-- Cauchy–Schwarz for arbitrary finite-dimensional real vectors
-- 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

-- Mathematical source: imandra-ai/imandra-100-theorems, src/cauchy_schwarz.iml.
-- Sequences restricted to n coordinates replace the source's equal-length lists.

module Grace.Cumulus.LinearAlgebra.CauchySchwarz
import FiniteDimensional

theorem dotProduct.zero_of_squaredNorm_zero (u v : Set, n : Nat)
  requires sequenceU : RealSequence(u)
  requires sequenceV : RealSequence(v)
  shows squaredNorm(v, n) = (0:R) → dotProduct(u, v, n) = (0:R)
proof
  induction n
  | zero =>
    intro zero
    exact dotProduct.zero at u, v
  | succ(k) using ih =>
    intro zero
    have uk := RealSequence.value at u, k with *
    have vk := RealSequence.value at v, k with *
    have nonnegative := squaredNorm.nonnegative at v, k with *
    have squareNonnegative := realMulValue.square_nonnegative at v[k] with *
    rw squaredNorm.succ at zero
    have prefixZero := realAddValue.nonnegative_eq_zero at squaredNorm(v, k), v[k] * v[k] with *
    have squareZero : (v[k]:R) * (v[k]:R) = 0
    . rw prefixZero at zero
      simpa [realAddValue.zero_left] using zero
    have valueZero := realMulValue.square_eq_zero_iff.forward at v[k] with *
    have dotZero := ih with prefixZero
    rw [dotProduct.succ, dotZero, valueZero]
    ring

theorem squaredNorm.scaled_difference (u v : Set, C D : Real, n : Nat)
  requires sequenceU : RealSequence(u)
  requires sequenceV : RealSequence(v)
  shows squaredNorm(subSequence(mulSequence(constantSequence(C), u), mulSequence(constantSequence(D), v)), n) =
    C * C * (squaredNorm(u, n):R) - (2:R) * C * D * (dotProduct(u, v, n):R) + D * D * (squaredNorm(v, n):R)
proof
  have constantC := constantSequence.sequence at C with *
  have constantD := constantSequence.sequence at D with *
  have scaledU := mulSequence.sequence at constantSequence(C), u with *
  have scaledV := mulSequence.sequence at constantSequence(D), v with *
  simp [squaredNorm.sub, squaredNorm.scale, dotProduct.scale_left, dotProduct.scale_right]
  ring

theorem dotProduct.cauchy_schwarz (u v : Set, n : Nat)
  requires sequenceU : RealSequence(u)
  requires sequenceV : RealSequence(v)
  shows (dotProduct(u, v, n):R) * (dotProduct(u, v, n):R) ≤
    (squaredNorm(u, n):R) * (squaredNorm(v, n):R)
proof
  have normNonnegative := squaredNorm.nonnegative at v, n with *
  have cases : squaredNorm(v, n) = (0:R) ∨ squaredNorm(v, n) ≠ (0:R)
  . tauto
  cases cases with zero | nonzero
  . have dotZero := dotProduct.zero_of_squaredNorm_zero at u, v, n with *
    rw [dotZero, zero]
    simpa [realMulValue.zero_left, realMulValue.zero_right] using realLe.reflexive at realZero with *
  . have positive := Real.positive_of_nonnegative_of_ne_zero with *
    have constantNorm := constantSequence.sequence at squaredNorm(v, n) with *
    have constantDot := constantSequence.sequence at dotProduct(u, v, n) with *
    have scaledU := mulSequence.sequence at constantSequence(squaredNorm(v, n)), u with *
    have scaledV := mulSequence.sequence at constantSequence(dotProduct(u, v, n)), v with *
    have difference := subSequence.sequence at mulSequence(constantSequence(squaredNorm(v, n)), u), mulSequence(constantSequence(dotProduct(u, v, n)), v) with *
    have nonnegative := squaredNorm.nonnegative at subSequence(mulSequence(constantSequence(squaredNorm(v, n)), u), mulSequence(constantSequence(dotProduct(u, v, n)), v)), n with *
    have expansion := squaredNorm.scaled_difference at u, v, squaredNorm(v, n), dotProduct(u, v, n), n with *
    rw expansion at nonnegative
    have rearrange : (squaredNorm(v, n):R) * (squaredNorm(v, n):R) * (squaredNorm(u, n):R) - (2:R) * (squaredNorm(v, n):R) * (dotProduct(u, v, n):R) * (dotProduct(u, v, n):R) + (dotProduct(u, v, n):R) * (dotProduct(u, v, n):R) * (squaredNorm(v, n):R) =
      ((squaredNorm(u, n):R) * (squaredNorm(v, n):R)) * (squaredNorm(v, n):R) - ((dotProduct(u, v, n):R) * (dotProduct(u, v, n):R)) * (squaredNorm(v, n):R)
    . ring
    rw rearrange at nonnegative
    have products := realLe.of_nonnegative_difference with *
    exact realLe.mul_right_cancel_positive with *

