読者です 読者をやめる 読者になる 読者になる

メモ

yukicoderで遊んでいる競プロゆるふわ勢

yukicoder No.181 A↑↑N mod M

問題はこちら
No.181 A↑↑N mod M - yukicoder

よく考えればそんなに難しくなくない…?と最初は思ったけれど、考えれば考えるほど泥沼になるタイプの問題だった

m=1のとき、及び、n≦2のときは明らか。そうでないとする。
a↑↑nを返す関数をtetra(a,n)、a↑↑n mod mを返す関数をmodtetra(a,n,m)、a^n mod mを返す関数をmodpow(a,n,m)とする。
modpowはa,n,mがintの範囲なら繰り返し二乗法で適当に実装できる
modtetra(a,n,m)が求めるもの。定義からmodtetra(a,n,m)=modpow(a,tetra(a,n-1),m)。
当然tetra(a,n-1)はそのままだとでかいので、ここも何かでいい感じに剰余を取る必要がある。
aとmが互いに素ならオイラーの定理からa^{\phi(m)}\equiv1\bmod mなので、modpow(a,n,m)=modpow(a,n%φ(m),m)が任意のnで成立するが、今回はそうとは限らない。
しかし実際には次の補題により、ほぼ同様のことが言える事がわかる

補題1. m\geqq 2 とする。任意の a に対し次が成立。
n\geqq\log_2 m \Rightarrow a^{n+\phi(m)}\equiv a^n\bmod m
証明. m素因数分解したものを \prod_i p_i^{e_i} とする。中国剰余定理から \bmod p_i^{e_i} ごとに合同式の成立が言えれば良い
ap_i を因数に持たない時
ap_i^{e_i} は互いに素なのでオイラーの定理から a^{\phi(p_i^{e_i})}\equiv 1\bmod p_i^{e_i}
\phi(p_i^{e_i})\phi(m) の約数なので a^{\phi(m)}\equiv 1\bmod p_i^{e_i} となり、
任意の n\geqq 0 に対し a^{n+\phi(m)}\equiv a^n\bmod p_i^{e_i} が成立。
ap_i を因数に持つ時
n\geqq {e_i} \Rightarrow a^n \equiv 0\bmod p_i^{e_i} より
n\geqq {e_i} \Rightarrow a^{n+\phi(m)}\equiv a^n\ mod \ p_i^{e_i}
2^{e_i}\leqq p_i^{e_i}\leqq m なので当然 {e_i}\leqq\log_2 m が成立。
よって n\geqq\log_2 m \Rightarrow n\geqq {e_i} となり求めるものが成立。
以上より示せた


ということでやるべき計算は「tetra(a,n-1)がlog2(m)以上ならmodをとってからφ(m)を足す。そうでないならそのまま」となる
そうすると今度はtetra(a,n-1)がlog2(m)以上かどうかを考えることになるが、
・制約からlog2(m)≦log2(2000)<11
・n=3かつa=2のときtetra(a,n-1)=4
・n=3かつa≧3のときtetra(a,n-1)≧tetra(3,2)=27
・n≧4のときtetra(a,n-1)≧tetra(2,3)=16
ということから、結局「n=3かつa=2」の場合を除き必ずφ(m)を足せば良い事がわかる

よって以上によりmodtetraを次のような再帰関数で書くことができた

modtetra(a,n,m){
	if(m==1)return 0;
	if(n==0)return 1;
	if(n==1)return a%m;
	if(n==2)return modpow(a%m,a,m);
	if(n==3&&a==2)return 16%m;
	return modpow(a%m,modtetra(a,n-1,phi(m))+phi(m),m);
}

ところでこれは間に合うのだろうか?

補題2.\phi(\phi(m))< m/2
証明.m素因数分解したものを \prod_i p_i^{e_i}とすれば、
\phi(m)=\prod(p_i-1)p_i^{e_i-1}=m\prod(1-\frac{1}{p_i})である
1つ目の形から、mが奇数ならば\phi(m)は偶数となることが分かり、
2つ目の形から、mが偶数ならば\phi(m)\leqq m/2となることが分かる。
これらから従う

ということで、modtetraは呼び出すたびに第3引数にφを適用しているので、高々2*log2(M)<22回の再帰で止まる。
これなら十分間に合いそう。
φ(m)はm以下の数全てとgcdを取るという愚直な計算で間に合う

