Theory Euler_Products

theory Euler_Products
imports Analysis Multiplicative_Function
(*
    File:      Euler_Products.thy
    Author:    Manuel Eberl, TU M√ľnchen
*)
section ‹Euler product expansions›
theory Euler_Products
imports 
  "HOL-Analysis.Analysis"
  Multiplicative_Function
begin

lemma prime_factors_power_subset:
  "prime_factors (x ^ n) ⊆ prime_factors x"
  by (cases "n = 0") (auto simp: prime_factors_power)

lemma prime_power_product_in_Pi:
  "(λg. ∏p∈{p. p ≤ (n::nat) ∧ prime p}. p ^ g p)
    ∈ ({p. p ≤ n ∧ prime p} →E UNIV) →
       {m. 0 < m ∧ prime_factors m ⊆ {..n}}"
proof (safe, goal_cases)
  case (2 f p)
  have "prime_factors (∏p∈{p. p ≤ n ∧ prime p}. p ^ f p) = 
          (⋃p∈{p. p ≤ n ∧ prime p}. prime_factors (p ^ f p))"
    by (subst prime_factors_prod) auto
  also have "… ⊆ (⋃p∈{p. p ≤ n ∧ prime p}. prime_factors p)"
    using prime_factors_power_subset by blast
  also have "… ⊆ (⋃p∈{p. p ≤ n ∧ prime p}. {p})"
    by (auto simp: prime_factors_dvd prime_gt_0_nat dest!: dvd_imp_le)
  also have "… ⊆ {..n}" by auto
  finally show ?case using 2 by auto
qed (auto simp: prime_gt_0_nat)

lemma inj_prime_power: "inj_on (λx. fst x ^ snd x :: nat) ({a. prime a} × {0<..})"
proof (intro inj_onI, clarify, goal_cases)
  case (1 p m q n)
  with prime_power_eq_imp_eq[of p q m n] and 1 
    have "p = q" by auto
  moreover from this have "m = n"
    using prime_gt_1_nat 1 by auto
  ultimately show ?case by simp
qed

lemma bij_betw_prime_powers:
  "bij_betw (λg. ∏p∈{p. p ≤ n ∧ prime p}. p ^ g p) ({p. p ≤ n ∧ prime p} →E UNIV)
     {m. 0 < m ∧ prime_factors m ⊆ {..(n::nat)}}"
proof (rule bij_betwI[of _ _ _ "(λm p. if p ≤ n ∧ prime p then multiplicity p m else undefined)"], 
         goal_cases)
  case 1
  show ?case by (rule prime_power_product_in_Pi)
next
  case 2
  show ?case
    by (auto split: if_splits)
next
  case (3 f)
  show ?case
  proof (rule ext, goal_cases)
    case (1 q)
    show ?case
    proof (cases "q ≤ n ∧ prime q")
      case True
      hence "multiplicity q (∏p∈{p. p ≤ n ∧ prime p}. p ^ f p) = 
               (∑x∈{p. p ≤ n ∧ prime p}. multiplicity q (x ^ f x))"
        by (subst prime_elem_multiplicity_prod_distrib) auto
      also have "… = (∑x∈{p. p ≤ n ∧ prime p}. if x = q then f q else 0)"
        using True by (intro sum.cong refl) (auto simp: multiplicity_distinct_prime_power)
      also have "… = f q" using True by auto
      finally show ?thesis using True by simp
    qed (insert 3, force+)
  qed
next
  case (4 m)
  have "(∏p | p ≤ n ∧ prime p. p ^ (if p ≤ n ∧ prime p then multiplicity p m else undefined)) =
          (∏p∈prime_factors m. p ^ multiplicity p m)"
  proof (rule prod.mono_neutral_cong)
    show "finite (prime_factors m)" by simp
  qed (insert 4, auto simp: prime_factors_multiplicity)
  also from 4 have "… = m"
    by (intro prime_factorization_nat [symmetric]) auto
  finally show ?case .
qed
      
lemma
  fixes f :: "nat ⇒ 'a :: {real_normed_field,banach,second_countable_topology}"
  assumes summable: "summable (λn. norm (f n))"
  assumes "multiplicative_function f"
  shows   abs_convergent_euler_product:
            "abs_convergent_prod (λp. if prime p then ∑n. f (p ^ n) else 1)"
    and   euler_product_LIMSEQ:
            "(λn. (∏p≤n. if prime p then ∑n. f (p ^ n) else 1)) ⇢ (∑n. f n)"
