AB fft

朱雀 2022-09-06 14:20 369阅读 0赞

# A\*B fft #

[P1919 【模板】A\*B Problem升级版（FFT快速傅里叶） - 洛谷 | 计算机科学教育新生态 (luogu.com.cn)][P1919 _A_B Problem_FFT_ - _ _ _ _luogu.com.cn]

參考代碼

#include<iostream>
    #include<algorithm>
    #include<vector>
    #include<cmath>
    #include<cstdio>
    #include<cstring>
    using namespace std;
    const int MAX_FFT = 2e7 + 10;
    const int MAX_AB = 2e7 + 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;
    string s1, s2;
    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 ans[MAX_FFT];
    int main()
    {
        ios::sync_with_stdio(false);
        cin >> s1 >> s2; n = s1.size() - 1; m = s2.size() - 1;
        for (int i = n, j = 0; i >= 0; i--, j++)a[j].x = s1[i] - '0';
        for (int i = m, j = 0; i >= 0; i--, j++)b[j].x = s2[i] - '0';
       int len = n + m+2;
        fft.init(len);
    
        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++)
        {
            ans[i] += (a[i].x / fft.limit + 0.49);
            ans[i + 1] += ans[i] / 10;
            ans[i] %= 10;
        }
      
        while (!ans[len] && len > 0)len--;
        for (int i = len; i >= 0; i--)cout << ans[i]; cout << '\n';
        //for (int i = 0; i <= n + m; i++)printf("%d ", (int)(a[i].x / fft.limit + 0.5));
        return 0;
    }

[P1919 _A_B Problem_FFT_ - _ _ _ _luogu.com.cn]: https://www.luogu.com.cn/problem/P1919