Mathlib.Analysis.Certified­Contour­Interval­Arithmetic.Complex

theorem piRaw_sound (fuel : ℕ) : (piRaw fuel).Contains pi := by have hfour := RatInterval.point_sound 4 have h5 := atanRaw_sound_of_nonneg_lt_one (1 / 5) (by norm_num) (by norm_num) fuel have h239 := atanRaw_sound_of_nonneg_lt_one (1 / 239) (by norm_num) (by norm_num) fuel have h := RatInterval.mul_sound hfour (RatInterval.sub_sound (RatInterval.mul_sound hfour h5) h239) have hm := Real.four_mul_arctan_inv_5_sub_arctan_inv_239 have heq : (4 : ℝ) * ((4 : ℝ) * Real.arctan (1 / 5 : ℝ) - Real.arctan (1 / 239 : ℝ)) = pi := by norm_num [div_eq_mul_inv] nlinarith [hm] rw [← heq] simpa [piRaw] using h

theorem piInterval_sound (fuel : ℕ) : (piInterval fuel).Contains pi := by induction fuel with | zero => exact piRaw_sound 0 | succ fuel ih => exact RatInterval.tighten_sound ih (piRaw_sound (fuel + 1))

theorem sinRaw_sound (q : ℚ) (fuel : ℕ) : (sinRaw q fuel).Contains (sin q) := by have hb := sin_taylor_error_bound (q : ℝ) fuel have hp : ((sinPartial q fuel : ℚ) : ℝ) = ∑ k ∈ range (fuel + 1), (-1 : ℝ) ^ k * (q : ℝ) ^ (2 * k + 1) / ((2 * k + 1).factorial : ℝ) := by simp [sinPartial] have he : ((sinError q fuel : ℚ) : ℝ) = |(q : ℝ)| ^ (2 * fuel + 3) / ((2 * fuel + 3).factorial : ℝ) := by simp [sinError] rw [← hp, ← he, abs_le] at hb exact ⟨by simpa [sinRaw] using hb.1, by simpa [sinRaw, add_comm] using hb.2⟩

theorem cosRaw_sound (q : ℚ) (fuel : ℕ) : (cosRaw q fuel).Contains (cos q) := by have hb := cos_taylor_error_bound (q : ℝ) fuel have hp : ((cosPartial q fuel : ℚ) : ℝ) = ∑ k ∈ range (fuel + 1), (-1 : ℝ) ^ k * (q : ℝ) ^ (2 * k) / ((2 * k).factorial : ℝ) := by simp [cosPartial] have he : ((cosError q fuel : ℚ) : ℝ) = |(q : ℝ)| ^ (2 * fuel + 2) / ((2 * fuel + 2).factorial : ℝ) := by simp [cosError] rw [← hp, ← he, abs_le] at hb exact ⟨by simpa [cosRaw] using hb.1, by simpa [cosRaw, add_comm] using hb.2⟩

theorem abs_sub_intervalMid_le_radius {I : RatInterval} {x : ℝ} (hx : I.Contains x) : |x - (intervalMid I : ℝ)| ≤ (intervalRadius I : ℝ) := by have hlo : (I.lo : ℝ) ≤ x := hx.1 have hhi : x ≤ (I.hi : ℝ) := hx.2 rw [abs_le] simp only [intervalMid, intervalRadius, RatInterval.width, Rat.cast_div, Rat.cast_add, Rat.cast_sub, Rat.cast_ofNat] constructor <;> linarith

theorem sinInterval_sound {I : RatInterval} {x : ℝ} (hx : I.Contains x) (fuel : ℕ) : (sinInterval I fuel).Contains (sin x) := by induction fuel with | zero => exact sinIntervalRaw_sound hx 0 | succ fuel ih => exact RatInterval.tighten_sound ih (sinIntervalRaw_sound hx (fuel + 1))

theorem cosInterval_sound {I : RatInterval} {x : ℝ} (hx : I.Contains x) (fuel : ℕ) : (cosInterval I fuel).Contains (cos x) := by induction fuel with | zero => exact cosIntervalRaw_sound hx 0 | succ fuel ih => exact RatInterval.tighten_sound ih (cosIntervalRaw_sound hx (fuel + 1))

theorem sinInterval_nested (I : RatInterval) (fuel : ℕ) : (sinInterval I (fuel + 1)).Subinterval (sinInterval I fuel) := by have hmid : I.Contains (intervalMid I : ℝ) := by constructor <;> exact_mod_cast (by dsimp [intervalMid] linarith [I.lo_le_hi]) exact RatInterval.tighten_subinterval_left (sinInterval_sound hmid fuel) (sinIntervalRaw_sound hmid (fuel + 1))

theorem cosInterval_nested (I : RatInterval) (fuel : ℕ) : (cosInterval I (fuel + 1)).Subinterval (cosInterval I fuel) := by have hmid : I.Contains (intervalMid I : ℝ) := by constructor <;> exact_mod_cast (by dsimp [intervalMid] linarith [I.lo_le_hi]) exact RatInterval.tighten_subinterval_left (cosInterval_sound hmid fuel) (cosIntervalRaw_sound hmid (fuel + 1))

theorem sinInterval_width (I : RatInterval) (fuel : ℕ) : (sinInterval I fuel).width ≤ I.width + 2 * sinError (intervalMid I) fuel := by have hraw : (sinIntervalRaw I fuel).width = I.width + 2 * sinError (intervalMid I) fuel := by rw [sinIntervalRaw] simp [RatInterval.expand, sinRaw, intervalRadius, RatInterval.width] ring cases fuel with | zero => exact hraw.le | succ fuel => have hmid : I.Contains (intervalMid I : ℝ) := by constructor <;> exact_mod_cast (by dsimp [intervalMid] linarith [I.lo_le_hi]) exact (RatInterval.width_mono (tighten_subinterval_right (sinInterval_sound hmid fuel) (sinIntervalRaw_sound hmid (fuel + 1)))).trans_eq hraw

theorem cosInterval_width (I : RatInterval) (fuel : ℕ) : (cosInterval I fuel).width ≤ I.width + 2 * cosError (intervalMid I) fuel := by have hraw : (cosIntervalRaw I fuel).width = I.width + 2 * cosError (intervalMid I) fuel := by rw [cosIntervalRaw] simp [RatInterval.expand, cosRaw, intervalRadius, RatInterval.width] ring cases fuel with | zero => exact hraw.le | succ fuel => have hmid : I.Contains (intervalMid I : ℝ) := by constructor <;> exact_mod_cast (by dsimp [intervalMid] linarith [I.lo_le_hi]) exact (RatInterval.width_mono (tighten_subinterval_right (cosInterval_sound hmid fuel) (cosIntervalRaw_sound hmid (fuel + 1)))).trans_eq hraw

theorem sq_sound {I : RatInterval} {x : ℝ} (hx : I.Contains x) : I.sq.Contains (x ^ 2) := by rcases hx with ⟨hxlo, hxhi⟩ unfold sq split_ifs with hneg hpos · have hneg' : (I.hi : ℝ) < 0 := by exact_mod_cast hneg constructor <;> simp only [Rat.cast_pow] · nlinarith [hxhi] · nlinarith [hxlo, hxhi] · have hpos' : (0 : ℝ) < I.lo := by exact_mod_cast hpos have hhi' : (0 : ℝ) < I.hi := hpos'.trans_le (by exact_mod_cast I.lo_le_hi) constructor <;> simp only [Rat.cast_pow] · nlinarith [hxlo] · nlinarith [hxhi] · constructor · simpa using sq_nonneg x · simp only [Rat.cast_max, Rat.cast_pow] by_cases hx0 : x ≤ 0 · exact le_max_of_le_left (by nlinarith [hxlo]) · exact le_max_of_le_right (by nlinarith [hxhi])

theorem sq_mono {I J : RatInterval} (hIJ : I.Subinterval J) : I.sq.Subinterval J.sq := by have Ilo : J.Contains (I.lo : ℝ) := by constructor · exact_mod_cast hIJ.1 · exact_mod_cast I.lo_le_hi.trans hIJ.2 have Ihi : J.Contains (I.hi : ℝ) := by constructor · exact_mod_cast hIJ.1.trans I.lo_le_hi · exact_mod_cast hIJ.2 have hlo := sq_sound Ilo have hhi := sq_sound Ihi rw [sq] split_ifs with Ihneg Ihpos · exact ⟨by exact_mod_cast hhi.1, by exact_mod_cast hlo.2⟩ · exact ⟨by exact_mod_cast hlo.1, by exact_mod_cast hhi.2⟩ · have hIlo : I.lo ≤ 0 := le_of_not_gt Ihpos have hIhi : 0 ≤ I.hi := le_of_not_gt Ihneg have hzero : J.Contains (0 : ℝ) := by constructor · exact_mod_cast hIJ.1.trans hIlo · exact_mod_cast hIhi.trans hIJ.2 have hz := sq_sound hzero constructor · exact_mod_cast hz.1 · exact max_le (by exact_mod_cast hlo.2) (by exact_mod_cast hhi.2)

theorem sqrt_gap_rate (q : ℚ) (hq : 0 ≤ q) (n : ℕ) : sqrtUpper q n - sqrtLower q n ≤ (q + 1) / (n + 1) := by induction n with | zero => rw [sqrtUpper, abs_of_nonneg hq] have hl : 0 ≤ sqrtLower q 0 := div_nonneg hq (sqrtUpper_invariants q hq 0).1.le norm_num linarith | succ n ih => have hstep : sqrtUpper q (n + 1) - sqrtLower q (n + 1) ≤ (sqrtUpper q n - sqrtLower q n) / 2 := by rw [sqrtUpper] rw [show q / sqrtUpper q n = sqrtLower q n by rfl] have hm := sqrtLower_mono_step q hq n linarith have hA : 0 ≤ q + 1 := by linarith calc sqrtUpper q (n + 1) - sqrtLower q (n + 1) ≤ (sqrtUpper q n - sqrtLower q n) / 2 := hstep _ ≤ ((q + 1) / (n + 1)) / 2 := div_le_div_of_nonneg_right ih (by norm_num) _ ≤ (q + 1) / (n + 1 + 1) := by rw [div_div] apply div_le_div_of_nonneg_left hA (by positivity) have hn0 : (0 : ℚ) ≤ n := by positivity push_cast nlinarith _ = (q + 1) / ((n + 1 : ℕ) + 1) := by norm_num [Nat.cast_add]

theorem sqrt_iterates_converge (q : ℚ) (hq : 0 ≤ q) (ε : PosRat) : ∃ n, sqrtUpper q n - sqrtLower q n ≤ ε.1 := by obtain ⟨n, hn⟩ := exists_nat_gt ((q + 1) / ε.1) refine ⟨n, (sqrt_gap_rate q hq n).trans ?_⟩ apply (div_le_iff₀ (by positivity : (0 : ℚ) < n + 1)).2 have hn' : q + 1 < ε.1 * n := by simpa [mul_comm] using (div_lt_iff₀ ε.2).mp hn push_cast nlinarith [ε.2]

theorem sqrtInterval_sound {I : RatInterval} (hI : 0 ≤ I.lo) {x : ℝ} (hx : I.Contains x) (fuel : ℕ) : (sqrtInterval I hI fuel).Contains (sqrt x) := by have hhi : 0 ≤ I.hi := hI.trans I.lo_le_hi have hlo_sound := sqrt_iterates_sound I.lo hI fuel have hhi_sound := sqrt_iterates_sound I.hi hhi fuel constructor · exact hlo_sound.1.trans (Real.sqrt_le_sqrt hx.1) · exact (Real.sqrt_le_sqrt hx.2).trans hhi_sound.2

