JZOJ 5970 space

Space

Description

一个四维空间中有$n^4$个点，每个点的每一维坐标都是整数，范围为$[1,n]$，接着给出四个$n$的排列$A，B，C，D$，一开始在$(1,1,1,1)$点，现在要遍历四维空间中所有的点，假设当前在点$(x,y,z,h)$，则从这个点去到$(A_x,B_y,C_z,D_h)$的代价为$1$，否则为$2$，求总代价最小的便利方案，求最小代价。

Data Constraints

$n\le 10^5$

Solution

让$(x,y,z,h)$向$(A_x,B_y,C_z,D_h)$连一条边，则每个点的出度入度都为1，则整个四维空间会形成若干个环，显然答案为$n^4$加上环的个数。

对每个排列中的环找出来，显然不同种类不会超过$O(\sqrt{n})$个，对于每个排列中的四个环$a,b,c,d$，他们的贡献显然是$\frac{abcd}{lcm(a,b,c,d)}$。
将$AB$的环暴力合并，$CD$的环也暴力合并，复杂度是$O(n)$的，大概的形式是可以得到这样的四个数组$r,r^{\prime },w,w^{\prime }$，现在要我们求${\sum }_i{\sum }_j\frac{r_i{w}_j}{lcm(r_i^{\prime },w_j^{\prime })}$
我们来推推式子
$$\sum _i\sum _j\frac{r_i{w}_j}{lcm(r_i^{\prime },w_j^{\prime })}$$

$$\sum _i\sum _j\frac{r_i{w}_j}{r_i^{\prime }{w}_j^{\prime }}\ast gcd(r_i^{\prime },w_j^{\prime })$$

令 $r_i=\frac{r_i}{r_i^{\prime }}{w}_i=\frac{w_i}{w_i^{\prime }}$

$$\sum _i\sum _j{r}_i{w}_{j}gcd(r_i^{\prime },w_j^{\prime })$$

$$\sum _i{r}_i\sum _{d|r_i^{\prime }}d\ast \sum _j[gcd(w_j^{\prime },r_i^{\prime })=d]w_j$$

$$\sum _i{r}_i\sum _{d|r_i^{\prime }}d\ast \sum _j[gcd(\frac{w_j^{\prime }}{d},\frac{r_i^{\prime }}{d})=1]w_j$$

$$\sum _i{r}_i\sum _{d|r_i^{\prime }}d\ast \sum _j{w}_j\sum _{k|\frac{w_j^{\prime }}{d}\&\&k|\frac{r_i^{\prime }}{d}}\mu (k)$$

$$\sum _i{r}_i\sum _{d|r_i^{\prime }}d\ast \sum _{k|\frac{r_i^{\prime }}{d}}\mu (k)\sum _j[kd|w_j^{\prime }]w_j$$

记$F(v)={\sum }_j[v|w_j^{\prime }]w_j$，则$F(v)$可以$O(\tau (n^2)n)$预处理出来。

$$\sum _i{r}_i\sum _{d|r_i^{\prime }}d\ast \sum _{k|\frac{r_i^{\prime }}{d}}\mu (k)F(kd)$$

记$t=kd$，更换主体，得

$$\sum _i{r}_i\sum _{t|r_i^{\prime }}F(t)\ast \sum _{d|r_i^{\prime }}d\ast \mu (\frac{r_i^{\prime }}{d})$$

$$\sum _i{r}_i\sum _{t|r_i^{\prime }}F(t)\ast \phi (t)$$

$O(\tau (n^2)n)$暴力枚举约数统计答案即可。

Code

#include<iostream>
#include<cstdio>
#include<algorithm>
#include<cstring>
#include<cmath>

#define fo(i,j,l) for(int i=j;i<=l;++i)
#define fd(i,j,l) for(int i=j;i>=l;--i)

using namespace std;
typedef long long ll;
const ll N=12e4,K=N*250,P=K-11,mo=998244353;

int a[N],b[N],c[N],d[N];
int d1[N],d2[N],d3[N],d4[N];
int e1[N],e2[N],e3[N],e4[N];
int g1[N],g2[N],g3[N],g4[N];
int bz[N],y[40][2],x[40][2];
int n,gh=0;
ll ans=0;

int R[N][20],gs[N],W[N][20];
ll ke[K][2];

inline int gcd(int a,int b)
{return (b==0)?a:(gcd(b,a%b));}

inline void add(ll x,ll del)
{
	int yu=x%P;
	while(ke[yu][0]&&ke[yu][0]!=x)yu=(yu==P-1)?0:(yu+1);
	ke[yu][0]=x;
	ke[yu][1]=(ke[yu][1]+del)%mo;
}

inline void dg1(int o,ll y,ll gx)
{
	if(o==gh+1){
		add(y,gx);
		return; 
	}
	dg1(o+1,y,gx);
	ll o1=gx*(x[o][0]-1)%mo;
	y=y*x[o][0];
	dg1(o+1,y,o1);
	fo(i,2,x[o][1]){
		o1=o1*(x[o][0])%mo;
		y=y*x[o][0];
		dg1(o+1,y,o1);
	}
}

inline ll ask(ll x)
{
	int yu=x%P;
	while(ke[yu][0]!=x&&ke[yu][0])yu=(yu==P-1)?0:(yu+1);
	return ke[yu][1];
}

