1. 应用场景
列生成算法是不断添加求解变量的算法，可参考论文。上一节我们介绍了约束条件存在分块结构的情形下可以使用单纯型法的判断规则进行列生成，下面我们介绍另一种场景，set-partitionning问题的列生成算法。 具体来说，场景如下：有
n
n
n个0-1变量
z
1
.
.
.
z
n
z_1...z_n
z1...zn，每个
z
i
z_i
zi带着很多中间变量
x
i
,
j
x_{i,j}
xi,j用来进行约束，这是个变量很多，约束也很多的模型。我们首先对问题模型进行转化，将所有
z
i
z_i
zi的组合枚举出来，在枚举的过程中，就把约束条件都过了一遍，只有满足所有约束条件的组合才会保留下来。最后的枚举结果分别对应
y
1
.
.
.
y
m
y_1...y_m
y1...ym。 不难看出，
m
m
m是
n
n
n的指数形式，使用set-partitioning模型之后，变量的数量大的可怕。好消息是，原问题的约束条件比较紧，只有少数组合成立，可以通过一定的规则产生，这样m的数量比较小，问题就比较好解决。
2. 例1：cutting stock problem
来看网上流传很多的cutting stock problem的例子，问题如下： 我们有一些纸筒，每个纸筒16m。顾客需要裁剪的长度为：25个3m，20个6m，15个7m。要求在满足顾客需求的情况下，裁剪的纸筒数最小。 我们可以枚举出每筒可能的裁剪方案： 方案1：5个3m，即3+3+3+3+3 方案2：2个6m，即6+6 方案3：2个7m，即7+7 我们可以继续枚举下去： 4：6+6+3 5：6+3+3+3 6：7+3+3+3 7：7+6+3 …
仅使用前三种方案，我们用23筒（5个方案1，10个方案2，8个方案3）是肯定可以满足要求的，这算问题的一个上界。下面我们探索一下其他的切筒方式，看能不能给出下界。
一般性的，我们假设
P
P
P是所有可行的裁剪方案的集合（其中前3个可以设置为我们前面的3个裁剪方案），里面方案的总数为
n
n
n(我们并不需要确切的知道这个值是多少，只需要知道它很大)。令
a
i
j
a_{ij}
aij表示第
j
j
j种方案里类别
i
i
i的个数，
y
j
y_j
yj表示第
j
j
j种方案的选择个数，原问题可以变为：
min
y
1
+
.
.
.
+
y
n
\min y_1+...+y_n
miny1+...+yn
a
1
,
1
y
1
+
a
1
,
2
y
2
+
.
.
.
+
a
1
,
n
y
n
≥
25
a_{1,1}y_1+a_{1,2}y_2+...+a_{1,n}y_n \geq 25
a1,1y1+a1,2y2+...+a1,nyn≥25
a
2
,
1
y
1
+
a
2
,
2
y
2
+
.
.
.
+
a
2
,
n
y
n
≥
20
a_{2,1}y_1+a_{2,2}y_2+...+a_{2,n}y_n \geq 20
a2,1y1+a2,2y2+...+a2,nyn≥20
a
3
,
1
y
1
+
a
3
,
2
y
2
+
.
.
.
+
a
3
,
n
y
n
≥
15
a_{3,1}y_1+a_{3,2}y_2+...+a_{3,n}y_n \geq 15
a3,1y1+a3,2y2+...+a3,nyn≥15 其中
a
1
,
1
=
5
a_{1,1}=5
a1,1=5，
a
2
,
1
=
0
a_{2,1}=0
a2,1=0，
a
3
,
1
=
0
a_{3,1}=0
a3,1=0
a
1
,
2
=
0
a_{1,2}=0
a1,2=0，
a
2
,
2
=
2
a_{2,2}=2
a2,2=2，
a
3
,
2
=
0
a_{3,2}=0
a3,2=0
a
1
,
3
=
0
a_{1,3}=0
a1,3=0，
a
2
,
3
=
0
a_{2,3}=0
a2,3=0，
a
3
,
3
=
2
a_{3,3}=2
a3,3=2 对应最初的3种方案。
a
i
,
j
a_{i,j}
ai,j满足条件：
3
∗
a
1
,
j
+
6
∗
a
2
,
j
+
7
∗
a
3
,
j
≤
16
3*a_{1,j}+6*a_{2,j}+7*a_{3,j}\leq 16
3∗a1,j+6∗a2,j+7∗a3,j≤16（前面所说的列生成的规则） 注意我们在使用列生成算法的时候，将整数变量松弛为了连续变量。
用矩阵形式重新表达如下：
min
Σ
j
y
j
\min \Sigma_{j} y_j
minΣjyj s.t.
A
y
≥
b
T
Ay ≥ b^T
Ay≥bT
A
A
A是一个矩阵，
b
=
[
b
1
,
.
.
.
,
b
m
]
b=[b_1,...,b_m]
b=[b1,...,bm]和
y
=
[
y
1
,
.
.
.
,
y
n
]
y=[y_1,...,y_n]
y=[y1,...,yn]是向量；并且有重要的一点：set-covering的变量
y
j
y_j
yj对应的列
a
j
=
[
a
1
j
,
.
.
.
,
a
m
j
]
T
a_j=[a_{1j},...,a_{mj}]^T
aj=[a1j,...,amj]T生成的规则，可以通过线性约束条件
B
∗
a
j
≤
D
B*a_j ≤ D
B∗aj≤D得到。这使得我们可以用线性规划的方法来求解子问题~ 注意由于列生成算法求解的是松弛线性规划问题，因此对整数规划模型需要结合branch and bound算法进行。
2.1 限制主问题
首先用启发式算法找出一部分
A
A
A，比如说选出了
k
k
k列。然后我们的线性规划问题就变成了： min
y
1
+
y
2
+
.
.
.
+
y
k
y_1+y_2+...+y_k
y1+y2+...+yk s.t.
a
i
,
1
∗
y
1
+
a
i
,
2
∗
y
2
+
.
.
.
+
a
i
,
k
∗
y
k
≥
b
i
a_{i,1}*y_1+a_{i,2}*y_2+...+a_{i,k}*y_k≥ b_i
ai,1∗y1+ai,2∗y2+...+ai,k∗yk≥bi，
i
=
1...
n
i = 1...n
i=1...n 相比原来的模型，相当于把
y
k
+
1
y_{k+1}
yk+1~
y
m
y_m
ym强制限制为非基变量了，称为限制主问题（Restricted Master Problem，RMP）。上面的限制主问题求解完成后，我们想使用单纯型法进行基变量的转换，看看
y
k
+
1
y_{k+1}
yk+1~
y
m
y_m
ym中，是否有可以转入基变量的列。检验数
σ
=
c
N
−
π
∗
a
j
\sigma = c_N - \pi * a_j
σ=cN−π∗aj，并且
π
=
c
B
B
−
1
π = c_BB^{-1}
π=cBB−1可以由求解器给出。我们要找出非基变量中最小的负数
σ
\sigma
σ，将其转入基变量。正如前面所述，我们这里使用线性规划来找这个新的列。
2.2 子问题
注意
c
=
[
1
,
.
.
.
,
1
]
c = [1,...,1]
c=[1,...,1]，因此 min
c
N
−
π
a
j
c_N - \pi a_j
cN−πaj等价于 min
1
−
π
∗
a
1 - \pi *a
1−π∗a s.t.
B
∗
a
≤
D
B*a ≤ D
B∗a≤D (列生成的规则) 称为子问题（sub problem）。如果目标函数最优值＜0，就将新生成的列yk+1转入基变量，生成新的限制主问题进行求解。如此往复，直至子问题的目标函数≥0。
2.3 求解质量分析
一般性的，我们假设原问题为
Z
=
max
c
x
Z=\max cx
Z=maxcx
A
x
=
b
,
x
≥
0
Ax=b,x\ge0
Ax=b,x≥0
RMP为
max
c
x
ˉ
\max c\bar{x}
maxcxˉ
A
x
ˉ
=
b
,
x
ˉ
≥
0
A\bar{x}=b,\bar{x}\ge0
Axˉ=b,xˉ≥0。 令
z
ˉ
\bar{z}
zˉ为RMP的最优解，
z
∗
z^*
z∗为原问题的最优解，我们有
z
∗
≤
z
ˉ
z^*\le \bar{z}
z∗≤zˉ
SP求得的最佳残差（负数）为
c
ˉ
∗
\bar{c}^*
cˉ∗，并令
k
≥
Σ
x
k\ge \Sigma x
k≥Σx。
z
ˉ
\bar{z}
zˉ同样也是对偶问题的最优解，因此我们有
z
ˉ
+
c
ˉ
∗
k
≤
z
∗
\bar{z}+ \bar{c}^*k\le z^*
zˉ+cˉ∗k≤z∗
在set-partitioning问题中，我们可以取
k
=
Σ
x
=
x
∗
k= \Sigma x=x^*
k=Σx=x∗，因此有
z
ˉ
/
(
1
−
c
ˉ
∗
)
≤
z
∗
≤
z
ˉ
\bar{z}/(1- \bar{c}^*)\le z^*\le \bar{z}
zˉ/(1−cˉ∗)≤z∗≤zˉ。同样，我们可以参考simplex的启发算法，令 则
2.4 迭代过程
2.4.1 第一轮
问题对应的初始解为
y
1
=
5
y_1=5
y1=5，
y
2
=
10
y_2=10
y2=10，
y
3
=
8
y_3=8
y3=8，
y
4
=
.
.
.
y
n
=
0
y_4=...y_n=0
y4=...yn=0 限制主问题（restricted master problem）如下：
min
y
1
+
y
2
+
y
3
\min y_1+y_2+y_3
miny1+y2+y3
5
y
1
+
0
y
2
+
0
y
3
≥
25
5 y_1+0 y_2+0y_3 \geq 25
5y1+0y2+0y3≥25
0
y
1
+
2
y
2
+
0
y
3
≥
20
0y_1+2y_2+0y_3 \geq 20
0y1+2y2+0y3≥20
0
y
1
+
0
y
2
+
2
y
3
≥
15
0y_1+0y_2+2y_3 \geq 15
0y1+0y2+2y3≥15
下面是pymprog的sensitivity输出： 对应
π
=
(
0.2
,
0.5
,
0.5
)
\pi = (0.2, 0.5, 0.5)
π=(0.2,0.5,0.5)，我们把新转入的变量记作
a
4
=
(
a
1
,
4
,
a
2
,
4
,
a
3
,
4
)
a_4 = (a_{1,4},a_{2,4},a_{3,4})
a4=(a1,4,a2,4,a3,4)，则子问题为： min
1
−
0.2
∗
a
1
,
4
−
0.5
∗
a
2
,
4
−
0.5
∗
a
3
,
4
1 - 0.2*a_{1,4}-0.5*a_{2,4}-0.5*a_{3,4}
1−0.2∗a1,4−0.5∗a2,4−0.5∗a3,4 s.t.
3
∗
a
1
,
4
+
6
∗
a
2
,
4
+
7
∗
a
3
,
4
≤
16
3*a_{1,4}+6*a_{2,4}+7*a_{3,4}\leq 16
3∗a1,4+6∗a2,4+7∗a3,4≤16 (列生成的规则)
a
i
,
j
∈
Z
a_{i,j}\in Z
ai,j∈Z 求解结果为
a
4
=
(
1
,
2
,
0
)
a_4=(1,2,0)
a4=(1,2,0)，检验数
σ
=
1
−
π
∗
a
4
=
−
0.2
≤
0
\sigma = 1 - \pi * a_4 =-0.2\le 0
σ=1−π∗a4=−0.2≤0。此时问题目标为22.5。 将
y
4
y_4
y4进基，把这个结果添加到主问题当中去，开始第二轮迭代。
2.4.2 第二轮
添加了
y
4
y_4
y4后，限制主问题为：
min
y
1
+
y
2
+
y
3
+
y
4
\min y_1+y_2+y_3+y_4
miny1+y2+y3+y4
5
y
1
+
0
y
2
+
0
y
3
+
1
y
4
≥
25
5 y_1+0 y_2+0y_3+1y_4 \geq 25
5y1+0y2+0y3+1y4≥25
0
y
1
+
2
y
2
+
0
y
3
+
2
y
4
≥
20
0y_1+2y_2+0y_3+2y_4 \geq 20
0y1+2y2+0y3+2y4≥20
0
y
1
+
0
y
2
+
2
y
3
+
0
y
4
≥
15
0y_1+0y_2+2y_3+0y_4 \geq 15
0y1+0y2+2y3+0y4≥15
直观上来看，是添加了一个新列，所以称为列生成算法。上面的主问题求解结果为： 对应
π
=
(
0.2
,
0.4
,
0.5
)
\pi = (0.2, 0.4, 0.5)
π=(0.2,0.4,0.5)，我们把新转入的变量记作
a
5
=
(
a
1
,
5
,
a
2
,
5
,
a
3
,
5
)
a_5 = (a_{1,5},a_{2,5},a_{3,5})
a5=(a1,5,a2,5,a3,5)，则子问题为： min
1
−
0.2
∗
a
1
,
5
−
0.4
∗
a
2
,
5
−
0.5
∗
a
3
,
5
1 - 0.2*a_{1,5}-0.4*a_{2,5}-0.5*a_{3,5}
1−0.2∗a1,5−0.4∗a2,5−0.5∗a3,5 s.t.
3
∗
a
1
,
5
+
6
∗
a
2
,
5
+
7
∗
a
3
,
5
≤
16
3*a_{1,5}+6*a_{2,5}+7*a_{3,5}\leq 16
3∗a1,5+6∗a2,5+7∗a3,5≤16 (列生成的规则)
a
i
,
j
∈
Z
a_{i,j}\in Z
ai,j∈Z 求解结果为
a
5
=
(
3
,
0
,
1
)
a_5=(3,0,1)
a5=(3,0,1)，检验数
σ
=
1
−
π
∗
a
5
=
−
0.1
<
0
\sigma = 1 - \pi * a_5=-0.1< 0
σ=1−π∗a5=−0.1<0。此时问题目标为20.5。 将
y
5
y_5
y5进基，把这个结果添加到主问题当中去，开始第三轮迭代。
2.4.3 第三轮
添加了
y
5
y_5
y5后，限制主问题为：
min
y
1
+
y
2
+
y
3
+
y
4
+
y
5
\min y_1+y_2+y_3+y_4+y_5
miny1+y2+y3+y4+y5
5
y
1
+
0
y
2
+
0
y
3
+
1
y
4
+
3
y
5
≥
25
5 y_1+0 y_2+0y_3+1y_4 +3y_5 \geq 25
5y1+0y2+0y3+1y4+3y5≥25
0
y
1
+
2
y
2
+
0
y
3
+
2
y
4
+
0
y
5
≥
20
0y_1+2y_2+0y_3+2y_4 +0y_5\geq 20
0y1+2y2+0y3+2y4+0y5≥20
0
y
1
+
0
y
2
+
2
y
3
+
0
y
4
+
1
y
5
≥
15
0y_1+0y_2+2y_3+0y_4 +1y_5\geq 15
0y1+0y2+2y3+0y4+1y5≥15
上面的主问题求解结果为： 对应
π
=
(
0.166667
,
0.416667
,
0.5
)
\pi = (0.166667, 0.416667, 0.5)
π=(0.166667,0.416667,0.5)，我们把新转入的变量记作
a
6
=
(
a
1
,
6
,
a
2
,
6
,
a
3
,
6
)
a_6 = (a_{1,6},a_{2,6},a_{3,6})
a6=(a1,6,a2,6,a3,6)，则子问题为： min
1
−
0.166667
∗
a
1
,
6
−
0.416667
∗
a
2
,
6
−
0.5
∗
a
3
,
6
1 - 0.166667*a_{1,6}-0.416667*a_{2,6}-0.5*a_{3,6}
1−0.166667∗a1,6−0.416667∗a2,6−0.5∗a3,6 s.t.
3
∗
a
1
,
6
+
6
∗
a
2
,
6
+
7
∗
a
3
,
6
≤
16
3*a_{1,6}+6*a_{2,6}+7*a_{3,6}\leq 16
3∗a1,6+6∗a2,6+7∗a3,6≤16 (列生成的规则)
a
i
,
j
∈
Z
a_{i,j}\in Z
ai,j∈Z 求解结果为
a
6
=
(
1
,
1
,
1
)
a_6=(1,1,1)
a6=(1,1,1)，检验数
σ
=
1
−
π
∗
a
6
=
−
0.083334
<
0
\sigma = 1 - \pi * a_6=-0.083334< 0
σ=1−π∗a6=−0.083334<0。此时问题目标为20。 将
y
6
y_6
y6进基，把这个结果添加到主问题当中去，开始第三轮迭代。
2.4.4 第四轮
添加了
y
6
y_6
y6后，限制主问题为：
min
y
1
+
y
2
+
y
3
+
y
4
+
y
5
+
y
6
\min y_1+y_2+y_3+y_4+y_5+y_6
miny1+y2+y3+y4+y5+y6
5
y
1
+
0
y
2
+
0
y
3
+
1
y
4
+
3
y
5
+
y
6
≥
25
5 y_1+0 y_2+0y_3+1y_4 +3y_5+y_6 \geq 25
5y1+0y2+0y3+1y4+3y5+y6≥25
0
y
1
+
2
y
2
+
0
y
3
+
2
y
4
+
0
y
5
+
y
6
≥
20
0y_1+2y_2+0y_3+2y_4 +0y_5+y_6\geq 20
0y1+2y2+0y3+2y4+0y5+y6≥20
0
y
1
+
0
y
2
+
2
y
3
+
0
y
4
+
1
y
5
+
y
6
≥
15
0y_1+0y_2+2y_3+0y_4 +1y_5+y_6\geq 15
0y1+0y2+2y3+0y4+1y5+y6≥15
上面的主问题求解结果为： 对应
π
=
(
0.2
,
0.4
,
0.4
)
\pi = (0.2, 0.4, 0.4)
π=(0.2,0.4,0.4)，我们把新转入的变量记作
a
7
=
(
a
1
,
7
,
a
2
,
7
,
a
3
,
7
)
a_7 = (a_{1,7},a_{2,7},a_{3,7})
a7=(a1,7,a2,7,a3,7)，则子问题为： min
1
−
0.2
∗
a
1
,
7
−
0.4
∗
a
2
,
7
−
0.4
∗
a
3
,
7
1 - 0.2*a_{1,7}-0.4*a_{2,7}-0.4*a_{3,7}
1−0.2∗a1,7−0.4∗a2,7−0.4∗a3,7 s.t.
3
∗
a
1
,
7
+
6
∗
a
2
,
7
+
7
∗
a
3
,
7
≤
16
3*a_{1,7}+6*a_{2,7}+7*a_{3,7}\leq 16
3∗a1,7+6∗a2,7+7∗a3,7≤16 (列生成的规则)
a
i
,
j
∈
Z
a_{i,j}\in Z
ai,j∈Z 求解结果为
a
7
=
(
3
,
1
,
0
)
a_7=(3,1,0)
a7=(3,1,0)，检验数
σ
=
1
−
π
∗
a
6
=
0
\sigma = 1 - \pi * a_6= 0
σ=1−π∗a6=0，结束迭代，此时问题目标为19。
2.4.5 LP和MIP的差别
我们把下面的模型称为Kantorovich模型：
LP问题可以给出一个下界，bound为 这个bound很直观，即假设纸可以拼接。 我们换一个例子
使用LP得到bound为： 和最优值156差别还是比较大的。
2.5 Julia调用GLPK代码示例
代码如下：
A = [[5,0,0],[0,2,0],[0,0,2]]
b = [25 20 15];
while true
# main model
model = Model(GLPK.Optimizer)
@variable(model, y[1:length(A[1])]>=0)
@objective(model, Min, sum(y))
r=@constraint(model,[i=1:3],sum(A[i].*y)>=b[i]);
optimize!(model)
# sub model
submodel = Model(GLPK.Optimizer)
@variable(submodel, a[1:3]>=0, Int)
@objective(submodel, Min, 1-sum(dual.(r).*a))
@constraint(submodel,sum([3,6,7].*a)<=16);
optimize!(submodel)
delta = objective_value(submodel)
if delta<-0.001
A = [[A[i];Int(JuMP.value(a[i]))] for i in 1:3]
@info objective_value(model),"add",Int.(JuMP.value.(a))
else
@info objective_value(model)
break
end
end
经过3轮迭代得到最优解
2.6 python调用cplex代码示例
使用下列代码求出最优解：
from pymprog import *
p = model('example')
y = p.var('y',6,int) # variables
p.minimize(sum(y), 'master')
r1 = 5*y[0]+ y[3]+3*y[4]+y[5]>= 25
r2 = 2*y[1]+2*y[3]+y[5] >= 20
r3 = 2*y[2] +y[4]+y[5]>= 15
p.solve()
for i in range(6):
print(y[i].primal)
求解结果为：
y
1
=
1
,
y
4
=
3
,
y
5
=
1
,
y
6
=
14
y_1=1,y_4=3,y_5=1,y_6=14
y1=1,y4=3,y5=1,y6=14，总计需要19个纸筒。 也可以将最后一轮的
y
y
y取整得到一个可行解：
y
1
=
2
,
y
4
=
3
,
y
6
=
15
y_1=2,y_4=3,y_6=15
y1=2,y4=3,y6=15，需要20个纸筒。 这里给出第四轮的限制主问题和子问题的代码供参考：
from pymprog import *
p = model('master problem')
y = p.var('y',6) # variables
p.minimize(sum(y), 'master')
r1 = 5*y[0]+ y[3]+3*y[4]+y[5]>= 25
r2 = 2*y[1]+2*y[3]+y[5] >= 20
r3 = 2*y[2] +y[4]+y[5]>= 15
p.solve()
p.sensitivity()
from pymprog import *
s = model('sub problem')
a = s.var('a',3) # variables
for i in range(3):
a[i].kind = int
s.minimize(1-(0.2*a[0] + 0.4*a[1] +0.4*a[2]), 'sub')
3*a[0] +6*a[1] +7*a[2]<= 16 # mountain bike limit
s.solve()
print(s.vobj())
for i in range(3):
print(a[i].primal)
然后下面给出直接求最优解的代码。pymprog用了一天都没有求出来，因此改用cplex进行求解：
import cplex
from cplex.exceptions import CplexError
my_obj = [0.0,0.0,0.0,1.0]*23
my_ub = [5.0,2.0,2.0,1.0]*23
my_lb = [0.0,0.0,0.0,0.0]*23
my_ctype = "IIII"*23
my_colnames = []
my_rhs = [25,20,15]+[0,16]*23
my_sense = "GGG"+"LL"*23
for i in range(23):
my_colnames = my_colnames + ["x0,"+str(i),"x1,"+str(i),"x2,"+str(i),"y"+str(i)]
def populatebyrow(prob):
prob.objective.set_sense(prob.objective.sense.minimize)
prob.variables.add(obj=my_obj, lb=my_lb, ub=my_ub, types=my_ctype,names=my_colnames)
row0 = [1.0,0.0,0.0,0.0]*23
row1 = [0.0,1.0,0.0,0.0]*23
row2 = [0.0,0.0,1.0,0.0]*23
rows = [[my_colnames, row0],
[my_colnames, row1],
[my_colnames, row2]]
for i in range(23):
row3 = [0.0,0.0,0.0,0.0]*23
row3[i*4:i*4+4]=[1.0,1.0,1.0,-100.0]
row4 = [0.0,0.0,0.0,0.0]*23
row4[i*4:i*4+4]=[3.0,6.0,7.0,0.0]
rows = rows+[[my_colnames, row3],
[my_colnames, row4]]
prob.linear_constraints.add(lin_expr=rows, senses=my_sense,rhs=my_rhs)
try:
my_prob = cplex.Cplex()
handle = populatebyrow(my_prob)
my_prob.solve()
except CplexError as exc:
print(exc)
print("Solution status = ", my_prob.solution.status[my_prob.solution.get_status()])
print("Solution value = ", my_prob.solution.get_objective_value())
x = my_prob.solution.get_values()
for i in range(23):
if x[i*4+3]>0:
print(x[i*4:i*4+3])
结果如下：
CPXPARAM_Read_DataCheck 1
Tried aggregator 1 time.
MIP Presolve modified 23 coefficients.
Reduced MIP has 49 rows, 92 columns, and 230 nonzeros.
Reduced MIP has 23 binaries, 69 generals, 0 SOSs, and 0 indicators.
Presolve time = 0.00 sec. (0.12 ticks)
Found incumbent of value 22.000000 after 0.00 sec. (0.48 ticks)
Probing time = 0.00 sec. (0.01 ticks)
Tried aggregator 1 time.
Reduced MIP has 49 rows, 92 columns, and 230 nonzeros.
Reduced MIP has 23 binaries, 69 generals, 0 SOSs, and 0 indicators.
Presolve time = 0.00 sec. (0.21 ticks)
Probing time = 0.00 sec. (0.01 ticks)
MIP emphasis: balance optimality and feasibility.
MIP search method: dynamic search.
Parallel mode: deterministic, using up to 4 threads.
Root relaxation solution time = 0.00 sec. (0.13 ticks)
Nodes Cuts/
Node Left Objective IInf Best Integer Best Bound ItCnt Gap
* 0+ 0 22.0000 0.0000 100.00%
0 0 6.6667 38 22.0000 6.6667 30 69.70%
0 0 10.2791 39 22.0000 Cuts: 62 143 53.28%
0 0 18.0972 40 22.0000 Cuts: 80 213 17.74%
* 0+ 0 21.0000 18.0972 13.82%
0 0 18.4688 35 21.0000 Cuts: 38 305 12.05%
0 0 18.7500 19 21.0000 Cuts: 41 395 10.71%
0 0 18.7500 12 21.0000 Cuts: 18 444 10.71%
0 0 18.8667 10 21.0000 Cuts: 12 496 10.16%
* 0+ 0 19.0000 18.8667 0.70%
0 0 cutoff 19.0000 18.8667 496 0.70%
Elapsed time = 0.18 sec. (7.56 ticks, tree = 0.01 MB, solutions = 3)
Implied bound cuts applied: 22
Mixed integer rounding cuts applied: 35
Zero-half cuts applied: 5
Gomory fractional cuts applied: 7
Root node processing (before b&c):
Real time = 0.19 sec. (7.57 ticks)
Parallel b&c, 4 threads:
Real time = 0.00 sec. (0.00 ticks)
Sync time (average) = 0.00 sec.
Wait time (average) = 0.00 sec.
------------
Total (root+branch&cut) = 0.19 sec. (7.57 ticks)
Solution status = MIP_optimal
Solution value = 19.0
[1.0, 1.0, 1.0]
[1.0, 1.0, 1.0]
[1.0, 2.0, 0.0]
[1.0, 1.0, 1.0]
[1.0, 1.0, 1.0]
[1.0, 1.0, 1.0]
[1.0, 1.0, 1.0]
[1.0, 1.0, 1.0]
[1.0, 1.0, 1.0]
[1.0, 1.0, 1.0]
[1.0, 2.0, 0.0]
[1.0, 1.0, 1.0]
[1.0, 2.0, 0.0]
[1.0, 1.0, 1.0]
[5.0, 0.0, 0.0]
[1.0, 1.0, 1.0]
[1.0, 1.0, 1.0]
[3.0, 0.0, 1.0]
[1.0, 1.0, 1.0]
最终是19个。
3. 例子2：多货物配送问题
多个货物在有路径容量限制的图上G=(V,E)进行配送。令
x
i
,
j
l
x^l_{i,j}
xi,jl为货物
l
l
l是否在路径
i
,
j
i,j
i,j上进行配送，
q
l
q^l
ql为货物的配送量，直观建模为：
b
i
l
=
1
b^l_i = 1
bil=1 if
i
i
i is the source node of commodity
l
l
l,let
b
i
l
=
−
1
b^l_i =−1
bil=−1if
i
i
i is the target node of commodity
l
l
l, and let
b
i
l
=
0
b^l_i = 0
bil=0 otherwise。 令P为path集合，我们使用set-partitioning模型重新建模：
其中
δ
i
j
p
δ^p_{ij}
δijp be a constant denoting whether or not path p visits edge
(
i
j
)
∈
E
(ij) \in E
(ij)∈E 模型中多了一个等式约束，并且子问题不能用线性规划简单求解。 首先我们来看列生成的规则，Let
π
i
j
π_{ij}
πij ≤0 be the dual of constraint (19) and
σ
l
∈
R
σ^l ∈R
σl∈R be the dual of constraint (20)， 则reduced cost for a column p 为：
我们有列生成的限制要求： 我们可以将距离矩阵从
c
i
j
c_{ij}
cij替换为
c
i
j
−
π
i
j
c_{ij}-\pi_{ij}
cij−πij，上面这个不等式的左边可以理解为寻找每个商品
l
l
l在新的距离矩阵下的最短路，然后看是否小于
σ
l
\sigma^l
σl，如果满足则将此路径p添加进限制主问题中。
4. 例子3：bin packing问题
bin的容量为
W
W
W，有
n
n
n个物品，大小为：
L
=
w
1
,
.
.
.
,
w
n
L={w_1,...,w_n}
L=w1,...,wn。将物品放在bin中，使bin的数量最少。 使用flow模型进行建模： 如上图，W=5，L={3,2}，下半部分是一条可行路径。目标函数是路径条数最小化，约束则是每个w都要有，并且中间的点需要流量平衡。 Min z s.t.