FFT&&FWT&&NTT

来源:互联网 发布:淘宝哪家男鞋店好最潮 编辑:程序博客网 时间:2024/04/26 22:33

FFT是计算卷积的,就是

这里写图片描述
这里写图片描述
FFT大数乘法模版:
把一个大数看成一个x为10的 多项式, 两个大数乘法,就转化成多项式乘法

#include<cstdio>#include<cmath>#include<cstring>#include<algorithm>using namespace std;const int N = 500005;const double pi = acos(-1.0);char s1[N],s2[N];int len,res[N];struct Complex{    double r,i;    Complex(double r=0,double i=0):r(r),i(i) {};    Complex operator+(const Complex &rhs)    {        return Complex(r + rhs.r,i + rhs.i);    }    Complex operator-(const Complex &rhs)    {        return Complex(r - rhs.r,i - rhs.i);    }    Complex operator*(const Complex &rhs)    {        return Complex(r*rhs.r - i*rhs.i,i*rhs.r + r*rhs.i);    }} va[N],vb[N];void rader(Complex F[],int len) //len = 2^M,reverse F[i] with  F[j] j为i二进制反转{    int j = len >> 1;    for(int i = 1;i < len - 1;++i)    {        if(i < j) swap(F[i],F[j]);  // reverse        int k = len>>1;         while(j>=k)        {            j -= k;            k >>= 1;        }        if(j < k) j += k;    }}void FFT(Complex F[],int len,int t){    rader(F,len);    for(int h=2;h<=len;h<<=1)    {        Complex wn(cos(-t*2*pi/h),sin(-t*2*pi/h));        for(int j=0;j<len;j+=h)        {            Complex E(1,0); //旋转因子            for(int k=j;k<j+h/2;++k)            {                Complex u = F[k];                Complex v = E*F[k+h/2];                F[k] = u+v;                F[k+h/2] = u-v;                E=E*wn;            }        }    }    if(t==-1)   //IDFT        for(int i=0;i<len;++i)            F[i].r/=len;}void Conv(Complex a[],Complex b[],int len) //求卷积{    FFT(a,len,1);    FFT(b,len,1);    for(int i=0;i<len;++i) a[i] = a[i]*b[i];     FFT(a,len,-1);}void init(char *s1,char *s2){    int n1 = strlen(s1),n2 = strlen(s2);    len = 1;    while(len < 2*n1 || len < 2*n2) len <<= 1;    int i;    for(i=0;i<n1;++i)    {        va[i].r = s1[n1-i-1]-'0';        //因为从字符串i=0为高位,在多项式中为最高项,所以应该置于后面        va[i].i = 0;    }    while(i<len)    {        va[i].r = va[i].i = 0;        ++i;    }    for(i=0;i<n2;++i)    {        vb[i].r = s2[n2-i-1]-'0';        vb[i].i = 0;    }    while(i<len)    {        vb[i].r = vb[i].i = 0;        ++i;    }}void gao(){    Conv(va,vb,len);    memset(res,0,sizeof res);    for(int i=0;i<len;++i)    {        res[i]=va[i].r + 0.5;    }    for(int i=0;i<len;++i)    {        res[i+1]+=res[i]/10;        res[i]%=10;    }    int high = 0;    for(int i=len-1;i>=0;--i)       {        if(res[i])          {            high = i;            break;        }    }    for(int i=high;i>=0;--i) putchar('0'+res[i]);    puts("");}int main(){    while(scanf("%s %s",s1,s2)==2)    {        init(s1,s2);        gao();    }    return 0;}

NTT 即FFT的数论版本(整数)
NTT 模版
hiho 1388

