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

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

Том 9, № 5, 2004

ДВУХПАРАМЕТРИЧЕСКИЙ ДВУЦИКЛИЧЕСКИЙ ИТЕРАЦИОННЫЙ МЕТОД РЕШЕНИЯ СИЛЬНО НЕСИММЕТРИЧНЫХ СИСТЕМ ЛИНЕЙНЫХ АЛГЕБРАИЧЕСКИХ УРАВНЕНИЙ*

Л. Г. ЧикинА, Б. Л. Крукиер Ростовский государственный университет, Ростов-на-Дону, Россия

e-mail: chikin@rsu.ru, bk@rsu.ru

An idea of including only skew-symmetric component of a matrix into the iterative operator of the triangular parts belongs to L.A. Krukier (1979). Further, this idea was developed by him and his collaborators for the product (PTKM) and the double cyclic (DTKM) triangular skew-symmetric methods, where only one relaxation parameter was included into the operator of the method. The novelty of the suggested methods consists in the generalization to the case of two different parameters (relaxation parameter is different from the parameter in the operator of the method) of the method DTKM and in the application of the new ideas for the proof of convergence.

Введение

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

В случае большой размерности системы для ее решения используются, как правило, итерационные методы. В настоящее время предложено большое количество методов решения систем линейных алгебраических уравнений (СЛАУ) [1, 2], большинство из которых эффективно работают со СЛАУ, матрицы их обладают некими дополнительными специальными свойствами, такими как симметричность, разреженность, знакопостоянство части или всех элементов, монотонность и др. Вместе с тем разработанный аппарат теории

* Работа выполнена при финансовой поддержке Российского фонда фундаментальных исследований (грант № 03-01-00005), гранта РФФИ и администрации Ростовской области (№ 04-01-96807) и программы Университеты России УР.03.01.024.

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

итерационных методов [2 - 4] не всегда удается эффективно использовать для решения задач с нарушением или отсутствием каких-то из вышеперечисленных свойств. Нас в первую очередь будут интересовать СЛАУ, матрицы которых потеряли свойство симметрии. Такие системы получаются, например, при наличии в дифференциальных уравнениях производных нечетных порядков. Ситуация еще более усложняется, когда несимметричная часть матрицы становится очень большой, т. е. норма кососимметричной части матрицы существенно превышает норму симметричной ее части. Такие СЛАУ называют сильно несимметричными. Матрицы такого класса появляются, например, в результате использования центрально-разностных схем при аппроксимации уравнения конвекции-диффузии с малым параметром при старшей производной.

Рассмотрим систему линейных алгебраических уравнений

Au = f. (1)

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

A = Alo + Ai,

Ao = (A + AT) /2 = AT,

Ai = (A - AT) /2 = -AT.

Для кососимметричной составляющей справедливо следующее разложение:

Ai = Ku + KL,

где Ku и Kl — строго нижняя и верхняя треугольные части матрицы A1, причем

Ku = -KT.

Определение 1. Матрица A называется диссипативной, если ее симметричная часть A0 положительно определена [2].

Большинство итерационных методов, которые применяются для решения линейных

систем, могут быть объединены общей формулой

yn+1 — yn

Bn--- + A Уп = f, (2)

где Bn — последовательность невырожденных матриц (операторов метода); A — исходная матрица; тп > 0 — последовательность итерационных параметров; y0 — начальное приближение; f — правая часть, f, y0 G H, H — конечномерное гильбертово пространство; yn — решение на n-й итерации. Этот итерационный процесс можно записать в эквивалентном виде

zn+1 = zn - ТпВ-1гп,

где rn = Ayn — f — вектор невязки на n-й итерации, а zn = yn — y — вектор ошибки этого метода (y-точное решение метода).

Определение 2. Итерационный метод (2) называется S-циклическим, если

Tn Bn Tn+S Bn+s

для любого п>0 и фиксированного s > 1 [4].

Определение 3. Итерационный метод (2) называется стационарным, если матрица Вп не зависит от номера итерации (оператор метода В = Вп является постоянной матрицей) [4].

