Docplex入门（1）——线性规划

BoJack_Horse

已于 2022-06-06 09:39:29 修改

阅读量1.9k

点赞数 1

分类专栏： (do)cplex 文章标签：线性代数

于 2022-06-05 23:28:44 首次发布

本文链接：https://blog.csdn.net/BoJack_Horse/article/details/125079788

版权

(do)cplex 专栏收录该内容

3 篇文章 0 订阅

订阅专栏

Docplex入门——线性规划

定义及前置知识

线性规划：解决最大/小化一个线性目标函数，所有决策变量都是连续，且要求目标和约束都由线性表达式组成
线性表达式（linear_expression）：
$\sum{a_ix_i} 亦可写为矩阵形式AX$
其中要求A为常数向量（已知），X为自变量（待求解且连续），对比一下非线性的例子：

1.两个/多个变量的乘法
2.指数
3.对数
4.绝对值
..

线性约束：
$linear_ -expression \quad symbol \quad linear_-expression$
其中symbol仅能取 $=$ 、 $\le$ 、 $\ge$
一般形式：
$\min C^t _x \\s.t. \quad Ax \ge B \\x \ge 0$
其中x为大小为n的向量，A为m*n的系数矩阵，B是大小为m的常向量。

线性规划样例

问题描述：
要生产台式和移动两类电话，要求每种类型电话至少生产100台，且最大化利润。但是公司生产能力有限，制造过程分为装配和涂装2个步骤，需要在不超过生产能力情况下计算每款手机最优生产数量。
建模：
决策变量： $A 、 B 各自的数量$
目标： $maximize: 12desk_-production + 20cell_-prodction$
约束条件：
$\ desk_-production \ge 100 \\ cell_-production \ge 100 \\ 0.2desk_-production + 0.4cell_-production \le 400 \\ 0.5desk_-production + 0.4cell_-production \le 490$
python实现：

# 导入模型包并创建实例：
from docplex.mp.model import Model
model = Model(name='telephone_production')

# 创建待求解变量desk_production和cell_production
desk_production = model.continuous_var(name='desk_production')
cell_production = model.continuous_var(name='cell_production')

# 目标函数
model.maximize(12*desk_production+20*cell_production)

# 4个约束条件，建议写约束条件都按照后两个的方式写等式
model.add_constraint(desk_production >= 100)
model.add_constraint(cell_production >= 100)
model.add_constraint(0.2*desk_production + 0.4*cell_production <= 400, ctname='ct_assembly')
model.add_constraint(0.5*desk_production + 0.4*cell_production <= 490, ctname='ct_painting')


# 打印模型信息
model.print_information()

# 对模型进行求解
sol = model.solve()

# 打印求解结果
print(sol)

输出结果：

Model: telephone_production
 - number of variables: 2
   - binary=0, integer=0, continuous=2
 - number of constraints: 4
   - linear=4
 - parameters: defaults
 - objective: maximize
 - problem type is: LP

solution for: telephone_production	#这个显然是上面模型初始化里面给的名称telephone_production
objective: 20600
desk_production=300.000
cell_production=850.000

特殊情况——无解及解决方案

注意事项

官网上的例子用的是desk_production >= 1000，这其实是不对的，不方便后面处理无解的方法，正确的应该是改cell_production，好承接后文修改。且即使按照官网做的话，本质上并没有修改一开始模型(m)的约束 $ct_-assembly = m.add_-constraint( 0.2 * desk + 0.4 * cell <= 400)$ ，并没有改成 $ct_-assembly = m.add_-constraint( 0.2 * desk + 0.4 * cell <= 1000)$ ——因为无解的情况修改的是模型m的副本im。这一点根据官网给出修改后的方案结果可以验证。
即使直接在原始模型m上的约束条件改成 $ct_-assembly = m.add_-constraint( 0.2 * desk + 0.4 * cell <= 1000)$ 也无法求得官方显示的正确结果。
因此，本样例中选取的无解情况改的是cell_production。

无解的情况

显然把上述例子中cell_production值改大，比如说cell_production >= 1000就得不到可行解，此时结果返回为None：

# 创建一个模型，直接复制上面那个
im = model.copy()

# 因为模型是复制过来的，已经没有上面那个desk_production变量了，现在就是要把这个变量给拿出来，重新换了个idesk的名字来使用
icell = im.get_var_by_name('cell_production')

# 加入新的约束——该约束实际上是覆盖了上面的那个desk_production >= 100的条件，但是desk_production >= 100并没有删除！！！
im.add_constraint(icell >= 1000)

# 输出一下此时的im模型信息
im.print_information()

# 如果没有得到可行解，就返回None进入if去输出
ims = im.solve()
if ims is None:
    print('- model is infeasible')

输出结果：
在这里插入图片描述
对比上面的输出结果可以发现constraint的数量变多了，说明原先的desk_production >= 100限制并没有删除，所以说这种方法没办法删除constraint。而输出“- model is infeasible”表明该条件下确实没有可行解。

解决方案

