Решение задачи Дирихле методом Монте-Карло

Введение

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

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

Цель данной работы рассмотреть решение задачи Дирихле для уравнения Лапласа и уравнения Пуассона методом Монте-Карло на основе метода сеток. 

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

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

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



















1. Метод Монте-Карло

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

а) речь идет о численных методах (и конкурировать они могут с классическими численными методами, а не с аналитическими методами решения задач);

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

Официальной датой рождения методов Монте-Карло считают 1949 год, когда появилась статья под заглавием «Метод Монте-Карло». Возникновение метода связывают обычно с именами Дж. Неймана, С. Улама, Н. Метрополиса, а также Г. Кана и Э. Ферми; все они в 40-х годах работали в Лос-Аламосе (США). 

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

Развитию методов Монте-Карло способствовало бурное развитие ЭВМ. Алгоритмы Монте-Карло (как правило, обладающие небольшой связностью) 

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

области их применения.

Важнейший прием построения методов Монте-Карло – сведение задачи к расчету математических ожиданий. Более подробно: для того чтобы приближенно вычислить некоторую скалярную величину а, надо придумать такую случайную величину , что  ; тогда, вычислив  независимых значений  величины , можно считать, что .

Пример. Требуется оценить объем  некоторой ограниченной пространственной фигуры .

Выберем параллелепипед , содержащий , объем которого  известен. Выберем  случайных точек, равномерно распределенных в , и обозначим через  количество точек, попавших в . Если  велико,

то, очевидно,  : , откуда получаем оценку 

	.

В этом примере случайная величина  равна , если случайная точка попадает в , и  равна нулю, если точка попадает в . Нетрудно проверить, что математическое ожидание , а среднее арифметическое 

	.

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

1) как выбрать удобную величину  для расчета той или иной задачи;

2) как находить значения произвольной случайной величины ?

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

Многие методы основаны на расчете математических ожиданий. Существуют методы случайного поиска (кроме простейшего) и стохастических приближений.

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



2. Задаче Дирихле для уравнения Лапласа

2.1 Основные сведения о задаче Дирихле для уравнения Лапласа

Определение. Функция , имеющая непрерывные частные второго порядка в области  и удовлетворяющая внутри  уравнению Лапласа, называется гармонической функцией: 

	.

Простейшим примером гармонической функции двух переменных является функция вида , где  (основное решение уравнения Лапласа).

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

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

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

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

Составим вспомогательную функцию 

	,

где  – диаметр области . Очевидно, имеем

	,

причем при  выполняется неравенство

	.

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

.

Из соотношения 



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

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

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

	,

где  на ,  на .

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

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

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

	.

Очевидно, что на  – гармоническая функция, обращающаяся в нуль на границе. По свойству I эта функция не может принимать внутри  значений больше или меньше нуля, следовательно,  внутри  и  .

Замечание. Из свойства II не следует, что задача Дирихле для ограниченной замкнутой области  имеет решение; это свойство лишь утверждает, что если существует решение задачи Дирихле для области , то оно единственно.

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

Свойство III (корректность задачи Дирихле). Решение задачи Дирихле для замкнутой и ограниченной области непрерывно зависит от граничных данных.

Доказательство. Допустим, что и  – решения задачи Дирихле, соответственно принимающее на границе значение  и .

Пусть всюду на  выполнено неравенство

	,

где  – произвольное малое положительное число.

Рассмотрим гармоническую функцию

	.

На границе  эта функция принимает значение 

	.

Так как  на , то по свойству I имеем

	  при ,

т. е.  или .

Таким образом, для задачи Дирихле требование корректности выполнено при .





2.2 Решение задачи Дирихле для уравнения Лапласа методом сеток

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

в плоскостной задаче, в которой разыскивается решение, строится сеточная область , состоящая из одинаковых ячеек (рис. 1, Приложение А) и приближающая данную область ;

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

на основании граничных условий устанавливаются значения искомого решения в граничных узлах области .

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

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

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

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

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

Покажем применение метода сеток для построения решения задачи Дирихле

	при и  при ,     (1)

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

