Применение схемы Неймана-Улама к решению первой краевой задачи для параболического уравнения Текст научной статьи по специальности «Математика»

УДК 519.71

Вестник СПбГУ. Сер. 1. Т. 1 (59). 2014. Вып. 1

ПРИМЕНЕНИЕ СХЕМЫ НЕЙМАНА—УЛАМА К РЕШЕНИЮ ПЕРВОЙ КРАЕВОЙ ЗАДАЧИ ДЛЯ ПАРАБОЛИЧЕСКОГО УРАВНЕНИЯ*

*

А. С. Сипин

Вологодский государственный педагогический университет, Российская Федерация, 160035, Вологда, ул. С.Орлова, 6

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

В данной работе построены несмещенные и малосмещенные оценки решения первой краевой задачи для уравнения теплопроводности с переменным коэффициентом при неизвестной функции на траекториях марковской цепи «блуждания по шароидам». Для исследования свойств цепи Маркова и свойств статистических оценок используется предложенное автором обобщение известной в теории методов Монте-Карло схемы Неймана—Улама на интегральные уравнения с субстохастическим ядром. Основой алгоритма является новое интегральное представление решения краевой задачи. Библиогр. 8 назв.

Ключевые слова: метод Монте-Карло, статистическое моделирование, уравнение теплопроводности, блуждание по шароидам.

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

Пусть Б — ограниченная область в Кп с границей Г. Рассмотрим первую краевую задачу

коэффициент ао(х, которого принадлежит классу (С^т), а < 1. Матрица коэф-

фициентов при старших производных предполагается постоянной и симметричной. В [1] показано, что если / (Е (<3т), € Функция Ф — непрерывна, а Г —

поверхность Ляпунова, то задача (1.1) имеет классическое, непрерывное вплоть до границы решение в области .

Несмещенные оценки для решения уравнения теплопроводности (при ао(х, £) = 0) построены как на траекториях блуждания по границе [2], так и на траекториях блуждания внутри пространственно-временного цилиндра [3]. В первом случае решение представляется в виде теплового потенциала двойного слоя, плотность которого удовлетворяет интегральному уравнению Вольтерра—Фредгольма. К этому уравнению применимы схема Неймана—Улама и общая теория статистических оценок для нее.

* Работа выполнена при финансовой поддержке РФФИ (грант № 11-01-00769-а).

Ьп = /, п\г = Ф(х, ¿), п\г=0 = у(х). в области Qт = Б х (0, Т), где Ь — параболический оператор,

(1.1)

(1.2)

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

Для уравнений с переменными коэффициентами фундаментальное решение неизвестно, что существенно осложняет методы построения несмещенных и малосмещен-ных статистических оценок. В данной работе построены несмещенные оценки решения u(x,t) задачи (1.1) на траекториях марковской цепи «блуждания по шароидам».

2. Схема Неймана—Улама для решения краевых задач. Стохастические методы решения краевых задач для уравнений в частных производных чаще всего основаны на интегральном представлениии решения краевой задачи. Это интегральное представление обычно рассматривается как интегральное уравнение, к которому применяется схема Неймана—Улама [4]. Особенности применения схемы к решению краевых задач подробно рассмотрены в [5]. Приведем необходимые определения и теоремы из данной работы.

Довольно часто полученное интегральное уравнение рассматривается в пространстве M(Q) ограниченных борелевских функций на некотором компакте Q в евклидовом пространстве Rn и имеет вид

u(x) = / u(y)P(x, dy) + F(x), x e Q, (2.1)

Jq

где P(x, dy) — субстохастическое ядро, а F(x) является интегралом от граничных условий или правой части уравнения. Известно [6, 7], что при F(x) > 0 решение уравнения (2.1) представимо в виде

u(x) = Y1 RiF(x) + K™u(x), x e Q, (2.2)

¿=0

где K — интегральный оператор в представлении (2.1), а

K^u(x) = lim Kiu(x). (2.3)

i—

Действительно, из уравнения (2.1) следует неравенство

n

Y,KiF(x) = u(x) - Kn+iu(x) < 2\\u\\.

i=0

Значит, ряд в (2.2) сходится и существует предел (2.3). Равенство (2.2) очевидно. Для функции v(x) = Kxu(x) в силу теоремы Лебега об ограниченной сходимости справедливо равенство v(x) = Kv(x). Такие функции будем называть инвариантными. Сумму ряда

^

GF (x) K¿F (x)

i=0

называют потенциалом функции F(x). Заметим, что Ku(x) < u(x). Функции, обладающие этим свойством, называются эксцессивными. Представление эксцессивной функции в виде суммы потенциала и инвариантной функции называется разложением Рисса. Оно единственное, так как K™GF(x) = 0. Суммируем все сказанное в следующей теореме:

Теорема 1. Пусть функции F(x) > 0, тогда ограниченное решение u(x) уравнения (2.1) единственным способом раскладывается в сумму потенциала и инвариантной функции.

Несмещенные оценки решения уравнения (2.1) обычно строят на траекториях цепи Маркова |xi}°=o, определяемой переходной вероятностью P(x, dy). Процесс обрывается в некоторый случайный марковский момент ti при переходе в поглощающее состояние Д, лежащее вне Q (при этом u(A) = 0).

2.1. Свойства траекторий цепи Маркова. Изучим некоторые свойства последовательности |xi}°=o , стартующей из точки xo = x. Пусть {Я^}ТО=о —последовательность а-алгебр, порожденная цепью до момента времени i, Xi —индикатор события {ti > i}. Пусть u(x) —ограниченная эксцессивная функция, удовлетворяющая уравнению (2.1). Определим стандартную последовательность несмещенных оценок

