1#include <assert.h>
  2#include <math.h>
  3#include <stdint.h>
  4
  5#include "q120_arithmetic.h"
  6#include "q120_common.h"
  7
  8EXPORT void q120_add_bbb_simple(uint64_t nn, q120b* const res, const q120b* const x, const q120b* const y) {
  9  const uint64_t* x_u64 = (uint64_t*)x;
 10  const uint64_t* y_u64 = (uint64_t*)y;
 11  uint64_t* res_u64 = (uint64_t*)res;
 12  for (uint64_t i = 0; i < 4 * nn; i += 4) {
 13    res_u64[i + 0] = x_u64[i + 0] % ((uint64_t)Q1 << 33) + y_u64[i + 0] % ((uint64_t)Q1 << 33);
 14    res_u64[i + 1] = x_u64[i + 1] % ((uint64_t)Q2 << 33) + y_u64[i + 1] % ((uint64_t)Q2 << 33);
 15    res_u64[i + 2] = x_u64[i + 2] % ((uint64_t)Q3 << 33) + y_u64[i + 2] % ((uint64_t)Q3 << 33);
 16    res_u64[i + 3] = x_u64[i + 3] % ((uint64_t)Q4 << 33) + y_u64[i + 3] % ((uint64_t)Q4 << 33);
 17  }
 18}
 19
 20EXPORT void q120_add_ccc_simple(uint64_t nn, q120c* const res, const q120c* const x, const q120c* const y) {
 21  const uint32_t* x_u32 = (uint32_t*)x;
 22  const uint32_t* y_u32 = (uint32_t*)y;
 23  uint32_t* res_u32 = (uint32_t*)res;
 24  for (uint64_t i = 0; i < 8 * nn; i += 8) {
 25    res_u32[i + 0] = (uint32_t)(((uint64_t)x_u32[i + 0] + (uint64_t)y_u32[i + 0]) % Q1);
 26    res_u32[i + 1] = (uint32_t)(((uint64_t)x_u32[i + 1] + (uint64_t)y_u32[i + 1]) % Q1);
 27    res_u32[i + 2] = (uint32_t)(((uint64_t)x_u32[i + 2] + (uint64_t)y_u32[i + 2]) % Q2);
 28    res_u32[i + 3] = (uint32_t)(((uint64_t)x_u32[i + 3] + (uint64_t)y_u32[i + 3]) % Q2);
 29    res_u32[i + 4] = (uint32_t)(((uint64_t)x_u32[i + 4] + (uint64_t)y_u32[i + 4]) % Q3);
 30    res_u32[i + 5] = (uint32_t)(((uint64_t)x_u32[i + 5] + (uint64_t)y_u32[i + 5]) % Q3);
 31    res_u32[i + 6] = (uint32_t)(((uint64_t)x_u32[i + 6] + (uint64_t)y_u32[i + 6]) % Q4);
 32    res_u32[i + 7] = (uint32_t)(((uint64_t)x_u32[i + 7] + (uint64_t)y_u32[i + 7]) % Q4);
 33  }
 34}
 35
 36EXPORT void q120_c_from_b_simple(uint64_t nn, q120c* const res, const q120b* const x) {
 37  const uint64_t* x_u64 = (uint64_t*)x;
 38  uint32_t* res_u32 = (uint32_t*)res;
 39  for (uint64_t i = 0, j = 0; i < 4 * nn; i += 4, j += 8) {
 40    res_u32[j + 0] = x_u64[i + 0] % Q1;
 41    res_u32[j + 1] = ((uint64_t)res_u32[j + 0] << 32) % Q1;
 42    res_u32[j + 2] = x_u64[i + 1] % Q2;
 43    res_u32[j + 3] = ((uint64_t)res_u32[j + 2] << 32) % Q2;
 44    res_u32[j + 4] = x_u64[i + 2] % Q3;
 45    res_u32[j + 5] = ((uint64_t)res_u32[j + 4] << 32) % Q3;
 46    res_u32[j + 6] = x_u64[i + 3] % Q4;
 47    res_u32[j + 7] = ((uint64_t)res_u32[j + 6] << 32) % Q4;
 48  }
 49}
 50
 51EXPORT void q120_b_from_znx64_simple(uint64_t nn, q120b* const res, const int64_t* const x) {
 52  static const int64_t MASK_HI = INT64_C(0x8000000000000000);
 53  static const int64_t MASK_LO = ~MASK_HI;
 54  static const uint64_t OQ[4] = {
 55      (Q1 - (UINT64_C(0x8000000000000000) % Q1)),
 56      (Q2 - (UINT64_C(0x8000000000000000) % Q2)),
 57      (Q3 - (UINT64_C(0x8000000000000000) % Q3)),
 58      (Q4 - (UINT64_C(0x8000000000000000) % Q4)),
 59  };
 60  uint64_t* res_u64 = (uint64_t*)res;
 61  for (uint64_t i = 0, j = 0; j < nn; i += 4, ++j) {
 62    uint64_t xj_lo = x[j] & MASK_LO;
 63    uint64_t xj_hi = x[j] & MASK_HI;
 64    res_u64[i + 0] = xj_lo + (xj_hi ? OQ[0] : 0);
 65    res_u64[i + 1] = xj_lo + (xj_hi ? OQ[1] : 0);
 66    res_u64[i + 2] = xj_lo + (xj_hi ? OQ[2] : 0);
 67    res_u64[i + 3] = xj_lo + (xj_hi ? OQ[3] : 0);
 68  }
 69}
 70
 71static int64_t posmod(int64_t x, int64_t q) {
 72  int64_t t = x % q;
 73  if (t < 0)
 74    return t + q;
 75  else
 76    return t;
 77}
 78
 79EXPORT void q120_c_from_znx64_simple(uint64_t nn, q120c* const res, const int64_t* const x) {
 80  uint32_t* res_u32 = (uint32_t*)res;
 81  for (uint64_t i = 0, j = 0; j < nn; i += 8, ++j) {
 82    res_u32[i + 0] = posmod(x[j], Q1);
 83    res_u32[i + 1] = ((uint64_t)res_u32[i + 0] << 32) % Q1;
 84    res_u32[i + 2] = posmod(x[j], Q2);
 85    res_u32[i + 3] = ((uint64_t)res_u32[i + 2] << 32) % Q2;
 86    res_u32[i + 4] = posmod(x[j], Q3);
 87    res_u32[i + 5] = ((uint64_t)res_u32[i + 4] << 32) % Q3;
 88    res_u32[i + 6] = posmod(x[j], Q4);
 89    res_u32[i + 7] = ((uint64_t)res_u32[i + 6] << 32) % Q4;
 90    ;
 91  }
 92}
 93
 94EXPORT void q120_b_to_znx128_simple(uint64_t nn, __int128_t* const res, const q120b* const x) {
 95  static const __int128_t Q = (__int128_t)Q1 * Q2 * Q3 * Q4;
 96  static const __int128_t Qm1 = (__int128_t)Q2 * Q3 * Q4;
 97  static const __int128_t Qm2 = (__int128_t)Q1 * Q3 * Q4;
 98  static const __int128_t Qm3 = (__int128_t)Q1 * Q2 * Q4;
 99  static const __int128_t Qm4 = (__int128_t)Q1 * Q2 * Q3;
100
101  const uint64_t* x_u64 = (uint64_t*)x;
102  for (uint64_t i = 0, j = 0; j < nn; i += 4, ++j) {
103    __int128_t tmp = 0;
104    tmp += (((x_u64[i + 0] % Q1) * Q1_CRT_CST) % Q1) * Qm1;
105    tmp += (((x_u64[i + 1] % Q2) * Q2_CRT_CST) % Q2) * Qm2;
106    tmp += (((x_u64[i + 2] % Q3) * Q3_CRT_CST) % Q3) * Qm3;
107    tmp += (((x_u64[i + 3] % Q4) * Q4_CRT_CST) % Q4) * Qm4;
108    tmp %= Q;
109    res[j] = (tmp >= (Q + 1) / 2) ? tmp - Q : tmp;
110  }
111}