theorem maxAbs_le_of_contains_width {I : RatInterval} {x : ℝ} {C w : ℚ} (hx : I.Contains x) (habs : |x| ≤ C) (hw : I.width ≤ w) : I.maxAbs ≤ C + w := by rw [abs_le] at habs rcases hx with ⟨hxlo, hxhi⟩ rcases habs with ⟨habslo, habshi⟩ have hw' : (I.hi : ℝ) - I.lo ≤ w := by exact_mod_cast hw apply max_le · rw [abs_le] constructor · exact_mod_cast (show (-(C + w) : ℝ) ≤ (I.lo : ℚ) by push_cast linarith) · exact_mod_cast (show (((I.lo : ℚ) : ℝ) ≤ (C + w : ℚ)) by push_cast linarith) · rw [abs_le] constructor · exact_mod_cast (show (-(C + w) : ℝ) ≤ (I.hi : ℚ) by push_cast linarith) · exact_mod_cast (show (((I.hi : ℚ) : ℝ) ≤ (C + w : ℚ)) by push_cast linarith)

theorem maxAbs_add (I J : RatInterval) : (I.add J).maxAbs ≤ I.maxAbs + J.maxAbs := by apply max_le · simpa [RatInterval.add, RatInterval.maxAbs] using (abs_add_le I.lo J.lo).trans (add_le_add (le_max_left _ _) (le_max_left _ _)) · simpa [RatInterval.add, RatInterval.maxAbs] using (abs_add_le I.hi J.hi).trans (add_le_add (le_max_right _ _) (le_max_right _ _))

theorem maxAbs_mul (I J : RatInterval) : (I.mul J).maxAbs ≤ I.maxAbs * J.maxAbs := by have hI0 : 0 ≤ I.maxAbs := (abs_nonneg I.lo).trans (le_max_left _ _) have hp (x y : ℚ) (hx : |x| ≤ I.maxAbs) (hy : |y| ≤ J.maxAbs) : |x * y| ≤ I.maxAbs * J.maxAbs := by rw [abs_mul] exact mul_le_mul hx hy (abs_nonneg y) hI0 have ha := hp I.lo J.lo (le_max_left _ _) (le_max_left _ _) have hb := hp I.lo J.hi (le_max_left _ _) (le_max_right _ _) have hc := hp I.hi J.lo (le_max_right _ _) (le_max_left _ _) have hd := hp I.hi J.hi (le_max_right _ _) (le_max_right _ _) simp only [RatInterval.maxAbs, RatInterval.mul] apply max_le · calc |min (min (I.lo * J.lo) (I.lo * J.hi)) (min (I.hi * J.lo) (I.hi * J.hi))| ≤ max |min (I.lo * J.lo) (I.lo * J.hi)| |min (I.hi * J.lo) (I.hi * J.hi)| := abs_min_le_max_abs_abs _ ≤ max (max |I.lo * J.lo| |I.lo * J.hi|) (max |I.hi * J.lo| |I.hi * J.hi|) := max_le_max abs_min_le_max_abs_abs abs_min_le_max_abs_abs _ ≤ I.maxAbs * J.maxAbs := max_le (max_le ha hb) (max_le hc hd) · calc |max (max (I.lo * J.lo) (I.lo * J.hi)) (max (I.hi * J.lo) (I.hi * J.hi))| ≤ max |max (I.lo * J.lo) (I.lo * J.hi)| |max (I.hi * J.lo) (I.hi * J.hi)| := abs_max_le_max_abs_abs _ ≤ max (max |I.lo * J.lo| |I.lo * J.hi|) (max |I.hi * J.lo| |I.hi * J.hi|) := max_le_max abs_max_le_max_abs_abs abs_max_le_max_abs_abs _ ≤ I.maxAbs * J.maxAbs := max_le (max_le ha hb) (max_le hc hd)

theorem point_sound (x y : ℚ) : (point x y).Contains ((x : ℝ) + (y : ℝ) * I) := by constructor · simpa [point] using RatInterval.point_sound x · simpa [point] using RatInterval.point_sound y

theorem maxAbs_le_of_contains_width {I : ComplexRatInterval} {z : ℂ} {C w : ℚ} (hz : I.Contains z) (hzabs : max |z.re| |z.im| ≤ (C : ℝ)) (hw : I.width ≤ w) : I.maxAbs ≤ C + w := by apply max_le · exact RatInterval.maxAbs_le_of_contains_width hz.1 ((le_max_left _ _).trans hzabs) ((le_max_left I.re.width I.im.width).trans hw) · exact RatInterval.maxAbs_le_of_contains_width hz.2 ((le_max_right _ _).trans hzabs) ((le_max_right I.re.width I.im.width).trans hw)

theorem mul_sound {I J : ComplexRatInterval} {z w : ℂ} (hz : I.Contains z) (hw : J.Contains w) : (I.mul J).Contains (z * w) := by constructor · simpa [mul, Complex.mul_re] using RatInterval.sub_sound (RatInterval.mul_sound hz.1 hw.1) (RatInterval.mul_sound hz.2 hw.2) · simpa [mul, Complex.mul_im] using RatInterval.add_sound (RatInterval.mul_sound hz.1 hw.2) (RatInterval.mul_sound hz.2 hw.1)

theorem normSq_sound {I : ComplexRatInterval} {z : ℂ} (hz : I.Contains z) : I.normSq.Contains (‖z‖ ^ 2) := by rw [Complex.sq_norm] simpa [normSq, Complex.normSq_apply, pow_two] using RatInterval.add_sound (RatInterval.sq_sound hz.1) (RatInterval.sq_sound hz.2)

theorem tighten_sound_left {I J : ComplexRatInterval} {z : ℂ} (hI : I.Contains z) (hJ : J.Contains z) : (I.tighten J).Contains z ∧ (I.tighten J).Subinterval I := by exact ⟨⟨RatInterval.tighten_sound hI.1 hJ.1, RatInterval.tighten_sound hI.2 hJ.2⟩, ⟨RatInterval.tighten_subinterval_left hI.1 hJ.1, RatInterval.tighten_subinterval_left hI.2 hJ.2⟩⟩

theorem mul_maxAbs (I J : ComplexRatInterval) : (I.mul J).maxAbs ≤ 2 * I.maxAbs * J.maxAbs := by have hIre : I.re.maxAbs ≤ I.maxAbs := le_max_left _ _ have hIim : I.im.maxAbs ≤ I.maxAbs := le_max_right _ _ have hJre : J.re.maxAbs ≤ J.maxAbs := le_max_left _ _ have hJim : J.im.maxAbs ≤ J.maxAbs := le_max_right _ _ have hI0 : 0 ≤ I.maxAbs := (abs_nonneg I.re.lo).trans ((le_max_left _ _).trans (le_max_left _ _)) have hJ0 : 0 ≤ J.maxAbs := (abs_nonneg J.re.lo).trans ((le_max_left _ _).trans (le_max_left _ _)) have hrr := RatInterval.maxAbs_mul I.re J.re have hii := RatInterval.maxAbs_mul I.im J.im have hri := RatInterval.maxAbs_mul I.re J.im have hir := RatInterval.maxAbs_mul I.im J.re have hrr' : (I.re.mul J.re).maxAbs ≤ I.maxAbs * J.maxAbs := hrr.trans (mul_le_mul hIre hJre ((abs_nonneg J.re.lo).trans (le_max_left _ _)) hI0) have hii' : (I.im.mul J.im).maxAbs ≤ I.maxAbs * J.maxAbs := hii.trans (mul_le_mul hIim hJim ((abs_nonneg J.im.lo).trans (le_max_left _ _)) hI0) have hri' : (I.re.mul J.im).maxAbs ≤ I.maxAbs * J.maxAbs := hri.trans (mul_le_mul hIre hJim ((abs_nonneg J.im.lo).trans (le_max_left _ _)) hI0) have hir' : (I.im.mul J.re).maxAbs ≤ I.maxAbs * J.maxAbs := hir.trans (mul_le_mul hIim hJre ((abs_nonneg J.re.lo).trans (le_max_left _ _)) hI0) apply max_le · exact (RatInterval.maxAbs_sub _ _).trans (by linarith) · exact (RatInterval.maxAbs_add _ _).trans (by linarith)

theorem rat_sq_hi_le (K : RatInterval) : K.sq.hi ≤ K.maxAbs ^ 2 := by have hlo : |K.lo| ≤ K.maxAbs := le_max_left _ _ have hhi : |K.hi| ≤ K.maxAbs := le_max_right _ _ have hM : 0 ≤ K.maxAbs := (abs_nonneg K.lo).trans hlo have hlo' : K.lo ^ 2 ≤ K.maxAbs ^ 2 := (sq_le_sq).2 (by simpa [abs_of_nonneg hM] using hlo) have hhi' : K.hi ^ 2 ≤ K.maxAbs ^ 2 := (sq_le_sq).2 (by simpa [abs_of_nonneg hM] using hhi) rw [RatInterval.sq] split_ifs · exact hlo' · exact hhi' · exact max_le hlo' hhi'

theorem rat_maxAbs_mono {K L : RatInterval} (hKL : K.Subinterval L) : K.maxAbs ≤ L.maxAbs := by apply max_le · exact abs_le_max_abs_abs hKL.1 (K.lo_le_hi.trans hKL.2) · exact abs_le_max_abs_abs (hKL.1.trans K.lo_le_hi) hKL.2

theorem mul_width (I J : ComplexRatInterval) : (I.mul J).width ≤ 2 * (I.maxAbs * J.width + J.maxAbs * I.width) := by have hIreA : I.re.maxAbs ≤ I.maxAbs := le_max_left _ _ have hIimA : I.im.maxAbs ≤ I.maxAbs := le_max_right _ _ have hJreA : J.re.maxAbs ≤ J.maxAbs := le_max_left _ _ have hJimA : J.im.maxAbs ≤ J.maxAbs := le_max_right _ _ have hIrew : I.re.width ≤ I.width := le_max_left _ _ have hIimw : I.im.width ≤ I.width := le_max_right _ _ have hJrew : J.re.width ≤ J.width := le_max_left _ _ have hJimw : J.im.width ≤ J.width := le_max_right _ _ have hIA0 : 0 ≤ I.maxAbs := (abs_nonneg I.re.lo).trans ((le_max_left _ _).trans (le_max_left _ _)) have hJA0 : 0 ≤ J.maxAbs := (abs_nonneg J.re.lo).trans ((le_max_left _ _).trans (le_max_left _ _)) have hrr := interval_mul_width_le I.re J.re have hii := interval_mul_width_le I.im J.im have hri := interval_mul_width_le I.re J.im have hir := interval_mul_width_le I.im J.re have hrr' : (I.re.mul J.re).width ≤ I.maxAbs * J.width + J.maxAbs * I.width := by nlinarith [mul_le_mul hIreA hJrew (RatInterval.width_nonneg J.re) hIA0, mul_le_mul hJreA hIrew (RatInterval.width_nonneg I.re) hJA0] have hii' : (I.im.mul J.im).width ≤ I.maxAbs * J.width + J.maxAbs * I.width := by nlinarith [mul_le_mul hIimA hJimw (RatInterval.width_nonneg J.im) hIA0, mul_le_mul hJimA hIimw (RatInterval.width_nonneg I.im) hJA0] have hri' : (I.re.mul J.im).width ≤ I.maxAbs * J.width + J.maxAbs * I.width := by nlinarith [mul_le_mul hIreA hJimw (RatInterval.width_nonneg J.im) hIA0, mul_le_mul hJimA hIrew (RatInterval.width_nonneg I.re) hJA0] have hir' : (I.im.mul J.re).width ≤ I.maxAbs * J.width + J.maxAbs * I.width := by nlinarith [mul_le_mul hIimA hJrew (RatInterval.width_nonneg J.re) hIA0, mul_le_mul hJreA hIimw (RatInterval.width_nonneg I.im) hJA0] change max ((I.re.mul J.re).sub (I.im.mul J.im)).width ((I.re.mul J.im).add (I.im.mul J.re)).width ≤ 2 * (I.maxAbs * J.width + J.maxAbs * I.width) rw [RatInterval.width_sub, RatInterval.width_add] exact max_le (by linarith) (by linarith)

