LGP5075题解

感觉并不是很难的说？

设 \(F(x)=\sum_{i=1}(Oi^2+Si+U)x^i\)

那么这题要求的就是 \(\sum_{i=1}^{n}F^i(x)=\frac{1-F^{n+1}(x)}{1-F(x)}-1\)。

然后来推一下 \(F(x)\) 是什么东西。

\[F(x)=Ox(x(\frac{1}{1-x})')'+Sx(\frac{1}{1-x})'+\frac{xU}{1-x} \]

\[F(x)=Ox(\frac{x}{(1-x)^2})'+\frac{Sx}{(1-x)^2}+\frac{xU}{1-x} \]

\[F(x)=Ox\times\frac{(1-x)^2+2x(1-x)}{(1-x^4)}+\frac{Sx-Sx^2}{(1-x)^3}+\frac{Ux-2Ux^2+Ux^3}{(1-x)^3} \]

\[F(x)=\frac{Ox-Ox^2+2Ox^2}{(1-x)^3}+\frac{Sx-Sx^2}{(1-x)^3}+\frac{Ux-2Ux^2+Ux^3}{(1-x)^3} \]

\[F(x)=\frac{Ux^3+(O-S-2U)x^2+(O+S+U)x}{(1-x)^3} \]

\[F(x)=\frac{x}{(1-x)^3}\times(Ux^2+(O-S-2U)x+(O+S+U)) \]

可以直接计算 \(F(x)\) 了（卷一次组合数），使用 FFT+快速幂+求逆就可以 \(O(m\log m\log A)\) 了。

#include<cstdio>
#include<cmath>
#define IMP(lim,act) for(int qwq=(lim),i=0;i^qwq;++i)act
typedef double db;
const db Pi=acos(-1);
const int M=1<<15|5;
int n,m,P,O,S,U,F[M],G[M];
struct Barrett{
	typedef unsigned long long ull;
	typedef __uint128_t LL;
	ull m,B;
	Barrett(const ull&m=2):m(m),B((LL(1)<<64)/m){}
	friend inline ull operator%(const ull&a,const Barrett&mod){
		ull r=a-mod.m*(LL(mod.B)*a>>64);return r>=mod.m?r-mod.m:r;
	}
}mod;
struct complex{
	db x,y;
	complex(const db&x=0,const db&y=0):x(x),y(y){}
	inline complex operator+(const complex&it)const{
		return complex(x+it.x,y+it.y);
	}
	inline complex operator-(const complex&it)const{
		return complex(x-it.x,y-it.y);
	}
	inline complex operator*(const complex&it)const{
		return complex(x*it.x-y*it.y,x*it.y+y*it.x);
	}
}buf[M<<2],*w[20];
inline int Getlen(const int&n){
	int len(0);while((1<<len)<n)++len;return len;
}
inline void swap(complex&a,complex&b){
	complex c=a;a=b;b=c;
}
inline void init(const int&n){
	const int&m=Getlen(n);complex*now=buf;w[m]=now;now+=1<<m;
	IMP(1<<m,w[m][i]=complex(std::cos(i*Pi/(1<<m)),std::sin(i*Pi/(1<<m))));
	for(int k=m-1;k>=0&&(w[k]=now,now+=1<<k);--k)IMP(1<<k,w[k][i]=w[k+1][i<<1]);
}
inline void DFT(complex*f,const int&M){
	const int&n=1<<M;
	for(int len=n>>1,d=M-1;d>=0;--d,len>>=1)for(int k=0;k^n;k+=len<<1){
		complex*W=w[d],*L=f+(k),*R=f+(k|len),x,y;IMP(len,(x=*L,y=*R)),*L++=(x+y),*R++=*W++*(x-y);
	}
}
inline void IDFT(complex*f,const int&M){
	const int&n=1<<M;
	for(int len=1,d=0;d^M;++d,len<<=1)for(int k=0;k^n;k+=len<<1){
		complex*W=w[d],*L=f+(k),*R=f+(k|len),x,y;IMP(len,(x=*L,y=*W++**R)),*L++=(x+y),*R++=(x-y);
	}
	IMP(n,(f[i].x/=n,f[i].y/=n));for(int i=1;(i<<1)<n;++i)swap(f[i],f[n-i]);
}
inline void MTT(int*f,int*g,const int&n){
	static complex T[M];const int&m=Getlen(n)+1;
	IMP(n,T[i]=complex(f[i],g[i]));DFT(T,m);IMP(1<<m,T[i]=T[i]*T[i]);IDFT(T,m);
	IMP(n,f[i]=((long long)(T[i].y+.5))/2%mod);IMP(1<<m,T[i]=complex());
}
inline void Inv(int*f,const int&n){
	static int b1[M],b2[M];static complex Q[M],P[M];const int&m=Getlen(n);b1[0]=1;
	for(int len=1;len<=m;++len){
		IMP(1<<len-1,b2[i]=2*b1[i]);IMP(1<<len,(Q[i].x=b1[i],P[i].x=f[i]));
		DFT(Q,len+1);DFT(P,len+1);IMP(1<<len+1,Q[i]=Q[i]*Q[i]*P[i]);IDFT(Q,len+1);
		IMP(1<<len,b1[i]=(::P+b2[i]-((long long)(Q[i].x+.5))%mod)%mod);IMP(1<<len+1,Q[i]=P[i]=complex());
	}
	IMP(n,f[i]=b1[i]);
}
inline void pow(int*f,const int&n,int b){
	static int g[M];g[0]=1;for(;b;b>>=1,MTT(f,f,n))if(b&1)MTT(g,f,n);IMP(n,f[i]=g[i]);
}
signed main(){
	scanf("%d%d%d%d%d%d",&n,&P,&m,&O,&S,&U);mod=Barrett(P);O=O%mod;S=S%mod;U=U%mod;++n;init(n<<1);
	const int&f0=(O+S+U)%mod,&f1=(3*P+O-S-2*U)%mod,&f2=U%mod;
	IMP(n,F[i]=(i+1)*(i+2)/2%mod);
	for(int i=n-1;i>=0;--i){
		F[i]=f0*F[i]%mod;if(i>=1)F[i]=(F[i]+f1*F[i-1])%mod;if(i>=2)F[i]=(F[i]+f2*F[i-2])%mod;
	}
	for(int i=n-1;i>=1;--i)G[i]=F[i]=F[i-1];F[0]=0;IMP(n,G[i]=(P-G[i])%mod);G[0]=1;
	if(m<n)pow(F,n,m+1);Inv(G,n);IMP(n,F[i]=(P-F[i])%mod);F[0]=1;
	if(m<n)MTT(G,F,n);printf("%d",G[n-1]);
}