-- Euclidean distance in the plane and the triangle inequality.
-- 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.Geometry.Distance
import Plane
import Real.Sqrt

namespace Euclidean
export open Euclidean

-- Section: Cauchy-Schwarz in the plane

-- The difference of the two sides is the square of a cross term.
theorem angleDot.cauchy_schwarz (a b c : Point)
  shows angleDot(a, b, c) * angleDot(a, b, c) ≤ Euclidean.distanceSquared(b, a) * Euclidean.distanceSquared(b, c)
proof
  have dotRule := angleDot.coordinates at a, b, c with *
  have firstRule := Euclidean.distanceSquared.coordinates at b, a with *
  have secondRule := Euclidean.distanceSquared.coordinates at b, c with *
  have crossNonnegative := realMulValue.square_nonnegative at (x(a) - x(b)) * (y(c) - y(b)) - (y(a) - y(b)) * (x(c) - x(b)) with *
  have difference : Euclidean.distanceSquared(b, a) * Euclidean.distanceSquared(b, c) - angleDot(a, b, c) * angleDot(a, b, c) =
      ((x(a) - x(b)) * (y(c) - y(b)) - (y(a) - y(b)) * (x(c) - x(b))) * ((x(a) - x(b)) * (y(c) - y(b)) - (y(a) - y(b)) * (x(c) - x(b)))
    rw [dotRule, firstRule, secondRule]
    ring
  have nonnegativeDifference : 0 ≤ Euclidean.distanceSquared(b, a) * Euclidean.distanceSquared(b, c) - angleDot(a, b, c) * angleDot(a, b, c)
    rwa difference
  exact realLe.of_nonnegative_difference with *

-- Section: Distance

theorem Euclidean.distanceSquared.nonnegative_real (a b : Point)
  shows NonnegativeReal(Euclidean.distanceSquared(a, b))
proof
  have le := Euclidean.distanceSquared.nonnegative at a, b with *
  exact NonnegativeReal.of_le with *

def distance (a b : Point) : Real := sqrt(Euclidean.distanceSquared(a, b))

theorem distance.real (a b : Point)
  shows Real(distance(a, b))
proof
  typecheck

theorem distance.nonnegative (a b : Point)
  shows 0 ≤ distance(a, b)
proof
  have nonnegative := Euclidean.distanceSquared.nonnegative_real at a, b with *
  have root := sqrt.nonnegative with nonnegative
  unfold distance
  exact root

theorem distance.square (a b : Point)
  shows distance(a, b) * distance(a, b) = Euclidean.distanceSquared(a, b)
proof
  have nonnegative := Euclidean.distanceSquared.nonnegative_real at a, b with *
  have square := sqrt.square with nonnegative
  unfold distance
  cc

theorem distance.symmetric (a b : Point)
  shows distance(a, b) = distance(b, a)
proof
  have swap := Euclidean.distanceSquared.symmetric at a, b with *
  unfold distance
  cc

theorem distance.eq_zero_iff (a b : Point)
  shows distance(a, b) = 0 ↔ a = b
proof
  have squared := Euclidean.distanceSquared.zero_iff at a, b with *
  have square := distance.square at a, b with *
  constructor
  . intro zero
    have squareZero : Euclidean.distanceSquared(a, b) = 0
      rw [← square, zero]
      exact realMulValue.zero_left with realZero.real
    exact squared.forward with squareZero
  . intro equal
    have squareZero := squared.backward with equal
    have root : distance(a, b) * distance(a, b) = 0 by rwa square
    exact realMulValue.square_eq_zero_iff.forward with root

theorem distance.positive_of_ne (a b : Point)
  requires distinct : a ≠ b
  shows 0 < distance(a, b)
proof
  have nonnegative := Euclidean.distanceSquared.nonnegative_real at a, b with *
  have positiveSquare := Euclidean.distanceSquared.positive with *
  have nonzero : Euclidean.distanceSquared(a, b) ≠ 0
    intro zero
    rw zero at positiveSquare
    exact realLess.not_self with realZero.real, positiveSquare
  have root := sqrt.positive with nonnegative, nonzero
  unfold distance
  exact root

-- Section: Triangle inequality

-- Cauchy-Schwarz as a lower bound on the dot product of two rays.
theorem angleDot.lower_bound (a b c : Point)
  shows 0 ≤ angleDot(a, b, c) + distance(b, a) * distance(b, c)
