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

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

КОРОТЧЕНКО Мария Андреевна

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

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

Автореферат

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

Новосибирск - 2008

003455685

Работа выполнена в Институте вычислительной математики и математической геофизики СО РАН.

Научный руководитель: член-корреспондент РАН, профессор

Михайлов Геннадий Алексеевич

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

профессор Григорьев Юрий Николаевич

кандидат физико-математических наук, старший научный сотрудник Плотников Михаил Юрьевич

Ведущая организация: Институт теоретической и прикладной

механики СО РАН (г. Новосибирск)

Защита состоится 23 декабря 2008 года в 15 часов на заседании диссертационного совета Д 003.061.02 при Институте вычислительной математики и математической геофизики СО РАН в конференц-зале Института вычислительной математики и математической геофизики СО РАН по адресу: 630090, Новосибирск, проспект ак. Лаврентьева, 6.

С диссертацией можно ознакомиться в библиотеке Института вычислительной математики и математической геофизики СО РАН.

Автореферат разослан 21 ноября 2008 года.

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

С.Б. Сорокин

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

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

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

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

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

Основные цели работы.

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

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

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

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

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

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

Апробация работы. Результаты, изложенные в диссертации, докладывались и обсуждались на семинаре Отдела статистического моделирования в физике ИВМиМГ СО РАН (2003 - 2008 гг.), на семинаре Отдела математических задач геофизики ИВМиМГ СО РАН (2008 г.), а также на следующих конференциях:

- Конференции молодых ученых ИВМиМГ СО РАН (2004-2006, 2008 гг.)

- Межвузовская научная студенческая конференция "Интеллектуальный потенциал Сибири" (г. Новосибирск, 2004 г.)

- XLII Международная научная студенческая конференция "Студент и научно-технический прогресс" (г. Новосибирск, 2004 г.)

- Пятый международный симпозиум по моделированию Simulation (г. Санкт-Петербург, 2005 г.)

- Всероссийская конференция по вычислительной математике КВМ-2007 (г. Новосибирск, 2007 г.)

- Международная Конференция по Математическим Методам в Геофизике ММГ-2008 (г. Новосибирск, 2008 г.)

Публикации. По тематике диссертации автором опубликовано 8 работ, среди которых 5 работ в изданиях из списка ВАК [1-5], в том числе одна работа без соавторов. Список опубликованных работ приведен в конце автореферата.

Структура и объём работы. Диссертация состоит из введения, четырёх глав, разбитых на параграфы, заключения и списка литературы. Диссертация изложена на 103 страницах, включает библиографический список из 42 наименований работ, 6 рисунков, 6 таблиц.

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

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

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

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

Уравнение Больцмана, записанное с использованием условной плотности вероятности w = u>(v', v'j |v, vi) перехода пары скоростей частиц от (у', v'j) к (v, vi), имеет вид:

|/(V, t) = J {/(v', t)f(y[,t) - /(v, t)f(yut)}w dv' dvj dvb

R9

Здесь /(v, t) — "одночастичная" функция плотности распределения по скоростям V в момент времени t € (О, Т].

Уравнение Смолуховского имеет вид:

i+J=l Ol

где n(l, t) — частота (числовая плотность) частиц размера I в момент времени t, причем размер частицы I принимает натуральные значения; Ki3 — коэффициенты коагуляции (предполагаются заданными) . После присоединения к уравнениям начальных данных, численное решение сформулированных задач Коши будем понимать в смысле нахождения линейных функционалов от функции /(v,i) и от функции n(l,t), соответственно.

Два следующих раздела содержат известные факты, приведенные в статье [3]. В разделе 1.1.1 описана общая схема построения весовых статистических алгоритмов для оценки линейных функционалов Ih = (<р, h) от решения интегральных уравнений второго рода

ц> = К<р +f.

Здесь же введено понятие "функция ценности" iр*{х) (точки х по отношению к функционалу Д). Функция </>*(■) является решением сопряжённого уравнения у* — К*<р* + h и используется в весовом моделировании для уменьшения трудоёмкости алгоритмов. В разделе 1.1.2 на основе свойств функции ценности введено частичное ценностное моделирование элементарных переходов в модифицированном фазовом пространстве.

В параграфе 1.2 описаны математические модели многочастичных систем [Кац М., 1965; Bird G., 1976; Бобылев A.B., 1975; Nanbu К., 1980; Лушников A.A., 1978; Marcus А.Н., 1968; Gillespie D.T., 1972] для изучаемых кинетических уравнений: уравнения Больцмана, которое описывает релаксацию однокомпонентно-

го газа (в разделе 1.2.1) и уравнения Смолуховского, описывающего коагуляцию частиц (в разделе 1.2.2). Фазовое состояние системы частиц в момент времени t определяется вектором скоростей всех частиц V = (vi, V2,..., vjv) € для уравнения Больцмана и вектором состояний X = (N, Ьдг) = (N, h,..., In), где 1г — размер (объём) i-ой частицы, N — количество частиц в системе, для уравнения Смолуховского.

Эти модели использованы в параграфе 1.3 для построения базового интегрального уравнения второго рода с помощью модификации фазового пространства состояний системы на основе введения новой неременной — номера взаимодействующей пары ж = (г, j) [Михайлов Г.А., Рогазинский C.B., 2002]:

t

F(Z, t) = J J F(Z',t')K(Z',t' -> Z,t) dZ' dt' + F0(Z,t)

о ъ

где Z = (я-, V) или Z — (тт,Х), a плотность F(Z,t) взаимодействий пары частиц в системе связана с плотностью столкновений <p{-,t) следующим образом: <p(-,t) = ^C-FXZ,£) = -),t). Там же за-

