洛谷 P6938 - [ICPC2017 WF]Son of Pipe Stream（网络流）

见过的最怪的网络流题，没有之一。

首先新建超级源点，向 $1,2$ 各连 $\infty$ 的边。设最大流为 $A$，那么显然最优方案中 flutter 和 water 流量之和为 $A$。

先分析一波答案函数。显然，最终答案关于 flutter 的流量 $x$ 的函数 $f(x)=x^a(A-x)^{1-a}$。求导得 $f'(x)=ax^{a-1}(A-x)^{1-a}-x^a(1-a)(A-x)^{-a}$，化简得 $f'(x)=x^{a-1}(A-x)^{-a}(aA-ax-x+ax)$，解得 $f(x)$ 在 $(0,aA)$ 上单调递增，$(aA,A)$ 上单调递减。当然，有时候 $x$ 并不能恰好取到 $aA$ —— 显然 flutter 流量上界有最大值 $F_{max}$，water 流量上界也有最大值 $W_{max}$，求出 $F_{max},W_{max}$ 之后找 $[A-W_{max},F_{max}]$ 中距离 $aA$ 最近的点的位置即可，下文中记这个位置为 $F'$。（当然这一步也可以三分，因为代价函数存在严格凸性）

然后很自然的想法是直接以 $1$ 为源点 $3$ 为汇点跑流量限制为 $F'$ 的最大流，但是这样并不能保证 F 和 W 流量方向相同。考虑调整：先在上面的残量网络中找出每条边的方向，然后只保留这个方向的边，双向改单向，流量为残量网络上这条边已经用掉的流量，然后再跑限流 $F'$ 的最大流，这样每条边的流量就是 F 的流量，剩余容量就是 W 的流量，感性理解一下保留最大流中的边不影响 $1$ 到 $3$ 的最大流，因此我们的构造是正确的。

const int MAXN=300;
const int MAXM=5e4;
const double INF=1e9;
int n,m;double v,a;
struct edge{int u,v;double w;}e[MAXM+5];
int S,T,hd[MAXN+5],to[MAXM*2+5],nxt[MAXM*2+5],ec=1;double cap[MAXM*2+5];
void clear(){memset(hd,0,sizeof(hd));ec=1;}
void adde(int u,int v,double f){
	to[++ec]=v;cap[ec]=f;nxt[ec]=hd[u];hd[u]=ec;
	to[++ec]=u;cap[ec]=f;nxt[ec]=hd[v];hd[v]=ec;
}
int dep[MAXN+5],now[MAXN+5];
bool getdep(){
	memset(dep,-1,sizeof(dep));queue<int>q;q.push(S);dep[S]=0;
	while(!q.empty()){
		int x=q.front();q.pop();now[x]=hd[x];
		for(int e=hd[x];e;e=nxt[e]){
			int y=to[e];double z=cap[e];
			if(!~dep[y]&&z>0)dep[y]=dep[x]+1,q.push(y);
		}
	}return ~dep[T];
}
double getflow(int x,double f){
	if(x==T)return f;double ret=0;
	for(int &e=now[x];e;e=nxt[e]){
		int y=to[e];double z=cap[e];
		if(dep[y]==dep[x]+1&&z>0){
			double w=getflow(y,min(f-ret,z));
			ret+=w;cap[e]-=w;cap[e^1]+=w;
			if(f==ret)return ret;
		}
	}return ret;
}
double dinic(){
	double ret=0;
	while(getdep())ret+=getflow(S,INF);
	return ret;
}
void readd_graph(){clear();for(int i=1;i<=m;i++)adde(e[i].u,e[i].v,e[i].w);}
double Fmax,Wmax,Z,ndF,ndW;int dir[MAXM+5];
int main(){
	freopen("transport.in","r",stdin);
	freopen("transport.out","w",stdout);
	scanf("%d%d%lf%lf",&n,&m,&v,&a);
	for(int i=1;i<=m;i++)scanf("%d%d%lf",&e[i].u,&e[i].v,&e[i].w);
	readd_graph();S=1;T=3;Fmax=dinic();
	readd_graph();S=2;T=3;Wmax=dinic();
	readd_graph();S=n+1;T=3;adde(S,1,INF);adde(S,2,INF);Z=dinic();
	if(Z-Wmax<=a*Z&&a*Z<=Fmax)ndF=a*Z;
	else if(Fmax<a*Z)ndF=Fmax;
	else ndF=Z-Wmax;
	ndW=Z-ndF;
	readd_graph();S=n+1;T=3;adde(S,1,ndF);adde(S,2,ndW);dinic();
	for(int i=2;i<=ec;i+=2){
		double mid=(cap[i]+cap[i^1])/2;
		if(cap[i]>mid)dir[i>>1]=-1,cap[i^1]=mid-cap[i^1],cap[i]=0;
		else dir[i>>1]=1,cap[i]=mid-cap[i],cap[i^1]=0;
	}
	for(int e=hd[n+1];e;e=nxt[e]){
		int y=to[e];
		if(y==1)cap[e]=ndF,cap[e^1]=0;
		else cap[e]=cap[e^1]=0;
	}
	dinic();
	for(int i=2;i<=ec-4;i+=2){
		if(dir[i>>1]==1)printf("%.10lf %.10lf\n",cap[i^1]/v,cap[i]);
		else printf("-%.10lf -%.10lf\n",cap[i]/v,cap[i^1]);
	}
	double res=pow(ndF/v,a)*pow(ndW,1-a);
	printf("%.10lf\n",res);
	return 0;
}