-- Greatest common divisors and Euclidean transformations
-- Part of Cumulus, the Grace17 mathematics library.
-- Copyright (c) 2026 Grant Passmore.
-- Released under the Apache License 2.0. See LICENSE in the repository root.
-- Authors: Grant Passmore

module Grace.Cumulus.NumberTheory.GCD
import Divisibility

-- Greatest in the divisibility order. This includes IsGCD(0, 0, 0).
-- Existence and an executable algorithm are separate from this specification.
predicate IsGCD (a b g : Nat) :=
  g ∣ a ∧ g ∣ b ∧ (∀ d, Nat(d) → d ∣ a → d ∣ b → d ∣ g)

theorem IsGCD.divides_left (a b g : Nat)
  requires gcd : IsGCD(a, b, g)
  shows g ∣ a
proof
  simp [IsGCD]

theorem IsGCD.divides_right (a b g : Nat)
  requires gcd : IsGCD(a, b, g)
  shows g ∣ b
proof
  simp [IsGCD]

theorem IsGCD.greatest (a b g d : Nat)
  requires gcd : IsGCD(a, b, g)
  requires left : d ∣ a
  requires right : d ∣ b
  shows d ∣ g
proof
  simp [IsGCD]

theorem IsGCD.intro (a b g : Nat)
  requires left : g ∣ a
  requires right : g ∣ b
  requires greatest : ∀ d, Nat(d) → d ∣ a → d ∣ b → d ∣ g
  shows IsGCD(a, b, g)
proof
  simp [IsGCD]

theorem IsGCD.unique (a b g h : Nat)
  requires first : IsGCD(a, b, g)
  requires second : IsGCD(a, b, h)
  shows g = h
proof
  have ga := IsGCD.divides_left at a, b, g with *
  have gb := IsGCD.divides_right at a, b, g with *
  have ha := IsGCD.divides_left at a, b, h with *
  have hb := IsGCD.divides_right at a, b, h with *
  have gh := IsGCD.greatest at a, b, h, g with *
  have hg := IsGCD.greatest at a, b, g, h with *
  exact divides.antisymmetric with *

theorem IsGCD.symmetric (a b g : Nat)
  requires gcd : IsGCD(a, b, g)
  shows IsGCD(b, a, g)
proof
  have left := IsGCD.divides_left with *
  have right := IsGCD.divides_right with *
  apply IsGCD.intro with right, left
  intros d dNat db da
  exact IsGCD.greatest with *

theorem IsGCD.zero_right (a : Nat)
  shows IsGCD(a, 0, a)
proof
  have self := divides.reflexive with *
  have zero := divides.zero with *
  apply IsGCD.intro with self, zero
  intros d dNat da dz
  exact da

theorem IsGCD.zero_left (a : Nat)
  shows IsGCD(0, a, a)
proof
  have right := IsGCD.zero_right with *
  exact IsGCD.symmetric with *

theorem IsGCD.of_divides (a b : Nat)
  requires divides : a ∣ b
  shows IsGCD(a, b, a)
proof
  have self := divides.reflexive at a with *
  apply IsGCD.intro with self, divides
  intros d dNat da db
  exact da

-- A single invariant covers subtraction and quotient/remainder algorithms.
-- No remainder bound is needed for preservation, only for termination.
theorem IsGCD.euclidean_forward (a b q r g : Nat)
  requires equation : a = b * q + r
  requires gcd : IsGCD(a, b, g)
  shows IsGCD(b, r, g)
proof
  have ga := IsGCD.divides_left with *
  have gb := IsGCD.divides_right with *
  have equivalence := divides.remainder_iff with *
  have gr : g ∣ r
  . tauto
  apply IsGCD.intro with gb, gr
  intros d dNat db dr
  have commonDivisor := divides.remainder_iff at d, a, b, q, r with *
  have da : d ∣ a
  . tauto
  exact IsGCD.greatest with *

theorem IsGCD.euclidean_backward (a b q r g : Nat)
  requires equation : a = b * q + r
  requires gcd : IsGCD(b, r, g)
  shows IsGCD(a, b, g)
proof
  have gb := IsGCD.divides_left with *
  have gr := IsGCD.divides_right with *
  have equivalence := divides.remainder_iff with *
  have ga : g ∣ a
  . tauto
  apply IsGCD.intro with ga, gb
  intros d dNat da db
  have commonDivisor := divides.remainder_iff at d, a, b, q, r with *
  have dr : d ∣ r
  . tauto
  exact IsGCD.greatest with *

theorem IsGCD.euclidean_iff (a b q r g : Nat)
  requires equation : a = b * q + r
  shows IsGCD(a, b, g) ↔ IsGCD(b, r, g)
proof
  constructor
  . intro gcd
    exact IsGCD.euclidean_forward with *
  . intro gcd
    exact IsGCD.euclidean_backward with *

theorem IsGCD.add_right_iff (a b g : Nat)
  shows IsGCD(a + b, b, g) ↔ IsGCD(a, b, g)
proof
  have equation : a + b = b * 1 + a
  . ring
  have step := IsGCD.euclidean_iff at a + b, b, 1, a, g with *
  constructor
  . intro gcd
    have reversed : IsGCD(b, a, g)
    . tauto
    exact IsGCD.symmetric with *
  . intro gcd
    have reversed := IsGCD.symmetric with *
    tauto