时间戳

首页 > 技术文章 > CF923E Perpetual Subtraction

CF923E Perpetual Subtraction

zzctommy 2021-01-11 20:21 原文

CF923E Perpetual Subtraction

这题超出了我这个初中生的数学能力（我目前只能掌握导数，积分有点生疏），没能自己做出来，后半部分参考了 rqy的题解 。

所以这篇题解主要是写给我自己看的。

首先有一个显然的 \(O(n^2m)\) 的dp做法。

\(f_{i,j}\) 表示第 \(i\) 轮 \(X=j\) 的概率，\(f_{0,i}=p_i\)

\[f_{i,j}=\sum_{i=j}^{n}\dfrac{f_{i-1,j}}{j+1} \]

当然可以直接后缀和优化成 \(O(nm)\) 不过半点用处没有。

注意到这个每次除 \(j+1\) 很像积分。事实上这步并不难想，沉下心来都能想到的，我都能想到。

然而对于积分这玩意我只会用初等方法推。事实上我初等方法都推出了类似的结果，但是根本不可能从dp式子看出这是一个多项式平移就GG了

设 \(k\) 轮之后答案的生成函数为（其实就是PGF）\(F_k(x)=\sum f_{k,i}x^i\) 。

\[F_{k+1}(x)=\sum_{i=0}^{n}\sum_{j=i}^{n}\dfrac{f_{k,j}}{j+1}x^i\\ =\sum_{j=0}^{n}\dfrac{f_{k,j}}{j+1}\sum_{i=0}^{j}x^i\\ =\sum_{j=0}^{n}\dfrac{f_{k,j}}{j+1}\dfrac{1-x^{j+1}}{1-x}\\ =\dfrac{1}{x-1}\sum_{j=0}^{n}f_{k,j}\sum_{i=0}^{j}\dfrac{x^{j+1}-1}{j+1}\\ =\dfrac{1}{x-1}\sum_{j=0}^{n}f_{k,j}\int_{1}^{x}t^j\mathrm{d}t\\ =\dfrac{1}{x-1}\int_{1}^{x}F_k(t)\mathrm{d}t\\ \]

我们下标是从 \(0\) 开始的，但是这个积分从 \(1\) 开始，并且除的是 \(x-1\)。如果除 \(x\) 就相当于系数平移了一位 ，会很方便。

考虑新开一个函数表示 \(F_k(x)\) 整体平移的结果。

令 \(G_k(x)=\sum\limits_{i=0}^{n}g_{k,i}x^i\)

\[G_{k+1}(x)=F_{k}(x+1)=\dfrac{\int_{1}^{x+1}F_{k}(t)\mathrm{d}t}{x+1-1}=\dfrac{\int_{0}^{x}G_{k}(t)\mathrm{d}t}{x} \]

所以 \(g_{k+1,i}=\dfrac{g_{k,i}}{i+1}\) ，所以 \(g_{m,i}=\dfrac{g_{0,i}}{(i+1)^{m}}\) 。

又有

\[G_{k}(x)=F_{k}(x+1)=\sum_{i=0}^{n}f_{k,i}(x+1)^{i}=\sum_{i=0}^{n}f_{k,i}\sum_{j=0}^{i}\binom{i}{j}x^j \]

对比两边系数可得

\[g_{k,i}=\sum\limits_{j=i}^{n}\dbinom{j}{i}f_{k,j} \]

二项式反演得

\[f_{k,i}=\sum_{j=i}^{n}(-1)^{j-i}\binom{j}{i}g_{k,j} \]

两个式子显然都可以NTT加速，目测必然可以卷积。

所以我们需要做的就是，把 \(f_{0,i}=p_i\) 转成 \(g_{0,i}\) ，然后快速幂把 \(g_{0,i}\) 推到 \(g_{m,i}\) ，然后再从 \(g_{m,i}\) 求出 \(f_{m,i}\) 。

过程中只需要两遍NTT就好了。

太强大了，数学太强大了！！！

不过式子有点多，一定要抄对啊qaq。

剩下部分为自己推卷积，可以直接跳过不看。

\[g_i=\sum_{j=i}^{n}\binom{j}{i}f_j\\ g_i=\sum_{j=i}^{n}\dfrac{j!}{i!(j-i)!}f_j\\ g_i=\dfrac{1}{i!}\sum_{j=i}^{n}\dfrac{j!f_j}{(j-i)!} \]

令 \(A_i=i!f_i,B_i=\dfrac{1}{(n-i)!}\) ，\([x^k](A*B)=\sum\limits_{i=0}^{k}A_iB_{k-i}=\sum\limits_{i=0}^{k}\dfrac{i!f_i}{(n-k+i)!}\) ，\(g_i=\dfrac{1}{i!}[x^{n+i}]\) 。

\[f_i=\sum_{j=i}^{n}(-1)^{j-i}\binom{j}{i}g_j\\ f_i=\dfrac{(-1)^{i}}{i!}\sum_{j=i}^{n}\dfrac{(-1)^{j}j!g_j}{(j-i)!} \]

