Метод Гаусса – Зейделя - Gauss–Seidel method

В числовая линейная алгебра, то Метод Гаусса – Зейделя, также известный как Метод Либмана или метод последовательного вытеснения, является итерационный метод используется для решения система линейных уравнений. Он назван в честь Немецкий математики Карл Фридрих Гаусс и Филипп Людвиг фон Зайдель, и аналогичен Метод Якоби. Хотя его можно применить к любой матрице с ненулевыми элементами на диагоналях, сходимость гарантируется только в том случае, если матрица либо строго по диагонали,^[1] или же симметричный и положительно определенный. Об этом упоминалось только в частном письме Гаусса своему ученику. Герлинг в 1823 г.^[2] До 1874 года Зайдель не выпустил ни одной публикации.

Описание

Метод Гаусса – Зейделя представляет собой итерационная техника для решения квадратной системы п линейные уравнения с неизвестными Икс:

{ Displaystyle А mathbf {x} = mathbf {b}}

.

Он определяется итерацией

{ Displaystyle L _ {*} mathbf {x} ^ {(k + 1)} = mathbf {b} -U mathbf {x} ^ {(k)},}

куда ${ Displaystyle mathbf {х} ^ {(к)}}$ это kое приближение или итерация ${ Displaystyle mathbf {x}, , mathbf {x} ^ {(к + 1)}}$ следующий или k + 1 итерация ${ displaystyle mathbf {x}}$ , а матрица А раскладывается на нижний треугольный компонент ${ displaystyle L _ {*}}$ , а строго верхнетреугольный компонент U: ${ displaystyle A = L _ {*} + U}$ .^[3]

Более подробно выпишите А, Икс и б в их составных частях:

{ displaystyle A = { begin {bmatrix} a_ {11} & a_ {12} & cdots & a_ {1n} a_ {21} & a_ {22} & cdots & a_ {2n} vdots & vdots & ddots & vdots a_ {n1} & a_ {n2} & cdots & a_ {nn} end {bmatrix}}, qquad mathbf {x} = { begin {bmatrix} x_ {1} x_ {2} vdots x_ {n} end {bmatrix}}, qquad mathbf {b} = { begin {bmatrix} b_ {1} b_ {2} vdots b_ {n} end {bmatrix}}.}

Тогда разложение А в его нижнюю треугольную составляющую и строго верхнюю треугольную составляющую:

{ displaystyle A = L _ {*} + U qquad { text {where}} qquad L _ {*} = { begin {bmatrix} a_ {11} & 0 & cdots & 0 a_ {21} & a_ {22 } & cdots & 0 vdots & vdots & ddots & vdots a_ {n1} & a_ {n2} & cdots & a_ {nn} end {bmatrix}}, quad U = { begin { bmatrix} 0 & a_ {12} & cdots & a_ {1n} 0 & 0 & cdots & a_ {2n} vdots & vdots & ddots & vdots 0 & 0 & cdots & 0 end {bmatrix}}.}

Систему линейных уравнений можно переписать как:

{ Displaystyle L _ {*} mathbf {x} = mathbf {b} -U mathbf {x}}

Теперь метод Гаусса – Зейделя решает левую часть этого выражения для Икс, используя предыдущее значение для Икс с правой стороны. Аналитически это можно записать как:

{ displaystyle mathbf {x} ^ {(k + 1)} = L _ {*} ^ {- 1} ( mathbf {b} -U mathbf {x} ^ {(k)}).}

Однако, пользуясь преимуществом треугольной формы ${ displaystyle L _ {*}}$ , элементы Икс^(k+1) можно вычислить последовательно, используя прямая замена:

{ displaystyle x_ {i} ^ {(k + 1)} = { frac {1} {a_ {ii}}} left (b_ {i} - sum _ {j = 1} ^ {i-1 } a_ {ij} x_ {j} ^ {(k + 1)} - sum _ {j = i + 1} ^ {n} a_ {ij} x_ {j} ^ {(k)} right), quad i = 1,2, dots, n.}

^[4]

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

Обсуждение

