在不使用递归的情况下将FFT应用于两个非常大的数的乘法

Question

最近学习了FFT算法

我按照这个伪代码把它应用到超大自然数的快速乘法问题中，

Let A be array of length m, w be primitive m-th root of unity.
Goal: produce DFT F(A): evaluation of A at 1, w, w^2,...,w^{m-1}.

FFT(A, m, w)
{
  if (m==1) return vector (a_0)
  else {
    A_even = (a_0, a_2, ..., a_{m-2})
    A_odd  = (a_1, a_3, ..., a_{m-1})
    F_even = FFT(A_even, m/2, w^2)    //w^2 is a primitive m/2-th root of unity
    F_odd = FFT(A_odd, m/2, w^2)
    F = new vector of length m
    x = 1
    for (j=0; j < m/2; ++j) {
      F[j] = F_even[j] + x*F_odd[j]
      F[j+m/2] = F_even[j] - x*F_odd[j]
      x = x * w
  }
  return F
}

效果很好，但我发现了一个更好的代码，无需递归即可完成相同的工作，而且运行速度更快。

我试图逐行弄清楚它是如何工作的，但是我失败了。

如果你能详细解释一下前两个 for 循环（不是数学部分）发生了什么，我将不胜感激

下面是新代码

typedef complex<double> base;

void fft(vector<base> &a, bool invert)
{
    int n = a.size();

    for (int i = 1, j = 0; i < n; i++){
        int bit = n >> 1;
        for (; j >= bit; bit >>= 1) j -= bit;
        j += bit;

        if (i < j) swap(a[i], a[j]);
    }
    for (int len = 2; len <= n; len <<= 1){

        double ang = 2 * M_PI / len * (invert ? -1 : 1);
        base wlen(cos(ang), sin(ang));

        for (int i = 0; i < n; i += len){
            base w(1);
            for (int j = 0; j < len / 2; j++){
                base u = a[i + j], v = a[i + j + len / 2] * w;
                a[i + j] = u + v;
                a[i + j + len / 2] = u - v;
                w *= wlen;
            }
        }
    }
    if (invert)
    {
        for (int i = 0; i < n; i++)
            a[i] /= n;
    }
}

Answer 1

Cooley–Tukey FFT 实现已被描述数百次。

Wiki page part with non-recursive method.

第一个循环是位反转部分——代码重新打包源数组，将第 i 个索引处的元素与 i 的反转位索引交换（因此对于长度=8，索引 6=110b 与索引 [= 交换11=]，索引 5=101b 保持在同一位置）。

这种重新排序允许就地处理数组，对由 1,2,4,8... 索引（len/2 此处的步骤）与相应的三角系数分隔的对进行计算。

P.S。您的答案包含 onlinejudge 标记，因此这种紧凑的实现非常适合您的目的。但对于实际工作，值得使用一些高度优化的库，如 fftw 等

在不使用递归的情况下将FFT应用于两个非常大的数的乘法

Applying FFT on multiplication of two very large number without using recursion

c++

algorithm

fft

dft