i-l

Vi = ^2 F(xj)xj + Xiu(xi), (2.4)

j=o

которая, очевидно, является равномерно интегрируемым мартингалом относительно потока {^¿}°=о и

Ex Xiu(xi) = Kiu(x).

По теореме сходимости мартингалов существует = limni Px п.н. и Ex= u(x).

Ряд E TO =0 F(xj )xj Px п.н. сходится в силу теоремы Леви. Поэтому по теореме Лебега

об ограниченной сходимости KTOu(x) = Ex limxiu(xi).

Теорема 2. 1. Если для всех x e Q вероятность Px(ti < ж) равна 1, то всякая ограниченная эксцессивная функция является потенциалом.

2. Если функция u(x) = 1 является потенциалом, то для всех x e Q вероятность Px(ti < ж) равна 1.

3. Если существует строго положительная инвариантная функция, то при всех x

Px (ti = ж) > 0.

4- Пусть Г —множество нулей неотрицательной непрерывной на Q функции F(x), потенциал которой всюду конечен. Пусть Px(ti = ж) > 0 и p(x, Г) — расстояние от x до Г, тогда

Px (lim p(xi, Г) = 0|ti = ж) = 1.

5. Пусть Г1 —множество нулей неотрицательной непрерывной на Q функции Fi(x), Г2 — множество нулей неотрицательной непрерывной на Q функции F2(x), потенциалы которых всюду конечны. Пусть Px(ti = ж) > 0, тогда

Px(limp(xi, Г1 П Г2) = 0|t1 = ж) = 1.

Как правило, сама последовательность {xi}°=o п.н. сходится. В частности, справедлива

Теорема 3. 1. Пусть при всех m = 1,... ,n координатные функции vm(x) таковы, что при некоторой ограниченной эксцессивной функции wm (x) сумма wm(x) +

ут(х) или разность шт(х) — ут(х) являются эксцессивными функциями; тогда последовательность {хг}?=о п.н. сходится на множестве {т1 = то}.

2. Пусть при всех т = 0,1,...,п координатные функции ут(х) таковы, что при некоторой постоянной сумма и>т + Ут(х) или разность — ут(х) являются эксцессивными функциями; тогда последовательность {х}°=0 п.н. сходится на множестве {т1 = то}.

3. Пусть при всех т = 0,1, .. . ,п координатные функции ут(х) таковы, что либо ут(х) или —ут(х) являются эксцессивными функциями, либо ут(х) — инвариантна; тогда последовательность {х}°=о п.н. сходится на множестве {т1 = то}.

2.2. Статистические оценки. Стандартная последовательность несмещенных оценок (2.4) нереализуема, так как содержит функции Г(х) и п(х), значения которых неизвестны. Чтобы получить реализуемую оценку, эти функции оценивают либо несмещенно, либо с «малым» смещением. Опишем соответствующую процедуру, следуя [3]. Будем предполагать, что выполнено условие

а) последовательность {х}°=о п.н. сходится на множестве {т1 = то} к некоторой точке хто € дQ.

