Научная статья на тему 'О РЕШЕНИИ МНОГОМЕРНЫХ РЕЛАКСАЦИОННЫХ УРАВНЕНИЙ ДИФФУЗИИ МЕТОДОМ МОНТЕ-КАРЛО'

О РЕШЕНИИ МНОГОМЕРНЫХ РЕЛАКСАЦИОННЫХ УРАВНЕНИЙ ДИФФУЗИИ МЕТОДОМ МОНТЕ-КАРЛО Текст научной статьи по специальности «Математика»

CC BY
35
7
i Надоели баннеры? Вы всегда можете отключить рекламу.
Ключевые слова
РЕЛАКСАЦИОННОЕ УРАВНЕНИЕ МАССОПЕ- РЕНОСА / УРАВНЕНИЕ ДИФФУЗИИ / DIFFUSION EQUATION / ВРЕМЯ РЕЛАКСАЦИИ / RELAXATION TIME / МЕТОД МОНТЕ-КАРЛО / MONTE CARLO METHOD / ФОРМУЛА КИРХГОФА / THE KIRCH- HOFF FORMULA / RELAXATION EQUATION OF MASS TRANSFER

Аннотация научной статьи по математике, автор научной работы — Сиренек Валерий Анатольевич

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

i Надоели баннеры? Вы всегда можете отключить рекламу.
iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.
i Надоели баннеры? Вы всегда можете отключить рекламу.

ON THE SOLUTIONS OF MULTIDIMENSIONAL DIFFUSION EQUATIONS RELAXATION MONTE-CARLO

The work is devoted to computational experiments (on model prima tures) associated with the decision by the Monte Carlo of the Cauchy problem for the three-dimensional telegraph equation. Application of the Monte Carlo method, but on the basis of a probabilistic representation of the solution in the form of mathematical expectation. Two variants of this rep- resentation. To verify the proposed method solutions found several functions satisfying-create the original equation. The calculation results are compared with analytic-cal results. The cases of both permanent and change-coefficients. The pro- posed method provides a high enough accuracy for practice.

Текст научной работы на тему «О РЕШЕНИИ МНОГОМЕРНЫХ РЕЛАКСАЦИОННЫХ УРАВНЕНИЙ ДИФФУЗИИ МЕТОДОМ МОНТЕ-КАРЛО»

УДК 541.12.012.2: 539.213.3

V.A. Sirenek

on the solutions of multidimensional diffusion equations relaxation monte carlo

St. Petersburg State Institute of Technology (Technical University), Moskovsky Pr., 26, St Petersburg, 190013, Russia e-mail: wasirenek@gmail.com

The work is devoted to computational experiments (on model prima tures) associated with the decision by the Monte Carlo of the Cauchy problem for the three-dimensional telegraph equation. Application of the Monte Carlo method, but on the basis of a probabilistic representation of the solution in the form of mathematical expectation. Two variants of this representation. To verify the proposed method solutions found several functions satisfying-create the original equation. The calculation results are compared with analytic-cal results. The cases of both permanent and change-coefficients. The proposed method provides a high enough accuracy for practice.

Keywords: relaxation equation of mass transfer, diffusion equation, relaxation time, the Monte Carlo method, the Kirchhoff formula.

й01 10.15217/1ззп1998984-9.2015.30.81

Введение

При моделировании диффузионных процессов с учетом релаксационных явлений (эффекта «запаздывания») массопереноса широко используются гиперболические уравнения-релаксационные модели [1]. Из них базовое (характерное) имеет форму телеграфного уравнения. Гиперболические уравнения, в отличие от параболических с постоянными коэффициентами (в том числе от классического уравнения диффузии Фика), учитывают конечность скорости распространения возмущений переносимой субстанции. Их решения позволяют адекватно определить как интегральные характеристики массопереноса (количество переносимого вещества, ширину диффузионной зоны), демонстрирующие «запаздывание» переносимых потоков, так и дифференциальные характеристики, т.е. профили концентрации с резко очерченным фронтом, проявляющим себя на начальной (релаксационной) стадии процесса. Классическое уравнение диффузии адекватно описывает стадию развитой диффузии с «размытым» фронтом профиля. Случаи успешного использования многомерных (по пространству) гиперболических уравнений в качестве релаксационных моделей тепломассопереноса, в отличие от использования одномерных уравнений этого типа, менее известны. Недостаточно развит и арсенал практических способов решения многомерных уравнений.

