Pyomo 中的双索引变量或多个单索引变量之间有区别吗？

100% 你想使用第二条路线并索引所有内容，而不是让一堆单例变量滚动。

• 跟踪所有内容要容易得多
• 您可以非常轻松地使用 Pyomo 的索引功能来创建“针对每个”的约束
• 减少错别字等风险。

在您的模型中，您有 2 个自然集： 和 ，我希望大多数变量都由这 2 个集合进行双索引。timevehicle

然而，首先（这可能是一篇很长的帖子来证明这一点......对不起！在我看来，你在如何对待模型构建中是在玩火。您隐式依赖实例变量来覆盖所有时间段，但它只是单个变量。它不能同时表示模型的所有周期中的 soc。你（有点）可以构建你正在做的表达式，但它很脆弱。含义：它正在为 SOC 制作正确的表达式以供模型使用，但它在您执行此操作时会损坏实例变量。下面显示了一个示例。请注意，SOC 的表达式是您所期望的，因为只要您以正确的顺序传递时间步长，它就会按顺序构建，但是当我再次使用该变量时，它在“坏”约束中已损坏，因为它仍然保留上次用法，并且它对时间一无所知。有意义？您可以（如果您小心的话）继续使用您生成的表达式（请参阅下一个约束），因为它在构建时捕获了正确的表达式。出于类似的原因，稍后必须再次使用该表达式来计算 SOC。socsocev.soc

易碎示例：

"""
EV charge model
"""

import pyomo.environ as pyo

class SimpleEV:
    idx = 1 # rolling counter
    def __init__(self):
        self.soc = 10  # initial charge
        self.idx = SimpleEV.idx
        SimpleEV.idx += 1

    def add_juice(self, charge):
        self.soc += charge

    def use_juice(self, charge):
        self.soc -= charge

    # for cleaner printing in pyomo, provide a __repr__
    def __repr__(self):
        return f'ev{self.idx}'

m = pyo.ConcreteModel('ev')

m.T = pyo.Set(initialize=range(3))      # time
m.E = pyo.Set(initialize=[SimpleEV()])  # 1 simple EV

m.charge = pyo.Var(m.E, m.T, domain=pyo.NonNegativeReals, doc='charge added to vehicle e at time t')

m.draws = pyo.Param(m.T, initialize={1:5}, default=0)

def expr_generator(m, e: SimpleEV, t):
    e.add_juice(m.charge[e, t])
    e.use_juice(m.draws[t])
    return e.soc

m.soc_expression = pyo.Expression(m.E, m.T, rule=expr_generator)

# this is NOT going to work using the .soc variable:
def bad_min_charge(m, e: SimpleEV, t):
    return e.soc >= 2
m.bad_min_charge_constraint = pyo.Constraint(m.E, m.T, rule=bad_min_charge)

def ok_min_charge(m, e: SimpleEV, t):
    return m.soc_expression[e, t] >= 2
m.ok_min_charge_constraint  = pyo.Constraint(m.E, m.T, rule=ok_min_charge)

m.pprint()

输出（截断）：

1 Expression Declarations
    soc_expression : Size=3, Index=soc_expression_index
        Key      : Expression
        (ev1, 0) : 10 + charge[ev1,0]
        (ev1, 1) : 10 + charge[ev1,0] + charge[ev1,1] - 5
        (ev1, 2) : 10 + charge[ev1,0] + charge[ev1,1] - 5 + charge[ev1,2]

2 Constraint Declarations
    bad_min_charge_constraint : Size=3, Index=bad_min_charge_constraint_index, Active=True
        Key      : Lower : Body                                                   : Upper : Active
        (ev1, 0) :   2.0 : 10 + charge[ev1,0] + charge[ev1,1] - 5 + charge[ev1,2] :  +Inf :   True
        (ev1, 1) :   2.0 : 10 + charge[ev1,0] + charge[ev1,1] - 5 + charge[ev1,2] :  +Inf :   True
        (ev1, 2) :   2.0 : 10 + charge[ev1,0] + charge[ev1,1] - 5 + charge[ev1,2] :  +Inf :   True
    ok_min_charge_constraint : Size=3, Index=ok_min_charge_constraint_index, Active=True
        Key      : Lower : Body                  : Upper : Active
        (ev1, 0) :   2.0 : soc_expression[ev1,0] :  +Inf :   True
        (ev1, 1) :   2.0 : soc_expression[ev1,1] :  +Inf :   True
        (ev1, 2) :   2.0 : soc_expression[ev1,2] :  +Inf :   True

