求
[g^{sum_{d|n}inom{n}{d}}
]
特判掉 (g,p) 不互素,我们目标是求 (sum_{d|n}inom{n}{d} mod 999911658)。
新模数不是质数,但他是 (2 imes 3 imes 4679 imes 35617),因子没有二次。这样我们就不用像普通的扩展 lucas 一样费死费活地求 (n! mod p_i^{c_i}) 和逆元了。
要求 (inom{n}{d}),只要将其分别对那四个素数用 lucas 求,然后用中国剩余定理合并就是了。
#include <iostream>
#include <cstdio>
#include <vector>
using namespace std;
typedef long long ll;
int n, g, ans, jie[5][40005], ni[5][40005];
const int p=999911659;
const int mod[]={0, 2, 3, 4679, 35617};
vector<int> vec;
int ksm(int a, int b, int p){
int re=1;
while(b){
if(b&1) re = ((ll)re * a) % p;
a = ((ll)a * a) % p;
b >>= 1;
}
return re;
}
void cal(int x){
jie[x][0] = jie[x][1] = ni[x][0] = ni[x][1] = 1;
for(int i=2; i<=mod[x]; i++)
jie[x][i] = jie[x][i-1] * i % mod[x];
for(int i=2; i<=mod[x]; i++)
ni[x][i] = (mod[x] - mod[x] / i) * ni[x][mod[x]%i] % mod[x];
for(int i=2; i<=mod[x]; i++)
ni[x][i] = ni[x][i] * ni[x][i-1] % mod[x];
}
int lucas(int a, int b, int x){
if(a<b) return 0;
else if(a<mod[x]) return (ll)jie[x][a]*ni[x][b]%mod[x]*ni[x][a-b]%mod[x];
else return (ll)lucas(a/mod[x],b/mod[x],x)*lucas(a%mod[x],b%mod[x],x)%mod[x];
}
int C(int n, int m){
int re=0;
for(int i=1; i<=4; i++)
re = (re + (ll)lucas(n,m,i)*((p-1)/mod[i])%(p-1)*ksm((p-1)/mod[i],mod[i]-2,mod[i])) % (p-1);
return re;
}
int main(){
cin>>n>>g;
if(g==p){
cout<<"0
";
return 0;
}
for(int i=1; i*i<=n; i++)
if(n%i==0){
vec.push_back(i);
if(i*i!=n) vec.push_back(n/i);
}
for(int i=1; i<=4; i++)
cal(i);
for(int i=0; i<vec.size(); i++)
ans = (ans + C(n, vec[i])) % (p - 1);
ans = ksm(g, ans, p);
cout<<ans<<endl;
return 0;
}