在Pyomo中定义一个变量在约束中相对于另一个变量的值而不使问题变为非线性

1 投票
1 回答
38 浏览
提问于 2025-04-14 15:26

背景参考:

[1] 在多个pyomo子集中生成带条件的约束

[2] 调试为什么一些请求意外未满足

我在解决上面提到的问题时,按照@Airsquid提供的方案得到了正确的结果,但当可用的能源供应量相对请求数量较少时,运行时间很长。

我添加了一些提前停止的设置,这在大多数情况下是可以的,因为我并不一定需要完全最优的解决方案,只要有一个足够好的就行。
solver.options = { 'sec': 600, 'threads':6, 'ratio': 5} #提前停止设置

我想尝试按照@Airsquid在[2]中给出的建议重构代码:
“在这种情况下,你从start变量就能知道请求的所有信息。所以,你可以进行大规模重构,可能会去掉一些辅助变量。如果它在第5期开始,请求持续10期,并且有一定的固定功率,你不需要一个变量来知道它是否在第8期运行——它就是在!同样,对于期内的分配,你仍然需要一个变量来跟踪每个请求分配了多少功率。”

由于我对线性编程和pyomo还很陌生,所以在语法的制定上遇到了一些问题。

我只想到两种方法来实现这个:

方法1:
定义m.dispatch[t,r]变量在m.start[t,r] == 1和m.start[t+request_length_periods,r]之间的值,下面是我对此的一个粗略语法尝试!我知道这不正确,也知道为什么不正确,但我不知道解决方案是什么。

    @m.Constraint(m.windows_flat)
    def request_length(m, t, r):
        request_length_periods = int(m.duration_periods[r])
        request_periods = set(list(range(t,t+request_length_periods+1)))
        return (m.dispatch[t,r] == m.request_power_size[r] for t in request_periods if m.start[t,r] == 1)

这显然给我带来了错误:

ERROR: Rule failed when generating expression for Constraint request_length
with index (910192, 88): ValueError: Constraint 'request_length[910192,88]'
does not have a proper value. Found '<generator object
generalOptimisation.<locals>.request_length.<locals>.<genexpr> at
0x168022b20>' Expecting a tuple or relational expression. Examples:
       sum(model.costs) == model.income (0, model.price[item], 50)
ERROR: Constructing component 'request_length' from data=None failed:
ValueError: Constraint 'request_length[910192,88]' does not have a proper
value. Found '<generator object

我知道我定义的request_periods是不正确的,但我不知道在不知道m.start[t,r] == 1的情况下,如何定义这个集合?

方法2:
m.start[t,r] *(m.dispatch[t,r]...m.dispatch[t+request_length_periods,r]) == 1
(我不知道在Pyomo中这个语法是什么,我知道这会使解决方案变得非线性)

如果有人能根据@Airsquid的建议给出正确的公式,我将非常感激。目前,我打算进一步放宽提前停止的设置来解决这个问题,但这并不是完全理想的。

关于提前停止设置的任何建议也欢迎——也许600秒的系统超时并不现实!

1 个回答

0

这是对你原始代码的一个改进版本,它使用了一个二进制变量来处理任务。这里面有一些地方可能还可以进一步优化,但这展示了一种解决问题的方法。它还假设每个请求所需的能量在不同时间段内是恒定的。

经过几次测试,得到了一个相当独特的结果(如图所示)

代码

# energy assignment
import sys
from collections import defaultdict
from random import randint

import pyomo.environ as pyo
from matplotlib import pyplot as plt

### DATA
num_periods = 12
periods = tuple(range(num_periods))


# (earliest, latest)
request_timespan = {
    "r1": (0, 7),
    "r2": (2, 3),
    "r3": (1, 5),
}
request_power_per_period = {"r1": 5, "r2": 4, "r3": 8}
request_duration = {"r1": 4, "r2": 8, "r3": 5}
power_avail = dict(zip(periods, (randint(8, 20) for _ in range(num_periods))))
# prevent mind-numbing errors:

assert request_timespan.keys() == request_power_per_period.keys()
assert request_timespan.keys() == request_duration.keys()

# we know that we are going to need to sum power across each period at some point, so one
# approach is to determine which request are eligible to be satisfied in each period.  So,
# you could probably either make a dictionary of requests: {timeslots} or perhaps better:
# timeslot: {eligible requests}... either way can work.

eiligible_starts = defaultdict(list)
for r, (start, end) in request_timespan.items():
    for t in [
        t for t in range(start, end + 1) if t + request_duration[r] <= num_periods
    ]:
        eiligible_starts[t].append(r)

