arnaucube
arnaucube/poulpyPublic Self-hosted repository
poulpy/poulpy-backend/src/cpu_spqlios/spqlios-arithmetic/spqlios/arithmetic/vec_rnx_svp_ref.c
poulpy-backend/src/cpu_spqlios/spqlios-arithmetic/spqlios/arithmetic/vec_rnx_svp_ref.cRaw
 1#include <string.h>
 2
 3#include "../coeffs/coeffs_arithmetic.h"
 4#include "vec_rnx_arithmetic_private.h"
 5
 6EXPORT uint64_t fft64_bytes_of_rnx_svp_ppol(const MOD_RNX* module) { return module->n * sizeof(double); }
 7
 8EXPORT RNX_SVP_PPOL* new_rnx_svp_ppol(const MOD_RNX* module) { return spqlios_alloc(bytes_of_rnx_svp_ppol(module)); }
 9
10EXPORT void delete_rnx_svp_ppol(RNX_SVP_PPOL* ppol) { spqlios_free(ppol); }
11
12/** @brief prepares a svp polynomial  */
13EXPORT void fft64_rnx_svp_prepare_ref(const MOD_RNX* module,  // N
14                                      RNX_SVP_PPOL* ppol,     // output
15                                      const double* pol       // a
16) {
17  double* const dppol = (double*)ppol;
18  rnx_divide_by_m_ref(module->n, module->m, dppol, pol);
19  reim_fft(module->precomp.fft64.p_fft, dppol);
20}
21
22EXPORT void fft64_rnx_svp_apply_ref(                  //
23    const MOD_RNX* module,                            // N
24    double* res, uint64_t res_size, uint64_t res_sl,  // output
25    const RNX_SVP_PPOL* ppol,                         // prepared pol
26    const double* a, uint64_t a_size, uint64_t a_sl   // a
27) {
28  const uint64_t nn = module->n;
29  double* const dppol = (double*)ppol;
30
31  const uint64_t auto_end_idx = res_size < a_size ? res_size : a_size;
32  for (uint64_t i = 0; i < auto_end_idx; ++i) {
33    const double* a_ptr = a + i * a_sl;
34    double* const res_ptr = res + i * res_sl;
35    // copy the polynomial to res, apply fft in place, call fftvec
36    // _mul, apply ifft in place.
37    memcpy(res_ptr, a_ptr, nn * sizeof(double));
38    reim_fft(module->precomp.fft64.p_fft, (double*)res_ptr);
39    reim_fftvec_mul(module->precomp.fft64.p_fftvec_mul, res_ptr, res_ptr, dppol);
40    reim_ifft(module->precomp.fft64.p_ifft, res_ptr);
41  }
42
43  // then extend with zeros
44  for (uint64_t i = auto_end_idx; i < res_size; ++i) {
45    memset(res + i * res_sl, 0, nn * sizeof(double));
46  }
47}