Pyomo：多 objective 优化
Pyomo: Multi-objective optimization
我正在使用此代码解决多 objective 优化模型（功率分配）并尝试在我的代码中调整示例。
example:
我试图跳过 'inefficient Pareto-front' 部分并直接绘制 'efficient Pareto-front'。
第一个标签可以运行正常生成Cost_min,Cost_max,Emission_min,Emission_max.
from pyomo.environ import *
import matplotlib.pyplot as plt
import numpy as np
import matplotlib.pyplot as plt
import random
# create a model
model = AbstractModel()
# declare decision variables
model.N = Param(mutable=True)
model.J = RangeSet(model.N)
model.A = Param(model.J)
model.B = Param(model.J)
model.C = Param(model.J)
model.D = Param(model.J)
model.E = Param(model.J)
model.F = Param(model.J)
model.P_min = Param(model.J, within=PositiveReals)
model.P_max = Param(model.J, within=PositiveReals)
model.demand = Param(mutable=True)
# declare constraints
def Pbounds(model, j):
return (model.P_min[j], model.P_max[j])
model.P = Var(model.J, bounds=Pbounds, domain=NonNegativeReals)
def P_LoadgenBalance(model):
return sum(model.P[j] for j in model.J) >= model.demand
model.P_LoadgenBalance = Constraint(rule=P_LoadgenBalance)
# declare objective_cost
def obj_cost(model):
return sum(model.A[j]* model.P[j] ** 2 + model.B[j] * model.P[j] + model.C[j] for j in model.J)
model.cost= Objective(rule=obj_cost, sense=minimize)
# declare objective_emission
def obj_emission(model):
return sum(model.E[j]* model.P[j] ** 2 + model.D[j] * model.P[j] + model.F[j] for j in model.J)
model.emission= Objective(rule=obj_emission, sense=minimize)
# deactivate model.emission calculate emission_max,cost_min
model.emission.deactivate()
instance = model.create_instance("E:\pycharm_project\PD\END-10units.dat")
opt = SolverFactory('Ipopt')
results = opt.solve(instance)
for i in instance.J:
print(i,value(instance.P[i]))
print( 'cost = ' + str(value(instance.cost)) )
print( 'emission = ' + str(value(instance.emission)) )
emission_max = value(instance.emission)
cost_min = value(instance.cost)
# ## max emission deactivate model.cost calculate emission_min,cost_max
model.emission.activate()
model.cost.deactivate()
instance = model.create_instance("E:\pycharm_project\PD\END-10units.dat")
results = opt.solve(instance)
for i in instance.J:
print(i,value(instance.P[i]))
print( 'cost = ' + str(value(instance.cost)) )
print( 'emission = ' + str(value(instance.emission)) )
emission_min = value(instance.emission)
cost_max = value(instance.cost)
运行在此选项卡中输入代码后，未生成任何错误。但是在输出Pareto-front的时候，这里只显示了一个点
# ## apply normal $\epsilon$-Constraint
model.emission.deactivate()
model.cost.activate()
model.emission_value = Param(initialize=0, mutable=True)
def c_epsilon(model):
return model.emission <= model.emission_value
model.C_epsilon = Constraint(rule=c_epsilon)
results = opt.solve(instance)
print('Each iteration will keep emission lower than some values between emission_min and emission_max, so [' + str(emission_min) + ', ' + str(emission_max) + ']')
n = 5
step = int((emission_max - emission_min) / n)
steps = list(range(int(emission_min), int(emission_max), step)) + [emission_max]
# ## apply augmented $\epsilon$-Constraint
# max emission + delta*epsilon <br>
# s.t. emission - s = emission_value
model.del_component(model.cost)
model.del_component(model.emission)
model.del_component(model.C_epsilon)
model.delta = Param(initialize=0.00001)
model.s = Var(within=NonNegativeReals)
def obj_cost_1(model):
return sum(model.cost+model.delta * model.s)
model.obj_cost_1 = Objective(rule=obj_cost_1, sense=maximize)
def C_e(model):
return model.emission-model.s==model.emission_value
model.C_e= Constraint(rule=C_e)
cost_l = []
emission_l = []
for i in steps:
model.emission_value = i
results = opt.solve(instance)
cost_l.append(value(instance.cost))
emission_l.append(value(instance.emission))
plt.plot(cost_l,emission_l,'o-.');
plt.title('efficient Pareto-front');
plt.grid(True);
plt.show()
结果如下图。不知道为什么这个不能输出正确的Paretochart.I不知道是哪一步代码错了
efficient Pareto-front
谁能帮我处理这段代码？
Thanks.Vivi
几件事.... :)
怎么了：
在您的循环内部，唯一影响模型的是您将新值分配给 model.e 那是什么？我认为这是一个打字错误，您错误地只是声明了一个名为 e 的新 和未使用的 模型组件实例变量。这就是为什么你没有得到不同的值。我想你想改成 model.emission.
此外，我不会一开始就尝试 1000 次，只尝试 5 次。
需要清理的内容：
您正在循环中实例化一个新求解器。没有必要。您不需要 1000 个不同的求解器，只需重新求解即可。您之前已经有一个求解器。
为了清晰起见，在您的代码中添加一些注释不会让您的手指着火，并且会有助于 T/S，以及一些重新组织。
此外，model.A model.B model.C ... 提供的信息不多。如果可以的话，我会建议更清晰的变量名。
我正在使用此代码解决多 objective 优化模型（功率分配）并尝试在我的代码中调整示例。
example: 我试图跳过 'inefficient Pareto-front' 部分并直接绘制 'efficient Pareto-front'。 第一个标签可以运行正常生成Cost_min,Cost_max,Emission_min,Emission_max. 运行在此选项卡中输入代码后，未生成任何错误。但是在输出Pareto-front的时候，这里只显示了一个点 结果如下图。不知道为什么这个不能输出正确的Paretochart.I不知道是哪一步代码错了 efficient Pareto-front
谁能帮我处理这段代码？
Thanks.Vivifrom pyomo.environ import *
import matplotlib.pyplot as plt
import numpy as np
import matplotlib.pyplot as plt
import random
# create a model
model = AbstractModel()
# declare decision variables
model.N = Param(mutable=True)
model.J = RangeSet(model.N)
model.A = Param(model.J)
model.B = Param(model.J)
model.C = Param(model.J)
model.D = Param(model.J)
model.E = Param(model.J)
model.F = Param(model.J)
model.P_min = Param(model.J, within=PositiveReals)
model.P_max = Param(model.J, within=PositiveReals)
model.demand = Param(mutable=True)
# declare constraints
def Pbounds(model, j):
return (model.P_min[j], model.P_max[j])
model.P = Var(model.J, bounds=Pbounds, domain=NonNegativeReals)
def P_LoadgenBalance(model):
return sum(model.P[j] for j in model.J) >= model.demand
model.P_LoadgenBalance = Constraint(rule=P_LoadgenBalance)
# declare objective_cost
def obj_cost(model):
return sum(model.A[j]* model.P[j] ** 2 + model.B[j] * model.P[j] + model.C[j] for j in model.J)
model.cost= Objective(rule=obj_cost, sense=minimize)
# declare objective_emission
def obj_emission(model):
return sum(model.E[j]* model.P[j] ** 2 + model.D[j] * model.P[j] + model.F[j] for j in model.J)
model.emission= Objective(rule=obj_emission, sense=minimize)
# deactivate model.emission calculate emission_max,cost_min
model.emission.deactivate()
instance = model.create_instance("E:\pycharm_project\PD\END-10units.dat")
opt = SolverFactory('Ipopt')
results = opt.solve(instance)
for i in instance.J:
print(i,value(instance.P[i]))
print( 'cost = ' + str(value(instance.cost)) )
print( 'emission = ' + str(value(instance.emission)) )
emission_max = value(instance.emission)
cost_min = value(instance.cost)
# ## max emission deactivate model.cost calculate emission_min,cost_max
model.emission.activate()
model.cost.deactivate()
instance = model.create_instance("E:\pycharm_project\PD\END-10units.dat")
results = opt.solve(instance)
for i in instance.J:
print(i,value(instance.P[i]))
print( 'cost = ' + str(value(instance.cost)) )
print( 'emission = ' + str(value(instance.emission)) )
emission_min = value(instance.emission)
cost_max = value(instance.cost)
# ## apply normal $\epsilon$-Constraint
model.emission.deactivate()
model.cost.activate()
model.emission_value = Param(initialize=0, mutable=True)
def c_epsilon(model):
return model.emission <= model.emission_value
model.C_epsilon = Constraint(rule=c_epsilon)
results = opt.solve(instance)
print('Each iteration will keep emission lower than some values between emission_min and emission_max, so [' + str(emission_min) + ', ' + str(emission_max) + ']')
n = 5
step = int((emission_max - emission_min) / n)
steps = list(range(int(emission_min), int(emission_max), step)) + [emission_max]
# ## apply augmented $\epsilon$-Constraint
# max emission + delta*epsilon <br>
# s.t. emission - s = emission_value
model.del_component(model.cost)
model.del_component(model.emission)
model.del_component(model.C_epsilon)
model.delta = Param(initialize=0.00001)
model.s = Var(within=NonNegativeReals)
def obj_cost_1(model):
return sum(model.cost+model.delta * model.s)
model.obj_cost_1 = Objective(rule=obj_cost_1, sense=maximize)
def C_e(model):
return model.emission-model.s==model.emission_value
model.C_e= Constraint(rule=C_e)
cost_l = []
emission_l = []
for i in steps:
model.emission_value = i
results = opt.solve(instance)
cost_l.append(value(instance.cost))
emission_l.append(value(instance.emission))
plt.plot(cost_l,emission_l,'o-.');
plt.title('efficient Pareto-front');
plt.grid(True);
plt.show()
几件事.... :)
怎么了：
在您的循环内部，唯一影响模型的是您将新值分配给 model.e 那是什么？我认为这是一个打字错误，您错误地只是声明了一个名为 e 的新 和未使用的 模型组件实例变量。这就是为什么你没有得到不同的值。我想你想改成 model.emission.
此外，我不会一开始就尝试 1000 次，只尝试 5 次。
需要清理的内容：
您正在循环中实例化一个新求解器。没有必要。您不需要 1000 个不同的求解器，只需重新求解即可。您之前已经有一个求解器。
为了清晰起见，在您的代码中添加一些注释不会让您的手指着火，并且会有助于 T/S，以及一些重新组织。
此外，model.A model.B model.C ... 提供的信息不多。如果可以的话，我会建议更清晰的变量名。