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