При решении гиперболических уравнений массопереноса универсальным методом конечных

В.А. Сиренек1

о решении многомерных

релаксационных уравнений диффузии

методом монте-карло

Санкт-Петербургский государственный технологический институт (технический университет), Московский пр., 26, Санкт-Петербург, 190013, Россия e-mail: wasirenek@gmail.com

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

Ключевые слова: релаксационное уравнение массопереноса, уравнение диффузии, время релаксации, метод Монте-Карло, формула Кирхгофа.

разностей происходит искажение крутого фронта профиля концентрации - эффекты диссипации и дисперсии, а также существуют трудности с выбором размеров конечно-разностной сетки [2]. Эффективным средством избежать подобных проблем служат основанные на методе Монте-Карло численно-вероятностные способы решения дифференциальных уравнений, не обладающие погрешностью, обусловленной дискретизацией задачи, и малочувствительные к её размерности. Решение в этом случае сводится к осреднению некоторого функционала от реализаций моделируемого случайного процесса. Обычно используются два типа случайных процессов - имитационные и фиктивные [3]. Представляющие собой модели «случайного блуждания» и лежащие в основе вывода дифференциальных уравнений массопереноса, случайные имитационные процессы позволяют не только решать эти уравнения методом Монте-Карло, но и изучать вероятностную природу характеристик массопереноса. Изучению связи одномерных гиперболических уравнений и их решений со случайными имитационными процессами посвящены работы [4-7]. Однако, для трехмерных гиперболических уравнений (в отличие от эллиптических и параболических уравнений) не существует порождающих их моделей «случайного блуждания» [8], и поэтому успех решения таких уравнений методом Монте-Карло определяется удачным выражением решения через математическое ожидание некоторой случайной величины. Так как выигрышным в вычислительном отношении при решении какой-либо задачи математической физики считается

1 Сиренек Валерий Анатольевич, канд. техн. наук, доцент каф. системного анализа, e-mail: wasirenek@gmail.com Sirenek Valeriy A., PhD (Eng.), Department of systems analysis, e-mail: wasirenek@gmail.com

Дата поступления - 17 февраля 2015 года Received February 17, 2015

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

д2С _ дС —- + 2a— = v dt2

2 д2С д2С д2С

Т +-Г +-Т

dt { дх2 ду2 dz2

v2 = D/т*, a =1/(2т*),

(1)

где ^х^) - концентрация переносимого вещества, х = {x,y,z}, т* - время релаксации среды, D - эффективный коэффициент диффузии, v - скорость распространения концентрационных возмущений. В случае решения задачи Коши к уравнению (1) следует добавить начальные условия:

