Многосеточные итерационные алгоритмы в методе конечных элементов с учетом численного интегрирования Текст научной статьи по специальности «Математика»

Вычислительные технологии

Том 14, № 1, 2009

Многосеточные итерационные алгоритмы в методе конечных элементов

*

с учетом численного интегрирования

Л. В. Гилева

Институт вычислительного моделирования СО РАН, Красноярск, Россия

е-таП^11еуа@1ст. кгаэп. ги

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

Ключевые слова: метод конечных элементов, многосеточный итерационный алгоритм, каскадный алгоритм, квадратурная формула, оценка числа операций.

Введение

Наиболее простая версия многосеточных методов — каскадный алгоритм. Он состоит в применении традиционных итерационных методов на последовательности вложенных сеток для достижения хорошего начального приближения на каждой сетке за счет интерполяции приближенного решения с предыдущей, более грубой, сетки. Эффективность такого подхода изучалась в работе В.П. Ильина и В.М. Свешникова [1], однако систематические исследования на эту тему начались только в 90-х годах. Свое название метод получил в статье П. Дойфльхарда [2], где продемонстрирована вычислительная эффективность алгоритма. Теоретическое обоснование дано в работах [3, 4], а затем повторено в [5].

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

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

© Институт вычислительных технологий Сибирского отделения Российской академии наук, 2009.

гладкости. В трехмерных краевых задачах алгоритм оказался даже более эффективным [9, 10].

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

Для построения дискретных задач на последовательности сеток используется следующий подход. Вначале строится система уравнений метода Бубнова—Галерина с использованием кусочно-линейных элементов на треугольниках на самой мелкой сетке. При этом для вычисления элементов матрицы и вектора правой части применяется квадратурная формула, точная для многочленов степени < 1. Далее элементы матрицы системы уравнений на более редкой сетке вычисляются как линейные комбинации элементов матрицы, соответствующей более мелкой сетке. Другими словами, матрицы систем уравнений на соседних сетках связаны с помощью операторов интерполяции и проектирования, как и в случае, когда элементы матриц вычисляются точно. Это свойство дает некоторые удобства при реализации метода на ЭВМ.

В данной работе рассматриваются каскадный алгоритм и полный многосеточный алгоритм на основе несимметричного Ш-цикла. Доказана оптимальная вычислительная сложность обоих алгоритмов, т.е. обоснована оценка Б < сЫ, где Б — число арифметических операций для достижения порядка точности, совпадающего с порядком точности решения дискретной задачи с учетом численного интегрирования, N — число неизвестных системы Бубнова—Галеркина, с — константа, не зависящая от N.

1. Формулировка дифференциальной задачи

В выпуклом открытом многоугольнике П С Я2 с границей Г рассмотрим задачу Дирихле

2

— Е дг(а.д.и) + аи = f в П, (1)

¿,.7 = 1

и = 0 на Г, (2)

где коэффициенты и правая часть уравнения (1) удовлетворяют условиям гладкости и эллиптичности:

дга. € Ья (П), д> 2, г] = 1, 2; а12 = а21 на П;

2 2 2 0 Е 62 < Е а. < на П V 6 € Я, (3)

г=1 ¿,7=1 ¿=1

V > 0 > 0; а, f € Ь2(П); а > 0 на П.

Тогда задача (1)-(2) имеет единственное решение в Ш|(П), удовлетворяющее оценке [11]:

||и||2,п < с^|о,п.

О

Перейдем к обобщенной формулировке для задачи (1)-(2): найти и € (П), такую, что

О

С(и,ь) = f (у) V V € Ш} (П), (4)

где билинейная форма С и функционал f соответственно определяются соотношениями

