FFT模板

港控/mmm° 2022-09-06 10:30 173阅读 0赞

# FFT模板 #

[P3803 【模板】多项式乘法（FFT） - 洛谷 | 计算机科学教育新生态 (luogu.com.cn)][P3803 _FFT_ - _ _ _ _luogu.com.cn]

## 系数表示法—>点值表示法—>系数表示法 ##

//后面的1表示要进行的变换是什么类型
        //1表示从系数变为点值
        //-1表示从点值变为系数 
    
    #include<iostream>
    #include<algorithm>
    #include<cmath>
    using namespace std;
    const int MAX_FFT = 1e7 + 10;
    const int MAX_AB = 4e6 + 10;//系数项的两倍
    
    inline int read()
    {
        char c = getchar(); int x = 0, f = 1;
        while (c < '0' || c>'9') { if (c == '-')f = -1; c = getchar(); }
        while (c >= '0' && c <= '9') { x = x * 10 + c - '0'; c = getchar(); }
        return x * f;
    }
    struct complex {//负数的加减乘
        double x, y;
        complex(double _x, double _y) { x = _x; y = _y; }
        complex() {}
        complex operator+(const complex& t) { return { x + t.x,y + t.y }; }
        complex operator-(const complex& t) { return { x - t.x,y - t.y }; }
        complex operator*(const complex& t) { return { x * t.x - y * t.y,x * t.y + y * t.x }; }
    };
    complex a[MAX_AB], b[MAX_AB];
    int n, m;
    struct  FFT//快速傅里叶变换
    {
        //存储
        double Pi = acos(-1.0);
        int l = 0, rev[MAX_FFT];
        int limit = 1;
        void init(int len)//len:为两个多项式乘积最高次幂为(degree(A)+degree(B))
        {
            while (limit <=len) limit <<= 1, l++;
            for (int i = 0; i < limit; i++)
                rev[i] = (rev[i >> 1] >> 1) | ((i & 1) << (l - 1));
        }
        //后面的type表示要进行的变换是什么类型
        //1表示从系数变为点值
        //-1表示从点值变为系数 
        void fast_fast_tle(complex* A, int type)
        {
         
    
            for (int i = 0; i < limit; i++)
                if (i < rev[i]) swap(A[i], A[rev[i]]);//求出要迭代的序列 
    
            for (int mid = 1; mid < limit; mid <<= 1)//待合并区间的中点
            {
                complex Wn(cos(Pi / mid), type * sin(Pi / mid)); //单位根 
                for (int R = mid << 1, j = 0; j < limit; j += R)//R是区间的右端点，j表示前已经到哪个位置了 
                {
                    complex w(1, 0);//幂 
                    for (int k = 0; k < mid; k++, w = w * Wn)//枚举左半部分 
                    {
                        complex x = A[j + k], y = w * A[j + mid + k];//蝴蝶效应 
                        A[j + k] = x + y;
                        A[j + mid + k] = x - y;
                    }
                }
            }
        }
    
    }fft;
    int main()
    {
        n = read(); m = read();
        for (int i = 0; i <= n; i++)a[i].x = read();
        for (int i = 0; i <= m; i++)b[i].x = read();
      
        fft.init(n+m);
    
        fft.fast_fast_tle(a, 1); //1表示从系数变为点值
        fft.fast_fast_tle(b, 1); //1表示从系数变为点值
    
        for (int i = 0; i <= fft.limit; i++)a[i] = a[i] * b[i];
        fft.fast_fast_tle(a, -1);  //-1表示从点值变为系数 
    
        for (int i = 0; i <= n + m; i++)printf("%d ", (int)(a[i].x / fft.limit + 0.5));
        return 0;
    }

[P3803 _FFT_ - _ _ _ _luogu.com.cn]: https://www.luogu.com.cn/problem/P3803