Forward Euler Method

We implement the explicit Euler method for ODEs and prove its convergence.

Generic infrastructure

• piecewiseLinear, piecewiseConst: Piecewise linear/constant interpolation on a regular grid.
• locallyFinite_Icc_grid: The regular grid is locally finite.
• ContinuousOn.of_Icc_grid: Cell-wise continuity implies continuity on [a, ∞).

Euler method

• ODE.EulerMethod.step, ODE.EulerMethod.point, ODE.EulerMethod.slope: The Euler iteration.
• ODE.EulerMethod.path, ODE.EulerMethod.deriv: Piecewise linear/constant interpolation of the Euler points.
• ODE.EulerMethod.dist_deriv_le: Global bound on the local truncation error.
• ODE.EulerMethod.dist_path_le: Error bound via Gronwall's inequality.
• ODE.EulerMethod.tendsto_path: Convergence as h → 0⁺.

Grid helpers

theorem Nat.floor_div_eq_of_mem_Ico {α : Type u_1} [Field α] [LinearOrder α] [FloorSemiring α] [IsStrictOrderedRing α] {h : α} (hh : 0 < h) {a : α} {n : ℕ} {t : α} (ht : t ∈ Set.Ico (a + ↑n * h) (a + (↑n + 1) * h)) : ⌊(t - a) / h⌋₊ = n

If t ∈ [a + n * h, a + (n + 1) * h) and 0 < h, then ⌊(t - a) / h⌋₊ = n.

theorem mem_Ico_Nat_floor_div {α : Type u_1} [Field α] [LinearOrder α] [FloorSemiring α] [IsStrictOrderedRing α] {h : α} (hh : 0 < h) {a t : α} (hat : a ≤ t) : t ∈ Set.Ico (a + ↑⌊(t - a) / h⌋₊ * h) (a + (↑⌊(t - a) / h⌋₊ + 1) * h)

If 0 < h and a ≤ t, then t lies in the floor interval [a + ⌊(t - a) / h⌋₊ * h, a + (⌊(t - a) / h⌋₊ + 1) * h).

Piecewise linear interpolation

noncomputable def piecewiseLinear {α : Type u_1} [Field α] [LinearOrder α] [FloorSemiring α] {E : Type u_2} [AddCommGroup E] [Module α E] (y c : ℕ → E) (h a t : α) : E

The piecewise linear interpolation of a sequence y with slopes c on a regular grid with step size h starting at a. On [a + n * h, a + (n + 1) * h), the value is y n + (t - (a + n * h)) • c n.

Equations

• piecewiseLinear y c h a t = y ⌊(t - a) / h⌋₊ + (t - (a + ↑⌊(t - a) / h⌋₊ * h)) • c ⌊(t - a) / h⌋₊

noncomputable def piecewiseConst {α : Type u_1} [Field α] [LinearOrder α] [FloorSemiring α] {E : Type u_2} (c : ℕ → E) (h a t : α) : E

The piecewise constant function taking value c n on [a + n * h, a + (n + 1) * h).

Equations

• piecewiseConst c h a t = c ⌊(t - a) / h⌋₊

theorem piecewiseConst_eq_on_Ico {α : Type u_1} [Field α] [LinearOrder α] [FloorSemiring α] [IsStrictOrderedRing α] {E : Type u_2} {c : ℕ → E} {h a : α} (hh : 0 < h) {n : ℕ} {t : α} (ht : t ∈ Set.Ico (a + ↑n * h) (a + (↑n + 1) * h)) : piecewiseConst c h a t = c n

The piecewise constant function equals c n on [a + n * h, a + (n + 1) * h).

theorem locallyFinite_Icc_grid {α : Type u_1} [Field α] [LinearOrder α] [FloorSemiring α] [IsStrictOrderedRing α] [TopologicalSpace α] [OrderTopology α] {h : α} (hh : 0 < h) (a : α) : LocallyFinite fun (n : ℕ) => Set.Icc (a + ↑n * h) (a + (↑n + 1) * h)

The regular grid of closed intervals [a + n * h, a + (n + 1) * h] is locally finite.