Поэлементная формула метода Гаусса – Зейделя чрезвычайно похожа на формулу метода Метод Якоби.

Расчет Икс^(k+1) использует элементы Икс^(k+1) которые уже были вычислены, и только элементы Икс^(k) которые не были вычислены в итерации k + 1. Это означает, что, в отличие от метода Якоби, требуется только один вектор хранения, поскольку элементы могут быть перезаписаны по мере их вычисления, что может быть выгодно для очень больших задач.

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

Гаусс-Зайдель такой же, как SOR (последовательное чрезмерное расслабление) с ${ displaystyle omega = 1}$ .

Конвергенция

Свойства сходимости метода Гаусса – Зейделя зависят от матрицы А. А именно, известно, что процедура сходится, если:

А симметричен положительно определенный,^[5] или же
А строго или несводимо диагонально доминирующий.^[6]

Метод Гаусса – Зейделя иногда сходится, даже если эти условия не выполняются.

Алгоритм

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

Входы: А, бВыход:  ${ displaystyle phi}$ Выберите первоначальное предположение  ${ displaystyle phi}$  к решениюповторение до схождения за я из 1 до того как п делать         ${ displaystyle sigma leftarrow 0}$         за j из 1 до того как п делать            если j ≠ я тогда                 ${ displaystyle sigma leftarrow sigma + a_ {ij} phi _ {j}}$             конец, если        конец (j-петля)  ${ displaystyle phi _ {i} leftarrow { frac {1} {a_ {ii}}} (b_ {i} - sigma)}$     конец (я-loop) проверить, достигнута ли сходимостьконец (повторение)

Примеры

Пример для матричной версии

Линейная система, показанная как ${ Displaystyle А mathbf {x} = mathbf {b}}$ дан кем-то:

{ displaystyle A = { begin {bmatrix} 16 & 3 7 & -11 end {bmatrix}}}

и

{ displaystyle b = { begin {bmatrix} 11 13 end {bmatrix}}.}

Мы хотим использовать уравнение

{ displaystyle mathbf {x} ^ {(k + 1)} = L _ {*} ^ {- 1} ( mathbf {b} -U mathbf {x} ^ {(k)})}

в виде

{ displaystyle mathbf {x} ^ {(k + 1)} = T mathbf {x} ^ {(k)} + C}

куда:

{ displaystyle T = -L _ {*} ^ {- 1} U}

и

{ displaystyle C = L _ {*} ^ {- 1} mathbf {b}.}

Мы должны разложить ${ displaystyle A _ {} ^ {}}$ в сумму нижней треугольной составляющей ${ displaystyle L _ {*} ^ {}}$ и строгий верхнетреугольный компонент ${ displaystyle U _ {} ^ {}}$ :

{ displaystyle L _ {*} = { begin {bmatrix} 16 & 0 7 & -11 end {bmatrix}}}

и

{ displaystyle U = { begin {bmatrix} 0 & 3 0 & 0 end {bmatrix}}.}

Обратное ${ displaystyle L _ {*} ^ {}}$ является:

{ displaystyle L _ {*} ^ {- 1} = { begin {bmatrix} 16 & 0 7 & -11 end {bmatrix}} ^ {- 1} = { begin {bmatrix} 0,0625 & 0,0000 0,0398 & -0.0909 конец {bmatrix}}}

.

Теперь мы можем найти:

{ displaystyle T = - { begin {bmatrix} 0,0625 & 0,0000 0,0398 & -0,0909 end {bmatrix}} times { begin {bmatrix} 0 & 3 0 & 0 end {bmatrix}} = { begin {bmatrix} 0,000 & -0,1875 0,000 & -0,1194 end {bmatrix}},}

{ displaystyle C = { begin {bmatrix} 0,0625 & 0,0000 0,0398 & -0,0909 end {bmatrix}} times { begin {bmatrix} 11 13 end {bmatrix}} = { begin { bmatrix} 0,6875 - 0,7439 end {bmatrix}}.}

