python - 使用神秘来解决不等式系统

Question

我正在尝试使用 Python 框架Mystic解决一个简单的不等式系统。我正在求解十个变量，这些变量应该都大于彼此（x1 > x2 > x3 > ... > x10）并且它们应该大于零。我希望这些变量仅限于整数，我还没有设法完成。

这十个变量用于具有大约 15-20 个方程的不等式系统。我已经用所有这些不等式创建了一个字符串，简化了它们，生成了求解器，从这些求解器中生成了约束，就像文档解释的那样。但是，我终其一生都无法理解如何从我的代码生成的约束对象中获得解决方案（即神秘分配给我想要解决的变量的数字）。我已经搜索了几个小时的文档，并检查了大多数关于神秘和不平等解决的 stackoverflow 帖子。大多数问题都比我的水平高得多，所以它们没有多大帮助。

这是我的代码的摘录。

inequalities = '''
x1 > x2
x2 > x3
# etc.
x9 > x10
x10 > 0
4*x1 + 4*x2 + 2*x3 + 2*x4 + 2*x6 + x7 + 2*x8 > x4 + x6 + x8 + x9
4*x1 + 4*x2 + 2*x3 + 2*x4 + 2*x6 + x7 + 2*x8 > 2*x1 + x2 + 2*x3 + x5 + 2*x6 + x7 + x8 + x10
4*x1 + 4*x2 + 2*x3 + 2*x4 + 2*x6 + x7 + 2*x8 > 2*x1 + 5*x3 + 3*x4 + x8 + x9 + x10
# etc.
'''

def solve_inequalities(inequalities):
    var = ms.get_variables(inequalities)
    eqns = ms.simplify(inequalities, variables=var)
    print("Simplified equations: ", eqns)

    solver = ms.generate_solvers(eqns, var)
    constraint = ms.generate_constraint(solver)
    # which of these objects do i need to call?
    solution = constraint([1,2,3,4,5,6,7,8,9,10])
    print("Solution: ", constraint)

solve_inequalities(inequalities)

哪个打印：

Simplified equations:  x6 > x7
x3 > x4
x2 > -x1 - x3/2 - x4/2 - x6/4 + 3*x9/4 + x10/4
x3 > -2*x1 - 3*x2/2 - x4 + x5 - x6/2 - x7/2 - x8/2 + x9 + 3*x10
x1 > -x2 - x3/4 - x4/2 - x6/2 - x8/4
x2 > -x1 - x3/4 - x4/4 - x6/4 - x8/2 + 3*x9/4 + x10/2
x5 < -2*x1 + x2 - x3 + x4 + x6 + x8 - x10
x5 > x6
x4 > -x1 - 3*x2/2 + x5/2 - x8/2 + x10/2
x4 > x5
x2 > x3
x8 > x9
x5 < 3*x1/2 + x2 + x3/2 + x4/2 - x6/2
x2 > -x1 + x4/4 + x5 - x6/4 - x7/4 - x8/4 + x9/4 + x10/4
x1 > -x2 - x3/2 - x4/4 - x6/2 - x7/4 - x8/2 + x9/2
x1 > -x2 - x3/4 - x4/4 + x5/4 + 3*x7/4
x1 > x2
x2 > -3*x1/4 - x3/2 - x4/2 - x6/2 - x7/4 - x8/2 + x9/2
x1 > -x2 - x3/2 - x4/2 - x6/2 - x7/4 - x8/4
x10 > 0
x7 > x8
x3 > -3*x1/2 - 3*x2/2 - x4/2 + 3*x5/2 - x6 + 3*x7/2 - x8/2 + x10
x3 > -3*x1/2 + x2 - x4/2 + 3*x5/2 - x6/2 + x7/2 - x8/2 + x10
x9 > x10
x5 < 3*x1/2 + 3*x2/2 + x3/2 + x8 - 3*x9/2 - x10/2
x2 > -x1 - x3/2 - x4/4 - x6/4 - x8/4 + x9/4 + x10/4
x2 > -x1/2 + 3*x3/4 + x4/4 - x6/2 - x7/4 - x8/4 + x9/4 + x10/4
x1 > -x2 - x3/2 - x4/4 - x6/4 - x7/4 - x8/4 + x9/4

Solution:  [20.875000000000057, 20.875000000000036, 20.875000000000014, 6.000000000000021, 6.000000000000014, 6.000000000000007, 8.000000000000018, 8.000000000000009, 9.00000000000001, 9, 10]

这显然是不正确的。（我想知道的其他事情：为什么这些不等式的原始顺序没有得到保留？）

任何帮助将不胜感激。

问候，伦纳特

标签： pythonoptimizationconstraintsmystic

Answer 1