Пусть 5 > 0, т2 — момент первого попадания цепи в ¿-окрестность границы дQ и т§ = шш(т1,т2). Для марковской цепи, удовлетворяющей условию а величина т§ п.н. конечна.

Последовательность несмещенных оценок {£г}°=о для решения п(х) задачи (2.1), будем называть допустимой, если существует последовательность а-алгебр {©¿}^=о таких, что С Bi и С ©¿+1, а имеет вид ^ = С,г + Хгп(хг), где С,г — ©¿-измерима. Для допустимой последовательности оценок определим случайную величину & равенством = Стг + хп(х*г), где х — индикатор события {т1 > т2}, а х—точка границы, отстоящая от точки хтг на расстояние не большее, чем 5. Справедлива

Теорема 4 ([3], теорема 2.3.2). Если допустимая последовательность оценок образует квадратично-интегрируемый мартингал относительно семейства а-алгебр {©¿}°=о, указанного в ее определении, то случайная величина является £(5)-смещенной для непрерывной функции п(х) (е(5) — модуль непрерывности функции п(х)), а ее дисперсия есть ограниченная функция параметра 5.

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

Лемма 1. Пусть уравнение (2.1) имеет ограниченные решения для правой части Г(х) и \Г(х)\, тогда стандартная последовательность несмещенных оценок образует квадратично-интегрируемый мартингал относительно потока а-алгебр

