Library Float.Expansions.FexpAdd
Require Export Fexp2.
Section FexpAdd.
Variable b : Fbound.
Variable precision : nat.
Let radix := 2%Z.
Let TMTO : (1 < radix)%Z := TwoMoreThanOne.
Hint Resolve TMTO: zarith.
Let FtoRradix := FtoR radix.
Coercion FtoRradix : float >-> R.
Hypothesis precisionGreaterThanOne : 1 < precision.
Hypothesis pGivesBound : Zpos (vNum b) = Zpower_nat radix precision.
Hypothesis Ngd : (1 <= pPred (vNum b) * (1 - / radix))%R.
Hypothesis Ngd2 : (6%nat <= pPred (vNum b) * (1 - / radix * / radix))%R.
Inductive NearEqual : list float -> list float -> Prop :=
| IsEqual : forall x : list float, NearEqual x x
| OneMoreR :
forall (x : list float) (e : float),
Fbounded b e -> NearEqual x (e :: x).
Hint Resolve IsEqual OneMoreR.
Definition Step (x y i x' y' : float) (input output output' : list float) :=
Fbounded b x' /\
Fbounded b y' /\
NearEqual output output' /\
(x + (y + (sum (i :: input) + sum output)))%R =
(x' + (y' + (sum input + sum output')))%R :>R /\
(Fexp i <= Fexp y')%Z /\
(Fexp y' <= Fexp x')%Z /\
IsExp b (i :: input) /\
(x' = 0%R :>R \/
y' = 0%R :>R \/
(Rabs y' <=
pPred (vNum b) * (Float 1%nat (Fexp x') - Float 1%nat (hdexp b input)))%R /\
(Rabs y' <=
pPred (vNum b) * (Float 1%nat (Fexp x') - Float 1%nat (Fexp i)))%R) /\
(output' = output \/
(Rabs (x' + y') <= 3%nat * radix * / pPred (vNum b) * Rabs (hd b output'))%R) /\
(output' = nil \/
(Rabs (hd b input) <=
3%nat * (radix * (Rabs (hd b output') * / (pPred (vNum b) - 1%nat))))%R /\
(input = nil \/
(Float (pPred (vNum b)) (Fexp (hd b input)) <=
3%nat * (radix * (Rabs (hd b output') * / (pPred (vNum b) - 1%nat))))%R /\
(Rabs (sum input) <=
length input *
(3%nat * (radix * (Rabs (hd b output') * / (pPred (vNum b) - 1%nat)))))%R)).
Theorem Rle_mult_pos :
forall r1 r2 : R, (0 <= r1)%R -> (0 <= r2)%R -> (0 <= r1 * r2)%R.
intros r1 r2 H H0; replace 0%R with (0 * 0)%R; auto with real.
Qed.
Theorem AddStep :
forall (x y i : float) (input output : list float),
Fbounded b x ->
Fbounded b y ->
IsExp b (i :: input) ->
(Fexp i <= Fexp y)%Z ->
(Fexp y <= Fexp x)%Z ->
(FtoR radix x = 0%R :>R \/ FtoR radix y = 0%R :>R \/ FtoR radix i = 0%R :>R) \/
(Rabs (FtoR radix y) <=
pPred (vNum b) *
(FtoR radix (Float 1%nat (Fexp x)) - FtoR radix (Float 1%nat (Fexp i))))%R ->
output = nil \/
(Rabs i <=
3%nat * (radix * (Rabs (hd b output) * / (pPred (vNum b) - 1%nat))))%R /\
(Float (pPred (vNum b)) (Fexp i) <=
3%nat * (radix * (Rabs (hd b output) * / (pPred (vNum b) - 1%nat))))%R /\
(Rabs (sum (i :: input)) <=
length (i :: input) *
(3%nat * (radix * (Rabs (hd b output) * / (pPred (vNum b) - 1%nat)))))%R ->
exists x' : float,
(exists y' : float,
(exists output' : list float, Step x y i x' y' input output output')).
cut (1 < pPred (vNum b))%Z;
[ intros Z1
| apply Zlt_trans with radix; auto with zarith;
apply (pPredMoreThanRadix b radix precision); auto with zarith ].
cut (0 < pPred (vNum b) - 1)%R;
[ intros Z2 | replace 1%R with (IZR 1); auto with real zarith ].
cut (0 < pPred (vNum b))%Z; [ intros Z3 | auto with real zarith ].
intros x y i input output H H0 H1 H2 H3 H4 G.
cut (Fbounded b i); [ intros K | inversion H1; auto ].
cut (TotalP (Closest b radix));
[ intros ExC | apply (ClosestTotal b radix precision); auto ].
cut (RoundedModeP b radix (Closest b radix));
[ intros Rc | apply (ClosestRoundedModeP b radix precision); auto ].
case (ExC (y + i)%R); intros u Hu.
cut (Fbounded b u);
[ intros Bu
| apply RoundedModeBounded with radix (Closest b radix) (y + i)%R; auto ].
case (ExC (x + u)%R); intros p' Hp'.
cut (Fbounded b p');
[ intros Bp'
| apply RoundedModeBounded with radix (Closest b radix) (x + u)%R; auto ].
case (ExC (x + u - p' + (y + i - u))%R); intros q' Hq'.
cut (Fbounded b q');
[ intros Bq'
| apply
RoundedModeBounded
with radix (Closest b radix) (x + u - p' + (y + i - u))%R;
auto ].
case
(errorBoundedPlus b radix precision TMTO precisionGreaterThanOne pGivesBound
y i u); auto.
intros v (H'1, (H'2, H'3)).
case
(errorBoundedPlus b radix precision TMTO precisionGreaterThanOne pGivesBound
x u p'); auto.
intros w (H'4, (H'5, H'6)).
case
(errorBoundedPlus b radix precision TMTO precisionGreaterThanOne pGivesBound
w v q'); auto.
replace (FtoR radix v) with (FtoR radix y + FtoR radix i - FtoR radix u)%R;
replace (FtoR radix w) with (FtoR radix x + FtoR radix u - FtoR radix p')%R;
auto with real.
intros r' (H'7, (H'8, H'9)).
case
(ThreeSumLoop b precision precisionGreaterThanOne pGivesBound x y i u v w p'
q' r'); auto.
fold radix in |- *;
replace (FtoR radix v) with (FtoR radix y + FtoR radix i - FtoR radix u)%R;
replace (FtoR radix w) with (FtoR radix x + FtoR radix u - FtoR radix p')%R;
auto with real.
intros p''
(q'',
(r'',
((H''0, (H''1, H''2)),
(((H''3, (H''4, H''5)), (H''6, ((H''7, H''8), H''9))),
[(Hr', H''10)| (Hr', H''10)])))).
exists p''; exists q''; exists output; repeat (split; auto).
replace (p'' + (q'' + (sum input + sum output)))%R with
(p'' + (q'' + r'') + (sum input + sum output))%R;
[ idtac | ring_simplify ].
replace (p'' + (q'' + r''))%R with (x + (y + i))%R;
replace (sum (i :: input)) with (i + sum input)%R.
ring.
simpl in |- *; auto with real.
replace (FtoRradix r'') with (FtoRradix r');
replace (FtoRradix r') with (w + v - q')%R; auto with real.
replace (FtoRradix q') with (FtoRradix q''); auto with real.
replace (FtoRradix p'') with (FtoRradix p'); auto with real.
replace (FtoRradix w) with (FtoR radix x + FtoR radix u - FtoR radix p')%R.
replace (FtoRradix v) with (y + i - u)%R; auto with real.
unfold FtoRradix in |- *; ring; ring.
simpl in |- *; auto with real.
replace (FtoRradix r'') with (FtoRradix r'); replace (FtoRradix r') with 0%R;
auto with real; ring.
apply Zle_trans with (Fexp r''); auto.
case (Req_dec p'' 0); intros; auto.
case H''10; auto; intros H''11.
right; right; split; auto.
apply
Rle_trans
with
(pPred (vNum b) *
(FtoR radix (Float 1%nat (Fexp p'')) - FtoR radix (Float 1%nat (Fexp i))))%R;
auto.
apply Rmult_le_compat_l; auto with real arith.
fold FtoRradix in |- *;
replace (Float 1%nat (Fexp p'') - Float 1%nat (hdexp b input))%R with
(Float 1%nat (Fexp p'') + - Float 1%nat (hdexp b input))%R;
replace (Float 1%nat (Fexp p'') - Float 1%nat (Fexp i))%R with
(Float 1%nat (Fexp p'') + - Float 1%nat (Fexp i))%R;
try ring.
apply Rplus_le_compat_l.
apply Ropp_le_contravar; unfold FtoRradix in |- *; unfold FtoR in |- *;
simpl in |- *.
repeat rewrite Rmult_1_l.
apply Rle_powerRZ; auto with real zarith.
apply IsExpZle; auto.
case G; clear G; auto.
generalize H1; clear H1; case input.
simpl in |- *; intros H1 (Hc1, (Hc2, Hc3)); right.
split; auto.
unfold FtoRradix, FtoR, Fzero in |- *; simpl in |- *; rewrite Rmult_0_l;
rewrite Rabs_R0; repeat apply Rle_mult_pos; auto with real zarith.
intros f l H1 (Hc1, (Hc2, Hc3)); right.
cut (Fbounded b f); [ intros Bf | idtac ].
cut
(Float (pPred (vNum b)) (Fexp f) <=
3%nat *
(radix * (Rabs (FtoRradix (hd b output)) * / (pPred (vNum b) - 1%nat))))%R;
[ intros HZ1 | idtac ].
cut
(FtoRradix (Float (pPred (vNum b)) (Fexp f)) <=
FtoRradix (Float (pPred (vNum b)) (Fexp i)))%R; [ intros HZ2 | idtac ].
split.
apply Rle_trans with (2 := HZ1); auto with real zarith.
simpl in |- *; (rewrite <- (Fabs_correct radix); auto with real zarith);
apply (maxMax1 radix); auto with real zarith.
right; split; auto.
apply sum_IsExp with b precision; auto.
inversion H1; auto with real zarith.
unfold FtoRradix, FtoR in |- *; simpl in |- *; apply Rmult_le_compat_l;
auto with real zarith.
apply Rle_powerRZ; auto with real zarith.
inversion H1; auto with real zarith.
apply Rle_trans with (2 := Hc2); auto with real zarith.
unfold FtoRradix, FtoR in |- *; simpl in |- *; apply Rmult_le_compat_l;
auto with real zarith.
apply Rle_powerRZ; auto with real zarith.
inversion H1; auto with real zarith.
inversion H1; auto.
exists q''; exists r''; exists (p'' :: output); repeat (split; auto).
simpl in |- *.
replace
(FtoRradix q'' + (FtoRradix r'' + (sum input + (FtoR 2 p'' + sum output))))%R
with (p'' + (q'' + r'') + (sum input + sum output))%R;
[ idtac | ring_simplify ].
replace (p'' + (q'' + r''))%R with (x + (y + i))%R; [ fold radix; fold FtoRradix; ring | idtac ].
replace (FtoRradix r'') with (FtoRradix r');
replace (FtoRradix r') with (w + v - q')%R; auto with real.
replace (FtoRradix w) with (FtoR radix x + FtoR radix u - FtoR radix p')%R.
replace (FtoRradix v) with (y + i - u)%R; auto with real.
replace (FtoRradix p'') with (FtoR 2 p');
replace (FtoRradix q'') with (FtoRradix q'); fold radix FtoRradix in |- *;
ring.
simpl in |- *; fold radix FtoRradix in |- *; ring.
right; right; split; auto.
apply Rle_trans with (1 := H''10).
apply Rmult_le_compat_l; auto with real zarith.
unfold Rminus in |- *; apply Rplus_le_compat_l; apply Ropp_le_contravar.
apply (oneExp_le radix); auto with real zarith.
apply IsExpZle; auto.
right; simpl in |- *.
apply Rlt_le; fold radix in |- *; unfold FtoRradix in |- *;
apply
bound3Sum
with
(precision := precision)
(w := w)
(v := v)
(p := x)
(q := y)
(r := i)
(u := u); auto; fold radix FtoRradix in |- *.
apply (ClosestCompatible b radix (x + u)%R (x + u)%R p'); auto with real.
replace (FtoRradix p'') with (FtoR radix p'); auto with real.
apply (ClosestCompatible b radix (w + v)%R (w + v)%R q'); auto with real.
replace (FtoRradix w) with (x + u - p')%R;
replace (FtoRradix v) with (y + i - u)%R; auto with real.
replace (FtoRradix r'') with (FtoRradix r');
replace (FtoRradix q'') with (FtoRradix q'); auto with real.
replace (FtoRradix r'') with (FtoRradix r'); auto with real.
right; simpl in |- *.
cut
(forall z : float,
(Fexp z <= Fexp i)%Z ->
Fbounded b z ->
(Rabs z <= 3%nat * (radix * (Rabs p'' * / (pPred (vNum b) - 1%nat))))%R);
[ intros Hz1 | idtac ].
split; [ apply Hz1; auto with zarith | idtac ].
replace (Fexp (hd b input)) with (hdexp b input).
apply IsExpZle; auto.
case input; simpl in |- *; auto with zarith.
generalize H1; case input; intros; simpl in |- *.
apply FboundedFzero.
inversion H5; auto.
generalize H1; clear H1; case input; auto.
intros f l H1; right.
cut
(Float (pPred (vNum b)) (Fexp f) <=
3%nat * (radix * (Rabs p'' * / (pPred (vNum b) - 1%nat))))%R;
[ intros Hz2 | idtac ].
split.
simpl in |- *; auto.
apply sum_IsExp with b precision; auto.
inversion H1; auto.
rewrite <- (fun x : float => Rabs_pos_eq x); auto with real zarith.
apply Hz1; auto with real zarith.
simpl in |- *; inversion H1; auto.
unfold pPred in |- *; apply maxFbounded; auto.
cut (Fbounded b f); [ intros tmp; case tmp | idtac ]; auto.
cut (IsExp b (f :: l)); [ intros tmp; inversion tmp | inversion H1 ]; auto.
auto with float zarith.
apply (LeFnumZERO radix); auto with real zarith.
intros z H5 H6;
apply Rle_trans with (FtoRradix (Float (pPred (vNum b)) (Fexp z))).
rewrite <- (Fabs_correct radix); auto with real zarith; apply (maxMax1 radix);
auto with real zarith.
replace (FtoRradix (Float (pPred (vNum b)) (Fexp z))) with
(pPred (vNum b) * powerRZ radix (Fexp z))%R;
[ idtac | auto with real zarith ].
replace (3%nat * (radix * (Rabs p'' * / (pPred (vNum b) - 1%nat))))%R with
(pPred (vNum b) *
(3%nat *
(radix *
(radix *
(Rabs p'' * / (pPred (vNum b) * (radix * pPred (vNum b) - radix)))))))%R.
apply Rmult_le_compat_l; auto with real zarith.
apply Rle_trans with (FtoRradix (Float 1%nat (Fexp i))).
unfold FtoRradix in |- *; unfold FtoR in |- *; simpl in |- *;
rewrite Rmult_1_l; auto with real zarith.
apply Rle_powerRZ; auto with real zarith.
apply Rlt_le;
apply
OutSum3
with
(precision := precision)
(p := x)
(q := y)
(r := i)
(u := u)
(v := v)
(w := w)
(p' := p'')
(q' := q'')
(r' := r''); auto with real arith; fold radix FtoRradix in |- *.
apply (ClosestCompatible b radix (x + u)%R (x + u)%R p'); auto with real.
replace (FtoRradix p'') with (FtoRradix p'); auto with real.
apply (ClosestCompatible b radix (w + v)%R (w + v)%R q'); auto with real.
replace (FtoRradix w) with (x + u - p')%R;
replace (FtoRradix v) with (y + i - u)%R; auto with real.
replace (FtoRradix r'') with (FtoRradix r');
replace (FtoRradix q'') with (FtoRradix q'); auto with real.
replace (FtoRradix r'') with (FtoRradix r'); auto with real.
replace (radix * pPred (vNum b) - radix)%R with
(radix * (pPred (vNum b) - 1%nat))%R; [ idtac | simpl; ring ].
repeat rewrite Rinv_mult_distr; auto with real zarith.
pattern (INR 3) at 2 in |- *;
replace (INR 3) with
(radix * / radix * (pPred (vNum b) * / pPred (vNum b) * 3%nat))%R;
[ simpl in |- *; ring | idtac ].
repeat rewrite Rinv_r; auto with real zarith; simpl in |- *; ring.
Qed.
Variable input : list float.
Inductive IsRleEpsExp : list float -> Prop :=
| IsRleEpsExpNil : IsRleEpsExp nil
| IsRleEpsExpSingle :
forall x : float, Fbounded b x -> IsRleEpsExp (x :: nil)
| IsRleEpsExpTop :
forall (x y : float) (L : list float),
Fbounded b x ->
Fbounded b y ->
(Rabs x <=
(6%nat * length input + 6%nat) *
/ (pPred (vNum b) - 1%nat - 6%nat * length input) *
Rabs y)%R -> IsRleEpsExp (y :: L) -> IsRleEpsExp (x :: y :: L).
Theorem endof_Rle_length :
forall (P Q : list float) (k : float),
endof input (P ++ k :: Q) -> (length P <= length input - 1)%R.
intros P Q k H; inversion H.
rewrite H0; repeat rewrite app_length; simpl in |- *.
replace 1%R with (INR 1); [ idtac | simpl in |- *; auto ].
replace (length x + (length P + S (length Q))) with
(S (length P) + (length Q + length x));
[ idtac
| repeat rewrite (fun x : list float => S_to_plus_one (length x)); ring ].
rewrite <- minus_INR; simpl in |- *; try apply Rle_INR; auto with arith.
rewrite <- minus_n_O; auto with arith.
Qed.
Theorem FexpAdd_aux :
(6%nat * length input * / (pPred (vNum b) - 1%nat) < 1)%R ->
forall (L output : list float) (x y : float),
endof input L ->
IsExp b L ->
IsRleEpsExp output ->
output = nil \/
(Rabs (x + y) <=
(pPred (vNum b) - 1%nat) * / pPred (vNum b) *
((6%nat * length input + 6%nat) *
/ (pPred (vNum b) - 1%nat - 6%nat * length input)) *
Rabs (hd b output))%R ->
output = nil \/
(exists L1 : list float,
(exists x1 : float,
(exists y1 : float,
IsExp b (L1 ++ L) /\
endof input (L1 ++ L) /\
(x1 + (y1 + sum L1))%R = (x + y)%R :>R /\
(Rabs (sum L1) <=
length L1 *
(3%nat *
(radix * (Rabs (hd b output) * / (pPred (vNum b) - 1%nat)))))%R /\
(Rabs (x1 + y1) <=
3%nat * radix * / pPred (vNum b) * Rabs (hd b output))%R))) ->
Fbounded b x ->
Fbounded b y ->
(hdexp b L <= Fexp y)%Z ->
(Fexp y <= Fexp x)%Z ->
(FtoR radix x = 0%R :>R \/ FtoR radix y = 0%R :>R) \/
(Rabs (FtoR radix y) <=
pPred (vNum b) *
(FtoR radix (Float 1%nat (Fexp x)) - FtoR radix (Float 1%nat (hdexp b L))))%R ->
output = nil \/
L = nil \/
(Rabs (hd b L) <=
3%nat * (radix * (Rabs (hd b output) * / (pPred (vNum b) - 1%nat))))%R /\
(Float (pPred (vNum b)) (Fexp (hd b L)) <=
3%nat * (radix * (Rabs (hd b output) * / (pPred (vNum b) - 1%nat))))%R /\
(Rabs (sum L) <=
length L *
(3%nat * (radix * (Rabs (hd b output) * / (pPred (vNum b) - 1%nat)))))%R ->
exists x' : float,
(exists y' : float,
(exists output' : list float,
Fbounded b x' /\
Fbounded b y' /\
(x + (y + (sum output + sum L)))%R = (x' + (y' + sum output'))%R :>R /\
length output' <= length L + length output /\
(Fexp y' <= Fexp x')%Z /\
IsRleEpsExp output' /\
(output = nil \/
(exists L1 : list float,
(exists x1 : float,
(exists y1 : float,
IsExp b (L1 ++ L) /\
endof input (L1 ++ L) /\
(x1 + (y1 + sum L1))%R = (x + y)%R :>R /\
(Rabs (sum L1) <=
length L1 *
(3%nat *
(radix *
(Rabs (hd b output) * / (pPred (vNum b) - 1%nat)))))%R /\
(Rabs (x1 + y1) <=
3%nat * radix * / pPred (vNum b) * Rabs (hd b output))%R)))) /\
endof output' output /\
(output' = nil \/
(Rabs (x' + y') <=
(pPred (vNum b) - 1%nat) * / pPred (vNum b) *
((6%nat * length input + 6%nat) *
/ (pPred (vNum b) - 1%nat - 6%nat * length input)) *
Rabs (hd b output'))%R))).
intros Z0.
cut (1 < pPred (vNum b))%Z;
[ intros Z1
| apply Zlt_trans with radix;
try apply (pPredMoreThanRadix b radix precision);
auto with zarith ].
cut (0 < pPred (vNum b))%Z;
[ intros Zl2 | apply Zlt_1_O; apply Zlt_le_weak; auto with zarith ].
cut (0 < pPred (vNum b) - 1%nat - 6%nat * length input)%R;
[ intros Z3
| apply Rlt_Rminus_ZERO;
apply Rmult_lt_reg_l with (/ (pPred (vNum b) - 1%nat))%R;
try rewrite Rinv_l; try rewrite (fun x => Rmult_comm (/ x));
auto with real zarith; replace (INR 1) with (IZR 1);
try rewrite Z_R_minus; auto with real zarith ].
cut
(3%nat * radix * / pPred (vNum b) <=
(pPred (vNum b) - 1%nat) * / pPred (vNum b) *
((6%nat * length input + 6%nat) *
/ (pPred (vNum b) - 1%nat - 6%nat * length input)))%R;
[ intros Z4 | auto with real zarith ].
2: replace (3%nat * radix)%R with (INR 6 * 1)%R;
[ idtac | rewrite Rmult_1_r; simpl in |- *; ring ].
2: repeat rewrite Rmult_assoc.
2: rewrite (Rmult_comm 1).
2: rewrite (Rmult_comm (pPred (vNum b) - 1%nat)).
2: replace (6%nat * length input + 6%nat)%R with
(INR 6 * (length input + 1%nat))%R; [ idtac | simpl in |- *; ring ].
2: repeat rewrite <- (Rmult_assoc (/ pPred (vNum b))).
2: repeat rewrite (Rmult_comm (/ pPred (vNum b)) 6%nat).
2: repeat rewrite Rmult_assoc.
2: repeat apply Rmult_le_compat_l; auto with real zarith.
2: rewrite (Rmult_comm (length input + 1%nat)).
2: apply
Rmult_le_reg_l with (pPred (vNum b) - 1%nat - 6%nat * length input)%R;
auto with real zarith.
2: repeat rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real zarith;
rewrite Rmult_1_r; rewrite Rmult_1_l.
2: rewrite (fun x => Rplus_comm x 1); rewrite Rmult_plus_distr_l;
rewrite Rmult_1_r.
2: unfold Rminus in |- *; apply Rplus_le_compat_l.
2: rewrite <- Ropp_mult_distr_l_reverse; apply Rmult_le_compat_r;
auto with real arith.
2: apply Rle_trans with (-0)%R; auto with real zarith arith.
2: rewrite Ropp_0; change (0 <= pPred (vNum b) - 1%nat)%R in |- *;
auto with real zarith.
intros L; elim L; clear L.
intros output x y H H0 H1 H2 H3 H4 H5 H6 H7 H8 H9; exists x; exists y;
exists output; repeat (split; simpl in |- *; auto);
[ ring | exists (nil (A:=float)) ]; auto.
intros a L HrecL output; case output; clear output.
intros x y H H0 H1 H2 H3 H4 H5 H6 H7 H8 H9.
cut (IsExp b L); [ intros Hp1 | inversion H0; auto; apply IsExpNil ].
cut (endof input L);
[ intros Hp2
| case H; intros LL H10; exists (LL ++ a :: nil); rewrite H10;
replace (a :: L) with ((a :: nil) ++ L); try apply ass_app;
simpl in |- *; auto ].
case (AddStep x y a L nil); auto with real zarith.
inversion H0; auto.
case H8; auto; intros tmp; case tmp; auto.
case H8; auto; intros tmp; case tmp; auto.
intros x' (y', (output', H'1)).
case H'1; intros T1 (T2, (T3, (T4, (T5, (T6, (T7, (T8, (T9, T10)))))))).
case (HrecL output' x' y'); clear H2 H3 H9 HrecL; auto with real zarith.
inversion T3; auto.
apply IsRleEpsExpSingle; auto.
CaseEq output'; auto; intros o output'' Eq1.
right; case T9; rewrite Eq1;
[ intros; discriminate | simpl in |- *; clear T9; intros T9 ].
apply
Rle_trans with (3%nat * radix * / pPred (vNum b) * Rabs (FtoRradix o))%R;
auto with real.
CaseEq output'; auto; intros o output'' Eq1.
right; exists (nil (A:=float)); exists x'; exists y';
repeat (split; simpl in |- *; auto with real).
rewrite Rabs_R0; rewrite Rmult_0_l; auto with real.
case T9; rewrite Eq1; [ intros; discriminate | simpl in |- *; auto ].
apply Zle_trans with (Fexp a); auto with zarith.
apply IsExpZle; auto.
repeat (case T8; auto; clear T8; intros T8).
case T10; auto; clear T10; intros (T10, T11); right; case T11; auto.
intros x'0
(y'0, (output'0, (V1, (V2, (V3, (V4, (V5, (V6, (V7, (V8, V9)))))))))).
exists x'0; exists y'0; exists output'0; repeat (split; auto).
rewrite (Rplus_comm (sum nil)); rewrite T4; rewrite <- V3; ring.
generalize V4; case T3; simpl in |- *; auto with arith; intros e Be V0;
rewrite plus_n_Sm; simpl in |- *; auto.
exists output'0; apply app_nil_end.
intros a0 output x y V' H H0 H1 M H2 H3 H4 H5 H6 H7.
cut (IsExp b L); [ intros Hp1 | inversion H; auto; apply IsExpNil ].
cut (endof input L);
[ intros Hp2
| case V'; intros LL H10; exists (LL ++ a :: nil); rewrite H10;
replace (a :: L) with ((a :: nil) ++ L); try apply ass_app;
simpl in |- *; auto ].
case H1; [ intros; discriminate | clear H1; intros H1 ].
case M; [ intros; discriminate | clear M; intros M ].
case H7; [ intros; discriminate | clear H7; intros H7 ].
case H7; [ intros; discriminate | clear H7; intros (H7, (H8, H9)) ].
case (AddStep x y a L (a0 :: output)); auto with real zarith.
repeat (case H6; auto; clear H6; intros H6).
intros x'
(y', (output', (V1, (V2, (V3, (V4, (V5, (V6, (V7, (V8, (V9, V10))))))))))).
cut (IsRleEpsExp output'); [ intros RR | idtac ].
case (HrecL output' x' y'); clear H2 H3 H9 HrecL; auto with real zarith.
right; auto.
4: repeat (case V8; auto; clear V8; intros V8).
3: apply Zle_trans with (Fexp a); auto with zarith.
3: apply IsExpZle; auto.
3: case V10; auto; try intros (tmp1, [tmp2| tmp3]); auto.
3: intros x'0
(y'0, (output'0, (W1, (W2, (W3, (W4, (W5, (W6, (W7, (W8, W9)))))))))).
3: exists x'0; exists y'0; exists output'0; repeat (split; auto).
3: rewrite <- W3; rewrite (Rplus_comm (sum output')); rewrite <- V4; ring.
3: apply le_trans with (1 := W4).
3: inversion V3; simpl in |- *; auto with arith; repeat rewrite plus_n_Sm;
simpl in |- *; auto with arith.
3: inversion V3; [ replace (a0 :: output) with output'; auto | idtac ].
3: rewrite <- H9 in W8; inversion W8; exists (x1 ++ e :: nil); rewrite H10;
replace (e :: a0 :: output) with ((e :: nil) ++ a0 :: output);
try (apply sym_equal; apply app_ass); simpl in |- *;
auto.
2: case M; intros L1 (x1, (y1, (X1, (X2, (X3, (X4, X5)))))).
2: inversion V3; auto.
2: right; exists (L1 ++ a :: nil); exists x1; exists y1; repeat (split; auto).
2: rewrite app_ass; simpl in |- *; auto.
2: rewrite app_ass; simpl in |- *; auto.
2: apply trans_eq with (FtoRradix x + FtoRradix y + a)%R;
fold FtoRradix in |- *;
[ rewrite <- X3; rewrite <- sum_app; fold radix; fold FtoRradix; ring | idtac ].
2: replace (FtoRradix x + FtoRradix y)%R with
(FtoRradix x + (FtoRradix y + (sum (a :: L) + sum (a0 :: output))) +
- (sum (a :: L) + sum (a0 :: output)))%R;
[ rewrite V4; rewrite H3; simpl in |- *; fold radix FtoRradix in |- *;
ring
| fold FtoRradix in |- *; ring ].
2: rewrite <- sum_app.
2: apply Rle_trans with (1 := Rabs_triang (FtoR 2 a) (sum L1)).
2: rewrite app_length; simpl in |- *.
2: rewrite plus_INR; simpl in |- *; rewrite Rmult_plus_distr_r;
rewrite Rmult_1_l.
2: rewrite (fun x => Rplus_comm (Rabs x)); apply Rplus_le_compat; auto.
2: right; exists (nil (A:=float)); exists x'; exists y';
repeat (split; simpl in |- *; auto).
2: rewrite Rplus_0_r; auto.
2: rewrite Rabs_R0; rewrite Rmult_0_l; auto with real.
2: case V9; rewrite <- H9; simpl in |- *; auto.
2: intros tmp; Contradict tmp; apply cons_neq; auto.
inversion V3.
case M; intros L1 (x1, (y1, (X1, (X2, (X3, (X4, X5)))))).
replace (FtoRradix x' + FtoRradix y')%R with
(a + (FtoRradix x1 + FtoRradix y1 + sum L1))%R.
apply
Rle_trans
with
(Rabs (FtoRradix a) +
(Rabs (FtoRradix x1 + FtoRradix y1) + Rabs (sum L1)))%R.
apply
Rle_trans
with (1 := Rabs_triang (FtoR 2 a) (FtoRradix x1 + FtoRradix y1 + sum L1));
auto with real.
apply Rplus_le_compat_l; apply Rabs_triang; auto with real.
apply
Rle_trans
with
(3%nat *
(2%nat *
(Rabs (FtoRradix (hd b (a0 :: output))) * / (pPred (vNum b) - 1%nat))) +
(Rabs (x1 + y1) + Rabs (sum L1)))%R;
[ apply Rplus_le_compat_r; auto | idtac ].
apply
Rle_trans
with
(3%nat *
(2%nat *
(Rabs (FtoRradix (hd b (a0 :: output))) * / (pPred (vNum b) - 1%nat))) +
(3%nat * 2%nat * / pPred (vNum b) *
Rabs (FtoRradix (hd b (a0 :: output))) + Rabs (sum L1)))%R;
[ apply Rplus_le_compat_l; apply Rplus_le_compat_r; auto | idtac ].
apply
Rle_trans
with
(3%nat *
(2%nat *
(Rabs (FtoRradix (hd b (a0 :: output))) * / (pPred (vNum b) - 1%nat))) +
(3%nat * 2%nat * / pPred (vNum b) *
Rabs (FtoRradix (hd b (a0 :: output))) +
length L1 *
(3%nat *
(2%nat *
(Rabs (FtoRradix (hd b (a0 :: output))) * / (pPred (vNum b) - 1%nat))))))%R;
[ apply Rplus_le_compat_l; apply Rplus_le_compat_l; auto | idtac ].
apply
Rle_trans
with
(3%nat *
(2%nat *
(Rabs (FtoRradix (hd b (a0 :: output))) * / (pPred (vNum b) - 1%nat))) +
(3%nat * 2%nat * / (pPred (vNum b) - 1%nat) *
Rabs (FtoRradix (hd b (a0 :: output))) +
length L1 *
(3%nat *
(2%nat *
(Rabs (FtoRradix (hd b (a0 :: output))) * / (pPred (vNum b) - 1%nat))))))%R;
[ apply Rplus_le_compat_l; apply Rplus_le_compat_r; simpl in |- *;
auto with real
| idtac ].
apply Rmult_le_compat_r; auto with real.
apply Rmult_le_compat_l; auto with real.
replace ((2 + 1) * 2)%R with (INR 6); auto with real arith; simpl in |- *;
ring.
apply Rle_Rinv; auto with real arith.
pattern (IZR (pPred (vNum b))) at 2 in |- *;
replace (IZR (pPred (vNum b))) with (IZR (pPred (vNum b)) + -0)%R;
[ unfold Rminus in |- *; auto with real | ring ].
apply
Rle_trans
with
((3%nat * 2%nat + (3%nat * 2%nat + length L1 * (3%nat * 2%nat))) *
/ (pPred (vNum b) - 1%nat) * Rabs (FtoRradix (hd b (a0 :: output))))%R;
[ right; ring; ring | idtac ].
apply Rmult_le_compat_r; auto with real.
replace (3%nat * 2%nat)%R with (INR 6); [ idtac | simpl in |- *; ring ].
apply
Rle_trans
with
((6%nat + (6%nat + (length input - 1) * 6%nat)) *
/ (pPred (vNum b) - 1%nat))%R.
apply Rmult_le_compat_r; auto with real.
apply Rplus_le_compat_l; apply Rplus_le_compat_l; apply Rmult_le_compat_r;
auto with arith zarith real.
apply endof_Rle_length with L a; auto.
replace (6%nat + (length input - 1) * 6%nat)%R with (6%nat * length input)%R;
[ idtac | ring; ring ].
apply
Rle_trans
with
((6%nat + 6%nat * length input) *
((pPred (vNum b) - 1%nat) * / pPred (vNum b) *
/ (pPred (vNum b) - 1%nat - 6%nat * length input)))%R;
[ idtac | right; ring; ring ].
apply Rmult_le_compat_l; auto with arith zarith real.
repeat rewrite <- plus_INR || rewrite <- mult_INR; apply pos_INR;
simpl in |- *; auto with arith.
apply Rmult_le_reg_l with (pPred (vNum b) - 1%nat)%R; auto with real.
apply Rmult_le_reg_l with (pPred (vNum b) - 1%nat - 6%nat * length input)%R;
auto with real.
apply Rmult_le_reg_l with (IZR (pPred (vNum b))); auto with real.
rewrite Rinv_r;
[ idtac
| simpl in |- *; replace 1%R with (IZR 1); try rewrite Z_R_minus;
auto with real zarith ].
apply
Rle_trans
with
((pPred (vNum b) - 1%nat) *
((pPred (vNum b) - 1%nat) *
(pPred (vNum b) * / pPred (vNum b) *
((pPred (vNum b) - 1%nat - 6%nat * length input) *
/ (pPred (vNum b) - 1%nat - 6%nat * length input)))))%R;
[ idtac | right; ring; ring ].
rewrite Rinv_r; [ idtac | auto with real zarith ].
rewrite Rinv_r; [ idtac | auto with real zarith ].
apply
Rle_trans
with
(pPred (vNum b) * pPred (vNum b) +
- (pPred (vNum b) * (6%nat * length input + 1)))%R;
[ right; simpl; ring | idtac ].
apply
Rle_trans
with (pPred (vNum b) * pPred (vNum b) + (- (pPred (vNum b) * 2) + 1))%R;
[ idtac | right; simpl; ring ].
apply Rplus_le_compat_l.
apply Rle_trans with (- (pPred (vNum b) * 2))%R.
apply Ropp_le_contravar; apply Rmult_le_compat_l; auto with real zarith.
apply Rle_trans with (INR 6 * INR 1)%R.
replace 2%R with (INR 2); repeat rewrite <- plus_INR || rewrite <- mult_INR;
try apply Rle_INR; simpl in |- *; auto with real arith.
apply Rle_trans with (6%nat * length input)%R; auto with real zarith.
apply Rmult_le_compat_l; auto with real arith.
case X2; intros LL HL; rewrite HL.
repeat rewrite app_length; simpl in |- *; repeat rewrite <- plus_n_Sm;
replace 1%R with (INR 1); auto with real arith.
rewrite <- (Rplus_0_r (- (pPred (vNum b) * 2))); auto with real.
repeat rewrite Rplus_assoc; rewrite X3;
replace (FtoRradix x' + FtoRradix y')%R with
(FtoRradix x' + (FtoRradix y' + (sum L + sum output')) +
- (sum L + sum output'))%R;
[ rewrite <- V4; rewrite <- H3; simpl in |- *; fold FtoRradix in |- *;
fold radix; fold FtoRradix; ring | ring ].
case V9;
[ rewrite <- H9; intros tmp; Contradict tmp; apply cons_neq
| rewrite H9; intros H'9 ].
apply Rle_trans with (1 := H'9); auto with real.
inversion V3; auto.
apply IsRleEpsExpTop; auto.
inversion H0; auto.
apply
Rmult_le_reg_l
with (1 - 6%nat * length (a :: L) * / (pPred (vNum b) - 1%nat))%R.
2: apply
Rle_trans
with
((1 -
(6%nat * / pPred (vNum b) +
6%nat * length L * / (pPred (vNum b) - 1%nat))) *
Rabs e)%R.
2: apply Rmult_le_compat_r; auto with real.
2: unfold Rminus in |- *; apply Rplus_le_compat_l; apply Ropp_le_contravar.
2: replace (6%nat * length (a :: L) * / (pPred (vNum b) + - 1%nat))%R with
(6%nat * / (pPred (vNum b) + -1) +
6%nat * length L * / (pPred (vNum b) + -1))%R;
[ apply Rplus_le_compat_r | idtac ].
2: apply Rmult_le_compat_l; apply Rlt_le; auto with real arith.
2: apply Rinv_1_lt_contravar; auto with real arith.
2: replace 1%R with (IZR 1);
[ rewrite <- Ropp_Ropp_IZR; rewrite <- plus_IZR; apply Rle_IZR
| simpl in |- *; auto ].
2: change (1 <= Zpred (pPred (vNum b)))%Z in |- *; auto with zarith.
2: pattern (IZR (pPred (vNum b))) at 2 in |- *;
replace (IZR (pPred (vNum b))) with (pPred (vNum b) + -0)%R;
auto with real; ring.
2: replace (INR (length (a :: L))) with (1 + INR (length L))%R;
[ simpl in |- *; ring
| simpl in |- *; case (length L); simpl in |- *; intros; ring ].
2: case V9;
[ rewrite <- H12; intros tmp; Contradict tmp; apply cons_neq
| rewrite <- H12; intros H'9 ].
2: replace
((1 -
(6%nat * / pPred (vNum b) +
6%nat * length L * / (pPred (vNum b) - 1%nat))) *
Rabs e)%R with
(Rabs e -
(6%nat * / pPred (vNum b) * Rabs e +
6%nat * length L * / (pPred (vNum b) - 1%nat) * Rabs e))%R;
[ idtac | ring; ring ].
2: apply
Rle_trans
with
(Rabs e -
(Rabs (x' + y') +
6%nat * length L * / (pPred (vNum b) - 1%nat) * Rabs e))%R.
2: unfold Rminus in |- *; apply Rplus_le_compat_l.
2: apply Ropp_le_contravar; apply Rplus_le_compat_r.
2: replace (INR 6) with (3%nat * 2%nat)%R; replace e with (hd b (e :: x0));
auto with real.
2: rewrite <- (mult_INR 3 2); auto with arith.
2: case V10;
[ rewrite <- H12; intros; discriminate
| rewrite <- H12; intros (H'10, H'11) ].
2: apply Rle_trans with (Rabs e - (Rabs (x' + y') + Rabs (sum L)))%R.
2: unfold Rminus in |- *; apply Rplus_le_compat_l; apply Ropp_le_contravar;
apply Rplus_le_compat_l.
2: case H'11;
[ intros H'12; rewrite H'12; simpl in |- *; rewrite Rabs_R0; right; ring
| intros (H'12, H'13); auto ].
2: apply Rle_trans with (1 := H'13); right; simpl in |- *.
2: ring_simplify (pPred (vNum b) +- 1)%R;ring.
elim H2; intros; clear H2.
2: apply Rle_trans with (Rabs (e + (x' + y' + sum L))).
2: rewrite <- (Rabs_Ropp (x' + y')); rewrite <- (Rabs_Ropp (sum L)).
2: replace (Rabs e - (Rabs (- (x' + y')) + Rabs (- sum L)))%R with
(Rabs e - Rabs (- (x' + y')) - Rabs (- sum L))%R;
[ idtac | ring; ring ].
2: apply Rle_trans with (Rabs (e - - (x' + y')) - Rabs (- sum L))%R.
2: unfold Rminus in |- *; apply Rplus_le_compat_r.
2: replace (Rabs e + - Rabs (- (x' + y')))%R with
(Rabs e - Rabs (- (x' + y')))%R;
replace (e + - - (x' + y'))%R with (e - - (x' + y'))%R;
try ring.
2: apply Rabs_triang_inv.
2: replace (e + (x' + y' + sum L))%R with (e - - (x' + y') - - sum L)%R;
try ring.
2: apply Rabs_triang_inv.
2: replace (e + (x' + y' + sum L))%R with (x + (y + sum (a :: L)))%R.
3: apply Rplus_eq_reg_l with (sum (a0 :: output)).
3: replace (sum (a0 :: output) + (x + (y + sum (a :: L))))%R with
(FtoRradix x + (FtoRradix y + (sum (a :: L) + sum (a0 :: output))))%R;
[ rewrite V4; rewrite <- H12; simpl in |- *; fold radix FtoRradix in |- *;
ring
| ring ].
2: case M; intros L1 (x1, (y1, (X1, (X2, (X3, (X4, X5)))))).
2: replace (x + (y + sum (a :: L)))%R with (x + y + sum (a :: L))%R;
[ rewrite <- X3 | ring ].
2: replace (x1 + (y1 + sum L1) + sum (a :: L))%R with
(x1 + y1 + (sum L1 + sum (a :: L)))%R; [ idtac | ring ].
2: apply Rle_trans with (Rabs (x1 + y1) + Rabs (sum L1 + sum (a :: L)))%R;
[ apply Rabs_triang | idtac ].
2: apply
Rle_trans
with
(6%nat * / pPred (vNum b) * Rabs a0 + Rabs (sum L1 + sum (a :: L)))%R.
2: apply Rplus_le_compat_r; apply Rle_trans with (1 := X5).
2: simpl in |- *; repeat apply Rmult_le_compat_r; auto with real.
2: replace ((2 + 1) * 2)%R with (INR 6);
replace (2 + 1 + 1 + 1 + 1)%R with (INR 6); auto with real arith;
simpl in |- *; ring.
2: apply
Rle_trans
with
(6%nat * / pPred (vNum b) * Rabs (FtoRradix a0) +
(Rabs (sum L1) + Rabs (sum (a :: L))))%R;
[ apply Rplus_le_compat_l; apply Rabs_triang | idtac ].
2: apply
Rle_trans
with
(6%nat * / pPred (vNum b) * Rabs (FtoRradix a0) +
(length L1 *
(3%nat *
(2%nat * (Rabs (hd b (a0 :: output)) * / (pPred (vNum b) - 1%nat)))) +
Rabs (sum (a :: L))))%R.
2: apply Rplus_le_compat_l; apply Rplus_le_compat_r; auto.
2: apply
Rle_trans
with
(6%nat * / pPred (vNum b) * Rabs (FtoRradix a0) +
(length L1 *
(3%nat *
(2%nat * (Rabs (hd b (a0 :: output)) * / (pPred (vNum b) - 1%nat)))) +
length (a :: L) *
(3%nat *
(2%nat * (Rabs (hd b (a0 :: output)) * / (pPred (vNum b) - 1%nat))))))%R.
2: apply Rplus_le_compat_l; apply Rplus_le_compat_l; auto.
2: replace (hd b (a0 :: output)) with a0; [ idtac | simpl in |- *; auto ].
2: replace
(length L1 * (3%nat * (2%nat * (Rabs a0 * / (pPred (vNum b) - 1%nat)))) +
length (a :: L) *
(3%nat * (2%nat * (Rabs a0 * / (pPred (vNum b) - 1%nat)))))%R with
((length L1 + length (a :: L)) *
(3%nat * 2%nat * (Rabs a0 * / (pPred (vNum b) - 1%nat))))%R;
[ idtac | ring; ring ].
2: replace (3%nat * 2%nat)%R with (INR 6); [ idtac | simpl in |- *; ring ].
2: replace
(6%nat * / pPred (vNum b) * Rabs a0 +
(length L1 + length (a :: L)) *
(6%nat * (Rabs a0 * / (pPred (vNum b) - 1%nat))))%R with
((6%nat * / pPred (vNum b) +
(length L1 + length (a :: L)) * (6%nat * / (pPred (vNum b) - 1%nat))) *
Rabs a0)%R; [ idtac | ring ].
2: apply
Rle_trans
with
((1 - 6%nat * length (a :: L) * / (pPred (vNum b) - 1%nat)) *
((6%nat * length input + 6%nat) *
/ (pPred (vNum b) - 1%nat - 6%nat * length input)) *
Rabs a0)%R; [ idtac | right; ring ].
2: apply Rmult_le_compat_r; auto with real.
2: apply
Rle_trans
with
(6%nat * / (pPred (vNum b) - 1%nat) +
(length L1 + length (a :: L)) * (6%nat * / (pPred (vNum b) - 1%nat)))%R.
2: apply Rplus_le_compat_r; apply Rmult_le_compat_l; auto with real arith.
2: apply Rlt_le; apply Rinv_1_lt_contravar; auto with real zarith.
2: simpl in |- *; replace 1%R with (IZR 1);
[ unfold Rminus in |- *; rewrite <- Ropp_Ropp_IZR; rewrite <- plus_IZR;
apply Rle_IZR
| simpl in |- *; auto ].
2: change (1 <= Zpred (pPred (vNum b)))%Z in |- *; auto with zarith.
2: simpl in |- *; unfold Rminus in |- *;
pattern (IZR (pPred (vNum b))) at 2 in |- *;
replace (IZR (pPred (vNum b))) with (pPred (vNum b) + -0)%R;
auto with real; ring.
2: replace
(6%nat * / (pPred (vNum b) - 1%nat) +
(length L1 + length (a :: L)) * (6%nat * / (pPred (vNum b) - 1%nat)))%R
with
((1 + (length L1 + length (a :: L))) *
(6%nat * / (pPred (vNum b) - 1%nat)))%R; [ idtac | ring ].
2: replace (length L1 + length (a :: L))%R with (INR (length (L1 ++ a :: L))).
2: cut (0 <= 6%nat * / (pPred (vNum b) - 1%nat))%R; [ intros T1 | idtac ].
2: apply
Rle_trans
with ((1 + length input) * (6%nat * / (pPred (vNum b) - 1%nat)))%R;
[ apply Rmult_le_compat_r; auto with real | idtac ].
2: apply Rplus_le_compat_l; apply Rle_INR.
2: apply endof_length; auto.
2: apply
Rle_trans
with (1 * / (pPred (vNum b) - 1%nat) * (6%nat * length input + 6%nat))%R;
[ right; ring; ring | idtac ].
2: apply
Rle_trans
with
(/ (pPred (vNum b) - 1%nat - 6%nat * length input) *
(1 - 6%nat * length (a :: L) * / (pPred (vNum b) - 1%nat)) *
(6%nat * length input + 6%nat))%R; [ idtac | right; ring ].
2: apply Rmult_le_compat_r; auto with real.
2: apply Rle_trans with (0 + 6%nat)%R;
[ apply Rlt_le | apply Rplus_le_compat_r ]; auto with real arith.
2: replace (0 + 6%nat)%R with (INR 6); auto with arith real.
2: replace 0%R with (0 * INR 0)%R;
[ auto with real arith | simpl in |- *; ring ].
2: replace (1 - 6%nat * length (a :: L) * / (pPred (vNum b) - 1%nat))%R with
((pPred (vNum b) - 1%nat - 6%nat * length (a :: L)) *
/ (pPred (vNum b) - 1%nat))%R.
2: rewrite <- Rmult_assoc; apply Rmult_le_compat_r.
2: apply Rlt_le; apply Rinv_0_lt_compat; auto.
2: simpl in |- *; replace 1%R with (IZR 1);
[ auto with real zarith | simpl in |- *; auto ].
2: apply
Rle_trans
with
(/ (pPred (vNum b) - 1%nat - 6%nat * length input) *
(pPred (vNum b) - 1%nat - 6%nat * length input))%R.
2: rewrite Rmult_comm; rewrite Rinv_r; auto with real.
2: apply Rmult_le_compat_l.
2: apply Rlt_le; apply Rinv_0_lt_compat; auto.
2: unfold Rminus in |- *; apply Rplus_le_compat_l; apply Ropp_le_contravar.
2: apply Rmult_le_compat_l; auto with real arith.
2: apply Rle_INR; apply endof_length; auto.
2: replace 1%R with ((pPred (vNum b) - 1%nat) * / (pPred (vNum b) - 1%nat))%R;
[ ring | idtac ].
2: apply Rinv_r; auto with real arith.
2: apply not_eq_sym; apply Rlt_dichotomy_converse; left; auto.
2: simpl in |- *; replace 1%R with (IZR 1);
[ auto with real zarith | simpl in |- *; auto ].
2: apply Rle_mult_pos; auto with real zarith arith.
2: rewrite app_length; rewrite plus_INR; auto.
replace (1 - 6%nat * length (a :: L) * / (pPred (vNum b) - 1%nat))%R with
((pPred (vNum b) - 1%nat - 6%nat * length (a :: L)) *
/ (pPred (vNum b) - 1%nat))%R.
apply Rmult_lt_0_compat; auto with real zarith arith.
2: replace 1%R with ((pPred (vNum b) - 1%nat) * / (pPred (vNum b) - 1%nat))%R;
[ ring | idtac ].
2: apply Rinv_r; auto with real zarith.
2: apply not_eq_sym; apply Rlt_dichotomy_converse; left; auto.
2: simpl in |- *; replace 1%R with (IZR 1);
[ auto with real zarith | simpl in |- *; auto ].
apply Rlt_le_trans with (pPred (vNum b) - 1%nat - 6%nat * length input)%R.
2: unfold Rminus in |- *; apply Rplus_le_compat_l; apply Ropp_le_contravar;
apply Rmult_le_compat_l; auto with real arith.
2: apply Rle_INR; apply endof_length; auto.
apply Rlt_Rminus_ZERO.
apply Rmult_lt_reg_l with (/ (pPred (vNum b) - 1%nat))%R;
auto with real zarith.
rewrite Rinv_l; try rewrite Rmult_comm; auto with real zarith.
apply not_eq_sym; apply Rlt_dichotomy_converse; left; auto with real zarith.
Qed.
Theorem FexpAdd_aux2 :
forall L : list float,
L = input ->
IsExp b input ->
(6%nat * length input * / (pPred (vNum b) - 1%nat) < 1)%R ->
exists output : list float,
(input <> nil -> length output <= S (length input)) /\
sum input = sum output :>R /\ IsRleEpsExp output.
cut (1 < pPred (vNum b))%Z;
[ intros Z1
| apply Zlt_trans with radix; auto with zarith;
apply (pPredMoreThanRadix b radix precision); auto with zarith ].
cut (0 < pPred (vNum b) - 1)%R;
[ intros Z2 | replace 1%R with (IZR 1); auto with real zarith ].
intros L; case L.
intros H H1 H2; rewrite <- H; exists (nil (A:=float)); repeat (split; auto);
apply IsRleEpsExpNil.
intros f l H H0 H1; rewrite <- H in H0.
case (FexpAdd_aux H1 l nil f (Fzero (Fexp f))); auto.
rewrite <- H; exists (f :: nil); auto.
inversion H0; auto.
apply IsExpNil.
apply IsRleEpsExpNil.
inversion H0; auto.
apply FboundedZeroSameExp; auto.
inversion H0; auto.
inversion H0; simpl in |- *; auto.
case H3; simpl in |- *; auto.
simpl in |- *; auto with zarith.
left; right; apply FzeroisReallyZero.
intros x' (y', (output', (E1, (E2, (E3, (E4, (E5, (E6, (E7, (E8, E9)))))))))).
cut (TotalP (Closest b radix));
[ intros CT | apply ClosestTotal with (precision := precision); auto ].
case (CT (x' + y')%R); intros p Ep.
case (errorBoundedPlus b radix precision) with (6 := Ep); auto with zarith;
intros q (Eq1, (Eq2, Eq3)).
exists (q :: p :: output'); repeat split.
intros H2; rewrite <- H; simpl in |- *.
apply le_trans with (S (S (length l + 0))); auto with arith zarith.
rewrite <- H; simpl in |- *; fold radix in |- *; rewrite Eq1.
replace
(FtoR radix x' + FtoR radix y' - FtoR radix p + (FtoR radix p + sum output'))%R
with (FtoRradix x' + (FtoRradix y' + sum output'))%R;
[ rewrite <- E3; simpl in |- * | fold FtoRradix in |- *; ring ].
rewrite (FzeroisReallyZero radix); fold FtoRradix; ring.
cut (Fbounded b p);
[ intros Bp
| apply (RoundedModeBounded b radix) with (P := Closest b radix) (2 := Ep);
apply ClosestRoundedModeP with precision; auto with float ].
apply IsRleEpsExpTop; auto.
case (Req_dec p 0); intros Hp0.
cut (FtoRradix q = 0%R); [ intros Hq0 | idtac ].
rewrite Hp0; rewrite Hq0; rewrite Rabs_R0; right; ring.
unfold FtoRradix in |- *; rewrite Eq1; apply TwoSumNul with b precision; auto.
apply Rle_trans with (Rabs p * / 2%nat * (radix * / Zpos (vNum b)))%R;
auto.
unfold FtoRradix in |- *; rewrite Eq1; rewrite <- Rabs_Ropp.
replace (- (FtoR radix x' + FtoR radix y' - FtoR radix p))%R with
(FtoR radix p - (FtoR radix x' + FtoR radix y'))%R;
[ idtac | ring ].
apply plusErrorBound1bis with precision; auto with real zarith.
Contradict Hp0; unfold FtoRradix in |- *; apply is_Fzero_rep1; auto.
rewrite Rmult_assoc; rewrite Rmult_comm; apply Rmult_le_compat_r;
auto with real.
rewrite <- Rmult_assoc.
replace (/ 2%nat * radix)%R with 1%R;
[ idtac | rewrite Rinv_l; auto with real zarith ].
apply
Rle_trans with ((6%nat * length input + 6%nat) * / Zpos (vNum b))%R;
[ apply Rmult_le_compat_r | apply Rmult_le_compat_l ];
auto with real arith zarith.
replace 1%R with (INR 1); repeat rewrite <- mult_INR;
repeat rewrite <- plus_INR; simpl in |- *; try rewrite <- plus_n_Sm;
auto with real zarith.
replace 1%R with (INR 1); repeat rewrite <- mult_INR;
repeat rewrite <- plus_INR; simpl in |- *; auto with real arith.
apply Rle_Rinv.
apply Rlt_Rminus_ZERO.
apply Rmult_lt_reg_l with (/ (pPred (vNum b) - 1%nat))%R;
auto with real zarith.
rewrite Rmult_comm; rewrite Rinv_l; auto with real.
apply Rle_trans with (Zpos (vNum b) - 6%nat * length input)%R;
auto with real arith zarith.
unfold Rminus in |- *; apply Rplus_le_compat_r; auto with real arith zarith.
replace (- 1%nat)%R with (IZR (- (1)));
[ rewrite <- plus_IZR; unfold pPred, Zpred in |- *; auto with real zarith
| simpl in |- *; auto ].
pattern (IZR (Zpos (vNum b))) at 2 in |- *;
replace (IZR (Zpos (vNum b))) with (Zpos (vNum b) - 0)%R;
unfold Rminus in |- *; try ring; apply Rplus_le_compat_l;
apply Ropp_le_contravar; auto with real arith zarith.
apply Rle_mult_pos; replace 0%R with (INR 0); auto with real arith zarith.
CaseEq output'.
intros; apply IsRleEpsExpSingle; auto.
intros f0 l0 Ho; apply IsRleEpsExpTop; auto.
3: rewrite <- Ho; auto.
rewrite Ho in E6; inversion E6; auto.
case E9.
rewrite Ho; intros; discriminate.
rewrite Ho; intros E10.
apply
Rle_trans
with (/ (pPred (vNum b) - 1%nat) * pPred (vNum b) * Rabs (x' + y'))%R.
apply
Rle_trans
with
(Rabs (FtoR 2%nat x' + FtoR 2%nat y') *
(2%nat * pPred (vNum b) * / (2%nat * pPred (vNum b) - radix)))%R.
unfold FtoRradix in |- *; apply RoundBound with precision; auto.
replace radix with (2%nat * 1)%Z; auto with zarith.
apply Zmult_gt_0_lt_compat_l; auto with zarith.
red in |- *; simpl in |- *; auto.
right.
replace (2%nat * pPred (vNum b) - radix)%R with
(2%nat * (pPred (vNum b) - 1%nat))%R; [ idtac | simpl in |- *; ring ].
rewrite Rinv_mult_distr; auto with real zarith arith.
pattern (IZR (pPred (vNum b))) at 4 in |- *;
replace (IZR (pPred (vNum b))) with
(2%nat * / 2%nat * IZR (pPred (vNum b)))%R.
simpl in |- *; fold radix FtoRradix in |- *; ring.
rewrite Rinv_r; auto with real zarith.
apply Rmult_le_reg_l with (/ pPred (vNum b))%R; auto with real zarith.
apply Rmult_le_reg_l with (pPred (vNum b) - 1%nat)%R; auto with real zarith.
apply
Rle_trans
with
(Rabs (x' + y') *
((pPred (vNum b) - 1%nat) * / (pPred (vNum b) - 1%nat) *
(pPred (vNum b) * / pPred (vNum b))))%R; [ right; ring; ring | idtac ].
repeat rewrite Rinv_r; auto with real zarith.
apply Rle_trans with (Rabs (x' + y')); [ right; ring | idtac ].
apply
Rle_trans
with
((pPred (vNum b) - 1%nat) * / pPred (vNum b) *
((6%nat * length input + 6%nat) *
/ (pPred (vNum b) - 1%nat - 6%nat * length input)) *
Rabs (hd b (f0 :: l0)))%R; [ auto | simpl in |- *; right; ring; ring ].
Qed.
Theorem FexpAdd_main :
IsExp b input ->
(6%nat * length input * / (pPred (vNum b) - 1%nat) < 1)%R ->
exists output : list float,
(input <> nil -> length output <= S (length input)) /\
sum input = sum output :>R /\ IsRleEpsExp output.
intros H H0; apply (FexpAdd_aux2 input); auto.
Qed.
End FexpAdd.
Section FexpAdd.
Variable b : Fbound.
Variable precision : nat.
Let radix := 2%Z.
Let TMTO : (1 < radix)%Z := TwoMoreThanOne.
Hint Resolve TMTO: zarith.
Let FtoRradix := FtoR radix.
Coercion FtoRradix : float >-> R.
Hypothesis precisionGreaterThanOne : 1 < precision.
Hypothesis pGivesBound : Zpos (vNum b) = Zpower_nat radix precision.
Hypothesis Ngd : (1 <= pPred (vNum b) * (1 - / radix))%R.
Hypothesis Ngd2 : (6%nat <= pPred (vNum b) * (1 - / radix * / radix))%R.
Inductive NearEqual : list float -> list float -> Prop :=
| IsEqual : forall x : list float, NearEqual x x
| OneMoreR :
forall (x : list float) (e : float),
Fbounded b e -> NearEqual x (e :: x).
Hint Resolve IsEqual OneMoreR.
Definition Step (x y i x' y' : float) (input output output' : list float) :=
Fbounded b x' /\
Fbounded b y' /\
NearEqual output output' /\
(x + (y + (sum (i :: input) + sum output)))%R =
(x' + (y' + (sum input + sum output')))%R :>R /\
(Fexp i <= Fexp y')%Z /\
(Fexp y' <= Fexp x')%Z /\
IsExp b (i :: input) /\
(x' = 0%R :>R \/
y' = 0%R :>R \/
(Rabs y' <=
pPred (vNum b) * (Float 1%nat (Fexp x') - Float 1%nat (hdexp b input)))%R /\
(Rabs y' <=
pPred (vNum b) * (Float 1%nat (Fexp x') - Float 1%nat (Fexp i)))%R) /\
(output' = output \/
(Rabs (x' + y') <= 3%nat * radix * / pPred (vNum b) * Rabs (hd b output'))%R) /\
(output' = nil \/
(Rabs (hd b input) <=
3%nat * (radix * (Rabs (hd b output') * / (pPred (vNum b) - 1%nat))))%R /\
(input = nil \/
(Float (pPred (vNum b)) (Fexp (hd b input)) <=
3%nat * (radix * (Rabs (hd b output') * / (pPred (vNum b) - 1%nat))))%R /\
(Rabs (sum input) <=
length input *
(3%nat * (radix * (Rabs (hd b output') * / (pPred (vNum b) - 1%nat)))))%R)).
Theorem Rle_mult_pos :
forall r1 r2 : R, (0 <= r1)%R -> (0 <= r2)%R -> (0 <= r1 * r2)%R.
intros r1 r2 H H0; replace 0%R with (0 * 0)%R; auto with real.
Qed.
Theorem AddStep :
forall (x y i : float) (input output : list float),
Fbounded b x ->
Fbounded b y ->
IsExp b (i :: input) ->
(Fexp i <= Fexp y)%Z ->
(Fexp y <= Fexp x)%Z ->
(FtoR radix x = 0%R :>R \/ FtoR radix y = 0%R :>R \/ FtoR radix i = 0%R :>R) \/
(Rabs (FtoR radix y) <=
pPred (vNum b) *
(FtoR radix (Float 1%nat (Fexp x)) - FtoR radix (Float 1%nat (Fexp i))))%R ->
output = nil \/
(Rabs i <=
3%nat * (radix * (Rabs (hd b output) * / (pPred (vNum b) - 1%nat))))%R /\
(Float (pPred (vNum b)) (Fexp i) <=
3%nat * (radix * (Rabs (hd b output) * / (pPred (vNum b) - 1%nat))))%R /\
(Rabs (sum (i :: input)) <=
length (i :: input) *
(3%nat * (radix * (Rabs (hd b output) * / (pPred (vNum b) - 1%nat)))))%R ->
exists x' : float,
(exists y' : float,
(exists output' : list float, Step x y i x' y' input output output')).
cut (1 < pPred (vNum b))%Z;
[ intros Z1
| apply Zlt_trans with radix; auto with zarith;
apply (pPredMoreThanRadix b radix precision); auto with zarith ].
cut (0 < pPred (vNum b) - 1)%R;
[ intros Z2 | replace 1%R with (IZR 1); auto with real zarith ].
cut (0 < pPred (vNum b))%Z; [ intros Z3 | auto with real zarith ].
intros x y i input output H H0 H1 H2 H3 H4 G.
cut (Fbounded b i); [ intros K | inversion H1; auto ].
cut (TotalP (Closest b radix));
[ intros ExC | apply (ClosestTotal b radix precision); auto ].
cut (RoundedModeP b radix (Closest b radix));
[ intros Rc | apply (ClosestRoundedModeP b radix precision); auto ].
case (ExC (y + i)%R); intros u Hu.
cut (Fbounded b u);
[ intros Bu
| apply RoundedModeBounded with radix (Closest b radix) (y + i)%R; auto ].
case (ExC (x + u)%R); intros p' Hp'.
cut (Fbounded b p');
[ intros Bp'
| apply RoundedModeBounded with radix (Closest b radix) (x + u)%R; auto ].
case (ExC (x + u - p' + (y + i - u))%R); intros q' Hq'.
cut (Fbounded b q');
[ intros Bq'
| apply
RoundedModeBounded
with radix (Closest b radix) (x + u - p' + (y + i - u))%R;
auto ].
case
(errorBoundedPlus b radix precision TMTO precisionGreaterThanOne pGivesBound
y i u); auto.
intros v (H'1, (H'2, H'3)).
case
(errorBoundedPlus b radix precision TMTO precisionGreaterThanOne pGivesBound
x u p'); auto.
intros w (H'4, (H'5, H'6)).
case
(errorBoundedPlus b radix precision TMTO precisionGreaterThanOne pGivesBound
w v q'); auto.
replace (FtoR radix v) with (FtoR radix y + FtoR radix i - FtoR radix u)%R;
replace (FtoR radix w) with (FtoR radix x + FtoR radix u - FtoR radix p')%R;
auto with real.
intros r' (H'7, (H'8, H'9)).
case
(ThreeSumLoop b precision precisionGreaterThanOne pGivesBound x y i u v w p'
q' r'); auto.
fold radix in |- *;
replace (FtoR radix v) with (FtoR radix y + FtoR radix i - FtoR radix u)%R;
replace (FtoR radix w) with (FtoR radix x + FtoR radix u - FtoR radix p')%R;
auto with real.
intros p''
(q'',
(r'',
((H''0, (H''1, H''2)),
(((H''3, (H''4, H''5)), (H''6, ((H''7, H''8), H''9))),
[(Hr', H''10)| (Hr', H''10)])))).
exists p''; exists q''; exists output; repeat (split; auto).
replace (p'' + (q'' + (sum input + sum output)))%R with
(p'' + (q'' + r'') + (sum input + sum output))%R;
[ idtac | ring_simplify ].
replace (p'' + (q'' + r''))%R with (x + (y + i))%R;
replace (sum (i :: input)) with (i + sum input)%R.
ring.
simpl in |- *; auto with real.
replace (FtoRradix r'') with (FtoRradix r');
replace (FtoRradix r') with (w + v - q')%R; auto with real.
replace (FtoRradix q') with (FtoRradix q''); auto with real.
replace (FtoRradix p'') with (FtoRradix p'); auto with real.
replace (FtoRradix w) with (FtoR radix x + FtoR radix u - FtoR radix p')%R.
replace (FtoRradix v) with (y + i - u)%R; auto with real.
unfold FtoRradix in |- *; ring; ring.
simpl in |- *; auto with real.
replace (FtoRradix r'') with (FtoRradix r'); replace (FtoRradix r') with 0%R;
auto with real; ring.
apply Zle_trans with (Fexp r''); auto.
case (Req_dec p'' 0); intros; auto.
case H''10; auto; intros H''11.
right; right; split; auto.
apply
Rle_trans
with
(pPred (vNum b) *
(FtoR radix (Float 1%nat (Fexp p'')) - FtoR radix (Float 1%nat (Fexp i))))%R;
auto.
apply Rmult_le_compat_l; auto with real arith.
fold FtoRradix in |- *;
replace (Float 1%nat (Fexp p'') - Float 1%nat (hdexp b input))%R with
(Float 1%nat (Fexp p'') + - Float 1%nat (hdexp b input))%R;
replace (Float 1%nat (Fexp p'') - Float 1%nat (Fexp i))%R with
(Float 1%nat (Fexp p'') + - Float 1%nat (Fexp i))%R;
try ring.
apply Rplus_le_compat_l.
apply Ropp_le_contravar; unfold FtoRradix in |- *; unfold FtoR in |- *;
simpl in |- *.
repeat rewrite Rmult_1_l.
apply Rle_powerRZ; auto with real zarith.
apply IsExpZle; auto.
case G; clear G; auto.
generalize H1; clear H1; case input.
simpl in |- *; intros H1 (Hc1, (Hc2, Hc3)); right.
split; auto.
unfold FtoRradix, FtoR, Fzero in |- *; simpl in |- *; rewrite Rmult_0_l;
rewrite Rabs_R0; repeat apply Rle_mult_pos; auto with real zarith.
intros f l H1 (Hc1, (Hc2, Hc3)); right.
cut (Fbounded b f); [ intros Bf | idtac ].
cut
(Float (pPred (vNum b)) (Fexp f) <=
3%nat *
(radix * (Rabs (FtoRradix (hd b output)) * / (pPred (vNum b) - 1%nat))))%R;
[ intros HZ1 | idtac ].
cut
(FtoRradix (Float (pPred (vNum b)) (Fexp f)) <=
FtoRradix (Float (pPred (vNum b)) (Fexp i)))%R; [ intros HZ2 | idtac ].
split.
apply Rle_trans with (2 := HZ1); auto with real zarith.
simpl in |- *; (rewrite <- (Fabs_correct radix); auto with real zarith);
apply (maxMax1 radix); auto with real zarith.
right; split; auto.
apply sum_IsExp with b precision; auto.
inversion H1; auto with real zarith.
unfold FtoRradix, FtoR in |- *; simpl in |- *; apply Rmult_le_compat_l;
auto with real zarith.
apply Rle_powerRZ; auto with real zarith.
inversion H1; auto with real zarith.
apply Rle_trans with (2 := Hc2); auto with real zarith.
unfold FtoRradix, FtoR in |- *; simpl in |- *; apply Rmult_le_compat_l;
auto with real zarith.
apply Rle_powerRZ; auto with real zarith.
inversion H1; auto with real zarith.
inversion H1; auto.
exists q''; exists r''; exists (p'' :: output); repeat (split; auto).
simpl in |- *.
replace
(FtoRradix q'' + (FtoRradix r'' + (sum input + (FtoR 2 p'' + sum output))))%R
with (p'' + (q'' + r'') + (sum input + sum output))%R;
[ idtac | ring_simplify ].
replace (p'' + (q'' + r''))%R with (x + (y + i))%R; [ fold radix; fold FtoRradix; ring | idtac ].
replace (FtoRradix r'') with (FtoRradix r');
replace (FtoRradix r') with (w + v - q')%R; auto with real.
replace (FtoRradix w) with (FtoR radix x + FtoR radix u - FtoR radix p')%R.
replace (FtoRradix v) with (y + i - u)%R; auto with real.
replace (FtoRradix p'') with (FtoR 2 p');
replace (FtoRradix q'') with (FtoRradix q'); fold radix FtoRradix in |- *;
ring.
simpl in |- *; fold radix FtoRradix in |- *; ring.
right; right; split; auto.
apply Rle_trans with (1 := H''10).
apply Rmult_le_compat_l; auto with real zarith.
unfold Rminus in |- *; apply Rplus_le_compat_l; apply Ropp_le_contravar.
apply (oneExp_le radix); auto with real zarith.
apply IsExpZle; auto.
right; simpl in |- *.
apply Rlt_le; fold radix in |- *; unfold FtoRradix in |- *;
apply
bound3Sum
with
(precision := precision)
(w := w)
(v := v)
(p := x)
(q := y)
(r := i)
(u := u); auto; fold radix FtoRradix in |- *.
apply (ClosestCompatible b radix (x + u)%R (x + u)%R p'); auto with real.
replace (FtoRradix p'') with (FtoR radix p'); auto with real.
apply (ClosestCompatible b radix (w + v)%R (w + v)%R q'); auto with real.
replace (FtoRradix w) with (x + u - p')%R;
replace (FtoRradix v) with (y + i - u)%R; auto with real.
replace (FtoRradix r'') with (FtoRradix r');
replace (FtoRradix q'') with (FtoRradix q'); auto with real.
replace (FtoRradix r'') with (FtoRradix r'); auto with real.
right; simpl in |- *.
cut
(forall z : float,
(Fexp z <= Fexp i)%Z ->
Fbounded b z ->
(Rabs z <= 3%nat * (radix * (Rabs p'' * / (pPred (vNum b) - 1%nat))))%R);
[ intros Hz1 | idtac ].
split; [ apply Hz1; auto with zarith | idtac ].
replace (Fexp (hd b input)) with (hdexp b input).
apply IsExpZle; auto.
case input; simpl in |- *; auto with zarith.
generalize H1; case input; intros; simpl in |- *.
apply FboundedFzero.
inversion H5; auto.
generalize H1; clear H1; case input; auto.
intros f l H1; right.
cut
(Float (pPred (vNum b)) (Fexp f) <=
3%nat * (radix * (Rabs p'' * / (pPred (vNum b) - 1%nat))))%R;
[ intros Hz2 | idtac ].
split.
simpl in |- *; auto.
apply sum_IsExp with b precision; auto.
inversion H1; auto.
rewrite <- (fun x : float => Rabs_pos_eq x); auto with real zarith.
apply Hz1; auto with real zarith.
simpl in |- *; inversion H1; auto.
unfold pPred in |- *; apply maxFbounded; auto.
cut (Fbounded b f); [ intros tmp; case tmp | idtac ]; auto.
cut (IsExp b (f :: l)); [ intros tmp; inversion tmp | inversion H1 ]; auto.
auto with float zarith.
apply (LeFnumZERO radix); auto with real zarith.
intros z H5 H6;
apply Rle_trans with (FtoRradix (Float (pPred (vNum b)) (Fexp z))).
rewrite <- (Fabs_correct radix); auto with real zarith; apply (maxMax1 radix);
auto with real zarith.
replace (FtoRradix (Float (pPred (vNum b)) (Fexp z))) with
(pPred (vNum b) * powerRZ radix (Fexp z))%R;
[ idtac | auto with real zarith ].
replace (3%nat * (radix * (Rabs p'' * / (pPred (vNum b) - 1%nat))))%R with
(pPred (vNum b) *
(3%nat *
(radix *
(radix *
(Rabs p'' * / (pPred (vNum b) * (radix * pPred (vNum b) - radix)))))))%R.
apply Rmult_le_compat_l; auto with real zarith.
apply Rle_trans with (FtoRradix (Float 1%nat (Fexp i))).
unfold FtoRradix in |- *; unfold FtoR in |- *; simpl in |- *;
rewrite Rmult_1_l; auto with real zarith.
apply Rle_powerRZ; auto with real zarith.
apply Rlt_le;
apply
OutSum3
with
(precision := precision)
(p := x)
(q := y)
(r := i)
(u := u)
(v := v)
(w := w)
(p' := p'')
(q' := q'')
(r' := r''); auto with real arith; fold radix FtoRradix in |- *.
apply (ClosestCompatible b radix (x + u)%R (x + u)%R p'); auto with real.
replace (FtoRradix p'') with (FtoRradix p'); auto with real.
apply (ClosestCompatible b radix (w + v)%R (w + v)%R q'); auto with real.
replace (FtoRradix w) with (x + u - p')%R;
replace (FtoRradix v) with (y + i - u)%R; auto with real.
replace (FtoRradix r'') with (FtoRradix r');
replace (FtoRradix q'') with (FtoRradix q'); auto with real.
replace (FtoRradix r'') with (FtoRradix r'); auto with real.
replace (radix * pPred (vNum b) - radix)%R with
(radix * (pPred (vNum b) - 1%nat))%R; [ idtac | simpl; ring ].
repeat rewrite Rinv_mult_distr; auto with real zarith.
pattern (INR 3) at 2 in |- *;
replace (INR 3) with
(radix * / radix * (pPred (vNum b) * / pPred (vNum b) * 3%nat))%R;
[ simpl in |- *; ring | idtac ].
repeat rewrite Rinv_r; auto with real zarith; simpl in |- *; ring.
Qed.
Variable input : list float.
Inductive IsRleEpsExp : list float -> Prop :=
| IsRleEpsExpNil : IsRleEpsExp nil
| IsRleEpsExpSingle :
forall x : float, Fbounded b x -> IsRleEpsExp (x :: nil)
| IsRleEpsExpTop :
forall (x y : float) (L : list float),
Fbounded b x ->
Fbounded b y ->
(Rabs x <=
(6%nat * length input + 6%nat) *
/ (pPred (vNum b) - 1%nat - 6%nat * length input) *
Rabs y)%R -> IsRleEpsExp (y :: L) -> IsRleEpsExp (x :: y :: L).
Theorem endof_Rle_length :
forall (P Q : list float) (k : float),
endof input (P ++ k :: Q) -> (length P <= length input - 1)%R.
intros P Q k H; inversion H.
rewrite H0; repeat rewrite app_length; simpl in |- *.
replace 1%R with (INR 1); [ idtac | simpl in |- *; auto ].
replace (length x + (length P + S (length Q))) with
(S (length P) + (length Q + length x));
[ idtac
| repeat rewrite (fun x : list float => S_to_plus_one (length x)); ring ].
rewrite <- minus_INR; simpl in |- *; try apply Rle_INR; auto with arith.
rewrite <- minus_n_O; auto with arith.
Qed.
Theorem FexpAdd_aux :
(6%nat * length input * / (pPred (vNum b) - 1%nat) < 1)%R ->
forall (L output : list float) (x y : float),
endof input L ->
IsExp b L ->
IsRleEpsExp output ->
output = nil \/
(Rabs (x + y) <=
(pPred (vNum b) - 1%nat) * / pPred (vNum b) *
((6%nat * length input + 6%nat) *
/ (pPred (vNum b) - 1%nat - 6%nat * length input)) *
Rabs (hd b output))%R ->
output = nil \/
(exists L1 : list float,
(exists x1 : float,
(exists y1 : float,
IsExp b (L1 ++ L) /\
endof input (L1 ++ L) /\
(x1 + (y1 + sum L1))%R = (x + y)%R :>R /\
(Rabs (sum L1) <=
length L1 *
(3%nat *
(radix * (Rabs (hd b output) * / (pPred (vNum b) - 1%nat)))))%R /\
(Rabs (x1 + y1) <=
3%nat * radix * / pPred (vNum b) * Rabs (hd b output))%R))) ->
Fbounded b x ->
Fbounded b y ->
(hdexp b L <= Fexp y)%Z ->
(Fexp y <= Fexp x)%Z ->
(FtoR radix x = 0%R :>R \/ FtoR radix y = 0%R :>R) \/
(Rabs (FtoR radix y) <=
pPred (vNum b) *
(FtoR radix (Float 1%nat (Fexp x)) - FtoR radix (Float 1%nat (hdexp b L))))%R ->
output = nil \/
L = nil \/
(Rabs (hd b L) <=
3%nat * (radix * (Rabs (hd b output) * / (pPred (vNum b) - 1%nat))))%R /\
(Float (pPred (vNum b)) (Fexp (hd b L)) <=
3%nat * (radix * (Rabs (hd b output) * / (pPred (vNum b) - 1%nat))))%R /\
(Rabs (sum L) <=
length L *
(3%nat * (radix * (Rabs (hd b output) * / (pPred (vNum b) - 1%nat)))))%R ->
exists x' : float,
(exists y' : float,
(exists output' : list float,
Fbounded b x' /\
Fbounded b y' /\
(x + (y + (sum output + sum L)))%R = (x' + (y' + sum output'))%R :>R /\
length output' <= length L + length output /\
(Fexp y' <= Fexp x')%Z /\
IsRleEpsExp output' /\
(output = nil \/
(exists L1 : list float,
(exists x1 : float,
(exists y1 : float,
IsExp b (L1 ++ L) /\
endof input (L1 ++ L) /\
(x1 + (y1 + sum L1))%R = (x + y)%R :>R /\
(Rabs (sum L1) <=
length L1 *
(3%nat *
(radix *
(Rabs (hd b output) * / (pPred (vNum b) - 1%nat)))))%R /\
(Rabs (x1 + y1) <=
3%nat * radix * / pPred (vNum b) * Rabs (hd b output))%R)))) /\
endof output' output /\
(output' = nil \/
(Rabs (x' + y') <=
(pPred (vNum b) - 1%nat) * / pPred (vNum b) *
((6%nat * length input + 6%nat) *
/ (pPred (vNum b) - 1%nat - 6%nat * length input)) *
Rabs (hd b output'))%R))).
intros Z0.
cut (1 < pPred (vNum b))%Z;
[ intros Z1
| apply Zlt_trans with radix;
try apply (pPredMoreThanRadix b radix precision);
auto with zarith ].
cut (0 < pPred (vNum b))%Z;
[ intros Zl2 | apply Zlt_1_O; apply Zlt_le_weak; auto with zarith ].
cut (0 < pPred (vNum b) - 1%nat - 6%nat * length input)%R;
[ intros Z3
| apply Rlt_Rminus_ZERO;
apply Rmult_lt_reg_l with (/ (pPred (vNum b) - 1%nat))%R;
try rewrite Rinv_l; try rewrite (fun x => Rmult_comm (/ x));
auto with real zarith; replace (INR 1) with (IZR 1);
try rewrite Z_R_minus; auto with real zarith ].
cut
(3%nat * radix * / pPred (vNum b) <=
(pPred (vNum b) - 1%nat) * / pPred (vNum b) *
((6%nat * length input + 6%nat) *
/ (pPred (vNum b) - 1%nat - 6%nat * length input)))%R;
[ intros Z4 | auto with real zarith ].
2: replace (3%nat * radix)%R with (INR 6 * 1)%R;
[ idtac | rewrite Rmult_1_r; simpl in |- *; ring ].
2: repeat rewrite Rmult_assoc.
2: rewrite (Rmult_comm 1).
2: rewrite (Rmult_comm (pPred (vNum b) - 1%nat)).
2: replace (6%nat * length input + 6%nat)%R with
(INR 6 * (length input + 1%nat))%R; [ idtac | simpl in |- *; ring ].
2: repeat rewrite <- (Rmult_assoc (/ pPred (vNum b))).
2: repeat rewrite (Rmult_comm (/ pPred (vNum b)) 6%nat).
2: repeat rewrite Rmult_assoc.
2: repeat apply Rmult_le_compat_l; auto with real zarith.
2: rewrite (Rmult_comm (length input + 1%nat)).
2: apply
Rmult_le_reg_l with (pPred (vNum b) - 1%nat - 6%nat * length input)%R;
auto with real zarith.
2: repeat rewrite <- Rmult_assoc; rewrite Rinv_r; auto with real zarith;
rewrite Rmult_1_r; rewrite Rmult_1_l.
2: rewrite (fun x => Rplus_comm x 1); rewrite Rmult_plus_distr_l;
rewrite Rmult_1_r.
2: unfold Rminus in |- *; apply Rplus_le_compat_l.
2: rewrite <- Ropp_mult_distr_l_reverse; apply Rmult_le_compat_r;
auto with real arith.
2: apply Rle_trans with (-0)%R; auto with real zarith arith.
2: rewrite Ropp_0; change (0 <= pPred (vNum b) - 1%nat)%R in |- *;
auto with real zarith.
intros L; elim L; clear L.
intros output x y H H0 H1 H2 H3 H4 H5 H6 H7 H8 H9; exists x; exists y;
exists output; repeat (split; simpl in |- *; auto);
[ ring | exists (nil (A:=float)) ]; auto.
intros a L HrecL output; case output; clear output.
intros x y H H0 H1 H2 H3 H4 H5 H6 H7 H8 H9.
cut (IsExp b L); [ intros Hp1 | inversion H0; auto; apply IsExpNil ].
cut (endof input L);
[ intros Hp2
| case H; intros LL H10; exists (LL ++ a :: nil); rewrite H10;
replace (a :: L) with ((a :: nil) ++ L); try apply ass_app;
simpl in |- *; auto ].
case (AddStep x y a L nil); auto with real zarith.
inversion H0; auto.
case H8; auto; intros tmp; case tmp; auto.
case H8; auto; intros tmp; case tmp; auto.
intros x' (y', (output', H'1)).
case H'1; intros T1 (T2, (T3, (T4, (T5, (T6, (T7, (T8, (T9, T10)))))))).
case (HrecL output' x' y'); clear H2 H3 H9 HrecL; auto with real zarith.
inversion T3; auto.
apply IsRleEpsExpSingle; auto.
CaseEq output'; auto; intros o output'' Eq1.
right; case T9; rewrite Eq1;
[ intros; discriminate | simpl in |- *; clear T9; intros T9 ].
apply
Rle_trans with (3%nat * radix * / pPred (vNum b) * Rabs (FtoRradix o))%R;
auto with real.
CaseEq output'; auto; intros o output'' Eq1.
right; exists (nil (A:=float)); exists x'; exists y';
repeat (split; simpl in |- *; auto with real).
rewrite Rabs_R0; rewrite Rmult_0_l; auto with real.
case T9; rewrite Eq1; [ intros; discriminate | simpl in |- *; auto ].
apply Zle_trans with (Fexp a); auto with zarith.
apply IsExpZle; auto.
repeat (case T8; auto; clear T8; intros T8).
case T10; auto; clear T10; intros (T10, T11); right; case T11; auto.
intros x'0
(y'0, (output'0, (V1, (V2, (V3, (V4, (V5, (V6, (V7, (V8, V9)))))))))).
exists x'0; exists y'0; exists output'0; repeat (split; auto).
rewrite (Rplus_comm (sum nil)); rewrite T4; rewrite <- V3; ring.
generalize V4; case T3; simpl in |- *; auto with arith; intros e Be V0;
rewrite plus_n_Sm; simpl in |- *; auto.
exists output'0; apply app_nil_end.
intros a0 output x y V' H H0 H1 M H2 H3 H4 H5 H6 H7.
cut (IsExp b L); [ intros Hp1 | inversion H; auto; apply IsExpNil ].
cut (endof input L);
[ intros Hp2
| case V'; intros LL H10; exists (LL ++ a :: nil); rewrite H10;
replace (a :: L) with ((a :: nil) ++ L); try apply ass_app;
simpl in |- *; auto ].
case H1; [ intros; discriminate | clear H1; intros H1 ].
case M; [ intros; discriminate | clear M; intros M ].
case H7; [ intros; discriminate | clear H7; intros H7 ].
case H7; [ intros; discriminate | clear H7; intros (H7, (H8, H9)) ].
case (AddStep x y a L (a0 :: output)); auto with real zarith.
repeat (case H6; auto; clear H6; intros H6).
intros x'
(y', (output', (V1, (V2, (V3, (V4, (V5, (V6, (V7, (V8, (V9, V10))))))))))).
cut (IsRleEpsExp output'); [ intros RR | idtac ].
case (HrecL output' x' y'); clear H2 H3 H9 HrecL; auto with real zarith.
right; auto.
4: repeat (case V8; auto; clear V8; intros V8).
3: apply Zle_trans with (Fexp a); auto with zarith.
3: apply IsExpZle; auto.
3: case V10; auto; try intros (tmp1, [tmp2| tmp3]); auto.
3: intros x'0
(y'0, (output'0, (W1, (W2, (W3, (W4, (W5, (W6, (W7, (W8, W9)))))))))).
3: exists x'0; exists y'0; exists output'0; repeat (split; auto).
3: rewrite <- W3; rewrite (Rplus_comm (sum output')); rewrite <- V4; ring.
3: apply le_trans with (1 := W4).
3: inversion V3; simpl in |- *; auto with arith; repeat rewrite plus_n_Sm;
simpl in |- *; auto with arith.
3: inversion V3; [ replace (a0 :: output) with output'; auto | idtac ].
3: rewrite <- H9 in W8; inversion W8; exists (x1 ++ e :: nil); rewrite H10;
replace (e :: a0 :: output) with ((e :: nil) ++ a0 :: output);
try (apply sym_equal; apply app_ass); simpl in |- *;
auto.
2: case M; intros L1 (x1, (y1, (X1, (X2, (X3, (X4, X5)))))).
2: inversion V3; auto.
2: right; exists (L1 ++ a :: nil); exists x1; exists y1; repeat (split; auto).
2: rewrite app_ass; simpl in |- *; auto.
2: rewrite app_ass; simpl in |- *; auto.
2: apply trans_eq with (FtoRradix x + FtoRradix y + a)%R;
fold FtoRradix in |- *;
[ rewrite <- X3; rewrite <- sum_app; fold radix; fold FtoRradix; ring | idtac ].
2: replace (FtoRradix x + FtoRradix y)%R with
(FtoRradix x + (FtoRradix y + (sum (a :: L) + sum (a0 :: output))) +
- (sum (a :: L) + sum (a0 :: output)))%R;
[ rewrite V4; rewrite H3; simpl in |- *; fold radix FtoRradix in |- *;
ring
| fold FtoRradix in |- *; ring ].
2: rewrite <- sum_app.
2: apply Rle_trans with (1 := Rabs_triang (FtoR 2 a) (sum L1)).
2: rewrite app_length; simpl in |- *.
2: rewrite plus_INR; simpl in |- *; rewrite Rmult_plus_distr_r;
rewrite Rmult_1_l.
2: rewrite (fun x => Rplus_comm (Rabs x)); apply Rplus_le_compat; auto.
2: right; exists (nil (A:=float)); exists x'; exists y';
repeat (split; simpl in |- *; auto).
2: rewrite Rplus_0_r; auto.
2: rewrite Rabs_R0; rewrite Rmult_0_l; auto with real.
2: case V9; rewrite <- H9; simpl in |- *; auto.
2: intros tmp; Contradict tmp; apply cons_neq; auto.
inversion V3.
case M; intros L1 (x1, (y1, (X1, (X2, (X3, (X4, X5)))))).
replace (FtoRradix x' + FtoRradix y')%R with
(a + (FtoRradix x1 + FtoRradix y1 + sum L1))%R.
apply
Rle_trans
with
(Rabs (FtoRradix a) +
(Rabs (FtoRradix x1 + FtoRradix y1) + Rabs (sum L1)))%R.
apply
Rle_trans
with (1 := Rabs_triang (FtoR 2 a) (FtoRradix x1 + FtoRradix y1 + sum L1));
auto with real.
apply Rplus_le_compat_l; apply Rabs_triang; auto with real.
apply
Rle_trans
with
(3%nat *
(2%nat *
(Rabs (FtoRradix (hd b (a0 :: output))) * / (pPred (vNum b) - 1%nat))) +
(Rabs (x1 + y1) + Rabs (sum L1)))%R;
[ apply Rplus_le_compat_r; auto | idtac ].
apply
Rle_trans
with
(3%nat *
(2%nat *
(Rabs (FtoRradix (hd b (a0 :: output))) * / (pPred (vNum b) - 1%nat))) +
(3%nat * 2%nat * / pPred (vNum b) *
Rabs (FtoRradix (hd b (a0 :: output))) + Rabs (sum L1)))%R;
[ apply Rplus_le_compat_l; apply Rplus_le_compat_r; auto | idtac ].
apply
Rle_trans
with
(3%nat *
(2%nat *
(Rabs (FtoRradix (hd b (a0 :: output))) * / (pPred (vNum b) - 1%nat))) +
(3%nat * 2%nat * / pPred (vNum b) *
Rabs (FtoRradix (hd b (a0 :: output))) +
length L1 *
(3%nat *
(2%nat *
(Rabs (FtoRradix (hd b (a0 :: output))) * / (pPred (vNum b) - 1%nat))))))%R;
[ apply Rplus_le_compat_l; apply Rplus_le_compat_l; auto | idtac ].
apply
Rle_trans
with
(3%nat *
(2%nat *
(Rabs (FtoRradix (hd b (a0 :: output))) * / (pPred (vNum b) - 1%nat))) +
(3%nat * 2%nat * / (pPred (vNum b) - 1%nat) *
Rabs (FtoRradix (hd b (a0 :: output))) +
length L1 *
(3%nat *
(2%nat *
(Rabs (FtoRradix (hd b (a0 :: output))) * / (pPred (vNum b) - 1%nat))))))%R;
[ apply Rplus_le_compat_l; apply Rplus_le_compat_r; simpl in |- *;
auto with real
| idtac ].
apply Rmult_le_compat_r; auto with real.
apply Rmult_le_compat_l; auto with real.
replace ((2 + 1) * 2)%R with (INR 6); auto with real arith; simpl in |- *;
ring.
apply Rle_Rinv; auto with real arith.
pattern (IZR (pPred (vNum b))) at 2 in |- *;
replace (IZR (pPred (vNum b))) with (IZR (pPred (vNum b)) + -0)%R;
[ unfold Rminus in |- *; auto with real | ring ].
apply
Rle_trans
with
((3%nat * 2%nat + (3%nat * 2%nat + length L1 * (3%nat * 2%nat))) *
/ (pPred (vNum b) - 1%nat) * Rabs (FtoRradix (hd b (a0 :: output))))%R;
[ right; ring; ring | idtac ].
apply Rmult_le_compat_r; auto with real.
replace (3%nat * 2%nat)%R with (INR 6); [ idtac | simpl in |- *; ring ].
apply
Rle_trans
with
((6%nat + (6%nat + (length input - 1) * 6%nat)) *
/ (pPred (vNum b) - 1%nat))%R.
apply Rmult_le_compat_r; auto with real.
apply Rplus_le_compat_l; apply Rplus_le_compat_l; apply Rmult_le_compat_r;
auto with arith zarith real.
apply endof_Rle_length with L a; auto.
replace (6%nat + (length input - 1) * 6%nat)%R with (6%nat * length input)%R;
[ idtac | ring; ring ].
apply
Rle_trans
with
((6%nat + 6%nat * length input) *
((pPred (vNum b) - 1%nat) * / pPred (vNum b) *
/ (pPred (vNum b) - 1%nat - 6%nat * length input)))%R;
[ idtac | right; ring; ring ].
apply Rmult_le_compat_l; auto with arith zarith real.
repeat rewrite <- plus_INR || rewrite <- mult_INR; apply pos_INR;
simpl in |- *; auto with arith.
apply Rmult_le_reg_l with (pPred (vNum b) - 1%nat)%R; auto with real.
apply Rmult_le_reg_l with (pPred (vNum b) - 1%nat - 6%nat * length input)%R;
auto with real.
apply Rmult_le_reg_l with (IZR (pPred (vNum b))); auto with real.
rewrite Rinv_r;
[ idtac
| simpl in |- *; replace 1%R with (IZR 1); try rewrite Z_R_minus;
auto with real zarith ].
apply
Rle_trans
with
((pPred (vNum b) - 1%nat) *
((pPred (vNum b) - 1%nat) *
(pPred (vNum b) * / pPred (vNum b) *
((pPred (vNum b) - 1%nat - 6%nat * length input) *
/ (pPred (vNum b) - 1%nat - 6%nat * length input)))))%R;
[ idtac | right; ring; ring ].
rewrite Rinv_r; [ idtac | auto with real zarith ].
rewrite Rinv_r; [ idtac | auto with real zarith ].
apply
Rle_trans
with
(pPred (vNum b) * pPred (vNum b) +
- (pPred (vNum b) * (6%nat * length input + 1)))%R;
[ right; simpl; ring | idtac ].
apply
Rle_trans
with (pPred (vNum b) * pPred (vNum b) + (- (pPred (vNum b) * 2) + 1))%R;
[ idtac | right; simpl; ring ].
apply Rplus_le_compat_l.
apply Rle_trans with (- (pPred (vNum b) * 2))%R.
apply Ropp_le_contravar; apply Rmult_le_compat_l; auto with real zarith.
apply Rle_trans with (INR 6 * INR 1)%R.
replace 2%R with (INR 2); repeat rewrite <- plus_INR || rewrite <- mult_INR;
try apply Rle_INR; simpl in |- *; auto with real arith.
apply Rle_trans with (6%nat * length input)%R; auto with real zarith.
apply Rmult_le_compat_l; auto with real arith.
case X2; intros LL HL; rewrite HL.
repeat rewrite app_length; simpl in |- *; repeat rewrite <- plus_n_Sm;
replace 1%R with (INR 1); auto with real arith.
rewrite <- (Rplus_0_r (- (pPred (vNum b) * 2))); auto with real.
repeat rewrite Rplus_assoc; rewrite X3;
replace (FtoRradix x' + FtoRradix y')%R with
(FtoRradix x' + (FtoRradix y' + (sum L + sum output')) +
- (sum L + sum output'))%R;
[ rewrite <- V4; rewrite <- H3; simpl in |- *; fold FtoRradix in |- *;
fold radix; fold FtoRradix; ring | ring ].
case V9;
[ rewrite <- H9; intros tmp; Contradict tmp; apply cons_neq
| rewrite H9; intros H'9 ].
apply Rle_trans with (1 := H'9); auto with real.
inversion V3; auto.
apply IsRleEpsExpTop; auto.
inversion H0; auto.
apply
Rmult_le_reg_l
with (1 - 6%nat * length (a :: L) * / (pPred (vNum b) - 1%nat))%R.
2: apply
Rle_trans
with
((1 -
(6%nat * / pPred (vNum b) +
6%nat * length L * / (pPred (vNum b) - 1%nat))) *
Rabs e)%R.
2: apply Rmult_le_compat_r; auto with real.
2: unfold Rminus in |- *; apply Rplus_le_compat_l; apply Ropp_le_contravar.
2: replace (6%nat * length (a :: L) * / (pPred (vNum b) + - 1%nat))%R with
(6%nat * / (pPred (vNum b) + -1) +
6%nat * length L * / (pPred (vNum b) + -1))%R;
[ apply Rplus_le_compat_r | idtac ].
2: apply Rmult_le_compat_l; apply Rlt_le; auto with real arith.
2: apply Rinv_1_lt_contravar; auto with real arith.
2: replace 1%R with (IZR 1);
[ rewrite <- Ropp_Ropp_IZR; rewrite <- plus_IZR; apply Rle_IZR
| simpl in |- *; auto ].
2: change (1 <= Zpred (pPred (vNum b)))%Z in |- *; auto with zarith.
2: pattern (IZR (pPred (vNum b))) at 2 in |- *;
replace (IZR (pPred (vNum b))) with (pPred (vNum b) + -0)%R;
auto with real; ring.
2: replace (INR (length (a :: L))) with (1 + INR (length L))%R;
[ simpl in |- *; ring
| simpl in |- *; case (length L); simpl in |- *; intros; ring ].
2: case V9;
[ rewrite <- H12; intros tmp; Contradict tmp; apply cons_neq
| rewrite <- H12; intros H'9 ].
2: replace
((1 -
(6%nat * / pPred (vNum b) +
6%nat * length L * / (pPred (vNum b) - 1%nat))) *
Rabs e)%R with
(Rabs e -
(6%nat * / pPred (vNum b) * Rabs e +
6%nat * length L * / (pPred (vNum b) - 1%nat) * Rabs e))%R;
[ idtac | ring; ring ].
2: apply
Rle_trans
with
(Rabs e -
(Rabs (x' + y') +
6%nat * length L * / (pPred (vNum b) - 1%nat) * Rabs e))%R.
2: unfold Rminus in |- *; apply Rplus_le_compat_l.
2: apply Ropp_le_contravar; apply Rplus_le_compat_r.
2: replace (INR 6) with (3%nat * 2%nat)%R; replace e with (hd b (e :: x0));
auto with real.
2: rewrite <- (mult_INR 3 2); auto with arith.
2: case V10;
[ rewrite <- H12; intros; discriminate
| rewrite <- H12; intros (H'10, H'11) ].
2: apply Rle_trans with (Rabs e - (Rabs (x' + y') + Rabs (sum L)))%R.
2: unfold Rminus in |- *; apply Rplus_le_compat_l; apply Ropp_le_contravar;
apply Rplus_le_compat_l.
2: case H'11;
[ intros H'12; rewrite H'12; simpl in |- *; rewrite Rabs_R0; right; ring
| intros (H'12, H'13); auto ].
2: apply Rle_trans with (1 := H'13); right; simpl in |- *.
2: ring_simplify (pPred (vNum b) +- 1)%R;ring.
elim H2; intros; clear H2.
2: apply Rle_trans with (Rabs (e + (x' + y' + sum L))).
2: rewrite <- (Rabs_Ropp (x' + y')); rewrite <- (Rabs_Ropp (sum L)).
2: replace (Rabs e - (Rabs (- (x' + y')) + Rabs (- sum L)))%R with
(Rabs e - Rabs (- (x' + y')) - Rabs (- sum L))%R;
[ idtac | ring; ring ].
2: apply Rle_trans with (Rabs (e - - (x' + y')) - Rabs (- sum L))%R.
2: unfold Rminus in |- *; apply Rplus_le_compat_r.
2: replace (Rabs e + - Rabs (- (x' + y')))%R with
(Rabs e - Rabs (- (x' + y')))%R;
replace (e + - - (x' + y'))%R with (e - - (x' + y'))%R;
try ring.
2: apply Rabs_triang_inv.
2: replace (e + (x' + y' + sum L))%R with (e - - (x' + y') - - sum L)%R;
try ring.
2: apply Rabs_triang_inv.
2: replace (e + (x' + y' + sum L))%R with (x + (y + sum (a :: L)))%R.
3: apply Rplus_eq_reg_l with (sum (a0 :: output)).
3: replace (sum (a0 :: output) + (x + (y + sum (a :: L))))%R with
(FtoRradix x + (FtoRradix y + (sum (a :: L) + sum (a0 :: output))))%R;
[ rewrite V4; rewrite <- H12; simpl in |- *; fold radix FtoRradix in |- *;
ring
| ring ].
2: case M; intros L1 (x1, (y1, (X1, (X2, (X3, (X4, X5)))))).
2: replace (x + (y + sum (a :: L)))%R with (x + y + sum (a :: L))%R;
[ rewrite <- X3 | ring ].
2: replace (x1 + (y1 + sum L1) + sum (a :: L))%R with
(x1 + y1 + (sum L1 + sum (a :: L)))%R; [ idtac | ring ].
2: apply Rle_trans with (Rabs (x1 + y1) + Rabs (sum L1 + sum (a :: L)))%R;
[ apply Rabs_triang | idtac ].
2: apply
Rle_trans
with
(6%nat * / pPred (vNum b) * Rabs a0 + Rabs (sum L1 + sum (a :: L)))%R.
2: apply Rplus_le_compat_r; apply Rle_trans with (1 := X5).
2: simpl in |- *; repeat apply Rmult_le_compat_r; auto with real.
2: replace ((2 + 1) * 2)%R with (INR 6);
replace (2 + 1 + 1 + 1 + 1)%R with (INR 6); auto with real arith;
simpl in |- *; ring.
2: apply
Rle_trans
with
(6%nat * / pPred (vNum b) * Rabs (FtoRradix a0) +
(Rabs (sum L1) + Rabs (sum (a :: L))))%R;
[ apply Rplus_le_compat_l; apply Rabs_triang | idtac ].
2: apply
Rle_trans
with
(6%nat * / pPred (vNum b) * Rabs (FtoRradix a0) +
(length L1 *
(3%nat *
(2%nat * (Rabs (hd b (a0 :: output)) * / (pPred (vNum b) - 1%nat)))) +
Rabs (sum (a :: L))))%R.
2: apply Rplus_le_compat_l; apply Rplus_le_compat_r; auto.
2: apply
Rle_trans
with
(6%nat * / pPred (vNum b) * Rabs (FtoRradix a0) +
(length L1 *
(3%nat *
(2%nat * (Rabs (hd b (a0 :: output)) * / (pPred (vNum b) - 1%nat)))) +
length (a :: L) *
(3%nat *
(2%nat * (Rabs (hd b (a0 :: output)) * / (pPred (vNum b) - 1%nat))))))%R.
2: apply Rplus_le_compat_l; apply Rplus_le_compat_l; auto.
2: replace (hd b (a0 :: output)) with a0; [ idtac | simpl in |- *; auto ].
2: replace
(length L1 * (3%nat * (2%nat * (Rabs a0 * / (pPred (vNum b) - 1%nat)))) +
length (a :: L) *
(3%nat * (2%nat * (Rabs a0 * / (pPred (vNum b) - 1%nat)))))%R with
((length L1 + length (a :: L)) *
(3%nat * 2%nat * (Rabs a0 * / (pPred (vNum b) - 1%nat))))%R;
[ idtac | ring; ring ].
2: replace (3%nat * 2%nat)%R with (INR 6); [ idtac | simpl in |- *; ring ].
2: replace
(6%nat * / pPred (vNum b) * Rabs a0 +
(length L1 + length (a :: L)) *
(6%nat * (Rabs a0 * / (pPred (vNum b) - 1%nat))))%R with
((6%nat * / pPred (vNum b) +
(length L1 + length (a :: L)) * (6%nat * / (pPred (vNum b) - 1%nat))) *
Rabs a0)%R; [ idtac | ring ].
2: apply
Rle_trans
with
((1 - 6%nat * length (a :: L) * / (pPred (vNum b) - 1%nat)) *
((6%nat * length input + 6%nat) *
/ (pPred (vNum b) - 1%nat - 6%nat * length input)) *
Rabs a0)%R; [ idtac | right; ring ].
2: apply Rmult_le_compat_r; auto with real.
2: apply
Rle_trans
with
(6%nat * / (pPred (vNum b) - 1%nat) +
(length L1 + length (a :: L)) * (6%nat * / (pPred (vNum b) - 1%nat)))%R.
2: apply Rplus_le_compat_r; apply Rmult_le_compat_l; auto with real arith.
2: apply Rlt_le; apply Rinv_1_lt_contravar; auto with real zarith.
2: simpl in |- *; replace 1%R with (IZR 1);
[ unfold Rminus in |- *; rewrite <- Ropp_Ropp_IZR; rewrite <- plus_IZR;
apply Rle_IZR
| simpl in |- *; auto ].
2: change (1 <= Zpred (pPred (vNum b)))%Z in |- *; auto with zarith.
2: simpl in |- *; unfold Rminus in |- *;
pattern (IZR (pPred (vNum b))) at 2 in |- *;
replace (IZR (pPred (vNum b))) with (pPred (vNum b) + -0)%R;
auto with real; ring.
2: replace
(6%nat * / (pPred (vNum b) - 1%nat) +
(length L1 + length (a :: L)) * (6%nat * / (pPred (vNum b) - 1%nat)))%R
with
((1 + (length L1 + length (a :: L))) *
(6%nat * / (pPred (vNum b) - 1%nat)))%R; [ idtac | ring ].
2: replace (length L1 + length (a :: L))%R with (INR (length (L1 ++ a :: L))).
2: cut (0 <= 6%nat * / (pPred (vNum b) - 1%nat))%R; [ intros T1 | idtac ].
2: apply
Rle_trans
with ((1 + length input) * (6%nat * / (pPred (vNum b) - 1%nat)))%R;
[ apply Rmult_le_compat_r; auto with real | idtac ].
2: apply Rplus_le_compat_l; apply Rle_INR.
2: apply endof_length; auto.
2: apply
Rle_trans
with (1 * / (pPred (vNum b) - 1%nat) * (6%nat * length input + 6%nat))%R;
[ right; ring; ring | idtac ].
2: apply
Rle_trans
with
(/ (pPred (vNum b) - 1%nat - 6%nat * length input) *
(1 - 6%nat * length (a :: L) * / (pPred (vNum b) - 1%nat)) *
(6%nat * length input + 6%nat))%R; [ idtac | right; ring ].
2: apply Rmult_le_compat_r; auto with real.
2: apply Rle_trans with (0 + 6%nat)%R;
[ apply Rlt_le | apply Rplus_le_compat_r ]; auto with real arith.
2: replace (0 + 6%nat)%R with (INR 6); auto with arith real.
2: replace 0%R with (0 * INR 0)%R;
[ auto with real arith | simpl in |- *; ring ].
2: replace (1 - 6%nat * length (a :: L) * / (pPred (vNum b) - 1%nat))%R with
((pPred (vNum b) - 1%nat - 6%nat * length (a :: L)) *
/ (pPred (vNum b) - 1%nat))%R.
2: rewrite <- Rmult_assoc; apply Rmult_le_compat_r.
2: apply Rlt_le; apply Rinv_0_lt_compat; auto.
2: simpl in |- *; replace 1%R with (IZR 1);
[ auto with real zarith | simpl in |- *; auto ].
2: apply
Rle_trans
with
(/ (pPred (vNum b) - 1%nat - 6%nat * length input) *
(pPred (vNum b) - 1%nat - 6%nat * length input))%R.
2: rewrite Rmult_comm; rewrite Rinv_r; auto with real.
2: apply Rmult_le_compat_l.
2: apply Rlt_le; apply Rinv_0_lt_compat; auto.
2: unfold Rminus in |- *; apply Rplus_le_compat_l; apply Ropp_le_contravar.
2: apply Rmult_le_compat_l; auto with real arith.
2: apply Rle_INR; apply endof_length; auto.
2: replace 1%R with ((pPred (vNum b) - 1%nat) * / (pPred (vNum b) - 1%nat))%R;
[ ring | idtac ].
2: apply Rinv_r; auto with real arith.
2: apply not_eq_sym; apply Rlt_dichotomy_converse; left; auto.
2: simpl in |- *; replace 1%R with (IZR 1);
[ auto with real zarith | simpl in |- *; auto ].
2: apply Rle_mult_pos; auto with real zarith arith.
2: rewrite app_length; rewrite plus_INR; auto.
replace (1 - 6%nat * length (a :: L) * / (pPred (vNum b) - 1%nat))%R with
((pPred (vNum b) - 1%nat - 6%nat * length (a :: L)) *
/ (pPred (vNum b) - 1%nat))%R.
apply Rmult_lt_0_compat; auto with real zarith arith.
2: replace 1%R with ((pPred (vNum b) - 1%nat) * / (pPred (vNum b) - 1%nat))%R;
[ ring | idtac ].
2: apply Rinv_r; auto with real zarith.
2: apply not_eq_sym; apply Rlt_dichotomy_converse; left; auto.
2: simpl in |- *; replace 1%R with (IZR 1);
[ auto with real zarith | simpl in |- *; auto ].
apply Rlt_le_trans with (pPred (vNum b) - 1%nat - 6%nat * length input)%R.
2: unfold Rminus in |- *; apply Rplus_le_compat_l; apply Ropp_le_contravar;
apply Rmult_le_compat_l; auto with real arith.
2: apply Rle_INR; apply endof_length; auto.
apply Rlt_Rminus_ZERO.
apply Rmult_lt_reg_l with (/ (pPred (vNum b) - 1%nat))%R;
auto with real zarith.
rewrite Rinv_l; try rewrite Rmult_comm; auto with real zarith.
apply not_eq_sym; apply Rlt_dichotomy_converse; left; auto with real zarith.
Qed.
Theorem FexpAdd_aux2 :
forall L : list float,
L = input ->
IsExp b input ->
(6%nat * length input * / (pPred (vNum b) - 1%nat) < 1)%R ->
exists output : list float,
(input <> nil -> length output <= S (length input)) /\
sum input = sum output :>R /\ IsRleEpsExp output.
cut (1 < pPred (vNum b))%Z;
[ intros Z1
| apply Zlt_trans with radix; auto with zarith;
apply (pPredMoreThanRadix b radix precision); auto with zarith ].
cut (0 < pPred (vNum b) - 1)%R;
[ intros Z2 | replace 1%R with (IZR 1); auto with real zarith ].
intros L; case L.
intros H H1 H2; rewrite <- H; exists (nil (A:=float)); repeat (split; auto);
apply IsRleEpsExpNil.
intros f l H H0 H1; rewrite <- H in H0.
case (FexpAdd_aux H1 l nil f (Fzero (Fexp f))); auto.
rewrite <- H; exists (f :: nil); auto.
inversion H0; auto.
apply IsExpNil.
apply IsRleEpsExpNil.
inversion H0; auto.
apply FboundedZeroSameExp; auto.
inversion H0; auto.
inversion H0; simpl in |- *; auto.
case H3; simpl in |- *; auto.
simpl in |- *; auto with zarith.
left; right; apply FzeroisReallyZero.
intros x' (y', (output', (E1, (E2, (E3, (E4, (E5, (E6, (E7, (E8, E9)))))))))).
cut (TotalP (Closest b radix));
[ intros CT | apply ClosestTotal with (precision := precision); auto ].
case (CT (x' + y')%R); intros p Ep.
case (errorBoundedPlus b radix precision) with (6 := Ep); auto with zarith;
intros q (Eq1, (Eq2, Eq3)).
exists (q :: p :: output'); repeat split.
intros H2; rewrite <- H; simpl in |- *.
apply le_trans with (S (S (length l + 0))); auto with arith zarith.
rewrite <- H; simpl in |- *; fold radix in |- *; rewrite Eq1.
replace
(FtoR radix x' + FtoR radix y' - FtoR radix p + (FtoR radix p + sum output'))%R
with (FtoRradix x' + (FtoRradix y' + sum output'))%R;
[ rewrite <- E3; simpl in |- * | fold FtoRradix in |- *; ring ].
rewrite (FzeroisReallyZero radix); fold FtoRradix; ring.
cut (Fbounded b p);
[ intros Bp
| apply (RoundedModeBounded b radix) with (P := Closest b radix) (2 := Ep);
apply ClosestRoundedModeP with precision; auto with float ].
apply IsRleEpsExpTop; auto.
case (Req_dec p 0); intros Hp0.
cut (FtoRradix q = 0%R); [ intros Hq0 | idtac ].
rewrite Hp0; rewrite Hq0; rewrite Rabs_R0; right; ring.
unfold FtoRradix in |- *; rewrite Eq1; apply TwoSumNul with b precision; auto.
apply Rle_trans with (Rabs p * / 2%nat * (radix * / Zpos (vNum b)))%R;
auto.
unfold FtoRradix in |- *; rewrite Eq1; rewrite <- Rabs_Ropp.
replace (- (FtoR radix x' + FtoR radix y' - FtoR radix p))%R with
(FtoR radix p - (FtoR radix x' + FtoR radix y'))%R;
[ idtac | ring ].
apply plusErrorBound1bis with precision; auto with real zarith.
Contradict Hp0; unfold FtoRradix in |- *; apply is_Fzero_rep1; auto.
rewrite Rmult_assoc; rewrite Rmult_comm; apply Rmult_le_compat_r;
auto with real.
rewrite <- Rmult_assoc.
replace (/ 2%nat * radix)%R with 1%R;
[ idtac | rewrite Rinv_l; auto with real zarith ].
apply
Rle_trans with ((6%nat * length input + 6%nat) * / Zpos (vNum b))%R;
[ apply Rmult_le_compat_r | apply Rmult_le_compat_l ];
auto with real arith zarith.
replace 1%R with (INR 1); repeat rewrite <- mult_INR;
repeat rewrite <- plus_INR; simpl in |- *; try rewrite <- plus_n_Sm;
auto with real zarith.
replace 1%R with (INR 1); repeat rewrite <- mult_INR;
repeat rewrite <- plus_INR; simpl in |- *; auto with real arith.
apply Rle_Rinv.
apply Rlt_Rminus_ZERO.
apply Rmult_lt_reg_l with (/ (pPred (vNum b) - 1%nat))%R;
auto with real zarith.
rewrite Rmult_comm; rewrite Rinv_l; auto with real.
apply Rle_trans with (Zpos (vNum b) - 6%nat * length input)%R;
auto with real arith zarith.
unfold Rminus in |- *; apply Rplus_le_compat_r; auto with real arith zarith.
replace (- 1%nat)%R with (IZR (- (1)));
[ rewrite <- plus_IZR; unfold pPred, Zpred in |- *; auto with real zarith
| simpl in |- *; auto ].
pattern (IZR (Zpos (vNum b))) at 2 in |- *;
replace (IZR (Zpos (vNum b))) with (Zpos (vNum b) - 0)%R;
unfold Rminus in |- *; try ring; apply Rplus_le_compat_l;
apply Ropp_le_contravar; auto with real arith zarith.
apply Rle_mult_pos; replace 0%R with (INR 0); auto with real arith zarith.
CaseEq output'.
intros; apply IsRleEpsExpSingle; auto.
intros f0 l0 Ho; apply IsRleEpsExpTop; auto.
3: rewrite <- Ho; auto.
rewrite Ho in E6; inversion E6; auto.
case E9.
rewrite Ho; intros; discriminate.
rewrite Ho; intros E10.
apply
Rle_trans
with (/ (pPred (vNum b) - 1%nat) * pPred (vNum b) * Rabs (x' + y'))%R.
apply
Rle_trans
with
(Rabs (FtoR 2%nat x' + FtoR 2%nat y') *
(2%nat * pPred (vNum b) * / (2%nat * pPred (vNum b) - radix)))%R.
unfold FtoRradix in |- *; apply RoundBound with precision; auto.
replace radix with (2%nat * 1)%Z; auto with zarith.
apply Zmult_gt_0_lt_compat_l; auto with zarith.
red in |- *; simpl in |- *; auto.
right.
replace (2%nat * pPred (vNum b) - radix)%R with
(2%nat * (pPred (vNum b) - 1%nat))%R; [ idtac | simpl in |- *; ring ].
rewrite Rinv_mult_distr; auto with real zarith arith.
pattern (IZR (pPred (vNum b))) at 4 in |- *;
replace (IZR (pPred (vNum b))) with
(2%nat * / 2%nat * IZR (pPred (vNum b)))%R.
simpl in |- *; fold radix FtoRradix in |- *; ring.
rewrite Rinv_r; auto with real zarith.
apply Rmult_le_reg_l with (/ pPred (vNum b))%R; auto with real zarith.
apply Rmult_le_reg_l with (pPred (vNum b) - 1%nat)%R; auto with real zarith.
apply
Rle_trans
with
(Rabs (x' + y') *
((pPred (vNum b) - 1%nat) * / (pPred (vNum b) - 1%nat) *
(pPred (vNum b) * / pPred (vNum b))))%R; [ right; ring; ring | idtac ].
repeat rewrite Rinv_r; auto with real zarith.
apply Rle_trans with (Rabs (x' + y')); [ right; ring | idtac ].
apply
Rle_trans
with
((pPred (vNum b) - 1%nat) * / pPred (vNum b) *
((6%nat * length input + 6%nat) *
/ (pPred (vNum b) - 1%nat - 6%nat * length input)) *
Rabs (hd b (f0 :: l0)))%R; [ auto | simpl in |- *; right; ring; ring ].
Qed.
Theorem FexpAdd_main :
IsExp b input ->
(6%nat * length input * / (pPred (vNum b) - 1%nat) < 1)%R ->
exists output : list float,
(input <> nil -> length output <= S (length input)) /\
sum input = sum output :>R /\ IsRleEpsExp output.
intros H H0; apply (FexpAdd_aux2 input); auto.
Qed.
End FexpAdd.