7Г 7Г

писано факторизованное ядро этого уравнения и сделано замечание о возможности представления его решения в виде ряда Неймана, что даёт возможность обоснованно применять весовые методы Монте-Карло [Ермаков С.М., Михайлов Г.А., 1976].

В параграфе 1.4 вводятся линейные функционалы от потока частиц Jh = (H, Ф), которые могут быть представлены в виде функционалов от плотности взаимодействий пары частиц: .7// = (H,F). Эти функционалы будут оцениваться с помощью весового моделирования, описанного в разделе 1.4.1. При этом рассматривается вспомогательная цепь Маркова (Zn,tn), п £ 0,1,..., с нормированной плотностью перехода, которую удобно представить в факторизован-ном виде

P(Z',t' ->Z,t) = Px{t' - t\Y')P2(n\Y')P3(Y' Y\n), и нормированной плотностью начального состояния

P0(Z,t) = P0(Z)5(t). Случайные веса вводятся по формулам

Со — ГЭ//7 \ ! Qn — Qn-lQ(Zn-l,tn~l',

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

Основные результаты диссертации представлены в главах 2, 3 и 4.

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

В параграфе 2.1 записано уравнение Больцмана для случая макс-велловских молекул

dt

= S9 (■т)[/(v''i)/(w''i} ~ /(v'i)/(w't)]

где u = v—w, u = |u|, g (cos -ü) = u-<r(u,cos îT} — произведение дифференциального сечения рассеяния на модуль относительной скорости; п — единичный вектор направления относительной скорости частиц после столкновения.

Начальные данные в тестовых задачах выбрались такими, для которых существует аналитическое решение [Бобылев A.B., 1975]:

где

T = r(t,a0) = 1 -ве'а°ьх; t > 0; 0 < 2/5, А = 1/6.

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

В параграфе 2.2 описан алгоритм моделирования процесса эволюции Лг-частичного ансамбля частиц. Алгоритм заключается в моделировании начального распределения частиц по скоростям с помощью метода суперпозиции и последовательном моделировании времени свободного пробега системы по плотности

номера пары взаимодействующих частиц и, наконец, новых скоростей частиц. При этом предложен новый алгоритм коррелированного выбора номера взаимодействующей пары дня одновременного моделирования ансамблей разных объёмов (Лг1 > N2). Предложенный алгоритм обладает следующим, полезным для коррелирования оценок, свойством: в случае ¿1,1*2,л,.72 ^ пары 7Г1 = (¿1,7*1) и л"2 = (¿2,7*2) совпадают.

В параграфе 2.3 приводится общий вид весовых оценок искомых функционалов. Здесь также предложены новые оценки решения уравнения Больцмана и его параметрических производных, которые позволяют изучать зависимость результатов от таких параметров задачи, как параметр начального распределения в (алгоритм описан в разделе 2.3.1) и параметр временнбго распределения яо (алгоритм описан в разделе 2.3.2).

В первом случае для оценки функционала Jн{T) при в = в\ в качестве моделируемой плотности начального распределения использована плотность с в = 62: /(V, 0) = ро (у, 62). В связи с этим вводится начальный вес

