метод Якоби

Не путать с алгоритмом собственных значений Якоби .

В числовой линейной алгебре метод Якоби (он же итерационный метод Якоби ) — итерационный алгоритм определения решений строго диагонально-доминантной системы линейных уравнений . Определяется каждый диагональный элемент и подставляется приблизительное значение. Затем процесс повторяется до тех пор, пока он не сойдётся. Этот алгоритм представляет собой урезанную версию метода диагонализации матрицы преобразования Якоби . Метод назван в честь Карла Густава Якоба Якоби .

Описание

Позволять $A\mathbf {x} =\mathbf {b}$ — квадратная система n линейных уравнений, где:

$$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}}.$$

Когда $A$ и $\mathbf {b}$ известны, и $\mathbf {x}$ неизвестно, мы можем использовать метод Якоби для аппроксимации $\mathbf {x}$. Вектор $\mathbf {x} ^{(0)}$ обозначает наше первоначальное предположение для $\mathbf {x}$ (часто $\mathbf {x} _{i}^{(0)}=0$ для $i=1,2,...,n$). Обозначим $\mathbf {x} ^{(k)}$ как k- е приближение или итерация $\mathbf {x}$, и $\mathbf {x} ^{(k+1)}$ это следующая (или k +1) итерация $\mathbf {x}$.

Формула на основе матрицы

Тогда A можно разложить на диагональную составляющую D , нижнюю треугольную часть L и верхнюю треугольную часть U :

$$A=D+L+U\qquad {\text{where}}\qquad D={\begin{bmatrix}a_{11}&0&\cdots &0\\0&a_{22}&\cdots &0\\\vdots &\vdots &\ddots &\vdots \\0&0&\cdots &a_{nn}\end{bmatrix}}{\text{ and }}L+U={\begin{bmatrix}0&a_{12}&\cdots &a_{1n}\\a_{21}&0&\cdots &a_{2n}\\\vdots &\vdots &\ddots &\vdots \\a_{n1}&a_{n2}&\cdots &0\end{bmatrix}}.$$

Затем решение получается итерационно через

$\mathbf {x} ^{(k+1)}=D^{-1}(\mathbf {b} -(L+U)\mathbf {x} ^{(k)}).$

Формула на основе элементов

Формула на основе элементов для каждой строки $i$ таким образом:

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

Вычисление $x_{i}^{(k+1)}$ требует каждого элемента в $\mathbf {x} ^{(k)}$кроме себя. В отличие от метода Гаусса–Зейделя , мы не можем перезаписать $x_{i}^{(k)}$ с $x_{i}^{(k+1)}$, так как это значение понадобится для остальной части вычислений. Минимальный объем памяти — два вектора размера n .

Алгоритм

Входные данные:   первоначальное предположение  x (0) к решению  , (диагональная доминанта) матрица  A  , вектор правой части  b  , критерий сходимости 
  Выходные данные:   решение при достижении сходимости. 
 Комментарии:  псевдокод, основанный на приведенной выше формуле на основе элементов. 

  k  = 0 
 , пока  сходимость не достигнута,  сделать 
     для   i  := 1  шаг, пока  n  сделать 
         σ  = 0 
         для   j  := 1  шаг, пока  n  сделать 
             , если   j  ≠  i   , то 
                 σ  =  σ  +  a  ij   x  j (  к  ) 
            конец 
         конец 
         х  я (  к  +1) знак равно (  б  я  -  σ  ) /  а  ii 
     конец 
      приращение  k 
 конец 

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

Стандартное условие сходимости (для любого итерационного метода) — это когда спектральный радиус итерационной матрицы меньше 1:

$\rho (D^{-1}(L+U))<1.$

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

$\left|a_{ii}\right|>\sum _{j\neq i}{\left|a_{ij}\right|}.$

Метод Якоби иногда сходится, даже если эти условия не выполняются.

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

$$A={\begin{pmatrix}29&2&1\\2&6&1\\1&1&{\frac {1}{5}}\end{pmatrix}}\quad \Rightarrow \quad D^{-1}(L+U)={\begin{pmatrix}0&{\frac {2}{29}}&{\frac {1}{29}}\\{\frac {1}{3}}&0&{\frac {1}{6}}\\5&5&0\end{pmatrix}}\quad \Rightarrow \quad \rho (D^{-1}(L+U))\approx 1.0661\,.$$

