ospf_rust_math/algebra/ordinary/
gcd.rs

1use std::mem::swap;
2use std::ops::{BitOr, SubAssign};
3
4use paste::paste;
5
6use crate::algebra::concept::{Arithmetic, Bits, SemiArithmetic};
7use crate::algebra::operator::TrailingZeros;
8
9pub fn gcd_stein<T: SemiArithmetic + Bits + TrailingZeros + for<'a> SubAssign<&'a T> + >(
10    mut x: T,
11    mut y: T,
12) -> T
13where
14    for<'a> &'a T: PartialEq + BitOr<Output = T>,
15{
16    debug_assert!(&x >= T::ZERO);
17    debug_assert!(&y >= T::ZERO);
18
19    if &x == T::ZERO {
20        return y;
21    }
22    if &y == T::ZERO {
23        return x;
24    }
25
26    let shift = T::trailing_zeros(&x | &y);
27    x >>= shift;
28    y >>= shift;
29    x >>= T::trailing_zeros(x.clone());
30
31    loop {
32        y >>= T::trailing_zeros(y.clone());
33        if x > y {
34            swap(&mut x, &mut y);
35        }
36        y -= &x;
37        if &y == T::ZERO {
38            break;
39        }
40    }
41    x << shift
42}
43
44pub fn gcd_euclid<T: SemiArithmetic + for<'a> SubAssign<&'a T>>(mut x: T, mut y: T) -> T
45where
46    for<'a> &'a T: PartialEq,
47{
48    debug_assert!(&x >= T::ZERO);
49    debug_assert!(&y >= T::ZERO);
50
51    if &x == T::ZERO {
52        return y;
53    }
54    if &y == T::ZERO {
55        return x;
56    }
57
58    loop {
59        if x > y {
60            swap(&mut x, &mut y);
61        }
62        y -= &x;
63        if &y == T::ZERO {
64            break;
65        }
66    }
67    x
68}
69
70macro_rules! gcd_stein_template {
71    ($($type:ident)*) => ($(
72        paste! {
73            pub fn [<gcd_stein_ $type>](mut x: $type, mut y: $type) -> $type {
74                debug_assert!(&x >= $type::ZERO);
75                debug_assert!(&y >= $type::ZERO);
76
77                if &x == $type::ZERO {
78                    return y;
79                }
80                if &y == $type::ZERO {
81                    return x;
82                }
83
84                let shift = $type::trailing_zeros(&x | &y);
85                x >>= shift;
86                y >>= shift;
87                x >>= $type::trailing_zeros(x.clone());
88
89                loop {
90                    y >>= $type::trailing_zeros(y.clone());
91                    if x > y {
92                        swap(&mut x, &mut y);
93                    }
94                    y -= &x;
95                    if &y == $type::ZERO {
96                        break;
97                    }
98                }
99                x << shift
100            }
101        }
102    )*)
103}
104gcd_stein_template! { u8 u16 u32 u64 u128 usize i8 i16 i32 i64 i128 isize }
105
106macro_rules! gcd_euclid_template {
107    ($($type:ident)*) => ($(
108        paste! {
109            pub fn [<gcd_euclid_ $type>](mut x: $type, mut y: $type) -> $type {
110                debug_assert!(&x >= $type::ZERO);
111                debug_assert!(&y >= $type::ZERO);
112
113                if &x == $type::ZERO {
114                    return y;
115                }
116                if &y == $type::ZERO {
117                    return x;
118                }
119
120                loop {
121                    if x > y {
122                        swap(&mut x, &mut y);
123                    }
124                    y -= &x;
125                    if &y == $type::ZERO {
126                        break;
127                    }
128                }
129                x
130            }
131        }
132    )*)
133}
134gcd_euclid_template! { u8 u16 u32 u64 u128 usize i8 i16 i32 i64 i128 isize f32 f64 }
135
136#[cfg(test)]
137mod tests {
138    use super::*;
139
140    #[test]
141    fn test_stein() {
142        assert_eq!(gcd_stein_i32(4, 6), 2);
143        assert_eq!(gcd_stein_i32(6, 9), 3);
144        assert_eq!(gcd_stein_i32(24, 30), 6);
145    }
146
147    #[test]
148    fn test_euclid() {
149        assert_eq!(gcd_euclid_i32(4, 6), 2);
150        assert_eq!(gcd_euclid_i32(6, 9), 3);
151        assert_eq!(gcd_euclid_i32(24, 30), 6);
152    }
153}