BZOJ.4514.[SDOI2016]数字配对(费用流SPFA 二分图)

\(Solution\)

很显然的建二分图后跑最大费用流，但有个问题是一个数是只能用一次的，这样二分图两部分都有这个数。
那么就用两倍的。如果\(i\)可以向\(j'\)连边，\(j\)也向\(i'\)连边，如果上一次走了\(i->j'\)，那么这一次一定走\(j->i'\)。
每次跑最大费用流，直至有一次费用变成负，然后加上当前正权值能抵消它的流量，最后总流量除以2就可以了。

\(Another Solution\)

两个数能匹配首先要能整除，其次它们所有质因子的次数和一定只相差1.
于是可以按这个次数的奇偶性建二分图，分别连源点汇点，一遍最大流就可以了。

~~ps：sb的我硬是把\(i\to j\)的费用设成\(0\)，然后源点汇点连\(\log b[i]\)的边，然后再取幂。。最后写完发现不对，不会处理上面的问题。。贪心也调不对。~~
数字匹配\(O(sqrt(n))\)直接枚举就可以，闲的写什么MillerRabin+set。。

~~考试的时候一个题写了三百多行还不算颓~~
时隔一年后的考试写出来了2333（虽然写过还是道水题）

第一种方法：

//3212kb	112ms
#include <set>
#include <queue>
#include <cmath>
#include <cstdio>
#include <cctype>
#include <algorithm>
#define gc() getchar()
typedef long long LL;
const int N=410,M=100005,INF=1e9;
const int P[9]={2,3,5,7,13,17,31};

int n,src,des,A[N],B[N],C[N],Enum,H[N],nxt[M],fr[M],to[M],cap[M],pre[N];
bool inq[N];
LL dis[N],cost[M];
std::queue<int> q;
std::set<int> is_p[2];

inline int read()
{
	int now=0,f=1;register char c=gc();
	for(;!isdigit(c);c=gc()) if(c=='-') f=-1;
	for(;isdigit(c);now=now*10+c-'0',c=gc());
	return now*f;
}
inline void AddEdge(int u,int v,int w,LL c)
{
	to[++Enum]=v, fr[Enum]=u, nxt[Enum]=H[u], H[u]=Enum, cap[Enum]=w, cost[Enum]=c;
	to[++Enum]=u, fr[Enum]=v, nxt[Enum]=H[v], H[v]=Enum, cap[Enum]=0, cost[Enum]=-c;
}
int FP(LL x,int k,int p)
{
	LL t=1;
	for(; k; k>>=1,x=x*x%p)
		if(k&1) t=t*x%p;
	return (int)t;
}
bool Miller_Rabin(int x)
{
	if(x==2) return 1;
	if(!(x&1)||x==1) return 0;
	for(int i=0; i<7; ++i)
		if(x==P[i]) return 1;
		else if(!(x%P[i])) return 0;
	int u=x-1,t=0;
	while(!(u&1)) u>>=1,++t;
	for(int i=0; i<7; ++i)
	{
		LL now=FP(P[i],u,x),las;
		for(int j=1; j<=t; ++j)
		{
			las=now, (now*=now)%=x;
			if(now==1&&las!=1&&las!=x-1) return 0;
		}
		if(now!=1) return 0;
	}
	return 1;
}
bool SPFA()
{
	for(int i=1; i<=des; ++i) dis[i]=-1e12;
	dis[src]=0, q.push(src);
	while(!q.empty())
	{
		int x=q.front(); q.pop();
		inq[x]=0;
		for(int v,i=H[x]; i; i=nxt[i])
			if(cap[i]&&dis[to[i]]<dis[x]+cost[i])
			{
				dis[v=to[i]]=dis[x]+cost[i], pre[v]=i;
				if(!inq[v]) q.push(v),inq[v]=1;
			}
	}
	return dis[des]>-1e12;
}
LL Solve()
{
	int flow=0,mn;
	LL sum_cost=0;
	while(SPFA())
	{//注意这只是一个二分图。。
		mn=1e9;
		for(int i=des; i!=src; i=fr[pre[i]]) mn=std::min(cap[pre[i]],mn);
		if(sum_cost+1ll*dis[des]*mn>=0){//dis[des]就是费用。
			for(int i=des; i!=src; i=fr[pre[i]]) cap[pre[i]]-=mn,cap[pre[i]^1]+=mn;
			sum_cost+=1ll*dis[des]*mn, flow+=mn;
		}
		else{
			flow+=(int)(sum_cost/std::abs(dis[des]));//加上剩余正价值能补多少负值。
			break;
		}
	}
	return flow>>1;
}