Пусть матрица А в (1) диссипативна. Для решения системы (1) рассмотрим следующий двуциклический итерационный метод:

уп+1/2 _ уп

^ --— + А Уп = /,

т (3)

уп+1 _ уп+1/2 +1/2 ( )

Т ----+ А уп+1/2 = /.

т

В системе (3) т > 0 — итерационный параметр; ^ и Т — обратимые операторы метода; у0 — начальное приближение, / — правая часть, /, у0 € Н, Н — конечномерное гильбертово пространство, уп — решение на п-й итерации.

На первом этапе находится значение уп+1/2 как решение первого уравнения системы (3), а на втором этапе решается второе уравнение системы (3), из которой находится значение

уп+1

Цикл вычислений состоит в поочередном применении двух итерационных методов с операторами ^ и Т.

Определим погрешности гп, гп+1/2 и гп+1 как разности

гп = уп _ у,

гп+1/2 = уп+1/2 _ у,

гп+1 = уп+1 _ у

между решениями уп,уп+1/2 и уп+1в (3) и точным решением у исходной системы (1). Перейдем к уравнениям для погрешностей

~п+1 _ -п+1/2

^ ^-г-+ А гп+1/2 = 0,

т (4)

Т -+ А гп = 0

т

и перепишем их в следующем виде:

= _ тА) гп+1/2, Тгп+1/2 = (т _ тА) гп .

Проведя алгебраические преобразования и избавясь от промежуточной неизвестной гп+1/2, метод (4) можно переписать в операторном виде

гп+1 = Сгп, (5)

где оператор перехода С имеет вид произведения

С = С^ От. (6)

В равенстве (6)

С^ = ^-1 (^ _ тА) ,

От = Т-1 (Т - тА).

Для сходимости итерационного метода (3) в энергетическом пространстве Ир достаточно потребовать [2]

\^п+1\\в <р \\гп\\в , 0 < р < 1,Б = > 0. (7)

Таким образом, в силу соотношении (7) сходимость итерационного метода (3) целиком определяется оператором перехода (5).

Представим оператор В = Вп в (2) в виде суммы симметричной и кососимметричной составляющих

В = Во + Вь

В [5] впервые было предложено при решении сильно несимметричных задач использовать треугольные части К^ и Ки кососимметричной составляющей матрицы системы (1) для посторения оператора В итерационного метода (2), причем таким образом, чтобы его кососимметричная составляющая удовлетворяла соотношению

Вг = тАь (8)

Условие (8), связывающее кососимметричные составляющие матрицы исходной системы (1) и оператора итерационного метода (2), обеспечивает самосопряженность оператора В — тА, входящего в оператор перехода:

О = В-1 (В —ТА) = (Во +тА1)-1 (Во —тАо). (9)

Равенство (8) достигается, когда оператор метода В имеет следующую треугольную структуру [5, 6]:

^ = Е + 2 тКь (10)

или

Т = Е + 2тКи. (11)

Итерационный метод (2) с операторами метода (10) или (11) называется треугольным кососимметрическим методом (ТКМ).

Требование диссипативности оператора В итерационного метода (2) и (8), т.е. В0 = Вот > 0, позволяет преобразовать оператор перехода (9) к виду

О = В-1/2(Е + тВ-1/2АгВ-1/2 )-1(Е — тВ-1/2АоВ-1/2)В-1/2

и доказать следующую теорему.

Теорема 1. Пусть диссипативный оператор В удовлетворяет соотношению (8) [5]. Тогда для сходимости итерационного метода (2) в энергетическом пространстве ИВо достаточно выполнения неравенства

\\(Е + тР1)-1 (Е — тРо)\\ < 1

или

\\ (Е — тРо)\\ < 1, (12)

где Ро = В-1/2АоВ-1/2, р = В-1/2 А1В-1/2.

Следствием условия (12) является достаточное условие сходимости (2), (8) в виде операторного неравенства

B0 > 0, 5тА0. (13)

В [6] было получено обобщение результата (13).