一个更好的主意：我会构建类似于下面的模型。几点说明：

• 如图所示，我会将路线与 EV 分开。我认为它更干净。
• 您可以安全地将 EV 直接放入套装中，只要您不使用 and 函数。__hash____eq__
• 您应该在代码中加入 a，以便输出更清晰。（参见 python dox 获取指导，我在这里偷工减料）__repr__()
• 你的语法对于制作组件来说非常笨拙，我认为我下面的语法对于制作集合、约束等来说更干净。
• 下面的结构将所有状态变量（如 SOC）与对象分开，我认为这更干净/更安全。您可以在优化后回填它们，正如我在最后所示。

更好的例子：

"""
improved EV charge model
"""

import pyomo.environ as pyo
class DriveProfile:
    def __init__(self, distances):
        self.distances = distances  # {time:km} dictionary

    def get_distance(self, time):
        return self.distances.get(time, None)  # NONE if no travels at this time block
class EV:
    idx = 1  # rolling counter

    def __init__(self, max_soc=100, efficiency=2, route: DriveProfile = None):
        self.initial_soc = 10  # initial charge
        self.max_soc = max_soc  # max charge
        self.discharge_rate = efficiency  # units/km
        self.route = route
        self.idx = EV.idx
        self.soc = None
        EV.idx += 1

    def load_soc(self, soc_dict: dict):
        self.soc = soc_dict

    # for cleaner printing in pyomo, provide a __repr__
    def __repr__(self):
        return f'ev{self.idx}'

# make some fake data...
times = list(range(4))

r1 = DriveProfile({1: 10, 3: 2.4})
r2 = DriveProfile({2: 20})

e1 = EV(route=r1, efficiency=3)
e2 = EV(max_soc=200, route=r2)

# make the model
m = pyo.ConcreteModel('ev')

# SETS
m.T = pyo.Set(initialize=times)     # time
m.E = pyo.Set(initialize=[e1, e2])  # ev's

# VARS
m.soc = pyo.Var(m.E, m.T, domain=pyo.NonNegativeReals, doc='state of charge of E at time T')
m.juice = pyo.Var(m.E, m.T, domain=pyo.NonNegativeReals, doc='juice added to E at time T')

# PARAMS
#   --- all params are internal to the EV objects.  We could capture them here, but not needed. ---

# OBJ:  Minimize the juice delivered
m.obj = pyo.Objective(expr=sum(m.juice[e, t] for e in m.E for t in m.T))


# CONSTRAINTS

# max charge
def max_charge(m, e: EV, t):
    return m.soc[e, t] <= e.max_soc
m.max_charge = pyo.Constraint(m.E, m.T, rule=max_charge)


def soc_balance(m, e: EV, t):
    # find any use (driving)
    use = e.route.get_distance(t)
    discharge = use * e.discharge_rate if use else 0
    if t == m.T.first():
        return m.soc[e, t] == e.initial_soc - discharge + m.juice[e, t]
    return m.soc[e, t] == m.soc[e, m.T.prev(t)] - discharge + m.juice[e, t]
m.soc_balance = pyo.Constraint(m.E, m.T, rule=soc_balance)

m.pprint()

# SOLVE
solver = pyo.SolverFactory('cbc')
soln = solver.solve(m)
print(soln)

# show the juice profile
m.juice.display()

# pass the soc info back to the ev's
print('The state of charge in the vehicles:')
for ev in m.E:
    # we can use items() to get tuples of ( (e, t), soc ) and screen them, then get the .value from the variable...
    res = {time: soc.value for (e, time), soc in m.soc.items() if e == ev}
    ev.load_soc(res)

    print(ev.__repr__(), ev.soc)

输出（截断第一部分）：

    soc_balance : Size=8, Index=soc_balance_index, Active=True
        Key      : Lower : Body                                                         : Upper : Active
        (ev1, 0) :   0.0 :                             soc[ev1,0] - (10 + juice[ev1,0]) :   0.0 :   True
        (ev1, 1) :   0.0 :                soc[ev1,1] - (soc[ev1,0] - 30 + juice[ev1,1]) :   0.0 :   True
        (ev1, 2) :   0.0 :                     soc[ev1,2] - (soc[ev1,1] + juice[ev1,2]) :   0.0 :   True
        (ev1, 3) :   0.0 : soc[ev1,3] - (soc[ev1,2] - 7.199999999999999 + juice[ev1,3]) :   0.0 :   True
        (ev2, 0) :   0.0 :                             soc[ev2,0] - (10 + juice[ev2,0]) :   0.0 :   True
        (ev2, 1) :   0.0 :                     soc[ev2,1] - (soc[ev2,0] + juice[ev2,1]) :   0.0 :   True
        (ev2, 2) :   0.0 :                soc[ev2,2] - (soc[ev2,1] - 40 + juice[ev2,2]) :   0.0 :   True
        (ev2, 3) :   0.0 :                     soc[ev2,3] - (soc[ev2,2] + juice[ev2,3]) :   0.0 :   True

