dit-fft的c++实现

以下是基于DIT(Divide in Time)算法的FFT(Fast Fourier Transform)的C++实现:

#include <bits/stdc++.h>
using namespace std;

const int N = 1 << 20;
const double PI = acos(-1.0);

struct Complex {
    double real, imag;
    Complex operator+(const Complex &rhs) const {
        return {real + rhs.real, imag + rhs.imag};
    }
    Complex operator-(const Complex &rhs) const {
        return {real - rhs.real, imag - rhs.imag};
    }
    Complex operator*(const Complex &rhs) const {
        return {real * rhs.real - imag * rhs.imag,
                real * rhs.imag + imag * rhs.real};
    }
} a[N], b[N], c[N], d[N];
int rev[N];

void FFT(Complex *a, int n, int inv) {
    for (int i = 0; i < n; ++i) {
        if (i < rev[i]) {
            swap(a[i], a[rev[i]]);
        }
    }
    for (int i = 2; i <= n; i <<= 1) {
        Complex wn = {cos(2 * PI / i), sin(2 * PI / i) * inv};
        for (int j = 0; j < n; j += i) {
            Complex w = {1, 0};
            for (int k = j; k < j + (i >> 1); ++k) {
                Complex x = a[k], y = w * a[k + (i >> 1)];
                a[k] = x + y, a[k + (i >> 1)] = x - y;
                w = w * wn;
            }
        }
    }
    if (inv == -1) {
        for (int i = 0; i < n; ++i) {
            a[i].real /= n, a[i].imag /= n;
        }
    }
}

int main() {
    int n = 1, m = 1;
    while (n <= 100000) {
        n <<= 1;
    }
    for (int i = 0; i < n; ++i) {
        rev[i] = (rev[i >> 1] >> 1) | ((i & 1) ? n >> 1 : 0);
    }
    for (int i = 0; i < m; ++i) {
        scanf("%lf", &a[i].real);
    }
    FFT(a, n, 1);
    for (int i = 0; i < n; ++i) {
        b[i] = a[i] * a[i];
    }
    FFT(b, n, -1);
    for (int i = 0; i < m; ++i) {
        printf("%.0f ", b[i].real);
    }
    return 0;
}

其中，`Complex`结构体表示复数，`FFT`函数为DIT-FFT算法的实现，`rev`数组为反转映射数组，`n`为最终计算的序列长度，`m`为输入序列的长度，`a`为输入序列，`b`为结果序列。该代码实现了对输入序列进行FFT变换后，计算其平方序列的逆变换，得到了平方序列。