section ‹Harmonic Numbers›
theory Harmonic_Numbers
imports
Complex_Transcendental
Summation
Integral_Test
begin
text ‹
The definition of the Harmonic Numbers and the Euler-Mascheroni constant.
Also provides a reasonably accurate approximation of @{term "ln 2 :: real"}
and the Euler-Mascheroni constant.
›
lemma ln_2_less_1: "ln 2 < (1::real)"
proof -
have "2 < 5/(2::real)" by simp
also have "5/2 ≤ exp (1::real)" using exp_lower_taylor_quadratic[of 1, simplified] by simp
finally have "exp (ln 2) < exp (1::real)" by simp
thus "ln 2 < (1::real)" by (subst (asm) exp_less_cancel_iff) simp
qed
lemma setsum_Suc_diff':
fixes f :: "nat ⇒ 'a::ab_group_add"
assumes "m ≤ n"
shows "(∑i = m..<n. f (Suc i) - f i) = f n - f m"
using assms by (induct n) (auto simp: le_Suc_eq)
subsection ‹The Harmonic numbers›
definition harm :: "nat ⇒ 'a :: real_normed_field" where
"harm n = (∑k=1..n. inverse (of_nat k))"
lemma harm_altdef: "harm n = (∑k<n. inverse (of_nat (Suc k)))"
unfolding harm_def by (induction n) simp_all
lemma harm_Suc: "harm (Suc n) = harm n + inverse (of_nat (Suc n))"
by (simp add: harm_def)
lemma harm_nonneg: "harm n ≥ (0 :: 'a :: {real_normed_field,linordered_field})"
unfolding harm_def by (intro setsum_nonneg) simp_all
lemma harm_pos: "n > 0 ⟹ harm n > (0 :: 'a :: {real_normed_field,linordered_field})"
unfolding harm_def by (intro setsum_pos) simp_all
lemma of_real_harm: "of_real (harm n) = harm n"
unfolding harm_def by simp
lemma norm_harm: "norm (harm n) = harm n"
by (subst of_real_harm [symmetric]) (simp add: harm_nonneg)
lemma harm_expand:
"harm 0 = 0"
"harm (Suc 0) = 1"
"harm (numeral n) = harm (pred_numeral n) + inverse (numeral n)"
proof -
have "numeral n = Suc (pred_numeral n)" by simp
also have "harm … = harm (pred_numeral n) + inverse (numeral n)"
by (subst harm_Suc, subst numeral_eq_Suc[symmetric]) simp
finally show "harm (numeral n) = harm (pred_numeral n) + inverse (numeral n)" .
qed (simp_all add: harm_def)
lemma not_convergent_harm: "¬convergent (harm :: nat ⇒ 'a :: real_normed_field)"
proof -
have "convergent (λn. norm (harm n :: 'a)) ⟷
convergent (harm :: nat ⇒ real)" by (simp add: norm_harm)
also have "… ⟷ convergent (λn. ∑k=Suc 0..Suc n. inverse (of_nat k) :: real)"
unfolding harm_def[abs_def] by (subst convergent_Suc_iff) simp_all
also have "... ⟷ convergent (λn. ∑k≤n. inverse (of_nat (Suc k)) :: real)"
by (subst setsum_shift_bounds_cl_Suc_ivl) (simp add: atLeast0AtMost)
also have "... ⟷ summable (λn. inverse (of_nat n) :: real)"
by (subst summable_Suc_iff [symmetric]) (simp add: summable_iff_convergent')
also have "¬..." by (rule not_summable_harmonic)
finally show ?thesis by (blast dest: convergent_norm)
qed
lemma harm_pos_iff [simp]: "harm n > (0 :: 'a :: {real_normed_field,linordered_field}) ⟷ n > 0"
by (rule iffI, cases n, simp add: harm_expand, simp, rule harm_pos)
lemma ln_diff_le_inverse:
assumes "x ≥ (1::real)"
shows "ln (x + 1) - ln x < 1 / x"
proof -
from assms have "∃z>x. z < x + 1 ∧ ln (x + 1) - ln x = (x + 1 - x) * inverse z"
by (intro MVT2) (auto intro!: derivative_eq_intros simp: field_simps)
then obtain z where z: "z > x" "z < x + 1" "ln (x + 1) - ln x = inverse z" by auto
have "ln (x + 1) - ln x = inverse z" by fact
also from z(1,2) assms have "… < 1 / x" by (simp add: field_simps)
finally show ?thesis .
qed
lemma ln_le_harm: "ln (real n + 1) ≤ (harm n :: real)"
proof (induction n)
fix n assume IH: "ln (real n + 1) ≤ harm n"
have "ln (real (Suc n) + 1) = ln (real n + 1) + (ln (real n + 2) - ln (real n + 1))" by simp
also have "(ln (real n + 2) - ln (real n + 1)) ≤ 1 / real (Suc n)"
using ln_diff_le_inverse[of "real n + 1"] by (simp add: add_ac)
also note IH
also have "harm n + 1 / real (Suc n) = harm (Suc n)" by (simp add: harm_Suc field_simps)
finally show "ln (real (Suc n) + 1) ≤ harm (Suc n)" by - simp
qed (simp_all add: harm_def)
subsection ‹The Euler--Mascheroni constant›
text ‹
The limit of the difference between the partial harmonic sum and the natural logarithm
(approximately 0.577216). This value occurs e.g. in the definition of the Gamma function.
›
definition euler_mascheroni :: "'a :: real_normed_algebra_1" where
"euler_mascheroni = of_real (lim (λn. harm n - ln (of_nat n)))"
lemma of_real_euler_mascheroni [simp]: "of_real euler_mascheroni = euler_mascheroni"
by (simp add: euler_mascheroni_def)
interpretation euler_mascheroni: antimono_fun_sum_integral_diff "λx. inverse (x + 1)"
by unfold_locales (auto intro!: continuous_intros)
lemma euler_mascheroni_sum_integral_diff_series:
"euler_mascheroni.sum_integral_diff_series n = harm (Suc n) - ln (of_nat (Suc n))"
proof -
have "harm (Suc n) = (∑k=0..n. inverse (of_nat k + 1) :: real)" unfolding harm_def
unfolding One_nat_def by (subst setsum_shift_bounds_cl_Suc_ivl) (simp add: add_ac)
moreover have "((λx. inverse (x + 1) :: real) has_integral ln (of_nat n + 1) - ln (0 + 1))
{0..of_nat n}"
by (intro fundamental_theorem_of_calculus)
(auto intro!: derivative_eq_intros simp: divide_inverse
has_field_derivative_iff_has_vector_derivative[symmetric])
hence "integral {0..of_nat n} (λx. inverse (x + 1) :: real) = ln (of_nat (Suc n))"
by (auto dest!: integral_unique)
ultimately show ?thesis
by (simp add: euler_mascheroni.sum_integral_diff_series_def atLeast0AtMost)
qed
lemma euler_mascheroni_sequence_decreasing:
"m > 0 ⟹ m ≤ n ⟹ harm n - ln (of_nat n) ≤ harm m - ln (of_nat m :: real)"
by (cases m, simp, cases n, simp, hypsubst,
subst (1 2) euler_mascheroni_sum_integral_diff_series [symmetric],
rule euler_mascheroni.sum_integral_diff_series_antimono, simp)
lemma euler_mascheroni_sequence_nonneg:
"n > 0 ⟹ harm n - ln (of_nat n) ≥ (0::real)"
by (cases n, simp, hypsubst, subst euler_mascheroni_sum_integral_diff_series [symmetric],
rule euler_mascheroni.sum_integral_diff_series_nonneg)
lemma euler_mascheroni_convergent: "convergent (λn. harm n - ln (of_nat n) :: real)"
proof -
have A: "(λn. harm (Suc n) - ln (of_nat (Suc n))) =
euler_mascheroni.sum_integral_diff_series"
by (subst euler_mascheroni_sum_integral_diff_series [symmetric]) (rule refl)
have "convergent (λn. harm (Suc n) - ln (of_nat (Suc n) :: real))"
by (subst A) (fact euler_mascheroni.sum_integral_diff_series_convergent)
thus ?thesis by (subst (asm) convergent_Suc_iff)
qed
lemma euler_mascheroni_LIMSEQ:
"(λn. harm n - ln (of_nat n) :: real) ⇢ euler_mascheroni"
unfolding euler_mascheroni_def
by (simp add: convergent_LIMSEQ_iff [symmetric] euler_mascheroni_convergent)
lemma euler_mascheroni_LIMSEQ_of_real:
"(λn. of_real (harm n - ln (of_nat n))) ⇢
(euler_mascheroni :: 'a :: {real_normed_algebra_1, topological_space})"
proof -
have "(λn. of_real (harm n - ln (of_nat n))) ⇢ (of_real (euler_mascheroni) :: 'a)"
by (intro tendsto_of_real euler_mascheroni_LIMSEQ)
thus ?thesis by simp
qed
lemma euler_mascheroni_sum:
"(λn. inverse (of_nat (n+1)) + ln (of_nat (n+1)) - ln (of_nat (n+2)) :: real)
sums euler_mascheroni"
using sums_add[OF telescope_sums[OF LIMSEQ_Suc[OF euler_mascheroni_LIMSEQ]]
telescope_sums'[OF LIMSEQ_inverse_real_of_nat]]
by (simp_all add: harm_def algebra_simps)
lemma alternating_harmonic_series_sums: "(λk. (-1)^k / real_of_nat (Suc k)) sums ln 2"
proof -
let ?f = "λn. harm n - ln (real_of_nat n)"
let ?g = "λn. if even n then 0 else (2::real)"
let ?em = "λn. harm n - ln (real_of_nat n)"
have "eventually (λn. ?em (2*n) - ?em n + ln 2 = (∑k<2*n. (-1)^k / real_of_nat (Suc k))) at_top"
using eventually_gt_at_top[of "0::nat"]
proof eventually_elim
fix n :: nat assume n: "n > 0"
have "(∑k<2*n. (-1)^k / real_of_nat (Suc k)) =
(∑k<2*n. ((-1)^k + ?g k) / of_nat (Suc k)) - (∑k<2*n. ?g k / of_nat (Suc k))"
by (simp add: setsum.distrib algebra_simps divide_inverse)
also have "(∑k<2*n. ((-1)^k + ?g k) / real_of_nat (Suc k)) = harm (2*n)"
unfolding harm_altdef by (intro setsum.cong) (auto simp: field_simps)
also have "(∑k<2*n. ?g k / real_of_nat (Suc k)) = (∑k|k<2*n ∧ odd k. ?g k / of_nat (Suc k))"
by (intro setsum.mono_neutral_right) auto
also have "… = (∑k|k<2*n ∧ odd k. 2 / (real_of_nat (Suc k)))"
by (intro setsum.cong) auto
also have "(∑k|k<2*n ∧ odd k. 2 / (real_of_nat (Suc k))) = harm n"
unfolding harm_altdef
by (intro setsum.reindex_cong[of "λn. 2*n+1"]) (auto simp: inj_on_def field_simps elim!: oddE)
also have "harm (2*n) - harm n = ?em (2*n) - ?em n + ln 2" using n
by (simp_all add: algebra_simps ln_mult)
finally show "?em (2*n) - ?em n + ln 2 = (∑k<2*n. (-1)^k / real_of_nat (Suc k))" ..
qed
moreover have "(λn. ?em (2*n) - ?em n + ln (2::real))
⇢ euler_mascheroni - euler_mascheroni + ln 2"
by (intro tendsto_intros euler_mascheroni_LIMSEQ filterlim_compose[OF euler_mascheroni_LIMSEQ]
filterlim_subseq) (auto simp: subseq_def)
hence "(λn. ?em (2*n) - ?em n + ln (2::real)) ⇢ ln 2" by simp
ultimately have "(λn. (∑k<2*n. (-1)^k / real_of_nat (Suc k))) ⇢ ln 2"
by (rule Lim_transform_eventually)
moreover have "summable (λk. (-1)^k * inverse (real_of_nat (Suc k)))"
using LIMSEQ_inverse_real_of_nat
by (intro summable_Leibniz(1) decseq_imp_monoseq decseq_SucI) simp_all
hence A: "(λn. ∑k<n. (-1)^k / real_of_nat (Suc k)) ⇢ (∑k. (-1)^k / real_of_nat (Suc k))"
by (simp add: summable_sums_iff divide_inverse sums_def)
from filterlim_compose[OF this filterlim_subseq[of "op * (2::nat)"]]
have "(λn. ∑k<2*n. (-1)^k / real_of_nat (Suc k)) ⇢ (∑k. (-1)^k / real_of_nat (Suc k))"
by (simp add: subseq_def)
ultimately have "(∑k. (- 1) ^ k / real_of_nat (Suc k)) = ln 2" by (intro LIMSEQ_unique)
with A show ?thesis by (simp add: sums_def)
qed
lemma alternating_harmonic_series_sums':
"(λk. inverse (real_of_nat (2*k+1)) - inverse (real_of_nat (2*k+2))) sums ln 2"
unfolding sums_def
proof (rule Lim_transform_eventually)
show "(λn. ∑k<2*n. (-1)^k / (real_of_nat (Suc k))) ⇢ ln 2"
using alternating_harmonic_series_sums unfolding sums_def
by (rule filterlim_compose) (rule mult_nat_left_at_top, simp)
show "eventually (λn. (∑k<2*n. (-1)^k / (real_of_nat (Suc k))) =
(∑k<n. inverse (real_of_nat (2*k+1)) - inverse (real_of_nat (2*k+2)))) sequentially"
proof (intro always_eventually allI)
fix n :: nat
show "(∑k<2*n. (-1)^k / (real_of_nat (Suc k))) =
(∑k<n. inverse (real_of_nat (2*k+1)) - inverse (real_of_nat (2*k+2)))"
by (induction n) (simp_all add: inverse_eq_divide)
qed
qed
subsection ‹Bounds on the Euler--Mascheroni constant›
lemma ln_inverse_approx_le:
assumes "(x::real) > 0" "a > 0"
shows "ln (x + a) - ln x ≤ a * (inverse x + inverse (x + a))/2" (is "_ ≤ ?A")
proof -
def f' ≡ "(inverse (x + a) - inverse x)/a"
have f'_nonpos: "f' ≤ 0" using assms by (simp add: f'_def divide_simps)
let ?f = "λt. (t - x) * f' + inverse x"
let ?F = "λt. (t - x)^2 * f' / 2 + t * inverse x"
have diff: "∀t∈{x..x+a}. (?F has_vector_derivative ?f t)
(at t within {x..x+a})" using assms
by (auto intro!: derivative_eq_intros
simp: has_field_derivative_iff_has_vector_derivative[symmetric])
from assms have "(?f has_integral (?F (x+a) - ?F x)) {x..x+a}"
by (intro fundamental_theorem_of_calculus[OF _ diff])
(auto simp: has_field_derivative_iff_has_vector_derivative[symmetric] field_simps
intro!: derivative_eq_intros)
also have "?F (x+a) - ?F x = (a*2 + f'*a⇧2*x) / (2*x)" using assms by (simp add: field_simps)
also have "f'*a^2 = - (a^2) / (x*(x + a))" using assms
by (simp add: divide_simps f'_def power2_eq_square)
also have "(a*2 + - a⇧2/(x*(x+a))*x) / (2*x) = ?A" using assms
by (simp add: divide_simps power2_eq_square) (simp add: algebra_simps)
finally have int1: "((λt. (t - x) * f' + inverse x) has_integral ?A) {x..x + a}" .
from assms have int2: "(inverse has_integral (ln (x + a) - ln x)) {x..x+a}"
by (intro fundamental_theorem_of_calculus)
(auto simp: has_field_derivative_iff_has_vector_derivative[symmetric] divide_simps
intro!: derivative_eq_intros)
hence "ln (x + a) - ln x = integral {x..x+a} inverse" by (simp add: integral_unique)
also have ineq: "∀xa∈{x..x + a}. inverse xa ≤ (xa - x) * f' + inverse x"
proof
fix t assume t': "t ∈ {x..x+a}"
with assms have t: "0 ≤ (t - x) / a" "(t - x) / a ≤ 1" by simp_all
have "inverse t = inverse ((1 - (t - x) / a) *⇩R x + ((t - x) / a) *⇩R (x + a))" (is "_ = ?A")
using assms t' by (simp add: field_simps)
also from assms have "convex_on {x..x+a} inverse" by (intro convex_on_inverse) auto
from convex_onD_Icc[OF this _ t] assms
have "?A ≤ (1 - (t - x) / a) * inverse x + (t - x) / a * inverse (x + a)" by simp
also have "… = (t - x) * f' + inverse x" using assms
by (simp add: f'_def divide_simps) (simp add: f'_def field_simps)
finally show "inverse t ≤ (t - x) * f' + inverse x" .
qed
hence "integral {x..x+a} inverse ≤ integral {x..x+a} ?f" using f'_nonpos assms
by (intro integral_le has_integral_integrable[OF int1] has_integral_integrable[OF int2] ineq)
also have "… = ?A" using int1 by (rule integral_unique)
finally show ?thesis .
qed
lemma ln_inverse_approx_ge:
assumes "(x::real) > 0" "x < y"
shows "ln y - ln x ≥ 2 * (y - x) / (x + y)" (is "_ ≥ ?A")
proof -
def m ≡ "(x+y)/2"
def f' ≡ "-inverse (m^2)"
from assms have m: "m > 0" by (simp add: m_def)
let ?F = "λt. (t - m)^2 * f' / 2 + t / m"
from assms have "((λt. (t - m) * f' + inverse m) has_integral (?F y - ?F x)) {x..y}"
by (intro fundamental_theorem_of_calculus)
(auto simp: has_field_derivative_iff_has_vector_derivative[symmetric] divide_simps
intro!: derivative_eq_intros)
also from m have "?F y - ?F x = ((y - m)^2 - (x - m)^2) * f' / 2 + (y - x) / m"
by (simp add: field_simps)
also have "((y - m)^2 - (x - m)^2) = 0" by (simp add: m_def power2_eq_square field_simps)
also have "0 * f' / 2 + (y - x) / m = ?A" by (simp add: m_def)
finally have int1: "((λt. (t - m) * f' + inverse m) has_integral ?A) {x..y}" .
from assms have int2: "(inverse has_integral (ln y - ln x)) {x..y}"
by (intro fundamental_theorem_of_calculus)
(auto simp: has_field_derivative_iff_has_vector_derivative[symmetric] divide_simps
intro!: derivative_eq_intros)
hence "ln y - ln x = integral {x..y} inverse" by (simp add: integral_unique)
also have ineq: "∀xa∈{x..y}. inverse xa ≥ (xa - m) * f' + inverse m"
proof
fix t assume t: "t ∈ {x..y}"
from t assms have "inverse t - inverse m ≥ f' * (t - m)"
by (intro convex_on_imp_above_tangent[of "{0<..}"] convex_on_inverse)
(auto simp: m_def interior_open f'_def power2_eq_square intro!: derivative_eq_intros)
thus "(t - m) * f' + inverse m ≤ inverse t" by (simp add: algebra_simps)
qed
hence "integral {x..y} inverse ≥ integral {x..y} (λt. (t - m) * f' + inverse m)"
using int1 int2 by (intro integral_le has_integral_integrable)
also have "integral {x..y} (λt. (t - m) * f' + inverse m) = ?A"
using integral_unique[OF int1] by simp
finally show ?thesis .
qed
lemma euler_mascheroni_lower:
"euler_mascheroni ≥ harm (Suc n) - ln (real_of_nat (n + 2)) + 1/real_of_nat (2 * (n + 2))"
and euler_mascheroni_upper:
"euler_mascheroni ≤ harm (Suc n) - ln (real_of_nat (n + 2)) + 1/real_of_nat (2 * (n + 1))"
proof -
def D ≡ "λn. inverse (of_nat (n+1)) + ln (of_nat (n+1)) - ln (of_nat (n+2)) :: real"
let ?g = "λn. ln (of_nat (n+2)) - ln (of_nat (n+1)) - inverse (of_nat (n+1)) :: real"
def inv ≡ "λn. inverse (real_of_nat n)"
fix n :: nat
note summable = sums_summable[OF euler_mascheroni_sum, folded D_def]
have sums: "(λk. (inv (Suc (k + (n+1))) - inv (Suc (Suc k + (n+1))))/2) sums ((inv (Suc (0 + (n+1))) - 0)/2)"
unfolding inv_def
by (intro sums_divide telescope_sums' LIMSEQ_ignore_initial_segment LIMSEQ_inverse_real_of_nat)
have sums': "(λk. (inv (Suc (k + n)) - inv (Suc (Suc k + n)))/2) sums ((inv (Suc (0 + n)) - 0)/2)"
unfolding inv_def
by (intro sums_divide telescope_sums' LIMSEQ_ignore_initial_segment LIMSEQ_inverse_real_of_nat)
from euler_mascheroni_sum have "euler_mascheroni = (∑k. D k)"
by (simp add: sums_iff D_def)
also have "… = (∑k. D (k + Suc n)) + (∑k≤n. D k)"
by (subst suminf_split_initial_segment[OF summable, of "Suc n"], subst lessThan_Suc_atMost) simp
finally have sum: "(∑k≤n. D k) - euler_mascheroni = -(∑k. D (k + Suc n))" by simp
note sum
also have "… ≤ -(∑k. (inv (k + Suc n + 1) - inv (k + Suc n + 2)) / 2)"
proof (intro le_imp_neg_le suminf_le allI summable_ignore_initial_segment[OF summable])
fix k' :: nat
def k ≡ "k' + Suc n"
hence k: "k > 0" by (simp add: k_def)
have "real_of_nat (k+1) > 0" by (simp add: k_def)
with ln_inverse_approx_le[OF this zero_less_one]
have "ln (of_nat k + 2) - ln (of_nat k + 1) ≤ (inv (k+1) + inv (k+2))/2"
by (simp add: inv_def add_ac)
hence "(inv (k+1) - inv (k+2))/2 ≤ inv (k+1) + ln (of_nat (k+1)) - ln (of_nat (k+2))"
by (simp add: field_simps)
also have "… = D k" unfolding D_def inv_def ..
finally show "D (k' + Suc n) ≥ (inv (k' + Suc n + 1) - inv (k' + Suc n + 2)) / 2"
by (simp add: k_def)
from sums_summable[OF sums]
show "summable (λk. (inv (k + Suc n + 1) - inv (k + Suc n + 2))/2)" by simp
qed
also from sums have "… = -inv (n+2) / 2" by (simp add: sums_iff)
finally have "euler_mascheroni ≥ (∑k≤n. D k) + 1 / (of_nat (2 * (n+2)))"
by (simp add: inv_def field_simps)
also have "(∑k≤n. D k) = harm (Suc n) - (∑k≤n. ln (real_of_nat (Suc k+1)) - ln (of_nat (k+1)))"
unfolding harm_altdef D_def by (subst lessThan_Suc_atMost) (simp add: setsum.distrib setsum_subtractf)
also have "(∑k≤n. ln (real_of_nat (Suc k+1)) - ln (of_nat (k+1))) = ln (of_nat (n+2))"
by (subst atLeast0AtMost [symmetric], subst setsum_Suc_diff) simp_all
finally show "euler_mascheroni ≥ harm (Suc n) - ln (real_of_nat (n + 2)) + 1/real_of_nat (2 * (n + 2))"
by simp
note sum
also have "-(∑k. D (k + Suc n)) ≥ -(∑k. (inv (Suc (k + n)) - inv (Suc (Suc k + n)))/2)"
proof (intro le_imp_neg_le suminf_le allI summable_ignore_initial_segment[OF summable])
fix k' :: nat
def k ≡ "k' + Suc n"
hence k: "k > 0" by (simp add: k_def)
have "real_of_nat (k+1) > 0" by (simp add: k_def)
from ln_inverse_approx_ge[of "of_nat k + 1" "of_nat k + 2"]
have "2 / (2 * real_of_nat k + 3) ≤ ln (of_nat (k+2)) - ln (real_of_nat (k+1))"
by (simp add: add_ac)
hence "D k ≤ 1 / real_of_nat (k+1) - 2 / (2 * real_of_nat k + 3)"
by (simp add: D_def inverse_eq_divide inv_def)
also have "… = inv ((k+1)*(2*k+3))" unfolding inv_def by (simp add: field_simps)
also have "… ≤ inv (2*k*(k+1))" unfolding inv_def using k
by (intro le_imp_inverse_le)
(simp add: algebra_simps, simp del: of_nat_add)
also have "… = (inv k - inv (k+1))/2" unfolding inv_def using k
by (simp add: divide_simps del: of_nat_mult) (simp add: algebra_simps)
finally show "D k ≤ (inv (Suc (k' + n)) - inv (Suc (Suc k' + n)))/2" unfolding k_def by simp
next
from sums_summable[OF sums']
show "summable (λk. (inv (Suc (k + n)) - inv (Suc (Suc k + n)))/2)" by simp
qed
also from sums' have "(∑k. (inv (Suc (k + n)) - inv (Suc (Suc k + n)))/2) = inv (n+1)/2"
by (simp add: sums_iff)
finally have "euler_mascheroni ≤ (∑k≤n. D k) + 1 / of_nat (2 * (n+1))"
by (simp add: inv_def field_simps)
also have "(∑k≤n. D k) = harm (Suc n) - (∑k≤n. ln (real_of_nat (Suc k+1)) - ln (of_nat (k+1)))"
unfolding harm_altdef D_def by (subst lessThan_Suc_atMost) (simp add: setsum.distrib setsum_subtractf)
also have "(∑k≤n. ln (real_of_nat (Suc k+1)) - ln (of_nat (k+1))) = ln (of_nat (n+2))"
by (subst atLeast0AtMost [symmetric], subst setsum_Suc_diff) simp_all
finally show "euler_mascheroni ≤ harm (Suc n) - ln (real_of_nat (n + 2)) + 1/real_of_nat (2 * (n + 1))"
by simp
qed
lemma euler_mascheroni_pos: "euler_mascheroni > (0::real)"
using euler_mascheroni_lower[of 0] ln_2_less_1 by (simp add: harm_def)
context
begin
private lemma ln_approx_aux:
fixes n :: nat and x :: real
defines "y ≡ (x-1)/(x+1)"
assumes x: "x > 0" "x ≠ 1"
shows "inverse (2*y^(2*n+1)) * (ln x - (∑k<n. 2*y^(2*k+1) / of_nat (2*k+1))) ∈
{0..(1 / (1 - y^2) / of_nat (2*n+1))}"
proof -
from x have norm_y: "norm y < 1" unfolding y_def by simp
from power_strict_mono[OF this, of 2] have norm_y': "norm y^2 < 1" by simp
let ?f = "λk. 2 * y ^ (2*k+1) / of_nat (2*k+1)"
note sums = ln_series_quadratic[OF x(1)]
def c ≡ "inverse (2*y^(2*n+1))"
let ?d = "c * (ln x - (∑k<n. ?f k))"
have "∀k. y⇧2^k / of_nat (2*(k+n)+1) ≤ y⇧2 ^ k / of_nat (2*n+1)"
by (intro allI divide_left_mono mult_right_mono mult_pos_pos zero_le_power[of "y^2"]) simp_all
moreover {
have "(λk. ?f (k + n)) sums (ln x - (∑k<n. ?f k))"
using sums_split_initial_segment[OF sums] by (simp add: y_def)
hence "(λk. c * ?f (k + n)) sums ?d" by (rule sums_mult)
also have "(λk. c * (2*y^(2*(k+n)+1) / of_nat (2*(k+n)+1))) =
(λk. (c * (2*y^(2*n+1))) * ((y^2)^k / of_nat (2*(k+n)+1)))"
by (simp only: ring_distribs power_add power_mult) (simp add: mult_ac)
also from x have "c * (2*y^(2*n+1)) = 1" by (simp add: c_def y_def)
finally have "(λk. (y^2)^k / of_nat (2*(k+n)+1)) sums ?d" by simp
} note sums' = this
moreover from norm_y' have "(λk. (y^2)^k / of_nat (2*n+1)) sums (1 / (1 - y^2) / of_nat (2*n+1))"
by (intro sums_divide geometric_sums) (simp_all add: norm_power)
ultimately have "?d ≤ (1 / (1 - y^2) / of_nat (2*n+1))" by (rule sums_le)
moreover have "c * (ln x - (∑k<n. 2 * y ^ (2 * k + 1) / real_of_nat (2 * k + 1))) ≥ 0"
by (intro sums_le[OF _ sums_zero sums']) simp_all
ultimately show ?thesis unfolding c_def by simp
qed
lemma
fixes n :: nat and x :: real
defines "y ≡ (x-1)/(x+1)"
defines "approx ≡ (∑k<n. 2*y^(2*k+1) / of_nat (2*k+1))"
defines "d ≡ y^(2*n+1) / (1 - y^2) / of_nat (2*n+1)"
assumes x: "x > 1"
shows ln_approx_bounds: "ln x ∈ {approx..approx + 2*d}"
and ln_approx_abs: "abs (ln x - (approx + d)) ≤ d"
proof -
def c ≡ "2*y^(2*n+1)"
from x have c_pos: "c > 0" unfolding c_def y_def
by (intro mult_pos_pos zero_less_power) simp_all
have A: "inverse c * (ln x - (∑k<n. 2*y^(2*k+1) / of_nat (2*k+1))) ∈
{0.. (1 / (1 - y^2) / of_nat (2*n+1))}" using assms unfolding y_def c_def
by (intro ln_approx_aux) simp_all
hence "inverse c * (ln x - (∑k<n. 2*y^(2*k+1)/of_nat (2*k+1))) ≤ (1 / (1-y^2) / of_nat (2*n+1))"
by simp
hence "(ln x - (∑k<n. 2*y^(2*k+1) / of_nat (2*k+1))) / c ≤ (1 / (1 - y^2) / of_nat (2*n+1))"
by (auto simp add: divide_simps)
with c_pos have "ln x ≤ c / (1 - y^2) / of_nat (2*n+1) + approx"
by (subst (asm) pos_divide_le_eq) (simp_all add: mult_ac approx_def)
moreover {
from A c_pos have "0 ≤ c * (inverse c * (ln x - (∑k<n. 2*y^(2*k+1) / of_nat (2*k+1))))"
by (intro mult_nonneg_nonneg[of c]) simp_all
also have "… = (c * inverse c) * (ln x - (∑k<n. 2*y^(2*k+1) / of_nat (2*k+1)))"
by (simp add: mult_ac)
also from c_pos have "c * inverse c = 1" by simp
finally have "ln x ≥ approx" by (simp add: approx_def)
}
ultimately show "ln x ∈ {approx..approx + 2*d}" by (simp add: c_def d_def)
thus "abs (ln x - (approx + d)) ≤ d" by auto
qed
end
lemma euler_mascheroni_bounds:
fixes n :: nat assumes "n ≥ 1" defines "t ≡ harm n - ln (of_nat (Suc n)) :: real"
shows "euler_mascheroni ∈ {t + inverse (of_nat (2*(n+1)))..t + inverse (of_nat (2*n))}"
using assms euler_mascheroni_upper[of "n-1"] euler_mascheroni_lower[of "n-1"]
unfolding t_def by (cases n) (simp_all add: harm_Suc t_def inverse_eq_divide)
lemma euler_mascheroni_bounds':
fixes n :: nat assumes "n ≥ 1" "ln (real_of_nat (Suc n)) ∈ {l<..<u}"
shows "euler_mascheroni ∈
{harm n - u + inverse (of_nat (2*(n+1)))<..<harm n - l + inverse (of_nat (2*n))}"
using euler_mascheroni_bounds[OF assms(1)] assms(2) by auto
text ‹
Approximation of @{term "ln 2"}. The lower bound is accurate to about 0.03; the upper
bound is accurate to about 0.0015.
›
lemma ln2_ge_two_thirds: "2/3 ≤ ln (2::real)"
and ln2_le_25_over_36: "ln (2::real) ≤ 25/36"
using ln_approx_bounds[of 2 1, simplified, simplified eval_nat_numeral, simplified] by simp_all
text ‹
Approximation of the Euler--Mascheroni constant. The lower bound is accurate to about 0.0015;
the upper bound is accurate to about 0.015.
›
lemma euler_mascheroni_gt_19_over_33: "(euler_mascheroni :: real) > 19/33" (is ?th1)
and euler_mascheroni_less_13_over_22: "(euler_mascheroni :: real) < 13/22" (is ?th2)
proof -
have "ln (real (Suc 7)) = 3 * ln 2" by (simp add: ln_powr [symmetric] powr_numeral)
also from ln_approx_bounds[of 2 3] have "… ∈ {3*307/443<..<3*4615/6658}"
by (simp add: eval_nat_numeral)
finally have "ln (real (Suc 7)) ∈ …" .
from euler_mascheroni_bounds'[OF _ this] have "?th1 ∧ ?th2" by (simp_all add: harm_expand)
thus ?th1 ?th2 by blast+
qed
end