Выбрав шаг , построим квадратную сетку

	

	

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

Точки (узлы)  называются соседними, если они удалены друг от друга в направлении оси или оси  на расстояние, равное шагу сетки . Узел  сетки  называется внутренним, если он принадлежит области , а все четыре соседних с ним узла – множеству  ; в противном случае он называется граничным (например, узлы сетки ) (рис. 1, Приложение А – внутренние узлы обозначены светлыми кружками, а граничные – темными кружками и темными треугольниками).

Граничный узел сетки  называется узлом первого рода, если он имеет соседний внутренний узел этой сетки (например, узел  – рис. 1, Приложение А); в противном случае граничный узел называется узлом второго рода (узел  – рис. 1, Приложение А). Внутренние узлы и граничные узлы первого рода сетки  называются расчетными точками. Граничные узлы второго рода не входят в вычисление и могут быть изъяты из сетки  (рис. 1, Приложение А –   граничные узлы второго рода обозначены темными треугольниками). 

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

Значение искомой функции  в точках  обозначим через . Для каждой внутренней точки  сетки  заменяем дифференциальное уравнение (1) конечно-разностным уравнением 

	,  (2)      

где  – расчетные точки.

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

	,     (3)

где  – ближайшая к точка границы .

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

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

	     (4)

для всех узлов сетки . В силу неравенства (4) будем иметь 

	.   (5)

На основании системы (2) получаем 

	.    (6)

Учитывая неравенство (4), заключаем, что 

	.

Ни одно из последних четырех неравенств не является строгим, так как если бы это имело место, то, складывая все четыре неравенства и учитывая формулу (6), мы получили бы .

Поэтому 

	.    (7)

Проводя аналогичные рассуждения для точек  сетки  с ближайшей точкой , где положено

	.

Таим образом, из цепи равенств (7) имеем  , что противоречит неравенству (5).

Так, однородная система (2) не может иметь положительных решений. Аналогично доказывается, что эта система не может иметь отрицательных решений. Следовательно,  для каждого решения, и, значит, неоднородная система (2) совместна и имеет единственное решение.

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





2.3 Понятие о решение задачи Дирихле для уравнения Лапласа методом Монте-Карло 

Пусть на плоскости дана область  с кусочно-гладкой границей . В области   построим квадратную сетку с шагом :

	,       (1)

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

Представим себе частицу , которая совершает равномерное случайное блуждание по узлам сетки (1). А именно, находясь во внутреннем узле  сетки , эта частица за один переход с одной и той же вероятностью, равной 1/4, может переместиться в один из четырех соседних узлов: или в   (шаг влево), или в  (шаг вправо), или в  (шаг вниз), или в  (шаг вверх), причем каждый такой единичный переход совершенно случаен и не зависит от положения частицы и ее прошлой истории. Будем считать, что блуждание частицы  заканчивается, как только эта частица попадет на границу ; в этом смысле граница  представляет собой «поглощающий экран». Можно доказать, что с вероятностью, равной 1, блуждание точки через конечное число шагов заканчивается на границе.

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

Равномерное случайное блуждание частицы на плоскости можно организовать с помощью равномерно распределенной последовательности одноразрядных случайных чисел, принимающих значения. Для этого, например, достаточно производить розыгрыш, т. е. случайную выборку из чисел , придерживаясь инструкции, указанной в таблице 1 (Приложение B); причем числа 8 и 9 переигрываются.

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

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

	.

Для краткости введем обозначение

	.

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

	,      (2)    

где суммирование распространяется на все точки  границы , причем 

	               (3)

где  –  граничный узел.

Составим сумму 

	,     (4)

где точка  пробегает всю границу . Если функцию  рассматривать как случайную величину, принимающую значения  на границе , то сумма (4) представляет собой математическое ожидание (среднее значение) функции  на границе  для траекторий, начинающихся в точке  («премия за выход на границу» из начальной точки ). Частица, начавшая свое случайное блуждание из внутреннего узла , после первого шага с вероятностью, равной 1/4, попадает в один из четырех соседних узлов. Поэтому случайные блуждания, начинающиеся в узле , в зависимости от вида траекторий распадаются на четыре категории новых случайных блужданий: 

	

