是否有可能在 pyomo 中获得二维集合的子集

Question

我想确定 BEV 沿电网的最佳充电。因此，我使用 pyomo 来构建描述优化的 LP 模型。我有一组时间（一天以一刻钟为单位）：

model.times = pe.Set(initialize=list(range(96)))

和一组节点：

model.buses = pe.Set(initialize=list(range(6)))

我使用这两个集合来索引变量和约束。例如，我有一个约束跟踪 BEV 的 SOC：

def track_socs_rule(model, t, b):
    if t >= model.bevs[b].t_start and t <= model.bevs[b].t_target:
        return (model.SOC[t, b] + model.I[t, b] * model.voltages[b] * model.resolution/60 / 1000
                / model.bevs[b].e_bat*100 - model.SOC[t+1, b]) == 0
    else:
        return pe.Constraint.Skip

所以我必须在定义约束的规则中小心，只包含有效的时间步长（即当相应节点的 BEV 正在充电时）。它有效，我得到了正确的结果——但必须有一种更优雅的方式。我想定义一个只包含有效时间步长的集合，比如：

model.valid_subset = pe.Set(within=model.times*model.charger_buses,
                            initialize=[[t for t in model.times if t >= model.bevs[b].t_start and
                                         t <= model.bevs[b].t_target] for b in model.charger_buses]

但这给了我一个错误：

ValueError: Cannot add value([bla, bla, ...], [bla, bla...], ...) to Set valid_subset.
The value is not in the domain valid_subset_domain

这是一个应该重现错误的最小工作示例：

import pyomo.environ as pe

class BEV:
    def __init__(self, home_bus, e_bat, soc_start, soc_target, t_start, t_target, resolution):
        self.resolution = resolution
        self.home_bus = home_bus
        self.e_bat = e_bat
        self.soc_start = soc_start
        self.soc_target = soc_target
        self.t_start = int(t_start * 60/self.resolution)
        self.t_target = int(t_target * 60/self.resolution)


RESOLUTION = 15 # minutes

home_buses = [0, 2, 3, 4]
e_bats = [50, 50, 50, 50] # kWh
soc_starts = [20, 25, 20, 30] # %
soc_targets = [80, 75, 90, 100] # %
t_starts = [10, 12, 9, 15] # hours
t_targets = [19, 18, 15, 23] # hours

bevs = []
for pos in range(4):
    bev = BEV(home_bus=home_buses[pos], e_bat=e_bats[pos], soc_start=soc_starts[pos],
              soc_target=soc_targets[pos], t_start=t_starts[pos],
              t_target=t_targets[pos], resolution=RESOLUTION)
    bevs.append(bev)

solver = pe.SolverFactory('glpk')
model = pe.ConcreteModel()

model.bevs = {bev.home_bus: bev for bev in bevs}
model.resolution = RESOLUTION
model.times = pe.Set(initialize=list(range(int(24 * 60/RESOLUTION))))
model.buses = pe.Set(initialize=list(range(6)))
model.charger_buses = pe.Set(within=model.buses, initialize=[bus for bus in model.bevs])
#model.valid_subset = pe.Set(within=model.times*model.charger_buses,
                            #initialize=[[t for t in model.times if t >= model.bevs[b].t_start
                                         #and t <= model.bevs[b].t_target] for b in model.charger_buses])

model.voltages = pe.Param(model.buses, initialize={i: 400-i/2 for i in model.buses})

model.SOC = pe.Var(model.times*model.charger_buses, domain=pe.PositiveReals)
model.I = pe.Var(model.times*model.charger_buses, domain=pe.PositiveReals)

def track_socs_rule(model, t, b):
    # instead leave if (once indexed in model.valid_subset)
    if t >= model.bevs[b].t_start and t <= model.bevs[b].t_target:
        return (model.SOC[t, b] + model.I[t, b] * model.voltages[b] * model.resolution/60 /1000
                /model.bevs[b].e_bat*100 - model.SOC[t+1, b]) == 0

    else:
        return pe.Constraint.Skip

# instead index in model.valid_subset (once it works...)
model.track_socs = pe.Constraint(model.times*model.charger_buses, rule=track_socs_rule)
model.track_socs.pprint()

一旦取消注释带有 model.valid_subset 的行，就会发生错误。我的希望是，一旦我使用了这个“聪明”的索引，我就可以不用 track_socs_rule 中的 if 了，因为这些组件已经被正确地索引了。这也将帮助我确定 I 和 SOC 的上限和下限（它们甚至没有包含在代码中——为了简单起见，我保留了它们）。

经过大量解释，我的问题：

有没有办法在 pyomo 中构建二维集合的子集？

非常感谢您的提前帮助！

Answer 1

这确实是构建这样一个所需子集的一种方法。其实你的提法还不错，有错误的地方还需要修改

在行 model.valid_subset = pe.Set(within=model.times*model.charger_buses,... 中，您告诉模型将子集限制在 model.times 和 model.charger_buses 的 cross product 范围内，但是你试图用 t 初始化子集，因此你不能添加这样的子集，因为域是 (t_i,b_i) 而你只是添加 t_i.
对于 Set 组件，initialize arg 应该是一个实际值的列表（我认为任何 np.array-like iterable 应该工作，但不是完全确定），但是您要传递一个列表列表，然后会出现错误。您需要将您的列表列表展平为一个列表。

为了扁平化您的列表列表，我将使用 itertools as in this answer，但您可以按照自己的方式进行。

import itertools
#All your model
...
model.valid_subset = pe.Set(within=model.times*model.charger_buses,
                            initialize=list(itertools.chain(*[[(t,b) for t in model.times if t >= model.bevs[b].t_start
                                         and t <= model.bevs[b].t_target] for b in model.charger_buses])))

您可能希望在实际 initialize arg 之外构造任何子集。它可能是相同的列表理解，但不是在 Set （或任何）组件的实际 arg 中。这将有助于使其更具可读性。

是否有可能在 pyomo 中获得二维集合的子集

is it possible to get a subset of a 2-dimensional Set in pyomo

python

optimization

mathematical-optimization

linear-programming

pyomo