Апостериорно адаптируемые (по градиенту решения) сетки в аппроксимации сингулярно возмущенных уравнений конвекции-диффузии Текст научной статьи по специальности «Математика»

Вычислительные технологии Том 6, № 1, 2001

АПОСТЕРИОРНО АДАПТИРУЕМЫЕ (ПО ГРАДИЕНТУ РЕШЕНИЯ) СЕТКИ В АППРОКСИМАЦИИ СИНГУЛЯРНО ВОЗМУЩЕННЫХ УРАВНЕНИЙ КОНВЕКЦИИ - ДИФФУЗИИ *

Г. И. Шишкин

Институт математики и механики УрО РАН, Екатеринбург, Россия

e-mail: Grigorii@shishkin.ural.ru

Dirichlet's problem for a parabolic equation of convection-diffusion with a small parameter £ of the highest derivative is considered on a segment. In order to improve the accuracy of the approximate solution consecutive a posteriori refinement of the grid is employed in subdomains determined by the gradients of intermediate discrete problems solutions. The grid solutions are corrected only in these subdomains, where uniform grids are used. Difference schemes are constructed to converge "almost £-uniformly" — with an error weakly depending on the parameter £.

Введение

Решения уравнений с частными производными, содержащих малый (возмущающий) параметр е при старших производных, имеют особенности типа пограничных и переходных слоев, что приводит к большим ошибкам при их численном решении классическими сеточными методами. Для сингулярно возмущенных задач хорошо известна проблема разработки специальных сеточных методов, погрешность которых слабо зависит от величины параметра е, в частности, методов, сходящихся е-равномерно [1-4]. Достаточно хорошо разработан метод построения е-равномерно сходящихся схем на специальных сетках (фиксированных в процессе решения сеточной задачи), априорно сгущающихся в пограничных слоях в случае экспоненциальных (см., например, [5-7]) и экспоненциальных, степенных и логарифмических пограничных слоев [8, 9]. Обзор прямых и итерационных методов построения сеток для произвольных пограничных и внутренних слоев дан в [9].

Заметим, что шаг специальной кусочно-равномерной сетки, описанной в [5], меняется достаточно резко в окрестности пограничного слоя, что, вообще говоря, может повлечь вычислительные трудности из-за сильной неравномерности сеток (соответствующий пример для регулярных краевых задач был построен А. А. Самарским [10]). С другой стороны, известно, что методы экспоненциальной подгонки [3, 4], преимущество которых

* Работа выполнена при финансовой поддержке Российского фонда фундаментальных исследований, грант №98-01-00362.

© Г. И. Шишкин, 2001.

состоит в использовании простейших равномерных сеток, неприменимы для построения е-равномерно сходящихся схем в случае широкого круга краевых задач с параболическими слоями [5, 11, 12]. Представляют интерес адаптивные методы, в частности, разностные схемы на сетках, переизмельчаемых по какому-либо закону в подобластях, в которых вычисленные решения оказываются недостаточно точными, — схемы на апостериорно сгущающихся сетках (см., например, [9, 13-15]). Привлекательны сгущающиеся сетки, являющиеся равномерными на подобластях, подвергающихся переизмельчению.

Особо отметим метод повышения точности приближенных решений регулярных уравнений, в котором применяется последовательное переизмельчение сетки в тех подобластях, где градиенты решения большие [14, 15]. В случае краевых задач, решения которых имеют локальные особенности, все больший интерес вызывает применение численных методов на апостериорно адаптирующихся сетках, позволяющих улучшить глобальную точность приближенного решения (см., например, работы [16, 17] и библиографию в них). В ряде методов процесс локального сгущения сетки определяется по градиентам решений промежуточных сеточных задач.

В работе рассматриваются сеточные аппроксимации задачи Дирихле для сингулярно возмущенного параболического уравнения конвекции — диффузии. При построении схем применяются классические разностные аппроксимации краевой задачи. Для повышения точности приближенного решения сетки переизмельчаются в окрестности пограничного слоя. Сгущение сетки проводится апостериорно — после анализа промежуточных решений, получаемых в процессе вычисления приближенного решения. На основе последовательно локально сгущающихся сеток (равномерных на подобластях, подвергающихся переизмельчению) строятся разностные схемы, позволяющие получать приближения с погрешностью, слабо зависящей от параметра е. Для таких "почти е-равномерно" сходящихся схем выполняется оценка (в равномерной сеточной норме): |и(х,£) — г(х,£)| < М[М-1/2 + е~иМ-1 + М0-1], где N1 + 1 и М0 + 1 — число узлов сетки по х и £ соответственно, V — произвольное число из (0,1].

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

Развитый в настоящей работе подход может быть использован при построении эффективных численных методов для достаточно представительных классов краевых задач с преобладающей конвекцией. Вызывают интерес схемы на адаптирующихся сетках для сингулярно возмущенных уравнений с конвективными членами, вырождающимися на границе по степенному закону, рассмотренных в [7].

1. Постановка задачи. Цель исследования

В области С с границей Б = С \ С, где

