Library Float.Expansions.Fast2Diff

Require Export EFast2Sum.
Section EDiff.
Variable b : Fbound.
Variable precision : nat.

Let radix := 2%Z.

Coercion Local FtoRradix := FtoR radix.
Hypothesis precisionGreaterThanOne : 1 < precision.
Hypothesis pGivesBound : Zpos (vNum b) = Zpower_nat radix precision.
Variable Iplus : floatfloatfloat.
Hypothesis
  IplusCorrect :
     p q : float,
    Fbounded b pFbounded b qClosest b radix (p + q) (Iplus p q).
Hypothesis
  IplusComp :
     p q r s : float,
    Fbounded b p
    Fbounded b q
    Fbounded b r
    Fbounded b sp = r :>Rq = s :>RIplus p q = Iplus r s :>R.
Hypothesis IplusSym : p q : float, Iplus p q = Iplus q p.
Hypothesis
  IplusOp : p q : float, Fopp (Iplus p q) = Iplus (Fopp p) (Fopp q).
Variable Iminus : floatfloatfloat.
Hypothesis IminusPlus : p q : float, Iminus p q = Iplus p (Fopp q).

Theorem MDekkerDiffAux1 :
  p q : float,
 Iminus p (Iminus p q) = (p - Iminus p q)%R :>R
 Fbounded b p
 Fbounded b q
 Iminus (Iminus p (Iminus p q)) q = (p - q - Iminus p q)%R :>R.
intros p q H' H'0 H'1.
elim
 (ErrorBoundedIplus b precision) with (Iplus := Iplus) (p := p) (q := Fopp q);
 fold radix in |- *; auto.
intros error H'2; elim H'2; intros H'3 H'4; clear H'2.
cut
 (Closest b radix (Iminus p (Iminus p q) - q)
    (Iminus (Iminus p (Iminus p q)) q)); auto.
rewrite H'.
replace (p - Iminus p q - q)%R with (p - q - Iminus p q)%R; [ idtac | ring ].
replace (p - q)%R with (p + - q)%R; [ idtac | ring ].
rewrite (IminusPlus p q).
unfold FtoRradix in |- *; rewrite <- Fopp_correct.
rewrite <- H'3.
intros H'2.
apply sym_eq; apply (ClosestIdem b radix); auto.
apply (IminusCorrect b Iplus); auto.
apply (IminusBounded b Iplus); auto.
apply (IminusBounded b Iplus); auto.
apply oppBounded; auto.
Qed.

Theorem MDekkerDiff :
  p q : float,
 Fbounded b p
 Fbounded b q
 (Rabs q Rabs p)%RIminus p (Iminus p q) = (p - Iminus p q)%R :>R.
intros p q H' H'0 H'1.
pattern (Iminus p q) at 2 in |- *; rewrite IminusPlus.
replace (p - Iplus p (Fopp q))%R with (- (Iplus p (Fopp q) - p))%R;
 [ idtac | ring ].
unfold FtoRradix in |- *;
 rewrite <- (MDekker b precision) with (Iminus := Iminus);
 auto.
rewrite <- Fopp_correct.
repeat rewrite IminusPlus || rewrite IplusOp || rewrite Fopp_Fopp.
rewrite IplusSym; auto.
apply oppBounded; auto.
rewrite Fopp_correct; auto.
rewrite Rabs_Ropp; auto.
Qed.

Theorem DekkerDiff :
  p q : float,
 Fbounded b p
 Fbounded b q
 (Rabs q Rabs p)%R
 Iminus (Iminus p (Iminus p q)) q = (p - q - Iminus p q)%R :>R.
intros p q H' H'0 H'1.
apply MDekkerDiffAux1; auto.
apply MDekkerDiff; auto.
Qed.

Theorem ExtMDekkerDiff :
  p q : float,
 Fbounded b p
 Fbounded b q
 (Fexp q Fexp p)%ZIminus p (Iminus p q) = (p - Iminus p q)%R :>R.
intros p q H' H'0 H'1.
pattern (Iminus p q) at 2 in |- *; rewrite IminusPlus.
replace (p - Iplus p (Fopp q))%R with (- (Iplus p (Fopp q) - p))%R;
 [ idtac | ring ].
unfold FtoRradix in |- *;
 rewrite <- (ExtMDekker b precision) with (Iminus := Iminus);
 auto.
rewrite <- Fopp_correct.
repeat rewrite IminusPlus || rewrite IplusOp || rewrite Fopp_Fopp.
rewrite IplusSym; auto.
apply oppBounded; auto.
Qed.

Theorem ExtDekkerDiff :
  p q : float,
 Fbounded b p
 Fbounded b q
 (Fexp q Fexp p)%Z
 Iminus (Iminus p (Iminus p q)) q = (p - q - Iminus p q)%R :>R.
intros p q H' H'0 H'1.
apply MDekkerDiffAux1; auto.
apply ExtMDekkerDiff; auto.
Qed.
End EDiff.