num_modular/
bigint.rs

1use crate::{ModularAbs, ModularCoreOps, ModularPow, ModularSymbols, ModularUnaryOps};
2use core::convert::TryInto;
3use num_integer::Integer;
4use num_traits::{One, ToPrimitive, Zero};
5
6// Efficient implementation for bigints can be found in "Handbook of Applied Cryptography"
7// Reference: https://cacr.uwaterloo.ca/hac/about/chap14.pdf
8
9// Forward modular operations to ref by ref
10macro_rules! impl_mod_ops_by_ref {
11    ($T:ty) => {
12        // core ops
13        impl ModularCoreOps<$T, &$T> for &$T {
14            type Output = $T;
15            #[inline]
16            fn addm(self, rhs: $T, m: &$T) -> $T {
17                self.addm(&rhs, &m)
18            }
19            #[inline]
20            fn subm(self, rhs: $T, m: &$T) -> $T {
21                self.subm(&rhs, &m)
22            }
23            #[inline]
24            fn mulm(self, rhs: $T, m: &$T) -> $T {
25                self.mulm(&rhs, &m)
26            }
27        }
28        impl ModularCoreOps<&$T, &$T> for $T {
29            type Output = $T;
30            #[inline]
31            fn addm(self, rhs: &$T, m: &$T) -> $T {
32                (&self).addm(rhs, &m)
33            }
34            #[inline]
35            fn subm(self, rhs: &$T, m: &$T) -> $T {
36                (&self).subm(rhs, &m)
37            }
38            #[inline]
39            fn mulm(self, rhs: &$T, m: &$T) -> $T {
40                (&self).mulm(rhs, &m)
41            }
42        }
43        impl ModularCoreOps<$T, &$T> for $T {
44            type Output = $T;
45            #[inline]
46            fn addm(self, rhs: $T, m: &$T) -> $T {
47                (&self).addm(&rhs, &m)
48            }
49            #[inline]
50            fn subm(self, rhs: $T, m: &$T) -> $T {
51                (&self).subm(&rhs, &m)
52            }
53            #[inline]
54            fn mulm(self, rhs: $T, m: &$T) -> $T {
55                (&self).mulm(&rhs, &m)
56            }
57        }
58
59        // pow
60        impl ModularPow<$T, &$T> for &$T {
61            type Output = $T;
62            #[inline]
63            fn powm(self, exp: $T, m: &$T) -> $T {
64                self.powm(&exp, &m)
65            }
66        }
67        impl ModularPow<&$T, &$T> for $T {
68            type Output = $T;
69            #[inline]
70            fn powm(self, exp: &$T, m: &$T) -> $T {
71                (&self).powm(exp, &m)
72            }
73        }
74        impl ModularPow<$T, &$T> for $T {
75            type Output = $T;
76            #[inline]
77            fn powm(self, exp: $T, m: &$T) -> $T {
78                (&self).powm(&exp, &m)
79            }
80        }
81
82        // unary ops and symbols
83        impl ModularUnaryOps<&$T> for $T {
84            type Output = $T;
85            #[inline]
86            fn negm(self, m: &$T) -> $T {
87                ModularUnaryOps::<&$T>::negm(&self, m)
88            }
89            #[inline]
90            fn invm(self, m: &$T) -> Option<$T> {
91                ModularUnaryOps::<&$T>::invm(&self, m)
92            }
93            #[inline]
94            fn dblm(self, m: &$T) -> $T {
95                ModularUnaryOps::<&$T>::dblm(&self, m)
96            }
97            #[inline]
98            fn sqm(self, m: &$T) -> $T {
99                ModularUnaryOps::<&$T>::sqm(&self, m)
100            }
101        }
102    };
103}
104
105#[cfg(feature = "num-bigint")]
106mod _num_bigint {
107    use super::*;
108    use num_bigint::{BigInt, BigUint};
109    use num_traits::Signed;
110
111    impl ModularCoreOps<&BigUint, &BigUint> for &BigUint {
112        type Output = BigUint;
113
114        #[inline]
115        fn addm(self, rhs: &BigUint, m: &BigUint) -> BigUint {
116            (self + rhs) % m
117        }
118        fn subm(self, rhs: &BigUint, m: &BigUint) -> BigUint {
119            let (lhs, rhs) = (self % m, rhs % m);
120            if lhs >= rhs {
121                lhs - rhs
122            } else {
123                m - (rhs - lhs)
124            }
125        }
126
127        fn mulm(self, rhs: &BigUint, m: &BigUint) -> BigUint {
128            let a = self % m;
129            let b = rhs % m;
130
131            if let Some(sm) = m.to_usize() {
132                let sself = a.to_usize().unwrap();
133                let srhs = b.to_usize().unwrap();
134                return BigUint::from(sself.mulm(srhs, &sm));
135            }
136
137            (a * b) % m
138        }
139    }
140
141    impl ModularUnaryOps<&BigUint> for &BigUint {
142        type Output = BigUint;
143        #[inline]
144        fn negm(self, m: &BigUint) -> BigUint {
145            let x = self % m;
146            if x.is_zero() {
147                BigUint::zero()
148            } else {
149                m - x
150            }
151        }
152
153        fn invm(self, m: &BigUint) -> Option<Self::Output> {
154            let x = if self >= m { self % m } else { self.clone() };
155
156            let (mut last_r, mut r) = (m.clone(), x);
157            let (mut last_t, mut t) = (BigUint::zero(), BigUint::one());
158
159            while r > BigUint::zero() {
160                let (quo, rem) = last_r.div_rem(&r);
161                last_r = r;
162                r = rem;
163
164                let new_t = last_t.subm(&quo.mulm(&t, m), m);
165                last_t = t;
166                t = new_t;
167            }
168
169            // if r = gcd(self, m) > 1, then inverse doesn't exist
170            if last_r > BigUint::one() {
171                None
172            } else {
173                Some(last_t)
174            }
175        }
176
177        #[inline]
178        fn dblm(self, m: &BigUint) -> BigUint {
179            let x = self % m;
180            let d = x << 1;
181            if &d > m {
182                d - m
183            } else {
184                d
185            }
186        }
187
188        #[inline]
189        fn sqm(self, m: &BigUint) -> BigUint {
190            self.modpow(&BigUint::from(2u8), m)
191        }
192    }
193
194    impl ModularPow<&BigUint, &BigUint> for &BigUint {
195        type Output = BigUint;
196        #[inline]
197        fn powm(self, exp: &BigUint, m: &BigUint) -> BigUint {
198            self.modpow(exp, m)
199        }
200    }
201
202    impl ModularSymbols<&BigUint> for BigUint {
203        #[inline]
204        fn checked_legendre(&self, n: &BigUint) -> Option<i8> {
205            let r = self.powm((n - 1u8) >> 1u8, n);
206            if r.is_zero() {
207                Some(0)
208            } else if r.is_one() {
209                Some(1)
210            } else if &(r + 1u8) == n {
211                Some(-1)
212            } else {
213                None
214            }
215        }
216
217        fn checked_jacobi(&self, n: &BigUint) -> Option<i8> {
218            if n.is_even() {
219                return None;
220            }
221            if self.is_zero() {
222                return Some(if n.is_one() { 1 } else { 0 });
223            }
224            if self.is_one() {
225                return Some(1);
226            }
227
228            let three = BigUint::from(3u8);
229            let five = BigUint::from(5u8);
230            let seven = BigUint::from(7u8);
231
232            let mut a = self % n;
233            let mut n = n.clone();
234            let mut t = 1;
235            while a > BigUint::zero() {
236                while a.is_even() {
237                    a >>= 1;
238                    if &n & &seven == three || &n & &seven == five {
239                        t *= -1;
240                    }
241                }
242                core::mem::swap(&mut a, &mut n);
243                if (&a & &three) == three && (&n & &three) == three {
244                    t *= -1;
245                }
246                a %= &n;
247            }
248            Some(if n.is_one() { t } else { 0 })
249        }
250
251        #[inline]
252        fn kronecker(&self, n: &BigUint) -> i8 {
253            if n.is_zero() {
254                return if self.is_one() { 1 } else { 0 };
255            }
256            if n.is_one() {
257                return 1;
258            }
259            if n == &BigUint::from(2u8) {
260                return if self.is_even() {
261                    0
262                } else {
263                    let seven = BigUint::from(7u8);
264                    if (self & &seven).is_one() || self & &seven == seven {
265                        1
266                    } else {
267                        -1
268                    }
269                };
270            }
271
272            let f = n.trailing_zeros().unwrap_or(0);
273            let n = n >> f;
274            let t1 = self.kronecker(&BigUint::from(2u8));
275            let t2 = self.jacobi(&n);
276            t1.pow(f.try_into().unwrap()) * t2
277        }
278    }
279
280    impl ModularSymbols<&BigInt> for BigInt {
281        #[inline]
282        fn checked_legendre(&self, n: &BigInt) -> Option<i8> {
283            if n < &BigInt::one() {
284                return None;
285            }
286            self.mod_floor(n)
287                .magnitude()
288                .checked_legendre(n.magnitude())
289        }
290
291        fn checked_jacobi(&self, n: &BigInt) -> Option<i8> {
292            if n < &BigInt::one() {
293                return None;
294            }
295            self.mod_floor(n).magnitude().checked_jacobi(n.magnitude())
296        }
297
298        #[inline]
299        fn kronecker(&self, n: &BigInt) -> i8 {
300            if n.is_negative() {
301                if n.magnitude().is_one() {
302                    return if self.is_negative() { -1 } else { 1 };
303                } else {
304                    return self.kronecker(&-BigInt::one()) * self.kronecker(&-n);
305                }
306            }
307
308            // n is positive from now on
309            let n = n.magnitude();
310            if n.is_zero() {
311                return if self.is_one() { 1 } else { 0 };
312            }
313            if n.is_one() {
314                return 1;
315            }
316            if n == &BigUint::from(2u8) {
317                return if self.is_even() {
318                    0
319                } else {
320                    let eight = BigInt::from(8u8);
321                    if (self.mod_floor(&eight)).is_one()
322                        || self.mod_floor(&eight) == BigInt::from(7u8)
323                    {
324                        1
325                    } else {
326                        -1
327                    }
328                };
329            }
330
331            let f = n.trailing_zeros().unwrap_or(0);
332            let n = n >> f;
333            let t1 = self.kronecker(&BigInt::from(2u8));
334            let t2 = self.jacobi(&n.into());
335            t1.pow(f.try_into().unwrap()) * t2
336        }
337    }
338
339    impl_mod_ops_by_ref!(BigUint);
340
341    impl ModularAbs<BigUint> for BigInt {
342        fn absm(self, m: &BigUint) -> BigUint {
343            if self.is_negative() {
344                self.magnitude().negm(m)
345            } else {
346                self.magnitude() % m
347            }
348        }
349    }
350
351    #[cfg(test)]
352    mod tests {
353        use super::*;
354        use rand::random;
355
356        const NRANDOM: u32 = 10; // number of random tests to run
357
358        #[test]
359        fn basic_tests() {
360            for _ in 0..NRANDOM {
361                let a = random::<u128>();
362                let ra = &BigUint::from(a);
363                let b = random::<u128>();
364                let rb = &BigUint::from(b);
365                let m = random::<u128>() | 1;
366                let rm = &BigUint::from(m);
367                assert_eq!(ra.addm(rb, rm), (ra + rb) % rm);
368                assert_eq!(ra.mulm(rb, rm), (ra * rb) % rm);
369
370                let a = random::<u8>();
371                let ra = &BigUint::from(a);
372                let e = random::<u8>();
373                let re = &BigUint::from(e);
374                let m = random::<u128>() | 1;
375                let rm = &BigUint::from(m);
376                assert_eq!(ra.powm(re, rm), ra.pow(e as u32) % rm);
377            }
378        }
379
380        #[test]
381        fn test_against_prim() {
382            for _ in 0..NRANDOM {
383                let a = random::<u128>();
384                let ra = &BigUint::from(a);
385                let b = random::<u128>();
386                let rb = &BigUint::from(b);
387                let m = random::<u128>();
388                let rm = &BigUint::from(m);
389                assert_eq!(ra.addm(rb, rm), a.addm(b, &m).into());
390                assert_eq!(ra.subm(rb, rm), a.subm(b, &m).into());
391                assert_eq!(ra.mulm(rb, rm), a.mulm(b, &m).into());
392                assert_eq!(ra.negm(rm), a.negm(&m).into());
393                assert_eq!(ra.invm(rm), a.invm(&m).map(|v| v.into()));
394                assert_eq!(ra.checked_legendre(rm), a.checked_legendre(&m));
395                assert_eq!(ra.checked_jacobi(rm), a.checked_jacobi(&m));
396                assert_eq!(ra.kronecker(rm), a.kronecker(&m));
397
398                let e = random::<u8>();
399                let re = &BigUint::from(e);
400                assert_eq!(ra.powm(re, rm), a.powm(e as u128, &m).into());
401
402                // signed integers
403                let a = random::<i128>();
404                let ra = &BigInt::from(a);
405                let m = random::<i128>();
406                let rm = &BigInt::from(m);
407                assert_eq!(ra.checked_legendre(rm), a.checked_legendre(&m));
408                assert_eq!(ra.checked_jacobi(rm), a.checked_jacobi(&m));
409                assert_eq!(ra.kronecker(rm), a.kronecker(&m));
410            }
411        }
412    }
413}