diff options
Diffstat (limited to 'flocq/Core/FLT.v')
-rw-r--r-- | flocq/Core/FLT.v | 406 |
1 files changed, 406 insertions, 0 deletions
diff --git a/flocq/Core/FLT.v b/flocq/Core/FLT.v new file mode 100644 index 00000000..bd48d4b7 --- /dev/null +++ b/flocq/Core/FLT.v @@ -0,0 +1,406 @@ +(** +This file is part of the Flocq formalization of floating-point +arithmetic in Coq: http://flocq.gforge.inria.fr/ + +Copyright (C) 2009-2018 Sylvie Boldo +#<br /># +Copyright (C) 2009-2018 Guillaume Melquiond + +This library is free software; you can redistribute it and/or +modify it under the terms of the GNU Lesser General Public +License as published by the Free Software Foundation; either +version 3 of the License, or (at your option) any later version. + +This library is distributed in the hope that it will be useful, +but WITHOUT ANY WARRANTY; without even the implied warranty of +MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +COPYING file for more details. +*) + +(** * Floating-point format with gradual underflow *) +Require Import Raux Defs Round_pred Generic_fmt Float_prop. +Require Import FLX FIX Ulp Round_NE. +Require Import Psatz. + +Section RND_FLT. + +Variable beta : radix. + +Notation bpow e := (bpow beta e). + +Variable emin prec : Z. + +Context { prec_gt_0_ : Prec_gt_0 prec }. + +Inductive FLT_format (x : R) : Prop := + FLT_spec (f : float beta) : + x = F2R f -> (Z.abs (Fnum f) < Zpower beta prec)%Z -> + (emin <= Fexp f)%Z -> FLT_format x. + +Definition FLT_exp e := Z.max (e - prec) emin. + +(** Properties of the FLT format *) +Global Instance FLT_exp_valid : Valid_exp FLT_exp. +Proof. +intros k. +unfold FLT_exp. +generalize (prec_gt_0 prec). +repeat split ; + intros ; zify ; omega. +Qed. + +Theorem generic_format_FLT : + forall x, FLT_format x -> generic_format beta FLT_exp x. +Proof. +clear prec_gt_0_. +intros x [[mx ex] H1 H2 H3]. +simpl in H2, H3. +rewrite H1. +apply generic_format_F2R. +intros Zmx. +unfold cexp, FLT_exp. +rewrite mag_F2R with (1 := Zmx). +apply Z.max_lub with (2 := H3). +apply Zplus_le_reg_r with (prec - ex)%Z. +ring_simplify. +now apply mag_le_Zpower. +Qed. + +Theorem FLT_format_generic : + forall x, generic_format beta FLT_exp x -> FLT_format x. +Proof. +intros x. +unfold generic_format. +set (ex := cexp beta FLT_exp x). +set (mx := Ztrunc (scaled_mantissa beta FLT_exp x)). +intros Hx. +rewrite Hx. +eexists ; repeat split ; simpl. +apply lt_IZR. +rewrite IZR_Zpower. 2: now apply Zlt_le_weak. +apply Rmult_lt_reg_r with (bpow ex). +apply bpow_gt_0. +rewrite <- bpow_plus. +change (F2R (Float beta (Z.abs mx) ex) < bpow (prec + ex))%R. +rewrite F2R_Zabs. +rewrite <- Hx. +destruct (Req_dec x 0) as [Hx0|Hx0]. +rewrite Hx0, Rabs_R0. +apply bpow_gt_0. +unfold cexp in ex. +destruct (mag beta x) as (ex', He). +simpl in ex. +specialize (He Hx0). +apply Rlt_le_trans with (1 := proj2 He). +apply bpow_le. +cut (ex' - prec <= ex)%Z. omega. +unfold ex, FLT_exp. +apply Z.le_max_l. +apply Z.le_max_r. +Qed. + + +Theorem FLT_format_bpow : + forall e, (emin <= e)%Z -> generic_format beta FLT_exp (bpow e). +Proof. +intros e He. +apply generic_format_bpow; unfold FLT_exp. +apply Z.max_case; try assumption. +unfold Prec_gt_0 in prec_gt_0_; omega. +Qed. + + + + +Theorem FLT_format_satisfies_any : + satisfies_any FLT_format. +Proof. +refine (satisfies_any_eq _ _ _ (generic_format_satisfies_any beta FLT_exp)). +intros x. +split. +apply FLT_format_generic. +apply generic_format_FLT. +Qed. + +Theorem cexp_FLT_FLX : + forall x, + (bpow (emin + prec - 1) <= Rabs x)%R -> + cexp beta FLT_exp x = cexp beta (FLX_exp prec) x. +Proof. +intros x Hx. +assert (Hx0: x <> 0%R). +intros H1; rewrite H1, Rabs_R0 in Hx. +contradict Hx; apply Rlt_not_le, bpow_gt_0. +unfold cexp. +apply Zmax_left. +destruct (mag beta x) as (ex, He). +unfold FLX_exp. simpl. +specialize (He Hx0). +cut (emin + prec - 1 < ex)%Z. omega. +apply (lt_bpow beta). +apply Rle_lt_trans with (1 := Hx). +apply He. +Qed. + +(** Links between FLT and FLX *) +Theorem generic_format_FLT_FLX : + forall x : R, + (bpow (emin + prec - 1) <= Rabs x)%R -> + generic_format beta (FLX_exp prec) x -> + generic_format beta FLT_exp x. +Proof. +intros x Hx H. +destruct (Req_dec x 0) as [Hx0|Hx0]. +rewrite Hx0. +apply generic_format_0. +unfold generic_format, scaled_mantissa. +now rewrite cexp_FLT_FLX. +Qed. + +Theorem generic_format_FLX_FLT : + forall x : R, + generic_format beta FLT_exp x -> generic_format beta (FLX_exp prec) x. +Proof. +clear prec_gt_0_. +intros x Hx. +unfold generic_format in Hx; rewrite Hx. +apply generic_format_F2R. +intros _. +rewrite <- Hx. +unfold cexp, FLX_exp, FLT_exp. +apply Z.le_max_l. +Qed. + +Theorem round_FLT_FLX : forall rnd x, + (bpow (emin + prec - 1) <= Rabs x)%R -> + round beta FLT_exp rnd x = round beta (FLX_exp prec) rnd x. +Proof. +intros rnd x Hx. +unfold round, scaled_mantissa. +rewrite cexp_FLT_FLX ; trivial. +Qed. + +(** Links between FLT and FIX (underflow) *) +Theorem cexp_FLT_FIX : + forall x, x <> 0%R -> + (Rabs x < bpow (emin + prec))%R -> + cexp beta FLT_exp x = cexp beta (FIX_exp emin) x. +Proof. +intros x Hx0 Hx. +unfold cexp. +apply Zmax_right. +unfold FIX_exp. +destruct (mag beta x) as (ex, Hex). +simpl. +cut (ex - 1 < emin + prec)%Z. omega. +apply (lt_bpow beta). +apply Rle_lt_trans with (2 := Hx). +now apply Hex. +Qed. + +Theorem generic_format_FIX_FLT : + forall x : R, + generic_format beta FLT_exp x -> + generic_format beta (FIX_exp emin) x. +Proof. +clear prec_gt_0_. +intros x Hx. +rewrite Hx. +apply generic_format_F2R. +intros _. +rewrite <- Hx. +apply Z.le_max_r. +Qed. + +Theorem generic_format_FLT_FIX : + forall x : R, + (Rabs x <= bpow (emin + prec))%R -> + generic_format beta (FIX_exp emin) x -> + generic_format beta FLT_exp x. +Proof with auto with typeclass_instances. +apply generic_inclusion_le... +intros e He. +unfold FIX_exp. +apply Z.max_lub. +omega. +apply Z.le_refl. +Qed. + +Lemma negligible_exp_FLT : + exists n, negligible_exp FLT_exp = Some n /\ (n <= emin)%Z. +Proof. +case (negligible_exp_spec FLT_exp). +{ intro H; exfalso; specialize (H emin); revert H. + apply Zle_not_lt, Z.le_max_r. } +intros n Hn; exists n; split; [now simpl|]. +destruct (Z.max_spec (n - prec) emin) as [(Hm, Hm')|(Hm, Hm')]. +{ now revert Hn; unfold FLT_exp; rewrite Hm'. } +revert Hn prec_gt_0_; unfold FLT_exp, Prec_gt_0; rewrite Hm'; lia. +Qed. + +Theorem generic_format_FLT_1 (Hemin : (emin <= 0)%Z) : + generic_format beta FLT_exp 1. +Proof. +unfold generic_format, scaled_mantissa, cexp, F2R; simpl. +rewrite Rmult_1_l, (mag_unique beta 1 1). +{ unfold FLT_exp. + destruct (Z.max_spec_le (1 - prec) emin) as [(H,Hm)|(H,Hm)]; rewrite Hm; + (rewrite <- IZR_Zpower; [|unfold Prec_gt_0 in prec_gt_0_; omega]); + (rewrite Ztrunc_IZR, IZR_Zpower, <-bpow_plus; + [|unfold Prec_gt_0 in prec_gt_0_; omega]); + now replace (_ + _)%Z with Z0 by ring. } +rewrite Rabs_R1; simpl; split; [now right|]. +rewrite IZR_Zpower_pos; simpl; rewrite Rmult_1_r; apply IZR_lt. +apply (Z.lt_le_trans _ 2); [omega|]; apply Zle_bool_imp_le, beta. +Qed. + +Theorem ulp_FLT_small: forall x, (Rabs x < bpow (emin+prec))%R -> + ulp beta FLT_exp x = bpow emin. +Proof with auto with typeclass_instances. +intros x Hx. +unfold ulp; case Req_bool_spec; intros Hx2. +(* x = 0 *) +case (negligible_exp_spec FLT_exp). +intros T; specialize (T (emin-1)%Z); contradict T. +apply Zle_not_lt; unfold FLT_exp. +apply Z.le_trans with (2:=Z.le_max_r _ _); omega. +assert (V:FLT_exp emin = emin). +unfold FLT_exp; apply Z.max_r. +unfold Prec_gt_0 in prec_gt_0_; omega. +intros n H2; rewrite <-V. +apply f_equal, fexp_negligible_exp_eq... +omega. +(* x <> 0 *) +apply f_equal; unfold cexp, FLT_exp. +apply Z.max_r. +assert (mag beta x-1 < emin+prec)%Z;[idtac|omega]. +destruct (mag beta x) as (e,He); simpl. +apply lt_bpow with beta. +apply Rle_lt_trans with (2:=Hx). +now apply He. +Qed. + +Theorem ulp_FLT_le : + forall x, (bpow (emin + prec - 1) <= Rabs x)%R -> + (ulp beta FLT_exp x <= Rabs x * bpow (1 - prec))%R. +Proof. +intros x Hx. +assert (Zx : (x <> 0)%R). + intros Z; contradict Hx; apply Rgt_not_le, Rlt_gt. + rewrite Z, Rabs_R0; apply bpow_gt_0. +rewrite ulp_neq_0 with (1 := Zx). +unfold cexp, FLT_exp. +destruct (mag beta x) as (e,He). +apply Rle_trans with (bpow (e-1)*bpow (1-prec))%R. +rewrite <- bpow_plus. +right; apply f_equal. +replace (e - 1 + (1 - prec))%Z with (e - prec)%Z by ring. +apply Z.max_l. +assert (emin+prec-1 < e)%Z; try omega. +apply lt_bpow with beta. +apply Rle_lt_trans with (1:=Hx). +now apply He. +apply Rmult_le_compat_r. +apply bpow_ge_0. +now apply He. +Qed. + +Theorem ulp_FLT_gt : + forall x, (Rabs x * bpow (-prec) < ulp beta FLT_exp x)%R. +Proof. +intros x; case (Req_dec x 0); intros Hx. +rewrite Hx, ulp_FLT_small, Rabs_R0, Rmult_0_l; try apply bpow_gt_0. +rewrite Rabs_R0; apply bpow_gt_0. +rewrite ulp_neq_0; try exact Hx. +unfold cexp, FLT_exp. +apply Rlt_le_trans with (bpow (mag beta x)*bpow (-prec))%R. +apply Rmult_lt_compat_r. +apply bpow_gt_0. +now apply bpow_mag_gt. +rewrite <- bpow_plus. +apply bpow_le. +apply Z.le_max_l. +Qed. + +Lemma ulp_FLT_exact_shift : + forall x e, + (x <> 0)%R -> + (emin + prec <= mag beta x)%Z -> + (emin + prec - mag beta x <= e)%Z -> + (ulp beta FLT_exp (x * bpow e) = ulp beta FLT_exp x * bpow e)%R. +Proof. +intros x e Nzx Hmx He. +unfold ulp; rewrite Req_bool_false; + [|now intro H; apply Nzx, (Rmult_eq_reg_r (bpow e)); + [rewrite Rmult_0_l|apply Rgt_not_eq, Rlt_gt, bpow_gt_0]]. +rewrite (Req_bool_false _ _ Nzx), <- bpow_plus; f_equal; unfold cexp, FLT_exp. +rewrite (mag_mult_bpow _ _ _ Nzx), !Z.max_l; omega. +Qed. + +Lemma succ_FLT_exact_shift_pos : + forall x e, + (0 < x)%R -> + (emin + prec <= mag beta x)%Z -> + (emin + prec - mag beta x <= e)%Z -> + (succ beta FLT_exp (x * bpow e) = succ beta FLT_exp x * bpow e)%R. +Proof. +intros x e Px Hmx He. +rewrite succ_eq_pos; [|now apply Rlt_le, Rmult_lt_0_compat, bpow_gt_0]. +rewrite (succ_eq_pos _ _ _ (Rlt_le _ _ Px)). +now rewrite Rmult_plus_distr_r; f_equal; apply ulp_FLT_exact_shift; [lra| |]. +Qed. + +Lemma succ_FLT_exact_shift : + forall x e, + (x <> 0)%R -> + (emin + prec + 1 <= mag beta x)%Z -> + (emin + prec - mag beta x + 1 <= e)%Z -> + (succ beta FLT_exp (x * bpow e) = succ beta FLT_exp x * bpow e)%R. +Proof. +intros x e Nzx Hmx He. +destruct (Rle_or_lt 0 x) as [Px|Nx]. +{ now apply succ_FLT_exact_shift_pos; [lra|lia|lia]. } +unfold succ. +rewrite Rle_bool_false; [|assert (H := bpow_gt_0 beta e); nra]. +rewrite Rle_bool_false; [|now simpl]. +rewrite Ropp_mult_distr_l_reverse, <-Ropp_mult_distr_l_reverse; f_equal. +unfold pred_pos. +rewrite mag_mult_bpow; [|lra]. +replace (_ - 1)%Z with (mag beta (- x) - 1 + e)%Z; [|ring]; rewrite bpow_plus. +unfold Req_bool; rewrite Rcompare_mult_r; [|now apply bpow_gt_0]. +fold (Req_bool (-x) (bpow (mag beta (-x) - 1))); case Req_bool. +{ rewrite mag_opp; unfold FLT_exp; do 2 (rewrite Z.max_l; [|lia]). + replace (_ - _)%Z with (mag beta x - 1 - prec + e)%Z; [|ring]. + rewrite bpow_plus; ring. } +rewrite ulp_FLT_exact_shift; [ring|lra| |]; rewrite mag_opp; lia. +Qed. + +(** FLT is a nice format: it has a monotone exponent... *) +Global Instance FLT_exp_monotone : Monotone_exp FLT_exp. +Proof. +intros ex ey. +unfold FLT_exp. +zify ; omega. +Qed. + +(** and it allows a rounding to nearest, ties to even. *) +Hypothesis NE_prop : Z.even beta = false \/ (1 < prec)%Z. + +Global Instance exists_NE_FLT : Exists_NE beta FLT_exp. +Proof. +destruct NE_prop as [H|H]. +now left. +right. +intros e. +unfold FLT_exp. +destruct (Zmax_spec (e - prec) emin) as [(H1,H2)|(H1,H2)] ; + rewrite H2 ; clear H2. +generalize (Zmax_spec (e + 1 - prec) emin). +generalize (Zmax_spec (e - prec + 1 - prec) emin). +omega. +generalize (Zmax_spec (e + 1 - prec) emin). +generalize (Zmax_spec (emin + 1 - prec) emin). +omega. +Qed. + +End RND_FLT. |