1. 1
2. 2
3. 3
4. 4
5. 5
6. 6
7. 7
8. 8
9. 9
10. 10
11. 11
12. 12
13. 13
14. 14
15. 15
16. 16
17. 17
18. 18
19. 19
20. 20
21. 21
22. 22
23. 23
24. 24
25. 25
26. 26
27. 27
28. 28
29. 29
30. 30
31. 31
32. 32
33. 33
34. 34
35. 35
36. 36
37. 37
38. 38
39. 39
40. 40
41. 41
42. 42
43. 43
44. 44
45. 45
46. 46
47. 47
48. 48
49. 49
50. 50
51. 51
52. 52
53. 53
54. 54
55. 55
56. 56
57. 57
58. 58
59. 59
60. 60
61. 61
62. 62
63. 63
64. 64
65. 65
66. 66
67. 67
68. 68
69. 69
70. 70
71. 71
72. 72
73. 73
74. 74
75. 75
76. 76
77. 77
78. 78
79. 79
80. 80
81. 81
82. 82
83. 83
84. 84
85. 85
86. 86
87. 87
88. 88
89. 89
90. 90
91. 91
92. 92
93. 93
94. 94
95. 95
96. 96
97. 97
98. 98
99. 99
100. 100
101. 101
102. 102
103. 103
104. 104
105. 105
106. 106
107. 107
108. 108
109. 109
110. 110
111. 111
112. 112
113. 113
114. 114
115. 115
116. 116
117. 117
118. 118
119. 119
120. 120
121. 121
122. 122
123. 123

import Mathlib.Algebra.Field.GeomSum
import Mathlib.Analysis.RCLike.Basic
import Mathlib.Tactic.IntervalCases
import Mathlib.Tactic.Rify
open Finset Filter Topology

-- Important: b should be in ℝ so that it's real div not nat div

lemma geom_sum (n : ℕ) (b : ℝ) (hn : 2 ≤ n) (hb : 2 ≤ b) :
    ∑ a ∈ Ico 2 n, 1 / (b ^ a) = 1 / (b - 1) - 1 / b - 1 / (b ^ (n - 1) * (b - 1)) :=
  by
  have : ∑ a ∈ Ico 2 n, 1 / (b ^ a) = ∑ a ∈ Ico 2 n, (1 / b) ^ a := by simp
  rw [this, geom_sum_Ico (by grind) hn]
  field_simp
  have : b ^ n ≠ 0 := by positivity
  grind [inv_pow, mul_pow_sub_one]

lemma telescope_sum (n : ℕ) (h : 2 ≤ n) :
    ∑ a ∈ Ico 2 n, ∑ b ∈ Ico 2 n, (1 : ℝ) / (b ^ a) =
      1 - (1 : ℝ) / (n - 1) - ∑ b ∈ Ico 2 n, (1 : ℝ) / (b ^ (n - 1) * (b - 1)) :=
  by
  calc
    ∑ a ∈ Ico 2 n, ∑ b ∈ Ico 2 n, (1 : ℝ) / (b ^ a) = ∑ b ∈ Ico 2 n, ∑ a ∈ Ico 2 n, (1 : ℝ) / (b ^ a) := sum_comm
    _ = ∑ b ∈ Ico 2 n, ((1 : ℝ) / (b - 1) - (1 : ℝ) / b - (1 : ℝ) / (b ^ (n - 1) * (b - 1))) :=
      by
      apply sum_congr rfl
      exact fun b bico ↦ geom_sum n b h <| Nat.ofNat_le_cast.mpr <| List.left_le_of_mem_range' bico
    _ = ∑ b ∈ Ico 2 n, ((1 : ℝ) / (b - 1) - (1 : ℝ) / b) - ∑ b ∈ Ico 2 n, (1 : ℝ) / (b ^ (n - 1) * (b - 1)) := by
      apply sum_sub_distrib
    _ = 1 - (1 : ℝ) / (n - 1) - ∑ b ∈ Ico 2 n, (1 : ℝ) / (b ^ (n - 1) * (b - 1)) :=
      by
      suffices ∑ b ∈ Ico 2 n, ((1 : ℝ) / (b - 1) - (1 : ℝ) / b) = 1 - (1 : ℝ) / (n - 1) by rw [this]
      induction h with
      | refl =>
        norm_num
        rfl
      | step h ih =>
        rw [sum_Ico_succ_top h, ih]
        simp

/-- Exponentials grow quickly -/
lemma exp_larger (n : ℕ) (h : 2 ≤ n) : (n - 2) * 2 ≤ 2 ^ (n - 1) := by induction h <;> grind [mul_pow_sub_one]

