Implemented taylor series expansion for approximation
authorhoelzl
Mon Jun 29 23:29:11 2009 +0200 (2009-06-29)
changeset 31863e391eee8bf14
parent 31861 1bb5fe96f61e
child 31864 90ec13cf7ab0
Implemented taylor series expansion for approximation
NEWS
src/HOL/Decision_Procs/Approximation.thy
src/HOL/Decision_Procs/ex/Approximation_Ex.thy
src/HOL/Library/Float.thy
     1.1 --- a/NEWS	Tue Jun 30 00:57:24 2009 +0200
     1.2 +++ b/NEWS	Mon Jun 29 23:29:11 2009 +0200
     1.3 @@ -73,7 +73,7 @@
     1.4  approximation method. 
     1.5  
     1.6  * "approximate" supports now arithmetic expressions as boundaries of intervals and implements
     1.7 -interval splitting.
     1.8 +interval splitting and taylor series expansion.
     1.9  
    1.10  
    1.11  *** ML ***
     2.1 --- a/src/HOL/Decision_Procs/Approximation.thy	Tue Jun 30 00:57:24 2009 +0200
     2.2 +++ b/src/HOL/Decision_Procs/Approximation.thy	Mon Jun 29 23:29:11 2009 +0200
     2.3 @@ -2069,8 +2069,7 @@
     2.4    | Atom nat
     2.5    | Num float
     2.6  
     2.7 -fun interpret_floatarith :: "floatarith \<Rightarrow> real list \<Rightarrow> real"
     2.8 -where
     2.9 +fun interpret_floatarith :: "floatarith \<Rightarrow> real list \<Rightarrow> real" where
    2.10  "interpret_floatarith (Add a b) vs   = (interpret_floatarith a vs) + (interpret_floatarith b vs)" |
    2.11  "interpret_floatarith (Minus a) vs    = - (interpret_floatarith a vs)" |
    2.12  "interpret_floatarith (Mult a b) vs   = (interpret_floatarith a vs) * (interpret_floatarith b vs)" |
    2.13 @@ -2117,7 +2116,6 @@
    2.14    and "interpret_floatarith (Num (Float 1 0)) vs = 1"
    2.15    and "interpret_floatarith (Num (Float (number_of a) 0)) vs = number_of a" by auto
    2.16  
    2.17 -
    2.18  subsection "Implement approximation function"
    2.19  
    2.20  fun lift_bin' :: "(float * float) option \<Rightarrow> (float * float) option \<Rightarrow> (float \<Rightarrow> float \<Rightarrow> float \<Rightarrow> float \<Rightarrow> (float * float)) \<Rightarrow> (float * float) option" where
    2.21 @@ -2632,6 +2630,576 @@
    2.22    shows "interpret_form f xs"
    2.23    using approx_form_aux[OF _ bounded_by_None] assms by auto
    2.24  
    2.25 +subsection {* Implementing Taylor series expansion *}
    2.26 +
    2.27 +fun isDERIV :: "nat \<Rightarrow> floatarith \<Rightarrow> real list \<Rightarrow> bool" where
    2.28 +"isDERIV x (Add a b) vs         = (isDERIV x a vs \<and> isDERIV x b vs)" |
    2.29 +"isDERIV x (Mult a b) vs        = (isDERIV x a vs \<and> isDERIV x b vs)" |
    2.30 +"isDERIV x (Minus a) vs         = isDERIV x a vs" |
    2.31 +"isDERIV x (Inverse a) vs       = (isDERIV x a vs \<and> interpret_floatarith a vs \<noteq> 0)" |
    2.32 +"isDERIV x (Cos a) vs           = isDERIV x a vs" |
    2.33 +"isDERIV x (Arctan a) vs        = isDERIV x a vs" |
    2.34 +"isDERIV x (Min a b) vs         = False" |
    2.35 +"isDERIV x (Max a b) vs         = False" |
    2.36 +"isDERIV x (Abs a) vs           = False" |
    2.37 +"isDERIV x Pi vs                = True" |
    2.38 +"isDERIV x (Sqrt a) vs          = (isDERIV x a vs \<and> interpret_floatarith a vs > 0)" |
    2.39 +"isDERIV x (Exp a) vs           = isDERIV x a vs" |
    2.40 +"isDERIV x (Ln a) vs            = (isDERIV x a vs \<and> interpret_floatarith a vs > 0)" |
    2.41 +"isDERIV x (Power a 0) vs       = True" |
    2.42 +"isDERIV x (Power a (Suc n)) vs = isDERIV x a vs" |
    2.43 +"isDERIV x (Num f) vs           = True" |
    2.44 +"isDERIV x (Atom n) vs          = True"
    2.45 +
    2.46 +fun DERIV_floatarith :: "nat \<Rightarrow> floatarith \<Rightarrow> floatarith" where
    2.47 +"DERIV_floatarith x (Add a b)         = Add (DERIV_floatarith x a) (DERIV_floatarith x b)" |
    2.48 +"DERIV_floatarith x (Mult a b)        = Add (Mult a (DERIV_floatarith x b)) (Mult (DERIV_floatarith x a) b)" |
    2.49 +"DERIV_floatarith x (Minus a)         = Minus (DERIV_floatarith x a)" |
    2.50 +"DERIV_floatarith x (Inverse a)       = Minus (Mult (DERIV_floatarith x a) (Inverse (Power a 2)))" |
    2.51 +"DERIV_floatarith x (Cos a)           = Minus (Mult (Cos (Add (Mult Pi (Num (Float 1 -1))) (Minus a))) (DERIV_floatarith x a))" |
    2.52 +"DERIV_floatarith x (Arctan a)        = Mult (Inverse (Add (Num 1) (Power a 2))) (DERIV_floatarith x a)" |
    2.53 +"DERIV_floatarith x (Min a b)         = Num 0" |
    2.54 +"DERIV_floatarith x (Max a b)         = Num 0" |
    2.55 +"DERIV_floatarith x (Abs a)           = Num 0" |
    2.56 +"DERIV_floatarith x Pi                = Num 0" |
    2.57 +"DERIV_floatarith x (Sqrt a)          = (Mult (Inverse (Mult (Sqrt a) (Num 2))) (DERIV_floatarith x a))" |
    2.58 +"DERIV_floatarith x (Exp a)           = Mult (Exp a) (DERIV_floatarith x a)" |
    2.59 +"DERIV_floatarith x (Ln a)            = Mult (Inverse a) (DERIV_floatarith x a)" |
    2.60 +"DERIV_floatarith x (Power a 0)       = Num 0" |
    2.61 +"DERIV_floatarith x (Power a (Suc n)) = Mult (Num (Float (int (Suc n)) 0)) (Mult (Power a n) (DERIV_floatarith x a))" |
    2.62 +"DERIV_floatarith x (Num f)           = Num 0" |
    2.63 +"DERIV_floatarith x (Atom n)          = (if x = n then Num 1 else Num 0)"
    2.64 +
    2.65 +lemma DERIV_chain'': "\<lbrakk>DERIV g (f x) :> E ; DERIV f x :> D; x' = E * D \<rbrakk> \<Longrightarrow>
    2.66 +  DERIV (\<lambda>x. g (f x)) x :> x'" using DERIV_chain' by auto
    2.67 +
    2.68 +lemma DERIV_cong: "\<lbrakk> DERIV f x :> X ; X = X' \<rbrakk> \<Longrightarrow> DERIV f x :> X'" by simp
    2.69 +
    2.70 +lemma DERIV_floatarith:
    2.71 +  assumes "n < length vs"
    2.72 +  assumes isDERIV: "isDERIV n f (vs[n := x])"
    2.73 +  shows "DERIV (\<lambda> x'. interpret_floatarith f (vs[n := x'])) x :>
    2.74 +               interpret_floatarith (DERIV_floatarith n f) (vs[n := x])"
    2.75 +   (is "DERIV (?i f) x :> _")
    2.76 +using isDERIV proof (induct f arbitrary: x)
    2.77 +  case (Add a b) thus ?case by (auto intro: DERIV_add)
    2.78 +next case (Mult a b) thus ?case by (auto intro!: DERIV_mult[THEN DERIV_cong])
    2.79 +next case (Minus a) thus ?case by (auto intro!: DERIV_minus[THEN DERIV_cong])
    2.80 +next case (Inverse a) thus ?case
    2.81 +    by (auto intro!: DERIV_inverse_fun[THEN DERIV_cong] DERIV_inverse[THEN DERIV_cong]
    2.82 +             simp add: algebra_simps power2_eq_square)
    2.83 +next case (Cos a) thus ?case
    2.84 +  by (auto intro!: DERIV_chain''[of cos "?i a"]
    2.85 +                   DERIV_cos[THEN DERIV_cong]
    2.86 +           simp del: interpret_floatarith.simps(5)
    2.87 +           simp add: interpret_floatarith_sin interpret_floatarith.simps(5)[of a])
    2.88 +next case (Arctan a) thus ?case
    2.89 +  by (auto intro!: DERIV_chain''[of arctan "?i a"] DERIV_arctan[THEN DERIV_cong])
    2.90 +next case (Exp a) thus ?case
    2.91 +  by (auto intro!: DERIV_chain''[of exp "?i a"] DERIV_exp[THEN DERIV_cong])
    2.92 +next case (Power a n) thus ?case
    2.93 +  by (cases n, auto intro!: DERIV_power_Suc[THEN DERIV_cong]
    2.94 +                    simp del: power_Suc simp add: real_eq_of_nat)
    2.95 +next case (Sqrt a) thus ?case
    2.96 +    by (auto intro!: DERIV_chain''[of sqrt "?i a"] DERIV_real_sqrt[THEN DERIV_cong])
    2.97 +next case (Ln a) thus ?case
    2.98 +    by (auto intro!: DERIV_chain''[of ln "?i a"] DERIV_ln[THEN DERIV_cong]
    2.99 +                     simp add: divide_inverse)
   2.100 +next case (Atom i) thus ?case using `n < length vs` by auto
   2.101 +qed auto
   2.102 +
   2.103 +declare approx.simps[simp del]
   2.104 +
   2.105 +fun isDERIV_approx :: "nat \<Rightarrow> nat \<Rightarrow> floatarith \<Rightarrow> (float * float) option list \<Rightarrow> bool" where
   2.106 +"isDERIV_approx prec x (Add a b) vs         = (isDERIV_approx prec x a vs \<and> isDERIV_approx prec x b vs)" |
   2.107 +"isDERIV_approx prec x (Mult a b) vs        = (isDERIV_approx prec x a vs \<and> isDERIV_approx prec x b vs)" |
   2.108 +"isDERIV_approx prec x (Minus a) vs         = isDERIV_approx prec x a vs" |
   2.109 +"isDERIV_approx prec x (Inverse a) vs       =
   2.110 +  (isDERIV_approx prec x a vs \<and> (case approx prec a vs of Some (l, u) \<Rightarrow> 0 < l \<or> u < 0 | None \<Rightarrow> False))" |
   2.111 +"isDERIV_approx prec x (Cos a) vs           = isDERIV_approx prec x a vs" |
   2.112 +"isDERIV_approx prec x (Arctan a) vs        = isDERIV_approx prec x a vs" |
   2.113 +"isDERIV_approx prec x (Min a b) vs         = False" |
   2.114 +"isDERIV_approx prec x (Max a b) vs         = False" |
   2.115 +"isDERIV_approx prec x (Abs a) vs           = False" |
   2.116 +"isDERIV_approx prec x Pi vs                = True" |
   2.117 +"isDERIV_approx prec x (Sqrt a) vs          =
   2.118 +  (isDERIV_approx prec x a vs \<and> (case approx prec a vs of Some (l, u) \<Rightarrow> 0 < l | None \<Rightarrow> False))" |
   2.119 +"isDERIV_approx prec x (Exp a) vs           = isDERIV_approx prec x a vs" |
   2.120 +"isDERIV_approx prec x (Ln a) vs            =
   2.121 +  (isDERIV_approx prec x a vs \<and> (case approx prec a vs of Some (l, u) \<Rightarrow> 0 < l | None \<Rightarrow> False))" |
   2.122 +"isDERIV_approx prec x (Power a 0) vs       = True" |
   2.123 +"isDERIV_approx prec x (Power a (Suc n)) vs = isDERIV_approx prec x a vs" |
   2.124 +"isDERIV_approx prec x (Num f) vs           = True" |
   2.125 +"isDERIV_approx prec x (Atom n) vs          = True"
   2.126 +
   2.127 +lemma isDERIV_approx:
   2.128 +  assumes "bounded_by xs vs"
   2.129 +  and isDERIV_approx: "isDERIV_approx prec x f vs"
   2.130 +  shows "isDERIV x f xs"
   2.131 +using isDERIV_approx proof (induct f)
   2.132 +  case (Inverse a)
   2.133 +  then obtain l u where approx_Some: "Some (l, u) = approx prec a vs"
   2.134 +    and *: "0 < l \<or> u < 0"
   2.135 +    by (cases "approx prec a vs", auto)
   2.136 +  with approx[OF `bounded_by xs vs` approx_Some]
   2.137 +  have "interpret_floatarith a xs \<noteq> 0" unfolding less_float_def by auto
   2.138 +  thus ?case using Inverse by auto
   2.139 +next
   2.140 +  case (Ln a)
   2.141 +  then obtain l u where approx_Some: "Some (l, u) = approx prec a vs"
   2.142 +    and *: "0 < l"
   2.143 +    by (cases "approx prec a vs", auto)
   2.144 +  with approx[OF `bounded_by xs vs` approx_Some]
   2.145 +  have "0 < interpret_floatarith a xs" unfolding less_float_def by auto
   2.146 +  thus ?case using Ln by auto
   2.147 +next
   2.148 +  case (Sqrt a)
   2.149 +  then obtain l u where approx_Some: "Some (l, u) = approx prec a vs"
   2.150 +    and *: "0 < l"
   2.151 +    by (cases "approx prec a vs", auto)
   2.152 +  with approx[OF `bounded_by xs vs` approx_Some]
   2.153 +  have "0 < interpret_floatarith a xs" unfolding less_float_def by auto
   2.154 +  thus ?case using Sqrt by auto
   2.155 +next
   2.156 +  case (Power a n) thus ?case by (cases n, auto)
   2.157 +qed auto
   2.158 +
   2.159 +lemma bounded_by_update_var:
   2.160 +  assumes "bounded_by xs vs" and "vs ! i = Some (l, u)"
   2.161 +  and bnd: "x \<in> { real l .. real u }"
   2.162 +  shows "bounded_by (xs[i := x]) vs"
   2.163 +proof (cases "i < length xs")
   2.164 +  case False thus ?thesis using `bounded_by xs vs` by auto
   2.165 +next
   2.166 +  let ?xs = "xs[i := x]"
   2.167 +  case True hence "i < length ?xs" by auto
   2.168 +{ fix j
   2.169 +  assume "j < length vs"
   2.170 +  have "case vs ! j of None \<Rightarrow> True | Some (l, u) \<Rightarrow> ?xs ! j \<in> { real l .. real u }"
   2.171 +  proof (cases "vs ! j")
   2.172 +    case (Some b)
   2.173 +    thus ?thesis
   2.174 +    proof (cases "i = j")
   2.175 +      case True
   2.176 +      thus ?thesis using `vs ! i = Some (l, u)` Some and bnd `i < length ?xs`
   2.177 +	by auto
   2.178 +    next
   2.179 +      case False
   2.180 +      thus ?thesis using `bounded_by xs vs`[THEN bounded_byE, OF `j < length vs`] Some
   2.181 +	by auto
   2.182 +    qed
   2.183 +  qed auto }
   2.184 +  thus ?thesis unfolding bounded_by_def by auto
   2.185 +qed
   2.186 +
   2.187 +lemma isDERIV_approx':
   2.188 +  assumes "bounded_by xs vs"
   2.189 +  and vs_x: "vs ! x = Some (l, u)" and X_in: "X \<in> { real l .. real u }"
   2.190 +  and approx: "isDERIV_approx prec x f vs"
   2.191 +  shows "isDERIV x f (xs[x := X])"
   2.192 +proof -
   2.193 +  note bounded_by_update_var[OF `bounded_by xs vs` vs_x X_in] approx
   2.194 +  thus ?thesis by (rule isDERIV_approx)
   2.195 +qed
   2.196 +
   2.197 +lemma DERIV_approx:
   2.198 +  assumes "n < length xs" and bnd: "bounded_by xs vs"
   2.199 +  and isD: "isDERIV_approx prec n f vs"
   2.200 +  and app: "Some (l, u) = approx prec (DERIV_floatarith n f) vs" (is "_ = approx _ ?D _")
   2.201 +  shows "\<exists>x. real l \<le> x \<and> x \<le> real u \<and>
   2.202 +             DERIV (\<lambda> x. interpret_floatarith f (xs[n := x])) (xs!n) :> x"
   2.203 +         (is "\<exists> x. _ \<and> _ \<and> DERIV (?i f) _ :> _")
   2.204 +proof (rule exI[of _ "?i ?D (xs!n)"], rule conjI[OF _ conjI])
   2.205 +  let "?i f x" = "interpret_floatarith f (xs[n := x])"
   2.206 +  from approx[OF bnd app]
   2.207 +  show "real l \<le> ?i ?D (xs!n)" and "?i ?D (xs!n) \<le> real u"
   2.208 +    using `n < length xs` by auto
   2.209 +  from DERIV_floatarith[OF `n < length xs`, of f "xs!n"] isDERIV_approx[OF bnd isD]
   2.210 +  show "DERIV (?i f) (xs!n) :> (?i ?D (xs!n))" by simp
   2.211 +qed
   2.212 +
   2.213 +fun lift_bin :: "(float * float) option \<Rightarrow> (float * float) option \<Rightarrow> (float \<Rightarrow> float \<Rightarrow> float \<Rightarrow> float \<Rightarrow> (float * float) option) \<Rightarrow> (float * float) option" where
   2.214 +"lift_bin (Some (l1, u1)) (Some (l2, u2)) f = f l1 u1 l2 u2" |
   2.215 +"lift_bin a b f = None"
   2.216 +
   2.217 +lemma lift_bin:
   2.218 +  assumes lift_bin_Some: "Some (l, u) = lift_bin a b f"
   2.219 +  obtains l1 u1 l2 u2
   2.220 +  where "a = Some (l1, u1)"
   2.221 +  and "b = Some (l2, u2)"
   2.222 +  and "f l1 u1 l2 u2 = Some (l, u)"
   2.223 +using assms by (cases a, simp, cases b, simp, auto)
   2.224 +
   2.225 +fun approx_tse where
   2.226 +"approx_tse prec n 0 c k f bs = approx prec f bs" |
   2.227 +"approx_tse prec n (Suc s) c k f bs =
   2.228 +  (if isDERIV_approx prec n f bs then
   2.229 +    lift_bin (approx prec f (bs[n := Some (c,c)]))
   2.230 +             (approx_tse prec n s c (Suc k) (DERIV_floatarith n f) bs)
   2.231 +             (\<lambda> l1 u1 l2 u2. approx prec
   2.232 +                 (Add (Atom 0)
   2.233 +                      (Mult (Inverse (Num (Float (int k) 0)))
   2.234 +                                 (Mult (Add (Atom (Suc (Suc 0))) (Minus (Num c)))
   2.235 +                                       (Atom (Suc 0))))) [Some (l1, u1), Some (l2, u2), bs!n])
   2.236 +  else approx prec f bs)"
   2.237 +
   2.238 +lemma bounded_by_Cons:
   2.239 +  assumes bnd: "bounded_by xs vs"
   2.240 +  and x: "x \<in> { real l .. real u }"
   2.241 +  shows "bounded_by (x#xs) ((Some (l, u))#vs)"
   2.242 +proof -
   2.243 +  { fix i assume *: "i < length ((Some (l, u))#vs)"
   2.244 +    have "case ((Some (l,u))#vs) ! i of Some (l, u) \<Rightarrow> (x#xs)!i \<in> { real l .. real u } | None \<Rightarrow> True"
   2.245 +    proof (cases i)
   2.246 +      case 0 with x show ?thesis by auto
   2.247 +    next
   2.248 +      case (Suc i) with * have "i < length vs" by auto
   2.249 +      from bnd[THEN bounded_byE, OF this]
   2.250 +      show ?thesis unfolding Suc nth_Cons_Suc .
   2.251 +    qed }
   2.252 +  thus ?thesis by (auto simp add: bounded_by_def)
   2.253 +qed
   2.254 +
   2.255 +lemma approx_tse_generic:
   2.256 +  assumes "bounded_by xs vs"
   2.257 +  and bnd_c: "bounded_by (xs[x := real c]) vs" and "x < length vs" and "x < length xs"
   2.258 +  and bnd_x: "vs ! x = Some (lx, ux)"
   2.259 +  and ate: "Some (l, u) = approx_tse prec x s c k f vs"
   2.260 +  shows "\<exists> n. (\<forall> m < n. \<forall> z \<in> {real lx .. real ux}.
   2.261 +      DERIV (\<lambda> y. interpret_floatarith ((DERIV_floatarith x ^^ m) f) (xs[x := y])) z :>
   2.262 +            (interpret_floatarith ((DERIV_floatarith x ^^ (Suc m)) f) (xs[x := z])))
   2.263 +   \<and> (\<forall> t \<in> {real lx .. real ux}.  (\<Sum> i = 0..<n. inverse (real (\<Prod> j \<in> {k..<k+i}. j)) *
   2.264 +                  interpret_floatarith ((DERIV_floatarith x ^^ i) f) (xs[x := real c]) *
   2.265 +                  (xs!x - real c)^i) +
   2.266 +      inverse (real (\<Prod> j \<in> {k..<k+n}. j)) *
   2.267 +      interpret_floatarith ((DERIV_floatarith x ^^ n) f) (xs[x := t]) *
   2.268 +      (xs!x - real c)^n \<in> {real l .. real u})" (is "\<exists> n. ?taylor f k l u n")
   2.269 +using ate proof (induct s arbitrary: k f l u)
   2.270 +  case 0
   2.271 +  { fix t assume "t \<in> {real lx .. real ux}"
   2.272 +    note bounded_by_update_var[OF `bounded_by xs vs` bnd_x this]
   2.273 +    from approx[OF this 0[unfolded approx_tse.simps]]
   2.274 +    have "(interpret_floatarith f (xs[x := t])) \<in> {real l .. real u}"
   2.275 +      by (auto simp add: algebra_simps)
   2.276 +  } thus ?case by (auto intro!: exI[of _ 0])
   2.277 +next
   2.278 +  case (Suc s)
   2.279 +  show ?case
   2.280 +  proof (cases "isDERIV_approx prec x f vs")
   2.281 +    case False
   2.282 +    note ap = Suc.prems[unfolded approx_tse.simps if_not_P[OF False]]
   2.283 +
   2.284 +    { fix t assume "t \<in> {real lx .. real ux}"
   2.285 +      note bounded_by_update_var[OF `bounded_by xs vs` bnd_x this]
   2.286 +      from approx[OF this ap]
   2.287 +      have "(interpret_floatarith f (xs[x := t])) \<in> {real l .. real u}"
   2.288 +	by (auto simp add: algebra_simps)
   2.289 +    } thus ?thesis by (auto intro!: exI[of _ 0])
   2.290 +  next
   2.291 +    case True
   2.292 +    with Suc.prems
   2.293 +    obtain l1 u1 l2 u2
   2.294 +      where a: "Some (l1, u1) = approx prec f (vs[x := Some (c,c)])"
   2.295 +      and ate: "Some (l2, u2) = approx_tse prec x s c (Suc k) (DERIV_floatarith x f) vs"
   2.296 +      and final: "Some (l, u) = approx prec
   2.297 +        (Add (Atom 0)
   2.298 +             (Mult (Inverse (Num (Float (int k) 0)))
   2.299 +                   (Mult (Add (Atom (Suc (Suc 0))) (Minus (Num c)))
   2.300 +                         (Atom (Suc 0))))) [Some (l1, u1), Some (l2, u2), vs!x]"
   2.301 +      by (auto elim!: lift_bin) blast
   2.302 +
   2.303 +    from bnd_c `x < length xs`
   2.304 +    have bnd: "bounded_by (xs[x:=real c]) (vs[x:= Some (c,c)])"
   2.305 +      by (auto intro!: bounded_by_update)
   2.306 +
   2.307 +    from approx[OF this a]
   2.308 +    have f_c: "interpret_floatarith ((DERIV_floatarith x ^^ 0) f) (xs[x := real c]) \<in> { real l1 .. real u1 }"
   2.309 +              (is "?f 0 (real c) \<in> _")
   2.310 +      by auto
   2.311 +
   2.312 +    { fix f :: "'a \<Rightarrow> 'a" fix n :: nat fix x :: 'a
   2.313 +      have "(f ^^ Suc n) x = (f ^^ n) (f x)"
   2.314 +	by (induct n, auto) }
   2.315 +    note funpow_Suc = this[symmetric]
   2.316 +    from Suc.hyps[OF ate, unfolded this]
   2.317 +    obtain n
   2.318 +      where DERIV_hyp: "\<And> m z. \<lbrakk> m < n ; z \<in> { real lx .. real ux } \<rbrakk> \<Longrightarrow> DERIV (?f (Suc m)) z :> ?f (Suc (Suc m)) z"
   2.319 +      and hyp: "\<forall> t \<in> {real lx .. real ux}. (\<Sum> i = 0..<n. inverse (real (\<Prod> j \<in> {Suc k..<Suc k + i}. j)) * ?f (Suc i) (real c) * (xs!x - real c)^i) +
   2.320 +           inverse (real (\<Prod> j \<in> {Suc k..<Suc k + n}. j)) * ?f (Suc n) t * (xs!x - real c)^n \<in> {real l2 .. real u2}"
   2.321 +          (is "\<forall> t \<in> _. ?X (Suc k) f n t \<in> _")
   2.322 +      by blast
   2.323 +
   2.324 +    { fix m z
   2.325 +      assume "m < Suc n" and bnd_z: "z \<in> { real lx .. real ux }"
   2.326 +      have "DERIV (?f m) z :> ?f (Suc m) z"
   2.327 +      proof (cases m)
   2.328 +	case 0
   2.329 +	with DERIV_floatarith[OF `x < length xs` isDERIV_approx'[OF `bounded_by xs vs` bnd_x bnd_z True]]
   2.330 +	show ?thesis by simp
   2.331 +      next
   2.332 +	case (Suc m')
   2.333 +	hence "m' < n" using `m < Suc n` by auto
   2.334 +	from DERIV_hyp[OF this bnd_z]
   2.335 +	show ?thesis using Suc by simp
   2.336 +      qed } note DERIV = this
   2.337 +
   2.338 +    have "\<And> k i. k < i \<Longrightarrow> {k ..< i} = insert k {Suc k ..< i}" by auto
   2.339 +    hence setprod_head_Suc: "\<And> k i. \<Prod> {k ..< k + Suc i} = k * \<Prod> {Suc k ..< Suc k + i}" by auto
   2.340 +    have setsum_move0: "\<And> k F. setsum F {0..<Suc k} = F 0 + setsum (\<lambda> k. F (Suc k)) {0..<k}"
   2.341 +      unfolding setsum_shift_bounds_Suc_ivl[symmetric]
   2.342 +      unfolding setsum_head_upt_Suc[OF zero_less_Suc] ..
   2.343 +    def C \<equiv> "xs!x - real c"
   2.344 +
   2.345 +    { fix t assume t: "t \<in> {real lx .. real ux}"
   2.346 +      hence "bounded_by [xs!x] [vs!x]"
   2.347 +	using `bounded_by xs vs`[THEN bounded_byE, OF `x < length vs`]
   2.348 +	by (cases "vs!x", auto simp add: bounded_by_def)
   2.349 +
   2.350 +      with hyp[THEN bspec, OF t] f_c
   2.351 +      have "bounded_by [?f 0 (real c), ?X (Suc k) f n t, xs!x] [Some (l1, u1), Some (l2, u2), vs!x]"
   2.352 +	by (auto intro!: bounded_by_Cons)
   2.353 +      from approx[OF this final, unfolded atLeastAtMost_iff[symmetric]]
   2.354 +      have "?X (Suc k) f n t * (xs!x - real c) * inverse (real k) + ?f 0 (real c) \<in> {real l .. real u}"
   2.355 +	by (auto simp add: algebra_simps)
   2.356 +      also have "?X (Suc k) f n t * (xs!x - real c) * inverse (real k) + ?f 0 (real c) =
   2.357 +               (\<Sum> i = 0..<Suc n. inverse (real (\<Prod> j \<in> {k..<k+i}. j)) * ?f i (real c) * (xs!x - real c)^i) +
   2.358 +               inverse (real (\<Prod> j \<in> {k..<k+Suc n}. j)) * ?f (Suc n) t * (xs!x - real c)^Suc n" (is "_ = ?T")
   2.359 +	unfolding funpow_Suc C_def[symmetric] setsum_move0 setprod_head_Suc
   2.360 +	by (auto simp add: algebra_simps setsum_right_distrib[symmetric])
   2.361 +      finally have "?T \<in> {real l .. real u}" . }
   2.362 +    thus ?thesis using DERIV by blast
   2.363 +  qed
   2.364 +qed
   2.365 +
   2.366 +lemma setprod_fact: "\<Prod> {1..<1 + k} = fact (k :: nat)"
   2.367 +proof (induct k)
   2.368 +  case (Suc k)
   2.369 +  have "{ 1 ..< Suc (Suc k) } = insert (Suc k) { 1 ..< Suc k }" by auto
   2.370 +  hence "\<Prod> { 1 ..< Suc (Suc k) } = (Suc k) * \<Prod> { 1 ..< Suc k }" by auto
   2.371 +  thus ?case using Suc by auto
   2.372 +qed simp
   2.373 +
   2.374 +lemma approx_tse:
   2.375 +  assumes "bounded_by xs vs"
   2.376 +  and bnd_x: "vs ! x = Some (lx, ux)" and bnd_c: "real c \<in> {real lx .. real ux}"
   2.377 +  and "x < length vs" and "x < length xs"
   2.378 +  and ate: "Some (l, u) = approx_tse prec x s c 1 f vs"
   2.379 +  shows "interpret_floatarith f xs \<in> { real l .. real u }"
   2.380 +proof -
   2.381 +  def F \<equiv> "\<lambda> n z. interpret_floatarith ((DERIV_floatarith x ^^ n) f) (xs[x := z])"
   2.382 +  hence F0: "F 0 = (\<lambda> z. interpret_floatarith f (xs[x := z]))" by auto
   2.383 +
   2.384 +  hence "bounded_by (xs[x := real c]) vs" and "x < length vs" "x < length xs"
   2.385 +    using `bounded_by xs vs` bnd_x bnd_c `x < length vs` `x < length xs`
   2.386 +    by (auto intro!: bounded_by_update_var)
   2.387 +
   2.388 +  from approx_tse_generic[OF `bounded_by xs vs` this bnd_x ate]
   2.389 +  obtain n
   2.390 +    where DERIV: "\<forall> m z. m < n \<and> real lx \<le> z \<and> z \<le> real ux \<longrightarrow> DERIV (F m) z :> F (Suc m) z"
   2.391 +    and hyp: "\<And> t. t \<in> {real lx .. real ux} \<Longrightarrow>
   2.392 +           (\<Sum> j = 0..<n. inverse (real (fact j)) * F j (real c) * (xs!x - real c)^j) +
   2.393 +             inverse (real (fact n)) * F n t * (xs!x - real c)^n
   2.394 +             \<in> {real l .. real u}" (is "\<And> t. _ \<Longrightarrow> ?taylor t \<in> _")
   2.395 +    unfolding F_def atLeastAtMost_iff[symmetric] setprod_fact by blast
   2.396 +
   2.397 +  have bnd_xs: "xs ! x \<in> { real lx .. real ux }"
   2.398 +    using `bounded_by xs vs`[THEN bounded_byE, OF `x < length vs`] bnd_x by auto
   2.399 +
   2.400 +  show ?thesis
   2.401 +  proof (cases n)
   2.402 +    case 0 thus ?thesis using hyp[OF bnd_xs] unfolding F_def by auto
   2.403 +  next
   2.404 +    case (Suc n')
   2.405 +    show ?thesis
   2.406 +    proof (cases "xs ! x = real c")
   2.407 +      case True
   2.408 +      from True[symmetric] hyp[OF bnd_xs] Suc show ?thesis
   2.409 +	unfolding F_def Suc setsum_head_upt_Suc[OF zero_less_Suc] setsum_shift_bounds_Suc_ivl by auto
   2.410 +    next
   2.411 +      case False
   2.412 +
   2.413 +      have "real lx \<le> real c" "real c \<le> real ux" "real lx \<le> xs!x" "xs!x \<le> real ux"
   2.414 +	using Suc bnd_c `bounded_by xs vs`[THEN bounded_byE, OF `x < length vs`] bnd_x by auto
   2.415 +      from Taylor.taylor[OF zero_less_Suc, of F, OF F0 DERIV[unfolded Suc] this False]
   2.416 +      obtain t where t_bnd: "if xs ! x < real c then xs ! x < t \<and> t < real c else real c < t \<and> t < xs ! x"
   2.417 +	and fl_eq: "interpret_floatarith f (xs[x := xs ! x]) =
   2.418 +	   (\<Sum>m = 0..<Suc n'. F m (real c) / real (fact m) * (xs ! x - real c) ^ m) +
   2.419 +           F (Suc n') t / real (fact (Suc n')) * (xs ! x - real c) ^ Suc n'"
   2.420 +	by blast
   2.421 +
   2.422 +      from t_bnd bnd_xs bnd_c have *: "t \<in> {real lx .. real ux}"
   2.423 +	by (cases "xs ! x < real c", auto)
   2.424 +
   2.425 +      have "interpret_floatarith f (xs[x := xs ! x]) = ?taylor t"
   2.426 +	unfolding fl_eq Suc by (auto simp add: algebra_simps divide_inverse)
   2.427 +      also have "\<dots> \<in> {real l .. real u}" using * by (rule hyp)
   2.428 +      finally show ?thesis by simp
   2.429 +    qed
   2.430 +  qed
   2.431 +qed
   2.432 +
   2.433 +fun approx_tse_form' where
   2.434 +"approx_tse_form' prec t f 0 l u cmp =
   2.435 +  (case approx_tse prec 0 t ((l + u) * Float 1 -1) 1 f [Some (l, u)]
   2.436 +     of Some (l, u) \<Rightarrow> cmp l u | None \<Rightarrow> False)" |
   2.437 +"approx_tse_form' prec t f (Suc s) l u cmp =
   2.438 +  (let m = (l + u) * Float 1 -1
   2.439 +   in approx_tse_form' prec t f s l m cmp \<and>
   2.440 +      approx_tse_form' prec t f s m u cmp)"
   2.441 +
   2.442 +lemma approx_tse_form':
   2.443 +  assumes "approx_tse_form' prec t f s l u cmp" and "x \<in> {real l .. real u}"
   2.444 +  shows "\<exists> l' u' ly uy. x \<in> { real l' .. real u' } \<and> real l \<le> real l' \<and> real u' \<le> real u \<and> cmp ly uy \<and>
   2.445 +                  approx_tse prec 0 t ((l' + u') * Float 1 -1) 1 f [Some (l', u')] = Some (ly, uy)"
   2.446 +using assms proof (induct s arbitrary: l u)
   2.447 +  case 0
   2.448 +  then obtain ly uy
   2.449 +    where *: "approx_tse prec 0 t ((l + u) * Float 1 -1) 1 f [Some (l, u)] = Some (ly, uy)"
   2.450 +    and **: "cmp ly uy" by (auto elim!: option_caseE)
   2.451 +  with 0 show ?case by (auto intro!: exI)
   2.452 +next
   2.453 +  case (Suc s)
   2.454 +  let ?m = "(l + u) * Float 1 -1"
   2.455 +  from Suc.prems
   2.456 +  have l: "approx_tse_form' prec t f s l ?m cmp"
   2.457 +    and u: "approx_tse_form' prec t f s ?m u cmp"
   2.458 +    by (auto simp add: Let_def)
   2.459 +
   2.460 +  have m_l: "real l \<le> real ?m" and m_u: "real ?m \<le> real u"
   2.461 +    unfolding le_float_def using Suc.prems by auto
   2.462 +
   2.463 +  with `x \<in> { real l .. real u }`
   2.464 +  have "x \<in> { real l .. real ?m} \<or> x \<in> { real ?m .. real u }" by auto
   2.465 +  thus ?case
   2.466 +  proof (rule disjE)
   2.467 +    assume "x \<in> { real l .. real ?m}"
   2.468 +    from Suc.hyps[OF l this]
   2.469 +    obtain l' u' ly uy
   2.470 +      where "x \<in> { real l' .. real u' } \<and> real l \<le> real l' \<and> real u' \<le> real ?m \<and> cmp ly uy \<and>
   2.471 +                  approx_tse prec 0 t ((l' + u') * Float 1 -1) 1 f [Some (l', u')] = Some (ly, uy)" by blast
   2.472 +    with m_u show ?thesis by (auto intro!: exI)
   2.473 +  next
   2.474 +    assume "x \<in> { real ?m .. real u }"
   2.475 +    from Suc.hyps[OF u this]
   2.476 +    obtain l' u' ly uy
   2.477 +      where "x \<in> { real l' .. real u' } \<and> real ?m \<le> real l' \<and> real u' \<le> real u \<and> cmp ly uy \<and>
   2.478 +                  approx_tse prec 0 t ((l' + u') * Float 1 -1) 1 f [Some (l', u')] = Some (ly, uy)" by blast
   2.479 +    with m_u show ?thesis by (auto intro!: exI)
   2.480 +  qed
   2.481 +qed
   2.482 +
   2.483 +lemma approx_tse_form'_less:
   2.484 +  assumes tse: "approx_tse_form' prec t (Add a (Minus b)) s l u (\<lambda> l u. 0 < l)"
   2.485 +  and x: "x \<in> {real l .. real u}"
   2.486 +  shows "interpret_floatarith b [x] < interpret_floatarith a [x]"
   2.487 +proof -
   2.488 +  from approx_tse_form'[OF tse x]
   2.489 +  obtain l' u' ly uy
   2.490 +    where x': "x \<in> { real l' .. real u' }" and "real l \<le> real l'"
   2.491 +    and "real u' \<le> real u" and "0 < ly"
   2.492 +    and tse: "approx_tse prec 0 t ((l' + u') * Float 1 -1) 1 (Add a (Minus b)) [Some (l', u')] = Some (ly, uy)"
   2.493 +    by blast
   2.494 +
   2.495 +  hence "bounded_by [x] [Some (l', u')]" by (auto simp add: bounded_by_def)
   2.496 +
   2.497 +  from approx_tse[OF this _ _ _ _ tse[symmetric], of l' u'] x'
   2.498 +  have "real ly \<le> interpret_floatarith a [x] - interpret_floatarith b [x]"
   2.499 +    by (auto simp add: diff_minus)
   2.500 +  from order_less_le_trans[OF `0 < ly`[unfolded less_float_def] this]
   2.501 +  show ?thesis by auto
   2.502 +qed
   2.503 +
   2.504 +lemma approx_tse_form'_le:
   2.505 +  assumes tse: "approx_tse_form' prec t (Add a (Minus b)) s l u (\<lambda> l u. 0 \<le> l)"
   2.506 +  and x: "x \<in> {real l .. real u}"
   2.507 +  shows "interpret_floatarith b [x] \<le> interpret_floatarith a [x]"
   2.508 +proof -
   2.509 +  from approx_tse_form'[OF tse x]
   2.510 +  obtain l' u' ly uy
   2.511 +    where x': "x \<in> { real l' .. real u' }" and "real l \<le> real l'"
   2.512 +    and "real u' \<le> real u" and "0 \<le> ly"
   2.513 +    and tse: "approx_tse prec 0 t ((l' + u') * Float 1 -1) 1 (Add a (Minus b)) [Some (l', u')] = Some (ly, uy)"
   2.514 +    by blast
   2.515 +
   2.516 +  hence "bounded_by [x] [Some (l', u')]" by (auto simp add: bounded_by_def)
   2.517 +
   2.518 +  from approx_tse[OF this _ _ _ _ tse[symmetric], of l' u'] x'
   2.519 +  have "real ly \<le> interpret_floatarith a [x] - interpret_floatarith b [x]"
   2.520 +    by (auto simp add: diff_minus)
   2.521 +  from order_trans[OF `0 \<le> ly`[unfolded le_float_def] this]
   2.522 +  show ?thesis by auto
   2.523 +qed
   2.524 +
   2.525 +definition
   2.526 +"approx_tse_form prec t s f =
   2.527 +  (case f
   2.528 +   of (Bound x a b f) \<Rightarrow> x = Atom 0 \<and>
   2.529 +     (case (approx prec a [None], approx prec b [None])
   2.530 +      of (Some (l, u), Some (l', u')) \<Rightarrow>
   2.531 +        (case f
   2.532 +         of Less lf rt \<Rightarrow> approx_tse_form' prec t (Add rt (Minus lf)) s l u' (\<lambda> l u. 0 < l)
   2.533 +          | LessEqual lf rt \<Rightarrow> approx_tse_form' prec t (Add rt (Minus lf)) s l u' (\<lambda> l u. 0 \<le> l)
   2.534 +          | AtLeastAtMost x lf rt \<Rightarrow>
   2.535 +            approx_tse_form' prec t (Add x (Minus lf)) s l u' (\<lambda> l u. 0 \<le> l) \<and>
   2.536 +            approx_tse_form' prec t (Add rt (Minus x)) s l u' (\<lambda> l u. 0 \<le> l)
   2.537 +          | _ \<Rightarrow> False)
   2.538 +       | _ \<Rightarrow> False)
   2.539 +   | _ \<Rightarrow> False)"
   2.540 +
   2.541 +lemma approx_tse_form:
   2.542 +  assumes "approx_tse_form prec t s f"
   2.543 +  shows "interpret_form f [x]"
   2.544 +proof (cases f)
   2.545 +  case (Bound i a b f') note f_def = this
   2.546 +  with assms obtain l u l' u'
   2.547 +    where a: "approx prec a [None] = Some (l, u)"
   2.548 +    and b: "approx prec b [None] = Some (l', u')"
   2.549 +    unfolding approx_tse_form_def by (auto elim!: option_caseE)
   2.550 +
   2.551 +  from Bound assms have "i = Atom 0" unfolding approx_tse_form_def by auto
   2.552 +  hence i: "interpret_floatarith i [x] = x" by auto
   2.553 +
   2.554 +  { let "?f z" = "interpret_floatarith z [x]"
   2.555 +    assume "?f i \<in> { ?f a .. ?f b }"
   2.556 +    with approx[OF _ a[symmetric], of "[x]"] approx[OF _ b[symmetric], of "[x]"]
   2.557 +    have bnd: "x \<in> { real l .. real u'}" unfolding bounded_by_def i by auto
   2.558 +
   2.559 +    have "interpret_form f' [x]"
   2.560 +    proof (cases f')
   2.561 +      case (Less lf rt)
   2.562 +      with Bound a b assms
   2.563 +      have "approx_tse_form' prec t (Add rt (Minus lf)) s l u' (\<lambda> l u. 0 < l)"
   2.564 +	unfolding approx_tse_form_def by auto
   2.565 +      from approx_tse_form'_less[OF this bnd]
   2.566 +      show ?thesis using Less by auto
   2.567 +    next
   2.568 +      case (LessEqual lf rt)
   2.569 +      with Bound a b assms
   2.570 +      have "approx_tse_form' prec t (Add rt (Minus lf)) s l u' (\<lambda> l u. 0 \<le> l)"
   2.571 +	unfolding approx_tse_form_def by auto
   2.572 +      from approx_tse_form'_le[OF this bnd]
   2.573 +      show ?thesis using LessEqual by auto
   2.574 +    next
   2.575 +      case (AtLeastAtMost x lf rt)
   2.576 +      with Bound a b assms
   2.577 +      have "approx_tse_form' prec t (Add rt (Minus x)) s l u' (\<lambda> l u. 0 \<le> l)"
   2.578 +	and "approx_tse_form' prec t (Add x (Minus lf)) s l u' (\<lambda> l u. 0 \<le> l)"
   2.579 +	unfolding approx_tse_form_def by auto
   2.580 +      from approx_tse_form'_le[OF this(1) bnd] approx_tse_form'_le[OF this(2) bnd]
   2.581 +      show ?thesis using AtLeastAtMost by auto
   2.582 +    next
   2.583 +      case (Bound x a b f') with assms
   2.584 +      show ?thesis by (auto elim!: option_caseE simp add: f_def approx_tse_form_def)
   2.585 +    next
   2.586 +      case (Assign x a f') with assms
   2.587 +      show ?thesis by (auto elim!: option_caseE simp add: f_def approx_tse_form_def)
   2.588 +    qed } thus ?thesis unfolding f_def by auto
   2.589 +next case Assign with assms show ?thesis by (auto simp add: approx_tse_form_def)
   2.590 +next case LessEqual with assms show ?thesis by (auto simp add: approx_tse_form_def)
   2.591 +next case Less with assms show ?thesis by (auto simp add: approx_tse_form_def)
   2.592 +next case AtLeastAtMost with assms show ?thesis by (auto simp add: approx_tse_form_def)
   2.593 +qed
   2.594 +
   2.595  subsection {* Implement proof method \texttt{approximation} *}
   2.596  
   2.597  lemmas interpret_form_equations = interpret_form.simps interpret_floatarith.simps interpret_floatarith_num
   2.598 @@ -2648,6 +3216,7 @@
   2.599  @{code_datatype form = Bound | Assign | Less | LessEqual | AtLeastAtMost}
   2.600  
   2.601  val approx_form = @{code approx_form}
   2.602 +val approx_tse_form = @{code approx_tse_form}
   2.603  val approx' = @{code approx'}
   2.604  
   2.605  end
   2.606 @@ -2675,6 +3244,7 @@
   2.607              "Float'_Arith.AtLeastAtMost/ (_,/ _,/ _)")
   2.608  
   2.609  code_const approx_form (Eval "Float'_Arith.approx'_form")
   2.610 +code_const approx_tse_form (Eval "Float'_Arith.approx'_tse'_form")
   2.611  code_const approx' (Eval "Float'_Arith.approx'")
   2.612  
   2.613  ML {*
   2.614 @@ -2712,30 +3282,49 @@
   2.615  
   2.616    val form_equations = PureThy.get_thms @{theory} "interpret_form_equations";
   2.617  
   2.618 -  fun rewrite_interpret_form_tac ctxt prec splitting i st = let
   2.619 +  fun rewrite_interpret_form_tac ctxt prec splitting taylor i st = let
   2.620 +      fun lookup_splitting (Free (name, typ))
   2.621 +        = case AList.lookup (op =) splitting name
   2.622 +          of SOME s => HOLogic.mk_number @{typ nat} s
   2.623 +           | NONE => @{term "0 :: nat"}
   2.624        val vs = nth (prems_of st) (i - 1)
   2.625                 |> Logic.strip_imp_concl
   2.626                 |> HOLogic.dest_Trueprop
   2.627                 |> Term.strip_comb |> snd |> List.last
   2.628                 |> HOLogic.dest_list
   2.629 -      val n = vs |> length
   2.630 -              |> HOLogic.mk_number @{typ nat}
   2.631 -              |> Thm.cterm_of (ProofContext.theory_of ctxt)
   2.632        val p = prec
   2.633                |> HOLogic.mk_number @{typ nat}
   2.634                |> Thm.cterm_of (ProofContext.theory_of ctxt)
   2.635 -      val s = vs
   2.636 -              |> map (fn Free (name, typ) =>
   2.637 -                      case AList.lookup (op =) splitting name of
   2.638 -                        SOME s => HOLogic.mk_number @{typ nat} s
   2.639 -                      | NONE => @{term "0 :: nat"})
   2.640 -              |> HOLogic.mk_list @{typ nat}
   2.641 +    in case taylor
   2.642 +    of NONE => let
   2.643 +         val n = vs |> length
   2.644 +                 |> HOLogic.mk_number @{typ nat}
   2.645 +                 |> Thm.cterm_of (ProofContext.theory_of ctxt)
   2.646 +         val s = vs
   2.647 +                 |> map lookup_splitting
   2.648 +                 |> HOLogic.mk_list @{typ nat}
   2.649 +                 |> Thm.cterm_of (ProofContext.theory_of ctxt)
   2.650 +       in
   2.651 +         (rtac (Thm.instantiate ([], [(@{cpat "?n::nat"}, n),
   2.652 +                                     (@{cpat "?prec::nat"}, p),
   2.653 +                                     (@{cpat "?ss::nat list"}, s)])
   2.654 +              @{thm "approx_form"}) i
   2.655 +          THEN simp_tac @{simpset} i) st
   2.656 +       end
   2.657 +
   2.658 +     | SOME t => if length vs <> 1 then raise (TERM ("More than one variable used for taylor series expansion", [prop_of st]))
   2.659 +       else let
   2.660 +         val t = t
   2.661 +              |> HOLogic.mk_number @{typ nat}
   2.662                |> Thm.cterm_of (ProofContext.theory_of ctxt)
   2.663 -    in
   2.664 -      rtac (Thm.instantiate ([], [(@{cpat "?n::nat"}, n),
   2.665 -                                  (@{cpat "?prec::nat"}, p),
   2.666 -                                  (@{cpat "?ss::nat list"}, s)])
   2.667 -           @{thm "approx_form"}) i st
   2.668 +         val s = vs |> map lookup_splitting |> hd
   2.669 +              |> Thm.cterm_of (ProofContext.theory_of ctxt)
   2.670 +       in
   2.671 +         rtac (Thm.instantiate ([], [(@{cpat "?s::nat"}, s),
   2.672 +                                     (@{cpat "?t::nat"}, t),
   2.673 +                                     (@{cpat "?prec::nat"}, p)])
   2.674 +              @{thm "approx_tse_form"}) i st
   2.675 +       end
   2.676      end
   2.677  
   2.678    (* copied from Tools/induct.ML should probably in args.ML *)
   2.679 @@ -2751,11 +3340,15 @@
   2.680    by auto
   2.681  
   2.682  method_setup approximation = {*
   2.683 -  Scan.lift (OuterParse.nat) --
   2.684 +  Scan.lift (OuterParse.nat)
   2.685 +  --
   2.686    Scan.optional (Scan.lift (Args.$$$ "splitting" |-- Args.colon)
   2.687      |-- OuterParse.and_list' (free --| Scan.lift (Args.$$$ "=") -- Scan.lift OuterParse.nat)) []
   2.688 +  --
   2.689 +  Scan.option (Scan.lift (Args.$$$ "taylor" |-- Args.colon)
   2.690 +    |-- (free |-- Scan.lift (Args.$$$ "=") |-- Scan.lift OuterParse.nat))
   2.691    >>
   2.692 -  (fn (prec, splitting) => fn ctxt =>
   2.693 +  (fn ((prec, splitting), taylor) => fn ctxt =>
   2.694      SIMPLE_METHOD' (fn i =>
   2.695        REPEAT (FIRST' [etac @{thm intervalE},
   2.696                        etac @{thm meta_eqE},
   2.697 @@ -2763,15 +3356,10 @@
   2.698        THEN METAHYPS (reorder_bounds_tac i) i
   2.699        THEN TRY (filter_prems_tac (K false) i)
   2.700        THEN DETERM (Reflection.genreify_tac ctxt form_equations NONE i)
   2.701 -      THEN print_tac "approximation"
   2.702 -      THEN rewrite_interpret_form_tac ctxt prec splitting i
   2.703 -      THEN simp_tac @{simpset} i
   2.704 +      THEN rewrite_interpret_form_tac ctxt prec splitting taylor i
   2.705        THEN gen_eval_tac eval_oracle ctxt i))
   2.706   *} "real number approximation"
   2.707  
   2.708 -lemma "\<phi> \<in> {0..1 :: real} \<Longrightarrow> \<phi> < \<phi> + 0.7"
   2.709 -  by (approximation 10 splitting: "\<phi>" = 1)
   2.710 -
   2.711  ML {*
   2.712    fun dest_interpret (@{const "interpret_floatarith"} $ b $ xs) = (b, xs)
   2.713    | dest_interpret t = raise TERM ("dest_interpret", [t])
     3.1 --- a/src/HOL/Decision_Procs/ex/Approximation_Ex.thy	Tue Jun 30 00:57:24 2009 +0200
     3.2 +++ b/src/HOL/Decision_Procs/ex/Approximation_Ex.thy	Mon Jun 29 23:29:11 2009 +0200
     3.3 @@ -8,13 +8,28 @@
     3.4  
     3.5  Here are some examples how to use the approximation method.
     3.6  
     3.7 -The parameter passed to the method specifies the precision used by the computations, it is specified
     3.8 -as number of bits to calculate. When a variable is used it needs to be bounded by an interval. This
     3.9 -interval is specified as a conjunction of the lower and upper bound. It must have the form
    3.10 -@{text "\<lbrakk> l\<^isub>1 \<le> x\<^isub>1 \<and> x\<^isub>1 \<le> u\<^isub>1 ; \<dots> ; l\<^isub>n \<le> x\<^isub>n \<and> x\<^isub>n \<le> u\<^isub>n \<rbrakk> \<Longrightarrow> F"} where @{term F} is the formula, and
    3.11 -@{text "x\<^isub>1, \<dots>, x\<^isub>n"} are the variables. The lower bounds @{text "l\<^isub>1, \<dots>, l\<^isub>n"} and the upper bounds
    3.12 -@{text "u\<^isub>1, \<dots>, u\<^isub>n"} must either be integer numerals, floating point numbers or of the form
    3.13 -@{term "m * pow2 e"} to specify a exact floating point value.
    3.14 +The approximation method has the following syntax:
    3.15 +
    3.16 +\verb| approximate "prec" (splitting: "x" = "depth" and "y" = "depth" \<dots>)? (taylor: "z" = "derivates")? |
    3.17 +
    3.18 +Here "prec" specifies the precision used in all computations, it is specified as
    3.19 +number of bits to calculate. In the proposition to prove an arbitrary amount of
    3.20 +variables can be used, but each one need to be bounded by an upper and lower
    3.21 +bound.
    3.22 +
    3.23 +To specify the bounds either @{term "l\<^isub>1 \<le> x \<and> x \<le> u\<^isub>1"},
    3.24 +@{term "x \<in> { l\<^isub>1 .. u\<^isub>1 }"} or @{term "x = bnd"} can be used. Where the
    3.25 +bound specification are again arithmetic formulas containing variables. They can
    3.26 +be connected using either meta level or HOL equivalence.
    3.27 +
    3.28 +To use interval splitting add for each variable whos interval should be splitted
    3.29 +to the "splitting:" parameter. The parameter specifies how often each interval
    3.30 +should be divided, e.g. when x = 16 is specified, there will be $65536 = 2^16$ intervals
    3.31 +to be calculated.
    3.32 +
    3.33 +To use taylor series expansion specify the variable to derive. You need to
    3.34 +specify the amount of derivations to compute. When using taylor series expansion
    3.35 +only one variable can be used.
    3.36  
    3.37  *}
    3.38  
    3.39 @@ -57,4 +72,7 @@
    3.40    shows "g / v * tan (35 * d) \<in> { 3 * d .. 3.1 * d }"
    3.41    using assms by (approximation 80)
    3.42  
    3.43 +lemma "\<phi> \<in> { 0 .. 1 :: real } \<longrightarrow> \<phi> ^ 2 \<le> \<phi>"
    3.44 +  by (approximation 30 splitting: \<phi>=1 taylor: \<phi> = 3)
    3.45 +
    3.46  end
     4.1 --- a/src/HOL/Library/Float.thy	Tue Jun 30 00:57:24 2009 +0200
     4.2 +++ b/src/HOL/Library/Float.thy	Mon Jun 29 23:29:11 2009 +0200
     4.3 @@ -59,6 +59,12 @@
     4.4     "real (Float -1 0) = -1" and "real (Float (number_of n) 0) = number_of n"
     4.5    by auto
     4.6  
     4.7 +lemma float_number_of[simp]: "real (number_of x :: float) = number_of x"
     4.8 +  by (simp only:number_of_float_def Float_num[unfolded number_of_is_id])
     4.9 +
    4.10 +lemma float_number_of_int[simp]: "real (Float n 0) = real n"
    4.11 +  by (simp add: Float_num[unfolded number_of_is_id] real_of_float_simp pow2_def)
    4.12 +
    4.13  lemma pow2_0[simp]: "pow2 0 = 1" by simp
    4.14  lemma pow2_1[simp]: "pow2 1 = 2" by simp
    4.15  lemma pow2_neg: "pow2 x = inverse (pow2 (-x))" by simp