arnaucube
arnaucube/poulpyPublic Self-hosted repository
poulpy/src/modulus/impl_u64/montgomery.rs
src/modulus/impl_u64/montgomery.rsRaw
  1
  2use crate::modulus::ReduceOnce;
  3use crate::modulus::montgomery::{MontgomeryPrecomp, Montgomery};
  4use crate::modulus::barrett::{BarrettPrecomp};
  5extern crate test;
  6
  7/// MontgomeryPrecomp is a set of methods implemented for MontgomeryPrecomp<u64>
  8/// enabling Montgomery arithmetic over u64 values.
  9impl MontgomeryPrecomp<u64>{
 10
 11    /// Returns an new instance of MontgomeryPrecomp<u64>.
 12    /// This method will fail if gcd(q, 2^64) != 1.
 13    #[inline(always)]
 14    pub fn new(q: u64) -> MontgomeryPrecomp<u64>{
 15        assert!(q & 1 != 0, "Invalid argument: gcd(q={}, radix=2^64) != 1", q);
 16        let mut q_inv: u64 = 1;
 17        let mut q_pow = q;
 18        for _i in 0..63{
 19            q_inv = q_inv.wrapping_mul(q_pow);
 20            q_pow = q_pow.wrapping_mul(q_pow);
 21        }
 22        let mut precomp = Self{ 
 23            q: q, 
 24            barrett: BarrettPrecomp::new(q), 
 25            q_inv: q_inv,
 26            one: Montgomery(0),
 27            minus_one: Montgomery(0),
 28        };
 29
 30        precomp.one = precomp.prepare(1);
 31        precomp.minus_one = Montgomery(q-precomp.one.value());
 32
 33        precomp
 34    }
 35
 36    /// Returns 2^64 mod q as a Montgomery<u64>.
 37    #[inline(always)]
 38    pub fn one(&self) -> Montgomery<u64>{
 39        self.one
 40    }
 41
 42    /// Returns (q-1) * 2^64 mod q as a Montgomery<u64>.
 43    #[inline(always)]
 44    pub fn minus_one(&self) -> Montgomery<u64>{
 45        self.minus_one
 46    }
 47
 48    /// Returns lhs * 2^64 mod q as a Montgomery<u64>.
 49    #[inline(always)]
 50    pub fn prepare(&self, lhs: u64) -> Montgomery<u64>{
 51        let mut rhs = Montgomery(0);
 52        self.prepare_assign(lhs, &mut rhs);
 53        rhs
 54    }
 55
 56    /// Assigns lhs * 2^64 mod q to rhs.
 57    #[inline(always)]
 58    pub fn prepare_assign(&self, lhs: u64, rhs: &mut Montgomery<u64>){
 59        self.prepare_lazy_assign(lhs, rhs);
 60        rhs.value_mut().reduce_once_assign(self.q);
 61    }
 62
 63    /// Returns lhs * 2^64 mod q in range [0, 2q-1] as a Montgomery<u64>.
 64    #[inline(always)]
 65    pub fn prepare_lazy(&self, lhs: u64) -> Montgomery<u64>{
 66        let mut rhs = Montgomery(0);
 67        self.prepare_lazy_assign(lhs, &mut rhs);
 68        rhs
 69    }
 70
 71    /// Assigns lhs * 2^64 mod q in range [0, 2q-1] to rhs.
 72    #[inline(always)]
 73    pub fn prepare_lazy_assign(&self, lhs: u64, rhs: &mut Montgomery<u64>){
 74        let (_, mhi) = lhs.widening_mul(*self.barrett.value_lo());
 75        *rhs = Montgomery((lhs.wrapping_mul(*self.barrett.value_hi()).wrapping_add(mhi)).wrapping_mul(self.q).wrapping_neg());
 76    }
 77
 78    /// Returns lhs * (2^64)^-1 mod q as a u64.
 79    #[inline(always)]
 80    pub fn unprepare(&self, lhs: Montgomery<u64>) -> u64{
 81        let mut rhs = 0u64;
 82        self.unprepare_assign(lhs, &mut rhs);
 83        rhs
 84    }
 85
 86    /// Assigns lhs * (2^64)^-1 mod q to rhs.
 87    #[inline(always)]
 88    pub fn unprepare_assign(&self, lhs: Montgomery<u64>, rhs: &mut u64){
 89        self.unprepare_lazy_assign(lhs, rhs);
 90        rhs.reduce_once_assign(self.q);
 91        
 92    }
 93
 94    /// Returns lhs * (2^64)^-1 mod q in range [0, 2q-1].
 95    #[inline(always)]
 96    pub fn unprepare_lazy(&self, lhs: Montgomery<u64>) -> u64{
 97        let mut rhs = 0u64;
 98        self.unprepare_lazy_assign(lhs, &mut rhs);
 99        rhs
100    }
101
102    /// Assigns lhs * (2^64)^-1 mod q in range [0, 2q-1] to rhs.
103    #[inline(always)]
104    pub fn unprepare_lazy_assign(&self, lhs: Montgomery<u64>, rhs: &mut u64){
105        let (_, r) = self.q.widening_mul(lhs.value().wrapping_mul(self.q_inv));
106        *rhs = self.q - r
107    }
108
109    /// Returns lhs * rhs * (2^{64})^-1 mod q.
110    #[inline(always)]
111    pub fn mul_external(&self, lhs: Montgomery<u64>, rhs: u64) -> u64{
112        let mut r = self.mul_external_lazy(lhs, rhs);
113        r.reduce_once_assign(self.q);
114        r
115    }
116
117    /// Assigns lhs * rhs * (2^{64})^-1 mod q to rhs.
118    #[inline(always)]
119    pub fn mul_external_assign(&self, lhs: Montgomery<u64>, rhs: &mut u64){
120        self.mul_external_lazy_assign(lhs, rhs);
121        rhs.reduce_once_assign(self.q);
122    }
123
124    /// Returns lhs * rhs * (2^{64})^-1 mod q in range [0, 2q-1].
125    #[inline(always)]
126    pub fn mul_external_lazy(&self, lhs: Montgomery<u64>, rhs: u64) -> u64{
127        let mut result: u64 = rhs;
128        self.mul_external_lazy_assign(lhs, &mut result);
129        result
130    }
131
132    /// Assigns lhs * rhs * (2^{64})^-1 mod q in range [0, 2q-1] to rhs.
133    #[inline(always)]
134    pub fn mul_external_lazy_assign(&self, lhs: Montgomery<u64>, rhs: &mut u64){
135        let (mlo, mhi) = lhs.value().widening_mul(*rhs);
136        let (_, hhi) = self.q.widening_mul(mlo.wrapping_mul(self.q_inv));
137        *rhs = mhi.wrapping_add(self.q - hhi)
138    }
139
140    /// Returns lhs * rhs * (2^{64})^-1 mod q in range [0, 2q-1].
141    #[inline(always)]
142    pub fn mul_internal(&self, lhs: Montgomery<u64>, rhs: Montgomery<u64>) -> Montgomery<u64>{
143        Montgomery(self.mul_external(lhs, *rhs.value()))
144    }
145
146    /// Assigns lhs * rhs * (2^{64})^-1 mod q to rhs.
147    #[inline(always)]
148    pub fn mul_internal_assign(&self, lhs: Montgomery<u64>, rhs: &mut Montgomery<u64>){
149        self.mul_external_assign(lhs, rhs.value_mut());
150    }
151
152    /// Returns lhs * rhs * (2^{64})^-1 mod q in range [0, 2q-1].
153    #[inline(always)]
154    pub fn mul_internal_lazy(&self, lhs: Montgomery<u64>, rhs: Montgomery<u64>) -> Montgomery<u64>{
155        Montgomery(self.mul_external_lazy(lhs, *rhs.value()))
156    }
157
158    /// Assigns lhs * rhs * (2^{64})^-1 mod q in range [0, 2q-1] to rhs.
159    #[inline(always)]
160    pub fn mul_internal_lazy_assign(&self, lhs: Montgomery<u64>, rhs: &mut Montgomery<u64>){
161        self.mul_external_lazy_assign(lhs, rhs.value_mut());
162    }
163
164    #[inline(always)]
165    pub fn add_internal(&self, lhs: Montgomery<u64>, rhs: Montgomery<u64>) -> Montgomery<u64>{
166        Montgomery(self.barrett.reduce(rhs.value() + lhs.value()))
167    }
168
169    /// Assigns lhs + rhs to rhs.
170    #[inline(always)]
171    pub fn add_internal_lazy_assign(&self, lhs: Montgomery<u64>, rhs: &mut Montgomery<u64>){
172        *rhs.value_mut() += lhs.value()
173    }
174
175    /// Assigns lhs + rhs - q if (lhs + rhs) >= q to rhs.
176    #[inline(always)]
177    pub fn add_internal_reduce_once_assign(&self, lhs: Montgomery<u64>, rhs: &mut Montgomery<u64>){
178        self.add_internal_lazy_assign(lhs, rhs);
179        rhs.value_mut().reduce_once_assign(self.q);
180    }
181
182    #[inline(always)]
183    pub fn reduce(&self, lhs: u64) -> u64{
184        self.barrett.reduce(lhs)
185    }
186
187    /// Returns lhs mod q in range [0, 2q-1].
188    #[inline(always)]
189    pub fn reduce_lazy(&self, lhs: u64) -> u64{
190        self.barrett.reduce_lazy(lhs)
191    }
192
193    #[inline(always)]
194    pub fn reduce_assign(&self, lhs: &mut u64){
195        self.barrett.reduce_assign(lhs)
196    }
197
198    /// Returns lhs mod q in range [0, 2q-1].
199    #[inline(always)]
200    pub fn reduce_lazy_assign(&self, lhs: &mut u64){
201        self.barrett.reduce_lazy_assign(lhs)
202    }
203
204    /// Returns (x^exponent) * 2^64 mod q.
205    #[inline(always)]
206    pub fn pow(&self, x: Montgomery<u64>, exponent:u64) -> Montgomery<u64>{
207        let mut y: Montgomery<u64> = self.one();
208        let mut x_mut: Montgomery<u64> = x;
209        let mut i: u64 = exponent;
210        while i > 0{
211            if i & 1 == 1{
212                self.mul_internal_assign(x_mut, &mut y);
213            }
214            self.mul_internal_assign(x_mut, &mut x_mut);
215            i >>= 1;
216        }
217
218        y.value_mut().reduce_once_assign(self.q);
219        y
220    }
221}
222
223/// Returns x^exponent mod q.
224/// This function internally instantiate a new MontgomeryPrecomp<u64>
225/// To be used when called only a few times and if there 
226/// is no Prime instantiated with q.
227fn pow(x:u64, exponent:u64, q:u64) -> u64{
228    let montgomery: MontgomeryPrecomp<u64> = MontgomeryPrecomp::<u64>::new(q);
229    let mut y_mont: Montgomery<u64> = montgomery.one();
230    let mut x_mont: Montgomery<u64> = montgomery.prepare(x);
231    while exponent > 0{
232        if exponent & 1 == 1{
233            montgomery.mul_internal_assign(x_mont, &mut y_mont);
234        }
235
236        montgomery.mul_internal_assign(x_mont, &mut x_mont);
237    }
238    
239    montgomery.unprepare(y_mont)
240}
241
242#[cfg(test)]
243mod tests {
244    use crate::modulus::montgomery;
245    use super::*;
246    use test::Bencher;
247
248    #[test]
249    fn test_mul_external() {
250        let q: u64 = 0x1fffffffffe00001;
251    	let m_precomp = montgomery::MontgomeryPrecomp::new(q);
252        let x: u64 = 0x5f876e514845cc8b;
253        let y: u64 = 0xad726f98f24a761a;
254        let y_mont = m_precomp.prepare(y);
255    	assert!(m_precomp.mul_external(y_mont, x) == (x as u128 * y as u128 % q as u128) as u64);
256    }
257
258    #[bench]
259    fn bench_mul_external(b: &mut Bencher){
260        let q: u64 = 0x1fffffffffe00001;
261    	let m_precomp = montgomery::MontgomeryPrecomp::new(q);
262        let mut x: u64 = 0x5f876e514845cc8b;
263        let y: u64 = 0xad726f98f24a761a;
264        let y_mont = m_precomp.prepare(y);
265        b.iter(|| m_precomp.mul_external_assign(y_mont, &mut x));
266    }
267}