1(*
2  File:     Real_Asymp_Approx.thy
3  Author:   Manuel Eberl, TU M��nchen
4  
5  Integrates the approximation method into the real_asymp framework as a sign oracle
6  to automatically prove positivity/negativity of real constants
7*)
8theory Real_Asymp_Approx
9  imports Real_Asymp "HOL-Decision_Procs.Approximation"
10begin
11
12text \<open>
13  For large enough constants (such as \<^term>\<open>exp (exp 10000 :: real)\<close>), the 
14  @{method approximation} method can require a huge amount of time and memory, effectively not
15  terminating and causing the entire prover environment to crash.
16
17  To avoid causing such situations, we first check the plausibility of the fact to prove using
18  floating-point arithmetic and only run the approximation method if the floating point evaluation
19  supports the claim. In particular, exorbitantly large numbers like the one above will lead to
20  floating point overflow, which makes the check fail.
21\<close>
22
23ML \<open>
24signature REAL_ASYMP_APPROX = 
25sig
26  val real_asymp_approx : bool Config.T
27  val sign_oracle_tac : Proof.context -> int -> tactic
28  val eval : term -> real
29end
30
31structure Real_Asymp_Approx : REAL_ASYMP_APPROX =
32struct
33
34val real_asymp_approx =
35  Attrib.setup_config_bool \<^binding>\<open>real_asymp_approx\<close> (K true)
36
37val nan = Real.fromInt 0 / Real.fromInt 0
38fun clamp n = if n < 0 then 0 else n
39
40fun eval_nat (\<^term>\<open>(+) :: nat => _\<close> $ a $ b) = bin (op +) (a, b)
41  | eval_nat (\<^term>\<open>(-) :: nat => _\<close> $ a $ b) = bin (clamp o op -) (a, b)
42  | eval_nat (\<^term>\<open>(*) :: nat => _\<close> $ a $ b) = bin (op *) (a, b)
43  | eval_nat (\<^term>\<open>(div) :: nat => _\<close> $ a $ b) = bin Int.div (a, b)
44  | eval_nat (\<^term>\<open>(^) :: nat => _\<close> $ a $ b) = bin (fn (a,b) => Integer.pow a b) (a, b)
45  | eval_nat (t as (\<^term>\<open>numeral :: _ => nat\<close> $ _)) = snd (HOLogic.dest_number t)
46  | eval_nat (\<^term>\<open>0 :: nat\<close>) = 0
47  | eval_nat (\<^term>\<open>1 :: nat\<close>) = 1
48  | eval_nat (\<^term>\<open>Nat.Suc\<close> $ a) = un (fn n => n + 1) a
49  | eval_nat _ = ~1
50and un f a =
51      let
52        val a = eval_nat a 
53      in
54        if a >= 0 then f a else ~1
55      end
56and bin f (a, b) =
57      let
58        val (a, b) = apply2 eval_nat (a, b) 
59      in
60        if a >= 0 andalso b >= 0 then f (a, b) else ~1
61      end
62
63fun sgn n =
64  if n < Real.fromInt 0 then Real.fromInt (~1) else Real.fromInt 1
65
66fun eval (\<^term>\<open>(+) :: real => _\<close> $ a $ b) = eval a + eval b
67  | eval (\<^term>\<open>(-) :: real => _\<close> $ a $ b) = eval a - eval b
68  | eval (\<^term>\<open>(*) :: real => _\<close> $ a $ b) = eval a * eval b
69  | eval (\<^term>\<open>(/) :: real => _\<close> $ a $ b) =
70      let val a = eval a; val b = eval b in
71        if Real.==(b, Real.fromInt 0) then nan else a / b
72      end
73  | eval (\<^term>\<open>inverse :: real => _\<close> $ a) = Real.fromInt 1 / eval a
74  | eval (\<^term>\<open>uminus :: real => _\<close> $ a) = Real.~ (eval a)
75  | eval (\<^term>\<open>exp :: real => _\<close> $ a) = Math.exp (eval a)
76  | eval (\<^term>\<open>ln :: real => _\<close> $ a) =
77      let val a = eval a in if a > Real.fromInt 0 then Math.ln a else nan end
78  | eval (\<^term>\<open>(powr) :: real => _\<close> $ a $ b) =
79      let
80        val a = eval a; val b = eval b
81      in
82        if a < Real.fromInt 0 orelse not (Real.isFinite a) orelse not (Real.isFinite b) then
83          nan
84        else if Real.==(a, Real.fromInt 0) then
85          Real.fromInt 0
86        else 
87          Math.pow (a, b)
88      end
89  | eval (\<^term>\<open>(^) :: real => _\<close> $ a $ b) =
90      let
91        fun powr x y = 
92          if not (Real.isFinite x) orelse y < 0 then
93            nan
94          else if x < Real.fromInt 0 andalso y mod 2 = 1 then 
95            Real.~ (Math.pow (Real.~ x, Real.fromInt y))
96          else
97            Math.pow (Real.abs x, Real.fromInt y)
98      in
99        powr (eval a) (eval_nat b)
100      end
101  | eval (\<^term>\<open>root :: nat => real => _\<close> $ n $ a) =
102      let val a = eval a; val n = eval_nat n in
103        if n = 0 then Real.fromInt 0
104          else sgn a * Math.pow (Real.abs a, Real.fromInt 1 / Real.fromInt n) end
105  | eval (\<^term>\<open>sqrt :: real => _\<close> $ a) =
106      let val a = eval a in sgn a * Math.sqrt (abs a) end
107  | eval (\<^term>\<open>log :: real => _\<close> $ a $ b) =
108      let
109        val (a, b) = apply2 eval (a, b)
110      in
111        if b > Real.fromInt 0 andalso a > Real.fromInt 0 andalso Real.!= (a, Real.fromInt 1) then
112          Math.ln b / Math.ln a
113        else
114          nan
115      end
116  | eval (t as (\<^term>\<open>numeral :: _ => real\<close> $ _)) =
117      Real.fromInt (snd (HOLogic.dest_number t))
118  | eval (\<^term>\<open>0 :: real\<close>) = Real.fromInt 0
119  | eval (\<^term>\<open>1 :: real\<close>) = Real.fromInt 1
120  | eval _ = nan
121
122fun sign_oracle_tac ctxt i =
123  let
124    fun tac {context = ctxt, concl, ...} =
125      let
126        val b =
127          case HOLogic.dest_Trueprop (Thm.term_of concl) of
128            \<^term>\<open>(<) :: real \<Rightarrow> _\<close> $ lhs $ rhs =>
129              let
130                val (x, y) = apply2 eval (lhs, rhs)
131              in
132                Real.isFinite x andalso Real.isFinite y andalso x < y
133              end
134          | \<^term>\<open>HOL.Not\<close> $ (\<^term>\<open>(=) :: real \<Rightarrow> _\<close> $ lhs $ rhs) =>
135              let
136                val (x, y) = apply2 eval (lhs, rhs)
137              in
138                Real.isFinite x andalso Real.isFinite y andalso Real.!= (x, y)
139              end
140          | _ => false
141     in
142       if b then HEADGOAL (Approximation.approximation_tac 10 [] NONE ctxt) else no_tac
143     end
144  in
145    if Config.get ctxt real_asymp_approx then
146      Subgoal.FOCUS tac ctxt i
147    else
148      no_tac
149  end
150
151end
152\<close>
153
154setup \<open>
155  Context.theory_map (
156    Multiseries_Expansion.register_sign_oracle
157      (\<^binding>\<open>approximation_tac\<close>, Real_Asymp_Approx.sign_oracle_tac))
158\<close>
159
160lemma "filterlim (\<lambda>n. (1 + (2 / 3 :: real) ^ (n + 1)) ^ 2 ^ n / 2 powr (4 / 3) ^ (n - 1))
161         at_top at_top"
162proof -
163  have [simp]: "ln 4 = 2 * ln (2 :: real)"
164    using ln_realpow[of 2 2] by simp
165  show ?thesis by (real_asymp simp: field_simps ln_div)
166qed
167
168end