theorem ContinuousOn.of_Icc_grid {α : Type u_1} [Field α] [LinearOrder α] [FloorSemiring α] [IsStrictOrderedRing α] [TopologicalSpace α] [OrderTopology α] {F : Type u_2} [TopologicalSpace F] {f : α → F} {h : α} (hh : 0 < h) {a : α} (hf : ∀ (n : ℕ), ContinuousOn f (Set.Icc (a + ↑n * h) (a + (↑n + 1) * h))) : ContinuousOn f (Set.Ici a)

A function continuous on each cell [a + n * h, a + (n + 1) * h] is continuous on [a, ∞).

theorem piecewiseLinear_apply_grid {E : Type u_2} [NormedAddCommGroup E] [NormedSpace ℝ E] {y c : ℕ → E} {h : ℝ} (hh : 0 < h) (a : ℝ) (n : ℕ) : piecewiseLinear y c h a (a + ↑n * h) = y n

The piecewise linear interpolation at a grid point a + n * h equals y n.

theorem piecewiseLinear_eq_on_Ico {E : Type u_2} [NormedAddCommGroup E] [NormedSpace ℝ E] {y c : ℕ → E} {h a : ℝ} (hh : 0 < h) {n : ℕ} {t : ℝ} (ht : t ∈ Set.Ico (a + ↑n * h) (a + (↑n + 1) * h)) : piecewiseLinear y c h a t = y n + (t - (a + ↑n * h)) • c n

The piecewise linear interpolation equals y n + (t - (a + n * h)) • c n on [a + n * h, a + (n + 1) * h).

theorem piecewiseLinear_continuousOn {E : Type u_2} [NormedAddCommGroup E] [NormedSpace ℝ E] {y c : ℕ → E} {h a : ℝ} (hh : 0 < h) (hstep : ∀ (n : ℕ), y (n + 1) = y n + h • c n) : ContinuousOn (piecewiseLinear y c h a) (Set.Ici a)

A piecewise linear function whose grid values satisfy y (n + 1) = y n + h • c n is continuous on [a, ∞).

theorem piecewiseLinear_hasDerivWithinAt {E : Type u_2} [NormedAddCommGroup E] [NormedSpace ℝ E] {y c : ℕ → E} {h a : ℝ} (hh : 0 < h) {t : ℝ} (hat : a ≤ t) : HasDerivWithinAt (piecewiseLinear y c h a) (piecewiseConst c h a t) (Set.Ici t) t

The right derivative of a piecewise linear function is the piecewise constant slope.

Euler method

def ODE.EulerMethod.step {𝕜 : Type u_3} {E : Type u_4} [Ring 𝕜] [AddCommGroup E] [Module 𝕜 E] (v : 𝕜 → E → E) (h t : 𝕜) (y : E) : E

A single step of the explicit Euler method: y + h • v(t, y).

Equations

• ODE.EulerMethod.step v h t y = y + h • v t y

def ODE.EulerMethod.point {𝕜 : Type u_3} {E : Type u_4} [Ring 𝕜] [AddCommGroup E] [Module 𝕜 E] (v : 𝕜 → E → E) (h t₀ : 𝕜) (y₀ : E) : ℕ → E

The sequence of Euler points, defined recursively: point v h t₀ y₀ 0 = y₀ and point v h t₀ y₀ (n+1) = step v h (t₀ + n*h) (point v h t₀ y₀ n).

Equations

• ODE.EulerMethod.point v h t₀ y₀ 0 = y₀
• ODE.EulerMethod.point v h t₀ y₀ n.succ = ODE.EulerMethod.step v h (t₀ + ↑n * h) (ODE.EulerMethod.point v h t₀ y₀ n)

noncomputable def ODE.EulerMethod.slope {E : Type u_2} [NormedAddCommGroup E] [NormedSpace ℝ E] (v : ℝ → E → E) (h t₀ : ℝ) (y₀ : E) (n : ℕ) : E

The slope of the Euler method on the n-th cell: v(t₀ + n * h, yₙ).

Equations

• ODE.EulerMethod.slope v h t₀ y₀ n = v (t₀ + ↑n * h) (ODE.EulerMethod.point v h t₀ y₀ n)

