BZOJ4836(CDQ分治+FFT（NTT会T））

题面
题目定义了一种奇怪的运算，是这样的

还给出了两个数组A,B，数组长度和元素大小都<=5*10^4。
有多个询问，给出一个k，问有多少对(i,j)使得
Ai opt Bj = k。

如果 x opt y = x+y，那么显然是一个桶状数组的卷积，FFT可以解决
若 x opt y = x-y，根据FFT题的套路，把一个数组反过来，发现它依然是个卷积。而题目把运算按x与y的大小关系分成了两种，两种都很类似，可以独立地用相同的方法做。

这题是一个不熟FFT的大佬搜CDQ分治搜到的，让我非常感兴趣，所以就来考虑分治。
其实也很简单，不过思路我没听过而已。
对于权值区间[l,r],我们二分一个mid，然后用FFT将[l,mid]和[mid+1,r]做个卷积，翻转一边，再做个卷积，统计答案。正确地清零，保证复杂度。再递归处理两个区间，然后就搞定了。
时间复杂度 O(N (log N)^2)
N=5*10^4。

为什么NTT在电脑上跑得飞快，交上去老是T。

#include <iostream>#include <fstream>#include <algorithm>#include <cmath>#include <ctime>#include <cstdio>#include <cstdlib>#include <cstring>#include <queue>using namespace std;#define mmst(a, b) memset(a, b, sizeof(a))#define mmcp(a, b) memcpy(a, b, sizeof(b))typedef long long LL;void read(int &hy){    hy=0;    char cc=getchar();    while(cc>'9'||cc<'0')    cc=getchar();    while(cc>='0'&&cc<='9')    {        hy=(hy<<3)+(hy<<1)+cc-'0';        cc=getchar();    }}const int N=400400;const double pi=acos(-1);struct yy{    double x,y;    yy(double a=0,double b=0):x(a),y(b){}};yy operator +(yy a,yy b) {return yy(a.x+b.x,a.y+b.y);}yy operator -(yy a,yy b) {return yy(a.x-b.x,a.y-b.y);}yy operator *(yy a,yy b) {return yy(a.x*b.x-a.y*b.y,a.x*b.y+a.y*b.x);}int n,rev[N];yy omega[N][2];void init(int lim){    int k=-1;    n=1;    while(n<=lim)    n<<=1,k++;    for(int i=0;i<n;i++)    rev[i]=(rev[i>>1]>>1) | ((i&1)<<k);    for(int i=0;i<n;i++)    {        omega[i][1]=yy(cos(2*pi*i/n),sin(2*pi*i/n));        omega[i][0]=yy(cos(2*pi*i/n),-sin(2*pi*i/n));    }}void fft(yy *a,int ops){    for(int i=0;i<n;i++)    if(i<rev[i])    swap(a[i],a[rev[i]]);    for(int l=2;l<=n;l<<=1)    {        int m=l>>1;        for(int i=0;i<n;i+=l)        for(int k=0;k<m;k++)        {            yy t=a[i+k+m]*omega[n/l*k][ops];            a[i+k+m]=a[i+k]-t;            a[i+k]=a[i+k]+t;        }    }    if(!ops)    for(int i=0;i<n;i++)    a[i].x=a[i].x/n;}int T,mx,v,m,q,by;yy aa[N],bb[N];LL c[N],d[N],ans[N];void work(int l1,int r1,int l2,int r2,int ops){    init(r1-l1+r2-l2+2);    for(int i=0;i<=r1-l1;i++)    aa[i].x=c[i+l1];    if(ops)    for(int i=0;i<=r2-l2;i++)    bb[i].x=d[i+l2];    else    for(int i=0;i<=r2-l2;i++)    bb[r2-l2-i].x=d[i+l2];    fft(aa,1);    fft(bb,1);    for(int i=0;i<n;i++)    aa[i]=aa[i]*bb[i];    fft(aa,0);    for(int i=0;i<n;i++)    if(ops)    ans[i+l1+l2]+=floor(aa[i].x+0.3);    else    ans[l1-r2+i]+=floor(aa[i].x+0.3);    for(int i=0;i<n;i++)    aa[i].x=aa[i].y=bb[i].x=bb[i].y=0;}void cdq(int l,int r){    if(l==r)    {        ans[0]+=c[l]*d[l];        return;    }    int mid=(l+r)/2;    cdq(l,mid);    cdq(mid+1,r);    work(l,mid,mid+1,r,1);    work(mid+1,r,l,mid,0);}int main(){    cin>>T;    while(T--)    {    mmst(c,0);    mmst(d,0);    mmst(ans,0);    mx=0;    cin>>v>>m>>q;    for(int i=1;i<=v;i++)    {        read(by);        c[by]++;        mx=max(mx,by);    }    for(int i=1;i<=m;i++)    {        read(by);        d[by]++;        mx=max(mx,by);    }    cdq(0,mx);    while(q--)    {        read(by);        printf("%lld\n",ans[by]);    }    }    return 0;}