1(* ------------------------------------------------------------------------ 2 Some basic properties of IEEE-754 (base 2) floating-point arithmetic 3 ------------------------------------------------------------------------ *) 4 5open HolKernel boolLib bossLib 6open binary_ieeeTheory realTheory wordsLib realLib 7 8val () = new_theory "lift_ieee" 9 10val () = Parse.temp_overload_on ("bias", ``words$INT_MAX``) 11 12(* ------------------------------------------------------------------------ *) 13 14val error_def = Define` 15 error (:'t # 'w) x = 16 float_to_real (round roundTiesToEven x : ('t, 'w) float) - x` 17 18val normalizes_def = Define` 19 normalizes (:'t # 'w) x = 20 1 < bias (:'w) /\ 21 inv (2 pow (bias (:'w) - 1)) <= abs x /\ abs x < threshold (:'t # 'w)` 22 23(* ------------------------------------------------------------------------ 24 Lifting comparison operations 25 ------------------------------------------------------------------------ *) 26 27val float_lt = Q.store_thm ("float_lt", 28 `!x y. float_is_finite x /\ float_is_finite y ==> 29 (float_less_than x y = float_to_real x < float_to_real y)`, 30 rw [float_less_than_def, float_compare_def, float_is_finite_def, 31 float_value_def] 32 \\ rw [] 33 ); 34 35val float_le = Q.store_thm ("float_le", 36 `!x y. float_is_finite x /\ float_is_finite y ==> 37 (float_less_equal x y = float_to_real x <= float_to_real y)`, 38 rw [float_less_equal_def, float_compare_def, float_is_finite_def, 39 float_value_def] 40 \\ rw [realTheory.REAL_LT_IMP_LE, 41 realLib.REAL_ARITH ``~(a < b : real) /\ a <> b ==> ~(a <= b)``] 42 ); 43 44val float_gt = Q.store_thm ("float_gt", 45 `!x y. float_is_finite x /\ float_is_finite y ==> 46 (float_greater_than x y = float_to_real x > float_to_real y)`, 47 rw [float_greater_than_def, float_compare_def, float_is_finite_def, 48 float_value_def] 49 \\ rw [realLib.REAL_ARITH ``a < b : real ==> ~(a > b)``, 50 realLib.REAL_ARITH ``~(a < b : real) /\ a <> b ==> a > b``, 51 realLib.REAL_ARITH ``~(a > a : real)``] 52 ); 53 54val float_ge = Q.store_thm ("float_ge", 55 `!x y. float_is_finite x /\ float_is_finite y ==> 56 (float_greater_equal x y = float_to_real x >= float_to_real y)`, 57 rw [float_greater_equal_def, float_compare_def, float_is_finite_def, 58 float_value_def] 59 \\ rw [realLib.REAL_ARITH ``a < b : real ==> ~(a >= b)``, 60 realLib.REAL_ARITH ``~(a < b : real) /\ a <> b ==> a >= b``, 61 realLib.REAL_ARITH ``a >= a : real``] 62 ); 63 64val float_eq = Q.store_thm ("float_eq", 65 `!x y. float_is_finite x /\ float_is_finite y ==> 66 (float_equal x y = (float_to_real x = float_to_real y))`, 67 rw [float_equal_def, float_compare_def, float_is_finite_def, 68 float_value_def] 69 \\ rw [realLib.REAL_ARITH ``a < b : real ==> a <> b``] 70 ); 71 72val float_eq_refl = Q.store_thm ("float_eq_refl", 73 `!x. float_equal x x = ~float_is_nan x`, 74 rw [float_equal_def, float_is_nan_def, float_compare_def, float_is_finite_def, 75 float_value_def] 76 ); 77 78 79(* ------------------------------------------------------------------------ 80 Closest 81 ------------------------------------------------------------------------ *) 82 83val is_closest_exits = Q.store_thm("is_closest_exits", 84 `!x s. FINITE s ==> s <> EMPTY ==> ?a. is_closest s x a`, 85 strip_tac 86 \\ HO_MATCH_MP_TAC pred_setTheory.FINITE_INDUCT 87 \\ simp [pred_setTheory.NOT_INSERT_EMPTY] 88 \\ strip_tac 89 \\ Cases_on `s = EMPTY` 90 >- simp [is_closest_def] 91 \\ Cases_on `FINITE s` 92 \\ rw [] 93 \\ Cases_on `abs (float_to_real a - x) <= abs (float_to_real e - x)` 94 \\ fs [is_closest_def] 95 >- (qexists_tac `a` \\ rw [] \\ simp []) 96 \\ qexists_tac `e` 97 \\ rw [] 98 >- simp [] 99 \\ qpat_x_assum `!b. P` (qspec_then `b` mp_tac) 100 \\ qpat_x_assum `~(abs (float_to_real a - x) <= 101 abs (float_to_real e - x))` mp_tac 102 \\ simp [] 103 \\ REAL_ARITH_TAC 104 ); 105 106val closest_is_everything = Q.store_thm("closest_is_everything", 107 `!p s x. FINITE s ==> s <> EMPTY ==> 108 is_closest s x (closest_such p s x) /\ 109 ((?b. is_closest s x b /\ p b) ==> p (closest_such p s x))`, 110 rw [closest_such_def] 111 \\ SELECT_ELIM_TAC 112 \\ metis_tac [is_closest_exits] 113 ); 114 115val closest_in_set = Q.store_thm("closest_in_set", 116 `!p s x. FINITE s ==> s <> EMPTY ==> closest_such p s x IN s`, 117 metis_tac [closest_is_everything, is_closest_def] 118 ); 119 120val closest_is_closest = Q.store_thm("closest_is_closest", 121 `!p s x. FINITE s ==> s <> EMPTY ==> is_closest s x (closest_such p s x)`, 122 metis_tac [closest_is_everything] 123 ); 124 125(* ------------------------------------------------------------------------ 126 127 ------------------------------------------------------------------------ *) 128 129val finite_word_cross = Q.prove( 130 `FINITE (univ (:word1) CROSS univ (:'t word) CROSS univ (:'w word))`, 131 metis_tac [pred_setTheory.FINITE_CROSS, wordsTheory.WORD_FINITE] 132 ); 133 134val inj_float_to_tuple = Q.prove( 135 `INJ (\x. ((x.Sign, x.Significand), x.Exponent)) 136 (univ (:('t, 'w) float)) 137 (univ (:word1) CROSS univ (:'t word) CROSS univ (:'w word))`, 138 rw [pred_setTheory.INJ_DEF, float_component_equality] 139 ); 140 141val float_finite = Q.store_thm("float_finite", 142 `FINITE (univ (:('t, 'w) float))`, 143 metis_tac [pred_setTheory.FINITE_INJ, inj_float_to_tuple, finite_word_cross] 144 ); 145 146val is_finite_finite = Q.store_thm("is_finite_finite", 147 `FINITE {a | float_is_finite a}`, 148 metis_tac [pred_setTheory.SUBSET_FINITE, float_finite, 149 pred_setTheory.SUBSET_UNIV] 150 ); 151 152val is_finite_nonempty = Q.store_thm("is_finite_nonempty", 153 `{a | float_is_finite a} <> EMPTY`, 154 rw [pred_setTheory.EXTENSION] 155 \\ qexists_tac `float_plus_zero (:'a # 'b)` 156 \\ simp [binary_ieeeTheory.zero_properties] 157 ); 158 159val is_finite_closest = Q.store_thm("is_finite_closest", 160 `!p x. float_is_finite (closest_such p {a | float_is_finite a} x)`, 161 rpt strip_tac 162 \\ `closest_such p {a | float_is_finite a} x IN {a | float_is_finite a}` 163 by metis_tac [closest_in_set, is_finite_finite, is_finite_nonempty] 164 \\ fs [] 165 ); 166 167(* ------------------------------------------------------------------------ 168 169 ------------------------------------------------------------------------ *) 170 171val REAL_ABS_INV = Q.prove( 172 `!x. abs (inv x) = inv (abs x)`, 173 gen_tac 174 \\ Cases_on `x = 0r` 175 \\ simp [REAL_INV_0, REAL_ABS_0, ABS_INV] 176 ); 177 178val REAL_ABS_DIV = Q.prove( 179 `!x y. abs (x / y) = abs x / abs y`, 180 rewrite_tac [real_div, REAL_ABS_INV, REAL_ABS_MUL]) 181 182val REAL_LE_LDIV2 = Q.prove( 183 `!x y z. 0r < z ==> (x / z <= y / z <=> x <= y)`, 184 simp [REAL_LE_LDIV_EQ, REAL_DIV_RMUL, REAL_POS_NZ] 185 ); 186 187val REAL_POW_LE_1 = Q.prove( 188 `!n x. 1r <= x ==> 1 <= x pow n`, 189 Induct 190 \\ rw [pow] 191 \\ Ho_Rewrite.GEN_REWRITE_TAC LAND_CONV [GSYM REAL_MUL_LID] 192 \\ match_mp_tac REAL_LE_MUL2 193 \\ simp [] 194 ); 195 196val REAL_POW_MONO = Q.prove( 197 `!m n x. 1r <= x /\ m <= n ==> x pow m <= x pow n`, 198 rw [arithmeticTheory.LESS_EQ_EXISTS] 199 \\ simp [REAL_POW_ADD] 200 \\ Ho_Rewrite.GEN_REWRITE_TAC LAND_CONV [GSYM REAL_MUL_RID] 201 \\ metis_tac [REAL_LE_LMUL_IMP, REAL_POW_LE_1, POW_POS, REAL_LE_TRANS, 202 REAL_LE_01] 203 ); 204 205val exponent_le = Q.prove( 206 `!e : 'a word. e <> -1w ==> (w2n e <= UINT_MAX (:'a) - 1)`, 207 simp_tac std_ss 208 [wordsTheory.WORD_NEG_1, wordsTheory.UINT_MAX_def, wordsTheory.word_T_def] 209 \\ Cases 210 \\ simp [] 211 ); 212 213val float_to_real_finite = Q.store_thm("float_to_real_finite", 214 `!x : ('t, 'w) float. 215 float_is_finite x ==> (abs (float_to_real x) <= largest (:'t # 'w))`, 216 Cases 217 \\ rename [`float s e f`] 218 \\ `float s e f = <| Sign := s; Exponent := e; Significand := f |>` 219 by simp [float_component_equality] 220 \\ `0 <= 1 + &w2n f / 2 pow dimindex (:'t)` 221 by (Cases_on `w2n f = 0` 222 \\ simp [REAL_ARITH ``0r < x ==> 0 <= 1 + x``, REAL_LT_RDIV_0]) 223 \\ `abs (1 + &w2n f / 2 pow dimindex (:'t)) = 224 1 + &w2n f / 2 pow dimindex (:'t)` 225 by simp [ABS_REFL] 226 \\ rw [float_tests, largest, float_to_real_def, GSYM POW_ABS, 227 REAL_ABS_MUL, REAL_ABS_DIV, ABS_NEG, ABS_N, POW_ONE, mult_rat, 228 GSYM REAL_OF_NUM_POW, REAL_LE_LDIV2, wordsTheory.dimword_def, 229 ONCE_REWRITE_RULE [REAL_MUL_COMM] mult_ratr] 230 >- ( 231 simp [largest_def, REAL_LE_RDIV_EQ, mult_ratl] 232 \\ CONV_TAC (PATH_CONV "lrr" (ONCE_REWRITE_CONV [REAL_MUL_COMM])) 233 \\ simp [REAL_DIV_RMUL_CANCEL, REAL_LE_LDIV_EQ, 234 REAL_ARITH ``a * (b - c) * d = a * (b * d - c * d) : real``, 235 GSYM REAL_OF_NUM_POW, REAL_DIV_RMUL] 236 \\ simp [REAL_SUB_LDISTRIB, REAL_LE_SUB_LADD] 237 \\ qabbrev_tac `u = UINT_MAX (:'w) - 1` 238 \\ `2 pow u * (2 * 2 pow dimindex (:'t)) = 239 2 pow u * (2 pow dimindex (:'t) - 1) + 240 2 pow u * (2 pow dimindex (:'t) - 1) + 241 2 * 2 pow u` 242 by (simp [REAL_SUB_LDISTRIB, 243 realLib.REAL_ARITH ``a * (2r * b) = a * b + a * b``, 244 realLib.REAL_ARITH ``2 * 2 pow x = 2 pow x + 2 pow x``] 245 \\ REAL_ARITH_TAC) 246 \\ pop_assum SUBST1_TAC 247 \\ once_rewrite_tac [GSYM REAL_MUL] 248 \\ simp_tac std_ss [realLib.REAL_ARITH ``2r * a = a + a``] 249 \\ match_mp_tac realTheory.REAL_LE_ADD2 250 \\ conj_tac 251 >- ( 252 match_mp_tac realTheory.REAL_LE_ADD2 253 \\ simp [] 254 \\ match_mp_tac 255 (realTheory.REAL_LE_MUL2 256 |> Q.SPEC `1` 257 |> SIMP_RULE (srw_ss()) []) 258 \\ simp 259 [REAL_OF_NUM_POW, REAL_SUB, GSYM wordsTheory.dimword_def, 260 DECIDE ``a <= 1n = ~(1 < a)``, wordsTheory.w2n_lt, 261 DECIDE ``a < n /\ 1n < n ==> a <= n - 1``, 262 REAL_POW_LE_1 263 |> Q.SPECL [`n`, `2`] 264 |> SIMP_RULE (srw_ss()) []] 265 ) 266 \\ simp [REAL_ARITH ``0r < a ==> a <= a + a``] 267 ) 268 \\ match_mp_tac REAL_LE_MUL2 269 \\ simp [REAL_POW_MONO, exponent_le, REAL_LE_SUB_LADD] 270 \\ simp_tac std_ss [GSYM REAL_ADD_ASSOC, REAL_DIV_ADD] 271 \\ match_mp_tac (realLib.REAL_ARITH ``(x:real) <= 1 ==> (1 + x <= 2)``) 272 \\ simp [REAL_LE_LDIV_EQ, REAL_OF_NUM_POW, GSYM wordsTheory.dimword_def, 273 wordsTheory.w2n_lt, DECIDE ``a < n ==> a + 1n <= n``] 274 ); 275 276val float_to_real_threshold = Q.store_thm("float_to_real_threshold", 277 `!x : ('t, 'w) float. 278 float_is_finite x ==> (abs (float_to_real x) < threshold (:'t # 'w))`, 279 metis_tac [REAL_LET_TRANS, float_to_real_finite, largest_lt_threshold] 280 ); 281 282(* ------------------------------------------------------------------------ 283 Lifting up of rounding to nearest 284 ------------------------------------------------------------------------ *) 285 286val bound_at_worst_lemma = Q.prove( 287 `!a : ('t, 'w) float x. 288 abs x < threshold (:'t # 'w) /\ float_is_finite a ==> 289 abs (float_to_real (round roundTiesToEven x : ('t, 'w) float) - x) <= 290 abs (float_to_real a - x)`, 291 rw [round_def, REAL_ARITH ``abs x < y = ~(x <= ~y) /\ ~(x >= y)``] 292 \\ match_mp_tac 293 (is_finite_finite 294 |> MATCH_MP closest_is_closest 295 |> Q.SPECL [`\a. ~word_lsb a.Significand`, `x`] 296 |> REWRITE_RULE [is_finite_nonempty, is_closest_def, 297 pred_setTheory.GSPEC_ETA] 298 |> CONJUNCT2) 299 \\ simp [pred_setTheory.SPECIFICATION] 300 ); 301 302val error_at_worst_lemma = Q.store_thm("error_at_worst_lemma", 303 `!a : ('t, 'w) float x. 304 abs x < threshold (:'t # 'w) /\ float_is_finite a ==> 305 abs (error (:'t # 'w) x) <= abs (float_to_real a - x)`, 306 simp [error_def, bound_at_worst_lemma]) 307 308val error_is_zero = Q.store_thm("error_is_zero", 309 `!a : ('t, 'w) float x. 310 float_is_finite a /\ (float_to_real a = x) ==> (error (:'t # 'w) x = 0)`, 311 rw [] 312 \\ match_mp_tac 313 (error_at_worst_lemma 314 |> Q.SPECL [`a : ('t, 'w) float`, `float_to_real (a : ('t, 'w) float)`] 315 |> SIMP_RULE (srw_ss()) 316 [REAL_ABS_0, REAL_ARITH ``abs x <= 0 = (x = 0r)``]) 317 \\ simp [float_to_real_threshold] 318 ); 319 320(* ------------------------------------------------------------------------ *) 321 322val lem = Q.prove( 323 `!a b. 0 < b /\ b < a ==> 1 < a / b : real`, 324 simp [realTheory.REAL_LT_RDIV_EQ] 325 ); 326 327val lem2 = Q.prove( 328 `!n x. 0r < n /\ n <= n * x ==> 1 <= x`, 329 rpt strip_tac 330 \\ spose_not_then assume_tac 331 \\ qpat_x_assum `n <= n * x` mp_tac 332 \\ fs [realTheory.REAL_NOT_LE, 333 GSYM (ONCE_REWRITE_RULE [REAL_MUL_COMM] realTheory.REAL_LT_RDIV_EQ), 334 realTheory.REAL_DIV_REFL, realTheory.REAL_POS_NZ] 335 ); 336 337val error_bound_lemma1 = Q.prove( 338 `!fracw x. 339 0r <= x /\ x < 1 /\ 0 < fracw ==> 340 ?n. n < 2n EXP fracw /\ &n / 2 pow fracw <= x /\ 341 x < &(SUC n) / 2 pow fracw`, 342 rpt strip_tac 343 \\ qspec_then `\n. &n / 2 pow fracw <= x` mp_tac 344 arithmeticTheory.EXISTS_GREATEST 345 \\ simp [] 346 \\ Lib.W (Lib.C SUBGOAL_THEN (fn th => rewrite_tac [th]) o lhs o lhand o snd) 347 >- (conj_tac 348 >- (qexists_tac `0n` \\ simp []) 349 \\ qexists_tac `2 ** fracw` 350 \\ rw [REAL_LE_LDIV_EQ] 351 \\ fs [realTheory.REAL_NOT_LE, realTheory.real_gt, 352 GSYM realTheory.REAL_LT_RDIV_EQ] 353 \\ `1 < &y / 2 pow fracw` 354 by (match_mp_tac lem \\ simp [realTheory.REAL_OF_NUM_POW]) 355 \\ metis_tac [realTheory.REAL_LT_TRANS] 356 ) 357 \\ disch_then (Q.X_CHOOSE_THEN `n` strip_assume_tac) 358 \\ pop_assum (qspec_then `SUC n` assume_tac) 359 \\ qexists_tac `n` 360 \\ fs [REAL_NOT_LE] 361 \\ fs [REAL_LE_LDIV_EQ] 362 \\ spose_not_then assume_tac 363 \\ `2 pow fracw <= &n` by simp [realTheory.REAL_OF_NUM_POW] 364 \\ `2 pow fracw <= x * 2 pow fracw` 365 by imp_res_tac realTheory.REAL_LE_TRANS 366 \\ metis_tac [binary_ieeeTheory.zero_lt_twopow, REAL_MUL_COMM, lem2, 367 realTheory.real_lt] 368 ); 369 370(* ------------------------------------------------------------------------ *) 371 372val error_bound_lemma2 = Q.prove( 373 `!fracw x. 374 0r <= x /\ x < 1 /\ 0 < fracw ==> 375 ?n. n <= 2 EXP fracw /\ 376 abs (x - &n / 2 pow fracw) <= inv (2 pow (fracw + 1))`, 377 ntac 2 gen_tac 378 \\ disch_then 379 (fn th => Q.X_CHOOSE_THEN `n` (CONJUNCTS_THEN2 ASSUME_TAC MP_TAC) 380 (MATCH_MP error_bound_lemma1 th) 381 \\ strip_assume_tac th) 382 \\ disch_then (mp_tac o Q.SPEC `inv (2 pow (fracw + 1))` o MATCH_MP 383 (REAL_ARITH ``!a:real b x d. a <= x /\ x < b ==> b <= a + 2 * d ==> 384 abs (x - a) <= d \/ abs (x - b) <= d``)) 385 \\ Lib.W (Lib.C SUBGOAL_THEN 386 (fn th => rewrite_tac [th]) o lhand o lhand o snd) 387 >- simp [realTheory.REAL_LE_LDIV_EQ, realTheory.REAL_RDISTRIB, 388 realTheory.REAL_DIV_RMUL, realTheory.REAL_MUL_LINV, 389 realTheory.REAL_POW_ADD, 390 REAL_ARITH ``a * inv (b * a) * b = inv (b * a) * (b * a)``] 391 \\ rw [] 392 >- (qexists_tac `n` \\ fs []) 393 \\ qexists_tac `SUC n` 394 \\ fs [] 395 ); 396 397(* ------------------------------------------------------------------------ *) 398 399val error_bound_lemma3 = Q.prove( 400 `!fracw x. 401 1r <= x /\ x < 2 /\ 0 < fracw ==> 402 ?n. n <= 2 EXP fracw /\ 403 abs ((1 + &n / 2 pow fracw) - x) <= inv (2 pow (fracw + 1))`, 404 rpt strip_tac 405 \\ Q.SUBGOAL_THEN `0r <= x - 1 /\ x - 1 < 1 /\ 0 < fracw` 406 (assume_tac o MATCH_MP error_bound_lemma2) 407 >- (simp [] 408 \\ pop_assum kall_tac 409 \\ ntac 2 (POP_ASSUM mp_tac) 410 \\ REAL_ARITH_TAC 411 ) 412 \\ metis_tac 413 [ABS_NEG, REAL_NEG_SUB, REAL_ARITH ``a - (b - c) = (c + a:real) - b``] 414 ); 415 416(* ------------------------------------------------------------------------ *) 417 418val two_pow_over_pre = Q.prove( 419 `!n. 0 < n ==> (2 pow n / 2 pow (n - 1) = 2)`, 420 rpt strip_tac 421 \\ imp_res_tac arithmeticTheory.LESS_ADD_1 422 \\ fs [POW_ADD, 423 realTheory.REAL_DIV_LMUL_CANCEL 424 |> Q.SPECL [`2 pow n`, `2`, `1`] 425 |> SIMP_RULE (srw_ss()) []] 426 ); 427 428val error_bound_lemma4 = Q.prove( 429 `!x. 1r <= x /\ x < 2 /\ 1 < dimindex (:'w) ==> 430 ?e f. 431 abs (float_to_real <| Sign := 0w; 432 Exponent := e : 'w word; 433 Significand := f : 't word |> - x) <= 434 inv (2 pow (dimindex (:'t) + 1)) /\ 435 ((e = n2w (bias (:'w))) \/ (e = n2w (INT_MIN (:'w))) /\ (f = 0w))`, 436 gen_tac 437 \\ DISCH_TAC 438 \\ Q.SUBGOAL_THEN `1 <= x /\ x < 2 /\ 0 < dimindex (:'t)` assume_tac 439 >- simp [] 440 \\ first_assum 441 (Q.X_CHOOSE_THEN `n` 442 (MP_TAC o REWRITE_RULE [arithmeticTheory.LESS_OR_EQ]) o 443 MATCH_MP error_bound_lemma3) 444 \\ strip_tac 445 >- (qexists_tac `n2w (bias (:'w))` 446 \\ qexists_tac `n2w n` 447 \\ fs [float_to_real_def, wordsTheory.INT_MAX_LT_DIMWORD, 448 GSYM wordsTheory.dimword_def, wordsTheory.ZERO_LT_INT_MAX, 449 realTheory.REAL_DIV_REFL, DECIDE ``0 < x ==> x <> 0n``] 450 ) 451 \\ qexists_tac `n2w (INT_MIN (:'w))` 452 \\ qexists_tac `0w` 453 \\ rfs [float_to_real_def, GSYM realTheory.REAL_OF_NUM_POW, two_pow_over_pre, 454 realTheory.REAL_DIV_REFL, wordsTheory.INT_MAX_def, 455 wordsTheory.INT_MIN_LT_DIMWORD, DECIDE ``0 < x ==> x <> 0n``] 456 ); 457 458(* ------------------------------------------------------------------------ *) 459 460val error_bound_lemma5 = Q.prove( 461 `!x. 1r <= abs x /\ abs x < 2 /\ 1 < dimindex (:'w) ==> 462 ?s e f. 463 abs (float_to_real <| Sign := s; 464 Exponent := e : 'w word; 465 Significand := f : 't word |> - x) <= 466 inv (2 pow (dimindex (:'t) + 1)) /\ 467 ((e = n2w (bias (:'w))) \/ (e = n2w (INT_MIN (:'w))) /\ (f = 0w))`, 468 gen_tac 469 \\ DISCH_TAC 470 \\ SUBGOAL_THEN ``1 <= (x:real) /\ x < 2 /\ 1 < dimindex (:'w) \/ 471 1 <= ~x /\ ~x < 2 /\ 1 < dimindex (:'w)`` 472 (DISJ_CASES_THEN 473 (Q.X_CHOOSE_THEN `e` (Q.X_CHOOSE_THEN `f` assume_tac) o 474 MATCH_MP error_bound_lemma4)) 475 >- (simp [] \\ pop_assum mp_tac \\ REAL_ARITH_TAC) 476 >| [qexists_tac `0w`, qexists_tac `1w`] 477 \\ qexists_tac `e` 478 \\ qexists_tac `f` 479 \\ ntac 2 (fs [float_to_real_def, wordsTheory.INT_MAX_LT_DIMWORD, 480 wordsTheory.INT_MIN_LT_DIMWORD, realTheory.REAL_DIV_REFL, 481 DECIDE ``0 < x ==> x <> 0n``, wordsTheory.ZERO_LT_INT_MAX, 482 REAL_ARITH ``abs (-2 - x) = abs (2 - -x)``, 483 REAL_ARITH ``abs (-1 * y - x) = abs (y - -x)``]) 484 \\ rfs [DECIDE ``0 < x ==> x <> 0n``, wordsTheory.ZERO_LT_INT_MAX] 485 ); 486 487(* ------------------------------------------------------------------------ *) 488 489val REAL_LE_LCANCEL_IMP = 490 METIS_PROVE [REAL_LE_LMUL] ``!x y z. 0r < x /\ x * y <= x * z ==> y <= z`` 491 492val lem = Q.prove( 493 `!a x. 494 1 < a ==> (2 pow (a - 2) * inv (2 pow (a - 1 + x)) = inv (2 pow (x + 1)))`, 495 rw [realTheory.REAL_INV_1OVER, realTheory.mult_ratr, realTheory.mult_ratl, 496 realTheory.REAL_EQ_RDIV_EQ, GSYM realTheory.POW_ADD, 497 realTheory.REAL_DIV_REFL] 498 ); 499 500val error_bound_lemma6 = Q.prove( 501 `!expw fracw x. 502 0 <= x /\ x < inv (2 pow (2 ** (expw - 1) - 2)) /\ 503 1 < expw /\ 0 < fracw ==> 504 ?n. n <= 2 EXP fracw /\ 505 abs (x - 2 / 2 pow (2 ** (expw - 1) - 1) * &n / 2 pow fracw) <= 506 inv (2 pow (2 ** (expw - 1) - 1 + fracw))`, 507 rpt strip_tac 508 \\ Q.SPECL_THEN [`fracw`, `2 pow (2 ** (expw - 1) - 2) * x`] mp_tac 509 error_bound_lemma2 510 \\ Lib.W (Lib.C SUBGOAL_THEN MP_TAC o lhand o lhand o snd) 511 >- (conj_tac 512 >- (match_mp_tac realTheory.REAL_LE_MUL \\ simp []) 513 \\ qpat_x_assum `x < _` mp_tac 514 \\ simp [realTheory.REAL_INV_1OVER, realTheory.REAL_LT_RDIV_EQ, 515 realTheory.lt_ratr] 516 \\ metis_tac [REAL_MUL_COMM] 517 ) 518 \\ DISCH_THEN (fn th => rewrite_tac [th]) 519 \\ DISCH_THEN (Q.X_CHOOSE_THEN `n` strip_assume_tac) 520 \\ qexists_tac `n` 521 \\ asm_rewrite_tac [] 522 \\ qspec_then `2 pow (2 ** (expw - 1) - 2)` match_mp_tac REAL_LE_LCANCEL_IMP 523 \\ conj_tac 524 >- simp [] 525 \\ rewrite_tac 526 [realTheory.ABS_MUL 527 |> GSYM 528 |> Q.SPEC `2 pow (2 ** (expw - 1) - 2)` 529 |> SIMP_RULE (srw_ss()) [GSYM realTheory.POW_ABS] 530 ] 531 \\ `1n < 2 ** (expw - 1)` 532 by metis_tac [EVAL ``2n ** 0``, bitTheory.TWOEXP_MONO, 533 DECIDE ``1n < n ==> 0 < n - 1``] 534 \\ fs [realTheory.REAL_SUB_LDISTRIB, realTheory.REAL_MUL_ASSOC, real_div, lem, 535 DECIDE ``1 < n ==> (SUC (n - 2) = n - 1)``, realTheory.REAL_MUL_RINV, 536 realTheory.pow 537 |> CONJUNCT2 538 |> GSYM 539 |> ONCE_REWRITE_RULE [REAL_MUL_COMM] 540 ] 541 ); 542 543(* ------------------------------------------------------------------------ *) 544 545val lem = Q.prove( 546 `!n. &(2 * 2 ** n) = 2 * 2 pow n`, 547 simp [realTheory.REAL_OF_NUM_POW] 548 ); 549 550val error_bound_lemma7 = Q.prove( 551 `!x. 0 <= x /\ x < inv (2 pow (bias (:'w) - 1)) /\ 1 < dimindex (:'w) ==> 552 ?e f. 553 abs (float_to_real <| Sign := 0w; 554 Exponent := e : 'w word; 555 Significand := f : 't word |> - x) <= 556 inv (2 pow (bias (:'w) + dimindex (:'t))) /\ 557 ((e = 0w) \/ (e = 1w) /\ (f = 0w))`, 558 gen_tac 559 \\ DISCH_TAC 560 \\ Q.SUBGOAL_THEN 561 `0 <= x /\ x < inv (2 pow (2 ** (dimindex (:'w) - 1) - 2)) /\ 562 1 < dimindex (:'w) /\ 0 < dimindex (:'t)` assume_tac 563 >- fs [wordsTheory.INT_MAX_def, wordsTheory.INT_MIN_def] 564 \\ FIRST_ASSUM (Q.X_CHOOSE_THEN `n` MP_TAC o MATCH_MP error_bound_lemma6) 565 \\ DISCH_THEN 566 (CONJUNCTS_THEN2 567 (strip_assume_tac o REWRITE_RULE [arithmeticTheory.LESS_OR_EQ]) 568 ASSUME_TAC) 569 >- ( 570 qexists_tac `0w` 571 \\ qexists_tac `n2w n` 572 \\ fs [float_to_real_def, GSYM wordsTheory.dimword_def, 573 wordsTheory.INT_MAX_def, wordsTheory.INT_MIN_def] 574 \\ simp [Once realTheory.ABS_SUB] 575 \\ fs [realTheory.mult_rat, realTheory.mult_ratl, 576 Once realTheory.div_ratl] 577 ) 578 \\ qexists_tac `1w` 579 \\ qexists_tac `0w` 580 \\ fs [float_to_real_def, wordsTheory.INT_MAX_def, wordsTheory.INT_MIN_def] 581 \\ simp [Once realTheory.ABS_SUB] 582 \\ rfs [realTheory.mult_rat, realTheory.mult_ratl, Once realTheory.div_ratl, 583 realTheory.REAL_DIV_RMUL_CANCEL, lem] 584 ); 585 586(* ------------------------------------------------------------------------ *) 587 588val error_bound_lemma8 = Q.prove( 589 `!x. abs x < inv (2 pow (bias (:'w) - 1)) /\ 1 < dimindex (:'w) ==> 590 ?s e f. 591 abs (float_to_real <| Sign := s; 592 Exponent := e : 'w word; 593 Significand := f : 't word |> - x) <= 594 inv (2 pow (bias (:'w) + dimindex (:'t))) /\ 595 ((e = 0w) \/ (e = 1w) /\ (f = 0w))`, 596 gen_tac 597 \\ DISCH_TAC 598 \\ SUBGOAL_THEN 599 ``0 <= x /\ x < inv (2 pow (bias (:'w) - 1)) /\ 1 < dimindex (:'w) \/ 600 0 <= ~x /\ ~x < inv (2 pow (bias (:'w) - 1)) /\ 1 < dimindex (:'w) `` 601 (DISJ_CASES_THEN 602 (Q.X_CHOOSE_THEN `e` (Q.X_CHOOSE_THEN `f` assume_tac) o 603 MATCH_MP error_bound_lemma7)) 604 >- (simp [] \\ pop_assum mp_tac \\ REAL_ARITH_TAC) 605 >| [qexists_tac `0w`, qexists_tac `1w`] 606 \\ qexists_tac `e` 607 \\ qexists_tac `f` 608 \\ ntac 2 (fs [float_to_real_def, wordsTheory.INT_MAX_LT_DIMWORD, 609 wordsTheory.INT_MIN_LT_DIMWORD, REAL_MUL_ASSOC, 610 DECIDE ``0 < x ==> x <> 0n``, wordsTheory.ZERO_LT_INT_MAX, 611 REAL_ARITH ``abs (y - -x) = abs (-1 * y - x)``]) 612 \\ rfs [DECIDE ``0 < x ==> x <> 0n``, wordsTheory.ZERO_LT_INT_MAX] 613 ); 614 615(* ------------------------------------------------------------------------ *) 616 617val float_to_real_scale_up = Q.prove( 618 `!s e f k. 619 e <> 0w /\ (e + n2w k <> 0w) /\ (w2n (e + n2w k) = w2n e + k) ==> 620 (float_to_real <| Sign := s; 621 Exponent := e + n2w k : 'w word; 622 Significand := f : 't word |> = 623 2 pow k * float_to_real <| Sign := s; 624 Exponent := e : 'w word; 625 Significand := f : 't word |>)`, 626 simp [float_to_real_def, REAL_POW_ADD, real_div, 627 AC REAL_MUL_ASSOC REAL_MUL_COMM] 628 ); 629 630val float_to_real_scale_down = Q.prove( 631 `!s e f k. 632 e <> 0w /\ n2w k <> e /\ 633 (w2n (e - n2w k + n2w k) = w2n (e - n2w k) + k) ==> 634 (float_to_real <| Sign := s; 635 Exponent := e - n2w k : 'w word; 636 Significand := f : 't word |> = 637 inv (2 pow k) * 638 float_to_real <| Sign := s; 639 Exponent := e : 'w word; 640 Significand := f : 't word |>)`, 641 rpt strip_tac 642 \\ `e + -n2w k <> 0w /\ (e = (e - n2w k) + n2w k)` 643 by metis_tac [wordsTheory.WORD_EQ_SUB_ZERO, wordsTheory.WORD_SUB_INTRO, 644 wordsTheory.WORD_LESS_NOT_EQ, wordsTheory.WORD_SUB_ADD] 645 \\ pop_assum (fn th => CONV_TAC (RAND_CONV (ONCE_REWRITE_CONV [th]))) 646 \\ asm_simp_tac (std_ss++simpLib.type_ssfrag ``:('a, 'b) float``) 647 [float_to_real_def, POW_ADD, AC REAL_MUL_ASSOC REAL_MUL_COMM] 648 \\ simp [SIMP_CONV (srw_ss()) [] ``a + b + -b : 'a word``, 649 REAL_MUL_ASSOC, realTheory.mult_ratr, REAL_MUL_LINV, POW_NZ, 650 REAL_ARITH ``inv a * b * c * a * d = (inv a * a) * b * c * d``] 651 ); 652 653(* ------------------------------------------------------------------------ *) 654 655val two_times_bias_lt = Q.prove( 656 `bias (:'a) + bias (:'a) < dimword (:'a) - 1`, 657 simp [wordsTheory.INT_MAX_def, wordsTheory.INT_MAX_def, 658 GSYM wordsTheory.dimword_IS_TWICE_INT_MIN, 659 DECIDE ``1n < n ==> 0 < n - 1``] 660 ); 661 662val int_min_bias_plus1 = Q.prove( 663 `INT_MIN (:'a) = INT_MAX (:'a) + 1`, 664 simp [wordsTheory.INT_MAX_def, DECIDE ``0n < n ==> (n - 1 + 1 = n)``]) 665 666 667val lem = Q.prove( 668 `1 < dimindex (:'a) ==> 669 2 pow (UINT_MAX (:'a) - 1) / 2 pow (INT_MAX (:'a)) <= 2 pow (INT_MAX (:'a))`, 670 rw [GSYM POW_ADD, realTheory.REAL_LE_LDIV_EQ] 671 \\ match_mp_tac REAL_POW_MONO 672 \\ simp [wordsTheory.UINT_MAX_def, wordsTheory.ZERO_LT_INT_MAX, 673 DECIDE ``0n < a ==> 0 < 2 * a``, 674 wordsTheory.dimword_IS_TWICE_INT_MIN] 675 \\ simp [wordsTheory.INT_MAX_def] 676 ); 677 678val lem2 = Q.prove( 679 `!a b c d : real. 0 < b /\ 0 <= c /\ a < b * c /\ b <= d ==> a < c * d`, 680 rw [] 681 \\ `?p. 0 <= p /\ (d = b + p)` 682 by (qexists_tac `d - b` 683 \\ simp [REAL_ARITH ``b <= d : real ==> (b + (d - b) = d) /\ 0 <= d - b``] 684 ) 685 \\ simp [realTheory.REAL_LDISTRIB, AC REAL_MUL_ASSOC REAL_MUL_COMM] 686 \\ `0 <= c * p` 687 by simp[ 688 realTheory.REAL_LE_MUL2 689 |> Q.SPECL [`0`, `c`, `0`, `p`] 690 |> SIMP_RULE (srw_ss()) []] 691 \\ simp [REAL_ARITH ``0 <= c /\ a < b ==> a < b + c : real``] 692 ); 693 694val error_bound_big1 = Q.prove( 695 `!x k. 2 pow k <= abs x /\ abs x < 2 pow SUC k /\ 696 abs x < threshold (:'t # 'w) /\ 1 < dimindex (:'w) ==> 697 ?a : ('t, 'w) float. 698 float_is_finite a /\ 699 abs (float_to_real a - x) <= 2 pow k / 2 pow (dimindex (:'t) + 1)`, 700 rpt strip_tac 701 \\ qspec_then `x / 2 pow k` mp_tac error_bound_lemma5 702 \\ Lib.W (Lib.C SUBGOAL_THEN mp_tac o lhand o lhand o snd) 703 >- (simp [ABS_DIV, GSYM realTheory.POW_ABS, ABS_N, POW_NZ, REAL_POW_LT, 704 REAL_LT_LDIV_EQ, GSYM (CONJUNCT2 pow)] 705 \\ match_mp_tac realTheory.REAL_LE_RDIV 706 \\ simp [realTheory.REAL_POW_LT]) 707 \\ DISCH_THEN (fn th => rewrite_tac [th]) 708 \\ `2 pow k < threshold (:'t # 'w)` by metis_tac [REAL_LET_TRANS] 709 \\ `k < bias (:'w) + 1` 710 by (spose_not_then (assume_tac o REWRITE_RULE [arithmeticTheory.NOT_LESS]) 711 \\ `2r pow (bias (:'w) + 1) <= 2 pow k` 712 by (match_mp_tac REAL_POW_MONO \\ simp []) 713 \\ `2r pow (bias (:'w) + 1) < threshold (:'t # 'w)` 714 by metis_tac [REAL_LET_TRANS] 715 \\ pop_assum mp_tac 716 \\ simp [threshold, realTheory.REAL_LT_RDIV_EQ, 717 GSYM realTheory.REAL_OF_NUM_POW, GSYM realTheory.POW_ADD, 718 wordsTheory.INT_MAX_def, wordsTheory.INT_MIN_def, 719 wordsTheory.UINT_MAX_def, 720 DECIDE ``0n < n ==> (n - 1 + 1 = n) /\ (SUC (n - 1) = n)``, 721 GSYM (CONJUNCT2 arithmeticTheory.EXP)] 722 \\ simp [realTheory.REAL_NOT_LT, GSYM wordsTheory.dimword_def, 723 realTheory.REAL_SUB_LDISTRIB, realTheory.mult_ratr, 724 DECIDE ``1n < n ==> (SUC (n - 2) = n - 1)``, 725 DECIDE ``1n < n ==> n <> 0``, 726 GSYM 727 (ONCE_REWRITE_RULE [REAL_MUL_COMM] (CONJUNCT2 realTheory.pow)) 728 ] 729 \\ match_mp_tac (REAL_ARITH ``0 < a /\ 0 <= b ==> (a - b <= a : real)``) 730 \\ simp [realTheory.REAL_LE_RDIV_EQ, DECIDE ``1n < n ==> 0 < 2 * n``] 731 ) 732 \\ strip_tac 733 >| [all_tac, 734 Cases_on `k = bias (:'w)` 735 >- ( 736 spose_not_then kall_tac 737 \\ qpat_x_assum `abs _ <= inv (2 pow _)` 738 (mp_tac o (MATCH_MP (REAL_ARITH 739 ``abs (a - b) <= c ==> abs(a) <= abs(b) + c``))) 740 \\ simp [realTheory.REAL_NOT_LE, REAL_ABS_MUL, GSYM POW_ABS, ABS_NEG, 741 ABS_DIV, ABS_N, float_to_real_def, 742 wordsTheory.INT_MIN_LT_DIMWORD, prim_recTheory.LESS_NOT_EQ] 743 \\ simp [int_min_bias_plus1, POW_ADD, realTheory.POW_ONE, 744 realTheory.REAL_LT_ADD_SUB, realTheory.REAL_LT_LDIV_EQ, 745 realTheory.REAL_DIV_LMUL_CANCEL 746 |> Q.SPECL [`2 pow n`, `2`, `1`] 747 |> SIMP_RULE (srw_ss()) []] 748 \\ match_mp_tac lem2 749 \\ qexists_tac `2 pow (UINT_MAX (:'w) - 1) / 2 pow bias (:'w)` 750 \\ fs [threshold_def, pow, lem, AC REAL_MUL_ASSOC REAL_MUL_COMM, 751 realTheory.REAL_LT_RDIV_0, realTheory.REAL_SUB_LE, 752 realTheory.REAL_INV_1OVER, realTheory.REAL_LE_LDIV_EQ, 753 realTheory.POW_2_LE1, REAL_ARITH ``0r < n ==> 0 < 2 * n``, 754 REAL_ARITH ``1r <= n ==> 1 <= 2 * (2 * n)``] 755 ) 756 \\ `k < bias (:'w)` by decide_tac 757 ] 758 \\ ( 759 qexists_tac `<| Sign := s; 760 Exponent := e + n2w k; 761 Significand := f |> : ('t, 'w) float` 762 \\ conj_tac 763 >- simp [float_tests, wordsTheory.word_add_n2w, int_min_bias_plus1, 764 wordsTheory.word_2comp_n2w, two_times_bias_lt, 765 DECIDE ``k < b + 1n /\ (b + b) < w - 1 ==> 766 k + b < w /\ k + b <> w - 1``, 767 DECIDE ``k < b /\ (b + b) < w - 1n ==> 768 k + (b + 1) < w /\ k + (b + 1) <> w - 1``] 769 \\ Q.SUBGOAL_THEN 770 `e <> 0w /\ e + n2w k <> 0w /\ (w2n (e + n2w k) = w2n e + k)` 771 (fn th => rewrite_tac [MATCH_MP float_to_real_scale_up th]) 772 >- ( 773 fs [wordsTheory.INT_MAX_LT_DIMWORD, prim_recTheory.LESS_NOT_EQ, 774 wordsTheory.INT_MIN_LT_DIMWORD, wordsTheory.ZERO_LT_INT_MAX, 775 wordsTheory.word_add_n2w, two_times_bias_lt, 776 DECIDE ``k < b + 1n /\ (b + b) < w - 1 ==> 777 k + b < w /\ k + b <> w - 1``] 778 \\ simp [int_min_bias_plus1, two_times_bias_lt, 779 DECIDE ``k < b /\ (b + b) < w - 1n ==> 780 k + (b + 1) < w /\ k + (b + 1) <> w - 1``] 781 ) 782 \\ match_mp_tac REAL_LE_LCANCEL_IMP 783 \\ qexists_tac `inv (2 pow k)` 784 \\ conj_tac 785 >- simp [REAL_LT_INV_EQ, REAL_POW_LT] 786 \\ `!x. inv (2 pow k) * abs x = abs (inv (2 pow k) * x)` 787 by rewrite_tac 788 [REAL_ABS_MUL, REAL_ABS_INV, GSYM realTheory.POW_ABS, ABS_N] 789 \\ qpat_x_assum `zz <= inv (2 pow _)` mp_tac 790 \\ simp [REAL_SUB_LDISTRIB, REAL_MUL_ASSOC, real_div, POW_NZ, 791 REAL_MUL_LINV, float_to_real_def] 792 \\ simp [AC REAL_MUL_COMM REAL_MUL_ASSOC, wordsTheory.ZERO_LT_INT_MAX, 793 wordsTheory.INT_MAX_LT_DIMWORD, prim_recTheory.LESS_NOT_EQ 794 ] 795 ) 796 ); 797 798val error_bound_big = Q.prove( 799 `!k x. 800 2 pow k <= abs x /\ abs x < 2 pow (SUC k) /\ 801 abs x < threshold (:'t # 'w) /\ 1 < dimindex (:'w) ==> 802 abs (error (:'t # 'w) x) <= 2 pow k / 2 pow (dimindex (:'t) + 1)`, 803 prove_tac [error_bound_big1, error_at_worst_lemma, REAL_LE_TRANS]) 804 805(* ------------------------------------------------------------------------ *) 806 807val suc_bias_lt_dimword = Q.prove( 808 `1 < dimindex (:'a) ==> bias (:'a) + 1 < dimword (:'a)`, 809 simp [wordsTheory.INT_MAX_def, wordsTheory.dimword_IS_TWICE_INT_MIN, 810 DECIDE ``0n < n ==> (n - 1 + 1 = n)``] 811 ); 812 813val error_bound_small1 = Q.prove( 814 `!x k. inv (2 pow SUC k) <= abs x /\ abs x < inv (2 pow k) /\ 815 k < bias (:'w) - 1 /\ 1 < dimindex (:'w) ==> 816 ?a : ('t, 'w) float. 817 float_is_finite a /\ 818 abs (float_to_real a - x) <= 819 inv (2 pow SUC k * 2 pow (dimindex (:'t) + 1))`, 820 rpt strip_tac 821 \\ qspec_then `x * 2 pow (SUC k)` mp_tac error_bound_lemma5 822 \\ Lib.W (Lib.C SUBGOAL_THEN mp_tac o lhand o lhand o snd) 823 >- (fs [ABS_MUL, GSYM POW_ABS, REAL_INV_1OVER, REAL_LE_LDIV_EQ, 824 REAL_LT_RDIV_EQ, REAL_POW_LT] 825 \\ simp [pow, REAL_ARITH ``a * (2r * b) < 2 = a * b < 1``]) 826 \\ DISCH_THEN (fn th => rewrite_tac [th]) 827 \\ DISCH_THEN 828 (Q.X_CHOOSE_THEN `s` 829 (Q.X_CHOOSE_THEN `e` 830 (Q.X_CHOOSE_THEN `f` (REPEAT_TCL CONJUNCTS_THEN assume_tac)))) 831 \\ qexists_tac `<| Sign := s; 832 Exponent := e - n2w (SUC k); 833 Significand := f |> : ('t, 'w) float` 834 \\ conj_tac 835 >- ( 836 fs [float_tests, wordsTheory.WORD_LITERAL_ADD, int_min_bias_plus1] 837 \\ `bias (:'w) - SUC k < dimword (:'w)` 838 by (match_mp_tac arithmeticTheory.LESS_TRANS 839 \\ qexists_tac `bias (:'w)` 840 \\ simp [wordsTheory.INT_MAX_LT_DIMWORD] 841 ) 842 \\ `bias (:'w) + 1 - SUC k < dimword (:'w)` 843 by (match_mp_tac arithmeticTheory.LESS_TRANS 844 \\ qexists_tac `bias (:'w) + 1` 845 \\ simp [wordsTheory.INT_MAX_def, 846 wordsTheory.dimword_IS_TWICE_INT_MIN, 847 DECIDE ``0n < n ==> (n - 1 + 1 = n)``] 848 ) 849 \\ simp_tac std_ss [wordsTheory.WORD_NEG_1, wordsTheory.word_T_def] 850 \\ simp [wordsTheory.BOUND_ORDER, wordsTheory.INT_MAX_LT_DIMWORD] 851 \\ simp [wordsTheory.INT_MAX_def, wordsTheory.UINT_MAX_def, 852 wordsTheory.dimword_IS_TWICE_INT_MIN, 853 DECIDE ``0 < a /\ 0 < b ==> a - b <> 2 * a - 1n`` 854 ] 855 ) 856 \\ `e <> 0w /\ n2w (SUC k) <> e /\ 857 (w2n (e - n2w (SUC k) + n2w (SUC k)) = w2n (e - n2w (SUC k)) + SUC k)` 858 by ( 859 `SUC k < dimword (:'w)` 860 by metis_tac [wordsTheory.ZERO_LT_INT_MAX, wordsTheory.INT_MAX_LT_DIMWORD, 861 arithmeticTheory.LESS_TRANS, 862 DECIDE ``0n < b /\ k < b - 1 ==> SUC k < b``] 863 \\ fs [wordsTheory.INT_MAX_LT_DIMWORD, wordsTheory.INT_MIN_LT_DIMWORD, 864 prim_recTheory.LESS_NOT_EQ, 865 int_min_bias_plus1, suc_bias_lt_dimword, 866 SIMP_CONV (srw_ss()) [] ``a + b + -b : 'a word``, 867 SIMP_CONV (srw_ss()) [] ``b + a + -b : 'a word``] 868 \\ simp [wordsTheory.WORD_LITERAL_ADD, wordsTheory.INT_MAX_LT_DIMWORD, 869 DECIDE ``0n < a /\ a < n ==> (a - SUC k < n) /\ 870 (a + 1 - SUC k < n)``] 871 ) 872 \\ NO_STRIP_FULL_SIMP_TAC std_ss [float_to_real_scale_down] 873 \\ match_mp_tac REAL_LE_LCANCEL_IMP 874 \\ qexists_tac `2 pow (SUC k)` 875 \\ `!x. 2 pow (SUC k) * abs x = abs (2 pow (SUC k) * x)` 876 by rewrite_tac [REAL_ABS_MUL, REAL_ABS_INV, GSYM POW_ABS, ABS_N] 877 \\ `!a b. 0 < a ==> (a * (inv a * b) = b)` 878 by simp [REAL_MUL_ASSOC, REAL_MUL_RINV, REAL_POS_NZ] 879 \\ simp [REAL_POW_LT, REAL_SUB_LDISTRIB, REAL_POS_NZ, REAL_INV_MUL] 880 \\ NO_STRIP_FULL_SIMP_TAC (srw_ss()) [AC REAL_MUL_ASSOC REAL_MUL_COMM] 881 ); 882 883val REAL_LE_INV2 = Q.prove( 884 `!x y. 0 < x /\ x <= y ==> inv y <= inv x`, 885 metis_tac [REAL_LE_LT, REAL_LT_INV]) 886 887val lem = Q.prove( 888 `!n m. 2n <= n /\ 0 < m ==> 889 2 pow (n - 1) < 2 pow (2 * n - 1) - 2 pow (2 * n - 2) / &(4 * m)`, 890 rw [realTheory.REAL_LT_SUB_LADD] 891 \\ `1 < 4 * m /\ 0 < 4 * m` by decide_tac 892 \\ `!x y:real. x < y = x * &(4 * m) < y * &(4 * m)` 893 by metis_tac [realTheory.REAL_LT_RMUL, 894 SIMP_CONV (srw_ss()) [] ``0r < &(4 * m)``] 895 \\ pop_assum (fn th => once_rewrite_tac [th]) 896 \\ simp [realTheory.REAL_ADD_RDISTRIB, realTheory.REAL_DIV_RMUL, 897 REAL_ARITH ``!n. x * (4 * n) = 2 * x * n + 2 * x * n : real`` 898 |> Q.SPEC `&n` 899 |> SIMP_RULE (srw_ss()) []] 900 \\ CONV_TAC (LAND_CONV (ONCE_REWRITE_CONV [REAL_ADD_COMM])) 901 \\ match_mp_tac realTheory.REAL_LTE_ADD2 902 \\ Q.SPECL_THEN [`2`, `n`] imp_res_tac arithmeticTheory.LESS_EQUAL_ADD 903 \\ fs [] 904 \\ rw [realTheory.REAL_DOUBLE] 905 >- ( 906 simp [realTheory.REAL_OF_NUM_POW, DECIDE ``x + 3 = SUC (x + 2)``, 907 arithmeticTheory.EXP, arithmeticTheory.RIGHT_ADD_DISTRIB, 908 arithmeticTheory.LEFT_ADD_DISTRIB] 909 \\ rewrite_tac [arithmeticTheory.MULT_ASSOC, 910 arithmeticTheory.LT_MULT_CANCEL_LBARE] 911 \\ simp [] 912 ) 913 \\ `m <> 0` by decide_tac 914 \\ asm_simp_tac std_ss 915 [realTheory.REAL_NZ_IMP_LT, realTheory.REAL_LE_RMUL, REAL_MUL_ASSOC] 916 \\ asm_simp_tac std_ss 917 [realTheory.REAL_LE_LMUL, GSYM REAL_MUL_ASSOC, REAL_ARITH ``0 < 2r``] 918 \\ simp [GSYM (CONJUNCT2 pow)] 919 \\ match_mp_tac REAL_POW_MONO 920 \\ simp [] 921 ); 922 923val threshold_gt1 = Q.prove( 924 `1 < dimindex (:'w) ==> 1 < threshold (:'t # 'w)`, 925 simp [threshold, realTheory.REAL_INV_1OVER, realTheory.REAL_LT_RDIV_EQ, 926 realTheory.mult_ratl, realTheory.mult_ratr, 927 GSYM realTheory.REAL_OF_NUM_POW, prim_recTheory.LESS_NOT_EQ, 928 wordsTheory.ZERO_LT_INT_MAX, two_pow_over_pre, 929 realTheory.REAL_SUB_LDISTRIB, DECIDE ``0n < n ==> (SUC (n - 1) = n)``, 930 GSYM (CONJUNCT2 arithmeticTheory.EXP)] 931 \\ simp [wordsTheory.UINT_MAX_def, wordsTheory.INT_MAX_def, 932 wordsTheory.dimword_IS_TWICE_INT_MIN] 933 \\ qabbrev_tac `n = INT_MIN (:'w)` 934 \\ qabbrev_tac `m = INT_MIN (:'t)` 935 \\ strip_tac 936 \\ `2n <= n` by simp [Abbr `n`, wordsTheory.INT_MIN_def] 937 \\ `0n < m` by simp [Abbr `m`, wordsTheory.INT_MIN_def] 938 \\ simp [lem] 939 ); 940 941val error_bound_small = Q.prove( 942 `!k x. 943 inv (2 pow (SUC k)) <= abs x /\ abs x < inv (2 pow k) /\ 944 k < bias (:'w) - 1 /\ 1 < dimindex (:'w) ==> 945 abs (error (:'t # 'w) x) <= 946 inv (2 pow (SUC k) * 2 pow (dimindex (:'t) + 1))`, 947 rpt strip_tac 948 \\ `?a : ('t, 'w) float. 949 float_is_finite a /\ 950 abs (float_to_real a - x) <= 951 inv (2 pow (SUC k) * 2 pow (dimindex (:'t) + 1))` 952 by metis_tac [error_bound_small1] 953 \\ match_mp_tac REAL_LE_TRANS 954 \\ qexists_tac `abs (float_to_real a - x)` 955 \\ simp [] 956 \\ match_mp_tac error_at_worst_lemma 957 \\ simp [] 958 \\ match_mp_tac REAL_LT_TRANS 959 \\ qexists_tac `inv (2 pow k)` 960 \\ simp [] 961 \\ match_mp_tac REAL_LET_TRANS 962 \\ qexists_tac `inv 1` 963 \\ conj_tac 964 >- (match_mp_tac REAL_LE_INV2 \\ simp [REAL_POW_LE_1]) 965 \\ simp [realTheory.REAL_INV_1OVER, threshold_gt1] 966 ); 967 968(* ------------------------------------------------------------------------ *) 969 970val lem = Q.prove( 971 `1 < dimindex (:'w) ==> -1w <> (1w : 'w word)`, 972 srw_tac [wordsLib.WORD_BIT_EQ_ss] [] 973 \\ qexists_tac `1` 974 \\ simp [wordsTheory.word_index] 975 ); 976 977val error_bound_tiny = Q.prove( 978 `!x. abs x < inv (2 pow (bias (:'w) - 1)) /\ 1 < dimindex (:'w) ==> 979 abs (error (:'t # 'w) x) <= inv (2 pow (bias (:'w) + dimindex (:'t)))`, 980 strip_tac 981 \\ DISCH_TAC 982 \\ `?a : ('t, 'w) float. 983 float_is_finite a /\ 984 abs (float_to_real a - x) <= inv (2 pow (bias (:'w) + dimindex (:'t)))` 985 by (pop_assum (fn th => mp_tac (MATCH_MP error_bound_lemma8 th) 986 \\ assume_tac th) 987 \\ DISCH_THEN 988 (Q.X_CHOOSE_THEN `s` 989 (Q.X_CHOOSE_THEN `e` 990 (Q.X_CHOOSE_THEN `f` (REPEAT_TCL CONJUNCTS_THEN assume_tac)))) 991 \\ qexists_tac `<|Sign := s; Exponent := e; Significand := f|>` 992 \\ fs [float_tests, wordsTheory.word_T_not_zero, lem] 993 ) 994 \\ match_mp_tac REAL_LE_TRANS 995 \\ qexists_tac `abs (float_to_real a - x)` 996 \\ simp [] 997 \\ match_mp_tac error_at_worst_lemma 998 \\ asm_rewrite_tac [] 999 \\ match_mp_tac REAL_LT_TRANS 1000 \\ qexists_tac `inv (2 pow (bias (:'w) - 1))` 1001 \\ asm_rewrite_tac [] 1002 \\ match_mp_tac realTheory.REAL_LET_TRANS 1003 \\ qexists_tac `1` 1004 \\ simp [realTheory.REAL_INV_1OVER, realTheory.REAL_LE_LDIV_EQ, threshold_gt1] 1005 \\ CONV_TAC (LAND_CONV (ONCE_REWRITE_CONV [GSYM (EVAL ``2r pow 0``)])) 1006 \\ match_mp_tac REAL_POW_MONO 1007 \\ simp [] 1008 ); 1009 1010(* ------------------------------------------------------------------------- 1011 Stronger versions not requiring exact location of the interval. 1012 ------------------------------------------------------------------------- *) 1013 1014val lem = Q.prove( 1015 `!n. 1 < n ==> (2 * inv (2 pow (n - 1)) = inv (2 pow (n - 2)))`, 1016 rw [realTheory.REAL_INV_1OVER, realTheory.REAL_EQ_RDIV_EQ, 1017 REAL_ARITH ``2 * (a:real) * b = a * (2 * b)``, GSYM (CONJUNCT2 pow), 1018 DECIDE ``1 < n ==> (SUC (n - 2) = n - 1)``, 1019 realTheory.REAL_DIV_RMUL 1020 ] 1021 ); 1022 1023val error_bound_norm_strong = Q.prove( 1024 `!x j. 1025 abs x < threshold (:'t # 'w) /\ 1026 abs x < 2 pow (SUC j) / 2 pow (bias (:'w) - 1) /\ 1 < bias (:'w) ==> 1027 abs (error (:'t # 'w) x) <= 2 pow j / 2 pow (bias (:'w) + dimindex (:'t))`, 1028 gen_tac 1029 \\ Induct 1030 >- ( 1031 rw_tac std_ss [pow, POW_1, real_div, REAL_MUL_LID, REAL_MUL_RID] 1032 \\ fs [lem] 1033 \\ `1 < dimindex (:'w)` 1034 by ( 1035 spose_not_then assume_tac 1036 \\ `(dimindex (:'w) = 0) \/ (dimindex (:'w) = 1)` by decide_tac 1037 \\ fs [wordsTheory.INT_MAX_def, wordsTheory.INT_MIN_def] 1038 ) 1039 \\ Cases_on `abs x < inv (2 pow (bias (:'w) - 1))` 1040 >- metis_tac [error_bound_tiny] 1041 \\ qspecl_then [`bias (:'w) - 2`, `x`] mp_tac error_bound_small 1042 \\ fs [GSYM REAL_POW_ADD, arithmeticTheory.ADD1, GSYM REAL_NOT_LT] 1043 ) 1044 \\ strip_tac 1045 \\ `1 < dimindex (:'w)` 1046 by ( 1047 spose_not_then assume_tac 1048 \\ `(dimindex (:'w) = 0) \/ (dimindex (:'w) = 1)` by decide_tac 1049 \\ fs [wordsTheory.INT_MAX_def, wordsTheory.INT_MIN_def] 1050 ) 1051 \\ Cases_on `abs x < 2 pow SUC j / 2 pow (bias (:'w) - 1)` 1052 >- (match_mp_tac REAL_LE_TRANS 1053 \\ qexists_tac `2 pow j / 2 pow (bias (:'w) + dimindex (:'t))` 1054 \\ asm_simp_tac std_ss [real_div, pow] 1055 \\ match_mp_tac realTheory.REAL_LE_RMUL_IMP 1056 \\ simp_tac std_ss [REAL_LE_INV_EQ, POW_POS, REAL_ARITH ``0 <= 2r``, 1057 REAL_ARITH ``0r <= a ==> a <= 2 * a``] 1058 ) 1059 \\ Cases_on `SUC j <= bias (:'w) - 2` 1060 >- ( 1061 `2 pow SUC j / 2 pow (bias (:'w) + dimindex (:'t)) = 1062 inv (2 pow SUC ((bias (:'w) - 2) - SUC j) * 2 pow (dimindex (:'t) + 1))` 1063 by simp [GSYM POW_ADD, realTheory.REAL_EQ_LDIV_EQ, 1064 realTheory.REAL_EQ_RDIV_EQ, 1065 arithmeticTheory.ADD1, REAL_INV_1OVER, realTheory.mult_ratl] 1066 \\ asm_rewrite_tac [] 1067 \\ match_mp_tac error_bound_small 1068 \\ `inv (2 pow (SUC (bias (:'w) - (SUC j + 2)))) = 1069 2 pow SUC j / 2 pow (bias (:'w) - 1)` 1070 by simp [GSYM POW_ADD, realTheory.REAL_EQ_LDIV_EQ, 1071 realTheory.REAL_EQ_RDIV_EQ, 1072 arithmeticTheory.ADD1, REAL_INV_1OVER, realTheory.mult_ratl] 1073 \\ `inv (2 pow (bias (:'w) - (SUC j + 2))) = 1074 2 pow SUC (SUC j) / 2 pow (bias (:'w) - 1)` 1075 by simp [GSYM POW_ADD, realTheory.REAL_EQ_LDIV_EQ, 1076 realTheory.REAL_EQ_RDIV_EQ, 1077 arithmeticTheory.ADD1, REAL_INV_1OVER, realTheory.mult_ratl] 1078 \\ fs [realTheory.REAL_NOT_LT] 1079 ) 1080 \\ `?i. j = (bias (:'w) - 2) + i` 1081 by (qexists_tac `j - (bias (:'w) - 2)` \\ decide_tac) 1082 \\ asm_simp_tac std_ss 1083 [DECIDE ``1n < b ==> (b + i = b - 1 + SUC i) /\ 1084 (SUC (b - 2 + i) = b - 1 + i)``] 1085 \\ simp_tac std_ss [POW_ADD] 1086 \\ simp [realTheory.REAL_DIV_LMUL_CANCEL, arithmeticTheory.ADD1] 1087 \\ match_mp_tac error_bound_big 1088 \\ qpat_x_assum `~(_ < _)` mp_tac 1089 \\ full_simp_tac bool_ss 1090 [realTheory.REAL_NOT_LT, POW_ADD, 1091 DECIDE ``1n < b ==> (SUC (b - 2 + i) = i + (b - 1))``, 1092 DECIDE ``SUC (i + (b - 1)) = SUC i + (b - 1)``, 1093 realTheory.REAL_DIV_RMUL_CANCEL 1094 |> Q.SPECL [`2 pow n`, `a`, `1`] 1095 |> SIMP_RULE (srw_ss()) [] 1096 ] 1097 ); 1098 1099(* ------------------------------------------------------------------------- 1100 "1 + Epsilon" property (relative error bounding). 1101 ------------------------------------------------------------------------- *) 1102 1103val bias_imp_dimindex = Q.prove( 1104 `1 < bias (:'a) ==> 1 < dimindex (:'a)`, 1105 rw [wordsTheory.INT_MAX_def, wordsTheory.INT_MIN_def] 1106 \\ spose_not_then assume_tac 1107 \\ `dimindex (:'a) = 1` by simp [DECIDE ``0n < n /\ ~(1 < n) ==> (n = 1)``] 1108 \\ fs [] 1109 ); 1110 1111val lem = Q.prove( 1112 `!n : num. n + SUC (n - 1) <= 2 ** n`, 1113 Induct \\ rw [arithmeticTheory.EXP]) 1114 1115val THRESHOLD_MUL_LT = Q.prove( 1116 `threshold (:'t # 'w) * 2 pow (bias (:'w) - 1) < 2 pow (2 EXP (bias (:'w)))`, 1117 `2 pow (UINT_MAX (:'w) - 1) * inv (2 pow (bias (:'w))) = 2 pow (bias (:'w))` 1118 by (simp [REAL_INV_1OVER, realTheory.mult_ratr, realTheory.REAL_EQ_LDIV_EQ, 1119 GSYM REAL_POW_ADD] 1120 \\ simp [realTheory.REAL_OF_NUM_POW, wordsTheory.UINT_MAX_def, 1121 wordsTheory.INT_MAX_def, wordsTheory.dimword_IS_TWICE_INT_MIN, 1122 arithmeticTheory.LEFT_SUB_DISTRIB]) 1123 \\ asm_simp_tac std_ss [threshold_def, real_div] 1124 \\ rewrite_tac 1125 [GSYM REAL_MUL_ASSOC, REAL_SUB_RDISTRIB, REAL_SUB_LDISTRIB, 1126 GSYM pow, GSYM POW_ADD] 1127 \\ match_mp_tac REAL_LTE_TRANS 1128 \\ qexists_tac `2 pow (bias (:'w) + SUC (bias (:'w) - 1))` 1129 \\ conj_tac 1130 >- ( 1131 match_mp_tac (REAL_ARITH ``0 < a /\ 0r < x ==> (a - x < a)``) 1132 \\ simp [realTheory.REAL_LT_RMUL_0, realTheory.REAL_LT_LMUL_0, 1133 realTheory.REAL_LT_INV_EQ] 1134 ) 1135 \\ match_mp_tac REAL_POW_MONO 1136 \\ simp [lem] 1137 ); 1138 1139val lem = Q.prove( 1140 `!a b c:real. 0 < b ==> ((a / b) * c = a * (c / b))`, 1141 rw [realTheory.mult_ratl, realTheory.mult_ratr] 1142 ); 1143 1144val LT_THRESHOLD_LT_POW_INV = Q.prove( 1145 `!x. 1 < dimindex (:'w) ==> x < threshold (:'t # 'w) ==> 1146 x < 2 pow (UINT_MAX (:'w) - 1) / 2 pow (bias (:'w) - 1)`, 1147 gen_tac 1148 \\ strip_tac 1149 \\ simp [threshold] 1150 \\ match_mp_tac (REAL_ARITH ``b < c ==> (a < b ==> a < c : real)``) 1151 \\ simp [realTheory.REAL_LT_LDIV_EQ, GSYM realTheory.REAL_OF_NUM_POW, lem, 1152 two_pow_over_pre, wordsTheory.ZERO_LT_INT_MAX, 1153 realTheory.REAL_LT_LMUL] 1154 \\ match_mp_tac (REAL_ARITH ``0r < a /\ 0r < b ==> a - b < a``) 1155 \\ `0r < &(2 * dimword (:'t))` by simp [DECIDE ``0n < n ==> 0 < 2 * n``] 1156 \\ simp [realTheory.REAL_LT_RDIV_0] 1157 ); 1158 1159val real_pos_in_binade = Q.prove( 1160 `!x. normalizes (:'t # 'w) x /\ 0 <= x ==> 1161 ?j. j <= UINT_MAX (:'w) - 2 /\ 2 pow j / 2 pow (bias (:'w) - 1) <= x /\ 1162 x < 2 pow (SUC j) / 2 pow (bias (:'w) - 1)`, 1163 rw_tac arith_ss [normalizes_def, abs] 1164 \\ imp_res_tac bias_imp_dimindex 1165 \\ qspec_then `\n. 2 pow n / 2 pow (bias (:'w) - 1) <= x` 1166 mp_tac arithmeticTheory.EXISTS_GREATEST 1167 \\ Lib.W (Lib.C SUBGOAL_THEN mp_tac o lhs o lhand o snd) 1168 >- ( 1169 conj_tac 1170 >- (qexists_tac `0` \\ asm_simp_tac std_ss [REAL_MUL_LID , pow, real_div]) 1171 \\ qexists_tac `2 EXP (bias (:'w))` 1172 \\ Q.X_GEN_TAC `n` 1173 \\ rw_tac bool_ss 1174 [GSYM real_lt, REAL_LT_RDIV_EQ, REAL_POW_LT, REAL_ARITH ``0 < 2r``] 1175 \\ match_mp_tac REAL_LT_TRANS 1176 \\ qexists_tac `2 pow (2 EXP (bias (:'w)))` 1177 \\ conj_tac 1178 >- ( 1179 match_mp_tac realTheory.REAL_LT_TRANS 1180 \\ qexists_tac `threshold (:'t # 'w) * 2 pow (bias (:'w) - 1)` 1181 \\ simp [REAL_LT_RMUL, THRESHOLD_MUL_LT] 1182 ) 1183 \\ match_mp_tac REAL_POW_MONO_LT 1184 \\ asm_simp_tac bool_ss 1185 [REAL_ARITH ``1 < 2r``, GSYM arithmeticTheory.GREATER_DEF] 1186 ) 1187 \\ DISCH_THEN (fn th => rewrite_tac [th]) 1188 \\ DISCH_THEN 1189 (X_CHOOSE_THEN ``n:num`` 1190 (CONJUNCTS_THEN2 ASSUME_TAC (MP_TAC o Q.SPEC `SUC n`))) 1191 \\ rw_tac arith_ss [REAL_NOT_LE] 1192 \\ qexists_tac `n` 1193 \\ full_simp_tac std_ss [] 1194 \\ imp_res_tac LT_THRESHOLD_LT_POW_INV 1195 \\ `2 pow n / 2 pow (bias (:'w) - 1) < 1196 2 pow (UINT_MAX (:'w) - 1) / 2 pow (bias (:'w) - 1)` 1197 by metis_tac [REAL_LET_TRANS] 1198 \\ spose_not_then assume_tac 1199 \\ `UINT_MAX (:'w) - 1 <= n` by decide_tac 1200 \\ `2 pow (UINT_MAX (:'w) - 1) <= 2 pow n` 1201 by metis_tac [REAL_POW_MONO, REAL_ARITH ``1 <= 2r``] 1202 \\ full_simp_tac std_ss 1203 [REAL_LT_RDIV, REAL_POW_LT, REAL_ARITH ``0 < 2r``, real_lte] 1204 ); 1205 1206val real_neg_in_binade = Q.prove( 1207 `!x. normalizes (:'t # 'w) x /\ 0 <= ~x ==> 1208 ?j. j <= UINT_MAX (:'w) - 2 /\ 2 pow j / 2 pow (bias (:'w) - 1) <= ~x /\ 1209 ~x < 2 pow (SUC j) / 2 pow (bias (:'w) - 1)`, 1210 metis_tac [normalizes_def, ABS_NEG, real_pos_in_binade]) 1211 1212val real_in_binade = Q.prove( 1213 `!x. normalizes (:'t # 'w) x ==> 1214 ?j. j <= UINT_MAX (:'w) - 2 /\ 1215 2 pow j / 2 pow (bias (:'w) - 1) <= abs x /\ 1216 abs x < 2 pow (SUC j) / 2 pow (bias (:'w) - 1)`, 1217 gen_tac 1218 \\ Cases_on `0 <= x` 1219 \\ asm_simp_tac arith_ss [abs, real_neg_in_binade, real_pos_in_binade, 1220 REAL_ARITH ``~(0r <= x) ==> 0 <= ~x``] 1221 ); 1222 1223(* ------------------------------------------------------------------------- *) 1224 1225val error_bound_norm_strong_normalize = Q.prove( 1226 `!x. normalizes (:'t # 'w) x ==> 1227 ?j. abs (error (:'t # 'w) x) <= 1228 2 pow j / 2 pow (bias (:'w) + dimindex (:'t))`, 1229 metis_tac [real_in_binade, error_bound_norm_strong, normalizes_def]) 1230 1231(* ------------------------------------------------------------------------- *) 1232 1233val inv_le = Q.prove( 1234 `!a b. 0 < a /\ 0 < b ==> (inv a <= inv b = b <= a)`, 1235 rw [realTheory.REAL_INV_1OVER, realTheory.REAL_LE_LDIV_EQ, 1236 realTheory.mult_ratl, realTheory.REAL_LE_RDIV_EQ] 1237 ); 1238 1239val relative_bound_lem = Q.prove( 1240 `!x j. x <> 0 ==> 1241 (2 pow j * inv (2 pow (bias (:'w) - 1)) <= abs x = 1242 inv (abs x) <= inv (2 pow j * inv (2 pow (bias (:'w) - 1))))`, 1243 REPEAT strip_tac 1244 \\ match_mp_tac (GSYM inv_le) 1245 \\ asm_simp_tac std_ss [REAL_ARITH ``x <> 0 ==> 0 < abs x``] 1246 \\ match_mp_tac realTheory.REAL_LT_MUL 1247 \\ simp_tac std_ss [realTheory.REAL_POW_LT, realTheory.REAL_LT_INV_EQ, 1248 REAL_ARITH ``0 < 2r``] 1249 ); 1250 1251val inv_mul = Q.prove( 1252 `!a b. a <> 0 /\ b <> 0 ==> (inv (a * inv b) = b / a)`, 1253 rw [realTheory.REAL_INV_MUL, realTheory.REAL_INV_NZ, realTheory.REAL_INV_INV] 1254 \\ simp [realTheory.REAL_INV_1OVER, realTheory.mult_ratl] 1255 ); 1256 1257val relative_error_zero = Q.prove( 1258 `!x. (x = 0) ==> 1259 ?e. abs e <= 1 / 2 pow (dimindex (:'t) + 1) /\ 1260 (float_to_real (round roundTiesToEven x : ('t, 'w) float) = 1261 x * (1 + e))`, 1262 rw [] 1263 \\ qexists_tac `0` 1264 \\ qspec_then `0` 1265 (fn th => simp [REWRITE_RULE [realTheory.REAL_SUB_RZERO] th]) 1266 (GSYM error_def) 1267 \\ match_mp_tac error_is_zero 1268 \\ qexists_tac `float_plus_zero (:'t # 'w)` 1269 \\ simp [binary_ieeeTheory.zero_to_real, binary_ieeeTheory.zero_properties] 1270 ); 1271 1272val relative_error = Q.store_thm ("relative_error", 1273 `!x. normalizes (:'t # 'w) x ==> 1274 ?e. abs e <= 1 / 2 pow (dimindex (:'t) + 1) /\ 1275 (float_to_real (round roundTiesToEven x : ('t, 'w) float) = 1276 x * (1 + e))`, 1277 rpt strip_tac 1278 \\ Cases_on `x = 0r` 1279 >- (match_mp_tac relative_error_zero \\ simp []) 1280 \\ imp_res_tac bias_imp_dimindex 1281 \\ `x < 0r \/ 0 < x` by (qpat_assum `x <> 0` MP_TAC \\ REAL_ARITH_TAC) 1282 \\ (qspec_then `x` mp_tac real_in_binade 1283 \\ rw_tac std_ss [] 1284 \\ full_simp_tac std_ss [normalizes_def] 1285 \\ qspecl_then [`x`,`j`] MP_TAC error_bound_norm_strong 1286 \\ rw_tac std_ss [] 1287 \\ `2 pow j * inv (2 pow (bias (:'w) - 1)) <= abs x = 1288 inv (abs x) <= inv (2 pow j * inv (2 pow (bias (:'w) - 1)))` 1289 by metis_tac [relative_bound_lem] 1290 \\ Q.UNDISCH_TAC `2 pow j / 2 pow (bias (:'w) - 1) <= abs x` 1291 \\ asm_simp_tac std_ss [real_div] 1292 \\ rw_tac std_ss [] 1293 \\ `0 <= inv (abs x)` by simp [REAL_LE_INV_EQ, ABS_POS] 1294 \\ qspec_then `error (:'t # 'w) x` assume_tac ABS_POS 1295 \\ qspecl_then 1296 [`abs (error (:'t # 'w) x)`, 1297 `2 pow j / 2 pow (bias (:'w) + dimindex (:'t))`, 1298 `inv (abs x)`, 1299 `inv (2 pow j * inv (2 pow (bias (:'w) - 1)))`] mp_tac REAL_LE_MUL2 1300 \\ asm_simp_tac std_ss [] 1301 \\ strip_tac 1302 \\ qexists_tac `error (:'t # 'w) x * inv x` 1303 \\ conj_tac 1304 >- (simp_tac std_ss [realTheory.ABS_MUL, realTheory.REAL_MUL_LID] 1305 \\ match_mp_tac realTheory.REAL_LE_TRANS 1306 \\ qexists_tac `2 pow j / 2 pow (bias (:'w) + dimindex (:'t)) * 1307 inv (2 pow j * inv (2 pow (bias (:'w) - 1)))` 1308 \\ asm_simp_tac std_ss [realTheory.ABS_INV] 1309 \\ simp_tac std_ss 1310 [inv_mul, realTheory.POW_NZ, REAL_ARITH ``2 <> 0r``, 1311 realTheory.REAL_POS_NZ, realTheory.REAL_INV_NZ, 1312 realTheory.REAL_DIV_OUTER_CANCEL] 1313 \\ simp [realTheory.REAL_INV_1OVER, realTheory.mult_ratl, 1314 realTheory.REAL_LE_LDIV_EQ, realTheory.REAL_LE_RDIV_EQ] 1315 \\ simp [GSYM POW_ADD] 1316 \\ EVAL_TAC 1317 ) 1318 \\ asm_simp_tac std_ss 1319 [error_def, REAL_LDISTRIB, REAL_MUL_RID, REAL_MUL_RINV, 1320 REAL_SUB_LDISTRIB, REAL_SUB_RDISTRIB, REAL_MUL_LID, REAL_SUB_ADD2, 1321 REAL_ARITH ``x * (float_to_real qq * inv x) = 1322 (x * inv x) * float_to_real qq``] 1323 ) 1324 ); 1325 1326(* ------------------------------------------------------------------------- 1327 Ensure that the result is actually finite. 1328 ------------------------------------------------------------------------- *) 1329 1330val float_round_finite = Q.store_thm ("float_round_finite", 1331 `!b x. abs x < threshold (:'t # 'w) ==> 1332 float_is_finite (float_round roundTiesToEven b x : ('t, 'w) float)`, 1333 rw [float_round_def, round_def, binary_ieeeTheory.zero_properties, 1334 REAL_ARITH ``abs x < y = ~(x <= ~y) /\ ~(x >= y)``, 1335 REWRITE_RULE [pred_setTheory.GSPEC_ETA] is_finite_closest] 1336 ); 1337 1338val float_value_finite = Q.prove( 1339 `!a. float_is_finite a ==> (float_value a = Float (float_to_real a))`, 1340 Cases 1341 \\ rename [`float s e f`] 1342 \\ `float s e f = <| Sign := s; Exponent := e; Significand := f |>` 1343 by simp [float_component_equality] 1344 \\ simp [binary_ieeeTheory.float_tests, float_value_def] 1345 ); 1346 1347(* ------------------------------------------------------------------------- 1348 Lifting of arithmetic operations. 1349 ------------------------------------------------------------------------- *) 1350 1351val finite_not = Q.prove( 1352 `!a. float_is_finite a ==> ~float_is_infinite a /\ ~float_is_nan a`, 1353 strip_tac 1354 \\ Cases_on `float_value a` 1355 \\ simp [float_is_finite_def, float_is_infinite_def, float_is_nan_def] 1356 ) 1357 1358val zero_le_ulp = Q.prove( 1359 `0 <= ulp (:'t # 'w)`, 1360 simp [ulp_def, ULP_def]) 1361 1362val round_zero = 1363 binary_ieeeTheory.round_roundTiesToEven_is_zero 1364 |> Q.SPEC `0` 1365 |> SIMP_RULE (srw_ss()) [zero_le_ulp] 1366 1367val lift_tac = 1368 rpt gen_tac 1369 \\ strip_tac 1370 \\ full_simp_tac (srw_ss()++realSimps.real_SS++boolSimps.LET_ss) 1371 [float_value_finite, error_def, float_round_finite, normalizes_def, 1372 float_add_def, float_sub_def, float_mul_def, float_div_def, 1373 float_sqrt_def, float_mul_add_def, float_mul_sub_def, 1374 binary_ieeeTheory.float_is_zero_to_real, float_round_with_flags_def] 1375 \\ rw [float_round_finite, finite_not, 1376 binary_ieeeTheory.float_is_zero_to_real, 1377 binary_ieeeTheory.zero_to_real, binary_ieeeTheory.zero_properties] 1378 \\ rw [float_round_def, finite_not, binary_ieeeTheory.float_is_zero_to_real, 1379 binary_ieeeTheory.zero_to_real, binary_ieeeTheory.zero_properties] 1380 1381val float_add = Q.store_thm ("float_add", 1382 `!a b : ('t, 'w) float. 1383 float_is_finite a /\ float_is_finite b /\ 1384 abs (float_to_real a + float_to_real b) < threshold (:'t # 'w) ==> 1385 float_is_finite (SND (float_add roundTiesToEven a b)) /\ 1386 (float_to_real (SND (float_add roundTiesToEven a b)) = 1387 float_to_real a + float_to_real b + 1388 error (:'t # 'w) (float_to_real a + float_to_real b))`, 1389 lift_tac 1390 ) 1391 1392val float_sub = Q.store_thm ("float_sub", 1393 `!a b : ('t, 'w) float. 1394 float_is_finite a /\ float_is_finite b /\ 1395 abs (float_to_real a - float_to_real b) < threshold (:'t # 'w) ==> 1396 float_is_finite (SND (float_sub roundTiesToEven a b)) /\ 1397 (float_to_real (SND (float_sub roundTiesToEven a b)) = 1398 float_to_real a - float_to_real b + 1399 error (:'t # 'w) (float_to_real a - float_to_real b))`, 1400 lift_tac 1401 ); 1402 1403val float_mul = Q.store_thm ("float_mul", 1404 `!a b : ('t, 'w) float. 1405 float_is_finite a /\ float_is_finite b /\ 1406 abs (float_to_real a * float_to_real b) < threshold (:'t # 'w) ==> 1407 float_is_finite (SND (float_mul roundTiesToEven a b)) /\ 1408 (float_to_real (SND (float_mul roundTiesToEven a b)) = 1409 float_to_real a * float_to_real b + 1410 error (:'t # 'w) (float_to_real a * float_to_real b))`, 1411 lift_tac) 1412 1413val float_div = Q.store_thm ("float_div", 1414 `!a b : ('t, 'w) float. 1415 float_is_finite a /\ float_is_finite b /\ ~float_is_zero b /\ 1416 abs (float_to_real a / float_to_real b) < threshold (:'t # 'w) ==> 1417 float_is_finite (SND (float_div roundTiesToEven a b)) /\ 1418 (float_to_real (SND (float_div roundTiesToEven a b)) = 1419 float_to_real a / float_to_real b + 1420 error (:'t # 'w) (float_to_real a / float_to_real b))`, 1421 lift_tac) 1422 1423val float_sqrt = Q.store_thm ("float_sqrt", 1424 `!a : ('t, 'w) float. 1425 float_is_finite a /\ (a.Sign = 0w) /\ 1426 abs (sqrt (float_to_real a)) < threshold (:'t # 'w) ==> 1427 float_is_finite (SND (float_sqrt roundTiesToEven a)) /\ 1428 (float_to_real (SND (float_sqrt roundTiesToEven a)) = 1429 sqrt (float_to_real a) + error (:'t # 'w) (sqrt (float_to_real a)))`, 1430 lift_tac) 1431 1432val float_mul_add = Q.store_thm ("float_mul_add", 1433 `!a b c : ('t, 'w) float. 1434 float_is_finite a /\ float_is_finite b /\ float_is_finite c /\ 1435 abs (float_to_real a * float_to_real b + float_to_real c) < 1436 threshold (:'t # 'w) ==> 1437 float_is_finite (SND (float_mul_add roundTiesToEven a b c)) /\ 1438 (float_to_real (SND (float_mul_add roundTiesToEven a b c)) = 1439 float_to_real a * float_to_real b + float_to_real c + 1440 error (:'t # 'w) (float_to_real a * float_to_real b + float_to_real c))`, 1441 lift_tac 1442 ) 1443 1444val float_mul_sub = Q.store_thm ("float_mul_add", 1445 `!a b c : ('t, 'w) float. 1446 float_is_finite a /\ float_is_finite b /\ float_is_finite c /\ 1447 abs (float_to_real a * float_to_real b - float_to_real c) < 1448 threshold (:'t # 'w) ==> 1449 float_is_finite (SND (float_mul_sub roundTiesToEven a b c)) /\ 1450 (float_to_real (SND (float_mul_sub roundTiesToEven a b c)) = 1451 float_to_real a * float_to_real b - float_to_real c + 1452 error (:'t # 'w) (float_to_real a * float_to_real b - float_to_real c))`, 1453 lift_tac) 1454 1455(*-----------------------*) 1456 1457fun try_gen q th = Q.GEN q th handle HOL_ERR _ => th 1458 1459val finite_rule = 1460 Q.GEN `a` o try_gen `b` o try_gen `c` o 1461 MATCH_MP (DECIDE ``(X ==> a /\ b) ==> (X ==> a)``) o 1462 Drule.SPEC_ALL 1463 1464val float_add_finite = save_thm ("float_add_finite", finite_rule float_add) 1465val float_sub_finite = save_thm ("float_sub_finite", finite_rule float_sub) 1466val float_mul_finite = save_thm ("float_mul_finite", finite_rule float_mul) 1467val float_div_finite = save_thm ("float_div_finite", finite_rule float_div) 1468val float_sqrt_finite = save_thm ("float_sqrt_finite", finite_rule float_sqrt) 1469 1470val float_mul_add_finite = save_thm ("float_mul_add_finite", 1471 finite_rule float_mul_add) 1472 1473val float_mul_sub_finite = save_thm ("float_mul_sub_finite", 1474 finite_rule float_mul_sub) 1475 1476(*-----------------------*) 1477 1478val relative_tac = 1479 rpt gen_tac 1480 \\ strip_tac 1481 \\ conj_tac 1482 >- fs [normalizes_def, float_add_finite, float_sub_finite, float_mul_finite, 1483 float_div_finite, float_sqrt_finite, float_mul_add_finite, 1484 float_mul_sub_finite] 1485 \\ imp_res_tac relative_error 1486 \\ qexists_tac `e` 1487 \\ full_simp_tac (srw_ss()++realSimps.real_SS++boolSimps.LET_ss) 1488 [float_value_finite, error_def, float_round_finite, normalizes_def, 1489 float_add_def, float_sub_def, float_mul_def, float_div_def, 1490 float_sqrt_def, float_mul_add_def, float_mul_sub_def, 1491 binary_ieeeTheory.float_is_zero_to_real, 1492 float_round_with_flags_def] 1493 \\ rw [float_round_def, binary_ieeeTheory.float_is_zero_to_real, finite_not, 1494 binary_ieeeTheory.zero_to_real, binary_ieeeTheory.zero_properties] 1495 \\ rw [real_to_float_def, float_round_def, finite_not, 1496 binary_ieeeTheory.float_is_zero_to_real, 1497 binary_ieeeTheory.zero_to_real, binary_ieeeTheory.zero_properties] 1498 1499val float_add_relative = Q.store_thm ("float_add_relative", 1500 `!a b : ('t, 'w) float. 1501 float_is_finite a /\ float_is_finite b /\ 1502 normalizes (:'t # 'w) (float_to_real a + float_to_real b) ==> 1503 float_is_finite (SND (float_add roundTiesToEven a b)) /\ 1504 ?e. abs e <= 1 / 2 pow (dimindex (:'t) + 1) /\ 1505 (float_to_real (SND (float_add roundTiesToEven a b)) = 1506 (float_to_real a + float_to_real b) * (1 + e))`, 1507 relative_tac 1508 ); 1509 1510val float_sub_relative = Q.store_thm ("float_sub_relative", 1511 `!a b : ('t, 'w) float. 1512 float_is_finite a /\ float_is_finite b /\ 1513 normalizes (:'t # 'w) (float_to_real a - float_to_real b) ==> 1514 float_is_finite (SND (float_sub roundTiesToEven a b)) /\ 1515 ?e. abs e <= 1 / 2 pow (dimindex (:'t) + 1) /\ 1516 (float_to_real (SND (float_sub roundTiesToEven a b)) = 1517 (float_to_real a - float_to_real b) * (1 + e))`, 1518 relative_tac 1519 ); 1520 1521val float_mul_relative = Q.store_thm ("float_mul_relative", 1522 `!a b : ('t, 'w) float. 1523 float_is_finite a /\ float_is_finite b /\ 1524 normalizes (:'t # 'w) (float_to_real a * float_to_real b) ==> 1525 float_is_finite (SND (float_mul roundTiesToEven a b)) /\ 1526 ?e. abs e <= 1 / 2 pow (dimindex (:'t) + 1) /\ 1527 (float_to_real (SND (float_mul roundTiesToEven a b)) = 1528 (float_to_real a * float_to_real b) * (1 + e))`, 1529 relative_tac 1530 ); 1531 1532val float_div_relative = Q.store_thm ("float_div_relative", 1533 `!a b : ('t, 'w) float. 1534 float_is_finite a /\ float_is_finite b /\ ~float_is_zero b /\ 1535 normalizes (:'t # 'w) (float_to_real a / float_to_real b) ==> 1536 float_is_finite (SND (float_div roundTiesToEven a b)) /\ 1537 ?e. abs e <= 1 / 2 pow (dimindex (:'t) + 1) /\ 1538 (float_to_real (SND (float_div roundTiesToEven a b)) = 1539 (float_to_real a / float_to_real b) * (1 + e))`, 1540 relative_tac 1541 ); 1542 1543val float_sqrt_relative = Q.store_thm ("float_sqrt_relative", 1544 `!a : ('t, 'w) float. 1545 float_is_finite a /\ (a.Sign = 0w) /\ 1546 normalizes (:'t # 'w) (sqrt (float_to_real a)) ==> 1547 float_is_finite (SND (float_sqrt roundTiesToEven a)) /\ 1548 ?e. abs e <= 1 / 2 pow (dimindex (:'t) + 1) /\ 1549 (float_to_real (SND (float_sqrt roundTiesToEven a)) = 1550 (sqrt (float_to_real a) * (1 + e)))`, 1551 relative_tac 1552 ); 1553 1554val float_mul_add_relative = Q.store_thm ("float_mul_add_relative", 1555 `!a b c : ('t, 'w) float. 1556 float_is_finite a /\ float_is_finite b /\ float_is_finite c /\ 1557 normalizes (:'t # 'w) 1558 (float_to_real a * float_to_real b + float_to_real c) ==> 1559 float_is_finite (SND (float_mul_add roundTiesToEven a b c)) /\ 1560 ?e. abs e <= 1 / 2 pow (dimindex (:'t) + 1) /\ 1561 (float_to_real (SND (float_mul_add roundTiesToEven a b c)) = 1562 (float_to_real a * float_to_real b + float_to_real c) * (1 + e))`, 1563 relative_tac 1564 ); 1565 1566val float_mul_sub_relative = Q.store_thm ("float_mul_sub_relative", 1567 `!a b c : ('t, 'w) float. 1568 float_is_finite a /\ float_is_finite b /\ float_is_finite c /\ 1569 normalizes (:'t # 'w) 1570 (float_to_real a * float_to_real b - float_to_real c) ==> 1571 float_is_finite (SND (float_mul_sub roundTiesToEven a b c)) /\ 1572 ?e. abs e <= 1 / 2 pow (dimindex (:'t) + 1) /\ 1573 (float_to_real (SND (float_mul_sub roundTiesToEven a b c)) = 1574 (float_to_real a * float_to_real b - float_to_real c) * (1 + e))`, 1575 relative_tac 1576 ); 1577 1578(* ------------------------------------------------------------------------- *) 1579 1580val finite_float_within_threshold = Q.store_thm ( 1581 "finite_float_within_threshold", 1582 `!f:('a , 'b) float. 1583 float_is_finite f ==> 1584 ~(float_to_real f <= -threshold (:'a # 'b)) /\ 1585 ~(float_to_real f >= threshold (:'a # 'b)) `, 1586 rpt strip_tac 1587 \\ Q.ISPECL_THEN [`f`] assume_tac float_to_real_threshold 1588 \\ fs[realTheory.abs] 1589 \\ BasicProvers.every_case_tac 1590 \\ res_tac 1591 \\ RealArith.REAL_ASM_ARITH_TAC); 1592 1593val round_finite_float_id = Q.store_thm( 1594"round_finite_normal_float_id", 1595 `!f. 1596 float_is_finite f /\ 1597 ~ float_is_zero f ==> 1598 (round roundTiesToEven (float_to_real f) = f)`, 1599 rw[] 1600 \\ qpat_assum `float_is_finite _` mp_tac 1601 \\ rewrite_tac [float_is_finite_def, float_value_def] 1602 \\ simp[] 1603 \\ strip_tac 1604 \\ once_rewrite_tac [round_def] 1605 \\ fs[finite_float_within_threshold] 1606 \\ once_rewrite_tac [closest_such_def] 1607 \\ SELECT_ELIM_TAC 1608 \\ rw[] 1609 >- (qexists_tac `f` 1610 \\ rw[is_closest_def, IN_DEF, realTheory.ABS_POS] 1611 \\ Cases_on `f = b` \\ fs[] 1612 \\ first_x_assum (qspec_then `f` mp_tac) 1613 \\ fs[realTheory.REAL_SUB_REFL] 1614 \\ strip_tac 1615 \\ `float_to_real b - float_to_real f = 0` 1616 by (RealArith.REAL_ASM_ARITH_TAC) 1617 \\ fs[float_to_real_eq] 1618 \\ rfs[]) 1619 \\ CCONTR_TAC 1620 \\ fs[is_closest_def, IN_DEF] 1621 \\ qpat_x_assum `!x._ ` mp_tac 1622 \\ first_x_assum (qspec_then `f` mp_tac) 1623 \\ fs[realTheory.REAL_SUB_REFL] 1624 \\ rpt strip_tac 1625 \\ `float_to_real x - float_to_real f = 0` 1626 by (RealArith.REAL_ASM_ARITH_TAC) 1627 \\ fs[float_to_real_eq] 1628 \\ rfs[]); 1629 1630val real_to_float_finite_id = Q.store_thm ( 1631 "real_to_float_finite_normal_id", 1632 `!f. 1633 float_is_finite f /\ 1634 ~ float_is_zero f ==> 1635 (real_to_float roundTiesToEven (float_to_real f) = f)`, 1636 rpt strip_tac 1637 \\ fs[real_to_float_def, float_round_def, round_finite_float_id]); 1638 1639val float_to_real_real_to_float_zero_id = Q.store_thm ( 1640 "float_to_real_real_to_float_zero_id", 1641 `float_to_real (real_to_float roundTiesToEven 0) = 0`, 1642 once_rewrite_tac[real_to_float_def] 1643 \\ `float_round roundTiesToEven F 0 = (float_plus_zero(:'a # 'b))` 1644 by (irule round_roundTiesToEven_is_plus_zero 1645 \\ fs[ulp_def, ULP_def]) 1646 \\ fs[float_to_real_def, float_plus_zero_def]); 1647 1648val non_representable_float_is_zero = store_thm ( 1649 "non_representable_float_is_zero", 1650 ``!ff P. 1651 2 * abs ff <= ulp ((:'a#'b) :('a#'b) itself) ==> 1652 (float_to_real ((float_round roundTiesToEven P ff):('a, 'b) float) = 0)``, 1653 rpt strip_tac \\ Cases_on `P` 1654 \\ fs [round_roundTiesToEven_is_plus_zero, 1655 round_roundTiesToEven_is_minus_zero, zero_to_real]); 1656 1657val () = export_theory () 1658