proof -
  interpret f: multiplicative_function f by fact
  define N where "N = (∑n. norm (f n))"

  have summable': "f abs_summable_on A" for A
    by (rule abs_summable_on_subset[of _ UNIV])
       (insert summable, auto simp: abs_summable_on_nat_iff')

  have summable'': "(λx. f (p ^ x)) abs_summable_on A" if "prime p" for A p
  proof (subst abs_summable_on_reindex_iff[of _ _ f])
    from ‹prime p› have "p > 1" 
      by (rule prime_gt_1_nat)
    thus "inj_on (λi. p ^ i) A"
      by (auto simp: inj_on_def)
  qed (intro summable')

  have "(λn. norm ((∑m. f m) - (∏p∈{p. p ≤ n ∧ prime p}. ∑i. f (p ^ i)))) ⇢ 0"
          (is "filterlim ?h _ _")
  proof (rule tendsto_sandwich)
    show "eventually (λn. ?h n ≤ N - (∑m≤n. norm (f m))) at_top"
    proof (intro always_eventually allI)
      fix n :: nat
      interpret product_sigma_finite "λ_::nat. count_space (UNIV :: nat set)"
        by (intro product_sigma_finite.intro sigma_finite_measure_count_space)

      have "(∏p | p ≤ n ∧ prime p. ∑i. f (p ^ i)) =
              (∏p | p ≤ n ∧ prime p. ∑a i∈UNIV. f (p ^ i))"
        by (intro prod.cong refl infsetsum_nat' [symmetric] summable'') auto
      also have "… = (∑ag∈{p. p ≤ n ∧ prime p} →E UNIV.
                         ∏x∈{p. p ≤ n ∧ prime p}. f (x ^ g x))"
        by (subst infsetsum_prod_PiE [symmetric])
           (auto simp: prime_gt_Suc_0_nat summable'')
      also have "… = (∑ag∈{p. p ≤ n ∧ prime p} →E UNIV.
                         f (∏x∈{p. p ≤ n ∧ prime p}. x ^ g x))"
        by (subst f.prod_coprime) (auto simp add: primes_coprime)
      also have "… = (∑am | m > 0 ∧ prime_factors m ⊆ {..n}. f m)"
        by (intro infsetsum_reindex_bij_betw bij_betw_prime_powers)
      also have "(∑am∈UNIV. f m) - … = (∑am∈UNIV - {m. m > 0 ∧ prime_factors m ⊆ {..n}}. f m)"
        by (intro infsetsum_Diff [symmetric] summable') auto
      also have "(∑am∈UNIV. f m) = (∑m. f m)"
        by (intro infsetsum_nat' summable')
      also have "UNIV - {m. m > 0 ∧ prime_factors m ⊆ {..n}} = 
                   insert 0 {m. ¬prime_factors m ⊆ {..n}}"
        by auto
      also have "(∑am∈…. f m) = (∑am | ¬prime_factors m ⊆ {..n}. f m)"
        by (intro infsetsum_cong_neutral) auto
      also have "norm … ≤ (∑am | ¬prime_factors m ⊆ {..n}. norm (f m))"
        by (rule norm_infsetsum_bound)
      also have "… ≤ (∑am∈{n<..}. norm (f m))" 
      proof (intro infsetsum_mono_neutral_left summable' abs_summable_on_normI)
        show "{m. ¬ prime_factors m ⊆ {..n}} ⊆ {n<..}"
        proof safe
          fix m k assume "¬m > n" and "k ∈ prime_factors m"
          thus "k ≤ n" by (cases "m = 0") (auto simp: prime_factors_dvd dest: dvd_imp_le)
        qed
      qed auto
      also have "{n<..} = UNIV - {..n}" 
        by auto
      also have "(∑am∈…. norm (f m)) = (∑am∈UNIV. norm (f m)) - (∑am∈{..n}. norm (f m))"
        using summable by (intro infsetsum_Diff) (auto simp: abs_summable_on_nat_iff')
      also have "(∑am∈UNIV. norm (f m)) = N"
        unfolding N_def using summable 
        by (intro infsetsum_nat') (auto simp: abs_summable_on_nat_iff')
      also have "(∑am∈{..n}. norm (f m)) = (∑m≤n. norm (f m))"
        by (simp add: suminf_finite)
      finally show "?h n ≤ N - (∑m≤n. norm (f m))" .
   qed
  next
    show "eventually (λn. ?h n ≥ 0) at_top" by simp
  next
    show "(λn. N - (∑m≤n. norm (f m))) ⇢ 0" unfolding N_def
      by (rule tendsto_eq_intros refl summable_LIMSEQ' summable)+ simp_all
  qed simp_all
  hence "(λn. (∑m. f m) - (∏p∈{p. p ≤ n ∧ prime p}. ∑i. f (p ^ i))) ⇢ 0"
    by (simp add: tendsto_norm_zero_iff)
  from tendsto_diff[OF tendsto_const[of "∑m. f m"] this]
    have "(λn. ∏p | p ≤ n ∧ prime p. ∑i. f (p ^ i)) ⇢ (∑m. f m)" by simp
  also have "(λn. ∏p | p ≤ n ∧ prime p. ∑i. f (p ^ i)) = 
                 (λn. ∏p≤n. if prime p then (∑i. f (p ^ i)) else 1)"
    by (intro ext prod.mono_neutral_cong_left) auto
  finally show "… ⇢ (∑m. f m)" .

  show "abs_convergent_prod (λp. if prime p then (∑i. f (p ^ i)) else 1)"
  proof (rule summable_imp_abs_convergent_prod)
    have "(λ(p,i). f (p ^ i)) abs_summable_on {p. prime p} × {0<..}"
      unfolding case_prod_unfold 
      by (subst abs_summable_on_reindex_iff[OF inj_prime_power]) fact
    hence "(λp. ∑ai∈{0<..}. f (p ^ i)) abs_summable_on {p. prime p}"
      by (rule abs_summable_on_Sigma_project1') simp_all
    also have "?this ⟷ (λp. (∑i. f (p ^ i)) - 1) abs_summable_on {p. prime p}"
    proof (intro abs_summable_on_cong refl)
      fix p :: nat assume p: "p ∈ {p. prime p}"
      have "{0<..} = UNIV - {0::nat}" by auto
      also have "(∑ai∈…. f (p ^ i)) = (∑i. f (p ^ i)) - 1"
        using p by (subst infsetsum_Diff) (simp_all add: infsetsum_nat' summable'')
      finally show "(∑ai∈{0<..}. f (p ^ i)) = (∑i. f (p ^ i)) - 1" .
    qed
    finally have "summable (λp. if prime p then norm ((∑i. f (p ^ i)) - 1) else 0)"
      (is "summable ?T") by (simp add: abs_summable_on_nat_iff)
    also have "?T = (λp. norm ((if prime p then ∑i. f (p ^ i) else 1) - 1))"
      by (rule ext) (simp add: if_splits)
    finally show "summable …" .
  qed
qed

lemma
  fixes f :: "nat ⇒ 'a :: {real_normed_field,banach,second_countable_topology}"
  assumes summable: "summable (λn. norm (f n))"
  assumes "completely_multiplicative_function f"
  shows   abs_convergent_euler_product':
            "abs_convergent_prod (λp. if prime p then inverse (1 - f p) else 1)"
    and   completely_multiplicative_summable_norm: 
            "⋀p. prime p ⟹ norm (f p) < 1"
    and   euler_product_LIMSEQ':
            "(λn. (∏p≤n. if prime p then inverse (1 - f p) else 1)) ⇢ (∑n. f n)"
proof -
  interpret f: completely_multiplicative_function f by fact
  {
    fix p :: nat assume "prime p"
    hence "inj (λi. p ^ i)" 
      by (auto simp: inj_on_def dest: prime_gt_1_nat)
    from summable_reindex[OF summable this] 
      have *: "summable (λi. norm (f (p ^ i)))" by (auto simp: o_def)
    also have "(λi. norm (f (p ^ i))) = (λi. norm (f p) ^ i)"
      by (simp add: f.power norm_power)
    finally show "norm (f p) < 1"
      by (subst (asm) summable_geometric_iff) simp_all
    note * and this
  } note summable' = this

  have eq: "(λp. if prime p then (∑i. f (p ^ i)) else 1) = 
              (λp. if prime p then inverse (1 - f p) else 1)"
  proof (rule ext, goal_cases)
    case (1 p)
    show ?case
    proof (cases "prime p")
      case True
      hence "norm (f p) < 1" by (rule summable')
      from suminf_geometric[OF this] and True show ?thesis
        by (simp add: field_simps f.power)
    qed simp_all
  qed
  hence eq': "(λn. ∏p≤n. if prime p then ∑n. f (p ^ n) else 1) =
                (λn. ∏p≤n. if prime p then inverse (1 - f p) else 1)"
    by (auto simp: fun_eq_iff)

  have f: "multiplicative_function f" ..
  from abs_convergent_euler_product[OF assms(1) f] and euler_product_LIMSEQ[OF assms(1) f]
    show "abs_convergent_prod (λp. if prime p then inverse (1 - f p) else 1)"
       and "(λn. ∏p≤n. if prime p then inverse (1 - f p) else 1) ⇢ (∑n. f n)"
      by (simp_all only: eq eq')
qed

end