Theory Prime_Number_Theorem

(*
  File:       Prime Number Theorem.thy
  Authors:    Manuel Eberl (TU München), Larry Paulson (University of Cambridge)

  A proof of the Prime Number Theorem and some related properties
*)
section ‹The Prime Number Theorem›
theory Prime_Number_Theorem
imports 
  Newman_Ingham_Tauberian
  Prime_Counting_Functions
begin

(*<*)
unbundle prime_counting_notation
(*>*)

subsection ‹Constructing Newman's function›

text ‹
  Starting from Mertens' first theorem, i.\,e.\ $\mathfrak M(x) = \ln x + O(1)$, we now 
  want to derive that $\mathfrak M(x) = \ln x + c + o(1)$. This result is considerably stronger
  and it implies the Prime Number Theorem quite directly.

  In order to do this, we define the Dirichlet series
  \[f(s) = \sum_{n=1}^\infty \frac{\mathfrak{M}(n)}{n^s}\ .\]
  We will prove that this series extends meromorphically to $\mathfrak{R}(s)\geq 1$ and
  apply Ingham's theorem to it (after we subtracted its pole at $s = 1$).
›
definition fds_newman where
  "fds_newman = fds (λn. complex_of_real (𝔐 n))"

lemma fds_nth_newman:
  "fds_nth fds_newman n = of_real (𝔐 n)"
  by (simp add: fds_newman_def fds_nth_fds)

lemma norm_fds_nth_newman:
  "norm (fds_nth fds_newman n) = 𝔐 n"
  unfolding fds_nth_newman norm_of_real
  by (intro abs_of_nonneg sum_nonneg divide_nonneg_pos) (auto dest: prime_ge_1_nat)

text ‹
  The Dirichlet series $f(s) + \zeta'(s)$ has the coefficients $\mathfrak{M}(n) - \ln n$,
  so by Mertens' first theorem, $f(s) + \zeta'(s)$ has bounded coefficients.
›
lemma bounded_coeffs_newman_minus_deriv_zeta:
  defines "f  fds_newman + fds_deriv fds_zeta"
  shows   "Bseq (λn. fds_nth f n)"
proof -
  have "(λn. 𝔐 (real n) - ln (real n))  O(λ_. 1)"
    using mertens_bounded by (rule landau_o.big.compose) real_asymp
  from natfun_bigo_1E[OF this, of 1]
    obtain c where c: "c  1" "n. ¦𝔐 (real n) - ln (real n)¦  c" by auto

  show ?thesis
  proof (intro BseqI[of c] allI)
    fix n :: nat
    show "norm (fds_nth f n)  c"
    proof (cases "n = 0")
      case False
      hence "fds_nth f n = of_real (𝔐 n - ln n)"
        by (simp add: f_def fds_nth_newman fds_nth_deriv fds_nth_zeta scaleR_conv_of_real)
      also from n  0 have "norm   c"
        using c(2)[of n] by (simp add: in_Reals_norm)
      finally show ?thesis .
    qed (insert c, auto)
  qed (insert c, auto)
qed

text ‹
  A Dirichlet series with bounded coefficients converges for all $s$ with
  $\mathfrak{R}(s)>1$ and so does $\zeta'(s)$, so we can conclude that $f(s)$ does as well.
›
lemma abs_conv_abscissa_newman: "abs_conv_abscissa fds_newman  1"
  and conv_abscissa_newman:     "conv_abscissa fds_newman  1"
proof -
  define f where "f = fds_newman + fds_deriv fds_zeta"
  have "abs_conv_abscissa f  1"
    using bounded_coeffs_newman_minus_deriv_zeta unfolding f_def
    by (rule bounded_coeffs_imp_abs_conv_abscissa_le_1)
  hence "abs_conv_abscissa (f - fds_deriv fds_zeta)  1"
    by (intro abs_conv_abscissa_diff_leI) (auto simp: abs_conv_abscissa_deriv)
  also have "f - fds_deriv fds_zeta = fds_newman" by (simp add: f_def)
  finally show "abs_conv_abscissa fds_newman  1" .
  from conv_le_abs_conv_abscissa and this show "conv_abscissa fds_newman  1"
    by (rule order.trans)
qed

text ‹
  We now change the order of summation to obtain an alternative form of $f(s)$ in terms of a 
  sum of Hurwitz $\zeta$ functions.
›
lemma eval_fds_newman_conv_infsetsum:
  assumes s: "Re s > 1"
  shows   "eval_fds fds_newman s = (ap | prime p. (ln (real p) / real p) * hurwitz_zeta p s)"
          "(λp. ln (real p) / real p * hurwitz_zeta p s) abs_summable_on {p. prime p}"