Лемма 1. Пусть А и B — диссипат,ивные операторы в (2) и оператор B удовлетворяет соотношению (8). Тогда выполнение операторных неравенств для симметричных составляющих операторов А и B

1 - р t 1 + р

-- Bo < Ао <-- Bo,

тт

где 0 < р < 1, достаточно для того, чтобы для любого z0 £ H и погрешности задачи (2) выполнялась оценка

llzn+1ll <Р \\zn\\ Ir IIbq < р \\z \\Bo •

Используя операторы (10) и (11), в работе [7] исследован на сходимость в энергетических нормах Bl0 и Bjj0 двуциклический треугольный кососимметрический метод (ДТКМ)

yn+1 — yn+1/2

(E + 2tKl) У-У-+ Ayn+1/2 = f,

t (14)

yn+1/2 - yn

(E + 2tKl) y-^ + Ayn = f.

т

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

H.И. Булеевым [8, 9], развитым и обоснованным в работах В.П. Ильина [10] и R.S. Varga

[11]. Для случая M-матриц обоснования неполного разложения Холецкого даны в работе

[12].

I. Достаточное условие сходимости двухпараметрического ДТКМ

В соответствии с идеей из [5] построим в (3) операторы F и T так, чтобы они включали треугольные части кососимметричной составляющей матрицы системы и зависели от второго итерационного параметра ш > 0:

F = Bl = Diag (Bl)+ uKl,

(15)

T = Bj = Diag (Bj ) + wKj . v у

Каждый из операторов Bl и Bj из (15) представим в виде суммы симметричной и кососимметричной составляющих:

Blo = Diag (Bl) + 0, 5ш (Kl + KTL) , Bl1 = 0, 5ш (Kl + Kj) = 0, 5шАь Bjo = Diag (Bj) + 0, 5ш (Kj + Kj) , BJl = 0,5ш (KL + Kj) = 0, 5шА1.

Обратим внимание, что кососимметричные составляющие операторов Бь и Би равны между собой и пропорциональны кососимметричной составляющей исходной матрицы. Для оценки нормы оператора перехода двухпараметрического ДТКМ (3), (15)

С = СьСи,

произведем некоторые преобразования операторов Сь и Си. Рассмотрим оператор

Сь = Б-1 (Бь - тА)

и подставим в выражение для оператора Бь его симметричную и кососимметричную составляющие:

Сь = (Бь0 +0, А1)"1 (Бь0 +0, -тАо -ТА1) .

Оператор Бь — тА не является симметричным. Добавляя и отнимая в обеих скобках оператор 0.5^Ао, получим

Сь = (БЬо —0,5^Ао +0,5^А1 +0,5^Ао) 1 х х (БЬо —0, 5^ Ао +0,5^ А1 —т Ао —т А1 +0, 5^ Ао) .

Введем оператор

^0 = Бьо —0, 5^Ао = N^0,

что позволит записать оператор перехода Сь в следующем виде:

Сь = (Мьо +0,5^ А)-1 (Мьо — (т — 0,5и) А) . Потребуем, чтобы оператор Мь0 (Ми0) был положительно определен, т.е.

^о = Бьо —0, 5^ Ао = N¿0 > 0 {Мио = Бщ —0,5^ Ао = > 0) . (16)

Тогда возможны следующие преобразования:

Сь = Мь2 [е +0,Мь2 А М^2) [е — (т — 0 5^) М^2 А М^2) МД . Введем оператор

_ 1 _ 1 / _ 1 _ 1

Рь = Мь 2 А ш 2 \ г> — дг 2 А ]\т 2

и получим

Рь = Мьо2 А Мьо2 Ри = Мио2 А МиоМ (17)

_ 1 1

Сь = Мь02 (Е +0,5^ Рь)_1 (Е — (т — 0, 5и) Рь) МД

или

1 1

2 Л дг2

Сь = Мьо2 Сь Мь2о

где

Сь = (Е +0, 5^ Рь)_1 (Е — (т — 0, 5и) Рь). Таким образом, получили оператор перехода Сь, подобный исходному оператору Сь.

