MaxLeafFormal.RSThetaRate

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

Statements of the paper formalised here: Proposition 7.7, Theorem 7.2.

import MaxLeafFormal.EncloseH
import MaxLeafFormal.RSThetaCert

/-!
# RS-θ: the rates of the two certificates

`prop:rst-certs` of the paper, the
value `0.7807 < ω*` of `tab:rst-coverage` and of the proof of `prop:rst-rate`, `2^{ω*} = 1.7203… <
1.7204` (proof of `thm:rstheta`), and the program's (R4) `H(0.2) < 0.7802`. The certificates are those
of `MaxLeafFormal.RSThetaCert`.
-/

namespace MaxLeafFormal.RSThetaRate

open RSThetaCert Enclose

/-- `ω_cert(σ) = log(1+y) - σ log e₁ - (1-σ) log e₂`. -/
noncomputable def omegaCert (y e₁ e₂ σ : ℝ) : ℝ :=
  Real.logb 2 (1 + y) - σ * Real.logb 2 e₁ - (1 - σ) * Real.logb 2 e₂

noncomputable def ω₀ (σ : ℝ) : ℝ := omegaCert y0 e1_0 e2_0 σ

noncomputable def ω₁ (σ : ℝ) : ℝ := omegaCert y1 e1_1 e2_1 σ

/-- The constant coefficient `a₀ = log(1+y) - log e₂` of `ω_cert`. -/
noncomputable def a0 (y e₂ : ℝ) : ℝ := Real.logb 2 (1 + y) - Real.logb 2 e₂

/-- The slope `a₁ = log e₂ - log e₁` of `ω_cert`. -/
noncomputable def a1 (e₁ e₂ : ℝ) : ℝ := Real.logb 2 e₂ - Real.logb 2 e₁

theorem omegaCert_eq (y e₁ e₂ σ : ℝ) : omegaCert y e₁ e₂ σ = a0 y e₂ + a1 e₁ e₂ * σ := by
  unfold omegaCert a0 a1
  ring

/-! ### Enclosures of the coefficients -/

/-- `a₀` enclosed from rational data. -/
def a0Q (y e₂ : ℚ) : QI := ((log2Q (1 + y) 60 80).sub (log2Q e₂ 60 80))

/-- `a₁` enclosed from rational data. -/
def a1Q (e₁ e₂ : ℚ) : QI := ((log2Q e₂ 60 80).sub (log2Q e₁ 60 80))

theorem a0Q_sound {y e₂ : ℚ} (hy : 0 ≤ y) (he₂ : 0 < e₂) : (a0Q y e₂).Contains (a0 y e₂) := by
  have h1 : (0:ℚ) < 1 + y := by linarith
  have h := QI.sub_sound (log2Q_sound h1 60 80) (log2Q_sound he₂ 60 80)
  have e : ((1 + y : ℚ) : ℝ) = 1 + (y:ℝ) := by push_cast; ring
  rw [e] at h
  unfold a0Q a0
  exact h

theorem a1Q_sound {e₁ e₂ : ℚ} (he₁ : 0 < e₁) (he₂ : 0 < e₂) : (a1Q e₁ e₂).Contains (a1 e₁ e₂) := by
  have h := QI.sub_sound (log2Q_sound he₂ 60 80) (log2Q_sound he₁ 60 80)
  unfold a1Q a1
  exact h

/-- An affine function with enclosed coefficients, on `[l, r]` with `0 ≤ l`, stays below `c` if the
upper corners do. -/
theorem affine_lt {A B : QI} {a b : ℝ} (ha : A.Contains a) (hb : B.Contains b) {l r c : ℚ}
    (hl : 0 ≤ l) (h1 : A.hi + B.hi * l < c) (h2 : A.hi + B.hi * r < c) :
    ∀ σ ∈ Set.Icc (l : ℝ) r, a + b * σ < c := by
  intro σ hσ
  obtain ⟨hl', hr'⟩ := hσ
  unfold QI.Contains at ha hb
  obtain ⟨_, ha2⟩ := ha
  obtain ⟨_, hb2⟩ := hb
  have hl0 : (0:ℝ) ≤ (l:ℝ) := by exact_mod_cast hl
  have hσ0 : 0 ≤ σ := le_trans hl0 hl'
  have h1' : (A.hi:ℝ) + (B.hi:ℝ) * (l:ℝ) < (c:ℝ) := by exact_mod_cast h1
  have h2' : (A.hi:ℝ) + (B.hi:ℝ) * (r:ℝ) < (c:ℝ) := by exact_mod_cast h2
  have hbσ : b * σ ≤ (B.hi:ℝ) * σ := mul_le_mul_of_nonneg_right hb2 hσ0
  rcases le_total 0 (B.hi:ℝ) with hB | hB
  · have k := mul_le_mul_of_nonneg_left hr' hB
    linarith
  · have k := mul_le_mul_of_nonpos_left hl' hB
    linarith

