[BZOJ 4555][Tjoi2016&Heoi2016]求和

题意

给定 $n$ , 求下式的值:

$$ f(n)= \sum_{i=0}^n\sum_{j=0}^i\begin{Bmatrix}i\\ j\end{Bmatrix}\times 2^j\times j!$$

题解

这题比较神仙...

那么我们可以思考如何来求一个比较简单的转移式.

首先我们发现, $f(n)$ 表达式中的第一重和式包含了 $f(n-1)$, 那么我们对 $f$ 的值做差分, 于是我们有 $f(n)-f(n-1)=\sum\limits_{i=0}^n\begin{Bmatrix}i\\ j\end{Bmatrix}\times 2^j\times j!$ .设 $g(n)=f(n)-f(n-1)$.

而我们知道第二类斯特林数的组合意义是在 $i$ 个物品中选 $j$ 个子集的方案数. 那么后面的 $j!$ 的组合意义相当于让划分的子集带标号, $2^j$ 的组合意义相当于每个子集贡献两种状态.

知道了表达式的组合意义, 我们可以尝试用朴素DP来解它. 我们枚举一下最后一个子集中的元素个数为 $i$ (因为子集带标号, 我们相当于枚举最后一个), 那么我们就有递推式:

$$g(n)=[n=0]+\sum_{i=1}^n {n\choose i} g(n-i)\times 2$$

最后乘 $2$ 是因为每个集合会贡献两种状态, $[n=0]$ 是边界条件.

观察这个式子, 发现是二项卷积的形式, 我们果断上EGF来解决. 设 $G(x)$ 是 $g(n)$ 的EGF, $H(x)$ 是数列 $\langle 0,2,2,2,\dots \rangle$ 的EGF, 那么我们有:

$$ G(x)=G(x)H(x)+1 $$

注意到 $H(x)$ 代表的数列的第 $0$ 项是 $0$, 因为 $g(n)$ 的递推式中求和下指标是 $1$.

那么我们移项之后就可以开心地得到:

$$ G(x)=\frac 1 {1-H(x)} $$

NTT爆算一发就可以了

代码实现

