CC BY
32
ВЕСТНИК ТОМСКОГО ГОСУДАРСТВЕННОГО УНИВЕРСИТЕТА
2011 Математика и механика № 2(14)
УДК 519.632.4
А. А. Фомин, Л.Н. Фомина
УСКОРЕНИЕ ПОЛИНЕЙНОГО РЕКУРРЕНТНОГО МЕТОДА В ПОДПРОСТРАНСТВАХ КРЫЛОВА
На примере алгоритма ЬЮ полинейного рекуррентного метода [1, 2] рассматриваются два механизма его ускорения в подпространствах Крылова.
В качестве ускоряющего метода используется алгоритм ВьСОБ1аЪ Р ван дер Ворста. Показано, что традиционный подход: построение предобуславлива-теля на базе алгоритма ЬЮ, не приводит к требуемому результату. В то время как прямое сочетание алгоритмов ЬЮ и ВьСОБ1аЬ Р позволяет значительно повысить скорость сходимости решения.
Ключевые слова: разностные эллиптические уравнения, итерационный метод решения, подпространства Крылова, полинейный рекуррентный метод.
Несмотря на бурный рост вычислительной техники и теоретические достижения в области решения систем линейных алгебраических уравнений [3], возникающих при моделировании многомерных явлений математической физики, проблема создания и развития быстродействующих алгоритмов не теряет своей актуальности. Происходит это потому, что с усложнением поставленных задач возрастают и требования к точности и надежности методов решения.
Современные методы решения задач математической физики зачастую сводятся к разностной аппроксимации многомерных дифференциальных уравнений [4], что, в свою очередь, приводит к построению системы линейных алгебраических уравнений (СЛАУ), матрица которой имеет большую размерность и разреженноупорядоченную структуру. Для многомерных эллиптических по пространству задач до сих пор не удалось разработать прямой, экономичный, устойчивый к ошибкам округления метод, способный решать СЛАУ с линейной трудоемкостью относительно числа неизвестных, наподобие того, как это было сделано для одномерного случая. В настоящее время наиболее перспективными можно считать те направления разработки новых методов решения СЛАУ, в которых на уровне алгоритма учитывается фундаментальное свойство краевых эллиптических задач обязательной чувствительности решения в каждой точке области определения задачи возмущения в любой иной, включая граничную, точке. Такая чувствительность обеспечивается современными градиентными методами [3, 5], но не напрямую от точки к точке, а опосредованно, через коэффициенты, полученные при минимизации норм соответствующих функционалов. Однако непосредственно градиентные методы не характеризуются высокими скоростями сходимости, и для их ускорения используются предобуславливатели, которые строятся, как правило, на основе известных релаксационных методов [3]. При этом, чем более эффективен исходный релаксационный метод, тем более эффективным получается сочетание соответствующего предобуславливателя с градиентным методом. Например, использование в методе бисопряженных градиентов со стабилизацией предо-буславливателя на базе явного метода Булеева позволяет реализовать более высокие скорости сходимости решения по сравнению с предобуславливателем на базе метода Гаусса - Зейделя [2, 6].
Постановка задачи
Пусть имеет место двумерная по пространству краевая задача в единичной области О = {(х,у): 0 < х < 1,0 < у < 1}. Внутри области О поведение искомой функции Ф(х,у) описывается дифференциальным уравнением
идФ+удФ=Нх ЗФ^^у дфи (1)
дх ду дх V дх) ду V ду)
где и, V - компоненты скорости, Vх, Vу - коэффициенты переноса, £ - источник. На границе области Г имеют место условия первого рода
Ф|г=ф. (2)
Здесь и, V, Vх, Vу - непрерывные дифференцируемые функции внутри и на границе О, причем Vх > 0, Vу > 0 ; ф - непрерывная функция, определенная на грани-
це О. Выбирая произвольные выражения для функций Ф, и, V, Vх, Vу, можно всегда, путем вычисления источника £ из (1) и определения ф из (2) по значениям Ф на границе О, сформулировать необходимую для тестирования численного метода задачу. В настоящей работе рассматриваются две такие тестовые задачи. Задача 1. Решение и коэффициенты задаются в виде следующих зависимостей:
Ф(х, у) = 256 [ (1 - х)(1 - у)]2; и=V=0;
Vх (х,у) = 1 + 2[(х-0,5)2 +(у-0,5)2 ], (3)
Vу (х,у) = 1 + 2[0,5-(х-0,5)2 -(у-0,5)2],
а на границе Ф|г = 0.
Задача 2. Решение и коэффициенты задаются в виде следующих зависимостей: Ф(х, у) = ехр(-101 2)соє(8 п 12);
V х (х, у ) = Vу (х, у ) = ехр(-512); (4)
12 = х2 + у2; V(х,у) = у3/(1 + х2); на левой границе и (0, у) = 0, а внутри области О и на оставшихся границах и(х,у)
ди дV
рассчитывается из соотношения------I---= 0. На границе О решение Ф опреде-
дх ду
ляется из соотношения (4) при соответствующем выборе 1:1(0,у) = у2, 1(1,у) = 1+у2, 1(х,0) = х2, 1(х, 1) = 1+х2.
Системы линейных алгебраических уравнений вида ^4Ф = Ь получаются путем разностной аппроксимации задач (3) и (4) методом контрольного объема (вариант экспоненциальной схемы) [7] на равномерной сетке, покрывающей расчетную область О. В качестве зависимости коэффициентов схемы от сеточного числа Пекле Рш используется рекомендованная С. Патанкаром функция /(Рш) = тах (0,(1--0,1Рш)5). Как отмечается в [7], получаемая при этом разностная схема
аРц Фу = аЕ. Фі+1 . + ашц Фі-1 . + аЫ. Фу+1 + а£. Фу-1 + Ьу , (5)
консервативна и имеет второй порядок аппроксимации. Причем в данном случае аР = аЕ + аш + а£ + аы , поскольку уравнение (1) стационарно и источник £ не за-
висит явным образом от Ф. Здесь 1 < i < n, 1 < j < m , где n, m - количество узлов
сеточного разбиения расчетной области по координатам x, у соответственно. Из вида уравнения (5) следует, что матрица A системы имеет пятидиагональную структуру (см. рис. 1), причем для задачи (3) матричный оператор A будет самосопряжен, а для задачи (4) - нет. Следует отметить, что в общем случае способ разностной аппроксимации исходной дифференциальной задачи не имеет принципиального значения при условии, что получаемая при этом матрица СЛАУ имеет пятидиагональную структуру и положительный тип [1].
Все расчеты проводились на вычислительной системе Intel Core i5 750 2.66GHz,
RAM 4Gb, Win32, IVF v.11. Итерационный процесс прекращался при выполнении условия |[^||2 /||R0|| < £, где R0 и Rk соответственно начальная и текущая невязки, а е - заданная точность решения. Во всех расчетах область Q покрывалась равномерной сеткой 1001x1001, точность е принималась равной 10-8, а начальное приближение решения бралось в виде единичного вектора Ф0 = 1.
Оценка числа обусловленности матрицы системы для первой задачи составила MA и 5,90x105, а для второй - MA и 9,40x106.
Кроме нормы невязки в итерационном процессе также контролировались отношение норм текущей погрешности к начальной||zk|| /||z0|| и средняя скорость сходимости Qk [4].
Градиентный метод решения СЛАУ
Среди многочисленных итерационных градиентных методов решения разностных СЛАУ наиболее эффективным по праву считается метод бисопряженных градиентов со стабилизацией ван дер Ворста Bi-CGStab, основанный на процедуре биортогонализации Ланцоша в подпространствах Крылова [3]. Алгоритм этого метода для случая использования предобуславливателя имеет следующий вид (модификация Bi-CGStab P) [5].
1. Ф0 -вектор начального приближения решения;
2. r0 = b - A Ф0;
3. r - произвольный вектор, такой, что (r0, r0) Ф 0 , например F0 = r0;
4. р0 = а = ю0 = 1;
5. и0 = Р0 =0;
6. для k = 1,2,3, ...
7. Pk = (ro,rk-1); Р = (Pk/Pk-1)(а/®k-1);
8. Pk = rk-1 + Р(Pk-1- ®k-1uk-1);
9. Определение вектора у из решения системы Ky = pk ;
10. uk = Ay;
11. a = Pk/(n^,uk);
Рис. 1
12- 5=гк-\ -аик;
13. Определение вектора г из решения системы Кг = 5 ;
14. ї = Аг;
15. юк = (ї,5)/(ї,ї);
16. фк = фк-і +ау+®к2;
17. Если Фк достигло требуемой точности - выход из цикла;
18. гк = 5 -юкї;
19. Конец цикла по к.
Здесь К - предобуславливатель, причем предполагается разложение К = К1К2, где К1 - нижняя, а К2 - верхняя треугольные матрицы. При наличии такого разложения строки 9 и 13 алгоритма Бі-Єв8іаЬ Р выполняются с линейной трудоемкостью относительно числа неизвестных.
Предобуславливатель на базе алгоритма ЬШ
Следуя общей идеологии ускорения градиентного метода с помощью соответствующего предобуславливателя, необходимо произвести неполную факторизацию матрицы А системы уравнений на сомножители К1 и К2 таким образом, чтобы в основе этого разложения лежал алгоритм ЬЯ1. Из описания полинейного рекуррентного метода следует, что сомножитель К2 в виде четырехдиагональной почти верхнетреугольной матрицы присутствует в самом алгоритме [8]. Таким образом остается только определиться с сомножителем К1. Структура этого сомножителя легко выводится из расчетных формул преобразования правых частей уравнения системы [1], аналогично тому, как это делается в прямом методе Гаусса решения СЛАУ, который, как известно, эквивалентен Ш-разложению матрицы системы с последующим последовательным решением двух систем с нижней и верхней треугольной матрицами. На рис. 2 представлены структуры матрицы К1 (слева) и матрицы К2 (справа).
Нетрудно видеть, что матрица К имеет две особенности: единичную главную диагональ и полностью заполненные клетки, прилегающие к главной диагонали снизу. Вторая особенность является критичной, поскольку в этом случае трудоемкость метода будет не линейной, а пропорциональной N3/2, где N = пхш - количество неизвестных в системе. Для преодоления этого осложнения можно восполь-
зоваться оценкой порядка величин, заполняющих клетки под главной диагональю. Анализ расчетных формул [1] показывает, что на главной диагонали «поддиаго-нальных» клеток стоят величины порядка 1/3, на прилегающих к ней сверху и снизу диагоналях - (1/3)2, на следующих диагоналях - (1/3)3 и так далее.
Аналогично тому, как это делается, например в [4, 9], в данном случае удобно в «поддиагональных» клетках выделить диагональную ленту шириной 2/0+1, элементы которой заметно отличны от нуля, а всеми остальными малыми элементами этих клеток пренебречь. Легко проверить, что уже при 10 = 10 порядок величин на краю ленты будет 10 -5, а с другой стороны, даже при обычном сеточном разбиении для двумерных задач порядка 100х 100 лента шириной в 21 элемент не будет серьезным обременением для вычислительного процесса.
В дальнейшем для удобства алгоритм метода бисопряженных градиентов с предобуславливателем на базе ЬЯ1 называется ЬЯШСвБ. Для более детального анализа влияния полуширины ленты 10 на скорость сходимости ЬЮЬСОб решались задачи 1 и 2. Результаты решения первой задачи показали, что минимум количества итераций (33 - 35), необходимых для достижения заданной точности, приходится на диапазон 10 < 10 < 18; а для второй задачи минимум количества итераций составил 100 - 110 при 5 < 10 < 17 . С дальнейшим ростом 10 количество итераций вновь начинало возрастать. Объясняется это очевидно тем, что вновь привлекаемые к расчету диагонали уже не влияют на сходимость в силу малости расположенных на них величин, а с другой стороны, дополнительные ошибки округления все более эффективно тормозят процесс сходимости. Примечательно, что увеличение числа диагоналей в области минимально постоянного числа итераций практически не увеличивает время расчета задачи (15 - 16 с для первой задачи и 40 - 47 с для второй), что говорит в пользу приема усечения полностью заполненных «поддиагональных» клеток до относительно не широкой ленты.
Прямое сочетание алгоритмов ЬИ1 и Б1-С081аЬ Р
Решение задач 1 и 2 при тех же условиях напрямую алгоритмом ЬЯ1 продемонстрировало сходимость решения в первой задаче за 26 итерация (9 с), а во второй - за 101 итерацию (33 с). С одной стороны, такой результат говорит о высокой эффективности собственно алгоритма ЬЯ1, а с другой стороны, что не удалось достичь поставленной цели: ускорить алгоритм ЬЯ1 в подпространствах Крылова. И дело не только в том, что алгоритм ЬЯШСвБ работает медленнее (в секундах) - в конце концов, его реализация предполагает двойное решение системы с помощью предобуславливателя на базе ЬЯ1 - но и, что важно, количество итераций требуется, вообще говоря, больше.
Выход из создавшегося положения заключается в особенностях алгоритма Б1-С081аЪ Р, в котором используется только правое предобуславливание. При этом в нем не применяется в явном виде разложение К на сомножители К и К2, но в неявном виде это разложение, естественно, предполагается для того, чтобы иметь возможность решить системы Ку = рк (строка 9 алгоритма Б1-С081аЪ Р) и Кг = 5 (строка 13 алгоритма Ы-С081аЪ Р) с линейной трудоемкостью относительно числа неизвестных. Действительно, путем последовательного решения двух соответствующих пар систем К1х = рк, К2у = % и К1| = 5, К2 г = это легко сделать в силу треугольных структур матриц К и К2.
По определению предобуславливатель - это легко обратимая матрица K, в некотором смысле близкая к исходной, то есть K и А . При этом считается, что качество предобуславливателя тем выше, чем ближе он к исходной матрице системы. Вполне естественно эту близость рассматривать как близость решений x и xA двух систем Kx = b и AxA = b . С этой точки зрения вектор у, который является точным решением системы Ky = pk (строка 9 алгоритма Bi-CGStab P), можно трактовать как приближенное решение системы AyA = pk . Соответственно вектор z (строка 13 алгоритма Bi-CGStab P) - как приближенное решение системы AzA = s. Получается, что вместо того, чтобы точно решать вспомогательную систему с матрицей предобуславливателя K, можно просто найти приближенное решение той же вспомогательной системы, но уже с основной матрицей A, применяя какой-либо дополнительный метод, который для удобства можно условно обозначить как метод «YZ». В этом случае строки 9 и 13 алгоритма Bi-CGStab P должны быть заменены на следующие:
9'. Определение с помощью метода «YZ» вектора приближенного решения у системы AyA = pk ;
13'. Определение с помощью метода «YZ» вектора приближенного решения z системы AzA = s;
Иными словами идеологически (но не математически!) нет особой разницы между тем, чтобы точно решить приближенную систему или тем, чтобы приближенно решить точную. На этом основана технология «внедрения» метода «YZ» в алгоритм Bi-CGStab P, которая и составляет суть прямого сочетания «YZ» с Bi-CGStab P, поскольку ни предобуславливатель, ни, тем более, его факторизация не используются. Понятно, что метод «YZ» должен быть итерационным (в противном случае им сразу можно было бы решить исходную систему АФ = b). При этом вопрос об оптимальном количестве итераций при использовании метода «YZ» остается открытым. Также остается открытым вопрос о начальных приближениях векторов у и z. В настоящей работе в целях максимального сближения результатов предлагаемого подхода и технологии предобуславливания с неполной факторизацией исходной матрицы в качестве начальных приближений у и z использовались нулевые вектора.
Поскольку эффективность прямого сочетания двух методов зависит от точности приближения векторов у и z при минимальных затратах машинного времени для их нахождения, то отсюда следует, что в качестве «YZ» лучше всего использовать такой метод, который за первую итерацию способен резко понизить начальную погрешность решения (соответственно и невязку). В этом смысле применение какого-либо из алгоритмов полинейного рекуррентного метода представляется наиболее оправданным. Дело в том, что данный метод характеризуется «фирменной» особенностью: он резко (на 2 - 4 порядка) понижает первоначальную невязку за первый итерационный проход [1, 2, 8].
В итоге на основании изложенных рассуждений можно определить новый метод, представляющий собой прямое сочетание алгоритмов LR1 и Bi-CGStab P, в дальнейшем для удобства называемый LR1sK. По отношению к предыдущему варианту LR1bCGs метод LR1sK обладает двумя очевидными преимуществами: во-
первых, он опосредованно учитывает все элементы предполагаемой нижней треугольной матрицы Кь сохраняя при этом линейную трудоемкость относительно числа неизвестных; во-вторых, в силу точности реализации алгоритма ЬЯ1 он наследует его свойство быть прямым методом относительно линейного по координатам решения. При этом понятно, что если вместо ЬЯ1 использовать алгоритм ЬЯ2 [1,2], то вновь образованный метод будет прямым для квадратичного по координатам решения.
Вычислительный эксперимент
Для оценки рассмотренных методов ЬЯШСвБ и ЬЮбК был проведен их сравнительный анализ путем решения тестовых задач 1 и 2. Их эффективность сравнивалась с эффективностью метода ЬЯ1 и трех алгоритмов метода бисопряжен-ных градиентов со стабилизацией: В1-Св81аЪ (без предобуславливания),
ВьСв81аЪ Р ЬИ (с предобуславливанием на базе ЬИ разложения матрицы А), ВьСв81аЪ Р В (с предобуславливанием на базе явного метода Булеева). Итерационные параметры (в тех алгоритмах, в которых они присутствуют) подбирались таким образом, чтобы реализовать максимальную скорость сходимости решения. Тем самым производилась своеобразная оценка сверху эффективности каждого из рассматриваемых алгоритмов.
Кривые сходимости в виде зависимостей Ня^1 Ь/11 Я0112 от номера итерации представлены на рис. 3, а (первая задача) и рис. 3, б (вторая задача). Видно, что отличительной особенностью графиков, так или иначе связанных с методом би-сопряженных градиентов (кривые 1, 3 - 6), является их резко немонотонное поведение, что качественно хорошо согласуется с результатами [5, 6]. Еще одной особенностью является отсутствие сходимости решения при использовании алгоритма В1-Св81аЪ для задачи 2 (несамосопряженный оператор; рис. 3, б, кривая 6).
Рис. 3. Кривые сходимости: а - первая задача; б - вторая задача; кр. 1 - ЬЮЖ.; кр. 2 -ЬЯ1; кр. 3 - ЬЯ1ЪС08; кр. 4 - В1-СОБ1аЪ Р В; кр. 5 - В1-СОБ1аЪ Р ЬИ; кр. 6 - В1-СОБ1аЪ
Итоговые результаты решения обеих задач всеми рассматриваемыми методами представлены в таблице. Кроме количества итераций и времени, необходимых для завершения итерационного процесса, здесь также представлены значения
ІНІ2/ІІ г 0І2 и & на момент его завершения. Поскольку алгоритмом Бі-Св8іаЬ за 2500 итераций не было достигнуто требуемой точности, то вместо количества итераций указано минимальное значение Н,Кк1І2/ІІ^0ІІ2, полученное в процессе итераций. Обращает на себя внимание, что отношение норм погрешностей всего на 1 - 3 порядка отличается от заявленной точности решения, что говорит об эффективных разрешающих возможностях всех методов, несмотря на достаточно высокие значения чисел обусловленности задач.
Метод
ЬЮїзК ЬЮ ЬЮЬС08 Бі-СОБІаЬ Р Б Бі-СОБІаЬ Р ЬИ Бі-СОБІаЬ
Итерации 9 26 33 62 430 1847
Время, с 5,5 8,9 14,9 4,8 30,3 70,5
1 / 2 г0 2 4,8-10-7 2,6-10-7 1,0-10-5 9,8-10-6 3,2-10-5 1,2-10-5
0і 1,620 0,580 0,350 0,190 0,024 0,006
Итерации 24 101 101 235 2431 1,86-10 -3
Время, с 13,9 33,0 42,7 16,7 170,0 -
2 ік / 2 г0 2 1,4-10-6 1,6-10-5 5,6-10-6 2,2-10-6 4,3-10-6 -
& 0,560 0,110 0,120 0,055 0,005 -
Из поведения графиков и данных таблицы хорошо видно, что по общему количеству итераций, необходимых для достижения заданной точности решения и итоговым значениям средней скорости сходимости, алгоритмы, в которых присутствует ЬЯ1 (кривые 1 - 3), имеют явное преимущество над остальными алгоритмами, представленными на рис. 3 (кривые 4 - 6). Однако по количеству времени, необходимому для завершения итерационного процесса не все так очевидно. Для задачи 1 (самосопряженный оператор) минимальное время реализуется в случае использования алгоритма В1-Св81аЪ Р В, хотя по количеству необходимых итераций он в семь раз превосходит алгоритм ЬЮбК. В случае более сложной задачи 2 (несамосопряженный оператор) ситуация меняется в пользу алгоритма ЬЮбК - теперь он явный лидер как по числу итераций, так и по затрачиваемому на сходимость решения времени. Отсюда можно сделать предположение, что для более сложных задач алгоритм ЬЮбК должен уверенно демонстрировать свои преимущества по отношению к алгоритму В1-Св81аЪ Р В (по отношению к остальным рассмотренным алгоритмам его превосходство сомнений не вызывает).
Для проверки этого предположения были решены алгоритмами ЬШбК и В1-Св81аЪ Р В модифицированные задачи 1 и 2 - вместо условий Дирихле на всех границах области были выставлены условия Неймана. Кривые сходимости всех четырех расчетов представлены на рис. 4. Нетрудно видеть, что для достижения заданной точности при использовании алгоритма ЬЮбК потребовалось 17 итераций (задача 1) и 40 итераций (задача 2), а при использовании алгоритма В1-Св81аЪ Р В - 212 и 998 итераций соответственно. При этом затраты по времени для первой и второй задачи соответственно составили 10,1 с и 22,7 с для алгоритма ЬЮбК и 15,2 с и 69,8 с для алгоритма В1-Св81аЪ Р В. Иными словами предположение о преимуществе ЬШбК перед В1-Св81аЪ Р В при решении более сложных задач в данном случае подтвердилось.
Рис. 4. Задача 1: кр. 1 - LRlsK, кр. 3 - Bi-CGStab P B; задача 2: кр. 2 - LRlsK, кр. 4 - Bi-CGStab P B
Выводы
На основании полученных результатов можно сделать следующие выводы.
1. Выявлена структура сомножителей неполной факторизации исходной матрицы системы на базе алгоритма LR1 полинейного рекуррентного метода решения эллиптических СЛАУ.
2. Показано, что использование предобуславливателя на базе LR1 в методе би-сопряженных градиентов со стабилизацией не приводит к ускорению вычислений по отношению к использованию самого алгоритма LR1.
3. Установлено, что в силу особенностей алгоритма Bi-CGStab P с ним напрямую можно сочетать любой итерационный метод решения СЛАУ.
4. Показано, что прямое сочетание алгоритмов LR1 и Bi-CGStab P позволяет значительно повысить скорость сходимости решения по отношению к использованию исходного алгоритма LR1.
ЛИТЕРАТУРА
1. Фомина Л.Н. Использование полинейного рекуррентного метода с переменным параметром компенсации для решения разностных эллиптических уравнений // Вычислительные технологии. ИВТ СО РАН. 2009. Т. 14. № 4. C. 108-120.
2. Фомин А.А., Фомина Л.Н. Сравнение эффективности высокоскоростных методов решения разностных эллиптических СЛАУ // Вестник Томского государственного университета. Математика и механика. 2009. № 2(6). C. 71-77.
3. Yousef Saad. Iterative Methods for Sparse Linear Systems. - N.Y.: PWS Publ., 1996. 460 p.
4. Ильин В.П. Методы конечных разностей и конечных объемов для эллиптических уравнений. Новосибирск: Изд-во Ин-та математики, 2000. 345 с.
5. Van der Vorst H.A. BI-CGSTAB: a fast and smoothly converging variant of BI-CG for the solution of nonsymmetric linear systems // SIAM J. Sci. Stat. Comput. 1992. V. 13. No. 2. P. 631-644.
6. Старченко А.В. Сравнительный анализ некоторых итерационных методов для численного решения пространственной краевой задачи для уравнений эллиптического типа // Вестник ТГУ. Бюллетень оперативной научной информации. Томск: ТГУ, 2003. № 10. C. 70-80.
7. Патанкар С. Численные методы решения задач теплообмена и динамики жидкости. М.: Энергоатомиздат, 1984. 152 с.
8. Фомин А.А., Фомина Л.Н. Об одном варианте полинейного рекуррентного метода решения разностных эллиптических уравнений // Вестник Томского государственного университета. Математика и механика. 2010. № 2(10). С. 20-27.
9. Вшивков В.А., Засыпкина О.А. Итерационный метод решения СЛАУ первого порядка сходимости с регулируемой матрицей перехода // Сибирский журнал индустриальной математики. 2008. Т. 11. № 2. C. 40-49.
Статья поступила 11.09.2010 г.
Fomin A.А., Fomina L.N. ACCELERATION OF THE LINE-BY-LINE RECURRETNT METHOD IN KRYLOV SUBSPACES. Two techniques of acceleration of line-by-line recurrent method in Krylov subspaces are considered by the example of the LR1 algorithm. The van der Vorst Bi-CGStab P algorithm is used as an accelerating method. It is shown that the traditional approach (generation of a preconditioner on the base of LR1 algorithm) doesn’t yield the required result. At the same time, the direct combination of LR1 and Bi-CGStab P algorithms allows to raise the convergence speed considerably.
Keywords: difference elliptic equations, iterative method, Krylov subspaces, line-by-line recurrent method.
FOMIN Alexander Arkadyevich E-mail: fomin_aa@mail.ru
FOMINA Lubov Nikolaevna (Kemerovo State University)
E-mail: lubafomina@mail.ru
