|
| 1 | +/- |
| 2 | +Copyright (c) 2021 Ashvni Narayanan. All rights reserved. |
| 3 | +Released under Apache 2.0 license as described in the file LICENSE. |
| 4 | +Author: Ashvni Narayanan |
| 5 | +-/ |
| 6 | + |
| 7 | +import number_theory.bernoulli |
| 8 | + |
| 9 | +/-! |
| 10 | +# Bernoulli polynomials |
| 11 | +
|
| 12 | +The Bernoulli polynomials (defined here : https://en.wikipedia.org/wiki/Bernoulli_polynomials) |
| 13 | +are an important tool obtained from Bernoulli numbers. |
| 14 | +
|
| 15 | +## Mathematical overview |
| 16 | +
|
| 17 | +The $n$-th Bernoulli polynomial is defined as |
| 18 | +$$ B_n(X) = ∑_{k = 0}^n {n \choose k} (-1)^k * B_k * X^{n - k} $$ |
| 19 | +where $B_k$ is the $k$-th Bernoulli number. The Bernoulli polynomials are generating functions, |
| 20 | +$$ t * e^{tX} / (e^t - 1) = ∑_{n = 0}^{\infty} B_n(X) * \frac{t^n}{n!} $$ |
| 21 | +
|
| 22 | +## Implementation detail |
| 23 | +
|
| 24 | +Bernoulli polynomials are defined using `bernoulli`, the Bernoulli numbers. |
| 25 | +
|
| 26 | +## Main theorems |
| 27 | +
|
| 28 | +- `sum_bernoulli_poly`: The sum of the $k^\mathrm{th}$ Bernoulli polynomial with binomial |
| 29 | + coefficients up to n is `(n + 1) * X^n`. |
| 30 | +- `exp_bernoulli_poly`: The Bernoulli polynomials act as generating functions for the exponential. |
| 31 | +
|
| 32 | +## TODO |
| 33 | +
|
| 34 | +- `bernoulli_poly_eval_one_neg` : $$ B_n(1 - x) = (-1)^n*B_n(x) $$ |
| 35 | +- ``bernoulli_poly_eval_one` : Follows as a consequence of `bernoulli_poly_eval_one_neg`. |
| 36 | +
|
| 37 | +-/ |
| 38 | + |
| 39 | +noncomputable theory |
| 40 | +open_locale big_operators |
| 41 | +open_locale nat |
| 42 | + |
| 43 | +open nat finset |
| 44 | + |
| 45 | +/-- The Bernoulli polynomials are defined in terms of the negative Bernoulli numbers. -/ |
| 46 | +def bernoulli_poly (n : ℕ) : polynomial ℚ := |
| 47 | + ∑ i in range (n + 1), polynomial.monomial (n - i) ((bernoulli i) * (choose n i)) |
| 48 | + |
| 49 | +lemma bernoulli_poly_def (n : ℕ) : bernoulli_poly n = |
| 50 | + ∑ i in range (n + 1), polynomial.monomial i ((bernoulli (n - i)) * (choose n i)) := |
| 51 | +begin |
| 52 | + rw [←sum_range_reflect, add_succ_sub_one, add_zero, bernoulli_poly], |
| 53 | + apply sum_congr rfl, |
| 54 | + rintros x hx, |
| 55 | + rw mem_range_succ_iff at hx, rw [choose_symm hx, nat.sub_sub_self hx], |
| 56 | +end |
| 57 | + |
| 58 | +namespace bernoulli_poly |
| 59 | + |
| 60 | +/- |
| 61 | +### examples |
| 62 | +-/ |
| 63 | + |
| 64 | +section examples |
| 65 | + |
| 66 | +@[simp] lemma bernoulli_poly_zero : bernoulli_poly 0 = 1 := |
| 67 | +by simp [bernoulli_poly] |
| 68 | + |
| 69 | +@[simp] lemma bernoulli_poly_eval_zero (n : ℕ) : (bernoulli_poly n).eval 0 = bernoulli n := |
| 70 | +begin |
| 71 | + rw [bernoulli_poly, polynomial.eval_finset_sum, sum_range_succ], |
| 72 | + have : ∑ (x : ℕ) in range n, bernoulli x * (n.choose x) * 0 ^ (n - x) = 0, |
| 73 | + { apply sum_eq_zero (λ x hx, _), |
| 74 | + have h : 0 < n - x := nat.sub_pos_of_lt (mem_range.1 hx), |
| 75 | + simp [h] }, |
| 76 | + simp [this], |
| 77 | +end |
| 78 | + |
| 79 | +end examples |
| 80 | + |
| 81 | +@[simp] theorem sum_bernoulli_poly (n : ℕ) : |
| 82 | + ∑ k in range (n + 1), ((n + 1).choose k : ℚ) • bernoulli_poly k = |
| 83 | + polynomial.monomial n (n + 1 : ℚ) := |
| 84 | +begin |
| 85 | + simp_rw [bernoulli_poly_def, finset.smul_sum, finset.range_eq_Ico, ←finset.sum_Ico_Ico_comm, |
| 86 | + finset.sum_Ico_eq_sum_range], |
| 87 | + simp only [cast_succ, nat.add_sub_cancel_left, nat.sub_zero, zero_add, linear_map.map_add], |
| 88 | + simp_rw [polynomial.smul_monomial, mul_comm (bernoulli _) _, smul_eq_mul, ←mul_assoc], |
| 89 | + conv_lhs { apply_congr, skip, conv |
| 90 | + { apply_congr, skip, |
| 91 | + rw [choose_mul ((nat.le_sub_left_iff_add_le (mem_range_le H)).1 (mem_range_le H_1)) |
| 92 | + (le.intro rfl), add_comm x x_1, nat.add_sub_cancel, mul_assoc, mul_comm, ←smul_eq_mul, |
| 93 | + ←polynomial.smul_monomial], }, |
| 94 | + rw [←sum_smul], }, |
| 95 | + rw sum_range_succ, |
| 96 | + simp only [add_right_eq_self, cast_succ, mul_one, cast_one, cast_add, nat.add_sub_cancel_left, |
| 97 | + choose_succ_self_right, one_smul, bernoulli_zero, sum_singleton, zero_add, |
| 98 | + linear_map.map_add, range_one], |
| 99 | + have f : ∀ x ∈ range n, 2 ≤ n + 1 - x, |
| 100 | + { rintros x H, rw [mem_range] at H, |
| 101 | + exact nat.le_sub_left_of_add_le (succ_le_succ H) }, |
| 102 | + apply sum_eq_zero (λ x hx, _), |
| 103 | + rw [sum_bernoulli _ (f x hx), zero_smul], |
| 104 | +end |
| 105 | + |
| 106 | +open power_series |
| 107 | +open polynomial (aeval) |
| 108 | +variables {A : Type*} [comm_ring A] [algebra ℚ A] |
| 109 | + |
| 110 | +-- TODO: define exponential generating functions, and use them here |
| 111 | +-- This name should probably be updated afterwards |
| 112 | + |
| 113 | +/-- The theorem that `∑ Bₙ(t)X^n/n!)(e^X-1)=Xe^{tX}` -/ |
| 114 | +theorem exp_bernoulli_poly' (t : A) : |
| 115 | + mk (λ n, aeval t ((1 / n! : ℚ) • bernoulli_poly n)) * (exp A - 1) = X * rescale t (exp A) := |
| 116 | +begin |
| 117 | + -- check equality of power series by checking coefficients of X^n |
| 118 | + ext n, |
| 119 | + -- n = 0 case solved by `simp` |
| 120 | + cases n, { simp }, |
| 121 | + -- n ≥ 1, the coefficients is a sum to n+2, so use `sum_range_succ` to write as |
| 122 | + -- last term plus sum to n+1 |
| 123 | + rw [coeff_succ_X_mul, coeff_rescale, coeff_exp, coeff_mul, |
| 124 | + nat.sum_antidiagonal_eq_sum_range_succ_mk, sum_range_succ], |
| 125 | + -- last term is zero so kill with `zero_add` |
| 126 | + simp only [ring_hom.map_sub, nat.sub_self, constant_coeff_one, constant_coeff_exp, |
| 127 | + coeff_zero_eq_constant_coeff, mul_zero, sub_self, zero_add], |
| 128 | + -- Let's multiply both sides by (n+1)! (OK because it's a unit) |
| 129 | + set u : units ℚ := ⟨(n+1)!, (n+1)!⁻¹, |
| 130 | + mul_inv_cancel (by exact_mod_cast factorial_ne_zero (n+1)), |
| 131 | + inv_mul_cancel (by exact_mod_cast factorial_ne_zero (n+1))⟩ with hu, |
| 132 | + rw ←units.mul_right_inj (units.map (algebra_map ℚ A).to_monoid_hom u), |
| 133 | + -- now tidy up unit mess and generally do trivial rearrangements |
| 134 | + -- to make RHS (n+1)*t^n |
| 135 | + rw [units.coe_map, mul_left_comm, ring_hom.to_monoid_hom_eq_coe, |
| 136 | + ring_hom.coe_monoid_hom, ←ring_hom.map_mul, hu, units.coe_mk], |
| 137 | + change _ = t^n * algebra_map ℚ A (((n+1)*n! : ℕ)*(1/n!)), |
| 138 | + rw [cast_mul, mul_assoc, mul_one_div_cancel |
| 139 | + (show (n! : ℚ) ≠ 0, from cast_ne_zero.2 (factorial_ne_zero n)), mul_one, mul_comm (t^n), |
| 140 | + ← polynomial.aeval_monomial, cast_add, cast_one], |
| 141 | + -- But this is the RHS of `sum_bernoulli_poly` |
| 142 | + rw [← sum_bernoulli_poly, finset.mul_sum, alg_hom.map_sum], |
| 143 | + -- and now we have to prove a sum is a sum, but all the terms are equal. |
| 144 | + apply finset.sum_congr rfl, |
| 145 | + -- The rest is just trivialities, hampered by the fact that we're coercing |
| 146 | + -- factorials and binomial coefficients between ℕ and ℚ and A. |
| 147 | + intros i hi, |
| 148 | + -- NB prime.choose_eq_factorial_div_factorial' is in the wrong namespace |
| 149 | + -- deal with coefficients of e^X-1 |
| 150 | + simp only [choose_eq_factorial_div_factorial' (mem_range_le hi), coeff_mk, |
| 151 | + if_neg (mem_range_sub_ne_zero hi), one_div, alg_hom.map_smul, coeff_one, units.coe_mk, |
| 152 | + coeff_exp, sub_zero, linear_map.map_sub, algebra.smul_mul_assoc, algebra.smul_def, |
| 153 | + mul_right_comm _ ((aeval t) _), ←mul_assoc, ← ring_hom.map_mul, succ_eq_add_one], |
| 154 | + -- finally cancel the Bernoulli polynomial and the algebra_map |
| 155 | + congr', |
| 156 | + apply congr_arg, |
| 157 | + rw [mul_assoc, div_eq_mul_inv, ← mul_inv'], |
| 158 | +end |
| 159 | + |
| 160 | +end bernoulli_poly |
0 commit comments