{»¿Шо.

Функцию Г(х), обычно, представляют в виде Г(х) = Н(х)Е/(У), где случайная величина У имеет распределение, зависящее от х, а функция /(у) является правой частью дифференциального уравнения, или значением его решения на границе. Пусть {Уз }^=о —последовательность случайных величин, таких что

Г(хг) = Н(хг)Е(/(уг) \ ») (2.5)

п.н. и ©г — минимальная а-алгебра, порожденная и последовательностью {уз-}^=о,

тогда последовательность несмещенных оценок

¿-1

л/ (у ,)Х.

& = Н(х3 )/(уз )Хз + X¿и(х) (2-6)

3=0

является допустимой. Такие оценки будем по традиции называть оценками по столкновениям. Если уравнение (2.1) имеет ограниченное решение й(х) для правой части Г(х) = Н(х), то из леммы 1 следует, что последовательность несмещенных оценок

¿-1

& = ^2 Ь(хз )Хз + X¿u(x¿) (2.7)

3=0

образует квадратично-интегрируемый мартингал относительно потока а-алгебр {^¿}°=0. Отсюда, очевидно, следует

Лемма 2. Пусть уравнение (2.1) имеет ограниченное решение для правой части Г(х) = Н(х). Тогда для правой части вида (2.5) с ограниченной функцией / (х) уравнение (2.1) также имеет ограниченное решение, и последовательность несмещенных оценок (2.6) образует квадратично-интегрируемый мартингал относительно потока а-алгебр {^¿}°=о.

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

Лемма 3. Пусть уравнение (2.1) имеет ограниченные решения для правой части Г(х) = Н(х) и существует постоянная о(6), такая что выполнено неравенство Н(х) > с(6) > 0 для всех х € (, удовлетворяющих условию р(х, д() > 6. Тогда справедливо неравенство Ехт§ < ОН(х)/е(6).

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

а2 (у,х) = (у - х)'А-1(у - х)

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

0/ 1 (а2 (у,х) г°{у -х,г,т) = -5-г-ехр - № '

[4тг(*-т)]т(Ле*А)* V 4(*-т)у '

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

Шароид (п определим равенством

(п = (к(х, I) = {(у, т)№0(у - х, 1, т) - (2пЕ2/п)-п/2(д.е1А)-1/2 > 0, т < 1} .

На сфероиде выполняется равенство

( Д2 \

*\у,х)=2п{1-г)ы(кщг—)у (3.1)

поэтому

Я2

—;-г > 1-

2п(Ь — т) -

Правую часть равенства (3.1) будем обозначать р2(т).

Из этих формул следует, что Ь — Я2/(2п) < т < Ь, причем гиперплоскости т = Ь и т = Ь — Я2/(2п) имеют пересечение со сфероидом в единственной точке, пространственные координаты которой равны х, а сечением сфероида гиперплоскостью т = с является эллипсоид, максимальный размер которого Д/л/ё.

Аналогично, при Ь — Я2/(2п) < 5 < Ь определяется усеченный шароид

QR,S = QR,s(x,t) = {(y,т)\2о(у — х,Ь,т) — (2пЯ2/п)-п/2(д.еЩх,Ь))-1/2 >0,5<т<ь} с нижней границей

1я,в = (х, Ь) = [(у,5)^ о(у — х,Ь,5) — (2пЯ2/п)-п/2 (йеЩх,Ь))-1/2 > 0} ,

имеющей проекцию I^ на координатное пространство. Определим функцию у(у, т) равенством

у(у,т )= 2 о(у — х,Ь,т) — (2пЯ2/п)-п/2(6.е1Л)-1/2.

Рассмотрим первую краевую задачу (1.1) для оператора с постоянными коэффициентами, содержащего только старшие производные. Зададимся функцией Я = Я(х,Ь), непрерывной в Qт и удовлетворяющей неравенствам

с^^х, Г) < Я(х,Ь) < с2^(х, Г)

при некоторых положительных постоянных С1 и С2 и одному из условий: QR С Qт или QRtо С Qт. Используя формулы Грина, нетрудно получить интегральное представление в шароиде (при Ь > Я2/(2п)):

I I (Р )

или интегральное уравнение в усеченном шароиде (при Ь < Я2/(2п)):

п(х,Ь) = у(у,т)/(у,т)<у<т + I у(у, 0)^(у)<у+

•'-Од,о

+ 1 д,тЕ(-У , ч-и(х + р(т)П,т), (3.3)

/п \2) Т(п 2 ¿-т Д™ ^ ' ' Ь У '

1о \2) Г(п/2)(Ь — т )Я где случайный вектор О распределен на единичном эллипсоиде

^ = {ш € Яп\штЛ-1 ш = 1}

с центром в нуле и имеет плотность распределения

р(х^,ш) =--. (3.4)

сг„ (д,еЛА(х, £))2 ЦА-^Ц

Величина

_ 2тг?

является площадью поверхности сферы радиуса 1 в Кп.

Таким образом, для решения и(х,Ь) первой краевой задачи получено интегральное уравнение с субстохастическим ядром. Формула (3.2) для уравнения теплопроводности получена другими методами Л. П. Купцовым [8].

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

р(т) = у/(2п(г-т) 1п(Д2/(2 п(г-т))) и выполняя замену переменной в = (п/2) \п(К2 / (2п(1 — т)), получаем равенства

рп(т)

<1тЕ (2 , ,Р}Т) и(х + р(т)П, г) = \2) Г(п/2)(£ — т)Яп у 'у ' ' у

Г ** (?) * - +

0 Л) Т(п/2)(г — т)

¿-(Д2/(2п)) ехр(-2в/п))Ж? =

Лп/2)1п(Д2/(2пе}) (п/2)г(п/2)

= Ех{-у>(п/2) 1п(ку(2ш))}и(х + Д,/7 ехр(-7/п)П, t - (Д2/(2п)) ехр(-27/п)),

где 7 = 7(1 + п/2) — случайная величина, имеющая гамма распределение с параметрами (1 + п/2, 1). Эта же оценка пригодна и для интеграла в (3.2), так как при Ь > К2/(2п) событие {7 > (п/2) Ы(Е,2/(2п))} является достоверным.

Процедура моделирования блуждания по сфероидам (хн^н) (к = 1, 2,...) на к-м шаге состоит из следующих действий:

1) моделируем случайную величину 7н с распределением гамма с параметрами (1+ п/2, 1);

2) проверяем условие ^н > (п/2)\п(КН_ 1 /(2пЬк~г));

3) если оно выполнено, то моделируем случайный вектор 0.н, распределенный с плотностью

р(ш) = -, -

на единичном эллипсоиде и вычисляем новые координаты по формулам

хк := + Дд;-1а/7й ехр(—7й/п)Пд;,

■= £н-1 — (Д|-1 /(2п)) ехр(—2чн/п);

4) в противном случае цепь обрывается.

0

оо

Из теорем 2 и 3 вытекает следующая лемма.

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

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

Р (х,Ь) = у(у,т)/(у,т )<у<т и Ро (х,Ь) = у(у,т)/(у,т )<у<т ¿Як ¿Як,о

