Library Flocq.Calc.Div

This file is part of the Flocq formalization of floating-point arithmetic in Coq:
Copyright (C) 2010-2018 Sylvie Boldo
Copyright (C) 2010-2018 Guillaume Melquiond
This library is free software; you can redistribute it and/or modify it under the terms of the GNU Lesser General Public License as published by the Free Software Foundation; either version 3 of the License, or (at your option) any later version.
This library is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the COPYING file for more details.

Helper function and theorem for computing the rounded quotient of two floating-point numbers.

From Coq Require Import ZArith Reals Lia.

Require Import Zaux Raux Defs Generic_fmt Float_prop Digits Bracket.

Set Implicit Arguments.

Section Fcalc_div.

Variable beta : radix.
Notation bpow e := (bpow beta e).

Variable fexp : Z Z.

Computes a mantissa of precision p, the corresponding exponent, and the position with respect to the real quotient of the input floating-point numbers.
The algorithm performs the following steps:
  • Shift dividend mantissa so that it has at least p2 + p digits.
  • Perform the Euclidean division.
  • Compute the position according to the division remainder.
Complexity is fine as long as p1 <= 2p and p2 <= p.

Lemma mag_div_F2R :
   m1 e1 m2 e2,
  (0 < m1)%Z (0 < m2)%Z
  let e := ((Zdigits beta m1 + e1) - (Zdigits beta m2 + e2))%Z in
  (e mag beta (F2R (Float beta m1 e1) / F2R (Float beta m2 e2)) e + 1)%Z.

Definition Fdiv_core m1 e1 m2 e2 e :=
  let (m1', m2') :=
    if Zle_bool e (e1 - e2)%Z
    then (m1 × Zpower beta (e1 - e2 - e), m2)%Z
    else (m1, m2 × Zpower beta (e - (e1 - e2)))%Z in
  let '(q, r) := Z.div_eucl m1' m2' in
  (q, new_location m2' r loc_Exact).

Theorem Fdiv_core_correct :
   m1 e1 m2 e2 e,
  (0 < m1)%Z (0 < m2)%Z
  let '(m, l) := Fdiv_core m1 e1 m2 e2 e in
  inbetween_float beta m e (F2R (Float beta m1 e1) / F2R (Float beta m2 e2)) l.

Definition Fdiv (x y : float beta) :=
  let (m1, e1) := x in
  let (m2, e2) := y in
  let e' := ((Zdigits beta m1 + e1) - (Zdigits beta m2 + e2))%Z in
  let e := Z.min (Z.min (fexp e') (fexp (e' + 1))) (e1 - e2) in
  let '(m, l) := Fdiv_core m1 e1 m2 e2 e in
  (m, e, l).

Theorem Fdiv_correct :
   x y,
  (0 < F2R x)%R (0 < F2R y)%R
  let '(m, e, l) := Fdiv x y in
  (e cexp beta fexp (F2R x / F2R y))%Z
  inbetween_float beta m e (F2R x / F2R y) l.

End Fcalc_div.