Аналогичными преобразованиями получим и оператор С и

_ 1 1 Си = С и ЛтД,, -1

Си = (Е +0, 5ц Ри) 1 (Е - (т - 0, 5ц) Ри) _ 1 1 Ри _

Жио2 АЖио2 .

Для сходимости итерационного метода достаточно [1, 2], чтобы норма оператора перехода была меньше 1, а в нашем случае — это норма произведения операторов:

1|С|| = ||Сь Си||<

Жь"о2 С ь N¿0

Ж^о2 С и N20

|С

|С

и

< 1,

что в свою очередь выполняется, если справедливы неравенства

С ь С и

Итак, оценим норму оператора Сь

|С

ь||мЬ|

< 1,

о

^о

|| Си |кг„ < 1.

(18)

С ь

вир

ж=0

((Е +0, 5ц Рь) 1 (Е - (т - 0, 5ц) Рь) х, (Е +0, 5ц Рь) 1 (Е - (т - 0,5ц) Рь) ж)

(ж, ж)

Операторы (Е + 0, 5цРь) 1 и (Е - (т - 0, 5ц) Рь) являются коммутативными, поэтому

С ь

вир

ж=0

((Е - (т - 0, 5ц) Рь) (Е +0, 5ц Рь)_ 1 ж, (Е - (т - 0, 5ц) Рь) (Е +0, 5ц Рь)-1)

(ж, ж)

Производя замену

у = (Е +0, 5ц Рь) 1 ж ^ ж = (Е +0, 5ц Рь) у,

получаем

С ь

= вир ■

у=0

1 - т т£

((Е - (т - 0, 5ц) Рь) у, (Е - (т - 0, 5ц) Рь) у)

((Е +0, 5ц Рь) у, (Е +0, 5ц Рь) у) 2 (Рь у, у) - (т - ц) (Рь у,Рь у)

(19)

У=0 ((у, у) + ц (Рь у, у) + (0,5ц)2 (Рь у, Рь у))' Из формулы (19) следует, что для выполнения первого из условий (18) достаточно

выполнения неравенства

2 (Рь у, у) - (т - ц) (Рь у,Рь у) > 0.

Произведя аналогичные выкладки, получаем, что для выполнения второго из условий (18) достаточно

2(Ри у, у) - (т - ц) (Ри у, Ри у) > 0.

Таким образом, получили достаточное условие сходимости метода (3), (15).

2

2

2

Теорема 2. Пусть выполняются условия (16). Тогда для сходимости, метода (3), (15) достаточно выполнения условий

2 (Рь у, у) _ (т _ и) (Рь у, Рь у) > 0,

2(Риу,у) _ (т _ и)(Риу,Риу) > 0, ( )

где Рь и Р^ определены формулами (17). Доказательство этой теоремы приведено выше.

Так как (Рь у, Рь у) > 0 и (Ри у, Ри у) > 0, условие (20) можно записать в виде

