Continuous representatives from Sobolev regularity

import PDETextbook.Base

open MeasureTheory
open scoped BigOperators ENNReal

noncomputable section

namespace PDETextbook.LarssonThomeeA5

/-- A bounded domain: an open, connected, bounded subset of `ℝ^d`. -/
def IsBoundedDomain {d : ℕ} (Ω : Set (Euclidean d)) : Prop :=
  IsOpen Ω ∧ IsConnected Ω ∧ Bornology.IsBounded Ω

/-- Smooth (C^∞) boundary, as a local level-set description. Near each boundary point `x`,
`Ω = {ρ < 0}` for a C^∞ function `ρ` with nonvanishing differential at `x`. -/
def HasSmoothBoundary {d : ℕ} (Ω : Set (Euclidean d)) : Prop :=
  ∀ x ∈ frontier Ω, ∃ U ∈ nhds x, ∃ ρ : Euclidean d → ℝ,
    ContDiff ℝ (⊤ : ℕ∞) ρ ∧ fderiv ℝ ρ x ≠ 0 ∧ Ω ∩ U = {y | ρ y < 0} ∩ U

/-- Lipschitz boundary: `Ω` lies locally on one side of its boundary, below a Lipschitz graph.
Near each boundary point `x` there is a unit direction `e` and a Lipschitz function `g`
which does not depend on the `e`-coordinate (`g (y + t • e) = g y`). Locally,
`Ω = {y | ⟪e, y⟫ < g y}`. This excludes slits and pinch points, e.g. two squares meeting
only at a corner. -/
def HasLipschitzBoundary {d : ℕ} (Ω : Set (Euclidean d)) : Prop :=
  ∀ x ∈ frontier Ω, ∃ e : Euclidean d, ‖e‖ = 1 ∧
    ∃ K : NNReal, ∃ g : Euclidean d → ℝ, LipschitzWith K g ∧
      (∀ (y : Euclidean d) (t : ℝ), g (y + t • e) = g y) ∧
      ∃ U ∈ nhds x, Ω ∩ U = {y | @inner ℝ _ _ e y < g y} ∩ U

/-- Polygonal (polyhedral, piecewise flat) boundary in the usual sense. `Ω` has a Lipschitz
boundary, so it lies on one side of its boundary with no slits or pinch points. The boundary
also lies in finitely many hyperplanes, i.e. it consists of finitely many flat faces. -/
def HasPolygonalBoundary {d : ℕ} (Ω : Set (Euclidean d)) : Prop :=
  HasLipschitzBoundary Ω ∧
  ∃ n : ℕ, ∃ a : Fin n → Euclidean d, ∃ b : Fin n → ℝ,
    (∀ j, a j ≠ 0) ∧ frontier Ω ⊆ ⋃ j, {x | @inner ℝ _ _ (a j) x = b j}

/-- `D` is a family of weak derivatives of `u` up to order `k`, all in `L²(Ω)`.
Here `D l` is the iterated weak derivative along the index list `l`, with `D [] = u`.
This is the base library's iterated weak-derivative characterization of `H^k(Ω)`. -/
def HasWeakDerivsUpTo {𝕜 : Type*} [RCLike 𝕜] {d : ℕ} (Ω : Set (Euclidean d)) (k : ℕ)
    (u : Euclidean d → 𝕜) (D : List (Fin d) → Euclidean d → 𝕜) : Prop :=
  D [] = u ∧
  (∀ l : List (Fin d), l.length ≤ k → InL2 Ω (D l)) ∧
  ∀ l : List (Fin d), l.length < k → ∀ i : Fin d, HasWeakPartial Ω (D l) i (D (i :: l))

/-- `v ∈ H^k(Ω)`. -/
def InHk {𝕜 : Type*} [RCLike 𝕜] {d : ℕ} (Ω : Set (Euclidean d)) (k : ℕ)
    (u : Euclidean d → 𝕜) : Prop :=
  ∃ D, HasWeakDerivsUpTo Ω k u D

/-- The canonical index list realizing the multi-index `α`, i.e. `D^α`. -/
def multiIndexList {d : ℕ} (α : Fin d → ℕ) : List (Fin d) :=
  (List.finRange d).flatMap (fun i => List.replicate (α i) i)

/-- `‖v‖_k = (∑_{|α| ≤ k} ‖D^α v‖²)^{1/2}`, computed from a weak-derivative family `D`.
Multi-indices with `|α| ≤ k` are encoded as `β : Fin d → Fin (k+1)` with `∑ β i ≤ k`. -/
def hkNorm {𝕜 : Type*} [RCLike 𝕜] {d : ℕ} (Ω : Set (Euclidean d)) (k : ℕ)
    (D : List (Fin d) → Euclidean d → 𝕜) : ℝ :=
  Real.sqrt (∑ β ∈ (Finset.univ.filter fun β : Fin d → Fin (k + 1) => ∑ i, (β i : ℕ) ≤ k),
    ∫ x in Ω, ‖D (multiIndexList (fun i => (β i : ℕ))) x‖ ^ 2)

/-- Theorem A.5 (Sobolev embedding). Let `Ω` be a bounded domain in `ℝ^d` with smooth or
polygonal boundary, and let `k > d/2`. Then `H^k(Ω) ⊂ C(Ω̄)`: every `v ∈ H^k` has a
representative that is continuous on `Ω̄`. Moreover there is a constant `C = C(Ω)` with
`‖v‖_C = max_{Ω̄} |v| ≤ C ‖v‖_k`. -/
theorem sobolev_embedding_continuous {𝕜 : Type*} [RCLike 𝕜] {d : ℕ}
    (Ω : Set (Euclidean d)) (hΩ : IsBoundedDomain Ω)
    (hbdry : HasSmoothBoundary Ω ∨ HasPolygonalBoundary Ω)
    (k : ℕ) (hk : (d : ℝ) / 2 < k) :
    ∃ C : ℝ, ∀ (v : Euclidean d → 𝕜) (D : List (Fin d) → Euclidean d → 𝕜),
      HasWeakDerivsUpTo Ω k v D →
      ∃ w : Euclidean d → 𝕜, ContinuousOn w (closure Ω) ∧ w =ᵐ[domainMeasure Ω] v ∧
        ∀ x ∈ closure Ω, ‖w x‖ ≤ C * hkNorm Ω k D := by
  sorry

end PDETextbook.LarssonThomeeA5