непрерывны и обращаются в ноль лишь на нижнем основании и боковой поверхности пространственно-временного цилиндра. По теореме 2 заключаем, что с вероятностью единица блуждание по шароидам либо обрывается, либо приближается к боковой поверхности или нижнему основанию пространственно-временного цилиндра. Рассматривая решение краевой задачи с /(х,Ь) = 0 нулевым граничным условием и согласованным с ним положительным начальным условием, исключаем п.н. приближение необрывающейся траектории блуждания к нижнему основанию пространственно-временного цилиндра, так как в этом случае Р(х, 0) = <^(х, 0) > 0. Очевидно, что субстохастическое ядро, определяемое равенствами (3.2)-(3.3) удовлетворяет пункту 2 теоремы 3, что и завершает доказательство леммы. ■

Построим мартингал несмещенных оценок решения и(х, Ь) задачи (1.1). Для этого преобразуем оставшиеся интегралы в формулах (3.2)-(3.3). Очевидно,

[ , , Гр(о) , 2гп-1 „ехр(—г2/(4Ь)) — ехр(—р2(0)/(4Ь)) ,

/ у(у,0Ыу)<1у= / ¿г—-—— - , ^ /о "<р(х + г9).

}ьа,0 К '} Уо Г(п/2) (44)»/2 ^ ^ ;

Записав разность экспонент в виде интеграла и изменив порядок интегрирования в повторных интегралах, получим под интегралом плотность случайной величины 7, имеющей гамма-распределение с параметром ^ + 1. Следовательно,

/ ъ(у,0)ср(у)с1у = Ех{7<(п/2)\п(ну(2гч))}Ах + 2Р^^), (3-5)

¿Ь к,о

п1

где р — случайная величина, имеющая плотность распределения пг 1 на отрезке [0,1].

Используя теорему Фубини, находим

/ ъ(у,т)/(у,т )<у<т = ¿Як,о

= Г *г Г - ,))) - ехр( Р2(г)/т - г)))

Уо Уо Г(п/2) (4(4-г))»/2 п '

Оценивая внутренний интеграл, получаем формулу

