Theory Multiseries_Expansion

```(*
File:    Multiseries_Expansion.thy
Author:  Manuel Eberl, TU München

Asymptotic bases and Multiseries expansions of real-valued functions.
Contains automation to prove limits and asymptotic relationships between such functions.
*)
section ‹Multiseries expansions›
theory Multiseries_Expansion
imports
"HOL-Library.BNF_Corec"
"HOL-Library.Landau_Symbols"
Lazy_Eval
Inst_Existentials
Eventuallize
begin

(* TODO Move *)
lemma real_powr_at_top:
assumes "(p::real) > 0"
shows   "filterlim (λx. x powr p) at_top at_top"
proof (subst filterlim_cong[OF refl refl])
show "LIM x at_top. exp (p * ln x) :> at_top"
by (rule filterlim_compose[OF exp_at_top filterlim_tendsto_pos_mult_at_top[OF tendsto_const]])
show "eventually (λx. x powr p = exp (p * ln x)) at_top"
using eventually_gt_at_top[of 0] by eventually_elim (simp add: powr_def)
qed

subsection ‹Streams and lazy lists›

codatatype 'a msstream = MSSCons 'a "'a msstream"

primrec mssnth :: "'a msstream ⇒ nat ⇒ 'a" where
"mssnth xs 0 = (case xs of MSSCons x xs ⇒ x)"
| "mssnth xs (Suc n) = (case xs of MSSCons x xs ⇒ mssnth xs n)"

codatatype 'a msllist = MSLNil | MSLCons 'a "'a msllist"
for map: mslmap

fun lcoeff where
"lcoeff MSLNil n = 0"
| "lcoeff (MSLCons x xs) 0 = x"
| "lcoeff (MSLCons x xs) (Suc n) = lcoeff xs n"

primcorec msllist_of_msstream :: "'a msstream ⇒ 'a msllist" where
"msllist_of_msstream xs = (case xs of MSSCons x xs ⇒ MSLCons x (msllist_of_msstream xs))"

lemma msllist_of_msstream_MSSCons [simp]:
"msllist_of_msstream (MSSCons x xs) = MSLCons x (msllist_of_msstream xs)"
by (subst msllist_of_msstream.code) simp

lemma lcoeff_llist_of_stream [simp]: "lcoeff (msllist_of_msstream xs) n = mssnth xs n"
by (induction "msllist_of_msstream xs" n arbitrary: xs rule: lcoeff.induct;
subst msllist_of_msstream.code) (auto split: msstream.splits)

subsection ‹Convergent power series›

subsubsection ‹Definition›

definition convergent_powser :: "real msllist ⇒ bool" where
"convergent_powser cs ⟷ (∃R>0. ∀x. abs x < R ⟶ summable (λn. lcoeff cs n * x ^ n))"

lemma convergent_powser_stl:
assumes "convergent_powser (MSLCons c cs)"
shows   "convergent_powser cs"
proof -
from assms obtain R where "R > 0" "⋀x. abs x < R ⟹ summable (λn. lcoeff (MSLCons c cs) n * x ^ n)"
unfolding convergent_powser_def by blast
hence "summable (λn. lcoeff cs n * x ^ n)" if "abs x < R" for x
using that by (subst (asm) summable_powser_split_head [symmetric]) simp
thus ?thesis using ‹R > 0› by (auto simp: convergent_powser_def)
qed

definition powser :: "real msllist ⇒ real ⇒ real" where
"powser cs x = suminf (λn. lcoeff cs n * x ^ n)"

lemma isCont_powser:
assumes "convergent_powser cs"
shows   "isCont (powser cs) 0"
proof -
from assms obtain R where R: "R > 0" "⋀x. abs x < R ⟹ summable (λn. lcoeff cs n * x ^ n)"
unfolding convergent_powser_def by blast
hence *: "summable (λn. lcoeff cs n * (R/2) ^ n)" by (intro R) simp_all
from ‹R > 0› show ?thesis unfolding powser_def
by (intro isCont_powser[OF *]) simp_all
qed

lemma powser_MSLNil [simp]: "powser MSLNil = (λ_. 0)"

lemma powser_MSLCons:
assumes "convergent_powser (MSLCons c cs)"
shows   "eventually (λx. powser (MSLCons c cs) x = x * powser cs x + c) (nhds 0)"
proof -
from assms obtain R where R: "R > 0" "⋀x. abs x < R ⟹ summable (λn. lcoeff (MSLCons c cs) n * x ^ n)"
unfolding convergent_powser_def by blast
from R have "powser (MSLCons c cs) x = x * powser cs x + c" if "abs x < R" for x
using that unfolding powser_def by (subst powser_split_head) simp_all
moreover have "eventually (λx. abs x < R) (nhds 0)"
using ‹R > 0› filterlim_ident[of "nhds (0::real)"]
by (subst (asm) tendsto_iff) (simp add: dist_real_def)
ultimately show ?thesis by (auto elim: eventually_mono)
qed

definition convergent_powser' :: "real msllist ⇒ (real ⇒ real) ⇒ bool" where
"convergent_powser' cs f ⟷ (∃R>0. ∀x. abs x < R ⟶ (λn. lcoeff cs n * x ^ n) sums f x)"

lemma convergent_powser'_imp_convergent_powser:
"convergent_powser' cs f ⟹ convergent_powser cs"
unfolding convergent_powser_def convergent_powser'_def by (auto simp: sums_iff)

lemma convergent_powser'_eventually:
assumes "convergent_powser' cs f"
shows   "eventually (λx. powser cs x = f x) (nhds 0)"
proof -
from assms obtain R where "R > 0" "⋀x. abs x < R ⟹ (λn. lcoeff cs n * x ^ n) sums f x"
unfolding convergent_powser'_def by blast
hence "powser cs x = f x" if "abs x < R" for x
using that by (simp add: powser_def sums_iff)
moreover have "eventually (λx. abs x < R) (nhds 0)"
using ‹R > 0› filterlim_ident[of "nhds (0::real)"]
by (subst (asm) tendsto_iff) (simp add: dist_real_def)
ultimately show ?thesis by (auto elim: eventually_mono)
qed

subsubsection ‹Geometric series›

primcorec mssalternate :: "'a ⇒ 'a ⇒ 'a msstream" where
"mssalternate a b = MSSCons a (mssalternate b a)"

lemma case_mssalternate [simp]:
"(case mssalternate a b of MSSCons c d ⇒ f c d) = f a (mssalternate b a)"
by (subst mssalternate.code) simp

lemma mssnth_mssalternate: "mssnth (mssalternate a b) n = (if even n then a else b)"
by (induction n arbitrary: a b; subst mssalternate.code) simp_all

lemma convergent_powser'_geometric:
"convergent_powser' (msllist_of_msstream (mssalternate 1 (-1))) (λx. inverse (1 + x))"
proof -
have "(λn. lcoeff (msllist_of_msstream (mssalternate 1 (-1))) n * x ^ n) sums (inverse (1 + x))"
if "abs x < 1" for x :: real
proof -
have "(λn. (-1) ^ n * x ^ n) sums (inverse (1 + x))"
using that geometric_sums[of "-x"] by (simp add: field_simps power_mult_distrib [symmetric])
also have "(λn. (-1) ^ n * x ^ n) =
(λn. lcoeff (msllist_of_msstream (mssalternate 1 (-1))) n * x ^ n)"
by (auto simp add: mssnth_mssalternate fun_eq_iff)
finally show ?thesis .
qed
thus ?thesis unfolding convergent_powser'_def by (force intro!: exI[of _ 1])
qed

subsubsection ‹Exponential series›

primcorec exp_series_stream_aux :: "real ⇒ real ⇒ real msstream" where
"exp_series_stream_aux acc n = MSSCons acc (exp_series_stream_aux (acc / n) (n + 1))"

lemma mssnth_exp_series_stream_aux:
"mssnth (exp_series_stream_aux (1 / fact n) (n + 1)) m = 1 / fact (n + m)"
proof (induction m arbitrary: n)
case (0 n)
thus ?case by (subst exp_series_stream_aux.code) simp
next
case (Suc m n)
from Suc.IH[of "n + 1"] show ?case
by (subst exp_series_stream_aux.code) (simp add: algebra_simps)
qed

definition exp_series_stream :: "real msstream" where
"exp_series_stream = exp_series_stream_aux 1 1"

lemma mssnth_exp_series_stream:
"mssnth exp_series_stream n = 1 / fact n"
unfolding exp_series_stream_def using mssnth_exp_series_stream_aux[of 0 n] by simp

lemma convergent_powser'_exp:
"convergent_powser' (msllist_of_msstream exp_series_stream) exp"
proof -
have "(λn. lcoeff (msllist_of_msstream exp_series_stream) n * x ^ n) sums exp x" for x :: real
using exp_converges[of x] by (simp add: mssnth_exp_series_stream field_split_simps)
thus ?thesis by (auto intro: exI[of _ "1::real"] simp: convergent_powser'_def)
qed

subsubsection ‹Logarithm series›

primcorec ln_series_stream_aux :: "bool ⇒ real ⇒ real msstream" where
"ln_series_stream_aux b n =
MSSCons (if b then -inverse n else inverse n) (ln_series_stream_aux (¬b) (n+1))"

lemma mssnth_ln_series_stream_aux:
"mssnth (ln_series_stream_aux b x) n =
(if b then - ((-1)^n) * inverse (x + real n) else ((-1)^n) * inverse (x + real n))"
by (induction n arbitrary: b x; subst ln_series_stream_aux.code) simp_all

definition ln_series_stream :: "real msstream" where
"ln_series_stream = MSSCons 0 (ln_series_stream_aux False 1)"

lemma mssnth_ln_series_stream:
"mssnth ln_series_stream n = (-1) ^ Suc n / real n"
unfolding ln_series_stream_def
by (cases n) (simp_all add: mssnth_ln_series_stream_aux field_simps)

lemma ln_series':
assumes "abs (x::real) < 1"
shows   "(λn. - ((-x)^n) / of_nat n) sums ln (1 + x)"
proof -
have "∀n≥1. norm (-((-x)^n)) / of_nat n ≤ norm x ^ n / 1"
by (intro exI[of _ "1::nat"] allI impI frac_le) (simp_all add: power_abs)
hence "∃N. ∀n≥N. norm (-((-x)^n) / of_nat n) ≤ norm x ^ n"
by (intro exI[of _ 1]) simp_all
moreover from assms have "summable (λn. norm x ^ n)"
by (intro summable_geometric) simp_all
ultimately have *: "summable (λn. - ((-x)^n) / of_nat n)"
by (rule summable_comparison_test)
moreover from assms have "0 < 1 + x" "1 + x < 2" by simp_all
from ln_series[OF this]
have "ln (1 + x) = (∑n. - ((-x) ^ Suc n) / real (Suc n))"
hence "ln (1 + x) = (∑n. - ((-x) ^ n / real n))"
by (subst (asm) suminf_split_head[OF *]) simp_all
ultimately show ?thesis by (simp add: sums_iff)
qed

lemma convergent_powser'_ln:
"convergent_powser' (msllist_of_msstream ln_series_stream) (λx. ln (1 + x))"
proof -
have "(λn. lcoeff (msllist_of_msstream ln_series_stream)n * x ^ n) sums ln (1 + x)"
if "abs x < 1" for x
proof -
from that have "(λn. - ((- x) ^ n) / real n) sums ln (1 + x)" by (rule ln_series')
also have "(λn. - ((- x) ^ n) / real n) =
(λn. lcoeff (msllist_of_msstream ln_series_stream) n * x ^ n)"
by (auto simp: fun_eq_iff mssnth_ln_series_stream power_mult_distrib [symmetric])
finally show ?thesis .
qed
thus ?thesis unfolding convergent_powser'_def by (force intro!: exI[of _ 1])
qed

subsubsection ‹Generalized binomial series›

primcorec gbinomial_series_aux :: "bool ⇒ real ⇒ real ⇒ real ⇒ real msllist" where
"gbinomial_series_aux abort x n acc =
(if abort ∧ acc = 0 then MSLNil else
MSLCons acc (gbinomial_series_aux abort x (n + 1) ((x - n) / (n + 1) * acc)))"

lemma gbinomial_series_abort [simp]: "gbinomial_series_aux True x n 0 = MSLNil"
by (subst gbinomial_series_aux.code) simp

definition gbinomial_series :: "bool ⇒ real ⇒ real msllist" where
"gbinomial_series abort x = gbinomial_series_aux abort x 0 1"

(* TODO Move *)
lemma gbinomial_Suc_rec:
"(x gchoose (Suc k)) = (x gchoose k) * ((x - of_nat k) / (of_nat k + 1))"
proof -
have "((x - 1) + 1) gchoose (Suc k) = x * (x - 1 gchoose k) / (1 + of_nat k)"
by (subst gbinomial_factors) simp
also have "x * (x - 1 gchoose k) = (x - of_nat k) * (x gchoose k)"
by (rule gbinomial_absorb_comp [symmetric])
finally show ?thesis by (simp add: algebra_simps)
qed

lemma lcoeff_gbinomial_series_aux:
"lcoeff (gbinomial_series_aux abort x m (x gchoose m)) n = (x gchoose (n + m))"
proof (induction n arbitrary: m)
case 0
show ?case by (subst gbinomial_series_aux.code) simp
next
case (Suc n m)
show ?case
proof (cases "abort ∧ (x gchoose m) = 0")
case True
with gbinomial_mult_fact[of m x] obtain k where "x = real k" "k < m"
by auto
hence "(x gchoose Suc (n + m)) = 0"
using gbinomial_mult_fact[of "Suc (n + m)" x]
with True show ?thesis by (subst gbinomial_series_aux.code) simp
next
case False
hence "lcoeff (gbinomial_series_aux abort x m (x gchoose m)) (Suc n) =
lcoeff (gbinomial_series_aux abort x (Suc m) ((x gchoose m) *
((x - real m) / (real m + 1)))) n"
by (subst gbinomial_series_aux.code) (auto simp: algebra_simps)
also have "((x gchoose m) * ((x - real m) / (real m + 1))) = x gchoose (Suc m)"
by (rule gbinomial_Suc_rec [symmetric])
also have "lcoeff (gbinomial_series_aux abort x (Suc m) …) n = x gchoose (n + Suc m)"
by (rule Suc.IH)
finally show ?thesis by simp
qed
qed

lemma lcoeff_gbinomial_series [simp]:
"lcoeff (gbinomial_series abort x) n = (x gchoose n)"
using lcoeff_gbinomial_series_aux[of abort x 0 n] by (simp add: gbinomial_series_def)

lemma gbinomial_ratio_limit:
fixes a :: "'a :: real_normed_field"
assumes "a ∉ ℕ"
shows "(λn. (a gchoose n) / (a gchoose Suc n)) ⇢ -1"
proof (rule Lim_transform_eventually)
let ?f = "λn. inverse (a / of_nat (Suc n) - of_nat n / of_nat (Suc n))"
from eventually_gt_at_top[of "0::nat"]
show "eventually (λn. ?f n = (a gchoose n) /(a gchoose Suc n)) sequentially"
proof eventually_elim
fix n :: nat assume n: "n > 0"
then obtain q where q: "n = Suc q" by (cases n) blast
let ?P = "∏i=0..<n. a - of_nat i"
from n have "(a gchoose n) / (a gchoose Suc n) = (of_nat (Suc n) :: 'a) *
(?P / (∏i=0..n. a - of_nat i))"
also from q have "(∏i=0..n. a - of_nat i) = ?P * (a - of_nat n)"
also have "?P / … = (?P / ?P) / (a - of_nat n)" by (rule divide_divide_eq_left[symmetric])
also from assms have "?P / ?P = 1" by auto
also have "of_nat (Suc n) * (1 / (a - of_nat n)) =
inverse (inverse (of_nat (Suc n)) * (a - of_nat n))" by (simp add: field_simps)
also have "inverse (of_nat (Suc n)) * (a - of_nat n) =
a / of_nat (Suc n) - of_nat n / of_nat (Suc n)"
by (simp add: field_simps del: of_nat_Suc)
finally show "?f n = (a gchoose n) / (a gchoose Suc n)" by simp
qed

have "(λn. norm a / (of_nat (Suc n))) ⇢ 0"
unfolding divide_inverse
by (intro tendsto_mult_right_zero LIMSEQ_inverse_real_of_nat)
hence "(λn. a / of_nat (Suc n)) ⇢ 0"
by (subst tendsto_norm_zero_iff[symmetric]) (simp add: norm_divide del: of_nat_Suc)
hence "?f ⇢ inverse (0 - 1)"
by (intro tendsto_inverse tendsto_diff LIMSEQ_n_over_Suc_n) simp_all
thus "?f ⇢ -1" by simp
qed

lemma summable_gbinomial:
fixes a z :: real
assumes "norm z < 1"
shows   "summable (λn. (a gchoose n) * z ^ n)"
proof (cases "z = 0 ∨ a ∈ ℕ")
case False
define r where "r = (norm z + 1) / 2"
from assms have r: "r > norm z" "r < 1" by (simp_all add: r_def field_simps)
from False have "abs z > 0" by auto
from False have "a ∉ ℕ" by (auto elim!: Nats_cases)
hence *: "(λn. norm (z / ((a gchoose n) / (a gchoose (Suc n))))) ⇢ norm (z / (-1))"
by (intro tendsto_norm tendsto_divide tendsto_const gbinomial_ratio_limit) simp_all
hence "∀⇩F x in at_top. norm (z / ((a gchoose x) / (a gchoose Suc x))) > 0"
using assms False by (intro order_tendstoD) auto
hence nz: "∀⇩F x in at_top. (a gchoose x) ≠ 0" by eventually_elim auto

have "∀⇩F x in at_top. norm (z / ((a gchoose x) / (a gchoose Suc x))) < r"
using assms r by (intro order_tendstoD(2)[OF *]) simp_all
with nz have "∀⇩F n in at_top. norm ((a gchoose (Suc n)) * z) ≤ r * norm (a gchoose n)"
by eventually_elim (simp add: field_simps abs_mult split: if_splits)
hence "∀⇩F n in at_top. norm (z ^ n) * norm ((a gchoose (Suc n)) * z) ≤
norm (z ^ n) * (r * norm (a gchoose n))"
by eventually_elim (insert False, intro mult_left_mono, auto)
hence "∀⇩F n in at_top. norm ((a gchoose (Suc n)) * z ^ (Suc n)) ≤
r * norm ((a gchoose n) * z ^ n)"
then obtain N where N: "⋀n. n ≥ N ⟹ norm ((a gchoose (Suc n)) * z ^ (Suc n)) ≤
r * norm ((a gchoose n) * z ^ n)"
unfolding eventually_at_top_linorder by blast
show "summable (λn. (a gchoose n) * z ^ n)"
by (intro summable_ratio_test[OF r(2), of N] N)
next
case True
thus ?thesis
proof
assume "a ∈ ℕ"
then obtain n where [simp]: "a = of_nat n" by (auto elim: Nats_cases)
from eventually_gt_at_top[of n]
have "eventually (λn. (a gchoose n) * z ^ n = 0) at_top"
by eventually_elim (simp add: binomial_gbinomial [symmetric])
from summable_cong[OF this] show ?thesis by simp
qed auto
qed

lemma gen_binomial_real:
fixes z :: real
assumes "norm z < 1"
shows   "(λn. (a gchoose n) * z^n) sums (1 + z) powr a"
proof -
define K where "K = 1 - (1 - norm z) / 2"
from assms have K: "K > 0" "K < 1" "norm z < K"
unfolding K_def by (auto simp: field_simps intro!: add_pos_nonneg)
let ?f = "λn. a gchoose n" and ?f' = "diffs (λn. a gchoose n)"
have summable_strong: "summable (λn. ?f n * z ^ n)" if "norm z < 1" for z using that
by (intro summable_gbinomial)
with K have summable: "summable (λn. ?f n * z ^ n)" if "norm z < K" for z using that by auto
hence summable': "summable (λn. ?f' n * z ^ n)" if "norm z < K" for z using that
by (intro termdiff_converges[of _ K]) simp_all

define f f' where [abs_def]: "f z = (∑n. ?f n * z ^ n)" "f' z = (∑n. ?f' n * z ^ n)" for z
{
fix z :: real assume z: "norm z < K"
from summable_mult2[OF summable'[OF z], of z]
have summable1: "summable (λn. ?f' n * z ^ Suc n)" by (simp add: mult_ac)
hence summable2: "summable (λn. of_nat n * ?f n * z^n)"
unfolding diffs_def by (subst (asm) summable_Suc_iff)

have "(1 + z) * f' z = (∑n. ?f' n * z^n) + (∑n. ?f' n * z^Suc n)"
unfolding f_f'_def using summable' z by (simp add: algebra_simps suminf_mult)
also have "(∑n. ?f' n * z^n) = (∑n. of_nat (Suc n) * ?f (Suc n) * z^n)"
by (intro suminf_cong) (simp add: diffs_def)
also have "(∑n. ?f' n * z^Suc n) = (∑n. of_nat n * ?f n * z ^ n)"
using summable1 suminf_split_initial_segment[OF summable1] unfolding diffs_def
by (subst suminf_split_head, subst (asm) summable_Suc_iff) simp_all
also have "(∑n. of_nat (Suc n) * ?f (Suc n) * z^n) + (∑n. of_nat n * ?f n * z^n) =
(∑n. a * ?f n * z^n)"
(insert summable'[OF z] summable2,
also have "… = a * f z" unfolding f_f'_def
by (subst suminf_mult[symmetric]) (simp_all add: summable[OF z] mult_ac)
finally have "a * f z = (1 + z) * f' z" by simp
} note deriv = this

have [derivative_intros]: "(f has_field_derivative f' z) (at z)" if "norm z < of_real K" for z
unfolding f_f'_def using K that
by (intro termdiffs_strong[of "?f" K z] summable_strong) simp_all
have "f 0 = (∑n. if n = 0 then 1 else 0)" unfolding f_f'_def by (intro suminf_cong) simp
also have "… = 1" using sums_single[of 0 "λ_. 1::real"] unfolding sums_iff by simp
finally have [simp]: "f 0 = 1" .

have "f z * (1 + z) powr (-a) = f 0 * (1 + 0) powr (-a)"
proof (rule DERIV_isconst3[where ?f = "λx. f x * (1 + x) powr (-a)"])
fix z :: real assume z': "z ∈ {-K<..<K}"
hence z: "norm z < K" using K by auto
with K have nz: "1 + z ≠ 0" by (auto dest!: minus_unique)
from z K have "norm z < 1" by simp
hence "((λz. f z * (1 + z) powr (-a)) has_field_derivative
f' z * (1 + z) powr (-a) - a * f z * (1 + z) powr (-a-1)) (at z)" using z
by (auto intro!: derivative_eq_intros)
also from z have "a * f z = (1 + z) * f' z" by (rule deriv)
also have "f' z * (1 + z) powr - a - (1 + z) * f' z * (1 + z) powr (- a - 1) = 0"
using ‹norm z < 1› by (auto simp add: field_simps powr_diff)
finally show "((λz. f z * (1 + z) powr (-a)) has_field_derivative 0) (at z)" .
qed (insert K, auto)
also have "f 0 * (1 + 0) powr -a = 1" by simp
finally have "f z = (1 + z) powr a" using K
by (simp add: field_simps dist_real_def powr_minus)
with summable K show ?thesis unfolding f_f'_def by (simp add: sums_iff)
qed

lemma convergent_powser'_gbinomial:
"convergent_powser' (gbinomial_series abort p) (λx. (1 + x) powr p)"
proof -
have "(λn. lcoeff (gbinomial_series abort p) n * x ^ n) sums (1 + x) powr p" if "abs x < 1" for x
using that gen_binomial_real[of x p] by simp
thus ?thesis unfolding convergent_powser'_def by (force intro!: exI[of _ 1])
qed

lemma convergent_powser'_gbinomial':
"convergent_powser' (gbinomial_series abort (real n)) (λx. (1 + x) ^ n)"
proof -
{
fix x :: real
have "(λk. if k ∈ {0..n} then real (n choose k) * x ^ k else 0) sums (x + 1) ^ n"
using sums_If_finite_set[of "{..n}" "λk. real (n choose k) * x ^ k"]
by (subst binomial_ring) simp
also have "(λk. if k ∈ {0..n} then real (n choose k) * x ^ k else 0) =
(λk. lcoeff (gbinomial_series abort (real n)) k * x ^ k)"
by (simp add: fun_eq_iff binomial_gbinomial [symmetric])
finally have "… sums (1 + x) ^ n" by (simp add: add_ac)
}
thus ?thesis unfolding convergent_powser'_def
by (auto intro!: exI[of _ 1])
qed

subsubsection ‹Sine and cosine›

primcorec sin_series_stream_aux :: "bool ⇒ real ⇒ real ⇒ real msstream" where
"sin_series_stream_aux b acc m =
MSSCons (if b then -inverse acc else inverse acc)
(sin_series_stream_aux (¬b) (acc * (m + 1) * (m + 2)) (m + 2))"

lemma mssnth_sin_series_stream_aux:
"mssnth (sin_series_stream_aux b (fact m) m) n =
(if b then -1 else 1) * (-1) ^ n / (fact (2 * n + m))"
proof (induction n arbitrary: b m)
case (0 b m)
thus ?case by (subst sin_series_stream_aux.code) (simp add: field_simps)
next
case (Suc n b m)
show ?case using Suc.IH[of "¬b" "m + 2"]
by (subst sin_series_stream_aux.code) (auto simp: field_simps)
qed

definition sin_series_stream where
"sin_series_stream = sin_series_stream_aux False 1 1"

lemma mssnth_sin_series_stream:
"mssnth sin_series_stream n = (-1) ^ n / fact (2 * n + 1)"
using mssnth_sin_series_stream_aux[of False 1 n] unfolding sin_series_stream_def by simp

definition cos_series_stream where
"cos_series_stream = sin_series_stream_aux False 1 0"

lemma mssnth_cos_series_stream:
"mssnth cos_series_stream n = (-1) ^ n / fact (2 * n)"
using mssnth_sin_series_stream_aux[of False 0 n] unfolding cos_series_stream_def by simp

lemma summable_pre_sin: "summable (λn. mssnth sin_series_stream n * (x::real) ^ n)"
proof -
have *: "summable (λn. 1 / fact n * abs x ^ n)"
using exp_converges[of "abs x"] by (simp add: sums_iff field_simps)
{
fix n :: nat
have "¦x¦ ^ n / fact (2 * n + 1) ≤ ¦x¦ ^ n / fact n"
by (intro divide_left_mono fact_mono) auto
hence "norm (mssnth sin_series_stream n * x ^ n) ≤ 1 / fact n * abs x ^ n"
by (simp add: mssnth_sin_series_stream abs_mult power_abs field_simps)
}
thus ?thesis
by (intro summable_comparison_test[OF _ *]) (auto intro!: exI[of _ 0])
qed

lemma summable_pre_cos: "summable (λn. mssnth cos_series_stream n * (x::real) ^ n)"
proof -
have *: "summable (λn. 1 / fact n * abs x ^ n)"
using exp_converges[of "abs x"] by (simp add: sums_iff field_simps)
{
fix n :: nat
have "¦x¦ ^ n / fact (2 * n) ≤ ¦x¦ ^ n / fact n"
by (intro divide_left_mono fact_mono) auto
hence "norm (mssnth cos_series_stream n * x ^ n) ≤ 1 / fact n * abs x ^ n"
by (simp add: mssnth_cos_series_stream abs_mult power_abs field_simps)
}
thus ?thesis
by (intro summable_comparison_test[OF _ *]) (auto intro!: exI[of _ 0])
qed

lemma cos_conv_pre_cos:
"cos x = powser (msllist_of_msstream cos_series_stream) (x ^ 2)"
proof -
have "(λn. cos_coeff (2 * n) * x ^ (2 * n)) sums cos x"
using cos_converges[of x]
by (subst sums_mono_reindex[of "λn. 2 * n"])
(auto simp: strict_mono_def cos_coeff_def elim!: evenE)
also have "(λn. cos_coeff (2 * n) * x ^ (2 * n)) =
(λn. mssnth cos_series_stream n * (x ^ 2) ^ n)"
by (simp add: fun_eq_iff mssnth_cos_series_stream cos_coeff_def power_mult)
finally have sums: "(λn. mssnth cos_series_stream n * x⇧2 ^ n) sums cos x" .
thus ?thesis by (simp add: powser_def sums_iff)
qed

lemma sin_conv_pre_sin:
"sin x = x * powser (msllist_of_msstream sin_series_stream) (x ^ 2)"
proof -
have "(λn. sin_coeff (2 * n + 1) * x ^ (2 * n + 1)) sums sin x"
using sin_converges[of x]
by (subst sums_mono_reindex[of "λn. 2 * n + 1"])
(auto simp: strict_mono_def sin_coeff_def elim!: oddE)
also have "(λn. sin_coeff (2 * n + 1) * x ^ (2 * n + 1)) =
(λn. x * mssnth sin_series_stream n * (x ^ 2) ^ n)"
by (simp add: fun_eq_iff mssnth_sin_series_stream sin_coeff_def power_mult [symmetric] algebra_simps)
finally have sums: "(λn. x * mssnth sin_series_stream n * x⇧2 ^ n) sums sin x" .
thus ?thesis using summable_pre_sin[of "x^2"]
by (simp add: powser_def sums_iff suminf_mult [symmetric] mult.assoc)
qed

lemma convergent_powser_sin_series_stream:
"convergent_powser (msllist_of_msstream sin_series_stream)"
(is "convergent_powser ?cs")
proof -
show ?thesis using summable_pre_sin unfolding convergent_powser_def
by (intro exI[of _ 1]) auto
qed

lemma convergent_powser_cos_series_stream:
"convergent_powser (msllist_of_msstream cos_series_stream)"
(is "convergent_powser ?cs")
proof -
show ?thesis using summable_pre_cos unfolding convergent_powser_def
by (intro exI[of _ 1]) auto
qed

subsubsection ‹Arc tangent›

definition arctan_coeffs :: "nat ⇒ 'a :: {real_normed_div_algebra, banach}" where
"arctan_coeffs n = (if odd n then (-1) ^ (n div 2) / of_nat n else 0)"

lemma arctan_converges:
assumes "norm x < 1"
shows   "(λn. arctan_coeffs n * x ^ n) sums arctan x"
proof -
have A: "(λn. diffs arctan_coeffs n * x ^ n) sums (1 / (1 + x^2))" if "norm x < 1" for x :: real
proof -
let ?f = "λn. (if odd n then 0 else (-1) ^ (n div 2)) * x ^ n"
from that have "norm x ^ 2 < 1 ^ 2" by (intro power_strict_mono) simp_all
hence "(λn. (-(x^2))^n) sums (1 / (1 - (-(x^2))))" by (intro geometric_sums) simp_all
also have "1 - (-(x^2)) = 1 + x^2" by simp
also have "(λn. (-(x^2))^n) = (λn. ?f (2*n))" by (auto simp: fun_eq_iff power_minus' power_mult)
also have "range ((*) (2::nat)) = {n. even n}"
by (auto elim!: evenE)
hence "(λn. ?f (2*n)) sums (1 / (1 + x^2)) ⟷ ?f sums (1 / (1 + x^2))"
by (intro sums_mono_reindex) (simp_all add: strict_mono_Suc_iff)
also have "?f = (λn. diffs arctan_coeffs n * x ^ n)"
by (simp add: fun_eq_iff arctan_coeffs_def diffs_def)
finally show "(λn. diffs arctan_coeffs n * x ^ n) sums (1 / (1 + x^2))" .
qed

have B: "summable (λn. arctan_coeffs n * x ^ n)" if "norm x < 1" for x :: real
proof (rule summable_comparison_test)
show "∃N. ∀n≥N. norm (arctan_coeffs n * x ^ n) ≤ 1 * norm x ^ n"
unfolding norm_mult norm_power
by (intro exI[of _ "0::nat"] allI impI mult_right_mono)
from that show "summable (λn. 1 * norm x ^ n)"
by (intro summable_mult summable_geometric) simp_all
qed

define F :: "real ⇒ real" where "F = (λx. ∑n. arctan_coeffs n * x ^ n)"
have [derivative_intros]:
"(F has_field_derivative (1 / (1 + x^2))) (at x)" if "norm x < 1" for x :: real
proof -
define K :: real where "K = (1 + norm x) / 2"
from that have K: "norm x < K" "K < 1" by (simp_all add: K_def field_simps)
have "(F has_field_derivative (∑n. diffs arctan_coeffs n * x ^ n)) (at x)"
unfolding F_def using K by (intro termdiffs_strong[where K = K] B) simp_all
also from A[OF that] have "(∑n. diffs arctan_coeffs n * x ^ n) = 1 / (1 + x^2)"
finally show ?thesis .
qed
from assms have "arctan x - F x = arctan 0 - F 0"
by (intro DERIV_isconst3[of "-1" 1 _ _ "λx. arctan x - F x"])
(auto intro!: derivative_eq_intros simp: field_simps)
hence "F x = arctan x" by (simp add: F_def arctan_coeffs_def)
with B[OF assms] show ?thesis by (simp add: sums_iff F_def)
qed

primcorec arctan_series_stream_aux :: "bool ⇒ real ⇒ real msstream" where
"arctan_series_stream_aux b n =
MSSCons (if b then -inverse n else inverse n) (arctan_series_stream_aux (¬b) (n + 2))"

lemma mssnth_arctan_series_stream_aux:
"mssnth (arctan_series_stream_aux b n) m = (-1) ^ (if b then Suc m else m) / (2*m + n)"
by (induction m arbitrary: b n; subst arctan_series_stream_aux.code)
(auto simp: field_simps split: if_splits)

definition arctan_series_stream where
"arctan_series_stream = arctan_series_stream_aux False 1"

lemma mssnth_arctan_series_stream:
"mssnth arctan_series_stream n = (-1) ^ n / (2 * n + 1)"

lemma summable_pre_arctan:
assumes "norm x < 1"
shows   "summable (λn. mssnth arctan_series_stream n * x ^ n)" (is "summable ?f")
proof (rule summable_comparison_test)
show "summable (λn. abs x ^ n)" using assms by (intro summable_geometric) auto
show "∃N. ∀n≥N. norm (?f n) ≤ abs x ^ n"
proof (intro exI[of _ 0] allI impI)
fix n :: nat
have "norm (?f n) = ¦x¦ ^ n / (2 * real n + 1)"
by (simp add: mssnth_arctan_series_stream abs_mult power_abs)
also have "… ≤ ¦x¦ ^ n / 1" by (intro divide_left_mono) auto
finally show "norm (?f n) ≤ abs x ^ n" by simp
qed
qed

lemma arctan_conv_pre_arctan:
assumes "norm x < 1"
shows   "arctan x = x * powser (msllist_of_msstream arctan_series_stream) (x ^ 2)"
proof -
from assms have "norm x * norm x < 1 * 1" by (intro mult_strict_mono) auto
hence norm_less: "norm (x ^ 2) < 1" by (simp add: power2_eq_square)
have "(λn. arctan_coeffs n * x ^ n) sums arctan x"
by (intro arctan_converges assms)
also have "?this ⟷ (λn. arctan_coeffs (2 * n + 1) * x ^ (2 * n + 1)) sums arctan x"
by (intro sums_mono_reindex [symmetric])
(auto simp: arctan_coeffs_def strict_mono_def elim!: oddE)
also have "(λn. arctan_coeffs (2 * n + 1) * x ^ (2 * n + 1)) =
(λn. x * mssnth arctan_series_stream n * (x ^ 2) ^ n)"
power_mult [symmetric] arctan_coeffs_def mult_ac)
finally have "(λn. x * mssnth arctan_series_stream n * x⇧2 ^ n) sums arctan x" .
thus ?thesis using summable_pre_arctan[OF norm_less] assms
by (simp add: powser_def sums_iff suminf_mult [symmetric] mult.assoc)
qed

lemma convergent_powser_arctan:
"convergent_powser (msllist_of_msstream arctan_series_stream)"
unfolding convergent_powser_def
using summable_pre_arctan by (intro exI[of _ 1]) auto

lemma arctan_inverse_pos: "x > 0 ⟹ arctan x = pi / 2 - arctan (inverse x)"
using arctan_inverse[of x] by (simp add: field_simps)

lemma arctan_inverse_neg: "x < 0 ⟹ arctan x = -pi / 2 - arctan (inverse x)"
using arctan_inverse[of x] by (simp add: field_simps)

subsection ‹Multiseries›

subsubsection ‹Asymptotic bases›

type_synonym basis = "(real ⇒ real) list"

lemma transp_ln_smallo_ln: "transp (λf g. (λx::real. ln (g x)) ∈ o(λx. ln (f x)))"
by (rule transpI, erule landau_o.small.trans)

definition basis_wf :: "basis ⇒ bool" where
"basis_wf basis ⟷ (∀f∈set basis. filterlim f at_top at_top) ∧
sorted_wrt (λf g. (λx. ln (g x)) ∈ o(λx. ln (f x))) basis"

lemma basis_wf_Nil [simp]: "basis_wf []"

lemma basis_wf_Cons:
"basis_wf (f # basis) ⟷
filterlim f at_top at_top ∧ (∀g∈set basis. (λx. ln (g x)) ∈ o(λx. ln (f x))) ∧ basis_wf basis"
unfolding basis_wf_def by auto

(* TODO: Move *)
lemma sorted_wrt_snoc:
assumes trans: "transp P" and "sorted_wrt P xs" "P (last xs) y"
shows   "sorted_wrt P (xs @ [y])"
proof (subst sorted_wrt_append, intro conjI ballI)
fix x y' assume x: "x ∈ set xs" and y': "y' ∈ set [y]"
from y' have [simp]: "y' = y" by simp
show "P x y'"
proof (cases "x = last xs")
case False
from x have eq: "xs = butlast xs @ [last xs]"
by (subst append_butlast_last_id) auto
from x and False have x': "x ∈ set (butlast xs)" by (subst (asm) eq) auto
have "sorted_wrt P xs" by fact
also note eq
finally have "P x (last xs)" using x'
by (subst (asm) sorted_wrt_append) auto
with ‹P (last xs) y› have "P x y" using transpD[OF trans] by blast
thus ?thesis by simp
qed (insert assms, auto)
qed (insert assms, auto)

lemma basis_wf_snoc:
assumes "bs ≠ []"
assumes "basis_wf bs" "filterlim b at_top at_top"
assumes "(λx. ln (b x)) ∈ o(λx. ln (last bs x))"
shows   "basis_wf (bs @ [b])"
proof -
have "transp (λf g. (λx. ln (g x)) ∈ o(λx. ln (f x)))"
by (auto simp: transp_def intro: landau_o.small.trans)
thus ?thesis using assms
by (auto simp add: basis_wf_def sorted_wrt_snoc[OF transp_ln_smallo_ln])
qed

lemma basis_wf_singleton: "basis_wf [b] ⟷ filterlim b at_top at_top"

lemma basis_wf_many: "basis_wf (b # b' # bs) ⟷
filterlim b at_top at_top ∧ (λx. ln (b' x)) ∈ o(λx. ln (b x)) ∧ basis_wf (b' # bs)"
unfolding basis_wf_def by (subst sorted_wrt2[OF transp_ln_smallo_ln]) auto

subsubsection ‹Monomials›

type_synonym monom = "real × real list"

definition eval_monom :: "monom ⇒ basis ⇒ real ⇒ real" where
"eval_monom = (λ(c,es) basis x. c * prod_list (map (λ(b,e). b x powr e) (zip basis es)))"

lemma eval_monom_Nil1 [simp]: "eval_monom (c, []) basis = (λ_. c)"
by (simp add: eval_monom_def split: prod.split)

lemma eval_monom_Nil2 [simp]: "eval_monom monom [] = (λ_. fst monom)"
by (simp add: eval_monom_def split: prod.split)

lemma eval_monom_Cons:
"eval_monom (c, e # es) (b # basis) = (λx. eval_monom (c, es) basis x * b x powr e)"
by (simp add: eval_monom_def fun_eq_iff split: prod.split)

definition inverse_monom :: "monom ⇒ monom" where
"inverse_monom monom = (case monom of (c, es) ⇒ (inverse c, map uminus es))"

lemma length_snd_inverse_monom [simp]:
"length (snd (inverse_monom monom)) = length (snd monom)"
by (simp add: inverse_monom_def split: prod.split)

lemma eval_monom_pos:
assumes "basis_wf basis" "fst monom > 0"
shows   "eventually (λx. eval_monom monom basis x > 0) at_top"
proof (cases monom)
case (Pair c es)
with assms have "c > 0" by simp
with assms(1) show ?thesis unfolding Pair
proof (induction es arbitrary: basis)
case (Cons e es)
note A = this
thus ?case
proof (cases basis)
case (Cons b basis')
with A(1)[of basis'] A(2,3)
have A: "∀⇩F x in at_top. eval_monom (c, es) basis' x > 0"
and B: "eventually (λx. b x > 0) at_top"
thus ?thesis by eventually_elim (simp add: Cons eval_monom_Cons)
qed simp
qed simp
qed

lemma eval_monom_uminus: "eval_monom (-c, es) basis x = -eval_monom (c, es) basis x"

lemma eval_monom_neg:
assumes "basis_wf basis" "fst monom < 0"
shows   "eventually (λx. eval_monom monom basis x < 0) at_top"
proof -
from assms have "eventually (λx. eval_monom (-fst monom, snd monom) basis x > 0) at_top"
by (intro eval_monom_pos) simp_all
thus ?thesis by (simp add: eval_monom_uminus)
qed

lemma eval_monom_nonzero:
assumes "basis_wf basis" "fst monom ≠ 0"
shows   "eventually (λx. eval_monom monom basis x ≠ 0) at_top"
proof (cases "fst monom" "0 :: real" rule: linorder_cases)
case greater
with assms have "eventually (λx. eval_monom monom basis x > 0) at_top"
thus ?thesis by eventually_elim simp
next
case less
with assms have "eventually (λx. eval_monom (-fst monom, snd monom) basis x > 0) at_top"
thus ?thesis by eventually_elim (simp add: eval_monom_uminus)
qed (insert assms, simp_all)

subsubsection ‹Typeclass for multiseries›

class multiseries = plus + minus + times + uminus + inverse +
fixes is_expansion :: "'a ⇒ basis ⇒ bool"
and expansion_level :: "'a itself ⇒ nat"
and eval :: "'a ⇒ real ⇒ real"
and zero_expansion :: 'a
and const_expansion :: "real ⇒ 'a"
and powr_expansion :: "bool ⇒ 'a ⇒ real ⇒ 'a"
and power_expansion :: "bool ⇒ 'a ⇒ nat ⇒ 'a"
and trimmed :: "'a ⇒ bool"
and dominant_term :: "'a ⇒ monom"

assumes is_expansion_length:
"is_expansion F basis ⟹ length basis = expansion_level (TYPE('a))"
assumes is_expansion_zero:
"basis_wf basis ⟹ length basis = expansion_level (TYPE('a)) ⟹
is_expansion zero_expansion basis"
assumes is_expansion_const:
"basis_wf basis ⟹ length basis = expansion_level (TYPE('a)) ⟹
is_expansion (const_expansion c) basis"
assumes is_expansion_uminus:
"basis_wf basis ⟹ is_expansion F basis ⟹ is_expansion (-F) basis"
"basis_wf basis ⟹ is_expansion F basis ⟹ is_expansion G basis ⟹
is_expansion (F + G) basis"
assumes is_expansion_minus:
"basis_wf basis ⟹ is_expansion F basis ⟹ is_expansion G basis ⟹
is_expansion (F - G) basis"
assumes is_expansion_mult:
"basis_wf basis ⟹ is_expansion F basis ⟹ is_expansion G basis ⟹
is_expansion (F * G) basis"
assumes is_expansion_inverse:
"basis_wf basis ⟹ trimmed F ⟹ is_expansion F basis ⟹
is_expansion (inverse F) basis"
assumes is_expansion_divide:
"basis_wf basis ⟹ trimmed G ⟹ is_expansion F basis ⟹ is_expansion G basis ⟹
is_expansion (F / G) basis"
assumes is_expansion_powr:
"basis_wf basis ⟹ trimmed F ⟹ fst (dominant_term F) > 0 ⟹ is_expansion F basis ⟹
is_expansion (powr_expansion abort F p) basis"
assumes is_expansion_power:
"basis_wf basis ⟹ trimmed F ⟹ is_expansion F basis ⟹
is_expansion (power_expansion abort F n) basis"

assumes is_expansion_imp_smallo:
"basis_wf basis ⟹ is_expansion F basis ⟹ filterlim b at_top at_top ⟹
(∀g∈set basis. (λx. ln (g x)) ∈ o(λx. ln (b x))) ⟹ e > 0 ⟹ eval F ∈ o(λx. b x powr e)"
assumes is_expansion_imp_smallomega:
"basis_wf basis ⟹ is_expansion F basis ⟹ filterlim b at_top at_top ⟹ trimmed F ⟹
(∀g∈set basis. (λx. ln (g x)) ∈ o(λx. ln (b x))) ⟹ e < 0 ⟹ eval F ∈ ω(λx. b x powr e)"
assumes trimmed_imp_eventually_sgn:
"basis_wf basis ⟹ is_expansion F basis ⟹ trimmed F ⟹
eventually (λx. sgn (eval F x) = sgn (fst (dominant_term F))) at_top"
assumes trimmed_imp_eventually_nz:
"basis_wf basis ⟹ is_expansion F basis ⟹ trimmed F ⟹
eventually (λx. eval F x ≠ 0) at_top"
assumes trimmed_imp_dominant_term_nz: "trimmed F ⟹ fst (dominant_term F) ≠ 0"

assumes dominant_term:
"basis_wf basis ⟹ is_expansion F basis ⟹ trimmed F ⟹
eval F ∼[at_top] eval_monom (dominant_term F) basis"
assumes dominant_term_bigo:
"basis_wf basis ⟹ is_expansion F basis ⟹
eval F ∈ O(eval_monom (1, snd (dominant_term F)) basis)"
assumes length_dominant_term:
"length (snd (dominant_term F)) = expansion_level TYPE('a)"
assumes fst_dominant_term_uminus [simp]: "fst (dominant_term (- F)) = -fst (dominant_term F)"
assumes trimmed_uminus_iff [simp]: "trimmed (-F) ⟷ trimmed F"

assumes add_zero_expansion_left [simp]: "zero_expansion + F = F"
assumes add_zero_expansion_right [simp]: "F + zero_expansion = F"

assumes eval_zero [simp]: "eval zero_expansion x = 0"
assumes eval_const [simp]: "eval (const_expansion c) x = c"
assumes eval_uminus [simp]: "eval (-F) = (λx. -eval F x)"
assumes eval_plus [simp]: "eval (F + G) = (λx. eval F x + eval G x)"
assumes eval_minus [simp]: "eval (F - G) = (λx. eval F x - eval G x)"
assumes eval_times [simp]: "eval (F * G) = (λx. eval F x * eval G x)"
assumes eval_inverse [simp]: "eval (inverse F) = (λx. inverse (eval F x))"
assumes eval_divide [simp]: "eval (F / G) = (λx. eval F x / eval G x)"
assumes eval_powr [simp]: "eval (powr_expansion abort F p) = (λx. eval F x powr p)"
assumes eval_power [simp]: "eval (power_expansion abort F n) = (λx. eval F x ^ n)"
assumes minus_eq_plus_uminus: "F - G = F + (-G)"
assumes times_const_expansion_1: "const_expansion 1 * F = F"
assumes trimmed_const_expansion: "trimmed (const_expansion c) ⟷ c ≠ 0"
begin

definition trimmed_pos where
"trimmed_pos F ⟷ trimmed F ∧ fst (dominant_term F) > 0"

definition trimmed_neg where
"trimmed_neg F ⟷ trimmed F ∧ fst (dominant_term F) < 0"

lemma trimmed_pos_uminus: "trimmed_neg F ⟹ trimmed_pos (-F)"

lemma trimmed_pos_imp_trimmed: "trimmed_pos x ⟹ trimmed x"

lemma trimmed_neg_imp_trimmed: "trimmed_neg x ⟹ trimmed x"

end

subsubsection ‹Zero-rank expansions›

instantiation real :: multiseries
begin

inductive is_expansion_real :: "real ⇒ basis ⇒ bool" where
"is_expansion_real c []"

definition expansion_level_real :: "real itself ⇒ nat" where
expansion_level_real_def [simp]: "expansion_level_real _ = 0"

definition eval_real :: "real ⇒ real ⇒ real" where
eval_real_def [simp]: "eval_real c = (λ_. c)"

definition zero_expansion_real :: "real" where
zero_expansion_real_def [simp]: "zero_expansion_real = 0"

definition const_expansion_real :: "real ⇒ real" where
const_expansion_real_def [simp]: "const_expansion_real x = x"

definition powr_expansion_real :: "bool ⇒ real ⇒ real ⇒ real" where
powr_expansion_real_def [simp]: "powr_expansion_real _ x p = x powr p"

definition power_expansion_real :: "bool ⇒ real ⇒ nat ⇒ real" where
power_expansion_real_def [simp]: "power_expansion_real _ x n = x ^ n"

definition trimmed_real :: "real ⇒ bool" where
trimmed_real_def [simp]: "trimmed_real x ⟷ x ≠ 0"

definition dominant_term_real :: "real ⇒ monom" where
dominant_term_real_def [simp]: "dominant_term_real c = (c, [])"

instance proof
fix basis :: basis and b :: "real ⇒ real" and e F :: real
assume *: "basis_wf basis" "filterlim b at_top at_top" "is_expansion F basis" "0 < e"
have "(λx. b x powr e) ∈ ω(λ_. 1)"
by (intro smallomegaI_filterlim_at_infinity filterlim_at_top_imp_at_infinity)
(auto intro!: filterlim_compose[OF real_powr_at_top] * )
with * show "eval F ∈ o(λx. b x powr e)"
by (cases "F = 0") (auto elim!: is_expansion_real.cases simp: smallomega_iff_smallo)
next
fix basis :: basis and b :: "real ⇒ real" and e F :: real
assume *: "basis_wf basis" "filterlim b at_top at_top" "is_expansion F basis"
"e < 0" "trimmed F"
from * have pos: "eventually (λx. b x > 0) at_top" by (simp add: filterlim_at_top_dense)
have "(λx. b x powr -e) ∈ ω(λ_. 1)"
by (intro smallomegaI_filterlim_at_infinity filterlim_at_top_imp_at_infinity)
(auto intro!: filterlim_compose[OF real_powr_at_top] * )
with pos have "(λx. b x powr e) ∈ o(λ_. 1)" unfolding powr_minus
by (subst (asm) landau_omega.small.inverse_eq2)
(auto elim: eventually_mono simp: smallomega_iff_smallo)
with * show "eval F ∈ ω(λx. b x powr e)"
by (auto elim!: is_expansion_real.cases simp: smallomega_iff_smallo)
qed (auto intro!: is_expansion_real.intros elim!: is_expansion_real.cases)

end

lemma eval_real: "eval (c :: real) x = c" by simp

subsubsection ‹Operations on multiseries›

lemma ms_aux_cases [case_names MSLNil MSLCons]:
fixes xs :: "('a × real) msllist"
obtains "xs = MSLNil" | c e xs' where "xs = MSLCons (c, e) xs'"
proof (cases xs)
case (MSLCons x xs')
with that(2)[of "fst x" "snd x" xs'] show ?thesis by auto
qed auto

definition ms_aux_hd_exp :: "('a × real) msllist ⇒ real option" where
"ms_aux_hd_exp xs = (case xs of MSLNil ⇒ None | MSLCons (_, e) _ ⇒ Some e)"

primrec ms_exp_gt :: "real ⇒ real option ⇒ bool" where
"ms_exp_gt x None = True"
| "ms_exp_gt x (Some y) ⟷ x > y"

lemma ms_aux_hd_exp_MSLNil [simp]: "ms_aux_hd_exp MSLNil = None"
by (simp add: ms_aux_hd_exp_def split: prod.split)

lemma ms_aux_hd_exp_MSLCons [simp]: "ms_aux_hd_exp (MSLCons x xs) = Some (snd x)"
by (simp add: ms_aux_hd_exp_def split: prod.split)

coinductive is_expansion_aux ::
"('a :: multiseries × real) msllist ⇒ (real ⇒ real) ⇒ basis ⇒ bool" where
is_expansion_MSLNil:
"eventually (λx. f x = 0) at_top ⟹ length basis = Suc (expansion_level TYPE('a)) ⟹
is_expansion_aux MSLNil f basis"
| is_expansion_MSLCons:
"is_expansion_aux xs (λx. f x - eval C x * b x powr e) (b # basis) ⟹
is_expansion C basis ⟹
(⋀e'. e' > e ⟹ f ∈ o(λx. b x powr e')) ⟹ ms_exp_gt e (ms_aux_hd_exp xs) ⟹
is_expansion_aux (MSLCons (C, e) xs) f (b # basis)"

inductive_cases is_expansion_aux_MSLCons: "is_expansion_aux (MSLCons (c, e) xs) f basis"

lemma is_expansion_aux_basis_nonempty: "is_expansion_aux F f basis ⟹ basis ≠ []"
by (erule is_expansion_aux.cases) auto

lemma is_expansion_aux_expansion_level:
assumes "is_expansion_aux (G :: ('a::multiseries × real) msllist) g basis"
shows   "length basis = Suc (expansion_level TYPE('a))"
using assms by (cases rule: is_expansion_aux.cases) (auto dest: is_expansion_length)

lemma is_expansion_aux_imp_smallo:
assumes "is_expansion_aux xs f basis" "ms_exp_gt p (ms_aux_hd_exp xs)"
shows   "f ∈ o(λx. hd basis x powr p)"
using assms
proof (cases rule: is_expansion_aux.cases)
case is_expansion_MSLNil
show ?thesis by (simp add: landau_o.small.in_cong[OF is_expansion_MSLNil(2)])
next
case (is_expansion_MSLCons xs C b e basis)
with assms have "f ∈ o(λx. b x powr p)"
by (intro is_expansion_MSLCons) (simp_all add: ms_aux_hd_exp_def)
thus ?thesis by (simp add: is_expansion_MSLCons)
qed

lemma is_expansion_aux_imp_smallo':
assumes "is_expansion_aux xs f basis"
obtains e where "f ∈ o(λx. hd basis x powr e)"
proof -
define e where "e = (case ms_aux_hd_exp xs of None ⇒ 0 | Some e ⇒ e + 1)"
have "ms_exp_gt e (ms_aux_hd_exp xs)"
by (auto simp add: e_def ms_aux_hd_exp_def split: msllist.splits)
from assms this have "f ∈ o(λx. hd basis x powr e)" by (rule is_expansion_aux_imp_smallo)
from that[OF this] show ?thesis .
qed

lemma is_expansion_aux_imp_smallo'':
assumes "is_expansion_aux xs f basis" "ms_exp_gt e' (ms_aux_hd_exp xs)"
obtains e where "e < e'" "f ∈ o(λx. hd basis x powr e)"
proof -
define e where "e = (case ms_aux_hd_exp xs of None ⇒ e' - 1 | Some e ⇒ (e' + e) / 2)"
from assms have "ms_exp_gt e (ms_aux_hd_exp xs)" "e < e'"
by (cases xs; simp add: e_def field_simps)+
from assms(1) this(1) have "f ∈ o(λx. hd basis x powr e)" by (rule is_expansion_aux_imp_smallo)
from that[OF ‹e < e'› this] show ?thesis .
qed

definition trimmed_ms_aux :: "('a :: multiseries × real) msllist ⇒ bool" where
"trimmed_ms_aux xs = (case xs of MSLNil ⇒ False | MSLCons (C, _) _ ⇒ trimmed C)"

lemma not_trimmed_ms_aux_MSLNil [simp]: "¬trimmed_ms_aux MSLNil"

lemma trimmed_ms_aux_MSLCons: "trimmed_ms_aux (MSLCons x xs) = trimmed (fst x)"
by (simp add: trimmed_ms_aux_def split: prod.split)

lemma trimmed_ms_aux_imp_nz:
assumes "basis_wf basis" "is_expansion_aux xs f basis" "trimmed_ms_aux xs"
shows   "eventually (λx. f x ≠ 0) at_top"
proof (cases xs rule: ms_aux_cases)
case (MSLCons C e xs')
from assms this obtain b basis' where *: "basis = b # basis'"
"is_expansion_aux xs' (λx. f x - eval C x * b x powr e) (b # basis')"
"ms_exp_gt e (ms_aux_hd_exp xs')" "is_expansion C basis'" "⋀e'. e' > e ⟹ f ∈ o(λx. b x powr e')"
by (auto elim!: is_expansion_aux_MSLCons)
from assms(1,3) * have nz: "eventually (λx. eval C x ≠ 0) at_top"
by (auto simp: trimmed_ms_aux_def MSLCons basis_wf_Cons
intro: trimmed_imp_eventually_nz[of basis'])
from assms(1) * have b_limit: "filterlim b at_top at_top" by (simp add: basis_wf_Cons)
hence b_nz: "eventually (λx. b x > 0) at_top" by (simp```