theorem normSq_width (I : ComplexRatInterval) : I.normSq.width ≤ 4 * I.maxAbs * I.width := by have hre := rat_sq_width_le I.re have him := rat_sq_width_le I.im have hA0 : 0 ≤ I.maxAbs := (abs_nonneg I.re.lo).trans ((le_max_left _ _).trans (le_max_left _ _)) have hreprod : I.re.maxAbs * I.re.width ≤ I.maxAbs * I.width := mul_le_mul (le_max_left _ _) (le_max_left _ _) (RatInterval.width_nonneg I.re) hA0 have himprod : I.im.maxAbs * I.im.width ≤ I.maxAbs * I.width := mul_le_mul (le_max_right _ _) (le_max_right _ _) (RatInterval.width_nonneg I.im) hA0 rw [normSq, RatInterval.width_add] nlinarith

theorem div_maxAbs {I J : ComplexRatInterval} (hJ : J.normSq.AwayFromZero) {δ : ℚ} (hδ : 0 < δ) (hsep : δ ≤ J.normSq.lo) : (I.div J hJ).maxAbs ≤ (I.mul J.conj).maxAbs / δ := by let N := I.mul J.conj let D := J.normSq have hlo : 0 < D.lo := hδ.trans_le hsep have hhi : 0 < D.hi := hlo.trans_le D.lo_le_hi have hinvA : (D.inv hJ).maxAbs ≤ 1 / δ := by have hinvord : D.hi⁻¹ ≤ D.lo⁻¹ := (inv_le_inv₀ hhi hlo).2 D.lo_le_hi simp only [RatInterval.inv, RatInterval.maxAbs, abs_of_pos (inv_pos.mpr hhi), abs_of_pos (inv_pos.mpr hlo)] rw [max_eq_right hinvord] simpa [one_div] using (inv_le_inv₀ hlo hδ).2 hsep have hN0 : 0 ≤ N.maxAbs := (abs_nonneg N.re.lo).trans ((le_max_left _ _).trans (le_max_left _ _)) have hinv0 : 0 ≤ (D.inv hJ).maxAbs := (abs_nonneg (D.inv hJ).lo).trans (le_max_left _ _) have hre := RatInterval.maxAbs_mul N.re (D.inv hJ) have him := RatInterval.maxAbs_mul N.im (D.inv hJ) have hre' : (N.re.mul (D.inv hJ)).maxAbs ≤ N.maxAbs / δ := by calc _ ≤ N.re.maxAbs * (D.inv hJ).maxAbs := hre _ ≤ N.maxAbs * (1 / δ) := mul_le_mul (le_max_left _ _) hinvA hinv0 hN0 _ = N.maxAbs / δ := by ring have him' : (N.im.mul (D.inv hJ)).maxAbs ≤ N.maxAbs / δ := by calc _ ≤ N.im.maxAbs * (D.inv hJ).maxAbs := him _ ≤ N.maxAbs * (1 / δ) := mul_le_mul (le_max_right _ _) hinvA hinv0 hN0 _ = N.maxAbs / δ := by ring change max (N.re.div D hJ).maxAbs (N.im.div D hJ).maxAbs ≤ N.maxAbs / δ change max (N.re.mul (D.inv hJ)).maxAbs (N.im.mul (D.inv hJ)).maxAbs ≤ N.maxAbs / δ exact max_le hre' him'

theorem div_width {I J : ComplexRatInterval} (hJ : J.normSq.AwayFromZero) {δ : ℚ} (hδ : 0 < δ) (hsep : δ ≤ J.normSq.lo) : (I.div J hJ).width ≤ (I.mul J.conj).maxAbs * J.normSq.width / (δ ^ 2) + (I.mul J.conj).width / δ := by let N := I.mul J.conj let D := J.normSq have hlo : 0 < D.lo := hδ.trans_le hsep have hhi : 0 < D.hi := hlo.trans_le D.lo_le_hi have hprod : δ ^ 2 ≤ D.lo * D.hi := by have := mul_le_mul hsep (hsep.trans D.lo_le_hi) hδ.le hlo.le simpa [pow_two] using this have hinvw : (D.inv hJ).width ≤ D.width / δ ^ 2 := by have heq : (D.inv hJ).width = D.width / (D.lo * D.hi) := by simp only [RatInterval.inv, RatInterval.width] field_simp rw [heq] exact div_le_div_of_nonneg_left (RatInterval.width_nonneg D) (sq_pos_of_pos hδ) hprod have hinvA : (D.inv hJ).maxAbs ≤ 1 / δ := by have hinvord : D.hi⁻¹ ≤ D.lo⁻¹ := (inv_le_inv₀ hhi hlo).2 D.lo_le_hi simp only [RatInterval.inv, RatInterval.maxAbs, abs_of_pos (inv_pos.mpr hhi), abs_of_pos (inv_pos.mpr hlo)] rw [max_eq_right hinvord] simpa [one_div] using (inv_le_inv₀ hlo hδ).2 hsep have hNA : 0 ≤ N.maxAbs := (abs_nonneg N.re.lo).trans ((le_max_left _ _).trans (le_max_left _ _)) have hNW : 0 ≤ N.width := (RatInterval.width_nonneg N.re).trans (le_max_left _ _) have hInvA0 : 0 ≤ (D.inv hJ).maxAbs := (abs_nonneg (D.inv hJ).lo).trans (le_max_left _ _) have hInvW0 : 0 ≤ (D.inv hJ).width := RatInterval.width_nonneg _ have hre := interval_mul_width_le N.re (D.inv hJ) have him := interval_mul_width_le N.im (D.inv hJ) have hreA : N.re.maxAbs ≤ N.maxAbs := le_max_left _ _ have himA : N.im.maxAbs ≤ N.maxAbs := le_max_right _ _ have hrew : N.re.width ≤ N.width := le_max_left _ _ have himw : N.im.width ≤ N.width := le_max_right _ _ have hre' : (N.re.mul (D.inv hJ)).width ≤ N.maxAbs * (D.width / δ ^ 2) + (1 / δ) * N.width := by nlinarith [mul_le_mul hreA hinvw hInvW0 hNA, mul_le_mul hinvA hrew (RatInterval.width_nonneg N.re) (by positivity : (0 : ℚ) ≤ 1 / δ)] have him' : (N.im.mul (D.inv hJ)).width ≤ N.maxAbs * (D.width / δ ^ 2) + (1 / δ) * N.width := by nlinarith [mul_le_mul himA hinvw hInvW0 hNA, mul_le_mul hinvA himw (RatInterval.width_nonneg N.im) (by positivity : (0 : ℚ) ≤ 1 / δ)] change max (N.re.div D hJ).width (N.im.div D hJ).width ≤ N.maxAbs * D.width / δ ^ 2 + N.width / δ change max (N.re.mul (D.inv hJ)).width (N.im.mul (D.inv hJ)).width ≤ _ apply max_le · exact hre'.trans_eq (by ring) · exact him'.trans_eq (by ring)

theorem normSq_lo_nonneg (I : ComplexRatInterval) : 0 ≤ I.normSq.lo := by have hsquare (K : RatInterval) : 0 ≤ K.sq.lo := by unfold RatInterval.sq split_ifs <;> positivity simp only [normSq, RatInterval.add] exact add_nonneg (hsquare I.re) (hsquare I.im)

theorem normInterval_sound {I : ComplexRatInterval} {z : ℂ} (hz : I.Contains z) (fuel : ℕ) : (I.normInterval fuel).Contains ‖z‖ := by have hs := RatInterval.sqrtInterval_sound (normSq_lo_nonneg I) (normSq_sound hz) fuel simpa [normInterval, Real.sqrt_sq (norm_nonneg z)] using hs

theorem expInterval_sound {I : RatInterval} {x : ℝ} (hx : I.Contains x) (fuel : ℕ) : (expInterval I fuel).Contains (exp x) := by induction fuel with | zero => exact expIntervalRaw_sound hx 0 | succ fuel ih => exact RatInterval.tighten_sound ih (expIntervalRaw_sound hx (fuel + 1))

theorem expInterval_nested (I : RatInterval) (fuel : ℕ) : (expInterval I (fuel + 1)).Subinterval (expInterval I fuel) := by have hmid : I.Contains (intervalMid I : ℝ) := by constructor <;> exact_mod_cast (by dsimp [intervalMid] linarith [I.lo_le_hi]) exact RatInterval.tighten_subinterval_left (expInterval_sound hmid fuel) (expIntervalRaw_sound hmid (fuel + 1))

theorem normApprox_spec (z : CertifiedComplex) (fuel : ℕ) : (normApprox z fuel).Contains ‖z.value‖ ∧ (normApprox z (fuel + 1)).Subinterval (normApprox z fuel) := by induction fuel with | zero => have h0 := ComplexRatInterval.normInterval_sound (z.contains 0) 0 have h1 := ComplexRatInterval.normInterval_sound (z.contains 1) 1 exact ⟨h0, RatInterval.tighten_subinterval_left h0 h1⟩ | succ fuel ih => have hnew := ComplexRatInterval.normInterval_sound (z.contains (fuel + 1)) (fuel + 1) have hcur : (normApprox z (fuel + 1)).Contains ‖z.value‖ := by rw [normApprox] exact RatInterval.tighten_sound ih.1 hnew have hnext := ComplexRatInterval.normInterval_sound (z.contains (fuel + 2)) (fuel + 2) constructor · exact hcur · rw [show fuel + 1 + 1 = fuel + 2 by omega, normApprox] exact RatInterval.tighten_subinterval_left hcur hnext

theorem complexExp_sound {I : ComplexRatInterval} {z : ℂ} (hz : I.Contains z) (fuel : ℕ) : (complexExp I fuel).Contains (exp z) := by induction fuel with | zero => exact complexExpRaw_sound hz 0 | succ fuel ih => exact (ComplexRatInterval.tighten_sound_left ih (complexExpRaw_sound hz (fuel + 1))).1

