Пользовательская реализация DFT возвращает ответ нечетного индекса в обратном порядке ⇐ C++

[Anonymous]

Это моя реализация DFT, тестируем с помощью {0, 1, 2, 3}.

Код: Выделить всё

#include 
#include 
#include 
#include 
#include 

std::vector DFT(std::vector&& P)
{
int n = P.size();
if (n == 1) {
return P;
}

std::vector Pe(n / 2);
std::vector Po(n / 2);
for (int i = 0; i < n / 2; ++i) {
Pe[i] = P[2 * i];
Po[i] = P[2 * i + 1];
}

auto ye = DFT(std::move(Pe)), yo = DFT(std::move(Po));

// in place algorithm, use input P to store output
auto wi = std::complex(1.0, 0.0);
auto wn = std::complex(std::cos(2.0 * std::numbers::pi_v / (double)n),
std::sin(2.0 * std::numbers::pi_v / (double)n));
for (int i = 0; i < n / 2; ++i) {
P[i]         = ye[i] + wi * yo[i];
P[i + n / 2] = ye[i] - wi * yo[i];
wi           = wi * wn;
}

return P;
}

int main()
{
int                               N = 4;
std::vector p_input(N);
for (int i = 0; i < N; ++i)
{
p_input[i] = {(double)i, 0.0};
}

auto p_output = DFT(std::move(p_input));

for (int i = 0; i < p_output.size(); ++i)
{
std::cout  fft(0:1:3)

ans =

6.0000 + 0.0000i  -2.0000 + 2.0000i  -2.0000 + 0.0000i  -2.0000 - 2.0000i

Я тестировал ввод большей длины, и мой результат всегда имеет обратный порядок в нечетной позиции индекса.
Я не знаю, какая часть неправильно? Моя справочная формула:
P(x): [p_0, p_1, ..., p_{n-1}]
w: [w^0, w^1, ..., w^{n-1}]
Pe(x^2): [p_0, p_2, ..., p_{ n-2}]
Po(x^2): [p_1, p_3, ..., p_{n-1}]
ye = [Pe(w^0), Pe(w^2), ..., Pe(w^{n-2})]
yo = [Po(w^ 0), Po(w^2), ..., Po(w^{n-2})]
P(w^j) = ye[j] + w^ j yo[j]
P(w^{j+n/2}) = ye[j] - w^j yo[j]
y = [P(w^0), P(w^1), ..., P(w^{n-1})]

Подробнее здесь: https://stackoverflow.com/questions/784 ... dex-answer