proof -
  from s have conv: "fds_abs_converges fds_newman s"
    by (intro fds_abs_converges le_less_trans[OF abs_conv_abscissa_newman]) auto
  define f where "f = (λn p. ln (real p) / real p / of_nat n powr s)"

  have eq: "(an{p..}. f n p) = ln (real p) / real p * hurwitz_zeta p s" if "prime p" for p
  proof -
    have "(an{p..}. f n p) = (ax{p..}. (ln (real p) / of_nat p) * (1 / of_nat x powr s))"
      by (simp add: f_def)
    also have " = (ln (real p) / of_nat p) * (ax{p..}. 1 / of_nat x powr s)"
      using abs_summable_hurwitz_zeta[of s 0 p] that s
      by (intro infsetsum_cmult_right) (auto dest: prime_gt_0_nat)
    also have "(ax{p..}. 1 / of_nat x powr s) = hurwitz_zeta p s"
      using s that by (subst hurwitz_zeta_nat_conv_infsetsum(2))
                      (auto dest: prime_gt_0_nat simp: field_simps powr_minus)
    finally show ?thesis .
  qed

  have norm_f: "norm (f n p) = ln p / p / n powr Re s" if "prime p" for n p :: nat
    by (auto simp: f_def norm_divide norm_mult norm_powr_real_powr)
  from conv have "(λn. norm (fds_nth fds_newman n / n powr s)) abs_summable_on UNIV"
    by (intro abs_summable_on_normI) (simp add: fds_abs_converges_altdef')
  also have "(λn. norm (fds_nth fds_newman n / n powr s)) =
               (λn. p | prime p  p  n. norm (f n p))"
    by (auto simp: norm_divide norm_fds_nth_newman sum_divide_distrib primes_M_def
                   prime_sum_upto_def norm_mult norm_f norm_powr_real_powr intro!: sum.cong)
  finally have summable1: "(λ(n,p). f n p) abs_summable_on (SIGMA n:UNIV. {p. prime p  p  n})"
    using conv by (subst abs_summable_on_Sigma_iff) auto
  also have "?this  (λ(p,n). f n p) abs_summable_on
                         (λ(n,p). (p,n)) ` (SIGMA n:UNIV. {p. prime p  p  n})"
    by (subst abs_summable_on_reindex_iff [symmetric]) (auto simp: case_prod_unfold inj_on_def)
  also have "(λ(n,p). (p,n)) ` (SIGMA n:UNIV. {p. prime p  p  n}) =
               (SIGMA p:{p. prime p}. {p..})" by auto
  finally have summable2: "(λ(p,n). f n p) abs_summable_on " .
  from abs_summable_on_Sigma_project1'[OF this]
    have "(λp. an{p..}. f n p) abs_summable_on {p. prime p}" by auto
  also have "?this  (λp. ln (real p) / real p * hurwitz_zeta p s) abs_summable_on {p. prime p}"
    by (intro abs_summable_on_cong eq) auto
  finally show  .

  have "eval_fds fds_newman s =
          (an. p | prime p  p  n. ln (real p) / real p / of_nat n powr s)"
    using conv by (simp add: eval_fds_altdef fds_nth_newman sum_divide_distrib
                             primes_M_def prime_sum_upto_def)
  also have " = (an. ap | prime p  p  n. f n p)"
    unfolding f_def by (subst infsetsum_finite) auto
  also have " = (a(n, p)  (SIGMA n:UNIV. {p. prime p  p  n}). f n p)"
    using summable1 by (subst infsetsum_Sigma) auto
  also have " = (a(p, n)  (λ(n,p). (p, n)) ` (SIGMA n:UNIV. {p. prime p  p  n}). f n p)"
    by (subst infsetsum_reindex) (auto simp: case_prod_unfold inj_on_def)
  also have "(λ(n,p). (p, n)) ` (SIGMA n:UNIV. {p. prime p  p  n}) =
               (SIGMA p:{p. prime p}. {p..})" by auto
  also have "(a(p,n). f n p) = (ap | prime p. an{p..}. f n p)"
    using summable2 by (subst infsetsum_Sigma) auto
  also have "(ap | prime p. an{p..}. f n p) =
               (ap | prime p. ln (real p) / real p * hurwitz_zeta p s)"
    by (intro infsetsum_cong eq) auto
  finally show "eval_fds fds_newman s =
                  (ap | prime p. (ln (real p) / real p) * hurwitz_zeta p s)" .
qed


text ‹
  We now define a meromorphic continuation of $f(s)$ on $\mathfrak{R}(s) > \frac{1}{2}$.

  To construct $f(s)$, we express it as
  \[f(s) = \frac{1}{z-1}\left(\bar f(s) - \frac{\zeta'(s)}{\zeta(s)}\right)\ ,\]
  where $\bar f(s)$ (which we shall call pre_newman›) is a function that is analytic on
  $\Re(s) > \frac{1}{2}$, which can be shown fairly easily using the Weierstra{\ss} M test.
  
  $\zeta'(s)/\zeta(s)$ is meromorphic except for a single pole at $s = 1$ and one $k$-th order
  pole for any $k$-th order zero of $\zeta$, but for the Prime Number Theorem, we are only
  concerned with the area $\mathfrak{R}(s) \geq 1$, where $\zeta$ does not have any zeros.

  Taken together, this means that $f(s)$ is analytic for $\mathfrak{R}(s)\geq 1$ except for a
  double pole at $s = 1$, which we will take care of later.
›

context
  fixes A :: "nat  complex  complex" and B :: "nat  complex  complex"
  defines "A  (λp s. (s - 1) * pre_zeta (real p) s - 
                         of_nat p / (of_nat p powr s * (of_nat p powr s - 1)))"
  defines "B  (λp s. of_real (ln (real p)) / of_nat p * A p s)"
begin

definition pre_newman :: "complex  complex" where
  "pre_newman s = (p. if prime p then B p s else 0)"

definition newman where "newman s = 1 / (s - 1) * (pre_newman s - deriv zeta s / zeta s)"

text ‹
  The sum used in the definition of pre_newman› converges uniformly on any disc within the
  half-space with $\mathfrak{R}(s) > \frac{1}{2}$ by the Weierstra{\ss} M test.
›
lemma uniform_limit_pre_newman:
  assumes r: "r  0" "Re s - r > 1 / 2"
  shows "uniform_limit (cball s r)
           (λn s. p<n. if prime p then B p s else 0) pre_newman at_top"
proof -
  from r have Re: "Re z > 1 / 2" if "dist s z  r" for z
    using abs_Re_le_cmod[of "s - z"] r that
    by (auto simp: dist_norm abs_if split: if_splits)

  define x where "x = Re s - r" ― ‹The lower bound for the real part in the disc›
  from r Re have "x > 1 / 2" by (auto simp: x_def)

  ― ‹The following sequence M› bounds the summand, and it is obviously $O(n^{-1-\epsilon})$
      and therefore summable›
  define C where "C = (norm s + r + 1) * (norm s + r) / x"
  define M where "M = (λp::nat. ln p * (C / p powr (x + 1) + 1 / (p powr x * (p powr x - 1))))"

  show ?thesis unfolding pre_newman_def
  proof (intro Weierstrass_m_test_ev[OF eventually_mono[OF eventually_gt_at_top[of 1]]] ballI)
    show "summable M"
    proof (rule summable_comparison_test_bigo)
      define ε where "ε = min (2 * x - 1) x / 2"
      from x > 1 / 2 have ε: "ε > 0" "1 + ε < 2 * x" "1 + ε < x + 1"
        by (auto simp: ε_def min_def field_simps)
      show "M  O(λn. n powr (- 1 - ε))" unfolding M_def distrib_left
        by (intro sum_in_bigo) (use ε in real_asymp)+
      from ε show "summable (λn. norm (n powr (- 1 - ε)))"
        by (simp add: summable_real_powr_iff)
    qed
  next
    fix p :: nat and z assume p: "p > 1" and z: "z  cball s r"
    from z r Re[of z] have x: "Re z  x" "x > 1 / 2" and "Re z > 1 / 2"
      using abs_Re_le_cmod[of "s - z"] by (auto simp: x_def algebra_simps dist_norm)
    have norm_z: "norm z  norm s + r"
      using z norm_triangle_ineq2[of z s] r by (auto simp: dist_norm norm_minus_commute)
    from p > 1 and x and r have "M p  0"
      by (auto simp: C_def M_def intro!: mult_nonneg_nonneg add_nonneg_nonneg divide_nonneg_pos)

    have bound: "norm ((z - 1) * pre_zeta p z)  
                   norm (z - 1) * (norm z / (Re z * p powr Re z))"
      using pre_zeta_bound'[of z p] p Re z > 1 / 2
      unfolding norm_mult by (intro mult_mono pre_zeta_bound) auto

    have "norm (B p z) = ln p / p * norm (A p z)"
      unfolding B_def using p > 1 by (simp add: B_def norm_mult norm_divide)
    also have "  ln p / p * (norm (z - 1) * norm z / Re z / p powr Re z + 
                                 p / (p powr Re z * (p powr Re z - 1)))"
      unfolding A_def using p > 1 and Re z > 1 / 2 and bound
      by (intro mult_left_mono order.trans[OF norm_triangle_ineq4 add_mono] mult_left_mono)
         (auto simp: norm_divide norm_mult norm_powr_real_powr
               intro!: divide_left_mono order.trans[OF _ norm_triangle_ineq2])
    also have " = ln p * (norm (z - 1) * norm z / Re z / p powr (Re z + 1) + 
                            1 / (p powr Re z * (p powr Re z - 1)))"
      using p > 1 by (simp add: field_simps powr_add powr_minus)
    also have "norm (z - 1) * norm z / Re z / p powr (Re z + 1)  C / p powr (x + 1)"
      unfolding C_def using r Re z > 1 / 2 norm_z p x
      by (intro mult_mono frac_le powr_mono order.trans[OF norm_triangle_ineq4]) auto
    also have "1 / (p powr Re z * (p powr Re z - 1)) 
                 1 / (p powr x * (p powr x - 1))" using p > 1 x
      by (intro divide_left_mono mult_mono powr_mono diff_right_mono mult_pos_pos)
         (auto simp: ge_one_powr_ge_zero)
    finally have "norm (B p z)  M p"
      using p > 1 by (simp add: mult_left_mono M_def)
    with M p  0 show "norm (if prime p then B p z else 0)  M p" by simp
  qed
qed

lemma sums_pre_newman: "Re s > 1 / 2  (λp. if prime p then B p s else 0) sums pre_newman s"
  using tendsto_uniform_limitI[OF uniform_limit_pre_newman[of 0 s]] by (auto simp: sums_def)

lemma analytic_pre_newman [THEN analytic_on_subset, analytic_intros]:
  "pre_newman analytic_on {s. Re s > 1 / 2}"
proof -
  have holo: "(λs::complex. if prime p then B p s else 0) holomorphic_on X"
    if "X  {s. Re s > 1 / 2}" for X and p :: nat using that
    by (cases "prime p")
       (auto intro!: holomorphic_intros simp: B_def A_def dest!: prime_gt_1_nat)
  have holo': "pre_newman holomorphic_on ball s r" if r: "r  0" "Re s - r > 1 / 2" for s r
  proof -
    from r have Re: "Re z > 1 / 2" if "dist s z  r" for z
      using abs_Re_le_cmod[of "s - z"] r that by (auto simp: dist_norm abs_if split: if_splits)
    show ?thesis
      by (rule holomorphic_uniform_limit[OF _ uniform_limit_pre_newman[of r s]])
         (insert that Re, auto intro!: always_eventually holomorphic_on_imp_continuous_on
                                       holomorphic_intros holo)
  qed
  show ?thesis unfolding analytic_on_def
  proof safe
    fix s assume "Re s > 1 / 2"
    thus "r>0. pre_newman holomorphic_on ball s r"
      by (intro exI[of _ "(Re s - 1 / 2) / 2"] conjI holo') (auto simp: field_simps)
  qed
qed

lemma holomorphic_pre_newman [holomorphic_intros]:
  "X  {s. Re s > 1 / 2}  pre_newman holomorphic_on X"
  using analytic_pre_newman by (rule analytic_imp_holomorphic)

lemma eval_fds_newman:
  assumes s: "Re s > 1"
  shows   "eval_fds fds_newman s = newman s"
proof -
  have eq: "(ln (real p) / real p) * hurwitz_zeta p s =
              1 / (s - 1) * (ln (real p) / (p powr s - 1) + B p s)"
    if p: "prime p" for p
  proof -
    have "(ln (real p) / real p) * hurwitz_zeta p s =
            ln (real p) / real p * (p powr (1 - s) / (s - 1) + pre_zeta p s)"
      using s by (auto simp add: hurwitz_zeta_def)
    also have " = 1 / (s - 1) * (ln (real p) / (p powr s - 1) + B p s)"
      using p s by (simp add: divide_simps powr_diff B_def)
                   (auto simp: A_def field_simps dest: prime_gt_1_nat)?
    finally show ?thesis .
  qed

  have "(λp. (ln (real p) / real p) * hurwitz_zeta p s) abs_summable_on {p. prime p}"
    using s by (intro eval_fds_newman_conv_infsetsum)
  hence "(λp. 1 / (s - 1) * (ln (real p) / (p powr s - 1) + B p s))
            abs_summable_on {p. prime p}"
    by (subst (asm) abs_summable_on_cong[OF eq refl]) auto
  hence summable:
    "(λp. ln (real p) / (p powr s - 1) + B p s) abs_summable_on {p. prime p}"
    using s by (subst (asm) abs_summable_on_cmult_right_iff) auto

  from s have [simp]: "s  1" by auto
  have "eval_fds fds_newman s =
          (ap | prime p. (ln (real p) / real p) * hurwitz_zeta p s)"
    using s by (rule eval_fds_newman_conv_infsetsum)
  also have " = (ap | prime p. 1 / (s - 1) * (ln (real p) / (p powr s - 1) + B p s))"
    by (intro infsetsum_cong eq) auto
  also have " = 1 / (s - 1) * (ap | prime p. ln (real p) / (p powr s - 1) + B p s)"
    (is "_ = _ * ?S") by (rule infsetsum_cmult_right[OF summable])
  also have "?S = (p. if prime p then 
                      ln (real p) / (p powr s - 1) + B p s else 0)"
    by (subst infsetsum_nat[OF summable]) auto
  also have " = (p. (if prime p then ln (real p) / (p powr s - 1) else 0) + 
                        (if prime p then B p s else 0))"
    by (intro suminf_cong) auto
  also have " = pre_newman s - deriv zeta s / zeta s"
    using sums_pre_newman[of s] sums_logderiv_zeta[of s] s
    by (subst suminf_add [symmetric]) (auto simp: sums_iff)
  finally show ?thesis by (simp add: newman_def)
