使用 pyomo 和 pao 优化电池容量和操作计划。尝试两阶段随机优化
Optimizing battery capacity AND operational schedule with pyomo and pao. Attempts at two-Stage Stochastic Optimization
更新：添加了可能的答案。
我正在尝试编写一个 multiobjective pyomo 优化脚本，以最小化成本为基础优化电池容量和运行计划。它应该回答通过在价格低时充电并在价格低时放电来减少电费的电池的最佳尺寸是多少，例如。在价格高时尽量减少电网采购。当给定电池容量时，我已经能够优化操作时间表。但是，我在实施容量优化时遇到问题。给出建筑物的每小时负荷曲线（用电量）和每小时购电价格（load[i] 和 price[i]）。成本与电池容量的每千瓦时价格是常数。
我正在使用 ipopt 求解器。
我尝试过不同的方法。第一个确实提供了最佳时间表和容量。但是，我不确定这是否是最佳容量的正确方法，因为它是在相同的 objective 函数中实现的，该函数最小化电池成本 + 购买电力的成本 - 从电池释放的能量成本。
m = AbstractModel()
m.cap = Var(domain=Reals, bounds=(10,1000)) #battery capacity
m.pe_c = Var(time,domain=NonNegativeReals,bounds=(0,1000)) #charge to battery
m.pe_d = Var(time,domain=NegativeReals,bounds=(-1000,0)) #discharge from battery
m.soc = Var(soc_time,bounds=(value(min_soc),value(max_soc))) #state of charge between 20% and 95%
def obj_expression(m):
return ((m.cap *cost)+sum([price[i]*(load[i]) for i in time])) + sum([price[i]*(m.pe_c[i]) for i in time]) + sum([price[i]*(m.pe_d[i]) for i in time])
m.OBJ = Objective(rule=obj_expression, sense=minimize)
def soc_constraint_rule(m,i):
return m.soc[i+1] - m.soc[i] <= float(100) * dt * (m.pe_d[i] + m.pe_c[i]*energy_efficiency + self_discharge_power) / m.cap
m.soc_constraints = Constraint(time, rule=soc_constraint_rule)
def soc_start_rule(m): #battery starts at initial value
return m.soc[0] == soc_init
m.soc_start = Constraint(rule=soc_start_rule)
def soc_end_rule(m): #battery returns to initial value at end of simulation
return m.soc[n] == m.soc[0]
m.soc_end = Constraint(rule=soc_end_rule)
def discharge_rule(m,i): # Cannot discharge more than what battery can hold
return m.pe_d[i] <= m.cap *C_rate #C_rate is an efficiency constant, in this case 0.9
m.discharge_rule = Constraint(time, rule=discharge_rule)
def charge_rule(m,i): #cannot charge more than battery can hold
return m.pe_c[i] <= m.cap *C_rate
m.charge_rule = Constraint(time, rule=charge_rule)
instance = m.create_instance()
solver = SolverFactory('ipopt')
solver.solve(instance,tee=True)
我的第二种方法无法正常工作以优化运营计划以最大限度地降低从电网购买能源的成本，但它确实 return 运营计划和容量。在这种方法中，我尝试使用能量状态 (SOE) 而不是 SOC，使容量成为 SOE 的限制，这应该避免 SOC 约束中的非凸性可能导致 ipopt 问题收敛。这个没有 return m.utility 值，我通常不确定 objective 函数是否正确，是否可以避免同时发生充电和放电，这在实际中是不可能的-生活场景
m = AbstractModel()
m.cap = Var(domain=Reals, bounds=(10,1000)) #battery capacity variable
m.utility = Var(time,domain=NonNegativeReals) #purchase from grid
m.pe_c = Var(time,domain=NonNegativeReals,bounds=(0,1000)) #charge to battery
m.pe_d = Var(time,domain=NegativeReals,bounds=(-1000,0)) #discharge from battery
m.SOE =Var(soe_time, domain=NonNegativeReals,bounds=(0,1000)) #state of energy inside battery
def obj_expression(m):
return (m.cap*cost) + sum([price[i]*(m.utility[i]) for i in time])
m.OBJ = Objective(rule=obj_expression, sense=minimize)
def battery_rule(m,i): #to connect the states of the amount of energy stored at, energy discharged from, energy charged to, the battery system
return m.SOE[i+1] == m.SOE[i] + m.pe_c[i] + m.pe_d[i]
m.battery_rule = Constraint(time, rule=battery_rule)
def capacity_rule(m,i): # a physical constraint that limits the amount of energy stored in the battery system
return m.SOE[i+1] <= m.cap
m.capacity_rule = Constraint(time, rule=capacity_rule)
def charging_rule(m,i): # a physical constraint that limits the amount of energy charged to the battery system
return m.pe_c[i] <= m.cap * C_rate
m.charging_rule = Constraint(time, rule=charging_rule)
def discharging_rule(m,i): # a physical constraint that limits the amount of energy discharged from the battery system
return m.pe_d[i] <= m.cap * C_rate
m.discharging_rule = Constraint(time, rule=discharging_rule)
def soe_start_rule(m):
return m.SOE[0] == soe_init
m.soe_start_rule = Constraint(rule=soe_start_rule)
def soe_end_rule(m):
return m.SOE[n] == m.SOE[0]
m.soe_end = Constraint(rule=soe_end_rule)
def demand_rule(m,i): #a physical constraint, representing the law of conservation of energy
return m.pe_d[i] + m.pe_c[i] + m.utility[i] == load[i]
m.demand_rule = Constraint(time, rule=demand_rule)
instance = m.create_instance()
solver = SolverFactory('ipopt')
solver.solve(instance,tee=True)
我的最终方法是使用 PAO 库访问双层优化，但根本不起作用。这可能是最好的方法，因为它清楚地定义了上层和下层 objective，但我在使用子模型构建迭代约束时遇到了问题。下面的代码 returns “TypeError: 'NoneType' object is not iterable”。但是，我确信由于人为错误，这种方法存在较大的逻辑错误。我还尝试通过交换变量 m 将容量 * 成本作为上层 objective。 m.sub;将 m.sub.cap 制作成 m.cap 并将 m.utility 制作成 m.sub.utility 等，但是 returns “AttributeError: 'NoneType' object has no属性 'transpose'" 当我 运行 求解器时。
m = AbstractModel()
m.sub = pao.SubModel()
m.sub.cap = Var(domain=Reals, bounds=(10,1000)) #battery capacity variable
m.utility = Var(time,domain=NonNegativeReals) #purchase from grid
m.pe_c = Var(time,domain=NonNegativeReals,bounds=(0,1000)) #charge to battery
m.pe_d = Var(time,domain=NegativeReals,bounds=(-1000,0)) #discharge from battery
m.SOE =Var(soe_time, domain=NonNegativeReals,bounds=(0,1000)) #state of energy inside battery
def obj_expression(m):
return sum([price[i]*(m.pe_c[i]) for i in time]) +sum([price[i]*(m.pe_d[i]) for i in time])
m.OBJ = Objective(rule=obj_expression, sense=minimize)
m.sub.obj = Objective(expr=m.sub.cap*cost, sense=minimize)
def battery_rule(m,i): #to connect the states of the amount of energy stored at, energy discharged from, energy charged to, the battery system
return m.SOE[i+1] == m.SOE[i] + m.pe_c[i] + m.pe_d[i]
m.battery_rule = Constraint(time, rule=battery_rule)
def capacity_rule(m,i): # a physical constraint that limits the amount of energy stored in the battery system
return m.SOE[i+1] <= m.sub.cap
m.capacity_rule = Constraint(time, rule=capacity_rule)
def charging_rule(m,i): # a physical constraint that limits the amount of energy charged to the battery system
return m.pe_c[i] <= m.sub.cap * C_rate
m.charging_rule = Constraint(time, rule=charging_rule)
def discharging_rule(m,i): # a physical constraint that limits the amount of energy discharged from the battery system
return m.pe_d[i] <= m.sub.cap * C_rate
m.discharging_rule = Constraint(time, rule=discharging_rule)
def soe_start_rule(m):
return m.SOE[0] == soe_init
m.soe_start_rule = Constraint(rule=soe_start_rule)
def soe_end_rule(m):
return m.SOE[n] == m.SOE[0]
m.soe_end = Constraint(rule=soe_end_rule)
def demand_rule(m,i): #a physical constraint, representing the law of conservation of energy
return m.pe_d[i] + m.pe_c[i] + m.utility[i] == load[i]
m.demand_rule = Constraint(time, rule=demand_rule)
nlp = pao.Solver('ipopt')
opt = pao.Solver('pao.pyomo.REG', nlp_solver=nlp)
results = opt.solve(m)
我写这篇文章 post 希望有人做过类似的多 objective 两阶段随机优化并且能够看到我的人为错误或给我反馈什么方法是更好的。如果有什么不清楚的地方请问我，我可以提供更多细节。
更新：我正在研究两个可能对上述模型有帮助也可能没有帮助的新约束。他们在下面
def cost_reduction(m,i): #making sure the without-battery cost of electricity is higher than the with-battery cost of electricity.
return (m.pe_c[i]*price[i]) + (m.pe_d[i]*price[i]) + (load[i]*price[i]) <= (load[i]*price[i])
m.cost_reduction = Constraint(time, rule=cost_reduction)
def battery_cost(m,i): #make sure recover battery cost M.pe_d is a negative number and m.year is bounded between 0,20.
return (m.cap*cost) + (m.pe_d[i]*price[i]) *m.year ==0 #I am aware that this is not the correct way but it shows the idea of a constraint that attempts to ensure a payback
m.battery_cost = Constraint(time, rule=battery_cost)
我相信我已经使用 gurobi 求解器解决了它。必须更改 demand_rule 并添加 binary_charge_constraint
m = ConcreteModel()
m.cap = Var(domain=Reals, bounds=(10,1000))#battery capacity variable
m.utility = Var(time,domain=NonNegativeReals) #purchase from grid
m.pe_c = Var(time,domain=NonNegativeReals,bounds=(0,1000)) #charge to battery
m.pe_d = Var(time,domain=NonNegativeReals,bounds=(0,1000)) #discharge from battery
m.SOE =Var(soe_time, domain=NonNegativeReals,bounds=(0,1000)) #state of energy inside battery
def obj_expression(m):
return (m.cap*cost) + sum([price[i]*(m.utility[i]) for i in time])
m.OBJ = Objective(rule=obj_expression, sense=minimize)
def battery_rule(m,i): #to connect the states of the amount of energy stored at, energy discharged from, energy charged to, the battery system
return m.SOE[i+1] == m.SOE[i] + m.pe_c[i] - m.pe_d[i]
m.battery_rule = Constraint(time, rule=battery_rule)
def capacity_rule(m,i): # a physical constraint that limits the amount of energy stored in the battery system
return m.SOE[i] <= m.cap
m.capacity_rule = Constraint(time, rule=capacity_rule)
def charging_rule(m,i): # a physical constraint that limits the amount of energy charged to the battery system
return m.pe_c[i] <= m.cap * C_rate
m.charging_rule = Constraint(time, rule=charging_rule)
def discharging_rule(m,i): # a physical constraint that limits the amount of energy discharged from the battery system
return m.pe_d[i] <= m.cap * C_rate
m.discharging_rule = Constraint(time, rule=discharging_rule)
def soe_start_rule(m):
return m.SOE[0] == soe_init
m.soe_start_rule = Constraint(rule=soe_start_rule)
def soe_end_rule(m):
return m.SOE[n] == m.SOE[0]
m.soe_end = Constraint(rule=soe_end_rule)
def demand_rule(m,i): #a physical constraint, representing the law of conservation of energy
return m.pe_d[i] - m.pe_c[i] + m.utility[i] == load[i]
m.demand_rule = Constraint(time, rule=demand_rule)
def binary_charge_constraint(m,i): #making sure it is not charging while discharging and viceversa
return m.pe_d[i]*m.pe_c[i]== 0
m.binary_charge_constraint = Constraint(time, rule=binary_charge_constraint)
instance = m.create_instance()
opt = SolverFactory("gurobi")
opt.options['NonConvex'] = 2
results = opt.solve(instance, tee=True)
更新：添加了可能的答案。
我正在尝试编写一个 multiobjective pyomo 优化脚本，以最小化成本为基础优化电池容量和运行计划。它应该回答通过在价格低时充电并在价格低时放电来减少电费的电池的最佳尺寸是多少，例如。在价格高时尽量减少电网采购。当给定电池容量时，我已经能够优化操作时间表。但是，我在实施容量优化时遇到问题。给出建筑物的每小时负荷曲线（用电量）和每小时购电价格（load[i] 和 price[i]）。成本与电池容量的每千瓦时价格是常数。
我正在使用 ipopt 求解器。
我尝试过不同的方法。第一个确实提供了最佳时间表和容量。但是，我不确定这是否是最佳容量的正确方法，因为它是在相同的 objective 函数中实现的，该函数最小化电池成本 + 购买电力的成本 - 从电池释放的能量成本。
m = AbstractModel()
m.cap = Var(domain=Reals, bounds=(10,1000)) #battery capacity
m.pe_c = Var(time,domain=NonNegativeReals,bounds=(0,1000)) #charge to battery
m.pe_d = Var(time,domain=NegativeReals,bounds=(-1000,0)) #discharge from battery
m.soc = Var(soc_time,bounds=(value(min_soc),value(max_soc))) #state of charge between 20% and 95%
def obj_expression(m):
return ((m.cap *cost)+sum([price[i]*(load[i]) for i in time])) + sum([price[i]*(m.pe_c[i]) for i in time]) + sum([price[i]*(m.pe_d[i]) for i in time])
m.OBJ = Objective(rule=obj_expression, sense=minimize)
def soc_constraint_rule(m,i):
return m.soc[i+1] - m.soc[i] <= float(100) * dt * (m.pe_d[i] + m.pe_c[i]*energy_efficiency + self_discharge_power) / m.cap
m.soc_constraints = Constraint(time, rule=soc_constraint_rule)
def soc_start_rule(m): #battery starts at initial value
return m.soc[0] == soc_init
m.soc_start = Constraint(rule=soc_start_rule)
def soc_end_rule(m): #battery returns to initial value at end of simulation
return m.soc[n] == m.soc[0]
m.soc_end = Constraint(rule=soc_end_rule)
def discharge_rule(m,i): # Cannot discharge more than what battery can hold
return m.pe_d[i] <= m.cap *C_rate #C_rate is an efficiency constant, in this case 0.9
m.discharge_rule = Constraint(time, rule=discharge_rule)
def charge_rule(m,i): #cannot charge more than battery can hold
return m.pe_c[i] <= m.cap *C_rate
m.charge_rule = Constraint(time, rule=charge_rule)
instance = m.create_instance()
solver = SolverFactory('ipopt')
solver.solve(instance,tee=True)
我的第二种方法无法正常工作以优化运营计划以最大限度地降低从电网购买能源的成本，但它确实 return 运营计划和容量。在这种方法中，我尝试使用能量状态 (SOE) 而不是 SOC，使容量成为 SOE 的限制，这应该避免 SOC 约束中的非凸性可能导致 ipopt 问题收敛。这个没有 return m.utility 值，我通常不确定 objective 函数是否正确，是否可以避免同时发生充电和放电，这在实际中是不可能的-生活场景
m = AbstractModel()
m.cap = Var(domain=Reals, bounds=(10,1000)) #battery capacity variable
m.utility = Var(time,domain=NonNegativeReals) #purchase from grid
m.pe_c = Var(time,domain=NonNegativeReals,bounds=(0,1000)) #charge to battery
m.pe_d = Var(time,domain=NegativeReals,bounds=(-1000,0)) #discharge from battery
m.SOE =Var(soe_time, domain=NonNegativeReals,bounds=(0,1000)) #state of energy inside battery
def obj_expression(m):
return (m.cap*cost) + sum([price[i]*(m.utility[i]) for i in time])
m.OBJ = Objective(rule=obj_expression, sense=minimize)
def battery_rule(m,i): #to connect the states of the amount of energy stored at, energy discharged from, energy charged to, the battery system
return m.SOE[i+1] == m.SOE[i] + m.pe_c[i] + m.pe_d[i]
m.battery_rule = Constraint(time, rule=battery_rule)
def capacity_rule(m,i): # a physical constraint that limits the amount of energy stored in the battery system
return m.SOE[i+1] <= m.cap
m.capacity_rule = Constraint(time, rule=capacity_rule)
def charging_rule(m,i): # a physical constraint that limits the amount of energy charged to the battery system
return m.pe_c[i] <= m.cap * C_rate
m.charging_rule = Constraint(time, rule=charging_rule)
def discharging_rule(m,i): # a physical constraint that limits the amount of energy discharged from the battery system
return m.pe_d[i] <= m.cap * C_rate
m.discharging_rule = Constraint(time, rule=discharging_rule)
def soe_start_rule(m):
return m.SOE[0] == soe_init
m.soe_start_rule = Constraint(rule=soe_start_rule)
def soe_end_rule(m):
return m.SOE[n] == m.SOE[0]
m.soe_end = Constraint(rule=soe_end_rule)
def demand_rule(m,i): #a physical constraint, representing the law of conservation of energy
return m.pe_d[i] + m.pe_c[i] + m.utility[i] == load[i]
m.demand_rule = Constraint(time, rule=demand_rule)
instance = m.create_instance()
solver = SolverFactory('ipopt')
solver.solve(instance,tee=True)
我的最终方法是使用 PAO 库访问双层优化，但根本不起作用。这可能是最好的方法，因为它清楚地定义了上层和下层 objective，但我在使用子模型构建迭代约束时遇到了问题。下面的代码 returns “TypeError: 'NoneType' object is not iterable”。但是，我确信由于人为错误，这种方法存在较大的逻辑错误。我还尝试通过交换变量 m 将容量 * 成本作为上层 objective。 m.sub;将 m.sub.cap 制作成 m.cap 并将 m.utility 制作成 m.sub.utility 等，但是 returns “AttributeError: 'NoneType' object has no属性 'transpose'" 当我 运行 求解器时。
m = AbstractModel()
m.sub = pao.SubModel()
m.sub.cap = Var(domain=Reals, bounds=(10,1000)) #battery capacity variable
m.utility = Var(time,domain=NonNegativeReals) #purchase from grid
m.pe_c = Var(time,domain=NonNegativeReals,bounds=(0,1000)) #charge to battery
m.pe_d = Var(time,domain=NegativeReals,bounds=(-1000,0)) #discharge from battery
m.SOE =Var(soe_time, domain=NonNegativeReals,bounds=(0,1000)) #state of energy inside battery
def obj_expression(m):
return sum([price[i]*(m.pe_c[i]) for i in time]) +sum([price[i]*(m.pe_d[i]) for i in time])
m.OBJ = Objective(rule=obj_expression, sense=minimize)
m.sub.obj = Objective(expr=m.sub.cap*cost, sense=minimize)
def battery_rule(m,i): #to connect the states of the amount of energy stored at, energy discharged from, energy charged to, the battery system
return m.SOE[i+1] == m.SOE[i] + m.pe_c[i] + m.pe_d[i]
m.battery_rule = Constraint(time, rule=battery_rule)
def capacity_rule(m,i): # a physical constraint that limits the amount of energy stored in the battery system
return m.SOE[i+1] <= m.sub.cap
m.capacity_rule = Constraint(time, rule=capacity_rule)
def charging_rule(m,i): # a physical constraint that limits the amount of energy charged to the battery system
return m.pe_c[i] <= m.sub.cap * C_rate
m.charging_rule = Constraint(time, rule=charging_rule)
def discharging_rule(m,i): # a physical constraint that limits the amount of energy discharged from the battery system
return m.pe_d[i] <= m.sub.cap * C_rate
m.discharging_rule = Constraint(time, rule=discharging_rule)
def soe_start_rule(m):
return m.SOE[0] == soe_init
m.soe_start_rule = Constraint(rule=soe_start_rule)
def soe_end_rule(m):
return m.SOE[n] == m.SOE[0]
m.soe_end = Constraint(rule=soe_end_rule)
def demand_rule(m,i): #a physical constraint, representing the law of conservation of energy
return m.pe_d[i] + m.pe_c[i] + m.utility[i] == load[i]
m.demand_rule = Constraint(time, rule=demand_rule)
nlp = pao.Solver('ipopt')
opt = pao.Solver('pao.pyomo.REG', nlp_solver=nlp)
results = opt.solve(m)
我写这篇文章 post 希望有人做过类似的多 objective 两阶段随机优化并且能够看到我的人为错误或给我反馈什么方法是更好的。如果有什么不清楚的地方请问我，我可以提供更多细节。
更新：我正在研究两个可能对上述模型有帮助也可能没有帮助的新约束。他们在下面
def cost_reduction(m,i): #making sure the without-battery cost of electricity is higher than the with-battery cost of electricity.
return (m.pe_c[i]*price[i]) + (m.pe_d[i]*price[i]) + (load[i]*price[i]) <= (load[i]*price[i])
m.cost_reduction = Constraint(time, rule=cost_reduction)
def battery_cost(m,i): #make sure recover battery cost M.pe_d is a negative number and m.year is bounded between 0,20.
return (m.cap*cost) + (m.pe_d[i]*price[i]) *m.year ==0 #I am aware that this is not the correct way but it shows the idea of a constraint that attempts to ensure a payback
m.battery_cost = Constraint(time, rule=battery_cost)
我相信我已经使用 gurobi 求解器解决了它。必须更改 demand_rule 并添加 binary_charge_constraint
m = ConcreteModel()
m.cap = Var(domain=Reals, bounds=(10,1000))#battery capacity variable
m.utility = Var(time,domain=NonNegativeReals) #purchase from grid
m.pe_c = Var(time,domain=NonNegativeReals,bounds=(0,1000)) #charge to battery
m.pe_d = Var(time,domain=NonNegativeReals,bounds=(0,1000)) #discharge from battery
m.SOE =Var(soe_time, domain=NonNegativeReals,bounds=(0,1000)) #state of energy inside battery
def obj_expression(m):
return (m.cap*cost) + sum([price[i]*(m.utility[i]) for i in time])
m.OBJ = Objective(rule=obj_expression, sense=minimize)
def battery_rule(m,i): #to connect the states of the amount of energy stored at, energy discharged from, energy charged to, the battery system
return m.SOE[i+1] == m.SOE[i] + m.pe_c[i] - m.pe_d[i]
m.battery_rule = Constraint(time, rule=battery_rule)
def capacity_rule(m,i): # a physical constraint that limits the amount of energy stored in the battery system
return m.SOE[i] <= m.cap
m.capacity_rule = Constraint(time, rule=capacity_rule)
def charging_rule(m,i): # a physical constraint that limits the amount of energy charged to the battery system
return m.pe_c[i] <= m.cap * C_rate
m.charging_rule = Constraint(time, rule=charging_rule)
def discharging_rule(m,i): # a physical constraint that limits the amount of energy discharged from the battery system
return m.pe_d[i] <= m.cap * C_rate
m.discharging_rule = Constraint(time, rule=discharging_rule)
def soe_start_rule(m):
return m.SOE[0] == soe_init
m.soe_start_rule = Constraint(rule=soe_start_rule)
def soe_end_rule(m):
return m.SOE[n] == m.SOE[0]
m.soe_end = Constraint(rule=soe_end_rule)
def demand_rule(m,i): #a physical constraint, representing the law of conservation of energy
return m.pe_d[i] - m.pe_c[i] + m.utility[i] == load[i]
m.demand_rule = Constraint(time, rule=demand_rule)
def binary_charge_constraint(m,i): #making sure it is not charging while discharging and viceversa
return m.pe_d[i]*m.pe_c[i]== 0
m.binary_charge_constraint = Constraint(time, rule=binary_charge_constraint)
instance = m.create_instance()
opt = SolverFactory("gurobi")
opt.options['NonConvex'] = 2
results = opt.solve(instance, tee=True)