1#ifndef SPQLIOS_CPLX_FFT_INTERNAL_H
  2#define SPQLIOS_CPLX_FFT_INTERNAL_H
  3
  4#include "cplx_fft.h"
  5
  6/** @brief a complex number contains two doubles real,imag */
  7typedef double CPLX[2];
  8
  9EXPORT void cplx_set(CPLX r, const CPLX a);
 10EXPORT void cplx_neg(CPLX r, const CPLX a);
 11EXPORT void cplx_add(CPLX r, const CPLX a, const CPLX b);
 12EXPORT void cplx_sub(CPLX r, const CPLX a, const CPLX b);
 13EXPORT void cplx_mul(CPLX r, const CPLX a, const CPLX b);
 14
 15/**
 16 * @brief splits 2h evaluations of one polynomials into 2 times h evaluations of even/odd polynomial
 17 * Input:  Q_0(y),...,Q_{h-1}(y),Q_0(-y),...,Q_{h-1}(-y)
 18 * Output: P_0(z),...,P_{h-1}(z),P_h(z),...,P_{2h-1}(z)
 19 * where Q_i(X)=P_i(X^2)+X.P_{h+i}(X^2) and y^2 = z
 20 * @param h number of "coefficients" h >= 1
 21 * @param data 2h complex coefficients interleaved and 256b aligned
 22 * @param powom y represented as (yre,yim)
 23 */
 24EXPORT void cplx_split_fft_ref(int32_t h, CPLX* data, const CPLX powom);
 25EXPORT void cplx_bisplit_fft_ref(int32_t h, CPLX* data, const CPLX powom[2]);
 26
 27/**
 28 * Input: Q(y),Q(-y)
 29 * Output: P_0(z),P_1(z)
 30 * where Q(X)=P_0(X^2)+X.P_1(X^2) and y^2 = z
 31 * @param data 2 complexes coefficients interleaved and 256b aligned
 32 * @param powom (z,-z) interleaved: (zre,zim,-zre,-zim)
 33 */
 34EXPORT void split_fft_last_ref(CPLX* data, const CPLX powom);
 35
 36EXPORT void cplx_ifft_naive(const uint32_t m, const double entry_pwr, CPLX* data);
 37EXPORT void cplx_ifft16_avx_fma(void* data, const void* omega);
 38EXPORT void cplx_ifft16_ref(void* data, const void* omega);
 39
 40/**
 41 * @brief compute the ifft evaluations of P in place
 42 * ifft(data) = ifft_rec(data, i);
 43 * function ifft_rec(data, omega) {
 44 *   if #data = 1: return data
 45 *   let s = sqrt(omega) w. re(s)>0
 46 *   let (u,v) = data
 47 *   return split_fft([ifft_rec(u, s), ifft_rec(v, -s)],s)
 48 * }
 49 * @param itables precomputed tables (contains all the powers of omega in the order they are used)
 50 * @param data vector of m complexes (coeffs as input, evals as output)
 51 */
 52EXPORT void cplx_ifft_ref(const CPLX_IFFT_PRECOMP* itables, void* data);
 53EXPORT void cplx_ifft_avx2_fma(const CPLX_IFFT_PRECOMP* itables, void* data);
 54EXPORT void cplx_fft_naive(const uint32_t m, const double entry_pwr, CPLX* data);
 55EXPORT void cplx_fft16_avx_fma(void* data, const void* omega);
 56EXPORT void cplx_fft16_ref(void* data, const void* omega);
 57
 58/**
 59 * @brief compute the fft evaluations of P in place
 60 * fft(data) = fft_rec(data, i);
 61 * function fft_rec(data, omega) {
 62 *   if #data = 1: return data
 63 *   let s = sqrt(omega) w. re(s)>0
 64 *   let (u,v) = merge_fft(data, s)
 65 *   return [fft_rec(u, s), fft_rec(v, -s)]
 66 * }
 67 * @param tables precomputed tables (contains all the powers of omega in the order they are used)
 68 * @param data vector of m complexes (coeffs as input, evals as output)
 69 */
 70EXPORT void cplx_fft_ref(const CPLX_FFT_PRECOMP* tables, void* data);
 71EXPORT void cplx_fft_avx2_fma(const CPLX_FFT_PRECOMP* tables, void* data);
 72
 73/**
 74 * @brief merges 2 times h evaluations of even/odd polynomials into 2h evaluations of a sigle polynomial
 75 * Input:  P_0(z),...,P_{h-1}(z),P_h(z),...,P_{2h-1}(z)
 76 * Output: Q_0(y),...,Q_{h-1}(y),Q_0(-y),...,Q_{h-1}(-y)
 77 * where Q_i(X)=P_i(X^2)+X.P_{h+i}(X^2) and y^2 = z
 78 * @param h number of "coefficients" h >= 1
 79 * @param data 2h complex coefficients interleaved and 256b aligned
 80 * @param powom y represented as (yre,yim)
 81 */
 82EXPORT void cplx_twiddle_fft_ref(int32_t h, CPLX* data, const CPLX powom);
 83
 84EXPORT void citwiddle(CPLX a, CPLX b, const CPLX om);
 85EXPORT void ctwiddle(CPLX a, CPLX b, const CPLX om);
 86EXPORT void invctwiddle(CPLX a, CPLX b, const CPLX ombar);
 87EXPORT void invcitwiddle(CPLX a, CPLX b, const CPLX ombar);
 88
 89// CONVERSIONS
 90
 91/** @brief r = x from ZnX (coeffs as signed int32_t's ) to double */
 92EXPORT void cplx_from_znx32_ref(const CPLX_FROM_ZNX32_PRECOMP* precomp, void* r, const int32_t* x);
 93EXPORT void cplx_from_znx32_avx2_fma(const CPLX_FROM_ZNX32_PRECOMP* precomp, void* r, const int32_t* x);
 94/** @brief r = x to ZnX (coeffs as signed int32_t's ) to double */
 95EXPORT void cplx_to_znx32_ref(const CPLX_TO_ZNX32_PRECOMP* precomp, int32_t* r, const void* x);
 96EXPORT void cplx_to_znx32_avx2_fma(const CPLX_TO_ZNX32_PRECOMP* precomp, int32_t* r, const void* x);
 97/** @brief r = x mod 1 from TnX (coeffs as signed int32_t's) to double */
 98EXPORT void cplx_from_tnx32_ref(const CPLX_FROM_TNX32_PRECOMP* precomp, void* r, const int32_t* x);
 99EXPORT void cplx_from_tnx32_avx2_fma(const CPLX_FROM_TNX32_PRECOMP* precomp, void* r, const int32_t* x);