记得EGF要除以阶乘...出答案的时候还得乘回来...

  1 #include <bits/stdc++.h>
  2 
  3 const int G=3;
  4 const int DFT=1;
  5 const int IDFT=-1;
  6 const int MAXN=550000;
  7 const int MOD=998244353;
  8 const int INV2=(MOD+1)>>1;
  9 const int PHI=MOD-1;
 10 
 11 typedef std::vector<int> Poly;
 12 
 13 Poly Sqrt(Poly);
 14 void Read(Poly&);
 15 Poly Inverse(Poly);
 16 Poly Ln(const Poly&);
 17 Poly Exp(const Poly&);
 18 void Print(const Poly&);
 19 void NTT(Poly&,int,int);
 20 Poly Pow(const Poly&,int);
 21 Poly Integral(const Poly&);
 22 Poly Derivative(const Poly&);
 23 Poly operator*(Poly,Poly);
 24 Poly operator/(Poly,Poly);
 25 Poly operator%(Poly,Poly);
 26 Poly operator+(const Poly&,const Poly&);
 27 Poly operator-(const Poly&,const Poly&);
 28 
 29 int rev[MAXN];
 30 
 31 int NTTPre(int);
 32 int Pow(int,int,int);
 33 
 34 int main(){
 35     int n;
 36     scanf("%d",&n);
 37     Poly h(n+1),one(1,1);
 38     h[1]=2;
 39     for(int i=2;i<=n;i++)
 40         h[i]=1ll*h[i-1]*Pow(i,MOD-2,MOD)%MOD;
 41     Poly g(Inverse(one-h));
 42     int ans=0;
 43     int fact=1;
 44     for(int i=0;i<=n;i++){
 45         (ans+=1ll*g[i]*fact%MOD)%=MOD;
 46         fact=1ll*fact*(i+1)%MOD;
 47     }
 48     printf("%d\n",ans);
 49     return 0;
 50 }
 51 
 52 void Read(Poly& a){
 53     for(auto& i:a)
 54         scanf("%d",&i);
 55 }
 56 
 57 void Print(const Poly& a){
 58     for(auto i:a)
 59         printf("%d ",i);
 60     puts("");
 61 }
 62 
 63 Poly Pow(const Poly& a,int k){
 64     Poly log=Ln(a);
 65     for(auto& i:log)
 66         i=1ll*i*k%MOD;
 67     return Exp(log);
 68 }
 69 
 70 Poly Sqrt(Poly a){
 71     int len=a.size();
 72     if(len==1)
 73         return Poly(1,int(sqrt(a[0])));
 74     else{
 75         Poly b=a;
 76         b.resize((len+1)>>1);
 77         b=Sqrt(b);
 78         b.resize(len);
 79         Poly inv=Inverse(b);
 80         int bln=NTTPre(inv.size()+a.size());
 81         NTT(a,bln,DFT);
 82         NTT(inv,bln,DFT);
 83         for(int i=0;i<bln;i++)
 84             a[i]=1ll*a[i]*INV2%MOD*inv[i]%MOD;
 85         NTT(a,bln,IDFT);
 86         for(int i=0;i<len;i++)
 87             b[i]=(1ll*b[i]*INV2%MOD+a[i])%MOD;
 88         return b;
 89     }
 90 }
 91 
 92 Poly Exp(const Poly& a){
 93     size_t len=1;
 94     Poly ans(1,1),one(1,1);
 95     while(len<(a.size()<<1)){
 96         len<<=1;
 97         Poly b=a;
 98         b.resize(len);
 99         ans=ans*(one-Ln(ans)+b);
100         ans.resize(len);
101     }
102     ans.resize(a.size());
103     return ans;
104 }
105 
106 Poly Ln(const Poly& a){
107     Poly ans=Integral(Derivative(a)*Inverse(a));
108     ans.resize(a.size());
109     return ans;
110 }
111 
112 Poly Integral(const Poly& a){
113     int len=a.size();
114     Poly ans(len+1);
115     for(int i=1;i<len;i++)
116         ans[i]=1ll*a[i-1]*Pow(i,MOD-2,MOD)%MOD;
117     return ans;
118 }
119 
120 Poly Derivative(const Poly& a){
121     int len=a.size();
122     Poly ans(len-1);
123     for(int i=1;i<len;i++)
124         ans[i-1]=1ll*a[i]*i%MOD;
125     return ans;
126 }
127 
128 Poly operator/(Poly a,Poly b){
129     int n=a.size()-1,m=b.size()-1;
130     Poly ans(1);
131     if(n>=m){
132         std::reverse(a.begin(),a.end());
133         std::reverse(b.begin(),b.end());
134         b.resize(n-m+1);
135         ans=Inverse(b)*a;
136         ans.resize(n-m+1);
137         std::reverse(ans.begin(),ans.end());
138     }
139     return ans;
140 }
141 
142 Poly operator%(Poly a,Poly b){
143     int n=a.size()-1,m=b.size()-1;
144     Poly ans;
145     if(n<m)
146         ans=a;
147     else
148         ans=a-(a/b)*b;
149     ans.resize(m);
150     return ans;
151 }
152 
153 Poly operator*(Poly a,Poly b){
154     int len=a.size()+b.size()-1;
155     int bln=NTTPre(len);
156     NTT(a,bln,DFT);
157     NTT(b,bln,DFT);
158     for(int i=0;i<bln;i++)
159         a[i]=1ll*a[i]*b[i]%MOD;
160     NTT(a,bln,IDFT);
161     a.resize(len);
162     return a;
163 }
164 
165 Poly operator+(const Poly& a,const Poly& b){
166     Poly ans(std::max(a.size(),b.size()));
167     std::copy(a.begin(),a.end(),ans.begin());
168     for(size_t i=0;i<b.size();i++)
169         ans[i]=(ans[i]+b[i])%MOD;
170     return ans;
171 }
172 
173 Poly operator-(const Poly& a,const Poly& b){
174     Poly ans(std::max(a.size(),b.size()));
175     std::copy(a.begin(),a.end(),ans.begin());
176     for(size_t i=0;i<b.size();i++)
177         ans[i]=(ans[i]+MOD-b[i])%MOD;
178     return ans;
179 }
180 
181 Poly Inverse(Poly a){
182     int len=a.size();
183     if(len==1)
184         return Poly(1,Pow(a[0],MOD-2,MOD));
185     else{
186         Poly b(a);
187         b.resize((len+1)>>1);
188         b=Inverse(b);
189         int bln=NTTPre(b.size()*2+a.size());
190         NTT(a,bln,DFT);
191         NTT(b,bln,DFT);
192         for(int i=0;i<bln;i++)
193             b[i]=(2ll*b[i]%MOD-1ll*b[i]*b[i]%MOD*a[i]%MOD+MOD)%MOD;
194         NTT(b,bln,IDFT);
195         b.resize(len);
196         return b;
197     }
198 }
199 
200 void NTT(Poly& a,int len,int opt){
201     a.resize(len);
202     for(int i=0;i<len;i++)
203         if(rev[i]>i)
204             std::swap(a[i],a[rev[i]]);
205     for(int i=1;i<len;i<<=1){
206         int step=i<<1;
207         int wn=Pow(G,(PHI+opt*PHI/step)%PHI,MOD);
208         for(int j=0;j<len;j+=step){
209             int w=1;
210             for(int k=0;k<i;k++,w=1ll*w*wn%MOD){
211                 int x=a[j+k];
212                 int y=1ll*a[j+k+i]*w%MOD;
213                 a[j+k]=(x+y)%MOD;
214                 a[j+k+i]=(x-y+MOD)%MOD;
215             }
216         }
217     }
218     if(opt==IDFT){
219         int inv=Pow(len,MOD-2,MOD);
220         for(int i=0;i<len;i++)
221             a[i]=1ll*a[i]*inv%MOD;
222     }
223 }
224 
225 int NTTPre(int n){
226     int bln=1,bct=0;
227     while(bln<n){
228         bln<<=1;
229         ++bct;
230     }
231     for(int i=0;i<bln;i++)
232         rev[i]=(rev[i>>1]>>1)|((i&1)<<(bct-1));
233     return bln;
234 }
235 
236 inline int Pow(int a,int n,int p){
237     int ans=1;
238     while(n>0){
239         if(n&1)
240             ans=1ll*a*ans%p;
241         a=1ll*a*a%p;
242         n>>=1;
243     }
244     return ans;
245 }
BZOJ 4555

 

日常图包

 

posted @ 2019-01-04 14:55  rvalue  阅读(302)  评论(0)  编辑  收藏  举报