POJ 3155 Hard Life（最大密度子图）

裸题。输入一个无向图，输出最大密度子图（输出子图结点数和升序编号）。

看了《最小割模型在信息学竞赛中的应用——胡伯涛》的一部分，感觉01分数规划问题又是个大坑。暂时还看不懂。

参考http://blog.csdn.net/power721/article/details/6781518

构图：

把原图中的无向边转换成两条有向边，容量为1。

设一源点，连接所有点，容量为U（取m）。

设一汇点，所有点连接汇点，容量为 U+2g-dv 。

二分枚举最大密度g，其中dv为v的度。

判断(U*n-MaxFlow)/2.>=0。

最后跳出的L就是最大密度。

拿这个L再重新建图，求最大流。

然后从s出发bfs或者dfs，走残留容量大于0的边，所有能到达的点就是答案。

具体分析过程见胡伯涛的论文《最小割模型在信息学竞赛中的应用》

#include <cstdio>
#include <cstring>
#include <iostream>
#include <cmath>
#include <algorithm>
#include <vector>
#include <string>
#include <set>
#include <queue>
using namespace std;

#define ll long long
#define MP make_pair

#define mxn 110
#define mxe 4500
#define inf 1e9
#define eps 1e-8

vector<int>ans;

struct SAP{
	int dis[mxn],pre[mxn],gap[mxn],arc[mxn];
	double f[mxe],cap[mxe];
	int head[mxn],nxt[mxe],vv[mxe],e;
	void init(){e=0;memset(head,-1,sizeof(head));}
	void add(int u,int v,double c){
		vv[e]=v,cap[e]=c,nxt[e]=head[u],head[u]=e++;
		vv[e]=u,cap[e]=0,nxt[e]=head[v],head[v]=e++;
	}
	void bfs(int s){
		int q[mxn];
		int ht=0,tl=0;
		q[tl++]=s;
		ans.clear();
		bool vis[mxn];
		memset(vis,false,sizeof(vis));
		vis[s]=true;
		while(ht<tl){
			int u = q[ht++];
			for(int i=head[u];i!=-1;i=nxt[i]){
				int v = vv[i];
				if(!vis[v] && cap[i]-f[i]){
					vis[v]=true;
					q[tl++]=v;
					if(v<s) ans.push_back(v);
				}
			}
		}
	}
	double max_flow(int s,int t,int n){
		int q[mxn],j,mindis;
		double ans=0;
		int ht=0,tl=1,u,v;
		double low;
		bool found,vis[mxn];
		memset(dis,0,sizeof(dis));
		memset(gap,0,sizeof(gap));
		memset(vis,0,sizeof(vis));
		memset(arc,-1,sizeof(arc));
		memset(f,0,sizeof(f));
		q[0]=t;vis[t]=true;dis[t]=0;gap[0]=1;
		while(ht<tl){
			u=q[ht++];
			for(int i=head[u];i!=-1;i=nxt[i]){
				v = vv[i];
				if(!vis[v]){
					vis[v]=true;
					dis[v]=dis[u]+1;
					q[tl++]=v;
					gap[dis[v]]++;
					arc[v]=head[v];
				}
			}
		}
		u=s;low=inf;pre[s]=s;
		while(dis[s]<n){
			found=false;
			for(int &i=arc[u];i!=-1;i=nxt[i]){
				if(dis[vv[i]]==dis[u]-1 && cap[i]>f[i]){
					found=true;v=vv[i];
					low=min(low,cap[i]-f[i]);
					pre[v]=u;u=v;
					if(u==t){
						while(u!=s){
							u=pre[u];
							f[arc[u]]+=low;
							f[arc[u]^1]-=low;
						}
						ans+=low;low=inf;
					}
					break;
				}
			}
			if(found) continue;
			mindis=n;
			for(int i=head[u];i!=-1;i=nxt[i]){
				if(mindis>dis[vv[i]] && cap[i]>f[i]){
					mindis=dis[vv[j=i]];
					arc[u]=i;
				}
			}
			if(--gap[dis[u]]==0) return ans;
			dis[u]=mindis+1;
			gap[dis[u]]++;
			u=pre[u];
		}
		return ans;
	}
}sap;
int e[1010][2];
int deg[110];
int n,m;
bool jud(double g){
	sap.init();
	for(int i=0;i<m;++i){
		sap.add(e[i][0],e[i][1],1);
		sap.add(e[i][1],e[i][0],1);
	}
	for(int i=1;i<=n;++i) sap.add(n+1,i,m);
	for(int i=1;i<=n;++i) sap.add(i,n+2,m+2*g-deg[i]);
	double mf = sap.max_flow(n+1,n+2,n+2);
	return (n*m-mf)/2>eps;
}
int main(){
	while(~scanf("%d%d",&n,&m)){
		if(m==0){
			puts("1");
			puts("1");
			continue;
		}
		memset(deg,0,sizeof(deg));
		for(int i=0;i<m;++i){
			scanf("%d%d",&e[i][0],&e[i][1]);
			deg[e[i][0]]++;
			deg[e[i][1]]++;
		}
		double l=0,r=m;
		while(r-l>eps){
			double mid = (l+r)/2;
			if(jud(mid)) l=mid;
			else r=mid;
		}
		//printf("%lf
",l);
		jud(l);
		sap.bfs(n+1);
		sort(ans.begin(),ans.end());
		printf("%d
",ans.size());
		for(int i=0;i<ans.size();++i) printf("%d
",ans[i]);
	}
	return 0;
}