FFT 学习笔记

#include <bits/stdc++.h>
using namespace std;
#define MAXN 100
const double PI = acos(-1);
/*快速傅立叶变换*/
struct FastFourierTransform {
    /*复数*/
    complex<double> omega[MAXN], omegaInverse[MAXN];
    inline void init(const int n) {
        for (register int i = 0; i < n; i++) {
            omega[i] = complex<double>(cos(2 * PI / n * i), sin(2 * PI / n * i));
            /*conj函数可以返回一个复数的共轭复数*/
            omegaInverse[i] = conj(omega[i]);
        }
    }
    inline void transform(complex<double> *a, const int n, const complex<double> *omega) {
        int k = 0;
        while ((1 << k) < n) k++;
        for (register int i = 0; i < n; i++) {
            int t = 0;
            for (register int j = 0; j < k; j++) if (i & (1 << j)) t |= (1 << (k - j - 1));
            if (i < t) swap(a[i], a[t]);
        }
        for (register int l = 2; l <= n; l *= 2) {
            int m = l / 2;
            for (complex<double> *p = a; p != a + n; p += l) {
                for (int i = 0; i < m; i++) {
                    complex<double> t = omega[n / l * i] * p[m + i];
                    p[m + i] = p[i] - t;
                    p[i] += t;
                }
            }
        }
    }
    inline void dft(complex<double> *a, const int n) {
        transform(a, n, omega);
    }
    inline void idft(complex<double> *a, const int n) {
        transform(a, n, omegaInverse);
        for (register int i = 0; i < n; i++) a[i] /= n;
    }
} fft;
inline void multiply(const int *a1, const int n1, const int *a2, const int n2, int *res) {
    int n = 1;
    while (n < n1 + n2) n <<= 1;
    static complex<double> c1[MAXN], c2[MAXN];
    for (register int i = 0; i < n1; i++) c1[i].real(a1[i]);
    for (register int i = 0; i < n2; i++) c2[i].real(a2[i]);
    fft.init(n);
    fft.dft(c1, n), fft.dft(c2, n);
    for (register int i = 0; i < n; i++) c1[i] *= c2[i];
    fft.idft(c1, n);
    for (register int i = 0; i < n1 + n2 - 1; i++) res[i] = static_cast<int>(floor(c1[i].real() + 0.5));
}
int main() {
    return 0;
}