Содержание к диссертации
Введение
ГЛАВА 1 Современные приложения мощных ультразвуковых пучков в медицине 14
1.1 Роль акустической нелинейности в режиме фокусировки мощного ультразвука 14
1.2 Слаборазреженные двумерные фазированные решетки для фокусировки мощного ультразвука 19
1.3 Акустические и температурные поля УЗ излучателей 24
1.4 Выводы 27
ГЛАВА 2 Нелинейное повышение эффективности нагрева биологической ткани мощным фокусированным ультразвуковым пучком 28
2.1 Постановка задачи 28
2.2 Акустическое поле 30
2.2.1 Теоретическая модель 30
2.2.2 Численный алгоритм 35
2.3 Температурное поле 39
2.3.1 Теоретическая модель 39
2.3.2 Численный алгоритм 42
2.4 Результаты 49
2.4.1 Влияние параметров ультразвука на процесс нагрева ткани 49
2.4.2 Эффекты акустической нелинейности при воздействии фокусированного ультразвука на ткань 51
2.4.3 Влияние преломления на границе вода-ткань 60
2.5 Выводы 67
ГЛАВА 3 Временная оптимизация режимов ультразвукового облучения ткани для повышения эффективности нагрева за счет акустической нелинейности 68
3.1 Описание эксперимента 69
3.1.1 Экспериментальная установка 69
3.1.2 Процедура калибровки излучателя 70
3.1.3 Подготовка установки и образцов печени 71
3.1.4 Исследуемые в эксперименте режимы 72
3.2 Численное моделирование 75
3.2.1 Акустическое поле 75
3.2.2 Трехмерное температурное поле 77
3.3 Результаты 81
3.4 Выводы 85
ГЛАВА 4 Двумерные фазированные решетки для пространственной оптимизации фокусировки мощного ультразвука в ткани 86
4.1 Постановка задачи 87
4.2 Теоретическая модель и численный алгоритм 89
4.3 Результаты 96
4.3.1 Влияние нерегулярного расположения элементов решетки на процесс нагрева ткани 96
4.3.2 Эффективность теплового воздействия на ткань при различной пространственной оптимизации ультразвукового облучения 103
4.4 Выводы 114
Заключение 115
Список литературы 117
Приложение 1 132
- Слаборазреженные двумерные фазированные решетки для фокусировки мощного ультразвука
- Температурное поле
- Численное моделирование
- Теоретическая модель и численный алгоритм
Введение к работе
Актуальность проблемы обусловлена важностью рассматриваемых задач о распространении нелинейных ультразвуковых волн в биологической ткани для многих современных медицинских приложений. Из наиболее перспективных новых направлений, интенсивно развивающихся в последнее время, можно выделить неинвазивное разрушение опухолей мощным фокусированным ультразвуком (высокотемпературная гипертермия или акустическая хирургия); остановку внутренних кровотечений (ультразвуковой гемостазис), направленный транспорт лекарств, разрушение почечных камней ударными импульсами (экстракорпоральная литотрипсия), визуализация биологических тканей по второй гармонике и с использованием контрастных агентов. Детальное понимание особенностей распространения мощного фокусированного ультразвука в биологической ткани и физических механизмов взаимодействия ультразвука с тканью представляет особую значимость для разработки, планирования и оптимизации протоколов облучения, а также для создания нового медицинского оборудования.
Эти задачи безусловно актуальны и с общефизической точки зрения, поскольку связаны с исследованием нелинейно-волновых явлений в фокусированных акустических пучках в средах с близким к линейному по частоте закону поглощения; пространственной структуры полей, создаваемых ультразвуковыми одиночными излучателями и фазированными решетками; теплового воздействия нелинейных волн на среду распространения. Построить аналитические решения для рассматриваемых в работе задач и развитых теоретических моделей удается лишь в исключительных случаях, поэтому большое значение имеет развитие методов численного моделирования. Это является самостоятельной проблемой диссертационной работы, актуальной для многих областей современной физики мощного ультразвука.
Способность мощного фокусированного ультразвукового пучка вызывать локальные разрушения в глубине органов без разрушения рядом лежащих тканей используется в системах ультразвуковой хирургии и гемостазиса, которые уже начали внедряться в клиническую практику. Основным механизмом разрушения ткани является ее быстрый наглев за
**ос национальная/
счет поглощения акустической энергии. Интенсивности ультразвука в фокусе таких систем достигают 103-104 Вт/см2. При таких интенсивностях большую роль играют эффекты акустической нелинейности, приводящие к каскадной генерации, гармоник исходной частоты и существенному увеличению эффективности нагрева, поскольку коэффициент поглощения в ткани увеличивается с частотой. Исследование роли акустической нелинейности при тепловом воздействии ультразвука на ткань в режимах акустической хирургии является, несомненно, важным.
Наряду с использованием одноэлементных пьезокерамических преобразователей, в последние годы.в УЗ хирургии все больший интерес привлекают к себе многоэлементные фазированные решетки. Они позволяют осуществлять электронное перемещение фокальной области по заданному объему ткани и создавать совокупности одиночных фокусов. Мощные акустические решетки обычно имеют большую апертуру (~10-13см) и слабую разреженность, размер элементов составляет несколько длин волн (~ 5Л, где Я ~ 1 мм). Дискретная структура решетки и относительно большой размер элементов способствуют появлению в акустическом поле при смещении основного фокуса в пространстве (или набора фокусов) сильно выраженных вторичных интерференционных максимумов, что является нежелательным для клинической практики. Существуют классические методы, позволяющие понизить уровень вторичных максимумов, например, уменьшение размера элементов решетки или использование разреженных решеток при отключении элементов решетки случайным образом. Недавно, в работах Гаврилова и Хэнда для мошлых ультразвуковых решеток был предложен метод, основанный на использовании случайного расположения элементов на поверхности. Этот метод был известен ранее, например, в радиолокации, однако, как оказалось, его использование не давало столь существенного выигрыша, как при конструировании двухмерных решеток для ультразвуковой хирургии. Задача выбора оптимального расположения элементов и режима облучения является многофакторной и требует детального исследования. Поэтому разработка теоретических моделей и численных алгоритмов, позволяющих описывать сложную структуру акустического и температурного поля в биологической ткани при
использовании мощных акустических решеток разной конфигурации, является важной- для: проектирования оптимальных терапевтических решеток и планирования эксперимента.
Цель работы:
Целью диссертационной работы являлось развитие общих подходов к исследованию различных режимов, теплового воздействия мощного ультразвука на биологическую ткань. В соответствии с поставленной целью было намечено решение следующих практически важных задач:
Разработка теоретических моделей и численных алгоритмов для описания акустических фокусированных полей с учетом нелинейных и дифракционных эффектов, поглощения и дисперсии скорости звука, слоистой неоднородности среды, а также температурных полей с учетом процесса диффузии тепла.
Экспериментальное и теоретическое исследование роли акустической нелинейности в процессе нагрева биологической ткани мощным фокусированным ультразвуком для различных режимов облучения.
Исследование пространственной структуры акустического и температурного полей различной геометрической конфигурации, создаваемых мощными ультразвуковыми решетками.
Научная новизна работы:
1. Развита новая численная модель, позволяющая решать задачу
теплового воздействия фокусированного мощного ультразвука на
биологическую ткань комплексно: с учетом дифракции акустической
волны, акустической нелинейности, степенного частотного закона
поглощения ткани, диффузии тепла, плоской слоистой неоднородности
ткани, пространственной и временной оптимизации воздействия.
Теоретически и экспериментально исследованы различные пространственные и временные режимы, формирования теплового разрушения биологической ткани в условиях высокотемпературной гипертермии.
Впервые выявлены основные особенности влияния акустической нелинейности на процесс нагрева и формирования теплового разрушения биологической ткани для различных режимов фокусировки ультразвука.
4. Впервые проведен анализ преимуществ использования двумерных многоэлементных слаборазреженных акустических решеток со случайным расположением элементов на поверхности для нагрева биологической ткани мощным фокусированным ультразвуком для различных режимов высокотемпературной гипертермии.
Научная и практическая значимость работы:
1. Разработаны оригинальные численные методы и соответствующие
пакеты программ, позволяющие моделировать процесс нагрева и
разрушения биологической ткани мощным фокусированным
ультразвуком.
Все результаты теоретических исследований, выполненных в диссертации, доведены до простых численных моделей, на основе которых с использованием современных пакетов программ (FORTRAN90, MATLAB) можно с наименьшими временными, затратами планировать протоколы облучения биологической ткани в клинической практике.
Развитые методики позволяют проводить анализ сложной пространственной структуры акустического: и температурного полей в биологической ткани для различной геометрической конфигурации и различного порядка пространственной дискретизации ультразвуковых излучателей.
Развитые в работе теоретические и численные модели дают хорошее количественное описание процесса нагрева и формирования теплового разрушения биологической ткани в режиме высокотемпературной гипертермии.
Знание особенностей процесса взаимодействия мощного ультразвука с биологической тканью необходимо для конструирования перспективных клинических, терапевтических и диагностических установок, где основной проблемой реализации является выбор оптимальных режимов облучения.
Положения, выносимые на защиту:
1. Развитие комплексной теоретической и численной модели,
описывающей тепловое воздействие мощного фокусированного
ультразвука на биологическую ткань с учетом дифракции,
акустической нелинейности, степенного частотного закона
поглощения биологической ткани, плоской слоистой неоднородности среды.
Экспериментальное и теоретическое исследование роли акустической нелинейности в процессе нагрева и формирования области теплового разрушения биологической ткани в режиме высокотемпературной гипертермии при. использовании одноэлементного поршневого фокусированного преобразователя.
Развитие теоретической модели и численного алгоритма для описания сложной структуры акустического и температурного поля двумерной многоэлементной слаборазреженной фазированной акустической решетки с регулярным или случайным расположением элементов на поверхности.
Обоснование преимуществ метода случайного расположения элементов на поверхности двумерной многоэлементной фазированной акустической решетки для подавления интерференционных вторичных максимумов при тепловом воздействии ультразвука на ткань. Анализ перспективных режимов. облучения, характерных для высокотемпературной гипертермии опухолей.
Апробация работы
Материалы диссертации докладывались на VI Международной конференции студентов и аспирантов по фундаментальным наукам "Ломоносов-99" (Физ. ф-т МГУ, Москва; 1999);.VI и VII Всероссийских школах - семинарах "Волновые явления в неоднородных средах" (Красновидово, Моск. обл., 1998 и 2000 г.г.); на Европейской Научной конференции (Франция, Лилль, 2000 г.); на 9-ом Межд. Конгрессе «Ультразвук в медицине и биологии» (Италия, Флоренция, 2000); на 2-й Международной конференции «Фундаментальные проблемы физики» (Саратов, 2000 г); на X и XI сессиях Российского Акустического Общества (Москва, 2000 и 2002 г.г.); на 1-ом Международном конгрессе «Применение фокусированного ультразвука высокой интенсивности в. медицине» (Китай, Чонкинг, 2001); на Международном симпозиуме по нелинейной акустике - ISNA16 (Россия, Москва, 2002); на 3-ем Международном симпозиуме по терапевтическому ультразвуку - ISTU3 (Франция, Лион, 2003); а так же обсуждались на семинарах в Институте исследования рака в Саттоне (Великобритания), на научном семинаре
С. А. Рыбака в Акустическом Институте им. Н. Н. Андреева, на семинаре в Научном Центре Волновых Исследований ИОФРАН, на научных семинарах кафедры акустики и кафедры общей физики и волновых процессов физического факультета МГУ.
Работа выполнена в рамках исследований,- проводимых по грантам РФФИ (№98-02-17318, 00-15-96530, 02-02-16999), гранту Президента РФ НШ-1575.003.2, по программам «Университеты России» (№1-5286), АФГИР (№RP2-2099, RP2_2384_MO_02), стипендий поддержки аспирантов Акустического Общества Америки (2000).
Публикации
Основные результаты диссертации изложены в 16 опубликованных печатных работах, список которых приводится в конце автореферата.
Структура и объем диссертации
Диссертация состоит из введения, четырех глав, заключения, двух приложений и списка цитируемой литературы из 151 наименования. Общий объем работы составляет 140 страниц, включающих 56 рисунков.
Слаборазреженные двумерные фазированные решетки для фокусировки мощного ультразвука
Интенсивное развитие современных методов воздействия мощного ультразвука на ткань направленно на получение эффективного, безболезненного, хорошо прогнозируемого быстрого лечения. С этой целью предлагается не только ряд новых режимов облучения, но разрабатываются и оптимизируются УЗ излучатели разной конфигурации [101]. В акустической хирургии для фокусировки мощного ультразвукового пучка используют либо одноэлементные преобразователи, излучающая поверхность которых по форме представляет собой сферическую вогнутую оболочку, изготовленную из пьезокерамического материала [6,10-10,72], либо многоэлементные акустические фазированные решетки [17,21,102,103]. Необходимость в перемещении фокального пятна в пространстве вызывает ряд трудностей при использовании одноэлементных преобразователей, фокус которых располагается вблизи центра кривизны оболочки. Это и сравнительно продолжительная по времени процедура лечения, и сложность регулирования положения зоны разрушения в ткани [6-Ю, 55]. Более универсальными УЗ устройствами являются двухмерные многоэлементные фазированные ультразвуковые решетки (пьезокерамические или пьезокомпозитные), которые позволяют формировать волновой фронт требуемой формы. А возможность одновременного формирования набора из одиночных фокусов и электронное смещение их в пространстве позволяет существенно облегчить и ускорить процедуру лечения [17,102-103]. Кроме того, фазированная решетка может использоваться одновременно и как излучатель мощного ультразвука, используемого для получения области разрушения ткани, и как визуализирующее устройство в режиме реального времени [103].
Однако наличие дискретной структуры решетки может способствовать появлению нежелательных вторичных интерференционных максимумов в акустическом поле, создаваемом решеткой [17,102,104]. Для того, чтобы разрушить приемлемый для практики объем опухолевой ткани ( 4 см3), двумерная фазированная решетка должна обеспечивать перемещение фокуса на требуемое расстояние в трех взаимно-перпендикулярных направлениях и при этом обладать акустической мощностью не менее 300-400 Вт [17]. Известно, что смещение фокусов в направлении перпендикулярном к оси регулярной решетки на расстояние 7-10 мм приводит к появлению вторичных интерференционных максимумов, уровень интенсивности в которых может достигать 50%-60% от максимума интенсивности в основных фокусах [21], а при более значительном смещении и превышать ее. Наличие подобных максимумов может привести к нежелательному перегреву и даже к разрушению ткани вне заданной области воздействия [41,53]. Безусловно, определяющим фактором при разработке фазированных решеток для применения в акустической хирургии является безопасность ультразвукового воздействия. Поэтому одним из основных критериев для акустических полей, создаваемых решеткой, является минимальный уровень вторичных интерференционных максимумов, находящихся вне фокальной области [21]. Реализация такой решетки представляет собой довольно сложную задачу, так как требует компромисса между рядом противоречивых требований. С одной стороны, чтобы увеличить расстояние, на которое можно перемещать фокальную область, а следовательно и расширить объем области разрушения, необходимо делать элементы менее направленным, т. е. уменьшать их размер [17]. С другой стороны, чтобы выполнить требования по излучаемой мощности при не слишком высоких значениях интенсивности на поверхности элементов ( 10 Вт/см), необходимо, чтобы активная площадь решетки была не менее 50 см2 (слаборазреженная) [17,21]. Известно, что чтобы существенно уменьшить влияние боковых лепестков и вторичных интерференционных максимумов в диаграмме направленности, расстояние между центрами элементов решетки должно быть Л/2 [22] {X - длина волны).
Это соответствует, например, 0.5 мм на частоте излучения 1.5 МГц, которая наиболее часто используется в эксперименте [17]. Однако, чтобы создать при столь малых размерах элементов решетку с достаточно большой апертурой и при этом реализовать требуемые для терапевтической решетки значения акустической мощности, необходимо использовать значительное число элементов и питающих их электронных каналов. Увеличение числа элементов с соответствующим уменьшением их размера на текущем уровне развития технологии является экономически нецелесообразным, а также вызывает ряд технических трудностей в процессе изготовления решетки [17]. Поэтому предлагается ряд других методов понижения уровня вторичных максимумов. Однако известные методы, основанные на уменьшении амплитуды на элементах решетки от ее центра к периферии (аподизация) [22] или расположении элементов с неравными расстояниями между ними [105] не всегда приемлемы в силу повышенных требований к акустической мощности решетки. Наиболее эффективный метод, используемый для понижения уровня вторичных интерференционных максимумов, основан на случайном пространственном расположении элементов на поверхности решетки [22]. Основанием для такого подхода является то обстоятельство, что уровень боковых лепестков в создаваемом решеткой поле существенно зависит от регулярности структуры решетки [22]. Используются либо решетки с «прореженными» элементами, случайным образом расположенными на поверхности решетки, т.е. используются не все элементы решетки, а лишь их определённая часть [106]. Однако повышение «прореженности» (активная площадь 40% от площади поверхности решетки) приводит к уменьшению излучаемой решеткой мощности и ухудшению качества распределений [17]. Либо используются решетки, у которых количество включенных элементов остается неизменным, однако их пространственное расположение на поверхности решетки носит случайный (рандомизированный) характер [17,23,107,108] (Рис, 5). Влияние случайного расположения элементов на структуру получаемого
Температурное поле
В настоящей работе для описания нестационарного теплового процесса используется неоднородное линейное уравнение теплопроводности, адаптированное для биологической ткани [35]: где Т=Т (r,t) -температура ткани, /-время, k,=K,/cvt -коэффициент температуропроводности, АГ, - коэффициент теплопроводности, Д - лапласиан по поперечным координатам, Ai /di +l/r д/дг в случае осесимметричного пучка [ИЗ], Т0=36.6 С- равновесная температура ткани. Первое слагаемое в уравнении (21) описывает охлаждение, обусловленное процессом диффузии тепла в среде. Второе - охлаждение, связанное с интенсивным теплообменом в кровеносных сосудах, находящихся как в самой нагреваемой области, так и вне ее [80-81, 128] (член объемного поглощения тепла [129]). Такой механизм охлаждения характеризуется временем перфузии ты= pbcvjfwcyb, где рь, сУъ, и w плотность, скорость течения и теплоемкость крови, соответственно, Cvi -теплоемкость единицы объема ткани. Индекс і указывает тип среды: 1 для воды и 2-для биологической ткани. Уравнение (21) позволяет рассчитывать температурные распределения для любой структуры тепловых источников gvfc г). Можно получить аналитическое решение уравнения (21) для точечного источника [130] где W = УЯАГя - мощность точечного излучателя, / = 2яр1/с0 р0 - интенсивность точечного излучателя, Е = exp(-r/Z), L=(kt)ll2 длина кровотока (L = (20-2) мм), Аналитическое решение (22) позволяет проводить качественное исследования влияния процесса перфузии и диффузии на нагрев ткани. В случае непрерывного распределения тепловых источников, температуру в любой точке пространства можно найти путем суммирования вкладов от всех тепловых источников.
В общем случае ведется интегрирование по всем направлениям, тогда аналитическое решение уравнения (21) записывается в виде: гдег = л]х +[z-z J. Решение (23) упрощается, если считать задачу цилиндрически симметричной Используя формулу (24), можно записать решение уравнения (21) в достаточно простой форме в случае нагрева ткани гауссовским пучком для выбранного времени воздействия. Однако в случае поршневого излучателя получить аналитическое решение уже не удается. В качестве «тестовой» задачи для отладки численного решения уравнения (21) (аппроксимации диффузионного оператора) используется аналитическое решение для начального сферического распределения температуры. Решение уравнения где g(4) - распределение температуры в начальный момент времени (t = 0), -пространственная координата. При подстановке соответствующего ядра K(t, , ) [132] в интеграл (26) аналитическое решение уравнения (25) записывается в виде: Akt+f2 В случае аксиально-симметричной задачи для отладки численного алгоритма можно также использовать аналитические решения уравнений, описывающих одномерные нестационарные диффузионные процессы вдоль оси пучка и в радиальном направлении (в фокальной плоскости): ) соответственно. Тогда при подстановке начального распределения g(4) = T0exp(-2/f2) и соответствующего для каждого из уравнений (27) ядра К (t, , ) [132] в интеграл (26) решения записываются в виде В случае слоисто-неоднородной среды типа вода-ткань, когда коэффициент температуропроводности в уравнение (21) не является постоянной величиной, условие теплопереноса на границе раздела может быть получено в предположении, что источников тепла на границе раздела нет [133].
Условие равенства тепловых потоков через границу записывается тогда в следующем виде: Конечный результат теплового ультразвукового воздействия обычно оценивается на основе развитой в работе [36] концепции тепловой дозы tdf математическая запись для которой имеет вид: Здесь („олнсе = tHaipta+tautened- полное время процедуры, Ttef- некоторая характерная температура, относительно которой рассчитывается доза, множитель R = 0.5 при Tt 43 С и 0.25 при Tt 43 С. За тепловую дозу, необходимую для создания теплового разрушения, обычно принимается доза, соответствующая поддержанию в ткани температуры 7 =43 С в течение 120-240 минут [36], Такое определение появилось в связи с работами по низкотемпературной гипертермии. Однако при фокусировке ультразвука высокой интенсивности используются значительно меньшие времена облучения (несколько секунд), в течение которых достигаются более высокие значения температуры 60 С в фокальной области, т.е. тепловая гибель клеток происходит почти мгновенно. Поэтому в данной работе удобнее выбрать в качестве точки отсчета Тге/ =56 С, что дает значение критической тепловой дозы tdjnj=\ с, соответствующее дозе tdTref= 136.6 минут относительно Тк/= 43 С,
Численное моделирование
Для численного моделирования процесса нагрева и формирования области разрушения ткани задача условно разбивается на два этапа. На первом этапе решается уравнение КЗК (1) при заданном граничном условии (6). Для режима формирования набора из восьми дискретных одиночных разрушений или «трека» задача не аксиально-симметричная, это приводит к необходимости расчета акустического поля в 3-мерном варианте. Однако такого рода расчеты требуют значительного компьютерного времени и большой памяти. Поэтому в настоящей работе получение 3-мерного поля тепловых источников проводится следующим образом. Вначале вычисляется полная интенсивность пучка с использованием алгоритма, описанного в главе 2. Полученное аксиально-симметричное поле интенсивности используется для вычисления двухмерного распределения тепловых источников (z, г) (12), которое конвертируется в трехмерное распределение (х, у, z) с помощью встроенных функций пакета MATLAB. Кратко опишем алгоритм процедуры конвертирования, которая условно разбивается на два шага. На первом - проводится конвертирование двухмерного массива поля тепловых источников с сетки с шагом вдоль оси пучка Az = 0.03 см и в радиальном направлении hr = 0.0042 см на сетку с более крупным шагом: hz = 0.1815 см и hr 0.00798 см, Более крупный шаг на новой сетке выбирается с целью сокращения затрат машинного времени при проведении расчета трехмерного поля температуры в дальнейшем. Пространственное окно интегрирования вдоль оси 0 z 27 см и в радиальном направлении 0 г 10.5 см не меняется. Для конвертирования массива, заданного на двухмерной сетке с меньшим пространственным шагом на двухмерную сетку с (a) (б) На втором шаге проводится конвертирование двумерного массива, заданного на новой двухмерной сетке с крупным пространственным шагом hz и hr в трехмерный массив на сетку с пространственным шагом: в продольном hz3D = hz = 0.1815 см и в поперечном Их = hy = hr = 0.00798 см направлениях. Процедуру перехода на 3-мерную сетку можно описать следующим образом. Рассмотрим двумерную сетку, показанную на Рис. 34о. Любому фиксированному значению z на этой сетке соответствует отрезок вдоль радиального направления г.
Длина этого отрезка равна, в нашей задаче, радиусу излучателя /. Построим квадрат со стороной 2/ и разобьем на нем сетку с шагом hr (Рис. 346), Получается сечение ху, для некоторого фиксированного значения z/. На следующем шаге каждому узлу (xi,yi) при фиксированном значении zt сопоставляется значение тепловых источников Q(ri) при том же значении z,-. Однако значения Q{ri) известны только в круге радиуса г , вписанного в квадрат со стороной 2/. Чтобы получить распределение поля тепловых источников на сетке во всем квадрате строится радиальное распределение для г гй V2 (Рис. 34е). В каждом узле r\ = rt V2 соответствующего радиального распределения экстраполируется новая функция поля тепловых источников Шаг круговой сетки совпадает с шагом квадратной сетки h г . Однако узлы этих сеток не совпадают. Так как отсчет в обоих сетках проводится от центра квадрата, то известно расстояние каждого узла от центра в обоих сетках. Зная значение функции Q\r ) в узле круговой сетки на расстоянии г/ значение функции Q (XJ, yi) на соответствующем расстоянии г, = у д-,2 + у] находится также с помощью процедуры линейной интерполяции. Описанный алгоритм получения распределения Q {х у]) в плоскости ху для некоторого фиксированного значения Zj, повторяется для каждого следующего значения z;+i- Полная процедура конвертирования двухмерного массива в трехмерный приводит к небольшому смещению (0,3 %) и понижению ( 4 %) пика в распределении поля тепловых источников по сравнению с имеющимся до процедуры конвертирования.
Для численного моделирования акустического поля значение интенсивности на источнике /о вычисляется на основе известной из эксперимента акустической мощности, подаваемой на излучатель (Таблица 4): /й = WacjSc, где Sc = 2n F Н = 56.52 см2 - площадь и Я- высота шарового сегмента. В Таблице 5 представлены значения исходной интенсивности Д , получаемые для различных режимов облучения. Расчет температурного поля для полученного 3-хмерного распределения тепловых источников проводится с использованием уравнения Пеннеса (21) в трехмерном представлении, когда лапласиан записывается в виде: Л = сР/дх2 + сР/ду2 + dVfe2-
Для описания температурного поля, получаемого при рассмотрении поршневого излучателя (6), задача решается численно с использованием конечно-разностного подхода. Для удобства вычислений уравнение (21) записывается в безразмерном виде Здесь 0 = t/t0 -время, нормированное na t0=\c, т= {Т-Т0)/Т0 -безразмерная температура, нормированная на равновесное значение температуры среды TQ = rconst, Z = zjF — продольная координата, нормированная на фокусное расстояние, Х=х/г0 и Y yfra — поперечные координаты, нормированные на радиус излучателя го, т6і =ты/{0 - безразмерное время перфузии. Индекс / указывает тип среды: 1 для воды и 2-для биологической ткани. Три безразмерных параметра at = kt to/F2 , Р,=к, t0/r и у і = A o/ci.j o F характеризуют диффузию тепла: вдоль оси пучка (а,-), в поперечных направлениях (fit), вклад тепловых источников в процесс нагрева (у;). Для решения задачи используется метод расщепления по физическим факторам Ь-Ь фф Л-и. На каждом временном шаге Ы при переходе от временного слоя $ к слою &+ht задача решается в два этапа, описанных в главе 2. Полученное для части, учитывающей влияние перфузии и тепловые источников, точное решение берется в качестве начального температурного распределения для решения диффузионной части Ь фуз. На втором этапе решается однородное уравнение диффузии
Теоретическая модель и численный алгоритм
Расчет акустического поля, создаваемого решеткой в ткани при заданном расположении элементов и при заданном наборе контрольных точек (фокусов) в пространстве, можно условно разбить на три этапа: расчет поля одиночного элемента решетки, определение оптимального набора фаз с последующим выравниванием абсолютных значений амплитуд на элементах и нахождение поля решетки путем суммирования полей всех элементов с найденным амплитудно-фазовым распределением. На первом этапе для описания акустического поля круглого поршневого элемента решетки используется интеграл Релея-Зоммерфельда [103, 150] где р - комплексная амплитуда давления, к = In -fQlcQ -волновое число, а-7.5 м 1 - коэффициент затухания на рабочей частоте решетки /о, и0 амплитуда колебательной скорости на поверхности элемента, -г, расстояние от точек г излучающей поверхности 5 круглого поршневого элемента до ї-ой точки в пространстве rf. Аналитическое решение выражения (46) удается найти только на оси излучателя и в фокальной плоскости. Для нахождения решения на оси, запишем выражение (46) в полярной системе координат где p- расстояние от центра источника до любой точки излучателя, р - полярный угол, d S = pdpd p - элемент площади. После процедуры интегрирования получаем выражение, описывающее поле давления круглого поршневого элемента вдоль оси которое используется в работе для отладки численного решения уравнения (46). Интеграл Рэлея-Зоммерфельда (46) рассчитывается численно в узлах достаточно частой пространственной сетки. Используются следующие значения параметров моделирования: область интегрирования по продольной координате 2 z 16 см, пространственное окно по поперечной координате 0 г 10 см, шаг сетки вдоль обоих направлений hz = hr= 0.2 мм, что составляет 1/5 длины волны. Для численного интегрирования поверхность элемента разбивается на 7849 квадратных участков со стороной 0.05 мм (-V20), что обеспечивает достаточную точность решения.
Полученное решение используется в дальнейшем для нахождения фаз на элементах и при расчете полного поля, создаваемого решеткой. На втором этапе решается обратная задача. Определяются значения комплексной амплитуды колебательной скорости на каждом из элементов, которые позволили бы сформировать в пространстве набор N=9 физических фокусов (контрольных точек) с заданными координатами. Накладываются дополнительные условия на выбор амплитуды и фазы давления в каждой контрольной точке: фазы равномерно вращаются по часовой стрелке относительно оси данной совокупности фокусов [151], с шагом ли. (Рис, 42), значения амплитуды давления в контрольных точках равны между собой. Основы метода расчета комплексной амплитуды колебательной скорости на элементах были предложены в работе [103], однако лаконичность изложения затрудняет его практическое использование. Поэтому ниже представлены основные составляющие данного подхода в двух его модификациях. В одной из них при вычислении амплитуды и фазы колебательной скорости на элементах используется акустическое поле (46) поршневого элемента [41], В другой модификации поле поршня заменяется полем точечного источника, расположенного в центре элемента [17]. Здесь рассмотрены обе модификации с целью сравнения получаемых наборов фаз и акустических полей, создаваемых решеткой.
В приближении линейного распространения акустической волны давление р(гп) в каждой и-ой контрольной точке (п = 1 N) является суперпозицией М парциальных давлений рт (т= 1,.,.,ЛҐ). создаваемых каждым из М элементов решетки р(7л )= Нпт ит , где ит — комплексная амплитуда колебательной скорости на поверхности /м-го элемента т = {1,..., М). Или, в матричном виде: Здесь U— вектор из М значений амплитуд колебательной скорости ит, Р - вектор давлений р(гп) в JV контрольных точках, матрица Н размерностью (NxM) определяет оператор прямого распространения от т-го элемента на решетке до я-й контрольной точки в пространстве. Если для решения (47) используется поле поршневого элемента с равномерным распределением колебательной скорости, то элементы матрицы Н„т определяются решением (46), рассчитанным для ти-го элемента решетки в л-й контрольной точке: Заметим, что тогда для нахождения значений элементов матрицы Нпт (48), необходимо рассчитать интеграл Рэлея-Зоммерфельда N-M раз (в N фокусах от М элементов решетки). Однако можно использовать рассчитанное для одного из элементов поле (46) для нахождения матричных элементов (48) от произвольного элемента т в произвольной контрольной точке п. Для этого определяется положение контрольной точки относительно цилиндрической системы координат, связанной с выбранным элементом (Рис.43), и находится взвешенное среднее значение поля в ближайших узлах, если положение точки попадает внутрь интервалов hz или кг сетки. Умножая далее полученный результат на комплексную амплитуду колебательной скорости элемента, можно получить значение поля давления, создаваемого элементом т в контрольной точке п. Если же использовать поле точечного источника с координатой fmt то элементы матрицы #„т имеют более простой вид [21]
