bzoj 4555 求和
题目大意:
求$\sum\limits_{i=0}^n \sum\limits_{j=0}^i S2(i,j) \times 2^j \times j!$
思路:
法1:
首先把斯特林数展开$S2(i,j)=\frac{1}{j!} \sum\limits_{k=0}^j (-1)^k * \binom{j}{k} * (j-k)^i$
然后先交换求和顺序提出只与$j$有关的项,同时把$j$的范围可以扩展到$n$
代入得到:$\sum\limits_{j=0}^n 2^j * j! \sum\limits_{i=0}^n \sum\limits_{k=0}^j \frac{(-1)^k}{j!} * \binom{j}{k} * (j-k)^i$
拆开组合数:$\sum\limits_{j=0}^n 2^j * j! \sum\limits_{i=0}^n \sum\limits_{k=0}^j \frac{(-1)^k}{k!} * \frac{(j-k)^i}{(j-k)!}$
把只于$k$有关的提出来:$\sum\limits_{j=0}^n 2^j * j! \sum\limits_{k=0}^j \frac{(-1)^k}{k!} * \frac{\sum\limits_{i=0}^n (j-k)^i}{(j-k)!}$
设$f(i)=\frac {(-1)^i}{i!},g(i)=\frac{\sum\limits_{j=0}^n i^j}{i!}$
则原式为$\sum\limits_{i=0}^n 2^i \times (i!) \times (f*g)(i)$
![](https://images.cnblogs.com/OutliningIndicators/ContractedBlock.gif)
1 #include<iostream> 2 #include<cstdio> 3 #include<cstring> 4 #include<cstdlib> 5 #include<cmath> 6 #include<algorithm> 7 #include<queue> 8 #include<vector> 9 #include<map> 10 #include<set> 11 #define ll long long 12 #define db double 13 #define inf 2139062143 14 #define MAXN 200100 15 #define MOD 998244353 16 #define rep(i,s,t) for(register int i=(s),i##__end=(t);i<=i##__end;++i) 17 #define dwn(i,s,t) for(register int i=(s),i##__end=(t);i>=i##__end;--i) 18 #define ren for(register int i=fst[x];i;i=nxt[i]) 19 #define pb(i,x) vec[i].push_back(x) 20 #define pls(a,b) ((a+b)%MOD+MOD)%MOD 21 #define mns(a,b) ((a%MOD-(b)%MOD)%MOD+MOD)%MOD 22 #define mul(a,b) (1LL*(a)*(b))%MOD 23 using namespace std; 24 inline int read() 25 { 26 int x=0,f=1;char ch=getchar(); 27 while(!isdigit(ch)) {if(ch=='-') f=-1;ch=getchar();} 28 while(isdigit(ch)) {x=x*10+ch-'0';ch=getchar();} 29 return x*f; 30 } 31 int n,a[MAXN<<1],b[MAXN<<1],fac[MAXN],ifac[MAXN]; 32 int rev[MAXN<<1],lim,lg,pw[30],ipw[30],ans; 33 int q_pow(int bas,int t,int res=1) 34 { 35 for(;t;t>>=1,bas=mul(bas,bas)) if(t&1) res=mul(res,bas);return res; 36 } 37 #define inv(x) q_pow(x,MOD-2) 38 void ntt(int *a,int n,int f) 39 { 40 rep(i,0,n-1) if(i<rev[i]) swap(a[i],a[rev[i]]); 41 for(int i=1,t=1;i<n;i<<=1,++t) 42 { 43 int wn= f>0?pw[t]:ipw[t];for(int j=0;j<n;j+=i<<1) 44 { 45 int w=1,x,y;for(int k=0;k<i;++k,w=mul(w,wn)) 46 x=a[j+k],y=mul(a[j+k+i],w),a[j+k]=pls(x,y),a[j+k+i]=mns(x,y); 47 } 48 } 49 if(f>0) return ;int nv=inv(n);rep(i,0,n-1) a[i]=mul(a[i],nv); 50 } 51 int main() 52 { 53 n=read();fac[0]=ifac[0]=1;rep(i,1,n) fac[i]=mul(fac[i-1],i),ifac[i]=inv(fac[i]); 54 int tmp=-1;rep(i,0,n) tmp=(MOD-tmp)%MOD,a[i]=mul(tmp,ifac[i]); 55 b[0]=1,b[1]=n+1;rep(i,2,n) b[i]=mul(ifac[i],mul(inv(i-1),(q_pow(i,n+1)-1))); 56 for(lim=1,lg=1;lim<=(n+1<<1);lim<<=1,lg++) 57 pw[lg]=q_pow(3,(MOD-1)/(1<<lg)),ipw[lg]=inv(pw[lg]); 58 rep(i,1,lim-1) rev[i]=(rev[i>>1]>>1)|((i&1)<<lg-2); 59 ntt(a,lim,1);ntt(b,lim,1);rep(i,0,lim-1) a[i]=mul(a[i],b[i]);ntt(a,lim,-1); 60 tmp=1;rep(i,0,n) ans=pls(ans,mul(mul(a[i],tmp),fac[i])),tmp=pls(tmp,tmp); 61 printf("%d\n",ans); 62 }
法2:
思考式子的组合意义,考虑$f(i)=\sum\limits_{j=0}^i S2(i,j) \times 2^j \times j!$
即表示将$i$个数分成$j$个有序集合且每个集合有两种状态的所有方案
那么我们可以枚举最后一个集合中数的数量即$f(n)=\sum\limits_{i=1}^n 2 * \binom{n}{i}* f(n-i)$
两边同除以阶乘,然后就可以分治或者求逆了
![](https://images.cnblogs.com/OutliningIndicators/ContractedBlock.gif)
1 #include<iostream> 2 #include<cstdio> 3 #include<cstring> 4 #include<cstdlib> 5 #include<cmath> 6 #include<algorithm> 7 #include<queue> 8 #include<vector> 9 #include<map> 10 #include<set> 11 #define ll long long 12 #define db double 13 #define inf 2139062143 14 #define MAXN 200100 15 #define MOD 998244353 16 #define rep(i,s,t) for(register int i=(s),i##__end=(t);i<=i##__end;++i) 17 #define dwn(i,s,t) for(register int i=(s),i##__end=(t);i>=i##__end;--i) 18 #define ren for(register int i=fst[x];i;i=nxt[i]) 19 #define pb(i,x) vec[i].push_back(x) 20 #define pls(a,b) ((a+b)%MOD+MOD)%MOD 21 #define mns(a,b) ((a%MOD-(b)%MOD)%MOD+MOD)%MOD 22 #define mul(a,b) (1LL*(a)*(b))%MOD 23 using namespace std; 24 inline int read() 25 { 26 int x=0,f=1;char ch=getchar(); 27 while(!isdigit(ch)) {if(ch=='-') f=-1;ch=getchar();} 28 while(isdigit(ch)) {x=x*10+ch-'0';ch=getchar();} 29 return x*f; 30 } 31 int n,a[MAXN<<1],b[MAXN<<1],fac[MAXN],ifac[MAXN],iv[MAXN<<1],tmp[MAXN<<1]; 32 int rev[MAXN<<1],lim,lg,pw[30],ipw[30],ans,l2[MAXN<<2]; 33 int q_pow(int bas,int t,int res=1) 34 { 35 for(;t;t>>=1,bas=mul(bas,bas)) if(t&1) res=mul(res,bas);return res; 36 } 37 #define inv(x) q_pow(x,MOD-2) 38 void ntt(int *a,int n,int f) 39 { 40 rep(i,0,n-1) if(i<rev[i]) swap(a[i],a[rev[i]]); 41 for(int i=1,t=1;i<n;i<<=1,++t) 42 { 43 int wn= f>0?pw[t]:ipw[t];for(int j=0;j<n;j+=i<<1) 44 { 45 int w=1,x,y;for(int k=0;k<i;++k,w=mul(w,wn)) 46 x=a[j+k],y=mul(a[j+k+i],w),a[j+k]=pls(x,y),a[j+k+i]=mns(x,y); 47 } 48 } 49 if(f>0) return ;int nv=inv(n);rep(i,0,n-1) a[i]=mul(a[i],nv); 50 } 51 void solve(int *a,int *b,int lmt) 52 { 53 rep(i,1,lmt-1) rev[i]=(rev[i>>1]>>1)|((i&1)<<l2[lmt]-1); 54 ntt(a,lmt,1);ntt(b,lmt,1); 55 rep(i,0,lmt-1) a[i]=mul(mns(2,mul(a[i],b[i])),a[i]);ntt(a,lmt,-1); 56 } 57 void Inv(int *a,int *f,int lmt) 58 { 59 f[0]=inv(a[0]);for(int t=2;t<=lmt;t<<=1) 60 {rep(i,0,t-1)tmp[i]=a[i];solve(f,tmp,t<<1);rep(i,t,(t<<1)-1) f[i]=0;} 61 } 62 int main() 63 { 64 n=read();fac[0]=ifac[0]=1;rep(i,1,n) fac[i]=mul(fac[i-1],i),ifac[i]=inv(fac[i]); 65 for(lim=1,lg=1;lim<=(n+1)<<1;lim<<=1,lg++) 66 l2[1<<lg]=lg,pw[lg]=q_pow(3,(MOD-1)/(1<<lg)),ipw[lg]=inv(pw[lg]); 67 lim>>=1;rep(i,1,n) b[i]=mns(MOD,mul(2,ifac[i]));b[0]=pls(b[0],1); 68 Inv(b,a,lim);rep(i,0,n) ans=pls(ans,mul(a[i],fac[i]));printf("%d\n",ans); 69 }