qed

end

text ‹
  Next, we shall attempt to get rid of the pole by subtracting suitable multiples of $\zeta(s)$
  and $\zeta'(s)$. To this end, we shall first prove the following alternative definition of 
  $\zeta'(s)$:
›
lemma deriv_zeta_eq':
  assumes "0 < Re s" "s  1"
  shows "deriv zeta s = deriv (λz. pre_zeta 1 z * (z - 1)) s / (s - 1) -
                          (pre_zeta 1 s * (s - 1) + 1) / (s - 1)2"
    (is "_ = ?rhs")
proof (rule DERIV_imp_deriv)
  have [derivative_intros]: "(pre_zeta 1 has_field_derivative deriv (pre_zeta 1) s) (at s)"
    by (intro holomorphic_derivI[of _ UNIV] holomorphic_intros) auto
  have *: "deriv (λz. pre_zeta 1 z * (z - 1)) s = deriv (pre_zeta 1) s * (s - 1) + pre_zeta 1 s"
    by (subst deriv_mult)
       (auto intro!: holomorphic_on_imp_differentiable_at[of _ UNIV] holomorphic_intros)
  hence "((λs. pre_zeta 1 s + 1 / (s - 1)) has_field_derivative
           deriv (pre_zeta 1) s - 1 / ((s - 1) * (s - 1))) (at s)"
    using assms by (auto intro!: derivative_eq_intros)
  also have "deriv (pre_zeta 1) s - 1 / ((s - 1) * (s - 1)) = ?rhs"
    using * assms by (simp add: divide_simps power2_eq_square, simp add: field_simps)
  also have "((λs. pre_zeta 1 s + 1 / (s - 1)) has_field_derivative ?rhs) (at s) 
               (zeta has_field_derivative ?rhs) (at s)"
    using assms
    by (intro has_field_derivative_cong_ev eventually_mono[OF t1_space_nhds[of _ 1]])
       (auto simp: zeta_def hurwitz_zeta_def)
  finally show  .
