首页 / 视频教程

python数学建模求解,pyomo声明式建模语言的使用,以一道规划问题建模为例

发布时间:2023-11-20 07:03:46

以一道线性规划问题为例,使用Pyomo构建简单的规划模型,调用求解器进行求解。

线性规划、数学建模、Pyomo、cbc求解器

 ✦ 

Pyomo简介...

Pyomo(Python Optimization Modeling Objects)是一个用于建模和解决优化问题的Python库。通过声明式建模语言定义规划问题的框架,调用求解器(Gurobi,GLPK,CBC,CPLEX,IPOPT)进行求解。模块化是它的优点,创建和维护起来会更轻松。

 ✦ 

过程演示...

from pyomo.environ import *
# 创建一个简单的线性规划模型model = ConcreteModel()
# 定义变量model.x = Var(within=NonNegativeReals)
# 定义目标函数和约束model.obj = Objective(expr=model.x, sense=minimize)model.con = Constraint(expr=model.x >= 5)
# 选择求解器(例如GLPK)并求解模型solver = SolverFactory('glpk')solver.solve(model)
# 打印结果print("Variable x =", value(model.x))

 ✦ 

例题...

在未来的六个月内, 公司计划生产四种元件, 型号分别为 X43-M1, X43-M2, Y54-N1,Y54-N2。这些元件的产量受到产能变化的影响,并且每次产能改变后都需要重新进行控制和调整, 因此会带来不可忽略的费用。因此公司希望最小化这些改变带来的费用,以及生产和库存的成本。 

在下表中列出了每种产品每个时期内的需求量,生产和库存成本,初始库存量,以及最后希望保留的库存量。

另外,当产量发生变化时, 需要对设备进行调整。由此带来的费用与当月产量较前一月的改变量(提高或降低)成正比。产量每提高一个产品单位, 则需要多支出 1 欧元;产量每降低一个产品单位, 只需要多支出 0.50 欧元。

目标是最小化由于产量改变引起的费用,以及生产和库存成本,给出总成本及每个月每种产品的生产量。

 ✦ 

建模设计一...

集合

月份:Month = [1,2,3,4,5,6]

产品:Products = ['X43-M1','X43-M2','Y54-N1','Y54-N2']

参数

Demands_{i,j}:产品 i 在 j 月的需求量

ProductCost_i:产品 i 的生产成本

StorageCost_i:产品 i 的存储成本

InitialStor_i:产品 i 的期初库存

FinalStor_i:产品 i 的期末库存

定义变量

Production_{i,j}:产品 i 在 j 月的生产量, 非负整数

Storage_{i,j}:产品 i 在 j 月的存储量, 非负整数

Increase_j:j 月比上一月增加的生产量, 非负整数

Decrease_j:j 月比上一月减少的生产量, 非正整数

定义目标函数

成本最小:总生产成本+总存储成本+总调整成本

定义约束

1、生产库存平衡约束

2、生产量调整约束

 ✦ 

代码求解...

