fft,ntt,fwt

Posted 2020-12-27 yinwuxiao

新搞了一个fft模板

比原来的快了4倍

1.自己写complex

2.利用mx的论文里的方法把另一半放在虚部

3.另外一个比较有点用的优化是

f[i+j+k]=a+k 当*a=f+i+j

#include <bits/stdc++.h>

using namespace std;

#define IL inline

#define rint register int

#define rep(i,h,t) for (rint i=h;i<=t;i++)

#define dep(i,t,h) for (rint i=t;i>=h;i--)

#define ll long long

const double pi = acos(-1.0);

char ss[1<<24],*A=ss,*B=ss;

IL char gc()

{

  return (A==B&&(B=(A=ss)+fread(ss,1,1<<24,stdin),A==B)?EOF:*A++);

}

template<class T>void read(T &x)

{

  rint f=1,c; while (c=gc(),c<48||c>57) if (c==‘-‘) f=-1; x=c^48;

  while (c=gc(),47<c&&c<58) x=(x<<3)+(x<<1)+(c^48); x*=f;

}
struct cp

{
    double a,b;
    cp operator +(const cp &o)const
    {
        return(cp){a+o.a,b+o.b};
    }
    cp operator -(const cp &o)const
   {
        return (cp){a-o.a,b-o.b};
    }
   cp operator *(const cp &o)const
    {
        return (cp){a*o.a-b*o.b,b*o.a+a*o.b};
    }
    cp operator *(const double &o)const
    {
        return (cp){a*o,b*o};
    }
    cp operator !()const
    {
        return (cp){a,-b};
    }
};
const int N = 3e6;
ll n, m, l, r[N+5], sum[N+5],ans[N+5];
cp a[N+5], b[N+5],a1[N+5],b1[N+5],w[N+5],z[N+5];
IL void fft(cp x[],int k,int o)
{
  cout<<k<<endl;
    for(int i=0,j=0; i<k; i++)
    {
        if(i>j)swap(x[i],x[j]);
        for(int l=k>>1; (j^=l)<l; l>>=1);
    }
    w[0]=(cp){1,0};
    for(int i=1;i<k;i<<=1)
    {
        cp g=(cp){cos(pi/i),sin(o*pi/i)};
       for (int j=1;j<i;j++) w[j]=w[j-1]*g;
        for(int j=0;j<k;j+=(i<<1))
        {
            cp *a=x+j,*b=a+i;
            for(int l=0;l<i;l++)
            {cp o=b[l]*w[l];b[l]=a[l]-o,a[l]=a[l]+o;}
        }
    }
    if(o==-1)for(int i=0;i<k;i++)x[i]=(cp){x[i].a/k,x[i].b/k};
}
char s1[N],s2[N];
IL void query()
{
    l=0;
    for (n = 1; n <= m/2; n <<= 1) l++;
      for (int i=0;i<n;i++) r[i]=(r[i/2]/2)|((i&1)<<(l-1));
    fft(a,n,1),
    fft(b,n,1);
    for(int i=0; i<n; i++)
    {
        int j=n-1&n-i;
        z[i]=(a[i]*b[i]*4-(a[i]-!a[j])*(b[i]-!b[j])*(((i&n>>1)?(cp){1,0}-w[i^n>>1]:w[i]+(cp){1,0})))*0.25;
    }
    fft(z,n,-1);
    for(int i=0; i<=m; i++)if(i&1)printf("%d ",((int)(z[i>>1].b+0.1)));
        else printf("%d ",((int)(z[i>>1].a+0.1)));
}
int main()
{
  freopen("1.in","r",stdin);
  freopen("1.out","w",stdout); 
  read(n); read(m);
  int x;
    for(int i=0; i<=n; i++)
    {
      read(x);
      i&1?a[i>>1].b=x:a[i>>1].a=x;
    }
    for(int i=0; i<=m; i++)
    {
      read(x);
      i&1?b[i>>1].b=x:b[i>>1].a=x;
    }
  m=n+m;
  query();
}