使用 Python 建模高斯赛德尔法
高斯赛德尔法是求解任何线性方程组的迭代法。虽然该方法与雅可比法非常相似，但在高斯赛德尔法中，迭代中获得的未知数 (x) 的值用于同一次迭代，而在雅可比法中，它们用于下一迭代级别。在同一步中更新x可以加快收敛速度。
线性方程组可以写成−
$$\mathrm{a_{1,1}x_{1} \: + \: a_{1,2}x_{2} \: + \: \dotso \: + \: a_{1,n}x_{n} \: = \: b_{1}}$$
$$\mathrm{a_{2,1}x_{1} \: + \: a_{2,2}x_{2} \: + \: \dotso \: + \: a_{2,n}x_{n} \: = \: b_{2}}$$
$$\mathrm{\vdots}$$
$$\mathrm{a_{n,1}x_{1} \: + \: a_{n,2}x_{2} \: + \: \dotso \: + \: a_{n,n}x_{n} \: = \: b_{n}}$$
因此，一般来说，高斯赛德尔法可以写成 −
$$\mathrm{x_{i_{new}} \: = \: − \frac{1}{a_{i,i}}(\sum_{j=1,j \neq i}^n a_{i,j}x_{j_{new}} \: − \: b_{i})}$$
其中，$\mathrm{x_{j_{new}}}$ 是刚刚在同一迭代循环中求解的方程式中 x 的新值。
高斯赛德尔算法可以总结如下 −
从一些未知数 (x) 的初始猜测数组开始。
通过在重新排列的方程式形式中用猜测值替换 x 来评估新的 x，如下所示 −
$$\mathrm{x_{i_{new}} \: = \: − \frac{1}{a_{i,i}}(\sum_{j=1,j \neq i}^n a_{i,j}x_{j_{new}} \: − \: b_{i})}$$
现在 $\mathrm{i_{new}}$ 将是从当前迭代获得的新 x 值。
下一步是评估新值和猜测值之间的误差，即 $\mathrm{\lvert x_{new} \: − \: x_{guess} \rvert}$。如果误差超过某个收敛标准(我们将其视为 $\mathrm{10^{−5}}$)，则将新值分配给旧猜测，即 $\mathrm{x_{guess} \: = \: x_{new}}$，然后开始下一次迭代。
否则，$\mathrm{x_{new}}$ 是最终答案。
重要的是要记住，方程的对角优势，通常称为 Scarborough 标准，适用于 Gauss-Seidel。如果 −
$$\mathrm{\frac{\sum_{i\neq j \lvert a_{i,j} \rvert}}{\lvert a_{i.i} \rvert}\lbrace \substack{\le \: 1 \: for \: all \: equations \\ \lt \: 1 \: for \: at \: least \: one \: equation}}$$
这就是为什么有时我们必须改变数组中方程的顺序来实现所谓的对角优势的原因。
让我们借助下面提到的例子来演示算法 −
$$\mathrm{20x \: + \: y \: − \: 2z \: = \: 17}$$
$$\mathrm{3x \: + \: 20y \: − \: z \: = \: −18}$$
$$\mathrm{2x \: − \: 3y \: + \: 20z \: = \: 25}$$
将上述方程重新排列如下 −
$$\mathrm{x_{new} \: = \: (−y_{guess} \: + \: 2z_{guess} \: + \: 17)/20}$$
$$\mathrm{y_{new} \: = \: (−3x_{guess} \: + \: z_{guess} \: − \: 18)/20}$$
$$\mathrm{z_{new} \: = \: (−2x_{guess} \: + \: 3y_{guess} \: + \: 25)/20}$$
现在将在 while 循环中求解这些方程，以根据猜测值和当前获得的值获得未知数的新值。
用于模拟高斯赛德尔方法的 Python 程序
用于实现高斯赛德尔方法(方程式实现)的程序如下所示 −
示例
# 导入模块
来自 pylab 导入 *
来自 numpy 导入 *
# 设置初始猜测
xg=0
yg=0
zg=0
#进入循环时出错
error=1
#迭代计数器
count=0
while error>1.E-5:
count+=1
x=(17-yg+2*zg)/20
y=(zg-18-3*x)/20 # 这里 x 是 x_new(仅在此步骤中获得)
z=(25-2*x+3*y)/20 # 这里 x 和 y 是 x_new 和 y_new(仅在此步骤中获得)
# 评估和绘制误差
error = abs(x-xg)+abs(y-yg)+abs(z-zg)
figure(1,dpi=300)
semilogy(count,error,'ko')
xlabel('iterations')
ylabel('error')
# 设置下一次迭代的猜测
xg=x
yg=y
zg=z
savefig('error_GS.jpg')
print(f'x={round(x,5)}, y={round(y,5)}, z={round(z,5)}')
输出
程序输出将是
$$\mathrm{x \: = \: 1.0 \: , \: y \: = \: −1.0 \: , \: z \: = \: 1.0}$$
但需要注意的关键是误差图，很明显，仅经过五步即可达到收敛解。
当方程数量变大时，我们很难在循环中编写方程并求解。因此，为了解决这个困难，系数被捆绑在一个矩阵中，可以轻松地从中插入到程序中。 Gauss Seidel 实现的矩阵方法如下 −
示例
# 导入模块
from pylab import *
from numpy import *
# 系数矩阵
a=array([[20,1,-2],[3,20,-1],[2,-3,20]])
# RHS 向量
b=array([17,-18,25])
# 矩阵中的行数和列数
n=len(b)
# 设置初始猜测
xn=zeros(len(b))
# 设置误差以在循环中移动
error=1
# 设置迭代计数器
count=0
# 将猜测数组复制到新的
xg=xn.copy()
while error>1.E-5:
count+=1
for i in range(n):
sum1=0
for j in range(n):
if i!=j:
sum1=sum1+a[i,j]*xn[j]
xn[i]=(-1/a[i,i])*(sum1-b[i])
# 评估并绘制误差
error = sum(abs(xn-xg))
figure(1,dpi=300)
semilogy(count,error,'ko')
xlabel('iterations')
ylabel('error')
# 更新猜测
xg=xn.copy()
savefig('error_GS1.jpg')
print('x: ',xn)
输出
程序输出将是
$$\mathrm{x \: : \: 0.99999999 \: , \: −1 \: , \: 1}$$
结论
在本教程中，我们尝试对著名的高斯赛德尔方法进行建模。我们介绍了深入的算法来对该方法进行编码。我们还解决了一个示例问题并展示了错误图。
高斯赛德尔方法的成功在于方程的编写方式。如果系统是对角占优的，那么请确保您将获得正确的答案，否则您必须重新排列方程式。
