看到这道题的一刻我是放弃想题的……
看了看暴力中的暴力,发现n和k不超过14的时候可以打个表……
下面是题解
首先可以发现,只要对于$a_i$而言, 二进制下$a_{i-1}$是$a_i$的一个真子集(将二进制位为1的位置当做集合),那么这个序列就是合法的
可以设$f_{i,j}$表示处理完了$a_1$到$a_i$,其中$a_1$按位或到$a_i$有$j$个二进制位是$0$的合法方案数
显然有$f_{i,j}=\sum\limits_{n=0}^{j-1}f_{i-1,n} \times 2^{n} \times C_{k-n}^{j-n}$
相当于枚举一下上一位放好了(包括继承$i-1$之前的)几个二进制位的$0$,这个是$f_{i-1,n}$
当然这些已经固定好的$0$是可以在$a_i$中选填的,所以是$2^{n}$
在填完这$n$个$0$后,还剩下$j-n$个$0$没有填,因此用总的位数减去填了的$0$的个数,得到可以填$0$的个数,即$k-n$,在这其中挑选$j-n$个位置来填上$0$,即$C_{k-n}^{j-n}$
边界是$f_{0,0}=1$,答案是$\sum_{j=1}^{k} f_{n,j}$
这样就得到了一个$O(nk^2)$的算法了 暂时还是没有什么用的
于是我们需要寻求一种更加有普适性的算法,通过上一个式子,可以看出近似于一个多项式的卷积形式了,那么就需要挖掘出一个可以卷积的式子了
定义$g_{i,j}$表示处理完了$a_1$到$a_i$,其中$a_1$按位或到$a_i$后有$j$位是$1$
有$g_{i,j}=\sum_{n=0}^{j-1}g_{i-1,n} \times 2^{n} \times C_{j}^{j-n}$,同时有$f_{i,j}=C_{k}^{j}g_{i,j}$
答案为$\sum_{j=1}^{k}C_{k}^{j}g_{n,j}$
发现这是一个卷积形式,可以用NTT进行优化,时间复杂度为$O(nk \log k)$
重新排列一下变量,得到$g_{i,j}=\sum_{n=1}^{j}g_{i-1,j-n} \times 2^{j-n} \times C_{j}^{n}$
即$$\frac{g_{i,j}}{j!}=\sum_{n=1}^{j}\frac{g_{i-1,j-n} \times 2^{j-n}}{(j-n)!} \times \frac{1}{n!}$$
设$$G_i(x)=\sum_{j=0}^{\infty}\frac{g_{i,j}}{j!}x^j$$
已知$e^x$的泰勒展开式为$e^x=\sum_{i=0}^{\infty}\frac{x^i}{i!}$
有$G_i(x)=G_{i-1}(2x) \times (e^x-1)$(手动展开一下即可)
再将$G_i(x)$手动展开一下,发现$G_i(x)=\prod\limits_{j=0}^{n-1}(e^{2^ix}-1)$,其中只需要保留$x^0$到$x^k$项即可,多余的可以直接抛弃
根据$e^x$的泰勒展开式,可以得知$e^{2^ix}-1=\sum_{n=0}^{\infty} \frac{(2^ix)^n}{n!}=\sum_{n=0}^{\infty} \frac{2^{i \cdot n}}{n!}x^n$
这个式子可以倍增进行求解,类似于快速幂
1 #include <bits/stdc++.h> 2 using namespace std; 3 typedef long long ll; 4 const int p = 998244353, g[2] = { 3, (p + 1) / 3 }, N = 7e4; 5 6 int n, k; 7 8 ll pw(ll a, ll b) { 9 ll r = 1; 10 for( ; b ; b >>= 1, a = a * a % p) if(b & 1) r = r * a % p; 11 return r; 12 } 13 14 ll fac[N], invfac[N]; 15 16 ll C(ll up, ll down) { 17 return up < down ? 0 : fac[up] * invfac[down] % p * invfac[up - down] % p; 18 } 19 20 ll ntt_tmp[N]; 21 22 int rev(int x, int n) { 23 int r = 0; 24 for(int i = 0 ; (1 << i) < n ; ++ i) r = (r << 1) | ((x >> i) & 1); 25 return r; 26 } 27 28 ll ntt(ll *a, int n, int ty) { 29 for(int i = 0 ; i < n ; ++ i) ntt_tmp[rev(i, n)] = a[i]; 30 for(int i = 2 ; i <= n ; i <<= 1) { 31 ll wn = pw(g[ty], (p - 1) / i); 32 for(int j = 0 ; j < n ; j += i) { 33 ll w = 1; 34 for(int k = j ; k < j + i / 2 ; ++ k) { 35 ll u = ntt_tmp[k], v = w * ntt_tmp[k + i / 2] % p; 36 ntt_tmp[k] = (u + v) % p; 37 ntt_tmp[k + i / 2] = ((u - v) % p + p) % p; 38 w = w * wn % p; 39 } 40 } 41 } 42 for(int i = 0 ; i < n ; ++ i) { 43 a[i] = ntt_tmp[i]; 44 if(ty == 1) a[i] = a[i] * pw(n, p - 2) % p; 45 } 46 } 47 48 ll F[N], x[N], y[N]; 49 50 void e2ix(int k, int n) { 51 memset(y, 0, sizeof y); 52 for(int i = 0 ; i < n ; ++ i) 53 y[i] = pw(2, i * k) * x[i] % p; 54 } 55 56 ll _x[N], _y[N]; 57 58 void mul(ll *a, ll *b, int _n) { 59 int n = 1; while(n < (2 * _n)) n <<= 1; 60 for(int i = 0 ; i < _n ; ++ i) _x[i] = a[i], _y[i] = b[i]; 61 for(int i = _n ; i < n ; ++ i) _x[i] = _y[i] = 0; 62 ntt(_x, n, 0), ntt(_y, n, 0); 63 for(int i = 0 ; i < n ; ++ i) _x[i] = (_x[i] * _y[i]) % p; 64 ntt(_x, n, 1); 65 for(int i = 0 ; i < n ; ++ i) a[i] = _x[i]; 66 } 67 68 void sol(int n) { 69 for(int i = 1 ; i <= k ; ++ i) x[i] = invfac[i]; 70 int now = 0, pw2 = 1; 71 F[0] = 1; 72 for( ; n ; n >>= 1, pw2 <<= 1) { 73 if(n & 1) { 74 e2ix(now, k + 1); 75 mul(F, y, k + 1); 76 now |= pw2; 77 } 78 e2ix(pw2, k + 1); 79 mul(x, y, k + 1); 80 } 81 } 82 83 int main() { 84 freopen("or.in", "r", stdin); 85 freopen("or.out", "w", stdout); 86 scanf("%d%d", &n, &k); 87 fac[0] = invfac[0] = 1; 88 for(int i = 1 ; i <= k ; ++ i) { 89 fac[i] = fac[i - 1] * i % p; 90 invfac[i] = pw(fac[i], p - 2); 91 } 92 sol(n); 93 ll ans = 0; 94 for(int i = 1 ; i <= k ; ++ i) { 95 ans = (ans + F[i] * C(k, i) % p * fac[i] % p) % p; 96 } 97 printf("%lld\n", ans); 98 }
参考文献:Orange Software [JZOJ 5785 Or]