CC BY
28
УЧЕНЫЕ ЗАПИСКИ КАЗАНСКОГО ГОСУДАРСТВЕННОГО УНИВЕРСИТЕТА
Том 14 7, ки. 3
Физико-математические пауки
2005
УДК 519.6
РЕШЕНИЕ ЗАДАЧИ О ПРЕПЯТСТВИИ МЕТОДОМ ДЕКОМПОЗИЦИИ ОБЛАСТИ
М.А. Игнатьева, A.B. Лапин
Аннотация
Метод декомпозиции области с пепалегающими подобластями применяется к решению задачи о препятствии. Препятствие, а значит и свободная граница, расположено в подобласти исходной области, следовательно, известно, где решение может потерять гладкость. Декомпозиция осуществляется с использованием информации о расположении препятствия.
Строятся сеточпые аппроксимации исходной задачи па несогласованных по подобластям сетках: при аппроксимации задачи в подобласти, содержащей свободную границу, используется более мелкая сетка.
Для решения сеточных задач предлагаются два итерационных метода, которые трактуются как нелинейные аналоги метода расщепления Дугласа Рэкфорда. Доказывается сходимость итераций, обсуждаются вопросы реализации методов, приводятся результаты вычислительных экспериментов.
Введение
Метод декомпозиции области с пепалегающими подобластями основан на разделении области, в которой нужно решить краевую задачу, на подобласти и аппроксимациях задач в подобластях. Это позволяет строить итерационные методы, на каждом шаге которых требуется решать приближенные задачи в подобластях, связанные между собой некоторыми условиями на линиях разрезов области. При сеточных аппроксимациях (методы конечных элементов, конечных разностей, конечных объемов) задачи в подобластях имеют меньшую алгебраическую размерность, чем исходная сеточная задача. Более того, в случае простой геометрии подобластей и при использовании структурированных сеток для решения задач в подобластях можно применять эффективные итерационные, а в ряде случаев и прямые методы. Отметим также, что метод декомпозиции области дает возможность применения различных аппроксимаций в подобластях и естественным образом приспособлен к построению параллельных алгоритмов.
В настоящее время методы декомпозиции области с пепалегающими подобластями наибольшее развитие получили для линейных эллиптических уравнений второго порядка (см.. например, монографии и обзорные статьи [1 7] и труды конференций по методам декомпозиции области [8. 9]). В то же время, построению и исследованию этих методов для нелинейных эллиптических задач со свободными границами посвящены лишь немногочисленные статьи (см. [10 13]).
В настоящей работе метод декомпозиции области с аппроксимацией исходной задачи по методу конечных элементов применяется для решения задачи о препятствии. когда препятствие расположено в некоторой подобласти исходной области. В подобласти, содержащей препятствие, т. е. и свободную границу, решение задачи в общем случае обладает малой гладкостью, поэтому в ней используется более
мелкая сотка. Сотки в подобластях строятся частично согласованными: на границе раздела подобластей множество узлов крупной сетки является подмножеством множества узлов мелкой сетки. На приближенное решение (кусочно полиномиальную функцию из пространства конечных элементов) накладывается условие непрерывности при переходе через границу раздела подобластей. Это условие является дополнительным ограничением в приближенной задаче и при построении итерационных методов либо присоединяется к имеющемуся ограничению в области. либо снимается с помощью множителей Лагранжа. В результате предлагаются два итерационных метода. Оба они могут рассматриваться как нелинейные аналоги метода расщепления Дугласа Рэкфорда. Доказывается сходимость итераций, обсуждаются вопросы реализации методов, приводятся результаты вычислительных экспериментов.
1. Постановка задачи, эквивалентные формулировки
Пусть О С М2 - ограниченная область с кусочно-гладкой границей дО, разделенная кусочно-гладкой кривой Б па две подобласти Ох и О2. Определим множество Кд = {и € Нд (0)| и(х) > 0 для п. в. х € 02} и рассмотрим следующую задачу о препятствии: найти функцию и € Кд, доставляющую минимум функционалу
1{и) = \ J \Чи\2<1х - J /иЛх, / € Ь2(0). (1)
п п
Известно [14], что эта задача имеет единственное решение и € Кд П Н2(О) и эквивалентна вариационному неравенству
и € Кд : J У и • У (у — и,)йх ^ J /(у — и,)йх Уу € Кд. (2)
п п
Приведем еще несколько эквивалентных формулировок задачи (1). которые будут использованы при аппроксимации по методу конечных элементов. Определим следующие пространства и множество:
V = {иг € Нх(0г) | щ(х) = 0 для п. в. х € д0г П дО, }, г =1, 2,
К = {(их, и2) € V х У2 | их(х) = и2(х) для п. в. х € Б, и2(х) ^ 0 для п. в. х € О2 }.
Будем искать пару (и1,и2) € К, доставляющую минимум функционалу
Л,2(«ь«<2) = -Ь^ч) + Мгч) = \ J \\/1ц\23,х- ! /ы.;Дг\ (3)
Нетрудно видеть, что задачи (3) и (1) эквивалентны в следующем смысле. Если (их, и2) - решение задачи (3), то функция и(х) = {и1(х),х € и2(х),х € 02}
и иг
сужения функции и на 0г, г = 1, 2, то (и\,и2) является решением (3).
Функционал дифференцируем по Гато на V х У2, и его дифференциал имеет вид
^12(и1,и2)(У1,У2) = 31 (их)у1 + 3'2(и2)У2, 3'г(иг)Уг = J УигУугйх — J /угё,х.
п, п,
Пусть Ik (u) = {0, если u G K; иначе} — индикаторная функция множества K и dIK - ее субдифференциал. Задача минимизации функционала (3) эквивалентна вариационному неравенству [14]
2 Í
найти (ui, U2) G Vi х V2 : 53 / Vu¿V(vi - ui)dx + ik (vi, V2)-
2 í
- ik(ui,u2) > 53 / f (v¿ - ui)dx V(vi,v2) G Vi х V2, (4)
которое в свою очередь может быть записано в виде включения
Ji(ui) + J2(u2) + BIk(ui,u2) Э 0. (5)
1.1. Аппроксимация. Построим сеточную аппроксимацию задачи (4). используя метод конечных элементов и квадратурные формулы. Для простоты будем считать, что Q = (0,1) х (0,1), Q = (0,0.5)_х (0,1), Q = (0.5,1) х (0,1). Построим в fi i квадратную сетку с шагом ff, а в 2 - квадратную сетку с шагом h = ff/m, где m > 0 - фиксированное целое число. Пусть wi и W2 — множества внутренних узлов этих сеток, - узлы, лежащие на границе области Q, и si = {x G dwi| x G S}. Множество элементов (квадратных ячеек), образованных линиями сетки и лежащих в ^обозначим через •
Определим конечномерные пространства и множество:
Vih = {uih G Vi| u,ih G Qi(x) для всех 5 G Th}, i = 1, 2,
Kh = {(uih,u2h) G Vih х V2h | uih(x) = u2h(x) для x G S, u2h(x) > 0 для x G Q2 },
где Qi(x) - пространство билинейных функций. 3аменив в (4) ui на uih и K на Kh
2 f
найти (uih,u2h) G Vih х V2h : 53 / Vu¿hV(vih - uih)dx + Ikh (vih,v2h)-
i=idi
2 í
- iKh (uih,u2h) > / f (vih - uih)dx V(vih,V2h) G Vih х V2h. (6)
i=io¿
Применим теперь квадратурные формулы для аппроксимации интегралов в (6). Пусть Sá - квадратурная формула трапеций па элементе 5:
i=i
где через {di} обозначены вершины элемента 5 = [xi, xi + h] х [x2, x2 + h]. Поставим в соответствие функциям uih G Vih векторы узловых параметров ui G RNí , i = 1, 2, такие, что uijj = uih(x¿) дая x¿ G wi U si, u2j- = u2h(x¿) дая x¿ G w2 U s2, a Ni равно числу точек множества wi U si. Положим
(A¿ui,vi) = 53 S¿(Vuih(x)Vvih(x)), (fi,vi) =53 Sá(f (x) vih(x)),
«eTÍ
где одним и тем же символом (•, •) обозначено евклидово скалярное произведение как в М^1, так и в М^2. Введем обозначения для множеств:
Mi = {u = (ui,u2) G RN1 x RN | U2,i = Пн(uí)(xi) для xi G S},
где Пн (ui)(xi) - значение в точке xi кусочно-линейной функции, построенной по значениям ui в лотках si,
M2 = {u2 G RN2 |u2,i > 0 для xi G ^2}, M2 = {u = (ui,u2) |u2 G M2}, M = М1ПМ2.
Сеточная схема для задачи (4) имеет вид: найти вектор u = (ui,u2) G RNl x RN2 такой, что для всех v = (vi, v2) G RNl x RN2
(Aiui,vi - ui) + (Á2u2,V2 - u2) + Im (v) - Im (u) > (/i,Vi - ui) + (/2,V2 - u2).
Запишем это вариационное неравенство в виде включения. Пусть
á=(Á /= (/:)■ C=dIм.
Тогда задачу (6) можно записать в виде
Áu + Cu э /. (7)
Матрица Á в (7) является симметричной и положительно определенной, C - максимально монотонный оператор. Отсюда следует однозначная разрешимость (7) [15].
M
равенство IM = IMlnM2 = IMl + Im2 ■ Более того, вектор (ui,u2) с положительными компонентами, удовлетворяющий условию ui,i = u2,i на S, принадлежит intM2 П Mi, поэтому (см. [16]) д(1м1 + 1м2) = д1м1 + д1м2 • Пусть Ci = д1м1, C2 = dIM2 и C2 = д1м[2 — субдифференциалы индикаторных функций соответствующих множеств. Тогда C = Ci + C2, и включение (7) можно написать в виде
Ái Á2)(:)+(¿2)+Ci{uOК/
1.2. Метод расщепления. Применим для решения включения (7) итерационный метод
D-iB(un+i - un) + Áun + C (B(un+i - un) + un) э /, n = 0,1, 2,... (8)
с некоторыми невырожденными матрицами B и D. Реализацию каждой итерации метода (8) будем осуществлять при помощи процедуры расщепления
ÍD-i(un+i/2 - un) + Áun + Cun+i/2 э /,
|b(u"+1 - un) = un+i/2 - un. ' ^
D
D = diag(Ti, ...,Ti,T2,.. .,T2), Ti,T2 > 0.
4-v-' 4-v-'
N1 N2
Выбор матрицы D в таком виде мотивирован тем, что в подобластях Qi имеем сетки с разными шагами, следовательно, каждая из задач в подобластях должна решаться со своим итерационным параметром.
Второе уравнение в (9) при выборе оператора В в блочно-диагональном виде расщепляется па две несвязанные системы линейных уравнений, соответствующие подобластям П1 и П2. Например, в случае В = Е + ВА (нелинейный вариант метода Дугласа Рэкфорда) получаем
' (Е1 + пА^м^1 - <) = иП+1/2 - < вш! и вь
_ (Е2 + Т2А2)(м1+1 - м1) = мЦ+1/2 - в Ш2 и 52,
где Е, — единичные матрицы соответствующих размеров.
Для включения вида (7) с положительно определенной матрицей А сходимость метода Дугласа Рэкфорда доказана в работе [10], где также получены оценки скорости сходимости и показано, что в случае блочно-диагональной матрицы А оптимальные итерационные параметры вычисляются по формуле т = = 1/ ^/А^шахА^щш, где Аг тах и Аг т;п максимальное и минимальное собственные числа матрицы А^.
Остановимся на вопросах решения включения в (9), которое при указанном выборе матрицы В приобретает вид
1 \ / о \ / .ы»+1/2 \ 1 ^l+h-Ayitl
T1 n+1/2 + C иП+1/2 + °Л
-U0 + 1/2 / V C2'"2
\ r2 2 / \ r2
T1 ' (10)
l —u2 + h ~ Aov'i
\ T2
Решение включения (10) для координат вектора мп+1/2 , соответствующих внутренним сеточным точкам подобластей П1 и П2, осуществляется по явным формулам:
1+1/2 = Т1(/1 - А1О+ < ВШ1,
1+1/2 = Рг^с(т2(/2 - А2«П)+ м1) В Ш2,
где Рг^о(х) = шах{х, 0} - операция проектирования на множество неотрицательных чисел.
Для координат вектора м"+1/2, соответствующих сеточным точкам на разрезе получим следующую задачу:
D-1un+1/2 + ö/m, (un+1/2) Э f — Aun + D
Обозначив правую часть включения через Ф = (Ф^Фо). Ф; = —и" + fi — А.;,ы" . i = 1, 2, запишем эквивалентную задачу минимизации:
min g(ui,u2), g(ui,U2) = 77—(Mb«<i) + -¡-(«,2,«,2) - - (Ф2,г/,2). (11)
(«i,«2)eMi 2T1 2t2
Введем следующие нумерации узлов S1 и S2:
S1 = {xfe| xfc = (0.5, (k - 1)H), k = 1,... ,n1},
s2 = {xj 1 = (0.5, (j - 1)h), j = 1 . . .,n2},
где n - число точек множества Sj. Тогда уз ел Xj G s2 с помер ом j = (k — 1)m + 1 будет совпадать с узлом G s 1. По определению множества M1 имеем
u2j = u1k, если j = (k — 1)m + 1,
i / i \ (12)
v>2j = u\ fc+i--1- 14 к 1--, если j = (к — l)m + i + 1, i = 1, .. ., m — 1.
mm
Здесь узлы с номерами (к — 1)т + г + 1, г = 1,... ,т — 1 являются «висячими» узлами множества (т. е. в этих точках нет соответствующих узлов множества ), принадлежащими интервалу ((к — 1)Н,кН) по Х2. Записывая функцию д в (11) более подробно, получим
1 П1 1 П1 1 щ —1 т—1
вЫ,«2) = — XX* + ^1<1(к-1)т+1 +
1 к=1 2 к=1 2 к=1 ¿=1
п1 п1 п1 —1 т—1
— У^ ф1,ки1,к Ф2,(к — 1)т+1 и2,(к —1)т+1 — ^ ^ Ф2,^к) и2,^к); к=1 к=1 к=1 ¿=1
где через (гк) обозначен индекс ] = (к — 1)т + г + 1. Использовав (12), исключим значения и2 и получим задачу безусловной минимизации функции
1 1 п1 2 /(«•!) = д{иии2{и{)) = + — ^ +
п1 —1 т—1 п1
+ + 1 + (1 - «¿)'ы1,й)2 - + Ф2,0-1)т + 1)«<1,й-
2 к=1 ¿=1 к=1
п1—1 т—1
_ ^2,{гк){аги1,к+1 + (1 ~ Я'г)м 1,/с ), аг =
к=1 ¿=1
Уравнение Эйлера для этой задачи система линейных алгебраических уравнений
ди = Ь, (13)
с трехдиагональной матрицей д = {дк,;} € МП1ХП1, элементы которой определены равенствами
( 1 1 \ 1 т—1 1 т—1
Чк,к =--1--Н--У^ (о? + (1 - о.;,)2), Як,к-1 = Чк,к+1 = — Ог(1 - Ог)-
\Т1 Т2 ) Т2 ^ Т2
4 ' ¿=1 ¿=1
Очевидно, что все ненулевые элементы матрицы д положительны и д имеет строгое диагональное преобладание. Для решения (13) может быть использован, например. метод прогонки или метод Холецкого.
2. Метод декомпозиции области с использованием функции Лагранжа
Пусть, как и выше, У^ = {щ € Н 1(^) | и^(х) = 0 для п. в. х € дQ¿ П }, г = 1, 2, и определено выпуклое замкнутое подмножество подпространства У2 :
К2 = {и2 € У | и2(х) > 0 для п. в. х € О2 }.
Для задачи (3) определим функцию Лагранжа
У1 х К2 х Ь2(Б) ^ М, ^?(иьи2,А) = 71 (и1) + + ^ А(и — (14)
я
и перейдем от задачи (3) к задаче поиска содловой точки лагранжиана (14). Если (и1, и2, А)
(и1; и2, А) = вир М у2, у) = М виру2, у),
то и(х) = {и1(х),х € П1;и2(х),х € П2} является решением исходной задачи [16]. Необходимыми и достаточными условиями того, что тройка (м1,м2,Л) есть седло-вая точка лагранжиана (14), являются следующие соотношения [16]:
^ («1,«2,Л)=0, («Ь«2,Л)+ д/к2 («2) Э 0, ^(М1,М2,Л)=0. (15)
В рассматриваемом случае система (15) принимает вид: найти тройку (м1,м2,Л) € € У х К2 х Ь2(£) такую, что
J + J = J € У1,
П1 я П1
У VM2V(V2 - Л(«2 - М2)^Г /2(^2 - У«2 € К2, (16)
я ^2
J(«1 - М2)м^Г = 0 € ).
Докажем существование и единственность решения задачи (16). Как уже отмечалось выше, решение и задачи (1) принадлежит Н2(^). Пусть и = и|п* - сужения функции и та подобласти ^. Из вариационного неравенства (2) с учетом гладкости решения и € Н2(П) следуют поточечные (в смысле «почти всюду») соотношения
-Ди1 = /1, х € О1,
-Дм2 - /2 > 0, М2 > 0, (-Дм2 - /2)и2 =0, х € ^2,
ди1 ди2 дп1 дп2
и1 = 0, и2 = 0,
ui = u2
x G S,
x G \ S, x G dQ2 \ S,
(17)
где ni и П2 — единичные векторы внешних нормалей к Qi и Q2 та границе S. Положив А = du/dn2 G ¿2(S), из поточечных постановок (17) стандартной техникой получаем обобщенные постановки интегральное тождество и вариационное неравенство в (16). Третье уравнение в (16) очевидно выполнено. Таким образом, тройка (u1,u2, А) - решение задачи (16). Из интегрального тождества и вариационного неравенства в (16) легко получить единственность ui и U2 и, как следствие, единственность А = du/dn2.
2.1. Аппроксимация. Так же, как в п. 1.1, построим в области Q сетку, определим конечномерные пространства
Vih = (ujh G V¿| ujh G Qi(x) для bcex 5 G Тк'}, Л^ = trV2h — пространство следов функций из V2h на S
и множество
K2h = {u2h G V2h | U2h(x) > 0 для x G ^2 }.
Аппроксимацией задачи (16) по методу конечных элементов назовем следующую задачу: найти тройку (и1^,и2^,Л^) € Ул х х Л^ такую, что для любых
(ъте Ущ х К2н х Ль
j + j \hVmdF = j
П1 Я П1
У Vu2hV(v2h - и-2н)йхЛь- и2ь)йГ ^У /- и-2н)йх, (18)
Я ^2
У (и1ь - и2н)^нйГ = 0.
Я
Для аппроксимации интегралов по подобластям г = 1, 2, снова будем использовать составную формулу трапеций по соответствующим квадратным элементам. Далее, обозначим через = {а} разбиение Б та отрезки длины Н, где Н - шаг мелкой сетки, и для аппроксимации интегралов по отрезкам а е Б применим формулу трапеций:
[ Н
/ г>скс « Еа(г>) = ~(г>(х2) + г>(х2 + /г)), а = [ж2, х2 + /г].
а
Поставим в соответствие функциям и^ е У^ векторы узловых параметров щ е М"*, г = 1, 2, где N равно числу точек множества и в^, в соответствие Ль — вектор Л е М"А , N равно числу точек в2, и положим
(А и^) = ^ Б<5 (Ущь(х)Чъ^ (х)), (/г,Ъг) = ^ Бй(/(х)ъ^(х)),
¿eT¿ йeT¿
(ДД^ )=^ Еа (Ль,^), г =1, 2. (19)
аеЗн
Как н в п. 1.1 определим множество
м2 = {и2 е М"21 и2 ^ > 0 для х е ^2} (20)
и максимально монотонный оператор С2 = .
Сеточной аппроксимацией (16) будет следующая система
А1и1 + Д1Л = /1,
А2и2 + С2и2 + Д2Л Э ¡2, (21)
ДТ и1 + ДТ и2 = 0.
Приведем также запись этой системы в терминах сеточных функций и сеточных операторов. Будем использовать обозначения Дни = их1г1 + их2г2, Дьи = их1г1 + +иЖ2г2 для сеточных операторов Лапласа на сетках с шагами Н в и Н в П2 , в тот время как для разностных производных ихi, иг¿ обозначения одинаковые для разных сеток. Тогда сеточная задача запишется в виде
-Дни1 = /1, х е
-^гЧх! ~ и 1х2х2 = + /ь -г' €
и1 = 0, х е 3^1 \ в1,
m— 1
где As(xj) = I ) (A(xi - jh) + A(x* + jh)) (1 - j/m) | /m,
j=i
+ <^2^2 Э /2, x G W2,
2 2
J/1-'- ~ u2x2x2 = ^A + /2, ж G s2,
U2 = 0,
x G \ S2,
(23)
U2j = uifc,
j = (k — 1)m + 1,
i I i = 'Mi,fc+1---h «<i,fc ( 1 — —— ) , 3 = {k - l)m + г + 1, i = 1,
.., m — 1.
(24)
Пусть
u =(Mi,M2)T, / =(/i,/2)T,
R
C =
00
A =
Ai 0
Ri
R2 У' ° l 0 (72 У ' л V 0 A2 /,u G ,Nn = Ni + N2; A,C : ^ ; R : RNi ^ .
Тогда систему (21) можно записать следующим образом:
f Au + Cu + RA э /,
—RT г
0.
Наконец, введя обозначения
A--
y =(u,A)T, F = (/, 0)T, C =
AR —RT 0
C0 00
получим
Ay + Cy Э F,
A
3a > 0 : (Ay, y) = (Au, u) > a||u||2. Теорема 1. Решение задачи (27) существует и единственно. Доказательство. Рассмотрим задачу
min, hAu,u) - (f,u), пемinм2 2
(25)
(26) (27)
(28)
где А = Ат > 0, А : ^ М1 = {и € | Дти = 0} с прямоугольной
матрицей Д, матрицы А и Д определены в (25), М2 определено в (20). Множество М1 П М2 те пусто, так как вектор и = 0 принадлежит этому множеству. Задача (28) имеет единственное решение, так как функционал квадратичный, а множество М1 П М2 — выпуклое, замкнутое и непустое.
и
щими узлам сетки о>2 и удовлетворяющими (24). Такой вектор принадлежит множеству ^М2 П М15 поэтому д(/л^1 + /м2) = + д/мд2. Следовательно, задача (28) равносильна следующей
Au + d/м (u) + ö/м (u) э /.
(29)
Иначе говоря, если и е М1 П М2 - решение (28), то
Зх е д1м2(и), З£ е д/^ (и) : Аи + х + С = /Последнее означает, что
-(Аи + х - /) = С е д1Й1 (и) ^ (Аи + х - /,у) = 0 Уу е Й1.
Множество М1 то определению является ядром оператора ЯТ. Таким образом, (Аи + х — /) е ( КегЯТ= 1т Я [17], и задача ЯХ = / - х - Аи имеет решение. Это означает, что пара (и, Х) является решением системы
Аи + х + ЯХ = /, х е Си,
Т (3°)
-ЯТ и = 0.
и
казать единственность Х. Пусть существуют Х1 и Х2, удовлетворяющие (30).
и
внутренних точках П2, то х® = 0 для г, соответствующих узлам в2. Подставив Х1 Х2
Я(Х1 - Х2) = 0 Я
диагональный блок размера N х , на диагонали которого стоят положительные Я
Х1 - Х2 = 0 . □
2.2. Метод расщепления. Запишем схему расщепления для задачи (27): ' Б-1 (Е + Бв)уп+1/2 э Б-1уп + Д - Ау
В(уп+1 - уп) = уп+1/2 - уп,
(31)
Здесь
Е1 + т1А1 0 -т1Я1 В = Е + Б А = | 0 Е2 + Т2А2 -Т2Я2 |. (32)
-тзЯТ -тзЯТ Ез
Б = diag(тl,... ,Т1,Т2,.. -,т2,тз, ...,тз), т® > 0, г = 1, 2, 3,
4-V-' 4-V-' 4-V-'
N1 N2 №а
является диагональной матрицей итерационных параметров, такой, что каждая из и1, и2 Х
Е®, г =1, 2, 3, — единичные матрицы соответствующих размеров. Первое включение в (31) распадается на следующие подзадачи:
хП+1/2 = т1(/1 - А1< - ДТХп)+ <,
хП+1/2 = Рг^с(т2(/2 - А2ип - ДТХп) + ип), (33)
,Хп+1/2 = тз(ЯТ ип + ДТ ип) + Хп,
решение которых осуществляется по явным формулам.
Второе уравнение, в отличие от случая п. 1.2. представляет собой связанную через Л"+1 систему в точках S
'(Е1 + Т1А1)(ип+1 - и") + Т1Д1(Л"+1 - Л") = и"+1/2 - и?,
(Е2 + Т2А2Ж+2 - и") + Т2Д2(Л"+1 - Л") = и"+1/2 - и?,
кЛ"+1 - Л" - тзДТ(и"+1 - и?) - тзДТ(и"+1 - и?) = Л"+1/2 - Л".
Приведем формулировку утверждений из [18]. которые будут использованы при доказательстве сходимости метода (31) (32).
Напомним, что оператор Т : К ^ К, где К С V — выпуклое, замкнутое множество гильбертова пространства V, называется нерастпягивающим, если
||Тх - Ту У < Ух - у|| У х,у € К
и жестко нерастягивающим. если
(Тх - Ту, х - у) > ||Тх - Ту||2 У х, у € К.
Предложение 1 [18]. Пусть Т1,Т2 : К ^ К - жестко нерастягивающие операторы, тогда Б = Т1(2Т2 - Е) + Е - Т2 - жестко нерастягивающий.
Предложение 2 [18]. Пусть оператор О : К ^ К - жестко нерастягивающий, и существует по крайней мере одна его неподвижная точка. Тогда итерации метода «"+1 = О«" слабо сходятся к V, где V = О« - неподвижная точка опера-О
Теорема 2. Итерации .метода (31) (32) сходятся к решению задачи (27) при любом выборе матрицы В = Вт > 0.
Доказательство. Введем обозначения
Ту = Су - Д, = (Е + ВА)-1, у = (Е + ВТ)-1у, V" = у". Тогда включение (27) примет вид
Ау + Ту э 0,
а итерационный метод (31) можно записать как
«"+1 = (27а - Е)«" + (Е - (34)
Пусть пространство векторов оснащено скалярным произведением (и, V) 2-1 = = (В-1и,«), где (•, •) - евклидово скалярное произведение. Обозначим через Н пространство векторов с введенным скалярным произведением и докажем, что оператор О = (27а - Е) + Е - является жестко нерастягнвающим в Н, т. е. для любых и, V € К" выполняется неравенство
(Ои - Gv, и - v)D-l ^ ||Ои - Оv||D-l.
Согласно предложению 1, для этого достаточно доказать, что и - жестко Н
Обозначив 7ам = х и V = у, получим
(7ам - и - V)2-1 = (х - у, (Е + ВА)х - (Е + ВА)у)2-1 =
= ||х - уНЪ-1 + (Ах - Ау,х - у)д-1 > ||х - у||2-1 = ||7ам - ^^,
т. е. За - жестко нерастягивающий. Аналогично.
(Зти — ЗтV, и — у)д-1 = (х — у, (Е + ПС)х — (Е + ПС)у)о-1 ^
=11 "Г" — о
> IIх — у\\0-1 = \\Зти — Зт•
Таким образом, О - жестко нерастягивающий в Н, и из предложения 2 следует, что уп ^ V при п ^ <х, где V = Оу - единственная неподвижная точка оператора О. Это равносильно существованию единственного решения задачи (27) у = (Е + ЛА)"Ч>. ' □
3. Численные результаты
Для сравнения предложенных методов решения двух полученных конечномерных задач (7) и (27) была проведена серия числовых расчетов. Методы тестировались па примере, когда известно точное решение иех дифференциальной задачи (1):
Гвт( ^)вт(лх2), х € (0, 0.75) х (0, 1), иех (х) = < V 3 /
[о, х е (0.75, 1) х (0, 1).
При расчетах шаг Н сетки в области ^ варьировался, в области П2 использовалась сетка с шагом Н = Н/ 2. Параметры т1 и т2 были взяты как теоретически оптимальные для задач в подобластях и ^2. В методе с множителями Лагран-жа параметр тз подбирался в ходе вычислительных экспериментов. Начальное приближение для А во всех тестах А0 = 0.
Табл. 1
Число итераций. Критерий остановки ||и — иех||с < £
£ И х 21. 21 х 41 21 х 41. 41 х 81 41 х 81. 81 х 161
I II I II I II
0.5 5 6 9 9 15 16
0.05 17 17 34 34 67 67
0.005 29 29 59 59 119 119
0.0005 - - 86 83 171 170
В табл. 1 показано число итераций метода Дугласа Рэкфорда для схемы без множителей Лагранжа (I) и с множителями (II) в зависимости от заданной точности на различных сетках. Начальное приближение и = 0. Знак « — » означает, что заданная точность не достигается.
Далее исследовалось поведение методов в зависимости от шага сетки и выбора начального приближения. В качестве критерия остановки использовалось условие ||и — иех||с < 0.005. В первой строке табл. 2 показаны варианты выбора начального
и1 = 0
ш 1 и в1, и2 = 0 в ш2 и в2 число итераций для обеих схем практически одинаково и при увеличении размерности задачи линейно зависит от числа точек по одному направлению.
Во втором случае (и1 = — 2 в ш1 и в1, и2 = 2 в ш2 и в2) число итераций для схем I и II одинаково, и оно значительно больше, чем в первом случае. В этом случае в процессе вычислений возникает локальное возмущение погрешности максимум погрешности сохраняется в течение итераций вблизи угловых точек области ^2 и медленно убывает в них, в то время как в остальных точках области
Табл. 2
Число итераций при различных начальных приближениях
иЧ = 0, и'о = 0 иЧ = -2, = 2 иЧ = 2, г^ = -2
I II I II I II
И х 21. 21 х 41 29 33 115 115 115 37
21 х 41. 41 х 81 59 63 224 224 232 66
41 х 81. 81 х 161 119 122 444 444 127 129
81 х 161, 161 х 321 238 241 885 885 979 251
решение уже фактически не меняется (рис. 3. 4). Отметим, что в случае схемы I возмущение погрешности наблюдается также около граничных точек, прилегающих к разрезу (рис. 3). В результате дополнительных расчетов, при варьировании начального приближения, было выяснено, что особенность в погрешности вызвана рассогласованием начального приближения с граничными условиями, которое особенно сказывается в окрестности коинцидентного множества.
В третьем случае выбор начального приближения = —2 не приводил к возмущению погрешности около угловых точек, что обеспечило меньшее число итераций в схеме II по сравнению с предыдущим случаем. Однако при решении схемы I снова наблюдался пик погрешности около разреза что привело к большому числу итераций. Можно было бы предположить, что это происходит из-за разрыва начального приближения на границе но выбор начального приближения в виде
п\ = —2, u2 = —2 не исправил ситуацию.
Основываясь на результатах вычислительных экспериментов можно утверждать. что в случае «хорошего» начального приближения метод Дугласа Рэкфорда для обеих схем ведет себя одинаково, но схема с множителями Лагранжа схема менее чувствительна к выбору начального приближения. Однако использование схемы Лагранжа осложняется отсутствием теоретических оценок параметра Т3 •
Ниже приведены некоторые иллюстрации к результатам расчетов на сетке с шагом H = 0.05 и h = 0.025. На рис. 1 показано численное решение сеточной задачи. Рис. 3. 4 иллюстрируют распределение погрешности в схемах I и II при остановке по критерию ||u — uex||c ^ 0.005 в случае начального приближения ui = = —2, u2 = 2.
В любом случае после достаточно большого количества итераций получено одинаковое распределение погрешности, как на рис. 2. что подтверждает совпадение сеточных решений обеих схем.
Работа выполнена при финансовой поддержке РФФИ (проект Х- 04-01-00484).
Summary
M.A. Ignatieva, A.V. Lapin. Solution of the obstacle problem by domain decomposition method.
Domain decomposition method with non-overlapping subdomains is applied for solving the obstacle problem. An obstacle being located in a known subdomain of the initial domain, partitioning of the domain is made by using this information, so, the information about possible lost of the solution regularity.
Finite element method with quadrature rules and non-matching grids in the subdomains is used to approximate corresponding variational inequality; finer grid is constructed in the subdomain containing the obstacle.
Two iterative methods are constructed to solve finite dimensional problems, they can be viewed as non-linear variants of Douglas - Rachford splitting iterative method. Convergence of the iterative algorithms is proved, their implementation is discussed, and numerical results are applied.
Литература
1. Лебедев В.И. Метод композиции. М.: ОВМ АН СССР, 1986. 192 с.
2. Лебедев В.И., Агогиков В.И. Операторы Пуанкаре Стеклова и их приложения в анализе. М.: ОВМ АН СССР, 1983. 184 с.
3. Агогиков В.И. Методы разделения области в задачах математической физики // Вычислительные процессы и системы. М.: Наука, 1991. Вып. 8. С. 4 51.
4. Quarteroni A. Domain decomposition and parallel processing for the numerical solution of partial differential equations // Survey on Mathematics for Industry. 1991. V. 1. P. 75 118.
5. Quarteroni A., Periaux J., Kuznetsov Yu., Widlund O.B. Domain decomposition methods in science and engineering. Providence: AMS, 1994.
6. Smith В., Bj0rstad P., Gropp W. Domain decomposition. Parallel multilevel method for elliptic partial differential equations. Cambridge University Press, 1996.
7. Le Tallec P. Domain decomposition methods in computational mechanics // Comput.. Mechanics Adv. 1994. V. 1. P. 121 220.
8. Herrera I., К eyes D., Witllund О., Yates R. (Eds.) Domain Decomposition Methods in Science and Engineering. Mexico City, Mexico: National Autonomous University of Mexico (UNAM), 2003. 466 p.
9. Kornhuber R., Huppe R., Periaux J., Pirunneau O., Witllund O., Xu J. (Eds.) Domain decomposition methods in science and engineering. Lecture Notes in Computational Science and Engineering. Springer, 2004. V. 40. 700 p.
10. Laitinen E., Lapin A. V., Pieskii J. Splitting iterative methods and parallel solution of variational inequalities//Lobaclievskii J. of Mathematics. 2001. V. 8. P. 167 184.
11. Dostal Z., Frietllantler A., Gomes F.A.M., Santos S.A. Preconditioning by projectors in the solution of contact problems: A parallel implementation // Ann. Oper. Res. 2002. V. 117. P. 117 129.
12. Dostal Z., Htrrak D. Scalability and FETI based algorithm for large discretized variational inequalities // Math. Comput. Simul. 2003. V. 61. P. 347 357.
13. Danek J., Hlavaeek I., Netloma J. Domain decomposition for generalized unilateral semi-coercive contact problem with given friction in elasticity // Math. Comput. Simul. 2005. V. 68. P. 271 300.
14. Лионе Ж.-Л. Некоторые методы решения нелинейных краевых задач. М.: Мир, 1972. 588 с.
15. Barbu V. Nonlinear semigroups and differential equations in Banacli spaces. Nordlioff Intern. Publ., 1976.
16. Эклаид И., Темам P. Выпуклый анализ и вариационные проблемы. М.: Мир, 1979. 399 с.
17. Воеводин В.В., Кузнецов Ю.А. Матрицы и вычисления. М.: Наука, 1984. 320 с.
18. Lions P. L., Mereier В. Splitting algorithms for the sum of two nonlinear operators // SIAM J. Kurner. Anal. 1979. V. 16. P. 964 979.
Поступила в редакцию 04.11.05
Лапин Александр Васильевич доктор физико-математических паук, профессор кафедры математической статистики Казанского государственного университета.
Е-шаП: alapinQksu.ru
Игнатьева Марина Александровна кандидат физико-математических паук, научный сотрудник НИИ математики и механики им. Н.Г. Чеботарева Казанского государственного университета.
Е-шаП: Marina.IgnatievaQksu.ru
