Ritz-projection approximation in L2 and H1

import PDETextbook.Base

open MeasureTheory
open scoped BigOperators ENNReal

noncomputable section

namespace PDETextbook.LarssonThomee.Thm55

abbrev E2 := Euclidean 2

/-- Polygonal domain in R^2: bounded open connected set whose boundary is a closed polygon. -/
def IsPolygonalDomain (Ω : Set E2) : Prop :=
  IsOpen Ω ∧ IsConnected Ω ∧ Bornology.IsBounded Ω ∧
  ∃ (n : ℕ) (p : Fin n → E2),
    frontier Ω = ⋃ i : Fin n, segment ℝ (p i) (p (finRotate n i))

/-- Closed triangle with vertices K 0, K 1, K 2. -/
def triangle (K : Fin 3 → E2) : Set E2 := convexHull ℝ (Set.range K)

/-- Triangulation of Ω with meshsize h = max_K diam K (conforming). -/
def IsTriangulation (Ω : Set E2) (T : Finset (Fin 3 → E2)) (h : ℝ) : Prop :=
  (∀ K ∈ T, AffineIndependent ℝ K) ∧
  closure Ω = (⋃ K ∈ T, triangle K) ∧
  (∀ K ∈ T, ∀ K' ∈ T, triangle K ≠ triangle K' →
    triangle K ∩ triangle K' = convexHull ℝ (Set.range K ∩ Set.range K')) ∧
  (∀ K ∈ T, Metric.diam (triangle K) ≤ h) ∧
  (∃ K ∈ T, Metric.diam (triangle K) = h)

/-- All angles of K are at least θ0. -/
def MinAngleAtLeast (K : Fin 3 → E2) (θ0 : ℝ) : Prop :=
  ∀ i j k : Fin 3, i ≠ j → j ≠ k → i ≠ k →
    θ0 ≤ EuclideanGeometry.angle (K i) (K j) (K k)

/-- S_h: continuous on closure Ω, linear on each K, vanishing on Γ. -/
def Sh (Ω : Set E2) (T : Finset (Fin 3 → E2)) : Set (E2 → ℝ) :=
  {v | ContinuousOn v (closure Ω) ∧
    (∀ K ∈ T, ∃ A : E2 →ᵃ[ℝ] ℝ, ∀ x ∈ triangle K, v x = A x) ∧
    (∀ x ∈ frontier Ω, v x = 0)}

/-- a(v,w) = ∫_Ω a ∇v · ∇w via weak gradients. -/
def energyForm (Ω : Set E2) (a : E2 → ℝ) (G1 G2 : Fin 2 → E2 → ℝ) : ℝ :=
  ∫ x in Ω, a x * ∑ i, G1 i x * G2 i x

/-- w = R_h v, (5.49): w ∈ S_h and a(w - v, χ) = 0 for all χ ∈ S_h. -/
def IsRitzProjection (Ω : Set E2) (a : E2 → ℝ) (T : Finset (Fin 3 → E2))
    (v w : E2 → ℝ) : Prop :=
  w ∈ Sh Ω T ∧
  ∀ Gv Gw : Fin 2 → E2 → ℝ, HasH1Gradient Ω v Gv → HasH1Gradient Ω w Gw →
    ∀ χ ∈ Sh Ω T, ∀ Gχ : Fin 2 → E2 → ℝ, HasH1Gradient Ω χ Gχ →
      energyForm Ω a (fun i x => Gw i x - Gv i x) Gχ = 0

def l2NormR (Ω : Set E2) (u : E2 → ℝ) : ℝ := Real.sqrt (∫ x in Ω, u x ^ 2)

/-- |u|_1 = ‖∇u‖ for a weak gradient G of u. -/
def h1Seminorm (Ω : Set E2) (G : Fin 2 → E2 → ℝ) : ℝ :=
  Real.sqrt (∫ x in Ω, ∑ i, G i x ^ 2)

def h1NormSq (Ω : Set E2) (u : E2 → ℝ) (G : Fin 2 → E2 → ℝ) : ℝ :=
  (∫ x in Ω, u x ^ 2) + ∫ x in Ω, ∑ i, G i x ^ 2

/-- ‖u‖_2^2 = Σ_{|α|≤2} ‖D^α u‖^2; H i j = ∂_j ∂_i u, each multi-index once. -/
def h2NormSq (Ω : Set E2) (u : E2 → ℝ) (G : Fin 2 → E2 → ℝ)
    (H : Fin 2 → Fin 2 → E2 → ℝ) : ℝ :=
  h1NormSq Ω u G + ∫ x in Ω, ∑ i, ∑ j, (if i ≤ j then H i j x ^ 2 else 0)

/-- Theorem 5.5: Ω convex polygonal, a smooth with a ≥ a0 > 0, minimum-angle family
of triangulations; Ritz projection error bounds for s = 1, 2. -/
theorem ritz_projection_error (Ω : Set E2) (hΩ : IsPolygonalDomain Ω)
    (hconv : Convex ℝ Ω) (a : E2 → ℝ) (ha : ContDiffOn ℝ (⊤ : ℕ∞) a (closure Ω))
    (a0 : ℝ) (ha0 : 0 < a0) (hapos : ∀ x ∈ closure Ω, a0 ≤ a x)
    (ℋ : Set ℝ) (hℋ : ℋ ⊆ Set.Ioo 0 1) (T : ℝ → Finset (Fin 3 → E2))
    (hT : ∀ h ∈ ℋ, IsTriangulation Ω (T h) h)
    (hangle : ∃ θ0 : ℝ, 0 < θ0 ∧ ∀ h ∈ ℋ, ∀ K ∈ T h, MinAngleAtLeast K θ0) :
    ∃ C : ℝ, ∀ h ∈ ℋ, ∀ v w : E2 → ℝ, ∀ Gv Gw : Fin 2 → E2 → ℝ,
      InH10Closure Ω v → HasH1Gradient Ω v Gv →
      IsRitzProjection Ω a (T h) v w → HasH1Gradient Ω w Gw →
      (l2NormR Ω (fun x => w x - v x) ≤ C * h * Real.sqrt (h1NormSq Ω v Gv) ∧
        h1Seminorm Ω (fun i x => Gw i x - Gv i x) ≤ C * Real.sqrt (h1NormSq Ω v Gv)) ∧
      (∀ H : Fin 2 → Fin 2 → E2 → ℝ, (∀ i, HasH1Gradient Ω (Gv i) (H i)) →
        l2NormR Ω (fun x => w x - v x) ≤ C * h ^ 2 * Real.sqrt (h2NormSq Ω v Gv H) ∧
        h1Seminorm Ω (fun i x => Gw i x - Gv i x) ≤
          C * h * Real.sqrt (h2NormSq Ω v Gv H)) := by
  sorry

end PDETextbook.LarssonThomee.Thm55