#include <bits/stdc++.h>using namespace std;typedef long long ll;const int N=3e5+10;const ll mod=( 1ll << 47 ) * 7 * 4451 + 1 ;const ll g=3;ll x1[N],x2[N];ll mul( ll x, ll y ){    return ( x * y - ( long long ) ( x / ( long double ) mod * y + 1e-3 ) * mod + mod ) % mod ;}ll power ( ll a, ll b ){    ll res = 1, tmp = a ;    while ( b )    {        if ( b & 1 ) res = mul ( res, tmp ) ;        tmp = mul ( tmp, tmp ) ;        b >>= 1 ;    }    return res ;}void DFT ( ll y[], int n, bool rev ){    for ( int i = 1, j, t, k ; i < n ; ++ i )    {        for ( k = n >> 1, t = i, j = 0 ; k ; k >>= 1, t >>= 1 )        {            j = j << 1 | t & 1 ;        }        if ( i < j ) swap ( y[i], y[j] ) ;    }    for ( int s = 2, ds = 1 ; s <= n ; ds = s, s <<= 1 )    {        ll wn = power ( g, ( mod - 1 ) / s ) ;        if ( !rev ) wn = power ( wn, mod - 2 ) ;        for ( int k = 0 ; k < n ; k += s )        {            ll w = 1, t ;            for ( int i = k ; i < k + ds ; ++ i, w = mul ( w, wn ) )            {                y[i + ds] = ( y[i] - ( t = mul ( y[i + ds], w ) ) + mod ) % mod ;                y[i] = ( y[i] + t ) % mod ;            }        }    }}void NTT ( ll x1[], ll x2[], int n ){    DFT ( x1, n, 1 ) ;    DFT ( x2, n, 1 ) ;    for ( int i = 0 ; i < n ; ++ i ) x1[i] = mul ( x1[i], x2[i] ) ;    DFT ( x1, n, 0 ) ;    ll vn = power ( n, mod - 2 ) ;    for ( int i = 0 ; i < n ; ++ i ) x1[i] = mul ( x1[i], vn ) ;}ll s1[N],s2[N];int main(){    int t;    scanf("%d",&t);    while(t--)    {        int n;        scanf("%d",&n);        ll ans=0;        for(int i=0;i<n;i++)        {            scanf("%lld",&s1[i]);            ans+=s1[i]*s1[i];        }        for(int i=n-1;i>=0;i--)        {            scanf("%lld",&s2[i]);            ans+=s2[i]*s2[i];        }        int len=1;        while(len<2*n)len<<=1;        for(int i=0; i<n; i++)            x1[i]=s1[i];        for(int i=n; i<len; i++)            x1[i]=0;        for(int i=0; i<n; i++)            x2[i]=s2[i];        for(int i=n; i<len; i++)            x2[i]=0;        NTT(x1,x2,len);        ll res=0;        for(int i=n;i<len;i++)            res=max(res,x1[i-n]+x1[i]);        printf("%lld\n",ans-2*res);    }}

FWT 异或,与,或模版
FWT 计算∑f(x)*g(y) ,但这里是计算所有的满足x^y = n,异或也可以换成&,|这些运算符。

//以下代码有mod,无mod的话直接去掉mid 和 转化rev 即可class FWT{public:    void fwt(int a[],int n)    {        for(int d=1; d<n; d<<=1)            for(int m=d<<1,i=0; i<n; i+=m)                for(int j=0; j<d; j++)                {                    int x=a[i+j],y=a[i+j+d];                    a[i+j]=(x+y)%mod,a[i+j+d]=(x-y+mod)%mod;                    //xor:a[i+j]=x+y,a[i+j+d]=(x-y+mod)%mod;                    //and:a[i+j]=x+y;                    //or:a[i+j+d]=x+y;                }    }    void ufwt(int a[],int n)    {        for(int d=1; d<n; d<<=1)            for(int m=d<<1,i=0; i<n; i+=m)                for(int j=0; j<d; j++)                {                    int x=a[i+j],y=a[i+j+d];                    a[i+j]=1LL*(x+y)*rev%mod,a[i+j+d]=(1LL*(x-y)*rev%mod+mod)%mod;                    // /上一个数  即 *rev这个数的逆元                     //xor:a[i+j]=(x+y)/2,a[i+j+d]=(x-y)/2;                    //and:a[i+j]=x-y;                    //or:a[i+j+d]=y-x;                }    }    void solve(int a[],int b[],int n)    {        fwt(a,n);        fwt(b,n);        for(int i=0; i<n; i++) a[i]=1LL*a[i]*b[i]%mod;        ufwt(a,n);    }} myfwt;