Hello, Solitude.

考虑任意时刻区间均为奇数的情况,即 \(n=2^k-1\)

观察到一个区间被取到的概率与区间外的具体状态无关,只与区间外取了多少次有关,因此考虑计算一个区间被选中后的出现次数,设 \(dp_{[l,r],t}\) 表示区间 \([l,r]\) 在剩余 \(t\) 次入座时被选中的所有状态的出现次数之和,贡献加入中点 \(ans_{mid}\),最终的答案就是 \(\cfrac{ans_i}{n^{\underline{m}}}\)

这个区间总长显然是 \(O(n)\) 的,因为在区间内选取一个位置分裂后的两边区间是本质相同的,只需要转移一边然后复制一下

对于一个区间到其子区间的转移,有:

\[dp_{[l,mid-1],i}=\sum\limits_{t>i}\binom{t-1}{i}(r-l+1)(r-mid)^{\underline{t-i-1}}dp_{[l,r],t} \]

其中 \(r-l+1\) 是选到区间内的方案数,\((r-mid)^{\underline{t-i-1}}\) 是另一个区间 \([mid+1,r]\) 选取的贡献,并且需要带一个组合系数

对于区间转移到答案,有:

\[ans_{mid}=\sum\limits_{t=1}^{r-l+1}(r-l+1)^{\underline t}dp_{[l,r],t} \]

将转移式改写一下,相关项放到一起:

\[\begin{aligned} dp_{[l,mid-1],i}&=(r-l+1)\sum\limits_{t>i}\frac{(t-1)!}{i!(t-i-1)!}\frac{(r-mid)!}{[r-mid-(t-i-1)]!}dp_{[l,r],t}\\ i!dp_{[l,mid-1],i}&=(r-l+1)\sum\limits_{t>i}\binom{r-mid}{t-i-1}\frac{t!dp_{[l,r],t}}{t} \end{aligned} \]

这是一个减法卷积的形式,可以用 NTT 优化

对于偶数长度的区间,只需要将贡献取半分别加入两个中点,分析一下状态数仍然是 \(O(n)\) 的,最终复杂度 \(O(n\log n)\),常数巨大

