Pyomo 访问/检索双变量 - 具有二元变量的影子价格



我对一般的优化,特别是pyomo很陌生,所以我提前为任何新手错误道歉。

我以 [2] 为起点定义了一个简单的单元承诺练习(示例 3.1 来自 [1])。我得到了正确的结果并且我的代码运行了,但我有几个关于如何访问东西的问题。

import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
import shutil
import sys
import os.path
import pyomo.environ as pyo
import pyomo.gdp as gdp #necessary if you use booleans to select active and innactive units
def bounds_rule(m, n, param='Cap_MW'):
# m because it pases the module
# n because it needs a variable from each set, in this case there was only m.N
return (0, Gen[n][param]) #returns lower and upper bounds.
def unit_commitment():
m = pyo.ConcreteModel()
m.dual = pyo.Suffix(direction=pyo.Suffix.IMPORT_EXPORT)
N=Gen.keys()
m.N = pyo.Set(initialize=N)    

m.Pgen = pyo.Var(m.N, bounds = bounds_rule) #amount of generation
m.Rgen = pyo.Var(m.N, bounds = bounds_rule) #amount of generation
# m.OnOff = pyo.Var(m.N, domain=pyo.Binary) #boolean on/off marker

# objective
m.cost = pyo.Objective(expr = sum( m.Pgen[n]*Gen[n]['energy_$MWh'] + m.Rgen[n]*Gen[n]['res_$MW'] for n in m.N), sense=pyo.minimize)

# demand
m.demandP = pyo.Constraint(rule=lambda m: sum(m.Pgen[n] for n in N) == Demand['ener_MWh'])
m.demandR = pyo.Constraint(rule=lambda m: sum(m.Rgen[n] for n in N) == Demand['res_MW'])

# machine production limits
# m.lb = pyo.Constraint(m.N, rule=lambda m, n: Gen[n]['Cap_min']*m.OnOff[n] <= m.Pgen[n]+m.Rgen[n] )
# m.ub = pyo.Constraint(m.N, rule=lambda m, n: Gen[n]['Cap_MW']*m.OnOff[n] >= m.Pgen[n]+m.Rgen[n])
m.lb = pyo.Constraint(m.N, rule=lambda m, n: Gen[n]['Cap_min'] <= m.Pgen[n]+m.Rgen[n] )
m.ub = pyo.Constraint(m.N, rule=lambda m, n: Gen[n]['Cap_MW'] >= m.Pgen[n]+m.Rgen[n])
m.rc = pyo.Suffix(direction=pyo.Suffix.IMPORT)

return m
Gen = {
'GenA' : {'Cap_MW': 100, 'energy_$MWh': 10, 'res_$MW': 0, 'Cap_min': 0},
'GenB' : {'Cap_MW': 100, 'energy_$MWh': 30, 'res_$MW': 25, 'Cap_min': 0},
} #starting data
Demand = {
'ener_MWh': 130, 'res_MW': 20,
} #starting data
m = unit_commitment()
pyo.SolverFactory('glpk').solve(m).write() 
m.display()

df = pd.DataFrame.from_dict([m.Pgen.extract_values(), m.Rgen.extract_values()]).T.rename(columns={0: "P", 1: "R"})
print(df)
print("Cost Function result: " + str(m.cost.expr()) + "$.")
print(m.rc.display())
print(m.dual.display())
print(m.dual[m.demandR])
da= {'duals': m.dual[m.demandP],
'uslack': m.demandP.uslack(),
'lslack': m.demandP.lslack(),
}
db= {'duals': m.dual[m.demandR],
'uslack': m.demandR.uslack(),
'lslack': m.demandR.lslack(),
}
duals = pd.DataFrame.from_dict([da, db]).T.rename(columns={0: "demandP", 1: "demandR"})
print(duals)

我的问题来了。

1-对偶/影子价格:根据定义,影子价格是约束的对偶变量(m.demandP和m.demandR)。有没有办法访问这些值并将它们放入数据帧中,而无需做我所做的"糟糕"的事情?我的意思是手动定义 da 和 db,然后在两个字典连接时创建数据帧?我想做一些更干净的事情,比如 df 保存系统中每个生成器的 P 和 R 的结果。

2-通常,在单位承诺问题中,二进制变量用于"标记"或"选择"活动和非活动单位。因此,"m.OnOff"变量(注释行)。对于我在 [3] 中发现的内容,对偶在包含二进制变量的模型中不存在。之后,我重写了问题,但不包括二进制文件。在这种所有单位都运行的简单练习中,这不是问题,但对于较大的单位来说,这不是问题。我需要能够让优化决定哪些单位会运行,哪些单位不会运行,我仍然需要影子价格。有没有办法在包含二进制变量的问题中获得影子价格/对偶? 我让基于二进制变量的约束定义也在那里,以防有人发现它有用。

注意:代码也使用二进制变量运行并获得正确的结果,但是我不知道如何获得影子价格。因此,我的问题。

[1] Morales, J. M.,Conejo, A. J., Madsen, H., Pinson, P., & Zugno, M. (2013).将可再生能源纳入电力市场:运营问题(第205卷)。施普林格科学与商业媒体。

[2] https://jckantor.github.io/ND-Pyomo-Cookbook/04.06-Unit-Commitment.html

[3] 对偶变量在 Pyomo 中不返回任何内容

要回答 1,您可以使用model.component_objects(pyo.Constraint)从模型中动态获取约束对象,这将返回约束的迭代器,这样就不必对约束名称进行硬编码。对于索引变量来说,这变得很棘手,因为您必须执行额外的步骤来获取每个索引的松弛,而不仅仅是约束对象。对于对偶,可以循环访问keys属性以检索这些值。

duals_dict = {str(key):m.dual[key] for key in m.dual.keys()}
u_slack_dict = {
# uslacks for non-indexed constraints
**{str(con):con.uslack() for con in m.component_objects(pyo.Constraint)
if not con.is_indexed()},
# indexed constraint uslack
# loop through the indexed constraints
# get all the indices then retrieve the slacks for each index of constraint
**{k:v for con in m.component_objects(pyo.Constraint) if con.is_indexed()
for k,v in {'{}[{}]'.format(str(con),key):con[key].uslack()
for key in con.keys()}.items()}
}
l_slack_dict = {
# lslacks for non-indexed constraints
**{str(con):con.lslack() for con in m.component_objects(pyo.Constraint)
if not con.is_indexed()},
# indexed constraint lslack
# loop through the indexed constraints
# get all the indices then retrieve the slacks for each index of constraint
**{k:v for con in m.component_objects(pyo.Constraint) if con.is_indexed()
for k,v in {'{}[{}]'.format(str(con),key):con[key].lslack()
for key in con.keys()}.items()}
}
# combine into a single df
df = pd.concat([pd.Series(d,name=name)
for name,d in {'duals':duals_dict,
'uslack':u_slack_dict,
'lslack':l_slack_dict}.items()],
axis='columns')

关于 2,我同意@Erwin关于使用二进制变量求解以获得最佳解决方案的评论,然后删除二进制限制,但将变量固定为最佳值以获得一些对偶变量值。

最新更新