( 2 \

С(п,у) = М а^д^пдгУ + апь I ¿х, (5)

П \гл=1 /

f (V) = I Мх. (6)

П

При условиях (3) задача (4) имеет также единственное решение. 2. Формулировка дискретной задачи

Для построения схемы метода конечных элементов проведем разбиение области П. Вначале разобьем ее на небольшое количество треугольников так, чтобы получившаяся триангуляция Т0 была согласованной, т. е. любые два треугольника должны иметь общую сторону либо общую вершину или не иметь общих точек. Обозначим через Л0 максимальную из длин сторон полученных треугольников. Для г = 1,... ,1 положим Ыг = 2г, Ы = Ло/Ы и разобьем каждый треугольник на Ыг2 равных треугольников.

Обозначим множество всех вершин полученной триангуляции Т через Пг и введем Пг = Пг Р| П. Через пг обозначим число точек множества Пг. Каждому узлу у £ Пг

о

поставим в соответствие базисную функцию ^у £ (П), которая линейна на каждом треугольнике из Т, равна единице в узле у и нулю во всех остальных узлах из Пг. Обозначим через Нг линейную оболочку этих функций, т.е. Нг = врап{^у}, у € Пг.

о

Рассматривая (4) на подпространстве Нг С (П), получим дискретную задачу: найти Vг £ Нг такую, что

С(^г,^) = f (V) V V £ Нг. (7)

Обозначим через Мг пространство размерности пг, состоящее из векторов w с компонентами и>(х), х £ Пг. Тогда задача (7) эквивалентна системе линейных алгебраических уравнений

Ь V = {г, (8)

где vг £ Мг — вектор неизвестных с компонентами vг(y), у £ Пг; ^г £ Мг — вектор правой части с компонентами /¿(х) = f (<^Х), х £ Пг; Ьг — матрица размером пг х пг с элементами Ьг(х,у) = ), х,у £ Пг.

Произвольному вектору V £ Мг поставим в соответствие восполнение в Нг:

= ^ v(y)pу(х), х £ П. (9)

уеП;

Из определения базисных функций следует, что

v(y) = ЭД, у £ Пг. (10)

Таким образом, мы определили изоморфизм между векторами V £ Мг и функциями V £ Нг.

Введем энергетическую норму для функций

||г|п = С(г,г)1/2, V £ Ж21(П).

Для функций из (П) она эквивалентна норме || • ||1;п [12]:

О

С2 11V11,п < ||у||п < сЗ|у|1,П, V € Ш1 (П). (11)

Введем также скалярное произведение и нормы для векторов:

(V, w)i =^2 У(х)и(х), ||у|г = (V, V)1/2, ¡V||г = (\,Ьг\)1/2, V, w € Мг.

На основании (9) и (10) для изоморфной пары V € Мг, V € Нг имеем

II V || г = || г || п. (12)

Норма |г|о,п эквивалентна норме ||V|г с множителем кг [13]:

С4^г | V |г < |г|о,П < С5 Л-г 11V 11 г. (13)

Подпространства Нг обладают свойством вложенности Нг-1 С Нг, т.е. любая базисная функция у1-1 € Нг-1 является линейной комбинацией базисных функций ^>гу € Нг:

^ = Е вху Ч>у, х € Пг-1. (14)

уеП;

Введем оператор интерполяции 1г-1 : Мг-1 ^ Мг. Пусть V — произвольный вектор из Мг-1. Тогда компоненты вектора w = 1г-\У € Мг определяются по формуле

и(У) = X у(х)в*У, У € Пг, (15)

с коэффициентами вху из (14). Отметим, что восполнения векторов V и w совпадают, т.е. V = и). Таким образом, оператор 1г-1 соответствует тождественному оператору на подпространстве Нг-1.

Введем также оператор проектирования Яг : Мг ^ Мг-1. Для вектора w € Мг компоненты вектора V = Я^ € Мг-1 вычисляются по формуле

у(х) = ^ вху и(у), X € Пг-1. (16)

уеП;

