我玩的游戏有一个谜题,涉及解决以下等式：
x*411 + y*295 + z*161 = 3200
不想以为我只是把它打成了同情,我还没有真正用到这一点：
>>> from sympy import *
>>> x,y,z = symbols('x y z',integer=True,positive=True)
>>> solve(x*411 + y*295 + z*161 - 3200,[x,z])
[{x: -295*y/411 - 161*z/411 + 3200/411}]
嗯,这只给了我一个依赖的解决方案,但我想在域中我所有可能的解决方案我将变量约束到,例如(假设没有其他解决方案)[{x：4,y：2,z：6}]或[(4,2,6)]
当然我现在可以在嵌套循环中手动替换两个变量,或者手动解决它(就像我上面的解决方案一样),但我想知道如何让sympy(或其他库)为我做这件事.
最佳答案
SymPy可以solve Diophantine equations,但没有内置的方法来生成积极的解决方案.使用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,6),)
在幕后,integral_points很可能只运行一个多循环;虽然当它似乎不起作用时,它会尝试使用史密斯普通形式.
我知道你想要积极的解决方案,但是(a)很容易从答案中排除任何零含有元组; (b)在解决之前,用x-1等替换x也很容易; (c)坚持使用“非负”可以很容易地使用Mixed Integer Linear Programming module创建一个多面体
如上.
根据文档,人们还可以直接从不等式系统(“Hrep”)构建Polyhedron object.这将允许一个人明确地说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,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,6].