OrientedEllipsoid.lean

1 import KomorebiLSM.EllipsoidGeometry
2 import KomorebiLSM.EllipsoidDualCS
3 import KomorebiLSM.PolynomialGram
4 
5 /-! Oriented quadratic certificates on the actual ellipsoid cone. -/
6 set_option autoImplicit false
7 namespace Komorebi.LSM.OrientedEllipsoid
8 open scoped BigOperators Matrix
9 variable {D : Type*} [Fintype D] [DecidableEq D]
10 
11 /-- The literal mass-first centering matrix from the monograph. -/
12 def centerMatrix (c : D → ℝ) : Matrix (Unit ⊕ D) (Unit ⊕ D) ℝ :=
13   Matrix.fromBlocks 1 0 (fun i _ => -c i) 1
14 
15 /-- The uncentered quadratic form has one positive mass coordinate. -/
16 def diagonalForm (H : Matrix D D ℝ) : Matrix (Unit ⊕ D) (Unit ⊕ D) ℝ :=
17   Matrix.fromBlocks 1 0 0 (-H)
18 
19 /-- This is literally S_cᵀ diag(1,-H) S_c. -/
20 def Q (c : D → ℝ) (H : Matrix D D ℝ) : Matrix (Unit ⊕ D) (Unit ⊕ D) ℝ :=
21   (centerMatrix c).transpose * diagonalForm H * centerMatrix c
22 
23 theorem centerMatrix_mulVec (c : D → ℝ) (z : State D) :
24     centerMatrix c *ᵥ z = point (mass z) (centered c z) := by
25   ext i
26   rcases i with i | i
27   · cases i
28     simp [Matrix.fromBlocks,centerMatrix,point,Matrix.mulVec,dotProduct,mass,Fintype.sum_sum_type]
29   · simp [Matrix.one_apply,Matrix.fromBlocks,centerMatrix,point,centered,spatial,mass,Matrix.mulVec,dotProduct,
30       Fintype.sum_sum_type,sub_eq_add_neg,mul_comm,add_comm]
31 
32 omit [DecidableEq D] in
33 theorem diagonalForm_mulVec (H : Matrix D D ℝ) (m : ℝ) (x : D → ℝ) :
34     diagonalForm H *ᵥ point m x = point m (-(H *ᵥ x)) := by
35   rw [diagonalForm,Matrix.fromBlocks_mulVec]
36   ext i
37   rcases i with i | i
38   · cases i
39     simp [point]
40   · simp [point,Matrix.neg_mulVec]
41 
42 omit [DecidableEq D] in
43 theorem point_dotProduct (m n : ℝ) (x y : D → ℝ) :
44     point m x ⬝ᵥ point n y = m*n + x ⬝ᵥ y := by
45   simp [dotProduct,Fintype.sum_sum_type,point]
46 
47 /-- The centered quadratic matrix represents exactly m²-(x-mc)ᵀH(x-mc). -/
48 theorem quadratic_form (c : D → ℝ) (H : Matrix D D ℝ) (z : State D) :
49     z ⬝ᵥ (Q c H *ᵥ z) = (mass z)^2 - centered c z ⬝ᵥ (H *ᵥ centered c z) := by
50   rw [Q,← Matrix.mulVec_mulVec,← Matrix.mulVec_mulVec,
51     Matrix.dotProduct_transpose_mulVec,centerMatrix_mulVec,
52     diagonalForm_mulVec,point_dotProduct,neg_dotProduct,dotProduct_comm (H *ᵥ centered c z)]
53   ring
54 
55 theorem cone_iff_quadratic (c : D → ℝ) (H : Matrix D D ℝ) (z : State D) :
56     z∈cone c H ↔ 0≤mass z ∧ 0≤z ⬝ᵥ (Q c H *ᵥ z) := by
57   rw [quadratic_form]
58   exact and_congr_right (fun _ => sub_nonneg.symm)
59 
60 def massConstant (B : Matrix (Unit ⊕ D) (Unit ⊕ D) ℝ) : ℝ := B (Sum.inl ()) (Sum.inl ())
61 def massSpatial (B : Matrix (Unit ⊕ D) (Unit ⊕ D) ℝ) (i : D) : ℝ := B (Sum.inl ()) (Sum.inr i)
62 def beta (c : D → ℝ) (B : Matrix (Unit ⊕ D) (Unit ⊕ D) ℝ) : ℝ :=
63   massConstant B + massSpatial B ⬝ᵥ c
64 
65 omit [DecidableEq D] in
66 theorem mass_mulVec (B : Matrix (Unit ⊕ D) (Unit ⊕ D) ℝ) (z : State D) :
67     mass (B *ᵥ z) = massConstant B * mass z + massSpatial B ⬝ᵥ spatial z := by
68   simp [mass,massConstant,massSpatial,spatial,Matrix.mulVec,dotProduct,Fintype.sum_sum_type]
69 
70 omit [DecidableEq D] in
71 theorem mass_mulVec_centered (c : D → ℝ)
72     (B : Matrix (Unit ⊕ D) (Unit ⊕ D) ℝ) (z : State D) :
73     mass (B *ᵥ z) = beta c B * mass z + massSpatial B ⬝ᵥ centered c z := by
74   rw [mass_mulVec,beta,centered_apply,dotProduct_sub,dotProduct_smul]
75   ring
76 
77 /-- The actual matrix identity propagates the quadratic half of the cone. -/
78 theorem quadratic_output_nonneg (c : D → ℝ) (H : Matrix D D ℝ)
79     (B R : Matrix (Unit ⊕ D) (Unit ⊕ D) ℝ) (tau : ℝ)
80     (htau : 0≤tau) (hR : R.PosSemidef)
81     (hidentity : B.transpose * Q c H * B - tau • Q c H = R)
82     {z : State D} (hz : 0≤z ⬝ᵥ (Q c H *ᵥ z)) :
83     0≤(B *ᵥ z) ⬝ᵥ (Q c H *ᵥ (B *ᵥ z)) := by
84   have hr : 0≤z ⬝ᵥ (R *ᵥ z) := by simpa only [star_trivial] using hR.dotProduct_mulVec_nonneg z
85   have he := congrArg (fun M : Matrix (Unit ⊕ D) (Unit ⊕ D) ℝ => z ⬝ᵥ (M *ᵥ z)) hidentity
86   simp only [Matrix.sub_mulVec,Matrix.smul_mulVec,dotProduct_sub,dotProduct_smul] at he
87   rw [← Matrix.mulVec_mulVec,← Matrix.mulVec_mulVec,Matrix.dotProduct_transpose_mulVec] at he
88   rw [dotProduct_comm (Q c H *ᵥ (B *ᵥ z))] at he
89   simp only [smul_eq_mul] at he
90   nlinarith [mul_nonneg htau hz]
91 
92 /-- Literal quadratic and orientation witnesses for one actual branch matrix. -/
93 structure BranchCertificate (c : D → ℝ) (H : Matrix D D ℝ)
94     (B : Matrix (Unit ⊕ D) (Unit ⊕ D) ℝ) where
95   tau : ℝ
96   remainder : Matrix (Unit ⊕ D) (Unit ⊕ D) ℝ
97   tau_nonneg : 0≤tau
98   remainder_positive : remainder.PosSemidef
99   identity : B.transpose * Q c H * B - tau • Q c H = remainder
100   beta_nonneg : 0≤beta c B
101   orientation : 0≤(beta c B)^2 - massSpatial B ⬝ᵥ (H⁻¹ *ᵥ massSpatial B)
102 
103 /-- The orientation witness selects the forward sheet, including zero-mass inputs. -/
104 theorem BranchCertificate.preserves {c : D → ℝ} {H : Matrix D D ℝ} (hH : H.PosDef)
105     {B : Matrix (Unit ⊕ D) (Unit ⊕ D) ℝ} (cert : BranchCertificate c H B)
106     {z : State D} (hz : z∈cone c H) : B *ᵥ z ∈ cone c H := by
107   apply (cone_iff_quadratic c H _).mpr
108   constructor
109   · rw [mass_mulVec_centered]
110     exact mass_orientation H hH (massSpatial B) (centered c z) (beta c B) (mass z)
111       hz.1 hz.2 cert.beta_nonneg cert.orientation
112   · exact quadratic_output_nonneg c H B cert.remainder cert.tau cert.tau_nonneg
113       cert.remainder_positive cert.identity ((cone_iff_quadratic c H z).mp hz).2
114 
115 end Komorebi.LSM.OrientedEllipsoid

Back to source controls