Python中是否有实现Gauss和Jacobi线性方程组求解算法的预定义函数？

高斯消元算法 #

• 有现成封装可以直接调用：numpy库的numpy.linalg.solve()底层就是经过数值优化的高斯消元实现，仅支持满秩方阵作为系数矩阵，调用示例：

import numpy as np
# 系数矩阵A，常数项向量b
A = np.array([[2, 1], [1, 2]], dtype=np.float64)
b = np.array([3, 3], dtype=np.float64)
x = np.linalg.solve(A, b)

雅可比迭代算法 #

• 主流数值计算库（numpy、scipy）没有单独封装雅可比迭代的预定义函数，原因是雅可比迭代收敛速度慢，工业场景更常用高斯-赛德尔、SOR、共轭梯度等性能更优的迭代解法，这类优化迭代方法在scipy.sparse.linalg模块有现成实现。
• 如需使用雅可比迭代，可自行编码实现，参考极简实现代码：

import numpy as np

def jacobi_solve(A, b, tol=1e-6, max_iter=1000):
    # 校验矩阵是否为方阵
    if A.shape[0] != A.shape[1] or A.shape[0] != len(b):
        raise ValueError("系数矩阵维度和常数项向量不匹配")
    n = A.shape[0]
    x = np.zeros(n, dtype=np.float64)
    # 提取对角矩阵的逆、L+U矩阵
    D_diag = np.diag(A)
    if np.any(D_diag == 0):
        raise ValueError("系数矩阵对角元素存在0，无法使用雅可比迭代")
    D_inv = np.diag(1 / D_diag)
    L_plus_U = A - np.diag(D_diag)
    # 迭代求解
    for _ in range(max_iter):
        x_new = D_inv @ (b - L_plus_U @ x)
        if np.linalg.norm(x_new - x, ord=np.inf) < tol:
            return x_new
        x = x_new
    raise RuntimeError(f"迭代{max_iter}次后仍未达到精度要求")

• 注意：雅可比迭代仅对严格对角占优或对称正定的系数矩阵收敛，调用前需确认矩阵属性。

内容的提问来源于stack exchange，提问作者Bogdan Doicin