我是mystic作者。抱歉，文档并不像预期的那样明显。我猜，文档中有几个地方有两个微妙的地方，包括generate_constraint. 一种是默认情况下迭代应用约束，用户有责任确保各个方程不冲突。这意味着，如果你有x1 = x2 + x3然后x1 = x2 + x4，它将首先设置x1为x2 + x3然后设置x1为x2 + x4... 从而替换 的初始值x1。要获得您所期望的，您需要在简化时更改顺序获取x1 = x2 + x3和x4 = x1 - x2...或者您需要设置join为and_. 前者可以使用一些关键字自动尝试simplify，例如cycle和target......但不能保证它会比你手动排序方程更好sort。我将在下面使用您的简化不等式演示后者（因为我没有您的原始方程式）。

inequalities = '''
x1 > x2
x2 > x3
x3 > x4
x4 > x5
x5 > x6
x6 > x7
x7 > x8
x8 > x9
x9 > x10
x10 > 0
x2 > -x1 - x3/2 - x4/2 - x6/4 + 3*x9/4 + x10/4
x3 > -2*x1 - 3*x2/2 - x4 + x5 - x6/2 - x7/2 - x8/2 + x9 + 3*x10
x1 > -x2 - x3/4 - x4/2 - x6/2 - x8/4
x2 > -x1 - x3/4 - x4/4 - x6/4 - x8/2 + 3*x9/4 + x10/2
x5 < -2*x1 + x2 - x3 + x4 + x6 + x8 - x10
x4 > -x1 - 3*x2/2 + x5/2 - x8/2 + x10/2
x5 < 3*x1/2 + x2 + x3/2 + x4/2 - x6/2
x2 > -x1 + x4/4 + x5 - x6/4 - x7/4 - x8/4 + x9/4 + x10/4
x1 > -x2 - x3/2 - x4/4 - x6/2 - x7/4 - x8/2 + x9/2
x1 > -x2 - x3/4 - x4/4 + x5/4 + 3*x7/4
x2 > -3*x1/4 - x3/2 - x4/2 - x6/2 - x7/4 - x8/2 + x9/2
x1 > -x2 - x3/2 - x4/2 - x6/2 - x7/4 - x8/4
x3 > -3*x1/2 - 3*x2/2 - x4/2 + 3*x5/2 - x6 + 3*x7/2 - x8/2 + x10
x3 > -3*x1/2 + x2 - x4/2 + 3*x5/2 - x6/2 + x7/2 - x8/2 + x10
x5 < 3*x1/2 + 3*x2/2 + x3/2 + x8 - 3*x9/2 - x10/2
x2 > -x1 - x3/2 - x4/4 - x6/4 - x8/4 + x9/4 + x10/4
x2 > -x1/2 + 3*x3/4 + x4/4 - x6/2 - x7/4 - x8/4 + x9/4 + x10/4
x1 > -x2 - x3/2 - x4/4 - x6/4 - x7/4 - x8/4 + x9/4
'''

所以，就像你之前做的那样，但这次，使用join=and_.

>>> import mystic as my
>>> var = my.symbolic.get_variables(inequalities)
>>> constraint = my.symbolic.generate_constraint(my.symbolic.generate_solvers(inequalities, var), join=my.constraints.and_)
>>> constraint([1,2,3,4,5,6,7,8,9,10])
[12.750000000000028, 12.750000000000014, 12.75, 10.000000000000064, 10.000000000000053, 10.000000000000043, 10.000000000000032, 10.000000000000021, 10.00000000000001, 10]

你可以看到它有效。

如果您想构建一个尝试使用整数的约束，那么您可以这样做：

>>> c = my.constraints.integers(float)(lambda x:x)
>>> c([1.1,2.3,3.7,4.3,5,6,7,8.932,9.0002,10])
[1.0, 2.0, 4.0, 4.0, 5.0, 6.0, 7.0, 9.0, 9.0, 10.0]

但是，将这两个约束耦合在一起很可能会失败，如您在此处看到的：

>>> con = my.constraints.and_(c, constraint, maxiter=1000)
>>> con([1,2,3,4,5,6,7,8,9,10])
[12.750000000000028, 12.750000000000014, 12.75, 10.000000000000064, 10.000000000000053, 10.000000000000043, 10.000000000000032, 10.000000000000021, 10.00000000000001, 10.0]

它失败的原因，这是另一个微妙的点，and_本质上是试图通过迭代应用约束的组合来循环......并希望它成功。通常它会这样做，但在你的情况下，它不会......所以它放弃并返回原来的x.

事实上，即使这样也不成功：

>>> ineq = '''
... x1 > x2
... x2 > x3
... x3 > x4
... x4 > x5
... x5 > x6
... x6 > x7
... x7 > x8
... x8 > x9
... x9 > x10
... x10 > 0
... '''
>>> 
>>> constraint = my.symbolic.generate_constraint(my.symbolic.generate_solvers(ineq, var), join=my.constraints.and_)
>>> con = my.constraints.and_(c, constraint, maxiter=1000)
>>> con([1,2,3,4,5,6,7,8,9,10])
[10.000000000000096, 10.000000000000085, 10.000000000000075, 10.000000000000064, 10.000000000000053, 10.000000000000043, 10.000000000000032, 10.000000000000021, 10.00000000000001, 10.0]
>>> 

