-- Coprime factors of a square are squares: a divisor whose square divides
-- a square divides its root, and when a * b is a square with gcd(a, b) = 1
-- both a and b are squares
-- 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.CoprimeSquares
import LCM
import PerfectSquares

-- A divisor of a is coprime to b whenever a is.
theorem coprime.of_divides_left (a b g : Nat)
  requires coprimeAB : gcd(a, b) = 1
  requires divides : g ∣ a
  shows gcd(g, b) = 1
proof
  have every := gcd.eq_one_iff at g, b
  apply every.backward
  intros d dNat dg db
  have da := divides.transitive at d, g, a with dg, divides
  exact coprime.divisor_one at a, b, d with coprimeAB, da, db

-- If d² divides c², then d divides c.
theorem divides.of_mul_self (d c : Nat)
  requires squares : d * d ∣ c * c
  shows d ∣ c
proof
  have split : c = 0 ∨ c ≠ 0 by tauto
  cases split with zero | nonzero
  . rw zero
    exact divides.zero at d
  . have gcdNonzero := gcd.nonzero_right at d, c with nonzero
    obtain ⟨p, pNat, dShape⟩ := gcd.divides_left at d, c
    obtain ⟨q, qNat, cShape⟩ := gcd.divides_right at d, c
    have cofactors := gcd.cofactor_coprime at d, c, p, q
      with gcdNonzero, dShape, cShape
    have squareNonzero := mul.nonzero at gcd(d, c), gcd(d, c)
      with gcdNonzero, gcdNonzero
    have cancel := divides.mul_cancel_left
      at gcd(d, c) * gcd(d, c), p * p, q * q with squareNonzero
    have regroup : gcd(d, c) * gcd(d, c) * (p * p) ∣
        gcd(d, c) * gcd(d, c) * (q * q)
      have leftSide : gcd(d, c) * gcd(d, c) * (p * p) = d * d by ring [dShape]
      have rightSide : gcd(d, c) * gcd(d, c) * (q * q) = c * c by ring [cShape]
      rw [leftSide, rightSide]
      exact squares
    have cofactorSquares := cancel.forward with regroup
    have pDivides := divides.of_mul_left at p, p, q * q with cofactorSquares
    have coprimeSquare := coprime.pow_right at p, q, 2 with cofactors
    have two := pow.two at q
    rw two at coprimeSquare
    have pSelf := divides.reflexive at p
    have pOne := coprime.divisor_one at p, q * q, p
      with coprimeSquare, pSelf, pDivides
    rw pOne at dShape
    have unit : gcd(d, c) * 1 = gcd(d, c) by ring
    rw unit at dShape
    rw ← dShape
    exact gcd.divides_right at d, c

-- A coprime factor of a square is a square: with g = gcd(a, c), a = g p and
-- c = g q, the equation a b = c² reduces to p b = g q², and the coprimality
-- of p with q² and of g with b force p ∣ g and g ∣ p, so a = g².
theorem IsSquare.left_of_coprime_mul (a b c : Nat)
  requires coprimeAB : gcd(a, b) = 1
  requires product : a * b = c * c
  shows IsSquare(a)
proof
  have split : a = 0 ∨ a ≠ 0 by tauto
  cases split with zero | nonzero
  . rw zero
    exact IsSquare.zero
  . have gcdNonzero := gcd.nonzero at a, c with nonzero
    obtain ⟨p, pNat, aShape⟩ := gcd.divides_left at a, c
    obtain ⟨q, qNat, cShape⟩ := gcd.divides_right at a, c
    have cofactors := gcd.cofactor_coprime at a, c, p, q
      with gcdNonzero, aShape, cShape
    have reduced : p * b = gcd(a, c) * (q * q)
      have expanded : gcd(a, c) * (p * b) = gcd(a, c) * (gcd(a, c) * (q * q))
        have leftSide : gcd(a, c) * (p * b) = a * b by ring [aShape]
        have rightSide : gcd(a, c) * (gcd(a, c) * (q * q)) = c * c
          ring [cShape]
        rw [leftSide, rightSide]
        exact product
      exact mul.left_cancel_positive at p * b, gcd(a, c) * (q * q), gcd(a, c)
        with gcdNonzero, expanded
    have coprimeSquare := coprime.pow_right at p, q, 2 with cofactors
    have two := pow.two at q
    rw two at coprimeSquare
    have pDividesProduct : p ∣ gcd(a, c) * (q * q)
      rw ← reduced
      exact divides.mul_left_self at p, b
    have coprimeSquareP := coprime.symmetric at p, q * q with coprimeSquare
    have pDivides := coprime.cancel_right at q * q, p, gcd(a, c)
      with coprimeSquareP, pDividesProduct
    have gDividesA := gcd.divides_left at a, c
    have coprimeGB := coprime.of_divides_left at a, b, gcd(a, c)
      with coprimeAB, gDividesA
    have gDividesProduct : gcd(a, c) ∣ p * b
      rw reduced
      exact divides.mul_left_self at gcd(a, c), q * q
    have coprimeBG := coprime.symmetric at gcd(a, c), b with coprimeGB
    have gDivides := coprime.cancel_right at b, gcd(a, c), p
      with coprimeBG, gDividesProduct
    have equal := divides.antisymmetric at p, gcd(a, c) with pDivides, gDivides
    rw equal at aShape
    witness gcd(a, c)
    symm
    exact aShape

theorem IsSquare.of_coprime_mul (a b c : Nat)
  requires coprimeAB : gcd(a, b) = 1
  requires product : a * b = c * c
  shows IsSquare(a) ∧ IsSquare(b)
proof
  constructor
  . exact IsSquare.left_of_coprime_mul at a, b, c with coprimeAB, product
  . have coprimeBA := coprime.symmetric at a, b with coprimeAB
    have commuted : b * a = c * c
      rw ← product
      ring
    exact IsSquare.left_of_coprime_mul at b, a, c with coprimeBA, commuted