mirror of
https://github.com/RGBCube/uutils-coreutils
synced 2025-07-28 19:47:45 +00:00
factor: Make the implementation of GCD more readable (#1576)
* factor::numeric::gcd: Switch variable names to be more consistent * factor::numeric::gcd: Improve comments * factor::numeric::gcd: Extend loop invariant to v
This commit is contained in:
parent
37f717f5e3
commit
d9be24e354
1 changed files with 24 additions and 20 deletions
|
@ -9,39 +9,43 @@
|
|||
use std::cmp::min;
|
||||
use std::mem::swap;
|
||||
|
||||
pub fn gcd(mut n: u64, mut m: u64) -> u64 {
|
||||
pub fn gcd(mut u: u64, mut v: u64) -> u64 {
|
||||
// Stein's binary GCD algorithm
|
||||
// Base cases: gcd(n, 0) = gcd(0, n) = n
|
||||
if n == 0 {
|
||||
return m;
|
||||
} else if m == 0 {
|
||||
return n;
|
||||
if u == 0 {
|
||||
return v;
|
||||
} else if v == 0 {
|
||||
return u;
|
||||
}
|
||||
|
||||
// Extract common factor-2: gcd(2ⁱ n, 2ⁱ m) = 2ⁱ gcd(n, m)
|
||||
// and reducing until odd gcd(2ⁱ n, m) = gcd(n, m) if m is odd
|
||||
// gcd(2ⁱ u, 2ʲ v) = 2ᵏ gcd(u, v) with u, v odd and k = min(i, j)
|
||||
// 2ᵏ is the greatest power of two that divides both u and v
|
||||
let k = {
|
||||
let k_n = n.trailing_zeros();
|
||||
let k_m = m.trailing_zeros();
|
||||
n >>= k_n;
|
||||
m >>= k_m;
|
||||
min(k_n, k_m)
|
||||
let i = u.trailing_zeros();
|
||||
let j = v.trailing_zeros();
|
||||
u >>= i;
|
||||
v >>= j;
|
||||
min(i, j)
|
||||
};
|
||||
|
||||
loop {
|
||||
// Invariant: n odd
|
||||
debug_assert!(n % 2 == 1, "n = {} is even", n);
|
||||
// Loop invariant: u and v are odd
|
||||
debug_assert!(u % 2 == 1, "u = {} is even", u);
|
||||
debug_assert!(v % 2 == 1, "v = {} is even", v);
|
||||
|
||||
if n > m {
|
||||
swap(&mut n, &mut m);
|
||||
// gcd(u, v) = gcd(|u - v|, min(u, v))
|
||||
if u > v {
|
||||
swap(&mut u, &mut v);
|
||||
}
|
||||
m -= n;
|
||||
v -= u;
|
||||
|
||||
if m == 0 {
|
||||
return n << k;
|
||||
if v == 0 {
|
||||
// Reached the base case; gcd is 2ᵏ u
|
||||
return u << k;
|
||||
}
|
||||
|
||||
m >>= m.trailing_zeros();
|
||||
// gcd(u, 2ʲ v) = gcd(u, v) as u is odd
|
||||
v >>= v.trailing_zeros();
|
||||
}
|
||||
}
|
||||
|
||||
|
|
Loading…
Add table
Add a link
Reference in a new issue