【CTSC2019】珍珠【二项式反演】【生成函数】

Description

  • 一个长度为(n),每一位都是值域([1,D])的整数随机变量的序列,问其中能选出至少(m)对相同元素(每个元素只能选一次)的概率。(Dle 10^5,nle 10^9)
  • 传送门

Solution

一个序列符合条件当且仅当其满足:

[sum_{i=1}^{D}lfloorfrac {c_i}{2} floorge m ]

其中(c_i)为元素(i)出现的次数。

变形得到:

[sum_{i=1}^{n}frac{c_i-(c_imod 2)}{2}ge m ]

[n-sum_{i=1}^{n}c_imod 2ge 2m ]

因此设(f_i)表示恰好有(i)种元素出现了奇数次的方案数,那么:

[ans=sum_{i=0}^{n-2m}f_i ]

直接求(f_i)比较麻烦,使用二项式反演经典套路,先求出(g_i)表示钦定(i)种元素出现了奇数次,其他元素随便取得方案数。那么:

[g_i=sum_{j=i}^{D}f_iinom{i}{j}Rightarrow f_i=sum_{j=i}^{D}g_jinom{i}{j}(-1)^{i-j} ]

因此如果求得(g),可以直接(mathcal O(Dlog(D)))卷积得到(f)

对于(g),因为数有顺序,所以写出(g)(EGF):考虑被钦定的那(k)个颜色,只在选择奇数次时作出贡献,因此贡献为([0,1,0,1,dots])(EGF)(dfrac{e^x-e^{-x}}{2})

[egin{aligned} g_k&=inom{D}{k}n![x^n](frac{e^x-e^{-x}}{2})^k(e^x)^{D-k}\ &=inom{D}{k}frac{n!}{2^k}[x^n](e^x-e^{-x})^k(e^x)^{D-k}\ &=inom{D}{k}frac{n!}{2^k}[x^n]sum_{j=0}^{k}inom{k}{j}(e^x)^j(-e^{-x})^{k-j}(e^x)^{D-k}\ &=inom{D}{k}frac{n!}{2^k}sum_{j=0}^{k}inom{k}{j}(-1)^{k-j}[x^n](e^x)^{D-2(k-j)}\ end{aligned} ]

注意到([x^n]e^{ax}=frac{a^n}{n!})

[egin{aligned} g_k&=inom{D}{k}frac{1}{2^k}sum_{j=0}^{k}inom{k}{j}(-1)^{k-j}[D-2(k-j)]^{n}\ &=inom{D}{k}frac{k!}{2^k}sum_{j=0}^{k}frac{(-1)^{k-j}[D-2(k-j)]^{n}}{j!(k-j)!}\ end{aligned} ]

于是令

[a_i=frac{1}{i!},b_i=frac{(-1)^{i}(D-i)^n}{i!} ]

即可卷积得到(g),总复杂度为(mathcal O(Dlog(D)))

Code

#include<bits/stdc++.h>
using namespace std;

const int N=(1<<21)+20;
const int mod=998244353;
const int gg=3;
const int giv=(mod+1)/3;
namespace Math{
	inline int add(int x,int y){return (x+y>=mod)?x+y-mod:x+y;}
	inline int dec(int x,int y){return (x-y<0)?x-y+mod:x-y;}
	inline void upd(int &x,int y){x=add(x,y);}
	inline int ksm(int x,int y){
		int ret=1;
		for(;y;y>>=1,x=1ll*x*x%mod) if(y&1) ret=1ll*ret*x%mod;
		return ret;
	}
}
using namespace Math;

