MaxLeafFormal.DetRate

Source file MaxLeafFormal/DetRate.lean of the Lean project (how it is built).

Statements of the paper formalised here: Lemma 10.1.

import Mathlib

/-!
# The deterministic rate certificate

`lem:det-cert` of the paper:
for every `u ∈ [0.15, 0.85]` there is `y ∈ {(1+√5)/8, 1/2, 5/16, 1/4, 27/128, 29/160}` with
`χ_y(u) ≤ log(√5 - 1/2)`. The proof follows the paper: `χ_y` is affine on each `[1/(a+1), 1/a]`; the
interval `a = 2` is algebraic, and the ten other endpoint inequalities are exact rational
inequalities with natural exponents.
-/

namespace MaxLeafFormal.DetRate

/-- `χ_y(u)` with `a = ⌊1/u⌋`. -/
noncomputable def chi (y u : ℝ) : ℝ :=
  let a : ℕ := ⌊1 / u⌋₊
  Real.logb 2 (1 + y) + 1 - ((a + 1) * u - 1) * Real.logb 2 (1 + y * 2 ^ a)
    - (1 - a * u) * Real.logb 2 (1 + y * 2 ^ (a + 1))

/-- The same formula with `a` fixed. -/
noncomputable def chiA (a : ℕ) (y u : ℝ) : ℝ :=
  Real.logb 2 (1 + y) + 1 - ((a + 1) * u - 1) * Real.logb 2 (1 + y * 2 ^ a)
    - (1 - a * u) * Real.logb 2 (1 + y * 2 ^ (a + 1))

/-! ### Reduction to endpoints (B.1) -/

theorem chi_eq_chiA_of_floor {y u : ℝ} {a : ℕ} (h : ⌊1 / u⌋₊ = a) : chi y u = chiA a y u := by
  unfold chi chiA
  simp only [h]

/-- At the left endpoint `u = 1/(a+1)` both formulas agree. -/
theorem chiA_succ_left (a : ℕ) (y : ℝ) :
    chiA (a + 1) y (1 / (a + 1)) = chiA a y (1 / (a + 1)) := by
  unfold chiA
  push_cast
  have hm : ((a:ℝ) + 1) * (1 / ((a:ℝ) + 1)) = 1 := mul_one_div_cancel (by positivity)
  linear_combination (Real.logb 2 (1 + y * 2 ^ a) + Real.logb 2 (1 + y * 2 ^ (a + 1 + 1))
    - 2 * Real.logb 2 (1 + y * 2 ^ (a + 1))) * hm

/-- On `[1/(a+1), 1/a]`, with `a ≥ 1`, `χ_y = χ_y^{(a)}`. -/
theorem chi_eq_chiA {a : ℕ} (ha : 1 ≤ a) {u : ℝ} (hu : u ∈ Set.Icc (1 / (a + 1 : ℝ)) (1 / a))
    (y : ℝ) : chi y u = chiA a y u := by
  obtain ⟨hu1, hu2⟩ := hu
  have hapos : (0:ℝ) < a := Nat.cast_pos.mpr (by omega)
  have hupos : 0 < u := lt_of_lt_of_le (by positivity) hu1
  rcases eq_or_lt_of_le hu1 with heq | hlt
  · have hfl : ⌊1 / u⌋₊ = a + 1 := by
      rw [← heq, one_div_one_div]
      rw [show ((a:ℝ) + 1) = ((a + 1 : ℕ) : ℝ) by push_cast; ring, Nat.floor_natCast]
    rw [chi_eq_chiA_of_floor hfl, ← heq]
    exact chiA_succ_left a y
  · have hfl : ⌊1 / u⌋₊ = a := by
      rw [Nat.floor_eq_iff (one_div_pos.mpr hupos).le]
      constructor
      · rw [le_div_iff₀ hupos]
        have := (le_div_iff₀ hapos).mp hu2
        linarith
      · rw [div_lt_iff₀ hupos]
        have := (div_lt_iff₀ (by positivity : (0:ℝ) < a + 1)).mp hlt
        linarith
    exact chi_eq_chiA_of_floor hfl

lemma aff_le {C D l r u c : ℝ} (hl : C + D * l ≤ c) (hr : C + D * r ≤ c) (hu1 : l ≤ u)
    (hu2 : u ≤ r) : C + D * u ≤ c := by
  rcases le_total 0 D with hD | hD
  · nlinarith [mul_le_mul_of_nonneg_left hu2 hD]
  · nlinarith [mul_le_mul_of_nonpos_left hu1 hD]

