Метод Якобі

Матеріал з Вікіпедії — вільної енциклопедії.
Перейти до: навігація, пошук

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

Опис методу[ред.ред. код]

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

A\mathbf x = \mathbf b

де:

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:

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

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

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

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

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

чи

 x^{(k+1)}_i  = \frac{1}{a_{ii}} \left(b_i -\sum_{j\ne i}a_{ij}x^{(k)}_j\right),\quad i=1,2,\ldots,n.


Збіжність[ред.ред. код]

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

\left | a_{ii} \right | > \sum_{i \ne j} {\left | a_{ij} \right |}.

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

\rho(D^{-1}(L+U)) < 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Введіть розмір сиcтеми:\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("Введіть розміри сиcтеми:\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);
}

Див. також[ред.ред. код]

Джерела[ред.ред. код]