Примеры

Пример вопроса

Линейная система вида $Ax=b$ с первоначальной оценкой $x^{(0)}$ дан кем-то

$A={\begin{bmatrix}2&1\\5&7\\\end{bmatrix}},\ b={\begin{bmatrix}11\\13\\\end{bmatrix}}\quad {\text{and}}\quad x^{(0)}={\begin{bmatrix}1\\1\\\end{bmatrix}}.$

Мы используем уравнение $x^{(k+1)}=D^{-1}(b-(L+U)x^{(k)})$, описанный выше, для оценки $x$. Сначала перепишем уравнение в более удобном виде $D^{-1}(b-(L+U)x^{(k)})=Tx^{(k)}+C$, где $T=-D^{-1}(L+U)$ и $C=D^{-1}b$. Из известных значений

$$D^{-1}={\begin{bmatrix}1/2&0\\0&1/7\\\end{bmatrix}},\ L={\begin{bmatrix}0&0\\5&0\\\end{bmatrix}}\quad {\text{and}}\quad U={\begin{bmatrix}0&1\\0&0\\\end{bmatrix}}.$$

мы определяем $T=-D^{-1}(L+U)$ как

$$T={\begin{bmatrix}1/2&0\\0&1/7\\\end{bmatrix}}\left\{{\begin{bmatrix}0&0\\-5&0\\\end{bmatrix}}+{\begin{bmatrix}0&-1\\0&0\\\end{bmatrix}}\right\}={\begin{bmatrix}0&-1/2\\-5/7&0\\\end{bmatrix}}.$$

Дальше, $C$ находится как

$$C={\begin{bmatrix}1/2&0\\0&1/7\\\end{bmatrix}}{\begin{bmatrix}11\\13\\\end{bmatrix}}={\begin{bmatrix}11/2\\13/7\\\end{bmatrix}}.$$

С $T$ и $C$ рассчитано, мы оцениваем $x$ как $x^{(1)}=Tx^{(0)}+C$:

$$x^{(1)}={\begin{bmatrix}0&-1/2\\-5/7&0\\\end{bmatrix}}{\begin{bmatrix}1\\1\\\end{bmatrix}}+{\begin{bmatrix}11/2\\13/7\\\end{bmatrix}}={\begin{bmatrix}5.0\\8/7\\\end{bmatrix}}\approx {\begin{bmatrix}5\\1.143\\\end{bmatrix}}.$$

Следующая итерация дает

$$x^{(2)}={\begin{bmatrix}0&-1/2\\-5/7&0\\\end{bmatrix}}{\begin{bmatrix}5.0\\8/7\\\end{bmatrix}}+{\begin{bmatrix}11/2\\13/7\\\end{bmatrix}}={\begin{bmatrix}69/14\\-12/7\\\end{bmatrix}}\approx {\begin{bmatrix}4.929\\-1.714\\\end{bmatrix}}.$$

Этот процесс повторяется до сходимости (т.е. до тех пор, пока $\|Ax^{(n)}-b\|$маленький). Решение после 25 итераций:

$x={\begin{bmatrix}7.111\\-3.222\end{bmatrix}}.$

Пример вопроса 2

Предположим, нам дана следующая линейная система:

${\begin{aligned}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{aligned}}$

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

$${\begin{aligned}x_{1}&=(6+0-(2*0))/10=0.6,\\x_{2}&=(25+0+0-(3*0))/11=25/11=2.2727,\\x_{3}&=(-11-(2*0)+0+0)/10=-1.1,\\x_{4}&=(15-(3*0)+0)/8=1.875.\end{aligned}}$$

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

$x_{1}$	$x_{2}$	$x_{3}$	$x_{4}$
0.6	2.27272	-1.1	1.875
1.04727	1.7159	-0.80522	0.88522
0.93263	2.05330	-1.0493	1.13088
1.01519	1.95369	-0.9681	0.97384
0.98899	2.0114	-1.0102	1.02135

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

Пример Python