from pyomo.environ import *#参数Demands = [[1500, 3000, 2000, 4000, 2000, 2500],            [1300, 800, 800, 1000, 1100, 900],            [2200, 1500, 2900, 1800, 1200, 2100],            [1400, 1600, 1500, 1000, 1100, 1200]]ProductCost = [20, 25, 10, 15]StorageCost = [0.4, 0.5, 0.3, 0.3]InitialStor = [10, 0, 0, 0]FinalStor = [50, 10, 10, 10]#变量维度model = ConcreteModel()model.I = Set(initialize = [i for i in range(len(Demands))])model.J = Set(initialize = [i for i in range(len(Demands[0]))])model.J1 = Set(initialize = [i for i in range(len(Demands[0]) - 1)])#变量model.Production = Var(model.I, model.J, within = NonNegativeIntegers, bounds = (0,10000), initialize=0)model.Storage = Var(model.I, model.J, within = NonNegativeIntegers, bounds = (0,10000), initialize=0)model.Increase = Var(model.J1, within = NonNegativeIntegers, bounds = (0,10000), initialize=0)model.Decrease = Var(model.J1, within = NonPositiveIntegers, bounds = (-10000,0), initialize=0)#目标函数表达式model.VarProductCost = Expression(expr=sum(ProductCost[i] * model.Production[i, j] for i in model.I for j in model.J))model.VarChangeCost = Expression(expr=sum(model.Increase[j] - 0.5 * model.Decrease[j] for j in model.J1))model.VarStorageCost = Expression(expr=sum(StorageCost[i] * model.Storage[i, j] for i in model.I for j in model.J))#目标函数model.cost = Objective(expr=model.VarProductCost + model.VarChangeCost + model.VarStorageCost, sense=minimize)#生产库存平衡约束model.StorageCost = ConstraintList()for i in model.I:    for j in model.J:        if j == 0:            model.StorageCost.add(model.Storage[i, j] ==  InitialStor[i] + model.Production[i, j] - Demands[i][j])        elif j == len(model.J) - 1:            model.StorageCost.add(model.Storage[i, j] == model.Storage[i, j-1] + model.Production[i, j] - Demands[i][j])            model.StorageCost.add(model.Storage[i, j] == FinalStor[i])        else:            model.StorageCost.add(model.Storage[i, j] == model.Storage[i, j-1] + model.Production[i, j] - Demands[i][j])#生产量调整约束def diff_rule(model,j):    return sum(model.Production[i, j + 1] - model.Production[i, j] for i in model.I) == model.Increase[j] + model.Decrease[j]model.conCapOver = Constraint(model.J1, rule=diff_rule)
#求解opt = SolverFactory('cbc')result = opt.solve(model)
print('++++++++++++++++++++++++++++++++++++++')model.Production.pprint()print('++++++++++++++++++++++++++++++++++++++')model.Increase.pprint()print('++++++++++++++++++++++++++++++++++++++')model.Decrease.pprint()print('++++++++++++++++++++++++++++++++++++++')model.Storage.pprint()print('Finally Cost: ', model.cost())
result.write()

 ✦ 

结果展示...

Production : Size=24, Index=Production_index    Key    : Lower : Value  : Upper : Fixed : Stale : Domain    (0, 0) :     0 : 1490.0 : 10000 : False : False : NonNegativeIntegers    (0, 1) :     0 : 3000.0 : 10000 : False : False : NonNegativeIntegers    (0, 2) :     0 : 2000.0 : 10000 : False : False : NonNegativeIntegers    (0, 3) :     0 : 4000.0 : 10000 : False : False : NonNegativeIntegers    (0, 4) :     0 : 2000.0 : 10000 : False : False : NonNegativeIntegers    (0, 5) :     0 : 2550.0 : 10000 : False : False : NonNegativeIntegers    (1, 0) :     0 : 1300.0 : 10000 : False : False : NonNegativeIntegers    (1, 1) :     0 :  800.0 : 10000 : False : False : NonNegativeIntegers    (1, 2) :     0 :  800.0 : 10000 : False : False : NonNegativeIntegers    (1, 3) :     0 : 1000.0 : 10000 : False : False : NonNegativeIntegers    (1, 4) :     0 : 1100.0 : 10000 : False : False : NonNegativeIntegers    (1, 5) :     0 :  910.0 : 10000 : False : False : NonNegativeIntegers    (2, 0) :     0 : 2882.0 : 10000 : False : False : NonNegativeIntegers    (2, 1) :     0 : 1672.0 : 10000 : False : False : NonNegativeIntegers    (2, 2) :     0 : 2046.0 : 10000 : False : False : NonNegativeIntegers    (2, 3) :     0 : 1800.0 : 10000 : False : False : NonNegativeIntegers    (2, 4) :     0 : 1890.0 : 10000 : False : False : NonNegativeIntegers    (2, 5) :     0 : 1420.0 : 10000 : False : False : NonNegativeIntegers    (3, 0) :     0 : 1400.0 : 10000 : False : False : NonNegativeIntegers    (3, 1) :     0 : 1600.0 : 10000 : False : False : NonNegativeIntegers    (3, 2) :     0 : 2227.0 : 10000 : False : False : NonNegativeIntegers    (3, 3) :     0 :  273.0 : 10000 : False : False : NonNegativeIntegers    (3, 4) :     0 : 1100.0 : 10000 : False : False : NonNegativeIntegers    (3, 5) :     0 : 1210.0 : 10000 : False : False : NonNegativeIntegers# ==========================================================Finally Cost:  684209.4# ==========================================================Time: 0.04673957824707031

 ✦ 

