1.前言
在之前的文章中,运筹优化工具ortools解读与实际-概览 曾经介绍过线性规划、整数布局、束缚布局的概念。
针对以上问题,ortools提供了多种解决方案:
- MPSolver接口:ortools提供一个名为MPSolver的接口,实现了对泛滥LP、MIP第三方求解器的封装,其中也蕴含了google自研的线性规划求解器Glop。
- CP-SAT求解器:google自研解决CP问题的求解器。为了进步求解器计算速度,CP-SAT被设计成只反对整数问题。因而,当束缚中含有非整数变量时,须要通过乘以一个较大的数转变为整数问题。
- Original CP Solver:解决CP问题的晚期求解器,官网倡议应用CP-SAT。
有意思的是,同样是自研求解器,Glop与第三方求解器并无差别,都须要MPSolver初始化之后应用;CP-SAT却是非凡的存在,不与MPSolver一起应用,编码格调也有较大变动。
这里以一个线性规划问题引出MPSolver及Glop求解器的应用办法(整数布局和混合整数布局同样实用)。
2.套路
针对线性规划问题,求解过程分为7步:
- 1.加载一个线性求解器的wrapper
- 2.实例化求解器
- 3.创立变量,为变量设置上上限
- 4.定义束缚
- 5.定义指标函数
- 6.调用求解办法
- 7.输入求解后果
在ortools官网文档中,主推google自研的GLOP求解器,从应用体验来讲,GLOP开发接口十分敌对,代码开发与数学建模完全一致,执行效率也十分高。
3.案例
依照运筹优化工具ortools解读与实际-概览中的“三段式”形容解决线性规划问题的过程。
定义问题
e.g.某企业有两种不同的投资渠道,A渠道有3倍潜在收益,B渠道有4倍的潜在收益。但思考到企业倒退和危险管制,制订了肯定的投资标准:
1)企业账户总金额为14W
2)每向B渠道减少一笔投资,须要留出雷同数额的准备金应答危险
3)投资B渠道的金额不应超过A渠道投资额的3倍
4)A渠道投资额与B渠道投资额之差,应该2W以内
定义问题须要将后续建模中优化指标、决策变量、约束条件以自然语言的模式明确分明。其次,作为算法人员,须要依据集体常识储备及过往教训大体明确问题类型和可能的解决方案。
在工作当中,定义问题是最重要的环节。
建模问题
优化指标:Maximize 3x + 4y
决策变量: x、y
约束条件:
- x + 2y ≤ 14
- 3x – y ≥ 0
- x – y ≤ 2
- x ≥ 0
- y ≥ 0
求解问题
#1.加载一个线性求解器的wrapperfrom ortools.linear_solver import pywraplp#2.实例化求解器,这里应用的是Google自研的GLOP线性规划求解器solver = pywraplp.Solver.CreateSolver('GLOP')# 3.创立变量,为变量设置上上限x = solver.NumVar(0, solver.infinity(), 'x')y = solver.NumVar(0, solver.infinity(), 'y')print('Number of variables =', solver.NumVariables())#4.定义束缚# Constraint 0: x + 2y <= 14.solver.Add(x + 2 * y <= 14.0)# Constraint 1: 3x - y >= 0.solver.Add(3 * x - y >= 0.0)# Constraint 2: x - y <= 2.solver.Add(x - y <= 2.0)print('Number of constraints =', solver.NumConstraints())#5.定义指标函数solver.Maximize(3 * x + 4 * y)#6.调用求解办法status = solver.Solve()#7.输入求解后果if status == pywraplp.Solver.OPTIMAL: print('Solution:') print('Objective value =', solver.Objective().Value()) print('x =', x.solution_value()) print('y =', y.solution_value())else: print('The problem does not have an optimal solution.')print('\nAdvanced usage:')print('Problem solved in %f milliseconds' % solver.wall_time())print('Problem solved in %d iterations' % solver.iterations())
4.MPSolver构造
除了GLOP,ortools同样反对其余求解器。
- CLP_LINEAR_PROGRAMMING or CLP
- CBC_MIXED_INTEGER_PROGRAMMING or CBC
- GLOP_LINEAR_PROGRAMMING or GLOP
- BOP_INTEGER_PROGRAMMING or BOP
- SAT_INTEGER_PROGRAMMING or SAT or CP_SAT
- SCIP_MIXED_INTEGER_PROGRAMMING or SCIP
- GUROBI_LINEAR_PROGRAMMING or GUROBI_LP
- GUROBI_MIXED_INTEGER_PROGRAMMING or GUROBI or GUROBI_MIP
- CPLEX_LINEAR_PROGRAMMING or CPLEX_LP
- CPLEX_MIXED_INTEGER_PROGRAMMING or CPLEX or CPLEX_MIP
- XPRESS_LINEAR_PROGRAMMING or XPRESS_LP
- XPRESS_MIXED_INTEGER_PROGRAMMING or XPRESS or XPRESS_MIP
- GLPK_LINEAR_PROGRAMMING or GLPK_LP
- GLPK_MIXED_INTEGER_PROGRAMMING or GLPK or GLPK_MIP
例如,将GLOP变更为SCIP或CP-SAT,同样能够失去最优解。
# solver = pywraplp.Solver.CreateSolver('GLOP')solver = pywraplp.Solver.CreateSolver('SCIP')
有几个重要的函数,别离是:
NumVar():创立数值型变量,当然也能够创立其余变量,例如整型变量IntVar、布尔型变量BoolVar
solver.Add():设置束缚
solver.Objective().Value():获取指标函数后果
solution_value:获取求解后果中变量的值
SetMaximization:指定求最大值还是求最小值
理解了以上次要内容,就能够照葫芦画瓢解决一些简略的线性规划问题了!遇到其余细节能够查阅ortools源代码,理解更多办法。
5.硬编码 VS 数组存储
能够看到,咱们将模型因素硬编码到代码中,这种形式可读性十分强。同样也有毛病,如果咱们的因素参数并不是固定的,会随着工夫的推移发生变化,此时,硬编码的模式将不利于工程化实现。
为此,ortools提供了一种数组存储模型因素的形式。
咱们将第3节案例,依照数组存储模型的形式进行编码,运行后果雷同。
#1.加载一个线性求解器的wrapperfrom ortools.linear_solver import pywraplpdef create_data_model(): """Stores the data for the problem.""" data = {} #束缚不等式左侧局部参数 data['constraint_coeffs'] = [ [1, 2], [-3, 1], [1, -1], ] #束缚不等式右侧局部参数 data['bounds'] = [14.0, 0.0, 2.0] #指标函数局部参数 data['obj_coeffs'] = [3,4] data['num_vars'] = 2 data['num_constraints'] = 3 return datadata = create_data_model()#2.实例化求解器,这里应用的是Google自研的GLOP线性规划求解器solver = pywraplp.Solver.CreateSolver('GLOP')infinity = solver.infinity()# 3.创立变量,为变量设置上上限x = {}for j in range(data['num_vars']): x[j] = solver.NumVar(0, infinity, 'x[%i]' % j)print('Number of variables =', solver.NumVariables())#4.定义束缚for i in range(data['num_constraints']): constraint = solver.RowConstraint(-infinity, data['bounds'][i], '') for j in range(data['num_vars']): constraint.SetCoefficient(x[j], data['constraint_coeffs'][i][j])print('Number of constraints =', solver.NumConstraints())#5.定义指标函数objective = solver.Objective()for j in range(data['num_vars']): objective.SetCoefficient(x[j], data['obj_coeffs'][j])objective.SetMaximization()#6.调用求解办法status = solver.Solve()#7.输入求解后果if status == pywraplp.Solver.OPTIMAL: print('Solution:') print('Objective value =', solver.Objective().Value()) for j in range(data['num_vars']): print(x[j].name(), ' = ', x[j].solution_value())else: print('The problem does not have an optimal solution.')print('\nAdvanced usage:')print('Problem solved in %f milliseconds' % solver.wall_time())print('Problem solved in %d iterations' % solver.iterations())
特地留神,为了束缚局部代码块编码不便,须要将束缚对立,我这里对立改成了≤的模式。对应的,create_data_model()办法中的编码与建模时存在差别。
原束缚表达形式:
- x + 2y ≤ 14
- -3x + y ≤ 0
- x – y ≤ 2
- x ≥ 0
- y ≥ 0
援用
调整后的表达形式:- x + 2y ≤ 14
- 3x – y ≥ 0
- x – y ≤ 2
- x ≥ 0
- y ≥ 0
6.结语
本篇文章次要解说了ortools利用MPSolver及其他求解器解决LP、IP、MIP问题的两种办法。这里引出两个问题:
- 不同求解器有哪些异同,理论我的项目中应该如何抉择?
- CP问题如何求解,google自研的CP-SAT有何特别之处?
以上两个问题将在后续文章中探讨。