11 Declarations: T E soc_index soc juice_index juice obj max_charge_index max_charge soc_balance_index soc_balance

Problem: 
- Name: unknown
  Lower bound: 57.2
  Upper bound: 57.2
  Number of objectives: 1
  Number of constraints: 16
  Number of variables: 16
  Number of nonzeros: 3
  Sense: minimize
Solver: 
- Status: ok
  User time: -1.0
  System time: 0.0
  Wallclock time: 0.0
  Termination condition: optimal
  Termination message: Model was solved to optimality (subject to tolerances), and an optimal solution is available.
  Statistics: 
    Branch and bound: 
      Number of bounded subproblems: None
      Number of created subproblems: None
    Black box: 
      Number of iterations: 2
  Error rc: 0
  Time: 0.007892131805419922
Solution: 
- number of solutions: 0
  number of solutions displayed: 0

juice : juice added to E at time T
    Size=8, Index=juice_index
    Key      : Lower : Value : Upper : Fixed : Stale : Domain
    (ev1, 0) :     0 :  27.2 :  None : False : False : NonNegativeReals
    (ev1, 1) :     0 :   0.0 :  None : False : False : NonNegativeReals
    (ev1, 2) :     0 :   0.0 :  None : False : False : NonNegativeReals
    (ev1, 3) :     0 :   0.0 :  None : False : False : NonNegativeReals
    (ev2, 0) :     0 :  30.0 :  None : False : False : NonNegativeReals
    (ev2, 1) :     0 :   0.0 :  None : False : False : NonNegativeReals
    (ev2, 2) :     0 :   0.0 :  None : False : False : NonNegativeReals
    (ev2, 3) :     0 :   0.0 :  None : False : False : NonNegativeReals
The state of charge in the vehicles:
ev1 {0: 37.2, 1: 7.2, 2: 7.2, 3: 0.0}
ev2 {0: 40.0, 1: 40.0, 2: 0.0, 3: 0.0}

除了建模风格之外，这两种情况之间几乎没有区别：

• 从求解器的角度来看，模型是相同的：相同的变量、相同的目标表达式和相同的约束表达式。（也就是说，对于一般问题，它们可能会导致不同的行/列顺序）。
• 从 Pyomo 的角度来看，案例 1 将使用更多的内存（由于 .. 多了 2 个容器，而案例 2 中只有一个 “”）。IndexedVarcharging_profile_ev1charging_profile_ev3IndexedVarcharging_profile

最后 - 这纯粹是风格上的 - 您只需要对案例 1 使用 /（因为组件名称是以编程方式生成的）。在案例 2 中，您可以使用以下任一方式声明它们add_componentcomponent

def expression_soc(m, ts, ev):
    charging_power_kw = m.charging_profiles[ts, ev]
    driving_distance_km = ev.driving_profile_km.loc[ts]
    ev.charge(charging_power_kw)
    ev.drive(driving_distance_km)
    return ev.soc
model.soc_expression = pyo.Expression(
    model.time_steps, model.electric_vehicles,
    rule=expression_soc,
)

def constraint_soc(m, ts, ev):
    soc_value = m.soc_expression[ts, ev]
    return 0.5, soc_value, 1.0
model.soc_constraint = pyo.Constraint(
    model.time_steps, model.electric_vehicles,
    rule=constraint_soc
)

或

@model.Expression(model.time_steps, model.electric_vehicles)
def soc_expression(m, ts, ev):
    charging_power_kw = m.charging_profiles[ts, ev]
    driving_distance_km = ev.driving_profile_km.loc[ts]
    ev.charge(charging_power_kw)
    ev.drive(driving_distance_km)
    return ev.soc

@model.Constraint(model.time_steps, model.electric_vehicles)
def soc_constraint(m, ts, ev):
    soc_value = m.soc_expression[ts, ev]
    return 0.5, soc_value, 1.0