Метод Якобі

Метод Якобі — класичний ітераційний метод розв'язку системи лінійних рівнянь

Опис методу

Для квадратної системи з n лінійних рівнянь:

${\displaystyle A\mathbf {x} =\mathbf {b} }$

де:

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

Матрицю A можна розкласти на два доданки: діагональну матрицю D, та все інше R:

${\displaystyle A=D+R,\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}},\qquad R={\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}}.}$

Систему лінійних рівнянь можна переписати в вигляді:

${\displaystyle D\mathbf {x} =\mathbf {b} -R\mathbf {x} }$

Ітераційний метод Якобі виражається формулою:

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

чи

${\displaystyle 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.}$

Збіжність

Метод є збіжним, коли матриця A має домінантну головну діагональ:

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

другою умовою збіжності є, те щоб не перевищував одиницю:

${\displaystyle \rho (D^{-1}R)<1.}$

Алгоритм

Реалізація на С++

#include <vector> #include <iostream> #include <cmath> using namespace std; void solve(const vector<float> a, //квадратна матриця  const vector<float> b, //вектор вільних елементів  vector<float>& x, //сюди буде записано розв'язок  const float allowed_error) //допустима похибка {  const unsigned n = x.size();  vector<float> tmp_x(n);    float error;    do  {  error = 0;    tmp_x = b;  for (unsigned i = 0; i < n; ++i)  {  for (unsigned j = 0; j < n; ++j)  {  if (i != j)  {  tmp_x[i] -= a[i * n + j] * x[j];  }  }    //оновити x[i] та порахувати похибку  const float x_updated = tmp_x[i] / a[i * (n + 1)];  const float e = fabs(x[i] - x_updated);  x[i] = x_updated;  if (e > error) { error = e; }  }  }  while (error > allowed_error); } //приклад використання. Користувач вводить вхідні дані. int main() {  unsigned n;    cout << "\nВведіть розмір системи:\nn = ";  cin >> n;    vector<float> x(n);  vector<float> a(n * n);  vector<float> b(n);    cout << "\nВведіть вектор вільних елементів: \n";  for (auto& b_elem : b)  {  cin >> b_elem;  }    cout << "\nВведіть коефіцієнти системи: \n";  for (auto& a_elem : a)  {  cin >> a_elem;  }    float allowed_error;  cout << "\nВведіть допустиму похибку: \n";  cin >> allowed_error;    solve(a, b, x, allowed_error);    cout << "\nРозв'язок системи:: \n";  for (unsigned i = 0; i < n; ++i)  {  cout << "\nx[" << i << "]= " << x[i];  } } 

Реалізація на С

#include<stdio.h> #include<stdlib.h> #include<conio.h> #include<math.h> void main(){ int n, i, j, count = 0; float **a, *b, *x, *tmp_x, exp, e; printf("Введіть розміри системи:\n"); printf("n = "); scanf("%i", &n); a = (float**)malloc(n*sizeof(float*)); for(i = 0; i < n; i++){ a[i] = (float*)malloc(n*sizeof(float)); } b = (float*)malloc(n*sizeof(float)); x = (float*)malloc(n*sizeof(float));  tmp_x = (float*)malloc(n*sizeof(float)); printf("Введіть коефіцієнти системи:\n"); for(i = 0; i < n; i++){ for(j = 0; j < n; j++){  scanf("%f", &a[i][j]); } } printf("Введіть вектор вільних елементів:\n"); for(i = 0; i < n; i++){ scanf("%f", &b[i]);  x[i] = 0; } printf("Введіть точність обчислення: "); scanf("%f", &e);   do{ count++; for(i = 0; i < n; i++){ tmp_x[i] = 0.0; for(j = 0; j < n; j++){ if(i != j){ tmp_x[i] = tmp_x[i] + (a[i][j] * x[j]); } } tmp_x[i] = (b[i] - tmp_x[i]) / a[i][i]; } exp = 0; for(i = 0; i < n; i++){ if(fabs(x[i] - tmp_x[i]) > exp){ exp = fabs(x[i] - tmp_x[i]); } x[i] = tmp_x[i]; } }while(exp > e);  free(tmp_x); for(i = 0; i < n; i++){ free(a[i]); } free(a); free(b);  printf("Розв'язок системи:\n"); for(i = 0; i < n; i++){ printf("x[%d] = %.6f\n", i+1, x[i]); } free(x); } 

Див. також

• Ітераційні методи розв'язування СЛАР
• Метод Гауса — Зейделя

Джерела

• Гантмахер Ф. Р. Теория матриц. — 5-е. — М: : Физматлит, 2010. — 559 с. — .(рос.)