单纯形法是计算“N维线性规划”的工具,在运筹学中占有重要的地位。本文从单纯形法的运行原理层面分析,意在让读者明白单纯形法表为什么要这么用,从而帮助OR初学者走出矩阵语言艰涩难懂的泥潭,为今后的OR学习打下坚实的基础,同时培养出对OR那朦胧而又深邃情感。
“N维线性规划”是单纯形法的理论基础,首先作以介绍。
与其他OR模型一样,线性规划模型由两部分组成:
1.目标函数(需要优化到max 或min)
2.约束函数(包括决策变量和约束条件)
约束函数类似“超女”的海选,满足约束条件的值都可以进入复试,也就是目标函数的检验。目标函数对这些复赛选手进行最优化选择(max或min),留下其中“最爷们”的那个,作为最终方案值。其中由约束函数确定的若干值称为“可行解”,可行解中是目标函数最优化的那个称为“最优解”。
在看到这种模型后,怎样确定他就是线性规划模型呢?依照【美】HamdyA .Taha(国际OR专家)的说法,线性规划模型具有以下三点特征:
1.比例性每个决策变量无论是在目标函数中,还是在约束函数中,其贡献与决策变量值直接成比例。
2.可加性所有决策变量在目标函数和约束函数中的总贡献,等于每个变量各自贡献的直接和。
3.确定性所有目标函数,约束函数中的都是确定的。(本质上,线性规划的各个系数都是概率分布的平均近似值)
满足以上三点特征的模型就是线性规划的模型。
当我们遇到线性规划模型时,该怎么应对呢?
图解法!
图解法的具体步骤如下:
1.确定可行域(可行解对应的点——可行点的集合)。
2.在可行域中的所有可行点中,选出最优点,最优点对应的解就是最优解。
同时给出两条对于图解法具有重大意义的定理,请读者务必自行证明!
定理一:若线性规划问题存在可行解,则可行域是是凸集。
定理二:线性规划问题的最优点在其可行域(凸集)的顶点上。
其中,每一个顶点都有成为最优点的可能,称这些顶点为“基本点(对应基本解)”;而基本点本身还是可行点,最终称之为“基本可行点(对应基本可行解)”。
随着实际问题的复杂化,决策变量和约束条件的数量大大增加,三维图解法很繁琐,三维以上的根本找不到对应的图形(虽然,爱因斯坦说六维空间更适合我们生存)!所以,务必创建新方法来解决多维线性规划问题!
单纯形法就此登场!
前文讲,当维数过多,导致图解法失效,而引入单纯形法解决了此问题。它的优越处在哪里呢?单纯形法的核心思想仍然是是前文所介绍的“线性规划图解法”,它尝试从空间的一个顶点移动到另一个顶点,直到找到最优点为止。但是,之所以单纯形法能解决多维问题,是因为它将图形法转化成了代数法,从而规避掉了多维空间的不可描述性。下面介绍一下这一从图形到代数的完美蜕变过程。这便是单纯形法的运作原理。
第一步:得到模型
为了使单纯形法标准化和简单化,对线性规划模型作两点要求:
1.所有变量非负。
2.除变量的非负限制外,其余约束函数都是等式且具有非负的最右端。
为满足上述两点要求,需要对原函数进行转化。
1.非负转化
a.若右端项为负,约束函数左右两边乘“-1”。
b.若变量未被约束,例如:x没有约束条件,则令x=x``-x`,其中,x``>=0,x`>=0。在实际问题中,x`只能为0,因为如果x`>=0,那么这个过程就做了无用功,降低了效率,这是我们不愿看到的。
2.等式转化
a.在(<=)约束中,右端项可视为对“资源可利用性”的描述,此时,左端项表示对这些有限 资源的使用量。所以,右端项与左端项的差就是“没被使用的资源量”,称为“松弛变量”。因此,左端加上松弛变量(非负)可以使约束函数转化为等式。
b.在(>=)约束中,右端项可视为对“资源达标性”的描述,此时,左端项同样表示对这些资源的使用量。所以,左端项与右端项的差就是“达标后剩余的资源量”,称为“剩余变量”。因此,左端减去剩余变量(非负)可以使约束函数转化为等式。
对约束函数转化后,把添加的变量同样写在目标函数里以作统一,其系数均为零。
对模型标准、简单化后就可以使用单纯形法了。
值得注意的是,加入的松弛变量,剩余变量并不会引起原函数维数的变化。因为这些变量代表的资源没有实际的物质意义,他们存在的意义是通过自己的调节,以达到调节决策变量的效果,因此,可以统一称为“虚拟变量”(笔者起)。
既然不改变维数,它们的几何意义是什么呢?它们是平行于原等式方程确定的直线的直线的集合,取得不同的值就会在不同的位置,变量为0时,直线与原约束方程直线重合。所以,当模型标准化后,所代表的图像发生了根本性的改变,由原来的“阴影覆盖区”变成了由若干直线围成的“栅栏”,而当有变量改变时,“栅栏”的形状也会改变。
第二步:确定一个基本可行解
我们知道,
若线性规划存在N维图像,则可行点有无数个,基本可行点(顶点)有限个。无限个可行点对应约束方程组无限个可行解,有限个基本可行点对应约束方程有限个基本可行解。而在标准约束函数中,变量的个数总是大于方程的个数,就是这个原因使得约束函数有无限组解,这不正好对应着N维图像中那无限多个可行点吗?所以,找到对应顶点的基本可行解的方法就是调整约束函数,让方程的个数等于变量的个数,得到的唯一解就是基本可行解!调整的方法很简单,让变量的个数等于方程的个数,其余的变量都得0(相当于不存在)。所以一个N维线性规划模型的图像,最多会有少个顶点,可以用排列组合算出来,你懂得。
选择的调整方案(哪些变量为0)不同,得到的基本可行解就会不同。从解的角度讲,参与求基本可行解的变量称为“基变量”,未参与求解过程直接0的变量称为“非基变量”;从变量本事角度讲,一些变量被强制为0称为“控零变量”(笔者起),参与求解的变量称为“自由变量”。
第一次得到的基本可行解称为单纯形法的“初始解”。说句实话,笔者很不情愿以具体的数字例子说事,那样不具有普遍性,不够严谨。但是,由于N维空间的不可描述性(怎么总是它),也为了使读者更好的理解,下面举例说明初始解的求法。
1.纯(<=)情况
Maxz=2x1+3x2
s.t{2x1+x2<=4
x1+2x2<=5
x1>=0,x2>=0
作图:图一
标准化后的模型:
Maxz=2x1+3x2+0x3+0x4
s.t{2x1+x2+x3=4
x1+2x2+x4=5
xi>=0,(i=1,2,3,4)
作图:图二
从图中可以很容易地得到一个基本可行解,你懂得,是原点。
推广:由于xi>=0,其他约束全是(<=)形式,原点一定是一个基本可行解!在这种情况下以原点所代表的解作为初始解即可。
2.含(>=)和(=)情况
将上例稍作修改
Maxz=2x1+3x2
s.t{2x1+x2>=4
x1+2x2<=5
x1>=0,x2>=0
作图:图三
标准化后的模型:
Maxz=2x1+3x2+0x3+0x4
s.t{2x1+x2-x3=4
x1+2x2+x4=5
xi>=0,(i=1,2,3,4)
作图:同图三
从图中可以看出,取值范围不包含原点,推广:由变量的非负限制,值直线无法移到原点。这种情况下,原点解不是方程的基本可行解。而对于N维函数来说找到一个解是很麻烦的,或者找到的解是不可行解(不满足所有约束条件)。所以还要使用原点解作为初始解,那么必须是原点包含在可行域里。我们所采用的方法就是OR江湖上的大M法(笔者以为这种叫法并没有抓住转化的本质,有待商榷)。
大M法有两种理解方式。
1.从等式左边理解——直线左移
正如前文所述,直线2x1+x2-x3=4中x1,x2是实变量x3是虚变量,虚变量只起调节作用,那么不妨再次引入一个虚变量x5,使得(x5-x3)成为一个可以为正的虚拟变量,它相当于x4的作用,这样在引入x5后的新图形中,原点就在它的取值范围中了。因为x5是认为引入的,所以称其为“人工变量”,得到的原点初始解成为“人工初始解”,而(x5-x4)成为“人工松弛变量”(笔者起)。其实,x5没有几何意义,它的作用是与x4组成人工松弛变量并调节人工松弛变量的平移度。
2.从等式右边理解——降低标准
如图二,观察x4=0的位置,如果x4=0过原点,那么直线能否左移就不重要了。那么,在2x1+x2-x3=4右侧减去一个变量x5(非负),当x5=4时,直线过原点就得以实现了。从实际意义上讲,右端项是资源使用量达标的标准,减去一个新变量实际上就降低了标准,所以此时可以称新引进的标量为“降标变量”。
含(=)的情况,也是直线不包含原点,原理与(<=)相同的。
无论是“降标变量”还是“人工变量”在实际问题中都是不存在的,所以在找到初始解后,要使问题回归到原始状态。采用的方法是加入严厉“惩罚”,强迫这些变量在最优解中为0,具体做法是,在目标函数中加入人工变量,并将其系数定位“-M”,M是一个充分大的数。
可以看出这种方法的核心是新变量的引入,而不是那个严厉的“惩罚”,所以,称其为“降低标准法”或是“剩余转松弛法”更恰当!
第三步 检验与优化
如果求得的基本可行解是最优解的话,直接导出最优方案。否则,回到第二步求另外的基本可行解。(这三步很像C程序中的循环结构吧,其实单纯形法就是一个运行程序,如果不是为了训练大脑,没什么学习价值,用应用软件就可以解决问题,不过我们还是训练一下大脑吧)
第三步是一个选择和循环的过程,依然用上面的初始例子说明。
Maxz=2x1+3x2
s.t{2x1+x2<=4
x1+2x2<=5
x1>=0,x2>=0
标准化后的模型:
Maxz=2x1+3x2+0x3+0x4
s.t{2x1+x2-x3=4
x1+2x2+x4=5
xi>=0,(i=1,2,3,4)
用第二步的方法得到初始解(0,0,4,5),此时目标函数的值是0,此时的解是否是最优解呢?检验方法是将当前目标函数值与可行域中的其他顶点的目标函数值进行比较,或用检验数的正负性作出判断。比较函数值实在是麻烦,在实际计算中我们是根据检验数的正负性来做取舍的。什么是检验数呢?由目标函数的系数可以看出,增加x1或x2,目标函数都会更优化,为了方便起见,单纯性法规定每次只选一个变量进行优化。那么,就有挑选那个能使目标函数优化程度最大的那个变量。目标函数中的各变量前的系数就描述了这种优化程度,所以称之为“检验数”。随着被选变量的增大,其余变量也会发生变化,产生新的图形,直到新解到达新图形的顶点增加结束。这种增加是渐进的,所以从一个顶点到另一个顶点的路径是沿着图形边缘的,是不能越过一个顶点找到下一个的。考虑到基本可行解只在图形顶点上,我们完全可以忽略中间的增加过程,直接“跳”到下一个顶点去。从基本可行解的确定原理可知,这种“跳跃”在代数中的表现形式就是基变量和非基变量的互换。控零变量不再为零,那么要就有一个自由变量变成控零变量,也就是说控零变量成了基变量,自由变量成了非基变量。这一过程称为“入基”和“出基”。
但是,出基变量时哪个呢?如上例,x1此时仍然为零,原约束可写成:
s.t{x2<=4
2x2<=5
x2>=0
得出0<=x2<=2.5
从几何角度看 4/1和5/2是非基变量关于入基变量轴的截距。入基变量增加到正向最小截距就停止了。
注意:x2的最终范围来自第二条约束。
这样引入虚拟变量后的约束:
x2-x3=4
2x2+x4=5
xi>=0,(i=2,3,4)
在第二条约束中,引入了松弛变量x4,只有x4满足上述x2的可行性,所以x4是那个出基变量。
入基变量和出基变量都找到之后就要进行入基和出基。用出基变量x4表示入基变量x2,得
x2=(5-x4)/2=2.5-x4/2,将其带回原模型,得
Maxz=2x1-3/2x4+0x3+0x4+7.5
s.t{2x1+(5-x4)/2-x3=4
x1+2(5-x4)/2+x4=5
xi>=0,(i=1,2,3,4)
接下来,求新解就可以了。
得到新解进入第三步。
单纯形法的原理到这里就结束了。
计算机程序按照上述过程编译即可。但是纸上作业要求规范和清晰。最后引入单纯形表。
单纯形表(含检验数行)的画法你懂的,笔者不啰嗦了。它的操作原理就是上面所写的,具体步骤你懂的,笔者不啰嗦了。需要强调的是由于采用表上作业,也就引入了矩阵,在求基本可行解时可以简便求解。入基和出基和求新解的过程可由高斯—诺尔当行运算搞定。它认为入基变量所在列为“轴枢列”,出基变量所在行为“轴枢行”,在轴枢行与轴枢列的交叉位置的元素成为“轴枢元素”。运算步骤如下:
1.轴枢行
a.在基列中,以入基变量替换出基变量。
b.新轴枢行=当前轴枢行/轴枢元素
2.所有其他行,包括z行
a.新行=当前行-当前轴枢列的系数*新轴枢列
以上就是笔者对单纯形法的理解,学术问题概无定论,欢迎有缘人指正。