С(X, у, z, 0) = Ф(х, у, — С(X, у, 0) = ¥(X, у, 7) (2)

&

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

О релаксационном характере диффузии примесей в водоемах

Показательным примером использования многомерных гиперболических уравнений может служить диффузионная модель распространения примесей в природных водоемах. Изучение влияния вредоносных стоков на качество воды является сегодня одной из важнейших проблем народного хозяйства. Так, нефть и продукты ее переработки, став в ХХ веке основным энергетическим ресурсом, превратились в самое распространенное загрязняющее вещество. Значительная роль в решении этой проблемы принадлежит диффузии, так как в основном диффузия обуславливает распределение промышленных и бытовых стоков в водоемах. Особенностью изучаемых процессов и явлений в воде является турбулентный характер миграции примесей. При этом использование традиционного (полуэмпирического) параболического уравнения турбулентной диффузии наталкивается на ряд принципиальных трудностей. Такие модели выводятся путем осреднения уравнения сохранения вещества в неподвижной среде для относительно больших времен (порядка десятков суток) и больших водных акваторий, т.е. по всем неоднородностям, существующим в водоеме. Информационная ценность такой модели мала в тех случаях, когда сброс стоков и забор воды происходит в прибрежной полосе водоемов, а потребители воды оценивают ее качество в момент употребления. Из опыта хорошо известно, что в окрестности источника воды область, занимаемая примесью, имеет четко выраженный фронт, который размывается лишь со временем. Следовательно, с точки зрения потребителей воды наиболее приемлемым методом описания распределения примеси будет тот, который дает концентрации примеси, осредненные по небольшому периоду времени и, соответственно, по относительно небольшому объему [9, 10]. Уменьшение периода осреднения уравнения сохранения вещества в трехмерном пространстве приводит, вообще говоря, к интегро-дифференциальному уравнению [11]. В частном случае (однородного и изотропного пространства) приходят к уравнению (1).

Вероятностное представление решения

Несмотря на то, что для трехмерного уравнения (1), в отличие от одномерного телеграфного уравнения, не существует порождающего его «случайного блуждания» [8], известны два способа вероятностного представления решения задачи Коши (1), (2) через математическое ожидание «рандомизированного по времени» точного решения волнового уравнения. При этом «рандомизация», имеющая в одномерном случае ясный физический смысл, инвариантна относительно размерности уравнения.

Способ 1. Теоретически доказано [12], что решение задачи (1), (2) имеет вероятностное представление, аналогичное представлению решения одномерной задачи (с той же «рандомизацией времени»):

С(x, t) = E {S(vt )Ф(x)} + E j J S(vr)¥(x)dт

(3)

где Е{-} - операция математического ожидания; 5(7)Ф(х) - операторная функция Кирхгофа, являющаяся аналитическим решением задачи Коши для трехмерного волнового уравнения (уравнения (1) с а = 0 и v = 1) с частными начальными условиями (Ф(х) = 0) и имеющая вид:

(4)

Величина t в (3) - «рандомизированное время» (случайная величина):

J = J (-1)"a (s)ds ,

(5)

где Na(s) - случайный (скачкообразный) процесс Пуассона интенсивности a. Надо заметить, что при аналогичном вероятностном представлении решения одномерного телеграфного уравнения величина Na(s) имеет конкретный физический смысл - она определяет число перемен направления движения у блуждающей на прямой со скоростью v фиктивной частицы за время s. Причем эта частица меняет направление движения с интенсивностью (частотой) a. Через замену переменной т = tr3 во втором слагаемом выражения (3) получаем интеграл по г от 0 до 1. Таким образом, (3) представляет собой сумму интегралов по единичной сфере и по единичному шару. Введем случайный вектор w, равномерно распределенный по единичной сфере: w = {sine • cosj, sine • sinj, cose}, где 0 < ф < 2п, 0 < 0 < п, и оператор дифференцирования Vw по направлению w. Формула (3) получает представление в виде двойного математического ожидания:

+t Y(x + Vtrсо) + vt У V й^(х + v/r со)}};

(6)

где ЕТ{...} и Е(г,ю){...} - математические ^жидания относительно независимых случайных величин t и (г,ю) .

Способ 2. В этом способе формула Кирхгофа используется непосредственно. Согласно [4, 5], решение задачи (1), (2) может быть представлено в виде математического ожидания от «рандомизированного по времени» аналитического решения той же задачи для волнового уравнения:

C(x, t) = E{C*(x, t )},

(7)

где С*(х,0 - точное решение волнового уравнения [уравнения (1) при а = 0] с начальными условиями общего вида (2) в виде формулы Кирхгофа [13]:

C*(x,t)=— } |[Ф(„) + v/(<Ь'х („) • sinG • coscp + Ф'у („) • sinG • siny + 471 о о

+Ф; („) • COS0) + f*F(„)]sine d 0 í/ф,

(8)

где (,,) обозначает (х + уйт0- ^ф, у + Уйт0- smф, 2 + ^^0). Численной реализации формулы (7) посвящены работы [14-16]. С учетом ш и Уш выражение (7) получает представление в виде двойного м атематического ожидания:

где Eт {...} и Eш{...} - операции внешнего и внутреннего матем атического ожидания относительно независимых случайных величин Т и ш.

Основные расчетные формулы

Вероятностный вид выражений (6) и (9) является предпосылкой для их численной реализации методом Монте-Карло:

(10)

+ 4>(x+vtmk))]

tj «I^ИГ1, /í=(-ln ф/a

при условии

I ti=t

(11)

(12)

где {ti} - экспоненциально распределенные независимые случайные величины с параметром a; {р,} - независимые случайные величины, равномерно распределенные на (0, 1); Wk ={sin 0k cos фк; sin 0k sin фк; cos 0k }. Известно [17], что значения Гк, 0к и фк следует вычислять по формулам: r = , Ф = 2па2, cos0 ■= 2 аз-1, где ai, а2, аз - независимые случайные величины, равномерно распределенные на (0,1); M и L - число реализаций соответственно t и (r, w) (в способе 1) или t и w (в способе 2).

Следует заметить, что в случае at < 1 большая часть реализаций процесса Пуассона на промежутке [0,t] не имеет скачков, что значительно понижает эффективность метода Монте-Карло. В связи с этим для повышения точности вычислений автором предложен прием, использующий условные математические ожидания. Выделим внутреннее математическое ожидание в формулах (6) и (9). Это соответственно n(t) = E(r, «){■■■} и n( t) = Ew{---}. Введем два против о полож ны х случайных события: A ={Na(t) = 0} и B = { Na(t) > 1}. Вероятности этих событий Р(Д)=ехр(- at), P(B) =1- exp(- at). Получаем:

C(x,í) = P (A) ■ ET{ it) | A) + P(B) ■ Ef{ H) IB) (13)

В первом слагаемом E T {n(f )| A} = n(t) -детерминированная величина, которую предлагается вычислять по формуле Симпсона. Во втором слагаемом - в случае события B - время до первого скачка процесса Пуассона tj имеет условную плотность вероятности: p(s)

= a exp(-as )/[1 - exp(-at )], 0 < s < t. В соответствие с этим значение tj определяется по формуле:

t1 = -(1 / a) ln(1 - r (1 - exp(-at))),

(14)

г - случайная величина, р. р. на (ОД). Другие г (/' > 1} вычисляются по формуле (12). Значение Е7ЦО\В} вычисляется по формуле (10) или (11).

В работе [16] доказано, что при ограниченных функциях Ф(х) и Ф(х) дисперсия Dп(t) конечна для конечных значений t. Такой результат дает основание для оценки погрешности в методе Монте-Карло [17]. Введем среднее значение величины п(^):

м

Пм = (1/ м )!л,,

}=1

тогда при достаточно больших значениях М имеет место оценка вероятности

Р{ | Пм - С(х,1) |< ^Щ(Т) /М}=Це) , (15)

где Це) = (2 />/2Л){0е ехр(-^2 / 2)&- интеграл вероятностей. Эта формула имеет целое семейство оценок, зависящих от параметра £. Задаваясь коэффициентом доверия р, можно найти £р - корень уравнения L(г) = р. Так, значению р = 0.997 соответствует £ = 3 (правило трех сигм). Несмещенную эмпирическую оценку величины Г)ц(t) можно получить в процессе вычисления С(хД) по формуле

1 M о 1

Dn(t ) = 777 7ГI (Пj)2 -

(M -1) р

M (M -1) j

M

I л,

(16)

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

Численный эксперимент

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

Случай постоянных коэффициентов

Пример 1. Аналитическим решением уравнения (1) является функция:

J (>£

2 „2

C(X, y, z, t) = (cos X + cos y + cos z)e'(Va v a) с начальными условиями типа (2):

(17)

<b(x,y,z) = COSX + COS J + cosz,

Ф(x, y, z) = (cos x + cos y + cos z)(Va2 - v2 - a) (18)

В уравнении (1) приняты значения коэффициентов a = 0,5, v = 0,1. Реализован первый способ вероятностного представления решения задачи Коши (1),(2), т.е. использованы формулы (13) и (10). Результаты расчета - в таблице 1. Время T = 2at - относительное. Время расчета каждого значения по методу Монте-Карло = 1,5 мин, M = 10000, L = 1000.

2

=1

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

Таблица 1. Результат решения уравнения (1) с начальными условиями (18).

Т х = у = z

0,1 0,5 1 5 10

0,01 2,984735 (2,984710) [8,3Е-4%] {1,24Е-5} 2,632447 (2,632481) [1,3Е-3%] {1,05Е-5} 1,62093 (1,62074) [0,011%] {8,09Е-6} 0,85106 (0,85090) [0,018%] {6,75Е-6} -2,517016 (-2,51696) [2,2Е-3%] {1,14Е-5}

0,1 2,981783 (2,98200) [7,3Е-3%] {4,05Е-4} 2,630042 (2,63009) [1,8Е-3%] {3,52Е-4} 1,6182 (1,6192) [0,062%] {2,47Е-4} 0,84984 (0,85012) [0,034%] {2,34Е-4} -2,51448 (-2,51467) [7,4Е-3%] {3,66Е-4}

0,5 2,97030 (2,96997) [0,011%] {4,1Е-3} 2,6206 (2,6195) [0,041%] {4,16Е-3} 1,6118 (1,6127) [0,063%] {2,8Е-3} 0,8413 (0,8467) [0,64%] {2,39Е-3} -2,5055 (-2,5045) [0,04%] {4,18Е-3}

1 2,95480 (2,95501) [7,2Е-2%] {1,2Е-2} 2,6098 (2,6063) [0,13%] {1,07Е-2} 1,60359 (1,60461) [0,064%] {6,68Е-3} 0,8535 (0,8424) [1,3%] {6,17Е-3} -2,49205 (-2,49191) [0,1%] {1,01Е-2}

5 2,8325 (2,8380) [0,19%] {0,136} 2,4821 (2,5031) [0,84%] {0,124} 1,5535 (1,5411) [0,81%] {7,13Е-2} 0,8048 (0,8091) [1,52%] {5,04Е-2} -2,3 (-2,27) [1,1%] {0,276}

Примечания:результат, полученный по формулам (13), (10), дан без скобок; в круглых скобках - аналитический результат; в квадратных скобках дано относительное отклонение результата, полученного по формулам (13), (10), от аналитического результата;

в фигурных скобках - среднеквадратическое отклонение

Пример 2. Аналитическим решением уравнения (1) является функция:

С(х,у, z,г) = ехр[(х + у + 7)а / V + аг] (19)

с начальными условиями типа (2):

Ф(х,у,г) = ехр|~(х + у + г)а! у],

Ч/(х,у,г) = аехр\\х + у + г)а/г~\ (20)

Были использованы формулы (13) и (10), а = 0,5, V = 1, T = 2а. Результаты расчета приведены в таблице 2. М = 10000, L = 1000.

Таблица 2. Результат решения уравнения (1) с начальными условиями (20)

Т х = у = z

0,1 0,5 1 5 10

0,01 1,16671 (1,16766) [0,081%] {2,42Е-4} 2,12592 (2,12760) [0,079%] {4,36Е-4} 4,5102 (4,5042) [0,13%] {9,62Е-4} 1817,168 (1817,110) [3,2Е-3%] {0,373} 3286550 (3285400) [0,035%] {702}

0,1 1,2115 (1,2214) [0,81%] {7,4Е-3} 2,2105 (2,2255) [0,68%] {1,34Е-2} 4,7435 (4,7115) [0,68%] {3,03Е-2} 1899,24 (1900,74) [0,079%] {12,3} 3,4147Е+6 (3,4366Е+6) [0,64%] {2,19Е+4}

0,5 1,448 (1,491) [2,9%] {7,46Е-2} 2,818 (2,718) [3,7%] {0,143} 5,682 (5,755) [1,3%] {0,296} 2198 (2321) [5,3%] {136} 4,355Е+6 (4,198Е+6) [3,8%] {2,25Е+5}

