洛谷 P4245 【模板】任意模数NTT

https://www.luogu.com.cn/problem/P4245

#include<bits/stdc++.h>
using namespace std;        //FFT模板
typedef long long ll;

const int maxn=1e5+5;

int idx=0;//次数
int n,m,realmod;
int limit,bit,len;
int wz[maxn<<2];
ll MOD[3]={167772161,998244353,1004535809};
int G[3]={3,3,3};
ll a[maxn<<2],b[maxn<<2],c[3][maxn<<2];
ll beg1[maxn<<2],beg2[maxn<<2];

inline ll fmul(ll a,ll b,ll p)
{
    if(a>=p)
        a%=p;
    if(b>=p)
        b%=p;
    return (a*b-(ll)((long double)a*b/p)*p+p)%p;
}

inline ll qpow(ll x,ll y,ll p)
{
    ll t1=x,t2=1;
    while(y)
    {
        if(y&1)
            t2=(t1*t2)%p;
        t1=(t1*t1)%p;
        y>>=1;
    }
    return t2;
}

inline void CRT()
{
    ll m1=MOD[0],m2=MOD[1],m3=MOD[2];
    for(int i=0;i<=len;i++)
    {
        ll M=m1*m2;
        ll ans=fmul(c[0][i]*m2,qpow(m2,m1-2,m1),M)+fmul(c[1][i]*m1,qpow(m1,m2-2,m2),M);
        ans%=M;
        ll t1=((c[2][i]-ans)%m3+m3)%m3;
        ll t2=t1*qpow(M%m3,m3-2,m3)%m3;
        ans%=realmod,M%=realmod,t2%=realmod;
        c[0][i]=(M*t2%realmod+ans)%realmod;
    }
}

void NTT(ll *A,int inv) //中国剩余定理
{
    for(int i=0;i<limit;i++)
        if(i<wz[i])
            swap(A[i],A[wz[i]]);
    ll gn,t1,t2;
    for(int mid=1;mid<limit;mid<<=1)
    {
        gn=qpow(G[idx],(MOD[idx]-1)/(mid<<1),MOD[idx]);
        if(inv==-1)
            gn=qpow(gn,MOD[idx]-2,MOD[idx]);
        for(int i=0;i<limit;i+=mid<<1)
        {
            ll g=1;
            for(int j=0;j<mid;j++,g=g*gn%MOD[idx])
            {
                t1=A[i+j];
                t2=g*A[i+mid+j]%MOD[idx];
                A[i+j]=(t1+t2)%MOD[idx];
                A[i+mid+j]=(t1-t2+MOD[idx])%MOD[idx];
            }
        }
    }
}

void work()
{
    limit=1,bit=0,len=n+m;
    while(limit<=len)
    {
        limit<<=1;
        bit++;
    }
    for(int i=0;i<limit;i++)
        wz[i]=(wz[i>>1]>>1)|((i&1)<<(bit-1));
    while(idx<3)
    {
        memcpy(a,beg1,sizeof(ll)*limit);
        memcpy(b,beg2,sizeof(ll)*limit);
        NTT(a,1);
        NTT(b,1);
        for(int i=0;i<limit;i++)
            a[i]=a[i]*b[i]%MOD[idx];
        NTT(a,-1);
        ll inv=qpow(limit,MOD[idx]-2,MOD[idx]); //长度的逆元
        for(int i=0;i<=len;i++)
            c[idx][i]=a[i]*inv%MOD[idx];
        ++idx;
    }
    CRT();
}

int main()
{
    scanf("%d %d %d",&n,&m,&realmod);
    for(int i=0;i<=n;i++) //从低位到高位
        scanf("%lld",&beg1[i]);
    for(int i=0;i<=m;i++) //从低位到高位
        scanf("%lld",&beg2[i]);
    work();
    printf("%lld",c[0][0]);
    for(int i=1;i<=len;i++)
        printf(" %lld",c[0][i]);
    return 0;
}