ads_rs/v1/math/
gcd.rs

1use std::mem::{replace, swap};
2
3/// Finds the GCD (Greatest Common Divisor) for an array of elements.
4/// # Examples
5/// ```
6/// use ads_rs::prelude::v1::math::gcd_many;
7///
8/// let res0 = gcd_many(&[42, 8, 144]);
9/// let res1 = gcd_many(&[89, 144, 233, 377, 610]);
10/// let res2 = gcd_many(&[25, 105, 235, 100]);
11///
12/// assert_eq!(2, res0);
13/// assert_eq!(1, res1);
14/// assert_eq!(5, res2);
15/// ```
16/// ## Corner cases
17/// - GCD of an empty array equals 0.
18/// - GCD of a single element array equals that element.
19/// ```
20/// use ads_rs::prelude::v1::math::gcd_many;
21///
22/// let res0 = gcd_many(&[]);
23/// let res1 = gcd_many(&[25]);
24///
25/// assert_eq!(0, res0);
26/// assert_eq!(25, res1);
27/// ```
28/// # Implementation details
29/// - Stein's algorithm is used.
30/// - Time complexity: O(K * N<sup>2</sup>) where:
31///     - N - bits count in the biggest number.
32///     - K - number's count
33pub fn gcd_many(elems: &[u64]) -> u64 {
34    if elems.is_empty() {
35        return 0;
36    }
37
38    if elems.len() == 1 {
39        return elems[0];
40    }
41
42    elems.iter().fold(0, |acc, e| {
43        let (mut lhs, mut rhs) = (acc, *e);
44
45        if lhs == 0 || rhs == 0 {
46            return lhs | rhs;
47        }
48
49        // find common factor of 2
50        let shift = (lhs | rhs).trailing_zeros();
51
52        // divide lhs and rhs by 2 until odd
53        rhs >>= rhs.trailing_zeros();
54        while lhs > 0 {
55            lhs >>= lhs.trailing_zeros();
56
57            if rhs > lhs {
58                swap(&mut lhs, &mut rhs);
59            }
60
61            lhs -= rhs
62        }
63
64        rhs << shift
65    })
66}
67
68/// Finds an extended GCD (Greatest Common Divisor) for a pair of numbers.  
69/// "Extended" means that algorithm will return not only GCD, but two coefficients `x` and `y` such that the equality
70///
71/// x * lhs + y * rhs = gcd(lhs, rhs)
72///
73/// holds.
74/// # Examples
75/// ```
76/// use ads_rs::prelude::v1::math::extended_gcd;
77///
78/// let res0 = extended_gcd(30, 20);
79/// let res1 = extended_gcd(15, 35);
80/// let res2 = extended_gcd(161, 28);
81///
82/// assert_eq!((10, 1, -1), res0);
83/// assert_eq!((5, -2, 1), res1);
84/// assert_eq!((7, -1, 6), res2);
85/// ```
86/// ## Corner case
87/// - Result of `extended_gcd(0, 0)` equals tuple `(0, 1, 0)`.
88/// - Negative numbers is not supported, but implementation allows it.
89/// ```
90/// use ads_rs::prelude::v1::math::extended_gcd;
91///
92/// let res = extended_gcd(0, 0);
93///
94/// assert_eq!((0, 1, 0), res);
95/// ```
96/// # Implementation details
97/// - Euclid's algorithm used, because its extended version is faster than Stein's algorithm
98/// - Time complexity is O(log<sub>2</sub>(min(lhs, rhs)))
99pub fn extended_gcd(lhs: u64, rhs: u64) -> (u64, i64, i64) {
100    let (mut x, mut y) = (1, 0);
101    let (mut x1, mut y1, mut lhs1, mut rhs1) = (0i64, 1i64, lhs, rhs);
102
103    while rhs1 > 0 {
104        let q = lhs1 / rhs1;
105
106        let new_x1 = x - (q as i64) * x1;
107        x = replace(&mut x1, new_x1);
108
109        let new_y1 = y - (q as i64) * y1;
110        y = replace(&mut y1, new_y1);
111
112        let new_rhs1 = lhs1 - q * rhs1;
113        lhs1 = replace(&mut rhs1, new_rhs1);
114    }
115
116    (lhs1, x, y)
117}
118
119/// Finds an GCD (Greatest Common Divisor) for a pair of numbers.
120/// # Examples
121/// ```
122/// use ads_rs::prelude::v1::math::gcd;
123///
124/// let res0 = gcd(42, 144);
125/// let res1 = gcd(377, 610);
126/// let res2 = gcd(105, 25);
127///
128/// assert_eq!(6, res0);
129/// assert_eq!(1, res1);
130/// assert_eq!(5, res2);
131/// ```
132/// ## Corner case
133/// GCD of both zero numbers equals 0.
134/// ```
135/// use ads_rs::prelude::v1::math::gcd;
136///
137/// let res = gcd(0, 0);
138///
139/// assert_eq!(0, res);
140/// ```
141/// # Implementation details
142/// - Stein's algorithm used (from [`gcd_many`]).
143/// - Time complexity: O(N<sup>2</sup>) where N - number of bits in the biggest number.
144#[inline]
145pub fn gcd(lhs: u64, rhs: u64) -> u64 {
146    gcd_many(&[lhs, rhs])
147}
148
149#[cfg(test)]
150mod tests {
151    use super::*;
152
153    #[test]
154    fn extended_gcd_works() {
155        // arrange
156        let test_suits = [
157            // Zero rhs
158            (123, 0, (123, 1, 0)),
159            // Zero lhs
160            (0, 123, (123, 0, 1)),
161            // Regular case
162            (2048, 48, (16, -1, 43)),
163            // Relative prime
164            (2052, 617, (1, 132, -439)),
165            // Zero lhs and rhs
166            (0, 0, (0, 1, 0)),
167        ];
168
169        // act
170        let result: Vec<(u64, i64, i64)> =
171            test_suits.iter().map(|t| extended_gcd(t.0, t.1)).collect();
172
173        // assert
174        for i in 0..test_suits.len() {
175            assert_eq!(test_suits[i].2, result[i]);
176        }
177    }
178    #[test]
179    fn gcd_many_works() {
180        // arrange
181        let test_suits = [
182            // Empty array
183            (vec![], 0),
184            // Single element
185            (vec![223], 223),
186            // Relative prime numbers
187            (vec![1, 2, 3, 4, 5], 1),
188            // Regular case
189            (vec![8, 24, 156, 36], 4),
190            // All zeros
191            (vec![0, 0, 0, 0], 0),
192        ];
193
194        // act
195        let result: Vec<u64> = test_suits.iter().map(|t| gcd_many(&t.0)).collect();
196
197        // assert
198        for i in 0..test_suits.len() {
199            assert_eq!(test_suits[i].1, result[i]);
200        }
201    }
202}