H1 error of semidiscrete heat finite elements

import PDETextbook.Base

open MeasureTheory
open scoped BigOperators ENNReal

noncomputable section

namespace PDETextbook.LarssonThomee.Thm10_2

abbrev E2 := Euclidean 2

/-- Ch. 10 setting: Ω ⊂ ℝ² a bounded convex domain with smooth (C^∞) boundary (level-set realization). -/
def SmoothBoundedConvexDomain (Ω : Set E2) : Prop :=
  Ω.Nonempty ∧ Bornology.IsBounded Ω ∧ Convex ℝ Ω ∧
  ∃ ρ : E2 → ℝ, ContDiff ℝ (⊤ : ℕ∞) ρ ∧ Ω = {x | ρ x < 0} ∧
    ∀ x, ρ x = 0 → fderiv ℝ ρ x ≠ 0

/-- A triangulation: finitely many vertex triples. -/
structure Triangulation where
  tri : Finset (Fin 3 → E2)

/-- The closed triangle K with vertices p. -/
def triK (p : Fin 3 → E2) : Set E2 := convexHull ℝ (Set.range p)

/-- Ω_h, the union of the triangles. -/
def Triangulation.region (T : Triangulation) : Set E2 := ⋃ p ∈ T.tri, triK p

/-- h = max_K h_K, h_K = diam K. -/
def Triangulation.meshsize (T : Triangulation) : ℝ := ⨆ p ∈ T.tri, Metric.diam (triK p)

/-- Sect. 5.2-type triangulation of Ω: nondegenerate, conforming (two triangles meet in a
common face), vertices in Ω̄, boundary vertices of Ω_h on Γ. -/
def IsTriangulationOf (Ω : Set E2) (T : Triangulation) : Prop :=
  (∀ p ∈ T.tri, AffineIndependent ℝ p) ∧
  (∀ p ∈ T.tri, ∀ q ∈ T.tri, triK p ∩ triK q = convexHull ℝ (Set.range p ∩ Set.range q)) ∧
  (∀ p ∈ T.tri, ∀ i, p i ∈ closure Ω) ∧
  (∀ p ∈ T.tri, ∀ i, p i ∈ frontier T.region → p i ∈ frontier Ω)

/-- All angles of all triangles are at least θ₀. -/
def MinAngleAtLeast (T : Triangulation) (θ₀ : ℝ) : Prop :=
  ∀ p ∈ T.tri, ∀ i j k : Fin 3, i ≠ j → j ≠ k → i ≠ k →
    θ₀ ≤ EuclideanGeometry.angle (p i) (p j) (p k)

/-- S_h: continuous, linear on each K, vanishing on Ω \ Ω_h and on ∂Ω. -/
def Sh (T : Triangulation) (Ω : Set E2) : Set (E2 → ℝ) :=
  {χ | Continuous χ ∧ (∀ p ∈ T.tri, ∃ A : E2 →ᵃ[ℝ] ℝ, Set.EqOn χ A (triK p)) ∧
    (∀ x, x ∉ T.region → χ x = 0) ∧ ∀ x ∈ frontier Ω, χ x = 0}

/-- a(v,w) = ∫_Ω ∇v·∇w. -/
def aForm (Ω : Set E2) (w χ : E2 → ℝ) : ℝ :=
  ∫ x in Ω, ∑ i, partialDeriv w i x * partialDeriv χ i x

def laplacian (w : E2 → ℝ) (x : E2) : ℝ :=
  ∑ i, partialDeriv (fun y => partialDeriv w i y) i x

/-- Classical solution of (10.1). -/
def IsHeatSolution (Ω : Set E2) (f : ℝ → E2 → ℝ) (v : E2 → ℝ) (u : ℝ → E2 → ℝ) : Prop :=
  ContinuousOn (fun p : ℝ × E2 => u p.1 p.2) (Set.Ici (0:ℝ) ×ˢ closure Ω) ∧
  (∀ t : ℝ, 0 < t → ContDiffOn ℝ 2 (u t) Ω) ∧
  (∀ t : ℝ, 0 < t → ∀ x ∈ Ω, HasDerivAt (fun τ => u τ x) (laplacian (u t) x + f t x) t) ∧
  (∀ t : ℝ, 0 < t → ∀ x ∈ frontier Ω, u t x = 0) ∧
  (∀ x ∈ Ω, u 0 x = v x)

