运筹学习笔记之列生成

列生成算法介绍

1. 什么是列生成

列生成算法是一种用于解决大规模线性规划问题的高效算法，它基于单纯形法的思想，通过求解子问题来找到可以进基的非基变量。在列生成算法中，每个变量都代表一列，因此称为列生成算法。该算法的优点在于其高效的计算性能和较好的收敛性，适用于处理大规模、复杂的线性规划问题。

在列生成算法的迭代过程中，因为会不断有变量入基，所以会导致限制主问题的列不断增加，所以叫做列生成算法。

2. 列生成的应用范围

列生成被广泛应用于调度问题、切割问题、车辆路径问题、选址问题等。 该算法的优点在于其高效的计算性能和较好的收敛性，适用于处理大规模、复杂的线性规划问题。对于变量数目很多的线性优化问题，单纯形法速度很慢，可以用到列生成方法来加快求解速度。

3. 列生成的原理

基本思路如下：

1、先把原问题限制到一个规模更小的限制主问题，在限制主问题的基础上用单纯形法求解，但此时的解并不是主问题的最优解。

2、通过一个子问题去检查那些未被考虑的变量中是否有使得reduced cost小于0的？如果有，就把这个变量的相关系数列加入到限制主问题的系数矩阵中，回到第一步。

经过反复迭代，知道子问题的reduced cost rate大于等于0，那么主问题就求到了最优解。

4. 基本概念

受限主问题

$$ min(y_1+y_2+\cdots+y_k) $$

$s.t.$ $$ R_1: a_{11}y_1+\cdots+a_{1k}y_k \geq b_1\ R_2: a_{21}y_1+\cdots+a_{2k}y_k \geq b_2\ \cdots\ R_m:a_{m1}y_1+\cdots+a_{mk}y_k \geq b_m $$ 就是从主问题中选取k个变量构成的松弛问题。

子问题

通过求解RMP问题或者RMP对偶问题后，得到了想要的$c_BB^{-1}$以后，子问题就是通过$\sigma_j = c_j - c_BB^{-1}a_j$

这条公式，在$y_{k+1}、y_m$中寻找检验数为负且最小的变量，将变量对应的那一列添加到RMP中。

列生成算法实现

1. 案例求解

问题描述

cutting stock problem是采用列生成算法求解的经典案例。该问题如下：有一些纸筒，每个纸筒长度为16m。顾客需要25个3m，20个6m，15个7m的纸筒。要求在满足顾客需求的情况下，裁剪的纸筒数最小。

建立模型

可以采用启发式算法求解到一个初始解如下：

5个纸筒采用切割方案1：3，3，3，3，3；

10个纸筒采用切割方案2：6，6；

8个纸筒采用切割方案3：7，7；

总计23个纸筒。可以看出，采用23个纸筒是一定可以满足要求的，这是问题的一个上界（upper bound）。

该问题可行的切割方案很多，我们可以采用P表示所有可行裁剪方案的集合，里面的方案总数为n（这里的n并不需要知道其确切取值，一般来说n是一个很大的数）。$a_{ij}$表示第$j$种方案里类别$i$的个数，$y_j$表示第$j$种方案的选择个数。建立数学模型如下： $$ min \quad y_1+y_2+\cdots+y_n $$ $s.t.$ $$ R1:a_{11}y_1 + \cdots + a_{1n}y_n \geq 25\ R2:a_{21}y_1 + \cdots + a_{2n}y_n \geq 20\ R3:a_{31}y_1 + \cdots + a_{3n}y_n \geq 15\ $$ 约束中的每一列对应的是一种切割方案，前三种方案已知，就是我们前面采用启发式算法求解的三种方案。其它的$n-3$种切割方案未知。

问题求解

第一轮循环

首先从上述模型中选出一些列，构成问题的限制主问题，这里我们选取前3列。则限制主问题如下 $$ min \quad y_1 + y_2 + y_3 \ 5y_1 + 0 y_2 + 0y_3 \geq 25\ 0y_1 + 2 y_2 + 0y_3 \geq 20\ 0y_1 + 0 y_2 + 2y_3 \geq 15\ $$ 使用python+gurobi求解上述线性规划问题

import gurobipy as gp
# 实例化模型
m = gp.Model('column generation')
y = m.addVars(3,vtype = gp.GRB.CONTINUOUS,name='y',lb=0)
m.addConstr(5*y[0] >= 25)
m.addConstr(2*y[1] >= 20)
m.addConstr(2*y[2] >= 15)
m.setObjective(y[0]+y[1]+y[2],gp.GRB.MINIMIZE)
m.setParam('OutputFlag',0)
m.optimize()
m.getAttr(gp.GRB.Attr.Pi)

可以求解出对偶变量的取值为$c_BB^{-1}=[0.2,0.5,0.5]$。现在要找一列加入RMP，但是并不知道其取值，我们将其记作$\alpha_4 = [a_{14},a_{24},a_{34}]^T$。得到非基变量的检验数为$\sigma_4 = c_4 - c_BB^{-1}\alpha_4 = 1-0.2\alpha_{14} - 0.5\alpha_{24} - 0.5\alpha_{34}$

