UOJ#449 喂鸽子

题意：有n个鸽子，你每秒随机喂一只鸽子，每只鸽子吃k次就饱了。求期望多少秒之后全饱了。n <= 50, k <= 1000。

解：有两种做法。一种直接DP的n²k做法在这。我用的是Min-Max容斥 + NTT优化DP。

先Min-Max容斥，由于鸽子是等价的，所以相当于求m只鸽子期望多少秒之后有一只饱了。

讲不清楚...看题解吧。

注意题解中m个盒子的答案为什么放了i + 1个球，这是因为我们把最后一个球拿出来考虑了。所以唯一放满k个球的盒子的贡献是1 / (k - 1)!

  1 #include <bits/stdc++.h>
  2 
  3 const int MO = 998244353;
  4 const int N = 50010;
  5 
  6 int n, K, fac[N], inv[N], invn[N];
  7 int f[55][N], g[55][N], A[N * 4], B[N * 4], r[N * 4];
  8 
  9 inline int C(int n, int m) {
 10     if(n < 0 || m < 0 || n < m) return 0;
 11     return 1ll * fac[n] * invn[m] % MO * invn[n - m] % MO;
 12 }
 13 
 14 inline int qpow(int a, int b) {
 15     int ans = 1;
 16     while(b) {
 17         if(b & 1) ans = 1ll * ans * a % MO;
 18         a = 1ll * a * a % MO;
 19         b = b >> 1;
 20     }
 21     return ans;
 22 }
 23 
 24 inline int Inv(int x) {
 25     return qpow(x, MO - 2);
 26 }
 27 
 28 inline void prework(int n) {
 29     static int R = 0;
 30     if(R == n) return;
 31     R = n;
 32     int lm = 1;
 33     while((1 << lm) < n) lm++;
 34     for(int i = 0; i < n; i++) {
 35         r[i] = (r[i >> 1] >> 1) | ((i & 1) << (lm - 1));
 36     }
 37     return;
 38 }
 39 
 40 inline void NTT(int *a, int n, int f) {
 41     prework(n);
 42     for(int i = 0; i < n; i++) {
 43         if(i < r[i]) {
 44             std::swap(a[i], a[r[i]]);
 45         }
 46     }
 47     for(int len = 1; len < n; len <<= 1) {
 48         int Wn = qpow(3, (MO - 1) / (len << 1));
 49         if(f == -1) Wn = Inv(Wn);
 50         for(int i = 0; i < n; i += (len << 1)) {
 51             int w = 1;
 52             for(int j = 0; j < len; j++) {
 53                 int t = 1ll * a[i + len + j] * w % MO;
 54                 a[i + len + j] = (a[i + j] - t) % MO;
 55                 a[i + j] = (a[i + j] + t) % MO;
 56                 w = 1ll * w * Wn % MO;
 57             }
 58         }
 59     }
 60     if(f == -1) {
 61         int inv = Inv(n);
 62         for(int i = 0; i < n; i++) {
 63             a[i] = 1ll * a[i] * inv % MO;
 64         }
 65     }
 66     return;
 67 }
 68 
 69 int main() {
 70     
 71     scanf("%d%d", &n, &K);
 72     
 73     fac[0] = inv[0] = invn[0] = 1;
 74     fac[1] = inv[1] = invn[1] = 1;
 75     for(int i = 2; i <= n * K; i++) {
 76         fac[i] = 1ll * fac[i - 1] * i % MO;
 77         inv[i] = 1ll * inv[MO % i] * (MO - MO / i) % MO;
 78         invn[i] = 1ll * invn[i - 1] * inv[i] % MO;
 79     }
 80     
 81     int len = 1;
 82     while(len <= n * K) len <<= 1;
 83     g[0][0] = 1;
 84     memcpy(B, invn, K * sizeof(int));
 85     NTT(B, len, 1);
 86     for(int i = 1; i <= n; i++) {
 87         /*for(int j = 0; j <= i * K; j++) {
 88             /// g[i][j] f[i][j] 
 89             for(int k = 0; k < K && k <= j; k++) {
 90                 g[i][j] = (g[i][j] + 1ll * g[i - 1][j - k] * invn[k] % MO) % MO;
 91                 f[i][j] = (f[i][j] + 1ll * f[i - 1][j - k] * invn[k] % MO) % MO;
 92             }
 93             if(j >= K) f[i][j] = (f[i][j] + 1ll * g[i - 1][j - K] * invn[K - 1] % MO) % MO;
 94         }*/
 95         memcpy(A, g[i - 1], n * K * sizeof(int));
 96         memset(A + n * K, 0, (len - n * K) * sizeof(int));
 97         NTT(A, len, 1);
 98         for(int j = 0; j < len; j++) {
 99             A[j] = 1ll * A[j] * B[j] % MO;
100         }
101         NTT(A, len, -1);
102         for(int j = 0; j <= i * K; j++) {
103             g[i][j] = A[j];
104         }
105         memcpy(A, f[i - 1], n * K * sizeof(int));
106         memset(A + n * K, 0, (len - n * K) * sizeof(int));
107         NTT(A, len, 1);
108         for(int j = 0; j < len; j++) {
109             A[j] = 1ll * A[j] * B[j] % MO;
110         }
111         NTT(A, len, -1);
112         for(int j = 0; j <= i * K; j++) {
113             f[i][j] = A[j];
114             if(j >= K) {
115                 f[i][j] = (f[i][j] + 1ll * g[i - 1][j - K] * invn[K - 1] % MO) % MO;
116             }
117         }
118     }
119     
120     int ans = 0;
121     for(int i = 1; i <= n; i++) {
122         int temp = 0;
123         for(int j = K - 1; j <= i * K; j++) {
124             temp = (temp + 1ll * fac[j - 1] * f[i][j] % MO * qpow(inv[i], j) % MO * j % MO) % MO;
125         }
126         temp = 1ll * temp * n % MO * inv[i] % MO * C(n, i) % MO;
127         if(i & 1) {
128             ans = (ans + temp) % MO;
129         }
130         else {
131             ans = (ans - temp) % MO;
132         }
133     }
134     printf("%d\n", (ans + MO) % MO);
135     return 0;
136 }

AC代码

 

转载于:https://www.cnblogs.com/huyufeifei/p/10942703.html