1 1,986 (1,916) [3,7%] {0,266} 3,447 (3,490) [4,2%] {0,472} 7,612 (7,389) [3%] {1,05} 3052 (2980) [6,4%] {365} 5,685 (5,390Е+6) [5,5%] {6,54Е+5}

См. примечания к таблице 1

Случай переменных коэффициентов

Если коэффициент а в уравнении (1) зависит от времени, т. е. а = а(Ь), то необходимо моделировать неоднородный процесс Пуассона. Таким образом, значения Ь в формулах (12), следует вычислять с использованием метода, предложенного в [18]. Этот метод требует, чтобы функция а(Ь) была неотрицательна, ограничена сверху на [0;да), ее наибольшее значение на [0;да) не было бы слишком большим и интеграл от нее был равен да. В противном случае объем вычислений сильно возрастает. Этим условиям удовлетворяет:

а(Ь)=1/(Ь+1) (21)

Решением уравнения (1) с учетом (20) является функция:

С (х, у, 7, г) = х2 + у2 + 72 + V2 (г +1)2 (22)

с начальными условиями:

Ф( х, у, 7) = х2 + у2 + 72 + V2, ^ (х, у, 7) = 2v2 (23)

В уравнении (1) принималось значение V = 1. Время Ь в таблице 3 - абсолютное. Время расчета каждого из значений - около минуты.

Таблица 3. Результат решения (1) с учетом (21) и с условиями (23).

Ь х=у=т

0,1 0,5 1 5 10

0.01 1,04995 (1,05010) [0,014%] {3,62Е-3} 1,76985 (1,77010) [0,014%] {4,25Е-3} 4,020155 (4,02010) [1.4Е-3%] {2,92Е-3} 76,0184 (76,0201) [2,2Е-3%] {1,21Е-2} 301,018 (301,020) [6,6Е-4%] {3,08Е-2}

0.1 1,251 (1,240) [0,9%] {4,51Е-2} 1,9569 (1,9600) [0,16%] {7,14Е-2} 4,203 (4,210) [0,18%] {8,68Е-2} 76,2186 (76,2100) [0,011%] {0,144} 301,178 (301,210) [0,011%] {0,266}

0.5 2,237 (2,280) [1,9%] {0,816} 3,1 (3,0) [1,4%] {0,735} 5,220 (5,250) [0,56%] {0,782} 77,13 (77,25) [0,15%] {0,951} 302,41 (302,25) [0,053%] {1,32}

1 4,039 (4,030) [2,23%] {2,34} 4,637 (4,750) [2,4%] {2,35} 7,022 (7,000) [0,31%] {2.26} 79,10 (79,00) [0,43%] {2,52} 304,36 (304,00) [0,12%] {3,3}

5 36,73 (36,03) [3,9%] {30,5} 42,4 (36,75) [15%] {30,2} 37,77 (39,00) [3,2%] {31,5} 112,4 (111,0) [1,3%] {31,5} 335,54 (336,00) [0,14%] {32,1}

См. примечания к таблице 1.

Замечание. Следует заметить, что формула (9) гораздо проще формулы (6). Фактически, в формуле (6), по сравнению с (9), функция заменяется на интеграл от ее производной по конечному промежутку времени и этот интеграл «рандомизируется». Известно, что такая операция может привести к увеличению дисперсии. Однако, как показали расчеты, использование более простой формулы (9) (в вычислительной реализации (11)) дает значения результатов и их среднеквадрати-ческих отклонений, близкие к значениям, рассчитанным по формуле (6) (соответственно (10)) (см. таблицу 4). Среднеквадратическое отклонение ап( г) зависит от значений М и L неявно. Например, при М = 100 и L = 200 получались результаты, сравнимые с результатами, полученными при М = 10000, L = 1000 и приведенными в таблицах.

Таблица 4. Результат решения уравнения (1) с условиями (18) двумя способами - по формулам (10) и (11), М=10000,1=1000.

x=y=z=1; T=0,5 x=y=z=1; T=1

По формуле (10) 1,6118; (1,6127); 1,6036; (1,6046);