lemma qle_r {a b : ℚ} (h : a ≤ b) : (a : ℝ) ≤ (b : ℝ) := Rat.cast_le.mpr h

lemma qlt_r {a b : ℚ} (h : a < b) : (a : ℝ) < (b : ℝ) := Rat.cast_lt.mpr h

lemma w0_lt {σ : ℝ} (hσ : σ ∈ Set.Icc (31/100 : ℝ) (4172629/10^7)) :
    ω₀ σ < 782701512/10^9 := by
  have hA := a0Q_sound (y := y0) (e₂ := e2_0) (by norm_num [y0, q]) (by norm_num [e2_0, q])
  have hB := a1Q_sound (e₁ := e1_0) (e₂ := e2_0) (by norm_num [e1_0, q]) (by norm_num [e2_0, q])
  have h := affine_lt hA hB (l := 31/100) (r := 4172629/10^7) (c := 782701512/10^9)
    (by norm_num) (by decide +kernel) (by decide +kernel) σ
    (by rw [Set.mem_Icc] at hσ ⊢; push_cast; exact hσ)
  have h' := h
  push_cast at h'
  unfold ω₀
  rw [omegaCert_eq]
  exact h'

lemma w1_lt {σ : ℝ} (hσ : σ ∈ Set.Icc (4172629/10^7 : ℝ) (54/100)) :
    ω₁ σ < 78270152/10^8 := by
  have hA := a0Q_sound (y := y1) (e₂ := e2_1) (by norm_num [y1, q]) (by norm_num [e2_1, q])
  have hB := a1Q_sound (e₁ := e1_1) (e₂ := e2_1) (by norm_num [e1_1, q]) (by norm_num [e2_1, q])
  have h := affine_lt hA hB (l := 4172629/10^7) (r := 54/100) (c := 78270152/10^8)
    (by norm_num) (by decide +kernel) (by decide +kernel) σ
    (by rw [Set.mem_Icc] at hσ ⊢; push_cast; exact hσ)
  have h' := h
  push_cast at h'
  unfold ω₁
  rw [omegaCert_eq]
  exact h'

lemma pow_enc :
    (17203/10^4 : ℝ) ≤ (2:ℝ) ^ (7827016/10^7 : ℝ) ∧ (2:ℝ) ^ (7827016/10^7 : ℝ) < 17204/10^4 := by
  have h := pow2Q_sound (7827016/10^7) 60 80
  have e : ((7827016/10^7 : ℚ) : ℝ) = 7827016/10^7 := by norm_num
  rw [e] at h
  have hl : (17203/10^4 : ℚ) ≤ (pow2Q (7827016/10^7) 60 80).lo := by decide +kernel
  have hu : (pow2Q (7827016/10^7) 60 80).hi < (17204/10^4 : ℚ) := by decide +kernel
  have hl' := qle_r hl
  have hu' := qlt_r hu
  push_cast at hl' hu'
  exact ⟨le_trans hl' h.1, lt_of_le_of_lt h.2 hu'⟩

