Greatest common divisor, using the Euclidean algorithm

Topics: Arithmetic / Divisibility

Tools: Why3

see also the index (by topic, by tool, by reference, by year)

```(* Greatest common divisor, using the Euclidean algorithm *)

module EuclideanAlgorithm

use mach.int.Int
use number.Gcd

let rec euclid (u v: int) : int
variant  { v }
requires { u >= 0 /\ v >= 0 }
ensures  { result = gcd u v }
=
if v = 0 then
u
else
euclid v (u % v)

end

module EuclideanAlgorithmIterative

use mach.int.Int
use ref.Ref
use number.Gcd

let euclid (u0 v0: int) : int
requires { u0 >= 0 /\ v0 >= 0 }
ensures  { result = gcd u0 v0 }
=
let ref u = u0 in
let ref v = v0 in
while v <> 0 do
invariant { u >= 0 /\ v >= 0 }
invariant { gcd u v = gcd u0 v0 }
variant   { v }
let tmp = v in
v <- u % v;
u <- tmp
done;
u

end

module BinaryGcd

use mach.int.Int
use number.Parity
use number.Gcd

lemma even1: forall n: int. 0 <= n -> even n <-> n = 2 * div n 2
lemma odd1: forall n: int. 0 <= n -> not (even n) <-> n = 2 * div n 2 + 1
lemma div_nonneg: forall n: int. 0 <= n -> 0 <= div n 2

use number.Coprime

lemma gcd_even_even: forall u v: int. 0 <= v -> 0 <= u ->
gcd (2 * u) (2 * v) = 2 * gcd u v
lemma gcd_even_odd: forall u v: int. 0 <= v -> 0 <= u ->
gcd (2 * u) (2 * v + 1) = gcd u (2 * v + 1)
lemma gcd_even_odd2: forall u v: int. 0 <= v -> 0 <= u ->
even u -> odd v -> gcd u v = gcd (div u 2) v
lemma odd_odd_div2: forall u v: int. 0 <= v -> 0 <= u ->
div ((2*u+1) - (2*v+1)) 2 = u - v

let lemma gcd_odd_odd (u v: int)
requires { 0 <= v <= u }
ensures { gcd (2 * u + 1) (2 * v + 1) = gcd (u - v) (2 * v + 1) }
= assert { gcd (2 * u + 1) (2 * v + 1) =
gcd ((2*u+1) - 1*(2*v+1)) (2 * v + 1) }

lemma gcd_odd_odd2: forall u v: int. 0 <= v <= u ->
odd u -> odd v -> gcd u v = gcd (div (u - v) 2) v

let rec binary_gcd (u v: int) : int
variant  { v, u }
requires { u >= 0 /\ v >= 0 }
ensures  { result = gcd u v }
=
if v > u then binary_gcd v u else
if v = 0 then u else
if mod u 2 = 0 then
if mod v 2 = 0 then 2 * binary_gcd (u / 2) (v / 2)
else binary_gcd (u / 2) v
else
if mod v 2 = 0 then binary_gcd u (v / 2)
else binary_gcd ((u - v) / 2) v

end

```

With machine integers. Note that we assume parameters u, v to be nonnegative. Otherwise, for u = v = min_int, the gcd could not be represented.

```(* does not work with extraction driver ocaml64
module EuclideanAlgorithm31

use mach.int.Int31
use number.Gcd

let rec euclid (u v: int31) : int31
variant  { to_int v }
requires { u >= 0 /\ v >= 0 }
ensures  { result = gcd u v }
=
if v = 0 then
u
else
euclid v (u % v)

end
*)

module EuclideanAlgorithm63

use mach.int.Int63
use number.Gcd

let rec euclid (u v: int63) : int63
variant  { to_int v }
requires { u >= 0 /\ v >= 0 }
ensures  { result = gcd u v }
=
if v = 0 then
u
else
euclid v (u % v)

end
```