inline void dg2(int o,ll y,ll gx)
{
	if(o==gh+1){
		ans=(ans+gx*ask(y))%mo;
		return; 
	}
	dg2(o+1,y,gx);
	fo(i,1,x[o][1]){
		y=y*x[o][0];
		dg2(o+1,y,gx);
	}
}

int main()
{
	scanf("%d",&n);
	fo(i,1,n)scanf("%d",&a[i]);
	fo(i,1,n)scanf("%d",&b[i]);
	fo(i,1,n)scanf("%d",&c[i]);
	fo(i,1,n)scanf("%d",&d[i]);
	fo(i,2,n){
		gs[i]=0;
		int k=i,p=(int)sqrt(i);
		fo(l,2,p)if(k%l==0){
			++gs[i];
			R[i][gs[i]]=l;
			while(k%l==0)++W[i][gs[i]],k/=l;
		}
		if(k!=1)R[i][++gs[i]]=k,W[i][gs[i]]=1;
	}
	
	fo(i,1,n)bz[i]=0;
	fo(i,1,n)if(!bz[i]){
		bz[i]=1;
		int k=1,w=i;
		while(!bz[a[w]])w=a[w],bz[w]=1,++k;
		d1[++d1[0]]=k;
	}
	sort(d1+1,d1+d1[0]+1);
	fo(i,1,d1[0])if(d1[i]!=d1[i-1]||i==1)
	e1[++e1[0]]=d1[i],g1[e1[0]]=1;
	else ++g1[e1[0]];
	
	fo(i,1,n)bz[i]=0;
	fo(i,1,n)if(!bz[i]){
		bz[i]=1;
		int k=1,w=i;
		while(!bz[b[w]])w=b[w],bz[w]=1,++k;
		d2[++d2[0]]=k;
	}
	sort(d2+1,d2+d2[0]+1);
	fo(i,1,d2[0])if(d2[i]!=d2[i-1]||i==1)
	e2[++e2[0]]=d2[i],g2[e2[0]]=1;
	else ++g2[e2[0]];
	
	fo(i,1,e1[0])fo(l,1,e2[0]){
		ll op=gcd(e1[i],e2[l]);
		ll ko=op*g1[i]%mo*g2[l]%mo;
		int j1=1,j2=1,a=e1[i],b=e2[l];
		gh=0;
		while(j1<=gs[a]||j2<=gs[b])
		if(j2>gs[b]||(R[a][j1]<R[b][j2]&&j1<=gs[a])){
			x[++gh][0]=R[a][j1];
			x[gh][1]=W[a][j1];
			++j1;
		}else
		if(j1>gs[a]||(R[b][j2]<R[a][j1]&&j2<=gs[b])){
			x[++gh][0]=R[b][j2];
			x[gh][1]=W[b][j2];
			++j2;
		}else{
			x[++gh][0]=R[a][j1];
			x[gh][1]=W[a][j1]+W[b][j2];
			++j1; ++j2;
		}
		fo(l,1,gh)while(op%x[l][0]==0)op/=x[l][0],--x[l][1];
		dg1(1,1,ko);
	}
	
	fo(i,1,n)bz[i]=0;
	fo(i,1,n)if(!bz[i]){
		bz[i]=1;
		int k=1,w=i;
		while(!bz[c[w]])w=c[w],bz[w]=1,++k;
		d3[++d3[0]]=k;
	}
	sort(d3+1,d3+d3[0]+1);
	fo(i,1,d3[0])if(d3[i]!=d3[i-1]||i==1)
	e3[++e3[0]]=d3[i],g3[e3[0]]=1;
	else ++g3[e3[0]];
	
	fo(i,1,n)bz[i]=0;
	fo(i,1,n)if(!bz[i]){
		bz[i]=1;
		int k=1,w=i;
		while(!bz[d[w]])w=d[w],bz[w]=1,++k;
		d4[++d4[0]]=k;
	}
	sort(d4+1,d4+d4[0]+1);
	fo(i,1,d4[0])if(d4[i]!=d4[i-1]||i==1)
	e4[++e4[0]]=d4[i],g4[e4[0]]=1;
	else ++g4[e4[0]];
	
	fo(i,1,e3[0])fo(l,1,e4[0]){
		ll op=gcd(e3[i],e4[l]);
		ll ko=op*g3[i]%mo*g4[l]%mo;
		int j1=1,j2=1,a=e3[i],b=e4[l];
		gh=0;
		while(j1<=gs[a]||j2<=gs[b])
		if(j2>gs[b]||(R[a][j1]<R[b][j2]&&j1<=gs[a])){
			x[++gh][0]=R[a][j1];
			x[gh][1]=W[a][j1];
			++j1;
		}else
		if(j1>gs[a]||(R[b][j2]<R[a][j1]&&j2<=gs[b])){
			x[++gh][0]=R[b][j2];
			x[gh][1]=W[b][j2];
			++j2;
		}else{
			x[++gh][0]=R[a][j1];
			x[gh][1]=W[a][j1]+W[b][j2];
			++j1; ++j2;
		}
		fo(l,1,gh)while(op%x[l][0]==0)op/=x[l][0],--x[l][1];
		dg2(1,1,ko);
	}
	ans=(ans+(ll)n*n%mo*n%mo*n)%mo;
	cout<<ans;
}