binary_ieeeScript.sml

1(* ------------------------------------------------------------------------
2   Theory of IEEE-754 (base 2) floating-point (basic) arithmetic
3   ------------------------------------------------------------------------ *)
4Theory binary_ieee
5Ancestors
6  words real intreal pred_set set_relation arithmetic
7Libs
8  wordsLib realLib RealArith
9
10val _ = diminish_srw_ss ["RMULCANON","RMULRELNORM","NORMEQ"]
11
12val _ = temp_delsimps ["lift_disj_eq", "lift_imp_disj"]
13
14local
15   open String
16   val mesg_to_string = !Feedback.MESG_to_string
17   fun f s = if isPrefix "mk_functional" s andalso isSubstring "completion" s
18                then ""
19             else mesg_to_string s
20in
21   val () = Feedback.set_trace "Theory.save_thm_reporting" 0
22   val () = Feedback.MESG_to_string := f
23end
24
25Overload tc[local] = “transitive_closure”
26
27(* ------------------------------------------------------------------------
28   Binary floating point representation
29   ------------------------------------------------------------------------ *)
30
31Datatype:
32   float = <| Sign : word1; Exponent : 'w word; Significand : 't word |>
33End
34Overload fsign = “λa. -1 pow w2n a.Sign”
35Overload sign[local] = “fsign”
36
37(* ------------------------------------------------------------------------
38   Maps to other representations
39   ------------------------------------------------------------------------ *)
40
41Datatype: float_value = Float real | Infinity | NaN
42End
43
44Overload precision[local] = “fcp$dimindex”
45Overload bias[local] = “words$INT_MAX”
46
47Definition float_to_real_def[nocompute]:
48   float_to_real (x: ('t, 'w) float) =
49      if x.Exponent = 0w
50         then -1r pow (w2n x.Sign) *
51              (2r / 2r pow (bias (:'w))) *
52              (&(w2n x.Significand) / 2r pow (precision (:'t)))
53      else -1r pow (w2n x.Sign) *
54           (2r pow (w2n x.Exponent) / 2r pow (bias (:'w))) *
55           (1r + &(w2n x.Significand) / 2r pow (precision (:'t)))
56End
57
58Definition float_value_def[nocompute]:
59   float_value (x: ('t, 'w) float) =
60      if x.Exponent = UINT_MAXw
61         then if x.Significand = 0w then Infinity else NaN
62      else Float (float_to_real x)
63End
64
65Theorem FINITE_floatsets[simp]:
66  !s : ('a,'b) float set. FINITE s
67Proof
68  gen_tac >>
69  irule SUBSET_FINITE_I >>
70  irule_at Any SUBSET_UNIV >>
71  qabbrev_tac ‘f = λa:('a,'b)float. (a.Sign, a.Significand, a.Exponent)’ >>
72  ‘!a1 a2. (f a1 = f a2) <=> (a1 = a2)’
73    by (simp[Abbr‘f’, theorem "float_component_equality"] >> metis_tac[]) >>
74  drule INJECTIVE_IMAGE_FINITE >>
75  disch_then (fn th => REWRITE_TAC [GSYM th]) >>
76  ‘IMAGE f UNIV = UNIV’
77    suffices_by (disch_then SUBST_ALL_TAC >> simp[]) >>
78  simp[EXTENSION, Abbr‘f’, pairTheory.FORALL_PROD] >>
79  qx_genl_tac[‘sn’, ‘m’, ‘e’] >>
80  qexists_tac ‘<| Sign := sn; Significand := m; Exponent := e|>’ >>
81  simp[]
82QED
83
84(* ------------------------------------------------------------------------
85   Tests
86   ------------------------------------------------------------------------ *)
87
88Definition float_is_nan_def[nocompute]:
89   float_is_nan (x: ('t, 'w) float) =
90      case float_value x of
91         NaN => T
92       | _ => F
93End
94
95Definition float_is_signalling_def[nocompute]:
96   float_is_signalling (x: ('t, 'w) float) <=>
97      float_is_nan x /\ ~word_msb x.Significand
98End
99
100Definition float_is_infinite_def[nocompute]:
101   float_is_infinite (x: ('t, 'w) float) =
102      case float_value x of
103         Infinity => T
104       | _ => F
105End
106
107Definition float_is_normal_def[nocompute]:
108   float_is_normal (x: ('t, 'w) float) <=>
109      x.Exponent <> 0w /\ x.Exponent <> UINT_MAXw
110End
111
112Definition float_is_subnormal_def[nocompute]:
113   float_is_subnormal (x: ('t, 'w) float) <=>
114      (x.Exponent = 0w) /\ x.Significand <> 0w
115End
116
117Definition float_is_zero_def[nocompute]:
118   float_is_zero (x: ('t, 'w) float) =
119      case float_value x of
120         Float r => r = 0
121       | _ => F
122End
123
124Definition float_is_finite_def[nocompute]:
125   float_is_finite (x: ('t, 'w) float) =
126      case float_value x of
127         Float _ => T
128       | _ => F
129End
130
131Theorem float_is_finite_thm:
132  float_is_finite f ⇔ ∃r. float_value f = Float r
133Proof
134  simp[float_is_finite_def] >> Cases_on ‘float_value f’ >> simp[]
135QED
136
137Definition is_integral_def[nocompute]:   is_integral r = ?n. abs r = &(n:num)
138End
139
140Definition float_is_integral_def[nocompute]:
141   float_is_integral (x: ('t, 'w) float) =
142      case float_value x of
143         Float r => is_integral r
144       | _ => F
145End
146
147(* ------------------------------------------------------------------------
148   Abs and Negate
149   (On some architectures the signalling behaviour changes from IEEE754:1985
150    and IEEE754:2008)
151   ------------------------------------------------------------------------ *)
152
153Definition float_negate_def[nocompute]:
154   float_negate (x: ('t, 'w) float) = x with Sign := ~x.Sign
155End
156
157Definition float_abs_def[nocompute]:
158   float_abs (x: ('t, 'w) float) = x with Sign := 0w
159End
160
161(* ------------------------------------------------------------------------
162   Some constants
163   ------------------------------------------------------------------------ *)
164
165Definition float_plus_infinity_def[nocompute]:
166   float_plus_infinity (:'t # 'w) =
167      <| Sign := 0w;
168         Exponent := UINT_MAXw: 'w word;
169         Significand := 0w: 't word |>
170End
171
172Definition float_plus_zero_def[nocompute]:
173   float_plus_zero (:'t # 'w) =
174      <| Sign := 0w;
175         Exponent := 0w: 'w word;
176         Significand := 0w: 't word |>
177End
178
179Definition float_top_def[nocompute]:
180   float_top (:'t # 'w) =
181      <| Sign := 0w;
182         Exponent := UINT_MAXw - 1w: 'w word;
183         Significand := UINT_MAXw: 't word |>
184End
185Overload FLT_MAX = “float_top(:'a # 'b)”
186
187Definition float_plus_min_def[nocompute]:
188   float_plus_min (:'t # 'w) =
189      <| Sign := 0w;
190         Exponent := 0w: 'w word;
191         Significand := 1w: 't word |>
192End
193
194Definition float_minus_infinity_def[nocompute]:
195   float_minus_infinity (:'t # 'w) =
196      float_negate (float_plus_infinity (:'t # 'w))
197End
198
199Definition float_minus_zero_def[nocompute]:
200   float_minus_zero (:'t # 'w) = float_negate (float_plus_zero (:'t # 'w))
201End
202
203Definition float_bottom_def[nocompute]:
204   float_bottom (:'t # 'w) = float_negate (float_top (:'t # 'w))
205End
206
207Definition float_minus_min_def[nocompute]:
208   float_minus_min (:'t # 'w) = float_negate (float_plus_min (:'t # 'w))
209End
210
211Overload POS0 = “float_plus_zero(:'a#'b)”
212Overload NEG0 = “float_minus_zero(:'a#'b)”
213
214(* ------------------------------------------------------------------------
215   Rounding reals to floating-point values
216   ------------------------------------------------------------------------ *)
217
218Datatype:
219   flags = <| DivideByZero : bool
220            ; InvalidOp : bool
221            ; Overflow : bool
222            ; Precision : bool
223            ; Underflow_BeforeRounding : bool
224            ; Underflow_AfterRounding : bool
225            |>
226End
227
228Definition clear_flags_def[nocompute]:
229  clear_flags = <| DivideByZero := F
230                 ; InvalidOp := F
231                 ; Overflow := F
232                 ; Precision := F
233                 ; Underflow_BeforeRounding := F
234                 ; Underflow_AfterRounding := F
235                 |>
236End
237
238Definition invalidop_flags_def[nocompute]:
239  invalidop_flags = clear_flags with InvalidOp := T
240End
241
242Definition dividezero_flags_def[nocompute]:
243  dividezero_flags = clear_flags with DivideByZero := T
244End
245
246Datatype:
247   rounding = roundTiesToEven
248            | roundTowardPositive
249            | roundTowardNegative
250            | roundTowardZero
251End
252
253Definition is_closest_def:
254   is_closest s x a <=>
255      a IN s /\
256      !b. b IN s ==> abs (float_to_real a - x) <= abs (float_to_real b - x)
257End
258
259Theorem is_closest_exists:
260  s <> {} ==> ?a. is_closest s r a
261Proof
262  simp[is_closest_def] >> rw[] >>
263  qabbrev_tac ‘dists = { d | ?b. b IN s /\ (d = abs (float_to_real b - r)) }’ >>
264  ‘FINITE dists’
265    by (‘?f. dists = IMAGE f s’ suffices_by
266          (rw[] >> irule IMAGE_FINITE >>
267           irule SUBSET_FINITE_I >> qexists_tac ‘UNIV’ >>
268           simp[]) >>
269        qexists_tac ‘λfl. abs (float_to_real fl - r)’ >>
270        simp[Abbr‘dists’, EXTENSION] >> metis_tac[]) >>
271  ‘dists <> {}’ by simp[Abbr‘dists’, EXTENSION, MEMBER_NOT_EMPTY] >>
272  ‘?d. d IN minimal_elements dists (UNCURRY $<)’
273    by (irule finite_acyclic_has_minimal >>
274        simp[acyclic_def] >>
275        ‘tc (UNCURRY ($< : real -> real -> bool)) = UNCURRY $<’
276          suffices_by simp[IN_DEF] >>
277        irule transitive_tc >>
278        simp[transitive_def, SF realSimps.REAL_ARITH_ss]) >>
279  pop_assum mp_tac >> simp[minimal_elements_def, Abbr‘dists’] >>
280  strip_tac >> rename [‘minfl IN s’, ‘d = abs (float_to_real minfl - r)’] >>
281  qexists_tac ‘minfl’ >> rw[] >>
282  metis_tac[REAL_NOT_LT, REAL_LT_REFL]
283QED
284
285Theorem zeroes_are_finite_floats[simp]:
286  float_is_finite (float_plus_zero (:'w # 't)) /\
287  float_is_finite (float_minus_zero (:'w # 't))
288Proof
289  simp[float_is_finite_def, float_plus_zero_def, float_minus_zero_def,
290       float_value_def, float_negate_def]
291QED
292
293Theorem float_to_real_zeroes[simp]:
294  (float_to_real (float_plus_zero (:'w # 't)) = 0) /\
295  (float_to_real (float_minus_zero (:'w # 't)) = 0)
296Proof
297  simp[float_to_real_def, float_plus_zero_def, float_negate_def,
298       float_minus_zero_def]
299QED
300
301Theorem float_to_real_EQ0:
302  (float_to_real (f : ('w,'t) float) = 0) <=>
303  (f = float_plus_zero (:'w # 't)) \/ (f = float_minus_zero (:'w # 't))
304Proof
305  simp[EQ_IMP_THM, DISJ_IMP_THM] >>
306  simp[float_to_real_def, AllCaseEqs(), REAL_DIV_ZERO] >> strip_tac
307  >- (simp[theorem "float_component_equality", float_plus_zero_def,
308           float_minus_zero_def, float_negate_def] >>
309      Cases_on ‘f.Sign’ using wordsTheory.ranged_word_nchotomy >>
310      gs[wordsTheory.word_eq_n2w, bitTheory.MOD_2EXP_MAX_def,
311         bitTheory.MOD_2EXP_def, bitTheory.MOD_2EXP_EQ_def]) >>
312  gs[REAL_ARITH “(1r + x = 0) <=> (x = -1)”,
313     SF realSimps.RMULRELNORM_ss, real_div] >>
314  Cases_on ‘f.Significand’ using wordsTheory.ranged_word_nchotomy >>
315  gs[REAL_OF_NUM_POW]
316QED
317
318Theorem is_closestP_finite_float_exists:
319  ?a : ('w,'t) float. is_closest float_is_finite r a /\
320                      !b. is_closest float_is_finite r b /\ P b ==> P a
321Proof
322  qabbrev_tac ‘cands = { a : ('w,'t) float | is_closest float_is_finite r a}’ >>
323  qabbrev_tac ‘candsP = { a | a IN cands /\ P a }’ >>
324  ‘cands <> {}’
325    by (simp[Abbr‘cands’, EXTENSION] >>
326        irule is_closest_exists >> simp[EXTENSION, IN_DEF] >>
327        irule_at Any (cj 1 zeroes_are_finite_floats)) >>
328  gs[GSYM MEMBER_NOT_EMPTY] >>
329  rename [‘c IN cands’] >>
330  Cases_on ‘candsP = {}’
331  >- (qexists_tac ‘c’ >> fs[Abbr‘cands’, Abbr‘candsP’] >>
332      fs[EXTENSION] >> metis_tac[]) >>
333  gs[GSYM MEMBER_NOT_EMPTY] >>
334  rename1 ‘cp IN candsP’ >> qexists_tac ‘cp’ >>
335  fs[Abbr‘cands’, Abbr‘candsP’]
336QED
337
338Theorem is_closest_float_is_finite_0:
339  is_closest float_is_finite 0 (f:('w,'t)float) <=>
340  (f = float_plus_zero (:'w#'t)) \/ (f = float_minus_zero(:'w#'t))
341Proof
342  eq_tac >> simp[is_closest_def, IN_DEF, DISJ_IMP_THM] >> rw[] >>
343  first_x_assum $ qspec_then ‘float_plus_zero (:'w # 't)’ mp_tac>>
344  simp[REAL_ABS_LE0, float_to_real_EQ0]
345QED
346
347Definition closest_such_def[nocompute]:
348   closest_such p s x =
349      @a. is_closest s x a /\ (!b. is_closest s x b /\ p b ==> p a)
350End
351
352Definition closest_def[nocompute]:   closest = closest_such (K T)
353End
354
355Definition largest_def[nocompute]:
356   largest (:'t # 'w) =
357      (2r pow (UINT_MAX (:'w) - 1) / 2r pow (INT_MAX (:'w))) *
358      (2r - inv (2r pow dimindex(:'t)))
359End
360
361Definition threshold_def[nocompute]:
362   threshold (:'t # 'w) =
363      (2r pow (UINT_MAX (:'w) - 1) / 2r pow (INT_MAX (:'w))) *
364      (2r - inv (2r pow SUC (dimindex(:'t))))
365End
366
367(* Unit in the Last Place (of least precision) *)
368
369(* For a given exponent (applies when significand is not zero) *)
370
371Definition ULP_def[nocompute]:
372   ULP (e:'w word, (:'t)) =
373   2 pow (if e = 0w then 1 else w2n e) / 2 pow (bias (:'w) + precision (:'t))
374End
375
376(* Smallest ULP *)
377
378Definition ulp_def[nocompute]:   ulp (:'t # 'w) = ULP (0w:'w word, (:'t))
379End
380
381Theorem ULP_positive[simp]:
382  0 < ULP (e, i) /\ 0 <= ULP (e, i) /\ ~(ULP (e,i) < 0) /\ ~(ULP (e,i) <= 0)
383Proof
384  csimp[REAL_LE_LT, REAL_NOT_LE, REAL_NOT_LT] >>
385  Cases_on ‘i’ >>
386  simp[ULP_def, REAL_LT_RDIV_0, REAL_OF_NUM_POW]
387QED
388
389Theorem ULP_nonzero[simp]:
390  ULP (e : 'w word, (:'t)) <> 0
391Proof
392  metis_tac[ULP_positive, REAL_LT_REFL]
393QED
394
395Theorem ulp_positive[simp]:
396  0 < ulp(:'t # 'w) /\ 0 <= ulp(:'t # 'w) /\ ~(ulp(:'t#'w) < 0) /\
397  ~(ulp(:'t # 'w) <= 0)
398Proof
399  simp[ulp_def]
400QED
401
402Theorem ulp_nonzero[simp]:
403  ulp (:'t # 'w) <> 0
404Proof
405  simp[ulp_def]
406QED
407
408
409(* rounding *)
410
411Definition round_def[nocompute]:
412   round mode (x: real) =
413   case mode of
414      roundTiesToEven =>
415        let t = threshold (:'t # 'w) in
416          if x <= -t
417             then float_minus_infinity (:'t # 'w)
418          else if x >= t
419             then float_plus_infinity (:'t # 'w)
420          else closest_such (\a. ~word_lsb a.Significand) float_is_finite x
421    | roundTowardZero =>
422        let t = largest (:'t # 'w) in
423          if x < -t
424             then float_bottom (:'t # 'w)
425          else if x > t
426             then float_top (:'t # 'w)
427          else closest
428                 {a | float_is_finite a /\ abs (float_to_real a) <= abs x} x
429    | roundTowardPositive =>
430        let t = largest (:'t # 'w) in
431          if x < -t
432             then float_bottom (:'t # 'w)
433          else if x > t
434             then float_plus_infinity (:'t # 'w)
435          else closest {a | float_is_finite a /\ float_to_real a >= x} x
436    | roundTowardNegative =>
437        let t = largest (:'t # 'w) in
438          if x < -t
439             then float_minus_infinity (:'t # 'w)
440          else if x > t
441             then float_top (:'t # 'w)
442          else closest {a | float_is_finite a /\ float_to_real a <= x} x
443End
444
445Definition integral_round_def[nocompute]:
446   integral_round mode (x: real) =
447   case mode of
448      roundTiesToEven =>
449        let t = threshold (:'t # 'w) in
450          if x <= -t
451             then float_minus_infinity (:'t # 'w)
452          else if x >= t
453             then float_plus_infinity (:'t # 'w)
454          else closest_such (\a. ?n. EVEN n /\ (abs (float_to_real a) = &n))
455                 float_is_integral x
456    | roundTowardZero =>
457        let t = largest (:'t # 'w) in
458          if x < -t
459             then float_bottom (:'t # 'w)
460          else if x > t
461             then float_top (:'t # 'w)
462          else closest
463                 {a | float_is_integral a /\ abs (float_to_real a) <= abs x} x
464    | roundTowardPositive =>
465        let t = largest (:'t # 'w) in
466          if x < -t
467             then float_bottom (:'t # 'w)
468          else if x > t
469             then float_plus_infinity (:'t # 'w)
470          else closest {a | float_is_integral a /\ float_to_real a >= x} x
471    | roundTowardNegative =>
472        let t = largest (:'t # 'w) in
473          if x < -t
474             then float_minus_infinity (:'t # 'w)
475          else if x > t
476             then float_top (:'t # 'w)
477          else closest {a | float_is_integral a /\ float_to_real a <= x} x
478End
479
480(* ------------------------------------------------------------------------
481   NaNs
482   ------------------------------------------------------------------------ *)
483
484Datatype:
485   fp_op =
486     FP_Sqrt rounding (('t, 'w) float)
487   | FP_Add rounding (('t, 'w) float) (('t, 'w) float)
488   | FP_Sub rounding (('t, 'w) float) (('t, 'w) float)
489   | FP_Mul rounding (('t, 'w) float) (('t, 'w) float)
490   | FP_Div rounding (('t, 'w) float) (('t, 'w) float)
491   | FP_MulAdd rounding (('t, 'w) float) (('t, 'w) float) (('t, 'w) float)
492   | FP_MulSub rounding (('t, 'w) float) (('t, 'w) float) (('t, 'w) float)
493End
494
495Definition float_some_qnan_def[nocompute]:
496   float_some_qnan (fp_op : ('t, 'w) fp_op) =
497   (@f. let qnan = f fp_op in float_is_nan qnan /\ ~float_is_signalling qnan)
498     fp_op : ('t, 'w) float
499End
500
501(* ------------------------------------------------------------------------
502   Some arithmetic operations
503   ------------------------------------------------------------------------ *)
504
505(* Round, choosing between -0.0 or +0.0 *)
506
507Definition float_round_def[nocompute]:
508   float_round mode toneg r =
509      let x = round mode r in
510         if float_is_zero x
511            then if toneg
512                    then float_minus_zero (:'t # 'w)
513                 else float_plus_zero (:'t # 'w)
514         else x
515End
516
517Definition float_round_with_flags_def[nocompute]:
518  float_round_with_flags mode to_neg r =
519  let x = float_round mode to_neg r : ('t, 'w) float and a = abs r in
520  let inexact = (float_value x <> Float r) in
521    ((clear_flags with
522        <| Overflow := (float_is_infinite x \/ 2 pow (INT_MIN (:'w)) <= a)
523         (* IEEE-754 permits a number of ways to detect underflow. Below
524            are two possible methods. *)
525         ; Underflow_BeforeRounding := (inexact /\ a < 2 / 2 pow (bias(:'w)))
526         ; Underflow_AfterRounding :=
527             (inexact /\
528              ((float_round mode to_neg r : ('t, 'w + 1) float).Exponent <=+
529               n2w (INT_MIN (:'w))))
530         ; Precision := inexact
531         |>), x)
532End
533
534Definition check_for_signalling_def[nocompute]:
535  check_for_signalling l =
536  clear_flags with InvalidOp := EXISTS float_is_signalling l
537End
538
539Definition real_to_float_def[nocompute]:
540   real_to_float m = float_round m (m = roundTowardNegative)
541End
542
543Definition real_to_float_with_flags_def[nocompute]:
544   real_to_float_with_flags m =
545   float_round_with_flags m (m = roundTowardNegative)
546End
547
548Definition float_round_to_integral_def[nocompute]:
549   float_round_to_integral mode (x: ('t, 'w) float) =
550      case float_value x of
551         Float r => integral_round mode r
552       | _ => x
553End
554
555Definition float_to_int_def[nocompute]:
556   float_to_int mode (x: ('t, 'w) float) =
557   case float_value x of
558      Float r =>
559       SOME (case mode of
560                roundTiesToEven =>
561                  let f = INT_FLOOR r in
562                  let df = abs (r - real_of_int f) in
563                  if (df < 1r / 2) \/ (df = 1r / 2) /\ EVEN (Num (ABS f)) then
564                    f
565                  else
566                    INT_CEILING r
567              | roundTowardPositive => INT_CEILING r
568              | roundTowardNegative => INT_FLOOR r
569              | roundTowardZero =>
570                  if x.Sign = 1w then INT_CEILING r else INT_FLOOR r)
571    | _ => NONE
572End
573
574Definition float_sqrt_def:
575   float_sqrt mode (x: ('t, 'w) float) =
576      if x.Sign = 0w then
577         case float_value x of
578            NaN => (check_for_signalling [x], float_some_qnan (FP_Sqrt mode x))
579          | Infinity => (clear_flags, float_plus_infinity (:'t # 'w))
580          | Float r => (float_round_with_flags mode F (sqrt r))
581      else if x = float_minus_zero (:'t # 'w) then
582        (clear_flags, float_minus_zero (:'t # 'w))
583      else
584        (invalidop_flags, float_some_qnan (FP_Sqrt mode x))
585End
586
587Definition float_add_def[nocompute]:
588   float_add mode (x: ('t, 'w) float) (y: ('t, 'w) float) =
589      case float_value x, float_value y of
590         NaN, _ => (check_for_signalling [x; y],
591                    float_some_qnan (FP_Add mode x y))
592       | _, NaN => (check_for_signalling [y],
593                    float_some_qnan (FP_Add mode x y))
594       | Infinity, Infinity =>
595            if x.Sign = y.Sign then
596               (clear_flags, x)
597            else
598               (invalidop_flags, float_some_qnan (FP_Add mode x y))
599       | Infinity, _ => (clear_flags, x)
600       | _, Infinity => (clear_flags, y)
601       | Float r1, Float r2 =>
602            float_round_with_flags mode
603               (if (r1 = 0) /\ (r2 = 0) /\ (x.Sign = y.Sign) then
604                  x.Sign = 1w
605                else mode = roundTowardNegative) (r1 + r2)
606End
607
608Definition float_sub_def[nocompute]:
609   float_sub mode (x: ('t, 'w) float) (y: ('t, 'w) float) =
610      case float_value x, float_value y of
611         NaN, _ => (check_for_signalling [x; y],
612                    float_some_qnan (FP_Sub mode x y))
613       | _, NaN => (check_for_signalling [y],
614                    float_some_qnan (FP_Sub mode x y))
615       | Infinity, Infinity =>
616            if x.Sign = y.Sign then
617               (invalidop_flags, float_some_qnan (FP_Sub mode x y))
618            else
619               (clear_flags, x)
620       | Infinity, _ => (clear_flags, x)
621       | _, Infinity => (clear_flags, float_negate y)
622       | Float r1, Float r2 =>
623            float_round_with_flags mode
624               (if (r1 = 0) /\ (r2 = 0) /\ x.Sign <> y.Sign then
625                  x.Sign = 1w
626                else mode = roundTowardNegative) (r1 - r2)
627End
628
629Definition float_mul_def[nocompute]:
630   float_mul mode (x: ('t, 'w) float) (y: ('t, 'w) float) =
631      case float_value x, float_value y of
632         NaN, _ => (check_for_signalling [x; y],
633                    float_some_qnan (FP_Mul mode x y))
634       | _, NaN => (check_for_signalling [y],
635                    float_some_qnan (FP_Mul mode x y))
636       | Infinity, Float r =>
637            if r = 0 then
638               (invalidop_flags, float_some_qnan (FP_Mul mode x y))
639            else
640               (clear_flags,
641                if x.Sign = y.Sign then
642                   float_plus_infinity (:'t # 'w)
643                else float_minus_infinity (:'t # 'w))
644       | Float r, Infinity =>
645            if r = 0 then
646               (invalidop_flags, float_some_qnan (FP_Mul mode x y))
647            else
648               (clear_flags,
649                if x.Sign = y.Sign then
650                   float_plus_infinity (:'t # 'w)
651                else float_minus_infinity (:'t # 'w))
652       | Infinity, Infinity =>
653            (clear_flags,
654             if x.Sign = y.Sign then
655                float_plus_infinity (:'t # 'w)
656             else float_minus_infinity (:'t # 'w))
657       | Float r1, Float r2 =>
658            float_round_with_flags mode (x.Sign <> y.Sign) (r1 * r2)
659End
660
661Definition float_div_def[nocompute]:
662   float_div mode (x: ('t, 'w) float) (y: ('t, 'w) float) =
663      case float_value x, float_value y of
664         NaN, _ => (check_for_signalling [x; y],
665                    float_some_qnan (FP_Div mode x y))
666       | _, NaN => (check_for_signalling [y],
667                    float_some_qnan (FP_Div mode x y))
668       | Infinity, Infinity =>
669            (invalidop_flags, float_some_qnan (FP_Div mode x y))
670       | Infinity, _ =>
671            (clear_flags,
672             if x.Sign = y.Sign then
673                float_plus_infinity (:'t # 'w)
674             else float_minus_infinity (:'t # 'w))
675       | _, Infinity =>
676            (clear_flags,
677             if x.Sign = y.Sign then
678                float_plus_zero (:'t # 'w)
679             else float_minus_zero (:'t # 'w))
680       | Float r1, Float r2 =>
681            if r2 = 0
682               then if r1 = 0 then
683                       (invalidop_flags, float_some_qnan (FP_Div mode x y))
684                    else
685                       (dividezero_flags,
686                        if x.Sign = y.Sign then
687                           float_plus_infinity (:'t # 'w)
688                        else float_minus_infinity (:'t # 'w))
689            else float_round_with_flags mode (x.Sign <> y.Sign) (r1 / r2)
690End
691
692Definition float_mul_add_def[nocompute]:
693   float_mul_add mode
694      (x: ('t, 'w) float) (y: ('t, 'w) float) (z: ('t, 'w) float) =
695      let signP = x.Sign ?? y.Sign in
696      let infP = (float_is_infinite x \/ float_is_infinite y)
697      in
698         if float_is_nan x \/ float_is_nan y \/ float_is_nan z then
699            (check_for_signalling [x; y; z],
700             float_some_qnan (FP_MulAdd mode x y z))
701         else if float_is_infinite x /\ float_is_zero y \/
702                 float_is_zero x /\ float_is_infinite y \/
703                 float_is_infinite z /\ infP /\ signP <> z.Sign then
704            (invalidop_flags, float_some_qnan (FP_MulAdd mode x y z))
705         else if float_is_infinite z /\ (z.Sign = 0w) \/ infP /\ (signP = 0w)
706            then (clear_flags, float_plus_infinity (:'t # 'w))
707         else if float_is_infinite z /\ (z.Sign = 1w) \/ infP /\ (signP = 1w)
708            then (clear_flags, float_minus_infinity (:'t # 'w))
709         else
710           let r1 = float_to_real x * float_to_real y ;
711               r2 = float_to_real z ;
712               r = r1 + r2 ;
713           in
714             float_round_with_flags
715               mode
716               ((r = 0) /\
717                (if (r1 = 0) /\ (r2 = 0) /\ (signP = z.Sign) then
718                   signP = 1w
719                 else mode = roundTowardNegative) \/
720                r < 0)
721               r
722End
723
724Definition float_mul_sub_def[nocompute]:
725   float_mul_sub mode
726      (x: ('t, 'w) float) (y: ('t, 'w) float) (z: ('t, 'w) float) =
727      let signP = x.Sign ?? y.Sign in
728      let infP = (float_is_infinite x \/ float_is_infinite y)
729      in
730         if float_is_nan x \/ float_is_nan y \/ float_is_nan z then
731            (check_for_signalling [x; y; z],
732             float_some_qnan (FP_MulSub mode x y z))
733         else if float_is_infinite x /\ float_is_zero y \/
734                 float_is_zero x /\ float_is_infinite y \/
735                 float_is_infinite z /\ infP /\ (signP = z.Sign) then
736            (invalidop_flags, float_some_qnan (FP_MulAdd mode x y z))
737         else if float_is_infinite z /\ (z.Sign = 1w) \/ infP /\ (signP = 0w)
738            then (clear_flags, float_plus_infinity (:'t # 'w))
739         else if float_is_infinite z /\ (z.Sign = 0w) \/ infP /\ (signP = 1w)
740            then (clear_flags, float_minus_infinity (:'t # 'w))
741         else
742            let r1 = float_to_real x * float_to_real y
743            and r2 = float_to_real z
744            in
745              float_round_with_flags mode
746                (if (r1 = 0) /\ (r2 = 0) /\ signP <> z.Sign then
747                   signP = 1w
748                 else mode = roundTowardNegative) (r1 - r2)
749End
750
751(* ------------------------------------------------------------------------
752   Some comparison operations
753   ------------------------------------------------------------------------ *)
754
755Datatype:  float_compare = LT | EQ | GT | UN
756End
757
758Definition float_compare_def[nocompute]:
759   float_compare (x: ('t, 'w) float) (y: ('t, 'w) float) =
760      case float_value x, float_value y of
761         NaN, _ => UN
762       | _, NaN => UN
763       | Infinity, Infinity =>
764            if x.Sign = y.Sign
765               then EQ
766            else if x.Sign = 1w
767               then LT
768            else GT
769       | Infinity, _ => if x.Sign = 1w then LT else GT
770       | _, Infinity => if y.Sign = 1w then GT else LT
771       | Float r1, Float r2 =>
772            if r1 < r2
773               then LT
774            else if r1 = r2
775               then EQ
776            else GT
777End
778
779Definition float_less_than_def[nocompute]:
780   float_less_than (x: ('t, 'w) float) y =
781      (float_compare x y = LT)
782End
783
784Definition float_less_equal_def[nocompute]:
785   float_less_equal (x: ('t, 'w) float) y =
786      case float_compare x y of
787         LT => T
788       | EQ => T
789       | _ => F
790End
791
792Definition float_greater_than_def[nocompute]:
793   float_greater_than (x: ('t, 'w) float) y =
794      (float_compare x y = GT)
795End
796
797Definition float_greater_equal_def[nocompute]:
798   float_greater_equal (x: ('t, 'w) float) y =
799      case float_compare x y of
800         GT => T
801       | EQ => T
802       | _ => F
803End
804
805Definition float_equal_def[nocompute]:
806   float_equal (x: ('t, 'w) float) y =
807      (float_compare x y = EQ)
808End
809
810Definition float_unordered_def[nocompute]:
811   float_unordered (x: ('t, 'w) float) y =
812      (float_compare x y = UN)
813End
814
815Definition exponent_boundary_def[nocompute]:
816   exponent_boundary (y: ('t, 'w) float) (x: ('t, 'w) float) <=>
817      (x.Sign = y.Sign) /\ (w2n x.Exponent = w2n y.Exponent + 1) /\
818      (x.Exponent <> 1w) /\ (y.Significand = -1w) /\ (x.Significand = 0w)
819End
820
821(* ------------------------------------------------------------------------
822   Some arithmetic theorems
823   ------------------------------------------------------------------------ *)
824
825val () = Feedback.set_trace "Theory.save_thm_reporting" 1
826
827val rrw = SRW_TAC [boolSimps.LET_ss, realSimps.REAL_ARITH_ss]
828
829(* |- !n. 0 < 2 pow n *)
830Theorem zero_lt_twopow[simp] =
831   REAL_POW_LT |> Q.SPEC `2` |> SIMP_RULE (srw_ss()) [];
832
833(* |- !n. 0 <= 2 pow n *)
834Theorem zero_le_twopow[simp] =
835   MATCH_MP REAL_LT_IMP_LE (Drule.SPEC_ALL zero_lt_twopow) |> GEN_ALL;
836
837
838Theorem zero_neq_twopow:   !n. 2 pow n <> 0
839Proof simp[]
840QED
841
842Theorem zero_le_pos_div_twopow[simp]:
843   !m n. 0 <= &m / 2 pow n
844Proof
845  rw [REAL_LE_DIV, REAL_LT_IMP_LE]
846QED
847
848Theorem div_eq0[simp]:
849   !a b. 0r < b ==> ((a / b = 0) = (a = 0))
850Proof
851   rw [REAL_EQ_LDIV_EQ]
852QED
853
854(* !b. 2 <= 2 ** b <=> 1 <= b *)
855Theorem exp_ge2[simp] =
856  logrootTheory.LE_EXP_ISO |> Q.SPECL [`2`, `1`] |> reduceLib.REDUCE_RULE
857  |> Conv.GSYM;
858
859
860(* !b. 4 <= 2 ** b <=> 2 <= b *)
861val exp_ge4 =
862  logrootTheory.LE_EXP_ISO
863  |> Q.SPECL [`2`, `2`]
864  |> reduceLib.REDUCE_RULE
865  |> Conv.GSYM
866
867Theorem exp_gt2[simp] =
868  logrootTheory.LT_EXP_ISO
869  |> Q.SPECL [`2`, `1`]
870  |> reduceLib.REDUCE_RULE
871  |> Conv.GSYM
872  ;
873
874(* !n x u. (x / 2 pow n = u / 2 pow n) <=> (x = u) *)
875val div_twopow =
876   eq_rat
877   |> Q.INST [`y` |-> `2 pow n`, `v` |-> `2 pow n`]
878   |> SIMP_RULE (srw_ss()) []
879   |> Q.GENL [`n`, `x`, `u`]
880
881Theorem div_le[local]:
882    !a b c. 0r < a ==> (b / a <= c / a <=> b <= c)
883Proof
884   metis_tac [REAL_LE_LMUL, REAL_DIV_RMUL,
885              REAL_POS_NZ, REAL_MUL_COMM]
886QED
887
888val tac =
889   NTAC 2 strip_tac
890   \\ Cases_on `n <= 2`
891   >- (`(n = 1) \/ (n = 2)` by decide_tac \\ rw [])
892   \\ `2 < n` by decide_tac
893   \\ lrw []
894
895Theorem two_mod_not_one[local]:
896    !n. 0 < n ==> 2 MOD n <> 1
897Proof tac
898QED
899
900Theorem two_mod_eq_zero[local]:
901    !n. 0 < n ==> (2 MOD n = 0 <=> n = 1 \/ n = 2)
902Proof
903   tac
904QED
905
906(* |- !c a. c <> 0 ==> (c * a / c = a) *)
907val mul_cancel =
908   REAL_DIV_LMUL_CANCEL
909   |> Drule.SPEC_ALL
910   |> Q.INST [`b` |-> `1`]
911   |> SIMP_RULE (srw_ss()) []
912   |> Q.GENL [`a`, `c`]
913
914Theorem ge2[local]:
915    dimindex(:'a) <> 1 ==> 2 <= dimindex (:'a)
916Proof
917   rw [DECIDE ``0 < a /\ a <> 1 ==> 2n <= a``]
918QED
919
920Theorem ge2b[local]:
921    !n. 2n <= n ==> 1 <= 2n ** (n - 1) - 1
922Proof
923   REPEAT strip_tac
924   \\ imp_res_tac arithmeticTheory.LESS_EQUAL_ADD
925   \\ simp [arithmeticTheory.EXP_ADD, DECIDE ``0 < n ==> 1n <= 2 * n - 1``]
926QED
927
928Theorem ge2c[local]:
929    !n m. 1r <= n /\ 2 < m ==> 2 < n * m
930Proof
931   rrw [GSYM REAL_LT_LDIV_EQ]
932   \\ `2r / m < 1` by (match_mp_tac REAL_LT_1 \\ simp [])
933   \\ METIS_TAC [REAL_LTE_TRANS]
934QED
935
936Theorem ge1_pow[local]:
937    !a b. 1 <= a /\ 1n <= b ==> a <= a pow b
938Proof
939   Induct_on `b`
940   \\ rw [pow]
941   \\ once_rewrite_tac [REAL_MUL_COMM]
942   \\ `0r < a /\ a <> 0`
943   by METIS_TAC [REAL_ARITH ``1 <= a ==> 0r < a``,
944                 REAL_POS_NZ]
945   \\ simp [GSYM REAL_LE_LDIV_EQ, REAL_DIV_REFL]
946   \\ Cases_on `b = 0`
947   \\ simp []
948   \\ `1 <= b` by decide_tac
949   \\ metis_tac [REAL_LE_TRANS]
950QED
951
952(* |- !n x. 1 < x /\ 1 < n ==> x < x pow n *)
953val gt1_pow =
954   REAL_POW_MONO_LT
955   |> Q.SPEC `1`
956   |> REWRITE_RULE [POW_1]
957
958(* |- !a b. 2 <= a /\ 1 <= b ==> 2 <= a * b *)
959val prod_ge2 =
960   REAL_LE_MUL2
961   |> Q.SPECL [`2`, `a`, `1`, `b`]
962   |> SIMP_RULE (srw_ss()) []
963   |> Q.GENL [`a`, `b`]
964
965Theorem le1[local]:
966    !x y. 0 < y /\ x <= y ==> x / y <= 1r
967Proof
968   REPEAT STRIP_TAC
969   \\ Cases_on `x = y`
970   \\ ASM_SIMP_TAC bool_ss
971        [REAL_LE_REFL, REAL_DIV_REFL,
972         REAL_POS_NZ]
973   \\ ASM_SIMP_TAC bool_ss [REAL_LE_LDIV_EQ, REAL_MUL_LID]
974QED
975
976Theorem le2: !n m. 2r <= n /\ 2 <= m ==> 2 <= n * m
977Proof rrw [prod_ge2]
978QED
979
980Theorem ge4[local]:
981    !n. n <> 0 ==> 4n <= 2 EXP n * 2
982Proof
983   Cases
984   \\ simp [arithmeticTheory.EXP]
985QED
986
987Theorem ge2d[local]:
988    !n m. 2r <= n /\ 2 <= m ==> 2 < n * m
989Proof
990   rrw [GSYM REAL_LT_LDIV_EQ]
991   \\ `2r / m <= 1`
992   by (match_mp_tac le1 \\ ASM_SIMP_TAC (srw_ss()++realSimps.REAL_ARITH_ss) [])
993   \\ imp_res_tac (REAL_ARITH ``a <= 1 ==> a < 2r``)
994   \\ METIS_TAC [REAL_LTE_TRANS]
995QED
996
997(* |- !b. 0 < w2n b <=> b <> 0w *)
998val word_lt0 =
999   wordsTheory.WORD_LO
1000   |> Q.SPEC `0w`
1001   |> REWRITE_RULE [wordsTheory.word_0_n2w, wordsTheory.WORD_LO_word_0]
1002   |> GSYM
1003
1004Theorem word_ge1[local]:
1005    !x. x <> 0w ==> 1 <= w2n x
1006Proof
1007   simp [GSYM word_lt0]
1008QED
1009
1010Theorem not_max_suc_lt_dimword[local]:
1011    !a:'a word. a <> -1w ==> w2n a + 1 < 2 ** dimindex(:'a)
1012Proof
1013   Cases
1014   \\ lrw [wordsTheory.word_eq_n2w, bitTheory.MOD_2EXP_MAX_def,
1015           bitTheory.MOD_2EXP_def, GSYM wordsTheory.dimword_def]
1016QED
1017
1018(* |- !a. a <> 0w ==> 2 <= 2 pow w2n a *)
1019val pow_ge2 =
1020   ge1_pow
1021   |> Q.SPECL [`2`, `w2n (a:'w word)`]
1022   |> SIMP_RULE (srw_ss()) [DECIDE ``1n <= n <=> 0 < n``, word_lt0]
1023   |> GEN_ALL
1024
1025Theorem mult_id[local]:
1026   !a b. 1 < a ==> ((a * b = a) = (b = 1n))
1027Proof
1028  Induct_on `b`
1029  \\ lrw [arithmeticTheory.MULT_CLAUSES]
1030QED
1031
1032(* |- !x y. 1 <= y /\ 0 < x ==> x <= x * y *)
1033val le_mult =
1034   REAL_LE_LDIV_EQ
1035   |> Q.SPECL [`x`, `y`, `x`]
1036   |> Q.DISCH `1 <= y`
1037   |> SIMP_RULE bool_ss [boolTheory.AND_IMP_INTRO, REAL_POS_NZ,
1038                         REAL_DIV_REFL]
1039   |> ONCE_REWRITE_RULE [REAL_MUL_COMM]
1040   |> Q.GENL [`x`, `y`]
1041
1042(* |- !x y. x < 1 /\ 0 < y ==> y * x < y *)
1043val lt_mult =
1044   REAL_LT_RDIV_EQ
1045   |> Q.SPECL [`x`, `y`, `y`]
1046   |> Q.DISCH `x < 1`
1047   |> SIMP_RULE bool_ss [boolTheory.AND_IMP_INTRO, REAL_POS_NZ,
1048                         REAL_DIV_REFL]
1049   |> ONCE_REWRITE_RULE [REAL_MUL_COMM]
1050   |> Q.GENL [`x`, `y`]
1051
1052(*  |- !x y. 1 < y /\ 0 < x ==> x < y * x  *)
1053val gt_mult =
1054   REAL_LT_LDIV_EQ
1055   |> Q.SPECL [`x`, `y`, `x`]
1056   |> Q.DISCH `1 < y`
1057   |> SIMP_RULE bool_ss [boolTheory.AND_IMP_INTRO, REAL_POS_NZ,
1058                         REAL_DIV_REFL]
1059   |> Q.GENL [`x`, `y`]
1060
1061Theorem exp_id[local]:
1062   !a b. 1 < a ==> ((a EXP b = a) = (b = 1))
1063Proof
1064  REPEAT strip_tac
1065  \\ Cases_on `b = 0`
1066  >- lrw [arithmeticTheory.EXP]
1067  \\ Cases_on `b = 1`
1068  >- lrw [arithmeticTheory.EXP]
1069  \\ `1 < b` by decide_tac
1070  \\ imp_res_tac arithmeticTheory.LESS_ADD
1071  \\ pop_assum kall_tac
1072  \\ pop_assum (SUBST1_TAC o REWRITE_RULE [GSYM arithmeticTheory.ADD1] o SYM)
1073  \\ lrw [arithmeticTheory.EXP, mult_id]
1074QED
1075
1076Theorem sub_rat_same_base[local]:
1077    !a b d. 0r < d ==> (a / d - b / d = (a - b) / d)
1078Proof
1079   rrw [REAL_EQ_RDIV_EQ, REAL_SUB_RDISTRIB,
1080        REAL_DIV_RMUL]
1081QED
1082
1083(* |- !x. 0 <= x ==> (abs x = x) *)
1084val gt0_abs =
1085   ABS_REFL
1086   |> Q.SPEC `x`
1087   |> Q.DISCH `0 <= x`
1088   |> SIMP_RULE bool_ss []
1089   |> Drule.GEN_ALL
1090
1091(*
1092(* !z x y. 0 <> z ==> ((x = y) <=> (x * z = y * z)) *)
1093val mul_into =
1094   REAL_EQ_RMUL
1095   |> Drule.SPEC_ALL
1096   |> Q.DISCH `z <> 0`
1097   |> SIMP_RULE std_ss []
1098   |> Conv.GSYM
1099   |> Q.GENL [`y`, `x`, `z`]
1100*)
1101
1102(* ------------------------------------------------------------------------
1103   Some basic theorems
1104   ------------------------------------------------------------------------ *)
1105
1106val rsimp = ASM_SIMP_TAC (srw_ss()++realSimps.REAL_ARITH_ss)
1107val rrw = SRW_TAC [boolSimps.LET_ss, realSimps.REAL_ARITH_ss]
1108val rlfs = full_simp_tac (srw_ss()++realSimps.REAL_ARITH_ss)
1109
1110val float_component_equality = DB.theorem "float_component_equality"
1111
1112val sign_inconsistent =
1113   Drule.GEN_ALL $ wordsLib.WORD_DECIDE “~(x <> -1w /\ x <> 0w: word1) /\
1114                                         ~(x <> 1w /\ x <> 0w)”
1115
1116
1117Theorem sign_neq[local]:
1118    !x. ~x <> x: word1
1119Proof
1120   wordsLib.Cases_word_value
1121   \\ simp []
1122QED
1123
1124Theorem lem[local]:
1125   (1w #>> 1 <> 0w : 'a word) /\ word_msb (1w : 'a word #>> 1)
1126Proof
1127  simp_tac (srw_ss()++wordsLib.WORD_BIT_EQ_ss) []
1128  \\ conj_tac
1129  >| [qexists_tac `dimindex(:'a) - 1`, all_tac]
1130  \\ simp [DECIDE ``0n < n ==> (n - 1 + 1 = n)``, wordsTheory.word_index]
1131QED
1132
1133Theorem some_nan_components[local]:
1134   !fp_op.
1135       ((float_some_qnan fp_op).Exponent = UINT_MAXw) /\
1136       ((float_some_qnan fp_op).Significand <> 0w)
1137Proof
1138   strip_tac \\ simp [float_some_qnan_def]
1139   \\ SELECT_ELIM_TAC
1140   \\ conj_tac
1141   >- (
1142       simp [float_is_nan_def, float_is_signalling_def]
1143       \\ EXISTS_TAC
1144             ``(K (<| Sign := 0w;
1145                      Exponent := UINT_MAXw: 'b word;
1146                      Significand := (1w #>> 1): 'a word |>)) :
1147               ('a, 'b) fp_op -> ('a, 'b) float``
1148       \\ simp [float_value_def, lem]
1149      )
1150   \\ strip_tac
1151   \\ Cases_on `float_value (x fp_op)`
1152   \\ simp [float_is_nan_def]
1153   \\ pop_assum mp_tac
1154   \\ rw [float_value_def]
1155QED
1156
1157Theorem float_components[simp]:
1158    ((float_plus_infinity (:'t # 'w)).Sign = 0w) /\
1159    ((float_plus_infinity (:'t # 'w)).Exponent = UINT_MAXw) /\
1160    ((float_plus_infinity (:'t # 'w)).Significand = 0w) /\
1161    ((float_minus_infinity (:'t # 'w)).Sign = 1w) /\
1162    ((float_minus_infinity (:'t # 'w)).Exponent = UINT_MAXw) /\
1163    ((float_minus_infinity (:'t # 'w)).Significand = 0w) /\
1164    ((float_plus_zero (:'t # 'w)).Sign = 0w) /\
1165    ((float_plus_zero (:'t # 'w)).Exponent = 0w) /\
1166    ((float_plus_zero (:'t # 'w)).Significand = 0w) /\
1167    ((float_minus_zero (:'t # 'w)).Sign = 1w) /\
1168    ((float_minus_zero (:'t # 'w)).Exponent = 0w) /\
1169    ((float_minus_zero (:'t # 'w)).Significand = 0w) /\
1170    ((float_plus_min (:'t # 'w)).Sign = 0w) /\
1171    ((float_plus_min (:'t # 'w)).Exponent = 0w) /\
1172    ((float_plus_min (:'t # 'w)).Significand = 1w) /\
1173    ((float_minus_min (:'t # 'w)).Sign = 1w) /\
1174    ((float_minus_min (:'t # 'w)).Exponent = 0w) /\
1175    ((float_minus_min (:'t # 'w)).Significand = 1w) /\
1176    ((float_top (:'t # 'w)).Sign = 0w) /\
1177    ((float_top (:'t # 'w)).Exponent = UINT_MAXw - 1w) /\
1178    ((float_top (:'t # 'w)).Significand = UINT_MAXw) /\
1179    ((float_bottom (:'t # 'w)).Sign = 1w) /\
1180    ((float_bottom (:'t # 'w)).Exponent = UINT_MAXw - 1w) /\
1181    ((float_bottom (:'t # 'w)).Significand = UINT_MAXw) /\
1182    (!fp_op. (float_some_qnan fp_op).Exponent = UINT_MAXw) /\
1183    (!fp_op. (float_some_qnan fp_op).Significand <> 0w) /\
1184    (!x. (float_negate x).Sign = ~x.Sign) /\
1185    (!x. (float_negate x).Exponent = x.Exponent) /\
1186    (!x. (float_negate x).Significand = x.Significand)
1187Proof
1188   rw [float_plus_infinity_def, float_minus_infinity_def,
1189       float_plus_zero_def, float_minus_zero_def, float_plus_min_def,
1190       float_minus_min_def, float_top_def, float_bottom_def,
1191       some_nan_components, float_negate_def]
1192QED
1193
1194Theorem float_distinct[simp]:
1195    (float_plus_infinity (:'t # 'w) <> float_minus_infinity (:'t # 'w)) /\
1196    (float_plus_infinity (:'t # 'w) <> float_plus_zero (:'t # 'w)) /\
1197    (float_plus_infinity (:'t # 'w) <> float_minus_zero (:'t # 'w)) /\
1198    (float_plus_infinity (:'t # 'w) <> float_top (:'t # 'w)) /\
1199    (float_plus_infinity (:'t # 'w) <> float_bottom (:'t # 'w)) /\
1200    (float_plus_infinity (:'t # 'w) <> float_plus_min (:'t # 'w)) /\
1201    (float_plus_infinity (:'t # 'w) <> float_minus_min (:'t # 'w)) /\
1202    (!fp_op. float_plus_infinity (:'t # 'w) <> float_some_qnan fp_op) /\
1203    (float_minus_infinity (:'t # 'w) <> float_plus_zero (:'t # 'w)) /\
1204    (float_minus_infinity (:'t # 'w) <> float_minus_zero (:'t # 'w)) /\
1205    (float_minus_infinity (:'t # 'w) <> float_top (:'t # 'w)) /\
1206    (float_minus_infinity (:'t # 'w) <> float_bottom (:'t # 'w)) /\
1207    (float_minus_infinity (:'t # 'w) <> float_plus_min (:'t # 'w)) /\
1208    (float_minus_infinity (:'t # 'w) <> float_minus_min (:'t # 'w)) /\
1209    (!fp_op. float_minus_infinity (:'t # 'w) <> float_some_qnan fp_op) /\
1210    (float_plus_zero (:'t # 'w) <> float_minus_zero (:'t # 'w)) /\
1211    (float_plus_zero (:'t # 'w) <> float_top (:'t # 'w)) /\
1212    (float_plus_zero (:'t # 'w) <> float_bottom (:'t # 'w)) /\
1213    (float_plus_zero (:'t # 'w) <> float_plus_min (:'t # 'w)) /\
1214    (float_plus_zero (:'t # 'w) <> float_minus_min (:'t # 'w)) /\
1215    (!fp_op. float_plus_zero (:'t # 'w) <> float_some_qnan fp_op) /\
1216    (float_minus_zero (:'t # 'w) <> float_top (:'t # 'w)) /\
1217    (float_minus_zero (:'t # 'w) <> float_bottom (:'t # 'w)) /\
1218    (float_minus_zero (:'t # 'w) <> float_plus_min (:'t # 'w)) /\
1219    (float_minus_zero (:'t # 'w) <> float_minus_min (:'t # 'w)) /\
1220    (!fp_op. float_minus_zero (:'t # 'w) <> float_some_qnan fp_op) /\
1221    (float_top (:'t # 'w) <> float_minus_min (:'t # 'w)) /\
1222    (float_top (:'t # 'w) <> float_bottom (:'t # 'w)) /\
1223    (!fp_op. float_top (:'t # 'w) <> float_some_qnan fp_op) /\
1224    (float_bottom (:'t # 'w) <> float_plus_min (:'t # 'w)) /\
1225    (!fp_op. float_bottom (:'t # 'w) <> float_some_qnan fp_op) /\
1226    (!fp_op. float_plus_min (:'t # 'w) <> float_some_qnan fp_op) /\
1227    (float_plus_min (:'t # 'w) <> float_minus_min (:'t # 'w)) /\
1228    (!fp_op. float_minus_min (:'t # 'w) <> float_some_qnan fp_op) /\
1229    (!x. float_negate x <> x)
1230Proof
1231   rw [float_component_equality, float_components, two_mod_not_one, sign_neq,
1232       wordsTheory.word_T_not_zero, wordsTheory.WORD_EQ_NEG]
1233QED
1234
1235Theorem float_values[simp]:
1236    (float_value (float_plus_infinity (:'t # 'w)) = Infinity) /\
1237    (float_value (float_minus_infinity (:'t # 'w)) = Infinity) /\
1238    (!fp_op. float_value (float_some_qnan fp_op) = NaN) /\
1239    (float_value (float_plus_zero (:'t # 'w)) = Float 0) /\
1240    (float_value (float_minus_zero (:'t # 'w)) = Float 0) /\
1241    (float_value (float_plus_min (:'t # 'w)) =
1242        Float (2r / (2r pow (bias (:'w) + precision (:'t))))) /\
1243    (float_value (float_minus_min (:'t # 'w)) =
1244        Float (-2r / (2r pow (bias (:'w) + precision (:'t)))))
1245Proof
1246   rw [float_plus_infinity_def, float_minus_infinity_def,
1247       float_plus_zero_def, float_minus_zero_def, float_plus_min_def,
1248       float_minus_min_def, some_nan_components, float_negate_def,
1249       float_value_def, float_to_real_def, wordsTheory.word_T_not_zero,
1250       mult_rat, POW_ADD, GSYM REAL_NEG_MINUS1,
1251       GSYM REAL_MUL_LNEG, neg_rat]
1252QED
1253
1254Theorem zero_to_real[simp]:
1255   (float_to_real (float_plus_zero (:'t # 'w)) = 0) /\
1256   (float_to_real (float_minus_zero (:'t # 'w)) = 0)
1257Proof
1258   rw [float_plus_zero_def, float_minus_zero_def,
1259       float_negate_def, float_to_real_def]
1260QED
1261
1262Theorem sign_not_zero: !s: word1. -1 pow w2n s <> 0
1263Proof wordsLib.Cases_word_value \\ EVAL_TAC
1264QED
1265
1266Theorem float_sets:
1267  (float_is_zero = {float_minus_zero (:'t # 'w);
1268                    float_plus_zero (:'t # 'w)}) /\
1269  (float_is_infinite = {float_minus_infinity (:'t # 'w);
1270                        float_plus_infinity (:'t # 'w)})
1271Proof
1272   rw [FUN_EQ_THM, float_is_infinite_def, float_is_zero_def, float_value_def,
1273       float_plus_infinity_def, float_minus_infinity_def,
1274       float_plus_zero_def, float_minus_zero_def, float_to_real_def,
1275       float_negate_def, float_component_equality]
1276   \\ rw [sign_not_zero, REAL_ARITH ``0r <= n ==> 1 + n <> 0``]
1277   \\ wordsLib.Cases_on_word_value `x.Sign`
1278   \\ simp []
1279QED
1280
1281val tac =
1282   rw [float_is_nan_def, float_is_normal_def, float_is_subnormal_def,
1283       float_is_finite_def, float_is_infinite_def, float_is_zero_def,
1284       float_is_integral_def, is_integral_def, float_values,
1285       some_nan_components]
1286   \\ rw [float_plus_infinity_def, float_minus_infinity_def,
1287          float_plus_zero_def, float_minus_zero_def, float_top_def,
1288          float_bottom_def, float_some_qnan_def, float_plus_min_def,
1289          float_minus_min_def, float_negate_def, float_value_def,
1290          wordsTheory.WORD_EQ_NEG, REAL_EQ_LDIV_EQ, two_mod_not_one]
1291
1292Theorem infinity_properties[simp]:
1293    ~float_is_zero (float_plus_infinity (:'t # 'w)) /\
1294    ~float_is_zero (float_minus_infinity (:'t # 'w)) /\
1295    ~float_is_finite (float_plus_infinity (:'t # 'w)) /\
1296    ~float_is_finite (float_minus_infinity (:'t # 'w)) /\
1297    ~float_is_integral (float_plus_infinity (:'t # 'w)) /\
1298    ~float_is_integral (float_minus_infinity (:'t # 'w)) /\
1299    ~float_is_nan (float_plus_infinity (:'t # 'w)) /\
1300    ~float_is_nan (float_minus_infinity (:'t # 'w)) /\
1301    ~float_is_normal (float_plus_infinity (:'t # 'w)) /\
1302    ~float_is_normal (float_minus_infinity (:'t # 'w)) /\
1303    ~float_is_subnormal (float_plus_infinity (:'t # 'w)) /\
1304    ~float_is_subnormal (float_minus_infinity (:'t # 'w)) /\
1305    float_is_infinite (float_plus_infinity (:'t # 'w)) /\
1306    float_is_infinite (float_minus_infinity (:'t # 'w))
1307Proof tac
1308QED
1309
1310Theorem zero_properties[simp]:
1311    float_is_zero (float_plus_zero (:'t # 'w)) /\
1312    float_is_zero (float_minus_zero (:'t # 'w)) /\
1313    float_is_finite (float_plus_zero (:'t # 'w)) /\
1314    float_is_finite (float_minus_zero (:'t # 'w)) /\
1315    float_is_integral (float_plus_zero (:'t # 'w)) /\
1316    float_is_integral (float_minus_zero (:'t # 'w)) /\
1317    ~float_is_nan (float_plus_zero (:'t # 'w)) /\
1318    ~float_is_nan (float_minus_zero (:'t # 'w)) /\
1319    ~float_is_normal (float_plus_zero (:'t # 'w)) /\
1320    ~float_is_normal (float_minus_zero (:'t # 'w)) /\
1321    ~float_is_subnormal (float_plus_zero (:'t # 'w)) /\
1322    ~float_is_subnormal (float_minus_zero (:'t # 'w)) /\
1323    ~float_is_infinite (float_plus_zero (:'t # 'w)) /\
1324    ~float_is_infinite (float_minus_zero (:'t # 'w))
1325Proof tac
1326QED
1327
1328
1329
1330Theorem some_nan_properties:
1331  !fp_op.
1332      ~float_is_zero (float_some_qnan fp_op) /\
1333      ~float_is_finite (float_some_qnan fp_op) /\
1334      ~float_is_integral (float_some_qnan fp_op) /\
1335      float_is_nan (float_some_qnan fp_op) /\
1336      ~float_is_signalling (float_some_qnan fp_op) /\
1337      ~float_is_normal (float_some_qnan fp_op) /\
1338      ~float_is_subnormal (float_some_qnan fp_op) /\
1339      ~float_is_infinite (float_some_qnan fp_op)
1340Proof
1341   tac
1342   \\ SELECT_ELIM_TAC
1343   \\ simp []
1344   \\ qexists_tac
1345        `(K (<| Sign := 0w;
1346                Exponent := UINT_MAXw: 'b word;
1347                Significand := (1w #>> 1): 'a word |>))`
1348   \\ simp [float_is_signalling_def]
1349   \\ tac
1350   \\ fs [lem]
1351QED
1352
1353Theorem NMUL_EQ_2:
1354  ((m:num) * n = 2) <=> (m = 1) /\ (n = 2) \/ (m = 2) /\ (n = 1)
1355Proof
1356  assume_tac dividesTheory.PRIME_2 >>
1357  gs[dividesTheory.prime_def, dividesTheory.divides_def, Excl "PRIME_2",
1358     PULL_EXISTS] >>
1359  eq_tac >> simp[DISJ_IMP_THM] >>
1360  disch_then (assume_tac o SYM) >> first_x_assum drule >>
1361  simp[DISJ_IMP_THM] >> strip_tac >> gs[]
1362QED
1363
1364Theorem min_properties[simp]:
1365    ~float_is_zero (float_plus_min (:'t # 'w)) /\
1366    float_is_finite (float_plus_min (:'t # 'w)) /\
1367    (float_is_integral (float_plus_min (:'t # 'w)) <=>
1368     (precision(:'w) = 1) /\ (precision(:'t) = 1)) /\
1369    ~float_is_nan (float_plus_min (:'t # 'w)) /\
1370    ~float_is_normal (float_plus_min (:'t # 'w)) /\
1371    float_is_subnormal (float_plus_min (:'t # 'w)) /\
1372    ~float_is_infinite (float_plus_min (:'t # 'w)) /\
1373    ~float_is_zero (float_minus_min (:'t # 'w)) /\
1374    float_is_finite (float_minus_min (:'t # 'w)) /\
1375    (float_is_integral (float_minus_min (:'t # 'w)) <=>
1376     (precision(:'w) = 1) /\ (precision(:'t) = 1)) /\
1377    ~float_is_nan (float_minus_min (:'t # 'w)) /\
1378    ~float_is_normal (float_minus_min (:'t # 'w)) /\
1379    float_is_subnormal (float_minus_min (:'t # 'w)) /\
1380    ~float_is_infinite (float_minus_min (:'t # 'w))
1381Proof
1382   tac (* after this the float_is_integral cases remain *) >>
1383   simp[ABS_DIV, iffRL ABS_REFL, real_div, REAL_POW_ADD,
1384        REAL_INV_MUL', SF realSimps.RMULRELNORM_ss] >>
1385   simp[REAL_OF_NUM_POW, wordsTheory.INT_MAX_def,
1386        wordsTheory.INT_MIN_def, NMUL_EQ_2] >>
1387   simp[DECIDE “x <= 1n <=> (x = 1) \/ (x = 0)”]
1388QED
1389
1390Theorem lem1[local]:
1391    0w <+ (-2w:'a word) <=> (dimindex(:'a) <> 1)
1392Proof
1393   once_rewrite_tac [wordsTheory.WORD_LESS_NEG_RIGHT]
1394   \\ simp [two_mod_eq_zero, wordsTheory.dimword_def, exp_id,
1395            DECIDE ``0 < n ==> n <> 0n``]
1396QED
1397
1398Theorem lem2[local]:
1399    dimindex(:'a) <> 1 ==> -2w <+ (-1w:'a word)
1400Proof
1401   once_rewrite_tac [wordsTheory.WORD_LESS_NEG_RIGHT]
1402   \\ simp [two_mod_eq_zero, wordsTheory.dimword_def, exp_id,
1403            DECIDE ``0 < n ==> n <> 0n``, wordsTheory.word_lo_n2w]
1404   \\ strip_tac
1405   \\ `1 < dimindex(:'a)` by simp [DECIDE ``0 < n /\ n <> 1 ==> 1n < n``]
1406   \\ imp_res_tac
1407        (bitTheory.TWOEXP_MONO
1408         |> Q.SPECL [`1`, `dimindex(:'a)`]
1409         |> numLib.REDUCE_RULE)
1410   \\ lrw []
1411QED
1412
1413val tac =
1414   tac
1415   \\ srw_tac[] [float_to_real_def, two_mod_eq_zero, wordsTheory.dimword_def,
1416                 REAL_ARITH ``0r <= n ==> 1 + n <> 0``, exp_id, lem1,
1417                 DECIDE ``0 < n ==> n <> 0n``]
1418   \\ Cases_on `dimindex(:'w) = 1`
1419   \\ lrw [lem2]
1420
1421Theorem top_properties[simp]:
1422    ~float_is_zero (float_top (:'t # 'w)) /\
1423    float_is_finite (float_top (:'t # 'w)) /\
1424    (* float_is_integral (float_top (:'t # 'w)) = ?? /\ *)
1425    ~float_is_nan (float_top (:'t # 'w)) /\
1426    (float_is_normal (float_top (:'t # 'w)) = (precision(:'w) <> 1)) /\
1427    (float_is_subnormal (float_top (:'t # 'w)) = (precision(:'w) = 1)) /\
1428    ~float_is_infinite (float_top (:'t # 'w))
1429Proof tac
1430QED
1431
1432Theorem bottom_properties[simp]:
1433    ~float_is_zero (float_bottom (:'t # 'w)) /\
1434    float_is_finite (float_bottom (:'t # 'w)) /\
1435    (* float_is_integral (float_bottom (:'t # 'w)) = ?? /\ *)
1436    ~float_is_nan (float_bottom (:'t # 'w)) /\
1437    (float_is_normal (float_bottom (:'t # 'w)) = (precision(:'w) <> 1)) /\
1438    (float_is_subnormal (float_bottom (:'t # 'w)) = (precision(:'w) = 1)) /\
1439    ~float_is_infinite (float_bottom (:'t # 'w))
1440Proof tac
1441QED
1442
1443Theorem float_is_zero:
1444  !x. float_is_zero x <=> (x.Exponent = 0w) /\ (x.Significand = 0w)
1445Proof
1446   rw [float_is_zero_def, float_value_def, float_to_real_def, sign_not_zero,
1447       REAL_ARITH ``0 <= x ==> 1 + x <> 0r``,
1448       REAL_LE_DIV, REAL_LT_IMP_LE]
1449QED
1450
1451Theorem float_is_finite:
1452  !x. float_is_finite x <=>
1453      float_is_normal x \/ float_is_subnormal x \/ float_is_zero x
1454Proof
1455  rw [float_is_finite_def, float_is_normal_def, float_is_subnormal_def,
1456      float_is_zero, float_value_def]
1457  \\ Cases_on `x.Exponent = 0w`
1458  \\ Cases_on `x.Significand = 0w`
1459  \\ fs []
1460QED
1461
1462Theorem float_cases_finite:
1463  !x. float_is_nan x \/ float_is_infinite x \/ float_is_finite x
1464Proof
1465  rw [float_is_nan_def, float_is_infinite_def, float_is_finite_def]
1466  \\ Cases_on `float_value x`
1467  \\ fs []
1468QED
1469
1470Theorem float_distinct_finite:
1471  !x. ~(float_is_nan x /\ float_is_infinite x) /\
1472      ~(float_is_nan x /\ float_is_finite x) /\
1473      ~(float_is_infinite x /\ float_is_finite x)
1474Proof
1475  rw [float_is_nan_def, float_is_infinite_def, float_is_finite_def]
1476  \\ Cases_on `float_value x`
1477  \\ fs []
1478QED
1479
1480Theorem float_cases:
1481  !x. float_is_nan x \/ float_is_infinite x \/ float_is_normal x \/
1482      float_is_subnormal x \/ float_is_zero x
1483Proof metis_tac [float_cases_finite, float_is_finite]
1484QED
1485
1486Theorem float_is_distinct:
1487  !x. ~(float_is_nan x /\ float_is_infinite x) /\
1488      ~(float_is_nan x /\ float_is_normal x) /\
1489      ~(float_is_nan x /\ float_is_subnormal x) /\
1490      ~(float_is_nan x /\ float_is_zero x) /\
1491      ~(float_is_infinite x /\ float_is_normal x) /\
1492      ~(float_is_infinite x /\ float_is_subnormal x) /\
1493      ~(float_is_infinite x /\ float_is_zero x) /\
1494      ~(float_is_normal x /\ float_is_subnormal x) /\
1495      ~(float_is_normal x /\ float_is_zero x) /\
1496      ~(float_is_subnormal x /\ float_is_zero x)
1497Proof
1498  rw []
1499  \\ TRY (metis_tac [float_is_finite, float_distinct_finite])
1500  \\ fs [float_is_normal_def, float_is_subnormal_def, float_is_zero]
1501  \\ Cases_on `x.Exponent = 0w`
1502  \\ Cases_on `x.Exponent = -1w`
1503  \\ Cases_on `x.Significand = 0w`
1504  \\ fs []
1505QED
1506
1507Theorem float_infinities:
1508  !x. float_is_infinite x <=>
1509       (x = float_plus_infinity (:'t # 'w)) \/
1510       (x = float_minus_infinity (:'t # 'w))
1511Proof
1512  strip_tac
1513  \\ Q.ISPEC_THEN `x : ('t, 'w) float` strip_assume_tac float_cases_finite
1514  \\ TRY (metis_tac [float_distinct_finite, infinity_properties])
1515  \\ Cases_on `float_value x`
1516  \\ Cases_on `x.Exponent = -1w`
1517  \\ Cases_on `x.Significand = 0w`
1518  \\ fs [float_is_infinite_def, float_value_def,
1519         float_plus_infinity_def, float_minus_infinity_def,
1520         float_negate_def, float_component_equality]
1521  \\ wordsLib.Cases_on_word_value `x.Sign`
1522  \\ simp []
1523QED
1524
1525Theorem float_infinities_distinct:
1526  !x. ~((x = float_plus_infinity (:'t # 'w)) /\
1527        (x = float_minus_infinity (:'t # 'w)))
1528Proof
1529  simp [float_plus_infinity_def, float_minus_infinity_def,
1530        float_negate_def, float_component_equality]
1531QED
1532(* ------------------------------------------------------------------------ *)
1533
1534Theorem float_to_real_negate:
1535   !x. float_to_real (float_negate x) = -float_to_real x
1536Proof
1537   rw [float_to_real_def, float_negate_def]
1538   \\ wordsLib.Cases_on_word_value `x.Sign`
1539   \\ rsimp []
1540QED
1541
1542Theorem float_negate_negate[simp]:
1543   !x. float_negate (float_negate x) = x
1544Proof
1545   simp [float_negate_def, float_component_equality]
1546QED
1547
1548(* ------------------------------------------------------------------------
1549   Lemma support for rounding theorems
1550   ------------------------------------------------------------------------ *)
1551
1552(*
1553val () = List.app Parse.clear_overloads_on ["bias", "precision"]
1554*)
1555
1556local
1557   val cnv =
1558      Conv.QCONV
1559        (REWRITE_CONV [REAL_LDISTRIB, REAL_RDISTRIB])
1560   val dec =
1561      METIS_PROVE
1562        [REAL_DIV_RMUL, REAL_MUL_COMM,
1563         REAL_MUL_ASSOC, zero_neq_twopow,
1564         mult_ratr
1565         |> Q.INST [`z` |-> `2 pow n`]
1566         |> REWRITE_RULE [zero_neq_twopow]
1567         |> GEN_ALL]
1568in
1569   fun CANCEL_PROVE tm =
1570      let
1571         val thm1 = cnv tm
1572         val tm1 = boolSyntax.rhs (Thm.concl thm1)
1573         val thm2 = Drule.EQT_INTRO (dec tm1)
1574      in
1575         Drule.EQT_ELIM (Thm.TRANS thm1 thm2)
1576      end
1577end
1578
1579Theorem cancel_rwts[local] = Q.prove(
1580   `(!a b. (2 pow a * b) / 2 pow a = b) /\
1581    (!a b c. 2 pow a * (b / 2 pow a * c) = b * c) /\
1582    (!a b c. a * (b / 2 pow c) * 2 pow c = a * b) /\
1583    (!a b c. a * (2 pow b * c) / 2 pow b = a * c) /\
1584    (!a b c. a / 2 pow b * (2 pow b * c) = a * c) /\
1585    (!a b c. a / 2 pow b * c * 2 pow b = a * c) /\
1586    (!a b c d. a / 2 pow b * c * (2 pow b * d) = a * c * d) /\
1587    (!a b c d. a * (2 pow b * c) * d / 2 pow b = a * c * d)`,
1588   metis_tac
1589     [REAL_DIV_RMUL, REAL_MUL_COMM,
1590      REAL_MUL_ASSOC, zero_neq_twopow,
1591      mult_ratr
1592      |> Q.INST [`z` |-> `2 pow n`]
1593      |> REWRITE_RULE [zero_neq_twopow]
1594      |> GEN_ALL]
1595   )
1596
1597val cancel_rwt =
1598   CANCEL_PROVE
1599     ``(!a b c d. a * (b + c / 2 pow d) * 2 pow d = a * b * 2 pow d + a * c)``
1600
1601Theorem ulp: ulp (:'t # 'w) = float_to_real (float_plus_min (:'t # 'w))
1602Proof
1603   simp [ulp_def, ULP_def, float_to_real_def, float_plus_min_def,
1604         mult_rat, GSYM POW_ADD]
1605QED
1606
1607Theorem neg_ulp:
1608   -ulp (:'t # 'w) = float_to_real (float_negate (float_plus_min (:'t # 'w)))
1609Proof simp [float_to_real_negate, ulp]
1610QED
1611
1612Theorem ULP_gt0[local]:
1613    !e. 0 < ULP (e:'w word, (:'t))
1614Proof
1615   rw [ULP_def, REAL_LT_RDIV_0]
1616QED
1617
1618val ulp_gt0 = (REWRITE_RULE [GSYM ulp_def] o Q.SPEC `0w`) ULP_gt0
1619
1620Theorem ULP_le_mono:
1621   !e1 e2. e2 <> 0w ==> (ULP (e1, (:'t)) <= ULP (e2, (:'t)) <=> e1 <=+ e2)
1622Proof
1623   Cases
1624   \\ Cases
1625   \\ lrw [ULP_def, wordsTheory.word_ls_n2w, div_le]
1626   \\ simp [REAL_OF_NUM_POW]
1627QED
1628
1629Theorem ulp_lt_ULP: !e: 'w word. ulp (:'t # 'w) <= ULP (e,(:'t))
1630Proof rw [ulp_def] \\ Cases_on `e = 0w` \\ simp [ULP_le_mono]
1631QED
1632
1633Theorem lem[local]:
1634   !n. 0 < n ==> 3 < 2 pow SUC n
1635Proof
1636   Induct
1637   \\ rw [Once pow]
1638   \\ Cases_on `0 < n`
1639   \\ simp [DECIDE ``~(0n < n) ==> (n = 0)``,
1640            REAL_ARITH ``3r < n ==> 3 < 2 * n``]
1641QED
1642
1643Theorem ulp_lt_largest: ulp (:'t # 'w) < largest (:'t # 'w)
1644Proof
1645   simp [ulp_def, ULP_def, largest_def, REAL_LT_RMUL_0, cancel_rwts,
1646         REAL_LT_LDIV_EQ, POW_ADD]
1647   \\ simp [REAL_SUB_RDISTRIB, REAL_SUB_LDISTRIB,
1648            REAL_MUL_LINV, GSYM realaxTheory.REAL_MUL_ASSOC,
1649            GSYM (CONJUNCT2 pow)]
1650   \\ simp [REAL_ARITH ``(a * b) - a = a * (b - 1): real``]
1651   \\ match_mp_tac ge2c
1652   \\ rw [GSYM REAL_LT_ADD_SUB, POW_2_LE1, lem]
1653QED
1654
1655Theorem ulp_lt_threshold:
1656   ulp (:'t # 'w) < threshold (:'t # 'w)
1657Proof
1658   simp [ulp_def, ULP_def, threshold_def, REAL_LT_RMUL_0,
1659         cancel_rwts, REAL_LT_LDIV_EQ, POW_ADD,
1660         pow]
1661   \\ simp [REAL_SUB_RDISTRIB, REAL_SUB_LDISTRIB,
1662            REAL_MUL_LINV, REAL_INV_MUL,
1663            GSYM realaxTheory.REAL_MUL_ASSOC]
1664   \\ simp [REAL_ARITH ``(a * b) - a * c = a * (b - c): real``]
1665   \\ match_mp_tac ge2c
1666   \\ rw [POW_2_LE1, GSYM REAL_LT_ADD_SUB,
1667          REAL_INV_1OVER, GSYM (CONJUNCT2 pow),
1668          REAL_LT_LDIV_EQ, lem,
1669          REAL_ARITH ``3r < n ==> 5 < n * 2``]
1670QED
1671
1672Theorem lt_ulp_not_infinity0[local] =
1673   MATCH_MP
1674      (REAL_ARITH ``u < l ==> abs x < u ==> ~(x < -l) /\ ~(x > l)``)
1675      ulp_lt_largest
1676   |> Drule.GEN_ALL
1677
1678Theorem lt_ulp_not_infinity1[local] =
1679   MATCH_MP
1680      (REAL_ARITH
1681         ``u < l ==> 2 * abs x <= u ==> ~(x <= -l) /\ ~(x >= l)``)
1682      ulp_lt_threshold
1683   |> Drule.GEN_ALL
1684
1685Theorem abs_float_value:
1686    (!b: word1 c d. abs (-1 pow w2n b * c * d) = abs (c * d)) /\
1687    (!b: word1 c. abs (-1 pow w2n b * c) = abs c)
1688Proof
1689   conj_tac
1690   \\ wordsLib.Cases_word_value
1691   \\ simp [ABS_MUL]
1692QED
1693
1694(* |- !x n. abs (1 + &n / 2 pow x) = 1 + &n / 2 pow x *)
1695Theorem abs_significand[local] =
1696   REAL_ARITH ``!a b. 0 <= a /\ 0 <= b ==> (abs (a + b) = a + b)``
1697   |> Q.SPECL [`1`, `&n / 2 pow x`]
1698   |> Conv.CONV_RULE
1699        (Conv.RATOR_CONV (SIMP_CONV (srw_ss()++realSimps.REAL_ARITH_ss)
1700            [REAL_LE_DIV, REAL_POS,
1701             REAL_LT_IMP_LE]))
1702   |> REWRITE_RULE []
1703   |> GEN_ALL
1704
1705Theorem less_than_ulp:
1706    !a: ('t, 'w) float.
1707       abs (float_to_real a) < ulp (:'t # 'w) <=>
1708       (a.Exponent = 0w) /\ (a.Significand = 0w)
1709Proof
1710   strip_tac
1711   \\ eq_tac
1712   \\ rw [ulp_def, ULP_def, float_to_real_def, abs_float_value, abs_significand,
1713          ABS_MUL, ABS_DIV, ABS_N,
1714          gt0_abs, mult_rat, REAL_LT_RDIV_0]
1715   >| [
1716      SPOSE_NOT_THEN strip_assume_tac
1717      \\ qpat_x_assum `x < y: real` mp_tac
1718      \\ simp [REAL_NOT_LT, GSYM POW_ADD,
1719               REAL_LE_RDIV_EQ, REAL_DIV_RMUL]
1720      \\ Cases_on `a.Significand`
1721      \\ lfs [],
1722      (* -- *)
1723      simp [REAL_NOT_LT, REAL_LDISTRIB,
1724            REAL_LE_LDIV_EQ, REAL_RDISTRIB,
1725            POW_ADD, cancel_rwts]
1726      \\ simp [GSYM realaxTheory.REAL_LDISTRIB]
1727      \\ match_mp_tac le2
1728      \\ conj_tac
1729      >- (match_mp_tac ge1_pow
1730          \\ Cases_on `a.Exponent`
1731          \\ lfs [])
1732      \\ match_mp_tac
1733           (REAL_ARITH ``2r <= a /\ 0 <= x ==> 2 <= a + x``)
1734      \\ simp [ge1_pow, DECIDE ``0n < a ==> 1 <= a``]
1735   ]
1736QED
1737
1738(* ------------------------------------------------------------------------ *)
1739
1740Theorem Float_is_finite[local]:
1741    !y: ('t, 'w) float r.
1742      (float_value y = Float r) ==> float_is_finite y
1743Proof
1744   rw [float_is_finite_def]
1745QED
1746
1747Theorem Float_float_to_real[local]:
1748    !y: ('t, 'w) float r.
1749      (float_value y = Float r) ==> (float_to_real y = r)
1750Proof rw [float_value_def]
1751QED
1752
1753Theorem float_is_zero_to_real:
1754   !x. float_is_zero x = (float_to_real x = 0)
1755Proof
1756   rw [float_is_zero_def, float_value_def, float_to_real_EQ0] >>
1757   simp[theorem "float_component_equality"]
1758QED
1759
1760(* !x. float_is_zero x ==> (float_to_real x = 0) *)
1761val float_is_zero_to_real_imp = iffLR float_is_zero_to_real
1762
1763Theorem pos_subnormal[local]:
1764    !a b n. 0 <= 2 / 2 pow a * (&n / 2 pow b)
1765Proof rrw [REAL_LE_MUL]
1766QED
1767
1768Theorem pos_normal[local]:
1769    !a b c n. 0 <= 2 pow a / 2 pow b * (1 + &n / 2 pow c)
1770Proof
1771   rw [REAL_LE_DIV, REAL_LE_MUL,
1772       REAL_ARITH ``0r <= n ==> 0 <= 1 + n``]
1773QED
1774
1775Theorem pos_normal2[local]:
1776    !a b c n. 0 <= 2 pow a / 2 pow b * (&n / 2 pow c)
1777Proof
1778   rw [REAL_LE_DIV, REAL_LE_MUL,
1779       REAL_ARITH ``0r <= n ==> 0 <= 1 + n``]
1780QED
1781
1782val thms =
1783   List.map REAL_ARITH
1784      [``a <= b /\ 0 <= c ==> a <= b + c: real``,
1785       ``0 <= b /\ a <= c ==> a <= b + c: real``,
1786       ``0 <= b /\ a <= c /\ 0 <= d ==> a <= b + (c + d): real``,
1787       ``a <= b /\ 0 <= c /\ 0 <= d ==> a <= b + c + d: real``,
1788       ``a <= b /\ 0 <= c /\ 0 <= d /\ 0 <= e ==> a <= b + c + (d + e): real``]
1789
1790Theorem diff_sign_ULP[local]:
1791    !x: ('t, 'w) float y: ('t, 'w) float.
1792        ~(float_is_zero x /\ float_is_zero y) /\ x.Sign <> y.Sign ==>
1793        ULP (x.Exponent,(:'t)) <= abs (float_to_real x - float_to_real y)
1794Proof
1795   NTAC 2 strip_tac
1796   \\ wordsLib.Cases_on_word_value `x.Sign`
1797   \\ wordsLib.Cases_on_word_value `y.Sign`
1798   \\ rw [ULP_def, float_to_real_def, float_is_zero, ABS_NEG,
1799          pos_normal, pos_subnormal,
1800          REAL_ARITH ``a - -1r * b * c = a + b * c``,
1801          REAL_ARITH ``-1r * b * c - a = -(b * c + a)``,
1802          REAL_ARITH ``0 <= a /\ 0 <= b ==> (abs (a + b) = a + b)``]
1803   \\ rw [REAL_LDISTRIB]
1804   >| (List.map match_mp_tac (thms @ thms))
1805   \\ rrw [pos_subnormal, pos_normal2, word_lt0, POW_ADD,
1806          REAL_LE_LDIV_EQ, le2, pow_ge2, ge1_pow, le_mult,
1807          fcpTheory.DIMINDEX_GE_1, REAL_LE_DIV, POW_2_LE1,
1808          cancel_rwts, DECIDE ``n <> 0n ==> 1 <= 2 * n``,
1809          REAL_ARITH ``2 <= a ==> 1r <= a``]
1810QED
1811
1812Theorem diff_sign_ULP_gt[local]:
1813   !x: ('t, 'w) float y: ('t, 'w) float.
1814     ~float_is_zero x /\ ~float_is_zero y /\ x.Sign <> y.Sign ==>
1815     ULP (x.Exponent,(:'t)) < abs (float_to_real x - float_to_real y)
1816Proof
1817   NTAC 2 strip_tac
1818   \\ wordsLib.Cases_on_word_value `x.Sign`
1819   \\ wordsLib.Cases_on_word_value `y.Sign`
1820   \\ rw [ULP_def, float_to_real_def, float_is_zero, ABS_NEG,
1821          pos_normal, pos_subnormal,
1822          REAL_ARITH ``a - -1r * b * c = a + b * c``,
1823          REAL_ARITH ``-1r * b * c - a = -(b * c + a)``,
1824          REAL_ARITH ``0 <= a /\ 0 <= b ==> (abs (a + b) = a + b)``]
1825   \\ rrw [REAL_LDISTRIB, REAL_RDISTRIB,
1826           REAL_LT_LDIV_EQ, REAL_LE_MUL,
1827           POW_2_LE1, POW_ADD,
1828           pos_subnormal, pos_normal2, word_lt0, cancel_rwts, word_lt0,
1829           prod_ge2, pow_ge2, le_mult,
1830           DECIDE ``0n < a /\ 0 < b ==> 2 < 2 * a + 2 * b``,
1831           DECIDE ``1n <= n <=> 0 < n``,
1832           DECIDE ``0n < x ==> 0 < 2 * x``,
1833           REAL_ARITH ``a <= b /\ 0r <= c /\ 1 <= d ==> a < b + c + d``,
1834           REAL_ARITH
1835              ``a <= b /\ 0r <= c /\ 2 <= d /\ 0 <= e ==> a < b + c + (d + e)``,
1836           REAL_ARITH
1837              ``1 <= a /\ 2r <= b /\ 0 <= c ==> 2 < 2 * a + (b + c)``
1838           |> Q.INST [`a` |-> `&n`]
1839           |> SIMP_RULE (srw_ss()) []
1840           ]
1841QED
1842
1843(* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -  *)
1844
1845(* |- !w. w2n w < 2 ** precision (:'a) *)
1846val w2n_lt_pow = REWRITE_RULE [wordsTheory.dimword_def] wordsTheory.w2n_lt
1847
1848Theorem w2n_lt_pow_sub1[local]:
1849  !x:'a word. x <> -1w ==> w2n x < 2 ** dimindex(:'a) - 1
1850Proof
1851   REPEAT strip_tac
1852   \\ match_mp_tac (DECIDE ``a < b /\ a <> b - 1 ==> a < b - 1n``)
1853   \\ simp [w2n_lt_pow]
1854   \\ Cases_on `x`
1855   \\ fs [wordsTheory.WORD_NEG_1, wordsTheory.word_T_def, wordsTheory.w2n_n2w,
1856          wordsTheory.UINT_MAX_def, wordsTheory.dimword_def]
1857QED
1858
1859Theorem nobias_denormal_lt_1[local]:
1860  !w:'t word. &w2n w / 2 pow precision (:'t) < 1
1861Proof
1862  rw [REAL_LT_LDIV_EQ, REAL_OF_NUM_POW, w2n_lt_pow]
1863QED
1864
1865Theorem nobias_denormal_lt_2[local]:
1866  !w:'t word. 2 * (&w2n w / 2 pow precision (:'t)) < 2
1867Proof
1868  rw [REAL_ARITH ``2r * n < 2 <=> n < 1``, nobias_denormal_lt_1]
1869QED
1870
1871Theorem subnormal_lt_normal[local]:
1872  !x y z.
1873    y <> 0w ==>
1874    2 / 2 pow bias (:'w) * (&w2n (x:'t word) / 2 pow precision (:'t)) <
1875    2 pow w2n (y:'w word) / 2 pow bias (:'w) *
1876    (1 + &w2n (z:'t word) / 2 pow precision (:'t))
1877Proof
1878   REPEAT strip_tac
1879   \\ once_rewrite_tac
1880        [REAL_LT_LMUL
1881         |> Q.SPEC `2 pow bias (:'w)`
1882         |> REWRITE_RULE [zero_lt_twopow]
1883         |> GSYM]
1884   \\ rewrite_tac [cancel_rwts, REAL_LDISTRIB]
1885   \\ match_mp_tac
1886        (REAL_ARITH ``a < 2r /\ 2 <= b /\ 0 <= c ==> a < b + c``)
1887   \\ rw [nobias_denormal_lt_2, pow_ge2, REAL_LE_MUL]
1888QED
1889
1890fun tac thm =
1891   REPEAT strip_tac
1892   \\ match_mp_tac (Q.SPECL [`a`, `b - c`, `x * b - c`] thm)
1893   \\ rsimp [REAL_LE_SUB_CANCEL2, GSYM REAL_LE_LDIV_EQ,
1894             REAL_DIV_REFL, REAL_SUB_ADD]
1895
1896Theorem weaken_leq[local]:
1897  !x a b c. 1r <= x /\ a <= b - c /\ 0 < b ==> a <= x * b - c
1898Proof tac REAL_LE_TRANS
1899QED
1900
1901Theorem weaken_lt[local]:
1902  !x a b c. 1r <= x /\ a < b - c /\ 0 < b ==> a < x * b - c
1903Proof tac REAL_LTE_TRANS
1904QED
1905
1906Theorem abs_diff1[local]:
1907    !s:word1 a b.
1908       a < b ==> (abs (-1 pow w2n s * a - -1 pow w2n s * b) = (b - a))
1909Proof
1910   wordsLib.Cases_word_value \\ rrw []
1911QED
1912
1913Theorem abs_diff2[local]:
1914    !s:word1 a b.
1915       b < a ==> (abs (-1 pow w2n s * a - -1 pow w2n s * b) = (a - b))
1916Proof
1917   wordsLib.Cases_word_value \\ rrw []
1918QED
1919
1920Theorem abs_diff1a[local] =
1921   abs_diff1
1922   |> Q.SPECL [`(y:('t,'w) float).Sign`,
1923               `(2 / 2 pow bias (:'w)) *
1924                (&w2n (x:('t,'w) float).Significand / 2 pow precision (:'t))`,
1925               `(2 pow w2n (y:('t,'w) float).Exponent / 2 pow bias (:'w)) *
1926                (1 + &w2n y.Significand / 2 pow precision (:'t))`]
1927
1928Theorem abs_diff1b[local] =
1929   abs_diff1
1930   |> Q.SPECL [`(y:('t,'w) float).Sign`,
1931               `(2 pow w2n (x:('t,'w) float).Exponent / 2 pow bias (:'w)) *
1932                (1 + &w2n x.Significand / 2 pow precision (:'t))`,
1933               `(2 pow (p + (w2n (x:('t,'w) float).Exponent + 1)) /
1934                2 pow bias (:'w)) *
1935                (1 + &w2n (y:('t,'w) float).Significand /
1936                2 pow precision (:'t))`]
1937
1938Theorem abs_diff1c[local] =
1939   abs_diff1
1940   |> Q.SPECL [`(y:('t,'w) float).Sign`,
1941               `(2 / 2 pow bias (:'w)) *
1942                (&w2n (x:('t,'w) float).Significand / 2 pow precision (:'t))`,
1943               `(2 / 2 pow bias (:'w)) *
1944                (&w2n (y:('t,'w) float).Significand / 2 pow precision (:'t))`]
1945
1946Theorem abs_diff1d[local] =
1947   abs_diff1
1948   |> Q.SPECL [`(y:('t,'w) float).Sign`,
1949               `(2 pow w2n (y:('t,'w) float).Exponent / 2 pow bias (:'w)) *
1950                (1 + &w2n (x:('t,'w) float).Significand /
1951                2 pow precision (:'t))`,
1952               `(2 pow w2n (y:('t,'w) float).Exponent / 2 pow bias (:'w)) *
1953                (1 + &w2n y.Significand / 2 pow precision (:'t))`]
1954
1955Theorem abs_diff1e[local] =
1956   abs_diff1
1957   |> Q.SPECL [`(y:('t,'w) float).Sign`,
1958               `(2 / 2 pow bias (:'w))`,
1959               `(2 pow w2n (y:('t,'w) float).Exponent / 2 pow bias (:'w)) *
1960                (1 + &w2n y.Significand / 2 pow precision (:'t))`]
1961
1962Theorem abs_diff1f[local] =
1963   abs_diff1
1964   |> Q.SPECL [`(y:('t,'w) float).Sign`,
1965               `(2 pow w2n (x:('t,'w) float).Exponent / 2 pow bias (:'w))`,
1966               `(2 pow w2n (y:('t,'w) float).Exponent / 2 pow bias (:'w)) *
1967                (1 + &w2n y.Significand / 2 pow precision (:'t))`]
1968
1969val abs_diff2a =
1970   abs_diff2
1971   |> Q.SPECL [`(y:('t,'w) float).Sign`,
1972               `(2 pow w2n (x:('t,'w) float).Exponent / 2 pow bias (:'w)) *
1973                (1 + &w2n x.Significand / 2 pow precision (:'t))`,
1974               `(2 / 2 pow bias (:'w)) *
1975                (&w2n (y:('t,'w) float).Significand / 2 pow precision (:'t))`]
1976
1977val abs_diff2b =
1978   abs_diff2
1979   |> Q.SPECL [`(y:('t,'w) float).Sign`,
1980               `(2 pow (p + (w2n (y:('t,'w) float).Exponent + 1)) /
1981                2 pow bias (:'w)) *
1982                (1 + &w2n (x:('t,'w) float).Significand /
1983                2 pow precision (:'t))`,
1984               `(2 pow w2n (y:('t,'w) float).Exponent / 2 pow bias (:'w)) *
1985                (1 + &w2n y.Significand / 2 pow precision (:'t))`]
1986
1987val abs_diff2c =
1988   abs_diff2
1989   |> Q.SPECL [`(y:('t,'w) float).Sign`,
1990               `(2 / 2 pow bias (:'w)) *
1991                (&w2n (x:('t,'w) float).Significand / 2 pow precision (:'t))`,
1992               `(2 / 2 pow bias (:'w)) *
1993                (&w2n (y:('t,'w) float).Significand / 2 pow precision (:'t))`]
1994
1995val abs_diff2d =
1996   abs_diff2
1997   |> Q.SPECL [`(y:('t,'w) float).Sign`,
1998               `(2 pow w2n (y:('t,'w) float).Exponent / 2 pow bias (:'w)) *
1999                (1 + &w2n (x:('t,'w) float).Significand /
2000                2 pow precision (:'t))`,
2001               `(2 pow w2n (y:('t,'w) float).Exponent / 2 pow bias (:'w)) *
2002                (1 + &w2n y.Significand / 2 pow precision (:'t))`]
2003
2004val abs_diff2e =
2005   abs_diff2
2006   |> Q.SPECL [`(y:('t,'w) float).Sign`,
2007               `(2 pow (w2n (y:('t,'w) float).Exponent + 1) /
2008                 2 pow bias (:'w))`,
2009               `(2 pow w2n (y:('t,'w) float).Exponent / 2 pow bias (:'w)) *
2010                (1 + &w2n (-1w: 't word) / 2 pow precision (:'t))`]
2011
2012fun abs_diff_tac thm =
2013   SUBGOAL_THEN
2014      (boolSyntax.rand (Thm.concl thm))
2015      (fn th => rewrite_tac [REWRITE_RULE [REAL_MUL_ASSOC] th])
2016
2017Theorem diff_exponent_boundary[local]:
2018  !x: ('t, 'w) float y: ('t, 'w) float.
2019    exponent_boundary y x ==>
2020    (abs (float_to_real x - float_to_real y) = ULP (y.Exponent, (:'t)))
2021Proof
2022   rw [ULP_def, exponent_boundary_def, float_to_real_def]
2023   >- (Cases_on `x.Exponent` \\ fs [])
2024   \\ abs_diff_tac abs_diff2e
2025   >- (match_mp_tac abs_diff2e
2026       \\ simp [REAL_LT_RDIV_EQ, REAL_LT_LMUL,
2027                REAL_ARITH ``2 * a = a * 2r``,
2028                REAL_ARITH ``1 + n < 2 <=> n < 1r``, cancel_rwts,
2029                REWRITE_RULE [arithmeticTheory.ADD1] pow]
2030       \\ simp [REAL_LT_LDIV_EQ, REAL_OF_NUM_POW,
2031                w2n_lt_pow]
2032      )
2033   \\ simp [REAL_EQ_RDIV_EQ, POW_ADD,
2034            REAL_SUB_RDISTRIB, cancel_rwts]
2035   \\ simp [REAL_ARITH
2036              ``(a * b * c - a * d * e) = a * (b * c - d * e: real)``,
2037            (GSYM o ONCE_REWRITE_RULE [REAL_MUL_COMM])
2038               REAL_EQ_RDIV_EQ,
2039            REAL_DIV_REFL]
2040   \\ simp [REAL_DIV_RMUL, REAL_EQ_SUB_RADD,
2041            REAL_ADD_RDISTRIB, wordsTheory.w2n_minus1,
2042            wordsTheory.dimword_def]
2043   \\ rsimp [REAL_OF_NUM_POW,
2044             DECIDE ``0 < n ==> (1 + (n + (n - 1)) = 2n * n)``]
2045QED
2046
2047val not_next_tac =
2048   REPEAT strip_tac
2049   \\ imp_res_tac arithmeticTheory.LESS_EQUAL_ADD
2050   \\ imp_res_tac arithmeticTheory.LESS_ADD_1
2051   \\ pop_assum kall_tac
2052   \\ REV_FULL_SIMP_TAC (srw_ss())
2053        [DECIDE ``1 < b ==> ((b - 1 = e) = (b = e + 1n))``]
2054   \\ simp [arithmeticTheory.LEFT_ADD_DISTRIB]
2055
2056local
2057   val lem1 = Q.prove(
2058      `!a b c. a < b /\ 1n < c ==> 2 * a + 2 <= b * c`,
2059      REPEAT strip_tac
2060      \\ imp_res_tac arithmeticTheory.LESS_ADD_1
2061      \\ simp []
2062      )
2063
2064   val lem1b = Q.prove(
2065      `!a b c. a + 1 < b /\ 1n < c ==> 2 * a + 2 < b * c`,
2066      REPEAT strip_tac
2067      \\ imp_res_tac arithmeticTheory.LESS_ADD_1
2068      \\ simp []
2069      )
2070
2071   val lem2 = Q.prove(
2072      `!x. x <> 0w ==> 1n < 2 EXP w2n x`,
2073      Cases
2074      \\ rw []
2075      \\ `0 < n` by decide_tac
2076      \\ imp_res_tac arithmeticTheory.LESS_ADD_1
2077      \\ simp [GSYM arithmeticTheory.ADD1, arithmeticTheory.EXP,
2078               DECIDE ``0n < n ==> 1 < 2 * n``]
2079      )
2080
2081   val lem3 = Q.prove(
2082      `!a b c. 2n <= b /\ 2 <= c /\ a < b ==>  2 * a + c <= b * c`,
2083      not_next_tac
2084      )
2085
2086   val lem3b = Q.prove(
2087      `!a b c d. 0 < b /\ 2n <= d /\ 2 <= c /\ a < d ==>
2088                 2 * a + c < b * c + d * c`,
2089      not_next_tac
2090      )
2091
2092   val lem3c = Q.prove(
2093      `!a b c. 2n <= b /\ 2 <= c /\ a + 1 < b ==> 2 * a + c < b * c`,
2094      not_next_tac
2095      )
2096
2097   val lem4 = Q.prove(
2098      `!a b c d. 1n < b /\ (4 <= a /\ c < b \/
2099                            2 <= a /\ c < b - 1 \/
2100                            2 <= a /\ c < b /\ 0 < d) ==>
2101                 a + (b + c) <= a * (b + d)`,
2102      not_next_tac
2103      )
2104
2105   val lem4b = Q.prove(
2106      `!a b c d. 2n <= a /\ 1n < b /\ 0 < d /\ c < b ==>
2107                 a + (b + c) < a * (b + d)`,
2108      not_next_tac
2109      )
2110
2111   (* |- 1 < 2 ** precision (:'a) *)
2112   val lem5 = REWRITE_RULE [wordsTheory.dimword_def] wordsTheory.ONE_LT_dimword
2113
2114   val t1 =
2115      simp [REAL_LE_LDIV_EQ, REAL_LT_LDIV_EQ,
2116            POW_ADD, REAL_SUB_RDISTRIB,
2117            REAL_LE_SUB_LADD, REAL_LT_SUB_LADD,
2118            GSYM realaxTheory.REAL_LDISTRIB, cancel_rwt, cancel_rwts]
2119      \\ simp [REAL_LDISTRIB]
2120   val t2 =
2121      once_rewrite_tac
2122          [REAL_LT_LMUL
2123           |> Q.SPEC `2 pow bias (:'w)`
2124           |> SIMP_RULE (srw_ss()) []
2125           |> GSYM]
2126      \\ rewrite_tac [cancel_rwts]
2127      \\ simp [POW_ADD]
2128      \\ once_rewrite_tac
2129           [REAL_ARITH ``a * (b * c) * d = b * (a * c * d): real``]
2130      \\ simp [REAL_LT_LMUL, REAL_LDISTRIB]
2131      \\ match_mp_tac
2132           (REAL_ARITH
2133               ``a < 1r /\ 1 <= b /\ 0 <= c ==> 1 + a < b * 2 + c``)
2134      \\ simp [nobias_denormal_lt_1, REAL_LE_MUL,
2135               POW_2_LE1]
2136   val t3 =
2137      simp [REAL_LE_LDIV_EQ, REAL_LT_LDIV_EQ]
2138      \\ simp [POW_ADD, cancel_rwts,
2139               REAL_SUB_RDISTRIB, REAL_DIV_RMUL]
2140      \\ simp [REAL_DIV_RMUL,
2141               REAL_ARITH
2142                   ``a * (b + c) * d = a * (b * d + c * d): real``]
2143      \\ simp [REAL_LE_SUB_LADD, REAL_LT_SUB_LADD,
2144               REAL_LDISTRIB]
2145   val t4 =
2146      match_mp_tac (REAL_ARITH ``a <= b /\ 0r <= c ==> a <= b + c``)
2147      \\ simp [REAL_LE_MUL, GSYM REAL_LE_SUB_LADD]
2148      \\ once_rewrite_tac
2149           [REAL_ARITH
2150               ``p * (a * 2r) * x - (a * x + a * y) =
2151                 a * ((p * (2 * x)) - (x + y))``]
2152      \\ match_mp_tac le_mult
2153      \\ simp []
2154      \\ match_mp_tac weaken_leq
2155      \\ simp [POW_2_LE1, REAL_LT_MUL,
2156               REAL_ARITH ``2r * a - (a + z) = a - z``]
2157   fun tac thm q =
2158      abs_diff_tac thm
2159      >- (match_mp_tac thm \\ t2)
2160      \\ Q.ABBREV_TAC `z:real = &w2n ^(Parse.Term q)`
2161      \\ t3
2162in
2163   fun tac1 thm =
2164      abs_diff_tac thm
2165      >- (match_mp_tac thm \\ simp [subnormal_lt_normal])
2166      \\ t1
2167      \\ match_mp_tac (REAL_ARITH ``0r <= c /\ a <= b ==> a <= b + c``)
2168      \\ simp [REAL_LE_MUL, REAL_OF_NUM_POW,
2169               fcpTheory.DIMINDEX_GE_1, lem1, lem2, lem3, lem5,
2170               word_ge1, w2n_lt_pow]
2171   val tac2 =
2172      tac abs_diff1b `(x: ('t, 'w) float).Significand`
2173      \\ t4
2174      \\ `?q. 1 <= q /\ (2 pow precision (:'t) = z + q)`
2175      by (ASSUME_TAC (Q.ISPEC `(x: ('t, 'w) float).Significand` w2n_lt_pow)
2176           \\ pop_assum
2177                (strip_assume_tac o MATCH_MP arithmeticTheory.LESS_ADD_1)
2178           \\ qexists_tac `&(p' + 1n)`
2179           \\ simp [REAL_OF_NUM_POW, Abbr `z`])
2180      \\ rsimp []
2181   val tac3 =
2182      tac abs_diff2b `(y: ('t, 'w) float).Significand`
2183      \\ once_rewrite_tac
2184           [div_le
2185            |> Q.SPEC `2r pow w2n (y:('t, 'w) float).Exponent`
2186            |> SIMP_RULE (srw_ss()) []
2187            |> GSYM]
2188      \\ simp [GSYM REAL_DIV_ADD, GSYM REAL_ADD_LDISTRIB,
2189               cancel_rwts]
2190      \\ rsimp [REAL_OF_NUM_POW, Abbr `z`]
2191      \\ match_mp_tac lem4
2192      \\ full_simp_tac (srw_ss()) [exponent_boundary_def]
2193      \\ REV_FULL_SIMP_TAC (srw_ss())
2194           [w2n_lt_pow, w2n_lt_pow_sub1, word_lt0, ge4, lem5]
2195      \\ `p <> 0` by (strip_tac \\
2196                      full_simp_tac (srw_ss())
2197                                [DECIDE ``(1 = x + 1) = (x = 0n)``])
2198      \\ full_simp_tac (srw_ss()) [ge4]
2199   val tac4 =
2200      abs_diff_tac abs_diff1a
2201      >- (match_mp_tac abs_diff1a \\ simp [subnormal_lt_normal])
2202      \\ t1
2203      \\ match_mp_tac (REAL_ARITH ``0r <= c /\ a < b ==> a < b + c``)
2204      \\ simp [REAL_LE_MUL, REAL_OF_NUM_POW,
2205               fcpTheory.DIMINDEX_GE_1, not_max_suc_lt_dimword, lem1b, lem2]
2206   val tac5 =
2207      abs_diff_tac abs_diff2a
2208      >- (match_mp_tac abs_diff2a \\ simp [subnormal_lt_normal])
2209      \\ t1
2210      \\ simp [REAL_OF_NUM_POW]
2211      \\ match_mp_tac lem3b
2212      \\ simp [REAL_LE_MUL, REAL_OF_NUM_POW,
2213               fcpTheory.DIMINDEX_GE_1, word_lt0, not_max_suc_lt_dimword,
2214               lem1b, lem2, word_ge1, w2n_lt_pow]
2215   val tac6 =
2216      tac abs_diff1b `(x: ('t, 'w) float).Significand`
2217      \\ match_mp_tac (REAL_ARITH ``a < b /\ 0r <= c ==> a < b + c``)
2218      \\ simp [REAL_LE_MUL, GSYM REAL_LT_SUB_LADD]
2219      \\ once_rewrite_tac
2220           [REAL_ARITH
2221               ``p * (a * 2r) * x - (a * x + a * y) =
2222                 a * ((p * (2 * x)) - (x + y))``]
2223      \\ match_mp_tac (ONCE_REWRITE_RULE [REAL_MUL_COMM] gt_mult)
2224      \\ simp []
2225      \\ match_mp_tac weaken_lt
2226      \\ simp [POW_2_LE1, REAL_LT_MUL,
2227               REAL_ARITH ``2r * a - (a + z) = a - z``]
2228      \\ simp [GSYM REAL_LT_ADD_SUB, not_max_suc_lt_dimword,
2229               REAL_OF_NUM_POW, Abbr `z`]
2230   val tac7 =
2231      tac abs_diff2b `(y: ('t, 'w) float).Significand`
2232      \\ once_rewrite_tac
2233           [REAL_LT_RDIV
2234            |> Q.SPECL [`x`, `y`, `2 pow w2n (y:('t, 'w) float).Exponent`]
2235            |> SIMP_RULE (srw_ss()) []
2236            |> GSYM]
2237      \\ simp [GSYM REAL_DIV_ADD, GSYM REAL_ADD_LDISTRIB,
2238               cancel_rwts]
2239      \\ rsimp [REAL_OF_NUM_POW, Abbr `z`]
2240      \\ match_mp_tac lem4b
2241      \\ REV_FULL_SIMP_TAC (srw_ss()) [w2n_lt_pow, word_lt0, lem5]
2242end
2243
2244Theorem diff_exponent_ULP[local]:
2245  !x: ('t, 'w) float y: ('t, 'w) float.
2246    (x.Sign = y.Sign) /\ x.Exponent <> y.Exponent /\
2247    ~exponent_boundary y x ==>
2248    ULP (x.Exponent, (:'t)) <= abs (float_to_real x - float_to_real y)
2249Proof
2250   rw [ULP_def, float_to_real_def]
2251   >- tac1 abs_diff1a
2252   >- tac1 abs_diff2a
2253   \\ `w2n x.Exponent < w2n y.Exponent \/ w2n y.Exponent < w2n x.Exponent`
2254   by metis_tac [arithmeticTheory.LESS_LESS_CASES, wordsTheory.w2n_11]
2255   \\ imp_res_tac arithmeticTheory.LESS_ADD_1
2256   \\ simp []
2257   >- tac2
2258   \\ fs []
2259   \\ tac3
2260QED
2261
2262Theorem diff_exponent_ULP_gt[local]:
2263  !x: ('t, 'w) float y: ('t, 'w) float.
2264    (x.Sign = y.Sign) /\ x.Exponent <> y.Exponent /\
2265    x.Significand NOTIN {0w; -1w} ==>
2266    ULP (x.Exponent, (:'t)) < abs (float_to_real x - float_to_real y)
2267Proof
2268   rw [ULP_def, float_to_real_def]
2269   >- tac4
2270   >- tac5
2271   \\ `w2n x.Exponent < w2n y.Exponent \/ w2n y.Exponent < w2n x.Exponent`
2272   by metis_tac [arithmeticTheory.LESS_LESS_CASES, wordsTheory.w2n_11]
2273   \\ imp_res_tac arithmeticTheory.LESS_ADD_1
2274   \\ simp []
2275   >- tac6
2276   \\ fs []
2277   \\ tac7
2278QED
2279
2280Theorem lem[local]:
2281    !a b m. 2n <= a /\ 2 <= b /\ 1 <= m ==> a * b + b < 2 * (m * a * b)
2282Proof
2283   REPEAT strip_tac
2284   \\ imp_res_tac arithmeticTheory.LESS_EQUAL_ADD
2285   \\ simp [arithmeticTheory.LEFT_ADD_DISTRIB]
2286QED
2287
2288Theorem diff_exponent_ULP_gt0[local]:
2289  !x: ('t, 'w) float y: ('t, 'w) float.
2290    (x.Sign = y.Sign) /\ x.Exponent <+ y.Exponent /\
2291    (x.Significand = 0w) /\ ~float_is_zero x ==>
2292    ULP (x.Exponent, (:'t)) < abs (float_to_real x - float_to_real y)
2293Proof
2294   rw [ULP_def, float_to_real_def, ABS_NEG, abs_float_value,
2295       abs_significand, ABS_MUL, ABS_DIV,
2296       ABS_N, gt0_abs, wordsTheory.WORD_LO]
2297   >- rfs [REAL_LT_LDIV_EQ, POW_ADD,
2298           REAL_SUB_LDISTRIB, REAL_SUB_RDISTRIB,
2299           REAL_LT_SUB_LADD, cancel_rwts, cancel_rwt, float_is_zero]
2300   \\ abs_diff_tac abs_diff1f
2301   >- (match_mp_tac abs_diff1f
2302       \\ simp [REAL_LT_LDIV_EQ, REAL_ADD_LDISTRIB,
2303                cancel_rwts]
2304       \\ simp [REAL_OF_NUM_POW, REAL_LE_MUL,
2305                REAL_ARITH ``a < b /\ 0 <= c ==> a < b + c: real``])
2306   \\ simp [REAL_LT_LDIV_EQ, POW_ADD,
2307            REAL_SUB_LDISTRIB, REAL_SUB_RDISTRIB,
2308            REAL_LT_SUB_LADD, cancel_rwts, cancel_rwt]
2309   \\ match_mp_tac (REAL_ARITH ``a < b /\ 0r <= c ==> a < b + c``)
2310   \\ simp [REAL_LE_MUL, REAL_OF_NUM_POW]
2311   \\ imp_res_tac arithmeticTheory.LESS_ADD_1
2312   \\ simp [arithmeticTheory.EXP_ADD, lem, fcpTheory.DIMINDEX_GE_1, exp_ge4,
2313            word_ge1]
2314QED
2315
2316Theorem lem[local]:
2317    !a b. 2 <= a /\ 4n <= b ==> 2 * a + 2 < a * b
2318Proof
2319   not_next_tac
2320QED
2321
2322Theorem diff_exponent_ULP_gt01[local]:
2323  !x: ('t, 'w) float y: ('t, 'w) float.
2324    (x.Sign = y.Sign) /\ x.Exponent <> y.Exponent /\
2325    y.Significand <> -1w  /\ (x.Significand = 0w) /\ (x.Exponent = 1w) ==>
2326    ULP (x.Exponent, (:'t)) < abs (float_to_real x - float_to_real y)
2327Proof
2328   rw [ULP_def, float_to_real_def, ABS_NEG, abs_float_value,
2329       abs_significand, ABS_MUL, ABS_DIV,
2330       ABS_N, gt0_abs, nobias_denormal_lt_1,
2331       REAL_ARITH ``a - a * b = a * (1 - b): real``,
2332       REAL_ARITH ``a < 1 ==> (abs (1 - a) = 1 - a)``]
2333   >- (simp [REAL_LT_LDIV_EQ, POW_ADD, cancel_rwts,
2334             REAL_SUB_LDISTRIB, REAL_SUB_RDISTRIB,
2335             REAL_LT_SUB_LADD]
2336       \\ rewrite_tac [simpLib.SIMP_PROVE (srw_ss()++ARITH_ss) []
2337                         ``&(2n * n + 2) = 2r * &(n + 1)``]
2338       \\ simp [REAL_LT_LMUL, REAL_OF_NUM_POW,
2339                not_max_suc_lt_dimword])
2340   \\ `1w <+ y.Exponent`
2341   by metis_tac
2342        [wordsLib.WORD_DECIDE ``a:'a word <> 0w /\ a <> 1w ==> 1w <+ a``]
2343   \\ fs [wordsTheory.WORD_LO]
2344   \\ abs_diff_tac abs_diff1e
2345   >- (match_mp_tac abs_diff1e
2346       \\ simp [REAL_LT_LDIV_EQ, REAL_LDISTRIB,
2347                cancel_rwts, cancel_rwt]
2348       \\ match_mp_tac (REAL_ARITH ``a < b /\ 0r <= c ==> a < b + c``)
2349       \\ simp [REAL_LE_MUL, REAL_OF_NUM_POW])
2350   \\ simp [REAL_LT_LDIV_EQ, POW_ADD,
2351            REAL_SUB_LDISTRIB, REAL_SUB_RDISTRIB,
2352            REAL_LT_SUB_LADD, cancel_rwts, cancel_rwt]
2353   \\ match_mp_tac (REAL_ARITH ``a < b /\ 0r <= c ==> a < b + c``)
2354   \\ simp [REAL_LE_MUL, REAL_OF_NUM_POW, lem,
2355            fcpTheory.DIMINDEX_GE_1, exp_ge4]
2356QED
2357
2358Theorem lem[local]:
2359    !a b c. a < b /\ 2n <= c ==> 2 * a < b * c
2360Proof
2361   not_next_tac
2362QED
2363
2364Theorem lem2[local]:
2365    !a b c. a < b /\ 1n <= c ==> a + b < 2 * (c * b)
2366Proof
2367   not_next_tac
2368QED
2369
2370Theorem float_to_real_lt_exponent_mono[local]:
2371  !x: ('t, 'w) float y: ('t, 'w) float.
2372    (x.Sign = y.Sign) /\ abs (float_to_real x) <= abs (float_to_real y) ==>
2373    x.Exponent <=+ y.Exponent
2374Proof
2375  rw [float_to_real_def, ABS_NEG, abs_float_value,
2376      abs_significand, ABS_MUL, ABS_DIV,
2377      ABS_N, gt0_abs, wordsTheory.WORD_LS]
2378   >| [
2379      Cases_on `x.Sign = y.Sign`
2380      \\ simp [REAL_NOT_LE]
2381      \\ once_rewrite_tac
2382           [REAL_LT_RMUL
2383            |> Q.SPECL [`x`, `y`, `2 pow bias (:'w)`]
2384            |> REWRITE_RULE [zero_lt_twopow]
2385            |> GSYM]
2386      \\ rewrite_tac [cancel_rwts, cancel_rwt]
2387      \\ simp [REAL_LDISTRIB, REAL_RDISTRIB,
2388               REAL_OF_NUM_POW, REAL_DIV_RMUL,
2389               REAL_LT_LDIV_EQ, mult_ratr,
2390               cancel_rwts, cancel_rwt, w2n_lt_pow,
2391               word_ge1, lem, DECIDE ``a < b ==> a < x + b: num``],
2392      (* --- *)
2393      pop_assum mp_tac
2394      \\ Cases_on `w2n x.Exponent <= w2n y.Exponent`
2395      \\ simp [REAL_NOT_LE]
2396      \\ fs [arithmeticTheory.NOT_LESS_EQUAL]
2397      \\ once_rewrite_tac
2398           [REAL_LT_RMUL
2399            |> Q.SPECL [`x`, `y`, `2 pow bias (:'w)`]
2400            |> REWRITE_RULE [zero_lt_twopow]
2401            |> GSYM]
2402      \\ rewrite_tac [cancel_rwts, cancel_rwt]
2403      \\ once_rewrite_tac
2404           [REAL_LT_RMUL
2405            |> Q.SPECL [`x`, `y`, `2 pow precision (:'t)`]
2406            |> REWRITE_RULE [zero_lt_twopow]
2407            |> GSYM]
2408      \\ rewrite_tac [cancel_rwts, cancel_rwt]
2409      \\ simp [REAL_OF_NUM_POW]
2410      \\ match_mp_tac (DECIDE ``a < b ==> a < x + b: num``)
2411      \\ imp_res_tac arithmeticTheory.LESS_ADD_1
2412      \\ asm_simp_tac bool_ss
2413             [arithmeticTheory.EXP_ADD, arithmeticTheory.LT_MULT_RCANCEL,
2414              GSYM arithmeticTheory.RIGHT_ADD_DISTRIB,
2415              DECIDE ``a * (b * (c * d)) = (a * c * d) * b: num``]
2416      \\ simp [lem2, w2n_lt_pow]
2417   ]
2418QED
2419
2420(* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -  *)
2421
2422fun tac thm =
2423   abs_diff_tac thm
2424   >- (match_mp_tac thm
2425       \\ simp [REAL_LT_LMUL, REAL_LT_DIV,
2426                REAL_LT_LDIV_EQ, REAL_DIV_RMUL])
2427   \\ simp [REAL_ARITH ``a < b ==> (abs (a - b) = b - a)``,
2428            REAL_ARITH ``b < a ==> (abs (a - b) = a - b)``,
2429            REAL_SUB_RDISTRIB, REAL_LDISTRIB,
2430            POW_ADD, mult_rat]
2431   \\ simp [mult_ratr]
2432
2433fun tac2 thm =
2434   abs_diff_tac thm
2435   >- (match_mp_tac thm
2436       \\ simp [REAL_LT_LMUL, REAL_LT_DIV,
2437                REAL_LT_LDIV_EQ, REAL_DIV_RMUL])
2438   \\ simp [REAL_ARITH ``a < b ==> (abs (a - b) = b - a)``,
2439            REAL_ARITH ``b < a ==> (abs (a - b) = a - b)``,
2440            REAL_ARITH ``1 + a - (1 + b) = a - b: real``,
2441            GSYM REAL_SUB_LDISTRIB, sub_rat_same_base]
2442   \\ simp [POW_ADD, mult_rat]
2443   \\ simp_tac (srw_ss()++realSimps.real_ac_SS) [mult_ratr]
2444
2445Theorem diff_significand_ULP_mul[local]:
2446    !x: ('t, 'w) float y: ('t, 'w) float.
2447        (x.Sign = y.Sign) /\ (x.Exponent = y.Exponent) ==>
2448        (abs (float_to_real x - float_to_real y) =
2449         abs (&w2n x.Significand - &w2n y.Significand) *
2450         ULP (x.Exponent, (:'t)))
2451Proof
2452   rw [ULP_def, float_to_real_def]
2453   \\ (Cases_on `x.Significand = y.Significand`
2454       >- rsimp [])
2455   \\ `w2n x.Significand < w2n y.Significand \/
2456       w2n y.Significand < w2n x.Significand`
2457      by metis_tac [arithmeticTheory.LESS_LESS_CASES, wordsTheory.w2n_11]
2458   >- tac abs_diff1c
2459   >- tac abs_diff2c
2460   >- tac2 abs_diff1d
2461   \\ tac2 abs_diff2d
2462QED
2463
2464Theorem diff_ge1[local]:
2465    !a b. 1 <= abs (&a - &b) <=> &a <> (&b: real)
2466Proof
2467   lrw [REAL_SUB, ABS_NEG, ABS_N]
2468QED
2469
2470Theorem diff_significand_ULP[local]:
2471    !x: ('t, 'w) float y: ('t, 'w) float.
2472        (x.Sign = y.Sign) /\ (x.Exponent = y.Exponent) /\
2473        x.Significand <> y.Significand ==>
2474        ULP (x.Exponent, (:'t)) <= abs (float_to_real x - float_to_real y)
2475Proof
2476   rw [diff_significand_ULP_mul, ULP_gt0, diff_ge1,
2477       ONCE_REWRITE_RULE [REAL_MUL_COMM] le_mult]
2478QED
2479
2480(* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -  *)
2481
2482Theorem ULP_same[local]:
2483    !x y.
2484      (x = y) ==>
2485      ~(ULP (x.Exponent, (:'t)) <= abs (float_to_real x - float_to_real y))
2486Proof
2487   rrw [ULP_gt0, REAL_NOT_LE]
2488QED
2489
2490Theorem diff_sign_neq[local]:
2491    !x: ('t, 'w) float y: ('t, 'w) float.
2492        ~(float_is_zero x /\ float_is_zero y) /\ x.Sign <> y.Sign ==>
2493        float_to_real x <> float_to_real y
2494Proof
2495   metis_tac [diff_sign_ULP, ULP_same]
2496QED
2497
2498Theorem diff_exponent_neq[local]:
2499    !x: ('t, 'w) float y: ('t, 'w) float.
2500        (x.Sign = y.Sign) /\ x.Exponent <> y.Exponent ==>
2501        float_to_real x <> float_to_real y
2502Proof
2503   REPEAT strip_tac
2504   \\ Cases_on `exponent_boundary y x`
2505   >- (fs []
2506       \\ imp_res_tac diff_exponent_boundary
2507       \\ rfs [ULP_gt0, REAL_POS_NZ])
2508   \\metis_tac [diff_exponent_ULP, ULP_same]
2509QED
2510
2511Theorem float_to_real_eq:
2512    !x: ('t, 'w) float y: ('t, 'w) float.
2513       (float_to_real x = float_to_real y) <=>
2514       (x = y) \/ (float_is_zero x /\ float_is_zero y)
2515Proof
2516   NTAC 2 strip_tac
2517   \\ Cases_on `x = y`
2518   \\ simp []
2519   \\ Cases_on `float_is_zero x /\ float_is_zero y`
2520   \\ simp [float_is_zero_to_real_imp]
2521   \\ Cases_on `x.Sign <> y.Sign`
2522   >- metis_tac [diff_sign_neq]
2523   \\ Cases_on `x.Exponent <> y.Exponent`
2524   >- metis_tac [diff_exponent_neq]
2525   \\ qpat_x_assum `~(p /\ q)` kall_tac
2526   \\ fs [float_component_equality]
2527   \\ rw [float_to_real_def, sign_not_zero, div_twopow]
2528QED
2529
2530Theorem diff_float_ULP:
2531    !x: ('t, 'w) float y: ('t, 'w) float.
2532       float_to_real x <> float_to_real y /\ ~exponent_boundary y x ==>
2533       ULP (x.Exponent, (:'t)) <= abs (float_to_real x - float_to_real y)
2534Proof
2535   rw [float_to_real_eq, float_component_equality]
2536   \\ metis_tac [diff_sign_ULP, diff_exponent_ULP, diff_significand_ULP]
2537QED
2538
2539(* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -  *)
2540
2541(* |- !x y. ~float_is_zero y ==>
2542            ((float_to_real x = float_to_real y) <=> (x = y)) *)
2543Theorem float_to_real_11_right[local] =
2544   float_to_real_eq
2545   |> Drule.SPEC_ALL
2546   |> Q.DISCH `~float_is_zero y`
2547   |> SIMP_RULE bool_ss []
2548   |> Q.GENL [`x`, `y`]
2549
2550(* |- !x y. ~float_is_zero x ==>
2551            ((float_to_real x = float_to_real y) <=> (x = y))
2552val float_to_real_11_left =
2553   float_to_real_eq
2554   |> Drule.SPEC_ALL
2555   |> Q.DISCH `~float_is_zero x`
2556   |> SIMP_RULE bool_ss []
2557   |> Drule.GEN_ALL
2558*)
2559
2560(* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -  *)
2561
2562Theorem diff1pos[local]:
2563  !a. a <> 0w ==> (&w2n a - &w2n (a + -1w) = 1r)
2564Proof
2565   Cases
2566   \\ Cases_on `n`
2567   \\ simp [wordsTheory.n2w_SUC]
2568   \\ rrw [REAL_SUB, bitTheory.SUC_SUB, DECIDE ``~(SUC n <= n)``]
2569QED
2570
2571Theorem diff1neg[local]:
2572  !a. a <> -1w ==> (&w2n a - &w2n (a + 1w) = -1r)
2573Proof
2574   rw [REAL_SUB, bitTheory.SUC_SUB, DECIDE ``~(SUC n <= n)``,
2575       GSYM wordsTheory.WORD_LS,
2576       ONCE_REWRITE_RULE [GSYM wordsTheory.WORD_ADD_COMM]
2577          wordsTheory.WORD_ADD_RIGHT_LS2]
2578   \\ lfs [wordsTheory.WORD_NOT_LOWER, wordsTheory.WORD_LS_word_T]
2579   \\ `a <=+ a + 1w` by wordsLib.WORD_DECIDE_TAC
2580   \\ simp [GSYM wordsTheory.word_sub_w2n]
2581QED
2582
2583Theorem must_be_1[local]:
2584  !a b: real. 0 < b ==> ((a * b = b) = (a = 1))
2585Proof
2586   REPEAT strip_tac
2587   \\ Cases_on `a = 1`
2588   >- simp []
2589   \\ Cases_on `a < 1`
2590   >- rsimp [REAL_LT_IMP_NE,
2591             ONCE_REWRITE_RULE [REAL_MUL_COMM] lt_mult]
2592   \\ `1 < a` by rsimp []
2593   \\ simp [REAL_LT_IMP_NE, gt_mult]
2594QED
2595
2596Theorem w2n_add1[local]:
2597  !a. a <> -1w ==> (w2n a + 1 = w2n (a + 1w))
2598Proof
2599   Cases
2600   \\ lrw [wordsTheory.word_eq_n2w, wordsTheory.word_add_n2w,
2601           bitTheory.MOD_2EXP_MAX_def, bitTheory.MOD_2EXP_def,
2602           GSYM wordsTheory.dimword_def]
2603QED
2604
2605Theorem diff_ulp_next_float[local]:
2606   !x y: ('t, 'w) float.
2607       ~float_is_zero x /\ y.Significand NOTIN {0w; -1w} ==>
2608       ((abs (float_to_real y - float_to_real x) = ULP (y.Exponent,(:'t))) <=>
2609        (x = y with Significand := y.Significand - 1w) \/
2610        (x = y with Significand := y.Significand + 1w))
2611Proof
2612   REPEAT strip_tac
2613   \\ eq_tac
2614   >| [
2615      `~float_is_zero y` by fs [float_is_zero]
2616      \\ Cases_on `x.Sign <> y.Sign`
2617      >- prove_tac [REAL_LT_IMP_NE, diff_sign_ULP_gt]
2618      \\ Cases_on `x.Exponent <> y.Exponent`
2619      >- prove_tac [REAL_LT_IMP_NE, diff_exponent_ULP_gt]
2620      \\ fs [diff_significand_ULP_mul, must_be_1, ULP_gt0,
2621             float_component_equality]
2622      \\ Cases_on `x.Significand = y.Significand + -1w`
2623      \\ simp []
2624      \\ Cases_on `x.Significand = y.Significand + 1w`
2625      \\ rsimp [REAL_ARITH ``(abs x = 1) <=> (x = 1) \/ (x = -1)``,
2626                REAL_ARITH ``(a = -1 + c) = (c = a + 1r)``,
2627                REAL_EQ_SUB_RADD, w2n_add1]
2628      \\ Cases_on `x.Significand = -1w`
2629      \\ simp [ONCE_REWRITE_RULE [arithmeticTheory.ADD_COMM] w2n_add1,
2630               wordsTheory.w2n_minus1, DECIDE ``0n < n ==> (1 + (n - 1) = n)``,
2631               wordsTheory.w2n_lt, prim_recTheory.LESS_NOT_EQ,
2632               wordsLib.WORD_ARITH_PROVE
2633                  ``a:'a word <> b + -1w ==> b <> a + 1w``],
2634      (* --- *)
2635      rw []
2636      \\ rw [float_to_real_def, abs_float_value, abs_significand,
2637             ABS_MUL, ABS_DIV, ABS_N,
2638             gt0_abs, GSYM REAL_SUB_LDISTRIB, sub_rat_same_base,
2639             REAL_ARITH ``1r + a - (1 + b) = a - b``]
2640      \\ fs [diff1pos, diff1neg, mult_rat, ULP_def,
2641             POW_ADD]
2642   ]
2643QED
2644
2645Theorem diff_ulp_next_float0[local]:
2646  !x y: ('t, 'w) float.
2647    ~float_is_zero x /\ ~float_is_zero y /\ (y.Significand = 0w) /\
2648    abs (float_to_real y) <= abs (float_to_real x) ==>
2649    ((abs (float_to_real y - float_to_real x) = ULP (y.Exponent,(:'t))) =
2650     (x = y with Significand := y.Significand + 1w))
2651Proof
2652   REPEAT strip_tac
2653   \\ eq_tac
2654   >| [
2655      Cases_on `x.Sign <> y.Sign`
2656      >- prove_tac [REAL_LT_IMP_NE, diff_sign_ULP_gt]
2657      \\ imp_res_tac float_to_real_lt_exponent_mono
2658      \\ Cases_on `x.Exponent <> y.Exponent`
2659      >- prove_tac
2660            [REAL_LT_IMP_NE, diff_exponent_ULP_gt0,
2661             wordsLib.WORD_DECIDE ``a <> b /\ a <=+ b ==> a <+ b:'a word``]
2662      \\ fs [diff_significand_ULP_mul, must_be_1, ULP_gt0,
2663             float_component_equality, ABS_NEG, ABS_N]
2664      \\ Cases_on `x.Significand`
2665      \\ simp [],
2666      (* --- *)
2667      rw []
2668      \\ rw [float_to_real_def, abs_float_value, abs_significand,
2669             ABS_MUL, ABS_DIV, ABS_N,
2670             ABS_NEG, gt0_abs, REAL_LDISTRIB,
2671             REAL_ARITH ``a - (a + c) = -c: real``]
2672      \\ fs [diff1pos, diff1neg, mult_rat, ULP_def,
2673             POW_ADD]
2674   ]
2675QED
2676
2677Theorem diff_ulp_next_float01[local]:
2678  !x y: ('t, 'w) float.
2679    ~float_is_zero x /\ ~float_is_zero y /\
2680    x.Significand <> -1w /\ (y.Significand = 0w) /\ (y.Exponent = 1w) ==>
2681    ((abs (float_to_real y - float_to_real x) = ULP (y.Exponent,(:'t))) =
2682     (x = y with Significand := y.Significand + 1w))
2683Proof
2684   REPEAT strip_tac
2685   \\ eq_tac
2686   >| [
2687      Cases_on `x.Sign <> y.Sign`
2688      >- prove_tac [REAL_LT_IMP_NE, diff_sign_ULP_gt]
2689      \\ Cases_on `x.Exponent <> y.Exponent`
2690      >- prove_tac [REAL_LT_IMP_NE, diff_exponent_ULP_gt01]
2691      \\ fs [diff_significand_ULP_mul, must_be_1, ULP_gt0,
2692             float_component_equality, ABS_NEG, ABS_N]
2693      \\ Cases_on `x.Significand`
2694      \\ simp [],
2695      (* --- *)
2696      rw []
2697      \\ rw [float_to_real_def, abs_float_value, abs_significand,
2698             ABS_MUL, ABS_DIV, ABS_N,
2699             ABS_NEG, gt0_abs, REAL_LDISTRIB,
2700             REAL_ARITH ``a - (a + c) = -c: real``]
2701      \\ fs [diff1pos, diff1neg, mult_rat, ULP_def,
2702             POW_ADD]
2703   ]
2704QED
2705
2706Theorem float_min_equiv_ULP_eq_float_to_real[local]:
2707  !y: ('t, 'w) float.
2708    (abs (float_to_real y) = ULP (y.Exponent,(:'t))) <=>
2709    y IN {float_plus_min (:'t # 'w); float_minus_min (:'t # 'w)}
2710Proof
2711   strip_tac
2712   \\ Cases_on `float_is_zero y`
2713   >- fs [float_sets, zero_to_real, float_components, float_distinct,
2714          GSYM float_distinct, ULP_gt0,
2715          REAL_ARITH ``0 < b ==> 0r <> b``]
2716   \\ Cases_on `(y = float_plus_min (:'t # 'w)) \/
2717                (y = float_minus_min (:'t # 'w))`
2718   >- rw [GSYM neg_ulp, GSYM ulp, float_minus_min_def, float_components,
2719          ulp_def, ULP_gt0, gt0_abs, REAL_LT_IMP_LE,
2720          ABS_NEG]
2721   \\ fs []
2722   \\ rw [float_to_real_def, ULP_def, abs_float_value, abs_significand,
2723          ABS_MUL, ABS_DIV, ABS_N, gt0_abs,
2724          REAL_EQ_RDIV_EQ]
2725   \\ simp [POW_ADD, GSYM REAL_LDISTRIB,
2726            cancel_rwts, cancel_rwt, REAL_DIV_REFL,
2727            REAL_EQ_RDIV_EQ
2728            |> ONCE_REWRITE_RULE [GSYM REAL_MUL_COMM]
2729            |> GSYM]
2730   >| [
2731      strip_tac
2732      \\ `y.Significand = 1w`
2733      by metis_tac [wordsTheory.w2n_11, wordsTheory.word_1_n2w]
2734      \\ fs [float_plus_min_def, float_minus_min_def, float_negate_def,
2735             float_component_equality]
2736      \\ metis_tac [sign_inconsistent],
2737      simp [REAL_OF_NUM_POW, GSYM wordsTheory.dimword_def,
2738            DECIDE ``1 < a ==> a + b <> 1n``]
2739   ]
2740QED
2741
2742(* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -  *)
2743
2744val tac =
2745   REPEAT strip_tac
2746   \\ spose_not_then assume_tac
2747   \\ `float_to_real a <> float_to_real b`
2748   by metis_tac [float_to_real_eq]
2749   \\ imp_res_tac diff_float_ULP
2750   \\ rlfs []
2751
2752Theorem diff_lt_ulp_eq0:
2753    !a: ('t, 'w) float b: ('t, 'w) float x.
2754       ~exponent_boundary b a /\
2755       abs (x - float_to_real a) < ULP (a.Exponent, (:'t)) /\
2756       abs (x - float_to_real b) < ULP (a.Exponent, (:'t)) /\
2757       abs (float_to_real a) <= abs x /\ abs (float_to_real b) <= abs x /\
2758       ~float_is_zero a ==>
2759       (b = a)
2760Proof tac
2761QED
2762
2763Theorem diff_lt_ulp_even:
2764    !a: ('t, 'w) float b: ('t, 'w) float x.
2765       ~exponent_boundary b a /\
2766       2 * abs (float_to_real a - x) < ULP (a.Exponent, (:'t)) /\
2767       2 * abs (float_to_real b - x) < ULP (a.Exponent, (:'t)) /\
2768       ~float_is_zero a ==>
2769       (b = a)
2770Proof
2771   REPEAT strip_tac
2772   \\ spose_not_then assume_tac
2773   \\ `float_to_real a <> float_to_real b`
2774   by metis_tac [float_to_real_eq]
2775   \\ imp_res_tac diff_float_ULP
2776   \\ rlfs []
2777QED
2778
2779Theorem diff_lt_ulp_even4:
2780    !a: ('t, 'w) float b: ('t, 'w) float x.
2781       ~exponent_boundary b a /\
2782       4 * abs (float_to_real a - x) <= ULP (a.Exponent, (:'t)) /\
2783       4 * abs (float_to_real b - x) <= ULP (a.Exponent, (:'t)) /\
2784       ~float_is_zero a ==>
2785       (b = a)
2786Proof
2787   REPEAT strip_tac
2788   \\ spose_not_then assume_tac
2789   \\ `float_to_real a <> float_to_real b`
2790   by metis_tac [float_to_real_eq]
2791   \\ imp_res_tac diff_float_ULP
2792   \\ rlfs []
2793QED
2794
2795(*
2796val diff_lt_ulp_eq_pos = Q.store_thm("diff_lt_ulp_eq_pos",
2797   `!a: ('t, 'w) float b: ('t, 'w) float x.
2798       ~exponent_boundary b a /\
2799       abs (x - float_to_real a) < ULP (a.Exponent, (:'t)) /\
2800       abs (x - float_to_real b) < ULP (a.Exponent, (:'t)) /\
2801       float_to_real a >= x /\ float_to_real b >= x /\
2802       ~float_is_zero b ==>
2803       (a = b)`,
2804   tac)
2805
2806val diff_lt_ulp_eq_neg = Q.store_thm("diff_lt_ulp_eq_neg",
2807   `!a: ('t, 'w) float b: ('t, 'w) float x.
2808       ~exponent_boundary b a /\
2809       abs (x - float_to_real a) < ULP (a.Exponent, (:'t)) /\
2810       abs (x - float_to_real b) < ULP (a.Exponent, (:'t)) /\
2811       float_to_real a <= x /\ float_to_real b <= x /\
2812       ~float_is_zero b ==>
2813       (a = b)`,
2814   tac)
2815*)
2816
2817Theorem exponent_boundary_lt[local]:
2818    !a b.
2819      exponent_boundary a b ==> abs (float_to_real a) < abs (float_to_real b)
2820Proof
2821   rrw [float_to_real_def, exponent_boundary_def, abs_float_value,
2822        abs_significand, ABS_MUL, ABS_DIV,
2823        ABS_N, gt0_abs]
2824   >- (match_mp_tac lt_mult
2825       \\ rsimp [nobias_denormal_lt_1, REAL_LT_DIV])
2826   \\ simp [REAL_LT_LMUL, REAL_LT_RDIV_EQ, cancel_rwts,
2827            POW_ADD, REAL_ARITH ``1 + x < 2 <=> x < 1r``,
2828            nobias_denormal_lt_1]
2829QED
2830
2831Theorem exponent_boundary_not_float_zero[local]:
2832    !x y. exponent_boundary x y ==> ~float_is_zero y
2833Proof
2834   rw [exponent_boundary_def, float_is_zero]
2835   \\ strip_tac
2836   \\ fs []
2837QED
2838
2839Theorem ULP_lt_float_to_real[local]:
2840    !y:('t,'w) float.
2841       ~float_is_zero y ==> ULP (y.Exponent,(:'t)) <= abs (float_to_real y)
2842Proof
2843   rw [ULP_def, float_to_real_def, abs_float_value, abs_significand,
2844       ABS_MUL, ABS_DIV, ABS_N,
2845       gt0_abs, REAL_LE_LDIV_EQ, float_is_zero, GSYM word_lt0]
2846   \\ simp [POW_ADD, cancel_rwt, cancel_rwts]
2847   \\ simp [GSYM REAL_LDISTRIB, POW_2_LE1,
2848            le_mult, REAL_ARITH ``1r <= x /\ 0 <= n ==> 1 <= x + n``]
2849QED
2850
2851(* |- !y. ~float_is_zero y ==> ulp (:'t # 'w) <= abs (float_to_real y) *)
2852val ulp_lt_float_to_real =
2853   diff_float_ULP
2854   |> Q.SPEC `float_plus_zero (:'t # 'w)`
2855   |> SIMP_RULE (srw_ss())
2856         [ABS_NEG, float_components, zero_to_real, zero_properties,
2857          exponent_boundary_def, GSYM ulp_def, GSYM float_is_zero_to_real]
2858
2859val abs_limits = REAL_ARITH ``!x l. abs x <= l <=> ~(x < -l) /\ ~(x > l)``
2860
2861val abs_limits2 =
2862   REAL_ARITH ``!x l. abs x < l <=> ~(x <= -l) /\ ~(x >= l)``
2863
2864(* ------------------------------------------------------------------------
2865   Rounding to regular value
2866   ------------------------------------------------------------------------ *)
2867
2868Theorem round_roundTowardZero:
2869    !y: ('t, 'w) float x r.
2870      (float_value y = Float r) /\
2871      abs (r - x) < ULP (y.Exponent, (:'t)) /\ abs r <= abs x /\
2872      ulp (:'t # 'w) <= abs x /\ abs x <= largest (:'t # 'w) ==>
2873      (round roundTowardZero x = y)
2874Proof
2875   lrw [round_def, closest_def, is_closest_def, closest_such_def]
2876   >- imp_res_tac abs_limits
2877   >- imp_res_tac abs_limits
2878   \\ SELECT_ELIM_TAC
2879   \\ rw []
2880   >| [
2881      qexists_tac `y`
2882      \\ imp_res_tac Float_float_to_real
2883      \\ rw [Float_is_finite]
2884      \\ Cases_on `float_to_real b = float_to_real y`
2885      >- simp []
2886      \\ Cases_on `exponent_boundary b y`
2887      >- (`ULP (y.Exponent,(:'t)) <= abs (float_to_real y)`
2888          by metis_tac [ULP_lt_float_to_real, exponent_boundary_not_float_zero]
2889          \\ match_mp_tac
2890               (REAL_ARITH
2891                   ``abs (a - x) < abs x /\ abs b < abs a /\ abs a <= abs x ==>
2892                     abs (a - x) <= abs (b - x)``)
2893          \\ rsimp [exponent_boundary_lt])
2894      \\ match_mp_tac
2895            (REAL_ARITH
2896               ``ULP ((y: ('t, 'w) float).Exponent, (:'t)) <= abs (ra - rb) /\
2897                 abs (ra - x) < ULP (y.Exponent, (:'t)) /\
2898                 abs ra <= abs x /\ abs rb <= abs x ==>
2899                 abs (ra - x) <= abs (rb - x)``)
2900      \\ simp []
2901      \\ metis_tac [diff_float_ULP],
2902      (* -- *)
2903      Cases_on `float_is_zero y`
2904      >- (
2905          `r = 0` by (pop_assum mp_tac \\ simp [float_is_zero_def])
2906          \\ `y.Exponent = 0w`
2907          by (qpat_x_assum `float_is_zero y` mp_tac \\ simp [float_is_zero])
2908          \\ rlfs [ulp_def]
2909          \\ metis_tac [REAL_ARITH ``~(x < b /\ b <= x: real)``]
2910         )
2911      \\ imp_res_tac Float_float_to_real
2912      \\ pop_assum (SUBST_ALL_TAC o SYM)
2913      \\ `abs (float_to_real x' - x) <= abs (float_to_real y - x)`
2914      by metis_tac [Float_is_finite]
2915      \\ `abs (x - float_to_real x') < ULP (y.Exponent, (:'t))`
2916      by metis_tac [REAL_LET_TRANS, ABS_SUB]
2917      \\ Cases_on `exponent_boundary x' y`
2918      >- (
2919          `ULP (y.Exponent,(:'t)) <= abs (float_to_real y)`
2920          by metis_tac [ULP_lt_float_to_real, exponent_boundary_not_float_zero]
2921          \\ `abs (float_to_real y - x) <= abs (float_to_real x' - x)`
2922          by (match_mp_tac
2923               (REAL_ARITH
2924                   ``abs (a - x) < abs x /\ abs b < abs a /\ abs a <= abs x ==>
2925                     abs (a - x) <= abs (b - x)``)
2926              \\ rsimp [exponent_boundary_lt]
2927             )
2928          \\ simp [GSYM float_to_real_11_right]
2929          \\ match_mp_tac
2930               (REAL_ARITH
2931                   ``abs a <= abs x /\ abs b <= abs x /\
2932                     (abs (a - x) = abs (b - x)) ==> (a = b)``)
2933          \\ rsimp []
2934         )
2935      \\ match_mp_tac diff_lt_ulp_eq0
2936      \\ qexists_tac `x`
2937      \\ rsimp []
2938   ]
2939QED
2940
2941(*
2942val ULP01 = Q.store_thm("ULP01",
2943   `ULP (0w:'w word, (:'t)) = ULP (1w:'w word, (:'t))`,
2944   rw [ULP_def]
2945   )
2946
2947val ULP_lt_mono = Q.store_thm("ULP_lt_mono",
2948   `!e1 e2. 1w <+ e2 ==> (ULP (e1, (:'t)) < ULP (e2, (:'t)) = e1 <+ e2)`,
2949   Cases
2950   \\ Cases
2951   \\ lrw [ULP_def, wordsTheory.word_lo_n2w, REAL_LT_RDIV]
2952   \\ simp [REAL_OF_NUM_POW]
2953   )
2954
2955val exponent_boundary_exp_gt1 = Q.prove(
2956   `!b y: ('t, 'w) float.
2957       exponent_boundary b y ==> b.Exponent <+ y.Exponent /\ 1w <+ y.Exponent`,
2958   rw [exponent_boundary_def, wordsTheory.WORD_LO]
2959   \\ Cases_on `b.Exponent`
2960   \\ Cases_on `y.Exponent`
2961   \\ lfs []
2962   )
2963*)
2964
2965Theorem word_lsb_plus_1[local]:
2966    !a. word_lsb (a + 1w) = ~word_lsb a
2967Proof
2968   Cases
2969   \\ simp [wordsTheory.word_add_n2w, arithmeticTheory.ODD,
2970            GSYM arithmeticTheory.ADD1]
2971QED
2972
2973Theorem word_lsb_minus_1[local]:
2974    !a. word_lsb (a - 1w) = ~word_lsb a
2975Proof
2976   Cases
2977   \\ Cases_on `n`
2978   \\ simp [GSYM wordsTheory.word_add_n2w, arithmeticTheory.ODD,
2979            arithmeticTheory.ADD1]
2980QED
2981
2982val tac =
2983   qpat_x_assum `!a. q \/ t` (qspec_then `y` strip_assume_tac)
2984   \\ fs [REAL_NOT_LE]
2985   \\ qpat_x_assum `!b. p` (qspec_then `b` imp_res_tac)
2986   \\ rlfs []
2987   \\ rfs []
2988
2989Theorem round_roundTiesToEven:
2990    !y: ('t, 'w) float x r.
2991      (float_value y = Float r) /\
2992      ((y.Significand = 0w) /\ y.Exponent <> 1w ==> abs r <= abs x) /\
2993      2 * abs (r - x) <= ULP (y.Exponent, (:'t)) /\
2994      ((2 * abs (r - x) = ULP (y.Exponent, (:'t))) ==>
2995       ~word_lsb (y.Significand)) /\
2996      ulp (:'t # 'w) < 2 * abs x /\ abs x < threshold (:'t # 'w) ==>
2997      (round roundTiesToEven x = y)
2998Proof
2999   lrw [round_def, closest_def, is_closest_def, closest_such_def,
3000        SPECIFICATION]
3001   >- imp_res_tac abs_limits2
3002   >- imp_res_tac abs_limits2
3003   \\ SELECT_ELIM_TAC
3004   \\ rw []
3005   >| [
3006      qexists_tac `y`
3007      \\ imp_res_tac Float_float_to_real
3008      \\ `float_is_finite y` by simp [Float_is_finite]
3009      \\ rw []
3010      >| [
3011         Cases_on `float_to_real y = float_to_real b`
3012         >- simp []
3013         \\ Cases_on `exponent_boundary b y`
3014         >- (
3015           `ULP (y.Exponent,(:'t)) <= abs (float_to_real y)`
3016           by metis_tac [ULP_lt_float_to_real, exponent_boundary_not_float_zero]
3017           \\ imp_res_tac exponent_boundary_lt
3018           \\ `2 * abs (float_to_real y - x) <= abs (float_to_real y)`
3019           by imp_res_tac REAL_LE_TRANS
3020           \\ match_mp_tac
3021                (REAL_ARITH
3022                    ``abs b < abs a /\ abs a <= abs x /\
3023                      2 * abs (a - x) <= abs a ==> abs (a - x) <= abs (b - x)``)
3024           \\ rlfs [exponent_boundary_def]
3025           )
3026         \\ metis_tac
3027               [diff_float_ULP,
3028                REAL_ARITH
3029                   ``2 * abs (r - x) <= u /\ u <= abs (r - b) ==>
3030                     abs (r - x) <= abs (b - x)``],
3031         (* -- *)
3032         strip_tac
3033         \\ fs []
3034         \\ `a' <> y` by metis_tac []
3035         \\ Cases_on `float_is_zero y`
3036         >- fs [float_is_zero]
3037         \\ `float_to_real y <> float_to_real a'`
3038         by simp [float_to_real_eq]
3039         \\ Cases_on `exponent_boundary a' y`
3040         >- fs [exponent_boundary_def]
3041         \\ imp_res_tac diff_float_ULP
3042         \\ `2 * abs (float_to_real y - x) <  ULP (y.Exponent,(:'t))`
3043         by rsimp []
3044         \\ metis_tac
3045              [REAL_ARITH
3046                  ``2 * abs (r - x) < u /\ u <= abs (r - a) ==>
3047                    ~(abs (a - x) <= abs (r - x))``]
3048      ],
3049      (* -- *)
3050      `float_is_finite y` by simp [Float_is_finite]
3051      \\ Cases_on `float_is_zero y`
3052      >- (`r = 0` by (pop_assum mp_tac \\ simp [float_is_zero_def])
3053          \\ `y.Exponent = 0w`
3054          by (qpat_x_assum `float_is_zero y` mp_tac \\ simp [float_is_zero])
3055          \\ rlfs [ulp_def]
3056          \\ metis_tac [ULP_gt0,
3057                REAL_ARITH ``~(0 < x /\ x < b /\ 2 * b <= x: real)``])
3058      \\ imp_res_tac Float_float_to_real
3059      \\ pop_assum (SUBST_ALL_TAC o SYM)
3060      \\ Cases_on `exponent_boundary x' y`
3061      >- (`abs (float_to_real y) <= abs x` by fs [exponent_boundary_def]
3062          \\ metis_tac
3063               [REAL_LE_TRANS, exponent_boundary_lt,
3064                ULP_lt_float_to_real,
3065                REAL_ARITH
3066                  ``~(2 * abs (a - x) <= abs a /\ abs a <= abs x /\
3067                      abs b < abs a /\ abs (b - x) <= abs (a - x))``])
3068      \\ Cases_on `2 * abs (float_to_real y - x) < ULP (y.Exponent,(:'t))`
3069      >- (`2 * abs (float_to_real x' - x) <= 2 * abs (float_to_real y - x)`
3070          by metis_tac [Float_is_finite,
3071               REAL_ARITH ``2 * abs a <= 2 * abs b <=> abs a <= abs b``]
3072          \\ metis_tac [REAL_LET_TRANS, diff_lt_ulp_even])
3073      \\ `2 * abs (float_to_real y - x) = ULP (y.Exponent,(:'t))` by rsimp []
3074      \\ fs []
3075      \\ Cases_on `float_to_real y = float_to_real x'`
3076      >- (fs [float_to_real_eq] \\ fs [])
3077      \\ imp_res_tac diff_float_ULP
3078      \\ `abs (float_to_real x' - x) <= abs (float_to_real y - x)`
3079      by metis_tac []
3080      \\ `abs (float_to_real y - float_to_real x') =
3081          ULP (y.Exponent,(:'t))`
3082      by metis_tac
3083           [REAL_ARITH
3084               ``(2 * abs (a - x) = u) /\ u <= abs (a - b) /\
3085                 abs (b - x) <= abs (a - x) ==> (abs (a - b) = u)``]
3086      \\ `y.Significand <> -1w` by (strip_tac \\ fs [])
3087      \\ `~float_is_zero x'`
3088      by (strip_tac
3089          \\ imp_res_tac float_is_zero_to_real
3090          \\ fs [float_min_equiv_ULP_eq_float_to_real]
3091          \\ fs [float_components])
3092      \\ Cases_on `y.Significand = 0w`
3093      >- (qpat_x_assum `~float_is_zero y` assume_tac
3094          \\ `x'.Significand <> -1w` by (strip_tac \\ fs [] \\ tac)
3095          \\ Cases_on `y.Exponent = 1w`
3096          >- (fs [diff_ulp_next_float01] \\ tac)
3097          \\ Cases_on `abs (float_to_real x') < abs (float_to_real y)`
3098          \\ fs [REAL_NOT_LT]
3099          >- metis_tac
3100               [REAL_LE_TRANS, ULP_lt_float_to_real,
3101                REAL_ARITH
3102                  ``~(2 * abs (a - x) <= abs a /\ abs a <= abs x /\
3103                      abs b < abs a /\ abs (b - x) <= abs (a - x))``]
3104          \\ `abs (float_to_real x' - x) = abs (float_to_real y - x)`
3105          by imp_res_tac
3106               (REAL_ARITH
3107                   ``(2 * abs (a - x) = u) /\ abs (b - x) <= abs (a - x) /\
3108                     (abs (a - b) = u) ==> (abs (b - x) = abs (a - x))``)
3109         \\ fs [diff_ulp_next_float0]
3110         \\ tac
3111         )
3112      \\ `y.Significand NOTIN {0w; -1w}` by simp []
3113      \\ fs [diff_ulp_next_float]
3114      \\ `word_lsb x'.Significand`
3115      by simp [word_lsb_plus_1, SIMP_RULE (srw_ss()) [] word_lsb_minus_1]
3116      \\ fs []
3117      \\ tac
3118   ]
3119QED
3120
3121val not_one_lem = wordsLib.WORD_DECIDE ``(x:'a word) <> 1w ==> w2n x <> 1``
3122val pow_add1 = REWRITE_RULE [arithmeticTheory.ADD1] pow
3123
3124Theorem exponent_boundary_ULPs[local]:
3125    !x y. exponent_boundary x y ==>
3126          (ULP (y.Exponent, (:'t)) = 2 * ULP (x.Exponent, (:'t)))
3127Proof
3128   srw_tac [] [exponent_boundary_def, ULP_def, pow_add1, mult_ratr]
3129   \\ fs [not_one_lem]
3130QED
3131
3132Theorem round_roundTiesToEven0:
3133    !y: ('t, 'w) float x r.
3134      (float_value y = Float r) /\
3135      ((y.Significand = 0w) /\ y.Exponent <> 1w /\ ~(abs r <= abs x)) /\
3136      4 * abs (r - x) <= ULP (y.Exponent, (:'t)) /\
3137      ulp (:'t # 'w) < 2 * abs x /\ abs x < threshold (:'t # 'w) ==>
3138      (round roundTiesToEven x = y)
3139Proof
3140   lrw [round_def, closest_def, is_closest_def, closest_such_def,
3141        SPECIFICATION]
3142   >- imp_res_tac abs_limits2
3143   >- imp_res_tac abs_limits2
3144   \\ SELECT_ELIM_TAC
3145   \\ rw []
3146   >| [
3147      qexists_tac `y`
3148      \\ imp_res_tac Float_float_to_real
3149      \\ `float_is_finite y` by simp [Float_is_finite]
3150      \\ rw []
3151      \\ Cases_on `float_to_real y = float_to_real b`
3152      >- simp []
3153      \\ Cases_on `exponent_boundary b y`
3154      >- (
3155        imp_res_tac diff_exponent_boundary
3156        \\ `2 * ULP (b.Exponent,(:'t)) = ULP (y.Exponent,(:'t))`
3157        by (fs [exponent_boundary_def, ULP_def]
3158            \\ rw [pow_add1, mult_ratr]
3159            \\ fs [not_one_lem])
3160        \\ match_mp_tac
3161             (REAL_ARITH
3162                ``~(abs a <= abs x) /\ 4 * abs (a - x) <= 2 * abs (a - b) ==>
3163                  abs (a - x) <= abs (b - x)``)
3164        \\ simp []
3165        )
3166      \\ metis_tac
3167            [diff_float_ULP,
3168             REAL_ARITH ``4 * abs (r - x) <= u /\ u <= abs (r - b) ==>
3169                                  abs (r - x) <= abs (b - x)``],
3170      (* -- *)
3171      `float_is_finite y` by simp [Float_is_finite]
3172      \\ Cases_on `float_is_zero y`
3173      >- (`r = 0` by (pop_assum mp_tac \\ simp [float_is_zero_def])
3174          \\ `y.Exponent = 0w`
3175          by (qpat_x_assum `float_is_zero y` mp_tac \\ simp [float_is_zero])
3176          \\ rlfs [ulp_def])
3177      \\ imp_res_tac Float_float_to_real
3178      \\ pop_assum (SUBST_ALL_TAC o SYM)
3179      \\ `abs (float_to_real x' - x) <= abs (float_to_real y - x)` by res_tac
3180      \\ `4 * abs (float_to_real x' - x) <= 4 * abs (float_to_real y - x)`
3181      by rsimp []
3182      \\ `4 * abs (float_to_real x' - x) <= ULP (y.Exponent,(:'t))`
3183      by metis_tac [REAL_LE_TRANS]
3184      \\ match_mp_tac diff_lt_ulp_even4
3185      \\ qexists_tac `x`
3186      \\ simp []
3187      \\ spose_not_then assume_tac
3188      \\ imp_res_tac exponent_boundary_ULPs
3189      \\ fs [REAL_ARITH ``4r * a <= 2 * b <=> 2 * a <= b``]
3190      \\ imp_res_tac diff_exponent_boundary
3191      \\ `abs (float_to_real x' - x) = abs (float_to_real y - x)`
3192      by metis_tac
3193           [REAL_ARITH
3194              ``2 * abs (a - x) <= u /\
3195                2 * abs (b - x) <= u /\
3196                (abs (a - b) = u) ==> (abs (a - x) = abs (b - x))``]
3197      \\ `~word_lsb y.Significand ==> ~word_lsb x'.Significand`
3198      by metis_tac []
3199      \\ rfs [exponent_boundary_def]
3200   ]
3201QED
3202
3203(*
3204
3205val round_roundTowardPositive = Q.store_thm("round_roundTowardPositive",
3206   `!y: ('t, 'w) float x r.
3207      (float_value y = Float r) /\
3208      abs (x - r) < ulp (:'t # 'w) /\ r >= x /\
3209      ulp (:'t # 'w) <= abs x /\ abs x <= largest (:'t # 'w) ==>
3210      (round roundTowardPositive x = y)`,
3211   tac (REAL_ARITH
3212          ``ulp (:'t # 'w) <= abs (ra - rb) /\
3213            abs (x - ra) < ulp (:'t # 'w) /\
3214            ra >= x /\ rb >= x ==>
3215            abs (ra - x) <= abs (rb - x)``)
3216       diff_lt_ulp_eq_pos)
3217
3218val round_roundTowardNegative = Q.store_thm("round_roundTowardNegative",
3219   `!y: ('t, 'w) float x r.
3220      (float_value y = Float r) /\
3221      abs (x - r) < ulp (:'t # 'w) /\ r <= x /\
3222      ulp (:'t # 'w) <= abs x /\ abs x <= largest (:'t # 'w) ==>
3223      (round roundTowardNegative x = y)`,
3224   tac (REAL_ARITH
3225          ``ulp (:'t # 'w) <= abs (ra - rb) /\
3226            abs (x - ra) < ulp (:'t # 'w) /\
3227            ra <= x /\ rb <= x ==>
3228            abs (ra - x) <= abs (rb - x)``)
3229       diff_lt_ulp_eq_neg)
3230
3231val tac =
3232   REPEAT strip_tac
3233   \\ Cases_on `float_is_zero y`
3234   >- (
3235       `r = 0` by (pop_assum mp_tac \\ simp [float_is_zero_def])
3236       \\ qpat_x_assum `float_is_zero y` mp_tac
3237       \\ rw [float_is_zero]
3238       \\ fs [ulp_def]
3239       \\ prove_tac [REAL_ARITH ``~(x < b /\ b <= x: real)``,
3240                     REAL_ARITH ``2 * abs (-x) <= u ==> abs x <= u``]
3241      )
3242   \\ lrw [float_round_def]
3243   \\ metis_tac [zero_properties, round_roundTowardZero, round_roundTiesToEven
3244                 (*, round_roundTowardPositive, round_roundTowardNegative *)]
3245
3246
3247val float_round_roundTowardZero = Q.store_thm(
3248   "float_round_roundTowardZero",
3249   `!b y: ('t, 'w) float x r.
3250      (float_value y = Float r) /\
3251      abs (x - r) < ULP (y.Exponent, (:'t)) /\ abs r <= abs x /\
3252      ulp (:'t # 'w) <= abs x /\ abs x <= largest (:'t # 'w) ==>
3253      (float_round roundTowardZero b x = y)`,
3254   tac
3255   )
3256
3257val float_round_roundTiesToEven = Q.store_thm("float_round_roundTiesToEven",
3258   `!b y: ('t, 'w) float x r.
3259      (float_value y = Float r) /\
3260      ((y.Significand = 0w) /\ y.Exponent <> 1w ==> abs r <= abs x) /\
3261      2 * abs (r - x) <= ULP (y.Exponent, (:'t)) /\
3262      ((2 * abs (r - x) = ULP (y.Exponent, (:'t))) ==>
3263       ~word_lsb (y.Significand)) /\
3264      ulp (:'t # 'w) < abs x /\ abs x < threshold (:'t # 'w) ==>
3265      (float_round roundTiesToEven b x = y)`,
3266   tac
3267   )
3268
3269val float_round_roundTowardPositive = Q.store_thm(
3270   "float_round_roundTowardPositive",
3271   `!b y: ('t, 'w) float x r.
3272      (float_value y = Float r) /\
3273      abs (x - r) < ulp (:'t # 'w) /\ r >= x /\
3274      ulp (:'t # 'w) <= abs x /\ abs x <= largest (:'t # 'w) ==>
3275      (float_round roundTowardPositive b x = y)`,
3276   tac)
3277
3278val float_round_roundTowardNegative = Q.store_thm(
3279   "float_round_roundTowardNegative",
3280   `!b y: ('t, 'w) float x r.
3281      (float_value y = Float r) /\
3282      abs (x - r) < ulp (:'t # 'w) /\ r <= x /\
3283      ulp (:'t # 'w) <= abs x /\ abs x <= largest (:'t # 'w) ==>
3284      (float_round roundTowardNegative b x = y)`,
3285   tac)
3286
3287*)
3288
3289(* ------------------------------------------------------------------------
3290   Rounding to +/- 0
3291   ------------------------------------------------------------------------ *)
3292
3293Theorem round_roundTowardZero_is_zero:
3294    !x. abs x < ulp (:'t # 'w) ==>
3295        (round roundTowardZero x = float_plus_zero (:'t # 'w)) \/
3296        (round roundTowardZero x = float_minus_zero (:'t # 'w))
3297Proof
3298   REPEAT strip_tac
3299   \\ qabbrev_tac `r: ('t, 'w) float = round roundTowardZero x`
3300   \\ pop_assum (mp_tac o SYM o REWRITE_RULE [markerTheory.Abbrev_def])
3301   \\ simp [round_def, lt_ulp_not_infinity0,
3302            closest_such_def, closest_def, is_closest_def]
3303   \\ SELECT_ELIM_TAC
3304   \\ rw []
3305   >| [
3306      qexists_tac `float_plus_zero (:'t # 'w)`
3307      \\ simp [zero_properties, zero_to_real, ABS_POS,
3308               ABS_NEG]
3309      \\ REPEAT strip_tac
3310      \\ imp_res_tac REAL_LET_TRANS
3311      \\ imp_res_tac less_than_ulp
3312      \\ Cases_on `b`
3313      \\ lfs [float_to_real_def, ABS_NEG],
3314      (* -- *)
3315      imp_res_tac REAL_LET_TRANS
3316      \\ imp_res_tac less_than_ulp
3317      \\ Cases_on `r`
3318      \\ lfs [float_plus_zero_def, float_minus_zero_def, float_negate_def,
3319              float_component_equality]
3320      \\ wordsLib.Cases_on_word_value `c`
3321      \\ simp []
3322   ]
3323QED
3324
3325Theorem is_closest_finite_AND:
3326  is_closest float_is_finite r f /\ Q f ==>
3327  is_closest { a | float_is_finite a /\ Q a} r f
3328Proof
3329  simp[is_closest_def, IN_DEF]
3330QED
3331
3332Theorem float_to_real_round0[simp]:
3333  float_to_real (round m 0) = 0
3334Proof
3335  Cases_on ‘m’ >>
3336  simp[round_def,
3337       SRULE [ulp_positive] (Q.SPEC ‘0’ lt_ulp_not_infinity0),
3338       SRULE [ulp_positive] (Q.SPEC ‘0’ lt_ulp_not_infinity1),
3339       closest_such_def, closest_def] >>
3340  SELECT_ELIM_TAC >>
3341  dsimp[is_closestP_finite_float_exists, is_closest_float_is_finite_0] >>
3342  qexists_tac ‘POS0’ >>
3343  simp[is_closest_finite_AND, is_closest_float_is_finite_0] >>
3344  simp[is_closest_def] >> rpt strip_tac >>
3345  first_x_assum $ qspec_then ‘POS0’ mp_tac >> gs[REAL_ABS_LE0]
3346QED
3347
3348Theorem float_to_real_min_pos[local]:
3349    !r: ('t, 'w) float.
3350       (abs (float_to_real r) = ulp (:'t # 'w)) <=>
3351       r IN {float_plus_min (:'t # 'w);
3352             float_negate (float_plus_min (:'t # 'w))}
3353Proof
3354   rw [float_plus_min_def, float_negate_def, ulp_def, ULP_def,
3355       float_to_real_def, float_component_equality, abs_float_value,
3356       abs_significand, ABS_MUL, ABS_DIV,
3357       ABS_N, gt0_abs]
3358   >| [
3359      wordsLib.Cases_on_word_value `r.Sign`
3360      \\ Cases_on `r.Significand = 1w`
3361      \\ simp [mult_rat, POW_ADD,
3362               div_twopow
3363               |> Q.SPEC `m + n`
3364               |> REWRITE_RULE [POW_ADD]
3365               |> Drule.GEN_ALL]
3366      \\ Cases_on `r.Significand`
3367      \\ fs [],
3368      simp [REAL_EQ_RDIV_EQ]
3369      \\ simp [POW_ADD, GSYM REAL_LDISTRIB,
3370               cancel_rwts, cancel_rwt]
3371      \\ match_mp_tac (REAL_ARITH ``2r < a * b ==> b * a <> 2``)
3372      \\ match_mp_tac ge2d
3373      \\ simp [REAL_OF_NUM_POW, pow_ge2, exp_ge2,
3374               DECIDE ``2n <= a ==> 2 <= b + a``,
3375               fcpTheory.DIMINDEX_GE_1 ]
3376   ]
3377QED
3378
3379val compare_with_zero_tac =
3380   qpat_x_assum `!b. float_is_finite b  ==> p`
3381      (fn th =>
3382         assume_tac
3383            (SIMP_RULE (srw_ss())
3384               [ABS_NEG, zero_to_real, zero_properties]
3385                  (Q.SPEC `float_plus_zero (:'t # 'w)` th))
3386         \\ assume_tac th)
3387
3388Theorem half_ulp[local]:
3389  !x r: ('t, 'w) float.
3390    (2 * abs x = ulp (:'t # 'w)) /\
3391    (!b: ('t, 'w) float.
3392       float_is_finite b ==>
3393       abs (float_to_real r - x) <= abs (float_to_real b - x)) ==>
3394    float_is_zero r \/
3395    r IN {float_plus_min (:'t # 'w);
3396          float_negate (float_plus_min (:'t # 'w))}
3397Proof
3398   REPEAT strip_tac
3399   \\ Cases_on `float_is_zero r`
3400   \\ simp []
3401   \\ compare_with_zero_tac
3402   \\ `abs (float_to_real r) = ulp (:'t # 'w)`
3403   by metis_tac
3404         [ulp_lt_float_to_real,
3405          REAL_ARITH
3406           ``(2 * abs x = u) /\ u <= abs r /\ abs (r - x) <= abs x ==>
3407             (abs r = u)``]
3408   \\ fs [float_to_real_min_pos]
3409QED
3410
3411Theorem min_pos_odd[local]:
3412  !r: ('t, 'w) float.
3413    r IN {float_plus_min (:'t # 'w);
3414          float_negate (float_plus_min (:'t # 'w))} ==>
3415    word_lsb r.Significand
3416Proof rw [float_plus_min_def, float_negate_def] \\ simp []
3417QED
3418
3419Theorem round_roundTiesToEven_is_zero:
3420    !x. 2 * abs x <= ulp (:'t # 'w) ==>
3421        (round roundTiesToEven x = float_plus_zero (:'t # 'w)) \/
3422        (round roundTiesToEven x = float_minus_zero (:'t # 'w))
3423Proof
3424   REPEAT strip_tac
3425   \\ qabbrev_tac `r: ('t, 'w) float = round roundTiesToEven x`
3426   \\ pop_assum (mp_tac o SYM o REWRITE_RULE [markerTheory.Abbrev_def])
3427   \\ simp [round_def, lt_ulp_not_infinity1, SPECIFICATION,
3428            closest_such_def, closest_def, is_closest_def]
3429   \\ SELECT_ELIM_TAC
3430   \\ rw []
3431   >| [
3432      qexists_tac `float_plus_zero (:'t # 'w)`
3433      \\ simp [zero_properties, zero_to_real, ABS_POS,
3434               ABS_NEG]
3435      \\ rw [float_plus_zero_def]
3436      \\ Cases_on `float_is_zero b`
3437      \\ rsimp [float_is_zero_to_real_imp]
3438      \\ metis_tac
3439           [ULP_lt_float_to_real, ulp_lt_ULP, REAL_LE_TRANS,
3440            REAL_ARITH ``2 * abs x <= abs c ==> abs x <= abs (c - x)``],
3441      (* -- *)
3442      Cases_on `float_is_zero r`
3443      >- fs [float_sets]
3444      \\ Cases_on `2 * abs x < ulp (:'t # 'w)`
3445      >| [
3446         imp_res_tac ulp_lt_float_to_real
3447         \\ compare_with_zero_tac
3448         \\ metis_tac
3449               [REAL_ARITH
3450                  ``~(2 * abs x < u /\ u <= abs r /\ abs (r - x) <= abs x)``],
3451         (* -- *)
3452         imp_res_tac
3453            (REAL_ARITH ``a <= b /\ ~(a < b) ==> (a = b: real)``)
3454         \\ imp_res_tac half_ulp
3455         \\ imp_res_tac min_pos_odd
3456         \\ compare_with_zero_tac
3457         \\ fs []
3458         \\ qpat_x_assum `!a. q \/ t`
3459               (qspec_then `float_plus_zero (:'t # 'w)` strip_assume_tac)
3460         \\ rfs [ABS_NEG, zero_properties, zero_to_real,
3461                 float_components, GSYM ulp, GSYM neg_ulp]
3462         \\ qpat_x_assum `!b. p` (qspec_then `b` imp_res_tac)
3463         \\ metis_tac
3464              [REAL_ARITH
3465                  ``~((2 * abs x = u) /\ abs (u - x) <= abs x /\
3466                      ~(abs x <= abs (b - x)) /\ abs (u - x) <= abs (b - x))``,
3467               REAL_ARITH
3468                  ``~((2 * abs x = u) /\ abs (-u - x) <= abs x /\
3469                      ~(abs x <= abs (b - x)) /\ abs (-u - x) <= abs (b - x))``]
3470      ]
3471   ]
3472QED
3473
3474val tac =
3475   lrw [float_round_def]
3476   \\ metis_tac [round_roundTowardZero_is_zero, round_roundTiesToEven_is_zero,
3477                 zero_properties]
3478
3479Theorem round_roundTowardZero_is_minus_zero:
3480    !x. abs x < ulp (:'t # 'w) ==>
3481        (float_round roundTowardZero T x = float_minus_zero (:'t # 'w))
3482Proof
3483   tac
3484QED
3485
3486Theorem round_roundTowardZero_is_plus_zero:
3487    !x. abs x < ulp (:'t # 'w) ==>
3488        (float_round roundTowardZero F x = float_plus_zero (:'t # 'w))
3489Proof
3490   tac
3491QED
3492
3493Theorem round_roundTiesToEven_is_minus_zero:
3494    !x. 2 * abs x <= ulp (:'t # 'w) ==>
3495        (float_round roundTiesToEven T x = float_minus_zero (:'t # 'w))
3496Proof
3497   tac
3498QED
3499
3500Theorem round_roundTiesToEven_is_plus_zero:
3501    !x. 2 * abs x <= ulp (:'t # 'w) ==>
3502        (float_round roundTiesToEven F x = float_plus_zero (:'t # 'w))
3503Proof
3504   tac
3505QED
3506
3507(* ------------------------------------------------------------------------
3508   Rounding to limits
3509   ------------------------------------------------------------------------ *)
3510
3511Theorem largest_is_positive[simp]:
3512   0 <= largest (:'t # 'w)
3513Proof
3514   simp [largest_def, REAL_LE_MUL, REAL_LE_DIV,
3515         REAL_SUB_LE, POW_2_LE1,
3516         REAL_INV_1OVER, REAL_LE_LDIV_EQ,
3517         REAL_ARITH ``1r <= n ==> 1 <= 2 * n``]
3518QED
3519
3520Theorem threshold_is_positive[simp]:
3521   0 < threshold (:'t # 'w)
3522Proof
3523   simp [threshold_def, REAL_LT_MUL, REAL_LT_DIV,
3524         REAL_SUB_LT, POW_2_LE1,
3525         REAL_INV_1OVER, REAL_LT_LDIV_EQ, pow,
3526         REAL_ARITH ``1r <= n ==> 1 < 2 * (2 * n)``]
3527QED
3528
3529val tac =
3530   rrw  [round_def]
3531   \\ rlfs [largest_is_positive,
3532           REAL_ARITH ``0 <= l /\ l < x ==> ~(x < -l: real)``]
3533   \\ metis_tac [threshold_is_positive, largest_is_positive,
3534                 REAL_ARITH “(0r < i /\ x <= -i ==> ~(i <= x)) /\
3535                                     (0r <= i /\ x < -i ==> ~(i < x))”]
3536
3537Theorem round_roundTiesToEven_plus_infinity:
3538    !y: ('t, 'w) float x.
3539      threshold (:'t # 'w) <= x ==>
3540      (round roundTiesToEven x = float_plus_infinity (:'t # 'w))
3541Proof
3542   tac
3543QED
3544
3545Theorem round_roundTiesToEven_minus_infinity:
3546    !y: ('t, 'w) float x.
3547      x <= -threshold (:'t # 'w) ==>
3548      (round roundTiesToEven x = float_minus_infinity (:'t # 'w))
3549Proof
3550   tac
3551QED
3552
3553Theorem round_roundTowardZero_top:
3554  !y: ('t, 'w) float x.
3555    largest (:'t # 'w) < x ==> (round roundTowardZero x = float_top (:'t # 'w))
3556Proof tac
3557QED
3558
3559Theorem round_roundTowardZero_bottom:
3560    !y: ('t, 'w) float x.
3561      x < -largest (:'t # 'w) ==>
3562      (round roundTowardZero x = float_bottom (:'t # 'w))
3563Proof
3564   tac
3565QED
3566
3567Theorem round_roundTowardPositive_plus_infinity:
3568    !y: ('t, 'w) float x.
3569      largest (:'t # 'w) < x ==>
3570      (round roundTowardPositive x = float_plus_infinity (:'t # 'w))
3571Proof
3572   tac
3573QED
3574
3575Theorem round_roundTowardPositive_bottom:
3576    !y: ('t, 'w) float x.
3577      x < -largest (:'t # 'w) ==>
3578      (round roundTowardPositive x = float_bottom (:'t # 'w))
3579Proof
3580   tac
3581QED
3582
3583Theorem round_roundTowardNegative_top:
3584    !y: ('t, 'w) float x.
3585      largest (:'t # 'w) < x ==>
3586      (round roundTowardNegative x = float_top (:'t # 'w))
3587Proof
3588   tac
3589QED
3590
3591Theorem round_roundTowardNegative_minus_infinity:
3592   !y: ('t, 'w) float x.
3593     x < -largest (:'t # 'w) ==>
3594     (round roundTowardNegative x = float_minus_infinity (:'t # 'w))
3595Proof tac
3596QED
3597
3598val tac =
3599   lrw [float_round_def, round_roundTowardZero_top,
3600        round_roundTowardZero_bottom, round_roundTowardPositive_plus_infinity,
3601        round_roundTowardPositive_bottom, round_roundTowardNegative_top,
3602        round_roundTowardNegative_minus_infinity,
3603        top_properties, bottom_properties, infinity_properties]
3604
3605Theorem float_round_roundTowardZero_top:
3606    !b y: ('t, 'w) float x.
3607      largest (:'t # 'w) < x ==>
3608      (float_round roundTowardZero b x = float_top (:'t # 'w))
3609Proof
3610   tac
3611QED
3612
3613Theorem float_round_roundTowardZero_bottom:
3614    !b y: ('t, 'w) float x.
3615      x < -largest (:'t # 'w) ==>
3616      (float_round roundTowardZero b x = float_bottom (:'t # 'w))
3617Proof
3618   tac
3619QED
3620
3621Theorem float_round_roundTowardPositive_plus_infinity:
3622    !b y: ('t, 'w) float x.
3623      largest (:'t # 'w) < x ==>
3624      (float_round roundTowardPositive b x = float_plus_infinity (:'t # 'w))
3625Proof
3626   tac
3627QED
3628
3629Theorem float_round_roundTowardPositive_bottom:
3630    !b y: ('t, 'w) float x.
3631      x < -largest (:'t # 'w) ==>
3632      (float_round roundTowardPositive b x = float_bottom (:'t # 'w))
3633Proof
3634   tac
3635QED
3636
3637Theorem float_round_roundTowardNegative_top:
3638    !b y: ('t, 'w) float x.
3639      largest (:'t # 'w) < x ==>
3640      (float_round roundTowardNegative b x = float_top (:'t # 'w))
3641Proof
3642   tac
3643QED
3644
3645Theorem float_round_roundTowardNegative_minus_infinity:
3646    !b y: ('t, 'w) float x.
3647      x < -largest (:'t # 'w) ==>
3648      (float_round roundTowardNegative b x = float_minus_infinity (:'t # 'w))
3649Proof
3650   tac
3651QED
3652
3653(* ------------------------------------------------------------------------
3654   Theorem support for evaluation
3655   ------------------------------------------------------------------------ *)
3656
3657Theorem float_minus_zero:
3658    float_minus_zero (:'t # 'w) =
3659       <| Sign := 1w; Exponent := 0w; Significand := 0w |>
3660Proof
3661   simp [float_minus_zero_def, float_plus_zero_def, float_negate_def]
3662QED
3663
3664Theorem float_minus_infinity:
3665    float_minus_infinity (:'t # 'w) =
3666       <| Sign := 1w; Exponent := UINT_MAXw; Significand := 0w |>
3667Proof
3668   simp [float_minus_infinity_def, float_plus_infinity_def, float_negate_def]
3669QED
3670
3671Theorem float_round_non_zero:
3672     !mode toneg r s e f.
3673        (round mode r = <| Sign := s; Exponent := e; Significand := f |>) /\
3674        (e <> 0w \/ f <> 0w) ==>
3675        (float_round mode toneg r =
3676         <| Sign := s; Exponent := e; Significand := f |>)
3677Proof
3678    lrw [float_round_def, float_is_zero]
3679QED
3680
3681Theorem float_round_plus_infinity:
3682     !mode toneg r.
3683        (round mode r = float_plus_infinity (:'t # 'w)) ==>
3684        (float_round mode toneg r = float_plus_infinity (:'t # 'w))
3685Proof
3686    lrw [float_round_def, infinity_properties]
3687QED
3688
3689Theorem float_round_minus_infinity:
3690     !mode toneg r.
3691        (round mode r = float_minus_infinity (:'t # 'w)) ==>
3692        (float_round mode toneg r = float_minus_infinity (:'t # 'w))
3693Proof
3694    lrw [float_round_def, infinity_properties]
3695QED
3696
3697Theorem float_round_top:
3698     !mode toneg r.
3699        (round mode r = float_top (:'t # 'w)) ==>
3700        (float_round mode toneg r = float_top (:'t # 'w))
3701Proof
3702    lrw [float_round_def, top_properties]
3703QED
3704
3705Theorem float_round_bottom:
3706     !mode toneg r.
3707        (round mode r = float_bottom (:'t # 'w)) ==>
3708        (float_round mode toneg r = float_bottom (:'t # 'w))
3709Proof
3710    lrw [float_round_def, bottom_properties]
3711QED
3712
3713fun tac thms =
3714   rrw ([largest_def, threshold_def, float_to_real_def, wordsTheory.dimword_def,
3715         GSYM REAL_NEG_MINUS1, REAL_OF_NUM_POW,
3716         wordsLib.WORD_DECIDE ``x <> 1w ==> (x = 0w: word1)``] @ thms)
3717
3718Theorem float_to_real:
3719    !s e:'w word f:'t word.
3720      float_to_real <| Sign := s; Exponent := e; Significand := f |> =
3721         let r = if e = 0w
3722                    then 2r / &(2 EXP INT_MAX (:'w)) * (&w2n f / &dimword (:'t))
3723                 else &(2 EXP (w2n e)) / &(2 EXP INT_MAX (:'w)) *
3724                      (1r + &w2n f / &dimword (:'t))
3725         in
3726            if s = 1w then -r else r
3727Proof
3728   tac []
3729QED
3730
3731Theorem largest:
3732    largest (:'t # 'w) =
3733       &(2 EXP (UINT_MAX (:'w) - 1)) * (2 - 1 / &dimword (:'t)) /
3734       &(2 EXP INT_MAX (:'w))
3735Proof
3736   tac [REAL_INV_1OVER, mult_ratl]
3737QED
3738
3739Theorem threshold:
3740    threshold (:'t # 'w) =
3741       &(2 EXP (UINT_MAX (:'w) - 1)) * (2 - 1 / &(2 * dimword (:'t))) /
3742       &(2 EXP INT_MAX (:'w))
3743Proof
3744   tac [REAL_INV_1OVER, mult_ratl, arithmeticTheory.EXP]
3745QED
3746
3747Theorem largest_top_lem[local]:
3748   w2n (n2w (UINT_MAX (:'w)) + -1w : 'w word) = UINT_MAX (:'w) - 1
3749Proof
3750  simp_tac arith_ss
3751     [wordsTheory.WORD_LITERAL_ADD
3752      |> CONJUNCT2
3753      |> Q.SPECL [`UINT_MAX (:'w)`, `1`]
3754      |> SIMP_RULE std_ss [wordsTheory.ZERO_LT_UINT_MAX,
3755                           DECIDE ``0n < x ==> 1 <= x``],
3756      wordsTheory.w2n_n2w, wordsTheory.BOUND_ORDER,
3757      DECIDE ``a < b ==> (a - 1n < b)``]
3758QED
3759
3760Theorem largest_top_lem2[local]:
3761   &UINT_MAX (:'t) + 1 = &dimword (:'t) : real
3762Proof
3763  simp [wordsTheory.UINT_MAX_def, DECIDE ``1n < n ==> (n - 1 + 1 = n)``]
3764QED
3765
3766Theorem largest_is_top:
3767   1 < dimindex(:'w) ==>
3768   (largest (:'t # 'w) = float_to_real (float_top (:'t # 'w)))
3769Proof
3770  strip_tac
3771  \\ `dimword(:'w) <> 2`
3772  by fs [wordsTheory.dimword_def,
3773         arithmeticTheory.EXP_BASE_INJECTIVE
3774         |> Q.SPEC `2`
3775         |> REWRITE_RULE [DECIDE ``1n < 2``]
3776         |> Q.SPECL [`n`, `1`]
3777         |> REWRITE_RULE [arithmeticTheory.EXP_1]]
3778  \\ `2 < dimword(:'w)` by simp [DECIDE ``1 < n /\ n <> 2 ==> 2n < n``]
3779  \\ `UINT_MAXw - 1w <> 0w : 'w word` by simp []
3780  \\ asm_simp_tac std_ss [largest, float_top_def, float_to_real]
3781  \\ simp_tac std_ss [wordsTheory.word_T_def]
3782  \\ simp [REAL_EQ_LDIV_EQ, DECIDE ``0n < n ==> n <> 0``,
3783           REAL_SUB_LDISTRIB, REAL_ADD_LDISTRIB,
3784           REAL_EQ_SUB_RADD, REAL_DIV_REFL,
3785           mult_ratr, mult_ratl, wordsTheory.BOUND_ORDER,
3786           ONCE_REWRITE_RULE [REAL_MUL_COMM] mul_cancel,
3787           largest_top_lem]
3788  \\ simp_tac std_ss
3789       [GSYM REAL_ADD_ASSOC, REAL_DIV_ADD,
3790        GSYM REAL_MUL, largest_top_lem2,
3791        mul_cancel |> Q.SPECL [`a`, `&(n : num)`] |> SIMP_RULE (srw_ss()) [],
3792        wordsTheory.ZERO_LT_dimword, DECIDE ``0 < n ==> n <> 0n``,
3793        REAL_ARITH ``a * b + b = (a + 1r) * b``, REAL_DOUBLE]
3794QED
3795
3796Theorem largest_lt_threshold:
3797  largest (:'t # 'w) < threshold (:'t # 'w)
3798Proof
3799  rw [largest, threshold, REAL_LT_RDIV, REAL_LT_LMUL,
3800      REAL_ARITH ``a - b < a - c <=> c < b : real``,
3801      REAL_LT_RDIV_EQ, REAL_LT_LDIV_EQ,
3802      mult_ratl] >>
3803  fs[wordsTheory.dimword_def]
3804QED
3805
3806Theorem float_tests:
3807    (!s e f.
3808       float_is_nan <| Sign := s; Exponent := e; Significand := f |> <=>
3809       (e = -1w) /\ (f <> 0w)) /\
3810    (!s e f.
3811       float_is_signalling <| Sign := s; Exponent := e; Significand := f |> <=>
3812       (e = -1w) /\ ~word_msb f /\ (f <> 0w)) /\
3813    (!s e f.
3814       float_is_infinite <| Sign := s; Exponent := e; Significand := f |> <=>
3815       (e = -1w) /\ (f = 0w)) /\
3816    (!s e f.
3817       float_is_normal <| Sign := s; Exponent := e; Significand := f |> <=>
3818       (e <> 0w) /\ (e <> -1w)) /\
3819    (!s e f.
3820       float_is_subnormal <| Sign := s; Exponent := e; Significand := f |> <=>
3821       (e = 0w) /\ (f <> 0w)) /\
3822    (!s e f.
3823       float_is_zero <| Sign := s; Exponent := e; Significand := f |> <=>
3824       (e = 0w) /\ (f = 0w)) /\
3825    (!s e f.
3826       float_is_finite <| Sign := s; Exponent := e; Significand := f |> <=>
3827       (e <> -1w))
3828Proof
3829   rw [float_is_nan_def, float_is_signalling_def, float_is_infinite_def,
3830       float_is_finite_def, float_is_normal_def, float_is_subnormal_def,
3831       float_value_def]
3832   \\ rw [float_sets, float_minus_zero_def, float_plus_zero_def,
3833          float_is_finite_def, float_negate_def]
3834   \\ wordsLib.Cases_on_word_value `s`
3835   \\ simp []
3836QED
3837
3838Theorem float_infinity_negate_abs:
3839    (float_negate (float_plus_infinity (:'t # 'w)) =
3840     float_minus_infinity (:'t # 'w)) /\
3841    (float_negate (float_minus_infinity (:'t # 'w)) =
3842     float_plus_infinity (:'t # 'w)) /\
3843    (float_abs (float_plus_infinity (:'t # 'w)) =
3844     float_plus_infinity (:'t # 'w)) /\
3845    (float_abs (float_minus_infinity (:'t # 'w)) =
3846     float_plus_infinity (:'t # 'w))
3847Proof
3848    rw [float_plus_infinity_def, float_minus_infinity_def,
3849        float_negate_def, float_abs_def]
3850QED
3851
3852(* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -  *)
3853
3854Theorem float_round_to_integral_compute:
3855    (!m. float_round_to_integral m (float_minus_infinity (:'t # 'w)) =
3856         float_minus_infinity (:'t # 'w)) /\
3857    (!m. float_round_to_integral m (float_plus_infinity (:'t # 'w)) =
3858         float_plus_infinity (:'t # 'w)) /\
3859    (!m fp_op.
3860         float_round_to_integral m (float_some_qnan fp_op) =
3861         float_some_qnan fp_op)
3862Proof
3863   simp [float_round_to_integral_def, float_values]
3864QED
3865
3866(* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -  *)
3867
3868Theorem float_add_compute:
3869    (!mode x fp_op.
3870       float_add mode (float_some_qnan fp_op) x =
3871       (check_for_signalling [x],
3872        float_some_qnan (FP_Add mode (float_some_qnan fp_op) x)))
3873       /\
3874    (!mode x fp_op.
3875       float_add mode x (float_some_qnan fp_op) =
3876       (check_for_signalling [x],
3877        float_some_qnan (FP_Add mode x (float_some_qnan fp_op))))
3878       /\
3879    (!mode.
3880       float_add mode (float_minus_infinity (:'t # 'w))
3881                      (float_minus_infinity (:'t # 'w)) =
3882       (clear_flags, float_minus_infinity (:'t # 'w))) /\
3883    (!mode.
3884       float_add mode (float_minus_infinity (:'t # 'w))
3885                      (float_plus_infinity (:'t # 'w)) =
3886       (invalidop_flags,
3887        float_some_qnan (FP_Add mode (float_minus_infinity (:'t # 'w))
3888                                     (float_plus_infinity (:'t # 'w))))) /\
3889    (!mode.
3890       float_add mode (float_plus_infinity (:'t # 'w))
3891                      (float_plus_infinity (:'t # 'w)) =
3892       (clear_flags, float_plus_infinity (:'t # 'w))) /\
3893    (!mode.
3894       float_add mode (float_plus_infinity (:'t # 'w))
3895                      (float_minus_infinity (:'t # 'w)) =
3896       (invalidop_flags,
3897        float_some_qnan (FP_Add mode (float_plus_infinity (:'t # 'w))
3898                                     (float_minus_infinity (:'t # 'w)))))
3899Proof
3900   simp [float_add_def, float_values, float_components, some_nan_properties,
3901         check_for_signalling_def]
3902   \\ strip_tac
3903   \\ strip_tac
3904   \\ Cases_on `float_value x`
3905   \\ simp [float_is_signalling_def, float_is_nan_def]
3906QED
3907
3908Theorem float_add_nan:
3909    !mode x y.
3910       (float_value x = NaN) \/ (float_value y = NaN) ==>
3911       (float_add mode x y =
3912        (check_for_signalling [x; y], float_some_qnan (FP_Add mode x y)))
3913Proof
3914   NTAC 3 strip_tac
3915   \\ Cases_on `float_value x`
3916   \\ Cases_on `float_value y`
3917   \\ simp [float_add_def, check_for_signalling_def,
3918            float_is_signalling_def, float_is_nan_def]
3919QED
3920
3921Theorem float_add_finite:
3922    !mode x y r1 r2.
3923       (float_value x = Float r1) /\ (float_value y = Float r2) ==>
3924       (float_add mode x y =
3925        float_round_with_flags mode
3926          (if (r1 = 0) /\ (r2 = 0) /\ (x.Sign = y.Sign) then
3927             x.Sign = 1w
3928           else mode = roundTowardNegative) (r1 + r2))
3929Proof
3930   simp [float_add_def]
3931QED
3932
3933Theorem float_add_finite_plus_infinity:
3934    !mode x r.
3935       (float_value x = Float r) ==>
3936       (float_add mode x (float_plus_infinity (:'t # 'w)) =
3937        (clear_flags, float_plus_infinity (:'t # 'w)))
3938Proof
3939   simp [float_add_def, float_values]
3940QED
3941
3942Theorem float_add_plus_infinity_finite:
3943    !mode x r.
3944       (float_value x = Float r) ==>
3945       (float_add mode (float_plus_infinity (:'t # 'w)) x =
3946        (clear_flags, float_plus_infinity (:'t # 'w)))
3947Proof
3948   simp [float_add_def, float_values]
3949QED
3950
3951Theorem float_add_finite_minus_infinity:
3952    !mode x r.
3953       (float_value x = Float r) ==>
3954       (float_add mode x (float_minus_infinity (:'t # 'w)) =
3955        (clear_flags, float_minus_infinity (:'t # 'w)))
3956Proof
3957   simp [float_add_def, float_values]
3958QED
3959
3960Theorem float_add_minus_infinity_finite:
3961    !mode x r.
3962       (float_value x = Float r) ==>
3963       (float_add mode (float_minus_infinity (:'t # 'w)) x =
3964        (clear_flags, float_minus_infinity (:'t # 'w)))
3965Proof
3966   simp [float_add_def, float_values]
3967QED
3968
3969(* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -  *)
3970
3971Theorem float_sub_compute:
3972    (!mode x fp_op.
3973       float_sub mode (float_some_qnan fp_op) x =
3974       (check_for_signalling [x],
3975        float_some_qnan (FP_Sub mode (float_some_qnan fp_op) x))) /\
3976    (!mode x fp_op.
3977       float_sub mode x (float_some_qnan fp_op) =
3978       (check_for_signalling [x],
3979        float_some_qnan (FP_Sub mode x (float_some_qnan fp_op)))) /\
3980    (!mode.
3981       float_sub mode (float_minus_infinity (:'t # 'w))
3982                      (float_minus_infinity (:'t # 'w)) =
3983       (invalidop_flags,
3984        float_some_qnan (FP_Sub mode (float_minus_infinity (:'t # 'w))
3985                                     (float_minus_infinity (:'t # 'w))))) /\
3986    (!mode.
3987       float_sub mode (float_minus_infinity (:'t # 'w))
3988                      (float_plus_infinity (:'t # 'w)) =
3989       (clear_flags, float_minus_infinity (:'t # 'w))) /\
3990    (!mode.
3991       float_sub mode (float_plus_infinity (:'t # 'w))
3992                      (float_plus_infinity (:'t # 'w)) =
3993       (invalidop_flags,
3994        float_some_qnan (FP_Sub mode (float_plus_infinity (:'t # 'w))
3995                                     (float_plus_infinity (:'t # 'w))))) /\
3996    (!mode.
3997       float_sub mode (float_plus_infinity (:'t # 'w))
3998                      (float_minus_infinity (:'t # 'w)) =
3999       (clear_flags, float_plus_infinity (:'t # 'w)))
4000Proof
4001   simp [float_sub_def, float_values, float_components, some_nan_properties,
4002         check_for_signalling_def]
4003   \\ strip_tac
4004   \\ strip_tac
4005   \\ Cases_on `float_value x`
4006   \\ simp [float_is_signalling_def, float_is_nan_def]
4007QED
4008
4009Theorem float_sub_nan:
4010    !mode x y.
4011       (float_value x = NaN) \/ (float_value y = NaN) ==>
4012       (float_sub mode x y =
4013        (check_for_signalling [x; y], float_some_qnan (FP_Sub mode x y)))
4014Proof
4015   NTAC 3 strip_tac
4016   \\ Cases_on `float_value x`
4017   \\ Cases_on `float_value y`
4018   \\ simp [float_sub_def, check_for_signalling_def,
4019            float_is_signalling_def, float_is_nan_def]
4020QED
4021
4022Theorem float_sub_finite:
4023    !mode x y r1 r2.
4024       (float_value x = Float r1) /\ (float_value y = Float r2) ==>
4025       (float_sub mode x y =
4026        float_round_with_flags mode
4027           (if (r1 = 0) /\ (r2 = 0) /\ x.Sign <> y.Sign then
4028              x.Sign = 1w
4029            else mode = roundTowardNegative) (r1 - r2))
4030Proof
4031   simp [float_sub_def]
4032QED
4033
4034Theorem float_sub_finite_plus_infinity:
4035    !mode x r.
4036       (float_value x = Float r) ==>
4037       (float_sub mode x (float_plus_infinity (:'t # 'w)) =
4038        (clear_flags, float_minus_infinity (:'t # 'w)))
4039Proof
4040   simp [float_sub_def, float_values, float_minus_infinity_def]
4041QED
4042
4043Theorem float_sub_plus_infinity_finite:
4044    !mode x r.
4045       (float_value x = Float r) ==>
4046       (float_sub mode (float_plus_infinity (:'t # 'w)) x =
4047        (clear_flags, float_plus_infinity (:'t # 'w)))
4048Proof
4049   simp [float_sub_def, float_values]
4050QED
4051
4052Theorem float_sub_finite_minus_infinity:
4053    !mode x r.
4054       (float_value x = Float r) ==>
4055       (float_sub mode x (float_minus_infinity (:'t # 'w)) =
4056        (clear_flags, float_plus_infinity (:'t # 'w)))
4057Proof
4058   simp [float_sub_def, float_values, float_negate_negate,
4059         float_minus_infinity_def]
4060QED
4061
4062Theorem float_sub_minus_infinity_finite:
4063    !mode x r.
4064       (float_value x = Float r) ==>
4065       (float_sub mode (float_minus_infinity (:'t # 'w)) x =
4066        (clear_flags, float_minus_infinity (:'t # 'w)))
4067Proof
4068   simp [float_sub_def, float_values]
4069QED
4070
4071(* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -  *)
4072
4073Theorem float_mul_compute:
4074    (!mode x fp_op.
4075       float_mul mode (float_some_qnan fp_op) x =
4076       (check_for_signalling [x],
4077        float_some_qnan (FP_Mul mode (float_some_qnan fp_op) x))) /\
4078    (!mode x fp_op.
4079       float_mul mode x (float_some_qnan fp_op) =
4080       (check_for_signalling [x],
4081        float_some_qnan (FP_Mul mode x (float_some_qnan fp_op)))) /\
4082    (!mode.
4083       float_mul mode (float_minus_infinity (:'t # 'w))
4084                      (float_minus_infinity (:'t # 'w)) =
4085       (clear_flags, float_plus_infinity (:'t # 'w))) /\
4086    (!mode.
4087       float_mul mode (float_minus_infinity (:'t # 'w))
4088                      (float_plus_infinity (:'t # 'w)) =
4089       (clear_flags, float_minus_infinity (:'t # 'w))) /\
4090    (!mode.
4091       float_mul mode (float_plus_infinity (:'t # 'w))
4092                      (float_plus_infinity (:'t # 'w)) =
4093       (clear_flags, float_plus_infinity (:'t # 'w))) /\
4094    (!mode.
4095       float_mul mode (float_plus_infinity (:'t # 'w))
4096                      (float_minus_infinity (:'t # 'w)) =
4097       (clear_flags, float_minus_infinity (:'t # 'w)))
4098Proof
4099   simp [float_mul_def, float_values, float_components, some_nan_properties,
4100         check_for_signalling_def]
4101   \\ strip_tac
4102   \\ strip_tac
4103   \\ Cases_on `float_value x`
4104   \\ simp [float_is_signalling_def, float_is_nan_def]
4105QED
4106
4107Theorem float_mul_nan:
4108    !mode x y.
4109       (float_value x = NaN) \/ (float_value y = NaN) ==>
4110       (float_mul mode x y =
4111        (check_for_signalling [x; y], float_some_qnan (FP_Mul mode x y)))
4112Proof
4113   NTAC 3 strip_tac
4114   \\ Cases_on `float_value x`
4115   \\ Cases_on `float_value y`
4116   \\ simp [float_mul_def, check_for_signalling_def,
4117            float_is_signalling_def, float_is_nan_def]
4118QED
4119
4120Theorem float_mul_finite:
4121    !mode x y r1 r2.
4122       (float_value x = Float r1) /\ (float_value y = Float r2) ==>
4123       (float_mul mode x y =
4124        float_round_with_flags mode (x.Sign <> y.Sign) (r1 * r2))
4125Proof
4126   simp [float_mul_def]
4127QED
4128
4129Theorem float_mul_finite_plus_infinity:
4130    !mode x r.
4131       (float_value x = Float r) ==>
4132       (float_mul mode x (float_plus_infinity (:'t # 'w)) =
4133        if r = 0 then
4134           (invalidop_flags,
4135            float_some_qnan (FP_Mul mode x (float_plus_infinity (:'t # 'w))))
4136        else (clear_flags,
4137              if x.Sign = 0w then
4138                float_plus_infinity (:'t # 'w)
4139              else float_minus_infinity (:'t # 'w)))
4140Proof
4141   rw [float_mul_def, float_values]
4142   \\ fs [float_plus_infinity_def]
4143QED
4144
4145Theorem float_mul_plus_infinity_finite:
4146    !mode x r.
4147       (float_value x = Float r) ==>
4148       (float_mul mode (float_plus_infinity (:'t # 'w)) x =
4149        if r = 0 then
4150           (invalidop_flags,
4151            float_some_qnan (FP_Mul mode (float_plus_infinity (:'t # 'w)) x))
4152        else (clear_flags,
4153              if x.Sign = 0w
4154                 then float_plus_infinity (:'t # 'w)
4155              else float_minus_infinity (:'t # 'w)))
4156Proof
4157   rw [float_mul_def, float_values]
4158   \\ fs [float_plus_infinity_def]
4159QED
4160
4161Theorem float_mul_finite_minus_infinity:
4162    !mode x r.
4163       (float_value x = Float r) ==>
4164       (float_mul mode x (float_minus_infinity (:'t # 'w)) =
4165        if r = 0 then
4166           (invalidop_flags,
4167            float_some_qnan (FP_Mul mode x (float_minus_infinity (:'t # 'w))))
4168        else (clear_flags,
4169              if x.Sign = 0w
4170                 then float_minus_infinity (:'t # 'w)
4171              else float_plus_infinity (:'t # 'w)))
4172Proof
4173   rw [float_mul_def, float_values]
4174   \\ fs [float_minus_infinity_def, float_plus_infinity_def, float_negate_def]
4175   \\ metis_tac [sign_inconsistent]
4176QED
4177
4178Theorem float_mul_minus_infinity_finite:
4179    !mode x r.
4180       (float_value x = Float r) ==>
4181       (float_mul mode (float_minus_infinity (:'t # 'w)) x =
4182        if r = 0 then
4183           (invalidop_flags,
4184            float_some_qnan (FP_Mul mode (float_minus_infinity (:'t # 'w)) x))
4185        else (clear_flags,
4186              if x.Sign = 0w
4187                 then float_minus_infinity (:'t # 'w)
4188              else float_plus_infinity (:'t # 'w)))
4189Proof
4190   rw [float_mul_def, float_values]
4191   \\ fs [float_minus_infinity_def, float_plus_infinity_def, float_negate_def]
4192   \\ metis_tac [sign_inconsistent]
4193QED
4194
4195(* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -  *)
4196
4197Theorem float_div_compute:
4198    (!mode x fp_op.
4199       float_div mode (float_some_qnan fp_op) x =
4200       (check_for_signalling [x],
4201        float_some_qnan (FP_Div mode (float_some_qnan fp_op) x))) /\
4202    (!mode x fp_op.
4203       float_div mode x (float_some_qnan fp_op) =
4204       (check_for_signalling [x],
4205        float_some_qnan (FP_Div mode x (float_some_qnan fp_op)))) /\
4206    (!mode.
4207       float_div mode (float_minus_infinity (:'t # 'w))
4208                      (float_minus_infinity (:'t # 'w)) =
4209       (invalidop_flags,
4210        float_some_qnan (FP_Div mode (float_minus_infinity (:'t # 'w))
4211                                     (float_minus_infinity (:'t # 'w))))) /\
4212    (!mode.
4213       float_div mode (float_minus_infinity (:'t # 'w))
4214                      (float_plus_infinity (:'t # 'w)) =
4215       (invalidop_flags,
4216        float_some_qnan (FP_Div mode (float_minus_infinity (:'t # 'w))
4217                                     (float_plus_infinity (:'t # 'w))))) /\
4218    (!mode.
4219       float_div mode (float_plus_infinity (:'t # 'w))
4220                      (float_plus_infinity (:'t # 'w)) =
4221       (invalidop_flags,
4222        float_some_qnan (FP_Div mode (float_plus_infinity (:'t # 'w))
4223                                     (float_plus_infinity (:'t # 'w))))) /\
4224    (!mode.
4225       float_div mode (float_plus_infinity (:'t # 'w))
4226                      (float_minus_infinity (:'t # 'w)) =
4227       (invalidop_flags,
4228        float_some_qnan (FP_Div mode (float_plus_infinity (:'t # 'w))
4229                                     (float_minus_infinity (:'t # 'w)))))
4230Proof
4231   simp [float_div_def, float_values, float_components, some_nan_properties,
4232         check_for_signalling_def]
4233   \\ strip_tac
4234   \\ strip_tac
4235   \\ Cases_on `float_value x`
4236   \\ simp [float_is_signalling_def, float_is_nan_def]
4237QED
4238
4239Theorem float_div_nan:
4240    !mode x y.
4241       (float_value x = NaN) \/ (float_value y = NaN) ==>
4242       (float_div mode x y =
4243        (check_for_signalling [x; y], float_some_qnan (FP_Div mode x y)))
4244Proof
4245   NTAC 3 strip_tac
4246   \\ Cases_on `float_value x`
4247   \\ Cases_on `float_value y`
4248   \\ simp [float_div_def, check_for_signalling_def,
4249            float_is_signalling_def, float_is_nan_def]
4250QED
4251
4252Theorem float_div_finite:
4253    !mode x y r1 r2.
4254       (float_value x = Float r1) /\ (float_value y = Float r2) ==>
4255       (float_div mode x y =
4256        if r2 = 0
4257           then if r1 = 0 then
4258                  (invalidop_flags, float_some_qnan (FP_Div mode x y))
4259                else
4260                  (dividezero_flags,
4261                   if x.Sign = y.Sign then float_plus_infinity (:'t # 'w)
4262                   else float_minus_infinity (:'t # 'w))
4263        else float_round_with_flags mode (x.Sign <> y.Sign) (r1 / r2))
4264Proof
4265   simp [float_div_def]
4266QED
4267
4268Theorem float_div_finite_plus_infinity:
4269    !mode x r.
4270       (float_value x = Float r) ==>
4271       (float_div mode x (float_plus_infinity (:'t # 'w)) =
4272        (clear_flags,
4273         if x.Sign = 0w then float_plus_zero (:'t # 'w)
4274         else float_minus_zero (:'t # 'w)))
4275Proof
4276   rw [float_div_def, float_values]
4277   \\ fs [float_plus_infinity_def]
4278QED
4279
4280Theorem float_div_plus_infinity_finite:
4281    !mode x r.
4282       (float_value x = Float r) ==>
4283       (float_div mode (float_plus_infinity (:'t # 'w)) x =
4284        (clear_flags,
4285         if x.Sign = 0w then float_plus_infinity (:'t # 'w)
4286         else float_minus_infinity (:'t # 'w)))
4287Proof
4288   rw [float_div_def, float_values]
4289   \\ fs [float_plus_infinity_def]
4290QED
4291
4292Theorem float_div_finite_minus_infinity:
4293    !mode x r.
4294       (float_value x = Float r) ==>
4295       (float_div mode x (float_minus_infinity (:'t # 'w)) =
4296        (clear_flags,
4297         if x.Sign = 0w then float_minus_zero (:'t # 'w)
4298         else float_plus_zero (:'t # 'w)))
4299Proof
4300   rw [float_div_def, float_values]
4301   \\ fs [float_minus_infinity_def, float_plus_infinity_def, float_negate_def]
4302   \\ metis_tac [sign_inconsistent]
4303QED
4304
4305Theorem float_div_minus_infinity_finite:
4306    !mode x r.
4307       (float_value x = Float r) ==>
4308       (float_div mode (float_minus_infinity (:'t # 'w)) x =
4309        (clear_flags,
4310         if x.Sign = 0w then float_minus_infinity (:'t # 'w)
4311         else float_plus_infinity (:'t # 'w)))
4312Proof
4313   rw [float_div_def, float_values]
4314   \\ fs [float_minus_infinity_def, float_plus_infinity_def, float_negate_def]
4315   \\ metis_tac [sign_inconsistent]
4316QED
4317
4318(* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -  *)
4319
4320Theorem float_is_nan_impl:
4321  !x. float_is_nan x <=> ~float_equal x x
4322Proof
4323   simp[float_is_nan_def, float_equal_def, float_compare_def]
4324   \\ strip_tac
4325   \\ Cases_on `float_value x`
4326   \\ simp[]
4327QED
4328
4329Theorem float_is_zero_impl:
4330  !x. float_is_zero x <=> float_equal (float_plus_zero (:'w # 't)) x
4331Proof
4332  simp[float_is_zero_def, float_equal_def, float_compare_def, AllCaseEqs(),
4333       SF CONJ_ss] >> qx_gen_tac ‘x’ >>
4334  Cases_on `float_value x` >> simp[EQ_SYM_EQ]
4335QED
4336
4337(* ------------------------------------------------------------------------ *)
4338
4339(* ----------------------------------------------------------------------
4340    operations working over ulps
4341   ---------------------------------------------------------------------- *)
4342
4343val _ = augment_srw_ss [realSimps.RMULRELNORM_ss, realSimps.RMULCANON_ss]
4344
4345Theorem abs_ULP[simp]:
4346  abs (ULP(x,y)) = ULP(x,y)
4347Proof
4348  Cases_on ‘y’ >>
4349  rw[ULP_def, ABS_REFL, real_div, REAL_POW_ADD, REAL_INV_MUL, REAL_LE_MUL,
4350     POW_POS]
4351QED
4352
4353Theorem abs_ulp[simp]:
4354  abs (ulp (:α # β)) = ulp (:α # β)
4355Proof
4356  simp[ulp_def, ULP_def]
4357QED
4358
4359Definition float_ulp_def:
4360  float_ulp (f : (α,β)float) = ULP(f.Exponent, (:α))
4361End
4362
4363Overload "ulpᶠ" = “float_ulp”
4364
4365Theorem float_ulp_negate[simp]:
4366  ulpᶠ (float_negate f) = ulpᶠ f
4367Proof
4368  simp[float_ulp_def, float_components]
4369QED
4370
4371Theorem float_ulp_updating_Significand[simp]:
4372  ulpᶠ (f with Significand := (s:α word)) = ulpᶠ (f:(α,β)float)
4373Proof
4374  simp[float_ulp_def]
4375QED
4376
4377Theorem float_ulp_updating_Sign[simp]:
4378  ulpᶠ (f with Sign := s) = ulpᶠ f
4379Proof
4380  simp[float_ulp_def]
4381QED
4382
4383Theorem ABS_REFL'[local]:
4384  0 ≤ x ⇒ abs x = x
4385Proof
4386  metis_tac[ABS_REFL]
4387QED
4388
4389fun NODP f ths = f (Excl "REAL_ARITH_DP"::ths)
4390
4391val ndps = NODP simp
4392val ndpf = NODP fs
4393val ndpr = NODP rw
4394val ndpg = NODP gs
4395Overload f2r[local] = “float_to_real”
4396
4397
4398val _ = augment_srw_ss [realSimps.REAL_ARITH_ss];
4399Theorem zero_le_two_pow_inv[simp]:
4400  0 ≤ inv (2 pow n)
4401Proof
4402  simp[REAL_LE_LT]
4403QED
4404
4405
4406Theorem abs_f2r_le_float_ulp_mono:
4407  abs (f2r (x:(α,β)float)) ≤ abs (f2r (y:(α,β)float)) ⇒
4408  ulpᶠ x ≤ ulpᶠ y
4409Proof
4410  simp[float_ulp_def, AllCaseEqs(), ULP_def] >> rw[] >>
4411  simp[REAL_LE_RDIV_CANCEL] >> simp[REAL_OF_NUM_POW]
4412  >- (Cases_on ‘y.Exponent’ >> gvs[dimword_def]) >>
4413  gvs[ABS_MUL, ABS_INV, ABS_REFL', REAL_LE_ADD, REAL_LE_MUL,
4414      float_to_real_def, real_div]
4415  >- (Cases_on ‘x.Exponent’ >> gvs[dimword_def] >>
4416      rename [‘x.Exponent = n2w xE’] >>
4417      CCONTR_TAC >> gvs[NOT_LE] >> qpat_x_assum ‘_:real ≤ _’ mp_tac >>
4418      simp[REAL_NOT_LE] >>
4419      irule REAL_LTE_TRANS >> qexists‘2 pow xE * 2 pow precision(:α)’ >>
4420      simp[REAL_LE_MUL] >>
4421      Cases_on ‘y.Significand’ >> gvs[dimword_def, REAL_OF_NUM_POW] >>
4422      irule LESS_TRANS >> qexists ‘2 * 2 ** precision(:α)’>> simp[]) >>
4423  map_every Cases_on [‘x.Exponent’, ‘y.Exponent’] >> gvs[dimword_def] >>
4424  CCONTR_TAC >> gvs[REAL_NOT_LE, NOT_LE] >> rename [‘yE :num < xE’] >>
4425  qpat_x_assum ‘_ : real ≤ _’ mp_tac >> simp[REAL_NOT_LE] >>
4426  irule REAL_LTE_TRANS >> qexists ‘2 pow xE’ >> simp[REAL_LE_MUL] >>
4427  dxrule (iffLR LT_EXISTS) >> rw[REAL_POW_ADD, pow] >>
4428  Cases_on ‘y.Significand’ >> gvs[dimword_def] >>
4429  rename [‘y.Significand = n2w yS’] >>
4430  ‘2 pow yE + &yS * 2 pow yE * inv (2 pow precision(:α)) =
4431   2 pow yE * (1 + &yS * inv (2 pow precision(:α)))’ by simp[] >>
4432  pop_assum SUBST1_TAC >> simp[REAL_POW_ADD] >>
4433  irule REAL_LTE_TRANS >> qexists ‘2’ >>
4434  simp[REAL_OF_NUM_POW, REAL_ARITH “1 + x < 2r ⇔ x < 1”]
4435QED
4436
4437Definition next_hi_def:
4438  next_hi (x:(τ, χ) float) =
4439    if x.Significand <₊ UINT_MAXw
4440    then x with Significand := (x.Significand + 1w)
4441    else <| Sign        := x.Sign
4442          ; Exponent    := x.Exponent + 1w
4443          ; Significand := 0w
4444          |>
4445End
4446
4447Definition next_lo_def:
4448  next_lo (x:(τ, χ) float) =
4449    if 0w <₊ x.Significand
4450    then x with Significand := (x.Significand - 1w)
4451    else <| Sign        := x.Sign
4452          ; Exponent    := x.Exponent - 1w
4453          ; Significand := UINT_MAXw
4454          |>
4455End
4456
4457Theorem next_lo_Sign[simp]:
4458  (next_lo f).Sign = f.Sign
4459Proof
4460  rw[next_lo_def]
4461QED
4462
4463Theorem next_hi_11[simp]:
4464  next_hi f = next_hi g ⇔ f = g
4465Proof
4466  simp[next_hi_def, EQ_IMP_THM, word_T_def, AllCaseEqs(), UINT_MAX_def,
4467       dimword_def] >>
4468  map_every Cases_on [‘f.Significand’, ‘g.Significand’, ‘f.Exponent’,
4469                      ‘g.Exponent’, ‘f.Sign’, ‘g.Sign’] >>
4470  gvs[dimword_def, word_lo_n2w, dimindex_1] >>
4471  simp[float_component_equality, dimword_def, dimindex_1] >>
4472  rw[] >> gvs[word_add_n2w, dimword_def]
4473QED
4474
4475Theorem zero_le_next_hi[simp]:
4476  float_is_finite f ⇒
4477  (0 ≤ f2r (next_hi f) ⇔ 0 ≤ f2r (f:(α,β)float) ∧ f ≠ NEG0)
4478Proof
4479  Cases_on ‘f.Significand’ >> Cases_on ‘f.Exponent’ >>
4480  gvs[dimword_def, float_to_real_def, next_hi_def, word_T_def, UINT_MAX_def,
4481      word_lo_n2w, float_minus_zero, word_add_n2w, float_is_finite_def,
4482      float_value_def] >> rw[] >>
4483  gvs[dimword_def, REAL_MUL_SIGN, GSYM REAL_LE_RNEG, REAL_OF_NUM_POW,
4484      REAL_LE_ADD, REAL_LE_MUL] >>
4485  Cases_on ‘EVEN (w2n f.Sign)’ >> gvs[] >>
4486  Cases_on ‘f.Sign’ >>
4487  gvs[dimword_def, dimindex_1, DECIDE “n < 2n ⇔ n = 0 ∨ n = 1”] >>~-
4488  ([‘f ≠ <| Sign := _; Exponent := _; Significand := _ |> (* g *)’],
4489   strip_tac >> gvs[]) >>
4490  Cases_on ‘n = 0’ >> gvs[] >>
4491  simp[float_component_equality]
4492QED
4493
4494
4495Theorem abs_sign_sub[local]:
4496  abs (sign f * x1 * x2 - sign f * y1 * y2) =
4497  abs (x1 * x2 - y1 * y2)
4498Proof
4499  qspec_then ‘f.Sign’ strip_assume_tac ranged_word_nchotomy >> simp[] >>
4500  ‘x1 * x2 * -1 pow n - y1 * y2 * -1 pow n =
4501   -1 pow n * (x1 * x2 - y1 * y2)’ by simp[] >> pop_assum SUBST1_TAC >>
4502  simp[POW_M1, REAL_ABS_MUL]
4503QED
4504
4505Theorem REAL_SUB'[local]:
4506  n ≤ m ⇒ &(m - n) : real = &m - &n
4507Proof
4508  simp[REAL_SUB] >> rw[]
4509QED
4510
4511Theorem abs_2pow[local,simp]:
4512  abs (2 pow n) = 2 pow n
4513Proof
4514  simp[ABS_REFL]
4515QED
4516
4517Theorem next_lo_difference:
4518  ¬float_is_zero (f:(α,β)float) ∧ float_is_finite f ⇒
4519  abs(float_to_real f - float_to_real (next_lo f)) = ulpᶠ (next_lo f)
4520Proof
4521  rw[next_lo_def, float_to_real_def] >> fs[] >>
4522  rfs[abs_sign_sub, word_lo_n2w] >>
4523  simp[REAL_ABS_MUL, GSYM REAL_SUB_LDISTRIB, real_div,
4524       GSYM REAL_SUB_RDISTRIB, GSYM POW_ABS, ABS_INV, POW_NZ]
4525  >- (simp[float_ulp_def, ULP_def] >>
4526      qspec_then ‘f.Significand’ strip_assume_tac ranged_word_nchotomy >>
4527      fs[word_lo_n2w, GSYM n2w_sub] >> simp[REAL_SUB] >>
4528      simp[real_div, REAL_POW_ADD, REAL_INV_MUL, POW_INV, WORD_LITERAL_ADD])
4529  >- gvs[float_is_zero, WORD_LO_word_0]
4530  >- (Cases_on ‘f.Significand’ >>
4531      Cases_on ‘f.Exponent’ >>
4532      rename [‘f.Significand = n2w s’, ‘f.Exponent = n2w e’] >>
4533      gvs[word_lo_n2w, GSYM n2w_sub, WORD_LITERAL_ADD] >>
4534      ‘e = 1’ by simp[] >> rw[float_ulp_def] >>
4535      qmatch_abbrev_tac ‘abs (2 * (SF * B) - 2 * (Y1 * SF * B * Y2)) = _’ >>
4536      ‘2 * (SF * B) - 2 * (Y1 * SF * B * Y2) = 2 * (SF * B) * (1 - Y1 * Y2)’
4537        by simp[] >> pop_assum SUBST1_TAC >>
4538      simp[REAL_ABS_MUL] >>
4539      simp[Abbr‘B’, Abbr‘SF’, REAL_ABS_MUL, ABS_INV, POW_NZ, GSYM POW_ABS] >>
4540      simp[ULP_def, real_div, REAL_POW_ADD, REAL_INV_MUL] >>
4541      map_every Q.UNABBREV_TAC [‘Y1’, ‘Y2’] >>
4542      simp[dimword_def, UINT_MAX_def, GSYM REAL_OF_NUM_POW,
4543           REAL_SUB', REAL_SUB_LDISTRIB, REAL_SUB_SUB2, w2n_minus1])
4544  >- (‘∀a b c:real. (a + b) - (a + c) = b - c’ by simp[] >> simp[] >>
4545      simp[REAL_ABS_MUL, GSYM REAL_SUB_LDISTRIB, real_div,
4546           GSYM REAL_SUB_RDISTRIB, ABS_INV, POW_NZ] >>
4547      simp[float_ulp_def, ULP_def] >>
4548      qspec_then ‘f.Significand’ strip_assume_tac ranged_word_nchotomy >>
4549      fs[word_lo_n2w, GSYM n2w_sub] >> simp[REAL_SUB] >>
4550      simp[real_div, REAL_POW_ADD, REAL_INV_MUL, WORD_LITERAL_ADD])
4551  >- (Cases_on ‘f.Significand’ >>
4552      Cases_on ‘f.Exponent’ >>
4553      rename [‘f.Significand = n2w s’, ‘f.Exponent = n2w e’] >>
4554      gvs[word_lo_n2w, GSYM n2w_sub, WORD_LITERAL_ADD, dimword_def] >>
4555      ‘1 < e’ by simp[] >>
4556      rw[w2n_minus1, float_ulp_def] >>
4557      fs[word_lo_n2w, GSYM n2w_sub] >> simp[REAL_SUB] >>
4558      simp[GSYM pow_inv_mul_invlt] >>
4559      qmatch_abbrev_tac ‘
4560       abs (SF * PF * BF - 1 / 2 * (SF * PF * BF * GROSS)) =
4561       ULP (n2w (e-1), (:α))
4562      ’ >>
4563      ‘SF * PF * BF - 1 / 2 * SF * PF * BF * GROSS =
4564       SF * PF * BF * (1 - 1 / 2 * GROSS)’ by simp[] >>
4565      pop_assum SUBST1_TAC >> simp[REAL_ABS_MUL] >>
4566      simp[Abbr‘BF’, Abbr‘SF’, Abbr‘PF’, ABS_INV, POW_NZ, GSYM POW_ABS] >>
4567      simp[Abbr‘GROSS’, word_T_def, UINT_MAX_def, ULP_def, dimword_def,
4568           REAL_SUB', GSYM REAL_OF_NUM_POW, REAL_SUB_LDISTRIB,
4569           REAL_LDISTRIB, REAL_POW_ADD] >>
4570      ‘∀x y. 1/2 * y + (1/2 * y - x) = y - x’
4571        by simp[REAL_INV_1OVER, REAL_DOUBLE,
4572                REAL_ARITH “(x:real) + (y - z) = x + y - z”] >>
4573      simp[REAL_SUB_SUB2, REAL_ABS_MUL, ABS_INV, GSYM POW_ABS] >>
4574      simp[GSYM pow_inv_mul_invlt] >> REWRITE_TAC [real_div])
4575QED
4576
4577Theorem next_hilo[simp]:
4578  next_hi (next_lo f) = f
4579Proof
4580  Cases_on ‘f.Significand’ >>
4581  Cases_on ‘f.Exponent’ >>
4582  simp[next_hi_def, next_lo_def, float_to_real_def] >>
4583  gvs[dimword_def, word_lo_n2w] >>
4584  rename [‘f.Significand = n2w fS’] >>
4585  Cases_on ‘0 < fS’ >> simp[WORD_ADD_LEFT_LO2, dimword_def, WORD_LO_word_T] >>
4586  simp[float_component_equality, dimword_def]
4587QED
4588
4589Theorem float_is_finite_Exponent:
4590  float_is_finite f ⇔ f.Exponent ≠ UINT_MAXw
4591Proof
4592  simp[float_is_finite_def, float_value_def] >> rw[]
4593QED
4594
4595Theorem next_lohi[simp]:
4596  next_lo (next_hi f) = f
4597Proof
4598  metis_tac[next_hilo, next_hi_11]
4599QED
4600
4601Theorem next_lo_11[simp]:
4602  next_lo f = next_lo g ⇔ f = g
4603Proof
4604  metis_tac[next_lohi, next_hilo]
4605QED
4606
4607Theorem next_hi_Sign[simp]:
4608  (next_hi f).Sign = f.Sign
4609Proof
4610  simp[next_hi_def] >> rw[]
4611QED
4612
4613Theorem Exponent_monotone:
4614  abs (f2r (f1:(α,β)float)) < abs (f2r (f2:(α,β)float)) ⇒
4615  f1.Exponent ≤₊ f2.Exponent
4616Proof
4617  rw[float_to_real_def] >>
4618  gvs[ABS_MUL, ABS_REFL', REAL_LE_ADD, REAL_LE_MUL]
4619  >- (simp[REAL_NOT_LT] >>
4620      map_every Cases_on [‘f1.Exponent’, ‘f1.Significand’, ‘f2.Significand’] >>
4621      gvs[dimword_def] >> irule REAL_LE_TRANS  >>
4622      rename [‘f1.Exponent = n2w E1’, ‘f1.Significand = n2w S1’,
4623              ‘f2.Significand = n2w S2’] >>
4624      qexists ‘2 pow E1 * 2 pow precision(:α)’ >> simp[] >>
4625      simp[REAL_OF_NUM_POW] >>
4626      irule LE_TRANS >> qexists ‘2 * 2 ** precision(:α)’ >> simp[]) >>
4627  map_every Cases_on [‘f1.Exponent’, ‘f2.Exponent’,
4628                      ‘f1.Significand’, ‘f2.Significand’] >>
4629  rename [‘n2w E1 ≤₊ n2w E2’,
4630          ‘f1.Exponent = n2w E1’, ‘f1.Significand = n2w S1’,
4631          ‘f2.Exponent = n2w E2’, ‘f2.Significand = n2w S2’] >>
4632  gvs[dimword_def, word_ls_n2w] >> CCONTR_TAC >> gvs[NOT_LE] >>
4633  qpat_x_assum ‘2 pow _ * _ < _’ mp_tac >> simp[REAL_NOT_LT] >>
4634  qabbrev_tac‘δ = E1 - E2’ >> ‘0 < δ ∧ E1 = E2 + δ’ by simp[Abbr‘δ’] >>
4635  rw[REAL_POW_ADD] >> irule REAL_LE_TRANS >>
4636  qexists ‘2 * (1 + &S1 / 2 pow precision(:α))’ >>
4637  irule_at Any REAL_LE_RMUL_IMP >>
4638  simp[REAL_OF_NUM_POW, REAL_LE_ADD] >> irule REAL_LE_TRANS >>
4639  qexists ‘2’ >> simp[] >>
4640  simp[REAL_ARITH “1 + x ≤ 2 ⇔ x ≤ 1”]
4641QED
4642
4643Theorem next_hi_discrete:
4644  abs (f2r (f0:(α,β)float)) < abs (f2r (f:(α,β)float)) ∧ float_is_finite f ⇒
4645  abs (f2r (next_hi f0)) ≤ abs (f2r f)
4646Proof
4647  rw[next_hi_def]
4648  >- ((* significand gets one larger *)
4649      ‘f0.Exponent ≤₊ f.Exponent’ by metis_tac[Exponent_monotone] >>
4650      Cases_on ‘f.Exponent = f0.Exponent’
4651      >- (‘f0.Significand <₊ f.Significand’
4652            by (qpat_x_assum ‘abs (f2r f0) < _’ mp_tac >>
4653                simp[float_to_real_def] >> rw[] >>
4654                gvs[ABS_MUL] >>
4655                map_every Cases_on [‘f0.Significand’, ‘f.Significand’] >>
4656                gvs[dimword_def, word_lo_n2w, ABS_REFL', REAL_LE_ADD]) >>
4657          simp[float_to_real_def] >> rw[] >>
4658          simp[ABS_MUL] >>
4659          map_every Cases_on [‘f0.Significand’, ‘f.Significand’] >>
4660          gvs[dimword_def, word_lo_n2w, ABS_REFL', REAL_LE_ADD, word_add_n2w]) >>
4661      fs[WORD_NOT_LOWER, WORD_LOWER_OR_EQ] >>
4662      ‘¬(f.Exponent ≤₊ (f0 with Significand := f0.Significand + 1w).Exponent)’
4663        by simp[WORD_NOT_LOWER_EQUAL] >>
4664      drule_at Concl Exponent_monotone >> simp[])
4665  >- (‘f0.Exponent ≤₊ f.Exponent’ by metis_tac[Exponent_monotone] >>
4666      Cases_on ‘f.Exponent = f0.Exponent’
4667      >- (‘f0.Significand <₊ f.Significand’
4668            by (qpat_x_assum ‘abs (f2r f0) < _’ mp_tac >>
4669                simp[float_to_real_def] >> rw[] >>
4670                gvs[ABS_MUL] >>
4671                map_every Cases_on [‘f0.Significand’, ‘f.Significand’] >>
4672                gvs[dimword_def, word_lo_n2w, ABS_REFL', REAL_LE_ADD,
4673                    w2n_minus1]) >>
4674          metis_tac[WORD_NOT_LOWER, WORD_LOWER_EQ_LOWER_TRANS, WORD_LO_word_T])>>
4675      map_every Cases_on [‘f0.Exponent’, ‘f.Exponent’] >>
4676      gs[word_add_n2w, word_lo_n2w, word_ls_n2w, dimword_def, word_T_def,
4677         UINT_MAX_def, NOT_LESS] >>
4678      rename [‘f.Exponent = n2w e’, ‘e0:num ≤ e’, ‘n2w (e0 + 1)’] >>
4679      ‘e0 + 1 ≤ e’ by simp[] >>
4680      Cases_on ‘e = e0 + 1’
4681      >- (qmatch_abbrev_tac ‘abs (float_to_real bump) ≤ abs _’ >>
4682          ‘bump.Exponent = f.Exponent’ by simp[Abbr‘bump’] >>
4683          simp[float_to_real_def, ABS_MUL, dimword_def, ABS_REFL',
4684               REAL_LE_ADD] >>
4685          simp[Abbr‘bump’]) >>
4686      qmatch_abbrev_tac ‘abs (f2r bump) ≤ _’ >>
4687      ‘¬(f.Exponent ≤₊ bump.Exponent)’
4688        by simp[Abbr‘bump’, word_ls_n2w, dimword_def] >>
4689      drule_at Concl Exponent_monotone >> simp[])
4690QED
4691
4692Theorem next_lo_smaller:
4693  ¬float_is_zero f ⇒ abs (f2r (next_lo f)) < abs (f2r f)
4694Proof
4695  simp[next_lo_def, float_is_zero, float_to_real_def] >> rw[] >>
4696  gvs[WORD_LO_word_0, ABS_MUL, ABS_REFL', REAL_LE_ADD]
4697  >- (Cases_on ‘f.Significand’ >> gvs[dimword_def] >>
4698      REWRITE_TAC[GSYM word_sub_def] >>
4699      simp[GSYM n2w_sub, dimword_def] )
4700  >- (simp[w2n_minus1, dimword_def, REAL_OF_NUM_POW] >>
4701      gvs[WORD_SUM_ZERO])
4702  >- (Cases_on ‘f.Significand’ >> gvs[dimword_def] >>
4703      REWRITE_TAC[GSYM word_sub_def] >>
4704      simp[GSYM n2w_sub, dimword_def]) >>
4705  simp[w2n_minus1, dimword_def] >> gvs[WORD_SUM_ZERO] >>
4706  Cases_on ‘f.Exponent’ >> gvs[dimword_def] >>
4707  REWRITE_TAC[GSYM word_sub_def] >>
4708  simp[GSYM n2w_sub, dimword_def] >>
4709  irule REAL_LTE_TRANS>> qexists ‘2 pow (n - 1) * 2’ >>
4710  simp[REAL_ARITH “1r + x < 2 ⇔ x < 1”] >>
4711  simp[REAL_OF_NUM_POW] >>
4712  simp[GSYM EXP]
4713QED
4714
4715Theorem float_is_zero_next_hi[simp]:
4716  float_is_finite a ⇒ ¬float_is_zero (next_hi a)
4717Proof
4718  rw[next_hi_def, float_is_zero] >> simp[WORD_ADD_EQ_SUB]
4719  >- (‘a.Significand ≠ -1w’ by (strip_tac >> gvs[]) >> simp[]) >>
4720  gvs[float_is_finite_Exponent]
4721QED
4722
4723Theorem next_hi_larger:
4724  float_is_finite f ⇒ abs (f2r f) < abs (f2r (next_hi f))
4725Proof
4726  qspec_then ‘next_hi f’ mp_tac (GEN_ALL next_lo_smaller) >> rw[] >>
4727  gvs[next_lohi] >> first_x_assum irule >> simp[]
4728QED
4729
4730Theorem next_hi_float_negate:
4731  next_hi (float_negate f) = float_negate (next_hi f)
4732Proof
4733  simp[float_negate_def, next_hi_def, float_is_finite_def, float_value_def,
4734       AllCaseEqs()] >> rw[] >> fs[]
4735QED
4736
4737Theorem float_negate :
4738  float_value (float_negate a) = Float r ⇔
4739  float_value (a:(α,β)float) = Float (-r)
4740Proof
4741  rw[float_negate_def, float_value_def, float_to_real_def] >>
4742  Cases_on ‘a.Sign’ >> gvs[dimword_def, DECIDE “n < 2n ⇔ n = 0 ∨ n = 1”]
4743QED
4744
4745Theorem float_is_finite_float_negate[simp]:
4746  float_is_finite (float_negate f) ⇔ float_is_finite f
4747Proof
4748  metis_tac[float_is_finite_thm, float_negate, float_negate_negate]
4749QED
4750
4751Theorem float_is_zero_float_value_EQ0:
4752  float_is_zero f ⇔ (float_value f = Float 0)
4753Proof
4754  simp[float_is_zero_to_real, float_value_def, CaseEq "bool"] >>
4755  ‘float_to_real f = 0 ⇒ f.Exponent ≠ -1w’ suffices_by csimp[GSYM WORD_NEG_1] >>
4756  rpt strip_tac >>
4757  fs[float_to_real_def,word_T_def, UINT_MAX_def, CaseEq "bool",
4758     dimword_def, add_ratr, REAL_OF_NUM_POW]
4759QED
4760
4761Theorem float_is_zero_float_is_finite:
4762  float_is_zero f ⇒ float_is_finite f
4763Proof
4764  simp[float_is_zero_float_value_EQ0, float_is_finite_thm]
4765QED
4766
4767Theorem float_is_zero_float_negate[simp]:
4768  float_is_zero (float_negate f) ⇔ float_is_zero f
4769Proof
4770  simp[float_is_zero_float_value_EQ0, float_negate]
4771QED
4772
4773val _ = augment_srw_ss [realSimps.REAL_ARITH_ss]
4774val _ = diminish_srw_ss [
4775    "word arith", "word ground", "word logic", "word shift",
4776    "word subtract", "words"
4777  ]
4778
4779val _ = augment_srw_ss [
4780    rewrites [w2n_n2w, WORD_AND_CLAUSES, n2w_11, WORD_ADD_0]
4781  ]
4782
4783
4784Theorem float_is_finite_next_hi:
4785  2 ≤ precision(:β) ∧
4786  float_is_finite (f:(α,β)float) ∧ abs (float_to_real f) < largest(:α#β) ⇒
4787  float_is_finite (next_hi f)
4788Proof
4789  qspec_then ‘f.Significand’ (qx_choose_then ‘fS’ strip_assume_tac)
4790             ranged_word_nchotomy >>
4791  qspec_then ‘f.Exponent’ (qx_choose_then ‘fE’ strip_assume_tac)
4792             ranged_word_nchotomy >>
4793  simp[float_is_finite_def] >>
4794  Cases_on ‘float_value f’ >> simp[] >>
4795  ‘float_to_real f = r’ by fs[float_value_def, AllCaseEqs()] >>
4796  simp[next_hi_def] >> Cases_on ‘n2w fS <₊ word_T’ >> simp[] >>
4797  gs[float_value_def, AllCaseEqs(), word_T_def, dimword_def, UINT_MAX_def,
4798     word_lo_n2w, WORD_LO_word_T, word_add_n2w, NOT_LESS] >> strip_tac >>
4799  rw[] >>
4800  rename [‘fE + 1 = 2 ** precision(:β) - 1’] >>
4801  qpat_x_assum ‘abs (f2r f) < largest _’ mp_tac >>
4802  simp[float_to_real_def, dimword_def] >> rw[] >> gvs[] >~
4803  [‘1 = 2 ** _ - 1’]
4804  >- gvs[DECIDE “1 ≤ x ⇒ (1n = x - 1 ⇔ x = 2)”] >>
4805  simp[ABS_MUL, REAL_NOT_LT, largest_is_top, float_top_def,
4806       float_to_real_def, GSYM n2w_sub, word_T_def, UINT_MAX_def,
4807       dimword_def, ABS_INV, GSYM POW_ABS, ABS_REFL', REAL_LE_ADD] >>
4808  ‘fE = 2 ** precision(:β) - 2’ by simp[] >> simp[] >>
4809  ‘fS = 2 ** precision(:α) - 1’ by simp[] >> simp[]
4810QED
4811
4812(*
4813    1 ≤ 2 * 2 pow precision(:α) ⇒
4814    (2 * 2 pow precision (:α) − 1 ≤
4815        2 pow maxExp * (2 * 2 pow precision (:α) − 1) ⇔ ??????) *)
4816
4817
4818Theorem abs_float_bounds:
4819  2 ≤ precision(:β) ∧ float_is_finite f ⇒
4820  f2r (float_abs (f:(α,β)float)) ≤ largest(:α#β)
4821Proof
4822  simp[float_abs_def, float_to_real_def, largest_def, UINT_MAX_def,
4823       dimword_def] >> rw[] >>
4824  qabbrev_tac ‘maxExp = 2n ** precision(:β) - 2’ >>
4825  simp[REAL_SUB_LDISTRIB] >>
4826  Cases_on ‘f.Significand’ >> gvs[dimword_def]
4827  >- (irule REAL_LE_TRANS >> qexists ‘2 * 2 pow precision(:α) - 1’ >>
4828      simp[] >>
4829      ‘2 * (2 pow maxExp * 2 pow precision(:α)) - 2 pow maxExp =
4830       2 pow maxExp * (2 * 2 pow precision(:α) - 1)’
4831        by simp[REAL_SUB_LDISTRIB] >>
4832      pop_assum SUBST1_TAC >>
4833      conj_tac >- simp[REAL_OF_NUM_POW, GSYM realaxTheory.REAL_OF_NUM_SUB] >>
4834      ‘1 ≤ 2 * 2 pow precision(:α)’ by simp[REAL_OF_NUM_POW] >>
4835      ‘2 * 2 pow precision(:α) - 1 = 1 * (2 * 2 pow precision(:α) - 1)’
4836        by simp[] >>
4837      pop_assum (CONV_TAC o LAND_CONV o REWR_CONV) >>
4838      irule REAL_LE_RMUL_IMP >> simp[POW_2_LE1]) >>
4839  Cases_on ‘f.Exponent’ >> gvs[dimword_def] >>
4840  rename [‘f.Significand = n2w fS’, ‘f.Exponent = n2w fE’] >>
4841  simp[RealArith.REAL_ARITH “x * y - y * z:real = y * (x - z)”] >>
4842  irule REAL_LE_TRANS >>
4843  qexists ‘2 pow maxExp * (1 + &fS / 2 pow precision(:α))’ >>
4844  simp[] >> conj_tac
4845  >- (irule REAL_LE_RMUL_IMP >>
4846      ‘fE ≤ maxExp’
4847        by (simp[Abbr‘maxExp’] >>
4848            gvs[float_is_finite_Exponent, word_T_def, UINT_MAX_def,
4849                dimword_def]) >>
4850      simp[REAL_POW_MONO, REAL_LE_ADD]) >>
4851  simp[real_div] >>
4852  simp[RealArith.REAL_ARITH “1 + x ≤ 2 - y ⇔ x + y ≤ 1r”,
4853       RealArith.REAL_ARITH “x * y + y = (x+1) * y:real”] >>
4854  simp[REAL_OF_NUM_POW]
4855QED
4856
4857Theorem float_to_real_float_abs[simp]:
4858  float_to_real (float_abs f) = abs (float_to_real f)
4859Proof
4860  simp[float_abs_def, float_to_real_def] >>
4861  rw[ABS_MUL, ABS_INV, GSYM POW_ABS, REAL_LE_ADD]
4862QED
4863
4864Theorem float_is_finite_float_value:
4865  float_is_finite f ⇒ float_value f = Float (f2r f)
4866Proof
4867  simp[float_value_def, float_is_finite_Exponent]
4868QED
4869
4870Theorem float_bounds:
4871  2 <= precision (:β) ∧ float_value (a:(α,β)float) = Float r ⇒
4872  -largest(:α # β) ≤ r ∧ r ≤ largest (:α # β)
4873Proof
4874  strip_tac >> ‘float_is_finite a’ by simp[float_is_finite_thm] >>
4875  drule_all_then strip_assume_tac abs_float_bounds >>
4876  gvs[float_to_real_float_abs] >>
4877  drule_then assume_tac float_is_finite_float_value >> gvs[]
4878QED
4879
4880Theorem next_hi_diff_lemma[local] =
4881  Q.INST [‘f’ |-> ‘next_hi f’] next_lo_difference |> SRULE [next_lohi]
4882
4883Theorem next_hi_difference:
4884  float_is_finite (f:(α,β)float) ⇒
4885  abs(float_to_real (next_hi f) - float_to_real f) = ulpᶠ f
4886Proof
4887  strip_tac >> Cases_on ‘float_is_finite (next_hi f)’
4888  >- metis_tac[next_hi_diff_lemma, float_is_zero_next_hi] >>
4889  gvs[next_hi_def, float_is_finite_Exponent, word_T_def, UINT_MAX_def,
4890      dimword_def] >>
4891  map_every Cases_on [‘f.Significand’, ‘f.Exponent’] >>
4892  gvs[dimword_def, word_lo_n2w] >> rw[] >> gvs[dimword_def, word_add_n2w] >>
4893  rename [‘f.Significand = n2w fS’, ‘f.Exponent = n2w fE’] >>
4894  ‘fE = 2 ** precision(:β) - 2’ by simp[] >> gvs[] >>
4895  simp[float_to_real_def, dimword_def] >>
4896  rw[] >~
4897  [‘precision(:β) ≤ 1’]
4898  >- (‘precision (:β) ≠ 0’ by simp[] >>
4899      ‘precision(:β) = 1’ by simp[] >> gvs[] >>
4900      simp[GSYM REAL_SUB_LDISTRIB, ABS_MUL] >>
4901      simp[float_ulp_def, ULP_def, REAL_POW_ADD] >>
4902      qabbrev_tac ‘B = 2 pow bias(:β)’ >>
4903      qabbrev_tac ‘AP = 2 pow precision(:α)’ >>
4904      qabbrev_tac ‘BP = 2n ** precision (:β)’ >>
4905      ‘B = abs B’ by simp[Abbr‘B’] >> pop_assum SUBST1_TAC >>
4906      REWRITE_TAC[GSYM ABS_MUL] >> simp[REAL_SUB_RDISTRIB] >>
4907      simp[Abbr‘B’] >>
4908      simp[REAL_ARITH “s - a * b * s = s * (1 - a * b)”] >>
4909      simp[ABS_MUL] >>
4910      ‘0 ≤ 1 - AP⁻¹ * &fS’
4911        by (simp[REAL_ARITH “(0r ≤ x - y ⇔ y ≤ x)”] >>
4912            ‘0 < AP’ by simp[Abbr‘AP’] >> simp[] >>
4913            simp[Abbr‘AP’, REAL_OF_NUM_POW]) >>
4914      simp[ABS_REFL'] >> simp[REAL_SUB_LDISTRIB, REAL_LDISTRIB] >>
4915      simp[Abbr‘AP’, REAL_OF_NUM_POW, REAL_SUB]) >>
4916  simp[float_ulp_def, ULP_def, dimword_def, REAL_POW_ADD] >>
4917  qabbrev_tac ‘B = 2 pow bias(:β)’ >>
4918  qabbrev_tac ‘AP = 2 pow precision(:α)’ >>
4919  qabbrev_tac ‘BP = 2n ** precision (:β)’ >>
4920  ‘B = abs B’ by simp[Abbr‘B’] >> pop_assum SUBST1_TAC >>
4921  REWRITE_TAC[GSYM ABS_MUL] >> simp[REAL_SUB_RDISTRIB] >>
4922  simp[Abbr‘B’] >>
4923  simp[REAL_ARITH “s * x - s * a * b = s * (x - a * b)”] >>
4924  simp[ABS_MUL] >>
4925  ‘2 pow (BP - 1) = 2 * 2 pow (BP - 2)’
4926    by (qpat_x_assum ‘BP - 2 + 1 = _’ (SUBST1_TAC o SYM)>>
4927        simp[REAL_POW_ADD]) >>
4928  pop_assum SUBST1_TAC >>
4929  simp[REAL_ARITH “x * b - b * y = b * (x - y):real”] >>
4930  simp[ABS_MUL] >> ‘0 < AP’ by simp[Abbr‘AP’] >>
4931  ‘0 ≤ 2 - (1 + &fS / AP)’
4932    by (simp[REAL_ARITH “(1 + x ≤ 2r ⇔ x ≤ 1) ∧ (0r ≤ x - y ⇔ y ≤ x)”] >>
4933        simp[Abbr‘AP’, REAL_OF_NUM_POW]) >>
4934  simp[ABS_REFL'] >> simp[REAL_SUB_LDISTRIB, REAL_LDISTRIB] >>
4935  ‘2 * AP - (AP + &fS) = AP - &fS’ by simp[] >> simp[] >>
4936  simp[Abbr‘AP’, REAL_OF_NUM_POW, REAL_SUB]
4937QED
4938
4939Theorem next_hi_idem[simp]:
4940  next_hi f ≠ f
4941Proof
4942  rw[next_hi_def, float_component_equality, WORD_ADD_RID_UNIQ]
4943QED
4944
4945Theorem next_lo_idem[simp]:
4946  next_lo f ≠ f
4947Proof
4948  metis_tac[next_hilo, next_hi_idem]
4949QED
4950
4951Theorem sign_float_abs[simp]:
4952  sign (float_abs f) = 1
4953Proof
4954  simp[float_abs_def]
4955QED
4956
4957Theorem float_abs_Exponent[simp]:
4958  (float_abs f).Exponent = f.Exponent
4959Proof
4960  simp[float_abs_def]
4961QED
4962
4963Theorem float_abs_Significand[simp]:
4964  (float_abs f).Significand = f.Significand
4965Proof
4966  simp[float_abs_def]
4967QED
4968
4969Theorem word_1comp_11[simp]:
4970  word_1comp w1 = word_1comp w2 <=> w1 = w2
4971Proof
4972  simp[WORD_NOT, WORD_LCANCEL_SUB, WORD_EQ_NEG]
4973QED
4974
4975Theorem float_negate_11[simp]:
4976  float_negate f1 = float_negate f2 <=> f1 = f2
4977Proof
4978  simp[float_negate_def, float_component_equality]
4979QED
4980
4981Theorem float_is_finite_next_lo:
4982  float_is_finite f ∧ ¬float_is_zero f ⇒ float_is_finite (next_lo f)
4983Proof
4984  map_every Cases_on [‘f.Exponent’, ‘f.Significand’] >>
4985  gvs[next_lo_def, float_is_finite_Exponent, float_is_zero, word_T_def,
4986      dimword_def, UINT_MAX_def, word_lo_n2w, GSYM n2w_sub] >> rw[] >>
4987  simp[dimword_def]
4988QED
4989
4990Theorem float_ulp_abs[simp]:
4991  ulpᶠ (float_abs a) = ulpᶠ a
4992Proof
4993  simp[float_ulp_def, float_abs_def]
4994QED
4995
4996Theorem float_is_finite_float_abs[simp]:
4997  float_is_finite (float_abs f) ⇔ float_is_finite f
4998Proof
4999  simp[float_abs_def, float_is_finite_Exponent]
5000QED
5001
5002Theorem next_hi_float_abs:
5003  next_hi (float_abs f) = float_abs (next_hi f)
5004Proof
5005  rw[next_hi_def, float_abs_def, AllCaseEqs(), float_component_equality]
5006QED
5007
5008Theorem float_ulp_EQ0[simp]:
5009  ulpᶠ f ≠ 0
5010Proof
5011  simp[float_ulp_def]
5012QED
5013
5014Theorem float_ulp_GT0[simp]:
5015  0 < ulpᶠ f
5016Proof
5017  simp[float_ulp_def]
5018QED
5019
5020Theorem abs_float_ulp[simp]:
5021  abs (ulpᶠ f) = ulpᶠ f
5022Proof
5023  simp[REAL_LE_LT]
5024QED