/-- `χ^{(a)}_y` is affine in `u`, so it is bounded by its endpoint values. -/
theorem chiA_le_of_endpoints {a : ℕ} {y l r c : ℝ} (hl : chiA a y l ≤ c) (hr : chiA a y r ≤ c)
    {u : ℝ} (hu : u ∈ Set.Icc l r) : chiA a y u ≤ c := by
  have hlin : ∀ v, chiA a y v = (Real.logb 2 (1 + y) + 1 + Real.logb 2 (1 + y * 2 ^ a)
      - Real.logb 2 (1 + y * 2 ^ (a + 1)))
      + (-((a:ℝ) + 1) * Real.logb 2 (1 + y * 2 ^ a) + (a:ℝ) * Real.logb 2 (1 + y * 2 ^ (a + 1))) * v := by
    intro v
    unfold chiA
    ring
  rw [hlin l] at hl
  rw [hlin r] at hr
  rw [hlin u]
  exact aff_le hl hr hu.1 hu.2

/-! ### Endpoint inequalities (B.2) -/

/-- The exponentiated form of an endpoint inequality. -/
theorem logb_endpoint {A X Y B : ℝ} (hA : 0 < A) (hX : 0 < X) (hY : 0 < Y) (hB : 0 < B)
    {h p₁ p₂ : ℕ} (hh : 0 < h) (hineq : A ^ h ≤ B ^ h * X ^ p₁ * Y ^ p₂) :
    Real.logb 2 A - (p₁ / h : ℝ) * Real.logb 2 X - (p₂ / h : ℝ) * Real.logb 2 Y ≤ Real.logb 2 B := by
  have hpos : (0:ℝ) < h := by exact_mod_cast hh
  have hh0 : (h:ℝ) ≠ 0 := hpos.ne'
  have h1 : Real.logb 2 (A ^ h) ≤ Real.logb 2 (B ^ h * X ^ p₁ * Y ^ p₂) :=
    (Real.logb_le_logb (by norm_num) (pow_pos hA h)
      (mul_pos (mul_pos (pow_pos hB h) (pow_pos hX p₁)) (pow_pos hY p₂))).2 hineq
  rw [Real.logb_mul (mul_pos (pow_pos hB h) (pow_pos hX p₁)).ne' (pow_pos hY p₂).ne',
    Real.logb_mul (pow_pos hB h).ne' (pow_pos hX p₁).ne'] at h1
  simp only [Real.logb_pow] at h1
  have ex1 : (h:ℝ) * ((p₁:ℝ) / h) = p₁ := by field_simp
  have ex2 : (h:ℝ) * ((p₂:ℝ) / h) = p₂ := by field_simp
  have key : (h:ℝ) * (Real.logb 2 A - (p₁ / h : ℝ) * Real.logb 2 X - (p₂ / h : ℝ) * Real.logb 2 Y)
      ≤ (h:ℝ) * Real.logb 2 B := by
    have e : (h:ℝ) * (Real.logb 2 A - (p₁ / h : ℝ) * Real.logb 2 X - (p₂ / h : ℝ) * Real.logb 2 Y)
        = h * Real.logb 2 A - ((h:ℝ) * ((p₁:ℝ) / h)) * Real.logb 2 X
          - ((h:ℝ) * ((p₂:ℝ) / h)) * Real.logb 2 Y := by ring
    rw [e, ex1, ex2]
    linarith [h1]
  exact le_of_mul_le_mul_left key hpos

theorem chiA_eq_logb_two_mul (a : ℕ) {y : ℝ} (hy : 0 ≤ y) (u : ℝ) :
    chiA a y u = Real.logb 2 (2 * (1 + y)) - ((a + 1) * u - 1) * Real.logb 2 (1 + y * 2 ^ a)
      - (1 - a * u) * Real.logb 2 (1 + y * 2 ^ (a + 1)) := by
  unfold chiA
  rw [Real.logb_mul (by norm_num) (by linarith : (0:ℝ) < 1 + y).ne',
    Real.logb_self_eq_one (by norm_num : (1:ℝ) < 2)]
  ring

/-- Endpoint inequality from the exact rational inequality. -/
theorem ep_gen {a : ℕ} {y u : ℝ} (hy : 0 ≤ y) {h p₁ p₂ : ℕ} (hh : 0 < h)
    (hc1 : ((a:ℝ) + 1) * u - 1 = (p₁ : ℝ) / h) (hc2 : 1 - (a:ℝ) * u = (p₂ : ℝ) / h)
    (hineq : (2 * (1 + y)) ^ h ≤ (1736/1000 : ℝ) ^ h * (1 + y * 2 ^ a) ^ p₁
      * (1 + y * 2 ^ (a + 1)) ^ p₂) :
    chiA a y u ≤ Real.logb 2 (1736/1000) := by
  rw [chiA_eq_logb_two_mul a hy u, hc1, hc2]
  have hX : (0:ℝ) < 1 + y * 2 ^ a := by
    have := mul_nonneg hy (by positivity : (0:ℝ) ≤ 2 ^ a)
    linarith
  have hY : (0:ℝ) < 1 + y * 2 ^ (a + 1) := by
    have := mul_nonneg hy (by positivity : (0:ℝ) ≤ 2 ^ (a + 1))
    linarith
  exact logb_endpoint (by linarith) hX hY (by norm_num) hh hineq

