Python——雅克比迭代求线性方程组的根

Python——高斯赛德尔迭代求线性方程组的根

雅克比迭代

计算代码

矩阵A——方程组的系数。代码中只需改变矩阵系数即可
矩阵B——方程组等号右侧的常数。

#雅克比迭代
#求解方程组：
#10*x1 - 2*x2 - 1*x3 = 3
#-2*x1 + 10*x2 - 1*x3 = 15
#-1*x1 - 2*x2 + 5*x3 = 10
#AX = B

#只用改变矩阵A、B的系数即可
A = [[10, -2, -1], [-2, 10, -1], [-1, -2, 5]]
B = [[3], [15], [10]]
#初始值
x = [0, 0, 0]
aij = None
n = 0
while True:
    x_new_list = []
    n += 1
    print(f'第{n}次迭代')
    for i in range(len(A)):
        a = A[i]
        ax = 0
        for j in range(len(a)):
            if j == i :
                aij = a[j]
            else:
                ax += a[j]*x[j]
        x_new = (B[i][0] - ax) / aij
        x_new_list.append(x_new)

    v_max = 0
    for m in range(len(x)):
        v = abs(x[m] - x_new_list[m])
        v_max = max(v_max, v)
    if v_max <= 10**-3:
        print(x)
        print(v_max)
        break
    else:
        x = x_new_list
        print(x)
        print(v_max)

结果：

第1次迭代
[0.3, 1.5, 2.0]
2.0
第2次迭代
[0.8, 1.7600000000000002, 2.66]
0.6600000000000001
第3次迭代
[0.9179999999999999, 1.9259999999999997, 2.864]
0.20399999999999974
第4次迭代
[0.9715999999999999, 1.97, 2.9539999999999997]
0.08999999999999986
第5次迭代
[0.9894000000000001, 1.9897199999999997, 2.98232]
0.028320000000000345
第6次迭代
[0.996176, 1.996112, 2.993768]
0.011448000000000125
第7次迭代
[0.9985991999999999, 1.998612, 2.99768]
0.0039119999999996935
第8次迭代
[0.9994904, 1.9994878399999998, 2.99916464]
0.001484640000000148
第9次迭代
[0.9994904, 1.9994878399999998, 2.99916464]
0.0005285759999997808

雅克比迭代收敛的充分必要条件

L——矩阵A的下三角矩阵
U——矩阵A的上三角矩阵
D——矩阵A的对角矩阵
例如：

from numpy import *
# numpy.linalg模块包含线性代数的函数。使用这个模块，可以计算逆矩阵、求特征值、解线性方程组以及求解行列式等
from numpy.linalg import *
A = mat([[25, 5, 1], [64, 8, 1], [144, 12, 1]])
B = mat([[106.8], [177.2], [279.2]])
L = tril(A, -1) #下三角矩阵
U = triu(A,1) #上三角矩阵
D = diag(diag(A))  #对角矩阵

J = eye(3) - dot(inv(D), A)
print(eigvals(J))  #特征值

结果：

[-3.33072758  2.4         0.93072758]

特征值中只要有一个绝对值＞1，用雅克比迭代求解就不收敛
所以上面的方程组用雅克比迭代不收敛

先判断是否收敛再迭代的综合python代码

#雅克比迭代
#求解方程组：
#10*x1 - 2*x2 - 1*x3 = 3
#-2*x1 + 10*x2 - 1*x3 = 15
#-1*x1 - 2*x2 + 5*x3 = 10
#AX = B
from numpy import *
# numpy.linalg模块包含线性代数的函数。使用这个模块，可以计算逆矩阵、求特征值、解线性方程组以及求解行列式等
from numpy.linalg import *
#只用改变矩阵A、B的系数即可
A = [[10, -2, -1], [-2, 10, -1], [-1, -2, 5]]
B = [[3], [15], [10]]
#先判断此方程组能否用雅克比迭代
L = tril(A, -1) #下三角矩阵
U = triu(A,1) #上三角矩阵
D = diag(diag(A))  #对角矩阵

J = eye(3) - dot(inv(D), A)
root = eigvals(J)   #特征值
flag = True
for r in root:
    if r > 1:
        print('雅克比迭代不收敛')
        flag = False
        break
    else:
        pass

#如果可以迭代，进行迭代求值
if flag:
    #初始值
    x = [0, 0, 0]
    aij = None
    n = 0
    while True:
        x_new_list = []
        n += 1
        print(f'第{n}次迭代')
        for i in range(len(A)):
            a = A[i]
            ax = 0
            for j in range(len(a)):
                if j == i :
                    aij = a[j]
                else:
                    ax += a[j]*x[j]
            x_new = (B[i][0] - ax) / aij
            x_new_list.append(x_new)

        v_max = 0
        for m in range(len(x)):
            v = abs(x[m] - x_new_list[m])
            v_max = max(v_max, v)
        if v_max <= 10**-3:
            print(x)
            print(v_max)
            break
        else:
            x = x_new_list
            print(x)
            print(v_max)

结果：

第1次迭代
[0.3, 1.5, 2.0]
2.0
第2次迭代
[0.8, 1.7600000000000002, 2.66]
0.6600000000000001
第3次迭代
[0.9179999999999999, 1.9259999999999997, 2.864]
0.20399999999999974
第4次迭代
[0.9715999999999999, 1.97, 2.9539999999999997]
0.08999999999999986
第5次迭代
[0.9894000000000001, 1.9897199999999997, 2.98232]
0.028320000000000345
第6次迭代
[0.996176, 1.996112, 2.993768]
0.011448000000000125
第7次迭代
[0.9985991999999999, 1.998612, 2.99768]
0.0039119999999996935
第8次迭代
[0.9994904, 1.9994878399999998, 2.99916464]
0.001484640000000148
第9次迭代
[0.9994904, 1.9994878399999998, 2.99916464]
0.0005285759999997808