1signature MULTISERIES_EXPANSION = sig
2
3type expansion_thm = thm
4type trimmed_thm = thm
5type expr = Exp_Log_Expression.expr
6type basis = Asymptotic_Basis.basis
7
8datatype trim_mode = Simple_Trim | Pos_Trim | Neg_Trim | Sgn_Trim
9
10datatype zeroness = IsZero | IsNonZero | IsPos | IsNeg
11
12datatype intyness = Nat of thm | Neg_Nat of thm | No_Nat
13datatype parity = Even of thm | Odd of thm | Unknown_Parity
14
15datatype limit =
16   Zero_Limit of bool option
17 | Finite_Limit of term
18 | Infinite_Limit of bool option
19
20datatype trim_result =
21    Trimmed of zeroness * trimmed_thm option
22  | Aborted of order
23
24val get_intyness : Proof.context -> cterm -> intyness
25val get_parity : cterm -> parity
26
27val get_expansion : thm -> term
28val get_coeff : term -> term
29val get_exponent : term -> term
30val get_expanded_fun : thm -> term
31val get_eval : term -> term
32val expands_to_hd : thm -> thm
33
34val mk_eval_ctxt : Proof.context -> Lazy_Eval.eval_ctxt
35val expand : Lazy_Eval.eval_ctxt -> expr -> basis -> expansion_thm * basis
36val expand_term : Lazy_Eval.eval_ctxt -> term -> basis -> expansion_thm * basis
37val expand_terms : Lazy_Eval.eval_ctxt -> term list -> basis -> expansion_thm list * basis
38
39val limit_of_expansion : bool * bool -> Lazy_Eval.eval_ctxt -> thm * basis -> limit * thm
40val compute_limit : Lazy_Eval.eval_ctxt -> term -> limit * thm
41val compare_expansions : 
42  Lazy_Eval.eval_ctxt -> expansion_thm * expansion_thm * basis -> 
43    order * thm * expansion_thm * expansion_thm
44
45(* TODO DEBUG *)
46datatype comparison_result =
47  Cmp_Dominated of order * thm list * zeroness * trimmed_thm * expansion_thm * expansion_thm 
48| Cmp_Asymp_Equiv of thm * thm
49val compare_expansions' :
50  Lazy_Eval.eval_ctxt ->
51      thm * thm * Asymptotic_Basis.basis ->
52        comparison_result
53
54val prove_at_infinity : Lazy_Eval.eval_ctxt -> thm * basis -> thm
55val prove_at_top : Lazy_Eval.eval_ctxt -> thm * basis -> thm
56val prove_at_bot : Lazy_Eval.eval_ctxt -> thm * basis -> thm
57val prove_nhds : Lazy_Eval.eval_ctxt -> thm * basis -> thm
58val prove_at_0 : Lazy_Eval.eval_ctxt -> thm * basis -> thm
59val prove_at_left_0 : Lazy_Eval.eval_ctxt -> thm * basis -> thm
60val prove_at_right_0 : Lazy_Eval.eval_ctxt -> thm * basis -> thm
61
62val prove_smallo : Lazy_Eval.eval_ctxt -> thm * thm * basis -> thm
63val prove_bigo : Lazy_Eval.eval_ctxt -> thm * thm * basis -> thm
64val prove_bigtheta : Lazy_Eval.eval_ctxt -> thm * thm * basis -> thm
65val prove_asymp_equiv : Lazy_Eval.eval_ctxt -> thm * thm * basis -> thm
66
67val prove_asymptotic_relation : Lazy_Eval.eval_ctxt -> thm * thm * basis -> order * thm
68val prove_eventually_less : Lazy_Eval.eval_ctxt -> thm * thm * basis -> thm
69val prove_eventually_greater : Lazy_Eval.eval_ctxt -> thm * thm * basis -> thm
70val prove_eventually_nonzero : Lazy_Eval.eval_ctxt -> thm * basis -> thm
71
72val extract_terms : int * bool -> Lazy_Eval.eval_ctxt -> basis -> term -> term * term option
73
74(* Internal functions *)
75val check_expansion : Exp_Log_Expression.expr -> expansion_thm -> expansion_thm
76
77val zero_expansion : basis -> expansion_thm
78val const_expansion : Lazy_Eval.eval_ctxt -> basis -> term -> expansion_thm
79val ln_expansion :
80  Lazy_Eval.eval_ctxt -> trimmed_thm -> expansion_thm -> basis -> expansion_thm * basis
81val exp_expansion : Lazy_Eval.eval_ctxt -> expansion_thm -> basis -> expansion_thm * basis
82val powr_expansion :
83  Lazy_Eval.eval_ctxt -> expansion_thm * expansion_thm * basis -> expansion_thm * basis
84val powr_const_expansion :
85  Lazy_Eval.eval_ctxt -> expansion_thm * term * basis -> expansion_thm
86val powr_nat_expansion :
87  Lazy_Eval.eval_ctxt -> expansion_thm * expansion_thm * basis -> expansion_thm * basis
88val power_expansion : Lazy_Eval.eval_ctxt -> expansion_thm * term * basis -> expansion_thm
89val root_expansion : Lazy_Eval.eval_ctxt -> expansion_thm * term * basis -> expansion_thm
90
91val sgn_expansion : Lazy_Eval.eval_ctxt -> expansion_thm * basis -> expansion_thm
92val min_expansion : Lazy_Eval.eval_ctxt -> expansion_thm * expansion_thm * basis -> expansion_thm
93val max_expansion : Lazy_Eval.eval_ctxt -> expansion_thm * expansion_thm * basis -> expansion_thm
94val arctan_expansion : Lazy_Eval.eval_ctxt -> basis -> expansion_thm -> expansion_thm
95
96val ev_zeroness_oracle : Lazy_Eval.eval_ctxt -> term -> thm option
97val zeroness_oracle : bool -> trim_mode option -> Lazy_Eval.eval_ctxt -> term -> zeroness * thm option
98
99val whnf_expansion : Lazy_Eval.eval_ctxt -> expansion_thm -> term option * expansion_thm * thm
100val simplify_expansion : Lazy_Eval.eval_ctxt -> expansion_thm -> expansion_thm 
101val simplify_term : Lazy_Eval.eval_ctxt -> term -> term
102
103val trim_expansion_while_greater :
104  bool -> term list option -> bool -> trim_mode option -> Lazy_Eval.eval_ctxt ->
105    thm * Asymptotic_Basis.basis -> thm * trim_result * (zeroness * thm) list
106val trim_expansion : bool -> trim_mode option -> Lazy_Eval.eval_ctxt -> expansion_thm * basis -> 
107  expansion_thm * zeroness * trimmed_thm option
108val try_drop_leading_term_ex : bool -> Lazy_Eval.eval_ctxt -> expansion_thm -> expansion_thm option
109
110val try_prove_real_eq : bool -> Lazy_Eval.eval_ctxt -> term * term -> thm option
111val try_prove_ev_eq : Lazy_Eval.eval_ctxt -> term * term -> thm option
112val prove_compare_expansions : order -> thm list -> thm
113
114val simplify_trimmed_expansion : Lazy_Eval.eval_ctxt -> expansion_thm * trimmed_thm -> 
115  expansion_thm * trimmed_thm
116val retrim_expansion : Lazy_Eval.eval_ctxt -> expansion_thm * basis -> expansion_thm * thm
117val retrim_pos_expansion : Lazy_Eval.eval_ctxt -> expansion_thm * basis * trimmed_thm ->
118  expansion_thm * thm * trimmed_thm
119
120val register_sign_oracle : 
121  binding * (Proof.context -> int -> tactic) -> Context.generic -> Context.generic
122val get_sign_oracles :
123  Context.generic -> (string * (Proof.context -> int -> tactic)) list
124
125val solve_eval_eq : thm -> thm
126
127end
128
129structure Multiseries_Expansion : MULTISERIES_EXPANSION = struct
130
131open Asymptotic_Basis
132open Exp_Log_Expression
133open Lazy_Eval
134
135structure Data = Generic_Data
136(
137  type T = (Proof.context -> int -> tactic) Name_Space.table;
138  val empty : T = Name_Space.empty_table "sign oracle tactics";
139  val extend = I;
140  fun merge (tactics1, tactics2) : T = Name_Space.merge_tables (tactics1, tactics2);
141);
142
143fun register_sign_oracle (s, tac) ctxt =
144  Data.map (Name_Space.define ctxt false (s, tac) #> snd) ctxt
145
146fun get_sign_oracles ctxt = Name_Space.fold_table cons (Data.get ctxt) []
147
148fun apply_sign_oracles ctxt tac =
149  let
150    val oracles = get_sign_oracles (Context.Proof ctxt)
151    fun tac' {context = ctxt, concl, ...} =
152      if Thm.term_of concl = \<^term>\<open>HOL.Trueprop HOL.False\<close> then
153        no_tac
154      else
155        FIRST (map (fn tac => HEADGOAL (snd tac ctxt)) oracles)
156  in
157    tac THEN_ALL_NEW (Subgoal.FOCUS_PREMS tac' ctxt)
158  end
159    
160
161type expansion_thm = thm
162type trimmed_thm = thm
163
164val dest_fun = dest_comb #> fst
165val dest_arg = dest_comb #> snd
166val concl_of' = Thm.concl_of #> HOLogic.dest_Trueprop
167
168fun get_expansion thm =
169  thm |> Thm.concl_of |> HOLogic.dest_Trueprop |> Term.dest_comb |> fst |> Term.dest_comb |> snd
170
171fun get_expanded_fun thm = thm |> concl_of' |> dest_fun |> dest_fun |> dest_arg
172
173(*
174  The following function is useful in order to detect whether a given real constant is
175  an integer, which allows us to use the "f(x) ^ n" operation instead of "f(x) powr n".
176  This usually leads to nicer results.
177*)
178datatype intyness = Nat of thm | Neg_Nat of thm | No_Nat
179
180fun get_intyness ctxt ct =
181  if Thm.typ_of_cterm ct = \<^typ>\<open>Real.real\<close> then
182    let
183      val ctxt' = put_simpset HOL_basic_ss ctxt addsimps @{thms intyness_simps}
184      val conv = 
185        Simplifier.rewrite ctxt then_conv Simplifier.rewrite ctxt'
186      fun flip (Nat thm) = Neg_Nat (thm RS @{thm intyness_uminus})
187        | flip _ = No_Nat
188      fun get_intyness' ct =
189        case Thm.term_of ct of
190          \<^term>\<open>0::real\<close> => Nat @{thm intyness_0}
191        | \<^term>\<open>1::real\<close> => Nat @{thm intyness_1}
192        | Const (\<^const_name>\<open>numeral\<close>, _) $ _ => 
193            Nat (Thm.reflexive (Thm.dest_arg ct) RS @{thm intyness_numeral})
194        | Const (\<^const_name>\<open>uminus\<close>, _) $ _ => flip (get_intyness' (Thm.dest_arg ct))
195        | Const (\<^const_name>\<open>of_nat\<close>, _) $ _ => 
196            Nat (Thm.reflexive (Thm.dest_arg ct) RS @{thm intyness_of_nat})
197        | _ => No_Nat
198      val thm = conv ct
199      val ct' = thm |> Thm.cprop_of |> Thm.dest_equals_rhs
200    in
201      case get_intyness' ct' of
202        Nat thm' => Nat (Thm.transitive thm thm' RS @{thm HOL.meta_eq_to_obj_eq})
203      | Neg_Nat thm' => Neg_Nat (Thm.transitive thm thm' RS @{thm HOL.meta_eq_to_obj_eq})
204      | No_Nat => No_Nat
205    end
206      handle CTERM _ => No_Nat
207  else
208    No_Nat
209
210datatype parity = Even of thm | Odd of thm | Unknown_Parity
211
212(* TODO: powers *)
213fun get_parity ct =
214  let
215    fun inst thm cts =
216      let
217        val tvs = Term.add_tvars (Thm.concl_of thm) []
218      in
219        case tvs of
220          [v] =>
221            let
222              val thm' = Thm.instantiate ([(v, Thm.ctyp_of_cterm ct)], []) thm
223              val vs = take (length cts) (rev (Term.add_vars (Thm.concl_of thm') []))
224            in
225              Thm.instantiate ([], vs ~~ cts) thm'
226            end
227        | _ => raise THM ("get_parity", 0, [thm])
228      end
229    val get_num = Thm.dest_arg o Thm.dest_arg
230  in
231    case Thm.term_of ct of
232      Const (\<^const_name>\<open>Groups.zero\<close>, _) => Even (inst @{thm even_zero} [])
233    | Const (\<^const_name>\<open>Groups.one\<close>, _) => Odd (inst @{thm odd_one} [])
234    | Const (\<^const_name>\<open>Num.numeral_class.numeral\<close>, _) $ \<^term>\<open>Num.One\<close> =>
235        Odd (inst @{thm odd_Numeral1} [])
236    | Const (\<^const_name>\<open>Num.numeral_class.numeral\<close>, _) $ (\<^term>\<open>Num.Bit0\<close> $ _) =>
237        Even (inst @{thm even_numeral} [get_num ct])
238    | Const (\<^const_name>\<open>Num.numeral_class.numeral\<close>, _) $ (\<^term>\<open>Num.Bit1\<close> $ _) =>
239        Odd (inst @{thm odd_numeral} [get_num ct])
240    | Const (\<^const_name>\<open>Groups.uminus\<close>, _) $ _ => (
241        case get_parity (Thm.dest_arg ct) of
242          Even thm => Even (@{thm even_uminusI} OF [thm])
243        | Odd thm => Odd (@{thm odd_uminusI} OF [thm])
244        | _ => Unknown_Parity)
245    | Const (\<^const_name>\<open>Groups.plus\<close>, _) $ _ $ _ => (
246        case apply2 get_parity (Thm.dest_binop ct) of
247          (Even thm1, Even thm2) => Even (@{thm even_addI(1)} OF [thm1, thm2])
248        | (Odd thm1, Odd thm2) => Even (@{thm even_addI(2)} OF [thm1, thm2])
249        | (Even thm1, Odd thm2) => Odd (@{thm odd_addI(1)} OF [thm1, thm2])
250        | (Odd thm1, Even thm2) => Odd (@{thm odd_addI(2)} OF [thm1, thm2])
251        | _ => Unknown_Parity)
252    | Const (\<^const_name>\<open>Groups.minus\<close>, _) $ _ $ _ => (
253        case apply2 get_parity (Thm.dest_binop ct) of
254          (Even thm1, Even thm2) => Even (@{thm even_diffI(1)} OF [thm1, thm2])
255        | (Odd thm1, Odd thm2) => Even (@{thm even_diffI(2)} OF [thm1, thm2])
256        | (Even thm1, Odd thm2) => Odd (@{thm odd_diffI(1)} OF [thm1, thm2])
257        | (Odd thm1, Even thm2) => Odd (@{thm odd_diffI(2)} OF [thm1, thm2])
258        | _ => Unknown_Parity)
259    | Const (\<^const_name>\<open>Groups.times\<close>, _) $ _ $ _ => (
260        case apply2 get_parity (Thm.dest_binop ct) of
261          (Even thm1, _) => Even (@{thm even_multI(1)} OF [thm1])
262        | (_, Even thm2) => Even (@{thm even_multI(2)} OF [thm2])
263        | (Odd thm1, Odd thm2) => Odd (@{thm odd_multI} OF [thm1, thm2])
264        | _ => Unknown_Parity)
265    | Const (\<^const_name>\<open>Power.power\<close>, _) $ _ $ _ =>
266        let
267          val (a, n) = Thm.dest_binop ct
268        in
269          case get_parity a of
270            Odd thm => Odd (inst @{thm odd_powerI} [a, n] OF [thm])
271          | _ => Unknown_Parity
272        end
273    | _ => Unknown_Parity
274  end
275
276fun simplify_term' facts ctxt =
277  let
278    val ctxt = Simplifier.add_prems facts ctxt
279  in
280    Thm.cterm_of ctxt #> Simplifier.rewrite ctxt #> 
281    Thm.concl_of #> Logic.dest_equals #> snd
282  end
283
284fun simplify_term ectxt = simplify_term' (get_facts ectxt) (get_ctxt ectxt)
285
286fun simplify_eval ctxt =
287  simplify_term' [] (put_simpset HOL_basic_ss ctxt addsimps @{thms eval_simps})
288
289datatype zeroness = IsZero | IsNonZero | IsPos | IsNeg
290
291
292(* Caution: The following functions assume that the given expansion is in normal form already
293   as far as needed. *)
294
295(* Returns the leading coefficient of the given expansion. This coefficient is a multiseries. *)
296fun try_get_coeff expr =
297  case expr of
298    Const (\<^const_name>\<open>MS\<close>, _) $ (
299      Const (\<^const_name>\<open>MSLCons\<close>, _) $ (
300        Const (\<^const_name>\<open>Pair\<close>, _) $ c $ _) $ _) $ _ =>
301      SOME c
302  | _ => NONE
303
304fun get_coeff expr = 
305  expr |> dest_comb |> fst |> dest_comb |> snd |> dest_comb |> fst |> dest_comb |> snd
306    |> dest_comb |> fst |> dest_comb |> snd
307
308(* Returns the coefficient of the leading term in the expansion (i.e. a real number) *)
309fun get_lead_coeff expr =
310  case try_get_coeff expr of
311    NONE => expr
312  | SOME c => get_lead_coeff c
313
314(* Returns the exponent (w.r.t. the fastest-growing basis element) of the leading term *)
315fun get_exponent expr = 
316  expr |> dest_comb |> fst |> dest_comb |> snd |> dest_comb |> fst |> dest_comb |> snd
317    |> dest_comb |> snd
318
319(* Returns the list of exponents of the leading term *)
320fun get_exponents exp =
321  if fastype_of exp = \<^typ>\<open>real\<close> then
322    []
323  else
324    get_exponent exp :: get_exponents (get_coeff exp)
325
326(* Returns the function that the expansion corresponds to *)
327fun get_eval expr =
328  if fastype_of expr = \<^typ>\<open>real\<close> then
329    Abs ("x", \<^typ>\<open>real\<close>, expr)
330  else
331    expr |> dest_comb |> snd
332
333val eval_simps = @{thms eval_simps [THEN eq_reflection]}
334
335(* Tries to prove that the given function is eventually zero *)
336fun ev_zeroness_oracle ectxt t = 
337  let
338    val ctxt = Lazy_Eval.get_ctxt ectxt
339    val goal = 
340      betapply (\<^term>\<open>\<lambda>f::real \<Rightarrow> real. eventually (\<lambda>x. f x = 0) at_top\<close>, t)
341      |> HOLogic.mk_Trueprop
342    fun tac {context = ctxt, ...} =
343      HEADGOAL (Method.insert_tac ctxt (get_facts ectxt))
344      THEN Local_Defs.unfold_tac ctxt eval_simps
345      THEN HEADGOAL (Simplifier.asm_full_simp_tac ctxt)
346  in
347    try (Goal.prove ctxt [] [] goal) tac
348  end
349
350(* 
351  Encodes the kind of trimming/zeroness checking operation to be performed.
352  Simple_Trim only checks for zeroness/non-zeroness. Pos_Trim/Neg_Trim try to prove either
353  zeroness or positivity (resp. negativity). Sgn_Trim tries all three possibilities (positive,
354  negative, zero). *)
355datatype trim_mode = Simple_Trim | Pos_Trim | Neg_Trim | Sgn_Trim
356
357(*
358  Checks (and proves) whether the given term (assumed to be a real number) is zero, positive,
359  or negative, depending on given flags. The "fail" flag determines whether an exception is
360  thrown if this fails.
361*)
362fun zeroness_oracle fail mode ectxt exp = 
363  let
364    val ctxt = Lazy_Eval.get_ctxt ectxt
365    val eq = (exp, \<^term>\<open>0::real\<close>) |> HOLogic.mk_eq
366    val goal1 = (IsZero, eq |> HOLogic.mk_Trueprop)
367    val goal2 = 
368      case mode of
369        SOME Pos_Trim => 
370          (IsPos, \<^term>\<open>(<) (0::real)\<close> $ exp |> HOLogic.mk_Trueprop)
371      | SOME Sgn_Trim => 
372          (IsPos, \<^term>\<open>(<) (0::real)\<close> $ exp |> HOLogic.mk_Trueprop)
373      | SOME Neg_Trim => 
374          (IsNeg, betapply (\<^term>\<open>\<lambda>x. x < (0::real)\<close>, exp) |> HOLogic.mk_Trueprop)
375      | _ =>
376          (IsNonZero, eq |> HOLogic.mk_not |> HOLogic.mk_Trueprop)
377    val goals =
378      (if mode = SOME Sgn_Trim then 
379         [(IsNeg, betapply (\<^term>\<open>\<lambda>x. x < (0::real)\<close>, exp) |> HOLogic.mk_Trueprop)] 
380       else 
381         [])
382    val goals = goal2 :: goals
383    fun tac {context = ctxt, ...} =
384      HEADGOAL (Method.insert_tac ctxt (get_facts ectxt))
385      THEN Local_Defs.unfold_tac ctxt eval_simps
386      THEN HEADGOAL (apply_sign_oracles ctxt (Simplifier.asm_full_simp_tac ctxt))
387    fun prove (res, goal) = try (fn goal => (res, SOME (Goal.prove ctxt [] [] goal tac))) goal
388    fun err () =
389      let
390        val mode_msg =
391          case mode of
392            SOME Simple_Trim => "whether the following constant is zero"
393          | SOME Pos_Trim => "whether the following constant is zero or positive"
394          | SOME Neg_Trim => "whether the following constant is zero or negative"
395          | SOME Sgn_Trim => "the sign of the following constant"
396          | _ => raise Match
397        val t = simplify_term' (get_facts ectxt) ctxt exp
398        val _ =
399          if #verbose (#ctxt ectxt) then
400            let
401              val p = Pretty.str ("real_asymp failed to determine " ^ mode_msg ^ ":")
402              val p = Pretty.chunks [p, Pretty.indent 2 (Syntax.pretty_term ctxt t)]
403            in
404              Pretty.writeln p
405            end else ()
406      in
407        raise TERM ("zeroness_oracle", [t])
408      end
409  in
410    case prove goal1 of
411      SOME res => res
412    | NONE => 
413        if mode = NONE then 
414          (IsNonZero, NONE)
415        else
416          case get_first prove (goal2 :: goals) of
417            NONE => if fail then err () else (IsNonZero, NONE)
418          | SOME res => res
419  end
420
421(* Tries to prove a given equality of real numbers. *)
422fun try_prove_real_eq fail ectxt (lhs, rhs) =
423  case zeroness_oracle false NONE ectxt (\<^term>\<open>(-) :: real => _\<close> $ lhs $ rhs) of
424    (IsZero, SOME thm) => SOME (thm RS @{thm real_eqI})
425  | _ => 
426    if not fail then NONE else
427      let
428        val ctxt = get_ctxt ectxt
429        val ts = map (simplify_term' (get_facts ectxt) ctxt) [lhs, rhs]
430        val _ =
431          if #verbose (#ctxt ectxt) then
432            let
433              val p = 
434                Pretty.str ("real_asymp failed to prove that the following two numbers are equal:")
435              val p = Pretty.chunks (p :: map (Pretty.indent 2 o Syntax.pretty_term ctxt) ts)
436            in
437              Pretty.writeln p
438            end else ()
439      in
440        raise TERM ("try_prove_real_eq", [lhs, rhs])
441      end
442
443(* Tries to prove a given eventual equality of real functions. *)
444fun try_prove_ev_eq ectxt (f, g) =
445  let
446    val t = Envir.beta_eta_contract (\<^term>\<open>\<lambda>(f::real=>real) g x. f x - g x\<close> $ f $ g)
447  in
448    Option.map (fn thm => thm RS @{thm eventually_diff_zero_imp_eq}) (ev_zeroness_oracle ectxt t)
449  end
450
451fun real_less a b = \<^term>\<open>(<) :: real \<Rightarrow> real \<Rightarrow> bool\<close> $ a $ b
452fun real_eq a b = \<^term>\<open>(=) :: real \<Rightarrow> real \<Rightarrow> bool\<close> $ a $ b
453fun real_neq a b = \<^term>\<open>(\<noteq>) :: real \<Rightarrow> real \<Rightarrow> bool\<close> $ a $ b
454
455(* The hook that is called by the Lazy_Eval module whenever two real numbers have to be compared *)
456fun real_sgn_hook ({pctxt = ctxt, facts, verbose, ...}) t =
457  let
458    val get_rhs = Thm.concl_of #> Logic.dest_equals #> snd
459    fun tac {context = ctxt, ...} = 
460      HEADGOAL (Method.insert_tac ctxt (Net.content facts) 
461        THEN' (apply_sign_oracles ctxt (Simplifier.asm_full_simp_tac ctxt)))
462    fun prove_first err [] [] =
463          if not verbose then raise TERM ("real_sgn_hook", [t])
464            else let val _ = err () in raise TERM ("real_sgn_hook", [t]) end
465      | prove_first err (goal :: goals) (thm :: thms) =
466          (case try (Goal.prove ctxt [] [] goal) tac of
467             SOME thm' => 
468               let val thm'' = thm' RS thm in SOME (get_rhs thm'', Conv.rewr_conv thm'') end
469           | NONE => prove_first err goals thms)
470      | prove_first _ _ _ = raise Match
471  in
472    case t of
473      \<^term>\<open>(=) :: real => _\<close> $ a $ \<^term>\<open>0 :: real\<close> =>
474        let
475          val goals =
476            map (fn c => HOLogic.mk_Trueprop (c a \<^term>\<open>0 :: real\<close>)) [real_neq, real_eq]
477          fun err () = 
478            let
479              val facts' = Net.content facts
480              val a' = simplify_term' facts' ctxt a
481              val p = Pretty.str ("real_asymp failed to determine whether the following " ^
482                                    "constant is zero: ")
483              val p = Pretty.chunks [p, Pretty.indent 2 (Syntax.pretty_term ctxt a')]
484            in
485              Pretty.writeln p
486            end
487        in
488          prove_first err goals @{thms Eq_FalseI Eq_TrueI}
489        end
490    | Const (\<^const_name>\<open>COMPARE\<close>, _) $ a $ b =>
491        let
492          val goals = map HOLogic.mk_Trueprop [real_less a b, real_less b a, real_eq a b]
493          fun err () = 
494            let
495              val facts' = Net.content facts
496              val (a', b') = apply2 (simplify_term' facts' ctxt) (a, b)
497              val p = Pretty.str ("real_asymp failed to compare" ^
498                        "the following two constants: ")
499              val p = Pretty.chunks (p :: map (Pretty.indent 2 o Syntax.pretty_term ctxt) [a', b'])
500            in
501              Pretty.writeln p
502            end
503        in
504          prove_first err goals @{thms COMPARE_intros}
505        end
506    | _ => NONE
507  end
508
509(* 
510  Returns the datatype constructors registered for use with the Lazy_Eval package.
511  All constructors on which pattern matching is performed need to be registered for evaluation
512  to work. It should be rare for users to add additional ones.
513*)
514fun get_constructors ctxt =
515  let
516    val thms = Named_Theorems.get ctxt \<^named_theorems>\<open>exp_log_eval_constructor\<close>
517    fun go _ [] acc = rev acc
518      | go f (x :: xs) acc =
519          case f x of
520            NONE => go f xs acc
521          | SOME y => go f xs (y :: acc)
522    fun map_option f xs = go f xs []
523    fun dest_constructor thm =
524      case Thm.concl_of thm of
525        Const (\<^const_name>\<open>HOL.Trueprop\<close>, _) $
526            (Const (\<^const_name>\<open>REAL_ASYMP_EVAL_CONSTRUCTOR\<close>, _) $ Const (c, T)) =>
527          SOME (c, length (fst (strip_type T)))
528     | _ => NONE
529  in
530    thms |> map_option dest_constructor
531  end
532
533(*
534  Creates an evaluation context with the correct setup of constructors,  equations, and hooks.
535*)
536fun mk_eval_ctxt ctxt =
537  let
538    val eval_eqs = (Named_Theorems.get ctxt \<^named_theorems>\<open>real_asymp_eval_eqs\<close>)
539    val constructors = get_constructors ctxt
540  in
541    Lazy_Eval.mk_eval_ctxt ctxt constructors eval_eqs
542    |> add_hook real_sgn_hook
543  end
544
545(* A pattern for determining the leading coefficient of a multiseries *)
546val exp_pat = 
547  let
548    val anypat = AnyPat ("_", 0)
549  in
550    ConsPat (\<^const_name>\<open>MS\<close>, 
551      [ConsPat (\<^const_name>\<open>MSLCons\<close>, 
552         [ConsPat (\<^const_name>\<open>Pair\<close>, [anypat, anypat]), anypat]), anypat])
553  end
554
555(*
556  Evaluates an expansion to (weak) head normal form, so that the leading coefficient and
557  exponent can be read off.
558*)
559fun whnf_expansion ectxt thm =
560  let
561    val ctxt = get_ctxt ectxt
562    val exp = get_expansion thm
563    val (_, _, conv) = match ectxt exp_pat exp (SOME [])
564    val eq_thm = conv (Thm.cterm_of ctxt exp)
565    val exp' = eq_thm |> Thm.concl_of |> Logic.dest_equals |> snd
566  in
567    case exp' of
568      Const (\<^const_name>\<open>MS\<close>, _) $ (Const (\<^const_name>\<open>MSLCons\<close>, _) $ 
569        (Const (\<^const_name>\<open>Pair\<close>, _) $ c $ _) $ _) $ _ =>
570          (SOME c, @{thm expands_to_meta_eq_cong} OF [thm, eq_thm], eq_thm)
571    | Const (\<^const_name>\<open>MS\<close>, _) $ Const (\<^const_name>\<open>MSLNil\<close>, _) $ _ =>
572        (NONE, @{thm expands_to_meta_eq_cong} OF [thm, eq_thm], eq_thm)
573    | _ => raise TERM ("whnf_expansion", [exp'])
574  end
575
576fun try_lift_function ectxt (thm, SEmpty) _ = (NONE, thm)
577  | try_lift_function ectxt (thm, basis) cont =
578  case whnf_expansion ectxt thm of
579    (SOME c, thm, _) =>
580      let
581        val f = get_expanded_fun thm
582        val T = fastype_of c
583        val t = Const (\<^const_name>\<open>eval\<close>, T --> \<^typ>\<open>real \<Rightarrow> real\<close>) $ c
584        val t = Term.betapply (Term.betapply (\<^term>\<open>\<lambda>(f::real\<Rightarrow>real) g x. f x - g x\<close>, f), t)
585      in
586        case ev_zeroness_oracle ectxt t of
587          NONE => (NONE, thm)
588        | SOME zero_thm =>
589            let
590              val thm' = cont ectxt (thm RS @{thm expands_to_hd''}, tl_basis basis)
591              val thm'' = @{thm expands_to_lift_function} OF [zero_thm, thm']
592            in
593              (SOME (lift basis thm''), thm)
594            end
595      end
596  | _ => (NONE, thm)
597
598(* Turns an expansion theorem into an expansion theorem for the leading coefficient. *)
599fun expands_to_hd thm = thm RS
600  (if fastype_of (get_expansion thm) = \<^typ>\<open>real ms\<close> then 
601     @{thm expands_to_hd'}
602   else 
603     @{thm expands_to_hd})
604
605fun simplify_expansion ectxt thm =
606  let
607    val exp = get_expansion thm
608    val ctxt = get_ctxt ectxt
609    val eq_thm = Simplifier.rewrite ctxt (Thm.cterm_of ctxt exp)
610  in
611    @{thm expands_to_meta_eq_cong} OF [thm, eq_thm]
612  end
613
614(*
615  Simplifies a trimmed expansion and returns the simplified expansion theorem and
616  the trimming theorem for that simplified expansion.
617*)
618fun simplify_trimmed_expansion ectxt (thm, trimmed_thm) =
619  let
620    val exp = get_expansion thm
621    val ctxt = get_ctxt ectxt
622    val eq_thm = Simplifier.rewrite ctxt (Thm.cterm_of ctxt exp)
623    val trimmed_cong_thm =
624      case trimmed_thm |> concl_of' |> dest_fun of
625        Const (\<^const_name>\<open>trimmed\<close>, _) => @{thm trimmed_eq_cong}
626      | Const (\<^const_name>\<open>trimmed_pos\<close>, _) => @{thm trimmed_pos_eq_cong}
627      | Const (\<^const_name>\<open>trimmed_neg\<close>, _) => @{thm trimmed_neg_eq_cong}
628      | _ => raise THM ("simplify_trimmed_expansion", 2, [thm, trimmed_thm])
629  in
630    (@{thm expands_to_meta_eq_cong} OF [thm, eq_thm], 
631      trimmed_cong_thm OF [trimmed_thm, eq_thm])     
632  end
633
634(*
635  Re-normalises a trimmed expansion (so that the leading term with its (real) coefficient and
636  all exponents can be read off. This may be necessary after lifting a trimmed expansion to
637  a larger basis.
638*)
639fun retrim_expansion ectxt (thm, basis) =
640  let
641    val (c, thm, eq_thm) = whnf_expansion ectxt thm
642  in
643    case c of
644      NONE => (thm, eq_thm)
645    | SOME c =>
646      if fastype_of c = \<^typ>\<open>real\<close> then 
647        (thm, eq_thm) 
648      else
649        let
650          val c_thm = thm RS @{thm expands_to_hd''}
651          val (c_thm', eq_thm') = retrim_expansion ectxt (c_thm, tl_basis basis)
652          val thm = @{thm expands_to_trim_cong} OF [thm, c_thm']
653        in
654          (thm, @{thm trim_lift_eq} OF [eq_thm, eq_thm'])
655        end
656  end
657
658fun retrim_pos_expansion ectxt (thm, basis, trimmed_thm) =
659  let
660    val (thm', eq_thm) = retrim_expansion ectxt (thm, basis)
661  in
662    (thm', eq_thm, @{thm trimmed_pos_eq_cong} OF [trimmed_thm, eq_thm])
663  end
664
665(*
666  Tries to determine whether the leading term is (identically) zero and drops it if it is.
667  If "fail" is set, an exception is thrown when that term is a real number and zeroness cannot
668  be determined. (Which typically indicates missing facts or case distinctions)
669*)
670fun try_drop_leading_term_ex fail ectxt thm =
671  let
672    val exp = get_expansion thm
673  in
674    if fastype_of exp = \<^typ>\<open>real\<close> then
675      NONE
676    else if fastype_of (get_coeff exp) = \<^typ>\<open>real\<close> then
677      case zeroness_oracle fail (SOME Simple_Trim) ectxt (get_coeff exp) of
678        (IsZero, SOME zero_thm) => SOME (@{thm drop_zero_ms'} OF [zero_thm, thm])
679      | _ => NONE
680    else
681      let
682        val c = get_coeff exp
683        val T = fastype_of c
684        val t = Const (\<^const_name>\<open>eval\<close>, T --> \<^typ>\<open>real \<Rightarrow> real\<close>) $ c
685      in
686        case ev_zeroness_oracle ectxt t of
687          SOME zero_thm => SOME (@{thm expands_to_drop_zero} OF [zero_thm, thm])
688        | _ => NONE
689      end
690  end
691
692(*
693  Tries to drop the leading term of an expansion. If this is not possible, an exception 
694  is thrown and an informative error message is printed.
695*)
696fun try_drop_leading_term ectxt thm =
697  let
698    fun err () =
699      let
700        val ctxt = get_ctxt ectxt
701        val exp = get_expansion thm
702        val c = get_coeff exp
703        val t = 
704          if fastype_of c = \<^typ>\<open>real\<close> then c else c |> dest_arg
705        val t = simplify_term' (get_facts ectxt) ctxt t
706        val _ =
707          if #verbose (#ctxt ectxt) then
708            let
709              val p = Pretty.str ("real_asymp failed to prove that the following term is zero: ")
710              val p = Pretty.chunks [p, Pretty.indent 2 (Syntax.pretty_term ctxt t)]
711            in
712              Pretty.writeln p
713            end else ()
714      in
715        raise TERM ("try_drop_leading_term", [t])
716      end
717  in
718    case try_drop_leading_term_ex true ectxt thm of
719      NONE => err ()
720    | SOME thm => thm
721  end
722
723
724datatype trim_result =
725    Trimmed of zeroness * trimmed_thm option
726  | Aborted of order
727
728fun cstrip_assms ct =
729  case Thm.term_of ct of
730    \<^term>\<open>(==>)\<close> $ _ $ _ => cstrip_assms (snd (Thm.dest_implies ct))
731  | _ => ct
732
733(*
734  Trims an expansion (i.e. drops leading zero terms) and provides a trimmedness theorem.
735  Optionally, a list of exponents can be given to instruct the function to only trim until
736  the exponents of the leading term are lexicographically less than (or less than or equal) than
737  the given ones. This is useful to avoid unnecessary trimming.
738
739  The "strict" flag indicates whether the trimming should already be aborted when the 
740  exponents are lexicographically equal or not.
741
742  The "fail" flag is passed on to the zeroness oracle and determines whether a failure to determine
743  the sign of a real number leads to an exception.
744
745  "mode" indicates what kind of trimmedness theorem will be returned: Simple_Trim only gives the
746  default trimmedness theorem, whereas Pos_Trim/Neg_Trim/Sgn_Trim will give trimmed_pos or
747  trimmed_neg. Giving "None" as mode will produce no trimmedness theorem; it will only drop 
748  leading zero terms until zeroness cannot be proven anymore, upon which it will stop.
749
750  The main result of the function is the trimmed expansion theorem.
751
752  The function returns whether the trimming has been aborted or not. If was aborted, either
753  LESS or EQUAL will be returned, indicating whether the exponents of the leading term are
754  now lexicographically smaller or equal to the given ones. In the other case, the zeroness
755  of the leading coefficient is returned (zero, non-zero, positive, negative) together with a
756  trimmedness theorem.
757
758  Lastly, a list of the exponent comparison results and associated theorems is also returned, so
759  that the caller can reconstruct the result of the lexicographic ordering without doing the
760  exponent comparisons again.
761*)
762fun trim_expansion_while_greater strict es fail mode ectxt (thm, basis) = 
763  let
764    val (_, thm, _) = whnf_expansion ectxt thm
765    val thm = simplify_expansion ectxt thm
766    val cexp = thm |> Thm.cprop_of |> cstrip_assms |> Thm.dest_arg |> Thm.dest_fun |> Thm.dest_arg
767    val c = try_get_coeff (get_expansion thm)
768    fun lift_trimmed_thm nz thm =
769      let
770        val cexp = thm |> Thm.cprop_of |> cstrip_assms |> Thm.dest_arg |> Thm.dest_fun |> Thm.dest_arg
771        val lift_thm =
772          case nz of
773            IsNonZero => @{thm trimmed_eq_cong[rotated, OF _ lift_trimmed]}
774          | IsPos => @{thm trimmed_pos_eq_cong[rotated, OF _ lift_trimmed_pos]}
775          | IsNeg => @{thm trimmed_neg_eq_cong[rotated, OF _ lift_trimmed_neg]}
776          | _ => raise TERM ("Unexpected zeroness result in trim_expansion", [])
777      in 
778        Thm.reflexive cexp RS lift_thm
779      end        
780    fun trimmed_real_thm nz = Thm.reflexive cexp RS (
781      case nz of
782        IsNonZero => @{thm trimmed_eq_cong[rotated, OF _ lift_trimmed[OF trimmed_realI]]}
783      | IsPos => @{thm trimmed_pos_eq_cong[rotated, OF _ lift_trimmed_pos[OF trimmed_pos_realI]]}
784      | IsNeg => @{thm trimmed_neg_eq_cong[rotated, OF _ lift_trimmed_neg[OF trimmed_neg_realI]]}
785      | _ => raise TERM ("Unexpected zeroness result in trim_expansion", []))
786    fun do_trim es =
787      let
788        val c = the c
789        val T = fastype_of c
790        val t = Const (\<^const_name>\<open>eval\<close>, T --> \<^typ>\<open>real \<Rightarrow> real\<close>) $ c
791      in
792        if T = \<^typ>\<open>real\<close> then (
793          case zeroness_oracle fail mode ectxt c of
794            (IsZero, SOME zero_thm) => 
795              trim_expansion_while_greater strict es fail mode ectxt
796                (@{thm drop_zero_ms'} OF [zero_thm, thm], basis)
797          | (nz, SOME nz_thm) => (thm, Trimmed (nz, SOME (nz_thm RS trimmed_real_thm nz)), [])
798          | (nz, NONE) => (thm, Trimmed (nz, NONE), []))
799        else
800          case trim_expansion_while_greater strict (Option.map tl es) fail mode ectxt
801                 (thm RS @{thm expands_to_hd''}, tl_basis basis) of
802            (c_thm', Aborted ord, thms) =>
803              (@{thm expands_to_trim_cong} OF [thm, c_thm'], Aborted ord, thms)
804          | (c_thm', Trimmed (nz, trimmed_thm), thms) =>
805              let
806                val thm = (@{thm expands_to_trim_cong} OF [thm, c_thm'])
807                fun err () =
808                  raise TERM ("trim_expansion: zero coefficient should have been trimmed", [c])
809              in
810                case (nz, trimmed_thm) of
811                  (IsZero, _) => 
812                    if #verbose (#ctxt ectxt) then
813                      let
814                        val ctxt = get_ctxt ectxt
815                        val t' = t |> simplify_eval ctxt |> simplify_term' (get_facts ectxt) ctxt
816                        val p = Pretty.str ("trim_expansion failed to recognise zeroness of " ^
817                          "the following term:")
818                        val p = Pretty.chunks [p, Pretty.indent 2 (Syntax.pretty_term ctxt t')]
819                        val _ = Pretty.writeln p
820                      in
821                        err ()
822                      end
823                    else err ()
824                | (_, SOME trimmed_thm) =>
825                      (thm, Trimmed (nz, SOME (trimmed_thm RS lift_trimmed_thm nz thm)), thms)
826                | (_, NONE) => (thm, Trimmed (nz, NONE), thms)
827              end
828      end
829    val minus = \<^term>\<open>(-) :: real => real => real\<close>
830  in
831    case (c, es) of
832      (NONE, _) => (thm, Trimmed (IsZero, NONE), [])
833    | (SOME c, SOME (e' :: _)) =>
834        let
835          val e = get_exponent (get_expansion thm)
836        in
837          case zeroness_oracle true (SOME Sgn_Trim) ectxt (minus $ e $ e') of
838            (IsPos, SOME pos_thm) => (
839              case try_drop_leading_term_ex false ectxt thm of
840                SOME thm =>
841                  trim_expansion_while_greater strict es fail mode ectxt (thm, basis)
842              | NONE => do_trim NONE |> @{apply 3(3)} (fn thms => (IsPos, pos_thm) :: thms))
843          | (IsNeg, SOME neg_thm) => (thm, Aborted LESS, [(IsNeg, neg_thm)])
844          | (IsZero, SOME zero_thm) =>
845              if not strict andalso fastype_of c = \<^typ>\<open>real\<close> then
846                (thm, Aborted EQUAL, [(IsZero, zero_thm)])
847              else (
848                case try_drop_leading_term_ex false ectxt thm of
849                  SOME thm => trim_expansion_while_greater strict es fail mode ectxt (thm, basis)
850                | NONE => (do_trim es |> @{apply 3(3)} (fn thms => (IsZero, zero_thm) :: thms)))
851          | _ => do_trim NONE
852        end
853    | _ => (
854      case try_drop_leading_term_ex false ectxt thm of
855          SOME thm => trim_expansion_while_greater strict es fail mode ectxt (thm, basis)
856        | NONE => do_trim NONE)
857  end
858
859(*
860  Trims an expansion without any stopping criterion.
861*)
862fun trim_expansion fail mode ectxt (thm, basis) = 
863  case trim_expansion_while_greater false NONE fail mode ectxt (thm, basis) of
864    (thm, Trimmed (zeroness, trimmed_thm), _) => (thm, zeroness, trimmed_thm)
865  | _ => raise Match
866
867(*
868  Determines the sign of an expansion that has already been trimmed.
869*)
870fun determine_trimmed_sgn ectxt exp =
871  if fastype_of exp = \<^typ>\<open>real\<close> then
872    (case zeroness_oracle true (SOME Sgn_Trim) ectxt exp of
873       (IsPos, SOME thm) => (IsPos, thm RS @{thm trimmed_pos_realI})
874     | (IsNeg, SOME thm) => (IsNeg, thm RS @{thm trimmed_neg_realI})
875     | _ => raise TERM ("determine_trimmed_sgn", []))
876  else
877    let
878      val ct = Thm.cterm_of (get_ctxt ectxt) exp
879    in
880      (case determine_trimmed_sgn ectxt (get_coeff exp) of
881         (IsPos, thm) => (IsPos, @{thm lift_trimmed_pos'} OF [thm, Thm.reflexive ct])
882       | (IsNeg, thm) => (IsNeg, @{thm lift_trimmed_neg'} OF [thm, Thm.reflexive ct])
883       | _ => raise TERM ("determine_trimmed_sgn", []))
884    end
885
886fun mk_compare_expansions_const T =
887      Const (\<^const_name>\<open>compare_expansions\<close>, 
888        T --> T --> \<^typ>\<open>cmp_result \<times> real \<times> real\<close>)
889
890datatype comparison_result =
891  Cmp_Dominated of order * thm list * zeroness * trimmed_thm * expansion_thm * expansion_thm 
892| Cmp_Asymp_Equiv of thm * thm
893
894fun compare_expansions' _ (thm1, thm2, SEmpty) = Cmp_Asymp_Equiv (thm1, thm2)
895  | compare_expansions' ectxt (thm1, thm2, basis) =
896  let
897    fun lift_trimmed_thm nz =
898      case nz of
899        IsPos => @{thm lift_trimmed_pos}
900      | IsNeg => @{thm lift_trimmed_neg}
901      | _ => raise TERM ("Unexpected zeroness result in compare_expansions'", [])
902    val (e1, e2) = apply2 (get_expansion #> get_exponent) (thm1, thm2)
903    val e = \<^term>\<open>(-) :: real => _\<close> $ e1 $ e2
904    fun trim thm = trim_expansion true (SOME Sgn_Trim) ectxt (thm, basis)
905    val try_drop = Option.map (whnf_expansion ectxt #> #2) o try_drop_leading_term_ex false ectxt
906    fun handle_result ord zeroness trimmed_thm thm1 thm2 =
907      let
908        val (e1, e2) = apply2 (get_expansion #> get_exponent) (thm1, thm2)
909        val e = \<^term>\<open>(-) :: real => _\<close> $ e1 $ e2
910        val mode = if ord = LESS then Neg_Trim else Pos_Trim
911      in
912        case zeroness_oracle true (SOME mode) ectxt e of
913          (_, SOME e_thm) => Cmp_Dominated (ord, [e_thm], zeroness, trimmed_thm, thm1, thm2)
914        | _ => raise Match
915      end
916    fun recurse e_zero_thm =
917      case basis of
918        SNE (SSng _) => Cmp_Asymp_Equiv (thm1, thm2)
919      | _ =>
920        let
921          val (thm1', thm2') = apply2 (fn thm => thm RS @{thm expands_to_hd''}) (thm1, thm2)
922          val (thm1', thm2') = apply2 (whnf_expansion ectxt #> #2) (thm1', thm2')
923        in
924          case compare_expansions' ectxt (thm1', thm2', tl_basis basis) of
925            Cmp_Dominated (order, e_thms, zeroness, trimmed_thm, thm1', thm2') =>
926              Cmp_Dominated (order, e_zero_thm :: e_thms, zeroness,
927                trimmed_thm RS lift_trimmed_thm zeroness,
928                @{thm expands_to_trim_cong} OF [thm1, thm1'],
929                @{thm expands_to_trim_cong} OF [thm2, thm2'])
930          | Cmp_Asymp_Equiv (thm1', thm2') => Cmp_Asymp_Equiv
931              (@{thm expands_to_trim_cong} OF [thm1, thm1'],
932                @{thm expands_to_trim_cong} OF [thm2, thm2'])
933        end
934  in
935    case zeroness_oracle false (SOME Sgn_Trim) ectxt e of
936      (IsPos, SOME _) => (
937        case try_drop thm1 of
938          SOME thm1 => compare_expansions' ectxt (thm1, thm2, basis)
939        | NONE => (
940            case trim thm1 of
941              (thm1, zeroness, SOME trimmed_thm) =>
942                handle_result GREATER zeroness trimmed_thm thm1 thm2
943            | _ => raise TERM ("compare_expansions", map get_expansion [thm1, thm2])))
944    | (IsNeg, SOME _) => (
945        case try_drop thm2 of
946          SOME thm2 => compare_expansions' ectxt (thm1, thm2, basis)
947        | NONE => (
948            case trim thm2 of
949              (thm2, zeroness, SOME trimmed_thm) =>
950                handle_result LESS zeroness trimmed_thm thm1 thm2
951            | _ => raise TERM ("compare_expansions", map get_expansion [thm1, thm2])))
952    | (IsZero, SOME e_zero_thm) => (
953        case try_drop thm1 of
954          SOME thm1 => compare_expansions' ectxt (thm1, thm2, basis)
955        | NONE => (
956            case try_drop thm2 of
957              SOME thm2 => compare_expansions' ectxt (thm1, thm2, basis)
958            | NONE => recurse e_zero_thm))
959    | _ =>
960        case try_drop thm1 of
961          SOME thm1 => compare_expansions' ectxt (thm1, thm2, basis)
962        | NONE => (
963            case try_drop thm2 of
964              SOME thm2 => compare_expansions' ectxt (thm1, thm2, basis)
965            | NONE => raise TERM ("compare_expansions", [e1, e2]))
966  end
967
968(* Uses a list of exponent comparison results to show that compare_expansions has a given result.*)
969fun prove_compare_expansions ord [thm] = (
970      case ord of
971        LESS => @{thm compare_expansions_LT_I} OF [thm]
972      | GREATER => @{thm compare_expansions_GT_I} OF [thm]
973      | EQUAL => @{thm compare_expansions_same_exp[OF _ compare_expansions_real]} OF [thm])
974  | prove_compare_expansions ord (thm :: thms) =
975      @{thm compare_expansions_same_exp} OF [thm, prove_compare_expansions ord thms]
976  | prove_compare_expansions _ [] = raise Match
977
978val ev_zero_pos_thm = Eventuallize.eventuallize \<^context>
979  @{lemma "\<forall>x::real. f x = 0 \<longrightarrow> g x > 0 \<longrightarrow> f x < g x" by auto} NONE
980  OF @{thms _ expands_to_imp_eventually_pos}
981
982val ev_zero_neg_thm = Eventuallize.eventuallize \<^context>
983  @{lemma "\<forall>x::real. f x = 0 \<longrightarrow> g x < 0 \<longrightarrow> f x > g x" by auto} NONE
984  OF @{thms _ expands_to_imp_eventually_neg}
985
986val ev_zero_zero_thm = Eventuallize.eventuallize \<^context>
987  @{lemma "\<forall>x::real. f x = 0 \<longrightarrow> g x = 0 \<longrightarrow> f x = g x" by auto} NONE
988
989fun compare_expansions_trivial ectxt (thm1, thm2, basis) =
990  case try_prove_ev_eq ectxt (apply2 get_expanded_fun (thm1, thm2)) of
991    SOME thm => SOME (EQUAL, thm, thm1, thm2)
992  | NONE =>
993      case apply2 (ev_zeroness_oracle ectxt o get_expanded_fun) (thm1, thm2) of
994        (NONE, NONE) => NONE
995      | (SOME zero1_thm, NONE) => (
996          case trim_expansion true (SOME Sgn_Trim) ectxt (thm2, basis) of
997            (thm2, IsPos, SOME trimmed2_thm) =>
998              SOME (LESS, ev_zero_pos_thm OF
999                [zero1_thm, get_basis_wf_thm basis, thm2, trimmed2_thm], thm1, thm2)
1000          | (thm2, IsNeg, SOME trimmed2_thm) =>
1001              SOME (GREATER, ev_zero_neg_thm OF
1002                [zero1_thm, get_basis_wf_thm basis, thm2, trimmed2_thm], thm1, thm2)
1003          | _ => raise TERM ("Unexpected zeroness result in compare_expansions", []))
1004      | (NONE, SOME zero2_thm) => (
1005          case trim_expansion true (SOME Sgn_Trim) ectxt (thm1, basis) of
1006            (thm1, IsPos, SOME trimmed1_thm) =>
1007              SOME (GREATER, ev_zero_pos_thm OF
1008                [zero2_thm, get_basis_wf_thm basis, thm1, trimmed1_thm], thm1, thm2)
1009          | (thm1, IsNeg, SOME trimmed1_thm) =>
1010              SOME (LESS, ev_zero_neg_thm OF
1011                [zero2_thm, get_basis_wf_thm basis, thm1, trimmed1_thm], thm1, thm2)
1012          | _ => raise TERM ("Unexpected zeroness result in compare_expansions", []))
1013      | (SOME zero1_thm, SOME zero2_thm) =>
1014          SOME (EQUAL, ev_zero_zero_thm OF [zero1_thm, zero2_thm] , thm1, thm2)
1015
1016fun compare_expansions ectxt (thm1, thm2, basis) =
1017  case compare_expansions_trivial ectxt (thm1, thm2, basis) of
1018    SOME res => res
1019  | NONE =>
1020    let
1021      val (_, thm1, _) = whnf_expansion ectxt thm1
1022      val (_, thm2, _) = whnf_expansion ectxt thm2
1023    in
1024      case compare_expansions' ectxt (thm1, thm2, basis) of
1025        Cmp_Dominated (order, e_thms, zeroness, trimmed_thm, thm1, thm2) =>
1026          let
1027            val wf_thm = get_basis_wf_thm basis
1028            val cmp_thm = prove_compare_expansions order e_thms
1029            val trimmed_thm' = trimmed_thm RS
1030              (if zeroness = IsPos then @{thm trimmed_pos_imp_trimmed}
1031                 else @{thm trimmed_neg_imp_trimmed})
1032            val smallo_thm = 
1033              (if order = LESS then @{thm compare_expansions_LT} else @{thm compare_expansions_GT}) OF
1034                [cmp_thm, trimmed_thm', thm1, thm2, wf_thm]
1035            val thm' = 
1036              if zeroness = IsPos then @{thm smallo_trimmed_imp_eventually_less} 
1037              else @{thm smallo_trimmed_imp_eventually_greater}
1038            val result_thm =
1039              thm' OF [smallo_thm, if order = LESS then thm2 else thm1, wf_thm, trimmed_thm]
1040          in
1041            (order, result_thm, thm1, thm2)
1042          end
1043       | Cmp_Asymp_Equiv (thm1, thm2) =>
1044          let
1045            val thm = @{thm expands_to_minus} OF [get_basis_wf_thm basis, thm1, thm2]
1046            val (order, result_thm) =
1047              case trim_expansion true (SOME Sgn_Trim) ectxt (thm, basis) of
1048                (thm, IsPos, SOME pos_thm) => (GREATER,
1049                  @{thm expands_to_imp_eventually_gt} OF [get_basis_wf_thm basis, thm, pos_thm])
1050              | (thm, IsNeg, SOME neg_thm) => (LESS,
1051                  @{thm expands_to_imp_eventually_lt} OF [get_basis_wf_thm basis, thm, neg_thm])
1052              | _ => raise TERM ("Unexpected zeroness result in prove_eventually_less", [])
1053          in
1054            (order, result_thm, thm1, thm2)
1055          end
1056    end
1057
1058
1059
1060(*
1061  Throws an exception and prints an error message indicating that the leading term could 
1062  not be determined to be either zero or non-zero.
1063*)
1064fun raise_trimming_error ectxt thm =
1065  let
1066    val ctxt = get_ctxt ectxt
1067    fun lead_coeff exp =
1068      if fastype_of exp = \<^typ>\<open>real\<close> then exp else lead_coeff (get_coeff exp)
1069    val c = lead_coeff (get_expansion thm)
1070    fun err () =
1071      let
1072        val t = simplify_term' (get_facts ectxt) ctxt c
1073        val _ =
1074          if #verbose (#ctxt ectxt) then
1075            let
1076              val p = Pretty.str 
1077                ("real_asymp failed to determine whether the following constant is zero:")
1078              val p = Pretty.chunks [p, Pretty.indent 2 (Syntax.pretty_term ctxt t)]
1079            in
1080              Pretty.writeln p
1081            end else ()
1082      in
1083        raise TERM ("zeroness_oracle", [t])
1084      end
1085  in
1086    err ()
1087  end
1088    
1089
1090(* TODO Here be dragons *)
1091fun solve_eval_eq thm =
1092  case try (fn _ => @{thm refl} RS thm) () of
1093    SOME thm' => thm'
1094  | NONE => 
1095      case try (fn _ => @{thm eval_real_def} RS thm) () of
1096        SOME thm' => thm'
1097      | NONE => @{thm eval_ms.simps} RS thm
1098
1099(*
1100  Returns an expansion theorem for the logarithm of the given expansion.
1101  May add one additional element to the basis at the end.
1102*)
1103fun ln_expansion _ _ _ SEmpty = raise TERM ("ln_expansion: empty basis", [])
1104  | ln_expansion ectxt trimmed_thm thm (SNE basis) =
1105  let
1106    fun trailing_exponent expr (SSng _) = get_exponent expr
1107      | trailing_exponent expr (SCons (_, _, tl)) = trailing_exponent (get_coeff expr) tl
1108    val e = trailing_exponent (get_expansion thm) basis
1109    fun ln_expansion_aux trimmed_thm zero_thm thm basis =
1110      let
1111        val t = betapply (\<^term>\<open>\<lambda>(f::real \<Rightarrow> real) x. f x - 1 :: real\<close>, get_expanded_fun thm)
1112      in
1113        case ev_zeroness_oracle ectxt t of
1114          NONE => ln_expansion_aux' trimmed_thm zero_thm thm basis
1115        | SOME zero_thm =>
1116            @{thm expands_to_ln_eventually_1} OF 
1117              [get_basis_wf_thm' basis, mk_expansion_level_eq_thm' basis, zero_thm]
1118      end
1119    and ln_expansion_aux' trimmed_thm zero_thm thm (SSng {wf_thm, ...}) =
1120          ( @{thm expands_to_ln} OF
1121            [trimmed_thm, wf_thm, thm, 
1122              @{thm expands_to_ln_aux_0} OF [zero_thm, @{thm expands_to_ln_const}]])
1123          |> solve_eval_eq
1124      | ln_expansion_aux' trimmed_thm zero_thm thm (SCons ({wf_thm, ...}, {ln_thm, ...}, basis')) =
1125          let
1126            val c_thm = 
1127              ln_expansion_aux (trimmed_thm RS @{thm trimmed_pos_hd_coeff}) zero_thm 
1128                (expands_to_hd thm) basis'
1129            val e = get_exponent (get_expansion thm)
1130            val c_thm' =
1131              case zeroness_oracle true NONE ectxt e of
1132                (IsZero, SOME thm) =>
1133                  @{thm expands_to_ln_to_expands_to_ln_eval [OF expands_to_ln_aux_0]} OF [thm,c_thm]
1134              | _ => 
1135                case try_prove_real_eq false ectxt (e, \<^term>\<open>1::real\<close>) of
1136                  SOME thm => 
1137                    @{thm expands_to_ln_to_expands_to_ln_eval [OF expands_to_ln_aux_1]}
1138                      OF [thm, wf_thm, c_thm, ln_thm]
1139                | NONE => 
1140                    @{thm expands_to_ln_to_expands_to_ln_eval [OF expands_to_ln_aux]} 
1141                      OF [wf_thm, c_thm, ln_thm]
1142          in
1143            (@{thm expands_to_ln} OF [trimmed_thm, wf_thm, thm, c_thm'])
1144            |> solve_eval_eq
1145          end
1146  in
1147    case zeroness_oracle true NONE ectxt e of
1148      (IsZero, SOME zero_thm) => (ln_expansion_aux trimmed_thm zero_thm thm basis, SNE basis)
1149    | _ => 
1150        let
1151          val basis' = insert_ln (SNE basis)
1152          val lifting = mk_lifting (get_basis_list' basis) basis'
1153          val thm' = lift_expands_to_thm lifting thm
1154          val trimmed_thm' = lift_trimmed_pos_thm lifting trimmed_thm
1155          val (thm'', eq_thm) = retrim_expansion ectxt (thm', basis')
1156          val trimmed_thm'' = @{thm trimmed_pos_eq_cong} OF [trimmed_thm', eq_thm]
1157        in
1158          ln_expansion ectxt trimmed_thm'' thm'' basis'
1159        end
1160  end
1161
1162(*
1163  Handles a possible basis change after expanding exp(c(x)) for an expansion of the form
1164  f(x) = c(x) + g(x). Expanding exp(c(x)) may have inserted an additional basis element. If the 
1165  old basis was b :: bs (i.e. c is an expansion w.r.t. bs) and the updated one is bs' (which
1166  agrees with bs except for one additional element b'), we need to argue that b :: bs' is still
1167  well-formed. This may require us to show that ln(b') is o(ln(b)), which the function takes
1168  as an argument.
1169*)
1170fun adjust_exp_basis basis basis' ln_smallo_thm =
1171  if length (get_basis_list basis) = length (get_basis_list basis') + 1 then
1172    basis
1173  else
1174    let
1175      val SNE (SCons (info, ln_info, tail)) = basis
1176      val SNE tail' = basis'
1177      val wf_thms = map get_basis_wf_thm [basis, basis']
1178      val wf_thm' = 
1179        case
1180          get_first (fn f => try f ())
1181            [fn _ => @{thm basis_wf_lift_modification} OF wf_thms,
1182             fn _ => @{thm basis_wf_insert_exp_near} OF (wf_thms @ [ln_smallo_thm]),
1183             fn _ => @{thm basis_wf_insert_exp_near} OF (wf_thms @ 
1184               [ln_smallo_thm RS @{thm basis_wf_insert_exp_uminus'}])] of
1185          SOME wf_thm => wf_thm
1186        | _ => raise TERM ("Lifting basis modification in exp_expansion failed.", map Thm.concl_of (wf_thms @ [ln_smallo_thm]))
1187      val info' = {wf_thm = wf_thm', head = #head info}
1188      val lifting = mk_lifting (get_basis_list' tail) basis'
1189      val ln_info' = 
1190        {trimmed_thm = lift_trimmed_pos_thm lifting (#trimmed_thm ln_info),
1191         ln_thm = lift_expands_to_thm lifting (#ln_thm ln_info)}
1192    in
1193      SNE (SCons (info', ln_info', tail'))
1194    end
1195
1196(* inserts the exponential of a given function at the beginning of the given basis *)
1197fun insert_exp _ _ _ _ _ SEmpty = raise TERM ("insert_exp", [])
1198  | insert_exp t ln_thm ln_smallo_thm ln_trimmed_thm lim_thm (SNE basis) =
1199      let
1200        val head = Envir.beta_eta_contract (\<^term>\<open>\<lambda>(f::real\<Rightarrow>real) x. exp (f x)\<close> $ t)
1201        val ln_smallo_thm = ln_smallo_thm RS @{thm ln_smallo_ln_exp}
1202        val wf_thm = @{thm basis_wf_manyI} OF [lim_thm, ln_smallo_thm, get_basis_wf_thm' basis]
1203        val basis' = SNE (SCons ({wf_thm = wf_thm, head = head}, 
1204          {ln_thm = ln_thm, trimmed_thm = ln_trimmed_thm} , basis))
1205      in
1206        check_basis basis'
1207      end 
1208
1209(* 
1210  Returns an expansion of the exponential of the given expansion. This may add several
1211  new basis elements at any position of the basis (except at the very end
1212*)
1213fun exp_expansion _ thm SEmpty = (thm RS @{thm expands_to_exp_real}, SEmpty)
1214  | exp_expansion ectxt thm basis =
1215    let
1216      val (_, thm, _) = whnf_expansion ectxt thm
1217    in
1218      case ev_zeroness_oracle ectxt (get_eval (get_expansion thm)) of
1219        SOME zero_thm => 
1220          (@{thm expands_to_exp_zero} OF 
1221             [thm, zero_thm, get_basis_wf_thm basis, mk_expansion_level_eq_thm basis], basis)
1222      | NONE =>
1223          let
1224            val ln =
1225              Option.map (fn x => (#ln_thm x, #trimmed_thm x)) (get_ln_info basis)
1226            val ln = Option.map (fn (x, y) => retrim_pos_expansion ectxt (x, basis, y)) ln
1227            val es' = \<^term>\<open>0::real\<close> :: (
1228              case ln of
1229                NONE => []
1230              | SOME (ln_thm, _, _) => get_exponents (get_expansion ln_thm))
1231            val trim_result =
1232              trim_expansion_while_greater true (SOME es') false (SOME Simple_Trim) ectxt (thm, basis)
1233          in
1234            exp_expansion' ectxt trim_result ln basis
1235          end
1236    end
1237and exp_expansion' _ (thm, _, _) _ SEmpty = (thm RS @{thm expands_to_exp_real}, SEmpty)
1238  | exp_expansion' ectxt (thm, trim_result, e_thms) ln basis =
1239  let
1240    val exp = get_expansion thm
1241    val wf_thm = get_basis_wf_thm basis
1242    val f = get_expanded_fun thm
1243    fun exp_expansion_insert ln_smallo_thm = (
1244      case determine_trimmed_sgn ectxt exp of
1245        (IsPos, trimmed_thm) =>
1246          let
1247            val [lim_thm, ln_thm', thm'] =
1248              @{thms expands_to_exp_insert_pos}
1249              |> map (fn thm' => thm' OF [thm, wf_thm, trimmed_thm, ln_smallo_thm])
1250            val basis' = insert_exp f ln_thm' ln_smallo_thm trimmed_thm lim_thm basis
1251          in
1252            (thm', basis')
1253          end
1254      | (IsNeg, trimmed_thm) =>
1255          let
1256            val [lim_thm, ln_thm', ln_trimmed_thm, thm'] = 
1257              @{thms expands_to_exp_insert_neg}
1258              |> map (fn thm' => thm' OF [thm, wf_thm, trimmed_thm, ln_smallo_thm])
1259            val ln_smallo_thm = ln_smallo_thm RS @{thm basis_wf_insert_exp_uminus}
1260            val f' = Envir.beta_eta_contract (\<^term>\<open>\<lambda>(f::real\<Rightarrow>real) x. -f x\<close> $ f)
1261            val basis' = insert_exp f' ln_thm' ln_smallo_thm ln_trimmed_thm lim_thm basis
1262          in
1263            (thm', basis')
1264          end
1265      | _ => raise TERM ("Unexpected zeroness result in exp_expansion", []))
1266    fun lexord (IsNeg :: _) = LESS
1267      | lexord (IsPos :: _) = GREATER
1268      | lexord (IsZero :: xs) = lexord xs
1269      | lexord [] = EQUAL
1270      | lexord _ = raise Match
1271    val compare_result = lexord (map fst e_thms)
1272  in
1273    case (trim_result, e_thms, compare_result) of
1274      (Aborted _, (IsNeg, e_neg_thm) :: _, _) =>
1275        (* leading exponent is negative; we can simply Taylor-expand exp(x) around 0 *)
1276        (@{thm expands_to_exp_neg} OF [thm, get_basis_wf_thm basis, e_neg_thm], basis)
1277    | (Trimmed (_, SOME trimmed_thm), (IsPos, e_pos_thm) :: _, GREATER) =>
1278        (* leading exponent is positive; exp(f(x)) or exp(-f(x)) is new basis element *)
1279        let
1280          val ln_smallo_thm =
1281            @{thm basis_wf_insert_exp_pos} OF [thm, get_basis_wf_thm basis, trimmed_thm, e_pos_thm]
1282        in
1283          exp_expansion_insert ln_smallo_thm
1284        end
1285    | (Trimmed (_, SOME trimmed_thm), _, GREATER) =>
1286        (* leading exponent is zero, but f(x) grows faster than ln(b(x)), so 
1287           exp(f(x)) or exp(-f(x)) must still be new basis elements *)
1288        let
1289          val ln_thm =
1290            case ln of
1291              SOME (ln_thm, _, _) => ln_thm
1292            | NONE => raise TERM ("TODO blubb", [])
1293          val ln_thm = @{thm expands_to_lift''} OF [get_basis_wf_thm basis, ln_thm]
1294          val ln_smallo_thm = 
1295             @{thm compare_expansions_GT} OF [prove_compare_expansions GREATER (map snd e_thms),
1296               trimmed_thm, thm, ln_thm, get_basis_wf_thm basis]
1297        in
1298          exp_expansion_insert ln_smallo_thm
1299        end
1300    | (Aborted LESS, (IsZero, e_zero_thm) :: e_thms', _) =>
1301        (* leading exponent is zero and f(x) grows more slowly than ln(b(x)), so 
1302           we can write f(x) = c(x) + g(x) and therefore exp(f(x)) = exp(c(x)) * exp(g(x)).
1303           The former is treated by a recursive call; the latter by Taylor expansion. *)
1304        let
1305          val (ln_thm, trimmed_thm) =
1306            case ln of
1307              SOME (ln_thm, _, trimmed_thm) =>
1308                (ln_thm, trimmed_thm RS @{thm trimmed_pos_imp_trimmed})
1309            | NONE => raise TERM ("TODO foo", [])
1310          val c_thm = expands_to_hd thm
1311          val ln_smallo_thm =
1312            @{thm compare_expansions_LT} OF [prove_compare_expansions LESS (map snd e_thms'),
1313              trimmed_thm, c_thm, ln_thm, get_basis_wf_thm (tl_basis basis)]
1314          val (c_thm, c_basis) = exp_expansion ectxt c_thm (tl_basis basis)
1315          val basis' = adjust_exp_basis basis c_basis ln_smallo_thm
1316          val wf_thm = get_basis_wf_thm basis'
1317          val thm' = lift basis' thm
1318          val (thm'', _) = retrim_expansion ectxt (thm', basis')
1319        in
1320          (@{thm expands_to_exp_0} OF [thm'', wf_thm, e_zero_thm, c_thm], basis')
1321        end
1322    | (Trimmed _, [(IsZero, e_zero_thm)], EQUAL) =>
1323        (* f(x) can be written as c + g(x) where c is just a real constant.
1324           We can therefore write exp(f(x)) = exp(c) * exp(g(x)), where the latter is
1325           a simple Taylor expansion. *)
1326        (@{thm expands_to_exp_0_real} OF [thm, wf_thm, e_zero_thm], basis)
1327    | (Trimmed _, (_, e_zero_thm) :: _, EQUAL) =>
1328        (* f(x) is asymptotically equivalent to c * ln(b(x)), so we can write f(x) as
1329           c * ln(b(x)) + g(x) and therefore exp(f(x)) = b(x)^c * exp(g(x)). The second
1330           factor is handled by a recursive call *)
1331        let
1332          val ln_thm =
1333            case ln of
1334              SOME (ln_thm, _, _) => ln_thm
1335            | NONE => raise TERM ("TODO blargh", [])
1336          val c =
1337            case (thm, ln_thm) |> apply2 (get_expansion #> get_lead_coeff) of
1338              (c1, c2) => \<^term>\<open>(/) :: real => _\<close> $ c1 $ c2
1339          val c = Thm.cterm_of (get_ctxt ectxt) c
1340          
1341          val thm' = 
1342            @{thm expands_to_exp_0_pull_out1} 
1343                OF [thm, ln_thm, wf_thm, e_zero_thm, Thm.reflexive c]
1344          val (thm'', basis') = exp_expansion ectxt thm' basis
1345          val pat = ConsPat ("MS", [AnyPat ("_", 0), AnyPat ("_", 0)])
1346          val (_, _, conv) = match ectxt pat (get_expansion thm'') (SOME [])
1347          val eq_thm = conv (Thm.cterm_of (get_ctxt ectxt) (get_expansion thm''))
1348          val thm''' = @{thm expands_to_meta_eq_cong} OF [thm'', eq_thm]
1349          val thm'''' = 
1350            case get_intyness (get_ctxt ectxt) c of
1351              No_Nat =>
1352                @{thm expands_to_exp_0_pull_out2} OF [thm''', get_basis_wf_thm basis']
1353             | Nat nat_thm =>
1354                @{thm expands_to_exp_0_pull_out2_nat} OF 
1355                  [thm''', get_basis_wf_thm basis', nat_thm]
1356             | Neg_Nat nat_thm =>
1357                @{thm expands_to_exp_0_pull_out2_neg_nat} OF 
1358                  [thm''', get_basis_wf_thm basis', nat_thm]
1359        in
1360          (thm'''', basis')
1361        end
1362    | (Trimmed (IsZero, _), [], _) =>
1363        (* Expansion is empty, i.e. f(x) is identically zero *)
1364        (@{thm expands_to_exp_MSLNil} OF [thm, get_basis_wf_thm basis], basis)
1365    | (Trimmed (_, NONE), _, GREATER) =>
1366        (* We could not determine whether f(x) grows faster than ln(b(x)) or not. *)
1367        raise_trimming_error ectxt thm
1368    | _ => raise Match
1369  end
1370
1371fun powr_expansion ectxt (thm1, thm2, basis) =
1372      case ev_zeroness_oracle ectxt (get_expanded_fun thm1) of
1373        SOME zero_thm =>
1374          (@{thm expands_to_powr_0} OF
1375             [zero_thm, Thm.reflexive (Thm.cterm_of (get_ctxt ectxt) (get_expanded_fun thm2)),
1376              get_basis_wf_thm basis, mk_expansion_level_eq_thm basis],
1377           basis)
1378      | NONE =>
1379          let
1380            val (thm1, _, SOME trimmed_thm) =
1381              trim_expansion true (SOME Pos_Trim) ectxt (thm1, basis)
1382            val (ln_thm, basis') = ln_expansion ectxt trimmed_thm thm1 basis
1383            val thm2' = lift basis' thm2 |> simplify_expansion ectxt
1384            val mult_thm = @{thm expands_to_mult} OF [get_basis_wf_thm basis', ln_thm, thm2']
1385            val (exp_thm, basis'') = exp_expansion ectxt mult_thm basis'
1386            val thm = @{thm expands_to_powr} OF 
1387              [trimmed_thm, get_basis_wf_thm basis, thm1, exp_thm]
1388          in  
1389            (thm, basis'')
1390          end
1391
1392fun powr_nat_expansion ectxt (thm1, thm2, basis) =
1393      case ev_zeroness_oracle ectxt (get_expanded_fun thm1) of
1394        SOME zero_thm => (
1395          case ev_zeroness_oracle ectxt (get_expanded_fun thm2) of
1396            SOME zero'_thm => (@{thm expands_to_powr_nat_0_0} OF
1397             [zero_thm, zero'_thm, get_basis_wf_thm basis, mk_expansion_level_eq_thm basis], basis)
1398          | NONE => (
1399              case trim_expansion true (SOME Simple_Trim) ectxt (thm2, basis) of
1400                (thm2, _, SOME trimmed_thm) =>
1401                  (@{thm expands_to_powr_nat_0} OF [zero_thm, thm2, trimmed_thm, 
1402                     get_basis_wf_thm basis, mk_expansion_level_eq_thm basis], basis)))
1403      | NONE =>
1404          let
1405            val (thm1, _, SOME trimmed_thm) =
1406              trim_expansion true (SOME Pos_Trim) ectxt (thm1, basis)
1407            val (ln_thm, basis') = ln_expansion ectxt trimmed_thm thm1 basis
1408            val thm2' = lift basis' thm2 |> simplify_expansion ectxt
1409            val mult_thm = @{thm expands_to_mult} OF [get_basis_wf_thm basis', ln_thm, thm2']
1410            val (exp_thm, basis'') = exp_expansion ectxt mult_thm basis'
1411            val thm = @{thm expands_to_powr_nat} OF 
1412              [trimmed_thm, get_basis_wf_thm basis, thm1, exp_thm]
1413          in  
1414            (thm, basis'')
1415          end
1416
1417fun is_numeral t =
1418  let
1419    val _ = HOLogic.dest_number t
1420  in
1421    true
1422  end
1423    handle TERM _ => false
1424
1425fun power_expansion ectxt (thm, n, basis) =
1426      case ev_zeroness_oracle ectxt (get_expanded_fun thm) of
1427        SOME zero_thm => @{thm expands_to_power_0} OF
1428          [zero_thm, get_basis_wf_thm basis, mk_expansion_level_eq_thm basis,
1429             Thm.reflexive (Thm.cterm_of (get_ctxt ectxt) n)]
1430      | NONE => (
1431          case trim_expansion true (SOME Simple_Trim) ectxt (thm, basis) of
1432            (thm', _, SOME trimmed_thm) =>
1433              let
1434                val ctxt = get_ctxt ectxt
1435                val thm =
1436                  if is_numeral n then @{thm expands_to_power[where abort = True]}
1437                    else @{thm expands_to_power[where abort = False]}
1438                val thm = 
1439                  Drule.infer_instantiate' ctxt [NONE, NONE, NONE, SOME (Thm.cterm_of ctxt n)] thm
1440              in                
1441                thm OF [trimmed_thm, get_basis_wf_thm basis, thm']
1442              end
1443          | _ => raise TERM ("Unexpected zeroness result in power_expansion", []))
1444
1445fun powr_const_expansion ectxt (thm, p, basis) =
1446  let
1447    val pthm = Thm.reflexive (Thm.cterm_of (get_ctxt ectxt) p)
1448  in
1449    case ev_zeroness_oracle ectxt (get_expanded_fun thm) of
1450      SOME zero_thm => @{thm expands_to_powr_const_0} OF 
1451        [zero_thm, get_basis_wf_thm basis, mk_expansion_level_eq_thm basis, pthm]
1452    | NONE =>
1453        case trim_expansion true (SOME Pos_Trim) ectxt (thm, basis) of
1454          (_, _, NONE) => raise TERM ("Unexpected zeroness result for powr", [])
1455        | (thm, _, SOME trimmed_thm) =>
1456            (if is_numeral p then @{thm expands_to_powr_const[where abort = True]}
1457                 else @{thm expands_to_powr_const[where abort = False]})
1458               OF [trimmed_thm, get_basis_wf_thm basis, thm, pthm]
1459  end
1460
1461fun sgn_expansion ectxt (thm, basis) =
1462  let
1463    val thms = [get_basis_wf_thm basis, mk_expansion_level_eq_thm basis]
1464  in
1465    case ev_zeroness_oracle ectxt (get_expanded_fun thm) of
1466      SOME zero_thm => @{thm expands_to_sgn_zero} OF (zero_thm :: thms)
1467    | NONE =>
1468        case trim_expansion true (SOME Sgn_Trim) ectxt (thm, basis) of
1469          (thm, IsPos, SOME trimmed_thm) =>
1470            @{thm expands_to_sgn_pos} OF ([trimmed_thm, thm] @ thms)
1471        | (thm, IsNeg, SOME trimmed_thm) =>
1472            @{thm expands_to_sgn_neg} OF ([trimmed_thm, thm] @ thms)
1473        | _ => raise TERM ("Unexpected zeroness result in sgn_expansion", [])
1474  end
1475
1476(*
1477  Returns an expansion of the sine and cosine of the given expansion. Fails if that function
1478  goes to infinity.
1479*)
1480fun sin_cos_expansion _ thm SEmpty =
1481      (thm RS @{thm expands_to_sin_real}, thm RS @{thm expands_to_cos_real})
1482  | sin_cos_expansion ectxt thm basis =
1483      let
1484        val exp = get_expansion thm
1485        val e = get_exponent exp
1486      in
1487        case zeroness_oracle true (SOME Sgn_Trim) ectxt e of
1488          (IsPos, _) => raise THM ("sin_cos_expansion", 0, [thm])
1489        | (IsNeg, SOME e_thm) =>
1490            let
1491              val [thm1, thm2] = 
1492                map (fn thm' => thm' OF [e_thm, get_basis_wf_thm basis, thm]) 
1493                  @{thms expands_to_sin_ms_neg_exp expands_to_cos_ms_neg_exp}
1494            in
1495              (thm1, thm2)
1496            end
1497        | (IsZero, SOME e_thm) =>
1498            let
1499              val (sin_thm, cos_thm) = (sin_cos_expansion ectxt (expands_to_hd thm) (tl_basis basis))
1500              fun mk_thm thm' = 
1501                (thm' OF [e_thm, get_basis_wf_thm basis, thm, sin_thm, cos_thm]) |> solve_eval_eq
1502              val [thm1, thm2] = 
1503                map mk_thm @{thms expands_to_sin_ms_zero_exp expands_to_cos_ms_zero_exp}
1504            in
1505              (thm1, thm2)
1506            end
1507        | _ => raise TERM ("Unexpected zeroness result in sin_exp_expansion", [])
1508      end
1509
1510fun abconv (t, t') = Envir.beta_eta_contract t aconv Envir.beta_eta_contract t'
1511
1512(*
1513  Makes sure that an expansion theorem really talks about the right function.
1514  This is basically a sanity check to make things fail early and in the right place.
1515*)
1516fun check_expansion e thm =
1517  if abconv (expr_to_term e, get_expanded_fun thm) then 
1518    thm 
1519  else
1520(* TODO Remove Debugging stuff *)
1521let val _ = \<^print> e
1522val _ = \<^print> (get_expanded_fun thm)
1523in
1524    raise TERM ("check_expansion", [Thm.concl_of thm, expr_to_term e])
1525end
1526
1527fun minmax_expansion max [less_thm, eq_thm, gt_thm] ectxt (thm1, thm2, basis) = (
1528      case compare_expansions ectxt (thm1, thm2, basis) of
1529        (LESS, less_thm', thm1, thm2) => less_thm OF [if max then thm2 else thm1, less_thm']
1530      | (GREATER, gt_thm', thm1, thm2) => gt_thm OF [if max then thm1 else thm2, gt_thm']
1531      | (EQUAL, eq_thm', thm1, _) => eq_thm OF [thm1, eq_thm'])
1532  | minmax_expansion _ _ _ _ = raise Match
1533
1534val min_expansion =
1535  minmax_expansion false @{thms expands_to_min_lt expands_to_min_eq expands_to_min_gt}
1536val max_expansion =
1537  minmax_expansion true @{thms expands_to_max_lt expands_to_max_eq expands_to_max_gt}
1538
1539fun zero_expansion basis =
1540  @{thm expands_to_zero} OF [get_basis_wf_thm basis, mk_expansion_level_eq_thm basis]
1541
1542fun const_expansion _ basis \<^term>\<open>0 :: real\<close> = zero_expansion basis
1543  | const_expansion ectxt basis t =
1544  let
1545    val ctxt = get_ctxt ectxt
1546    val thm = Drule.infer_instantiate' ctxt [NONE, SOME (Thm.cterm_of ctxt t)] 
1547                @{thm expands_to_const}
1548  in
1549    thm OF [get_basis_wf_thm basis, mk_expansion_level_eq_thm basis]
1550  end
1551
1552fun root_expansion ectxt (thm, n, basis) =
1553  let
1554    val ctxt = get_ctxt ectxt
1555    fun tac {context = ctxt, ...} =
1556      HEADGOAL (Method.insert_tac ctxt (get_facts ectxt))
1557      THEN Local_Defs.unfold_tac ctxt eval_simps
1558      THEN HEADGOAL (Simplifier.asm_full_simp_tac ctxt)
1559    fun prove goal =
1560      try (Goal.prove ctxt [] [] (HOLogic.mk_Trueprop (Term.betapply (goal, n)))) tac
1561    fun err () =
1562      let
1563        val t = simplify_term' (get_facts ectxt) ctxt n
1564        val _ =
1565          if #verbose (#ctxt ectxt) then
1566            let
1567              val p = Pretty.str ("real_asymp failed to determine whether the following constant " ^
1568                "is zero or not:")
1569              val p = Pretty.chunks [p, Pretty.indent 2 (Syntax.pretty_term ctxt t)]
1570            in
1571              Pretty.writeln p
1572            end else ()
1573      in
1574        raise TERM ("zeroness_oracle", [n])
1575      end
1576    fun aux nz_thm =
1577      case trim_expansion true (SOME Sgn_Trim) ectxt (thm, basis) of
1578        (thm, IsPos, SOME trimmed_thm) =>
1579          @{thm expands_to_root} OF [nz_thm, trimmed_thm, get_basis_wf_thm basis, thm]
1580      | (thm, IsNeg, SOME trimmed_thm) =>
1581          @{thm expands_to_root_neg} OF [nz_thm, trimmed_thm, get_basis_wf_thm basis, thm]
1582      | _ => raise TERM ("Unexpected zeroness result in root_expansion", [])
1583  in
1584    case prove \<^term>\<open>\<lambda>n::nat. n = 0\<close> of
1585      SOME zero_thm =>
1586        @{thm expands_to_0th_root} OF
1587          [zero_thm, get_basis_wf_thm basis, mk_expansion_level_eq_thm basis, 
1588             Thm.reflexive (Thm.cterm_of ctxt (get_expanded_fun thm))]
1589    | NONE =>
1590        case prove \<^term>\<open>\<lambda>n::nat. n > 0\<close> of
1591          NONE => err ()
1592        | SOME nz_thm =>
1593            case ev_zeroness_oracle ectxt (get_expanded_fun thm) of
1594              SOME zero_thm => @{thm expands_to_root_0} OF
1595                [nz_thm, zero_thm, get_basis_wf_thm basis, mk_expansion_level_eq_thm basis]
1596            | NONE => aux nz_thm
1597  end
1598
1599
1600fun arctan_expansion _ SEmpty thm =
1601      @{thm expands_to_real_compose[where g = arctan]} OF [thm]
1602  | arctan_expansion ectxt basis thm =
1603      case ev_zeroness_oracle ectxt (get_expanded_fun thm) of
1604        SOME zero_thm => @{thm expands_to_arctan_zero} OF [zero_expansion basis, zero_thm]
1605      | NONE =>
1606          let
1607            val (thm, _, _) = trim_expansion true (SOME Simple_Trim) ectxt (thm, basis)
1608            val e = get_exponent (get_expansion thm)
1609            fun cont ectxt (thm, basis) = arctan_expansion ectxt basis thm
1610          in
1611            case zeroness_oracle true (SOME Sgn_Trim) ectxt e of
1612                (IsNeg, SOME neg_thm) =>
1613                  @{thm expands_to_arctan_ms_neg_exp} OF [neg_thm, get_basis_wf_thm basis, thm]
1614              | (IsPos, SOME e_pos_thm) => (
1615                  case determine_trimmed_sgn ectxt (get_expansion thm) of
1616                    (IsPos, trimmed_thm) =>
1617                      @{thm expands_to_arctan_ms_pos_exp_pos} OF 
1618                        [e_pos_thm, trimmed_thm, get_basis_wf_thm basis, thm]
1619                  | (IsNeg, trimmed_thm) =>
1620                      @{thm expands_to_arctan_ms_pos_exp_neg} OF 
1621                        [e_pos_thm, trimmed_thm, get_basis_wf_thm basis, thm]
1622                  | _ => raise TERM ("Unexpected trim result during expansion of arctan", []))
1623              | (IsZero, _) => (
1624                  case try_lift_function ectxt (thm, basis) cont of
1625                    (SOME thm', _) => thm'
1626                  | _ =>
1627                      let
1628                        val _ = if get_verbose ectxt then 
1629                          writeln "Unsupported occurrence of arctan" else ()
1630                      in
1631                        raise TERM ("Unsupported occurence of arctan", [])
1632                      end)
1633              | _ => raise TERM ("Unexpected trim result during expansion of arctan", [])
1634          end
1635
1636(* Returns an expansion theorem for a function that is already a basis element *)
1637fun expand_basic _ t SEmpty = raise TERM ("expand_basic", [t])
1638  | expand_basic thm t basis =
1639      if abconv (get_basis_head basis, t) then
1640        thm (get_basis_wf_thm basis) (mk_expansion_level_eq_thm (tl_basis basis))
1641      else
1642        @{thm expands_to_lift'} OF [get_basis_wf_thm basis, expand_basic thm t (tl_basis basis)]
1643  
1644fun expand_unary ectxt thm e basis =
1645      let
1646        val (thm', basis') = expand' ectxt e basis |> apfst (simplify_expansion ectxt)
1647      in
1648        (thm OF [get_basis_wf_thm basis', thm'], basis')
1649      end
1650and expand_binary ectxt thm (e1, e2) basis =
1651      let
1652        val (thm1, basis') = expand' ectxt e1 basis |> apfst (simplify_expansion ectxt)
1653        val (thm2, basis'') = expand' ectxt e2 basis' |> apfst (simplify_expansion ectxt)
1654        val thm1 = lift basis'' thm1 |> simplify_expansion ectxt
1655      in
1656        (thm OF [get_basis_wf_thm basis'', thm1, thm2], basis'')
1657      end
1658and trim_nz mode ectxt e basis =
1659      let
1660        val (thm, basis') = expand' ectxt e basis |> apfst (simplify_expansion ectxt)
1661        val (thm', nz, trimmed_thm) = trim_expansion true (SOME mode) ectxt (thm, basis')
1662      in
1663        case trimmed_thm of
1664          NONE => raise TERM ("expand: zero denominator", [get_expansion thm])
1665        | SOME trimmed_thm => (thm', basis', nz, trimmed_thm)
1666      end
1667and expand'' ectxt (ConstExpr c) basis = (const_expansion ectxt basis c, basis)
1668  | expand'' _ X basis = (lift basis @{thm expands_to_X}, basis)
1669  | expand'' ectxt (Uminus e) basis = expand_unary ectxt @{thm expands_to_uminus} e basis
1670  | expand'' ectxt (Add e12) basis = expand_binary ectxt @{thm expands_to_add} e12 basis
1671  | expand'' ectxt (Minus e12) basis = expand_binary ectxt @{thm expands_to_minus} e12 basis
1672  | expand'' ectxt (Mult e12) basis = expand_binary ectxt @{thm expands_to_mult} e12 basis
1673  | expand'' ectxt (Powr' (e, p)) basis = (* TODO zero basis *)
1674      let
1675        val (thm, basis') = expand' ectxt e basis |> apfst (simplify_expansion ectxt)
1676      in
1677        (powr_const_expansion ectxt (thm, p, basis'), basis')
1678      end
1679  | expand'' ectxt (Powr (e1, e2)) basis =
1680      let
1681        val (thm2, basis1) = expand' ectxt e2 basis |> apfst (simplify_expansion ectxt)
1682        val (thm1, basis2) = expand' ectxt e1 basis1 |> apfst (simplify_expansion ectxt)
1683      in
1684        powr_expansion ectxt (thm1, thm2, basis2)
1685      end
1686  | expand'' ectxt (Powr_Nat (e1, e2)) basis =
1687      let
1688        val (thm2, basis1) = expand' ectxt e2 basis |> apfst (simplify_expansion ectxt)
1689        val (thm1, basis2) = expand' ectxt e1 basis1 |> apfst (simplify_expansion ectxt)
1690      in
1691        powr_nat_expansion ectxt (thm1, thm2, basis2)
1692      end
1693  | expand'' ectxt (LnPowr (e1, e2)) basis =
1694      let (* TODO zero base *)
1695        val (thm2, basis1) = expand' ectxt e2 basis |> apfst (simplify_expansion ectxt)
1696        val (thm1, basis2, _, trimmed_thm) = trim_nz Pos_Trim ectxt e1 basis1
1697        val (ln_thm, basis3) = ln_expansion ectxt trimmed_thm thm1 basis2
1698        val thm2' = lift basis3 thm2 |> simplify_expansion ectxt
1699        val mult_thm = @{thm expands_to_mult} OF [get_basis_wf_thm basis3, ln_thm, thm2']
1700        val thm = @{thm expands_to_ln_powr} OF 
1701          [trimmed_thm, get_basis_wf_thm basis2, thm1, mult_thm]
1702      in  
1703        (thm, basis3)
1704      end
1705  | expand'' ectxt (ExpLn e) basis =
1706      let
1707        val (thm, basis', _, trimmed_thm) = trim_nz Pos_Trim ectxt e basis
1708        val thm = @{thm expands_to_exp_ln} OF [trimmed_thm, get_basis_wf_thm basis', thm]
1709      in  
1710        (thm, basis')
1711      end
1712  | expand'' ectxt (Power (e, n)) basis =
1713      let
1714        val (thm, basis') = expand' ectxt e basis |> apfst (simplify_expansion ectxt)
1715      in
1716        (power_expansion ectxt (thm, n, basis'), basis')
1717      end
1718  | expand'' ectxt (Root (e, n)) basis =
1719      let
1720        val (thm, basis') = expand' ectxt e basis |> apfst (simplify_expansion ectxt)
1721      in
1722        (root_expansion ectxt (thm, n, basis'), basis')
1723      end
1724  | expand'' ectxt (Inverse e) basis = 
1725      (case trim_nz Simple_Trim ectxt e basis of
1726         (thm, basis', _, trimmed_thm) => 
1727           (@{thm expands_to_inverse} OF [trimmed_thm, get_basis_wf_thm basis', thm], basis'))
1728  | expand'' ectxt (Div (e1, e2)) basis =
1729      let
1730        val (thm1, basis') = expand' ectxt e1 basis
1731        val (thm2, basis'', _, trimmed_thm) = trim_nz Simple_Trim ectxt e2 basis'
1732        val thm1 = lift basis'' thm1
1733      in
1734        (@{thm expands_to_divide} OF [trimmed_thm, get_basis_wf_thm basis'', thm1, thm2], basis'')
1735      end
1736  | expand'' ectxt (Ln e) basis =
1737      let
1738        val (thm, basis', _, trimmed_thm) = trim_nz Pos_Trim ectxt e basis
1739      in
1740        ln_expansion ectxt trimmed_thm thm basis'
1741      end
1742  | expand'' ectxt (Exp e) basis =
1743      let
1744        val (thm, basis') = expand' ectxt e basis
1745      in
1746        exp_expansion ectxt thm basis'
1747      end
1748  | expand'' ectxt (Absolute e) basis =
1749      let
1750        val (thm, basis', nz, trimmed_thm) = trim_nz Sgn_Trim ectxt e basis
1751        val thm' =
1752          case nz of 
1753            IsPos => @{thm expands_to_abs_pos} 
1754          | IsNeg => @{thm expands_to_abs_neg}
1755          | _ => raise TERM ("Unexpected trim result during expansion of abs", [])
1756      in
1757        (thm' OF [trimmed_thm, get_basis_wf_thm basis', thm], basis')
1758      end
1759  | expand'' ectxt (Sgn e) basis =
1760      let
1761        val (thm, basis') = expand' ectxt e basis
1762      in
1763        (sgn_expansion ectxt (thm, basis'), basis')
1764      end
1765  | expand'' ectxt (Min (e1, e2)) basis = (
1766      case try_prove_ev_eq ectxt (apply2 expr_to_term (e1, e2)) of
1767        SOME eq_thm =>
1768          expand' ectxt e1 basis
1769          |> apfst (fn thm => @{thm expands_to_min_eq} OF [thm, eq_thm])
1770      | NONE =>
1771          let
1772            val (thm1, basis') = expand' ectxt e1 basis
1773            val (thm2, basis'') = expand' ectxt e2 basis'
1774            val thm1' = lift basis'' thm1
1775          in
1776            (min_expansion ectxt (thm1', thm2, basis''), basis'')
1777          end)
1778  | expand'' ectxt (Max (e1, e2)) basis = (
1779      case try_prove_ev_eq ectxt (apply2 expr_to_term (e1, e2)) of
1780        SOME eq_thm =>
1781          expand' ectxt e1 basis
1782          |> apfst (fn thm => @{thm expands_to_max_eq} OF [thm, eq_thm])
1783      | NONE =>
1784          let
1785            val (thm1, basis') = expand' ectxt e1 basis
1786            val (thm2, basis'') = expand' ectxt e2 basis'
1787            val thm1' = lift basis'' thm1
1788          in
1789            (max_expansion ectxt (thm1', thm2, basis''), basis'')
1790          end)
1791  | expand'' ectxt (Sin e) basis =
1792      let
1793        val (thm, basis', _, _) = trim_nz Simple_Trim ectxt e basis (* TODO could be relaxed *)
1794      in
1795        (sin_cos_expansion ectxt thm basis' |> fst, basis')
1796      end
1797  | expand'' ectxt (Cos e) basis =
1798      let
1799        val (thm, basis', _, _) = trim_nz Simple_Trim ectxt e basis (* TODO could be relaxed *)
1800      in
1801        (sin_cos_expansion ectxt thm basis' |> snd, basis')
1802      end
1803  | expand'' _ (Floor _) _ =
1804      raise TERM ("floor not supported.", [])
1805  | expand'' _ (Ceiling _) _ =
1806      raise TERM ("ceiling not supported.", [])
1807  | expand'' _ (Frac _) _ =
1808      raise TERM ("frac not supported.", [])
1809  | expand'' _ (NatMod _) _ =
1810      raise TERM ("mod not supported.", [])
1811  | expand'' ectxt (ArcTan e) basis =
1812      let
1813        (* TODO: what if it's zero *)
1814        val (thm, basis') = expand' ectxt e basis
1815      in
1816        (arctan_expansion ectxt basis' thm, basis')
1817      end
1818  | expand'' ectxt (Custom (name, t, args)) basis =
1819      let
1820        fun expand_args acc basis [] = (rev acc, basis)
1821          | expand_args acc basis (arg :: args) =
1822              case expand' ectxt arg basis of
1823                (thm, basis') => expand_args (thm :: acc) basis' args          
1824      in
1825        case expand_custom (get_ctxt ectxt) name of
1826          NONE => raise TERM ("Unsupported custom function: " ^ name, [t])
1827        | SOME e => e ectxt t (expand_args [] basis args)
1828      end
1829
1830and expand' ectxt (e' as (Inverse e)) basis =
1831      let
1832        val t = expr_to_term e
1833        fun thm wf_thm len_thm =
1834          @{thm expands_to_basic_inverse} OF [wf_thm, len_thm]
1835      in
1836        if member abconv (get_basis_list basis) t then
1837          (expand_basic thm t basis, basis) |> apfst (check_expansion e')
1838        else
1839          expand'' ectxt e' basis |> apfst (check_expansion e')
1840      end
1841  | expand' ectxt (Div (e1, e2)) basis =
1842      let
1843        val (thm1, basis') = expand' ectxt e1 basis
1844        val t = expr_to_term e2
1845        fun thm wf_thm len_thm =
1846          @{thm expands_to_basic_inverse} OF [wf_thm, len_thm]
1847      in
1848        if member abconv (get_basis_list basis') t then
1849          (@{thm expands_to_div'} OF [get_basis_wf_thm basis', thm1, expand_basic thm t basis'], 
1850             basis')
1851        else
1852          let
1853            val (thm2, basis'', _, trimmed_thm) = trim_nz Simple_Trim ectxt e2 basis'
1854            val thm1 = lift basis'' thm1
1855          in
1856            (@{thm expands_to_divide} OF [trimmed_thm, get_basis_wf_thm basis'', thm1, thm2], 
1857               basis'')
1858          end
1859      end
1860  | expand' ectxt (e' as (Powr' (e, p))) basis =
1861      let
1862        val t = expr_to_term e
1863        val ctxt = get_ctxt ectxt
1864        fun thm wf_thm len_thm =
1865          (Drule.infer_instantiate' ctxt [NONE, NONE, SOME (Thm.cterm_of ctxt p)]
1866            @{thm expands_to_basic_powr}) OF [wf_thm, len_thm]
1867      in
1868        if member abconv (get_basis_list basis) t then
1869          (expand_basic thm t basis, basis) |> apfst (check_expansion e')
1870        else
1871          expand'' ectxt e' basis |> apfst (check_expansion e')
1872      end
1873  | expand' ectxt e basis =
1874      let
1875        val t = expr_to_term e
1876        fun thm wf_thm len_thm = @{thm expands_to_basic} OF [wf_thm, len_thm]
1877      in
1878        if member abconv (get_basis_list basis) t then
1879          (expand_basic thm t basis, basis) |> apfst (check_expansion e)
1880        else
1881          expand'' ectxt e basis |> apfst (check_expansion e)
1882      end
1883
1884fun expand ectxt e basis = 
1885  expand' ectxt e basis
1886  |> apfst (simplify_expansion ectxt) 
1887  |> apfst (check_expansion e)
1888
1889fun expand_term ectxt t basis =
1890  let
1891    val ctxt = get_ctxt ectxt
1892    val (e, eq_thm) = reify ctxt t
1893    val (thm,  basis) = expand ectxt e basis
1894  in
1895    (@{thm expands_to_meta_eq_cong'} OF [thm, eq_thm], basis)
1896  end
1897
1898fun expand_terms ectxt ts basis =
1899  let
1900    val ctxt = get_ctxt ectxt
1901    val e_eq_thms = map (reify ctxt) ts
1902    fun step (e, eq_thm) (thms, basis) =
1903      let
1904        val (thm, basis) = expand' ectxt e basis
1905        val thm = @{thm expands_to_meta_eq_cong'} OF [simplify_expansion ectxt thm, eq_thm]
1906      in
1907        (thm :: thms, basis)
1908      end
1909    val (thms, basis) = fold step e_eq_thms ([], basis)
1910    fun lift thm = lift_expands_to_thm (mk_lifting (extract_basis_list thm) basis) thm
1911  in
1912    (map lift (rev thms), basis)
1913  end
1914
1915datatype limit =
1916   Zero_Limit of bool option
1917 | Finite_Limit of term
1918 | Infinite_Limit of bool option
1919
1920fun is_empty_expansion (Const (\<^const_name>\<open>MS\<close>, _) $ Const (\<^const_name>\<open>MSLNil\<close>, _) $ _) = true
1921  | is_empty_expansion _ = false
1922
1923fun limit_of_expansion_aux ectxt basis thm =
1924  let
1925    val n = length (get_basis_list basis)
1926    val (thm, res, e_thms) =
1927      trim_expansion_while_greater false (SOME (replicate n \<^term>\<open>0::real\<close>)) true
1928        (SOME Simple_Trim) ectxt (thm, basis)
1929    val trimmed_thm = case res of Trimmed (_, trimmed_thm) => trimmed_thm | _ => NONE
1930    val res = case res of Trimmed _ => GREATER | Aborted res => res
1931    val exp = get_expansion thm
1932    val _ = if res = GREATER andalso is_none trimmed_thm andalso not (is_empty_expansion exp) then
1933              raise TERM ("limit_of_expansion", [get_expansion thm]) else ()
1934    fun go thm _ _ [] = (
1935          case zeroness_oracle false (SOME Simple_Trim) ectxt (get_expansion thm) of
1936            (IsZero, _) => (Zero_Limit NONE, @{thm expands_to_real_imp_filterlim} OF [thm])
1937          | _ => (Finite_Limit \<^term>\<open>0::real\<close>, @{thm expands_to_real_imp_filterlim} OF [thm]))
1938      | go thm _ basis ((IsNeg, neg_thm) :: _) = (Zero_Limit NONE,
1939          @{thm expands_to_neg_exponent_imp_filterlim} OF
1940            [thm, get_basis_wf_thm basis, neg_thm RS @{thm compare_reals_diff_sgnD(1)}])
1941      | go thm trimmed_thm basis ((IsPos, pos_thm) :: _) = (Infinite_Limit NONE,
1942          @{thm expands_to_pos_exponent_imp_filterlim} OF
1943            [thm, the trimmed_thm, get_basis_wf_thm basis,
1944             pos_thm RS @{thm compare_reals_diff_sgnD(3)}])
1945      | go thm trimmed_thm basis ((IsZero, zero_thm) :: e_thms) =
1946          let
1947             val thm' = thm RS @{thm expands_to_hd''}
1948             val trimmed_thm' = Option.map (fn thm => thm RS @{thm trimmed_hd}) trimmed_thm
1949             val (lim, lim_thm) = go thm' trimmed_thm' (tl_basis basis) e_thms
1950             val lim_lift_thm =
1951                case lim of
1952                  Infinite_Limit _ => @{thm expands_to_zero_exponent_imp_filterlim(1)}
1953                | _ => @{thm expands_to_zero_exponent_imp_filterlim(2)}
1954             val lim_thm' = 
1955               lim_lift_thm OF [thm, get_basis_wf_thm basis, 
1956                 zero_thm RS @{thm compare_reals_diff_sgnD(2)}, lim_thm]
1957          in
1958              (lim, lim_thm')
1959          end
1960      | go _ _ _ _ = raise Match
1961  in
1962    if is_empty_expansion exp then
1963      (Zero_Limit NONE, thm RS @{thm expands_to_MSLNil_imp_filterlim}, thm)
1964    else
1965      case go thm trimmed_thm basis e_thms of
1966        (lim, lim_thm) => (lim, lim_thm, thm)
1967  end
1968
1969(* 
1970  Determines the limit of a function from its expansion. The two flags control whether the
1971  the sign of the approach should be determined for the infinite case (i.e. at_top/at_bot instead
1972  of just at_infinity) and the zero case (i.e. at_right 0/at_left 0 instead of just nhds 0)
1973*)
1974fun limit_of_expansion (sgn_zero, sgn_inf) ectxt (thm, basis) =
1975  let
1976    val (lim, lim_thm, thm) = limit_of_expansion_aux ectxt basis thm
1977  in
1978    case lim of
1979      Zero_Limit _ => (
1980        if sgn_zero then
1981          case trim_expansion false (SOME Sgn_Trim) ectxt (thm, basis) of
1982            (thm, IsPos, SOME pos_thm) => (Zero_Limit (SOME true),
1983              @{thm tendsto_imp_filterlim_at_right[OF _ expands_to_imp_eventually_pos]} OF
1984                [lim_thm, get_basis_wf_thm basis, thm, pos_thm])
1985          | (thm, IsNeg, SOME neg_thm) => (Zero_Limit (SOME false),
1986              @{thm tendsto_imp_filterlim_at_left[OF _ expands_to_imp_eventually_neg]} OF
1987                [lim_thm, get_basis_wf_thm basis, thm, neg_thm])
1988          | _ => (Zero_Limit NONE, lim_thm)
1989        else (Zero_Limit NONE, lim_thm))
1990    | Infinite_Limit _ => (
1991        if sgn_inf then
1992          case trim_expansion false (SOME Sgn_Trim) ectxt (thm, basis) of
1993            (thm, IsPos, SOME pos_thm) => (Infinite_Limit (SOME true),
1994              (@{thm filterlim_at_infinity_imp_filterlim_at_top[OF _ expands_to_imp_eventually_pos]} OF
1995                 [lim_thm, get_basis_wf_thm basis, thm, pos_thm]))
1996          | (thm, IsNeg, SOME neg_thm) => (Infinite_Limit (SOME false),
1997              @{thm filterlim_at_infinity_imp_filterlim_at_bot[OF _ expands_to_imp_eventually_neg]} OF
1998                [lim_thm, get_basis_wf_thm basis, thm, neg_thm])
1999          | _ => (Infinite_Limit NONE, lim_thm)
2000        else (Infinite_Limit NONE, lim_thm))
2001    | Finite_Limit c => (Finite_Limit c, lim_thm)
2002  end
2003
2004fun compute_limit ectxt t =
2005  case expand_term ectxt t default_basis of
2006    (thm, basis) => limit_of_expansion (true, true) ectxt (thm, basis)
2007
2008fun prove_at_infinity ectxt (thm, basis) =
2009  let
2010    fun err () = raise TERM ("prove_at_infinity", [get_expanded_fun thm])
2011    val (thm, _, SOME trimmed_thm) = trim_expansion true (SOME Simple_Trim) ectxt (thm, basis)
2012    fun go basis thm trimmed_thm =
2013      if fastype_of (get_expansion thm) = \<^typ>\<open>real\<close> then
2014        err ()
2015      else
2016        case zeroness_oracle true (SOME Pos_Trim) ectxt (get_exponent (get_expansion thm)) of
2017          (IsPos, SOME pos_thm) =>
2018            @{thm expands_to_pos_exponent_imp_filterlim} OF
2019              [thm, trimmed_thm, get_basis_wf_thm basis, pos_thm]
2020        | (IsZero, SOME zero_thm) =>
2021            @{thm expands_to_zero_exponent_imp_filterlim(1)} OF
2022              [thm, get_basis_wf_thm basis, zero_thm,
2023                 go (tl_basis basis) (thm RS @{thm expands_to_hd''})
2024                   (trimmed_thm RS @{thm trimmed_hd})]
2025        | _ => err ()
2026  in
2027    go basis thm trimmed_thm
2028  end
2029
2030fun prove_at_top_at_bot mode ectxt (thm, basis) =
2031  let
2032    val s = if mode = Pos_Trim then "prove_at_top" else "prove_at_bot"
2033    fun err () = raise TERM (s, [get_expanded_fun thm])
2034    val (thm, _, SOME trimmed_thm) = trim_expansion true (SOME mode) ectxt (thm, basis)
2035    val trimmed_thm' = trimmed_thm RS
2036      (if mode = Pos_Trim then @{thm trimmed_pos_imp_trimmed} else @{thm trimmed_neg_imp_trimmed})
2037    fun go basis thm trimmed_thm =
2038      if fastype_of (get_expansion thm) = \<^typ>\<open>real\<close> then
2039        err ()
2040      else
2041        case zeroness_oracle true (SOME Pos_Trim) ectxt (get_exponent (get_expansion thm)) of
2042          (IsPos, SOME pos_thm) =>
2043            @{thm expands_to_pos_exponent_imp_filterlim} OF
2044              [thm, trimmed_thm, get_basis_wf_thm basis, pos_thm]
2045        | (IsZero, SOME zero_thm) =>
2046            @{thm expands_to_zero_exponent_imp_filterlim(1)} OF
2047              [thm, get_basis_wf_thm basis, zero_thm,
2048                 go (tl_basis basis) (thm RS @{thm expands_to_hd''})
2049                   (trimmed_thm RS @{thm trimmed_hd})]
2050        | _ => err ()
2051    val lim_thm = go basis thm trimmed_thm'
2052    val add_sign_thm =
2053      if mode = Pos_Trim then
2054        @{thm filterlim_at_infinity_imp_filterlim_at_top[OF _ expands_to_imp_eventually_pos]}
2055      else
2056        @{thm filterlim_at_infinity_imp_filterlim_at_bot[OF _ expands_to_imp_eventually_neg]}
2057  in
2058    add_sign_thm OF [lim_thm, get_basis_wf_thm basis, thm, trimmed_thm]
2059  end
2060
2061val prove_at_top = prove_at_top_at_bot Pos_Trim
2062val prove_at_bot = prove_at_top_at_bot Neg_Trim
2063
2064
2065fun prove_at_aux mode ectxt (thm, basis) =
2066  let
2067    val (s, add_sign_thm) =
2068      case mode of
2069        Simple_Trim =>
2070          ("prove_at_0", @{thm Topological_Spaces.filterlim_atI[OF _ expands_to_imp_eventually_nz]})
2071      | Pos_Trim =>
2072          ("prove_at_right_0",
2073             @{thm tendsto_imp_filterlim_at_right[OF _ expands_to_imp_eventually_pos]})
2074      | Neg_Trim =>
2075          ("prove_at_left_0",
2076             @{thm tendsto_imp_filterlim_at_left[OF _ expands_to_imp_eventually_neg]})
2077    fun err () = raise TERM (s, [get_expanded_fun thm])
2078    val (thm, _, SOME trimmed_thm) = trim_expansion true (SOME mode) ectxt (thm, basis)
2079    fun go basis thm =
2080      if fastype_of (get_expansion thm) = \<^typ>\<open>real\<close> then
2081        err ()
2082      else
2083        case zeroness_oracle true (SOME Neg_Trim) ectxt (get_exponent (get_expansion thm)) of
2084          (IsNeg, SOME neg_thm) =>
2085            @{thm expands_to_neg_exponent_imp_filterlim} OF
2086              [thm, get_basis_wf_thm basis, neg_thm]
2087        | (IsZero, SOME zero_thm) =>
2088            @{thm expands_to_zero_exponent_imp_filterlim(2)} OF
2089              [thm, get_basis_wf_thm basis, zero_thm,
2090                 go (tl_basis basis) (thm RS @{thm expands_to_hd''})]
2091        | _ => err ()
2092    val lim_thm = go basis thm
2093  in
2094    add_sign_thm OF [lim_thm, get_basis_wf_thm basis, thm, trimmed_thm]
2095  end
2096
2097val prove_at_0 = prove_at_aux Simple_Trim
2098val prove_at_left_0 = prove_at_aux Neg_Trim
2099val prove_at_right_0 = prove_at_aux Pos_Trim
2100
2101
2102fun prove_nhds ectxt (thm, basis) =
2103  let
2104    fun simplify (a, b, c) = (a, simplify_expansion ectxt b, c)
2105    fun go thm basis =
2106      if fastype_of (get_expansion thm) = \<^typ>\<open>real\<close> then
2107        @{thm expands_to_real_imp_filterlim} OF [thm]
2108      else
2109        case whnf_expansion ectxt thm |> simplify of
2110          (NONE, thm, _) => @{thm expands_to_MSLNil_imp_filterlim} OF [thm]
2111        | (SOME _, thm, _) => (
2112            case zeroness_oracle true (SOME Sgn_Trim) ectxt (get_exponent (get_expansion thm)) of
2113              (IsZero, SOME zero_thm) =>
2114                @{thm expands_to_zero_exponent_imp_filterlim(2)} OF
2115                  [thm, get_basis_wf_thm basis, zero_thm,
2116                    go (thm RS @{thm expands_to_hd''}) (tl_basis basis)]
2117            | (IsNeg, SOME neg_thm) =>
2118                @{thm expands_to_neg_exponent_imp_filterlim} OF
2119                  [thm, get_basis_wf_thm basis, neg_thm]
2120            | (IsPos, _) =>
2121                go (try_drop_leading_term ectxt thm) basis
2122            | _ => raise TERM ("Unexpected zeroness result in prove_nhds",
2123                     [get_exponent (get_expansion thm)]))
2124  in
2125    go thm basis
2126  end
2127
2128fun prove_equivalent theta ectxt (thm1, thm2, basis) =
2129  let
2130    val ((thm1, _, SOME trimmed_thm1), (thm2, _, SOME trimmed_thm2)) =
2131      apply2 (trim_expansion true (SOME Simple_Trim) ectxt) ((thm1, basis), (thm2, basis))
2132    val pat = ConsPat (\<^const_name>\<open>Pair\<close>, [ConsPat (\<^const_name>\<open>Lazy_Eval.cmp_result.EQ\<close>, []), 
2133                ConsPat (\<^const_name>\<open>Pair\<close>, [AnyPat ("_", 0), AnyPat ("_", 0)])])
2134    val (exp1, exp2) = apply2 get_expansion (thm1, thm2)
2135    val T = fastype_of exp1
2136    val t = mk_compare_expansions_const T $ exp1 $ exp2
2137    fun eq_thm conv = HOLogic.mk_obj_eq (conv (Thm.cterm_of (get_ctxt ectxt) t))
2138    val imp_thm =
2139      if theta then @{thm compare_expansions_EQ_imp_bigtheta}
2140      else @{thm compare_expansions_EQ_same}
2141  in
2142    case match ectxt pat t (SOME []) of
2143      (SOME _, t, conv) =>
2144        let
2145          val [_, c1, c2] = HOLogic.strip_tuple t
2146          val c12_thm = if theta then [] else [the (try_prove_real_eq true ectxt (c1, c2))]
2147        in
2148          imp_thm OF ([eq_thm conv, trimmed_thm1, trimmed_thm2, thm1, thm2, get_basis_wf_thm basis] 
2149            @ c12_thm)
2150        end
2151    | _ => raise TERM ("prove_equivalent", map get_expanded_fun [thm1, thm2])
2152  end
2153
2154val prove_bigtheta = prove_equivalent true
2155val prove_asymp_equiv = prove_equivalent false
2156
2157fun print_trimming_error s ectxt exp =
2158  let
2159    val c = get_coeff exp
2160    val t = if fastype_of c = \<^typ>\<open>real\<close> then c else get_eval c
2161  in
2162    if #verbose (#ctxt ectxt) then
2163      let
2164        val ctxt = get_ctxt ectxt
2165        val p = Pretty.str "real_asymp failed to show zeroness of the following expression:"
2166        val p = Pretty.chunks [p, Pretty.indent 2 (Syntax.pretty_term ctxt t)]
2167        val _ = Pretty.writeln p
2168      in
2169        raise TERM (s, [t])
2170      end
2171    else
2172      raise TERM (s, [t])
2173  end
2174
2175fun prove_smallo ectxt (thm1, thm2, basis) =
2176  let
2177    val (thm2, _, SOME trimmed_thm) = trim_expansion true (SOME Simple_Trim) ectxt (thm2, basis)
2178    val es = get_exponents (get_expansion thm2)
2179  in
2180    case trim_expansion_while_greater true (SOME es) false NONE ectxt (thm1, basis) of
2181      (thm1, Aborted LESS, thms) =>
2182        @{thm compare_expansions_LT} OF [prove_compare_expansions LESS (map snd thms),
2183          trimmed_thm, thm1, thm2, get_basis_wf_thm basis]
2184    | (thm1, Aborted _, _) =>
2185        print_trimming_error "prove_smallo" ectxt (get_expansion thm1)
2186    | (thm1, Trimmed _, _) =>
2187        print_trimming_error "prove_smallo" ectxt (get_expansion thm1)
2188  end
2189
2190fun prove_bigo ectxt (thm1, thm2, basis) =
2191  let
2192    val (thm2, _, SOME trimmed_thm) = trim_expansion true (SOME Simple_Trim) ectxt (thm2, basis)
2193    val es = get_exponents (get_expansion thm2)
2194  in
2195    case trim_expansion_while_greater false (SOME es) false NONE ectxt (thm1, basis) of
2196      (thm1, Aborted LESS, thms) =>
2197        @{thm landau_o.small_imp_big[OF compare_expansions_LT]} OF
2198          [prove_compare_expansions LESS (map snd thms), trimmed_thm, thm1, thm2,
2199           get_basis_wf_thm basis]
2200    | (thm1, Aborted EQ, thms) =>
2201        @{thm compare_expansions_EQ_imp_bigo} OF [prove_compare_expansions EQ (map snd thms),
2202          trimmed_thm, thm1, thm2, get_basis_wf_thm basis]
2203    | (thm1, Trimmed _, _) =>
2204        print_trimming_error "prove_bigo" ectxt (get_expansion thm1)
2205  end
2206
2207
2208fun prove_asymptotic_relation_aux mode f ectxt (thm1, thm2, basis) = f (
2209  let
2210    val thm = @{thm expands_to_minus} OF [get_basis_wf_thm basis, thm1, thm2]
2211  in
2212    case ev_zeroness_oracle ectxt (get_expanded_fun thm) of
2213      SOME zero_thm => (EQUAL, zero_thm RS @{thm eventually_diff_zero_imp_eq})
2214    | _ => (
2215      case trim_expansion true (SOME mode) ectxt (thm, basis) of
2216        (thm, IsPos, SOME pos_thm) =>
2217          (GREATER, @{thm expands_to_imp_eventually_gt} OF [get_basis_wf_thm basis, thm, pos_thm])
2218      | (thm, IsNeg, SOME neg_thm) =>
2219          (LESS, @{thm expands_to_imp_eventually_lt} OF [get_basis_wf_thm basis, thm, neg_thm])
2220      | _ => raise TERM ("Unexpected zeroness result in prove_asymptotic_relation", []))
2221  end)
2222
2223val prove_eventually_greater = prove_asymptotic_relation_aux Pos_Trim snd
2224val prove_eventually_less = prove_asymptotic_relation_aux Neg_Trim snd
2225val prove_asymptotic_relation = prove_asymptotic_relation_aux Sgn_Trim I
2226
2227fun prove_eventually_nonzero ectxt (thm, basis) =
2228  case trim_expansion true (SOME Simple_Trim) ectxt (thm, basis) of
2229    (thm, _, SOME trimmed_thm) =>
2230      @{thm expands_to_imp_eventually_nz} OF [get_basis_wf_thm basis, thm, trimmed_thm]
2231  | _ => raise TERM ("prove_eventually_nonzero", [get_expanded_fun thm])
2232
2233fun extract_terms (n, strict) ectxt basis t =
2234  let
2235    val bs = get_basis_list basis
2236    fun mk_constfun c = (Abs ("x", \<^typ>\<open>real\<close>, c))
2237    val const_0 = mk_constfun \<^term>\<open>0 :: real\<close>
2238    val const_1 = mk_constfun \<^term>\<open>1 :: real\<close>
2239    fun uminus t = Term.betapply (\<^term>\<open>\<lambda>(f::real\<Rightarrow>real) x. -f x\<close>, t)
2240    fun betapply2 a b c = Term.betapply (Term.betapply (a, b), c)
2241
2242    fun mk_sum' [] acc = acc
2243      | mk_sum' ((t, sgn) :: ts) acc = mk_sum' ts (
2244          if sgn then
2245            betapply2 \<^term>\<open>%(f::real=>real) g x. f x - g x\<close> acc t
2246          else
2247            betapply2 \<^term>\<open>%(f::real=>real) g x. f x + g x\<close> acc t)
2248    fun mk_sum [] = const_0
2249      | mk_sum ((t, sgn) :: ts) = mk_sum' ts (if sgn then uminus t else t) 
2250
2251    fun mk_mult a b =
2252      if a aconv const_0 then
2253        const_0
2254      else if b aconv const_0 then
2255        const_0
2256      else if a aconv \<^term>\<open>\<lambda>_::real. 1 :: real\<close> then
2257        b
2258      else if b aconv \<^term>\<open>\<lambda>_::real. 1 :: real\<close> then
2259        a
2260      else if a aconv \<^term>\<open>\<lambda>_::real. -1 :: real\<close> then
2261        Term.betapply (\<^term>\<open>\<lambda>(f::real\<Rightarrow>real) x. -f x\<close>, b)
2262      else if b aconv \<^term>\<open>\<lambda>_::real. -1 :: real\<close> then
2263        Term.betapply (\<^term>\<open>\<lambda>(f::real\<Rightarrow>real) x. -f x\<close>, a)
2264      else
2265        Abs ("x", \<^typ>\<open>real\<close>, \<^term>\<open>(*) :: real => _\<close> $
2266          (Term.betapply (a, Bound 0)) $ (Term.betapply (b, Bound 0)))
2267
2268    fun mk_powr b e =
2269      if e = \<^term>\<open>0 :: real\<close> then
2270        const_1
2271      else
2272        let
2273          val n = HOLogic.dest_number e |> snd
2274        in
2275          if n >= 0 then
2276            Term.betapply (Term.betapply (\<^term>\<open>%(b::real=>real) e x. b x ^ e\<close>, b),
2277              HOLogic.mk_number \<^typ>\<open>nat\<close> n)
2278          else
2279            Term.betapply (Term.betapply (\<^term>\<open>%(b::real=>real) e x. b x powr e\<close>, b), e)
2280        end
2281      handle TERM _ =>
2282        Term.betapply (Term.betapply (\<^term>\<open>%(b::real=>real) e x. b x powr e\<close>, b), e)
2283
2284    fun mk_scale_elem b e acc =
2285      let
2286        val e = simplify_term ectxt e
2287      in
2288        if e = \<^term>\<open>0 :: real\<close> then
2289          acc
2290        else if e = \<^term>\<open>1 :: real\<close> then
2291          mk_mult acc b
2292        else
2293          mk_mult acc (mk_powr b e)
2294      end
2295
2296    fun mk_scale_elems [] _ acc = acc
2297      | mk_scale_elems (b :: bs) (e :: es) acc =
2298          mk_scale_elems bs es (mk_scale_elem b e acc)
2299      | mk_scale_elems _ _ _ = raise Match
2300
2301    fun mk_summand c es =
2302      let
2303        val es = mk_scale_elems bs es \<^term>\<open>\<lambda>_::real. 1 :: real\<close>
2304      in
2305        case c of
2306          Const (\<^const_name>\<open>uminus\<close>, _) $ c => ((c, true), es)
2307        | _ => ((c, false), es)
2308      end
2309
2310    fun go _ _ _ acc 0 = (acc, 0)
2311      | go 0 es t acc n =
2312          let
2313            val c = simplify_term ectxt t
2314          in
2315            if strict andalso c = \<^term>\<open>0 :: real\<close> then
2316              (acc, n)
2317            else
2318              (mk_summand c (rev es) :: acc, n - 1)
2319          end
2320      | go m es t acc n =
2321          case Lazy_Eval.whnf ectxt t |> fst of
2322            Const (\<^const_name>\<open>MS\<close>, _) $ cs $ _ =>
2323              go' m es (simplify_term ectxt cs) acc n
2324          | _ => raise TERM("extract_terms", [t])
2325    and go' _ _ _ acc 0 = (acc, 0)
2326      | go' m es cs acc n =
2327          case Lazy_Eval.whnf ectxt cs |> fst of
2328            Const (\<^const_name>\<open>MSLNil\<close>, _) => (acc, n)
2329          | Const (\<^const_name>\<open>MSLCons\<close>, _) $ c $ cs => (
2330              case Lazy_Eval.whnf ectxt c |> fst |> HOLogic.dest_prod of
2331                (c, e) =>
2332                  case go (m - 1) (e :: es) c acc n of
2333                   (acc, n) => go' m es (simplify_term ectxt cs) acc n)
2334          | _ => raise TERM("extract_terms", [t])
2335    val (summands, remaining) = go (basis_size basis) [] t [] (n + 1)
2336    val (summands, error) =
2337      if remaining = 0 then (rev (tl summands), SOME (snd (hd summands))) else (rev summands, NONE)
2338    val summands = map (fn ((c, sgn), es) => (mk_mult (mk_constfun c) es, sgn)) summands
2339    val error = Option.map (fn err => Term.betapply (\<^term>\<open>\<lambda>f::real\<Rightarrow>real. O(f)\<close>, err)) error
2340    val expansion = mk_sum summands 
2341  in
2342    (expansion, error)
2343  end
2344
2345end
2346
2347
2348structure Multiseries_Expansion_Basic : EXPANSION_INTERFACE =
2349struct
2350open Multiseries_Expansion;
2351
2352type T = expansion_thm
2353
2354val expand_term = expand_term
2355val expand_terms = expand_terms
2356
2357val prove_nhds = prove_nhds
2358val prove_at_infinity = prove_at_infinity
2359val prove_at_top = prove_at_top
2360val prove_at_bot = prove_at_bot
2361val prove_at_0 = prove_at_0
2362val prove_at_right_0 = prove_at_right_0
2363val prove_at_left_0 = prove_at_left_0
2364val prove_eventually_nonzero = prove_eventually_nonzero
2365
2366val prove_eventually_less = prove_eventually_less
2367val prove_eventually_greater = prove_eventually_greater
2368
2369val prove_smallo = prove_smallo
2370val prove_bigo = prove_bigo
2371val prove_bigtheta = prove_bigtheta
2372val prove_asymp_equiv = prove_asymp_equiv
2373
2374end
2375