一句话总结：放宽限制约束，硬约束转软约束

硬约束：在任何情况都不能违反的约束，上述的所有约束条件都是硬约束
软约束：在某些情况可以违反的约束

如何转化：

比如说想要把约束 $0.2desk_-production + 0.4cell_-production \le 400$ 改为 $0.2desk_-production + 0.4cell_-production \le 440$

原始硬约束：
$\\ 0.2desk_-production + 0.4cell_-production \le 400$
转化成软约束：
$\\ 0.2desk_-production + 0.4cell_-production \le 400 + overtime$
同时要限制overtime硬约束： $\le 40$

由于变成了软约束，需要引入惩罚项。假定超过400后的额外费用为$2/h，优化目标就需要修改为：
$maximize: 12*desk_-production + 20*cell_-prodction - 2*overtime$

对应代码（直接在上述代码后面加）：

# 增加硬约束overtime
overtime = im.continuous_var(name='overtime', ub=40)

# 根据条件等式名称获得ct_assembly，并对其等式右项进行修改
# 显然ct_assembly.rhs代表等式右项，那么ct_assembly.lhs代表等式左项
ct_assembly = im.get_constraint_by_name('ct_assembly')
ct_assembly.rhs = 400 + overtime

# 增加惩罚项，修改优化目标
im.maximize(12*desk_production + 20*cell_production - 2*overtime)

imsol = im.solve()
print(imsol)

输出结果：
在这里插入图片描述

显然这是改了cell_production的取值范围，与上文注意事项相互印证。

无界变量vs无界模型

定义与关系

无界变量：变量没有上界||没有下界
无界模型：模型的目标没有上界||没有下界，说明限制条件少了
关系：相互不影响

Docplex中的默认设置：

变量：无上界，下确界为0

如果不对上述例子中的3、4项约束，解空间会是如下：
在这里插入图片描述

解决线性规划问题（LP）的算法（3种）

单纯形法（Simplex Optimizer）

证明不细究，主要是给一些结论，证明可以参考这一篇，具体怎么算可以参考这一篇

LP问题的解集是凸集
凸集的顶点一定是可行解

单纯形法核心思想：找到每一个可行解，带入目标函数后，按照要求取最大/小值即可
几何解释：由于解集是凸集，只要找到顶点处也就是最极端的条件，单纯形法本质是遍历凸集的每一个顶点

对偶单纯形法（Dual-simplex Optimizer）

对偶：每个线性规划问题（假设叫P问题）都有一个关联的LP问题（假设叫D问题），那么P称之为原问题，D称之为P的对偶
特点：如果P是一个最小化问题，则D是一个最大化问题，反之亦然
对偶价格（dual prices）/影子价格：对偶变量的值

怎么从原问题P推出对偶：
原问题：
$\quad c^Tx \\subject \quad to \quad \quad Ax \ge b \\x \ge 0$
对偶：
$\quad y^Tb \\subject \quad to \quad \quad y^TA \le c^T \\y \ge 0$

示例：拿的这篇的
在这里插入图片描述

特性：

P问题中的约束都有一个关联的对偶变量 $y_i$
D的任何可行解都是P的上限，P的任何可行解都是D的下界
LP问题中D和P的最终目标本质是等价的，都是在边界处
如果原始问题很困难不便于解决，可以考虑从对偶问题来解

在这里插入图片描述

障碍法（Barrier Optimizer）

证明暂时放过，以后看了再补充
几何解释：从可行区域的某个地方开始，采用了预测-校正算法，通过可行区域的中间路径不断调整路径。相比单纯形法沿着边缘走，这个是直接从中间穿

在这里插入图片描述

调用方法

代码示例

# 修改lpmethod的值，不同的值采用不同的算法，lpmethod=4时采用的是用的障碍法求解
im.parameters.lpmethod = 4

# 求解并输出log
im.solve(log_output=True)

输出结果：
在这里插入图片描述

参数取值：
鉴于docplex是直接封装的cplex，底层的值应该没改，对应于docplex里面的lpmethod应该同cplex（但是我在docplex里面没找到直接查找的参数说明，要是有大佬知道麻烦踢我一脚，不甚感激）

>>> print(c.parameters.lpmethod.help())
method for linear optimization  :
  0 = automatic
  1 = primal simplex
  2 = dual simplex
  3 = network simplex
  4 = barrier
  5 = sifting
  6 = concurrent optimizers

参考链接：

线性规划官方Docplex示例
 Python 中 CPLEX 参数的帮助

BoJack_Horse

关注

1
点赞
踩
9

收藏

觉得还不错? 一键收藏
0
评论
Docplex入门（1）——线性规划

线性约束：linear−expressionsymbollinear−expressionlinear_ -expression \quad symbol \quad linear_-expressionlinear−expressionsymbollinear−expression其中symbol仅能取===、≤\le≤、≥\ge≥一般形式：min⁡Cxts.t.Ax≥Bx≥0\min C^t _x\\s.t. \quad Ax \ge B\\x \ge 0minCxts.
复制链接

扫一扫