int main()
{
	n=read(),Enum=1,src=0,des=n<<1|1;
	for(int i=1; i<=n; ++i) A[i]=read();
	for(int i=1; i<=n; ++i) B[i]=read();
	for(int i=1; i<=n; ++i) C[i]=read();
	for(int i=1; i<=n; ++i) AddEdge(src,i,B[i],0),AddEdge(i+n,des,B[i],0);
	for(int i=1; i<=n; ++i)
		for(int j=1; j<=n; ++j)
			if(A[i]>A[j]&&!(A[i]%A[j]))
			{
				int t=A[i]/A[j];
				if(is_p[0].count(t)) continue;
				if(is_p[1].count(t)) AddEdge(i,j+n,INF,1ll*C[i]*C[j]),AddEdge(j,i+n,INF,1ll*C[i]*C[j]);
				else{
					bool f=Miller_Rabin(t);
					if(f) is_p[1].insert(t),AddEdge(i,j+n,INF,1ll*C[i]*C[j]),AddEdge(j,i+n,INF,1ll*C[i]*C[j]);
					else is_p[0].insert(t);
				}
			}
	printf("%lld",Solve());

	return 0;
}

第二种方法（19.2.17）：

/*
2840kb	28ms
常规做建二分图跑费用流很容易，但因为拆了点每对数之间可能选两次。但是既然选了一次那就一定会选第二次。最后费用即将变成负数时判一下还能选多少个，然后流量除以二即可。
因为匹配条件是$i$可整除$j$且结果是质数，那么$i,j$之间质因子次数一定相差一，可以直接判，同时根据质因子次数奇偶性可以分别连源点汇点，就没有之前的问题了。
另外这张二分图里费用就是$dis[T]$。
*/
#include <queue>
#include <cstdio>
#include <cctype>
#include <algorithm>
#define gc() getchar()
typedef long long LL;
const int N=205,M=(N*N+N)*2,INF=1<<30;
const LL INFll=1ll<<60;

int T,Enum,H[N],nxt[M],fr[M],to[M],cap[M],pre[N];
LL dis[N],cost[M];

inline int read()
{
	int now=0,f=1;register char c=gc();
	for(;!isdigit(c);c=='-'&&(f=-1),c=gc());
	for(;isdigit(c);now=now*10+c-48,c=gc());
	return now*f;
}
inline void AE(int u,int v,int w,LL c)
{
	to[++Enum]=v, fr[Enum]=u, nxt[Enum]=H[u], H[u]=Enum, cap[Enum]=w, cost[Enum]=c;
	to[++Enum]=u, fr[Enum]=v, nxt[Enum]=H[v], H[v]=Enum, cap[Enum]=0, cost[Enum]=-c;
}
int Calc(int x)
{
	int t=0;
	for(int i=2; 1ll*i*i<=x; ++i)
		if(!(x%i))
		{
			x/=i, ++t;
			while(!(x%i)) x/=i, ++t;
		}
	if(x!=1) ++t;
	return t;
}
bool SPFA()
{
	static bool inq[N];
	static std::queue<int> q;
	for(int i=1; i<=T; ++i) dis[i]=-INFll;
	q.push(0), dis[0]=0;
	while(!q.empty())
	{
		int x=q.front(); q.pop();
		inq[x]=0;
		for(int i=H[x],v; i; i=nxt[i])
			if(dis[v=to[i]]<dis[x]+cost[i] && cap[i])
				pre[v]=i, dis[v]=dis[x]+cost[i], !inq[v]&&(q.push(v),inq[v]=1);
	}
	return dis[T]>-INFll;
}
int Solve()
{
	int res=0; LL now=0;
	while(SPFA())
	{
		int mn=INF;
		for(int i=T; i; i=fr[pre[i]]) mn=std::min(mn,cap[pre[i]]);
		if(now+dis[T]*mn>=0)
		{
			now+=dis[T]*mn, res+=mn;
			for(int i=T; i; i=fr[pre[i]]) cap[pre[i]]-=mn, cap[pre[i]^1]+=mn;
		}
		else
		{
			res+=now/std::abs(dis[T]);//abs!!
			break;
		}
	}
	return res;
}