Теперь у нас есть ${ displaystyle T _ {} ^ {}}$ и ${ displaystyle C _ {} ^ {}}$ и мы можем использовать их для получения векторов ${ displaystyle mathbf {x}}$ итеративно.

Прежде всего, мы должны выбрать ${ Displaystyle mathbf {х} ^ {(0)}}$ : мы можем только догадываться. Чем точнее предположение, тем быстрее будет работать алгоритм.

Мы предполагаем:

{ displaystyle x ^ {(0)} = { begin {bmatrix} 1.0 1.0 end {bmatrix}}.}

Затем мы можем вычислить:

{ displaystyle x ^ {(1)} = { begin {bmatrix} 0,000 & -0,1875 0,000 & -0,1193 end {bmatrix}} times { begin {bmatrix} 1.0 1.0 end {bmatrix} } + { begin {bmatrix} 0,6875 - 0,7443 end {bmatrix}} = { begin {bmatrix} 0,5000 - 0,8636 end {bmatrix}}.}

{ displaystyle x ^ {(2)} = { begin {bmatrix} 0,000 & -0,1875 0,000 & -0,1193 end {bmatrix}} times { begin {bmatrix} 0,5000 - 0,8636 end {bmatrix }} + { begin {bmatrix} 0,6875 - 0,7443 end {bmatrix}} = { begin {bmatrix} 0,8494 - 0,6413 end {bmatrix}}.}

{ displaystyle x ^ {(3)} = { begin {bmatrix} 0,000 & -0,1875 0,000 & -0,1193 end {bmatrix}} times { begin {bmatrix} 0,8494 - 0,6413 конец {bmatrix}} + { begin {bmatrix} 0,6875 - 0,7443 end {bmatrix}} = { begin {bmatrix} 0,8077 - 0,6678 end {bmatrix}}.}

{ displaystyle x ^ {(4)} = { begin {bmatrix} 0,000 & -0,1875 0,000 & -0,1193 end {bmatrix}} times { begin {bmatrix} 0.8077 - 0,6678 end {bmatrix }} + { begin {bmatrix} 0,6875 - 0,7443 end {bmatrix}} = { begin {bmatrix} 0,8127 - 0,6646 end {bmatrix}}.}

{ displaystyle x ^ {(5)} = { begin {bmatrix} 0,000 & -0,1875 0,000 & -0,1193 end {bmatrix}} times { begin {bmatrix} 0,8127 - 0,6646 end {bmatrix }} + { begin {bmatrix} 0,6875 - 0,7443 end {bmatrix}} = { begin {bmatrix} 0,8121 - 0,6650 end {bmatrix}}.}

{ displaystyle x ^ {(6)} = { begin {bmatrix} 0,000 & -0,1875 0,000 & -0,1193 end {bmatrix}} times { begin {bmatrix} 0,8121 - 0,6650 end {bmatrix }} + { begin {bmatrix} 0,6875 - 0,7443 end {bmatrix}} = { begin {bmatrix} 0,8122 - 0,6650 end {bmatrix}}.}

{ displaystyle x ^ {(7)} = { begin {bmatrix} 0,000 & -0,1875 0,000 & -0,1193 end {bmatrix}} times { begin {bmatrix} 0,8122 - 0,6650 end {bmatrix }} + { begin {bmatrix} 0,6875 - 0,7443 end {bmatrix}} = { begin {bmatrix} 0,8122 - 0,6650 end {bmatrix}}.}

Как и ожидалось, алгоритм сходится к точному решению:

{ displaystyle mathbf {x} = A ^ {- 1} mathbf {b} приблизительно { begin {bmatrix} 0.8122 - 0.6650 end {bmatrix}}.}

На самом деле матрица A строго диагонально доминирует (но не положительно определена).

Другой пример для матричной версии

Другая линейная система, обозначенная как ${ Displaystyle А mathbf {x} = mathbf {b}}$ дан кем-то:

{ displaystyle A = { begin {bmatrix} 2 & 3 5 & 7 end {bmatrix}}}

и

{ displaystyle b = { begin {bmatrix} 11 13 end {bmatrix}}.}

Мы хотим использовать уравнение

