num_modular/
lib.rs

1//! This crate provides efficient Modular arithmetic operations for various integer types,
2//! including primitive integers and `num-bigint`. The latter option is enabled optionally.
3//!
4//! To achieve fast modular arithmetics, convert integers to any [ModularInteger] implementation
5//! use static `new()` or associated [ModularInteger::convert()] functions. Some builtin implementations
6//! of [ModularInteger] includes [MontgomeryInt], [FixedMersenneInt] and [FixedSolinasInt].
7//!
8//! Example code:
9//! ```rust
10//! use num_modular::{ModularCoreOps, ModularInteger, MontgomeryInt};
11//!
12//! // directly using methods in ModularCoreOps
13//! let (x, y, m) = (12u8, 13u8, 5u8);
14//! assert_eq!(x.mulm(y, &m), x * y % m);
15//!
16//! // convert integers into ModularInteger
17//! let mx = MontgomeryInt::new(x, &m);
18//! let my = mx.convert(y); // faster than static MontgomeryInt::new(y, m)
19//! assert_eq!((mx * my).residue(), x * y % m);
20//! ```
21//!
22//! # Comparison of fast division / modular arithmetics
23//! Several fast division / modulo tricks are provided in these crate, the difference of them are listed below:
24//! - [PreModInv]: pre-compute modular inverse of the divisor, only applicable to exact division
25//! - Barrett (to be implemented): pre-compute (rational approximation of) the reciprocal of the divisor,
26//!   applicable to fast division and modulo
27//! - [Montgomery]: Convert the dividend into a special form by shifting and pre-compute a modular inverse,
28//!   only applicable to fast modulo, but faster than Barrett reduction
29//! - [FixedMersenne]: Specialization of modulo in form `2^P-K` under 2^127.
30//! - [FixedTrinomialSolinas]: Specialization of modulo in form `2^P1-2^P2+K` under 2^127.
31//!
32
33// XXX: Other fast modular arithmetic tricks
34// REF: https://github.com/lemire/fastmod & https://arxiv.org/pdf/1902.01961.pdf
35// REF: https://eprint.iacr.org/2014/040.pdf
36// REF: https://github.com/ridiculousfish/libdivide/
37// REF: Faster Interleaved Modular Multiplication Based on Barrett and Montgomery Reduction Methods (work for modulus in certain form)
38
39#![no_std]
40#[cfg(any(feature = "std", test))]
41extern crate std;
42
43use core::ops::{Add, Mul, Neg, Sub};
44
45/// Core modular arithmetic operations.
46///
47/// Note that all functions will panic if the modulus is zero.
48pub trait ModularCoreOps<Rhs = Self, Modulus = Self> {
49    type Output;
50
51    /// Return (self + rhs) % m
52    fn addm(self, rhs: Rhs, m: Modulus) -> Self::Output;
53
54    /// Return (self - rhs) % m
55    fn subm(self, rhs: Rhs, m: Modulus) -> Self::Output;
56
57    /// Return (self * rhs) % m
58    fn mulm(self, rhs: Rhs, m: Modulus) -> Self::Output;
59}
60
61/// Core unary modular arithmetics
62///
63/// Note that all functions will panic if the modulus is zero.
64pub trait ModularUnaryOps<Modulus = Self> {
65    type Output;
66
67    /// Return (-self) % m and make sure the result is normalized in range [0,m)
68    fn negm(self, m: Modulus) -> Self::Output;
69
70    /// Calculate modular inverse (x such that self*x = 1 mod m).
71    ///
72    /// This operation is only available for integer that is coprime to `m`. If not,
73    /// the result will be [None].
74    fn invm(self, m: Modulus) -> Option<Self::Output>;
75
76    /// Calculate modular double ( x+x mod m)
77    fn dblm(self, m: Modulus) -> Self::Output;
78
79    /// Calculate modular square ( x*x mod m )
80    fn sqm(self, m: Modulus) -> Self::Output;
81
82    // TODO: Modular sqrt aka Quadratic residue, follow the behavior of FLINT `n_sqrtmod`
83    // fn sqrtm(self, m: Modulus) -> Option<Self::Output>;
84    // REF: https://stackoverflow.com/questions/6752374/cube-root-modulo-p-how-do-i-do-this
85}
86
87/// Modular power functions
88pub trait ModularPow<Exp = Self, Modulus = Self> {
89    type Output;
90
91    /// Return (self ^ exp) % m
92    fn powm(self, exp: Exp, m: Modulus) -> Self::Output;
93}
94
95/// Math symbols related to modular arithmetics
96pub trait ModularSymbols<Modulus = Self> {
97    /// Calculate Legendre Symbol (a|n), where a is `self`.
98    ///
99    /// Note that this function doesn't perform a full primality check, since
100    /// is costly. So if n is not a prime, the result can be not reasonable.
101    ///
102    /// # Panics
103    /// Only if n is not prime
104    #[inline]
105    fn legendre(&self, n: Modulus) -> i8 {
106        self.checked_legendre(n).expect("n shoud be a prime")
107    }
108
109    /// Calculate Legendre Symbol (a|n), where a is `self`. Returns [None] only if n is
110    /// not a prime.
111    ///
112    /// Note that this function doesn't perform a full primality check, since
113    /// is costly. So if n is not a prime, the result can be not reasonable.
114    ///
115    /// # Panics
116    /// Only if n is not prime
117    fn checked_legendre(&self, n: Modulus) -> Option<i8>;
118
119    /// Calculate Jacobi Symbol (a|n), where a is `self`
120    ///
121    /// # Panics
122    /// if n is negative or even
123    #[inline]
124    fn jacobi(&self, n: Modulus) -> i8 {
125        self.checked_jacobi(n)
126            .expect("the Jacobi symbol is only defined for non-negative odd integers")
127    }
128
129    /// Calculate Jacobi Symbol (a|n), where a is `self`. Returns [None] if n is negative or even.
130    fn checked_jacobi(&self, n: Modulus) -> Option<i8>;
131
132    /// Calculate Kronecker Symbol (a|n), where a is `self`
133    fn kronecker(&self, n: Modulus) -> i8;
134}
135
136// TODO: Discrete log aka index, follow the behavior of FLINT `n_discrete_log_bsgs`
137// REF: https://github.com/vks/discrete-log
138// fn logm(self, base: Modulus, m: Modulus);
139
140/// Collection of common modular arithmetic operations
141pub trait ModularOps<Rhs = Self, Modulus = Self, Output = Self>:
142    ModularCoreOps<Rhs, Modulus, Output = Output>
143    + ModularUnaryOps<Modulus, Output = Output>
144    + ModularPow<Rhs, Modulus, Output = Output>
145    + ModularSymbols<Modulus>
146{
147}
148impl<T, Rhs, Modulus> ModularOps<Rhs, Modulus> for T where
149    T: ModularCoreOps<Rhs, Modulus, Output = T>
150        + ModularUnaryOps<Modulus, Output = T>
151        + ModularPow<Rhs, Modulus, Output = T>
152        + ModularSymbols<Modulus>
153{
154}
155
156/// Collection of operations similar to [ModularOps], but takes operands with references
157pub trait ModularRefOps: for<'r> ModularOps<&'r Self, &'r Self> + Sized {}
158impl<T> ModularRefOps for T where T: for<'r> ModularOps<&'r T, &'r T> {}
159
160/// Provides a utility function to convert signed integers into unsigned modular form
161pub trait ModularAbs<Modulus> {
162    /// Return self % m, but accepting signed integers
163    fn absm(self, m: &Modulus) -> Modulus;
164}
165
166/// Represents an number defined in a modulo ring ℤ/nℤ
167///
168/// The operators should panic if the modulus of two number
169/// are not the same.
170pub trait ModularInteger:
171    Sized
172    + PartialEq
173    + Add<Self, Output = Self>
174    + Sub<Self, Output = Self>
175    + Neg<Output = Self>
176    + Mul<Self, Output = Self>
177{
178    /// The underlying representation type of the integer
179    type Base;
180
181    /// Return the modulus of the ring
182    fn modulus(&self) -> Self::Base;
183
184    /// Return the normalized residue of this integer in the ring
185    fn residue(&self) -> Self::Base;
186
187    /// Check if the integer is zero
188    fn is_zero(&self) -> bool;
189
190    /// Convert an normal integer into the same ring.
191    ///
192    /// This method should be perferred over the static
193    /// constructor to prevent unnecessary overhead of pre-computation.
194    fn convert(&self, n: Self::Base) -> Self;
195
196    /// Calculate the value of self + self
197    fn double(self) -> Self;
198
199    /// Calculate the value of self * self
200    fn square(self) -> Self;
201}
202
203// XXX: implement ModularInteger for ff::PrimeField?
204// TODO: implement invm_range (Modular inverse in certain range) and crt (Chinese Remainder Theorem), REF: bubblemath crate
205
206/// Utility function for exact division, with precomputed helper values
207///
208/// # Available Pre-computation types:
209/// - `()`: No pre-computation, the implementation relies on native integer division
210/// - [PreModInv]: With Pre-computed modular inverse
211pub trait DivExact<Rhs, Precompute>: Sized {
212    type Output;
213
214    /// Check if d divides self with the help of the precomputation. If d divides self,
215    /// then the quotient is returned.
216    fn div_exact(self, d: Rhs, pre: &Precompute) -> Option<Self::Output>;
217}
218
219/// A modular reducer that can ensure that the operations on integers are all performed
220/// in a modular ring.
221///
222/// Essential information for performing the modulo operation will be stored in the reducer.
223pub trait Reducer<T> {
224    /// Create a reducer for a modulus m
225    fn new(m: &T) -> Self;
226
227    /// Transform a normal integer into reduced form
228    fn transform(&self, target: T) -> T;
229
230    /// Check whether target is a valid reduced form
231    fn check(&self, target: &T) -> bool;
232
233    /// Get the modulus in original integer type
234    fn modulus(&self) -> T;
235
236    /// Transform a reduced form back to normal integer
237    fn residue(&self, target: T) -> T;
238
239    /// Test if the residue() == 0
240    fn is_zero(&self, target: &T) -> bool;
241
242    /// Calculate (lhs + rhs) mod m in reduced form
243    fn add(&self, lhs: &T, rhs: &T) -> T;
244
245    #[inline]
246    fn add_in_place(&self, lhs: &mut T, rhs: &T) {
247        *lhs = self.add(lhs, rhs)
248    }
249
250    /// Calculate 2*target mod m
251    fn dbl(&self, target: T) -> T;
252
253    /// Calculate (lhs - rhs) mod m in reduced form
254    fn sub(&self, lhs: &T, rhs: &T) -> T;
255
256    #[inline]
257    fn sub_in_place(&self, lhs: &mut T, rhs: &T) {
258        *lhs = self.sub(lhs, rhs);
259    }
260
261    /// Calculate -monty mod m in reduced form
262    fn neg(&self, target: T) -> T;
263
264    /// Calculate (lhs * rhs) mod m in reduced form
265    fn mul(&self, lhs: &T, rhs: &T) -> T;
266
267    #[inline]
268    fn mul_in_place(&self, lhs: &mut T, rhs: &T) {
269        *lhs = self.mul(lhs, rhs);
270    }
271
272    /// Calculate target^-1 mod m in reduced form,
273    /// it may return None when there is no modular inverse.
274    fn inv(&self, target: T) -> Option<T>;
275
276    /// Calculate target^2 mod m in reduced form
277    fn sqr(&self, target: T) -> T;
278
279    /// Calculate base ^ exp mod m in reduced form
280    fn pow(&self, base: T, exp: &T) -> T;
281}
282
283mod barrett;
284mod double;
285mod mersenne;
286mod monty;
287mod preinv;
288mod prim;
289mod reduced;
290mod solinas;
291mod word;
292
293pub use barrett::{
294    Normalized2by1Divisor, Normalized3by2Divisor, PreMulInv1by1, PreMulInv2by1, PreMulInv3by2,
295};
296pub use double::{imax, udouble, umax};
297pub use mersenne::{FixedMersenne, FixedMersenne32, FixedMersenne64};
298pub use monty::Montgomery;
299pub use preinv::PreModInv;
300pub use reduced::{ReducedInt, Vanilla, VanillaInt};
301pub use solinas::{FixedTrinomialSolinas, FixedTrinomialSolinas32, FixedTrinomialSolinas64};
302
303/// An integer in modulo ring based on [Montgomery form](https://en.wikipedia.org/wiki/Montgomery_modular_multiplication#Montgomery_form)
304pub type MontgomeryInt<T> = ReducedInt<T, Montgomery<T>>;
305
306/// An integer in modulo ring with a fixed (pseudo) Mersenne number as modulus
307pub type FixedMersenneInt<const P: u8, const K: umax> = ReducedInt<umax, FixedMersenne<P, K>>;
308
309/// An integer in modulo ring with a fixed (pseudo) Mersenne number as modulus (32-bit)
310pub type FixedMersenneInt32<const P: u8, const K: u32> = ReducedInt<u32, FixedMersenne32<P, K>>;
311
312/// An integer in modulo ring with a fixed (pseudo) Mersenne number as modulus (64-bit)
313pub type FixedMersenneInt64<const P: u8, const K: u64> = ReducedInt<u64, FixedMersenne64<P, K>>;
314
315/// An integer in modulo ring with a fixed Solinas number as modulus
316pub type FixedSolinasInt<const P1: u8, const P2: u8, const K: imax> =
317    ReducedInt<umax, FixedTrinomialSolinas<P1, P2, K>>;
318
319/// An integer in modulo ring with a fixed Solinas number as modulus (32-bit)
320pub type FixedSolinasInt32<const P1: u8, const P2: u8, const K: i32> =
321    ReducedInt<u32, FixedTrinomialSolinas32<P1, P2, K>>;
322
323/// An integer in modulo ring with a fixed Solinas number as modulus (64-bit)
324pub type FixedSolinasInt64<const P1: u8, const P2: u8, const K: i64> =
325    ReducedInt<u64, FixedTrinomialSolinas64<P1, P2, K>>;
326
327// pub type BarrettInt<T> = ReducedInt<T, BarrettInt<T>>;
328
329#[cfg(feature = "num-bigint")]
330mod bigint;
num_modular/lib.rs

num_modular/
lib.rs
