HUNAN 11562 The Triangle Division of the Convex Polygon(大卡特兰数)
http://acm.hunnu.edu.cn/online/?action=problem&type=show&id=11562&courseid=0
求n边形分解成三角形的方案数。
就是求n-2个卡特兰数,从大神那盗取了一份模板,效率极高.同时也很复杂.
1 #include <cstdio> 2 #include <cmath> 3 #include <stdlib.h> 4 #include <memory.h> 5 typedef int typec; 6 typec GCD(typec a, typec b) 7 { 8 return b? GCD(b, a % b) : a; 9 } 10 typec extendGCD(typec a, typec b, typec& x, typec& y) 11 { 12 if(!b) return x = 1, y = 0, a; 13 typec res = extendGCD(b, a % b, x, y), tmp = x; 14 x = y, y = tmp -(a/b)*y; 15 return res; 16 } 17 typec power(typec x, typec k) 18 { 19 typec res = 1; 20 while(k) 21 { 22 if(k&1) res *= x; 23 x *= x, k >>= 1; 24 } 25 return res; 26 } 27 typec powerMod(typec x, typec k, typec m) 28 { 29 typec res = 1; 30 while(x %= m, k) 31 { 32 if(k&1) res *= x, res %= m; 33 x *= x, k >>= 1; 34 } 35 return res; 36 } 37 typec inverse(typec a, typec p, typec t = 1) 38 { 39 typec pt = power(p, t); 40 typec x, y; 41 y = extendGCD(a, pt, x, y); 42 return x < 0? x += pt : x; 43 } 44 typec linearCongruence(typec a, typec b, typec p, typec q) 45 { 46 typec x, y; 47 y = extendGCD(p, q, x, y); 48 x *= b - a, x = p * x + a, x %= p * q; 49 if(x < 0) x += p * q; 50 return x; 51 } 52 const int PRIMEMAX = 1000; 53 int prime[PRIMEMAX + 1]; 54 int getPrime() 55 { 56 memset(prime, 0, sizeof(int) * (PRIMEMAX + 1)); 57 for(int i = 2; i <= PRIMEMAX; i++) 58 { 59 if(!prime[i]) prime[++prime[0]] = i; 60 for(int j = 1; j <= prime[0] && prime[j] <= PRIMEMAX/i; j++) 61 { 62 prime[prime[j]*i] = 1; 63 if(i % prime[j] == 0) break; 64 } 65 } 66 return prime[0]; 67 } 68 int factor[100][3], facCnt; 69 int getFactors(int x) 70 { 71 facCnt = 0; 72 int tmp = x; 73 for(int i = 1; prime[i] <= tmp / prime[i]; i++) 74 { 75 factor[facCnt][1] = 1, factor[facCnt][2] = 0; 76 if(tmp % prime[i] == 0) 77 factor[facCnt][0] = prime[i]; 78 while(tmp % prime[i] == 0) 79 factor[facCnt][2]++, factor[facCnt][1] *= prime[i], tmp /= prime[i]; 80 if(factor[facCnt][2]) facCnt++; 81 } 82 if(tmp != 1) factor[facCnt][0] = tmp, factor[facCnt][1] = tmp, factor[facCnt++][2] = 1; 83 return facCnt; 84 } 85 typec combinationModPt(typec n, typec k, typec p, typec t = 1) 86 { 87 if(k > n) return 0; 88 if(n - k < k) k = n - k; 89 typec pt = power(p, t); 90 typec a = 1, b = k + 1, x, y; 91 int pcnt = 0; 92 while(b % p == 0) pcnt--, b /= p; 93 b %= pt; 94 for(int i = 1; i <= k; i++) 95 { 96 x = n - i + 1, y = i; 97 while(x % p == 0) pcnt++, x /= p; 98 while(y % p == 0) pcnt--, y /= p; 99 x %= pt, y %= pt, a *= x, b *= y; 100 a %= pt, b %= pt; 101 } 102 if(pcnt >= t) return 0; 103 extendGCD(b, pt, x, y); 104 if(x < 0) x += pt; 105 a *= x, a %= pt; 106 return a * power(p, pcnt) % pt; 107 } 108 const typec PTMAX = 45000; 109 typec facmod[PTMAX]; 110 void initFacMod(typec p, typec t = 1) 111 { 112 typec pt = power(p, t); 113 facmod[0] = 1 % pt; 114 for(int i = 1; i < pt; i++) 115 { 116 if(i % p) facmod[i] = facmod[i - 1] * i % pt; 117 else facmod[i] = facmod[i - 1]; 118 } 119 } 120 typec factorialMod(typec n, typec &pcnt, typec p, typec t = 1) 121 { 122 typec pt = power(p, t), res = 1; 123 typec stepCnt = 0; 124 while(n) 125 { 126 res *= facmod[n % pt], res %= pt; 127 stepCnt += n /pt, n /= p, pcnt += n; 128 } 129 res *= powerMod(facmod[pt - 1], stepCnt, pt); 130 return res %= pt; 131 } 132 typec combinationModPtFac(typec n, typec k, typec p, typec t = 1) 133 { 134 if(k > n || p == 1) return 0; 135 if(n - k < k) k = n - k; 136 typec pt = power(p, t), pcnt = 0, pmcnt = 0; 137 if(k < pt) return combinationModPt(n, k, p, t); 138 initFacMod(p, t); 139 typec a = factorialMod(n, pcnt, p, t); 140 typec b = factorialMod(k, pmcnt, p, t); 141 b *= b, pmcnt <<= 1, b %= pt; 142 typec tmp = k + 1; 143 while(tmp % p == 0) tmp /= p, pmcnt++; 144 b *= tmp % pt, b %= pt; 145 pcnt -= pmcnt; 146 if(pcnt >= t) return 0; 147 a *= inverse(b, p, t), a %= pt; 148 return a * power(p, pcnt) % pt; 149 } 150 typec combinationModFac(typec n, typec k, typec m) 151 { 152 getFactors(m); 153 typec a, b, p, q; 154 for(int i = 0; i < facCnt; i++) 155 { 156 if(!i) a = combinationModPtFac(n, k, factor[i][0], factor[i][2]), p = factor[i][1]; 157 else b = combinationModPtFac(n, k, factor[i][0], factor[i][2]), q = factor[i][1]; 158 if(!i) continue; 159 a = linearCongruence(a, b, p, q), p *= q; 160 } 161 return a; 162 } 163 int main() 164 { 165 getPrime(); 166 typec n, k; 167 while(scanf("%d %d", &n, &k) != EOF) 168 printf("%d\n", combinationModFac(2 * (n-2), n-2, k)); 169 return 0; 170 }