src/modinteger/cl_MI_pow2p1.h

   1 // m > 0, m = 2^m1 + 1 (m1 > 1)
   2
   3 namespace cln {
   4
   5 class cl_heap_modint_ring_pow2p1 : public cl_heap_modint_ring {
   6         SUBCLASS_cl_heap_modint_ring()
   7 public:
   8         // Constructor.
   9         cl_heap_modint_ring_pow2p1 (const cl_I& m, uintL m1); // m = 2^m1 + 1
  10         // Virtual destructor.
  11         ~cl_heap_modint_ring_pow2p1 () {}
  12         // Additional information.
  13         uintL m1;
  14 };
  15
  16 static inline const cl_I pow2p1_reduce_modulo (cl_heap_modint_ring* _R, const cl_I& x)
  17 {
  18         var cl_heap_modint_ring_pow2p1* R = (cl_heap_modint_ring_pow2p1*)_R;
  19         // Method:
  20         // If x>=0, split x into pieces of m1 bits and sum them up.
  21         //   x = x0 + 2^m1*x1 + 2^(2*m1)*x2 + ... ==>
  22         //   mod(x,m) = mod(x0-x1+x2-+...,m).
  23         // If x<0, apply this to -1-x, and use mod(x,m) = m-1-mod(-1-x,m).
  24  {      Mutable(cl_I,x);
  25         var bool sign = minusp(x);
  26         if (sign) { x = lognot(x); }
  27         var const uintL m1 = R->m1;
  28         while (x >= R->modulus) {
  29                 var uintL xlen = integer_length(x);
  30                 var cl_I y = ldb(x,cl_byte(m1,0));
  31                 for (var uintL i = m1; ; ) {
  32                         y = y - ldb(x,cl_byte(m1,i));
  33                         i += m1;
  34                         if (i >= xlen)
  35                                 break;
  36                         y = y + ldb(x,cl_byte(m1,i));
  37                         i += m1;
  38                         if (i >= xlen)
  39                                 break;
  40                 }
  41                 if (minusp(y))
  42                         { sign = !sign; x = lognot(y); }
  43                 else
  44                         x = y;
  45         }
  46         // Now 0 <= x < m.
  47         if (sign) { x = R->modulus - 1 - x; }
  48         return x;
  49 }}
  50
  51 static const _cl_MI pow2p1_canonhom (cl_heap_modint_ring* R, const cl_I& x)
  52 {
  53         return _cl_MI(R, pow2p1_reduce_modulo(R,x));
  54 }
  55
  56 static const _cl_MI pow2p1_mul (cl_heap_modint_ring* _R, const _cl_MI& x, const _cl_MI& y)
  57 {
  58         var cl_heap_modint_ring_pow2p1* R = (cl_heap_modint_ring_pow2p1*)_R;
  59         var const uintL m1 = R->m1;
  60         var cl_I zr = x.rep * y.rep;
  61         // Now 0 <= zr <= 2^(2*m1).
  62         zr = ldb(zr,cl_byte(1,2*m1)) - ldb(zr,cl_byte(m1,m1)) + ldb(zr,cl_byte(m1,0));
  63         // Now -(2^m1-1) <= zr <= 2^m1.
  64         return _cl_MI(R, minusp(zr) ? zr + R->modulus : zr);
  65 }
  66
  67 static const _cl_MI pow2p1_square (cl_heap_modint_ring* _R, const _cl_MI& x)
  68 {
  69         var cl_heap_modint_ring_pow2p1* R = (cl_heap_modint_ring_pow2p1*)_R;
  70         var const uintL m1 = R->m1;
  71         var cl_I zr = square(x.rep);
  72         // Now 0 <= zr <= 2^(2*m1).
  73         zr = ldb(zr,cl_byte(1,2*m1)) - ldb(zr,cl_byte(m1,m1)) + ldb(zr,cl_byte(m1,0));
  74         // Now -(2^m1-1) <= zr <= 2^m1.
  75         return _cl_MI(R, minusp(zr) ? zr + R->modulus : zr);
  76 }
  77
  78 #define pow2p1_addops std_addops
  79 static cl_modint_mulops pow2p1_mulops = {
  80         std_one,
  81         pow2p1_canonhom,
  82         pow2p1_mul,
  83         pow2p1_square,
  84         std_expt_pos,
  85         std_recip,
  86         std_div,
  87         std_expt,
  88         pow2p1_reduce_modulo,
  89         std_retract
  90 };
  91
  92 // Constructor.
  93 inline cl_heap_modint_ring_pow2p1::cl_heap_modint_ring_pow2p1 (const cl_I& m, uintL _m1)
  94         : cl_heap_modint_ring (m, &std_setops, &pow2p1_addops, &pow2p1_mulops), m1 (_m1) {}
  95
  96 }  // namespace cln