在不使用递归的情况下将 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
}

它工作得很好，但我找到了一个更好的代码，它可以在没有递归的情况下完成相同的工作，并且运行速度也快得多。

我

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

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

下面是新代码

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实现已经被描述了数百次。

具有非递归方法的维基页面部分。

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

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

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