-- a = 1, y = 1/2, u ∈ [1/2, 17/20]

theorem ep_1_left : chiA 1 (1/2) (1/2) ≤ Real.logb 2 (1736/1000) :=
  ep_gen (by norm_num) (h := 2) (p₁ := 0) (p₂ := 1) (by norm_num) (by norm_num) (by norm_num)
    (by norm_num)

theorem ep_1_right : chiA 1 (1/2) (17/20) ≤ Real.logb 2 (1736/1000) :=
  ep_gen (by norm_num) (h := 20) (p₁ := 14) (p₂ := 3) (by norm_num) (by norm_num) (by norm_num)
    (by norm_num)
-- a = 3, y = 5/16, u ∈ [1/4, 1/3]

theorem ep_3_left : chiA 3 (5/16) (1/4) ≤ Real.logb 2 (1736/1000) :=
  ep_gen (by norm_num) (h := 4) (p₁ := 0) (p₂ := 1) (by norm_num) (by norm_num) (by norm_num)
    (by norm_num)

theorem ep_3_right : chiA 3 (5/16) (1/3) ≤ Real.logb 2 (1736/1000) :=
  ep_gen (by norm_num) (h := 3) (p₁ := 1) (p₂ := 0) (by norm_num) (by norm_num) (by norm_num)
    (by norm_num)
-- a = 4, y = 1/4, u ∈ [1/5, 1/4]

theorem ep_4_left : chiA 4 (1/4) (1/5) ≤ Real.logb 2 (1736/1000) :=
  ep_gen (by norm_num) (h := 5) (p₁ := 0) (p₂ := 1) (by norm_num) (by norm_num) (by norm_num)
    (by norm_num)

theorem ep_4_right : chiA 4 (1/4) (1/4) ≤ Real.logb 2 (1736/1000) :=
  ep_gen (by norm_num) (h := 4) (p₁ := 1) (p₂ := 0) (by norm_num) (by norm_num) (by norm_num)
    (by norm_num)
-- a = 5, y = 27/128, u ∈ [1/6, 1/5]

theorem ep_5_left : chiA 5 (27/128) (1/6) ≤ Real.logb 2 (1736/1000) :=
  ep_gen (by norm_num) (h := 6) (p₁ := 0) (p₂ := 1) (by norm_num) (by norm_num) (by norm_num)
    (by norm_num)

theorem ep_5_right : chiA 5 (27/128) (1/5) ≤ Real.logb 2 (1736/1000) :=
  ep_gen (by norm_num) (h := 5) (p₁ := 1) (p₂ := 0) (by norm_num) (by norm_num) (by norm_num)
    (by norm_num)
-- a = 6, y = 29/160, u ∈ [3/20, 1/6]

theorem ep_6_left : chiA 6 (29/160) (3/20) ≤ Real.logb 2 (1736/1000) :=
  ep_gen (by norm_num) (h := 20) (p₁ := 1) (p₂ := 2) (by norm_num) (by norm_num) (by norm_num)
    (by norm_num)

theorem ep_6_right : chiA 6 (29/160) (1/6) ≤ Real.logb 2 (1736/1000) :=
  ep_gen (by norm_num) (h := 6) (p₁ := 1) (p₂ := 0) (by norm_num) (by norm_num) (by norm_num)
    (by norm_num)

theorem bound_lt : (1736/1000 : ℝ) < Real.sqrt 5 - 1/2 := by
  have h1 : (2236/1000 : ℝ) < Real.sqrt 5 := (Real.lt_sqrt (by norm_num)).2 (by norm_num)
  linarith

theorem logb_mono_bound : Real.logb 2 (1736/1000) ≤ Real.logb 2 (Real.sqrt 5 - 1/2) :=
  (Real.logb_le_logb (by norm_num) (by norm_num) (by linarith [bound_lt])).2 bound_lt.le

/-! ### The interval `a = 2` (B.3) -/

