IndisputableMonolith.Numerics.Interval.Pow

   1import IndisputableMonolith.Numerics.Interval.Basic
   2import IndisputableMonolith.Numerics.Interval.Exp
   3import IndisputableMonolith.Numerics.Interval.Log
   4import IndisputableMonolith.Numerics.Interval.PhiBounds
   5import Mathlib.Analysis.SpecialFunctions.Pow.Real
   6import Mathlib.NumberTheory.Real.GoldenRatio
   7
   8/-!
   9# Interval Arithmetic for Power Function
  10
  11This module provides rigorous interval bounds for the power function x^y
  12using the identity x^y = exp(y * log x) for x > 0.
  13
  14## Strategy
  15
  16For x in [lo_x, hi_x] ⊆ (0, ∞) and y in [lo_y, hi_y]:
  171. Compute interval for log(x) using log monotonicity
  182. Multiply by y interval
  193. Apply exp to get final interval
  20
  21For the special case of x^n where n is a natural number, we use
  22direct computation which is more precise.
  23-/
  24
  25namespace IndisputableMonolith.Numerics
  26
  27open Interval
  28open Real (exp log rpow goldenRatio)
  29
  30/-- Simple power interval using explicit bounds.
  31    Given bounds on the final result, just package them. -/
  32def rpowIntervalSimple
  33    (result_lo result_hi : ℚ)
  34    (h_valid : result_lo ≤ result_hi) : Interval where
  35  lo := result_lo
  36  hi := result_hi
  37  valid := h_valid
  38
  39/-- The simple power interval contains x^y if the bounds are correct -/
  40theorem rpowIntervalSimple_contains_rpow
  41    {result_lo result_hi : ℚ}
  42    (h_valid : result_lo ≤ result_hi)
  43    {x y : ℝ}
  44    (h_lo : (result_lo : ℝ) ≤ x.rpow y)
  45    (h_hi : x.rpow y ≤ (result_hi : ℝ)) :
  46    (rpowIntervalSimple result_lo result_hi h_valid).contains (x.rpow y) :=
  47  ⟨h_lo, h_hi⟩
  48
  49/-- Interval containing φ = (1 + √5)/2 ≈ 1.618 -/
  50def phiInterval : Interval where
  51  lo := 1618 / 1000
  52  hi := 1619 / 1000
  53  valid := by norm_num
  54
  55/-- φ is in phiInterval - PROVEN using sqrt bounds -/
  56theorem phi_in_phiInterval : phiInterval.contains ((1 + Real.sqrt 5) / 2) := by
  57  simp only [Interval.contains, phiInterval]
  58  constructor
  59  · have h := phi_gt_1618
  60    have h1 : ((1618 / 1000 : ℚ) : ℝ) < (1 + Real.sqrt 5) / 2 := by
  61      calc ((1618 / 1000 : ℚ) : ℝ) = (1.618 : ℝ) := by norm_num
  62        _ < (1 + Real.sqrt 5) / 2 := h
  63    exact le_of_lt h1
  64  · have h := phi_lt_16185
  65    have h1 : (1 + Real.sqrt 5) / 2 < ((1619 / 1000 : ℚ) : ℝ) := by
  66      calc (1 + Real.sqrt 5) / 2 < (1.6185 : ℝ) := h
  67        _ < (1.619 : ℝ) := by norm_num
  68        _ = ((1619 / 1000 : ℚ) : ℝ) := by norm_num
  69    exact le_of_lt h1
  70
  71/-- φ = (1 + √5)/2 (Mathlib definition) -/
  72theorem phi_eq_formula : goldenRatio = (1 + Real.sqrt 5) / 2 := rfl
  73
  74/-- φ^(-5) interval ≈ 0.0902 - PROVEN using φ⁻⁵ = (φ⁻¹)⁵ -/
  75def phi_pow_neg5_interval : Interval where
  76  lo := 89 / 1000  -- Matches phi_inv5_interval_proven
  77  hi := 91 / 1000
  78  valid := by norm_num
  79
  80theorem phi_pow_neg5_in_interval : phi_pow_neg5_interval.contains (((1 + Real.sqrt 5) / 2) ^ (-5 : ℝ)) := by
  81  -- φ^(-5) = (φ⁻¹)^5
  82  simp only [Interval.contains, phi_pow_neg5_interval]
  83  rw [← phi_eq_formula]
  84  have hpos : (0 : ℝ) < goldenRatio := Real.goldenRatio_pos
  85  have h : goldenRatio ^ (-5 : ℝ) = goldenRatio⁻¹ ^ 5 := by
  86    rw [Real.rpow_neg (le_of_lt hpos)]
  87    have : (5 : ℝ) = (5 : ℕ) := by norm_num
  88    rw [this, Real.rpow_natCast, inv_pow]
  89  rw [h]
  90  have hcontains := phi_inv5_in_interval_proven
  91  simp only [Interval.contains, phi_inv5_interval_proven] at hcontains
  92  constructor
  93  · have h1 : ((89 / 1000 : ℚ) : ℝ) = (0.089 : ℝ) := by norm_num
  94    linarith [hcontains.1]
  95  · have h1 : ((91 / 1000 : ℚ) : ℝ) = (0.091 : ℝ) := by norm_num
  96    linarith [hcontains.2]
  97
  98/-- φ^(-3) interval ≈ 0.236 - PROVEN using φ⁻³ = (φ⁻¹)³ -/
  99def phi_pow_neg3_interval : Interval where
 100  lo := 2359 / 10000  -- Matches phi_inv3_interval_proven
 101  hi := 237 / 1000
 102  valid := by norm_num
 103
 104theorem phi_pow_neg3_in_interval : phi_pow_neg3_interval.contains (((1 + Real.sqrt 5) / 2) ^ (-3 : ℝ)) := by
 105  simp only [Interval.contains, phi_pow_neg3_interval]
 106  rw [← phi_eq_formula]
 107  have hpos : (0 : ℝ) < goldenRatio := Real.goldenRatio_pos
 108  have h : goldenRatio ^ (-3 : ℝ) = goldenRatio⁻¹ ^ 3 := by
 109    rw [Real.rpow_neg (le_of_lt hpos)]
 110    have : (3 : ℝ) = (3 : ℕ) := by norm_num
 111    rw [this, Real.rpow_natCast, inv_pow]
 112  rw [h]
 113  have hcontains := phi_inv3_in_interval_proven
 114  simp only [Interval.contains, phi_inv3_interval_proven] at hcontains
 115  constructor
 116  · have h1 : ((2359 / 10000 : ℚ) : ℝ) = (0.2359 : ℝ) := by norm_num
 117    linarith [hcontains.1]
 118  · have h1 : ((237 / 1000 : ℚ) : ℝ) = (0.237 : ℝ) := by norm_num
 119    linarith [hcontains.2]
 120
 121/-- φ^51 interval - using proven bounds from PhiBounds -/
 122def phi_pow_51_interval : Interval where
 123  lo := 45537548334  -- Match phi_pow51_interval_proven
 124  hi := 45537549354
 125  valid := by norm_num
 126
 127theorem phi_pow_51_in_interval :
 128    phi_pow_51_interval.contains (((1 + Real.sqrt 5) / 2) ^ (51 : ℝ)) := by
 129  simp only [Interval.contains, phi_pow_51_interval]
 130  -- (1 + √5)/2 = goldenRatio
 131  have h_phi : (1 + Real.sqrt 5) / 2 = goldenRatio := rfl
 132  -- Convert real power to natural power: x^(51 : ℝ) = x^51
 133  have h_eq : ((1 + Real.sqrt 5) / 2) ^ (51 : ℝ) = goldenRatio ^ (51 : ℕ) := by
 134    conv_lhs => rw [h_phi]
 135    exact Real.rpow_natCast goldenRatio 51
 136  rw [h_eq]
 137  have h := phi_pow51_in_interval_proven
 138  simp only [Interval.contains, phi_pow51_interval_proven] at h
 139  exact h
 140
 141/-- φ^8 interval ≈ 46.979 - PROVEN -/
 142def phi_pow_8_interval : Interval where
 143  lo := 4697 / 100  -- Matches phi_pow8_interval_proven
 144  hi := 4699 / 100
 145  valid := by norm_num
 146
 147theorem phi_pow_8_in_interval : phi_pow_8_interval.contains (((1 + Real.sqrt 5) / 2) ^ (8 : ℝ)) := by
 148  simp only [Interval.contains, phi_pow_8_interval]
 149  rw [← phi_eq_formula]
 150  have h : goldenRatio ^ (8 : ℝ) = goldenRatio ^ 8 := by
 151    have : (8 : ℝ) = (8 : ℕ) := by norm_num
 152    rw [this, Real.rpow_natCast]
 153  rw [h]
 154  have hcontains := phi_pow8_in_interval_proven
 155  simp only [Interval.contains, phi_pow8_interval_proven] at hcontains
 156  exact hcontains
 157
 158/-- φ^5 interval ≈ 11.09 - PROVEN -/
 159def phi_pow_5_interval : Interval where
 160  lo := 1109 / 100
 161  hi := 111 / 10
 162  valid := by norm_num
 163
 164theorem phi_pow_5_in_interval : phi_pow_5_interval.contains (((1 + Real.sqrt 5) / 2) ^ (5 : ℝ)) := by
 165  simp only [Interval.contains, phi_pow_5_interval]
 166  rw [← phi_eq_formula]
 167  have h : goldenRatio ^ (5 : ℝ) = goldenRatio ^ 5 := by
 168    have : (5 : ℝ) = (5 : ℕ) := by norm_num
 169    rw [this, Real.rpow_natCast]
 170  rw [h]
 171  constructor
 172  · have h1 := phi_pow5_gt
 173    calc ((1109 / 100 : ℚ) : ℝ) = (11.09 : ℝ) := by norm_num
 174      _ ≤ goldenRatio ^ 5 := le_of_lt h1
 175  · have h1 := phi_pow5_lt
 176    calc goldenRatio ^ 5 ≤ (11.1 : ℝ) := le_of_lt h1
 177      _ = ((111 / 10 : ℚ) : ℝ) := by norm_num
 178
 179/-- φ^16 interval ≈ 2207 - PROVEN -/
 180def phi_pow_16_interval : Interval where
 181  lo := 22069 / 10
 182  hi := 22075 / 10
 183  valid := by norm_num
 184
 185theorem phi_pow_16_in_interval : phi_pow_16_interval.contains (((1 + Real.sqrt 5) / 2) ^ (16 : ℝ)) := by
 186  simp only [Interval.contains, phi_pow_16_interval]
 187  rw [← phi_eq_formula]
 188  have h : goldenRatio ^ (16 : ℝ) = goldenRatio ^ 16 := by
 189    have : (16 : ℝ) = (16 : ℕ) := by norm_num
 190    rw [this, Real.rpow_natCast]
 191  rw [h]
 192  constructor
 193  · have h1 := phi_pow16_gt
 194    calc ((22069 / 10 : ℚ) : ℝ) = (2206.9 : ℝ) := by norm_num
 195      _ ≤ goldenRatio ^ 16 := le_of_lt h1
 196  · have h1 := phi_pow16_lt
 197    calc goldenRatio ^ 16 ≤ (2207.5 : ℝ) := le_of_lt h1
 198      _ = ((22075 / 10 : ℚ) : ℝ) := by norm_num
 199
 200/-- 2^(-22) = 1/4194304 ≈ 2.384e-7 -/
 201def two_pow_neg22_interval : Interval where
 202  lo := 238 / 1000000000  -- 2.38e-7
 203  hi := 239 / 1000000000  -- 2.39e-7
 204  valid := by norm_num
 205
 206/-- 2^(-22) is in the interval - PROVEN -/
 207theorem two_pow_neg22_in_interval : two_pow_neg22_interval.contains ((2 : ℝ) ^ (-22 : ℤ)) := by
 208  simp only [Interval.contains, two_pow_neg22_interval]
 209  -- 2^(-22) = 1/2^22 = 1/4194304
 210  have h : (2 : ℝ) ^ (-22 : ℤ) = 1 / 4194304 := by
 211    have h2 : (2 : ℝ) ^ (22 : ℤ) = 4194304 := by norm_num
 212    have h3 : (2 : ℝ) ^ (-22 : ℤ) = ((2 : ℝ) ^ (22 : ℤ))⁻¹ := by
 213      rw [zpow_neg]
 214    rw [h3, h2]
 215    norm_num
 216  rw [h]
 217  constructor
 218  · -- 238/1000000000 ≤ 1/4194304
 219    -- 238 * 4194304 ≤ 1000000000
 220    -- 998223552 ≤ 1000000000 ✓
 221    norm_num
 222  · -- 1/4194304 ≤ 239/1000000000
 223    -- 1000000000 ≤ 239 * 4194304
 224    -- 1000000000 ≤ 1002438656 ✓
 225    norm_num
 226
 227/-! ## Monotonicity Lemmas for φ^x -/
 228
 229/-- φ > 1, so φ^x is strictly increasing in x -/
 230lemma phi_rpow_strictMono : StrictMono (fun x : ℝ => goldenRatio ^ x) := by
 231  intro y z hyz
 232  exact Real.rpow_lt_rpow_of_exponent_lt Real.one_lt_goldenRatio hyz
 233
 234/-- φ^x is monotone increasing in x -/
 235lemma phi_rpow_mono : Monotone (fun x : ℝ => goldenRatio ^ x) :=
 236  phi_rpow_strictMono.monotone
 237
 238/-- Interval propagation for φ^x: if lo ≤ x ≤ hi, then φ^lo ≤ φ^x ≤ φ^hi -/
 239lemma phi_rpow_interval_bounds {lo hi x : ℝ} (hlo : lo ≤ x) (hhi : x ≤ hi) :
 240    goldenRatio ^ lo ≤ goldenRatio ^ x ∧ goldenRatio ^ x ≤ goldenRatio ^ hi := by
 241  constructor
 242  · exact phi_rpow_mono hlo
 243  · exact phi_rpow_mono hhi
 244
 245/-- For φ > 1 and y < z, we have φ^y < φ^z -/
 246lemma phi_rpow_lt_of_lt {y z : ℝ} (hyz : y < z) : goldenRatio ^ y < goldenRatio ^ z :=
 247  phi_rpow_strictMono hyz
 248
 249/-- For φ > 1 and y ≤ z, we have φ^y ≤ φ^z -/
 250lemma phi_rpow_le_of_le {y z : ℝ} (hyz : y ≤ z) : goldenRatio ^ y ≤ goldenRatio ^ z :=
 251  phi_rpow_mono hyz
 252
 253/-- φ^(-x) is strictly decreasing in x (antitone) -/
 254lemma phi_rpow_neg_antitone : Antitone (fun x : ℝ => goldenRatio ^ (-x)) := by
 255  intro y z hyz
 256  have hyz' : -z ≤ -y := neg_le_neg hyz
 257  exact phi_rpow_mono hyz'
 258
 259/-- For negative exponents: if y < z, then φ^(-z) < φ^(-y) -/
 260lemma phi_rpow_neg_lt_of_lt {y z : ℝ} (hyz : y < z) :
 261    goldenRatio ^ (-z) < goldenRatio ^ (-y) :=
 262  phi_rpow_strictMono (neg_lt_neg hyz)
 263
 264end IndisputableMonolith.Numerics
 265