0 <т<и + 2(Рь«>

0 < т < и +

(Рь у, Рь у)' 2(Ри у, у) (Ри y, Ри у)'

Следствие 1. Пусть выполняются условия (16). Если оператор А в (1) диссипативен, то для сходимости метода (3), (15) достаточно выполнения условия

0 < т < и.

Доказательство. В условиях (20)

(Рь у, Рь у) > 0, (Ри у, Ри у) > 0.

Так как матрица А диссипативна, операторы Рь и Ри, заданные равенством (17), также диссипативны в силу свойств диссипативных операторов [2]. Значит, для выполнения условий (20) достаточно потребовать и _ т> 0 ^ 0 <т<и.

Следствие доказано.

2. Ускорение двуциклических методов

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

= Diag (Вь) + 0, 5и (К + ) _ 0, 5иАо > 0, = Diag (Ви) _ 0, 5и (К + КТ) _ 0, 5иАо > 0

или

= Diag (Вь) _ 0, 5иDiag (А0) + 0, 5и (Кь + КТ) _ 0, 5иА0 > 0, = Diag (В^) _ 0, 5иDiag (А0) _ 0,5и (Кь + ) _ 0, 5иА > 0,

N

где А0 = {А^.. , г = .7, А0 = К0ь + К0и.

Так как матрицы и симметричны, их собственные числа действительны и по теореме Гершгорина [1] для положительной определенности операторов и достаточно выполнения следующих условий:

{N¿0}« > 0,

гг >Е { N^0}

г.?

{N^0}гг > 0,

гг >£ {

г=7

г.?

причем хотя бы для одной строки в каждой системе неравенство должно быть строгим. Запишем условия (21) в развернутом виде:

Г (Вь) - 0, (АоШ > 0,

(Ви) - 0,(Ао)}гг > 0.

{^ (Вь)}« - 0, (Ао) ^ 0,5^£ ( | {К^- + {Ко^| + {к!}.. + {к^}

¿=7 ^

(Ви- 0,(Ао) ^ 0, 5^ {Ки^ + + + {К^. ) •

Тогда, очевидно, что если в качестве диагоналей операторов метода взять следующие выражения:

{^ (Вь)}« = 0,5^ ( |{кь^.| + |{А^| + |{к! }.|) +0, 5^ {Diag (Д,)}«,

¿=.7

{^ (Ви )}« = 0, 5^ ( |{Ки + |{Ао}г^.| + |{ки +0, 5^ {Diag (Д,)}«,

¿=7

то они обеспечат выполнение условий (16) и тем самым дадут сходимость двухпарамет-рического ДТКМ с диссипативным оператором системы, сохраняя при этом в операторах В^ и Ви информацию об изменениях в строках и столбцах исходной матрицы. Кроме того, построение диагоналей матриц В^и Ви не требует существенных вычислительных затрат, что не понижает эффективность метода.

3. Численные исследования ДТКМ с различными ускорителями

Численное исследование итерационных методов проводилось на следующей модельной задаче: в замкнутой области П = [0,1] х [0,1] рассматривалось стационарное уравнение конвекции-диффузии

1 д +1/ дм + д (^м) + дм + д (м2м)\ . ) Ре 2 \ дж дх ду ду /

(ймУ = 0, V = {м1,м2},

(22)

конвективные члены которого записаны в симметричном виде, т. е. как полусумма дивергентной и недивергентной форм записи. Если есть условия (22), все три формы записи уравнения конвекции-диффузии (дивергентная, недивергентная и симметричная) эквивалентны [13]. На границе области расчета ставились нулевые краевые условия 1-го рода.

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

Итерационный процесс прекращался, если

,(fc)|

-(0)1

< £, £

10

-6

где r(k) и r(0) - невязки соответственно на k-й и нулевой итерациях. В качестве точного решения бралась функция

U(x,y) = exy sinпх sinny,

обращающаяся в нуль на границе.

При проведении численного исследования рассмотрены четыре варианта задания коэффициентов при конвективных членах (табл. 1).

Коэффициенты скоростей были подобраны таким образом, чтобы удовлетворить условию (22) для каждой задачи. Расчеты проводились при числах Ре = 103,104,105. Исследовалось влияние числа Ре и итерационного параметра т на число итераций п.

Было проведено численное исследование трех методов решения системы (1) с диссипа-тивной матрицей А, получаемой в результате описанной конечно-разностной аппроксимации уравнения конвекции-диффузии:

— стандартного двуциклического треугольного кососимметрического метода с одним параметром т, где В^ = Е + 2тК^, Ви = Е + 2тКи, обозначаемого в табл. 2 ДТКМ(т);

Т аб л и ц а 1

Коэффициенты при конвективных членах

Задача 1 Задача 2 Задача 3 Задача 4

v1 =1, V2 =- 1 v1 =1 - 2x, V2 =2y-1 Vi = x + y, V2 = x - y v1 =sin2nx, v2= -2ny cos2nx

Таблица2

Количество итераций ДТКМ и SSOR

Pe ДТКМ(т) ДТКМ0(т) ДТКМ(2,т) SSOR(w) SSOR(w)

ДТКМ(2,т)

Задача 1

103 103 103 68 101 1,48

104 753 750 517 747 1,44

105 5725 5773 4126 5816 1,4

Задача 2

103 470 127 34 102 3

104 611 442 205 375 1,82

105 4733 3355 1201 2940 2,44

Задача 3

103 118 74 50 101 2

104 629 423 215 604 2,8

105 4733 3211 1851 4655 2,51

Задача 4

103 225 83 58 147 2,53

104 1601 476 319 1067 3,34

105 13714 3534 1590 7990 5,02

2

2

— двуциклического треугольного кососимметрического метода с одним параметром т и ускорителем D0 = 0.5Diag(KLKUKL), где BL = E-wD0 + 2tKl, B^ = E-wD0 + 2tKl, обозначаемого в таблице ДТКМ0(т);

— двуциклического треугольного кососимметрического метода с двумя параметрами — т и и, операторы метода первого и второго цикла имеют вид

Bl = Diag (Bl) + uKL, Bu = Diag (Bu) + uKи,

обозначаемого в таблице ДТКМ(и,т) с и = 2.

Достаточные условия сходимости методов ДТКМ(т) и ДТКМ0(т) исследовались в работе [7];

Для сходимости ДТКМ(2,т) достаточно взять диагональные матрицы в виде

Diag (Bl) = Diag (Bu ) =\ = J] I ^ | + ^ | kU] | + ^ | fcuji , (23)

l j j j J

так как по теореме Гершгорина [1] операторы

Blo - Ao > 0, Buo - Ao > 0

будут положительно определены.

Количество итераций каждого метода сравнивалось с эталонным методом SSOR [13]. Численно подтвердилось существование таких значений Topt, для которых при достижении заданной точности число итераций минимально.

Выводы

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

2. Получены достаточные условия сходимости метода в виде легко проверяемых операторных неравенств.

3. Предложен способ ускорения сходимости метода, основанный на построении диагонали обратимого оператора по специальной формуле (23).

4. Среди рассмотренных двуциклических методов ДТКМ(2,т) — самый эффективный метод решения сильно несимметричных систем линейных уравнений.

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

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

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

[3] Хейгеман Л., Янг Д. Прикладные итерационные методы. М.: Мир, 1986.

[4] Марчук Г.И. Методы вычислительной математики. М.: Наука, 1989.

[5] Крукиер Л.А. Итерационный метод решения неявных конечно-разностных схем аппроксимирующих один класс квазилинейных систем уравнений // Изв. вузов. Математика. 1979. № 7. C. 41-52.

[6] Крукиер Л.А., ЧикинА Л.Г. Кососимметрический итерационный метод решения стационарного уравнения конвекции-диффузии // Изв. вузов. Математика. 2000. № 11. C. 62-75.

[7] Крукиер Л.А., Чикина Л.Г. Двуциклический треугольный кососимметрический метод решения сильно несимметричных систем линейных алгебраических уравнений // Изв. вузов. Математика. 2001. № 5. C. 36-42.

[8] Булеев Н.И. Численный метод решения двумерных и трехмерных уравненний диффузии // Мат. сб. 1960. Т. 51, № 2.

[9] Булеев Н.И. Пространственная модель турбулентного обмена. М.: Наука, 1989. 344 с.

[10] Il'in V.P. Iterative incoplite factorization methods. Singapure, World Sci., Publ. Co., 1992.

[11] VARGA R.S. Matrix iterative analysis. Prentince Hall, Englewood Cliffs, NJ, 1962. 322 p.

[12] Meijerink J.A. Van der vorst an iterative solution methods for linear equation systems of which coefficient matrix is a symmetric M-matrix // Math. Comp. 1977. Vol. 31. P. 148-162.

[13] Крукиер Л.А., Мартынова Т.С. Влияние формы записи уравнения конвекции-диффузии на сходимость метода верхней релаксации // Журн. вычисл. математики и мат. физики. 1999. Т. 39, № 11. C. 1821-1827

Поступила в редакцию 9 апреля 2004 г., в переработанном виде — 8 июня 2004 г.