"searchsorted" 在 pyomo 优化中查找

Question

我有一个可以简化为以下内容的系统：使用发电和存储单元来满足需求。 objective 函数是发电成本乘以发电量。然而，所产生的电力被分成不同成本的箱子，发电的“清算价格”是每小时生产的最高箱子的成本：

T = np.arange(5, dtype=int)
produce_cap = 90 # MW
store_cap = 100  # MWh
store_init = 0   # MWh

m = pyo.ConcreteModel()

m.T = pyo.Set(initialize=T) # time, hourly

m.produce = pyo.Var(m.T, within=pyo.NonNegativeReals, initialize=0) # generation
m.store = pyo.Var(m.T, within=pyo.Reals, initialize=0)              # storage
stack = np.arange(10, 91, 20)             # cumulative sum of generation subidivisions
price = np.arange(0.9, 0.01, -0.2)        # marginal cost for subdivision of generation
demand = np.asarray([35, 5, 75, 110, 15]) # load to meet

m.produce_cap = pyo.Constraint(m.T, rule=lambda m, t: m.produce[t] <= produce_cap)
m.store_max = pyo.Constraint(m.T, rule=lambda m, t: m.store[t] <= store_cap)
m.store_min = pyo.Constraint(m.T, rule=lambda m, t: m.store[t] >= -store_cap)
rule = lambda m, t: m.produce[t] + m.store[t] == demand[t] # conservation rule
m.consv = pyo.Constraint(m.T, rule=rule)

# objective
def obj(stack, price, demand, m):
  cost = 0
  for t in m.T:
    load = m.produce[t]
    idx = np.searchsorted(stack, m.produce[t])
    p = price[idx] if idx < len(price) else 1000 # penalty for exceeding production capability
    cost += m.produce[t] * p
  return cost
rule = functools.partial(obj, stack, price, demand)
m.objective = pyo.Objective(rule=rule, sense=pyo.minimize)

# more constraints added below ...

问题似乎出在 objective 函数定义中，使用了 np.searchsorted 算法。具体错误是

Cannot create a compound inequality with identical upper and lower
    bounds using strict inequalities: constraint infeasible:
    produce[0]  <  produce[0] and 50.0 < produce[0]

如果我尝试实现自己的类似搜索排序的算法，我会遇到类似的错误。我收集了 Pyomo 试图创建的 objective 函数的表达式无法处理这种 table 查找，至少我是如何实现它的。我可以考虑另一种方法或重新制定吗？

Answer 1

这里发生了很多事情。

根本原因是对 Pyomo 工作原理的概念性误解：形成约束和 objective 的规则是 而不是 在优化期间调用的回调函数。相反，规则是 Pyomo 调用来生成模型的函数，这些规则应该是 return expression 对象。然后，Pyomo 通过几种标准中间格式（例如 LP、NL、BAR、GMS 格式）之一将这些表达式传递给底层求解器。因此，作为一般规则，您应该 而不是 具有以变量值为条件的逻辑的规则（规则可以运行，但结果将是初始变量值的函数，在优化过程中不会是updated/changed。

对于您的具体示例，挑战在于 searchsorted 遍历 m.produce 变量并将其与分割点进行比较。这导致 Pyomo 开始生成表达式对象（通过运算符重载）。然后，您运行与一个（已弃用的）功能发生冲突，其中 Pyomo 允许使用类似“lower <= m.x <= upper”的语法生成复合（范围）不等式表达式。

解决方案是您需要将 objective 重新表述为 return 一个 表达式 以获得 objective 成本。有几种方法可以做到这一点，“最佳”方法取决于模型的平衡和成本曲线的实际形状。从你的例子来看，成本曲线看起来是分段线性的，所以我会考虑直接重新制定表达式（使用中间变量和一组约束），或者使用 Pyomo 的“Piecewise”组件用于生成分段线性表达式。

"searchsorted" 在 pyomo 优化中查找

"searchsorted" lookup in pyomo optimization

python

numpy

pyomo