皮莫 | Couenne 求解器 |将索引变量域限制为两个整数值

Posted 2023-02-18

【中文标题】皮莫 | Couenne 求解器 |将索引变量域限制为两个整数值

【英文标题】：Pyomo | Couenne solver | restrict indexed variable domain to two integer values

【发布时间】：2021-07-08 08:24:39

【问题描述】：

我是 Pyomo 的新手，正在努力理解 Pyomo 语法背后的直觉以及它如何构建模型。这可能就是为什么我无法弄清楚如何定义和解决“二元”问题，其中 N 变量必须使用 Pyomo 和 Couenne 求解器仅取 ±1 值。

首先我用bounds=(-1, 1) 尝试了Integers 域并尝试添加严格的不等式：

import pyomo.environ as pyo
import numpy as np
N = 5
w = np.ones(range(N))
pyoModel = pyo.ConcreteModel('binary model')
pyoModel.o = pyo.Var(range(N), bounds=(-1, 1), within=pyo.Integers, initialize=1)
pyoModel.binaryConstraintP = pyo.Constraint(range(N), rule=strictlyPositive)
pyoModel.binaryConstraintN = pyo.Constraint(range(N), rule=strictlyNegative)
pyoModel.objective = pyo.Objective(expr=pyo.sum_product(pyoModel.o, w, index=range(N)), sense=pyo.maximize)
def strictlyPositive(model, i):
    return model.o[i] > 0
def strictlyNegative(model, i):
    return model.o[i] < 0

最后是：

ValueError：约束 'binaryConstraintP[0]' 遇到严格的不等式表达式（'>' 或 '=' 或 '==' 来表述。

好的，不允许严格的不等式（不知道为什么！），我尝试切换到 Binary 域并通过操作目标中的变量来解决问题，使其位于 -1, 1 - 即，如果 o ∈ 0, 1 那么 2 xo - 1 ∈ -1, 1：

import pyomo.environ as pyo
import numpy as np
N = 5
w = np.ones(range(N))
pyoModel = pyo.ConcreteModel('binary model')
pyoModel.o = pyo.Var(range(N), within=pyo.Binary, initialize=1)
pyoModel.objective = pyo.Objective(expr=pyo.sum_product(2 * pyoModel.o - 1, w, index=range(N)), sense=pyo.maximize)

得到：

TypeError: *: 'int' 和 'IndexedVar' 的操作数类型不受支持

所以我使用了一个二和一的数组，而不是 2 和 1，但在形状广播方面出现了另一个错误。我确定我在这里遗漏了一些东西，因为在目标中构造一个线性方程应该很容易吧？

我也尝试将域更改为用户定义的域：

...
pyoModel.domain = pyo.Set(initialize=[-1, 1])
...
pyoModel.o = pyo.Var(range(N), domain=pyoModel.domain, initialize=1)
...
with SolverFactory('couenne') as opt:
    results = opt.solve(pyoModel, load_solutions=False)
    ...

并以 Couenne 错误结束：

TypeError：名称为“%s”的变量的域类型无效。变量不是连续的、整数的或二进制的。

我也想过使用 SOS，但更难理解它们的工作原理！

同样，我必须在每种方法中遗漏一些东西。如有任何帮助，我们将不胜感激。

旁注：我尽可能简化了原始代码以使其更易于阅读。

【参考方案1】：

不允许严格不等式的原因是它可能不受底层解决方案理论的支持。我相信这适用于 MIP。您是否尝试在您定义的约束范围内将限制设置为非常小的值？例如

def strictlyPositive(model, i):
    return model.o[i] >= 0.000000000000001
def strictlyNegative(model, i):
    return model.o[i] <= -0.000000000000001

以这种方式模拟严格的不等式并不理想，因为它可能会导致数值问题，但值得一试。

【讨论】：

你认为这适用于整数吗？

老实说，我不太确定这是一个疯狂的猜测，抱歉。

【参考方案2】：

由于使用了严格的不等式，您的第一次尝试失败了，这是一个禁忌。这背后有理论，因为这类问题的解决者在问题空间的“凸包”上工作。有关更多信息，请参阅有关线性规划的文本——它超出了堆栈溢出答案的范围。

您的第二次尝试是正确的。如果您正在寻找 ±1 的数学等效值，则让 x 为二进制变量并使用线性转换 2x-1 是完全合适的。您在那里的尝试失败了，因为您通过在 Var() 构造中提供索引来将 x 变量声明为索引变量，但您没有在目标中使用该索引。

这是一个使用索引变量的示例。如果您只有一个单例，则只需删除索引集引用。

from pyomo.environ import *

some_constants =   0: 100,
                    1: 200,
                    2: -50,
                    3: 300,
                    4: 50

m = ConcreteModel('plus & minus ones project')

m.S = Set(initialize=range(5))
m.c = Param(m.S, initialize=some_constants)
m.x = Var(m.S, within=Binary)    # creates m.x[0], m.x[1], ... , m.x[4]

# try to maximize the sum of x*c
m.obj = Objective(expr=sum((2*m.x[s] - 1)*m.c[s] for s in m.S), sense=maximize)

# some constraint to limit the number of "+1 picks" to 2...easy with binary vars.
m.C1 = Constraint(expr=sum(m.x[s] for s in m.S) <= 2)

m.pprint()

产量：

1 Set Declarations
    S : Size=1, Index=None, Ordered=Insertion
        Key  : Dimen : Domain : Size : Members
        None :     1 :    Any :    5 : 0, 1, 2, 3, 4

1 Param Declarations
    c : Size=5, Index=S, Domain=Any, Default=None, Mutable=False
        Key : Value
          0 :   100
          1 :   200
          2 :   -50
          3 :   300
          4 :    50

1 Var Declarations
    x : Size=5, Index=S
        Key : Lower : Value : Upper : Fixed : Stale : Domain
          0 :     0 :  None :     1 : False :  True : Binary
          1 :     0 :  None :     1 : False :  True : Binary
          2 :     0 :  None :     1 : False :  True : Binary
          3 :     0 :  None :     1 : False :  True : Binary
          4 :     0 :  None :     1 : False :  True : Binary

1 Objective Declarations
    obj : Size=1, Index=None, Active=True
        Key  : Active : Sense    : Expression
        None :   True : maximize : (2*x[0] - 1)*100 + (2*x[1] - 1)*200 + (2*x[2] - 1)*-50 + (2*x[3] - 1)*300 + (2*x[4] - 1)*50

1 Constraint Declarations
    C1 : Size=1, Index=None, Active=True
        Key  : Lower : Body                             : Upper : Active
        None :  -Inf : x[0] + x[1] + x[2] + x[3] + x[4] :   2.0 :   True

5 Declarations: S c x obj C1

【讨论】：

感谢您的快速回答！所以我错过的地方是我宣布目标的方式。事实上，我想避免 for 循环，因为我的“真实”案例更复杂：xt * M * x。顺便说一句，您的解决方案是否扩展到参数m.c 的矩阵形式？