但是，您会注意到，它似乎mystic是在约束的末尾键入10，然后确保它之前的每个条目都稍大一些。因此，知道了这一点，我们可以像这样重写约束：

>>> ineq = '''
... x1 > 0
... x1 < x2
... x2 < x3
... x3 < x4
... x4 < x5
... x5 < x6
... x6 < x7
... x7 < x8
... x8 < x9
... x9 < x10
... '''
>>> constraint = my.symbolic.generate_constraint(my.symbolic.generate_solvers(ineq, var), join=my.constraints.and_)
>>> con = my.constraints.and_(c, constraint, maxiter=10)
>>> con([1,2,3,4,5,6,7,8,9,10])
[1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0, 9.0, 10.0]
>>> con([1.1,2.3,3.7,4.3,5,6,7,8.932,9.0002,10])
[1.0, 2.0, 3.999999999999995, 4.0, 5.0, 6.0, 7.0, 8.99999999999999, 9.0, 10.0]

你可以看到它做得更好，这次几乎可以工作......但在第二种情况下仍然失败。

所以，在这种情况下......如果你要在优化问题中使用这些约束，我会使用constraint硬约束（如上所示），然后应用整数约束作为惩罚（例如软约束）... 或相反亦然。当我有两个失败的约束时and_，我通常会应用一个作为惩罚。

如果能够mystic更好地应用来自不公平系统的约束与附加约束（如整数约束）相结合，那就太好了，但它在幕后主要是蛮力，所以它可能会失败并迫使你依赖使用penalty.

编辑：啊哈......我有第二个想法。如果您有合理的界限，对于x，您可以执行以下操作......这可以更好地工作。

>>> cu = my.constraints.and_(my.constraints.discrete(range(100))(lambda x:x), my.constraints.impose_unique(range(100))(lambda x:x))
>>> 
>>> cu([1.1,2.3,3.9,4.3,5,6,7,8.9,9.0002,10])
[1.0, 2.0, 4.0, 38.0, 5.0, 6.0, 7.0, 9.0, 75.0, 10.0]
>>> 
>>> con = my.constraints.and_(c, constraint, cu, maxiter=1000)
>>> con([1.1,2.3,3.9,4.3,5,6,7,8.9,9.0002,10])
[1.0, 2.0, 4.0, 5.0, 6.0, 7.0, 9.0, 10.0, 32.0, 47.0]
>>> 

因此，这就是说，从 [0,99] 中的离散整数集中选择解决方案，其中所有项目都是唯一的......并确保它们满足您的不等式约束。

所以，现在，回到你更难的原始案例......它可能会起作用：

>>> constraint = my.symbolic.generate_constraint(my.symbolic.generate_solvers(inequalities, var), join=my.constraints.and_)
>>> con = my.constraints.and_(c, constraint, cu, maxiter=1000)
>>> con([1.1,2.3,3.9,4.3,5,6,7,8.9,9.0002,10])
[99.0, 98.0, 97.0, 96.0, 67.0, 54.0, 29.0, 20.0, 6.0, 1.0]

并检查...到目前为止看起来不错...

>>> x1, x2, x3, x4, x5, x6, x7, x8, x9, x10 = [99.0, 98.0, 97.0, 96.0, 67.0, 54.0, 29.0, 20.0, 6.0, 1.0]
>>> x1 > -x2 - x3/2 - x4/4 - x6/4 - x7/4 - x8/4 + x9/4
True
>>> x2 > -x1/2 + 3*x3/4 + x4/4 - x6/2 - x7/4 - x8/4 + x9/4 + x10/4
True
>>> x2 > -x1 - x3/2 - x4/4 - x6/4 - x8/4 + x9/4 + x10/4
True

耶！

编辑：有趣的是，通过构建一个小函数来检查哪些不等式失败......我注意到了一些有趣的事情。

def failures(x):
  return tuple(eval(i, dict(zip(var,x))) for i in inequalities.strip().split('\n'))

从本质上讲，即使使用新的约束求解器，解决方案有时也会失败——但失败只出现在一个方程中，而且只是有时。

这显然是困难的：

x5 < -2*x1 + x2 - x3 + x4 + x6 + x8 - x10

当包含这条线时，不等式求解器会运行最多 1000 次尝试......有时成功，有时失败。但是，如果我删除了这个约束，那么成功的几率是 100%，而且非常快（求解器使用了几次尝试）。

所以，这里的结论是，这条线应该被删除，并作为惩罚重铸……而其他一切都应该作为约束来处理。