重修 MCMF 最小费用最大流

Primal-Dual 原始对偶算法

Upd 2023年10月24日

我们尝试将 f 次 spfa 改成 dijkstra 来优化复杂度。

但是边权不一定是非负的，考虑搬用全源最短路（Johnson）的做法（这里能找到），给每个点设一个势能 hi（hS=0）使得每条边加上势能差后都非负。

但是网络流的图是会变的，会有反向边出现，我们只需要每次增广之后将所有点势能 hi 加上 disi，其中 disi 为新边权跑出来的最短路。

证明略，可看 OI-wiki

Code

由于初始求势能要一次 spfa，复杂度 O(nm+mflog)。

稠密图可以用 O(n2) dijkstra，复杂度 O(nm+n2f)，具体应用有 二分图最大权匹配。

用途

求在满足最大流的情况下最小的费用。

算法思路

每次找到一条费用之和最小的增广路径增广（反向边的费用是正向边的相反数）。

易证明贪心正确，或参考 OI Wiki 的证明。

具体实现

每次在残量网络中跑 Bellman-Ford 或 SPFA（按费用跑最短路），取路径上流量最小值增广，修改答案并改动网络上的这条路径。

注意我们不像 Dinic 一样按树形增广，而是像 EK 一样单条路径。树形也行，不过时间一样，何必呢。

时间复杂度

O(nmf) 点数、边数、最大流的积。

最坏可以到达 O(n32n/2)，其实不刻意卡的话是伪多项式复杂度的。

由于 SPFA 已经玄学了，MCMF 时间复杂度更是玄学。

例题 & Code

P3381 【模板】最小费用最大流

#include<bits/stdc++.h>
using namespace std;
#define IOS ios::sync_with_stdio(0);cin.tie(0);cout.tie(0)
#define For(i,j,k) for(int i=j;i<=k;i++)
#define Rof(i,j,k) for(int i=j;i>=k;i--)
#define N 5010
#define M 50010
#define inf 0x7fffffff
struct edge{int to,flow,cost,nxt;}e[2*M];//2 times of edges!!!
int head[N],tot=1;
int n,m,S,T;
inline void adde(int x,int y,int flow,int cost){
	e[++tot]=(edge){y,flow, cost,head[x]}; head[x]=tot;
	e[++tot]=(edge){x,0   ,-cost,head[y]}; head[y]=tot;
}
int dis[N];//min toll to i (which SPFA works)
int flow[N];//the flow of the path with the min toll from S to i (just ONE path)
int pre[N],lst[N];//the last node and edge of the path with the min toll
queue<int> q;//for SPFA
bool ins[N];//in the queue or not to reduce the useless 'q.push()'
bool spfa(){
	while(!q.empty()) q.pop();
	For(i,1,n) ins[i]=pre[i]=0;
	For(i,1,n) dis[i]=inf;
	dis[S]=0;
	ins[S]=1;
	flow[S]=inf;
	q.push(S);
	int x;
	while(!q.empty()){
		x=q.front();
		q.pop();
		ins[x]=0;
		for(int i=head[x];i;i=e[i].nxt){
			int to=e[i].to;
			assert(e[i].flow>=0);
			if(e[i].flow==0 || dis[to]<=dis[x]+e[i].cost) continue;
			dis[to]=dis[x]+e[i].cost;
			pre[to]=x;
			lst[to]=i;
			flow[to]=min(flow[x],e[i].flow);
			if(!ins[to]){
				ins[to]=1;
				q.push(to);
			}
		}
	}
	return dis[T]<inf;
}
void mcmf(){
	int mf=0;
	int mc=0;
	while(spfa()){//find ONE augmentation path at a time to augment
		mf+=flow[T];
		mc+=flow[T]*dis[T];
		int x=T;
		while(x!=S){
			e[lst[x]  ].flow-=flow[T];
			e[lst[x]^1].flow+=flow[T];
			x=pre[x];
		}
	}
	cout<<mf<<" "<<mc<<endl;
}
signed main(){IOS;
	cin>>n>>m>>S>>T;
	int x,y,w,c;
	while(m--){
		cin>>x>>y>>w>>c;
		adde(x,y,w,c);
	}
	mcmf();
return 0;}

去注释版：

#include<bits/stdc++.h>
using namespace std;
#define IOS ios::sync_with_stdio(0);cin.tie(0);cout.tie(0)
#define For(i,j,k) for(int i=j;i<=k;i++)
#define Rof(i,j,k) for(int i=j;i>=k;i--)
#define N 5010
#define M 50010
#define inf 0x7fffffff
struct edge{int to,flow,cost,nxt;}e[2*M];
int head[N],tot=1;
int n,m,S,T;
inline void adde(int x,int y,int flow,int cost){
	e[++tot]=(edge){y,flow, cost,head[x]}; head[x]=tot;
	e[++tot]=(edge){x,0   ,-cost,head[y]}; head[y]=tot;
}
int dis[N];
int flow[N];
int pre[N],lst[N];
queue<int> q;
bool ins[N];
bool spfa(){
	while(!q.empty()) q.pop();
	For(i,1,n) ins[i]=pre[i]=0;
	For(i,1,n) dis[i]=inf;
	dis[S]=0;
	ins[S]=1;
	flow[S]=inf;
	q.push(S);
	int x;
	while(!q.empty()){
		x=q.front();
		q.pop();
		ins[x]=0;
		for(int i=head[x];i;i=e[i].nxt){
			int to=e[i].to;
			assert(e[i].flow>=0);
			if(e[i].flow==0 || dis[to]<=dis[x]+e[i].cost) continue;
			dis[to]=dis[x]+e[i].cost;
			pre[to]=x;
			lst[to]=i;
			flow[to]=min(flow[x],e[i].flow);
			if(!ins[to]){
				ins[to]=1;
				q.push(to);
			}
		}
	}
	return dis[T]<inf;
}
void mcmf(){
	int mf=0;
	int mc=0;
	while(spfa()){
		mf+=flow[T];
		mc+=flow[T]*dis[T];
		int x=T;
		while(x!=S){
			e[lst[x]  ].flow-=flow[T];
			e[lst[x]^1].flow+=flow[T];
			x=pre[x];
		}
	}
	cout<<mf<<" "<<mc<<endl;
}
signed main(){IOS;
	cin>>n>>m>>S>>T;
	int x,y,w,c;
	while(m--){
		cin>>x>>y>>w>>c;
		adde(x,y,w,c);
	}
	mcmf();
return 0;}

