fft.cpp

#include "benchmark_base/fft.h"

namespace comparison_benchmarks {
namespace base {
namespace fft {

void fill_input(fft::complex_vector &data) {
  for (size_t i = 0; i < data.size(); i++) {
    data[i] = std::complex<double>(sin(i), 0.0);
  }
}

void divide(complex_vector::iterator data, complex_vector::iterator tmp_odd_elements, int n) {
  for (int i = 0; i < n / 2; i++) {
    tmp_odd_elements[i] = data[i * 2 + 1];
  }
  for (int i = 0; i < n / 2; i++) {
    data[i] = data[i * 2];
  }
  for (int i = 0; i < n / 2; i++) {
    data[i + n / 2] = tmp_odd_elements[i];
  }
}

void combine(complex_vector::iterator data, int n) {
  for (int i = 0; i < n / 2; i++) {
    std::complex<double> even = data[i];
    std::complex<double> odd = data[i + n / 2];

    // w is the "twiddle-factor".
    // this could be cached, but we run the same 'base' algorithm parallel/serial,
    // so it won't impact the performance comparison.
    std::complex<double> w = exp(std::complex<double>(0, -2. * M_PI * i / n));

    data[i] = even + w * odd;
    data[i + n / 2] = even - w * odd;
  }
}

void conquer(complex_vector::iterator data, complex_vector::iterator swap_array, int n) {
  if (n < 2) {
    return;
  }

  divide(data, swap_array, n);
  conquer(data, swap_array, n / 2);
  conquer(data + n / 2, swap_array + n / 2, n / 2);
  combine(data, n);
}

}
}
}