Higher-order semidiscrete heat approximation

import PDETextbook.Base

open MeasureTheory
open scoped BigOperators ENNReal

noncomputable section

namespace PDETextbook.LarssonThomee.Thm104

open PDETextbook

/-- Energy form a(v,w) = ∫_Ω ∇v·∇w dx, evaluated on given weak gradients. -/
def energyForm {d : ℕ} (Ω : Set (Euclidean d)) (Gv Gw : Fin d → Euclidean d → ℝ) : ℝ :=
  ∫ x in Ω, ∑ i, Gv i x * Gw i x

/-- `D α` is the weak derivative D^α w on Ω for every multi-index |α| ≤ r, each in L²
(built iteratively from D 0 = w a.e.). -/
def IsWeakDerivFamily {d : ℕ} (Ω : Set (Euclidean d)) (r : ℕ) (w : Euclidean d → ℝ)
    (D : (Fin d → ℕ) → Euclidean d → ℝ) : Prop :=
  D 0 =ᵐ[domainMeasure Ω] w ∧
  (∀ α : Fin d → ℕ, ∑ i, α i ≤ r → InL2 Ω (D α)) ∧
  ∀ α : Fin d → ℕ, ∀ i : Fin d, ∑ j, α j < r →
    HasWeakPartial Ω (D α) i (D (α + Pi.single i 1))

/-- ‖w‖_r = (∑_{|α| ≤ r} ‖D^α w‖²)^{1/2}, computed from a weak-derivative family. -/
def HrNorm {d : ℕ} (Ω : Set (Euclidean d)) (r : ℕ)
    (D : (Fin d → ℕ) → Euclidean d → ℝ) : ℝ≥0∞ :=
  (∑ β ∈ (Finset.univ.filter fun β : Fin d → Fin (r + 1) => ∑ i, (β i : ℕ) ≤ r),
      eLpNorm (D fun i => (β i : ℕ)) 2 (domainMeasure Ω) ^ 2) ^ ((1 : ℝ) / 2)

/-- (5.49): Rw ∈ S and a(Rw - w, χ) = 0 for all χ ∈ S (Ritz projection). -/
def IsRitzProjection {d : ℕ} (Ω : Set (Euclidean d)) (S : Submodule ℝ (L2 ℝ Ω))
    (w Rw : L2 ℝ Ω) : Prop :=
  Rw ∈ S ∧ ∀ χ ∈ S, ∀ GR Gw Gχ : Fin d → Euclidean d → ℝ,
    HasH1Gradient Ω (⇑Rw) GR → HasH1Gradient Ω (⇑w) Gw → HasH1Gradient Ω (⇑χ) Gχ →
    energyForm Ω GR Gχ = energyForm Ω Gw Gχ

/-- Solution of (10.1) in the variational form (8.29)-(8.30): u(t) ∈ H₀¹, u_t = ut is the
L²-time derivative on [0,∞), (u_t,φ)+a(u,φ)=(f,φ) for φ ∈ H₀¹, t>0, and u(0)=v. -/
def IsHeatSolution {d : ℕ} (Ω : Set (Euclidean d)) (f : ℝ → L2 ℝ Ω) (v : L2 ℝ Ω)
    (u ut : ℝ → L2 ℝ Ω) : Prop :=
  (∀ t, 0 ≤ t → HasDerivWithinAt u (ut t) (Set.Ici 0) t) ∧
  (∀ t, 0 ≤ t → InH10Closure Ω (⇑(u t))) ∧
  (∀ t, 0 < t → ∀ φ : L2 ℝ Ω, InH10Closure Ω (⇑φ) →
    ∀ Gu Gφ : Fin d → Euclidean d → ℝ, HasH1Gradient Ω (⇑(u t)) Gu →
      HasH1Gradient Ω (⇑φ) Gφ →
      l2Inner (ut t) φ + energyForm Ω Gu Gφ = l2Inner (f t) φ) ∧
  u 0 = v

/-- Semidiscrete problem (10.3): u_h(t) ∈ S_h, (u_{h,t},χ)+a(u_h,χ)=(f,χ) ∀χ ∈ S_h, t>0,
u_h(0)=v_h. -/
def IsSemidiscreteSolution {d : ℕ} (Ω : Set (Euclidean d)) (S : Submodule ℝ (L2 ℝ Ω))
    (f : ℝ → L2 ℝ Ω) (vh : L2 ℝ Ω) (uh uht : ℝ → L2 ℝ Ω) : Prop :=
  (∀ t, 0 ≤ t → uh t ∈ S) ∧
  (∀ t, 0 ≤ t → HasDerivWithinAt uh (uht t) (Set.Ici 0) t) ∧
  (∀ t, 0 < t → ∀ χ ∈ S, ∀ Gu Gχ : Fin d → Euclidean d → ℝ,
    HasH1Gradient Ω (⇑(uh t)) Gu → HasH1Gradient Ω (⇑χ) Gχ →
      l2Inner (uht t) χ + energyForm Ω Gu Gχ = l2Inner (f t) χ) ∧
  uh 0 = vh

/-- Theorem 10.4: under (10.19), for suitably chosen v_h,
‖u_h(t) - u(t)‖ ≤ C h^r (‖v‖_r + ∫_0^t ‖u_t‖_r ds), t ≥ 0. -/
theorem thm_10_4 {d : ℕ} (Ω : Set (Euclidean d)) (hΩ : IsOpen Ω) (r : ℕ)
    (hs : Set ℝ) (hhs : hs ⊆ Set.Ioi 0) (S : ℝ → Submodule ℝ (L2 ℝ Ω))
    (hSfin : ∀ h ∈ hs, FiniteDimensional ℝ (S h))
    (hSH : ∀ h ∈ hs, ∀ χ ∈ S h, InH10Closure Ω (⇑χ))
    (h1019 : ∃ C₀ : ℝ, ∀ h ∈ hs, ∀ w : L2 ℝ Ω, InH10Closure Ω (⇑w) →
      ∀ D, IsWeakDerivFamily Ω r (⇑w) D → ∀ Rw, IsRitzProjection Ω (S h) w Rw →
        ENNReal.ofReal ‖Rw - w‖ ≤ ENNReal.ofReal (C₀ * h ^ r) * HrNorm Ω r D) :
    ∃ C : ℝ, 0 < C ∧ ∀ h ∈ hs, ∀ v : L2 ℝ Ω, ∃ vh ∈ S h,
      ∀ f u ut : ℝ → L2 ℝ Ω, IsHeatSolution Ω f v u ut →
      ∀ uh uht : ℝ → L2 ℝ Ω, IsSemidiscreteSolution Ω (S h) f vh uh uht →
      ∀ t : ℝ, 0 ≤ t →
      ∀ Dv, IsWeakDerivFamily Ω r (⇑v) Dv →
      ∀ Dut : ℝ → (Fin d → ℕ) → Euclidean d → ℝ,
        (∀ s ∈ Set.Icc 0 t, IsWeakDerivFamily Ω r (⇑(ut s)) (Dut s)) →
        ENNReal.ofReal ‖uh t - u t‖ ≤ ENNReal.ofReal (C * h ^ r) *
          (HrNorm Ω r Dv + ∫⁻ s in Set.Icc 0 t, HrNorm Ω r (Dut s)) := by
  sorry

end PDETextbook.LarssonThomee.Thm104