7.3. Pyomo 的建模与优化

Pyomo 是基于 Python 开发的第三方开源建模语言,它支持对线性规划、混合整数规划、非线性规划等问题建模和分析,并调用其他商用或开源的求解器进行求解。

目前,MindOpt 可对在 Windows/Linux/OSX 平台上通过 Pyomo 建立的 线性规划模型 进行求解。关于 Pyomo 的更多详细内容,请参考 Pyomo 官方文档库

在本节中,我们将介绍如何使用 Pyomo API 来建立 线性规划问题示例 中的优化问题,并调用 MindOpt 求解。

7.3.1. 安装 Pyomo

用户需首先安装 MindOpt。关于 MindOpt 的安装与配置请参考 单机版安装MindOpt 安装后,用户可以通过以下 pip 的命令方式来安装 Pyomo:

pip install pyomo

关于 Pyomo 的详细安装方式,请参考 Pyomo Installation

7.3.2. Pyomo 调用接口

MindOpt 为Pyomo提供两种不同的插件接口,分别是mindopt_direct和mindopt_persistent。

mindopt_direct是一种简单的接口,它在每次调用求解器时都会重新建立一个Mindopt模型,并通过将pyomo模型中的约束和变量转化为Mindopt模型的约束和变量来进行求解。这种接口适用于简单的优化问题,对于复杂的问题可能会导致一些求解器性能上的损失。

mindopt_persistent是一种更高级的接口,它将模型对象存于内存中,允许多次修改和重复求解。对于需要多次求解相同模型的问题非常有用,可以避免重复创建和销毁模型对象产生的开销,从而提高求解效率。

7.3.2.1. 调用 Pyomo Direct 接口

MindOpt 在 Pyomo Direct接口文件( mindopt_pyomo.py )内定义了 Pyomo 调用 MindOpt 所需的相关接口。此接口继承自 Pyomo 的 DirectSolver 类别,实现细节见安装包内的接口文件:

<MDOHOME>/<VERSION>/<PLATFORM>/lib/pyomo/mindopt_persistent.py
<MDOHOME>/<VERSION>/<PLATFORM>/lib/pyomo/mindopt_pyomo.py

用户在使用时需首先将该接口文件移到当前工作目录中,并在 Python 代码中加载该模块中定义的 MindoDirect 类:

30from mindopt_pyomo import MindoDirect

接着,我们调用 Pyomo API 来建立 线性规划问题示例 中的优化问题。关于 Pyomo API 的详细说明,请参考 Pyomo 官方文档API

41    # Define variables and constraints.
42    variable_names = ['x0', 'x1', 'x2', 'x3']
43    var_lb = {'x0':0, 'x1':0, 'x2':0, 'x3':0}
44    var_ub = {'x0':10, 'x1':None, 'x2':None, 'x3':None}
45    objective_coeff = {'x0': 1, 'x1': 1, 'x2': 1, 'x3': 1}
46    constraint_names = ['c0', 'c1']
47    constraint_bound = {'c0': 1, 'c1': 1}
48    constraint_coeff = {('x0', 'c0'): 1, ('x1', 'c0'): 1, ('x2', 'c0'): 2, ('x3', 'c0'): 3,
49                        ('x0', 'c1'): 1, ('x1', 'c1'): -1, ('x2', 'c1'): 0, ('x3', 'c1'): 6}
50
51    # Create model.
52    model = ConcreteModel(name="ex1")
53
54    # Build decision variables.
55    model.Set = Set(initialize=variable_names)
56    model.Variable = Var(model.Set, within=NonNegativeReals, bounds=fb)
57
58    # Objective.
59    model.obj = Objective(expr=sum(objective_coeff[var_name] * model.Variable[var_name] for var_name in variable_names), sense=minimize)
60
61    # Constraints.
62    model.dual = Suffix(direction=Suffix.IMPORT)
63    model.cons1 = Constraint(expr = ConstraintsRule(model, 'c0') >= constraint_bound['c0'])
64    model.cons2 = Constraint(expr = ConstraintsRule(model, 'c1') == constraint_bound['c1'])

求解前,我们指定使用 MindOpt 求解器,并对求解的相关参数进行设置(求解器参数请查阅 参数):

69    # Solve problem by MindOpt solver.
70    opt = SolverFactory("mindo_direct")
71    
72    # Set options.
73    opt.options['Method'] = -1
74    opt.options['IPM/PrimalTolerance'] = 1e-10

最后,调用 Pyomo 的求解函数 solve() 进行求解,并获取相关的结果:

79    # Solve.
80    results = opt.solve(model)
81
82    # Summary of result.
83    results.write()
84
85    if (results.solver.status == SolverStatus.ok) and (results.solver.termination_condition == TerminationCondition.optimal):
86        print("The solution is optimal and feasible.")
87        model.Variable.display()
88        model.dual.display()
89        model.obj.display()
90    elif (results.solver.termination_condition == TerminationCondition.infeasible):
91        print("The model is infeasible.")
92        print("Solver Status: ",  results.solver.status)
93    else:
94        print("Something else is wrong.")
95        print("Solver Status: ",  results.solver.status)

7.3.2.2. 调用 Pyomo Persistent 接口文件

MindOpt 在 Pyomo Persistent接口文件( mindopt_persistent.py )内定义了 Pyomo 调用 MindOpt 所需的相关persistent接口。此接口继承自 Pyomo 的 PersistentSolver 类别,实现细节见安装包内的接口文件:

<MDOHOME>/<VERSION>/<PLATFORM>/sh/pyomo/mindopt_persistent.py

用户在使用时需首先将该接口文件移到当前工作目录中,并在 Python 代码中加载该模块中定义的 MindoptPersistent 类:

from mindopt_persistent import MindoptPersistent

接着,我们调用 Pyomo API 来建立 线性规划问题示例 中的优化问题。关于 Pyomo API 的详细说明,请参考 Pyomo 官方文档API

model = ConcreteModel()

products = ['A', 'B']
materials = ['X', 'Y']

model.x = Var(products, domain=NonNegativeReals)

model.profit = Objective(expr=10 * model.x['A'] + 15 * model.x['B'], sense=maximize)

model.material_x = Constraint(expr=2 * model.x['A'] + model.x['B'] <= 10)
model.material_y = Constraint(expr=model.x['A'] + 2 * model.x['B'] <= 6)

solver = SolverFactory('mindopt_persistent')

在求解前必须使用set_instance将模型保存在内存中。

solver.set_instance(model)
solver.solve(model)

求解后也可以使用Persistent提供的函数获取模型属性/参数,也可以修改原有模型。

var_attr = solver.get_var_attr(model.x['A'], 'lb')
print(var_attr)

7.3.3. 使用 SOS 功能

SOS(Special Ordered Sets)约束是一种特殊的约束类型,用于限制一组变量中最多只能有N个变量取非零值。SOS约束在某些优化问题中非常有用,比如在调度问题中,可以用SOS约束来确保每个时间点只能选择一个任务执行。

无论是使用mindopt_direct还是mindopt_persistent接口,都可以使用 MindOpt 的SOS功能。

用户在通过Pyomo建模代码使用SOS时,需要首先定义一个SOS集合,将特定变量分组到有序集合中。例如,假设有两个变量x和y属于同一个SOS组,则可实现如下代码:

from pyomo.environ import *

model = ConcreteModel()
model.x = Var(within=NonNegativeReals)
model.y = Var(within=NonNegativeReals)

model.sos = SOSConstraint(
   rule=lambda m: [(m.x, m.y)]
)

在上述示例中, SOSConstraint 的参数rule用于指定变量属于SOS集合的规则。在此示例中,我们将变量x和y分组到同一个SOS集合中。

7.3.4. 使用 Callback 功能

无论是使用mindopt_direct还是mindopt_persistent接口,用户都可以使用 MindOpt 的回调功能,提供一些启发式决策用以优化求解速度。

但是考虑到需要使用回调功能的模型一般比较复杂,推荐用户使用mindopt_persistent接口。

在Pyomo中使用Mindopt的回调函数功能,首先需要定义一个回调函数来处理特定的事件。

from pyomo.environ import *

def my_callback(model, where):
   where_names = {
         3: "MIPSTART",
         4: "MIPSOL",
         5: "MIPNODE"
   }

   print("where = {}".format(where_names[where]))

   # Counter for calls
   model._calls += 1

   # Print the constraint matrix of the presolved model. The getA() function returns a scipy.sparse.csc_matrix.
   # Note that the model passed to the callback function is a presolved model in Python SDK.
   # Any modifications to the constraint matrix should not be allowed as it may lead to unforeseen issues.
   print(model.getA())

   bounds = []
   for v in model.getVars():
         bounds.append((model.cbGet(MDO.Callback.PRE_LB, v), model.cbGet(MDO.Callback.PRE_LB, v)))

   # Print the upper and lower bounds for the original problem
   print(bounds)

   bounds.clear()
   for v in model.getVars():
         bounds.append((model.cbGet(MDO.Callback.MIP_LB, v), model.cbGet(MDO.Callback.MIP_LB, v)))