/-- Solution of the semidiscrete problem (10.3). -/
def IsSemidiscreteSolution (Ω : Set E2) (T : Triangulation) (f : ℝ → E2 → ℝ)
    (vh : E2 → ℝ) (uh : ℝ → E2 → ℝ) : Prop :=
  (∀ t : ℝ, 0 ≤ t → uh t ∈ Sh T Ω) ∧
  (∀ x, ContinuousOn (fun τ => uh τ x) (Set.Ici 0)) ∧
  (∀ t : ℝ, 0 < t → ∀ x, DifferentiableAt ℝ (fun τ => uh τ x) t) ∧
  (∀ t : ℝ, 0 < t → ∀ χ ∈ Sh T Ω,
    (∫ x in Ω, deriv (fun τ => uh τ x) t * χ x) + aForm Ω (uh t) χ = ∫ x in Ω, f t x * χ x) ∧
  uh 0 = vh

/-- |w|_1 (⊤ if w ∉ H¹). -/
def h1SeminormE {d : ℕ} (Ω : Set (Euclidean d)) (w : Euclidean d → ℝ) : ℝ≥0∞ :=
  ⨅ (G : Fin d → Euclidean d → ℝ) (_ : HasH1Gradient Ω w G),
    (∑ i, eLpNorm (G i) 2 (domainMeasure Ω) ^ 2) ^ (1/2 : ℝ)

/-- ‖w‖_1 (⊤ if w ∉ H¹). -/
def h1NormE {d : ℕ} (Ω : Set (Euclidean d)) (w : Euclidean d → ℝ) : ℝ≥0∞ :=
  ⨅ (G : Fin d → Euclidean d → ℝ) (_ : HasH1Gradient Ω w G),
    (eLpNorm w 2 (domainMeasure Ω) ^ 2 + ∑ i, eLpNorm (G i) 2 (domainMeasure Ω) ^ 2) ^ (1/2 : ℝ)

/-- ‖w‖_2 = (Σ_{|α|≤2} ‖D^α w‖²)^{1/2}, each second-order multi-index once (i ≤ j). -/
def h2NormE {d : ℕ} (Ω : Set (Euclidean d)) (w : Euclidean d → ℝ) : ℝ≥0∞ :=
  ⨅ (G : Fin d → Euclidean d → ℝ) (H : Fin d → Fin d → Euclidean d → ℝ)
    (_ : HasH1Gradient Ω w G ∧ ∀ i, HasH1Gradient Ω (G i) (H i)),
    (eLpNorm w 2 (domainMeasure Ω) ^ 2 + ∑ i, eLpNorm (G i) 2 (domainMeasure Ω) ^ 2 +
      ∑ i, ∑ j, (if i ≤ j then eLpNorm (H i j) 2 (domainMeasure Ω) ^ 2 else 0)) ^ (1/2 : ℝ)

/-- Theorem 10.2 (H¹ error of the semidiscrete piecewise-linear method). -/
theorem semidiscrete_H1_error (Ω : Set E2) (hΩ : SmoothBoundedConvexDomain Ω)
    (𝒯 : Set Triangulation) (θ₀ c₀ : ℝ) (hθ₀ : 0 < θ₀)
    (hfam : ∀ T ∈ 𝒯, IsTriangulationOf Ω T ∧ 0 < T.meshsize ∧ T.meshsize < 1 ∧
      MinAngleAtLeast T θ₀ ∧ ∀ x ∈ Ω, Metric.infDist x T.region ≤ c₀ * T.meshsize ^ 2) :
    ∃ C : ℝ, 0 < C ∧ ∀ T ∈ 𝒯, ∀ (f : ℝ → E2 → ℝ) (v vh : E2 → ℝ) (u uh : ℝ → E2 → ℝ),
      IsHeatSolution Ω f v u → IsSemidiscreteSolution Ω T f vh uh → ∀ t : ℝ, 0 ≤ t →
      h1SeminormE Ω (fun x => uh t x - u t x) ≤
        h1SeminormE Ω (fun x => vh x - v x) + ENNReal.ofReal (C * T.meshsize) *
          (h2NormE Ω v + h2NormE Ω (u t) +
            (∫⁻ s in Set.Ioc 0 t, h1NormE Ω (fun x => deriv (fun τ => u τ x) s) ^ 2) ^ (1/2 : ℝ)) := by
  sorry

end PDETextbook.LarssonThomee.Thm10_2