jacobi迭代-python_python实现jacobi迭代

相关公式：（老师的ppt）

代码：

import numpy as np
#输入A、b
n=3
A=np.array([[10,-2,-1],
            [-2,10,-1],
            [-1,-2,5]],dtype=np.float_)
b=np.array([[3,15,10]]).T
 
#初始化xk+1，xk
xk=np.zeros((3,1),dtype=np.float_)
xk1=np.zeros((3,1),dtype=np.float_)
 
 
 
#生成D、L、U
D=np.diag(np.diag(A))
#  D=[[10,0,0],
#     [0,10,0],
#     [0,0,5]]
 
L=-np.tril(A,k=-1)
#  L=[[0,0,0],
#     [-1,0,0],
#     [-1,-1,0]]
 
U=-np.triu(A,k=1)
#  U=[[0,-1,-2],
#     [0,0,-2],
#     [0,0,0]]
 
 
#计算M、g
M=np.linalg.inv(D)@(L+U)
g=np.linalg.inv(D)@b
 
#迭代
for i in range(10):
    print(xk)
    xk1=M@xk+g
    xk=xk1
 
#结果
print(xk.T)