CC BY
32
Том 154, kii. 4
УЧЕНЫЕ ЗАПИСКИ КАЗАНСКОГО УНИВЕРСИТЕТА
Физико-математические пауки
2012
УДК 519.63
ДВУХСЕТОЧНЫЙ МЕТОД ДЛЯ ЭЛЛИПТИЧЕСКОГО УРАВНЕНИЯ С ПОГРАНИЧНЫМИ СЛОЯМИ НА СЕТКЕ ШИШКИНА
C.B. Тиховская
Аннотация
Рассмотрено лилейное эллиптическое уравнение с регулярным и параболическими пограничными слоями. Для его решения использована схема направленных разностей па сетке Шишкина, обладающая свойством сходимости, равномерной по малому параметру. Исследован вопрос уменьшения необходимого количества арифметических действий для реализации разностной схемы па основе двухсеточпого метода. Исследован метод Ричардсона для повышения точности разностной схемы. Приведены результаты численных экспериментов.
Ключевые слова: эллиптическое уравнение, сингулярное возмущение, сетка Шишкина, разностная схема, итерационный метод, двухсеточпый метод, экстраполяция Ричардсона.
Введение
Известно, что в случае линейной эллиптической задачи с пограничными слоями равномерная сходимость разностных схем может быть обеспечена сгущением сетки в пограничных слоях [1]. Разностная схема является пятиточечной и представляет собой систему линейных уравнений, которая обычно решается на основе итераций.
Количество итераций можно уменьшить, применяя двухсеточпый метод, когда краевая задача предварительно решается на более редкой вспомогательной сетке с последующей интерполяцией найденного решения на исходную сетку. Найденное на основе интерполяции сеточное решение далее принимается за начальное приближение для итераций на исходной сетке, что и приводит к уменьшению числа арифметических действий.
В работе исследован двухсеточпый метод решения эллиптического уравнения с пограничными слоями на сетке Шишкина [2. 3]. При использовании двухсеточпого метода решение разностной схемы известно на двух сетках, поэтому можно применить экстраполяцию Ричардсона для повышения точности разностной схемы.
1. Двухсеточный метод
Рассмотрим краевую задачу
епхх + епуу + а(х)пх - с(х,у)и = /(х,у), (х,у) € О, и(х,у) = д(х,у), (х, у) € Г,
где О = (0,1) х (0,1), Г = О \ О, а, с, /, д достаточно гладкие функции,
а(х) ^ а, с(х,у) ^ 0, е > 0. (2)
Известно [2]. что при выполнении условий (2) решение задачи (1) является равномерно ограниченным и имеет регулярный пограничный слой у границы х = 0 и параболические пограничные слои у границ у = 0, у = 1. Зададим в области О кусочно-равномерную сетку [2]:
^N,aN = {(Хг,у,-), 1,3 = 0, Н = Хг - Хг_1, Т} = У} - у—},
Н 1 ^ • ^ N Н 2(1 - дх) N ^ А Нг = -, 1 ^ г ^ —; Н = -, — < г < N.
г N ' ^ ^ 2 ' г N '2 '
4ау N _ 2(1 - 2ау) N 3N
т} = —-, 1 < 3 < —, - <3 < N; т} = —-^, — < 3 < -, (3)
} N ' 7 4 4 ^ 7 ^ ' } N 4 4
ах =ш1^|11, — 1п , ау =тш|1, 1п .
На заданной сетке QN,o■N выпишем схему с направленными разностями: 0 ( N N N N \
- иг} иг} - иг-1А +
hi + hi+1 \ hi+i hi
0 f N _ N N _ N
2£ I Ui,j+1 Ui,j Ui,j Ui,j-1
Tj + Tj+i V Tj+i
u.
N - u N
i+1,j ui,j , \ N
a(xi)-h-:--c(xi,yj)u4,j = f(xi,yj), (xi,yj) G , (4)
hi+1
i,j = g(xi,yj), (xi,yj) G Fn.^n = Г n Qn,vn. Определим норму непрерывной функции я и норму сеточной функции zN со-
ответственно:
||z|| = max_|z(x,y)|, \\zN||n = max |zN
(x,y)en i,j
Пусть [u]nNj^N - проекция функции непрерывного аргумента u(x, y) на сетку
Отметим, что матрица схемы (4) является M-матрицей. В соответствии с [2, 4]
||uN - [u]fiN,,N ||n < C1An, An = lnN/N, (5)
где под Ci понимаем положительные постоянные, не зависящие от е и N. Разностную схему (4) запишем в общем виде
LNuN = fN. (6)
Схему (4) можно разрешить на основе итераций
u(m+1) = G( u(m)), (7)
где u(0) задано. Матрица системы (6) является M-матрицей, что обеспечивает сходимость многих итерационных методов [5, 6]. Пусть для решения системы (6) применяется сходящийся итерационный метод и справедлива оценка сходимости итераций
||u(m+1) - u N||n < 4N ||u(m) - u N ||n, 4N < 1. (8)
Чтобы итерационная погрешность не доминировала над погрешностью разностной схемы, итерации необходимо продолжать, пока не выполнится условие
||и(т) - им< Дм. (9)
Пусть
5м = ||и(0) - им||м,
тогда, учитывая (8). несложно заключить, что для выполнения (9) потребуется
Им « ¿К2 ^(Дм/6м) (10)
арифметических действий в предположении, что для реализации одной итерации необходимо выполнить ¿К2 арифметических действий (пропорциональное числу неизвестных).
Теперь рассмотрим двухсеточный метод для нахождения решения разностной схемы (4). Введём сетку Оп,аы такую же то структуре, как только с на-
много меньшим количеством узлов п : п ^ N. Предварительно с использованием итерационного метода (7) решаем задачу (1) НЕ С6ТК6 ' Итвр<ЩИИ НЕ С6ТК6
Qn,aN продолжаем, пока не выполнится условие
||и(тп) - ип||п < Дп, (И)
где тп - необходимое количество итераций.
Далее необходимо найденное сеточное решение и(тп) интерполировать на узлы исходной сетки . При этом точность интерполяционной формулы должна
быть не ниже точности используемой разностной схемы.
Пусть I(уп, х, у) - пнтерполянт сеточной функции «п в областп О и справедлива оценка погрешности
||1 ([и]п„><7„,х,у) - и|| < С2^,п, (12)
где и(х, у) - решение задачи (1), причем Дтг,п ^ Дп-
Пнтерполянт I(уп,х,у) должен быть устойчивым к возмущению интерполируемой сеточной функции уп:
||1 (*п,х,у) - I(йп,х,у)|| < Сэ||*п - 5п|п. (13)
Определившись с формулой интерполяции, можно найденное сеточное решение и(тп) интерполировать та узлы исходной сетки Пусть
ип = [1 (и(тп),х,у)]^.
Покажем, что для некоторой постоянной С4 имеет место оценка
||ип - им||м < С4Дп. С учетом неравенств (5), (11) (13) имеем
||< - им< ||1 (и(тп), х, у) - I(ип, х, у)|| + ||1 (ип, х, у) - I([и]^, х, у)|| +
+ У^([иК,^,х - yN + У[иЬ^ - иМУм <
^ СзДп + СзС\Дп + С2Д^п(,п + С\Дм ^ С4Дп.
Итак, с помощью решения задачи (1) на вспомогательной сетке Пп,^ получено приближение к решению схемы (4) с точностью 0(Дп). Используя это приближение как начальную итерацию в методе (7), уменьшаем количество итераций на исходной сетке ■ Вычислим необходимое количество арифметических дей-
ствий двухсеточного метода:
Ипн « dn2 ^ Дп + ¿М2 ^ Д^ + (14)
^п Дп
где Jn — количество арифметических действий, необходимых для интерполяции.
Сравнивая (10), (14), оцепим выигрыш по количеству арифметических операций при использовании двухсеточного метода:
мк - ИпК > ¿(М2 - п2) ^ (Д") - Зп-
2. Экстраполяция Ричардсона
Исследуем экстраполяцию Ричардсона для повышения точности разностной схемы при использовании двухсеточного метода. Метод экстраполяции Ричардсона для сингулярно возмущенных эллиптических задач исследовался, например, в работах [3, 7]. Рассмотрим случай, когда вспомогательная сетка Шишкина имеет вид (3) и содержит в к раз меньше по каждому направлению сеточных интервалов, чем исходная сетка, где к > 1 - целое число.
Решение схемы с направленными разностями (4) вычисляется па сетке Q.N,aN ■ Для повышения точности будем также использовать решение этой схемы на сетке Пп,^, которая содержит п сеточных интервалов и имеет те же параметры ах, ау, что и сетка П^^ • Эти сетки вложены так, что Qn,o■N — {(Х;,Ут)} С С ^N,ON — {(xi,Уj)}• Очевидно, что сетку П^^ можно получить из Qn,o■N деле-
к
Обозначим решение схемы (4) НЕ Пп,а_^ 3 Пп , ПуСТЬ
п _ 1 _ М к
кп — Гт — - 7; kN —
М - п к - 1' М - п к - 1
В соответствии с методом экстраполяции Ричардсона зададим функцию пп-^ на сетке П^^, приближающую решение п(х, у) с более высоким порядком точности, чем п^ ^ ^^^ ^^^^^ ^^^^^^ ^ ^^ах вспомогательной сетки Пn,ON определим сеточную функцию nnN следующим образом:
^ (Х;,Ут) = кпип(Х1 ,Ут) + kN пN (Х1,Ут), (Х; ,Ут) е П
п,0N •
В узлах исходной сетки ПN,ON, те совпадающих с узлам и сетки Пn,ON, зададим сеточную функцию пп-^(xi,yj), используя интерполяцию. Тогда для каждого узла (xi,yj) е ПN,ON , принадлежащего некоторой ячейке — [Х;_1,Х;+1]х х[Ут-Ь Ут+1], определим
п
^(х^-)— I([nnN]п„^ ,Xi,Уj), (15)
где интерполяционную формулу для ячейки £;,т зададим следующим образом [8]:
I ([п^ ]о ^^ ) — ^ ^^ - ^^ (yj - Ут) +
+
Ут+1 Ут
\хн , Ут+1) - 2пфУ (^, Ут ) + пф^^ (^, Ут_ 1)
©т+1 - 2©т + 0т_1
(©(Уз ) - е(Ут)-
в(Ут) !Ут-1) (Уз - У»)) + (з^Ут), (16)
У т+1 У т '
ф
где
пЮ V „, \ „,гМ
Г ) = ^) + """(Х"уХ -и\(Х1-1,у])(х —№
Лг+1 — XI
-_____-^ — Ф(Хг) —
и
Фг+1 — 2Ф1 + Фг_1
Ф(Хг) — Ф(Х,_1)
Хг+1 — Х\
(Хг — Х,)) .
Интерполяция (16) является точной на функциях Ф(х) и ©(у) и будет квадра-тической при задании Ф(х) = х2, ©(у) = у2. Функции Ф(х) и ©(у) можно выбирать в соответствии с пограничными слоями задачи (1).
3. Результаты численных экспериментов
Рассмотрим краевую задачу
еихх + еиуу + и,х — и = /(х, у), и(х,у) = д(х,у), (17)
где /, д соответствуют решению
и(х, у) = в-х/е + + е-(1-у)/^ + еов(пх).
Решение задачи (17) находим на основе схемы (4). Начальное приближение для используемых итерационных методов задаем следующим образом:
(0)( ^ |0, (хг ) е \ т
)(хг, у2 ) =
{д(хг,уз ), (хг ,уз ) е
оы ■
Итерационный метод (7) на сетке ^N,ON завершаем, если выполнено условие
тах 1Ь?„ и(ты< аАм, г,]
тогда будет выполнена оценка ||и(ты) — им||N ^ А^.
Использование формул линейной и квадратнческой интерполяций сеточной функции с тетки Ип,оы на исходную сетку приводило к тому, что
погрешность метода Ричардсона становилась значительно выше, чем в узлах вспомогательной сетки Пп,оы , поэтому функции Ф, © были заданы следующим образом:
Ф(х)= х2, ©(у) = + е-(1-у)/^.
Исследуем реализацию схемы (4) на основе явного метода Зейделя [6]. Пятиточечную схему (6) можно представить как
аг,] ^-1] + Ьг,] и4]-1 + Сг,] и|+1,2 + ¿г,]и^'+1 — ег,]^Ч,] = /, 0 < "т ] < N. Тогда векторно-матричная запись метода Зейделя будет иметь вид
и(т) = В-1 / + Ьи(т) + ии(т-1)) ,
(Ьу)г1 = а] Уг-1] + Ьг,] у^-^ (Ву)-; = , = Сц «¿+1,] + «¿¿+1.
Табл. 1
Количество итераций односеточного и двухсеточного методов Зейделя, е = 10-4
n N
8 16 32 64 128 256 512
4 39(11) 118(11) 335(10) 970(10) 2936(10) 9331(9) 30847(9)
8 110(37) 320(33) 940(30) 2869(28) 9164(27) 30393(26)
16 292(104) 882(91) 2738(82) 8845(75) 29540(69)
32 789(291) 2525(252) 8305(223) 28087(202)
64 2218(845) 7502(731) 25883(646)
128 6537(2592) 23009(2255)
256 21191(8357)
43 123 347 1002 3045 9694 32068
Табл. 2
Норма погрешности односеточного метода Зейделя и двухсеточного метода Зейделя с экстраполяцией Ричардсона, е = 10-4
n N
8 16 32 64 128 256 512
4 2.67e-1 2.40e-1 2.57e-1 3.15e-1 3.66e- 1 4.09e- 1 4.47e-1
8 7.43e-2 9.43e-2 1.28e- 1 1.60e- 1 1.94e- 1 2.26e- 1
16 3.03e-2 3.74e- 2 4.80e- 2 6.05e-2 7.48e- 2
32 9.82e- 3 1.15e- 2 1.38e- 2 1.73e- 2
64 3.58e-3 2.84e- 3 3.28e- 3
128 1.15e- 3 8.83e-4
256 3.74e- 4
2.30e-1 1.37e-1 7.23e- 2 3.70e-2 2.21e- 2 1.39e- 2 8.20e- 3
Результаты численных экспериментов для метода Зейделя при е = 10-4 приведены в табл. 1 и 2.
В табл. 1 при различных значениях N и n указано количество итераций двухсеточного метода на исходной сетке , при этом в скобках приведено количество итераций на более редкой вспомогательной сетке . В нижней строке
N
Nn
сеточного метода с экстраполяцией Ричардсона (15). В нижней строке для срав-
N
Из таблицы следует, что экстраполяция Ричардсона (15) повышает точность схемы до порядка O(ln2 N/N2).
Помимо метода Зейделя в качестве итерационного метода использовался метод последовательной верхней релаксации [6]
u(m) = wD-1 f + Lu(m) + Uu(m-1)) + (1 - w)u(m-1),
где ш - итерационный параметр, а также метод Писмапа- Речфорда (продольно-поперечных прогонок) [6]
u(m-1/2) = u(m-1) + т (Aiu(m-1/2) + A2U(m-1) - f) ,
u(m) = M(m-1/2) + т (A1U(m-1/2) + A2U(m) - f) ,
т
Табл. 3
Сравнение количества итераций, е = 10-4
Метод N
8 16 32 64 128 256
Зейделя 39(11) 43 110(37) 123 292(104) 347 789(291) 1002 2218(845) 3045 6537(2592) 9694
верхней релаксации 32(7) 35 92(30) 103 245(87) 291 663(244) 841 1862(709) 2553 5471(2174) 8123
Писмана Речфорда 11(2) 12 26(13) 30 57(31) 68 122(70) 152 258(155) 333 559(340) 744
(A1v)i,j = ai,j vi-1,j - e'i,j vi,j + ci,j vi+1,j > (A2v)i,j = bi,jvi,j-1 - ei,jvi,j + di,jvi,j+1■
Здесь ei}j = eij + ei'j , ei j соответствует аппроксимации производных по направлению х, a ei ^j - то направлению y.
В табл. 3 для е = 10-4 проведено сравнение исследуемых итерационных методов. Двухсеточный метод используем в случае n = N/2. Для метода последовательной верхней релаксации итерационный параметр ш = 2/(1 + %/0Т7), для метода Писмана-Речфорда итерационный параметр т = 8/N.
Применение методов последовательной верхней релаксации и Писмана Реч-форда привело к уменьшению количества итераций по сравнению с методом Зей-деля. Повышение точности на основе экстраполяции Ричардсона не зависит от конкретного итерационного метода.
Итак, применение двухсеточного метода приводит к выигрышу в количестве арифметических действий, а использование экстраполяции Ричардсона повышает точность разностной схемы на порядок.
Автор благодарит профессора А.И. Задорина за полезные обсуждения при подготовке настоящей работы.
Работа выполнена при финансовой поддержке РФФИ (проект Х- 11-01-00875).
Summary
S. V. Tikhovskaya. A Two-Grid Method for an Elliptic Equation with Boundary Layers on a Sliishkin Mesh.
In this article, we consider a linear elliptic equation with regular and parabolic boundary layers. To solve this equation, we used an upwind difference scheme on a Sliishkin mesh, possessing the property of the uniform convergence with respect to a small parameter. We investigated the problem of decreasing the required number of arithmetic operations for implementation of the difference scheme on the basis of a two-grid method, and studied the Richardson method aimed at improving the accuracy of this scheme. Here, we present results of some numerical experiments.
Key words: elliptic equation, singular perturbation, Sliishkin mesh, difference scheme, iterative method, two-grid method, Richardson extrapolation.
Литература
1. Miller J.J.H., O'Riordan E., Shishkin G.I. Fitted Numerical Methods for Singular Perturbation Problems. Singapore: World Scientific, 1996. 163 p.
2. Шишкин Г.И. Сеточные аппроксимации сингулярно возмущенных эллиптических и параболических уравнений. Екатеринбург: Изд-во УрО РАН. 1992. 232 с.
3. Shishkin G.I., Shishkina L.P. Difference Methods for Singular Perturbation Problems. Boca Raton: Cliapmanfe Hall/CRC, 2009. 408 p.
4. Rous E.G., Stynes M., Tubiska L. Robust Numerical Methods for Singularly Perturbed Differential Equations. Berlin: Springer-Verlag, 2008. 604 p.
5. Воеводин B.B., Кузнецов Ю.А. Матрицы и вычисления. М.: Наука, 1984. 320 с.
6. Ильин В. П. Методы конечных разностей и конечных объемов для эллиптических уравнений. Новосибирск: Изд-во ИВМиМГ, 2001. 318 с.
7. Shishkin G.I., Shishkina L.P. A Higher-Order Richardson Method for a Quasilinear Singularly Perturbed Elliptic Reaction-Diffusion Equation // Differ. Equat. 2005. V. 41, No 7 P. 1030 1039.
8. Задорин А.И., Задорин H.A. Интерполяция функций с пограпелойпыми составляющими и ее применение в двухсеточпом методе // Сиб. электр. матем. изв. 2011. Т. 8. С. 247 267.
Поступила в редакцию 01.10.12
Тиховская Светлана Валерьевна аспирант Омского филиала Института математики им. С.Л. Соболева СО РАН. Е-шаП: аЛИюьакауа вуапйех. ги
