COGS 2294. [HZOI 2015] 释迦

额，其实就是裸的三模数NTT，上一篇已经说过了

哦，还有一个就是对乘起来炸long long的数取模，用long double之类的搞一下就好，精度什么的，，（看出题人心情？？）

#include<cstdio>
#include<cstdlib>
#include<algorithm>
#include<cstring>
#define LL long long
#define N 300005
using namespace std;
inline int ra()
{
    int x=0; char ch=getchar();
    while (ch<'0' || ch>'9') ch=getchar();
    while (ch>='0' && ch<='9') {x=x*10+ch-'0'; ch=getchar();}
    return x;
}

const int mod=23333333;
const int M[]={998244353,1004535809,469762049};
const int G[]={3,3,3};
const LL _M=(LL)M[0]*M[1];

LL ksm(LL a, int b, int P)
{
    LL sum=1;
    for (;b;b>>=1,a=a*a%P)
        if (b&1) sum=sum*a%P;
    return sum;
}
/*LL mul(LL a, LL b, LL p)
{
    a%=p; b%=p;
    return ((a*b-(LL)((LL)((long double)a/p*b+1e-3)*p))%p+p)%p;
}
const int m1=M[0],m2=M[1],m3=M[2];
const int inv1=ksm(m1%m2,m2-2,m2),inv2=ksm(m2%m1,m1-2,m1),inv12=ksm(_M%m3,m3-2,m3);
int CRT(int a1, int a2, int a3)
{
    LL A=(mul((LL)a1*m2%_M,inv2,_M)+mul((LL)a2*m1%_M,inv1,_M))%_M;
    LL k=((LL)a3+m3-A%m3)*inv12%m3;
    return (k*(_M%mod)+A)%mod;
}*/

inline LL mul(LL a,LL b,LL p){
  a%=p; b%=p;
  return ((a*b-(LL)((LL)((long double)a/p*b+1e-3)*p))%p+p)%p;
}

const int m1=M[0],m2=M[1],m3=M[2];
const int inv1=ksm(m1%m2,m2-2,m2),inv2=ksm(m2%m1,m1-2,m1),inv12=ksm(_M%m3,m3-2,m3);
inline int CRT(int a1,int a2,int a3){
  LL A=(mul((LL)a1*m2%_M,inv2,_M)+mul((LL)a2*m1%_M,inv1,_M))%_M;
  LL k=((LL)a3+m3-A%m3)*inv12%m3;
  return (k*(_M%mod)+A)%mod;
}
int rev[N];
struct NTT{
    int num,P,G;
    int w[2][N];
    void pre(int _P, int _G, int n)
    {
        num=n; P=_P; G=_G;
        int g=ksm(G,(P-1)/num,P);
        w[1][0]=1; for (int i=1; i<n; i++) w[1][i]=(LL)w[1][i-1]*g%P;
        w[0][0]=1; for (int i=1; i<n; i++) w[0][i]=w[1][n-i];
    }
    void FFT(int *a, int n, int f)
    {
        for (int i=0; i<n; i++) if (i<rev[i]) swap(a[i],a[rev[i]]);
        for (int i=1; i<n; i<<=1)
            for (int j=0; j<n; j+=(i<<1))
                for (int k=0; k<i; k++)
                {
                    int x=a[j+k],y=(LL)w[f][num/(i<<1)*k]*a[i+j+k]%P;
                    a[j+k]=(x+y)%P; a[i+j+k]=(x-y+P)%P;
                }
        if (!f) for (int i=0,inv=ksm(n,P-2,P); i<n; i++) a[i]=(LL)a[i]*inv%P;
    }
}ntt[3];

int n,m;
int ans[3][N];
int a[N],b[N],c[N],d[N];

int main()
{
      freopen("annona_squamosa.in","r",stdin); freopen("annona_squamosa.out","w",stdout);
    n=ra();
    for (int i=0; i<n; i++) a[i]=ra();
    for (int i=0; i<n; i++) b[i]=ra();
    for (m=1; m<n<<1; m<<=1);
    int L=0,x=m>>1; while (!(x&1)) x>>=1,L++;
    for (int i=0; i<m; i++) rev[i]=(rev[i>>1]>>1)|((i&1)<<L);
//      for (m=1;m<=2*(n-1);m<<=1);
//    int L=0,x=m; while (x>>=1) L++;
//    for (int i=1;i<=m;i++) rev[i]=(rev[i>>1]>>1)|((i&1)<<(L-1));
//    for (int i=0; i<m; i++) printf("%d ",rev[i]); cout<<endl;
    for (int i=0; i<3; i++) ntt[i].pre(M[i],G[i],m);
    for (int i=0; i<3; i++)
    {
        memcpy(c,a,sizeof(int)*(m+5)); memcpy(d,b,sizeof(int)*(m+5));
        ntt[i].FFT(c,m,1); ntt[i].FFT(d,m,1);
        for (int j=0; j<m; j++) c[j]=(LL)c[j]*d[j]%ntt[i].P;
        ntt[i].FFT(c,m,0);
        for (int j=0; j<m; j++) ans[i][j]=c[j];
    }
    for (int i=0; i<n; i++) printf("%d ",CRT(ans[0][i],ans[1][i],ans[2][i]));
    return 0;
}