theorem dotProduct.cauchy_schwarz_scalar_equality (v : Set, C : Real, n : Nat)
  requires sequence : RealSequence(v)
  shows (dotProduct(mulSequence(constantSequence(C), v), v, n):R) *
      (dotProduct(mulSequence(constantSequence(C), v), v, n):R) =
    (squaredNorm(mulSequence(constantSequence(C), v), n):R) * (squaredNorm(v, n):R)
proof
  have constants := constantSequence.sequence with *
  have scaled := mulSequence.sequence at constantSequence(C), v with *
  simp [dotProduct.scale_left, squaredNorm.scale]
  unfold squaredNorm
  ring

-- A zero norm is exactly coordinatewise zero on the observed prefix.
theorem squaredNorm.zero_prefix (u : Set, n : Nat)
  requires sequence : RealSequence(u)
  shows squaredNorm(u, n) = (0:R) → ∀ i ∈ ℕ, i < n → u[i] = (0:R)
proof
  induction n
  | zero =>
    intros zero i member before
    have naturalI := nat.omega.backward with member
    have impossible : i < i
    . arith
    exfalso
    exact Less.irreflexive with naturalI, impossible
  | succ(k) using ih =>
    intros zero i member before
    have naturalI := nat.omega.backward with member
    have uk := RealSequence.value at u, k with *
    have nonnegative := squaredNorm.nonnegative at u, k with *
    have squareNonnegative := realMulValue.square_nonnegative at u[k] with *
    rw squaredNorm.succ at zero
    have prefixZero := realAddValue.nonnegative_eq_zero at squaredNorm(u, k), u[k] * u[k] with *
    have squareZero : (u[k]:R) * (u[k]:R) = 0
    . rw prefixZero at zero
      simpa [realAddValue.zero_left] using zero
    have valueZero := realMulValue.square_eq_zero_iff.forward at u[k] with *
    have included := membership.of_Less at i, succ(k) with *, before
    have alternatives := membership.successor_forward{n := k} with included
    cases alternatives with earlier | equal
    . have smaller := Less.of_membership at i, k with *, earlier
      exact ih with prefixZero, i, member, smaller
    . rw equal
      exact valueZero

theorem squaredNorm.of_zero_prefix (u : Set, n : Nat)
  requires sequence : RealSequence(u)
  shows (∀ i ∈ ℕ, i < n → u[i] = (0:R)) → squaredNorm(u, n) = (0:R)
proof
  induction n
  | zero =>
    intro coordinates
    unfold squaredNorm
    exact dotProduct.zero at u, u
  | succ(k) using ih =>
    intro coordinates
    have prefix : ∀ i ∈ ℕ, i < k → u[i] = (0:R)
    . intros i member before
      have naturalI := nat.omega.backward with member
      have smaller : i < succ(k)
      . arith
      exact coordinates with i, member, smaller
    have zeroPrefix := ih with prefix
    have memberK := nat.omega.forward at k with *
    have beforeK : k < succ(k)
    . arith
    have zeroLast := coordinates with k, memberK, beforeK
    rw [squaredNorm.succ, zeroPrefix, zeroLast]
    ring

theorem squaredNorm.eq_zero_iff (u : Set, n : Nat)
  requires sequence : RealSequence(u)
  shows squaredNorm(u, n) = (0:R) ↔ (∀ i ∈ ℕ, i < n → u[i] = (0:R))
proof
  constructor
  . exact squaredNorm.zero_prefix with *
  . exact squaredNorm.of_zero_prefix with *

theorem squaredNorm.eq_prefix (u v : Set, n : Nat)
  requires sequenceU : RealSequence(u)
  requires sequenceV : RealSequence(v)
  requires coordinates : ∀ i ∈ ℕ, i < n → u[i] = v[i]
  shows squaredNorm(u, n) = squaredNorm(v, n)
proof
  unfold squaredNorm
  exact dotProduct.eq_prefix with *