С = В х (0,Т], В = {х : 0 < х <

(1.1)

рассмотрим краевую задачу для сингулярно возмущенного параболического уравнения

f д2 д д 1

L(1 2)u(x, t) = < ea(x, t)^^ + b(x, t) —--c(x, t) — p(x, t) — > u(x, t) = f (x, t), (x, t) G G,

I dx2 dx dt I

(1.2)

u(x,t) = ip(x,t), (x,t) G S.

Здесь a(x,t), b(x,t), c(x,t), p(x,t), f (x,t), (x,t) G G, p(x,t), (x,t) G S — достаточно гладкие функции, причем ao < a(x,t) < a0, b0 < b(x,t) < b0, c(x,t) > 0, p0 < p(x,t) < p0, (x,t) G G, a0, bo, p0 > 0; e — параметр, принимающий произвольные значения из полуинтервала (0,1]. Запись L(jk) (M(jk), f(j.k)(x,t)) означает, что эти операторы, постоянные и функции введены в формуле (j.k). Предполагаем, что на множестве y0 = Г х {t = 0}, где Г = D \ D, в угловых точках (0, 0), (d, 0) выполнены условия согласования ([10] или [18]), обеспечивающие гладкость решения. Считаем, что граница S состоит из множеств SL и S0: S = S0 U SL, где SL — боковая граница, SL = (J S2"; и SLl — соответственно левая и правая границы множества G; S0 — нижнее основание множества G (S0 = S0).

При e ^ 0 появляется пограничный слой в окрестности боковой границы Sf", на которую направлен конвективный поток.

Отметим некоторые проблемы, возникающие при численном решении задачи (1.2), (1.1) с использованием классической разностной схемы (монотонной схемы А. А. Самарского). На множестве G введем прямоугольную сетку

Gh = Ui х й0, (1.3)

где U1 и U0 — сетки на отрезках [0, d] и [0,T] соответственно; U1 — сетка с произвольным распределением узлов, удовлетворяющим лишь условию h < MN-1 (h = maxi h%, h% = xl+l — xi, xi, xi+1 G U1); U0 — равномерная сетка с шагом ht = TN—1. Здесь N +1 и N0 + 1 — число узлов сеток U1 и U0 соответственно. Здесь и ниже через M, Mi (m, mi) обозначаем достаточно большие (малые) положительные постоянные, не зависящие от e. В случае сеточных задач эти постоянные не зависят и от параметров разностных схем. Наибольший интерес для нас будут представлять равномерные по x и t сетки

Gh = Gh(1.3), (1.4)

где и1 — равномерная сетка с шагом h = dN-1.

Для решения задачи используем неявную схему [10]

Л(15) z(x, t) = \^ea(x, t)5xx + b(x, t)8x — c(x, t) — p(x, t)&t} z(x, t) = f (x, t), (x, t) G Gh,

z(x,t) = ip(x,t), (x,t) G Sh, (1.5)

где Gh = GQGh; Sh = SQGh; Sxxxz(x,t) и 5xz(x,t), 5xz(x,t) — вторая и первые (вперед и назад) разностные производные по x, 5xsz(x,t) = 2(h% + hi-1) 1 [$x — 8x] z(x,t), x = xi; &j;z(x,t) — разностная производная (назад) по t. В случае равномерных сеток вида (1.4) с шагом h будем использовать (вместо 8x$z(x, t)) вторую разностную производную 5xxz(x, t) = (z(x + h,t) — 2z(x, t) — z (x — h, t))/h2.

Для решений схем (1.5), (1.3) и (1.5), (1.4) справедлива оценка

lu(x,t) — z(x, t)| < M[(e + N-1)-1N-1 + N0-1], (x,t) G Gh. (1.6)

Из оценки (1.6) вытекает, что решения схем (1.5), (1.3) и (1.5), (1.4) сходятся при условии к = о(е) или

е-1 = о(Ж). (1.7)

Если условие (1.7) нарушается, например, при е-1 = 0(М), то, вообще говоря, решения разностных схем (1.5), (1.3) и (1.5), (1.4) при N N ^ то не сходятся к решению задачи (1.2).

В связи с таким поведением сеточных решений возникает интерес к построению специальных разностных схем, погрешность решений которых не зависит от величины параметра е, в частности схем, сходящихся при более слабом условии, чем (1.7) — условие сходимости решений сеточных задач (1.5), (1.3) и (1.5), (1.4).

Сформулируем цель работы.

Пусть ¿(ж,*), (ж,*) € — решение некоторой разностной схемы. По определению эта схема сходится равномерно по параметру е (или е-равномерно), если для функции ¿(ж,*) выполняется оценка |и(ж,*) — ¿(ж,*)| < Л(Ж-1,Ж0-1), (ж,*) € где величина Л(Ж-1,Ж0"1) стремится к нулю при N N ^ то равномерно по параметру. Будем говорить, что решение схемы сходится почти е-равномерно, если для любого сколь угодно малого числа V > 0 найдется функция Л(е-^N-1,Ж0-1), такая, что для сеточной функции ¿(ж,*) выполняется оценка

|и(ж,*) — г(ж,*)|< МЛ(е-^Ж-1,Ж0-1), (ж,*) € Си, (1.8)

где Л(Ж-1,Ж-1) ^ 0 при Ж, N ^ то равномерно относительно параметра е. В случае оценки (1.8) будем говорить, что решение схемы сходится почти е-равномерно с дефектом V. Дефект схем (1.5), (1.4) и (1.5), (1.3) равен единице.

Таким образом, в случае задачи (1.2), (1.1) приходим к проблеме построения схем с дефектом, меньшим единицы, в частности, к проблеме построения почти е-равномерно сходящихся схем.

2. Сеточные аппроксимации на локально переизмельчаемых сетках

Приведем алгоритм построения локально переизмельчаемой (в пограничном слое) сетки. На областях, подвергающихся переизмельчению, этот алгоритм использует равномерные сетки по пространству и времени (сетка по времени не переизмельчается).

Предварительно опишем формальный итерационный алгоритм построения приближенных решений краевой задачи (1.2), (1.1), который будем использовать при конструировании разностных схем. Пусть на множестве С найдена функция м1(ж,^) — приближение решения краевой задачи, причем

|и(ж,*) — и1 (ж,*)|< М8, (ж,*) € С, ^ < ж, (2.1а)

где 8 > 0 — произвольное малое число; постоянная М не зависит от 8; € (0,^). Пусть М(2)(ж,*), (ж,*) € С(2), где С(2) = ^(2) х (0,Т], Д(2) = (0,^1), есть решение задачи

¿(1.2) и(2)(М) = /(М), (М) € С(2),

(2.1Ь)

( *) =) м1(ж,*), (ж,*) € ^(2) \

^ ) =1 ¥>М), (ж,*) € 5(2) П

