尝试对 FFT 进行递归编码的问题

Question

我最近一直在尝试对 FFT 进行编码，但未能获得有关示例问题的正确结果。由于我已经递归编码，调试一直很困难。这也是我第一次使用复杂的向量，所以我的问题很可能在于我对 valarrays 的声明。下面是我的代码（我使用模板的原因是因为我计划测量单精度和双精度来验证算法的稳定性）。

#include <iostream>
#include <complex>
#include "dftfft.h"
#include <valarray>

int main()
{
    std::valarray<std::complex<double>> V{ {-1,0}, {2,0}, {3,0}, {0,0} };
    std::valarray<std::complex<double>> VFFT(FFT<double>(V));
    for (size_t i = 0; i < VFFT.size(); i++)
    {
        std::cout << VFFT[i] << std::endl;
    }
}

在“dftfft.h”中：

template <typename T>
std::valarray<std::complex<T>> FFT(std::valarray<std::complex<T>> P)
{
    double pi = 2 * acos(0.0);
    size_t n = P.size();
    if (n == 1)
    {
        return P;
    }
    std::valarray<std::complex<T>> P_even = P[std::slice(0, n/2, 2)];
    std::valarray<std::complex<T>> P_odd = P[std::slice(1, n/2, 2)];
    std::valarray<std::complex<T>> y_even = FFT(P_even);
    std::valarray<std::complex<T>> y_odd = FFT(P_odd);
    std::valarray<std::complex<T>> y(n);
    for (size_t i = 0; i < n / 2; i++)
    {
        std::complex<T> omega = std::polar(1.0, 2.0 * pi / n);
        y[i] = P_even[i] + pow(omega, i) * P_odd[i];
        y[i + n / 2] = P_even[i] - pow(omega, i) * P_odd[i];
    }
    return y;
}

对于这个特定的向量 V，输出是：

(1,0) (3,0) (-3,0) (3,0)

我知道这是不正确的，但不知道为什么。正确答案是： (4,0) (-4,-2) (0,0) (-4,2)

我们将不胜感激。

Answer 1

解决了这个问题。在我的组合循环中，我需要使用 y_even/odd 而不是 P_even/odd。