автореферат диссертации по информатике, вычислительной технике и управлению, 05.13.18, диссертация на тему:Математическое моделирование процессов переноса излучения в многослойных средах с подвижными рассеивателями

На правах рукописи

СТАРУХИН Павел Юрьевич

МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ ПРОЦЕССОВ ПЕРЕНОСА ИЗЛУЧЕНИЯ В МНОГОСЛОЙНЫХ СРЕДАХ С ПОДВИЖНЫМИ РАССЕИВАТЕЛЯМИ

Специальность 05.13.18 - Математическое моделирование, численные методы и комплексы программ

Автореферат диссертации на соискание ученой степени кандидата физико-математических наук

2 НЮН 2011

Саратов-2011

4848987

Работа выполнена в ГОУ технический университет»

ВПО «Саратовский государственный

Научный руководитель: доктор физико-математических наук,

профессор Клинаев Юрий Васильевич

Официальные оппоненты: доктор физико-математических наук,

Заслуженный деятель науки РФ, Лауреат государственной премии РФ, профессор Синицын Николай Иванович

доктор физико-математических наук, профессор Губатенко Валерий Петрович

Ведущая организация: Саратовский государственный

университет им. Н.Г. Чернышевского

Защита состоится 15 июня 2011 года в 13-00 часов на заседании диссертационного совета Д 212.242.08 при ГОУ ВПО «Саратовский государственный технический университет» по адресу: 410054, г. Саратов, ул. Политехническая, д.77, корп. 1, ауд. 319.

С диссертацией можно ознакомиться в научно-технической библиотеке ГОУ ВПО «Саратовский государственный технический университет».

Автореферат размещён на сайте Саратовского государственного технического университета www.sstu.ru 13 мая 2011 г.

Автореферат разослан « 13 » мая 2011г.

Ученый секретарь диссертационного совета

Терентьев А. А.

ОБЩАЯ ХАРАКТЕРИСТИКА РАБОТЫ

Актуальность исследований. В настоящее время одно из актуальных научных направлений связано с анализом влияния физиологических и структурных свойств биологических объектов на их физические и, в частности, оптические свойства [Приезжее А.В., Тучин В.В., Ульянов С.С.]. В этой связи важное значение приобретает математическое моделирование физических процессов, протекающих в объекте исследования в условиях лазерного облучения, при этом возникает задача решения интегродифференциального уравнения теории переноса излучения (ТПИ) (Чандрасекар, Исимару).

Одним из часто используемых методов численного решения задач ТПИ является метод статистического моделирования Монте-Карло (М. Кейзер, С. Жак, И.В. Ярославский, И.В. Меглинский, Ф. де Мул, М. Кулинк и др.).

К настоящему моменту разработаны программные средства для расчета распределения рассеянного излучения в неоднородных средах, в том числе с учетом локальных неоднородностей, например опухолей. Существует коммерческий программный комплекс ТгасеРго фирмы Lambda Research Corporation для расчета оптических систем, позволяющий проводить трассировку лучей методом Монте-Карло. Изучены особенности динамики импульсов рассеянного излучения (A. Popov, A. Priezzhev, and R. Myllyla) и возможности диагностики скорости потоков жидкостей доплеровским методом (A.V. Bykov, M.Yu. Kirillin, A.V. Priezzhev, R.A. Myllyla).

Несмотря на значительные достижения исследователей, ряд важных как в теоретическом, так и в практическом плане задач рассмотрен недостаточно.

В частности, не достаточно изучены возможности моделирования многослойных сред с различными видами циркуляции биожидкостей. Не достаточно полно исследовано влияние подвижных рассеивателей на характеристики рассеянного излучения. Не изучены некоторые возможности применения доплеровских методов исследования для задач томографии и диагностики. Также актуальным является вопрос разработки новых эффективных вычислительных методик решения задач ТПИ с использованием современных средств массивно-параллельных вычислений. Перечисленные факты позволили сформулировать цель настоящей работы.

Целью работы является разработка математических моделей и высокопроизводительных численных алгоритмов исследования процессов распространения лазерного излучения в сильно рассеивающих и поглощающих многослойных биологических средах, содержащих подвижные рассеиватели, а также разработка моделей доплеровской томографии биологических объектов. Для достижения поставленной цели сформулированы следующие задачи:

• Разработка математических моделей распространения оптического излучения в многослойных, рассеивающих и поглощающих средах, учитывающих возможное уширение спектра рассеянного излучения вследствие эффекта Доплера, обусловленного рассеянием излучения подвижными частицами среды.

• Разработка численных алгоритмов решения стационарного и нестационарного уравнений переноса излучения в непрозрачных средах на основе метода моделирования Монте-Карло, позволяющего учитывать чрезвычайно большие объёмы (~1010) рассеиваемых пакетов фотонов.

• Реализация разработанных алгоритмов моделирования в виде комплекса программ для многопроцессорных вычислительных систем.

• Анализ возможности применения современных многопоточных центральных процессоров и многопоточных графических процессоров для производства параллельных вычислений по методу Монте-Карло.

Научная новизна работы:

• Построена математическая модель рассеивающей и поглощающей многослойной среды, отличающаяся учётом оптической анизотропии, присутствия подвижных рассеивателей в форме сети микроциркуляции и уединенного крупного сосуда. Разработана модель учёта уширения спектра излучения рассеянного подвижными частицами жидкости.

• Разработаны численные алгоритмы, реализующие процесс моделирования переноса излучения на основе метода Монте-Карло.

• Предложены эффективные подходы к реализации моделирования переноса излучения с использованием многопоточных графических процессоров, позволившие сократить время моделирования на 2 порядка.

• Разработан проблемно-ориентированный комплекс программ для моделирования распространения излучения в рассеивающих средах.

• Проведен ряд вычислительных экспериментов, в результате которых установлены следующие функциональные зависимости:

-при непрерывном характере облучения соотношение интенсивности доплеровской фракции излучения к полной интенсивности рассеянного излучения квадратично растет с увеличением объемной концентрации подвижных рассеивателей в поверхностных слоях среды;

-эллиптичность контура равной интенсивности доплеровской фракции рассеянного излучения линейно убывает с увеличением глубины залегания уединенного сосуда с жидкостью в поверхностных слоях биоткани;

-при импульсном режиме облучения двухслойной среды со статичным верхним слоем и вторым слоем, содержащим подвижные рассеиватели в форме сети микроциркуляции, временная задержка между появлением на детекторе доплеровской и несмещенной компонент рассеянного излучения линейно возрастает с увеличением толщины верхнего статичного слоя. Научная и практическая значимость.

Научная значимость работы заключается в применимости разработанного на основе методов Монте-Карло проблемно-ориентированного комплекса программ моделирования распространения излучения в мутных средах для исследования распределения интенсивности рассеянного излучения и уширения спектра рассеянного излучения вследствие эффекта Доплера.

Практическая значимость заключается в полученных результатах моделирования и установленных закономерностях поведения доплеровской компоненты рассеянного излучения при различных условиях.

Достоверность предоставленных научных результатов подтверждается тем, что полученные результаты отличаются непротиворечивостью и находятся в хорошем соответствии с результатами исследований других авторов и экспериментальными данными, опубликованными в научной литературе.

Положения и результаты, выносимые на защиту:

1. Предложенная математическая модель процесса распространения излучения в неоднородных средах учитывает многослойное строение биоткани, оптическую анизотропию среды, явления многократного рассеяния, поглощения, отражения и преломления излучения, наличие подвижных рассеивателей в форме сети микроциркуляции и в форме уединенного крупного сосуда.

2. Программный комплекс «Mcml-doppler» позволяет решать задачи стационарной и нестационарной ТПИ, рассчитывать пространственно-временное распределение интенсивности рассеянного излучения и характеристики доплеровского уширения спектра излучения рассеянного подвижными частицами.

3. Увеличение объемной концентрации подвижных рассеивателей в многослойной среде с микроциркуляцией приводит к увеличению соотношения интенсивности доплеровской фракции рассеянного излучения и значения полной интенсивности рассеянного излучения.

4. Глубина залегания протяженного сосуда с жидкостью в среде определяет форму распределения интенсивности доплеровской фракции рассеянного излучения, а именно эллиптичность контура равной интенсивности доплеровской фракции рассеянного излучения убывает с увеличением глубины залегания сосуда.

5. При распространении короткого импульса излучения в многослойной среде со статичным поверхностным слоем и с подвижными рассеивателями в нижележащих слоях длительность задержки между поступлением смещенной по частоте и несмещенной компоненты рассеянного излучения возрастает с увеличением толщины поверхностного слоя.

6. Разработанный алгоритм моделирования переноса излучения позволяет распараллелить вычислительный процесс в требуемом масштабе для эффективного использования многоядерных/многопоточных CPU, и многократно сократить время вычислительного эксперимента при использовании современных графических процессоров, либо увеличить точность вычислений при заданном времени, затрачиваемом на процесс моделирования.

Апробация работы. Основные результаты диссертационной работы докладывались на российских и международных конференциях: Молодежная научная школа по оптике, лазерной физике и оптоэлектронике «Проблемы оптической физики» (Саратов, Россия, 1997), Международная научно-техническая конференция «Актуальные проблемы электронного приборостроения» (Саратов, Россия, 2010), Всероссийская научно-практическая конференция молодых ученых «Инновации и актуальные проблемы техники и технологий 2010» (Саратов, Россия, 2010), Ежегодная Всероссийская научная школа-семинар «Методы компьютерной диагностики в биологии и медицине 2010» (Саратов, Россия, 2010), XX Ежегодная международная конференция-выставка «Информационные технологии в образовании «ИТО - 2010» (Москва, Россия, 2010).

Публикации. По теме диссертации опубликовано 20 печатных работ, в том числе пять статей - в изданиях, рекомендованных ВАК РФ. Программный комплекс моделирования зарегистрирован ИНИМ РАО РФ / Объединённый фонд электронных ресурсов «Наука и образование» 02.12.2010г. Список 13 публикаций приведен в конце автореферата.

Jltvtep Hü-N¿ 633hm

üclcKíup спсктроашлншор

Структура и объем работы. Диссертационная работа состоит из введения, 4 глав, заключения, списка цитируемой литературы из 128 наименований и 4 приложений. Диссертация содержит 7 таблиц, иллюстрирована 43 рисунками. Общий объем диссертационной работы составляет 162 страницы.

КРАТКОЕ СОДЕРЖАНИЕ РАБОТЫ

Во введении дана общая характеристика работы: обоснована актуальность, отмечена научная новизна и практическая значимость работы. Формулируются цели и задачи исследования. Проанализированы и обоснованы применяемые в работе методы исследования.

В первой главе изложены принципы применения методов светорассеяния для задач

биомедицинской диагностики