theorem rst_rates :
    (∀ σ ∈ Set.Icc (31/100 : ℝ) (4172629/10^7), ω₀ σ < 78270152/10^8) ∧
    (∀ σ ∈ Set.Icc (4172629/10^7 : ℝ) (54/100), ω₁ σ < 78270152/10^8) ∧
    (78270152/10^8 : ℝ) < 7827016/10^7 ∧ (7807/10000 : ℝ) < 7827016/10^7 ∧
    MLS.H2 (1/5) < 7802/10000 ∧
    (17203/10^4 : ℝ) ≤ (2:ℝ) ^ (7827016/10^7 : ℝ) ∧ (2:ℝ) ^ (7827016/10^7 : ℝ) < 17204/10^4 := by
  refine ⟨fun σ hσ => ?_, fun σ hσ => w1_lt hσ, ?_, ?_, ?_, pow_enc.1, pow_enc.2⟩
  · exact lt_trans (w0_lt hσ) (by norm_num)
  · norm_num
  · norm_num
  · have h := H2Q_sound (x := 1/5) (by norm_num) (by norm_num) 60 80
    have hh : (H2Q (1/5) 60 80).hi < 7802/10000 := by decide +kernel
    have hh' := qlt_r hh
    push_cast at hh'
    have h2 := h.2
    have e : ((1/5 : ℚ) : ℝ) = 1/5 := by norm_num
    rw [e] at h2
    linarith

/-- The enclosure of `ω₁(r)` at `r = 4172629/10^7`. -/
def ω1Q : QI := (a0Q y1 e2_1).add ((a1Q e1_1 e2_1).smul (4172629/10^7))

lemma omega1_r_bounds :
    (78270151357/10^11 : ℝ) ≤ ω₁ (4172629/10^7) ∧ ω₁ (4172629/10^7) < 78270151358/10^11 := by
  have hA := a0Q_sound (y := y1) (e₂ := e2_1) (by norm_num [y1, q]) (by norm_num [e2_1, q])
  have hB := a1Q_sound (e₁ := e1_1) (e₂ := e2_1) (by norm_num [e1_1, q]) (by norm_num [e2_1, q])
  have h : ω1Q.Contains (a0 y1 e2_1 + ((4172629/10^7 : ℚ) : ℝ) * a1 e1_1 e2_1) :=
    QI.add_sound hA (QI.smul_sound (4172629/10^7 : ℚ) hB)
  have hlo : (78270151357/10^11 : ℚ) ≤ ω1Q.lo := by decide +kernel
  have hhi : ω1Q.hi < (78270151358/10^11 : ℚ) := by decide +kernel
  have hlo' := qle_r hlo
  have hhi' := qlt_r hhi
  push_cast at hlo' hhi'
  have e : ω₁ (4172629/10^7) = a0 y1 e2_1 + ((4172629/10^7 : ℚ) : ℝ) * a1 e1_1 e2_1 := by
    unfold ω₁
    rw [omegaCert_eq]
    push_cast
    ring
  have c1 : ((4172629/10^7 : ℚ) : ℝ) = 4172629/10^7 := by norm_num
  rw [e]
  unfold QI.Contains at h
  obtain ⟨h1, h2⟩ := h
  rw [c1] at h1 h2 ⊢
  constructor <;> linarith [h1, h2]

/-- Numerator and denominator of the crossing point, as enclosures. -/
def numQ : QI := (a0Q y1 e2_1).sub (a0Q y0 e2_0)

def denQ : QI := (a1Q e1_0 e2_0).sub (a1Q e1_1 e2_1)

def ratQ : QI := numQ.div denQ

lemma abs_of_enc {I : QI} {x : ℝ} (h : I.Contains x) {c d : ℚ} (hl : c - d ≤ I.lo)
    (hh : I.hi ≤ c + d) : |x - c| ≤ d := by
  obtain ⟨h1, h2⟩ := h
  have hl' : ((c - d : ℚ) : ℝ) ≤ (I.lo : ℝ) := by exact_mod_cast hl
  have hh' : (I.hi : ℝ) ≤ ((c + d : ℚ) : ℝ) := by exact_mod_cast hh
  rw [abs_sub_le_iff]
  push_cast at hl' hh'
  constructor <;> linarith

