pyomo 中的逻辑或强制连续赋值

Question

我想优化我的模型，找到项目 (probes) 与序列 (sequence) 的最大可能重叠 (assignment)。我有所有项目的开始和结束位置，因此可以按如下方式构建我的模型：

import pyomo
import pyomo.environ as pe
import pyomo.opt as po

sequence = [0, 1, 2, 3, 4, 5]
probes = ["a", "b", "c"]
probe_starts = {"a": 0, "b": 2, "c": 3}
probe_ends = {"a": 2, "b": 4, "c": 5}

# Model definition
model = pe.ConcreteModel()
model.sequence = pe.Set(initialize=sequence)
model.probes = pe.Set(initialize=probes)
model.starts = pe.Param(model.probes, initialize=probe_starts)
model.ends = pe.Param(model.probes, initialize=probe_ends)
model.assignment = pe.Var(model.sequence, model.probes, domain=pe.Binary)

# Objective
expr = sum([model.assignment[j, i] for j in model.sequence for i in model.probes])
model.objective = pe.Objective(expr=expr, sense=pe.maximize)

我现在有以下三个约束：

项目一次只能绑定一个（没有重叠的项目）
因为项目有开始和结束位置，所以它们只能在各自的开始和结束后分配
如果项目被分配，它们必须从开始到结束全部绑定

# One probe per sequence position
model.one_probe_bound = pe.ConstraintList()
for s in model.sequence:
    model.one_probe_bound.add(sum(model.assignment[s, p] for p in model.probes) <= 1)

# No assignment before/after start/end
model.define_length = pe.ConstraintList()
for s in model.sequence:
    for p in model.probes:
        if s < model.starts[p]:
            model.define_length.add(model.assignment[s, p] == 0)
        if s > model.ends[p]:
            model.define_length.add(model.assignment[s, p] == 0)

上面的两个约束都没有问题，但我找不到从我的第三个条件输入逻辑或的方法。我尝试使用 :

中描述的析取

# Only allow full assignment or none
def disjunct_rule(b, p, i):
    m = b.model()
    expr = sum(m.assignment[s, p] for s in m.sequence)
    if i:
        return expr == m.ends[p] - m.starts[p]
    else:
        return expr == 0

def disjunction_rule(m, p):
    return [m.disjunct[p, i] for i in [True, False]]

def xfrm(m):
    pe.TransformationFactory("gdp.bigm").apply_to(m)

model.disjunct = pyomo.gdp.Disjunct(model.probes, [True, False], rule=disjunct_rule)
model.disjunction = pyomo.gdp.Disjunction(model.probes, rule=disjunction_rule)
model.xfrm = pe.BuildAction(rule=xfrm)

查看 model.assignment 的矩阵表示，其中 sequence 沿列，probes 沿行，我得到以下内容：

array([[1, 1, 1, 0, 0, 0],
       [0, 0, 0, 1, 1, 0],
       [0, 0, 0, 0, 0, 1]])

如上所示，我得到的分配没有跨越项目的整个长度（例如，c / 3rd 项目只分配在最后一个位置，而它也应该绑定在前两个位置. 我可以在这个玩具示例中看到的唯一有效解决方案如下：

array([[1, 1, 1, 0, 0, 0],
      [0, 0, 0, 0, 0, 0],
      [0, 0, 0, 1, 1, 1]])

其中项目 a 和 c 被全部选中，而 b 根本没有被选中。这样我们就匹配了所有约束。我使用的求解器是 glpk。在此先感谢您的任何建议。

Answer 1

这是一个剪辑....

我提到的另一种方法是为序列中的每个点引入另一个二进制变量，并通过分配探针（以某种方式）控制它。下面的方法应该比那个更有效。

尝试使用更大的数据集...当前只有 1 个可行的解决方案。另外，我假设更好的解决方案会使用更少的探针，所以我 re-wrote objective

此解决方案的基础是您必须连接到起点（约束 1）一次和终点（约束 2）一次的约束。并且任何中间连接必须一致（约束3）。

我在一些需要的地方使用了一些 sub-setting。

代码：

# model to make contiguous connections across a sequence
# using as few connections (probes) as possible

import pyomo
import pyomo.environ as pe

sequence = [0, 1, 2, 3, 4, 5]
probes = ["a", "b", "c"]
probe_starts = {"a": 0, "b": 2, "c": 3}
probe_ends = {"a": 2, "b": 4, "c": 5}

# Model definition
model = pe.ConcreteModel()
model.sequence = pe.Set(initialize=sequence)
model.probes = pe.Set(initialize=probes)
model.starts = pe.Param(model.probes, initialize=probe_starts)
model.ends = pe.Param(model.probes, initialize=probe_ends)

model.assign = pe.Var(model.probes, domain=pe.Binary)  # 1 if probe p is used...

# Objective
obj = sum(model.assign[p] for p in model.probes) # use as few as possible...?
model.objective = pe.Objective(expr=obj, sense=pe.minimize)

# Constraints

# must connect start once
model.C1 = pe.Constraint(expr=sum(model.assign[p] for p in model.probes 
    if model.starts[p] == sequence[0]) == 1)

# must connect end once
model.C2 = pe.Constraint(expr=sum(model.assign[p] for p in model.probes 
    if model.ends[p] == sequence[-1]) == 1)

# must connect any intermediate connections...
# if probe p1 is selected, must select an eligible p2 follow on
def connect(model, p1):
    # create subset on the fly of legal follow-on connections
    # assumption here that sequence is a sequential list of ints by using the "+1"
    p2s = [p for p in model.probes if model.starts[p] == model.ends[p1] + 1]
    if not p2s:
        return pe.Constraint.Skip
    return sum(model.assign[p2] for p2 in p2s) == model.assign[p1]

non_completing_probes = [p for p in model.probes if model.ends[p] != sequence[-1]]
model.C3 = pe.Constraint(non_completing_probes, rule=connect)

solver = pe.SolverFactory('glpk')
result = solver.solve(model)
print(result)

model.display()

产量：

Problem: 
- Name: unknown
  Lower bound: 2.0
  Upper bound: 2.0
  Number of objectives: 1
  Number of constraints: 4
  Number of variables: 4
  Number of nonzeros: 5
  Sense: minimize
Solver: 
- Status: ok
  Termination condition: optimal
  Statistics: 
    Branch and bound: 
      Number of bounded subproblems: 0
      Number of created subproblems: 0
  Error rc: 0
  Time: 0.006771087646484375
Solution: 
- number of solutions: 0
  number of solutions displayed: 0

Model unknown

  Variables:
    assign : Size=3, Index=probes
        Key : Lower : Value : Upper : Fixed : Stale : Domain
          a :     0 :   1.0 :     1 : False : False : Binary
          b :     0 :   0.0 :     1 : False : False : Binary
          c :     0 :   1.0 :     1 : False : False : Binary

  Objectives:
    objective : Size=1, Index=None, Active=True
        Key  : Active : Value
        None :   True :   2.0

  Constraints:
    C1 : Size=1
        Key  : Lower : Body : Upper
        None :   1.0 :  1.0 :   1.0
    C2 : Size=1
        Key  : Lower : Body : Upper
        None :   1.0 :  1.0 :   1.0
    C3 : Size=1
        Key : Lower : Body : Upper
          a :   0.0 :  0.0 :   0.0

pyomo 中的逻辑或强制连续赋值

Logical OR in pyomo to force consecutive assignment

python

optimization

pyomo

代码：

产量：