noncomputable def ODE.EulerMethod.path {E : Type u_2} [NormedAddCommGroup E] [NormedSpace ℝ E] (v : ℝ → E → E) (h t₀ : ℝ) (y₀ : E) : ℝ → E

The piecewise linear Euler path, interpolating the Euler points with Euler slopes.

Equations

• ODE.EulerMethod.path v h t₀ y₀ = piecewiseLinear (ODE.EulerMethod.point v h t₀ y₀) (ODE.EulerMethod.slope v h t₀ y₀) h t₀

noncomputable def ODE.EulerMethod.deriv {E : Type u_2} [NormedAddCommGroup E] [NormedSpace ℝ E] (v : ℝ → E → E) (h t₀ : ℝ) (y₀ : E) : ℝ → E

The piecewise constant right derivative of the Euler path.

Equations

• ODE.EulerMethod.deriv v h t₀ y₀ = piecewiseConst (ODE.EulerMethod.slope v h t₀ y₀) h t₀

theorem ODE.EulerMethod.dist_deriv_le {E : Type u_2} [NormedAddCommGroup E] [NormedSpace ℝ E] {h : ℝ} {v : ℝ → E → E} {K L : NNReal} {M : ℝ} (hv : ∀ (t : ℝ), LipschitzWith K (v t)) (hvt : ∀ (y : E), LipschitzWith L fun (t : ℝ) => v t y) (hM : ∀ (t : ℝ) (y : E), ‖v t y‖ ≤ M) {t₀ : ℝ} {y₀ : E} (hh : 0 < h) {t : ℝ} (ht₀ : t₀ ≤ t) : dist (deriv v h t₀ y₀ t) (v t (path v h t₀ y₀ t)) ≤ h * (↑L + ↑K * M)

Global bound on the difference between the Euler derivative and the vector field along the Euler path.

theorem ODE.EulerMethod.dist_path_le {E : Type u_2} [NormedAddCommGroup E] [NormedSpace ℝ E] {h : ℝ} {v : ℝ → E → E} {K L : NNReal} {M : ℝ} (hv : ∀ (t : ℝ), LipschitzWith K (v t)) (hvt : ∀ (y : E), LipschitzWith L fun (t : ℝ) => v t y) (hM : ∀ (t : ℝ) (y : E), ‖v t y‖ ≤ M) {t₀ : ℝ} {y₀ : E} (hh : 0 < h) {T : ℝ} {sol : ℝ → E} (hsol : ContinuousOn sol (Set.Icc t₀ T)) (hsol' : ∀ t ∈ Set.Ico t₀ T, HasDerivWithinAt sol (v t (sol t)) (Set.Ici t) t) (hsol₀ : sol t₀ = y₀) (t : ℝ) : t ∈ Set.Icc t₀ T → dist (path v h t₀ y₀ t) (sol t) ≤ gronwallBound 0 (↑K) (h * (↑L + ↑K * M)) (t - t₀)

Error bound for the Euler method via Gronwall's inequality.

theorem ODE.EulerMethod.tendsto_path {E : Type u_2} [NormedAddCommGroup E] [NormedSpace ℝ E] {v : ℝ → E → E} {K L : NNReal} {M : ℝ} (hv : ∀ (t : ℝ), LipschitzWith K (v t)) (hvt : ∀ (y : E), LipschitzWith L fun (t : ℝ) => v t y) (hM : ∀ (t : ℝ) (y : E), ‖v t y‖ ≤ M) {t₀ : ℝ} {y₀ : E} {T : ℝ} {sol : ℝ → E} (hsol : ContinuousOn sol (Set.Icc t₀ T)) (hsol' : ∀ t ∈ Set.Ico t₀ T, HasDerivWithinAt sol (v t (sol t)) (Set.Ici t) t) (hsol₀ : sol t₀ = y₀) (t : ℝ) : t ∈ Set.Icc t₀ T → Filter.Tendsto (fun (δ : ℝ) => path v δ t₀ y₀ t) (nhdsWithin 0 (Set.Ioi 0)) (nhds (sol t))

The Euler method converges to the true solution as h → 0⁺.