theorem dotProduct.cauchy_schwarz_equality_residual (u v : Set, n : Nat)
  requires sequenceU : RealSequence(u)
  requires sequenceV : RealSequence(v)
  requires equality : (dotProduct(u, v, n):R) * (dotProduct(u, v, n):R) =
    (squaredNorm(u, n):R) * (squaredNorm(v, n):R)
  shows squaredNorm(subSequence(mulSequence(constantSequence(squaredNorm(v, n)), u),
    mulSequence(constantSequence(dotProduct(u, v, n)), v)), n) = (0:R)
proof
  rw squaredNorm.scaled_difference
  linear_combination -(squaredNorm(v, n):R) * equality

theorem dotProduct.cauchy_schwarz_equality_iff (u v : Set, n : Nat)
  requires sequenceU : RealSequence(u)
  requires sequenceV : RealSequence(v)
  shows (dotProduct(u, v, n):R) * (dotProduct(u, v, n):R) =
      (squaredNorm(u, n):R) * (squaredNorm(v, n):R) ↔
    (squaredNorm(v, n) = (0:R) ∨
      (∃ C ∈ ℝ, (∀ i ∈ ℕ, i < n → u[i] = C * (v[i]:R))))
proof
  constructor
  . intro equality
    have alternatives : squaredNorm(v, n) = (0:R) ∨ squaredNorm(v, n) ≠ (0:R)
    . tauto
    cases alternatives with zero | nonzero
    . left
      exact zero
    . right
      have inverseLaw := realInvValue.mul_right with *
      have coefficient := realMulValue.real at dotProduct(u, v, n), realInvValue(squaredNorm(v, n)) with *
      have coefficientMember : dotProduct(u, v, n) * realInvValue(squaredNorm(v, n)) ∈ ℝ
      . settle with Real
      refine ⟨dotProduct(u, v, n) * realInvValue(squaredNorm(v, n)), coefficientMember, _⟩
      intros i member before
      have naturalI := nat.omega.backward with member
      have ui := RealSequence.value at u, i with *
      have vi := RealSequence.value at v, i with *
      have constantNorm := constantSequence.sequence at squaredNorm(v, n) with *
      have constantDot := constantSequence.sequence at dotProduct(u, v, n) with *
      have scaledU := mulSequence.sequence at constantSequence(squaredNorm(v, n)), u with *
      have scaledV := mulSequence.sequence at constantSequence(dotProduct(u, v, n)), v with *
      have difference := subSequence.sequence at mulSequence(constantSequence(squaredNorm(v, n)), u), mulSequence(constantSequence(dotProduct(u, v, n)), v) with *
      have normZero := dotProduct.cauchy_schwarz_equality_residual with *
      have pointZero := squaredNorm.zero_prefix with *, normZero, i, member, before
      have firstValue := RealSequence.value at mulSequence(constantSequence(squaredNorm(v, n)), u), i with *
      have secondValue := RealSequence.value at mulSequence(constantSequence(dotProduct(u, v, n)), v), i with *
      have normValue := RealSequence.value at constantSequence(squaredNorm(v, n)), i with *
      have dotValue := RealSequence.value at constantSequence(dotProduct(u, v, n)), i with *
      simp [subSequence.value, mulSequence.value, constantSequence.value] at pointZero
      linear_combination realInvValue(squaredNorm(v, n)) * pointZero - (u[i]:R) * inverseLaw
  . intro alternatives
    cases alternatives with zero | proportional
    . have dotZero := dotProduct.zero_of_squaredNorm_zero at u, v, n with *, zero
      rw [dotZero, zero]
      ring
    . obtain ⟨C, memberC, coordinates⟩ := proportional
      have realC : Real(C)
      . settle with Real
      have constants := constantSequence.sequence with *
      have scaled := mulSequence.sequence at constantSequence(C), v with *
      have equalPrefix : ∀ i ∈ ℕ, i < n → u[i] = mulSequence(constantSequence(C), v)[i]
      . intros i member before
        have naturalI := nat.omega.backward with member
        have vi := RealSequence.value at v, i with *
        have ci := RealSequence.value at constantSequence(C), i with *
        simp [mulSequence.value, constantSequence.value]
        exact coordinates with i, member, before
      have reflexive : ∀ i ∈ ℕ, i < n → v[i] = v[i]
      . intros i member before
        rfl
      have equalDot := dotProduct.eq_prefix at u, v, mulSequence(constantSequence(C), v), v, n with *
      have equalNorm := squaredNorm.eq_prefix at u, mulSequence(constantSequence(C), v), n with *
      rw [equalDot, equalNorm]
      exact dotProduct.cauchy_schwarz_scalar_equality with *