bzoj1061--线性规划

线性规划裸题。。。

根据题目很容易可以得到线性规划方程(以样例为例)：

Min(2*x1+5*x2+2*x3)

x1+ 0+ 0>=2

x1+x2+ 0>=3

0+x2+x3>=4

x1,x2,x3>=0

再将方程对偶，得到：

Max(2*x1+3*x2+4*x3)

x1+x2+ 0<=2

0+x2+x3<=5

0+ 0+x3<=2

x1,x2,x3>=0

这就是线性规划的标准型了。

为了方便单纯型算法，加入变量x4,x5,x6:

Max(2*x1+3*x2+4*x3)

x4+x1+x2+ 0=2

x5+ 0+x2+x3=5

x6+ 0+ 0+x3=2

x1,x2,x3,x4,x5,x6>=0

这就是松弛型。显然此时最优解不变。

将松弛型写成矩阵的形式：

x1 x2 x3

x4  1  1  0  2

x5  0  1  1  5

x6  0  1  1  2

2   3  4  0(k)

当x1,x2,x3取0时，显然满足条件，此时答案为右下角的常数k

我们只需不断增大k，当k达到最大值时最优解就是k了。

那么怎么增大k呢？显然如果我们增大x1，答案会更优。

但x1不能无限制地增大，对于前3个方程，我们得到x1的限制：

1、x1<=2

2、x1无限制

3、x1无限制

我们选择最紧的一个限制1，将x1增大到它，再交换x1,x4。

交换之后再将某些系数改变，使其满足方程就可以了。

于是我们可以不断交换，直到矩阵最后一行的系数都不为正就可以了。最优解就是k。

具体看代码。

 1 #include<iostream>
 2 #include<cstdio>
 3 #include<cstring>
 4 #include<cmath>
 5 using namespace std;
 6 #define N 1001
 7 #define M 10001
 8 #define DB double
 9 #define Eps 1e-7
10 #define INF 0x3f3f3f3f3f3f3f3f
11 DB a[M][N],c[N],b[M],Ans,Tmp;
12 int i,j,n,m,l,r,x;
13 inline void Pivot(int x,int y){                                //转轴操作，使矩阵满足方程
14     b[x]/=a[x][y];
15     for(int i=1;i<=n;i++)if(i!=y)a[x][i]/=a[x][y];
16     a[x][y]=1/a[x][y];
17     for(int i=1;i<=m;i++)
18     if(i!=x&&fabs(a[i][y])>Eps){
19         b[i]-=a[i][y]*b[x];
20         for(int j=1;j<=n;j++)if(j!=y)a[i][j]-=a[i][y]*a[x][j];
21         a[i][y]*=-a[x][y];
22     }
23     Ans+=c[y]*b[x];
24     for(int i=1;i<=n;i++)if(i!=y)c[i]-=c[y]*a[x][i];
25     c[y]*=-a[x][y];
26 }
27 inline DB Simplex(){
28     while(1){                                                   //不断交换
29         for(i=1;i<=n;i++)if(c[i]>Eps)break;
30         if(i>n)return Ans;
31         Tmp=INF;
32         for(j=1;j<=m;j++)
33         if(a[j][i]>Eps&&b[j]/a[j][i]<Tmp)Tmp=b[j]/a[j][i],x=j;
34         if(Tmp==INF)return INF;
35         Pivot(x,i);                                         //交换第x行，第i列
36     }
37 }
38 int main()
39 {
40     scanf("%d%d",&n,&m);
41     for(i=1;i<=n;i++)scanf("%lf",&c[i]);
42     for(i=1;i<=m;i++){
43         scanf("%d%d%lf",&l,&r,&b[i]);
44         for(j=l;j<=r;j++)a[i][j]=1;
45     }
46     printf("%d",(int)(Simplex()+0.5));
47 }

bzoj1061