从而构造子问题如下： $$ min(1-0.2a_{14} - 0.5a_{24} - 0.5a_{34}) $$ $s.t.$ $$ 3a_{14}+6a_{24}+7a_{34} \leq 16\ a_{ij} \in Z $$ 同样采用python+gurobi求解上述问题

import gurobipy as gp
m1 = gp.Model('subproblem')
a = m1.addVars(3,vtype=gp.GRB.INTEGER,lb=0,name='a')
m1.setObjective(1-0.2*a[0]-0.5*a[1]-0.5*a[2])
m1.addConstr(3*a[0]+6*a[1]+7*a[2] <= 16)  # 列生成规则
m1.setParam('OutputFlag',0)
m1.optimize()
print(a[0].x,a[1].x,a[2].x)

求解出$\alpha_4 = [1,2,0]^T$，检验数为$\sigma_4 = c_4 - c_BB^{-1}\alpha_4 = 1-0.2\times1 - 0.5\times2 - 0.5\times0 = -0.2 < 0$，因为检验数小于0，所以需要将$y_4$入基，添加到主问题中，开始第二轮迭代。

第二轮循环

加入$y_4$后，限制主问题变为 $$ min \quad y_1 + y_2 + y_3 + y_4 \ 5y_1 + 0 y_2 + 0y_3 + 1y_4\geq 25\ 0y_1 + 2 y_2 + 0y_3 + 2y_4 \geq 20\ 0y_1 + 0 y_2 + 2y_3 + 0y_4 \geq 15\ $$ 限制主问题与之前相比，多了一列，所以被称为列生成算法。求解该限制主问题，其对偶变量取值为$c_BB^{-1}=[0.2,0.4,0.5]$，下一个需要入基的变量记为$\alpha_5 = [a_{15},a_{25},a_{35}]^T$，构造子问题如下 $$ min(1-0.2a_{15} - 0.4a_{25} - 0.5a_{35}) $$ $s.t.$ $$ 3a_{15}+6a_{25}+7a_{35} \leq 16\ a_{ij} \in Z $$ 采用同样方法求解该子问题，得到结果$\alpha_5 = [3,0,1]^T$，检验数为$\sigma_5 = c_5 - c_BB^{-1}\alpha_5 = 1-0.2\times3 - 0.4\times0 - 0.5\times1 = -0.1 < 0$，因为检验数小于0，所以需要将$y_5$入基，添加到主问题中，开始下一轮迭代。

第三轮循环

加入$y_5$后，限制主问题变为 $$ min \quad y_1 + y_2 + y_3 + y_4 + y_5 \ 5y_1 + 0 y_2 + 0y_3 + 1y_4 + 3y_5\geq 25\ 0y_1 + 2 y_2 + 0y_3 + 2y_4 + 0y_5 \geq 20\ 0y_1 + 0 y_2 + 2y_3 + 0y_4 + 1y_5 \geq 15\ $$ 还是求解限制主问题，其对偶变量取值为$c_BB^{-1}=[0.1667,0.4167,0.5]$，下一个需要入基的变量记为$\alpha_6 = [a_{16},a_{26},a_{36}]^T$，构造子问题如下 $$ min(1-0.1667a_{16} - 0.4167a_{26} - 0.5a_{36}) $$ $s.t.$ $$ 3a_{16}+6a_{26}+7a_{36} \leq 16\ a_{ij} \in Z $$ 求得结果$\alpha_5 = [1,1,1]^T$，检验数$\sigma_6= c_6 - c_BB^{-1}\alpha_6 = 1-0.1667\times1- 0.4167\times1 - 0.5\times1 = -0.08333 < 0$，所以将$y_6$代入

第四轮循环

加入$y_6$后，限制主问题变为 $$ min \quad y_1 + y_2 + y_3 + y_4 + y_5 + y_6 \ 5y_1 + 0 y_2 + 0y_3 + 1y_4 + 3y_5 + y_6\geq 25\ 0y_1 + 2 y_2 + 0y_3 + 2y_4 + 0y_5 + y_6 \geq 20\ 0y_1 + 0 y_2 + 2y_3 + 0y_4 + 1y_5 + y_6 \geq 15\ $$ 还是求解限制主问题，其对偶变量取值为$c_BB^{-1}=[0.2,0.4,0.4]$，下一个需要入基的变量记为$\alpha_7 = [a_{17},a_{27},a_{37}]^T$，构造子问题如下 $$ min(1-0.2a_{17} - 0.4a_{27} - 0.4a_{37}) $$ $s.t.$ $$ 3a_{17}+6a_{27}+7a_{37} \leq 16\ a_{ij} \in Z $$ 求得结果$\alpha_5 = [5,0,0]^T$，检验数$\sigma_6= c_6 - c_BB^{-1}\alpha_6 = 1-0.2\times5- 0.4\times0 - 0.4\times0 = 0$，结束迭代，此时列生成算法结束。

此时，我们将$y_7$代入模型，此时求解结果就是最优解。

求解出最优解为$y=[1,0,0,3,1,14,0]$，此时目标函数值为19。

所以我们得到的最终切割方案为