/ т)/(у, т)(1у(1т = / (1тЕх{1<{п/2) 1п(д2/(2п(4-т)))}/(ж + 2ру/(г - т)г) =

¿Як,о ¿о

= tEx{1<{n/^2)\n{в?/(2п40))>/(х + -гв).

Точно так же оцениваем последний интеграл:

с д2

/ у(у, т)/(у, т)(1уд,т = —ЕХ{1<{п/2) 1п(1 /9))}1(х + 2рКу/01/(2п)П, t - Н2в/(2п)). ■¡О.ъ 2п

Здесь случайная величина в распределена равномерно на отрезке [0,1].

Пусть N1 — момент обрыва траектории блуждания, х, — индикатор события {N1 > г}. Положим = шт(^, Д2/(2п)) и обозначим ф, —индикатор события {7®+1 < (п/2)\п(К2/(2пт,в,+1))}. Рассмотрим последовательность несмещенных оценок (к = 0,1,...) для решения и(х,г) задачи (1.1), определяемую равенством

н1

i=0

+ Хки(хк,гк) + х{^=к}(р(хк-1 + Ярк^к-Пк^к)-

Пусть = а(в,, р,, 7,, г = 1, 2,..., к) — а-алгебра, порожденная независимыми случайными элементами, распределение которых определено выше. Из определения последовательности (к = 0,1, 2,...) следует, что она является мартингалом относительно потока а-алгебр , (к = 0,1, 2,...), порожденного блужданием. Из леммы 2 вытекает теорема

Теорема 5. Пусть функции /(х,г), <р(х) и Ф(х,г) ограничены. Тогда мартингал £н(к = 0,1, 2,. ..) квадратично интегрируемый.

Пусть 6 > 0, N2 — момент первого попадания траектории в ¿-окрестность границы области В, а Щ = тогда справелива следующая лемма.

Лемма 5. Случайная величина Щ имеет конечное математическое ожидание. Случайная величина является несмещенной оценкой и(х,г) и имеет конечную дисперсию.

Доказательство. Будем следить за изменением времени Ьн до момента времени N. Заметим, что на шаге к происходит уменьшение времени Ьн-1 на величину (Д|_1 /(2п)) вхр(—2^н/п). При этом Дн-1 > С16. Рассмотрим последовательность независимых случайных величин тн = (С16)2 вхр(—2^н/п)/(2п)(к = 1, 2,...). Эти случайные величины одинаково распределены и ограничены. Очевидно, что

N — 1 < N = штх{к ■ т1 + т2 + ... + тн < г}.

По теореме восстановления EN = г/Ет1 + о(г).

Второе утверждение является очевидным следствием равномерной интегрируемости мартингала = 0,1, 2,...). ■

4. Блуждание по шароидам для уравнения с переменным коэффициентом при неизвестной функции. Пусть Ах > 0 и Л2 > 0 удовлетворяют в С}т условиям А1 > —ао(х,г) и А2 > А1 + ао(х,г). Функция и,1(х,г) = е-^1*и(х,г) удовлетворяет уравнению Ьи1(х,г) + А1П1(х,г) = в-Хгг/(х,г), в котором коэффициент при неизвестной функции ао1(х,г) = А1 + ао(х,г) > 0. Аналогично, функция и2(х,г) = в^2*и1(х,г) удовлетворяет уравнению Ьи2(х,Ь) + (А1 — А2)и2(х,Ь) = в(Х2-Хг)г/(х,Ь), в котором

коэффициент при неизвестной функции ао2(х,Ь) = —Х2 + ао1(х,Ь) < 0. Перенося выражение ао2(х,Ь)и2(х,Ь) в правую часть уравнения и применяя формулы среднего значения (3.2)—(3.3), получаем интегральное представление для функции и1(х,Ь) в шароиде:

и1(х,Ь)= у(у,т )в-Х2(г-т )в-Х1 т / (у,т)<1у<1т+ ■'Як

+ I у(у, т)(\2 — (Л1 + ао(у, т)))е-Хз(—)и1(у, т)<у<т+

Як

Як

+ <1тЕ (§) * г(^ад-пд.+ '>• <">

или интегральное представление в усеченном шароиде (при Ь < К2/(2п)):

и1(х,Ь)= у(уу,т)в-Х2 (г-т )в-Х1Т / (у,т)<1у<1т+

Як, о

+ у(у,т)(Л2 — (Л1 + ао(у,т)))в-х2(г-т)и1(у,т)<у<т + у(у, 0)е-х*<р(у)<у+

¿Як,о •'С к,о

+С "тЕ (§) * шт^'"'''^+<4-2»

Очевидно, что ядро в каждом из этих представлений субстохастическое. Следовательно, для построения процесса блуждания и оценок на его траекториях снова можно использовать результаты § 2, которые их однозначно определяют. Для этого достаточно несмещенно оценить интегралы в формулах (4.1-4.2). Используя несмещенные оценки на траекториях блуждания по сфероидам, находим

= Ех{1>{п/2)\п{И2/(2пг))} ехр (—Л2(К2/(2п)) ехр(—2^/п)) х

х иг(х + Дл/7ехр(-7/п)П,4 - (Д2/(2п)) ехр(-27/п)),

где 7 = 1 + п/2). Эта же оценка пригодна и для аналогичного интеграла в (4.1), так как при Ь > К2/(2п) событие {7 > (п/2) \п(К2/(2пЬ))} является достоверным. Далее,

/ 0)е~Х2г1р(уу)(1у = Ех{~,<(п/2) 1п(пу(2ш))}е~Х2^{х + 2ру/руП), (4.3)

■' ¿к,о

где величины 7 и р имеют такие же распределения, как в формуле (3.5). Точно так же оцениваем интеграл

[ у(у, т)е-х2(г-т)е-х1т/(у, т)<у<т = Як,о

Наконец,

[ ю(у,т )в-Х2 (г-т )в-Х1Т / (у,т )йуйт

¿О.п

'О

т> 2

Простой заменой функции из уже полученных интегралов, находим / у(у, т)(Л2 - (Ах + ао(у, т)))в-Х2(г-Т^щ(у, т)<1у<1т =

I у(у, т)(Л2 - (Л1 + ао(у, т)))в-х2(г-т^щ(у, т)йуйт = ■'Ян

= ^-БХ{7<(п/2)1п(1/0))}(А2 - (А1 +а0(ж + 2рлД&уП, t ~Щ)х

х е-А2л2®/(2п)и1(х + 2рКл/Щ2п)П,г - К2в/(2п)).

Заметим, что {7 < (п/2)Ы(1/в))} = {в < ехр(—27/п)}, поэтому процедура моделирования блуждания по сфероидам и шароидам (хк) (к = 1, 2, .. .) на к-м шаге состоит из следующих действий:

1) моделируем случайную величину 7 к, имеющую гамма-распределение с параметрами (1 + п/2,1) и случайную величину щ = — \п(вк), имеющую показательное распределение;

2) моделируем случайный вектор Ок, распределенный с плотностью

р(ш) =

на единичном эллипсоиде;

3) проверяем неравенство щ > (2пехр(27/п));

4) если оно не выполнено, то переходим к (8);

5) проверяем условие ^к > (п/2)\п(В^_1 /(2пЬ/-1));

6) если оно выполнено, то вычисляем новые координаты по формулам

хк := + Дд;-1а/7й ехр(—

^ := ^-1 — (Я1-1 /(2п)) ехр(—2^к/п);

полученная точка лежит на сфероиде;

7) в противном случае цепь обрывается;

8) вычисляем новые координаты по формулам

Хк ■= хк-1 + 2рк ^щ^к/М^-к,

1к := ¿к-1 — Пк/Л2;

полученная точка лежит в шароиде;

9) процесс обрывается с вероятностью (Л1 + ао(хк))/Л2 в полученной точке шароида.

Несмещенные оценки для ui(x, t) на траекториях построенного блуждания легко строятся по аналогии с оценками для блуждания по сфероидам и обладают аналогичными свойствами.

Отметим, что возможны и другие замены неизвестной функции, приводящие к тождеству с неотрицательным субстохастическим ядром. Например, если ао(у, т) > 0, то можно использовать замену и\(у,т) = (1 — c(t — т))и(у,т), которая при выполнении условий c(t — т) < 1, c > m'AX(y Tао(у,т) приводит к интегральному уравнению с субстохастическим ядром в шароиде, который следует выбирать усеченным, если cR2/(2n) > 1.

Литература

1. Ладыженская О. А., Солонников В. А., Уральцева Н.Н. Линейные и квазилинейные уравнения параболического типа. М.: Наука, 1967.

2. Курбанмурадов О. А., Сабельфельд К. К., Симонов Н. А. Алгоритмы случайного блуждания по границе. Новосибирск, 1989.

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

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

5. Сипин А. С. Об особенностях применения схемы Неймана—Улама к решению краевых задач // Математические модели. Теория и приложения. Сборник научных статей под ред. Проф. М. К. Чиркова. Вып. 13. СПб.: Изд-во ВВМ НИИХ СПбГУ, 2012. C. 37-47.

6. Меер П. А. Вероятность и потенциалы. М.: Мир, 1973.

7. Doob J.L. Discrete potential theory and boundaries // J. Math. Mech. Vol. 8,3. 1959. P. 433-458.

8. Купцов Л. П. Свойства среднего и принцип максимума для параболических уравнений второго порядка. ДАН СССР. 242, №3. 1978. C. 56-59.

Статья поступила в редакцию 24 октября 2013 г.

Сведения об авторе: Сипин Александр Степанович — кандидат физико-математических наук, доцент; [email protected]

APPLICATION VON-NEUMANN—ULAM SCHEME TO THE SOLUTION OF THE FIRST BOUNDARY VALUE PROBLEM FOR A PARABOLIC EQUATION

Aleksandr S. Sipin

Vologda State Pedagogical University, ul. S. Orlova, 6, Vologda, 160035, Russian Federation; [email protected]

Statistical estimators for the solutions of boundary value problems for parabolic equations with constant coefficients are usually constructed on the random walk trajectories.The phase space of these random walk is the domain for the boundary value problem or the boundary of this domain. It is necessary to know the exact fundamental solution to construct the random walk. So, this algorithm cannot be used for heat equations when its coefficients are not constant.

In this paper we construct unbiased estimators and small biased estimators for the solution of the first boundary value problem for parabolic equations in case when the coefficient of the unknown function is sufficiently smooth and other coefficients are constant. We use the random walk on balloids and von-Neumann—Ulam scheme generalization for integral equation with sub-stochastic kernel to investigate the Markov chain properties and statistical estimators. The algorithm is based on a new integral representation for the solution of the boundary value problem. Refs 0. Figs 0. Tables 0.

Keywords: Monte Carlo methods,simulation method, heat equation, random walk on balloids.