UOJ #179. 线性规划

#179. 线性规划

分析：
　　单纯形算法。
代码：
 1 #include<bits/stdc++.h>
 2 using namespace std;
 3 typedef long long LL;
 4 
 5 inline int read() {
 6     int x=0,f=1;char ch=getchar();for(;!isdigit(ch);ch=getchar())if(ch=='-')f=-1;
 7     for (;isdigit(ch);ch=getchar())x=x*10+ch-'0';return x*f;
 8 }
 9 
10 const int N = 25;
11 const double eps = 1e-8;
12 double a[N][N], ans[N];
13 int B[N<<1], n, m;
14 /*
15 a[0][i] -> ci 目标函数中第i个元素系数
16 a[i][0] -> bi 第i条约束中小于等于的数 
17 a[i][j] -> Aij 第i条约束中第j个元素的系数
18 最大化 sigma(ci*xi),i∈N
19 约束 xj=bj-sigma(aji*xi) ,j∈B
20 */
21 void Pivot(int l,int e) {
22     swap(B[n+l], B[e]); // 交换基本变量（n+l）与非基本变量（e） 
23     double t = a[l][e];
24     a[l][e] = 1;
25     for (int j=0; j<=n; ++j) a[l][j] /= t;
26     for (int i=0; i<=m; ++i) 
27         if (i !=l && abs(a[i][e]) > eps) {
28             t = a[i][e];
29             a[i][e] = 0;
30             for (int j=0; j<=n; ++j) a[i][j] -= t * a[l][j];
31         }
32 }
33 bool Simplex() {
34     while (true) {
35         int l = 0, e = 0;
36         double mn = 1e18;
37         for (int j=1; j<=n; ++j) 
38             if (a[0][j] > eps) {e = j; break;}
39         if (!e) return true;
40         for (int i=1; i<=m; ++i) 
41             if (a[i][e] > eps && a[i][0]/a[i][e] < mn) 
42                 mn = a[i][0] / a[i][e], l = i; // 找到约束最紧的 
43         if (!l) {
44             puts("Unbounded");
45             return false;
46         }
47         Pivot(l,e);
48     }
49 }
50 /*
51 初始化
52 方法一：引入一个辅助线性规划 要求最大化-x0
53 约束为 xj=bj-sigma(aji*xi)+x0 ,j∈B然后用x0替换bj为负的约束
54 下面的是方法二：
55 在bi为负的时候，把所有基变量设为0不是一组合法的初始解
56 所以选择一个bi为负的基变量x[i+n]
57 然后在该约束右边找一个系数为正(即原系数为负)的非基变量进行Pivot操作
58 如果没有系数为正显然就无解了
59 */
60 bool init() {
61     while (true) {
62         int e = 0, l = 0;
63         for (int i=1; i<=m; ++i) 
64             if (a[i][0] < -eps && (!l || (rand() & 1))) l = i;
65         if (!l) return true;
66         for (int j=1; j<=n; ++j) 
67             if (a[l][j] < -eps && (!e || (rand() & 1))) e = j;
68         if (!e) {
69             puts("Infeasible");
70             return false;
71         }
72         Pivot(l,e);
73     }
74 }
75 int main() {
76     n = read(), m = read();
77     int Type = read();
78     for (int i=1; i<=n; ++i) a[0][i] = read();
79     for (int i=1; i<=m; ++i) {
80         for (int j=1; j<=n; ++j) a[i][j] = read();
81         a[i][0] = read();
82     }
83     for (int i=1; i<=n; ++i) B[i] = i;
84     if (init() && Simplex()) {
85         printf("%.8lf
",-a[0][0]);
86         if (Type) {
87             for (int i=1; i<=m; ++i) ans[B[i+n]] = a[i][0];
88             for (int i=1; i<=n; ++i) printf("%.8lf ",ans[i]);
89         }
90     }
91     return 0;
92 }