По формуле полной вероятности имеем 

	  (5)

Отсюда, умножая обе части равенства (5) на граничные значения  и суммируя по всем возможным значениям  и , на основании формулы (4) получим 

	.           (6)

Кроме того, в силу формулы (3) имеем

	,           (7)

если точка .

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

	        (8)

Сравнивая формулы (8) с формулами (6), (7), мы усматриваем, что они совпадают с точностью до обозначений. Следовательно, искомые неизвестные  можно рассматривать как математические ожидания . Величины   допускают экспериментальное определение. Рассмотрим достаточно большое число  равномерных случайных блужданий частицы по узлам сетки , исходящих из фиксированного узла  и заканчивающихся на границе . Пусть   соответствующие точки выхода частицы на границу . Заменяя математическое ожидание  эмпирическим математическим ожиданием, будем иметь 

	.    (9)

Формула (9) дает статистическую оценку величины  и может быть применена для приближенного решения задачи Дирихле. Метод решения задач, основанный на использовании случайных величин, получил общее название метода Монте-Карло. 

Заметим, что с помощью формулы (9) можно непосредственно найти приближенное значение  решения задачи Дирихле в единственной фиксированной точке  сетки , не зная решения задачи для остальных точек сетки. Этим обстоятельством метод Монте-Карло для задачи Дирихле резко отличается от обычных стандартных способов решения этой задачи.

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

	.

Построив такую функцию Грина, мы получаем возможность, применяя формулу (9), просто

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

Недостатком рассмотренного варианта метода Монте-Карло для задачи Дирихле является слабая сходимость по вероятности при  эмпирического математического ожидания 

	

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

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





2.4 Метод  «блуждания»  по сферам

Укажем другой метод Монте-Карло для решения задачи Дирихле для уравнения Лапласа, не связанный с разностными уравнениями. Пусть задана ограниченная связная область  и точка . Определим случайную траекторию  следующим образом: положим ; далее, если точка  известна, то построим окружность произвольного радиуса , расположенную внутри , и на этой окружности выберем случайную точку  (рис. 2, Приложение C). 

Таким образом, , где , и угол  равномерно распределен в интервале . 

Приведем теорему: если функция удовлетворяет в области уравнению Лапласа

	,               (1)

то при каждом и при любых математическое ожидание  равно значению  в начале траектории. 

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

	.

По теореме о среднем значении гармонической функции 

	.

Поэтому 

	.

При  точка  и . Применяя индукцию, получим утверждение теоремы.

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

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

Фиксируем достаточно малую окрестность  границы (рис. 3, Приложение D). Чтобы вычислить , будем строить траектории вида до тех пор, пока случайная точка  не попадет в . Пусть  – ближайшая к  точка границы . Можем считать, что значение случайной величины  приближенно равно . Построив  траекторий такого типа, получим значения , по которым оценивается искомое решение 

	.     (2)

Замети, что сходимость по вероятности 

	,    (3)

когда   не вытекает из теоремы Хинчина, говорящей о том, что последовательность одинаково распределенных независимых величин, у которых существуют математические ожидания, подчиняется закону больших чисел, так как в сумме (3) фигурируют  различных случайных величин, различающихся правилами выбора  Можно, однако воспользоваться другой формой закона больших чисел – теоремой Чебышева:

Если величины  независимы и существует  и , то при 

	

(Доказательство этой теоремы легко получить, применяя к величине  неравенство Чебышева – ).

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

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

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







Задача Дирихле для уравнения Пуассона и ее решение методом Монте-Карло с использованием метода сеток

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

	   (1)

где  – искомая функция. Уравнение (1) при  называется уравнением Лапласа, а при  – уравнением Пуассона.

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

Предположим, что на границе  задана некоторая функция (часто пишут g(S), где  — длина дуги границы, отсчитываемая от какой-нибудь фиксированной точки). Требуется найти такое решение  уравнения (1), которое на границе совпадает с :

          .      (2)

Задачу об отыскании решения уравнения 

    (3)                                    

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