theorem complexExp_nested (I : ComplexRatInterval) (fuel : ℕ) : (complexExp I (fuel + 1)).Subinterval (complexExp I fuel) := by have hz : I.Contains ((intervalMid I.re : ℝ) + (intervalMid I.im : ℝ) * I) := by constructor · simpa using (show I.re.Contains (intervalMid I.re : ℝ) from by constructor <;> exact_mod_cast (by dsimp [intervalMid] linarith [I.re.lo_le_hi])) · simpa using (show I.im.Contains (intervalMid I.im : ℝ) from by constructor <;> exact_mod_cast (by dsimp [intervalMid] linarith [I.im.lo_le_hi])) exact (ComplexRatInterval.tighten_sound_left (complexExp_sound hz fuel) (complexExpRaw_sound hz (fuel + 1))).2

theorem complexExp_width (I : ComplexRatInterval) (fuel : ℕ) : (complexExp I fuel).width ≤ 2 * (max (expInterval I.re fuel).maxAbs 0 * max (cosInterval I.im fuel).width (sinInterval I.im fuel).width + max (cosInterval I.im fuel).maxAbs (sinInterval I.im fuel).maxAbs * max (expInterval I.re fuel).width 0) := by -- View `complexExpRaw` as multiplication of `⟨E, 0⟩` by `⟨cos, sin⟩`. -- First show the recursively tightened result is a subrectangle of the raw -- result at this fuel, then apply `ComplexRatInterval.mul_width` and simplify -- the point-zero coordinate widths and magnitudes. let E := expInterval I.re fuel let T : ComplexRatInterval := ⟨cosInterval I.im fuel, sinInterval I.im fuel⟩ let A : ComplexRatInterval := ⟨E, RatInterval.point 0⟩ have hsub : (complexExp I fuel).Subinterval (complexExpRaw I fuel) := by cases fuel with | zero => exact ⟨RatInterval.subinterval_refl _, RatInterval.subinterval_refl _⟩ | succ fuel => rw [complexExp] have hmid : I.Contains ((intervalMid I.re : ℝ) + (intervalMid I.im : ℝ) * I) := by constructor · simpa using (show I.re.Contains (intervalMid I.re : ℝ) from by constructor <;> exact_mod_cast (by dsimp [intervalMid] linarith [I.re.lo_le_hi])) · simpa using (show I.im.Contains (intervalMid I.im : ℝ) from by constructor <;> exact_mod_cast (by dsimp [intervalMid] linarith [I.im.lo_le_hi])) have hprev := complexExp_sound hmid fuel have hraw := complexExpRaw_sound hmid (fuel + 1) exact ⟨ Causalean.Mathlib.Analysis.CertifiedContourIntervalArithmetic.Transcendental.tighten_subinterval_right hprev.1 hraw.1, Causalean.Mathlib.Analysis.CertifiedContourIntervalArithmetic.Transcendental.tighten_subinterval_right hprev.2 hraw.2⟩ have hwidth : (complexExp I fuel).width ≤ (complexExpRaw I fuel).width := max_le_max (RatInterval.width_mono hsub.1) (RatInterval.width_mono hsub.2) have hraw_eq : complexExpRaw I fuel = A.mul T := by ext <;> simp [complexExpRaw, A, T, E, ComplexRatInterval.mul, RatInterval.point, RatInterval.mul, RatInterval.sub, RatInterval.add, RatInterval.neg] rw [hraw_eq] at hwidth refine hwidth.trans ((ComplexRatInterval.mul_width A T).trans_eq ?_) dsimp [A, T, E, ComplexRatInterval.maxAbs, ComplexRatInterval.width] simp only [RatInterval.point, RatInterval.maxAbs, RatInterval.width, abs_zero, max_self, sub_self]

theorem complexExpNameApprox_spec (z : CertifiedComplex) (n : ℕ) : (complexExpNameApprox z n).Contains (exp z.value) ∧ (complexExpNameApprox z (n + 1)).Subinterval (complexExpNameApprox z n) := by induction n with | zero => have h0 := complexExpStage_sound z 0 have h1 := complexExpStage_sound z 1 exact ⟨h0, (ComplexRatInterval.tighten_sound_left h0 h1).2⟩ | succ n ih => have hnew := complexExpStage_sound z (n + 1) have hcur : (complexExpNameApprox z (n + 1)).Contains (exp z.value) := by rw [complexExpNameApprox] exact (ComplexRatInterval.tighten_sound_left ih.1 hnew).1 have hnext := complexExpStage_sound z (n + 2) constructor · exact hcur · rw [show n + 1 + 1 = n + 2 by omega, complexExpNameApprox] exact (ComplexRatInterval.tighten_sound_left hcur hnext).2