Здесь $(2) = С(2) \ С(2). Пусть и2(ж,С), (ж, С) € С(2) — приближение решения Ч(2)(ж,С), удовлетворяющее неравенству | Ч(2)(ж,С) — и2(ж,С) | < М^, (ж, С) € С(2), < ж. Полагаем

(ж,С) € С(2),

Ч2(ж, С) = ^ — —

их(х, ¿), (ж, С) € С \ С(2).

Тогда |и(ж,С) — и2(х, С)| < М#, (ж, С) € С, < ж.

Пусть при к > 3 построена функция Чк_1(ж,С), (ж, С) € С, имеющая удобное для вычислений представление при ¿к_1 < ж, ¿к_1 € (0,^), причем

|и(ж, С) — ик-1(ж,С)| < М5, (ж, С) € С, < ж. (2.1с)

Постоянная М зависит от к. Через Ч(к)(ж,С), (ж, С) € С(к), где С(к) = Д(к) х (0,Т], = (0,^^_1), обозначим решение задачи

¿(1.2) Ч(к)(ж,С) = /(ж,С), (ж,С) € С(к^

(2.1а)

( ,) = I мл_1(ж, С), (ж, С) € $П) \ Ч(к) ) [ ^(ж,С), (ж, С) € % П $,

где $(к) = С(к) \ С(к). Пусть Ч^(ж,С), (ж,С) € С(к) — приближение функции Ч(к)(ж,С) и

| Ч(к)(ж,С) — ик(ж, С) | < М5, (ж, С) € С(к), < ж. Полагаем

^^ (ж,С) € С(^

Чк (ж, с) = ^ — —

Чк_1(ж,С), (ж, С) € С \ С(к).

Для функции ик(ж, С), (ж, С) € С имеем оценку |и(ж,С) — ик(ж, С)| < М#, (ж, С) € С, ¿к < ж. Если при некотором значении к = К0 оказалось, что ¿К0 = 0, то ¿к = 0 при к > К0. При к > К0 + 1 множества С(к) считаем пустыми и функции Ч(к)(ж, С) не определяем. Например, при к > К0 имеем ик(ж,С) = иКо(ж,С), (ж,С) € С.

При к = К, где К — заданное фиксированное число (К > 1), полагаем

иК(ж,С) = иК(ж,С), (ж,С) € С. (2.1е)

Функции чК(ж, С) и ик(ж, С) (к = 1, ... , К, (ж, С) € С) назовем решением и соответственно компонентами решения итерационного процесса (2.1).

Функции чК(ж, С) имеют удобное для вычислений представление на подобластях, расширяющихся с ростом К. Для функций Чк (ж, С) и чК (ж, С) = Чк (ж, С) справедлива оценка

|и(ж,С) — Чк(ж,С)| < М5, (ж, С) € С, ¿к < ж, к =1, ... ,К. (2.2)

Лемма 1. Для функций чК(ж,С) = ик(ж,С) и чк(ж,С), (ж,С) € С, к = 1, ...,К (решения итерационного процесса (2.1) и его компонент) выполняется оценка (2.2).

Опишем сеточную конструкцию, аппроксимирующую итерационный процесс (2.1). На множестве С введем грубую (исходную) сетку

Сш = х (2.3а)

где и w0 — равномерные сетки, w0 = ^0(1.3); шаг сетки w есть hi = dN-1. Решение задачи (1.5), (2.3a) обозначим через z1(x,t), (x,t) G G1h, где G1h = G1h(2.3) = Gh(1.4).

Пусть каким-то образом найдена величина d1 G W1, такая, что при d1 < x сеточное решение z1(x,t), (x,t) G G1h хорошо приближает решение задачи (1.2), (1.1). Если окажется, что d1 > 0, то определим подобласть, на которой будем переизмельчать сетку:

G(2) = G(2)(d1), G(2) = D(2) X (0,T], D(2) = (0,d1). (2.3b)

На подобласти G(2) введем сетку G(2)h = W(2) X w0, где W(2) — равномерная сетка с числом узлов N +1. На множестве G(2)h найдем решение z(2)(x,t) сеточной задачи

Л(1.5) Z(2)(x, t) = f (x, t), (x, t) G G(2)h,

( t) i Z1(x,t), (x,t) G S(2)h\ S, Z(2) ) { ^(x,t), (x,t) G S(2)h П S,

где G(2)h = G(2) f| G(2)h, S(2)h = S(2) f| G(2)h, S(2) = G(2) \ G(2). Сеточное множество G2h и функцию z2(x,t), (x, t) G G2h определим соотношениями

G = G Mir \ G 1 z (xt)=Z Z(2)(x,t), (x,t) G G(2)h, G2h = G(2)hU |G1h \ G(2)j , z2(x, t) Л

[ Z1(x,t), (x, t) G G1h \ G(2).

Пусть при k > 3 уже построены сеточные множество Gfc_1;h и функция Zk-1(x,t) на этом множестве, а также каким-то образом найдена величина dk_1 G такая, что при

dk_1 < x сеточное решение Zk_1(x,t), (x,t) G Gfc_1;h хорошо приближает решение задачи (1.2), (1.1). Здесь — сетка, порождающая сетку Gfc_1;h, Gfc_1;h = X w0; N + 1 — число узлов сетки , k > 1; N1 = N. Если окажется, что dk_1 > 0, то определим область

G(k) = G(fc)(dfc_1), G(fc) = D(fc) x (0,T], D(fc) = (0,dfc_1). (2.3c)

На множестве G(k) введем сетку

G(k)h = w(fc) X wo, (2.3d)

где W(k) — равномерная сетка с числом узлов N + 1; h(k) — шаг сетки W(k). Пусть z(k)(x, t), (x, t) G G(k)h — решение сеточной задачи

Л(1.5) z(k)(x,t) = f (x,t), (x,t) G G(k)h

r (2.3e)

( t) I zfc_1(x,t), (x,t) G S(k)h \ S, Z(k) ) [ ^(x,t), (x,t) G S(k)h n S.

Полагаем = G(fc)h U{Gfc_1;h \ G(fc)},

z (x t) = f z(k)(x,t), (x,t) G G(fc)h,

' I zfc_1(x,t), (x, t) G Gfc_1;h \ G(fc).

Если при некотором значении k = K0 оказалось, что dx0 = 0, то полагаем dk = 0 при k > K0. При k > K0 + 1 множества G(k) считаем пустыми и функции z(k)(x, t) не вычисляем. Например, при k > K0 имеем zk(x,t) = zKo (x,t), Gkh = GKoh.

При к = К, где К — заданное фиксированное число, К > 1, полагаем

= Окн = Он, гк (ж,*) = гк (ж,*) = г(х,г). (2.3£)

Функцию г(2.з)(х,*), (х,*) £ Сн(2.з) назовем решением схемы (1.5), (2.3), а функции гк(ж,*), (ж,*) £ Окь, к =1, ..., К — компонентами решения разностной схемы.

Приведенный алгоритм (назовем его А(2 3)) позволяет строить сетки, сгущающиеся в пограничных слоях. Величина Мк + 1 — число узлов сетки йк = и к, используемой при построении функции гк (ж,*), — не превосходит величины N (К) = К (М + 1).

В схемах (1.5), (2.3) при решении промежуточных задач (2.3е) не требуется интерполяция для определения значений функций г»(ж,*) на границе £(к)н.

Сетки Окь (к = 1, ... , К) получаемые по алгоритму А(2 3), определяются законом выбора величин ¿к (к = 1, 2, ... ,К), а также величинами К и М, М0. Таким образом, алгоритм А(2.з) определяет класс разностных схем (1.5), (2.3). В этом классе схем граница подобласти, в которой проводится переизмельчение сетки, проходит через узлы более грубой измельчаемой сетки. Заметим, что самый маленький шаг сетки ик = йк не меньше величины ¿М-к. В сетках, получаемых по алгоритму А(2 3), величины ¿к будем находить на основе промежуточных результатов вычислений. Таким образом, сетки Окн относятся к апостериорно сгущающимся сеткам.

Для разностных схем из класса (1.5), (2.3) справедлив принцип максимума. Заметим, что в этом классе не существует схем, решения которых сходятся е-равномерно к решению краевой задачи (1.2), (1.1).

3. Адаптивная схема на основе оценки градиента решения

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

Предварительно приведем оценки решения краевой задачи и их производных [5]. Решение задачи (1.2) представим в виде суммы функций

и(ж,*) = и(ж,*) + V(ж,*), (ж,*) £ О, (3.1а)

где и (ж,*), V (ж,*) — регулярная и сингулярная части решения;

и (ж,*) = ио(ж,*) + еи^ж,*) + Уи (ж,*),

(ж,*) £ О.

(3.1Ь)

Функции и0(х,1), и1(х,1), уи(х,Ь), V(х,Ь) находятся из решения задач

Г д д

¿(3.1) ио(х,г) = 1 Ь(х,г)— - с(х,г) - р(х,г) — \ио(х,г) = f (х,г), (х,г) е си^,

щ(х,г) = (р(х,г), (х,г) е Б \ Б^; д2

¿(3.1) и1(х,г) = -а(х,г) — и0(х,г), (х,г) е СиБ^, и1(х,г) = 0, (х,г) е Б \ Б[; 2 д2

Ь(1ф2) уи(х,Ь) = —£2 а(х,Ь) —2 и1(х,Ь), (х,Ь) е С, уи(х,Ь) = 0, (х,Ь) е Б;

ь(1.2) V(х,г) = 0, (х,г) е С, V(х,г) = ^(х,г) - и(х,г), (х,г) е Б. Для функций п(х,Ь), и(х,Ь), V(х,Ь) справедливы оценки

(3.1с)

дк+к° * \

< Ые-к,

дк+ко

дхк дЬко

и(х,г) < м[1 + е2-к]

дк+к

дхк дЬко

(3.2)

V(х,г) < Ме-к ехр(-ш£-1т(х,Г1)), (х,г) е С, к + 2 ко < 4, к < 3,

где ш — произвольное число из интервала (0,шо), ш0 = (а0)-1Ь0; постоянная М = М(ш), вообще говоря, увеличивается с ростом ш; г(х,Г1) — расстояние от точки х до левой границы Г1 множества О.

Теорема 1. Пусть для данных краевой задачи (1.2), (1.1) выполняется условие а, Ь, с, р, f е С5+а,1+а(С). Тогда для решения краевой задачи и его компонент из представления (3.1а) справедливы оценки (3.2).

Замечание. Для функции V(х,Ь) справедлива также оценка IV(х,Ь)1 < МШ(х,Ь), (х,Ь) е С, где Ш(х,Ь) — решение задачи

¿(1.2) Ш(х,г) = 0, (х,г) е С, Ш(х,г) = 1, (х,г) е Б[, Ш(х,г) = 0, (х,г) е Бо и

Приведем конструкцию индикатора на основе оценки градиента решения.

Определим ширину пограничного слоя в случае краевой задачи (1.2), (1.1). Пусть для компоненты и(х,Ь) из представления (3.1а) выполняется оценка

д

тхи (хЛ)

< мъ (х,г) е С.

(3.3а)

Предположим, что значения параметра е достаточно малы, е < е0. Будем говорить, что а0 = а0(е, М0) есть ширина пограничного слоя в окрестности стороны Б^, если а0 есть минимум величины а, для которой справедлива оценка

д

шу (х л)

< Мо, (х,г) е С, г(х,б£) > а,

(3.3Ь)

где М0 — произвольное достаточно большое число (М0 > М1).

В качестве модельного примера рассмотрим краевую задачу

еи"(ж)+ п'(ж) = 1, ж € П = (0,1), п(0) = п(1) = 1. (3.4)

Решение задачи (3.4) представим в виде суммы функций: п(ж) = и (ж) + V (ж), ж € П, где и (ж) и V (ж) — регулярная и сингулярная части решения; пограничный слой возникает в окрестности точки ж = 0.

Для ширины слоя 0о(з.з) имеем асимптотическое поведение а0 ~ е 1пе-1 при е = о(1). Справедлива также оценка а0 < Ме 1п(е-1М-1), е € (0,е0], е0 = е0(М0), е0 < тМ0-1, где М0, т — произвольные постоянные, удовлетворяющие условиям М0 > 1, т < 1, М> 1.

Определим ширину пограничного слоя в случае разностной схемы (1.5), (1.3). Пусть ¿V(ж,*), (ж,*) € есть решение разностной задачи

Л(1.5) ¿(ж,*) = ¿(12) г>(ж,*), (ж,*) € г (ж, *) = г>(ж, *), (ж,*) €

где ^(ж,*) — произвольная достаточно гладкая функция, V € С2,1 (С) Р| С (С). Решение задачи (1.5), (1.3) представим в виде суммы функций

¿(ж,*) = ¿и(ж,*) + ¿у(ж,*), (ж,*) € (3.5)

Здесь ¿и (ж, *) и ¿у (ж, *) — сеточные функции, приближающие компоненты и (ж, *) и V (ж, *) из представления (3.1а). Пусть для компоненты ¿и(ж,*) выполняется оценка

| 8х ¿и(ж,*) |< М1, (ж,*) € Си. (3.6а)

Будем говорить, что а0 = а0(М0) = а0(М0; Д,е,к) (е € (0,е0], М0 и е0 — достаточно большая и малая постоянные, М0 > М1, е0 = е0(М0)) есть ширина сеточного пограничного слоя

оценка

в окрестности стороны если а0 есть минимум величины а, для которой выполняется

|8^у(ж,*)| < М0, (ж,*) € г(ж,5^) > а. (3.6Ь)

Таким образом, функция а0(М0) = а0(3.6)(М0; Д,е,к) построена. В случае разностной схемы

е8х$ ^(ж) + 8Х ¿(ж) = 1, ж € Пи, ¿(0) = ¿(1) = 1,

аппроксимирующей задачу (3.4), для ширины пограничного слоя на равномерных (с шагом к) сетках имеем асимптотику

е 1п е-1, к < Ме,

к 1п-1(1 + е-1к)1пк-1, к > те; е, к = о(1).

Справедлива оценка

а0 < М [е 1пе-1 + к 1пк-1] , е € (0,е0], к < к0,

где е0 = е0(М0), к0 = к0(М0) — достаточно малые величины.

Для того чтобы формальная сеточная структура (1.5), (2.3) стала конструктивной, требуется задать величины К и ^, к = 1, 2, ..., К. Пусть К > 1. Определим величины ^(2.3). Полагаем

= а1, (3.7а)

где 01 = 0о(з.6)(М0; О,е,Н); О = О^л); Н = Н(^4); М0 — достаточно большое число. Величину ¿к_1 считаем уже найденной. Далее находим величину

0к = ао(кМо; О(к),е, Н(к)), к > 2, (3.7Ь)

где 0о(з.7) (М; О,е,Н) = 0о(з.6)(м; О,е,Н), О(к) = О(к)(2.з), Н(к) = Н(к)(2.з). Если выполняется соотношение 0к < т0 0к_ 1, то полагаем

¿к = 0к. (3.7с)

Здесь т0 — достаточно малое число. Если же при некотором значении к = к0 оказалось, что 0ко > т0 0ко_1, то полагаем = ¿ко при к > к0.

Разностная схема (1.5), (2.3), (3.7) — схема на адаптирующихся сетках, которые строятся на основе оценки градиента сеточных решений, получаемых в процессе промежуточных вычислений. Переизмельчение сетки проводится лишь в окрестности пограничных слоев; диаметр такой окрестности (ширина пограничного слоя) сужается с ростом величины к.

4. Исследование схемы (1.5), (2.3), (3.7)

Приведем некоторые оценки для решения разностной схемы (1.5), (1.4). Через (ж, С) обозначим решение задачи

Л г(ж,С) = {еа(1)&кг + 6(1)4} г(ж,С) = 0, (ж, С) € С \

г(ж,С) = 1, (ж, С) € г(ж, С) = 0, (ж, С) €

где С = С^(1.з); а(1) = шах^ а(ж,С); 6(1) = шт^ 6(ж,С). Для функции (ж, С) = (ж) справедлива оценка

^(ж) < , (ж, С) € (4.1а)

где г1 = г(ж, $"); д = 1 + а_1) 6(1) е_1 Н. Таким образом, имеем оценку

!Мехр(—те_1 г1), Н < Ме, 1 _

, , >, ж € Г1 = г(ж,$Г). (4.1Ь)

М(еН ) - 1г1, Н>те / ' 1 17

Функция (ж) является мажорантой для компоненты гу (ж, С) из представления (3.5)

|гу(ж, С)| < М^(ж), (ж, С) € (4.2)

Для решения разностной схемы (1.5), (1.4) получаем оценки

|и(ж,С) — г(ж,С)| < М [(е + N_1)_1 N_1 + Ж0_1] , (4.3а)

|и(ж,С) — г(ж,С)| < М [г^(ж) + N_1 + N,7^ , (ж, С) € (4.3Ь)

Из оценок (4.2), (4.3) следует, что при условии (1.7) схема сходится на С^, а также сходится е-равномерно вне 00-окрестности множества

|и(ж,С) — г(ж,С)|< М [N_1/2 + N0_1] , (ж, С) € г(ж, > 00, (4.4а)

где 00 = 00(3.6)(М0; О, е, Н),

00 < М [е 1п е_1 + N_11п Щ . (4.4Ь)

Окрестность, вне которой справедлива оценка (4.4а), сужается при е ^ 0, N ^ то.

Теорема 2. Пусть для решения краевой задачи (2.2), (2.1) выполняются оценки теоремы 1. Тогда решение разностной схемы (1.5), (1.4) сходится на С к решению краевой задачи при условии (1.7), а также е-равномерно (со скоростью -1/2 + Л^-1)) вне а0-окрестности множества Б[. Для сеточного решения справедливы оценки (4.2) - (4.4).

Рассмотрим разностную схему (1.5), (2.3), (3.7). Для компоненты ^(ж,*) = г(ж,*) решения этой схемы выполняются оценки (4.3). С учетом оценки (4.1) для функции (ж) находим оценку величины а1 — ширины пограничного слоя:

а1 < М [е 1пе-1 + N-11пN] , е € (0,ео], Н < Н0. Для величины Н2 — шага сетки ^(2) — имеем

Н2 < МЛ-1 [е 1п е-1 + N-1 1п.

Принимая во внимание неравенство (4.3Ь), оценим и(ж,*) — ^(ж,*) на границе множества С(2)й, а также и(ж, *) — г(2)(ж,*) — на самом множестве С(2)^. Для и(ж, *) — г2(ж,*) получается оценка

|и(ж,*) — *2(М)| < М [Л-1/2 + е-1Л-21пN + Л0-1] , (ж,*) €

а вне а2-окрестности множества Б^ — оценка

|и(ж,*) — г2(ж,*)|< М [N-1/2 + ^-1] , (ж,*) € ^ г(ж,Б1ь) > а2.

Для величины а2 имеем

а2 < М N-1 1п N.

Подобным образом находим

|и(ж,*) — гк(ж,*)| < М ^-1/2 + е-1N-к 1пк-1 N + ^-1] , (ж,*) € С^,

|и(ж,*) — ^(ж,*)| < М [N-1/2 + ^-1] , (ж,*) € С^, г(ж,Б1) > ак,

М [е 1п е-1 + N-11п N1, к = 1, ак < ^ ^ , к = 1,2, ..., К; (4.5)

MN-к+11пк-1 N, к > 2

|и(ж,*) — г(ж,*)| < М [N-1/2 + е-1N-к 1пк-1 N + , (ж,*) € С, |и(ж,*) — г(ж,*)|< М [N-1/2 + ^-1] , (ж,*) € СН} г(ж,Б1) > ак,

Г М [е 1п е-1 + N-11п N1, К = 1, 1 ак < < > , (4.6)

[ MN-к+11пк-1 N, К > 2 ]

где г(ж,*) = г(2.з)(ж, *), С = С^(2.з).

Функции г(ж,*) и (ж,*) при N, N0 ^ то сходятся (к решению краевой задачи (1.2), (1.1)) е-равномерно вне ак- и а^-окрестностей множества Б^, а также на множествах и Сдл при достаточно малых (но не слишком малых) значениях параметра е, а именно, при условии

е > ео(М), е-1(М) = о(Мк 1п-к+1 N),

е > ек(М), е-1(М) = о(Мк 1п-к+1 N), к = 1, 2, ..., К. (4.7)

Таким образом, разностная схема (1.5), (2.3), (3.7) — схема на адаптирующихся сетках — сходится почти е-равномерно. Для того чтобы обеспечить дефект сходимости функции г(ж, *) не выше величины ^(1.8), требуется выбрать величину К удовлетворяющей условию

К > К(V), К(V) = V-1. (4.8)

Теорема 3. Пусть выполняется предположение теоремы 2. Тогда функции г(ж,*), (ж,*) £ и гк(ж,*), (ж,*) £ Окн, к = 1, ... , К — решение разностной схемы (1.5), (2.3), (3.7) и его компоненты — сходятся на О к решению краевой задачи (1.2), (1.1) при условии (4.7), а также е-равномерно (со скоростью 0(М-1/2 + М—1)) вне ак-окрестности множества ; решение схемы (1.5), (2.3), (3.7), (4.8) сходится к решению краевой задачи почти е-равномерно с дефектом V. Для сеточных решений справедливы оценки (4.5), (4.6).

Для иллюстрации сеток и сеточных решений, получаемых при измельчении сетки на основе оценки градиента решения, рассмотрим модельный пример для обыкновенного дифференциального уравнения

( й2 й \

Ьм(ж) = + и(ж) = /(x), ж £ Д

u(0) = u(1) = 0, (4.9)

где D = (0,1), f (x) = — 2(x + е). Точное решение этой краевой задачи имеет вид

1 — exp(—x/е) 2 U(X) = 1 — exp(—1/е) — Х •

На рисунке приведены функции Zk(x) на сетках Dkh при k = 1, 2, решение u(x) краевой задачи (4.9) для е = 2-4, N = 8 при М0(з.бь) = 4, а также конструкции сеток D^ и D2h. Здесь di = 0.25.

Функции Zk (x), x G Dkh можно представить в виде суммы функций

Zk(x) = Zku(x) + Zkv(x), x G Dkh, k = 1, 2,

где компоненты Zku (x) и Zkv (x) отвечают регулярной U(x) и сингулярной V(x) частям решения u(x) задачи (4.9). Заметим, что функции Zku(x) при N ^ то сходятся к компоненте U(x) е-равномерно, где U(x) = [1 — exp(—1/е)]-1 — x2. Ошибки для сингулярных компонент z1v(x), x G D1h и z2 v(x), x G D2h: E(N, е) = max |ZkV(x) — V(x)|, k =1, 2 приведены

xeDfch

в таблице, откуда следует, что компонента z1v(x) (а значит, и функция Z1(x)) сходится при N-1 ^ е, а компонента z2v(x) (и функция z2(x)) — уже при N-1 ^ е1/2, т.е. при более слабом условии. Таким образом, приведенные численные результаты согласуются с теоретическими, соответствующими стационарным задачам.

О 0.2 0.4 0.6 0.8 1

Рис. 1. Графики функций ^(х), х Е — Ш, ^(х), х Е — х и точное решение задачи (4.9) — сплошная линия (а); конструкции сеток (б).

Ошибки для сингулярных компонент

£ Число узлов сетки (Ж)

4 16 64 256

к = 1

2 -0 1.2703Е- -02 3.6236Е- - 03 9.3435Е- -04 2.3537Е- - 04

2 -2 1.1058Е- -01 3.6415Е- - 02 9.9301Е- -03 2.5405Е- - 03

2 -4 1.8040Е- -01 1.3211Е- - 01 4.1720Е- -02 1.1206Е- - 02

2 -6 5.8812Е- -02 1.8168Е- - 01 1.3212Е- -01 4.1721Е- - 02

2 -8 1.5385Е- -02 5.8823Е- - 02 1.8168Е- -01 1.3212Е- - 01

2- -10 3.8911Е- -03 1.5385Е- - 02 5.8823Е- -02 1.8168Е- - 01

2- -12 9.7561Е- -04 3.8911Е- - 03 1.5385Е- -02 5.8823Е- - 02

2- -14 2.4408Е -04 9.7561Е- - 04 3.8911Е- -03 1.5385Е- - 02

к =2

2 -0 1.2703Е -02 3.6236Е- - 03 9.3435Е- -04 2.3537Е- - 04

2 -2 1.1058Е- -01 3.6415Е- - 02 9.9301Е- -03 2.5405Е- - 03

2 -4 2.2364Е -01 7.9490Е- - 02 3.5941Е- -02 1.0469Е- - 02

2 -6 2.2753Е -01 1.0028Е- - 01 4.5345Е- -02 1.6201Е- - 02

2 -8 7.3292Е -02 2.0030Е- - 01 3.7098Е- -02 9.8053Е- - 03

2- -10 1.9216Е -02 1.1099Е- - 01 7.8487Е- -02 1.2224Е- - 02

2- -12 4.8629Е -03 3.0318Е- - 02 1.9816Е- -01 3.2092Е- - 02

2- -14 1.2195Е -03 7.7529Е- - 03 1.1079Е- -01 7.6696Е- - 02

5. Замечания и обобщения

1. Утверждение теоремы 3 сохраняется и в том случае, когда величина а0(3.6) в (3.7Ь) заменяется на а0 — минимум величины а, для которой выполняется оценка

|Мж,*)| < М0, (ж,*) € Си, г(ж,^) > а,

где М0 = М0(з.б); а0 = а0(Ма) = а^(М0; Д,е, к).

2. Пусть ¿(ж,*), (ж,*) € С — функция, построенная линейной интерполяцией по ж,* по значениям решения схемы (1.5), (2.3), (3.7). При условиях теоремы 3 справедлива оценка

|п(ж,*) — ¿(ж,*)| < М[Ж-1/2 + е-1Ж-к 1пк-1 N + Ж0-1], (ж,*) € С; |п(ж, *) — ¿(ж,*)|< М [Ж-1/2 + Ж0-1], (ж,*) € С, г(ж,^) > а к.

3. Приведенные апостериорные сетки применимы для построения почти е-равномерно сходящихся схем и в случае квазилинейных уравнений. На множестве С(1.1) рассмотрим квазилинейное параболическое уравнение

¿(1.2) п(ж, *) = #(ж, п(ж, *)), (ж, *) € С,

(5.1)

п(ж,*) = <^(ж, *), (ж,*) € 5.

Коэффициенты уравнения и функция <^(ж,*) удовлетворяют условиям, указанным в постановке задачи (1.2), (1.1); функция $(ж,*, п) — достаточно гладкая на множестве С х Л,

5 —

причем —М < — д(ж,*,п) < то, (ж,*,п) € С х Л. дп

Для решения задачи (5.1) справедливы оценки, подобные полученным для задачи (1.2) [19] (стационарный случай). Задаче (5.1) сопоставим разностную схему

¿(1.5)^ (ж, *)= #(ж,*, £ (ж,*)), (ж,*) € Си, (5.2)

¿(ж,*) = <^(ж,*), (ж,*) € Си = Си(1.з).

Решения разностных схем (5.2), (2.3), (3.7) сходятся к решению краевой задачи (5.1), (1.1) почти е-равномерно с дефектом V (е-равномерно вне -окрестности множества 5^). Для решений этих схем справедливы утверждения, подобные утверждениям теоремы 3.

4. Использование в качестве индикаторов разностных производных порядка выше первого позволяет повысить скорость сходимости схемы. Рассмотрим индикатор на основе разностных производных второго порядка.

Пусть | 8Хх ¿и(ж, *) | < М2, (ж, *) € Си. Будем говорить, что а0 = а0(М0) = а0(М0; Д е, к) есть ширина пограничного слоя по производной второго порядка, если а0 есть минимум величины а, для которой выполняется оценка

| 8хг¿у(ж,*) | < М0, (ж,*) € Си, г(ж,5^) > а, (5.3а)

где М0 > М2. В алгоритме А(2.3) величины а0(3.7) определим соотношением

а0(з.7)(М; Де,к) = а0(5.з«)(М; Де,к).

(5.3Ь)

Для решений разностной схемы (1.5), (2.3), (3.7с), (5.3) получаем оценки

|и(ж,С) — г(ж,С)| < М [ N_2/3 + _к 1пк_1 N + N0_1 ] , (ж, С) € С^;

|и(ж,С) — г(ж,С)|< М [ N_2/3 + N0_1 ] , (ж, С) € СН} г(ж, > 0к,

Г М [ е 1п е_1 + N_11п N ] , К =1,1 0К < | М N_к+11пк_1 N, К > 2 у

Величина К>^_1 обеспечивает дефект сходимости функции г(ж,С) не выше величины V.

Автор выражает признательность участникам IV Сибирского конгресса по индустриальной и прикладной математике ИНПРИМ-2000 (Новосибирск, 26 июня - 1 июля 2000 г.) за плодотворные обсуждения вычислительных методов для задач с негладкими решениями, а также Л. П. Шишкиной за проведение расчетов.

Список литературы

[1] Allen D. N., Southwell R. V. Relaxation methods applied to determine the motion, in two dimensions, of viscous fluid past a fixed cylinder // Quart. J. Mech. and Appl. Math. 1955. V. 8. P. 129-145.

[2] Бахвалов Н. C. К оптимизации методов решения краевых задач при наличии пограничного слоя // Журн. вычисл. математики и мат. физики. 1969. Т. 9, №4. С. 841-859.

[3] Ильин A. M. Разностная схема для дифференциального уравнения с малым параметром при старшей производной // Мат. заметки. 1969. Т. 6, вып. 2. С. 237-248.

[4] Дулан Э., Миллер Дж., Шилдерс У. Равномерные численные методы решения задач с пограничным слоем. M.: Мир, 1983.

[5] Шишкин Г. И. Сеточные аппроксимации сингулярно возмущенных эллиптических и параболических уравнений. Екатеринбург: Изд-во УрО РАН, 1992.

[6] Miller J. J. H., O'Riordan E., Shishkin G. I. Fitted Numerical Methods for Singular Perturbation Problems. Singapore: World Sci., 1996.

[7] Farrell P. A., Hegarty A. F., Miller J.J. H. ET AL. Robust Computational Techniques for Boundary Layers. Boca Raton: CRC Press, 2000.

[8] ЛисЕйкин В. Д., Петренко В.Е. Адаптивно-инвариантный метод численного решения задач с пограничными и внутренними слоями. Новосибирск: ВЦ СО АН СССР, 1989.

[9] Liseikin V. D. Grid Generation Methods. Berlin: Springer-Verlag, 1999.

[10] Самарский А. А. Теория разносных схем. М.: Наука, 1989.

[11] Шишкин Г. И. Аппроксимация решений сингулярно возмущенных краевых задач с параболическим пограничным слоем // Журн. вычисл. математики и мат. физики. 1989. Т. 29, №7. С. 963-977.

[12] SHISHKIN G. I. On finite difference fitted schemes for singularly perturbed boundary value problems with a parabolic boundary layer //J. Math. Anal. and Appl. 1997. V. 208, No. 1. P. 181-204.

[13] Roos H.-G., Stynes M., TüBISKA L. Numerical Methods for Singularly Perturbed Differential Equations. Heidelberg: Springer, 1996.

[14] МАРЧУК Г. И., ШАйДУРОВ В. В. Повышение точности решений разностных схем. М.: Наука, 1979.

[15] Birkhoff G., Lynch R. E. Numerical Solution of Elliptic Problems. Philadelphia: SIAM, 1984.

[16] Adaptive Methods for Partial Differential Equations / J. E. Flaherty, P. J. Paslow, M.S. Shephard, J. D. Vasilakis (Eds). Philadelphia: SIAM, 1989.

[17] MCCORMICK S.F. Multilevel Adaptive Methods for Partial Differential Equations. Philadelphia: SIAM, 1989.

[18] ЛАДЫЖЕНСКАЯ О. А., Солонников В. А., УРАЛЬЦЕВА Н. Н. Линейные и квазилинейные уравнения параболического типа. М.: Наука, 1967.

[19] Шишкин Г. И. Разностная аппроксимация сингулярно возмущенной краевой задачи для квазилинейных эллиптических уравнений, вырождающихся в уравнение первого порядка // Журн. вычисл. матем. и мат. физики. 1992. Т. 32, №4. С. 550-566.

Поступила в редакцию 21 марта 2000 г. в переработанном виде — 24 октября 2000 г.