UOJ #179. 线性规划

#179. 线性规划

http://uoj.ac/problem/179

分析：

　　单纯形算法。

代码：

#include<bits/stdc++.h>
using namespace std;
typedef long long LL;

inline int read() {
    int x=0,f=1;char ch=getchar();for(;!isdigit(ch);ch=getchar())if(ch=='-')f=-1;
    for (;isdigit(ch);ch=getchar())x=x*10+ch-'0';return x*f;
}

const int N = 25;
const double eps = 1e-8;
double a[N][N], ans[N];
int B[N<<1], n, m;
/*
a[0][i] -> ci 目标函数中第i个元素系数
a[i][0] -> bi 第i条约束中小于等于的数 
a[i][j] -> Aij 第i条约束中第j个元素的系数
最大化 sigma(ci*xi),i∈N
约束 xj=bj-sigma(aji*xi) ,j∈B
*/
void Pivot(int l,int e) {
    swap(B[n+l], B[e]); // 交换基本变量（n+l）与非基本变量（e） 
    double t = a[l][e];
    a[l][e] = 1;
    for (int j=0; j<=n; ++j) a[l][j] /= t;
    for (int i=0; i<=m; ++i) 
        if (i !=l && abs(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 Simplex() {
    while (true) {
        int l = 0, e = 0;
        double mn = 1e18;
        for (int j=1; j<=n; ++j) 
            if (a[0][j] > eps) {e = j; break;}
        if (!e) return true;
        for (int i=1; i<=m; ++i) 
            if (a[i][e] > eps && a[i][0]/a[i][e] < mn) 
                mn = a[i][0] / a[i][e], l = i; // 找到约束最紧的 
        if (!l) {
            puts("Unbounded");
            return false;
        }
        Pivot(l,e);
    }
}
/*
初始化
方法一：引入一个辅助线性规划 要求最大化-x0
约束为 xj=bj-sigma(aji*xi)+x0 ,j∈B然后用x0替换bj为负的约束
下面的是方法二：
在bi为负的时候，把所有基变量设为0不是一组合法的初始解
所以选择一个bi为负的基变量x[i+n]
然后在该约束右边找一个系数为正(即原系数为负)的非基变量进行Pivot操作
如果没有系数为正显然就无解了
*/
bool init() {
    while (true) {
        int e = 0, l = 0;
        for (int i=1; i<=m; ++i) 
            if (a[i][0] < -eps && (!l || (rand() & 1))) l = i;
        if (!l) return true;
        for (int j=1; j<=n; ++j) 
            if (a[l][j] < -eps && (!e || (rand() & 1))) e = j;
        if (!e) {
            puts("Infeasible");
            return false;
        }
        Pivot(l,e);
    }
}
int main() {
    n = read(), m = read();
    int Type = read();
    for (int i=1; i<=n; ++i) a[0][i] = read();
    for (int i=1; i<=m; ++i) {
        for (int j=1; j<=n; ++j) a[i][j] = read();
        a[i][0] = read();
    }
    for (int i=1; i<=n; ++i) B[i] = i;
    if (init() && Simplex()) {
        printf("%.8lf
",-a[0][0]);
        if (Type) {
            for (int i=1; i<=m; ++i) ans[B[i+n]] = a[i][0];
            for (int i=1; i<=n; ++i) printf("%.8lf ",ans[i]);
        }
    }
    return 0;
}