kmjp's blog

競技プログラミング参加記です

yukicoder : No.590 Replacement

これは本番解けて良かった。
https://yukicoder.me/problems/no/590

問題

1~Nの順列を成す2つの数列A,Bが与えられる。

今1~Nの範囲の2つの整数x,yがあるとする。
1回処理を行うと、x,yがA[x],B[y]に置き換わるとしよう。
f(x,y) := 両者の数字が一致するまでの最小処理回数 (ただしいくら処理しても一致しない場合は0)
とする。f(x,y)の総和を求めよ。

解法

2つのN頂点のグラフを考える。
それぞれ、x→A[x]という有向辺を張ったグラフと、y→B[y]という有向辺を張ったグラフとすると、それぞれいくつかの閉路の集合となる。
f(x,y)は初期状態で前者のグラフで頂点x,後者のグラフで頂点yから始め、ともに何度か辺に沿って移動したところ同じ数字の頂点にたどり着く最小回数となる。

ここで、前者から1つ、後者から1つの閉路を選んだ場合を考える。
両者が同じ数字aを含むとすると、両閉路に含まれるx,yから開始し、何回か処理を行うと両者が共にa,aになる可能性がある。
よって、逆に最初にaで一致するような(x,y)の数を求めることを考えよう。
前者の閉路長をP、後者をQとすると、処理を繰り返すとLCM(P,Q)回で両方の数字の遷移が1順する。
その過程で、両者の数字が一致するタイミングがいくつかあったとする。
例えば適当な状態から開始して3回目で(a,a)、10回目で(b,b)と一致するタイミングがあったとする。
この場合、その適当な状態から4~10回目相当の位置に相当する頂点対(x,y)を開始位置とすると、それぞれ6~0回で(b,b)に到達する。
よって4~10回目相当の位置に相当する頂点対(x,y)におけるf(x,y)の和は0+1+2+...+6であることがわかる。
このように、前LCM(P,Q)回の状態遷移中、条件を満たす位置が何か所かわかると、それ以外の位置から状態遷移を始めた際条件を満たすまでの回数がわかる。
そのためそこからそれらの総和を求めることができる。

条件を満たす位置の算出については以下のコードでは拡張ユークリッドの互除法を用いた。
また、P,Qが互いに素でない場合、同じ2つの閉路内であってもいくら巡回しても到達できない状態が生じるので、それらは個別に考える。

int N;

int A[2][101010];
int uf[2][101010];
int id[2][101010];
int size[2][101010];
ll mo=1000000007;

map<pair<int,int>,vector<int>> M;

ll ext_gcd(ll p,ll q,ll& x, ll& y) { // get px+qy=gcd(p,q)
	if(q==0) return x=1,y=0,p;
	ll g=ext_gcd(q,p%q,y,x);
	y-=p/q*x;
	return g;
}

pair<ll,ll> crt(ll a1,ll mo1,ll a2,ll mo2) { // return (x,y) y=lcm(a1,a2),x%mo1=a1,x%mo2=a2
	ll g,x,y,z;
	g=ext_gcd(mo1,mo2,x,y);
	a1=(a1%mo1+mo1)%mo1;a2=(a2%mo2+mo2)%mo2;
	if(a1%g != a2%g) return pair<ll,ll>(-1,0); // N/A
	ll lcm=mo1*(mo2/g);
	if(lcm<mo1) return pair<ll,ll>(-2,0); // overflow
	ll v=a1+((a2-a1)%lcm+lcm)*x%lcm*(mo1/g);
	return make_pair(((v%lcm)+lcm) % lcm,lcm);
}

ll hoge(int X,int Y,vector<int> R) {
	int i;
	int g=__gcd(X,Y);
	
	map<int,vector<ll>> V;
	FORR(r,R) {
		int a=id[0][r];
		int b=id[1][r];
		int dif=(a-b+Y)%g;
		a=(a+X-dif)%X;
		pair<ll,ll> tmp=crt(a,X,b,Y);
		V[dif].push_back(tmp.first);
	}
	ll ret=0;
	FORR(r,V) {
		vector<ll> vec=r.second;
		sort(ALL(vec));
		vec.push_back(vec[0]+1LL*X*Y/g);
		FOR(i,vec.size()-1) {
			ll dif=vec[i+1]-vec[i];
			(ret+=(dif%mo)*((dif+mo-1)%mo)%mo*((mo+1)/2))%=mo;
		}
	}
	return ret;
	
}

void solve() {
	int i,j,k,l,r,x,y; string s;
	
	cin>>N;
	FOR(j,2) {
		FOR(i,N) {
			cin>>A[j][i], A[j][i]--;
			uf[j][i]=-1;
		}
		FOR(i,N) if(uf[j][i]==-1) {
			int pre=i;
			size[j][i]=1;
			uf[j][i]=i;
			while(A[j][pre]!=i) {
				pre=A[j][pre];
				uf[j][pre]=i;
				id[j][pre]=size[j][i];
				size[j][i]++;
			}
		}
	}
	FOR(i,N) M[{uf[0][i],uf[1][i]}].push_back(i);
	ll ret=0;
	FORR(r,M) {
		(ret+=hoge(size[0][r.first.first],size[1][r.first.second],r.second))%=mo;
	}
	cout<<ret<<endl;
	
}

まとめ

これ文字で説明するの非常に面倒臭い。
でもEditorialもあんまり文章書いてないな…。