{0,0028} {0,0067}

По формуле (11) 1,6139; (1,6127); 1,6028; (1,6046);

{0,0011} {0,0041}

Заключение

Как показал вычислительный эксперимент предложенный способ реализации методом Монте-Карло вероятностных представлений решения трехмерного релаксационного уравнения диффузии типа телеграфного дает достаточно высокую для практики точность. С ростом значения времени увеличивается погрешность метода. При этом следует иметь в виду, что при относительно больших временах (T > 10) диффузионные процессы могут быть достаточно хорошо описаны параболическими уравнениями. Показано также, что данный способ решения уравнения (1) можно использовать и для случая a = a(t). Расчеты проводились на контрольных примерах, не связанных с конкретными диффузионными процессами. Исследование лишь задачи Коши не ограничивает применимость данного способа при решении краевых задач. Так, хорошо известны приемы сведения краевых задач для волнового уравнения к задаче Коши [метод распространяющихся волн - см. в [13]). Развитие рассмотренного способа решения на случай краевых задач для одномерного телеграфного уравнения успешно проведено автором в [6, 7]. В случае краевых задач для многомерного телеграфного уравнения такой подход может быть реализован лишь при наличии аналитических решений краевых задач для многомерного волнового уравнения.

Литература

1. Таганов И.Н, СиренекВ.А. Волновая диффузия СПб.: НИИХ СПбГУ, 2000. С. 107-209.

2. Абиев Р.Ш. Вычислительная гидромеханика и тепломассообмен. СПб.: 2002. 576 с.

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

4. Кац М. Несколько вероятностных задач физики и математики. М.: Наука, 1967. 176 с.

5. Hersh R. Stochastic solution of hyperbolic equations // Lect. Notes in Math. 1975. N446. P.283-300.

6. Сиренек В.А., Крючков А.Ф. Вероятностное решение краевой задачи для гиперболического уравнения массопереноса // Математическое моделирование. 1998. Т. 10. № 6. С. 107-117.

7. Сиренек В.А. Связь релаксационных уравнений массопереноса с моделями случайного блуждания и способы решения диффузионных задач на их основе // Известия СПбГТИ(ТУ). 2014. № 23(49). С. 93-96.

8. Polya G. Uber eine Aufgabe der Wahrschein-lichkeitrechung betreffend die Irrfahrt im Straennetz // Math. Ann. 1921. Bd. 84. S. 149-160.

9. Галкин Л.М. Решение диффузионных задач методом Монте-Карло. М.: Наука, 1975. 96 с.

10. Монин А.С., Яглом А.М. Статистическая гидромеханика. М.: Наука, 1965. 320 с.

11. Лыков А.В. Тепломассообмен. Справочник. М.: Энергия, 1972. 560 с.

12. Kisynski J. On M. Kac's probabilistic formula for the solution of the telegraphist's equations // Ann. Polonici Math. 1974. V. 29. P. 259-272.

13. Тихонов А.Н., Самарский А.А. Уравнения математической физики. М.: Наука, 1978. 735 с.

14. Кабанихина И.И. Численная реализация вероятностного представления решения задачи Коши для телеграфного уравнения: сб. науч. тр. Методы статистического моделирования. Новосибирск: ВЦ СО АН СССР, 1986. С. 86-90.

15. Веселовская А.З. О решении первой краевой задачи для волнового уравнения методом Монте-Карло: сб. науч. тр. Методы Монте-Карло в вычислительной математике и математической физике / под ред. Г.И Марчука Новосибирск: ВЦ СО АН СССР, 1974. С. 128135.

16. Крючков А.Ф., Сиренек В.А. Численно-вероятностные методы решения трехмерного гиперболического уравнения диффузии // Журн. вычисл. матем. и матем. Физики. 1998. Т. 38. N 1. С. 103-110.

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

18. Coleman, W.A., Mathematical verification of certain Monte Carlo sampling technique and applications of the technique to radiation transport problems // Nucl. Sci. and Engng. 1968. V. 32. P. 78-81.

i Надоели баннеры? Вы всегда можете отключить рекламу.