proof
  have dotReal := angleDot.real at a, b, c with *
  have firstNonnegative := distance.nonnegative at b, a with *
  have secondNonnegative := distance.nonnegative at b, c with *
  have productNonnegative := realMulValue.nonnegative at distance(b, a), distance(b, c) with *
  have firstSquare := distance.square at b, a with *
  have secondSquare := distance.square at b, c with *
  have productSquare : (distance(b, a) * distance(b, c)) * (distance(b, a) * distance(b, c)) = Euclidean.distanceSquared(b, a) * Euclidean.distanceSquared(b, c)
    rw [← firstSquare, ← secondSquare]
    ring
  have schwarz := angleDot.cauchy_schwarz at a, b, c with *
  rw ← productSquare at schwarz
  have sign := Real.sign_cases with dotReal
  cases sign with nonpositive | nonnegative
  . have negReal := realNegValue.real at angleDot(a, b, c) with dotReal
    have negNonnegative := realNegValue.nonnegative_of_nonpositive with nonpositive
    have negLe := realLe.of_subset with realZero.real, negNonnegative
    have negSquare := realMulValue.neg_both at angleDot(a, b, c), angleDot(a, b, c) with *
    have negSchwarz : (-angleDot(a, b, c)) * (-angleDot(a, b, c)) ≤ (distance(b, a) * distance(b, c)) * (distance(b, a) * distance(b, c))
      rwa negSquare
    have negLeProduct := realMulValue.le_of_square_le at -angleDot(a, b, c), distance(b, a) * distance(b, c) with *
    have shifted := realLe.add_right at -angleDot(a, b, c), distance(b, a) * distance(b, c), angleDot(a, b, c) with *
    have cancel : -angleDot(a, b, c) + angleDot(a, b, c) = 0 by ring
    rw cancel at shifted
    have swap : distance(b, a) * distance(b, c) + angleDot(a, b, c) = angleDot(a, b, c) + distance(b, a) * distance(b, c)
      ring
    rwa swap at shifted
  . have dotLe := realLe.of_subset at realZero, angleDot(a, b, c) with realZero.real, nonnegative
    have sum := realLe.add at realZero, angleDot(a, b, c), realZero, distance(b, a) * distance(b, c) with *
    have zeroZero := realAddValue.zero_right with realZero.real
    rwa zeroZero at sum

theorem distance.triangle (a b c : Point)
  shows distance(a, c) ≤ distance(a, b) + distance(b, c)
proof
  have acNonnegative := distance.nonnegative at a, c with *
  have abNonnegative := distance.nonnegative at a, b with *
  have bcNonnegative := distance.nonnegative at b, c with *
  have sumNonnegative : 0 ≤ distance(a, b) + distance(b, c)
    have step := realLe.add at realZero, distance(a, b), realZero, distance(b, c) with *
    have zeroZero := realAddValue.zero_right with realZero.real
    rwa zeroZero at step
  have acSquare := distance.square at a, c with *
  have abSquare := distance.square at a, b with *
  have bcSquare := distance.square at b, c with *
  have symmetric := distance.symmetric at b, a with *
  have bound := angleDot.lower_bound at a, b, c with *
  rw symmetric at bound
  have identity := cosine_identity at a, b, c with *
  have doubled : 0 ≤ (angleDot(a, b, c) + distance(a, b) * distance(b, c)) + (angleDot(a, b, c) + distance(a, b) * distance(b, c))
    have step := realLe.add at realZero, angleDot(a, b, c) + distance(a, b) * distance(b, c), realZero, angleDot(a, b, c) + distance(a, b) * distance(b, c) with *
    have zeroZero := realAddValue.zero_right with realZero.real
    rwa zeroZero at step
  have difference : (distance(a, b) + distance(b, c)) * (distance(a, b) + distance(b, c)) - distance(a, c) * distance(a, c) =
      (angleDot(a, b, c) + distance(a, b) * distance(b, c)) + (angleDot(a, b, c) + distance(a, b) * distance(b, c))
    linear_combination abSquare + bcSquare - acSquare - identity
  have nonnegativeDifference : 0 ≤ (distance(a, b) + distance(b, c)) * (distance(a, b) + distance(b, c)) - distance(a, c) * distance(a, c)
    rwa difference
  have squares := realLe.of_nonnegative_difference with *
  exact realMulValue.le_of_square_le with *

end Euclidean