состояния биологических объектов на основе анализа их оптических свойств, проявляющихся при облучении. Приведен обзор работ, посвященных этой проблематике. Изложена содержательная и математическая постановка задачи. На рис 1 представлена традиционная схема организации оптического зондирования объекта исследования - кожного покрова человека in vivo. Рассеивающая среда облучается направленным пучком

монохроматического излучения, и

Рис. 1. Схема многослойной модели кожи и моделируемого эксперимента

отклик среды на исходное возмущение анализируется с помощью комплекса детекторов диффузно отраженного излучения. Каждый плоскопараллельный слой среды на рис. 1 описывается следующими величинами, зависящими от длины волны падающего излучения:

- коэффициент поглощения (величина, обратная расстоянию, на котором интенсивность уменьшается за счет поглощения в е раз) [1/см]; ^ - коэффициент рассеяния [1/см]; (1)

g- фактор анизотропии или средний косинус угла рассеяния; п - показатель преломления; ¿/-толщинаслоя [см]).

Процессы распространения излучения в подобных средах описываются в рамках теории переноса излучения (ТПИ). Интегро-дифференциальное уравнение нестационарной ТПИ записывается следующим образом:

ds dt 4л-,

\l(r,s\Of(t,t'W

p{s,s')d£í+S(r,s,t)

(2)

где 1{г, {) - лучевая интенсивность в точке г в направлении я, в момент времени ([Вт'М~2-ср-1] (ср - стерадиан), (?=(/;, с)"' - среднее время между взаимодействиями,

fi,=fis+fia - полный коэффициент ослабления или экстинкции (в биологических средах fJs»/./„), f(t,t') описывает временную деформацию ¿-образного импульса после акта рассеяния, p(s,s') - фазовая функция рассеяния, dO! - единичный телесный угол в направлении s', S(r,s,t) описывает функцию источника излучения в среде, с - скорость света в среде. Уравнение (2) дополняется начальными условиями:

l(r,~s)\„0=4>(7,~s). <2')

Исключение из (2) переменной t, приводит к уравнению стационарной ТПИ.

^н! = + . (3)

as 4/г

Если процесс переноса излучения исследуется в области G с R3, a dG - граница области, то граничные условия на dG могут быть записаны в виде

(,,„<о = Sir's) + RI (г, (J,0>0, (3*:

где reôG,n- внешняя нормаль к ÔG, S(r,s) - лучевая интенсивность падающего излучения, R- оператор отражения. При наличии в области G отражающих и преломляющих поверхностей условия, аналогичные (3') должны быть записаны для каждой из них. Данное граничное условие отражает тот факт, что входящее в G излучение складывается из излучения внешних источников и той части вышедшего из G излучения, которое снова попало в G за счет отражения, описываемого оператором R (например, Френелевское отражение).

Уравнения (2-3) носят эвристический характер и показывают, что при отсутствии в среде источников интенсивность потока излучения, проходящего через элементарную площадку s, ослабевает за счет поглощения и рассеяния (первый член в правой части (2)), а также может возрасти за счет многократного рассеяния (второй член в правой части (2)).

ТПИ предусматривает постановку прямой и обратной задачи. Прямая задача ТПИ состоит в нахождении распределения отраженного и поглощенного излучения по известным параметрам рассеивающей среды. Обратная задача заключается в определении характеристик сред по известному распределению интенсивности возбуждающего импульса, что представляет значительный интерес в областях диагностики и томографии объектов.

Недостатки известных методов решения прямой задачи ТПИ связаны, в частности, с невозможностью учета подвижных рассеивателей в среде и возникающим в связи с этим уширением спектра рассеянного излучения за счет эффекта Доплера, а также с ограничениями по учету граничных условий.

Поэтому для решения поставленной прямой задачи и учета эффектов, обусловленных многослойностыо среды и присутствием подвижных рассеивателей, целесообразно использование методов статистического моделирования распространения излучения, например, методов Монте-Карло.

Поскольку ТПИ оперирует терминами лучевой интенсивности, а не напряженности поля излучения, можно рассматривать интенсивность как поток дискретных фотонов через единицу площади в единицу времени.

Тогда суть применения метода Монте-Карло для нахождения решения прямой задачи теории переноса излучения выразится в следующем. Непрерывное излучение рассматривается как поток N дискретных частиц (фотонов), каждая из которых обладает энергией. Это дает возможность представить процесс переноса излучения в виде последовательности событий (рассеяние, поглощение и т.д.), происходящих с каждым из N фотонов и в результате статистического моделирования большого числа траекторий рассчитать вероятное распределение энергии, а также другие необходимые функционалы от решения уравнения ТПИ.

Метод статистического моделирования Монте-Карло включает следующие основные этапы: формирование модели фотона на границе среды, расчет траектории, учет поглощения, прекращение отслеживания траектории фотона, регистрация фотона детектором на границе среды.

В простейшем виде в процессе моделирования можно отслеживать траектории N фотонов, распространяющихся в среде. Тогда на каждом шаге моделирования требуется решать, будет ли фотон рассеян или поглощен средой. В случае поглощения информация о рассчитанной траектории теряется и осуществляется переход к расчету следующей траектории фотона. Для того чтобы не терять информацию о промежуточных траекториях, используется технология «пакетов фотонов»: считается, что по траектории движется не один, а множество фотонов, часть которых поглощается на очередном шаге моделирования, а оставшиеся продолжают движение.

Тогда на каждом шаге алгоритма моделирования выполняются действия:

1. Формирование очередного пакета фотонов на границе среды.

Плоский возбуждающий импульс длительности Т, соответствующий члену 5,('".-г>') уравнения (2), моделируется путем облучения поверхности среды N пакетами фотонов, для каждого из которых время влета ¡0 в среду равномерно распределено в диапазоне от 0 до Т и определяется

г0=(1-^)-Г. №

Здесь и далее случайное число, равномерно распределенное на отрезке [0,1]. Пространственное и угловое распределение пакетов фотонов в импульсе соответствует профилю падающего излучения. Например, для точечного источника, излучающего вдоль оси 2 по нормали к поверхности среды, начальные координаты нового пакета фотонов:

*о=Л=г,=0. (5>

Исходные направляющие косинусы волнового вектора пакета фотонов:

Безразмерный параметр, отражающий энергию поступающего в неоднородную среду пакета фотонов, называется весом пакета и имеет начальное значение:

Г0= 1.0. (7)

2. Расчет траектории.

Поглощение излучения в однородном слое среды описывается законом Ламберта-Бугера-Бера. Ослабление интенсивности за счет рассеяния также носит

экспоненциальный характер. Поэтому ослабление интенсивности излучения в заданном направлении с учетом поглощения и рассеяния будет определяться коэффициентом полного ослабления (экстинкции) //, = м„ + ns и при распространении излучения в направлении оси z можно записать

/ = /0ехр(-Л2). (8)

Из (8) можно получить выражение для моделирования длины свободного пробега пакета фотонов /„ до очередного взаимодействия с частицей среды с использованием случайного числа:

/ =JïiÎL, (9)

/Л

что отражает ослабление интенсивности излучения в соответствии с первым членом из правой части (2).

Дальнейшее направление движения будет определяться фазовой функцией рассеяния pis,s') из второго члена правой части (2).

При переходе к полярным координатам направление определяется двумя углами рассеяния - полярным 9 и азимутальным <р. Тогда можно записать

/;(?,?) =р(в)р(<р). 0°)

Полярный угол рассеяния в определяет отклонение направляющего вектора рассеянного пакета относительно исходного направления. В однородном слое биологической среды полярный угол фазовой функции однократного рассеяния хорошо описывается функцией Хени-Гринштейна:

р{в) =--(11)

Расчет полярного угла в в соответствии с (11), осуществляется по формуле:

cos(É>) = —

ч

1

l-g+2^:

з J

(12)

В случае, если рассеивающая среда характеризуется изотропным рассеянием, т.е. g=0, угол в можно определить на основе выражения

cos(0) = 2£-l. (13)

Азимутальный угол рассеяния (р в однородном слое среды считается произвольным и равновероятным в диапазоне от 0 до 2л:

На основе (9)-(14) формируется модель учета многократного рассеянного излучения в соответствии со вторым членом из правой части (2).

Направляющие косинусы //,', fty', единичного вектора скорости пакета фотонов ks после рассеяния могут быть определены из соотношений:

Мг= {f, f. c°s <Р ~ Му sin <p)+ M, eos в

v-м,

M,= i , \MyM, cos <р + цх sin <р)+ My cos в

VI

M,'= -sinecosqhjl-M¡* + M, cosí?,

где //,, // , цг - направляющие косинусы единичного вектора скорости пакета фотонов Л/ до рассеяния.

Координаты пакета после рассеяния определяются по формулам:

х'=х+(гхЧсм

у=у+м;-1с, <16)

Время, затрачиваемое пакетом фотонов на перемещение между последовательными рассеяниями, учитывается следующим образом:

.ЬЛ. (17)

с

где /с> - длина пробега между рассеяниями, п - показатель преломления среды. Учет доплеровского ушмрения спектра рассеянного излучения При наличии в среде подвижных частиц и в случае однократного рассеяния фотона (пакета фотонов) на частице, движущейся со скоростью v, возникает эффект Доплера, учитывать который можно в форме сдвига частоты излучения следующим образом:

^ = О»)

где 8/ - сдвиг частоты излучения, к/, - единичные векторы скорости фотона «до» и «после» рассеяния, V - вектор скорости движущейся частицы.

V/

/ / /

- / .

е

к,

Рис. 2. Эффект Доплера при рассеянии излучения подвижным рассеивателем

Если в процессе распространения фотон испытывает и рассеяний движущимися частицами, его суммарный сдвиг частоты вычисляется как:

б/ = . В результате расчета распространения в среде заданного числа У = 1

пакетов фотонов (106 - Ю10), детектора заданной геометрии достигнут т пакетов

т

фотонов с общим весом: ^полное ~ У, Щ . Из них к < т испытавших

»=1

доплеровский сдвиг частоты. Тогда средний доплеровский сдвиг частоты

компоненты излучения, рассеянного движущимися частицами: № =-¡=1-.

/.1

В средах биологических подвижные рассеиватели могут встречаться либо в крупных кровеносных сосудах (артерии, вены и т.д.), либо в виде распределенной сети микроциркуляции. В первом случае неоднородность моделируется 10

цилиндром заданного диаметра, и все рассеяния пакета фотонов в пределах неоднородности считаются доплеровскими. Во втором случае каждый слой многослойной среды характеризуется значением относительной объемной концентрации подвижных рассеивателей:

где Уподв - объем, занимаемый подвижными частицами в слое среды, Vсреды -общий объем частиц.

Пространственное распределение подвижных рассеивателей в каждом слое среды считалось случайным и равномерным. Поэтому, если в процессе моделирования случайная величина равномерно распределенная в интервале [0,1], £ < Сслоя, пакет фотонов считался рассеянным движущейся частицей и приобретал доплеровский сдвиг частоты в соответствии с (18).

Учет граничных условий на разделе сред с различными значениями показателя преломления

Интенсивности отраженной и преломленной компонент излучения на границе показателя преломления рассчитываются на основе закона преломления и формул Френеля:

где и,,и, - показатели преломления среды, из которой поступает излучение, и среды в которую проникает излучение, соответственно, а.,а, - углы падения и преломления, Л, - коэффициент внутреннего отражения.

Если при расчете новой координаты местоположения оказалось, что пакет фотонов пересекает границу однородного слоя среды и попадает в область (другой слой или локальную неоднородность, например, кровеносный сосуд) с другими значениями оптических параметров (показателя преломления, коэффициентов рассеяния, поглощения), то рассчитывается вероятность отражения или преломления и координаты пакета фотонов с учетом (20).

Если случайное число считается, что пакет фотонов отражен на

границе показателя преломления, иначе происходит преломление и пакет попадает в соседний слой многослойной среды или другую область неоднородности. 3. Поглощение.

На входе в непрозрачную среду вес каждого пакета фотонов задается в соответствии с (7). За счет поглощения и рассеяния после каждого /-го акта взаимодействия со средой вес пакета уменьшается в соответствии с поглощающими и рассеивающими свойствами среды по формуле

среды

(19)

я- бш а, = п, 51п а,, 5т2(а, - а,) ^ Гап2(а, - а,) 2 [_51П2(£Г,. +аг,) \т1(а1 +а,)

(20)

(21)

Поглощенная доля энергии пакета, равная учитывается при расчете

Л

пространственного распределения поглощенного излучения. Таким образом, рассчитывается ослабление потока излучения за счет поглощения и рассеяния, что соответствует первому члену из правой части уравнения (2).

4. Прекращение отслеживания траектории пакета.

Если после 1-го рассеяния вес пакета фотонов становится меньше некоторого минимума, считается, что он полностью поглощен средой и осуществляется возврат к п. 1. алгоритма, иначе производится переход к п. 2.

5. Регистрация пакета фотонов детектором излучения на границе среды.

Если в результате преломления на границе верхнего слоя (или нижнего при ограниченной толщине среды) пакет фотонов пересек границу среда-воздух, его текущие параметры (координата, направление движения, вес, сдвиг частоты) регистрируются в памяти компьютера.

После повторения шагов 1-5 для достаточного количества N траекторий пакетов фотонов карта распределения энергии отраженного, пропущенного и поглощенного излучения рассчитывается и накапливается в компьютере.

С учетом нормировки на суммарную энергию исходного возбуждающего импульса Ы1¥о можно записать

• А«+т* + к<1 =1 (22) Ш,

где Аш -суммарная доля поглощенной средой энергии из всех N построенных траекторий пакетов, и - суммарные доли диффузно-отраженной и диффузно-пропущенной энергии соответственно.

Выражение (22) отражает суть закона сохранения энергии. В результате проведения ряда статистических испытаний, формируется пространственно-временное распределение интенсивности и среднего значения доплеровского сдвига частоты рассеянного излучения.

Разработанный алгоритм моделирования методом Монте-Карло предназначен для численного решения прямых задач ТПИ в многослойных рассеивающих и поглощающих средах. Решением классической задачи ТПИ будет являться распределение интенсивности падающего излучения (например, расчет параметров , называемых также коэффициентами диффузного

пропускания и отражения соответственно, или пространственное распределение интенсивности отраженного и пропущенного излучения на границах среды). Учет эффекта Доплера при рассеянии излучения подвижными частицами будет выражаться в расчете пространственного распределения среднего значения сдвига частоты рассеянного излучения, или в расчете более информативной гистограммы частот (спектра) рассеянного излучения.

При исследовании прохождения импульсов излучения конечной длительности в рассеивающих средах задача ТПИ приобретает выраженную зависимость от времени и решается в рамках нестационарной ТПИ. Для установившихся непрерывных режимов облучения интенсивность в точке не зависит от времени и допустимо использование стационарного подхода.

С прикладной точки зрения особый интерес представляют методы восстановления физических свойств среды по результатам измерения

характеристик рассеянного излучения. Для нахождения решения обратных задач ТПИ использовались результаты решения прямой задачи. В общем случае постановка обратной задачи ТПИ является некорректной, поскольку приводит к целому классу решений. Однако для узкопоставленной задачи со всеми известными, кроме одной переменной, можно свести задачу к однопараметрической. Тогда метод решения обратной задачи, например определения объемной концентрации подвижных рассеивателей С в среде по результатам измерения соотношения интенсивности доплеровской фракции рассеянного излучения /д к полной интенсивности I излучения достигшего детектора, можно представить в виде последовательности шагов:

1. Определение основных характеристик объекта согласно (1);

2. Численный поиск зависимости соотношения /д/Г от концентрации подвижных рассеивателей С методом решения последовательности прямых задач ТПИ в виде V/ =ДС) в заданном интервале С;

3. Получение обратной зависимости С =/ ''(1[/1) в аналитической или численной форме на основе аппроксимации функции/;

4. Использование / для определения концентрации подвижных рассеивателей С по экспериментально измеряемому значению соотношения 1р/1.

Во второй главе на основе разработанного алгоритма для решения прямых задач стационарной и нестационарной ТПИ реализован проблемно-ориентированный комплекс программ ¡\fcml-doppler.

В отличие от известных программ топ! и gcmcm¡, предназначенных для расчета распределения интенсивности рассеянного излучения, предложенный комплекс позволяет: учитывать подвижность рассеивателей, рассчитывать доплеровское уширение спектра рассеянного излучения, использовать генераторы случайных чисел с большим периодом (~109 чисел), распараллеливать процесс численного моделирования с применением многопоточных центральных и графических процессоров, управлять временем решения на основе выбора технических средств вычислительного комплекса.

Обобщённая архитектура программного комплекса представлена на рис. 3.

В главе 3 для иллюстрации возможностей разработанного программного комплекса рассмотрен ряд моделей рассеивающих сред.

3.1. В целях верификации программного комплекса решена прямая задача стационарной ТПИ по расчету коэффициентов отражения и пропускания однослойной рассеивающей среды со следующими параметрами: //в=10см~', ^5=90 см"1, ¿=0.75, и=1, толщина среды с/=0.02 см, показатели преломления для окружающей среды и„,=1. Сравнение полученных результатов с данными более ранних исследований представлено в табл. 1. В столбцах среднее и Ошибка табл. 1 указаны среднее значение и стандартная ошибка для коэффициента отражения соответственно, в столбцах Та, среднее и Ошибка Та для коэффициента пропускания соответственно.

Из табл. 1 видно, что результаты, полученные с использованием комплекса Мст1-</орр1ег, незначительно отличаются от результатов, полученных другими исследователями, и их можно считать достоверными.

Рис. 3. Схема программного комплекса МстЫорр1ег

Таблица 1

Результаты верификации комплекса Мст1-Лорр1ег

Источник Метод расчета средн. Ошибка Я(| Та, средн. Ошибка Т()

деХюльст, 1980 Теория Ми 0.09739 0.66096

Пралл, 1989 Монте-Карло 0.09711 0.00033 0.66159 0.00049

Мст1-йорр1ег Монте-Карло 0.09750 0.00163 0.66091 0.00041

3.2. Решена задача по моделированию распространения непрерывного излучения в двухслойной среде с расположенным во втором слое сосудом цилиндрической формы с подвижными рассеивателями. Цилиндр диаметром 0.06-0.10 см располагался на глубине 0.05-0.25 см параллельно внешней поверхности среды, ось пучка падающего излучения проходила через центр цилиндра. Значения параметров слоев среды выбраны близкими к значениям характерным для нормальной кожи человека на длине волны 633 нм.: //„(=0.8 см" ^/==300 см"1, я,=0.8, и/=1.4, с/,=0.02 см, С,=0.0, ца2=\.0 см"1, ^=200 см"1, #,=0.9, /72=1-3, с/г= 1.0 см. Параметры для подвижных частиц выбраны близкими к значениям обычным для крови человека: ца=2.0см"', ц5=200 см"1, §=0.95, п=1.3. В рамках решения прямой задачи стационарной ТПИ получено, что контур равной интенсивности доплеровской компоненты диффузно отраженного излучения имеет форму эллипса, параметры которого описываются значениями полуосей (а -большая полуось, Ь - малая полуось). Коэффициент сжатия эллипса

(эллиптичность) определяется как: К = —. Для всех рассмотренных значений

а

диаметра с увеличением глубины залегания цилиндра эллиптичность контура равной интенсивности снижается в зависимости, близкой к линейной. С точки зрения постановки обратной задачи ТПИ, измерение Кттт может быть положено в основу методики томографии кровеносных сосудов.

3.3. Численно исследованы возможности диагностики состояния кровеносной системы человека на примере семислойной модели кожного покрова, все слои которой (кроме двух верхних, моделирующих роговой слой и эпидермис) содержали сеть микроциркуляции крови. Значения параметров слоев среды выбраны близкими к значениям, характерным для нормальной кожи человека на длине волны 633 нм.

Интегральным параметром, отражающим локальное состояние кровеносной системы, будет являться, например, объемное кровенаполнение микроциркуляторного русла Укр, определяемое как отношение объема крови к общему объему локальной области среды. При разработке диагностических комплексов Укр обычно измеряют по соотношению энергии излучения поглощенного в кровеносной системе (гемоглобином) к обшей энергии поглощенной средой.

Нормальные значения объемной концентрации крови в каждом отдельном слое кожного покрова составляют от 4% для верхней дермы с капиллярными петлями до 14% в нижнем слое дермы с глубинными сплетениями сосудов кровеносной системы. При таких значениях концентрации крови можно говорить о нормальном состоянии объемного кровенаполнения микроциркуляторного русла Укр. При изменениях состояния кровенаполнения, соответственно будут изменяться и значения объемной концентрации крови в каждом слое.

Задача состоит в том, чтобы связать измеряемые характеристики рассеянного излучения со значением Укр для решения задач диагностики кровеносной системы.

Использовав нормальные значения объемной концентрации крови в каждом слое как базовые (Укр= 1.0), численно исследовано влияние объемного кровенаполнения микроциркуляторного русла кожи человека Укр на соотношение интенсивности доплеровской компоненты к полной интенсивности излучения достигшего детектора в диапазоне Укр от 0.5 до 1.5.

На рис. 4 приведено построение изолиний соотношения интенсивности доплеровской компоненты к полной интенсивности излучения, попавшего на детектор, при различных значениях Укр и расстояния между источником и детектором излучения.

С использованием методики решения обратной задачи, описанной в гл. 1, и полученного решения прямой задачи ТПИ найдена зависимость объемного кровенаполнения микроциркуляторного русла Укр от измеряемого соотношения интенсивностей доплеровской и полной интенсивности рассеянного излучения: Укр = 104,61Ц10 /1)2 - 21,3(/д / /) +1,513, (23)

где 10 и / — интенсивность доплеровской компоненты и полная интенсивность рассеянного излучения, соответственно. Результат решения обратной задачи для одного значения расстояния источник-детектор излучения приведены на рис. 5. Полученный результат может быть полезен при разработке методов диагностики биологических тканей.

0.01 С,!]; 0.0.1 0,04 0.05 0.06 РйССгошдм псточипк-деиго]). си

Рис. 4. Изолинии соотношения интенсивности доплеровской компоненты /в к суммарной интенсивности / рассеянного излучения при различных Укр и расстояния источник-детектор

/

/

0.» Ф.И . 0.17 0.1« «,10 0,70

4х"

Рис.5. Зависимость Угр от измеряемого отношения интенсивности доплеровской компоненты к суммарной интенсивности рассеянного излучения

3.4. Возможности разработанного комплекса по решению прямых задач нестационарной ТЛИ изучены на примере задачи о прохождении короткого импульса излучения в многослойной среде, содержащей подвижные рассеиватели. В качестве модели выбрана двухслойная среда, со статическим верхним слоем и вторым слоем с сетью микроциркуляции. Моделировалось распространение в среде прямоугольных импульсов монохроматического излучения, длительностью от 0.1 до 10 пс. На рис. 6 приведен пример формы импульса рассеянного излучения и его доплеровской составляющей, сформированной излучением рассеянным подвижными частицами. Наблюдается существенное "расплывание" импульса во времени, обусловленное многократным рассеянием в среде. Доплеровская компонента импульса запаздывает относительно общего сигнала на детекторе излучения и имеет меньшую интенсивность. При анализе заднего фронта импульса установлено, что интенсивность доплеровской компоненты вплотную приближается к значению полной интенсивности рассеянного излучения. На рис. 6 этот момент соответствует задержке в 5-6 пс для импульса длительностью 5 пс. Такое поведение сигнала свидетельствует о том, что, при заданной концентрации подвижных рассеивателей, практически все фотоны, достигающие детектора излучения, успевают приобрести доплеровский сдвиг частоты вследствие одного или более рассеяния на подвижной частице.

Общий смгиап Аоплсроккм «оняомит*

3 А & 4

■ Момте-Клрло иодешфоинив - Аллроксимацич лолнидмсм 2 степени

Рис. 6. Импульс длительностью 5пс и его доплеровская компонента

0.00 О.М 0,08 0.12 0.16 0.20 Толщина статического слоя, см

Рис. 7. Влияние толщины статического слоя на задержку доплеровской компоненты импульса

3.4а) Задержка между появлением несмещенной и доплеровской компонентами сигнала возникает вследствие прохождения излучением толщи

первого статического слоя среды без подвижных рассеивателей. В результате исследования для модели двухслойной среды выявлена квадратичная зависимость временной задержки доплеровской компоненты рассеянного излучения от толщины статического слоя среды (рис. 7).

3.4.6) Верхние слои кожного покрова содержат сеть микроциркуляции биожидкости, в то время как нижележащие слои содержат крупные сосуды с более скоростными частицами.

Тогда более реалистичиая математическая модель биоткани будет представлена тремя слоями: верхний слой статический толщиной 0,02 см, средний с 10% содержанием частиц со скоростью 0,1 мм/с, и нижний толщиной 0,5 см с 10% содержанием частиц со скоростью 1,0 мм/с. На рис. 8 представлены результаты

моделирования для трех случаев: в первом рассматривалась 2-слойная модель со вторым слоем толщиной 5мм и абсолютным значением скорости

рассеивателей 0,1 мм/с, во 2 и 3 случаях поведение доплеровского уширения спектра трехслойная среда, для которой рассеянного сигнала

толщина 2-го слоя (со скоростью рассеивателей 0,1 мм/с) составляла 0,09 и 0,11 см соответственно, толщина 3-го слоя со скоростью рассеивателей 1,0 мм/с 5 мм. При исследовании зависимости среднего доплеровского сдвига частоты рассеянного излучения от времени наблюдения обнаружен излом данной зависимости (рис. 8), связанный с проникновением излучения в слой с высокоскоростными рассеивателями. Положение излома (обозначено окружностями на рис. 8) определяется толщиной второго слоя с низкоскоростными частицами. Область Г на рис. 8 соответствует толщине второго слоя 0,09 см, область II - 0,11см. В диапазоне толщин второго слоя от 0,08 до 0,15 см, (соответствующем поверхностным слоям кожного покрова) зависимость положения излома от толщины второго слоя носит характер, близкий к линейному.

Полученные результаты могут быть полезны при разработке новых импульсных методов диагностики и томографии биологических сред.

В главе 4 проведено исследование производительности многопроцессорных вычислительных систем в задачах моделирования переноса излучения методом Монте-Карло. Рассмотрены как многоядерные/многопоточные центральные процессоры общего назначения - CPU, так и многоядерные графические процессоры (GPU - Graphics Processing Units).

Алгоритмы реализации метода Монте-Карло позволяют распараллелить процесс поиска решения прямой задачи ТПИ на этапе моделирования случайного блуждания пакетов фотонов в среде. При проведении статистического моделирования с использованием многоядерных CPU общего назначения сокращение времени вычислений составило до 40-45% на каждое дополнительное процессорное ядро и до 20% на дополнительный поток по технологии Intel Hyper-Threading. Современные графические процессоры Nvidia позволяют проводить массивно-параллельные математические вычисления для решения неграфических задач с использованием интерфейса прикладного программирования CUDA

17

А*

2 ОПОЯ

3 слоя 3 слоя 2-й слои 0,1 1см

3 слоя 2-и слои O.OScu ^ / Л-у"^*.

0 5 10 15 20 2S 30 Задержка, пс

Рис. 8. Влияние строения среды на

(Compute Unified Device Architecture). Для оценки эффективности использования технологии CUDA проведена серия численных экспериментов по моделированию распространения излучения в трехслойной модельной среде, основные оптические параметры которой близки к характеристикам биологических сред. В задаче использовались программные комплексы Mcml-Doppler, MCML (http://omlc.ogi.edu/software/mc/) а также GCMCML (http://bmp.hust.edu.cny GPU_Cluster/GPU_Cluster_MCML.HTM). Количество моделируемых траекторий фотонов составило 3-107. Время моделирования заданного количества траекторий с использованием CPU Intel Core i3 ЗЗОМ (2133 МГц) составила 1700с. Характеристики использованных GPU и результаты измерения времени моделирования представлены в табл. 2. Переход к использованию GPU сократил время моделирования на 2 порядка (в -140 раз) по сравнению с CPU общего назначения, для поставленной задачи. Учитывая сильную взаимосвязь между производительностью и двумя параметрами - шириной шины данных и числом

потоковых процессоров GPU, предложен оценочный параметр производительности:

р = (24)

где Np - число потоковых процессоров, Qd - ширина шины данных в битах. Полученная зависимость времени моделирования от характеристического параметра р аппроксимируется экспоненциальной функцией:

Г = 13,71 + 146,43е-°'0394 р. (25)

Присутствие в параметре р двух характеристик, влияющих на производительность GPU, позволяет управлять временем моделирования, подбирая устройства с различным числом потоковых процессоров и шириной шины данных (24), обеспечивающих либо значение параметра р по требуемому времени моделирования (25), либо сохраняя значение р как инвариант.

Таблица 2

Характеристики использованных GPU и измеренное время моделирования тестовой _задачи_

Тип GPU Потоковых процессоров, шт. Шина данных, бит. Параметр р Время расчета, сек.

1 Nvidia GT210 16 64 1024 63,28

2 Nvidia GT9500 32 128 4096 28,05

3 Nvidia GT240M 48 128 6144 22,11

4 Nvidia GT240 96 128 12288 20,86

5 Nvidia GT9800 112 256 28672 12,12

При реализации методов Монте-Карло для моделирования рассеяния излучения с применением CPU общего назначения, как правило, используется представление чисел с плавающей точкой с двойной точностью, в соответствии со стандартом IEEE 754. Архитектура современных GPU позволяет эффективно выполнять обработку только чисел с одинарной точностью. Исследование сходимости результатов моделирования с использованием GPU показало, что для достижения дисперсии вычисляемых величин (например, коэффициентов отражения Rd и пропускания Td) на уровне результатов получаемых при моделировании на CPU с двойной точностью, необходимо увеличение числа экспериментов в методе Монте-Карло приблизительно в 4 раза.

В целом при использовании современных GPU время моделирования методами Монте-Карло сокращается многократно.

Приложения 1-4 содержат: техническое описание комплекса Mcntl-doppler, программный код для моделирования переноса излучения (при стационарном и нестационарном режимах облучения), свидетельство о регистрации программного комплекса, перечень из 20 публикаций автора по теме работы.

ЗАКЛЮЧЕНИЕ

1. На основе метода Монте-Карло разработаны модели рассеивающих сред и алгоритмы решения прямых задач стационарной и нестационарной ТПИ для сложных, оптически неоднородных объектов.

2. Учтено уширение спектра рассеянного излучения, вследствие эффекта Доплера при рассеянии подвижными частицами. Разработанные модели сред позволяют рассматривать как крупные уединенные сосуды, так и распределенную сеть микроциркуляции биологических жидкостей.

3. Построенная математическая модель позволяет распараллелить процесс вычислений и при использовании графических процессоров, с поддержкой технологии Nvidia CUDA, многократно сократить время моделирования методами Монте-Карло.

4. Разработан комплекс программ Монте-Карло моделирования для использования на компьютерах с ОС Microsoft Windows и Linux.

5. Увеличение объемной концентрации подвижных рассеивателей в многослойной среде с микроциркуляцией приводит к увеличению соотношения интенсивности доплеровской фракции рассеянного излучения и значения полной интенсивности рассеянного излучения. Глубина залегания протяженного сосуда с жидкостью в среде определяет форму распределения интенсивности доплеровской фракции рассеянного излучения, а именно - эллиптичность контура равной интенсивности доплеровской фракции рассеянного излучения убывает с увеличением глубины залегания сосуда. При распространении короткого импульса излучения в многослойной среде со статичным поверхностным слоем и с подвижными рассеивателями в нижележащих слоях длительность задержки между поступлением смещенной по частоте и несмещенной компоненты рассеянного излучения возрастает с увеличением толщины поверхностного слоя.

6. К числу перспективных направлений будущих исследований относится развитие использования технологии CUDA, и реализация эффективных генераторов случайных чисел и методов квази-Монте-Карло.

СПИСОК ОСНОВНЫХ ПУБЛИКАЦИЙ

1. Публикации в изданиях, рекомендованных ВАК РФ:

1. Старухин П.Ю. Моделирование доплеровского уширения спектра рассеянного лазерного излучения и диагностика кровотока в биологических тканях / П.Ю. Старухин, Ю.В. Клинаев // Вестник Саратовского государственного технического университета. 2010. № 4 (49). С. 28-35.

2. Старухин П.Ю. Математическое моделирование уширения спектров рассеянного лазерного излучения при распространении в биоткани методом Монте-Карло / П.Ю. Старухин, Ю.В. Клинаев, Д.В. Терин // Вестник Саратовского государственного технического университета. 2010. №3(48). С. 24-32.

3. Старухин П.Ю. Эффективность использования параллельных вычислений по технологии CUDA в решении задач биомедицинской физики / П.Ю. Старухин, Ю.В. Клинаев // Вестник Саратовского государственного технического университета. 2010. №3(48). С. 126-130.

4. Старухин П.Ю. Применение метода Монте-Карло для моделирования прохождения сверхкоротких лазерных импульсов через неоднородную среду с подвижными рассеивателями / П.Ю. Старухин, Ю.В. Клинаев // Журнал прикладной спектроскопии. 2011. Т. 78. №2.С. 277-281.

5. Starukhin P. Blood-flow measurements with a small number of scattering events / P. Starukhin', S. Ulyanov, E. Galanzha, V. Tuchin // Applied Optics. 2000. Vol.39. № 16. P.2823-2830.

2. Другие публикации:

6. Старухин П.Ю. Расчет доплеровского сдвига методом Монте-Карло / П.Ю. Старухин // Проблемы оптической физики - 1997 : материалы молодежной научной школы по оптике, лазерной физике и оптоэлектронике, Саратов, ноябрь 1997 г. Саратов: Изд-во Сарат. ун-та, 1997. С. 170-171. ISBN 5-292-01579-2

7. Старухин П.Ю. Исследование влияния геометрии детекторов излучения в задачах лазерного доплеровского измерения скорости потока крови в биоткани методом математического моделирования Монте-Карло / П.Ю. Старухин, Ю.В. Клинаев // Актуальные проблемы электронного приборостроения АПЭП-2010 : материалы Междунар. науч.-техн. конф., Саратов, 22-23 сентября 2010 г. Саратов, 2010. С.163-169.

8. Старухин П.Ю. Моделирование применения эффекта Доплера для диагностики параметров микроциркуляции крови в биотканях / П.Ю. Старухин, Ю.В. Клинаев // Методы компьютерной диагностики в биологии и медицине - 2010 : материалы ежегодной Всерос. школы-семинара, Саратов, 9-11 ноября 2010 г. Саратов: Изд-во Сарат. ун-та, 2010. С. 27-35.

9. Старухин П.Ю. О сходимости результатов моделирования методами Монте-Карло при вычислениях по технологии CUDA / П.Ю. Старухин, Ю.В. Клинаев, Д.В. Терин // Инновации и актуальные проблемы техники и технологий.: материалы Всероссийской научно-практической конференции молодых ученых. Саратов, 26-29 октября 2010 г. Саратов: Сарат. гос. техн. ун-т, 2010. Т. 1. С.19-22.

10. Старухин П.Ю. Эффективность использования параллельных вычислений по технологии CUDA в задачах моделирования методами Монте-Карло / П.Ю. Старухин, Ю.В. Клинаев, Д.В. Терин // Информационные технологии в образовании. XX Международная конференция-выставка: сб. трудов. Ч.Ш. М.: МИФИ, 2010. С.68-70. ISBN 978-5-7262-1363-7

11. Старухин П.Ю., Клинаев Ю.В. // Программа моделирования процесса распространения лазерного излучения в мутных средах Mcml-doppler // ХРОНИКИ ОБЪЕДИНЕННОГО ФОНДА ЭЛЕКТРОННЫХ РЕСУРСОВ «НАУКА И ОБРАЗОВАНИЕ». №12 (19), декабрь 2010. С. 1.

12. Старухин П.Ю. Моделирование методом Монте-Карло задач теории переноса излучения через биологические многослойные среды с подвижными рассеивателями / П.Ю Старухин., Ю.В. Клинаев, O.A. Монахова // Исследования в области естественных наук и методики их преподавания: сб. науч. статей. Саратов : ООО Изд. Центр «Наука», 2011. С. 102-1 И. ISBN 978-5-9999-0773-8

3. Патентные документы:

13. Свидетельство о регистрации электронного ресурса № 16458 ИНИМ РАО РФ. Объединённый фонд электронных ресурсов «Наука и образование»: Программа моделирования процесса распространения лазерного излучения в мутных средах Mcml-doppler / П.Ю. Старухин, Ю.В. Клинаев. Инв. № 50201050204, дата регистрации 02.12.2010.

Подписано в печать 10.05.11 Формат 60x84 1/16

Бум. офсет. Усл. печ. л. 1,16 (1,25) Уч.-изд. л. 1.0

Тираж 100 экз. Заказ 85 Бесплатно Саратовский государственный технический университет

410054, Саратов, Политехническая ул., 77 Отпечатано в Издательстве СГТУ. 410054, Саратов, Политехническая ул., 77 Тел.: 24-95-70; 99-87-39, e-mail: izdat@sstu.ru

Обозначения и сокращения.

Введение.

ГЛАВА 1. МАТЕМАТИЧЕСКИЕ МЕТОДЫ РЕШЕНИЯ ЗАДАЧ ТЕОРИИ

ПЕРЕНОСА ИЗЛУЧЕНИЯ В БИОЛОГИЧЕСКИХ СРЕДАХ

1 Л. Кожный покров человека, как пример сильно рассеивающей и поглощающей электромагнитное излучение, среды.

1.2. Аналитический обзор исследований биологических сред методами математического моделирования.

1.3. Математическая постановка задачи.

1.4. Использование теории переноса для описания процесса переноса излучения в непрозрачных средах.

1.5. Аналитические методы решения интегро-дифференциального уравнения теории переноса излучения.

1.5.1. Теория Кубелки-Мунка.

1.5.2. Метод сферических гармоник.

1.5.3. Диффузионное приближение.

1.6. Использование методов статистического моделирования Монте-Карло для решения задач теории переноса излучения.

1.6.1. Метод статистического моделирования Монте-Карло.

1.6.2. Применение метода Монте-Карло для вычисления интегралов.

1.6.3. Генераторы случайных чисел.

1.6.4. Алгоритм переноса излучения в рассеивающих средах.

1.6.5 Учет доплеровского уширения спектра излучения при рассеянии подвижными рассеивателями.

1.7. Выводы по главе 1.

ГЛАВА 2. ПРОГРАММНЫЙ КОМПЛЕКС MCML-DOPPLER.

2.1. Описание функциональных возможностей комплекса Mcml-doppler.

2.2. Схема программного комплекса.

2.3. Верификация построенной модели и разработанного алгоритма.

2.4. Геометрические аспекты некоторых прямых задач стационарной ТПИ.

2.4.1. Расстояние между осями источника и приемника излучения.

2.4.2. Угловая апертура детектора.

2.4.3. Влияние угла падения излучения.

2.4.4. Гистограмма частот диффузного отражения.

2.5. Выводы по главе 2.

ГЛАВА 3. МОДЕЛИРОВАНИЕ УШИРЕНИЯ СПЕКТРОВ РАССЕЯННОГО ИЗЛУЧЕНИЯ ДЛЯ ПОСТРОЕНИЯ МЕТОДОВ ТОМОГРАФИИ И

ДИАГНОСТИКИ БИОЛОГИЧЕСКИХ ОБЪЕКТОВ.

3.1. Решение некоторых задач стационарной теории переноса излучения методом Монте-Карло.

3.1.1. Изучение возможностей томографии кровеносного сосуда на основе анализа пространственного распределения доплеровской фракции рассеянного излучения.

3.1.2. Исследование объемного кровенаполнения микроциркуляторного русла кожи человека.

3.2. Применение метода статистического моделирования к решению задач нестационарной ТПИ.

3.2.1. Моделирование прохождения коротких импульсов излучения в среде с подвижными рассеивателями.

3.2.2. Изучение возможностей томографии многослойных биологических тканей на основе анализа уширения спектра импульса рассеянного излучения

3.3. Выводы по главе 3.

ГЛАВА 4. ПРИМЕНЕНИЕ МНОГОПРОЦЕССОРНЫХ ВЫЧИСЛИТЕЛЬНЫХ СИСТЕМ В МОДЕЛИРОВАНИИ МЕТОДОМ МОНТЕ-КАРЛО.

4.1. Использование многоядерных CPU, технология HyperThreading.

4.2. Современные графические процессоры. Технология CUDA.

4.3. Особенности использования GPU с одинарной точностью представления чисел.

4.4. Выводы по главе 4.

Актуальность исследований. В настоящее время, одним из актуальных направлений оптики является разработка методов диагностики биологических объектов. Так как физиологические и структурные изменения биологических объектов ведут к изменению оптических свойств рассеивающей среды, то изменение характеристик рассеянного излучения будет свидетельством изменения состояния биологического объекта, что отражено в работах В.В. Тучина, A.B., Приезжева и др.[1-5].

В силу присутствия в живых тканях подвижных частиц биологических жидкостей, при решении значительного круга биомедицинских задач широкое применение находят доплеровские методы исследований, позволяющие диагностировать состояние кровеносной системы in vivo [6-10].

При этом важное значение приобретает математическое моделирование физических процессов, протекающих в объекте исследования, в основе которого лежит интегро-дифференциальное уравнение теории переноса излучения (ТПИ), предложенной А. Шустером и развитой С. Чандрасекаром, де Хюлстом и А. Исимару [11-14].

Для поиска решений уравнения ТПИ в задачах с граничными условиями, соответствующими реальным биологическим объектам, приходится применять численные методы. Одним из часто используемых методов численного решения задач ТПИ для рассеивающих и поглощающих сред является метод статистического моделирования Монте-Карло, примененный, например С. Жаком, И.В. Ярославским[15-18].

Несмотря на значительные достижения исследователей (решен широкий класс задач стационарной ТПИ для расчетов распределения поглощенной энергии и исследованы возможности доплеровских методов для измерения скоростей потоков жидкости в биологических тканях), ряд важных, как в теоретическом, так и в практическом плане, задач рассмотрен недостаточно.

Не исследовано влияние характеристик подвижных рассеивателей на распределение рассеянного излучения, а также не изучены возможности применения доплеровских методов для задач диагностики и томографии биологических объектов.

Учёт и анализ этих эффектов позволил бы, развить использование доплеровских методов в целях медицинской томографии. При этом возникает необходимость проведения значительных объемов вычислений.

В связи с появлением эффективных и общедоступных реализаций многопоточных вычислительных систем, этот актуальный вопрос может стать разрешимым, если поставить задачу о параллелизации процесса моделирования методом Монте-Карло.

Перечисленные факты позволили сформулировать цель настоящей работы.

Целью диссертационной работы является разработка математических моделей и высокопроизводительных численных алгоритмов исследования процессов распространения лазерного излучения в сильно рассеивающих и поглощающих многослойных биологических средах, содержащих подвижные рассеиватели, а также разработка моделей доплеровской томографии биологических объектов и методик диагностики ряда физиологических параметров живых тканей.

Для достижения поставленной цели сформулированы следующие задачи: • разработать математическую модель распространения лазерного излучения в многослойных, сильно рассеивающих и поглощающих средах, учитывающую возможное уширение спектра рассеянного излучения вследствие эффекта Доплера, обусловленного рассеянием излучения подвижными частицами среды.

• разработать численные алгоритмы решения стационарного и нестационарного уравнений переноса излучения в мутных средах на основе метода статистического моделирования Монте-Карло, позволяющие учитывать чрезвычайно большие объёмы (~Ю10) рассеиваемых пакетов фотонов.

• реализовать разработанные алгоритмы моделирования в виде комплекса программ для ЭВМ под управлением операционных систем Microsoft и Linux.

• на основе численного эксперимента исследовать распределение доплеровской фракции излучения, рассеянного многослойной непрозрачной средой с цилиндрическим сосудом.

• изучить влияние концентрации подвижных рассеивателей на распределение доплеровской фракции излучения, рассеянного многослойной непрозрачной средой с равномерным распределением подвижных рассеивателей.

• исследовать динамику деформации короткого импульса излучения обратно рассеянного, рассеянного непрозрачной средой с цилиндрическим сосудом.

• исследовать распределение доплеровской фракции излучения, рассеянного непрозрачной средой с цилиндрическим сосудом.

• Применение комплекса программ для моделирования новых методик доплеровской томографии биологических объектов.

• Анализ возможности применения современных многопоточных центральных процессоров и многопоточных графических процессоров для производства параллельных вычислений по методу Монте-Карло.

Научная новизна работы

Построена математическая модель рассеивающей и поглощающей многослойной среды, отличающаяся учётом оптической анизотропии, присутствия подвижных рассеивателей в форме сети микроциркуляции и уединенного крупного сосуда. Разработана модель учёта уширения спектра излучения рассеянного подвижными частицами жидкости.

• Разработаны численные алгоритмы реализующие процесс моделирования переноса излучения на основе метода Монте-Карло.

• Предложены эффективные подходы к реализации моделирования переноса излучения с использованием многопоточных графических процессоров, позволившие сократить время моделирования на 2 порядка.

• Разработан проблемно-ориентированный комплекс программ для моделирования распространения излучения в рассеивающих средах.

• Проведен ряд вычислительных экспериментов в результате которых установлены следующие функциональные зависимости:

- при непрерывном характере облучения соотношение интенсивности доплеровской фракции излучения к полной интенсивности рассеянного излучения квадратично растет с увеличением объемной концентрации подвижных рассеивателей в поверхностных слоях среды.

- эллиптичность контура равной интенсивности доплеровской фракции рассеянного излучения линейно убывает с увеличением глубины залегания уединенного сосуда с жидкостью в поверхностных слоях биоткани.

- при импульсном режиме облучения двухслойной среды со статичным верхним слоем и вторым слоем содержащем подвижные рассеиватели в форме сети микроциркуляции, временная задержка между появлением на детекторе доплеровской и несмещенной компонент рассеянного излучения линейно возрастает с увеличением толщины верхнего статичного слоя. Научная и практическая значимость

Научная значимость работы заключается в применимости разработанного на основе методов Монте-Карло проблемно-ориентированного комплекса программ моделирования распространения излучения в мутных средах для исследования распределения интенсивности рассеянного излучения и уширения спектра рассеянного излучения, вследствие эффекта Доплера.

Практическая значимость заключается в полученных результатах моделирования и установленных закономерностях поведения доплеровской компоненты рассеянного излучения при различных начальных и граничных условиях. Что позволило предложить новые методики диагностики и томографии биологических объектов.

Достоверность предоставленных научных результатов подтверждается тем, что полученные результаты отличаются непротиворечивостью и находятся в хорошем соответствии с результатами исследований других авторов и экспериментальными данными, опубликованными в научной литературе.

Положения и результаты, выносимые на защиту:

1. Предложенная математическая модель процесса распространения излучения в неоднородных средах учитывает многослойное строение биоткани, оптическую анизотропию среды, явления многократного рассеяния, поглощения, отражения и преломления излучения, наличие подвижных рассеивателей в форме сети микроциркуляции и в форме уединенного крупного сосуда.

2. Программный комплекс «Мcml—doppler» позволяет решать задачи стационарной и нестационарной ТПИ, рассчитывать пространственно-временное распределение интенсивности рассеянного излучения и характеристики доплеровского уширения спектра излучения рассеянного подвижными частицами.

3. Увеличение объемной концентрации подвижных рассеивателей в многослойной среде с микроциркуляцией приводит к увеличению соотношения интенсивности доплеровской фракции рассеянного излучения и значения полной интенсивности рассеянного излучения.

4. Глубина залегания протяженного сосуда с жидкостью в среде определяет форму распределения интенсивности доплеровской фракции рассеянного излучения, а именно — эллиптичность контура равной интенсивности доплеровской фракции рассеянного излучения убывает с увеличением глубины залегания сосуда.

5. При распространении короткого импульса излучения в многослойной среде со статичным поверхностным слоем и с подвижными рассеивателями в нижележащих слоях длительность задержки между поступлением смещенной по частоте и несмещенной компоненты рассеянного излучения возрастает с увеличением толщины поверхностного слоя.

6. Разработанный алгоритм моделирования переноса излучения позволяет распараллелить вычислительный процесс в требуемом масштабе для эффективного использования многоядерных/многопоточных CPU, и многократно сократить время вычислительного эксперимента при использовании современных графических процессоров, либо увеличить точность вычислений при заданном времени, затрачиваемом на процесс моделирования.

Диссертация включает в себя введение, четыре главы, заключение, список использованных библиографических источников и четыре приложения. По каждой главе сделаны выводы.

4.4. Выводы по главе 4

Метод статистического моделирования Монте-Карло позволяет распараллелить процесс моделирования процессов переноса излучения для эффективного использования многоядерных процессоров общего назначения.

Технология Intel HyperThreading, позволяющая одновременно выполнять два потока вычислений на одном процессорном ядре позволяет достичь повышения производительности моделирования до 80%. Применение многопоточных графических процессоров с поддержкой интерфейса прикладного программирования CUDA позволяет сократить время вычислений на два порядка, несмотря на использование одинарной точности представления чисел с плавающей запятой.

ЗАКЛЮЧЕНИЕ

1. Разработаны модели рассеивающих сред и алгоритмы решения основных уравнений стационарной и нестационарной ТПИ для сложных оптически неоднородных объектов.

2. Построенный алгоритм моделирования распространения излучения учитывает основные физические характеристики исследуемых объектов (многослойность биологических тканей, геометрические параметры источников и детекторов излучения) и происходящие в среде физические процессы (рассеяние, поглощение, преломление, отражение излучения).

3. Учтено уширение спектра рассеянного излучения, вследствие эффекта Доплера при рассеянии подвижными частицами. Разработанные модели сред позволяют рассматривать как крупные уединенные сосуды, так и распределенную сеть микроциркуляции биологических жидкостей.

4. Разработан комплекс программ Монте-Карло моделирования для использования на компьютерах с архитектурой х86, под управлением операционных систем Microsoft Windows и Linux.

5. Увеличение объемной концентрации подвижных рассеивателей в многослойной среде с микроциркуляцией приводит к увеличению соотношения интенсивности доплеровской фракции рассеянного излучения и значения полной интенсивности рассеянного излучения.

6. Глубина залегания протяженного сосуда с жидкостью в среде определяет форму распределения интенсивности доплеровской фракции рассеянного излучения, а именно — эллиптичность контура равной интенсивности доплеровской фракции рассеянного излучения убывает с увеличением глубины залегания сосуда.

7. При распространении короткого импульса излучения в многослойной среде со статичным поверхностным слоем и с подвижными рассеивателями в нижележащих слоях длительность задержки между поступлением смещенной по частоте и несмещенной компоненты рассеянного излучения возрастает с увеличением толщины поверхностного слоя.

8. Построенная математическая модель позволяет распараллелить процесс вычислений в требуемом масштабе. Применение графических процессоров, с поддержкой технологии Nvidia CUDA, позволяет многократно сократить время моделирования методами Монте-Карло.

9. К числу перспективных направлений будущих исследований относится развитие использования технологии CUDA, и реализация эффективных генераторов случайных чисел и методов квази-Монте-Карло.

В заключение хочу выразить глубокую благодарность д.ф.-м.н., заведующему кафедрой оптики и биофотоники СГУ профессору Валерию Викторовичу Тучину, д.ф.-м.н., профессору кафедры оптики и биофотоники СГУ Сергею Сергеевичу Ульянову, заложившим фундамент этой работы, д.ф.-м.н., заведующему кафедрой программного обеспечения вычислительной техники и автоматизированных систем СГТУ профессору Байбурину Вилу Бариевичу, д.ф.-м.н., заведующему кафедрой технической физики и информационных технологий ЭТИ СГТУ профессору Альберту Марковичу Кацу за научные дискуссии и моему научному руководителю, д.ф.-м.н., профессору Юрию Васильевичу Клинаеву, за постановку задачи и чуткое руководство моими исследованиями.

1. Приезжев, A.B. Лазерная диагностика в биологии и медицине / A.B. Приезжев, В.В. Тучин, Л.П. Шубочкин. М.: Наука, 1989. - 240 с.

2. Тучин, В.В. Исследование биотканей методами светорассеяния / В.В. Тучин // Успехи физических наук. — 1997. — 167(5). С. 517-539.

3. Тучин, В. В. Лазеры и волоконная оптика в биомедицинских исследованиях / В.В. Тучин. Саратов : Изд-во Сарат. ун-та, 1998. - 384 с. - ISBN 5-292-01779-5.

4. Тучин, В.В. Оптическая медицинская диагностика / В.В. Тучин // Известия Саратовского университета. — 2005. 5(1). - С. 39-53.

5. Tuchin, V.V. Selected Papers on Tissue Optics / Valery Tuchin. Bellingham: SPIE, 1994. ( Applications in Medical Diagnostics and Therapy. Vol. MS 102).

6. Riva, С. Laser Doppler measurements of blood flow in capillary tubes and retinal arteries / C. Riva, В. Ross, GB. Benedek // Invest Ophthalmol. 1972. - 11(11). - P. 936-44.

7. Model for Photon Migration in Turbid Biological Media / R. F. Bonner, R. Nossal, S. Havlin, G. H. Weiss // J. Opt. Soc. Am. A. 1987. - 4/ - P. 423-432.

8. In-vivo evaluation of signal processors for laser Doppler tissue flowmeters / H. Ahn, К. Johansson, О. Lundgren, G.E. Nilsson // Medical and Biological Engineering and Computing. 1987. - 25. - P. 207-211.

9. Oberg, P.A. Laser Doppler Flowmetry / P.A. Oberg // Critical Review. 1990. -18(2). - P. 125-162.

10. Ю.Ульянов, C.C. Динамика спеклов и эффект Доплера / С.С. Ульянов // Соросовский образовательный журнал. 2001. - 7(10). - С. 109-114.

11. Schuster, А. // Astrophys. J. 1905. - 21. - P. 1.

12. Чандрасекар, С. Перенос лучистой энергии / С. Чандрасекар. М.: Изд-во иностр. лит., 1953. - 432 с.

13. Исимару, А. Распространение и рассеяние волн в случайно-неоднородных средах: в 2 т. Т.1. Однократное рассеяние и теория переноса / А. Исимару. М.: Мир, 1981. - 281 с.

14. Хюлст, Г. Рассеяние света малыми частицами / Г. ван де Хюлст. М.: Изд-во иностр. лит., 1961. - 536 с.

15. Соболь, И.М. Численные методы Монте-Карло / И.М. Соболь. М.: Наука, 1973. -312 с.

16. Cashwell, E.D. A Practical Manual on the Monte Carlo Method for Random Walk Problems / E.D. Cashwell, C.J. Everett. New York: Pergamon Press, 1959. - 103 p.

17. Wang, L.H. Mcml Monte-Carlo Modeling of Light Transport in Multilayered Tissues / L.H. Wang, S.L. Jacques, L.Q. Zheng // Computer Methods and Programs in Biomedicine. - 1995. - 47(2). - P. 131-146.

18. New approach to Monte Carlo simulation of photon transport in the frequency-domain / Yaroslavsky I.V., A.N. Yaroslavsky, H.J. Schwarzmaier, G.G. Akchurin, V.V. Tuchin // SPIE Proc. 1995. - 2626. - P. 45-55.

19. Anderson, R.R. Optical properties of the human skin / R.R. Anderson, J.A.Parrish // J.D. Regan & J.A. Parrish (Eds.), The Science of Photomedicine. New York: Plenum, 1982. - P. 147-194.

20. Tuchin, V.V. Tissue Optics / V.V. Tuchin. -Bellingham: SPIE Press, 2000.21.3имняков Д.А. Оптическая томография тканей / Д.А. Зимняков, В.В. Тучин // Квантовая электроника. 2002. - 32(10). - С. 849-867.

21. Muller, G. Medical Optical Tomography: Functional Imaging and Monitoring / G. Muller, ed. Vol. IS11. - 1993. - Bellingham: SPIE.

22. Chance,B. Optical method / B. Chance // Arm. Rev. Biophys. Biophys. Chem. — 1991. 20. - P. 1-28.

23. Соболь, И.М. Метод Монте-Карло / И.М. Соболь // Популярные лекции по математике. — Вып.46. М.: Наука, 1985. — 64 с.

24. Бусленко, Н.П. Метод статистических испытаний (Монте-Карло) и его реализация на цифровых вычислительных машинах / Н.П.Бусленко, Ю.А.Шрейдер. -М.: Физматгиз, 1961. 226 с.

25. Войтишек, A.B. Дополнительные сведения о численном моделировании случайных элементов: учеб. пособие / A.B. Войтишек. Новосибирск: НГУ, 2007.

26. Войтишек, A.B. Функциональные оценки метода Монте-Карло: учеб. пособие / A.B. Войтишек. Новосибирск: НГУ, 2007.

27. Войтишек, A.B., Основы метода Монте-Карло в алгоритмах и задачах: в 6 ч. / A.B. Войтишек. Новосибирск: НГУ, 1997 - 2004.

28. Методы и алгоритмы статистического моделирования: сб. науч. тр. / под ред. Г.А. Михайлова. Новосибирск: ВЦ СО АН СССР, 1983.

29. Методы статистического моделирования: сб. науч. тр. / под ред. Г.А. Михайлова. Новосибирск: ВЦ СО АН СССР, 1986.

30. Теория и приложения статистического моделирования: сб. науч. ст. / под ред. Г.А.Михайлова. Новосибирск: ВЦ СО АН СССР, 1989.

31. Михайлов, Г.А. Алгоритмы метода Монте-Карло для восстановления индикатрисы рассеяния с учетом поляризации / Г.А. Михайлов, С.А. Ухтинов, A.C. Чимаева // Доклады Академии наук. 2008. - 423(2). - С. 161-164.

32. Fishman, G.S. Monte Carlo: concepts, algorithms, and applications / G.S. Fishman. New York: Springer, 1996.

33. Briesmcister, J.F. MCNP A General Monte Carlo N-Particle Transport Code, Version 4A / J.F. Briesmeister. - Los Alamos: Los Alamos National Laboratory, 1993.

34. Brown, F. A Tutorial on Using MCNP for 1-Group Transport Calculations / F. Brown, N. Barnett. Los Alamos: Los Alamos National Laboratory, 2007.

35. Метод Монте-Карло в атмосферной оптике / под общ. ред. Г.И. Марчука; Академия наук СССР. Сибирское отделение. Новосибирск: Вычислительный центр. Наука. Сиб. отд-ние, 1976.

36. Решение прямых и некоторых обратных задач атмосферной оптики методом Монте-Карло: Методика, алгоритмы, программа, результаты расчетов / под ред. Г.И. Марчука. Новосибирск: Наука, Сиб. отд-ние, 1968.

37. Ермаков, С.М. Метод Монте-Карло и смежные вопросы / С.М. Ермаков. — М.: Наука, 1975.

38. Избранные алгоритмы метода Монте-Карло / С.М. Ермаков, A.C. Расулов, М.Т. Бакоев, А.З. Веселовская. Ташкент: Изд-во Ташкент, гос. ун-та, 1992.

39. Апресян, JI.A. Теория переноса излучения. Статистические и волновые аспекты / JI.A. Апресян, Ю.А. Кравцов. М.: Наука, 1983.

40. Будак, В.П. Методы решения уравнения переноса излучения / В.П. Будак. -М.: Изд. дом МЭИ, 2007. 52 с.

41. Гермогенова, Т.А. Локальные свойства решений уравнения переноса / Т.А. Гермогенова. М.: Наука, 1986. - 272 с.

42. Соболь, И.М. Точки равномерно заполняющие многомерный куб / И.М. Соболь. -М.: Знание, 1985. 32 с.

43. Jacques S.L. Light Distributions from Point, Line and Plane Sources for Photochemical Reactions and Fluorescence in Turbid Biological Tissues / S.L. Jacques // Photochemistry and Photobiology. 1998. - 67(1). - P. 23-32.

44. Wang, L. Monte Carlo Modeling of Light Transport in Multi-layered Tissues in Standard С / L. Wang, S.L. Jacques. University of Texas M. D. Anderson Cancer, Center, 1992.

45. Jacques S.L. Time-resolved reflectance spectroscopy in turbid tissues / S.L. Jacques // IEEE Trans. Biomed. Eng. 1989. - 36. - P. 1155-1161.

46. Polarization Memory of Multiply / F.C. MacKintosh, J.X. Zhu, D.J. Pine, D.A. Weitz // Phys. Rev. Rapid Commun. 1989. - B40. - P. 9342-9345.

47. Depolarization of multiply scattering waves by spherical diffusers: Influence of size parameter / D. Bicout, C. Brosseau, A.S. Martinez, J.M. Schmitt // Phys.Rev.E. -1994. 49. - P. 1767-1770.

48. Signal attenuation and localization in optical coherence tomography studied by Monte Carlo simulation / Derek J. Smithies, T. Lindmo, Zhongping Chen, J. Stuart Nelson, Thomas E. Milner // Phys. Med. Biol. 1998. - 43. - P. 3025-3044.

49. Optical Doppler Tomography / Zhongping Chen, Y. Zhao, Shyam M. Srinivas, J. Stuart Nelson, Neal Prakash, Ron D. Frostig // IEEE journal of selected topics in quantum electronics. 1999. - 5(4). - P. 1134-1142.

50. Accuracy and noise in optical Doppler tomography studied by Monte Carlo simulation / Tore Lindmo, D.J. Smithies, Zhongping Chen, J. Stuart, Nelson Thomas, E. Milner // Phys. Med. Biol. -1998. 43. - P. 3045-3064.

51. Kienle, A. Non-invasive determination of muscle blood flow in the extremities from laser Doppler spectra / A.Kienle // Phys. Med. Biol. 2001. - 46. - P. 12311244.

52. Меглинский, И.В. Моделирование спектров отражения оптического излучения от случайно-неоднородных многослойных сильнорассеивающих и поглощающих свет сред методом Монте-Карло / И.В.Меглинский // Квантовая электроника. 2001. -31. - С. 1101-1107.

53. Jentink, H.W. Monte-Carlo Simulations of Laser Doppler Blood-Flow Measurements in Tissue / H.W. Jentink, F.F. de Mul et al. // Applied Optics. 1990. 29(16). P. 2371-2381.

54. Моделирование распространения оптического излучения в фантоме биологической ткани на суперЭВМ МВСЮОО/М / Л.П. Басс, О.В. Николаева, B.C. Кузнецов, А.В. Быков, А.В. Приезжев, А.А. Дергачев // Математическое моделирование. — 2006. — 18(1). С. 29-42.

55. Binzoni, Т. Full-field laser-Doppler imaging and its physiological significance for tissue blood perfusion / T. Binzoni, D.Van De Ville // Phys. Med. Biol. 2008. - 53. - P. 6673-6694.

56. Non-invasive laser Doppler perfusion measurements of large tissue volumes and human skeletal muscle blood RMS velocity / T. Binzoni, T.S. Leung, D. Boggett, D. Delpy // Physics in Medicine and Biology. 2003. - 48(15). - P. 2527-2549.

57. Binzoni, T. The photo-electric current in laser-Doppler flowmetry by Monte Carlo simulations / T. Binzoni, T.S. Luing, D. Van De Ville // Phys. Med. Biol. 2009. - 54. -P. 303-318.

58. Fischer, D.G. Monte Carlo modeling of spatial coherence:free-space diffraction / D.G. Fischer, S.A.Prahl, Donald D. Duncan // J. Opt. Soc. Am. A. 2008. - 25(10). -P. 2571-2581.

59. Исследование возможности повышения эффективности лазерного удаления татуировок с помощью оптического просветления кожи / Э.А. Генина, А.Н. Башкатов, В.В. Тучин и др. // Квантовая электроника. 2008. - 38(6). - С. 580-587.

60. Fang, Q. Mesh-based Monte Carlo method using fast ray-tracing in Pliicker coordinates. / Q. Fang // Biomedical optics express. 2010. - 1(1). — P. 165-175.

61. Введение в рассеяние света биологическими объектами / Н.Г. Хлебцов, И.Л. Максимова, В.В. Тучин, J1.B. Ванг // Оптическая медицинская диагностика: сб. / под ред. В.В. Тучина. М.: Физматлит, 2007. - С. 36-156.

62. Can the light scattering depolarization ratio of small particles be greater than 1/3? / N.G. Khlebtsov, A.G. Melnikov, V.A. Bogatyrev, L.A. Dykman, A.V. Alekseeva, L.A. Trachuk, B.N. Khlebtsov // J. Phys. Chem. B. 2005. - 109(28). - P. 13578-13584.

63. Chimaeva, A.S. Study of polarization estimates variance by the Monte Carlo method / A.S. Chimaeva, G.A. Mikhailov // Russian Journal of Numerical analysis and Mathematical Modelling. 2006. - 20(3). - P. 305-317.

64. Михайлов, Г.А. Дисперсия стандартной векторной оценки метода Монте-Карло в теории переноса поляризованного излучения / Г.А. Михайлов, С.А. Ухинов, А.С. Чимаева // Журн. вычисл. математики и мат. физики. 2006. - 46(11). - С. 2199-2212.

65. Niemz, M.II. Laser-Tissue Interactions: Fundamentals and Applications / M.H. Niemz. Berlin: Springer, 1996. - 305 p.

66. Patterson, M.S. Time resolved reflectance and transmittance for the noninvasive measurement of tissue optical properties / M.S. Patterson, B. Chance, В. C. Wilson // Applied optics. 1986. - 28(12). - P. 2331-2336.

67. Keijzer, M. Optical diffusion in layered media / M. Keijzer, W.M. Star, P. Storchi // Applied Optics. 1988. - 27. - P. 1820-1824.

68. The optical properties of blood in the near infrared spectral range / A.N.Yaroslavsky, Y.l.Yaroslavsky, T. Goldbach, H.J.Schwarzmaier // SPIE Proc. 1996.- 2678. P. 314-332.

69. Борн, М. Основы оптики / М. Борн, Э. Вольф. М.: Наука, 1970. - 856 с.

70. Mie, G. Beitrage zur Optik triiber Medien, speziell kolloidaler Metallosungen / G.Mie // Ann. Phys. Leipzig. 1908. - 25. - P. 377-445.

71. Henyey, L.G. Diffuse Radiation in the Galaxy / L.G.Henyey, J.L. Greenstein // Astrophys. J. 1941. - 93. - P. 70-83.

72. Star, W.M. Diffusion Theory of Light Transport. Optical-Thermal Response of Laser-Irradiated Tissue / ed. by A.J. Welch and M.J.C.van Gemert. 1995. - P. 131206.

73. Рогаткин, Д.А. Об одном подходе крешению многомерных задач теории рассеяния света в мутных средах / Д.А.Рогаткин // Квантовая электроника. — 2001. — 31(5). С. 279-281.

74. Lahaye, C.T.W. Optimal Laser Parameters for Port Wine Stain Therapy: a Theoretical Approach / C.T.W. Lahaye, M.J.C. van Gemert // Physics in Medicine and Biology. 1985. - 30(6). - P. 573-588.

75. Gemert, van M.J.C. Is There an Optimal Laser Treatment for Port Wine Stains? / M.J.C. van Gemert, A.J. Welch, P.A. Alpesh // Lasers in Surgery and Medicine. 1986. -6(1). - P. 76-83.

76. Скипетров, C.E. Диффузионно-волновая спектроскопия в средах с пространственно неоднородной динамикой рассеивателей :дис. канд. физ.-мат. наук / С.Е. Скипетров. -М.: МГУ, 1998.

77. Михайлов, Г.А. Оптимизация весовых методов Монте-Карло / Г.А. Михайлов. -М.: Наука, 1987. 240 с.

78. Минин, И.Н. Теория переноса излучения в атмосфере планет / И.Н. Минин. -М.: Наука, 1988. 264 с.

79. Light Distributions in a Port Wine Stain Model Containing Multiple Cylindrical and Curved Blood Vessels / G.W. Lucassen, W. Verkruysse, M.Keijzer, M.J.C. van Gemert // Lasers in Surgery and Medicine. 1996. - 18(4). - P. 345-357.

80. Wavelengths for Port Wine Stain Laser Treatment: Influence of Vessel Radius and Skin Anatomy / M.J.C. van Gemert, D.J.Smithies, W.Verkruysse, T.E.Milner, J.S.Nelson // Physics in Medicine and Biology. 1997. - 42(1). - P. 41-50.

81. Modeling of Temperature Distribution in the Skin Irradiated by Visible Laser Light / Y.N.Scherbakov, A.N. Yakunin, I.V. Yaroslavsky, V.V.Tuchin // Proc. SPIE. 1994.- 2082(3). P. 268-275.

82. Smithies, D.J. Modelling the Distribution of Laser Light in Port-Wine Stains with the Monte Carlo Method / D.J. Smithies, P.H.Butler // Physics in Medicine and Biology.- 1995. 40. - P. 701-733.

83. Hall, A. On an experiment determination of л / A. Hall // Messeng. Math. -1873(2).

84. Metropolis N. The Monte Carlo Method / N. Metropolis, S.Ulam // J. Amer. statistical assoc. 1949. - 44(247). - P. 335-341.

85. Владимиров, B.C. Расчёт наименьшего характеристического числа уравнения Пайерса методом Монте-Карло / В.С.Владимиров, И.М.Соболь // Вычислит, математика. 1958. — 3.

86. Кнут Д.Э. Искусство программирования. Т. 2. Получисленные алгоритмы / Д.Э. Кнут. - М.: Вильяме, 2000. - 832 с.

87. Старухин, П.Ю. Расчет доплеровского сдвига методом Монте-Карло // Проблемы оптической физики 1997: материалы молодежной научной школы по оптике, лазерной физике и оптоэлектронике. - Саратов: Изд-во Сарат. ун-та, 1997. -С. 170-171.

88. Starukhin, P. Blood-flow measurements with a small number of scattering events / P. Starukhin, U.S., E. Galanzha, V. Tuchin // Applied Optics. 2000. - 39(16): - P. 2823-2830.

89. Wilson, B.C. Thephysicsofphotodynamictherapy / B.C. Wilson, M.S. Patterson // Phys. Med. Biol. 1986. - 31. - P. 327-360.

90. De Mul F.F. Doppler Monte Carlo simulations of light scattering in tissue to support laser-Doppler perfusion measurements / F.F. De Mul, W.Steenbergen, J. Greve // Technology and Health Care. 1999. - 7(2-3). - P. 171-183.

91. Программа моделирования процесса распространения лазерного излучения в мутных средах Mcml-doppler / П.Ю. Старухин, Ю.В. Клинаев // Хроники объединенного фонда электронных ресурсов «Наука и образование», 2010. 12(19). - С. 1.

92. Monte-Carlo Simulations // Сайт Монте-Карло моделирования Oregon Medical Laser Center. 2007. - Электронный ресурс. URL: http:// omlc.ogi.edu/software/mc (дата обращения: 11.02.2011).

93. Jacques, S.L. Angular dependence of HeNe Laser light scattering by human dermis / S.L. Jacques, C.A. Alter., S.A. Prahl // Las. Life Sci. 1987. - 1. - P. 309333.

94. A Monte Carlo Model of Light Propagation in Tissue / S.A. Prahl, M. Keijzer, S.L. Jacques, and A.J. Welch // SPIE Institute Series. 1989. - 5. - P. 102-111.

95. Doppler techniques for blood microcirculation monitoring in dentistry / Pavel Y. Starukhin, N.A. Kharish, Sergey S. Ulyanov, Alexander V. Lepilin, Valery V. Tuchin // International Conference on Correlation Optics. 1997. - Proc. SPIE.

96. Monte-Carlo simulation of Doppler shift for laser light propagation in human teeth / P.Y. Starukhin, N.A. Kharish, A. Karpovitch, A.V. Lepilin, S.S. Ulyanov, V.V. Tuchin // Photon Propagation in Tissues: Proc. SPIE— 1998. —III.

97. Recovery of optical parameters in multiple-layered diffusive media: theory and experiments / J. Ripoll, J.P. Culver, Deva N. Pattanayak, Arjun G. Yodh // J. Opt. Soc. Am. 2001. - 18(4). - P. 821-830.

98. Kolinko V.G. Monte-Carlo simulations of transit time dependence of Doppler spectra in dynamic time-resolved tomography / V.G. Kolinko, A.V. Priezzhev, F.F.M. de Mul // Proc. SPIE. 1996. - 2925. - P. 160-168.

99. Intel® Hyper-Threading Technology (Intel® HT Technology) // Сайт технологии Hyper-Threading фирмы Intel. 2011. — Электронный ресурс. URL: http://vvww.intel.com/technology/platform-technology/hyper-threading/index.htm (дата обращения: 11.02.2011).

100. Боресков, A.B. Основы работы с технологией CUDA / Боресков A.B. А.В.Харламов. М.: ДМК Пресс, 2010. - 232 с.

101. Nvidia Developer Zone // Сайт для разработчиков программного обеспечения устройств фирмы Nvidia. 2011. - Электронный ресурс. URL: http://developer.nvidia.com. (дата обращения: 10.02.2011).

102. Сайт GPU-accelerated Monte Carlo Simulation for Light Propagation in Multi-Layered Tissue Geometry — 2011. — Электронный ресурс. URL: http://c0de.g00gle.c0m/p/gpumcml (дата обращения: 17.02.2011).

103. Сайт Atomfysik: GPU Monte Carlo 2009. - Электронный ресурс. URL: http://www.atomic.physics.lu.se/biophotonics/ourresearch/montecarlosimulations/gpumontecarlo (дата обращения: 03.12.2010).

104. Старухин, П.Ю. Эффективность использования параллельных вычислений по технологии CUDA в решении задач биомедицинской физики / П.Ю. Старухин, Ю.В. Клинаев // Вестник Саратовского государственного технического университета. — 2010. №3(48). - С. 126-130.