theorem complexExpStage_width (z : CertifiedComplex) (n : ℕ) : (complexExpStage z n).width ≤ (complexExpStageTolerance n).1 := by -- Let η be the stage tolerance, B the initial exponential magnitude scale, -- τ = η/(16B), α = τ/8, and K the scheduled input. Nestedness of `z` -- gives `K ⊆ z.approx 0`, while `z.width_modulus` gives `width K ≤ τ/(4B)`. -- Use `power_div_factorial_le` for the sine/cosine errors and -- `expScalar_width` plus scalar nestedness for both exponential calls. -- Bound each factor's `maxAbs` from one contained exact value plus its width, -- feed those estimates to `complexExp_width`, and finish the deliberately -- loose constants by `nlinarith`/`field_simp` using η ≤ 1 and B ≥ 2. let η : ℚ := (complexExpStageTolerance n).1 let B : ℚ := complexExpMagnitude z let τ : ℚ := (complexExpInnerTolerance z n).1 let α : PosRat := ⟨τ / 8, by change 0 < (complexExpInnerTolerance z n).1 / 8 exact div_pos (complexExpInnerTolerance z n).2 (by norm_num)⟩ let K : ComplexRatInterval := complexExpStageInput z n let N : ℕ := complexExpStageFuel z n let d : ℕ := α.1.den + 1 let A : ℕ := (intervalMid K.im).num.natAbs + 2 let S0 := expScalar (z.approx 0).re.maxAbs 0 let E := expInterval K.re N let C := cosInterval K.im N let S := sinInterval K.im N have hηpos : 0 < η := by dsimp [η, complexExpStageTolerance] positivity have hηone : η ≤ 1 := by dsimp [η, complexExpStageTolerance] apply (div_le_iff₀ (by positivity : (0 : ℚ) < (n : ℚ) + 1)).2 norm_num have hB : 2 ≤ B := by dsimp [B, complexExpMagnitude, S0] linarith [abs_nonneg ((expScalar (z.approx 0).re.maxAbs 0).hi)] have hBpos : 0 < B := lt_of_lt_of_le (by norm_num) hB have hτeq : τ = η / (16 * B) := by rfl have hαeq : α.1 = τ / 8 := by rfl have hτpos : 0 < τ := by rw [hτeq]; positivity have hτone : τ ≤ 1 := by rw [hτeq] apply (div_le_iff₀ (by positivity : (0 : ℚ) < 16 * B)).2 nlinarith [hηone] have hαpos : 0 < α.1 := α.2 have hαone : α.1 ≤ 1 := by rw [hαeq]; nlinarith [hτone] have hK0 : K.Subinterval (z.approx 0) := by dsimp [K, complexExpStageInput] exact complexApprox_subinterval_of_le z (zero_le _) have hKw : K.width ≤ τ / (4 * B) := by have hm := z.width_modulus (complexExpInputTolerance z n) simpa [K, complexExpStageInput, complexExpInputTolerance, τ] using hm have hKrew : K.re.width ≤ τ / (4 * B) := (le_max_left _ _).trans hKw have hKimw : K.im.width ≤ τ / (4 * B) := (le_max_right _ _).trans hKw have hNtrig : 32 * d * A ^ 2 ≤ N := by dsimp [N, complexExpStageFuel, d, A, α, K] exact le_max_left _ _ have hNexpMid : expPrecision (intervalMid K.re) α ≤ N := by dsimp [N, complexExpStageFuel, K, α] exact (le_max_left _ _).trans (le_max_right _ _) have hNexpAbs : expPrecision K.re.maxAbs α ≤ N := by dsimp [N, complexExpStageFuel, K, α] exact (le_max_right _ _).trans (le_max_right _ _) have hAone : 1 ≤ A := by dsimp [A]; omega have hdone : 1 ≤ d := by dsimp [d]; omega have hmidA : |intervalMid K.im| ≤ (A : ℚ) := by have hdq : (0 : ℚ) < (intervalMid K.im).den := by exact_mod_cast (intervalMid K.im).den_pos have hnum : |intervalMid K.im| ≤ ((intervalMid K.im).num.natAbs : ℚ) := by calc |intervalMid K.im| = |((intervalMid K.im).num : ℚ) / ((intervalMid K.im).den : ℚ)| := by rw [Rat.num_div_den] _ = ((intervalMid K.im).num.natAbs : ℚ) / (intervalMid K.im).den := by rw [abs_div, abs_of_pos hdq] norm_num _ ≤ ((intervalMid K.im).num.natAbs : ℕ) := by apply (div_le_iff₀ hdq).2 have hd1 : (1 : ℚ) ≤ (intervalMid K.im).den := by exact_mod_cast (intervalMid K.im).den_pos nlinarith calc |intervalMid K.im| ≤ ((intervalMid K.im).num.natAbs : ℕ) := hnum _ ≤ A := by exact_mod_cast le_add_right (intervalMid K.im).num.natAbs 2 have hsinErr : sinError (intervalMid K.im) N ≤ α.1 / 4 := by have hp := power_div_factorial_le (intervalMid K.im) A d N 3 hmidA hAone hdone (by omega) hNtrig have hinv : 1 / (d : ℚ) ≤ α.1 := by exact (div_le_div_of_nonneg_left (by norm_num) (by positivity) (by dsimp [d]; exact_mod_cast Nat.le_succ α.1.den)).trans (inv_den_le_of_pos α.1 α.2) calc _ ≤ 1 / ((4 * d : ℕ) : ℚ) := by simpa [sinError] using hp _ = (1 / (d : ℚ)) / 4 := by push_cast; field_simp _ ≤ α.1 / 4 := div_le_div_of_nonneg_right hinv (by norm_num) have hcosErr : cosError (intervalMid K.im) N ≤ α.1 / 4 := by have hp := power_div_factorial_le (intervalMid K.im) A d N 2 hmidA hAone hdone (by omega) hNtrig have hinv : 1 / (d : ℚ) ≤ α.1 := by exact (div_le_div_of_nonneg_left (by norm_num) (by positivity) (by dsimp [d]; exact_mod_cast Nat.le_succ α.1.den)).trans (inv_den_le_of_pos α.1 α.2) calc _ ≤ 1 / ((4 * d : ℕ) : ℚ) := by simpa [cosError] using hp _ = (1 / (d : ℚ)) / 4 := by push_cast; field_simp _ ≤ α.1 / 4 := div_le_div_of_nonneg_right hinv (by norm_num) have hCw : C.width ≤ τ := by have h := cosInterval_width K.im N dsimp [C] rw [hαeq] at hcosErr have hdiv : τ / (4 * B) ≤ τ / 8 := by apply div_le_div_of_nonneg_left hτpos.le (by norm_num) nlinarith nlinarith have hSw : S.width ≤ τ := by have h := sinInterval_width K.im N dsimp [S] rw [hαeq] at hsinErr have hdiv : τ / (4 * B) ≤ τ / 8 := by apply div_le_div_of_nonneg_left hτpos.le (by norm_num) nlinarith nlinarith have hscalarMid : (expScalar (intervalMid K.re) N).width ≤ α.1 := by exact (RatInterval.width_mono (expScalar_subinterval_of_le (intervalMid K.re) hNexpMid)).trans (expScalar_width (intervalMid K.re) α) have hscalarAbs : (expScalar K.re.maxAbs N).width ≤ α.1 := by exact (RatInterval.width_mono (expScalar_subinterval_of_le K.re.maxAbs hNexpAbs)).trans (expScalar_width K.re.maxAbs α) have hKabs : K.re.maxAbs ≤ (z.approx 0).re.maxAbs := ratMaxAbs_mono hK0.1 have hmidK : K.re.Contains (intervalMid K.re : ℝ) := by constructor <;> exact_mod_cast (by dsimp [intervalMid] linarith [K.re.lo_le_hi]) have hmidM : (intervalMid K.re : ℝ) ≤ (z.approx 0).re.maxAbs := by exact_mod_cast (le_trans (le_abs_self _) ((abs_le_max_abs_abs (by exact_mod_cast hmidK.1) (by exact_mod_cast hmidK.2)).trans hKabs)) have hvalM : z.value.re ≤ ((z.approx 0).re.maxAbs : ℚ) := by have hz0 := (z.contains 0).1 exact (le_abs_self _).trans (by exact_mod_cast abs_le_max_abs_abs hz0.1 hz0.2) have hS0 := Causalean.Mathlib.Analysis.CertifiedContourIntervalArithmetic.Transcendental.expScalar_sound (z.approx 0).re.maxAbs 0 have hS0hi0 : 0 ≤ S0.hi := by dsimp [S0] exact_mod_cast (Real.exp_pos ((z.approx 0).re.maxAbs : ℝ)).le.trans hS0.2 have hUpper : (expScalar K.re.maxAbs N).hi ≤ B := by change (expScalar K.re.maxAbs N).hi ≤ |S0.hi| + 2 have hUs := Causalean.Mathlib.Analysis.CertifiedContourIntervalArithmetic.Transcendental.expScalar_sound K.re.maxAbs N have hexp : exp (K.re.maxAbs : ℝ) ≤ exp ((z.approx 0).re.maxAbs : ℝ) := Real.exp_le_exp.mpr (by exact_mod_cast hKabs) have hUhi : ((expScalar K.re.maxAbs N).hi : ℝ) ≤ exp (K.re.maxAbs : ℝ) + α.1 := by have hw : (((expScalar K.re.maxAbs N).hi : ℝ) - (expScalar K.re.maxAbs N).lo) ≤ α.1 := by exact_mod_cast hscalarAbs nlinarith [hUs.1] have hUS0 : ((expScalar K.re.maxAbs N).hi : ℝ) ≤ (S0.hi : ℝ) + α.1 := by nlinarith [hS0.2] have hUS0q : (expScalar K.re.maxAbs N).hi ≤ S0.hi + α.1 := by exact_mod_cast hUS0 exact hUS0q.trans (by rw [abs_of_nonneg hS0hi0] nlinarith [hαone]) have hEw : E.width ≤ τ := by have h := expInterval_width K.re N dsimp [E] rw [show intervalRadius K.re = K.re.width / 2 by rfl] at h rw [hαeq] at hscalarMid have hprod : (expScalar K.re.maxAbs N).hi * K.re.width ≤ τ / 4 := by have hu0 : 0 ≤ (expScalar K.re.maxAbs N).hi := by exact_mod_cast (Real.exp_pos (K.re.maxAbs : ℝ)).le.trans (Causalean.Mathlib.Analysis.CertifiedContourIntervalArithmetic.Transcendental.expScalar_sound K.re.maxAbs N).2 have := mul_le_mul hUpper hKrew (RatInterval.width_nonneg K.re) hBpos.le calc _ ≤ B * (τ / (4 * B)) := this _ = τ / 4 := by field_simp nlinarith have hEcontains : E.Contains (exp z.value.re) := by exact expInterval_sound (z.contains _).1 N have hEabs : E.maxAbs ≤ B := by have hexp : exp z.value.re ≤ (S0.hi : ℝ) := (Real.exp_le_exp.mpr hvalM).trans hS0.2 have h := ratMaxAbs_le_of_contains hEcontains (by rw [Real.abs_exp]; exact hexp.trans (by exact_mod_cast le_abs_self S0.hi)) hEw (abs_nonneg S0.hi) hτpos.le change E.maxAbs ≤ |S0.hi| + 2 exact h.trans (by nlinarith [hτone]) have hCcontains : C.Contains (cos z.value.im) := cosInterval_sound (z.contains _).2 N have hScontains : S.Contains (sin z.value.im) := sinInterval_sound (z.contains _).2 N have hCabs : C.maxAbs ≤ 2 := by exact (ratMaxAbs_le_of_contains (C := (1 : ℚ)) (w := τ) hCcontains (by simpa using Real.abs_cos_le_one z.value.im) hCw (by norm_num) hτpos.le).trans (by nlinarith [hτone]) have hSabs : S.maxAbs ≤ 2 := by exact (ratMaxAbs_le_of_contains (C := (1 : ℚ)) (w := τ) hScontains (by simpa using Real.abs_sin_le_one z.value.im) hSw (by norm_num) hτpos.le).trans (by nlinarith [hτone]) have hout := complexExp_width K N change (complexExp K N).width ≤ η have hEA0 : 0 ≤ E.maxAbs := (abs_nonneg E.lo).trans (le_max_left _ _) have hCA0 : 0 ≤ C.maxAbs := (abs_nonneg C.lo).trans (le_max_left _ _) have hSA0 : 0 ≤ S.maxAbs := (abs_nonneg S.lo).trans (le_max_left _ _) have hEw0 := RatInterval.width_nonneg E have hCw0 := RatInterval.width_nonneg C have hSw0 := RatInterval.width_nonneg S have hout' : (complexExp K N).width ≤ 2 * (E.maxAbs * max C.width S.width + -- … truncated; follow the source link for the rest …

theorem circleNode_width_at_selected_precision (radius : ℚ) (hradius : 0 ≤ radius) (schedule : Schedule) (target : PosRat) (hprecision : schedule.inputPrecision = circleInputPrecision radius target) (hfuel : circleExpFuel radius schedule.mesh target ≤ schedule.fuel) {k : ℕ} (hk : k ≤ schedule.mesh) : (circleNode radius schedule k).width ≤ target.1 := by let τ := circleInnerTolerance radius target let c : ℚ := 2 * k / schedule.mesh let P := Transcendental.piInterval schedule.inputPrecision let K := (circleAngle schedule k).im let E := Transcendental.expInterval (circleAngle schedule k).re schedule.fuel let C := Transcendental.cosInterval K schedule.fuel let S := Transcendental.sinInterval K schedule.fuel have hτ0 : 0 < τ.1 := τ.2 have hτone : τ.1 ≤ 1 / 1024 := by exact min_le_right _ _ have hτtarget : τ.1 ≤ target.1 / (256 * (|radius| + 1)) := by exact min_le_left _ _ have hc0 : 0 ≤ c := by dsimp [c] positivity have hcle : c ≤ 2 := by apply (div_le_iff₀ (show (0 : ℚ) < schedule.mesh by exact_mod_cast schedule.mesh_pos)).2 push_cast exact_mod_cast (Nat.mul_le_mul_left 2 hk) have hPw : P.width ≤ τ.1 := by dsimp [P] rw [hprecision] exact Transcendental.piInterval_width τ have hKw : K.width ≤ 2 * τ.1 := by change ((RatInterval.point c).mul P).width ≤ 2 * τ.1 rw [point_mul_width c hc0 P] calc c * P.width ≤ 2 * P.width := mul_le_mul_of_nonneg_right hcle (RatInterval.width_nonneg P) _ ≤ 2 * τ.1 := mul_le_mul_of_nonneg_left hPw (by norm_num) have hP0 : P.Subinterval (Transcendental.piInterval 0) := by exact CertifiedReal.approx_mono Transcendental.piName (zero_le _) have hP0w : (Transcendental.piInterval 0).width ≤ 1 := by norm_num [Transcendental.piInterval, Transcendental.piRaw, Transcendental.atanRaw, Transcendental.atanPartial, Transcendental.atanError, RatInterval.mul, RatInterval.sub, RatInterval.add, RatInterval.neg, RatInterval.point, RatInterval.width] have hPwide : P.width ≤ 1 := (RatInterval.width_mono hP0).trans hP0w have hPabs : P.maxAbs ≤ 5 := by have h := interval_maxAbs_le_of_contains_width (C := 4) (w := 1) (Transcendental.piInterval_sound schedule.inputPrecision) (by simpa [abs_of_pos Real.pi_pos] using Real.pi_le_four) hPwide convert h using 1 <;> norm_num [P] have hKabs : K.maxAbs ≤ 10 := by change ((RatInterval.point c).mul P).maxAbs ≤ 10 rw [point_mul_maxAbs c hc0 P] nlinarith have hmidK : K.Contains (Transcendental.intervalMid K : ℝ) := by constructor <;> exact_mod_cast (by dsimp [Transcendental.intervalMid] linarith [K.lo_le_hi]) have hmidA : |Transcendental.intervalMid K| ≤ (12 : ℚ) := by exact (abs_le_max_abs_abs (by exact_mod_cast hmidK.1) (by exact_mod_cast hmidK.2)).trans (hKabs.trans (by norm_num)) have hfuelExp : expPrecision 0 τ ≤ schedule.fuel := by exact (le_max_left _ _).trans (by simpa [circleExpFuel, τ] using hfuel) have hfuelTrig : 32 * (τ.1.den + 1) * 12 ^ 2 ≤ schedule.fuel := by exact (le_max_right _ _).trans (by simpa [circleExpFuel, τ] using hfuel) have hEw : E.width ≤ τ.1 := by have hraw := Transcendental.expInterval_width (RatInterval.point 0) schedule.fuel have hscalar : (expScalar 0 schedule.fuel).width ≤ τ.1 := by exact (RatInterval.width_mono (CertifiedReal.approx_mono (expName 0) hfuelExp)).trans (expScalar_width 0 τ) change (Transcendental.expInterval (RatInterval.point 0) schedule.fuel).width ≤ τ.1 exact hraw.trans (by simpa [RatInterval.point, Transcendental.intervalMid, Transcendental.intervalRadius, RatInterval.width] using hscalar) have hsinErr : Transcendental.sinError (Transcendental.intervalMid K) schedule.fuel ≤ τ.1 / 4 := by have herr := Transcendental.power_div_factorial_le (Transcendental.intervalMid K) 12 (τ.1.den + 1) schedule.fuel 3 hmidA (by norm_num) (by omega) (by omega) hfuelTrig have hinv : 1 / ((τ.1.den + 1 : ℕ) : ℚ) ≤ τ.1 := by exact (div_le_div_of_nonneg_left (by norm_num) (by positivity) (by exact_mod_cast Nat.le_succ τ.1.den)).trans (inv_den_le_of_pos τ.1 τ.2) calc _ ≤ 1 / ((4 * (τ.1.den + 1) : ℕ) : ℚ) := by simpa [Transcendental.sinError] using herr _ = (1 / ((τ.1.den + 1 : ℕ) : ℚ)) / 4 := by push_cast; field_simp _ ≤ τ.1 / 4 := div_le_div_of_nonneg_right hinv (by norm_num) have hcosErr : Transcendental.cosError (Transcendental.intervalMid K) schedule.fuel ≤ τ.1 / 4 := by have herr := Transcendental.power_div_factorial_le (Transcendental.intervalMid K) 12 (τ.1.den + 1) schedule.fuel 2 hmidA (by norm_num) (by omega) (by omega) hfuelTrig have hinv : 1 / ((τ.1.den + 1 : ℕ) : ℚ) ≤ τ.1 := by exact (div_le_div_of_nonneg_left (by norm_num) (by positivity) (by exact_mod_cast Nat.le_succ τ.1.den)).trans (inv_den_le_of_pos τ.1 τ.2) calc _ ≤ 1 / ((4 * (τ.1.den + 1) : ℕ) : ℚ) := by simpa [Transcendental.cosError] using herr _ = (1 / ((τ.1.den + 1 : ℕ) : ℚ)) / 4 := by push_cast; field_simp _ ≤ τ.1 / 4 := div_le_div_of_nonneg_right hinv (by norm_num) have hCw : C.width ≤ 5 * τ.1 / 2 := (Transcendental.cosInterval_width K schedule.fuel).trans (by nlinarith) have hSw : S.width ≤ 5 * τ.1 / 2 := (Transcendental.sinInterval_width K schedule.fuel).trans (by nlinarith) have hEcontains : E.Contains (1 : ℝ) := by change (Transcendental.expInterval (RatInterval.point 0) schedule.fuel).Contains 1 simpa using Transcendental.expInterval_sound (RatInterval.point_sound (0 : ℚ)) schedule.fuel have hEabs : E.maxAbs ≤ 2 := by exact (interval_maxAbs_le_of_contains_width (C := 1) (w := τ.1) hEcontains (by norm_num) hEw).trans (by nlinarith) have hCcontains : C.Contains (cos (Transcendental.intervalMid K : ℝ)) := Transcendental.cosInterval_sound hmidK schedule.fuel have hScontains : S.Contains (sin (Transcendental.intervalMid K : ℝ)) := Transcendental.sinInterval_sound hmidK schedule.fuel have hCabs : C.maxAbs ≤ 2 := by exact (interval_maxAbs_le_of_contains_width (C := 1) (w := 5 * τ.1 / 2) hCcontains (by simpa only [Rat.cast_one] using Real.abs_cos_le_one (Transcendental.intervalMid K : ℝ)) hCw).trans (by nlinarith) have hSabs : S.maxAbs ≤ 2 := by exact (interval_maxAbs_le_of_contains_width (C := 1) (w := 5 * τ.1 / 2) hScontains (by simpa only [Rat.cast_one] using Real.abs_sin_le_one (Transcendental.intervalMid K : ℝ)) hSw).trans (by nlinarith) have hout := Transcendental.complexExp_width (circleAngle schedule k) schedule.fuel have hout' : (Transcendental.complexExp (circleAngle schedule k) schedule.fuel).width ≤ 2 * (max E.maxAbs 0 * max C.width S.width + max C.maxAbs S.maxAbs * max E.width 0) := by simpa only [E, C, S, K] using hout have hcomplex : (Transcendental.complexExp (circleAngle schedule k) schedule.fuel).width ≤ 14 * τ.1 := by change _ ≤ 14 * τ.1 have hTW : max C.width S.width ≤ 5 * τ.1 / 2 := max_le hCw hSw have hTA : max C.maxAbs S.maxAbs ≤ 2 := max_le hCabs hSabs have hEA : max E.maxAbs 0 ≤ 2 := max_le hEabs (by norm_num) have hEW : max E.width 0 ≤ τ.1 := max_le hEw hτ0.le have hTW0 : 0 ≤ max C.width S.width := (RatInterval.width_nonneg C).trans (le_max_left _ _) have hEW0 : 0 ≤ max E.width 0 := le_max_right _ _ have hprod1 := mul_le_mul hEA hTW hTW0 (by norm_num : (0 : ℚ) ≤ 2) have hprod2 := mul_le_mul hTA hEW hEW0 (by norm_num : (0 : ℚ) ≤ 2) exact hout'.trans (by nlinarith) rw [circleNode, complex_smul_width radius hradius] rw [abs_of_nonneg hradius] at hτtarget have hden : 0 < 256 * (radius + 1) := by positivity have hprod : 256 * (radius + 1) * τ.1 ≤ target.1 := by simpa [mul_comm, mul_left_comm] using (le_div_iff₀ hden).mp hτtarget calc radius * (Transcendental.complexExp (circleAngle schedule k) schedule.fuel).width ≤ radius * (14 * τ.1) := mul_le_mul_of_nonneg_left hcomplex hradius _ ≤ target.1 := by nlinarith

theorem circleAngle_sound (schedule : Schedule) {k : ℕ} (hk : k ≤ schedule.mesh) : (circleAngle schedule k).Contains (((2 : ℝ) * pi * ((k : ℝ) / schedule.mesh)) * I) := by constructor · simpa [circleAngle] using RatInterval.point_sound (0 : ℚ) · have h := RatInterval.mul_sound (RatInterval.point_sound (2 * k / schedule.mesh : ℚ)) (Transcendental.piInterval_sound schedule.inputPrecision) have heq : (((2 * k / schedule.mesh : ℚ) : ℝ) * pi) = ((((2 : ℝ) * pi * ((k : ℝ) / schedule.mesh)) * I).im) := by push_cast simp ring rw [← heq] exact h

theorem circleNode_endpoint (radius : ℚ) (schedule : Schedule) : (circleNode radius schedule schedule.mesh).Contains (radius : ℂ) ∧ exactCircleNode radius schedule schedule.mesh = (radius : ℂ) := by have hratio : ((schedule.mesh : ℝ) : ℂ) / schedule.mesh = 1 := by apply div_self exact_mod_cast schedule.mesh_pos.ne' have hexact : exactCircleNode radius schedule schedule.mesh = (radius : ℂ) := by rw [exactCircleNode, hratio] simp [Complex.exp_two_pi_mul_I] constructor · rw [← hexact] exact circleNode_sound radius schedule (le_refl schedule.mesh) · exact hexact

theorem circleNode_width (radius : ℚ) (hradius : 0 ≤ radius) (schedule : Schedule) (k : ℕ) : (circleNode radius schedule k).width ≤ radius * (Transcendental.complexExp (circleAngle schedule k) schedule.fuel).width := by let I := Transcendental.complexExp (circleAngle schedule k) schedule.fuel have hre : ((RatInterval.point radius).mul I.re).width = radius * I.re.width := by have h := mul_le_mul_of_nonneg_left I.re.lo_le_hi hradius simp [RatInterval.mul, RatInterval.point, RatInterval.width, min_eq_left h, max_eq_right h] ring have him : ((RatInterval.point radius).mul I.im).width = radius * I.im.width := by have h := mul_le_mul_of_nonneg_left I.im.lo_le_hi hradius simp [RatInterval.mul, RatInterval.point, RatInterval.width, min_eq_left h, max_eq_right h] ring change max ((RatInterval.point radius).mul I.re).width ((RatInterval.point radius).mul I.im).width ≤ radius * max I.re.width I.im.width rw [hre, him, ← mul_max_of_nonneg I.re.width I.im.width hradius]

theorem circleNode_schedule_correspondence (radius : ℚ) (schedule : Schedule) (k operation : ℕ) : (circleNodeEvent schedule k operation).schedule.operationCount = schedule.operationCount ∧ (circleNodeEvent schedule k operation).schedule.inputPrecision = schedule.inputPrecision ∧ (circleNodeEvent schedule k operation).schedule.fuel = schedule.fuel ∧ (circleNodeEvent schedule k operation).schedule.mesh = schedule.mesh ∧ circleNode radius (circleNodeEvent schedule k operation).schedule k = circleNode radius schedule k := by exact ⟨rfl, rfl, rfl, rfl, rfl⟩

theorem normalizedContourIntegral_eq (program : ContourProgram) (schedule : Schedule) (hden : ∀ u ∈ Icc (0 : ℝ) 1, program.denominator.value ((program.radius : ℂ) * exp (((2 : ℝ) * Real.pi * u) * I)) ≠ 0) : program.normalizedContourIntegral = ∫ u in (0 : ℝ)..1, program.normalizedIntegrand schedule u := by rw [ContourProgram.normalizedContourIntegral, CircleMesh.circleContourIntegral] calc (∫ u in (0 : ℝ)..1, CircleMesh.circleIntegrand (fun z => program.numerator.value z / program.denominator.value z) 0 program.radius u) / (((2 : ℝ) * Real.pi : ℂ) * I) = ∫ u in (0 : ℝ)..1, CircleMesh.circleIntegrand (fun z => program.numerator.value z / program.denominator.value z) 0 program.radius u / (((2 : ℝ) * Real.pi : ℂ) * I) := (intervalIntegral.integral_div _ _).symm _ = ∫ u in (0 : ℝ)..1, program.normalizedIntegrand schedule u := by apply intervalIntegral.integral_congr intro u hu simp only [CircleMesh.circleIntegrand, CircleMesh.circleMap, CircleMesh.circleTangent, ContourProgram.normalizedIntegrand, zero_add] have htwoPi : (((2 : ℝ) * Real.pi : ℂ) * I) ≠ 0 := by exact mul_ne_zero (mul_ne_zero (by norm_num) (Complex.ofReal_ne_zero.mpr Real.pi_ne_zero)) Complex.I_ne_zero have huIcc : u ∈ Icc (0 : ℝ) 1 := by simpa [uIcc_of_le (by norm_num : (0 : ℝ) ≤ 1)] using hu have hdenu := hden u huIcc field_simp [htwoPi, hdenu] norm_num ring_nf

theorem integrandNodes_sound (program : ContourProgram) (schedule : Schedule) (certificate : DenominatorCertificate program schedule) {k : ℕ} (hk : k ≤ schedule.mesh) : (program.integrandNodes schedule certificate k).Contains (program.normalizedIntegrand schedule (CircleMesh.meshPoint schedule.mesh k)) := by have hnode : (program.nodeBox schedule k).Contains (exactCircleNode program.radius schedule k) := by exact circleNode_sound program.radius schedule hk have hnum := program.numerator.sound hnode schedule.fuel have hden := program.denominator.sound hnode schedule.fuel have hquot := ComplexRatInterval.div_sound (certificate.away k hk) hnum hden have hmul := ComplexRatInterval.mul_sound hquot hnode simpa [ContourProgram.integrandNodes, hk, ContourProgram.nodeBox, ContourProgram.normalizedIntegrand, exactCircleNode, CircleMesh.meshPoint] using hmul

theorem integrandNodes_width_le_propagation (program : ContourProgram) (bounds : ContourValueBounds program) (separation : PosRat) (scheduled : CertifiedProgramSchedule program bounds separation) (certificate : DenominatorCertificate program scheduled.schedule) (hseparation : separation.1 ≤ certificate.separation) {k : ℕ} (hk : k ≤ scheduled.schedule.mesh) : (program.integrandNodes scheduled.schedule certificate k).width ≤ program.nodePropagationBound bounds separation.1 scheduled.target.1 := by let schedule := scheduled.schedule let node := program.nodeBox schedule k let numerator := program.numerator.eval node schedule.fuel let denominator := program.denominator.eval node schedule.fuel let numeratorWidth := scheduled.target.1 + program.numerator.amplification * scheduled.target.1 let denominatorWidth := scheduled.target.1 + program.denominator.amplification * scheduled.target.1 let numeratorMax := bounds.numerator + numeratorWidth let denominatorMax := bounds.denominator + denominatorWidth let divisionMax := 2 * numeratorMax * denominatorMax let divisionWidth := 2 * (numeratorMax * denominatorWidth + denominatorMax * numeratorWidth) let normWidth := 4 * denominatorMax * denominatorWidth let quotientWidth := divisionMax * normWidth / separation.1 ^ 2 + divisionWidth / separation.1 have hnodeWidth : node.width ≤ scheduled.target.1 := by exact circleNode_width_at_selected_precision program.radius program.radius_pos.le schedule scheduled.target scheduled.inputPrecision_eq scheduled.circleFuel_le hk have hnodeContains : node.Contains (exactCircleNode program.radius schedule k) := by exact circleNode_sound program.radius schedule hk have hnumContains : numerator.Contains (program.numerator.value (exactCircleNode program.radius schedule k)) := by exact program.numerator.sound hnodeContains schedule.fuel have hdenContains : denominator.Contains (program.denominator.value (exactCircleNode program.radius schedule k)) := by exact program.denominator.sound hnodeContains schedule.fuel have hnumAlgorithm : program.numerator.algorithmError schedule.fuel ≤ scheduled.target.1 := by exact CertifiedComplexMap.algorithmError_le_of_modulus scheduled.numeratorFuel_le have hdenAlgorithm : program.denominator.algorithmError schedule.fuel ≤ scheduled.target.1 := by exact CertifiedComplexMap.algorithmError_le_of_modulus scheduled.denominatorFuel_le have hnumWidth : numerator.width ≤ numeratorWidth := by calc numerator.width ≤ program.numerator.algorithmError schedule.fuel + program.numerator.amplification * node.width := program.numerator.width_le node schedule.fuel _ ≤ scheduled.target.1 + program.numerator.amplification * scheduled.target.1 := by exact add_le_add hnumAlgorithm (mul_le_mul_of_nonneg_left hnodeWidth program.numerator.amplification_nonneg) _ = numeratorWidth := rfl have hdenWidth : denominator.width ≤ denominatorWidth := by calc denominator.width ≤ program.denominator.algorithmError schedule.fuel + program.denominator.amplification * node.width := program.denominator.width_le node schedule.fuel _ ≤ scheduled.target.1 + program.denominator.amplification * scheduled.target.1 := by exact add_le_add hdenAlgorithm (mul_le_mul_of_nonneg_left hnodeWidth program.denominator.amplification_nonneg) _ = denominatorWidth := rfl have hmesh : CircleMesh.meshPoint schedule.mesh k ∈ Icc (0 : ℝ) 1 := CircleMesh.meshPoint_mem schedule.mesh_pos hk have hnumValueBound : max |(program.numerator.value (exactCircleNode program.radius schedule k)).re| |(program.numerator.value (exactCircleNode program.radius schedule k)).im| ≤ (bounds.numerator : ℝ) := by simpa [exactCircleNode, CircleMesh.meshPoint] using bounds.numerator_bound (CircleMesh.meshPoint schedule.mesh k) hmesh have hdenValueBound : max |(program.denominator.value (exactCircleNode program.radius schedule k)).re| |(program.denominator.value (exactCircleNode program.radius schedule k)).im| ≤ (bounds.denominator : ℝ) := by simpa [exactCircleNode, CircleMesh.meshPoint] using bounds.denominator_bound (CircleMesh.meshPoint schedule.mesh k) hmesh have hnumMax : numerator.maxAbs ≤ numeratorMax := by exact (ComplexRatInterval.maxAbs_le_of_contains_width hnumContains hnumValueBound hnumWidth) have hdenMax : denominator.maxAbs ≤ denominatorMax := by exact (ComplexRatInterval.maxAbs_le_of_contains_width hdenContains hdenValueBound hdenWidth) have htarget0 : 0 ≤ scheduled.target.1 := scheduled.target.2.le have hnWidth0 : 0 ≤ numeratorWidth := by dsimp [numeratorWidth] nlinarith [mul_nonneg program.numerator.amplification_nonneg htarget0] have hdWidth0 : 0 ≤ denominatorWidth := by dsimp [denominatorWidth] nlinarith [mul_nonneg program.denominator.amplification_nonneg htarget0] have hnMax0 : 0 ≤ numeratorMax := by dsimp [numeratorMax] linarith [bounds.numerator_nonneg] have hdMax0 : 0 ≤ denominatorMax := by dsimp [denominatorMax] linarith [bounds.denominator_nonneg] have hdivisionMax0 : 0 ≤ divisionMax := by dsimp [divisionMax] positivity have hdivisionWidth0 : 0 ≤ divisionWidth := by dsimp [divisionWidth] positivity have hnormWidth0 : 0 ≤ normWidth := by dsimp [normWidth] positivity have hnumDenMax : (numerator.mul denominator.conj).maxAbs ≤ divisionMax := by calc (numerator.mul denominator.conj).maxAbs ≤ 2 * numerator.maxAbs * denominator.conj.maxAbs := ComplexRatInterval.mul_maxAbs numerator denominator.conj _ = 2 * numerator.maxAbs * denominator.maxAbs := by rw [ComplexRatInterval.maxAbs_conj] _ ≤ 2 * numeratorMax * denominatorMax := by nlinarith [mul_le_mul hnumMax hdenMax ((abs_nonneg denominator.re.lo).trans ((le_max_left _ _).trans (le_max_left _ _))) hnMax0] _ = divisionMax := rfl have hnumDenWidth : (numerator.mul denominator.conj).width ≤ divisionWidth := by calc (numerator.mul denominator.conj).width ≤ 2 * (numerator.maxAbs * denominator.conj.width + denominator.conj.maxAbs * numerator.width) := ComplexRatInterval.mul_width numerator denominator.conj _ = 2 * (numerator.maxAbs * denominator.width + denominator.maxAbs * numerator.width) := by rw [ComplexRatInterval.width_conj, ComplexRatInterval.maxAbs_conj] _ ≤ 2 * (numeratorMax * denominatorWidth + denominatorMax * numeratorWidth) := by have hnumA0 : 0 ≤ numerator.maxAbs := (abs_nonneg numerator.re.lo).trans ((le_max_left _ _).trans (le_max_left _ _)) have hdenA0 : 0 ≤ denominator.maxAbs := (abs_nonneg denominator.re.lo).trans ((le_max_left _ _).trans (le_max_left _ _)) nlinarith [mul_le_mul hnumMax hdenWidth (show 0 ≤ denominator.width from (RatInterval.width_nonneg denominator.re).trans (le_max_left _ _)) hnMax0, mul_le_mul hdenMax hnumWidth (show 0 ≤ numerator.width from (RatInterval.width_nonneg numerator.re).trans (le_max_left _ _)) hdMax0] _ = divisionWidth := rfl have hnormWidth : denominator.normSq.width ≤ normWidth := by calc denominator.normSq.width ≤ 4 * denominator.maxAbs * denominator.width := ComplexRatInterval.normSq_width denominator _ ≤ 4 * denominatorMax * denominatorWidth := by nlinarith [mul_le_mul hdenMax hdenWidth (show 0 ≤ denominator.width from (RatInterval.width_nonneg denominator.re).trans (le_max_left _ _)) hdMax0] _ = normWidth := rfl have hsep : separation.1 ≤ denominator.normSq.lo := hseparation.trans (certificate.lower k hk) have hquotientWidth : (numerator.div denominator (certificate.away k hk)).width ≤ quotientWidth := by calc _ ≤ (numerator.mul denominator.conj).maxAbs * denominator.normSq.width / separation.1 ^ 2 + (numerator.mul denominator.conj).width / separation.1 := ComplexRatInterval.div_width (certificate.away k hk) separation.2 hsep _ ≤ divisionMax * normWidth / separation.1 ^ 2 + divisionWidth / separation.1 := by have hproduct := mul_le_mul hnumDenMax hnormWidth (RatInterval.width_nonneg denominator.normSq) hdivisionMax0 exact add_le_add (div_le_div_of_nonneg_right hproduct (sq_nonneg separation.1)) (div_le_div_of_nonneg_right hnumDenWidth separation.2.le) _ = quotientWidth := rfl have hquotientMax : (numerator.div denominator (certificate.away k hk)).maxAbs ≤ divisionMax / separation.1 := by exact (ComplexRatInterval.div_maxAbs (certificate.away k hk) separation.2 hsep).trans (div_le_div_of_nonneg_right hnumDenMax separation.2.le) have hexactNorm : ‖exactCircleNode program.radius schedule k‖ = program.radius := by have hradius0 : (0 : ℝ) ≤ program.radius := by exact_mod_cast program.radius_pos.le simp [exactCircleNode, Complex.norm_exp, abs_of_nonneg hradius0] have hnodeValueBound : max |(exactCircleNode program.radius schedule k).re| |(exactCircleNode program.radius schedule k).im| ≤ (program.radius : ℝ) := by exact max_le ((Complex.abs_re_le_norm _).trans_eq hexactNorm) ((Complex.abs_im_le_norm _).trans_eq hexactNorm) have hnodeMax : node.maxAbs ≤ program.radius + scheduled.target.1 := by exact ComplexRatInterval.maxAbs_le_of_contains_width hnodeContains hnodeValueBound hnodeWidth have hquotientWidth0 : 0 ≤ quotientWidth := by dsimp [quotientWidth] exact add_nonneg (div_nonneg (mul_nonneg hdivisionMax0 hnormWidth0) (sq_nonneg separation.1)) (div_nonneg hdivisionWidth0 separation.2.le) rw [ContourProgram.integrandNodes, dif_pos hk] change ((numerator.div denominator (certificate.away k hk)).mul node).width ≤ _ calc ((numerator.div denominator (certificate.away k hk)).mul node).width ≤ 2 * ((numerator.div denominator (certificate.away k hk)).maxAbs * node.width + node.maxAbs * (numerator.div denominator (certificate.away k hk)).width) := ComplexRatInterval.mul_width _ _ _ ≤ 2 * ((divisionMax / separation.1) * scheduled.target.1 + (program.radius + scheduled.target.1) * quotientWidth) := by have hquotientA0 : 0 ≤ (numerator.div denominator (certificate.away k hk)).maxAbs := (abs_nonneg (numerator.div denominator (certificate.away k hk)).re.lo).trans ((le_max_left _ _).trans (le_max_left _ _)) have hnodeA0 : 0 ≤ node.maxAbs := (abs_nonneg node.re.lo).trans ((le_max_left _ _).trans (le_max_left _ _)) nlinarith [mul_le_mul hquotientMax hnodeWidth (show 0 ≤ node.width from (RatInterval.width_nonneg node.re).trans (le_max_left _ _)) (div_nonneg hdivisionMax0 separation.2.le), mul_le_mul hnodeMax hquotientWidth (show 0 ≤ (numerator.div denominator (certificate.away k hk)).width from (RatInterval.width_nonneg (numerator.div denominator (certificate.away k hk)).re).trans (le_max_left _ _)) (by linarith [program.radius_pos])] _ = program.nodePropagationBound bounds separation.1 scheduled.target.1 := by rfl

theorem denominator_ne_zero_at_nodes (program : ContourProgram) (schedule : Schedule) (certificate : DenominatorCertificate program schedule) {k : ℕ} (hk : k ≤ schedule.mesh) : program.denominator.value (exactCircleNode program.radius schedule k) ≠ 0 := by have hnode : (program.nodeBox schedule k).Contains (exactCircleNode program.radius schedule k) := by exact circleNode_sound program.radius schedule hk have hden := program.denominator.sound hnode schedule.fuel have hnorm := ComplexRatInterval.normSq_sound hden have hnorm_ne : ‖program.denominator.value (exactCircleNode program.radius schedule k)‖ ^ 2 ≠ 0 := RatInterval.ne_zero_of_contains (certificate.away k hk) hnorm intro hz apply hnorm_ne simp [hz]

theorem evaluate_contains (program : ContourProgram) (schedule : Schedule) (certificate : DenominatorCertificate program schedule) (hden : ∀ u ∈ Icc (0 : ℝ) 1, program.denominator.value ((program.radius : ℂ) * exp (((2 : ℝ) * pi * u) * I)) ≠ 0) (hLip : ∀ s ∈ Icc (0 : ℝ) 1, ∀ t ∈ Icc (0 : ℝ) 1, ‖program.normalizedIntegrand schedule s - program.normalizedIntegrand schedule t‖ ≤ (schedule.magnitude : ℝ) * |s - t|) : (program.evaluate schedule certificate).Contains program.normalizedContourIntegral := by rw [normalizedContourIntegral_eq program schedule hden] apply ComplexRatInterval.expand_contains exact CircleMesh.integralEnclosure_sound schedule.magnitude_nonneg schedule.mesh_pos hLip (fun k hk => integrandNodes_sound program schedule certificate hk)

theorem evaluate_width (program : ContourProgram) (bounds : ContourValueBounds program) (separation : PosRat) (scheduled : CertifiedProgramSchedule program bounds separation) (certificate : DenominatorCertificate program scheduled.schedule) (hseparation : separation.1 ≤ certificate.separation) : (program.evaluate scheduled.schedule certificate).width ≤ scheduled.schedule.tolerance.1 := by let schedule := scheduled.schedule have hcoordinate : ∀ k ≤ schedule.mesh, (program.integrandNodes schedule certificate k).re.width ≤ schedule.nodeBudget ∧ (program.integrandNodes schedule certificate k).im.width ≤ schedule.nodeBudget := by intro k hk have h := (integrandNodes_width_le_propagation program bounds separation scheduled certificate hseparation hk).trans scheduled.propagation_le exact ⟨(le_max_left _ _).trans h, (le_max_right _ _).trans h⟩ have hbase := CircleMesh.integralEnclosure_width schedule.nodeBudget_nonneg schedule.magnitude_nonneg schedule.mesh_pos hcoordinate rw [ContourProgram.evaluate, ComplexRatInterval.width, ComplexRatInterval.expand, RatInterval.width_expand, RatInterval.width_expand] apply max_le · linarith [hbase.1, schedule.mesh_error_le, schedule.budget_sum_le] · linarith [hbase.2, schedule.mesh_error_le, schedule.budget_sum_le]

theorem certified_contour_evaluation_inverseMax (n : ℕ) (program : ContourProgram) (bounds : ContourValueBounds program) (separation : PosRat) (scheduled : CertifiedProgramSchedule program bounds separation) (hschedule : scheduled.schedule.tolerance = inverseMaxTolerance n) (certificate : DenominatorCertificate program scheduled.schedule) (hseparation : separation.1 ≤ certificate.separation) (hden : ∀ u ∈ Icc (0 : ℝ) 1, program.denominator.value ((program.radius : ℂ) * exp (((2 : ℝ) * pi * u) * I)) ≠ 0) (hLip : ∀ s ∈ Icc (0 : ℝ) 1, ∀ t ∈ Icc (0 : ℝ) 1, ‖program.normalizedIntegrand scheduled.schedule s - program.normalizedIntegrand scheduled.schedule t‖ ≤ (scheduled.schedule.magnitude : ℝ) * |s - t|) : (program.evaluate scheduled.schedule certificate).re.width ≤ 1 / (max n 1 : ℚ) := by have h := (certified_contour_evaluation program bounds separation scheduled certificate hseparation hden hLip).2 rw [hschedule] at h exact h

theorem midpoint_error_le_half_width {I : RatInterval} {x : ℝ} (hx : I.Contains x) : |(((I.lo + I.hi) / 2 : ℚ) : ℝ) - x| ≤ (I.width : ℝ) / 2 := by simp only [RatInterval.Contains] at hx simp only [RatInterval.width] rw [abs_le] push_cast constructor <;> linarith [hx.1, hx.2]

theorem sinNameApprox_spec (x : CertifiedReal) (fuel : ℕ) : (sinNameApprox x fuel).Contains (sin x.value) ∧ (sinNameApprox x (fuel + 1)).Subinterval (sinNameApprox x fuel) := by induction fuel with | zero => have h0 := sinInterval_sound (x.contains 0) 0 have h1 := sinInterval_sound (x.contains 1) 1 exact ⟨h0, RatInterval.tighten_subinterval_left h0 h1⟩ | succ fuel ih => have hnew := sinInterval_sound (x.contains (fuel + 1)) (fuel + 1) have hcur : (sinNameApprox x (fuel + 1)).Contains (sin x.value) := by rw [sinNameApprox] exact RatInterval.tighten_sound ih.1 hnew have hnext := sinInterval_sound (x.contains (fuel + 2)) (fuel + 2) constructor · exact hcur · rw [show fuel + 1 + 1 = fuel + 2 by omega, sinNameApprox] exact RatInterval.tighten_subinterval_left hcur hnext

theorem cosNameApprox_spec (x : CertifiedReal) (fuel : ℕ) : (cosNameApprox x fuel).Contains (cos x.value) ∧ (cosNameApprox x (fuel + 1)).Subinterval (cosNameApprox x fuel) := by induction fuel with | zero => have h0 := cosInterval_sound (x.contains 0) 0 have h1 := cosInterval_sound (x.contains 1) 1 exact ⟨h0, RatInterval.tighten_subinterval_left h0 h1⟩ | succ fuel ih => have hnew := cosInterval_sound (x.contains (fuel + 1)) (fuel + 1) have hcur : (cosNameApprox x (fuel + 1)).Contains (cos x.value) := by rw [cosNameApprox] exact RatInterval.tighten_sound ih.1 hnew have hnext := cosInterval_sound (x.contains (fuel + 2)) (fuel + 2) constructor · exact hcur · rw [show fuel + 1 + 1 = fuel + 2 by omega, cosNameApprox] exact RatInterval.tighten_subinterval_left hcur hnext

theorem power_div_factorial_le (q : ℚ) (A d N r : ℕ) (hq : |q| ≤ A) (hA : 1 ≤ A) (hd : 1 ≤ d) (hr : 1 ≤ r) (hN : 32 * d * A ^ 2 ≤ N) : |q| ^ (2 * N + r) / ((2 * N + r).factorial : ℚ) ≤ 1 / (4 * d : ℕ) := by have hp : |q| ^ (2 * N + r) ≤ (A : ℚ) ^ (2 * N + r) := pow_le_pow_left₀ (abs_nonneg q) hq _ have hf : (4 * d : ℚ) * (A : ℚ) ^ (2 * N + r) ≤ ((2 * N + r).factorial : ℚ) := by exact_mod_cast factorial_power_bound A d N r hA hd hr hN apply (div_le_div_iff₀ (by positivity : (0 : ℚ) < (2 * N + r).factorial) (by positivity : (0 : ℚ) < (4 * d : ℕ))).2 calc |q| ^ (2 * N + r) * (4 * d : ℕ) ≤ (A : ℚ) ^ (2 * N + r) * (4 * d : ℕ) := mul_le_mul_of_nonneg_right hp (by positivity) _ = (4 * d : ℚ) * (A : ℚ) ^ (2 * N + r) := by push_cast ring _ ≤ ((2 * N + r).factorial : ℚ) := hf _ = 1 * ((2 * N + r).factorial : ℚ) := by ring

theorem cosName_width_at_precision (x : CertifiedReal) (ε : PosRat) : (cosNameApprox x (trigNamePrecision x ε)).width ≤ ε.1 := by let δ : PosRat := ⟨ε.1 / 2, div_pos ε.2 (by norm_num)⟩ let d : ℕ := ε.1.den + 1 let A : ℕ := (x.approx 0).maxAbs.num.natAbs + 2 let N : ℕ := trigNamePrecision x ε let K := x.approx N have hmod : x.modulus δ ≤ N := by exact le_max_left _ _ have hN : 32 * d * A ^ 2 ≤ N := by exact le_max_right _ _ have hKw : K.width ≤ δ.1 := (RatInterval.width_mono (CertifiedReal.approx_mono x hmod)).trans (x.width_modulus δ) have hKK0 : K.Subinterval (x.approx 0) := CertifiedReal.approx_mono x (zero_le N) have hmidK : K.Contains (intervalMid K : ℝ) := by constructor <;> exact_mod_cast (by dsimp [intervalMid] linarith [K.lo_le_hi]) have hmid0 := RatInterval.Contains.mono hKK0 hmidK have hM0 : 0 ≤ (x.approx 0).maxAbs := (abs_nonneg (x.approx 0).lo).trans (le_max_left _ _) have hMle : (x.approx 0).maxAbs ≤ ((x.approx 0).maxAbs.num.natAbs : ℕ) := by have hn0 : 0 ≤ (x.approx 0).maxAbs.num := Rat.num_nonneg.mpr hM0 calc (x.approx 0).maxAbs = ((x.approx 0).maxAbs.num : ℚ) / ((x.approx 0).maxAbs.den : ℕ) := (Rat.num_div_den _).symm _ ≤ ((x.approx 0).maxAbs.num : ℚ) := div_le_self (by exact_mod_cast hn0) (by exact_mod_cast Rat.den_pos _) _ = ((x.approx 0).maxAbs.num.natAbs : ℕ) := by have hi : ((x.approx 0).maxAbs.num.natAbs : ℤ) = (x.approx 0).maxAbs.num := Int.natAbs_of_nonneg hn0 exact (congrArg (fun z : ℤ => (z : ℚ)) hi).symm have hmidA : |intervalMid K| ≤ (A : ℚ) := by exact (abs_le_max_abs_abs (by exact_mod_cast hmid0.1) (by exact_mod_cast hmid0.2)).trans (hMle.trans (by exact_mod_cast le_add_right (x.approx 0).maxAbs.num.natAbs 2)) have herr : cosError (intervalMid K) N ≤ 1 / (4 * d : ℕ) := by simpa [cosError] using power_div_factorial_le (intervalMid K) A d N 2 hmidA (by dsimp [A]; omega) (by dsimp [d]; omega) (by omega) hN have hinv : 1 / (d : ℚ) ≤ ε.1 := by exact (div_le_div_of_nonneg_left (by norm_num) (by positivity) (by dsimp [d]; exact_mod_cast Nat.le_succ ε.1.den)).trans (inv_den_le_of_pos ε.1 ε.2) have hcur := cosInterval_width K N have hout := RatInterval.width_mono (cosNameApprox_subinterval_current x N) have herr' : cosError (intervalMid K) N ≤ ε.1 / 4 := by calc _ ≤ (1 : ℚ) / ((4 * d : ℕ) : ℚ) := herr _ = (1 / (d : ℚ)) / 4 := by push_cast; field_simp _ ≤ _ := div_le_div_of_nonneg_right hinv (by norm_num) dsimp [δ] at hKw exact hout.trans (hcur.trans (by nlinarith [herr']))