在上述示例中,我们定义了一个名为my_callback的回调函数。在此回调函数中,我们使用cbGet方法获取当前预处理模型的一些信息,比如变量的lower bound等。 然后,您需要将回调函数与Mindopt求解器连接起来,并为特定的事件注册回调函数。

opt = SolverFactory('mindopt_persistent')
opt.set_instance(model)

# Register the callback function
opt.set_callback(my_callback)

results = opt.solve()

在上述示例中,我们使用SolverFactory创建了一个mindopt_persistent求解器实例,并将回调函数my_callback注册到该求解器中。然后,我们使用solve方法解决模型,并触发相应的回调事件。

7.3.5. 建模示例: mdo_pyomo_lo_ex1

文件链接 mdo_pyomo_lo_ex1.py 中提供了完整代码:

 1"""
 2/**
 3 *  Description
 4 *  -----------
 5 *
 6 *  Linear optimization.
 7 *
 8 *  Formulation
 9 *  -----------
10 *
11 *  Minimize
12 *    obj: 1 x0 + 1 x1 + 1 x2 + 1 x3
13 *  Subject To
14 *   c1 : 1 x0 + 1 x1 + 2 x2 + 3 x3 >= 1
15 *   c2 : 1 x0 - 1 x2        + 6 x3 == 1
16 *  Bounds
17 *    0 <= x0 <= 10
18 *    0 <= x1
19 *    0 <= x2
20 *    0 <= x3
21 *  End
22 */
23"""
24
25from pyomo.environ import *
26from pyomo.core.base.util import Initializer
27from pyomo.opt import SolverFactory
28from pyomo.opt import SolverStatus, TerminationCondition
29
30from mindopt_pyomo import MindoDirect
31from mindoptpy import *
32
33def ConstraintsRule(model, p):
34    return sum(constraint_coeff[i, p] * model.Variable[i] for i in variable_names)
35
36def fb(model, i):
37    return (var_lb[i], var_ub[i])
38
39if __name__ == '__main__':
40
41    # Define variables and constraints.
42    variable_names = ['x0', 'x1', 'x2', 'x3']
43    var_lb = {'x0':0, 'x1':0, 'x2':0, 'x3':0}
44    var_ub = {'x0':10, 'x1':None, 'x2':None, 'x3':None}
45    objective_coeff = {'x0': 1, 'x1': 1, 'x2': 1, 'x3': 1}
46    constraint_names = ['c0', 'c1']
47    constraint_bound = {'c0': 1, 'c1': 1}
48    constraint_coeff = {('x0', 'c0'): 1, ('x1', 'c0'): 1, ('x2', 'c0'): 2, ('x3', 'c0'): 3,
49                        ('x0', 'c1'): 1, ('x1', 'c1'): -1, ('x2', 'c1'): 0, ('x3', 'c1'): 6}
50
51    # Create model.
52    model = ConcreteModel(name="ex1")
53
54    # Build decision variables.
55    model.Set = Set(initialize=variable_names)
56    model.Variable = Var(model.Set, within=NonNegativeReals, bounds=fb)
57
58    # Objective.
59    model.obj = Objective(expr=sum(objective_coeff[var_name] * model.Variable[var_name] for var_name in variable_names), sense=minimize)
60
61    # Constraints.
62    model.dual = Suffix(direction=Suffix.IMPORT)
63    model.cons1 = Constraint(expr = ConstraintsRule(model, 'c0') >= constraint_bound['c0'])
64    model.cons2 = Constraint(expr = ConstraintsRule(model, 'c1') == constraint_bound['c1'])
65
66    # Print formulation of model.
67    model.pprint()
68
69    # Solve problem by MindOpt solver.
70    opt = SolverFactory("mindo_direct")
71    
72    # Set options.
73    opt.options['Method'] = -1
74    opt.options['IPM/PrimalTolerance'] = 1e-10
75    
76    # Solve.
77    results = opt.solve(model)
78
79    # Summary of result.
80    results.write()
81
82    if (results.solver.status == SolverStatus.ok) and (results.solver.termination_condition == TerminationCondition.optimal):
83        print("The solution is optimal and feasible.")
84        model.Variable.display()
85        model.dual.display()
86        model.obj.display()
87    elif (results.solver.termination_condition == TerminationCondition.infeasible):
88        print("The model is infeasible.")
89        print("Solver Status: ",  results.solver.status)
90    else:
91        print("Something else is wrong.")
92        print("Solver Status: ",  results.solver.status)

其他 Pyomo 的示例请参考 Pyomo Gallery