qed

text ‹
  From this, it follows that $(s - 1) \zeta'(s) - \zeta'(s) / \zeta(s)$ is analytic 
  for $\mathfrak{R}(s) \geq 1$:
›
lemma analytic_zeta_derivdiff:
  obtains a where
    "(λz. if z = 1 then a else (z - 1) * deriv zeta z - deriv zeta z / zeta z)
          analytic_on {s. Re s  1}" 
proof 
  have neq: "pre_zeta 1 z * (z - 1) + 1  0" if "Re z  1" for z
    using zeta_Re_ge_1_nonzero[of z] that
    by (cases "z = 1") (auto simp: zeta_def hurwitz_zeta_def divide_simps)
  let ?g = "λz. (1 - inverse (pre_zeta 1 z * (z - 1) + 1)) * ((z - 1) *
                deriv ((λu. pre_zeta 1 u * (u - 1))) z - (pre_zeta 1 z * (z - 1) + 1))"
  show "(λz. if z = 1 then deriv ?g 1 else (z - 1) * deriv zeta z - deriv zeta z / zeta z)
          analytic_on {s. Re s  1}" (is "?f analytic_on _")
  proof (rule pole_theorem_analytic_0)
    show "?g analytic_on {s. 1  Re s}" using neq
      by (auto intro!: analytic_intros)
  next
    show "d>0. wball z d - {1}. ?g w = (w - 1) * ?f w"
      if z: "z  {s. 1  Re s}" for z
    proof -
      have *: "isCont (λz. pre_zeta 1 z * (z - 1) + 1) z"
        by (auto intro!: continuous_intros)
       obtain e where "e > 0" and e: "y. dist z y < e  pre_zeta (Suc 0) y * (y-1) + 1  0"
         using continuous_at_avoid [OF * neq[of z]] z by auto
      show ?thesis
      proof (intro exI ballI conjI)
        fix w
        assume w: "w  ball z (min e 1) - {1}"
        then have "Re w > 0"
          using complex_Re_le_cmod [of "z-w"] z by (simp add: dist_norm)
        with w show "?g w = (w - 1) * (if w = 1 then deriv ?g 1 else
                        (w - 1) * deriv zeta w - deriv zeta w / zeta w)"
          by (subst (1 2) deriv_zeta_eq', 
              simp_all add: zeta_def hurwitz_zeta_def divide_simps e power2_eq_square)
             (simp_all add: algebra_simps)?
      qed (use e > 0 in auto)
    qed
  qed auto