100/** @brief r = x mod 1 from TnX (coeffs as signed int32_t's) */
101EXPORT void cplx_to_tnx32_ref(const CPLX_TO_TNX32_PRECOMP* precomp, int32_t* x, const void* c);
102EXPORT void cplx_to_tnx32_avx2_fma(const CPLX_TO_TNX32_PRECOMP* precomp, int32_t* x, const void* c);
103/** @brief r = x from RnX (coeffs as doubles ) to double */
104EXPORT void cplx_from_rnx64_ref(const CPLX_FROM_RNX64_PRECOMP* precomp, void* r, const double* x);
105EXPORT void cplx_from_rnx64_avx2_fma(const CPLX_FROM_RNX64_PRECOMP* precomp, void* r, const double* x);
106/** @brief r = x to RnX (coeffs as doubles ) to double */
107EXPORT void cplx_to_rnx64_ref(const CPLX_TO_RNX64_PRECOMP* precomp, double* r, const void* x);
108EXPORT void cplx_to_rnx64_avx2_fma(const CPLX_TO_RNX64_PRECOMP* precomp, double* r, const void* x);
109/** @brief r = x to integers in RnX (coeffs as doubles ) to double */
110EXPORT void cplx_round_to_rnx64_ref(const CPLX_ROUND_TO_RNX64_PRECOMP* precomp, double* r, const void* x);
111EXPORT void cplx_round_to_rnx64_avx2_fma(const CPLX_ROUND_TO_RNX64_PRECOMP* precomp, double* r, const void* x);
112
113// fftvec operations
114/** @brief element-wise addmul r += ab */
115EXPORT void cplx_fftvec_addmul_ref(const CPLX_FFTVEC_ADDMUL_PRECOMP* precomp, void* r, const void* a, const void* b);
116EXPORT void cplx_fftvec_addmul_fma(const CPLX_FFTVEC_ADDMUL_PRECOMP* tables, void* r, const void* a, const void* b);
117EXPORT void cplx_fftvec_addmul_sse(const CPLX_FFTVEC_ADDMUL_PRECOMP* precomp, void* r, const void* a, const void* b);
118EXPORT void cplx_fftvec_addmul_avx512(const CPLX_FFTVEC_ADDMUL_PRECOMP* precomp, void* r, const void* a, const void* b);
119/** @brief element-wise mul r = ab */
120EXPORT void cplx_fftvec_mul_ref(const CPLX_FFTVEC_MUL_PRECOMP* tables, void* r, const void* a, const void* b);
121EXPORT void cplx_fftvec_mul_fma(const CPLX_FFTVEC_MUL_PRECOMP* tables, void* r, const void* a, const void* b);
122
123#endif  // SPQLIOS_CPLX_FFT_INTERNAL_H