int gcd(int p,int q){return q?gcd(q,p%q):p;}
int phi(int m){
	int i,s=0;
	for(i=1;i<m;i++)s+=gcd(m,i)==1;
	return s;
}
int modpow(int a,int n,int m){
	if(n==0)return 1;
	return modpow(a*a%m,n/2,m)*(n%2?a:1)%m;
}
int modtetra(int a,int n,int m){
	if(m==1)return 0;
	if(n==0)return 1;
	if(n==1)return a%m;
	if(n==2)return modpow(a%m,a,m);
	if(n==3&&a==2)return 16%m;
	return modpow(a%m,modtetra(a,n-1,phi(m))+phi(m),m);
}

int main(){
	int a,n,m;
	scanf("%d%d%d",&a,&n,&m);
	printf("%d",modtetra(a,n,m));
	return 0;
}


========
計算量の話(間違ってたら教えてください)
modtetra(a,n,m)を呼び出す毎に行われる計算はmodpowとphiの2つ。
phi(m)を求めるのは、以下で述べるように、素因数分解と同等のことをすれば良いので、O(√m)。
powmodの計算は、指数が高々2*φ(m)<2mなのでO(log(m))以下。再帰ごとに重いのはphiの計算。
ここで2回再帰する毎にmは半分未満になっていくので全体では
2*√m*(1+1/2+…)という感じとなり、トータルでO(√m)
========


トーシェント関数の計算は、補題の2つ目の変形を用いてgcdを使わずに書くことが出来る

phi(m){i=1;for(s=m;m/++i;)for(s-=m%i?0:s/i;m%i<1;m/=i);return s;}

またmodtetraもぐっと睨むと、ほとんどすべての場合はmodpowにまとめられる事がわかるので、トータルでは以下のよう

s,i;p(m){i=1;for(s=m;m/++i;)for(s-=m%i?0:s/i;m%i<1;m/=i);return s;}
g(a,n,m){return n?n%2?g(a*a%m,n/2,m)*a%m:g(a*a%m,n/2,m):1;}
f(a,n,m){return~-m?n?g(a%m,--n?n-2|a-2?n-1?p(m)+f(a,n,p(m)):a:4:1,m):1:0;}
n,m;main(a){scanf("%d%d%d",&a,&n,&m);n=!printf("%d",f(a,n,m));}

p(m)とf(a,n,m)の計算をよく睨むと、グローバル変数sにφ(m)の値が入るタイミングをみて

s,i;p(m){i=1;for(s=m;m/++i;)for(s-=m%i?0:s/i;m%i<1;m/=i);return s;}
g(a,n,m){return n?n%2?g(a*a%m,n/2,m)*a%m:g(a*a%m,n/2,m):1;}
f(a,n,m){return~-m?n?g(a%m,--n?n-2|a-2?n-1?p(m)+f(a,n,s):a:4:1,m):1:0;}
n,m;main(a){scanf("%d%d%d",&a,&n,&m);n=!printf("%d",f(a,n,m));}

とできる。
gの計算は再帰を使わずに書けば

r;g(a,n,m){for(r=1;n;n/=2)n&1?r=r*a%m:0,a=a*a%m;return r;}

となることから、rに値が代入されるタイミングをみてreturnを省略でき

s,i;p(m){i=1;for(s=m;m/++i;)for(s-=m%i?0:s/i;m%i<1;m/=i);return s;}
r;g(a,n,m){for(r=1;n;n/=2)n&1?r=r*a%m:0,a=a*a%m;}
f(a,n,m){return~-m?n?g(a%m,--n?n-2|a-2?n-1?p(m)+f(a,n,s):a:4:1,m),r:1:0;}
n,m;main(a){scanf("%d%d%d",&a,&n,&m);n=!printf("%d",f(a,n,m));}

252B
だったらpのreturnも省略できるんじゃね?と思い

s,i;p(m){i=1;for(s=m;m/++i;)for(s-=m%i?0:s/i;m%i<1;m/=i);}
r;g(a,n,m){for(r=1;n;n/=2)n&1?r=r*a%m:0,a=a*a%m;}
f(a,n,m){return~-m?n?g(a%m,--n?n-2|a-2?n-1?p(m),s+f(a,n,s):a:4:1,m),r:1:0;}
n,m;main(a){scanf("%d%d%d",&a,&n,&m);n=!printf("%d",f(a,n,m));}

とするとWA。うーむ