gcd/
lib.rs

1#![no_std]
2use core::num::{NonZeroU128, NonZeroU16, NonZeroU32, NonZeroU64, NonZeroU8, NonZeroUsize};
3
4pub trait Gcd {
5    /// Determine [greatest common divisor](https://en.wikipedia.org/wiki/Greatest_common_divisor)
6    /// using [`gcd_binary`].
7    ///
8    /// [`gcd_binary`]: #method.gcd_binary
9    ///
10    /// # Examples
11    ///
12    /// ```
13    /// use gcd::Gcd;
14    ///
15    /// assert_eq!(0, 0u8.gcd(0));
16    /// assert_eq!(10, 10u8.gcd(0));
17    /// assert_eq!(10, 0u8.gcd(10));
18    /// assert_eq!(10, 10u8.gcd(20));
19    /// assert_eq!(44, 2024u32.gcd(748));
20    /// ```
21    fn gcd(self, other: Self) -> Self;
22
23    /// Determine [greatest common divisor](https://en.wikipedia.org/wiki/Greatest_common_divisor)
24    /// using the [Binary GCD algorithm](https://en.wikipedia.org/wiki/Binary_GCD_algorithm).
25    fn gcd_binary(self, other: Self) -> Self;
26
27    /// Determine [greatest common divisor](https://en.wikipedia.org/wiki/Greatest_common_divisor)
28    /// using the [Euclidean algorithm](https://en.wikipedia.org/wiki/Euclidean_algorithm).
29    fn gcd_euclid(self, other: Self) -> Self;
30}
31
32macro_rules! gcd_impl {
33    ($(($T:ty) $binary:ident $euclid:ident),*) => {$(
34        #[doc = concat!("Const binary GCD implementation for `", stringify!($T), "`.")]
35        pub const fn $binary(mut u: $T, mut v: $T) -> $T
36        {
37            if u == 0 { return v; }
38            if v == 0 { return u; }
39
40            let shift = (u | v).trailing_zeros();
41            u >>= shift;
42            v >>= shift;
43            u >>= u.trailing_zeros();
44
45            loop {
46                v >>= v.trailing_zeros();
47
48                #[allow(clippy::manual_swap)]
49                if u > v {
50                    // mem::swap(&mut u, &mut v);
51                    let temp = u;
52                    u = v;
53                    v = temp;
54                }
55
56                v -= u; // here v >= u
57
58                if v == 0 { break; }
59            }
60
61            u << shift
62        }
63
64        #[doc = concat!("Const euclid GCD implementation for `", stringify!($T), "`.")]
65        pub const fn $euclid(a: $T, b: $T) -> $T
66        {
67            // variable names based off euclidean division equation: a = b · q + r
68            let (mut a, mut b) = if a > b {
69                (a, b)
70            } else {
71                (b, a)
72            };
73
74            #[allow(clippy::manual_swap)]
75            while b != 0 {
76                // mem::swap(&mut a, &mut b);
77                let temp = a;
78                a = b;
79                b = temp;
80
81                b %= a;
82            }
83
84            a
85        }
86
87        impl Gcd for $T {
88            #[inline]
89            fn gcd(self, other: $T) -> $T {
90                self.gcd_binary(other)
91            }
92
93            #[inline]
94            fn gcd_binary(self, v: $T) -> $T {
95                $binary(self, v)
96            }
97
98            #[inline]
99            fn gcd_euclid(self, other: $T) -> $T {
100                $euclid(self, other)
101            }
102        }
103    )*};
104}
105
106gcd_impl! {
107    (u8) binary_u8 euclid_u8,
108    (u16) binary_u16 euclid_u16,
109    (u32) binary_u32 euclid_u32,
110    (u64) binary_u64 euclid_u64,
111    (u128) binary_u128 euclid_u128,
112    (usize) binary_usize euclid_usize
113}
114
115macro_rules! gcd_impl_nonzero {
116    ($(($T:ty) $binary_nonzero:ident/$binary:ident $euclid_nonzero:ident/$euclid:ident),*) => {$(
117        #[doc = concat!("Const binary GCD implementation for `", stringify!($T), "`.")]
118        pub const fn $binary_nonzero(u: $T, v: $T) -> $T
119        {
120            match <$T>::new($binary(u.get(), v.get())) {
121                Some(x) => x,
122                None => unreachable!(),
123            }
124        }
125
126        #[doc = concat!("Const euclid GCD implementation for `", stringify!($T), "`.")]
127        pub const fn $euclid_nonzero(a: $T, b: $T) -> $T
128        {
129            match <$T>::new($euclid(a.get(), b.get())) {
130                Some(x) => x,
131                None => unreachable!(),
132            }
133        }
134
135        impl Gcd for $T {
136            #[inline]
137            fn gcd(self, other: $T) -> $T {
138                self.gcd_binary(other)
139            }
140
141            #[inline]
142            fn gcd_binary(self, v: $T) -> $T {
143                $binary_nonzero(self, v)
144            }
145
146            #[inline]
147            fn gcd_euclid(self, other: $T) -> $T {
148                $euclid_nonzero(self, other)
149            }
150        }
151    )*}
152}
153
154gcd_impl_nonzero! {
155    (NonZeroU8) binary_nonzero_u8/binary_u8 euclid_nonzero_u8/euclid_u8,
156    (NonZeroU16) binary_nonzero_u16/binary_u16 euclid_nonzero_u16/euclid_u16,
157    (NonZeroU32) binary_nonzero_u32/binary_u32 euclid_nonzero_u32/euclid_u32,
158    (NonZeroU64) binary_nonzero_u64/binary_u64 euclid_nonzero_u64/euclid_u64,
159    (NonZeroU128) binary_nonzero_u128/binary_u128 euclid_nonzero_u128/euclid_u128,
160    (NonZeroUsize) binary_nonzero_usize/binary_usize euclid_nonzero_usize/euclid_usize
161}
162
163#[cfg(test)]
164mod test {
165    use super::*;
166    use core::fmt::Debug;
167
168    const U8_GCD_A: [u8; 5] = [140, 1, 140, 33, 225];
169    const U8_GCD_B: [u8; 5] = [136, 123, 203, 252, 153];
170    const U8_GCD_R: [u8; 5] = [4, 1, 7, 3, 9];
171
172    const U16_GCD_A: [u16; 5] = [53144, 44062, 65054, 60568, 11932];
173    const U16_GCD_B: [u16; 5] = [41105, 5088, 35332, 19184, 54004];
174    const U16_GCD_R: [u16; 5] = [1, 2, 22, 8, 4];
175
176    const U32_GCD_A: [u32; 5] = [3392079986, 273672341, 1353048788, 1491301950, 3569727686];
177    const U32_GCD_B: [u32; 5] = [2080089626, 3912533700, 1969135932, 1356732645, 58056677];
178    const U32_GCD_R: [u32; 5] = [2, 1, 4, 15, 7];
179
180    const U64_GCD_A: [u64; 5] = [
181        190266297176832000,
182        2040134905096275968,
183        16611311494648745984,
184        14863931409971066880,
185        11777713923171739648,
186    ];
187    const U64_GCD_B: [u64; 5] = [
188        10430732356495263744,
189        5701159354248194048,
190        7514969329383038976,
191        7911906750992527360,
192        1994469765110767616,
193    ];
194    const U64_GCD_R: [u64; 5] = [6144, 2048, 4096, 10240, 14336];
195
196    const U128_GCD_A: [u128; 5] = [
197        183222947567111613556380400704880115712,
198        115621006611964852903362423926779019264,
199        50724538437787115589243518273596686336,
200        18298803717624646317403958239767298048,
201        196929845599653749349770751890136498176,
202    ];
203    const U128_GCD_B: [u128; 5] = [
204        283620717889381409474181015983148236800,
205        152390035351551984363917166384150216704,
206        74996138554240857099554660445327458304,
207        245604784002268488089190010796573196288,
208        194671916188106984823441978656659865600,
209    ];
210    const U128_GCD_R: [u128; 5] = [
211        37778931862957161709568,
212        75557863725914323419136,
213        113336795588871485128704,
214        151115727451828646838272,
215        302231454903657293676544,
216    ];
217
218    const USIZE_GCD_A: [usize; 5] = [335286345, 3125888386, 3550412466, 924335944, 2870209473];
219    const USIZE_GCD_B: [usize; 5] = [1843742025, 2080426243, 16052620, 1587387560, 24708111];
220    const USIZE_GCD_R: [usize; 5] = [15, 1, 2, 8, 3];
221
222    #[test]
223    fn test_gcd_basic() {
224        // some base cases
225        assert_eq!(0, 0u8.gcd(0));
226        assert_eq!(10, 10u8.gcd(0));
227        assert_eq!(10, 0u8.gcd(10));
228    }
229
230    fn verify_gcd<T>(a: T, b: T, r: T)
231    where
232        T: Gcd + Copy + PartialEq + Debug,
233    {
234        let gcd = a.gcd(b);
235        let egcd = a.gcd_euclid(b);
236        let bgcd = a.gcd_binary(b);
237        assert_eq!(r, gcd, "{:?}.gcd({:?})", a, b);
238        assert_eq!(r, egcd, "{:?}.gcd_euclid({:?})", a, b);
239        assert_eq!(r, bgcd, "{:?}.gcd_binary({:?})", a, b);
240    }
241
242    fn test_gcd_ty<T, NZ, const N: usize>(
243        new: impl Fn(T) -> Option<NZ>,
244        zero: T,
245        a: &[T; N],
246        b: &[T; N],
247        r: &[T; N],
248    ) where
249        T: Gcd + Copy + PartialEq + Debug,
250        NZ: Gcd + Copy + PartialEq + Debug,
251    {
252        for ind in 0..N {
253            let a = new(a[ind]).unwrap();
254            let b = new(b[ind]).unwrap();
255            let r = new(r[ind]).unwrap();
256            verify_gcd(a, b, r);
257        }
258
259        let a = a[0];
260        verify_gcd(zero, a, a);
261        verify_gcd(a, zero, a);
262    }
263
264    #[test]
265    fn test_gcd() {
266        test_gcd_ty(NonZeroU8::new, 0, &U8_GCD_A, &U8_GCD_B, &U8_GCD_R);
267        test_gcd_ty(NonZeroU16::new, 0, &U16_GCD_A, &U16_GCD_B, &U16_GCD_R);
268        test_gcd_ty(NonZeroU32::new, 0, &U32_GCD_A, &U32_GCD_B, &U32_GCD_R);
269        test_gcd_ty(NonZeroU64::new, 0, &U64_GCD_A, &U64_GCD_B, &U64_GCD_R);
270        test_gcd_ty(NonZeroU128::new, 0, &U128_GCD_A, &U128_GCD_B, &U128_GCD_R);
271        test_gcd_ty(
272            NonZeroUsize::new,
273            0,
274            &USIZE_GCD_A,
275            &USIZE_GCD_B,
276            &USIZE_GCD_R,
277        );
278    }
279
280    const U32_GCD_R_0: u32 = binary_u32(U32_GCD_A[0], U32_GCD_B[0]);
281    const U32_GCD_R_1: u32 = euclid_u32(U32_GCD_A[1], U32_GCD_B[1]);
282    const U32_GCD_R_2: u32 = binary_u32(U32_GCD_A[2], U32_GCD_B[2]);
283    const U32_GCD_R_3: u32 = euclid_u32(U32_GCD_A[3], U32_GCD_B[3]);
284    const U32_GCD_R_4: u32 = binary_u32(U32_GCD_A[4], U32_GCD_B[4]);
285
286    #[test]
287    fn test_const_gcd() {
288        assert_eq!(U32_GCD_R[0], U32_GCD_R_0);
289        assert_eq!(U32_GCD_R[1], U32_GCD_R_1);
290        assert_eq!(U32_GCD_R[2], U32_GCD_R_2);
291        assert_eq!(U32_GCD_R[3], U32_GCD_R_3);
292        assert_eq!(U32_GCD_R[4], U32_GCD_R_4);
293    }
294}