点击查看代码
#include<bits/stdc++.h>
#define ll long long
#define N 1048577
using namespace std;
int read(){
	int w=0,h=1;char ch=getchar();
	while(ch<'0'||ch>'9'){if(ch=='-')h=-h;ch=getchar();}
	while(ch>='0'&&ch<='9'){w=w*10+ch-'0';ch=getchar();}
	return w*h;
}
int n,m,ans[N];
namespace Polynomial{
	const int mod=998244353,g=3,ig=332748118,B=25000;
	int pwg[B+1],PWG[B+1],pwig[B+1],PWIG[B+1],inv[N<<1],lg[N<<1],pw[N<<1],trans[N<<1];
	inline void add(int&x,int y){(x+=y)>=mod?x-=mod:x;}
	inline void Init(){
		for(int i=pwg[0]=1;i<=B;i++)pwg[i]=1ll*pwg[i-1]*g%mod;
		for(int i=PWG[0]=1;i<=B;i++)PWG[i]=1ll*PWG[i-1]*pwg[B]%mod;
		for(int i=pwig[0]=1;i<=B;i++)pwig[i]=1ll*pwig[i-1]*ig%mod;
		for(int i=PWIG[0]=1;i<=B;i++)PWIG[i]=1ll*PWIG[i-1]*pwig[B]%mod;
		for(int i=pw[0]=inv[1]=1;i<N;i++)pw[i]=pw[i-1]*2ll%mod;
		for(int i=2;i<N;i++)inv[i]=1ll*(mod-mod/i)*inv[mod%i]%mod,lg[i]=lg[i>>1]+1;
	}
	inline int qpow(int k){
		if(k==0)return 1;
		if(k>0)return 1ll*PWG[k/B]*pwg[k%B]%mod;
		if(k<0)return 1ll*PWIG[(-k)/B]*pwig[(-k)%B]%mod;
		return 114514;
	}
	inline int qpow(int b,int k){
		int s=1;
		while(k){
			if(k&1)s=1ll*s*b%mod;
			b=1ll*b*b%mod;k>>=1;
		}
		return s;
	}
	struct Poly{
#define ite vector<int>::iterator
#define rite vector<int>::reverse_iterator
		vector<int>f;
		Poly(){f=vector<int>(0,0);}
		Poly(int x=0,int y=0){f=vector<int>(x,y);}
		Poly(vector<int>x=vector<int>(0)){f=x;}
		inline int size(){return f.size();}
		inline int back(){return f.back();}
		inline bool empty(){return f.empty();}
		inline int&operator[](int x){return f[x];}
		inline ite begin(){return f.begin();}
		inline ite end(){return f.end();}
		inline rite rbegin(){return f.rbegin();}
		inline rite rend(){return f.rend();}
		inline void swap(Poly&b){vector<int>tmp=f;f=b.f;b.f=tmp;}
		inline void sort(ite BEGIN,ite END){sort(BEGIN,END);}
		inline void sort(rite BEGIN,rite END){sort(BEGIN,END);}
		inline void erase(ite BEGIN,ite END){f.erase(BEGIN,END);}
		inline void reverse(ite BEGIN,ite END){reverse(BEGIN,END);}
		inline void emplace_back(int x){f.emplace_back(x);}
		inline void clear(){f.clear();}
		inline void erase(ite x){f.erase(x);}
		inline void pop_back(){f.pop_back();}
		inline void resize(int x){f.resize(x);}
		inline void resize(int x,int y){f.resize(x,y);}
		inline void shrink_to_fit(){f.shrink_to_fit();}
		inline bool operator==(Poly b){return f==b.f;}
		inline bool operator!=(Poly b){return f!=b.f;}
		inline bool operator<(Poly b){return f<b.f;}
		inline bool operator>(Poly b){return f>b.f;}
		inline bool operator<=(Poly b){return f<=b.f;}
		inline bool operator>=(Poly b){return f>=b.f;}
	};
	inline void Complete(Poly&f,int len){while(f.size()<len)f.emplace_back(0);}
	inline void Init(Poly&x,int len){
		int p=pw[lg[len]];
		while(p<len)p<<=1;
		if(len<p)Complete(x,p);
	}// Complete to 2^x
	inline void NTT(Poly&f,int t){
		/*
		if(f.size()==1)return;
		int n=f.size(),buf=1,G=qpow(t*(mod-1)/n);
		Poly l(n>>1),r(n>>1);
		for(int i=0;i<n>>1;++i)l[i]=f[i<<1],r[i]=f[i<<1|1];
		NTT(l,t);NTT(r,t);
		for(int i=0;i<n>>1;++i,buf=1ll*buf*G%mod){
			int p=1ll*buf*r[i]%mod;
			f[i]=(l[i]+p)%mod;
			f[i+(n>>1)]=(l[i]+mod-p)%mod;
		}
		*/
		int n=f.size();
		for(int i=0;i<n;i++)trans[i]=(trans[i>>1]>>1)|((i&1)?n>>1:0);
		for(int i=0;i<n;i++)if(i<trans[i])swap(f[i],f[trans[i]]);
		for(int p=2;p<=n;p<<=1){
			int len=p>>1,G=qpow(t*(mod-1)/p);
			for(int i=0;i<n;i+=p)
				for(int j=i,buf=1;j<i+len;j++,buf=1ll*buf*G%mod){
					int d=1ll*buf*f[j+len]%mod;
					f[j+len]=(f[j]+mod-d)%mod;
					f[j]=(f[j]+d)%mod;
				}
		}
	}// NTT
	inline Poly Convol(Poly a,Poly b){
		int len=a.size()+b.size(),p=pw[lg[len]];
		while(p<len)p<<=1;
		Complete(a,p);NTT(a,1);
		Complete(b,p);NTT(b,1);
		for(int i=0;i<p;i++)a[i]=a[i]*b[i]%mod;
		NTT(a,-1);
		for(int i=0;i<p;i++)a[i]=a[i]*inv[p]%mod;
		while(a.back()==0)a.pop_back();
		return a;
	}// Convolution
	inline Poly Inverse(Poly f){
		Init(f,f.size());
		if(f.size()==1)return Poly(1,qpow(f[0],mod-2));
		int n=f.size(),p=pw[lg[n]];
		Poly hlf(n>>1);
		for(int i=0;i<n>>1;i++)hlf[i]=f[i];
		Poly g=Inverse(hlf);
		while(p<n<<1)p<<=1;
		Complete(f,p);NTT(f,1);
		Complete(g,p);NTT(g,1);
		for(int i=0;i<p;i++)g[i]=(2+mod-f[i]*g[i]%mod)%mod*g[i]%mod;
		NTT(g,-1);
		for(int i=n;i<p;i++)g[i]=0;
		for(int i=0;i<p;i++)g[i]=g[i]*inv[p]%mod;
		return g;
	}// Inverse
	inline Poly operator+(Poly a,Poly b){
		Poly ret(max(a.size(),b.size()));
		for(int i=0;i<a.size();i++)add(ret[i],a.f[i]);
		for(int i=0;i<b.size();i++)add(ret[i],b.f[i]);
		return ret;
	}
	inline Poly operator-(Poly a,Poly b){
		Poly ret(max(a.size(),b.size()));
		for(int i=0;i<a.size();i++)add(ret[i],a.f[i]);
		for(int i=0;i<b.size();i++)add(ret[i],mod-b.f[i]);
		while(ret.f.back()==0)ret.f.pop_back();
		return ret;
	}
	inline Poly operator*(Poly a,Poly b){return Poly(Convol(a.f,b.f));}
	inline Poly operator*(Poly x,int v){return Poly(x*Poly(1,v));}
	inline Poly operator/(Poly a,Poly b){return Poly(Convol(a.f,Inverse(b.f)));}
	inline Poly Derivative(Poly x){
		if(x.size()==0)return x;
		Poly ret(x.size()-1);
		for(int i=0;i<ret.size();i++)ret[i]=x[i+1]*(i+1)%mod;
		return ret;
	}// F'(x)
	inline Poly Integral(Poly x){
		Poly ret(x.size()+1);
		for(int i=1;i<ret.size();i++)ret[i]=x[i-1]*inv[i]%mod;
		return ret;
	}// ∫ F(x)dx
	inline Poly Ln(Poly x){return Integral(Derivative(x)*Inverse(x));}// ln(F(x))
	inline Poly Exp(Poly f){
		Init(f,f.size());
		if(f.size()==1)return Poly(1,1);
		int n=f.size(),p=pw[lg[n]];
		Poly hlf(n>>1);
		for(int i=0;i<n>>1;i++)hlf[i]=f[i];
		Poly g=Exp(hlf);
		while(p<n<<1)p<<=1;
		Complete(f,p);
		Complete(g,p);
		Poly lng=Ln(g);
		Complete(lng,p);f[0]++;
		for(int i=n;i<p;i++)f[i]=lng[i]=0;
		for(int i=0;i<n;i++)f[i]=(f[i]-lng[i]+mod)%mod;
		NTT(f,1);NTT(g,1);
		for(int i=0;i<p;i++)g[i]=g[i]*f[i]%mod;
		NTT(g,-1);
		for(int i=n;i<p;i++)g[i]=0;
		for(int i=0;i<p;i++)g[i]=g[i]*inv[p]%mod;
		return g;
	}// exp(F(x))
	inline Poly operator^(Poly x,int k){return Exp(Ln(x)*k);}
	inline Poly operator^(Poly x,string k){
		int len=k.size(),k1=0,k2=0,k3=0,shift=0;
		for(int i=0;i<len;i++)k1=(k1*10ll%mod+k[i]-'0')%mod;
		for(int i=0;i<len;i++)k2=(k2*10ll%(mod-1)+k[i]-'0')%(mod-1);
		for(int i=0;i<len;i++)if(k3*10ll+k[i]-'0'<=mod)k3=k3*10ll+k[i]-'0';
		while(shift<x.size()&&x[shift]==0)shift++;
		if(shift>=x.size()||shift*k3>=x.size())return Poly(0);
		int inver=qpow(x[shift],mod-2),coef=qpow(x[shift],k2);
		Poly tmp=Poly(x.size()-shift);
		for(int i=0;i+shift<x.size();i++)tmp[i]=x[i+shift]*inver%mod;
		tmp=tmp^(k1%mod);shift=shift*k1;x=Poly(tmp.size()+shift);
		for(int i=0;i<shift;i++)x[i]=0;
		for(int i=0;i<tmp.size();i++)x[i+shift]=tmp[i]*coef%mod;
		return x;
	}// F^k(x)
}
using namespace Polynomial;// call Init() before use
namespace Math{
	int Fac[N],IFac[N];
	void Init(int n){
		for(int i=Fac[0]=1;i<=n;i++)Fac[i]=1ll*Fac[i-1]*i%mod;
		IFac[n]=qpow(Fac[n],mod-2);
		for(int i=n-1;~i;i--)IFac[i]=1ll*IFac[i+1]*(i+1)%mod;
	}
	inline int C(int n,int m){return n<m?0:1ll*Fac[n]*IFac[m]%mod*IFac[n-m]%mod;}
}
using namespace Math;
inline Poly Transform(Poly dp,int n,int m){
	Poly F(m+1),G(dp.size());
	for(int i=0;i<=m;++i)F[i]=C(m,i);
	for(int i=1;i<dp.size();++i)G[i-1]=1ll*dp[i]*inv[i]%mod;
	int len=F.size()+G.size()-1,p=pw[lg[len]];
	while(p<len)p<<=1;
	Complete(F,p);NTT(F,1);
	Complete(G,p);NTT(G,1);
	for(int i=0;i<p;++i)F[i]=1ll*F[i]*G[i]%mod;
	NTT(F,-1);
	for(int i=0;i<p;++i)F[i]=1ll*F[i]*inv[p]%mod;
	for(int i=m;i<n+m+1;++i)F[i-m]=1ll*F[i]*(n+m+1)%mod;
	F.resize(n+1);
	return F;
}
Poly Solve(int len,Poly poly){
	if(len==0)return Poly(0);
	Poly dp(len);
	int mid=(len-1)>>1;
	for(int i=1;i<=len;++i)add(dp[mid],1ll*poly[i]*C(len,i)%mod);
	if(len&1){
		Poly f=Solve(mid,Transform(poly,mid,len-mid-1));
		for(int i=0;i<f.size();++i)add(dp[i],f[i]);
	}
	else{
		Poly f=Solve(mid,Transform(poly,mid,len-mid-1));
		Poly g=Solve(mid+1,Transform(poly,mid+1,len-mid-2));
		for(int i=0;i<f.size();++i)add(dp[i],f[i]);
		for(int i=0;i<g.size();++i)add(dp[i],g[i]);
		for(int i=0;i<g.size();++i)dp[i]=dp[i]*499122177ll%mod;
	}
	for(int i=0;i<len>>1;++i)dp[len-i-1]=dp[i];
	return dp;
}
signed main(){
	Init(n=read());m=read();Init();
	Poly dp(n+1);dp[m]=Fac[m];dp=Solve(n,dp);
	for(int i=0;i<n;i++)printf("%lld\n",1ll*dp[i]*IFac[n]%mod*Fac[n-m]%mod);
	return 0;
}
posted @ 2023-10-11 16:33  pidan007  阅读(27)  评论(0编辑  收藏  举报