На правах рукописи
Прохоров Игорь Васильевич Математические задачи теории переноса излучения 05.13.18 - математическое моделирование, численные методы и комплексы программ
АВТОРЕФЕРАТ
диссертации на соискание ученой степени доктора физико-математических наук
ВЛАДИВОСТОК - 2007
Работа выполнена в лаборатории вычислительных методов математической физики Института прикладной математики Дальневосточного Отделения Российской Академии наук
Научный консультант: доктор физико-математических наук, профессор Аниконов Дмитрий Сергеевич
Официальные оппоненты: доктор физико-математических наук, профессор Прилепко Алексей Иванович доктор физико-математических наук, член-корреспондент РАН, профессор Смагин Сергей Иванович доктор физико-математических наук, старший научный сотрудник Ярощук Игорь Олегович
Ведущая организация: Институт математики имени С.Л. Соболева Сибирского Отделения Российской Академии наук, г. Новосибирск
Защита состоится 2 ноября 2007г. в 10 часов на заседании диссертационного совета Д005.007.01 в Институте автоматики и процессов управления ДВО РАН по адресу: 690041, г.Владивосток, ул. Радио, 5.
С диссертацией можно ознакомиться в библиотеке Института автоматики и процессов управления ДВО РАН.
Автореферат разослан
Ученый секретарь диссертационного совета Д005.007.01 А.В. Лебедев
Актуальность темы и объект исследования. Широкий круг физических процессов и явлений в природе описываются кинетическими уравнениями переноса. С точки зрения практики представляет интерес исследование линейных интегродифференциальных уравнений в частных производных первого порядка, моделирующих рассеяние света в атмосфере и биологических тканях, распространение рентгеновских лучей и гамма-квантов в неоднородных средах. Эти уравнения используются также при описании диффузии нейтронов в веществе, в теории плазмы, в газо-кинетических моделях и в теории распространения звука.
Большой вклад в исследование краевых задач для интегродифференциальных уравнений переноса и в создание математически обоснованных вычислительных методов их решения внесли работы В.С. Владимирова, Е.С.
Кузнецова, Г.И. Марчука, С.Б. Шихова, В.И. Лебедева, М.В. Масленникова, Т.А. Гермогеновой, Г.А. Михайлова, В.И. Агошкова, Т.А. Сушкевич, В.В.
Смелова, A. Исимару, К. Кейза, С. Чандрасекхара и др.
При описании процесса распространения излучения в веществе приходится сталкиваться с тем, что среда составлена из разнородных по своим радиационным свойствам материалов. В терминах уравнения переноса это означает разрывность его коэффициентов. На поверхностях разрыва ставятся условия сопряжения, которые накладывают дополнительные ограничения на функцию распределения плотности излучения. Наиболее полно исследованы задачи уравнения переноса с условиями сопряжения типа непрерывной склейки решения на границе раздела сред. Даже в этом одном из самых простых случаев локальная структура решения оказывается весьма сложной и на этапе получения общих качественных утверждений многие трудности обходятся рассмотрением достаточно широких классов решений.
Рассмотрение более общих условий сопряжения, позволяет учитывать отражение и преломление света на контактных границах, а не только его поглощение и рассеяние на случайно распределенных микро-неоднородностях среды. Модель, основанная на уравнения переноса с обобщенными условиями сопряжения, весьма актуальна и находит свое применение в атмосферной оптике, в лазерной медицине и в трехмерной компьютерной графике. Несмотря на ее более чем полувековую историю и повышенный интерес в последнее время, теоретические результаты в этой области развиты в незначительной степени. По всей видимости, это объясняется увеличивающейся сложностью поставленных задач.
С еще большими проблемами приходится сталкиваться при решении обратных и экстремальных задач теории переноса излучения, имеющих многочисленные приложения в медицине, физике и технике. На сегодняшний день теория обратных задач весьма обширна и развита в работах Г.И.Марчука, М.В.Масленникова, А.И. Прилепко, Д.С. Аниконова, Г.А. Михайлова, Ю.Е.
Аниконова, В.Г. Романова, В.А.Шарафутдинова, В.Р. Кирейтова, У.Н. Султангазина и других. Теория обратных задач отличается большим разнообразием постановок задач, различными ограничениями и методами.
Экстремальные задачи теории переноса излучения в основном рассматривались в рамках волновых моделей, применительно к теории синтеза оптических многослойных систем (А.Н. Тихонов, А.Г. Свешников, В.Б. Гласко, А.В.
Тихонравов и др.).
Работа посвящена исследованию прямых, обратных и экстремальных задач для стационарного линейного уравнения переноса излучения с различными условиями сопряжения на границах раздела сред.
Цель работы. Исследовать качественные свойства решения прямых задач для стационарного уравнения переноса. Сформулировать необходимые и достаточные условия единственности решения обратных задач для уравнения переноса и разработать соответствующие численные алгоритмы. Рассмотреть проблему постановки и решения экстремальных задач в рамках кинетического подхода к описанию процесса переноса излучения.
Методы исследования. При изучении качественных свойств решения краевых задач применялись методы теории интегральных уравнений и дифференциальных уравнений в частных производных. Исследование обратных задач проводилось в рамках теории условно-корректных задач и существенно основывалось на изучении локальных свойств решения прямой задачи.
При построение численных алгоритмов и проведении экспериментов использовались методы компьютерного моделирования и некоторые разделы теории методов Монте-Карло.
Научная новизна. В диссертационной работе получены следующие новые результаты.
Построена качественная теория краевых задач для стационарного уравнения переноса излучения с обобщенными условиями сопряжения на границе раздела сред.
В рамках построенной теории, разработан и обоснован метод решения задачи томографии, позволяющий определять отражающую и преломляющую границу раздела сред по потоку прошедшему через рассеивающий объект неизвестной структуры.
Для задачи томографии в модели, не учитывающей эффекты отражения и преломления, получены необходимые условия единственности решения.
Предложены новые подходы для нахождения оптимального решения задачи томографии, используя в качестве управляющей переменной энергетическую зависимость внешних источников излучения.
Для полихроматического уравнения переноса разработан метод нахождения коэффициента ослабления по излучению известному на границе, основанный на использовании специального типа внешнего излучателя.
Предложен новый подход к исследованию теоретических проблем оптики просветляющих и маскировочных покрытий случайно-неоднородных сред, рассматривающий указанные проблемы с позиции решения экстремальных задач для уравнения переноса с френелевскими условиями сопряжения.
Разработаны и программно реализованы следующие численные алгоритмы: для определения границ неоднородностей по излучению прошедшему через среду; для решения оптимизационных задач рентгеновской томографии и оптической маскировки; для нахождения коэффициента ослабления уравнения переноса при разрывном типе внешнего источника.
Теоретическая и практическая ценность. При исследовании краевых задач для стационарного уравнения переноса излучения получены теоретически значимые новые результаты. Изученные свойства решений прямых задач будут полезны при разработке и обосновании эффективности численных алгоритмов в атмосферной оптике, астрофизике, компьютерной графике. Предложенные методы решения обратных и экстремальных задач для уравнения переноса допускают практическое применение в неразрушающем контроле промышленных изделий, в медицинской томографии и биооптике.
Исследования диссертации, позволили создать "Базу данных радиационных характеристик веществ, представляющих интерес в рентгенодиагностике"( и запатентовать способ маскировки изделий [21], которые будут полезны специалистам при интерпретации данных, полученных при рентгеновском облучении веществ с целью определения их внутренней структуры или в связи с созданием материалов с заданными радиационными свойствами.
Исследования автора были поддержаны РФФИ (проекты №№97-01-00078а, 01-01-00128-а, 04-01-00126, 05-07-90055-в), программами ведущих научных школ Российской Федерации (проект НШ-9004.2006.1), "Университеты России" (проект УР.03.01.002) и "Минобразования России" (проект Е02-1.0-128), грантами: 6-го конкурса-экспертизы научных проектов молодых ученых РАН (проект №4), Конкурса проектов ДВО РАН (проекты №№ 05-III-А-01-101, 06III-А-01-014), Конкурса интеграционных проектов Дальневосточного, Сибирского и Уральского отделений РАН (проекты №№04-2-1-00-006, 06-II-СУ-01001).
Апробация работы. Основные результаты диссертации были представлены: на Всесоюзной конференции "Условно-корректные задачи" (Новосибирск, 1992), на международном симпозиуме по компьютерной томографии (Новосибирск, 1993), на Российско-японском семинаре "Дифференциальные уравнения в прикладной математики" (Хабаровск, 1994), на международных конференциях "Inverse Problems in Engineering Sciences" (IPES-94, Osaka, Japan, 1994), "Mathematical Modeling and Cryptography" (PMMC-95, Владивосток, 1995), "Computational Radiology and Imaging: Therapy and Diagnostics" (University of Minnesota, USA, 1997), "Recent Developments in Theories and Numerics" (Hong Kong, China, 2002), "Тихонов и современная математика" (Москва, МГУ, 2006), на Дальневосточных математических школах-семинарах им. академика Е.В. Золотова (Владивосток, Хабаровск, 1998Ц2006).
По материалам диссертации автор выступал на семинарах лаборатории обратных задач математической физики Института математики им. С.Л. Соболева СО РАН (1993, рук. д.ф.-м.н. Ю.Е. Аниконов); отдела статистического моделирования в физике Института вычислительной математики и математической геофизики СО РАН (1993, рук. чл.-корр. РАН Г.А. Михайлов); отдела условно-корректных задач ИМ им. С.Л. Соболева СО РАН (2004, рук.
академик М.М. Лаврентьев, чл.-корр. РАН В.Г. Романов); лаборатории математического моделирования фазовых переходов Института гидродинамики им. М.А. Лаврентьева СО РАН (2004, рук. чл.-корр. РАН П.И. Плотников);
на семинаре "Обратные задачи анализа, математической физики и естествознания" кафедры математического анализа МГУ (2005, рук. д.ф.-м.н. А.И.
Прилепко); на семинаре по "Математической физике" ИПМ им. М.В. Келдыша РАН (2005, рук. д.ф.-м.н. М.В. Масленников); на расширенном семинаре лаборатории механики деформируемого твердого тела ИАПУ ДВО РАН (2007, рук. д.ф.-м.н. А.А. Буренин), на общеинститутских семинарах ИПМ ДВО РАН (2004, 2007, рук. чл.-корр. РАН Н.В. Кузнецов).
Публикации. Всего по теме диссертации опубликовано 45 работ, том числе три монографии [1-3], патент [21] и 13 статей [4,8,10-14,16-18,22-24] в журналах, рекомендованных ВАК РФ для публикации результатов докторских диссертаций.
ичный вклад автора. Из результатов совместных работ в диссертацию включены только те, которые получены непосредственно автором, либо те, в которых степень участия соавторов в исследовании равная. Результаты, представленные во второй главе диссертации, получены совместно с Д.С.
Аниконовым на паритетных началах.
Структура и объем работы. Диссертация состоит из введения, 27-ми параграфов, структурно разделенных на пять глав, заключения и списка литературы из 266 наименования. Работа изложена на 256 страницах и содержит 18 рисунков.
Содержание работы Во введении обоснована актуальность темы диссертационной работы, приведен обзор литературы по вопросам связанным с диссертацией и коротко перечислены результаты, изложенные в остальных главах.
В главе I исследуются краевые задачи для полихроматического уравнения переноса:
rf(r, , E) + (r, E)f(r, , E) = Emax = K(r, , E, E )f(r, , E )d dE + J(r, , E), (1) Emin где f(r, , E) трактуется как плотность потока частиц (плотность потока излучения) в точке r, r = (r1, r2, r3) G R3, движущихся в направлении единичного вектора , = (1, 2, 3) = { R3 : || = 1} и имеющих энергию E, E [Emin, Emax]; (r, E) - коэффициент ослабления в точке r частиц с энергией E; K(r, , E, E ) - дифференциальное сечение рассеяния в точке r, характеризующее вероятность перехода частиц из состояния (, E ) в (, E); J(r, , E) - плотность внутренних источников в точке r, в направлении и на энергетическом уровне E.
Процесс миграции фотонов происходит в среде, заполняющей ограниченную выпуклую трехмерную область G. Область G разбивается поверхностью на конечное число подобластей G1, G2,..., Gp (рис. 1), объединение которых обозначим через G0. Множества Gi интерпретируются как некоторые части неоднородной среды G, заполненные i-м веществом. Поверхность является границей раздела сред G1,..., Gp и называется внутренней границей G0, а Рис. 1.
G - внешней. Множество G0 предполагается обобщенно выпуклым (любая прямая, имеющая общую точку с G0, пересекает границу G0 в конечном числе точек).
Аналогично вводятся конечные разбиения 1, 2,... и I1, I2,... попарнонепересекающихся открытых областей, объединения которых 0 и I0, определяют множества непрерывного изменения параметров излучения по угловой и спектральной переменной E, соответственно.
Прямая задача для уравнения переноса заключается в нахождении функции f из уравнения (1) и граничного условия f(z, , E) = h(z, , E), (z, , E) -, (2) в которых известны , K, J и функция h, интерпретируемая как плотность потока излучения проникающего в среду G через поверхность G. На множествах + и + заданы выходящий и выходящий потоки, и они определяются через расстояние d(r, ) от точки r до G в направлении вектора соотношением = {(z, , E) G I : z = r d(r, ), r G0} (рис. 1).
Подмножество есть пересечение и G 0 I0.
Условия на границе раздела записываются с помощью оператора сопря жения B:
f-(z, , E) = (Bf+)(z, , E), (z, , E) [Emin, Emax], (3) где f(z, , E) = lim f(z, , E) предельные значения функции f(r, , E) +при соответствующем стремлении точки r к z . В случае непрерывности решения при переходе через оператор сопряжения имеет вид:
(Bf+)(z, , E) = f+(z, , E). (4) Условия (3) с оператором сопряжения вида (4) достаточно традиционны и адекватно описывают процесс на границах раздела сред при распространении рентгеновского излучения и гамма-квантов. В первых двух главах диссертации рассматриваются задачи с оператором сопряжения вида (4).
Решение прямой задачи ищется в пространстве Dc(G0 0 I0), которому принадлежат функции f(r, , E), удовлетворяющие условиям: при любых (r, , E) G0 0 I0 функция f(r + t, , E) абсолютно непрерывна по переменной t, t [-d(r, -), d(r, )]; функции f и rf принадлежат Cb(G0 0 I0), где Cb(X) - пространство ограниченных и непрерывных функции на множестве X и символ rf(r, , E) понимается как производная по переменной r в направлении вектора .
Из определения пространства решений Dc, вытекает, что условие (3) на поверхности с оператором сопряжения в форме (4) (условие непрерывной склейки) автоматически удовлетворяется. Это существенно облегчает исследование прямой задачи и позволяет свести исходную дифференциальную задачу (1),(2) к эквивалентной интегральной постановке, причем в более широком классе, не заботясь о поведении следов функций из Dc на поверхности .
В предположении, что функции , K, J, h неотрицательные и Cb(G0 I0), J Cb(G0 I0), h Cb(-), и некоторых естественных ограничениях относительно функции K(r, , E, E ), в зз 1 - 4 доказаны теоремы существования и единственности решения и проведен анализ его непрерывных свойств, которые впоследствии используются при изучении обратной задачи в зз 5, 6.
В з 5 ставится задача об определении коэффициента (r, E) из соотношений (1), (2) и условия f(z, , E) = H(z, , E), (z, , E) +, (5) в которых известны только функции h и H. Функцию H, по аналогии с h, иногда называют плотностью потока выходящего из среды G излучения. В обратной задаче функции K, J не предполагаются известными и в то же время не подлежат определению в этой задаче.
В з 6 для решения обратной задачи применяется специальный тип внешнего источника, который в терминах модели описывается функцией h. В предположении, что она имеет разрыв первого рода по угловой переменной при = (1, 2, 0) доказана формула, связывающая преобразование Радона от функции в горизонтальном сечении с функциями h и H.
Теорема 1. Пусть 0 = 12, где 1 = { : sgn(3) = 1}, 2 = { : sgn(3) = -1}, функция h(z, , E) Cb(-) и для всех (z, , E) -\0 удовлетворяет условию h(z, , E) = lim h(z, , E) - lim h(z, , E) = 0. (6) 1, 2, тогда для всех (r, , E) G0 \ 0 I0 справедливо равенство d(r,) h(r - d(r, -), , E) (r + t, E)dt = ln. (7) H(r + d(r, ), , E) -d(r,-) Используя формулу (7), доказана теорема единственности решения обратной задачи. Из (7) видно, что неизвестные функции K и J не влияют на процедуру определения функции . В физических терминах это означает, что влияние рассеяния и наличие радиоактивных источников в среде подавляется за счет выбора внешнего источника. В диссертации обсуждаются возможности использования формулы (7) в спектротомографии - сравнительно новой области в радиодиагностике сред.
(a) (b) (c) (d) Рис. 2: Восстановление фантома Кормака в горизонтальном сечении с использованием параллельной схемы сканирования: (a) - оригинал; (b),(c),(d) - реконструкция среды при ухудшающейся угловой дискретизации выходящего излучения.
В з 7 обосновывается принципиальная возможность создания требуемого внешнего источника. Приводятся достаточные условия на структуру среды такие, что при просвечивании среды непрерывно-распределенным потоком выходящее излучение будет иметь заданный разрыв в горизонтальном направлении.
Завершает первую главу з 8, посвященный проведению компьютерных экспериментов, тестирующих предложенный численный алгоритм определения функции . Численные эксперименты проводились в три этапа. Сначала производилось вычисление функции H. Для этого использовалась одна из версий метода Монте-Карло, называемая методом сопряженных блужданий. На практике прошедшее через среду излучение измеряется детекторами и этот трудоемкий этап отсутствует. На втором этапе, используя формулу (7), находилось преобразование Радона от функции в некотором горизонтальном сечении. Затем, применяя для обращения преобразования Радона алгоритм свертки и обратной проекции, определялась функция в рассматриваемом сечении. Показана эффективность предложенного алгоритма для сильно рассеивающих сред и продемонстрировано влияние ошибок дискретизации и неполноты данных на качество реконструкции функции . На рис. приведены типичные томограммы, полученные в результате тестирования алгоритма.
Вторая глава посвящена моноэнергетическому варианту уравнения переноса:
rf(r, , E) + (r, E)f(r, , E) = = K(r, , E)f(r, , E)d + J(r, , E). (1 ) В этом случае переменная E в уравнении (1 ) играет роль параметра и везде, кроме з 13 и частично з 11, зависимость от него опущена.
Дифференциальное сечение рассеяния K в (1 ) имеет вид M K(r, ) = i(r)Si( ), Si( )d = 1, i = 1,..., M, (8) i= где каждый член в сумме соответствует определенному типу рассеяния, а M коэффициенты рассеяния и поглощения даются формулами s(r) = i(r) i=и a(r) = (r) - s(r).
В з 9 утверждения, относящиеся к прямой задаче и доказанные в главе I, конкретизируются для моноэнергетического случая.
В з 10 ставится задача томографии, заключающаяся в нахождение поверхности разрывов коэффициентов уравнения (1 ), если известна только плотность выходящего потока H на границе среды G.
Рассмотрим функцию d(z,) - (z + t)dt m(z, ), mv(z, ) = exp (9) M m(z, ) = [[(z)]]f(z, )- [[i(r)]] Si( )f(z, )d +[[J(r, )]]. (10) i= где через [[]] обозначен скачок функции по пространственной переменной r в точке контакта двух областей Gj и Gi.
Величина |mv(z, )|, называется мерой видимости среды G в точке z G0 в направлении при ее радиодиагностике рентгеновским излучением.
Д.С. Аниконовым в [1] показано, что условие видимости |mv(z, )| > 0, z G0, гарантирует единственность решения задачи томографии.
В з 11 исследуется необходимость условия видимости для задач томографии вообще и, в частности, для рассматриваемой задачи. Имеет место Теорема 2. Пусть существуют функции , K, J, h и соответствующее решение прямой задачи f такие, что:
a) [[(z)]] = 0 или [[s(z)]] = 0 для всех z Gj, где j некоторый фиксированный номер 1 j p;
b)J(r, ) = 0 для всех (r, ) Gj и [[J(z, )]] = 0, (z, ) Gj ;
c) [[i(z)]][[k(z)]] 0, i, k = 1,..., M, z Gj;
d) mv(z, ) = 0 для всех (z, ) Gj .
Тогда: [[a(z)]] = [[(z)]] - [[s(z)]] = 0, z Gj и существует бесконечно много различных функций , s и областей Gj заменяющих , s и Gj, со ответственно, для которых решение прямых задач f(r, ) одно и то же везде в G и совпадает с f(r, ).
Предположение a) теоремы 2, означает, что существует по крайней мере одно включение (Gj). Пункт b) накладывает ограничения на внутренние источники. В частности, случай J 0 подходит под условие b). Неравенство c) можно трактовать как преобладание рассеяния в одной среде над другим, например каждое сечение i(r), i = 1,..., M в Gl больше чем i(r) в Gj. Условие c) также выполняется в часто встречающемся на практике случае M = 1, то есть когда рассматривается только один тип рассеяния.
Предположение d) наиболее важное в теореме, поскольку именно оно приводит к невозможности решения томографической проблемы ни каким методом в рамках изучаемой модели. Среду G, имеющую свойство mv(z, ) = 0, z G0, (11) будем называть невидимой при ее диагностике радиационным потоком. Последнее определение зависит от многих параметров процесса и, следовательно, его трудно использовать на практике. В работе показано, что [[a(z)]] = является аппроксимацией условия (11) при некоторых предположениях, которые часто имеют место в томографии.
В связи с этим, среда G, имеющая свойство [[a(z)]] = 0, z G0 называется плоховидимой при ее томографическом просвечивании. Равенство [[a(z)]] = 0 означает непрерывность коэффициента поглощения на границе включения Gj. Два вещества называются плохо-видимой парой на уровне энергии E, если они имеют одно и тоже значение коэффициента поглощения a(E). В монографии [3] приведены таблицы плоховидимых пар для 92 химических элементов и 48 веществ. Эти таблицы были сформированы на основе спектральных данных о коэффициентах ослабления и поглощения, предоставленных Дж. Хабблом и С. Сельтцзером. Расширенный вариант этих таблиц содержится в базе данных "sxray.iam.dvo.ru".
В этом же параграфе дан алгоритм решения задачи томографии, основанный на построение индикатора неоднородности Ind(r) и предназначенного для послойного восстановления внутренней структуры среды (в сечении P = {r3 = const}). Он представляет собой интегродифференциальный оператор вида:
f(r + d(r, ), )d, Ind(r) = (12) r где единичная окружность 1 состоит из векторов = (1, 2, 0), а символ означает градиент по r1, r2 при r3 = const.
r Д.С Аниконовым и В.Г. Назаровым в [3] доказано, что в достаточно малой окрестности контактной точки z Gl Gj P имеет место представление:
Ind(r) = 2M(z)| ln |r - z|| + O(1) (13) где O(1) - ограниченная функция, а M(z, E) определена равенством:
M(z) = |mv(z, ) + mv(z, -)|. (14) Здесь функция mv определена в (9), а вектор 1 является ортогональным к нормали nl(z) к кривой Gl P. Поскольку таких векторов всего два и оба симметрично представлены в (14), то левая часть в последнем равенстве не зависит от .
Если M(z) > 0, то из (13) видно, что функция Ind(r) неограничена только вблизи G0. При численной реализации алгоритма, построенного по формуле (12), поверхность G0 определяется как место, где значения индикатора значительно больше других.
Результаты, проведены численных экспериментов, подтвердили полезность введения понятий меры видимости, невидимых и плоховидимых сред. Кроме того, проделан вычислительный эксперимент, показывающий, что плоховидимые и тем более невидимые среды - это особые состояния, которые не всегда подходят для использования некоторых упрощенных моделей процесса распространения излучения.
На рис. 3a изображен оригинал структуры среды в горизонтальном сечении, на рис. 3b - реконструкция, сделанная с помощью индикатора неоднородности (12). Улучшение качества реконструкции происходит при увеличении скачка коэффициента поглощения на границах включений. Причем, граница плоховидимого включения G1, где наблюдается непрерывность коэффициента поглощения, вообще "не видна". На рис. 3c представлена реконструкция среды, если выходящее излучение вычисляется по упрощенной модели, в которой не учитывается наличие 4-х процентного рассеяния в среде. В этом случае, границы неоднородностей G1 и G3 восстановились отчетливо, а G2 - нет. Это показывает, что использование такой упрощенной модели может приводить как к потери информации о структуре среды, так и к появлению изображений, не имеющих место при реконструкции с использованием моделей, более адекватно описывающих реальный процесс переноса излучения.
(a) (b) (c) Рис. 3.
В з 12 построены примеры неединственности решения задачи томографии.
Показано, что для них выполняется условие невидимости (11).
з 13 посвящен постановке и изучению задачи оптимизации в рентгеновской томографии. Рассматривается задача определения внутренних границ между различными материалами в произвольной среде. Исследования предыдущих параграфов позволяют не только решить поставленную задачу, но и предложить подход к проблеме выбора наилучших условий рентгенодиагностики.
В данном параграфе изучается экстремальная задача, заключающаяся в определении энергии E источников h(z, , E), которая обеспечивает наилучшее качество реконструкции искомых границ.
Все характеристики процесса имеют параметрическую зависимость от E, поэтому от E зависит и функция M(z, E), определенная в (14). Хотя теоретическое обоснование алгоритма решения задачи томографии, основанного на построение оператора (12), справедливо при любом положительном M(z, E), даже сколь угодно малом, для численной реализации величина M(z, E) довольно существенна. Из теоретических и численных результатов следует вывод: чем больше значение M(z, E), тем выше качество (отчетливость) реконструкции.
В связи с этим возникает задача выбора условий, обеспечивающих наилучшее качество восстановления искомых границ. Предварительный анализ показал, что ее конкретные постановки могут иметь много разновидностей.
В диссертации рассматривается следующий вариант: найти значение энергии E из заданного интервала [Emin, Emax], доставляющее наибольшее значение целевой функции M(z, E) при фиксированном z.
Для случая табличного задания величин , i, Si, J, h, что обычно принято на практике, решение этой задачи легко реализуемо путем численного перебора. Более экономичному варианту решения препятствует прежде всего недостаточный уровень исследования свойств функции f(r, , E). Знание M(z, E) подразумевает известность всех других характеристик процесса переноса, что для дефектоскопии неестественно. Эти обстоятельства ставят вопрос о такой аппроксимации целевой функции, которая зависела бы от меньшего числа параметров и имела бы более простой вид, удобный для теоретического анализа. В связи с этим, в этом параграфе предлагается и исследуется упрощенная оптимизационная задача поиска энергии, близкой к оптимальной в случае малости дефекта. Результаты численных экспериментов полностью подтвердили ожидания, основанные на теоретическом анализе. Томограммы, полученные при оптимальных значениях параметра E, соответствовали наи лучшей реконструкции среды. Целевая функция M(z, E), введенная в задаче квазиоптимизации, достаточно хорошо аппроксимировала функцию M(z, E) в случае малого дефекта. Алгоритм решения задачи квазиоптимизации программно реализован в базе данных "sxray.iam.dvo.ru".
Следующие три главы посвящены исследованию краевых задач для монохроматического уравнения переноса с обобщенными условиями сопряжения.
В главе III рассматривается прямая задача для уравнения переноса rf(r, ) + (r)f(r, ) = s(r) S(r, )f(r, )d + J(r, ), (15) с граничными условиями f-(z, ) = h(z, ), (z, ) -, (16) и условиями сопряжения f-(z, ) = (Bf+)(z, ), (z, ) . (17) В прямой задаче функции , s, J, S, h и оператор B заданы, а f требуется найти. Фазовая функция S(r, ) для всех r удовлетворяет условию нор мировки (8). Оператор сопряжения предполагается френелевским B = Bfr и имеет вид:
(Bfrf+)(z, ) = R(z, )f+(z, R) + T (z, )f+(z, T ). (18) Здесь R(или T ) направление потока излучения падающего на поверхность и в результате зеркального отражения (или соответственно преломления по закону Снеллиуса, меняющего свое направление на :
R = - 2n, T = (z, )n + (z, )( - n), = n(z), (19) i/j, если z Gi Gj, 0 < (z) 1, (z, ) = j/i, если z Gi Gj, -1 (z) < 0, sgn() 1 - 2(z, )(1 - 2), если 1 - 2(z, )(1 - 2) 0, (z, ) = 0, иначе, (20) где константы i и j обозначают показатели преломления сред Gi и Gj, а n(z) единичный вектор нормали к Gj, i < j.
Коэффициенты R и T в (18) характеризуют отражательную и пропускательную способность границы раздела сред Gi и Gj при френелевском отражении для неполяризованного излучения:
1 1 (z, ) 2 2 2 R(z, ) = (R (z, ) + R(z, )), T (z, ) = (T (z, ) + T(z, )), 2 2 (z, ) (21) (z, )(z, ) - () - (z, ) R (z, ) =, R(z, ) =, (22) (z, )(z, ) + (z, ) + (z, ) 2(z, ) 2(z, ) T (z, ) =, T(z, ) =. (23) (z, )(z, ) + (z, ) + (z, ) В зз 14,15 вводятся основные ограничения и функциональные пространства, необходимые для определения решения краевой задачи. Устанавлива ются некоторые свойства оператора сопряжения Bfr.
Функция f(r, ) принадлежит к классу D(G0 ) если почти при всех (r, ) G0 функция f(r+t, ) абсолютно непрерывна по t на интервале (t1(r, -), t1(r, )], f и rf принадлежат L(G ).
Величина t1(r, ) в определении класса D, в отличие от d(r, ), обозначает ближайшее расстояние от точки r в направлении до G0 (рис. 1).
Показано, что для любой функции f(r, ) D(G0 ) предельные значения f(z, ) принадлежат пространству L(). Под решением прямой задачи (15)-(17) понимается функция f(r, ) из пространства D(G0 ) с нормой f D = max f- L, rf + f, L удовлетворяющая почти всюду на G уравнению (15), почти всюду на и условиям (16) и (17), соответственно.
В з 16 доказываются некоторые вспомогательные утверждения с доста точно произвольным оператором B. В параграфе з 17 приводится пример условий сопряжения, для которых имеет место неединственность решения краевой задачи, и там же доказаны теоремы, касающиеся единственности граничной задачи с условиями сопряжения на основе формул Френеля.
Сформулируем основной результат главы.
Теорема 3. Пусть оператор сопряжения B = Bfr, , s L(G), J L(G ), S L(G ), h L(-) и = s/ < 1. Тогда решение краевой задачи f существует, единственно и почти всюду на G удовлетворяет неравенству J |f(r, )| max h L (-),, 1 - L(G) а при J 0 неравенству |f(r, )| h L.
В з 18 рассматриваются приложения прямой задачи к проблеме компьютерной визуализации трехмерных объектов. Предложен численный алгоритм нахождения решения, основанный на применении метода Монте-Карло с ветвлением траекторий. Проведены компьютерные эксперименты, показывающие эффективность предложенного алгоритма.
В главе IV изучается обратная задача для уравнения переноса излучения с обобщенными условиями сопряжения, заключающаяся в определении по данным радиационного просвечивания.
В з 19 формулируется постановка обратной задачи и приводятся основные ограничения, используемые в главе IV. В этой главе представлен случай, когда в условиях на границе раздела используется оператор сопряжения, являющийся линейной комбинацией френелевского и диффузного операторов (Bf+)(z, ) = fr(z, )(Bfrf+)(z, ) + d(z, )(Bdf+)(z, ). (24) В соотношении (24) оператор Bfr задается формулой (18), a Bd, описывающий диффузное отражение, имеет вид (Bdf+)(z, ) = Sd(z, , )f+(z, )d, (25) где функция Sd(z, , ), называемая индикатрисой отражения, ограничена и непрерывна на . Неотрицательные функции fr, d принадлежат Cb( ) и fr + d < 1. Они определяют френелевскую и диффузную составляющую в отраженном свете.
В з 20 вводятся функциональные пространства, и исследуются их свойства. В з 21 изучается непрерывность решения прямой задачи при условии, что h Cb(-) и коэффициенты уравнения непрерывные и ограниченные функции на G0 . Показано, что френелевская составляющая в операторе сопряжения приводит к достаточно сложной структуре областей непрерывности решения прямой задачи.
Основываясь на результатах параграфов з 20-21, в з 22 исследуется единственность решения обратной задачи. Рассматривается случай когда разбиение области G, состоит из двух областей G1 и G2 (G0 = G1 G2), причем G1 строго выпуклая область с гладкой границей и G1 G. Под выражени ем [f(r + d(r, ), )] понимается величина разрыва функции f(r, ) = f(r + d(r, ), ) в точке (r, ) по переменной r в направлении . В этих обозначениях и предположениях доказана Теорема 4. Пусть вектора и ортогональны друг к другу и лежат в горизонтальной плоскости. Если функция f является решением прямой задачи (15)-(17), тогда [f(r + d(r, ), )] = 0 для любой точки r G0 P, P = {r3 = const}, и любого горизонтального направления такого, что прямая {r + t : - < t < } не является касательной к кривой P, если же она касается ее в точке z = z(r, ) P, то |[f(r + d(r, ), )] | + |[f(r - d(r, -), -)] | = M(z), (26) где M(z) = lim (|M1(z, z )| + |M1(z, z )|), (27) z,z z M1(z, z ) = ((Bf+)(z, ) - f+(z, )) d(z,) z - z exp (z + t)dt, z, z , =, |z - z | При выполнения условия видимости M(z) > 0, показана единственность решения обратной задачи. В этом же параграфе приводятся ограничения на оператор сопряжения, при которых условие видимости выполняется. Для положительности функции M(z) достаточно, чтобы 0 < fr(z, ) < 1 для всех (z, ) , d(z, ) 0 при n(z) 0 и функция h была положительной. В частности, в случае чисто френелевского оператора при fr(z, ) < условие M(z) > 0 выполняется.
В конце параграфа приводятся численные эксперименты, проверяющие возможности алгоритма, указанного в доказательстве теоремы о единственности решения. Он заключается в построении некоторого оператора, переводящего плотность выходящего излучения в характеристическую функцию области G1. Для расчета функции H, необходимой для тестирования алгоритма решения обратной задачи, использовался численный алгоритм разработанный в з 18. В случае, когда в среде находится более одного включения предложена модификация метода, которая позволяет обнаруживать несколько выпуклых включений.
В главе V рассматриваются экстремальные задачи для уравнения переноса с френелевскими условиями сопряжения. Кроме последнего параграфа исследования проводятся для случая плоскопараллельной симметрии среды.
В з 23 исследуется прямая задача для уравнения переноса в плоско-параллельном неоднородном слое c оператором сопряжения общего вида. При ограничениях общего характера показана разрешимость краевой задачи на некотором подмножестве пространства кусочно-непрерывных функций. Для единственности решения прямой задачи в плоско-параллельном слое достаточно, чтобы норма оператора сопряжения не превосходила единицы. В трехмерной ограниченной среды единственность решения при этом же условии доказана только для френелевского оператора сопряжения.
В з 24 рассматриваются постановки экстремальных задач для уравнения переноса излучения и обсуждаются их физические аспекты.
Исследуется случай, когда неоднородная среда, имеющая плоскую геометрию, заполняет область G = {z : z (z0, z3)} и состоит всего лишь из трех веществ Gi = {z : z (zi-1, zi)}, i = 1, 2, 3, отделенных друг от друга плоскостями z = z1 и z = z2, и для каждого вещества Gi задан коэффициент преломления i (рис. 4a).
(a) (b) Рис. 4.
В этом случае уравнение переноса и граничные условия имеют вид:
fz(z, ) + (z)f(z, ) = s(z) S(z, , )f(z, )d + J(z, ). (28) -f-(z0, ) = h(), > 0, f-(z3, ) = h(), < 0, (29) f-(zi, ) = (Bfrf+)(zi, ), i = 1, 2. (30) Величина , [-1, 1] в (28)-(30) есть косинус угла между направлением распространения излучения и положительным направлением оси симметрии среды.
Рассматривается задача, заключающаяся в определении коэффициент преломления 2 из (28),(29),(30) и условия f+(z3, 1) max, (31) если коэффициенты 1,3 и функции , s, S, J, h известны.
Физическое содержание поставленной задачи заключается в выборе показателя преломления промежуточного слоя из условия наилучшего просветления среды. В последнее время интерес к этой задаче возрос не только в связи с изготовлением просветляющих покрытий, но и благодаря развитию оптических исследований биотканей. Добавление промежуточного слоя с оптимальным показателем преломления, который покрывает основную ткань пациента, позволяет добиваться достижения наибольшего терапевтического эффекта в лазерной медицине.
В з 25,26 проводится исследование задачи просветления оптики для некоторых упрощенных вариантов уравнения переноса. В трехслойной слаборассеивающей или рассеивающей преимущественно вперед среде получены аналитические решения, которые совпадают с классическим решением 2 = 13. Наличие рассеяния в среде может приводить как к увеличению, так и к уменьшению оптимального показателя преломления промежуточного слоя.
В з 27 рассматриваются экстремальные задачи для уравнения переноса в трехмерном случае.
С физической точки зрения суть задачи заключается в выборе оптимального коэффициента преломления среды G2, покрывающей включение G3 так, чтобы минимизировать влияние составного тела G2 G3 на отраженное и рассеянное от него излучение в некоторой точке (r, ). В этом случае, покрытие G2 является маскирующим для объекта G3 при визуализации его вблизи точки r в направлении (рис. 4b).
Сформулируем математическую постановку экстремальной задачи. Пред полагается, что оператор сопряжения B = Bfr, и внутренние источники отсутствуют. Обозначим через f решение прямой задачи (15)-(17) с заданными функциями , s, S, h и показателем преломления 1 в случае, когда среда G не содержит включений, то есть G0 = G. Поместим внутри области G два включения G2, G3. Обозначим через G1 = G \ G2 G3 и предположим, что G2 полностью "покрывает" область G3 и является разделяющим множеством между областями G1 и G3. Через f обозначим решение прямой задачи (15)(17) для разбиения G1, G2, G3 с известными коэффициентами (r), s(r), (r) и функциями S, h, причем (r) = (r), s(r) = s(r), (r) = 1 для всех r G1, а функция (r) равна 2 и 3 в областях G2 и G3, соответственно.
Экстремальная задача заключается в нахождении коэффициента преломления 2 среды G2 так, чтобы в фиксированной точки (r, ) он сообщал минимум функции F (2) = (f(r, ) - f(r, , 2))2, (32) если для областей G1, G3 заданы коэффициенты 1, 3 и функции , s,, s,S,h известны.
Для решения поставленной задачи предложен численный метод, заключа ющийся в последовательном вычислении функции f(r, , 2) на некотором дискретном множестве из заданного интервала [2,min, 2,max] с последующим нахождением соответствующего значения , доставляющего минимум целевой функции (32). Вычислительные эксперименты показали, что при компьютерной визуализации объектов покрытых маскирующим слоем действительно происходит ухудшение их видимости.
В заключении приводятся основные результаты диссертации:
1. Доказана единственность решения обратной задачи для полихроматического уравнения переноса, заключающаяся в определении неизвестного коэффициента ослабления по заданному излучению на внешней границе просвечиваемой среды. Для ее решения предложен специальный тип внешнего источника излучения, имеющий разрыв первого рода по угловой переменной. Разработан численный алгоритм решения обратной задачи и проведены вычислительные эксперименты.
2. Для задачи определения границ разрывов коэффициентов монохроматического уравнения переноса показано, что равенство нулю меры видимости на границе разрывов приводит к существованию бесконечного множества решений. Построены примеры неединственности. Теоретически и численно показано, что в оптически плотных средах непрерывность коэффициента поглощения на границе раздела является хорошей аппроксимацией условия невидимости среды.
3. Рассмотрена задача поиска оптимальной энергии излучения внешних источников, с целью получения наилучшего качества реконструкции среды по данным томографического просвечивания. Ее решение базируется на обеспечении условий, увеличивающих меру видимости. Предложена упрощенная постановка задачи оптимизации и обосновано, что в дефектоскопии использование ее решения оправдано.
4. Доказана разрешимость прямой задачи для уравнения переноса с обобщенными условиями сопряжения на границе раздела сред в ограниченной трехмерной области и для пространственно одномерного случая плоскопараллельной симметрии. Исследованы непрерывные свойства решения задачи.
Получены оценки типа принципа максимума.
5. Показана единственность решения задачи определения границы раздела сред при условии, что оператор сопряжения является линейной комбинацией френелевского и диффузного операторов, и среда содержит одно выпуклое включение. Предложен численный метод нахождения решения.
6. Сформулированы новые экстремальные задачи для уравнения переноса, применительно к оптике просветляющих и маскирующих покрытий. Получены частные аналитические решения и построены численные алгоритмы нахождения решения в общем случае. Показано, что наличие рассеяния в среде существенно сказывается на решение экстремальной задачи.
Основные публикации автора по теме диссертации 1. Аниконов Д.С., Ковтанюк А.Е., Прохоров И.В., Использование уравнения переноса в томографии. М.: Логос, 2000, 224с.
2. Anikonov D.S., Kovtanyuk A.E., and Prokhorov I.V., Transport Equation and Tomography. Utrecht-Boston: VSP, 2002, 216p.
3. Anikonov D.S., Nazarov V.G., and Prokhorov I.V., Poorly Visible Media in X-ray Tomography. Utrecht-Boston: VSP, 2002. 302p.
4. Аниконов Д.С., Прохоров И.В., Определения коэффициента уравнения переноса при энергетических и угловых особенностях внешнего излучения. //Доклады АН. 1992. Т. 327. №2. С. 205-207.
5. Anikonov D.S., Prokhorov I.V., and Kovtanyuk A.E., Investigation of Scattering and Absorbing Media by the Method of X-ray Tomography. //Journal of Inverse and Ill-Posed Problems. 1993. V. 1. №4. P. 259-281.
6. Anikonov D.S., Kovtanyuk A.E., and Prokhorov I.V., Some Numerical Experiments in Tomography in of Scattering Medium. //In book: Computerized Tomography, by editor M.M. LavrentТev. Inverse and Ill-Posed Problems Series. Boston-Utrecht. VSP. 1995. P. 25-31.
7. Прохоров И.В., Некоторые свойства решений уравнения переноса. //ДВ матем. сборник. 1996. T. 2. С. 161-172.
8. Аниконов Д.С., Назаров В.Г., Прохоров И.В., Видимые и невидимые среды в томографии. //Доклады АН. 1997. Т. 357. №5, С. 599-603.
9. Anikonov D.S., Kovtanyuk A.E., and Prokhorov I.V., Tomography through the transport equation. //Proc. IMA Volumes in Mathematics and its applications. In book: Computational Radiology and Imaging: Therapy and Diagnostics. Springer-Verlag. New York. 1999. Vol. 110. P. 33-44.
10. Аниконов Д.С., Прохоров И.В., Значение коэффициента поглощения в диагностике рассеивающих и поглощающих сред. //Доклады АН. 1999.
T. 368. №1. С. 24-26.
11. Аниконов Д.С., Прохоров И.В., Некоторые математические модели томографии для особых состояний сред. //Доклады АН. 2000. T. 371. №4.
С. 452-456.
12. Прохоров И.В., Краевая задача теории переноса излучения в неоднородной среде с условиями отражения на границе. //Дифференциальные уравнения. 2000. Т. 36. №6. С. 848-851.
13. Аниконов Д.С., Прохоров И.В., Необходимые и достаточные условия единственности одной задачи томографии. // Журнал вычислительной математики и математической физики. 2002. Т. 42. №3. С. 370-379.
14. Прохоров И.В., Определение поверхности раздела сред по данным томографического просвечивания. // Журнал вычислительной математики и математической физики. 2002. Т. 42. №10. С. 1542-1555.
15. Anikonov D.S., Konovalova D.S., Kovtanyuk A.E., Nazarov V.G., and Prokhorov I.V., Investigation of Inverse and Other Non-classical Problems in the Russian Far East. //In book: "Resent Developments in Theories and Numerics". Word Scientific Publisher. 2003. P. 13-26.
16. Прохоров И.В., О разрешимости краевой задачи для уравнения переноса излучения с обобщенными условиями сопряжения на границе раздела сред. //Известия РАН. Серия математическая. 2003. Т. 67. №6. C. 169192.
17. Прохоров И.В., Яровенко И.П., Краевая задача теории переноса в многослойной среде с обобщенными условиями сопряжения. //Сибирский журнал индустриальной математики, 2003. T. 6. №1. С. 93-107.
18. Аниконов Д.С. Прохоров И.В., Формальная оценка качества реконструкции в рентгеновской томографии. // Доклады АН. 2005. Т. 401. №3.
С. 312-315.
19. Прохоров И.В., Яровенко И.П., Численное решение дифракционных задач для уравнения переноса излучения. // Сибирские электронные математические известия. 2005. Т. 2. С. 88-101.
20. Prokhorov I.V., Yarovenko I.P., and Krasnikova T.V., An extremum problem for the radiation transfer equation. // Journal of Inverse and Ill-Posed Problems. 2005. V. 13. №4. P. 365-382.
21. Аниконов Д.С., Прохоров И.В., Назаров В.Г., Солнышко Н.В., Способ маскировки изделий. //Патент Российской Федерации №2264424. Бюллетень №32, 20.11.2005.
22. Аниконов Д.С., Прохоров И.В., Постановка и численное решение задачи оптимизации в рентгеновской томографии. // Журнал вычислительной математики и математической физики. 2006. Т. 46. №1. С. 18-25.
23. Аниконов Д.С., Прохоров И.В., Выбор оптимальной энергии излучения в рентгеновской дефектоскопии. //Доклады АН. 2006. Т. 408. №4. С. 455459.
24. Прохоров И.В., Яровенко И.П., Исследование задач оптической томографии методами теории переноса излучения. //Оптика и спектроскопия.
2006. T. 101. №5. С. 817-824.
25. Kovtanyuk A.E., Prokhorov I.V., Tomography problem for the polarizedradiation transfer equation. // Journal of Inverse and Ill-Posed Problems.
2006. V. 14. №6. P. 609-620.
Прохоров Игорь Васильевич Математические задачи теории переноса излучения Автореферат Лицензия ЛР N 040118 от 15.10.96 г.
Подписано к печати 9.07.2007 Формат 60х84/Усл.п.л. 1.8 Уч.-изд.л. 1.Тираж 100 экз. Заказ 2Отпечатано в типографии издательства "Дальнаука" 690041, Владивосток, Радио, Авторефераты по всем темам >> Авторефераты по техническим специальностям