int main()
{
	static int A[N],B[N],C[N],tm[N];

	freopen("pair.in","r",stdin);
	freopen("pair.out","w",stdout);

	int n=read(); Enum=1, T=n+1;
	for(int i=1; i<=n; ++i) tm[i]=Calc(A[i]=read());
	for(int i=1; i<=n; ++i) B[i]=read();
	for(int i=1; i<=n; ++i) C[i]=read();
	for(int i=1; i<=n; ++i)
	{
		if(tm[i]&1)
		{
			AE(0,i,B[i],0);
			for(int j=1; j<=n; ++j)
				if(!(tm[j]&1) && tm[i]-tm[j]==1 && !(A[i]%A[j])) AE(i,j,INF,1ll*C[i]*C[j]);
		}
		else
		{
			AE(i,T,B[i],0);
			for(int j=1; j<=n; ++j)
				if(tm[j]&1 && tm[i]-tm[j]==1 && !(A[i]%A[j])) AE(j,i,INF,1ll*C[i]*C[j]);
		}
	}
	printf("%d\n",Solve());

	return 0;
}

弃疗的暴力+最大流（18.4.1）：

#include <set>
#include <queue>
#include <cstdio>
#include <cctype>
#include <algorithm>
#define gc() getchar()
typedef long long LL;
const int N=410,M=100005,INF=1e9;
const int P[9]={2,3,5,7,13,17,31};

int n,range,src,des,A[N],B[N],C[N],Enum,cur[N],H[N],nxt[M],fr[M],to[M],cap[M],pre[N],lev[N],num[N],que[N];
std::queue<int> q;
std::set<int> is_p[2];

inline int read()
{
	int now=0,f=1;register char c=gc();
	for(;!isdigit(c);c=gc()) if(c=='-') f=-1;
	for(;isdigit(c);now=now*10+c-'0',c=gc());
	return now*f;
}
inline void AddEdge(int u,int v,int w)
{
	printf("%d->%d %d\n",u,v,w);
	to[++Enum]=v, fr[Enum]=u, nxt[Enum]=H[u], H[u]=Enum, cap[Enum]=w;
	to[++Enum]=u, fr[Enum]=v, nxt[Enum]=H[v], H[v]=Enum, cap[Enum]=0;
}
int FP(LL x,int k,int p)
{
	LL t=1;
	for(; k; k>>=1,x=x*x%p)
		if(k&1) t=t*x%p;
	return (int)t;
}
bool Miller_Rabin(int x)
{
	if(x==2) return 1;
	if(!(x&1)||x==1) return 0;
	for(int i=0; i<7; ++i)
		if(x==P[i]) return 1;
		else if(!(x%P[i])) return 0;
	int u=x-1,t=0;
	while(!(u&1)) u>>=1,++t;
	for(int i=0; i<7; ++i)
	{
		LL now=FP(P[i],u,x),las;
		for(int j=1; j<=t; ++j)
		{
			las=now, (now*=now)%=x;
			if(now==1&&las!=1&&las!=x-1) return 0;
		}
		if(now!=1) return 0;
	}
	return 1;
}
namespace Subtask1
{
	int Ans;
	bool can[13][13],vis[13];

