Метод Якоби

Материал из Википедии — свободной энциклопедии
Перейти к: навигация, поиск

Метод Якоби — разновидность метода простой итерации для решения системы линейных алгебраических уравнений. Назван в честь Карла Густава Якоби.

Постановка задачи[править | править вики-текст]

Возьмём систему линейных уравнений:

A\vec{x}=\vec{b}, где A=\left(
\begin{array}{ccc}
a_{11} & \ldots & a_{1n} \\
\vdots & \ddots & \vdots \\
a_{n1} & \ldots & a_{nn} 
\end{array} \right),\quad \vec{b}=\left(
\begin{array}{c}
b_1 \\
\vdots \\
b_n 
\end{array} \right)

Или \left\{
\begin{array}{rcl}
a_{11}x_1 + \ldots + a_{1n}x_n& = & b_{1} \\
\ldots ~~~~~~~~~~~~~~~~~~~~~ \\
a_{n1}x_1 + \ldots + a_{nn}x_n & = & b_{n} 
\end{array} \right.

Описание метода[править | править вики-текст]

Для того, чтобы построить итеративную процедуру метода Якоби, необходимо провести предварительное преобразование системы уравнений A\vec{x}=\vec{b} к итерационному виду \vec{x}=B\vec{x}+\vec{g}. Оно может быть осуществлено по одному из следующих правил:

  • B = E-D^{-1}A = D^{-1}(D - A),\quad \vec{g}=D^{-1}\vec{b};
  • B = -D^{-1}(L + U) = -D^{-1}(A - D),\quad \vec{g}=D^{-1}\vec{b}
  • D^{-1}_{ii} = 1 / D_{ii}, D_{ii} \neq 0,\, i = 1,2, ..., n\quad;


где в принятых обозначениях D означает матрицу, у которой на главной диагонали стоят соответствующие элементы матрицы A, а все остальные нули; тогда как матрицы U и L содержат верхнюю и нижнюю треугольные части A, на главной диагонали которых нули, E — единичная матрица.

Тогда процедура нахождения решения имеет вид:

 \vec{x}^{(k+1)}  = B \vec{x}^{(k)}+\vec{g},

Или в виде поэлементной формулы:

 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.

где k счётчик итерации.

В отличие от метода Гаусса-Зейделя мы не можем заменять x^{(k)}_i\, на x^{(k+1)}_i в процессе итерационной процедуры, так как эти значения понадобятся для остальных вычислений. Это наиболее значимое различие между методом Якоби и методом Гаусса-Зейделя решения СЛАУ. Таким образом на каждой итерации придётся хранить оба вектора приближений: старый и новый.

Условие сходимости[править | править вики-текст]

Приведем достаточное условие сходимости метода.

Logo arte.jpg Теорема.
Пусть \| B \| < 1\!. Тогда при любом выборе начального приближения \vec{x}^{(0)}\!:
  1. метод сходится;
  2. скорость сходимости метода равна скорости сходимости геометрической прогрессии со знаменателем q= \|B\|\!;
  3. верна оценка погрешности: \|\vec{x}^{(k)}-\vec{x}\| < q^k \, \|\vec{x}^{(0)}-\vec{x}\|\!.

Условие остановки[править | править вики-текст]

Условие окончания итерационного процесса при достижении точности \varepsilon в упрощённой форме имеет вид:

max_j\|x_j^{(k+1)}-x_j^{(k)}\|/(1-q) < \varepsilon\!

(Существует более точное условие окончания итерационного процесса, которое более сложно и требует дополнительных вычислений)[источник не указан 1341 день]


Алгоритм[править | править вики-текст]

Ниже приведён алгоритм реализации на C++

#include <math.h>
const double eps = 0.001; ///< желаемая точность 
 
..........................
 
/// N - размерность матрицы; A[N][N] - матрица коэффициентов, F[N] - столбец свободных членов,
/// X[N] - начальное приближение, ответ записывается также в X[N];
void Jacobi (int N, double** A, double* F, double* X)
{
	double* TempX = new double[N];
	double norm; // норма, определяемая как наибольшая разность компонент столбца иксов соседних итераций.
 
	do {
		for (int i = 0; i < N; i++) {
			TempX[i] = F[i];
			for (int g = 0; g < N; g++) {
				if (i != g)
					TempX[i] -= A[i][g] * X[g];
			}
			TempX[i] /= A[i][i];
		}
                norm = fabs(X[0] - TempX[0]);
		for (int h = 0; h < N; h++) {
			if (fabs(X[h] - TempX[h]) > norm)
				norm = fabs(X[h] - TempX[h]);
			X[h] = TempX[h];
		}
	} while (norm > eps);
	delete[] TempX;
}

См. также[править | править вики-текст]