src/HOL/Real_Asymp/Real_Asymp_Approx.thy
changeset 69597 ff784d5a5bfb
parent 69064 5840724b1d71
equal deleted inserted replaced
69596:c8a2755bf220 69597:ff784d5a5bfb
     8 theory Real_Asymp_Approx
     8 theory Real_Asymp_Approx
     9   imports Real_Asymp "HOL-Decision_Procs.Approximation"
     9   imports Real_Asymp "HOL-Decision_Procs.Approximation"
    10 begin
    10 begin
    11 
    11 
    12 text \<open>
    12 text \<open>
    13   For large enough constants (such as @{term "exp (exp 10000 :: real)"}), the 
    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
    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.
    15   terminating and causing the entire prover environment to crash.
    16 
    16 
    17   To avoid causing such situations, we first check the plausibility of the fact to prove using
    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
    18   floating-point arithmetic and only run the approximation method if the floating point evaluation
    30 
    30 
    31 structure Real_Asymp_Approx : REAL_ASYMP_APPROX =
    31 structure Real_Asymp_Approx : REAL_ASYMP_APPROX =
    32 struct
    32 struct
    33 
    33 
    34 val real_asymp_approx =
    34 val real_asymp_approx =
    35   Attrib.setup_config_bool @{binding real_asymp_approx} (K true)
    35   Attrib.setup_config_bool \<^binding>\<open>real_asymp_approx\<close> (K true)
    36 
    36 
    37 val nan = Real.fromInt 0 / Real.fromInt 0
    37 val nan = Real.fromInt 0 / Real.fromInt 0
    38 fun clamp n = if n < 0 then 0 else n
    38 fun clamp n = if n < 0 then 0 else n
    39 
    39 
    40 fun eval_nat (@{term "(+) :: nat => _"} $ a $ b) = bin (op +) (a, b)
    40 fun eval_nat (\<^term>\<open>(+) :: nat => _\<close> $ a $ b) = bin (op +) (a, b)
    41   | eval_nat (@{term "(-) :: nat => _"} $ a $ b) = bin (clamp o op -) (a, b)
    41   | eval_nat (\<^term>\<open>(-) :: nat => _\<close> $ a $ b) = bin (clamp o op -) (a, b)
    42   | eval_nat (@{term "(*) :: nat => _"} $ a $ b) = bin (op *) (a, b)
    42   | eval_nat (\<^term>\<open>(*) :: nat => _\<close> $ a $ b) = bin (op *) (a, b)
    43   | eval_nat (@{term "(div) :: nat => _"} $ a $ b) = bin Int.div (a, b)
    43   | eval_nat (\<^term>\<open>(div) :: nat => _\<close> $ a $ b) = bin Int.div (a, b)
    44   | eval_nat (@{term "(^) :: nat => _"} $ a $ b) = bin (fn (a,b) => Integer.pow a b) (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 "numeral :: _ => nat"} $ _)) = snd (HOLogic.dest_number t)
    45   | eval_nat (t as (\<^term>\<open>numeral :: _ => nat\<close> $ _)) = snd (HOLogic.dest_number t)
    46   | eval_nat (@{term "0 :: nat"}) = 0
    46   | eval_nat (\<^term>\<open>0 :: nat\<close>) = 0
    47   | eval_nat (@{term "1 :: nat"}) = 1
    47   | eval_nat (\<^term>\<open>1 :: nat\<close>) = 1
    48   | eval_nat (@{term "Nat.Suc"} $ a) = un (fn n => n + 1) a
    48   | eval_nat (\<^term>\<open>Nat.Suc\<close> $ a) = un (fn n => n + 1) a
    49   | eval_nat _ = ~1
    49   | eval_nat _ = ~1
    50 and un f a =
    50 and un f a =
    51       let
    51       let
    52         val a = eval_nat a 
    52         val a = eval_nat a 
    53       in
    53       in
    61       end
    61       end
    62 
    62 
    63 fun sgn n =
    63 fun sgn n =
    64   if n < Real.fromInt 0 then Real.fromInt (~1) else Real.fromInt 1
    64   if n < Real.fromInt 0 then Real.fromInt (~1) else Real.fromInt 1
    65 
    65 
    66 fun eval (@{term "(+) :: real => _"} $ a $ b) = eval a + eval b
    66 fun eval (\<^term>\<open>(+) :: real => _\<close> $ a $ b) = eval a + eval b
    67   | eval (@{term "(-) :: real => _"} $ a $ b) = eval a - eval b
    67   | eval (\<^term>\<open>(-) :: real => _\<close> $ a $ b) = eval a - eval b
    68   | eval (@{term "(*) :: real => _"} $ a $ b) = eval a * eval b
    68   | eval (\<^term>\<open>(*) :: real => _\<close> $ a $ b) = eval a * eval b
    69   | eval (@{term "(/) :: real => _"} $ a $ b) =
    69   | eval (\<^term>\<open>(/) :: real => _\<close> $ a $ b) =
    70       let val a = eval a; val b = eval b in
    70       let val a = eval a; val b = eval b in
    71         if Real.==(b, Real.fromInt 0) then nan else a / b
    71         if Real.==(b, Real.fromInt 0) then nan else a / b
    72       end
    72       end
    73   | eval (@{term "inverse :: real => _"} $ a) = Real.fromInt 1 / eval a
    73   | eval (\<^term>\<open>inverse :: real => _\<close> $ a) = Real.fromInt 1 / eval a
    74   | eval (@{term "uminus :: real => _"} $ a) = Real.~ (eval a)
    74   | eval (\<^term>\<open>uminus :: real => _\<close> $ a) = Real.~ (eval a)
    75   | eval (@{term "exp :: real => _"} $ a) = Math.exp (eval a)
    75   | eval (\<^term>\<open>exp :: real => _\<close> $ a) = Math.exp (eval a)
    76   | eval (@{term "ln :: real => _"} $ 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
    77       let val a = eval a in if a > Real.fromInt 0 then Math.ln a else nan end
    78   | eval (@{term "(powr) :: real => _"} $ a $ b) =
    78   | eval (\<^term>\<open>(powr) :: real => _\<close> $ a $ b) =
    79       let
    79       let
    80         val a = eval a; val b = eval b
    80         val a = eval a; val b = eval b
    81       in
    81       in
    82         if a < Real.fromInt 0 orelse not (Real.isFinite a) orelse not (Real.isFinite b) then
    82         if a < Real.fromInt 0 orelse not (Real.isFinite a) orelse not (Real.isFinite b) then
    83           nan
    83           nan
    84         else if Real.==(a, Real.fromInt 0) then
    84         else if Real.==(a, Real.fromInt 0) then
    85           Real.fromInt 0
    85           Real.fromInt 0
    86         else 
    86         else 
    87           Math.pow (a, b)
    87           Math.pow (a, b)
    88       end
    88       end
    89   | eval (@{term "(^) :: real => _"} $ a $ b) =
    89   | eval (\<^term>\<open>(^) :: real => _\<close> $ a $ b) =
    90       let
    90       let
    91         fun powr x y = 
    91         fun powr x y = 
    92           if not (Real.isFinite x) orelse y < 0 then
    92           if not (Real.isFinite x) orelse y < 0 then
    93             nan
    93             nan
    94           else if x < Real.fromInt 0 andalso y mod 2 = 1 then 
    94           else if x < Real.fromInt 0 andalso y mod 2 = 1 then 
    96           else
    96           else
    97             Math.pow (Real.abs x, Real.fromInt y)
    97             Math.pow (Real.abs x, Real.fromInt y)
    98       in
    98       in
    99         powr (eval a) (eval_nat b)
    99         powr (eval a) (eval_nat b)
   100       end
   100       end
   101   | eval (@{term "root :: nat => real => _"} $ n $ a) =
   101   | eval (\<^term>\<open>root :: nat => real => _\<close> $ n $ a) =
   102       let val a = eval a; val n = eval_nat n in
   102       let val a = eval a; val n = eval_nat n in
   103         if n = 0 then Real.fromInt 0
   103         if n = 0 then Real.fromInt 0
   104           else sgn a * Math.pow (Real.abs a, Real.fromInt 1 / Real.fromInt n) end
   104           else sgn a * Math.pow (Real.abs a, Real.fromInt 1 / Real.fromInt n) end
   105   | eval (@{term "sqrt :: real => _"} $ a) =
   105   | eval (\<^term>\<open>sqrt :: real => _\<close> $ a) =
   106       let val a = eval a in sgn a * Math.sqrt (abs a) end
   106       let val a = eval a in sgn a * Math.sqrt (abs a) end
   107   | eval (@{term "log :: real => _"} $ a $ b) =
   107   | eval (\<^term>\<open>log :: real => _\<close> $ a $ b) =
   108       let
   108       let
   109         val (a, b) = apply2 eval (a, b)
   109         val (a, b) = apply2 eval (a, b)
   110       in
   110       in
   111         if b > Real.fromInt 0 andalso a > Real.fromInt 0 andalso Real.!= (a, Real.fromInt 1) then
   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
   112           Math.ln b / Math.ln a
   113         else
   113         else
   114           nan
   114           nan
   115       end
   115       end
   116   | eval (t as (@{term "numeral :: _ => real"} $ _)) =
   116   | eval (t as (\<^term>\<open>numeral :: _ => real\<close> $ _)) =
   117       Real.fromInt (snd (HOLogic.dest_number t))
   117       Real.fromInt (snd (HOLogic.dest_number t))
   118   | eval (@{term "0 :: real"}) = Real.fromInt 0
   118   | eval (\<^term>\<open>0 :: real\<close>) = Real.fromInt 0
   119   | eval (@{term "1 :: real"}) = Real.fromInt 1
   119   | eval (\<^term>\<open>1 :: real\<close>) = Real.fromInt 1
   120   | eval _ = nan
   120   | eval _ = nan
   121 
   121 
   122 fun sign_oracle_tac ctxt i =
   122 fun sign_oracle_tac ctxt i =
   123   let
   123   let
   124     fun tac {context = ctxt, concl, ...} =
   124     fun tac {context = ctxt, concl, ...} =
   125       let
   125       let
   126         val b =
   126         val b =
   127           case HOLogic.dest_Trueprop (Thm.term_of concl) of
   127           case HOLogic.dest_Trueprop (Thm.term_of concl) of
   128             @{term "(<) :: real \<Rightarrow> _"} $ lhs $ rhs =>
   128             \<^term>\<open>(<) :: real \<Rightarrow> _\<close> $ lhs $ rhs =>
   129               let
   129               let
   130                 val (x, y) = apply2 eval (lhs, rhs)
   130                 val (x, y) = apply2 eval (lhs, rhs)
   131               in
   131               in
   132                 Real.isFinite x andalso Real.isFinite y andalso x < y
   132                 Real.isFinite x andalso Real.isFinite y andalso x < y
   133               end
   133               end
   134           | @{term "HOL.Not"} $ (@{term "(=) :: real \<Rightarrow> _"} $ lhs $ rhs) =>
   134           | \<^term>\<open>HOL.Not\<close> $ (\<^term>\<open>(=) :: real \<Rightarrow> _\<close> $ lhs $ rhs) =>
   135               let
   135               let
   136                 val (x, y) = apply2 eval (lhs, rhs)
   136                 val (x, y) = apply2 eval (lhs, rhs)
   137               in
   137               in
   138                 Real.isFinite x andalso Real.isFinite y andalso Real.!= (x, y)
   138                 Real.isFinite x andalso Real.isFinite y andalso Real.!= (x, y)
   139               end
   139               end
   152 \<close>
   152 \<close>
   153 
   153 
   154 setup \<open>
   154 setup \<open>
   155   Context.theory_map (
   155   Context.theory_map (
   156     Multiseries_Expansion.register_sign_oracle
   156     Multiseries_Expansion.register_sign_oracle
   157       (@{binding approximation_tac}, Real_Asymp_Approx.sign_oracle_tac))
   157       (\<^binding>\<open>approximation_tac\<close>, Real_Asymp_Approx.sign_oracle_tac))
   158 \<close>
   158 \<close>
   159 
   159 
   160 lemma "filterlim (\<lambda>n. (1 + (2 / 3 :: real) ^ (n + 1)) ^ 2 ^ n / 2 powr (4 / 3) ^ (n - 1))
   160 lemma "filterlim (\<lambda>n. (1 + (2 / 3 :: real) ^ (n + 1)) ^ 2 ^ n / 2 powr (4 / 3) ^ (n - 1))
   161          at_top at_top"
   161          at_top at_top"
   162 proof -
   162 proof -