任意模数NTT（二）

这一版是mx发明的MTT

速度极快，精度基本有保证，在奇技淫巧无效时可以考虑这个东西...

（但是无论如何我都不想用真正的任意模数NTT，那种东西简直毒瘤而且常常数巨大...）

原理：并不关心

#include <cstdio>
#include <cmath>
#include <cstring>
#include <cstdlib>
#include <iostream>
#include <algorithm>
#include <queue>
#include <stack>
#define ld long double
#define ll long long
using namespace std;
const ld pi=acos(-1.0);
const int Siz=(1<<21)+5;
const int siz=(1<<15)-1;
struct cp
{
    ld x,y;
    cp (){}
    cp (ld a,ld b):x(a),y(b){}
    friend cp operator + (cp a,cp b)
    {
        return cp(a.x+b.x,a.y+b.y);
    }
    friend cp operator - (cp a,cp b)
    {
        return cp(a.x-b.x,a.y-b.y);
    }
    friend cp operator * (cp a,cp b)
    {
        return cp(a.x*b.x-a.y*b.y,a.x*b.y+a.y*b.x);
    }
};
int to[Siz];
ll p;
int n,m;
ll a[Siz],b[Siz];
ll ret[Siz] ;
void FFT(cp *a,int len,int k)
{
    for(int i=0;i<len;i++)if(i<to[i])swap(a[i],a[to[i]]);
    for(int i=1;i<len;i<<=1)
    {
        cp w0=cp(cos(pi/i),k*sin(pi/i));
        for(int j=0;j<len;j+=(i<<1))
        {
            cp w=cp(1,0);
            for(int o=0;o<i;o++,w=w*w0)
            {
                cp w1=a[j+o],w2=a[j+o+i]*w;
                a[j+o]=w1+w2,a[j+o+i]=w1-w2;
            }
        }
    }
}
cp A[Siz],B[Siz],C[Siz],D[Siz],E[Siz],F[Siz],Q[Siz],G[Siz],H[Siz];
inline void MTT(ll *f,ll *g,ll *re,int len)
{
    int lim=1,l=0;
    while(lim<=len)lim<<=1,l++;
    for(int i=0;i<lim;i++)to[i]=((to[i>>1]>>1)|((i&1)<<(l-1)));
    for(int i=0;i<=len;i++)A[i]=cp((ld)(f[i]&siz),(ld)(f[i]>>15));
    for(int i=0;i<=len;i++)B[i]=cp((ld)(g[i]&siz),(ld)(g[i]>>15));
    FFT(A,lim,1),FFT(B,lim,1);
    for(int i=0;i<lim;i++)
    {
        int j=(lim-i)&(lim-1);
        C[j]=cp(0.5*(A[i].x+A[j].x),0.5*(A[i].y-A[j].y))*B[i];
        D[j]=cp(0.5*(A[i].y+A[j].y),0.5*(A[j].x-A[i].x))*B[i];
    }
    FFT(C,lim,1),FFT(D,lim,1);
    for(int i=0;i<lim;i++)
    {
        ll tempA=(ll)(C[i].x/lim+0.5)%p;
        ll tempa=(ll)(C[i].y/lim+0.5)%p;
        ll tempB=(ll)(D[i].x/lim+0.5)%p;
        ll tempb=(ll)(D[i].y/lim+0.5)%p;
        re[i]=(ll)(tempA+p+(ll)((ll)(tempa+tempB)<<15)%p+p+(ll)((ll)tempb<<30)%p+p)%p;
    }
}
int main()
{
    scanf("%d%d%lld",&n,&m,&p);
    for(int i=0;i<=n;i++)scanf("%lld",&a[i]);
    for(int i=0;i<=m;i++)scanf("%lld",&b[i]);
    MTT(a,b,ret,n+m-1);
    for(int i=0;i<=n+m;i++)printf("%lld ",ret[i]);
    printf("
");
    return 0;
}