御坂网络(枚举基准,二分图)

御坂网络(枚举基准,二分图)

现在有n个A点,m个B点(n,m<=200),和最大半径rmax。可以选定一个半径r<rmax,以r为半径,A点或B点为圆心,同时保证A点构造出的圆和B点构造出的圆不相交。问最大面积和。

首先ckr大佬考场上写三分,其实是不对的,但结果把数据卡过去了qwq,orz。

那正解是什么呢?首先有一个显然的性质:可能作为答案的半径一定会使得一些类别不同的圆相切,不然还可以扩张半径直到相切为止,并且圆的个数仍然不变。也有可能半径就是rmax。

所以首先要枚举半径。正常人可能都会认为,半径作自变量,最大面积和作因变量,这个函数是一个凸函数,但其实不是的,因此不能三分。枚举完半径,剩下的工作就是找到最大独立集,这个可以通过网络流跑二分图实现。

枚举半径是\(O(n^2)\)的,dinic跑二分图是\(的O(n^2m)的\),那岂不是\(O(n^6)\)了吗?具体分析一下,如果将半径从小到大排序,那么就只有加边操作,也就是说网络流不用完全重新跑,只要更新一下就行了。每次枚举半径都要至少bfs一次,因此bfs的总时间复杂度是\(O(n^4)\)。由于最多会增广n次,并且每bfs一遍就至少dfs一遍,因此dfs的总时间复杂度也是\(O(n^4)\)

但是\(O(n^4)\)也过不了啊。观察一波方法,我们可以发现时间复杂度主要用在了判断源点和汇点的连通性上。如果我们能提前判断好S和T是否联通,那么最多只会进行n次bfs,时间复杂度就被降到\(O(n^3)\)了。因此我们不用bfs判断连通性而是用bitset传递闭包来做。如果当前只是加边,那么只需要更新与当前点相连点的传递闭包即可,时间复杂度是\(O(n^2)*O(n^2/32)=O(n^4/32)\)。如果在bfs更新了网络流,那就必须跑一遍全新的传递闭包,时间复杂度依然是\(O(n)*O(n^3/32)=O(n^4/32)\)(O记号这样用其实是不严谨的,不过理解万岁)。因此我们就达到了\(O(n^4/32)\)的复杂度。

做完这道题我的收获很大啊。虽然用了整整一天,但是我学到了二分图网络流+bitset传递闭包,就是后面半天耗费了太多时间在一个感叹号上qwq。

#include <bitset>
#include <cstdio>
#include <cstring>
#include <algorithm>
using namespace std;

const int maxn=205*2, INF=1e9;
const double PI=3.1415926535898;
double sqr(double x){ return x*x; }

struct Edge{
    int to, next, cap;
}e[maxn*maxn];
bitset<maxn> f[maxn];  //连通性
int fir[maxn], cnte, src, dst;
void addedge(int x, int y, int v){
    Edge &ed=e[++cnte];
    ed.to=y; ed.next=fir[x]; ed.cap=v; fir[x]=cnte;
    if (!v) return;
    f[x]|=f[y];
    for (int i=src; i<=dst; ++i)
        if (f[i][x]) f[i]|=f[x];
}

int dep[maxn], q[maxn], cur[maxn];

int bfs(){
    int head=0, tail=0;
    memset(dep, 0, sizeof(dep)); memset(cur, 0, sizeof(cur));
    dep[src]=1; q[tail++]=src; int tmp;
    while (head<tail){
        tmp=q[head++];
        for (int i=fir[tmp]; i; i=e[i].next)
            if (e[i].cap>0&&!dep[e[i].to]){
                dep[e[i].to]=dep[tmp]+1;
                q[tail++]=e[i].to;
            }
    }
    return dep[dst]?1:0;
}

int dfs(int now, int flow){
    if (now==dst) return flow;
    if (cur[now]==-1) return 0;
    for (int i=(cur[now]?cur[now]:fir[now]); i; i=e[i].next){
        cur[now]=i;
        if (dep[e[i].to]==dep[now]+1&&e[i].cap){
            int minm=dfs(e[i].to, min(flow, e[i].cap));
            if (minm>0){
                e[i].cap-=minm; e[(i-1^1)+1].cap+=minm;
                return minm;
            }
        }
    } cur[now]=-1;
    return 0;
}

int n, m, x[maxn], y[maxn], totc, maxr;
struct Couple{
    int x, y, d;
}c[maxn*maxn];  //表示所有可能的半径(存的是直径的平方),同时可以用用来加边
bool cmp(const Couple &a, const Couple &b){ return a.d<b.d; }

double ans;

inline void update(){
    for (int i=src; i<=dst; ++i){
        f[i].reset();
        for (int j=fir[i]; j; j=e[j].next)
            if (e[j].cap) f[i][e[j].to]=1;
    }
    for (int i=src; i<=dst; ++i)
        for (int j=src; j<=dst; ++j)
            if (f[i][j]) f[i]|=f[j];
}

int main(){
    scanf("%d%d%d", &n, &m, &maxr); src=0, dst=n+m+1;
    for (int i=src; i<=dst; ++i) f[i][i]=1;
    for (int i=1; i<=n; ++i){ addedge(0, i, 1); addedge(i, 0, 0); }
    for (int i=n+1; i<=n+m; ++i){ addedge(i, dst, 1); addedge(dst, i, 0); }
    for (int i=1; i<=n; ++i) scanf("%d", &x[i]);
    for (int i=1; i<=n; ++i) scanf("%d", &y[i]);
    for (int i=n+1; i<=n+m; ++i) scanf("%d", &x[i]);
    for (int i=n+1; i<=n+m; ++i) scanf("%d", &y[i]);
    for (int i=1; i<=n; ++i)
        for (int j=n+1; j<=n+m; ++j){
            c[++totc].x=i; c[totc].y=j;
            c[totc].d=sqr(x[j]-x[i])+sqr(y[j]-y[i]);
            if (c[totc].d>=maxr*maxr*4){ --totc; continue; }
        }
    sort(c+1, c+totc+1, cmp); int t, tot=n+m;
    for (int i=1; i<=totc; ++i){
        if (i==1||c[i].d!=c[i-1].d){  //不要统计多种情况
            if (f[src][dst]){ //如果源汇不联通就说明找无可找,不用更新最大独立集大小
                while (bfs()) while (t=dfs(src, INF)) tot-=t;
                update(); }
            ans=max(ans, tot*PI*c[i].d/4);  //记录当前枚举的半径的圆面积*最大独立集中点的个数
        }
        addedge(c[i].x, c[i].y, 1);  //后面半径继续扩大,这个半径就是非法的了。
        addedge(c[i].y, c[i].x, 0);
    }
    while (bfs()) while (t=dfs(src, INF)) tot-=t;  //半径扩张到rmax的话,又会有之前的半径非法
    printf("%.4lf\n", max(ans, tot*maxr*maxr*PI));  //懒得特判rmax=最后半径的情况了
    return 0;
}
posted @ 2018-05-22 10:25  pechpo  阅读(390)  评论(0编辑  收藏  举报