/-- Another similar lemma -/
lemma exp_larger' (n : ℕ) (h : 2 ≤ n) : (n - 2) * (n * (n - 1)) < (2 * n - 3) * 2 ^ (n - 1) :=
  by
  by_cases h : 4 < n
  · have : n * (n - 1) * (n - 2) < 2 ^ (n + 1) := by
      induction h
      · decide
      · rename_i m _ _
        suffices (m + 1) * (m * (m - 1)) ≤ 2 * (m - 2) * (m * (m - 1)) by grind
        apply Nat.mul_le_mul_right (m * (m - 1))
        grind
    suffices 2 ^ (n + 1) ≤ (2 * n - 3) * 2 ^ (n - 1) by grind
    have : 2 ^ (n + 1) = 4 * 2 ^ (n - 1) := by grind [mul_pow_sub_one]
    rw [this, mul_le_mul_iff_left₀ (by positivity)]
    grind
  · interval_cases n <;> decide

theorem double_sum : Tendsto (fun n : ℕ ↦ ∑ a ∈ Ico 2 n, ∑ b ∈ Ico 2 n, (1 : ℝ) / (b ^ a)) atTop (𝓝 (1 : ℝ)) :=
  by
  rw [atTop_basis.tendsto_iff (nhds_basis_Ioo_pos 1)]
  intro ε εpos
  simp only [true_and, Set.mem_Ici, Set.mem_Ioo]
  use max (Nat.ceil (3 / ε)) 2
  intro n nlarge
  have h₁ := le_of_max_le_right nlarge
  rw [telescope_sum n h₁]
  constructor
  · -- Greater than 1 - ε

    suffices
      (∑ b ∈ Ico 2 n, (1 : ℝ) / (b ^ (n - 1) * (b - 1)) ≤ (n - 2) / 2 ^ (n - 1)) ∧
        1 / (n - 1) + (n - 2) / 2 ^ (n - 1) < ε
      by grind
    constructor
    · have h₂ b (bico : b ∈ Ico 2 n) : (1 : ℝ) / (b ^ (n - 1) * (b - 1)) ≤ 1 / 2 ^ (n - 1) :=
        by
        rw [← one_div_mul_one_div ((b : ℝ) ^ (n - 1))]
        have h₃ : 2 ≤ (b : ℝ) := by
          norm_cast
          exact (List.left_le_of_mem_range' bico)
        have h₄ : 1 / ((b : ℝ) - 1) ≤ 1 := div_le_one₀ (by grind) |>.mpr (by grind)
        have h₅ : 1 / (b : ℝ) ^ (n - 1) ≤ 1 / 2 ^ (n - 1) :=
          by
          simp_rw [one_div]
          exact inv_anti₀ (by norm_num) <| pow_le_pow_left₀ (by norm_num) h₃ (n - 1)
        grw [h₄, h₅]
        simp
      grw [sum_le_sum h₂]
      simp
      norm_cast
    · by_cases 3 / 2 < ε
      · have : 2 ≤ (n : ℝ) := Nat.ofNat_le_cast.mpr h₁
        have : 1 / ((n : ℝ) - 1) ≤ 1 := div_le_one₀ (by grind) |>.mpr (by grind)
        suffices ((n : ℝ) - 2) / 2 ^ (n - 1) ≤ 1 / 2 by grind
        have := exp_larger n h₁
        rw [← one_mul (2 ^ (n - 1))] at this
        exact div_le_div_iff₀ (by norm_num) (by norm_num) |>.mpr (by norm_cast)
      · have := (div_le_comm₀ (by positivity) εpos).mpr <| Nat.ceil_le.mp (le_of_max_le_left nlarge)
        grw [← this]
        have : 3 / (n : ℝ) - 1 / (n - 1) = (2 * n - 3) / (n * (n - 1)) :=
          by
          field_simp
          grind
        rw [← lt_sub_iff_add_lt', this]
        have h₂ : ((n : ℝ) - 2) * (n * (n - 1)) < (2 * n - 3) * 2 ^ (n - 1) :=
          by
          have := exp_larger' n h₁
          rify at this
          grind [Nat.cast_sub]
        have h₃ : 0 < (2 : ℝ) ^ (n - 1) := by norm_num
        exact div_lt_div_iff₀ h₃ (by nlinarith) |>.mpr h₂
  · -- Less than 1

    have h₂ b (bico : b ∈ Ico 2 n) : 0 ≤ (1 : ℝ) / (b ^ (n - 1) * (b - 1)) :=
      by
      rw [one_div_nonneg]
      have : 1 < (b : ℝ) := Nat.one_lt_cast.mpr (List.left_le_of_mem_range' bico)
      exact mul_nonneg (pow_nonneg (by grind) (n - 1)) (by grind)
    have : 0 < (1 : ℝ) / (n - 1) := by simp [Nat.one_lt_cast.mpr h₁]
    grind [sum_nonneg h₂]