古代猪文

数论大杂烩
题目简意：给定n, g，p表示n的所有因子,求$$g^{\sum_{p|n} C_n^p} \mod 999911659$$
（markdown看不清楚可以放大看）
注意到999911659是质数，所以使用欧拉定理变成$$g^{\sum_{p|n} C_n^p \mod \phi(999911659)} \mod 999911659$$
即$$g^{\sum_{p|n} C_n^p \mod 999911658} \mod 999911659$$
显然重点要算的是指数，把指数算出来之后用快速幂直接算即可；
不过由于999911658非常大，直接用lucas显然炸了，所以分解一下变成：\(2 \times 3 \times 4679 \times 35617\)
枚举n的因子，分解计算在模2，3，4679，35617下的结果，记为\(a_1\), \(a_2\), \(a_3\), \(a_4\), 用excrt 解出x：

\[\begin{cases} x \equiv a_1 \pmod{2} \\ x \equiv a_2 \pmod{3} \\ x \equiv a_3 \pmod{4679} \\ x \equiv a_4 \pmod{35617} \end{cases} \]

最后快速幂计算 \(g^x\);
代码如下：

#include<bits/stdc++.h>
#define int long long
using namespace std;
const int mod = 999911659;
int n, g, a[5], b[5] = {0, 2, 3, 4679, 35617};
int qpow(int a, int b, int p){
	int m = 1;
	for(; b; b >>= 1){
		if(b & 1) m = (m * a) % p;
		a = (a * a) % p;
	}
	return m;
}
int C(int n, int m, int p){
	if(n > m) return 0;
	int a, b, c; a = b = c = 1;
	for(int i = 1; i <= m; i++) a = (a * i) % p;
	for(int i = 1; i <= n; i++) b = (b * i) % p;
	for(int i = 1; i <= m - n; i++) c = (c * i) % p;
	return (a * qpow(b * c % p, p - 2, p) % p + p) % p;
}
int lucas(int n, int m, int p){
	if(n == 0) return 1;
	int c = lucas(n / p, m / p, p);
	return (c * C(n % p, m % p, p) % p + p) % p;
}
int exgcd(int a, int b, int &x, int &y){
	if(b == 0){
		x = 1; y = 0;
		return a;
	}
	int d = exgcd(b, a % b, x, y);
	int z = x;
	x = y;
	y = z - (a / b) * x;
	return d;
}
void excrt(int &a1, int &a2, int &b1, int &b2){
	int yo = 0, ys = 0;
	int a = b1, b = b2;
	int d = exgcd(a, b, yo, ys);
	b1 = a * b / d;
	ys = ys * (a1 - a2) / d;
	a1 = (a2 + b2 * ys % b1 + b1) % b1;
}
signed main(){
	cin >> n >> g;
	if(g % mod == 0){
		cout << 0;
		return 0;
	}
	for(int j = 1; j <= 4; j++){
		for(int i = 1; i * i <= n; i++){
			if(n % i == 0){
				a[j] = (a[j] + lucas(i, n, b[j])) % b[j];
				if(i * i != n){
					a[j] = (a[j] + lucas(n / i, n, b[j])) % b[j];
				}
			}
		}
	}
	for(int i = 2; i <= 4; i++){
		excrt(a[i - 1], a[i], b[i - 1], b[i]);
		a[i] = a[i - 1]; b[i] = b[i - 1];
	}
	cout << qpow(g, a[4], mod);
	return 0;
}