[sdoi2010]古代猪文

2019-04-13 20:48发布

站内文章 / 模拟电子
16648 0 1145
写这道题首先srO qw。。。。题目就不贴了。。。大意：求 (G^(Σd|n C(n, d))) mod M。 首先发现模是一个质数，所以根据费马小定理可知指数只需模模减一即可。 然后又发现模减一太大，而且不是质数。。。然后发现质因子分解后都是比较小的数字。。果断中国剩余定理合并啊。。。 然后记得lyp讲过什么n的正约数个数均摊是O(logn)级的。。。果断暴枚约数啊。。。 然后是组合数取模。。。n太大，实在想不出什么好办法直接处理阶乘，好像可以用勒让德定理，但是觉得太麻烦了。。。于是用zfone讲的lucas定理做吧。。。 然后组合数就线性预处理阶乘的模和逆元，每次直接用好了。。。 然后时间就被虐爆了。。。不知道好多版去了。。。 代码写了一百多行。。。各种繁琐。。。脑残错误有一开始搞错了模。。。 #include #include #include #include #include #include #include #include #include #include #include #include #include using namespace std; #ifdef WIN32 #define fmt64 "%I64d" #else #define fmt64 "%lld" #endif #define PI M_PI #define oo 0x13131313 #define iter iterator #define PB push_back #define PO pop_back #define MP make_pair #define fst first #define snd second #define cstr(a) (a).c_str() #define FOR(i, j, k) for(i = (j); i <= (k); ++i) #define ROF(i, j, k) for(i = (j); i >= (k); --i) #define FER(e, d, u) for(e = d[u]; e; e = e->n) #define FRE(i, a) for(i = (a).begin(); i != (a).end(); ++i) typedef unsigned int uint; typedef long long int64; typedef unsigned long long uint64; typedef long double real; template inline bool minim(T &a, const T &b) {return b < a ? a = b, 1 : 0;} template inline bool maxim(T &a, const T &b) {return b > a ? a = b, 1 : 0;} template inline T sqr(const T &a) {return a * a;} #define maxn 100005 #define maxm 35700 int pt, pm[maxn / 5]; bool np[maxn + 1]; int ft, fac[20], cnt[20]; int dt, d[maxn]; void get_prime() { int i, j, k; FOR(i, 2, maxn) { if (!np[i]) pm[++pt] = i; FOR(j, 1, pt) { if ((k = i * pm[j]) > maxn) break; np[k] = 1; if (!(i % pm[j])) break; } } } void divide(int k) { int i, j; j = ft = 0; FOR(i, 1, pt) { if (sqr(pm[i]) > k) break; if (!(k % pm[i])) { fac[++ft] = pm[i], cnt[ft] = 1; for (k /= pm[i]; !(k % pm[i]); k /= pm[i]) ++cnt[ft]; } } if (k != 1) fac[++ft] = k, cnt[ft] = 1; } void dfs(int t, int rest, int now) { if (!rest) {d[++dt] = now; return;} if (t > ft) return; dfs(t + 1, rest, now); for (int i = 1; i <= cnt[t]; ++i) { if (i > rest) break; dfs(t + 1, rest - i, now *= fac[t]); } } const int mm = 2 * 3 * 4679 * 35617; const int m[4] = {2, 3, 4679, 35617}; const int M[4] = {499955829, 333303886, 213702, 28074}; int inv[4]; int f[4][maxm + 1], v[4][maxm + 1]; int ex_gcd(int a, int b, int &x, int &y) { if (!b) return x = 1, y = 0, a; int ret = ex_gcd(b, a % b, x, y), t; return t = x, x = y, y = t - (a / b) * y, ret; } int inverse(int a, int b) { int x, y; ex_gcd(a, b, x, y); if ((x %= b) < 0) x += b; return x; } void prepare() { int i, j, k, l; FOR(l, 0, 3) { inv[l] = inverse(M[l], m[l]); f[l][0] = f[l][1] = v[l][0] = v[l][1] = 1; } FOR(l, 0, 3) FOR(i, 2, m[l] - 1) { f[l][i] = (f[l][i - 1] * i) % m[l]; if (!np[i]) v[l][i] = inverse(i, m[l]); FOR(j, 1, pt) { if ((k = pm[j] * i) > m[l]) break; (v[l][k] = v[l][i] * v[l][pm[j]]) %= m[l]; if (!(i % pm[j])) break; } } FOR(l, 0, 3) FOR(i, 2, m[l] - 1) (v[l][i] *= v[l][i - 1]) %= m[l]; } int C(int a, int b, int l) { return a < b ? 0 : (int64(f[l][a] * v[l][b]) * v[l][a - b]) % m[l]; } int lucas(int a, int b, int l) { int64 ret = 1; for (; b; a /= m[l], b /= m[l]) (ret *= C(a % m[l], b % m[l], l)) %= m[l]; return ret; } int main() { freopen("ancient.in", "r", stdin); freopen("ancient.out", "w", stdout); int i, j; int64 k; j = k = 0; int n; int64 G; cin >> n >> G; if (G == mm + 1) puts("0"), exit(0); get_prime(), divide(n); FOR(i, 1, ft) j += cnt[i]; FOR(i, 0, j) dfs(1, i, 1); prepare(); FOR(i, 1, dt) { int64 sum = 0; FOR(j, 0, 3) { int64 t = lucas(n, d[i], j); (sum += inv[j] * M[j] * t) %= mm; } (k += sum) %= mm; } int64 ans = 1; for (i = k; i; (G *= G) %= mm + 1, i >>= 1) if (i & 1) (ans *= G) %= mm + 1; cout << ans << endl; return 0; }

热门文章