bzoj 3993 星际战争 - 二分答案 - 最大流

3333年，在银河系的某星球上，X军团和Y军团正在激烈地作战。在战斗的某一阶段，Y军团一共派遣了N个巨型机器人进攻X军团的阵地，其中第i个巨型机器人的装甲值为Ai。当一个巨型机器人的装甲值减少到0或者以下时，这个巨型机器人就被摧毁了。X军团有M个激光武器，其中第i个激光武器每秒可以削减一个巨型机器人Bi的装甲值。激光武器的攻击是连续的。这种激光武器非常奇怪，一个激光武器只能攻击一些特定的敌人。Y军团看到自己的巨型机器人被X军团一个一个消灭，他们急需下达更多的指令。为了这个目标，Y军团需要知道X军团最少需要用多长时间才能将Y军团的所有巨型机器人摧毁。但是他们不会计算这个问题，因此向你求助。

Input

第一行，两个整数，N、M。

第二行，N个整数，A1、A2…AN。

第三行，M个整数，B1、B2…BM。

接下来的M行，每行N个整数，这些整数均为0或者1。这部分中的第i行的第j个整数为0表示第i个激光武器不可以攻击第j个巨型机器人，为1表示第i个激光武器可以攻击第j个巨型机器人。

Output

 一行，一个实数，表示X军团要摧毁Y军团的所有巨型机器人最少需要的时间。输出结果与标准答案的绝对误差不超过10-3即视为正确。

Sample Input

2 2
3 10
4 6
0 1
1 1

Sample Output

1.300000

Hint

 【样例说明1】

战斗开始后的前0.5秒，激光武器1攻击2号巨型机器人，激光武器2攻击1号巨型机器人。1号巨型机器人被完全摧毁，2号巨型机器人还剩余8的装甲值；

接下来的0.8秒，激光武器1、2同时攻击2号巨型机器人。2号巨型机器人被完全摧毁。

对于全部的数据，1<=N, M<=50，1<=Ai<=105，1<=Bi<=1000，输入数据保证X军团一定能摧毁Y军团的所有巨型机器人

---

　　二分时间，然后建图，激光武器和源点相连，容量为这个激光武器在这个时间内能够造成的伤害，机器人和汇点连边，容量为机器人的装甲值，激光武器和它能够攻击的目标连一条边，容量为无限大。

Code

/**
 * bzoj
 * Problem#3993
 * Accepted
 * Time:48ms
 * Memory:1688k
 */
#include <iostream>
#include <cstdio>
#include <ctime>
#include <cmath>
#include <cctype>
#include <cstring>
#include <cstdlib>
#include <fstream>
#include <sstream>
#include <algorithm>
#include <map>
#include <set>
#include <stack>
#include <queue>
#include <vector>
#include <stack>
#ifndef WIN32
#define Auto "%lld"
#else
#define Auto "%I64d"
#endif
using namespace std;
typedef bool boolean;
const signed int inf = (signed)((1u << 31) - 1);
const double eps = 1e-6;
#define smin(a, b) a = min(a, b)
#define smax(a, b) a = max(a, b)
#define max3(a, b, c) max(a, max(b, c))
#define min3(a, b, c) min(a, min(b, c))
template<typename T>
inline boolean readInteger(T& u){
    char x;
    int aFlag = 1;
    while(!isdigit((x = getchar())) && x != '-' && x != -1);
    if(x == -1) {
        ungetc(x, stdin);
        return false;
    }
    if(x == '-'){
        x = getchar();
        aFlag = -1;
    }
    for(u = x - '0'; isdigit((x = getchar())); u = (u << 1) + (u << 3) + x - '0');
    ungetc(x, stdin);
    u *= aFlag;
    return true;
}

typedef class Edge {
    public:
        int end;
        int next;
        double flow;
        double cap;
        Edge(int end = 0, int next = -1, double flow = 0, double cap = 0):end(end), next(next), flow(flow), cap(cap) {    }
}Edge;

