Дискретно-стохастические численные методы

Оглавление
Введение................................................................... 3
Глава 1. Функциональная сходимость численных моделей случайных
полей.................................................................. 26
1.1. Специальные условия функциональной сходимости последовательностей случайных функций в пространствах С(Г) и О(Г).......................... 26
1.2. Сходимость численных рандомизированных спектральных моделей случайных полей........................................................ 49
1.3. Слабая сходимость вычислительных моделей случайных полей, связанных
с точечными потоками Пальма........................................ 64
Глава 2. Дискретно-стохастические методы вычисления многократных интегралов............................................................. 71
2.1. Геометрический метод Монте-Карло и его модификации................ 71
2.2. Построение выборки по важности с использованием аппроксимации Стренга-Фикса.......................................................... 77
2.3. Дискретно-стохастические методы уменьшения дисперсии.............. 84
2.4. Использование существенной выборки................................ 97
Глава 3. Дискретно—стохастические методы аппроксимации интегралов, зависящих от параметра................................................ 101
3.1. Стохастические оценки в узлах сетки. Достаточные условия сходимости метода зависимых испытаний............................................ 101
3.2. Построение верхних границ погрешностей функциональных алгоритмов 108
3.3. Условная оптимизация функциональных алгоритмов................... 120
3.4. Многоуровневый метод зависимых испытаний......................... 125
Глава 4. Дискретно-стохастические методы аппроксимации решения
интегрального уравнения второго рода.................................. 131
4.1. Несмещенные стохастические оценки в узлах сетки. Верхние границы компонент погрешности................................................. 131
4.2. Общие утверждения о сходимости и условно-оптимальные параметры функциональных алгоритмов с несмещенными оценками в узлах сетки 143
4.3. Многомерный аналог метода полигона частот........................ 151
Глава 5. Применение дискретно-стохастических численных методов ... 166
5.1. Об одной стохастической задаче теории переноса излучения......... 166
5.2. Исследование смешанных алгоритмов решения интегро-дифференциальных систем................................................................ 179
5.2.1. Модельная каталитическая задача............................ 179
5.2.2. Одномерная задача радиационно-кондуктивного теплопереноса 183
5.3. Об одном итерационном методе решения нелинейных интегральных уравнений............................................................. 204
5.4. Тестирование дискретно-стохастических численных методов.......... 213
5.4.1. Использование спектральных моделей при тестировании
алгоритмов численного интегрирования........................ 213
1
5.4.2. Тестирование алгоритмов численного интегрирования из Разделов
2.2 и 2.3.................................................. 219
5.4.3. Вычисление констант в выражениях для условно-оптимальных параметров из Разделов 3.3, 4.2 и 4.3............................ 222
5.4.4. Тестирование алгоритмов аппроксимации интеграла, зависящего от параметра, из Раздела 3.1........................................ 224
5.4.5. Тестирование многоуровневого метода зависимых испытаний из Раздела 3.4...................................................... 227
5.4.6. Тестирование алгоритмов аппроксимации решения интегрального уравнения второго рода из Разделов 4.1 и 4.3..................... 230
Приложение 1. Реализация случайных величин и случайных векторов
с кусочно-полиномиальными плотностями........................... 234
Приложение 2. Методологические основы электронного учебника по
методам Монте-Карло............................................. 240
Заключение............................................................. 245
Литература............................................................. 247
2
Введение
С развитием вычислительной техники возрастает интерес к численным методам решения прикладных задач, в частности, к статистическому моделированию (или методу Монте-Карло) [1 - 15, 16 - 20, 21 - 23] (см. также обзоры литературы в этих книгах); здесь и далее ссылки на работы автора диссертации выделены жирным шрифтом. В теории методов Монте-Карло можно выделить пять основных разделов [4 - 6, 10]:
а) численное моделирование случайных величин, векторов и функций;
б) вычисление многократных интегралов;
в) приближение интегралов, зависящих от параметра;
г) решение интегральных уравнений второго рода;
д) приложения методов Монте-Карло к задачам вычислительной математики и математической физики.
Традиционно методы Монте-Карло рассматриваются в качестве альтернативных "детерминированным” численным методам (в частности, конечно-разностным и конечно-элементным схемам) (см., например, [24 - 27]). Однако во многих случаях эффективными оказываются алгоритмы, содержащие в себе элементы детерминированных и стохастических численных схем. Такие комбинированные алгоритмы будем называть дискретно-стохастическими численными методами. Исследованию этих методов посвящена данная работа, которая состоит из Введения, пяти глав, двух приложений, Заключения и списка литературы из 267 наименований.
Следует сразу отметить, что спектр дискретно-стохастических численных методов достаточно широк. Комбинированные алгоритмы возникают во всех перечисленных выше разделах теории методов Монте-Карло. Для каждого из разделов мы приведем важнейшие примеры таких алгоритмов и изучим вопросы сходимости и оптимизации соответствующих численных схем.
Начнем с моделирования случайных величин. Как хорошо известно [1 - б, 10, 16], численная реализация случайных величин состоит из двух этапов:
1) реализуются значения осу,... , а* стандартного случайного числа а, равномерно распределенного па полуинтервале [0,1), с помощью специального устройства или программы, которое называется генератором случайных (псевдослучайных) чисел;
2) с помощью некоторых преобразований полученных чисел {ау} вычисляются значения случайных величин с более сложными законами распределения.
Большинство расчетов по методу Монте-Карло произведено и производится с помощью генераторов псевдослучайных чисел, представляющих собой некоторые вычислительные программы, которые, в свою очередь, реализуют рекуррентные последовательности вида
ап+1 = 5(ап), (0.1)
или их модификации; здесь 5 - неслучайная функция, определенная на отрезке [0; 1]. То есть по существу речь идет о детерминированной последовательности, обладающей свойствами последовательности независимых случайных величин. Наиболее часто используется функция $ вида 6т(г) = {М х] для достаточно большого натурального множителя А/, здесь {Л} обозначает дробную часть числа А. В этом случае алгоритм (0.1) называется мультипликативным методом вычетов [4, 5, 10, 16, 23]. Таким образом, в самой основе алгоритмов численного статистического моделирования, на уровне генераторов псевдослучайных чисел, используются по сути комбинированные методы вида (0.1).
3
При моделировании дискретной случайной величины Ç с конечным числом значений ДГ|,... ,хп и распределением вероятностей Р(£ = а;,-) = /?,, i = 1,..., п используется следующая стандартная процедура.
АЛГОРИТМ 0.1 [1 - б, 10, 17]. Реализуем значение Q := а и полагаем т := 1. Производим переприсваивапие
Q'=Q-Pm (0.2)
(то есть заносим новое значение (a—pi) в ячейку Q). Если повое Q не положительно, то о качестве т выбираем текущее его значение и считаем, что для данного испытания £ — хт. Если Q > 0 в m = п - 1, то £ = хп. Если же (J>0«m<n-1, то в производим персприсваивания т := т + 1 и (0.2) и вновь производим проверку Q на положительность и т.д.
Алгоритм 0.1 реализуем также в случае моделирования дискретной случайной величины со счетным числом значений; при этом должны иметься формулы вычисления Pi — ш ИЛИ рекуррентного пересчета р,+1 = /2(1»Pi) вероятностей, и перед очередным переприсваиванием (0.2) должно производиться вычисление вероятности рт. Реализация случайной величины ( с конечным числом значений заметно упрощается, когда все значения xl5. . .,хп равновероятны, то есть все р,- равны 1 /гг (такое распределение вероятностей называется дискретным равпо.иериым [28]):
m=|anj + l, £ = хт, (0.3)
где [AJ - целая часть числа Л [1 - б, 10, 17].
Следующий прием, который мы назовем квантильным методом [17, 29, 30], позволяет существенно снизить трудоемкость стандартного Алгоритма 0.1 для случая большого или даже счетного числа значений {х,} дискретной случайной величины f, и здесь существенно используется формула (0.3).
Зададим целое число К и введем на полуинтервале [0,1) равномерную сетку {wk = k/K\ k = 0,1,____, /С} (число К, как правило, в два-три раза меньше п). Да-
лее построим массив целых чисел Xk = min{j : Rj — pi +. •. +Pj > u>k-1}, к = 1,..., который называется массивом нижних квантилей. Этот массив задает номер j элемента массива {Я,- ; i = 1,... ,п}, с которого следует начинать поиск ”вверх” (то есть, как и в Алгоритме 0.1, вычитать величины Rq , q = j, j + 1,... из a до получения первого отрицательного значения) при w*-! < а < Wk■ Окончательно моделирование дискретной случайной величины выглядит следующим образом.
АЛГОРИТМ 0.2 [17, 29, 30]. I. Генерируем значение а равномерно распределенной па полуинтервале [0,1) случайной величины.
2. Вычисляем номер к полуинтервала [iOfc_i,tü*), в который попадает а по формуле (0.3):fc= [Ка\ +1.
3. Реализуем последовательный поиск ”снизу вверх” начиная с Rxk-
В связи с введением сетки {к;*} Алгоритм 0.2 можно назвать дискретно-стохастическим. Заметим также, что хотя квантильный метод существенно сокращает вычислительные затраты, он требует дополнительной оперативной памяти ЭВМ.
Во многих ситуациях в данной работе рассматривается следующая конструкция. Пусть в R1 определена неотрицательная функция <р(й) и задано компактное множество
4
ИСК1. Для простоты полагаем, что и ~ выпуклая ограниченная область с границей в И/. Рассмотрим интерполяцию функции основанную на ее разложении по неотрицательным базисным функциям {Хт}> построенным па дискретной сетке {хт}:
<р(й) (й), (0.4)
т
коэффициенты {и>т} выражаются через значения функции у и ее производных в узлах сетки. Сделаем следующее
ПРЕДПОЛОЖЕНИЕ 0.1. Число узлов хт, в которых гЬтхт(й) ф 0 для всех й = (и*1),..., нО)) ^ и, конечно и равно М.
Иереобозначим через {т1?..., хд/} и {гох,... , те узлы и коэффициенты из (0.4), которые удовлетворяют Предположению 0.1. Тогда в правой части (0.4) получим конечную сумму
м
¥■(«) и Ь{м)ч>(й) = £ ш;х.(м). (0.5)
1 = 1
В данной работе в основном будет рассматриваться случай, когда ш,- = <^(т,-), при этом
м
<р(й) « Ь{м)ч>(й) = Х^(50х<(«)- (°-6)
1=1
В качестве основного примера базиса {хт}> удовлетворяющего Предположению 0.1, рассмотрим далее аппроксимацию Стренга-Фикса [26, 27], которая строится следующим образом.
Будем предполагать, что {ят} является равномерной прямоугольной сеткой, то есть каждому узлу хт можно сопоставить мультииндекс = (д[т),... ,<7/Ш^) так, что
хт = (д{ш^1,... ,д[т^0> здесь шаг сетки Н по всем переменным один и тот же (несложно
построить аппроксимацию Стренга-Фикса и для неравномерных сеток, однако при этом нужно использовать более громоздкие обозначения). В этом случае
М ~ /г~' и И ~ ЛГ1'', (0.7)
где знак заменяет слово ” пропорционально”. Рассмотрим базис
Ый) = *,,<->.„Г,,(«(■>..«<'>) = х(х - Ч'п)) х • • • х х(х - 9'(т))’ (°-8>
где х ~ финитная производящая функция [27]. В силу ограниченности множества I] и
финитности функций (0.8) базис {хш} удовлетворяет Предположению 0.1 и условию (0.6). Как и в [27], в качестве производящей функции х возьмем /7-сплайн порядка г:
х = /3(г). (0.9)
Как известно (см., например, [27]), В-сплайн можно определить как (г + 1)-й член рекуррентной последовательности /?(,+1)(и) = 0^ * /9^(п), где
^(0,(«) =
1 для —1/2 < и < 1/2; 0 иначе,
5
а знак ” *” обозначает свертку. Выбор аппроксимации Стренга-Фикса (0.5), (0.6), (0.8), (0.9) объясняется тем, что эта интерполяция обладает рядом свойств, позволяющих эффективно использовать ее как при моделировании случайных величин и случайных векторов (см. далее Задачу 0.1), так и при построении дискретно-стохастических алгоритмов численного интегрирования (см. Главу 2) и глобальной аппроксимации функций, представленных в интегральной форме (см. Главы 3 и 4). Важными частными случаями аппроксимации (0.5), (0.8), (0.9) являются г = 0 - это кусочно-постоянная интерполяции функции V?, а также / = 1, г = 1 - это кусочно-линейная интерполяция функции V? одного переменного. Кроме того, при / > 1 и г = 1 будем называть (0.5), (0.8), (0.9) мулътилинейной интерполяцией функции ср (здесь линейность функции Ь^мур получается только вдоль координатных осей, а в целом линии уровня этой функции получаются ’’гиперболического” типа, так как значения координат в (0.8) перемножаются). Для всех этих случаев приближение функции <р имеет вид (0.6) (то есть можно положить т.- = у?(х;)) [27]. Рассмотрим следующую
ЗАДАЧА 0.1. Построить алгоритм численного моделирования случайного вектора (для I = 1 - случайной величины) Н, имеющего совместную плотность распределения
В Приложении 1 предсталено решение этой задачи (здесь использованы работы автора [17, 18, 31 - 34]).
Отметим три ситуации, в которых возникает необходимость моделирования случайных векторов по плотности (0.10).
Пусть нам известно, что параметр 0 физического процесса V является случайной величиной, реализации которой можно получать с помощью проведения (достаточно дорогих) экспериментов. Пусть также имеется математическая модель IУу, описывающая процесс V с помощью уравнения или алгоритма, причем для получения требуемых величин необходимо численно решать уравнение или реализовывать алгоритм и использовать при этом достаточно много значений случайной величины в (такая ситуация имеет место, например, при решении стохастических задач методами статистического моделирования).
Здесь можно использовать следующий прием. С помощью физических экспериментов получаем значения 0^\ 0$ случайной величины 0 и построим гистограмму
случайной величины в [35]. Для этого разделим весь диапазон наблюденных значений
р$(й) = С £(л/)<^(й), С = 1/с,
(0.10)
на разряды
а = 1,..., М; а = щ < щ < • • • < йм = Ь
и подсчитаем частоты и* = пц/№у (здесь т* - число значений {0^}, приходящихся на г-й разряд). Гистограмму
, м,--! < й < «,•; х = 1,...,М (0.11)
щ - д,-1 Иу (щ - йі-і)
можно переписать в виде (0.10) для I = 1, г = 0, и = [а, 6], полагая
її,- - й,_1 = (6 - а)/М = Л; х,- = (и, + й,-_і)/2; у(х{) =
Иу (щ - Д|'—і)
6
*(«)={ п ПРИ й€Ій;Г’й;)’,и (012)
4 [0 иначе, г = 1,...,Л/. 4 }
Функцию (0.11) можно теперь использовать в качестве моделирующей плотности для получения новых значений 0^\в[*У\... у0^ при численной реализации математической модели при этом нужно решать Задачу 0.1.
Интересен вопрос об оптимальном соотношении параметров ^ и М, определяющих функцию (0.11). Учитывая то, что базисные функции (0.12) являются ортогональными на (о, 6], то есть
(Хп,Хі2) = / ДЙ1 («)»»(“)<й = ( о’
если і! = :2; иначе
(при г > 0 базисные функции (0.8), (0.9) не являются ортогональными), то здесь применим следующий подход, предложенный в [36] (см. также [5, б, 10]).
Пусть 0 - случайная величина, распределенная на конечном промежутке [а, 6] и имеющая на этом промежутке достаточно гладкую плотность распределения р$у и 01,, $цу - ее независимые реализации. Если ф\,...уфм ~ ортогональные и нормированные на [а, Ь] функции:
{Фіі,Фі,)= І Ф:Лй)Фі7(й)<ій =
^а
1, если ?1 = І2; 0 иначе,
то (ре,ф,•) являются коэффициентами Фурье функции ре. Кроме того, (р$уф{) = БФ0)у и тогда функция
М т V
(0-13)
1=1 у ;=1
служит несмещенной оценкой ряда Фурье £?=,(р0,^,)0,(й), то есть при достаточно большом М сумма (0.13) будет служить хорошим приближением в среднем ДЛЯ Рв. Если
(см. (0.12)), то приближение (0.13) вновь представляет собой гистограмму (0.11).
ЛЕММА 0.1 [10, 36]. Если производная сРро/Ии2 ограничена на [а, 6], то уклонение гистограммы (0.11) от графика плотности при числе ступеней М ~ Лу3 о лучшем случае имеет по вероятности порядок Ыу .
Имеется также теория разложения плотности по системе ортонормированных с весом функций и использования приближений типа (0.13). Из этой теории следует, что для гладких плотностей ре целесообразно использовать не гистограммы, а болсс сложные ортогональные разложения (см., например, [10]). Для этого же случая можно использовать неортогональные представления вида (0.10). Например, при для / = 1, г =
1, р = м,
щ - йг-1 = (Ь - а)/М = Л; Х{ = (и,- + «,*-1)/2; х0 = а — Л/2; хм+\ = Ы- /1/2;
^(*0 = V —\; = = 0;
(ц,- -
7
х.'Н =
(ц-ж,-_і)/Л при и Є [х,_ьіі),
(х|+1-й)/^ при й Є [х,-,жі+1), (0.14)
0 иначе
(суммирование в (0.6) производится от г = 0 до г = М + 1) получаем кусочно-линейное представление плотности, которое носит название полигон частот.
Вторая ситуация, в которой требуется решать Задачу 0.1, возникает в случае, когда совместная плотность р§ случайного вектора 0, распределенного в области Я С И1, имеет сложный вид (то есть не удается построить простых моделирующих формул или алгоритмов для компонент 0), и строится приближение (0.10) ДЛЯ ф = р$ И моделирование осуществляется по построенному приближению
р*(й) = Я Ь{м)рр(й);
нормировочная константа Н введена здесь для того, чтобы функция рд осталась плотностью: /Н1 рв{й) с1й = 1. В связи с тем, что метод Монте-Карло используется в основном для оценки линейных функционалов вида [к,р$) = JRl /г(й)р^(й) дй, для оценки близости плотностей р§ и ро можно использовать, как и в классических работах по теории интерполяции плотностей (см., например, [39]), норму пространства 1ч (Я) и исследовать погрешность
^'! = 1п‘
В частности, для / = 1, Я = [а, Ь) справедлива следующяя
ЛЕММЛ 0.2 [37]. Пусть А - множество всех борслсвских множеств интервала
(а, Ь). Тогда
$(^1) = 2 5ир| [ р§(й)с1й— I Ро(й)с1йI.
, *Ча )л '
Л€Л
В качестве следствия Леммы 0.2 можно получить оцепку погрешности вычисления функционалов:
|(/г,рд) — (/г,р^)| < уга1 Бир |/г(х)| (2Бир!/ рд(й)(1й- / рд(й)<1йI), К € Ь^Я).
*€[а,Ь] 4 а'->А ->Л
АСА
Однако, например, для исследования вопроса об изменении дисперсии соответствующего метола Монте-Карло, получаемого при замене плотности р$ на р$, требуются более сильные метрики, чем (например, метрики пространств Ьг(Я) и даже С(Я)) - см. Раздел 2.2. В связи с этим приведем утверждения об аппроксимирующих свойствах приближения Стренга-Фикса (0.5), (0.6), (0.8) в соответствующих функциональных пространствах. Обозначим
^ } Чм)?)= ( £ [„№т Мй) “ ьшМи)))2 ^)1 >
2
здесь
^1С(Г)) = Рсм(<Рі ^(Д/)^) = £ ™Р\Пт{ч>{й) - Ь(А/)У?(Й))|;
т:|т|<гй€У
*га»'
8
где т = (ті,...,га/), пи ■ целые неотрицательные числа, |т| = гаї +... + тп/. Заметим,
найдется набор коэффициентов {и?,-} из (0.5), для которого справедливо неравенство
ЛЕММА 0.4 [26, 27]. Если {хі,... ,хд/} являются точками непрерывности функции
найдется набор коэффициентов {д\} из (0.5), для которого справедливо неравенство
и что для г = 0 и г = 1 указанный в Леммах 0.3 и 0.4 набор коэффициентов {м/,-} представляет собой значения функции р в узлах сетки, то есть выполнено соотношение (0.6) (см. [27], Глава 2, §5).
Чтобы обеспечить условие непрерывности функции <р и ее производных в точках {*!,... у хм} (см. Лемму 0.4), используем следующий результат.
ЛЕММА 0.5 (теорема вложения) [38]. Если Ц - ограниченная область, звездная относительно некоторого шара, / Є \У2Г+1^(£/) и
то существует функция / 6 0(11), совпадающая почти всюду на II с функцией /, такая, что
полнения упомяпутого выше условия Леммы 0.4 разумно сделать
ПРЕДПОЛОЖЕНИЕ 0.2. Для всех рассматриваемых далее функций / из ■уу(г+1)(£/) справедливы соотношения (0.20) и / = /€ С((/).
Потребность решения Задачи 0.1 может также возникнуть при реализациии мажорантного метода исключения [4 - 6, 10, 18, 29].
что *<ь’> = *Г2'0)) и *‘с) =
ЛЕММА 0.3 [26, 27]. Для х <Е и
V € С(г>([/)
(0.15)
(0.17)
(0.18)
Заметим, что
X = Р{г) € \У<г)(К.) и х = Р[г) Є С^(К),
(0.19)
І < 2(г + 1)
(0.20)
где Н - константа, не зависящая от выбора функции /.
Из Леммы 0.5 следует, что при выполнении условия (0.20) функцию / € ^2+1\и) можно отождествлять с эквивалентной ей непрерывной функцией /. Поэтому для вы-
9
Пусть в области IV С К/ задана положительная функция (/|(й>) > 0 при й> £ IV и 1(Ы) = 0 при XV $ И7), причем 0\ = /{(го) дй) < оо (при этом сама константа 61
может быть неизвестной) и требуется построить алгоритм моделирования случайного вектора | с совместной плотностью
р^ги) = /((й)/6г. (0.21)
Если стандартные приемы моделирования случайного вектора £ (например, представление совместной плотности в виде произведения условных плотностей и последовательное моделирование компонент вектора ( методом обратной функции распределения или методом суперпозиции [4 - б, 10, 18, 29]) не реализуемы, то можно использовать следующий прием.
Построим функцию /ир, заданную на той же области ТУ £ Я1 и такую, что /ир(й) > /((хо) > 0 для всех гд в И7, при этом бир = /ир(хЪ)с1й) < оо, и имеются эффективные моделирующие формулы (в дальнейшем - ФОРМУЛЫ (0.22)) реализации случайного вектора г} по совместной плотности р7-,(й>) ~ /иР(й>)/ёир.
АЛГОРИТМ 0.3 [5, 10, 18]. 1. Моделируем вектор щ по формулам (0.22). Реализуем также случайную величину 7 = а /ир(^о) (здесь а по-прежнему стандартное случайное число).
2. Если 7 > /{(чо), то испытание считаем несостоявхиимся и повторяем п. 1 алгоритма, если же 7 < /*(т}0), то принимаем £ = г}0.
Алгоритм 0.3 является частным случаем общей схемы метода исключения.
АЛГОРИТМ 0.4. Пусть случайный вектор ф распределен в некотором множестве X и дано подмножество А С X. Проводим статистическое испытание Т, при этом считаем, что Т состоялось, если ф € А, и Т не состоялось, если ф <£ А.
ЛЕММА 0.6 [4, 18]. Трудоемкость реализации испытания Т (среднее число попыток до получения Т) в Алгоритме 0.4 пропорциональна величине в = 1 /р, где р = Р(ф € А).
Из Леммы 0.6 следует, что трудоемкость Алгоритма 0.3 пропорциональна величине
________1________Оир
Р((чьл)бС) ~ <5 '
Таким образом, мажоранту /ир функции Д следует подбирать так, чтобы площади (7 и Сир были близки; при этом, очевидно, близкими должны быть и сами функции /ир и /£. В качестве /ир целесообразно взять кусочно-постоянные ”окаймления” функции (особенно в многомерном случае при / > 1) и для реализации значений г) (см. пункт 1 Алгоритма 0.3) требуется решать Задачу 0.1.
Рассмотрим также модификацию Алгоритма 0.3, которая может оказаться эффективной для случая, когда вычисление значения функции Д(г)0) является трудоемким, а потому велики затраты на проверку неравенства 7 > /Аг}о) в пункте 2 Алгоритма 0.3.
Построим функцию //ои>> заданную па той же области И7 € И* и такую, что /<(й>) > ЛпДй) > 0 для всех гй £ И7, при этом вычисление значения построенной функции /(ои. является значительно менее трудоемким, чем вычисление значения функции Тогда можно рассмотреть следующую модификацию Алгоритма 0.3.
10
АЛГОРИТМ 0.5 (двусторонний метод исключения) [29, 39, 40, 41]. 1. Как и в Алгоритме 0.3, моделируем оектор т/о по формулам (0.22). Реализуем также случайную величину 7 = а /ир(т/о)-
2. Производим сначала проверку неравенства 7 < fi0w(^o)‘ Если оно верно, то принимаем Ç =Jjq. В противном случае проверяем далее неравенство 7 < если оно
верно, то £ = т/о, иначе испытание считаем песостоявшимся и повторяем пункт 1 алгоритма.
Если функции /ир, и fiow близки, то, с одной стороны, вероятность р = G/GUp близка к единице и среднее число испытаний до получения ” состоявшегося1' испытания близко к единице, а с другой стороны, трудоемкая проверка 7 < h {fjo) будет происходить достаточно редко. В качестве миноранты //ош также целесообразно использовать кусочно-постоянные приближения ’’снизу” функции /<- [40].
Отметим также, что метод исключения имеет аналогию с геометрическим методом Монте-Карло для вычисления интегралов, и здесь также возможно использование двусторонних технологий (см. далее описание Раздела 2.1).
Таким образом, дискретно-стохастические методы (в частности, связанные с решением Задачи 0.1) возникают уже на стадии моделирования случайных величин и случайных векторов.
В Главе 1 исследуются комбинированные алгоритмы, возникающие при реализации траекторий случайных функций.
В Разделе 1.1 представлены результаты по теории слабой (функциональной) сходимости последовательностей случайных функций в пространствах С(Т) (непрерывных на Т С R1 функций с метрикой ’’супремум модуля разности”) и D(T) (функций без разрывов второго рода с обобщенной метрикой Скорохода) [15, 21, 42 - 60], которые использованы далее в Разделах 1.2, 1.3, 3.1 для обоснованного применения численных моделей случайных функций и метода зависимых испытаний. Показано, как можпо преобразовать общие необходимые и достаточные условия слабой сходимости (условие сходимости на алгебре, образующей о-алгебру выборочного пространства, и условие слабой компактности мер) к более жестким достаточным, но легко проверяемым условиям. Так, для пространств С(Т) и D(T) условие сходимости на алгебре сводится к условию сходимости конечномерных распределений (и здесь уже можно использовать предельные теоремы теории случайных функций, в частности, центральную предельную теорему - см., например, [50, 51]). В свою очередь, теорема Арцела [61] и ее аналоги позволяют свести условие слабой компактности мер для пространств С(7’) и D(7’) к условиям малости модуля непрерывности. Дальнейшие упрощения утверждений о слабой сходимости последовательностей случайных полей связаны с переходом к условиям слабой компактности в терминах приращений, а затем - к дифференциальным и мо-ментным условиям. В Разделе 1.1 также рассмотрен вопрос о том, какие функционалы непрерывны в метриках пространств С(Г) и D(T) (в частности, доказана непрерывность функционала F(y) = supf€Tt/(t), у 6 С(Т) V D(T)). Результаты этого раздела опубликованы в работах автора [19, 62 — 70].
В Разделе 1.2 рассмотрены различные виды сходимости численных рандомизированных спектральных моделей гауссовских случайных полей: в среднеквадратическом, слабая, равномерная по вероятности [10 - 12, 15, 21, 71 - 76]. Численные методы построения таких моделей можно трактовать как дискретно-стохастические, так как они основаны на разбиении спектрального пространства на непересекающиеся множества и на построении допредельных интегральных сумм для спектральных представлений мо-
11
делируемых случайных функций. При этом можно достичь совпадения корреляционных характеристик модели и моделируемого поля за счет рандомизации используемых точек спектра [71]. Такие конструкции применяются в основном для моделирования гауссовских случайных полей, так как на практике удается моделировать только независимые слагаемые интегральной суммы спектрального представления, а в этом случае, в силу центральной предельной теоремы для случайных функций [50], предел интегральной суммы представляет собой гауссовскую случайную функцию. В свою очередь, для реализации пегауссовских полей следует моделировать зависимые слагаемые интегральной суммы, что, вообще говоря, непросто. Кроме того, для получения информации о конечномерных распределениях предельного поля требуется доказывать более топкие (по сравнению с центральной) предельные теоремы, которые будут к тому же обладать малой степенью общности. Для получения негауссовских траекторий более перспективным является рассмотрение различных преобразований спектральных моделей гауссовских случайных полей. В Разделе 1.2, в частности, исследованы преобразования, основанные на рассмотрении функций на гауссовских траекториях, а также на комбинировании гауссовских спектральных моделей со случайными величинами [21, 78 - 83]; исследован вопрос о функциональной сходимости моделей, полученных с помощью этих преобразований. Результаты данного раздела опубликованы в работах автора [19, 63 - 70, 84 - 87].
В Разделе 1.3 исследованы вопросы слабой сходимости моделей случайных процессов и полей на точечных потоках Пальма [10, 11, 21, 72, 88]. Эти модели представляют собой кусочно-постоянные аппроксимации некоторых предельных случайных функций на стохастических сетках. Особо отметим, что автором для специального класса корреляционных функций предложена новая численная модель такого типа, более экономичная (по сравнению с ранее предложенными) как в смысле численной реализации, так и простоты обоснования функциональной сходимости. Результаты этого раздела опубликованы в работах [89 - 91].
В Главе 2 исследованы способы решения следующей
ЗАДАЧА 0.2. Построить алгоритм вычисления интеграла
Рассмотрены различные модификации стандартного метода Монте-Карло. Этот алгоритм основан на представлении
и /П1 р^(й)би = 1, а /-мерный случайный вектор 77 распределен по плотности рГ).
I = [ Дй)«Й, и = {й Є Б.' : Дй) ^ 0}. •>и
(0.23)
АЛГОРИТМ 0.0 (Стандартный метод Монте-Карло вычисления многократного интеграла) [1 - 6, 10, 16, 20]. 1. Реализуем п значений случайного вектора г}: г}і,..., уп. 2. Вычисляем приближенное значение интеграла (0.23);
Если в (0.25) трактовать т/ь...,7/п как набор независимых одинаково распределенных случайных векторов с совместной плотностью распределения то случайные величины = Ф(*д),. ..,£п = Ф(*?п) будут также независимыми одинаково распределенными с математическим ожиданием I (см. соотношение (0.24)) и дисперсией
сг2 = УФ(*7) = [ Ф2(й) р5(и) йй - /2 = / <Й - 1\ (0.26)
^и Ju Рч\и)
Если величина (0.26) конечна, то в силу закона больших чисел (см., например, [28, 92]) формула (0.25) верна для достаточно больших п. В силу центральной предельной теоремы (см., например, [28, 92]) для достаточно больших п имеет место соотношение
р(|/-|п| <с(е) А=) = 1 - £, (0.27)
из которого следует, что стандартный Алгоритм 0.6 имеет по вероятности погрешность порядка 1/у/п. Поэтому для приближенного вычисления интеграла с гладкой подынтегральной функцией и небольшой размерности I области I/ метод Монте-Карло заведомо проиграет детерминированным квадратурным и кубатурным формулам [24, 25, 93 - 95]. Однако метод Монте-Карло имеет и ряд несомненных преимуществ. Прежде всего, скорость сходимости метода не зависит от увеличения размерности и возрастания сложности задачи. Поэтому метод Монте-Карло применяется при решении Задачи 0.2 для многомерного случая и (или) ” сложной” подынтегральной функции / (и иногда -для ’’сложной” области интегрирования Ц).
Важным достоинством метода Монте-Карло является также возможность оценки погрешности по результатам вычислений, которую дает соотношение (0.27). Здесь уместно заметить, что если даже величина <т заранее неизвестна, то несложно получить ее приближенное значение, используя ту же выборку т) 1,...,77п, что используется при подсчете математического ожидания Е£ [4 - б, 10, 16]:
\
Исторически некоторое время (в 50-е годы прошлого столетия) активно использовался так называемый геометрический метод Монте-Карло, основанный на построении мажоранты положительной подыптегральпой функции / и розыгрыше равномерно распределенной точки в ’’подграфике” этой мажоранты. В монографии [4] этот метод рассмотрен в качестве альтернативы Алгоритму 0.6. Геометрический метод может быть эффективным (сравнительно с Алгоритмом 0.6) в достаточно экзотическом случае, когда функция Ф“1 существует и ее значения вычисляются значительно быстрее, чем значения Ф.
В Разделе 2.1 показано, что геометрический метод является частным случаем стандартного Алгоритма 0.6, а также установлена аналогия геометрического метода с методами исключения (Алгоритмы 0.3 и 0.5). В частности, использование двусторонних технологий позволяет существенно расширить сферу применимости геометрического метода. Результаты данного раздела опубликованы в [20, 39, 40, 96, 97].
Количество вычислительной работы при реализации Алгоритма 0.6 определяется величиной 5 = п £, где * - среднее время ЭВМ, затрачиваемое на вычисление одного значения Из соотношения (0.27) следует, что при фиксированном е количество
13
испытаний п пропорционально дисперсии V £, поэтому в качестве величины, характеризующей затраты метода Монте-Карло, можно выбрать
S = ta\ (0.28)
Величина S называется трудоемкостью метода Монте-Карло [4 - 6, 10, 16, 20].
Таким образом, плотность следует выбирать так, чтобы минимизировать трудоемкость (0.28).
ЛЕММА 0.7 [4 - 6,10, 20]. Минимальная дисперсия о^ы реализуется о случае, когда плотность р$ пропорциональна модулю подынтегральной функции:
_ „л _ 1/(*)1
) " Ы/(5)|<Й
и равна <r2min = (/<, |/(ü)| dû)2 - I2.
Заметим, что если подынтегральная функция / не меняет знак на U, то = 0.
Выборка т/l,... ,т/п по плотности (0.29) называется существенной. Из Леммы 0.7 можно сделать вывод о том, что в ряде случаев можно добиться уменьшения трудоемкости S Алгоритма 0.6, выбирая плотность р^ близкой к функции (0.29). Алгоритм 0.6 в этом случае называется выборкой по важности [4 - 6, 10, 20], что соответствует английскому термину ”important sampling”. Такое название объясняется тем, что если рп пропорциональна /, то в тех частях области U, в которых |/| больше и вклад которых в интеграл I более существенен, будет выбираться больше случайных точек.
В Разделе 2.2 исследован Алгоритм 0.6, в котором в качестве плотности выбрана нормированиная аппроксимация Стренга-Фикса (0.5), (0.6), (0.8), (0.9) для = С F, где F = |/| и С - нормирующая константа:
C=l/fRiL{M)F(w)dw.
При таком выборе плотности рв первом пункте Алгоритма 0.6 требуется решать Задачу 0.1 для случая / > 1. Получено выражение для условно-оптимального числа узлов, используемых при построении аппроксимации. Предложен ряд модификаций исследуемого алгоритма (в частности, концепция ’’отделения плотности от нуля”). Результаты данного раздела опубликованы в [20, 32 - 34].
В Разделе 2.3 аппроксимация Стренга-Фикса подынтегральной функции использована уже в алгоритме выделения главной части [4 - 6,10, 20]. Показано, что получаемый метод по скорости сходимости и простоте реализации не уступает алгоритму выборки по важности из Раздела 2.2. Предложена также смешанная численная схема, использующая идеи построения алгоритмов выборки по важности и выделения главной части. В этом же разделе дастся обзор возможностей применения ’’детерминированных” технологий в других методах понижения дисперсии стандартного алгоритма метода Монте-Карло. Рассмотрены, в частности, интегрирование по части области, метод условных математических ожиданий, выборка по важности по части переменных, сложная симметризация переменных и метод расслоенной выборки [4 - 6, 10, 20]. Далее в Разделе 2.3 приведен обзор работ, посвященных построению оптимальных (в смысле теории сложности) алгоритмов численного интегрирования [1,3 - 6, 10, 24, 25, 93 - 95, 98 - 106]. Обсуждается, в частности, соотношение между ’’детерминированными” кубатурными
14
(0.29)
формулами [24, 25, 93 - 95] и методами Монте-Карло [4, 6, 10, 24, 25]. Здесь важно отмстить, что явной границы между детерминированными и стохастическими численными алгоритмами интегрирования не существует, и оптимальные алгоритмы этого типа являются по сути дискретно-стохастическими [24, 25, 99 - 106]. К сожалению, последние алгоритмы, как правило, трудно реализуемы на практике. В конце Раздела
2.3 приведены примеры комбинированных дискретно-стохастических алгоритмов вычисления континуальных (в частности, винеровских) интегралов. Результаты данного раздела опубликованы в [20, 34].
В Разделе 2.4 исследован вопрос о возможности использования заданной существенной выборки при решении Задачи 0.2. Предложен алгоритм, основанный на введении функции, интеграл от которой известен. Эта функция должна быть достаточно близкой к подынтегральной, и для ее построения вновь можно использовать конечноэлементные аппроксимации подынтегральной функции. Кроме того, отмечено, что существенную выборку можно получать мажорантным методом исключения. Показано, что в случае, когда существенная выборка не задана, при вычислении интеграла не следует стремиться ее воспроизводить, а лучше использовать стандартный метод Монте-Карло со случайными узлами, близкими (но в точности не совпадающими) по распределению к существенной выборке. Результаты данного раздела опубликованы в работе [39].
Переходя к обзору Глав 3 и 4, отметим, что помимо задач вычисления многократных интегралов "классические” методы статистического моделирования используются при вычислении функционалов от решений интегральных уравнений; при этом используются представления таких функционалов в виде интегралов счетной размерности. Указанные алгоритмы позволяют получать отдельные значения интегралов, в том числе, значения решений интегральных уравнений в выбранпых точках. Если же речь идет о параметрическом интегрировании или об оценке решений интегральных уравнений в целом, то здесь требуются специальные подходы. Одним из первых таких подходов явился так называемый метод зависимых испытаний, предложенный в [107, 108]. В последние годы активно разрабатываются новые функциональные алгоритмы метода Монте-Карло (см. обзоры литературы в [13,14, 22, 23]). В Главах 3 и 4 исследуются комбинированные алгоритмы, которые мы будем называть дискретно-стохастическими (или функциональными) численными методами глобальной аппроксимации функций <р, представленных в интегральной форме. В качестве основных примеров таких функций рассматриваются интеграл, зависящий от параметра
<,£>,(«) = й£ и С И', й€КсГ, (0.30)
и решение <^г(й) интегрального уравнения второго рода
<^2(й) = ^ &(й, й)1р2(й)дй + 0(й), й £ и С V С П/. (0.31)
Предполагается, что существуют достаточно эффективные алгоритмы метода Монте-
Карло для оценки функции <р(й) в одной или нескольких точках й € и. Дискретно-
стохастические методы имеют следующую общую структуру.
АЛГОРИТМ 0.7 По аналогии с (0.4) рассмотрим интерполяцию функции <р, основанную па разложении по базисным функциям {Хт)> построенным на дискретной сетке {хт};
¥>(«) ~ 1^тХт(“)- (°-32)
т
15
Полагаем, что для этого разложения выполнено Предположение 0.1. Переобозиачим через vf} и через {xi»• • • >Хм) и {wi,...,u>m} я* с ?/злы из {xm} и соответ-
ствующие им базисные функции и коэффициенты, которые удовлетворяют этому предположению. Тогда, по аналогии с (0.5), в правой части (0.32) получим конечную сумму
м
¥>(й) я: L(M)V»(ü) = Y. «’iXj(ü)- (0.33)
»=1
Более подробно будет изучаться случай
м
<р(й) и Ь(Щvj(u) = Y. ¥>(*«)Х;(“)- (О-34)
1=1
В этом случае для оценки значений {^(ж,-)} реализуем алгоритмы метода Монте-Карло с числами реализаций {п,};
<p{xi) w Фт{хх) = — £fjf\ (0-35)
j=i
г<?е - j-я реализация - стохастической оценки метода Монте-Карло для значения г = 1,..., М. Окончательно приближение вида (0.34) функции ip имеет
вид
м
V»(w) « L(lU)V{ü) = £6ii(*i)X;(ä). (0.36)
«=1
Для более общего случая (0.33) окончательное приближение имеет вид
м
<р(й) « L(M)<fi(ü) = £ Wi Х;(“)> (0-37)
t=l
где W; явлются соответствующими комбинациями оценок вида (0.35) для значений функции (р в узлах сетки (включающими конечно разностные приближения производных).
В качестве основного примера системы функций {Хт} рассмотрим конечно-элементный базис (0.8) аппроксимации Стренга-Фикса с производящей функцией (0.9). При этом для г = 0 и г = 1 можно использовать приближение (0.36). Выбор аппроксимации Стренга-Фикса определяет некоторое отличие исследуемых здесь численных схем от алгоритмов
из [13, 14, 22, 23] (см. также списки литературы в этих монографиях), в которых ис-
пользуется специальное линейное восполнение. Одной из основных проблем, изучаемых в Главах 3 и 4, будет вопрос о скорости сходимости Алгоритма 0.7 при согласованном увеличении параметров М и п = min(ni,... ,п*/).
Особенностью исследуемых численных процедур является наличие детерминированной и стохастической составляющих погрешности Алгоритма 0.7, что приводит к некоторым сложностям при исследовании сходимости и при оптимизации упомянутых процедур. В частности, возможны различные подходы к выбору вероятностной меры оценки погрешности и критерия оптимальности алгоритма. С точки зрепия теории методов Монте-Карло при исследовании погрешностей комбинированных дискретно-стохастических численных методов глобальной аппроксимации функций ip = ip\ V v?2 наиболее естественным является так называемый L2-подход. Здесь в качестве меры
16
погрешности выбирается Ь2-мсра. Погрешность дискретно-стохастической численной процедуры являтся случайной величиной, и поэтому требуется выбрать вид вероятностной сходимости этой величины к нулю при увеличении числа узлов сетки и числа испытаний при реализации алгоритмов метода Монте-Карло в узлах сетки. Для Ь2-подхода выбирается сходимость в среднем. Тогда стохастическая компонента погрешности определяется дисперсиями оценок значений решения в узлах сетки, которые, в свою очередь являются монте-карловскими критериями эффективности этих оценок. В ряде задач требуется вычислять значения ’’жестких” функционалов от решения 9 (в частности, супремум на компактном множестве). В этом случае в качестве меры погрешности следует выбирать С-меру (супремум модуля разности двух непрерывных функций). В монографиях [13, 14, 22, 23] (см. также обзоры литературы в этих изданиях) для перехода к этой мере предлагается использовать теоремы вложения [38]. При этом нужно изучать сходимость в среднем к нулю производных решения. Кроме того, следует отметить, что верхние границы погрешности, получаемые из теорем вложения, являются достаточно грубыми и реализуемыми лишь для небольших размерностей задач. Более точные верхние границы погрешности в метрике С можно получить, если взять более слабую, чем в среднем, сходимость погрешности к нулю по вероятности. Такой метод исследования погрешности дискретно-стохастических численных алгоритмов мы будем называть С-подходом. Сразу отметим, что С-подход позволяет применить предельные теоремы теории вероятностей при получении верхних границ для стохастической компоненты погрешности Алгоритма 0.7.
При построении верхних границ для детерминированных компонент погрешности используется теория соответствующего восполнения (для восполнения (0.5), (0.6), (0.8) применяются Леммы 0.3, 0.4). При этом важным является вопрос о том, какому функциональному пространству принадлежит исследуемая функция = ч>\ V у>2 - см. соотношения (0.15), (0.17).
Наибольший интерес (и часто - трудность) составляет построение верхних границ для стохастических компонент погрешностей функциональных численных методов для того или иного подхода. Здесь играют роль как свойства соответствующего восполнения, так и особенности используемых стохастических оценок {£^} значений решения в узлах. Важно смещенные они или несмещенные, зависимые или независимые, векторные или скалярпые. Естественно, имеет значение и то, для приближения какой функции - или <?2 ~ используется та или иная комбинированная численная процедура. Для интеграла зависящего от параметра, в качестве стохастических оценок в узлах сетки в данной работе рассматриваются независимые оценки (в каждом узле подбирается своя вероятностная плотность), оценки по методу зависимых испытаний (здесь и плотность, и реализуемые по ней случайные векторы одни и те же для всех узлов), а также смешанные оценки. Предлагается также концепция зависимых оценок с оптимальными плотностями. Все эти оценки являются несмещенными. Для решения интегрального уравнения (0.31) в работе изучаются несмещенные оценки по методу сопряженных блужданий, локальные, векторные и локально-векторные оцепки, а также смещенные оценки по многомерному аналогу метода полигона частот. Здесь уместно заметить, что оценки по методу сопряженных блужданий И векторные оценки ДЛЯ функции 1р2 из (0.31) являются соответственно скалярным и векторным аналогами независимых оценок для функции </?! из (0.30). В свою очередь, локальные и локально-векторные оценки являются скалярным и векторным аналогами оценок по методу зависимых испытаний.
17
Сформулируем следующее
ПРЕДПОЛОЖЕНИЕ 0.3. Коэффициенты {ш,} из (0.33) и их стохастические оценки {И^} из (0.37) таковы, что одновременно выполнены Леммы 0.3, 0.4 и соотношения
Отметим, что в левых частях неравенств (0.38), (0.39) стоят стохастические компоненты погрешностей Алгоритма 0.7 для С- и Ь2-подходов соотвествепно. Соотношения (0.38), (0.39) определяют свойство ’’сноса погрешности в узлы” аппроксимации Стренга-Фикса (0.5), (0.8), (0.9) (это аналог свойства устойчивости приближения функции). Отметим, что для г = 0, г = 1 и приближения (0.36) соотношения (0.38), (0.39) несложно доказать (см. Раздел 3.2). Соотношения (0.38), (0.39), очевидно, выполняются и в случае, когда коэффициенты ги,- из (0.33) выбираются с помощью решения задачи интерполяции Ь\мур{х{) — 4>[х{) (в силу линейности этой задачи), однако при этом могут не выполняться Леммы 0.3 и 0.4 (а значит, и Предположение 0.3). Класс функций, удовлетворяющих Предположению 0.3 при г > 1 требует отдельного изучения.
Для максимума дисперсий из (0.39) и для рассматриваемых функций у>\ и на основании Лемм 0.3, 0.4 удается построить верхние границы. При построении верхних границ для максимума из (0.38) для независимых оценок £0) в узлах (скалярных и векторных) в данной работе используются результаты теории порядковых статистик из [109], а для методов с зависимыми оценками в узлах - теория метода зависимых испытаний.
Для многих используемых в приложениях интегральных уравнений второго рода (в том числе, нелинейных и стохастических) можно с успехом использовать аппроксимацию решения по многомерному аналогу метода полигона частот. Здесь предполагается, что узлы {т,} являются центрами одинаковых малых кубов {{/;}, на которые разбивается область 11. С помощью основной оценки метода Монте-Карло (см., например, [10]) вычисляются функционалы (<^2> А») = /од (^2 (й) /гоев [/{) дй (здесь для всех функционалов используются одни и те же траектории цепи Маркова), которые и принимаются за приближенные значения {^„.(х,)} в узлах сетки {*»■}; затем реализуется восполнение (0.36) или (0.37) с базисными функциями (0.8), (0.9). Соответствующие стохастические оценки {£^} являются смещенными и ’’слабо” зависимыми, причем и смещение, и зависимость уменьшаются с ростом числа узлов. При построении верхних границ для стохастических компонент погрешности этого алгоритма в рамках С-подхода потребовалось доказывать специальные аналоги утверждений теории порядковых статистик для случая зависимых случайных величии, корреляции которых убывают с ростом числа величин.
Отметим, что все условия, обеспечивающие возможность построения верхних границ погрешностей дискретно-стохастических алгоритмов глобальной аппроксимации функций, выражены в терминах подынтегральной функции / (для функции ср1 из (0.30)), ядра к и свободного члена ф интегрального уравнения (для функции из
Особое внимание в работе уделено проблеме выбора параметров функциональных численных методов. Выбору подлежат число узлов сетки А/, вероятностные плотности р, определяющие стохастические оценки в узлах, а также минимум числа реализаций
вир |£(Д#)¥>(«) - £(М)£(Й)| < Сх тах \<р{зц) -
и£(/ 1—1,... ,т
(0.38)
(0.39)
(0.31)).
18
оценок в узлах сетки п. Используется подход из [13, 14]. Фиксируется допустимый уровень погрешности 7, которому приравнивается верхняя граница Х(М, п;р) погрешности Алгоритма 0.7. При этом условии минимизируется функция трудоемкости, аргументами которой являются выбираемые параметры. Практически для всех рассматриваемых в данной работе функциональных алгоритмов такую задачу удалось решить при условии, что вероятностные плотности р выбраны и фиксированы, то есть минимизация трудоемкости производилась по М и п, соответствующие критические значения выбираемых параметров названы далее условно-оптимальными.
В Главе 3 рассмотрены общие подходы к оптимизации дискретно-стохастических численных методов и их реализация для функции (0.30).
В Разделе 3.1 описаны независимые и зависимые стохастические оценки в узлах сетки для интеграла, зависящего от параметра. Рассмотрен вопрос об ослаблении условий сходимости метода зависимых испытаний (по сравнению с работами [107, 108]). Здесь существенно использованы полученные в Разделе 1.1 достаточные условия функциональной сходимости последовательностей случайных функций в пространствах С((/) и 0(6/). Предложены концепции смешанных оценок и зависимых оценок с оптимальными плотностями. Результаты данного раздела опубликованы в [20, 63, 67 - 70, 110-117].
В разделе 3.2 определены Ьг- и С-подходы к исследованию сходимости функциональных алгоритмов. Показано, что для этих подходов погрешность распадается на детерминированную и стохастическую составляющие (для смещенных оценок появляется также компонента смещения). В связи с получением верхних границ для детерминированных компонент погрешностей (с использованием Лемм 0.3, 0.4) получены условия, обеспечивающие соотношения (0.15), (0.17) для функции (р\. Доказаны неравенства (0.38). (0.39) для аппроксимации Стренга-Фикса (0.6), (0.8), (0.9) при г = 0 и г = 1. Эти неравенства позволяют строить верхние границы для стохастических компонент погрешностей и смещения функциональных алгоритмов. В Разделе 3.2, в частности, получены верхние границы стохастических компонент погрешностей функциональных алгоритмов аппроксимации интеграла <р\, зависящего от параметра. Для процедур с независимыми несмещенными оценками {£^} в узлах {я,} для этого использована теория порядковых статистик, а для зависимых оценок - рассмотренная в Разделе 3.1 теория метода зависимых испытаний. Приведены общие утверждения о погрешностях дискретно-стохастических процедур аппроксимации функции из (0.30) для Ь2- и С-подходов. Результаты данного раздела опубликованы в [20, 33, 111 — 125].
В Разделе 3.3 сформулирована задача условной оптимизации. Замечено, что в большинстве случаев верхняя граница Т имеет вид Т(М,п;р) = (М,п;р) +
л;р), при этом порядки возрастания условно-оптимальных параметров при 7 —> 0 можно получить из соотношений Ф\{М, п\р) ~ 7 и Ф2(М,п;р) ~ 7. Минимизация трудоемкости позволяет уточнить константы в выражениях для условно-оптимальных параметров. На основе верхних границ погрешностей, полученных в Разделе 3.2, выписаны выражения для условно-оптимальных параметров функционального Алгоритма 0.7 с зависимыми и независимыми оценками в узлах сетки. Результаты данного раздела опубликованы в [20, 112, 117, 120, 123, 124, 126].
В Разделе 3.4 описан многоуровневый метод зависимых испытаний из [127, 128], оптимальный в смысле теории информационной сложности [101, 102]. Этот алгоритм основан на технологии выделения главной части. Отмечены обстоятельства, ограничивающие возможности применения этого алгоритма на практике. Сформулированы выводы из численных экспериментов подраздела П.3.4. В частности, отмечено, что па-
19
раметры многоуровневого алгоритма, выбираемые на основании простой стоимостной модели, на деле не являются оптимальными. Поэтому требуются более тонкие исследования этого метода на основе весовых стоимостных моделей. Результаты данного раздела опубликованы в [129].
В Главе 4 рассматриваются те же вопросы, что и в Разделах 3.1 - 3.3, но па этот раз - для более сложной функции у?2 из (0.31).
В Разделе 4.1, в частности, описаны несмещенные стохастические оценки в узлах сетки для решения интегрального уравнения второго рода: независимые (полученные по методу сопряженных блужданий) и зависимые (или локальные) оценки, а также их векторные аналоги. Затем рассмотрены вопросы построения верхних границ для детерминированных и стохастических компонент погрешностей рассматриваемых функциональных алгоритмов. При использовании Лемм 0.3 и 0.4 для построения верхних границ детерминированных компонент погрешностей пришлось доказывать утверждения о принадлежности функции <^2 используемым функциональным пространствам в = с<г>(1/) V \у<г+1)(1/) (см. условия (0.15), (0.17)), причем условия этих утверждений удалось выразить в терминах ядра к и свободного члена ф интегрального уравнения (0.31). В свою очередь, при получении верхних границ для стохастических компонент погрешности потребовались аналогичные (то есть выраженные в терминах к и ф) утверждения об ограниченности дисперсий оценок функции (£>2 в узлах сетки. Соответствующие формулировки и доказательства приведены в Разделе 4.1. Результаты этого раздела представлены в [20, 67 - 70, 110 - 113, 115, 117 — 120, 123 - 126, 130 -139].
В Разделе 4.2 получены общие теоремы о погрешностях комбинированных алгоритмов аппроксимации решения интегрального уравнения (0.31) (всего б утверждений; не удалось пока получить верхних границ для погрешностей локально-векторных оценок) и выписаны выражения для условно-оптимальных параметров рассматриваемых дискретно-стохастических процедур. Приведены соображения о выборе плотностей р, определяющих цепи Маркова, используемые для построения стохастических оценок в узлах. Результаты данного раздела опубликованы в [20, 112, 113, 115, 117 - 120, 123 - 126] (см. также работу [140], выполненную под научным руководством автора диссертации).
В Разделе 4.3 рассматриваются вопросы сходимости и условной оптимизации многомерного аналога метода полигона частот. Погрешность этого функционального алгоритма распадается (как для 1/2 -> так и для С-подхода) на три части - помимо детерминированной и стохастической компонент возникает смещение. Предложены подходы к оценке смещения на основе разложения решения в ряд Тейлора. Приведены соображения о нецелесообразности использования гладких восполнений при реализации многомерного аналога метода полигона частот (то есть достаточно использовать базисные функции , (0.8), (0.9) сг = 1 и приближение (0.36)). Получены также верхние границы для дисперсий оценок метода полигона частот, что даст возможность ограничить стохастическую компоненту погрешности для 1/2-подхода. С учетом того, что корреляции оценок в узлах убывают с ростом числа узлов, удалось, по аналогии с рассуждениями из теории экстремумов нормальных последовательностей, получить асимптотические верхние границы стохастической компоненты погрешности метода полигона частот для С-подхода. На этой основе сформулированы окончательные утверждения о погрешностях многомерного аналога метода полигона частот и получены выражения для условно-оптимальных параметров рассматриваемой процедуры. Результаты данного раздела опубликованы в 112, 113, 117, 120, 125] (см. также работы [141, 142], выполненные под научным
20
руководством автора диссертации).
В Главе 5 рассмотрены различные приложения алгоритмов из Глав 1 - 4 и их аналогов. Приведены также результаты тестирования этих алгоритмов.
В Разделе 5.1 рассмотрен ряд вопросов решения стохастических задач математической физики. В частности, построенная в Разделе 1.3 экономичная вычислительная модель случайного поля использована при получении асимптотических (при росте оптической толщины слоя) формул для проходящего излучения через плоские стохастические слои вещества (реализуется методология работ [143, 144]). Здесь уместно заметить, что предположение о стохастической структуре облучаемого слоя является весьма естественным для ряда задач переноса излучения, в которых требуется оценить суммарный поток проходящего излучения [7, И, 15, 21 - 23]. Аналитически и численно для используемой стохастической модели среды показало, что оптически толстый стохастический слой пропускает больше излучения, чем ’’детерминированный” (усредненный). Продемонстрирована возможность эффективного использования метода полигона частот для получения профилей (вдоль толщины слоя) плотностей столкновения и их производных, что позволило определить участки слоя, в которых происходит наиболее быстрое расхождение интенсивностей излучения для стохастической и усреденной сред. В Разделе 5.1 также отмечено, что введение стохастических элементов (случайных коэффициентов, случайных функций и т.п.) часто эффективно осуществляется за счет дополнительной рандомизации вероятностных представлений решений задач математической физики (рассматриваемая в разделе стохастическая задача теории переноса излучения может служить примером такой рандомизации). Это является одной из причин интенсивного развития теории таких вероятностных представлений (см. [9 - 14, 22, 23, 145 - 150], а также списки литературы в этих монографиях). Приведено вероятностное представление решения нелинейного разностного уравнения из [151]. Отмечено также, что в последнее время интенсивно развиваются численные методы, связанные с решением стохастических дифференциальных уравнений [152- 154]. Результаты данного подраздела опубликованы в [16, 90, 91, 117, 121, 151].
В Разделе 5.2 исследуются смешанные дискретно-стохастические алгоритмы решения шггегро-дифференциальных систем из двух уравнений. Первое из них, так называемое ’’полевое” уравнение, содержит неизвестную функцию 0 и интеграл по некоторой неизвестной функции распределения Ф. Для функции Ф, в свою очередь, справедливо соотношение типа уравнения переноса излучения; задаются также соответствующие начальные и граничные условия. Идея построения алгоритмов численного решения таких систем достаточно естественна и проста: ”полевое” уравнение решать детерминированным методом (конечно-разностным, конечных элементов, коллокаций и т.д.), а источники, являющиеся линейными функционалами от решения уравнения переноса, оценивать с помощью метода Монте-Карло. В подразделе 5.2.1 приведены примеры практически важных задач из [155, 156 - 159], для числепного решения которых применяются такие численные схемы. Приведены результаты тестирования смешанного алгоритма рассматриваемого класса на примере решения модельной каталитической задачи из [155]. Предложен подход к построению верхней границы погрешности этого алгоритма. В подразделе 5.2.2 рассмотрен алгоритм ’’конечных элементов - Монте-Карло” для решения одномерной нелинейной задачи радиационно-кондуктивиого теплоперено-са (РКТ) [160 - 167]. Построены различные стохастические оценки конечно-элементных вкладов в источник уравнения теплопроводности (оценки по столкновениям, по поглощениям; локальные, векторные и локально-вскторпыс оценки); проведен их сравнительный анализ и численное тестирование. Для оптически тонких слоев наилучшими оказа-
21
лись оценки по пробегу, а для оптически толстых слоев - локально-векторные оценки. Для случая "слабо линейной” задачи РКТ предложен подход к построению верхней границы погрешности для рассматриваемого комбинированного алгоритма. Результаты данного подраздела опубликованы в [31, 68, 130 - 133, 136, 137, 139, 168].
В Разделе 5.3 рассмотрено применение многомерного аналога метода полигона частот при решении нелинейных интегральных уравнений
V? = G(v?), (0.40)
где ip принадлежит пространству С^(Х), а множество X и ее образ У = <р(Х) совпадают и являются ограниченными выпуклыми областями с границей в R;. Предполагается, что решение этого уравнения может быть найдено методом простой итерации, начиная с некоторого приближения то есть
где = G(^m~^), (0*41)
причем
l|G(y?H) - y>w]|c(Y) < gl<f№ - <рЩс(у) =дыр I0 < 9 < 1 (°*42)
x&r
(здесь и далее индекс в квадратных скобках обозначает номер итерации) и на каждой итерации решается линейное интегральное уравнение второго рода
$sH(*) = J *М(у,*; ^Im-1)) ^Н(у) dy + у*"-'!) (0.43)
методом полигона частот и затем в качестве берется соответствующее приближение
L{M)$mKx)> то есть
G(V[m-4) = 1(М) (О.44)
Алгоритм такого типа впервые предложен в [169] (см. также [13, 170]); в этой же работе построены верхние границы погрешности и выражения для условно-оптимальных параметров в L2- норме. Под научным руководством автора диссертации в [171, 172] Е.В.Шкарупа удалось решить эти же вопросы для более сложного С-подхода. Проведено тестирование алгоритма на примере модельных нелинейных уравнений из [173 - 175].
В Разделе 5.4 приведены результаты тестирования алгоритмов, рассмотренных в Главах 2 - 4. В подразделе 5.4.1 приведены соображения о целесообразности использования траекторий спектральных моделей случайных полей из Раздела 1.2 в качестве подынтегральных функций при тестировании алгоритмов численного (в том числе, параметрического) интегрирования. Такой выбор подынтегральных функций дает воз-I можность соблюсти "независимость” тестирования, получить (хотя бы в простейших случаях) точные значения решения и порядки верхних границ оценки погрешности, контролировать гладкость и трудоемкость вычисления подынтегральной функции. Результаты данного подраздела опубликованы в [87]. В подразделе 5.4.2 исследуются адаптивные алгоритмы численного интегрирования из Разделов 2.2, 2.3 (соответствующие результаты опубликованы в [32, 176]). В подразделе 5.4.3 предложены алгоритмы вычисления констант, имеющихся в выражениях для верхних границ погрешностей и условно-оптимальных параметров исследуемых в Главах 3 и 4 численных процедур [117]. В подразделе 5.4.4 представлены результаты тестирования алгоритмов аппроксимации функций (0.30) из [111]. В подразделе 5.4.5 приведены результаты расчетов по многоуровневому алгоритму, позволившие сформулировать основные выводы Раздела 3.4 [129]. Наконец, в подразделе 5.4.6 рассмотрено специальное тестовое
22