namespace Container{
	struct poly{
		vector<int> v;
		inline int& operator[](int x){
			while(x>=v.size()) v.push_back(0);
			return v[x];
		}	
		inline void resize(int x){v.resize(x);}
	};
	int Wn[2][N],p1,p2,lg[N];
	inline void init_poly(int n){
		p1=1;
		while((p1<<1)<=n) p1<<=1;
		p2=p1<<1;
		for(int i=2;i<p2;++i) lg[i]=lg[i>>1]+1;
		int wn0=ksm(gg,(mod-1)/p2),wn1=ksm(giv,(mod-1)/p2);
		Wn[0][p1]=Wn[1][p1]=1;
		for(int i=p1+1;i<=p2;++i) Wn[0][i]=1ll*Wn[0][i-1]*wn0%mod,Wn[1][i]=1ll*Wn[1][i-1]*wn1%mod;
		for(int i=p1-1;i>=1;--i) Wn[0][i]=Wn[0][i<<1],Wn[1][i]=Wn[1][i<<1];
	}
}
using namespace Container;

namespace basic{
	int r[N],fr[N];
	inline void init(int lim,int len){
		for(int i=0;i<lim;++i) r[i]=(r[i>>1]>>1)|((i&1)<<len);
	}
	inline void NTT(poly& f,int lim,int tp){
		for(int i=0;i<lim;++i) fr[i]=f[r[i]];
		for(int mid=1;mid<lim;mid<<=1){
			for(int len=mid<<1,l=0;l+len-1<lim;l+=len){
				for(int k=l;k<l+mid;++k){
					int w1=fr[k],w2=1ll*fr[k+mid]*Wn[tp][k-l+mid]%mod;
					fr[k]=add(w1,w2),fr[k+mid]=dec(w1,w2);
				}
			}
		}
		for(int i=0;i<lim;++i) f[i]=fr[i];
	}
	inline poly mul(int n,int m,poly f,poly g){
		f.resize(n);g.resize(m);
		int len=lg[n+m-1],lim=1<<(len+1);
		init(lim,len);
		NTT(f,lim,0);NTT(g,lim,0);
		for(int i=0;i<lim;++i) f[i]=1ll*f[i]*g[i]%mod;
		NTT(f,lim,1);
		int iv=ksm(lim,mod-2);
		for(int i=0;i<lim;++i) f[i]=1ll*f[i]*iv%mod;
		return f;
	}
} 
using namespace basic;

int n,D,m,fac[N],inv[N];
inline void init_math(int n){
	fac[0]=1;for(int i=1;i<=n;++i) fac[i]=1ll*fac[i-1]*i%mod;
	inv[n]=ksm(fac[n],mod-2);
	for(int i=n-1;i>=0;--i) inv[i]=1ll*inv[i+1]*(i+1)%mod;
}
inline int binom(int n,int m){return 1ll*fac[n]*inv[m]%mod*inv[n-m]%mod;}
inline int sgn(int x){return (x&1)?mod-1:1;}
int main(){
	scanf("%d%d%d",&D,&n,&m);
	if(D<=n-2*m){printf("%d
",ksm(D,n));return 0;}
	if(n<2*m){puts("0");return 0;}
	init_poly(D<<1);init_math(D<<1);
	poly f,g;
	for(int i=0;i<=D;++i) f[i]=1ll*inv[i]*sgn(i)%mod*ksm((D-2*i+mod)%mod,n)%mod,g[i]=inv[i];
	f=mul(D+1,D+1,f,g);
	int iv=ksm(2,mod-2);
	for(int i=0;i<=D;++i){
		f[i]=1ll*f[i]*binom(D,i)%mod*ksm(iv,i)%mod;
		f[i]=1ll*f[i]*fac[i]%mod*fac[i]%mod;
	}
	for(int i=0;i<=D;++i) g[i]=1ll*sgn(D-i)*inv[D-i]%mod;
	f=mul(D+1,D+1,f,g);
	int ans=0;
	for(int i=0;i<=n-2*m;++i) upd(ans,1ll*inv[i]*f[i+D]%mod);
	printf("%d
",ans);
	return 0;
}
原文地址:https://www.cnblogs.com/tqxboomzero/p/14504437.html