Во внутреннем узле  уравнение заменим разностным 

                                          (4)

которое перепишем в виде

                (5)

В граничных узлах  

                       .                        (6)

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

Перенумеруем все узлы, принадлежащие (в произвольном порядке), и перепишем в том же порядке уравнения (3), (4). Тогда получим систему вида

                                                                   (7)

Матрица этой системы имеет следующую структуру: внутреннему узлу с номером  отвечает строка , в котором четыре элемента равны ¼, а остальные – нули; граничному узлу с номером  отвечает строка; все диагональные элементы=0. (Все собственные значения такой матрицы по абсолютной величине меньше единицы.) Свободные члены этой системы, если узел  внутренний, и , если узел  граничный.

Один из методов решения системы (5) является метод Монте-Карло. Построим данный метод для расчета  - значения решения в одном заранее заданном узле. Выберем матрицу переходов 



Здесь – символ Кронекера:.

Далее строим следующую цепь:

1)  

2) если узел  внутренний, то с одинаковой вероятностью ¼ выбираем в качестве   номер одного из соседних с ним узлов;

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

                                  ,                                     (8)   

где – номер первого выхода цепи на границу.

В формуле (6) все  вычисляются по формуле  и лишь последнее равно значению .

Замечание. Если вместо граничных условий (2) заданы более сложные условия, например:

,

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

Если количество цепей достаточно велико, то решение задачи Дирихле в узле определяется по формуле 

                                        .                                       (9)







4. Прикладное применение метода Монте-Карло с использованием метода сеток для решения задачи Дирихле для уравнения Лапласа

Пусть  – решение уравнения Лапласа  в единичном квадрате , удовлетворяющее граничным условиям . Вычислить значение .

Выберем в квадрате сетку с шагом  и перенумеруем узлы (рис. (4), Приложение Е). Для уравнения Лапласа формула (8) все более упрощается: , так что  равно значению  в том узле, в котором цепь попадает на границу. Возле каждого граничного узла на рис. (4) (Приложение Е) проставлено значение для данного примера.

Для построения цепей необходимо воспользоваться таблицей случайных цифр (таблица 1, Приложение В). 

Если случайная цифра  окажется 0 или 4, то будем перемещаться в соседний узел справа, если  окажется 1 или 5, то будем перемещаться влево,  окажется 2 или 6, то перемещаться вверх, если  окажется 3 или 7, то перемещаться вниз; значения , равные 8 или 9, опускаем. 

В таблице 2 (Приложение F) приведены 16 случайных цепей. В первой строке записаны использованные случайные цифры, а в третьей – сама цепь (номера ). Соответствующие этим цепям значения  равны . Среднее арифметическое этих величин дает нам приближенное значение решения в точке :

	.

   Из эмпирической оценки дисперсии 

	

следует, что вероятная ошибка .

Точное решение рассмотренной задачи , так что , и фактическая ошибка расчета равна 0,08.

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

































5. Математическое обоснование решения задачи Дирихле для уравнения Пуассона

Найдем решение задачи Дирихле для уравнения Пуассона: 

	,

где  .

Выберем в квадрате  сетку с шагом . Для построения цепей используем таблицу случайных цифр (таблица 1, Приложение В). Если случайная цифра  окажется 0 или 4, то будем перемещаться в соседний узел справа, если  окажется 1 или 5, то будем перемещаться влево,  окажется 2 или 6, то перемещаться вверх, если  окажется 3 или 7, то перемещаться вниз; значение , равные 8 или 9, опускаем. 

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

                                  ,                                     (8)   

где – номер первого выхода цепи на границу.

В формуле (6) все  вычисляются по формуле , где , и лишь последнее  равно значению : .

Итоговое значение функции  получаем по формуле , где .





















Заключение

В данной работе были рассмотрены основные сведения, связанные с задачей Дирихле для уравнений Лапласа и Пуассона – определения, свойства и методы решения. Было приведено два метода решения данной задачи с помощью метода Монте-Карло – метод сеток и метод «блуждания» по сферам для уравнения Лапласа и метод сеток для уравнения Пуассона. Приведено математическое обоснование решения задачи Дирихле для уравнения Пуассона методом Монте-Карло с использованием метода сеток.

