pyomo的全称是Python Optimization Modeling Objects(python优化建模对象)
模型类型
抽象模型
例如,下面的方程代表一个线性程序(LP),以找到(参数n和b,以及参数向量a和c的)向量x的最优值。
我们称其为抽象的或符号的数学模型,因为它依赖于未指定的参数值。AbstractModel
类为Pyomo中定义和初始化抽象优化模型提供了一个接口。
具体模型
下面的LP模型是前面抽象模型的一个具体实例。
ConcreteModel
类用于定义Pyomo中的具体优化模型。
Python程序员可能更喜欢写具体的模型,而其他一些代数建模语言的用户可能更倾向于写抽象的模型。选择在很大程度上是一个品味的问题;一些应用可能使用其中一种更直接一些。
个人感觉,抽象模型更像是具体模型的类。
模型代码样例
具体的模型:
import pyomo.environ as pyo
#建立具体模型
model = pyo.ConcreteModel()
# 模型第三行
model.x = pyo.Var([1,2], domain=pyo.NonNegativeReals)
# 模型第一行
model.OBJ = pyo.Objective(expr = 2*model.x[1] + 3*model.x[2])
# 代表模型第二行
model.Constraint1 = pyo.Constraint(expr = 3*model.x[1] + 4*model.x[2] >= 1)
抽象的模型:
在Pyomo中实现这一点的一种方法如下所示:
from __future__ import division
import pyomo.environ as pyo
model = pyo.AbstractModel()
model.m = pyo.Param(within=pyo.NonNegativeIntegers)
model.n = pyo.Param(within=pyo.NonNegativeIntegers)
model.I = pyo.RangeSet(1, model.m)
model.J = pyo.RangeSet(1, model.n)
model.a = pyo.Param(model.I, model.J)
model.b = pyo.Param(model.I)
model.c = pyo.Param(model.J)
# the next line declares a variable indexed by the set J
model.x = pyo.Var(model.J, domain=pyo.NonNegativeReals)
def obj_expression(m):
return pyo.summation(m.c, m.x)
model.OBJ = pyo.Objective(rule=obj_expression)
def ax_constraint_rule(m, i):
# return the expression for the constraint for i
return sum(m.a[i,j] * m.x[j] for j in m.J) >= m.b[i]
# the next line creates one constraint for each member of the set model.I
model.AxbConstraint = pyo.Constraint(model.I, rule=ax_constraint_rule)
这里的代码比较复杂,下面对其进行解释:
from __future__ import division
导入python未来支持的语言特征division(精确除法),当我们没有在程序中导入该特征时,"/“操作符执行的是截断除法(Truncating Division),当我们导入精确除法之后,”/"执行的是精确除法,如果有兴趣的可以参考这篇文章:python from future import division
model = pyo.AbstractModel()
模型的声明也是必需的。使用模型这个名字不是必须的。几乎任何名字都可以使用,但是我们将在大多数的例子中使用模型这个名字。在这个例子中,我们声明它将是一个抽象的模型。
model.m = pyo.Param(within=pyo.NonNegativeIntegers)
model.n = pyo.Param(within=pyo.NonNegativeIntegers)
m与n是模型中的Param
,这里within
是给参数设置范围为非负。如果不设置,那么Pyomo就不会反对任何类型的数据被分配给这些参数。
model.I = pyo.RangeSet(1, model.m)
model.J = pyo.RangeSet(1, model.n)
尽管不是必须的,但定义索引集是很方便的。在这个例子中,我们使用RangeSet
组件来声明这些集合将是一个从1开始到由model.m
和model.n
参数指定的值结束的整数序列。这其实就是相当于range(1,m)。
model.a = pyo.Param(model.I, model.J)
model.b = pyo.Param(model.I)
model.c = pyo.Param(model.J)
系数和右边的数据被定义为索引参数。当集合作为参数给Param组件的参数时,表示该集合将对参数进行索引。其实就是根据式子里的标号来,比如式子里是bi
那b的Param
里面就写I
。这其实相当于把a变成了一个矩阵,把b和c变成了向量
# the next line declares a variable indexed by the set J
model.x = pyo.Var(model.J, domain=pyo.NonNegativeReals)
Var 组件的第一个参数是一个集合,所以它被定义为该变量的索引集。在这种情况下,变量只有一个索引集,但是可以使用多个索引集,就像参数model.a
的声明那样。第二个参数为变量指定一个域。这个信息是模型的一部分,在提供数据和求解模型时将传递给求解器。NonNegativeReals
域的指定实现了变量大于或等于零的要求。
def obj_expression(m):
return pyo.summation(m.c, m.x)
当给定两个参数时,summation()函数返回两个参数在其索引上的乘积之和的表达式。当然,这只在两个参数有相同的索引时有效。如果只给了一个参数,则返回该参数所有索引上的和的表达式。
在本式中:summation(m.c,m.x)
相当于表达式
∑
j
=
1
n
=
c
j
x
j
\sum_{j=1}^n=c_jx_j
∑j=1n=cjxj
model.OBJ = pyo.Objective(rule=obj_expression)
为了声明一个目标函数,使用了Objective
的Pyomo
类。规则参数给出了一个返回目标表达式的函数的名称。默认的意义是最小化。对于最大化,必须使用sense=pyo.maximum
参数。
def ax_constraint_rule(m, i):
# return the expression for the constraint for i
return sum(m.a[i,j] * m.x[j] for j in m.J) >= m.b[i]
约束的声明是类似的。一个函数被声明用来生成约束表达式。在这种情况下可以有多个相同形式的约束,因为我们在表达式中用i来索引约束 ∑ j = 1 n = a i j x j ≥ b i ∀ i = 1 ⋅ ⋅ ⋅ m \sum_{j=1}^n=a_ijx_j≥b_i \forall i=1···m ∑j=1n=aijxj≥bi ∀i=1⋅⋅⋅m
# the next line creates one constraint for each member of the set model.I
model.AxbConstraint = pyo.Constraint(model.I, rule=ax_constraint_rule)
为了声明使用这种表达方式的约束,我们使用PyomoConstraint
组件,该组件接受各种参数。在这种情况下,我们的模型规定我们可以有一个以上的相同形式的约束,我们已经创建了一个集合,model.I
,这些约束可以被索引,所以这就是约束声明的第一个参数。下一个参数给出了将用于生成约束条件表达式的规则。从整体上看,这个约束声明说,将创建一个由集合model.I
索引的约束列表,对于model.I
的每个成员,将调用函数ax_constraint_rule
,它将被传递给model
对象以及model.I
的成员
为了使用这个模型,必须提供参数值的数据。下面是一个提供数据的文件 :
(in AMPL “.dat” format).
# one way to input the data in AMPL format
# for indexed parameters, the indexes are given before the value
param m := 1 ;
param n := 2 ;
param a :=
1 1 3
1 2 4
;
param c:=
1 2
2 3
;
param b := 1 1 ;
只有一个约束,所以model.a
只需要两个值。当用AMPL格式给数组和向量赋值时,一种方法是给出索引和值。这一行1 2 4
使model.a[1,2]
得到4
的值。由于model.c
只有一个索引,所以只需要一个索引值,例如,1 2
使model.c[1]
得到值2
。在AMPL格式的数据文件中,换行通常并不重要,所以model.b的单一索引值的赋值是在一行中给出的,因为这很容易阅读。
有多种格式可以用来向Pyomo模型提供数据,但AMPL格式对我们的目的来说很好,因为它包含了数据元素的名称和数据。在AMPL数据文件中,#之后的文本被当作注释。行数一般不重要,但语句必须用分号来结束。
处理符号索引集
当使用Pyomo(或任何其他AML)工作时,使用包含字符串的索引集,而不是由1…m或从1到n的总和所暗示的索引集,以某种更抽象的方式编写抽象模型是很方便的。当这样做时,集合的大小是由输入暗示的,而不是直接指定。此外,索引条目可能没有真正的顺序。通常,在同一个模型中需要混合使用整数和索引以及字符串作为索引。为了开始说明一般的索引,请考虑我们刚才介绍的模型的一个稍微不同的Pyomo实现。
# abstract2.py
from __future__ import division
from pyomo.environ import *
model = AbstractModel()
model.I = Set()
model.J = Set()
model.a = Param(model.I, model.J)
model.b = Param(model.I)
model.c = Param(model.J)
# the next line declares a variable indexed by the set J
model.x = Var(model.J, domain=NonNegativeReals)
def obj_expression(model):
return summation(model.c, model.x)
model.OBJ = Objective(rule=obj_expression)
def ax_constraint_rule(model, i):
# return the expression for the constraint for i
return sum(model.a[i,j] * model.x[j] for j in model.J) >= model.b[i]
# the next line creates one constraint for each member of the set model.I
model.AxbConstraint = Constraint(model.I, rule=ax_constraint_rule)
其实就是不定义m,n,在文件里直接设置I和J
而且,这个模型也可以使用有意义的索引为相同的一般形式的问题提供不同的数据。我们可以对比一下,左边是原来的AMPL,右边是新的AMPL
param m := 1 ; set I := TV Film ;
param n := 2 ; set J := Graham John Carol ;
param a := param a :=
1 1 3 TV Graham 3
1 2 4 TV John 4.4
TV Carol 4.9
Film Graham 1
Film John 2.4
Film Carol 1.1
; ;
param c:= param c := [*]
1 2 Graham 2.2
2 3 John 3.1416
Carol 3
; ;
param b := 1 1 ; param b := TV 1 Film 1 ;
便于复制粘贴:
# abstract2.dat AMPL data format
set I := TV Film ;
set J := Graham John Carol ;
param a :=
TV Graham 3
TV John 4.4
TV Carol 4.9
Film Graham 1
Film John 2.4
Film Carol 1.1
;
param c := [*]
Graham 2.2
John 3.1416
Carol 3
;
param b := TV 1 Film 1 ;
模型求解
具体模型的求解
如果你有一个ConcreteModel
,在你的Python脚本的底部添加这几行来求解它
opt = pyo.SolverFactory('glpk')
opt.solve(model)
抽象模型的求解
如果你有一个AbstractModel
,你必须在求解你的模型之前用上述同样的行来创建一个具体实例。
instance = model.create_instance()
opt = pyo.SolverFactory('glpk')
opt.solve(instance)
pyomo solve
命令行求解
Pyomo支持建模和脚本,但不会自动安装求解器。为了解决一个模型,必须在要使用的计算机上安装一个求解器。如果有一个求解器,那么就可以使用pyomo
命令来解决一个问题实例。
假设计算机上安装了名为glpk
(也叫glpsol
)的求解器。假设抽象模型在名为 abstract1.py
的文件中,其数据文件在名为 abstract1.dat
的文件中。在命令提示符下,如果这两个文件都在当前目录下,可以通过以下命令获得一个解决方案。
pyomo solve abstract1.py abstract1.dat --solver=glpk
由于glpk是默认的求解器,确实没有必要指定它,所以可以放弃--solver
如果安装了CPLEX
:
pyomo solve abstract1.py abstract1.dat --solver=cplex
显示结果:
[ 0.00] Setting up Pyomo environment
[ 0.00] Applying Pyomo preprocessing actions
[ 0.07] Creating model
[ 0.15] Applying solver
[ 0.37] Processing results
Number of solutions: 1
Solution Information
Gap: 0.0
Status: optimal
Function Value: 0.666666666667
Solver results file: results.json
[ 0.39] Applying Pyomo postprocessing actions
[ 0.39] Pyomo Finished
我的是--solver=ipopt
,如果是在ros路径里面,打开放python的文件,命令行输入pwd
打印当前路径,复制到文件前面就行。
方括号中的数字表示每个步骤需要多少时间。结果被写入名为results.json
的文件中,该文件具有特殊的结构,对后期处理很有帮助。要看到写在屏幕上的结果摘要,请使用--summary
。
pyomo solve abstract1.py abstract1.dat --solver=cplex --summary
查看pyomo命令行:
pyomo solve --help