计算机数值方法-雅可比迭代和高斯赛德尔迭代

算法流程

其实迭代法前面已经学习过啦，这里的迭代是在前面迭代的基础上的高阶形式——即解决线性方程组的问题。

下面简单介绍雅克比迭代的基本流程。

雅可比迭代

有一线性方程组，$Ax=b$,其中：

我们可以将其化为以下形式：
$x_i=B{x}_j+f,(i=1,2,3......n,j=1,2,3,\neg i.....n)$
则迭代形式可化为：
$x^i=B{x}^{i+1}+f$
$jacobi$迭代法的流程是：
若系数矩阵$A$是非奇异矩阵且$a_{ii}/\ne 0$，则可以将$A$分裂成:
$A=D+L+U$
其中$D$为对角矩阵，$L$为下三角矩阵，$U$为上三角矩阵
则迭代公式可以转换为：
$x^i=-D^{-1}(L+U)x^{i+1}+f$
整理得：

具此求解.

高斯-赛德尔迭代

在雅可比迭代的流程中我们不难发现

前一步计算出来的$x_i^{k+1}$在下一步中并没有利用到，而新计算出来的值必定比前置更为精确，故为了使计算更为精确，我们将下一步中的$x_i^k$替换为上一步中计算出来的$x_i^{k+1}$进行计算，这种算法就叫做高斯-赛德尔迭代（Gauss-Seidel）
化简得到：

C++代码

雅可比迭代：

#include <bits/stdc++.h>
using namespace std;
typedef pair<int, int> PII;
#define int long long
const int N = 1e3 + 10;
double A[N][N], B[N], X[N];
int n;
void jacobi()
{int k = N;while (k--){double X2[N];for (int i = 0; i < n; i++){double cnt = 0;for (int j = 0; j < n; j++){if (j == i)continue;elsecnt += A[i][j] * X[j];}X2[i] = (B[i] - cnt) / A[i][i];}for(int i= 0; i < n; i++) X[i]=X2[i];}for (int i = 0; i < n; i++)printf("X[%d]=%lf%c", i + 1, X2[i], i == n - 1 ? '\n' : ' ');
}signed main()
{cin >> n;for (int i = 0; i < n; i++)for (int j = 0; j < n; j++)cin >> A[i][j];for (int i = 0; i < n; i++)cin >> B[i];jacobi();return 0;
}

高斯赛德尔迭代：

#include <bits/stdc++.h>
using namespace std;
typedef pair<int, int> PII;
#define int long long
const int N = 1e3 + 10;
double A[N][N], B[N], X[N];
int n;
void gauss_seidel()
{int k = N;while (k--){for (int i = 0; i < n; i++){double cnt = 0;for (int j = 0; j < n; j++){if (j == i)continue;elsecnt += A[i][j] * X[j];}X[i] = (B[i] - cnt) / A[i][i];}}for (int i = 0; i < n; i++)printf("X[%d]=%lf%c", i + 1, X[i], i == n - 1 ? '\n' : ' ');
}signed main()
{cin >> n;for (int i = 0; i < n; i++)for (int j = 0; j < n; j++)cin >> A[i][j];for (int i = 0; i < n; i++)cin >> B[i];gauss_seidel();return 0;
}

python代码

雅可比迭代：

在这里插入代码片

高斯-赛德尔迭代：

在这里插入代码片