Диссертация (1145371), страница 16
Текст из файла (страница 16)
Проверяем условие γk > (n/2) ln Rk−16. Если оно выполнено, то вычисляем новые координаты по формулам:√xk := xk−1 + Rk−1 γk exp(−γk /n)Ωk ,2tk := tk−1 − Rk−1/(2n) exp(−2γk /n),Полученная точка лежит на сфероиде.7. В противном случае цепь обрывается.8. Вычисляем новые координаты по формулам:xk := xk−1 + 2ρkpηk γk /λ2 Ωk ,tk := tk−1 − ηk /λ2 .Полученная точка лежит в шароиде.9. Процесс обрывается с вероятностью λ1 + a0 (xk , tk ) /λ2 в полученной точкешароида.Несмещенные оценки для u1 (x, t) на траекториях построенного блуждания легко строятся по аналогии с оценками для блуждания по сфероидам иобладают аналогичными свойствами.Отметим, что возможны и другие замены неизвестной функции, приводящие к тождеству с неотрицательным субстохастическим ядром.
Например,если коэффициент a0 (y, τ ) > 0, то можно использовать замену u1 (y, τ ) == 1 − c(t − τ ) u(y, τ ), которая при выполнении условий c(t − τ ) 6 1,123c > max(y,τ )∈QR a0 (y, τ ) приводит к интегральному уравнению с субстохастическис ядром в шароиде, который следует выбирать усеченным, еслиcR2 /(2n) > 1.2.3.5Алгоритмы, связанные с дискретизацией времениЗаменяя производную по времени в параболическом уравнении ее раз-ностной аппроксимацией, можно свести решение первой краевой (2.3.1) задачик решению последовательности краевых задач для уравнения эллиптическоготипа. Рассмотрим данный метод на примере неявной (по времени) разностнойсхемы Эйлера.
Не умаляя общности, можно считать, что c0 > a0 (x, t) > c > 1.Выполнения этого условия всегда можно добиться экспоненциальным преобразованием неизвестной функции u(x, t).Пусть τ — шаг по времени, tk = kτ , k = 0, 1, 2, . . . , N и T = N τ . ПустьMk — эллиптическая часть параболического оператора при t = tk , uk (x) == u(x, tk ), fk (x) = f (x, tk ), Φk (x) = Φ(x, tk ), и, наконец, uek (x) — приближениедля u(x, tk ), удовлетворяющее системе разностных уравненийMk uek (x) = fk (x) −uek (x) − uek−1 (x),τx ∈ D,k = 1, 2, . . . , N,(2.3.58)краевым условиямuek (x) = Φk (x),x ∈ ∂D,k = 1, 2, .
. . , N(2.3.59)и начальному условиюue0 (x) = ϕ(x),x ∈ D.(2.3.60)Погрешность аппроксимации разностной схемы (2.3.58) определяется равенствомψk (x) = fk (x) −uk (x) − uk−1 (x)− Mk uk (x),τx ∈ D,k = 1, 2, . . . , N, (2.3.61)124и имеет порядок hα/2 , если решение краевой задачи (2.3.1) является элеменαтом пространства H 2+α,1+ 2 (QT ), (0 6 α < 1). Достаточные для этого условияможно найти в [28](гл. IV, теорема 5.2). Погрешностьzj (x) = uej (x) − uj (x),удовлетворяет разностным уравнениямMk zk (x) = ψk (x) −zk (x) − zk−1 (x),τx ∈ D,k = 1, 2, .
. . , N,(2.3.62)с однородными граничными и начальным условиями.Единственность решения разностно-дифференциальных уравнений с заданными граничными и начальным условиями следует из принципа максимумадля эллиптических уравнений, который справедлив в силу неравенстваa0 (x, t) + 1/τ > 0,(x, t) ∈ QT .Определим оператор Mkτ равенствомMkτnXnX∂2∂=−aij (x, tk )+ai (x, tk )+ (c0 + 1/τ ),∂x∂x∂xijii,j=1i=1(2.3.63)Запишем уравнение (2.3.58) в виде1ek−1 (x) + (c0 − a0 (x, tk ))euk (x),Mkτ uek (x) = fk (x) + uτx ∈ D,(2.3.64)k = 1, 2, . .
. , N,Разностно-дифференциальные уравнения легко преобразуются в разностноинтегральные в результате применения какой-либо теоремы о среднем значениидля эллиптического оператора. Воспользуемся алгоритмом, решающим первуюкраевую задачу для эллиптического уравнения, основанном на функции Левииз главы 3. Пусть Lk (y, x) — функция Леви для оператора Mkτ , построенная для125области Tk = Tk (x) ⊂ D, и удовлетворяющая условиямLk (y, x) = 0,∂Lk (y, x)= 0,∂yiNy L(y, x) > 0,i = 1, 2, . .
. , n,y ∈ ∂Tk (x);(2.3.65)Lk (y, x) > 0 y ∈ Tk (x) \ ∂Tk (x),тогдаZuek (x) =Tk1ek−1 (y) + c0 − a0 (y, tk ) uek (y) +Lk (y, x) fk (y) + uτ!τ+uek (y)Nk,yLk (y, x) dy. (2.3.66)Здесь, как обычно, Nkτ — оператор, формально сопряженный к оператору Mkτ .Выражениеτpk (x, y) = (c0 + 1/τ )Lk (y, x) + Nk,yLk (y, x)является по переменной y плотностью вероятностей в области Tk (x). Запишемего в виде смеси плотностейpk (x, y) = αk (x)pk,1 (x, y) + 1 − αk (x) pk,2 (x, y),τпропорциональных Lk (y, x) и Nk,yLk (y, x).
При этомZαk (x) = (c0 + 1/τ )Lk (y, x)dy.Tk (x)Используя представление (2.3.66), получим несмещенную оценку ξk (x)для uek (x). Пусть Yk,1 имеет плотность распределения pk,1 (x, y), а Yk,2 имеетплотность распределения pk,2 (x, y). Пусть α и β — случайные величины, распределенные равномерно на отрезке [0, 1] и независимые в совокупности с Yk,1и Yk,2 . Несмещенную оценку ξk (x) проще всего определить выписав реализующий её условный оператор126If α > αk (x) Then ξk (x) := uek (Yk,2 )ElseIf β(c0 + 1/τ ) < 1/τ Then ξk (x) := uek−1 (Yk,1 )ElseIf β(c0 + 1/τ ) < 1/τ + a0 (Yk,1 , tk ) Then ξk (x) := fk (Yk,1 )/a0 (Yk,1 , tk )Else ξk (x) := uek (Yk,1 )Рассмотрим теперь множество Q = {0, 1, . . . , N + 1} × D и определим нанем однородную цепь Маркова, для которой плотность вероятности переходаp (k, x) → (m, y) относительно прямого произведения “считывающей” меры νи меры Лебега λ определяется предыдущими построениями.
Таким образом, заодин шаг по времени возможны следующие изменения в состоянии цепи.1. С вероятностью 1 − αk (x) переходим в точку (k, Yk,2 ).2. С вероятностью αk (x) переходим в точку y = Yk,1 , а величина m принимаетзначенияk,k − 1,N +1с вероятностями(c0 − a0 (x, tk ))/(1/τ + c0 ),(1/τ )/(1/τ + c0 ),a0 (x, tk )/(1/τ + c0 ).Все состояния вида(0, x),(N + 1, x),x ∈ D;(k, x),x ∈ ∂Dявляются поглощающими.Введя обозначение ue(k, x) = uek (x), получаем из (2.3.66) интегральноеуравнениеZue(k, x) =Qp (k, x) → (m, y) ue(m, y)dµdλ + F (k, x),(2.3.67)127гдеZLk (y, x)fk (y)dλ.F (k, x) =TkВыбирая в качестве областей Tk (x) систему эллипсоидов, построенную вглаве 3 для решения первой краевой задачи для эллиптического уравнения иfk (x) ≡ 1, Φ(x) ≡ 0 и ϕ(x) ≡ 0, мы можем применить к уравнению (2.3.67)теорему 1.3.2. В результате получим следующее утверждениеЛемма 2.3.8 С вероятность 1 цепь Маркова {(kj , xj )}∞j=0 с переходной плотностью p (k, x) → (m, y) либо обрывается за конечное число шагов в некоторой точке (m, y), либо последовательность {kj }∞j=0 становится стационарной(kj = m,j > j0m 6= 0,m 6= N + 1), а расстояние ρ(xj , ∂D) → 0.Пусть δ = min{j : kj ∈ {0, N + 1}} — момент обрыва цепи, стартующейиз точки (k, x), а χj — индикатор события {δ > j}.
Используя оценку типа ξk (x)на каждом шаге марковской цепи за j шагов получим случайную величинуηj = χ j ue(kj , xj ) + (1 − χj ) χ{kδ =0} ϕ(xkδ ) + χ{kδ =N +1} fkδ−1 (xkδ )/a0 (xkδ , tkδ−1 ) ,которая является несмещенной оценкой для ue(k, x). Очевидно, что последовательность {ηj }∞j=0 образует ограниченный мартингал из которого стандартнымспособом строится малосмещенная оценка для ue(k, x), учитывающая граничныеусловия и не содержащая неизвестной функции.Естественно назвать построенную статистическую оценку — оценкой попоглощению.
Соответствующую выбранной цепи оценку по столкновениям также нетрудно записать. Конечность ее дисперсии вытекает из общих теорем илемм, доказанных в разделе 1.Ограниченный мартингал сходится с вероятностью 1 к некоторой случайной величине η∞ , которая несмещенно оценивает ue(k, x).Можно показать, что при нулевых граничных и начальных условиях этавеличина либо равна нулю, либо совпадает с отношением fkδ−1 (xkδ )/a0 (xkδ , tkδ−1 ).128Заменяя в оценке η∞ , функцию f на погрешность аппроксимации ψ, получаем несмещенную оценку погрешности zk (x).
Следовательно, порядок точности разностно-дифференциальной схемы совпадает с порядком погрешностиаппроксимации.2.4Одна нелинейная краевая задачаВ данном параграфе рассматривается начально-краевая задача для урав-нения теплопроводности с нелинейным условием Стефана-Больцмана, связывающим тепловой поток на поверхности абсолютно черного тела с его температурой. Задача рассматривается в ограниченной выпуклой области D ⊂ R3 ,граница которой S — достаточно гладкая поверхность с внешней нормалью ν,например, поверхность Ляпунова.Температура тела u(x, t) как функция точки x и времени t является решением уравнения теплопроводности∂u= a∆u + f (x, t),∂tx ∈ D,t > 0,(2.4.1)удовлетворяет начальному условиюu(x, 0) = u0 (x),x∈D(2.4.2)и граничному условию Стефана-Больцмана∂u= β(x, t) − κu4 (x, t),∂νx ∈ S,t > 0,(2.4.3)где a и κ являются постоянными.Решение (2.4.1–2.4.3) существует и единственно на достаточно малом интервале времени [59].