Вычислительные технологии
Том 5, № 2, 2000
ТРЕХТОЧЕЧНАЯ РАЗНОСТНАЯ СХЕМА НА ПОЛУБЕСКОНЕЧНОМ ИНТЕРВАЛЕ
А.И.Задорин
Омский филиал института математики им. С. Л. Соболева СО РАН Омск, Россия e-mail: [email protected]
A three-point difference scheme with the infinite number of mesh points is considered.
A method of reducing this scheme to the scheme with a finite number of mesh points is proposed. The transformation error is estimated. The case of degenerate difference scheme is considered. The asymptotic method for reduction of that scheme is investigated.
При численном решении краевых задач на полубесконечном интервале может возникнуть разностная схема с предельным условием на бесконечности. Для нахождения решения такой схемы необходимо ее редуцировать к конечному числу узлов. В данной работе это делается на основе выделения многообразия решений схемы, удовлетворяющих предельному условию на бесконечности. Такое многообразие будет задаваться в виде двухточечной разностной схемы, что позволит свести исходную схему к схеме с конечным числом узлов. Аналогичный подход в случае дифференциальной задачи принимался, например, в [1, 2, 4]. В случае, когда при стремлении некоторого параметра к нулю разностная схема вырождается, используется асимптотический подход для решения вспомогательной начальной задачи.
Определим норму ограниченной сеточной функции:
||zh 11 = max |zQ|.
n
Итак, рассмотрим трехточечную разностную схему:
LhnUh = AnUhn-1 - CnUhn + BnUhn+1 = Fn, n> 0, (1)
uq = G, uh — 0 при n —— to. (2)
Предполагаем, что при всех n
An > 0, Bn > 0, Ch > An + Bn + A, A > 0, (3)
An — A0, Bn — B0, Cn — C0, Fn — 0 при n — to. (4)
Дополнительные ограничения на коэффициенты разностной схемы будем делать отдельно.
Перейдем к анализу свойств схемы (1)-(2). Для оператора схемы (1)-(2) справедлив принцип максимума [5, с. 40], в соответствии с которым в случае бесконечного числа узлов из условий:
> 0, < 0, n > 0, lim ФП > 0 (5)
n—
© А. И. Задорин, 2000.
следует, что при всех п > 0 ФП — 0-
На основании принципа максимума нетрудно показать, что справедлива оценка
|И|< Л-1||^|| + |С|-
Аналогичная оценка в случае конечного числа узлов получена в [5, с. 42].
Выделим многообразие решений схемы (1), удовлетворяющих предельному условию на бесконечности. Для этого определим соотношение
ип апип— 1 + п > (6)
где коэффициенты ап и вп задаем как решения двухточечных разностных схем с предельным условием на бесконечности:
А 9 4°
ап = —-----ат--------, ап ^ а°, п ^ то, а° =------------. , (7)
Сп - Вп«п+^ С° + VѰѰ - 4А°В°
вп = Впвп+в - Рп • вп ^ о, п ^то. (8)
Сп Впап+1
Учитывая (7), (8), можно показать, что если сеточная функция ин удовлетворяет соотношениям (6), то она удовлетворяет (1).
Покажем, что при всех п > 0 ап < 1. Учитывая условие диагонального преобладания
(3), получим:
° < 2А° < А° <
а < С° + |А° - В°| + А < А° + А < 1-
Для достаточно больших п (п > М) в силу предельного условия ап < 1. В случае п < М
ап < 1 в силу условий (3) и соотношения (7). Итак, при всех п > 0 0 < ап < 1.
Учитывая (7), можно показать, что при всех п > 0
Ап
0 < ап < ! I” ^ атах < 1.
Ап + А
Покажем, что соотношение (6) задает многообразие решений схемы (1), удовлетворяющих предельному условию на бесконечности. Из (6) следует, что при всех п
п
к! < атах |и0 1 + Х ^-Х|ві1'
і=1
Учитывая, что в этом неравенстве 0 < атах < 1, ві ^ 0, і ^ то, можно показать, что «П ^ 0, п ^ то.
Зададим некоторое N > 0 и, учитывая соотношение (6), перейдем от (1)-(2) к разностной схеме с конечным числом узлов:
ЬПи" = Ап«П-1 - Сп«П + Вп«П+1 = Рп, П =1, 2, ..., N - 1,
ио = ^, им — ам «N-1 = вм. (9)
Коэффициенты ам и вм из задач (7), (8) могут быть найдены приближенно. Оценим влияние погрешности при их вычислении на решение схемы (9). Анализ влияния возмущения коэффициентов схемы на решение проводился, например, в [3].
Итак, перейдем от (9) к разностной схеме с возмущенными и /м :
АпиП—1 — СпиП + ВпиП+1 = п =1, 9, . . . , N — 1
ио = ^, им — им—1 = вМ. (10)
Лемма 1. Пусть
|ам — | < Аъ 0 < < 1, |вм — /Зм| < А2.
Тогда В случае
выполнится
\un — Un\ < (1 — aN) {А1 \UN —1 \ Т А2}, О < n < N (11)
An
q = min —h > 1 (12)
n<N Bn
\uh - Uh\ < —^—{АіК-і\ Т А2^п^, О < n < N (13)
q — aN
Доказательство. Определим = «О — «о. Тогда является решением задачи
= 0, п =1, 2, ..., N — 1, =0, — ^°-1 = (ам — )«м-1 + (/м — вм). (14)
На основании [5, с. 44] можно убедиться, что для оператора задачи (14) справедлив принцип максимума: из условий
ФО > 0, ^пФ° < 0, п =1, ...,N — 1, £°Ф° > 0 (15)
следует ФП У О при всех n = 0,1, ... , N nh
о
Определим сеточную функцию Фо:
ФП = (1 — ) 1{А1 |им—1| + А2}±
Тогда выполнены условия (15) и в силу принципа максимума при всех п ФП > 0, что доказывает (11).
Рассмотрим случай условий (12). Определим сеточную функцию фо: фП = 9П—М. Тогда при всех п
Lh^h < 0hq 1{Bnq2 — (An т Bn)q т An} = Bn^hq 1(q—1Hq — b~ ^ <0. (1б)
Определим сеточную функцию Ф
h
ФП = {Ai|uN_,| + А2}ФП ± zj,
q —
Учитывая (16), получим, что для функции Ф^ выполнены условия (15). В силу принципа максимума при всех n ФП > 0. Это доказывает оценку (13). Лемма доказана.
Лемма 2. Пусть при всех n > M An > Bn, Тогда справедлива оценка
max |en| < A-1 max |Fn|.
n>M n>M
Доказательство. Перепишем уравнение (8) в виде
ВПв = Вп[вп+1 — вп] — {Сп — Вп — ВПаП+1}вП = ^
Нетрудно убедиться, что если для какой-либо сеточной функции Фо выполнены условия
ДПФ° < 0, п > М, Иш ФП > 0, (17)
то при всех n > M ФП > 0.
Определим сеточную функцию Ф^ :
ФП = A-1 max |Fn| ± А,
n>M
Для такой функции Ф^ выполнятся условия (17). Следовательно, при всех n > M ФП > 0. Это доказывает лемму.
Исследуем, как влияет возмущение коэффициентов схемы (1)-(2) на а и в, которые соответствуют задачам (7) и (8).
Лемма 3. Пусть коэффициенты схемы (1)-(2) возмущены: при всех n > M
| An — An | < O, | Bn — Bn | < O, | Cn — Cn | < O,
An у Bn, An у Bn > о, Cn у An + Bn + А, /А > о An — A0, Bn — B0, Cn — C0 при n — то.
Пусть
Aan
aan
Cn Bnan+1
an — a0, n — то, a
0
2A0
Аг =
-E^n/^n+l Fn
C0 + VC0C0 - 4A0B0:
/5n — 0, n — то.
СП ВПаП+1
Тогда для некоторой постоянной С, не зависящей от а, при всех п > М
|ап — (5п| < |вп — Аг| < Са.
(18)
(19)
Доказательство. Начнем с первой оценки в (19). В силу условий (18), 0 < аП < 1. Пусть = а — а. Тогда
ВП = апВп(^+1 — гП) — (Сп — Впап+1 — апВп^П = С*
где
Gh = An — An + an(Cn — (Cn) + anan+1 (Bn — Bn).
Можно показать, что
\a0 - a0\ <
2a
(C0 + А)(С0 + А)
A0 + 2C0 +
C0 + C0 + 4A0 + 4BB0 A0 - B0 + А + A0 — B0 + А
Определим
Ф« = Са ± г«.
Тогда для некоторой постоянной — выполнятся условия (17) и поэтому при всех п > М Ф« > 0. Это доказывает первую оценку в (19).
Далее получим вторую оценку в (19). Пусть гн = в — /3- Тогда
= ВВга(гга+1 — ^га) — (—« — Вгаага+1 — ^т = 0
где
^га = (Вп — Вга)вга+1 + Рга — ^га + (—п — —га)в«. + (Вп — Вга)вгаага+1 + Вгавга(ага+1 — ап+1)-
В соответствии с условиями леммы,
—п — В„ага+1 — В* > Д-
Теперь вторую оценку в (19) можно получить на основании принципа максимума. Лемма доказана.
Перейдем к вопросу приближенного нахождения и вм из задач (7) и (8). Восполь-
зуемся леммой 3, согласно которой малым изменениям коэффициентов схемы (1) - (2) при п > М соответствуют малые изменения ап и вп при п > М.
Исследуем два подхода для приближенного нахождения ам и вм -
При первом подходе исходим из предельных условий (4). В случае дифференциальной задачи такой подход применялся, например, в [2]. Пусть при достаточно больших п (п > М) для некоторого г справедливы разложения:
г 4(0 /1 \ Г+1 г С(0 /1 \г+1
А* = Х -г + О(-) , С« = ^ -г + О(-) , -(0) = А0, С(0) = С0,
' пг \п п \п
г=0 4 7 г=0 4 7
Г / 1 \ Г+1 Г Р(0 / 1 \ Г+1
в„ = £В- + о(-) , ^ + о(-) , в(0) = в». р(0) = 0' пг \п п \п
г=0 4 7 г=0 4 7
Тогда при п > М ап и вп из(7)-(8) можно искать в виде:
p а(0 p
' « = £ ^ • Р < r. (20)
ap
П , y „ y .
n‘ z—' n
i=0 i=0
Подставляя эти разложения в (7), (8), получим рекуррентную формулу относительно коэффициентов а*г) и вПг)-
Например, в случае р = 1 имеем:
0 -(1) — а0(С(1) — В(1)а0) ~ Р(1)
ам = а + 77^0) п 0 п0\ лт-------------------, = —"
(C0 - 2a0B0)N ’ tN (C0 - а0В0 - B0)N'
Пусть при n > N An > Bn. В соответствии с леммой 3 для достаточно больших N, обеспечивающих выполнение условий (18), для некоторой постоянной C
|«N — <5n |, |вм — @N |< CN-2.
С учетом леммы 1 найдется постоянная — такая, что при всех п < N
|и« — и«1 < —N-2.
Точность данного подхода повышается с увеличением N.
Остановимся на втором подходе. При сеточной аппроксимации уравнений второго порядка с малым параметром при старшей производной возникает вырождающаяся разностная схема. Например, рассмотрим краевую задачу:
ем" — а(х)и/ — 6(х)м = f (х), м(0) = С, м(х) ^ 0, х ^то (21)
в предположении
а(х) > ат;п > 0, е € (0,1], Ь(х) > 6т;п > 0, а(х) ^ а0, Ь(х) ^ Ь0, f (х) ^ 0, х ^ то.
Применяя к этой дифференциальной задаче схему направленных разностей, получим схему (1) - (2) при
—п = ^12 + ^, —п = |е + + 6(х„), Вп = ^, Рп = f (х„). (22)
Из (22) следует, что Вп ^ 0 при е ^ 0. Аналогичным образом может быть рассмотрена другая монотонная разностная схема. Таким образом, разностная схема (1)-(2) в этом случае вырождается в двухточечную при стремлении параметра к нулю. Это обстоятельство можно использовать при переходе к схеме на конечном интервале.
Итак, пусть коэффициенты схемы (1) - (2) зависят от малого положительного параметра е и для некоторого г > 0 справедливы разложения
Г Г
-п = X-V + 0(ег+1), Сп = X-IV + 0(ег+1),
г=0 г=0
г г
Вп = XВ«е* + 0(ег+1), Рп = XР!г)ег + 0(ег+1). (23)
г=1 г=0
Рассмотрим вопрос приближенного нахождения ап из задачи (7). Пусть
р
< = Ха(г)е'’ р <г. (24)
г=0
Учитывая разложения (23) и (24) в (7), собирая члены при одинаковых степенях параметра е, для г > 1 получим рекуррентную формулу:
а(г) = 1 п (0)
С
Сп
п
г—1 г г—]
ап(
.7=0 .7=1 к=0
4(0) —п
(
С
Итак, ап можно искать на основе асимптотической формулы (24), учитывая определенное число членов этого разложения.
Теперь получим асимптотические формулы для вп. Пусть
р
Др = £ в(°е\ р < г. (26)
г=0
Учитывая разложения (23), (24) и (26) в (8), собирая члены при одинаковых степенях параметра е, для г > 1 получим рекуррентную формулу:
А‘> 1
C(0)
с(0)
i—1 i i—j i
-f<->-Ej)+EEB'j)e'k)a'!;-;—k)+Ej) , вп0)=-j=0 j=1 fc=0 j=1
(27)
Итак, значение вп может быть найдено на основе асимптотического разложения (26) и рекуррентной формулы (27) для членов этого разложения. Под C ниже будем понимать положительные постоянные, не зависящие от параметра е.
Лемма 4. Пусть Ап > Вп при всех n > M, а параметр е достаточно мал. Тогда найдется постоянная C такая, что при всех n > M
|an - <|< Cep+1, |/t - вП|< ep+1.
Доказательство. При построении и Дп мы ограничились p членами разложения коэффициентов (23), поэтому ap и (Зр являются решением задач (7), (8) в случае возмущенных коэффициентов схемы (1) - (2) и в условиях (18) выполнится а < Cep+1. Параметр е должен быть мал настолько, чтобы для схемы с возмущенными коэффициентами выполнились условия (18). Теперь утверждение леммы следует из леммы 3.
Точность расчета коэффициентов ап и вп повышается с увеличением числа членов асимптотического разложения и с уменьшением параметра е. Согласно лемме 1 схема (9) устойчива к погрешностям, возникающим при расчете этих коэффициентов. Решение схемы (9) может быть найдено методом прогонки [5], который в силу диагонального преобладания устойчив.
Используя соотношение (6), исходную схему (1)-(2) можно свести к двухточечной разностной схеме:
Un = 1 + вп, n > 1, Uh = G. (28)
Выше было показано, что предельное условие на бесконечности (2) для решения схемы
(28) выполнится.
Лемма 5. Для решения схемы (28) справедлива оценка
А + A
l|uh|| < |G| + max А п||вII.
п А
Доказательство. Перепишем схему (28) в виде
Rn= ап(иП— иП-1) +(1 — ап)иП = вп-
Нетрудно показать, что если для сеточной функции Ф^ выполнены условия
ф£ > 0, > 0, n > 1, (29)
то при всех n > 0 Фп > 0. Определим Ф^ :
А + A
фп =|G| + max—а-п ||в|1 ± ^
пА
Для такой функции Ф^ выполнены условия (29), и поэтому при всех п > 0 Фп > 0. Это доказывает лемму.
Коэффициенты схемы (28) с учетом малости параметра е могут быть найдены приближенно, согласно формулам (24)-(27). Покажем, что схема (28) устойчива к возмущению этих коэффициентов.
Лемма 6. Пусть йп — решение схемы (28) в случае возмущенных коэффициентов ап и /Зп. Пусть при всех п > 1
|ап «п| — ^, |вп 15п\ — ^, 0 < ап — атах < 1.
Тогда при всех п > 0
!«п - й!1 — (1 + 1!«п!!К
1 атах
Доказательство. Пусть гп = мп — 5п. Тогда гп является решением задачи
^п^Ш = ап(^Ш — ^Ш-1) + (1 — ап)^Ш = (ап — ап)йп-1 + вп — /Зп ^0 = °.
Определим Фп :
фп=, а, (1+!1иш1!)^ ±
1 атах
Для такой функции Фп выполнены условия (29), и поэтому при всех п > 0 Фп > 0. Это доказывает лемму.
Остановимся на результатах численных экспериментов. Рассмотрим схему (1)-(2) с коэффициентами из (22). Из разложений (24) и (26) с учетом формул (25) и (27) можно получить:
~0 ап Я0 /пЬ ~1 ~0 Ьпап+1 апЬп+1 + ЬпЬп+1Ь
а = ---------- в =---------------а = а + е
п \ Т Т') Уп | 1 и ? ^п ^п 1 °
ап + ЬпЬ ап + Ьп'
ап + ^h п йп + 6п^ п п (йп + 6п^)2(ап+1 + бп+1^)
51 50 /пап+1 ап/п+1 + 2/п Ьп+1Ь /п+1ЬпЬ
п п (ап + Ьп')2(ап+1 + Ьп+1Ь)
где ап а(хп), Ьп ^(^п^ /п /(хп).
Если в соотношении (6) принять ап = ап, вп = вп, то (6) будет соответствовать вырожденной разностной схеме (1)-(2). В численных экспериментах в качестве решения схемы (1)-(2) принималось решение схемы на достаточно длинном (Ь0 = 100) интервале, когда решение схемы не зависит от способа задания краевого условия.
Рассматривалось четыре способа задания краевого условия при переходе к схеме с конечным числом узлов:
1) < = 0;
lN
2) uN = uN-1;
3) согласно (10) с aN = aN, /N = /N;
4) согласно (10) с aN = aN, /N = /?Nv.
При проведении экспериментов предполагалось G = 1, h = 0.1, L =1, где L — длина интервала, к которому сводилась исходная разностная схема, h — шаг сетки. Остановимся на примере:
йп = 1 + ехр(-хп), Ьп = 1 + ехр(-хп), /п = ехр(-хп), Хп = nh.
В табл. 1 приведена норма погрешности
s = max |«п -
в зависимости от способа задания краевого условия, где — решение схемы (1)-(2) на достаточно длинном интервале [0,Ьо], — решение схемы на интервале [0,Ь]. Анало-
гичным образом в табл. 1 приведена норма погрешности, возникающая при переходе от (1)-(2) к начальной задаче (28). При этом при расчете погрешности в соответствует решению схемы (1)-(2) на достаточно длинном интервале, соответствует решению задачи (28) на интервале [0,Ь].
Теперь рассмотрим пример:
йт(ж„)
2 +
Хп + 1
— 3 +
0О8(жга)
Хп + 1
/п
1
Хп + 1
В табл. 2 приведена норма погрешности для данного случая. Результаты вычислений подтверждают полученные оценки.
Таблица 1
п
Способ переноса условия £
1 10-1 10-2 10-3
и% = 0 0.38 0.20 0.17 0.17
-1 0.19 0.69Е-1 0.47Е-1 0.44Е-1
в 1 ) р* "вм 0.86Е-1 0.54Е-2 0.40Е-3 0.38Е-4
С0 т—Н в 0.11 0.62Е-3 0.46Е-5 0.60Е-7
в (28) й°, РП 0.21 0.32Е-1 0.34Е-2 0.34Е-3
в 2 ) а , 0.23 0.30Е-2 0.31Е-4 0.30Е-6
Таблица 2
Способ переноса условия £
1 10-1 10-2 10-3
и% = 0 0.26 0.13 0.10 0.10
-1 0.14 0.48Е-1 0.36Е-1 0.34Е-1
в 1 ) р* , "вм 0.42Е-1 0.24Е-2 0.19Е-3 0.19Е-4
0) (1 в 0.42Е-1 0.24Е-3 0.19Е-5 0.60Е-7
в (28) «п,вп 0.18 0.31Е-1 0.34Е-2 0.34Е-3
в 2 ) р* , ^р' 0.28 0.45Е-2 0.50Е-4 0.39Е-6
Список литературы
[1] АБРАМОВ А. А., Балла К., Конюховл Н. Б. Перенос граничных условий из особых точек для систем обыкновенных дифференциальных уравнений. В “Сообщ. по вычисл. матем.” ВЦ АН СССР, М., 1981.
[2] Биргер Е. С., ЛяликовА Н. Б. О нахождении для некоторых систем обыкновенных дифференциальных уравнений решений с заданным условием на бесконечности. Журн. вычисл. матем. и матем. физики. 5, №6, 1965, 979-990.
[3] Годунов С. К., Рябенький В. С. Разностные схемы. Наука, М., 1977.
[4] Задорин А. И. Перенос краевого условия из бесконечности при численном решении уравнений второго порядка с малым параметром. Сиб. журн. вычисл. матем., 2, №1, 1999, 21-35.
[5] Самарский А. А. Теория разностных схем. Наука, М., 1989.
Поступила в редакцию 6 января 1999 г., в переработанном виде 26 марта 1999 г.