洛谷 P4245 [模板]任意模数NTT —— 三模数NTT

#include<iostream>
#include<cstdio>
#include<cstring>
#include<algorithm>
using namespace std;
typedef long long ll;
int const xn=(1<<19);
int n,m,lim,rev[xn],a[5][xn],b[5][xn],p[5]={0,469762049,998244353,1004535809};
ll mod;
int rd()
{
  int ret=0,f=1; char ch=getchar();
  while(ch<‘0‘||ch>‘9‘){if(ch==‘-‘)f=0; ch=getchar();}
  while(ch>=‘0‘&&ch<=‘9‘)ret=ret*10+ch-‘0‘,ch=getchar();
  return f?ret:-ret;
}
ll upt(ll x,int md){while(x>=md)x-=md; while(x<0)x+=md; return x;}
ll mul(ll a,ll b,int md)
{
  ll ret=0; a=a%md; b=b%md;
  if(a<0)a+=md; if(b<0)b+=md;//
  for(;b;b>>=1ll,a=(a+a)%md)if(b&1)ret=(ret+a)%md;
  return ret;
}
ll pw(ll a,int b,int md)
{
  ll ret=1; a=a%md;
  for(;b;b>>=1,a=mul(a,a,md))if(b&1)ret=mul(ret,a,md);//mul!!
  return ret;
}
void ntt(int *a,int tp,int md)
{
  for(int i=0;i<lim;i++)
    if(i<rev[i])swap(a[i],a[rev[i]]);
  for(int mid=1;mid<lim;mid<<=1)
    {
      int wn=pw(3,(md-1)/(mid<<1),md);
      if(tp==-1)wn=pw(wn,md-2,md);
      for(int j=0,len=(mid<<1);j<lim;j+=len)
    {
      int w=1;
      for(int k=0;k<mid;k++,w=(ll)w*wn%md)
        {
          int x=a[j+k],y=(ll)w*a[j+mid+k]%md;
          a[j+k]=upt(x+y,md); a[j+mid+k]=upt(x-y,md);
        }
    }
    }
  if(tp==1)return; int inv=pw(lim,md-2,md);
  for(int i=0;i<lim;i++)a[i]=(ll)a[i]*inv%md;
}
ll uni(ll r1,ll r2,ll m1,int m2,int tp,int inv)
{
  ll k=mul(r2-r1,inv,m2);
  if(!tp)return (r1+k*m1)%(m1*m2);
  return upt((r1+mul(k,m1,mod))%mod,mod);//%mod!!
}
int main()
{
  n=rd(); m=rd(); mod=rd();
  for(int i=0;i<=n;i++)a[1][i]=a[2][i]=a[3][i]=rd();
  for(int i=0;i<=m;i++)b[1][i]=b[2][i]=b[3][i]=rd();
  lim=1; int l=0;
  while(lim<=n+m+2)lim<<=1,l++;//+2!
  for(int i=0;i<lim;i++)
    rev[i]=((rev[i>>1]>>1)|((i&1)<<(l-1)));
  for(int i=1;i<=3;i++)
    {
      ntt(a[i],1,p[i]); ntt(b[i],1,p[i]);
      for(int j=0;j<lim;j++)a[i][j]=(ll)a[i][j]*b[i][j]%p[i];
      ntt(a[i],-1,p[i]);
    }
  int inv1=pw(p[1],p[2]-2,p[2]);
  int inv2=pw((ll)p[1]*p[2],p[3]-2,p[3]);//inv!!
  for(int i=0;i<=n+m;i++)
    {
      ll ans=uni(a[1][i],a[2][i],p[1],p[2],0,inv1);
      ans=uni(ans,a[3][i],(ll)p[1]*p[2],p[3],1,inv2);//(ll)!!!
      printf("%lld ",ans);
    }
  puts("");
  return 0;
}