Введение
История вопроса
Развитие методов Монте-Карло связано с применением стохастических алгоритмов в области теории переноса излучения ([8], [13]). В частности, схема Неймана-Улама ([10]) возникла в связи с проблемой точности моделирования решения задачи переноса излучения. При этом схема Неймайй-*УлаМ&'с^стоит в моделировании цепей Маркова и, с другой стороны, тесно связана с итерационным процессом
Xn+1 = АХп + F.
Развитые в связи с этим алгоритмы расходовали небольшой объем машинной памяти, однако для их применения требовалось выполнение условия сходимости мажорантного процесса, что резко ограничивало область применимости схемы Неймана-Улама.
Э го обстоятельство потребовало развития повьгх методов, применимых к более широкому классу уравнений. Существенный шаг в этом направлении был сделан в работах [5], [3]. Во второй нз них исследование было проведено в применении к разностному аналогу волнового уравнения, сформулированы условия стохастической устойчивости, приведен пример неустойчивой схемы.
Понятие устойчивости играет большую роль в исследовании свойств стохастических алгоримов. Во многих случаях неустойчивость, проявляющаяся в экспоненциальном росте дисперсии оценок, делает применение схемы Неймана-Улама невозможных!.
Альтернативой схеме Неймана-Улама могут служить, в частности, методы, используюшне аппарат стохастических дифференциальных уравнений. Для решения эллиптических уравнений по-
2
строены эффективные алгоритмы такого рода ([18]). Строго говоря, такие алгоритмы не укладываются в схему Неймана-Улама. Безусловно, представляет интерес разработка общего языка, позволяющего с единой точки зрения взглянуть на эти два типа стохастических алгоритмов. Определенный прогресс в установлении такой общности был достигнут в работе [4].
Однако, задачи этого типа достаточно сложны. Их анализ осуществляется сравнительно просто после введения дискретного времени. Заметим, что при решении практических задач дискретизация по времени является стандартным приемом и служит достаточно общим инструментом исследования. В особенности она важна, когда применяются метод расщепления или метод дробных шагов ([21]), позволяющие заменить исходную задачу последовательностью более простых задач.
Если исходное уравнение имеет вид
(0.1)
где - линейный оператор, а и - искомая (вектор-) функция, то после дискретизации с шагом Д/ по времени получаем приближенные соотношения:
^Д+1 м уЯ
---------= Ьпии + Vй (явная схема) или (0.2)
О С мП+1 _ ип
— = 1п+1нп+) +г,п (неявная схема). (0.3)
Если же /,{ представить в виде суммы /,< = то можно
построить явно-неявную схему:
= £(.)„» + £Й1“”+1+*'”- С-4)
Полученное рекуррентное соотношение определяет рекуррентный характер связанных с ним стохастических алгоритмов.
Язык рекуррентного алгоритма дает возможность связать разные по своей природе стохастические алгоритмы и является инструментом для исследования асимптотического поведения их трудоемкости.
— = Ь1и + ь{1),
3
Возникает также вопрос о предельном поведении серий последо-вательностей оценок рекуррентного алгоритма п связанных с ними марковских процессов прп Дг —► 0. В предлагаемой работе получены условия сходимости конечномерных распределений таких процессов к конечномерным распределениям некоторого диффузионного процесса. Полученные результаты позволяют проследить нюансы взаимосвязи предельных переходов по пространственным и временной переменных! с основными классами оценок Монте-Карло.
В диссертационной работе получены некоторые обобщения результатов, опубликованных ранее в работах [3], [4]. Эти обобщения потребовали, в частности, введения некоторого нового класса цепей Маркова, с помощью которого строится и обосновывается рекуррентный алгоритм.
Структура диссертации
В главе 1 описан рекуррентный алгоритм построения несмещенных оценок г,л > 0 величин ии,п > 0, удовлетворяющих рекуррентному соотношению
ип+| = £0 )ад«+1 + /£)„« + „»+> (0.5)
с заданными значениями и°,ип,л > 0 и известными операторами
Г (О т(2)
Дп •
Для этого сначала вводится специальный класс цепей Маркова, множество состояний которых - целочисленная решетка, ограниченная слева, справа и снизу. Траектории цепей устроены так: сначала движение по горизонтали, затем переход на верхний слой, после чего цепь обрывается. Движение вниз запрещено.
При этом параметры марковской цепи связаны с операторами условиями согласования.
По аналогии с классической схемой Неймана-Улама выделены четыре основные оценки: "сверху вниз”, "снизу вверх”, по поглощению п по столкновениям. Доказана их несмещенность, вычислены вторые моменты.
Рекуррентный характер алгоритма позволяет понизить дисперсию, вводя дополнительное усреднение на каждом шаге итераций.
4
В главе 1 также подробно обсуждаются примеры применения рекуррентного алгоритма: для решения уравнений в частных производных, интегральных уравнении, а также систем алгебраических уравнений. В последнем примере применение рекуррентного алгоритма позволяет расширить по сравнению со схемой Неймана-Улама границы применимости методов Монте-Карло.
Глава 2 посвящена вопросам устойчивости рекуррентного алгоритма. При этом мы ограничиваемся исследованием спектральной устойчивости. Так, слабая стохастическая устойчивость определяется как свойство алгоритма, заключающееся в линейном росте нормы матрицы ковариации при росте п. Сильная стохастическая устойчивость определяется как ограниченность последовательности норм матриц ковариации оценок при всех п.
Таким образом, задача сводится к отысканию линейного рекуррентного соотношения, связывающего элементы матриц ковариаций оценок £",п > 0:
со^а+1 = ЯпсоуС + рп. (0.6)
Для этого используется следующий прием: матрица ковариаций "растягивается” в вектор, составленный из ее элементов.
Устойчивость рекуррентного алгоритма, для оценок которого выполняется соотношение типа (0.6), определяется значением нормы операторов Кп и ограниченностью последовательности ри. Получены достаточные, а также в некоторых случаях необходимые и достаточные условия слабой и сильной стохастической устойчивости.
В случае интегральных уравнений все соотношения понимаются в слабом смысле.
Операторы, входящие в соотношения (0.6), получены в явном виде для общего случая и для некоторых частных, более удобных в практических задачах.
Завершает главу серия примеров исследования устойчивости рекуррентных алгоритмов при помощи доказанных теорем.
В главе 3 рекуррентный алгоритм рассматривается в применении к системам линейных дифференциальых уравнений. После дискретизации по времени с шагом Д(. они сводятся к (явно-неявной) рекуррентной процедуре.
Рассмотрим последовательность серий случайных величин, где каждая серия есть носледовательность оценок £п величин, удо-
влетворяюшнх этой рекуррентной процедуре, полученная С ПОМОЩЬЮ рекуррентного алгоритма при фиксированном Д*. Рассмотрим кусочно-линейный процесс с вершинами в точках
(£п,пД«),п > 0.
В главе 3 изучаются условия сходимости конечномерных распределений процессов и распределений функционалов от к соответствующим распределениям марковского процесса являющегося решением некоторого стохастического дифференциального уравнения.
При некоторых условиях на параметры алгоритма получены достаточные условия сходимости и явный вид предельного стохастического дифференциального уравнения.
Полученные результаты служат теоретической основой для построения алгоритмов приближенного (в смысле сходимости конечномерных распределений) решения стохастических дифференциальных уравнений произвольной размерности. Пример такого типа завершает главу 3.
В главе 4 предложен метод частиц 1>ешення уравнений параболического типа с граничными условиями. Под методом частиц понимают, вообще говоря, любой алгоритм моделирования траектории случайного процесса. В диссертации рассматривается более узкий класс методов моделирования марковских процессов с дискретным временем. При этом учет границ осуществляется за счет искусственного приема поддержания постоянного числа частиц.
В этом смысле результаты главы 3 являются обобщением результатов, полученных в работах [11], [10], [12].
Далее, показано, что метод части и может быть определен на языке рекуррентного алгоритма, и для него, следовательно, применим разработанный в предыдущих главах аппарат.
Показано, что опенка метода частиц является асимптотически несмещенной при N —» сю, где N число частиц. Получена оценка порядка трудоемкости метода. Это позволяет сравнивать алгоритм метода частиц с другими стохастическими алгоритмами.
Об обозначениях
Перечислим обозначения, принятые в диссертационной работе.
б
1) Под точкой подразумеваем вектор-столбец. Для вся-
кого вектора х = (xi,..., ха) = (яі)<=і:
1. Xі - транспонированный вектор (строка);
2. |я| - модуль вектора:
И = У*2 + ... + х2.
3. Обозначим символом hk к-й базисный вектор:
h\k) = l \,1 = к] при і = I,... ,d,k = I,... >d.
• 0 иначе
2) Для всякой матрицы .4 порядка dх d:
1. Ат - транспонированная матрица;
2. ding А - диагональ матрицы A, diag А Є Rrf;
3. А(.4) - максимальное по модулю собственное число матрицы А.
Единичную матрицу будем обозначать буквой I.
3) Обозначения, связанные с интегральными уравнениями.
1. Для всякой ограниченной функции f(x) и заряда //((/#)
</* /*) = J /(*)/*(<**)•
2. Дельта-меру, сосредоточенную в точке X, будем обозначать 6x{dx).
3. Символом |/<| будем обозначать полную вариацию заряда д.
4. Для двух зарядов (мер) //, v символом /іхі/ будем обозначать их прямое произведение. Будем также использовать обозначение
/і2 = fiXfl.
5. Символом X будем обозначать тождественный оператор.
7
6. Для всякого интегрального оператора А символом Л(Л} будем обозначать наибольшее по модулю собственное число оператора А.
4) В выражениях для вторых моментов оценок Монте-Карло обычно встречаются матрицы и векторы, составленные из произведений элементов двух соответствующих матриц или векторов. Во избежание путаницы для таких операций введем специальные обозначения.
Для любых двух векторов х, у € К :
1. [х,у]с - вектор, составленный из произведений соответствующих компонент векторов X и у:
[х,!/!* = (£іУі,...,гдо);
Для векторов, составленных из квадратов компонент вектора, будем прпмепять обозначение
х[2] = [х,х)е = (х?,...,х5).
И вообще, для любой степени п будем обозначать
*1"' = м *3).
2. р - вектор, составленный из частных соответствующих компонент векторов х и у:
X _ ҐХ\ Х(і
І 1 ' 1 »
У \У і Ул
3. <т(х, у) - матрица (ІХ(іу составленная из попарных произведений компонент векторов х и у:
<г(*,У) = \\*іУі\\Ь=і •
Вместо сг(х,х) будем для краткости писать а(х).
5) Введем обозначения для операций поэлементного умножения и деления матриц. Для любых гіх^-матриц А = 1М.* = 1М
8
1. [А,В]е - матрица, составленная из произведений соответствующих элементов матриц А и В:
И,в]е = ||«уМ1?л=1;
2. ^ - матрица, составленная из частных соответствующих элементов матриц А а В:
Л
В
а
'ч
і, 3=1
Для матриц, составленных из квадратов элементов матрицы, будем примепять обозначение
дй = [А, А], = ||4||у=1 ■
И вообще, для любой степени п будем обозначать
№ = 1К11и=. •
3. Для всяких матрицы А и вектора х будем обозначать
А
А
х
хз
10=1
6) Пусть А = ||ау|| - матрица с/х</. Тогда будем обозначать
1. Вектор длины <12
1 = (ву) 6 К«5;
2. Матрица dy.fi
0(Л) = 11«(А),Л.,•
При этом порядок, в котором элементы матрицы выписываются
в строку, не имеет значения, важно только, чтобы Л и в(А) были согласованы с точки зрения этого порядка. Для исследования устойчивости (глава 2) нам понадобятся специфические обозначения.
7) Пусть А, В - матрицы dxd, х € К*.
9
- Київ+380960830922