qed

text ‹
  Finally, $f(s) + \zeta'(s) + c\zeta(s)$ is analytic.
›
lemma analytic_newman_variant:
  obtains c a where
     "(λz. if z = 1 then a else newman z + deriv zeta z + c * zeta z) analytic_on {s. Re s  1}"
proof -
  obtain c where (* -euler_mascheroni *)
    c: "(λz. if z = 1 then c else (z - 1) * deriv zeta z - deriv zeta z / zeta z)
        analytic_on {s. Re s  1}"
    using analytic_zeta_derivdiff by blast
  let ?g = "λz. pre_newman z +
          (if z = 1 then c
           else (z - 1) * deriv zeta z -
                deriv zeta z / zeta z) - (c + pre_newman 1) * (pre_zeta 1 z * (z - 1) + 1)"
  have "(λz. if z = 1 then deriv ?g 1 else newman z + deriv zeta z + (-(c + pre_newman 1)) * zeta z)
        analytic_on {s. Re s  1}"  (is "?f analytic_on _")
  proof (rule pole_theorem_analytic_0)
    show "?g analytic_on {s. 1  Re s}"
      by (intro c analytic_intros) auto
  next
    show "d>0. wball z d - {1}. ?g w = (w - 1) * ?f w"
      if "z  {s. 1  Re s}" for z using that
      by (intro exI[of _ 1], simp_all add: newman_def divide_simps zeta_def hurwitz_zeta_def)
         (auto simp: field_simps)?
  qed auto
  with that show ?thesis by blast
qed


subsection ‹The asymptotic expansion of 𝔐›

text ‹
  Our next goal is to show the key result that $\mathfrak{M}(x) = \ln n + c + o(1)$.

  As a first step, we invoke Ingham's Tauberian theorem on the function we have
  just defined and obtain that the sum
  \[\sum\limits_{n=1}^\infty \frac{\mathfrak{M}(n) - \ln n + c}{n}\]
  exists.
›
lemma mertens_summable:
  obtains c :: real where "summable (λn. (𝔐 n - ln n + c) / n)"
