arnaucube
arnaucube/poulpyPublic Self-hosted repository
poulpy/spqlios/lib/spqlios/arithmetic/zn_vmp_int32_ref.c
spqlios/lib/spqlios/arithmetic/zn_vmp_int32_ref.cRaw
 1// This file is actually a template: it will be compiled multiple times with
 2// different INTTYPES
 3#ifndef INTTYPE
 4#define INTTYPE int32_t
 5#define INTSN i32
 6#endif
 7
 8#include <memory.h>
 9
10#include "zn_arithmetic_private.h"
11
12#define concat_inner(aa, bb, cc) aa##_##bb##_##cc
13#define concat(aa, bb, cc) concat_inner(aa, bb, cc)
14#define zn32_vec_fn(cc) concat(zn32_vec, INTSN, cc)
15
16// the ref version shares the same implementation for each fixed column size
17// optimized implementations may do something different.
18static __always_inline void IMPL_zn32_vec_matcols_ref(
19    const uint64_t NCOLS,            // fixed number of columns
20    uint64_t nrows,                  // nrows of b
21    int32_t* res,                    // result: size NCOLS, only the first min(b_sl, NCOLS) are relevant
22    const INTTYPE* a,                // a: nrows-sized vector
23    const int32_t* b, uint64_t b_sl  // b: nrows * min(b_sl, NCOLS) matrix
24) {
25  memset(res, 0, NCOLS * sizeof(int32_t));
26  for (uint64_t row = 0; row < nrows; ++row) {
27    int32_t ai = a[row];
28    const int32_t* bb = b + row * b_sl;
29    for (uint64_t i = 0; i < NCOLS; ++i) {
30      res[i] += ai * bb[i];
31    }
32  }
33}
34
35void zn32_vec_fn(mat32cols_ref)(uint64_t nrows, int32_t* res, const INTTYPE* a, const int32_t* b, uint64_t b_sl) {
36  IMPL_zn32_vec_matcols_ref(32, nrows, res, a, b, b_sl);
37}
38void zn32_vec_fn(mat24cols_ref)(uint64_t nrows, int32_t* res, const INTTYPE* a, const int32_t* b, uint64_t b_sl) {
39  IMPL_zn32_vec_matcols_ref(24, nrows, res, a, b, b_sl);
40}
41void zn32_vec_fn(mat16cols_ref)(uint64_t nrows, int32_t* res, const INTTYPE* a, const int32_t* b, uint64_t b_sl) {
42  IMPL_zn32_vec_matcols_ref(16, nrows, res, a, b, b_sl);
43}
44void zn32_vec_fn(mat8cols_ref)(uint64_t nrows, int32_t* res, const INTTYPE* a, const int32_t* b, uint64_t b_sl) {
45  IMPL_zn32_vec_matcols_ref(8, nrows, res, a, b, b_sl);
46}
47
48typedef void (*vm_f)(uint64_t nrows,                  //
49                     int32_t* res,                    //
50                     const INTTYPE* a,                //
51                     const int32_t* b, uint64_t b_sl  //
52);
53static const vm_f zn32_vec_mat8kcols_ref[4] = {  //
54    zn32_vec_fn(mat8cols_ref),                   //
55    zn32_vec_fn(mat16cols_ref),                  //
56    zn32_vec_fn(mat24cols_ref),                  //
57    zn32_vec_fn(mat32cols_ref)};
58
59/** @brief applies a vmp product (int32_t* input) */
60EXPORT void concat(default_zn32_vmp_apply, INTSN, ref)(  //
61    const MOD_Z* module,                                 //
62    int32_t* res, uint64_t res_size,                     //
63    const INTTYPE* a, uint64_t a_size,                   //
64    const ZN32_VMP_PMAT* pmat, uint64_t nrows, uint64_t ncols) {
65  const uint64_t rows = a_size < nrows ? a_size : nrows;
66  const uint64_t cols = res_size < ncols ? res_size : ncols;
67  const uint64_t ncolblk = cols >> 5;
68  const uint64_t ncolrem = cols & 31;
69  // copy the first full blocks
70  const uint32_t full_blk_size = nrows * 32;
71  const int32_t* mat = (int32_t*)pmat;
72  int32_t* rr = res;
73  for (uint64_t blk = 0;  //
74       blk < ncolblk;     //
75       ++blk, mat += full_blk_size, rr += 32) {
76    zn32_vec_fn(mat32cols_ref)(rows, rr, a, mat, 32);
77  }
78  // last block
79  if (ncolrem) {
80    uint64_t orig_rem = ncols - (ncolblk << 5);
81    uint64_t b_sl = orig_rem >= 32 ? 32 : orig_rem;
82    int32_t tmp[32];
83    zn32_vec_mat8kcols_ref[(ncolrem - 1) >> 3](rows, tmp, a, mat, b_sl);
84    memcpy(rr, tmp, ncolrem * sizeof(int32_t));
85  }
86  // trailing bytes
87  memset(res + cols, 0, (res_size - cols) * sizeof(int32_t));
88}