/-- `χ^{(2)}_{y₀}(u) = log(√5 - 1/2)` for every `u`: it is constant in `u`. -/
theorem chiA_two_y0 (u : ℝ) :
    chiA 2 ((1 + Real.sqrt 5) / 8) u = Real.logb 2 (Real.sqrt 5 - 1/2) := by
  have hs : Real.sqrt 5 ^ 2 = 5 := Real.sq_sqrt (by norm_num)
  have hgt : (2236/1000 : ℝ) < Real.sqrt 5 := (Real.lt_sqrt (by norm_num)).2 (by norm_num)
  have hA : Real.sqrt 5 - 1/2 ≠ 0 := by
    intro h
    linarith
  have hB : (1 + Real.sqrt 5) / 2 ≠ 0 := by positivity
  rw [chiA_eq_logb_two_mul 2 (by positivity)]
  have e1 : (1:ℝ) + (1 + Real.sqrt 5) / 8 * 2 ^ 2 = ((1 + Real.sqrt 5) / 2) ^ 2 := by
    linear_combination (-1/4 : ℝ) * hs
  have e2 : (1:ℝ) + (1 + Real.sqrt 5) / 8 * 2 ^ (2 + 1) = ((1 + Real.sqrt 5) / 2) ^ 3 := by
    linear_combination ((-Real.sqrt 5 - 3) / 8) * hs
  have e3 : (2:ℝ) * (1 + (1 + Real.sqrt 5) / 8) = (Real.sqrt 5 - 1/2) * ((1 + Real.sqrt 5) / 2) := by
    linear_combination (-1/2 : ℝ) * hs
  rw [e1, e2, e3, Real.logb_mul hA hB]
  simp only [Real.logb_pow]
  push_cast
  ring

/-! ### The lemma (B.4) -/

/-- lem:det-cert. -/
theorem det_cert : ∀ u ∈ Set.Icc (15/100 : ℝ) (85/100),
    ∃ y ∈ ({(1 + Real.sqrt 5) / 8, 1/2, 5/16, 1/4, 27/128, 29/160} : Set ℝ),
      chi y u ≤ Real.logb 2 (Real.sqrt 5 - 1/2) := by
  intro u hu
  obtain ⟨hu1, hu2⟩ := hu
  rcases le_total u (1/6) with h1 | h1
  · refine ⟨29/160, by simp, ?_⟩
    rw [chi_eq_chiA (a := 6) (u := u) (by norm_num)
      ⟨by norm_num; linarith, by norm_num; linarith⟩ (29/160)]
    exact le_trans (chiA_le_of_endpoints (u := u) ep_6_left ep_6_right ⟨by linarith, h1⟩)
      logb_mono_bound
  · rcases le_total u (1/5) with h2 | h2
    · refine ⟨27/128, by simp, ?_⟩
      rw [chi_eq_chiA (a := 5) (u := u) (by norm_num)
        ⟨by norm_num; linarith, by norm_num; linarith⟩ (27/128)]
      exact le_trans (chiA_le_of_endpoints (u := u) ep_5_left ep_5_right ⟨h1, h2⟩)
        logb_mono_bound
    · rcases le_total u (1/4) with h3 | h3
      · refine ⟨1/4, by simp, ?_⟩
        rw [chi_eq_chiA (a := 4) (u := u) (by norm_num)
          ⟨by norm_num; linarith, by norm_num; linarith⟩ (1/4)]
        exact le_trans (chiA_le_of_endpoints (u := u) ep_4_left ep_4_right ⟨h2, h3⟩)
          logb_mono_bound
      · rcases le_total u (1/3) with h4 | h4
        · refine ⟨5/16, by simp, ?_⟩
          rw [chi_eq_chiA (a := 3) (u := u) (by norm_num)
            ⟨by norm_num; linarith, by norm_num; linarith⟩ (5/16)]
          exact le_trans (chiA_le_of_endpoints (u := u) ep_3_left ep_3_right ⟨h3, h4⟩)
            logb_mono_bound
        · rcases le_total u (1/2) with h5 | h5
          · refine ⟨(1 + Real.sqrt 5) / 8, by simp, ?_⟩
            rw [chi_eq_chiA (a := 2) (u := u) (by norm_num)
              ⟨by norm_num; linarith, by norm_num; linarith⟩ ((1 + Real.sqrt 5) / 8)]
            exact (chiA_two_y0 u).le
          · refine ⟨1/2, by simp, ?_⟩
            rw [chi_eq_chiA (a := 1) (u := u) (by norm_num)
              ⟨by norm_num; linarith, by norm_num; linarith⟩ (1/2)]
            exact le_trans (chiA_le_of_endpoints (u := u) ep_1_left ep_1_right
              ⟨h5, by linarith⟩) logb_mono_bound

end MaxLeafFormal.DetRate