Метод Якобі

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

Опис методу

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

де:

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

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

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

чи


Збіжність

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

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

Алгоритм

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

 
#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);
}

Див. також

Джерела
This article is issued from Wikipedia. The text is licensed under Creative Commons - Attribution - Sharealike. Additional terms may apply for the media files.