	void DFS(int x,int tm,LL sum)
	{//O(n!)
		if(x>n)
			if(sum>=0 && tm>Ans) Ans=tm;
			else ;
		else{
			if(vis[x]) DFS(x+1,tm,sum);
			else{
				vis[x]=1;
				for(int i=1; i<=n; ++i)
					if(can[x][i]&&!vis[i])
						vis[i]=1, DFS(x+1,tm+1,sum+1ll*C[x]*C[i]), vis[i]=0;
				vis[x]=0;
				DFS(x+1,tm,sum);
			}
		}
	}
	void Solve()
	{
		for(int i=1; i<=n; ++i) A[i]=read();
		for(int i=1; i<=n; ++i) B[i]=read();
		for(int i=1; i<=n; ++i) C[i]=read();
		for(int i=1; i<=n; ++i)
			for(int j=1; j<=n; ++j)
				if(!(A[i]%A[j])&&i!=j)
				{
					int t=A[i]/A[j];
					if(is_p[0].count(t)) continue;
					if(is_p[1].count(t)) can[i][j]=1;
					else{
						bool f=Miller_Rabin(t);
						if(f) is_p[1].insert(t),can[i][j]=1;
						else is_p[0].insert(t);
					}
				}
		DFS(1,0,0);
		printf("%d",Ans);
	}
}
bool BFS()
{
	for(int i=src; i<des; ++i) lev[i]=des+1;
	lev[des]=0, que[0]=des; int h=0,t=1;
	while(h<t)
	{
		int x=que[h++];
		for(int i=H[x]; i; i=nxt[i])
			if(lev[to[i]]==des+1 && cap[i^1])
				lev[to[i]]=lev[x]+1, que[t++]=to[i];
	}
	return lev[src]<=des;
}
void Augment(){
	for(int i=des; i!=src; i=fr[pre[i]])
		--cap[pre[i]], ++cap[pre[i]^1];
}
int ISAP()
{
	if(!BFS()) return 0;
	for(int i=src; i<=des; ++i) ++num[lev[i]],cur[i]=H[i];
	int x=src,res=0; bool f=0;
	while(lev[src]<=des)
	{
		if(x==des) x=src,f=0,++res,Augment();
		bool can=0;
		for(int i=cur[x]; i; i=nxt[i])
			if((to[i]==des&&f)||(to[i]!=des &&  cap[i]))
			{
				printf("%d->%d lev:%d f:%d ",x,to[i],lev[x],f);
				if(x>n&&to[i]<=n) f=1;
				can=1, cur[x]=i, pre[x=to[i]]=i;
				printf("now:%d\n",f);
				break;
			}
		if(!can)
		{
			int mn=des;
			for(int i=H[x]; i; i=nxt[i])
				if((to[i]==des&&f)||(to[i]!=des&&cap[i])) mn=std::min(mn,lev[to[i]]);
			if(!--num[lev[x]]) break;
			printf("Change:%d %d->",x,lev[x]);
			++num[lev[x]=mn+1];
			printf("%d\n",lev[x]);
			cur[x]=H[x];
			if(x!=src)
			{
				if(x<=n&&fr[pre[x]]>n) f=0;
				x=fr[pre[x]];
			}
		}
	}
	return res;
}

int main()
{
	freopen("pair.in","r",stdin);
	freopen("pair.out","w",stdout);

	n=read();
	if(n<=10) {Subtask1::Solve(); return 0;}
	Enum=1,src=0,des=n<<1|1;
	for(int i=1; i<=n; ++i) A[i]=read();
	for(int i=1; i<=n; ++i) B[i]=read();
	for(int i=1; i<=n; ++i) C[i]=read();
	for(int i=1; i<=n; ++i) AddEdge(src,i,INF),AddEdge(i+n,des,INF),AddEdge(i,i+n,B[i]);
	for(int i=1; i<=n; ++i)
		for(int j=1; j<=n; ++j)
			if(!(A[i]%A[j])&&i!=j)
			{
				int t=A[i]/A[j];
				if(is_p[0].count(t)) continue;
				if(is_p[1].count(t)) AddEdge(i+n,j,INF);
				else{
					bool f=Miller_Rabin(t);
					if(f) is_p[1].insert(t),AddEdge(i+n,j,INF);
					else is_p[0].insert(t);
				}
			}
	printf("%d",ISAP());

	return 0;
}