Для весовой "оценки по пересечениям" (, которая в случае прямого моделирования свободного пробега имеет вид С = <5 оН(у„), где у = и{Т) — номер столкновения, непосредственно предшествующего выходу системы за рассматриваемую временную границу Т, имеем

Для оценки производной функционала ,1н{Т) по в необходимо оценить математическое ожидание следующей величины:

= Ро(УЛ) Ро(У,в2)'

ЕС = Jн{T).

дв (1-0) ^1.2 2(1 - в) + Уг20 - 30(1 - в) + 2(1 - в)2

В разделе 2.3.2 в качестве моделируемой плотности свободного пробега системы рассмотрена плотность экспоненциального распределения с параметром N/2, соответствующим значению ао = 1. Для вычисления функционала Jн{T) при а ~ а0 + со использовалась оценка С = (а)Н(Уи), где результирующий вес имеет вид:

Для оценки производной функционала Jц(T) по а в точке а — а0 используется величина

В параграфе 2.4 для исследования зависимости функционалов от решения уравнения Больцмана от числа тестовых частиц N рассматривается возможность одновременного построения траекторий /^-частичного случайного процесса при различных N (N1 > N2) с положительной корреляционной зависимостью.

В разделе 2.4.1 приведены два способа реализации такого построения с использованием "глобально-весового" метода. Рекомендовано свободный пробег моделировать для разных ансамблей соответственно одной плотности сехр{— аЬ}, причем аЫг/2 ^ <т < аМ^/ 2. В этом случае итоговые веса в момент времени Т равны

Показано, что оптимальный (в смысле минимизации дисперсий) параметр имеет вид

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

В параграфе 2.5 для оценки "хвоста" скоростного распределения в момент времени Т, то есть величины

<2{Т)(а) = а^ехр{-£0АгГ/2}.

а( + N2)

J^)(T) = V (|v(i)| > vg\t = T) = I /(v, T) dv,

va

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

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

Результаты главы 2 опубликованы в статьях [1-4,6-8].

В третьей главе рассматриваются "ценностные" весовые модификации статистического моделирования для численного решения уравнения Смолуховского в однородном по пространству случае.

В параграфе 3.1 представлены два вида коэффициентов и начальных данных уравнения Смолуховского, для которых известно аналитическое решение: при Ki3 = 1 и п(1,0) = <5;д, I > 1 имеем [Во-лощук В.М., 1984]

