Skip to main content

improv_verified/
gcd.rs

1//! GCD and LCM (the spreadsheet `GCD` / `LCM` functions).
2
3#[allow(unused_imports)]
4use vstd::arithmetic::div_mod::*;
5#[allow(unused_imports)]
6use vstd::prelude::*;
7
8verus! {
9
10/// `d` divides `n`: `n` is an integer multiple of `d`.
11pub open spec fn divides(d: int, n: int) -> bool {
12    exists|k: int| #[trigger] (d * k) == n
13}
14
15/// Euclid's algorithm, as the specification of the greatest common divisor.
16pub open spec fn gcd_spec(a: nat, b: nat) -> nat
17    decreases b,
18{
19    if b == 0 {
20        a
21    } else {
22        gcd_spec(b, a % b)
23    }
24}
25
26/// `gcd_spec(a, b)` divides both `a` and `b`.
27pub proof fn lemma_gcd_divides(a: nat, b: nat)
28    ensures
29        divides(gcd_spec(a, b) as int, a as int),
30        divides(gcd_spec(a, b) as int, b as int),
31    decreases b,
32{
33    if b == 0 {
34        let g = gcd_spec(a, b) as int;
35        assert(g * 1 == a as int);
36        assert(g * 0 == b as int);
37    } else {
38        lemma_gcd_divides(b, a % b);
39        let g = gcd_spec(b, a % b) as int;
40        let k1 = choose|k: int| #[trigger] (g * k) == b as int;
41        let k2 = choose|k: int| #[trigger] (g * k) == (a % b) as int;
42        lemma_fundamental_div_mod(a as int, b as int);
43        let q = a as int / b as int;
44        assert(g * (q * k1 + k2) == a as int) by (nonlinear_arith)
45            requires
46                g * k1 == b as int,
47                g * k2 == (a % b) as int,
48                a as int == b as int * q + (a % b) as int,
49        ;
50    }
51}
52
53/// Every common divisor of `a` and `b` divides `gcd_spec(a, b)`.
54pub proof fn lemma_gcd_greatest(a: nat, b: nat, d: int)
55    requires
56        divides(d, a as int),
57        divides(d, b as int),
58    ensures
59        divides(d, gcd_spec(a, b) as int),
60    decreases b,
61{
62    if b != 0 {
63        let ka = choose|k: int| #[trigger] (d * k) == a as int;
64        let kb = choose|k: int| #[trigger] (d * k) == b as int;
65        lemma_fundamental_div_mod(a as int, b as int);
66        let q = a as int / b as int;
67        assert(d * (ka - kb * q) == (a % b) as int) by (nonlinear_arith)
68            requires
69                d * ka == a as int,
70                d * kb == b as int,
71                a as int == b as int * q + (a % b) as int,
72        ;
73        lemma_gcd_greatest(b, a % b, d);
74    }
75}
76
77/// The greatest common divisor of `a` and `b` (Euclid): it divides both, and
78/// every common divisor divides it. `gcd(a, 0) == a`, so folding from 0 gives
79/// the gcd of a list.
80pub fn gcd(a: u64, b: u64) -> (g: u64)
81    ensures
82        g == gcd_spec(a as nat, b as nat),
83        divides(g as int, a as int),
84        divides(g as int, b as int),
85        forall|d: int| divides(d, a as int) && divides(d, b as int) ==> #[trigger] divides(d, g as int),
86        b == 0 ==> g == a,
87{
88    let mut x = a;
89    let mut y = b;
90    while y != 0
91        invariant
92            gcd_spec(x as nat, y as nat) == gcd_spec(a as nat, b as nat),
93        decreases y,
94    {
95        let t = x % y;
96        x = y;
97        y = t;
98    }
99    proof {
100        lemma_gcd_divides(a as nat, b as nat);
101        assert forall|d: int| divides(d, a as int) && divides(d, b as int) implies #[trigger] divides(
102            d,
103            x as int,
104        ) by {
105            lemma_gcd_greatest(a as nat, b as nat, d);
106        }
107    }
108    x
109}
110
111/// The least common multiple as the engine computes it: `a / gcd(a, b) * b`.
112pub open spec fn lcm_spec(a: nat, b: nat) -> int {
113    (a / gcd_spec(a, b)) as int * b as int
114}
115
116/// `lcm(a, b)` if it is below `cap`, else `None`. Positive `a` and `b`; the
117/// product is formed in `u128`, so the overflow check is exact: `None` exactly
118/// when the true lcm is at least `cap`. A result is a positive common multiple
119/// of `a` and `b`.
120pub fn lcm_below(a: u64, b: u64, cap: u64) -> (r: Option<u64>)
121    requires
122        a > 0,
123        b > 0,
124    ensures
125        r is None <==> lcm_spec(a as nat, b as nat) >= cap,
126        r matches Some(m) ==> m as int == lcm_spec(a as nat, b as nat) && m < cap && m > 0
127            && divides(a as int, m as int) && divides(b as int, m as int),
128{
129    let g = gcd(a, b);
130    proof {
131        // g divides a > 0, so g > 0 and a == g * (a / g).
132        let k = choose|k: int| #[trigger] (g as int * k) == a as int;
133        assert(g > 0) by (nonlinear_arith)
134            requires
135                g as int * k == a as int,
136                a > 0,
137                g >= 0,
138        ;
139        lemma_div_multiples_vanish(k, g as int);
140        assert(a as int / g as int == k);
141    }
142    let q = a / g;
143    proof {
144        assert(q as int * b as int <= u64::MAX as int * u64::MAX as int) by (nonlinear_arith)
145            requires
146                q <= u64::MAX,
147                b <= u64::MAX,
148        ;
149    }
150    let m = (q as u128) * (b as u128);
151    if m < cap as u128 {
152        proof {
153            let k = choose|k: int| #[trigger] (g as int * k) == a as int;
154            let kb = choose|k: int| #[trigger] (g as int * k) == b as int;
155            assert(a as int * kb == m as int) by (nonlinear_arith)
156                requires
157                    g as int * k == a as int,
158                    g as int * kb == b as int,
159                    q as int == k,
160                    m as int == q as int * b as int,
161            ;
162            assert(b as int * q as int == m as int) by (nonlinear_arith)
163                requires
164                    m as int == q as int * b as int,
165            ;
166            assert(m > 0) by (nonlinear_arith)
167                requires
168                    q as int == k,
169                    g as int * k == a as int,
170                    a > 0,
171                    g > 0,
172                    b > 0,
173                    m as int == q as int * b as int,
174            ;
175        }
176        Some(m as u64)
177    } else {
178        None
179    }
180}
181
182} // verus!
183
184#[cfg(test)]
185#[cfg_attr(coverage_nightly, coverage(off))]
186mod tests {
187    use super::*;
188
189    #[test]
190    fn gcd_matches_euclid() {
191        assert_eq!(gcd(0, 0), 0);
192        assert_eq!(gcd(12, 0), 12);
193        assert_eq!(gcd(0, 12), 12);
194        assert_eq!(gcd(12, 18), 6);
195        assert_eq!(gcd(17, 5), 1);
196    }
197
198    #[test]
199    fn lcm_caps_exactly() {
200        assert_eq!(lcm_below(4, 6, 100), Some(12));
201        assert_eq!(lcm_below(4, 6, 12), None);
202        assert_eq!(lcm_below(4, 6, 13), Some(12));
203        assert_eq!(lcm_below(u64::MAX, u64::MAX - 1, u64::MAX), None);
204    }
205}