Содержание к диссертации
Введение
1 Многомасштабная модель криогенного теплопереноса в биотканях 16
1.1 Макроскопическая модель Пенне для описания процессов теплопере-носа в биотканях 16
1.2 Энтальпийный метод расчета фазового перехода первого рода 18
1.3 Альтернативная постановка задачи о фазовом переходе 22
1.4 Основные модели внутриклеточной кристаллизации и трансмембранного переноса. Пределы детализации макроскопической модели 23
1.5 Дифференциальная модель и разностная схема для расчета динамики кристаллизации внутриклеточной жидкости 26
1.6 Результаты численного моделирования 27
1.7 Криодеструкция злокачественного образования в предположении сферической симметрии наконечника
1.7.1 Разностная схема для сферически симметричной задачи 36
1.7.2 Результаты численного моделирования криогенного разрушения биоткани сферически симметричным наконечником 38
1.8 Выводы к Главе 1 40
2 Метод усредненного фазового распределения для макроскопического расчета процессов теплопереноса в биотканях 42
2.1 Параметрическое представление тепловой функции 42
2.2 Макроскопическая постановка и разностная схема для сферически симметричной задачи 45
2.3 Оптимальные параметры фазового распределения для зондовой криохирургии 46
2.4 Расчет глубины проникновения зоны заморозки. Сравнение с экспериментальными данными 50
2.5 Трехмерная математическая модель зондовой криохирургии 51
2.5.1 Разностная схема для трехмерной задачи зондовой криохирургии.
2.5.2 Метод релаксации потоков для повышения устойчивости явной схемы 58
2.5.3 Исследование сходимости разностной схемы по сетке. 59
2.5.4 Результаты численного моделирования трехзондовой криохирургии на предстательной железе с опухолью
2.6 Результаты трехмерного численного моделирования зондовой криохирургии на злокачественном образовании округлой формы 65
2.7 Выводы к Главе 2 66
3 Оптимизация криохирургических операций 67
3.1 Задача многокритериальной оптимизации в вычислительной медицине 67
3.2 Основные положения теории многокритериальной оптимизации по Па-рето 69
3.3 Простейшие задачи многокритериальной оптимизации при моделировании криогенного воздействия на биоткань 71
3.4 Частный квазиградиентный метод 74
3.5 Общий квазиградиентный метод
3.5.1 Общий квазиградиентный метод для двумерной задачи 78
3.5.2 Квазиградиентные схемы для трехмерных отображений 81
3.5.3 Сужение двумерной квазиградиентной схемы на фронт Парето
3.6 Оптимизация зондовой криохирургии 85
3.7 Численное моделирование и оптимизация операций кожной криохирургии 85
3.8 Выводы к Главе 3 91
Заключение 93
Приложения:
A Реализация алгоритма двухмасштабного моделирования воздей ствия сферического крионаконечника 94
B Реализация основных алгоритмов для трехмерного моделирования и оптимизации зондовой криохирургии 97
Список литературы
- Основные модели внутриклеточной кристаллизации и трансмембранного переноса. Пределы детализации макроскопической модели
- Оптимальные параметры фазового распределения для зондовой криохирургии
- Результаты численного моделирования трехзондовой криохирургии на предстательной железе с опухолью
- Простейшие задачи многокритериальной оптимизации при моделировании криогенного воздействия на биоткань
Введение к работе
Актуальность работы. Рак предстательной железы – одно из самых распространенных злокачественных заболеваний, диагностируемых у мужчин. Ежегодно в мире выявляется свыше 400 тысяч новых случаев рака предстательной железы и около 200 тысяч человек каждый год умирает от этой болезни. При этом прирост заболеваемости опережает рост заболеваемости всеми другими видами рака. В этой связи актуальной задачей является поиск новых эффективных методов лечения онкологических заболеваний. В настоящее время наряду с химиотерапией и лучевыми методами активно развиваются и внедряются методы криохирургического воздействия для лечения широкого спектра воспалительных и онкологических заболеваний.
Впервые использование низких температур для лечения рака молочной железы и шейки матки предложил англичанин James Arnott в середине XIX века. Впечатляющие клинические результаты и доступность криохирургических методов способствовали тому, что последние получили прогрессивное развитие. Совершенствование техники привело к расширению области применения криохирургии: онкология, дерматология, гастроэнтерология, оториноларингология и др.
Известно, что при охлаждении биоткани происходит дегидратация клетки из-за резкого повышения осмотического давления внутриклеточной жидкости за счет гипергликемии и разницы в концентрациях электролитов во вне- и внутриклеточном пространствах. В этой связи понижается температура кристаллизации внутриклеточной жидкости, однако нормальное функционирование клеточной системы оказывается уже невозможным, что приводит к гибели клетки.
Снижение температуры до -10С приводит к началу кристаллообразования во внеклеточном пространстве, при этом также в интервале температур до -15C начинается образование льда внутри клетки, что увеличивает поражение биологических тканей. Максимальный эффект достигается при охлаждении биоткани до -50С. В то же время интенсивность клеточного некроза в значительной степени зависит от скорости заморозки биоткани. Экспериментально установлено, что скорость охлаждения, обеспечивающая наименьшую выживаемость клеток ткани, составляет порядка 50K в минуту.
Для количественного описания биологических процессов, возникающих в биотканях в ходе операций криохирургии, последние несколько десятилетий активно развиваются методы математического моделирования биологических систем. Появилась возможность разработки программных комплексов, позволяющих на основе
численного моделирования криогенного разрушения образований сложной геометрической формы давать в разумное время рекомендации практикующим хирургам по выбору оптимальных режимов проведения операций криохирургии.
Задачи с фазовыми переходами повсеместно встречаются в широком спектре приложений, поэтому методы их численного решения достаточно хорошо развиты. Однако, для корректного описания фазового превращения жидкости в биологической ткани необходимо учесть, что последнее оказывается размытым на некотором температурном диапазоне, вообще говоря, зависящем от конкретного температурного режима. В этой связи значительную часть работ по математическому моделированию криогенных процессов в биотканях можно разделить на два класса по используемым подходам: макроскопические и микроскопические модели.
В рамках макроскопического подхода, как правило, проводится численное моделирование процессов теплопроводности в биоткани с учетом фазовых переходов. При этом используются как разностные методы, так и конечно-элементный подход, а особенности кристаллизации биологических растворов как правило не учитываются.
В рамках микроскопического описания криогенных процессов в биологических тканях рассматривается система «клетка – внеклеточное пространство». При этом детально исследуются возникающие процессы дегидратации и внутриклеточной кристаллазции рассматриваемой системы. Актуальную и по сей день модель дегидратации клетки при ее охлаждении в 1963 г. предложил и численно исследовал Mazur. Для описания процесса внутриклеточной кристаллизации в 1990 г. в работе Toner et al. была предложена термодинамическая модель, получившая развитие в более поздних работах.
Несмотря на развитость методов макроскопического численного моделирования процессов теплопереноса и математических моделей для описания криогенных процессов на клеточном масштабе, работ, посвященных изучению взаимного влияния микро- и макроскопических процессов в биотканях друг на друга и, главным образом, непосредственно на свойства тканей, крайне мало.
Другим важным аспектом вычислительной медицины является задача отыскания оптимальных режимов того или иного рода воздействий. Так в приложении к операциям криохирургии разумна постановка задачи о нахождении оптимального способа расположения крионаконечников в целевой области для минимизации ущерба прилегающей здоровой ткани и максимального разрушения тканей злокачественного образования.
Интерес представляют работы, в которых оптимизация и прогнозирование операций рассматриваются как неразрывно связанные задачи и решаются на основе непосредственного математического моделирования тепловых процессов в биотканях. Для решения задачи минимизации разрушения здоровых тканей в свое время
предлагались симплексные методы, статистические методы муравьиных колоний и методы силовых полей.
В зависимости от целевой области (почка, печень, кожный покров, предстательная железа и т.д.) требования к режимам проведения операций могут отличаться. Например, даже малое повреждение тонких стенок кишечника может привести к необратимым последствиям. Ввиду этого проведение повторного точечного воздействия оказывается более предпочтительным, чем разовое разрушение всего объема злокачественного образования с сопутствующим ущербом здоровой биоткани. При предоставлении практикующим специалистам дифференцированного выбора между разрушением целевой области и сохранением здоровой биоткани естественным образом возникают задачи многокритериальной оптимизации.
Целью диссертационной работы является разработка эффективных численных алгоритмов прогнозирования и оптимизации режимов операций криохирургии и создание соответствующих программных комплексов.
Методы исследования. Исследование процессов теплопереноса в биотканях в ходе операций зондовой и кожной криохирургии строилось на основе методов математического моделирования. В основу предложенных математических моделей были положены законы сохранения и термодинамического равновесия. В значительной степени для микроскопического описания процессов кристаллизации биологических растворов использовались экспериментально установленные соотношения. При разработке алгоритмов численного моделирования исследуемых процессов применялась теория разностных схем и метод конечных объемов. При разработке численных методов построения решений задачи отыскания оптимальных режимов проведения операций использовались методы теории оптимизации. Реализация численных алгоритмов проводилась в среде Microsoft Visual C++ Express на языке C++.
В диссертационной работе решены следующие задачи
Разработана математическая модель для описания процессов кристаллизации и дегидратации внутриклеточной жидкости при криогенном воздействии на биоткань;
Проведено численное исследование влияния температурных режимов на характер кристаллизации внутриклеточной жидкости;
Разработана двухмасштабная математическая модель для описания криогенного воздействия, учитывающая влияние особенностей клеточных процессов на макроскопические свойства биоткани;
Проведено численное моделирование и построены температурные профили для операции зондовой криохирургии с наконечником сферической формы;
Для прогнозирования операций зондовой криохирургии на областях сложной геометрической формы предложено двухпараметрическое семейство макроскопических моделей, основанных на фазовом осреднении макроскопических параметров биоткани;
На основе численных экспериментов определены оптимальные параметры фазового осреднения, гарантирующие согласованность макроскопической модели с клеточными процессами;
Разработаны алгоритмы и реализованы программные комплексы для трехмерного численного моделирования операций зондовой и кожной криохирургии;
Разработан экономичный алгоритм численного решения задач многокритериальной оптимизации;
Для повышения эффективности операций криохирургии реализован алгоритм численного решения задачи отыскания оптимальных режимов их проведения.
Научная новизна работы
Предложена новая двухмасштабная математическая модель для описания криогенного воздействия, учитывающая влияние особенностей клеточных процессов на макроскопические свойства биоткани;
Для двухмашстабного моделирования операции зондовой криохирургии разработан и реализован алгоритм раздельного счета по макроскопическим и микроскопическим параметрам биоткани;
Предложена модификация энтальпийного метода макроскопического расчета фазовых переходов, согласованная с внутриклеточными процессами;
Показано, что параметры фазового осреднения являются функциями температурного режима и должны определяться для каждого типа операции (расположение и размеры опухоли, количество и спецификация зондов) отдельно;
Разработан экономичный алгоритм для трехмерного численного моделирования операций криохирургии на злокачественных образованиях сложной геометрической формы;
Для повышения устойчивости явных схем численного решения задач теплопе-реноса с фазовым переходом впервые использован метод релаксации потоков;
Для уменьшения ущерба здоровым тканям при операциях кожной криохирургии предложен метод теплового удержания фронта клеточного некроза,
основанный на введении в биоткань дополнительных нагревательных элементов;
Предложены квазиградиентные схемы численного решения задач многокритериальной оптимизации, учитывающие неединственность решения;
Показано, что применение методов многокритериальной оптимизации позволяет значительно увеличить эффективность и безопасность операций криохирургии.
Обоснованность и достоверность результатов работы подтверждается согласованием полученных результатов с экспериментально установленными закономерностями и результатами других авторов.
Апробация работы. Основные результаты и положения диссертации докладывались и обсуждались на научных конференциях и семинарах:
Научный семинар кафедры «Прикладная математика» НИЯУ МИФИ «Проблемы современной математики» под руководством Н.А. Кудряшова в 2013– 2016 годах;
Научная сессия НИЯУ МИФИ, Москва, 2014–2015 гг.;
13-th International Conference of Numerical Analysis and Applied Mathematics, Greece, Rhodes, September 2015;
V международная конференция «Проблемы математической и теоретической физики и математическое моделирование», Москва, апрель 2016 года;
Практическая значимость работы состоит в следующем:
Созданный в рамках диссертационной работы программный комплекс может быть использован в медицинских учреждениях на этапе планирования проведения операций криохирургии и позволяет давать рекомендации по выбору оптимальных режимов их проведения;
Разработанная двухмасштабная модель и созданный на ее основе программный комплекс могут быть использованы в исследовательских учреждениях для детального изучения особенностей биологических процессов в биотканях и разработки методов управления последними;
Разработанные квазиградиентные схемы решения задач многокритериальной оптимизации могут быть использованы при решении широкого круга прикладных задач связанных с отысканием оптимальных значений некоторых управляющих параметров а также при исследовании свойств сюрьективных отображений;
На защиту выносятся:
Двухмасштабная модель для описания криогенных процессов теплопереноса в биотканях;
Программный комплекс, позволяющий проводить двухмасштабное математическое моделирование процессов теплопереноса вблизи поверхности крионако-нечника;
Результаты численного моделирования криогенного разрушения биоткани вблизи поверхности сферического криозонда;
Согласованная макроскопическая модель теплопереноса в биотканях для крупномасштабного моделирования зондовой и кожной криохирургии;
Квазиградиентные схемы численного решения задач многокритериальной оптимизации с целевыми функциями, не имеющими аналитической формы записи;
Метод тепловой локализации фронта клеточного некроза при проведении операций кожной криохирургии;
Комплексы программ для трехмерного численного моделирования операций зондовой и кожной криохирургии и отыскания оптимальных режимов их проведения;
Результаты серии вычислительных экспериментов по моделированию и оптимизации зондовой криохирургии;
Результаты трехмерного численного моделирования теплового удержания фронта клеточного некроза с помощью дополнительных нагревательных элементов при проведении операции кожной криохирургии.
Структура и объем диссертации. Диссертация состоит из введения, трех глав, заключения, двух приложений и списка литературы. Диссертация содержит 111 машинописных страниц, 30 рисунков и 5 таблиц. В список литературы включено 105 наименований.
Основные модели внутриклеточной кристаллизации и трансмембранного переноса. Пределы детализации макроскопической модели
Диссертационная работа состоит из введения, трех глав, заключения, списка литературы и двух приложений.
Первая глава диссертационной работы посвящена разработке двухмасштабной модели для описания криогенных процессов в биотканях. Для учета особенностей процессов дегидратации клетки и внутриклеточной кристаллизации предлагается дифференциальная модель, основанная на законах сохранения и описывающая динамику поведения жидкой и кристаллической компонент во внутриклеточном пространстве. Численно исследуется влияние температурных режимов на фазовый состав внутриклеточной среды. Получены количественные характеристики для температурных диапазонов размытия фазового превращения и долей кристаллизованной и дегидратированной внутриклеточной жидкости при различных скоростях охлаждения.
Также обсуждаются особенности применения энтальпийного метода для макроскопического расчета задач с размытым фазовым переходом. С учетом поправок на осмотически неактивную и связанную жидкость предлагается модель двухмасштаб-ного описания криогенных процессов в биотканях. При этом на масштабах, много бoльших размеров клетки, макроскопические параметры среды, такие как теплопроводность и теплоемкость описываются как усредненные функции фазового состава системы «клетка-внеклеточное просранство».
В разделе 1.7 рассматривается задача о криогенном разрушении биоткани криона-конечником в предположении его сферической симметрии. Для предложенной двух-масштабной модели строится явная конечно-объемная аппроксимация. При этом для перехода на новый временной слой в каждой ячейке расчетной области численно решается задача Коши для эволюции параметров фазового состава среды. По результатам численного эксперимента строятся температурные профили для криогенного разрушения биоткани вблизи поверхности крионаконечника.
Вторая глава диссертационной работы посвящена методам трехмерного макроскопического моделирования операций зондовой криохирургии на областях сложной геометрической формы. Для согласования макроскопического описания криогенных процессов в биоткани с особенностями клеточных процессов вводится двухпарамет-рическое семейство энтальпийных моделей отвечающих наборам усредненных параметров фазового распределения. В качестве параметров фазового осреднения используются нижняя пороговая температура фазового превращения и величина, имеющая смысл приведенной температуры фазового перехода. С помощью численных экспериментов на задаче о сферическом криозонде показано, что с помощью макро скопической модели с фазовым осреднением можно добиться 98% согласования с результатами, полученными на основе двухмасштабной модели. Показано, что отличие значений параметров фазового осреднения от их априорных оценок носит локальный характер и определяется так называемой прицельной скоростью охлаждения.
В разделе 2.5 рассматривается задача трехзондовой криохирургии на модельном злокачественном образовании округлой формы. Ввиду того, что высокая скорость охлаждения достигается в малых по сравнению с размерами образования областях, типичные протоколы криохирургии глобально характеризуются достаточно мягкими температурными режимами. Таким образом в пределах приемлемой погрешности в 5-7% значения параметров фазового распределения при малом количестве криона-конечников отвечают простейшим априорным оценкам.
На основе макроскопической модели для рассматриваемой задачи получена явная трехмерная конечно-объемная аппроксимация и разработан программный комплекс, допускающий использование неортогональных индексных сеток. Для повышения устойчивости используемой явной схемы применяется метод релаксации потоков.
В третьей главе диссертационной работы рассматриваются способы повышения эффективности и безопасности операций криохирургии и численные методы решения возникающих при этом оптимизационных задач. В начале даются основные положения теории многокритериальной оптимизации по Парето, и на примере модельной задачи зондовой криохирургии с образованием кубической формы рассматриваются примеры возникновения неоднозначности решения. В качестве набора целевых переменных воздействия рассматривается конфигурация системы крионаконечников, а целевыми функциями, подлежащими минимизации, являются объемы неразрушенной ткани злокачественного образования и поврежденной здоровой биоткани.
Рассматриваемые в данной главе задачи оптимизации обладают той особенностью, что целевые функции не имеют аналитической записи и вычисляются непосредственно для каждого заданного набора целевых переменных. При этом оптимальное решение по Парето не единственно, а представляет собой кусочно-связное множество точек, так называемый фронт Парето. Для эффективного решения таких задач предлагается квазиградиентный подход, позволяющий итерационно строить сегменты фронта Парето. Отличительной особенностью предлагаемого подхода является то, что задачи моделирования операции криохирургии и ее оптимизация рассматриваются неразрывно друг от друга.
Раздел 3.5 посвящен рассмотрению возможных обобщений квазиградиентных схем. Для отыскания всех сегментов фронта Парето задач многокритериальной оптимизации производится обобщение постановки, приводящее к задаче отыскания границы множества образов сюрьективных отображений. Вводятся понятия локального выпуклого продолжения и дополнения, позволяющие произвести обобщение двумерной задачи на случай трех целевых функций. По аналогии с теорией разностных схем строятся явные и неявные квазиградиентные схемы, реализующие численное построение искомой поверхности.
В разделе 3.6 рассматривается применение частного квазиградиентного метода для отыскания оптимального расположения крионаконечников внутри модельного образования рассмотренного ранее. Показано, что оптимальный выбор конфигурации криозондов может значительно повысить степень криогенной деструкции нежелательной области и снизить вред прилегающей здоровой биоткани.
В разделе 3.7 численно решается задача моделирования операции кожной криохирургии на модельном образовании сложной геометрической формы. В отличие от зондовой криохирургии в этом случае хладогент доставляется в целевую область путем напыления LN2-спрея на кожный покров. Для прогнозирования кожной криохирургии используется модификация программного комплекса для моделирования зондовых операций. Ввиду того, что распространение зоны промерзания и, следовательно, зоны клеточного некроза осуществляется анизотропно, предлагается и исследуется методика теплового удержания фронта некроза с помощью дополнительных нагревательных элементов, вводимых в прилегающую здоровую биоткань. Рассматривается задача оптимизации такой операции, при этом целевыми переменными являются конфигурация сдерживающих элементов и время криохирургического воздействия. Для численного решения экстремальной задачи используется частный квазиградиентный метод с расщеплением по времени. Показано, что введение нагревателей слабой мощности позволяет также значительно уменьшить ущерб здоровым тканям кожного покрова.
В Приложении A представлены фрагменты кода программного комплекса для двухмасштабного моделирования температурного распределения вблизи поверхности крионаконечника. Приводятся реализации алгоритмов вычисления микроскопических и макроскопических параметров биоткани и расчета эволюции сеточной функции температурного распределения.
Приложение B содержит основные фрагменты кода разработанного программного комплекса для трехмерного моделирования операций криохирургии и реализацию этапа получения Парето оптимального решения в рамках квазиградиентного метода отыскания оптимальной конфигурации крионаконечников.
Оптимальные параметры фазового распределения для зондовой криохирургии
Как правило, практические задачи прогнозирования результатов криохирургического воздействия связаны с численным расчетом пространственно неоднородных процессов теплопереноса в областях сложной геометрической формы. В этой связи экономичность расчетного алгоритма оказывается одним из важнейших критериев применимости того или иного расчетного комплекса в медицинской практике.
Несмотря на согласованность макроскопических тепловых процессов с особенностями внутриклеточной кристаллизации в предложенной в Главе 1 многомасштабной модели, ресурсоемкость последней значительно увеличивается при измельчении расчетной сетки и отсутствии каких-либо групп симметрии. Таким образом, метод введения эффективной теплоемкости С подобной (1.10), (1.16) и (1.17), на основе некоторой сглаженной аппроксимации тепловой функции Я вида, например, (1.11) или (1.15) представляется более перспективным с точки зрения экономичности расчетного алгоритма.
С другой стороны, температурные распределения и характеристики фазового состава клеточного пространства в рамках многомасштабной модели хорошо согласуются с экспериментальными данными, в то время как при макроскопическом энталь-пийном подходе в интервале размытого фазового перехода зависимость соотношения между жидкой и кристаллической компонентами в единице объема от температурного режима либо учитывается грубо, либо не учитывается вовсе. В связи с этим возникает вопрос о возможности модификации энтальпийного метода таким образом, чтобы полученные на ее основе результаты численного моделирования с хорошей точностью приближали расчеты для многомасштабной модели. В рамках задачи численного моделирования криогенного воздействия на биоткань априорная оценка усредненных параметров фазового перехода оказывается возможной.
Известно, что в биотканях процесс фазового перехода начинается при температуре Тс 272.62Х. Введем первый параметр фазового перехода - величину Та, суть температуру, при которой для заданного температурного режима массовая доля жидкой компоненты в биоткани становится много меньше соответствующей величины для кристаллической составляющей. Иначе говоря Та - это характерная температура завершения процесса фазового превращения, а температурный интервал (Та,Тс) - характерный интервал, на котором размыт фазовый переход.
Далее, введем второй безразмерный параметр 6 = (Т - Та)/(ТС - Та), где Т есть некоторое число, причем Т Є (Та,Тс). Экстраполируем тепловую функцию Н{Т : Т Тс) на интервал (1 - 8){ТС - Та) влево. Тогда на интервале (Та,Тс) для приращения тепловой функции имеем: т AH = pL+ f Cfr(T )dT + Cunfr(Tc - Т ), Т = Т {8). (2.1) та С помощью линейной интерполяции сгладим Н при Т Є (Та, Тс), тогда для усредненной теплоемкости с учетом малости поправки, связанной с нелинейностью Н при Т Т , по сравнению с pL получаем: С = Ci(T), Т Та т + 5Сг(Та) + (1 - 5)С2 Та Т Тс (2.2) С2, т тс Для случая 5 {0.7-0.9} (2.2) можно обобщить. Используя оценку f Т в (1.12), или разложив C(Taipha) по степеням {Таіфа — Т) и отбросив члены, много меньшие слагаемого, пропорционального pL, можно ввести эффективную теплоемкость в ви де: Сі(Т), т та с = ф + 6С\(Т) + (1 - 5)С2 Та Т ТС с2, т тс (2.3) Отметим один предельный случай, когда Сг = Const, Тс - Т„ = 2Д, А/Тс 1 и 6 = 1/2. Из последнего предположения и связи между величинами Т и 6 очевидно, что Т = Тс - А = Та + А. Таким образом для эффективной теплоемкости (2.2) получаем: Си С = Т Та Ci+C2 U+ 2 та т тс с2, т тс После замены Та с = Т ± А окончательно получаем: (2.4) Си Т Т - А С = 2f + + , Т -А Т Т + А (2.5) Со Т Т + А Легко видеть, что с точностью до переобозначений Т Тс, А 6 выражение (2.5) совпадает с (1.10) для классической задачи Стефана. Таким образом, величина Т имеет физический смысл приведенной температуры фазового перехода, а параметр 6 = 5(Та,Т ) отвечает усредненному распределению жидкой и твердой фаз на характерном интервале фазового перехода (Та,Тс). Ясно, что при энталь-пийном подходе к решению задачи Стефана в чистой среде нет видимых причин для неравноправного распределения фаз на интервале размытия фазового перехода, этим объясняется слагаемое вида (Сг + С2)/2 в (2.5) и (1.10).
Отметим, что хотя из тройки величин {Т ,6,Та} только любые две являются независимыми, использование их всех помогает сократить записи и делает их более наглядными. 2.2 Макроскопическая постановка и разностная схема для сферически симметричной задачи
Рассмотрим задачу о криогенном разрушении биоткани сферически симметричным наконечником, численное двухмасштабное решение которой было рассмотрено в разделах 1.7.1 и 1.7.2. В рамках макроскопического подхода естественно ввести бипараметрическое семейство моделей теплопереноса: С(Т; Та,Т ) Т = \ г2к(Т;Та,Т ) Т + рьиьСь(Ть - Т) + qmet, v dt r2dr v dr v (2.6) t 0,L r R0, где коэффициент теплопроводности линейно интерполируется на интревале (Та,Т ): к(Т;Та,Т ) h(T), Т Та h(Ta) + k l(Ta) Та Т Т (2.7) к2, Т Т При этом функция теплоемкости определяется в соответствии с (2.2). Как показали численные эксперименты, при рассмотрении случая (2.3) влияние малой нелинейной поправки по сравнению с выделяемой теплотой кристаллизации оказывается много слабее выбора значений параметров макроскопического осреднения.
Результаты численного моделирования трехзондовой криохирургии на предстательной железе с опухолью
В связи с развитием методов математического моделирования и увеличения мощностей современных вычислительных систем для изучения широчайшего круга приложений в последние несколько лет часто ставится задача нахождения таких параметров исследуемой системы, при которых результаты численного эксперимента в определенной степени соответствуют некоторым ожиданиям. Как правило, в этом случае наряду с задачей численного моделирования изучаемого физического процесса ставится задача оптимизации одной или нескольких целевых функций, вычисляемых на основе данных численного эксперимента.
Задачи подобного рода возникают не только в связи с численным моделированием физических процессов, но и во многих других областях, где необходимо отыскание оптимальных значений некоторых величин. Так методы многокритериальной оптимизации находят свои приложения в экономических задачах пакетного инвестирования, логистики, военной стратегии и государственного управления. В качестве конкретного примера рассмотрим круг вычислительных задач, связанных с прогнозированием операций криохрурургии. Для более эффективного и безопасного лечения разумно иметь представление о рисках и перспективах применения криохирургического воздействия. В этой связи находят приложение методы математического моделирования. Можно также наблюдать расширение области научных исследований от непосредственного математического моделирования процессов биотеплопереноса в тканях в ходе криохирургического вмешательства [58,61] к разработке методов оптимизации параметров проведения операций [73, 77]. В этом случае задача состоит в нахождении таких параметров проведения операции криохирургии, при которых ее результаты удовлетворяют определенным критериям успеха и вероятного ущерба от криогенного воздействия.
Особенностью рассматриваемого круга задач является то, что подлежащие оптимизации целевые функции будут определяться для заданного набора аргументов из области определения без явной аналитической зависимости своих компонент от значений аргументов. Таким образом целевые функции могут мыслиться как некоторый «черный ящик» или автомат, возвращающий некоторый набор значений по введенным аргументам. Как следствие для любого набора целевых переменных может быть определена только аппроксимация якобиана перехода от пространства целевых переменных к пространству целевых функций, естественно в предположении достаточной гладкости последних. В этой связи классические градиентные методы оказываются неприменимы.
Вторая проблема связана с нарушением в привычном смысле единственности решения задачи многокритериальной оптимизации. В рамках подхода Парето [101–103] вводится отношение предпочтительности, и задача сводится к нахождению наиболее предпочтительных наборов целевых переменных. При этом все их наборы, относительно которых данное отношение не определено, называются решениями по Парето или Парето-оптимальными векторами, а их совокупность – фронтом Парето.
Наибольшей популярностью при решении задач многокритериальной оптимизации пользуются методы, связанные с различного вида скаляризацией вектора целевых функций, что в свою очередь нарушает общность решения как некоторого множества. Для случая, когда целевые функции вычисляются на основе тех или иных расчетов, активно применяются генетические алгоритмы. Именно такой подход использовался для оптимизации расположения крионаконечников при проведении криохирургических операций, путем двумерного моделирования теплофизических процессов в биотканях в работе [77], где для нахождения оптимального расположения наконечников ставилась двухкритериальная задача, и проводилась ее скаляризация.
Ясно, что жесткая скалярная фиксация целевых функций приводит к потере фронта Парето, и, как результат, потере альтернативных решений, хотя в различных приложениях могут допускаться различные вариации в предпочтении какому-либо из параметров. Выбор же конкретной конфигурации лежит непосредственно на хирурге и может быть основан на его опыте и понимании особенностей конкретной операции. Так в одном случае требуется тотальная деструкция образования, и допускается определенный уровень ущерба здоровым тканям. В других ситуациях даже минимальный ущерб здоровым тканям может повлечь непоправимые последствия и даже летальный исход, как, например, в случае со стенками пищевода [104].
Под решением задачи многокритериальной оптимизации понимают вектор целевых переменных X Є Df, который в некотором смысле доставляет минимум целевой функции F(X), формально: X = arg min F(X). (3.3) При этом необходимо определение смысла выражения minF(X). Заметим, что в случае, если существует такой номер что соответствующую компоненту целевой функции fio требуется максимизировать, после умножения её на -1, постановка (3.3) будет корректной. Наиболее общим подходом является оптимизация по Парето.
На множестве Df определим соотношение предпочтительности «У » по векторной функции F. Вектор целевых переменных Xi будем называть более предпочтитель ным по F, чем вектор Х2, и обозначать Хі У Х2, если выполняются следующие соотношения: /i(Xi) /i(X2), /2(Xi) /2(Х2), (3.4) /n(Xi) /„(Х2), и при этом одно из неравенств в (3.4) является строгим. Оптимальным по Парето называется такой вектор целевых переменных X , что X У X, для всех X Є Df, относительно которых соотношение предпочтительности вообще может быть определено. Введенное описанным образом соотношение предпочтительности определяет отношение частичного порядка на множестве Df. Нетрудно видеть, что оптимальное по Парето решение задачи (3.3) неединственно и представляет собой множество векторов Р С Df, между которыми отношение « » не может быть определено, тем не менее, для любого Х« Є Р и для любого Х/з Є Df\P справедливо соотношение Ха У Х .
Простейшие задачи многокритериальной оптимизации при моделировании криогенного воздействия на биоткань
Сужение двумерной квазиградиентной схемы на фронт Парето. При рассмотрении задач оптимизации операций криохирургии, как правило, рассматривают сюрьективные отображения F : R3n R2, т.е. V = F(9t), где Ж = (Ri, R2,... , Rn), Rj Є і?3 - вектор целевых переменных, характеризующих расположение крионаконечников, а V = {VT,VH} - целевые функции отвечающие объемам неразрушенной ткани образования и поврежденной здоровой биоткани. Очевидно, что по физическому смыслу целевой функции, подлежащей минимизации по Паре-то, ее компоненты VT иУЯ неотрицательны.
Далее, в силу определения отношения предпочтительности всякое Парето-оптимальное решение Н удовлетворяет неравенствам: VWarg min(Vff)) VT№) IC1 (3.33) VWarg тіп(Уг)) VH( R) Vn, т.е. образ каждого Парето-оптимального решения находится внутри прямоугольника, построенного на вершинах, порожденных решениями однокритериальных задач. В предположении гладкости компонент Уу и Уя как функций целевых переменных Н для решения задачи оптимизации квазиградиентную схему можно сузить, взяв в качестве первой Парето-оптимальной точки решение однокритериальной задачи, например Н = arg min VT. Тогда для построения локальных выпуклых продолже-ний следует определить направление убывания VH. В случае если VWarg min (VH)) ф
Vfnm, когда фронт Парето состоит из одной точки, и задача уже решена, локально это направление можно определить. При достижении точки второго минимума, построение решения Парето можно считать завершенным. Очевидно, что все сегменты последнего будут лежать на кривой K(t) : VT{SR{t0)) = Vf1 1, VH{?R{ti)) = V m. В предположении гладкости целевых функций из существования кривой на dVf, соединяющей точки (FT(arg min(VH), V7?ira) и (V m,VWarg min(VT))) по построению кривая в R3n, порождающая ее, конечно же будет существовать.
Заметим также, что особое значение имеют конфигурации крионаконечников, при которых злокачественное образование разрушается полностью. Ясно, что в общем случае добиться этого удается не всегда. Однако при увеличении числа степеней свободы системы за счет введения параметра времени воздействия соответствующее калибровочное условие часто можно наложить. Такая возможность при достаточной тепловой мощности крионаконечников обусловлена продвижением зоны криогенного Рис. 3.6. Эволюция зоны клеточного некроза. Парето-оптимальная конфигурация крионаконечников некроза от поверхностей крионаконечников при увеличении времени воздействия.
Вернемся теперь к задаче о криогенном разрушении злокачественного образования (Рис. 2.4) тремя крионаконечниками, рассмотренной в 2.5. В разделе 2.6 уже упоминалось о том, что для эффективного проведения такой операции возникает необходимость в отыскании оптимального расположения криозондов.
С помощью частного квазиградиентного метода удалось отыскать ряд Парето-оптимальных способов расположения трех крионаконечников. На Рис. 3.6 представлена эволюция зоны клеточного некроза для одной из компромиссных конфигураций.
По сравнению с результатами, полученными в разделе 2.6, проведение операции с такой конфигурацией крионаконечников приводит к увеличению доли разрушенной злокачественной биоткани на 23% и снижению ущебра здоровой ее части в 3 раза.
Двумерное численное моделирование такой задачи проводилось, например, в [105]. Рассмотрим трехмерную постановку задачи для модельного злокачественного образования в кожном покрове.
На Рис. 3.7 изображена форма рассматриваемого злокачественного образования в кожной ткани, аппроксимированная выпуклыми шестигранниками.
Поскольку при поверхностном нанесении хладогента на кожный покров фронт клеточного некроза от зоны контакта будет распространяться во всех направлениях, он, в общем случае, может нанести серьезный вред здоровой биоткани. Для снижения риска повреждения здоровых тканей, прилегающих к пораженной области, предлагается методика тепловой локализации фронта, при которой для контроля распространения зоны клеточного некроза в здоровой ткани размещаются дополнительные нагревательные элементы.
Ясно, что как и конфигурация крионаконечников в зондовой криохирургии, расположение сдерживающих нагревательных элементов должно влиять на результирующее распределение зоны клеточного некроза в биоткани. Для увеличения сдерживающего эффекта следует рассматривать задачу оптимизации расположения нагревательных элементов, родственную задаче оптимизации зондной криохирургии.
Для описания процессов теплопереноса в кожной биоткани будем использовать модель Пенне: дТ рС dt (3.34) -dzvQ + pbujbCb(Tb), Q = -KVT, где T, р, С, к - температура кожной ткани, плотность, удельная теплоемкость и теплопроводность соответственно. Для учета фазового перехода, размытого в температурном диапазоне от верхней границы Ти = — 1С до нижней границы 7) = —8С осредненные термодинамические характеристики среды определяются следующим образом: pC(T Tl) = pCf(T), РС{Тг Т TU) = - {pCfiTt) + рСи) + , А=\Ти-Т1\/2 (3.35) рС(Т Ти) = рСи. Здесь Си и С/ - удельные теплоемкости незамороженной и замороженной тканей соответственно. Коэффициент шь в (3.34) предполагается линейно убывающим на интервале температур от шъ = ш0 при Т = Ти до шъ = 0 при Т = Т%.
Предлагаемая математическая модель дополняется нелинейной зависимостью коэффициента теплопроводности кожной ткани от температуры в области ниже температуры полной кристаллизации жидкой компоненты биоткани. Учитывается также линейная температурная зависимость теплоемкости среды в той же области температур. Численные значения основных характеристик среды представлены в Таблице 3.1.
Термодинамические параметры кожной биоткани плотность р 1116 теплоемкость(для замороженной ткани) Сиcf 2721 4к89.7 + 5.04Т теплопроводность(для замороженной ткани) ки 0 39 №. т К1.46 + 3.88-10"3(273 Т)1-156 \ удельная теплота кристаллизации L 217,100 1 я Начальное температурное распределение в биоткани перед криохирургической операцией предполагается равномерным: T(t = 0,r) = Tbody (3.36) Здесь Tbody - равновесная температура тела (37 C). Разумно изолировать от окружающей среды поверхность кожного покрова, не подвергающегося нанесению хла-догента. Область нанесения LN2 спрея обозначим как П.
Без ограничения общности будем рассматриваеть случай, когда область нанесе ния спрея Q находится на границе расчетной области Г. Предполагая границы расчетной области Г, не являющиеся поверхностью кожи, достаточно удаленными от области криогенного воздействия считаем: (Q,n)nn = 0. (3.37) Здесь п - вектор нормали к границе расчетной области. Нанесение LN2 спрея на поверхность контакта с кожей считается достаточно интенсивным для мгновенного образования на последней слоя кипящего азота, поддерживающего поверхность контакта при постоянной температуре: T(t, г Є ft) = ТLN2 (3.38)
Поскольку температура нагревателей предполагается недостаточно высокой для нанесения ощутимого вреда здоровой ткани, будем считать, что элементы нагрева помещены в биоткань до начала операции, и в малой их окрестности установилась и поддерживается постоянная температура Th. Таким образом при численном моделировании операции из расчетной области исключаются указанные окрестности, и вводятся дополнительные граничные условия, описывающие теплообмен на границе по закону Ньютона с коэффициентом обмена, равным коэффициенту теплопроводности среды.
Аналогично разделу 2.5 для численного решения поставленной задачи будем использовать явную конечно-объемную аппроксимацию (2.28) с расчетом потоков методом контрольных объемов (2.34). Для повышения устойчивости схемы при численном решения задачи так же использовался метод релаксации потоков (2.36).