diff --git a/Mathlib/NumberTheory/Bernoulli.lean b/Mathlib/NumberTheory/Bernoulli.lean index bd0e027a4f49bc..3cec04d70d916d 100644 --- a/Mathlib/NumberTheory/Bernoulli.lean +++ b/Mathlib/NumberTheory/Bernoulli.lean @@ -1,13 +1,17 @@ /- Copyright (c) 2020 Johan Commelin. All rights reserved. Released under Apache 2.0 license as described in the file LICENSE. -Authors: Johan Commelin, Kevin Buzzard +Authors: Johan Commelin, Kevin Buzzard, Seewoo Lee -/ module public import Mathlib.Algebra.BigOperators.Field -public import Mathlib.RingTheory.PowerSeries.Inverse +public import Mathlib.Algebra.Field.GeomSum +public import Mathlib.Data.Nat.Choose.Bounds public import Mathlib.RingTheory.PowerSeries.Exp +public import Mathlib.FieldTheory.Finite.Basic +public import Mathlib.RingTheory.ZMod.UnitsCyclic +public import Mathlib.NumberTheory.Padics.PadicNumbers /-! # Bernoulli numbers @@ -46,11 +50,23 @@ This formula is true for all $n$ and in particular $B_0=1$. Note that this is th for positive Bernoulli numbers, which we call `bernoulli'`. The negative Bernoulli numbers are then defined as `bernoulli := (-1)^n * bernoulli'`. +The proof of von Staudt-Clausen's theorem follows Rado's JLMS 1934 paper +"A New Proof of a Theorem of v. Staudt" + ## Main theorems -`sum_bernoulli : ∑ k ∈ Finset.range n, (n.choose k : ℚ) * bernoulli k = if n = 1 then 1 else 0` +* `sum_bernoulli : ∑ k ∈ Finset.range n, (n.choose k : ℚ) * bernoulli k = + if n = 1 then 1 else 0` +* `bernoulli.vonStaudt_clausen : bernoulli (2 * k) + ∑ p ∈ Finset.range (2 * k + 2) + with p.Prime ∧ (p - 1) ∣ 2 * k, (1 : ℚ) / p ∈ Set.range Int.cast` + +## References + +* https://en.wikipedia.org/wiki/Bernoulli_number +* [R. Rado, *A New Proof of a Theorem of v. Staudt*][Rado1934] -/ + @[expose] public section @@ -235,7 +251,6 @@ theorem bernoulli_spec' (n : ℕ) : rw [if_neg (succ_ne_zero _)] -- algebra facts have h₁ : (1, n) ∈ antidiagonal n.succ := by simp [mem_antidiagonal, add_comm] - have h₂ : (n : ℚ) + 1 ≠ 0 := by norm_cast have h₃ : (1 + n).choose n = n + 1 := by simp [add_comm] -- key equation: the corresponding fact for `bernoulli'` have H := bernoulli'_spec' n.succ @@ -383,3 +398,357 @@ theorem sum_Ico_pow (n p : ℕ) : _ = ∑ i ∈ range p.succ.succ, f' i := by simp_rw [sum_range_succ'] end Faulhaber + +section vonStaudtClausen + +/-! +### The von Staudt-Clausen Theorem + +Here we formalize Rado's proof of von Staudt-Clausen's theorem, which states that for any $k \ge 0$, +$$B_{2k} + \sum_{p \text{ prime}, (p - 1) \mid 2k} \frac{1}{p} \in \mathbb{Z}.$$ +Rado's proof is based on Faulhaber's theorem and induction on $k$. +-/ + +namespace bernoulli + +/- Indicator function that is `1` if `(p - 1) ∣ k` and `0` otherwise. -/ +private noncomputable def vonStaudtIndicator (k p : ℕ) : ℚ := + if (p - 1) ∣ k then 1 else 0 + +/-- Over `ZMod p`, the nonzero `l`-th power sum equals the negative indicator of `(p - 1) ∣ l`. -/ +private lemma sum_pow_add_indicator_eq_zero (p l : ℕ) [Fact p.Prime] : + (∑ v ∈ Finset.range p with v ≠ 0, (v : ZMod p) ^ l) + + (if (p - 1) ∣ l then (1 : ZMod p) else 0) = 0 := by + have hbij : (∑ v ∈ Finset.range p with v ≠ 0, (v : ZMod p) ^ l) = + ∑ u : (ZMod p)ˣ, (u : ZMod p) ^ l := + Finset.sum_bij' + (fun v hv ↦ Units.mk0 (v : ZMod p) (mt (ZMod.natCast_eq_zero_iff v p).mp (by + obtain ⟨hlt, hne⟩ := Finset.mem_filter.mp hv + exact Nat.not_dvd_of_pos_of_lt (Nat.pos_of_ne_zero hne) (Finset.mem_range.mp hlt)))) + (fun u _ ↦ (u : ZMod p).val) + (fun _ _ ↦ Finset.mem_univ _) + (fun u _ ↦ by simp [ZMod.val_lt, u.ne_zero]) + (fun v hv ↦ by simp [ZMod.val_cast_of_lt (Finset.mem_range.mp (Finset.mem_filter.mp hv).1)]) + (fun u _ ↦ Units.ext (ZMod.natCast_zmod_val _)) + (fun _ _ ↦ rfl) + rw [hbij, FiniteField.sum_pow_units, ZMod.card] + grind + +/-- A rational number `x` is `p`-integral if `p` does not divide its denominator. -/ +private def pIntegral (p : ℕ) (x : ℚ) [Fact p.Prime] : Prop := Rat.padicValuation p x ≤ 1 + +private lemma pIntegral_iff_not_dvd (p : ℕ) (x : ℚ) [Fact p.Prime] : + pIntegral p x ↔ ¬ p ∣ x.den := by + simp only [pIntegral, Rat.padicValuation_le_one_iff] + +private lemma sum_den_dvd_prod_den {ι : Type*} (s : Finset ι) (f : ι → ℚ) : + (∑ i ∈ s, f i).den ∣ ∏ i ∈ s, (f i).den := by + classical + induction s using Finset.induction_on with + | empty => simp + | @insert a s has ih => + rw [Finset.sum_insert has, Finset.prod_insert has] + exact dvd_trans (Rat.add_den_dvd (f a) (∑ x ∈ s, f x)) (mul_dvd_mul_left _ ih) + +private lemma pIntegral_sub (p : ℕ) [Fact p.Prime] (x y : ℚ) + (hx : pIntegral p x) (hy : pIntegral p y) : pIntegral p (x - y) := by + simpa [pIntegral] using (Rat.padicValuation p).map_sub_le hx hy + +private lemma pIntegral_sum {ι : Type*} (p : ℕ) [Fact p.Prime] (s : Finset ι) (f : ι → ℚ) + (hf : ∀ i ∈ s, pIntegral p (f i)) : pIntegral p (∑ i ∈ s, f i) := + (Rat.padicValuation p).map_sum_le hf + +private lemma pIntegral_ofInt (p : ℕ) [Fact p.Prime] (z : ℤ) : pIntegral p z := + (pIntegral_iff_not_dvd p _).2 (by simpa using (Nat.Prime.not_dvd_one Fact.out)) + +private lemma pIntegral_mul (p : ℕ) [Fact p.Prime] (x y : ℚ) + (hx : pIntegral p x) (hy : pIntegral p y) : pIntegral p (x * y) := by + simp only [pIntegral, map_mul] at * + exact mul_le_one' hx hy + +/-- Denominators of the "other primes" part of the indicator sum +stay coprime to a fixed prime `p`. -/ +private lemma prod_one_div_prime_den_coprime (k p : ℕ) [Fact p.Prime] : + (∏ q ∈ Finset.range (2 * k + 2) with q.Prime ∧ (q - 1) ∣ 2 * k ∧ q ≠ p, + ((1 : ℚ) / q).den).Coprime p := by + apply Nat.Coprime.prod_left + intro q hq + have h1 : q.Prime := (Finset.mem_filter.mp hq).2.1 + have h2 : ((1 : ℚ) / q).den = q := by simp [h1.ne_zero] + rw [h2] + exact (Nat.coprime_primes h1 Fact.out).mpr (Finset.mem_filter.mp hq).2.2.2 + +/-- Splits the prime-indexed correction sum into the `p`-term (`vonStaudtIndicator / p`) +plus the rest. -/ +private lemma sum_one_div_prime_eq_indicator_div_add (k p : ℕ) (hk : k > 0) [Fact p.Prime] : + (∑ q ∈ Finset.range (2 * k + 2) with q.Prime ∧ (q - 1) ∣ 2 * k, (1 : ℚ) / q) = + vonStaudtIndicator (2 * k) p / p + ∑ q ∈ Finset.range (2 * k + 2) with + q.Prime ∧ (q - 1) ∣ 2 * k ∧ q ≠ p, (1 : ℚ) / q := by + by_cases hdvd : (p - 1) ∣ 2 * k + · -- p is in the filtered set; extract its term + have hp_mem : p ∈ (Finset.range (2 * k + 2)).filter (fun q ↦ q.Prime ∧ (q - 1) ∣ 2 * k) := by + simp only [Finset.mem_filter, Finset.mem_range] + exact ⟨by have := Nat.le_of_dvd (by lia) hdvd; lia, Fact.out, hdvd⟩ + rw [← Finset.add_sum_erase _ _ hp_mem] + simp only [vonStaudtIndicator, if_pos hdvd] + congr 1 + apply Finset.sum_congr _ (fun _ _ ↦ rfl) + grind + · -- p is not in the filtered set; indicator is 0, filter sets are equal + simp only [vonStaudtIndicator, if_neg hdvd, zero_div, zero_add] + exact Finset.sum_congr (Finset.filter_congr fun q _ ↦ + ⟨fun ⟨hpr, hd⟩ ↦ ⟨hpr, hd, fun h ↦ hdvd (h ▸ hd)⟩, + fun ⟨hpr, hd, _⟩ ↦ ⟨hpr, hd⟩⟩) fun _ _ ↦ rfl + +/-- If the `p`-adic valuation of `M` is at most `N`, then `p^N / M` is `p`-integral. -/ +private lemma pIntegral_pow_div (p M N : ℕ) [Fact p.Prime] (hM : M ≠ 0) + (hv : M.factorization p ≤ N) : pIntegral p ((p : ℚ) ^ N / M) := by + set e := M.factorization p + set M' := M / p ^ e + have hM'_cop : M'.Coprime p := (Nat.coprime_ordCompl Fact.out hM).symm + have hp_ne : (p : ℚ) ≠ 0 := Nat.cast_ne_zero.mpr (Nat.Prime.ne_zero Fact.out) + -- Rewrite p^N / M as p^(N-e) / M' where M' = M / p^e is coprime to p + have hdecomp : p ^ e * M' = M := Nat.ordProj_mul_ordCompl_eq_self M p + have hM_eq : (M : ℚ) = ↑(p ^ e) * ↑M' := by rw [← hdecomp]; simp + have hrw : (p : ℚ) ^ N / M = (p : ℚ) ^ (N - e) / M' := by + rw [hM_eq, Nat.cast_pow, div_mul_eq_div_div] + congr 1 + rw [div_eq_iff (pow_ne_zero e hp_ne), ← pow_add, Nat.sub_add_cancel hv] + have hM'_eq : ((p : ℚ) ^ (N - e) / (M' : ℚ)) = Rat.divInt (p ^ (N - e) : ℤ) (M' : ℤ) := by + norm_cast + simp + rw [hrw] + exact (pIntegral_iff_not_dvd p _).2 ((Nat.Prime.coprime_iff_not_dvd Fact.out).1 + (hM'_cop.coprime_dvd_left (by + rw [hM'_eq]; exact Int.natCast_dvd_natCast.mp (Rat.den_dvd _ _))).symm) + +/-- The `i = 1` Faulhaber term is `p`-integral (handled separately for `p = 2` and odd `p`). -/ +private lemma pIntegral_bernoulli_one_term (k p : ℕ) (hk : k > 0) [Fact p.Prime] : + pIntegral p (bernoulli 1 * (2 * k) * (p : ℚ) ^ (2 * k - 1) / (2 * k)) := by + rw [bernoulli_one] + obtain rfl | hp2 := eq_or_ne p 2 + · have hsub : 2 * k - 1 = (2 * k - 2) + 1 := by lia + have h : ((-1 / 2 : ℚ) * (2 * k) * (2 : ℚ) ^ (2 * k - 1) / (2 * k)) = + -(2 : ℤ) ^ (2 * k - 2) := by + rw [hsub, pow_succ]; push_cast; field_simp + simpa [h] using pIntegral_ofInt _ (-(2 : ℤ) ^ (2 * k - 2)) + · have hk_ne : (2 * k : ℚ) ≠ 0 := by positivity + have hrw : (-1 / 2 : ℚ) * (2 * k) * (p : ℚ) ^ (2 * k - 1) / (2 * k) = + (-1 : ℤ) * ((p : ℚ) ^ (2 * k - 1) / 2) := by + field_simp [hk_ne]; push_cast; ring + rw [hrw] + exact pIntegral_mul p _ _ (by exact_mod_cast pIntegral_ofInt p (-1)) + (pIntegral_pow_div p 2 (2 * k - 1) two_ne_zero (by + rw [Nat.factorization_eq_zero_of_not_dvd (fun h ↦ by + have := Nat.le_of_dvd two_pos h; have := (Fact.out : p.Prime).two_le; lia)] + exact Nat.zero_le _)) + +/-- Main valuation estimate behind the contradiction step for even-index summands. -/ +private lemma factorization_succ_le_sub_one (p d : ℕ) [Fact p.Prime] (hd : d ≥ 2) : + (d + 1).factorization p ≤ d - 1 := by + -- Special case: p = 2, d = 2 (since 3 is not divisible by 2) + obtain ⟨rfl, rfl⟩ | hne : (p = 2 ∧ d = 2) ∨ ¬(p = 2 ∧ d = 2) := by tauto + · simp [Nat.factorization_eq_zero_of_not_dvd (by decide : ¬(2 ∣ 3))] + · -- For all other prime p and d ≥ 2, we have d + 1 ≤ p^(d-1) + apply Nat.factorization_le_of_le_pow + have hp2 := (Fact.out : p.Prime).two_le + suffices ∀ n : ℕ, n ≥ 2 → ¬(p = 2 ∧ n = 2) → n + 1 ≤ p ^ (n - 1) from this d hd hne + intro n hn hne' + induction hn with + | refl => norm_num at hne' ⊢; lia + | @step m hm IH => + by_cases hm2 : p = 2 ∧ m = 2 + · obtain ⟨rfl, rfl⟩ := hm2 + norm_num + · have hm_eq : m = m - 1 + 1 := by lia + calc m + 1 + 1 ≤ p ^ (m - 1) + 1 := by linarith [IH hm2] + _ ≤ p ^ (m - 1) * p := by nlinarith [Nat.one_le_pow (m - 1) p (by lia)] + _ = p ^ m := by rw [hm_eq]; exact pow_succ .. + +/-- Multiplicative variant of the binomial coefficient denominator rewrite +as in Rado's summand. -/ +private lemma choose_two_mul_succ_mul_div_eq (k m : ℕ) (x : ℚ) (hm_lt : m < k) : + ((2 * k + 1).choose (2 * m) : ℚ) * x / (2 * k + 1) = + ((2 * k).choose (2 * m) : ℚ) * x / (2 * k - 2 * m + 1) := by + have h : ((2 * k + 1).choose (2 * m) : ℚ) / (2 * k + 1) = + ((2 * k).choose (2 * m) : ℚ) / (2 * k - 2 * m + 1) := by + have hm_le : 2 * m ≤ 2 * k + 1 := by lia + rw [div_eq_div_iff (by norm_cast) (by norm_cast; lia)] + conv_rhs => norm_cast; rw [Nat.choose_mul_succ_eq] + push_cast [Nat.cast_sub hm_le]; ring + rw [mul_comm ((2 * k + 1).choose (2 * m) : ℚ) x, mul_div_assoc, + mul_comm ((2 * k).choose (2 * m) : ℚ) x, mul_div_assoc, h] + +/-- `p`-integrality of the core even-index summand after denominator normalization. -/ +private lemma pIntegral_choose_mul_pow_div (k m p : ℕ) (hm_lt : m < k) [Fact p.Prime] + (hd : 2 * k - 2 * m ≥ 2) : + pIntegral p (((2 * k).choose (2 * m) : ℚ) * (p : ℚ) ^ (2 * k - 2 * m - 1) / + (2 * k - 2 * m + 1)) := by + set d := 2 * k - 2 * m with hd_def + have ⟨hd_ne_zero, hd_plus_one_ne_zero, h_exp, hkm⟩ : + d ≠ 0 ∧ d + 1 ≠ 0 ∧ 2 * k - 2 * m - 1 = d - 1 ∧ 2 * m ≤ 2 * k := by lia + have h_denom_rat : (2 * (k : ℚ) - 2 * m + 1) = ((d + 1 : ℕ) : ℚ) := by + simp only [hd_def]; push_cast [Nat.cast_sub hkm]; ring + rw [h_exp, h_denom_rat, mul_div_assoc] + exact pIntegral_mul p _ _ (by exact_mod_cast pIntegral_ofInt p ((2 * k).choose (2 * m))) + (pIntegral_pow_div p (d + 1) (d - 1) hd_plus_one_ne_zero (factorization_succ_le_sub_one p d hd)) + +/-- Uses the induction hypothesis on `B_{2m} + e_{2m}(p)/p` +to prove `p`-integrality of the even term. -/ +private lemma pIntegral_bernoulli_even_term (k m p : ℕ) (hm_lt : m < k) [Fact p.Prime] + (ih : pIntegral p (bernoulli (2 * m) + vonStaudtIndicator (2 * m) p / p)) : + pIntegral p (bernoulli (2 * m) * ((2 * k + 1).choose (2 * m)) * + (p : ℚ) ^ (2 * k - 2 * m) / (2 * k + 1)) := by + have hp_ne : (p : ℚ) ≠ 0 := Nat.cast_ne_zero.mpr (Nat.Prime.ne_zero Fact.out) + set P := (p : ℚ) ^ (2 * k - 2 * m - 1) + have hpow : (p : ℚ) ^ (2 * k - 2 * m) = P * p := by + rw [(by lia : 2 * k - 2 * m = (2 * k - 2 * m - 1) + 1), pow_succ] + have hdecomp : bernoulli (2 * m) * ((2 * k + 1).choose (2 * m)) * + (p : ℚ) ^ (2 * k - 2 * m) / (2 * k + 1) = + (bernoulli (2 * m) + vonStaudtIndicator (2 * m) p / p) * + ((2 * k + 1).choose (2 * m)) * (p : ℚ) ^ (2 * k - 2 * m) / (2 * k + 1) - + vonStaudtIndicator (2 * m) p * ((2 * k + 1).choose (2 * m)) * + P / (2 * k + 1) := by rw [hpow]; field_simp [hp_ne]; ring + rw [hdecomp] + have hcmp := pIntegral_choose_mul_pow_div k m p hm_lt (by lia) + apply pIntegral_sub + · rw [mul_assoc, mul_div_assoc] + apply pIntegral_mul _ _ _ ih + have hpow_mul : ((2 * k).choose (2 * m) : ℚ) * (p : ℚ) ^ (2 * k - 2 * m) / + (2 * k - 2 * m + 1) = + (p : ℚ) * (((2 * k).choose (2 * m) : ℚ) * P / (2 * k - 2 * m + 1)) := by + rw [hpow]; ring + rw [choose_two_mul_succ_mul_div_eq k m _ hm_lt, hpow_mul] + exact pIntegral_mul _ _ _ (pIntegral_ofInt p p) hcmp + · unfold vonStaudtIndicator; split_ifs + · simp only [one_mul]; rw [choose_two_mul_succ_mul_div_eq k m _ hm_lt]; exact hcmp + · simp only [zero_mul, zero_div, pIntegral, map_zero]; exact bot_le + +/-- The full remainder sum in Faulhaber's formula is `p`-integral. -/ +private lemma pIntegral_faulhaber_sum (k p : ℕ) (hk : k > 0) [Fact p.Prime] + (ih : ∀ m, 0 < m → m < k → pIntegral p (bernoulli (2 * m) + vonStaudtIndicator (2 * m) p / p)) : + pIntegral p (∑ i ∈ Finset.range (2 * k), + bernoulli i * ((2 * k + 1).choose i) * (p : ℚ) ^ (2 * k - i) / (2 * k + 1)) := by + apply pIntegral_sum + intro i hi + rw [Finset.mem_range] at hi + rcases i with _ | _ | i + · simp only [bernoulli_zero, one_mul, Nat.choose_zero_right, Nat.cast_one, Nat.sub_zero] + exact_mod_cast pIntegral_pow_div p (2 * k + 1) (2 * k) (by lia) + (Nat.factorization_le_of_le_pow <| calc 2 * k + 1 = (2 * k + 1).choose 1 := by simp + _ ≤ 2 ^ (2 * k) := Nat.choose_succ_le_two_pow _ 1 + _ ≤ p ^ (2 * k) := Nat.pow_le_pow_left (Fact.out : p.Prime).two_le _) + · simp only [zero_add, Nat.choose_one_right] + convert pIntegral_bernoulli_one_term k p hk using 1 + push_cast; field_simp + · rcases Nat.even_or_odd (i + 2) with ⟨m, hm⟩ | hodd + · have ⟨hm_pos, hm_lt, hi_eq⟩ : 0 < m ∧ m < k ∧ i + 2 = 2 * m := by lia + simp only [hi_eq] + exact pIntegral_bernoulli_even_term k m p hm_lt (ih m hm_pos hm_lt) + · simp only [bernoulli_eq_zero_of_odd hodd (by lia), + zero_mul, zero_div, pIntegral, map_zero]; exact bot_le + +private lemma exists_int_sum_pow_add_indicator_eq (k p : ℕ) [Fact p.Prime] : + ∃ T : ℤ, (∑ v ∈ Finset.range p with v ≠ 0, (v : ℚ) ^ (2 * k)) + + vonStaudtIndicator (2 * k) p = p * T := by + -- Get divisibility in ℤ from the ZMod result + have ⟨T, hT⟩ : ∃ T : ℤ, (∑ v ∈ Finset.range p with v ≠ 0, (v : ℤ) ^ (2 * k)) + + (if (p - 1) ∣ 2 * k then 1 else 0) = p * T := by + have : (↑((∑ v ∈ Finset.range p with v ≠ 0, (v : ℤ) ^ (2 * k)) + + (if (p - 1) ∣ 2 * k then 1 else 0)) : ZMod p) = 0 := by + push_cast; exact sum_pow_add_indicator_eq_zero p (2 * k) + exact (ZMod.intCast_zmod_eq_zero_iff_dvd _ _).mp this + refine ⟨T, ?_⟩ + have : vonStaudtIndicator (2 * k) p = + ((if (p - 1) ∣ 2 * k then (1 : ℤ) else 0) : ℚ) := by + unfold vonStaudtIndicator; split_ifs <;> simp + rw [this]; exact_mod_cast hT + +private lemma sum_pow_filter_eq_faulhaber (k p : ℕ) (hk : 0 < k) : + (∑ v ∈ Finset.range p with v ≠ 0, (v : ℚ) ^ (2 * k)) = + (∑ i ∈ Finset.range (2 * k), bernoulli i * ((2 * k + 1).choose i) * + (p : ℚ) ^ (2 * k + 1 - i) / (2 * k + 1)) + p * bernoulli (2 * k) := by + have hfilter : (∑ v ∈ Finset.range p with v ≠ 0, (v : ℚ) ^ (2 * k)) = + ∑ v ∈ Finset.range p, (v : ℚ) ^ (2 * k) := by + rw [Finset.sum_filter]; exact Finset.sum_congr rfl fun v _ ↦ by + by_cases hv : v = 0 <;> simp [hv, (by lia : 2 * k ≠ 0)] + have hsub : 2 * k + 1 - 2 * k = 1 := by lia + have hne : (2 * k + 1 : ℚ) ≠ 0 := by positivity + rw [hfilter, sum_range_pow, Finset.sum_range_succ] + push_cast + congr 1 + rw [Nat.choose_succ_self_right, hsub, pow_one] + push_cast; field_simp [hne] + +private lemma faulhaber_sum_div_prime_eq (k p : ℕ) [Fact p.Prime] : + (∑ i ∈ Finset.range (2 * k), bernoulli i * ((2 * k + 1).choose i : ℚ) * + (p : ℚ) ^ (2 * k + 1 - i) / (2 * k + 1 : ℚ)) / (p : ℚ) = + ∑ i ∈ Finset.range (2 * k), bernoulli i * ((2 * k + 1).choose i : ℚ) * + (p : ℚ) ^ (2 * k - i) / (2 * k + 1 : ℚ) := by + have hp_ne : (p : ℚ) ≠ 0 := Nat.cast_ne_zero.mpr (Fact.out : p.Prime).ne_zero + rw [Finset.sum_div] + exact Finset.sum_congr rfl fun i hi ↦ by + have hi_lt := Finset.mem_range.mp hi + have hsub : 2 * k + 1 - i = (2 * k - i) + 1 := by lia + rw [hsub, pow_succ] + field_simp [hp_ne] + +/-- Rearranges the Faulhaber identity and power-sum congruence to isolate +`bernoulli (2*k) + vonStaudtIndicator (2*k) p / p`. -/ +private lemma bernoulli_add_indicator_eq_sub (k p : ℕ) (hk : k > 0) [Fact p.Prime] : + ∃ T : ℤ, bernoulli (2 * k) + vonStaudtIndicator (2 * k) p / p = + T - (∑ i ∈ Finset.range (2 * k), + bernoulli i * ((2 * k + 1).choose i) * (p : ℚ) ^ (2 * k - i) / (2 * k + 1)) := by + obtain ⟨T, hT⟩ := exists_int_sum_pow_add_indicator_eq k p + use T + have hp_ne : (p : ℚ) ≠ 0 := Nat.cast_ne_zero.mpr (Fact.out : p.Prime).ne_zero + have hAlg : bernoulli (2 * k) + vonStaudtIndicator (2 * k) p / p = + T - (∑ i ∈ Finset.range (2 * k), bernoulli i * ((2 * k + 1).choose i) * + (p : ℚ) ^ (2 * k + 1 - i) / (2 * k + 1)) / p := by + field_simp [hp_ne]; linarith [hT, sum_pow_filter_eq_faulhaber k p hk] + rw [hAlg]; congr 1; simpa using faulhaber_sum_div_prime_eq k p + +/-- For fixed prime `p`, the denominator of `B_{2k} + e_{2k}(p)/p` is not divisible by `p`. -/ +private lemma bernoulli_add_indicator_den_not_dvd (k p : ℕ) (hk : k > 0) [Fact p.Prime] : + ¬ p ∣ (bernoulli (2 * k) + vonStaudtIndicator (2 * k) p / p).den := by + induction k using Nat.strong_induction_on with + | _ k ih => + obtain ⟨T, hT⟩ := bernoulli_add_indicator_eq_sub k p hk + rw [hT] + have hT_int : pIntegral p (T : ℚ) := pIntegral_ofInt p T + have hR : pIntegral p (∑ i ∈ Finset.range (2 * k), + bernoulli i * ((2 * k + 1).choose i) * (p : ℚ) ^ (2 * k - i) / (2 * k + 1)) := by + apply pIntegral_faulhaber_sum k p hk + intro m hm_pos hm_lt + exact (pIntegral_iff_not_dvd p _).2 (ih m hm_lt hm_pos) + exact (pIntegral_iff_not_dvd p _).1 (pIntegral_sub p T _ hT_int hR) + +/-- Extends the fixed-prime nondivisibility result to the full prime correction sum. -/ +private lemma vonStaudt_sum_den_not_dvd (k p : ℕ) (hk : k > 0) [Fact p.Prime] : + ¬ p ∣ (bernoulli (2 * k) + ∑ q ∈ Finset.range (2 * k + 2) with + q.Prime ∧ (q - 1) ∣ 2 * k, (1 : ℚ) / q).den := by + rw [sum_one_div_prime_eq_indicator_div_add k p hk, ← add_assoc] + have hcop_ind := ((Nat.Prime.coprime_iff_not_dvd Fact.out).2 + (bernoulli_add_indicator_den_not_dvd k p hk)).symm + have hcop_rest := Nat.Coprime.of_dvd_left (sum_den_dvd_prod_den _ _) + (prod_one_div_prime_den_coprime k p) + have hcop := (Nat.Coprime.of_dvd_left (Rat.add_den_dvd _ _) (hcop_ind.mul_left hcop_rest)).symm + exact (Nat.Prime.coprime_iff_not_dvd Fact.out).1 hcop + +/-- **von Staudt-Clausen theorem:** For any natural number $k$, the sum +$$B_{2k} + \sum_{p - 1 \mid 2k} \frac{1}{p}$$ is an integer. +-/ +theorem vonStaudt_clausen (k : ℕ) : + bernoulli (2 * k) + ∑ p ∈ Finset.range (2 * k + 2) with p.Prime ∧ (p - 1) ∣ 2 * k, + (1 : ℚ) / p ∈ Set.range Int.cast := by + rcases Nat.eq_zero_or_pos k with rfl | hk + · exact ⟨1, by decide +kernel⟩ + · rw [Set.mem_range] + refine ⟨_, Rat.coe_int_num_of_den_eq_one ?_⟩ + by_contra h + obtain ⟨p, hp, hdvd⟩ := ne_one_iff_exists_prime_dvd.mp h + exact absurd hdvd (by let : Fact p.Prime := ⟨hp⟩; exact vonStaudt_sum_den_not_dvd k p hk) + +end bernoulli + +end vonStaudtClausen diff --git a/docs/references.bib b/docs/references.bib index 10df6668f614c2..e990cc99da7ffc 100644 --- a/docs/references.bib +++ b/docs/references.bib @@ -4673,6 +4673,17 @@ @Book{ Quillen1967 mrnumber = {223432} } +@Article{ Rado1934, + title = {A new proof of a theorem of V. Staudt}, + author = {Rado, R}, + journal = {Journal of the London Mathematical Society}, + volume = {1}, + number = {2}, + pages = {85--88}, + year = {1934}, + publisher = {Oxford University Press} +} + @InCollection{ ribenboim1971, author = {Ribenboim, Paulo}, title = {\'{E}pimorphismes de modules qui sont n\'{e}cessairement