Library Float.FnElem.Axpy
Require Export MinOrMax.
Section AxpyMisc.
Let radix := 2%Z.
Let FtoRradix := FtoR radix.
Coercion FtoRradix : float >-> R.
Variable b : Fbound.
Variable precision : nat.
Hypothesis precisionGreaterThanOne : 1 < precision.
Hypothesis pGivesBound : Zpos (vNum b) = Zpower_nat radix precision.
Theorem TwoMoreThanOne : (1 < radix)%Z.
unfold radix in |- *; red in |- *; simpl in |- *; auto.
Qed.
Hint Resolve TwoMoreThanOne.
Theorem TwoMoreThanOneR : (1 < radix)%R.
replace 1%R with (IZR 1); auto with real.
Qed.
Hint Resolve TwoMoreThanOneR.
Theorem FulpLeGeneral :
∀ p : float,
Fbounded b p →
(Fulp b radix precision p ≤
Rabs (FtoRradix p) × powerRZ radix (Zsucc (- precision)) +
powerRZ radix (- dExp b))%R.
intros p Hp.
cut (Fcanonic radix b (Fnormalize radix b precision p));
[ intros H | apply FnormalizeCanonic; auto with arith zarith ].
case H; intros H1.
apply
Rle_trans with (Rabs (FtoR radix p) × powerRZ radix (Zsucc (- precision)))%R.
apply FulpLe2; auto with zarith.
apply Rle_trans with (Rabs p × powerRZ radix (Zsucc (- precision)) + 0)%R;
[ right; fold FtoRradix; ring | apply Rplus_le_compat_l; auto with real zarith ].
apply Rle_trans with (powerRZ radix (- dExp b)).
right; unfold Fulp in |- *; simpl in |- ×.
replace (Fexp (Fnormalize radix b precision p)) with (- dExp b)%Z;
[ auto with real zarith | elim H1; intuition ].
apply Rle_trans with (0 + powerRZ radix (- dExp b))%R;
[ right; ring | apply Rplus_le_compat_r ].
apply Rmult_le_pos; auto with real zarith.
Qed.
Theorem RoundLeGeneral :
∀ (p : float) (z : R),
Fbounded b p →
Closest b radix z p →
(Rabs p ≤
Rabs z × / (1 - powerRZ radix (- precision)) +
powerRZ radix (Zpred (- dExp b)) × / (1 - powerRZ radix (- precision)))%R.
intros p z Hp H.
cut (0 < 1 - powerRZ radix (- precision))%R; [ intros H1 | idtac ].
2: apply Rplus_lt_reg_r with (powerRZ radix (- precision)).
2: ring_simplify.
2: replace 1%R with (powerRZ radix 0);
[ auto with real zarith | simpl in |- *; auto ].
apply Rmult_le_reg_l with (1 - powerRZ radix (- precision))%R; auto.
ring_simplify ((1 - powerRZ radix (- precision)) × Rabs p)%R.
apply
Rle_trans
with
(Rabs z ×
((1 - powerRZ radix (- precision)) × / (1 - powerRZ radix (- precision))) +
powerRZ radix (Zpred (- dExp b)) ×
((1 - powerRZ radix (- precision)) × / (1 - powerRZ radix (- precision))))%R;
[ idtac | right; ring; ring ].
repeat rewrite Rinv_r; auto with real.
ring_simplify.
apply Rplus_le_reg_l with (- powerRZ radix (Zpred (- dExp b)))%R.
ring_simplify.
apply
Rle_trans
with
(Rabs p +
-
(Rabs p × powerRZ radix (- precision) + powerRZ radix (Zpred (- dExp b))))%R;
[ right; ring | idtac ].
apply Rle_trans with (Rabs p + - (/ radix × Fulp b radix precision p))%R.
apply Rplus_le_compat_l; apply Ropp_le_contravar.
apply Rmult_le_reg_l with (IZR radix); auto with real zarith.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real zarith.
ring_simplify (1 × Fulp b radix precision p)%R.
apply
Rle_trans
with
(Rabs p × (powerRZ radix 1 × powerRZ radix (- precision)) +
powerRZ radix 1 × powerRZ radix (Zpred (- dExp b)))%R;
[ idtac | right; simpl in |- *; ring ].
repeat rewrite <- powerRZ_add; auto with zarith real.
replace (1 + - precision)%Z with (Zsucc (- precision)); auto with zarith.
replace (1 + Zpred (- dExp b))%Z with (- dExp b)%Z; auto with zarith.
apply FulpLeGeneral; auto.
unfold Zpred in |- *; auto with zarith.
apply Rplus_le_reg_l with (- Rabs z + / radix × Fulp b radix precision p)%R.
ring_simplify.
apply Rle_trans with (Rabs (p - z)).
apply Rle_trans with (Rabs p - Rabs z)%R;
[ right; ring | apply Rabs_triang_inv ].
rewrite <- Rabs_Ropp.
replace (- (p - z))%R with (z - p)%R; [ idtac | ring ].
apply Rmult_le_reg_l with (IZR radix); auto with real zarith.
apply Rle_trans with (Fulp b radix precision p × (radix × / radix))%R;
[ rewrite Rinv_r; auto with real zarith | right; ring ].
ring_simplify (Fulp b radix precision p × 1)%R.
unfold FtoRradix in |- *; apply ClosestUlp; auto.
Qed.
Theorem RoundGeGeneral :
∀ (p : float) (z : R),
Fbounded b p →
Closest b radix z p →
(Rabs z × / (1 + powerRZ radix (- precision)) -
powerRZ radix (Zpred (- dExp b)) × / (1 + powerRZ radix (- precision)) ≤
Rabs p)%R.
intros p z Hp H.
cut (0 < 1 + powerRZ radix (- precision))%R; [ intros H1 | idtac ].
2: apply Rlt_le_trans with 1%R; auto with real.
2: apply Rle_trans with (1 + 0)%R; auto with real zarith.
apply Rmult_le_reg_l with (1 + powerRZ radix (- precision))%R; auto.
apply
Rle_trans
with
(Rabs z ×
((1 + powerRZ radix (- precision)) × / (1 + powerRZ radix (- precision))) -
powerRZ radix (Zpred (- dExp b)) ×
((1 + powerRZ radix (- precision)) × / (1 + powerRZ radix (- precision))))%R;
[ right; ring; ring | idtac ].
repeat rewrite Rinv_r; auto with real.
ring_simplify.
apply Rplus_le_reg_l with (powerRZ radix (Zpred (- dExp b))).
ring_simplify
(powerRZ radix (Zpred (- dExp b)) +
(Rabs z - powerRZ radix (Zpred (- dExp b))))%R.
apply
Rle_trans
with
(Rabs p +
(Rabs p × powerRZ radix (- precision) + powerRZ radix (Zpred (- dExp b))))%R;
[ idtac | right; ring ].
apply Rle_trans with (Rabs p + / radix × Fulp b radix precision p)%R.
apply Rplus_le_reg_l with (- Rabs p)%R.
ring_simplify (- Rabs p + (Rabs p + / radix × Fulp b radix precision p))%R.
apply Rle_trans with (Rabs (z - p)).
apply Rle_trans with (Rabs z - Rabs p)%R;
[ right; ring | apply Rabs_triang_inv ].
apply Rmult_le_reg_l with (IZR radix); auto with real zarith.
apply Rle_trans with (Fulp b radix precision p × (radix × / radix))%R;
[ rewrite Rinv_r; auto with real zarith | right; ring ].
ring_simplify (Fulp b radix precision p × 1)%R.
unfold FtoRradix in |- *; apply ClosestUlp; auto.
apply Rplus_le_compat_l.
apply Rmult_le_reg_l with (IZR radix); auto with real zarith.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real zarith.
ring_simplify (1 × Fulp b radix precision p)%R.
apply
Rle_trans
with
(Rabs p × (powerRZ radix 1 × powerRZ radix (- precision)) +
powerRZ radix 1 × powerRZ radix (Zpred (- dExp b)))%R;
[ idtac | right; simpl in |- *; ring ].
repeat rewrite <- powerRZ_add; auto with zarith real.
replace (1 + - precision)%Z with (Zsucc (- precision)); auto with zarith.
replace (1 + Zpred (- dExp b))%Z with (- dExp b)%Z; auto with zarith.
apply FulpLeGeneral; auto.
Qed.
Theorem ExactSum_Near :
∀ p q f : float,
Fbounded b p →
Fbounded b q →
Fbounded b f →
Closest b radix (p + q) f →
Fexp p = (- dExp b)%Z →
(Rabs (p + q - f) < powerRZ radix (- dExp b))%R → (p + q)%R = f.
intros p q f Hp Hq Hf H H12 H1.
case
errorBoundedPlus
with
(b := b)
(radix := radix)
(precision := precision)
(p := p)
(q := q)
(pq := f); auto with zarith.
intros x0 H'; elim H'; intros H2 H'1; elim H'1; intros H3 H4; clear H' H'1.
apply Rplus_eq_reg_l with (- FtoRradix f)%R; ring_simplify (-f+f)%R.
apply trans_eq with (FtoR radix p + FtoR radix q - FtoR radix f)%R;
[ fold FtoRradix; ring | rewrite <- H2 ].
generalize H1; unfold FtoRradix in |- *; rewrite <- H2; intros H5.
cut (∀ r : R, Rabs r = 0%R → r = 0%R);
[ intros V; apply V | auto with real ].
rewrite <- Fabs_correct; auto with zarith.
apply is_Fzero_rep1; unfold is_Fzero in |- ×.
cut (∀ z : Z, (0 ≤ z)%Z → (z < 1)%Z → z = 0%Z);
[ intros W; apply W | auto with zarith ].
unfold Fabs in |- *; simpl in |- *; apply Zle_ZERO_Zabs.
replace 1%Z with (Fnum (Float 1 (- dExp b))); [ idtac | simpl in |- *; auto ].
apply Rlt_Fexp_eq_Zlt with (radix := radix); auto with zarith real float.
rewrite Fabs_correct; auto with zarith; apply Rlt_le_trans with (1 := H5).
right; unfold FtoR in |- *; simpl in |- *; ring.
unfold Fabs in |- *; simpl in |- ×.
rewrite H4; rewrite H12; apply Zmin_le1; auto with zarith.
elim Hq; auto with zarith.
intros m; case (Rle_or_lt 0 m); intros H'.
rewrite Rabs_right; auto with real.
rewrite Faux.Rabsolu_left1; auto with real; intros H'1.
rewrite <- (Ropp_involutive m); rewrite H'1; auto with real.
Qed.
End AxpyMisc.
Section AxpyAux.
Add Field RField : Rfield (completeness Zeq_bool_complete).
Let radix := 2%Z.
Let FtoRradix := FtoR radix.
Coercion FtoRradix : float >-> R.
Variable b : Fbound.
Variable precision : nat.
Hypothesis precisionGreaterThanOne : 1 < precision.
Hypothesis pGivesBound : Zpos (vNum b) = Zpower_nat radix precision.
Variables a1 x1 y1 : R.
Variables a x y t u r : float.
Hypothesis Fa : Fbounded b a.
Hypothesis Fx : Fbounded b x.
Hypothesis Fy : Fbounded b y.
Hypothesis Ft : Fbounded b t.
Hypothesis Fu : Fbounded b u.
Hypothesis tDef : Closest b radix (a × x) t.
Hypothesis uDef : Closest b radix (t + y) u.
Hypothesis rDef : isMin b radix (a1 × x1 + y1) r.
Theorem Axpy_aux1 :
Fcanonic radix b u →
(Rabs (FtoRradix a × FtoRradix x - FtoRradix t) ≤
/ 4%nat × Fulp b radix precision (FPred b radix precision u))%R →
(0 < u)%R →
(4%nat × Rabs t ≤ Rabs u)%R →
(Rabs (y1 - y) + Rabs (a1 × x1 - a × x) <
/ 4%nat × Fulp b radix precision (FPred b radix precision u))%R →
MinOrMax radix b (a1 × x1 + y1) u.
intros H Hblop H1 H2 H3.
cut
(Rabs (a1 × x1 + y1 - u) <
/ 2%nat × Fulp b radix precision (FPred b radix precision u) +
Rabs (t + y - u))%R; [ intros H4 | idtac ].
case (Rle_or_lt (t + y) u); intros H5.
apply MinOrMax1 with precision; auto with zarith arith.
apply Rlt_le_trans with (1 := H4).
apply
Rplus_le_reg_l
with (- (/ 2%nat × Fulp b radix precision (FPred b radix precision u)))%R.
ring_simplify
(- (/ 2%nat × Fulp b radix precision (FPred b radix precision u)) +
(/ 2%nat × Fulp b radix precision (FPred b radix precision u) +
Rabs (FtoRradix t + FtoRradix y - FtoRradix u)))%R.
apply
Rle_trans
with (/ 2%nat × Fulp b radix precision (FPred b radix precision u))%R;
[ idtac | right ].
generalize uDef; unfold Closest in |- *; intros H6; elim H6; intros H7 H8;
clear H6 H7.
case
(Rle_or_lt (Rabs (FtoRradix t + FtoRradix y - FtoRradix u))
(/ 2%nat × Fulp b radix precision (FPred b radix precision u)));
auto; intros H6.
cut
(Rabs (FtoR radix u - (FtoRradix t + FtoRradix y)) ≤
Rabs (FtoR radix (FPred b radix precision u) - (FtoRradix t + FtoRradix y)))%R;
[ intros H7 | apply H8 ].
2: apply FBoundedPred; auto with zarith arith.
Contradict H7; apply Rlt_not_le.
apply Rlt_le_trans with (Rabs (FtoRradix t + FtoRradix y - FtoRradix u)).
2: right; rewrite <- Rabs_Ropp.
2: replace (- (FtoRradix t + FtoRradix y - FtoRradix u))%R with
(FtoR radix u - (FtoRradix t + FtoRradix y))%R;
auto with real; ring.
apply Rle_lt_trans with (2 := H6).
rewrite Faux.Rabsolu_left1.
apply
Rplus_le_reg_l
with
(u - (FtoRradix t + FtoRradix y) -
/ 2%nat × Fulp b radix precision (FPred b radix precision u))%R.
ring_simplify.
fold FtoRradix in |- *; pattern (FtoRradix u) at 1 in |- *;
replace (FtoRradix u) with
(FtoRradix (FPred b radix precision u) +
Fulp b radix precision (FPred b radix precision u))%R;
[ idtac
| unfold FtoRradix in |- *; apply FpredUlpPos; auto with zarith arith ].
apply
Rle_trans
with
(Fulp b radix precision (FPred b radix precision u) -
Fulp b radix precision (FPred b radix precision u) × / 2%nat)%R;
[ right; ring | idtac ].
apply
Rle_trans
with (/ 2%nat × Fulp b radix precision (FPred b radix precision u))%R.
right;
apply
trans_eq
with
((1 - / 2%nat) × Fulp b radix precision (FPred b radix precision u))%R;
[ ring | auto with real arith ].
replace (1 - / 2%nat)%R with (/ 2%nat)%R; auto with real arith.
rewrite <- (Rinv_r (INR 2)); auto with real arith.
simpl; field; auto with real.
apply Rlt_le; apply Rlt_le_trans with (1 := H6).
rewrite Faux.Rabsolu_left1; [ right; ring | idtac ].
apply Rplus_le_reg_l with (FtoRradix u).
ring_simplify; auto.
fold FtoRradix in |- *;
apply
Rplus_le_reg_l with (Fulp b radix precision (FPred b radix precision u)).
apply
Rle_trans
with
(FtoRradix (FPred b radix precision u) +
Fulp b radix precision (FPred b radix precision u) +
- (FtoRradix t + FtoRradix y))%R;
[ right; ring
| unfold FtoRradix in |- *; rewrite FpredUlpPos; auto with zarith arith ].
ring_simplify (Fulp b radix precision (FPred b radix precision u) + 0)%R.
apply Rle_trans with (Rabs (FtoRradix u + - (FtoRradix t + FtoRradix y))).
apply RRle_abs.
rewrite <- Rabs_Ropp.
replace (- (FtoRradix u + - (FtoRradix t + FtoRradix y)))%R with
(FtoRradix t + FtoRradix y - u)%R; [ idtac | ring ].
apply Rmult_le_reg_l with (INR 2); auto with arith real.
apply Rle_trans with (Fulp b radix precision u).
unfold FtoRradix in |- *; apply ClosestUlp; auto with arith zarith.
replace (INR 2) with (IZR radix); auto with zarith arith real.
apply FulpFPredLe; auto with zarith.
apply
trans_eq
with ((1 - / 2%nat) × Fulp b radix precision (FPred b radix precision u))%R;
[ idtac | ring ].
replace (1 - / 2%nat)%R with (/ 2%nat)%R; auto with real arith.
rewrite <- (Rinv_r (INR 2)); auto with real arith.
simpl; field; auto with real.
case
(Rle_or_lt (t + y)
(u + / 2%nat × Fulp b radix precision (FPred b radix precision u)));
intros H6.
apply MinOrMax1 with precision; auto with arith zarith real float.
apply Rlt_le_trans with (1 := H4); rewrite Rabs_right.
apply
Rplus_le_reg_l
with
(u + - (/ 2%nat × Fulp b radix precision (FPred b radix precision u)))%R.
ring_simplify.
apply Rle_trans with (1 := H6).
apply
Rle_trans
with
(FtoRradix u +
(Fulp b radix precision (FPred b radix precision u) +
- (Fulp b radix precision (FPred b radix precision u) × / 2%nat)))%R;
[ apply Rplus_le_compat_l; right | right; ring ].
apply
trans_eq
with ((1 - / 2%nat) × Fulp b radix precision (FPred b radix precision u))%R;
[ idtac | ring ].
replace (1 - / 2%nat)%R with (/ 2%nat)%R; auto with real arith.
rewrite <- (Rinv_r (INR 2)); auto with real arith.
simpl; field; auto with real.
apply Rle_ge; apply Rplus_le_reg_l with (FtoRradix u); auto with real.
apply MinOrMax2 with precision; auto with arith zarith float real.
apply Rlt_le_trans with (1 := H4).
apply
Rle_trans
with
(/ 2%nat × Fulp b radix precision u +
Rabs (FtoRradix t + FtoRradix y - FtoRradix u))%R;
[ apply Rplus_le_compat_r | idtac ].
apply Rmult_le_compat_l; auto with real arith.
apply FulpFPredGePos; auto with arith zarith.
apply
Rle_trans
with
(/ 2%nat × Fulp b radix precision u + / 2%nat × Fulp b radix precision u)%R;
[ apply Rplus_le_compat_l | right ].
apply Rmult_le_reg_l with (INR 2); auto with real arith.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real arith.
ring_simplify (1 × Fulp b radix precision u)%R.
unfold FtoRradix in |- *; apply ClosestUlp; auto with real arith zarith float.
apply trans_eq with ((/ 2%nat + / 2%nat) × Fulp b radix precision u)%R;
[ ring | auto with real arith ].
replace (/ 2%nat + / 2%nat)%R with 1%R; [ ring | auto with real arith ].
rewrite <- (Rinv_r (INR 2)); auto with real arith; simpl in |- *; ring.
cut (∀ z z' : R, (Rabs z ≤ z')%R → (- z' ≤ z)%R);
[ intros V | idtac ].
replace (a1 × x1 + y1)%R with
(a1 × x1 - a × x + (y1 - y) + (a × x - t + (t + y)))%R;
[ fold FtoRradix in |- × | ring ].
replace (FtoRradix u) with
(- (/ 4%nat × Fulp b radix precision (FPred b radix precision u)) +
(- (/ 4%nat × Fulp b radix precision (FPred b radix precision u)) +
(FtoRradix u +
/ 2%nat × Fulp b radix precision (FPred b radix precision u))))%R.
apply Rplus_le_compat.
apply V; auto with real.
apply
Rle_trans
with
(Rabs (y1 - FtoRradix y) + Rabs (a1 × x1 - FtoRradix a × FtoRradix x))%R;
auto with real.
rewrite Rplus_comm; apply Rabs_triang.
apply Rplus_le_compat.
apply V; auto with real.
auto with real.
apply
trans_eq
with
(FtoRradix u +
(- / 4%nat + (- / 4%nat + / 2%nat)) ×
Fulp b radix precision (FPred b radix precision u))%R;
[ ring | idtac ].
replace (- / 4%nat + (- / 4%nat + / 2%nat))%R with 0%R; [ ring | idtac ].
replace (/ 2%nat)%R with (2%nat × / 4%nat)%R; [ simpl in |- *; ring | idtac ].
replace (INR 4) with (2%nat × 2%nat)%R; auto with arith real zarith.
rewrite Rinv_mult_distr; auto with real arith.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real arith.
rewrite <- mult_INR; auto with real arith.
intros z z'; case (Rle_or_lt 0 z); intros P1 P2.
apply Rle_trans with (-0)%R; [ apply Ropp_le_contravar | simpl in |- × ];
auto with real.
apply Rle_trans with (Rabs z); auto with real.
replace (-0)%R with 0%R; auto with real.
apply Ropp_le_cancel; rewrite Ropp_involutive; rewrite <- Rabs_left;
auto with real.
replace (a1 × x1 + y1 - FtoRradix u)%R with
(y1 - y + (a1 × x1 - a × x) + (a × x - t + (t + y - u)))%R;
[ idtac | ring ].
apply
Rle_lt_trans
with
(Rabs (y1 - FtoRradix y + (a1 × x1 - FtoRradix a × FtoRradix x)) +
Rabs
(FtoRradix a × FtoRradix x - FtoRradix t +
(FtoRradix t + FtoRradix y - FtoRradix u)))%R;
[ apply Rabs_triang | idtac ].
apply
Rle_lt_trans
with
(Rabs (y1 - FtoRradix y + (a1 × x1 - FtoRradix a × FtoRradix x)) +
(Rabs (FtoRradix a × FtoRradix x - FtoRradix t) +
Rabs (FtoRradix t + FtoRradix y - FtoRradix u)))%R;
[ apply Rplus_le_compat_l; apply Rabs_triang | idtac ].
apply
Rle_lt_trans
with
(Rabs (y1 - FtoRradix y) + Rabs (a1 × x1 - FtoRradix a × FtoRradix x) +
(Rabs (FtoRradix a × FtoRradix x - FtoRradix t) +
Rabs (FtoRradix t + FtoRradix y - FtoRradix u)))%R.
apply Rplus_le_compat_r; apply Rabs_triang.
apply
Rlt_le_trans
with
(/ 4%nat × Fulp b radix precision (FPred b radix precision u) +
(Rabs (FtoRradix a × FtoRradix x - FtoRradix t) +
Rabs (FtoRradix t + FtoRradix y - FtoRradix u)))%R;
auto with real.
apply
Rle_trans
with
(/ 4%nat × Fulp b radix precision (FPred b radix precision u) +
(/ 4%nat × Fulp b radix precision (FPred b radix precision u) +
Rabs (FtoRradix t + FtoRradix y - FtoRradix u)))%R;
auto with real.
apply
Rle_trans
with
((/ 4%nat + / 4%nat) × Fulp b radix precision (FPred b radix precision u) +
Rabs (FtoRradix t + FtoRradix y - FtoRradix u))%R;
[ right; ring | apply Rplus_le_compat_r ].
replace (/ 4%nat + / 4%nat)%R with (/ 2%nat)%R; auto with real.
replace (/ 2%nat)%R with (2%nat × / 4%nat)%R; [ simpl in |- *; ring | idtac ].
replace (INR 4) with (2%nat × 2%nat)%R; auto with arith real zarith.
rewrite Rinv_mult_distr; auto with real arith.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real arith.
rewrite <- mult_INR; auto with real arith.
Qed.
Theorem Axpy_aux1_aux1 :
Fnormal radix b t →
Fcanonic radix b u →
(0 < u)%R →
(4%nat × Rabs t ≤ Rabs u)%R →
(Rabs (FtoRradix a × FtoRradix x - FtoRradix t) ≤
/ 4%nat × Fulp b radix precision (FPred b radix precision u))%R.
intros H1 H2 H3 H4.
cut (Fcanonic radix b t); [ intros H5 | left; auto ].
apply Rle_trans with (/ 2%nat × Fulp b radix precision t)%R.
apply Rmult_le_reg_l with (INR 2); auto with real arith.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real arith.
ring_simplify (1 × Fulp b radix precision t)%R.
unfold FtoRradix in |- *; apply ClosestUlp; auto with real arith zarith.
apply Rle_trans with (/ 2%nat × (/ 4%nat × Fulp b radix precision u))%R.
apply Rmult_le_compat_l; auto with real arith.
apply Rmult_le_reg_l with (INR 4); auto with real arith.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real arith.
ring_simplify (1 × Fulp b radix precision u)%R.
cut (Fbounded b (Float (Fnum t) (Zsucc (Zsucc (Fexp t)))));
[ intros H6 | idtac ].
2: elim H1; intros H7 H8; elim H7; intros H9 H10.
2: repeat (split; simpl in |- *; auto with zarith).
apply
Rle_trans
with (Fulp b radix precision (Float (Fnum t) (Zsucc (Zsucc (Fexp t))))).
right; unfold Fulp in |- ×.
replace (Fnormalize radix b precision t) with t;
[ idtac
| apply
FcanonicUnique with (radix := radix) (b := b) (precision := precision);
auto with real float zarith ].
2: apply sym_eq; apply FnormalizeCorrect; auto with real zarith.
replace
(Fnormalize radix b precision (Float (Fnum t) (Zsucc (Zsucc (Fexp t)))))
with (Float (Fnum t) (Zsucc (Zsucc (Fexp t)))).
replace (Fexp (Float (Fnum t) (Zsucc (Zsucc (Fexp t))))) with
(Zsucc (Zsucc (Fexp t))); [ idtac | simpl in |- *; auto ].
replace (Zsucc (Zsucc (Fexp t))) with (2 + Fexp t)%Z;
[ rewrite powerRZ_add | unfold Zsucc in |- × ]; auto with zarith real.
replace (powerRZ radix 2) with (INR 4); [ idtac | simpl in |- × ];
auto with real zarith; ring.
apply FcanonicUnique with (radix := radix) (b := b) (precision := precision);
auto with real float arith zarith.
2: apply sym_eq; apply FnormalizeCorrect; auto with real zarith.
elim H1; intros H7 H8.
left; repeat (split; simpl in |- *; auto with zarith).
rewrite
FulpFabs with b radix precision (Float (Fnum t) (Zsucc (Zsucc (Fexp t))));
auto with zarith.
rewrite FulpFabs with b radix precision u; auto with zarith.
apply LeFulpPos; auto with zarith real float.
unfold FtoRradix in |- *; rewrite Fabs_correct; auto with real zarith.
replace (FtoR radix (Fabs (Float (Fnum t) (Zsucc (Zsucc (Fexp t)))))) with
(Zabs (Fnum t) × powerRZ radix (Zsucc (Zsucc (Fexp t))))%R;
[ idtac | unfold FtoR in |- *; simpl in |- *; auto ].
replace (Zsucc (Zsucc (Fexp t))) with (2 + Fexp t)%Z;
[ rewrite powerRZ_add | unfold Zsucc in |- × ]; auto with zarith real.
replace (powerRZ radix 2) with (INR 4);
[ idtac | simpl in |- *; auto with real zarith; ring ].
rewrite Rmult_comm; rewrite Rmult_assoc.
replace (powerRZ radix (Fexp t) × Zabs (Fnum t))%R with (Rabs (FtoRradix t));
[ unfold FtoRradix in |- *; repeat rewrite Fabs_correct;
auto with real zarith
| idtac ].
unfold FtoRradix in |- *; rewrite <- Fabs_correct; auto with zarith;
unfold FtoR in |- *; simpl in |- *; ring.
apply Rmult_le_reg_l with (INR 2); auto with arith real.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with arith real.
ring_simplify (1 × (/ 4%nat × Fulp b radix precision u))%R.
apply
Rle_trans
with
(/ 4%nat × (2%nat × Fulp b radix precision (FPred b radix precision u)))%R;
[ apply Rmult_le_compat_l; auto with real arith | right; ring ].
replace (INR 2) with (IZR radix); auto with arith zarith real.
apply FulpFPredLe; auto with arith zarith real.
Qed.
Theorem Axpy_aux2 :
Fcanonic radix b u →
Fsubnormal radix b t →
(0 < u)%R →
FtoRradix u = (t + y)%R →
(Rabs (y1 - y) + Rabs (a1 × x1 - a × x) <
/ 4%nat × Fulp b radix precision (FPred b radix precision u))%R →
MinOrMax radix b (a1 × x1 + y1) u.
intros H1 H2 H3 H4 H5.
apply MinOrMax1 with precision; auto with zarith; fold FtoRradix in |- ×.
replace (a1 × x1 + y1 - FtoRradix u)%R with
(y1 - y + (a1 × x1 - a × x) + (a × x - t + (t + y - u)))%R;
[ idtac | ring ].
apply
Rle_lt_trans
with
(Rabs (y1 - FtoRradix y + (a1 × x1 - FtoRradix a × FtoRradix x)) +
Rabs
(FtoRradix a × FtoRradix x - FtoRradix t +
(FtoRradix t + FtoRradix y - FtoRradix u)))%R;
[ apply Rabs_triang | idtac ].
apply
Rle_lt_trans
with
(Rabs (y1 - FtoRradix y) + Rabs (a1 × x1 - FtoRradix a × FtoRradix x) +
Rabs
(FtoRradix a × FtoRradix x - FtoRradix t +
(FtoRradix t + FtoRradix y - FtoRradix u)))%R;
[ apply Rplus_le_compat_r; apply Rabs_triang | idtac ].
apply
Rle_lt_trans
with
(Rabs (y1 - FtoRradix y) + Rabs (a1 × x1 - FtoRradix a × FtoRradix x) +
(Rabs (FtoRradix a × FtoRradix x - FtoRradix t) +
Rabs (FtoRradix t + FtoRradix y - FtoRradix u)))%R;
[ apply Rplus_le_compat_l; apply Rabs_triang | idtac ].
apply
Rlt_le_trans
with
(/ 4%nat × Fulp b radix precision (FPred b radix precision u) +
(Rabs (FtoRradix a × FtoRradix x - FtoRradix t) +
Rabs (FtoRradix t + FtoRradix y - FtoRradix u)))%R;
[ apply Rplus_lt_compat_r; auto | idtac ].
apply
Rle_trans
with
(/ 4%nat × Fulp b radix precision (FPred b radix precision u) +
(/ 2%nat × powerRZ radix (- dExp b) + 0))%R.
apply Rplus_le_compat; auto with real.
apply Rplus_le_compat.
apply Rle_trans with (/ 2%nat × Fulp b radix precision t)%R;
[ apply Rmult_le_reg_l with (INR 2) | apply Rmult_le_compat_l ];
auto with real arith.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real arith.
ring_simplify (1 × Fulp b radix precision t)%R; unfold FtoRradix in |- *;
apply ClosestUlp; auto with zarith.
unfold Fulp in |- *; replace (Fnormalize radix b precision t) with t.
elim H2; intros H6 H7; elim H7; intros H8 H9; rewrite H8;
auto with zarith real.
apply FcanonicUnique with (radix := radix) (b := b) (precision := precision);
auto with real zarith float.
right; auto.
apply sym_eq; apply FnormalizeCorrect; auto with zarith.
rewrite H4; right.
ring_simplify (FtoRradix t + FtoRradix y - (FtoRradix t + FtoRradix y))%R;
apply Rabs_R0.
replace (/ 2%nat × powerRZ radix (- dExp b) + 0)%R with
(/ 2%nat × powerRZ radix (- dExp b))%R; [ idtac | ring ].
apply
Rle_trans
with
(/ 4%nat × Fulp b radix precision (FPred b radix precision u) +
/ 2%nat × Fulp b radix precision (FPred b radix precision u))%R;
[ apply Rplus_le_compat_l; apply Rmult_le_compat_l; auto with real arith
| idtac ].
unfold Fulp in |- *; apply Rle_powerRZ; auto with zarith real.
cut (Fbounded b (Fnormalize radix b precision (FPred b radix precision u)));
[ intros H6; elim H6; auto | idtac ].
apply FnormalizeBounded; auto with zarith arith.
apply FBoundedPred; auto with zarith arith.
apply
Rle_trans
with
((/ 4%nat + / 2%nat) × Fulp b radix precision (FPred b radix precision u))%R;
[ right; ring | idtac ].
apply
Rle_trans with (1 × Fulp b radix precision (FPred b radix precision u))%R;
[ apply Rmult_le_compat_r; auto with real zarith | right; ring ].
unfold Fulp in |- *; auto with real zarith.
apply Rmult_le_reg_l with (INR 4); auto with real arith.
apply Rle_trans with (4%nat × / 4%nat + 2%nat × (2%nat × / 2%nat))%R;
[ right; simpl; ring | idtac ].
replace (2%nat × 2%nat)%R with (INR (2 × 2)); auto with arith real.
repeat rewrite Rinv_r; auto with real arith.
simpl; ring_simplify; auto with real.
apply Rle_trans with (3+1)%R; auto with real.
right; ring.
Qed.
Theorem Axpy_aux1_aux2 :
Fsubnormal radix b t →
Fcanonic radix b u →
(0 < u)%R →
(Zsucc (- dExp b) ≤ Fexp (FPred b radix precision u))%Z →
(Rabs (FtoRradix a × FtoRradix x - FtoRradix t) ≤
/ 4%nat × Fulp b radix precision (FPred b radix precision u))%R.
intros H1 H2 H3 H4.
apply Rle_trans with (/ 2%nat × Fulp b radix precision t)%R;
[ apply Rmult_le_reg_l with (INR 2) | idtac ]; auto with real arith.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real arith.
ring_simplify (1 × Fulp b radix precision t)%R; unfold FtoRradix in |- *;
apply ClosestUlp; auto with zarith.
unfold Fulp in |- *; replace (Fnormalize radix b precision t) with t.
2: apply
FcanonicUnique with (radix := radix) (b := b) (precision := precision);
auto with real zarith float.
2: right; auto.
2: apply sym_eq; apply FnormalizeCorrect; auto with zarith.
replace (Fnormalize radix b precision (FPred b radix precision u)) with
(FPred b radix precision u).
2: apply
FcanonicUnique with (radix := radix) (b := b) (precision := precision);
auto with real zarith float.
2: apply sym_eq; apply FnormalizeCorrect; auto with zarith.
elim H1; intros H5 H6; elim H6; intros H7 H8; rewrite H7.
apply Rmult_le_reg_l with (INR 4); auto with real arith.
repeat rewrite <- Rmult_assoc.
pattern (INR 4) at 1 in |- *; replace (INR 4) with (2%nat × 2%nat)%R.
rewrite Rmult_comm; rewrite Rmult_assoc; repeat rewrite Rinv_r;
auto with real arith.
ring_simplify.
replace (INR 2) with (powerRZ radix 1);
[ rewrite <- powerRZ_add | simpl in |- × ]; auto with zarith real.
replace (- dExp b + 1)%Z with (Zsucc (- dExp b));
[ apply Rle_powerRZ | unfold Zsucc in |- × ]; auto with zarith real.
rewrite <- mult_INR; auto with arith.
Qed.
Theorem Axpy_aux1_aux3 :
Fsubnormal radix b t →
Fcanonic radix b u →
(0 < u)%R →
(Zsucc (- dExp b) ≤ Fexp (FPred b radix precision u))%Z →
(Rabs (FtoRradix a × FtoRradix x - FtoRradix t) ≤
/ 4%nat × Fulp b radix precision (FPred b radix precision u))%R.
intros H1 H2 H3 H4.
apply Rle_trans with (/ 2%nat × Fulp b radix precision t)%R.
apply Rmult_le_reg_l with (INR 2); auto with real arith;
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with arith real.
ring_simplify (1 × Fulp b radix precision t)%R; unfold FtoRradix in |- *;
apply ClosestUlp; auto with zarith.
apply Rmult_le_reg_l with (INR 4); auto with real arith;
repeat rewrite <- Rmult_assoc.
replace (INR 4) with (2%nat × 2%nat)%R;
[ rewrite (Rmult_assoc 2%nat 2%nat (/ 2%nat)); repeat rewrite Rinv_r;
auto with real arith
| rewrite <- mult_INR; auto with real arith ].
ring_simplify; unfold Fulp in |- ×.
replace (Fnormalize radix b precision t) with t;
[ elim H1; intros H9 H10; elim H10; intros H11 H12; rewrite H11 | idtac ].
2: apply
FcanonicUnique with (radix := radix) (b := b) (precision := precision);
auto with real zarith float.
2: right; auto.
2: apply sym_eq; apply FnormalizeCorrect; auto with zarith.
replace (Fnormalize radix b precision (FPred b radix precision u)) with
(FPred b radix precision u).
2: apply
FcanonicUnique with (radix := radix) (b := b) (precision := precision);
auto with real zarith float.
2: apply sym_eq; apply FnormalizeCorrect; auto with zarith.
replace (2%nat × powerRZ radix (- dExp b))%R with
(powerRZ radix (Zsucc (- dExp b)));
[ apply Rle_powerRZ | unfold Zsucc in |- × ]; auto with zarith real float.
rewrite powerRZ_add; auto with real zarith; simpl in |- *; ring.
Qed.
Theorem Axpy_aux3 :
Fcanonic radix b u →
Fsubnormal radix b t →
(0 < u)%R →
Fexp (FPred b radix precision u) = (- dExp b)%Z →
(Zsucc (- dExp b) ≤ Fexp u)%Z →
(Rabs (y1 - y) + Rabs (a1 × x1 - a × x) <
/ 4%nat × Fulp b radix precision (FPred b radix precision u))%R →
MinOrMax radix b (a1 × x1 + y1) u.
intros H1 H2 H3 H4 H5 H6.
case (Rle_or_lt u (a1 × x1 + y1)); intros H7.
apply MinOrMax2 with precision; auto with zarith float real;
fold FtoRradix in |- ×.
replace (a1 × x1 + y1 - FtoRradix u)%R with
(y1 - y + (a1 × x1 - a × x) + (a × x - t + (t + y - u)))%R;
[ idtac | ring ].
apply
Rlt_le_trans
with
(/ 4%nat × powerRZ radix (- dExp b) +
(/ 2%nat × powerRZ radix (- dExp b) + / 2%nat × Fulp b radix precision u))%R.
apply
Rle_lt_trans
with
(Rabs (y1 - FtoRradix y + (a1 × x1 - FtoRradix a × FtoRradix x)) +
Rabs
(FtoRradix a × FtoRradix x - FtoRradix t +
(FtoRradix t + FtoRradix y - FtoRradix u)))%R;
[ apply Rabs_triang | idtac ].
apply
Rlt_le_trans
with
(/ 4%nat × powerRZ radix (- dExp b) +
Rabs
(FtoRradix a × FtoRradix x - FtoRradix t +
(FtoRradix t + FtoRradix y - FtoRradix u)))%R;
[ apply Rplus_lt_compat_r | apply Rplus_le_compat_l ].
apply
Rle_lt_trans
with
(Rabs (y1 - FtoRradix y) + Rabs (a1 × x1 - FtoRradix a × FtoRradix x))%R;
[ apply Rabs_triang | idtac ].
apply Rlt_le_trans with (1 := H6); unfold Fulp in |- ×.
replace (Fexp (Fnormalize radix b precision (FPred b radix precision u)))
with (- dExp b)%Z; auto with real zarith.
replace (Fnormalize radix b precision (FPred b radix precision u)) with
(FPred b radix precision u); auto with zarith.
apply FcanonicUnique with (radix := radix) (b := b) (precision := precision);
auto with real zarith float.
apply sym_eq; apply FnormalizeCorrect; auto with zarith.
apply
Rle_trans
with
(Rabs (FtoRradix a × FtoRradix x - FtoRradix t) +
Rabs (FtoRradix t + FtoRradix y - FtoRradix u))%R;
[ apply Rabs_triang | apply Rplus_le_compat ].
apply Rmult_le_reg_l with (INR 2);
[ idtac | rewrite <- Rmult_assoc; rewrite Rinv_r ];
auto with arith real.
ring_simplify; apply Rle_trans with (Fulp b radix precision t);
[ unfold FtoRradix in |- *; apply ClosestUlp; auto with zarith
| unfold Fulp in |- × ].
replace (Fexp (Fnormalize radix b precision t)) with (- dExp b)%Z;
auto with real zarith.
replace (Fnormalize radix b precision t) with t; elim H2; auto with zarith.
intros;
apply FcanonicUnique with (radix := radix) (b := b) (precision := precision);
auto with real zarith float.
right; auto.
apply sym_eq; apply FnormalizeCorrect; auto with zarith.
apply Rmult_le_reg_l with (INR 2);
[ idtac | rewrite <- Rmult_assoc; rewrite Rinv_r ];
auto with arith real.
ring_simplify (1 × Fulp b radix precision u)%R; unfold FtoRradix in |- *;
apply ClosestUlp; auto with zarith.
apply Rmult_le_reg_l with (INR 4); auto with arith real.
apply Rle_trans with (powerRZ radix (-dExp b) × 3%nat
+ Fulp b radix precision u × 2%nat)%R.
right; simpl; field; auto with real.
apply
Rle_trans
with
(Fulp b radix precision u × 2%nat + Fulp b radix precision u × 2%nat)%R;
[ apply Rplus_le_compat_r | idtac ]; auto with real arith.
apply Rle_trans with (3%nat × powerRZ radix (- dExp b))%R; [ right;ring | idtac ].
apply Rle_trans with (4%nat × powerRZ radix (- dExp b))%R;
auto with arith zarith real.
unfold Fulp in |- *; replace (INR 2) with (powerRZ radix 1);
[ idtac | simpl in |- *; auto with zarith real ].
replace (INR 4) with (powerRZ radix 2);
[ idtac | simpl in |- *; auto with zarith real; ring ].
repeat rewrite <- powerRZ_add; auto with zarith real.
replace (2 + - dExp b)%Z with (Zsucc (- dExp b) + 1)%Z;
[ apply Rle_powerRZ | unfold Zsucc in |- × ]; auto with zarith real.
replace (Fnormalize radix b precision u) with u; auto with zarith arith.
apply FcanonicUnique with (radix := radix) (b := b) (precision := precision);
auto with real zarith float.
apply sym_eq; apply FnormalizeCorrect; auto with zarith.
right; replace (INR 4) with (2%nat + 2%nat)%R;
[ ring | rewrite <- plus_INR; auto with arith real ].
case (Rle_or_lt (t + y) u); intros H8.
apply Axpy_aux2; auto.
apply sym_eq; unfold FtoRradix in |- *; apply ExactSum_Near with b precision;
auto with zarith.
elim H2; intuition.
cut (∀ v w : R, (v ≤ w)%R → v ≠ w → (v < w)%R); [ intros V | idtac ].
2: intros V W H' H'1; case H'; auto with real.
2: intros H'2; absurd (V = W); auto with real.
apply V.
apply Rmult_le_reg_l with (INR 2); auto with real arith.
apply Rle_trans with (Fulp b radix precision u);
[ unfold FtoRradix in |- *; apply ClosestUlp; auto with zarith | idtac ].
replace (powerRZ 2%Z (- dExp b)) with
(Fulp b radix precision (FPred b radix precision u)).
replace (INR 2) with (IZR radix); auto with zarith real.
apply FulpFPredLe; auto with zarith.
unfold Fulp in |- *;
replace (Fnormalize radix b precision (FPred b radix precision u)) with
(FPred b radix precision u).
rewrite H4; auto with zarith real.
apply FcanonicUnique with (radix := radix) (b := b) (precision := precision);
auto with real zarith float.
apply sym_eq; apply FnormalizeCorrect; auto with zarith.
replace (powerRZ 2%Z (- dExp b)) with
(Fulp b radix precision (FPred b radix precision u)).
2: unfold Fulp in |- *;
replace (Fnormalize radix b precision (FPred b radix precision u)) with
(FPred b radix precision u).
2: rewrite H4; auto with zarith real.
rewrite Faux.Rabsolu_left1; auto with real.
2: apply Rplus_le_reg_l with (FtoRradix u); unfold FtoRradix, radix in |- ×.
2: ring_simplify; auto with real zarith.
case
(Req_dec (- (FtoR 2 t + FtoR 2 y - FtoR 2 u))
(Fulp b radix precision (FPred b radix precision u)));
auto.
intros W; Contradict uDef.
replace (FtoRradix t + FtoRradix y)%R with
(FtoRradix (FPred b radix precision u)); auto with float real.
cut (FtoRradix u ≠ FPred b radix precision u); [ intros V' | idtac ].
Contradict V'.
cut (ProjectorP b radix (Closest b radix));
[ unfold ProjectorP in |- *; intros W' | idtac ].
apply sym_eq; apply W'; auto with zarith float real.
apply RoundedProjector; apply ClosestRoundedModeP with precision;
auto with zarith.
apply not_eq_sym; apply Rlt_dichotomy_converse; left.
unfold FtoRradix in |- *; apply FPredLt; auto with zarith.
apply
Rplus_eq_reg_l with (Fulp b radix precision (FPred b radix precision u)).
rewrite Rplus_comm; unfold FtoRradix in |- *; rewrite FpredUlpPos;
auto with zarith.
rewrite <- W; unfold FtoRradix, radix in |- *; ring.
apply FcanonicUnique with (radix := radix) (b := b) (precision := precision);
auto with real zarith float.
apply sym_eq; apply FnormalizeCorrect; auto with zarith.
apply MinOrMax1 with precision; auto with zarith; fold FtoRradix in |- ×.
rewrite Faux.Rabsolu_left1.
2: apply Rplus_le_reg_l with (FtoRradix u).
2: ring_simplify; auto with real.
apply Rle_lt_trans with (y - y1 + (t - a1 × x1) - (t + y - u))%R;
[ right; ring | idtac ].
apply Rlt_le_trans with (Rabs (y - y1 + (t - a1 × x1))).
cut (∀ u v : R, (0 < v)%R → (u - v < Rabs u)%R);
[ intros V; apply V | idtac ].
apply Rplus_lt_reg_r with (FtoRradix u).
ring_simplify; auto with real.
intros v w H'1.
apply Rlt_le_trans with v; auto with real.
unfold Rminus in |- *; apply Rlt_le_trans with (v + -0)%R;
[ auto with real | right; ring ].
apply RRle_abs.
replace (FtoRradix y - y1 + (FtoRradix t - a1 × x1))%R with
(- (y1 - FtoRradix y + (a1 × x1 - a × x) + (a × x - t)))%R;
[ idtac | ring ].
rewrite Rabs_Ropp.
apply
Rle_trans
with
(Rabs (y1 - FtoRradix y + (a1 × x1 - FtoRradix a × FtoRradix x)) +
Rabs (FtoRradix a × FtoRradix x - FtoRradix t))%R;
[ apply Rabs_triang | idtac ].
apply Rmult_le_reg_l with (INR 2); auto with real arith.
apply
Rle_trans
with
(2%nat × Rabs (y1 - FtoRradix y + (a1 × x1 - FtoRradix a × FtoRradix x)) +
2%nat × Rabs (FtoRradix a × FtoRradix x - FtoRradix t))%R;
[ right; ring | idtac ].
apply
Rle_trans
with
(Fulp b radix precision (FPred b radix precision u) +
Fulp b radix precision (FPred b radix precision u))%R;
[ idtac | right; simpl; ring ].
apply Rplus_le_compat.
apply
Rle_trans
with
(2%nat × (/ 4%nat × Fulp b radix precision (FPred b radix precision u)))%R;
auto with real.
apply
Rle_trans
with
(2%nat ×
(Rabs (y1 - FtoRradix y) + Rabs (a1 × x1 - FtoRradix a × FtoRradix x)))%R;
[ apply Rmult_le_compat_l; auto with real arith; apply Rabs_triang | idtac ].
apply Rmult_le_compat_l; auto with real arith.
apply
Rle_trans with (1 × Fulp b radix precision (FPred b radix precision u))%R;
[ rewrite <- Rmult_assoc; apply Rmult_le_compat_r | right; ring ];
auto with real arith zarith.
unfold Fulp in |- *; auto with real float zarith.
rewrite Rmult_comm; apply Rmult_le_reg_l with (INR 4); auto with arith real;
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real arith.
simpl; ring_simplify; auto with real arith.
apply Rle_trans with (Fulp b radix precision t);
[ unfold FtoRradix in |- *; apply ClosestUlp; auto with zarith
| unfold Fulp in |- × ].
apply Rle_powerRZ; auto with zarith real.
replace (Fnormalize radix b precision t) with t;
[ elim H2; intros H9 H10; elim H10; intros H11 H12; rewrite H11 | idtac ].
2: apply sym_eq; apply FcanonicFnormalizeEq; auto with zarith; right; auto.
replace (Fnormalize radix b precision (FPred b radix precision u)) with
(FPred b radix precision u); [ rewrite H4; auto with zarith | idtac ].
apply sym_eq; apply FcanonicFnormalizeEq; auto with zarith float.
Qed.
Theorem AxpyPos :
Fcanonic radix b u →
Fcanonic radix b t →
(0 < u)%R →
(4%nat × Rabs t ≤ Rabs u)%R →
(Rabs (y1 - y) + Rabs (a1 × x1 - a × x) <
/ 4%nat × Fulp b radix precision (FPred b radix precision u))%R →
MinOrMax radix b (a1 × x1 + y1) u.
intros H1 H2 H3 H4 H5.
case H2; intros H6.
apply Axpy_aux1; auto; apply Axpy_aux1_aux1; auto.
cut (∀ z1 z2 : Z, (z1 ≤ z2)%Z → z1 = z2 ∨ (Zsucc z1 ≤ z2)%Z);
[ intros V | idtac ].
2: intros z1 z2 V; omega.
case (V (- dExp b)%Z (Fexp u));
[ elim Fu; intuition | intros H7 | intros H7 ].
apply Axpy_aux2; auto.
unfold FtoRradix in |- *; apply plusExact1 with b precision;
auto with zarith float.
rewrite <- H7; elim H6; intros H8 H9; elim H9; intros H10 H11; rewrite H10.
rewrite Zmin_le1; auto with zarith; elim Fy; auto.
case (V (- dExp b)%Z (Fexp (FPred b radix precision u))).
cut (Fbounded b (FPred b radix precision u)); [ intros H8 | idtac ];
auto with float.
apply FBoundedPred; auto with zarith.
intros H8; apply Axpy_aux3; auto.
intros H8; apply Axpy_aux1; auto.
apply Axpy_aux1_aux3; auto.
Qed.
Definition FLess (p : float) :=
match Rcase_abs p with
| left _ ⇒ FSucc b radix precision p
| right _ ⇒ FPred b radix precision p
end.
Theorem UlpFlessuGe_aux :
∀ p : float,
Fbounded b p →
Fcanonic radix b p →
(Rabs p - Fulp b radix precision p ≤ Rabs (FLess p))%R.
intros p H H1.
cut
(∀ p : float,
Fbounded b p →
Fcanonic radix b p →
(0 < p)%R →
(Rabs p - Fulp b radix precision p ≤ Rabs (FPred b radix precision p))%R);
[ intros V | idtac ].
unfold FLess in |- *; case (Rcase_abs p); intros H2.
rewrite <- Rabs_Ropp; unfold FtoRradix in |- *; rewrite <- Fopp_correct.
pattern p at 3 in |- *; rewrite <- Fopp_Fopp.
rewrite <- (Rabs_Ropp (FtoR radix (FSucc b radix precision (Fopp (Fopp p)))));
rewrite <- (Fopp_correct radix (FSucc b radix precision (Fopp (Fopp p)))).
rewrite <- FPredFopFSucc with b radix precision (Fopp p); auto with zarith.
replace (Fulp b radix precision p) with (Fulp b radix precision (Fopp p));
auto with float.
fold FtoRradix in |- *; apply V; auto with float.
unfold FtoRradix in |- *; rewrite Fopp_correct; auto with real.
unfold Fulp in |- *; rewrite Fnormalize_Fopp; auto with arith.
unfold Fopp in |- *; simpl in |- *; auto with real.
cut (0 ≤ p)%R; [ intros H3; case H3; intros H4 | auto with real ].
apply V; auto.
rewrite <- H4; rewrite Rabs_R0.
ring_simplify (0 - Fulp b radix precision p)%R; apply Rle_trans with 0%R;
auto with real zarith.
unfold Fulp in |- *; auto with real zarith.
intros q H2 H3 H4.
unfold FtoRradix in |- *; rewrite <- FpredUlpPos with b radix precision q;
auto with zarith; fold FtoRradix in |- ×.
apply
Rle_trans
with
(Rabs (FtoRradix (FPred b radix precision q)) +
Rabs (Fulp b radix precision (FPred b radix precision q)) -
Fulp b radix precision q)%R;
[ unfold Rminus in |- *; apply Rplus_le_compat_r; apply Rabs_triang
| idtac ].
rewrite (Rabs_right (Fulp b radix precision (FPred b radix precision q)));
[ idtac | apply Rle_ge; unfold Fulp in |- *; auto with real zarith ].
apply Rplus_le_reg_l with (Fulp b radix precision q).
ring_simplify
(Fulp b radix precision q +
(Rabs (FtoRradix (FPred b radix precision q)) +
Fulp b radix precision (FPred b radix precision q) -
Fulp b radix precision q))%R.
rewrite Rplus_comm; apply Rplus_le_compat_r.
apply FulpFPredGePos; auto with zarith real float.
Qed.
Theorem UlpFlessuGe :
Fcanonic radix b u →
(/
(4%nat × (powerRZ radix precision - 1) × (1 + powerRZ radix (- precision))) ×
((1 - powerRZ radix (Zsucc (- precision))) × Rabs y) +
-
(/
(4%nat × (powerRZ radix precision - 1) × (1 + powerRZ radix (- precision)) ×
(1 - powerRZ radix (- precision))) ×
((1 - powerRZ radix (Zsucc (- precision))) × Rabs (a × x))) +
-
(powerRZ radix (Zpred (- dExp b)) ×
(/ (2%nat × (powerRZ radix precision - 1)) +
/
(4%nat × (powerRZ radix precision - 1) ×
(1 + powerRZ radix (- precision)) × (1 - powerRZ radix (- precision))) ×
(1 - powerRZ radix (Zsucc (- precision))))) ≤
/ 4%nat × Fulp b radix precision (FLess u))%R.
intros Cu.
cut (0 < 1 - powerRZ radix (- precision))%R; [ intros H'1 | idtac ].
2: apply Rplus_lt_reg_r with (powerRZ radix (- precision)).
2: ring_simplify.
2: replace 1%R with (powerRZ radix 0);
[ auto with real zarith | simpl in |- *; auto ].
cut (0 < 1 + powerRZ radix (- precision))%R; [ intros H'2 | idtac ].
2: apply Rlt_le_trans with 1%R; auto with real.
2: apply Rle_trans with (1 + 0)%R; auto with real zarith.
cut (0 < powerRZ radix precision - 1)%R; [ intros H'3 | idtac ].
2: apply Rplus_lt_reg_r with 1%R.
2: ring_simplify.
2: replace 1%R with (powerRZ radix 0);
[ auto with real zarith | simpl in |- *; auto ].
cut (0 < 1 - powerRZ radix (Zsucc (- precision)))%R; [ intros H'4 | idtac ].
2: apply Rplus_lt_reg_r with (powerRZ radix (Zsucc (- precision))).
2: ring_simplify.
2: replace 1%R with (powerRZ radix 0);
[ auto with real zarith | simpl in |- *; auto ].
cut (0 < INR 4)%R; [ intros H'5 | auto with arith real ].
rewrite
(Rinv_mult_distr
(4%nat × (powerRZ radix precision - 1) ×
(1 + powerRZ radix (- precision))) (1 - powerRZ radix (- precision)))
; auto with real.
2: apply Rmult_integral_contrapositive; split; auto with real.
apply
Rle_trans
with
(/
(4%nat × (powerRZ radix precision - 1) ×
(1 + powerRZ radix (- precision))) ×
(1 - powerRZ radix (Zsucc (- precision))) ×
(Rabs y -
(Rabs (FtoRradix a × FtoRradix x) × / (1 - powerRZ radix (- precision)) +
powerRZ radix (Zpred (- dExp b)) × / (1 - powerRZ radix (- precision)))) +
-
(powerRZ radix (Zpred (- dExp b)) ×
/ (2%nat × (powerRZ radix precision - 1))))%R;
[ right; ring; ring | idtac ].
apply
Rle_trans
with
(/
(4%nat × (powerRZ radix precision - 1) ×
(1 + powerRZ radix (- precision))) ×
(1 - powerRZ radix (Zsucc (- precision))) ×
(Rabs (FtoRradix y) - Rabs t) +
-
(powerRZ radix (Zpred (- dExp b)) ×
/ (2%nat × (powerRZ radix precision - 1))))%R;
[ apply Rplus_le_compat_r | idtac ].
apply Rmult_le_compat_l.
apply Rlt_le; apply Rmult_lt_0_compat; auto with real.
apply Rinv_0_lt_compat; apply Rmult_lt_0_compat; auto with real;
apply Rmult_lt_0_compat; auto with real.
unfold Rminus in |- *; apply Rplus_le_compat_l; apply Ropp_le_contravar.
apply
Rle_trans
with
(Rabs (FtoRradix a × FtoRradix x) × / (1 - powerRZ 2%Z (- precision)) +
powerRZ 2%Z (Zpred (- dExp b)) × / (1 - powerRZ 2%Z (- precision)))%R;
[ idtac | right; unfold FtoRradix, radix, Rminus in |- *; auto with real ].
apply RoundLeGeneral; auto with zarith.
rewrite
(Rinv_mult_distr (4%nat × (powerRZ radix precision - 1))
(1 + powerRZ radix (- precision))); auto with real.
apply
Rle_trans
with
(/ (4%nat × (powerRZ radix precision - 1)) ×
(1 - powerRZ radix (Zsucc (- precision))) ×
(/ (1 + powerRZ radix (- precision)) ×
(Rabs (FtoRradix y) - Rabs (FtoRradix t))) +
-
(powerRZ radix (Zpred (- dExp b)) ×
/ (2%nat × (powerRZ radix precision - 1))))%R;
[ right; ring | idtac ].
apply
Rle_trans
with
(/ (4%nat × (powerRZ radix precision - 1)) ×
(1 - powerRZ radix (Zsucc (- precision))) × Rabs u +
-
(powerRZ radix (Zpred (- dExp b)) ×
/ (2%nat × (powerRZ radix precision - 1))))%R.
apply Rplus_le_compat_r; apply Rmult_le_compat_l.
apply Rlt_le; apply Rmult_lt_0_compat; auto with real.
apply Rinv_0_lt_compat; apply Rmult_lt_0_compat; auto with real.
apply
Rle_trans
with
(/ (1 + powerRZ radix (- precision)) × Rabs (FtoRradix y + FtoRradix t))%R.
apply Rmult_le_compat_l; auto with real.
rewrite <- (Rabs_Ropp (FtoRradix t)).
replace (FtoRradix y + FtoRradix t)%R with (FtoRradix y - - FtoRradix t)%R;
[ apply Rabs_triang_inv | ring ].
case Cu; intros H1.
apply Rmult_le_reg_l with (1 + powerRZ radix (- precision))%R; auto.
apply
Rle_trans
with
(Rabs (FtoRradix y + FtoRradix t) ×
((1 + powerRZ radix (- precision)) × / (1 + powerRZ radix (- precision))))%R;
[ right; ring; ring | idtac ].
rewrite Rinv_r; auto with real.
ring_simplify.
apply Rle_trans with (Rabs u + / radix × Fulp b radix precision u)%R.
apply Rplus_le_reg_l with (- Rabs u)%R.
ring_simplify.
apply Rle_trans with (Rabs (FtoRradix y + FtoRradix t - FtoRradix u)).
apply
Rle_trans with (Rabs (FtoRradix y + FtoRradix t) - Rabs (FtoRradix u))%R;
[ right; ring | apply Rabs_triang_inv ].
apply Rmult_le_reg_l with (IZR radix); auto with real zarith.
apply Rle_trans with (Fulp b radix precision u × (radix × / radix))%R;
[ rewrite Rinv_r; auto with real zarith | right; ring ].
ring_simplify (Fulp b radix precision u × 1)%R.
unfold FtoRradix in |- *; apply ClosestUlp; auto with real zarith.
rewrite Rplus_comm; auto.
rewrite Rplus_comm; apply Rplus_le_compat_r.
apply Rmult_le_reg_l with (IZR radix); auto with real zarith.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real zarith.
ring_simplify (1 × Fulp b radix precision u)%R.
apply
Rle_trans with (Rabs (FtoRradix u) × powerRZ radix (Zsucc (- precision)))%R.
unfold FtoRradix in |- *; apply FulpLe2; auto with zarith.
replace (Fnormalize radix b precision u) with u;
[ idtac
| apply
FcanonicUnique with (radix := radix) (precision := precision) (b := b) ];
auto with zarith real.
apply FnormalizeCanonic; auto with zarith.
apply sym_eq; apply FnormalizeCorrect; auto with zarith.
unfold Zsucc in |- *; rewrite powerRZ_add; auto with zarith real;
simpl in |- *; auto with zarith real; right; ring.
replace (FtoRradix y + FtoRradix t)%R with (FtoRradix u); auto with real.
apply Rle_trans with (1 × Rabs (FtoRradix u))%R;
[ apply Rmult_le_compat_r | right; ring ]; auto with real.
pattern 1%R at 2 in |- *; replace 1%R with (/ 1)%R; auto with real zarith.
apply Rle_Rinv; auto with real zarith.
pattern 1%R at 1 in |- *; replace 1%R with (1 + 0)%R; auto with real zarith.
unfold FtoRradix in |- *; apply plusExact1 with b precision;
auto with real zarith.
rewrite Rplus_comm; auto.
elim H1; intros H5 H6; elim H6; intros H7 H8; rewrite H7; apply Zmin_Zle;
elim Fy; elim Ft; auto with zarith.
apply
Rle_trans
with
(/ (4%nat × (powerRZ radix precision - 1)) ×
(Rabs (FtoRradix u) -
(Rabs u × powerRZ radix (Zsucc (- precision)) +
powerRZ radix (- dExp b))))%R.
right; simpl; ring_simplify.
unfold Zpred in |- *; rewrite powerRZ_add; auto with zarith real.
simpl; field; auto with real.
apply
Rle_trans
with
(/ (4%nat × (powerRZ radix precision - 1)) ×
(Rabs (FtoRradix u) - Fulp b radix precision u))%R.
apply Rmult_le_compat_l; auto with real.
apply Rlt_le; apply Rinv_0_lt_compat; apply Rmult_lt_0_compat; auto with real.
unfold Rminus in |- *; apply Rplus_le_compat_l; apply Ropp_le_contravar.
apply FulpLeGeneral; auto.
apply
Rle_trans
with (/ (4%nat × (powerRZ radix precision - 1)) × Rabs (FLess u))%R.
apply Rmult_le_compat_l; auto with real.
apply Rlt_le; apply Rinv_0_lt_compat; apply Rmult_lt_0_compat; auto with real.
apply UlpFlessuGe_aux; auto.
rewrite Rinv_mult_distr; auto with real.
rewrite Rmult_assoc; apply Rmult_le_compat_l; auto with real.
apply Rmult_le_reg_l with (powerRZ radix precision - 1)%R; auto.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real.
ring_simplify (1 × Rabs (FtoRradix (FLess u)))%R; unfold FtoRradix in |- *;
apply FulpGe; auto with zarith.
unfold FLess in |- *; case (Rcase_abs u); intros H; auto with float.
apply FBoundedSuc; auto with zarith.
apply FBoundedPred; auto with zarith.
Qed.
Theorem UlpFlessuGe2 :
Fcanonic radix b u →
(powerRZ radix (Zpred (Zpred (- precision))) ×
(1 - powerRZ radix (Zsucc (- precision))) × Rabs y +
- (powerRZ radix (Zpred (Zpred (- precision))) × Rabs (a × x)) +
- powerRZ radix (Zpred (Zpred (- dExp b))) <
/ 4%nat × Fulp b radix precision (FLess u))%R.
intros H.
apply
Rlt_le_trans
with
(/
(4%nat × (powerRZ radix precision - 1) ×
(1 + powerRZ radix (- precision))) ×
((1 - powerRZ radix (Zsucc (- precision))) × Rabs (FtoRradix y)) +
-
(/
(4%nat × (powerRZ radix precision - 1) ×
(1 + powerRZ radix (- precision)) × (1 - powerRZ radix (- precision))) ×
((1 - powerRZ radix (Zsucc (- precision))) ×
Rabs (FtoRradix a × FtoRradix x))) +
-
(powerRZ radix (Zpred (- dExp b)) ×
(/ (2%nat × (powerRZ radix precision - 1)) +
/
(4%nat × (powerRZ radix precision - 1) ×
(1 + powerRZ radix (- precision)) × (1 - powerRZ radix (- precision))) ×
(1 - powerRZ radix (Zsucc (- precision))))))%R;
[ idtac | apply UlpFlessuGe; auto ].
cut (0 < 1 - powerRZ radix (- precision))%R; [ intros H'1 | idtac ].
2: apply Rplus_lt_reg_r with (powerRZ radix (- precision)).
2: ring_simplify.
2: replace 1%R with (powerRZ radix 0);
[ auto with real zarith | simpl in |- *; auto ].
cut (0 < 1 + powerRZ radix (- precision))%R; [ intros H'2 | idtac ].
2: apply Rlt_le_trans with 1%R; auto with real.
2: apply Rle_trans with (1 + 0)%R; auto with real zarith.
cut (0 < powerRZ radix precision - 1)%R; [ intros H'3 | idtac ].
2: apply Rplus_lt_reg_r with 1%R.
2: ring_simplify.
2: replace 1%R with (powerRZ radix 0);
[ auto with real zarith | simpl in |- *; auto ].
cut (0 < 1 - powerRZ radix (Zsucc (- precision)))%R; [ intros H'4 | idtac ].
2: apply Rplus_lt_reg_r with (powerRZ radix (Zsucc (- precision))).
2: ring_simplify.
2: replace 1%R with (powerRZ radix 0);
[ auto with real zarith | simpl in |- *; auto ].
cut (0 < INR 4)%R; [ intros H'5 | auto with arith real ].
apply
Rle_lt_trans
with
(/
(4%nat × (powerRZ radix precision - 1) ×
(1 + powerRZ radix (- precision))) ×
((1 - powerRZ radix (Zsucc (- precision))) × Rabs (FtoRradix y)) +
-
(/
(4%nat × (powerRZ radix precision - 1) ×
(1 + powerRZ radix (- precision)) × (1 - powerRZ radix (- precision))) ×
((1 - powerRZ radix (Zsucc (- precision))) ×
Rabs (FtoRradix a × FtoRradix x))) +
- powerRZ radix (Zpred (Zpred (- dExp b))))%R;
[ apply Rplus_le_compat_r | apply Rplus_lt_compat_l ].
apply Rplus_le_compat.
repeat rewrite <- Rmult_assoc.
apply Rmult_le_compat_r; auto with real.
rewrite Rmult_assoc.
ring_simplify ((powerRZ radix precision - 1) × (1 + powerRZ radix (- precision)))%R.
rewrite <- powerRZ_add; auto with real zarith.
replace (powerRZ radix (- precision + precision)) with 1%R;
[ idtac
| ring_simplify (- precision + precision)%Z; simpl in |- *; auto with real zarith ].
ring_simplify (-1 + (- powerRZ radix (- precision) + (powerRZ radix precision + 1)))%R.
apply Rmult_le_compat_r; auto with real.
apply Rle_trans with (/ (4%nat × powerRZ radix precision))%R;
[ right | idtac ].
unfold Zpred in |- ×.
replace (- precision + -1 + -1)%Z with (- 2%nat + - precision)%Z;
[ rewrite powerRZ_add | ring ]; auto with zarith real.
rewrite Rinv_mult_distr; auto with real zarith.
replace (powerRZ radix (- 2%nat)) with (/ 4%nat)%R.
replace (/ powerRZ radix precision)%R with (powerRZ radix (- precision));
auto with real.
apply powerRZ_Zopp; auto with zarith real.
rewrite powerRZ_Zopp; auto with zarith real.
replace (powerRZ radix 2%nat) with (INR 4);
[ auto with zarith real | simpl in |- *; ring ].
apply Rle_Rinv; auto with real.
apply
Rlt_le_trans
with
(4%nat ×
(powerRZ radix precision - 1 + (1 - powerRZ radix (- precision))))%R.
apply Rle_lt_trans with (4%nat × (0 + 0))%R; [ right; ring | auto with real ].
ring_simplify (precision+-precision)%Z; simpl; right;ring.
apply Rmult_le_compat_l; auto with real zarith.
ring_simplify (precision+-precision)%Z.
apply
Rle_trans with (powerRZ radix precision + - powerRZ radix (- precision))%R;
[ right; simpl; ring | idtac ]; auto with real zarith.
apply Rle_trans with (powerRZ radix precision + -0)%R;
[ idtac | right; ring ]; auto with real zarith.
apply Ropp_le_contravar.
repeat rewrite <- Rmult_assoc.
apply Rmult_le_compat_r; auto with real.
unfold Zpred in |- ×.
replace (- precision + -1 + -1)%Z with (- 2%nat + - precision)%Z;
[ rewrite powerRZ_add | ring ]; auto with zarith real.
repeat rewrite Rmult_assoc; rewrite Rinv_mult_distr; auto with real zarith.
2: apply Rmult_integral_contrapositive; split; auto with real.
replace (powerRZ radix (- 2%nat)) with (/ 4%nat)%R.
2: rewrite powerRZ_Zopp; auto with zarith real.
2: replace (powerRZ radix 2%nat) with (INR 4);
[ auto with zarith real | simpl in |- *; ring ].
repeat rewrite Rmult_assoc; apply Rmult_le_compat_l; auto with real.
apply
Rmult_le_reg_l
with
((powerRZ radix precision - 1) ×
((1 + powerRZ radix (- precision)) × (1 - powerRZ radix (- precision))))%R;
auto with real.
apply Rmult_lt_0_compat; auto with real.
apply Rmult_lt_0_compat; auto with real.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real.
2: apply Rmult_integral_contrapositive; split; auto with real.
ring_simplify.
repeat rewrite Ropp_mult_distr_l_reverse.
repeat rewrite <- powerRZ_add; auto with real zarith.
replace (powerRZ radix (precision + - precision)) with 1%R;
[ idtac | ring_simplify ( precision +- precision)%Z; simpl in |- *; auto ].
ring_simplify ( precision + - precision + - precision + -precision)%Z.
apply
Rle_trans
with
(1 +
(- powerRZ radix (- precision) + (- powerRZ radix (- precision) + 0)))%R;
[ right | idtac ].
rewrite powerRZ_add; auto with real zarith; simpl; ring.
ring_simplify (- precision + - precision + - precision)%Z.
apply
Rle_trans
with
(1 +
(- powerRZ radix (- precision) +
(- powerRZ radix (-2 × precision) + powerRZ radix (-3 × precision))))%R;
[ idtac | right; ring ].
apply Rplus_le_compat_l; apply Rplus_le_compat_l.
apply Rplus_le_compat; auto with real zarith.
apply Ropp_lt_contravar.
apply Rmult_lt_reg_l with (powerRZ radix (Z_of_nat 2 + dExp b));
auto with real zarith.
rewrite <- Rmult_assoc; unfold Zpred in |- ×.
repeat rewrite <- powerRZ_add; auto with real zarith.
replace (2%nat + dExp b + (- dExp b + -1 + -1))%Z with 0%Z by ring.
replace (2%nat + dExp b + (- dExp b + -1))%Z with 1%Z by ring.
replace (powerRZ radix 0) with 1%R;
[ idtac | simpl in |- *; auto with real zarith ].
replace (powerRZ radix 1) with 2%R;
[ idtac | simpl in |- *; ring; auto with real zarith ].
rewrite Rmult_plus_distr_l.
apply Rlt_le_trans with (/ 3%nat + 2%nat × / 3%nat)%R.
2: simpl; right;field; auto with real.
apply
Rle_lt_trans
with
(/ 3%nat +
2 ×
(/
(4%nat × (powerRZ radix precision - 1) × (1 + powerRZ radix (- precision)) ×
(1 - powerRZ radix (- precision))) × (1 - powerRZ radix (- precision + 1))))%R;
[ apply Rplus_le_compat_r | apply Rplus_lt_compat_l ].
rewrite Rinv_mult_distr; auto with arith real.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with arith real.
rewrite Rmult_1_l; replace (INR 3) with (powerRZ radix 2%nat - 1)%R;
auto with real zarith arith.
apply Rle_Rinv; auto with real arith.
replace 0%R with (powerRZ radix 0 - 1)%R;
[ unfold Rminus in |- *; auto with real arith zarith | simpl in |- *; ring ].
unfold Rminus in |- *; auto with real arith zarith.
apply Rplus_le_compat_r; apply Rle_powerRZ; auto with real arith zarith.
simpl in |- *; ring; auto with arith zarith real.
apply Rmult_lt_compat_l; auto with real arith.
repeat rewrite Rinv_mult_distr; auto with real.
2: apply Rmult_integral_contrapositive; split; auto with real.
apply
Rle_lt_trans
with
(/ 4%nat × / (powerRZ radix 2%nat - 1) × / (1 + 0) ×
/ (1 - powerRZ radix (- 2%nat)) × (1 - 0))%R.
repeat rewrite Rmult_assoc.
apply Rmult_le_compat_l; auto with real arith.
apply Rmult_le_compat; auto with real arith zarith.
apply Rmult_le_pos; auto with real.
apply Rmult_le_pos; auto with real.
apply Rle_Rinv; auto with real arith zarith.
apply Rle_lt_trans with (powerRZ radix 0 - 1)%R;
[ right; simpl in |- *; ring | auto with real arith zarith ].
unfold Rminus in |- *; apply Rplus_lt_compat_r; auto with real arith zarith.
unfold Rminus in |- *; apply Rplus_le_compat_r; apply Rle_powerRZ;
auto with real arith zarith.
apply Rmult_le_compat; auto with real arith zarith.
apply Rmult_le_pos; auto with real.
apply Rmult_le_compat; auto with real arith zarith.
apply Rle_Rinv; auto with real arith zarith.
apply Rle_lt_trans with (1 - powerRZ radix 0)%R;
[ simpl in |- *; right; ring | auto with real arith zarith ].
unfold Rminus in |- *; apply Rplus_lt_compat_l; apply Ropp_lt_contravar;
auto with real arith zarith.
unfold Rminus in |- *; apply Rplus_le_compat_l; apply Ropp_le_contravar;
apply Rle_powerRZ; auto with real arith zarith.
unfold Rminus in |- *; apply Rplus_le_compat_l; apply Ropp_le_contravar;
auto with real arith zarith.
ring_simplify (1 - 0)%R; rewrite Rmult_1_r.
ring_simplify (1 + 0)%R; rewrite Rinv_1; rewrite Rmult_1_r.
rewrite powerRZ_Zopp; auto with real arith zarith.
replace (powerRZ radix 2%nat) with (INR 4);
[ idtac | simpl in |- *; ring; auto with arith real zarith ].
replace (4%nat - 1)%R with (INR 3);
[ idtac | simpl in |- *; ring; auto with arith real zarith ].
replace (1 - / 4%nat)%R with (3%nat × / 4%nat)%R.
rewrite Rinv_mult_distr; auto with real arith zarith.
rewrite Rinv_involutive; auto with real arith zarith.
apply Rle_lt_trans with (4%nat × / 4%nat × (/ 3%nat × / 3%nat))%R;
[ right; ring | idtac ].
rewrite Rinv_r; auto with real arith zarith; rewrite Rmult_1_l.
apply Rlt_le_trans with (/ 3%nat × 1)%R;
[ apply Rmult_lt_compat_l; auto with arith real zarith | right; ring ].
apply Rmult_lt_reg_l with (INR 3); auto with arith real.
rewrite Rinv_r; auto with arith real; rewrite Rmult_1_r; auto with arith real.
simpl; field; auto with real.
Qed.
End AxpyAux.
Section Axpy.
Let radix := 2%Z.
Let FtoRradix := FtoR radix.
Coercion FtoRradix : float >-> R.
Variable b : Fbound.
Variable precision : nat.
Hypothesis precisionGreaterThanOne : 1 < precision.
Hypothesis pGivesBound : Zpos (vNum b) = Zpower_nat radix precision.
Theorem Axpy_tFlessu :
∀ (a1 x1 y1 : R) (a x y t u : float),
Fbounded b a →
Fbounded b x →
Fbounded b y →
Fbounded b t →
Fbounded b u →
Closest b radix (a × x) t →
Closest b radix (t + y) u →
Fcanonic radix b u →
Fcanonic radix b t →
(4%nat × Rabs t ≤ Rabs u)%R →
(Rabs (y1 - y) + Rabs (a1 × x1 - a × x) <
/ 4%nat × Fulp b radix precision (FLess b precision u))%R →
MinOrMax radix b (a1 × x1 + y1) u.
intros a1 x1 y1 a x y t u Fa Fx Fy Ft Fu tDef uDef Cu Ct H1.
unfold FLess in |- *; case (Rcase_abs (FtoR 2 u)); fold radix in |- *;
intros H3 H2.
2: cut (0 ≤ u)%R;
[ intros H'; case H'; intros H4; clear H' | auto with real ].
2: apply AxpyPos with precision a x y t; auto.
apply MinOrMax_Fopp.
replace (- (a1 × x1 + y1))%R with (- a1 × x1 + - y1)%R; [ idtac | ring ].
apply AxpyPos with precision (Fopp a) x (Fopp y) (Fopp t);
fold radix FtoRradix in |- *; auto with real float.
replace (FtoRradix (Fopp a) × FtoRradix x)%R with (- (a × x))%R;
[ apply ClosestOpp; auto
| unfold FtoRradix in |- *; rewrite Fopp_correct; ring ].
replace (FtoRradix (Fopp t) + FtoRradix (Fopp y))%R with (- (t + y))%R;
[ apply ClosestOpp; auto
| unfold FtoRradix in |- *; repeat rewrite Fopp_correct; ring ].
unfold FtoRradix in |- *; rewrite Fopp_correct; auto with real.
unfold FtoRradix in |- *; repeat rewrite Fopp_correct;
repeat rewrite Rabs_Ropp; auto with real.
unfold FtoRradix in |- *; repeat rewrite Fopp_correct;
rewrite <- (Rabs_Ropp (- a1 × x1 - - FtoR radix a × FtoR radix x));
auto with real.
replace (- (- a1 × x1 - - FtoR radix a × FtoR radix x))%R with
(a1 × x1 - FtoRradix a × FtoRradix x)%R;
[ auto | unfold FtoRradix in |- *; ring ].
rewrite <- Rabs_Ropp.
replace (- (- y1 - - FtoR radix y))%R with (y1 - FtoRradix y)%R;
[ idtac | unfold FtoRradix in |- *; ring ].
rewrite FPredFopFSucc; auto with zarith; rewrite Fopp_Fopp; auto with zarith.
replace (Fulp b radix precision (Fopp (FSucc b radix precision u))) with
(Fulp b radix precision (FSucc b radix precision u));
auto.
unfold Fulp in |- *; rewrite Fnormalize_Fopp; simpl in |- *;
auto with real zarith.
apply MinOrMax3 with precision; auto with zarith; fold FtoRradix in |- ×.
cut (FtoRradix t = 0%R); [ intros H' | idtac ].
cut (FtoRradix y = 0%R); [ intros H5 | idtac ].
replace (a1 × x1 + y1 - FtoRradix u)%R with
(y1 - y + (a1 × x1 - a × x) + a × x)%R.
2: rewrite <- H4; rewrite H5; ring.
apply
Rle_lt_trans
with
(Rabs (y1 - FtoRradix y + (a1 × x1 - FtoRradix a × FtoRradix x)) +
Rabs (FtoRradix a × FtoRradix x))%R; [ apply Rabs_triang | idtac ].
apply
Rlt_le_trans
with
(/ 4%nat × Fulp b radix precision (FPred b radix precision u) +
Rabs (FtoRradix a × FtoRradix x))%R; auto with real.
apply Rplus_lt_compat_r; apply Rle_lt_trans with (2 := H2); apply Rabs_triang.
apply
Rle_trans
with
(/ 4%nat × Fulp b radix precision (FPred b radix precision u) +
/ 2%nat × Fulp b radix precision (FPred b radix precision u))%R;
[ apply Rplus_le_compat_l | idtac ].
replace (FtoRradix a × FtoRradix x)%R with (FtoRradix a × FtoRradix x - t)%R;
[ idtac | rewrite H'; ring ].
apply Rmult_le_reg_l with (INR 2); auto with real arith.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real arith.
ring_simplify (1 × Fulp b radix precision (FPred b radix precision u))%R;
apply Rle_trans with (Fulp b radix precision t).
unfold FtoRradix in |- *; apply ClosestUlp; auto with zarith.
unfold Fulp in |- *; apply Rle_powerRZ; auto with real zarith.
replace (Fexp (Fnormalize radix b precision t)) with (- dExp b)%Z.
cut (Fbounded b (Fnormalize radix b precision (FPred b radix precision u)));
[ intros H6; elim H6; auto | apply FnormalizeBounded; auto with zarith ].
apply FBoundedPred; auto with zarith.
replace (Fnormalize radix b precision t) with t.
replace t with (Float 0 (- dExp b)); [ simpl in |- *; auto | idtac ].
apply FcanonicUnique with (radix := radix) (precision := precision) (b := b);
auto with zarith.
right; repeat (split; simpl in |- *; auto with zarith).
fold FtoRradix in |- *; rewrite H'; unfold FtoRradix, FtoR in |- *;
simpl in |- *; ring.
apply FcanonicUnique with (radix := radix) (precision := precision) (b := b);
auto with zarith.
apply FnormalizeCanonic; auto with zarith.
apply sym_eq; apply FnormalizeCorrect; auto with zarith.
apply
Rle_trans
with
((/ 4%nat + / 2%nat) × Fulp b radix precision (FPred b radix precision u))%R;
[ right; ring | idtac ].
apply
Rle_trans with (1 × Fulp b radix precision (FPred b radix precision u))%R;
[ apply Rmult_le_compat_r | right; ring ].
unfold Fulp in |- *; auto with real zarith.
apply Rmult_le_reg_l with (INR 4); auto with real arith.
apply Rle_trans with 3%R;simpl;[right; field|idtac].
apply Rle_trans with (3+1)%R; auto with real; right; ring.
cut (ProjectorP b radix (Closest b radix));
[ unfold ProjectorP in |- *; intros V | idtac ].
replace 0%R with (FtoRradix (Float 0 (- dExp b)));
[ idtac | unfold FtoRradix, FtoR in |- *; simpl in |- *; ring ].
unfold FtoRradix in |- *; apply V; auto.
replace (Float 0 (- dExp b)) with u; fold FtoRradix in |- ×.
replace (FtoRradix y) with (FtoRradix t + FtoRradix y)%R; auto.
rewrite H'; ring.
apply FcanonicUnique with (radix := radix) (precision := precision) (b := b);
auto with zarith.
right; repeat (split; simpl in |- *; auto with zarith).
fold FtoRradix in |- *; rewrite <- H4; unfold FtoRradix, FtoR in |- *;
simpl in |- *; ring.
apply RoundedProjector; auto.
apply ClosestRoundedModeP with precision; auto with zarith.
cut (∀ z : R, (Rabs z ≤ 0)%R → z = 0%R); [ intros V; apply V | idtac ].
apply Rmult_le_reg_l with (INR 4); auto with real arith.
apply Rle_trans with (1 := H1); right; rewrite <- H4; auto with real.
ring_simplify (4%nat×0)%R; apply Rabs_R0.
intros z V; case (Req_dec 0 z); auto with real.
intros V'; Contradict V.
apply Rlt_not_le; apply Rabs_pos_lt; auto.
Qed.
Theorem Axpy_opt :
∀ (a1 x1 y1 : R) (a x y t u : float),
(Fbounded b a) → (Fbounded b x) → (Fbounded b y) →
(Fbounded b t) → (Fbounded b u) →
(Closest b radix (a × x) t) →
(Closest b radix (t + y) u) →
(Fcanonic radix b u) → (Fcanonic radix b t) →
((5%nat + 4%nat × (powerRZ radix (- precision))) ×
/ (1 - powerRZ radix (- precision)) ×
(Rabs (a × x) + (powerRZ radix (Zpred (- dExp b))))
≤ Rabs y)%R →
(Rabs (y1 - y) + Rabs (a1 × x1 - a × x) ≤
(powerRZ radix (Zpred (Zpred (- precision)))) ×
(1 - powerRZ radix (Zsucc (- precision))) × Rabs y +
- (powerRZ radix (Zpred (Zpred (- precision))) × Rabs (a × x)) +
- powerRZ radix (Zpred (Zpred (- dExp b))))%R →
(MinOrMax radix b (a1 × x1 + y1) u).
intros a1 x1 y1 a x y t u Fa Fx Fy Ft Fu tDef uDef Cu Ct H1 H2.
apply Axpy_tFlessu with a x y t; auto.
cut (0 < 1 - powerRZ radix (- precision))%R; [ intros H'1 | idtac ].
2: apply Rplus_lt_reg_r with (powerRZ radix (- precision)).
2: ring_simplify (powerRZ radix (- precision) + 0)%R.
2: ring_simplify (powerRZ radix (- precision) + (1 - powerRZ radix (- precision)))%R.
2: replace 1%R with (powerRZ radix 0);
[ auto with real zarith | simpl in |- *; auto ].
cut (0 < 1 + powerRZ radix (- precision))%R; [ intros H'2 | idtac ].
2: apply Rlt_le_trans with 1%R; auto with real.
2: apply Rle_trans with (1 + 0)%R; auto with real zarith.
cut ((5%nat + 4%nat × powerRZ radix (- precision)) × Rabs t ≤ Rabs y)%R;
[ intros H3 | idtac ].
apply Rplus_le_reg_l with (- Rabs (FtoRradix u))%R.
ring_simplify (- Rabs (FtoRradix u) + Rabs (FtoRradix u))%R.
apply
Rle_trans
with
(- ((Rabs y + - Rabs t) × / (1 + powerRZ radix (- precision))) +
4%nat × Rabs (FtoRradix t))%R;
[ apply Rplus_le_compat_r; apply Ropp_le_contravar | idtac ].
apply
Rle_trans
with
(Rabs (FtoRradix y + FtoRradix t) × / (1 + powerRZ radix (- precision)))%R.
apply Rmult_le_compat_r; auto with real zarith.
rewrite <- (Rabs_Ropp (FtoRradix t)).
apply Rle_trans with (Rabs (FtoRradix y) - Rabs (- FtoRradix t))%R;
[ right; ring | idtac ].
replace (FtoRradix y + FtoRradix t)%R with (FtoRradix y - - FtoRradix t)%R;
[ apply Rabs_triang_inv | ring ].
case Cu; intros H4.
apply Rmult_le_reg_l with (1 + powerRZ radix (- precision))%R; auto.
apply
Rle_trans
with
(Rabs (FtoRradix y + FtoRradix t) ×
((1 + powerRZ radix (- precision)) × / (1 + powerRZ radix (- precision))))%R;
[ right; ring; ring | idtac ].
rewrite Rinv_r; auto with real.
ring_simplify (Rabs (FtoRradix y + FtoRradix t) × 1)%R.
rewrite Rmult_plus_distr_r.
apply Rle_trans with (Rabs u + / radix × Fulp b radix precision u)%R.
apply Rplus_le_reg_l with (- Rabs u)%R.
apply Rle_trans with (/ radix × Fulp b radix precision u)%R;
[idtac|right; ring].
apply Rle_trans with (Rabs (FtoRradix y + FtoRradix t - FtoRradix u)).
apply
Rle_trans with (Rabs (FtoRradix y + FtoRradix t) - Rabs (FtoRradix u))%R;
[ right; ring | apply Rabs_triang_inv ].
apply Rmult_le_reg_l with (IZR radix); auto with real zarith.
apply Rle_trans with (Fulp b radix precision u × (radix × / radix))%R;
[ rewrite Rinv_r; auto with real zarith | right; ring ].
ring_simplify (Fulp b radix precision u × 1)%R.
unfold FtoRradix in |- *; apply ClosestUlp; auto with real zarith.
rewrite Rplus_comm; auto.
ring_simplify (1×Rabs u)%R; apply Rplus_le_compat_l.
apply Rmult_le_reg_l with (IZR radix); auto with real zarith.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real zarith.
ring_simplify (1 × Fulp b radix precision u)%R.
apply
Rle_trans with (Rabs (FtoRradix u) × powerRZ radix (Zsucc (- precision)))%R.
unfold FtoRradix in |- *; apply FulpLe2; auto with zarith.
replace (Fnormalize radix b precision u) with u;
[ idtac
| apply
FcanonicUnique with (radix := radix) (precision := precision) (b := b) ];
auto with zarith real.
apply FnormalizeCanonic; auto with zarith.
apply sym_eq; apply FnormalizeCorrect; auto with zarith.
unfold Zsucc in |- *; rewrite powerRZ_add; auto with zarith real;
simpl in |- *; auto with zarith real; right; ring.
replace (FtoRradix y + FtoRradix t)%R with (FtoRradix u); auto with real.
apply Rle_trans with (Rabs (FtoRradix u) × 1)%R;
[ apply Rmult_le_compat_l | right; ring ]; auto with real.
pattern 1%R at 2 in |- *; replace 1%R with (/ 1)%R; auto with real zarith.
apply Rle_Rinv; auto with real zarith.
pattern 1%R at 1 in |- *; replace 1%R with (1 + 0)%R; auto with real zarith.
unfold FtoRradix in |- *; apply plusExact1 with b precision;
auto with real zarith.
rewrite Rplus_comm; auto.
elim H4; intros H5 H6; elim H6; intros H7 H8; rewrite H7; apply Zmin_Zle;
elim Fy; elim Ft; auto with zarith.
apply
Rle_trans
with
(((5%nat + 4%nat × powerRZ radix (- precision)) × Rabs (FtoRradix t) +
- Rabs (FtoRradix y)) × / (1 + powerRZ radix (- precision)))%R.
right; replace (INR 5) with (1 + 4%nat)%R;
[ idtac
| replace 1%R with (INR 1); [ rewrite <- plus_INR | idtac ];
auto with real arith ].
replace (4%nat × Rabs (FtoRradix t))%R with
(4%nat × Rabs (FtoRradix t) ×
((1 + powerRZ radix (- precision)) × / (1 + powerRZ radix (- precision))))%R.
ring; ring.
rewrite Rinv_r; auto with real; ring.
apply Rle_trans with (0 × / (1 + powerRZ radix (- precision)))%R;
[ apply Rmult_le_compat_r; auto with real | right; ring ].
apply Rplus_le_reg_l with (Rabs (FtoRradix y)).
apply
Rle_trans
with ((5%nat + 4%nat × powerRZ radix (- precision)) × Rabs (FtoRradix t))%R;
[ right; ring | ring_simplify (Rabs (FtoRradix y) + 0)%R ];
auto.
apply Rle_trans with (2 := H1).
rewrite Rmult_assoc.
apply Rmult_le_compat_l; auto with real arith.
apply Rle_trans with (5%nat + 4%nat × 0)%R; auto with real arith.
ring_simplify (5%nat + 4%nat × 0)%R; auto with real arith.
apply
Rle_trans
with
(Rabs (FtoRradix a × FtoRradix x) × / (1 - powerRZ radix (- precision)) +
powerRZ radix (Zpred (- dExp b)) × / (1 - powerRZ 2%Z (- precision)))%R;
[ idtac | right; simpl; ring ].
unfold FtoRradix in |- *; apply RoundLeGeneral; auto with real zarith.
apply Rle_lt_trans with (1 := H2).
apply UlpFlessuGe2 with t; auto.
Qed.
Theorem Axpy_Simpl1 :
∀ (a1 x1 y1 : R) (a x y t u : float),
4 ≤ precision →
Fbounded b a →
Fbounded b x →
Fbounded b y →
Fbounded b t →
Fbounded b u →
Closest b radix (a × x) t →
Closest b radix (t + y) u →
Fcanonic radix b u →
Fcanonic radix b t →
(6%nat × (Rabs (a × x) + powerRZ radix (Zpred (- dExp b))) ≤ Rabs y)%R →
(Rabs (y1 - y) + Rabs (a1 × x1 - a × x) ≤
powerRZ radix (Zpred (Zpred (- precision))) ×
(1 - powerRZ radix (Zsucc (- precision))) × Rabs y +
- (powerRZ radix (Zpred (Zpred (- precision))) × Rabs (a × x)) +
- powerRZ radix (Zpred (Zpred (- dExp b))))%R →
MinOrMax radix b (a1 × x1 + y1) u.
intros a1 x1 y1 a x y t u Enoughp Fa Fx Fy Ft Fu tDef uDef Cu Ct H1 H2.
apply Axpy_opt with a x y t; auto.
cut (0 < 1 - powerRZ radix (- precision))%R; [ intros H'1 | idtac ].
2: apply Rplus_lt_reg_r with (powerRZ radix (- precision)).
2: ring_simplify (powerRZ radix (- precision) + 0)%R.
2: ring_simplify (powerRZ radix (- precision) + (1 - powerRZ radix (- precision)))%R.
2: replace 1%R with (powerRZ radix 0);
[ auto with real zarith | simpl in |- *; auto ].
apply Rle_trans with (2 := H1).
apply Rmult_le_compat_r.
apply Rle_trans with (Rabs (FtoRradix a × FtoRradix x) + 0)%R;
auto with real zarith.
ring_simplify (Rabs (FtoRradix a × FtoRradix x) + 0)%R; auto with real.
apply Rmult_le_reg_l with (1 - powerRZ radix (- precision))%R; auto.
apply
Rle_trans
with
((5%nat + 4%nat × powerRZ radix (- precision)) ×
((1 - powerRZ radix (- precision)) × / (1 - powerRZ radix (- precision))))%R;
[ right; ring; ring | rewrite Rinv_r; auto with real ].
apply Rplus_le_reg_l with (- 5%nat + 6%nat × powerRZ radix (- precision))%R.
simpl; ring_simplify.
apply Rle_trans with (powerRZ radix 4%nat × powerRZ radix (- precision))%R;
[ apply Rmult_le_compat_r; auto with real arith zarith
| rewrite <- powerRZ_add; auto with real zarith ].
apply Rle_trans with (INR 10); [ right; simpl in |- *; ring | idtac ].
apply Rle_trans with (INR 16); [ idtac | right; simpl in |- *; ring ].
apply Rle_INR; auto with arith.
apply Rle_trans with (powerRZ radix 0);
[ idtac | right; simpl in |- *; ring ].
apply Rle_powerRZ; auto with zarith arith real.
apply Zle_trans with (4-precision)%Z; auto with zarith.
Qed.
Theorem Axpy_Simpl1bis :
∀ (a1 x1 y1 : R) (a x y t u : float),
4 ≤ precision →
Fbounded b a →
Fbounded b x →
Fbounded b y →
Fbounded b t →
Fbounded b u →
Closest b radix (a × x) t →
Closest b radix (t + y) u →
Fcanonic radix b u →
Fcanonic radix b t →
(6%nat × (Rabs (a × x) + powerRZ radix (Zpred (- dExp b))) ≤ Rabs y)%R →
(Rabs (y1 - y) + Rabs (a1 × x1 - a × x) ≤
powerRZ radix (Zpred (Zpred (- precision))) ×
(5%nat × / 6%nat - powerRZ radix (Zsucc (- precision))) ×
Rabs y + - powerRZ radix (Zpred (Zpred (- dExp b))))%R →
MinOrMax radix b (a1 × x1 + y1) u.
intros a1 x1 y1 a x y t u Enoughp Fa Fx Fy Ft Fu tDef uDef Cu Ct H1 H2.
apply Axpy_Simpl1 with a x y t; auto.
apply Rle_trans with (1 := H2).
replace (5%nat × / 6%nat)%R with (1 - / 6%nat)%R.
2: apply Rmult_eq_reg_l with (INR 6); auto with arith real.
2: apply trans_eq with (6%nat - 6%nat × / 6%nat)%R;
[ ring; ring | rewrite Rinv_r; auto with arith real ].
2: apply trans_eq with (5%nat × (6%nat × / 6%nat))%R;
[ rewrite Rinv_r; auto with arith real | ring ];
simpl in |- *; ring.
apply
Rplus_le_reg_l
with
(-
(powerRZ radix (Zpred (Zpred (- precision))) ×
(1 - powerRZ radix (Zsucc (- precision))) × Rabs (FtoRradix y)))%R.
ring_simplify.
unfold Rminus;apply Rplus_le_compat_r.
repeat rewrite Ropp_mult_distr_l_reverse; apply Ropp_le_contravar.
rewrite Rmult_assoc; apply Rmult_le_compat_l; auto with real zarith.
apply Rmult_le_reg_l with (INR 6); auto with arith real.
apply Rle_trans with (Rabs (FtoRradix y) × (6%nat × / 6%nat))%R;
[ rewrite Rinv_r; auto with arith real | right; ring ].
ring_simplify (Rabs (FtoRradix y) × 1)%R; apply Rle_trans with (2 := H1).
apply Rle_trans with (6%nat × (Rabs (FtoRradix a × FtoRradix x) + 0))%R;
auto with real zarith.
Qed.
Theorem Axpy_Simpl2 :
∀ (a1 x1 y1 : R) (a x y t u : float),
4 ≤ precision →
Fbounded b a →
Fbounded b x →
Fbounded b y →
Fbounded b t →
Fbounded b u →
Closest b radix (a × x) t →
Closest b radix (t + y) u →
Fcanonic radix b u →
Fcanonic radix b t →
(6%nat × (Rabs (a × x) + powerRZ radix (Zpred (- dExp b))) ≤ Rabs y)%R →
(Rabs (y1 - y) + Rabs (a1 × x1 - a × x) ≤
powerRZ radix (Zpred (Zpred (- precision))) × (2%nat × / 3%nat) × Rabs y +
- powerRZ radix (Zpred (Zpred (- dExp b))))%R →
MinOrMax radix b (a1 × x1 + y1) u.
intros a1 x1 y1 a x y t u Enoughp Fa Fx Fy Ft Fu tDef uDef Cu Ct H1 H2.
apply Axpy_Simpl1bis with a x y t; auto.
apply Rle_trans with (1 := H2).
apply Rplus_le_compat_r; apply Rmult_le_compat_r; auto with real.
apply Rmult_le_compat_l; auto with real zarith.
apply Rle_trans with (5%nat × / 6%nat - / 8%nat)%R.
apply Rmult_le_reg_l with (6%nat × 8%nat)%R;
[ apply Rmult_lt_0_compat; auto with real arith | idtac ].
apply
Rle_trans
with (6%nat × / 6%nat × (8%nat × 5%nat) - 8%nat × / 8%nat × 6%nat)%R;
[ idtac | right; ring ].
pattern (INR 6) at 1 in |- *; replace (INR 6) with (2%nat × 3%nat)%R;
[ idtac | rewrite <- mult_INR; auto with real arith ].
apply Rle_trans with (3%nat × / 3%nat × 8%nat × (2%nat × 2%nat))%R;
[ right; ring | idtac ].
repeat rewrite Rinv_r; auto with real arith.
apply Rle_trans with (INR 32);[simpl; right; ring|idtac].
apply Rle_trans with (INR 34);[auto with zarith real|simpl; right; ring].
unfold Rminus in |- *; apply Rplus_le_compat_l; apply Ropp_le_contravar.
apply Rle_trans with (powerRZ radix (Zsucc (- 4%nat)));
[ apply Rle_powerRZ; auto with real zarith | idtac ].
simpl; right; field; auto with real.
Qed.
End Axpy.
Section AxpyFmac.
Let radix := 2%Z.
Let FtoRradix := FtoR radix.
Coercion FtoRradix : float >-> R.
Variable b : Fbound.
Variable precision : nat.
Hypothesis precisionGreaterThanOne : 1 < precision.
Hypothesis pGivesBound : Zpos (vNum b) = Zpower_nat radix precision.
Theorem AxpyPos_Fmac :
∀ (a1 x1 y1 : R) (a x y u : float),
Fbounded b a →
Fbounded b x →
Fbounded b y →
Fbounded b u →
Closest b radix (a × x + y) u →
Fcanonic radix b u →
(0 < u)%R →
(Rabs (y1 - y) + Rabs (a1 × x1 - a × x) <
/ 2%nat × Fulp b radix precision (FPred b radix precision u))%R →
MinOrMax radix b (a1 × x1 + y1) u.
intros a1 x1 y1 a x y u Fa Fx Fy Fu uDef Cu Pu H.
cut
(Rabs (a1 × x1 + y1 - u) <
/ 2%nat × Fulp b radix precision (FPred b radix precision u) +
Rabs (a × x + y - u))%R; [ intros H4 | idtac ].
case (Rle_or_lt (a × x + y) u); intros H5.
apply MinOrMax1 with precision; auto with zarith arith.
apply Rlt_le_trans with (1 := H4).
replace 2%Z with radix; auto.
apply
Rplus_le_reg_l
with (- (/ 2%nat × Fulp b radix precision (FPred b radix precision u)))%R.
ring_simplify.
apply
Rle_trans
with (/ 2%nat × Fulp b radix precision (FPred b radix precision u))%R;
[ idtac | right ].
generalize uDef; unfold Closest in |- *; intros H6; elim H6; intros H7 H8;
clear H6 H7.
case
(Rle_or_lt (Rabs (FtoRradix a × FtoRradix x + FtoRradix y - FtoRradix u))
(/ 2%nat × Fulp b radix precision (FPred b radix precision u)));
auto; intros H6.
cut
(Rabs (FtoR radix u - (FtoRradix a × FtoRradix x + FtoRradix y)) ≤
Rabs
(FtoR radix (FPred b radix precision u) -
(FtoRradix a × FtoRradix x + FtoRradix y)))%R;
[ intros H7 | apply H8 ].
2: apply FBoundedPred; auto with zarith arith.
Contradict H7; apply Rlt_not_le.
apply
Rlt_le_trans
with (Rabs (FtoRradix a × FtoRradix x + FtoRradix y - FtoRradix u)).
2: right; rewrite <- Rabs_Ropp.
2: replace (- (FtoRradix a × FtoRradix x + FtoRradix y - FtoRradix u))%R with
(FtoR radix u - (FtoRradix a × FtoRradix x + FtoRradix y))%R;
auto with real; ring.
apply Rle_lt_trans with (2 := H6).
rewrite Faux.Rabsolu_left1.
apply
Rplus_le_reg_l
with
(u - (FtoRradix a × FtoRradix x + FtoRradix y) -
/ 2%nat × Fulp b radix precision (FPred b radix precision u))%R.
ring_simplify.
pattern (FtoRradix u) at 1 in |- *;
replace (FtoRradix u) with
(FtoRradix (FPred b radix precision u) +
Fulp b radix precision (FPred b radix precision u))%R;
[ idtac
| unfold FtoRradix in |- *; apply FpredUlpPos; auto with zarith arith ].
apply
Rle_trans
with
(Fulp b radix precision (FPred b radix precision u) -
Fulp b radix precision (FPred b radix precision u) × / 2%nat)%R;
[ right; fold FtoRradix; ring | idtac ].
apply
Rle_trans
with (/ 2%nat × Fulp b radix precision (FPred b radix precision u))%R.
right; simpl; field.
apply Rlt_le; apply Rlt_le_trans with (1 := H6).
rewrite Faux.Rabsolu_left1; [ right; ring | idtac ].
apply Rplus_le_reg_l with (FtoRradix u).
ring_simplify; auto with real.
fold FtoRradix in |- *;
apply
Rplus_le_reg_l with (Fulp b radix precision (FPred b radix precision u)).
apply
Rle_trans
with
(FtoRradix (FPred b radix precision u) +
Fulp b radix precision (FPred b radix precision u) +
- (FtoRradix a × FtoRradix x + FtoRradix y))%R;
[ right; ring
| unfold FtoRradix in |- *; rewrite FpredUlpPos; auto with zarith arith ].
ring_simplify (Fulp b radix precision (FPred b radix precision u) + 0)%R.
apply
Rle_trans
with (Rabs (FtoRradix u + - (FtoRradix a × FtoRradix x + FtoRradix y))).
apply RRle_abs.
rewrite <- Rabs_Ropp.
replace (- (FtoRradix u + - (FtoRradix a × FtoRradix x + FtoRradix y)))%R
with (FtoRradix a × FtoRradix x + FtoRradix y - FtoRradix u)%R;
[ idtac | ring ].
apply Rmult_le_reg_l with (INR 2); auto with arith real.
apply Rle_trans with (Fulp b radix precision u).
unfold FtoRradix in |- *; apply ClosestUlp; auto with arith zarith.
replace (INR 2) with (IZR radix); auto with zarith arith real.
apply FulpFPredLe; auto with zarith.
simpl; field.
case
(Rle_or_lt (a × x + y)
(u + / 2%nat × Fulp b radix precision (FPred b radix precision u)));
intros H6.
apply MinOrMax1 with precision; auto with arith zarith real float.
apply Rlt_le_trans with (1 := H4); rewrite Rabs_right.
apply
Rplus_le_reg_l
with
(u + - (/ 2%nat × Fulp b radix precision (FPred b radix precision u)))%R.
ring_simplify.
apply Rle_trans with (1 := H6).
right; simpl; field.
apply Rle_ge; apply Rplus_le_reg_l with (FtoRradix u); auto with real.
apply MinOrMax2 with precision; auto with arith zarith float real.
apply Rlt_le_trans with (1 := H4).
apply
Rle_trans
with
(/ 2%nat × Fulp b radix precision u +
Rabs (FtoRradix a × FtoRradix x + FtoRradix y - FtoRradix u))%R;
[ apply Rplus_le_compat_r | idtac ].
apply Rmult_le_compat_l; auto with real arith.
apply FulpFPredGePos; auto with arith zarith.
apply
Rle_trans
with
(/ 2%nat × Fulp b radix precision u + / 2%nat × Fulp b radix precision u)%R;
[ apply Rplus_le_compat_l | right ].
apply Rmult_le_reg_l with (INR 2); auto with real arith.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real arith.
ring_simplify (1 × Fulp b radix precision u)%R.
unfold FtoRradix in |- *; apply ClosestUlp; auto with real arith zarith float.
apply trans_eq with ((/ 2%nat + / 2%nat) × Fulp b radix precision u)%R;
[ ring | auto with real arith ].
replace (/ 2%nat + / 2%nat)%R with 1%R; [ ring | auto with real arith ].
rewrite <- (Rinv_r (INR 2)); auto with real arith; simpl in |- *; ring.
replace 2%Z with radix; auto; fold FtoRradix in |- *; auto with real.
apply Rlt_le;
apply
Rplus_lt_reg_r
with (/ 2%nat × Fulp b radix precision (FPred b radix precision u))%R.
rewrite Rplus_comm; apply Rlt_le_trans with (1 := H6).
apply Rplus_le_reg_l with (- (a1 × x1 + y1))%R.
ring_simplify
(- (a1 × x1 + y1) +
(/ 2%nat × Fulp b radix precision (FPred b radix precision u) +
(a1 × x1 + y1)))%R.
apply
Rle_trans
with (Rabs (- (a1 × x1 + y1) + (FtoRradix a × FtoRradix x + FtoRradix y)));
[ apply RRle_abs | idtac ].
rewrite <- Rabs_Ropp.
replace (- (- (a1 × x1 + y1) + (FtoRradix a × FtoRradix x + FtoRradix y)))%R
with (y1 - FtoRradix y + (a1 × x1 - FtoRradix a × FtoRradix x))%R;
[ auto with real | ring ].
apply
Rle_trans
with
(Rabs (y1 - FtoRradix y) + Rabs (a1 × x1 - FtoRradix a × FtoRradix x))%R;
[ apply Rabs_triang | auto with real ].
replace (a1 × x1 + y1 - FtoRradix u)%R with
(y1 - y + (a1 × x1 - a × x) + (a × x + y - u))%R;
[ idtac | ring ].
apply
Rle_lt_trans
with
(Rabs (y1 - FtoRradix y + (a1 × x1 - FtoRradix a × FtoRradix x)) +
Rabs (FtoRradix a × FtoRradix x + FtoRradix y - FtoRradix u))%R;
[ apply Rabs_triang | idtac ].
apply Rplus_lt_compat_r; auto.
apply Rle_lt_trans with (2 := H); apply Rabs_triang.
Qed.
Theorem AxpyFLessu_Fmac :
∀ (a1 x1 y1 : R) (a x y u : float),
Fbounded b a →
Fbounded b x →
Fbounded b y →
Fbounded b u →
Closest b radix (a × x + y) u →
Fcanonic radix b u →
(Rabs (y1 - y) + Rabs (a1 × x1 - a × x) <
/ 2%nat × Fulp b radix precision (FLess b precision u))%R →
MinOrMax radix b (a1 × x1 + y1) u.
intros a1 x1 y1 a x y u Fa Fx Fy Fu uDef Cu.
unfold FLess in |- *; case (Rcase_abs (FtoR 2 u)); fold radix in |- *;
intros H3 H2.
2: cut (0 ≤ u)%R;
[ intros H'; case H'; intros H4; clear H' | auto with real ].
2: apply AxpyPos_Fmac with a x y; auto.
apply MinOrMax_Fopp.
replace (- (a1 × x1 + y1))%R with (- a1 × x1 + - y1)%R; [ idtac | ring ].
apply AxpyPos_Fmac with (Fopp a) x (Fopp y); fold radix FtoRradix in |- *;
auto with real float.
replace (FtoRradix (Fopp a) × FtoRradix x + FtoRradix (Fopp y))%R with
(- (a × x + y))%R;
[ apply ClosestOpp; auto
| unfold FtoRradix in |- *; repeat rewrite Fopp_correct;ring ].
unfold FtoRradix in |- *; rewrite Fopp_correct; auto with real.
unfold FtoRradix in |- *; repeat rewrite Fopp_correct; auto with real.
rewrite <- Rabs_Ropp;
rewrite <- (Rabs_Ropp (- a1 × x1 - - FtoR radix a × FtoR radix x)).
replace (- (- y1 - - FtoR radix y))%R with (y1 - FtoRradix y)%R;
[ idtac | fold FtoRradix; ring ].
replace (- (- a1 × x1 - - FtoR radix a × FtoR radix x))%R with
(a1 × x1 - FtoRradix a × FtoRradix x)%R;
[ auto | unfold FtoRradix in |- *; ring ].
rewrite FPredFopFSucc; auto with zarith; rewrite Fopp_Fopp; auto with zarith.
replace (Fulp b radix precision (Fopp (FSucc b radix precision u))) with
(Fulp b radix precision (FSucc b radix precision u));
auto.
unfold Fulp in |- *; rewrite Fnormalize_Fopp; simpl in |- *;
auto with real zarith.
apply MinOrMax3 with precision; auto with zarith.
cut (u = Float 0 (- dExp b)); [ intros H1 | idtac ].
2: apply
FcanonicUnique with (radix := radix) (precision := precision) (b := b);
auto with zarith.
2: right; repeat (split; simpl in |- *; auto with zarith).
2: fold FtoRradix in |- *; rewrite <- H4; unfold FtoRradix, FtoR in |- *;
simpl in |- *; ring.
replace 2%Z with radix; auto; fold FtoRradix in |- ×.
replace (a1 × x1 + y1 - FtoRradix u)%R with
(y1 - y + (a1 × x1 - a × x) + (FtoRradix a × FtoRradix x + FtoRradix y - u))%R;
[ idtac | ring ].
apply
Rle_lt_trans
with
(Rabs (y1 - FtoRradix y + (a1 × x1 - FtoRradix a × FtoRradix x)) +
Rabs (FtoRradix a × FtoRradix x + FtoRradix y - FtoRradix u))%R;
[ apply Rabs_triang | idtac ].
apply
Rlt_le_trans
with
(/ 2%nat × Fulp b radix precision (FPred b radix precision u) +
Rabs (FtoRradix a × FtoRradix x + FtoRradix y - FtoRradix u))%R.
apply Rplus_lt_compat_r; auto.
apply Rle_lt_trans with (2 := H2); apply Rabs_triang.
apply
Rle_trans
with
(/ 2%nat × Fulp b radix precision (FPred b radix precision u) +
/ 2%nat × Fulp b radix precision u)%R.
apply Rplus_le_compat_l.
apply Rmult_le_reg_l with (INR 2); auto with real arith.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real arith.
ring_simplify (1 × Fulp b radix precision u)%R.
unfold FtoRradix in |- *; apply ClosestUlp; auto with zarith.
unfold Fulp in |- ×.
repeat rewrite FcanonicFnormalizeEq; auto with zarith float.
rewrite H1; rewrite FPredSimpl4; auto with zarith; simpl in |- ×.
2: cut (0 < pPred (vNum b))%Z;
[ idtac | apply pPredMoreThanOne with radix precision ];
auto with zarith.
2: cut (0 < nNormMin radix precision)%Z; [ idtac | apply nNormPos ];
auto with zarith.
right; apply trans_eq with ((/ 2 + / 2) × powerRZ 2 (- dExp b))%R;
[ ring | idtac ].
replace (/ 2 + / 2)%R with 1%R; [ ring | idtac ].
field; auto with real.
Qed.
Theorem Axpy_opt_Fmac :
∀ (a1 x1 y1 : R) (a x y u : float),
(Fbounded b a) → (Fbounded b x) → (Fbounded b y) → (Fbounded b u) →
(Closest b radix (a × x + y) u) →
(Fcanonic radix b u) →
(Rabs (y1 - y) + Rabs (a1 × x1 - a × x) <
(Rabs (a × x + y)) × (powerRZ radix (Zpred (- precision)) ×
(1 - powerRZ radix (Zsucc (- precision)))) -
powerRZ radix (Zpred (Zpred (- dExp b))) ×
(3 × / (powerRZ radix precision - powerRZ radix (- precision))))%R →
(MinOrMax radix b (a1 × x1 + y1) u).
intros a1 x1 y1 a x y u Fa Fx Fy Fu uDef Cu H.
apply AxpyFLessu_Fmac with a x y; auto.
cut (0 < powerRZ radix precision - 1)%R; [ intros H'1 | idtac ].
2: apply Rplus_lt_reg_r with 1%R.
2: ring_simplify.
2: replace 1%R with (powerRZ radix 0); [ idtac | simpl in |- × ];
auto with real zarith.
cut (0 < powerRZ radix precision - powerRZ radix (- precision))%R;
[ intros H'2 | idtac ].
2: apply Rplus_lt_reg_r with (powerRZ radix (- precision)).
2: ring_simplify; auto with real zarith.
apply Rlt_le_trans with (1 := H).
apply
Rle_trans
with
(/ 2%nat ×
(/ (powerRZ radix precision - 1) ×
(Rabs u × (1 - powerRZ radix (Zsucc (- precision))) -
powerRZ radix (- dExp b))))%R.
apply
Rle_trans
with
(/ 2%nat ×
(/ (powerRZ radix precision - 1) ×
((Rabs (a × x + y) × / (1 + powerRZ radix (- precision)) -
powerRZ radix (Zpred (- dExp b)) ×
/ (1 + powerRZ radix (- precision))) ×
(1 - powerRZ radix (Zsucc (- precision))) -
powerRZ radix (- dExp b))))%R.
apply
Rle_trans
with
(/ 2%nat ×
(/ (powerRZ radix precision - 1) ×
(Rabs (FtoRradix a × FtoRradix x + FtoRradix y) ×
/ (1 + powerRZ radix (- precision)) ×
(1 - powerRZ radix (Zsucc (- precision))))) -
(/ 2%nat ×
(/ (powerRZ radix precision - 1) ×
((1 - powerRZ radix (Zsucc (- precision))) ×
(powerRZ radix (Zpred (- dExp b)) ×
/ (1 + powerRZ radix (- precision))))) +
/ 2%nat × (/ (powerRZ radix precision - 1) × powerRZ radix (- dExp b))))%R;
[ idtac | right; ring; ring ].
apply
Rle_trans
with
(/ 2%nat ×
(/ (powerRZ radix precision - 1) ×
(Rabs (FtoRradix a × FtoRradix x + FtoRradix y) ×
/ (1 + powerRZ radix (- precision)) ×
(1 - powerRZ radix (Zsucc (- precision))))) -
powerRZ radix (Zpred (Zpred (- dExp b))) ×
(3 × / (powerRZ radix precision - powerRZ radix (- precision))))%R;
unfold Rminus in |- *; [ apply Rplus_le_compat_r | apply Rplus_le_compat_l ].
apply
Rle_trans
with
(Rabs (FtoRradix a × FtoRradix x + FtoRradix y) ×
(/ 2%nat ×
(/ (powerRZ radix precision + -1) ×
(/ (1 + powerRZ radix (- precision)) ×
(1 + - powerRZ radix (Zsucc (- precision)))))))%R;
[ idtac | right; ring ].
apply Rmult_le_compat_l; auto with real.
apply
Rle_trans
with
(/ 2%nat ×
(/ (powerRZ radix precision + -1) × / (1 + powerRZ radix (- precision))) ×
(1 + - powerRZ radix (Zsucc (- precision))))%R;
[ idtac | right; ring ].
apply Rmult_le_compat_r; auto with real zarith.
apply Rplus_le_reg_l with (powerRZ radix (Zsucc (- precision))).
ring_simplify.
replace 1%R with (powerRZ radix 0); auto with real zarith.
replace (Zpred (- precision)) with (- Zsucc precision)%Z;
[ idtac | unfold Zsucc, Zpred in |- *; ring ].
rewrite <- Rinv_powerRZ; auto with real zarith.
rewrite <- Rinv_mult_distr; auto with real zarith.
rewrite <- Rinv_mult_distr; auto with real zarith.
apply Rle_Rinv; auto with real zarith.
apply Rmult_lt_0_compat; auto with real arith zarith.
apply Rmult_lt_0_compat; auto with real arith zarith.
apply Rlt_le_trans with (1 + 0)%R; auto with real arith zarith.
ring_simplify (1 + 0)%R; auto with real.
replace (INR 2) with (powerRZ radix 1); auto with real zarith.
ring_simplify
(powerRZ radix 1 ×
((powerRZ radix precision + -1) × (1 + powerRZ radix (- precision))))%R.
repeat rewrite <- powerRZ_add; auto with zarith real.
ring_simplify (1+ precision + -precision)%Z.
apply
Rle_trans
with (powerRZ radix (precision + 1) + - powerRZ radix (- precision + 1))%R;
[ right | idtac ]; auto with real zarith.
unfold Zminus; rewrite Zplus_comm; rewrite Zplus_comm
with (n:=(-precision)%Z); ring.
apply Rle_trans with (powerRZ radix (precision + 1) + -0)%R;
auto with real zarith.
right;ring.
apply Rmult_integral_contrapositive; split; auto with real zarith.
apply not_eq_sym; apply Rlt_not_eq; auto with real zarith.
apply Rlt_le_trans with (1 + 0)%R;
[ ring_simplify (1 + 0)%R; auto with real | auto with real zarith ].
apply not_eq_sym; apply Rlt_not_eq; auto with real zarith.
apply Rlt_le_trans with (1 + 0)%R;
[ ring_simplify (1 + 0)%R; auto with real | auto with real zarith ].
apply Ropp_le_contravar.
replace (powerRZ radix precision + -1)%R with (powerRZ radix precision - 1)%R;
[ idtac | ring ].
replace (powerRZ radix precision + - powerRZ radix (- precision))%R with
(powerRZ radix precision - powerRZ radix (- precision))%R;
[ idtac | ring ].
cut (0 < 1 + powerRZ radix (- precision))%R; [ intros H'3 | idtac ].
2: apply Rlt_le_trans with (1 + 0)%R;
[ ring_simplify (1 + 0)%R; auto with real | auto with real zarith ].
apply Rmult_le_reg_l with (powerRZ radix precision - 1)%R; auto with real.
apply
Rle_trans
with
(/ 2%nat ×
((powerRZ radix precision - 1) × / (powerRZ radix precision - 1) ×
((1 + - powerRZ radix (Zsucc (- precision))) ×
(powerRZ radix (Zpred (- dExp b)) ×
/ (1 + powerRZ radix (- precision))))) +
/ 2%nat ×
((powerRZ radix precision - 1) × / (powerRZ radix precision - 1) ×
powerRZ radix (- dExp b)))%R; [ right; ring | idtac ].
repeat rewrite Rinv_r; auto with real.
apply Rmult_le_reg_l with (1 + powerRZ radix (- precision))%R; auto with real.
apply
Rle_trans
with
(/ 2%nat ×
((1 + - powerRZ radix (Zsucc (- precision))) ×
(powerRZ radix (Zpred (- dExp b)) ×
((1 + powerRZ radix (- precision)) ×
/ (1 + powerRZ radix (- precision))))) +
/ 2%nat × ((1 + powerRZ radix (- precision)) × powerRZ radix (- dExp b)))%R;
[ right; ring; ring | idtac ].
rewrite Rinv_r; auto with real.
apply
Rle_trans
with
(powerRZ radix (Zpred (Zpred (- dExp b))) ×
(3 ×
((1 + powerRZ radix (- precision)) × (powerRZ radix precision - 1) ×
/ (powerRZ radix precision - powerRZ radix (- precision)))))%R;
[ idtac | right; ring ].
replace ((1 + powerRZ radix (- precision)) × (powerRZ radix precision - 1))%R
with (powerRZ radix precision - powerRZ radix (- precision))%R.
rewrite Rinv_r; auto with real.
right; unfold Zsucc, Zpred in |- ×.
repeat rewrite powerRZ_add; auto with real zarith.
replace (powerRZ radix (-1)) with (/ 2%nat)%R;
[ idtac | simpl in |- *; auto with real zarith ].
simpl; field.
ring_simplify; repeat rewrite <- powerRZ_add; auto with real zarith.
ring_simplify (precision + - precision)%Z; simpl in |- *; ring.
apply Rmult_le_compat_l; auto with real arith.
apply Rmult_le_compat_l; auto with real arith.
unfold Rminus in |- *; apply Rplus_le_compat_r.
apply Rmult_le_compat_r; auto with real.
apply Rplus_le_reg_l with (powerRZ radix (Zsucc (- precision))).
ring_simplify.
replace 1%R with (powerRZ radix 0); auto with real zarith.
fold
(Rabs (FtoRradix a × FtoRradix x + FtoRradix y) ×
/ (1 + powerRZ radix (- precision)) -
powerRZ radix (Zpred (- dExp b)) × / (1 + powerRZ radix (- precision)))%R
in |- ×.
unfold FtoRradix in |- *; apply RoundGeGeneral; auto.
apply Rmult_le_compat_l; auto with real arith.
apply
Rle_trans
with
(/ (powerRZ radix precision - 1) ×
(Rabs (FtoRradix u) - Fulp b radix precision u))%R.
apply Rmult_le_compat_l; auto with real arith.
apply
Rplus_le_reg_l
with
(Fulp b radix precision u +
-
(Rabs (FtoRradix u) × (1 - powerRZ radix (Zsucc (- precision))) -
powerRZ radix (- dExp b)))%R.
ring_simplify.
apply FulpLeGeneral; auto.
apply
Rle_trans
with (/ (powerRZ radix precision - 1) × Rabs (FLess b precision u))%R.
apply Rmult_le_compat_l; auto with real arith.
unfold FtoRradix in |- *; apply UlpFlessuGe_aux; auto.
apply Rmult_le_reg_l with (powerRZ radix precision - 1)%R; auto with real.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real.
ring_simplify (1 × Rabs (FtoRradix (FLess b precision u)))%R.
unfold FtoRradix in |- *; apply FulpGe; auto with zarith.
unfold FLess in |- *; case Rcase_abs; intros H3; auto with float zarith.
Qed.
End AxpyFmac.
Section AxpyMisc.
Let radix := 2%Z.
Let FtoRradix := FtoR radix.
Coercion FtoRradix : float >-> R.
Variable b : Fbound.
Variable precision : nat.
Hypothesis precisionGreaterThanOne : 1 < precision.
Hypothesis pGivesBound : Zpos (vNum b) = Zpower_nat radix precision.
Theorem TwoMoreThanOne : (1 < radix)%Z.
unfold radix in |- *; red in |- *; simpl in |- *; auto.
Qed.
Hint Resolve TwoMoreThanOne.
Theorem TwoMoreThanOneR : (1 < radix)%R.
replace 1%R with (IZR 1); auto with real.
Qed.
Hint Resolve TwoMoreThanOneR.
Theorem FulpLeGeneral :
∀ p : float,
Fbounded b p →
(Fulp b radix precision p ≤
Rabs (FtoRradix p) × powerRZ radix (Zsucc (- precision)) +
powerRZ radix (- dExp b))%R.
intros p Hp.
cut (Fcanonic radix b (Fnormalize radix b precision p));
[ intros H | apply FnormalizeCanonic; auto with arith zarith ].
case H; intros H1.
apply
Rle_trans with (Rabs (FtoR radix p) × powerRZ radix (Zsucc (- precision)))%R.
apply FulpLe2; auto with zarith.
apply Rle_trans with (Rabs p × powerRZ radix (Zsucc (- precision)) + 0)%R;
[ right; fold FtoRradix; ring | apply Rplus_le_compat_l; auto with real zarith ].
apply Rle_trans with (powerRZ radix (- dExp b)).
right; unfold Fulp in |- *; simpl in |- ×.
replace (Fexp (Fnormalize radix b precision p)) with (- dExp b)%Z;
[ auto with real zarith | elim H1; intuition ].
apply Rle_trans with (0 + powerRZ radix (- dExp b))%R;
[ right; ring | apply Rplus_le_compat_r ].
apply Rmult_le_pos; auto with real zarith.
Qed.
Theorem RoundLeGeneral :
∀ (p : float) (z : R),
Fbounded b p →
Closest b radix z p →
(Rabs p ≤
Rabs z × / (1 - powerRZ radix (- precision)) +
powerRZ radix (Zpred (- dExp b)) × / (1 - powerRZ radix (- precision)))%R.
intros p z Hp H.
cut (0 < 1 - powerRZ radix (- precision))%R; [ intros H1 | idtac ].
2: apply Rplus_lt_reg_r with (powerRZ radix (- precision)).
2: ring_simplify.
2: replace 1%R with (powerRZ radix 0);
[ auto with real zarith | simpl in |- *; auto ].
apply Rmult_le_reg_l with (1 - powerRZ radix (- precision))%R; auto.
ring_simplify ((1 - powerRZ radix (- precision)) × Rabs p)%R.
apply
Rle_trans
with
(Rabs z ×
((1 - powerRZ radix (- precision)) × / (1 - powerRZ radix (- precision))) +
powerRZ radix (Zpred (- dExp b)) ×
((1 - powerRZ radix (- precision)) × / (1 - powerRZ radix (- precision))))%R;
[ idtac | right; ring; ring ].
repeat rewrite Rinv_r; auto with real.
ring_simplify.
apply Rplus_le_reg_l with (- powerRZ radix (Zpred (- dExp b)))%R.
ring_simplify.
apply
Rle_trans
with
(Rabs p +
-
(Rabs p × powerRZ radix (- precision) + powerRZ radix (Zpred (- dExp b))))%R;
[ right; ring | idtac ].
apply Rle_trans with (Rabs p + - (/ radix × Fulp b radix precision p))%R.
apply Rplus_le_compat_l; apply Ropp_le_contravar.
apply Rmult_le_reg_l with (IZR radix); auto with real zarith.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real zarith.
ring_simplify (1 × Fulp b radix precision p)%R.
apply
Rle_trans
with
(Rabs p × (powerRZ radix 1 × powerRZ radix (- precision)) +
powerRZ radix 1 × powerRZ radix (Zpred (- dExp b)))%R;
[ idtac | right; simpl in |- *; ring ].
repeat rewrite <- powerRZ_add; auto with zarith real.
replace (1 + - precision)%Z with (Zsucc (- precision)); auto with zarith.
replace (1 + Zpred (- dExp b))%Z with (- dExp b)%Z; auto with zarith.
apply FulpLeGeneral; auto.
unfold Zpred in |- *; auto with zarith.
apply Rplus_le_reg_l with (- Rabs z + / radix × Fulp b radix precision p)%R.
ring_simplify.
apply Rle_trans with (Rabs (p - z)).
apply Rle_trans with (Rabs p - Rabs z)%R;
[ right; ring | apply Rabs_triang_inv ].
rewrite <- Rabs_Ropp.
replace (- (p - z))%R with (z - p)%R; [ idtac | ring ].
apply Rmult_le_reg_l with (IZR radix); auto with real zarith.
apply Rle_trans with (Fulp b radix precision p × (radix × / radix))%R;
[ rewrite Rinv_r; auto with real zarith | right; ring ].
ring_simplify (Fulp b radix precision p × 1)%R.
unfold FtoRradix in |- *; apply ClosestUlp; auto.
Qed.
Theorem RoundGeGeneral :
∀ (p : float) (z : R),
Fbounded b p →
Closest b radix z p →
(Rabs z × / (1 + powerRZ radix (- precision)) -
powerRZ radix (Zpred (- dExp b)) × / (1 + powerRZ radix (- precision)) ≤
Rabs p)%R.
intros p z Hp H.
cut (0 < 1 + powerRZ radix (- precision))%R; [ intros H1 | idtac ].
2: apply Rlt_le_trans with 1%R; auto with real.
2: apply Rle_trans with (1 + 0)%R; auto with real zarith.
apply Rmult_le_reg_l with (1 + powerRZ radix (- precision))%R; auto.
apply
Rle_trans
with
(Rabs z ×
((1 + powerRZ radix (- precision)) × / (1 + powerRZ radix (- precision))) -
powerRZ radix (Zpred (- dExp b)) ×
((1 + powerRZ radix (- precision)) × / (1 + powerRZ radix (- precision))))%R;
[ right; ring; ring | idtac ].
repeat rewrite Rinv_r; auto with real.
ring_simplify.
apply Rplus_le_reg_l with (powerRZ radix (Zpred (- dExp b))).
ring_simplify
(powerRZ radix (Zpred (- dExp b)) +
(Rabs z - powerRZ radix (Zpred (- dExp b))))%R.
apply
Rle_trans
with
(Rabs p +
(Rabs p × powerRZ radix (- precision) + powerRZ radix (Zpred (- dExp b))))%R;
[ idtac | right; ring ].
apply Rle_trans with (Rabs p + / radix × Fulp b radix precision p)%R.
apply Rplus_le_reg_l with (- Rabs p)%R.
ring_simplify (- Rabs p + (Rabs p + / radix × Fulp b radix precision p))%R.
apply Rle_trans with (Rabs (z - p)).
apply Rle_trans with (Rabs z - Rabs p)%R;
[ right; ring | apply Rabs_triang_inv ].
apply Rmult_le_reg_l with (IZR radix); auto with real zarith.
apply Rle_trans with (Fulp b radix precision p × (radix × / radix))%R;
[ rewrite Rinv_r; auto with real zarith | right; ring ].
ring_simplify (Fulp b radix precision p × 1)%R.
unfold FtoRradix in |- *; apply ClosestUlp; auto.
apply Rplus_le_compat_l.
apply Rmult_le_reg_l with (IZR radix); auto with real zarith.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real zarith.
ring_simplify (1 × Fulp b radix precision p)%R.
apply
Rle_trans
with
(Rabs p × (powerRZ radix 1 × powerRZ radix (- precision)) +
powerRZ radix 1 × powerRZ radix (Zpred (- dExp b)))%R;
[ idtac | right; simpl in |- *; ring ].
repeat rewrite <- powerRZ_add; auto with zarith real.
replace (1 + - precision)%Z with (Zsucc (- precision)); auto with zarith.
replace (1 + Zpred (- dExp b))%Z with (- dExp b)%Z; auto with zarith.
apply FulpLeGeneral; auto.
Qed.
Theorem ExactSum_Near :
∀ p q f : float,
Fbounded b p →
Fbounded b q →
Fbounded b f →
Closest b radix (p + q) f →
Fexp p = (- dExp b)%Z →
(Rabs (p + q - f) < powerRZ radix (- dExp b))%R → (p + q)%R = f.
intros p q f Hp Hq Hf H H12 H1.
case
errorBoundedPlus
with
(b := b)
(radix := radix)
(precision := precision)
(p := p)
(q := q)
(pq := f); auto with zarith.
intros x0 H'; elim H'; intros H2 H'1; elim H'1; intros H3 H4; clear H' H'1.
apply Rplus_eq_reg_l with (- FtoRradix f)%R; ring_simplify (-f+f)%R.
apply trans_eq with (FtoR radix p + FtoR radix q - FtoR radix f)%R;
[ fold FtoRradix; ring | rewrite <- H2 ].
generalize H1; unfold FtoRradix in |- *; rewrite <- H2; intros H5.
cut (∀ r : R, Rabs r = 0%R → r = 0%R);
[ intros V; apply V | auto with real ].
rewrite <- Fabs_correct; auto with zarith.
apply is_Fzero_rep1; unfold is_Fzero in |- ×.
cut (∀ z : Z, (0 ≤ z)%Z → (z < 1)%Z → z = 0%Z);
[ intros W; apply W | auto with zarith ].
unfold Fabs in |- *; simpl in |- *; apply Zle_ZERO_Zabs.
replace 1%Z with (Fnum (Float 1 (- dExp b))); [ idtac | simpl in |- *; auto ].
apply Rlt_Fexp_eq_Zlt with (radix := radix); auto with zarith real float.
rewrite Fabs_correct; auto with zarith; apply Rlt_le_trans with (1 := H5).
right; unfold FtoR in |- *; simpl in |- *; ring.
unfold Fabs in |- *; simpl in |- ×.
rewrite H4; rewrite H12; apply Zmin_le1; auto with zarith.
elim Hq; auto with zarith.
intros m; case (Rle_or_lt 0 m); intros H'.
rewrite Rabs_right; auto with real.
rewrite Faux.Rabsolu_left1; auto with real; intros H'1.
rewrite <- (Ropp_involutive m); rewrite H'1; auto with real.
Qed.
End AxpyMisc.
Section AxpyAux.
Add Field RField : Rfield (completeness Zeq_bool_complete).
Let radix := 2%Z.
Let FtoRradix := FtoR radix.
Coercion FtoRradix : float >-> R.
Variable b : Fbound.
Variable precision : nat.
Hypothesis precisionGreaterThanOne : 1 < precision.
Hypothesis pGivesBound : Zpos (vNum b) = Zpower_nat radix precision.
Variables a1 x1 y1 : R.
Variables a x y t u r : float.
Hypothesis Fa : Fbounded b a.
Hypothesis Fx : Fbounded b x.
Hypothesis Fy : Fbounded b y.
Hypothesis Ft : Fbounded b t.
Hypothesis Fu : Fbounded b u.
Hypothesis tDef : Closest b radix (a × x) t.
Hypothesis uDef : Closest b radix (t + y) u.
Hypothesis rDef : isMin b radix (a1 × x1 + y1) r.
Theorem Axpy_aux1 :
Fcanonic radix b u →
(Rabs (FtoRradix a × FtoRradix x - FtoRradix t) ≤
/ 4%nat × Fulp b radix precision (FPred b radix precision u))%R →
(0 < u)%R →
(4%nat × Rabs t ≤ Rabs u)%R →
(Rabs (y1 - y) + Rabs (a1 × x1 - a × x) <
/ 4%nat × Fulp b radix precision (FPred b radix precision u))%R →
MinOrMax radix b (a1 × x1 + y1) u.
intros H Hblop H1 H2 H3.
cut
(Rabs (a1 × x1 + y1 - u) <
/ 2%nat × Fulp b radix precision (FPred b radix precision u) +
Rabs (t + y - u))%R; [ intros H4 | idtac ].
case (Rle_or_lt (t + y) u); intros H5.
apply MinOrMax1 with precision; auto with zarith arith.
apply Rlt_le_trans with (1 := H4).
apply
Rplus_le_reg_l
with (- (/ 2%nat × Fulp b radix precision (FPred b radix precision u)))%R.
ring_simplify
(- (/ 2%nat × Fulp b radix precision (FPred b radix precision u)) +
(/ 2%nat × Fulp b radix precision (FPred b radix precision u) +
Rabs (FtoRradix t + FtoRradix y - FtoRradix u)))%R.
apply
Rle_trans
with (/ 2%nat × Fulp b radix precision (FPred b radix precision u))%R;
[ idtac | right ].
generalize uDef; unfold Closest in |- *; intros H6; elim H6; intros H7 H8;
clear H6 H7.
case
(Rle_or_lt (Rabs (FtoRradix t + FtoRradix y - FtoRradix u))
(/ 2%nat × Fulp b radix precision (FPred b radix precision u)));
auto; intros H6.
cut
(Rabs (FtoR radix u - (FtoRradix t + FtoRradix y)) ≤
Rabs (FtoR radix (FPred b radix precision u) - (FtoRradix t + FtoRradix y)))%R;
[ intros H7 | apply H8 ].
2: apply FBoundedPred; auto with zarith arith.
Contradict H7; apply Rlt_not_le.
apply Rlt_le_trans with (Rabs (FtoRradix t + FtoRradix y - FtoRradix u)).
2: right; rewrite <- Rabs_Ropp.
2: replace (- (FtoRradix t + FtoRradix y - FtoRradix u))%R with
(FtoR radix u - (FtoRradix t + FtoRradix y))%R;
auto with real; ring.
apply Rle_lt_trans with (2 := H6).
rewrite Faux.Rabsolu_left1.
apply
Rplus_le_reg_l
with
(u - (FtoRradix t + FtoRradix y) -
/ 2%nat × Fulp b radix precision (FPred b radix precision u))%R.
ring_simplify.
fold FtoRradix in |- *; pattern (FtoRradix u) at 1 in |- *;
replace (FtoRradix u) with
(FtoRradix (FPred b radix precision u) +
Fulp b radix precision (FPred b radix precision u))%R;
[ idtac
| unfold FtoRradix in |- *; apply FpredUlpPos; auto with zarith arith ].
apply
Rle_trans
with
(Fulp b radix precision (FPred b radix precision u) -
Fulp b radix precision (FPred b radix precision u) × / 2%nat)%R;
[ right; ring | idtac ].
apply
Rle_trans
with (/ 2%nat × Fulp b radix precision (FPred b radix precision u))%R.
right;
apply
trans_eq
with
((1 - / 2%nat) × Fulp b radix precision (FPred b radix precision u))%R;
[ ring | auto with real arith ].
replace (1 - / 2%nat)%R with (/ 2%nat)%R; auto with real arith.
rewrite <- (Rinv_r (INR 2)); auto with real arith.
simpl; field; auto with real.
apply Rlt_le; apply Rlt_le_trans with (1 := H6).
rewrite Faux.Rabsolu_left1; [ right; ring | idtac ].
apply Rplus_le_reg_l with (FtoRradix u).
ring_simplify; auto.
fold FtoRradix in |- *;
apply
Rplus_le_reg_l with (Fulp b radix precision (FPred b radix precision u)).
apply
Rle_trans
with
(FtoRradix (FPred b radix precision u) +
Fulp b radix precision (FPred b radix precision u) +
- (FtoRradix t + FtoRradix y))%R;
[ right; ring
| unfold FtoRradix in |- *; rewrite FpredUlpPos; auto with zarith arith ].
ring_simplify (Fulp b radix precision (FPred b radix precision u) + 0)%R.
apply Rle_trans with (Rabs (FtoRradix u + - (FtoRradix t + FtoRradix y))).
apply RRle_abs.
rewrite <- Rabs_Ropp.
replace (- (FtoRradix u + - (FtoRradix t + FtoRradix y)))%R with
(FtoRradix t + FtoRradix y - u)%R; [ idtac | ring ].
apply Rmult_le_reg_l with (INR 2); auto with arith real.
apply Rle_trans with (Fulp b radix precision u).
unfold FtoRradix in |- *; apply ClosestUlp; auto with arith zarith.
replace (INR 2) with (IZR radix); auto with zarith arith real.
apply FulpFPredLe; auto with zarith.
apply
trans_eq
with ((1 - / 2%nat) × Fulp b radix precision (FPred b radix precision u))%R;
[ idtac | ring ].
replace (1 - / 2%nat)%R with (/ 2%nat)%R; auto with real arith.
rewrite <- (Rinv_r (INR 2)); auto with real arith.
simpl; field; auto with real.
case
(Rle_or_lt (t + y)
(u + / 2%nat × Fulp b radix precision (FPred b radix precision u)));
intros H6.
apply MinOrMax1 with precision; auto with arith zarith real float.
apply Rlt_le_trans with (1 := H4); rewrite Rabs_right.
apply
Rplus_le_reg_l
with
(u + - (/ 2%nat × Fulp b radix precision (FPred b radix precision u)))%R.
ring_simplify.
apply Rle_trans with (1 := H6).
apply
Rle_trans
with
(FtoRradix u +
(Fulp b radix precision (FPred b radix precision u) +
- (Fulp b radix precision (FPred b radix precision u) × / 2%nat)))%R;
[ apply Rplus_le_compat_l; right | right; ring ].
apply
trans_eq
with ((1 - / 2%nat) × Fulp b radix precision (FPred b radix precision u))%R;
[ idtac | ring ].
replace (1 - / 2%nat)%R with (/ 2%nat)%R; auto with real arith.
rewrite <- (Rinv_r (INR 2)); auto with real arith.
simpl; field; auto with real.
apply Rle_ge; apply Rplus_le_reg_l with (FtoRradix u); auto with real.
apply MinOrMax2 with precision; auto with arith zarith float real.
apply Rlt_le_trans with (1 := H4).
apply
Rle_trans
with
(/ 2%nat × Fulp b radix precision u +
Rabs (FtoRradix t + FtoRradix y - FtoRradix u))%R;
[ apply Rplus_le_compat_r | idtac ].
apply Rmult_le_compat_l; auto with real arith.
apply FulpFPredGePos; auto with arith zarith.
apply
Rle_trans
with
(/ 2%nat × Fulp b radix precision u + / 2%nat × Fulp b radix precision u)%R;
[ apply Rplus_le_compat_l | right ].
apply Rmult_le_reg_l with (INR 2); auto with real arith.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real arith.
ring_simplify (1 × Fulp b radix precision u)%R.
unfold FtoRradix in |- *; apply ClosestUlp; auto with real arith zarith float.
apply trans_eq with ((/ 2%nat + / 2%nat) × Fulp b radix precision u)%R;
[ ring | auto with real arith ].
replace (/ 2%nat + / 2%nat)%R with 1%R; [ ring | auto with real arith ].
rewrite <- (Rinv_r (INR 2)); auto with real arith; simpl in |- *; ring.
cut (∀ z z' : R, (Rabs z ≤ z')%R → (- z' ≤ z)%R);
[ intros V | idtac ].
replace (a1 × x1 + y1)%R with
(a1 × x1 - a × x + (y1 - y) + (a × x - t + (t + y)))%R;
[ fold FtoRradix in |- × | ring ].
replace (FtoRradix u) with
(- (/ 4%nat × Fulp b radix precision (FPred b radix precision u)) +
(- (/ 4%nat × Fulp b radix precision (FPred b radix precision u)) +
(FtoRradix u +
/ 2%nat × Fulp b radix precision (FPred b radix precision u))))%R.
apply Rplus_le_compat.
apply V; auto with real.
apply
Rle_trans
with
(Rabs (y1 - FtoRradix y) + Rabs (a1 × x1 - FtoRradix a × FtoRradix x))%R;
auto with real.
rewrite Rplus_comm; apply Rabs_triang.
apply Rplus_le_compat.
apply V; auto with real.
auto with real.
apply
trans_eq
with
(FtoRradix u +
(- / 4%nat + (- / 4%nat + / 2%nat)) ×
Fulp b radix precision (FPred b radix precision u))%R;
[ ring | idtac ].
replace (- / 4%nat + (- / 4%nat + / 2%nat))%R with 0%R; [ ring | idtac ].
replace (/ 2%nat)%R with (2%nat × / 4%nat)%R; [ simpl in |- *; ring | idtac ].
replace (INR 4) with (2%nat × 2%nat)%R; auto with arith real zarith.
rewrite Rinv_mult_distr; auto with real arith.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real arith.
rewrite <- mult_INR; auto with real arith.
intros z z'; case (Rle_or_lt 0 z); intros P1 P2.
apply Rle_trans with (-0)%R; [ apply Ropp_le_contravar | simpl in |- × ];
auto with real.
apply Rle_trans with (Rabs z); auto with real.
replace (-0)%R with 0%R; auto with real.
apply Ropp_le_cancel; rewrite Ropp_involutive; rewrite <- Rabs_left;
auto with real.
replace (a1 × x1 + y1 - FtoRradix u)%R with
(y1 - y + (a1 × x1 - a × x) + (a × x - t + (t + y - u)))%R;
[ idtac | ring ].
apply
Rle_lt_trans
with
(Rabs (y1 - FtoRradix y + (a1 × x1 - FtoRradix a × FtoRradix x)) +
Rabs
(FtoRradix a × FtoRradix x - FtoRradix t +
(FtoRradix t + FtoRradix y - FtoRradix u)))%R;
[ apply Rabs_triang | idtac ].
apply
Rle_lt_trans
with
(Rabs (y1 - FtoRradix y + (a1 × x1 - FtoRradix a × FtoRradix x)) +
(Rabs (FtoRradix a × FtoRradix x - FtoRradix t) +
Rabs (FtoRradix t + FtoRradix y - FtoRradix u)))%R;
[ apply Rplus_le_compat_l; apply Rabs_triang | idtac ].
apply
Rle_lt_trans
with
(Rabs (y1 - FtoRradix y) + Rabs (a1 × x1 - FtoRradix a × FtoRradix x) +
(Rabs (FtoRradix a × FtoRradix x - FtoRradix t) +
Rabs (FtoRradix t + FtoRradix y - FtoRradix u)))%R.
apply Rplus_le_compat_r; apply Rabs_triang.
apply
Rlt_le_trans
with
(/ 4%nat × Fulp b radix precision (FPred b radix precision u) +
(Rabs (FtoRradix a × FtoRradix x - FtoRradix t) +
Rabs (FtoRradix t + FtoRradix y - FtoRradix u)))%R;
auto with real.
apply
Rle_trans
with
(/ 4%nat × Fulp b radix precision (FPred b radix precision u) +
(/ 4%nat × Fulp b radix precision (FPred b radix precision u) +
Rabs (FtoRradix t + FtoRradix y - FtoRradix u)))%R;
auto with real.
apply
Rle_trans
with
((/ 4%nat + / 4%nat) × Fulp b radix precision (FPred b radix precision u) +
Rabs (FtoRradix t + FtoRradix y - FtoRradix u))%R;
[ right; ring | apply Rplus_le_compat_r ].
replace (/ 4%nat + / 4%nat)%R with (/ 2%nat)%R; auto with real.
replace (/ 2%nat)%R with (2%nat × / 4%nat)%R; [ simpl in |- *; ring | idtac ].
replace (INR 4) with (2%nat × 2%nat)%R; auto with arith real zarith.
rewrite Rinv_mult_distr; auto with real arith.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real arith.
rewrite <- mult_INR; auto with real arith.
Qed.
Theorem Axpy_aux1_aux1 :
Fnormal radix b t →
Fcanonic radix b u →
(0 < u)%R →
(4%nat × Rabs t ≤ Rabs u)%R →
(Rabs (FtoRradix a × FtoRradix x - FtoRradix t) ≤
/ 4%nat × Fulp b radix precision (FPred b radix precision u))%R.
intros H1 H2 H3 H4.
cut (Fcanonic radix b t); [ intros H5 | left; auto ].
apply Rle_trans with (/ 2%nat × Fulp b radix precision t)%R.
apply Rmult_le_reg_l with (INR 2); auto with real arith.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real arith.
ring_simplify (1 × Fulp b radix precision t)%R.
unfold FtoRradix in |- *; apply ClosestUlp; auto with real arith zarith.
apply Rle_trans with (/ 2%nat × (/ 4%nat × Fulp b radix precision u))%R.
apply Rmult_le_compat_l; auto with real arith.
apply Rmult_le_reg_l with (INR 4); auto with real arith.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real arith.
ring_simplify (1 × Fulp b radix precision u)%R.
cut (Fbounded b (Float (Fnum t) (Zsucc (Zsucc (Fexp t)))));
[ intros H6 | idtac ].
2: elim H1; intros H7 H8; elim H7; intros H9 H10.
2: repeat (split; simpl in |- *; auto with zarith).
apply
Rle_trans
with (Fulp b radix precision (Float (Fnum t) (Zsucc (Zsucc (Fexp t))))).
right; unfold Fulp in |- ×.
replace (Fnormalize radix b precision t) with t;
[ idtac
| apply
FcanonicUnique with (radix := radix) (b := b) (precision := precision);
auto with real float zarith ].
2: apply sym_eq; apply FnormalizeCorrect; auto with real zarith.
replace
(Fnormalize radix b precision (Float (Fnum t) (Zsucc (Zsucc (Fexp t)))))
with (Float (Fnum t) (Zsucc (Zsucc (Fexp t)))).
replace (Fexp (Float (Fnum t) (Zsucc (Zsucc (Fexp t))))) with
(Zsucc (Zsucc (Fexp t))); [ idtac | simpl in |- *; auto ].
replace (Zsucc (Zsucc (Fexp t))) with (2 + Fexp t)%Z;
[ rewrite powerRZ_add | unfold Zsucc in |- × ]; auto with zarith real.
replace (powerRZ radix 2) with (INR 4); [ idtac | simpl in |- × ];
auto with real zarith; ring.
apply FcanonicUnique with (radix := radix) (b := b) (precision := precision);
auto with real float arith zarith.
2: apply sym_eq; apply FnormalizeCorrect; auto with real zarith.
elim H1; intros H7 H8.
left; repeat (split; simpl in |- *; auto with zarith).
rewrite
FulpFabs with b radix precision (Float (Fnum t) (Zsucc (Zsucc (Fexp t))));
auto with zarith.
rewrite FulpFabs with b radix precision u; auto with zarith.
apply LeFulpPos; auto with zarith real float.
unfold FtoRradix in |- *; rewrite Fabs_correct; auto with real zarith.
replace (FtoR radix (Fabs (Float (Fnum t) (Zsucc (Zsucc (Fexp t)))))) with
(Zabs (Fnum t) × powerRZ radix (Zsucc (Zsucc (Fexp t))))%R;
[ idtac | unfold FtoR in |- *; simpl in |- *; auto ].
replace (Zsucc (Zsucc (Fexp t))) with (2 + Fexp t)%Z;
[ rewrite powerRZ_add | unfold Zsucc in |- × ]; auto with zarith real.
replace (powerRZ radix 2) with (INR 4);
[ idtac | simpl in |- *; auto with real zarith; ring ].
rewrite Rmult_comm; rewrite Rmult_assoc.
replace (powerRZ radix (Fexp t) × Zabs (Fnum t))%R with (Rabs (FtoRradix t));
[ unfold FtoRradix in |- *; repeat rewrite Fabs_correct;
auto with real zarith
| idtac ].
unfold FtoRradix in |- *; rewrite <- Fabs_correct; auto with zarith;
unfold FtoR in |- *; simpl in |- *; ring.
apply Rmult_le_reg_l with (INR 2); auto with arith real.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with arith real.
ring_simplify (1 × (/ 4%nat × Fulp b radix precision u))%R.
apply
Rle_trans
with
(/ 4%nat × (2%nat × Fulp b radix precision (FPred b radix precision u)))%R;
[ apply Rmult_le_compat_l; auto with real arith | right; ring ].
replace (INR 2) with (IZR radix); auto with arith zarith real.
apply FulpFPredLe; auto with arith zarith real.
Qed.
Theorem Axpy_aux2 :
Fcanonic radix b u →
Fsubnormal radix b t →
(0 < u)%R →
FtoRradix u = (t + y)%R →
(Rabs (y1 - y) + Rabs (a1 × x1 - a × x) <
/ 4%nat × Fulp b radix precision (FPred b radix precision u))%R →
MinOrMax radix b (a1 × x1 + y1) u.
intros H1 H2 H3 H4 H5.
apply MinOrMax1 with precision; auto with zarith; fold FtoRradix in |- ×.
replace (a1 × x1 + y1 - FtoRradix u)%R with
(y1 - y + (a1 × x1 - a × x) + (a × x - t + (t + y - u)))%R;
[ idtac | ring ].
apply
Rle_lt_trans
with
(Rabs (y1 - FtoRradix y + (a1 × x1 - FtoRradix a × FtoRradix x)) +
Rabs
(FtoRradix a × FtoRradix x - FtoRradix t +
(FtoRradix t + FtoRradix y - FtoRradix u)))%R;
[ apply Rabs_triang | idtac ].
apply
Rle_lt_trans
with
(Rabs (y1 - FtoRradix y) + Rabs (a1 × x1 - FtoRradix a × FtoRradix x) +
Rabs
(FtoRradix a × FtoRradix x - FtoRradix t +
(FtoRradix t + FtoRradix y - FtoRradix u)))%R;
[ apply Rplus_le_compat_r; apply Rabs_triang | idtac ].
apply
Rle_lt_trans
with
(Rabs (y1 - FtoRradix y) + Rabs (a1 × x1 - FtoRradix a × FtoRradix x) +
(Rabs (FtoRradix a × FtoRradix x - FtoRradix t) +
Rabs (FtoRradix t + FtoRradix y - FtoRradix u)))%R;
[ apply Rplus_le_compat_l; apply Rabs_triang | idtac ].
apply
Rlt_le_trans
with
(/ 4%nat × Fulp b radix precision (FPred b radix precision u) +
(Rabs (FtoRradix a × FtoRradix x - FtoRradix t) +
Rabs (FtoRradix t + FtoRradix y - FtoRradix u)))%R;
[ apply Rplus_lt_compat_r; auto | idtac ].
apply
Rle_trans
with
(/ 4%nat × Fulp b radix precision (FPred b radix precision u) +
(/ 2%nat × powerRZ radix (- dExp b) + 0))%R.
apply Rplus_le_compat; auto with real.
apply Rplus_le_compat.
apply Rle_trans with (/ 2%nat × Fulp b radix precision t)%R;
[ apply Rmult_le_reg_l with (INR 2) | apply Rmult_le_compat_l ];
auto with real arith.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real arith.
ring_simplify (1 × Fulp b radix precision t)%R; unfold FtoRradix in |- *;
apply ClosestUlp; auto with zarith.
unfold Fulp in |- *; replace (Fnormalize radix b precision t) with t.
elim H2; intros H6 H7; elim H7; intros H8 H9; rewrite H8;
auto with zarith real.
apply FcanonicUnique with (radix := radix) (b := b) (precision := precision);
auto with real zarith float.
right; auto.
apply sym_eq; apply FnormalizeCorrect; auto with zarith.
rewrite H4; right.
ring_simplify (FtoRradix t + FtoRradix y - (FtoRradix t + FtoRradix y))%R;
apply Rabs_R0.
replace (/ 2%nat × powerRZ radix (- dExp b) + 0)%R with
(/ 2%nat × powerRZ radix (- dExp b))%R; [ idtac | ring ].
apply
Rle_trans
with
(/ 4%nat × Fulp b radix precision (FPred b radix precision u) +
/ 2%nat × Fulp b radix precision (FPred b radix precision u))%R;
[ apply Rplus_le_compat_l; apply Rmult_le_compat_l; auto with real arith
| idtac ].
unfold Fulp in |- *; apply Rle_powerRZ; auto with zarith real.
cut (Fbounded b (Fnormalize radix b precision (FPred b radix precision u)));
[ intros H6; elim H6; auto | idtac ].
apply FnormalizeBounded; auto with zarith arith.
apply FBoundedPred; auto with zarith arith.
apply
Rle_trans
with
((/ 4%nat + / 2%nat) × Fulp b radix precision (FPred b radix precision u))%R;
[ right; ring | idtac ].
apply
Rle_trans with (1 × Fulp b radix precision (FPred b radix precision u))%R;
[ apply Rmult_le_compat_r; auto with real zarith | right; ring ].
unfold Fulp in |- *; auto with real zarith.
apply Rmult_le_reg_l with (INR 4); auto with real arith.
apply Rle_trans with (4%nat × / 4%nat + 2%nat × (2%nat × / 2%nat))%R;
[ right; simpl; ring | idtac ].
replace (2%nat × 2%nat)%R with (INR (2 × 2)); auto with arith real.
repeat rewrite Rinv_r; auto with real arith.
simpl; ring_simplify; auto with real.
apply Rle_trans with (3+1)%R; auto with real.
right; ring.
Qed.
Theorem Axpy_aux1_aux2 :
Fsubnormal radix b t →
Fcanonic radix b u →
(0 < u)%R →
(Zsucc (- dExp b) ≤ Fexp (FPred b radix precision u))%Z →
(Rabs (FtoRradix a × FtoRradix x - FtoRradix t) ≤
/ 4%nat × Fulp b radix precision (FPred b radix precision u))%R.
intros H1 H2 H3 H4.
apply Rle_trans with (/ 2%nat × Fulp b radix precision t)%R;
[ apply Rmult_le_reg_l with (INR 2) | idtac ]; auto with real arith.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real arith.
ring_simplify (1 × Fulp b radix precision t)%R; unfold FtoRradix in |- *;
apply ClosestUlp; auto with zarith.
unfold Fulp in |- *; replace (Fnormalize radix b precision t) with t.
2: apply
FcanonicUnique with (radix := radix) (b := b) (precision := precision);
auto with real zarith float.
2: right; auto.
2: apply sym_eq; apply FnormalizeCorrect; auto with zarith.
replace (Fnormalize radix b precision (FPred b radix precision u)) with
(FPred b radix precision u).
2: apply
FcanonicUnique with (radix := radix) (b := b) (precision := precision);
auto with real zarith float.
2: apply sym_eq; apply FnormalizeCorrect; auto with zarith.
elim H1; intros H5 H6; elim H6; intros H7 H8; rewrite H7.
apply Rmult_le_reg_l with (INR 4); auto with real arith.
repeat rewrite <- Rmult_assoc.
pattern (INR 4) at 1 in |- *; replace (INR 4) with (2%nat × 2%nat)%R.
rewrite Rmult_comm; rewrite Rmult_assoc; repeat rewrite Rinv_r;
auto with real arith.
ring_simplify.
replace (INR 2) with (powerRZ radix 1);
[ rewrite <- powerRZ_add | simpl in |- × ]; auto with zarith real.
replace (- dExp b + 1)%Z with (Zsucc (- dExp b));
[ apply Rle_powerRZ | unfold Zsucc in |- × ]; auto with zarith real.
rewrite <- mult_INR; auto with arith.
Qed.
Theorem Axpy_aux1_aux3 :
Fsubnormal radix b t →
Fcanonic radix b u →
(0 < u)%R →
(Zsucc (- dExp b) ≤ Fexp (FPred b radix precision u))%Z →
(Rabs (FtoRradix a × FtoRradix x - FtoRradix t) ≤
/ 4%nat × Fulp b radix precision (FPred b radix precision u))%R.
intros H1 H2 H3 H4.
apply Rle_trans with (/ 2%nat × Fulp b radix precision t)%R.
apply Rmult_le_reg_l with (INR 2); auto with real arith;
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with arith real.
ring_simplify (1 × Fulp b radix precision t)%R; unfold FtoRradix in |- *;
apply ClosestUlp; auto with zarith.
apply Rmult_le_reg_l with (INR 4); auto with real arith;
repeat rewrite <- Rmult_assoc.
replace (INR 4) with (2%nat × 2%nat)%R;
[ rewrite (Rmult_assoc 2%nat 2%nat (/ 2%nat)); repeat rewrite Rinv_r;
auto with real arith
| rewrite <- mult_INR; auto with real arith ].
ring_simplify; unfold Fulp in |- ×.
replace (Fnormalize radix b precision t) with t;
[ elim H1; intros H9 H10; elim H10; intros H11 H12; rewrite H11 | idtac ].
2: apply
FcanonicUnique with (radix := radix) (b := b) (precision := precision);
auto with real zarith float.
2: right; auto.
2: apply sym_eq; apply FnormalizeCorrect; auto with zarith.
replace (Fnormalize radix b precision (FPred b radix precision u)) with
(FPred b radix precision u).
2: apply
FcanonicUnique with (radix := radix) (b := b) (precision := precision);
auto with real zarith float.
2: apply sym_eq; apply FnormalizeCorrect; auto with zarith.
replace (2%nat × powerRZ radix (- dExp b))%R with
(powerRZ radix (Zsucc (- dExp b)));
[ apply Rle_powerRZ | unfold Zsucc in |- × ]; auto with zarith real float.
rewrite powerRZ_add; auto with real zarith; simpl in |- *; ring.
Qed.
Theorem Axpy_aux3 :
Fcanonic radix b u →
Fsubnormal radix b t →
(0 < u)%R →
Fexp (FPred b radix precision u) = (- dExp b)%Z →
(Zsucc (- dExp b) ≤ Fexp u)%Z →
(Rabs (y1 - y) + Rabs (a1 × x1 - a × x) <
/ 4%nat × Fulp b radix precision (FPred b radix precision u))%R →
MinOrMax radix b (a1 × x1 + y1) u.
intros H1 H2 H3 H4 H5 H6.
case (Rle_or_lt u (a1 × x1 + y1)); intros H7.
apply MinOrMax2 with precision; auto with zarith float real;
fold FtoRradix in |- ×.
replace (a1 × x1 + y1 - FtoRradix u)%R with
(y1 - y + (a1 × x1 - a × x) + (a × x - t + (t + y - u)))%R;
[ idtac | ring ].
apply
Rlt_le_trans
with
(/ 4%nat × powerRZ radix (- dExp b) +
(/ 2%nat × powerRZ radix (- dExp b) + / 2%nat × Fulp b radix precision u))%R.
apply
Rle_lt_trans
with
(Rabs (y1 - FtoRradix y + (a1 × x1 - FtoRradix a × FtoRradix x)) +
Rabs
(FtoRradix a × FtoRradix x - FtoRradix t +
(FtoRradix t + FtoRradix y - FtoRradix u)))%R;
[ apply Rabs_triang | idtac ].
apply
Rlt_le_trans
with
(/ 4%nat × powerRZ radix (- dExp b) +
Rabs
(FtoRradix a × FtoRradix x - FtoRradix t +
(FtoRradix t + FtoRradix y - FtoRradix u)))%R;
[ apply Rplus_lt_compat_r | apply Rplus_le_compat_l ].
apply
Rle_lt_trans
with
(Rabs (y1 - FtoRradix y) + Rabs (a1 × x1 - FtoRradix a × FtoRradix x))%R;
[ apply Rabs_triang | idtac ].
apply Rlt_le_trans with (1 := H6); unfold Fulp in |- ×.
replace (Fexp (Fnormalize radix b precision (FPred b radix precision u)))
with (- dExp b)%Z; auto with real zarith.
replace (Fnormalize radix b precision (FPred b radix precision u)) with
(FPred b radix precision u); auto with zarith.
apply FcanonicUnique with (radix := radix) (b := b) (precision := precision);
auto with real zarith float.
apply sym_eq; apply FnormalizeCorrect; auto with zarith.
apply
Rle_trans
with
(Rabs (FtoRradix a × FtoRradix x - FtoRradix t) +
Rabs (FtoRradix t + FtoRradix y - FtoRradix u))%R;
[ apply Rabs_triang | apply Rplus_le_compat ].
apply Rmult_le_reg_l with (INR 2);
[ idtac | rewrite <- Rmult_assoc; rewrite Rinv_r ];
auto with arith real.
ring_simplify; apply Rle_trans with (Fulp b radix precision t);
[ unfold FtoRradix in |- *; apply ClosestUlp; auto with zarith
| unfold Fulp in |- × ].
replace (Fexp (Fnormalize radix b precision t)) with (- dExp b)%Z;
auto with real zarith.
replace (Fnormalize radix b precision t) with t; elim H2; auto with zarith.
intros;
apply FcanonicUnique with (radix := radix) (b := b) (precision := precision);
auto with real zarith float.
right; auto.
apply sym_eq; apply FnormalizeCorrect; auto with zarith.
apply Rmult_le_reg_l with (INR 2);
[ idtac | rewrite <- Rmult_assoc; rewrite Rinv_r ];
auto with arith real.
ring_simplify (1 × Fulp b radix precision u)%R; unfold FtoRradix in |- *;
apply ClosestUlp; auto with zarith.
apply Rmult_le_reg_l with (INR 4); auto with arith real.
apply Rle_trans with (powerRZ radix (-dExp b) × 3%nat
+ Fulp b radix precision u × 2%nat)%R.
right; simpl; field; auto with real.
apply
Rle_trans
with
(Fulp b radix precision u × 2%nat + Fulp b radix precision u × 2%nat)%R;
[ apply Rplus_le_compat_r | idtac ]; auto with real arith.
apply Rle_trans with (3%nat × powerRZ radix (- dExp b))%R; [ right;ring | idtac ].
apply Rle_trans with (4%nat × powerRZ radix (- dExp b))%R;
auto with arith zarith real.
unfold Fulp in |- *; replace (INR 2) with (powerRZ radix 1);
[ idtac | simpl in |- *; auto with zarith real ].
replace (INR 4) with (powerRZ radix 2);
[ idtac | simpl in |- *; auto with zarith real; ring ].
repeat rewrite <- powerRZ_add; auto with zarith real.
replace (2 + - dExp b)%Z with (Zsucc (- dExp b) + 1)%Z;
[ apply Rle_powerRZ | unfold Zsucc in |- × ]; auto with zarith real.
replace (Fnormalize radix b precision u) with u; auto with zarith arith.
apply FcanonicUnique with (radix := radix) (b := b) (precision := precision);
auto with real zarith float.
apply sym_eq; apply FnormalizeCorrect; auto with zarith.
right; replace (INR 4) with (2%nat + 2%nat)%R;
[ ring | rewrite <- plus_INR; auto with arith real ].
case (Rle_or_lt (t + y) u); intros H8.
apply Axpy_aux2; auto.
apply sym_eq; unfold FtoRradix in |- *; apply ExactSum_Near with b precision;
auto with zarith.
elim H2; intuition.
cut (∀ v w : R, (v ≤ w)%R → v ≠ w → (v < w)%R); [ intros V | idtac ].
2: intros V W H' H'1; case H'; auto with real.
2: intros H'2; absurd (V = W); auto with real.
apply V.
apply Rmult_le_reg_l with (INR 2); auto with real arith.
apply Rle_trans with (Fulp b radix precision u);
[ unfold FtoRradix in |- *; apply ClosestUlp; auto with zarith | idtac ].
replace (powerRZ 2%Z (- dExp b)) with
(Fulp b radix precision (FPred b radix precision u)).
replace (INR 2) with (IZR radix); auto with zarith real.
apply FulpFPredLe; auto with zarith.
unfold Fulp in |- *;
replace (Fnormalize radix b precision (FPred b radix precision u)) with
(FPred b radix precision u).
rewrite H4; auto with zarith real.
apply FcanonicUnique with (radix := radix) (b := b) (precision := precision);
auto with real zarith float.
apply sym_eq; apply FnormalizeCorrect; auto with zarith.
replace (powerRZ 2%Z (- dExp b)) with
(Fulp b radix precision (FPred b radix precision u)).
2: unfold Fulp in |- *;
replace (Fnormalize radix b precision (FPred b radix precision u)) with
(FPred b radix precision u).
2: rewrite H4; auto with zarith real.
rewrite Faux.Rabsolu_left1; auto with real.
2: apply Rplus_le_reg_l with (FtoRradix u); unfold FtoRradix, radix in |- ×.
2: ring_simplify; auto with real zarith.
case
(Req_dec (- (FtoR 2 t + FtoR 2 y - FtoR 2 u))
(Fulp b radix precision (FPred b radix precision u)));
auto.
intros W; Contradict uDef.
replace (FtoRradix t + FtoRradix y)%R with
(FtoRradix (FPred b radix precision u)); auto with float real.
cut (FtoRradix u ≠ FPred b radix precision u); [ intros V' | idtac ].
Contradict V'.
cut (ProjectorP b radix (Closest b radix));
[ unfold ProjectorP in |- *; intros W' | idtac ].
apply sym_eq; apply W'; auto with zarith float real.
apply RoundedProjector; apply ClosestRoundedModeP with precision;
auto with zarith.
apply not_eq_sym; apply Rlt_dichotomy_converse; left.
unfold FtoRradix in |- *; apply FPredLt; auto with zarith.
apply
Rplus_eq_reg_l with (Fulp b radix precision (FPred b radix precision u)).
rewrite Rplus_comm; unfold FtoRradix in |- *; rewrite FpredUlpPos;
auto with zarith.
rewrite <- W; unfold FtoRradix, radix in |- *; ring.
apply FcanonicUnique with (radix := radix) (b := b) (precision := precision);
auto with real zarith float.
apply sym_eq; apply FnormalizeCorrect; auto with zarith.
apply MinOrMax1 with precision; auto with zarith; fold FtoRradix in |- ×.
rewrite Faux.Rabsolu_left1.
2: apply Rplus_le_reg_l with (FtoRradix u).
2: ring_simplify; auto with real.
apply Rle_lt_trans with (y - y1 + (t - a1 × x1) - (t + y - u))%R;
[ right; ring | idtac ].
apply Rlt_le_trans with (Rabs (y - y1 + (t - a1 × x1))).
cut (∀ u v : R, (0 < v)%R → (u - v < Rabs u)%R);
[ intros V; apply V | idtac ].
apply Rplus_lt_reg_r with (FtoRradix u).
ring_simplify; auto with real.
intros v w H'1.
apply Rlt_le_trans with v; auto with real.
unfold Rminus in |- *; apply Rlt_le_trans with (v + -0)%R;
[ auto with real | right; ring ].
apply RRle_abs.
replace (FtoRradix y - y1 + (FtoRradix t - a1 × x1))%R with
(- (y1 - FtoRradix y + (a1 × x1 - a × x) + (a × x - t)))%R;
[ idtac | ring ].
rewrite Rabs_Ropp.
apply
Rle_trans
with
(Rabs (y1 - FtoRradix y + (a1 × x1 - FtoRradix a × FtoRradix x)) +
Rabs (FtoRradix a × FtoRradix x - FtoRradix t))%R;
[ apply Rabs_triang | idtac ].
apply Rmult_le_reg_l with (INR 2); auto with real arith.
apply
Rle_trans
with
(2%nat × Rabs (y1 - FtoRradix y + (a1 × x1 - FtoRradix a × FtoRradix x)) +
2%nat × Rabs (FtoRradix a × FtoRradix x - FtoRradix t))%R;
[ right; ring | idtac ].
apply
Rle_trans
with
(Fulp b radix precision (FPred b radix precision u) +
Fulp b radix precision (FPred b radix precision u))%R;
[ idtac | right; simpl; ring ].
apply Rplus_le_compat.
apply
Rle_trans
with
(2%nat × (/ 4%nat × Fulp b radix precision (FPred b radix precision u)))%R;
auto with real.
apply
Rle_trans
with
(2%nat ×
(Rabs (y1 - FtoRradix y) + Rabs (a1 × x1 - FtoRradix a × FtoRradix x)))%R;
[ apply Rmult_le_compat_l; auto with real arith; apply Rabs_triang | idtac ].
apply Rmult_le_compat_l; auto with real arith.
apply
Rle_trans with (1 × Fulp b radix precision (FPred b radix precision u))%R;
[ rewrite <- Rmult_assoc; apply Rmult_le_compat_r | right; ring ];
auto with real arith zarith.
unfold Fulp in |- *; auto with real float zarith.
rewrite Rmult_comm; apply Rmult_le_reg_l with (INR 4); auto with arith real;
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real arith.
simpl; ring_simplify; auto with real arith.
apply Rle_trans with (Fulp b radix precision t);
[ unfold FtoRradix in |- *; apply ClosestUlp; auto with zarith
| unfold Fulp in |- × ].
apply Rle_powerRZ; auto with zarith real.
replace (Fnormalize radix b precision t) with t;
[ elim H2; intros H9 H10; elim H10; intros H11 H12; rewrite H11 | idtac ].
2: apply sym_eq; apply FcanonicFnormalizeEq; auto with zarith; right; auto.
replace (Fnormalize radix b precision (FPred b radix precision u)) with
(FPred b radix precision u); [ rewrite H4; auto with zarith | idtac ].
apply sym_eq; apply FcanonicFnormalizeEq; auto with zarith float.
Qed.
Theorem AxpyPos :
Fcanonic radix b u →
Fcanonic radix b t →
(0 < u)%R →
(4%nat × Rabs t ≤ Rabs u)%R →
(Rabs (y1 - y) + Rabs (a1 × x1 - a × x) <
/ 4%nat × Fulp b radix precision (FPred b radix precision u))%R →
MinOrMax radix b (a1 × x1 + y1) u.
intros H1 H2 H3 H4 H5.
case H2; intros H6.
apply Axpy_aux1; auto; apply Axpy_aux1_aux1; auto.
cut (∀ z1 z2 : Z, (z1 ≤ z2)%Z → z1 = z2 ∨ (Zsucc z1 ≤ z2)%Z);
[ intros V | idtac ].
2: intros z1 z2 V; omega.
case (V (- dExp b)%Z (Fexp u));
[ elim Fu; intuition | intros H7 | intros H7 ].
apply Axpy_aux2; auto.
unfold FtoRradix in |- *; apply plusExact1 with b precision;
auto with zarith float.
rewrite <- H7; elim H6; intros H8 H9; elim H9; intros H10 H11; rewrite H10.
rewrite Zmin_le1; auto with zarith; elim Fy; auto.
case (V (- dExp b)%Z (Fexp (FPred b radix precision u))).
cut (Fbounded b (FPred b radix precision u)); [ intros H8 | idtac ];
auto with float.
apply FBoundedPred; auto with zarith.
intros H8; apply Axpy_aux3; auto.
intros H8; apply Axpy_aux1; auto.
apply Axpy_aux1_aux3; auto.
Qed.
Definition FLess (p : float) :=
match Rcase_abs p with
| left _ ⇒ FSucc b radix precision p
| right _ ⇒ FPred b radix precision p
end.
Theorem UlpFlessuGe_aux :
∀ p : float,
Fbounded b p →
Fcanonic radix b p →
(Rabs p - Fulp b radix precision p ≤ Rabs (FLess p))%R.
intros p H H1.
cut
(∀ p : float,
Fbounded b p →
Fcanonic radix b p →
(0 < p)%R →
(Rabs p - Fulp b radix precision p ≤ Rabs (FPred b radix precision p))%R);
[ intros V | idtac ].
unfold FLess in |- *; case (Rcase_abs p); intros H2.
rewrite <- Rabs_Ropp; unfold FtoRradix in |- *; rewrite <- Fopp_correct.
pattern p at 3 in |- *; rewrite <- Fopp_Fopp.
rewrite <- (Rabs_Ropp (FtoR radix (FSucc b radix precision (Fopp (Fopp p)))));
rewrite <- (Fopp_correct radix (FSucc b radix precision (Fopp (Fopp p)))).
rewrite <- FPredFopFSucc with b radix precision (Fopp p); auto with zarith.
replace (Fulp b radix precision p) with (Fulp b radix precision (Fopp p));
auto with float.
fold FtoRradix in |- *; apply V; auto with float.
unfold FtoRradix in |- *; rewrite Fopp_correct; auto with real.
unfold Fulp in |- *; rewrite Fnormalize_Fopp; auto with arith.
unfold Fopp in |- *; simpl in |- *; auto with real.
cut (0 ≤ p)%R; [ intros H3; case H3; intros H4 | auto with real ].
apply V; auto.
rewrite <- H4; rewrite Rabs_R0.
ring_simplify (0 - Fulp b radix precision p)%R; apply Rle_trans with 0%R;
auto with real zarith.
unfold Fulp in |- *; auto with real zarith.
intros q H2 H3 H4.
unfold FtoRradix in |- *; rewrite <- FpredUlpPos with b radix precision q;
auto with zarith; fold FtoRradix in |- ×.
apply
Rle_trans
with
(Rabs (FtoRradix (FPred b radix precision q)) +
Rabs (Fulp b radix precision (FPred b radix precision q)) -
Fulp b radix precision q)%R;
[ unfold Rminus in |- *; apply Rplus_le_compat_r; apply Rabs_triang
| idtac ].
rewrite (Rabs_right (Fulp b radix precision (FPred b radix precision q)));
[ idtac | apply Rle_ge; unfold Fulp in |- *; auto with real zarith ].
apply Rplus_le_reg_l with (Fulp b radix precision q).
ring_simplify
(Fulp b radix precision q +
(Rabs (FtoRradix (FPred b radix precision q)) +
Fulp b radix precision (FPred b radix precision q) -
Fulp b radix precision q))%R.
rewrite Rplus_comm; apply Rplus_le_compat_r.
apply FulpFPredGePos; auto with zarith real float.
Qed.
Theorem UlpFlessuGe :
Fcanonic radix b u →
(/
(4%nat × (powerRZ radix precision - 1) × (1 + powerRZ radix (- precision))) ×
((1 - powerRZ radix (Zsucc (- precision))) × Rabs y) +
-
(/
(4%nat × (powerRZ radix precision - 1) × (1 + powerRZ radix (- precision)) ×
(1 - powerRZ radix (- precision))) ×
((1 - powerRZ radix (Zsucc (- precision))) × Rabs (a × x))) +
-
(powerRZ radix (Zpred (- dExp b)) ×
(/ (2%nat × (powerRZ radix precision - 1)) +
/
(4%nat × (powerRZ radix precision - 1) ×
(1 + powerRZ radix (- precision)) × (1 - powerRZ radix (- precision))) ×
(1 - powerRZ radix (Zsucc (- precision))))) ≤
/ 4%nat × Fulp b radix precision (FLess u))%R.
intros Cu.
cut (0 < 1 - powerRZ radix (- precision))%R; [ intros H'1 | idtac ].
2: apply Rplus_lt_reg_r with (powerRZ radix (- precision)).
2: ring_simplify.
2: replace 1%R with (powerRZ radix 0);
[ auto with real zarith | simpl in |- *; auto ].
cut (0 < 1 + powerRZ radix (- precision))%R; [ intros H'2 | idtac ].
2: apply Rlt_le_trans with 1%R; auto with real.
2: apply Rle_trans with (1 + 0)%R; auto with real zarith.
cut (0 < powerRZ radix precision - 1)%R; [ intros H'3 | idtac ].
2: apply Rplus_lt_reg_r with 1%R.
2: ring_simplify.
2: replace 1%R with (powerRZ radix 0);
[ auto with real zarith | simpl in |- *; auto ].
cut (0 < 1 - powerRZ radix (Zsucc (- precision)))%R; [ intros H'4 | idtac ].
2: apply Rplus_lt_reg_r with (powerRZ radix (Zsucc (- precision))).
2: ring_simplify.
2: replace 1%R with (powerRZ radix 0);
[ auto with real zarith | simpl in |- *; auto ].
cut (0 < INR 4)%R; [ intros H'5 | auto with arith real ].
rewrite
(Rinv_mult_distr
(4%nat × (powerRZ radix precision - 1) ×
(1 + powerRZ radix (- precision))) (1 - powerRZ radix (- precision)))
; auto with real.
2: apply Rmult_integral_contrapositive; split; auto with real.
apply
Rle_trans
with
(/
(4%nat × (powerRZ radix precision - 1) ×
(1 + powerRZ radix (- precision))) ×
(1 - powerRZ radix (Zsucc (- precision))) ×
(Rabs y -
(Rabs (FtoRradix a × FtoRradix x) × / (1 - powerRZ radix (- precision)) +
powerRZ radix (Zpred (- dExp b)) × / (1 - powerRZ radix (- precision)))) +
-
(powerRZ radix (Zpred (- dExp b)) ×
/ (2%nat × (powerRZ radix precision - 1))))%R;
[ right; ring; ring | idtac ].
apply
Rle_trans
with
(/
(4%nat × (powerRZ radix precision - 1) ×
(1 + powerRZ radix (- precision))) ×
(1 - powerRZ radix (Zsucc (- precision))) ×
(Rabs (FtoRradix y) - Rabs t) +
-
(powerRZ radix (Zpred (- dExp b)) ×
/ (2%nat × (powerRZ radix precision - 1))))%R;
[ apply Rplus_le_compat_r | idtac ].
apply Rmult_le_compat_l.
apply Rlt_le; apply Rmult_lt_0_compat; auto with real.
apply Rinv_0_lt_compat; apply Rmult_lt_0_compat; auto with real;
apply Rmult_lt_0_compat; auto with real.
unfold Rminus in |- *; apply Rplus_le_compat_l; apply Ropp_le_contravar.
apply
Rle_trans
with
(Rabs (FtoRradix a × FtoRradix x) × / (1 - powerRZ 2%Z (- precision)) +
powerRZ 2%Z (Zpred (- dExp b)) × / (1 - powerRZ 2%Z (- precision)))%R;
[ idtac | right; unfold FtoRradix, radix, Rminus in |- *; auto with real ].
apply RoundLeGeneral; auto with zarith.
rewrite
(Rinv_mult_distr (4%nat × (powerRZ radix precision - 1))
(1 + powerRZ radix (- precision))); auto with real.
apply
Rle_trans
with
(/ (4%nat × (powerRZ radix precision - 1)) ×
(1 - powerRZ radix (Zsucc (- precision))) ×
(/ (1 + powerRZ radix (- precision)) ×
(Rabs (FtoRradix y) - Rabs (FtoRradix t))) +
-
(powerRZ radix (Zpred (- dExp b)) ×
/ (2%nat × (powerRZ radix precision - 1))))%R;
[ right; ring | idtac ].
apply
Rle_trans
with
(/ (4%nat × (powerRZ radix precision - 1)) ×
(1 - powerRZ radix (Zsucc (- precision))) × Rabs u +
-
(powerRZ radix (Zpred (- dExp b)) ×
/ (2%nat × (powerRZ radix precision - 1))))%R.
apply Rplus_le_compat_r; apply Rmult_le_compat_l.
apply Rlt_le; apply Rmult_lt_0_compat; auto with real.
apply Rinv_0_lt_compat; apply Rmult_lt_0_compat; auto with real.
apply
Rle_trans
with
(/ (1 + powerRZ radix (- precision)) × Rabs (FtoRradix y + FtoRradix t))%R.
apply Rmult_le_compat_l; auto with real.
rewrite <- (Rabs_Ropp (FtoRradix t)).
replace (FtoRradix y + FtoRradix t)%R with (FtoRradix y - - FtoRradix t)%R;
[ apply Rabs_triang_inv | ring ].
case Cu; intros H1.
apply Rmult_le_reg_l with (1 + powerRZ radix (- precision))%R; auto.
apply
Rle_trans
with
(Rabs (FtoRradix y + FtoRradix t) ×
((1 + powerRZ radix (- precision)) × / (1 + powerRZ radix (- precision))))%R;
[ right; ring; ring | idtac ].
rewrite Rinv_r; auto with real.
ring_simplify.
apply Rle_trans with (Rabs u + / radix × Fulp b radix precision u)%R.
apply Rplus_le_reg_l with (- Rabs u)%R.
ring_simplify.
apply Rle_trans with (Rabs (FtoRradix y + FtoRradix t - FtoRradix u)).
apply
Rle_trans with (Rabs (FtoRradix y + FtoRradix t) - Rabs (FtoRradix u))%R;
[ right; ring | apply Rabs_triang_inv ].
apply Rmult_le_reg_l with (IZR radix); auto with real zarith.
apply Rle_trans with (Fulp b radix precision u × (radix × / radix))%R;
[ rewrite Rinv_r; auto with real zarith | right; ring ].
ring_simplify (Fulp b radix precision u × 1)%R.
unfold FtoRradix in |- *; apply ClosestUlp; auto with real zarith.
rewrite Rplus_comm; auto.
rewrite Rplus_comm; apply Rplus_le_compat_r.
apply Rmult_le_reg_l with (IZR radix); auto with real zarith.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real zarith.
ring_simplify (1 × Fulp b radix precision u)%R.
apply
Rle_trans with (Rabs (FtoRradix u) × powerRZ radix (Zsucc (- precision)))%R.
unfold FtoRradix in |- *; apply FulpLe2; auto with zarith.
replace (Fnormalize radix b precision u) with u;
[ idtac
| apply
FcanonicUnique with (radix := radix) (precision := precision) (b := b) ];
auto with zarith real.
apply FnormalizeCanonic; auto with zarith.
apply sym_eq; apply FnormalizeCorrect; auto with zarith.
unfold Zsucc in |- *; rewrite powerRZ_add; auto with zarith real;
simpl in |- *; auto with zarith real; right; ring.
replace (FtoRradix y + FtoRradix t)%R with (FtoRradix u); auto with real.
apply Rle_trans with (1 × Rabs (FtoRradix u))%R;
[ apply Rmult_le_compat_r | right; ring ]; auto with real.
pattern 1%R at 2 in |- *; replace 1%R with (/ 1)%R; auto with real zarith.
apply Rle_Rinv; auto with real zarith.
pattern 1%R at 1 in |- *; replace 1%R with (1 + 0)%R; auto with real zarith.
unfold FtoRradix in |- *; apply plusExact1 with b precision;
auto with real zarith.
rewrite Rplus_comm; auto.
elim H1; intros H5 H6; elim H6; intros H7 H8; rewrite H7; apply Zmin_Zle;
elim Fy; elim Ft; auto with zarith.
apply
Rle_trans
with
(/ (4%nat × (powerRZ radix precision - 1)) ×
(Rabs (FtoRradix u) -
(Rabs u × powerRZ radix (Zsucc (- precision)) +
powerRZ radix (- dExp b))))%R.
right; simpl; ring_simplify.
unfold Zpred in |- *; rewrite powerRZ_add; auto with zarith real.
simpl; field; auto with real.
apply
Rle_trans
with
(/ (4%nat × (powerRZ radix precision - 1)) ×
(Rabs (FtoRradix u) - Fulp b radix precision u))%R.
apply Rmult_le_compat_l; auto with real.
apply Rlt_le; apply Rinv_0_lt_compat; apply Rmult_lt_0_compat; auto with real.
unfold Rminus in |- *; apply Rplus_le_compat_l; apply Ropp_le_contravar.
apply FulpLeGeneral; auto.
apply
Rle_trans
with (/ (4%nat × (powerRZ radix precision - 1)) × Rabs (FLess u))%R.
apply Rmult_le_compat_l; auto with real.
apply Rlt_le; apply Rinv_0_lt_compat; apply Rmult_lt_0_compat; auto with real.
apply UlpFlessuGe_aux; auto.
rewrite Rinv_mult_distr; auto with real.
rewrite Rmult_assoc; apply Rmult_le_compat_l; auto with real.
apply Rmult_le_reg_l with (powerRZ radix precision - 1)%R; auto.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real.
ring_simplify (1 × Rabs (FtoRradix (FLess u)))%R; unfold FtoRradix in |- *;
apply FulpGe; auto with zarith.
unfold FLess in |- *; case (Rcase_abs u); intros H; auto with float.
apply FBoundedSuc; auto with zarith.
apply FBoundedPred; auto with zarith.
Qed.
Theorem UlpFlessuGe2 :
Fcanonic radix b u →
(powerRZ radix (Zpred (Zpred (- precision))) ×
(1 - powerRZ radix (Zsucc (- precision))) × Rabs y +
- (powerRZ radix (Zpred (Zpred (- precision))) × Rabs (a × x)) +
- powerRZ radix (Zpred (Zpred (- dExp b))) <
/ 4%nat × Fulp b radix precision (FLess u))%R.
intros H.
apply
Rlt_le_trans
with
(/
(4%nat × (powerRZ radix precision - 1) ×
(1 + powerRZ radix (- precision))) ×
((1 - powerRZ radix (Zsucc (- precision))) × Rabs (FtoRradix y)) +
-
(/
(4%nat × (powerRZ radix precision - 1) ×
(1 + powerRZ radix (- precision)) × (1 - powerRZ radix (- precision))) ×
((1 - powerRZ radix (Zsucc (- precision))) ×
Rabs (FtoRradix a × FtoRradix x))) +
-
(powerRZ radix (Zpred (- dExp b)) ×
(/ (2%nat × (powerRZ radix precision - 1)) +
/
(4%nat × (powerRZ radix precision - 1) ×
(1 + powerRZ radix (- precision)) × (1 - powerRZ radix (- precision))) ×
(1 - powerRZ radix (Zsucc (- precision))))))%R;
[ idtac | apply UlpFlessuGe; auto ].
cut (0 < 1 - powerRZ radix (- precision))%R; [ intros H'1 | idtac ].
2: apply Rplus_lt_reg_r with (powerRZ radix (- precision)).
2: ring_simplify.
2: replace 1%R with (powerRZ radix 0);
[ auto with real zarith | simpl in |- *; auto ].
cut (0 < 1 + powerRZ radix (- precision))%R; [ intros H'2 | idtac ].
2: apply Rlt_le_trans with 1%R; auto with real.
2: apply Rle_trans with (1 + 0)%R; auto with real zarith.
cut (0 < powerRZ radix precision - 1)%R; [ intros H'3 | idtac ].
2: apply Rplus_lt_reg_r with 1%R.
2: ring_simplify.
2: replace 1%R with (powerRZ radix 0);
[ auto with real zarith | simpl in |- *; auto ].
cut (0 < 1 - powerRZ radix (Zsucc (- precision)))%R; [ intros H'4 | idtac ].
2: apply Rplus_lt_reg_r with (powerRZ radix (Zsucc (- precision))).
2: ring_simplify.
2: replace 1%R with (powerRZ radix 0);
[ auto with real zarith | simpl in |- *; auto ].
cut (0 < INR 4)%R; [ intros H'5 | auto with arith real ].
apply
Rle_lt_trans
with
(/
(4%nat × (powerRZ radix precision - 1) ×
(1 + powerRZ radix (- precision))) ×
((1 - powerRZ radix (Zsucc (- precision))) × Rabs (FtoRradix y)) +
-
(/
(4%nat × (powerRZ radix precision - 1) ×
(1 + powerRZ radix (- precision)) × (1 - powerRZ radix (- precision))) ×
((1 - powerRZ radix (Zsucc (- precision))) ×
Rabs (FtoRradix a × FtoRradix x))) +
- powerRZ radix (Zpred (Zpred (- dExp b))))%R;
[ apply Rplus_le_compat_r | apply Rplus_lt_compat_l ].
apply Rplus_le_compat.
repeat rewrite <- Rmult_assoc.
apply Rmult_le_compat_r; auto with real.
rewrite Rmult_assoc.
ring_simplify ((powerRZ radix precision - 1) × (1 + powerRZ radix (- precision)))%R.
rewrite <- powerRZ_add; auto with real zarith.
replace (powerRZ radix (- precision + precision)) with 1%R;
[ idtac
| ring_simplify (- precision + precision)%Z; simpl in |- *; auto with real zarith ].
ring_simplify (-1 + (- powerRZ radix (- precision) + (powerRZ radix precision + 1)))%R.
apply Rmult_le_compat_r; auto with real.
apply Rle_trans with (/ (4%nat × powerRZ radix precision))%R;
[ right | idtac ].
unfold Zpred in |- ×.
replace (- precision + -1 + -1)%Z with (- 2%nat + - precision)%Z;
[ rewrite powerRZ_add | ring ]; auto with zarith real.
rewrite Rinv_mult_distr; auto with real zarith.
replace (powerRZ radix (- 2%nat)) with (/ 4%nat)%R.
replace (/ powerRZ radix precision)%R with (powerRZ radix (- precision));
auto with real.
apply powerRZ_Zopp; auto with zarith real.
rewrite powerRZ_Zopp; auto with zarith real.
replace (powerRZ radix 2%nat) with (INR 4);
[ auto with zarith real | simpl in |- *; ring ].
apply Rle_Rinv; auto with real.
apply
Rlt_le_trans
with
(4%nat ×
(powerRZ radix precision - 1 + (1 - powerRZ radix (- precision))))%R.
apply Rle_lt_trans with (4%nat × (0 + 0))%R; [ right; ring | auto with real ].
ring_simplify (precision+-precision)%Z; simpl; right;ring.
apply Rmult_le_compat_l; auto with real zarith.
ring_simplify (precision+-precision)%Z.
apply
Rle_trans with (powerRZ radix precision + - powerRZ radix (- precision))%R;
[ right; simpl; ring | idtac ]; auto with real zarith.
apply Rle_trans with (powerRZ radix precision + -0)%R;
[ idtac | right; ring ]; auto with real zarith.
apply Ropp_le_contravar.
repeat rewrite <- Rmult_assoc.
apply Rmult_le_compat_r; auto with real.
unfold Zpred in |- ×.
replace (- precision + -1 + -1)%Z with (- 2%nat + - precision)%Z;
[ rewrite powerRZ_add | ring ]; auto with zarith real.
repeat rewrite Rmult_assoc; rewrite Rinv_mult_distr; auto with real zarith.
2: apply Rmult_integral_contrapositive; split; auto with real.
replace (powerRZ radix (- 2%nat)) with (/ 4%nat)%R.
2: rewrite powerRZ_Zopp; auto with zarith real.
2: replace (powerRZ radix 2%nat) with (INR 4);
[ auto with zarith real | simpl in |- *; ring ].
repeat rewrite Rmult_assoc; apply Rmult_le_compat_l; auto with real.
apply
Rmult_le_reg_l
with
((powerRZ radix precision - 1) ×
((1 + powerRZ radix (- precision)) × (1 - powerRZ radix (- precision))))%R;
auto with real.
apply Rmult_lt_0_compat; auto with real.
apply Rmult_lt_0_compat; auto with real.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real.
2: apply Rmult_integral_contrapositive; split; auto with real.
ring_simplify.
repeat rewrite Ropp_mult_distr_l_reverse.
repeat rewrite <- powerRZ_add; auto with real zarith.
replace (powerRZ radix (precision + - precision)) with 1%R;
[ idtac | ring_simplify ( precision +- precision)%Z; simpl in |- *; auto ].
ring_simplify ( precision + - precision + - precision + -precision)%Z.
apply
Rle_trans
with
(1 +
(- powerRZ radix (- precision) + (- powerRZ radix (- precision) + 0)))%R;
[ right | idtac ].
rewrite powerRZ_add; auto with real zarith; simpl; ring.
ring_simplify (- precision + - precision + - precision)%Z.
apply
Rle_trans
with
(1 +
(- powerRZ radix (- precision) +
(- powerRZ radix (-2 × precision) + powerRZ radix (-3 × precision))))%R;
[ idtac | right; ring ].
apply Rplus_le_compat_l; apply Rplus_le_compat_l.
apply Rplus_le_compat; auto with real zarith.
apply Ropp_lt_contravar.
apply Rmult_lt_reg_l with (powerRZ radix (Z_of_nat 2 + dExp b));
auto with real zarith.
rewrite <- Rmult_assoc; unfold Zpred in |- ×.
repeat rewrite <- powerRZ_add; auto with real zarith.
replace (2%nat + dExp b + (- dExp b + -1 + -1))%Z with 0%Z by ring.
replace (2%nat + dExp b + (- dExp b + -1))%Z with 1%Z by ring.
replace (powerRZ radix 0) with 1%R;
[ idtac | simpl in |- *; auto with real zarith ].
replace (powerRZ radix 1) with 2%R;
[ idtac | simpl in |- *; ring; auto with real zarith ].
rewrite Rmult_plus_distr_l.
apply Rlt_le_trans with (/ 3%nat + 2%nat × / 3%nat)%R.
2: simpl; right;field; auto with real.
apply
Rle_lt_trans
with
(/ 3%nat +
2 ×
(/
(4%nat × (powerRZ radix precision - 1) × (1 + powerRZ radix (- precision)) ×
(1 - powerRZ radix (- precision))) × (1 - powerRZ radix (- precision + 1))))%R;
[ apply Rplus_le_compat_r | apply Rplus_lt_compat_l ].
rewrite Rinv_mult_distr; auto with arith real.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with arith real.
rewrite Rmult_1_l; replace (INR 3) with (powerRZ radix 2%nat - 1)%R;
auto with real zarith arith.
apply Rle_Rinv; auto with real arith.
replace 0%R with (powerRZ radix 0 - 1)%R;
[ unfold Rminus in |- *; auto with real arith zarith | simpl in |- *; ring ].
unfold Rminus in |- *; auto with real arith zarith.
apply Rplus_le_compat_r; apply Rle_powerRZ; auto with real arith zarith.
simpl in |- *; ring; auto with arith zarith real.
apply Rmult_lt_compat_l; auto with real arith.
repeat rewrite Rinv_mult_distr; auto with real.
2: apply Rmult_integral_contrapositive; split; auto with real.
apply
Rle_lt_trans
with
(/ 4%nat × / (powerRZ radix 2%nat - 1) × / (1 + 0) ×
/ (1 - powerRZ radix (- 2%nat)) × (1 - 0))%R.
repeat rewrite Rmult_assoc.
apply Rmult_le_compat_l; auto with real arith.
apply Rmult_le_compat; auto with real arith zarith.
apply Rmult_le_pos; auto with real.
apply Rmult_le_pos; auto with real.
apply Rle_Rinv; auto with real arith zarith.
apply Rle_lt_trans with (powerRZ radix 0 - 1)%R;
[ right; simpl in |- *; ring | auto with real arith zarith ].
unfold Rminus in |- *; apply Rplus_lt_compat_r; auto with real arith zarith.
unfold Rminus in |- *; apply Rplus_le_compat_r; apply Rle_powerRZ;
auto with real arith zarith.
apply Rmult_le_compat; auto with real arith zarith.
apply Rmult_le_pos; auto with real.
apply Rmult_le_compat; auto with real arith zarith.
apply Rle_Rinv; auto with real arith zarith.
apply Rle_lt_trans with (1 - powerRZ radix 0)%R;
[ simpl in |- *; right; ring | auto with real arith zarith ].
unfold Rminus in |- *; apply Rplus_lt_compat_l; apply Ropp_lt_contravar;
auto with real arith zarith.
unfold Rminus in |- *; apply Rplus_le_compat_l; apply Ropp_le_contravar;
apply Rle_powerRZ; auto with real arith zarith.
unfold Rminus in |- *; apply Rplus_le_compat_l; apply Ropp_le_contravar;
auto with real arith zarith.
ring_simplify (1 - 0)%R; rewrite Rmult_1_r.
ring_simplify (1 + 0)%R; rewrite Rinv_1; rewrite Rmult_1_r.
rewrite powerRZ_Zopp; auto with real arith zarith.
replace (powerRZ radix 2%nat) with (INR 4);
[ idtac | simpl in |- *; ring; auto with arith real zarith ].
replace (4%nat - 1)%R with (INR 3);
[ idtac | simpl in |- *; ring; auto with arith real zarith ].
replace (1 - / 4%nat)%R with (3%nat × / 4%nat)%R.
rewrite Rinv_mult_distr; auto with real arith zarith.
rewrite Rinv_involutive; auto with real arith zarith.
apply Rle_lt_trans with (4%nat × / 4%nat × (/ 3%nat × / 3%nat))%R;
[ right; ring | idtac ].
rewrite Rinv_r; auto with real arith zarith; rewrite Rmult_1_l.
apply Rlt_le_trans with (/ 3%nat × 1)%R;
[ apply Rmult_lt_compat_l; auto with arith real zarith | right; ring ].
apply Rmult_lt_reg_l with (INR 3); auto with arith real.
rewrite Rinv_r; auto with arith real; rewrite Rmult_1_r; auto with arith real.
simpl; field; auto with real.
Qed.
End AxpyAux.
Section Axpy.
Let radix := 2%Z.
Let FtoRradix := FtoR radix.
Coercion FtoRradix : float >-> R.
Variable b : Fbound.
Variable precision : nat.
Hypothesis precisionGreaterThanOne : 1 < precision.
Hypothesis pGivesBound : Zpos (vNum b) = Zpower_nat radix precision.
Theorem Axpy_tFlessu :
∀ (a1 x1 y1 : R) (a x y t u : float),
Fbounded b a →
Fbounded b x →
Fbounded b y →
Fbounded b t →
Fbounded b u →
Closest b radix (a × x) t →
Closest b radix (t + y) u →
Fcanonic radix b u →
Fcanonic radix b t →
(4%nat × Rabs t ≤ Rabs u)%R →
(Rabs (y1 - y) + Rabs (a1 × x1 - a × x) <
/ 4%nat × Fulp b radix precision (FLess b precision u))%R →
MinOrMax radix b (a1 × x1 + y1) u.
intros a1 x1 y1 a x y t u Fa Fx Fy Ft Fu tDef uDef Cu Ct H1.
unfold FLess in |- *; case (Rcase_abs (FtoR 2 u)); fold radix in |- *;
intros H3 H2.
2: cut (0 ≤ u)%R;
[ intros H'; case H'; intros H4; clear H' | auto with real ].
2: apply AxpyPos with precision a x y t; auto.
apply MinOrMax_Fopp.
replace (- (a1 × x1 + y1))%R with (- a1 × x1 + - y1)%R; [ idtac | ring ].
apply AxpyPos with precision (Fopp a) x (Fopp y) (Fopp t);
fold radix FtoRradix in |- *; auto with real float.
replace (FtoRradix (Fopp a) × FtoRradix x)%R with (- (a × x))%R;
[ apply ClosestOpp; auto
| unfold FtoRradix in |- *; rewrite Fopp_correct; ring ].
replace (FtoRradix (Fopp t) + FtoRradix (Fopp y))%R with (- (t + y))%R;
[ apply ClosestOpp; auto
| unfold FtoRradix in |- *; repeat rewrite Fopp_correct; ring ].
unfold FtoRradix in |- *; rewrite Fopp_correct; auto with real.
unfold FtoRradix in |- *; repeat rewrite Fopp_correct;
repeat rewrite Rabs_Ropp; auto with real.
unfold FtoRradix in |- *; repeat rewrite Fopp_correct;
rewrite <- (Rabs_Ropp (- a1 × x1 - - FtoR radix a × FtoR radix x));
auto with real.
replace (- (- a1 × x1 - - FtoR radix a × FtoR radix x))%R with
(a1 × x1 - FtoRradix a × FtoRradix x)%R;
[ auto | unfold FtoRradix in |- *; ring ].
rewrite <- Rabs_Ropp.
replace (- (- y1 - - FtoR radix y))%R with (y1 - FtoRradix y)%R;
[ idtac | unfold FtoRradix in |- *; ring ].
rewrite FPredFopFSucc; auto with zarith; rewrite Fopp_Fopp; auto with zarith.
replace (Fulp b radix precision (Fopp (FSucc b radix precision u))) with
(Fulp b radix precision (FSucc b radix precision u));
auto.
unfold Fulp in |- *; rewrite Fnormalize_Fopp; simpl in |- *;
auto with real zarith.
apply MinOrMax3 with precision; auto with zarith; fold FtoRradix in |- ×.
cut (FtoRradix t = 0%R); [ intros H' | idtac ].
cut (FtoRradix y = 0%R); [ intros H5 | idtac ].
replace (a1 × x1 + y1 - FtoRradix u)%R with
(y1 - y + (a1 × x1 - a × x) + a × x)%R.
2: rewrite <- H4; rewrite H5; ring.
apply
Rle_lt_trans
with
(Rabs (y1 - FtoRradix y + (a1 × x1 - FtoRradix a × FtoRradix x)) +
Rabs (FtoRradix a × FtoRradix x))%R; [ apply Rabs_triang | idtac ].
apply
Rlt_le_trans
with
(/ 4%nat × Fulp b radix precision (FPred b radix precision u) +
Rabs (FtoRradix a × FtoRradix x))%R; auto with real.
apply Rplus_lt_compat_r; apply Rle_lt_trans with (2 := H2); apply Rabs_triang.
apply
Rle_trans
with
(/ 4%nat × Fulp b radix precision (FPred b radix precision u) +
/ 2%nat × Fulp b radix precision (FPred b radix precision u))%R;
[ apply Rplus_le_compat_l | idtac ].
replace (FtoRradix a × FtoRradix x)%R with (FtoRradix a × FtoRradix x - t)%R;
[ idtac | rewrite H'; ring ].
apply Rmult_le_reg_l with (INR 2); auto with real arith.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real arith.
ring_simplify (1 × Fulp b radix precision (FPred b radix precision u))%R;
apply Rle_trans with (Fulp b radix precision t).
unfold FtoRradix in |- *; apply ClosestUlp; auto with zarith.
unfold Fulp in |- *; apply Rle_powerRZ; auto with real zarith.
replace (Fexp (Fnormalize radix b precision t)) with (- dExp b)%Z.
cut (Fbounded b (Fnormalize radix b precision (FPred b radix precision u)));
[ intros H6; elim H6; auto | apply FnormalizeBounded; auto with zarith ].
apply FBoundedPred; auto with zarith.
replace (Fnormalize radix b precision t) with t.
replace t with (Float 0 (- dExp b)); [ simpl in |- *; auto | idtac ].
apply FcanonicUnique with (radix := radix) (precision := precision) (b := b);
auto with zarith.
right; repeat (split; simpl in |- *; auto with zarith).
fold FtoRradix in |- *; rewrite H'; unfold FtoRradix, FtoR in |- *;
simpl in |- *; ring.
apply FcanonicUnique with (radix := radix) (precision := precision) (b := b);
auto with zarith.
apply FnormalizeCanonic; auto with zarith.
apply sym_eq; apply FnormalizeCorrect; auto with zarith.
apply
Rle_trans
with
((/ 4%nat + / 2%nat) × Fulp b radix precision (FPred b radix precision u))%R;
[ right; ring | idtac ].
apply
Rle_trans with (1 × Fulp b radix precision (FPred b radix precision u))%R;
[ apply Rmult_le_compat_r | right; ring ].
unfold Fulp in |- *; auto with real zarith.
apply Rmult_le_reg_l with (INR 4); auto with real arith.
apply Rle_trans with 3%R;simpl;[right; field|idtac].
apply Rle_trans with (3+1)%R; auto with real; right; ring.
cut (ProjectorP b radix (Closest b radix));
[ unfold ProjectorP in |- *; intros V | idtac ].
replace 0%R with (FtoRradix (Float 0 (- dExp b)));
[ idtac | unfold FtoRradix, FtoR in |- *; simpl in |- *; ring ].
unfold FtoRradix in |- *; apply V; auto.
replace (Float 0 (- dExp b)) with u; fold FtoRradix in |- ×.
replace (FtoRradix y) with (FtoRradix t + FtoRradix y)%R; auto.
rewrite H'; ring.
apply FcanonicUnique with (radix := radix) (precision := precision) (b := b);
auto with zarith.
right; repeat (split; simpl in |- *; auto with zarith).
fold FtoRradix in |- *; rewrite <- H4; unfold FtoRradix, FtoR in |- *;
simpl in |- *; ring.
apply RoundedProjector; auto.
apply ClosestRoundedModeP with precision; auto with zarith.
cut (∀ z : R, (Rabs z ≤ 0)%R → z = 0%R); [ intros V; apply V | idtac ].
apply Rmult_le_reg_l with (INR 4); auto with real arith.
apply Rle_trans with (1 := H1); right; rewrite <- H4; auto with real.
ring_simplify (4%nat×0)%R; apply Rabs_R0.
intros z V; case (Req_dec 0 z); auto with real.
intros V'; Contradict V.
apply Rlt_not_le; apply Rabs_pos_lt; auto.
Qed.
Theorem Axpy_opt :
∀ (a1 x1 y1 : R) (a x y t u : float),
(Fbounded b a) → (Fbounded b x) → (Fbounded b y) →
(Fbounded b t) → (Fbounded b u) →
(Closest b radix (a × x) t) →
(Closest b radix (t + y) u) →
(Fcanonic radix b u) → (Fcanonic radix b t) →
((5%nat + 4%nat × (powerRZ radix (- precision))) ×
/ (1 - powerRZ radix (- precision)) ×
(Rabs (a × x) + (powerRZ radix (Zpred (- dExp b))))
≤ Rabs y)%R →
(Rabs (y1 - y) + Rabs (a1 × x1 - a × x) ≤
(powerRZ radix (Zpred (Zpred (- precision)))) ×
(1 - powerRZ radix (Zsucc (- precision))) × Rabs y +
- (powerRZ radix (Zpred (Zpred (- precision))) × Rabs (a × x)) +
- powerRZ radix (Zpred (Zpred (- dExp b))))%R →
(MinOrMax radix b (a1 × x1 + y1) u).
intros a1 x1 y1 a x y t u Fa Fx Fy Ft Fu tDef uDef Cu Ct H1 H2.
apply Axpy_tFlessu with a x y t; auto.
cut (0 < 1 - powerRZ radix (- precision))%R; [ intros H'1 | idtac ].
2: apply Rplus_lt_reg_r with (powerRZ radix (- precision)).
2: ring_simplify (powerRZ radix (- precision) + 0)%R.
2: ring_simplify (powerRZ radix (- precision) + (1 - powerRZ radix (- precision)))%R.
2: replace 1%R with (powerRZ radix 0);
[ auto with real zarith | simpl in |- *; auto ].
cut (0 < 1 + powerRZ radix (- precision))%R; [ intros H'2 | idtac ].
2: apply Rlt_le_trans with 1%R; auto with real.
2: apply Rle_trans with (1 + 0)%R; auto with real zarith.
cut ((5%nat + 4%nat × powerRZ radix (- precision)) × Rabs t ≤ Rabs y)%R;
[ intros H3 | idtac ].
apply Rplus_le_reg_l with (- Rabs (FtoRradix u))%R.
ring_simplify (- Rabs (FtoRradix u) + Rabs (FtoRradix u))%R.
apply
Rle_trans
with
(- ((Rabs y + - Rabs t) × / (1 + powerRZ radix (- precision))) +
4%nat × Rabs (FtoRradix t))%R;
[ apply Rplus_le_compat_r; apply Ropp_le_contravar | idtac ].
apply
Rle_trans
with
(Rabs (FtoRradix y + FtoRradix t) × / (1 + powerRZ radix (- precision)))%R.
apply Rmult_le_compat_r; auto with real zarith.
rewrite <- (Rabs_Ropp (FtoRradix t)).
apply Rle_trans with (Rabs (FtoRradix y) - Rabs (- FtoRradix t))%R;
[ right; ring | idtac ].
replace (FtoRradix y + FtoRradix t)%R with (FtoRradix y - - FtoRradix t)%R;
[ apply Rabs_triang_inv | ring ].
case Cu; intros H4.
apply Rmult_le_reg_l with (1 + powerRZ radix (- precision))%R; auto.
apply
Rle_trans
with
(Rabs (FtoRradix y + FtoRradix t) ×
((1 + powerRZ radix (- precision)) × / (1 + powerRZ radix (- precision))))%R;
[ right; ring; ring | idtac ].
rewrite Rinv_r; auto with real.
ring_simplify (Rabs (FtoRradix y + FtoRradix t) × 1)%R.
rewrite Rmult_plus_distr_r.
apply Rle_trans with (Rabs u + / radix × Fulp b radix precision u)%R.
apply Rplus_le_reg_l with (- Rabs u)%R.
apply Rle_trans with (/ radix × Fulp b radix precision u)%R;
[idtac|right; ring].
apply Rle_trans with (Rabs (FtoRradix y + FtoRradix t - FtoRradix u)).
apply
Rle_trans with (Rabs (FtoRradix y + FtoRradix t) - Rabs (FtoRradix u))%R;
[ right; ring | apply Rabs_triang_inv ].
apply Rmult_le_reg_l with (IZR radix); auto with real zarith.
apply Rle_trans with (Fulp b radix precision u × (radix × / radix))%R;
[ rewrite Rinv_r; auto with real zarith | right; ring ].
ring_simplify (Fulp b radix precision u × 1)%R.
unfold FtoRradix in |- *; apply ClosestUlp; auto with real zarith.
rewrite Rplus_comm; auto.
ring_simplify (1×Rabs u)%R; apply Rplus_le_compat_l.
apply Rmult_le_reg_l with (IZR radix); auto with real zarith.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real zarith.
ring_simplify (1 × Fulp b radix precision u)%R.
apply
Rle_trans with (Rabs (FtoRradix u) × powerRZ radix (Zsucc (- precision)))%R.
unfold FtoRradix in |- *; apply FulpLe2; auto with zarith.
replace (Fnormalize radix b precision u) with u;
[ idtac
| apply
FcanonicUnique with (radix := radix) (precision := precision) (b := b) ];
auto with zarith real.
apply FnormalizeCanonic; auto with zarith.
apply sym_eq; apply FnormalizeCorrect; auto with zarith.
unfold Zsucc in |- *; rewrite powerRZ_add; auto with zarith real;
simpl in |- *; auto with zarith real; right; ring.
replace (FtoRradix y + FtoRradix t)%R with (FtoRradix u); auto with real.
apply Rle_trans with (Rabs (FtoRradix u) × 1)%R;
[ apply Rmult_le_compat_l | right; ring ]; auto with real.
pattern 1%R at 2 in |- *; replace 1%R with (/ 1)%R; auto with real zarith.
apply Rle_Rinv; auto with real zarith.
pattern 1%R at 1 in |- *; replace 1%R with (1 + 0)%R; auto with real zarith.
unfold FtoRradix in |- *; apply plusExact1 with b precision;
auto with real zarith.
rewrite Rplus_comm; auto.
elim H4; intros H5 H6; elim H6; intros H7 H8; rewrite H7; apply Zmin_Zle;
elim Fy; elim Ft; auto with zarith.
apply
Rle_trans
with
(((5%nat + 4%nat × powerRZ radix (- precision)) × Rabs (FtoRradix t) +
- Rabs (FtoRradix y)) × / (1 + powerRZ radix (- precision)))%R.
right; replace (INR 5) with (1 + 4%nat)%R;
[ idtac
| replace 1%R with (INR 1); [ rewrite <- plus_INR | idtac ];
auto with real arith ].
replace (4%nat × Rabs (FtoRradix t))%R with
(4%nat × Rabs (FtoRradix t) ×
((1 + powerRZ radix (- precision)) × / (1 + powerRZ radix (- precision))))%R.
ring; ring.
rewrite Rinv_r; auto with real; ring.
apply Rle_trans with (0 × / (1 + powerRZ radix (- precision)))%R;
[ apply Rmult_le_compat_r; auto with real | right; ring ].
apply Rplus_le_reg_l with (Rabs (FtoRradix y)).
apply
Rle_trans
with ((5%nat + 4%nat × powerRZ radix (- precision)) × Rabs (FtoRradix t))%R;
[ right; ring | ring_simplify (Rabs (FtoRradix y) + 0)%R ];
auto.
apply Rle_trans with (2 := H1).
rewrite Rmult_assoc.
apply Rmult_le_compat_l; auto with real arith.
apply Rle_trans with (5%nat + 4%nat × 0)%R; auto with real arith.
ring_simplify (5%nat + 4%nat × 0)%R; auto with real arith.
apply
Rle_trans
with
(Rabs (FtoRradix a × FtoRradix x) × / (1 - powerRZ radix (- precision)) +
powerRZ radix (Zpred (- dExp b)) × / (1 - powerRZ 2%Z (- precision)))%R;
[ idtac | right; simpl; ring ].
unfold FtoRradix in |- *; apply RoundLeGeneral; auto with real zarith.
apply Rle_lt_trans with (1 := H2).
apply UlpFlessuGe2 with t; auto.
Qed.
Theorem Axpy_Simpl1 :
∀ (a1 x1 y1 : R) (a x y t u : float),
4 ≤ precision →
Fbounded b a →
Fbounded b x →
Fbounded b y →
Fbounded b t →
Fbounded b u →
Closest b radix (a × x) t →
Closest b radix (t + y) u →
Fcanonic radix b u →
Fcanonic radix b t →
(6%nat × (Rabs (a × x) + powerRZ radix (Zpred (- dExp b))) ≤ Rabs y)%R →
(Rabs (y1 - y) + Rabs (a1 × x1 - a × x) ≤
powerRZ radix (Zpred (Zpred (- precision))) ×
(1 - powerRZ radix (Zsucc (- precision))) × Rabs y +
- (powerRZ radix (Zpred (Zpred (- precision))) × Rabs (a × x)) +
- powerRZ radix (Zpred (Zpred (- dExp b))))%R →
MinOrMax radix b (a1 × x1 + y1) u.
intros a1 x1 y1 a x y t u Enoughp Fa Fx Fy Ft Fu tDef uDef Cu Ct H1 H2.
apply Axpy_opt with a x y t; auto.
cut (0 < 1 - powerRZ radix (- precision))%R; [ intros H'1 | idtac ].
2: apply Rplus_lt_reg_r with (powerRZ radix (- precision)).
2: ring_simplify (powerRZ radix (- precision) + 0)%R.
2: ring_simplify (powerRZ radix (- precision) + (1 - powerRZ radix (- precision)))%R.
2: replace 1%R with (powerRZ radix 0);
[ auto with real zarith | simpl in |- *; auto ].
apply Rle_trans with (2 := H1).
apply Rmult_le_compat_r.
apply Rle_trans with (Rabs (FtoRradix a × FtoRradix x) + 0)%R;
auto with real zarith.
ring_simplify (Rabs (FtoRradix a × FtoRradix x) + 0)%R; auto with real.
apply Rmult_le_reg_l with (1 - powerRZ radix (- precision))%R; auto.
apply
Rle_trans
with
((5%nat + 4%nat × powerRZ radix (- precision)) ×
((1 - powerRZ radix (- precision)) × / (1 - powerRZ radix (- precision))))%R;
[ right; ring; ring | rewrite Rinv_r; auto with real ].
apply Rplus_le_reg_l with (- 5%nat + 6%nat × powerRZ radix (- precision))%R.
simpl; ring_simplify.
apply Rle_trans with (powerRZ radix 4%nat × powerRZ radix (- precision))%R;
[ apply Rmult_le_compat_r; auto with real arith zarith
| rewrite <- powerRZ_add; auto with real zarith ].
apply Rle_trans with (INR 10); [ right; simpl in |- *; ring | idtac ].
apply Rle_trans with (INR 16); [ idtac | right; simpl in |- *; ring ].
apply Rle_INR; auto with arith.
apply Rle_trans with (powerRZ radix 0);
[ idtac | right; simpl in |- *; ring ].
apply Rle_powerRZ; auto with zarith arith real.
apply Zle_trans with (4-precision)%Z; auto with zarith.
Qed.
Theorem Axpy_Simpl1bis :
∀ (a1 x1 y1 : R) (a x y t u : float),
4 ≤ precision →
Fbounded b a →
Fbounded b x →
Fbounded b y →
Fbounded b t →
Fbounded b u →
Closest b radix (a × x) t →
Closest b radix (t + y) u →
Fcanonic radix b u →
Fcanonic radix b t →
(6%nat × (Rabs (a × x) + powerRZ radix (Zpred (- dExp b))) ≤ Rabs y)%R →
(Rabs (y1 - y) + Rabs (a1 × x1 - a × x) ≤
powerRZ radix (Zpred (Zpred (- precision))) ×
(5%nat × / 6%nat - powerRZ radix (Zsucc (- precision))) ×
Rabs y + - powerRZ radix (Zpred (Zpred (- dExp b))))%R →
MinOrMax radix b (a1 × x1 + y1) u.
intros a1 x1 y1 a x y t u Enoughp Fa Fx Fy Ft Fu tDef uDef Cu Ct H1 H2.
apply Axpy_Simpl1 with a x y t; auto.
apply Rle_trans with (1 := H2).
replace (5%nat × / 6%nat)%R with (1 - / 6%nat)%R.
2: apply Rmult_eq_reg_l with (INR 6); auto with arith real.
2: apply trans_eq with (6%nat - 6%nat × / 6%nat)%R;
[ ring; ring | rewrite Rinv_r; auto with arith real ].
2: apply trans_eq with (5%nat × (6%nat × / 6%nat))%R;
[ rewrite Rinv_r; auto with arith real | ring ];
simpl in |- *; ring.
apply
Rplus_le_reg_l
with
(-
(powerRZ radix (Zpred (Zpred (- precision))) ×
(1 - powerRZ radix (Zsucc (- precision))) × Rabs (FtoRradix y)))%R.
ring_simplify.
unfold Rminus;apply Rplus_le_compat_r.
repeat rewrite Ropp_mult_distr_l_reverse; apply Ropp_le_contravar.
rewrite Rmult_assoc; apply Rmult_le_compat_l; auto with real zarith.
apply Rmult_le_reg_l with (INR 6); auto with arith real.
apply Rle_trans with (Rabs (FtoRradix y) × (6%nat × / 6%nat))%R;
[ rewrite Rinv_r; auto with arith real | right; ring ].
ring_simplify (Rabs (FtoRradix y) × 1)%R; apply Rle_trans with (2 := H1).
apply Rle_trans with (6%nat × (Rabs (FtoRradix a × FtoRradix x) + 0))%R;
auto with real zarith.
Qed.
Theorem Axpy_Simpl2 :
∀ (a1 x1 y1 : R) (a x y t u : float),
4 ≤ precision →
Fbounded b a →
Fbounded b x →
Fbounded b y →
Fbounded b t →
Fbounded b u →
Closest b radix (a × x) t →
Closest b radix (t + y) u →
Fcanonic radix b u →
Fcanonic radix b t →
(6%nat × (Rabs (a × x) + powerRZ radix (Zpred (- dExp b))) ≤ Rabs y)%R →
(Rabs (y1 - y) + Rabs (a1 × x1 - a × x) ≤
powerRZ radix (Zpred (Zpred (- precision))) × (2%nat × / 3%nat) × Rabs y +
- powerRZ radix (Zpred (Zpred (- dExp b))))%R →
MinOrMax radix b (a1 × x1 + y1) u.
intros a1 x1 y1 a x y t u Enoughp Fa Fx Fy Ft Fu tDef uDef Cu Ct H1 H2.
apply Axpy_Simpl1bis with a x y t; auto.
apply Rle_trans with (1 := H2).
apply Rplus_le_compat_r; apply Rmult_le_compat_r; auto with real.
apply Rmult_le_compat_l; auto with real zarith.
apply Rle_trans with (5%nat × / 6%nat - / 8%nat)%R.
apply Rmult_le_reg_l with (6%nat × 8%nat)%R;
[ apply Rmult_lt_0_compat; auto with real arith | idtac ].
apply
Rle_trans
with (6%nat × / 6%nat × (8%nat × 5%nat) - 8%nat × / 8%nat × 6%nat)%R;
[ idtac | right; ring ].
pattern (INR 6) at 1 in |- *; replace (INR 6) with (2%nat × 3%nat)%R;
[ idtac | rewrite <- mult_INR; auto with real arith ].
apply Rle_trans with (3%nat × / 3%nat × 8%nat × (2%nat × 2%nat))%R;
[ right; ring | idtac ].
repeat rewrite Rinv_r; auto with real arith.
apply Rle_trans with (INR 32);[simpl; right; ring|idtac].
apply Rle_trans with (INR 34);[auto with zarith real|simpl; right; ring].
unfold Rminus in |- *; apply Rplus_le_compat_l; apply Ropp_le_contravar.
apply Rle_trans with (powerRZ radix (Zsucc (- 4%nat)));
[ apply Rle_powerRZ; auto with real zarith | idtac ].
simpl; right; field; auto with real.
Qed.
End Axpy.
Section AxpyFmac.
Let radix := 2%Z.
Let FtoRradix := FtoR radix.
Coercion FtoRradix : float >-> R.
Variable b : Fbound.
Variable precision : nat.
Hypothesis precisionGreaterThanOne : 1 < precision.
Hypothesis pGivesBound : Zpos (vNum b) = Zpower_nat radix precision.
Theorem AxpyPos_Fmac :
∀ (a1 x1 y1 : R) (a x y u : float),
Fbounded b a →
Fbounded b x →
Fbounded b y →
Fbounded b u →
Closest b radix (a × x + y) u →
Fcanonic radix b u →
(0 < u)%R →
(Rabs (y1 - y) + Rabs (a1 × x1 - a × x) <
/ 2%nat × Fulp b radix precision (FPred b radix precision u))%R →
MinOrMax radix b (a1 × x1 + y1) u.
intros a1 x1 y1 a x y u Fa Fx Fy Fu uDef Cu Pu H.
cut
(Rabs (a1 × x1 + y1 - u) <
/ 2%nat × Fulp b radix precision (FPred b radix precision u) +
Rabs (a × x + y - u))%R; [ intros H4 | idtac ].
case (Rle_or_lt (a × x + y) u); intros H5.
apply MinOrMax1 with precision; auto with zarith arith.
apply Rlt_le_trans with (1 := H4).
replace 2%Z with radix; auto.
apply
Rplus_le_reg_l
with (- (/ 2%nat × Fulp b radix precision (FPred b radix precision u)))%R.
ring_simplify.
apply
Rle_trans
with (/ 2%nat × Fulp b radix precision (FPred b radix precision u))%R;
[ idtac | right ].
generalize uDef; unfold Closest in |- *; intros H6; elim H6; intros H7 H8;
clear H6 H7.
case
(Rle_or_lt (Rabs (FtoRradix a × FtoRradix x + FtoRradix y - FtoRradix u))
(/ 2%nat × Fulp b radix precision (FPred b radix precision u)));
auto; intros H6.
cut
(Rabs (FtoR radix u - (FtoRradix a × FtoRradix x + FtoRradix y)) ≤
Rabs
(FtoR radix (FPred b radix precision u) -
(FtoRradix a × FtoRradix x + FtoRradix y)))%R;
[ intros H7 | apply H8 ].
2: apply FBoundedPred; auto with zarith arith.
Contradict H7; apply Rlt_not_le.
apply
Rlt_le_trans
with (Rabs (FtoRradix a × FtoRradix x + FtoRradix y - FtoRradix u)).
2: right; rewrite <- Rabs_Ropp.
2: replace (- (FtoRradix a × FtoRradix x + FtoRradix y - FtoRradix u))%R with
(FtoR radix u - (FtoRradix a × FtoRradix x + FtoRradix y))%R;
auto with real; ring.
apply Rle_lt_trans with (2 := H6).
rewrite Faux.Rabsolu_left1.
apply
Rplus_le_reg_l
with
(u - (FtoRradix a × FtoRradix x + FtoRradix y) -
/ 2%nat × Fulp b radix precision (FPred b radix precision u))%R.
ring_simplify.
pattern (FtoRradix u) at 1 in |- *;
replace (FtoRradix u) with
(FtoRradix (FPred b radix precision u) +
Fulp b radix precision (FPred b radix precision u))%R;
[ idtac
| unfold FtoRradix in |- *; apply FpredUlpPos; auto with zarith arith ].
apply
Rle_trans
with
(Fulp b radix precision (FPred b radix precision u) -
Fulp b radix precision (FPred b radix precision u) × / 2%nat)%R;
[ right; fold FtoRradix; ring | idtac ].
apply
Rle_trans
with (/ 2%nat × Fulp b radix precision (FPred b radix precision u))%R.
right; simpl; field.
apply Rlt_le; apply Rlt_le_trans with (1 := H6).
rewrite Faux.Rabsolu_left1; [ right; ring | idtac ].
apply Rplus_le_reg_l with (FtoRradix u).
ring_simplify; auto with real.
fold FtoRradix in |- *;
apply
Rplus_le_reg_l with (Fulp b radix precision (FPred b radix precision u)).
apply
Rle_trans
with
(FtoRradix (FPred b radix precision u) +
Fulp b radix precision (FPred b radix precision u) +
- (FtoRradix a × FtoRradix x + FtoRradix y))%R;
[ right; ring
| unfold FtoRradix in |- *; rewrite FpredUlpPos; auto with zarith arith ].
ring_simplify (Fulp b radix precision (FPred b radix precision u) + 0)%R.
apply
Rle_trans
with (Rabs (FtoRradix u + - (FtoRradix a × FtoRradix x + FtoRradix y))).
apply RRle_abs.
rewrite <- Rabs_Ropp.
replace (- (FtoRradix u + - (FtoRradix a × FtoRradix x + FtoRradix y)))%R
with (FtoRradix a × FtoRradix x + FtoRradix y - FtoRradix u)%R;
[ idtac | ring ].
apply Rmult_le_reg_l with (INR 2); auto with arith real.
apply Rle_trans with (Fulp b radix precision u).
unfold FtoRradix in |- *; apply ClosestUlp; auto with arith zarith.
replace (INR 2) with (IZR radix); auto with zarith arith real.
apply FulpFPredLe; auto with zarith.
simpl; field.
case
(Rle_or_lt (a × x + y)
(u + / 2%nat × Fulp b radix precision (FPred b radix precision u)));
intros H6.
apply MinOrMax1 with precision; auto with arith zarith real float.
apply Rlt_le_trans with (1 := H4); rewrite Rabs_right.
apply
Rplus_le_reg_l
with
(u + - (/ 2%nat × Fulp b radix precision (FPred b radix precision u)))%R.
ring_simplify.
apply Rle_trans with (1 := H6).
right; simpl; field.
apply Rle_ge; apply Rplus_le_reg_l with (FtoRradix u); auto with real.
apply MinOrMax2 with precision; auto with arith zarith float real.
apply Rlt_le_trans with (1 := H4).
apply
Rle_trans
with
(/ 2%nat × Fulp b radix precision u +
Rabs (FtoRradix a × FtoRradix x + FtoRradix y - FtoRradix u))%R;
[ apply Rplus_le_compat_r | idtac ].
apply Rmult_le_compat_l; auto with real arith.
apply FulpFPredGePos; auto with arith zarith.
apply
Rle_trans
with
(/ 2%nat × Fulp b radix precision u + / 2%nat × Fulp b radix precision u)%R;
[ apply Rplus_le_compat_l | right ].
apply Rmult_le_reg_l with (INR 2); auto with real arith.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real arith.
ring_simplify (1 × Fulp b radix precision u)%R.
unfold FtoRradix in |- *; apply ClosestUlp; auto with real arith zarith float.
apply trans_eq with ((/ 2%nat + / 2%nat) × Fulp b radix precision u)%R;
[ ring | auto with real arith ].
replace (/ 2%nat + / 2%nat)%R with 1%R; [ ring | auto with real arith ].
rewrite <- (Rinv_r (INR 2)); auto with real arith; simpl in |- *; ring.
replace 2%Z with radix; auto; fold FtoRradix in |- *; auto with real.
apply Rlt_le;
apply
Rplus_lt_reg_r
with (/ 2%nat × Fulp b radix precision (FPred b radix precision u))%R.
rewrite Rplus_comm; apply Rlt_le_trans with (1 := H6).
apply Rplus_le_reg_l with (- (a1 × x1 + y1))%R.
ring_simplify
(- (a1 × x1 + y1) +
(/ 2%nat × Fulp b radix precision (FPred b radix precision u) +
(a1 × x1 + y1)))%R.
apply
Rle_trans
with (Rabs (- (a1 × x1 + y1) + (FtoRradix a × FtoRradix x + FtoRradix y)));
[ apply RRle_abs | idtac ].
rewrite <- Rabs_Ropp.
replace (- (- (a1 × x1 + y1) + (FtoRradix a × FtoRradix x + FtoRradix y)))%R
with (y1 - FtoRradix y + (a1 × x1 - FtoRradix a × FtoRradix x))%R;
[ auto with real | ring ].
apply
Rle_trans
with
(Rabs (y1 - FtoRradix y) + Rabs (a1 × x1 - FtoRradix a × FtoRradix x))%R;
[ apply Rabs_triang | auto with real ].
replace (a1 × x1 + y1 - FtoRradix u)%R with
(y1 - y + (a1 × x1 - a × x) + (a × x + y - u))%R;
[ idtac | ring ].
apply
Rle_lt_trans
with
(Rabs (y1 - FtoRradix y + (a1 × x1 - FtoRradix a × FtoRradix x)) +
Rabs (FtoRradix a × FtoRradix x + FtoRradix y - FtoRradix u))%R;
[ apply Rabs_triang | idtac ].
apply Rplus_lt_compat_r; auto.
apply Rle_lt_trans with (2 := H); apply Rabs_triang.
Qed.
Theorem AxpyFLessu_Fmac :
∀ (a1 x1 y1 : R) (a x y u : float),
Fbounded b a →
Fbounded b x →
Fbounded b y →
Fbounded b u →
Closest b radix (a × x + y) u →
Fcanonic radix b u →
(Rabs (y1 - y) + Rabs (a1 × x1 - a × x) <
/ 2%nat × Fulp b radix precision (FLess b precision u))%R →
MinOrMax radix b (a1 × x1 + y1) u.
intros a1 x1 y1 a x y u Fa Fx Fy Fu uDef Cu.
unfold FLess in |- *; case (Rcase_abs (FtoR 2 u)); fold radix in |- *;
intros H3 H2.
2: cut (0 ≤ u)%R;
[ intros H'; case H'; intros H4; clear H' | auto with real ].
2: apply AxpyPos_Fmac with a x y; auto.
apply MinOrMax_Fopp.
replace (- (a1 × x1 + y1))%R with (- a1 × x1 + - y1)%R; [ idtac | ring ].
apply AxpyPos_Fmac with (Fopp a) x (Fopp y); fold radix FtoRradix in |- *;
auto with real float.
replace (FtoRradix (Fopp a) × FtoRradix x + FtoRradix (Fopp y))%R with
(- (a × x + y))%R;
[ apply ClosestOpp; auto
| unfold FtoRradix in |- *; repeat rewrite Fopp_correct;ring ].
unfold FtoRradix in |- *; rewrite Fopp_correct; auto with real.
unfold FtoRradix in |- *; repeat rewrite Fopp_correct; auto with real.
rewrite <- Rabs_Ropp;
rewrite <- (Rabs_Ropp (- a1 × x1 - - FtoR radix a × FtoR radix x)).
replace (- (- y1 - - FtoR radix y))%R with (y1 - FtoRradix y)%R;
[ idtac | fold FtoRradix; ring ].
replace (- (- a1 × x1 - - FtoR radix a × FtoR radix x))%R with
(a1 × x1 - FtoRradix a × FtoRradix x)%R;
[ auto | unfold FtoRradix in |- *; ring ].
rewrite FPredFopFSucc; auto with zarith; rewrite Fopp_Fopp; auto with zarith.
replace (Fulp b radix precision (Fopp (FSucc b radix precision u))) with
(Fulp b radix precision (FSucc b radix precision u));
auto.
unfold Fulp in |- *; rewrite Fnormalize_Fopp; simpl in |- *;
auto with real zarith.
apply MinOrMax3 with precision; auto with zarith.
cut (u = Float 0 (- dExp b)); [ intros H1 | idtac ].
2: apply
FcanonicUnique with (radix := radix) (precision := precision) (b := b);
auto with zarith.
2: right; repeat (split; simpl in |- *; auto with zarith).
2: fold FtoRradix in |- *; rewrite <- H4; unfold FtoRradix, FtoR in |- *;
simpl in |- *; ring.
replace 2%Z with radix; auto; fold FtoRradix in |- ×.
replace (a1 × x1 + y1 - FtoRradix u)%R with
(y1 - y + (a1 × x1 - a × x) + (FtoRradix a × FtoRradix x + FtoRradix y - u))%R;
[ idtac | ring ].
apply
Rle_lt_trans
with
(Rabs (y1 - FtoRradix y + (a1 × x1 - FtoRradix a × FtoRradix x)) +
Rabs (FtoRradix a × FtoRradix x + FtoRradix y - FtoRradix u))%R;
[ apply Rabs_triang | idtac ].
apply
Rlt_le_trans
with
(/ 2%nat × Fulp b radix precision (FPred b radix precision u) +
Rabs (FtoRradix a × FtoRradix x + FtoRradix y - FtoRradix u))%R.
apply Rplus_lt_compat_r; auto.
apply Rle_lt_trans with (2 := H2); apply Rabs_triang.
apply
Rle_trans
with
(/ 2%nat × Fulp b radix precision (FPred b radix precision u) +
/ 2%nat × Fulp b radix precision u)%R.
apply Rplus_le_compat_l.
apply Rmult_le_reg_l with (INR 2); auto with real arith.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real arith.
ring_simplify (1 × Fulp b radix precision u)%R.
unfold FtoRradix in |- *; apply ClosestUlp; auto with zarith.
unfold Fulp in |- ×.
repeat rewrite FcanonicFnormalizeEq; auto with zarith float.
rewrite H1; rewrite FPredSimpl4; auto with zarith; simpl in |- ×.
2: cut (0 < pPred (vNum b))%Z;
[ idtac | apply pPredMoreThanOne with radix precision ];
auto with zarith.
2: cut (0 < nNormMin radix precision)%Z; [ idtac | apply nNormPos ];
auto with zarith.
right; apply trans_eq with ((/ 2 + / 2) × powerRZ 2 (- dExp b))%R;
[ ring | idtac ].
replace (/ 2 + / 2)%R with 1%R; [ ring | idtac ].
field; auto with real.
Qed.
Theorem Axpy_opt_Fmac :
∀ (a1 x1 y1 : R) (a x y u : float),
(Fbounded b a) → (Fbounded b x) → (Fbounded b y) → (Fbounded b u) →
(Closest b radix (a × x + y) u) →
(Fcanonic radix b u) →
(Rabs (y1 - y) + Rabs (a1 × x1 - a × x) <
(Rabs (a × x + y)) × (powerRZ radix (Zpred (- precision)) ×
(1 - powerRZ radix (Zsucc (- precision)))) -
powerRZ radix (Zpred (Zpred (- dExp b))) ×
(3 × / (powerRZ radix precision - powerRZ radix (- precision))))%R →
(MinOrMax radix b (a1 × x1 + y1) u).
intros a1 x1 y1 a x y u Fa Fx Fy Fu uDef Cu H.
apply AxpyFLessu_Fmac with a x y; auto.
cut (0 < powerRZ radix precision - 1)%R; [ intros H'1 | idtac ].
2: apply Rplus_lt_reg_r with 1%R.
2: ring_simplify.
2: replace 1%R with (powerRZ radix 0); [ idtac | simpl in |- × ];
auto with real zarith.
cut (0 < powerRZ radix precision - powerRZ radix (- precision))%R;
[ intros H'2 | idtac ].
2: apply Rplus_lt_reg_r with (powerRZ radix (- precision)).
2: ring_simplify; auto with real zarith.
apply Rlt_le_trans with (1 := H).
apply
Rle_trans
with
(/ 2%nat ×
(/ (powerRZ radix precision - 1) ×
(Rabs u × (1 - powerRZ radix (Zsucc (- precision))) -
powerRZ radix (- dExp b))))%R.
apply
Rle_trans
with
(/ 2%nat ×
(/ (powerRZ radix precision - 1) ×
((Rabs (a × x + y) × / (1 + powerRZ radix (- precision)) -
powerRZ radix (Zpred (- dExp b)) ×
/ (1 + powerRZ radix (- precision))) ×
(1 - powerRZ radix (Zsucc (- precision))) -
powerRZ radix (- dExp b))))%R.
apply
Rle_trans
with
(/ 2%nat ×
(/ (powerRZ radix precision - 1) ×
(Rabs (FtoRradix a × FtoRradix x + FtoRradix y) ×
/ (1 + powerRZ radix (- precision)) ×
(1 - powerRZ radix (Zsucc (- precision))))) -
(/ 2%nat ×
(/ (powerRZ radix precision - 1) ×
((1 - powerRZ radix (Zsucc (- precision))) ×
(powerRZ radix (Zpred (- dExp b)) ×
/ (1 + powerRZ radix (- precision))))) +
/ 2%nat × (/ (powerRZ radix precision - 1) × powerRZ radix (- dExp b))))%R;
[ idtac | right; ring; ring ].
apply
Rle_trans
with
(/ 2%nat ×
(/ (powerRZ radix precision - 1) ×
(Rabs (FtoRradix a × FtoRradix x + FtoRradix y) ×
/ (1 + powerRZ radix (- precision)) ×
(1 - powerRZ radix (Zsucc (- precision))))) -
powerRZ radix (Zpred (Zpred (- dExp b))) ×
(3 × / (powerRZ radix precision - powerRZ radix (- precision))))%R;
unfold Rminus in |- *; [ apply Rplus_le_compat_r | apply Rplus_le_compat_l ].
apply
Rle_trans
with
(Rabs (FtoRradix a × FtoRradix x + FtoRradix y) ×
(/ 2%nat ×
(/ (powerRZ radix precision + -1) ×
(/ (1 + powerRZ radix (- precision)) ×
(1 + - powerRZ radix (Zsucc (- precision)))))))%R;
[ idtac | right; ring ].
apply Rmult_le_compat_l; auto with real.
apply
Rle_trans
with
(/ 2%nat ×
(/ (powerRZ radix precision + -1) × / (1 + powerRZ radix (- precision))) ×
(1 + - powerRZ radix (Zsucc (- precision))))%R;
[ idtac | right; ring ].
apply Rmult_le_compat_r; auto with real zarith.
apply Rplus_le_reg_l with (powerRZ radix (Zsucc (- precision))).
ring_simplify.
replace 1%R with (powerRZ radix 0); auto with real zarith.
replace (Zpred (- precision)) with (- Zsucc precision)%Z;
[ idtac | unfold Zsucc, Zpred in |- *; ring ].
rewrite <- Rinv_powerRZ; auto with real zarith.
rewrite <- Rinv_mult_distr; auto with real zarith.
rewrite <- Rinv_mult_distr; auto with real zarith.
apply Rle_Rinv; auto with real zarith.
apply Rmult_lt_0_compat; auto with real arith zarith.
apply Rmult_lt_0_compat; auto with real arith zarith.
apply Rlt_le_trans with (1 + 0)%R; auto with real arith zarith.
ring_simplify (1 + 0)%R; auto with real.
replace (INR 2) with (powerRZ radix 1); auto with real zarith.
ring_simplify
(powerRZ radix 1 ×
((powerRZ radix precision + -1) × (1 + powerRZ radix (- precision))))%R.
repeat rewrite <- powerRZ_add; auto with zarith real.
ring_simplify (1+ precision + -precision)%Z.
apply
Rle_trans
with (powerRZ radix (precision + 1) + - powerRZ radix (- precision + 1))%R;
[ right | idtac ]; auto with real zarith.
unfold Zminus; rewrite Zplus_comm; rewrite Zplus_comm
with (n:=(-precision)%Z); ring.
apply Rle_trans with (powerRZ radix (precision + 1) + -0)%R;
auto with real zarith.
right;ring.
apply Rmult_integral_contrapositive; split; auto with real zarith.
apply not_eq_sym; apply Rlt_not_eq; auto with real zarith.
apply Rlt_le_trans with (1 + 0)%R;
[ ring_simplify (1 + 0)%R; auto with real | auto with real zarith ].
apply not_eq_sym; apply Rlt_not_eq; auto with real zarith.
apply Rlt_le_trans with (1 + 0)%R;
[ ring_simplify (1 + 0)%R; auto with real | auto with real zarith ].
apply Ropp_le_contravar.
replace (powerRZ radix precision + -1)%R with (powerRZ radix precision - 1)%R;
[ idtac | ring ].
replace (powerRZ radix precision + - powerRZ radix (- precision))%R with
(powerRZ radix precision - powerRZ radix (- precision))%R;
[ idtac | ring ].
cut (0 < 1 + powerRZ radix (- precision))%R; [ intros H'3 | idtac ].
2: apply Rlt_le_trans with (1 + 0)%R;
[ ring_simplify (1 + 0)%R; auto with real | auto with real zarith ].
apply Rmult_le_reg_l with (powerRZ radix precision - 1)%R; auto with real.
apply
Rle_trans
with
(/ 2%nat ×
((powerRZ radix precision - 1) × / (powerRZ radix precision - 1) ×
((1 + - powerRZ radix (Zsucc (- precision))) ×
(powerRZ radix (Zpred (- dExp b)) ×
/ (1 + powerRZ radix (- precision))))) +
/ 2%nat ×
((powerRZ radix precision - 1) × / (powerRZ radix precision - 1) ×
powerRZ radix (- dExp b)))%R; [ right; ring | idtac ].
repeat rewrite Rinv_r; auto with real.
apply Rmult_le_reg_l with (1 + powerRZ radix (- precision))%R; auto with real.
apply
Rle_trans
with
(/ 2%nat ×
((1 + - powerRZ radix (Zsucc (- precision))) ×
(powerRZ radix (Zpred (- dExp b)) ×
((1 + powerRZ radix (- precision)) ×
/ (1 + powerRZ radix (- precision))))) +
/ 2%nat × ((1 + powerRZ radix (- precision)) × powerRZ radix (- dExp b)))%R;
[ right; ring; ring | idtac ].
rewrite Rinv_r; auto with real.
apply
Rle_trans
with
(powerRZ radix (Zpred (Zpred (- dExp b))) ×
(3 ×
((1 + powerRZ radix (- precision)) × (powerRZ radix precision - 1) ×
/ (powerRZ radix precision - powerRZ radix (- precision)))))%R;
[ idtac | right; ring ].
replace ((1 + powerRZ radix (- precision)) × (powerRZ radix precision - 1))%R
with (powerRZ radix precision - powerRZ radix (- precision))%R.
rewrite Rinv_r; auto with real.
right; unfold Zsucc, Zpred in |- ×.
repeat rewrite powerRZ_add; auto with real zarith.
replace (powerRZ radix (-1)) with (/ 2%nat)%R;
[ idtac | simpl in |- *; auto with real zarith ].
simpl; field.
ring_simplify; repeat rewrite <- powerRZ_add; auto with real zarith.
ring_simplify (precision + - precision)%Z; simpl in |- *; ring.
apply Rmult_le_compat_l; auto with real arith.
apply Rmult_le_compat_l; auto with real arith.
unfold Rminus in |- *; apply Rplus_le_compat_r.
apply Rmult_le_compat_r; auto with real.
apply Rplus_le_reg_l with (powerRZ radix (Zsucc (- precision))).
ring_simplify.
replace 1%R with (powerRZ radix 0); auto with real zarith.
fold
(Rabs (FtoRradix a × FtoRradix x + FtoRradix y) ×
/ (1 + powerRZ radix (- precision)) -
powerRZ radix (Zpred (- dExp b)) × / (1 + powerRZ radix (- precision)))%R
in |- ×.
unfold FtoRradix in |- *; apply RoundGeGeneral; auto.
apply Rmult_le_compat_l; auto with real arith.
apply
Rle_trans
with
(/ (powerRZ radix precision - 1) ×
(Rabs (FtoRradix u) - Fulp b radix precision u))%R.
apply Rmult_le_compat_l; auto with real arith.
apply
Rplus_le_reg_l
with
(Fulp b radix precision u +
-
(Rabs (FtoRradix u) × (1 - powerRZ radix (Zsucc (- precision))) -
powerRZ radix (- dExp b)))%R.
ring_simplify.
apply FulpLeGeneral; auto.
apply
Rle_trans
with (/ (powerRZ radix precision - 1) × Rabs (FLess b precision u))%R.
apply Rmult_le_compat_l; auto with real arith.
unfold FtoRradix in |- *; apply UlpFlessuGe_aux; auto.
apply Rmult_le_reg_l with (powerRZ radix precision - 1)%R; auto with real.
rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real.
ring_simplify (1 × Rabs (FtoRradix (FLess b precision u)))%R.
unfold FtoRradix in |- *; apply FulpGe; auto with zarith.
unfold FLess in |- *; case Rcase_abs; intros H3; auto with float zarith.
Qed.
End AxpyFmac.