Из (15) и (16) следует, что (ЯгW, v)г-l = ^,/г-^)г, V € Мг-1, w € Мг, т.е. Яг = /*-1.

При условиях (3) задача (7) имеет единственное решение, которое удовлетворяет оценкам [12, 13]

IIи - Уг||п < с6кг^|о,п, (17)

|и - Уг|о,п < с7к2|f |о,п. (18)

3. Использование численного интегрирования

На практике, даже если функции а?, а и / имеют простые аналитические выражения, элементы матрицы Ь» и вектора ^ в (8) редко вычисляются точно. Вместо этого они аппроксимируются с помощью численного интегрирования.

Мы начнем с построения дискретной задачи с учетом численного интегрирования на самой мелкой сетке. Используя (5) и (6), можем записать

¿г(ж, у) = ^) = Е Л Е а?<9?^Д^ + а^У ) (19)

ш€Т ш У

fг(ж) = f (у£) = Е / М^ж, ж, у € Пг. (20)

шеТ ш

Пусть ш — произвольный треугольник из Т. Квадратурная схема состоит в замене

м

интеграла / д(ж)^ж конечной суммой вида У] ат,ш д(Ьт,ш), где ат,ш — веса, а точки

ш т=1

— узлы квадратурной формулы. Будем использовать следующую квадратурную формулу:

/д(ж)^ж - д(Ьт,ш), (21)

т=1

где твав(ш) — площадь треугольника ш, а Ьт,ш, т = 1, 2, 3, — вершины ш. Формула (21) точна для многочленов степени < 1.

Введем аппроксимирующую билинейную форму

I \ 3 2

- теаз(ш)

£ (и, V) = -3-¿^ (а?д?идг^ + амг>)(от,ш) (22)

ш€Т т=1 =1

и аппроксимирующий функционал

f (V) = Е -^ Е^Х^ш). (23)

ш€Т т=1

Подставляя (21) в (19) и (20), получим вместо матрицы ¿г матрицу ¿г с элементами ¿г(ж, у) = £ , ^У), ж, у € Пг, а вместо вектора правой части íг вектор ^ с компонентами Л (ж) = f (<^Х), ж € Пг. В итоге мы получаем систему уравнений ¿г ^г = , которая эквивалентна следующей задаче: найти гг € И1 такую, что

£ (гг, V) = /(V) V V € Иг,

где билинейная форма £ и функционал / соответственно заданы равенствами (22) и (23).

Перейдем к построению систем уравнений на более крупных сетках. Матрицы систем уравнений Бубнова—Галеркина на соседних сетках связаны соотношением [13] = ЛгЬ»/г-1, что позволяет по матрице ¿г, соответствующей самой мелкой сетке, строить матрицы г < /, с помощью простых операций суммирования строк и столбцов. Мы

потребуем, чтобы аналогичное свойство выполнялось для матриц Ьг. Предположим, что на сетке Пг мы имеем систему уравнений

^г = Ъ, (24)

где Ьг — матрица размером пг х пг с элементами Ьг(х,у) = ), х,у € Пг, а

fг € Мг — вектор с компонентами fг(x) = f (^Х), X € Пг. Положим

Lг-1 = Яг Lгh-1, (25)

f г—1 = Яг (26)

Систему уравнений на сетке Пг-1 определим следующим образом:

Ьг-1 Wг-l = f г-1. (27)

Используя равенство (26), определение оператора проектирования (16), линейность f и соотношение (14), для компонент вектора правой части (27) получим

/г-1(х) = /(^Т1), X € Пг-1. (28)

Используя (14), (22), определения операторов интерполяции и проектирования, равенство (25) и учитывая вид элементов матрицы Ьг, получим, что элементы матрицы Ьг-1 определяются соотношением

¿г-1(х,у) = £(^Г1,^-1), х,у € Пг-1. (29)

Из (28) и (29) следует, что для любого г = 0,... ,1 система уравнений (24) эквивалентна следующей задаче: найти иг г такую, что

£(и)г,у) = 7 (у) V V € Нг, (30)

где билинейная форма £ и функционал f соответственно заданы соотношениями (22) и (23).

Отметим, что задача (30) может быть получена из (7) с использованием составной квадратурной формулы:

д(х)йх « ^ -ъТ. 9(Ъш,Ш1), (31)

,, т=1

где шг — произвольный треугольник триангуляции Тг. Очевидно, что квадратурная схема (31) точна для многочленов степени < 1 .

Так как квадратурные формулы (21) и (31) точны для многочленов степени 0, то решение задачи (30) удовлетворяет оценке [12]

||и - и)г11,п < с8кг. (32)

Кроме того, поскольку квадратурные формулы (21) и (31) точны и для многочленов степени 1, имеет место оценка в норме Ь2 [12]

|и - и)г|о,п < с9Н2. (33)

В квадратурной формуле (21) все веса строго положительны, и она точна для многочленов степени 0. Поэтому существуют константы ci0 и Сц, не зависящие от hi, такие, что [14]

L (v,v) > Cio|M|?>n, (34)

L (v,w) < cii||v||i;n ||w||i;n, v,w G Hl. (35)

В силу вложенности подпространств H% оценки (34) и (35) справедливы для любых v,w G H\ Введем норму для функций из Hг:

[v]n = L(v,v)i/2, v G Hi.

На основании (34) и (35) она эквивалентна норме || • ||i;n:

Ci2 ||v|| i;Q < [v]n < Ci3 ||v|i>n, v G H1. (36)

Отсюда с учетом (11) следует, что

См|М|п < [v]n < Ci5 || v ||n, v G H \ (37)

т.е. для функций из H% нормы || • ||q и [^]q эквивалентны.

На основании (31) матрица L i положительно определена, поэтому мы можем ввести норму для векторов

[v]i = (v,LiV)i/2, v G Mi. Используя (9) и (10), легко показать, что для изоморфной пары v G Mi, v G Hi

[v]i = [v]n. (38)

Из (37), (12) и (38) следует эквивалентность норм || • ||i и [• ]i, т.е.

ci41|v || i < [v]i < ci51|v || i, v G Mi. (39)

Нам потребуется оценка собственных чисел матрицы Li. Лемма 1. Собственные числа матрицы Li удовлетворяют неравенству

0 < Ai < cie, i = 0,...,l. (40)

Доказательство. Для собственных чисел матрицы Li мы имеем оценку [13]

0 < Ai < ci7, i = 0,...,l. (41)

Обозначим через A* и A* максимальные собственные числа матриц Li и Li соответственно. На основании соотношения Рэлея и оценок (39) и (41) можем записать неравенства

i* (v,L iv)i [v]2 A* = sup —---— = sup --— <

vEMi, (v, v)i vEMi, (v, v)i

v=0 v=0

^ С1б |v|2 2 (v,Lj v)i 2 л* ^ 2 < sup —-— = С15 SUp —--— = С15Лг < C15C17,

vEMi, (V, v)j v6Mi, (V, V)

v=0 v=0

т. е. верхняя оценка в (40) справедлива с константой c16 = c25c17. Нижняя оценка следует из положительной определенности матрицы Lj. □

Итак, на последовательности сеток Q мы получили последовательность задач: для заданного fj G Mj найти вектор wj G Mj такой, что

LjWj = fj. (42)

Наша цель состоит в решении задачи (42) при i = /. Для этого мы используем два варианта многосеточных методов (каскадный алгоритм и многосеточный алгоритм с использованием несимметричного W-цикла) и докажем их сходимость.

4. Формулировка каскадного алгоритма

Каскадный алгоритм состоит в последовательном решении задач (42). При i = 0 число уравнений в (42) невелико и задача решается прямым методом. На более мелких сетках приближенные решения получают итерационным методом. Сформулируем каскадный алгоритм с некоторым абстрактным итерационным процессом Sj (сглаживателем). Каскадный алгоритм:

1) uo = L-1f о;

2) для i = 1,2,...,/ цикл { (43)

2.1. Zj = ij_1Uj_1;

2.2. полагаем uj = Sj(Lj, zj, fj); }.

В качестве сглаживателя рассмотрим два итерационных процесса.

Метод сопряженных градиентов (mj итераций); процедура Sj(Lj, zj,fj):

3) yo = Zj; po = ro = f j - Ljyo; ao = (ro, ro)j;

4) для i = 1, 2,..., mj цикл {

если = 0, то ymi = y^_1 и перейти на 5; afc-1 = afc-1 /(Pfc_b LjPfc_1)j;

yfc = yfc_1 + afc_1pfc_1; (44)

rfc = rfc_1 — afc_1LjPfc_1;

afc = (rfc, rfc)j; ^fc = afc/afc_1;

Pfc = rfc + ^fc Pfc_1;};

5) полагаем Sj = ym.

Метод простых итераций (mj итераций); процедура Sj(Lj, zj,fj):

3) yo = Zj;

4) для k = 1, 2,... , mj цикл {

rfc_1 = (Л*)_1 cos_2 П(2к + 1) ; (45)

k У j> 2(2mj + 1)' V '

yfc = yfc_1 - Tfc_1(Ljyfc_1 - fj)j };

5) полагаем Sj = ymi.

Здесь Л* — верхняя оценка собственных чисел Л оператора Li в пространстве Mi, т.е. L№ = В (45) она нужна в явной форме, поэтому она находится из леммы Гершго-рина [15] и удовлетворяет неравенству

Л* = max Л < Л* < c18 max Л = c18Л*. (46)

Aes^Li)

В (44) в явном виде она не используется, поэтому в теоретических выкладках ее можно положить равной Л*, т.е. (46) тогда выполняется с константой c18 = 1.

Зафиксируем целое i G [1, /] и определим линейный оператор Bi : Mi ^ Mi, ставящий в соответствие погрешности начального приближения = wi — Z для решения задачи (42) погрешность конечной аппроксимации = wi — ui, полученной на i-м уровне каскадного алгоритма (43) с помощью метода сопряженных градиентов (44), т.е. = Bi бо.

Введем норму для оператора B : Mi ^ Mi:

[Bli = sup |Bv|i

v6 Mi, [V]

v=0

Оператор Bi является матричным полиномом от Li, т.е.

mi

Bi = I + X Yk L k. (47)

fc=i

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

Аналогично определим оператор Qi : Ыi ^ Ыi для метода простых итераций (45), т.е. = Qiд0. Оператор Qi является матричным полиномом и имеет вид

mi-1

Qi =Ц(1 — тк L i). (48)

k=0

Для него выполняются следующие оценки [3]:

iw]i < v i llwHi, (49)

i 2mi + 111 Hi' v '

iw]i < [w]i V w G Mi. (50)

5. Сходимость каскадного алгоритма

Теорема 1. Пусть оператор ^ задачи (42) является самосопряженным, положительно определенным и удовлетворяет условию (25). Тогда для приближенного решения и^ полученного на 1-м уровне каскадного алгоритма (43) с одним из итерационных процессов (44) или (45), справедлива оценка

. к

^- uili £ 2ттг■ (51)

3=1 ■'

Доказательство. Обозначим через ei погрешность на г-м уровне каскадного алгоритма, т. е.

ei = wi — ui, г = 0,1,... ,1.

Для метода сопряженных градиентов (44) имеем ei = Bi(wi — zi), а для метода простых итераций (45) — ei = Qi(wi — zi). Оператор Qi имеет вид (47), поэтому на основании приведенного выше свойства оператора В. получаем

[ЫЬ < [Qi(Wi — Zi)]i.

На основании (43) можем записать

N1. < [Qi(wi — )] + [Qi/i_l(wi_l — Ui_l)]i =

= [Qi(wi — Wi_l)]i + [Qi/i_lei_l]i. (52)

Оценим первое слагаемое в правой части (52). Для этого воспользуемся (49), (13), (40) и (33):

' \ * . / \ *

_1 V Л i _1

¿(Wi — Wi_l)]i < у ||Wi — 1._^_1||. < с_ у \щ — 1У._1\о,П <

2т. + 1 2т. + 1

Л/ с !СоС1/2

< С_! 2т. + 1 ^^ — и|0>П + |и — ^М < 24т. К_1(Кг2 + К2_1)-

Из построения сеток следует, что К. = К._1/2. Поэтому

5 _1 1/2

[Qi(Wi — /i_lWi_l)]i < 2С24тСэ+11 К._1. (53)

Для оценки второго слагаемого в правой части (52) используем (50):

/ч 1/2

[Qi/i_lei_l] < [/i_lei_l] = (LiIi_lei_l,Ii_lei_l)i/ = = (Я^^i_lei_l, ei_l)г1^21 = фi_lei_l, ei_l)г1_2 = [ei_l]i_l. Объединим эту оценку с (52) и (53):

5 _1 1/2

Ы < К._1 + [^1 ]._!• (54)

Далее для доказательства используем индукцию по г. При г = 0 имеем [e0]1 = 0. Из (54) при г = 1 следует, что

5 _1 1/2

Ы. < 2С4тСэ+ГКо, 2 2т1 + 1

т.е. (51) справедливо при г = 1 с константой с1э = 5с_1сэс1/2. Предположим, что (51) справедливо на уровне г — 1 :

._1 и -1]._1 " " "

[ei_l]i_l < с1э X

2т, — 1 7=1 3

Объединяя это неравенство с (54), получаем, что оценка (51) выполняется с константой

_ 5 -1 1/2

с19 — ^ С4 С9С14 •

□

Путем простого анализа последовательности вычислений с учетом разреженности матриц устанавливается верхняя оценка числа арифметических операций в каскадном алгоритме:

г

5 < С20 ^(т. + 021)% (55)

.7 = 1

с константами с20, с21, не зависящими от п. и т., но разными для итерационных процессов (44) и (45). Ясно, что для последнего процесса эти константы существенно меньше.

Число итераций т1,... ,т^-1 выберем так, чтобы минимизировать 5г как функцию от т1, • • • , тг_1 при фиксированной величине правой части (51) при г — I [3]. В результате тг выбирается как наименьшее целое, удовлетворяющее неравенству

пЛ-Л 1/2

2тг + 1 > (2тг + 1) -V"1 • (56)

Теорема 2. Погрешность каскадного алгоритма (43) с тг итерациями метода сопряженных градиентов (44) или метода простых итераций (45), где тг выбирается из условия (56) при фиксированном тг, оценивается как

|иг - то'г || I < С20 • (57)

2тг + 1

Для восполнения иг € И1 вектора иг справедлива оценка

( С22 1

||йг - и|п < ( С3С8 + 2тТ+Г ) ^ (58)

Число арифметических операций оценивается сверху величиной

5 < (с2зтг + 024)^ (59)

Доказательство. Для задачи Дирихле имеем оценку

П < 4г_гпь (60)

Из построения сеток следует равенство

— 2'_гйг • (61)

С учетом этим соотношений можно переписать (56) в виде

2тг + 1 > (2тг + l)23(г_i)/2• (62)

Используя это неравенство вместе с (61) в (51) при г — I и оценивая сумму геометрической прогрессии, получим

- иг]г < 019 V 23(7_г)/2 — 2в19—^^— V 2(7_г)/2 < С19—^^— •

[ ]г < 19 2тг + 1 19 2тг + 1 < уД - 119 2тг + 1

7=1 7=1

Привлекая эквивалентность норм (39), мы приходим к оценке (57) с константой

_ 2^2 _1

с22 = _ 1 с14 с1Э.

На основании неравенства треугольника, эквивалентности норм (12) и (11) и оценки (32) мы получаем цепочку неравенств

\\\ щ — и \\\ П < \\\ щ — \\\ п + \\\ — и \\\ п < \\\ и — wг \\\ 1 + сз||гУг — и||1,п < \\\ и — wг \\\ г + сэсзКг.

Вместе с уже доказанной оценкой (57) это приводит к (58).

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

т. < 1(2т, + 1)23(1_.)/2 + 1.

Используя это неравенство и соотношение (60) и оценивая сумму геометрической прогрессии, из (55) получим

' /2(2тг + 1)23(1_3)/2 + 1 + си) V-1 щ

Я < с2о £ (1(2т1 + 1)23(1_3)/2 + 1 + с^ 43"

(2т + 1)) £ 2(3_1)/2 + + с2^ £ 43_0 <

3 = 1

с*, 2(2т + 1)) £ 2(3_1)/2 + (± 1 с ^ £ 44~1'

/ 3=1 V / 3=1

< (с*-^т1 +с20 + 3 (с21 + \\ I ^

т. е. справедлива оценка (59) с константами

У2 1 , 4 ( +14

с23 =с20 72—1и с24 =с20 2—72 + 31с21 + 2).

□

Из оценки (58) вытекает, что число итераций на самом верхнем уровне следует выбирать из условия с3с8 ~ с22/(2тг + 1). Тогда погрешность итерационного процесса приобретает такой же порядок малости, как и погрешность дискретной задачи с учетом численного интегрирования. Хотя константы с3,с8 и с22 неизвестны, видна независимость т^ от числа уровней и количества неизвестных.

6. Формулировка многосеточного алгоритма с использованием несимметричного ^-цикла

Начнем с формулировки алгоритма несимметричного Ш-цикла для решения задачи (42) (МС.-алгоритма). Предположим, что мы имеем некоторое начальное приближение z0 € М. решения задачи (42). В результате получим приближенное решение zг1 = МС.^0, í.). Алгоритм МОг

Если i = 0, то полагаем

z0 = MGo(z0, f0) = L-1 fо

и алгоритм MGi завершен.

Если i > 0, то алгоритм состоит из нескольких шагов.

1. Положим uQ = z0 и вычислим

uk = uk_i - Tfc_i(LiUk_i - fi), k = 1,... ,m,

где

=_1 + cos a__= П (63)

Tk-1 = A*(cos a - cos(2k + 1)a)' a = 2m + 2. ( )

2. Положим gi_i = R^L^m — fi), где R : Mi ^ Mi_1 — оператор проектирования.

3. Возьмем zQ = 0 G Mi_1 и повторим алгоритм MGi_1 2 раза:

zS = MGi_1(zQ_1, gi_1), s = 1, 2.

4. Положим

z1 = MGi(z0, fi) = um - Ii_1 z2.

Теперь сформулируем полный многосеточный алгоритм для решения последовательности задач (42). Алгоритм FMG

1) uo = L^1 f о;

2) для i = 1, 2,... , l цикл {

2.1. Ui = 1i_1Ui_1;

2.2. повторить алгоритм MGi t раз : Ui := MGi(ui, fi).

В результате мы получим u — приближенное решение задачи (42) при i = l. Рассмотрим матричный полином (48), где параметры Tk заданы соотношениями (63). Для него выполняется оценка [13]

,, f „ ,, . с Q tg a/2 n

IILiQili < AQ^f-, a = ——. (64)

m +1 2m + 2

7. Сходимость алгоритма FMG

Лемма 2. Существует константа c25 такая, что имеет место оценка

IIL"1 - Ii_ 1 L__11 Ri|i < I5. (65)

AQ

Доказательство. Пусть gi G Mi — произвольный вектор. Рассмотрим систему уравнений

LiVi = gi. (66)

Любому вектору gi G Mi соответствует функция g G Hi такая, что (66) является системой Бубнова—Галеркина для задачи [13]:

о

L(v,w) = (g,w)n v w G W (П),

т. е. система (66) эквивалентна задаче

L(Vi,w) = (g,w)n V w е H\ При этом выполняется оценка (18), т.е.

|v - Vi|o,n < crh2|g|o,Q. (67)

Рассмотрим также систему уравнений

L iw* = gi,

которая эквивалентна задаче

L(w*, w) = (g, w)n V w е Hi. Очевидно, что для w* справедлива оценка вида (33), т.е.

|v - w*|o,n < cgh2. (68)

С помощью неравенства треугольника и оценок (67) и (68) получаем, что |Vi - wD*|o,n < |Vi - v|o,n + |v - wD*|o,n < (c7|g|o,n + Сд)Л| < С2бЛ|. С учетом эквивалентности норм (13) для изоморфных векторов получим оценку

IIVi - w*||i < С-1С26hi. (69)

Используя неравенство треугольника, можем записать неравенство

||L- 1 - Ii-i-Li_1iRi< II-1 - 1|i + + ||L-1 - Ii_ 1L"11Ri|i + ||Ii_ 1(L"11 - L-_11)Ri|i. (70)

Оценим первое слагаемое в правой части (70). На основании определения матричной нормы имеем

||L_1 - L_1||i = SUp = ||(L-1 -L-1)gi"i.

Si^Mi, 11 gi N i

Si=0

Пусть gi е Mi — вектор, на котором достигается супремум в правой части последнего равенства. Ему соответствуют векторы Vi и W*, являющиеся решениями систем уравнений LiVi = gi и Li\v* = gi соответственно. Тогда, используя оценку (69), получим неравенства

-1 Г-1Ц _ |(L_-1 - L_-1 ^Ji _ ||Vi - W*^ / c—1 c26hi

||_^_1 __. _ II У .__ 11 . . 11. <

. П. Уё.У.

Отсюда следует, что при фиксированном существует константа с27 такая, что

П^1 — < с27. (71)

Наконец, с помощью (40) получим, что

1 - < ^ (72)

А*

Для второго слагаемого в правой части (70) имеем оценку [13]

с28 А*

При доказательстве леммы 1 мы, в частности, получили оценку А* < с15А*• Используя ее, можем записать:

[Ь-1 - /г-Л-Д^г < ^• (73)

А*

Нам осталось оценить третье слагаемое в правой части (70). Пусть V £ Мг-1 — произвольный вектор. На основании эквивалентности норм (13) можем получить неравенство

'^г-^Нг ^ « -1

1М|г-1

Норма оператора интерполяции определяется соотношением

ИГ II _ '' 1г— 1V ' г

Н1г-1 "Мг — ёир

< 2е4-1С5. (74)

vGMi_l, " v " г-1

С учетом (74) получаем оценку

"1г-1"Мг^Мг_1 < 2с- 1с5. (75)

Из определения оператора проектирования следует равенство

(/г-^, — )г — (V, Яг-^г-^ V £ Мг-1, W £ Мг. (76)

Положим V — Яг—. Тогда с помощью неравенства Коши и оценки (75) для левой части (76) получим оценку

| (1г-1V, -)г| < " 1г-1V " г "-"г < " /¿-1 " М;^М;_1 НvНг-lНwНг < 2с-1С5 " Яг—" ¿-1 " — " (77) В результате из (76), (77) имеем

Н Яг— Н г-1

I- 4

< 2с-1С5. (78)

На основании определения нормы оператора проектирования

Н Яг— Н г-1

Н Яг Н Мг-1^Мг — 8Ир

w=0

из (78) в силу произвольности — следует оценка

НЯг НМг_1^Мг < 2с-Ч. (79)

Используя оценку (71), которая справедлива для любого г, и оценку (40), получим

А

А"11 - А-1Нг-1 < ^• (80)

В результате с помощью (75), (79) и (80) получаем оценку

л л 4с 2с2с1бС27 Н /г-1 -1 - -1)ЯгНг < Н/г-1 Н Мг^Мг_1 НЬг-11 - Ьг--1 Нг-1 Н Яг Н М;_1^М; < -5-• (81)

А*

Наконец, привлекая (72), (73) и (81) в (70), мы приходим к оценке (65) с константой

2 -2 2 с25 — с16с27 + с15с28 + 4с4 с5с16с27^

□

Объединяя (64) и (65), получим, что выполняется критерий сходимости

НА-1 - /г-Л--11ЯгНгНЬг^гНг < п(т), (82)

где

П(т) — • (83)

т + 1 4(т +1)2

Зафиксируем целое г £ [1,1 ] и определим линейный оператор Вг : Мг ^ Мг, ставящий в соответствие погрешности начального приближения ¿0 — —г - для решения задачи (42) погрешность конечной аппроксимации — —г - , полученной алгоритмом МСг, т.е. ¿1 —

Теперь мы можем сформулировать результат, характеризующий сходимость алгоритма МСг [13].

Теорема 3. Пусть операторы Ьг являются самосопряженными и положительно определенными. Пусть также выполняется критерий сходимости (82), где функция п(т) не зависит от г и стремится к нулю при т ^ то. Тогда для любого £ £ (0,1) существует т0 такое, что для всех т > т0:

НВНг < £, г — (84)

Из оценки (84) вытекает, что для любого е £ (0,1) существует т0, не зависящее от Л,г, такое, что для всех т > т0 погрешность приближенного решения задачи (42), полученного алгоритмом МСг, уменьшается с множителем е по сравнению с погрешностью начального приближения , т.е.

Н—г - < еН—г - z0Нг• (85)

Отметим, что на основании (83) множитель е в (85) зависит от числа итераций т следующим образом [13]:

е — 0((т + 1)-2)

Рассмотрим реализацию алгоритма ^МС. Выберем т так, чтобы выполнялось неравенство е < 1 для е из (85). Тогда мы можем выбрать целое ¿, такое что

е4 < С4С-1/4^ (86)

Приведем результат, характеризующий сходимость алгоритма ^МС [13, теорема 4.15].

Теорема 4. Пусть m выбрано так, что (85) выполняется с £ < 1 и пусть t удовлетворяет условию (86). Тогда при использовании алгоритма FMG мы получим последовательность приближенных решений щ, i = 1,... , I, задач (42), которые удовлетворяют оценке

||Wj - Uj||j < Cg hj—--, p = 4£tC-1C5 < 1,

c5 (1 - P)

а для их восполнений щ £ H справедлива оценка

|u - Ui|ü,n < Cgh2 ■ 2/(1 - p).

Наконец, сформулируем результат, содержащий оценку числа арифметических операций в алгоритме FMG (см. [13, лемма 5.1.4]).

Теорема 5. В алгоритме FMG для достижения точности, обусловленной погрешностью дискретной задачи с учетом численного интегрирования, т. е. для выполнения оценок

|W - щ|ü,Q < cgh2, |u - Ul|ü,n < 2cgh2

требуется O(n) арифметических операций. Список литературы

[1] Ильин В.П., Свешников В.М. О разностных методах на последовательности сеток // Численные методы механики сплошной среды: Информ. бюл. / АН СССР. Сиб. отд-ние. Отд-ние механики и процессов управления. 1971. Т. 2, № 1. С. 43-55.

[2] Deuflhard P. Cascadic conjugate-gradient methods for elliptic partial differential equations. I. Algorithm and numerical results // Technical Report SC 93-23. Berlin: Konrad-Zuse-Zentrum (ZIB), 1993.

[3] Shaidurov V.V. Some estimates of the rate of convergence for the cascadic conjugate-gradient method // Comp. Math. Applic. 1996. Vol. 31, N 4-5. P. 161-171.

[4] Shaidurov V.V. Some estimates of the rate of convergence for the cascadic conjugate-gradient method. Magdeburg, 1994 (Preprint Otto-von-Guericke Universität. № 4).

[5] Bornemann F.A., Deuflhard P. The cascadic multigrid method for elliptic problems // Numer. Math. 1996. Vol. 75. P. 135-152.

[6] Shaidurov V.V., Timmermann G. A cascadic multigrid algorithm for semilinear indefinite elliptic problems // Computing. 2000. Vol. 64. P. 349-366.

[7] Shaidurov V.V. Cascadic algorithm with nested subspaces in domains with curvilinear boundary / / Advanced Mathematics: Computations and Applications. Novosibirsk: Computing Center Publisher, 1995. P. 588-595.

[8] Shaidurov V.V., Tobiska L. The convergence of the cascadic conjugate-gradient method applied to elliptic problems in domain with re-entrant corners // Math. Comput. 2000. Vol. 69, N 230. P. 501-520.

[9] ГилЕВА Л.В. Каскадный многосеточный алгоритм в методе конечных элементов для трехмерной задачи Дирихле // Сиб. журн. вычисл. математики. 1998. Т. 1, № 3. C. 217-226.

[10] Гилева Л.В., ШАйдуров В.В. Каскадный многосеточный алгоритм в методе конечных элементов для трехмерной задачи Дирихле в области с криволинейной границей // Сиб. журн. вычисл. математики. 2002. Т. 5, № 2. C. 127-147.

[11] Ладыженская О.А., Уральцева Н.Н. Линейные и квазилинейные уравнения эллиптического типа. М.: Наука, 1973. 578 с.

[12] Сьярле Ф. Метод конечных элементов для эллиптических задач: пер. с англ. М.: Мир, 1980.

[13] Shaidürov V.V. Multigrid Methods for Finite Elements. Dordrecht: Kluwer Academic Publishers, 1995.

[14] Даутов Р.З., Карчевский М.М. Введение в теорию метода конечных элементов: учеб. пособие. Казань: Изд-во Каз. ГУ, 2004.

[15] Воеводин В.В., Кузнецов Ю.А. Матрицы и вычисления. М.: Наука, 1984. 320 с.

[16] Самарский А.А., Николаев Е.С. Методы решения сеточных уравнений. М.: Наука, 1978. 592 с.

Поступила в редакцию 11 ноября 2007 г., в переработанном виде —13 мая 2008 г.