洛谷P2480 古代猪文

这道题把我坑了好久......

原因竟是CRT忘了取正数!

题意:求

指数太大了,首先用欧拉定理取模。

由于模数是质数所以不用加上phi(p)

然后发现phi(p)过大,不能lucas,但是它是个square free,可以分解质因数然后lucas然后CRT。

然后就没有然后了......模板套来套去......

注意CRT的结果可能是负数,要取正。

  1 #include <cstdio>
  2 #include <algorithm>
  3 #define say(a) printf(#a); printf(" = %lld \n", a)
  4 
  5 typedef long long LL;
  6 const int N = 200010;
  7 const LL MO = 999911659;
  8 const LL prime[4] = {2, 3, 4679, 35617};
  9 
 10 int T;
 11 LL n, G;
 12 LL nn[4][N], inv[4][N], invn[4][N], aa[4];
 13 
 14 inline LL qpow(LL a, LL b, LL c) {
 15     LL ans = 1;
 16     while(b) {
 17         if(b & 1) {
 18             ans = ans * a % c;
 19         }
 20         a = a * a % c;
 21         b = b >> 1;
 22     }
 23     return ans;
 24 }
 25 LL C(int a, int b) {
 26     return nn[T][a] * invn[T][b] % prime[T] * invn[T][a - b] % prime[T];
 27 }
 28 LL lucas(int a, int b) {
 29     if(!b) {
 30         return 1;
 31     }
 32     return C(a % prime[T], b % prime[T]) * lucas(a / prime[T], b / prime[T]) % prime[T];
 33 }
 34 LL exgcd(LL a, LL b, LL &x, LL &y) {
 35     if(!b) {
 36         x = 1;
 37         y = 0;
 38         return a;
 39     }
 40     LL g = exgcd(b, a % b, x, y);
 41     std::swap(x, y);
 42     y -= (a / b) * x;
 43     return g;
 44 }
 45 
 46 inline void prework(int h) {
 47     inv[h][1] = 1;
 48     LL p = prime[h];
 49     for(int i = 2; i < p; i++) {
 50         inv[h][i] = (p - p / i) * inv[h][p % i] % p;
 51     }
 52     nn[h][0] = 1;
 53     invn[h][0] = 1;
 54     for(int i = 1; i < p; i++) {
 55         nn[h][i] = nn[h][i - 1] * i % p;
 56         invn[h][i] = invn[h][i - 1] * inv[h][i] % p;
 57     }
 58     return;
 59 }
 60 
 61 inline LL CRT() {
 62 
 63     LL M = MO - 1, ans = 0, b, t;
 64     for(int i = 0; i < 4; i++) {
 65         LL m = M / prime[i];
 66         exgcd(m, prime[i], t, b);
 67         t = (t % M + M) % M;
 68         ans += aa[i] * m % M * t % M;
 69         ans %= M;
 70     }
 71 
 72     return ans;
 73 }
 74 
 75 int main() {
 76 
 77     scanf("%lld%lld", &n, &G);
 78     if(G == MO) {
 79         printf("0");
 80         return 0;
 81     }
 82     for(int i = 0; i < 4; i++) {
 83         prework(i);
 84     }
 85 
 86     for(T = 0; T < 4; T++) {
 87         for(int i = 1; i * i <= n; i++) {
 88             if(n % i) {
 89                 continue;
 90             }
 91             aa[T] += lucas(n, i);
 92             aa[T] %= prime[T];
 93             if(i * i < n) {
 94                 aa[T] += lucas(n, n / i);
 95                 aa[T] %= prime[T];
 96             }
 97         }
 98         //printf("aa[%d] = %lld \n", T, aa[T]);
 99     }
100 
101     LL ans = CRT();
102     //printf("ans = %lld \n", ans);
103 
104 
105     ans = qpow(G, ans, MO);
106     printf("%lld\n", ans);
107     return 0;
108 }
AC代码

 

posted @ 2018-11-29 20:17  huyufeifei  阅读(172)  评论(0编辑  收藏  举报
试着放一个广告栏(虽然没有一分钱广告费)

『Flyable Heart 応援中!』 HHG 高苗京铃 闪十PSS 双六 電動伝奇堂 章鱼罐头制作组 はきか 祝姬 星降夜