В приложении приведена программа, написанная на BorlandPascal 7.0, реализующая данный метод с заданными исходными данными: 

	,

.

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































Список использованной литературы 

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

2. Демидович Б. П., Марон И. А., Шувалова Э. З. Численные методы анализа. – М.:Наука, 1967. – 368 с. 

3. Березин И. С., Жидков Н. П. Методы вычислений. – М.:Государственной литературы, 1959. – 602 с.

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



















































Приложения

А. Сеточная область  





(рис. 1)































































B. Таблица 1

(блуждание частицы на плоскости)

Случайное число

Характер перемещения

0 или 4

 (шаг вправо)

1 » 5

 (шаг вверх)

2 » 6

 (шаг влево)

3 » 7

 (шаг вниз0











































С. Ограниченная область  



(рис. 2)









































D. Ограниченная область  с границей 



(рис. 3)







































E. Единичный квадрат



(рис. 4)

Единичный квадрат: в нем сетка с шагом , краевые условия: 













































F. Таблица 2 

(случайные цепи)

6   5   1

   

13–18–17–16

  5   0   7   5    6  6   1

       

13–12–13–8–7–12–17–16

5  5



13–12–11

6  6

       

13–18–23

4   3   4

            

13–14–9–10

5   6   5

  

13–12–17–16

5  1



13–12–11

2   3   3   2   4  3   7

        

13–18–13–8–13–12–17–16

7   5   7

                   

13–8–7–2

0   2   6

               

13–14–19–24

1   6   0   3   3  3  

          

13–12–17–18–13–8–3

4   2   5   0   2   

        

13–14–19–18–19–8

2  2

 

13–18–23

            4   5     5    5      

        

13–14–13–12–11

3  7

 

13–18–3

          5  1

         

13–12–11



























G. Программа, реализующая решение задачи Дирихле методом Монте-Карло

program kp;

Uses crt;

Const  x=3/4; y=1/2; N=10; h=1/4;

Var s1,f,s,u,u1,g,x_1,y_1: real;

    i : integer;



procedure Loop;

var f1,f2,f3,f4,f5:real;

begin



 y_1:=y;  x_1:=x; f:=0;

repeat

g:=random(10);

writeln('g=',g);



f1:=f;

if (g=0) or (g=4) then begin

 x_1:= x_1 + h;

 f1:=f1-0.25*h*h*(x_1*x_1+y_1*y_1);

 f:=f1;

 end;



f2:=f;

if (g=1) or (g=5) then begin

  y_1:= y_1 + h;

 f2:=f2-0.25*h*h*(x_1*x_1+y_1*y_1);

  f:=f2;

  end;



f3:=f;

if (g=2) or (g=6) then begin

  x_1:= x_1 - h;

  f3:=f3-0.25*h*h*(x_1*x_1+y_1*y_1);

  f:=f3;

  end;



f4:=f;

if (g=3) or (g=7) then begin

  y_1:= y_1 - h;

  f4:=f4-0.25*h*h*(x_1*x_1+y_1*y_1);

  f:=f4;

   end;



f5:=f;

if (g=8) or (g=9) then begin

   y_1:=y_1;

   x_1:=x_1;

   f5:=f5-0.25*h*h*(x_1*x_1+y_1*y_1);

   f:=f5;

     end;

until (x_1=0) or (x_1=1) or (y_1=0) or (y_1=1);



if y_1=1 then

s:=x_1*x_1;

if (y_1=0) or (x_1=1) then

s:=0;

if x_1=0 then

s:=y_1*y_1;

s1:=s1+s;

writeln('s=',s);

writeln('s1=',s1);

end;



begin

 clrscr;

 randomize;

 u:=0; u1:=0;

 for i:=1 to N do begin

 Loop;

 writeln('f=',f);

 u:=u+f+s;

 writeln('u=',u);

 end;



u1:=u/N;

 writeln('u1=',u1);

writeln('press any key');

 readkey;

END.

5