Bzoj1951: [Sdoi2010]古代豬文

#題面 傳送門

#Sol 超級數學板子題！！！

費馬小定理，擴展歐幾里德定理，中國剩餘定理，盧卡斯定理等

題意就是求 $$G^{\sum_{k|N}C_N^k}\ mod\ 999911659$$

首先根據擴展歐拉定理或者費馬小定理 $a^x\equiv a^{x% (p-1)}(mod\ p)$ $p$爲質數

那麼直接求那個 $\sum_{k|N}C_N^k (mod\ 999911658)$ 就行了

模數不是質數$lucas+CRT$合併便可

而後$CRT$用上$ExGcd$

# include <bits/stdc++.h>
# define RG register
# define IL inline
# define Fill(a, b) memset(a, b, sizeof(a))
using namespace std;
typedef long long ll;

template <class Int>
IL void Input(RG Int &x){
    RG int z = 1; RG char c = getchar(); x = 0;
    for(; c < '0' || c > '9'; c = getchar()) z = c == '-' ? -1 : 1;
    for(; c >= '0' && c <= '9'; c = getchar()) x = (x << 1) + (x << 3) + (c ^ 48);
    x *= z;
}

const int mod(999911659);
const int phi(999911658);
const int maxn(40000);

int prime[5] = {0, 2, 3, 4679, 35617};
int n, g, fac[5][maxn], inv[5][maxn], ans;

IL int Pow(RG ll x, RG ll y, RG int p){
	RG ll ret = 1;
	for(; y; y >>= 1, x = x * x % p)
		if(y & 1) ret = ret * x % p;
	return ret;
}

IL int ExGcd(RG int a, RG int b, RG int &x, RG int &y){
	if(!b){
		x = 1, y = 0;
		return a;
	}
	RG ll d = ExGcd(b, a % b, y, x);
	y -= 1LL * a / b * x;
	return d;
}

IL int C(RG int x, RG int y, RG int tp){
	return 1LL * fac[tp][x] * inv[tp][y] % prime[tp] * inv[tp][x - y] % prime[tp];
}

IL int Lucas(RG int x, RG int y, RG int tp){
	if(y > x) return 0;
	if(!y) return 1;
	if(x < prime[tp] && y < prime[tp]) return C(x, y, tp);
	RG ll ret = Lucas(x % prime[tp], y % prime[tp], tp);
	return 1LL * ret * Lucas(x / prime[tp], y / prime[tp], tp) % prime[tp];
}

IL int CRT(RG int k){
	RG int x, y, a[5], ret = 0, tmp;
	for(RG int i = 1; i <= 4; ++i){
		a[i] = Lucas(n, k, i);
		ExGcd(tmp = phi / prime[i], prime[i], x, y);
		ret = (ret + 1LL * x * a[i] % phi * tmp % phi) % phi;
	}
	return (ret + phi) % phi;
}

IL void Init(){
	for(RG int j = 1; j <= 4; ++j) fac[j][0] = inv[j][0] = 1;
	for(RG int j = 1; j <= 4; ++j)
		for(RG int i = 1; i < prime[j]; ++i)
			fac[j][i] = 1LL * fac[j][i - 1] * i % prime[j];
	for(RG int j = 1; j <= 4; ++j)
		inv[j][prime[j] - 1] = Pow(fac[j][prime[j] - 1], prime[j] - 2, prime[j]);
	for(RG int j = 1; j <= 4; ++j)
		for(RG int i = prime[j] - 2; i; --i)
			inv[j][i] = 1LL * inv[j][i + 1] * (i + 1) % prime[j];
}

int main(RG int argc, RG char* argv[]){
	Input(n), Input(g);
	if(g % mod == 0) return puts("0"), 0;
	Init();
	for(RG int i = 1; i * i <= n; ++i)
		if(n % i == 0){
			ans = (ans + CRT(i)) % phi;
			RG int j = n / i;
			if(i != j) ans = (ans + CRT(j)) % phi;
		}
	printf("%d\n", Pow(g, ans, mod));
    return 0;
}