获得线性方程的所有正整数解

Question

我玩的游戏有一个谜题,涉及解决以下等式:

x*411 + y*295 + z*161 = 3200

不想以为我只是打了它sympy,我还没有真正用到那一点:

>>> from sympy import *
>>> x, y, z = symbols('x y z', integer=True, positive=True)
>>> solve(x*411 + y*295 + z*161 - 3200, [x, y, z])
[{x: -295*y/411 - 161*z/411 + 3200/411}]

嗯,这只给了我一个依赖的解决方案,但我想在域中所有可能的解决方案我将变量约束到,例如(假设没有其他解决方案)[{x: 4, y: 2, z:6}]或[(4, 2, 6)]

当然我现在可以在嵌套循环中手动替换两个变量,或者手动解决它(就像我上面的解决方案一样),但我想知道如何让sympy(或其他库)为我做这件事.

Answer 1

SymPy 可以求解丢番图方程，但没有内置方法来生成正解。使用Sage，我们可以轻松做到这一点：这里有四行代码，可以生成方程的所有非负整数解。

p = MixedIntegerLinearProgram()
w = p.new_variable(integer=True, nonnegative=True)
p.add_constraint(411*w[0] + 295*w[1] + 161*w[2] == 3200)
p.polyhedron().integral_points()

输出是((4, 2, 6),)

在幕后，integral_points很可能只是运行多个循环；尽管当这似乎不起作用时，它会尝试使用史密斯范式。

我知道您想要正解，但是（a）很容易从答案中排除任何包含零的元组；(b) 在求解之前，也很容易将 x 替换为 x-1 等； (c) 坚持“非负”可以很容易地使用上述 混合整数线性规划模块创建多面体。

根据文档，人们还可以直接从不等式系统（“Hrep”）构建多面体对象。这将允许人们明确地说 x >= 1 等，但我在这条路线上没有成功。

与 SymPy

SymPy 的 Diophantine 模块的输出是参数解，例如

(t_0, 2627*t_0 + 161*t_1 - 19200, -4816*t_0 - 295*t_1 + 35200)

在你的例子中。这可以在循环中使用，以非常有效的方式生成解决方案。棘手的问题是寻找参数 t_0 和 t_1 的界限。由于这只是一个示例，因此我查看了上面的最后一个表达式，并将限制 35200/4816 和 35200/295 直接插入下面的循环中。

from sympy import *
x, y, z = symbols('x y z')
[s] = diophantine(x*411 + y*295 + z*161 - 3200)
print(s)
t_0, t_1 = s[2].free_symbols
for t0 in range(int(35200/4816)+1):
    for t1 in range(int(35200/295)+1):
        sol = [expr.subs({t_0: t0, t_1: t1}) for expr in s]
        if min(sol) > 0:
            print(sol)

输出是[4, 2, 6].