• 模板


    #include <bits/stdc++.h>
    using namespace std;
    typedef long long ll;
    
    const int MAXN = 4e6;
    const double PI = acos(-1.0);
    
    struct Complex {
        double x, y;
        Complex() {}
        Complex(double x, double y): x(x), y(y) {}
        friend Complex operator+(const Complex &a, const Complex &b) {
            return Complex(a.x + b.x, a.y + b.y);
        }
        friend Complex operator-(const Complex &a, const Complex &b) {
            return Complex(a.x - b.x, a.y - b.y);
        }
        friend Complex operator*(const Complex &a, const Complex &b) {
            return Complex(a.x * b.x - a.y * b.y, a.x * b.y + a.y * b.x);
        }
    } A[MAXN + 5], B[MAXN + 5];
    
    void FFT(Complex a[], int n, int op) {
        for(int i = 1, j = n >> 1; i < n - 1; ++i) {
            if(i < j)
                swap(a[i], a[j]);
            int k = n >> 1;
            while(k <= j) {
                j -= k;
                k >>= 1;
            }
            j += k;
        }
        for(int len = 2; len <= n; len <<= 1) {
            Complex wn(cos(2.0 * PI / len), sin(2.0 * PI / len)*op);
            for(int i = 0; i < n; i += len) {
                Complex w(1.0, 0.0);
                for(int j = i; j < i + (len >> 1); ++j) {
                    Complex u = a[j], t = a[j + (len >> 1)] * w ;
                    a[j] = u + t, a[j + (len >> 1)] = u - t;
                    w = w * wn;
                }
            }
        }
        if(op == -1) {
            for(int i = 0; i < n; ++i)
                a[i].x = (int)(a[i].x / n + 0.5);
        }
    }
    
    int pow2(int x) {
        int res = 1;
        while(res < x)
            res <<= 1;
        return res;
    }
    
    void convolution(Complex A[], Complex B[], int Asize, int Bsize) {
        int n = pow2(Asize + Bsize - 1);
        for(int i = 0; i < n; ++i) {
            A[i].y = 0.0;
            B[i].y = 0.0;
        }
        for(int i = Asize; i < n; ++i)
            A[i].x = 0;
        for(int i = Bsize; i < n; ++i)
            B[i].x = 0;
        FFT(A, n, 1);
        FFT(B, n, 1);
        for(int i = 0; i < n; ++i)
            A[i] = A[i] * B[i];
        FFT(A, n, -1);
        return;
    }
    
    int main() {
    #ifdef Yinku
        freopen("Yinku.in", "r", stdin);
    #endif // Yinku
        int n, m;
        scanf("%d%d", &n, &m);
        for(int i = 0; i <= n; ++i) {
            scanf("%lf", &A[i].x);
        }
        for(int i = 0; i <= m; ++i) {
            scanf("%lf", &B[i].x);
        }
        convolution(A, B, n + 1, m + 1);
        for(int i = 0; i <= n + m; i++) {
            printf("%d%c", (int)A[i].x, " 
    "[i == n + m]);
        }
        return 0;
    }
    

    有一个精度更高(在多次调用时速度也更快)的版本,需要预处理单位根,花费多一点空间:

    #include <bits/stdc++.h>
    using namespace std;
    typedef long long ll;
    
    const double PI = acos(-1.0);
    const int MAXN = 1e6;
    
    struct Complex {
        double x, y;
        Complex(): x(0), y(0) {}
        Complex(double x, double y): x(x), y(y) {}
        friend Complex operator+(const Complex &a, const Complex &b) {
            return Complex(a.x + b.x, a.y + b.y);
        }
        friend Complex operator-(const Complex &a, const Complex &b) {
            return Complex(a.x - b.x, a.y - b.y);
        }
        friend Complex operator*(const Complex &a, const Complex &b) {
            return Complex(a.x * b.x - a.y * b.y, a.x * b.y + a.y * b.x);
        }
    };
    
    typedef vector<Complex> Poly;
    
    Poly w[2 * MAXN + 5][2];
    int rev[4 * MAXN + 5];
    
    inline void FFT(Poly &a, int n, int op) {
        for(int i = 0; i < n ; ++i) {
            if(i < rev[i])
                swap(a[i], a[rev[i]]);
        }
        for(int len = 2; len <= n; len <<= 1) {
            register int m = len >> 1;
            for(int i = 0; i < m; ++i) {
                Complex &tw = w[m][op == 1][i];
                for(int j = i; j < n; j += len) {
                    Complex u = a[j], t = a[j + m] * tw ;
                    a[j] = u + t, a[j + m] = u - t;
                }
            }
        }
        if(op == -1) {
            for(int i = 0; i < n; ++i) {
                a[i].x = (ll)(a[i].x / n + 0.5);
                a[i].y = 0;
            }
        }
    }
    
    inline int pow2(int x, int &lgn) {
        int res = 1;
        lgn = 0;
        while(res < x) {
            if(!w[res][0].size()) {
                w[res][0].resize(res);
                w[res][1].resize(res);
                for(int i = 0; i < res; ++i) {
                    //0是逆变换需要的
                    w[res][0][i] = Complex(cos(-PI * i  / res), sin(-PI * i / res));
                    w[res][1][i] = Complex(cos(PI * i  / res), sin(PI * i  / res));
                }
            }
            res <<= 1;
            ++lgn;
        }
        return res;
    }
    
    inline int convolution(Poly &A, Poly &B, int Asize, int Bsize) {
        int lgn, n = pow2(Asize + Bsize - 1, lgn);
        for(int i = 0; i < n; ++i)
            rev[i] = (rev[i >> 1] >> 1) | ((i & 1) << lgn - 1);
        A.resize(n);
        B.resize(n);
        FFT(A, n, 1);
        FFT(B, n, 1);
        for(int i = 0; i < n; ++i)
            A[i] = A[i] * B[i];
        FFT(A, n, -1);
        while(n && (A[n - 1].x == 0))
            --n;
        return n;
    }
    
    Poly poly[MAXN + 5];
    
    
    int main() {
    #ifdef Yinku
        freopen("Yinku.in", "r", stdin);
    #endif // Yinku
        int n, m;
        scanf("%d%d", &n, &m);
        Poly A, B;
        for(int i = 0; i <= n; ++i) {
            int x;
            scanf("%d", &x);
            A.push_back(Complex(x, 0));
        }
        for(int i = 0; i <= m; ++i) {
            int x;
            scanf("%d", &x);
            B.push_back(Complex(x, 0));
        }
    
        int C = convolution(A, B, n + 1, m + 1);
        for(int i = 0; i <= n + m; ++i) {
            printf("%lld%c", (ll)A[i].x, " 
    "[i == n + m]);
        }
        return 0;
    }
    
  • 相关阅读:
    读写锁操作(ReaderWriterLockSlim)
    VirtualBox的小秘密:命令行
    云的始祖概念,认识Linux瘦客户机
    Flash ActionScript 3.0 通过asp.net 访问 数据库
    js刷新iframe框架的几种情况分析
    Mozilla两款火狐插件包含恶意代码被紧急喊停
    asp.net中DataBinder.Eval的用法总结
    实现firebird的Embedded模式(.net 3.5)
    Flash Player 9 支持H.264视频和aac音频(附官方代码)
    右下角浮动广告代码DEMO
  • 原文地址:https://www.cnblogs.com/Inko/p/11503377.html
Copyright © 2020-2023  润新知