Deterministic rangefinder error bound

import RNLAPhaseB.B3Base
import RNLAPhaseB.Complexity

open MeasureTheory
open scoped Matrix
open scoped Matrix.Norms.L2Operator
open scoped BigOperators
open scoped ENNReal

namespace RNLAPhaseB.HMT2011MainErrorBd

/-- Trailing index embedding `Fin (n - k) → Fin n`, `i ↦ k + i`, selecting the last `n - k` indices. -/
def trailingFin (k n : ℕ) : Fin (n - k) → Fin n :=
  fun i => ⟨k + i.val, by have := i.isLt; omega⟩

/-- `A = U Σ Vᴴ` is a singular value decomposition. `U` is `m × p` with orthonormal columns
(full SVD `p = m`, or thin SVD). `Σ` is `p × n` diagonal with real nonnegative, nonincreasing
diagonal entries. `V` is `n × n` unitary. -/
def IsSVD {𝕜 : Type*} [RCLike 𝕜] {m n p : ℕ} (A : Matrix (Fin m) (Fin n) 𝕜)
    (U : Matrix (Fin m) (Fin p) 𝕜) (S : Matrix (Fin p) (Fin n) 𝕜)
    (V : Matrix (Fin n) (Fin n) 𝕜) : Prop :=
  min m n ≤ p ∧ Uᴴ * U = 1 ∧ Vᴴ * V = 1 ∧ V * Vᴴ = 1 ∧ A = U * S * Vᴴ ∧
  (∀ (i : Fin p) (j : Fin n), (i : ℕ) ≠ j → S i j = 0) ∧
  (∀ (i : Fin p) (j : Fin n), (i : ℕ) = j → ∃ s : ℝ, 0 ≤ s ∧ S i j = (s : 𝕜)) ∧
  (∀ (i : Fin p) (j : Fin n) (i' : Fin p) (j' : Fin n),
    (i : ℕ) = j → (i' : ℕ) = j' → (i : ℕ) ≤ i' → ‖S i' j'‖ ≤ ‖S i j‖)

/-- `V₁`: the first `k` columns of `V` (partition [eq:mugpart]). -/
def V1 {𝕜 : Type*} {n : ℕ} (k : ℕ) (hk : k ≤ n) (V : Matrix (Fin n) (Fin n) 𝕜) :
    Matrix (Fin n) (Fin k) 𝕜 :=
  V.submatrix id (Fin.castLE hk)

/-- `V₂`: the last `n - k` columns of `V` (partition [eq:mugpart]). -/
def V2 {𝕜 : Type*} {n : ℕ} (k : ℕ) (V : Matrix (Fin n) (Fin n) 𝕜) :
    Matrix (Fin n) (Fin (n - k)) 𝕜 :=
  V.submatrix id (trailingFin k n)

/-- `Σ₂`: the trailing diagonal block of `Σ` (partition [eq:mugpart]). -/
def Sigma2 {𝕜 : Type*} {p n : ℕ} (k : ℕ) (S : Matrix (Fin p) (Fin n) 𝕜) :
    Matrix (Fin (p - k)) (Fin (n - k)) 𝕜 :=
  S.submatrix (trailingFin k p) (trailingFin k n)

/-- Full row rank: rank equals the number of rows. -/
def FullRowRank {𝕜 : Type*} [Field 𝕜] {r c : ℕ} (M : Matrix (Fin r) (Fin c) 𝕜) : Prop :=
  M.rank = r

/-- `X` is the Moore–Penrose pseudoinverse `M⁺` of `M` (the four Penrose conditions). -/
def IsMoorePenroseInverse {𝕜 : Type*} [RCLike 𝕜] {r c : ℕ}
    (M : Matrix (Fin r) (Fin c) 𝕜) (X : Matrix (Fin c) (Fin r) 𝕜) : Prop :=
  M * X * M = M ∧ X * M * X = X ∧ (M * X)ᴴ = M * X ∧ (X * M)ᴴ = X * M

/-- Frobenius norm `(∑ᵢⱼ |Mᵢⱼ|²)^{1/2}`. -/
noncomputable def frobeniusNorm {𝕜 : Type*} [RCLike 𝕜] {r c : ℕ}
    (M : Matrix (Fin r) (Fin c) 𝕜) : ℝ :=
  Real.sqrt (∑ i, ∑ j, ‖M i j‖ ^ 2)

/-- HMT 2011, Theorem 9.1 (thm:main-error-bd). `A = UΣVᴴ` is an SVD, `k ≥ 0` with `k ≤ n`,
which is required for the partition [eq:mugpart] to make sense. `Ω` is any test matrix and
`Y = AΩ`. Define `Ω₁ = V₁ᴴΩ` and `Ω₂ = V₂ᴴΩ` [eq:def_Omegaj]. If `Ω₁` has full row rank, then
`|||(I - P_Y)A|||² ≤ |||Σ₂|||² + |||Σ₂Ω₂Ω₁⁺|||²`. This holds both for the spectral norm (first
conjunct) and for the Frobenius norm (second conjunct). `Ω₁⁺` is any matrix satisfying the
Penrose conditions for `Ω₁`, i.e. its unique pseudoinverse. -/
theorem main_error_bd {𝕜 : Type*} [RCLike 𝕜] {m n p ℓ : ℕ} (k : ℕ) (hk : k ≤ n)
    (A : Matrix (Fin m) (Fin n) 𝕜) (U : Matrix (Fin m) (Fin p) 𝕜)
    (S : Matrix (Fin p) (Fin n) 𝕜) (V : Matrix (Fin n) (Fin n) 𝕜)
    (hSVD : IsSVD A U S V)
    (Ω : Matrix (Fin n) (Fin ℓ) 𝕜)
    (hrank : FullRowRank ((V1 k hk V)ᴴ * Ω))
    (Ω1pinv : Matrix (Fin ℓ) (Fin k) 𝕜)
    (hpinv : IsMoorePenroseInverse ((V1 k hk V)ᴴ * Ω) Ω1pinv) :
    ‖(1 - orthogonalProjector (A * Ω)) * A‖ ^ 2 ≤
        ‖Sigma2 k S‖ ^ 2 + ‖Sigma2 k S * ((V2 k V)ᴴ * Ω) * Ω1pinv‖ ^ 2 ∧
    frobeniusNorm ((1 - orthogonalProjector (A * Ω)) * A) ^ 2 ≤
        frobeniusNorm (Sigma2 k S) ^ 2 +
          frobeniusNorm (Sigma2 k S * ((V2 k V)ᴴ * Ω) * Ω1pinv) ^ 2 := by
  sorry

end RNLAPhaseB.HMT2011MainErrorBd