代码

#include<bits/stdc++.h>
#define LL long long
#define LLL __int128
#define uint unsigned
#define ldb long double
#define uLL unsigned long long
using namespace std;
typedef vector<int> poly;
typedef vector<int> Vec;
typedef vector<Vec> Mat;
typedef tuple<poly,poly> Vec2;
typedef tuple<poly,poly,poly,poly> Mat2;
mt19937 rng(chrono::system_clock::now().time_since_epoch().count());
const int Mod=998244353,G=3,img=86583718;
const int winv=ULLONG_MAX%Mod+1,rinv=2-Mod;
const LL Mod2=(LL)Mod*Mod;
poly Grt({1}),Grt0,Grt1,frc({1,1}),inv({0,1}),ivf({1,1});vector<poly>Tr,Ts;
inline int fadd(const int&x,const int&y){int z=x+y;return z<0?z+Mod:z-Mod;}
inline int fsub(const int&x,const int&y){int z=x-y;return z<0?z+Mod:z-Mod;}
inline int reduce(const LL&x){return (x>>32)-(((LL)((int)(x*rinv))*Mod)>>32);}
inline int qpow(int x,int y,int z=1){
	for(;y;(y>>=1)&&(x=(LL)x*x%Mod))if(y&1)z=(LL)z*x%Mod;return z;
}
inline void dft(poly&P){
	const int n=P.size(),n2=n/2;
    int L=Grt.size();
	if(L*2<n)for(Grt.resize(n2);L*2<n;L<<=1){
		const int cr=qpow(G,Mod/(L<<2));
		for(int i=0;i<L;++i)Grt[i+L]=(LL)cr*Grt[i]%Mod;
	}
	if(Grt0.size()<Grt.size()){
		const int pr=Grt0.size(),ed=Grt.size();
		Grt0.resize(ed),Grt1.resize(ed);
		for(int i=pr;i<ed;++i)
			Grt0[i]=reduce((LL)winv*Grt[i]),Grt1[i]=(uint)rinv*Grt0[i];
	}
    int m=1,k=n;
	if(__lg(n)&1){
		const int R2=n2;
		for(int i=0;i<R2;++i){
			const int j=i+R2,v=fsub(P[i],P[j]);
			P[i]=fadd(P[i],P[j]),P[j]=v;
		}
        m*=2,k/=2;
	}
	for(;k>1;k/=4,m*=4)
        for(int i=0;i<m;++i){
            const int R=k,R4=R/4,g0=Grt1[i+i],g1=Grt0[i+i],g2=Grt1[i],g3=Grt0[i],g4=Grt1[i+i+1],g5=Grt0[i+i+1];
            int*P0=&P[i*R],*P1=P0+R4,*P2=P1+R4,*P3=P2+R4;
#pragma GCC unroll 8
            for(int j=0,z=0,m=0,z0=0,z1=0,z2=0,z3=0,z4=0,z5=0,z6=0,z7=0;j<R4;++j){
                m=(uint)g2*P2[j],z0=((LL)P2[j]*g3-(LL)m*Mod)>>32;
                m=(uint)g2*P3[j],z1=((LL)P3[j]*g3-(LL)m*Mod)>>32;
                z4=P0[j]+z0,z5=P0[j]-z0,z6=P1[j]+z1,z7=P1[j]-z1;
                m=(uint)g0*z6,z2=((LL)z6*g1-(LL)m*Mod)>>32;
                m=(uint)g4*z7,z3=((LL)z7*g5-(LL)m*Mod)>>32;
                P0[j]=fadd(z4,z2),P1[j]=fsub(z4,z2),P2[j]=fadd(z5,z3),P3[j]=fsub(z5,z3);
            }
        }
}
inline void idft(poly&P){
	const int m=P.size(),n=m,ni=reduce((LL)qpow(n,Mod-2)*winv),n2=n/2;
	poly rev(n);
	for(int i=1;i<n;++i)rev[i]=(rev[i>>1]>>1)+(i&1)*n2;
	for(int&i:P)i=reduce((LL)i*ni);
	for(int i=1;i<n;++i)if(i<rev[i])swap(P[i],P[rev[i]]);
	reverse(P.begin()+1,P.end());
	dft(P);
    for(int&i:P)(i<0&&(i+=Mod));
	for(int i=1;i<n;++i)if(i<rev[i])swap(P[i],P[rev[i]]);
}
inline poly Mul(poly P,poly Q){
    if(P.empty()||Q.empty())return poly();
	const int pn=P.size(),qn=Q.size(),rn=pn+qn-1;
	const int m=2<<__lg(max(1,rn-1));
	P.resize(m),Q.resize(m);
	dft(P),dft(Q);
	for(int i=m;i--;)P[i]=(LL)P[i]*Q[i]%Mod;
	idft(P);
	return P.resize(rn),P;
}
void poly_multiply(unsigned *a, int n, unsigned *b, int m, unsigned *c)
{
	poly P=Mul(poly(a,a+n+1),poly(b,b+m+1));
copy(P.begin(),P.end(),c);
}

评测结果