Added new Float theory and moved old Library/Float.thy to ComputeFloat
authorhoelzl
Thu Feb 05 11:45:15 2009 +0100 (2009-02-05)
changeset 29804e15b74577368
parent 29803 c56a5571f60a
child 29805 a5da150bd0ab
Added new Float theory and moved old Library/Float.thy to ComputeFloat
src/HOL/IsaMakefile
src/HOL/Library/Float.thy
src/HOL/Matrix/cplex/Cplex.thy
src/HOL/Matrix/cplex/matrixlp.ML
src/HOL/Tools/ComputeFloat.thy
src/HOL/Tools/ComputeNumeral.thy
src/HOL/Tools/float_arith.ML
     1.1 --- a/src/HOL/IsaMakefile	Thu Feb 05 11:34:42 2009 +0100
     1.2 +++ b/src/HOL/IsaMakefile	Thu Feb 05 11:45:15 2009 +0100
     1.3 @@ -335,7 +335,7 @@
     1.4    Library/Mapping.thy	Library/Numeral_Type.thy	Library/Reflection.thy		\
     1.5    Library/Boolean_Algebra.thy Library/Countable.thy	\
     1.6    Library/RBT.thy	Library/Univ_Poly.thy	\
     1.7 -  Library/Enum.thy Library/Float.thy $(SRC)/Tools/float.ML $(SRC)/HOL/Tools/float_arith.ML \
     1.8 +  Library/Enum.thy Library/Float.thy $(SRC)/Tools/float.ML \
     1.9    Library/reify_data.ML Library/reflection.ML
    1.10  	@cd Library; $(ISABELLE_TOOL) usedir $(OUT)/HOL Library
    1.11  
    1.12 @@ -883,6 +883,7 @@
    1.13    $(SRC)/Tools/Compute_Oracle/am_ghc.ML					\
    1.14    $(SRC)/Tools/Compute_Oracle/am_sml.ML					\
    1.15    $(SRC)/Tools/Compute_Oracle/compute.ML	\
    1.16 +  Tools/ComputeFloat.thy Tools/float_arith.ML \
    1.17    Matrix/Matrix.thy Matrix/SparseMatrix.thy Matrix/LP.thy		\
    1.18    Matrix/document/root.tex Matrix/ROOT.ML Matrix/cplex/Cplex.thy	\
    1.19    Matrix/cplex/CplexMatrixConverter.ML Matrix/cplex/Cplex_tools.ML	\
     2.1 --- a/src/HOL/Library/Float.thy	Thu Feb 05 11:34:42 2009 +0100
     2.2 +++ b/src/HOL/Library/Float.thy	Thu Feb 05 11:45:15 2009 +0100
     2.3 @@ -1,30 +1,56 @@
     2.4 -(*  Title:  HOL/Real/Float.thy
     2.5 -    Author: Steven Obua
     2.6 -*)
     2.7 -
     2.8 -header {* Floating Point Representation of the Reals *}
     2.9 -
    2.10 +(* Title:    HOL/Library/Float.thy
    2.11 + * Author:   Steven Obua 2008
    2.12 + *           Johannes Hölzl, TU Muenchen <hoelzl@in.tum.de> 2008 / 2009
    2.13 + *)
    2.14  theory Float
    2.15  imports Complex_Main
    2.16 -uses "~~/src/Tools/float.ML" ("~~/src/HOL/Tools/float_arith.ML")
    2.17  begin
    2.18  
    2.19  definition
    2.20    pow2 :: "int \<Rightarrow> real" where
    2.21 -  "pow2 a = (if (0 <= a) then (2^(nat a)) else (inverse (2^(nat (-a)))))"
    2.22 +  [simp]: "pow2 a = (if (0 <= a) then (2^(nat a)) else (inverse (2^(nat (-a)))))"
    2.23 +
    2.24 +datatype float = Float int int
    2.25 +
    2.26 +fun Ifloat :: "float \<Rightarrow> real" where
    2.27 +"Ifloat (Float a b) = real a * pow2 b"
    2.28  
    2.29 -definition
    2.30 -  float :: "int * int \<Rightarrow> real" where
    2.31 -  "float x = real (fst x) * pow2 (snd x)"
    2.32 +instantiation float :: zero begin
    2.33 +definition zero_float where "0 = Float 0 0" 
    2.34 +instance ..
    2.35 +end
    2.36 +
    2.37 +instantiation float :: one begin
    2.38 +definition one_float where "1 = Float 1 0"
    2.39 +instance ..
    2.40 +end
    2.41 +
    2.42 +instantiation float :: number begin
    2.43 +definition number_of_float where "number_of n = Float n 0"
    2.44 +instance ..
    2.45 +end
    2.46  
    2.47 -lemma pow2_0[simp]: "pow2 0 = 1"
    2.48 -by (simp add: pow2_def)
    2.49 +fun mantissa :: "float \<Rightarrow> int" where
    2.50 +"mantissa (Float a b) = a"
    2.51 +
    2.52 +fun scale :: "float \<Rightarrow> int" where
    2.53 +"scale (Float a b) = b"
    2.54 +
    2.55 +lemma Ifloat_neg_exp: "e < 0 \<Longrightarrow> Ifloat (Float m e) = real m * inverse (2^nat (-e))" by auto
    2.56 +lemma Ifloat_nge0_exp: "\<not> 0 \<le> e \<Longrightarrow> Ifloat (Float m e) = real m * inverse (2^nat (-e))" by auto
    2.57 +lemma Ifloat_ge0_exp: "0 \<le> e \<Longrightarrow> Ifloat (Float m e) = real m * (2^nat e)" by auto
    2.58  
    2.59 -lemma pow2_1[simp]: "pow2 1 = 2"
    2.60 -by (simp add: pow2_def)
    2.61 +lemma Float_num[simp]: shows
    2.62 +   "Ifloat (Float 1 0) = 1" and "Ifloat (Float 1 1) = 2" and "Ifloat (Float 1 2) = 4" and 
    2.63 +   "Ifloat (Float 1 -1) = 1/2" and "Ifloat (Float 1 -2) = 1/4" and "Ifloat (Float 1 -3) = 1/8" and
    2.64 +   "Ifloat (Float -1 0) = -1" and "Ifloat (Float (number_of n) 0) = number_of n"
    2.65 +  by auto
    2.66  
    2.67 -lemma pow2_neg: "pow2 x = inverse (pow2 (-x))"
    2.68 -by (simp add: pow2_def)
    2.69 +lemma pow2_0[simp]: "pow2 0 = 1" by simp
    2.70 +lemma pow2_1[simp]: "pow2 1 = 2" by simp
    2.71 +lemma pow2_neg: "pow2 x = inverse (pow2 (-x))" by simp
    2.72 +
    2.73 +declare pow2_def[simp del]
    2.74  
    2.75  lemma pow2_add1: "pow2 (1 + a) = 2 * (pow2 a)"
    2.76  proof -
    2.77 @@ -96,281 +122,43 @@
    2.78    qed
    2.79  qed
    2.80  
    2.81 -lemma "float (a, e) + float (b, e) = float (a + b, e)"
    2.82 -by (simp add: float_def algebra_simps)
    2.83 +lemma float_components[simp]: "Float (mantissa f) (scale f) = f" by (cases f, auto)
    2.84 +
    2.85 +lemma float_split: "\<exists> a b. x = Float a b" by (cases x, auto)
    2.86  
    2.87 -definition
    2.88 -  int_of_real :: "real \<Rightarrow> int" where
    2.89 -  "int_of_real x = (SOME y. real y = x)"
    2.90 +lemma float_split2: "(\<forall> a b. x \<noteq> Float a b) = False" by (auto simp add: float_split)
    2.91 +
    2.92 +lemma float_zero[simp]: "Ifloat (Float 0 e) = 0" by simp
    2.93 +
    2.94 +lemma abs_div_2_less: "a \<noteq> 0 \<Longrightarrow> a \<noteq> -1 \<Longrightarrow> abs((a::int) div 2) < abs a"
    2.95 +by arith
    2.96  
    2.97 -definition
    2.98 -  real_is_int :: "real \<Rightarrow> bool" where
    2.99 -  "real_is_int x = (EX (u::int). x = real u)"
   2.100 +function normfloat :: "float \<Rightarrow> float" where
   2.101 +"normfloat (Float a b) = (if a \<noteq> 0 \<and> even a then normfloat (Float (a div 2) (b+1)) else if a=0 then Float 0 0 else Float a b)"
   2.102 +by pat_completeness auto
   2.103 +termination by (relation "measure (nat o abs o mantissa)") (auto intro: abs_div_2_less)
   2.104 +declare normfloat.simps[simp del]
   2.105  
   2.106 -lemma real_is_int_def2: "real_is_int x = (x = real (int_of_real x))"
   2.107 -by (auto simp add: real_is_int_def int_of_real_def)
   2.108 +theorem normfloat[symmetric, simp]: "Ifloat f = Ifloat (normfloat f)"
   2.109 +proof (induct f rule: normfloat.induct)
   2.110 +  case (1 a b)
   2.111 +  have real2: "2 = real (2::int)"
   2.112 +    by auto
   2.113 +  show ?case
   2.114 +    apply (subst normfloat.simps)
   2.115 +    apply (auto simp add: float_zero)
   2.116 +    apply (subst 1[symmetric])
   2.117 +    apply (auto simp add: pow2_add even_def)
   2.118 +    done
   2.119 +qed
   2.120  
   2.121 -lemma float_transfer: "real_is_int ((real a)*(pow2 c)) \<Longrightarrow> float (a, b) = float (int_of_real ((real a)*(pow2 c)), b - c)"
   2.122 -by (simp add: float_def real_is_int_def2 pow2_add[symmetric])
   2.123 +lemma pow2_neq_zero[simp]: "pow2 x \<noteq> 0"
   2.124 +  by (auto simp add: pow2_def)
   2.125  
   2.126  lemma pow2_int: "pow2 (int c) = 2^c"
   2.127  by (simp add: pow2_def)
   2.128  
   2.129 -lemma float_transfer_nat: "float (a, b) = float (a * 2^c, b - int c)"
   2.130 -by (simp add: float_def pow2_int[symmetric] pow2_add[symmetric])
   2.131 -
   2.132 -lemma real_is_int_real[simp]: "real_is_int (real (x::int))"
   2.133 -by (auto simp add: real_is_int_def int_of_real_def)
   2.134 -
   2.135 -lemma int_of_real_real[simp]: "int_of_real (real x) = x"
   2.136 -by (simp add: int_of_real_def)
   2.137 -
   2.138 -lemma real_int_of_real[simp]: "real_is_int x \<Longrightarrow> real (int_of_real x) = x"
   2.139 -by (auto simp add: int_of_real_def real_is_int_def)
   2.140 -
   2.141 -lemma real_is_int_add_int_of_real: "real_is_int a \<Longrightarrow> real_is_int b \<Longrightarrow> (int_of_real (a+b)) = (int_of_real a) + (int_of_real b)"
   2.142 -by (auto simp add: int_of_real_def real_is_int_def)
   2.143 -
   2.144 -lemma real_is_int_add[simp]: "real_is_int a \<Longrightarrow> real_is_int b \<Longrightarrow> real_is_int (a+b)"
   2.145 -apply (subst real_is_int_def2)
   2.146 -apply (simp add: real_is_int_add_int_of_real real_int_of_real)
   2.147 -done
   2.148 -
   2.149 -lemma int_of_real_sub: "real_is_int a \<Longrightarrow> real_is_int b \<Longrightarrow> (int_of_real (a-b)) = (int_of_real a) - (int_of_real b)"
   2.150 -by (auto simp add: int_of_real_def real_is_int_def)
   2.151 -
   2.152 -lemma real_is_int_sub[simp]: "real_is_int a \<Longrightarrow> real_is_int b \<Longrightarrow> real_is_int (a-b)"
   2.153 -apply (subst real_is_int_def2)
   2.154 -apply (simp add: int_of_real_sub real_int_of_real)
   2.155 -done
   2.156 -
   2.157 -lemma real_is_int_rep: "real_is_int x \<Longrightarrow> ?! (a::int). real a = x"
   2.158 -by (auto simp add: real_is_int_def)
   2.159 -
   2.160 -lemma int_of_real_mult:
   2.161 -  assumes "real_is_int a" "real_is_int b"
   2.162 -  shows "(int_of_real (a*b)) = (int_of_real a) * (int_of_real b)"
   2.163 -proof -
   2.164 -  from prems have a: "?! (a'::int). real a' = a" by (rule_tac real_is_int_rep, auto)
   2.165 -  from prems have b: "?! (b'::int). real b' = b" by (rule_tac real_is_int_rep, auto)
   2.166 -  from a obtain a'::int where a':"a = real a'" by auto
   2.167 -  from b obtain b'::int where b':"b = real b'" by auto
   2.168 -  have r: "real a' * real b' = real (a' * b')" by auto
   2.169 -  show ?thesis
   2.170 -    apply (simp add: a' b')
   2.171 -    apply (subst r)
   2.172 -    apply (simp only: int_of_real_real)
   2.173 -    done
   2.174 -qed
   2.175 -
   2.176 -lemma real_is_int_mult[simp]: "real_is_int a \<Longrightarrow> real_is_int b \<Longrightarrow> real_is_int (a*b)"
   2.177 -apply (subst real_is_int_def2)
   2.178 -apply (simp add: int_of_real_mult)
   2.179 -done
   2.180 -
   2.181 -lemma real_is_int_0[simp]: "real_is_int (0::real)"
   2.182 -by (simp add: real_is_int_def int_of_real_def)
   2.183 -
   2.184 -lemma real_is_int_1[simp]: "real_is_int (1::real)"
   2.185 -proof -
   2.186 -  have "real_is_int (1::real) = real_is_int(real (1::int))" by auto
   2.187 -  also have "\<dots> = True" by (simp only: real_is_int_real)
   2.188 -  ultimately show ?thesis by auto
   2.189 -qed
   2.190 -
   2.191 -lemma real_is_int_n1: "real_is_int (-1::real)"
   2.192 -proof -
   2.193 -  have "real_is_int (-1::real) = real_is_int(real (-1::int))" by auto
   2.194 -  also have "\<dots> = True" by (simp only: real_is_int_real)
   2.195 -  ultimately show ?thesis by auto
   2.196 -qed
   2.197 -
   2.198 -lemma real_is_int_number_of[simp]: "real_is_int ((number_of \<Colon> int \<Rightarrow> real) x)"
   2.199 -proof -
   2.200 -  have neg1: "real_is_int (-1::real)"
   2.201 -  proof -
   2.202 -    have "real_is_int (-1::real) = real_is_int(real (-1::int))" by auto
   2.203 -    also have "\<dots> = True" by (simp only: real_is_int_real)
   2.204 -    ultimately show ?thesis by auto
   2.205 -  qed
   2.206 -
   2.207 -  {
   2.208 -    fix x :: int
   2.209 -    have "real_is_int ((number_of \<Colon> int \<Rightarrow> real) x)"
   2.210 -      unfolding number_of_eq
   2.211 -      apply (induct x)
   2.212 -      apply (induct_tac n)
   2.213 -      apply (simp)
   2.214 -      apply (simp)
   2.215 -      apply (induct_tac n)
   2.216 -      apply (simp add: neg1)
   2.217 -    proof -
   2.218 -      fix n :: nat
   2.219 -      assume rn: "(real_is_int (of_int (- (int (Suc n)))))"
   2.220 -      have s: "-(int (Suc (Suc n))) = -1 + - (int (Suc n))" by simp
   2.221 -      show "real_is_int (of_int (- (int (Suc (Suc n)))))"
   2.222 -        apply (simp only: s of_int_add)
   2.223 -        apply (rule real_is_int_add)
   2.224 -        apply (simp add: neg1)
   2.225 -        apply (simp only: rn)
   2.226 -        done
   2.227 -    qed
   2.228 -  }
   2.229 -  note Abs_Bin = this
   2.230 -  {
   2.231 -    fix x :: int
   2.232 -    have "? u. x = u"
   2.233 -      apply (rule exI[where x = "x"])
   2.234 -      apply (simp)
   2.235 -      done
   2.236 -  }
   2.237 -  then obtain u::int where "x = u" by auto
   2.238 -  with Abs_Bin show ?thesis by auto
   2.239 -qed
   2.240 -
   2.241 -lemma int_of_real_0[simp]: "int_of_real (0::real) = (0::int)"
   2.242 -by (simp add: int_of_real_def)
   2.243 -
   2.244 -lemma int_of_real_1[simp]: "int_of_real (1::real) = (1::int)"
   2.245 -proof -
   2.246 -  have 1: "(1::real) = real (1::int)" by auto
   2.247 -  show ?thesis by (simp only: 1 int_of_real_real)
   2.248 -qed
   2.249 -
   2.250 -lemma int_of_real_number_of[simp]: "int_of_real (number_of b) = number_of b"
   2.251 -proof -
   2.252 -  have "real_is_int (number_of b)" by simp
   2.253 -  then have uu: "?! u::int. number_of b = real u" by (auto simp add: real_is_int_rep)
   2.254 -  then obtain u::int where u:"number_of b = real u" by auto
   2.255 -  have "number_of b = real ((number_of b)::int)"
   2.256 -    by (simp add: number_of_eq real_of_int_def)
   2.257 -  have ub: "number_of b = real ((number_of b)::int)"
   2.258 -    by (simp add: number_of_eq real_of_int_def)
   2.259 -  from uu u ub have unb: "u = number_of b"
   2.260 -    by blast
   2.261 -  have "int_of_real (number_of b) = u" by (simp add: u)
   2.262 -  with unb show ?thesis by simp
   2.263 -qed
   2.264 -
   2.265 -lemma float_transfer_even: "even a \<Longrightarrow> float (a, b) = float (a div 2, b+1)"
   2.266 -  apply (subst float_transfer[where a="a" and b="b" and c="-1", simplified])
   2.267 -  apply (simp_all add: pow2_def even_def real_is_int_def algebra_simps)
   2.268 -  apply (auto)
   2.269 -proof -
   2.270 -  fix q::int
   2.271 -  have a:"b - (-1\<Colon>int) = (1\<Colon>int) + b" by arith
   2.272 -  show "(float (q, (b - (-1\<Colon>int)))) = (float (q, ((1\<Colon>int) + b)))"
   2.273 -    by (simp add: a)
   2.274 -qed
   2.275 -
   2.276 -lemma int_div_zdiv: "int (a div b) = (int a) div (int b)"
   2.277 -by (rule zdiv_int)
   2.278 -
   2.279 -lemma int_mod_zmod: "int (a mod b) = (int a) mod (int b)"
   2.280 -by (rule zmod_int)
   2.281 -
   2.282 -lemma abs_div_2_less: "a \<noteq> 0 \<Longrightarrow> a \<noteq> -1 \<Longrightarrow> abs((a::int) div 2) < abs a"
   2.283 -by arith
   2.284 -
   2.285 -function norm_float :: "int \<Rightarrow> int \<Rightarrow> int \<times> int" where
   2.286 -  "norm_float a b = (if a \<noteq> 0 \<and> even a then norm_float (a div 2) (b + 1)
   2.287 -    else if a = 0 then (0, 0) else (a, b))"
   2.288 -by auto
   2.289 -
   2.290 -termination by (relation "measure (nat o abs o fst)")
   2.291 -  (auto intro: abs_div_2_less)
   2.292 -
   2.293 -lemma norm_float: "float x = float (split norm_float x)"
   2.294 -proof -
   2.295 -  {
   2.296 -    fix a b :: int
   2.297 -    have norm_float_pair: "float (a, b) = float (norm_float a b)"
   2.298 -    proof (induct a b rule: norm_float.induct)
   2.299 -      case (1 u v)
   2.300 -      show ?case
   2.301 -      proof cases
   2.302 -        assume u: "u \<noteq> 0 \<and> even u"
   2.303 -        with prems have ind: "float (u div 2, v + 1) = float (norm_float (u div 2) (v + 1))" by auto
   2.304 -        with u have "float (u,v) = float (u div 2, v+1)" by (simp add: float_transfer_even)
   2.305 -        then show ?thesis
   2.306 -          apply (subst norm_float.simps)
   2.307 -          apply (simp add: ind)
   2.308 -          done
   2.309 -      next
   2.310 -        assume "~(u \<noteq> 0 \<and> even u)"
   2.311 -        then show ?thesis
   2.312 -          by (simp add: prems float_def)
   2.313 -      qed
   2.314 -    qed
   2.315 -  }
   2.316 -  note helper = this
   2.317 -  have "? a b. x = (a,b)" by auto
   2.318 -  then obtain a b where "x = (a, b)" by blast
   2.319 -  then show ?thesis by (simp add: helper)
   2.320 -qed
   2.321 -
   2.322 -lemma float_add_l0: "float (0, e) + x = x"
   2.323 -  by (simp add: float_def)
   2.324 -
   2.325 -lemma float_add_r0: "x + float (0, e) = x"
   2.326 -  by (simp add: float_def)
   2.327 -
   2.328 -lemma float_add:
   2.329 -  "float (a1, e1) + float (a2, e2) =
   2.330 -  (if e1<=e2 then float (a1+a2*2^(nat(e2-e1)), e1)
   2.331 -  else float (a1*2^(nat (e1-e2))+a2, e2))"
   2.332 -  apply (simp add: float_def algebra_simps)
   2.333 -  apply (auto simp add: pow2_int[symmetric] pow2_add[symmetric])
   2.334 -  done
   2.335 -
   2.336 -lemma float_add_assoc1:
   2.337 -  "(x + float (y1, e1)) + float (y2, e2) = (float (y1, e1) + float (y2, e2)) + x"
   2.338 -  by simp
   2.339 -
   2.340 -lemma float_add_assoc2:
   2.341 -  "(float (y1, e1) + x) + float (y2, e2) = (float (y1, e1) + float (y2, e2)) + x"
   2.342 -  by simp
   2.343 -
   2.344 -lemma float_add_assoc3:
   2.345 -  "float (y1, e1) + (x + float (y2, e2)) = (float (y1, e1) + float (y2, e2)) + x"
   2.346 -  by simp
   2.347 -
   2.348 -lemma float_add_assoc4:
   2.349 -  "float (y1, e1) + (float (y2, e2) + x) = (float (y1, e1) + float (y2, e2)) + x"
   2.350 -  by simp
   2.351 -
   2.352 -lemma float_mult_l0: "float (0, e) * x = float (0, 0)"
   2.353 -  by (simp add: float_def)
   2.354 -
   2.355 -lemma float_mult_r0: "x * float (0, e) = float (0, 0)"
   2.356 -  by (simp add: float_def)
   2.357 -
   2.358 -definition 
   2.359 -  lbound :: "real \<Rightarrow> real"
   2.360 -where
   2.361 -  "lbound x = min 0 x"
   2.362 -
   2.363 -definition
   2.364 -  ubound :: "real \<Rightarrow> real"
   2.365 -where
   2.366 -  "ubound x = max 0 x"
   2.367 -
   2.368 -lemma lbound: "lbound x \<le> x"   
   2.369 -  by (simp add: lbound_def)
   2.370 -
   2.371 -lemma ubound: "x \<le> ubound x"
   2.372 -  by (simp add: ubound_def)
   2.373 -
   2.374 -lemma float_mult:
   2.375 -  "float (a1, e1) * float (a2, e2) =
   2.376 -  (float (a1 * a2, e1 + e2))"
   2.377 -  by (simp add: float_def pow2_add)
   2.378 -
   2.379 -lemma float_minus:
   2.380 -  "- (float (a,b)) = float (-a, b)"
   2.381 -  by (simp add: float_def)
   2.382 -
   2.383 -lemma zero_less_pow2:
   2.384 +lemma zero_less_pow2[simp]:
   2.385    "0 < pow2 x"
   2.386  proof -
   2.387    {
   2.388 @@ -387,182 +175,1258 @@
   2.389      done
   2.390  qed
   2.391  
   2.392 +lemma normfloat_imp_odd_or_zero: "normfloat f = Float a b \<Longrightarrow> odd a \<or> (a = 0 \<and> b = 0)"
   2.393 +proof (induct f rule: normfloat.induct)
   2.394 +  case (1 u v)
   2.395 +  from 1 have ab: "normfloat (Float u v) = Float a b" by auto
   2.396 +  {
   2.397 +    assume eu: "even u"
   2.398 +    assume z: "u \<noteq> 0"
   2.399 +    have "normfloat (Float u v) = normfloat (Float (u div 2) (v + 1))"
   2.400 +      apply (subst normfloat.simps)
   2.401 +      by (simp add: eu z)
   2.402 +    with ab have "normfloat (Float (u div 2) (v + 1)) = Float a b" by simp
   2.403 +    with 1 eu z have ?case by auto
   2.404 +  }
   2.405 +  note case1 = this
   2.406 +  {
   2.407 +    assume "odd u \<or> u = 0"
   2.408 +    then have ou: "\<not> (u \<noteq> 0 \<and> even u)" by auto
   2.409 +    have "normfloat (Float u v) = (if u = 0 then Float 0 0 else Float u v)"
   2.410 +      apply (subst normfloat.simps)
   2.411 +      apply (simp add: ou)
   2.412 +      done
   2.413 +    with ab have "Float a b = (if u = 0 then Float 0 0 else Float u v)" by auto
   2.414 +    then have ?case
   2.415 +      apply (case_tac "u=0")
   2.416 +      apply (auto)
   2.417 +      by (insert ou, auto)
   2.418 +  }
   2.419 +  note case2 = this
   2.420 +  show ?case
   2.421 +    apply (case_tac "odd u \<or> u = 0")
   2.422 +    apply (rule case2)
   2.423 +    apply simp
   2.424 +    apply (rule case1)
   2.425 +    apply auto
   2.426 +    done
   2.427 +qed
   2.428 +
   2.429 +lemma float_eq_odd_helper: 
   2.430 +  assumes odd: "odd a'"
   2.431 +  and floateq: "Ifloat (Float a b) = Ifloat (Float a' b')"
   2.432 +  shows "b \<le> b'"
   2.433 +proof - 
   2.434 +  {
   2.435 +    assume bcmp: "b > b'"
   2.436 +    from floateq have eq: "real a * pow2 b = real a' * pow2 b'" by simp
   2.437 +    {
   2.438 +      fix x y z :: real
   2.439 +      assume "y \<noteq> 0"
   2.440 +      then have "(x * inverse y = z) = (x = z * y)"
   2.441 +	by auto
   2.442 +    }
   2.443 +    note inverse = this
   2.444 +    have eq': "real a * (pow2 (b - b')) = real a'"
   2.445 +      apply (subst diff_int_def)
   2.446 +      apply (subst pow2_add)
   2.447 +      apply (subst pow2_neg[where x = "-b'"])
   2.448 +      apply simp
   2.449 +      apply (subst mult_assoc[symmetric])
   2.450 +      apply (subst inverse)
   2.451 +      apply (simp_all add: eq)
   2.452 +      done
   2.453 +    have "\<exists> z > 0. pow2 (b-b') = 2^z"
   2.454 +      apply (rule exI[where x="nat (b - b')"])
   2.455 +      apply (auto)
   2.456 +      apply (insert bcmp)
   2.457 +      apply simp
   2.458 +      apply (subst pow2_int[symmetric])
   2.459 +      apply auto
   2.460 +      done
   2.461 +    then obtain z where z: "z > 0 \<and> pow2 (b-b') = 2^z" by auto
   2.462 +    with eq' have "real a * 2^z = real a'"
   2.463 +      by auto
   2.464 +    then have "real a * real ((2::int)^z) = real a'"
   2.465 +      by auto
   2.466 +    then have "real (a * 2^z) = real a'"
   2.467 +      apply (subst real_of_int_mult)
   2.468 +      apply simp
   2.469 +      done
   2.470 +    then have a'_rep: "a * 2^z = a'" by arith
   2.471 +    then have "a' = a*2^z" by simp
   2.472 +    with z have "even a'" by simp
   2.473 +    with odd have False by auto
   2.474 +  }
   2.475 +  then show ?thesis by arith
   2.476 +qed
   2.477 +
   2.478 +lemma float_eq_odd: 
   2.479 +  assumes odd1: "odd a"
   2.480 +  and odd2: "odd a'"
   2.481 +  and floateq: "Ifloat (Float a b) = Ifloat (Float a' b')"
   2.482 +  shows "a = a' \<and> b = b'"
   2.483 +proof -
   2.484 +  from 
   2.485 +     float_eq_odd_helper[OF odd2 floateq] 
   2.486 +     float_eq_odd_helper[OF odd1 floateq[symmetric]]
   2.487 +  have beq: "b = b'"  by arith
   2.488 +  with floateq show ?thesis by auto
   2.489 +qed
   2.490 +
   2.491 +theorem normfloat_unique:
   2.492 +  assumes Ifloat_eq: "Ifloat f = Ifloat g"
   2.493 +  shows "normfloat f = normfloat g"
   2.494 +proof - 
   2.495 +  from float_split[of "normfloat f"] obtain a b where normf:"normfloat f = Float a b" by auto
   2.496 +  from float_split[of "normfloat g"] obtain a' b' where normg:"normfloat g = Float a' b'" by auto
   2.497 +  have "Ifloat (normfloat f) = Ifloat (normfloat g)"
   2.498 +    by (simp add: Ifloat_eq)
   2.499 +  then have float_eq: "Ifloat (Float a b) = Ifloat (Float a' b')"
   2.500 +    by (simp add: normf normg)
   2.501 +  have ab: "odd a \<or> (a = 0 \<and> b = 0)" by (rule normfloat_imp_odd_or_zero[OF normf])
   2.502 +  have ab': "odd a' \<or> (a' = 0 \<and> b' = 0)" by (rule normfloat_imp_odd_or_zero[OF normg])
   2.503 +  {
   2.504 +    assume odd: "odd a"
   2.505 +    then have "a \<noteq> 0" by (simp add: even_def, arith)
   2.506 +    with float_eq have "a' \<noteq> 0" by auto
   2.507 +    with ab' have "odd a'" by simp
   2.508 +    from odd this float_eq have "a = a' \<and> b = b'" by (rule float_eq_odd)
   2.509 +  }
   2.510 +  note odd_case = this
   2.511 +  {
   2.512 +    assume even: "even a"
   2.513 +    with ab have a0: "a = 0" by simp
   2.514 +    with float_eq have a0': "a' = 0" by auto 
   2.515 +    from a0 a0' ab ab' have "a = a' \<and> b = b'" by auto
   2.516 +  }
   2.517 +  note even_case = this
   2.518 +  from odd_case even_case show ?thesis
   2.519 +    apply (simp add: normf normg)
   2.520 +    apply (case_tac "even a")
   2.521 +    apply auto
   2.522 +    done
   2.523 +qed
   2.524 +
   2.525 +instantiation float :: plus begin
   2.526 +fun plus_float where
   2.527 +[simp del]: "(Float a_m a_e) + (Float b_m b_e) = 
   2.528 +     (if a_e \<le> b_e then Float (a_m + b_m * 2^(nat(b_e - a_e))) a_e 
   2.529 +                   else Float (a_m * 2^(nat (a_e - b_e)) + b_m) b_e)"
   2.530 +instance ..
   2.531 +end
   2.532 +
   2.533 +instantiation float :: uminus begin
   2.534 +fun uminus_float where [simp del]: "uminus_float (Float m e) = Float (-m) e"
   2.535 +instance ..
   2.536 +end
   2.537 +
   2.538 +instantiation float :: minus begin
   2.539 +fun minus_float where [simp del]: "(z::float) - w = z + (- w)"
   2.540 +instance ..
   2.541 +end
   2.542 +
   2.543 +instantiation float :: times begin
   2.544 +fun times_float where [simp del]: "(Float a_m a_e) * (Float b_m b_e) = Float (a_m * b_m) (a_e + b_e)"
   2.545 +instance ..
   2.546 +end
   2.547 +
   2.548 +fun float_pprt :: "float \<Rightarrow> float" where
   2.549 +"float_pprt (Float a e) = (if 0 <= a then (Float a e) else 0)"
   2.550 +
   2.551 +fun float_nprt :: "float \<Rightarrow> float" where
   2.552 +"float_nprt (Float a e) = (if 0 <= a then 0 else (Float a e))" 
   2.553 +
   2.554 +instantiation float :: ord begin
   2.555 +definition le_float_def: "z \<le> w \<equiv> Ifloat z \<le> Ifloat w"
   2.556 +definition less_float_def: "z < w \<equiv> Ifloat z < Ifloat w"
   2.557 +instance ..
   2.558 +end
   2.559 +
   2.560 +lemma Ifloat_add[simp]: "Ifloat (a + b) = Ifloat a + Ifloat b"
   2.561 +  by (cases a, cases b, simp add: algebra_simps plus_float.simps, 
   2.562 +      auto simp add: pow2_int[symmetric] pow2_add[symmetric])
   2.563 +
   2.564 +lemma Ifloat_minus[simp]: "Ifloat (- a) = - Ifloat a"
   2.565 +  by (cases a, simp add: uminus_float.simps)
   2.566 +
   2.567 +lemma Ifloat_sub[simp]: "Ifloat (a - b) = Ifloat a - Ifloat b" 
   2.568 +  by (cases a, cases b, simp add: minus_float.simps)
   2.569 +
   2.570 +lemma Ifloat_mult[simp]: "Ifloat (a*b) = Ifloat a * Ifloat b"
   2.571 +  by (cases a, cases b, simp add: times_float.simps pow2_add)
   2.572 +
   2.573 +lemma Ifloat_0[simp]: "Ifloat 0 = 0"
   2.574 +  by (auto simp add: zero_float_def float_zero)
   2.575 +
   2.576 +lemma Ifloat_1[simp]: "Ifloat 1 = 1"
   2.577 +  by (auto simp add: one_float_def)
   2.578 +
   2.579  lemma zero_le_float:
   2.580 -  "(0 <= float (a,b)) = (0 <= a)"
   2.581 -  apply (auto simp add: float_def)
   2.582 -  apply (auto simp add: zero_le_mult_iff zero_less_pow2)
   2.583 +  "(0 <= Ifloat (Float a b)) = (0 <= a)"
   2.584 +  apply auto
   2.585 +  apply (auto simp add: zero_le_mult_iff)
   2.586    apply (insert zero_less_pow2[of b])
   2.587    apply (simp_all)
   2.588    done
   2.589  
   2.590  lemma float_le_zero:
   2.591 -  "(float (a,b) <= 0) = (a <= 0)"
   2.592 -  apply (auto simp add: float_def)
   2.593 +  "(Ifloat (Float a b) <= 0) = (a <= 0)"
   2.594 +  apply auto
   2.595    apply (auto simp add: mult_le_0_iff)
   2.596    apply (insert zero_less_pow2[of b])
   2.597    apply auto
   2.598    done
   2.599  
   2.600 -lemma float_abs:
   2.601 -  "abs (float (a,b)) = (if 0 <= a then (float (a,b)) else (float (-a,b)))"
   2.602 -  apply (auto simp add: abs_if)
   2.603 -  apply (simp_all add: zero_le_float[symmetric, of a b] float_minus)
   2.604 +declare Ifloat.simps[simp del]
   2.605 +
   2.606 +lemma Ifloat_pprt[simp]: "Ifloat (float_pprt a) = pprt (Ifloat a)"
   2.607 +  by (cases a, auto simp add: float_pprt.simps zero_le_float float_le_zero float_zero)
   2.608 +
   2.609 +lemma Ifloat_nprt[simp]: "Ifloat (float_nprt a) = nprt (Ifloat a)"
   2.610 +  by (cases a,  auto simp add: float_nprt.simps zero_le_float float_le_zero float_zero)
   2.611 +
   2.612 +instance float :: ab_semigroup_add
   2.613 +proof (intro_classes)
   2.614 +  fix a b c :: float
   2.615 +  show "a + b + c = a + (b + c)"
   2.616 +    by (cases a, cases b, cases c, auto simp add: algebra_simps power_add[symmetric] plus_float.simps)
   2.617 +next
   2.618 +  fix a b :: float
   2.619 +  show "a + b = b + a"
   2.620 +    by (cases a, cases b, simp add: plus_float.simps)
   2.621 +qed
   2.622 +
   2.623 +instance float :: comm_monoid_mult
   2.624 +proof (intro_classes)
   2.625 +  fix a b c :: float
   2.626 +  show "a * b * c = a * (b * c)"
   2.627 +    by (cases a, cases b, cases c, simp add: times_float.simps)
   2.628 +next
   2.629 +  fix a b :: float
   2.630 +  show "a * b = b * a"
   2.631 +    by (cases a, cases b, simp add: times_float.simps)
   2.632 +next
   2.633 +  fix a :: float
   2.634 +  show "1 * a = a"
   2.635 +    by (cases a, simp add: times_float.simps one_float_def)
   2.636 +qed
   2.637 +
   2.638 +(* Floats do NOT form a cancel_semigroup_add: *)
   2.639 +lemma "0 + Float 0 1 = 0 + Float 0 2"
   2.640 +  by (simp add: plus_float.simps zero_float_def)
   2.641 +
   2.642 +instance float :: comm_semiring
   2.643 +proof (intro_classes)
   2.644 +  fix a b c :: float
   2.645 +  show "(a + b) * c = a * c + b * c"
   2.646 +    by (cases a, cases b, cases c, simp, simp add: plus_float.simps times_float.simps algebra_simps)
   2.647 +qed
   2.648 +
   2.649 +(* Floats do NOT form an order, because "(x < y) = (x <= y & x <> y)" does NOT hold *)
   2.650 +
   2.651 +instance float :: zero_neq_one
   2.652 +proof (intro_classes)
   2.653 +  show "(0::float) \<noteq> 1"
   2.654 +    by (simp add: zero_float_def one_float_def)
   2.655 +qed
   2.656 +
   2.657 +lemma float_le_simp: "((x::float) \<le> y) = (0 \<le> y - x)"
   2.658 +  by (auto simp add: le_float_def)
   2.659 +
   2.660 +lemma float_less_simp: "((x::float) < y) = (0 < y - x)"
   2.661 +  by (auto simp add: less_float_def)
   2.662 +
   2.663 +lemma Ifloat_min: "Ifloat (min x y) = min (Ifloat x) (Ifloat y)" unfolding min_def le_float_def by auto
   2.664 +lemma Ifloat_max: "Ifloat (max a b) = max (Ifloat a) (Ifloat b)" unfolding max_def le_float_def by auto
   2.665 +
   2.666 +instantiation float :: power begin 
   2.667 +fun power_float where [simp del]: "(Float m e) ^ n = Float (m ^ n) (e * int n)"
   2.668 +instance ..
   2.669 +end
   2.670 +
   2.671 +instance float :: recpower
   2.672 +proof (intro_classes)
   2.673 +  fix a :: float show "a ^ 0 = 1" by (cases a, auto simp add: power_float.simps one_float_def)
   2.674 +next
   2.675 +  fix a :: float and n :: nat show "a ^ (Suc n) = a * a ^ n" 
   2.676 +  by (cases a, auto simp add: power_float.simps times_float.simps algebra_simps)
   2.677 +qed
   2.678 +
   2.679 +lemma float_power: "Ifloat (x ^ n) = (Ifloat x) ^ n"
   2.680 +proof (cases x)
   2.681 +  case (Float m e)
   2.682 +  
   2.683 +  have "pow2 e ^ n = pow2 (e * int n)"
   2.684 +  proof (cases "e >= 0")
   2.685 +    case True hence e_nat: "e = int (nat e)" by auto
   2.686 +    hence "pow2 e ^ n = (2 ^ nat e) ^ n" using pow2_int[of "nat e"] by auto
   2.687 +    thus ?thesis unfolding power_mult[symmetric] unfolding pow2_int[symmetric] int_mult e_nat[symmetric] .
   2.688 +  next
   2.689 +    case False hence e_minus: "-e = int (nat (-e))" by auto
   2.690 +    hence "pow2 (-e) ^ n = (2 ^ nat (-e)) ^ n" using pow2_int[of "nat (-e)"] by auto
   2.691 +    hence "pow2 (-e) ^ n = pow2 ((-e) * int n)" unfolding power_mult[symmetric] unfolding pow2_int[symmetric] int_mult e_minus[symmetric] zmult_zminus .
   2.692 +    thus ?thesis unfolding pow2_neg[of "-e"] pow2_neg[of "-e * int n"] unfolding zmult_zminus zminus_zminus nonzero_power_inverse[OF pow2_neq_zero, symmetric]
   2.693 +      using nonzero_inverse_eq_imp_eq[OF _ pow2_neq_zero pow2_neq_zero] by auto
   2.694 +  qed
   2.695 +  thus ?thesis by (auto simp add: Float power_mult_distrib Ifloat.simps power_float.simps)
   2.696 +qed
   2.697 +
   2.698 +lemma zero_le_pow2[simp]: "0 \<le> pow2 s"
   2.699 +  apply (subgoal_tac "0 < pow2 s")
   2.700 +  apply (auto simp only:)
   2.701 +  apply auto
   2.702    done
   2.703  
   2.704 -lemma float_zero:
   2.705 -  "float (0, b) = 0"
   2.706 -  by (simp add: float_def)
   2.707 -
   2.708 -lemma float_pprt:
   2.709 -  "pprt (float (a, b)) = (if 0 <= a then (float (a,b)) else (float (0, b)))"
   2.710 -  by (auto simp add: zero_le_float float_le_zero float_zero)
   2.711 -
   2.712 -lemma pprt_lbound: "pprt (lbound x) = float (0, 0)"
   2.713 -  apply (simp add: float_def)
   2.714 -  apply (rule pprt_eq_0)
   2.715 -  apply (simp add: lbound_def)
   2.716 +lemma pow2_less_0_eq_False[simp]: "(pow2 s < 0) = False"
   2.717 +  apply auto
   2.718 +  apply (subgoal_tac "0 \<le> pow2 s")
   2.719 +  apply simp
   2.720 +  apply simp
   2.721    done
   2.722  
   2.723 -lemma nprt_ubound: "nprt (ubound x) = float (0, 0)"
   2.724 -  apply (simp add: float_def)
   2.725 -  apply (rule nprt_eq_0)
   2.726 -  apply (simp add: ubound_def)
   2.727 +lemma pow2_le_0_eq_False[simp]: "(pow2 s \<le> 0) = False"
   2.728 +  apply auto
   2.729 +  apply (subgoal_tac "0 < pow2 s")
   2.730 +  apply simp
   2.731 +  apply simp
   2.732    done
   2.733  
   2.734 -lemma float_nprt:
   2.735 -  "nprt (float (a, b)) = (if 0 <= a then (float (0,b)) else (float (a, b)))"
   2.736 -  by (auto simp add: zero_le_float float_le_zero float_zero)
   2.737 -
   2.738 -lemma norm_0_1: "(0::_::number_ring) = Numeral0 & (1::_::number_ring) = Numeral1"
   2.739 +lemma float_pos_m_pos: "0 < Float m e \<Longrightarrow> 0 < m"
   2.740 +  unfolding less_float_def Ifloat.simps Ifloat_0 zero_less_mult_iff
   2.741    by auto
   2.742  
   2.743 -lemma add_left_zero: "0 + a = (a::'a::comm_monoid_add)"
   2.744 -  by simp
   2.745 +lemma float_pos_less1_e_neg: assumes "0 < Float m e" and "Float m e < 1" shows "e < 0"
   2.746 +proof -
   2.747 +  have "0 < m" using float_pos_m_pos `0 < Float m e` by auto
   2.748 +  hence "0 \<le> real m" and "1 \<le> real m" by auto
   2.749 +  
   2.750 +  show "e < 0"
   2.751 +  proof (rule ccontr)
   2.752 +    assume "\<not> e < 0" hence "0 \<le> e" by auto
   2.753 +    hence "1 \<le> pow2 e" unfolding pow2_def by auto
   2.754 +    from mult_mono[OF `1 \<le> real m` this `0 \<le> real m`]
   2.755 +    have "1 \<le> Float m e" by (simp add: le_float_def Ifloat.simps)
   2.756 +    thus False using `Float m e < 1` unfolding less_float_def le_float_def by auto
   2.757 +  qed
   2.758 +qed
   2.759 +
   2.760 +lemma float_less1_mantissa_bound: assumes "0 < Float m e" "Float m e < 1" shows "m < 2^(nat (-e))"
   2.761 +proof -
   2.762 +  have "e < 0" using float_pos_less1_e_neg assms by auto
   2.763 +  have "\<And>x. (0::real) < 2^x" by auto
   2.764 +  have "real m < 2^(nat (-e))" using `Float m e < 1`
   2.765 +    unfolding less_float_def Ifloat_neg_exp[OF `e < 0`] Ifloat_1
   2.766 +          real_mult_less_iff1[of _ _ 1, OF `0 < 2^(nat (-e))`, symmetric] 
   2.767 +          real_mult_assoc by auto
   2.768 +  thus ?thesis unfolding real_of_int_less_iff[symmetric] by auto
   2.769 +qed
   2.770 +
   2.771 +function bitlen :: "int \<Rightarrow> int" where
   2.772 +"bitlen 0 = 0" | 
   2.773 +"bitlen -1 = 1" | 
   2.774 +"0 < x \<Longrightarrow> bitlen x = 1 + (bitlen (x div 2))" | 
   2.775 +"x < -1 \<Longrightarrow> bitlen x = 1 + (bitlen (x div 2))"
   2.776 +  apply (case_tac "x = 0 \<or> x = -1 \<or> x < -1 \<or> x > 0")
   2.777 +  apply auto
   2.778 +  done
   2.779 +termination by (relation "measure (nat o abs)", auto)
   2.780 +
   2.781 +lemma bitlen_ge0: "0 \<le> bitlen x" by (induct x rule: bitlen.induct, auto)
   2.782 +lemma bitlen_ge1: "x \<noteq> 0 \<Longrightarrow> 1 \<le> bitlen x" by (induct x rule: bitlen.induct, auto simp add: bitlen_ge0)
   2.783 +
   2.784 +lemma bitlen_bounds': assumes "0 < x" shows "2^nat (bitlen x - 1) \<le> x \<and> x + 1 \<le> 2^nat (bitlen x)" (is "?P x")
   2.785 +  using `0 < x`
   2.786 +proof (induct x rule: bitlen.induct)
   2.787 +  fix x
   2.788 +  assume "0 < x" and hyp: "0 < x div 2 \<Longrightarrow> ?P (x div 2)" hence "0 \<le> x" and "x \<noteq> 0" by auto
   2.789 +  { fix x have "0 \<le> 1 + bitlen x" using bitlen_ge0[of x] by auto } note gt0_pls1 = this
   2.790 +
   2.791 +  have "0 < (2::int)" by auto
   2.792  
   2.793 -lemma add_right_zero: "a + 0 = (a::'a::comm_monoid_add)"
   2.794 -  by simp
   2.795 +  show "?P x"
   2.796 +  proof (cases "x = 1")
   2.797 +    case True show "?P x" unfolding True by auto
   2.798 +  next
   2.799 +    case False hence "2 \<le> x" using `0 < x` `x \<noteq> 1` by auto
   2.800 +    hence "2 div 2 \<le> x div 2" by (rule zdiv_mono1, auto)
   2.801 +    hence "0 < x div 2" and "x div 2 \<noteq> 0" by auto
   2.802 +    hence bitlen_s1_ge0: "0 \<le> bitlen (x div 2) - 1" using bitlen_ge1[OF `x div 2 \<noteq> 0`] by auto
   2.803  
   2.804 -lemma mult_left_one: "1 * a = (a::'a::semiring_1)"
   2.805 -  by simp
   2.806 +    { from hyp[OF `0 < x div 2`]
   2.807 +      have "2 ^ nat (bitlen (x div 2) - 1) \<le> x div 2" by auto
   2.808 +      hence "2 ^ nat (bitlen (x div 2) - 1) * 2 \<le> x div 2 * 2" by (rule mult_right_mono, auto)
   2.809 +      also have "\<dots> \<le> x" using `0 < x` by auto
   2.810 +      finally have "2^nat (1 + bitlen (x div 2) - 1) \<le> x" unfolding power_Suc2[symmetric] Suc_nat_eq_nat_zadd1[OF bitlen_s1_ge0] by auto
   2.811 +    } moreover
   2.812 +    { have "x + 1 \<le> x - x mod 2 + 2"
   2.813 +      proof -
   2.814 +	have "x mod 2 < 2" using `0 < x` by auto
   2.815 + 	hence "x < x - x mod 2 +  2" unfolding algebra_simps by auto
   2.816 +	thus ?thesis by auto
   2.817 +      qed
   2.818 +      also have "x - x mod 2 + 2 = (x div 2 + 1) * 2" unfolding algebra_simps using `0 < x` zdiv_zmod_equality2[of x 2 0] by auto
   2.819 +      also have "\<dots> \<le> 2^nat (bitlen (x div 2)) * 2" using hyp[OF `0 < x div 2`, THEN conjunct2] by (rule mult_right_mono, auto)
   2.820 +      also have "\<dots> = 2^(1 + nat (bitlen (x div 2)))" unfolding power_Suc2[symmetric] by auto
   2.821 +      finally have "x + 1 \<le> 2^(1 + nat (bitlen (x div 2)))" .
   2.822 +    }
   2.823 +    ultimately show ?thesis
   2.824 +      unfolding bitlen.simps(3)[OF `0 < x`] nat_add_distrib[OF zero_le_one bitlen_ge0]
   2.825 +      unfolding add_commute nat_add_distrib[OF zero_le_one gt0_pls1]
   2.826 +      by auto
   2.827 +  qed
   2.828 +next
   2.829 +  fix x :: int assume "x < -1" and "0 < x" hence False by auto
   2.830 +  thus "?P x" by auto
   2.831 +qed auto
   2.832 +
   2.833 +lemma bitlen_bounds: assumes "0 < x" shows "2^nat (bitlen x - 1) \<le> x \<and> x < 2^nat (bitlen x)"
   2.834 +  using bitlen_bounds'[OF `0<x`] by auto
   2.835 +
   2.836 +lemma bitlen_div: assumes "0 < m" shows "1 \<le> real m / 2^nat (bitlen m - 1)" and "real m / 2^nat (bitlen m - 1) < 2"
   2.837 +proof -
   2.838 +  let ?B = "2^nat(bitlen m - 1)"
   2.839 +
   2.840 +  have "?B \<le> m" using bitlen_bounds[OF `0 <m`] ..
   2.841 +  hence "1 * ?B \<le> real m" unfolding real_of_int_le_iff[symmetric] by auto
   2.842 +  thus "1 \<le> real m / ?B" by auto
   2.843 +
   2.844 +  have "m \<noteq> 0" using assms by auto
   2.845 +  have "0 \<le> bitlen m - 1" using bitlen_ge1[OF `m \<noteq> 0`] by auto
   2.846  
   2.847 -lemma mult_right_one: "a * 1 = (a::'a::semiring_1)"
   2.848 -  by simp
   2.849 +  have "m < 2^nat(bitlen m)" using bitlen_bounds[OF `0 <m`] ..
   2.850 +  also have "\<dots> = 2^nat(bitlen m - 1 + 1)" using bitlen_ge1[OF `m \<noteq> 0`] by auto
   2.851 +  also have "\<dots> = ?B * 2" unfolding nat_add_distrib[OF `0 \<le> bitlen m - 1` zero_le_one] by auto
   2.852 +  finally have "real m < 2 * ?B" unfolding real_of_int_less_iff[symmetric] by auto
   2.853 +  hence "real m / ?B < 2 * ?B / ?B" by (rule divide_strict_right_mono, auto)
   2.854 +  thus "real m / ?B < 2" by auto
   2.855 +qed
   2.856 +
   2.857 +lemma float_gt1_scale: assumes "1 \<le> Float m e"
   2.858 +  shows "0 \<le> e + (bitlen m - 1)"
   2.859 +proof (cases "0 \<le> e")
   2.860 +  have "0 < Float m e" using assms unfolding less_float_def le_float_def by auto
   2.861 +  hence "0 < m" using float_pos_m_pos by auto
   2.862 +  hence "m \<noteq> 0" by auto
   2.863 +  case True with bitlen_ge1[OF `m \<noteq> 0`] show ?thesis by auto
   2.864 +next
   2.865 +  have "0 < Float m e" using assms unfolding less_float_def le_float_def by auto
   2.866 +  hence "0 < m" using float_pos_m_pos by auto
   2.867 +  hence "m \<noteq> 0" and "1 < (2::int)" by auto
   2.868 +  case False let ?S = "2^(nat (-e))"
   2.869 +  have "1 \<le> real m * inverse ?S" using assms unfolding le_float_def Ifloat_nge0_exp[OF False] by auto
   2.870 +  hence "1 * ?S \<le> real m * inverse ?S * ?S" by (rule mult_right_mono, auto)
   2.871 +  hence "?S \<le> real m" unfolding mult_assoc by auto
   2.872 +  hence "?S \<le> m" unfolding real_of_int_le_iff[symmetric] by auto
   2.873 +  from this bitlen_bounds[OF `0 < m`, THEN conjunct2]
   2.874 +  have "nat (-e) < (nat (bitlen m))" unfolding power_strict_increasing_iff[OF `1 < 2`, symmetric] by (rule order_le_less_trans)
   2.875 +  hence "-e < bitlen m" using False bitlen_ge0 by auto
   2.876 +  thus ?thesis by auto
   2.877 +qed
   2.878 +
   2.879 +lemma normalized_float: assumes "m \<noteq> 0" shows "Ifloat (Float m (- (bitlen m - 1))) = real m / 2^nat (bitlen m - 1)"
   2.880 +proof (cases "- (bitlen m - 1) = 0")
   2.881 +  case True show ?thesis unfolding Ifloat.simps pow2_def using True by auto
   2.882 +next
   2.883 +  case False hence P: "\<not> 0 \<le> - (bitlen m - 1)" using bitlen_ge1[OF `m \<noteq> 0`] by auto
   2.884 +  show ?thesis unfolding Ifloat_nge0_exp[OF P] real_divide_def by auto
   2.885 +qed
   2.886 +
   2.887 +lemma bitlen_Pls: "bitlen (Int.Pls) = Int.Pls" by (subst Pls_def, subst Pls_def, simp)
   2.888 +
   2.889 +lemma bitlen_Min: "bitlen (Int.Min) = Int.Bit1 Int.Pls" by (subst Min_def, simp add: Bit1_def) 
   2.890 +
   2.891 +lemma bitlen_B0: "bitlen (Int.Bit0 b) = (if iszero b then Int.Pls else Int.succ (bitlen b))"
   2.892 +  apply (auto simp add: iszero_def succ_def)
   2.893 +  apply (simp add: Bit0_def Pls_def)
   2.894 +  apply (subst Bit0_def)
   2.895 +  apply simp
   2.896 +  apply (subgoal_tac "0 < 2 * b \<or> 2 * b < -1")
   2.897 +  apply auto
   2.898 +  done
   2.899  
   2.900 -lemma int_pow_0: "(a::int)^(Numeral0) = 1"
   2.901 -  by simp
   2.902 +lemma bitlen_B1: "bitlen (Int.Bit1 b) = (if iszero (Int.succ b) then Int.Bit1 Int.Pls else Int.succ (bitlen b))"
   2.903 +proof -
   2.904 +  have h: "! x. (2*x + 1) div 2 = (x::int)"
   2.905 +    by arith    
   2.906 +  show ?thesis
   2.907 +    apply (auto simp add: iszero_def succ_def)
   2.908 +    apply (subst Bit1_def)+
   2.909 +    apply simp
   2.910 +    apply (subgoal_tac "2 * b + 1 = -1")
   2.911 +    apply (simp only:)
   2.912 +    apply simp_all
   2.913 +    apply (subst Bit1_def)
   2.914 +    apply simp
   2.915 +    apply (subgoal_tac "0 < 2 * b + 1 \<or> 2 * b + 1 < -1")
   2.916 +    apply (auto simp add: h)
   2.917 +    done
   2.918 +qed
   2.919 +
   2.920 +lemma bitlen_number_of: "bitlen (number_of w) = number_of (bitlen w)"
   2.921 +  by (simp add: number_of_is_id)
   2.922  
   2.923 -lemma int_pow_1: "(a::int)^(Numeral1) = a"
   2.924 -  by simp
   2.925 +lemma [code]: "bitlen x = 
   2.926 +     (if x = 0  then 0 
   2.927 + else if x = -1 then 1 
   2.928 +                else (1 + (bitlen (x div 2))))"
   2.929 +  by (cases "x = 0 \<or> x = -1 \<or> 0 < x") auto
   2.930 +
   2.931 +definition lapprox_posrat :: "nat \<Rightarrow> int \<Rightarrow> int \<Rightarrow> float"
   2.932 +where
   2.933 +  "lapprox_posrat prec x y = 
   2.934 +   (let 
   2.935 +       l = nat (int prec + bitlen y - bitlen x) ;
   2.936 +       d = (x * 2^l) div y
   2.937 +    in normfloat (Float d (- (int l))))"
   2.938 +
   2.939 +lemma pow2_minus: "pow2 (-x) = inverse (pow2 x)"
   2.940 +  unfolding pow2_neg[of "-x"] by auto
   2.941 +
   2.942 +lemma lapprox_posrat: 
   2.943 +  assumes x: "0 \<le> x"
   2.944 +  and y: "0 < y"
   2.945 +  shows "Ifloat (lapprox_posrat prec x y) \<le> real x / real y"
   2.946 +proof -
   2.947 +  let ?l = "nat (int prec + bitlen y - bitlen x)"
   2.948 +  
   2.949 +  have "real (x * 2^?l div y) * inverse (2^?l) \<le> (real (x * 2^?l) / real y) * inverse (2^?l)" 
   2.950 +    by (rule mult_right_mono, fact real_of_int_div4, simp)
   2.951 +  also have "\<dots> \<le> (real x / real y) * 2^?l * inverse (2^?l)" by auto
   2.952 +  finally have "real (x * 2^?l div y) * inverse (2^?l) \<le> real x / real y" unfolding real_mult_assoc by auto
   2.953 +  thus ?thesis unfolding lapprox_posrat_def Let_def normfloat Ifloat.simps
   2.954 +    unfolding pow2_minus pow2_int minus_minus .
   2.955 +qed
   2.956  
   2.957 -lemma zero_eq_Numeral0_nring: "(0::'a::number_ring) = Numeral0"
   2.958 -  by simp
   2.959 +lemma real_of_int_div_mult: 
   2.960 +  fixes x y c :: int assumes "0 < y" and "0 < c"
   2.961 +  shows "real (x div y) \<le> real (x * c div y) * inverse (real c)"
   2.962 +proof -
   2.963 +  have "c * (x div y) + 0 \<le> c * x div y" unfolding zdiv_zmult1_eq[of c x y]
   2.964 +    by (rule zadd_left_mono, 
   2.965 +        auto intro!: mult_nonneg_nonneg 
   2.966 +             simp add: pos_imp_zdiv_nonneg_iff[OF `0 < y`] `0 < c`[THEN less_imp_le] pos_mod_sign[OF `0 < y`])
   2.967 +  hence "real (x div y) * real c \<le> real (x * c div y)" 
   2.968 +    unfolding real_of_int_mult[symmetric] real_of_int_le_iff zmult_commute by auto
   2.969 +  hence "real (x div y) * real c * inverse (real c) \<le> real (x * c div y) * inverse (real c)"
   2.970 +    using `0 < c` by auto
   2.971 +  thus ?thesis unfolding real_mult_assoc using `0 < c` by auto
   2.972 +qed
   2.973 +
   2.974 +lemma lapprox_posrat_bottom: assumes "0 < y"
   2.975 +  shows "real (x div y) \<le> Ifloat (lapprox_posrat n x y)" 
   2.976 +proof -
   2.977 +  have pow: "\<And>x. (0::int) < 2^x" by auto
   2.978 +  show ?thesis
   2.979 +    unfolding lapprox_posrat_def Let_def Ifloat_add normfloat Ifloat.simps pow2_minus pow2_int
   2.980 +    using real_of_int_div_mult[OF `0 < y` pow] by auto
   2.981 +qed
   2.982 +
   2.983 +lemma lapprox_posrat_nonneg: assumes "0 \<le> x" and "0 < y"
   2.984 +  shows "0 \<le> Ifloat (lapprox_posrat n x y)" 
   2.985 +proof -
   2.986 +  show ?thesis
   2.987 +    unfolding lapprox_posrat_def Let_def Ifloat_add normfloat Ifloat.simps pow2_minus pow2_int
   2.988 +    using pos_imp_zdiv_nonneg_iff[OF `0 < y`] assms by (auto intro!: mult_nonneg_nonneg)
   2.989 +qed
   2.990 +
   2.991 +definition rapprox_posrat :: "nat \<Rightarrow> int \<Rightarrow> int \<Rightarrow> float"
   2.992 +where
   2.993 +  "rapprox_posrat prec x y = (let
   2.994 +     l = nat (int prec + bitlen y - bitlen x) ;
   2.995 +     X = x * 2^l ;
   2.996 +     d = X div y ;
   2.997 +     m = X mod y
   2.998 +   in normfloat (Float (d + (if m = 0 then 0 else 1)) (- (int l))))"
   2.999  
  2.1000 -lemma one_eq_Numeral1_nring: "(1::'a::number_ring) = Numeral1"
  2.1001 -  by simp
  2.1002 +lemma rapprox_posrat:
  2.1003 +  assumes x: "0 \<le> x"
  2.1004 +  and y: "0 < y"
  2.1005 +  shows "real x / real y \<le> Ifloat (rapprox_posrat prec x y)"
  2.1006 +proof -
  2.1007 +  let ?l = "nat (int prec + bitlen y - bitlen x)" let ?X = "x * 2^?l"
  2.1008 +  show ?thesis 
  2.1009 +  proof (cases "?X mod y = 0")
  2.1010 +    case True hence "y \<noteq> 0" and "y dvd ?X" using `0 < y` by auto
  2.1011 +    from real_of_int_div[OF this]
  2.1012 +    have "real (?X div y) * inverse (2 ^ ?l) = real ?X / real y * inverse (2 ^ ?l)" by auto
  2.1013 +    also have "\<dots> = real x / real y * (2^?l * inverse (2^?l))" by auto
  2.1014 +    finally have "real (?X div y) * inverse (2^?l) = real x / real y" by auto
  2.1015 +    thus ?thesis unfolding rapprox_posrat_def Let_def normfloat if_P[OF True] 
  2.1016 +      unfolding Ifloat.simps pow2_minus pow2_int minus_minus by auto
  2.1017 +  next
  2.1018 +    case False
  2.1019 +    have "0 \<le> real y" and "real y \<noteq> 0" using `0 < y` by auto
  2.1020 +    have "0 \<le> real y * 2^?l" by (rule mult_nonneg_nonneg, rule `0 \<le> real y`, auto)
  2.1021  
  2.1022 -lemma zero_eq_Numeral0_nat: "(0::nat) = Numeral0"
  2.1023 -  by simp
  2.1024 +    have "?X = y * (?X div y) + ?X mod y" by auto
  2.1025 +    also have "\<dots> \<le> y * (?X div y) + y" by (rule add_mono, auto simp add: pos_mod_bound[OF `0 < y`, THEN less_imp_le])
  2.1026 +    also have "\<dots> = y * (?X div y + 1)" unfolding zadd_zmult_distrib2 by auto
  2.1027 +    finally have "real ?X \<le> real y * real (?X div y + 1)" unfolding real_of_int_le_iff real_of_int_mult[symmetric] .
  2.1028 +    hence "real ?X / (real y * 2^?l) \<le> real y * real (?X div y + 1) / (real y * 2^?l)" 
  2.1029 +      by (rule divide_right_mono, simp only: `0 \<le> real y * 2^?l`)
  2.1030 +    also have "\<dots> = real y * real (?X div y + 1) / real y / 2^?l" by auto
  2.1031 +    also have "\<dots> = real (?X div y + 1) * inverse (2^?l)" unfolding nonzero_mult_divide_cancel_left[OF `real y \<noteq> 0`] 
  2.1032 +      unfolding real_divide_def ..
  2.1033 +    finally show ?thesis unfolding rapprox_posrat_def Let_def normfloat Ifloat.simps if_not_P[OF False]
  2.1034 +      unfolding pow2_minus pow2_int minus_minus by auto
  2.1035 +  qed
  2.1036 +qed
  2.1037 +
  2.1038 +lemma rapprox_posrat_le1: assumes "0 \<le> x" and "0 < y" and "x \<le> y"
  2.1039 +  shows "Ifloat (rapprox_posrat n x y) \<le> 1"
  2.1040 +proof -
  2.1041 +  let ?l = "nat (int n + bitlen y - bitlen x)" let ?X = "x * 2^?l"
  2.1042 +  show ?thesis
  2.1043 +  proof (cases "?X mod y = 0")
  2.1044 +    case True hence "y \<noteq> 0" and "y dvd ?X" using `0 < y` by auto
  2.1045 +    from real_of_int_div[OF this]
  2.1046 +    have "real (?X div y) * inverse (2 ^ ?l) = real ?X / real y * inverse (2 ^ ?l)" by auto
  2.1047 +    also have "\<dots> = real x / real y * (2^?l * inverse (2^?l))" by auto
  2.1048 +    finally have "real (?X div y) * inverse (2^?l) = real x / real y" by auto
  2.1049 +    also have "real x / real y \<le> 1" using `0 \<le> x` and `0 < y` and `x \<le> y` by auto
  2.1050 +    finally show ?thesis unfolding rapprox_posrat_def Let_def normfloat if_P[OF True]
  2.1051 +      unfolding Ifloat.simps pow2_minus pow2_int minus_minus by auto
  2.1052 +  next
  2.1053 +    case False
  2.1054 +    have "x \<noteq> y"
  2.1055 +    proof (rule ccontr)
  2.1056 +      assume "\<not> x \<noteq> y" hence "x = y" by auto
  2.1057 +      have "?X mod y = 0" unfolding `x = y` using zmod_zmult_self2 by auto
  2.1058 +      thus False using False by auto
  2.1059 +    qed
  2.1060 +    hence "x < y" using `x \<le> y` by auto
  2.1061 +    hence "real x / real y < 1" using `0 < y` and `0 \<le> x` by auto
  2.1062  
  2.1063 -lemma one_eq_Numeral1_nat: "(1::nat) = Numeral1"
  2.1064 -  by simp
  2.1065 +    from real_of_int_div4[of "?X" y]
  2.1066 +    have "real (?X div y) \<le> (real x / real y) * 2^?l" unfolding real_of_int_mult times_divide_eq_left real_of_int_power[symmetric] real_number_of .
  2.1067 +    also have "\<dots> < 1 * 2^?l" using `real x / real y < 1` by (rule mult_strict_right_mono, auto)
  2.1068 +    finally have "?X div y < 2^?l" unfolding real_of_int_less_iff[of _ "2^?l", symmetric] by auto
  2.1069 +    hence "?X div y + 1 \<le> 2^?l" by auto
  2.1070 +    hence "real (?X div y + 1) * inverse (2^?l) \<le> 2^?l * inverse (2^?l)"
  2.1071 +      unfolding real_of_int_le_iff[of _ "2^?l", symmetric] real_of_int_power[symmetric] real_number_of
  2.1072 +      by (rule mult_right_mono, auto)
  2.1073 +    hence "real (?X div y + 1) * inverse (2^?l) \<le> 1" by auto
  2.1074 +    thus ?thesis unfolding rapprox_posrat_def Let_def normfloat Ifloat.simps if_not_P[OF False]
  2.1075 +      unfolding pow2_minus pow2_int minus_minus by auto
  2.1076 +  qed
  2.1077 +qed
  2.1078  
  2.1079 -lemma zpower_Pls: "(z::int)^Numeral0 = Numeral1"
  2.1080 -  by simp
  2.1081 +lemma zdiv_greater_zero: fixes a b :: int assumes "0 < a" and "a \<le> b"
  2.1082 +  shows "0 < b div a"
  2.1083 +proof (rule ccontr)
  2.1084 +  have "0 \<le> b" using assms by auto
  2.1085 +  assume "\<not> 0 < b div a" hence "b div a = 0" using `0 \<le> b`[unfolded pos_imp_zdiv_nonneg_iff[OF `0<a`, of b, symmetric]] by auto
  2.1086 +  have "b = a * (b div a) + b mod a" by auto
  2.1087 +  hence "b = b mod a" unfolding `b div a = 0` by auto
  2.1088 +  hence "b < a" using `0 < a`[THEN pos_mod_bound, of b] by auto
  2.1089 +  thus False using `a \<le> b` by auto
  2.1090 +qed
  2.1091 +
  2.1092 +lemma rapprox_posrat_less1: assumes "0 \<le> x" and "0 < y" and "2 * x < y" and "0 < n"
  2.1093 +  shows "Ifloat (rapprox_posrat n x y) < 1"
  2.1094 +proof (cases "x = 0")
  2.1095 +  case True thus ?thesis unfolding rapprox_posrat_def True Let_def normfloat Ifloat.simps by auto
  2.1096 +next
  2.1097 +  case False hence "0 < x" using `0 \<le> x` by auto
  2.1098 +  hence "x < y" using assms by auto
  2.1099 +  
  2.1100 +  let ?l = "nat (int n + bitlen y - bitlen x)" let ?X = "x * 2^?l"
  2.1101 +  show ?thesis
  2.1102 +  proof (cases "?X mod y = 0")
  2.1103 +    case True hence "y \<noteq> 0" and "y dvd ?X" using `0 < y` by auto
  2.1104 +    from real_of_int_div[OF this]
  2.1105 +    have "real (?X div y) * inverse (2 ^ ?l) = real ?X / real y * inverse (2 ^ ?l)" by auto
  2.1106 +    also have "\<dots> = real x / real y * (2^?l * inverse (2^?l))" by auto
  2.1107 +    finally have "real (?X div y) * inverse (2^?l) = real x / real y" by auto
  2.1108 +    also have "real x / real y < 1" using `0 \<le> x` and `0 < y` and `x < y` by auto
  2.1109 +    finally show ?thesis unfolding rapprox_posrat_def Let_def normfloat Ifloat.simps if_P[OF True]
  2.1110 +      unfolding pow2_minus pow2_int minus_minus by auto
  2.1111 +  next
  2.1112 +    case False
  2.1113 +    hence "(real x / real y) < 1 / 2" using `0 < y` and `0 \<le> x` `2 * x < y` by auto
  2.1114  
  2.1115 -lemma zpower_Min: "(z::int)^((-1)::nat) = Numeral1"
  2.1116 +    have "0 < ?X div y"
  2.1117 +    proof -
  2.1118 +      have "2^nat (bitlen x - 1) \<le> y" and "y < 2^nat (bitlen y)"
  2.1119 +	using bitlen_bounds[OF `0 < x`, THEN conjunct1] bitlen_bounds[OF `0 < y`, THEN conjunct2] `x < y` by auto
  2.1120 +      hence "(2::int)^nat (bitlen x - 1) < 2^nat (bitlen y)" by (rule order_le_less_trans)
  2.1121 +      hence "bitlen x \<le> bitlen y" by auto
  2.1122 +      hence len_less: "nat (bitlen x - 1) \<le> nat (int (n - 1) + bitlen y)" by auto
  2.1123 +
  2.1124 +      have "x \<noteq> 0" and "y \<noteq> 0" using `0 < x` `0 < y` by auto
  2.1125 +
  2.1126 +      have exp_eq: "nat (int (n - 1) + bitlen y) - nat (bitlen x - 1) = ?l"
  2.1127 +	using `bitlen x \<le> bitlen y` bitlen_ge1[OF `x \<noteq> 0`] bitlen_ge1[OF `y \<noteq> 0`] `0 < n` by auto
  2.1128 +
  2.1129 +      have "y * 2^nat (bitlen x - 1) \<le> y * x" 
  2.1130 +	using bitlen_bounds[OF `0 < x`, THEN conjunct1] `0 < y`[THEN less_imp_le] by (rule mult_left_mono)
  2.1131 +      also have "\<dots> \<le> 2^nat (bitlen y) * x" using bitlen_bounds[OF `0 < y`, THEN conjunct2, THEN less_imp_le] `0 \<le> x` by (rule mult_right_mono)
  2.1132 +      also have "\<dots> \<le> x * 2^nat (int (n - 1) + bitlen y)" unfolding mult_commute[of x] by (rule mult_right_mono, auto simp add: `0 \<le> x`)
  2.1133 +      finally have "real y * 2^nat (bitlen x - 1) * inverse (2^nat (bitlen x - 1)) \<le> real x * 2^nat (int (n - 1) + bitlen y) * inverse (2^nat (bitlen x - 1))"
  2.1134 +	unfolding real_of_int_le_iff[symmetric] by auto
  2.1135 +      hence "real y \<le> real x * (2^nat (int (n - 1) + bitlen y) / (2^nat (bitlen x - 1)))" 
  2.1136 +	unfolding real_mult_assoc real_divide_def by auto
  2.1137 +      also have "\<dots> = real x * (2^(nat (int (n - 1) + bitlen y) - nat (bitlen x - 1)))" using power_diff[of "2::real", OF _ len_less] by auto
  2.1138 +      finally have "y \<le> x * 2^?l" unfolding exp_eq unfolding real_of_int_le_iff[symmetric] by auto
  2.1139 +      thus ?thesis using zdiv_greater_zero[OF `0 < y`] by auto
  2.1140 +    qed
  2.1141 +
  2.1142 +    from real_of_int_div4[of "?X" y]
  2.1143 +    have "real (?X div y) \<le> (real x / real y) * 2^?l" unfolding real_of_int_mult times_divide_eq_left real_of_int_power[symmetric] real_number_of .
  2.1144 +    also have "\<dots> < 1/2 * 2^?l" using `real x / real y < 1/2` by (rule mult_strict_right_mono, auto)
  2.1145 +    finally have "?X div y * 2 < 2^?l" unfolding real_of_int_less_iff[of _ "2^?l", symmetric] by auto
  2.1146 +    hence "?X div y + 1 < 2^?l" using `0 < ?X div y` by auto
  2.1147 +    hence "real (?X div y + 1) * inverse (2^?l) < 2^?l * inverse (2^?l)"
  2.1148 +      unfolding real_of_int_less_iff[of _ "2^?l", symmetric] real_of_int_power[symmetric] real_number_of
  2.1149 +      by (rule mult_strict_right_mono, auto)
  2.1150 +    hence "real (?X div y + 1) * inverse (2^?l) < 1" by auto
  2.1151 +    thus ?thesis unfolding rapprox_posrat_def Let_def normfloat Ifloat.simps if_not_P[OF False]
  2.1152 +      unfolding pow2_minus pow2_int minus_minus by auto
  2.1153 +  qed
  2.1154 +qed
  2.1155 +
  2.1156 +lemma approx_rat_pattern: fixes P and ps :: "nat * int * int"
  2.1157 +  assumes Y: "\<And>y prec x. \<lbrakk>y = 0; ps = (prec, x, 0)\<rbrakk> \<Longrightarrow> P" 
  2.1158 +  and A: "\<And>x y prec. \<lbrakk>0 \<le> x; 0 < y; ps = (prec, x, y)\<rbrakk> \<Longrightarrow> P"
  2.1159 +  and B: "\<And>x y prec. \<lbrakk>x < 0; 0 < y; ps = (prec, x, y)\<rbrakk> \<Longrightarrow> P"
  2.1160 +  and C: "\<And>x y prec. \<lbrakk>x < 0; y < 0; ps = (prec, x, y)\<rbrakk> \<Longrightarrow> P"
  2.1161 +  and D: "\<And>x y prec. \<lbrakk>0 \<le> x; y < 0; ps = (prec, x, y)\<rbrakk> \<Longrightarrow> P"
  2.1162 +  shows P
  2.1163  proof -
  2.1164 -  have 1:"((-1)::nat) = 0"
  2.1165 -    by simp
  2.1166 -  show ?thesis by (simp add: 1)
  2.1167 +  obtain prec x y where [simp]: "ps = (prec, x, y)" by (cases ps, auto)
  2.1168 +  from Y have "y = 0 \<Longrightarrow> P" by auto
  2.1169 +  moreover { assume "0 < y" have P proof (cases "0 \<le> x") case True with A and `0 < y` show P by auto next case False with B and `0 < y` show P by auto qed } 
  2.1170 +  moreover { assume "y < 0" have P proof (cases "0 \<le> x") case True with D and `y < 0` show P by auto next case False with C and `y < 0` show P by auto qed }
  2.1171 +  ultimately show P by (cases "y = 0 \<or> 0 < y \<or> y < 0", auto)
  2.1172  qed
  2.1173  
  2.1174 -lemma fst_cong: "a=a' \<Longrightarrow> fst (a,b) = fst (a',b)"
  2.1175 -  by simp
  2.1176 +function lapprox_rat :: "nat \<Rightarrow> int \<Rightarrow> int \<Rightarrow> float"
  2.1177 +where
  2.1178 +  "y = 0 \<Longrightarrow> lapprox_rat prec x y = 0"
  2.1179 +| "0 \<le> x \<Longrightarrow> 0 < y \<Longrightarrow> lapprox_rat prec x y = lapprox_posrat prec x y"
  2.1180 +| "x < 0 \<Longrightarrow> 0 < y \<Longrightarrow> lapprox_rat prec x y = - (rapprox_posrat prec (-x) y)"
  2.1181 +| "x < 0 \<Longrightarrow> y < 0 \<Longrightarrow> lapprox_rat prec x y = lapprox_posrat prec (-x) (-y)"
  2.1182 +| "0 \<le> x \<Longrightarrow> y < 0 \<Longrightarrow> lapprox_rat prec x y = - (rapprox_posrat prec x (-y))"
  2.1183 +apply simp_all by (rule approx_rat_pattern)
  2.1184 +termination by lexicographic_order
  2.1185  
  2.1186 -lemma snd_cong: "b=b' \<Longrightarrow> snd (a,b) = snd (a,b')"
  2.1187 -  by simp
  2.1188 +lemma compute_lapprox_rat[code]:
  2.1189 +      "lapprox_rat prec x y = (if y = 0 then 0 else if 0 \<le> x then (if 0 < y then lapprox_posrat prec x y else - (rapprox_posrat prec x (-y))) 
  2.1190 +                                                             else (if 0 < y then - (rapprox_posrat prec (-x) y) else lapprox_posrat prec (-x) (-y)))"
  2.1191 +  by auto
  2.1192 +            
  2.1193 +lemma lapprox_rat: "Ifloat (lapprox_rat prec x y) \<le> real x / real y"
  2.1194 +proof -      
  2.1195 +  have h[rule_format]: "! a b b'. b' \<le> b \<longrightarrow> a \<le> b' \<longrightarrow> a \<le> (b::real)" by auto
  2.1196 +  show ?thesis
  2.1197 +    apply (case_tac "y = 0")
  2.1198 +    apply simp
  2.1199 +    apply (case_tac "0 \<le> x \<and> 0 < y")
  2.1200 +    apply (simp add: lapprox_posrat)
  2.1201 +    apply (case_tac "x < 0 \<and> 0 < y")
  2.1202 +    apply simp
  2.1203 +    apply (subst minus_le_iff)   
  2.1204 +    apply (rule h[OF rapprox_posrat])
  2.1205 +    apply (simp_all)
  2.1206 +    apply (case_tac "x < 0 \<and> y < 0")
  2.1207 +    apply simp
  2.1208 +    apply (rule h[OF _ lapprox_posrat])
  2.1209 +    apply (simp_all)
  2.1210 +    apply (case_tac "0 \<le> x \<and> y < 0")
  2.1211 +    apply (simp)
  2.1212 +    apply (subst minus_le_iff)   
  2.1213 +    apply (rule h[OF rapprox_posrat])
  2.1214 +    apply simp_all
  2.1215 +    apply arith
  2.1216 +    done
  2.1217 +qed
  2.1218  
  2.1219 -lemma lift_bool: "x \<Longrightarrow> x=True"
  2.1220 -  by simp
  2.1221 +lemma lapprox_rat_bottom: assumes "0 \<le> x" and "0 < y"
  2.1222 +  shows "real (x div y) \<le> Ifloat (lapprox_rat n x y)" 
  2.1223 +  unfolding lapprox_rat.simps(2)[OF assms]  using lapprox_posrat_bottom[OF `0<y`] .
  2.1224 +
  2.1225 +function rapprox_rat :: "nat \<Rightarrow> int \<Rightarrow> int \<Rightarrow> float"
  2.1226 +where
  2.1227 +  "y = 0 \<Longrightarrow> rapprox_rat prec x y = 0"
  2.1228 +| "0 \<le> x \<Longrightarrow> 0 < y \<Longrightarrow> rapprox_rat prec x y = rapprox_posrat prec x y"
  2.1229 +| "x < 0 \<Longrightarrow> 0 < y \<Longrightarrow> rapprox_rat prec x y = - (lapprox_posrat prec (-x) y)"
  2.1230 +| "x < 0 \<Longrightarrow> y < 0 \<Longrightarrow> rapprox_rat prec x y = rapprox_posrat prec (-x) (-y)"
  2.1231 +| "0 \<le> x \<Longrightarrow> y < 0 \<Longrightarrow> rapprox_rat prec x y = - (lapprox_posrat prec x (-y))"
  2.1232 +apply simp_all by (rule approx_rat_pattern)
  2.1233 +termination by lexicographic_order
  2.1234  
  2.1235 -lemma nlift_bool: "~x \<Longrightarrow> x=False"
  2.1236 -  by simp
  2.1237 +lemma compute_rapprox_rat[code]:
  2.1238 +      "rapprox_rat prec x y = (if y = 0 then 0 else if 0 \<le> x then (if 0 < y then rapprox_posrat prec x y else - (lapprox_posrat prec x (-y))) else 
  2.1239 +                                                                  (if 0 < y then - (lapprox_posrat prec (-x) y) else rapprox_posrat prec (-x) (-y)))"
  2.1240 +  by auto
  2.1241  
  2.1242 -lemma not_false_eq_true: "(~ False) = True" by simp
  2.1243 +lemma rapprox_rat: "real x / real y \<le> Ifloat (rapprox_rat prec x y)"
  2.1244 +proof -      
  2.1245 +  have h[rule_format]: "! a b b'. b' \<le> b \<longrightarrow> a \<le> b' \<longrightarrow> a \<le> (b::real)" by auto
  2.1246 +  show ?thesis
  2.1247 +    apply (case_tac "y = 0")
  2.1248 +    apply simp
  2.1249 +    apply (case_tac "0 \<le> x \<and> 0 < y")
  2.1250 +    apply (simp add: rapprox_posrat)
  2.1251 +    apply (case_tac "x < 0 \<and> 0 < y")
  2.1252 +    apply simp
  2.1253 +    apply (subst le_minus_iff)   
  2.1254 +    apply (rule h[OF _ lapprox_posrat])
  2.1255 +    apply (simp_all)
  2.1256 +    apply (case_tac "x < 0 \<and> y < 0")
  2.1257 +    apply simp
  2.1258 +    apply (rule h[OF rapprox_posrat])
  2.1259 +    apply (simp_all)
  2.1260 +    apply (case_tac "0 \<le> x \<and> y < 0")
  2.1261 +    apply (simp)
  2.1262 +    apply (subst le_minus_iff)   
  2.1263 +    apply (rule h[OF _ lapprox_posrat])
  2.1264 +    apply simp_all
  2.1265 +    apply arith
  2.1266 +    done
  2.1267 +qed
  2.1268  
  2.1269 -lemma not_true_eq_false: "(~ True) = False" by simp
  2.1270 +lemma rapprox_rat_le1: assumes "0 \<le> x" and "0 < y" and "x \<le> y"
  2.1271 +  shows "Ifloat (rapprox_rat n x y) \<le> 1"
  2.1272 +  unfolding rapprox_rat.simps(2)[OF `0 \<le> x` `0 < y`] using rapprox_posrat_le1[OF assms] .
  2.1273 +
  2.1274 +lemma rapprox_rat_neg: assumes "x < 0" and "0 < y"
  2.1275 +  shows "Ifloat (rapprox_rat n x y) \<le> 0"
  2.1276 +  unfolding rapprox_rat.simps(3)[OF assms] using lapprox_posrat_nonneg[of "-x" y n] assms by auto
  2.1277 +
  2.1278 +lemma rapprox_rat_nonneg_neg: assumes "0 \<le> x" and "y < 0"
  2.1279 +  shows "Ifloat (rapprox_rat n x y) \<le> 0"
  2.1280 +  unfolding rapprox_rat.simps(5)[OF assms] using lapprox_posrat_nonneg[of x "-y" n] assms by auto
  2.1281  
  2.1282 -lemmas binarith =
  2.1283 -  normalize_bin_simps
  2.1284 -  pred_bin_simps succ_bin_simps
  2.1285 -  add_bin_simps minus_bin_simps mult_bin_simps
  2.1286 +lemma rapprox_rat_nonpos_pos: assumes "x \<le> 0" and "0 < y"
  2.1287 +  shows "Ifloat (rapprox_rat n x y) \<le> 0"
  2.1288 +proof (cases "x = 0") 
  2.1289 +  case True hence "0 \<le> x" by auto show ?thesis unfolding rapprox_rat.simps(2)[OF `0 \<le> x` `0 < y`]
  2.1290 +    unfolding True rapprox_posrat_def Let_def by auto
  2.1291 +next
  2.1292 +  case False hence "x < 0" using assms by auto
  2.1293 +  show ?thesis using rapprox_rat_neg[OF `x < 0` `0 < y`] .
  2.1294 +qed
  2.1295 +
  2.1296 +fun float_divl :: "nat \<Rightarrow> float \<Rightarrow> float \<Rightarrow> float"
  2.1297 +where
  2.1298 +  "float_divl prec (Float m1 s1) (Float m2 s2) = 
  2.1299 +    (let
  2.1300 +       l = lapprox_rat prec m1 m2;
  2.1301 +       f = Float 1 (s1 - s2)
  2.1302 +     in
  2.1303 +       f * l)"     
  2.1304  
  2.1305 -lemma int_eq_number_of_eq:
  2.1306 -  "(((number_of v)::int)=(number_of w)) = iszero ((number_of (v + uminus w))::int)"
  2.1307 -  by (rule eq_number_of_eq)
  2.1308 +lemma float_divl: "Ifloat (float_divl prec x y) \<le> Ifloat x / Ifloat y"
  2.1309 +proof - 
  2.1310 +  from float_split[of x] obtain mx sx where x: "x = Float mx sx" by auto
  2.1311 +  from float_split[of y] obtain my sy where y: "y = Float my sy" by auto
  2.1312 +  have "real mx / real my \<le> (real mx * pow2 sx / (real my * pow2 sy)) / (pow2 (sx - sy))"
  2.1313 +    apply (case_tac "my = 0")
  2.1314 +    apply simp
  2.1315 +    apply (case_tac "my > 0")       
  2.1316 +    apply (subst pos_le_divide_eq)
  2.1317 +    apply simp
  2.1318 +    apply (subst pos_le_divide_eq)
  2.1319 +    apply (simp add: mult_pos_pos)
  2.1320 +    apply simp
  2.1321 +    apply (subst pow2_add[symmetric])
  2.1322 +    apply simp
  2.1323 +    apply (subgoal_tac "my < 0")
  2.1324 +    apply auto
  2.1325 +    apply (simp add: field_simps)
  2.1326 +    apply (subst pow2_add[symmetric])
  2.1327 +    apply (simp add: field_simps)
  2.1328 +    done
  2.1329 +  then have "Ifloat (lapprox_rat prec mx my) \<le> (real mx * pow2 sx / (real my * pow2 sy)) / (pow2 (sx - sy))"
  2.1330 +    by (rule order_trans[OF lapprox_rat])
  2.1331 +  then have "Ifloat (lapprox_rat prec mx my) * pow2 (sx - sy) \<le> real mx * pow2 sx / (real my * pow2 sy)"
  2.1332 +    apply (subst pos_le_divide_eq[symmetric])
  2.1333 +    apply simp_all
  2.1334 +    done
  2.1335 +  then have "pow2 (sx - sy) * Ifloat (lapprox_rat prec mx my) \<le> real mx * pow2 sx / (real my * pow2 sy)"
  2.1336 +    by (simp add: algebra_simps)
  2.1337 +  then show ?thesis
  2.1338 +    by (simp add: x y Let_def Ifloat.simps)
  2.1339 +qed
  2.1340  
  2.1341 -lemma int_iszero_number_of_Pls: "iszero (Numeral0::int)"
  2.1342 -  by (simp only: iszero_number_of_Pls)
  2.1343 +lemma float_divl_lower_bound: assumes "0 \<le> x" and "0 < y" shows "0 \<le> float_divl prec x y"
  2.1344 +proof (cases x, cases y)
  2.1345 +  fix xm xe ym ye :: int
  2.1346 +  assume x_eq: "x = Float xm xe" and y_eq: "y = Float ym ye"
  2.1347 +  have "0 \<le> xm" using `0 \<le> x`[unfolded x_eq le_float_def Ifloat.simps Ifloat_0 zero_le_mult_iff] by auto
  2.1348 +  have "0 < ym" using `0 < y`[unfolded y_eq less_float_def Ifloat.simps Ifloat_0 zero_less_mult_iff] by auto
  2.1349  
  2.1350 -lemma int_nonzero_number_of_Min: "~(iszero ((-1)::int))"
  2.1351 -  by simp
  2.1352 +  have "\<And>n. 0 \<le> Ifloat (Float 1 n)" unfolding Ifloat.simps using zero_le_pow2 by auto
  2.1353 +  moreover have "0 \<le> Ifloat (lapprox_rat prec xm ym)" by (rule order_trans[OF _ lapprox_rat_bottom[OF `0 \<le> xm` `0 < ym`]], auto simp add: `0 \<le> xm` pos_imp_zdiv_nonneg_iff[OF `0 < ym`])
  2.1354 +  ultimately show "0 \<le> float_divl prec x y"
  2.1355 +    unfolding x_eq y_eq float_divl.simps Let_def le_float_def Ifloat_0 by (auto intro!: mult_nonneg_nonneg)
  2.1356 +qed
  2.1357 +
  2.1358 +lemma float_divl_pos_less1_bound: assumes "0 < x" and "x < 1" and "0 < prec" shows "1 \<le> float_divl prec 1 x"
  2.1359 +proof (cases x)
  2.1360 +  case (Float m e)
  2.1361 +  from `0 < x` `x < 1` have "0 < m" "e < 0" using float_pos_m_pos float_pos_less1_e_neg unfolding Float by auto
  2.1362 +  let ?b = "nat (bitlen m)" and ?e = "nat (-e)"
  2.1363 +  have "1 \<le> m" and "m \<noteq> 0" using `0 < m` by auto
  2.1364 +  with bitlen_bounds[OF `0 < m`] have "m < 2^?b" and "(2::int) \<le> 2^?b" by auto
  2.1365 +  hence "1 \<le> bitlen m" using power_le_imp_le_exp[of "2::int" 1 ?b] by auto
  2.1366 +  hence pow_split: "nat (int prec + bitlen m - 1) = (prec - 1) + ?b" using `0 < prec` by auto
  2.1367 +  
  2.1368 +  have pow_not0: "\<And>x. (2::real)^x \<noteq> 0" by auto
  2.1369  
  2.1370 -lemma int_iszero_number_of_Bit0: "iszero ((number_of (Int.Bit0 w))::int) = iszero ((number_of w)::int)"
  2.1371 -  by simp
  2.1372 +  from float_less1_mantissa_bound `0 < x` `x < 1` Float 
  2.1373 +  have "m < 2^?e" by auto
  2.1374 +  with bitlen_bounds[OF `0 < m`, THEN conjunct1]
  2.1375 +  have "(2::int)^nat (bitlen m - 1) < 2^?e" by (rule order_le_less_trans)
  2.1376 +  from power_less_imp_less_exp[OF _ this]
  2.1377 +  have "bitlen m \<le> - e" by auto
  2.1378 +  hence "(2::real)^?b \<le> 2^?e" by auto
  2.1379 +  hence "(2::real)^?b * inverse (2^?b) \<le> 2^?e * inverse (2^?b)" by (rule mult_right_mono, auto)
  2.1380 +  hence "(1::real) \<le> 2^?e * inverse (2^?b)" by auto
  2.1381 +  also
  2.1382 +  let ?d = "real (2 ^ nat (int prec + bitlen m - 1) div m) * inverse (2 ^ nat (int prec + bitlen m - 1))"
  2.1383 +  { have "2^(prec - 1) * m \<le> 2^(prec - 1) * 2^?b" using `m < 2^?b`[THEN less_imp_le] by (rule mult_left_mono, auto)
  2.1384 +    also have "\<dots> = 2 ^ nat (int prec + bitlen m - 1)" unfolding pow_split zpower_zadd_distrib by auto
  2.1385 +    finally have "2^(prec - 1) * m div m \<le> 2 ^ nat (int prec + bitlen m - 1) div m" using `0 < m` by (rule zdiv_mono1)
  2.1386 +    hence "2^(prec - 1) \<le> 2 ^ nat (int prec + bitlen m - 1) div m" unfolding zdiv_zmult_self1[OF `m \<noteq> 0`] .
  2.1387 +    hence "2^(prec - 1) * inverse (2 ^ nat (int prec + bitlen m - 1)) \<le> ?d"
  2.1388 +      unfolding real_of_int_le_iff[of "2^(prec - 1)", symmetric] by auto }
  2.1389 +  from mult_left_mono[OF this[unfolded pow_split power_add inverse_mult_distrib real_mult_assoc[symmetric] right_inverse[OF pow_not0] real_mult_1], of "2^?e"]
  2.1390 +  have "2^?e * inverse (2^?b) \<le> 2^?e * ?d" unfolding pow_split power_add by auto
  2.1391 +  finally have "1 \<le> 2^?e * ?d" .
  2.1392 +  
  2.1393 +  have e_nat: "0 - e = int (nat (-e))" using `e < 0` by auto
  2.1394 +  have "bitlen 1 = 1" using bitlen.simps by auto
  2.1395 +  
  2.1396 +  show ?thesis 
  2.1397 +    unfolding one_float_def Float float_divl.simps Let_def lapprox_rat.simps(2)[OF zero_le_one `0 < m`] lapprox_posrat_def `bitlen 1 = 1`
  2.1398 +    unfolding le_float_def Ifloat_mult normfloat Ifloat.simps pow2_minus pow2_int e_nat
  2.1399 +    using `1 \<le> 2^?e * ?d` by (auto simp add: pow2_def)
  2.1400 +qed
  2.1401  
  2.1402 -lemma int_iszero_number_of_Bit1: "\<not> iszero ((number_of (Int.Bit1 w))::int)"
  2.1403 -  by simp
  2.1404 +fun float_divr :: "nat \<Rightarrow> float \<Rightarrow> float \<Rightarrow> float"
  2.1405 +where
  2.1406 +  "float_divr prec (Float m1 s1) (Float m2 s2) = 
  2.1407 +    (let
  2.1408 +       r = rapprox_rat prec m1 m2;
  2.1409 +       f = Float 1 (s1 - s2)
  2.1410 +     in
  2.1411 +       f * r)"  
  2.1412  
  2.1413 -lemma int_less_number_of_eq_neg: "(((number_of x)::int) < number_of y) = neg ((number_of (x + (uminus y)))::int)"
  2.1414 -  unfolding neg_def number_of_is_id by simp
  2.1415 +lemma float_divr: "Ifloat x / Ifloat y \<le> Ifloat (float_divr prec x y)"
  2.1416 +proof - 
  2.1417 +  from float_split[of x] obtain mx sx where x: "x = Float mx sx" by auto
  2.1418 +  from float_split[of y] obtain my sy where y: "y = Float my sy" by auto
  2.1419 +  have "real mx / real my \<ge> (real mx * pow2 sx / (real my * pow2 sy)) / (pow2 (sx - sy))"
  2.1420 +    apply (case_tac "my = 0")
  2.1421 +    apply simp
  2.1422 +    apply (case_tac "my > 0")
  2.1423 +    apply auto
  2.1424 +    apply (subst pos_divide_le_eq)
  2.1425 +    apply (rule mult_pos_pos)+
  2.1426 +    apply simp_all
  2.1427 +    apply (subst pow2_add[symmetric])
  2.1428 +    apply simp
  2.1429 +    apply (subgoal_tac "my < 0")
  2.1430 +    apply auto
  2.1431 +    apply (simp add: field_simps)
  2.1432 +    apply (subst pow2_add[symmetric])
  2.1433 +    apply (simp add: field_simps)
  2.1434 +    done
  2.1435 +  then have "Ifloat (rapprox_rat prec mx my) \<ge> (real mx * pow2 sx / (real my * pow2 sy)) / (pow2 (sx - sy))"
  2.1436 +    by (rule order_trans[OF _ rapprox_rat])
  2.1437 +  then have "Ifloat (rapprox_rat prec mx my) * pow2 (sx - sy) \<ge> real mx * pow2 sx / (real my * pow2 sy)"
  2.1438 +    apply (subst pos_divide_le_eq[symmetric])
  2.1439 +    apply simp_all
  2.1440 +    done
  2.1441 +  then show ?thesis
  2.1442 +    by (simp add: x y Let_def algebra_simps Ifloat.simps)
  2.1443 +qed
  2.1444  
  2.1445 -lemma int_not_neg_number_of_Pls: "\<not> (neg (Numeral0::int))"
  2.1446 -  by simp
  2.1447 +lemma float_divr_pos_less1_lower_bound: assumes "0 < x" and "x < 1" shows "1 \<le> float_divr prec 1 x"
  2.1448 +proof -
  2.1449 +  have "1 \<le> 1 / Ifloat x" using `0 < x` and `x < 1` unfolding less_float_def by auto
  2.1450 +  also have "\<dots> \<le> Ifloat (float_divr prec 1 x)" using float_divr[where x=1 and y=x] by auto
  2.1451 +  finally show ?thesis unfolding le_float_def by auto
  2.1452 +qed
  2.1453 +
  2.1454 +lemma float_divr_nonpos_pos_upper_bound: assumes "x \<le> 0" and "0 < y" shows "float_divr prec x y \<le> 0"
  2.1455 +proof (cases x, cases y)
  2.1456 +  fix xm xe ym ye :: int
  2.1457 +  assume x_eq: "x = Float xm xe" and y_eq: "y = Float ym ye"
  2.1458 +  have "xm \<le> 0" using `x \<le> 0`[unfolded x_eq le_float_def Ifloat.simps Ifloat_0 mult_le_0_iff] by auto
  2.1459 +  have "0 < ym" using `0 < y`[unfolded y_eq less_float_def Ifloat.simps Ifloat_0 zero_less_mult_iff] by auto
  2.1460 +
  2.1461 +  have "\<And>n. 0 \<le> Ifloat (Float 1 n)" unfolding Ifloat.simps using zero_le_pow2 by auto
  2.1462 +  moreover have "Ifloat (rapprox_rat prec xm ym) \<le> 0" using rapprox_rat_nonpos_pos[OF `xm \<le> 0` `0 < ym`] .
  2.1463 +  ultimately show "float_divr prec x y \<le> 0"
  2.1464 +    unfolding x_eq y_eq float_divr.simps Let_def le_float_def Ifloat_0 Ifloat_mult by (auto intro!: mult_nonneg_nonpos)
  2.1465 +qed
  2.1466  
  2.1467 -lemma int_neg_number_of_Min: "neg (-1::int)"
  2.1468 -  by simp
  2.1469 +lemma float_divr_nonneg_neg_upper_bound: assumes "0 \<le> x" and "y < 0" shows "float_divr prec x y \<le> 0"
  2.1470 +proof (cases x, cases y)
  2.1471 +  fix xm xe ym ye :: int
  2.1472 +  assume x_eq: "x = Float xm xe" and y_eq: "y = Float ym ye"
  2.1473 +  have "0 \<le> xm" using `0 \<le> x`[unfolded x_eq le_float_def Ifloat.simps Ifloat_0 zero_le_mult_iff] by auto
  2.1474 +  have "ym < 0" using `y < 0`[unfolded y_eq less_float_def Ifloat.simps Ifloat_0 mult_less_0_iff] by auto
  2.1475 +  hence "0 < - ym" by auto
  2.1476 +
  2.1477 +  have "\<And>n. 0 \<le> Ifloat (Float 1 n)" unfolding Ifloat.simps using zero_le_pow2 by auto
  2.1478 +  moreover have "Ifloat (rapprox_rat prec xm ym) \<le> 0" using rapprox_rat_nonneg_neg[OF `0 \<le> xm` `ym < 0`] .
  2.1479 +  ultimately show "float_divr prec x y \<le> 0"
  2.1480 +    unfolding x_eq y_eq float_divr.simps Let_def le_float_def Ifloat_0 Ifloat_mult by (auto intro!: mult_nonneg_nonpos)
  2.1481 +qed
  2.1482 +
  2.1483 +fun round_down :: "nat \<Rightarrow> float \<Rightarrow> float" where
  2.1484 +"round_down prec (Float m e) = (let d = bitlen m - int prec in
  2.1485 +     if 0 < d then let P = 2^nat d ; n = m div P in Float n (e + d)
  2.1486 +              else Float m e)"
  2.1487 +
  2.1488 +fun round_up :: "nat \<Rightarrow> float \<Rightarrow> float" where
  2.1489 +"round_up prec (Float m e) = (let d = bitlen m - int prec in
  2.1490 +  if 0 < d then let P = 2^nat d ; n = m div P ; r = m mod P in Float (n + (if r = 0 then 0 else 1)) (e + d) 
  2.1491 +           else Float m e)"
  2.1492  
  2.1493 -lemma int_neg_number_of_Bit0: "neg ((number_of (Int.Bit0 w))::int) = neg ((number_of w)::int)"
  2.1494 -  by simp
  2.1495 -
  2.1496 -lemma int_neg_number_of_Bit1: "neg ((number_of (Int.Bit1 w))::int) = neg ((number_of w)::int)"
  2.1497 -  by simp
  2.1498 +lemma round_up: "Ifloat x \<le> Ifloat (round_up prec x)"
  2.1499 +proof (cases x)
  2.1500 +  case (Float m e)
  2.1501 +  let ?d = "bitlen m - int prec"
  2.1502 +  let ?p = "(2::int)^nat ?d"
  2.1503 +  have "0 < ?p" by auto
  2.1504 +  show "?thesis"
  2.1505 +  proof (cases "0 < ?d")
  2.1506 +    case True
  2.1507 +    hence pow_d: "pow2 ?d = real ?p" unfolding pow2_int[symmetric] power_real_number_of[symmetric] by auto
  2.1508 +    show ?thesis
  2.1509 +    proof (cases "m mod ?p = 0")
  2.1510 +      case True
  2.1511 +      have m: "m = m div ?p * ?p + 0" unfolding True[symmetric] using zdiv_zmod_equality2[where k=0, unfolded monoid_add_class.add_0_right, symmetric] .
  2.1512 +      have "Ifloat (Float m e) = Ifloat (Float (m div ?p) (e + ?d))" unfolding Ifloat.simps arg_cong[OF m, of real]
  2.1513 +	by (auto simp add: pow2_add `0 < ?d` pow_d)
  2.1514 +      thus ?thesis
  2.1515 +	unfolding Float round_up.simps Let_def if_P[OF `m mod ?p = 0`] if_P[OF `0 < ?d`]
  2.1516 +	by auto
  2.1517 +    next
  2.1518 +      case False
  2.1519 +      have "m = m div ?p * ?p + m mod ?p" unfolding zdiv_zmod_equality2[where k=0, unfolded monoid_add_class.add_0_right] ..
  2.1520 +      also have "\<dots> \<le> (m div ?p + 1) * ?p" unfolding left_distrib zmult_1 by (rule add_left_mono, rule pos_mod_bound[OF `0 < ?p`, THEN less_imp_le])
  2.1521 +      finally have "Ifloat (Float m e) \<le> Ifloat (Float (m div ?p + 1) (e + ?d))" unfolding Ifloat.simps add_commute[of e]
  2.1522 +	unfolding pow2_add mult_assoc[symmetric] real_of_int_le_iff[of m, symmetric]
  2.1523 +	by (auto intro!: mult_mono simp add: pow2_add `0 < ?d` pow_d)
  2.1524 +      thus ?thesis
  2.1525 +	unfolding Float round_up.simps Let_def if_not_P[OF `\<not> m mod ?p = 0`] if_P[OF `0 < ?d`] .
  2.1526 +    qed
  2.1527 +  next
  2.1528 +    case False
  2.1529 +    show ?thesis
  2.1530 +      unfolding Float round_up.simps Let_def if_not_P[OF False] .. 
  2.1531 +  qed
  2.1532 +qed
  2.1533  
  2.1534 -lemma int_le_number_of_eq: "(((number_of x)::int) \<le> number_of y) = (\<not> neg ((number_of (y + (uminus x)))::int))"
  2.1535 -  unfolding neg_def number_of_is_id by (simp add: not_less)
  2.1536 +lemma round_down: "Ifloat (round_down prec x) \<le> Ifloat x"
  2.1537 +proof (cases x)
  2.1538 +  case (Float m e)
  2.1539 +  let ?d = "bitlen m - int prec"
  2.1540 +  let ?p = "(2::int)^nat ?d"
  2.1541 +  have "0 < ?p" by auto
  2.1542 +  show "?thesis"
  2.1543 +  proof (cases "0 < ?d")
  2.1544 +    case True
  2.1545 +    hence pow_d: "pow2 ?d = real ?p" unfolding pow2_int[symmetric] power_real_number_of[symmetric] by auto
  2.1546 +    have "m div ?p * ?p \<le> m div ?p * ?p + m mod ?p" by (auto simp add: pos_mod_bound[OF `0 < ?p`, THEN less_imp_le])
  2.1547 +    also have "\<dots> \<le> m" unfolding zdiv_zmod_equality2[where k=0, unfolded monoid_add_class.add_0_right] ..
  2.1548 +    finally have "Ifloat (Float (m div ?p) (e + ?d)) \<le> Ifloat (Float m e)" unfolding Ifloat.simps add_commute[of e]
  2.1549 +      unfolding pow2_add mult_assoc[symmetric] real_of_int_le_iff[of _ m, symmetric]
  2.1550 +      by (auto intro!: mult_mono simp add: pow2_add `0 < ?d` pow_d)
  2.1551 +    thus ?thesis
  2.1552 +      unfolding Float round_down.simps Let_def if_P[OF `0 < ?d`] .
  2.1553 +  next
  2.1554 +    case False
  2.1555 +    show ?thesis
  2.1556 +      unfolding Float round_down.simps Let_def if_not_P[OF False] .. 
  2.1557 +  qed
  2.1558 +qed
  2.1559 +
  2.1560 +definition lb_mult :: "nat \<Rightarrow> float \<Rightarrow> float \<Rightarrow> float" where
  2.1561 +"lb_mult prec x y = (case normfloat (x * y) of Float m e \<Rightarrow> let
  2.1562 +    l = bitlen m - int prec
  2.1563 +  in if l > 0 then Float (m div (2^nat l)) (e + l)
  2.1564 +              else Float m e)"
  2.1565  
  2.1566 -lemmas intarithrel =
  2.1567 -  int_eq_number_of_eq
  2.1568 -  lift_bool[OF int_iszero_number_of_Pls] nlift_bool[OF int_nonzero_number_of_Min] int_iszero_number_of_Bit0
  2.1569 -  lift_bool[OF int_iszero_number_of_Bit1] int_less_number_of_eq_neg nlift_bool[OF int_not_neg_number_of_Pls] lift_bool[OF int_neg_number_of_Min]
  2.1570 -  int_neg_number_of_Bit0 int_neg_number_of_Bit1 int_le_number_of_eq
  2.1571 +definition ub_mult :: "nat \<Rightarrow> float \<Rightarrow> float \<Rightarrow> float" where
  2.1572 +"ub_mult prec x y = (case normfloat (x * y) of Float m e \<Rightarrow> let
  2.1573 +    l = bitlen m - int prec
  2.1574 +  in if l > 0 then Float (m div (2^nat l) + 1) (e + l)
  2.1575 +              else Float m e)"
  2.1576  
  2.1577 -lemma int_number_of_add_sym: "((number_of v)::int) + number_of w = number_of (v + w)"
  2.1578 -  by simp
  2.1579 +lemma lb_mult: "Ifloat (lb_mult prec x y) \<le> Ifloat (x * y)"
  2.1580 +proof (cases "normfloat (x * y)")
  2.1581 +  case (Float m e)
  2.1582 +  hence "odd m \<or> (m = 0 \<and> e = 0)" by (rule normfloat_imp_odd_or_zero)
  2.1583 +  let ?l = "bitlen m - int prec"
  2.1584 +  have "Ifloat (lb_mult prec x y) \<le> Ifloat (normfloat (x * y))"
  2.1585 +  proof (cases "?l > 0")
  2.1586 +    case False thus ?thesis unfolding lb_mult_def Float Let_def float.cases by auto
  2.1587 +  next
  2.1588 +    case True
  2.1589 +    have "real (m div 2^(nat ?l)) * pow2 ?l \<le> real m"
  2.1590 +    proof -
  2.1591 +      have "real (m div 2^(nat ?l)) * pow2 ?l = real (2^(nat ?l) * (m div 2^(nat ?l)))" unfolding real_of_int_mult real_of_int_power[symmetric] real_number_of unfolding pow2_int[symmetric] 
  2.1592 +	using `?l > 0` by auto
  2.1593 +      also have "\<dots> \<le> real (2^(nat ?l) * (m div 2^(nat ?l)) + m mod 2^(nat ?l))" unfolding real_of_int_add by auto
  2.1594 +      also have "\<dots> = real m" unfolding zmod_zdiv_equality[symmetric] ..
  2.1595 +      finally show ?thesis by auto
  2.1596 +    qed
  2.1597 +    thus ?thesis unfolding lb_mult_def Float Let_def float.cases if_P[OF True] Ifloat.simps pow2_add real_mult_commute real_mult_assoc by auto
  2.1598 +  qed
  2.1599 +  also have "\<dots> = Ifloat (x * y)" unfolding normfloat ..
  2.1600 +  finally show ?thesis .
  2.1601 +qed
  2.1602  
  2.1603 -lemma int_number_of_diff_sym: "((number_of v)::int) - number_of w = number_of (v + (uminus w))"
  2.1604 -  by simp
  2.1605 +lemma ub_mult: "Ifloat (x * y) \<le> Ifloat (ub_mult prec x y)"
  2.1606 +proof (cases "normfloat (x * y)")
  2.1607 +  case (Float m e)
  2.1608 +  hence "odd m \<or> (m = 0 \<and> e = 0)" by (rule normfloat_imp_odd_or_zero)
  2.1609 +  let ?l = "bitlen m - int prec"
  2.1610 +  have "Ifloat (x * y) = Ifloat (normfloat (x * y))" unfolding normfloat ..
  2.1611 +  also have "\<dots> \<le> Ifloat (ub_mult prec x y)"
  2.1612 +  proof (cases "?l > 0")
  2.1613 +    case False thus ?thesis unfolding ub_mult_def Float Let_def float.cases by auto
  2.1614 +  next
  2.1615 +    case True
  2.1616 +    have "real m \<le> real (m div 2^(nat ?l) + 1) * pow2 ?l"
  2.1617 +    proof -
  2.1618 +      have "m mod 2^(nat ?l) < 2^(nat ?l)" by (rule pos_mod_bound) auto
  2.1619 +      hence mod_uneq: "real (m mod 2^(nat ?l)) \<le> 1 * 2^(nat ?l)" unfolding zmult_1 real_of_int_less_iff[symmetric] by auto
  2.1620 +      
  2.1621 +      have "real m = real (2^(nat ?l) * (m div 2^(nat ?l)) + m mod 2^(nat ?l))" unfolding zmod_zdiv_equality[symmetric] ..
  2.1622 +      also have "\<dots> = real (m div 2^(nat ?l)) * 2^(nat ?l) + real (m mod 2^(nat ?l))" unfolding real_of_int_add by auto
  2.1623 +      also have "\<dots> \<le> (real (m div 2^(nat ?l)) + 1) * 2^(nat ?l)" unfolding real_add_mult_distrib using mod_uneq by auto
  2.1624 +      finally show ?thesis unfolding pow2_int[symmetric] using True by auto
  2.1625 +    qed
  2.1626 +    thus ?thesis unfolding ub_mult_def Float Let_def float.cases if_P[OF True] Ifloat.simps pow2_add real_mult_commute real_mult_assoc by auto
  2.1627 +  qed
  2.1628 +  finally show ?thesis .
  2.1629 +qed
  2.1630 +
  2.1631 +fun float_abs :: "float \<Rightarrow> float" where
  2.1632 +"float_abs (Float m e) = Float \<bar>m\<bar> e"
  2.1633 +
  2.1634 +instantiation float :: abs begin
  2.1635 +definition abs_float_def: "\<bar>x\<bar> = float_abs x"
  2.1636 +instance ..
  2.1637 +end
  2.1638  
  2.1639 -lemma int_number_of_mult_sym: "((number_of v)::int) * number_of w = number_of (v * w)"
  2.1640 -  by simp
  2.1641 +lemma Ifloat_abs: "Ifloat \<bar>x\<bar> = \<bar>Ifloat x\<bar>" 
  2.1642 +proof (cases x)
  2.1643 +  case (Float m e)
  2.1644 +  have "\<bar>real m\<bar> * pow2 e = \<bar>real m * pow2 e\<bar>" unfolding abs_mult by auto
  2.1645 +  thus ?thesis unfolding Float abs_float_def float_abs.simps Ifloat.simps by auto
  2.1646 +qed
  2.1647 +
  2.1648 +fun floor_fl :: "float \<Rightarrow> float" where
  2.1649 +"floor_fl (Float m e) = (if 0 \<le> e then Float m e
  2.1650 +                                  else Float (m div (2 ^ (nat (-e)))) 0)"
  2.1651  
  2.1652 -lemma int_number_of_minus_sym: "- ((number_of v)::int) = number_of (uminus v)"
  2.1653 -  by simp
  2.1654 +lemma floor_fl: "Ifloat (floor_fl x) \<le> Ifloat x"
  2.1655 +proof (cases x)
  2.1656 +  case (Float m e)
  2.1657 +  show ?thesis
  2.1658 +  proof (cases "0 \<le> e")
  2.1659 +    case False
  2.1660 +    hence me_eq: "pow2 (-e) = pow2 (int (nat (-e)))" by auto
  2.1661 +    have "Ifloat (Float (m div (2 ^ (nat (-e)))) 0) = real (m div 2 ^ (nat (-e)))" unfolding Ifloat.simps by auto
  2.1662 +    also have "\<dots> \<le> real m / real ((2::int) ^ (nat (-e)))" using real_of_int_div4 .
  2.1663 +    also have "\<dots> = real m * inverse (2 ^ (nat (-e)))" unfolding power_real_number_of[symmetric] real_divide_def ..
  2.1664 +    also have "\<dots> = Ifloat (Float m e)" unfolding Ifloat.simps me_eq pow2_int pow2_neg[of e] ..
  2.1665 +    finally show ?thesis unfolding Float floor_fl.simps if_not_P[OF `\<not> 0 \<le> e`] .
  2.1666 +  next
  2.1667 +    case True thus ?thesis unfolding Float by auto
  2.1668 +  qed
  2.1669 +qed
  2.1670  
  2.1671 -lemmas intarith = int_number_of_add_sym int_number_of_minus_sym int_number_of_diff_sym int_number_of_mult_sym
  2.1672 +lemma floor_pos_exp: assumes floor: "Float m e = floor_fl x" shows "0 \<le> e"
  2.1673 +proof (cases x)
  2.1674 +  case (Float mx me)
  2.1675 +  from floor[unfolded Float floor_fl.simps] show ?thesis by (cases "0 \<le> me", auto)
  2.1676 +qed
  2.1677 +
  2.1678 +declare floor_fl.simps[simp del]
  2.1679  
  2.1680 -lemmas natarith = add_nat_number_of diff_nat_number_of mult_nat_number_of eq_nat_number_of less_nat_number_of
  2.1681 +fun ceiling_fl :: "float \<Rightarrow> float" where
  2.1682 +"ceiling_fl (Float m e) = (if 0 \<le> e then Float m e
  2.1683 +                                    else Float (m div (2 ^ (nat (-e))) + 1) 0)"
  2.1684  
  2.1685 -lemmas powerarith = nat_number_of zpower_number_of_even
  2.1686 -  zpower_number_of_odd[simplified zero_eq_Numeral0_nring one_eq_Numeral1_nring]
  2.1687 -  zpower_Pls zpower_Min
  2.1688 +lemma ceiling_fl: "Ifloat x \<le> Ifloat (ceiling_fl x)"
  2.1689 +proof (cases x)
  2.1690 +  case (Float m e)
  2.1691 +  show ?thesis
  2.1692 +  proof (cases "0 \<le> e")
  2.1693 +    case False
  2.1694 +    hence me_eq: "pow2 (-e) = pow2 (int (nat (-e)))" by auto
  2.1695 +    have "Ifloat (Float m e) = real m * inverse (2 ^ (nat (-e)))" unfolding Ifloat.simps me_eq pow2_int pow2_neg[of e] ..
  2.1696 +    also have "\<dots> = real m / real ((2::int) ^ (nat (-e)))" unfolding power_real_number_of[symmetric] real_divide_def ..
  2.1697 +    also have "\<dots> \<le> 1 + real (m div 2 ^ (nat (-e)))" using real_of_int_div3[unfolded diff_le_eq] .
  2.1698 +    also have "\<dots> = Ifloat (Float (m div (2 ^ (nat (-e))) + 1) 0)" unfolding Ifloat.simps by auto
  2.1699 +    finally show ?thesis unfolding Float ceiling_fl.simps if_not_P[OF `\<not> 0 \<le> e`] .
  2.1700 +  next
  2.1701 +    case True thus ?thesis unfolding Float by auto
  2.1702 +  qed
  2.1703 +qed
  2.1704 +
  2.1705 +declare ceiling_fl.simps[simp del]
  2.1706 +
  2.1707 +definition lb_mod :: "nat \<Rightarrow> float \<Rightarrow> float \<Rightarrow> float \<Rightarrow> float" where
  2.1708 +"lb_mod prec x ub lb = x - ceiling_fl (float_divr prec x lb) * ub"
  2.1709 +
  2.1710 +definition ub_mod :: "nat \<Rightarrow> float \<Rightarrow> float \<Rightarrow> float \<Rightarrow> float" where
  2.1711 +"ub_mod prec x ub lb = x - floor_fl (float_divl prec x ub) * lb"
  2.1712  
  2.1713 -lemmas floatarith[simplified norm_0_1] = float_add float_add_l0 float_add_r0 float_mult float_mult_l0 float_mult_r0 
  2.1714 -          float_minus float_abs zero_le_float float_pprt float_nprt pprt_lbound nprt_ubound
  2.1715 +lemma lb_mod: fixes k :: int assumes "0 \<le> Ifloat x" and "real k * y \<le> Ifloat x" (is "?k * y \<le> ?x")
  2.1716 +  assumes "0 < Ifloat lb" "Ifloat lb \<le> y" (is "?lb \<le> y") "y \<le> Ifloat ub" (is "y \<le> ?ub")
  2.1717 +  shows "Ifloat (lb_mod prec x ub lb) \<le> ?x - ?k * y"
  2.1718 +proof -
  2.1719 +  have "?lb \<le> ?ub" by (auto!)
  2.1720 +  have "0 \<le> ?lb" and "?lb \<noteq> 0" by (auto!)
  2.1721 +  have "?k * y \<le> ?x" using assms by auto
  2.1722 +  also have "\<dots> \<le> ?x / ?lb * ?ub" by (metis mult_left_mono[OF `?lb \<le> ?ub` `0 \<le> ?x`] divide_right_mono[OF _ `0 \<le> ?lb` ] times_divide_eq_left nonzero_mult_divide_cancel_right[OF `?lb \<noteq> 0`])
  2.1723 +  also have "\<dots> \<le> Ifloat (ceiling_fl (float_divr prec x lb)) * ?ub" by (metis mult_right_mono order_trans `0 \<le> ?lb` `?lb \<le> ?ub` float_divr ceiling_fl)
  2.1724 +  finally show ?thesis unfolding lb_mod_def Ifloat_sub Ifloat_mult by auto
  2.1725 +qed
  2.1726  
  2.1727 -(* for use with the compute oracle *)
  2.1728 -lemmas arith = binarith intarith intarithrel natarith powerarith floatarith not_false_eq_true not_true_eq_false
  2.1729 +lemma ub_mod: fixes k :: int assumes "0 \<le> Ifloat x" and "Ifloat x \<le> real k * y" (is "?x \<le> ?k * y")
  2.1730 +  assumes "0 < Ifloat lb" "Ifloat lb \<le> y" (is "?lb \<le> y") "y \<le> Ifloat ub" (is "y \<le> ?ub")
  2.1731 +  shows "?x - ?k * y \<le> Ifloat (ub_mod prec x ub lb)"
  2.1732 +proof -
  2.1733 +  have "?lb \<le> ?ub" by (auto!)
  2.1734 +  hence "0 \<le> ?lb" and "0 \<le> ?ub" and "?ub \<noteq> 0" by (auto!)
  2.1735 +  have "Ifloat (floor_fl (float_divl prec x ub)) * ?lb \<le> ?x / ?ub * ?lb" by (metis mult_right_mono order_trans `0 \<le> ?lb` `?lb \<le> ?ub` float_divl floor_fl)
  2.1736 +  also have "\<dots> \<le> ?x" by (metis mult_left_mono[OF `?lb \<le> ?ub` `0 \<le> ?x`] divide_right_mono[OF _ `0 \<le> ?ub` ] times_divide_eq_left nonzero_mult_divide_cancel_right[OF `?ub \<noteq> 0`])
  2.1737 +  also have "\<dots> \<le> ?k * y" using assms by auto
  2.1738 +  finally show ?thesis unfolding ub_mod_def Ifloat_sub Ifloat_mult by auto
  2.1739 +qed
  2.1740  
  2.1741 -use "~~/src/HOL/Tools/float_arith.ML"
  2.1742 +lemma le_float_def': "f \<le> g = (case f - g of Float a b \<Rightarrow> a \<le> 0)"
  2.1743 +proof -
  2.1744 +  have le_transfer: "(f \<le> g) = (Ifloat (f - g) \<le> 0)" by (auto simp add: le_float_def)
  2.1745 +  from float_split[of "f - g"] obtain a b where f_diff_g: "f - g = Float a b" by auto
  2.1746 +  with le_transfer have le_transfer': "f \<le> g = (Ifloat (Float a b) \<le> 0)" by simp
  2.1747 +  show ?thesis by (simp add: le_transfer' f_diff_g float_le_zero)
  2.1748 +qed
  2.1749 +
  2.1750 +lemma float_less_zero:
  2.1751 +  "(Ifloat (Float a b) < 0) = (a < 0)"
  2.1752 +  apply (auto simp add: mult_less_0_iff Ifloat.simps)
  2.1753 +  done
  2.1754 +
  2.1755 +lemma less_float_def': "f < g = (case f - g of Float a b \<Rightarrow> a < 0)"
  2.1756 +proof -
  2.1757 +  have less_transfer: "(f < g) = (Ifloat (f - g) < 0)" by (auto simp add: less_float_def)
  2.1758 +  from float_split[of "f - g"] obtain a b where f_diff_g: "f - g = Float a b" by auto
  2.1759 +  with less_transfer have less_transfer': "f < g = (Ifloat (Float a b) < 0)" by simp
  2.1760 +  show ?thesis by (simp add: less_transfer' f_diff_g float_less_zero)
  2.1761 +qed
  2.1762  
  2.1763  end
     3.1 --- a/src/HOL/Matrix/cplex/Cplex.thy	Thu Feb 05 11:34:42 2009 +0100
     3.2 +++ b/src/HOL/Matrix/cplex/Cplex.thy	Thu Feb 05 11:45:15 2009 +0100
     3.3 @@ -4,7 +4,7 @@
     3.4  *)
     3.5  
     3.6  theory Cplex 
     3.7 -imports SparseMatrix LP "~~/src/HOL/Real/Float" "~~/src/HOL/Tools/ComputeNumeral"
     3.8 +imports SparseMatrix LP "~~/src/HOL/Tools/ComputeFloat" "~~/src/HOL/Tools/ComputeNumeral"
     3.9  uses "Cplex_tools.ML" "CplexMatrixConverter.ML" "FloatSparseMatrixBuilder.ML"
    3.10    "fspmlp.ML" ("matrixlp.ML")
    3.11  begin
     4.1 --- a/src/HOL/Matrix/cplex/matrixlp.ML	Thu Feb 05 11:34:42 2009 +0100
     4.2 +++ b/src/HOL/Matrix/cplex/matrixlp.ML	Thu Feb 05 11:45:15 2009 +0100
     4.3 @@ -67,7 +67,7 @@
     4.4  
     4.5  local
     4.6      val ths = ComputeHOL.prep_thms @{thms "ComputeHOL.compute_list_case" "ComputeHOL.compute_let"
     4.7 -      "ComputeHOL.compute_if" "Float.arith" "SparseMatrix.sparse_row_matrix_arith_simps"
     4.8 +      "ComputeHOL.compute_if" "ComputeFloat.arith" "SparseMatrix.sparse_row_matrix_arith_simps"
     4.9        "ComputeHOL.compute_bool" "ComputeHOL.compute_pair"
    4.10        "SparseMatrix.sorted_sp_simps" "ComputeNumeral.number_norm"
    4.11        "ComputeNumeral.natnorm"};
     5.1 --- /dev/null	Thu Jan 01 00:00:00 1970 +0000
     5.2 +++ b/src/HOL/Tools/ComputeFloat.thy	Thu Feb 05 11:45:15 2009 +0100
     5.3 @@ -0,0 +1,568 @@
     5.4 +(*  Title:  HOL/Tools/ComputeFloat.thy
     5.5 +    Author: Steven Obua
     5.6 +*)
     5.7 +
     5.8 +header {* Floating Point Representation of the Reals *}
     5.9 +
    5.10 +theory ComputeFloat
    5.11 +imports Complex_Main
    5.12 +uses "~~/src/Tools/float.ML" ("~~/src/HOL/Tools/float_arith.ML")
    5.13 +begin
    5.14 +
    5.15 +definition
    5.16 +  pow2 :: "int \<Rightarrow> real" where
    5.17 +  "pow2 a = (if (0 <= a) then (2^(nat a)) else (inverse (2^(nat (-a)))))"
    5.18 +
    5.19 +definition
    5.20 +  float :: "int * int \<Rightarrow> real" where
    5.21 +  "float x = real (fst x) * pow2 (snd x)"
    5.22 +
    5.23 +lemma pow2_0[simp]: "pow2 0 = 1"
    5.24 +by (simp add: pow2_def)
    5.25 +
    5.26 +lemma pow2_1[simp]: "pow2 1 = 2"
    5.27 +by (simp add: pow2_def)
    5.28 +
    5.29 +lemma pow2_neg: "pow2 x = inverse (pow2 (-x))"
    5.30 +by (simp add: pow2_def)
    5.31 +
    5.32 +lemma pow2_add1: "pow2 (1 + a) = 2 * (pow2 a)"
    5.33 +proof -
    5.34 +  have h: "! n. nat (2 + int n) - Suc 0 = nat (1 + int n)" by arith
    5.35 +  have g: "! a b. a - -1 = a + (1::int)" by arith
    5.36 +  have pos: "! n. pow2 (int n + 1) = 2 * pow2 (int n)"
    5.37 +    apply (auto, induct_tac n)
    5.38 +    apply (simp_all add: pow2_def)
    5.39 +    apply (rule_tac m1="2" and n1="nat (2 + int na)" in ssubst[OF realpow_num_eq_if])
    5.40 +    by (auto simp add: h)
    5.41 +  show ?thesis
    5.42 +  proof (induct a)
    5.43 +    case (1 n)
    5.44 +    from pos show ?case by (simp add: algebra_simps)
    5.45 +  next
    5.46 +    case (2 n)
    5.47 +    show ?case
    5.48 +      apply (auto)
    5.49 +      apply (subst pow2_neg[of "- int n"])
    5.50 +      apply (subst pow2_neg[of "-1 - int n"])
    5.51 +      apply (auto simp add: g pos)
    5.52 +      done
    5.53 +  qed
    5.54 +qed
    5.55 +
    5.56 +lemma pow2_add: "pow2 (a+b) = (pow2 a) * (pow2 b)"
    5.57 +proof (induct b)
    5.58 +  case (1 n)
    5.59 +  show ?case
    5.60 +  proof (induct n)
    5.61 +    case 0
    5.62 +    show ?case by simp
    5.63 +  next
    5.64 +    case (Suc m)
    5.65 +    show ?case by (auto simp add: algebra_simps pow2_add1 prems)
    5.66 +  qed
    5.67 +next
    5.68 +  case (2 n)
    5.69 +  show ?case
    5.70 +  proof (induct n)
    5.71 +    case 0
    5.72 +    show ?case
    5.73 +      apply (auto)
    5.74 +      apply (subst pow2_neg[of "a + -1"])
    5.75 +      apply (subst pow2_neg[of "-1"])
    5.76 +      apply (simp)
    5.77 +      apply (insert pow2_add1[of "-a"])
    5.78 +      apply (simp add: algebra_simps)
    5.79 +      apply (subst pow2_neg[of "-a"])
    5.80 +      apply (simp)
    5.81 +      done
    5.82 +    case (Suc m)
    5.83 +    have a: "int m - (a + -2) =  1 + (int m - a + 1)" by arith
    5.84 +    have b: "int m - -2 = 1 + (int m + 1)" by arith
    5.85 +    show ?case
    5.86 +      apply (auto)
    5.87 +      apply (subst pow2_neg[of "a + (-2 - int m)"])
    5.88 +      apply (subst pow2_neg[of "-2 - int m"])
    5.89 +      apply (auto simp add: algebra_simps)
    5.90 +      apply (subst a)
    5.91 +      apply (subst b)
    5.92 +      apply (simp only: pow2_add1)
    5.93 +      apply (subst pow2_neg[of "int m - a + 1"])
    5.94 +      apply (subst pow2_neg[of "int m + 1"])
    5.95 +      apply auto
    5.96 +      apply (insert prems)
    5.97 +      apply (auto simp add: algebra_simps)
    5.98 +      done
    5.99 +  qed
   5.100 +qed
   5.101 +
   5.102 +lemma "float (a, e) + float (b, e) = float (a + b, e)"
   5.103 +by (simp add: float_def algebra_simps)
   5.104 +
   5.105 +definition
   5.106 +  int_of_real :: "real \<Rightarrow> int" where
   5.107 +  "int_of_real x = (SOME y. real y = x)"
   5.108 +
   5.109 +definition
   5.110 +  real_is_int :: "real \<Rightarrow> bool" where
   5.111 +  "real_is_int x = (EX (u::int). x = real u)"
   5.112 +
   5.113 +lemma real_is_int_def2: "real_is_int x = (x = real (int_of_real x))"
   5.114 +by (auto simp add: real_is_int_def int_of_real_def)
   5.115 +
   5.116 +lemma float_transfer: "real_is_int ((real a)*(pow2 c)) \<Longrightarrow> float (a, b) = float (int_of_real ((real a)*(pow2 c)), b - c)"
   5.117 +by (simp add: float_def real_is_int_def2 pow2_add[symmetric])
   5.118 +
   5.119 +lemma pow2_int: "pow2 (int c) = 2^c"
   5.120 +by (simp add: pow2_def)
   5.121 +
   5.122 +lemma float_transfer_nat: "float (a, b) = float (a * 2^c, b - int c)"
   5.123 +by (simp add: float_def pow2_int[symmetric] pow2_add[symmetric])
   5.124 +
   5.125 +lemma real_is_int_real[simp]: "real_is_int (real (x::int))"
   5.126 +by (auto simp add: real_is_int_def int_of_real_def)
   5.127 +
   5.128 +lemma int_of_real_real[simp]: "int_of_real (real x) = x"
   5.129 +by (simp add: int_of_real_def)
   5.130 +
   5.131 +lemma real_int_of_real[simp]: "real_is_int x \<Longrightarrow> real (int_of_real x) = x"
   5.132 +by (auto simp add: int_of_real_def real_is_int_def)
   5.133 +
   5.134 +lemma real_is_int_add_int_of_real: "real_is_int a \<Longrightarrow> real_is_int b \<Longrightarrow> (int_of_real (a+b)) = (int_of_real a) + (int_of_real b)"
   5.135 +by (auto simp add: int_of_real_def real_is_int_def)
   5.136 +
   5.137 +lemma real_is_int_add[simp]: "real_is_int a \<Longrightarrow> real_is_int b \<Longrightarrow> real_is_int (a+b)"
   5.138 +apply (subst real_is_int_def2)
   5.139 +apply (simp add: real_is_int_add_int_of_real real_int_of_real)
   5.140 +done
   5.141 +
   5.142 +lemma int_of_real_sub: "real_is_int a \<Longrightarrow> real_is_int b \<Longrightarrow> (int_of_real (a-b)) = (int_of_real a) - (int_of_real b)"
   5.143 +by (auto simp add: int_of_real_def real_is_int_def)
   5.144 +
   5.145 +lemma real_is_int_sub[simp]: "real_is_int a \<Longrightarrow> real_is_int b \<Longrightarrow> real_is_int (a-b)"
   5.146 +apply (subst real_is_int_def2)
   5.147 +apply (simp add: int_of_real_sub real_int_of_real)
   5.148 +done
   5.149 +
   5.150 +lemma real_is_int_rep: "real_is_int x \<Longrightarrow> ?! (a::int). real a = x"
   5.151 +by (auto simp add: real_is_int_def)
   5.152 +
   5.153 +lemma int_of_real_mult:
   5.154 +  assumes "real_is_int a" "real_is_int b"
   5.155 +  shows "(int_of_real (a*b)) = (int_of_real a) * (int_of_real b)"
   5.156 +proof -
   5.157 +  from prems have a: "?! (a'::int). real a' = a" by (rule_tac real_is_int_rep, auto)
   5.158 +  from prems have b: "?! (b'::int). real b' = b" by (rule_tac real_is_int_rep, auto)
   5.159 +  from a obtain a'::int where a':"a = real a'" by auto
   5.160 +  from b obtain b'::int where b':"b = real b'" by auto
   5.161 +  have r: "real a' * real b' = real (a' * b')" by auto
   5.162 +  show ?thesis
   5.163 +    apply (simp add: a' b')
   5.164 +    apply (subst r)
   5.165 +    apply (simp only: int_of_real_real)
   5.166 +    done
   5.167 +qed
   5.168 +
   5.169 +lemma real_is_int_mult[simp]: "real_is_int a \<Longrightarrow> real_is_int b \<Longrightarrow> real_is_int (a*b)"
   5.170 +apply (subst real_is_int_def2)
   5.171 +apply (simp add: int_of_real_mult)
   5.172 +done
   5.173 +
   5.174 +lemma real_is_int_0[simp]: "real_is_int (0::real)"
   5.175 +by (simp add: real_is_int_def int_of_real_def)
   5.176 +
   5.177 +lemma real_is_int_1[simp]: "real_is_int (1::real)"
   5.178 +proof -
   5.179 +  have "real_is_int (1::real) = real_is_int(real (1::int))" by auto
   5.180 +  also have "\<dots> = True" by (simp only: real_is_int_real)
   5.181 +  ultimately show ?thesis by auto
   5.182 +qed
   5.183 +
   5.184 +lemma real_is_int_n1: "real_is_int (-1::real)"
   5.185 +proof -
   5.186 +  have "real_is_int (-1::real) = real_is_int(real (-1::int))" by auto
   5.187 +  also have "\<dots> = True" by (simp only: real_is_int_real)
   5.188 +  ultimately show ?thesis by auto
   5.189 +qed
   5.190 +
   5.191 +lemma real_is_int_number_of[simp]: "real_is_int ((number_of \<Colon> int \<Rightarrow> real) x)"
   5.192 +proof -
   5.193 +  have neg1: "real_is_int (-1::real)"
   5.194 +  proof -
   5.195 +    have "real_is_int (-1::real) = real_is_int(real (-1::int))" by auto
   5.196 +    also have "\<dots> = True" by (simp only: real_is_int_real)
   5.197 +    ultimately show ?thesis by auto
   5.198 +  qed
   5.199 +
   5.200 +  {
   5.201 +    fix x :: int
   5.202 +    have "real_is_int ((number_of \<Colon> int \<Rightarrow> real) x)"
   5.203 +      unfolding number_of_eq
   5.204 +      apply (induct x)
   5.205 +      apply (induct_tac n)
   5.206 +      apply (simp)
   5.207 +      apply (simp)
   5.208 +      apply (induct_tac n)
   5.209 +      apply (simp add: neg1)
   5.210 +    proof -
   5.211 +      fix n :: nat
   5.212 +      assume rn: "(real_is_int (of_int (- (int (Suc n)))))"
   5.213 +      have s: "-(int (Suc (Suc n))) = -1 + - (int (Suc n))" by simp
   5.214 +      show "real_is_int (of_int (- (int (Suc (Suc n)))))"
   5.215 +        apply (simp only: s of_int_add)
   5.216 +        apply (rule real_is_int_add)
   5.217 +        apply (simp add: neg1)
   5.218 +        apply (simp only: rn)
   5.219 +        done
   5.220 +    qed
   5.221 +  }
   5.222 +  note Abs_Bin = this
   5.223 +  {
   5.224 +    fix x :: int
   5.225 +    have "? u. x = u"
   5.226 +      apply (rule exI[where x = "x"])
   5.227 +      apply (simp)
   5.228 +      done
   5.229 +  }
   5.230 +  then obtain u::int where "x = u" by auto
   5.231 +  with Abs_Bin show ?thesis by auto
   5.232 +qed
   5.233 +
   5.234 +lemma int_of_real_0[simp]: "int_of_real (0::real) = (0::int)"
   5.235 +by (simp add: int_of_real_def)
   5.236 +
   5.237 +lemma int_of_real_1[simp]: "int_of_real (1::real) = (1::int)"
   5.238 +proof -
   5.239 +  have 1: "(1::real) = real (1::int)" by auto
   5.240 +  show ?thesis by (simp only: 1 int_of_real_real)
   5.241 +qed
   5.242 +
   5.243 +lemma int_of_real_number_of[simp]: "int_of_real (number_of b) = number_of b"
   5.244 +proof -
   5.245 +  have "real_is_int (number_of b)" by simp
   5.246 +  then have uu: "?! u::int. number_of b = real u" by (auto simp add: real_is_int_rep)
   5.247 +  then obtain u::int where u:"number_of b = real u" by auto
   5.248 +  have "number_of b = real ((number_of b)::int)"
   5.249 +    by (simp add: number_of_eq real_of_int_def)
   5.250 +  have ub: "number_of b = real ((number_of b)::int)"
   5.251 +    by (simp add: number_of_eq real_of_int_def)
   5.252 +  from uu u ub have unb: "u = number_of b"
   5.253 +    by blast
   5.254 +  have "int_of_real (number_of b) = u" by (simp add: u)
   5.255 +  with unb show ?thesis by simp
   5.256 +qed
   5.257 +
   5.258 +lemma float_transfer_even: "even a \<Longrightarrow> float (a, b) = float (a div 2, b+1)"
   5.259 +  apply (subst float_transfer[where a="a" and b="b" and c="-1", simplified])
   5.260 +  apply (simp_all add: pow2_def even_def real_is_int_def algebra_simps)
   5.261 +  apply (auto)
   5.262 +proof -
   5.263 +  fix q::int
   5.264 +  have a:"b - (-1\<Colon>int) = (1\<Colon>int) + b" by arith
   5.265 +  show "(float (q, (b - (-1\<Colon>int)))) = (float (q, ((1\<Colon>int) + b)))"
   5.266 +    by (simp add: a)
   5.267 +qed
   5.268 +
   5.269 +lemma int_div_zdiv: "int (a div b) = (int a) div (int b)"
   5.270 +by (rule zdiv_int)
   5.271 +
   5.272 +lemma int_mod_zmod: "int (a mod b) = (int a) mod (int b)"
   5.273 +by (rule zmod_int)
   5.274 +
   5.275 +lemma abs_div_2_less: "a \<noteq> 0 \<Longrightarrow> a \<noteq> -1 \<Longrightarrow> abs((a::int) div 2) < abs a"
   5.276 +by arith
   5.277 +
   5.278 +function norm_float :: "int \<Rightarrow> int \<Rightarrow> int \<times> int" where
   5.279 +  "norm_float a b = (if a \<noteq> 0 \<and> even a then norm_float (a div 2) (b + 1)
   5.280 +    else if a = 0 then (0, 0) else (a, b))"
   5.281 +by auto
   5.282 +
   5.283 +termination by (relation "measure (nat o abs o fst)")
   5.284 +  (auto intro: abs_div_2_less)
   5.285 +
   5.286 +lemma norm_float: "float x = float (split norm_float x)"
   5.287 +proof -
   5.288 +  {
   5.289 +    fix a b :: int
   5.290 +    have norm_float_pair: "float (a, b) = float (norm_float a b)"
   5.291 +    proof (induct a b rule: norm_float.induct)
   5.292 +      case (1 u v)
   5.293 +      show ?case
   5.294 +      proof cases
   5.295 +        assume u: "u \<noteq> 0 \<and> even u"
   5.296 +        with prems have ind: "float (u div 2, v + 1) = float (norm_float (u div 2) (v + 1))" by auto
   5.297 +        with u have "float (u,v) = float (u div 2, v+1)" by (simp add: float_transfer_even)
   5.298 +        then show ?thesis
   5.299 +          apply (subst norm_float.simps)
   5.300 +          apply (simp add: ind)
   5.301 +          done
   5.302 +      next
   5.303 +        assume "~(u \<noteq> 0 \<and> even u)"
   5.304 +        then show ?thesis
   5.305 +          by (simp add: prems float_def)
   5.306 +      qed
   5.307 +    qed
   5.308 +  }
   5.309 +  note helper = this
   5.310 +  have "? a b. x = (a,b)" by auto
   5.311 +  then obtain a b where "x = (a, b)" by blast
   5.312 +  then show ?thesis by (simp add: helper)
   5.313 +qed
   5.314 +
   5.315 +lemma float_add_l0: "float (0, e) + x = x"
   5.316 +  by (simp add: float_def)
   5.317 +
   5.318 +lemma float_add_r0: "x + float (0, e) = x"
   5.319 +  by (simp add: float_def)
   5.320 +
   5.321 +lemma float_add:
   5.322 +  "float (a1, e1) + float (a2, e2) =
   5.323 +  (if e1<=e2 then float (a1+a2*2^(nat(e2-e1)), e1)
   5.324 +  else float (a1*2^(nat (e1-e2))+a2, e2))"
   5.325 +  apply (simp add: float_def algebra_simps)
   5.326 +  apply (auto simp add: pow2_int[symmetric] pow2_add[symmetric])
   5.327 +  done
   5.328 +
   5.329 +lemma float_add_assoc1:
   5.330 +  "(x + float (y1, e1)) + float (y2, e2) = (float (y1, e1) + float (y2, e2)) + x"
   5.331 +  by simp
   5.332 +
   5.333 +lemma float_add_assoc2:
   5.334 +  "(float (y1, e1) + x) + float (y2, e2) = (float (y1, e1) + float (y2, e2)) + x"
   5.335 +  by simp
   5.336 +
   5.337 +lemma float_add_assoc3:
   5.338 +  "float (y1, e1) + (x + float (y2, e2)) = (float (y1, e1) + float (y2, e2)) + x"
   5.339 +  by simp
   5.340 +
   5.341 +lemma float_add_assoc4:
   5.342 +  "float (y1, e1) + (float (y2, e2) + x) = (float (y1, e1) + float (y2, e2)) + x"
   5.343 +  by simp
   5.344 +
   5.345 +lemma float_mult_l0: "float (0, e) * x = float (0, 0)"
   5.346 +  by (simp add: float_def)
   5.347 +
   5.348 +lemma float_mult_r0: "x * float (0, e) = float (0, 0)"
   5.349 +  by (simp add: float_def)
   5.350 +
   5.351 +definition 
   5.352 +  lbound :: "real \<Rightarrow> real"
   5.353 +where
   5.354 +  "lbound x = min 0 x"
   5.355 +
   5.356 +definition
   5.357 +  ubound :: "real \<Rightarrow> real"
   5.358 +where
   5.359 +  "ubound x = max 0 x"
   5.360 +
   5.361 +lemma lbound: "lbound x \<le> x"   
   5.362 +  by (simp add: lbound_def)
   5.363 +
   5.364 +lemma ubound: "x \<le> ubound x"
   5.365 +  by (simp add: ubound_def)
   5.366 +
   5.367 +lemma float_mult:
   5.368 +  "float (a1, e1) * float (a2, e2) =
   5.369 +  (float (a1 * a2, e1 + e2))"
   5.370 +  by (simp add: float_def pow2_add)
   5.371 +
   5.372 +lemma float_minus:
   5.373 +  "- (float (a,b)) = float (-a, b)"
   5.374 +  by (simp add: float_def)
   5.375 +
   5.376 +lemma zero_less_pow2:
   5.377 +  "0 < pow2 x"
   5.378 +proof -
   5.379 +  {
   5.380 +    fix y
   5.381 +    have "0 <= y \<Longrightarrow> 0 < pow2 y"
   5.382 +      by (induct y, induct_tac n, simp_all add: pow2_add)
   5.383 +  }
   5.384 +  note helper=this
   5.385 +  show ?thesis
   5.386 +    apply (case_tac "0 <= x")
   5.387 +    apply (simp add: helper)
   5.388 +    apply (subst pow2_neg)
   5.389 +    apply (simp add: helper)
   5.390 +    done
   5.391 +qed
   5.392 +
   5.393 +lemma zero_le_float:
   5.394 +  "(0 <= float (a,b)) = (0 <= a)"
   5.395 +  apply (auto simp add: float_def)
   5.396 +  apply (auto simp add: zero_le_mult_iff zero_less_pow2)
   5.397 +  apply (insert zero_less_pow2[of b])
   5.398 +  apply (simp_all)
   5.399 +  done
   5.400 +
   5.401 +lemma float_le_zero:
   5.402 +  "(float (a,b) <= 0) = (a <= 0)"
   5.403 +  apply (auto simp add: float_def)
   5.404 +  apply (auto simp add: mult_le_0_iff)
   5.405 +  apply (insert zero_less_pow2[of b])
   5.406 +  apply auto
   5.407 +  done
   5.408 +
   5.409 +lemma float_abs:
   5.410 +  "abs (float (a,b)) = (if 0 <= a then (float (a,b)) else (float (-a,b)))"
   5.411 +  apply (auto simp add: abs_if)
   5.412 +  apply (simp_all add: zero_le_float[symmetric, of a b] float_minus)
   5.413 +  done
   5.414 +
   5.415 +lemma float_zero:
   5.416 +  "float (0, b) = 0"
   5.417 +  by (simp add: float_def)
   5.418 +
   5.419 +lemma float_pprt:
   5.420 +  "pprt (float (a, b)) = (if 0 <= a then (float (a,b)) else (float (0, b)))"
   5.421 +  by (auto simp add: zero_le_float float_le_zero float_zero)
   5.422 +
   5.423 +lemma pprt_lbound: "pprt (lbound x) = float (0, 0)"
   5.424 +  apply (simp add: float_def)
   5.425 +  apply (rule pprt_eq_0)
   5.426 +  apply (simp add: lbound_def)
   5.427 +  done
   5.428 +
   5.429 +lemma nprt_ubound: "nprt (ubound x) = float (0, 0)"
   5.430 +  apply (simp add: float_def)
   5.431 +  apply (rule nprt_eq_0)
   5.432 +  apply (simp add: ubound_def)
   5.433 +  done
   5.434 +
   5.435 +lemma float_nprt:
   5.436 +  "nprt (float (a, b)) = (if 0 <= a then (float (0,b)) else (float (a, b)))"
   5.437 +  by (auto simp add: zero_le_float float_le_zero float_zero)
   5.438 +
   5.439 +lemma norm_0_1: "(0::_::number_ring) = Numeral0 & (1::_::number_ring) = Numeral1"
   5.440 +  by auto
   5.441 +
   5.442 +lemma add_left_zero: "0 + a = (a::'a::comm_monoid_add)"
   5.443 +  by simp
   5.444 +
   5.445 +lemma add_right_zero: "a + 0 = (a::'a::comm_monoid_add)"
   5.446 +  by simp
   5.447 +
   5.448 +lemma mult_left_one: "1 * a = (a::'a::semiring_1)"
   5.449 +  by simp
   5.450 +
   5.451 +lemma mult_right_one: "a * 1 = (a::'a::semiring_1)"
   5.452 +  by simp
   5.453 +
   5.454 +lemma int_pow_0: "(a::int)^(Numeral0) = 1"
   5.455 +  by simp
   5.456 +
   5.457 +lemma int_pow_1: "(a::int)^(Numeral1) = a"
   5.458 +  by simp
   5.459 +
   5.460 +lemma zero_eq_Numeral0_nring: "(0::'a::number_ring) = Numeral0"
   5.461 +  by simp
   5.462 +
   5.463 +lemma one_eq_Numeral1_nring: "(1::'a::number_ring) = Numeral1"
   5.464 +  by simp
   5.465 +
   5.466 +lemma zero_eq_Numeral0_nat: "(0::nat) = Numeral0"
   5.467 +  by simp
   5.468 +
   5.469 +lemma one_eq_Numeral1_nat: "(1::nat) = Numeral1"
   5.470 +  by simp
   5.471 +
   5.472 +lemma zpower_Pls: "(z::int)^Numeral0 = Numeral1"
   5.473 +  by simp
   5.474 +
   5.475 +lemma zpower_Min: "(z::int)^((-1)::nat) = Numeral1"
   5.476 +proof -
   5.477 +  have 1:"((-1)::nat) = 0"
   5.478 +    by simp
   5.479 +  show ?thesis by (simp add: 1)
   5.480 +qed
   5.481 +
   5.482 +lemma fst_cong: "a=a' \<Longrightarrow> fst (a,b) = fst (a',b)"
   5.483 +  by simp
   5.484 +
   5.485 +lemma snd_cong: "b=b' \<Longrightarrow> snd (a,b) = snd (a,b')"
   5.486 +  by simp
   5.487 +
   5.488 +lemma lift_bool: "x \<Longrightarrow> x=True"
   5.489 +  by simp
   5.490 +
   5.491 +lemma nlift_bool: "~x \<Longrightarrow> x=False"
   5.492 +  by simp
   5.493 +
   5.494 +lemma not_false_eq_true: "(~ False) = True" by simp
   5.495 +
   5.496 +lemma not_true_eq_false: "(~ True) = False" by simp
   5.497 +
   5.498 +lemmas binarith =
   5.499 +  normalize_bin_simps
   5.500 +  pred_bin_simps succ_bin_simps
   5.501 +  add_bin_simps minus_bin_simps mult_bin_simps
   5.502 +
   5.503 +lemma int_eq_number_of_eq:
   5.504 +  "(((number_of v)::int)=(number_of w)) = iszero ((number_of (v + uminus w))::int)"
   5.505 +  by (rule eq_number_of_eq)
   5.506 +
   5.507 +lemma int_iszero_number_of_Pls: "iszero (Numeral0::int)"
   5.508 +  by (simp only: iszero_number_of_Pls)
   5.509 +
   5.510 +lemma int_nonzero_number_of_Min: "~(iszero ((-1)::int))"
   5.511 +  by simp
   5.512 +
   5.513 +lemma int_iszero_number_of_Bit0: "iszero ((number_of (Int.Bit0 w))::int) = iszero ((number_of w)::int)"
   5.514 +  by simp
   5.515 +
   5.516 +lemma int_iszero_number_of_Bit1: "\<not> iszero ((number_of (Int.Bit1 w))::int)"
   5.517 +  by simp
   5.518 +
   5.519 +lemma int_less_number_of_eq_neg: "(((number_of x)::int) < number_of y) = neg ((number_of (x + (uminus y)))::int)"
   5.520 +  unfolding neg_def number_of_is_id by simp
   5.521 +
   5.522 +lemma int_not_neg_number_of_Pls: "\<not> (neg (Numeral0::int))"
   5.523 +  by simp
   5.524 +
   5.525 +lemma int_neg_number_of_Min: "neg (-1::int)"
   5.526 +  by simp
   5.527 +
   5.528 +lemma int_neg_number_of_Bit0: "neg ((number_of (Int.Bit0 w))::int) = neg ((number_of w)::int)"
   5.529 +  by simp
   5.530 +
   5.531 +lemma int_neg_number_of_Bit1: "neg ((number_of (Int.Bit1 w))::int) = neg ((number_of w)::int)"
   5.532 +  by simp
   5.533 +
   5.534 +lemma int_le_number_of_eq: "(((number_of x)::int) \<le> number_of y) = (\<not> neg ((number_of (y + (uminus x)))::int))"
   5.535 +  unfolding neg_def number_of_is_id by (simp add: not_less)
   5.536 +
   5.537 +lemmas intarithrel =
   5.538 +  int_eq_number_of_eq
   5.539 +  lift_bool[OF int_iszero_number_of_Pls] nlift_bool[OF int_nonzero_number_of_Min] int_iszero_number_of_Bit0
   5.540 +  lift_bool[OF int_iszero_number_of_Bit1] int_less_number_of_eq_neg nlift_bool[OF int_not_neg_number_of_Pls] lift_bool[OF int_neg_number_of_Min]
   5.541 +  int_neg_number_of_Bit0 int_neg_number_of_Bit1 int_le_number_of_eq
   5.542 +
   5.543 +lemma int_number_of_add_sym: "((number_of v)::int) + number_of w = number_of (v + w)"
   5.544 +  by simp
   5.545 +
   5.546 +lemma int_number_of_diff_sym: "((number_of v)::int) - number_of w = number_of (v + (uminus w))"
   5.547 +  by simp
   5.548 +
   5.549 +lemma int_number_of_mult_sym: "((number_of v)::int) * number_of w = number_of (v * w)"
   5.550 +  by simp
   5.551 +
   5.552 +lemma int_number_of_minus_sym: "- ((number_of v)::int) = number_of (uminus v)"
   5.553 +  by simp
   5.554 +
   5.555 +lemmas intarith = int_number_of_add_sym int_number_of_minus_sym int_number_of_diff_sym int_number_of_mult_sym
   5.556 +
   5.557 +lemmas natarith = add_nat_number_of diff_nat_number_of mult_nat_number_of eq_nat_number_of less_nat_number_of
   5.558 +
   5.559 +lemmas powerarith = nat_number_of zpower_number_of_even
   5.560 +  zpower_number_of_odd[simplified zero_eq_Numeral0_nring one_eq_Numeral1_nring]
   5.561 +  zpower_Pls zpower_Min
   5.562 +
   5.563 +lemmas floatarith[simplified norm_0_1] = float_add float_add_l0 float_add_r0 float_mult float_mult_l0 float_mult_r0 
   5.564 +          float_minus float_abs zero_le_float float_pprt float_nprt pprt_lbound nprt_ubound
   5.565 +
   5.566 +(* for use with the compute oracle *)
   5.567 +lemmas arith = binarith intarith intarithrel natarith powerarith floatarith not_false_eq_true not_true_eq_false
   5.568 +
   5.569 +use "~~/src/HOL/Tools/float_arith.ML"
   5.570 +
   5.571 +end
     6.1 --- a/src/HOL/Tools/ComputeNumeral.thy	Thu Feb 05 11:34:42 2009 +0100
     6.2 +++ b/src/HOL/Tools/ComputeNumeral.thy	Thu Feb 05 11:45:15 2009 +0100
     6.3 @@ -1,5 +1,5 @@
     6.4  theory ComputeNumeral
     6.5 -imports ComputeHOL Float
     6.6 +imports ComputeHOL ComputeFloat
     6.7  begin
     6.8  
     6.9  (* normalization of bit strings *)
     7.1 --- a/src/HOL/Tools/float_arith.ML	Thu Feb 05 11:34:42 2009 +0100
     7.2 +++ b/src/HOL/Tools/float_arith.ML	Thu Feb 05 11:45:15 2009 +0100
     7.3 @@ -1,5 +1,4 @@
     7.4 -(*  Title:  HOL/Real/float_arith.ML
     7.5 -    ID:     $Id$
     7.6 +(*  Title:  HOL/Tools/float_arith.ML
     7.7      Author: Steven Obua
     7.8  *)
     7.9