线性规划与单纯形——学习笔记

什么是单纯形？

就是用于解决线性规划问题的一般算法。
时间复杂度比较玄学，不是多项式算法，但是实际表现不错。

前置知识

线性规划

在给定有限资源和竞争约束的情况下，要最大化或最小化某个目标，如果可以把目标描述为某些目标的线性函数，且约束为某些变量的不等式或等式，那我们就可以得到一个线性规划问题,如网络流问题就是特殊的线性规划。

几个经典问题的线性规划表达：

最短路：
最大化:dt
满足约束:dv<=du+w(u,v)
ds=0,di≥0

最大流：
最大化:∑v∈Vflowsv−∑v∈Vflowvs
满足约束:flowuv≤cap(u,v)
∑v∈Vflowvu=∑v∈Vflowuv
fuv≥0
最小费用最大流：
最小化:∑(u,v)cost(u,v)flowuv
满足约束:flowuv≤cap(u,v)
∑v∈Vflowsv−∑v∈Vflowvs=MaxFlow
∑v∈Vflowvu=∑v∈Vflowuv
fuv≥0

线性规划表示方法

标准型：

最大化:∑i=1ncixi

满足约束：∑j=1naij∗xj≤bii=1,2...,m

xi≥0i=1,2,...,n

要把普通的式子转化为标准型怎么搞呢？很简单：
若目标是要最小化，取个负号就变最大化了。
若有变量没有非负约束，拆成xi=x′i−x′′i,x′i,x′′i≥0。
若有等式，就变成大于等于和小于等于两个约束。若有大于等于约束，就取反。
松弛型：
为了用单纯形算法，我们更喜欢将约束表示成等式,

对于 m 个标准型中的约束,我们加 m 个变量，改为

满足约束：xi+n=bi−∑j=1naijxi

xi≥0

其中左边的 x 被称为基本变量，右边为非基本变量。我们记基本变量的集合为 B，非基本变量集合为 N 。这样我们就能用一个元组 (N,B,A,b,c,v) 表示一个松弛形：

最大化:z=v+∑i∈Ncixi

满足约束:xi=bi−∑j∈Naijxj,i∈B

xi≥0,i∈(B∪N)

单纯形算法

主要思想：大概可以理解成”不等式的高斯消元”。单纯形通过不断改变基本变量与非基本变量的集合，不断重写松弛形，直到变成最优解显而易见的形式。
松弛形的基本解: 所有非基本变量xi=0，得到基本变量 xi+n=bi,这样一组解，此时目标函数z=v+∑i∈Ncixi=v 。
主要步骤:
1.先找到一组初始的基本可行解。(下面会解释)
2.找到一个在目标函数中系数为正的非基本变量 e，我们希望使他变大。若找不到说明已经到最优解了，结束算法。
3.在m个约束中找到一个限制e的值最紧的约束 l 。
4.交换 e 与 l中的基本变量 的位置，即重写约束l,使 e 变为基本变量，原本 l 中的基本变量变为非基本变量。然后把其他的约束和目标函数中的e替换掉。这个过程称为一次转动。
5.回到步骤2。

关于第一个步骤，找基本解可行解：
由于一开始的 bi 不一定非负，所以基本解不合法，我们需要进行一些转动使得bi 为正。
如何搞呢？
有一种玄学的简单搞法是这样的：每次瞎jb随机选一个负的bi,在对应的约束中瞎jb随机选一个系数负的x, 然后转动一下，这样这个约束就负负得正了。不断重复上述过程，直到找不到负的bi的为止。
这个搞法不能保证正确性，但是简单一些。
正常的方法应该是建立一个辅助线性规划，不过好像用上面的玄学算法的人更多。
辅助线性规划：
(留坑…)

正确性证明什么的不会。
下面是模板(uoj176)，有很多细节，但理解步骤后都不是问题，各人有各人的写法。

#include<cstdio>#include<cstring>#include<algorithm>#include<cmath>using namespace std;const int maxn=55;const double eps=1e-8;int n,m,Type;double a[maxn][maxn],ans[maxn];int id[maxn*2];void Pivot(int l,int e){    swap(id[l+n],id[e]);    double t=a[l][e]; a[l][e]=1;     for(int i=0;i<=n;i++) a[l][i]/=t;    for(int i=0;i<=m;i++) if(i!=l&&fabs(a[i][e])>eps){        t=a[i][e]; a[i][e]=0;         for(int j=0;j<=n;j++) a[i][j]-=t*a[l][j];     }}bool Init(){    while(1){        int l=0,e=0;        for(int i=1;i<=m;i++) if(a[i][0]<-eps&&(!l||(rand()&1))) l=i;        if(!l) return true;        for(int i=1;i<=n;i++) if(a[l][i]<-eps&&(!e||(rand()&1))) e=i;        if(!e) return printf("Infeasible\n"),false;         Pivot(l,e);    }    return true;}bool Simplex(){    while(1){        int l=0,e=0; double _min=1e+50;        for(int i=1;i<=n;i++) if(a[0][i]>eps){ e=i; break; }        if(!e) break;        for(int i=1;i<=m;i++)          if(-a[i][e]<-eps&&a[i][0]/a[i][e]<_min) _min=a[i][0]/a[i][e], l=i;         if(!l) return printf("Unbounded\n"),false;         Pivot(l,e);    }    return true;}int main(){    freopen("uoj179.in","r",stdin);    freopen("uoj179.out","w",stdout);    scanf("%d%d%d",&n,&m,&Type);    for(int i=1;i<=n;i++) scanf("%lf",&a[0][i]);    for(int i=1;i<=m;i++){        for(int j=1;j<=n;j++) scanf("%lf",&a[i][j]);        scanf("%lf",&a[i][0]);    }    for(int i=1;i<=n;i++) id[i]=i;     if(Init()&&Simplex()){        printf("%.8lf\n",-a[0][0]);         if(Type){            for(int i=1;i<=m;i++) ans[id[i+n]]=a[i][0];             for(int i=1;i<=n;i++) printf("%.8lf ",ans[i]);        }    }    return 0;}