Диссертация (1145371), страница 19
Текст из файла (страница 19)
Следовательно, кблужданию по сферам применима теорема 3.2.2. Несмещенная оценка решениязадачи Дирихле на траекториях блуждания по сферам определяется формулойξδ =τXδ −1j=0R2 (xj )f (yj ) + u(xτδ ),2nв которую входят величины, определяемые равенствами:x0 = x,xj+1 = xj + R(xj )ωj ,yj = xj + R(xj )ρj ωj ,j = 0, 1, .
. .Отметим, что величина u(xτδ ) = ϕ(xτδ ), если τδ = τ1 . Если же τδ < τ1 , то u(xτδ )заменяют на ϕ(x∗τδ ), где x∗τδ — ближайшая к xτδ точка границы Γ. Алгоритммоделирования величины ρ представлен в § A.4.1.Пример 3.2.2. Блуждание по сферам и шарамРешается задача Дирихле для уравнения Гельмгольцаx ∈ D ⊂ Rn ,∆u(x) − a(x)u(x) = −f (x),n > 2.Предполагается, что коэффициент a(x) непрерывен и удовлетворяетнеравенству 0 6 a(x) 6 c2 , где c — положительная постоянная.В качестве области T (x) ⊂ D, берется шар радиуса R(x) с центромв точке x. Функция R(x) удовлетворяет неравенствам c1 dist(x, Γ) 6 R(x) 66 dist(x, Γ) и непрерывна.Функция Леви определяется равенством −cree−cR1,−L(y, x) =σn (n − 2) rn−2 Rn−2где r = |x−y|. Асимптотика (3.1.7) для данной функции справедлива при λ = 1.Вычисления показывают, что∆y L(y, x) = e2r + c(n − 3).σn (n − 2)rn−1−cr c146Следовательно,Ny L(y, x) =e−cR a(y)e−cr2(c−a(y))r+c(n−3)+.σn (n − 2)rn−1σn (n − 2)Rn−2Нормальная производная функции L(y, x) на поверхности шара T (x) совпадаетс производной по радиусу r и равна −cR∂LcRe=− 1+.∂νn − 2 σn Rn−1Значит, для выбранной функции Леви выполнены условия (3.1.12), а для решения задачи Дирихле справедливо уравнение (3.1.13).
Функция h(x) удовлетворяет неравенствамZh(x) >2(cR)R2 e−cRe−cRcRe − 1 − cR −>L(y, x)dy = 2>c (n − 2)n2nT (x)(c1 δ)2 e−cd,> c(δ) =2nпри dist(x, Γ) > δ. Здесь d — диаметр области D.В данном примере удобно строить несмещенные оценки с весами (1.4.4)на траекториях необрывающейся внутри области цепи, переходная вероятность для которой совпадает с ядром уравнения при a(y) ≡ 0. Она является смесью изотропного распределения в шаре и равномерного распределенияна окружности. Опишем процедуру моделирования следующей точки цепи потекущей.
Пусть ω — изотропный единичный вектор, γ1 , γ2 — экспоненциально(с параметром λ = 1) распределенные независимые случайные величины, α —распределена равномерно и независима от уже определенных случайных элементов. Пусть χ(n) — индикатор интервала [0, 1/(n − 2)], случайная величинаρe = γ1 + χ(n) (α)γ2 /c, и, наконец, ρ = min(eρ, R(x)), тогда x1 = x + ρω. Припереходе из x в x1 , весовой множитель умножается на величинуa(x1 ) ρ (n−2) −c(R−ρ)p(x, x1 ) = 1 − χ (α) 2e.cR(n)147ФункциюZL(y, x)f (y)dyF (x) =Tможно несмещенно оценить нулем, при ρ1 > cR, и величиной ρ (n−2)f (y)111−,e−c(R−ρ1 )n−2Rc2при ρ1 6 cR.
Где ρ1 = (γ1 + γ2 )/c, а точка y = x + ρ1 ω.3.2.1Блуждание по эллипсоидамСледуя работе [42], опишем алгоритм блуждания по эллипсоидам длярешения задачи Дирихле (3.2.1). В качестве областей T (x) выберем областиограниченные эллипсоидамиT (x) = TR (x) = {y : σ(y, x) = (A−1 (x)(y − x), (y − x))1/2 6 R} ⊆ D.(3.2.2)Справедлива лемма.Лемма 3.2.1 Существует непрерывная функция R = R(x) удовлетворяющая в области D, при некоторых положительных постоянных c1 и c2 , неравенствам c1 dist(x, Γ) 6 R(x) 6 c2 dist(x, Γ) и такая, что выполнено включение TR (x) ⊂ D.Доказательство.
Из условия равномерной эллиптичности оператораследует, что в качестве множества T (x) можно взять максимальный эллипсоидTR (x), лежащий в D. Действительно, пусть y ∈ Γ — точка касания максимального эллипсоида с границей Γ. Тогда, из условия равномерной эллиптичностиоператора M следуют неравенства1111µ− 2 dist(x, y) = µ− 2 |x − y| 6 σ(y, x) = R 6 ν − 2 |x − y| = ν − 2 dist(x, y).Функция R(x) = inf y∈Γ σ(y, x) является непрерывной. Действительно, из условий на коэффициенты оператора следует, что σ(y, x) является равномерно148непрерывной функцией на компакте Γ × D.
По заданному ε > 0 определим δ > 0 так, что при всех y ∈ Γ и |x − x1 | < δ выполнено неравенство|σ(y, x) − σ(y, x1 )| < ε. Тогда, R(x) 6 σ(y, x) 6 σ(y, x1 ) + ε, и, следовательно,R(x) 6 R(x1 ) + ε. Меняя местами x и x1 , получаем обратное неравенство.Для произвольной положительной функции p(ρ) определим функциюL(y, x) равенствомL(y, x) =ZR1pq(R)σn (n − 2) det(A(x))σ 2−n − ρ2−n p(ρ)dρ,(3.2.3)σгдеZRσ = σ(y, x),R = R(x),q(R) =p(ρ)dρ.0Лемма 3.2.2 Пусть p(ρ) — ограниченная (p(ρ) 6 C) и непрерывная функцияна [0, ∞), тогда функция L(y, x), определенная формулой (3.2.3), являетсяфункцией Леви.Доказательство. Проверим справедливость асимптотики (3.1.8).
Положим−1pµ(x) = q(R)σn (n − 2) det(A(x)).Тогда,L(y, x) − H(y, x) = −µ(x) σ 2−nZσZRp(ρ)dρ +0ρ2−n p(ρ)dρ .σЕсли x изменяется в замкнутой подобласти T ⊂ D, то функция µ(x) ограниченасверху некоторой постоянной C1 .При n = 3, получаем неравенство|L(y, x) − H(y, x)| 6 C · C1 [1 + | ln R| + | ln σ|] 6 C · C1 [const +| ln r|] .149Следовательно, при 0 < λ < 1, верно неравенство |L(y, x) − H(y, x)| 66 const ·rλ−1 .При n > 3,"|L(y, x) − H(y, x)| 6 C · C1 σ 3−n +# 3−nσ− R3−n n−36 const ·r3−n .Значит, первое из условий (3.1.8) выполнено при λ = 1.Вычисляя производные, находим∂(L − H)∂σ= µ(x)(n − 2)σ 1−n∂yi∂yiZσp(ρ)dρ,0∂ 2 (L − H)1 ∂ 2H=−∂yi ∂yjq(R) ∂yi ∂yjZσp(ρ)dρ −1 ∂H∂σp(σ).q(R) ∂yi∂yj0Оценивая производные, получаем неравенства: ∂(L − H) ∂H Cσ6 ∂yi q(R) ∂yi , 2 2 ∂ (L − H) ∂ H ∂H ∂σ C ∂yi ∂yj 6 q(R) ∂yi ∂yj σ + ∂yi ∂yj .Если x изменяется в замкнутой подобласти T ⊂ D, то функция q R(x) ограничена снизу положительной постоянной.
Используя асимптотику для производных параметрикса H(y, x), получаем неравенство ∂(L − H) 6 const ·r2−n . ∂yiПоскольку, n µ−1 nr ∂σ 1 Xν 1/2 n==,A(y−x)6jk kk ∂yj σ ν −1/2 rµk=1справедливо неравенство 2 ∂ (L − H) 1−n ∂yi ∂yj 6 const ·rСледовательно, второе и третье из условий (3.1.8) выполнены при λ = 1.150Подставляя выражения для производных1 ∂H∂L=∂yiq(R) ∂yiZRp(ρ)dρ,σ∂ 2L1 ∂ 2H=∂yi ∂yjq(R) ∂yi ∂yjZRp(ρ)dρ −1 ∂H∂σp(σ)q(R) ∂yi∂yjσв формулу (3.1.3) для сопряженного оператора, получаем1Ny L(y, x) =q(R)ZRσnnp(σ) X X∂H ∂σp(ρ)dρ · Ny H −−aij (y)q(R) i=1 j=1∂yi ∂yjZR− µ(x)d(y)ρ2−n p(ρ)dρ, (3.2.4)σгдеnnn XXX∂ 2 aij∂bi+ c(y) −.d(y) =∂y∂y∂yijii=1i=1 j=1Выберем вес p(ρ) так, чтобы функция Ny L(y, x) была неотрицательнойв T (x).
Оценим снизу каждое слагаемое в формуле (3.2.4). Заметим, чтоnn XXaij (y)i=1 j=1∂σ ∂σ1= 2 (A(y)A−1 (x)(y − x), A−1 (x)(y − x)) >∂yi ∂yjσ2 ν 2 −12 ν 2ν −1> 2 A (x)(y − x) > 2 A (x)(y − x) > 2 = κ1 . (3.2.5)σrµПоэтому,nnp(σ) X X∂H ∂σ−=aij (y)q(R) i=1 j=1∂yi ∂yj= µ(x)(n − 2)p(σ)σ1−nn XnXi=1 j=1aij (y)∂σ ∂σ> µ(x)(n − 2)p(σ)σ 1−n κ1 . (3.2.6)∂yi ∂yjЕсли коэффициенты оператора M достаточно гладкие, например, aij ∈ C (2) (D),bi ∈ C (1) (D) и c ∈ C(D), то Ny H 1−n6 κ2 µ(x)ν (1−n)/2 σ 1−n q(R) 6 κ2 µ(x)r(3.2.7)151Пусть κ3 = maxx∈D |d(x)|.
Из формул (3.2.4–3.2.7) следует, что функцию p(ρ),обеспечивающую неотрицательность Ny L на T (x), достаточно подчинить условию(n − 2)κ1 p(σ) > κ2 ν (1−n)/2ZRZRp(ρ)dρ + κ3 σσρ2−np(ρ)dρ,σ 2−nσкоторое заведомо будет выполнено, еслиZRp(σ) > κp(ρ)dρ,σприκ2 ν (1−n)/2 + κ3 R.(n − 2)κ1является искомой функцией.
Таким образом, докаκ>Очевидно, что p(σ) = e−κσзана следующая теорема.Теорема 3.2.3 Пусть коэффициенты оператора M являются гладкимифункциями: aij ∈ C (2) (D), bi ∈ C (1) (D), c ∈ C(D). Тогда существуют непрерывная функция R(x) и функция Леви L(y, x), удовлетворяющие условиям:1) для x ∈ D и некоторых положительных постоянных c1 , c2 справедливынеравенства c1 dist(x, Γ) 6 R(x) 6 c2 dist(x, Γ);2)L(y, x) =для y ∈ ∂T (x),∂L(y, x)=0∂yii = 1, 2 . .
. , n;3) Ny L(y, x) > 0 для y ∈ T (x);4) для решения первой краевой задачи справедливо интегральное представлениеZZNy L(y, x)u(y)dy +u(x) =T (x)L(y, x)f (y)dy,T (x)в котором ядро оператора является субстохастическим, если выполненоусловие c(x) 6 0.152Определим переходную вероятность цепи Маркова x0 = x, x1 , . . . , равенствомZNy L(y, x)dy.P (x, B) =T (x)∩BПрименим к уравнениюZZNy L(y, x)u(y)dy +u(x) =T (x)L(y, x)f (y)dy(3.2.8)T (x)теорему 3.2.2. Функцию h(x) легко вычислить, применяя лемму 2.1.2 к семейству эллипсоидов σ(y, x) = r, 0 6 r 6 R(x):Zh(x) =1L(y, x)dy =q(R)(n − 2)ZRdrrn−11=q(R)(n − 2)ZρZRdρ p(ρ)0(r2−n − ρ2−n )p(ρ)dρ =r0T (x)ZR1(r − rn−1 ρ2−n )dr =2nq(R)0ZRρ2 p(ρ)dρ.0Рассмотрим функцию1h1 (r) =2nq(r)Zrρ2 p(ρ)dρ.0Ee производная положительна2dh1=drRrr p(r)q(r) − p(r) ρ2 p(ρ)dρ02nq 2 (r)> 0.Следовательно, при dist(x, Γ) > δ выполнено неравенство R(x) > c1 δ и неравенство h(x) > h1 (c1 δ) > 0.Таким образом, все условия теоремы 3.2.2 выполнены.
Следовательно,построенная цепь Маркова с вероятностью 1 либо обрывается внутри областив некоторый конечный момент τ1 , либо сходится к точке на границе Γ. Пустьδ > 0, τ2 — момент первого попадания цепи в δ-окрестность границы Γ и τδ == min(τ1 , τ2 ). Математическое ожидание Eτδ < ∞.153Последовательность оценок ξi =Pi−1j=0 h(xj )f (yj )χj+ χi u(xi ), определен-ная в (1.4.2) образует квадратично интегрируемый мартингал. Случайная величина ξτδ является несмещенной оценкой u(x) и имеет конечную дисперсию.Случайные векторы yj распределены в T (xj ) с плотностью L(y, xj )/h(xj ).Алгоритмы моделирования этой плотности и самого блуждания по эллипсоидам содержатся в §A.4.2.Результаты вычислительных экспериментовАлгоритм блуждания по эллипсоидам был применен для решения первойкраевой задачи в кубе [0, 1]3 .Матрица A(x) была выбрана диагональной. Элементами её главной диагонали равны:a11 (x) = x21 + 1,a22 (x) = x22 + 2,a33 (x) = 1.Остальные коэффициенты эллиптического оператора M равны нулю.
Праваячасть уравнения f (x) равна −2x21 − 2x22 − 6. В качестве точного решения и краевых условий выбрана функция u(x) = x21 + x22 + 2x3 . При попадании блужданияв δ - окрестность границы в статистической оценке использовалось точное значение решения, а не его приближенное значение в ближайшей точке границы.Таким образом, использовалась несмещенная статистическая оценка решенияуравнения. Параметр k, обеспечивающий неотрицательность ядра интегрального уравнения был выбран равным 30. Отметим, что этот параметр влияет насреднюю величину перемещения за один шаг и, следовательно, на длину траектории блуждания и время работы алгоритма. Ниже приведены результатывычислений для десяти точек куба, выполненных для параметра δ = 0.001 ивыборки объёмом 10000 .