УДК 519.615.5
О. В. Белякова
ВАРИАНТ МНОГОСЕТОЧНОГО МЕТОДА С ПОЛУУКРУПНЕНИЕМ
Для решения систем линейных уравнений с блочно-трехдиагональной матрицей предложен вариант многосеточного метода с полуукруп-нением. Представлены результаты численных экспериментов, исследующих эффективность предложенного алгоритма.
Вестник Балтийского федерального университета им. И. Канта. 2011. Вып. 10. С. 102 — 109.
For solving systems of linear equations with block tridiagonal matrices there is presenting a variant of multi-grid method with semi-coarse. There are giving the results of numerical experiments, that show high efficiency of presented method.
Ключевые слова: многосеточный метод, сетки с полуукрупнением, системы линейных уравнений, блочно-трехдиагональные матрицы.
Key words: multi-grid method, semi-coarse grids, systems of linear equation, block tridiagonal matrices.
Введение
При численном решении краевых задач для эллиптических уравнений методом сеток или Галёркина решаются системы линейных алгебраических уравнений с блочно-трехдиагональной матрицей
Au = f, (1)
A = blocktridiag(Li, Di, Ri), Li, Di, Ri <z ^mxm, i = 1,..., n. (2)
Одним из эффективных итерационных методов решения задачи (1 — 2) является многосеточный метод (MGM, MultiGrid Method). Главное его преимущество — независимая от числа узлов сетки скорость сходимости, что обеспечивает получение решения системы за время, пропорциональное числу неизвестных.
Основная идея MGM состоит в уменьшении низкочастотной части невязки с помощью решения вспомогательных систем, соответствующих грубым сеткам. Для этого вводится последовательность сеток с убывающим количеством узлов, с каждой из которых связана своя система. В большинстве работ [1 — 4] для построения матриц вспомогательных систем используется метод Галеркина, в двумерном случае при этом число узлов на каждой более крупной сетке уменьшается примерно в четыре раза, а в качестве сглаживающих итераций применяются простые итерации или вариант блочного метода Гаусса — Зейделя с чередованием по строкам (black-red).
В работе [5] для численного решения задачи Дирихле для двумерного уравнения диффузии методом сеток был предложен и исследован вариант MGM с полуукрупнением сетки MGMSC (SC — semi-coarse grids). В нем для решения каждой системы уравнений, соответствующей сетке с параметрами (hx, hy), используется система уравнений, соответствующая сетке с параметрами (2hx, hy), так что количество неизвестных во вспомогательной системе уменьшается в два раза, MGMSC. Такой способ уменьшения неизвестных позволил эффективно применить в качестве сглаживающих итераций блочный метод Гаусса — Зейделя с шахматным (черно-белым) упорядочением (zebra-line), рассмотреть дополнительные варианты построения вспомогательных матриц, в том числе вариант [6], основанный на исключении отдельных строк системы, подобно исключению в известном методе четно-нечетной редукции [7].
104
Для численного решения системы (1 — 2) с симметричной положительно определенной матрицей предлагается вариант многосеточного метода с полуукрупнением МОМ^С, использующий композицию МЄМ^с и МЄМ^с, исследуется его эффективность на численных примерах при различных способах построения матриц вспомогательных систем. Предложенный метод достаточно просто реализуется и во многих случаях имеет явное преимущество по сравнению с МСМ^
1^с.
1. Описание варианта многосеточного метода с полуукрупнением прямоугольной сетки
Пусть решается система (1), в которой А = ЫоскМ.^^._1, Б,, Ь,), Б,, Ь, е Жп‘хп-, А е ЖПхПХПхП, А >0. (3)
Векторы и и £ — блочные вектор-столбцы из пространства Пу:
£ = (Ц, Б2,^, ТПг )т, и = (и1, и2,^, Цх )т, Б, Ц еШПу,. = 1,..., пж.
Одна итерация предлагаемого варианта многосеточного метода МвМ^С состоит в следующем. Для решения (1 и 3) применим МСМ^с, использующий последовательность систем A( u( = £ , 1Х = 1,..., кх, в каждой из которых количество блоков Ц е ЭТПу , . = 1, ., пХ неизвестного вектора и, уменьшается, но размерность каждого Ц остается неизменной. Затем каждую из этих систем будем решать аналогично, применяя МвМ^с: А(!у)и|,у) = £{ХУ), 1у = 1,..., ку. Теперь в каждой системе количество блоков неизвестного вектора остается без изменения, но размерность блоков меняется. В результате получаются (кх _ 1) • (ку _ 1)
вспомогательных систем, каждая из которых определяется из предыдущей вычеркиванием нечетных столбцов или нечетных строк.
Матрицы А, вспомогательных систем строятся рекуррентно. При
I = 1 А, совпадает с матрицей заданной системы. Матрица А,+1 строится по выбранному правилу, например методом Галёркина: А;+1 = гх А;РХ, где рХ — оператор продолжения в пространстве ШПх, кото-
рый соответствует шаблону р =
, Г = (р )т — оператор сужения. По-
лучается матрица исходной структуры: Л!+1 = blocktridiag(L(l,+11), 0*,+1), Ь(1,+1)).
Такая же структура будет и у Л,+1, если она строится следующими методами:
• исключением строк [6]: Л,+1 = гегїисе(А,);
• модифицированным методом Галёркина [5]: Л,+1 = г*Л,рХ + Q,+1, где Q^+1 — некоторая вспомогательная матрица;
• Л,+1 = ге^се(Л,)+^+1.
Опишем одну итерацию МЄМ^с и МЄМ^С на псевдокоде (табл. 1).
Таблица 1
Псевдокоды многосеточных итераций MGMSC, MGMSC
Я>
MGMXc (ly), ly = 1 MGM5C
Smoothlter = 1,5; {количество сглаживающих итераций} for lx = 1 to kx - 1 do begin BlockGS_x (lx, ly, Smoothlter); {сглаживающие итерации Гаусса-Зейделя по строкам} Restriction_x (lx + 1, ly); {сужение на сетку с номером (lx + 1)} Coarse_f_x (lx, ly); {пересчет столбца свободных членов для следующей вспомогательной системы} end; lx = kx; { kx — номер самой грубой сетки} SolveInLine_x (lx, ly); {решение системы с трехдиагональной матрицей} While lx > 1 do begin Plolongation_x (lx, ly); {операция продолжения} lx = lx — 1; BlockGS_x (lx, ly, Smoothlter); end Smoothlter = 1,5; ly = 1; for lx = 1 to kx - 1 do begin BlockGS_x (lx, ly, Smoothlter); Restriction_x (lx + 1, ly); Coarse_f_x (lx, ly); end; lx = kx; SolveInLine_x (lx, ly); While lx > 1 do begin Plolongation_x (lx, ly); lx = lx — 1; HalfBlockGS_x (lx, ly); {половина итерации Гаусса — Зейделя} MGM|c (lx); end
105
2. Оценка вычислительной сложности алгоритма
Оценим и сравним количество арифметических операций, необходимых для выполнения каждой из итераций МСМ^с и МСМ^С. Поскольку элементы матриц вспомогательных систем вычисляются один раз, то не будем включать эти вычисления в общий объем работы. В многосеточной итерации наиболее емкими по количеству арифметических действий являются сглаживающие итерации Б1оскСБ_х, Б1оскСБ_у. Поэтому для оценки вычислительной сложности алгоритма будем учитывать только сглаживающие итерации.
Обозначим через ССБ • пх •пу количество арифметических операций
на выполнение сглаживающей итерации блочным методом Гаусса — Зейделя для приближенного решения системы размером пх -пу, где
ССБ — некоторая положительная константа, не зависящая от пх,пу . Подсчитаем число арифметических операций, необходимых для выполнения итерации МОМ^С с 1,5 сглаживающей итерации. При переходе от исходной системы уравнений к системе с наименьшим количеством неизвестных выполняются арифметические действия в количестве:
• 1,5Сся •пх-пу — для исходной системы (1х = 1 в таблице 1);
106
1,5Ссб • —^'пу — для следующей вспомогательной системы (1х = 2);
• 1,5С • —• п — для предпоследней вспомогательной системы
2 х у
(1х = кх — 1) (последняя вспомогательная система — трехдиагональная порядка пу, которая решается точно методом прогонки).
Всего операций
3 11 2^2кх-1 -1 3
-С“ п п(1+- + ...+-1Т) = 2 2. , 1 •-С“ • пх п « 3Сс5 п п .
2 х у 2 2 .х 2 .х-1 2 х у х у
Из таблицы 1 видно, что итерация МОМ^С реализована по схеме симметричного У-цикла. Поэтому такое же количество сглаживающих итераций затрачивается при обратном переходе: от системы с номером (кх - 1) к исходной системе с номером 1. Следовательно, общее количе-
ство арифметических операций будет примерно 6 С -пх-п.
у
Проанализируем теперь . Как видно из таблицы 1, количе-
ство арифметических действий на выполнение сглаживающих итераций при переходе от исходной системы, соответствующей мелкой сетке, к системе, соответствующей крупной сетке, такое же, как в МОМХС, 3ССБпх • пу. При обратном переходе для каждой системы выполняются:
1) половина итерации блочного метода Гаусса — Зейделя;
2) симметричный У-цикл итерации МОМ^С, всего
13 ^ сб „ „Г 1 , 1 , -Л п
С пхп\ +—-Г +... +1 |« 13С •пх-пу
2 х у і 2*х 2кх ~ I х у
арифметических действий. Итак, получим, что на выполнение одной итерации МвМ^С по схеме симметричного У-цикла необходимо затратить число арифметических действий примерно 16 • СОБ •пх-пу .
Если для реализации МОМ^С-итерации применить схему несимметричного У-цикла, а именно при переходах от системы, соответствующей мелкой сетке, к решению системы, соответствующей крупной, выполнять только половину сглаживающей итерации, то общее число арифметических действий уменьшится так, что получим 10СОБ • пх •пу .
3. Численные результаты
Исследуем на численных примерах эффективность вариантов МОМХС и МОМ^С многосеточного метода при различных способах построения матриц вспомогательных систем и схем выполнения У-цикла итерации. Для этого рассмотрим однородную задачу Дирихле в единичном квадрате О = (0,1) х (0,1) для уравнения
ах I( р(х'у > § ^ду ( <’(х'у) I )=-; •
Применение сеточного метода для численного решения этой задачи приводит к системе (1, 2). Предположим, что пх = пу = п, обозначим
||г< . 11 = £и|, нии^-1! Р=11т V,
,•,Л'"* "Г ГГ
где Гк) — невязка после к-й многосеточной итерации. Значение скорости сходимости р вычислялось в момент, когда первоначальная невязка уменьшилась не менее чем в 1010 раз; V],, VI — значения скоростей сходимости итераций МСМ;!С и МСМ^С соответственно.
Таблица 2
Скорость сходимости и количество итераций МСМ^е , МСМдС по схеме симметричного У-цикла в зависимости от коэффициентов системы и мелкости сетки при условии Лг+1 = ге^се(Лг ) + Qг+1
Р<х, у) = <?(х,у) п 99 256 512 777
1 Р<х, у ) =1 V1 (к) 0,05096(8) 0,05309(8) 0,05355(8) 0,05310(8)
V2 <к) 0,00071(3) 0,00067(3) 0,00062(3) 0,00060(3)
2. р(х, у) = 1 - е- х у V1 (к) 0,05240(8) 0,05376(8) 0,05409(8) 0,05405(8)
V2 (к) 0,00090(3) 0,00082(3) 0,00078(3) 0,00076(3)
3 Р< х, у) =1 + + 0,85т(14лх)5т(14лу) V1 (к) 0,07810(10) 0,08124(10) 0,08181(10) 0,08283(10)
V2 (к) 0,00184(4) 0,00139(4) 0,00070(3) 0,00046(3)
Видно, что для монотонных коэффициентов в уравнении МСМ£С -метод малоэффективен по сравнению с МСМ^С, так как для уменьшения первоначальной невязки в 1010 раз требуется примерно одинаковое число действий — 48Сп2. Строка 3 показывает, что МСМ^С может иметь некоторое преимущество при решении задач с быстро изменяющимися коэффициентами со значительной амплитудой колебаний. Если итерация МСМХС реализуется по схеме несимметричного У-цикла, то для уменьшения первоначальной невязки в 1010 раз требуется примерно 40Сп2 действий против 48Сп2 МСМ^С -методом. Следовательно, этот вариант МСМ£С -итерации более эффективен по сравнению с МСМ;!С.
Здесь и далее через к1, к2 обозначено количество многосеточных итераций МСМХС, МСМ^С соответственно. Из таблицы 3 видно, что в большинстве случаев есть явное преимущество МСМ^С-метода, особенно выразительно оно в случае 2: для уменьшения первоначальной невязки в 1010 раз требуется примерно 50Сп2 арифметических действий против 72Сп2 МСМ^С -методом.
107
108
Таблица 3
Количество итераций МСМ^с , МСМдС по схеме несимметричного У-цикла в зависимости от коэффициентов системы и мелкости сетки при условии
Л,+1 =геаисе(Л, ) + Qг+l
Р( х , у ), <?(х, у ) п 99 256 512 777
1. р(х, у) = 1 - е- х' у ц( х, у) = 1 + 0,8вт(14ях )вт(14яу) к1 11 11 11 11
к2 7 6 5 5
2. ц(х, у) = 1 - е- х' у р(х, у) = 1 + 0,8вт(14ях )вт(14яу) к1 12 13 13 13
к2 5 5 5 4
Из таблицы 4 видно, что при указанном способе построения матриц вспомогательных систем, более экономичном по сравнению с предыдущим способом, также имеются преимущества МСМ^С по сравнению с МСМХС. В большинстве приведенных примеров количество арифметических операций для уменьшения невязки в 1010 раз МСМХС -методом примерно в 1,5 раза меньше по сравнению с МСМ^С -методом.
Таблица 4
Количество итераций МСМ^с , МСМ^С по схеме несимметричного У-цикла в зависимости от коэффициентов системы и мелкости сетки при условии
Лг+1 = ге^се(Аг)
р(х, у ), ф, у) п 99 256 512 777
1 p(x, у) = q(x, у) =1 к1 7 8 8 8
к2 3 3 3 3
2. р( х, у ) = Ч(х, у ) =1 ех у к1 8 8 8 8
к2 3 3 3 3
3. р(х, у) = q(x, у) = 1 + 0,8вт(14ях )вт(14яу) к1 9 9 9 10
к2 5 4 4 4
4. р(х, у) = 1 - е-х' у q(x, у) = 1 + 0,8вт(14ях )вт(14яу) к1 10 10 10 10
к2 6 5 4 4
5. ц(х, у) = 1 - е-х' у р(х, у) = 1 + 0,8вт(14ях)вт(14яу) к1 12 13 13 13
к2 4 4 4 4
Из таблицы 5 видно, что МСМ^С -метод имеет явное преимущество по сравнению с МСМХС -методом.
-©
Таблица З
Скорость сходимости и количество итераций МСМ^е , МСМдС по схеме несимметричного У-цикла в зависимости от коэффициентов системы и мелкости сетки при условии А,+1 = Г А, р,
p(x, y), q(x, y) n 128 25б 512
L p(x, y) - q(x, y) -1 vk(k) 0,04958(8) 0,05043(8) 0,05087(8)
vl(k) 0,00132(3) 0,00133(3) 0,00133(3)
2. p(x, y) - q( x y) -1 ex y vl(k) 0,2бЗЗ9(17) 0,32б41(20) 0,38538(23)
vl (k) 0,01бб1(б) 0,01807(б) 0,0191б(б)
3. p(x y) - q(x, y) - 1 + + 0,8sin(14nx)sin(14ny) vl (k) 0,08387(10) 0,085бЗ(10) 0,08б07(10)
vl(k) 0,00419(5) 0,00248(4) 0,00178(4)
4. p(x, y) -1 - е-x' y q(x,y) - 1 + 0,8sin(14nx)sin(14ny) v1 (k) 0,11729(11) 0,1239б(11) 0,1530б(12)
vl(k) 0,03480(б) 0,0148б(5) 0,010б0(5)
5. q(x, y) -1 - е-^ p(x,y) - 1 + 0,8sin(14nx)sin(14ny) v1 (k) 0,1б512(13) 0,1б941(13) 0,17058(13)
v2k (k) 0,00250(4) 0,00204(4) 0,00180(4)
lO9
Список литературы
1. Mandel J., McCormick S., Ruge J. An algebraic theory for the multi-grid method including V-cycle // SIAM J. Numer. Anal. 1983. Vol. 20.
2. Hackbusch W. Multi-grid method and applications. Springer, 1985.
3. Maitre J.-F., Musy F. Multigrid methods for symmetric variational problems: a general theory and convergence estimates for usual smoothers // Appl. Math. and Comp. 1987. Vol. 21.
4. Ольшанский М. А. Лекции и упражнения по многосеточным методам. М., 2005.
5. Белякова О. В., Буздин А. А. Многосеточный метод с полуукрупнением для решения систем с блочной трехдиагональной матрицей // Методы вычислений. 2005. № 21. С. 5-19.
6. Буздин А. А., Дедух С. С. Вариант .многосеточного метода с полуукрупнением сетки / / Вестник Российского государственного университета им. И. Канта. Калининград, 2009. Вып. 10. С. 74 — 81.
7. Самарский А. А., Николаев Е. С. Методы решения сеточных уравнений. М., 1978.
Об авторе
Ольга Владиславовна Белякова — канд. физ.-мат. наук, доц., Балтийский федеральный университет им. И. Канта, e-mail: [email protected].
Author
Dr Olga Belyakova — assistant professor, I. Kant Baltic Federal University, e-mail: [email protected].