建模设计二(大M法)...

参数

Demands_{i,j}:产品 i 在 j 月的需求量

ProductCost_i:产品 i 的生产成本

StorageCost_i:产品 i 的存储成本

InitialStor_i:产品 i 的期初库存

FinalStor_i:产品 i 的期末库存

M:极大正数

定义变量

Production_{i,j}:产品 i 在 j 月的生产量, 非负整数

Storage_{i,j}:产品 i 在 j 月的存储量, 非负整数

Change_j:j 月相对上一月调整成本, 非负整数

Sita_j:0-1变量,为1代表 j 月的下一月产量增加,为0减少

定义目标函数

成本最小:总生产成本+总存储成本+总调整成本

定义约束

1、生产库存平衡约束

2、生产量调整约束

 ✦ 

代码求解...

from pyomo.environ import *#参数Demands = [[1500, 3000, 2000, 4000, 2000, 2500],            [1300, 800, 800, 1000, 1100, 900],            [2200, 1500, 2900, 1800, 1200, 2100],            [1400, 1600, 1500, 1000, 1100, 1200]]ProductCost = [20, 25, 10, 15]StorageCost = [0.4, 0.5, 0.3, 0.3]InitialStor = [10, 0, 0, 0]FinalStor = [50, 10, 10, 10]#变量维度model = ConcreteModel()model.I = Set(initialize = [i for i in range(len(Demands))])model.J = Set(initialize = [i for i in range(len(Demands[0]))])model.J1 = Set(initialize = [i for i in range(len(Demands[0]) - 1)])#变量model.Production = Var(model.I, model.J, within = NonNegativeIntegers, bounds = (0,10000), initialize=0)model.Storage = Var(model.I, model.J, within = NonNegativeIntegers, bounds = (0,10000), initialize=0)model.Change = Var(model.J1, within = NonNegativeIntegers, bounds = (0,10000), initialize=0)model.sita = Var(model.J1, within=Binary, initialize=1)#目标函数表达式model.VarProductCost = Expression(expr=sum(ProductCost[i] * model.Production[i, j] for i in model.I for j in model.J))model.VarChangeCost = Expression(expr=sum(model.Change[j] for j in model.J1))model.VarStorageCost = Expression(expr=sum(StorageCost[i] * model.Storage[i, j] for i in model.I for j in model.J))#目标函数model.cost = Objective(expr=model.VarProductCost + model.VarChangeCost + model.VarStorageCost, sense=minimize)#生产库存平衡约束model.StorageCost = ConstraintList()for i in model.I:    for j in model.J:        if j == 0:            model.StorageCost.add(model.Storage[i, j] ==  InitialStor[i] + model.Production[i, j] - Demands[i][j])        elif j == len(model.J) - 1:            model.StorageCost.add(model.Storage[i, j] == model.Storage[i, j-1] + model.Production[i, j] - Demands[i][j])            model.StorageCost.add(model.Storage[i, j] == FinalStor[i])        else:            model.StorageCost.add(model.Storage[i, j] == model.Storage[i, j-1] + model.Production[i, j] - Demands[i][j])#生产量调整约束model.ChangeCost = ConstraintList()for j in model.J1:    model.ChangeCost.add(sum(model.Production[i, j] - model.Production[i, j + 1] for i in model.I) <= (1 - model.sita[j]) * 9999999)    model.ChangeCost.add(sum(model.Production[i, j] - model.Production[i, j + 1] for i in model.I) >= - model.sita[j] * 9999999)    model.ChangeCost.add(model.Change[j] >= sum(model.Production[i, j + 1] - model.Production[i, j] for i in model.I) - (1 - model.sita[j]) * 9999999)    model.ChangeCost.add(model.Change[j] <= sum(model.Production[i, j + 1] - model.Production[i, j] for i in model.I) + (1 - model.sita[j]) * 9999999)    model.ChangeCost.add(model.Change[j] >= -0.5 * sum(model.Production[i, j + 1] - model.Production[i, j] for i in model.I) - model.sita[j] * 9999999)    model.ChangeCost.add(model.Change[j] <= -0.5 * sum(model.Production[i, j + 1] - model.Production[i, j] for i in model.I) + model.sita[j] * 9999999)#求解opt = SolverFactory('cbc')result = opt.solve(model)
print('++++++++++++++++++++++++++++++++++++++')model.Production.pprint()print('++++++++++++++++++++++++++++++++++++++')model.Change.pprint()print('++++++++++++++++++++++++++++++++++++++')model.Storage.pprint()print('++++++++++++++++++++++++++++++++++++++')model.sita.pprint()print('Finally Cost: ', model.cost())result.write()

 ✦ 

