CC BY
19УДК 535.33/.34 535.3
DOI: 10.21779/2542-0321-2018-33-2-28-39
1 2 1 2 1 2 1 Э.Х. Исрапов ' , К.М. Гираев ' , М.А. Магомедов ' , Н.А. Ашурбеков
Численное моделирование пространственного распределения температуры в биотканях по мере их лазерного нагрева
1 Дагестанский государственный университет; Россия, 367001, Махачкала, ул. Гаджи-ева, 43а; Ed-gar1993@yandex.ru;
2 Институт физики им. Х.И. Амирханова ДНЦ РАН; Россия, 367003, г. Махачкала, ул. М. Ярагского, 94
В работе посредством решения нестационарного уравнения теплопроводности рассчитано пространственное распределение температуры в тканях серого вещества головного мозга по мере развития процессов лазерной гипертермии. Обнаружено, что для данных условий лазерного облучения по достижению температуры 60.0 0С наибольший объем фототермического повреждения 1.0 см3 и характеристического расстояния 0.5 мм при экспозиции 6.0 мин. наблюдается на длине волны 1064 нм. На других длинах волн (300 и 1460 нм) эти параметры могут быть в 2-4 раза меньше. Показано, что развитие лазерной термокоагуляции приводит к дополнительному сокращению объема, характеристического рассеяния и времени экспозиции в 1.5-2 раза на всех исследуемых длинах волн.
Ключевые слова: биоткани, уравнение теплопроводности, численное моделирование, лазерная термокоагуляция, экспозиция, объем повреждения.
Введение
Изучение процессов фотонагрева и тепломобмена в биологических тканях является важной задачей современной биологической и медицинской физики [1, 2]. Наиболее актуальными эти исследования стали при развитии экспериментальных и теоретических основ различных методов лазерной терапии и диагностики онкологических заболеваний. В частности, для корректного определения параметров дозиметрии лазерного излучения (длины волны и мощности излучения, длительности экспозиции) при проведении процедур радиометрической диагностики и фотодинамической или лазерно-индуцированной термотерапии злокачественных новообразований необходимы сведения о пространственном распределении температуры в биотканях и ее динамике в процессе лазерного нагрева (см., например, [1-5]).
В настоящее время теоретические методы исследования процессов теплообмена базируются преимущественно на математическом моделировании, которое включает в себя различные аналитические и численные методы (см., например, [1-5]). При этом наиболее распространённым и эффективным методом исследования процессов тепло-переноса в биотканях является метод конечно-разностных схем, применение которого стало возможным благодаря значительному прогрессу в развитии вычислительных методов и увеличению мощности современных компьютеров [4-7]. В то же время при анализе эффектов лазерной гипертермии необходимо принимать во внимание ряд фи-
зиологических и морфофункциональных характеристик, определяемых индивидуальными особенностями тканевой и внутриклеточной структуры, распределением соединительной и жировой ткани и др., что, в свою очередь, может быть учтено посредством исследования оптических параметров среды - коэффициентов поглощения (да), рассеяния и фактора анизотропии ($).
Данная работа посвящена исследованию процессов переноса тепла и их динамике в ходе лазерного нагрева биологических тканей. С этой целью был разработан алгоритм численного решения уравнения теплопроводности и определено пространственное распределение температуры в биотканях интактного (неизменного) состояния и по мере воздействия мощного лазерного излучения с использованием данных их оптических свойств.
Материалы и методы
Расчет пространственного распределения температуры в биотканях в процессе лазерного зондирования производился при помощи решения нестационарного уравнения теплопроводности (уравнение Фурье-Кирхгофа) (1), записанного в декартовой системе координат, методом конечных разностей (МКР) [8].
РС—^Т~ = Л I ^ + ду2 + У, г, О, (1)
где р - плотность ткани (кг/м3); с - удельная теплоемкость ткани (Дж/кг°С); Я - коэффициент теплопроводности ткани (Вт/м°С); Т - температура (°С) и Qabs - мощность внутренних источников тепла (Вт/м3). Последний параметр напрямую зависит от оптических коэффициентов исследуемой среды - да, и д, как:
Qabs=^re-^ff, (2)
где Р - мощность лазерного излучения; Б= ————-—- - коэффициент оптической
диффузии; = + (1 — д)1*5) - коэффициент эффективной экстинкции.
Уравнение (1) описывает множество вариантов развития процесса теплопроводности и устанавливает связь между временным и пространственным изменением температуры в любой точке исследуемой среды. Для того чтобы из бесчисленного количества этих вариантов выбрать требуемый и дать его полное математическое описание, к уравнению (1) необходимо добавить условия однозначности, которые содержат геометрические, физические, начальные и граничные условия.
Геометрические условия определяют форму и размеры тела, в котором протекает исследуемый процесс. Физические условия определяются набором теплофизических характеристик тела (Я, р, с). Начальные условия описывают распределение температуры в начальный момент времени: & = = /(х,у,г) - в общем случае. При этом, если считать распределение температуры в теле равномерным, то начальное условие упрощается: t = 0; Т = Т0 = сопз1. В качестве граничных условий использовалось условие второго рода, которое задает значение теплового потока для каждой точки границы
среды в любой момент времени: = 0уу(х,У,О, где п - нормаль к поверхности
тела. В данной работе считалось, что qw = сопз1.
Дифференциальное уравнение (1) вместе с граничными условиями, описанными выше, дает полную математическую формулировку краевой задачи теплопроводности. Для решения подобного уравнения наиболее часто используют метод конечных разностей, что обусловлено его точностью и простотой реализации. Идея метода конечных
разностей решения краевых задач видна уже из самого названия: вместо производных в дифференциальном уравнении используются их конечно-разностные аппроксимации. При использовании МКР для задач теплопроводности тело представляют в виде совокупности узлов. Заменяя частные производные уравнения (1) конечными разностями, получают систему линейных алгебраических уравнений для определения температуры как локальной характеристики в каждом узле пространственной сетки:
тп + 1_тп /тП+1 _-?тП+1,тП + 1 тП+1 _-?тП + 1,тП+1
рс—-— = А I----— Н—--;----\-
н т V ^ ^у
'М,к+1 А11,}Л | /-лп
где , ку, к2 - шаги интегрирования по пространственным координатам, т - шаг интегрирования по временной координате. Такая же конечно-разностная аппроксимация применяется для начальных и граничных условий. Подобные системы линейных алгебраических уравнений легко решаются методом прогонки [9].
Таким образом, располагая данными коэффициентов поглощения (да), рассеяния (д5) и фактора анизотропии рассеяния (д) для биотканей в норме и при лазерном нагреве, можно, используя описанную выше численную модель, сформировать картину пространственного распределения температуры и оценить ее динамику по мере развития процессов лазерно-индуцированной гипертермии.
Для определения оптических коэффициентов биотканей в работе проводилась процедура лазерно-индуцированной термокоагуляции, заключающаяся в зондировании мощным лазерным излучением образцов тканей серого вещества головного мозга (далее биоткани) в течение 5 и 10 минут. В последующем методом гистоморфологической микроскопии установлено, что термический эффект, вызванный лазерным нагревом с экспозицией 5 минут, при достижении температуры
70.0±2.0 0С классифицировался как коагуляция средней стадии, в то время как эффект от последовательного воздействия двух аналогичных актов - как коагуляция высокой стадии. При этом под коагулированным образованием в настоящей работе понимается процесс фототермического повреждения биоткани, при котором наблюдаются необратимые повреждения морфологии ткани, гибель клеток или денатурация молекулярных комплексов.
В качестве источника зондирующего излучения в работе использовался твердотельный Кё:УЛО лазер LQ529 (Солар Лазер Систем, Республика Беларусь) с длиной волны излучения 1064 нм, длительностью импульса 10 нс, частотой следования 100 Гц и средней мощностью 3 Вт.
Определение оптических коэффициентов (да, и д) биотканей проводилось при помощи метода спектрофотометрии интегрирующих сфер в спектральном диапазоне 250-2000 нм и инверсного метода Монте-Карло, подробно описанного в работе [10, 11]. При этом значения теплофизических величин для исследуемых биотканей принимались близкими к значениям таковых для воды (ткани головного мозга на 80 % состоят из воды), и составили согласно [6, 7]: р = 1,038 — 1,041 кг/м3, Я = 0,511 Вт/м°С и с = 4200 Дж/кг Т.
Результаты и обсуждение
Результаты численного моделирования пространственного распределения температуры для тканей серого вещества головного мозга по мере развития процессов лазерной термокоагуляции на длинах волн 300, 1064 и 1460 нм показаны на рисунках 1-3,
каждый из которых состоит из трех групп спектральных данных: А - спектр пространственного распределения температуры в биоткани; В - радиальное распределение температуры на различной глубине; С - азимутальное распределение температуры по профилю лазерного пучка. Значения оптических коэффициентов, длительности экспозиции лазерного излучения, необходимые для достижения фототермического повреждения исследуемых биотканей при t = 65 — 70 °С, а также в их динамика по мере развития процессов лазерно-индуцированной термокоагуляции на длинах волн 300, 1064 и 1460 нм приведены в таблице 1.
Таблица 1. Данные оптических коэффициентов поглощения - и транспортного рассеяния - = — , а также длительности экспозиции* - ^ объема повреждения** - V и характеристического расстояния*** - 1тах на различных длинах волн для тканей серого вещества головного мозга в интактном состоянии и при различных стадиях лазерно-индуцированной термокоагуляции
Состояние биоткани Длина волны, Я, нм ßa, мм 1 ,,' -1 ßs, мм Экспозиции t, мин Объем повреждения V, см3 Характ-ое расстояние Imax, мм
Интакт-ное состояние 300 0.59 6.1 1.5 0.3±0.05 0.1±0.05
1064 0.009 1.24 6.1 1.0±0.2 0.5±0.1
1460 0.31 0.9 3.26 0.54±0.05 0.3±0.05
Средняя стадия коагуляции 300 0.67 6.83 1.31 0.21±0.03 0.07±0.01
1064 0.015 1.94 5.2 0.67±0.05 0.4±0.05
1460 0.29 1.62 2.5 0.32±0.03 0.2±0.05
Высокая стадия коагуляции 300 1.1 9.12 1.18 0.16±0.02 0.05±0.01
1064 0.04 2.86 4.45 0.56±0.05 0.3±0.05
1460 0.28 2.24 2.27 0.28±0.03 0.15±0.05
* Длительность экспозиции лазерного излучения, необходимая для осуществления нагрева поверхности биоткани до температуры 70 °С.
** Объем фототермического повреждения биоткани при достижении температуры 60.0±5.0 °С. *** Глубина биоткани, на которую приходится максимальная температура 73.0 °С
Анализ результатов моделирования пространственного распределения температуры позволил установить, что для формирования области коагуляционного некроза объемом К~1.0 см3 при температурном пороге ~60.0 + 5.0 °С в интактных тканях серого вещества головного мозга при лазерном облучении с длиной волны 1064 нм и мощностью 3.0 Вт требуется экспозиция длительностью примерно минут. Вместе с этим зондирование биотканей лазерным излучением других длин волн при прочих в равных условиях приводит к уменьшению длительности экспозиции и объема фототермического повреждения, что для длины волны 300 нм составляет - 0,3 см3 и Ь~1.5 минуты и У~0,6 см и ?~3.3 минуты - для длины волны 1460 нм
/ о. 3\ В
17 1.4 \\
Х=1460 им
Рис. 1. Пространственное распределение температуры для тканей серого вещества головного мозга в интактном состоянии на длинах волн 300, 1064 и 1460 нм: А - карта пространственного распределения температуры в биоткани; В - радиальное распределение температуры на различной глубине; С - азимутальное распределение температуры по профилю
лазерного пучка
Л=1064 им
°с
70 ВО 50 40 30
г, mm -so
То X в
...........................А/.-.н
Л=1460 HM
Рис. 2. Пространственное распределение температуры для тканей серого вещества головного мозга при средней стадии термокоагуляции на длинах волн 300, 1064 и 1460 нм: А - карта пространственного распределения температуры в биоткани; В - радиальное распределение температуры на различной глубине; С - азимутальное распределение температуры по профилю лазерного пучка
°С
г, mm -so -40 -зо -го -ю о ю го эо 40 so S о ё 3 S 3
°С
70 60
¡0.15\ В
____________// Л:0 ' _________
1 \
..................... / 7 а \
—■—1—г—I— ^^ \ —-1-'-г--•-
Л=1460 нм
-12 0 7
Рис. 3. Пространственное распределение температуры для тканей серого вещества головного мозга при высокой стадии термокоагуляции на длинах волн 300, 1064 и 1460 нм: А - карта пространственного распределения температуры в биоткани; В - радиальное распределение температуры на различной глубине; С - азимутальное распределение температуры по профилю лазерного пучка
Кроме того, по результатам моделирования было выявлено (см. рис. 1-3), что максимум температуры (73.0 + 2.0 °С) приходится не на поверхность среды, как следовало ожидать, а на приповерхностную область на некотором характеристическом расстоянии 1тах под поверхностью биоткани, и далее монотонно снижается вглубь среды. Наиболее ярко это выражено при воздействии лазерным излучением на длине волны 1064 нм и составляет /т0^ = 0.5 + 0.1 мм, а менее всего - при зондировании на длине
волны 300 нм - I
зоо
0.1 + 0.05 мм. Данный факт, по-видимому, связан как с тепло-
обменными процессами на поверхности биоткани с окружающей средой, так и с эффектом усиления интенсивности лазерного излучения в приповерхностном слое среды вследствие влияния процессов светорассеяния, что полностью согласуется с результатами работы [12, 13].
Из рисунков видно, что при зондировании лазерным излучением инфракрасного диапазона (1064 и 1460 нм) наблюдается как максимальная глубина (порядка 4-6 мм), так и объем эффективного прогревания биоткани, что является следствием наибольшей проникающей способности лазерного излучения, обусловленной оптическими свойствами среды (дя, и д) на данных длинах волн. В то же время для УФ области спектра характерны минимальные глубина и объем прогрева, что объясняется высоким уровнем поглощения биотканей, ослабляющем интенсивность лазерного излучения и ограничивающем его распространение.
Сравнительный анализ полученных результатов (см. рисунки и данные таблицы) позволил установить, что по мере развития процессов гипертермии длительность экспозиции лазерного облучения t сокращается примерно на 15-20 % для средней и до 30 % для высокой стадии термокоагуляции, тогда как объем фототермического повреждения V уменьшается соответственно на 30-40 % и 50-60 %. Вследствие этого происходит уменьшение как характеристического расстояния 1тах, так и области эффективного прогревания в 1.5-2 раза. Однако дальнейшее увеличение экспозиции лазерного облучения приводит не столько к прогрессированию термокоагуляции биотканей, сколько к их обугливанию.
Объяснить подобную зависимость можно следующим образом. Как было показано в работах [6, 7, 10, 14], развитие процессов лазерно-индуцированной термокоагуляции в биотканях сопровождается дегидратацией клеток и денатурацией протеинов, приводящих к уплотнению и гранулированию гистоструктуры среды. Такие нарушения нативной структуры сводятся к увеличению плотности поглощающих и рассеивающих частиц, образующих в биотканях «экранирующий слой» с еще более высоким значением коэффициентов поглощения и рассеяния и низкой теплопроводностью. В свою очередь это приводит к дополнительной аккумуляции тепла в биотканях с образованием слоя нагара (карбонизации), препятствующего распространению лазерного излучения и нагреву более глубоких областей.
Заключение
Была проделана следующая работа:
1. Разработан алгоритм численного моделирования уравнения теплопроводности и определена динамика пространственного распределения температуры и параметры, характеризующие температурный и временной порог, а также объем фототермического повреждения для тканей серого вещества головного мозга по мере развития лазерно-индуцированной термокоагуляции.
2. Установлено, что для формирования области коагуляционного некроза объемом К1064~1.0 см3 при температурном пороге ~60.0 + 5.0 °С в биотканях при лазерном облучении на длине волны 1064 нм и мощностью 3.0 Вт требуется экспозиция длительностью примерно ^0б4~6 мин. На других длинах волн (300 и 1460 нм) значения этих параметров сокращаются и составляют соответственно К300~0,3 см3 и ^00~1.5 мин. и К14б0~0,6 см3 и *;1460~3.3 мин.
3. Показано, что максимум температуры (73.0 + 2.0 °С) приходится на приповерхностную область на некотором характеристическом расстоянии под поверхностью биоткани, равном /^д4 = 0.5 + 0.1 мм. В то же время на других длинах волн это расстояние сокращается и составляет - = 0.1 + 0.05 мм и = 0.3 + 0.1 мм, что объясняется наличием теплообменных процессов на поверхности биоткани с окружающей средой и усилением интенсивности лазерного излучения в приповерхностной области по причине рассеяния света.
4. Обнаружено, что по мере развития процессов гипертермии длительность экспозиции лазерного облучения t сокращается примерно на 15-20 % для средней и до 30 % для высокой стадии термокоагуляции, объем фототермического повреждения V уменьшается, соответственно, на 30-40 % и 50-60 %, а характеристическое расстояние 1тах и область эффективного прогревания в 1.5-2 раза. При этом дополнительное увеличение экспозиции приводит к обугливанию биотканей.
5. Анализ характера пространственного распределения температуры позволил предположить, что по мере развития процессов фотогипертермии за счет уплотнения структуры в биотканях формируется «экранирующий слой» с высоким уровнем поглощения и рассеяния, приводящий к аккумуляции тепла с образованием слоя нагара, препятствующего распространению лазерного излучения и нагреву более глубоких областей.
Литература
1. Optical-Thermal Response of Laser Irradiation Tissues / Ed. A.J. Welch and M.J.C. Van Gemert // Springer Science+Business Media B.V. - 2011. - 958 p.
2. Тучин В.В. Оптика биологических тканей. Методы рассеяния света в медицинской диагностике. - М.: Физматлит, 2013. - 811 с.
3. Ghosh S., Hanson W., Abdollahzadeh N. and Han B. Effects of light-tissue interaction on quantum dot mediated fluorescence thermometry // Meas. Sci. Technol. - 2012. -V. 23 045704. - P. 13.
4. Saccomandi P., Schena E., Giurazza F., Riccardo Del Vescovo et al. Temperature monitoring and lesion volume estimation during double-applicator laser-induced thermother-apy in ex vivo swine pancreas: a preliminary study // Lasers Med Sci. - 2013. - V. 29 (2). -607 p.
5. Sazgarnia A., Naghavi N., Mehdizadeh H. et al. Investigation of thermal distribution for pulsed laser radiation in cancer treatment with nanoparticle-mediated hyperthermia// J. of Thermal Biology. - 2015. - V. 47. - P. 32-41.
6. Marqa M.-F., Colin P., Nevoux P., Mordon S.R., Betrouni N. Focal laser ablation of prostate cancer: numerical simulation of temperature and damage distribution // Biomedical Engineering Online. - 2011. - V. 10. - P. 45.
7. Marqa M.-F., Mordon S.R., Betrouni N. Laser interstitial thermotherapy of small breast fibroadenomas: Numerical simulations // Lasers in Surgery and Medicine. - 2012. -V. 44. - P. 832-839.
8. Кузнецов Г.В., Шеремет М.А. Разностные методы решения задач теплопроводности. - Томск: Изд-во Томского политехнического университета, 2007. - 172 с.
9. Вержбицкий В.М. Основы численных методов. - М.: Высшая школа, 2002. -
840 с.
10. Гираев К.М., Ашурбеков Н.А., Магомедов М.А. и др. Влияние гипертермии на динамику спектров поглощения и светорассеяния нормальных и опухолевых биологических тканей // Вестник ДГУ. Естественные науки. - 2013. - Вып. 6. - С. 56-64.
11. Гираев К.М., Ашурбеков Н.А., Магомедов М.А. и др. Влияние патологических процессов на оптические спектры поглощения и рассеяния проб желчи и панкреатического сока // Оптика и спектроскопия. - 2015. - Т. 119. - С. 162-160.
12. Исрапов Э.Х., Магомедов М.А., Гираев К.М. и др. Численное моделирование пространственного распределения интенсивности лазерного излучения в биотканях в зависимости от степени их лазерного нагрева методом Монте-Карло // Вестник ДГУ. Сер.: Естественные науки. - 2017. - Т. 32. - С. 7-19.
13. Li Y., Sun Ch., Ji Ch., Liang X., Huang L., Kuang J., Wu J. Comparison of Photothermal Effects on In-Vivo Tissue between Moxibustion Therapy and Laser Irradiation // IEEE Transactions on Biomedical Engineering. - 2018. - V. 65, № 5. - P. 1181-1192.
14. Гираев К.М., Ашурбеков Н.А., Расулов М.Т. и др. Исследование влияния тер-мо-индуцированной коагуляции на структурные и спектральные свойства нормальных
и опухолевых биотканей методом светооптической и флуоресцентной микроскопии // Известия вузов. Сев.-Кав. регион. Сер.: Естественные науки. - 2011. - № 6. - С. 92-97.
Поступила в редакцию 15 февраля 2018 г.
UDC 535.33/.34 535.3
DOI: 10.21779/2542-0321-2018-33-2-28-39
Numerical modeling of the temperature spatial distribution in biological tissues throughout the process of laser heating
Israpov Ed.Kh.12, Giraev K.M.12, Magomedov M.A.12, Ashurbekov N.A.1
1 Dagestan State University; Russia, 367001, Makhachkala, M. Gadzhiev st., 43a; ed-gar1993@yandex.ru;
2 Institute of Physics, Dagestan Scientific Center, Russian Academy of Sciences; Russia, 367003, Makhachkala, Yaragskiy st., 94
The spatial distribution of temperature in the tissues of the brain gray matter as the laser hyperthermia processes develop is calculated by solving the non-stationary heat conduction equation. It was found that for the given conditions of laser irradiation upon reaching a temperature of 60.0 degree C the largest volume of photothermal damage 1.0 cm3 and the characteristic distance 0.5 mm at exposure 6.0 min are observed at a wavelength of 1064 nm. At other wavelengths (300 and 1460 nm), these parameters can be 2-4 times smaller. It is shown that the development of laser thermocoagulation leads to an additional reduction of the damage volume, distance and exposure 1.5-2 times at all the investigated wave lengths.
Keywords: biotissue, heat conduction equation, numerical simulation, laser thermocoagulation, exposure, damage volume.
Received 15 February, 2018
