formal_ineqs/arith/float.hl

   1 (* =========================================================== *)
   2 (* Formal floating point arithmetic                            *)
   3 (* Author: Alexey Solovyev                                     *)
   4 (* Date: 2012-10-27                                            *)
   5 (* =========================================================== *)
   6
   7 (* Dependencies *)
   8 needs "arith/nat.hl";;
   9 needs "arith/num_exp_theory.hl";;
  10 needs "arith/float_theory.hl";;
  11 needs "arith/interval_arith.hl";;
  12 needs "misc/vars.hl";;
  13
  14 (* FLOOR_DIV_DIV *)
  15 needs "Library/floor.ml";;
  16 (* sqrt and its properties *)
  17 needs "Multivariate/vectors.ml";;
  18
  19 prioritize_real();;
  20
  21 module type Arith_float_sig =
  22   sig
  23     val mk_num_exp : term -> term -> term
  24     val dest_num_exp : term -> (term * term)
  25     val dest_float : term -> (string * term * term)
  26     val float_lt0 : term -> thm
  27     val float_gt0 : term -> thm
  28     val float_lt : term -> term -> thm
  29     val float_le0 : term -> thm
  30     val float_ge0 : term -> thm
  31     val float_le : term -> term -> thm
  32     val float_min : term -> term -> thm
  33     val float_max : term -> term -> thm
  34     val float_min_max : term -> term -> (thm * thm)
  35     val float_mul_eq : term -> term -> thm
  36     val float_mul_lo : int -> term -> term -> thm
  37     val float_mul_hi : int -> term -> term -> thm
  38     val float_div_lo : int -> term -> term -> thm
  39     val float_div_hi : int -> term -> term -> thm
  40     val float_add_lo : int -> term -> term -> thm
  41     val float_add_hi : int -> term -> term -> thm
  42     val float_sub_lo : int -> term -> term -> thm
  43     val float_sub_hi : int -> term -> term -> thm
  44     val float_sqrt_lo : int -> term -> thm
  45     val float_sqrt_hi : int -> term -> thm
  46
  47     val reset_stat : unit -> unit
  48     val reset_cache : unit -> unit
  49     val print_stat : unit -> unit
  50
  51     val dest_float_interval : term -> term * term * term
  52     val mk_float_interval_small_num : int -> thm
  53     val mk_float_interval_num : num -> thm
  54
  55     val float_lo : int -> term -> thm
  56     val float_hi : int -> term -> thm
  57
  58     val float_interval_round : int -> thm -> thm
  59
  60     val float_interval_neg : thm -> thm
  61     val float_interval_mul : int -> thm -> thm -> thm
  62     val float_interval_div : int -> thm -> thm -> thm
  63     val float_interval_add : int -> thm -> thm -> thm
  64     val float_interval_sub : int -> thm -> thm -> thm
  65     val float_interval_sqrt : int -> thm -> thm
  66
  67     val float_abs : term -> thm
  68
  69     val FLOAT_TO_NUM_CONV : term -> thm
  70 end;;
  71
  72
  73 module Arith_float : Arith_float_sig = struct
  74
  75 open Big_int;;
  76 open Arith_misc;;
  77 open Arith_nat;;
  78 open Num_exp_theory;;
  79 open Float_theory;;
  80 open Interval_arith;;
  81 open Misc_vars;;
  82
  83 (* interval *)
  84
  85 let APPROX_INTERVAL' = (UNDISCH_ALL o REWRITE_RULE[GSYM IMP_IMP] o SPEC_ALL) APPROX_INTERVAL;;
  86
  87 let interval_const = `interval_arith` and
  88     num_exp_const = `num_exp`;;
  89
  90 let b0_const = (fst o dest_comb o lhand o concl) (Arith_hash.def_array.(0));;
  91 let b0_name = (fst o dest_const) b0_const;;
  92 let base_const = mk_small_numeral Arith_hash.arith_base;;
  93
  94 let NUM_REMOVE = prove(mk_eq(mk_comb(Arith_hash.num_const, n_var_num), n_var_num),
  95                        REWRITE_TAC[Arith_hash.num_def; NUMERAL]);;
  96
  97 (* B0 n = base * n *)
  98 let b0_thm = prove(mk_eq(mk_comb(b0_const, n_var_num),
  99                          mk_binop mul_op_num base_const n_var_num),
 100                    REWRITE_TAC[Arith_hash.def_array.(0)] THEN
 101                      TRY ARITH_TAC THEN
 102                      GEN_REWRITE_TAC (LAND_CONV o DEPTH_CONV) [BIT0] THEN
 103                      ARITH_TAC);;
 104
 105 let dest_num_exp tm =
 106   let ltm, e_tm = dest_comb tm in
 107     rand ltm, e_tm;;
 108
 109 let num_exp_const = `num_exp`;;
 110 let mk_num_exp n_tm e_tm = mk_binop num_exp_const n_tm e_tm;;
 111
 112 (* float_num s n e -> "s", n, e *)
 113 let dest_float tm =
 114   let ltm, e_tm = dest_comb tm in
 115   let ltm, n_tm = dest_comb ltm in
 116   let float_tm, s_tm = dest_comb ltm in
 117     if (fst o dest_const) float_tm <> "float_num" then
 118       failwith "dest_float: not float"
 119     else
 120       (fst o dest_const) s_tm, n_tm, e_tm;;
 121
 122
 123 (************************************)
 124
 125 let NUM_EXP_EXP' = SPEC_ALL NUM_EXP_EXP;;
 126 let NUM_EXP_0' = (SPEC_ALL o REWRITE_RULE[NUMERAL]) NUM_EXP_0;;
 127 let NUM_EXP_LE' = (UNDISCH_ALL o SPEC_ALL) NUM_EXP_LE;;
 128 let NUM_EXP_LT' = (UNDISCH_ALL o SPEC_ALL) NUM_EXP_LT;;
 129
 130 (* B0 n = num_exp n bits *)
 131 let normal_lemma1 = prove(mk_eq(mk_comb(b0_const, n_var_num), `num_exp n 1`),
 132    REWRITE_TAC[Arith_hash.def_array.(0); num_exp] THEN
 133      TRY ARITH_TAC THEN
 134      GEN_REWRITE_TAC (LAND_CONV o DEPTH_CONV) [BIT0] THEN
 135      ARITH_TAC);;
 136
 137 let NORMAL_LEMMA1 = NUMERALS_TO_NUM normal_lemma1;;
 138
 139 let normal_lemma2 = prove(mk_eq (mk_comb (b0_const, `num_exp n e`), `num_exp n (SUC e)`),
 140    REWRITE_TAC[normal_lemma1; NUM_EXP_EXP] THEN ARITH_TAC);;
 141
 142 let rec normalize tm =
 143   if (is_comb tm) then
 144     let ltm, rtm = dest_comb tm in
 145     let lname = (fst o dest_const) ltm in
 146       if (lname = b0_name) then
 147         let lth = INST[rtm, n_var_num] NORMAL_LEMMA1 in
 148         let rth, flag = normalize rtm in
 149           if flag then
 150             let ltm, lexp = (dest_comb o snd o dest_eq o concl) lth in
 151             let ltm, rtm = dest_comb ltm in
 152             let rn, rexp = (dest_comb o snd o dest_eq o concl) rth in
 153             let rn = rand rn in
 154             let th1 = AP_THM (AP_TERM ltm rth) lexp in
 155             let th2 = INST[rexp, e1_var_num; lexp, e2_var_num; rn, n_var_num] NUM_EXP_EXP' in
 156             let th3 = TRANS lth (TRANS th1 th2) in
 157             let ltm, rtm = (dest_comb o snd o dest_eq o concl) th3 in
 158             let add_th = raw_add_conv_hash rtm in
 159             let th4 = AP_TERM ltm add_th in
 160               (TRANS th3 th4, true)
 161           else
 162             (lth, true)
 163       else
 164         (REFL tm, false)
 165   else
 166     (REFL tm, false);;
 167
 168 (* Converts a raw numeral to a num_exp expression *)
 169 let to_num_exp tm =
 170   let x, flag = normalize tm in
 171     if flag then x
 172     else
 173       INST[tm, n_var_num] NUM_EXP_0';;
 174
 175 (************************************)
 176
 177 let SYM_NUM_EXP_0' = SYM NUM_EXP_0';;
 178 let NUM_EXP_n0 = prove(`!e. num_exp 0 e = 0`, REWRITE_TAC[num_exp; MULT_CLAUSES]);;
 179 let NUM_EXP_n0' = (REWRITE_RULE[NUMERAL] o SPEC_ALL) NUM_EXP_n0;;
 180
 181 let NUM_EXP_DENORM = (UNDISCH_ALL o prove)
 182   (mk_imp(`e = _0 <=> F`, mk_eq(`num_exp n e`, mk_comb (b0_const, `num_exp n (PRE e)`))),
 183    REWRITE_TAC[] THEN ONCE_REWRITE_TAC[SYM (REWRITE_CONV[NUMERAL] `0`)] THEN
 184      REWRITE_TAC[num_exp; b0_thm] THEN
 185      REWRITE_TAC[ARITH_RULE (mk_eq(mk_binop mul_op_num base_const `n * a:num`,
 186                                    mk_binop mul_op_num `n:num` (mk_binop mul_op_num base_const `a:num`)))] THEN
 187      REWRITE_TAC[GSYM EXP] THEN
 188      SIMP_TAC[ARITH_RULE `~(e = 0) ==> SUC (PRE e) = e`]);;
 189
 190 (* Converts num_exp n e to a numeral by adding e B0's *)
 191 let rec denormalize tm =
 192   let ltm, etm = dest_comb tm in
 193   let ntm = rand ltm in
 194     if (etm = zero_const) then
 195       INST[ntm, n_var_num] SYM_NUM_EXP_0'
 196     else
 197       if ntm = zero_const then
 198         INST[etm, e_var_num] NUM_EXP_n0'
 199       else
 200         let e_th = raw_eq0_hash_conv etm in
 201         let th0' = INST[etm, e_var_num; ntm, n_var_num] NUM_EXP_DENORM in
 202         let th0 = MY_PROVE_HYP e_th th0' in
 203         let b0_tm, rtm = dest_comb(rand(concl th0)) in
 204         let ltm, pre_tm = dest_comb rtm in
 205         let pre_th = raw_pre_hash_conv pre_tm in
 206         let th1 = AP_TERM ltm pre_th in
 207         let th2 = denormalize (rand(concl th1)) in
 208           TRANS th0 (AP_TERM b0_tm (TRANS th1 th2));;
 209
 210 (***************************************)
 211
 212 let rec comb_number tm n =
 213   if (is_comb tm) then comb_number ((snd o dest_comb) tm) (n + 1) else n;;
 214
 215 let make_lo_thm i =
 216   let th_concl = mk_binop `(<=):num->num->bool`
 217     (mk_comb (Arith_hash.const_array.(0), n_var_num))
 218     (mk_comb (Arith_hash.const_array.(i), n_var_num)) in
 219     prove(th_concl,
 220           REWRITE_TAC[Arith_hash.def_array.(i); Arith_hash.def_array.(0)] THEN
 221             REWRITE_TAC[ARITH_LE; LE_REFL] THEN
 222             ARITH_TAC);;
 223
 224 let lo_thm_array = Array.init Arith_hash.arith_base make_lo_thm;;
 225 let lo_thm_table = Hashtbl.create Arith_hash.arith_base;;
 226
 227 for i = 0 to Arith_hash.arith_base - 1 do
 228   Hashtbl.add lo_thm_table Arith_hash.const_array.(i) lo_thm_array.(i);
 229 done;;
 230
 231 let make_lo_thm2 i =
 232   let th_concl = mk_imp (`n <= m:num`,
 233                          mk_binop `(<=):num->num->bool`
 234                           (mk_comb (Arith_hash.const_array.(0), n_var_num))
 235                           (mk_comb (Arith_hash.const_array.(i), m_var_num))) in
 236     (UNDISCH_ALL o prove) (th_concl,
 237           REWRITE_TAC[Arith_hash.def_array.(i); Arith_hash.def_array.(0); ARITH_LE] THEN
 238                           ARITH_TAC);;
 239
 240 let lo_thm2_array = Array.init Arith_hash.arith_base make_lo_thm2;;
 241 let lo_thm2_table = Hashtbl.create Arith_hash.arith_base;;
 242
 243 for i = 0 to Arith_hash.arith_base - 1 do
 244   Hashtbl.add lo_thm2_table Arith_hash.const_array.(i) lo_thm2_array.(i);
 245 done;;
 246
 247
 248 let make_hi_thm i =
 249   let th_concl = mk_imp (`n < m:num`,
 250                          mk_binop `(<):num->num->bool`
 251                           (mk_comb (Arith_hash.const_array.(i), n_var_num))
 252                           (mk_comb (Arith_hash.const_array.(0), m_var_num))) in
 253     (UNDISCH_ALL o prove) (th_concl,
 254                            REWRITE_TAC[Arith_hash.def_array.(i); Arith_hash.def_array.(0); ARITH_LT] THEN
 255                              ARITH_TAC);;
 256
 257 let hi_thm_array = Array.init Arith_hash.arith_base make_hi_thm;;
 258 let hi_thm_table = Hashtbl.create Arith_hash.arith_base;;
 259
 260 for i = 0 to Arith_hash.arith_base - 1 do
 261   Hashtbl.add hi_thm_table Arith_hash.const_array.(i) hi_thm_array.(i);
 262 done;;
 263
 264 (***************************************)
 265
 266 let LE_REFL' = SPEC_ALL LE_REFL;;
 267 let LE_TRANS' = (UNDISCH_ALL o SPEC_ALL o REWRITE_RULE[GSYM IMP_IMP]) LE_TRANS;;
 268
 269 let lo_num_conv p tm =
 270   let n = comb_number tm 0 in
 271     if (n <= p) then
 272       INST[tm, n_var_num] LE_REFL'
 273     else
 274       let rec lo_bound n tm =
 275         let btm, rtm = dest_comb tm in
 276         let th0 = INST[rtm, n_var_num] (Hashtbl.find lo_thm_table btm) in
 277           if n > 1 then
 278             let rth = lo_bound (n - 1) rtm in
 279             let xtm = rand (rator (concl rth)) in
 280             let th1' = INST[xtm, n_var_num; rtm, m_var_num] (Hashtbl.find lo_thm2_table btm) in
 281             let th1 = MY_PROVE_HYP rth th1' in
 282               th1
 283           else
 284             th0 in
 285
 286         lo_bound (n - p) tm;;
 287
 288 let N_LT_SUC = ARITH_RULE `n < SUC n`;;
 289 let LT_IMP_LE' = (UNDISCH_ALL o SPEC_ALL) LT_IMP_LE;;
 290 let N_LT_SUC = ARITH_RULE `n < SUC n`;;
 291 let LT_LE_TRANS = (UNDISCH_ALL o ARITH_RULE) `n < e ==> e <= m ==> n < m:num`;;
 292
 293 (* Generates a theorem |- n <= m such that m contains at most p non-zero digits *)
 294 let hi_num_conv p tm =
 295   let n = comb_number tm 0 in
 296     if (n <= p) then
 297       INST[tm, n_var_num] LE_REFL'
 298     else
 299       let k = n - p in
 300
 301       let rec check_b0s n tm =
 302         let btm, rtm = dest_comb tm in
 303           if ((fst o dest_const) btm = b0_name) then
 304             if n > 1 then check_b0s (n - 1) rtm else true
 305           else
 306             false in
 307
 308         if (check_b0s k tm) then
 309           INST[tm, n_var_num] LE_REFL'
 310         else
 311           let rec hi_bound n tm =
 312             if n > 0 then
 313               let btm, rtm = dest_comb tm in
 314               let r_th = hi_bound (n - 1) rtm in
 315               let xtm = rand (concl r_th) in
 316               let th0 = INST[rtm, n_var_num; xtm, m_var_num] (Hashtbl.find hi_thm_table btm) in
 317                 MY_PROVE_HYP r_th th0
 318             else
 319               let th0 = INST[tm, n_var_num] N_LT_SUC in
 320               let ltm, suc_tm = dest_comb (concl th0) in
 321               let suc_th = raw_suc_conv_hash suc_tm in
 322                 EQ_MP (AP_TERM ltm suc_th) th0 in
 323
 324           let th = hi_bound k tm in
 325           let m_tm, l_tm = dest_comb (concl th) in
 326             MY_PROVE_HYP th (INST[rand m_tm, m_var_num; l_tm, n_var_num] LT_IMP_LE');;
 327
 328 (* Generates a theorem |- n < m such that m contains at most p non-zero digits *)
 329 let hi_lt_num_conv p tm =
 330   let n = comb_number tm 0 in
 331     if (n <= p) then
 332       let th0 = INST[tm, n_var_num] N_LT_SUC in
 333       let ltm, rtm = dest_comb(concl th0) in
 334       let suc_th = raw_suc_conv_hash rtm in
 335         EQ_MP (AP_TERM ltm suc_th) th0
 336     else
 337       let k = n - p in
 338
 339       let rec check_b0s n tm =
 340         let btm, rtm = dest_comb tm in
 341           if ((fst o dest_const) btm = b0_name) then
 342             if n > 1 then check_b0s (n - 1) rtm else true
 343           else
 344             false in
 345
 346         if (check_b0s k tm) then
 347           let th0 = INST[tm, n_var_num] N_LT_SUC in
 348           let ltm, rtm = dest_comb (concl th0) in
 349           let suc_th = raw_suc_conv_hash rtm in
 350           let suc_tm = rand(concl suc_th) in
 351           let th1 = hi_num_conv p suc_tm in
 352           let th2 = EQ_MP (AP_TERM ltm suc_th) th0 in
 353           let th = INST[tm, n_var_num; suc_tm, e_var_num; rand(concl th1), m_var_num] LT_LE_TRANS in
 354             MY_PROVE_HYP th1 (MY_PROVE_HYP th2 th)
 355
 356         else
 357           let rec hi_bound n tm =
 358             if n > 0 then
 359               let btm, rtm = dest_comb tm in
 360               let r_th = hi_bound (n - 1) rtm in
 361               let xtm = rand (concl r_th) in
 362               let th0 = INST[rtm, n_var_num; xtm, m_var_num] (Hashtbl.find hi_thm_table btm) in
 363                 MY_PROVE_HYP r_th th0
 364             else
 365               let th0 = INST[tm, n_var_num] N_LT_SUC in
 366               let ltm, suc_tm = dest_comb (concl th0) in
 367               let suc_th = raw_suc_conv_hash suc_tm in
 368                 EQ_MP (AP_TERM ltm suc_th) th0 in
 369             hi_bound k tm;;
 370
 371 (*****************************************)
 372
 373 let num_exp_lo p tm =
 374   let ltm, e_tm = dest_comb tm in
 375   let n_tm = rand ltm in
 376   let n_th = lo_num_conv p n_tm in
 377   let m_tm = rand (rator (concl n_th)) in
 378   let m_norm, flag = normalize m_tm in
 379
 380   let th0' = INST[m_tm, m_var_num; n_tm, n_var_num; e_tm, e_var_num] NUM_EXP_LE' in
 381   let th0 = MY_PROVE_HYP n_th th0' in
 382     if flag then
 383       let th1 = AP_THM (AP_TERM (rator ltm) m_norm) e_tm in
 384       let m_tm, me_tm = (dest_comb o rand o concl) m_norm in
 385       let th2 = INST[me_tm, e1_var_num; e_tm, e2_var_num; rand m_tm, n_var_num] NUM_EXP_EXP' in
 386       let th3 = TRANS th1 th2 in
 387       let ltm, rtm = (dest_comb o rand o concl) th3 in
 388       let th_add = raw_add_conv_hash rtm in
 389       let th4 = TRANS th3 (AP_TERM ltm th_add) in
 390         EQ_MP (AP_THM (AP_TERM le_op_num th4) tm) th0
 391     else
 392       th0;;
 393
 394 let num_exp_hi p tm =
 395   let ltm, e_tm = dest_comb tm in
 396   let n_tm = rand ltm in
 397   let n_th = hi_num_conv p n_tm in
 398   let m_tm = rand (concl n_th) in
 399   let m_norm, flag = normalize m_tm in
 400
 401   let th0' = INST[m_tm, n_var_num; n_tm, m_var_num; e_tm, e_var_num] NUM_EXP_LE' in
 402   let th0 = MY_PROVE_HYP n_th th0' in
 403     if flag then
 404       let th1 = AP_THM (AP_TERM (rator ltm) m_norm) e_tm in
 405       let m_tm, me_tm = (dest_comb o rand o concl) m_norm in
 406       let th2 = INST[me_tm, e1_var_num; e_tm, e2_var_num; rand m_tm, n_var_num] NUM_EXP_EXP' in
 407       let th3 = TRANS th1 th2 in
 408       let ltm, rtm = (dest_comb o rand o concl) th3 in
 409       let th_add = raw_add_conv_hash rtm in
 410       let th4 = TRANS th3 (AP_TERM ltm th_add) in
 411         EQ_MP (AP_TERM (rator (concl th0)) th4) th0
 412     else
 413       th0;;
 414
 415 let num_exp_hi_lt p tm =
 416   let ltm, e_tm = dest_comb tm in
 417   let n_tm = rand ltm in
 418   let n_th = hi_lt_num_conv p n_tm in
 419   let m_tm = rand (concl n_th) in
 420   let m_norm, flag = normalize m_tm in
 421
 422   let th0' = INST[m_tm, n_var_num; n_tm, m_var_num; e_tm, e_var_num] NUM_EXP_LT' in
 423   let th0 = MY_PROVE_HYP n_th th0' in
 424     if flag then
 425       let th1 = AP_THM (AP_TERM (rator ltm) m_norm) e_tm in
 426       let m_tm, me_tm = (dest_comb o rand o concl) m_norm in
 427       let th2 = INST[me_tm, e1_var_num; e_tm, e2_var_num; rand m_tm, n_var_num] NUM_EXP_EXP' in
 428       let th3 = TRANS th1 th2 in
 429       let ltm, rtm = (dest_comb o rand o concl) th3 in
 430       let th_add = raw_add_conv_hash rtm in
 431       let th4 = TRANS th3 (AP_TERM ltm th_add) in
 432         EQ_MP (AP_TERM (rator (concl th0)) th4) th0
 433     else
 434       th0;;
 435
 436 (***************************************)
 437 (* num_exp_lt, num_exp_le *)
 438
 439 let transform = UNDISCH_ALL o REWRITE_RULE[GSYM IMP_IMP] o SPEC_ALL;;
 440
 441 let NUM_EXP_LT1_EQ' = transform NUM_EXP_LT1_EQ and
 442     NUM_EXP_LT2_EQ' = transform NUM_EXP_LT2_EQ;;
 443
 444 let NUM_EXP_LE1_EQ' = transform NUM_EXP_LE1_EQ and
 445     NUM_EXP_LE2_EQ' = transform NUM_EXP_LE2_EQ;;
 446
 447 let num_exp_lt tm1 tm2 =
 448   let n1_tm, e1_tm = dest_num_exp tm1 in
 449   let n2_tm, e2_tm = dest_num_exp tm2 in
 450   let sub_th, le_th = raw_sub_and_le_hash_conv e1_tm e2_tm in
 451   let r_tm = rand(concl sub_th) in
 452     if (rand(concl le_th) = e1_tm) then
 453       let x_expr = mk_num_exp n1_tm r_tm in
 454       let x_th = denormalize x_expr in
 455       let x_tm = rand(concl x_th) in
 456
 457       let th0 = INST[e2_tm, e2_var_num; e1_tm, e1_var_num;
 458                      r_tm, r_var_num; x_tm, x_var_num;
 459                      n1_tm, n1_var_num; n2_tm, n2_var_num] NUM_EXP_LT1_EQ' in
 460       let th1 = MY_PROVE_HYP x_th (MY_PROVE_HYP sub_th (MY_PROVE_HYP le_th th0)) in
 461       let lt_th = raw_lt_hash_conv (rand(concl th1)) in
 462         TRANS th1 lt_th
 463     else
 464       let x_expr = mk_num_exp n2_tm r_tm in
 465       let x_th = denormalize x_expr in
 466       let x_tm = rand(concl x_th) in
 467
 468       let th0 = INST[e2_tm, e2_var_num; e1_tm, e1_var_num;
 469                      r_tm, r_var_num; x_tm, x_var_num;
 470                      n1_tm, n1_var_num; n2_tm, n2_var_num] NUM_EXP_LT2_EQ' in
 471       let th1 = MY_PROVE_HYP x_th (MY_PROVE_HYP sub_th (MY_PROVE_HYP le_th th0)) in
 472       let lt_th = raw_lt_hash_conv (rand(concl th1)) in
 473         TRANS th1 lt_th;;
 474
 475
 476 let num_exp_le tm1 tm2 =
 477   let n1_tm, e1_tm = dest_num_exp tm1 in
 478   let n2_tm, e2_tm = dest_num_exp tm2 in
 479   let sub_th, le_th = raw_sub_and_le_hash_conv e1_tm e2_tm in
 480   let r_tm = rand(concl sub_th) in
 481     if (rand(concl le_th) = e1_tm) then
 482       let x_expr = mk_num_exp n1_tm r_tm in
 483       let x_th = denormalize x_expr in
 484       let x_tm = rand(concl x_th) in
 485
 486       let th0 = INST[e2_tm, e2_var_num; e1_tm, e1_var_num;
 487                      r_tm, r_var_num; x_tm, x_var_num;
 488                      n1_tm, n1_var_num; n2_tm, n2_var_num] NUM_EXP_LE1_EQ' in
 489       let th1 = MY_PROVE_HYP x_th (MY_PROVE_HYP sub_th (MY_PROVE_HYP le_th th0)) in
 490       let le_th = raw_le_hash_conv (rand(concl th1)) in
 491         TRANS th1 le_th
 492     else
 493       let x_expr = mk_num_exp n2_tm r_tm in
 494       let x_th = denormalize x_expr in
 495       let x_tm = rand(concl x_th) in
 496
 497       let th0 = INST[e2_tm, e2_var_num; e1_tm, e1_var_num;
 498                      r_tm, r_var_num; x_tm, x_var_num;
 499                      n1_tm, n1_var_num; n2_tm, n2_var_num] NUM_EXP_LE2_EQ' in
 500       let th1 = MY_PROVE_HYP x_th (MY_PROVE_HYP sub_th (MY_PROVE_HYP le_th th0)) in
 501       let le_th = raw_le_hash_conv (rand(concl th1)) in
 502         TRANS th1 le_th;;
 503
 504 (***************************************)
 505 (* num_exp_mul *)
 506
 507 let NUM_EXP_MUL' = SPEC_ALL NUM_EXP_MUL;;
 508
 509 let num_exp_mul tm1 tm2 =
 510   let n1_tm, e1_tm = dest_comb tm1 in
 511   let n1_tm = rand n1_tm in
 512   let n2_tm, e2_tm = dest_comb tm2 in
 513   let n2_tm = rand n2_tm in
 514   let th0 = INST[n1_tm, n1_var_num; e1_tm, e1_var_num;
 515                  n2_tm, n2_var_num; e2_tm, e2_var_num] NUM_EXP_MUL' in
 516   let ltm, tm_add = dest_comb (rand (concl th0)) in
 517   let tm_mul = rand ltm in
 518   let th_mul = raw_mul_conv_hash tm_mul in
 519   let th_add = raw_add_conv_hash tm_add in
 520     TRANS th0 (MK_COMB (AP_TERM (rator ltm) th_mul, th_add));;
 521
 522
 523 (**********************************)
 524 (* num_exp_add *)
 525
 526 let NUM_EXP_ADD' = (UNDISCH_ALL o SPEC_ALL) NUM_EXP_ADD;;
 527 let ADD_COMM = ARITH_RULE `m + n = n + m:num`;;
 528
 529 let num_exp_add tm1 tm2 =
 530   let n1_tm, e1_tm = dest_comb tm1 in
 531   let n1_tm = rand n1_tm in
 532   let n2_tm, e2_tm = dest_comb tm2 in
 533   let n2_tm = rand n2_tm in
 534   let e_sub, e_le = raw_sub_and_le_hash_conv e1_tm e2_tm in
 535
 536   let flag = (rand(concl e_le) = e2_tm) in
 537
 538   let th0' =
 539     if flag then
 540       INST[n1_tm, n1_var_num; e1_tm, e1_var_num;
 541            n2_tm, n2_var_num; e2_tm, e2_var_num] NUM_EXP_ADD'
 542     else
 543       INST[n2_tm, n1_var_num; e2_tm, e1_var_num;
 544            n1_tm, n2_var_num; e1_tm, e2_var_num] NUM_EXP_ADD' in
 545
 546   let th0 = MY_PROVE_HYP e_le th0' in
 547   let ltm, e0_tm = dest_comb(rand(concl th0)) in
 548   let exp_tm, add_tm = dest_comb ltm in
 549   let ltm, d_tm = dest_comb add_tm in
 550   let th1 = AP_TERM (rator d_tm) e_sub in
 551   let th2 = denormalize (rand(concl th1)) in
 552   let th3 = AP_TERM ltm (TRANS th1 th2) in
 553   let th4 = raw_add_conv_hash (rand(concl th3)) in
 554   let th5 = AP_THM (AP_TERM exp_tm (TRANS th3 th4)) e0_tm in
 555   let th = TRANS th0 th5 in
 556     if flag then th else
 557       TRANS (INST[tm1, m_var_num; tm2, n_var_num] ADD_COMM) th;;
 558
 559 (****************************************)
 560 (* num_exp_sub *)
 561
 562 let NUM_EXP_SUB1' = (UNDISCH_ALL o REWRITE_RULE[GSYM IMP_IMP] o SPEC_ALL) NUM_EXP_SUB1 and
 563     NUM_EXP_SUB2' = (UNDISCH_ALL o REWRITE_RULE[GSYM IMP_IMP] o SPEC_ALL) NUM_EXP_SUB2 and
 564     NUM_EXP_LE1' = (UNDISCH_ALL o REWRITE_RULE[GSYM IMP_IMP] o SPEC_ALL) NUM_EXP_LE1 and
 565     NUM_EXP_LE2' = (UNDISCH_ALL o REWRITE_RULE[GSYM IMP_IMP] o SPEC_ALL) NUM_EXP_LE2;;
 566
 567 (* Returns two theorems: |- tm1 - tm2 = tm, |- tm2 <= tm1 or
 568    |- tm2 - tm1 = tm, |- tm1 <= tm2 *)
 569 let num_exp_sub tm1 tm2 =
 570   let n1_tm, e1_tm = dest_num_exp tm1 in
 571   let n2_tm, e2_tm = dest_num_exp tm2 in
 572   let e_sub, e_le = raw_sub_and_le_hash_conv e1_tm e2_tm in
 573
 574     if rand(concl e_le) = e1_tm then
 575       (* e2 <= e1 *)
 576       let e1_sub_e2 = rand(concl e_sub) in
 577       let a0 = mk_num_exp n1_tm e1_sub_e2 in
 578       let b = n2_tm in
 579       let a_th = denormalize a0 in
 580       let a = rand(concl a_th) in
 581
 582       let th_sub, th_le = raw_sub_and_le_hash_conv a b in
 583         if rand(concl th_le) = a then
 584           (* b <= a *)
 585           let a_sub_b = TRANS (AP_THM (AP_TERM sub_op_num a_th) b) th_sub in
 586           let b_le_a = EQ_MP (SYM (AP_TERM (rator(concl th_le)) a_th)) th_le in
 587           let th0 = AP_THM (AP_TERM num_exp_const a_sub_b) e2_tm in
 588
 589           let inst = INST[n1_tm, n1_var_num; e1_tm, e1_var_num;
 590                           n2_tm, n2_var_num; e2_tm, e2_var_num; e1_sub_e2, r_var_num] in
 591           let th1_sub = inst NUM_EXP_SUB1' in
 592           let th1_le = inst NUM_EXP_LE1' in
 593           let th2_sub = MY_PROVE_HYP e_sub (MY_PROVE_HYP e_le th1_sub) in
 594           let th2_le = MY_PROVE_HYP e_sub (MY_PROVE_HYP b_le_a (MY_PROVE_HYP e_le th1_le)) in
 595             TRANS th2_sub th0, th2_le
 596
 597         else
 598           (* a <= b *)
 599           let b_sub_a = TRANS (AP_TERM (rator(lhand(concl th_sub))) a_th) th_sub in
 600           let a_le_b = EQ_MP (SYM (AP_THM (AP_TERM le_op_num a_th) b)) th_le in
 601           let th0 = AP_THM (AP_TERM num_exp_const b_sub_a) e2_tm in
 602           let inst = INST[n2_tm, n1_var_num; e2_tm, e1_var_num;
 603                           n1_tm, n2_var_num; e1_tm, e2_var_num; e1_sub_e2, r_var_num] in
 604           let th1_sub = inst NUM_EXP_SUB2' in
 605           let th1_le = inst NUM_EXP_LE2' in
 606           let th2_sub = MY_PROVE_HYP e_sub (MY_PROVE_HYP e_le th1_sub) in
 607           let th2_le = MY_PROVE_HYP e_sub (MY_PROVE_HYP a_le_b (MY_PROVE_HYP e_le th1_le)) in
 608             TRANS th2_sub th0, th2_le
 609
 610     else
 611       (* e1 <= e2 *)
 612       let e2_sub_e1 = rand(concl e_sub) in
 613       let b0 = mk_num_exp n2_tm e2_sub_e1 in
 614       let a = n1_tm in
 615       let b_th = denormalize b0 in
 616       let b = rand(concl b_th) in
 617
 618       let th_sub, th_le = raw_sub_and_le_hash_conv a b in
 619         if rand(concl th_le) = a then
 620           (* b <= a *)
 621           let a_sub_b = TRANS (AP_TERM (rator(lhand(concl th_sub))) b_th) th_sub in
 622           let b_le_a = EQ_MP (SYM (AP_THM (AP_TERM le_op_num b_th) a)) th_le in
 623           let th0 = AP_THM (AP_TERM num_exp_const a_sub_b) e1_tm in
 624           let inst = INST[n1_tm, n1_var_num; e1_tm, e1_var_num;
 625                           n2_tm, n2_var_num; e2_tm, e2_var_num; e2_sub_e1, r_var_num] in
 626           let th1_sub = inst NUM_EXP_SUB2' in
 627           let th1_le = inst NUM_EXP_LE2' in
 628           let th2_sub = MY_PROVE_HYP e_sub (MY_PROVE_HYP e_le th1_sub) in
 629           let th2_le = MY_PROVE_HYP e_sub (MY_PROVE_HYP b_le_a (MY_PROVE_HYP e_le th1_le)) in
 630             TRANS th2_sub th0, th2_le
 631
 632         else
 633           (* a <= b *)
 634           let b_sub_a = TRANS (AP_THM (AP_TERM sub_op_num b_th) a) th_sub in
 635           let a_le_b = EQ_MP (SYM (AP_TERM (rator(concl th_le)) b_th)) th_le in
 636           let th0 = AP_THM (AP_TERM num_exp_const b_sub_a) e1_tm in
 637           let inst = INST[n2_tm, n1_var_num; e2_tm, e1_var_num;
 638                           n1_tm, n2_var_num; e1_tm, e2_var_num; e2_sub_e1, r_var_num] in
 639           let th1_sub = inst NUM_EXP_SUB1' in
 640           let th1_le = inst NUM_EXP_LE1' in
 641           let th2_sub = MY_PROVE_HYP e_sub (MY_PROVE_HYP e_le th1_sub) in
 642           let th2_le = MY_PROVE_HYP e_sub (MY_PROVE_HYP a_le_b (MY_PROVE_HYP e_le th1_le)) in
 643             TRANS th2_sub th0, th2_le;;
 644
 645
 646 (*************************************)
 647 (* division *)
 648
 649 let NUM_EXP_DIV1' = (UNDISCH_ALL o PURE_REWRITE_RULE[NUMERAL] o
 650                        PURE_ONCE_REWRITE_RULE[ARITH_RULE `~(x = 0) <=> (x = 0 <=> F)`] o
 651                        REWRITE_RULE[GSYM IMP_IMP]) NUM_EXP_DIV1;;
 652 let NUM_EXP_DIV2' = (UNDISCH_ALL o PURE_REWRITE_RULE[NUMERAL] o
 653                        PURE_ONCE_REWRITE_RULE[ARITH_RULE `~(x = 0) <=> (x = 0 <=> F)`] o
 654                        REWRITE_RULE[GSYM IMP_IMP]) NUM_EXP_DIV2;;
 655
 656 let num_exp_div tm1 tm2 =
 657   let n1_tm, e1_tm = dest_comb tm1 in
 658   let n1_tm = rand n1_tm in
 659   let n2_tm, e2_tm = dest_comb tm2 in
 660   let n2_tm = rand n2_tm in
 661   let e_sub, e_le = raw_sub_and_le_hash_conv e1_tm e2_tm in
 662
 663   let inst = INST[n1_tm, n1_var_num; e1_tm, e1_var_num;
 664                   n2_tm, n2_var_num; e2_tm, e2_var_num] in
 665
 666   let n2_not_0 = raw_eq0_hash_conv n2_tm in
 667     if ((fst o dest_const o rand o concl) n2_not_0 = "T") then
 668       failwith "num_exp_div: n2 = 0"
 669     else
 670       if (rand(concl e_le) = e1_tm) then
 671          let th0' = inst NUM_EXP_DIV1' in
 672         let th0 = MY_PROVE_HYP n2_not_0 (MY_PROVE_HYP e_le th0') in
 673
 674         let ltm, rtm = dest_comb(rand(concl th0)) in
 675         let div_tm, rtm2 = dest_comb ltm in
 676         let num_exp_tm = rator rtm2 in
 677
 678         let th1 = AP_THM (AP_TERM div_tm (AP_TERM num_exp_tm e_sub)) rtm in
 679         let ltm, rtm = dest_comb(rand(concl th1)) in
 680         let tm1 = rand ltm in
 681
 682         let th2 = AP_THM (AP_TERM div_tm (denormalize tm1)) rtm in
 683         let th3 = raw_div_hash_conv (rand(concl th2)) in
 684         let th = TRANS th0 (TRANS th1 (TRANS th2 th3)) in
 685           TRANS th (INST[rand(concl th), n_var_num] NUM_EXP_0')
 686
 687       else
 688         let th0' = inst NUM_EXP_DIV2' in
 689         let th0 = MY_PROVE_HYP n2_not_0 (MY_PROVE_HYP e_le th0') in
 690
 691         let ltm, rtm = dest_comb(rand(concl th0)) in
 692         let num_exp_tm = rator rtm in
 693         let th1 = AP_TERM ltm (AP_TERM num_exp_tm e_sub) in
 694
 695         let ltm, rtm = dest_comb(rand(concl th1)) in
 696         let th2 = AP_TERM ltm (denormalize rtm) in
 697         let th3 = raw_div_hash_conv (rand(concl th2)) in
 698         let th = TRANS th0 (TRANS th1 (TRANS th2 th3)) in
 699           TRANS th (INST[rand(concl th), n_var_num] NUM_EXP_0');;
 700
 701
 702 (*****************************)
 703 (* Computes a lower bound for (op tm1 tm2) with p significant digits *)
 704 let num_exp_op_lo p op tm1 tm2 =
 705   let op_th = op tm1 tm2 in
 706   let rtm = rand (concl op_th) in
 707   let lo_th = num_exp_lo p rtm in
 708   let ltm = rator (concl lo_th) in
 709   let th0 = AP_TERM ltm op_th in
 710     EQ_MP (SYM th0) lo_th;;
 711
 712 (* Computes an upper bound for (op tm1 tm2) with p significant digits *)
 713 let num_exp_op_hi p op tm1 tm2 =
 714   let op_th = op tm1 tm2 in
 715   let rtm = rand (concl op_th) in
 716   let hi_th = num_exp_hi p rtm in
 717   let tm = rand (concl hi_th) in
 718   let th0 = AP_THM (AP_TERM le_op_num op_th) tm in
 719     EQ_MP (SYM th0) hi_th;;
 720
 721 (* Computes a strict upper bound for (op tm1 tm2) with p significant digits *)
 722 let num_exp_op_hi_lt p op tm1 tm2 =
 723   let op_th = op tm1 tm2 in
 724   let rtm = rand (concl op_th) in
 725   let hi_lt_th = num_exp_hi_lt p rtm in
 726   let tm = rand (concl hi_lt_th) in
 727   let th0 = AP_THM (AP_TERM lt_op_num op_th) tm in
 728     EQ_MP (SYM th0) hi_lt_th;;
 729
 730 (******************************************)
 731 (* float *)
 732
 733 let mod_plus = new_definition `mod_plus s1 s2 = (~(s1 /\ s2) /\ (s1 \/ s2))`;;
 734
 735
 736 (********************)
 737 (* Float operations *)
 738 (********************)
 739
 740 module Float_ops = struct
 741
 742
 743 (**********************************)
 744
 745 (* FLOAT_LT *)
 746
 747 let FLOAT_LT_FF = prove(`float_num F n1 e1 < float_num F n2 e2 <=> num_exp n1 e1 < num_exp n2 e2`,
 748                         REWRITE_TAC[float; GSYM REAL_OF_NUM_LT; REAL_MUL_LID; real_div] THEN
 749                           MATCH_MP_TAC REAL_LT_RMUL_EQ THEN
 750                           MATCH_MP_TAC REAL_LT_INV THEN
 751                           REWRITE_TAC[REAL_OF_NUM_LT; LT_NZ; NUM_EXP_EQ_0] THEN
 752                           ARITH_TAC);;
 753
 754 let FLOAT_LT_TT = prove(`float_num T n1 e1 < float_num T n2 e2 <=> num_exp n2 e2 < num_exp n1 e1`,
 755                         REWRITE_TAC[FLOAT_NEG_T; REAL_ARITH `--a < --b <=> b < a`] THEN
 756                           REWRITE_TAC[FLOAT_LT_FF]);;
 757
 758 let FLOAT_LT_FT = prove(`float_num F n1 e1 < float_num T n2 e2 <=> F`,
 759                         MP_TAC (SPECL [`n2:num`; `e2:num`] FLOAT_T_NEG) THEN
 760                           MP_TAC (SPECL [`n1:num`; `e1:num`] FLOAT_F_POS) THEN
 761                           REAL_ARITH_TAC);;
 762
 763 let FLOAT_LT_TF_00 = (PURE_REWRITE_RULE[NUMERAL] o prove)
 764   (`float_num T 0 e1 < float_num F 0 e2 <=> F`,
 765    MP_TAC (SPECL [`T`; `0`; `e1:num`] FLOAT_EQ_0) THEN
 766      MP_TAC (SPECL [`F`; `0`; `e2:num`] FLOAT_EQ_0) THEN
 767      REWRITE_TAC[] THEN
 768      REPLICATE_TAC 2 (DISCH_THEN (fun th -> REWRITE_TAC[th])) THEN
 769      REAL_ARITH_TAC);;
 770
 771 let FLOAT_LT_TF_1 = (UNDISCH_ALL o PURE_REWRITE_RULE[NUMERAL] o prove)
 772   (`(n1 = 0 <=> F) ==> (float_num T n1 e1 < float_num F n2 e2 <=> T)`,
 773    DISCH_TAC THEN
 774      MATCH_MP_TAC (REAL_ARITH `a < &0 /\ &0 <= b ==> (a < b <=> T)`) THEN
 775      REWRITE_TAC[FLOAT_F_POS] THEN
 776      MATCH_MP_TAC (REAL_ARITH `~(a = &0) /\ a <= &0 ==> a < &0`) THEN
 777      ASM_REWRITE_TAC[FLOAT_T_NEG; FLOAT_EQ_0]);;
 778
 779 let FLOAT_LT_TF_2 = (UNDISCH_ALL o PURE_REWRITE_RULE[NUMERAL] o prove)
 780   (`(n2 = 0 <=> F) ==> (float_num T n1 e1 < float_num F n2 e2 <=> T)`,
 781    DISCH_TAC THEN
 782      MATCH_MP_TAC (REAL_ARITH `a <= &0 /\ &0 < b ==> (a < b <=> T)`) THEN
 783      REWRITE_TAC[FLOAT_T_NEG] THEN
 784      MATCH_MP_TAC (REAL_ARITH `~(a = &0) /\ &0 <= a ==> &0 < a`) THEN
 785      ASM_REWRITE_TAC[FLOAT_F_POS; FLOAT_EQ_0]);;
 786
 787 let FLOAT_F_LT_0 = prove(`float_num F n e < &0 <=> F`,
 788                          MP_TAC (SPEC_ALL FLOAT_F_POS) THEN
 789                            REAL_ARITH_TAC);;
 790
 791 let FLOAT_T_LT_0 = (CONV_RULE (RAND_CONV (REWRITE_CONV[NUMERAL])) o prove)
 792   (`float_num T n e < &0 <=> (0 < n)`,
 793    REWRITE_TAC[REAL_ARITH `a < &0 <=> a <= &0 /\ ~(a = &0)`] THEN
 794      REWRITE_TAC[FLOAT_T_NEG; FLOAT_EQ_0] THEN
 795      ARITH_TAC);;
 796
 797 let FLOAT_F_GT_0 = (CONV_RULE (RAND_CONV (REWRITE_CONV[NUMERAL])) o prove)
 798   (`&0 < float_num F n e <=> 0 < n`,
 799    REWRITE_TAC[REAL_ARITH `&0 < a <=> &0 <= a /\ ~(a = &0)`] THEN
 800      REWRITE_TAC[FLOAT_F_POS; FLOAT_EQ_0] THEN
 801      ARITH_TAC);;
 802
 803 let FLOAT_T_GT_0 = prove(`&0 < float_num T n e <=> F`,
 804                          MP_TAC (SPEC_ALL FLOAT_T_NEG) THEN
 805                            REAL_ARITH_TAC);;
 806
 807
 808 (* float_lt0, float_gt0 *)
 809
 810 let float_lt0 f1 =
 811   let s, n_tm, e_tm = dest_float f1 in
 812   let inst = INST[n_tm, n_var_num; e_tm, e_var_num] in
 813     if s = "F" then
 814       inst FLOAT_F_LT_0
 815     else
 816       let gt_th = raw_gt0_hash_conv n_tm in
 817         TRANS (inst FLOAT_T_LT_0) gt_th;;
 818
 819 let float_gt0 f1 =
 820   let s, n_tm, e_tm = dest_float f1 in
 821   let inst = INST[n_tm, n_var_num; e_tm, e_var_num] in
 822     if s = "F" then
 823       let gt_th = raw_gt0_hash_conv n_tm in
 824         TRANS (inst FLOAT_F_GT_0) gt_th
 825     else
 826       inst FLOAT_T_GT_0;;
 827
 828 (* float_lt *)
 829 let float_lt f1 f2 =
 830   let s1, n1, e1 = dest_float f1 in
 831   let s2, n2, e2 = dest_float f2 in
 832   let inst = INST[n1, n1_var_num; e1, e1_var_num;
 833                   n2, n2_var_num; e2, e2_var_num] in
 834     if s1 = "F" then
 835       if s2 = "F" then
 836         (* FF *)
 837         let th0 = inst FLOAT_LT_FF in
 838         let ltm, tm2 = dest_comb (rand (concl th0)) in
 839         let lt_th = num_exp_lt (rand ltm) tm2 in
 840           TRANS th0 lt_th
 841       else
 842         (* FT *)
 843         inst FLOAT_LT_FT
 844     else
 845       if s2 = "F" then
 846         (* TF *)
 847         if (is_const n1 && is_const n2) then
 848           (* n1 = _0 and n2 = _0 *)
 849           inst FLOAT_LT_TF_00
 850         else
 851           let n1_0 = raw_eq0_hash_conv n1 in
 852             if (fst o dest_const o rand o concl) n1_0 = "F" then
 853               (* n1 <> _0 *)
 854               MY_PROVE_HYP n1_0 (inst FLOAT_LT_TF_1)
 855             else
 856               let n2_0 = raw_eq0_hash_conv n2 in
 857                 if (fst o dest_const o rand o concl) n2_0 = "F" then
 858                   (* n2 <> _0 *)
 859                   MY_PROVE_HYP n2_0 (inst FLOAT_LT_TF_2)
 860                 else
 861                   failwith "float_lt: D0 _0 exception"
 862       else
 863         (* TT *)
 864         let th0 = inst FLOAT_LT_TT in
 865         let ltm, tm2 = dest_comb (rand (concl th0)) in
 866         let lt_th = num_exp_lt (rand ltm) tm2 in
 867           TRANS th0 lt_th;;
 868
 869 (**********************************)
 870 (* FLOAT_LE *)
 871
 872 let FLOAT_LE_FF = prove(`float_num F n1 e1 <= float_num F n2 e2 <=> num_exp n1 e1 <= num_exp n2 e2`,
 873                         REWRITE_TAC[float; GSYM REAL_OF_NUM_LE; REAL_MUL_LID; real_div] THEN
 874                           MATCH_MP_TAC REAL_LE_RMUL_EQ THEN
 875                           MATCH_MP_TAC REAL_LT_INV THEN
 876                           REWRITE_TAC[REAL_OF_NUM_LT; LT_NZ; NUM_EXP_EQ_0] THEN
 877                           ARITH_TAC);;
 878
 879 let FLOAT_LE_TT = prove(`float_num T n1 e1 <= float_num T n2 e2 <=> num_exp n2 e2 <= num_exp n1 e1`,
 880                         REWRITE_TAC[FLOAT_NEG_T; REAL_ARITH `--a <= --b <=> b <= a`] THEN
 881                           REWRITE_TAC[FLOAT_LE_FF]);;
 882
 883 let FLOAT_LE_TF = prove(`float_num T n1 e1 <= float_num F n2 e2 <=> T`,
 884                         MP_TAC (SPECL [`n1:num`; `e1:num`] FLOAT_T_NEG) THEN
 885                           MP_TAC (SPECL [`n2:num`; `e2:num`] FLOAT_F_POS) THEN
 886                           REAL_ARITH_TAC);;
 887
 888 let FLOAT_LE_FT = prove(`float_num F n1 e1 <= float_num T n2 e2 <=> n1 = 0 /\ n2 = 0`,
 889                         REWRITE_TAC[REAL_LE_LT; FLOAT_LT_FT] THEN EQ_TAC THENL
 890                           [
 891                             DISCH_TAC THEN SUBGOAL_THEN `float_num F n1 e1 = &0 /\ float_num T n2 e2 = &0` MP_TAC THENL
 892                               [
 893                                 MP_TAC (SPECL [`n2:num`; `e2:num`] FLOAT_T_NEG) THEN
 894                                   MP_TAC (SPECL [`n1:num`; `e1:num`] FLOAT_F_POS) THEN
 895                                   ASM_REWRITE_TAC[] THEN REAL_ARITH_TAC;
 896                                 ALL_TAC
 897                               ] THEN
 898                               REWRITE_TAC[FLOAT_EQ_0];
 899                             DISCH_THEN (fun th -> REWRITE_TAC[th]) THEN
 900                               REWRITE_TAC[float; NUM_EXP_n0; real_div; REAL_MUL_LZERO; REAL_MUL_RZERO]
 901                           ]);;
 902
 903 let FLOAT_LE_FT_00 = (PURE_REWRITE_RULE[NUMERAL] o prove)
 904   (`float_num F 0 e1 <= float_num T 0 e2 <=> T`, REWRITE_TAC[FLOAT_LE_FT]);;
 905
 906 let FLOAT_LE_FT_1 = (UNDISCH_ALL o PURE_REWRITE_RULE[NUMERAL] o prove)
 907   (`(n1 = 0 <=> F) ==> (float_num F n1 e1 <= float_num T n2 e2 <=> F)`,
 908    DISCH_TAC THEN ASM_REWRITE_TAC[FLOAT_LE_FT]);;
 909
 910 let FLOAT_LE_FT_2 = (UNDISCH_ALL o PURE_REWRITE_RULE[NUMERAL] o prove)
 911   (`(n2 = 0 <=> F) ==> (float_num F n1 e1 <= float_num T n2 e2 <=> F)`,
 912    DISCH_TAC THEN ASM_REWRITE_TAC[FLOAT_LE_FT]);;
 913
 914 let FLOAT_F_LE_0 = (CONV_RULE (RAND_CONV (REWRITE_CONV[NUMERAL])) o prove)
 915   (`float_num F n e <= &0 <=> n = 0`,
 916    REWRITE_TAC[GSYM (SPEC `F` FLOAT_EQ_0)] THEN
 917      MP_TAC (SPEC_ALL FLOAT_F_POS) THEN
 918      REAL_ARITH_TAC);;
 919
 920 let FLOAT_T_LE_0 = prove(`float_num T n e <= &0 <=> T`, REWRITE_TAC[FLOAT_T_NEG]);;
 921
 922 let FLOAT_F_GE_0 = prove(`&0 <= float_num F n e <=> T`, REWRITE_TAC[FLOAT_F_POS]);;
 923
 924 let FLOAT_T_GE_0 = (CONV_RULE (RAND_CONV (REWRITE_CONV[NUMERAL])) o prove)
 925   (`&0 <= float_num T n e <=> n = 0`,
 926    REWRITE_TAC[GSYM (SPEC `T` FLOAT_EQ_0)] THEN
 927      MP_TAC (SPEC_ALL FLOAT_T_NEG) THEN
 928      REAL_ARITH_TAC);;
 929
 930 (* float_le0, float_ge0 *)
 931 let float_le0 f1 =
 932   let s, n_tm, e_tm = dest_float f1 in
 933   let inst = INST[n_tm, n_var_num; e_tm, e_var_num] in
 934     if s = "T" then
 935       inst FLOAT_T_LE_0
 936     else
 937       let eq_th = raw_eq0_hash_conv n_tm in
 938         TRANS (inst FLOAT_F_LE_0) eq_th;;
 939
 940 let float_ge0 f1 =
 941   let s, n_tm, e_tm = dest_float f1 in
 942   let inst = INST[n_tm, n_var_num; e_tm, e_var_num] in
 943     if s = "T" then
 944       let eq_th = raw_eq0_hash_conv n_tm in
 945         TRANS (inst FLOAT_T_GE_0) eq_th
 946     else
 947       inst FLOAT_F_GE_0;;
 948
 949 (* float_le *)
 950 let float_le f1 f2 =
 951   let s1, n1, e1 = dest_float f1 in
 952   let s2, n2, e2 = dest_float f2 in
 953   let inst = INST[n1, n1_var_num; e1, e1_var_num;
 954                   n2, n2_var_num; e2, e2_var_num] in
 955     if s2 = "F" then
 956       if s1 = "F" then
 957         (* FF *)
 958         let th0 = inst FLOAT_LE_FF in
 959         let ltm, tm2 = dest_comb (rand (concl th0)) in
 960         let le_th = num_exp_le (rand ltm) tm2 in
 961           TRANS th0 le_th
 962       else
 963         (* TF *)
 964         inst FLOAT_LE_TF
 965     else
 966       if s1 = "F" then
 967         (* FT *)
 968         if (is_const n1 && is_const n2) then
 969           (* n1 = _0 and n2 = _0 *)
 970           inst FLOAT_LE_FT_00
 971         else
 972           let n1_0 = raw_eq0_hash_conv n1 in
 973             if (fst o dest_const o rand o concl) n1_0 = "F" then
 974               (* n1 <> _0 *)
 975               MY_PROVE_HYP n1_0 (inst FLOAT_LE_FT_1)
 976             else
 977               let n2_0 = raw_eq0_hash_conv n2 in
 978                 if (fst o dest_const o rand o concl) n2_0 = "F" then
 979                   (* n2 <> _0 *)
 980                   MY_PROVE_HYP n2_0 (inst FLOAT_LE_FT_2)
 981                 else
 982                   failwith "float_lt: D0 _0 exception"
 983       else
 984         (* TT *)
 985         let th0 = inst FLOAT_LE_TT in
 986         let ltm, tm2 = dest_comb (rand (concl th0)) in
 987         let le_th = num_exp_le (rand ltm) tm2 in
 988           TRANS th0 le_th;;
 989
 990
 991 (*************************************)
 992 (* float_max, float_min *)
 993
 994 let FLOAT_MIN_1 = (UNDISCH_ALL o prove)(`(f1 <= f2 <=> T) ==> min f1 f2 = f1`, REAL_ARITH_TAC);;
 995 let FLOAT_MIN_2 = (UNDISCH_ALL o prove)(`(f1 <= f2 <=> F) ==> min f1 f2 = f2`, REAL_ARITH_TAC);;
 996
 997 let FLOAT_MAX_1 = (UNDISCH_ALL o prove)(`(f1 <= f2 <=> T) ==> max f1 f2 = f2`, REAL_ARITH_TAC);;
 998 let FLOAT_MAX_2 = (UNDISCH_ALL o prove)(`(f1 <= f2 <=> F) ==> max f1 f2 = f1`, REAL_ARITH_TAC);;
 999
1000 let float_min f1 f2 =
1001   let inst = INST[f1, f1_var_real; f2, f2_var_real] in
1002   let le_th = float_le f1 f2 in
1003   let th0 =
1004     if (fst o dest_const o rand o concl) le_th = "T" then
1005       inst FLOAT_MIN_1
1006     else
1007       inst FLOAT_MIN_2 in
1008     MY_PROVE_HYP le_th th0;;
1009
1010 let float_max f1 f2 =
1011   let inst = INST[f1, f1_var_real; f2, f2_var_real] in
1012   let le_th = float_le f1 f2 in
1013   let th0 =
1014     if (fst o dest_const o rand o concl) le_th = "T" then
1015       inst FLOAT_MAX_1
1016     else
1017       inst FLOAT_MAX_2 in
1018     MY_PROVE_HYP le_th th0;;
1019
1020 let float_min_max f1 f2 =
1021   let inst = INST[f1, f1_var_real; f2, f2_var_real] in
1022   let le_th = float_le f1 f2 in
1023   let th_min, th_max =
1024     if (fst o dest_const o rand o concl) le_th = "T" then
1025       inst FLOAT_MIN_1, inst FLOAT_MAX_1
1026     else
1027       inst FLOAT_MIN_2, inst FLOAT_MAX_2 in
1028     MY_PROVE_HYP le_th th_min, MY_PROVE_HYP le_th th_max;;
1029
1030
1031 (*************************************)
1032 (* FLOAT_MUL *)
1033
1034 let FLOAT_MUL = prove(`!s1 s2. min_exp <= e /\ num_exp n1 e1 * num_exp n2 e2 = num_exp n e
1035                           ==> float_num s1 n1 e1 * float_num s2 n2 e2 =
1036                               float_num (mod_plus s1 s2) n (e - min_exp)`,
1037    REPEAT STRIP_TAC THEN
1038      REWRITE_TAC[float] THEN
1039      ONCE_REWRITE_TAC[REAL_ARITH `(a * b / c) * (d * e / f) = (a * d) * (b * e) / c / f`] THEN
1040
1041      SUBGOAL_THEN `(if s1 then -- &1 else &1) * (if s2 then -- &1 else &1) = if mod_plus s1 s2 then -- &1 else &1` MP_TAC THENL
1042      [
1043        REWRITE_TAC[mod_plus] THEN
1044          COND_CASES_TAC THEN COND_CASES_TAC THEN
1045          REWRITE_TAC[REAL_ARITH `-- &1 * -- &1 = &1`; REAL_MUL_LID; REAL_MUL_RID];
1046        ALL_TAC
1047      ] THEN
1048
1049      DISCH_THEN (fun th -> REWRITE_TAC[th]) THEN
1050      REWRITE_TAC[real_div] THEN
1051      REWRITE_TAC[GSYM REAL_MUL_ASSOC] THEN
1052      REWRITE_TAC[REAL_EQ_MUL_LCANCEL] THEN
1053      DISJ2_TAC THEN
1054
1055      MP_TAC (SPECL[`n:num`; `e:num`; `min_exp`] NUM_EXP_SUB_lemma) THEN
1056      ASM_REWRITE_TAC[] THEN
1057      DISCH_THEN (fun th -> REWRITE_TAC[th]) THEN
1058      REWRITE_TAC[REAL_MUL_ASSOC] THEN
1059      ASM_REWRITE_TAC[REAL_OF_NUM_MUL]);;
1060
1061 let FLOAT_MUL_FF = prove(`min_exp <= e /\ e - min_exp = r /\ num_exp n1 e1 * num_exp n2 e2 = num_exp n e ==>
1062                              float_num F n1 e1 * float_num F n2 e2 = float_num F n r`,
1063                          SIMP_TAC[FLOAT_MUL; mod_plus]);;
1064 let FLOAT_MUL_FT = prove(`min_exp <= e /\ e - min_exp = r /\ num_exp n1 e1 * num_exp n2 e2 = num_exp n e ==>
1065                              float_num F n1 e1 * float_num T n2 e2 = float_num T n r`,
1066                          SIMP_TAC[FLOAT_MUL; mod_plus]);;
1067 let FLOAT_MUL_TF = prove(`min_exp <= e /\ e - min_exp = r /\ num_exp n1 e1 * num_exp n2 e2 = num_exp n e ==>
1068                              float_num T n1 e1 * float_num F n2 e2 = float_num T n r`,
1069                          SIMP_TAC[FLOAT_MUL; mod_plus]);;
1070 let FLOAT_MUL_TT = prove(`min_exp <= e /\ e - min_exp = r /\ num_exp n1 e1 * num_exp n2 e2 = num_exp n e ==>
1071                              float_num T n1 e1 * float_num T n2 e2 = float_num F n r`,
1072                          SIMP_TAC[FLOAT_MUL; mod_plus]);;
1073
1074
1075 let FLOAT_MUL_0x_hi, FLOAT_MUL_0x_lo, FLOAT_MUL_x0_hi, FLOAT_MUL_x0_lo =
1076   let mul_0x_hi = `(n1 = 0 <=> T) ==> float_num s1 n1 e1 * f2 <= float_num F 0 min_exp` in
1077   let mul_0x_lo = `(n1 = 0 <=> T) ==> float_num F 0 min_exp <= float_num s1 n1 e1 * f2` in
1078   let mul_x0_hi = `(n2 = 0 <=> T) ==> f1 * float_num s2 n2 e2 <= float_num F 0 min_exp` in
1079   let mul_x0_lo = `(n2 = 0 <=> T) ==> float_num F 0 min_exp <= f1 * float_num s2 n2 e2` in
1080   let proof = MP_TAC (GEN_ALL (SPECL [`s:bool`; `0`] FLOAT_EQ_0)) THEN
1081     SIMP_TAC[REAL_MUL_LZERO; REAL_MUL_RZERO; REAL_LE_REFL] in
1082     prove(mul_0x_hi, proof), prove(mul_0x_lo, proof),
1083   prove(mul_x0_hi, proof), prove(mul_x0_lo, proof);;
1084
1085
1086 let FLOAT_MUL_FF_hi, FLOAT_MUL_FF_lo =
1087   let ff_hi = `min_exp <= e /\ e - min_exp = r /\ num_exp n1 e1 * num_exp n2 e2 <= num_exp n e
1088                       ==> float_num F n1 e1 * float_num F n2 e2 <= float_num F n r` in
1089   let ff_lo = `min_exp <= e /\ e - min_exp = r /\ num_exp n e <= num_exp n1 e1 * num_exp n2 e2
1090                       ==> float_num F n r <= float_num F n1 e1 * float_num F n2 e2` in
1091   let proof =
1092     REPEAT STRIP_TAC THEN
1093       POP_ASSUM MP_TAC THEN
1094       POP_ASSUM (fun th -> REWRITE_TAC[SYM th]) THEN
1095       REWRITE_TAC[GSYM REAL_OF_NUM_LE; GSYM REAL_OF_NUM_MUL] THEN
1096       DISCH_TAC THEN
1097       MAP_EVERY ABBREV_TAC [`z = &(num_exp n e)`; `x = &(num_exp n1 e1)`; `y = &(num_exp n2 e2)`] THEN
1098       ASM_REWRITE_TAC[float; REAL_MUL_LID] THEN
1099       REWRITE_TAC[REAL_ARITH `a / b * c / d = (a * c) / b / d`] THEN
1100       REWRITE_TAC[real_div] THEN
1101       REWRITE_TAC[REAL_MUL_ASSOC] THEN
1102       MATCH_MP_TAC REAL_LE_RMUL THEN
1103       REWRITE_TAC[REAL_LE_INV_EQ; REAL_POS] THEN
1104
1105       MP_TAC (SPECL [`n:num`; `e:num`; `min_exp`] NUM_EXP_SUB_lemma) THEN
1106       ASM_REWRITE_TAC[] THEN
1107       DISCH_THEN (fun th -> REWRITE_TAC[th]) THEN
1108       MATCH_MP_TAC REAL_LE_RMUL THEN
1109       ASM_REWRITE_TAC[REAL_LE_INV_EQ; REAL_POS] in
1110     prove(ff_hi, proof), prove(ff_lo, proof);;
1111
1112 let FLOAT_MUL_TT_hi, FLOAT_MUL_TT_lo =
1113   let tt_hi = `min_exp <= e /\ e - min_exp = r /\ num_exp n1 e1 * num_exp n2 e2 <= num_exp n e
1114                      ==> float_num T n1 e1 * float_num T n2 e2 <= float_num F n r` in
1115   let tt_lo = `min_exp <= e /\ e - min_exp = r /\ num_exp n e <= num_exp n1 e1 * num_exp n2 e2
1116                      ==> float_num F n r <= float_num T n1 e1 * float_num T n2 e2` in
1117   let proof =
1118     REWRITE_TAC[FLOAT_NEG_T] THEN
1119       REWRITE_TAC[REAL_ARITH `--a * --b = a * b`] THEN
1120       REWRITE_TAC[FLOAT_MUL_FF_hi; FLOAT_MUL_FF_lo] in
1121     prove(tt_hi, proof), prove(tt_lo, proof);;
1122
1123 let FLOAT_MUL_FT_hi, FLOAT_MUL_FT_lo =
1124   let ft_hi = `min_exp <= e /\ e - min_exp = r /\ num_exp n e <= num_exp n1 e1 * num_exp n2 e2
1125                      ==> float_num F n1 e1 * float_num T n2 e2 <= float_num T n r` in
1126   let ft_lo = `min_exp <= e /\ e - min_exp = r /\ num_exp n1 e1 * num_exp n2 e2 <= num_exp n e
1127                      ==> float_num T n r <= float_num F n1 e1 * float_num T n2 e2` in
1128   let proof =
1129     REWRITE_TAC[FLOAT_NEG_T] THEN
1130       REWRITE_TAC[REAL_ARITH `a * --b <= --c <=> c <= a * b`] THEN
1131       REWRITE_TAC[REAL_ARITH `--c <= a * --b <=> a * b <= c`] THEN
1132       REWRITE_TAC[FLOAT_MUL_FF_hi; FLOAT_MUL_FF_lo] in
1133     prove(ft_hi, proof), prove(ft_lo, proof);;
1134
1135 let FLOAT_MUL_TF_hi, FLOAT_MUL_TF_lo =
1136   let ft_hi = `min_exp <= e /\ e - min_exp = r /\ num_exp n e <= num_exp n1 e1 * num_exp n2 e2
1137                      ==> float_num T n1 e1 * float_num F n2 e2 <= float_num T n r` in
1138   let ft_lo = `min_exp <= e /\ e - min_exp = r /\ num_exp n1 e1 * num_exp n2 e2 <= num_exp n e
1139                      ==> float_num T n r <= float_num T n1 e1 * float_num F n2 e2` in
1140   let proof =
1141     REWRITE_TAC[FLOAT_NEG_T] THEN
1142       REWRITE_TAC[REAL_ARITH `--a * b <= --c <=> c <= a * b`] THEN
1143       REWRITE_TAC[REAL_ARITH `--c <= --a * b <=> a * b <= c`] THEN
1144       REWRITE_TAC[FLOAT_MUL_FF_hi; FLOAT_MUL_FF_lo] in
1145     prove(ft_hi, proof), prove(ft_lo, proof);;
1146
1147
1148 (*********************************************)
1149 (* float_mul_lo, float_mul_hi *)
1150
1151 let transform = UNDISCH_ALL o NUMERALS_TO_NUM o PURE_REWRITE_RULE[min_exp_def; GSYM IMP_IMP];;
1152 let FLOAT_MUL_FF_hi' = transform FLOAT_MUL_FF_hi and
1153     FLOAT_MUL_FF_lo' = transform FLOAT_MUL_FF_lo and
1154     FLOAT_MUL_TT_hi' = transform FLOAT_MUL_TT_hi and
1155     FLOAT_MUL_TT_lo' = transform FLOAT_MUL_TT_lo and
1156     FLOAT_MUL_FT_hi' = transform FLOAT_MUL_FT_hi and
1157     FLOAT_MUL_FT_lo' = transform FLOAT_MUL_FT_lo and
1158     FLOAT_MUL_TF_hi' = transform FLOAT_MUL_TF_hi and
1159     FLOAT_MUL_TF_lo' = transform FLOAT_MUL_TF_lo and
1160     FLOAT_MUL_0x_hi' = transform FLOAT_MUL_0x_hi and
1161     FLOAT_MUL_0x_lo' = transform FLOAT_MUL_0x_lo and
1162     FLOAT_MUL_x0_hi' = transform FLOAT_MUL_x0_hi and
1163     FLOAT_MUL_x0_lo' = transform FLOAT_MUL_x0_lo;;
1164
1165 let FLOAT_MUL_FF' = transform FLOAT_MUL_FF and
1166     FLOAT_MUL_TT' = transform FLOAT_MUL_TT and
1167     FLOAT_MUL_FT' = transform FLOAT_MUL_FT and
1168     FLOAT_MUL_TF' = transform FLOAT_MUL_TF;;
1169
1170 let float_mul_eq f1 f2 =
1171   let s1, n1, e1 = dest_float f1 and
1172       s2, n2, e2 = dest_float f2 in
1173   let flag = s1 = s2 in
1174   let num_exp1 = mk_num_exp n1 e1 and
1175       num_exp2 = mk_num_exp n2 e2 in
1176
1177   let mul_th = num_exp_mul num_exp1 num_exp2 in
1178   let n_tm, e_tm = dest_num_exp (rand (concl mul_th)) in
1179
1180   let sub_th, le_th = raw_sub_and_le_hash_conv e_tm min_exp_num_const in
1181     if (rand(concl le_th) <> e_tm) then
1182       failwith "float_mul_eq: underflow"
1183     else
1184       let r_tm = rand(concl sub_th) in
1185       let inst = INST[e_tm, e_var_num; r_tm, r_var_num; n_tm, n_var_num;
1186                       n1, n1_var_num; e1, e1_var_num; n2, n2_var_num; e2, e2_var_num] in
1187       let th0 = inst
1188         (if flag then
1189            if s1 = "F" then FLOAT_MUL_FF' else FLOAT_MUL_TT'
1190          else
1191            if s1 = "F" then FLOAT_MUL_FT' else FLOAT_MUL_TF') in
1192         MY_PROVE_HYP sub_th (MY_PROVE_HYP mul_th (MY_PROVE_HYP le_th th0));;
1193
1194
1195 let float_mul_lo pp f1 f2 =
1196   let s1, n1, e1 = dest_float f1 and
1197       s2, n2, e2 = dest_float f2 in
1198     (* Multiplication by zero *)
1199   let n1_eq0_th = raw_eq0_hash_conv n1 in
1200     if (rand o concl) n1_eq0_th = t_const then
1201       (MY_PROVE_HYP n1_eq0_th o
1202          INST[e1, e1_var_num; f2, f2_var_real; n1, n1_var_num;
1203               (if s1 = "T" then t_const else f_const), s1_var_bool]) FLOAT_MUL_0x_lo'
1204     else
1205       let n2_eq0_th = raw_eq0_hash_conv n2 in
1206         if (rand o concl) n2_eq0_th = t_const then
1207           (MY_PROVE_HYP n2_eq0_th o
1208              INST[e2, e2_var_num; f1, f1_var_real; n2, n2_var_num;
1209                   (if s2 = "T" then t_const else f_const), s2_var_bool]) FLOAT_MUL_x0_lo'
1210         else
1211           let flag = s1 = s2 in
1212           let num_exp1 = mk_num_exp n1 e1 and
1213               num_exp2 = mk_num_exp n2 e2 in
1214
1215           let mul_th, n_tm, e_tm =
1216             if flag then
1217               let th = num_exp_op_lo pp num_exp_mul num_exp1 num_exp2 in
1218               let n_tm, e_tm = dest_num_exp (lhand (concl th)) in
1219                 th, n_tm, e_tm
1220             else
1221               let th = num_exp_op_hi pp num_exp_mul num_exp1 num_exp2 in
1222               let n_tm, e_tm = dest_num_exp (rand (concl th)) in
1223                 th, n_tm, e_tm in
1224
1225           let sub_th, le_th = raw_sub_and_le_hash_conv e_tm min_exp_num_const in
1226             if (rand(concl le_th) <> e_tm) then
1227               failwith "float_mul_lo: underflow"
1228             else
1229               let r_tm = rand(concl sub_th) in
1230               let inst = INST[e_tm, e_var_num; r_tm, r_var_num; n_tm, n_var_num;
1231                               n1, n1_var_num; e1, e1_var_num; n2, n2_var_num; e2, e2_var_num] in
1232               let th0 = inst
1233                 (if flag then
1234                    if s1 = "F" then FLOAT_MUL_FF_lo' else FLOAT_MUL_TT_lo'
1235                  else
1236                    if s1 = "F" then FLOAT_MUL_FT_lo' else FLOAT_MUL_TF_lo') in
1237                 MY_PROVE_HYP sub_th (MY_PROVE_HYP mul_th (MY_PROVE_HYP le_th th0));;
1238
1239
1240 let float_mul_hi pp f1 f2 =
1241   let s1, n1, e1 = dest_float f1 and
1242       s2, n2, e2 = dest_float f2 in
1243     (* Multiplication by zero *)
1244   let n1_eq0_th = raw_eq0_hash_conv n1 in
1245     if (rand o concl) n1_eq0_th = t_const then
1246       (MY_PROVE_HYP n1_eq0_th o
1247          INST[e1, e1_var_num; f2, f2_var_real; n1, n1_var_num;
1248               (if s1 = "T" then t_const else f_const), s1_var_bool]) FLOAT_MUL_0x_hi'
1249     else
1250       let n2_eq0_th = raw_eq0_hash_conv n2 in
1251         if (rand o concl) n2_eq0_th = t_const then
1252           (MY_PROVE_HYP n2_eq0_th o
1253              INST[e2, e2_var_num; f1, f1_var_real; n2, n2_var_num;
1254                   (if s2 = "T" then t_const else f_const), s2_var_bool]) FLOAT_MUL_x0_hi'
1255         else
1256           let flag = s1 = s2 in
1257           let num_exp1 = mk_num_exp n1 e1 and
1258               num_exp2 = mk_num_exp n2 e2 in
1259
1260           let mul_th, n_tm, e_tm =
1261             if flag then
1262               let th = num_exp_op_hi pp num_exp_mul num_exp1 num_exp2 in
1263               let n_tm, e_tm = dest_num_exp (rand (concl th)) in
1264                 th, n_tm, e_tm
1265             else
1266               let th = num_exp_op_lo pp num_exp_mul num_exp1 num_exp2 in
1267               let n_tm, e_tm = dest_num_exp (lhand (concl th)) in
1268                 th, n_tm, e_tm in
1269
1270           let sub_th, le_th = raw_sub_and_le_hash_conv e_tm min_exp_num_const in
1271             if (rand(concl le_th) <> e_tm) then
1272               failwith "float_mul_hi: underflow"
1273             else
1274               let r_tm = rand(concl sub_th) in
1275               let inst = INST[e_tm, e_var_num; r_tm, r_var_num; n_tm, n_var_num;
1276                               n1, n1_var_num; e1, e1_var_num; n2, n2_var_num; e2, e2_var_num] in
1277               let th0 = inst
1278                 (if flag then
1279                    if s1 = "F" then FLOAT_MUL_FF_hi' else FLOAT_MUL_TT_hi'
1280                  else
1281                    if s1 = "F" then FLOAT_MUL_FT_hi' else FLOAT_MUL_TF_hi') in
1282                 MY_PROVE_HYP sub_th (MY_PROVE_HYP mul_th (MY_PROVE_HYP le_th th0));;
1283
1284
1285
1286 (*********************************************)
1287 (* FLOAT_DIV *)
1288
1289 let DIV_lemma = prove(`!x y. ~(y = 0) ==> &(x DIV y) <= &x / &y /\ &x / &y <= &(x DIV y + 1)`,
1290    REPEAT GEN_TAC THEN DISCH_TAC THEN
1291      MP_TAC (SPECL [`y:num`; `x:num`] FLOOR_DIV_DIV) THEN
1292      ASM_REWRITE_TAC[GSYM REAL_OF_NUM_ADD] THEN
1293      DISCH_THEN (fun th -> REWRITE_TAC[SYM th]) THEN
1294      SIMP_TAC[FLOOR; REAL_LT_IMP_LE]);;
1295
1296 let FLOAT_DIV_FF = prove(`e2 + k <= min_exp + e + e1 /\ ~(n2 = 0) /\
1297                            num_exp n1 k DIV num_exp n2 0 = num_exp n e
1298                              ==> float_num F n ((min_exp + e + e1) - (e2 + k)) <= float_num F n1 e1 / float_num F n2 e2`,
1299    MAP_EVERY ABBREV_TAC [`z = num_exp n e`; `x = num_exp n1 k`; `y = num_exp n2 0`] THEN
1300      REPEAT STRIP_TAC THEN
1301      REWRITE_TAC[float; REAL_MUL_LID] THEN
1302      REWRITE_TAC[real_div; REAL_INV_MUL; REAL_INV_INV] THEN
1303      REWRITE_TAC[REAL_ARITH `(a * b) * c * d = (b * d) * (a * c)`] THEN
1304      SUBGOAL_THEN `~(&(num_exp 1 min_exp) = &0)` ASSUME_TAC THENL
1305      [
1306        REWRITE_TAC[num_exp; REAL_OF_NUM_EQ; MULT_CLAUSES; EXP_EQ_0] THEN
1307          ARITH_TAC;
1308        ALL_TAC
1309      ] THEN
1310
1311      ASM_SIMP_TAC[REAL_MUL_LINV; REAL_MUL_LID] THEN
1312
1313      ASM_SIMP_TAC[NUM_EXP_SUB_lemma] THEN
1314      SUBGOAL_THEN `&(num_exp n1 e1) * inv(&(num_exp n2 e2)) = (&x / &y) * &(num_exp 1 e1) * inv(&(num_exp 1 (e2 + k)))` MP_TAC THENL
1315      [
1316        EXPAND_TAC "x" THEN EXPAND_TAC "y" THEN
1317          REWRITE_TAC[real_div] THEN
1318          REWRITE_TAC[num_exp; GSYM REAL_OF_NUM_MUL; GSYM REAL_OF_NUM_POW] THEN
1319          REWRITE_TAC[REAL_MUL_LID; REAL_INV_MUL; REAL_INV_1; real_pow; REAL_MUL_RID] THEN
1320          REWRITE_TAC[REAL_POW_ADD; REAL_INV_MUL] THEN
1321          REWRITE_TAC[REAL_ARITH `((a * b) * c) * d * e * f = (b * f) * (a * c * d * e)`] THEN
1322
1323          SUBGOAL_THEN (mk_comb(`(~)`, mk_eq(mk_binop `pow` (mk_comb (`&`, base_const)) `k:num`, `&0`))) ASSUME_TAC THENL
1324          [
1325            REWRITE_TAC[REAL_POW_EQ_0] THEN
1326              REAL_ARITH_TAC;
1327            ALL_TAC
1328          ] THEN
1329
1330          ASM_SIMP_TAC[REAL_MUL_RINV; REAL_MUL_LID] THEN
1331          REAL_ARITH_TAC;
1332        ALL_TAC
1333      ] THEN
1334
1335      DISCH_THEN (fun th -> REWRITE_TAC[th]) THEN
1336      ONCE_REWRITE_TAC[REAL_ARITH `(a * b) * c = (a * c) * b`] THEN
1337      REWRITE_TAC[REAL_MUL_ASSOC] THEN
1338      MATCH_MP_TAC REAL_LE_RMUL THEN
1339      REWRITE_TAC[REAL_LE_INV_EQ; REAL_POS] THEN
1340      ONCE_REWRITE_TAC[NUM_EXP_SUM1] THEN
1341      REWRITE_TAC[NUM_EXP_SUM] THEN
1342      REWRITE_TAC[GSYM REAL_OF_NUM_MUL] THEN
1343      ASM_REWRITE_TAC[REAL_ARITH `(a * b * c) * d = (d * a) * b * c`] THEN
1344      ASM_SIMP_TAC[REAL_MUL_LINV; REAL_MUL_LID] THEN
1345      MATCH_MP_TAC REAL_LE_RMUL THEN
1346      REWRITE_TAC[REAL_LE_INV_EQ; REAL_POS] THEN
1347      MP_TAC (SPEC_ALL DIV_lemma) THEN
1348      ANTS_TAC THENL
1349      [
1350        EXPAND_TAC "y" THEN
1351          REWRITE_TAC[num_exp; MULT_EQ_0; DE_MORGAN_THM] THEN
1352          ASM_REWRITE_TAC[EXP] THEN
1353          ARITH_TAC;
1354        ALL_TAC
1355      ] THEN
1356
1357      ASM_SIMP_TAC[]);;
1358
1359 let FLOAT_DIV_0x_lo = prove(`(n1 = 0 <=> T) ==> float_num F 0 min_exp <= float_num s1 n1 e1 / f2`,
1360                             SIMP_TAC[real_div; FLOAT_MUL_0x_lo]);;
1361
1362 let FLOAT_DIV_0x_hi = prove(`(n1 = 0 <=> T) ==> float_num s1 n1 e1 / f2 <= float_num F 0 min_exp`,
1363                             SIMP_TAC[real_div; FLOAT_MUL_0x_hi]);;
1364
1365 let FLOAT_DIV_FF_lo = prove(`e2 + k = r1 /\ min_exp + e + e1 = r2 /\ r2 - r1 = r /\
1366                              r1 <= r2 /\ ~(n2 = 0) /\
1367                            num_exp n e <= num_exp n1 k DIV num_exp n2 0
1368                              ==> float_num F n r <= float_num F n1 e1 / float_num F n2 e2`,
1369    MAP_EVERY ABBREV_TAC [`z = num_exp n e`; `x = num_exp n1 k`; `y = num_exp n2 0`] THEN
1370      REPEAT STRIP_TAC THEN
1371      REPLICATE_TAC 3 (POP_ASSUM MP_TAC) THEN
1372      REPLICATE_TAC 3 (POP_ASSUM (fun th -> REWRITE_TAC[SYM th])) THEN
1373      REPEAT STRIP_TAC THEN
1374      REWRITE_TAC[float; REAL_MUL_LID] THEN
1375      REWRITE_TAC[real_div; REAL_INV_MUL; REAL_INV_INV] THEN
1376      REWRITE_TAC[REAL_ARITH `(a * b) * c * d = (b * d) * (a * c)`] THEN
1377      SUBGOAL_THEN `~(&(num_exp 1 min_exp) = &0)` ASSUME_TAC THENL
1378      [
1379        REWRITE_TAC[num_exp; REAL_OF_NUM_EQ; MULT_CLAUSES; EXP_EQ_0] THEN
1380          ARITH_TAC;
1381        ALL_TAC
1382      ] THEN
1383
1384      ASM_SIMP_TAC[REAL_MUL_LINV; REAL_MUL_LID] THEN
1385
1386      ASM_SIMP_TAC[NUM_EXP_SUB_lemma] THEN
1387      SUBGOAL_THEN `&(num_exp n1 e1) * inv(&(num_exp n2 e2)) = (&x / &y) * &(num_exp 1 e1) * inv(&(num_exp 1 (e2 + k)))` MP_TAC THENL
1388      [
1389        EXPAND_TAC "x" THEN EXPAND_TAC "y" THEN
1390          REWRITE_TAC[real_div] THEN
1391          REWRITE_TAC[num_exp; GSYM REAL_OF_NUM_MUL; GSYM REAL_OF_NUM_POW] THEN
1392          REWRITE_TAC[REAL_MUL_LID; REAL_INV_MUL; REAL_INV_1; real_pow; REAL_MUL_RID] THEN
1393          REWRITE_TAC[REAL_POW_ADD; REAL_INV_MUL] THEN
1394          REWRITE_TAC[REAL_ARITH `((a * b) * c) * d * e * f = (b * f) * (a * c * d * e)`] THEN
1395          SUBGOAL_THEN
1396          (mk_comb(`(~)`, mk_eq(mk_binop `pow` (mk_comb (`&`, base_const)) `k:num`, `&0`))) ASSUME_TAC THENL
1397          [
1398            REWRITE_TAC[REAL_POW_EQ_0] THEN
1399              REAL_ARITH_TAC;
1400            ALL_TAC
1401          ] THEN
1402
1403          ASM_SIMP_TAC[REAL_MUL_RINV; REAL_MUL_LID] THEN
1404          REAL_ARITH_TAC;
1405        ALL_TAC
1406      ] THEN
1407
1408      DISCH_THEN (fun th -> REWRITE_TAC[th]) THEN
1409      ONCE_REWRITE_TAC[REAL_ARITH `(a * b) * c = (a * c) * b`] THEN
1410      REWRITE_TAC[REAL_MUL_ASSOC] THEN
1411      MATCH_MP_TAC REAL_LE_RMUL THEN
1412      REWRITE_TAC[REAL_LE_INV_EQ; REAL_POS] THEN
1413      ONCE_REWRITE_TAC[NUM_EXP_SUM1] THEN
1414      REWRITE_TAC[NUM_EXP_SUM] THEN
1415      REWRITE_TAC[GSYM REAL_OF_NUM_MUL] THEN
1416      ASM_REWRITE_TAC[REAL_ARITH `(a * b * c) * d = (d * a) * b * c`] THEN
1417      ASM_SIMP_TAC[REAL_MUL_LINV; REAL_MUL_LID] THEN
1418      MATCH_MP_TAC REAL_LE_RMUL THEN
1419      REWRITE_TAC[REAL_LE_INV_EQ; REAL_POS] THEN
1420      MP_TAC (SPEC_ALL DIV_lemma) THEN
1421      ANTS_TAC THENL
1422      [
1423        EXPAND_TAC "y" THEN
1424          REWRITE_TAC[num_exp; MULT_EQ_0; DE_MORGAN_THM] THEN
1425          ASM_REWRITE_TAC[EXP] THEN
1426          ARITH_TAC;
1427        ALL_TAC
1428      ] THEN
1429
1430      STRIP_TAC THEN
1431      MATCH_MP_TAC REAL_LE_TRANS THEN
1432      EXISTS_TAC `&(x DIV y)` THEN
1433      ASM_REWRITE_TAC[REAL_OF_NUM_LE]);;
1434
1435 let FLOAT_DIV_FF_hi = prove(`e2 + k = r1 /\ min_exp + e + e1 = r2 /\ r2 - r1 = r /\
1436                               r1 <= r2 /\ ~(n2 = 0) /\
1437                               num_exp n1 k DIV num_exp n2 0 < num_exp n e
1438                               ==> float_num F n1 e1 / float_num F n2 e2 <= float_num F n r`,
1439    MAP_EVERY ABBREV_TAC [`z = num_exp n e`; `x = num_exp n1 k`; `y = num_exp n2 0`] THEN
1440      REPEAT STRIP_TAC THEN
1441      REPLICATE_TAC 3 (POP_ASSUM MP_TAC) THEN
1442      REPLICATE_TAC 3 (POP_ASSUM (fun th -> REWRITE_TAC[SYM th])) THEN
1443      REPEAT STRIP_TAC THEN
1444      REWRITE_TAC[float; REAL_MUL_LID] THEN
1445      REWRITE_TAC[real_div; REAL_INV_MUL; REAL_INV_INV] THEN
1446      REWRITE_TAC[REAL_ARITH `(a * b) * c * d = (b * d) * (a * c)`] THEN
1447      SUBGOAL_THEN `~(&(num_exp 1 min_exp) = &0)` ASSUME_TAC THENL
1448      [
1449        REWRITE_TAC[num_exp; REAL_OF_NUM_EQ; MULT_CLAUSES; EXP_EQ_0] THEN
1450          ARITH_TAC;
1451        ALL_TAC
1452      ] THEN
1453
1454      ASM_SIMP_TAC[REAL_MUL_LINV; REAL_MUL_LID] THEN
1455
1456      ASM_SIMP_TAC[NUM_EXP_SUB_lemma] THEN
1457      SUBGOAL_THEN `&(num_exp n1 e1) * inv(&(num_exp n2 e2)) = (&x / &y) * &(num_exp 1 e1) * inv(&(num_exp 1 (e2 + k)))` MP_TAC THENL
1458      [
1459        EXPAND_TAC "x" THEN EXPAND_TAC "y" THEN
1460          REWRITE_TAC[real_div] THEN
1461          REWRITE_TAC[num_exp; GSYM REAL_OF_NUM_MUL; GSYM REAL_OF_NUM_POW] THEN
1462          REWRITE_TAC[REAL_MUL_LID; REAL_INV_MUL; REAL_INV_1; real_pow; REAL_MUL_RID] THEN
1463          REWRITE_TAC[REAL_POW_ADD; REAL_INV_MUL] THEN
1464          REWRITE_TAC[REAL_ARITH `((a * b) * c) * d * e * f = (b * f) * (a * c * d * e)`] THEN
1465          SUBGOAL_THEN
1466          (mk_comb(`(~)`, mk_eq(mk_binop `pow` (mk_comb (`&`, base_const)) `k:num`, `&0`))) ASSUME_TAC THENL
1467          [
1468            REWRITE_TAC[REAL_POW_EQ_0] THEN
1469              REAL_ARITH_TAC;
1470            ALL_TAC
1471          ] THEN
1472
1473          ASM_SIMP_TAC[REAL_MUL_RINV; REAL_MUL_LID] THEN
1474          REAL_ARITH_TAC;
1475        ALL_TAC
1476      ] THEN
1477
1478      DISCH_THEN (fun th -> REWRITE_TAC[th]) THEN
1479      ONCE_REWRITE_TAC[REAL_ARITH `(a * b) * c = (a * c) * b`] THEN
1480      REWRITE_TAC[REAL_MUL_ASSOC] THEN
1481      MATCH_MP_TAC REAL_LE_RMUL THEN
1482      REWRITE_TAC[REAL_LE_INV_EQ; REAL_POS] THEN
1483      ONCE_REWRITE_TAC[NUM_EXP_SUM1] THEN
1484      REWRITE_TAC[NUM_EXP_SUM] THEN
1485      REWRITE_TAC[GSYM REAL_OF_NUM_MUL] THEN
1486      ASM_REWRITE_TAC[REAL_ARITH `(a * b * c) * d = (d * a) * b * c`] THEN
1487      ASM_SIMP_TAC[REAL_MUL_LINV; REAL_MUL_LID] THEN
1488      MATCH_MP_TAC REAL_LE_RMUL THEN
1489      REWRITE_TAC[REAL_LE_INV_EQ; REAL_POS] THEN
1490      MP_TAC (SPEC_ALL DIV_lemma) THEN
1491      ANTS_TAC THENL
1492      [
1493        EXPAND_TAC "y" THEN
1494          REWRITE_TAC[num_exp; MULT_EQ_0; DE_MORGAN_THM] THEN
1495          ASM_REWRITE_TAC[EXP] THEN
1496          ARITH_TAC;
1497        ALL_TAC
1498      ] THEN
1499
1500      STRIP_TAC THEN
1501      MATCH_MP_TAC REAL_LE_TRANS THEN
1502      EXISTS_TAC `&(x DIV y + 1)` THEN
1503      ASM_REWRITE_TAC[REAL_OF_NUM_LE] THEN
1504      UNDISCH_TAC `x DIV y < z` THEN
1505      ARITH_TAC);;
1506
1507 let FLOAT_DIV_TT_lo = prove(`e2 + k = r1 /\ min_exp + e + e1 = r2 /\ r2 - r1 = r /\
1508                               r1 <= r2 /\ ~(n2 = 0) /\
1509                               num_exp n e <= num_exp n1 k DIV num_exp n2 0
1510              ==> float_num F n r <= float_num T n1 e1 / float_num T n2 e2`,
1511    REWRITE_TAC[FLOAT_NEG_T] THEN
1512      REWRITE_TAC[real_div; REAL_INV_NEG; REAL_NEG_MUL2] THEN
1513      REWRITE_TAC[GSYM real_div] THEN
1514      REWRITE_TAC[FLOAT_DIV_FF_lo]);;
1515
1516 let FLOAT_DIV_TT_hi = prove(`e2 + k = r1 /\ min_exp + e + e1 = r2 /\
1517                               r2 - r1 = r /\ r1 <= r2 /\ ~(n2 = 0) /\
1518                               num_exp n1 k DIV num_exp n2 0 < num_exp n e
1519              ==> float_num T n1 e1 / float_num T n2 e2 <= float_num F n r`,
1520    REWRITE_TAC[FLOAT_NEG_T] THEN
1521      REWRITE_TAC[real_div; REAL_INV_NEG; REAL_NEG_MUL2] THEN
1522      REWRITE_TAC[GSYM real_div] THEN
1523      REWRITE_TAC[FLOAT_DIV_FF_hi]);;
1524
1525
1526 let FLOAT_DIV_FT_lo = prove(`e2 + k = r1 /\ min_exp + e + e1 = r2 /\
1527                               r2 - r1 = r /\ r1 <= r2 /\ ~(n2 = 0) /\
1528                               num_exp n1 k DIV num_exp n2 0 < num_exp n e
1529              ==> float_num T n r <= float_num F n1 e1 / float_num T n2 e2`,
1530    REWRITE_TAC[FLOAT_NEG_T] THEN
1531      REWRITE_TAC[real_div; REAL_INV_NEG] THEN
1532      REWRITE_TAC[REAL_ARITH `--a <= b * --c <=> b * c <= a`] THEN
1533      REWRITE_TAC[GSYM real_div] THEN
1534      REWRITE_TAC[FLOAT_DIV_FF_hi]);;
1535
1536 let FLOAT_DIV_FT_hi = prove(`e2 + k = r1 /\ min_exp + e + e1 = r2 /\
1537                               r2 - r1 = r /\ r1 <= r2 /\ ~(n2 = 0) /\
1538                               num_exp n e <= num_exp n1 k DIV num_exp n2 0
1539              ==> float_num F n1 e1 / float_num T n2 e2 <= float_num T n r`,
1540    REWRITE_TAC[FLOAT_NEG_T] THEN
1541      REWRITE_TAC[real_div; REAL_INV_NEG] THEN
1542      REWRITE_TAC[REAL_ARITH `a * --b <= --c <=> c <= a * b`] THEN
1543      REWRITE_TAC[GSYM real_div] THEN
1544      REWRITE_TAC[FLOAT_DIV_FF_lo]);;
1545
1546
1547 let FLOAT_DIV_TF_lo = prove(`e2 + k = r1 /\ min_exp + e + e1 = r2 /\
1548                               r2 - r1 = r /\ r1 <= r2 /\ ~(n2 = 0) /\
1549                               num_exp n1 k DIV num_exp n2 0 < num_exp n e
1550              ==> float_num T n r <= float_num T n1 e1 / float_num F n2 e2`,
1551    REWRITE_TAC[FLOAT_NEG_T] THEN
1552      REWRITE_TAC[real_div; REAL_INV_NEG] THEN
1553      REWRITE_TAC[REAL_ARITH `--a <= --b * c <=> b * c <= a`] THEN
1554      REWRITE_TAC[GSYM real_div] THEN
1555      REWRITE_TAC[FLOAT_DIV_FF_hi]);;
1556
1557 let FLOAT_DIV_TF_hi = prove(`e2 + k = r1 /\ min_exp + e + e1 = r2 /\
1558                               r2 - r1 = r /\ r1 <= r2 /\ ~(n2 = 0) /\
1559                               num_exp n e <= num_exp n1 k DIV num_exp n2 0
1560              ==> float_num T n1 e1 / float_num F n2 e2 <= float_num T n r`,
1561    REWRITE_TAC[FLOAT_NEG_T] THEN
1562      REWRITE_TAC[real_div; REAL_INV_NEG] THEN
1563      REWRITE_TAC[REAL_ARITH `--a * b <= --c <=> c <= a * b`] THEN
1564      REWRITE_TAC[GSYM real_div] THEN
1565      REWRITE_TAC[FLOAT_DIV_FF_lo]);;
1566
1567
1568 (******************************************)
1569 (* float_div_lo, float_div_hi *)
1570
1571 let transform = UNDISCH_ALL o PURE_REWRITE_RULE[TAUT `~P <=> (P <=> F)`] o
1572   NUMERALS_TO_NUM o PURE_REWRITE_RULE[GSYM IMP_IMP; min_exp_def];;
1573
1574 let FLOAT_DIV_FF_hi' = transform FLOAT_DIV_FF_hi and
1575     FLOAT_DIV_FF_lo' = transform FLOAT_DIV_FF_lo and
1576     FLOAT_DIV_TT_hi' = transform FLOAT_DIV_TT_hi and
1577     FLOAT_DIV_TT_lo' = transform FLOAT_DIV_TT_lo and
1578     FLOAT_DIV_FT_hi' = transform FLOAT_DIV_FT_hi and
1579     FLOAT_DIV_FT_lo' = transform FLOAT_DIV_FT_lo and
1580     FLOAT_DIV_TF_hi' = transform FLOAT_DIV_TF_hi and
1581     FLOAT_DIV_TF_lo' = transform FLOAT_DIV_TF_lo and
1582     FLOAT_DIV_0x_hi' = transform FLOAT_DIV_0x_hi and
1583     FLOAT_DIV_0x_lo' = transform FLOAT_DIV_0x_lo;;
1584
1585 let float_div_lo pp f1 f2 =
1586   let s1, n1, e1 = dest_float f1 and
1587       s2, n2, e2 = dest_float f2 in
1588   let n1_eq0_th = raw_eq0_hash_conv n1 in
1589     if (rand o concl) n1_eq0_th = t_const then
1590       (MY_PROVE_HYP n1_eq0_th o
1591          INST[e1, e1_var_num; f2, f2_var_real; n1, n1_var_num;
1592               (if s1 = "T" then t_const else f_const), s1_var_bool]) FLOAT_DIV_0x_lo'
1593     else
1594       let flag = s1 = s2 in
1595
1596       let k_tm = rand (mk_small_numeral_array (2 * pp)) in
1597       let num_exp1 = mk_num_exp n1 k_tm and
1598           num_exp2 = mk_num_exp n2 zero_const in
1599       let div_th, n_tm, e_tm =
1600         if flag then
1601           let th = num_exp_op_lo pp num_exp_div num_exp1 num_exp2 in
1602           let n_tm, e_tm = dest_num_exp (lhand(concl th)) in
1603             th, n_tm, e_tm
1604         else
1605           let th = num_exp_op_hi_lt pp num_exp_div num_exp1 num_exp2 in
1606           let n_tm, e_tm = dest_num_exp (rand(concl th)) in
1607             th, n_tm, e_tm in
1608
1609       let r1_th = raw_add_conv_hash (mk_binop add_op_num e2 k_tm) in
1610       let r1_tm = rand(concl r1_th) in
1611       let e_plus_e1 = raw_add_conv_hash (mk_binop add_op_num e_tm e1) in
1612       let ltm, rtm = dest_comb(concl e_plus_e1) in
1613       let r2_th' = raw_add_conv_hash (mk_binop add_op_num min_exp_num_const rtm) in
1614       let r2_th = TRANS (AP_TERM (mk_comb (add_op_num, min_exp_num_const)) e_plus_e1) r2_th' in
1615       let r2_tm = rand(concl r2_th) in
1616       let sub_th, le_th = raw_sub_and_le_hash_conv r2_tm r1_tm in
1617         if rand(concl le_th) <> r2_tm then
1618           failwith "float_div_lo: underflow"
1619         else
1620           let r_tm = rand(concl sub_th) in
1621           let n2_not_zero = raw_eq0_hash_conv n2 in
1622           let inst = INST[r1_tm, r1_var_num; r2_tm, r2_var_num;
1623                           n1, n1_var_num; e1, e1_var_num;
1624                           e_tm, e_var_num; k_tm, k_var_num;
1625                           n2, n2_var_num; e2, e2_var_num;
1626                           n_tm, n_var_num; r_tm, r_var_num] in
1627           let th0 = inst
1628             (if flag then
1629                if s1 = "F" then FLOAT_DIV_FF_lo' else FLOAT_DIV_TT_lo'
1630              else
1631                if s1 = "F" then FLOAT_DIV_FT_lo' else FLOAT_DIV_TF_lo') in
1632           let th1 = MY_PROVE_HYP n2_not_zero (MY_PROVE_HYP div_th (MY_PROVE_HYP le_th th0)) in
1633             MY_PROVE_HYP sub_th (MY_PROVE_HYP r2_th (MY_PROVE_HYP r1_th th1));;
1634
1635
1636 let float_div_hi pp f1 f2 =
1637   let s1, n1, e1 = dest_float f1 and
1638       s2, n2, e2 = dest_float f2 in
1639   let n1_eq0_th = raw_eq0_hash_conv n1 in
1640     if (rand o concl) n1_eq0_th = t_const then
1641       (MY_PROVE_HYP n1_eq0_th o
1642          INST[e1, e1_var_num; f2, f2_var_real; n1, n1_var_num;
1643               (if s1 = "T" then t_const else f_const), s1_var_bool]) FLOAT_DIV_0x_hi'
1644     else
1645       let flag = s1 = s2 in
1646
1647       let k_tm = rand (mk_small_numeral_array (2 * pp)) in
1648       let num_exp1 = mk_num_exp n1 k_tm and
1649           num_exp2 = mk_num_exp n2 zero_const in
1650       let div_th, n_tm, e_tm =
1651         if flag then
1652           let th = num_exp_op_hi_lt pp num_exp_div num_exp1 num_exp2 in
1653           let n_tm, e_tm = dest_num_exp (rand(concl th)) in
1654             th, n_tm, e_tm
1655         else
1656           let th = num_exp_op_lo pp num_exp_div num_exp1 num_exp2 in
1657           let n_tm, e_tm = dest_num_exp (lhand(concl th)) in
1658             th, n_tm, e_tm in
1659
1660       let r1_th = raw_add_conv_hash (mk_binop add_op_num e2 k_tm) in
1661       let r1_tm = rand(concl r1_th) in
1662       let e_plus_e1 = raw_add_conv_hash (mk_binop add_op_num e_tm e1) in
1663       let ltm, rtm = dest_comb(concl e_plus_e1) in
1664       let r2_th' = raw_add_conv_hash (mk_binop add_op_num min_exp_num_const rtm) in
1665       let r2_th = TRANS (AP_TERM (mk_comb (add_op_num, min_exp_num_const)) e_plus_e1) r2_th' in
1666       let r2_tm = rand(concl r2_th) in
1667       let sub_th, le_th = raw_sub_and_le_hash_conv r2_tm r1_tm in
1668         if rand(concl le_th) <> r2_tm then
1669           failwith "float_div_hi: underflow"
1670         else
1671           let r_tm = rand(concl sub_th) in
1672           let n2_not_zero = raw_eq0_hash_conv n2 in
1673           let inst = INST[r1_tm, r1_var_num; r2_tm, r2_var_num;
1674                           n1, n1_var_num; e1, e1_var_num;
1675                           e_tm, e_var_num; k_tm, k_var_num;
1676                           n2, n2_var_num; e2, e2_var_num;
1677                           n_tm, n_var_num; r_tm, r_var_num] in
1678           let th0 = inst
1679             (if flag then
1680                if s1 = "F" then FLOAT_DIV_FF_hi' else FLOAT_DIV_TT_hi'
1681              else
1682                if s1 = "F" then FLOAT_DIV_FT_hi' else FLOAT_DIV_TF_hi') in
1683           let th1 = MY_PROVE_HYP n2_not_zero (MY_PROVE_HYP div_th (MY_PROVE_HYP le_th th0)) in
1684             MY_PROVE_HYP sub_th (MY_PROVE_HYP r2_th (MY_PROVE_HYP r1_th th1));;
1685
1686
1687 (***********************************)
1688 (* FLOAT_ADD *)
1689
1690 let FLOAT_ADD_FF = prove(`num_exp n1 e1 + num_exp n2 e2 = num_exp n e
1691     ==> float_num F n1 e1 + float_num F n2 e2 = float_num F n e`,
1692    REPEAT STRIP_TAC THEN
1693      REWRITE_TAC[float; REAL_MUL_LID] THEN
1694      REWRITE_TAC[REAL_ARITH `a / b + c / b = (a + c) / b`] THEN
1695      ASM_REWRITE_TAC[REAL_OF_NUM_ADD]);;
1696
1697 let FLOAT_ADD_TT = prove(`num_exp n1 e1 + num_exp n2 e2 = num_exp n e
1698     ==> float_num T n1 e1 + float_num T n2 e2 = float_num T n e`,
1699    REWRITE_TAC[FLOAT_NEG_T; REAL_ARITH `--a + --b = --c <=> a + b = c`] THEN
1700      REWRITE_TAC[FLOAT_ADD_FF]);;
1701
1702 let FLOAT_ADD_FF_lo = prove(`num_exp n e <= num_exp n1 e1 + num_exp n2 e2
1703                               ==> float_num F n e <= float_num F n1 e1 + float_num F n2 e2`,
1704    REWRITE_TAC[GSYM REAL_OF_NUM_LE; GSYM REAL_OF_NUM_ADD] THEN
1705      REPEAT STRIP_TAC THEN
1706      MAP_EVERY ABBREV_TAC [`z = &(num_exp n e)`; `x = &(num_exp n1 e1)`; `y = &(num_exp n2 e2)`] THEN
1707      ASM_REWRITE_TAC[float; REAL_MUL_LID] THEN
1708      REWRITE_TAC[REAL_ARITH `a / b + c / b = (a + c) / b`] THEN
1709      REWRITE_TAC[real_div] THEN
1710      MATCH_MP_TAC REAL_LE_RMUL THEN
1711      ASM_REWRITE_TAC[REAL_LE_INV_EQ; REAL_POS]);;
1712
1713 let FLOAT_ADD_FF_hi = prove(`num_exp n1 e1 + num_exp n2 e2 <= num_exp n e
1714                               ==> float_num F n1 e1 + float_num F n2 e2 <= float_num F n e`,
1715    REWRITE_TAC[GSYM REAL_OF_NUM_LE; GSYM REAL_OF_NUM_ADD] THEN
1716      REPEAT STRIP_TAC THEN
1717      MAP_EVERY ABBREV_TAC [`z = &(num_exp n e)`; `x = &(num_exp n1 e1)`; `y = &(num_exp n2 e2)`] THEN
1718      ASM_REWRITE_TAC[float; REAL_MUL_LID] THEN
1719      REWRITE_TAC[REAL_ARITH `a / b + c / b = (a + c) / b`] THEN
1720      REWRITE_TAC[real_div] THEN
1721      MATCH_MP_TAC REAL_LE_RMUL THEN
1722      ASM_REWRITE_TAC[REAL_LE_INV_EQ; REAL_POS]);;
1723
1724 let FLOAT_ADD_TT_lo = prove(`num_exp n1 e1 + num_exp n2 e2 <= num_exp n e
1725                               ==> float_num T n e <= float_num T n1 e1 + float_num T n2 e2`,
1726    REWRITE_TAC[FLOAT_NEG_T; REAL_ARITH `--a <= --b + --c <=> b + c <= a`] THEN
1727      REWRITE_TAC[FLOAT_ADD_FF_hi]);;
1728
1729 let FLOAT_ADD_TT_hi = prove(`num_exp n e <= num_exp n1 e1 + num_exp n2 e2
1730                               ==> float_num T n1 e1 + float_num T n2 e2 <= float_num T n e`,
1731    REWRITE_TAC[FLOAT_NEG_T; REAL_ARITH `--b + --c <= --a <=> a <= b + c`] THEN
1732      REWRITE_TAC[FLOAT_ADD_FF_lo]);;
1733
1734 let FLOAT_ADD_FT_F_lo = prove(`num_exp n2 e2 <= num_exp n1 e1 ==>
1735                                 num_exp n e <= num_exp n1 e1 - num_exp n2 e2
1736                                 ==> float_num F n e <= float_num F n1 e1 + float_num T n2 e2`,
1737    MAP_EVERY ABBREV_TAC[`z = num_exp n e`; `x = num_exp n1 e1`; `y = num_exp n2 e2`] THEN
1738      ASM_REWRITE_TAC[FLOAT_NEG_T; float; REAL_MUL_LID] THEN
1739      DISCH_TAC THEN
1740      ASM_SIMP_TAC[GSYM REAL_OF_NUM_LE; GSYM REAL_OF_NUM_SUB] THEN
1741      REWRITE_TAC[num_exp; min_exp_def; MULT_CLAUSES; GSYM REAL_OF_NUM_POW] THEN
1742      REAL_ARITH_TAC);;
1743
1744 let FLOAT_ADD_FT_T_lo = prove(`num_exp n1 e1 <= num_exp n2 e2 ==>
1745                                 num_exp n2 e2 - num_exp n1 e1 <= num_exp n e
1746                                 ==> float_num T n e <= float_num F n1 e1 + float_num T n2 e2`,
1747    MAP_EVERY ABBREV_TAC[`z = num_exp n e`; `x = num_exp n1 e1`; `y = num_exp n2 e2`] THEN
1748      ASM_REWRITE_TAC[FLOAT_NEG_T; float; REAL_MUL_LID] THEN
1749      DISCH_TAC THEN
1750      ASM_SIMP_TAC[GSYM REAL_OF_NUM_LE; GSYM REAL_OF_NUM_SUB] THEN
1751      REWRITE_TAC[num_exp; min_exp_def; MULT_CLAUSES; GSYM REAL_OF_NUM_POW] THEN
1752      REAL_ARITH_TAC);;
1753
1754 let FLOAT_ADD_FT_F_hi = prove(`num_exp n2 e2 <= num_exp n1 e1 ==>
1755                                 num_exp n1 e1 - num_exp n2 e2 <= num_exp n e
1756                                 ==> float_num F n1 e1 + float_num T n2 e2 <= float_num F n e`,
1757    REWRITE_TAC[FLOAT_NEG_T; REAL_ARITH `a + --b <= c <=> --c <= b + --a`] THEN
1758      REWRITE_TAC[GSYM FLOAT_NEG_T; FLOAT_ADD_FT_T_lo]);;
1759
1760 let FLOAT_ADD_FT_T_hi = prove(`num_exp n1 e1 <= num_exp n2 e2 ==>
1761                                 num_exp n e <= num_exp n2 e2 - num_exp n1 e1
1762                                 ==> float_num F n1 e1 + float_num T n2 e2 <= float_num T n e`,
1763    REWRITE_TAC[FLOAT_NEG_T; REAL_ARITH `a + --b <= --c <=> c <= b + --a`] THEN
1764      REWRITE_TAC[GSYM FLOAT_NEG_T; FLOAT_ADD_FT_F_lo]);;
1765
1766
1767 (******************************************)
1768 (* float_add_lo, float_add_hi *)
1769
1770 let REAL_ADD_COMM = CONJUNCT1 REAL_ADD_AC;;
1771
1772 let transform = UNDISCH_ALL o REWRITE_RULE[GSYM IMP_IMP] o NUMERALS_TO_NUM;;
1773 let FLOAT_ADD_FF_hi' = transform FLOAT_ADD_FF_hi and
1774     FLOAT_ADD_FF_lo' = transform FLOAT_ADD_FF_lo and
1775     FLOAT_ADD_TT_hi' = transform FLOAT_ADD_TT_hi and
1776     FLOAT_ADD_TT_lo' = transform FLOAT_ADD_TT_lo and
1777     FLOAT_ADD_FT_F_lo' = transform FLOAT_ADD_FT_F_lo and
1778     FLOAT_ADD_FT_T_lo' = transform FLOAT_ADD_FT_T_lo and
1779     FLOAT_ADD_FT_F_hi' = transform FLOAT_ADD_FT_F_hi and
1780     FLOAT_ADD_FT_T_hi' = transform FLOAT_ADD_FT_T_hi;;
1781
1782 let float_add_lo pp f1 f2 =
1783   let s1, n1, e1 = dest_float f1 in
1784   let s2, n2, e2 = dest_float f2 in
1785     if s1 = s2 then
1786       let num_exp1 = mk_num_exp n1 e1 in
1787       let num_exp2 = mk_num_exp n2 e2 in
1788
1789         if s1 = "F" then
1790           (* F + F *)
1791           let add_th = num_exp_op_lo pp num_exp_add num_exp1 num_exp2 in
1792           let n_tm, e_tm = dest_num_exp (lhand(concl add_th)) in
1793           let th0 = INST[e_tm, e_var_num; n_tm, n_var_num; n1, n1_var_num;
1794                          e1, e1_var_num; n2, n2_var_num; e2, e2_var_num] FLOAT_ADD_FF_lo' in
1795             MY_PROVE_HYP add_th th0
1796         else
1797           (* T + T *)
1798           let add_th = num_exp_op_hi pp num_exp_add num_exp1 num_exp2 in
1799           let n_tm, e_tm = dest_num_exp (rand(concl add_th)) in
1800           let th0 = INST[e_tm, e_var_num; n_tm, n_var_num; n1, n1_var_num;
1801                          e1, e1_var_num; n2, n2_var_num; e2, e2_var_num] FLOAT_ADD_TT_lo' in
1802             MY_PROVE_HYP add_th th0
1803     else
1804       (* F + T or T + F *)
1805       let th0, n1, e1, n2, e2 =
1806         if s1 = "T" then
1807           INST[f2, m_var_real; f1, n_var_real] REAL_ADD_COMM, n2, e2, n1, e1
1808         else
1809           REFL(mk_binop add_op_real f1 f2), n1, e1, n2, e2 in
1810
1811       let num_exp1 = mk_num_exp n1 e1 in
1812       let num_exp2 = mk_num_exp n2 e2 in
1813
1814       let sub_th, le_th = num_exp_sub num_exp1 num_exp2 in
1815       let sub_tm = rand(concl sub_th) in
1816
1817         if rand(concl le_th) = num_exp1 then
1818           let lo_th = num_exp_lo pp sub_tm in
1819           let n_tm, e_tm = dest_num_exp (lhand(concl lo_th)) in
1820           let lo_sub_th = EQ_MP (AP_TERM (rator(concl lo_th)) (SYM sub_th)) lo_th in
1821
1822           let th1 = INST[n1, n1_var_num; e1, e1_var_num; n2, n2_var_num; e2, e2_var_num;
1823                          n_tm, n_var_num; e_tm, e_var_num] FLOAT_ADD_FT_F_lo' in
1824           let th2 = MY_PROVE_HYP lo_sub_th (MY_PROVE_HYP le_th th1) in
1825             EQ_MP (AP_TERM (rator(concl th2)) th0) th2
1826
1827         else
1828           let hi_th = num_exp_hi pp sub_tm in
1829           let n_tm, e_tm = dest_num_exp(rand(concl hi_th)) in
1830           let hi_sub_th = EQ_MP (SYM (AP_THM (AP_TERM le_op_num sub_th) (rand(concl hi_th)))) hi_th in
1831
1832           let th1 = INST[n1, n1_var_num; e1, e1_var_num; n2, n2_var_num; e2, e2_var_num;
1833                          n_tm, n_var_num; e_tm, e_var_num] FLOAT_ADD_FT_T_lo' in
1834           let th2 = MY_PROVE_HYP hi_sub_th (MY_PROVE_HYP le_th th1) in
1835             EQ_MP (AP_TERM (rator(concl th2)) th0) th2;;
1836
1837
1838 let float_add_hi pp f1 f2 =
1839   let s1, n1, e1 = dest_float f1 in
1840   let s2, n2, e2 = dest_float f2 in
1841     if s1 = s2 then
1842       let num_exp1 = mk_num_exp n1 e1 in
1843       let num_exp2 = mk_num_exp n2 e2 in
1844
1845         if s1 = "F" then
1846           (* F + F *)
1847           let add_th = num_exp_op_hi pp num_exp_add num_exp1 num_exp2 in
1848           let n_tm, e_tm = dest_num_exp (rand(concl add_th)) in
1849           let th0 = INST[e_tm, e_var_num; n_tm, n_var_num; n1, n1_var_num;
1850                          e1, e1_var_num; n2, n2_var_num; e2, e2_var_num] FLOAT_ADD_FF_hi' in
1851             MY_PROVE_HYP add_th th0
1852         else
1853           (* T + T *)
1854           let add_th = num_exp_op_lo pp num_exp_add num_exp1 num_exp2 in
1855           let n_tm, e_tm = dest_num_exp (lhand(concl add_th)) in
1856           let th0 = INST[e_tm, e_var_num; n_tm, n_var_num; n1, n1_var_num;
1857                          e1, e1_var_num; n2, n2_var_num; e2, e2_var_num] FLOAT_ADD_TT_hi' in
1858             MY_PROVE_HYP add_th th0
1859     else
1860       (* F + T or T + F *)
1861       let th0, n1, e1, n2, e2 =
1862         if s1 = "T" then
1863           INST[f2, m_var_real; f1, n_var_real] REAL_ADD_COMM, n2, e2, n1, e1
1864         else
1865           REFL(mk_binop add_op_real f1 f2), n1, e1, n2, e2 in
1866
1867       let num_exp1 = mk_num_exp n1 e1 in
1868       let num_exp2 = mk_num_exp n2 e2 in
1869
1870       let sub_th, le_th = num_exp_sub num_exp1 num_exp2 in
1871       let sub_tm = rand(concl sub_th) in
1872
1873         if rand(concl le_th) = num_exp1 then
1874           let hi_th = num_exp_hi pp sub_tm in
1875           let n_tm, e_tm = dest_num_exp (rand(concl hi_th)) in
1876           let hi_sub_th = EQ_MP (SYM (AP_THM (AP_TERM le_op_num sub_th) (rand(concl hi_th)))) hi_th in
1877
1878           let th1 = INST[n1, n1_var_num; e1, e1_var_num; n2, n2_var_num; e2, e2_var_num;
1879                          n_tm, n_var_num; e_tm, e_var_num] FLOAT_ADD_FT_F_hi' in
1880           let th2 = MY_PROVE_HYP hi_sub_th (MY_PROVE_HYP le_th th1) in
1881             EQ_MP (AP_THM (AP_TERM le_op_real th0) (rand(concl th2))) th2
1882
1883         else
1884           let lo_th = num_exp_lo pp sub_tm in
1885           let n_tm, e_tm = dest_num_exp(lhand(concl lo_th)) in
1886           let lo_sub_th = EQ_MP (AP_TERM (rator(concl lo_th)) (SYM sub_th)) lo_th in
1887
1888           let th1 = INST[n1, n1_var_num; e1, e1_var_num; n2, n2_var_num; e2, e2_var_num;
1889                          n_tm, n_var_num; e_tm, e_var_num] FLOAT_ADD_FT_T_hi' in
1890           let th2 = MY_PROVE_HYP lo_sub_th (MY_PROVE_HYP le_th th1) in
1891             EQ_MP (AP_THM (AP_TERM le_op_real th0) (rand(concl th2))) th2;;
1892
1893
1894 (******************************************)
1895 (* float_sub_lo, float_sub_hi *)
1896
1897 let FLOAT_SUB_F_EQ_ADD = (SYM o prove)(`f1 - float_num F n2 e2 = f1 + float_num T n2 e2`,
1898                                        REWRITE_TAC[FLOAT_NEG_T] THEN REAL_ARITH_TAC);;
1899
1900 let FLOAT_SUB_T_EQ_ADD = (SYM o prove)(`f1 - float_num T n2 e2 = f1 + float_num F n2 e2`,
1901                                        REWRITE_TAC[FLOAT_NEG_T] THEN REAL_ARITH_TAC);;
1902
1903 let float_sub_lo pp f1 f2 =
1904   let s2, n2, e2 = dest_float f2 in
1905   let th0 =
1906     INST[f1, f1_var_real; n2, n2_var_num; e2, e2_var_num]
1907       (if s2 = "F" then FLOAT_SUB_F_EQ_ADD else FLOAT_SUB_T_EQ_ADD) in
1908   let ltm,f2_tm = dest_comb(lhand(concl th0)) in
1909   let f1_tm = rand ltm in
1910   let lo_th = float_add_lo pp f1_tm f2_tm in
1911     EQ_MP (AP_TERM (rator(concl lo_th)) th0) lo_th;;
1912
1913
1914 let float_sub_hi pp f1 f2 =
1915   let s2, n2, e2 = dest_float f2 in
1916   let th0 =
1917     INST[f1, f1_var_real; n2, n2_var_num; e2, e2_var_num]
1918       (if s2 = "F" then FLOAT_SUB_F_EQ_ADD else FLOAT_SUB_T_EQ_ADD) in
1919   let ltm, f2_tm = dest_comb(lhand(concl th0)) in
1920   let f1_tm = rand ltm in
1921   let hi_th = float_add_hi pp f1_tm f2_tm in
1922     EQ_MP (AP_THM (AP_TERM le_op_real th0) (rand(concl hi_th))) hi_th;;
1923
1924
1925 (*******************************************)
1926 (* FLOAT_SQRT *)
1927
1928 (* float_num F m e = float_num F (B0 m) (PRE e) *)
1929 let FLOAT_PRE_EXP = prove(mk_imp(`~(e = 0) /\ PRE e = e1`,
1930                                  mk_eq(`float_num F m e`,
1931                                        mk_comb(mk_comb(`float_num F`, mk_comb(b0_const, m_var_num)), `e1:num`))),
1932    STRIP_TAC THEN POP_ASSUM (fun th -> REWRITE_TAC[SYM th]) THEN
1933      REWRITE_TAC[float; REAL_MUL_LID; real_div; REAL_EQ_MUL_RCANCEL] THEN
1934      DISJ1_TAC THEN
1935      REWRITE_TAC[num_exp; b0_thm; REAL_OF_NUM_EQ] THEN
1936      SUBGOAL_THEN `e = SUC (PRE e)` MP_TAC THENL
1937      [
1938        POP_ASSUM MP_TAC THEN ARITH_TAC;
1939        ALL_TAC
1940      ] THEN
1941      DISCH_THEN (fun th -> CONV_TAC (LAND_CONV (ONCE_REWRITE_CONV[th]))) THEN
1942      REWRITE_TAC[EXP] THEN
1943      ARITH_TAC);;
1944
1945 let DIV2_EVEN_lemma = prove(`!n. EVEN n ==> 2 * (n DIV 2) = n`,
1946    GEN_TAC THEN
1947      REWRITE_TAC[EVEN_EXISTS] THEN
1948      STRIP_TAC THEN
1949      ASM_REWRITE_TAC[] THEN
1950      MATCH_MP_TAC (ARITH_RULE `x = y ==> 2 * x = 2 * y`) THEN
1951      MATCH_MP_TAC DIV_MULT THEN
1952      ARITH_TAC);;
1953
1954 let FLOAT_SQRT_EVEN_lo = prove(`f1 * f1 = f2 /\ f2 <= x /\
1955                                 num_exp m (2 * p) = x /\ f1 = num_exp n1 e1
1956                                    /\ EVEN e /\ e DIV 2 = e2 /\
1957                                    e1 + e2 + (min_exp DIV 2) = r /\
1958                                    p <= r /\ r - p = r2
1959                                        ==> float_num F n1 r2 <= sqrt (float_num F m e)`,
1960    STRIP_TAC THEN
1961      UNDISCH_TAC `f2 <= x:num` THEN
1962      UNDISCH_TAC `num_exp m (2 * p) = x` THEN
1963      DISCH_THEN (fun th -> REWRITE_TAC[SYM th]) THEN
1964      UNDISCH_TAC `f1 * f1 = f2:num` THEN
1965      DISCH_THEN (fun th -> REWRITE_TAC[SYM th]) THEN
1966      UNDISCH_TAC `e1 + e2 + min_exp DIV 2 = r` THEN
1967      UNDISCH_TAC `e DIV 2 = e2` THEN
1968      DISCH_THEN (fun th -> REWRITE_TAC[SYM th]) THEN
1969      REPEAT (POP_ASSUM MP_TAC) THEN
1970      REWRITE_TAC[num_exp; float; REAL_MUL_LID; GSYM REAL_OF_NUM_MUL] THEN
1971      REPEAT STRIP_TAC THEN
1972      UNDISCH_TAC `r - p = r2:num` THEN
1973      DISCH_THEN (fun th -> REWRITE_TAC[SYM th]) THEN
1974
1975      MATCH_MP_TAC REAL_LE_RSQRT THEN
1976      REWRITE_TAC[GSYM REAL_OF_NUM_POW; REAL_POW_DIV] THEN
1977      REWRITE_TAC[REAL_POW_2; real_div; REAL_INV_MUL] THEN
1978      REWRITE_TAC[REAL_MUL_ASSOC] THEN
1979      MATCH_MP_TAC REAL_LE_RMUL THEN
1980      CONJ_TAC THENL
1981      [
1982        REWRITE_TAC[REAL_ARITH `(((a * b) * a) * b) * c = (a * a) * (b * b) * c:real`] THEN
1983          REWRITE_TAC[GSYM REAL_POW_ADD] THEN
1984          REWRITE_TAC[ARITH_RULE `r - p + r - p = 2 * r - 2 * p`] THEN
1985          MP_TAC (SPECL[mk_comb(amp_op_real, base_const); `2 * r`; `2 * p`] REAL_DIV_POW2) THEN
1986          ANTS_TAC THENL [ REAL_ARITH_TAC; ALL_TAC ] THEN
1987          ASM_SIMP_TAC[ARITH_RULE `p <= r ==> 2 * p <= 2 * r`] THEN
1988          DISCH_THEN (fun th -> REWRITE_TAC[SYM th; real_div]) THEN
1989
1990          SUBGOAL_THEN `2 * r = (e1 + e1) + min_exp + e` (fun th -> REWRITE_TAC[th]) THENL
1991          [
1992            EXPAND_TAC "r" THEN
1993              REWRITE_TAC[ARITH_RULE `2 * (e1 + b + c) = (e1 + e1) + 2 * c + 2 * b`] THEN
1994              MATCH_MP_TAC (ARITH_RULE `b1 = b2 /\ c1 = c2 ==> a + b1 + c1 = a + b2 + c2:num`) THEN
1995              SUBGOAL_THEN `EVEN min_exp` ASSUME_TAC THENL
1996              [
1997                REWRITE_TAC[min_exp_def] THEN ARITH_TAC;
1998                ALL_TAC
1999              ] THEN
2000              ASM_SIMP_TAC[DIV2_EVEN_lemma];
2001            ALL_TAC
2002          ] THEN
2003
2004          REWRITE_TAC[REAL_POW_ADD] THEN
2005          REWRITE_TAC[REAL_ARITH `(n * n) * (((e * e) * x * y) * z) * u = (n * e) * (n * e) * (x * u) * z * y:real`] THEN
2006          SUBGOAL_THEN `~(&(num_exp 1 min_exp) = &0)` MP_TAC THENL
2007          [
2008            REWRITE_TAC[REAL_OF_NUM_EQ; NUM_EXP_EQ_0] THEN ARITH_TAC;
2009            ALL_TAC
2010          ] THEN
2011          REWRITE_TAC[num_exp; REAL_MUL_LID; GSYM REAL_OF_NUM_POW; GSYM REAL_OF_NUM_MUL] THEN
2012          DISCH_THEN (fun th -> SIMP_TAC[th; REAL_MUL_RINV; REAL_MUL_LID]) THEN
2013          FIRST_X_ASSUM (MP_TAC o check(fun th -> (fst o dest_var o lhand o concl) th = "f1")) THEN
2014          REWRITE_TAC[GSYM REAL_OF_NUM_EQ; GSYM REAL_OF_NUM_MUL; GSYM REAL_OF_NUM_POW] THEN
2015          DISCH_THEN (fun th -> REWRITE_TAC[SYM th]) THEN
2016
2017          REWRITE_TAC[REAL_MUL_ASSOC] THEN
2018          MATCH_MP_TAC REAL_LE_RMUL THEN
2019          CONJ_TAC THENL
2020          [
2021            SUBGOAL_THEN `!x y z. &0 < x /\ y <= z * x ==> y * inv x <= z` MP_TAC THENL
2022              [
2023                REPEAT STRIP_TAC THEN
2024                  MP_TAC (SPECL [`y * inv x`; `z:real`; `x:real`] REAL_LE_RMUL_EQ) THEN
2025                  ASM_REWRITE_TAC[] THEN
2026                  DISCH_THEN (fun th -> REWRITE_TAC[SYM th; GSYM REAL_MUL_ASSOC]) THEN
2027                  ASM_SIMP_TAC[REAL_ARITH `&0 < x ==> ~(x = &0)`; REAL_MUL_LINV; REAL_MUL_RID];
2028                ALL_TAC
2029              ] THEN
2030
2031              DISCH_THEN (MP_TAC o SPECL[`&(num_exp 1 (2 * p))`; `&(f1 * f1)`; `&m`]) THEN
2032              REWRITE_TAC[num_exp; GSYM REAL_OF_NUM_MUL; GSYM REAL_OF_NUM_POW; REAL_MUL_LID] THEN
2033              DISCH_THEN MATCH_MP_TAC THEN
2034              ASM_REWRITE_TAC[REAL_OF_NUM_MUL; REAL_OF_NUM_POW; REAL_OF_NUM_LE] THEN
2035              REWRITE_TAC[REAL_OF_NUM_LT; EXP_LT_0] THEN
2036              ARITH_TAC;
2037            ALL_TAC
2038          ] THEN
2039
2040          MATCH_MP_TAC REAL_POW_LE THEN
2041          ARITH_TAC;
2042        ALL_TAC
2043      ] THEN
2044
2045      MATCH_MP_TAC REAL_LE_INV THEN
2046      MATCH_MP_TAC REAL_POW_LE THEN
2047      ARITH_TAC);;
2048
2049
2050 let FLOAT_SQRT_EVEN_hi = prove(`f1 * f1 = f2 /\ x <= f2 /\
2051                                 num_exp m (2 * p) = x /\ f1 = num_exp n1 e1
2052                                    /\ EVEN e /\ e DIV 2 = e2 /\
2053                                    e1 + e2 + (min_exp DIV 2) = r /\
2054                                    p <= r /\ r - p = r2
2055                                        ==> sqrt (float_num F m e) <= float_num F n1 r2`,
2056    STRIP_TAC THEN
2057      UNDISCH_TAC `x <= f2:num` THEN
2058      UNDISCH_TAC `num_exp m (2 * p) = x` THEN
2059      DISCH_THEN (fun th -> REWRITE_TAC[SYM th]) THEN
2060      UNDISCH_TAC `f1 * f1 = f2:num` THEN
2061      DISCH_THEN (fun th -> REWRITE_TAC[SYM th]) THEN
2062      UNDISCH_TAC `e1 + e2 + min_exp DIV 2 = r` THEN
2063      UNDISCH_TAC `e DIV 2 = e2` THEN
2064      DISCH_THEN (fun th -> REWRITE_TAC[SYM th]) THEN
2065      REPEAT (POP_ASSUM MP_TAC) THEN
2066      REWRITE_TAC[num_exp; float; REAL_MUL_LID; GSYM REAL_OF_NUM_MUL] THEN
2067      REPEAT STRIP_TAC THEN
2068      UNDISCH_TAC `r - p = r2:num` THEN
2069      DISCH_THEN (fun th -> REWRITE_TAC[SYM th]) THEN
2070
2071      MATCH_MP_TAC REAL_LE_LSQRT THEN
2072      REPEAT CONJ_TAC THENL
2073      [
2074        REWRITE_TAC[REAL_OF_NUM_MUL; real_div] THEN
2075          MATCH_MP_TAC REAL_LE_MUL THEN
2076          REWRITE_TAC[REAL_POS] THEN
2077          MATCH_MP_TAC REAL_LE_INV THEN
2078          REWRITE_TAC[REAL_POS];
2079        REWRITE_TAC[REAL_OF_NUM_MUL; real_div] THEN
2080          MATCH_MP_TAC REAL_LE_MUL THEN
2081          REWRITE_TAC[REAL_POS] THEN
2082          MATCH_MP_TAC REAL_LE_INV THEN
2083          REWRITE_TAC[REAL_POS];
2084        ALL_TAC
2085      ] THEN
2086
2087      REWRITE_TAC[GSYM REAL_OF_NUM_POW; REAL_POW_DIV] THEN
2088      REWRITE_TAC[REAL_POW_2; real_div; REAL_INV_MUL] THEN
2089      REWRITE_TAC[REAL_MUL_ASSOC] THEN
2090      MATCH_MP_TAC REAL_LE_RMUL THEN
2091      CONJ_TAC THENL
2092      [
2093        REWRITE_TAC[REAL_ARITH `(((a * b) * a) * b) * c = (a * a) * (b * b) * c:real`] THEN
2094          REWRITE_TAC[GSYM REAL_POW_ADD] THEN
2095          REWRITE_TAC[ARITH_RULE `r - p + r - p = 2 * r - 2 * p`] THEN
2096          MP_TAC (SPECL[mk_comb(amp_op_real, base_const); `2 * r`; `2 * p`] REAL_DIV_POW2) THEN
2097          ANTS_TAC THENL [ REAL_ARITH_TAC; ALL_TAC ] THEN
2098          ASM_SIMP_TAC[ARITH_RULE `p <= r ==> 2 * p <= 2 * r`] THEN
2099          DISCH_THEN (fun th -> REWRITE_TAC[SYM th; real_div]) THEN
2100
2101          SUBGOAL_THEN `2 * r = (e1 + e1) + min_exp + e` (fun th -> REWRITE_TAC[th]) THENL
2102          [
2103            EXPAND_TAC "r" THEN
2104              REWRITE_TAC[ARITH_RULE `2 * (e1 + b + c) = (e1 + e1) + 2 * c + 2 * b`] THEN
2105              MATCH_MP_TAC (ARITH_RULE `b1 = b2 /\ c1 = c2 ==> a + b1 + c1 = a + b2 + c2:num`) THEN
2106              SUBGOAL_THEN `EVEN min_exp` ASSUME_TAC THENL
2107              [
2108                REWRITE_TAC[min_exp_def] THEN ARITH_TAC;
2109                ALL_TAC
2110              ] THEN
2111              ASM_SIMP_TAC[DIV2_EVEN_lemma];
2112            ALL_TAC
2113          ] THEN
2114
2115          REWRITE_TAC[REAL_POW_ADD] THEN
2116          REWRITE_TAC[REAL_ARITH `(n * n) * (((e * e) * x * y) * z) * u = (n * e) * (n * e) * (x * u) * z * y:real`] THEN
2117          SUBGOAL_THEN `~(&(num_exp 1 min_exp) = &0)` MP_TAC THENL
2118          [
2119            REWRITE_TAC[REAL_OF_NUM_EQ; NUM_EXP_EQ_0] THEN ARITH_TAC;
2120            ALL_TAC
2121          ] THEN
2122          REWRITE_TAC[num_exp; REAL_MUL_LID; GSYM REAL_OF_NUM_POW; GSYM REAL_OF_NUM_MUL] THEN
2123          DISCH_THEN (fun th -> SIMP_TAC[th; REAL_MUL_RINV; REAL_MUL_LID]) THEN
2124          FIRST_X_ASSUM (MP_TAC o check(fun th -> (fst o dest_var o lhand o concl) th = "f1")) THEN
2125          REWRITE_TAC[GSYM REAL_OF_NUM_EQ; GSYM REAL_OF_NUM_MUL; GSYM REAL_OF_NUM_POW] THEN
2126          DISCH_THEN (fun th -> REWRITE_TAC[SYM th]) THEN
2127
2128          REWRITE_TAC[REAL_MUL_ASSOC] THEN
2129          MATCH_MP_TAC REAL_LE_RMUL THEN
2130          CONJ_TAC THENL
2131          [
2132            SUBGOAL_THEN `!x y z. &0 < x /\ z * x <= y ==> z <= y * inv x` MP_TAC THENL
2133              [
2134                REPEAT STRIP_TAC THEN
2135                  MP_TAC (SPECL [`z:real`; `y * inv x`; `x:real`] REAL_LE_RMUL_EQ) THEN
2136                  ASM_REWRITE_TAC[] THEN
2137                  DISCH_THEN (fun th -> REWRITE_TAC[SYM th; GSYM REAL_MUL_ASSOC]) THEN
2138                  ASM_SIMP_TAC[REAL_ARITH `&0 < x ==> ~(x = &0)`; REAL_MUL_LINV; REAL_MUL_RID];
2139                ALL_TAC
2140              ] THEN
2141
2142              DISCH_THEN (MP_TAC o SPECL[`&(num_exp 1 (2 * p))`; `&(f1 * f1)`; `&m`]) THEN
2143              REWRITE_TAC[num_exp; GSYM REAL_OF_NUM_MUL; GSYM REAL_OF_NUM_POW; REAL_MUL_LID] THEN
2144              DISCH_THEN MATCH_MP_TAC THEN
2145              ASM_REWRITE_TAC[REAL_OF_NUM_MUL; REAL_OF_NUM_POW; REAL_OF_NUM_LE] THEN
2146              REWRITE_TAC[REAL_OF_NUM_LT; EXP_LT_0] THEN
2147              ARITH_TAC;
2148            ALL_TAC
2149          ] THEN
2150
2151          MATCH_MP_TAC REAL_POW_LE THEN
2152          ARITH_TAC;
2153        ALL_TAC
2154      ] THEN
2155
2156      MATCH_MP_TAC REAL_LE_INV THEN
2157      MATCH_MP_TAC REAL_POW_LE THEN
2158      ARITH_TAC);;
2159
2160
2161 (******************)
2162 let transform = UNDISCH_ALL o
2163   PURE_ONCE_REWRITE_RULE[TAUT `EVEN e <=> (EVEN e <=> T)`] o
2164   NUMERALS_TO_NUM o
2165   CONV_RULE (DEPTH_CONV NUM_DIV_CONV) o
2166   REWRITE_RULE[GSYM IMP_IMP; min_exp_def];;
2167
2168 let FLOAT_SQRT_EVEN_lo' = transform FLOAT_SQRT_EVEN_lo and
2169     FLOAT_SQRT_EVEN_hi' = transform FLOAT_SQRT_EVEN_hi and
2170     FLOAT_PRE_EXP' = (UNDISCH_ALL o
2171                         PURE_ONCE_REWRITE_RULE[TAUT `~(e = _0) <=> ((e = _0) <=> F)`] o
2172                         REWRITE_RULE[GSYM IMP_IMP; NUMERAL]) FLOAT_PRE_EXP;;
2173
2174 let even_const = `EVEN` and
2175     pre_const = `PRE` and
2176     two_num = rand(mk_small_numeral_array 2) and
2177     min_exp_div2 = rand(mk_small_numeral_array (min_exp / 2)) and
2178     f2_var_num = `f2:num` and
2179     f1_var_num = `f1:num` and
2180     p_var_num = `p:num`;;
2181
2182 (* Returns the list of digits of the given Big_int n in the base b *)
2183 let rec get_big_int_digits b n =
2184   let bb = big_int_of_int b in
2185   if le_big_int n zero_big_int then []
2186   else
2187     let q, r = quomod_big_int n bb in
2188       r :: get_big_int_digits b q;;
2189
2190 (* [1;2;3] -> 123 (base = 10) *)
2191 let rec big_int_from_list b list =
2192   let rec proc acc list =
2193     match list with
2194         [] -> acc
2195       | h::t -> proc (add_big_int h (mult_int_big_int b acc)) t in
2196     proc zero_big_int list;;
2197
2198 (* Returns n first elements of the list *)
2199 let rec take n list =
2200   match list with
2201       x :: xs -> if n > 0 then x :: take (n - 1) xs else []
2202     |   [] -> [];;
2203
2204 (* Returns an integer number that contains at most pp significant digits
2205    in the given base b *)
2206 let big_int_round_lo base pp n =
2207   let digits = rev (get_big_int_digits base n) in
2208   let n_digits = length digits in
2209     if n_digits <= pp then
2210       n
2211     else
2212       let m = big_int_from_list base (take pp digits) in
2213         mult_big_int (power_int_positive_int base (n_digits - pp)) m;;
2214
2215 let big_int_round_hi base pp n =
2216   let digits = rev (get_big_int_digits base n) in
2217   let n_digits = length digits in
2218     if n_digits <= pp then n
2219     else
2220       let l1, l2 = chop_list pp digits in
2221         if forall (eq_big_int zero_big_int) l2 then n
2222         else
2223           let m = succ_big_int (big_int_from_list base l1) in
2224             mult_big_int (power_int_positive_int base (n_digits - pp)) m;;
2225
2226
2227 (******************)
2228 let rec float_sqrt_lo pp tm =
2229   let s, m_tm, e_tm = dest_float tm in
2230   let p_tm = rand (mk_small_numeral_array pp) in
2231     if s <> "F" then
2232       failwith "float_sqrt_lo: negative argument"
2233     else
2234       let even_th = raw_even_hash_conv (mk_comb (even_const, e_tm)) in
2235         if (fst o dest_const o rand o concl) even_th <> "T" then
2236           (* ODD e *)
2237           let pre_e = raw_pre_hash_conv (mk_comb (pre_const, e_tm)) in
2238           let e_neq_0 = raw_eq0_hash_conv e_tm in
2239           let e1_tm = rand (concl pre_e) in
2240           let th0 = INST[e1_tm, e1_var_num; e_tm, e_var_num; m_tm, m_var_num] FLOAT_PRE_EXP' in
2241           let th1 = MY_PROVE_HYP pre_e (MY_PROVE_HYP e_neq_0 th0) in
2242           let th2 = float_sqrt_lo pp (rand(concl th1)) in
2243           let ltm, rtm = dest_comb (concl th2) in
2244             EQ_MP (SYM (AP_TERM ltm (AP_TERM (rator rtm) th1))) th2
2245         else
2246           (* EVEN e *)
2247           let p2_tm = mk_binop mul_op_num two_num p_tm in
2248           let p2_th = raw_mul_conv_hash p2_tm in
2249           let f1_1 = AP_TERM (mk_comb(num_exp_const, m_tm)) p2_th in
2250           let f1_2 = TRANS f1_1 (denormalize (rand (concl f1_1))) in
2251
2252           let x_tm = rand(concl f1_2) in
2253           let x = raw_dest_hash x_tm in
2254           let f1' = Big_int.sqrt_big_int (big_int_of_num x) in
2255           let f1 = num_of_big_int (big_int_round_lo Arith_hash.arith_base pp f1') in
2256           let f1_tm = rand(mk_numeral_array f1) in
2257           let f1_num_exp = to_num_exp f1_tm in
2258
2259           let n1_tm, e1_tm = dest_num_exp (rand (concl f1_num_exp)) in
2260           let f1f1_eq_f2 = raw_mul_conv_hash (mk_binop mul_op_num f1_tm f1_tm) in
2261           let f2_tm = rand(concl f1f1_eq_f2) in
2262           let f2_le_x = EQT_ELIM (raw_le_hash_conv (mk_binop le_op_num f2_tm x_tm)) in
2263
2264           let e_div2_eq_e2 = raw_div_hash_conv (mk_binop div_op_num e_tm two_num) in
2265           let e2_tm = rand(concl e_div2_eq_e2) in
2266           let r_th1 = raw_add_conv_hash (mk_binop add_op_num e2_tm min_exp_div2) in
2267           let r_th2 = AP_TERM (mk_comb(add_op_num, e1_tm)) r_th1 in
2268           let r_th = TRANS r_th2 (raw_add_conv_hash (rand (concl r_th2))) in
2269
2270           let r_tm = rand(concl r_th) in
2271           let r_sub_p, p_le_r = raw_sub_and_le_hash_conv p_tm r_tm in
2272           let r2_tm = rand(concl r_sub_p) in
2273             if (rand(concl p_le_r) <> r_tm) then
2274               failwith "float_sqrt_lo: underflow"
2275             else
2276               let th0 = INST[f2_tm, f2_var_num; x_tm, x_var_num; p_tm, p_var_num; r_tm, r_var_num;
2277                              f1_tm, f1_var_num; n1_tm, n1_var_num; e1_tm, e1_var_num; e2_tm, e2_var_num;
2278                              e_tm, e_var_num; m_tm, m_var_num; r2_tm, r2_var_num]
2279                         FLOAT_SQRT_EVEN_lo' in
2280                 MY_PROVE_HYP f1_2 (
2281                   MY_PROVE_HYP e_div2_eq_e2 (
2282                     MY_PROVE_HYP r_sub_p (
2283                       MY_PROVE_HYP r_th (
2284                         MY_PROVE_HYP f1f1_eq_f2 (
2285                           MY_PROVE_HYP f1_num_exp (
2286                             MY_PROVE_HYP even_th (
2287                               MY_PROVE_HYP f2_le_x (
2288                                 MY_PROVE_HYP p_le_r th0
2289                               ))))))));;
2290
2291
2292 let rec float_sqrt_hi pp tm =
2293   let s, m_tm, e_tm = dest_float tm in
2294   let p_tm = rand (mk_small_numeral_array pp) in
2295     if s <> "F" then
2296       failwith "float_sqrt_lo: negative argument"
2297     else
2298       let even_th = raw_even_hash_conv (mk_comb (even_const, e_tm)) in
2299         if (fst o dest_const o rand o concl) even_th <> "T" then
2300           (* ODD e *)
2301           let pre_e = raw_pre_hash_conv (mk_comb (pre_const, e_tm)) in
2302           let e_neq_0 = raw_eq0_hash_conv e_tm in
2303           let e1_tm = rand (concl pre_e) in
2304           let th0 = INST[e1_tm, e1_var_num; e_tm, e_var_num; m_tm, m_var_num] FLOAT_PRE_EXP' in
2305           let th1 = MY_PROVE_HYP pre_e (MY_PROVE_HYP e_neq_0 th0) in
2306           let th2 = float_sqrt_hi pp (rand(concl th1)) in
2307           let ltm, rtm = dest_comb (concl th2) in
2308           let ltm2, rtm2 = dest_comb ltm in
2309           let th3 = AP_THM (AP_TERM ltm2 (AP_TERM (rator rtm2) th1)) rtm in
2310             EQ_MP (SYM th3) th2
2311         else
2312           (* EVEN e *)
2313           let p2_tm = mk_binop mul_op_num two_num p_tm in
2314           let p2_th = raw_mul_conv_hash p2_tm in
2315           let f1_1 = AP_TERM (mk_comb(num_exp_const, m_tm)) p2_th in
2316           let f1_2 = TRANS f1_1 (denormalize (rand (concl f1_1))) in
2317
2318           let x_tm = rand(concl f1_2) in
2319           let x = raw_dest_hash x_tm in
2320           let x' = big_int_of_num x in
2321           let f1' = sqrt_big_int x' in
2322           let f1 = (num_of_big_int o big_int_round_hi Arith_hash.arith_base pp)
2323             (if eq_big_int (mult_big_int f1' f1') x' then f1' else succ_big_int f1') in
2324
2325           let f1_tm = rand(mk_numeral_array f1) in
2326           let f1_num_exp = to_num_exp f1_tm in
2327
2328           let n1_tm, e1_tm = dest_num_exp (rand (concl f1_num_exp)) in
2329           let f1f1_eq_f2 = raw_mul_conv_hash (mk_binop mul_op_num f1_tm f1_tm) in
2330           let f2_tm = rand(concl f1f1_eq_f2) in
2331           let x_le_f2 = EQT_ELIM (raw_le_hash_conv (mk_binop le_op_num x_tm f2_tm)) in
2332
2333           let e_div2_eq_e2 = raw_div_hash_conv (mk_binop div_op_num e_tm two_num) in
2334           let e2_tm = rand(concl e_div2_eq_e2) in
2335           let r_th1 = raw_add_conv_hash (mk_binop add_op_num e2_tm min_exp_div2) in
2336           let r_th2 = AP_TERM (mk_comb(add_op_num, e1_tm)) r_th1 in
2337           let r_th = TRANS r_th2 (raw_add_conv_hash (rand (concl r_th2))) in
2338
2339           let r_tm = rand(concl r_th) in
2340           let r_sub_p, p_le_r = raw_sub_and_le_hash_conv p_tm r_tm in
2341           let r2_tm = rand(concl r_sub_p) in
2342             if (rand(concl p_le_r) <> r_tm) then
2343               failwith "float_sqrt_lo: underflow"
2344             else
2345               let th0 = INST[f2_tm, f2_var_num; x_tm, x_var_num; p_tm, p_var_num; r_tm, r_var_num;
2346                              f1_tm, f1_var_num; n1_tm, n1_var_num; e1_tm, e1_var_num; e2_tm, e2_var_num;
2347                              e_tm, e_var_num; m_tm, m_var_num; r2_tm, r2_var_num]
2348                         FLOAT_SQRT_EVEN_hi' in
2349                 MY_PROVE_HYP f1_2 (
2350                   MY_PROVE_HYP e_div2_eq_e2 (
2351                     MY_PROVE_HYP r_sub_p (
2352                       MY_PROVE_HYP r_th (
2353                         MY_PROVE_HYP f1f1_eq_f2 (
2354                           MY_PROVE_HYP f1_num_exp (
2355                             MY_PROVE_HYP even_th (
2356                               MY_PROVE_HYP x_le_f2 (
2357                                 MY_PROVE_HYP p_le_r th0
2358                               ))))))));;
2359
2360 end;; (* Float_ops module *)
2361
2362 (************************************)
2363 (* Cached floating point operations *)
2364 (************************************)
2365
2366 (* Counters for collecting stats *)
2367 let lt0_c = ref 0 and
2368     gt0_c = ref 0 and
2369     lt_c = ref 0 and
2370     le0_c = ref 0 and
2371     ge0_c = ref 0 and
2372     le_c = ref 0 and
2373     min_c = ref 0 and
2374     max_c = ref 0 and
2375     min_max_c = ref 0 and
2376     mul_lo_c = ref 0 and
2377     mul_hi_c = ref 0 and
2378     div_lo_c = ref 0 and
2379     div_hi_c = ref 0 and
2380     add_lo_c = ref 0 and
2381     add_hi_c = ref 0 and
2382     sub_lo_c = ref 0 and
2383     sub_hi_c = ref 0 and
2384     sqrt_lo_c = ref 0 and
2385     sqrt_hi_c = ref 0;;
2386
2387 (* Hash tables *)
2388 let cache_size = if !Arith_options.float_cached then !Arith_options.init_cache_size else 1;;
2389
2390 let my_add h key v =
2391   if Hashtbl.length h >= !Arith_options.max_cache_size then
2392         Hashtbl.clear h
2393 (*    let _ = Hashtbl.clear h in
2394       print_string "Clearing a float hash table" *)
2395   else
2396     ();
2397   Hashtbl.add h key v;;
2398
2399 let mul_table = Hashtbl.create cache_size and
2400     div_table = Hashtbl.create cache_size and
2401     add_table = Hashtbl.create cache_size and
2402     sub_table = Hashtbl.create cache_size and
2403     sqrt_table = Hashtbl.create cache_size and
2404     le_table = Hashtbl.create cache_size and
2405     max_table = Hashtbl.create cache_size;;
2406
2407 let reset_cache () =
2408   Hashtbl.clear mul_table;
2409   Hashtbl.clear div_table;
2410   Hashtbl.clear add_table;
2411   Hashtbl.clear sub_table;
2412   Hashtbl.clear sqrt_table;
2413   Hashtbl.clear le_table;
2414   Hashtbl.clear max_table;;
2415
2416 let reset_stat () =
2417   lt0_c := 0;
2418   gt0_c := 0;
2419   lt_c := 0;
2420   le0_c := 0;
2421   ge0_c := 0;
2422   le_c := 0;
2423   min_c := 0;
2424   max_c := 0;
2425   min_max_c := 0;
2426   mul_lo_c := 0;
2427   mul_hi_c := 0;
2428   div_lo_c := 0;
2429   div_hi_c := 0;
2430   add_lo_c := 0;
2431   add_hi_c := 0;
2432   sub_lo_c := 0;
2433   sub_hi_c := 0;
2434   sqrt_lo_c := 0;
2435   sqrt_hi_c := 0;;
2436
2437 let print_stat () =
2438   let len = Hashtbl.length in
2439   let cmp_str1 = sprintf "lt0 = %d\ngt0 = %d\nlt = %d\n" !lt0_c !gt0_c !lt_c and
2440       cmp_str2 = sprintf "le0 = %d\nge0 = %d\n" !le0_c !ge0_c and
2441       cmp_str3 = sprintf "min = %d\nmin_max = %d\n" !min_c !min_max_c and
2442       le_str = sprintf "le = %d (le_hash = %d)\n" !le_c (len le_table) and
2443       max_str = sprintf "max = %d (max_hash = %d)\n" !max_c (len max_table) and
2444       mul_str = sprintf "mul_lo = %d, mul_hi = %d (mul_hash = %d)\n" !mul_lo_c !mul_hi_c (len mul_table) and
2445       div_str = sprintf "div_lo = %d, div_hi = %d (div_hash = %d)\n" !div_lo_c !div_hi_c (len div_table) and
2446       add_str = sprintf "add_lo = %d, add_hi = %d (add_hash = %d)\n" !add_lo_c !add_hi_c (len add_table) and
2447       sub_str = sprintf "sub_lo = %d, sub_hi = %d (sub_hash = %d)\n" !sub_lo_c !sub_hi_c (len sub_table) and
2448       sqrt_str = sprintf "sqrt_lo = %d, sqrt_hi = %d (sqrt_hash = %d)\n" !sqrt_lo_c !sqrt_hi_c (len sqrt_table) in
2449     print_string (cmp_str1 ^ cmp_str2 ^ cmp_str3 ^
2450                     le_str ^ max_str ^ mul_str ^ div_str ^ add_str ^ sub_str ^ sqrt_str);;
2451
2452
2453 (* lt0 *)
2454 let float_lt0 =
2455   if !Arith_options.float_cached then
2456     fun tm ->
2457       let _ = lt0_c := !lt0_c + 1 in
2458         Float_ops.float_lt0 tm
2459   else
2460     Float_ops.float_lt0;;
2461
2462 (* gt0 *)
2463 let float_gt0 =
2464   if !Arith_options.float_cached then
2465     fun tm ->
2466       let _ = gt0_c := !gt0_c + 1 in
2467         Float_ops.float_gt0 tm
2468   else
2469     Float_ops.float_gt0;;
2470
2471 (* lt *)
2472 let float_lt =
2473   if !Arith_options.float_cached then
2474     fun tm1 tm2 ->
2475       let _ = lt_c := !lt_c + 1 in
2476         Float_ops.float_lt tm1 tm2
2477   else
2478     Float_ops.float_lt;;
2479
2480 (* le0 *)
2481 let float_le0 =
2482   if !Arith_options.float_cached then
2483     fun tm ->
2484       let _ = le0_c := !le0_c + 1 in
2485         Float_ops.float_le0 tm
2486   else
2487     Float_ops.float_le0;;
2488
2489 (* ge0 *)
2490 let float_ge0 =
2491   if !Arith_options.float_cached then
2492     fun tm ->
2493       let _ = ge0_c := !ge0_c + 1 in
2494         Float_ops.float_ge0 tm
2495   else
2496     Float_ops.float_ge0;;
2497
2498 (* min *)
2499 let float_min =
2500   if !Arith_options.float_cached then
2501     fun tm1 tm2 ->
2502       let _ = min_c := !min_c + 1 in
2503         Float_ops.float_min tm1 tm2
2504   else
2505     Float_ops.float_min;;
2506
2507 (* min_max *)
2508 let float_min_max =
2509   if !Arith_options.float_cached then
2510     fun tm1 tm2 ->
2511       let _ = min_max_c := !min_max_c + 1 in
2512         Float_ops.float_min_max tm1 tm2
2513   else
2514     Float_ops.float_min_max;;
2515
2516
2517 (***************)
2518 let float_hash tm =
2519   let s, n_tm, e_tm = dest_float tm in
2520     s ^ (Arith_cache.num_tm_hash n_tm) ^ "e" ^ (Arith_cache.num_tm_hash e_tm);;
2521
2522 let float_op_hash pp tm1 tm2 =
2523   string_of_int pp ^ float_hash tm1 ^ "x" ^ float_hash tm2;;
2524
2525 let float_op_hash1 pp tm =
2526         string_of_int pp ^ float_hash tm;;
2527
2528
2529 (* le *)
2530 let float_le =
2531   if !Arith_options.float_cached then
2532     fun tm1 tm2 ->
2533       let _ = le_c := !le_c + 1 in
2534       let hash = float_op_hash 0 tm1 tm2 in
2535         try
2536           Hashtbl.find le_table hash
2537         with Not_found ->
2538           let result = Float_ops.float_le tm1 tm2 in
2539           let _ = my_add le_table hash result in
2540             result
2541   else
2542     Float_ops.float_le;;
2543
2544 (* max *)
2545 let float_max =
2546   if !Arith_options.float_cached then
2547     fun tm1 tm2 ->
2548       let _ = max_c := !max_c + 1 in
2549       let hash = float_op_hash 0 tm1 tm2 in
2550         try
2551           Hashtbl.find max_table hash
2552         with Not_found ->
2553           let result = Float_ops.float_max tm1 tm2 in
2554           let _ = my_add max_table hash result in
2555             result
2556   else
2557     Float_ops.float_max;;
2558
2559 (* mul_eq *)
2560 let float_mul_eq = Float_ops.float_mul_eq;;
2561
2562 (* mul_lo *)
2563 let float_mul_lo =
2564   if !Arith_options.float_cached then
2565     fun pp tm1 tm2 ->
2566       let _ = mul_lo_c := !mul_lo_c + 1 in
2567       let hash = "lo" ^ float_op_hash pp tm1 tm2 in
2568         try
2569           Hashtbl.find mul_table hash
2570         with Not_found ->
2571           let result = Float_ops.float_mul_lo pp tm1 tm2 in
2572           let _ = my_add mul_table hash result in
2573             result
2574   else
2575     Float_ops.float_mul_lo;;
2576
2577 (* mul_hi *)
2578 let float_mul_hi =
2579   if !Arith_options.float_cached then
2580     fun pp tm1 tm2 ->
2581       let _ = mul_hi_c := !mul_hi_c + 1 in
2582       let hash = "hi" ^ float_op_hash pp tm1 tm2 in
2583         try
2584           Hashtbl.find mul_table hash
2585         with Not_found ->
2586           let result = Float_ops.float_mul_hi pp tm1 tm2 in
2587           let _ = my_add mul_table hash result in
2588             result
2589   else
2590     Float_ops.float_mul_hi;;
2591
2592 (* div_lo *)
2593 let float_div_lo =
2594   if !Arith_options.float_cached then
2595     fun pp tm1 tm2 ->
2596       let _ = div_lo_c := !div_lo_c + 1 in
2597       let hash = "lo" ^ float_op_hash pp tm1 tm2 in
2598         try
2599           Hashtbl.find div_table hash
2600         with Not_found ->
2601           let result = Float_ops.float_div_lo pp tm1 tm2 in
2602           let _ = my_add div_table hash result in
2603             result
2604   else
2605     Float_ops.float_div_lo;;
2606
2607 (* div_hi *)
2608 let float_div_hi =
2609   if !Arith_options.float_cached then
2610     fun pp tm1 tm2 ->
2611       let _ = div_hi_c := !div_hi_c + 1 in
2612       let hash = "hi" ^ float_op_hash pp tm1 tm2 in
2613         try
2614           Hashtbl.find div_table hash
2615         with Not_found ->
2616           let result = Float_ops.float_div_hi pp tm1 tm2 in
2617           let _ = my_add div_table hash result in
2618             result
2619   else
2620     Float_ops.float_div_hi;;
2621
2622 (* add_lo *)
2623 let float_add_lo =
2624   if !Arith_options.float_cached then
2625     fun pp tm1 tm2 ->
2626       let _ = add_lo_c := !add_lo_c + 1 in
2627       let hash = "lo" ^ float_op_hash pp tm1 tm2 in
2628         try
2629           Hashtbl.find add_table hash
2630         with Not_found ->
2631           let result = Float_ops.float_add_lo pp tm1 tm2 in
2632           let _ = my_add add_table hash result in
2633             result
2634   else
2635     Float_ops.float_add_lo;;
2636
2637 (* add_hi *)
2638 let float_add_hi =
2639   if !Arith_options.float_cached then
2640     fun pp tm1 tm2 ->
2641       let _ = add_hi_c := !add_hi_c + 1 in
2642       let hash = "hi" ^ float_op_hash pp tm1 tm2 in
2643         try
2644           Hashtbl.find add_table hash
2645         with Not_found ->
2646           let result = Float_ops.float_add_hi pp tm1 tm2 in
2647           let _ = my_add add_table hash result in
2648             result
2649   else
2650     Float_ops.float_add_hi;;
2651
2652 (* sub_lo *)
2653 let float_sub_lo =
2654   if !Arith_options.float_cached then
2655     fun pp tm1 tm2 ->
2656       let _ = sub_lo_c := !sub_lo_c + 1 in
2657       let hash = "lo" ^ float_op_hash pp tm1 tm2 in
2658         try
2659           Hashtbl.find sub_table hash
2660         with Not_found ->
2661           let result = Float_ops.float_sub_lo pp tm1 tm2 in
2662           let _ = my_add sub_table hash result in
2663             result
2664   else
2665     Float_ops.float_sub_lo;;
2666
2667 (* sub_hi *)
2668 let float_sub_hi =
2669   if !Arith_options.float_cached then
2670     fun pp tm1 tm2 ->
2671       let _ = sub_hi_c := !sub_hi_c + 1 in
2672       let hash = "hi" ^ float_op_hash pp tm1 tm2 in
2673         try
2674           Hashtbl.find sub_table hash
2675         with Not_found ->
2676           let result = Float_ops.float_sub_hi pp tm1 tm2 in
2677           let _ = my_add sub_table hash result in
2678             result
2679   else
2680     Float_ops.float_sub_hi;;
2681
2682 (* sqrt_lo *)
2683 let float_sqrt_lo =
2684   if !Arith_options.float_cached then
2685     fun pp tm ->
2686       let _ = sqrt_lo_c := !sqrt_lo_c + 1 in
2687       let hash = "lo" ^ float_op_hash1 pp tm in
2688         try
2689           Hashtbl.find sqrt_table hash
2690         with Not_found ->
2691           let result = Float_ops.float_sqrt_lo pp tm in
2692           let _ = my_add sqrt_table hash result in
2693             result
2694   else
2695     Float_ops.float_sqrt_lo;;
2696
2697 (* sqrt_hi *)
2698 let float_sqrt_hi =
2699   if !Arith_options.float_cached then
2700     fun pp tm ->
2701       let _ = sqrt_hi_c := !sqrt_hi_c + 1 in
2702       let hash = "hi" ^ float_op_hash1 pp tm in
2703         try
2704           Hashtbl.find sqrt_table hash
2705         with Not_found ->
2706           let result = Float_ops.float_sqrt_hi pp tm in
2707           let _ = my_add sqrt_table hash result in
2708             result
2709   else
2710     Float_ops.float_sqrt_hi;;
2711
2712
2713 (******************************************)
2714 (* float intervals *)
2715
2716 let FLOAT_OF_NUM' = (SPEC_ALL o REWRITE_RULE[min_exp_def]) FLOAT_OF_NUM;;
2717
2718 let FLOAT_INTERVAL_OF_NUM = (NUMERALS_TO_NUM o REWRITE_RULE[min_exp_def] o prove)(`interval_arith (&n) (float_num F n min_exp, float_num F n min_exp)`,
2719    REWRITE_TAC[FLOAT_OF_NUM; CONST_INTERVAL]);;
2720
2721 let FLOAT_F_bound' = (UNDISCH_ALL o SPEC_ALL) FLOAT_F_bound;;
2722
2723 let FLOAT_T_bound' = (UNDISCH_ALL o SPEC_ALL) FLOAT_T_bound;;
2724
2725 (* interval_arith x (float_num s1 n1 e1, float_num s2 n2 e2) -> x, float_num s1 n1 e1, float_num s2 n2 e2 *)
2726 let dest_float_interval tm =
2727   let ltm, rtm = dest_comb tm in
2728   let f1, f2 = dest_pair rtm in
2729     rand ltm, f1, f2;;
2730
2731 let mk_float_interval_small_num n =
2732   let n_tm0 = mk_small_numeral n in
2733   let n_th = NUMERAL_TO_NUM_CONV n_tm0 in
2734   let n_tm = rand(rand(concl n_th)) in
2735   let n_th1 = TRANS n_th (INST[n_tm, n_var_num] NUM_REMOVE) in
2736   let th1 = AP_TERM amp_op_real n_th1 in
2737   let int_th = INST[n_tm, n_var_num] FLOAT_INTERVAL_OF_NUM in
2738   let rtm = rand(concl int_th) in
2739     EQ_MP (SYM (AP_THM (AP_TERM interval_const th1) rtm)) int_th;;
2740
2741 let mk_float_interval_num n =
2742   let n_tm0 = mk_numeral n in
2743   let n_th = NUMERAL_TO_NUM_CONV n_tm0 in
2744   let n_tm = rand(rand(concl n_th)) in
2745   let n_th1 = TRANS n_th (INST[n_tm, n_var_num] NUM_REMOVE) in
2746   let th1 = AP_TERM amp_op_real n_th1 in
2747   let int_th = INST[n_tm, n_var_num] FLOAT_INTERVAL_OF_NUM in
2748   let rtm = rand(concl int_th) in
2749     EQ_MP (SYM (AP_THM (AP_TERM interval_const th1) rtm)) int_th;;
2750
2751
2752 (* Returns the lower bound for the given float *)
2753 let float_lo p tm =
2754   let s, n_tm, e_tm = dest_float tm in
2755     if s = "F" then
2756       let num_exp_tm = mk_num_exp n_tm e_tm in
2757       let th0 = num_exp_lo p num_exp_tm in
2758       let ltm, e1_tm = dest_comb(lhand(concl th0)) in
2759       let n1_tm = rand ltm in
2760       let th1 = INST[n1_tm, n1_var_num; e1_tm, e1_var_num; n_tm, n2_var_num; e_tm, e2_var_num] FLOAT_F_bound' in
2761         MY_PROVE_HYP th0 th1
2762     else
2763       let num_exp_tm = mk_num_exp n_tm e_tm in
2764       let th0 = num_exp_hi p num_exp_tm in
2765       let ltm, e1_tm = dest_comb(rand(concl th0)) in
2766       let n1_tm = rand ltm in
2767       let th1 = INST[n_tm, n1_var_num; e_tm, e1_var_num; n1_tm, n2_var_num; e1_tm, e2_var_num] FLOAT_T_bound' in
2768         MY_PROVE_HYP th0 th1;;
2769
2770 (* Returns the upper bound for the given float *)
2771 let float_hi p tm =
2772   let s, n_tm, e_tm = dest_float tm in
2773     if s = "F" then
2774       let num_exp_tm = mk_num_exp n_tm e_tm in
2775       let th0 = num_exp_hi p num_exp_tm in
2776       let ltm, e2_tm = dest_comb(rand(concl th0)) in
2777       let n2_tm = rand ltm in
2778       let th1 = INST[n_tm, n1_var_num; e_tm, e1_var_num; n2_tm, n2_var_num; e2_tm, e2_var_num] FLOAT_F_bound' in
2779         MY_PROVE_HYP th0 th1
2780     else
2781       let num_exp_tm = mk_num_exp n_tm e_tm in
2782       let th0 = num_exp_lo p num_exp_tm in
2783       let ltm, e1_tm = dest_comb(lhand(concl th0)) in
2784       let n1_tm = rand ltm in
2785       let th1 = INST[n1_tm, n1_var_num; e1_tm, e1_var_num; n_tm, n2_var_num; e_tm, e2_var_num] FLOAT_T_bound' in
2786         MY_PROVE_HYP th0 th1;;
2787
2788
2789 (* Approximates the given interval with p-digits floating point numbers *)
2790 let float_interval_round p th =
2791   let x_tm, f1, f2 = dest_float_interval (concl th) in
2792   let lo_th = float_lo p f1 in
2793   let hi_th = float_hi p f2 in
2794   let lo_tm = lhand(concl lo_th) in
2795   let hi_tm = rand(concl hi_th) in
2796   let th0 = INST[x_tm, x_var_real; f1, lo_var_real; f2, hi_var_real; lo_tm, a_var_real; hi_tm, b_var_real] APPROX_INTERVAL' in
2797     MY_PROVE_HYP lo_th (MY_PROVE_HYP hi_th (MY_PROVE_HYP th th0));;
2798
2799
2800 (****************************************)
2801 (* float_interval_lt *)
2802
2803 let FLOAT_INTERVAL_LT = prove(`interval_arith x (lo1, hi1) /\ interval_arith y (lo2, hi2) /\ hi1 < lo2
2804                                 ==> x < y`,
2805                               REWRITE_TAC[interval_arith] THEN
2806                                 REAL_ARITH_TAC);;
2807
2808
2809 (****************************************)
2810 (* float_interval_neg *)
2811
2812 let FLOAT_INTERVAL_NEG = prove(`!s1 s2. interval_arith x (float_num s1 n1 e1, float_num s2 n2 e2)
2813                                  ==> interval_arith (--x) (float_num (~s2) n2 e2, float_num (~s1) n1 e1)`,
2814    REPEAT GEN_TAC THEN
2815      DISCH_THEN (fun th -> MP_TAC (MATCH_MP INTERVAL_NEG th)) THEN
2816      SIMP_TAC[FLOAT_NEG]);;
2817
2818 let FLOAT_INTERVAL_NEG_FF = (UNDISCH_ALL o REWRITE_RULE[] o SPECL[`F`; `F`]) FLOAT_INTERVAL_NEG;;
2819 let FLOAT_INTERVAL_NEG_FT = (UNDISCH_ALL o REWRITE_RULE[] o SPECL[`F`; `T`]) FLOAT_INTERVAL_NEG;;
2820 let FLOAT_INTERVAL_NEG_TF = (UNDISCH_ALL o REWRITE_RULE[] o SPECL[`T`; `F`]) FLOAT_INTERVAL_NEG;;
2821 let FLOAT_INTERVAL_NEG_TT = (UNDISCH_ALL o REWRITE_RULE[] o SPECL[`T`; `T`]) FLOAT_INTERVAL_NEG;;
2822
2823 (* |- interval x (float s1 n1 e1, float s2 n2 e2) ->
2824    |- interval (--x) (float ~s2 n2 e2, float ~s1 n1 e1 *)
2825 let float_interval_neg th =
2826   let x_tm, f1, f2 = dest_float_interval (concl th) in
2827   let s1, n1_tm, e1_tm = dest_float f1 in
2828   let s2, n2_tm, e2_tm = dest_float f2 in
2829   let inst = INST[x_tm, x_var_real; n1_tm, n1_var_num; e1_tm, e1_var_num;
2830                   n2_tm, n2_var_num; e2_tm, e2_var_num] in
2831   let th0 =
2832     if s1 = "F" then
2833       if s2 = "F" then
2834         inst FLOAT_INTERVAL_NEG_FF
2835       else
2836         inst FLOAT_INTERVAL_NEG_FT
2837     else
2838       if s2 = "F" then
2839         inst FLOAT_INTERVAL_NEG_TF
2840       else
2841         inst FLOAT_INTERVAL_NEG_TT in
2842     MY_PROVE_HYP th th0;;
2843
2844
2845 (***********************************************)
2846 (* float_interval_mul *)
2847
2848 let f1_1_var = `f1_1:real` and
2849     f1_2_var = `f1_2:real` and
2850     f2_1_var = `f2_1:real` and
2851     f2_2_var = `f2_2:real`;;
2852
2853
2854 let FLOAT_INTERVAL_FT_IMP_0 = prove(`interval_arith x (float_num F n1 e1, float_num T n2 e2) ==> x = &0`,
2855                                     REWRITE_TAC[interval_arith] THEN STRIP_TAC THEN
2856                                       REWRITE_TAC[GSYM REAL_LE_ANTISYM] THEN
2857                                       CONJ_TAC THEN MATCH_MP_TAC REAL_LE_TRANS THENL [
2858                                         EXISTS_TAC `float_num T n2 e2` THEN ASM_REWRITE_TAC[FLOAT_T_NEG];
2859                                         EXISTS_TAC `float_num F n1 e1` THEN ASM_REWRITE_TAC[FLOAT_F_POS]
2860                                       ]);;
2861
2862 (* FT_xx *)
2863 let FLOAT_INTERVAL_MUL_FT_xx = (UNDISCH_ALL o NUMERALS_TO_NUM o REWRITE_RULE[GSYM IMP_IMP; min_exp_def] o prove)(
2864   `interval_arith x (float_num F n1 e1, float_num T n2 e2)
2865   ==> interval_arith (x * y) (float_num F 0 min_exp, float_num F 0 min_exp)`,
2866   STRIP_TAC THEN
2867     FIRST_X_ASSUM (fun th -> REWRITE_TAC[MATCH_MP FLOAT_INTERVAL_FT_IMP_0 th]) THEN
2868     REWRITE_TAC[REAL_MUL_LZERO; interval_arith] THEN
2869     MP_TAC (GEN_ALL (SPECL [`s:bool`; `0`] FLOAT_EQ_0)) THEN SIMP_TAC[REAL_LE_REFL]);;
2870
2871
2872 (* xx_FT *)
2873 let FLOAT_INTERVAL_MUL_xx_FT = (UNDISCH_ALL o NUMERALS_TO_NUM o REWRITE_RULE[GSYM IMP_IMP; min_exp_def] o prove)(
2874   `interval_arith y (float_num F m1 r1, float_num T m2 r2)
2875     ==> interval_arith (x * y) (float_num F 0 min_exp, float_num F 0 min_exp)`,
2876   STRIP_TAC THEN
2877     FIRST_X_ASSUM (fun th -> REWRITE_TAC[MATCH_MP FLOAT_INTERVAL_FT_IMP_0 th]) THEN
2878     REWRITE_TAC[REAL_MUL_RZERO; interval_arith] THEN
2879     MP_TAC (GEN_ALL (SPECL [`s:bool`; `0`] FLOAT_EQ_0)) THEN SIMP_TAC[REAL_LE_REFL]);;
2880
2881
2882
2883 (* FF_FF *)
2884 let FLOAT_INTERVAL_MUL_FF_FF = (UNDISCH_ALL o REWRITE_RULE[GSYM IMP_IMP] o prove)(
2885   `interval_arith x (float_num F n1 e1, float_num F n2 e2) /\
2886     interval_arith y (float_num F m1 r1, float_num F m2 r2) /\
2887     f1 <= float_num F n1 e1 * float_num F m1 r1 /\
2888     float_num F n2 e2 * float_num F m2 r2 <= f2
2889     ==> interval_arith (x * y) (f1, f2)`,
2890    MAP_EVERY ABBREV_TAC [`a = float_num F n1 e1`; `b = float_num F n2 e2`; `c = float_num F m1 r1`; `d = float_num F m2 r2`] THEN
2891      SUBGOAL_THEN `&0 <= a /\ &0 <= b /\ &0 <= c /\ &0 <= d` MP_TAC THENL
2892      [
2893        MAP_EVERY EXPAND_TAC ["a"; "b"; "c"; "d"] THEN
2894          REWRITE_TAC[FLOAT_F_POS];
2895        ALL_TAC
2896      ] THEN
2897      REPEAT (POP_ASSUM (fun th -> ALL_TAC)) THEN
2898      REWRITE_TAC[interval_arith] THEN
2899      REPEAT STRIP_TAC THENL
2900      [
2901        MATCH_MP_TAC REAL_LE_TRANS THEN
2902          EXISTS_TAC `a * c:real` THEN
2903          ASM_REWRITE_TAC[] THEN
2904          MATCH_MP_TAC REAL_LE_MUL2 THEN
2905          ASM_REWRITE_TAC[];
2906        MATCH_MP_TAC REAL_LE_TRANS THEN
2907          EXISTS_TAC `b * d:real` THEN
2908          ASM_REWRITE_TAC[] THEN
2909          MATCH_MP_TAC REAL_LE_MUL2 THEN
2910          ASM_REWRITE_TAC[] THEN
2911          CONJ_TAC THEN MATCH_MP_TAC REAL_LE_TRANS THENL
2912          [
2913            EXISTS_TAC `a:real` THEN ASM_REWRITE_TAC[];
2914            EXISTS_TAC `c:real` THEN ASM_REWRITE_TAC[]
2915          ]
2916      ]);;
2917
2918 (* TT_TT *)
2919 let FLOAT_INTERVAL_MUL_TT_TT = (UNDISCH_ALL o REWRITE_RULE[GSYM IMP_IMP] o prove)(
2920   `interval_arith x (float_num T n1 e1, float_num T n2 e2) /\
2921     interval_arith y (float_num T m1 r1, float_num T m2 r2) /\
2922     f1 <= float_num T n2 e2 * float_num T m2 r2 /\
2923     float_num T n1 e1 * float_num T m1 r1 <= f2
2924     ==> interval_arith (x * y) (f1, f2)`,
2925   REWRITE_TAC[FLOAT_NEG_T] THEN
2926     REWRITE_TAC[interval_arith] THEN
2927     MAP_EVERY ABBREV_TAC [`a = float_num F n1 e1`; `b = float_num F n2 e2`; `c = float_num F m1 r1`; `d = float_num F m2 r2`] THEN
2928     SUBGOAL_THEN `&0 <= a /\ &0 <= b /\ &0 <= c /\ &0 <= d` MP_TAC THENL
2929     [
2930       MAP_EVERY EXPAND_TAC ["a"; "b"; "c"; "d"] THEN
2931         REWRITE_TAC[FLOAT_F_POS];
2932       ALL_TAC
2933     ] THEN
2934     REWRITE_TAC[REAL_NEG_MUL2] THEN
2935     REPEAT STRIP_TAC THENL
2936     [
2937       MATCH_MP_TAC REAL_LE_TRANS THEN
2938         EXISTS_TAC `b * d:real` THEN
2939         ASM_REWRITE_TAC[] THEN
2940         ONCE_REWRITE_TAC[REAL_ARITH `a <= x * y <=> a <= --x * --y`] THEN
2941         MATCH_MP_TAC REAL_LE_MUL2 THEN
2942         ONCE_REWRITE_TAC[REAL_ARITH `b <= --x <=> x <= --b`] THEN
2943         ASM_REWRITE_TAC[];
2944       MATCH_MP_TAC REAL_LE_TRANS THEN
2945         EXISTS_TAC `a * c:real` THEN
2946         ASM_REWRITE_TAC[] THEN
2947         ONCE_REWRITE_TAC[REAL_ARITH `x * y <= a <=> --x * --y <= a`] THEN
2948         MATCH_MP_TAC REAL_LE_MUL2 THEN
2949         ONCE_REWRITE_TAC[REAL_ARITH `--x <= c <=> --c <= x`] THEN
2950         ASM_REWRITE_TAC[] THEN
2951         ASSUME_TAC (REAL_ARITH `!b x. &0 <= b /\ x <= --b ==> &0 <= --x`) THEN
2952         CONJ_TAC THEN FIRST_X_ASSUM MATCH_MP_TAC THENL
2953         [
2954           EXISTS_TAC `b:real` THEN ASM_REWRITE_TAC[];
2955           EXISTS_TAC `d:real` THEN ASM_REWRITE_TAC[]
2956         ]
2957     ]);;
2958
2959 (* FF_TT *)
2960 let FLOAT_INTERVAL_MUL_FF_TT = (UNDISCH_ALL o REWRITE_RULE[GSYM IMP_IMP] o prove)(
2961   `interval_arith x (float_num F n1 e1, float_num F n2 e2) /\
2962     interval_arith y (float_num T m1 r1, float_num T m2 r2) /\
2963     f1 <= float_num F n2 e2 * float_num T m1 r1 /\
2964     float_num F n1 e1 * float_num T m2 r2 <= f2
2965     ==> interval_arith (x * y) (f1, f2)`,
2966   REWRITE_TAC[FLOAT_NEG_T] THEN
2967     REWRITE_TAC[interval_arith] THEN
2968     MAP_EVERY ABBREV_TAC [`a = float_num F n1 e1`; `b = float_num F n2 e2`; `c = float_num F m1 r1`; `d = float_num F m2 r2`] THEN
2969     SUBGOAL_THEN `&0 <= a /\ &0 <= b /\ &0 <= c /\ &0 <= d` MP_TAC THENL
2970     [
2971       MAP_EVERY EXPAND_TAC ["a"; "b"; "c"; "d"] THEN
2972         REWRITE_TAC[FLOAT_F_POS];
2973       ALL_TAC
2974     ] THEN
2975
2976     REPEAT STRIP_TAC THENL
2977     [
2978       MATCH_MP_TAC REAL_LE_TRANS THEN
2979         EXISTS_TAC `b * --c` THEN
2980         ASM_REWRITE_TAC[] THEN
2981         ONCE_REWRITE_TAC[REAL_ARITH `b * --c <= x * y <=> x * --y <= b * c`] THEN
2982         MATCH_MP_TAC REAL_LE_MUL2 THEN
2983         ONCE_REWRITE_TAC[REAL_ARITH `--y <= c <=> --c <= y`] THEN
2984         ASM_REWRITE_TAC[] THEN
2985         CONJ_TAC THEN MATCH_MP_TAC REAL_LE_TRANS THENL
2986         [
2987           EXISTS_TAC `a:real` THEN ASM_REWRITE_TAC[];
2988           EXISTS_TAC `d:real` THEN
2989             ONCE_REWRITE_TAC[REAL_ARITH `d <= --y <=> y <= --d`] THEN
2990             ASM_REWRITE_TAC[]
2991         ];
2992
2993       MATCH_MP_TAC REAL_LE_TRANS THEN
2994         EXISTS_TAC `a * --d` THEN
2995         ASM_REWRITE_TAC[] THEN
2996         ONCE_REWRITE_TAC[REAL_ARITH `x * y <= a * --d <=> a * d <= x * --y`] THEN
2997         MATCH_MP_TAC REAL_LE_MUL2 THEN
2998         ONCE_REWRITE_TAC[REAL_ARITH `d <= --y <=> y <= --d`] THEN
2999         ASM_REWRITE_TAC[]
3000     ]);;
3001
3002 (* TT_FF *)
3003 let FLOAT_INTERVAL_MUL_TT_FF = (UNDISCH_ALL o REWRITE_RULE[GSYM IMP_IMP] o prove)(
3004   `interval_arith x (float_num T n1 e1, float_num T n2 e2) /\
3005     interval_arith y (float_num F m1 r1, float_num F m2 r2) /\
3006     f1 <= float_num T n1 e1 * float_num F m2 r2 /\
3007     float_num T n2 e2 * float_num F m1 r1 <= f2
3008     ==> interval_arith (x * y) (f1, f2)`,
3009    STRIP_TAC THEN
3010      MP_TAC ((GEN_ALL o DISCH_ALL) FLOAT_INTERVAL_MUL_FF_TT) THEN
3011      DISCH_THEN (MP_TAC o SPECL[`n1:num`; `e1:num`; `n2:num`; `e2:num`; `m1:num`; `r1:num`; `m2:num`; `r2:num`]) THEN
3012      DISCH_THEN (MP_TAC o SPECL[`y:real`; `x:real`; `f1:real`; `f2:real`]) THEN
3013      REPEAT (POP_ASSUM MP_TAC) THEN
3014      SIMP_TAC[REAL_MUL_AC]);;
3015
3016 (* TF_FF *)
3017 let FLOAT_INTERVAL_MUL_TF_FF = (UNDISCH_ALL o REWRITE_RULE[GSYM IMP_IMP] o prove)(
3018   `interval_arith x (float_num T n1 e1, float_num F n2 e2) /\
3019     interval_arith y (float_num F m1 r1, float_num F m2 r2) /\
3020     f1 <= float_num T n1 e1 * float_num F m2 r2 /\
3021     float_num F n2 e2 * float_num F m2 r2 <= f2
3022     ==> interval_arith (x * y) (f1, f2)`,
3023   REWRITE_TAC[FLOAT_NEG_T] THEN
3024     REWRITE_TAC[interval_arith] THEN
3025     MAP_EVERY ABBREV_TAC [`a = float_num F n1 e1`; `b = float_num F n2 e2`; `c = float_num F m1 r1`; `d = float_num F m2 r2`] THEN
3026     STRIP_TAC THEN
3027     SUBGOAL_THEN `&0 <= a /\ &0 <= b /\ &0 <= c /\ &0 <= d` ASSUME_TAC THENL
3028     [
3029       MAP_EVERY EXPAND_TAC ["a"; "b"; "c"; "d"] THEN
3030         REWRITE_TAC[FLOAT_F_POS];
3031       ALL_TAC
3032     ] THEN
3033
3034     SUBGOAL_THEN `&0 <= y` ASSUME_TAC THENL
3035     [
3036       MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `c:real` THEN ASM_REWRITE_TAC[];
3037       ALL_TAC
3038     ] THEN
3039
3040     CONJ_TAC THENL
3041     [
3042       DISJ_CASES_TAC (REAL_ARITH `&0 <= x \/ &0 <= --x`) THENL
3043         [
3044           MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `&0` THEN
3045             CONJ_TAC THENL
3046             [
3047               MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `--a * d` THEN
3048                 ASM_REWRITE_TAC[REAL_ARITH `--a * d <= &0 <=> &0 <= a * d`] THEN
3049                 MATCH_MP_TAC REAL_LE_MUL THEN ASM_REWRITE_TAC[];
3050               MATCH_MP_TAC REAL_LE_MUL THEN ASM_REWRITE_TAC[]
3051             ];
3052           ALL_TAC
3053         ] THEN
3054
3055         MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `--a * d` THEN ASM_REWRITE_TAC[] THEN
3056         ONCE_REWRITE_TAC[REAL_ARITH `--a * d <= x * y <=> --x * y <= a * d`] THEN
3057         MATCH_MP_TAC REAL_LE_MUL2 THEN
3058         ONCE_REWRITE_TAC[REAL_ARITH `--x <= a <=> --a <= x`] THEN ASM_REWRITE_TAC[];
3059
3060       DISJ_CASES_TAC (REAL_ARITH `&0 <= --x \/ &0 <= x`) THENL
3061         [
3062           MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `&0` THEN
3063             CONJ_TAC THENL
3064             [
3065               REWRITE_TAC[REAL_ARITH `x * y <= &0 <=> &0 <= --x * y`] THEN
3066                 MATCH_MP_TAC REAL_LE_MUL THEN ASM_REWRITE_TAC[];
3067               MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `b * d` THEN ASM_REWRITE_TAC[] THEN
3068                 MATCH_MP_TAC REAL_LE_MUL THEN ASM_REWRITE_TAC[]
3069             ];
3070           ALL_TAC
3071         ] THEN
3072
3073         MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `b * d` THEN ASM_REWRITE_TAC[] THEN
3074         MATCH_MP_TAC REAL_LE_MUL2 THEN ASM_REWRITE_TAC[]
3075     ]);;
3076
3077 (* TF_TT *)
3078 let FLOAT_INTERVAL_MUL_TF_TT = (UNDISCH_ALL o REWRITE_RULE[GSYM IMP_IMP] o prove)(
3079   `interval_arith x (float_num T n1 e1, float_num F n2 e2) /\
3080     interval_arith y (float_num T m1 r1, float_num T m2 r2) /\
3081     f1 <= float_num F n2 e2 * float_num T m1 r1 /\
3082     float_num T n1 e1 * float_num T m1 r1 <= f2
3083     ==> interval_arith (x * y) (f1, f2)`,
3084   REWRITE_TAC[FLOAT_NEG_T] THEN
3085     REWRITE_TAC[interval_arith] THEN
3086     MAP_EVERY ABBREV_TAC [`a = float_num F n1 e1`; `b = float_num F n2 e2`; `c = float_num F m1 r1`; `d = float_num F m2 r2`] THEN
3087     STRIP_TAC THEN
3088     SUBGOAL_THEN `&0 <= a /\ &0 <= b /\ &0 <= c /\ &0 <= d` ASSUME_TAC THENL
3089     [
3090       MAP_EVERY EXPAND_TAC ["a"; "b"; "c"; "d"] THEN
3091         REWRITE_TAC[FLOAT_F_POS];
3092       ALL_TAC
3093     ] THEN
3094
3095     SUBGOAL_THEN `&0 <= --y` ASSUME_TAC THENL
3096     [
3097       MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `d:real` THEN
3098         ONCE_REWRITE_TAC[REAL_ARITH `d <= --y <=> y <= --d`] THEN ASM_REWRITE_TAC[];
3099       ALL_TAC
3100     ] THEN
3101
3102     CONJ_TAC THENL
3103     [
3104       DISJ_CASES_TAC (REAL_ARITH `&0 <= --x \/ &0 <= x`) THENL
3105         [
3106           MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `&0` THEN
3107             CONJ_TAC THENL
3108             [
3109               MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `b * --c` THEN
3110                 ASM_REWRITE_TAC[REAL_ARITH `b * --c <= &0 <=> &0 <= b * c`] THEN
3111                 MATCH_MP_TAC REAL_LE_MUL THEN ASM_REWRITE_TAC[];
3112               ONCE_REWRITE_TAC[REAL_ARITH `x * y = --x * --y`] THEN
3113               MATCH_MP_TAC REAL_LE_MUL THEN ASM_REWRITE_TAC[]
3114             ];
3115           ALL_TAC
3116         ] THEN
3117
3118         MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `b * --c` THEN ASM_REWRITE_TAC[] THEN
3119         ONCE_REWRITE_TAC[REAL_ARITH `b * --c <= x * y <=> x * --y <= b * c`] THEN
3120         MATCH_MP_TAC REAL_LE_MUL2 THEN ONCE_REWRITE_TAC[REAL_ARITH `--y <= c <=> --c <= y`] THEN
3121         ASM_REWRITE_TAC[];
3122
3123       DISJ_CASES_TAC (REAL_ARITH `&0 <= x \/ &0 <= --x`) THENL
3124         [
3125           MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `&0` THEN
3126             CONJ_TAC THENL
3127             [
3128               REWRITE_TAC[REAL_ARITH `x * y <= &0 <=> &0 <= x * --y`] THEN
3129                 MATCH_MP_TAC REAL_LE_MUL THEN ASM_REWRITE_TAC[];
3130               MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `--a * --c` THEN ASM_REWRITE_TAC[REAL_NEG_MUL2] THEN
3131                 MATCH_MP_TAC REAL_LE_MUL THEN ASM_REWRITE_TAC[]
3132             ];
3133           ALL_TAC
3134         ] THEN
3135
3136         MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `--a * --c` THEN ASM_REWRITE_TAC[REAL_NEG_MUL2] THEN
3137         ONCE_REWRITE_TAC[REAL_ARITH `x * y <= a * c <=> --x * --y <= a * c`] THEN
3138         MATCH_MP_TAC REAL_LE_MUL2 THEN ONCE_REWRITE_TAC[REAL_ARITH `--x <= a <=> --a <= x`] THEN
3139         ASM_REWRITE_TAC[]
3140     ]);;
3141
3142 (* FF_TF *)
3143 let FLOAT_INTERVAL_MUL_FF_TF = (UNDISCH_ALL o REWRITE_RULE[GSYM IMP_IMP] o prove)(
3144   `interval_arith x (float_num F n1 e1, float_num F n2 e2) /\
3145     interval_arith y (float_num T m1 r1, float_num F m2 r2) /\
3146     f1 <= float_num F n2 e2 * float_num T m1 r1 /\
3147     float_num F n2 e2 * float_num F m2 r2 <= f2
3148     ==> interval_arith (x * y) (f1, f2)`,
3149   STRIP_TAC THEN
3150     MP_TAC ((SPECL [`n1:num`; `e1:num`; `n2:num`; `e2:num`; `m1:num`; `r1:num`; `m2:num`; `r2:num`; `y:real`; `x:real`; `f1:real`; `f2:real`] o GEN_ALL o DISCH_ALL) FLOAT_INTERVAL_MUL_TF_FF) THEN
3151     ASM_REWRITE_TAC[REAL_MUL_SYM] THEN DISCH_THEN MATCH_MP_TAC THEN
3152     POP_ASSUM MP_TAC THEN REAL_ARITH_TAC);;
3153
3154 (* TT_TF *)
3155 let FLOAT_INTERVAL_MUL_TT_TF = (UNDISCH_ALL o REWRITE_RULE[GSYM IMP_IMP] o prove)(
3156   `interval_arith x (float_num T n1 e1, float_num T n2 e2) /\
3157     interval_arith y (float_num T m1 r1, float_num F m2 r2) /\
3158     f1 <= float_num T n1 e1 * float_num F m2 r2 /\
3159     float_num T n1 e1 * float_num T m1 r1 <= f2
3160     ==> interval_arith (x * y) (f1, f2)`,
3161   STRIP_TAC THEN
3162     MP_TAC ((SPECL [`n1:num`; `e1:num`; `n2:num`; `e2:num`; `m1:num`; `r1:num`; `m2:num`; `r2:num`; `y:real`; `x:real`; `f1:real`; `f2:real`] o GEN_ALL o DISCH_ALL) FLOAT_INTERVAL_MUL_TF_TT) THEN
3163     ASM_REWRITE_TAC[REAL_MUL_SYM] THEN POP_ASSUM MP_TAC THEN POP_ASSUM MP_TAC THEN
3164     SIMP_TAC[REAL_MUL_SYM]);;
3165
3166 (* TF_TF *)
3167 let FLOAT_INTERVAL_MUL_TF_TF = (UNDISCH_ALL o REWRITE_RULE[GSYM IMP_IMP] o prove)(
3168   `interval_arith x (float_num T n1 e1, float_num F n2 e2) /\
3169     interval_arith y (float_num T m1 r1, float_num F m2 r2) /\
3170     f1_1 <= float_num T n1 e1 * float_num F m2 r2 /\
3171     f1_2 <= float_num F n2 e2 * float_num T m1 r1 /\
3172     min f1_1 f1_2 = f1 /\
3173     float_num T n1 e1 * float_num T m1 r1 <= f2_1 /\
3174     float_num F n2 e2 * float_num F m2 r2 <= f2_2 /\
3175     max f2_1 f2_2 = f2
3176     ==> interval_arith (x * y) (f1, f2)`,
3177   REWRITE_TAC[EQ_SYM_EQ; FLOAT_NEG_T] THEN
3178     REWRITE_TAC[interval_arith] THEN
3179     MAP_EVERY ABBREV_TAC [`a = float_num F n1 e1`; `b = float_num F n2 e2`; `c = float_num F m1 r1`; `d = float_num F m2 r2`] THEN
3180     STRIP_TAC THEN
3181     SUBGOAL_THEN `&0 <= a /\ &0 <= b /\ &0 <= c /\ &0 <= d` ASSUME_TAC THENL
3182     [
3183       MAP_EVERY EXPAND_TAC ["a"; "b"; "c"; "d"] THEN
3184         REWRITE_TAC[FLOAT_F_POS];
3185       ALL_TAC
3186     ] THEN
3187
3188     DISJ_CASES_TAC (REAL_ARITH `&0 <= x \/ &0 <= --x`) THENL
3189     [
3190       DISJ_CASES_TAC (REAL_ARITH `&0 <= y \/ &0 <= --y`) THENL
3191         [
3192           CONJ_TAC THENL
3193             [
3194               MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `&0` THEN
3195                 CONJ_TAC THENL
3196                 [
3197                   MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `f1_1:real` THEN
3198                     ASM_REWRITE_TAC[REAL_MIN_MIN] THEN
3199                     MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `--a * d` THEN
3200                     ASM_REWRITE_TAC[REAL_ARITH `--a * d <= &0 <=> &0 <= a * d`] THEN
3201                     MATCH_MP_TAC REAL_LE_MUL THEN ASM_REWRITE_TAC[];
3202                   ALL_TAC
3203                 ] THEN
3204                 MATCH_MP_TAC REAL_LE_MUL THEN ASM_REWRITE_TAC[];
3205               ALL_TAC
3206             ] THEN
3207
3208             MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `b * d` THEN
3209             CONJ_TAC THENL
3210             [
3211               MATCH_MP_TAC REAL_LE_MUL2 THEN ASM_REWRITE_TAC[];
3212               ALL_TAC
3213             ] THEN
3214             MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `f2_2:real` THEN
3215             ASM_REWRITE_TAC[REAL_MAX_MAX];
3216           ALL_TAC
3217         ] THEN
3218
3219         CONJ_TAC THENL
3220         [
3221           MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `f1_2:real` THEN
3222             ASM_REWRITE_TAC[REAL_MIN_MIN] THEN
3223             MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `b * --c` THEN
3224             ONCE_REWRITE_TAC[REAL_ARITH `b * --c <= x * y <=> x * --y <= b * c`] THEN ASM_REWRITE_TAC[] THEN
3225             MATCH_MP_TAC REAL_LE_MUL2 THEN ONCE_REWRITE_TAC[REAL_ARITH `--y <= c <=> --c <= y`] THEN
3226             ASM_REWRITE_TAC[];
3227           ALL_TAC
3228         ] THEN
3229         MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `&0` THEN
3230         CONJ_TAC THENL
3231         [
3232           REWRITE_TAC[REAL_ARITH `x * y <= &0 <=> &0 <= x * --y`] THEN
3233             MATCH_MP_TAC REAL_LE_MUL THEN ASM_REWRITE_TAC[];
3234           ALL_TAC
3235         ] THEN
3236         MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `f2_2:real` THEN
3237         ASM_REWRITE_TAC[REAL_MAX_MAX] THEN
3238         MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `b * d` THEN ASM_REWRITE_TAC[] THEN
3239         MATCH_MP_TAC REAL_LE_MUL THEN ASM_REWRITE_TAC[];
3240       ALL_TAC
3241     ] THEN
3242
3243     DISJ_CASES_TAC (REAL_ARITH `&0 <= y \/ &0 <= --y`) THENL
3244     [
3245       CONJ_TAC THENL
3246         [
3247           MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `f1_1:real` THEN
3248             ASM_REWRITE_TAC[REAL_MIN_MIN] THEN
3249             MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `--a * d` THEN
3250             ASM_REWRITE_TAC[] THEN ONCE_REWRITE_TAC[REAL_ARITH `--a * d <= x * y <=> --x * y <= a * d`] THEN
3251             MATCH_MP_TAC REAL_LE_MUL2 THEN ONCE_REWRITE_TAC[REAL_ARITH `--x <= a <=> --a <= x`] THEN
3252             ASM_REWRITE_TAC[];
3253           ALL_TAC
3254         ] THEN
3255         MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `&0` THEN
3256         CONJ_TAC THENL
3257         [
3258           REWRITE_TAC[REAL_ARITH `x * y <= &0 <=> &0 <= --x * y`] THEN
3259             MATCH_MP_TAC REAL_LE_MUL THEN ASM_REWRITE_TAC[];
3260           ALL_TAC
3261         ] THEN
3262         MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `f2_2:real` THEN
3263         ASM_REWRITE_TAC[REAL_MAX_MAX] THEN
3264         MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `b * d` THEN ASM_REWRITE_TAC[] THEN
3265         MATCH_MP_TAC REAL_LE_MUL THEN ASM_REWRITE_TAC[];
3266       ALL_TAC
3267     ] THEN
3268
3269     CONJ_TAC THENL
3270     [
3271       MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `&0` THEN
3272         CONJ_TAC THENL
3273         [
3274           MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `f1_1:real` THEN
3275             ASM_REWRITE_TAC[REAL_MIN_MIN] THEN
3276             MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `--a * d` THEN
3277             ASM_REWRITE_TAC[REAL_ARITH `--a * d <= &0 <=> &0 <= a * d`] THEN
3278             MATCH_MP_TAC REAL_LE_MUL THEN ASM_REWRITE_TAC[];
3279           ALL_TAC
3280         ] THEN
3281         ONCE_REWRITE_TAC[GSYM REAL_NEG_MUL2] THEN
3282         MATCH_MP_TAC REAL_LE_MUL THEN ASM_REWRITE_TAC[];
3283       ALL_TAC
3284     ] THEN
3285
3286     MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `--a * --c` THEN
3287     CONJ_TAC THENL
3288     [
3289       ONCE_REWRITE_TAC[GSYM REAL_NEG_MUL2] THEN REWRITE_TAC[REAL_NEG_NEG] THEN
3290         MATCH_MP_TAC REAL_LE_MUL2 THEN ONCE_REWRITE_TAC[REAL_ARITH `--x <= a <=> --a <= x`] THEN
3291         ASM_REWRITE_TAC[];
3292       ALL_TAC
3293     ] THEN
3294     MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `f2_1:real` THEN
3295     ASM_REWRITE_TAC[REAL_MAX_MAX]);;
3296
3297
3298 (****************************)
3299 let float_interval_mul =
3300   let mul_ft_xx th1 x y n1 e1 n2 e2 =
3301     let th0 = INST[x, x_var_real; y, y_var_real;
3302                    n1, n1_var_num; e1, e1_var_num;
3303                    n2, n2_var_num; e2, e2_var_num] FLOAT_INTERVAL_MUL_FT_xx in
3304       MY_PROVE_HYP th1 th0 in
3305   let mul_xx_ft th2 x y m1 r1 m2 r2 =
3306     let th0 = INST[x, x_var_real; y, y_var_real;
3307                    m1, m1_var_num; r1, r1_var_num;
3308                    m2, m2_var_num; r2, r2_var_num] FLOAT_INTERVAL_MUL_xx_FT in
3309       MY_PROVE_HYP th2 th0 in
3310
3311     fun pp th1 th2 ->
3312       let x, l_lo, l_hi = dest_float_interval (concl th1) and
3313           y, r_lo, r_hi = dest_float_interval (concl th2) in
3314       let s1, n1, e1 = dest_float l_lo and
3315           s2, n2, e2 = dest_float l_hi and
3316           s3, m1, r1 = dest_float r_lo and
3317           s4, m2, r2 = dest_float r_hi in
3318
3319         (* Special case 1 *)
3320         if s1 <> s2 && s1 = "F" then
3321           mul_ft_xx th1 x y n1 e1 n2 e2
3322         else if s3 <> s4 && s3 = "F" then
3323           mul_xx_ft th2 x y m1 r1 m2 r2
3324         else
3325           (* Special case 2 *)
3326           if s1 <> s2 && s3 <> s4 then
3327             let lo1, lo2 = float_mul_lo pp l_lo r_hi, float_mul_lo pp l_hi r_lo and
3328                 hi1, hi2 = float_mul_hi pp l_lo r_lo, float_mul_hi pp l_hi r_hi in
3329             let f1_1 = (lhand o concl) lo1 and
3330                 f1_2 = (lhand o concl) lo2 and
3331                 f2_1 = (rand o concl) hi1 and
3332                 f2_2 = (rand o concl) hi2 in
3333             let min_th = float_min f1_1 f1_2 and
3334                 max_th = float_max f2_1 f2_2 in
3335             let f1_tm = (rand o concl) min_th and
3336                 f2_tm = (rand o concl) max_th in
3337             let th0 = INST[x, x_var_real; n1, n1_var_num; e1, e1_var_num;
3338                            y, y_var_real; n2, n2_var_num; e2, e2_var_num;
3339                            m1, m1_var_num; r1, r1_var_num;
3340                            m2, m2_var_num; r2, r2_var_num;
3341                            f1_tm, f1_var_real; f2_tm, f2_var_real;
3342                            f1_1, f1_1_var; f1_2, f1_2_var;
3343                            f2_1, f2_1_var; f2_2, f2_2_var] FLOAT_INTERVAL_MUL_TF_TF in
3344               (MY_PROVE_HYP min_th o MY_PROVE_HYP max_th o MY_PROVE_HYP lo1 o MY_PROVE_HYP lo2 o
3345                  MY_PROVE_HYP hi1 o MY_PROVE_HYP hi2 o MY_PROVE_HYP th1 o MY_PROVE_HYP th2) th0
3346           else
3347             let lo_th, hi_th, th0 =
3348               if s1 <> s2 then
3349                 if s3 = "F" then
3350                   float_mul_lo pp l_lo r_hi, float_mul_hi pp l_hi r_hi, FLOAT_INTERVAL_MUL_TF_FF
3351                 else
3352                   float_mul_lo pp l_hi r_lo, float_mul_hi pp l_lo r_lo, FLOAT_INTERVAL_MUL_TF_TT
3353               else
3354                 if s3 <> s4 then
3355                   if s1 = "F" then
3356                     float_mul_lo pp l_hi r_lo, float_mul_hi pp l_hi r_hi, FLOAT_INTERVAL_MUL_FF_TF
3357                   else
3358                     float_mul_lo pp l_lo r_hi, float_mul_hi pp l_lo r_lo, FLOAT_INTERVAL_MUL_TT_TF
3359                 else
3360                   if s1 = "F" then
3361                     if s3 = "F" then
3362                       float_mul_lo pp l_lo r_lo, float_mul_hi pp l_hi r_hi, FLOAT_INTERVAL_MUL_FF_FF
3363                     else
3364                       float_mul_lo pp l_hi r_lo, float_mul_hi pp l_lo r_hi, FLOAT_INTERVAL_MUL_FF_TT
3365                   else
3366                     if s3 = "F" then
3367                       float_mul_lo pp l_lo r_hi, float_mul_hi pp l_hi r_lo, FLOAT_INTERVAL_MUL_TT_FF
3368                     else
3369                       float_mul_lo pp l_hi r_hi, float_mul_hi pp l_lo r_lo, FLOAT_INTERVAL_MUL_TT_TT in
3370
3371             let f1_tm = lhand(concl lo_th) and
3372                 f2_tm = rand(concl hi_th) in
3373
3374             let th = INST[x, x_var_real; n1, n1_var_num; e1, e1_var_num;
3375                           y, y_var_real; n2, n2_var_num; e2, e2_var_num;
3376                           m1, m1_var_num; r1, r1_var_num;
3377                           m2, m2_var_num; r2, r2_var_num;
3378                           f1_tm, f1_var_real; f2_tm, f2_var_real] th0 in
3379               MY_PROVE_HYP lo_th (MY_PROVE_HYP hi_th (MY_PROVE_HYP th1 (MY_PROVE_HYP th2 th)));;
3380
3381
3382 (*************************************)
3383 (* float_interval_div *)
3384
3385 (* FT_xx *)
3386 let FLOAT_INTERVAL_DIV_FT_xx = prove(
3387   `interval_arith x (float_num F n1 e1, float_num T n2 e2)
3388   ==> interval_arith (x / y) (float_num F 0 min_exp, float_num F 0 min_exp)`,
3389   REWRITE_TAC[real_div] THEN DISCH_THEN (MP_TAC o MATCH_MP FLOAT_INTERVAL_FT_IMP_0) THEN
3390     SIMP_TAC[REAL_MUL_LZERO; interval_arith] THEN
3391     MP_TAC (GEN_ALL (SPECL [`s:bool`; `0`] FLOAT_EQ_0)) THEN SIMP_TAC[REAL_LE_REFL]);;
3392
3393 (* FF_FF *)
3394 let FLOAT_INTERVAL_DIV_FF_FF = prove(
3395   `~(m1 = 0) /\
3396     interval_arith x (float_num F n1 e1, float_num F n2 e2) /\
3397     interval_arith y (float_num F m1 r1, float_num F m2 r2) /\
3398     f1 <= float_num F n1 e1 / float_num F m2 r2 /\
3399     float_num F n2 e2 / float_num F m1 r1 <= f2
3400     ==> interval_arith (x / y) (f1, f2)`,
3401    MAP_EVERY ABBREV_TAC [`a = float_num F n1 e1`; `b = float_num F n2 e2`; `c = float_num F m1 r1`; `d = float_num F m2 r2`] THEN
3402      REWRITE_TAC[real_div] THEN
3403      STRIP_TAC THEN
3404      SUBGOAL_THEN `&0 <= a /\ &0 <= b /\ &0 <= c /\ &0 <= d` MP_TAC THENL
3405      [
3406        MAP_EVERY EXPAND_TAC ["a"; "b"; "c"; "d"] THEN
3407          REWRITE_TAC[FLOAT_F_POS];
3408        ALL_TAC
3409      ] THEN
3410      SUBGOAL_THEN `~(c = &0)` ASSUME_TAC THENL
3411      [
3412        EXPAND_TAC "c" THEN ASM_REWRITE_TAC[FLOAT_EQ_0];
3413        ALL_TAC
3414      ] THEN
3415      STRIP_TAC THEN
3416      SUBGOAL_THEN `~(d = &0)` MP_TAC THENL
3417      [
3418        MATCH_MP_TAC (REAL_ARITH `~(c = &0) /\ &0 <= c /\ c <= d ==> ~(d = &0)`) THEN
3419          ASM_REWRITE_TAC[] THEN
3420          UNDISCH_TAC `interval_arith y (c,d)` THEN
3421          REWRITE_TAC[interval_arith] THEN
3422          REAL_ARITH_TAC;
3423        ALL_TAC
3424      ] THEN
3425
3426      REPLICATE_TAC 10 (POP_ASSUM MP_TAC) THEN
3427      REPEAT (POP_ASSUM (fun th -> ALL_TAC)) THEN
3428      REWRITE_TAC[interval_arith] THEN
3429      REPEAT STRIP_TAC THENL
3430      [
3431        MATCH_MP_TAC REAL_LE_TRANS THEN
3432          EXISTS_TAC `a * inv d` THEN
3433          ASM_REWRITE_TAC[] THEN
3434          MATCH_MP_TAC REAL_LE_MUL2 THEN
3435          ASM_REWRITE_TAC[REAL_LE_INV_EQ] THEN
3436          MATCH_MP_TAC REAL_LE_INV2 THEN
3437          ASM_REWRITE_TAC[] THEN
3438          MATCH_MP_TAC REAL_LTE_TRANS THEN
3439          EXISTS_TAC `c:real` THEN
3440          ASM_REWRITE_TAC[REAL_ARITH `&0 < c <=> ~(c = &0) /\ &0 <= c`];
3441
3442        MATCH_MP_TAC REAL_LE_TRANS THEN
3443          EXISTS_TAC `b * inv c` THEN
3444          ASM_REWRITE_TAC[] THEN
3445          MATCH_MP_TAC REAL_LE_MUL2 THEN
3446          ASM_REWRITE_TAC[REAL_LE_INV_EQ] THEN
3447          REPEAT CONJ_TAC THENL
3448          [
3449            MATCH_MP_TAC REAL_LE_TRANS THEN
3450              EXISTS_TAC `a:real` THEN ASM_REWRITE_TAC[];
3451            MATCH_MP_TAC REAL_LE_TRANS THEN
3452              EXISTS_TAC `c:real` THEN ASM_REWRITE_TAC[];
3453            ALL_TAC
3454          ] THEN
3455          MATCH_MP_TAC REAL_LE_INV2 THEN
3456          ASM_REWRITE_TAC[REAL_ARITH `&0 < c <=> ~(c = &0) /\ &0 <= c`]
3457      ]);;
3458
3459 (* TT_TT *)
3460 let FLOAT_INTERVAL_DIV_TT_TT = prove(
3461   `~(m2 = 0) /\
3462     interval_arith x (float_num T n1 e1, float_num T n2 e2) /\
3463     interval_arith y (float_num T m1 r1, float_num T m2 r2) /\
3464     f1 <= float_num T n2 e2 / float_num T m1 r1 /\
3465     float_num T n1 e1 / float_num T m2 r2 <= f2
3466     ==> interval_arith (x / y) (f1,f2)`,
3467    REWRITE_TAC[FLOAT_NEG_T] THEN
3468      REWRITE_TAC[real_div; REAL_INV_NEG; REAL_NEG_MUL2] THEN
3469      GEN_REWRITE_TAC (LAND_CONV o DEPTH_CONV)[interval_arith] THEN
3470      REWRITE_TAC[REAL_ARITH `--a <= x /\ x <= --b <=> b <= --x /\ --x <= a`] THEN
3471      GEN_REWRITE_TAC (LAND_CONV o DEPTH_CONV)[GSYM interval_arith] THEN
3472      STRIP_TAC THEN
3473      MP_TAC ((SPECL[`n2:num`; `e2:num`; `m1:num`; `r1:num`; `n1:num`; `e1:num`; `m2:num`; `r2:num`; `--x`; `--y`] o GEN_ALL) FLOAT_INTERVAL_DIV_FF_FF) THEN
3474      REWRITE_TAC[real_div; REAL_INV_NEG; REAL_NEG_MUL2] THEN
3475      DISCH_THEN MATCH_MP_TAC THEN
3476      ASM_REWRITE_TAC[]);;
3477
3478 (* FF_TT *)
3479 let FLOAT_INTERVAL_DIV_FF_TT = prove(
3480   `~(m2 = 0) /\
3481     interval_arith x (float_num F n1 e1, float_num F n2 e2) /\
3482     interval_arith y (float_num T m1 r1, float_num T m2 r2) /\
3483     f1 <= float_num F n2 e2 / float_num T m2 r2 /\
3484     float_num F n1 e1 / float_num T m1 r1 <= f2
3485     ==> interval_arith (x / y) (f1,f2)`,
3486    REWRITE_TAC[FLOAT_NEG_T] THEN
3487      REWRITE_TAC[real_div; REAL_INV_NEG] THEN
3488      REWRITE_TAC[REAL_ARITH `a * --b <= c <=> --c <= a * b`] THEN
3489      REWRITE_TAC[REAL_ARITH `c <= a * --b <=> a * b <= --c`] THEN
3490      GEN_REWRITE_TAC (LAND_CONV o DEPTH_CONV)[interval_arith] THEN
3491      REWRITE_TAC[REAL_ARITH `--a <= x /\ x <= --b <=> b <= --x /\ --x <= a`] THEN
3492      GEN_REWRITE_TAC (LAND_CONV o DEPTH_CONV)[GSYM interval_arith] THEN
3493      STRIP_TAC THEN
3494      MP_TAC ((SPECL[`n1:num`; `e1:num`; `m1:num`; `r1:num`; `n2:num`; `e2:num`; `m2:num`; `r2:num`; `x:real`; `--y`; `--f2`; `--f1`] o GEN_ALL) FLOAT_INTERVAL_DIV_FF_FF) THEN
3495      ANTS_TAC THENL
3496      [
3497        ASM_REWRITE_TAC[real_div];
3498        ALL_TAC
3499      ] THEN
3500      REWRITE_TAC[real_div; REAL_INV_NEG; interval_arith] THEN
3501      REAL_ARITH_TAC);;
3502
3503 (* TT_FF *)
3504 let FLOAT_INTERVAL_DIV_TT_FF = prove(
3505   `~(m1 = 0) /\
3506     interval_arith x (float_num T n1 e1, float_num T n2 e2) /\
3507     interval_arith y (float_num F m1 r1, float_num F m2 r2) /\
3508     f1 <= float_num T n1 e1 / float_num F m1 r1 /\
3509     float_num T n2 e2 / float_num F m2 r2 <= f2
3510     ==> interval_arith (x / y) (f1,f2)`,
3511    REWRITE_TAC[FLOAT_NEG_T] THEN
3512      REWRITE_TAC[real_div; REAL_INV_NEG] THEN
3513      REWRITE_TAC[REAL_ARITH `--a * b <= c <=> --c <= a * b`] THEN
3514      REWRITE_TAC[REAL_ARITH `c <= --a * b <=> a * b <= --c`] THEN
3515      GEN_REWRITE_TAC (LAND_CONV o DEPTH_CONV)[interval_arith] THEN
3516      REWRITE_TAC[REAL_ARITH `--a <= x /\ x <= --b <=> b <= --x /\ --x <= a`] THEN
3517      GEN_REWRITE_TAC (LAND_CONV o DEPTH_CONV)[GSYM interval_arith] THEN
3518      STRIP_TAC THEN
3519      MP_TAC ((SPECL[`n2:num`; `e2:num`; `m2:num`; `r2:num`; `n1:num`; `e1:num`; `m1:num`; `r1:num`; `--x:real`; `y:real`; `--f2`; `--f1`] o GEN_ALL) FLOAT_INTERVAL_DIV_FF_FF) THEN
3520      ANTS_TAC THENL
3521      [
3522        ASM_REWRITE_TAC[real_div];
3523        ALL_TAC
3524      ] THEN
3525      REWRITE_TAC[real_div; REAL_INV_NEG; interval_arith] THEN
3526      REAL_ARITH_TAC);;
3527
3528 let FLOAT_0 = prove(`!s e. float_num s 0 e = &0`, REWRITE_TAC[FLOAT_EQ_0]);;
3529
3530 (* TF_FF *)
3531 let FLOAT_INTERVAL_DIV_TF_FF = prove(
3532   `~(m1 = 0) /\
3533     interval_arith x (float_num T n1 e1, float_num F n2 e2) /\
3534     interval_arith y (float_num F m1 r1, float_num F m2 r2) /\
3535     f1 <= float_num T n1 e1 / float_num F m1 r1 /\
3536     float_num F n2 e2 / float_num F m1 r1 <= f2
3537     ==> interval_arith (x / y) (f1,f2)`,
3538    REWRITE_TAC[FLOAT_NEG_T] THEN STRIP_TAC THEN
3539      MAP_EVERY ABBREV_TAC [`a = float_num F n1 e1`; `b = float_num F n2 e2`; `c = float_num F m1 r1`; `d = float_num F m2 r2`] THEN
3540      SUBGOAL_THEN `&0 <= a /\ &0 <= b /\ &0 <= c /\ &0 <= d` ASSUME_TAC THENL
3541      [
3542        MAP_EVERY EXPAND_TAC ["a"; "b"; "c"; "d"] THEN REWRITE_TAC[FLOAT_F_POS];
3543        ALL_TAC
3544      ] THEN
3545
3546      DISJ_CASES_TAC (REAL_ARITH `&0 <= x \/ x <= &0`) THENL
3547      [
3548        SUBGOAL_THEN `interval_arith x (float_num F 0 0, b)` ASSUME_TAC THENL
3549          [
3550            UNDISCH_TAC `interval_arith x (--a, b)` THEN POP_ASSUM MP_TAC THEN
3551              REWRITE_TAC[interval_arith; FLOAT_0] THEN REAL_ARITH_TAC;
3552            ALL_TAC
3553          ] THEN
3554
3555          MP_TAC ((SPEC_ALL o SPECL [`0`; `0`] o GEN_ALL) FLOAT_INTERVAL_DIV_FF_FF) THEN
3556          ASM_REWRITE_TAC[] THEN DISCH_THEN MATCH_MP_TAC THEN
3557          MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `&0` THEN
3558          REWRITE_TAC[real_div] THEN CONJ_TAC THENL
3559          [
3560            MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `--a / c` THEN
3561              ASM_REWRITE_TAC[real_div; ARITH_RULE `--a * b <= &0 <=> &0 <= a * b`] THEN
3562              MATCH_MP_TAC REAL_LE_MUL THEN ASM_REWRITE_TAC[REAL_LE_INV_EQ];
3563            ALL_TAC
3564          ] THEN
3565
3566          MATCH_MP_TAC REAL_LE_MUL THEN ASM_REWRITE_TAC[REAL_LE_INV_EQ; FLOAT_F_POS];
3567        ALL_TAC
3568      ] THEN
3569
3570      SUBGOAL_THEN `interval_arith x (--a, float_num T 0 0)` ASSUME_TAC THENL
3571      [
3572        UNDISCH_TAC `interval_arith x (--a, b)` THEN POP_ASSUM MP_TAC THEN
3573          REWRITE_TAC[interval_arith; FLOAT_0] THEN REAL_ARITH_TAC;
3574        ALL_TAC
3575      ] THEN
3576
3577      MP_TAC ((INST[`0`, `n2:num`; `0`, `e2:num`]) FLOAT_INTERVAL_DIV_TT_FF) THEN
3578      ASM_REWRITE_TAC[FLOAT_NEG_T] THEN DISCH_THEN MATCH_MP_TAC THEN
3579      ASM_REWRITE_TAC[GSYM FLOAT_NEG_T] THEN MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `&0` THEN
3580      REWRITE_TAC[real_div] THEN CONJ_TAC THENL
3581      [
3582        REWRITE_TAC[FLOAT_0; REAL_MUL_LZERO; REAL_LE_REFL];
3583        ALL_TAC
3584      ] THEN
3585
3586      MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `b / c` THEN ASM_REWRITE_TAC[] THEN
3587      MATCH_MP_TAC REAL_LE_DIV THEN ASM_REWRITE_TAC[]);;
3588
3589 (* TF_TT *)
3590 let FLOAT_INTERVAL_DIV_TF_TT = prove(
3591   `~(m2 = 0) /\
3592     interval_arith x (float_num T n1 e1, float_num F n2 e2) /\
3593     interval_arith y (float_num T m1 r1, float_num T m2 r2) /\
3594     f1 <= float_num F n2 e2 / float_num T m2 r2 /\
3595     float_num T n1 e1 / float_num T m2 r2 <= f2
3596     ==> interval_arith (x / y) (f1,f2)`,
3597    REWRITE_TAC[FLOAT_NEG_T] THEN STRIP_TAC THEN
3598      MAP_EVERY ABBREV_TAC [`a = float_num F n1 e1`; `b = float_num F n2 e2`; `c = float_num F m1 r1`; `d = float_num F m2 r2`] THEN
3599      SUBGOAL_THEN `&0 <= a /\ &0 <= b /\ &0 <= c /\ &0 <= d` ASSUME_TAC THENL
3600      [
3601        MAP_EVERY EXPAND_TAC ["a"; "b"; "c"; "d"] THEN REWRITE_TAC[FLOAT_F_POS];
3602        ALL_TAC
3603      ] THEN
3604
3605      DISJ_CASES_TAC (REAL_ARITH `x <= &0 \/ &0 <= x`) THENL
3606      [
3607        SUBGOAL_THEN `interval_arith x (--a, float_num T 0 0)` ASSUME_TAC THENL
3608          [
3609            UNDISCH_TAC `interval_arith x (--a, b)` THEN POP_ASSUM MP_TAC THEN
3610              REWRITE_TAC[interval_arith; FLOAT_0] THEN REAL_ARITH_TAC;
3611            ALL_TAC
3612          ] THEN
3613
3614          MP_TAC ((SPEC_ALL o SPECL [`0`; `0`] o GEN_ALL) FLOAT_INTERVAL_DIV_TT_TT) THEN
3615          ASM_REWRITE_TAC[] THEN DISCH_THEN MATCH_MP_TAC THEN ASM_REWRITE_TAC[FLOAT_NEG_T] THEN
3616          ASM_REWRITE_TAC[GSYM FLOAT_NEG_T] THEN
3617          MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `&0` THEN
3618          REWRITE_TAC[real_div] THEN CONJ_TAC THENL
3619          [
3620            MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `b / --d` THEN
3621              ASM_REWRITE_TAC[real_div; REAL_INV_NEG; REAL_ARITH `b * --d <= &0 <=> &0 <= b * d`] THEN
3622              MATCH_MP_TAC REAL_LE_MUL THEN ASM_REWRITE_TAC[REAL_LE_INV_EQ];
3623            ALL_TAC
3624          ] THEN
3625
3626          REWRITE_TAC[FLOAT_0; REAL_MUL_LZERO; REAL_LE_REFL];
3627        ALL_TAC
3628      ] THEN
3629
3630      SUBGOAL_THEN `interval_arith x (float_num F 0 0, b)` ASSUME_TAC THENL
3631      [
3632        UNDISCH_TAC `interval_arith x (--a, b)` THEN POP_ASSUM MP_TAC THEN
3633          REWRITE_TAC[interval_arith; FLOAT_0] THEN REAL_ARITH_TAC;
3634        ALL_TAC
3635      ] THEN
3636
3637      MP_TAC ((INST[`0`, `n1:num`; `0`, `e1:num`]) FLOAT_INTERVAL_DIV_FF_TT) THEN
3638      ASM_REWRITE_TAC[FLOAT_NEG_T] THEN DISCH_THEN MATCH_MP_TAC THEN
3639      MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `&0` THEN
3640      REWRITE_TAC[real_div] THEN CONJ_TAC THENL
3641      [
3642        REWRITE_TAC[FLOAT_0; REAL_MUL_LZERO; REAL_LE_REFL];
3643        ALL_TAC
3644      ] THEN
3645
3646      MATCH_MP_TAC REAL_LE_TRANS THEN EXISTS_TAC `--a / --d` THEN ASM_REWRITE_TAC[] THEN
3647      REWRITE_TAC[real_div; REAL_INV_NEG; REAL_NEG_MUL2] THEN MATCH_MP_TAC REAL_LE_MUL THEN
3648      ASM_REWRITE_TAC[REAL_LE_INV_EQ]);;
3649
3650 let transform = UNDISCH_ALL o
3651   PURE_REWRITE_RULE[TAUT `~P <=> (P <=> F)`] o
3652   NUMERALS_TO_NUM o
3653   REWRITE_RULE[GSYM IMP_IMP; min_exp_def];;
3654
3655 let FLOAT_INTERVAL_DIV_FT_xx' = transform FLOAT_INTERVAL_DIV_FT_xx and
3656     FLOAT_INTERVAL_DIV_FF_FF' = transform FLOAT_INTERVAL_DIV_FF_FF and
3657     FLOAT_INTERVAL_DIV_TT_TT' = transform FLOAT_INTERVAL_DIV_TT_TT and
3658     FLOAT_INTERVAL_DIV_FF_TT' = transform FLOAT_INTERVAL_DIV_FF_TT and
3659     FLOAT_INTERVAL_DIV_TT_FF' = transform FLOAT_INTERVAL_DIV_TT_FF and
3660     FLOAT_INTERVAL_DIV_TF_FF' = transform FLOAT_INTERVAL_DIV_TF_FF and
3661     FLOAT_INTERVAL_DIV_TF_TT' = transform FLOAT_INTERVAL_DIV_TF_TT;;
3662
3663
3664 let float_interval_div pp th1 th2 =
3665   let x, l_lo, l_hi = dest_float_interval (concl th1) and
3666       y, r_lo, r_hi = dest_float_interval (concl th2) in
3667   let s1, n1, e1 = dest_float l_lo and
3668       s2, n2, e2 = dest_float l_hi and
3669       s3, m1, r1 = dest_float r_lo and
3670       s4, m2, r2 = dest_float r_hi in
3671
3672     if s1 <> s2 && s1 = "F" then
3673       let th0 = INST[x, x_var_real; y, y_var_real;
3674                      n1, n1_var_num; e1, e1_var_num;
3675                      n2, n2_var_num; e2, e2_var_num] FLOAT_INTERVAL_DIV_FT_xx' in
3676         MY_PROVE_HYP th1 th0
3677     else if s3 <> s4 then
3678       failwith "float_interval_div: division by an interval containing 0"
3679     else
3680       let lo_th, hi_th, th0, zero_th =
3681         if s1 = s2 then
3682           if s1 = "F" then
3683             if s3 = "F" then
3684               float_div_lo pp l_lo r_hi, float_div_hi pp l_hi r_lo,
3685         FLOAT_INTERVAL_DIV_FF_FF', raw_eq0_hash_conv m1
3686             else
3687               float_div_lo pp l_hi r_hi, float_div_hi pp l_lo r_lo,
3688         FLOAT_INTERVAL_DIV_FF_TT', raw_eq0_hash_conv m2
3689           else
3690             if s3 = "F" then
3691               float_div_lo pp l_lo r_lo, float_div_hi pp l_hi r_hi,
3692         FLOAT_INTERVAL_DIV_TT_FF', raw_eq0_hash_conv m1
3693             else
3694               float_div_lo pp l_hi r_lo, float_div_hi pp l_lo r_hi,
3695         FLOAT_INTERVAL_DIV_TT_TT', raw_eq0_hash_conv m2
3696         else
3697           if s3 = "F" then
3698             float_div_lo pp l_lo r_lo, float_div_hi pp l_hi r_lo,
3699         FLOAT_INTERVAL_DIV_TF_FF', raw_eq0_hash_conv m1
3700           else
3701             float_div_lo pp l_hi r_hi, float_div_hi pp l_lo r_hi,
3702         FLOAT_INTERVAL_DIV_TF_TT', raw_eq0_hash_conv m2 in
3703
3704       let f1_tm = lhand(concl lo_th) and
3705           f2_tm = rand(concl hi_th) in
3706
3707       let th = INST[x, x_var_real; n1, n1_var_num; e1, e1_var_num;
3708                     y, y_var_real; n2, n2_var_num; e2, e2_var_num;
3709                     m1, m1_var_num; r1, r1_var_num;
3710                     m2, m2_var_num; r2, r2_var_num;
3711                     f1_tm, f1_var_real; f2_tm, f2_var_real] th0 in
3712         (MY_PROVE_HYP lo_th o MY_PROVE_HYP hi_th o
3713            MY_PROVE_HYP th1 o MY_PROVE_HYP th2 o MY_PROVE_HYP zero_th) th;;
3714
3715
3716 (*****************************************)
3717 (* float_interval_add, float_interval_sub *)
3718
3719 let n1_var_real = `n1:real` and
3720     n2_var_real = `n2:real` and
3721     m1_var_real = `m1:real` and
3722     m2_var_real = `m2:real` and
3723     n_var_real = `n:real` and
3724     m_var_real = `m:real`;;
3725
3726 let INTERVAL_ADD = (UNDISCH_ALL o REWRITE_RULE[GSYM IMP_IMP] o prove)(
3727   `interval_arith x (n1, m1) /\
3728     interval_arith y (n2, m2) /\
3729     n <= n1 + n2 /\ m1 + m2 <= m
3730     ==> interval_arith (x + y) (n, m)`,
3731    REWRITE_TAC[interval_arith] THEN REAL_ARITH_TAC);;
3732
3733 let INTERVAL_SUB = (UNDISCH_ALL o REWRITE_RULE[GSYM IMP_IMP] o prove)(
3734   `interval_arith x (n1, m1) /\
3735     interval_arith y (n2, m2) /\
3736     n <= n1 - m2 /\ m1 - n2 <= m
3737     ==> interval_arith (x - y) (n, m)`,
3738   REWRITE_TAC[interval_arith] THEN REAL_ARITH_TAC);;
3739
3740
3741  let float_interval_add pp th1 th2 =
3742   let x, n1, m1 = dest_float_interval (concl th1) in
3743   let y, n2, m2 = dest_float_interval (concl th2) in
3744   let lo_th = float_add_lo pp n1 n2 in
3745   let hi_th = float_add_hi pp m1 m2 in
3746   let n_tm = lhand (concl lo_th) in
3747   let m_tm = rand (concl hi_th) in
3748   let th0 = INST[x, x_var_real; n1, n1_var_real; m1, m1_var_real;
3749                  y, y_var_real; n2, n2_var_real; m2, m2_var_real;
3750                  n_tm, n_var_real; m_tm, m_var_real] INTERVAL_ADD in
3751     MY_PROVE_HYP lo_th (MY_PROVE_HYP hi_th (MY_PROVE_HYP th2 (MY_PROVE_HYP th1 th0)));;
3752
3753 let float_interval_sub pp th1 th2 =
3754   let x, n1, m1 = dest_float_interval (concl th1) in
3755   let y, n2, m2 = dest_float_interval (concl th2) in
3756   let lo_th = float_sub_lo pp n1 m2 in
3757   let hi_th = float_sub_hi pp m1 n2 in
3758   let n_tm = lhand(concl lo_th) in
3759   let m_tm = rand(concl hi_th) in
3760   let th0 = INST[x, x_var_real; n1, n1_var_real; m1, m1_var_real;
3761                  y, y_var_real; n2, n2_var_real; m2, m2_var_real;
3762                  n_tm, n_var_real; m_tm, m_var_real] INTERVAL_SUB in
3763     MY_PROVE_HYP lo_th (MY_PROVE_HYP hi_th (MY_PROVE_HYP th2 (MY_PROVE_HYP th1 th0)));;
3764
3765
3766 (********************************************)
3767 (* FLOAT_ABS *)
3768
3769 let s_var_bool = `s:bool`;;
3770
3771 let FLOAT_ABS = prove(`abs (float_num s n e) = float_num F n e`,
3772    BOOL_CASES_TAC `s:bool` THEN
3773      REWRITE_TAC[FLOAT_NEG_T; REAL_ABS_NEG; REAL_ABS_REFL; FLOAT_F_POS]);;
3774
3775 let float_abs tm =
3776   let ltm, rtm = dest_comb tm in
3777     if ((fst o dest_const) ltm <> "real_abs") then
3778       failwith "float_abs: no abs"
3779     else
3780       let ltm, e_tm = dest_comb rtm in
3781       let ltm, n_tm = dest_comb ltm in
3782       let s_tm = rand ltm in
3783         INST[s_tm, s_var_bool; n_tm, n_var_num; e_tm, e_var_num] FLOAT_ABS;;
3784
3785
3786 (*******************************)
3787 (* float_interval_sqrt *)
3788
3789 let FLOAT_INTERVAL_SQRT = prove(`interval_arith x (float_num F n1 e1, hi) /\
3790                                   f1 <= sqrt (float_num F n1 e1) /\ sqrt hi <= f2
3791                                   ==> interval_arith (sqrt x) (f1, f2)`,
3792    ABBREV_TAC `lo = float_num F n1 e1` THEN
3793      REWRITE_TAC[interval_arith] THEN
3794      STRIP_TAC THEN
3795      SUBGOAL_THEN `&0 <= lo /\ &0 <= hi` ASSUME_TAC THENL
3796      [
3797        EXPAND_TAC "lo" THEN
3798          REWRITE_TAC[FLOAT_F_POS] THEN
3799          MATCH_MP_TAC REAL_LE_TRANS THEN
3800          EXISTS_TAC `x:real` THEN
3801          ASM_REWRITE_TAC[] THEN
3802          MATCH_MP_TAC REAL_LE_TRANS THEN
3803          EXISTS_TAC `lo:real` THEN
3804          ASM_REWRITE_TAC[] THEN
3805          EXPAND_TAC "lo" THEN
3806          REWRITE_TAC[FLOAT_F_POS];
3807        ALL_TAC
3808      ] THEN
3809      SUBGOAL_THEN `&0 <= x` ASSUME_TAC THENL
3810      [
3811        MATCH_MP_TAC REAL_LE_TRANS THEN
3812          EXISTS_TAC `lo:real` THEN
3813          ASM_REWRITE_TAC[];
3814        ALL_TAC
3815      ] THEN
3816      CONJ_TAC THENL
3817      [
3818        MATCH_MP_TAC REAL_LE_TRANS THEN
3819          EXISTS_TAC `sqrt lo` THEN
3820          ASM_REWRITE_TAC[] THEN
3821          MATCH_MP_TAC SQRT_MONO_LE THEN
3822          ASM_REWRITE_TAC[];
3823        MATCH_MP_TAC REAL_LE_TRANS THEN
3824          EXISTS_TAC `sqrt hi` THEN
3825          ASM_REWRITE_TAC[] THEN
3826          MATCH_MP_TAC SQRT_MONO_LE THEN
3827          ASM_REWRITE_TAC[]
3828      ]);;
3829
3830
3831 let FLOAT_INTERVAL_SQRT' = (UNDISCH_ALL o REWRITE_RULE[GSYM IMP_IMP]) FLOAT_INTERVAL_SQRT;;
3832
3833 let float_interval_sqrt pp th =
3834   let x_tm, lo_tm, hi_tm = dest_float_interval (concl th) in
3835   let s1, n1_tm, e1_tm = dest_float lo_tm in
3836     if s1 <> "F" then
3837       failwith "float_interval_sqrt: negative low bound"
3838     else
3839       let lo_th = float_sqrt_lo pp lo_tm in
3840       let hi_th = float_sqrt_hi pp hi_tm in
3841       let f1_tm = lhand (concl lo_th) in
3842       let f2_tm = rand (concl hi_th) in
3843       let th0 = INST[x_tm, x_var_real; n1_tm, n1_var_num; e1_tm, e1_var_num;
3844                      hi_tm, hi_var_real; f1_tm, f1_var_real;
3845                      f2_tm, f2_var_real] FLOAT_INTERVAL_SQRT' in
3846         MY_PROVE_HYP lo_th (MY_PROVE_HYP hi_th (MY_PROVE_HYP th th0));;
3847
3848
3849 (******************************************)
3850 (* FLOAT_TO_NUM_CONV *)
3851
3852 let FLOAT_TO_NUM_CONV tm =
3853   let ltm, e_tm = dest_comb tm in
3854   let f_tm, n_tm = dest_comb ltm in
3855     if (fst o dest_const o rator) f_tm <> "float_num" then
3856       failwith "FLOAT_TO_NUM_CONV"
3857     else
3858       let n_th' = SYM (INST[n_tm, n_var_num] Arith_hash.NUM_THM) in
3859       let e_th' = SYM (INST[e_tm, n_var_num] Arith_hash.NUM_THM) in
3860       let n_th = TRANS n_th' (NUM_TO_NUMERAL_CONV (mk_comb(Arith_hash.num_const, n_tm))) in
3861       let e_th = TRANS e_th' (NUM_TO_NUMERAL_CONV (mk_comb(Arith_hash.num_const, e_tm))) in
3862       let th0 = MK_COMB (AP_TERM f_tm n_th, e_th) in
3863       let tm0 = rand(concl th0) in
3864       let th1 = REWRITE_CONV[float; num_exp; REAL_MUL_LID; GSYM REAL_OF_NUM_MUL; GSYM REAL_OF_NUM_POW; min_exp_def] tm0 in
3865       let th2 = REAL_RAT_REDUCE_CONV (rand(concl th1)) in
3866         TRANS th0 (TRANS th1 th2);;
3867
3868 end;;
3869
3870
3871 (**************************************)
3872 (* Printer for floating-point numbers *)
3873 (**************************************)
3874
3875 let print_float tm=
3876   try
3877     let s, m_tm, e_tm = Arith_float.dest_float tm in
3878     let m = Arith_hash.raw_dest_hash m_tm and
3879         e = Arith_hash.raw_dest_hash e_tm -/ Num.num_of_int Float_theory.min_exp in
3880     let s_str = if s = "T" then "-" else "" in
3881     let m_str = Num.string_of_num m in
3882     let e_str = if e = num_0 then ""
3883     else "*" ^ string_of_int Arith_hash.arith_base ^ "^" ^ Num.string_of_num e in
3884     let str = "##" ^ s_str ^ m_str ^ e_str in
3885       Format.print_string str
3886   with _ -> failwith "print_float";;
3887
3888
3889 install_user_printer ("float_num", print_float);;