LG P2480 [SDOI2010] 古代猪文 Solution
Description
给定 n,gn,gn,g，定义 f(n)=∑k∣n(nk)f(n)=\sum\limits_{k\mid n}\binom{n}{k}f(n)=k∣n∑(kn)，求 gf(n) mod 999911659g^{f(n)}\bmod 999911659gf(n)mod999911659 的值.
Limitations
1≤n,g≤1091\le n,g\le 10^91≤n,g≤109
1s,125MB1\text{s},125\text{MB}1s,125MB
Solution
首先，如果 g=999911659g=999911659g=999911659，直接特判掉输出 000.
接下来，根据欧拉定理，我们需要求 f(n) mod 999911658f(n)\bmod 999911658f(n)mod999911658.
虽然可以暴力枚举 n\sqrt nn 个因数求答案，但是 999911658∉P999911658\notin \mathbb{P}999911658∈/P，在求组合数时无法直接使用 Lucas 定理.
注意到 999911658=2×3×4679×35617999911658=2\times 3\times 4679\times 35617999911658=2×3×4679×35617，无平方因子.
现在我们就可以用 Lucas 分别求出 f(n)f(n)f(n) 模这四个数的结果 a1,a2,a3,a4a_1,a_2,a_3,a_4a1,a2,a3,a4，然后用 CRT 解如下方程：
{x≡a1(mod2)x≡a2(mod3)x≡a3(mod4679)x≡a4(mod35617)
\begin{equation*}
\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}
\end{equation*}
⎩⎨⎧x≡a1(mod2)x≡a2(mod3)x≡a3(mod4679)x≡a4(mod35617)
即可得到 f(n) mod 999911658f(n)\bmod 999911658f(n)mod999911658，再用一次快速幂即可求得答案.
Code
2.14KB,94ms,864KB (C++20 with O2)2.14\text{KB},94\text{ms},864\text{KB}\;\texttt{(C++20 with O2)}2.14KB,94ms,864KB(C++20 with O2)
#include <bits/stdc++.h>
using namespace std;using i64 = long long;
using ui64 = unsigned long long;
using i128 = __int128;
using ui128 = unsigned __int128;
using f4 = float;
using f8 = double;
using f16 = long double;template<class T>
bool chmax(T &a, const T &b){if(a < b){ a = b; return true; }return false;
}template<class T>
bool chmin(T &a, const T &b){if(a > b){ a = b; return true; }return false;
}constexpr int P = 999911659;inline int add(int a, int b, int mod) { return (a + b >= mod) ? (a + b - mod) : (a + b); }
inline int mul(int a, int b, int mod) { return 1LL * a * b % mod; }
inline int qpow(int a, int b, int mod) {int res = 1;for (; b; b >>= 1, a = mul(a, a, mod))if (b & 1) res = mul(res, a, mod);return res;
}struct Lucas {int P;vector<int> fac, ifac;inline Lucas() {}inline void init(int p) {P = p;fac.resize(P);ifac.resize(P);fac[0] = 1;for (int i = 1; i < P; ++i) {fac[i] = mul(fac[i - 1], i, P);}ifac[P - 1] = qpow(fac[P - 1], P - 2, P);for (int i = P - 2; i >= 0; --i) {ifac[i] = mul(ifac[i + 1], i + 1, P);}}inline int comb(int n, int m) {if (m < 0 || m > n) return 0;return mul(fac[n], mul(ifac[m], ifac[n - m], P), P);}inline int lucas(int n, int m) {if (m == 0) return 1;return mul(comb(n % P, m % P), lucas(n / P, m / P), P);}
};signed main() {ios::sync_with_stdio(0);cin.tie(0), cout.tie(0);int n, g; cin >> n >> g;if (g == P) {cout << 0;return 0;}array<Lucas, 4> L;L[0].init(2), L[1].init(3);L[2].init(4679), L[3].init(35617);array<int, 4> sum; sum.fill(0);for (int i = 1; i <= n / i; i++) if (n % i == 0)for (int j = 0; j < 4; j++) {sum[j] = add(sum[j], L[j].lucas(n, i), L[j].P); if (i != n / i) sum[j] = add(sum[j], L[j].lucas(n, n / i), L[j].P);}int ans = 0;for (int i = 0; i < 4; i++) {int M = (P - 1) / L[i].P, Minv = qpow(M, L[i].P - 2, L[i].P);ans = add(ans, mul(mul(M, Minv, P - 1), sum[i], P - 1), P - 1);}cout << qpow(g, ans, P);return 0;
}