AcWing382. K取方格数

算是 传纸条 的加强版吧。

#include<bits/stdc++.h>
using namespace std;
#define IOS ios::sync_with_stdio(0);cin.tie(0);cout.tie(0)
#define fir first
#define sec second
#define mkp make_pair
#define pb emplace_back
#define mem(x,y) memset(x,y,sizeof(x))
#define For(i,j,k) for(int i=j;i<=k;i++)
#define Rof(i,j,k) for(int i=j;i>=k;i--)
#define ckmx(a,b) a=max(a,b)
#define ckmn(a,b) a=min(a,b)
#define fin(s) freopen(s,"r",stdin)
#define fout(s) freopen(s,"w",stdout)
#define file(s) fin(s".in");fout(s".out")
#define cerr cerr<<'_'
#define debug cerr<<"Passed line #"<<__LINE__<<endl
template<typename T>T ov(T x){cerr<<"Value: "<<x<<endl;return x;}
#define N 52
#define inf 0x3f3f3f3f
//Since the algorithm template is the minimum cost flow, 
//and the topic requires the maximum value, 
//the edge weights take the opposite number. 
struct edge{int to,flow,cost,nxt;}e[8*N*N];
int head[2*N*N],tot=1;
int n,k,S,T,lim;
inline int num(int x,int y,bool z){return (((x-1)*n+y)<<1)+z;}
inline void adde(int x,int y,int flow,int cost){
	e[++tot]=(edge){y,flow, cost,head[x]}; head[x]=tot;
	e[++tot]=(edge){x,0   ,-cost,head[y]}; head[y]=tot;
}
int dis[2*N*N];//min toll to i (which SPFA works)
int flow[2*N*N];//the flow of the path with the min toll from S to i (just ONE path)
int pre[2*N*N],lst[2*N*N];//the last node and edge of the path with the min toll
queue<int> q;//for SPFA
bool ins[2*N*N];//in the queue or not to reduce the useless 'q.push()'
bool spfa(){
	while(!q.empty()) q.pop();
	For(i,1,lim) ins[i]=pre[i]=0;
	For(i,1,lim) dis[i]=inf;
	dis[S]=0;
	ins[S]=1;
	flow[S]=inf;
	q.push(S);
	int x;
	while(!q.empty()){
		x=q.front();
		q.pop();
		ins[x]=0;
		for(int i=head[x];i;i=e[i].nxt){
			int to=e[i].to;
			assert(e[i].flow>=0);
			if(e[i].flow==0 || dis[to]<=dis[x]+e[i].cost) continue;
			dis[to]=dis[x]+e[i].cost;
			pre[to]=x;
			lst[to]=i;
			flow[to]=min(flow[x],e[i].flow);
			if(!ins[to]){
				ins[to]=1;
				q.push(to);
			}
		}
	}
	return dis[T]<inf;
}
int mcmf(){
	int mf=0;//useless in this problem
	int mc=0;
	while(spfa()){//find ONE augmentation path at a time to augment
		mf+=flow[T];
		mc+=flow[T]*dis[T];
		int x=T;
		while(x!=S){
			e[lst[x]  ].flow-=flow[T];
			e[lst[x]^1].flow+=flow[T];
			x=pre[x];
		}
	}
	return mc;
}
signed main(){IOS;
	cin>>n>>k;
	int x;
	For(i,1,n) For(j,1,n){
		cin>>x;
		adde(num(i,j,0),num(i,j,1),1,-x);
		adde(num(i,j,0),num(i,j,1),k,0);//'k' means inf 
	}
	For(i,1,n) For(j,1,n-1) adde(num(i,j,1),num(i,j+1,0),k,0);
	For(i,1,n-1) For(j,1,n) adde(num(i,j,1),num(i+1,j,0),k,0);
	S=1,lim=T=num(n,n,1)+1;
	adde(S,num(1,1,0),k,0);//but only here not,JUST k 
	adde(num(n,n,1),T,k,0);
	cout<<(-mcmf())<<endl;//DONT forget to *(-1)
return 0;}
posted @   ShaoJia  阅读(248)  评论(0)  编辑  收藏  举报
相关博文：
·  重修 网络最大流
·  最小斯坦纳树
·  图论 —— 网络流 —— 费用流 —— MCMF 算法
·  「Note」最小费用最大流 MCMF
·  Note -「Dijkstra 求解 MCMF」
阅读排行：
· TypeScript + Deepseek 打造卜卦网站：技术与玄学的结合
· 阿里巴巴 QwQ-32B真的超越了 DeepSeek R-1吗？
· 【译】Visual Studio 中新的强大生产力特性
· 【设计模式】告别冗长if-else语句：使用策略模式优化代码结构
· 10年+ .NET Coder 心语 ── 封装的思维：从隐藏、稳定开始理解其本质意义
点击右上角即可分享
微信分享提示
SQL AI 助手