\(A_i=(-1)^{i}i!g_i,B_i=\dfrac{1}{(n-i)!}\) ，\([x^k](A*B)=\sum\limits_{i=0}^{k}\dfrac{(-1)^{i}i!g_i}{(n-k+i)!}\) ，\(f_i=\dfrac{(-1)^i}{i!}[x^{n+i}]\)

CP-Editor这玩意挺好玩的。

// Problem: CF923E Perpetual Subtraction
// Contest: Luogu
// URL: https://www.luogu.com.cn/problem/CF923E
// Memory Limit: 250 MB
// Time Limit: 2000 ms
// 
// Powered by CP Editor (https://cpeditor.org)

#include<bits/stdc++.h>
using namespace std;
#define fi first
#define se second
#define mkp(x,y) make_pair(x,y)
#define pb(x) push_back(x)
#define sz(v) (int)v.size()
typedef long long LL;
typedef double db;
template<class T>bool ckmax(T&x,T y){return x<y?x=y,1:0;}
template<class T>bool ckmin(T&x,T y){return x>y?x=y,1:0;}
#define rep(i,x,y) for(int i=x,i##end=y;i<=i##end;++i)
#define per(i,x,y) for(int i=x,i##end=y;i>=i##end;--i)
inline int read(){
	int x=0,f=1;char ch=getchar();
	while(!isdigit(ch)){if(ch=='-')f=0;ch=getchar();}
	while(isdigit(ch))x=x*10+ch-'0',ch=getchar();
	return f?x:-x;
}
const int N=100005;
const int M=N<<2;
#define mod 998244353

namespace math{
inline int qpow(int n,int k){int res=1;for(;k;k>>=1,n=1ll*n*n%mod)if(k&1)res=1ll*n*res%mod;return res;}
int fac[N],ifc[N],inv[N];
void initmath(const int&n=N-1){
	fac[0]=1;for(int i=1;i<=n;++i)fac[i]=1ll*i*fac[i-1]%mod;
	ifc[n]=qpow(fac[n],mod-2);for(int i=n-1;i>=0;--i)ifc[i]=1ll*(i+1)*ifc[i+1]%mod;
	inv[1]=1;for(int i=2;i<=n;++i)inv[i]=1ll*inv[mod%i]*(mod-mod/i)%mod;
}

}

using math::qpow;

namespace poly{
	
int rev[M],lg,lim;
void init_poly(const int&n){
	for(lg=0,lim=1;lim<n;lim<<=1,++lg);
	for(int i=0;i<lim;++i)rev[i]=(rev[i>>1]>>1)|((i&1)<<(lg-1));
}
void NTT(int*a,int op){
	for(int i=0;i<lim;++i)
		if(i>rev[i])swap(a[i],a[rev[i]]);
	const int g=op?3:math::inv[3];
	for(int i=1;i<lim;i<<=1){
		const int wn=qpow(g,(mod-1)/(i<<1));
		for(int j=0;j<lim;j+=i<<1){
			int w0=1;
			for(int k=0;k<i;++k,w0=1ll*w0*wn%mod){
				const int X=a[j+k],Y=1ll*a[i+j+k]*w0%mod;
				a[j+k]=(X+Y)%mod,a[i+j+k]=(X+mod-Y)%mod;
			}
		}
	}
	if(op)return;const int ilim=qpow(lim,mod-2);
	for(int i=0;i<lim;++i)a[i]=1ll*a[i]*ilim%mod;
}
#define clr(a,n) memset(a,0,sizeof(int)*(n))
#define cpy(a,b,n) memcpy(a,b,sizeof(int)*(n))
void poly_mul(int*f,int*g,int*ans,int n,int m){
	static int A[M],B[M];init_poly(n+m);
	cpy(A,f,n),clr(A+n,lim-n),NTT(A,1);
	cpy(B,g,m),clr(B+m,lim-m),NTT(B,1);
	for(int i=0;i<lim;++i)ans[i]=1ll*A[i]*B[i]%mod;
	NTT(ans,0);
}

}

int n,p[N],A[M],B[M],g[M],f[M];
LL m;
signed main(){
	math::initmath();scanf("%d%lld",&n,&m);
	for(int i=0;i<=n;++i)p[i]=read();
	for(int i=0;i<=n;++i)A[i]=1ll*math::fac[i]*p[i]%mod;
	for(int i=0;i<=n;++i)B[i]=math::ifc[n-i];
	poly::poly_mul(A,B,A,n+1,n+1);
	for(int i=0;i<=n;++i)g[i]=1ll*A[n+i]*math::ifc[i]%mod;
	
	m%=mod-1;
	for(int i=0;i<=n;++i)g[i]=1ll*g[i]*qpow(math::inv[i+1],m)%mod;
	
	for(int i=0;i<=n;++i)A[i]=1ll*math::fac[i]*g[i]%mod*(i&1?mod-1:1)%mod;
	for(int i=0;i<=n;++i)B[i]=math::ifc[n-i];
	poly::poly_mul(A,B,A,n+1,n+1);
	for(int i=0;i<=n;++i)f[i]=1ll*A[n+i]*(i&1?mod-math::ifc[i]:math::ifc[i])%mod;
	
	for(int i=0;i<=n;++i)printf("%d ",f[i]);
	return 0;
}

推荐阅读