使用自定义约束实现Python解算器

我有两个相互关联的变量,我想找到一个最佳解决方案,在这种情况下,这是其总和的最小值.现在,我们将它们称为X和Y,以及预定义的常量,它们加起来就是一组“变量” s1和s2(稍后将它们提供给约束)：

105896649.59 + X = s1
    -6738.82 + Y = s2

在搜索SciPy文档时,我遇到了linear programming解决方案,其中具有最小化函数(在这种情况下为X Y)以及一组与变量绑定的不等式和相等约束.就我而言,它们如下：

> X> = 0,Y> = 0
> s1> = 1,s2> = 1
> s2 /(s1 s2)= 0.0001％

对于这种特定情况,该代码很容易实现：

from scipy.optimize import linprog

lstConst = [105896649.59, -6738.82]

# function to minimise: X + Y
c= [1, 1]

# left-hand side of the equation for s2 / (s1 + s2) = 0.0001%
# i.e., -0.000001 * X + 0.999999 * Y
Aeq = [[-0.000001, 0.999999]]

# right-hand side of the equation
beq = [0.000001 * (lstConst[0] + lstConst[1]) - lstConst[1]]

# ensures minimum can't be a negative number
minX = max(1, max(1 -lstConst[0], 0))
minY = max(1, max(1 -lstConst[1], 0))

X_bounds = (minX, None)

Y_bounds = (minY, None)

res = linprog(c, A_eq=Aeq, b_eq=beq, bounds=[X_bounds, Y_bounds])

因此,我们具有X和Y的值以最小化x参数上的函数：

In [1]: res.x
Out[1]: array([1.00000000e+00, 6.84471676e+03])

我想以此方法为基础：

>实际上,还有另一组限制：s1和s2也必须是整数(请注意,X和Y浮点数没有问题).
>我没有为s1和s2之间的比率定义单个值,而是提供了不同可能比率的列表.

本质上,给定s1和s2之间的几个不同比率,我想找到X Y函数的最小值.这可以通过在列表上进行迭代以在每次迭代中定义Aeq和beq或定义其他限制(如果可能)来实现.

但是,对于整数限制以及如何使线性编程算法考虑在内,我一无所知.

如果有人建议使用SciPy和linprog以外的其他库/优化器,也欢迎您.

解决方法:

首先,重述问题：

minimize x + y, subject to:

    k1 + x = s1
    k2 + y = s2
    x >= 0
    y >= 0
    s1 >= 1
    s2 >= 1
    s2 / (s1 + s2) = k3

Where:

    k1 = 105896649.59
    k2 = -6738.82
    k3 = 0.000001

注意,您不需要s1和s2变量来在linprog中编码问题.没有s1和s2辅助变量,问题是：

minimize x + y, subject to:

  x >= 0
  y >= 0
  x + k1 >= 1,
  y + k2 >= 1,
  (1-k3)y - k3x = (k1 + k2)k3 - k2

在linprog中阅读和编写代码比较容易：

import numpy as np
from scipy.optimize import linprog
k1, k2, k3 = 105896649.59, -6738.82, 0.000001
A_ub = -np.eye(2)
b_ub = [k1-1, k2-1]
A_eq = [[-k3, (1-k3)]]
b_eq = (k1 + k2)*k3 -k2
res = linprog([1,1], A_ub=A_ub, b_ub=b_ub, A_eq=A_eq, b_eq=b_eq, bounds=[[0,None], [0, None]])
print(res)

给出[0.,6844.71675549]的x = 1,因为实际上您已经将x和y的下限设置为1(我认为这是一个错字…),但是在问题的背景：

关于这个问题：

… I’m clueless as to the integer restriction and how to make the linear programming algorithm take it into account.

If anyone has an alternative suggestion that uses a library/optimizer other than SciPy and linprog, that’s also welcome.

您的要求是mixed integer linear programming (MILP).MILP和线性编程(LP)通常是用不同的算法解决的,而MILP问题通常很难精确解决. SciPy Optimize不支持MILP.有许多开放源代码工具,例如OrTools和PySCIPOpt,它们是SCIP的Python包装器.

PySCIPOpt中的示例：

PySCIPOpt很不错,因为它具有约束编程类型API.在PySCIPOpt中,您的问题很容易以可读的形式陈述.重新引入辅助变量,我们可以逐字逐句地输入约束条件：

from pyscipopt import Model

k1, k2, k3 = 105896649.59, -6738.82, 0.000001
model = Model()
x = model.addVar(vtype="CONTINUOUS", name="x", lb=0)
y = model.addVar(vtype="CONTINUOUS", name="y", lb=0)
s1 = model.addVar(vtype="CONTINUOUS", name="s1", lb=None, ub=None)
s2 = model.addVar(vtype="CONTINUOUS" name="s2", lb=None, ub=None)
o = model.addVar(vtype="CONTINUOUS", name="Objective Value", lb=0, ub=None)
model.addCons(k1 + x == s1)
model.addCons(k2 + y == s2)
model.addCons(s1 >= 1)
model.addCons(s2 >= 1)
model.addCons(s2/(s1+s2) == k3)
model.addCons(x + y == o)
model.setObjective(o, "minimize")
model.optimize()
print('x + y = o -> (%.4f + %.4f = %.4f)' % (model.getVal(x), model.getVal(y), model.getVal(o)))

与linprog给出的答案相同,因为它只是一个线性程序.但是,由于SCIP支持MILP,因此我们可以引入整数变量.要处理您的情况#1,只需将s1和s2更改为整数：

...
s1 = model.addVar(vtype="INTEGER", name="s1", lb=None, ub=None)
s2 = model.addVar(vtype="INTEGER", name="s2", lb=None, ub=None)

给出：

...
SCIP Status        : problem is solved [optimal solution found]
Solving Time (sec) : 0.00
Solving Nodes      : 1
Primal Bound       : +1.10089229999989e+05 (1 solutions)
Dual Bound         : +1.10089229999989e+05
Gap                : 0.00 %
x + y = o -> (103244.4100 + 6844.8200 = 110089.2300)

这是完全不同的解决方案…但这就是为什么MILP不是LP的原因.

从上面的示例中,并通过阅读docs,您应该能够弄清楚如何编写#2格-基本上像1 / k3这样的内容成为模型中的另一个整数变量.

标签：scipy,linear-programming,python
来源： https://codeday.me/bug/20191210/2104349.html