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
124. 124
125. 125
126. 126
127. 127
128. 128
129. 129
130. 130
131. 131
132. 132
133. 133
134. 134
135. 135
136. 136
137. 137
138. 138
139. 139

import Mathlib
open Finset BigOperators Filter Topology

-- Important: b should be in ℝ so that it's real div not nat div
lemma geom_sum (b : ℝ) (N : ℕ) (h₁ : 2 ≤ N) (h₂ : 2 ≤ b) : ∑ a ∈ Ico 2 N, 1 / (b ^ a) = 1 / (b - 1) - 1 / b - 1 / (b ^ (N - 1) * (b - 1)) :=
  calc
    ∑ a ∈ Ico 2 N, 1 / (b ^ a) = ∑ a ∈ Ico 2 N, (1 / b) ^ a := by simp
    _ = ((1 / b) ^ 2 - (1 / b) ^ N) / (1 - (1 / b)) := by
      rw [geom_sum_Ico' (by grind) h₁]
    _ = 1 / (b - 1) - 1 / b - b / (b ^ N * (b - 1)) := by
      have : b - 1 ≠ 0 := by linarith
      field_simp
      ring
    _ = 1 / (b - 1) - 1 / b - 1 / (b ^ (N - 1) * (b - 1)) := by
      have h₃ : b ≠ 0 := by linarith
      have h₄ : b ^ N * (b - 1) = b * (b ^ (N - 1) * (b - 1)) := by
        rw [← mul_assoc, ← pow_succ' b (N - 1), Nat.sub_one_add_one]
        linarith
      simp [h₄, div_mul_cancel_left₀ h₃]

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 Finset.sum_congr rfl
      intro b bico
      have h₁ : 2 ≤ (b : ℝ) := Nat.ofNat_le_cast.mpr (List.left_le_of_mem_range' bico)
      exact geom_sum b N h h₁
    _ = ∑ 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 =>
        field_simp
        norm_num
      | step h ih =>
        rw [sum_Ico_succ_top, ih]
        field_simp
        exact h

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

-- 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
      · trivial
      · 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 linarith
    have : 2 ^ (n + 1) = 4 * 2 ^ (n - 1) := by grind [mul_pow_sub_one]
    simp [this]
    omega
  · interval_cases n <;> simp

theorem double_sum : Tendsto (fun N : ℕ => ∑ a ∈ Ico 2 N, ∑ b ∈ Ico 2 N, (1 : ℝ) / (b ^ a)) atTop (𝓝 (1 : ℝ)) := by
  have : atTop.HasBasis (fun _ : ℕ ↦ True) Set.Ici := atTop_basis
  rw [this.tendsto_iff (nhds_basis_Ioo_pos 1)]
  intro ε εpos
  simp
  use max (Nat.ceil (3 / ε)) 2
  intro n nlarge
  have h₁ : 2 ≤ n := le_of_max_le_right nlarge
  constructor <;> field_simp <;> rw [telescope_sum n h₁] <;> ring_nf
  · -- Greater than 1 - ε
    have h₂ : ∑ b ∈ Ico 2 n, (1 : ℝ) / (b ^ (n - 1) * (b - 1)) ≤ (n - 2) / 2 ^ (n - 1) := by
      have h₃ b (bico : b ∈ Ico 2 n) : (1 : ℝ) / (b ^ (n - 1) * (b - 1)) ≤ 1 / 2 ^ (n - 1) := by
        have h₄ : 2 ≤ (b : ℝ) := by
          norm_num
          exact (List.left_le_of_mem_range' bico)
        have h₅ : 1 / ((b : ℝ) - 1) ≤ 1 := (div_le_one₀ (by linarith)).mpr (by linarith)
        have h₆ : 1 / (b : ℝ) ^ (n - 1) ≤ 1 / 2 ^ (n - 1) := by
          rw [← one_div_pow, ← one_div_pow]
          apply pow_le_pow_left₀ (by norm_num)
          simp
          exact inv_anti₀ (by norm_num) h₄
        rw [← one_div_mul_one_div, mul_comm]
        nth_rw 5 [← one_mul 1]
        rw [mul_div_assoc]
        exact mul_le_mul h₅ h₆ (by norm_num) (by norm_num)
      have h₄ : ∑ b ∈ Ico 2 n, (1 : ℝ) / (b ^ (n - 1) * (b - 1)) ≤ ∑ b ∈ Ico 2 n, 1 / 2 ^ (n - 1) := sum_le_sum h₃
      field_simp at h₄
      field_simp
      exact h₄
    ring_nf at h₂
    have h₃ : 1 < 1 + ε - 1 / (n - 1) - (n - 2) / 2 ^ (n - 1) := by
      suffices 1 / (n - 1) + (n - 2) / 2 ^ (n - 1) < ε by linarith
      by_cases hε : 3 / 2 < ε
      · have : 2 ≤ (n : ℝ) := Nat.ofNat_le_cast.mpr h₁
        have : 1 / ((n : ℝ) - 1) ≤ 1 := div_le_one₀ (by linarith) |>.mpr (by linarith)
        suffices ((n : ℝ) - 2) / 2 ^ (n - 1) ≤ 1 / 2 by linarith
        have h₄ := exp_larger n h₁
        rw [mul_comm, ← one_mul (2 ^ (n - 1))] at h₄
        exact (div_le_div_iff₀ (by norm_num) (by norm_num)).mpr (by norm_cast)
      · simp at hε
        have : 3 / n ≤ ε := by
          have h₃ : 3 / ε ≤ n := Nat.ceil_le.mp (le_of_max_le_left nlarge)
          have h₄ : 0 < (n : ℝ) := by
            norm_num
            omega
          exact (div_le_comm₀ h₄ εpos).mpr h₃
        suffices 1 / ((n : ℝ) - 1) + ((n : ℝ) - 2) / 2 ^ (n - 1) < 3 / n by linarith
        rw [← lt_sub_iff_add_lt']
        have h₃ : 3 / (n : ℝ) - 1 / (n - 1) = (2 * n - 3) / (n * (n - 1)) := by
          have : (n : ℝ) ≠ 0 := by grind [Nat.cast_ne_zero]
          have : (n : ℝ) - 1 ≠ 0 := by grind
          field_simp
          ring
        rw [h₃]
        have h₄ : ((n : ℝ) - 2) * (n * (n - 1)) < (2 * n - 3) * 2 ^ (n - 1) := by
          have h₅ := exp_larger' n h₁
          have h₆ : ↑(2 * n - 3) = 2 * (n : ℝ) - 3 := by
            rw [Nat.cast_sub (by linarith)]
            simp
          rify [h₁, h₆] at h₅
          grind [Nat.cast_sub]
        have h₅ : 0 < (2 : ℝ) ^ (n - 1) := by norm_num
        exact (div_lt_div_iff₀ h₅ (by nlinarith)).mpr h₄
    ring_nf at h₃
    linarith
  · -- Less than 1
    have h₂ b (bico : b ∈ Ico 2 n) : 0 < (1 : ℝ) / (b ^ (n - 1) * (b - 1)) := by
      field_simp
      have : 1 < (b : ℝ) := Nat.one_lt_cast.mpr (List.left_le_of_mem_range' bico)
      exact mul_pos (pow_pos (by linarith) (n - 1)) (by linarith)
    have : 0 ≤ ∑ b ∈ Ico 2 n, (1 : ℝ) / (b ^ (n - 1) * (b - 1)) := sum_nonneg (fun b a => le_of_lt (h₂ b a))
    have : 0 < (1 : ℝ) / (n - 1) := by
      field_simp
      omega
    grind