typedef class MapManager {
    public:
        int ce;
        vector<Edge> edge;
        int* h;
        
        MapManager():ce(0), h(NULL) {        }
        MapManager(int nodes):ce(0) {
            h = new int[(const int)(nodes + 1)];
            memset(h, -1, sizeof(int) * (nodes + 1));
        }
        
        inline void addEdge(int from, int end, double flow, double cap) {
            edge.push_back(Edge(end, h[from], flow, cap));
            h[from] = ce++;
        }
        
        inline void addDoubleEdge(int from, int end, double cap) {
            if(cap == 0)    return;
            addEdge(from, end, 0, cap);
            addEdge(end, from, cap, cap);
        }
        
        Edge& operator [] (int pos) {
            return edge[pos];
        }
        
        inline void clear() {
            delete[] h;
            edge.clear();
        }
}MapManager;
#define m_begin(g, i) (g).h[(i)]
#define m_endpos -1

inline boolean dcmp(double a, double b) {
    return fabs(a - b) < eps;
}

template<typename T>class Matrix{
    public:
        T *p;
        int lines;
        int rows;
        Matrix():p(NULL){    }
        Matrix(int rows, int lines):lines(lines), rows(rows){
            p = new T[(lines * rows)];
        }
        T* operator [](int pos){
            return (p + pos * lines);
        }
};
#define matset(m, i, s) memset((m).p, (i), (s) * (m).lines * (m).rows)

int n, m;
int *A, *B;
int sA = 0;
int s, t;
Matrix<boolean> atable;

inline void init() {
    readInteger(n);
    readInteger(m);
    A = new int[(const int)(n + 1)];
    B = new int[(const int)(m + 1)];
    atable = Matrix<boolean>(m + 1, n + 1);
    for(int i = 1; i <= n; i++)
        readInteger(A[i]), sA += A[i];
    for(int i = 1; i <= m; i++)
        readInteger(B[i]);
    for(int i = 1; i <= m; i++)
        for(int j = 1; j <= n; j++)
            readInteger(atable[i][j]);
    s = 0, t = n + m + 1;
}

MapManager g;
inline void mkmap(double mid) {
    g = MapManager(n + m + 1);
    for(int i = 1; i <= m; i++)
        g.addDoubleEdge(s, i, B[i] * mid); 
    for(int i = 1; i <= n; i++)
        g.addDoubleEdge(i + m, t, A[i]); 
    for(int i = 1; i <= m; i++)
        for(int j = 1; j <= n; j++)
            if(atable[i][j])
                g.addDoubleEdge(i, j + m, inf); 
}

int* dis;
boolean* vis;
queue<int> que;
inline boolean bfs() {
    memset(vis, false, sizeof(boolean) * (t + 1));
    que.push(s);
    vis[s] = true;
    dis[s] = 0;
    while(!que.empty()) {
        int e = que.front();
        que.pop();
        for(int i = m_begin(g, e); i != m_endpos; i = g[i].next) {
            if(dcmp(g[i].cap, g[i].flow))    continue;
            int eu = g[i].end;
            if(vis[eu])    continue;
            vis[eu] = true;
            dis[eu] = dis[e] + 1;
            que.push(eu);
        }
    }
    return vis[t];
}

int *cur;
inline double blockedflow(int node, double minf) {
    if((node == t) || (minf < eps))    return minf;
    double f, flow = 0;
    for(int& i = cur[node]; i != m_endpos; i = g[i].next) {
        int& eu = g[i].end;
        if(dis[eu] == (dis[node] + 1) && (f = blockedflow(eu, min(minf, g[i].cap - g[i].flow))) >= eps) {
            minf -= f;
            flow += f;
            g[i].flow += f;
            g[i ^ 1].flow -= f;
            if(minf < eps)    return flow;
        }
    }
    return flow;
}

inline void init_dinic() {
    vis = new boolean[(const int)(t + 1)];
    dis = new int[(const int)(t + 1)];
    cur = new int[(const int)(t + 1)];
}

inline boolean dinic(double mid) {
    mkmap(mid);
    double maxflow = 0.0;
    while(bfs()) {
        for(int i = s; i <= t; i++)
            cur[i] = m_begin(g, i);
        maxflow += blockedflow(s, inf);
    }
    g.clear();
    return dcmp(maxflow, sA);
} 

inline void solve() {
    init_dinic();
    double l = 0.0, r = sA;
    while(l + eps <= r) {
        double mid = (l + r) / 2;
        if(dinic(mid))    r = mid;
        else l = mid;
    }
    printf("%.6lf", r);
}

int main() {
    init();
    solve();
    return 0;
}