import   numpy   as   np 

 ITERATION_LIMIT   =   1000 

 # инициализировать матрицу 
 A   =   np  .   массив  ([[  10.  ,   -  1.  ,   2.  ,   0.  ], 
               [  -  1.  ,   11.  ,   -  1.  ,   3.  ], 
               [  2.  ,   -  1.  ,   10.  ,   -  1.  ], 
               [  0.0  ,   3.  ,   -  1.  ,   8.  ]]) 
 # инициализируем вектор RHS 
 b   =   np  .   array  ([  6.  ,   25.  ,   -  11.  ,   15.  ]) 

 # печатает системную 
 печать  (  "System:"  ) 
 для   i   в   диапазоне  (  A  .  shape  [  0  ]): 
     row   =   [  f  "  {  A  [  i  ,   j  ]  }  *x  {  j   +   1  }  "   for   j   в   диапазоне  (  A  .  shape  [  1  ])] 
     print  (  f  '  {  " + "  .  join  (  row  )  }  =  {  b  [  i  ]  }  '  ) 
 print  ( ) 

 Икс знак   равно   np  .   нули_подобно  (  b  ) 
 для   it_count   в   диапазоне  (  ITERATION_LIMIT  ): 
     если   it_count   !=   0  : 
         напечатайте  (  f  "Итерация  {  it_count  }  :  {  x  }  "  ) 
     x_new   =   np  .   Zeros_like  (  x  ) 

     для   i   в   диапазоне  (  A  .  shape  [  0  ]): 
         s1   =   np  .   точка  (  А  [  я  ,   :  я  ],   Икс  [:  я  ]) 
         s2   знак равно   np  .   точка  (  A  [  i  ,   i   +   1  :],   x  [  i   +   1  :]) 
         x_new  [  i  ]   =   (  b  [  i  ]   -   s1   -   s2  )   /   A  [  i  ,   i  ] 
         if   x_new  [  i  ]   ==   x_new  [  i  -  1  ]: 
           сломать, 

     если   np  .   allclose  (  x  ,   x_new  ,   atol  =  1e-10  ,   rtol  =  0.  ): 
         перерыв 

     x   =   x_new 

 print  (  «Решение:»  ) 
 print  (  x  ) 
 error   =   np  .   точка  (  A  ,   x  )   -   b 
 печать  (  «Ошибка:»  ) 
 печать  (  ошибка  ) 

Взвешенный метод Якоби

Взвешенная итерация Якоби использует параметр $\omega$ вычислить итерацию как

$\mathbf {x} ^{(k+1)}=\omega D^{-1}(\mathbf {b} -(L+U)\mathbf {x} ^{(k)})+\left(1-\omega \right)\mathbf {x} ^{(k)}$

с $\omega =2/3$ это обычный выбор. Из отношения $L+U=A-D$, это также может быть выражено как

$\mathbf {x} ^{(k+1)}=\omega D^{-1}\mathbf {b} +\left(I-\omega D^{-1}A\right)\mathbf {x} ^{(k)}$.

в симметричном положительно

В случае, если системная матрица $A$ имеет симметричный положительно определенный тип, можно показать сходимость.

Позволять $C=C_{\omega }=I-\omega D^{-1}A$быть итерационной матрицей. Тогда сходимость гарантирована для

$\rho (C_{\omega })<1\quad \Longleftrightarrow \quad 0<\omega <{\frac {2}{\lambda _{\text{max}}(D^{-1}A)}}\,,$

где $\lambda _{\text{max}}$ является максимальным собственным значением.

Спектральный радиус может быть минимизирован при конкретном выборе $\omega =\omega _{\text{opt}}$ следующее

$$\min _{\omega }\rho (C_{\omega })=\rho (C_{\omega _{\text{opt}}})=1-{\frac {2}{\kappa (D^{-1}A)+1}}\quad {\text{for}}\quad \omega _{\text{opt}}:={\frac {2}{\lambda _{\text{min}}(D^{-1}A)+\lambda _{\text{max}}(D^{-1}A)}}\,,$$

где $\kappa$ – число обусловленности матрицы .

См. также

• Метод Гаусса – Зейделя
• Последовательное чрезмерное расслабление
• Итерационный метод § Линейные системы
• Распространение гауссовских убеждений
• Расщепление матрицы

Внешние ссылки

• Эта статья включает текст из статьи Jacobi_method на CFD-Wiki , которая находится под лицензией GFDL .

• Блэк, Ноэль; Мур, Ширли и Вайсштейн, Эрик В. «Метод Якоби» . Математический мир .
• Метод Якоби с сайта www.math-linux.com.