1 ( 0.5* У"1 ,

= (Г+ому* (ттобг)

При Ktj ~ (i+j)/2 и п(1,0) = Ju, I ^ 1 имеем [Aldous D.J., 1999]

n(l, t) = е~° 5tB (l — е-0 5i, /) I B{x,l)={lx)l~^'l\ 1^1.

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

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

В параграфе 3.2 приводятся выражения статистических весовых оценок искомых функционалов.

Далее, в параграфе 3.3 рассматривается ценностное моделирование ансамбля тестовых частиц с целью оценки искомых функционалов от уравнения Смолуховского с постоянными коэффициентами. В разделе 3.3.1 доказана теорема, утверждающая, что при определённых условиях для задачи коагуляции имеет место равенство

<р'{Х,Ь) = Ъ (Х)С(М0,1),

где N1 (X) — количество мономеров в ансамбле. В частности теорема верна для постоянных и линейных коэффициентов коагуляции. Таким образом, это уникальные случаи, когда точно известна функция ценности. Эту функция была использована для ценностного моделировании номера пары взаимодействующих частиц. С целью моделирования свободного пробега, на основе теоремы, можно строить приближения к известной функции ценности, удобные для моделирования.

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

функции ценности:

где Те -Т + £, 1е(£) — индикатор отрезка [0,Те]. Величина е продолжения временного промежутка, на котором строится цепь Маркова, вводится с целью искусственного обрыва этой цепи. Для задачи оценки общего числа мономеров и димеров использовано одно приближение к функции ценности следующего вида:

Здесь же сделано замечание о начальном размере ансамбля частиц, необходимого для достижения требуемого уровня точности.

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

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

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

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

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

В численных экспериментах для модельных задач коагуляции полное ценностное моделирование дало выигрыш в трудоёмкости по сравнению с прямым моделированием порядка 10-100 раз.

Результаты главы 3 опубликованы в статьях [2-5,8].

В четвёртой главе представлены результаты численных расчётов. Решения тестовых задач для уравнений Больцмана и Смо-

луховского с помощью алгоритмов, предложенных в предыдущих двух главах, представлены в параграфах 4.1 и 4.2, соответственно. Численные результаты были проверены на тестовых задачах и подтвердили эффективность разработанных алгоритмов. В частности, с помощью достаточно точных расчётов показано, что относительная детерминированная погрешность имеет порядок С(1/Аг) по числу модельных частиц N.

Заключение содержит перечень основных результатов диссертационной работы.

ЗАКЛЮЧЕНИЕ

Сформулируем основные результаты диссертационной работы.

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

2. Построены новые алгоритмы коррелированного и весового моделирования для изучения зависимости решения кинетического уравнения Больцмана от начального и временнбго параметров.

3. Предложены и апробированы алгоритмы для частичного ценностного моделирования элементарных переходов при решении кинетических уравнений.

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

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

4. Разработанные весовые методы, в сочетании с методом коррелированной выборки, применены для исследования важной параметрической зависимости приближенного решения уравнения Больцмана ог числа модельных частиц. Для модельных ЛГ-частичных кинетических уравнений с помощью высокоточных тестовых расчётов впервые получен порядок относительной погрешности 1/N.

5. Для оценки функционалов от решения модельных задач осуществлен ряд численных экспериментов. В частности, для уравнения Больцмана в заданный момент времени оценивались

• 10-й момент скоростей частиц ансамбля;

• производные от этого момента по начальному и временному параметрам;

• частота частиц со скоростями, превышающими заданный уровень.

Для уравнения Смолуховского в заданный момент времени оценивались частота мономеров и суммарная частота мономеров и димеров

• для случая постоянных коэффициентов коагуляции;

• для случая линейных коэффициентов коагуляции.

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

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

СПИСОК РАБОТ ПО ТЕМЕ ДИССЕРТАЦИИ

Публикации в журналах списка ВАК

[1] Иванов М.С., Михайлов Г.А., Коротченко М.А., Рогазин-ский C.B. Глобально-весовой метод Монте-Карло для нелиней-

ного уравнения Больцмана // Журн. вычисл. мат. и мат. физики. 2005. Т. 45, № 10, С. 1860-1870.

[2] Михайлов Г.А., Коротченко М.А., Рогазинский С.В. Ценностные алгоритмы статистического моделирования для нелинейных кинетических уравнений // Доклады РАН. 2007. Т. 415, № 1. С. 26-30.

[3] М.А. Korotchenko, G.A. Mikhailov and S.V. Rogasinsky Value modifications of weighted statistical modeling for solving nonlinear kinetic equations // Russian J. Numer. Anal. Math. Modelling. 2007. Vol. 22, No. 5. P. 471-486.

[4] Коротченко M.A., Михайлов Г.А., Рогазинский С.В. Модификации весовых алгоритмов метода Монте-Карло для решения нелинейных кинетических уравнений // Журн. вычисл. мат. и мат. физики. 2007. Т. 47, № 12. С. 2110-2121.

[5] Коротченко М.А. Статистические алгоритмы ценностного моделирования для решения уравнения Смолуховского // Вычисл. технологии. 2008. Т. 13, Специальный выпуск 4. С. 68-74.

Прочие публикации

[6] Коротченко М.А. Новый глобально-весовой метод Монте-Карло для изучения параметрической зависимости решения нелинейного кинетического уравнения Больцмана // Труды конференции молодых ученых ИВМиМГ СО ГАН. - Новосибирск, 2004. С. 87-95.

[7] Коротченко М.А. Исследование зависимости приближенного решения уравнения Больцмана от числа модельных частиц // Труды конференции молодых ученых ИВМиМГ СО РАН. - Новосибирск, 2005. С. 79-87.

[8] Коротченко М.А., Урева Н.М. Алгоритм "ценностного" моделирования для решения нелинейных уравнений Больцмана и Смолуховского методом Монте-Карло // Труды конференции молодых ученых ИВМиМГ СО РАН. - Новосибирск, 2006. С. 119-129.

КОРОТЧЕНКО Мария Андреевна

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

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

Автореферат

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

Лицензия ИД n0 02202 от 30 июня 2000 г. Подписано в печать 16 М.2008 г.

Формат бумаги 60 х 841/16 Объем 1,0 п. л. 0,9 уч.-изд.л. Тираж 100 экз. Заказ №2

ООО «Омега Принт», Новосибирск-90, пр. Лаврентьева, 6

Введение

1. Элементы теории методов Монте-Карло для решения кинетических уравнений

1.1. Постановка задач.

1.1.1. Интегральные уравнении 2-го рода. Весовые оценки.

1.1.2. Модификации фазового пространства. Ценностное моделирование.

1.2. Математические модели многочастичных систем.

1.2.1. Пространственно-однородная релаксация простого од-нокомпонентного газа.

1.2.2. Процесс пространственно-однородной коагуляции

1.3. Построение базового интегрального уравнения.

1.4. Численная оценка функционалов.

1.4.1. Весовые оценки для кинетических уравнений.

1.4.2. Алгоритмы ценностного моделирования для кинетических уравнений

2. Весовые алгоритмы метода Монте-Карло для решения уравнения Больцмана

2.1. Тестовая задача для уравнения Больцмана.

2.2. Алгоритм моделирования для оценки решения задачи Коши

2.3. Весовые алгоритмы и параметрические оценки.

2.3.1. Изучение зависимости решения от параметра начального распределения в.

2.3.2. Исследование зависимости решения от параметра временного распределения ао

2.4. О выборе количества частиц N.

2.4.1. Коррелированное моделирование ансамблей частиц при различных N

2.5. Весовые модификации алгоритмов для оценки "хвоста" скоростного распределения.

2.5.1. Весовое моделирование начального распределения скоростей

2.5.2. Алгоритм ценностного моделирования номера взаимодействующей пары.

2.5.3. Ценностное моделирование длины пробега.

3. Весовые алгоритмы метода Монте-Карло для решения уравнения Смолуховского

3.1. Тестовые задачи для уравнения Смолуховского.

3.2. Весовые оценки

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

3.3.1. Алгоритмы ценностного моделирования длины свободного пробега

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

3.4. Алгоритмы ценностного моделирования для задачи с линейными коэффициентами коагуляции.

4. Численные результаты

4.1. Решение тестовых задач для кинетического уравнения Больцмана.

4.1.1. Использование глобального весового метода.

4.1.2. Изучение зависимости оценок от количества частиц в ансамбле

4.1.3. Весовые оценки "хвоста" скоростного распределения

4.2. Решение тестовых задач для кинетического уравнения Смолуховского

4.2.1. Ценностное моделирование для решения задачи коагуляции с постоянными коэффициентами.

4.2.2. Ценностное моделирование для решения задачи коагуляции с линейными коэффициентами.

Методы статистического моделирования (методы Монте-Карло) — это "численные методы решения математических задач при помощи моделирования случайных величин и статистической оценки их характеристик" [32]. Их применение для решения разнообразных задач в различных областях физики, математики и техники имеет богатую историю [9, 10, 32]. В 50-х годах состоялись несколько первых симпозиумов по методам Монте-Карло, которые продемонстрировали широкие перспективы их применения. Хотя методы Монте-Карло были известны и до 40-х годов, тогда они не получили широкого распространения из-за больших объемов вычислений. Появление компьютеров сделало их практически применимыми и за последние 50 лет вместе с развитием компьютерных технологий методы Монте-Карло всё более активно используются во многих научных исследованиях.

Резко возросший в мире за последние два десятилетия интерес к исследованиям в области статистического моделирования обусловлен двумя основными причинами. Первая связана с тем, что методы Монте-Карло оказались удачно адаптируемыми к современной вычислительной технике с высоким уровнем распараллеливания. Вторая обусловлена новым уровнем описания задачи — статистический подход стал чрезвычайно часто использоваться уже на этапе постановки задачи и выбора модели. Это в полной мере относится и к задачам теории кинетических уравнений, которые описывают эволюцию функции распределения F(t, v, х) молекул или других объектов (электронов, ионов, звезд, галактик, пор и др.) по скоростям v и пространству х в момент времени Ь. Это означает, что число частиц в элементе фазового объема dv drc есть F(t: v, х) di; da;. Простейшее такое уравнение — уравнение свободного движения:

В настоящей диссертации будут рассмотрены два кинетических уравнения в пространственно-однородном случае: уравнение Больцмана и уравнение Смолуховского.

Уравнение Больцмана — интегро-дифференциальное уравнение, описывающее поведение разреженного газа, которое было выведено Людвигом Больцманом в 1872 году. Оно использовалось для обоснования молекулярно-кинетической теории, второго закона термодинамики о возрастании энтропии, выводе уравнений гидродинамики и до сих пор остаётся основой кинетической теории газов.

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

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

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

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

Расщепление непрерывного процесса движения молекул и их столкновений в разреженном газе на два последовательных этапа в течение временного шага At лежит в основе метода прямого статистического моделирования (ПСМ, Direct Simulation Monte Carlo Method), в котором используются вероятностные представления для моделирования течений разреженного газа. Он впервые был предложен Г. Бёрдом (см., например, [2]). Дальнейшее развитие метода шло по пути построений различных модификаций моделирования столкновительной релаксации в модельном газе. В работе [8] для расчёта пространственно-однородной релаксации была использована УУ-частичная модель Каца [14]. Более простой алгоритм приближённого моделирования столкновений в N-частичной модели газа, основанный на вероятностной схеме испытаний Бернулли, был построен в [1]. Принципиально другой подход, основанный на получении статистического алгоритма моделирования однородной релаксации из разностной по времени формы уравнения Больцмана, был предложен Нанбу в [41]. Отличие этой схемы заключается в том, что законы сохранения импульса и энергии в системе выполняются только в среднем, а не при каждом столкновении. Также строились весовые модификации метода ПСМ, в которых отдельным частицам присваивались "веса", связанные с искусственным перераспределением типов частиц в начальный момент времени. В работе [31] рассмотрена наиболее общая весовая модификация метода ПСМ такого типа под названием метода "дополнительной переменной".

Наиболее простой способ моделирования процесса коагуляции это прямое моделирование, называемое также процессом Маркуса

Лушникова. В процессе прямого моделирования на каждом временном шаге согласно коэффициентам коагуляции Кц выбираются частицы размера г и j, а в результате появляется новая частица размера г + j. Такие вычислительные модели по своей простоте подобны модели Бёрда [2] для уравнения Больцмана; они были описаны в [21, 35, 39]. Позже были предложены различные модификации алгоритмов, увеличивающие эффективность, такие как "метод мажорантной частоты" [13, 42], а также его обобщение — "метод фиктивных прыжков" [34]. Некоторые алгоритмы используют дополнительный параметр, множитель к шагу по времени, который позволяет учитывать несколько парных взаимодействий между частицами за один временной шаг. Такие алгоритмы (см., например, [42]) аналогичны алгоритму Нанбу [41] для уравнения Больцмана.

Настоящая диссертация посвящена построению, обоснованию и применению весовых модификаций алгоритмов для оценки функционалов от решений пространственно-однородных кинетических уравнений. Для этого использовано лилейное интегральное уравнение Фредгольма второго рода, которое, в случае уравнения Больцмана, эквивалентно N-частичному уравнению Каца, чья связь с нелинейным кинетическим уравнением известна [14]. В случае уравнения Смолуховского интегральное уравнение связано с уравнением Колмогорова, которое представляет собой вероятностное описание эволюции системы из N частиц и даёт приближённое решение соответствующего нелинейного кинетического уравнения (см. [30, 31]).

Таким образом решение кинетических уравнений оценивается с помощью статистического моделирования марковских процессов эволюции соответствующих iV-частичных ансамблей, фазовые состояния которых меняются вследствие парных взаимодействий частиц (см., например, [14, 29]). Однако, связанные с этим процессом стандартные базовые интегральные уравнения использовать непосредственно для построения весовых модификаций статистического моделирования затруднительно, так как их ядра представляет собой сумму взаимосингулярных слагаемых. Это затруднение было преодолено в [29] посредством модификации фазового пространства системы путем введения номера взаимодействующей пары частиц в число координат фазового пространства системы. Это привело к "расслоению" распределения взаимодействий в системе по номеру пары взаимодействующих частиц 7г = [г. j). Данный прием позволил в [29] сформулировать новое интегральное уравнение, на основе которого удобно строить стандартные весовые модификации статистического моделирования рассматриваемой многочастичной системы, поскольку ядро уравнения содержит сингулярности лишь в виде сомножителей.

В настоящей диссертации для такого ядра разрабатываются алгоритмы с "глобальным" весом, который после каждого элементарного перехода в моделируемой цепи Маркова домножается на стандартный для теории методов Монте-Карло [25, 40] весовой множитель. Это позволяет распространить хорошо разработанную теорию весовых методов [25, 40] на рассматриваемый класс задач. Более того, это даёт возможность оценивать параметрические производные от решения, что особенно важно при численном исследовании влияния различных параметров на решение нелинейных кинетических уравнений.

Величиной, определяющей эффективность алгоритма метода Мон-тс-Карло при его практическом использовании, является трудоёмкость S(C) соответствующей случайной оценки С (см., например, [9, 27]). В диссертационной работе под трудоёмкостью S(Q будем понимать среднее количество вычислительной работы, необходимой для достижения заданной погрешности, которая пропорциональна величине

S(C)=VarC-i(C), где £(£) — среднее время, затрачиваемое на моделирование одного выборочного значения случайной величины С (например, для кинетических уравнений величина t(Q определяется средним числом переходов в цепи до момента выхода за пределы заданного временнбго интервала). Как известно, среднее арифметическое выборочных значений Sm = ^ ^^ Сг5 г=1 которое приближает математическое ожидание оценки при достаточно большом М имеет нормальное распределение и

V - SM\ где о; — уровень доверия, Р(а) — константа, определяющаяся стандартным нормальным распределением по значению а. Из этого равенства видно, что вероятностная погрешность статистического алгоритма пропорциональна величине \j . Погрешность оценки обычно состоит из двух частей: детерминированной и вероятностной, поэтому целесообразно выбирать параметры алгоритма таким образом, чтобы обе эти погрешности имели один и тот же порядок. Детерминированная погрешность возникает, например, при замене системы физических частиц модельным ансамблем из малого количества частиц.

Диссертация состоит из введения, четырёх глав, разбитых на параграфы, заключения и списка литературы из 42 наименований.

Заключение

Сформулируем основные результаты диссертационной работы.

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

2. Построены новые алгоритмы коррелированного и весового моделирования для изучения зависимости решения кинетического уравнения Больцмана от начального и временного параметров.

3. Предложены и апробированы алгоритмы для частичного ценностного моделирования элементарных переходов при решении кинетических уравнений. о Для моделирования длины свободного пробега iV-частим ной системы использовано несколько приближений к известной функции ценности, которые могут быть использованы в реальных задачах с неизвестными решениями. о Для ценностного моделирования номера пары взаимодействующих частиц, в зависимости от оцениваемого функционала, ценность ансамбля полагали пропорциональной количеству мономеров или суммарному числу мономеров и димеров (для уравнения Смолуховского), а для уравнения Больцмана — количеству частиц со скоростями выше заданного уровня.

4. Разработанные весовые методы, в сочетании с методом коррелированной выборки, применены для исследования важной параметрической зависимости приближенного решения уравнения Больцмана от числа модельных частиц. Для модельных TV-частичных кинетических уравнений с помощью высокоточных тестовых расчётов впервые получен порядок относительной погрешности 1/N.

5. Для оценки функционалов от решения модельных задач осуществлен ряд численных экспериментов. В частности, для уравнения Больцмана в заданный момент времени оценивались

• 10-й момент скоростей частиц ансамбля;

• производные от этого момента по начальному и временному параметрам;

• частота частиц со скоростями, превышающими заданный уровень.

Для уравнения Смолуховского в заданный момент времени оценивались частота мономеров и суммарная частота мономеров и димеров

• для случая постоянных коэффициентов коагуляции;

• для случая линейных коэффициентов коагуляции.

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

1. Велоцерковский О.М., Япицкий В.Е. Статистический метод "частиц в ячейках" для решения задач динамики разреженного газа // Журн. вычисл. мат. и мат. физики. 1975. Т. 15, № 5. С. 1195-1208.

2. Бёрд Г. Молекулярная газовая динамика. М.: Мир, 1981.

3. Бобылев А.В. О точных решениях уравнения Больцмана // Доклады АН СССР. 1975. Т. 225, № 6. С. 1296-1299.

4. Бобылев А.В. Точные решения нелинейного уравнения Больцмана и теория релаксации максвелловского газа // Теоретическая и математическая физика. 1984. Т. 60, № 2. С. 280-310.

5. Владимиров B.C., Марчук Г.И. Об определении сопряженного оператора для нелинейных задач // Доклады РАН. 2000. Т. 372, № 2. С. 165-168.

6. Волощук В.М. Кинетическая теория коагуляции. Ленинград: Гидрометеоиз-дат, 1984.

7. Денисик С.А., Лебедев С.Н., Малама Ю.Г. Об одной проверке нелинейной схемы метода Монте-Карло // Журн. вычисл. мат. и мат. физики. 1971. Т. 11, N2 3. С. 783-785.

8. Денисик С.А., Малама Ю.Г., Лебедев С.Н., Полак Л.С. Решение задач физической и химической кинетики методом Монте-Карло // Применение вычисл. матем. в хим. и физ. кинетике. М.: Наука, 1969. С. 179-231.

9. Ермаков С.М., Михайлов Г.А. Статистическое моделирование. М.: Наука, 1982.

10. Ермаков С.М., Некруткин В.В., Сипин А. С. Случайные процессы для решения классических уравнений математической физики. М.: Наука, 1984.

11. Иванов М.С., Михайлов Г.А., Коротченко М.А., Рогазипский С.В. Глобально-весовой метод Монте-Карло для нелинейного уравнения Больцмана // Журн. вычисл. мат. и мат. физики. 2005. Т. 45, № 10, С. 1860-1870.

12. Иванов М.С., Рогазипский С.В. Метод прямого статистического моделирования в динамике разреженного газа. Новосибирск: ВЦ СО АН СССР, 1988.

13. Иванов М.С., Рогазипский С.В. Экономичные схемы статистического моделирования течений разреженного газа // Математическое моделирование. 1989. Т. 1, № 7. С. 130-145.

14. Кац М. Вероятность и смежные вопросы в физике. М.: Мир, 1965.

15. Коротченко М.А. Исследование зависимости приближенного решения уравнения Больцмана от числа модельных частиц // Труды конференции молодых ученых ИВМиМГ СО РАН. Новосибирск, 2005. С. 79-87.

16. Коротченко М.А. Новый глобально-весовой метод Монте-Карло для изучения параметрической зависимости решения нелинейного кинетического уравнения Больцмана // Труды конференции молодых ученых ИВМиМГ СО РАН. -Новосибирск, 2004. С. 87-95.

17. Коротченко М.А. Статистические алгоритмы ценностного моделирования для решения уравнения Смолуховского // Вычисл. технологии. 2008. Т. 13, Специальный выпуск 4. С. 68-74.

18. Коротченко М.А., Михайлов Г.А., Рогазипский С.В. Модификации весовых алгоритмов метода Монте-Карло для решения нелинейных кинетических уравнений // Журн. вычисл. мат. и мат. физики. 2007. Т. 47, № 12. С. 2110-2121.

19. Коротченко М.А., Урева Н.М. Алгоритм "ценностного" моделирования для решения нелинейных уравнений Больцмана и Смолуховского методом Монте-Карло // Труды конференции молодых ученых ИВМиМГ СО РАН. Новосибирск, 2006. С. 119-129.

20. Леонтович М.А. Основные уравнения кинетической теории газов с точки зрения теории случайных процессов // Журнал экспериментальной и теоретической физики. 1935. Т. 5. С. 75-79.

21. Лушников А.А. Некоторые новые аспекты теории коагуляции // Изв. АН СССР, Физ. атмосферы и океана. 1978. Т. 14, № 10. С. 1048-1055.

22. МарченкоМ.А., Михайлов Г.А. Распределенные вычисления по методу Монте-Карло // Автоматика и Телемеханика. 2007. № 5. С. 157-170.

23. Марчук Г.И., Агошков В.И., Шутяев В.П. Сопряженные уравнения и методы возмущений в нелинейных задачах математической физики. М.: Наука, 1993.

24. Михайлов Г.А. Весовые алгоритмы статистического моделирования. Новосибирск: ИВМиМГ СО РАН, 2003.

25. Михайлов Г. А. Весовые методы Монте-Карло. Новосибирск: Изд-во СО РАН, 2000.

26. Михайлов Г. А. Построение весовых методов Монте-Карло на основе увеличения размерности фазового пространства // Доклады РАН. 2003. Т. 389, № 4. С. 461464.

27. Михайлов Г.А., Войтишек А.В. Численное статистическое моделирование (Метод Монте-Карло). М.: Издательский центр "Академия", 2006.

28. Михайлов Г.А., Коротченко М.А., Рогазинский С.В. Ценностные алгоритмы статистического моделирования для нелинейных кинетических уравнений // Доклады РАН. 2007. Т. 415, № 1. С. 26-30.

29. Михайлов Г.А., Рогазинский С.В. Весовые методы Монте-Карло для приближенного решения нелинейного уравнения Больцмана // Сиб. мат. журн. 2002. Т. 43, № 3. С. 620-628.

30. Михайлов Г.А., Рогазинский С.В. Поправки к статье "Весовые методы Монте-Карло для приближенного решения нелинейного уравнения Больцмана" // Сиб. мат. журн. 2003. Т. 44, № 2. С. 473-474.

31. Михайлов Г.А., Рогазинский С.В., Урева Н.М. Весовой мегод Монте-Карло для приближенного решения нелинейного уравнения коагуляции // Журн. вычисл. мат. и мат. физики. 2006. Т. 46, № 4. С. 714-725.

32. Рогазинский С.В. Об одном подходе к решению однородного уравнения Больцмана // Журн. вычисл. мат. и мат. физики. 1987. Т. 27, № 4. С. 564-574.

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

34. Aldous D.J. Deterministic and Stochastic Models for Coalescence (Agregation and Coagulation): a Review of the Mean-Field Theory for Probabilities // Bernulli. 1999. V. 5, № 1. P. 3-48.

35. Eibeck A., Wagner W. An efficient stochastic algorithm for studing coagulation dynamics and gelation phenomena // SIAM J. Sci. Comput. 2000. V. 22, № 3. P 802-821.

36. Gillespie D. T. An exact method for numerically simulating the stochastic coalescence process in a cloud // J. Atmospheric Sci. 1975. V. 32, № 10. P 1977-1989.

37. Gillespie D. T. The stochastic coalescence model for cloud droplet growth // J. Atmospheric Sci. 1972. V. 29, № 8. P 1496-1510.

38. Ivanov M.S., Rogasinsky S.V. Analysis of numerical techniques of the direct simulation Monte Carlo method in the rarefied gas dynamics // Sov. J. Numer. Anal. Math. Modelling. 1988. Vol. 3, No. 6. P. 453-465.

39. M.A. Korotchenko, G.A. Mikhailov and S.V. Rogasinsky Value modifications of weighted statistical modeling for solving nonlinear kinetic equations // Russian J. Numer. Anal. Math. Modelling. 2007. Vol. 22, No. 5. P. 471-486.

40. Marcus A.H. Stochastic coalescence // Technometrics. 1968. V. 10, № 1. P. 133-143.

41. Mikhailov G.A. Parametric estimates by the Monte Carlo method. Utrecht: VSP,

42. Nanbu K. Direct simulation scheme derived from Boltzmann equation. I. Monocom-ponent gases // J. Phys. Soc. of Japan. 1980. V. 49, № 5. P. 2042-2049.

43. Sabelfeld K.K., Rogasinsky S.V., Kolodko A.A., Levykin A.I. Stochastic algorithms for solving Smoluchovsky coagulation equation and applications to aerosol growth simulation // Monte Carlo Methods and Appl. 1996. V. 2, № 1. P. 41-87.1999.