丧心病狂的费用流~~(其实就是sb了)~~（18.4.1）

#include <set>
#include <queue>
#include <cmath>
#include <cstdio>
#include <cctype>
#include <algorithm>
#define gc() getchar()
typedef long long LL;
const int N=410,M=100005,INF=1e9;
const int P[9]={2,3,5,7,13,17,31};

int n,range,src,des,A[N],B[N],C[N],Enum,H[N],nxt[M],fr[M],to[M],cap[M],Cap[M],pre[N];
double dis[N],cost[M],neg[M];
bool inq[N],ignore[M];
std::queue<int> q;
std::set<int> is_p[2];

inline int read()
{
	int now=0,f=1;register char c=gc();
	for(;!isdigit(c);c=gc()) if(c=='-') f=-1;
	for(;isdigit(c);now=now*10+c-'0',c=gc());
	return now*f;
}
inline void AddEdge(int u,int v,int w,int c)
{
	if(!c){
		to[++Enum]=v, fr[Enum]=u, nxt[Enum]=H[u], H[u]=Enum, Cap[Enum]=cap[Enum]=w, cost[Enum]=0, ignore[Enum]=1;
		to[++Enum]=u, fr[Enum]=v, nxt[Enum]=H[v], H[v]=Enum, Cap[Enum]=cap[Enum]=0, cost[Enum]=0, ignore[Enum]=1;
	}
	else if(c<0){
		to[++Enum]=v, fr[Enum]=u, nxt[Enum]=H[u], H[u]=Enum, Cap[Enum]=cap[Enum]=w, cost[Enum]=-log2((double)(-c)), neg[Enum]=1;
		to[++Enum]=u, fr[Enum]=v, nxt[Enum]=H[v], H[v]=Enum, Cap[Enum]=cap[Enum]=0, cost[Enum]=log2((double)(-c)), neg[Enum]=1;
	}
	else{
		to[++Enum]=v, fr[Enum]=u, nxt[Enum]=H[u], H[u]=Enum, Cap[Enum]=cap[Enum]=w, cost[Enum]=log2((double)c);
		to[++Enum]=u, fr[Enum]=v, nxt[Enum]=H[v], H[v]=Enum, Cap[Enum]=cap[Enum]=0, cost[Enum]=-log2((double)c);
	}
	printf("%d->%d %d %d %.3lf\n",u,v,w,c,cost[Enum-1]);
}
int FP(LL x,int k,int p)
{
	LL t=1;
	for(; k; k>>=1,x=x*x%p)
		if(k&1) t=t*x%p;
	return (int)t;
}
bool Miller_Rabin(int x)
{
	if(x==2) return 1;
	if(!(x&1)||x==1) return 0;
	for(int i=0; i<7; ++i)
		if(x==P[i]) return 1;
		else if(!(x%P[i])) return 0;
	int u=x-1,t=0;
	while(!(u&1)) u>>=1,++t;
	for(int i=0; i<7; ++i)
	{
		LL now=FP(P[i],u,x),las;
		for(int j=1; j<=t; ++j)
		{
			las=now, (now*=now)%=x;
			if(now==1&&las!=1&&las!=x-1) return 0;
		}
		if(now!=1) return 0;
	}
	return 1;
}
bool SPFA()
{
	for(int i=1; i<=des; ++i) dis[i]=1e14;
	dis[src]=0, q.push(src);
	while(!q.empty())
	{
		int x=q.front(); q.pop();
		inq[x]=0;
		for(int v,i=H[x]; i; i=nxt[i])
			if(cap[i]&&(neg[i]||dis[to[i]]>dis[x]+cost[i]))
			{
				dis[v=to[i]]=dis[x]+cost[i], pre[v]=i;
				if(!inq[v]) q.push(v);
			}
	}
	return dis[des]<1e14;
}
bool debug=1;
double MCMF(LL &flow)
{
	int mn=1e9; double res=0,tmp1=0;
	bool f=0;
	for(int i=des; i!=src; i=fr[pre[i]])
	{
		mn=std::min(cap[pre[i]],mn);
		if(ignore[i]&&i<=range) f=1;
	}
//	printf("flow:%I64d ",flow);
	if(mn>=flow) mn=flow,flow=0;
	else flow-=mn;
	if(f){
		for(int i=des; i!=src; i=fr[pre[i]])
		{
			if(debug) printf("%d<-",i),
			cap[pre[i]]-=mn, cap[pre[i]^1]+=mn;
			if(neg[i]) ;
		}
		if(debug) printf("%d\n",src);
		if(debug) printf("then:%I64d %d\n",flow,mn);
		return 0;
	}
	for(int i=des; i!=src; i=fr[pre[i]])
	{
		if(debug) printf("%d<-",i),
		cap[pre[i]]-=mn, cap[pre[i]^1]+=mn, res+=cost[pre[i]];
		if(neg[pre[i]]&&!cost[pre[i]]) tmp1+=1.0;
	}
	if(debug) printf("%d\n",src);
	if(debug) printf("then:%I64d res:%.3lf 2:%.3lf\n",flow,res,res<0 ? -1.0*mn*pow(2.0,-res) : 1.0*mn*pow(2.0,res));
	return res<0 ? -1.0*mn*pow(2.0,-res)-tmp1 : 1.0*mn*pow(2.0,res)-tmp1;
}
bool Check(LL flow)
{
	double sigma=0;
	LL tmp=flow;
	for(int i=2; i<=Enum; ++i) cap[i]=Cap[i];
	while(SPFA())
	{
		sigma+=MCMF(flow);
		if(!flow||sigma>0) break;
	}
	if(debug) printf("End:flow:%I64d bef:%I64d sigma:%.3lf\n\n",flow,tmp,sigma);
	return !flow&&sigma<=0;
}

