quadratic/
lib.rs

1//! Helper functions related to [quadratic reciprocity](https://en.wikipedia.org/wiki/Quadratic_reciprocity).
2//!
3//! This library calculates the [Legendre](http://mathworld.wolfram.com/LegendreSymbol.html)
4//! and [Jacobi symbol](http://mathworld.wolfram.com/JacobiSymbol.html) using the
5//! [standard algorithm](https://en.wikipedia.org/wiki/Jacobi_symbol#Calculating_the_Jacobi_symbol).
6//!
7//! Since the Jacobi symbol is a generalization of the Legendre symbol,
8//! only a `jacobi` function is provided.
9//!
10
11extern crate num;
12use num::Integer;
13use num::traits::ToPrimitive;
14
15
16/// Supplementary law
17fn two_over(n: usize) -> (i8) {
18    if n % 8 == 1 || n % 8 == 7 { 1 } else { -1 }
19}
20
21/// Legendre's version of quadratic reciprocity law
22/// Returns the sign change needed to keep track of if flipping
23fn reciprocity(num: usize, den: usize) -> (i8) {
24    if num % 4 == 3 && den % 4 == 3 { -1 } else { 1 }
25}
26
27/// Returns the value of the Jacobi symbol for `(a \ n)`
28/// Where `a` is an integer and `n` is an odd positive integer
29///
30/// If the intergers passed in do not match these assumptions,
31/// the return value will be `None`. Otherwise the algorithm can be assumed to
32/// always succeed.
33///
34/// # Example
35///
36/// ```rust
37/// # extern crate quadratic;
38/// // compute the Jacobi symbol of (2 \ 5)
39/// let symb = quadratic::jacobi(2, 5);
40/// assert_eq!(Some(-1), symb);
41/// ```
42
43pub fn jacobi(a: isize, n: isize) -> Option<i8> {
44    // jacobi symbol is only defined for odd positive moduli
45    if n.is_even() || n <= 0 {
46        return None;
47    }
48
49    // Raise a mod n, then start the unsigned algorithm
50    let mut acc = 1;
51    let mut num = a.mod_floor(&n).to_usize().unwrap();
52    let mut den = n as usize;
53    loop {
54        // reduce numerator
55        num = num % den;
56        if num == 0 {
57            return Some(0);
58        }
59
60        // extract factors of two from numerator
61        while num.is_even() {
62            acc *= two_over(den);
63            num /= 2;
64        }
65        // if numerator is 1 => this sub-symbol is 1
66        if num == 1 {
67            return Some(acc);
68        }
69        // shared factors => one sub-symbol is zero
70        if num.gcd(&den) > 1 {
71            return Some(0);
72        }
73        // num and den are now odd co-prime, use reciprocity law:
74        acc *= reciprocity(num, den);
75        let tmp = num;
76        num = den;
77        den = tmp;
78    }
79}
80
81#[cfg(test)]
82mod tests {
83    use jacobi;
84
85    // legendre tests
86    #[test]
87    fn minus_one_over_p() {
88        // 1 mod 4 => 1
89        assert_eq!(Some(1), jacobi(-1, 5));
90        assert_eq!(Some(1), jacobi(-1, 13));
91        // 3 mod 4 => -1
92        assert_eq!(Some(-1), jacobi(-1, 3));
93        assert_eq!(Some(-1), jacobi(-1, 7));
94    }
95    #[test]
96    fn two_over_p() {
97        assert_eq!(Some(-1), jacobi(2, 3));
98        assert_eq!(Some(-1), jacobi(2, 5));
99        assert_eq!(Some(1), jacobi(2, 7));
100        assert_eq!(Some(1), jacobi(2, 17)); // 17 = 1 mod 8
101    }
102    #[test]
103    fn three_over_p() {
104        assert_eq!(Some(0), jacobi(3, 3));
105        assert_eq!(Some(-1), jacobi(3, 5));
106        assert_eq!(Some(-1), jacobi(3, 7));
107    }
108    #[test]
109    fn periodicity() {
110        assert_eq!(jacobi(3,5), jacobi(-2,5));
111        assert_eq!(jacobi(-1,5), jacobi(4,5));
112        assert_eq!(jacobi(11,7), jacobi(4,7));
113        assert_eq!(jacobi(-3,7), jacobi(4,7));
114        assert_eq!(jacobi(10,7), jacobi(3,7));
115    }
116
117    #[test]
118    fn jacobi_simple() { // 45 = 3*3*5
119        assert_eq!(Some(-1), jacobi(2, 45)); // -1 * -1 * -1
120        assert_eq!(Some(0), jacobi(3, 45)); // 0 * 0 * -1
121        assert_eq!(Some(-1), jacobi(7, 45)); // -1 * -1 * -1
122        assert_eq!(Some(1), jacobi(2, 15)); // (2\5) * (2\3)
123        assert_eq!(Some(-1), jacobi(1001, 9907)); // wikepedia example
124    }
125
126    #[test]
127    fn even_moduli_fails() {
128        assert_eq!(None, jacobi(2, 4));
129    }
130    #[test]
131    fn negative_moduli_fails() {
132        assert_eq!(None, jacobi(2, -3));
133    }
134}
quadratic/lib.rs

quadratic/
lib.rs
