首页 模拟笔记

# Jacobi 迭代和 Gauss-Seidel 迭代

# Jacobi 迭代

先是线性方程组Ax=bAx=b, 其解我们可以写成方程组a11x1+a12x2+…a1nxn=b1a_{11}x_1+a_{12}x_2+\ldots a_{1n}x_n=b_1, 其中xx 是未知数，A 和 b 是已知的数，A 和 b 的行列数是一样的。

# 例题

解方程组

{8x1−3x2+2x3=204x1+11x2−x3=336x1+3x2+12x3=36\begin{cases}{8x_1-3x_2+2x_3=20}\\{4x_1+11x_2-x_3=33}\\{6x_1+3x_2+12x_3=36}&\end{cases}

其中

A=[8−32411−16312]x=[x1x2x3]b=[203336]A=\begin{bmatrix}8&-3&2\\4&11&-1\\6&3&12\end{bmatrix}\quad x=\begin{bmatrix}x_1\\x_2\\x_3\end{bmatrix}\quad b=\begin{bmatrix}20\\33\\36\end{bmatrix}

精准解为

x∗=(3,2,1)T\mathbf{x}^{*}=(\mathbf{3},\mathbf{2},\mathbf{1})^{\mathrm{T}}

可以将原方程写为

{x1=(3x2−2x3+20)/8x2=(−4x1+x3+33)/11x3=(−6x1−3x2+36)/12\begin{cases}x_1=(3x_2-2x_3+20)/8\\x_2=(-4x_1+x_3+33)/11\\x_3=(-6x_1-3x_2+36)/12&\end{cases}

那么他的迭代公式是

{x1(k+1)=(3x2(k)−2x3(k)+20)/8x2(k+1)=(−4x1(k)+x3(k)+33)/11x3(k+1)=(−6x1(k)−3x2(k)+36)/12\begin{cases}x_1^{(k+1)}=(3x_2^{(k)}-2x_3^{(k)}+20)/8\\x_2^{(k+1)}=(-4x_1^{(k)}+x_3^{(k)}+33)/11\\x_3^{(k+1)}=(-6x_1^{(k)}-3x_2^{(k)}+36)/12&\end{cases}

那么 jocobi 迭代方程组一般形式为

{a11x1+a12x2+...+a1nxn=b1a21x1+a22x2+...+a2nxn=b2...an1x1+an2x2+...+annxn=bn\begin{cases}{a_{11}x_1+a_{12}x_2+...+a_{1n}x_n=b_1}\\{a_{21}x_1+a_{22}x_2+...+a_{2n}x_n=b_2}\\...\\{a_{n1}x_1+a_{n2}x_2+...+a_{nn}x_n=b_n}&\end{cases}

变形为

{x1=(−a12x2−...−a1nxn+b1)/a11x2=(−a21x1−...−a2nxn+b2)/a22...xn=(−an1x1−an2x2−...−an,n−1xn+bn)/ann\begin{cases}x_1=\quad(-a_{12}x_2-...-a_{1n}x_n+b_1)/a_{11}\\x_2=(-a_{21}x_1\quad-...-a_{2n}x_n+b_2)/a_{22}\\...\\x_n=(-a_{n1}x_1-a_{n2}x_2-...-a_{n,n-1}x_n+b_n)/a_{nn}&\end{cases}

那么迭代方程就是

{x(0)=(x1(0),x2(0),...,xn(0))T(初始向量)xi(k+1)=1aii(bi−∑i=1naijxj(k))\color{red}\begin{cases}x^{(0)}=(x_1^{(0)},x_2^{(0)},...,x_n^{(0)})^T(\text{初始向量})\\x_i^{(k+1)}=\frac1{a_{ii}}(b_i-\sum_{i=1}^na_{ij}x_j^{(k)})&\end{cases}

对于矩阵来说，A 可以分解为A=D+L+UA = D+L+U，那么方程就是

Dx=−(L+U)x+bx=−D−1(L+U)x+D−1bx=Bx+d\begin{gathered}\mathbf{Dx=-(L+U)x+b} \\\mathbf{x=-D^{-1}(L+U)x+D^{-1}b} \\\mathbf{x}=\mathbf{B}\mathbf{x}+\mathbf{d} \end{gathered}

""
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
#include<bits/stdc++.h>
#define MAXSIZE 100
using namespace std;
int main() {
	double A[MAXSIZE][MAXSIZE], x[MAXSIZE], b[MAXSIZE],re[MAXSIZE];
	int n;
	double e;
	cout << " 请输入原方程的阶数：";
	cin >> n;
	cout << "请输入原方程的增广矩阵：";
	for (int i = 0; i < n; ++i) {
		for (int j = 0; j < n ; ++j) {
			cin >> A[i][j];
		}
		cin >> b[i];
	}
	cout << "请输入初始迭代向量值";
	for (int i = 0; i < n; ++i) {
		cin >> x[i];
	}
	cout << "请输入误差上限";
	cin >> e;
	int count = 0;
	while (true) {
		cout << "第" << ++count << "次迭代：";
		int flag = 0;
		for (int i = 0; i < n; ++i) {
			re[i] = 0;
			for (int j = 0; j < n; ++j) {
				if (i != j) {
					re[i] += - A[i][j] * x[j];
				}
			}
			re[i] = (re[i] + b[i]) / A[i][i];
			if (fabs(x[i] - re[i]) < e)
				++flag;
			cout << re[i]<<" ";
		}
		cout << endl;
		if (flag == n) break;
		for (int i = 0; i < n; ++i) {
			x[i] = re[i];
		}
	}
	for (int i = 0; i < n; ++i) {
		cout << re[i] << " ";
	}
}

# 高斯 - 赛德尔迭代法

G-S 迭代就是在 Jocobi 上作为改进，公式是

{x(0)=(x1(0),x2(0),...,xn(0))T(初始向量)xi(k+1)=1aii(bi−∑j=1i−1aijxj(k+1)−∑j=i+1naijxj(k))\color{red}\begin{cases}x^{(0)}=(x_1^{(0)},x_2^{(0)},...,x_n^{(0)})^T(\text{初始向量})\\x_i^{(k+1)}=\dfrac{1}{a_{ii}}(b_i-\sum_{j=1}^{i-1}a_{ij}x_j^{(k+1)}-\sum_{j=i+1}^{n}a_{ij}x_j^{(k)})\end{cases}

而矩阵形式就是

Dx=−(L+U)x+bDx(k+1)=−Lx(k+1)−Ux(k)+bx(k+1)=−(D+L)−1UX(k)+(D+L)−1b\begin{gathered}\mathbf{Dx=-(L+U)x+b} \\\mathbf{Dx^{(k+1)}=-Lx^{(k+1)}-Ux^{(k)}+b} \\{x^{(k+1)}=-(D+L)^{-1}U_X^{(k)}+(D+L)^{-1}b}\end{gathered}

{x(0)=(x1(0),x2(0),...,xn(0))T(初始向量)x(k+1)=−(D+L)−1Ux(k)+(D+L)−1b\begin{cases}\mathbf{x^{(0)}=(x_1^{(0)},x_2^{(0)},...,x_n^{(0)})^T(\text{初始向量})}\\\mathbf{x^{(k+1)}=-(D+L)^{-1}Ux^{(k)}+(D+L)^{-1}b}&\end{cases}

""
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
#include<bits/stdc++.h>
#define MAXSIZE 100
using namespace std;
int main() {
	double A[MAXSIZE][MAXSIZE], x[MAXSIZE], b[MAXSIZE];
	int n;
	double e;
	cout << " 请输入原方程的阶数：";
	cin >> n;
	cout << "请输入原方程的增广矩阵：";
	for (int i = 0; i < n; ++i) {
		for (int j = 0; j < n ; ++j) {
			cin >> A[i][j];
		}
		cin >> b[i];
	}
	cout << "请输入初始迭代向量值";
	for (int i = 0; i < n; ++i) {
		cin >> x[i];
	}
	cout << "请输入误差上限";
	cin >> e;
	int count = 0;
	while (true) {
		cout << "第" << ++count << "次迭代：";
		int flag = 0;
		for (int i = 0; i < n; ++i) {
			double tmp = x[i];
			x[i] = 0;
			for (int j = 0; j < n; ++j) {
				if (i != j) {
					x[i] += - A[i][j] * x[j];
				}
			}
			x[i] = (x[i] + b[i]) / A[i][i];
			if (fabs(x[i] - tmp) < e)
				++flag;
			cout << x[i]<<" ";
		}
		cout << endl;
		if (flag == n) break;
	}
	for (int i = 0; i < n; ++i) {
		cout << x[i] << " ";
	}
}

# 牛顿迭代法

Untitled