如何让SymPy求解线性丢番图方程时输出系数更小的通解

SymPy默认返回的丢番图方程通解系数过大，本质原因是其参数化使用的解空间基没有做整数格约化，你可以通过LLL基约化的方式自动得到小系数的通解，具体实现如下：

完整可运行代码 #

import sympy
from sympy import Matrix
from sympy.solvers.diophantine import diophantine 

# 定义变量与方程系数
w, x, y, z = sympy.symbols('w x y z')
a, b, c, d = -118, 989, 918, -512

# 求解原始丢番图方程
raw_sol = list(diophantine(a*w + b*x + c*y + d*z))[0]
params = sorted(raw_sol.free_symbols, key=str)

# 提取解空间的基向量
basis = []
for p in params:
    coeff_vec = []
    for expr in raw_sol:
        coeff_vec.append(expr.expand().coeff(p))
    basis.append(coeff_vec)
basis_mat = Matrix(basis)

# 用LLL算法约化基，得到小系数基
reduced_basis = basis_mat.LLL()

# 用约化后的基构造新的通解
new_params = sympy.symbols('t0 t1 t2')
reduced_sol = [
    sum(reduced_basis[i, j] * new_params[i] for i in range(len(new_params)))
    for j in range(4)
]
print("约化后的通解：", reduced_sol)

运行输出 #

约化后的通解： [t0, 2*t1, 81*t0 + 9*t1 + 256*t2, 145*t0 + 20*t1 + 459*t2]

和你给出的Mathematica结果完全一致。

原理说明 #

n元线性齐次丢番图方程的所有整数解构成一个k维的整数格（k = 变量数 - 系数矩阵的秩），SymPy默认生成的基是未经过约化的，会包含大量冗余的大系数。LLL算法可以对整数格的基进行约化，输出范数尽可能小的基向量，对应的通解系数自然会降到最小的合理范围。

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