Диссертация (1145332), страница 11
Текст из файла (страница 11)
Несмотряна большое количество методов, задача определения коэффициента ослабления на практике сопряжена с рядом трудностей. Одна из них заключается вналичии в большинстве случаев сильного шума, обусловленного эффектамирассеяния и внутренними источниками.В Главе 1 предлагается подход, основанный на использовании внешнихисточников излучения специального типа.
Предполагается, что внешнее излучение имеет серийный (многократный) характер. Для его описания в качестве дополнительного аргумента функции, описывающей интенсивность излучения, вводится параметр. Предполагается, что производная интенсивности входящего излучения по параметру испытывает неограниченный рост приприближении к некоторому его значению. Указанное свойство и определяетспециальный характер внешних источников. В основе алгоритма определениякоэффициента ослабления положено уравнение переноса с параметрическойзависимостью. Соответственно, исследование прямой задачи является важным этапом решения задачи компьютерной томографии.Несмотря на то, что восстановление коэффициента ослабления в отдельном сечении исследуемого объекта происходит относительно быстро, восстановление трехмерной структуры по набору сечений может потребовать значительных вычислительных ресурсов.
К их увеличению может привести требование высокого качества разрешения. Таким образом, использование технологий параллельных вычислений для ряда задач компьютерной томографииможет являться весьма актуальным.К наиболее значимым результатам Главы 1 можно отнести следующие.1. Исследовано параметризованное уравнение переноса излучения, описывающее серийный (многократный) характер облучения. Изучены свойстварешения прямой задачи для параметризованного уравнения переноса при специальном типе внешнего источника.2. Разработан метод многократного облучения восстановления структу-70ры трехмерного объекта, основанный на использовании внешних источников специального типа. Доказана однозначная разрешимость задачи томографии.
Эффективность алгоритма реконструкции продемонстрирована натестовых примерах, содержащих неблагоприятные с вычислительной точкизрения факторы, к которым относится высокий уровень рассеяния (более90%) и сильное ослабление нерассеянного потока (более чем в 450000 раз).3. Разработаны алгоритмы параллельных вычислений, основанные на технологиях MPI и CUDA, восстановления трехмерной структуры объекта. Исследована эффективность параллелизации, которую обеспечивают предложенные алгоритмы.4.
Разработано специализированное программное обеспечение для решения прямых задач для уравнения переноса на основе метода Монте-Карлои решения задач компьютерной томографии на основе метода многократного облучения и алгоритма свертки и обратной проекции. Осуществлена программная реализация алгоритмов параллельных вычислений решения задачикомпьютерной томографии.Основные результаты Главы 1 представлены в работах [14, 37, 38, 43, 47,87, 101–103, 105, 139]. Результаты в [37, 38, 47, 139] получены непосредственноавтором. В остальных (совместных) работах творческий вклад по отношениюк соавторам как минимум равный.
Так, в монографии [14] из 19 параграфовавтором полностью написаны пять (§1.1, 1.6, 1.7, 1.8 и 2.9), три параграфа написаны совместно (§1.2, 1.3, 1.9). В англоязычной монографии [103],близкой по структуре к [14], творческий вклад аналогичен. В совместныхработах [101, 102] автором написаны разделы, в которых обосновывается метод многократного облучения и проводится тестирование алгоритма. В работе [87] автором разработан параллельный алгоритм решения трехмернойзадачи компьютерной томографии.71ГЛАВА 2. ПЕРЕНОС ПОЛЯРИЗОВАННОГО ИЗЛУЧЕНИЯ ВТРЕХМЕРНОЙ НЕОДНОРОДНОЙ СРЕДЕ§3.Прямые и обратные задачи для уравнения переносаполяризованного излучения3.1.Постановка и исследование прямой задачи для векторногоуравнения переноса излученияОсновной характеристикой поляризованного излучения является четырехмерный вектор-параметр Стокса f =(f1 , f2 , f3 , f4 ).
Соответствующее уравнение переноса в изотропной среде для этого вектора имеет вид [73, 77, 91]:Zω · ∇r f (r, ω) + µ(r)f (r, ω) =P (r, ω, ω 0 )f (r, ω 0 )dω 0 + J(r, ω),(2.1)Ωгде r = (r1 , r2 , r3 ) ∈ G, G – выпуклая ограниченная область в трехмерномевклидовом пространстве R3 , ω ∈ Ω = {ω ∈ R3 : |ω| = 1}. Здесь J(r, ω) – четырехкомпонентный вектор, описывающий внутренние источники излучения,µ(r) – коэффициент ослабления, P (r, ω, ω 0 ) – матрица рассеяния размерности4 × 4.Для характеристики неоднородности среды G будем использовать введенное в §1.1 разбиение G0 области G.
Множество G0 предполагается открытыми плотным в G, то есть G0 = G. Кроме того, оно является объединениемконечного числа непересекающихся областей Gi , которые можно интерпретировать как некоторые части неоднородной среды G, заполненные i-м веществом. Будем полагать, что множество G0 обладает свойством обобщеннойвыпуклости [25], то есть любой луч Lr,ω , исходящий из точки r ∈ G0 в направлением ω ∈ Ω, пересекает ∂G0 в конечном числе точек.72Как и в предыдущей главе через Cb (X), X ∈ Rm будем обозначать банахово пространство функций, определенных на X, ограниченных и непрерывныхна X, с нормой||f || = sup |f (x)|.x∈X(4)Аналогично определим пространство Cb (X), образованное вектор-функциями f =(f1 , f2 , f3 , f4 ), каждая компонента которой принадлежит Cb (X) исоответствующая норма определяется равенством||f ||4 = max ||fi ||.1≤i≤4Относительно коэффициентов уравнения (2.1) будем предполагать следующее.
Функция µ(r) – неотрицательная и принадлежит пространству Cb (G0 ),(4)а вектор-функция J(r, ω) ∈ Cb (G0 ×Ω). Все компоненты матрицы P (r, ω, ω 0 )принадлежат Cb (G0 × Ω × Ω).В дальнейшем будем использовать некоторые обозначения, введенные в§1.1. Напомним, что функция d(r, ω) описывает расстояние от точки r ∈ Gдо границы ∂G = G \ G в направлении ω, множества Γ− и Γ+ являютсяобластями определения входящего и выходящего излучения соответственно.К уравнению (2.1) присоединим граничное условиеf (ξ, ω) = h(ξ, ω),(ξ, ω) ∈ Γ− .(2.2)Вектор-функция h(ξ, ω) описывает характеристики входящего в среду G излучения. Учитывая определения множества Γ− и функции d(r, ω), краевоеусловие (2.2) можно переписать в видеf (r − d(r, −ω)ω, ω) = h(r − d(r, −ω)ω, ω),(r, ω) ∈ G0 × Ω.(2.20 )При решении задач томографии мы будем иметь дело с функциями h, которые могут иметь разрывы в области своего определения.
Для этого введем73конечное разбиениеΩ0 =q[Ωi ∩ Ωj = ∅,Ωi ,i 6= j,i=1где Ωi – области на единичной сфере Ω, Ω0 = Ω.Пусть функцияh(ξ, ω)такова, что ее продолжениеeh(r, ω) =(4)h(r − d(r, −ω)ω, ω) ∈ Cb (G0 × Ω0 ). Простейший пример функции h, удовлетворяющий вышеприведенным ограничениям, может иметь следующий вид:h(ω) = (h1 (ω), 0, 0, 0), гдеh1 (ω) = sgn(ω3 ) + 1,sgn(x) = 1,x > 0,0, x = 0, −1, x < 0.В этом примере Ω0 = Ω1 ∪ Ω2 , где Ωk = {ω ∈ Ω : sgn(ω3 ) = (−1)k }, и(4)eh(r, ω) = h(ω) ∈ C (Ω0 ).bПрямая задача (2.1), (2.2). Задачу определения функции f из (2.1),(2.2) при известных µ, P, J, h будем называть прямой задачей (2.1), (2.2).Обозначим(lf )(r, ω) = ω · ∇r f (r, ω) + µ(r)f (r, ω),ZN (r, ω) = (Sf )(r, ω) = P (r, ω, ω 0 )f (r, ω 0 )dω 0 .(2.3)(2.4)ΩПоскольку все компоненты матрицы P (r, ω, ω 0 ) есть ограниченные и непрерывные функции по переменным ω, ω 0 , то согласно известным свойствам интегральных операторов [23] заключаем, что оператор S, определенный через(4)(4)(2.4), переводит пространство Cb (G0 × Ω0 ) в Cb (G0 × Ω).
Таким образом,векторный интеграл столкновений N (r, ω) обладает сглаживающими свойствами по угловой переменной ω. Это свойство типично в теории переноса [23, 25] и нередко используется при решении ряда задач компьютернойтомографии [99, 101, 103, 104].Пусть везде далее Ω0 обозначает открытое подмножество единичной сфе-74ры плотное в Ω. Определим класс D, в котором будем искать решение уравнения (2.1).Определение 2.1. Вектор-функция f (r, ω) ∈ D(G0 × Ω0 ), если для любыхточек (r, ω) ∈ G0 × Ω0 функция f (r + tω, ω) абсолютно непрерывна по переменной t ∈ [−d(r, −ω), d(r, ω)] и функции f (r, ω), ω · ∇r f (r, ω) принадлежат(4)пространству Cb (G0 × Ω0 ).Отметим следующее.
Так как функция µ(r) ∈ Cb (G0 ), то оператор l, опре(4)деленный равенством (2.3), переводит D(G0 × Ω0 ) в Cb (G0 × Ω0 ).Определение 2.2. Решением прямой задачи (2.1), (2.2) будем называтьфункцию f (r, ω) ∈ D(G0 × Ω0 ), которая для всех (r, ω) ∈ G0 × Ω0 удовлетворяет соотношениям(lf )(r, ω) = (Sf )(r, ω) + J(r, ω),f (r − d(r, −ω)ω, ω) = h(r − d(r, −ω)ω, ω).Приведем необходимое для дальнейших рассуждений следующее утверждение [103, 104].Лемма 2.1. Пусть множество G0 удовлетворяет условию обобщенной выпуклости, тогда если функция ϕ(r, ω) ∈ Cb (G0 × Ω0 ), тоd(r,−ω)Zϕ(r − ωt, ω)dt ∈ Cb (G0 × Ω0 ),Φ(r, ω) =0ZtΨ(r, ω, t) =ϕ(r − ωt0 , ω)dt0 ∈ Cb (G0 × Ω0 × [0, d(r, −ω)]).0Следствие 2.1.
Функции оптического расстояния, определенные формуламиZtτ (r, ω, t) =0µ(r − t0 ω)dt0 ,d(r,−ω)Zµ(r − tω)dt,τ (r, ω) =075принадлежат соответственно пространствам Cb (G0 × Ω × [0, d(r, −ω)])и Cb (G0 × Ω).Лемма 2.2. Пусть функция f (r, ω) ∈ D(G0 × Ω0 ), тогда(4)f ± (r, ω) = f (r ± d(r, ±ω)ω, ω) ∈ Cb (G0 × Ω0 ).Доказательство. Так как f (r, ω) ∈ D(G0 × Ω0 ), то функции fi (r + tω, ω),i = 1, 4, абсолютно непрерывны по t ∈ [−d(r, −ω), d(r, ω)] для любых (r, ω) ∈G0 × Ω0 , следовательноd(r,−ω)Zω · ∇r fi (r − ωt, ω)dt = fi (r − d(r, −ω)ω, ω) − fi (r, ω).0Учитывая, что функция в левой части этого равенства согласно лемме 2.1.принадлежит Cb (G0 × Ω0 ) и fi (r, ω) принадлежат Cb (G0 × Ω0 ), заключаем, что(4)fi (r − d(r, −ω)ω, ω) ∈ Cb (G0 × Ω0 ).
Следовательно, f − (r, ω) ∈ Cb (G0 × Ω0 ).(4)Аналогично показывается, что f + (r, ω) ∈ Cb (G0 × Ω0 ).Из леммы 2.2. вытекает очевидное следствие.Следствие 2.2. Пусть функции ϕ ∈ D(G0 × Ω), ψ ∈ D(G0 × Ω0 ), тогда(4)ϕ± (r, ω) = ϕ(r ± d(r, ±ω)ω, ω) ∈ Cb (G0 × Ω),(4)ψ ± (r, ω) = ψ(r ± d(r, ±ω)ω, ω) ∈ Cb (G0 × Ω0 ).Приведем еще одно утверждение из [103, 104].Лемма 2.3. Пусть функция ϕ(r, ω) ∈ Cb (G0 × Ω0 ), тогда функцияd(r,−ω)Zexp(−τ (r, ω, t))ϕ(r − tω, ω)dt ∈ Cb (G0 × Ω0 ),ψ(r, ω) =0причемω · ∇r ψ(r, ω) ∈ Cb (G0 × Ω0 ),76ω · ∇r ψ(r, ω) + µ(r)ψ(r, ω) = ϕ(r, ω), ϕk ψ k≤ 1 − e−µd k k,µгде µ = sup µ(r), d - диаметр области G.r∈G0В дальнейшем мы будем использовать условия, вызванные физическими ограничениями, на компоненты функции f = (f1 , f2 , f3 , f4 ) – решенияуравнения (2.1). Функции fi (r, ω) являются параметрами Стокса и должныудовлетворять неравенствам [24, 91]f1 ≥ 0,f12 ≥ f22 + f32 + f42 .(2.5)(4)Обозначим через K – конус в пространстве Cb (G0 × Ω0 ), образованный(4)функциями f (r, ω) = (f1 (r, ω), f2 (r, ω), f3 (r, ω), f4 (r, ω)) ∈ Cb (G0 × Ω0 ), удовлетворяющими условиям (2.5).