proof -
  (* c = euler_mascheroni - pre_newman 1 *)
  from analytic_newman_variant obtain c a where
    analytic: "(λz. if z = 1 then a else newman z + deriv zeta z + c * zeta z)
                 analytic_on {s. Re s  1}" .
  define f where "f = (λz. if z = 1 then a else newman z + deriv zeta z + c * zeta z)"
  have analytic: "f analytic_on {s. Re s  1}" using analytic by (simp add: f_def)
  define F where "F = fds_newman + fds_deriv fds_zeta + fds_const c * fds_zeta"

  note le = conv_abscissa_add_leI conv_abscissa_deriv_le conv_abscissa_newman conv_abscissa_mult_const_left
  note intros = le le[THEN le_less_trans] le[THEN order.trans] fds_converges
  have eval_F: "eval_fds F s = f s" if s: "Re s > 1" for s
  proof -
    have "eval_fds F s = eval_fds (fds_newman + fds_deriv fds_zeta) s +
                           eval_fds (fds_const c * fds_zeta) s"
      unfolding F_def using s by (subst eval_fds_add) (auto intro!: intros)
    also have " = f s" using s unfolding f_def
      by (subst eval_fds_add)
         (auto intro!: intros simp: eval_fds_newman eval_fds_deriv_zeta eval_fds_mult eval_fds_zeta)
    finally show ?thesis .
  qed

  have conv: "fds_converges F s" if "Re s  1" for s
  proof (rule Newman_Ingham_1)
    have "(λn. 𝔐 (real n) - ln (real n))  O(λ_. 1)"
      using mertens_bounded by (rule landau_o.big.compose) real_asymp
    from natfun_bigo_1E[OF this, of 1]
      obtain c' where c': "c'  1" "n. ¦𝔐 (real n) - ln (real n)¦  c'" by auto
    have "Bseq (fds_nth F)"
    proof (intro BseqI allI)
      fix n :: nat
      show "norm (fds_nth F n)  (c' + norm c)" unfolding F_def using c'
        by (auto simp: fds_nth_zeta fds_nth_deriv fds_nth_newman scaleR_conv_of_real in_Reals_norm
                 intro!: order.trans[OF norm_triangle_ineq] add_mono)
    qed (insert c', auto intro: add_pos_nonneg)
    thus "fds_nth F  O(λ_. 1)" by (simp add: natfun_bigo_iff_Bseq)
  next
    show "f analytic_on {s. Re s  1}" by fact
  next
    show "eval_fds F s = f s" if "Re s > 1" for s using that by (rule eval_F)
  qed (insert that, auto simp: F_def intro!: intros)
  from conv[of 1] have "summable (λn. fds_nth F n / of_nat n)"
    unfolding fds_converges_def by auto
  also have "?this  summable (λn. (𝔐 n - Ln n + c) / n)"
    by (intro summable_cong eventually_mono[OF eventually_gt_at_top[of 0]])
       (auto simp: F_def fds_nth_newman fds_nth_deriv fds_nth_zeta scaleR_conv_of_real
             intro!: sum.cong dest: prime_gt_0_nat)
  finally have "summable (λn. (𝔐 n - Re (Ln (of_nat n)) + Re c) / n)"
    by (auto dest: summable_Re)
  also have "?this  summable (λn. (𝔐 n - ln n + Re c) / n)"
    by (intro summable_cong eventually_mono[OF eventually_gt_at_top[of 0]]) (auto intro!: sum.cong)
  finally show ?thesis using that[of "Re c"] by blast
qed

text ‹
  Next, we prove a lemma given by Newman stating that if the sum $\sum a_n / n$ exists and
  $a_n + \ln n$ is nondecreasing, then $a_n$ must tend to 0. Unfortunately, the proof is
  rather tedious, but so is the paper version by Newman.
›
lemma sum_goestozero_lemma:
  fixes d::real
  assumes d: "¦i = M..N. a i / i¦ < d" and le: "n. a n + ln n  a (Suc n) + ln (Suc n)"
      and "0 < M" "M < N"
    shows "a M  d * N / (real N - real M) + (real N - real M) / M 
          -a N  d * N / (real N - real M) + (real N - real M) / M"
proof -
  have "0  d"
    using assms by linarith+
  then have "0  d * N / (N - M + 1)" by simp
  then have le_dN: "0  x  x  d * N / (N - M + 1)  x  d * N / (N - M + 1)" for x::real
    by linarith
  have le_a_ln: "a m + ln m  a n + ln n" if "n  m" for n m
    by (rule transitive_stepwise_le) (use le that in auto)
  have *: "x  b  y  b" if "a  b" "x  a" "y  a" for a b x y::real
    using that by linarith
  show ?thesis
  proof (rule *)
    show "d * N / (N - M) + ln (N / M)  d * N / (real N - real M) + (real N - real M) / M"
      using 0 < M M < N ln_le_minus_one [of "N / M"]
      by (simp add: of_nat_diff) (simp add: divide_simps)
  next
    have "a M - ln (N / M)  (d * N) / (N - M + 1)"
    proof (rule le_dN)
      assume 0: "0  a M - ln (N / M)"
      have "(Suc N - M) * (a M - ln (N / M)) / N = (i = M..N. (a M - ln (N / M)) / N)"
        by simp
      also have "  (i = M..N. a i / i)"
      proof (rule sum_mono)
        fix i
        assume i: "i  {M..N}"
        with 0 < M have "0 < i" by auto
        have "(a M - ln (N / M)) / N  (a M - ln (N / M)) / i"
          using 0 using i 0 < M by (simp add: frac_le_eq divide_simps mult_left_mono)
        also have "a M + ln (real M)  a i + ln (real N)"
          by (rule order.trans[OF le_a_ln[of M i]]) (use i assms in auto)
        hence "(a M - ln (N / M)) / i  a i / real i"
          using assms i by (intro divide_right_mono) (auto simp: ln_div field_simps)
        finally show "(a M - ln (N / M)) / real N  a i / real i" .
      qed
      finally have "((Suc N) - M) * (a M - ln (N / M)) / N  ¦i = M..N. a i / i¦"
        by simp
      also have "  d" using d by simp
      finally have "((Suc N) - M) * (a M - ln (N / M)) / N  d" .
      then show ?thesis
        using M < N  by (simp add: of_nat_diff field_simps)
    qed
    also have "  d * N / (N - M)"
      using assms(1,4) by (simp add: field_simps)
    finally show "a M  d * N / (N - M) + ln (N / M)" by simp
  next
    have "- a N - ln (N / M)  (d * N) / (N - M + 1)"
    proof (rule le_dN)
      assume 0: "0  - a N - ln (N / M)"
      have "(i = M..N. a i / i)  (i = M..N. (a N + ln (N / M)) / N)"
      proof (rule sum_mono)
        fix i
        assume i: "i  {M..N}"
        with 0 < M have "0 < i" by auto
        have "a i + ln (real M)  a N + ln (real N)"
          by (rule order.trans[OF _ le_a_ln[of i N]]) (use i assms in auto)
        hence "a i / i  (a N + ln (N / M)) / i"
          using assms(3,4) by (intro divide_right_mono) (auto simp: field_simps ln_div)
        also have "  (a N + ln (N / M)) / N"
          using i i > 0 0 by (intro divide_left_mono_neg) auto
        finally show "a i / i  (a N + ln (N / M)) / N" .
      qed
      also have " = ((Suc N) - M) * (a N + ln (N / M)) / N"
        by simp
      finally have "(i = M..N. a i / i)  (real (Suc N) - real M) * (a N + ln (N / M)) / N"
        using M < N by (simp add: of_nat_diff)
      then have "-((real (Suc N) - real M) * (a N + ln (N / M)) / N)  ¦i = M..N. a i / i¦"
        by linarith
      also have "  d" using d by simp
      finally have "- ((real (Suc N) - real M) * (a N + ln (N / M)) / N)  d" .
      then show ?thesis
        using M < N  by (simp add: of_nat_diff field_simps)
    qed
    also have "  d * N / real (N - M)"
      using 0 < M M < N 0  d by (simp add: field_simps)
    finally show "-a N  d * N / real (N - M) + ln (N / M)" by simp
  qed
qed

proposition sum_goestozero_theorem:
  assumes summ: "summable (λi. a i / i)"
      and le:   "n. a n + ln n  a (Suc n) + ln (Suc n)"
    shows "a  0"
proof (clarsimp simp: lim_sequentially)
  fix r::real
  assume "r > 0"
  have *: "n0. nn0. ¦a n¦ < ε" if ε: "0 < ε" "ε < 1" for ε
  proof -
    have "0 < (ε / 8)2" using 0 < ε  by simp
    then obtain N0 where N0: "m n. m  N0  norm (k=m..n. (λi. a i / i) k) < (ε / 8)2"
      by (metis summable_partial_sum_bound summ)
    obtain N1 where "real N1 > 4 / ε"
      using reals_Archimedean2[of "4 / ε"] ε by auto
    hence "N1  0" and N1: "1 / real N1 < ε / 4" using ε
      by (auto simp: divide_simps mult_ac intro: Nat.gr0I)

    have "¦a n¦ < ε" if n: "n  2 * N0 + N1 + 7" for n
    proof -
      define k where "k = n * ε/4"
      have "n * ε / 4 > 1" and "n * ε / 4  n / 4" and "n / 4 < n"
        using less_le_trans[OF N1, of "n / N1 * ε / 4"] N1  0 ε n by (auto simp: field_simps)
      hence k: "k > 0" "4 * k  n" "nat k < n" "(n * ε / 4) - 1 < k" "k  (n * ε / 4)"
        unfolding k_def by linarith+

      have "-a n < ε"
      proof -
        have "N0  n - nat k"
          using n k by linarith
        then have *: "¦k = n - nat k .. n. a k / k¦ < (ε / 8)2"
          using N0 [of "n - nat k" n] by simp
        have "-a n  (ε / 8)2 * n / n * ε / 4 + n * ε / 4 / (n - k)"
          using sum_goestozero_lemma [OF * le, THEN conjunct2] k by (simp add: of_nat_diff k_def)
        also have "< ε"
        proof -
          have "ε / 16 * n / k < 2"
            using k by (auto simp: field_simps)
          then have "ε * (ε / 16 * n / k) < ε * 2"
            using ε mult_less_cancel_left_pos by blast
          then have "(ε / 8)2 * n / k < ε / 2"
            by (simp add: field_simps power2_eq_square)
          moreover have "k / (n - k) < ε / 2"
          proof -
            have "(ε + 2) * k < 4 * k" using k ε by simp
            also have "  ε * real n" using k by (auto simp: field_simps)
            finally show ?thesis using k by (auto simp: field_simps)
          qed
          ultimately show ?thesis unfolding k_def by linarith
        qed
        finally show ?thesis .
      qed
      moreover have "a n < ε"
      proof -
        have "N0  n" using n k by linarith
        then have *: "¦k = n .. n + nat k. a k / k¦ < (ε/8)2"
          using N0 [of n "n + nat k"] by simp
        have "a n  (ε/8)2 * (n + nat k) / k + k / n"
          using sum_goestozero_lemma [OF * le, THEN conjunct1] k by (simp add: of_nat_diff)
        also have "< ε"
        proof -
          have "4  28 * real_of_int k" using k by linarith
          then have "ε/16 * n / k < 2" using k by (auto simp: field_simps)
          have "ε * (real n + k) < 32 * k"
          proof -
            have "ε * n / 4 < k + 1" by (simp add: mult.commute k_def)
            then have "ε * n < 4 * k + 4" by (simp add: divide_simps)
            also have "  8 * k" using k by auto
            finally have 1: "ε * real n < 8 * k" .
            have 2: "ε * k < k" using k ε by simp
            show ?thesis using k add_strict_mono [OF 1 2] by (simp add: algebra_simps)
        qed
          then have "(ε / 8)2 * real (n + nat k) / k < ε / 2"
            using ε k by (simp add: divide_simps mult_less_0_iff power2_eq_square)
          moreover have "k / n < ε / 2"
            using k ε by (auto simp: k_def field_simps)
          ultimately show ?thesis by linarith
        qed
        finally show ?thesis .
      qed
      ultimately show ?thesis by force
    qed
    then show ?thesis by blast
  qed
  show "n0. nn0. ¦a n¦ < r"
    using * [of "min r (1/5)"] 0 < r by force
qed


text ‹
  This leads us to the main intermediate result:
›
lemma Mertens_convergent: "convergent (λn::nat. 𝔐 n - ln n)"
proof -
  obtain c where c: "summable (λn. (𝔐 n - ln n + c) / n)"
    by (blast intro: mertens_summable)
  then obtain l where l: "(λn. (𝔐 n - ln n + c) / n) sums l"
    by (auto simp: summable_def)
  have *: "(λn. 𝔐 n - ln n + c)  0"
    by (rule sum_goestozero_theorem[OF c]) auto
  hence "(λn. 𝔐 n - ln n)  -c"
    by (simp add: tendsto_iff dist_norm)
  thus ?thesis by (rule convergentI)
qed

corollary 𝔐_minus_ln_limit:
  obtains c where "((λx::real. 𝔐 x - ln x)  c) at_top"
proof -
  from Mertens_convergent obtain c where "(λn. 𝔐 n - ln n)  c"
    by (auto simp: convergent_def)
  hence 1: "((λx::real. 𝔐 (nat x) - ln (nat x))  c) at_top" 
    by (rule filterlim_compose) real_asymp
  have 2: "((λx::real. ln (nat x) - ln x)  0) at_top"
    by real_asymp
  have 3: "((λx. 𝔐 x - ln x)  c) at_top"
    using tendsto_add[OF 1 2] by simp
  with that show ?thesis by blast
qed


subsection ‹The asymptotics of the prime-counting functions›

text ‹
  We will now use the above result to prove the asymptotics of the prime-counting functions
  $\vartheta(x) \sim x$, $\psi(x) \sim x$, and $\pi(x) \sim x / \ln x$. The last of these is 
  typically called the Prime Number Theorem, but since these functions can be expressed in terms 
  of one another quite easily, knowing the asymptotics of any of them immediately gives the 
  asymptotics of the other ones.

  In this sense, all of the above are equivalent formulations of the Prime Number Theorem.
  The one we shall tackle first, due to its strong connection to the $\mathfrak{M}$ function, is
  $\vartheta(x) \sim x$.

  We know that $\mathfrak{M}(x)$ has the asymptotic expansion
  $\mathfrak{M}(x) = \ln x + c + o(1)$. We also know that
  \[\vartheta(x) = x\mathfrak{M}(x) - \int\nolimits_2^x \mathfrak{M}(t) \,\mathrm{d}t\ .\]
  Substituting in the above asymptotic equation, we obtain:
  \begin{align*}
  \vartheta(x) &= x\ln x + cx + o(x) - \int\nolimits_2^x \ln t + c + o(1) \,\mathrm{d}t\\
            &= x\ln x + cx + o(x) - (x\ln x - x + cx + o(x))\\
            &= x + o(x)
  \end{align*}
  In conclusion, $\vartheta(x) \sim x$.
›
theorem θ_asymptotics: "θ ∼[at_top] (λx. x)"
proof -
  from 𝔐_minus_ln_limit obtain c where c: "((λx. 𝔐 x - ln x)  c) at_top"
    by auto
  define r where "r = (λx. 𝔐 x - ln x - c)"
  have 𝔐_expand: "𝔐 = (λx. ln x + c + r x)"
    by (simp add: r_def)
  have r: "r  o(λ_. 1)" unfolding r_def
    using tendsto_add[OF c tendsto_const[of "-c"]] by (intro smalloI_tendsto) auto

  define r' where "r' = (λx. integral {2..x} r)"
  have integrable_r: "r integrable_on {x..y}"
    if "2  x" for x y :: real using that unfolding r_def
    by (intro integrable_diff integrable_primes_M)
       (auto intro!: integrable_continuous_real continuous_intros)
  hence integral: "(r has_integral r' x) {2..x}" if "x  2" for x
    by (auto simp: has_integral_iff r'_def)
  have r': "r'  o(λx. x)" using integrable_r unfolding r'_def
    by (intro integral_smallo[OF r]) (auto simp: filterlim_ident)

  define C where "C = 2 * (c + ln 2 - 1)"
  have "θ ∼[at_top] (λx. x + (r x * x + C - r' x))"
  proof (intro asymp_equiv_refl_ev eventually_mono[OF eventually_gt_at_top])
    fix x :: real assume x: "x > 2"
    have "(𝔐 has_integral ((x * ln x - x + c * x) - (2 * ln 2 - 2 + c * 2) + r' x)) {2..x}"
      unfolding 𝔐_expand using x
      by (intro has_integral_add[OF fundamental_theorem_of_calculus integral])
         (auto simp flip: has_real_derivative_iff_has_vector_derivative
               intro!: derivative_eq_intros continuous_intros)
    from has_integral_unique[OF θ_conv_𝔐_integral this]
      show "θ x = x + (r x * x + C - r' x)" using x
      by (simp add: field_simps 𝔐_expand C_def)
  qed
  also have "(λx. r x * x + C - r' x)  o(λx. x)"
  proof (intro sum_in_smallo r)
    show "(λ_. C)  o(λx. x)" by real_asymp
  qed (insert landau_o.small_big_mult[OF r, of "λx. x"] r', simp_all)
  hence "(λx. x + (r x * x + C - r' x)) ∼[at_top] (λx. x)"
    by (subst asymp_equiv_add_right) auto
  finally show ?thesis by auto
qed

text ‹
  The various other forms of the Prime Number Theorem follow as simple corollaries.
›
corollary ψ_asymptotics: "ψ ∼[at_top] (λx. x)"
  using θ_asymptotics PNT4_imp_PNT5 by simp
  
corollary prime_number_theorem: "π ∼[at_top] (λx. x / ln x)"
  using θ_asymptotics PNT4_imp_PNT1 by simp

corollary ln_π_asymptotics: "(λx. ln (π x)) ∼[at_top] ln"
  using prime_number_theorem PNT1_imp_PNT1' by simp

corollary π_ln_π_asymptotics: "(λx. π x * ln (π x)) ∼[at_top] (λx. x)"
  using prime_number_theorem PNT1_imp_PNT2 by simp

corollary nth_prime_asymptotics: "(λn. real (nth_prime n)) ∼[at_top] (λn. real n * ln (real n))"
  using π_ln_π_asymptotics PNT2_imp_PNT3 by simp


text ‹
  The following versions use a little less notation.
›
corollary prime_number_theorem': "((λx. π x / (x / ln x))  1) at_top"
  using prime_number_theorem
  by (rule asymp_equivD_strong[OF _ eventually_mono[OF eventually_gt_at_top[of 1]]]) auto

corollary prime_number_theorem'':
  "(λx. card {p. prime p  real p  x}) ∼[at_top] (λx. x / ln x)"
proof -
  have "π = (λx. card {p. prime p  real p  x})"
    by (intro ext) (simp add: π_def prime_sum_upto_def)
  with prime_number_theorem show ?thesis by simp
qed

corollary prime_number_theorem''':
  "(λn. card {p. prime p  p  n}) ∼[at_top] (λn. real n / ln (real n))"
proof -
  have "(λn. card {p. prime p  real p  real n}) ∼[at_top] (λn. real n / ln (real n))"
    using prime_number_theorem''
    by (rule asymp_equiv_compose') (simp add: filterlim_real_sequentially)
  thus ?thesis by simp
qed

(*<*)
unbundle no_prime_counting_notation
(*>*)

end