{ displaystyle mathbf {x} ^ {(k + 1)} = L _ {*} ^ {- 1} ( mathbf {b} -U mathbf {x} ^ {(k)})}

в виде

{ displaystyle mathbf {x} ^ {(k + 1)} = T mathbf {x} ^ {(k)} + C}

куда:

{ displaystyle T = -L _ {*} ^ {- 1} U}

и

{ displaystyle C = L _ {*} ^ {- 1} mathbf {b}.}

Мы должны разложить ${ displaystyle A _ {} ^ {}}$ в сумму нижней треугольной составляющей ${ displaystyle L _ {*} ^ {}}$ и строгий верхнетреугольный компонент ${ displaystyle U _ {} ^ {}}$ :

{ displaystyle L _ {*} = { begin {bmatrix} 2 & 0 5 & 7 end {bmatrix}}}

и

{ displaystyle U = { begin {bmatrix} 0 & 3 0 & 0 end {bmatrix}}.}

Обратное ${ displaystyle L _ {*} ^ {}}$ является:

{ displaystyle L _ {*} ^ {- 1} = { begin {bmatrix} 2 & 0 5 & 7 end {bmatrix}} ^ {- 1} = { begin {bmatrix} 0,500 и 0,000 - 0,357 и 0,143 end {bmatrix}}}

.

Теперь мы можем найти:

{ displaystyle T = - { begin {bmatrix} 0,500 и 0,000 - 0,357 и 0,143 end {bmatrix}} times { begin {bmatrix} 0 & 3 0 & 0 end {bmatrix} } = { begin {bmatrix} 0,000 & -1,500 0,000 & 1.071 end {bmatrix}},}

{ displaystyle C = { begin {bmatrix} 0,500 и 0,000 - 0,357 и 0,143 end {bmatrix}} times { begin {bmatrix} 11 13 end {bmatrix}} = { begin {bmatrix} 5.500 - 2.071 end {bmatrix}}.}

Теперь у нас есть ${ displaystyle T _ {} ^ {}}$ и ${ displaystyle C _ {} ^ {}}$ и мы можем использовать их для получения векторов ${ displaystyle mathbf {x}}$ итеративно.

Прежде всего, мы должны выбрать ${ Displaystyle mathbf {х} ^ {(0)}}$ : мы можем только догадываться. Чем точнее предположение, тем быстрее будет выполнен алгоритм.

Мы предполагаем:

{ displaystyle x ^ {(0)} = { begin {bmatrix} 1.1 2.3 end {bmatrix}}.}

Затем мы можем вычислить:

{ displaystyle x ^ {(1)} = { begin {bmatrix} 0 & -1,500 0 & 1.071 end {bmatrix}} times { begin {bmatrix} 1.1 2.3 end { bmatrix}} + { begin {bmatrix} 5.500 - 2.071 end {bmatrix}} = { begin {bmatrix} 2.050 0.393 end {bmatrix}}.}

{ displaystyle x ^ {(2)} = { begin {bmatrix} 0 & -1.500 0 & 1.071 end {bmatrix}} times { begin {bmatrix} 2.050 0.393 end { bmatrix}} + { begin {bmatrix} 5.500 - 2.071 end {bmatrix}} = { begin {bmatrix} 4.911 - 1.651 end {bmatrix}}.}

{ Displaystyle х ^ {(3)} = cdots. ,}

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

{ displaystyle mathbf {x} = A ^ {- 1} mathbf {b} = { begin {bmatrix} -38 29 end {bmatrix}}}

не гарантируется и в этом случае не произойдет.

Пример для версии уравнения

Предположим, что k уравнения, где Икс_п - векторы этих уравнений и отправная точка Икс₀. Из первого уравнения решите для Икс₁ с точки зрения ${ displaystyle x_ {n + 1}, x_ {n + 2}, dots, x_ {n}.}$ Для следующих уравнений подставьте предыдущие значенияИксс.

Чтобы было понятно, рассмотрим пример.