/-- The detail of `prop:rst-certs`: the ten-decimal coefficients of
the table, the crossing point `σ× = 0.4172629270…`, and the maximum of the display, which lies in
`[0.78270151357, 0.78270151358)`. -/
theorem rst_rates_detail :
    |a0 y0 e2_0 - 7776061742/10^10| ≤ 1/(2*10^10) ∧
    |a1 e1_0 e2_0 - 122113367/10^10| ≤ 1/(2*10^10) ∧
    |a0 y1 e2_1 - 8022330021/10^10| ≤ 1/(2*10^10) ∧
    |a1 e1_1 e2_1 - (-468085913/10^10)| ≤ 1/(2*10^10) ∧
    (4172629270/10^10 : ℝ) ≤ (a0 y1 e2_1 - a0 y0 e2_0) / (a1 e1_0 e2_0 - a1 e1_1 e2_1) ∧
    (a0 y1 e2_1 - a0 y0 e2_0) / (a1 e1_0 e2_0 - a1 e1_1 e2_1) < 4172629271/10^10 ∧
    ∃ m : ℝ, IsGreatest (ω₀ '' Set.Icc (31/100) (4172629/10^7) ∪
        ω₁ '' Set.Icc (4172629/10^7) (54/100)) m ∧
      78270151357/10^11 ≤ m ∧ m < 78270151358/10^11 := by
  have hA0 := a0Q_sound (y := y0) (e₂ := e2_0) (by norm_num [y0, q]) (by norm_num [e2_0, q])
  have hB0 := a1Q_sound (e₁ := e1_0) (e₂ := e2_0) (by norm_num [e1_0, q]) (by norm_num [e2_0, q])
  have hA1 := a0Q_sound (y := y1) (e₂ := e2_1) (by norm_num [y1, q]) (by norm_num [e2_1, q])
  have hB1 := a1Q_sound (e₁ := e1_1) (e₂ := e2_1) (by norm_num [e1_1, q]) (by norm_num [e2_1, q])
  refine ⟨?_, ?_, ?_, ?_, ?_, ?_, ?_⟩
  · have := abs_of_enc hA0 (c := 7776061742/10^10) (d := 1/(2*10^10))
      (by decide +kernel) (by decide +kernel)
    push_cast at this
    exact this
  · have := abs_of_enc hB0 (c := 122113367/10^10) (d := 1/(2*10^10))
      (by decide +kernel) (by decide +kernel)
    push_cast at this
    exact this
  · have := abs_of_enc hA1 (c := 8022330021/10^10) (d := 1/(2*10^10))
      (by decide +kernel) (by decide +kernel)
    push_cast at this
    exact this
  · have := abs_of_enc hB1 (c := -468085913/10^10) (d := 1/(2*10^10))
      (by decide +kernel) (by decide +kernel)
    push_cast at this
    exact this
  · have hJ0 : (0:ℚ) < denQ.lo := by decide +kernel
    have hR : ratQ.Contains ((a0 y1 e2_1 - a0 y0 e2_0) / (a1 e1_0 e2_0 - a1 e1_1 e2_1)) :=
      QI.div_sound (QI.sub_sound hA1 hA0) (QI.sub_sound hB0 hB1) hJ0
    have dl : (4172629270/10^10 : ℚ) ≤ ratQ.lo := by decide +kernel
    have dl' := qle_r dl
    push_cast at dl'
    unfold QI.Contains at hR
    linarith [hR.1]
  · have hJ0 : (0:ℚ) < denQ.lo := by decide +kernel
    have hR : ratQ.Contains ((a0 y1 e2_1 - a0 y0 e2_0) / (a1 e1_0 e2_0 - a1 e1_1 e2_1)) :=
      QI.div_sound (QI.sub_sound hA1 hA0) (QI.sub_sound hB0 hB1) hJ0
    have du : ratQ.hi < (4172629271/10^10 : ℚ) := by decide +kernel
    have du' := qlt_r du
    push_cast at du'
    unfold QI.Contains at hR
    linarith [hR.2]
  · have hb : a1 e1_1 e2_1 < 0 := by
      have h2 := hB1.2
      have hneg : ((a1Q e1_1 e2_1).hi : ℝ) < 0 := by
        exact_mod_cast (by decide +kernel : (a1Q e1_1 e2_1).hi < 0)
      exact lt_of_le_of_lt h2 hneg
    have hm := omega1_r_bounds
    refine ⟨ω₁ (4172629/10^7), ⟨?_, ?_⟩, hm.1, hm.2⟩
    · exact Set.mem_union_right _ ⟨4172629/10^7, ⟨le_rfl, by norm_num⟩, rfl⟩
    · rintro x (⟨σ, hσ, rfl⟩ | ⟨σ, hσ, rfl⟩)
      · have := w0_lt hσ
        linarith [hm.1]
      · unfold ω₁
        rw [omegaCert_eq, omegaCert_eq]
        have k := mul_le_mul_of_nonpos_left hσ.1 hb.le
        linarith

end MaxLeafFormal.RSThetaRate