Содержание к диссертации
Введение
1 Задачи доплеровской томографии для моделей вихревого поля . 10
1.1 Задача определения параметров моделей вихревого поля 11
1.1.1 Постановка задачи определения параметров моделей вихревого поля. 11
1.1.2 Метод определения параметров пятипараметрической модели вихревого поля 17
1.1.3 Метод определения параметров шестипараметрической модели вихревого поля 23
1.2 Задача определения параметров комбинированных моделей вихревого поля. 30
1.2.1 Постановка задачи определения параметров комбинированных моделей вихревого поля 30
1.2.2 Метод определения параметров одиннадцатипараметрической комбинированной модели вихревого поля 31
1.2.3 Метод определения параметров двенадцатипараметрической комбинированной модели вихревого поля 39
1.3 Алгоритмы численного решения задач доплеровской томографии для моделей вихревого поля и реализующий их программный комплекс. Результаты вычислительного экперимента 45
2 Задача доплеровской томографии для модели кусочно-постоянного поля . 62
2.1 Задачи доплеровской томографии для модели кусочно-постоянного поля с различными вариантами задания информации 63
2.1.1 Постановка задач для модели кусочно-постоянного поля 63
2.1.2 Методы решения поставленных задач 67
2.1.3 Анализ существования решения рассматриваемых задач 72
2.2 Задача доплеровской томографии в случае нескольких измерений вдоль каждой прямой 74
2.3 Численное решение задач доплеровской томографии для модели кусочно-постоянного поля. Программная реализация. Результаты вычислительного эксперимента 82
Заключение. 94
Список литературы.
- Метод определения параметров пятипараметрической модели вихревого поля
- Постановка задачи определения параметров комбинированных моделей вихревого поля
- Постановка задач для модели кусочно-постоянного поля
- Задача доплеровской томографии в случае нескольких измерений вдоль каждой прямой
Введение к работе
Вычислительная томография (или просто томография) — область математики, занимающаяся разработкой математических методов и алгоритмов восстановления внутренней структуры объектов по проекционным данным. Математические основы томографии были заложены в 1917 году Радоном, который в своей работе [31] рассматривал оператор интегрирования по гиперплоскостям в евклидовом пространстве, названный позже преобразованием Радона, и предложил явный метод решения задачи восстановления функции по ее проекционным данным.
Однако интенсивное исследование задач томографии началось с середины 20-го века в связи с появлением компьютерных томографов. Внедрение методов компьютерной томографии в медицину позволило существенно повысить эффективность диагностики и обеспечило создание новых методов лечения. Со временем, методы компьютерной томографии стали также широко использоваться в электронной и рентгеновской микроскопии — для получения структур кристаллов и макромолекул, в геофизике — для поиска месторождений полезных ископаемых, на производственных предприятиях — для анализа качества продукции и поиска дефектов и в других областях науки и техники.
Радоном была решена задача восстановления функции f(x,y) по функции TZf(p, р), заданной для всех возможных р, (р.
Реальные математические постановки и методы решения задач компьютерной томографией обусловлены конкретными схемами их практического применения. Набор исходных данных для решения задачи восстановления зависит от применяемой схемы скани 4
рования.
На практике (рассматриваем задачу на плоскости, задача восстановления трехмерной функции обычно решается как задача восстановления двумерной на множестве сечений) наиболее распространены параллельная [8], [9], [10] и веерная [8] схемы. При этом существует два типа задач реконструкции: задача с полными данными и задача с неполными данными [2], [8]. В задаче с полными данными предполагается, что интегральные проекции известны для всех прямых, проходящих через исследуемый объект, т.е. для всех необходимых р, if, а задачи с неполными данными могут быть следующих типов:
• задача с ограниченным диапазоном углов,
• внутренняя задача, когда 7lf(p, (р) известна для р а, а 0,
• внешняя задача, когда TZf(p, ф) известна для р а 0,
• задача с ограниченным числом источников и приемников излучения.
Существуют различные численные методы решения задач компьютерной томографии, такие как метод преобразования Фурье, метод свертки и обратной проекции, метод итераций и методы регуляризации, которые достаточно подробно исследованы (напри-мер, [8], [9], [10], [13], [16]).
Одним из важных классов задач томографии являются задачи доплеровской томографии ([1, [21], [19]), отличительной особенностью которых является восстановление векторных функций по проекционным данным. Такие задачи возникают, например, в радиометеорологии, в связи с необходимостью определения поля скорости ветра в атмосфере [1], в магнитно-резонансной томографии (29], при измерении потоков плазмы [15], [18, [17], в астрофизике [23] и других областях современной науки.
Многие задачи доплеровской томографии являются нелинейными [3], [11]. Поэтому их решение требует поиска особых подходов в каждом конкретном случае и сопряжено со значительными трудностями.
В то же время большое количество работ по доплеровской томографии посвящено линейным постановкам задачи, которые можно интерпретировать как линеаризацию задачи в нелинейной постановке путем отбрасывания части исходных данных, оставляя только линейную часть измеренного спектра [3]. Известно [35], что векторное поле можно представить в виде суммы двух составляющих: соленоидалыюго и потенциального поля. В ряде работ (например [21], [12], [38], [24]) показано, что по линейной части спектра можно однозначно восстановить только соленоидальную часть поля, в то время как потенциальнаю часть не может быть найдена [20]. Чтобы полностью восстановить поле, используется различная дополнительная информация. В работе [12] авторы используют данные дополнительных измерений, дающие информацию о компоненте поля, поперечном линии интегрирования. В работах [21], [24] и [28] дивергенция поля полагается равной нулю внутри исследуемой области, и налагаются граничные условия Неймана.
Процесс восстановления соленоидального и потенциального полей в трехмерном пространстве по измерениям интегралов скалярного произведения поля на направляющие вектора измерительных плоскостей описан в работах [29], [30]. Аналогичный метод для ограниченных полей рассматривается и [25], [26], [27[ и [40].
Ряд работ посвящен численному решению задач векторной томографии. В работах [14], [36], [45], [43] разработаны алгоритмы итеративного восстановления соленоидалыюго поля, основанные на методе алгебраического восстановления (ART). В работах [34], [35] изучены способы детектирования потенциальной составляющей поля, возникающей при решении задачи восстановления соленоидального поля из-за ошибок измерений, а также способы увеличения точности восстановления.
Установлено, что эта задача в общем случае, несмотря на существенную переопределенность, имеет неединственное решение ([3], [4], [11]), и что по заданной функции М(р, (р,и ) могут быть определены функции [3]: Fn(p, ip) — f ° Wn(p, ip, s)ds, n = 1,2,..., p со, (p Є [0,2тг], называемые доплеровскими моментами порядка п. Также в работе [3] показано, что эту задачу можно свести к системе уравнений в частных производных
где /„(ж, у) однозначно определяются функциями Fn(p,ip), а значит и М(р,(р,ш).
В связи с неединственностью решения задачи доплеровской томографии в общей постановке возникает вопрос об использовании некоторой априорной информации о виде неизвестного поля, использование которой позволило бы достичь единственности решения. Одним из возможных подходов здесь является переход к конечномерным моделям поля.
Диссертационная работа посвящена исследованию задач доплеровской томографии для некоторых конечномерных моделей векторного поля, разработке методов их решения и программной реализации предложенных алгоритмов.
В первой главе поставлены и изучены задачи доплеровской томографии для конечно-параметрических моделей вихревого поля. Предложены и программно реализованы численные алгоритмы их решения. Под конечно-параметрическим вихревым полем подразумевается поле вида с неизвестными параметрами А, П, Гг, R, х0, уо Для этих моделей найдено число прямых, достаточное для однозначного восстановления вихревого поля, предложен метод определения неизвестных параметров.
Во втором параграфе исследованы задачи доплеровской томографии для комбинированных параметрических моделей вихревого ноля, представляющих собой суперпозицию вихревого поля с функцией амплитуды заданной (1) или (2) и поля вида {е\Х + в2У + е-л, fix + fiV + /з}5 гДе еъ е2, Єз, /і, /я, /з — неизвестные параметры. Для полей такого вида также определено число прямых L(pi, (pi), достаточное для того, чтобы задача имела единственное решение, представлены методы решения.
В третьем параграфе предложены вычислительные алгоритмы решения задач доплеровской томографии для конечно-параметрических моделей векторного поля. Эти алгоритмы реализованы в разработанном программном комплексе, написанном на языке программирования C++ и имеющем модульную архитектуру, что позволяет включать в него новые модули решения задачи для иных параметрических моделей векторного поля. С использованием разработанного программного обеспечения проведены вычислительные эксперименты, в ходе которых изучены свойства и границы применимости предложенных вычислительных алгоритмов.
Во второй главе исследуются задачи доплеровской томографии для модели кусочно-постоянного поля. Кусочно-постоянное поле определяется следующим образом. На области D = [0, а] х [0, а], за пределами которой V(x,y) = 0, введена прямоугольная сетка
Для задачи в такой постановке проведен анализ существования и единственности решения, предложен метод его определения.
Во втором параграфе исследуется задача доплеровской томографии для модели кусочно-постоянного поля, в которой известны результаты нескольких измерений вдоль каждой прямой. А именно, для каждой прямой Ь{р,ф) считаем известной М (р, р,си) — меру Лебега множества точек прямой L(p, ср) лежащих внутри области D, для которых выполнено W(p, ср, s) -Jr , где А = A«i,..., /.
Для данной задачи доказаны теорема о единственности решения и теорема об условиях существования решения.
В третьем параграфе разработаны вычислительные алгоритмы решения задач до-плеровской томографии для модели кусочно-постоянного поля, которые реализованы в компьютерной программе, написанной на языке для математических вычислений Octave. Программа использована для проведения расчетов, в ходе которых изучены свойства предложенных методов. В конце параграфа представлены примеры, иллюстрирующие работу этих методов.
Основные рещультаты, полученные в диссертационной работе, были представлены на нескольких научных конференциях и опубликованы в ряде изданий. Ссылки на публикации автора по теме диссертации представлены в общем библиографическом списке в конце работы.
Метод определения параметров пятипараметрической модели вихревого поля
Определить, какой из трех вышеописанных случаев имеет место для каждой конкретной прямой семейства L(p, /3) можно, анализируя значения функции М(р, (3,ш). В случае 1, М(р, /3,о ) принимает всего два различных значения: Мр и 0. В случае 2, согласно (1.1.21) и (1.1.22), М(р,/3,ш) непрерывна при всех значениях и. А в случае 3, согласно (1.1.26) и (1.1.27), функция М(р,(3,ш) будет иметь разрыв в точке cu = — Arfc m".
Условия теоремы на Д/? гарантируют, что не менее трех прямых из семейства {L}c пересекают множество (х — XQ)2 + (у — уо)2 R2 ПРИ этом как минимум две из них пересекают облать (а; — х0)2 4- (у — уя)2 г2, и как минимум одна нет.
Получив значения функций М(р3,[33,ш), для всех j33 — j-A/З и соответствующих им pj, находим знак выражения — Asgna и sgna для каждой прямой, используя информацию о взаимном расположении прямых.
Найдем центр вихревого поля. Для этого выберем две прямые L(pi,/5i), L(p2)(52) из семейства {L]ci и две прямые Ь(р3,Рз), (Р4 АО из семейства {Ь}с2 так, чтобы все они пересекали область (х — XQ)2 + (у — т/о)2 г2. (Согласно условиям на Д/9 такие прямые существуют.)
Поясним, что здесь и далее обозначение L(pi,/?i) (и аналогичные ему) использовано вместо обозначения L(pJ1,PJl), чтобы не перегружать текст индексами. И его следует понимать, не как "прямая со значением / = .7 Д/3 при j = 1", а как "прямая с некоторым значением /?л".
Обозначим di = y/(xQ - xCl)2 + (z/o - Усі)2 и d2 = у/(х0 - хС2)2 + {уо - Ус2)2 Рассмотрим сначала прямые L(pi,0i), L(p2,/5 2)- Для них, используя, (1.1.26) и (1.1.27), находим такие Ш{, і = 1,2, для которых выполнено: f Arginfw(M(pi,/3i,w) = MKiA) при - Л sin «і О, и і = s (1.1.28) l Argsupw(M(pi,$,a/) = 0) при - Ashing 0. Эти значения ш соответствуют экстремальным значениям функции W(pi, Pi,s), и в соответствии с (1.1.25), ШІ = — disinai. Получаем систему уравнений -jdiRinoii =u)u -fdlsma2 = w2, а\ — а2 = А, от которой переходим к системе Ґ ш2 sin «і = uj\ sin а2, I ai - a2 — A. Здесь A — известный угол между прямыми L(pi,/3i), L(p2,P2), равный /Зі — /. Из системы получаем уравнение ш2 sin(a:2 + А) — wisina . На интервале (—7г/2;7г/2) это уравнение имеет единственное решение, которое можно найти численно. Докажем единственность решения. Пусть это не так и есть два решения: а и а; А Є (0;тг), «,ає(-7г/2;7г/2) sin(a + A) = sina, 8ІП(Й + Д) = 8ІПЙ, sin(a: + A) sine sin(a + A) sin a sin(a + A) sin a = sin(a + A) sin a. Следовательно, sin a sin a cos A + cos a sin a sin A = sin a sin a cos A + sin a cos a sin A, sin A(sin(a — a)) — 0. И т.к. А Є (0; тг), то sin(a — a) = 0 и a — a = 7rn, n Є Z. Получили противоречие. Таким образом, «і и cv2 определены однозначно. Для прямых Ь(рз,Рз), L(p4, /З4) аналогичным образом получаем систему уравнений - d2sina3 =шз, — c/2sintt4 = W4 аз — ОІА = А ) где A = / — / и для г = 3,4, а; задается (1.1.28). Решая эту систему, находим а3 и «4.
Для прямых L(pi,(3i) и Ь(рз,Рз) верны следующие соотношения: жсі + rfi cos (/ - ai) = жс2 + d2 cos(/33 - a3) yCl + rfi sin( ! -ai) = yc2+ d2 sin(/33 - a3) где ai,/?i, а?з,/?з уже известны. Это линейная система уравнений относительно d\ и d2, с определителем, отличным от нуля (Заметим, что если для пары прямых L(pi,/3i) и Ь(рз,Рз) эта система имеет определитель, равный нулю, то тогда прямую 1/(р35Дз) заменяем прямой L(p4, А) и система становится разрешимой единственным образом.). Из нее находим d\ и d2. Зная di и «і, в соответствии с (1.1.7), находим центр вихревого поля. Перейдем к определению остальных параметров модели. Для рассмотренной выше прямой L(pi,/3i), пересекающей область (х — XQ)2 + (у — Уо)2 г2, согласно (1.1.25), (1.1.26), (1.1.27), можно найти wo такое, что со0 = — Ad smai _ Получаем уравнение А = Кгг, (1.1.29) где К\ — (wor)/(c/isina;i) — вычисленная константа. Рассмотрим прямую L(p5, /35), пересекающую область (х — XQ)2 + (у — у0)2 R2, но не пересекающую (х — Хо)2 + (у — уо)2 г2, т.е. такую, для которой имеет место описанный выше случай 2. Для нее, согласно (1.1.21), (1.1.22), вычисляем Arg infw(M(р5, Аь и) = MP6tPi) при - A sin а5 О, U)5 = , Argsupw(M(p5,As,w) = 0) при - Asina5 0, и из (1.1.18) получим уравнение А -(dismay — -Rsgna5) = OJ5. (1.1.30) R-r Возьмем теперь два значения to: wfi, uij такие, что для г = 6, 7, - (di sin а Rsgn a5) Wj 0, если — Asma5 0, и 0 0 37( 1 sina5 — .Rsgiia ), если —.4sina5 0. Тогда из (1.1.21), (1.1.22) получим два уравнения (г = 6,7) МІ о /« 2 f Л2Д2 Т 2ч / af sin af5 -r-r-;—: т тттг — 1 , у l J\(Adlsinar, uJi{R-r))2 / 23 где МІ — вычисленные константы. Преобразуя эти уравнения, получаем ( AR = K3{Adi sin а5 - (R - г)ш6), \ AR — K4(Adislna5 — (R — г)ш7), и, следовательно, {R - г){К3ш6 - К4ш7) = Adx sina;(i-G - К4) (1.1.31) 4df sinz Q!5 у 4df smz a5 Преобразуя (1.1.31), а также используя (1.1.29) и (1.1.30), получаем систему уравнений относительно A, R, г: -jj {.dx sin a5 - Rsgn a5) = w6 A = Kir, r _ і Лгії sin 05(/ 3-/ 4)
Эта система имеет единственное решение, которое находится в явном виде. Все параметры модели найдены. Из описанного метода следует справедливость утверждения теоремы.
Метод определения параметров шестипараметрической модели вихревого поля.
Теперь рассмотрим шестипараметрическую модель вихревого поля, в которой функция амплитуды а(г) задается формулой (1.1.5). Параметрами данной модели являются координаты положения центра вихревого поля (х0,уо), а также параметры A, R, г і, r i его функции амплитуды.
Пусть известно, что внешний радиус вихревого поля R не меньше некоторого значения Rmin, а также что rimin гх гХтах и г2тш г2 г2тах- Тогда верна следующая теорема:
Теорема 1.2. Если Д/3 2 min{arcsin Кті пц с arcsin r,m J max, arcsin г Ш}, то все параметры шестипараметрической модели вихревого поля могут быть однозначно определяются результатами измерений вдоль прямых двух семейств {L}Cl и {L}c2, где С\ и Сі — две различные произвольные точки границы области D.
Доказательство. Рассмотрим семейство прямых, проходящих через произвольную точку С границы области D. Для каждой прямой семейства {L}c возможны четыре взаимоисключающих случая пересечения с областью D.
Постановка задачи определения параметров комбинированных моделей вихревого поля
Постановка задачи определения параметров комбинированных моделей вихревого поля. Рассмотрим модель векторного поля специального вида. Пусть V(x,y) = R{x,y) + E(x,y), (1.2.1) где R(x, у) = \ а(х, у) У ===== , а(х, у) Х [ \ V \х - хо) +{у- УоУ J \ V(x - хо) + (У - УоУ и Ё(х, у) = {ехх + е2у + е3, fxx + f2y + /з}
Пусть V(x, у) тождественно равна нулю вне области D = [0,Ь] х [0,Ь]. Функция а(х,у) = а(г) — радиально-симметрична относительно точки (хо,уо) и a(r) = 0 при f R, где f = х/(:г - х0)2 + (у - т/о)2, а ег, е2, е3, /ь /2, /з — числовые константы, причем еь е2, /і, /2 не обращаются в нуль одновременно. На поле наложены ограничения (1.1.3). Такая модель представляет собой суперпозицию вихревого поля, рассмотренного в 1.1, и линейного поля с неизвестными параметрами ел, е2, 3,/1,/2,/3 Будем исследовать две комбинированные модели вихревого поля, для которых функция а(г), задается формулами (1.1.4), (1.1.5). Множество D пересекают семейства прямых {L}Ci (как и в 1.1, см. рис. 1.2). Каждая прямая г-го семейства проходит через точку СІ, принадлежащую границе D. Координаты точек j-той прямой г-го семейства задаются формулами X = Хс{ + S COS /3j, У = УСг +ssin/3j, \s\ 00, где /3j = jA/3 — угол между прямой и координатной осью Ох, а хсл и уа — координаты точки Cj.
Для каждой прямой Ь{р,ф) нам известна функция М(р,ір,ш) — мера Лебега множества, на котором W(p,cp,s) и , где W(p,ip,s) определяется (1.1.6). Заметим, что Р = (р + тг/2 и А/3 = Аїр.
Сформулируем задачу: Найти максимально возможное значение А/3, при котором можно определить все параметры модели и найти эти параметры.
Метод определения параметров одиннадцатипараметриче-ской комбинированной модели вихревого поля.
Рассмотрим комбинированную модель вихревого поля (1.2.1), где функция амплитуды а{г) поля R(x,y) задается формулой (1.1.4). Начнем решение задачи с определения параметров линейной составляющей поля Е{х,у). Для произвольной прямой Ь(р,(р,ш), проходящей через область D, представим функцию W(p, ср, s) в виде суммы двух функций W(p, if, s) = WR(p, p, s) + WE(p, tp, s), где WR(P, p,s) = R- k, WE(p, if, s) — E k, k= {—sirup, cos p}. Запишем функцию Ё(х,у) в параметрической системе координат (p,cp,s), связанной с прямой L(p,(p): Е(р, р, s) = {ei(pcos(p — s sin ip) + ea(psin cp + s cos ер) + e3, /i (p cos /5 — .s sin ( ) + /2 (p sin cp + s cos cp) 4- /3 }. Тогда функция WE{P, ip,s) = E к задается формулой WE(P, cp, ) = s(ei sin2 p — e2 sin ip cos p — fi sin cos + /2 cos2 /?)+ +p(fi cos2 cp + f2 sin y? cos /? — Єї sin (/? cos p — Є2 sin2 ) + /3 cos 93 — e3 sin уз. (1.2.2)
Выделим два случая прохождения прямой L(p, cp) через область D. 1) Прямая не пересекает множество (х — х0)2 + (у — уо)2 R2. Тогда W(p,cp,s) = WE(p, р, s), т.к. R(p, р, s) = 0. В этом случае результат измерений — функция М(р, ip, со) линейна, т.к. линейна функция W(p,cp,s): М{р,ср,со) = kisu + k2, где sw - такое значение s, при котором W(p, cp, su) — со (s линейно зависит от ш), а и &2 — некоторые константы. 2) Прямая пересекает множество {x — x0)2 + (y — yQ)2 R2. В этом случае R(p, cp, ,ч) ф 0, и, как было показано в 1.1, для рассматриваемых моделей R(p, ср,со) функция М(р, ср,со) нелинейна по со. Таким образом, анализируя значения М(р,ср,со), для любой прямой L(p,cp), можно определить, пересекает она область (х — х0)2 + (у — у0)2 R2 или нет. В [3] было доказано, что по известной функции М(р, ср, со) могут быть определены следующие функции: оо РП(Р, Р)= f Wn(p,cp,s)ds, п = 1,2,..., Н оо, ре[0,27г]. — оо Значит для каждой прямой L(p, ip) нам известна функция оо Fi(p,cp)= J W(P,cp,s)ds. (1.2.3) Для нахождения линейной составляющей поля достаточно выбрать шесть прямых, для которых имеет место случай 1, попарно непараллельных между собой. Для этих прямых, согласно (1.2.2), (1.2.3) имеем «і Fi(P, Ч ) = / (ais + a2)ds = -(sj- sg)oi + (вг - s0)a2, so где so, б і — значения переменной s, соответствующие точкам пересечения прямой Ь(р, ф) с границами области D, Ч\ = а\(р,ф) = еі sin2 у? — е2 sin (р cos ip — /і simp cos (p + /2 cos2 p, (1.2.4) «2 = аг(р, ) =p(/iCos2 v3+/2sin cos -e1 sin /?cos —e2 sin2 99)+/3 cosy?—e3 sin (p. (1.2.5) Обозначив Fu — 1(7 ,1 ,-), получаем систему из 6 линейных уравнений относительно еь е2,е3,/1,/2,/3: Ь = 2 (S 50г) 1і + (»1 - S0i)O2i, « = 1, —, 6, где an = аі(рг,(рі), а2г = a2(p;, (/). Решив эту систему, определяем все параметры функции Е(х,у).
Сформулируем условия, достаточные для нахождения линейной составляющей поля. Для этого необходимо, чтобы по меньшей мере шесть прямых не пересекали область (х — х0)2 + (у — у0)2 В2. При этом не имеет значения, будут ли все эти прямые принадлежать одному семейству или разным.
Пусть известно, что внешний радиус вихревой составляющей поля 7?, Rmax. Чтобы определить линейную составляющую поля, используя только одно семейство прямых {Ь}с, проходящих через точку С с координатами (хс, ус), достаточно потребовать, чтобы для точки выполнялось условие \/(хс — жо)2 + (ус З/о)2 dmin и угол между прямыми семейства удовлетворял соотношению А/3 ДДв, где f 1 ( - 2 arcsin f—) , если С Є {(0,0),(0,6), (6,0),(6,6)}, А/?Е = , „""" (1.2.6) і 7Г-2arcsin-f - , иначе. k 6 у «mm / При использовании измерений вдоль двух семейев прямых {L}ci и {L}c2i для которых выполнено
Постановка задач для модели кусочно-постоянного поля
Поскольку формулы (2.1.12) верны для всех / = 1,2,..., и А; = 1, 2,..., п, то полученные таким образом решения удовлетворяют требованиям теоремы. И других решений система иметь не может. Теорема доказана. Как видно из теоремы, задача 3 в общем случае может иметь 2 ? решений. Однако заметим, что суіцествуют такие классы wik, в которых задача имеет единственное решение. Таким классом являются, например, Wik, монотонные по /. (Для таких w в формулах (2.1.11) можно однозначно определить w2i-itk и W2i,k из условия монотонности.)
Перейдем к анализу задачи 4. Напомним, что задача 3 в общем случае может иметь 2- решений. И, чтобы уменьшить это множество неединственности, в задаче 4 используются данные измерений по дополнительным п2 прямым. Найдем ограничения на ги , при которых задача 4 имеет единственное решение. Теорема 2.2. Пусть решение wik задачи 4 существует и для всех I = 1,2,..,, к — 1,2,.., Ц выполняется хотя бы одно из условий 21-\,2к-1 Ф W2l-l,2k, (2.1.13) W2l,2k-1 Ф V 2l,2k, (2.1.14) 21-1,2к-1 = 21,2к-1, (2.1.15) W2l l,2k = U)2l,2k, (2.1.16) 21-\,2к — W2l,2k-l- (2.1.17) Тогда это решение единственно.
Доказател ьство. Рассмотрим соотношения (2.1.8) и (2.1.9) как систему уравнений ДЛЯ W2l-l,2k И W2l,2k-l (I, к = 1,2,..., f): 2fc-l Ik Cv)2i-i,2k + Cw2iflk-i = Pi„2+(A:_l)n+/ + Y Fq+i - Y Щ+1-1, 2k-і 2k C w2i-i,2k + Cw2i,2k-i = in2+(fc-i)f+; + X) Єщ+і - Y GVi+l-\ 1=1 2 = 1 Из этой системы находим _ F\n2+(k-l) +l , і / 1 0/ , 1 fi2 W2l-l,2k - 2C ± 2y C/(jn2+(fc-l)+/ - c I„2+(A._1)n+;, — F2"2+(fc-Df-+ -r- 1 /lor 1 rv2 " 2/,2 -1 - 2C + 2ycZL, 2+(fc-l)f+l с - 1,,2+( )-+,, где соответственно (2.1.18) (2.1.19) 2/fc-l 2k = 1 2=1 2fc-l 2fc (2.1.20) 1=1 1=1 Также будем принимать во внимание систему уравнений (2.1.10) и ее решение (2.1.11), рассмотренные в доказательстве теоремы 2.1.
Теперь разделим область определения wik на группы из четырех точек: (2/ -1,2/:-1), (21 — 1, 2к), (21, 2к — 1), (21, 2к) и будем рассматривать каждую такую группу независимо при / = 1,2,..,, Л=1,2,..,. Из (2.1.11) и (2.1.19) следует, что 2l-l,2k-l = Rlk ± П/с, W2l,2k-l = Rik =F Пл, (2.1.21) W2l-l,2k = Sik ± SjA;, W2f,2fc = - T Sjfc, (2.1.22) W2l,2k-1 = ї/fc ± hk, U 2l-l,2k — Flk T Jfc- (2.1.23) Здесь fc = 2C (2fc-l)+b Hfe = \ 2G(2k-l)%+l - (22Ц)+Ц 5«fc = 2C 2 t+ s fc = 2у c2 fcf+i - C Alfcf+l (2.1.24) Tjfc = 2С 1„2+(,._1)й+л ifc = i2Gin2 + (fc_1)+j - 2- n2 + (fc_])a+r
Опишем способ выбора значения в каждой точке. Из (2.1.21)-(2.1.23) видно, что значения в точках (21,2к — 1) и (21 — 1, 2к) посчитаны 2 раза как решения систем уравнений (2.1.10) и (2.1.18).
Рассмотрим сначала W2i,2k-i- Согласно (2.1.21) W2ipk-\ — Rik Т rik, и в то же время согласно (2.1.23) w2i,2k-i = Tlk ± fc. Рассмотрим множество f%2fc-i СОСТОЯ1Цее из ДВУХ элементов {Rlk + rtk,Rik — Г[к}, МНОЖеСТВО «4J) 2A._] = ( + 6 ,-%-} И МНОЖеСТВО u4j,2Jfe-l = {Tik + tik,Tlklk}. Видно, ЧТО W2/j2fc-l Є W .! И W2 ,2fc-1 Є й 2,2 -1 Отдельно рассмотрим случаи, когда множество w2[2k_ состоит из одного элемента. Это возможно, когда Гц- = 0. При этом w2i-i,2k-i = 2i,2fc-i = (выполняется условие (2.1.15) теоремы). Тогда ги2і-і,2 = «S.at-i \ W2i,2 -i " u 2i,2fc = й 2і,2 —і \ w2j-i,2fc- Таким образом, значения искомой функции в четырех точках определяются однозначно.
Такой же результат получается, если из одного элемента состоит множество w2j2fc-i или Щі2к-]_- Тогда u)2i-i,2k — и 2і,2к — Sfc (выполняется условие (2.1.16) теоремы) или соответственно W2i,2k-i — W2i-i,2k — Тік (выполняется условие (2.1.17) теоремы). тт с- -(1) -(2) -(3) Далее будем проводить рассуждения, предполагая, что w2l2k_l, w2i2k-v w2i2k-\ состоят из двух элементов. Возможны три случая: 1) Щі- к-і ГЇ "4/2fc-i — 0- В этом случае решение задачи не существует, что противоречит условиям теоремы; 2) {2%к-1 2%к-1 Ф 0 НО 21,2к-\ Ф «;M,2fc-l- В ЭТОМ Случае W2l,2k-l = Й 2 2А:-1П S 2l,2fc-1 -Тогда W2l-l,2k-l = 2/,k-l\W2«,2fc-l И W2l-l,2k = Шя,Й-1 \W2l,2k-l, Я 7/ ,2 = 4?,Ll\W2(-l,2t; 3) Й2І,2А-1 = 2?,2fe-l ЭТО ВОЗМОЖНО, КОГДа / = Ttk И 7 = t/jfc, а ЗНачИТ И W2l-l,2k-l = w2i,2k-i- (Это означает невыполнение условия (2.1.13) теоремы, иначе имеет место случай 2.) В этом случае нельзя однозначно определить значения в рассматриваемых точках, и НеобхОДИМО рассмотреть W2/-l,2fc Согласно (2.1.22) w2j-i,2fc = Slk±sik, и в то же время согласно (2.1.23) %-ід = Tlk tlk. Определим множества Щі-\,2к — {Rlk + rlk,Rik - rlk}, Щ?_і2к = {Sik + slk,Sik - slk} и 21-1,2Jfc = {Tlk + tlk,Tik Uk}- ВИДНО, ЧТО U)-2l-l,2k Є WS-l,2k И w2l-l,2k Є 4/-1,2 Снова возможны трн случая: 1) Щі-і 2к L-i 2к = 0- В этом случае решение задачи не существует, что противо речит условиям теоремы; 2) "4t-l,2fcn"4t-l,2k Ф 0 НО "4l-l,2fc Ф 2l-l,2fc- В ЭТ0М Случае Ш2І_1ій = 2?,2Л-1Пїї"4?-1,2А; Тогда , = Шя-1,2 \ад2/-1,2 И 2 ,2/:-1 = Щі-\,2к \w2L-l,2k, & W 2J-l,2fc-l = « 2I,2fc-l \ 2J,2fc-l5 3) 2 -1,2 = 2 -1,2 - T0 возможно, когда Sik = 7Ъ и 6- = tik, а значит и w2i,2k = w2i,2k-i- (Это означает невыполнения условия (2.1.14) теоремы, иначе имеет место случай 2.) В этом случае нельзя однозначно определить значения в рассматриваемых точках (2/ - 1,2к - 1), (21 - 1,2к), (21,2к-1), (21,2к). Теорема доказана. Вернемся к задаче 1 и задаче 2 определения гф, vfj и переформулируем для нес основные результаты.
Задача доплеровской томографии в случае нескольких измерений вдоль каждой прямой
Для численного решения задач доплеровской томографии в случае кусочно-постоянного поля на компьютере была разработана программа на языке программтрования Oclave (http://www.gnu.0rg/.H0ftware/0ctave/), совместимом с MATLAB кроссплатформенном языке для математических вычислений высокого уровня.
Программа состоит из следующих взаимодействующих между собой блоков: Функция, возвращающая значения v , vfj в узлах заданной сетки. Данная функция нужна для представления исходных данных прямой задачи. Ее использование позволяет решать прямую задачу для произвольного кусочно-постоянного поля
Блок решения прямой задачи. Этот блок необходим для получения исходных данных задачи восстановления кусочно-постоянного поля. Блок реализован в виде отдельной функции, и для каждой из задач (задача 1, задача 2 из параграфа 2.1 и задача 5 из параграфа 2.2) используется своя функция. Функция решения прямой задачи принимает на вход кусочно-постоянную векторную функцию V(x, у),а на выходе выдает набор данных Fm, Gm (для задачи 1 и задачи 2) или Fm , Gm(w) (для задачи 5), где \ц — параметр функции. Во все полученные данные вносится некоторая погрешность, величина которой является параметром вычислительного эксперимента.
Основная функция программы — восстановление векторного поля по данным измерений. Для решения каждой из задач (задачи 1, задачи 2 и задачи 5) реализована специальная функция восстановления кусочно-постоянного поля. Но общая схема работы всех трех функций одинакова. На вход принимаются данные измерений Fm, Gm (или Fm(/n), С?т( г),- (/12)і Ст(р2) в случае задачи 5), обрабатываются программой, и на выходе получаем восстановленное поле в виде 2п2 значений.
Блок вывода на экран результатов вычислений. Для вывода на экран используются встроенные функции для работы с графикой языка Octave. По завершении работы программы демонстрируются исходные данные прямой задачи и восстановленное в результате работы программы векторное ноле. На экран векторное поле выводится покомпонентно: отдельно выводится Vi(x, у) и отдельно v2(x, у). Также есть возможность вывести величину ошибки восстановления ноля в узлах сетки.
Далее рассмотрим алгоритмы численного восстановления поля, реализованные в описанной программе. 2.3.2 Алгоритм численного решения задачи 1. Для задачи 1, постановка которой приведена в 2.1, были доказаны теорема 2.3, описывающая ее множество неединственности, и теорема 2.5 об условиях существования решения. Опишем вычислительный алгоритм решения задачи 1, реализованный в программ Исходные данные для программы решения задачи восстановления поля есть набор п2 чисел F(pm, ipk), G(pm, ірк), к = 1, 2, m = 1,..., n2 — результаты измерений вдоль со ветствующих прямых, определенных в 2.1. Также считаем известными размеры облас D = [0, а] х [0, а]. Процесс восстановления поля проходит в три этапа: 1-й этап. По значениям F(pm, ipi), G(pm,(pi), т — 1,..., \п2 восстанавливаем ги- , i,j 1,..., гг. Алгоритм работает пошагово. На каждом шаге определяем два значе (і) :i) ш W, 21-1,к ш21,к , делая в общей сложности п шагов, где / — 1,2,...,, к = 1,2,.. Решаем при систему уравнений (2.1.11), получаем два решение (2.1.12). При отс: ствие каких либо априорных данных о поле V{x, у), случайным образом выбир одно из решений. Если же есть информация, например, о монотонности восста: ливаемого поля, то выбираем правильное решение, исходя из условий монотонно И так для всех / = 1,2,..., , к = 1, 2,..., п.
2-й этап. По значениям F(pm, p2), G(pm,ip2)i т = 1,..., п2 восстанавливаем іщ \ г, 1,...,гг. Делаем это точно точно таким же способом, какой приведен в оииса, 1-ого этапа.
3-й эт,ап. После первых двух этапов нам известны w\j\ ги- , i,j = 1,...,п. На тре: этапе решаем п2 систем уравнений