# check it...
print(eiligible_starts)

### MODEL BUILD

m = pyo.ConcreteModel("dispatch")

### SETS
m.T = pyo.Set(initialize=periods)
m.R = pyo.Set(initialize=tuple(request_timespan.keys()))
# Note:  This will be an "indexed" set, where we have sets indexed by some index, in this case, m.T
m.windows = pyo.Set(
    m.T,
    initialize=eiligible_starts,
    within=m.R,
    doc="requests eligible in each timeslot",
)
# We can make a flat-set here which will be a good basis for the dispatch variable
m.windows_flat = pyo.Set(
    initialize={(t, r) for t in eiligible_starts for r in eiligible_starts[t]},
    within=m.T * m.R,
)

### PARAMS
m.request_size = pyo.Param(m.R, initialize=request_power_per_period)
m.power_limit = pyo.Param(m.T, initialize=power_avail)

### VARS
m.dispatch = pyo.Var(
    m.windows_flat, domain=pyo.Binary, doc="dispatch power in timeslot t to request r"
)

### OBJ
m.obj = pyo.Objective(
    expr=sum(m.dispatch[t, r] for (t, r) in m.windows_flat), sense=pyo.maximize
)


### CONSTRAINTS
@m.Constraint(m.R)
def satisfy_only_once(m, r):
    return sum(m.dispatch[t, rr] for (t, rr) in m.windows_flat if r == rr) <= 1


@m.Constraint(m.T)
def supply_limit(m, t):
    # we need to sum across all dispatches that could be running in this period
    possible_dispatches = {
        (tt, r) for (tt, r) in m.windows_flat if 0 <= t - tt < request_duration[r]
    }
    if not possible_dispatches:
        return pyo.Constraint.Skip
    return (
        sum(m.dispatch[tt, r] * m.request_size[r] for tt, r in possible_dispatches)
        <= m.power_limit[t]
    )


# check it
m.pprint()

# solve it
solver = pyo.SolverFactory("cbc")
res = solver.solve(m)
print(res)

m.dispatch.display()


plt.step(periods, [power_avail[p] for p in periods], color="g")
assigned_periods = {}
for t, r in m.dispatch:
    if pyo.value(m.dispatch[t, r]) > 0.5:  # request was assigned
        assigned_periods[r] = list(range(t, t + request_duration[r]))
        print("hit", t, r)
total_period_power_assigned = []
for t in m.T:
    assigned = 0
    for r in m.R:
        if t in assigned_periods.get(r, set()):
            assigned += request_power_per_period[r]
    total_period_power_assigned.append(assigned)

print(total_period_power_assigned)

plt.step(periods, total_period_power_assigned)
plt.show()

输出 [截断]

dispatch : dispatch power in timeslot t to request r
    Size=15, Index=windows_flat
    Key       : Lower : Value : Upper : Fixed : Stale : Domain
    (0, 'r1') :     0 :   0.0 :     1 : False : False : Binary
    (1, 'r1') :     0 :   1.0 :     1 : False : False : Binary
    (1, 'r3') :     0 :   0.0 :     1 : False : False : Binary
    (2, 'r1') :     0 :   0.0 :     1 : False : False : Binary
    (2, 'r2') :     0 :   0.0 :     1 : False : False : Binary
    (2, 'r3') :     0 :   0.0 :     1 : False : False : Binary
    (3, 'r1') :     0 :   0.0 :     1 : False : False : Binary
    (3, 'r2') :     0 :   1.0 :     1 : False : False : Binary
    (3, 'r3') :     0 :   0.0 :     1 : False : False : Binary
    (4, 'r1') :     0 :   0.0 :     1 : False : False : Binary
    (4, 'r3') :     0 :   1.0 :     1 : False : False : Binary
    (5, 'r1') :     0 :   0.0 :     1 : False : False : Binary
    (5, 'r3') :     0 :   0.0 :     1 : False : False : Binary
    (6, 'r1') :     0 :   0.0 :     1 : False : False : Binary
    (7, 'r1') :     0 :   0.0 :     1 : False : False : Binary
hit 1 r1
hit 4 r3
hit 3 r2
{'r1': [1, 2, 3, 4], 'r3': [4, 5, 6, 7, 8], 'r2': [3, 4, 5, 6, 7, 8, 9, 10]}
[0, 5, 5, 9, 17, 12, 12, 12, 12, 4, 4, 0]

在这里输入图片描述

撰写回答