题解 CF1349F2 【Slime and Sequences (Hard Version)】

官方题解比较简略,这里我来补充一个详细一点的。 此题 : 国际计数水平.jpg **题意** : 定义好序列为 : 都是正整数 , 如果某个数 $k$ 出现过 , 那么一定有 $k-1$ 出现在其前面。 现在对 $k=1...n$ ,分别统计其在长度为 $n$ 的所有好序列中,出现次数的和。 - `Easy` $\ n\leq 5000$ , 时限$\texttt{2s}$ - `Hard` $\ n\leq 10^5$ , 时限$\texttt{3s}$ ------------ 中 国 特 色 计数工业总产值跃居世界第一! 这位老哥您先别走,能不能告诉我 GF 是什么意思, 拉格朗日反演又是啥? ## 0.前置芝士 - 分式域 - (扩展)拉格朗日反演 - 欧拉数 - 二元 `GF` ## 1.转化 这个好序列的限制十分刻意,首先考虑转化一下。 事实上,可以建立起长度为 $n$ 的排列和长度为 $n$ 的好序列之间的**对应关系**。 在排列 $p$ 相邻数之间填写不等号,序列中 $a_{p_i}$ 的值就是 $p_{1...i}$ 之间的小于号个数$+1$。 - 证明其满足好序列的条件 : 考虑一个 $<$ ,以此为分界线,前面为 $p_{i-1}$ ,后面为 $p_{i}$ 。有 $a_{p_{i-1}}+1=a_{p_i}$ ,因为中间隔了一个小于号。 同时又有 $p_{i-1}<p_i$(位置) ,也是因为小于号。实际上小于号就是把两个限制结合在了一起。 通过这个构造,也容易从好序列还原排列。具体证明略。 ~~学长说这可以入选“十大炫酷映射”~~ 以本人的憨憨水平,考场上根本想不出来这种精妙的构造。甚至还根本没去想等价转化,直接对着原问题分析,连正解的边都没摸着,直接被拒之门外。这可能也就是我正式比赛中没切过任何一道计数工业的原因吧…… 目前大概想到了几点。 首先如果原问题复杂(条件信息量太大),感觉不可做,那么大概率是出题人做了神仙转化,要去寻找构造或者充要条件。 其次,可以统计需要计数的东西的总数,比如这个题,如果发现总数是 $n!$ ,就更容易想到构造成排列。 ## 2.初步分析 现在就是要对排列统计 $<$ 的个数。不难发现其实就是 [P5825 排列计数](https://www.luogu.com.cn/problem/P5825) 中的**欧拉数**。 设 $\left\langle\begin{matrix}n \\ k\end{matrix}\right\rangle$ 为 : 长度为 $n$ 的排列,小于号的个数为 $k$ 的方案数。 我们分别考虑每个位置的贡献。 $${\rm Ans}[k]=\sum\limits_{i=\max(k,1)}^n\left\langle\begin{matrix}i \\ k\end{matrix}\right\rangle\dbinom{n}{i}(n-i)!$$ 意思是 : 枚举在好序列中产生 $k$ 的位置。此时整个排列的方案数,就是对局部排列进行标号分配。后面的部分没有顺序限制,还能随意排列。 其实只有 ${\rm Ans}[0]$ 会被 $\max(k,1)$ 影响,我们可以特判,下界就会变得更漂亮。 存在一个$O(n^2)$递推欧拉数的做法,然后暴力计算上式,即可通过`Easy`版。 ```cpp #include<algorithm> #include<cstdio> #define MaxN 5050 #define ll long long using namespace std; const int mod=998244353; ll powM(ll a,int t=mod-2){ ll ret=1; while(t){ if (t&1)ret=ret*a%mod; a=a*a%mod;t>>=1; }return ret; } ll fac[MaxN],ifac[MaxN],s[MaxN],s2[MaxN],ans[MaxN]; int n; int main() { scanf("%d",&n); fac[0]=1; for (int i=1;i<=n;i++) fac[i]=fac[i-1]*i%mod; ifac[n]=powM(fac[n]); for (int i=n;i;i--) ifac[i-1]=ifac[i]*i%mod; s2[0]=1; for (int i=1;i<=n;i++){ s[0]=1;ans[0]=(ans[0]+ifac[i])%mod; for (int j=1;j<=i;j++){ s[j]=(s2[j]*(j+1)+s2[j-1]*(i-j))%mod; ans[j]=(ans[j]+s[j]*ifac[i])%mod; }for (int j=0;j<=i;j++)s2[j]=s[j]; }for (int i=0;i<n;i++) printf("%lld ",ans[i]*fac[n]%mod); return 0; } ``` ## 3.加速计算 这个式子涉及到了$O(n^2)$个的大量欧拉数,显然分离考虑是不行的。我们需要把 “黑盒” 欧拉数换成更易于批量分析的形式。 先设 $\left|\begin{matrix}n \\ k\end{matrix}\right|$ 为 : 长度为 $n$ 的排列,钦定 $k$ 个小于号,其余随意的方案数(可能重复计数)。 那么,根据**二项式反演**,非常经典地有 $$\left|\begin{matrix}n \\ k\end{matrix}\right|=\sum\limits_{i=k}^n\dbinom{i}{k}\left\langle\begin{matrix}n \\ i\end{matrix}\right\rangle\quad\Longrightarrow\quad \left\langle\begin{matrix}n \\ k\end{matrix}\right\rangle=\sum\limits_{i=k}^n\dbinom{i}{k}(-1)^{i-k}\left|\begin{matrix}n \\ i\end{matrix}\right|$$ 我们考虑 $\left|\begin{matrix}n \\ k\end{matrix}\right|$ 如何计算,注意到一串**连续**的 $<$ 可以视作一块,然后块中放置数的顺序就是固定的。 考虑单个块的`EGF`,为 $\{0,1,1,1,\dots\}\longrightarrow e^x-1$ (不考虑空集) 如果恰好需要 $k$ 个 $<$ ,则有 $n-k$ 个块。 这是个**有序组合并分配标号**,能得到 $\left|\begin{matrix}n \\ k\end{matrix}\right|=n^{n-k}$ 则 $\left\langle\begin{matrix}n \\ k\end{matrix}\right\rangle=n!\sum\limits_{i=k}^n\dbinom{i}{k}(-1)^{i-k}[x^n](e^x-1)^{n-i}$ 代入`Easy`的式子中,得到 $${\rm Ans}[k]=\sum\limits_{i=k}^n\dbinom{n}{i}(n-i)!i!\sum\limits_{j=k}^i\dbinom{j}{k}(-1)^{j-k}[x^i](e^x-1)^{i-j}$$ $$=n!\sum\limits_{i=k}^n\sum\limits_{j=k}^i\dbinom{j}{k}(-1)^{j-k}[x^i](e^x-1)^{i-j}$$ $$=\dfrac{n!}{k!}\sum\limits_{j=k}^n\dfrac{(-1)^{j-k}j!}{(j-k)!}\sum\limits_{i=j}^n[x^i](e^x-1)^{i-j}$$ 设$S[k]=\sum\limits_{i=k}^n[x^i](e^x-1)^{i-k}$ $$=\dfrac{n!}{k!}\sum\limits_{j=k}^n\dfrac{(-1)^{j-k}}{(j-k)!}S[j]j!$$ 得知 $S[0...k]$ 之后显然可以**差卷积**计算。 - 观察 $S[k]=\sum\limits_{i=k}^n[x^i](e^x-1)^{i-k}$ $=\sum\limits_{i=k}^n[x^k]\Big(\dfrac{e^x-1}{x}\Big)^{i-k}$ $=\sum\limits_{i=0}^{n-k}[x^k]\Big(\dfrac{e^x-1}{x}\Big)^{i}$ 形式统一了许多。 设 $F(x)=\dfrac{e^x-1}{x}$ 。 我们的目标就是对 $k=0...(n-1)$ 求出 $S[k]=\sum\limits_{i=0}^{n-k}[x^k]F(x)^{i}$ $=[x^k]\bigg(\dfrac{1-F(x)^{n-k+1}}{1-F(x)}\bigg)$ $=[x^k]\dfrac{1}{1-F(x)}-[x^k]\dfrac{F(x)^{n-k+1}}{1-F(x)}$ 第一部分可以多项式求逆,但是注意到分母的常数项为 $0$ ,不能直接上。 - 考虑修改求逆的定义 : $G(x)=x^kF(x)$, 使得 $F(x)$ 有常数项，那么 $\dfrac{1}{G(x)} = x^{-k}\dfrac{1}{F(x)}$ **附** : 这需要引入**分式域**(负数次数),拉格朗日反演在这里仍然成立,所以不用担心干扰推导。 通俗地讲,就是先把整个式子乘以 $x^k$ ,相当于对分母除 $x^k$ ,然后可以经典求逆。提取系数的时候要考虑多余的 $x^k$ ,提取高$k$位,所以多项式的界要大一点。 这里实际上最低只有 $[x^{-1}]$ 项,所以也不是很麻烦。 接下来考虑第二部分。(分母的逆同样拥有负指数,不过没关系) 把次数补齐 : $[x^k]\dfrac{F(x)^{n-k+1}}{1-F(x)}=[x^{n+1}]\dfrac{\big(xF(x)\big)^{n-k+1}}{1-F(x)}$ 考虑**加入一个元以区分信息**。 $=[x^{n+1}y^{n-k+1}]\sum\limits_{i=0}\dfrac{\big(xF(x)y\big)^i}{1-F(x)}$ $=[x^{n+1}y^{n-k+1}]\dfrac{1}{1-F(x)}\dfrac{1}{1-xF(x)y}$ 注意到右边是与 $k$ 无关的,我们只需要求出 $[x^{n+1}]$ 项关于 $y$ 的多项式就得到了答案。 考虑使用拉格朗日反演,多项式复合系列理论。 但注意到 $\dfrac{1}{1-F(x)}\dfrac{1}{1-xF(x)y}$ 除了 $F(x)$ 还有额外的 $x$ ,不容易用经典的多项式复合直接描述。 先设$W(x)=xF(x)$ , $H(x)$ 满足 $\dfrac{W(x)}{H(W(x))}=x$, 即 $\dfrac{xF(x)}{H(W(x))}=x$ , 则 $F(x)=H(W(x))$ (我们使用 $H(x)$ 给 $xF(x)$ 消去了 $x$ ) 又设 $G(x)=\dfrac{1}{1-H(x)}\dfrac{1}{1-xy}$ - **注** : 手动爆算可求得 $G'(x)=\dfrac{y+H'(x)-yH(x)-xyH'(x)}{\big(1-H(x)\big)^2\big(1-xy\big)^2}$ 那么有 $G(W(x))=\dfrac{1}{1-H(W(x))}\dfrac{1}{1-W(x)y}=\dfrac{1}{1-F(x)}\dfrac{1}{1-xF(x)y}$ (现在我们以经典多项式复合描述了问题) 现在的问题就是 $[x^{n+1}]G(W(x))$。 设 $P(x)$ 为 $W(x)$ 的复合逆 ,不难发现 $P(x)=\dfrac{x}{H(x)}$ , 则 $\dfrac{P(x)}{x}=\dfrac{1}{H(x)}$。 使用**扩展拉格朗日反演**可得 $$[x^{n+1}]H(W(x))=\dfrac{1}{n+1}[x^n]\dfrac{y+H'(x)-yH(x)-xyH'(x)}{\big(1-H(x)\big)^2\big(1-xy\big)^2}\left(\dfrac{P(x)}{x}\right)^{-(n+1)}$$ 忽略常数 $\frac{1}{n+1}$ 并代换 $$=[x^n]\dfrac{y\big(1-H(x)\big)+\big(1-xy\big)H'(x)}{\big(1-H(x)\big)^2\big(1-xy\big)^2}H(x)^{n+1}$$ $$=[x^n]\left(\dfrac{y}{\big(1-H(x)\big)\big(1-xy\big)^2}+\dfrac{H'(x)}{\big(1-H(x)\big)^2\big(1-xy\big)}\right)H(x)^{n+1}$$ 展开关于 $y$ 的封闭形式 $$=[x^n]\left(\dfrac{1}{\big(1-H(x)\big)}\sum\limits_{i=0}(i+1)x^iy^{i+1}+\dfrac{H'(x)}{\big(1-H(x)\big)^2}\sum\limits_{i=0}(xy)^i\right)H(x)^{n+1}$$ 考虑提取 $[y^m]$ 系数 $$[x^ny^m]\left(\dfrac{1}{\big(1-H(x)\big)}\sum\limits_{i=0}(i+1)x^iy^{i+1}+\dfrac{H'(x)}{\big(1-H(x)\big)^2}\sum\limits_{i=0}(xy)^i\right)H(x)^{n+1}$$ $$=[x^n]\left(\dfrac{mx^{m-1}}{\big(1-H(x)\big)}+\dfrac{H'(x)x^m}{\big(1-H(x)\big)^2}\right)H(x)^{n+1}$$ $$=[x^{n-m+1}]m\dfrac{H(x)^{n+1}}{\big(1-H(x)\big)}+[x^{n-m}]\dfrac{H'(x)H(x)^{n+1}}{\big(1-H(x)\big)^2}$$ 分别求出 $\dfrac{H(x)^{n+1}}{\big(1-H(x)\big)}$ 和 $\dfrac{H'(x)H(x)^{n+1}}{\big(1-H(x)\big)^2}$ 即可,问题变为了求 $H(x)$。 常规的多项式复合逆肯定是不行的,注意到问题的特殊性 : $W(x)=xF(x)=e^x-1$ 容易构造 $\ln(e^x-1+1)=x$ , 那么有 $P(x)=\ln(x+1)$ ,则 $H(x)=\dfrac{x}{\ln(x+1)}$。 ( 此处有 $1-H(x)$ 的常数项为 $0$ ,仍需要分式域,记得偏移次数 ) 综上,我们以 $O(n\log n)$ 或者 $O(n\log^2n)$ 的复杂度解决了本题。 实现中使用了大力倍增快速幂,时限比较松就卡过去了。 ```cpp #include<algorithm> #include<cstring> #include<cstdio> #define ll long long #define clr(f,n) memset(f,0,sizeof(ll)*(n)) #define cpy(f,g,n) memcpy(f,g,sizeof(ll)*(n)) using namespace std; const int _G=3,mod=998244353,Maxn=135000; ll powM(ll a,ll t=mod-2){ ll ans=1; while(t){ if(t&1)ans=ans*a%mod; a=a*a%mod;t>>=1; }return ans; } const int invG=powM(_G); int tr[Maxn<<1],tf; void tpre(int n){ if (tf==n)return ;tf=n; for(int i=0;i<n;i++) tr[i]=(tr[i>>1]>>1)|((i&1)?n>>1:0); } void NTT(ll *g,bool op,int n) { #define ull unsigned long long tpre(n); static ull f[Maxn<<1],w[Maxn]={1}; for (int i=0;i<n;i++)f[i]=(((ll)mod<<5)+g[tr[i]])%mod; for(int l=1;l<n;l<<=1){ ull tG=powM(op?_G:invG,(mod-1)/(l+l)); for (int i=1;i<l;i++)w[i]=w[i-1]*tG%mod; for(int k=0;k<n;k+=l+l) for(int p=0;p<l;p++){ int tt=w[p]*f[k|l|p]%mod; f[k|l|p]=f[k|p]+mod-tt; f[k|p]+=tt; } if (l==(1<<10)) for (int i=0;i<n;i++)f[i]%=mod; }if (!op){ ull invn=powM(n); for(int i=0;i<n;++i) g[i]=f[i]%mod*invn%mod; }else for(int i=0;i<n;++i)g[i]=f[i]%mod; } void px(ll *f,ll *g,int n) {for(int i=0;i<n;++i)f[i]=f[i]*g[i]%mod;} void times(ll *f,ll *g,int m) { int n=1;while(n<m+m)n<<=1; NTT(f,1,n);NTT(g,1,n); px(f,g,n);NTT(f,0,n);clr(f+m,n-m); } void invp(ll *f,int m) { static ll sav[Maxn<<1],w[Maxn<<1],r[Maxn<<1]; int n;for (n=1;n<m;n<<=1); w[0]=powM(f[0]); for (int len=2;len<=n;len<<=1){ for (int i=0;i<(len>>1);i++) r[i]=(w[i]<<1)%mod; cpy(sav,f,len); NTT(w,1,len<<1);px(w,w,len<<1); NTT(sav,1,len<<1);px(w,sav,len<<1); NTT(w,0,len<<1);clr(w+len,len); for (int i=0;i<len;i++) w[i]=(r[i]-w[i]+mod)%mod; }cpy(f,w,m);clr(sav,n+n);clr(w,n+n);clr(r,n+n); } void dao(ll *f,int m){ for (int i=1;i<m;i++) f[i-1]=f[i]*i%mod; f[m-1]=0; } void powp(ll *f,int m,int t) { static ll g[Maxn<<1]; int n=1;while(n<m+m)n<<=1;g[0]=1; while(t){ NTT(f,1,n); if (t&1){ NTT(g,1,n); px(g,f,n);NTT(g,0,n); clr(g+m,n-m); }px(f,f,n);NTT(f,0,n); clr(f+m,n-m);t>>=1; }cpy(f,g,m);clr(f+m,n-m);clr(g,n); } ll fac[Maxn],ifac[Maxn]; void Init(int m) { fac[0]=1; for (int i=1;i<=m;i++) fac[i]=fac[i-1]*i%mod; ifac[m]=powM(fac[m]); for (int i=m;i;i--) ifac[i-1]=ifac[i]*i%mod; } int n,m; ll F[Maxn<<1],G[Maxn<<1],S[Maxn<<1],H[Maxn<<1],H1[Maxn<<1],H2[Maxn<<1]; int main() { scanf("%d",&n);Init(n+6); int len=1;while(len<n+n+4)len<<=1; m=n+2; for (int i=0;i<m+1;i++)G[i]=mod-ifac[i+2]; invp(G,m+1);cpy(S,G+1,m+1); for (int i=0;i<m+1;i++){ H[i]=fac[i]*ifac[i+1]%mod; if (i&1)H[i]=mod-H[i]; }invp(H,m+1);cpy(H1,H,m); powp(H1,m,n+1); cpy(G,H,m+1);dao(G,m+1); cpy(H2,H1,m);times(H2,G,m);clr(G,len); for (int i=0;i<m;i++)G[i]=mod-H[i+1]; invp(G,m); NTT(G,1,len);NTT(H1,1,len);NTT(H2,1,len); px(H1,G,len);NTT(H1,0,len);clr(H1+m,len-m); px(H2,G,len);NTT(H2,0,len);clr(H2+m,len-m);NTT(H2,1,len); px(H2,G,len);NTT(H2,0,len);clr(H2+m,len-m); for (int i=1,sav=powM(n+1);i<n+2;i++){ int m=n-i+1; S[i]=(S[i]-(H1[n-m+2]*m+H2[n-m+2])%mod*sav)%mod; }S[0]=n; for (int i=0;i<n;i++)S[i]=(S[i]+mod)%mod*fac[i]%mod; for (int i=0;i<n;i++)F[i]=(i&1) ? mod-ifac[i] : ifac[i]; clr(F+n,len-n);clr(S+n,len-n);reverse(S,S+n); NTT(F,1,len);NTT(S,1,len); px(S,F,len);NTT(S,0,len); clr(S+n,len-n);reverse(S,S+n); for (int i=0;i<n;i++) printf("%lld ",S[i]*fac[n]%mod*ifac[i]%mod); return 0; } ```