本番つめが甘くて正答できず。
アプローチまではよかったんだけどなぁ。
http://community.topcoder.com/stat?c=problem_statement&pm=13016
問題
数N,Kが与えられる。
low～highのうちのいずれかの整数をN個並べた数列を作る。
この時、数列の全要素の最大公約数がKとなる数列の作り方は何通りあるか。
解法
まず数列中の要素はすべてKの倍数でなければならないので、lowより大きくhighより小さいKの倍数のみが解の候補になりうる。
よって、low=(low+(K-1))/K、high=high/Kで置き換えよう。
以後はKを無視でき、数列の全要素の最大公約数を1にすることを考えればよい。
i=1～(high-low)の範囲で総当たりし、「選んだ数値がいずれもiの倍数になる」ような数の組み合わせを考える。
この数の計算は容易で(low+(i-1))/i*i～hi/i*iの間のiの倍数である。
この数の個数をKとすると、それらで構成される数列はK^N-K通りである。
K引いているのは全部が同じ数字だと、各数値がiの倍数であることが意味を持たない(iが何だろうとなんかしらの倍数になるため)である。
あとはこれを包除原理の要領で足したり引いたりすればよい。
素因数の数の偶奇に合わせて(K^N-K)を足し引きする。
ただし素因数が2乗以上になっている場合は、足しも引きもしない。
9の倍数の数の集合は3の倍数の集合に含まれているので、さらに足し引きする必要がないためである。
最後に注意点として、low=1の場合、lowだけで構成される数列は題意を満たすので1加算する。
本番は(K^N-K)の-Kの下りを導きだせず正答できなかった。
上記包除原理や平方数の下りだが、どうもメビウス関数・メビウスの反転公式と関係があるようだ。
後でもう一度蟻本を確認してみよう。
ll mo=1000000007; int NP,prime[100000]; const int prime_max = 100000; void cprime(int hoge) { static int p[prime_max]; int i,j; NP=0; for(i=2;i<hoge;i++) { if(p[i]) continue; prime[NP++]=i; for(j=i*2;j<hoge;j+=i) p[j]=i; } } ll modpow(ll a, ll n) { ll r=1; while(n) r=r*((n%2)?a:1)%mo,a=a*a%mo,n>>=1; return r; } int vis[100001]; class RandomGCD { public: int countTuples(int N, int K, int low, int high) { low=(low+(K-1))/K; high=high/K; if(low>high) return 0; cprime(high-low+1); int i,j; ll ret=0; for(i=1;i<=high-low;i++) { int nd=1,xx=i; int lo=(low+(i-1))/i; int hi=high/i; if(hi<lo) continue; int ng=0; FOR(j,NP) { if(xx<prime[j] || ng) break; if(xx%prime[j]==0) { int ii=0; nd=-nd; while(xx%prime[j]==0) xx/=prime[j],ii++; if(ii>1) ng++; } } if(ng) continue; if(xx>1) nd=-nd; ret=(ret+mo+nd*(modpow(hi-lo+1,N)-(hi-lo+1)))%mo; } if(low==1) ret++; ret %= mo; return (int)ret; }