int main()
{
	freopen("pair.in","r",stdin);
//	freopen("pair.out","w",stdout);

	n=read(),Enum=1,src=0,des=n<<1|1;
	for(int i=1; i<=n; ++i) A[i]=read();
	for(int i=1; i<=n; ++i) B[i]=read();
	for(int i=1; i<=n; ++i) C[i]=read();
	for(int i=1; i<=n; ++i) AddEdge(src,i,INF,0),AddEdge(i+n,des,INF,0),AddEdge(i,i+n,B[i],-C[i]);
//	range=Enum;
	for(int i=1; i<=n; ++i)
		for(int j=1; j<=n; ++j)
			if(!(A[i]%A[j])&&i!=j)
			{
				int t=A[i]/A[j];
				if(is_p[0].count(t)) continue;
				if(is_p[1].count(t)) AddEdge(i+n,j,INF,0);
				else{
					bool f=Miller_Rabin(t);
					if(f) is_p[1].insert(t),AddEdge(i+n,j,INF,0);
					else is_p[0].insert(t);
				}
			}
	
//	long long l=0,r=1e10,mid,ans=0;
	long long l=1,r=15,mid,ans=0;
	while(l<=r)
		if(Check(mid=l+r>>1)) ans=mid,l=mid+1;
		else r=mid-1;
	printf("%I64d",ans);

	return 0;
}

posted @ 2018-04-01 14:19 SovietPower 阅读(227) 评论(0) 编辑收藏举报

刷新页面返回顶部

SovietPower

因你而起的一泓喜悲权当年轻留个纪念

BZOJ.4514.[SDOI2016]数字配对(费用流SPFA 二分图)

\(Solution\)

\(Another Solution\)

公告

SovietPower

因你而起的一泓喜悲 权当年轻 留个纪念

BZOJ.4514.[SDOI2016]数字配对(费用流SPFA 二分图)

\(Solution\)

\(Another Solution\)

公告

因你而起的一泓喜悲权当年轻留个纪念