Skip to content

Commit

Permalink
feat(archive/100-theorems-list): Stirling, part 2 (#14875)
Browse files Browse the repository at this point in the history
feat(archive/100-theorems-list): Stirling

Part 2

- [x] depends on: #14874 [part 1]

Co-authored-by: fabian-kruse <fabian_kruse@gmx.de>
Co-authored-by: nick-kuhn <46423253+nick-kuhn@users.noreply.github.com>



Co-authored-by: Moritz Firsching <firsching@google.com>
Co-authored-by: nick-kuhn <46423253+nick-kuhn@users.noreply.github.com>
  • Loading branch information
3 people committed Sep 16, 2022
1 parent 104f141 commit fc35e1e
Show file tree
Hide file tree
Showing 2 changed files with 102 additions and 5 deletions.
3 changes: 3 additions & 0 deletions docs/100.yaml
Expand Up @@ -368,6 +368,9 @@
author : Chris Hughes
90:
title : Stirling’s Formula
decls :
- stirling.tendsto_stirling_seq_sqrt_pi
author : mathlib (Moritz Firsching, Fabian Kruse, Nikolas Kuhn, Heather Macbeth)
91:
title : The Triangle Inequality
decl : norm_add_le
Expand Down
104 changes: 99 additions & 5 deletions src/analysis/special_functions/stirling.lean
Expand Up @@ -6,29 +6,36 @@ Authors: Moritz Firsching, Fabian Kruse, Nikolas Kuhn
import analysis.p_series
import analysis.special_functions.log.deriv
import tactic.positivity
import data.real.pi.wallis

/-!
# Stirling's formula
This file proves Stirling's formula for the factorial.
It states that $n!$ grows asymptotically like $\sqrt{2\pi n}(\frac{n}{e})^n$.
TODO: Add Part 2 to complete the proof
## Proof outline
The proof follows: <https://proofwiki.org/wiki/Stirling%27s_Formula>.
We proceed in two parts.
### Part 1
We consider the fraction sequence $a_n$ of fractions $\frac{n!}{\sqrt{2n}(\frac{n}{e})^n}$ and
prove that this sequence converges against a real, positive number $a$. For this the two main
ingredients are
- taking the logarithm of the sequence and
- use the series expansion of $\log(1 + x)$.
### Part 2
We use the fact that the series defined in part 1 converges againt a real number $a$ and prove that
$a = \sqrt{\pi}$. Here the main ingredient is the convergence of the Wallis product.
-/

open_locale topological_space big_operators
open_locale topological_space real big_operators nat
open finset filter nat real

namespace stirling
/-!
### Part 1
https://proofwiki.org/wiki/Stirling%27s_Formula#Part_1
Expand All @@ -39,7 +46,7 @@ Define `stirling_seq n` as $\frac{n!}{\sqrt{2n}/(\frac{n}{e})^n$.
Stirling's formula states that this sequence has limit $\sqrt(π)$.
-/
noncomputable def stirling_seq (n : ℕ) : ℝ :=
n.factorial / (sqrt (2 * n) * (n / exp 1) ^ n)
n! / (sqrt (2 * n) * (n / exp 1) ^ n)

@[simp] lemma stirling_seq_zero : stirling_seq 0 = 0 :=
by rw [stirling_seq, cast_zero, mul_zero, real.sqrt_zero, zero_mul, div_zero]
Expand All @@ -52,9 +59,9 @@ We have the expression
`log (stirling_seq (n + 1)) = log(n + 1)! - 1 / 2 * log(2 * n) - n * log ((n + 1) / e)`.
-/
lemma log_stirling_seq_formula (n : ℕ) : log (stirling_seq n.succ) =
log n.succ.factorial - 1 / 2 * log (2 * n.succ) - n.succ * log (n.succ / exp 1) :=
log n.succ!- 1 / 2 * log (2 * n.succ) - n.succ * log (n.succ / exp 1) :=
begin
have h1 : (0 : ℝ) < n.succ.factorial := cast_pos.mpr n.succ.factorial_pos,
have h1 : (0 : ℝ) < n.succ!:= cast_pos.mpr n.succ.factorial_pos,
have h2 : (0 : ℝ) < (2 * n.succ) := mul_pos two_pos (cast_pos.mpr (succ_pos n)),
have h3 := real.sqrt_pos.mpr h2,
have h4 := pow_pos (div_pos (cast_pos.mpr n.succ_pos ) (exp_pos 1)) n.succ,
Expand Down Expand Up @@ -203,3 +210,90 @@ begin
rw ←filter.tendsto_add_at_top_iff_nat 1,
exact tendsto_at_top_cinfi stirling_seq'_antitone ⟨x, hx'⟩,
end

/-!
### Part 2
https://proofwiki.org/wiki/Stirling%27s_Formula#Part_2
-/

/-- For `n : ℕ`, define `w n` as `2^(4*n) * n!^4 / ((2*n)!^2 * (2*n + 1))` -/
noncomputable def w (n : ℕ) : ℝ :=
(2 ^ (4 * n) * n! ^ 4) / ((2 * n)!^ 2 * (2 * n + 1))

/-- The sequence `w n` converges to `π/2` -/
lemma tendsto_w_at_top: tendsto (λ (n : ℕ), w n) at_top (𝓝 (π/2)) :=
begin
convert tendsto_prod_pi_div_two,
funext n,
induction n with n ih,
{ rw [w, prod_range_zero, cast_zero, mul_zero, pow_zero, one_mul, mul_zero, factorial_zero,
cast_one, one_pow, one_pow, one_mul, mul_zero, zero_add, div_one] },
rw [w, prod_range_succ, ←ih, w, _root_.div_mul_div_comm, _root_.div_mul_div_comm],
refine (div_eq_div_iff (ne_of_gt _) (ne_of_gt _)).mpr _,
{ exact mul_pos (pow_pos (cast_pos.mpr (2 * n.succ).factorial_pos) 2)
(add_pos (mul_pos two_pos (cast_pos.mpr n.succ_pos)) one_pos) },
{ have h : 0 < 2 * (n : ℝ) + 1 :=
add_pos_of_nonneg_of_pos (mul_nonneg zero_le_two n.cast_nonneg) one_pos,
exact mul_pos (mul_pos (pow_pos (cast_pos.mpr (2 * n).factorial_pos) 2) h)
(mul_pos h (add_pos_of_nonneg_of_pos (mul_nonneg zero_le_two n.cast_nonneg) three_pos)) },
{ simp_rw [nat.mul_succ, factorial_succ, succ_eq_add_one, pow_succ],
push_cast,
ring_nf },
end

/-- The sequence `n / (2 * n + 1)` tends to `1/2` -/
lemma tendsto_self_div_two_mul_self_add_one :
tendsto (λ (n : ℕ), (n : ℝ) / (2 * n + 1)) at_top (𝓝 (1 / 2)) :=
begin
conv { congr, skip, skip, rw [one_div, ←add_zero (2 : ℝ)] },
refine (((tendsto_const_div_at_top_nhds_0_nat 1).const_add (2 : ℝ)).inv₀
((add_zero (2 : ℝ)).symm ▸ two_ne_zero)).congr' (eventually_at_top.mpr ⟨1, λ n hn, _⟩),
rw [add_div' (1 : ℝ) (2 : ℝ) (n : ℝ) (cast_ne_zero.mpr (one_le_iff_ne_zero.mp hn)), inv_div],
end


/-- For any `n ≠ 0`, we have the identity
`(stirling_seq n)^4/(stirling_seq (2*n))^2 * (n / (2 * n + 1)) = w n`. -/
lemma stirling_seq_pow_four_div_stirling_seq_pow_two_eq (n : ℕ) (hn : n ≠ 0) :
((stirling_seq n) ^ 4 / (stirling_seq (2 * n)) ^ 2) * (n / (2 * n + 1)) = w n :=
begin
rw [bit0_eq_two_mul, stirling_seq, pow_mul, stirling_seq, w],
simp_rw [div_pow, mul_pow],
rw [sq_sqrt (mul_nonneg two_pos.le n.cast_nonneg),
sq_sqrt (mul_nonneg two_pos.le (2 * n).cast_nonneg)],
have : (n : ℝ) ≠ 0, from cast_ne_zero.mpr hn,
have : (exp 1) ≠ 0, from exp_ne_zero 1,
have : ((2 * n)!: ℝ) ≠ 0, from cast_ne_zero.mpr (factorial_ne_zero (2 * n)),
have : 2 * (n : ℝ) + 10, by {norm_cast, exact succ_ne_zero (2*n)},
field_simp,
simp only [mul_pow, mul_comm 2 n, mul_comm 4 n, pow_mul],
ring,
end

/--
Suppose the sequence `stirling_seq` (defined above) has the limit `a ≠ 0`.
Then the sequence `w` has limit `a^2/2`
-/
lemma second_wallis_limit (a : ℝ) (hane : a ≠ 0) (ha : tendsto stirling_seq at_top (𝓝 a)) :
tendsto w at_top (𝓝 (a ^ 2 / 2)):=
begin
refine tendsto.congr' (eventually_at_top.mpr ⟨1, λ n hn,
stirling_seq_pow_four_div_stirling_seq_pow_two_eq n (one_le_iff_ne_zero.mp hn)⟩) _,
have h : a ^ 2 / 2 = (a ^ 4 / a ^ 2) * (1 / 2),
{ rw [mul_one_div, ←mul_one_div (a ^ 4) (a ^ 2), one_div, ←pow_sub_of_lt a],
norm_num },
rw h,
exact ((ha.pow 4).div ((ha.comp (tendsto_id.const_mul_at_top' two_pos)).pow 2)
(pow_ne_zero 2 hane)).mul tendsto_self_div_two_mul_self_add_one,
end

/-- **Stirling's Formula** -/
theorem tendsto_stirling_seq_sqrt_pi : tendsto (λ (n : ℕ), stirling_seq n) at_top (𝓝 (sqrt π)) :=
begin
obtain ⟨a, hapos, halimit⟩ := stirling_seq_has_pos_limit_a,
have hπ : π / 2 = a ^ 2 / 2 := tendsto_nhds_unique tendsto_w_at_top
(second_wallis_limit a (ne_of_gt hapos) halimit),
rwa [(div_left_inj' (show (2 : ℝ) ≠ 0, from two_ne_zero)).mp hπ, sqrt_sq hapos.le],
end

end stirling

0 comments on commit fc35e1e

Please sign in to comment.