1. 程式人生 > >洛谷P2480 古代豬文

洛谷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程式碼