线性规划专题——SIMPLEX 单纯形算法(一)
线性规划专题——SIMPLEX 单纯形算法(二)
前面两篇博文已经把单纯形算法里面的核心思想给解释清楚了,主要是要认识到在线性规划里面的以下几点:
-
目标函数的最优值一定在可行域的顶点取得。
-
可行域的顶点对应这系数矩阵的一组基;系数矩阵的一组基也对应这一个可行域上的顶点
-
顶点的转移是通过在旧的基本列里面加入新的列,同时为了保持rank一致,再从基本列里面删去一列。在转移的时候,重点就是要求出那个 λ ⃗ \vec \lambda λ 来,它其实是使得 A λ ⃗ = 0 A\vec\lambda=0 Aλ=0 的 λ ⃗ \vec \lambda λ的解,只不过在解这个方程的时候,选择 A A A 的那组旧的基本列来求解。
-
单纯形法的终止条件是,添加任意的非基本列都不能改善目标函数,此时目标函数到达最小值。
OK,本博客来看看到底如何来用这个单纯形法来解线性规划。
一般,单纯形法会使用一个表格。使用表格来记录。我们来举几个例子。
再次使用如下记号来表示线性规划的松弛型:
几个例子
例1
假设存在如下的线性规划:
这是一个标准型的线性规划。我们添加松弛变量,得到松弛型:
这个线性规划的右边的所有 b i b_i bi 都是非负的,所以: X = [ 0 , 0 , 0 , 4 , 2 , 3 , 6 ] X=[0,0,0,4,2,3,6] X=[0,0,0,4,2,3,6] 就是满足条件的一个顶点。
我们画出下面这个表格出来:
这个表格一共有5部分组成。
第1部分,表示各个变量。
第2部分,目标函数的各个系数,这些系数是与第一部分的变量是对应起来的。 c ‾ i \overline c_i ci 与 x i x_i xi对应。
第3部分,当前得到的目标函数值的相反数。
第4部分,对应于 A X = b AX=b AX=b 的 b,它其实表示了选定基本列后,基本列对应的 x i x_i xi的值,那些非基本列的 x j x_j xj 全部为0。上面的表格说明 基本列是第4,5,6,7列,这组基对应的顶点是 X = [ 0 , 0 , 0 , 4 , 2 , 3 , 6 ] X=[0,0,0,4,2,3,6] X=[0,0,0,4,2,3,6]
第5部分,系数矩阵。 每一列与变量也是对应的,第i列表示第i个变量的系数列。
注意,我们需要始终保持基本列都是 e i e_i ei , e i e_i ei是单位阵的第 i i i列。化成这种形式是为了方便的解方程和求 λ ⃗ \vec \lambda λ。
怎么操作呢?
每次从第2部分中,选择一个负的
c
‾
i
\overline c_{i}
ci ,负的意味着把
c
i
c_i
ci对应的列添加进来以后,目标函数是会减少的量的多少。
例如,现在 x 4 , x 5 , x 6 , x 7 x_4,x_5,x_6,x_7 x4,x5,x6,x7 对应的列都是基本列,现在 c ‾ 1 \overline c_1 c1 是小于0的。那么,我们就想把这一列换到基本列去,然后把旧的基本列某一列给丢掉。
既然要换入第1列,那么我们就要解出 λ ⃗ \vec \lambda λ来,也就是要用第4,5,6,7列去表示第1列。
我们发现第1列是 [ 1 , 1 , 0 , 0 ] [1,1,0,0] [1,1,0,0],而基本列刚好是单位阵的列,所以第1列其实也就表明了如何用基本列去表出它。即我们可以很快的写出: λ ⃗ = [ − 1 , 0 , 0 , 1 , 1 , 0 , 0 ] \vec \lambda=[-1,0,0,1,1,0,0] λ=[−1,0,0,1,1,0,0].后面4个元素 就是第1列的4个元素,这得益于基本列是单位阵的列。
找到了
λ
⃗
\vec \lambda
λ,我们选择合适非负的
θ
\theta
θ,使得
X
′
=
X
−
θ
λ
⃗
X'=X-\theta \vec \lambda
X′=X−θλ 也是个顶点。为啥要选择非负的
θ
\theta
θ 呢。
这是因为:
A
λ
⃗
=
0
⃗
A\vec\lambda=\vec 0
Aλ=0,而
A
X
=
b
AX=b
AX=b,于是
A
X
′
=
A
X
−
θ
A
λ
⃗
=
b
AX'=AX-\theta A\vec\lambda=b
AX′=AX−θAλ=b 是满足第一个约束。
考虑到 X ′ = X − θ λ ⃗ X'=X- \theta \vec \lambda X′=X−θλ 里面, X = [ 0 , 0 , 0 , 4 , 2 , 3 , 6 ] X=[0,0,0,4,2,3,6] X=[0,0,0,4,2,3,6] 的每一项都是非负,而且是非基本列对应的项为0,第1列不是基本列,所以 x 1 x_1 x1也是0。而 λ ⃗ \vec \lambda λ 出了第1项为-1,其他也全部非负。当 θ ≤ 0 \theta \leq 0 θ≤0,那 X ′ X' X′的第1项 x 1 ′ = − θ λ 1 = θ < 0 x'_1=-\theta\lambda_1=\theta < 0 x1′=−θλ1=θ<0 就不满足条件了。
所以得考虑非负的
θ
\theta
θ 。
这个
θ
\theta
θ 的取法是,用表格的第4部分的每一行的元素去对应除以第1列的元素,得到他们的最小非负的商。对于
3
0
\frac{3}{0}
03 和
6
0
\frac{6}{0}
06 的情况,记他们的商为正无穷。
这样除下去,我们得到所有可能的
θ
⃗
\vec \theta
θ为
[
4
,
2
,
+
∞
,
+
∞
]
[4,2,+\infty,+\infty]
[4,2,+∞,+∞],其中最小的是
θ
2
=
2
\theta_2=2
θ2=2,于是选择
θ
=
2
\theta=2
θ=2,然后 现在 关键来了:以第1列的第2行为主元,通过高斯消元法,把第1列第2行的其他元都消掉,包括
c
‾
\overline c
c那一行,消元的时候带着第4部分一起消。
经过处理后,我们得到:
变化后:
变化前:
注意变化前后的对比。
变化后,我们获得一组新的基本列:第1列,第4列,第6列,第7列,而且每一列都是正交单位向量(这意味这把原来的第5列给踢掉了)。而且
c
‾
1
\overline c_1
c1 的值是0。这组基下的顶点的坐标是:
X
=
[
2
,
0
,
0
,
2
,
0
,
3
,
6
]
X=[2,0,0,2,0,3,6]
X=[2,0,0,2,0,3,6],此时的目标函数值为:-2。这个值通过看第3部分可得到。
接下来,我们再看新状态下,第2部分是否有小于0的
c
‾
i
\overline c_i
ci ,发现
c
‾
3
\overline c_3
c3 小于0。
于是,我们选择第3列,准备把它换到基本列去。同样的,计算
θ
⃗
=
[
2
,
+
∞
,
3
,
6
]
\vec \theta=[2,+\infty,3,6]
θ=[2,+∞,3,6],选择非负最小的那个,得到
θ
1
=
2
\theta_1=2
θ1=2。于是以第3列第1行为主元,把主元所在的列的其他元素全部消掉,包括
c
‾
3
\overline c_3
c3.于是得到:
再看第2部分是否有小于0的,发现
c
‾
2
,
c
‾
5
\overline c_2,\overline c_5
c2,c5 都小于0,任选一个。我们选择
c
‾
2
\overline c_2
c2所在的列。
计算
θ
⃗
=
[
2
,
+
∞
,
−
1
,
2
]
\vec \theta =[2,+\infty,-1,2]
θ=[2,+∞,−1,2],最小的非负为2,有两个2,那么任选一个。假设选择第2列的第1行为主元,消去主元上下的所有元素。得到:
接下来,再看第二部分是否有小于0的
c
‾
i
\overline c_i
ci ,发现
c
‾
5
\overline c_5
c5小于0。
选择第5列,计算
θ
⃗
=
[
−
2
,
2
,
+
∞
,
0
]
\vec \theta=[-2,2,+\infty,0]
θ=[−2,2,+∞,0],最小的非负数为0,于是选择0所在的位置作为主元。
然而,此时,就有问题了:
θ
=
0
\theta=0
θ=0,这说明 转换后的
X
′
=
X
−
0
λ
⃗
=
X
X'=X-0\vec \lambda=X
X′=X−0λ=X.顶点就不变化了,顶点不再变化,那么最优值当然也就不再变化。此时即使更换了基本列,但是旧的基本列和新的基本列对应的顶点是一样的了。
但是,我们依然可以进行消元,把主元所在的列的非主元位置全部消掉,只不过目标函数值不变而已。
然后,我们选第3列 入基,方法与先前类似:
计算
θ
⃗
\vec \theta
θ 然后选择非负最小的那个所在位置的元作为主元,完了再高斯消元。
OK,此时,第2部分所有元素都是非负了,目标函数收敛。
收敛的最小值为-32,对应的基本列是第2,3,5,6列,这组基对应的顶点为
X
=
[
0
,
1
,
3
,
0
,
2
,
0
,
0
]
X=[0,1,3,0,2,0,0]
X=[0,1,3,0,2,0,0]。
通常,如果每一步 θ = 0 \theta=0 θ=0,那么很容易诱发死循环,因为出现了不同基本列对应于同一个顶点的情况。
单纯形的死循环是指:存在某个状态,它选的基本列与后面某个状态的基本列完全一致,那么就会出现死循环。
处理死循环的方法:
- 给系数矩阵或增广矩阵添加一些微小的扰动,多次计算取平均即可。
- 当第二部分存在多个小于0的候选 c ‾ i \overline c_i ci 时,选择那个最小的列。同时,在同一个列里面,存在多个 θ j \theta_j θj 为0 的时候,选择那个j最小的作为主元。
例2
考虑如下的线性规划:
列出,它们的松弛型:
显然,
X
=
[
0
,
0
,
1
,
1
]
X=[0,0,1,1]
X=[0,0,1,1]是一个可行顶点,因为
b
1
=
b
2
=
1
>
0
b_1=b_2 =1 >0
b1=b2=1>0。
于是列出下面的单纯形表格:
按照上面一样的方法,选择第2部分小于0的,选最小的
c
‾
1
=
−
1
\overline c_1=-1
c1=−1。
计算
θ
⃗
=
[
1
,
−
1
]
\vec \theta=[1,-1]
θ=[1,−1],选择非负最小的那个位置第1行为主元位置,消元后得到。
再看第二部分只剩下第2列可选,计算
θ
=
[
−
1
,
+
∞
]
\theta=[-1,+\infty]
θ=[−1,+∞]。最小的非负
θ
\theta
θ为正无穷。这意味这啥呢?
我们要知道: X ′ = X − θ λ ⃗ = [ 1 , 0 , 0 , 2 ] − θ [ − 1 , 0 , 0 , 0 ] = [ 1 + θ , 0 , 0 , 2 ] X'=X-\theta \vec \lambda=[1,0,0,2]-\theta[-1,0,0,0]=[1+\theta,0,0,2] X′=X−θλ=[1,0,0,2]−θ[−1,0,0,0]=[1+θ,0,0,2],当 θ \theta θ可以取任何非负值,得到的 X ′ X' X′都是可行解,但是它一定不会是顶点。任何得到的 − z = 1 + 2 θ -z=1+2\theta −z=1+2θ 是无上界的,那么 z z z无下界。
因此,当选择的 θ \theta θ 为正无穷,那么目标函数无极小值。 ===>只要存在一个这样的列,那目标函数就没有极小值了。
例3
我们来看下面这个线性规划:
转换为标准型:
此时,不等号的右边有小于0的了,不能一眼看出答案,解辅助线性规划。
先解这个线性规划,注意到它的约束就是松弛型配上一个
x
0
x_0
x0。
看
x
0
x_0
x0是否最小值为0,如果为0,那么把
x
0
x_0
x0除掉,得到的顶点也就是原线性规划的松弛型下的顶点。
列出辅助线性规划的表格:
此时,选择最后两列为基本列,得到的顶点是
X
=
[
0
,
0
,
0
,
−
2
,
−
1
]
X=[0,0,0,-2,-1]
X=[0,0,0,−2,−1] ,不满足条件所以这组基不行。
使用上一篇博客结束的方法,我们做这样的变化。
取所有的 b i b_i bi 的最小值 b j = b 1 = − 2 b_j=b_1=-2 bj=b1=−2,最小值所在行号 l l l 为1。
以x0所在的列里面的第 l l l 个元素为主元,执行一次高斯消元,把x0所在的列里面的所有非主元元素化为0。
c那一行的变化是,选择
x
0
x_0
x0所在的列换入。然后选择这一列的第1行为主元位置消元得到的。
然后再选择第2列换入:
此时,辅助线性规划的最小值为
3
2
\frac{3}{2}
23 ,说明原线性方程可行解为空。其实画出他们的图来也比较直观,没有重叠的。
编程实现
按照上面的表格的方式,实现即可。实现后续下篇博客给出。
线性规划专题——SIMPLEX 单纯形算法(四)——实现