Содержание к диссертации
Введение
Глава 1. Вероятностное представление для решения метаэллиптического уравнения 18
1.1. Основные определения 18
1.2. Вероятностное представление для решения первой краевой задачи для уравнения (L + С)Р+1ІІ д 19
1.3. Дискретный метод Эйлера для решения бигармонического уравнения ААи — д .' 22
Глава 2. Алгоритмы "блуждания по решётке" 26
2.1. Переход от дифференциальных уравнений к конечно-разностным 26
2.2. Оценка решения метагармонического уравнения 27
2.3. Оценка решения уравнения 31
2.4. Глобальная оценка решения уравнения и задача минимизации трудоёмкости 38
2.5. Оценки решений уравнения, уравнения со слабой нелинейностью и задач со смешанными краевыми условиями, включая условие Неймана 41
2.6. Оценки решения метагармонического уравнения и первого собственного числа многомерного оператора Лапласа 44
2.7. Численные результаты 47
2.7.1. Метод, основанный на "блуждании по решётке" 47
2.7.2. Сравнение алгоритмов,' соответствующих "блужданию по решётке" и дискретному методу Эйлера 49
Глава 3. Алгоритмы "блуждания по сферам" 51
3.1. Основные обозначения и определения 51
3.2. Оценки решения метагармонического уравнения 53
3.3. Оценки для решения бигармонического уравнения при п = 2, 3 58
3.4. Оценка ковариационной функции решения уравнения А?и ~ д со случайными параметрами при п ~ 2 59
3.5. Оценка ковариационной функции решения уравнения А2и -\- си = д со случайными параметрами 64
3.6. Численные примеры 69
3.6.1. Нахождение решения трёхмерного бигармонического уравнения 69
3.6.2. Сравнение двух методов для вычисления ко вариации бигармонического уравнения 70
3.6.3. Реализация двух алгоритмов для нахождения ковариационной функции решения уравнения А2и + си = д .71
Заключение 78
Список литературы 79
- Вероятностное представление для решения первой краевой задачи для уравнения (L + С)Р+1ІІ д
- Глобальная оценка решения уравнения и задача минимизации трудоёмкости
- Оценки решения метагармонического уравнения и первого собственного числа многомерного оператора Лапласа
- Оценки для решения бигармонического уравнения при п = 2, 3
Введение к работе
Официальной датой рождения метода Монте-Карло принято считать 1949 год, когда была опубликована статья [46] С. Улама и Н. Мет-рололиса. Сам термин был предложен еще во время Второй мировой войны выдающимися учеными XX века математиком Дж. фон Нейманом и физиком Энрико Ферми в Лос-Аламосе (США) в процессе работ по ядерной тематике. Хотя методы Монте-Карло были известны и до 40-х годов, интенсивное развитие статистическое моделирование получило несколько позже в связи с появлением компьютеров, что позволило проводить вычисления больших объемов. С другой стороны, всё большее распространение получает статистическое описание тех или иных сложных физических процессов в связи с чем методы Монте-Карло всё более активно используются во многих научных областях (теория переноса, теория массового обслуживания, теория надежности, статистическая физика и др.). -
Основными преимуществами данных методов являются физическая наглядность и простота реализации, малая зависимость трудоемкости задачи от её размерности, возможность решения задач со сложной геометрией, оценивание отдельных функционалов от решения без запоминания значений решения во всей области, вероятностные представления позволяют строить обобщённые решения уравнений, одновременное оценивание вероятностной погрешности оценки искомого функционала, простое распараллеливание методов.
Одна из схем решения краевых задач методом Монте-Карло заключается в сведении исходной дифференциальной задачи к некоторому интеграл ьно'му уравнению, что даёт возможность использовать развитой аппарат методов Монте-Карло [10, 11] для решения интегральных уравнений второго рода. Среди подходов такого рода выделим следующие два. .
Первый подход связан с использованием различных теорем о среднем и формул Грина для стандартных областей, содержащихся в исходной области, например, для шара, сферы, эллипсоида и т.д. При этом локальное уравнение записывается на само решение исходной задачи, а ядро этого уравнения является, как правило, обобщенным. Решение локальных интегральных уравнений с обобщенными ядрами практически невозможно традиционными численными методами, но такой вид ядер допускает естественную реализацию методом статистического моделирования. Необходимо отметить, что алгоритмы статистического моделирования строятся в рамках предположения о сходимости соответствующего ряда Неймана.
Второй подход включает методы, основанные на использовании не локальных, как в первом подходе, а глобальных интегральных уравнений. То-есть интегральное уравнение рассматривается либо на границе, либо по области, либо одновременно по области и границе, причём уравнение может быть записано как на само решение исходной задачи, так и на некоторую вспомогательную функцию. Последний случай реализуется, например, при построении случайного блуждания по границе на основе граничных интегральных уравнений теории потенциала [33]. Данное направление потребовало разработки алгоритмов метода Монте-Карло для решения интегральных уравнений, спектральный радиус которых не обязательно меньше единицы.
Другой способ состоит в использовании тесной связи эллиптических операторов с диффузионными процессами (см., например, [5, 40]) и основан на вероятностном представлении решения исходной дифференциальной задачи в виде функционала от решения соответствующей системы .стохастических дифференциальных уравнений (СДУ). В качестве вспомогательной задачи в диссертации рассматривается n-мерная задача Ди- рихле для эллиптического уравнения lm + cu+g = -^aij(r)^^+^bi(r)—+cu+g = 0, и^=<р (1) в ограниченной области D с границей, которая предполагается односвяз-ной и регулярной. Дифференциальному оператору L из (1) можно поставить в соответствие систему СДУ в смысле Ито [38]: 'b = b + jb{Qn,cr(r) - нижняя треугольная матрица диффузии, определяемая разложением Холесского [13]:
А~а{г)- ат{г), a w(t) - стандартный винеровский процесс.
Настоящая работа посвящена построению, обоснованию и применению методов Монте-Карло в рамках последнего подхода с целью приближенного решения метагармонического уравнения, являющегося част-ным случаем мегпаэллиптического уравнения {L + с)р+1и = ~9t (L + c)ku\T =
Впервые эвристическое описание алгоритма "блуждания по сферам" было дано Брауном [4] для приближенного решения задачи Дирихле для уравнения Лапласа. Строгое обоснование алгоритма было предложено Мюллером в [45], где изучался также вопрос об эффективности такого алгоритма. Дальнейшее развитие алгоритмы "блуждания по сферам" получили в работах Михайлова Г.А. и Елепова Б.С. В работе [9] был предложен алгоритм для решения уравнения Гельмгольца, использующий вероятностное представление решения. Позднее в [8] был построен алгоритм для решения этой задачи, исходя из специального интегрального уравнения второго рода. Данный подход позволил разработать алгоритмы статистического моделирования для решения более общих уравнений и систем (см., например, [12, 32, 33, 41]). Алгоритм "блуждания по сферам" основан на моделировании точек последовательного выхода винеровского процесса из максимальных сфер, целиком лежащих в рассматриваемой области. В работах [7, 8, 12, 33] для широкого класса границ Г была получена логарифмическая оценка для среднего числа переходов q(ro, є) в цепи "блуждания по сферам" до момента попадания в е-окрестность границы: -Eq(r0,e)
Кроме того, для больших значений п размерности пространства верна [45] линейная зависимость Сп ~ п.
В первой главе диссертации путём итерирования вероятностного представления решения задачи Дирихле для эллиптического уравнения (1) построено вероятностное представление для решения метаэллиптиче-ского уравнения (3). Это даёт возможность, в рамках третьего подхода, применить для нахождения функционалов от решения хорошо развитую теорию оценок метода Монте-Карло, использующих сферический процесс. Разработанные алгоритмы "блуждания по сферам" для решения как детерминированных, так и стохастических задач приводятся в третьей главе диссертации. Следует отметить, что в рамках первого и второго подходов векторные оценки для решения метагармонического уравнения были получены в [33]. Однако в этой работе было рассмотрено только однородное уравнение. Кроме того, полученные с помощью "двойной рандомизации" [41] оценки "блуждания по малым сферам" для ковариационной функции решения задачи А2и + ku = f со случайными параметрами имеют существенно большую дисперсию, сравнительно с новыми оценками из п. 3.4 и п. 3.5, полученными частичным осреднением скалярных оценок.
Для приближенного статистического моделирования решения системы СДУ (2) можно, в частности, использовать простой в реализации и физически наглядный дискретный метод Эйлера с постоянным шагом At [35]: + Ь(гя) At + ff (rjifo л/ДЇ, (4) где rn - численная оценка решения (2) в узлах равномерной сетки по времени {n-At}, а {г}п} - последовательность независимых между собой случайных векторов с независимыми стандартными гауссовскими компонентами. Однако из-за невысокой скорости сходимости метода (4) приходится использовать методы более высоких порядков с более сложной реализацией (см., например, [43, 35]). В диссертации предложен дискретный метод Эйлера для решения б и гармонического уравнения, который естественным образом распространяется для решения метаэллиптического уравнения.
Во второй главе диссертации использован ещё один подход к решению краевых задач, широко применяемый в практических расчётах в силу своей простоты и сравнительной универсальности. Речь идёт о методе "блуждания по решётке", которое может интерпретироваться как дискретный вариант упомянутых выше подходов. Этот подход прост в реализации, не требует большой памяти ЭВМ, но является сравнительно трудоёмким из-за необходимости моделирования длинных случайных траекторий [32, 34]. Для глобальных оценок решения рассматриваемых задач приходится решать задачу выбора оптимального распределения начальной точки траектории, минимизирующую трудоёмкость метода.
Величиной, определяющей эффективность алгоритмов метода Монте-Карло при их практическом использовании, является трудоёмкость 5(C) случайной оценки ( (см., например, [11])- В диссертационной работе под трудоёмкостью S(Q будем понимать среднее количество вычислительной работы, необходимой для достижения заданной погрешности:- S(0 = DC-^(0, где t(Q - среднее время, затрачиваемое на моделирование одного выборочного значения случайной величины С (например, для процесса "блуждания по сферам" величина t() определяется средним числом переходов в цепи до момента попадания в є - окрестность границы Г). Поскольку общая погрешность оценки, полученной методом Монте-Карло, обычно состоит из двух частей: детерминированной и вероятностной, целесообразно выбирать параметры алгоритма таким образом, чтобы обе эти погрешности имели один и тот же порядок. Вероятностная погрешность возникает в связи с заменой математического ожидания Е ариф-
1 N метическим средним выборочных значений -— у Q и пропорциональна величине — ) где N - количество моделируемых выборочных значений. Детерминированная погрешность появляется, например, при введении є - окрестности границы (в алгоритмах "блуждания по сферам"), при замене диффузионного процесса кусочно-линейным процессом с шагом по времени At (в методе Эйлера) или при замене исходной дифференциальной задачи разностной (в алгоритмах "блуждания по решётке").
В диссертационной работе предполагаются выполненными все условия регулярности данных задачи (коэффициенты, граничная поверхность области, краевые условия и т.п.) для существования и единственности некоторого обобщенного решения задачи (1), а также его вероятностного и интегрального представление с помощью шаровой функции Грина.
Далее следует краткое содержание диссертационной работы по главам и параграфам.
В первой главе исходя из стандартного вероятностного представления решения задачи Дирихле для уравнения эллиптического построено вероятностное представление решения первой краевой задачи для ме-таэллиптического уравнения. Это представление показывает, что рас- сматриваемое решение может быть получено путём параметрического дифференцирования решения специальной краевой задачи для эллиптического уравнения . Таким образом, открывается возможность построения и обоснования требуемых алгоритмов метода Монте-Карло путём параметрического дифференцирования известных стандартных оценок, связанных с "блужданием по решётке" и "блужданием по сферам". Результаты главы 1 опубликованы в [17, 20, 24].
В параграфе 1.1 вводится вспомогательная задача для эллиптического уравнения в области DcR" (L + с 4- Х)и ~ —д, и T = v, формулируются условия, обеспечивающие существование и единственность решения данной задачи; выписывается конкретный вид вероятностного представления для этого решения.
В параграфе 1,2 рассматривается краевая задача для метаэллипти-ческого уравнения (I + с)р+1и ^д, (L + с)ки = >*, fc = 0,...,p.
Посредством итерирования вероятностного представления для решения эллиптического уравнения строится вероятностное представление решения и = ЕС, = Е - Р' 1=0 для предварительно расщеплённой задачи настоящего параграфа. На основе данного представления показывается, что решение метаэллипти-ческого уравнения может быть получено путём параметрического дифференцирования стандартного решения специальной задачи для эллиптического уравнения. Доказывается конечность дисперсии случайной величины Ср-
В параграфе 1.3 для полученного вероятностного представления реализована схема Эйлера с постоянным шагом по времени с целью нахождения решения бигармонического уравнения ЛДи = д. В конце параграфа 1.3 приведены результаты численных расчётов и сделано замечание об улучшении порядка погрешности в методе Эйлера.
Во второй главе построены и обоснованы новые весовые методы Монте-Карло для оценки решения задачи Дирихле для многомерного разностного бигармонического уравнения на основе моделирования "блуждания по решетке". Векторные варианты построенных алгоритмов непосредственно распространяются на разностные метагармонические уравнения с сохранением вида условий несмещенности оценок и ограниченности их дисперсий. В связи с этим построен простой алгоритм для оценки первого собственного числа многомерного разностного оператора Лапласа. Кроме того, построены специальные алгоритмы "блуждания по решетке", позволяющие при определенных условиях оценивать решения задачи Дирихле для бигармонического уравнения со слабой нелинейностью и для задач со смешанными краевыми условиями, включающими условие Неймана. Результаты главы 2 опубликованы в [16, 18, 20, 22, 23].
В параграфе 2.1 формулируется разностная задача, аппроксимирующая вспомогательное уравнение Гельгольца; вводятся соответствующие обозначения; определяется случайный процесс "блуждания по решётке".
В параграфе 2.2 с помощью приема описанного в параграфе 1.2 строится оценка решения разностного метагарионического уравнения (Ал + сл)р+1іі'* — gh. Доказывается теорема о равномерной ограниченности полученной оценки.
В параграфе 2.3 рассматривается задача вида: (Д + с)(Д+ &)« = -#, Аи + Ъи и эквивалентная ей система уравнений = ср и Г г = *
Аи + Ьи = v, Av + си = —з,
В соответствующую оценку решения разностной задачи для первого уравнения данной системы входят неизвестные значения г, но можно использовать их случайные оценки на этой же траектории. Таким образом, строится оценка для решения исходной задачи. На основе векторного представления доказываются конечность дисперсии данной оценки. Из рекуррентного представления оценки выводится равномерная по h ограниченность дисперсии.
В параграфе 2.4 строится глобальная оценка решения задачи из параграфа 2.3. Применительно к таким задачам известно несколько подходов к моделированию цепи [32, 30]. В этом параграфе было выбрано равномерное по сеточной границе области начальное распределение цепи, для которого среднее число шагов есть величина порядка О (Л-1). Для глобальной оценки, которая получается линейным восполнением оценок решения в узлах равномерной сетки, решается задача минимизации трудоёмкости.
В параграфе 2.5 рассматривается задача вида:
ДА« + сАи + Ьи = —ду = <>.
С помощью рандомизации соответствующей разностной задачи получается скалярная оценка решения. Доказывается конечность дисперсии представленной оценки.
Также из рекуррентного представления полученной оценки с помощью ветвления строится оценка решения аналогичной задачи для слабо нелинейного уравнения:
Л Аи 4- сАи + Ьит - -д.
В конце параграфа рассматривается рассматривается модификация схемы блуждания, сформулированной в параграфе 2.1, позволяющая строить оценки решения для задач со смешанными краевыми условиями, включая условие Неймана.
В параграфе 2.6 оценка "блуждания по решётке", построенная в п. 2.3, распространяется на следующую задачу Дирихле для метагар-монического уравнения: (Л + с^)(А + с<2>)... (Л + с^)и = -д,
П(Д + с<*))« =^k'l\ Ar = 2,...,m,
1=к г с комплексными коэффициентами с^к\ При этом обосновывается конечность и равномерная ограниченность по h дисперсии полученной оценки. На основе соотношения (при с = const) * ch - с> где uh - решение задачи Ahuh + cuh - О, изучается метод Монте-Карло для оценивания с* по формуле р[1 - ch?/{2n)] XN+p-l)[h2/{2n)l Здесь —с* - первое собственное число многомерного оператора Лапласа для области D.
В параграфе 2.7 приводятся результаты численных расчётов. В пп. 2.7.1 с помощью оценки, полученной в п. 2.3 численно решаются би-гармоническое и метагармоническое уравнения. В пп. 2.7.2 проводится сравнение двух алгоритмов, соответствующих "блужданию по решётке" и методу Эйлера.
В третьей главе при построении оценок параметрических производных в случае "блуждания по сферам" преодолены существенные трудности. Здесь дополнительно необходимо изучать смещение "- оценок" и обосновывать равномерную ограниченность их дисперсий. Это сделано в настоящей работе на основе использования "мартингального свойства" стандартного вероятностного представления, которое состоит в следующем: представление сохраняет свойство несмещенности, если траекторию процесса оборвать в некоторый момент времени, добавив к оценке значение стандартного решения с соответствующим весом. Таким образом построены и обоснованы оценки "блуждания по сферам" для решения первой краевой задачи для метагармонического уравнения.
В случае бигарионического уравнения- А2и = —д полученные оценки записаны в явном виде. Они представляют собой линейные комбинации значений функции д и граничных функций. Это дало возможность построения просто реализуемых алгоритмов оценки ковариационных функций решений уравнений Д2и — д> А2и + си — д со случайными функциональными параметрами. Показано, что неограниченность числа шагов в случайном блуждании не является существенной. Результаты главы 3 опубликованы в [17,19-21,24]
В параграфе 3.1 вводятся обозначения и определения, связанные с процессом "блуждания по сферам".
В параграфе 3.2 на основе полученного в п. 1.2 вероятностного представления строятся оценки решения метагармонического уравнения, связанные с процессом "блуждания по сферам". Доказывается конечность дисперсий и обосновывается є-смещенность построенных оценок.
В параграфе 3.3 выписывается конкретный вид оценок решения би-гарионического уравнения для случая двумерной и трёхмерной области D.
В параграфе 3.4 рассматривается задача Дирихле для бигармониче-ского уравнения
АЛи — —5) Аи Vi> и со случайными функциональными параметрами д, <^о, Уь Путём условного (для фиксированных траекторий) осреднения произведения соот- ветствующих оценок решения в точках г и п1' строится оценка ковариационной функции v(r,г^) решения задачи, записанной выше. Для реализации полученного метода предлагается при превышении количества точек {гі} процесса "блуждания по окружностям" некоторого достаточно большого уровня М не сохранять информацию о следующих точках с номерами і — М+1,..., N — 1, но, при попадании в є-окрестность границы, использовать полученные вес Qn и координаты г^. При этом Qi = О для і = М + 1,..,, JV—1, и можно заменить Qm на Qu{N—М+1). Показывается, что порядок (по е) величины М можно эвристически оценить по формуле
М х|1пе|2.
В параграфе 3.5 рассматривается задача вида:
Аи + си~д, и\ — (ро, Ди = 9?ъ со случайными функциональными параметрами gt v?o, >ь не зависящими от случайной константы с такой, что
Ес = 0, Dc Решение и(г,с) аппроксимируется суммой двух первых членов ряда Лорана, по формуле и (г, с) и и(г) 0) -|- и^ (г, 0)с, с погрешностью приблизительно равной ^г/2'(г, 0)с2. Выписывается конкретный вид оценок параметрических производных и^(г, 0), уг-'{г} 0). Путём условного осреднения (для фиксированных траекторий г і, г\ и значений случайных параметров д, (ро, (pi) произведения полученных оценок решения в точках Ті и г\ строится оценка ковариационной функции v(ri,r\ ) решения исходной задачи этого параграфа. Оценивается погрешность предложенного метода Монте-Карло. В параграфе 3.6 приводятся результаты числовых расчётов для оценок, полученных в третьей главе. В пп. 3.6.1 решается неоднородное бигармоническое уравнение с помощью оценки, выписанной в п. 3.3. В пп. 3.6.2 тестируются два алгоритма, предложенных в п. 3.4. В алгоритме А, после того как количество центров сфер из соответствующей оценки достигало уровня М, случайное "блуждание по кругам", исходящее из точки г, прекращалось. В алгоритме В в подобном случае процесс не обрывался, но соответствующие веса Qi и случайные точки рі не запоминались, а последний вес Qjv0 домножался на величину (JV — М + 1). Численно показывается, что эвристическая оценка для М (М х |гп|2) является несколько завышенной. В пп. 3.6.3 сравниваются алгоритм, полученный в п. 3.5, и рандомизированный алгоритм "блуждания по окружностям и кругам" [27]. Результаты расчётов подтверждают, что предложенный в настоящей работе алгоритм имеет значительно меньшую (приблизительно в 5 раз) дисперсию, На защиту выносятся следующие основные результаты. Используя вероятностное представление (1.3) для решения второй задачи из (1.6) и теорему Фубини, перепишем соотношение (1.7) в следующем виде: Дифференцируя равенство (1.3), можно вносить производные по параметру А под знак математического ожидания (см., например, [14]), так как в интервале А Л Ао имеет место следующее мажорирование Соответствующее этим мажорантам математическое ожидание конечно, так как оно представляет собой решение задачи Дирихле со следующими значениями функциональных параметров: Предполагая справедливость перестановки оператора L с производной по параметру А, получаем следующее соотношение для р-кратной параметрической производной от решения эллиптического уравнения (1.1): Следовательно, можно предполагать, что р-кратная производная по параметру А при А = 0 от решения уравнения (1.1) удовлетворяет следующей краевой задаче для метаэллиптического уравнения Теорема 1. В условиях леммы 2 параметрическая производная р-ой степени от решения задачи (1-1) с функциональными параметрами Доказательство. Для заданных равенствами (1-Ю) функций (р, д зада ча (1.9) совпадает с задачей (1.4) и, следовательно, выражение (1.5) сов падает с выражением, находящимся под знаком математического ожи дания в (1.8). Это доказывает справедливость утверждения теоремы и тем самым допустимость перестановки оператора L с производной по параметру А. Теорема 2. Если с CQ с /2, то дисперсии случайных величин Qp конечны. Здесь использовано мажорирование решением уравнения Гельм гольца со следующими со следующими значениями функциональных параметров: Один из подходов построения оценки математического ожидания функционалов от траекторий случайных процессов вида (1.11) состоит в следующем [6]. Используя метод Эйлера [25] с постоянным шагом по времени At для численного решения СДУ (1.2), статистически моделируются приближения гп к траектории случайного процесса & в моменты времени п -At: где {%} - последовательность независимых между собой случайных векторов с независимыми стандартными гауссовскими компонентами. Далее с помощью различных численных методов приближённо находится значение функционала к от кусочно-линейной траектории д$, полученной по значениям {гп}. Соответствующая оценка имеет следующий вид: Согласно рассмотренной выше схеме Эйлера построим на основе вероятностного представления (1.4) оценку решения для би гармони чес кого уравнения со следующими значениями входных параметров в области D = [0; 1] X [0; 1] х [0; 1]. Точное решение задачи имеет вид и — exeyez. Координаты г, г, г приближения д4 к винеровскому процессу & моделируются по формулам [11]: где а\, с 2, «зі а4 - независимые случайные величины равномерно распределённые в (0; 1). Значение функционала h(t,t) аппроксимируется следующим образом [1]: где в - момент времени выхода моделируемой траектории на границу области. Результаты расчётов приведены в таблице 1. В таблице 1 приняты обозначения: г - координаты точки, At шаг по времени, JV количество моделируемых траекторий, а7 - дисперсия случайной оценки разности и(г) — и(г). Замечание 1. Известно, что без учета выхода процесса на границу верно неравенство (см., например, [35]): то есть порядок погрешности дискретного метода Эйлера равен 1, а при появлении границы он уменьшается до [At)q — (At)1/2 [39]: Особо отметим, что в рассмотренных оценках порядка погрешности функционал ft(&) зависит только от траектории & и не зависит от времени Приведённые численные оценки показывают, что порядок О {(At)1 2) погрешности схемы Эйлера для решения краевых задач, по-видимому, имеет место и для оценок вида (1.14). Это можно объяснить детерминированностью дополнительного временного множителя. Для увеличения порядка погрешности q можно перейти к непрерывному методу Эйлера в котором учитывается вероятность невыхода траектории t из области за время At. Другой способ увеличения порядка погрешности - использование блуждания по малым сферам [43] Здесь п - размерность пространства, vn - последовательность независимых равномерно распределённых по поверхности единичной сферы случайных векторов, поэтому приращения процесса ограничены, тогда как в дискретном методе Эйлера приращения процесса на произвольном временном интервале могут принимать сколь угодно большие значения. Глава 2. Алгоритмы блуждания по решётке" 2.1. Переход от дифференциальных уравнений к конечно-разностным В вычислительной математике для нахождения приближенного решения классическим подходом является замена дифференциальных уравнений соответствующей разностной задачей. Полученную систему линейных уравнений возможно решить методом Монте-Карло, используя случайные "блуждания по решётке", после приведения её к специальному виду: где р{А) - спектральный радиус матрицы А (см. [И, 37]). Воспользуемся данным подходом применительно к рассматриваемым ниже задачам. В этой главе допускается, что функции с, д, р и и могут быть комплексными, причём М = Re(c). В области D строится равномерная сетка с шагом Лив качестве оценки решения задачи (1.1) для L = А в узлах сетки г = (ііh,..., inh) рассматривается решение разностной задачи Для глобальной оценки решения задачи (2.11) в области D, можно оценить значения решения системы (2.13), (2.14) в узлах сетки и построить линейное восполнение й(г) [16]. При этом трудоемкость метода определяется величиной где t - время моделирования одного перехода, N - число моделируемых траекторий, необходимых для достижения заданной погрешности, Е1 - средняя длина траекторий (асимптотически линейно по п растущая функция). Для оценки решения в разных точках можно использовать одни и те же траектории, полагая, что траектория началась в состоянии г, если она попала в него впервые на каком-то шаге блуждания. Для глобальной оценки решения систем уравнения (2.13) и (2.14) можно предложить следующий алгоритм [23]. Начальное состояние цепи Маркова выбирается из области Dft U 1\ согласно заданному начальному распределению с вероятностями При моделировании процесса "блуждания по решётке" запоминаем узлы сетки, через которые проходит траектория. Цепь обрывается, когда траектория попадает на границу области. На обратном проходе вдоль траектории, осуществляем суммирование с соответствующими весами. Глобальная оценка решения имеет вид 1-1 I где среднее берется по траекториям, прошедшим через i-ый узел, тгц -момент первого в него попадания. Справедливы следующие леммы. Лемма 4. [32] Если 7Г» = тг = 1/L, т.е. начальное состояние цепи "блуждания по решётке" выбирается равновероятно по области D/t U Гл, Лемма 5. [32J Пусть начальное состояние цепи Маркова выбирается равновероятно из множества точек границы, т.е. 7Г = 0 в D , и щ = L\ = C hn на Гд. Тогда среднее число посещений заданного внутреннего узла равно C$hn l т.е. одинаково для всех узлов, и средняя длина траектории Поэтому предлагается следующий алгоритм [23]. Начальное состояние цепи Маркова выбирается равномерно по границе, затем с вероятностью единица траектория попадает в ближайший внутренний узел, после чего с одинаковой вероятностью она может перейти в любой идлина траекторий (асимптотически линейно по п растущая функция). Для оценки решения в разных точках можно использовать одни и те же траектории, полагая, что траектория началась в состоянии г, если она попала в него впервые на каком-то шаге блуждания. Для глобальной оценки решения систем уравнения (2.13) и (2.14) можно предложить следующий алгоритм [23]. Начальное состояние цепи Маркова выбирается из области Dft U 1\ согласно заданному начальному распределению с вероятностями При моделировании процесса "блуждания по решётке" запоминаем узлы сетки, через которые проходит траектория. Цепь обрывается, когда траектория попадает на границу области. На обратном проходе вдоль траектории, осуществляем суммирование с соответствующими весами. Глобальная оценка решения имеет вид 1-1 I где среднее берется по траекториям, прошедшим через i-ый узел, тгц -момент первого в него попадания. Справедливы следующие леммы. Лемма 4. [32] Если 7Г» = тг = 1/L, т.е. начальное состояние цепи "блуждания по решётке" выбирается равновероятно по области D/t U Гл, Лемма 5. [32J Пусть начальное состояние цепи Маркова выбирается равновероятно из множества точек границы, т.е. 7Г = 0 в D , и щ = L\ = C hn на Гд. Тогда среднее число посещений заданного внутреннего узла равно C$hn l т.е. одинаково для всех узлов, и средняя длина траектории Поэтому предлагается следующий алгоритм [23]. Начальное состояние цепи Маркова выбирается равномерно по границе, затем с вероятностью единица траектория попадает в ближайший внутренний узел, после чего с одинаковой вероятностью она может перейти в любой из соседних узлов и т.д.. Цепь обрывается, когда траектория вновь окажется на границез соседних узлов и т.д.. Цепь обрывается, когда траектория вновь окажется на границе. Известно, что, если и{т) Є C4(D), то верно неравенство и, следовательно, где й(г) - линейное восполнение оценки разностного решения по N траекториям. Отсюда в случае вещественных функций и и й получается [41, 32] следующая оценка погрешности в метрике пространства 2 (D): Достаточно ясным способом, рассмотренная в п. 2.3, векторная оценка распространяется на следующую задачу Дирихле для метагармонического уравнения: с комплексными коэффициентами с -. При этом к ая строка блочной матрицы системы типа (2.1), а тем самым и матрицы весовых множите лей (см. п. 2.3) строится рекуррентно путем умножения элементов предыдущей (к — 1)-й строки на величину с добавлением блока q K на диагонали. Следовательно, элементы нижней треугольной матрицы весовых множителей выражаются формулой: г—1 Свободные элементы системы здесь равны Используя бл очно-треугольный вид системы также как в разделе 2.3, получаем, что если то где j0 - последняя компонента векторной оценки Если же Mh с/2, то EJ2 +00 и, более того, величина Е02 ограничена равномерно по г о и по Л. Отметим, что в случае постоянных коэффициентов с№ строго по аналогии с вышеизложенным можно построить оценку на "блуждании по сферам" [11, 41] для решения задачи (2.27) непосредственно, без использования разностного приближения. Обоснование такой оценки осуществлены на основе результатов работы [31] и треугольного вида соответствующих базовых операторов К и Кр согласно изложенному в п. 2.3. Трудоемкость оценки на "блуждании по сферам" ограничена величиной порядка c„lm! 5-2, причем Сп — асимптотически линейно по п растущая функция. Нетрудно видеть, что при с = const р-кратные производные по этому параметру от решения задачи удовлетворяют соотношениям то есть фактически реализуют итерации резольвенты. Следовательно, Прямое дифференцирование весовой оценки = [1 — ch2/2n] N решения задачи (2.28) дает соотношение Согласно изложенному в п. 2.3 настоящего раздела, величина Df ограничена равномерно по h при с c h/2. Ясно, что рассмотренный здесь алгоритм автоматически распространяется на оценку первого собственного числа для более общего оператора Ли + с(г), соответствующего уравнению Шреди J2 +00 и, более того, величина Е02 ограничена равномерно по г о и по Л. Отметим, что в случае постоянных коэффициентов с№ строго по аналогии с вышеизложенным можно построить оценку на "блуждании по сферам" [11, 41] для решения задачи (2.27) непосредственно, без использования разностного приближения. Обоснование такой оценки осуществлены на основе результатов работы [31] и треугольного вида соответствующих базовых операторов К и Кр согласно изложенному в п. 2.3. Трудоемкость оценки на "блуждании по сферам" ограничена величиной порядка c„lm! 5-2, причем Сп — асимптотически линейно по п растущая функция. Нетрудно видеть, что при с = const р-кратные производные по этому параметру от решения задачи удовлетворяют соотношениям то есть фактически реализуют итерации резольвенты. Следовательно, Прямое дифференцирование весовой оценки = [1 — ch2/2n] N решения задачи (2.28) дает соотношение Согласно изложенному в п. 2.3 настоящего раздела, величина Df ограничена равномерно по h при с c h/2. Ясно, что рассмотренный здесь алгоритм автоматически распространяется на оценку первого собственного числа для более общего оператора Ли + с(г), соответствующего уравнению Шредингера. Отметим, что метод Монте-Карло дает возможность одновременно оценивать производные по параметрам задачи и вероятностные моменты решения в задачах со случайными (в допустимых пределах) параметрами. Он допускает также идеальное распараллеливание вычислений путем простого рапределения независимых испытаний по вычислительным процессорам. 2.2) в единичном квадрате D — [0,1] х [0,1]. Решениями для задач (2.1), (2.2) являются нгера. Отметим, что метод Монте-Карло дает возможность одновременно оценивать производные по параметрам задачи и вероятностные моменты решения в задачах со случайными (в допустимых пределах) параметрами. Он допускает также идеальное распараллеливание вычислений путем простого рапределения независимых испытаний по вычислительным процессорам. 2.2) в единичном квадрате D — [0,1] х [0,1]. Решениями для задач (2.1), (2.2) являются соответственно функции и = sin(;r) sin(y) и и = ехеу. Оценка решений проводилась в одной точке области D согласно формуле полученной в теореме 5. j=o Замечание 3. Если сумму из (2.17) записать в виде ft- зі то спРа_ ведливо рекуррентное соотношение Рассмотрим теперь в области D С R" первую краевую задачу для неоднородного бигармонического уравнения А2и = -д, и\г = fo, Аи\т = р\. (3.15) Для построения соответствующей оценки щ j используем равенства Случайная величина і/,, распределенная в интервале (0,() с плотностью бл;(1 — x/di)d и единичный изотропный вектор ШІ моделируются при помощи известных формул [11]. Замечание 5. Для случая 5 = 0 выражение (3.17) можно получить также из известной векторной оценки [33]. Колебания пластины в ограниченной области D С it!2 под действием случайного поля нагрузок cr(r) = —д(г) описываются уравнением вида (3.15) [2], при этом можно учитывать также случайность граничных функций ро(г) и ipi(r). Решение этой задачи также является случайным полем. Требуется определить его ковариационную функцию г? (г, г ) = Щи(г)и(г )}, Предполагается, что E#(r) = E o(r) = E j(r) = 0 и, следовательно, EU(T-) = 0. Согласно формулам (3.4) и (3.5) имеем K(r, г 1)) = EbWff t1))], /f0(r, r «) = Е[5(г) (гСі))], /Ti rW) = EbM frW)], "«,(»-, г) = E[V0(r) (r{1))], #10(7-, ) = EfoW r)], ЛГц(г,г(") = E[Vl(r) (rC1})]. Таким образом, при помощи выражений (3.22), (3.23) можно оценивать ковариационную функцию решения v(r}г ), используя только ковариационные функции случайного поля нагрузок и случайных граничных условий. Замечание 6. Построенный метод реализуется и исследуется существенно проще логически сложного алгоритма рассмотренного в [28]. Соответствующая оценка имеет дисперсию, заведомо меньшую дисперсии использованной в [33] "двойной рандомизации" [11], так как здесь осуществлено частичное аналитическое осреднение. Однако, при реализации полученного метода возникает проблема, связанная с ограниченностью машинной памяти. В отличии от реализации оценок вида (3.17) здесь необходимо сохранять все веса и координаты центров случайных кругов в (3.23) хотя бы для одной траектории. Известно, что процесс "блуждания по окружностям" весьма быстро сходится к границе области. Следовательно, можно при превышении количества точек {г{} некоторого достаточно большого уровня М не сохранять информацию о следующих точках с номерами і — JW +1,..., ЛГ — 1, но в Гє использовать полученные вес QN И координаты TN- При этом Qi — О для і = М + lt... ,N — 1, и можно заменить QM на QM{N — М + \). Если такая замена практически не влияет на результат, то оценка удовлетворительна. Порядок (по є) величины М можно эвристически оценить с использованием асимптотической теории восстановления по аналогии с [11]. В качестве "эталонной" рассмотрим задачу вычисления величины EiV2. Введем обозначения Использование асимптотической теории восстановления (см. [11]) дает возможность эвристически предположить, что при є — О Последнее соотношение может быть практически вполне удовлетворительным для достаточно больших є, но оно весьма завышает величину М при существенно малом є в связи с тем, что величина асимптотически стандартно нормальна [3]. При этом можно использовать известную асимптотику гауссовской функции распределения: использовать, например, значение М х [1пє2. Поскольку неравенства N М и t(N) t(M) эквивалентны, то соответственно (3.24) следует оценивать величину Таким образом, значение М порядка 1пє)2 при малом є, обеспечивающем достаточную точность нормального приближения, может быть вполне удовлетворительным. Это подтверждается также следующими соображениями. При использовании г-смещённой оценки "блуждания по сферам" целесообразно моделировать Сє 2 случайных траекторий [11]. Проведённые выкладки показывают, что для М х ]1пг2 имеем P(N М) = о(є2) и, следовательно, вероятность появления события N М в ансамбле обьёма Се 2 есть величина о(1) при е — 0. Отметим, что полученные оценки являются асимптотически точными для плоской и, следовательно, выпуклой границы [11]. Эвристически ясно, что они справедливы и в случае гладкой границы ввиду быстрой сходимости к ней процесса "блуждания по сферам (кругам)" [11]. Окончательный выбор М должен определятся тестовыми расчётами с использованием замены веса Qu на QM{N — М +1), как указано выше (см. далее раздел 3.6.2). 3.5. Оценка ковариационной функции решения уравнения А2и 4- си = д со случайными параметрами Если пластина лежит на упругом основании, то нормальный прогиб пластины и{т) удовлетворяет уравнению В задачах оценки показателей надёжности различных сооружений, подверженных случайным вибрациям, функциональные параметры задачи имеют следующий физический смысл [2]: с = к/В - коэффициент постели, д = f(-)/B - случайное поле нагрузки. Здесь /() - интенсивность нормальной нагрузки; В Eh3/12(1 — а2), где Е - модуль Юнга; т -постоянная Пуассона материала пластины; 2h - её толщина. Далее, рассмотрим случай, когда константа с является случайной величиной, причём Ее = 0, Dc 1, с с . Функции д, щ и р\ случайны и не зависят от параметра с. Пренебрегая в разложении решения и(г7с) в ряд Лорана слагаемым о(с), имеем и(г, с, д, с/70» pi) « и{г, 0,5, о, ?i) + cu(1) (г, 0, gt (р0, (р{). (3.26) Следовательно, для ковариационной функции решения задачи (3.25) выполняется соотношение Cov(u(rh с,д, р0,(pi)tu(r2,с,д, р0, р{)) ft v(rbг2) + г 1 (3.27) Погрешность разложения (3.26) можно приблизительно считать равной 2w(2)(r S, Ab i)c2. Таким образом, при вычислении ковариационной функции по формуле (3.27) погрешность имеет вид 25 = Е[-ц(2) (г, 0,д, о v?i)u(r, 0,5, о, Pi)\Dc + Е[« (г, 0, д, р0, i)u(1) (г, 0,5j о, i)]E[c3] + Щи{2\г, О, с,, Vo, i)u(2)(r, 0, д, р0, Ы]Е[с4] Замечание 7. Решение и(г, с, д, 9 о, i) можно приближённо заменять соответствующей частичной суммой ряда Лорана с любым количеством слагаемых. Результаты численных расчётов, приведённые в пункте 3.6.2, показывают, что при достаточно малых отклонениях случайной величины с от 0 замена (3.26) приемлема. . Для построения оценки неизвестных параметрических производных и (г, О, д1 о, y?i), u ir О, д, (/?о, і) воспользуемся следующим приёмом. Дифференцируя по параметру с уравнение (3.25), получаем при с—О: А2и = и, и =Д« = 0. Следовательно, параметрическая производная и \г, 0, ?, о, рі)/2 является решением следующей задачи: (3.33) t-1 где Si(c) = Yl s(c,dj). Из представления оценки видно, что необходимо 7=0 найти значения производных 5 (0) для к 0,3. Дифференцируя по параметру с логарифмическую производную 5 1) — 5(1п5) \ нетрудно получить следующие соотношения 5(2) = 5 (1п5)« + 5(1п5) , 5(3) = S 2)lnS 4 + 25W(lnS)W + 5(ln5) 3), 5 4 = W1) + 35 2 (ln5) 2 + 35 (hi5) + 5(Ь5) 4\ 5(5) = 1115 + 45(ln5)W + 65 2 (ln5) 3) + 45 )(1п5)(4) + 5(1д5) ). Производные (In 5) ), к равенств 1,5 при с = 0 легко находятся с помощью . ,. , cP 3 f , 13d6 з 211 , 1217 , . e,Вероятностное представление для решения первой краевой задачи для уравнения (L + С)Р+1ІІ д
Глобальная оценка решения уравнения и задача минимизации трудоёмкости
Оценки решения метагармонического уравнения и первого собственного числа многомерного оператора Лапласа
Оценки для решения бигармонического уравнения при п = 2, 3
Похожие диссертации на Скалярные алгоритмы метода Монте-Карло для решения метагармонических уравнений