{ displaystyle { begin {array} {rrrrl} 10x_ {1} & - x_ {2} & + 2x_ {3} && = 6, - x_ {1} & + 11x_ {2} & - x_ {3 } & + 3x_ {4} & = 25, 2x_ {1} & - x_ {2} & + 10x_ {3} & - x_ {4} & = - 11, & 3x_ {2} & - x_ { 3} & + 8x_ {4} & = 15. end {array}}}

Решение для ${ displaystyle x_ {1}, x_ {2}, x_ {3}}$ и ${ displaystyle x_ {4}}$ дает:

{ displaystyle { begin {align} x_ {1} & = x_ {2} / 10-x_ {3} / 5 + 3/5, x_ {2} & = x_ {1} / 11 + x_ { 3} / 11-3x_ {4} / 11 + 25/11, x_ {3} & = - x_ {1} / 5 + x_ {2} / 10 + x_ {4} / 10-11 / 10, x_ {4} & = - 3x_ {2} / 8 + x_ {3} /8+15/8.end {выровнено}}}

Предположим, мы выбираем $(0, 0, 0, 0)$ в качестве начального приближения первое приближенное решение дается выражением

{ displaystyle { begin {align} x_ {1} & = 3/5 = 0,6, x_ {2} & = (3/5) /11+25/11=3/55+25/11=2.3272 , x_ {3} & = - (3/5) / 5 + (2.3272) /10-11/10=-3/25+0.23272-1.1=-0.9873, x_ {4} & = - 3 (2.3272) / 8 + (- 0.9873) /8+15/8=0.8789.end {выровнено}}}

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

{ displaystyle { begin {array} {llll} x_ {1} & x_ {2} & x_ {3} & x_ {4} hline 0.6 & 2.32727 & -0.987273 & 0.878864 1.03018 & 2.03694 & -1.01446 & 0 .984341 1.00659 & 2.00356 & -1.00253 & 0.998351 1.00086 & 2.0003 & -1.00031 & 0.99985 end {array}}}

Точное решение системы $(1, 2, -1, 1)$ .

Пример использования Python и NumPy

Следующая численная процедура просто выполняет итерацию для получения вектора решения.

импорт тупой в качестве нпITERATION_LIMIT = 1000# инициализируем матрицуА = нп.множество([[10., -1., 2., 0.],              [-1., 11., -1., 3.],              [2., -1., 10., -1.],              [0., 3., -1., 8.]])# инициализировать вектор RHSб = нп.множество([6.0, 25.0, -11.0, 15.0])Распечатать(«Система уравнений:»)за я в классифицировать(А.форма[0]):    ряд = ["{0: 3g}*Икс{1}".формат(А[я, j], j + 1) за j в классифицировать(А.форма[1])]    Распечатать("[{0}] = [{1: 3g}]".формат(" + ".присоединиться(ряд), б[я]))Икс = нп.zeros_like(б)за it_count в классифицировать(1, ITERATION_LIMIT):    x_new = нп.zeros_like(Икс)    Распечатать("Итерация {0}: {1}".формат(it_count, Икс))    за я в классифицировать(А.форма[0]):        s1 = нп.точка(А[я, :я], x_new[:я])        s2 = нп.точка(А[я, я + 1 :], Икс[я + 1 :])        x_new[я] = (б[я] - s1 - s2) / А[я, я]    если нп.всезакрыть(Икс, x_new, rtol=1e-8):        перемена    Икс = x_newРаспечатать("Решение: {0}".формат(Икс))ошибка = нп.точка(А, Икс) - бРаспечатать("Ошибка: {0}".формат(ошибка))

Производит вывод:

Система из уравнения:[ 10*x1 +  -1*x2 +   2*x3 +   0*x4] = [  6][ -1*x1 +  11*x2 +  -1*x3 +   3*x4] = [ 25][  2*x1 +  -1*x2 +  10*x3 +  -1*x4] = [-11][  0*x1 +   3*x2 +  -1*x3 +   8*x4] = [ 15]Итерация 1: [ 0.  0.  0.  0.]Итерация 2: [ 0.6         2.32727273 -0.98727273  0.87886364]Итерация 3: [ 1.03018182  2.03693802 -1.0144562   0.98434122]Итерация 4: [ 1.00658504  2.00355502 -1.00252738  0.99835095]Итерация 5: [ 1.00086098  2.00029825 -1.00030728  0.99984975]Итерация 6: [ 1.00009128  2.00002134 -1.00003115  0.9999881 ]Итерация 7: [ 1.00000836  2.00000117 -1.00000275  0.99999922]Итерация 8: [ 1.00000067  2.00000002 -1.00000021  0.99999996]Итерация 9: [ 1.00000004  1.99999999 -1.00000001  1.        ]Итерация 10: [ 1.  2. -1.  1.]Решение: [ 1.  2. -1.  1.]Ошибка: [  2.06480930e-08  -1.25551054e-08   3.61417563e-11   0,00000000e + 00]

Программы для решения произвольной нет. уравнений с помощью Matlab

В следующем коде используется формула ${ displaystyle x_ {i} ^ {(k + 1)} = { frac {1} {a_ {ii}}} left (b_ {i} - sum _ {j i} a_ {ij} x_ {j} ^ {(k)} right), quad i, j = 1,2, ldots , n}$

функцияИкс =gauss_seidel(А, б, х, итерс)за я = 1:iters        за j = 1:размер(А,1)            Икс(j) = (1/А(j,j)) * (б(j) - А(j,:)*Икс + А(j,j)*Икс(j));        конец    конецконец

Смотрите также

Последовательное чрезмерное расслабление
Качмарц метод ("ориентированный на строку" метод, тогда как метод Гаусса-Зейделя "ориентированный на столбец". См., например, Эта бумага.)
Итерационный метод. Линейные системы
Распространение гауссовских убеждений
Расщепление матрицы
Итерация Ричардсона

Примечания

^ Зауэр, Тимоти (2006). Числовой анализ (2-е изд.). Pearson Education, Inc. стр. 109. ISBN 978-0-321-78367-7.
^ Гаусс 1903, п. 279; Прямая ссылка.
^ Голуб и Ван Лоан 1996, п. 511.
^ Голуб и Ван Лоан 1996, уравнение (10.1.3).
^ Голуб и Ван Лоан 1996, Thm 10.1.2.
^ Багнара, Роберто (март 1995). «Единое доказательство сходимости методов Якоби и Гаусса-Зейделя». SIAM Обзор. 37 (1): 93–97. Дои:10.1137/1037008. JSTOR 2132758.

внешняя ссылка

«Метод Зейделя», Энциклопедия математики, EMS Press, 2001 [1994]
Гаусс – Зайдель с сайта www.math-linux.com
Гаусс-Зейдель Из Института целостных численных методов
Gauss Siedel Iteration с сайта www.geocities.com
Метод Гаусса-Зейделя
Биксон
Код Matlab
Пример кода C

[1] Зауэр, Тимоти (2006). Числовой анализ (2-е изд.). Pearson Education, Inc. стр. 109. ISBN 978-0-321-78367-7.

[2] Гаусс 1903, п. 279; Прямая ссылка.

[3] Голуб и Ван Лоан 1996, п. 511.

[4] Голуб и Ван Лоан 1996, уравнение (10.1.3).

[5] Голуб и Ван Лоан 1996, Thm 10.1.2.

[6] Багнара, Роберто (март 1995). «Единое доказательство сходимости методов Якоби и Гаусса-Зейделя». SIAM Обзор. 37 (1): 93–97. Дои:10.1137/1037008. JSTOR 2132758.

[1]

[2]

[3]

[4]

[5]

[6]

Числовая линейная алгебра
Ключевые идеи	Плавающая точка Численная стабильность
Проблемы	Система линейных уравнений Матричные разложения Умножение матриц (алгоритмы ) Расщепление матрицы Редкие проблемы
Аппаратное обеспечение	Кэш процессора TLB Алгоритм без кеширования SIMD Многопроцессорность
Программного обеспечения	MATLAB Подпрограммы базовой линейной алгебры (BLAS) ЛАПАК Специализированные библиотеки Программное обеспечение общего назначения