结果展示...

Production : Size=24, Index=Production_index    Key    : Lower : Value  : Upper : Fixed : Stale : Domain    (0, 0) :     0 : 1490.0 : 10000 : False : False : NonNegativeIntegers    (0, 1) :     0 : 3000.0 : 10000 : False : False : NonNegativeIntegers    (0, 2) :     0 : 2000.0 : 10000 : False : False : NonNegativeIntegers    (0, 3) :     0 : 4000.0 : 10000 : False : False : NonNegativeIntegers    (0, 4) :     0 : 2000.0 : 10000 : False : False : NonNegativeIntegers    (0, 5) :     0 : 2550.0 : 10000 : False : False : NonNegativeIntegers    (1, 0) :     0 : 1300.0 : 10000 : False : False : NonNegativeIntegers    (1, 1) :     0 :  800.0 : 10000 : False : False : NonNegativeIntegers    (1, 2) :     0 :  800.0 : 10000 : False : False : NonNegativeIntegers    (1, 3) :     0 : 1000.0 : 10000 : False : False : NonNegativeIntegers    (1, 4) :     0 : 1100.0 : 10000 : False : False : NonNegativeIntegers    (1, 5) :     0 :  910.0 : 10000 : False : False : NonNegativeIntegers    (2, 0) :     0 : 2882.0 : 10000 : False : False : NonNegativeIntegers    (2, 1) :     0 : 1672.0 : 10000 : False : False : NonNegativeIntegers    (2, 2) :     0 : 2773.0 : 10000 : False : False : NonNegativeIntegers    (2, 3) :     0 : 1073.0 : 10000 : False : False : NonNegativeIntegers    (2, 4) :     0 : 1891.0 : 10000 : False : False : NonNegativeIntegers    (2, 5) :     0 : 1419.0 : 10000 : False : False : NonNegativeIntegers    (3, 0) :     0 : 1400.0 : 10000 : False : False : NonNegativeIntegers    (3, 1) :     0 : 1600.0 : 10000 : False : False : NonNegativeIntegers    (3, 2) :     0 : 1500.0 : 10000 : False : False : NonNegativeIntegers    (3, 3) :     0 : 1000.0 : 10000 : False : False : NonNegativeIntegers    (3, 4) :     0 : 1100.0 : 10000 : False : False : NonNegativeIntegers    (3, 5) :     0 : 1210.0 : 10000 : False : False : NonNegativeIntegers# ==========================================================Finally Cost:  684210.2# ==========================================================Time: 1.277825117111206

参考资料:

未完待续

相关推荐