Library Flocq.Core.Ulp

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

Unit in the Last Place: our definition using fexp and its properties, successor and predecessor

Require Import Reals Psatz.
Require Import Raux Defs Round_pred Generic_fmt Float_prop.

Section Fcore_ulp.

Variable beta : radix.

Notation bpow e := (bpow beta e).

Variable fexp : Z Z.

Definition and basic properties about the minimal exponent, when it exists

Lemma Z_le_dec_aux: x y : Z, (x y)%Z ¬ (x y)%Z.

negligible_exp is either none (as in FLX) or Some n such that n <= fexp n.
Definition negligible_exp: option Z :=
  match (LPO_Z _ (fun zZ_le_dec_aux z (fexp z))) with
   | inleft NSome (proj1_sig N)
   | inright _None
 end.

Inductive negligible_exp_prop: option Z Prop :=
  | negligible_None: ( n, (fexp n < n)%Z) negligible_exp_prop None
  | negligible_Some: n, (n fexp n)%Z negligible_exp_prop (Some n).

Lemma negligible_exp_spec: negligible_exp_prop negligible_exp.

Lemma negligible_exp_spec': (negligible_exp = None n, (fexp n < n)%Z)
            n, (negligible_exp = Some n (n fexp n)%Z).

Context { valid_exp : Valid_exp fexp }.

Lemma fexp_negligible_exp_eq: n m, (n fexp n)%Z (m fexp m)%Z fexp n = fexp m.

Definition and basic properties about the ulp Now includes a nice ulp(0): ulp(0) is now 0 when there is no minimal exponent, such as in FLX, and beta^(fexp n) when there is a n such that n <= fexp n. For instance, the value of ulp(O) is then beta^emin in FIX and FLT. The main lemma to use is ulp_neq_0 that is equivalent to the previous "unfold ulp" provided the value is not zero.

Definition ulp x := match Req_bool x 0 with
  | truematch negligible_exp with
            | Some nbpow (fexp n)
            | None ⇒ 0%R
            end
  | falsebpow (cexp beta fexp x)
 end.

Lemma ulp_neq_0 :
   x, x 0%R
  ulp x = bpow (cexp beta fexp x).

Notation F := (generic_format beta fexp).

Theorem ulp_opp :
   x, ulp (- x) = ulp x.

Theorem ulp_abs :
   x, ulp (Rabs x) = ulp x.

Theorem ulp_ge_0:
   x, (0 ulp x)%R.

Theorem ulp_le_id:
   x,
    (0 < x)%R
    F x
    (ulp x x)%R.

Theorem ulp_le_abs:
   x,
    (x 0)%R
    F x
    (ulp x Rabs x)%R.

Theorem round_UP_DN_ulp :
   x, ¬ F x
  round beta fexp Zceil x = (round beta fexp Zfloor x + ulp x)%R.

Theorem ulp_bpow :
   e, ulp (bpow e) = bpow (fexp (e + 1)).

Lemma generic_format_ulp_0 :
  F (ulp 0).

Lemma generic_format_bpow_ge_ulp_0 :
   e, (ulp 0 bpow e)%R
  F (bpow e).

The three following properties are equivalent: Exp_not_FTZ ; forall x, F (ulp x) ; forall x, ulp 0 <= ulp x

Lemma generic_format_ulp :
  Exp_not_FTZ fexp
   x, F (ulp x).

Lemma not_FTZ_generic_format_ulp :
  ( x, F (ulp x))
  Exp_not_FTZ fexp.

Lemma ulp_ge_ulp_0 :
  Exp_not_FTZ fexp
   x, (ulp 0 ulp x)%R.

Lemma not_FTZ_ulp_ge_ulp_0:
  ( x, (ulp 0 ulp x)%R)
  Exp_not_FTZ fexp.

Lemma ulp_le_pos :
   { Hm : Monotone_exp fexp },
   x y: R,
  (0 x)%R (x y)%R
  (ulp x ulp y)%R.

Theorem ulp_le :
   { Hm : Monotone_exp fexp },
   x y: R,
  (Rabs x Rabs y)%R
  (ulp x ulp y)%R.

Properties when there is no minimal exponent
Theorem eq_0_round_0_negligible_exp :
   negligible_exp = None rnd {Vr: Valid_rnd rnd} x,
     round beta fexp rnd x = 0%R x = 0%R.

Definition and properties of pred and succ

Definition pred_pos x :=
  if Req_bool x (bpow (mag beta x - 1)) then
    (x - bpow (fexp (mag beta x - 1)))%R
  else
    (x - ulp x)%R.

Definition succ x :=
  if (Rle_bool 0 x) then
    (x+ulp x)%R
  else
    (- pred_pos (-x))%R.

Definition pred x := (- succ (-x))%R.

Theorem pred_eq_pos :
   x, (0 x)%R
  pred x = pred_pos x.

Theorem succ_eq_pos :
   x, (0 x)%R
  succ x = (x + ulp x)%R.

Theorem succ_opp :
   x, succ (-x) = (- pred x)%R.

Theorem pred_opp :
   x, pred (-x) = (- succ x)%R.

pred and succ are in the format

Theorem id_m_ulp_ge_bpow :
   x e, F x
  x ulp x
  (bpow e < x)%R
  (bpow e x - ulp x)%R.

Theorem id_p_ulp_le_bpow :
   x e, (0 < x)%R F x
  (x < bpow e)%R
  (x + ulp x bpow e)%R.

Lemma generic_format_pred_aux1:
   x, (0 < x)%R F x
  x bpow (mag beta x - 1)
  F (x - ulp x).

Lemma generic_format_pred_aux2 :
   x, (0 < x)%R F x
  let e := mag_val beta x (mag beta x) in
  x = bpow (e - 1)
  F (x - bpow (fexp (e - 1))).

Lemma generic_format_succ_aux1 :
   x, (0 < x)%R F x
  F (x + ulp x).

Lemma generic_format_pred_pos :
   x, F x (0 < x)%R
  F (pred_pos x).

Theorem generic_format_succ :
   x, F x
  F (succ x).

Theorem generic_format_pred :
   x, F x
  F (pred x).

Lemma pred_pos_lt_id :
   x, (x 0)%R
  (pred_pos x < x)%R.

Theorem succ_gt_id :
   x, (x 0)%R
  (x < succ x)%R.

Theorem pred_lt_id :
   x, (x 0)%R
  (pred x < x)%R.

Theorem succ_ge_id :
   x, (x succ x)%R.

Theorem pred_le_id :
   x, (pred x x)%R.

Lemma pred_pos_ge_0 :
   x,
  (0 < x)%R F x (0 pred_pos x)%R.

Theorem pred_ge_0 :
   x,
  (0 < x)%R F x (0 pred x)%R.

Lemma pred_pos_plus_ulp_aux1 :
   x, (0 < x)%R F x
  x bpow (mag beta x - 1)
  ((x - ulp x) + ulp (x-ulp x) = x)%R.

Lemma pred_pos_plus_ulp_aux2 :
   x, (0 < x)%R F x
  let e := mag_val beta x (mag beta x) in
  x = bpow (e - 1)
  (x - bpow (fexp (e-1)) 0)%R
  ((x - bpow (fexp (e-1))) + ulp (x - bpow (fexp (e-1))) = x)%R.

Lemma pred_pos_plus_ulp_aux3 :
   x, (0 < x)%R F x
  let e := mag_val beta x (mag beta x) in
  x = bpow (e - 1)
  (x - bpow (fexp (e-1)) = 0)%R
  (ulp 0 = x)%R.

The following one is false for x = 0 in FTZ

Lemma pred_pos_plus_ulp :
   x, (0 < x)%R F x
  (pred_pos x + ulp (pred_pos x) = x)%R.

Theorem pred_plus_ulp :
   x, (0 < x)%R F x
  (pred x + ulp (pred x))%R = x.

Rounding x + small epsilon

Theorem mag_plus_eps :
   x, (0 < x)%R F x
   eps, (0 eps < ulp x)%R
  mag beta (x + eps) = mag beta x :> Z.

Theorem round_DN_plus_eps_pos :
   x, (0 x)%R F x
   eps, (0 eps < ulp x)%R
  round beta fexp Zfloor (x + eps) = x.

Theorem round_UP_plus_eps_pos :
   x, (0 x)%R F x
   eps, (0 < eps ulp x)%R
  round beta fexp Zceil (x + eps) = (x + ulp x)%R.

Theorem round_UP_pred_plus_eps_pos :
   x, (0 < x)%R F x
   eps, (0 < eps ulp (pred x) )%R
  round beta fexp Zceil (pred x + eps) = x.

Theorem round_DN_minus_eps_pos :
   x, (0 < x)%R F x
   eps, (0 < eps ulp (pred x))%R
  round beta fexp Zfloor (x - eps) = pred x.

Theorem round_DN_plus_eps:
   x, F x
   eps, (0 eps < if (Rle_bool 0 x) then (ulp x)
                                     else (ulp (pred (-x))))%R
  round beta fexp Zfloor (x + eps) = x.

Theorem round_UP_plus_eps :
   x, F x
   eps, (0 < eps if (Rle_bool 0 x) then (ulp x)
                                     else (ulp (pred (-x))))%R
  round beta fexp Zceil (x + eps) = (succ x)%R.

Lemma le_pred_pos_lt :
   x y,
  F x F y
  (0 x < y)%R
  (x pred_pos y)%R.

Lemma succ_le_lt_aux:
   x y,
  F x F y
  (0 x)%R (x < y)%R
  (succ x y)%R.

Theorem succ_le_lt:
   x y,
  F x F y
  (x < y)%R
  (succ x y)%R.

Theorem pred_ge_gt :
   x y,
  F x F y
  (x < y)%R
  (x pred y)%R.

Theorem succ_gt_ge :
   x y,
  (y 0)%R
  (x y)%R
  (x < succ y)%R.

Theorem pred_lt_le :
   x y,
  (x 0)%R
  (x y)%R
  (pred x < y)%R.

Lemma succ_pred_pos :
   x, F x (0 < x)%R succ (pred x) = x.

Theorem pred_ulp_0 :
  pred (ulp 0) = 0%R.

Theorem succ_0 :
  succ 0 = ulp 0.

Theorem pred_0 :
  pred 0 = Ropp (ulp 0).

Lemma pred_succ_pos :
   x, F x (0 < x)%R
  pred (succ x) = x.

Theorem succ_pred :
   x, F x
  succ (pred x) = x.

Theorem pred_succ :
   x, F x
  pred (succ x) = x.

Theorem round_UP_pred_plus_eps :
   x, F x
   eps, (0 < eps if Rle_bool x 0 then ulp x
                          else ulp (pred x))%R
  round beta fexp Zceil (pred x + eps) = x.

Theorem round_DN_minus_eps:
   x, F x
   eps, (0 < eps if (Rle_bool x 0) then (ulp x)
                                     else (ulp (pred x)))%R
  round beta fexp Zfloor (x - eps) = pred x.

Error of a rounding, expressed in number of ulps false for x=0 in the FLX format
Theorem error_lt_ulp :
   rnd { Zrnd : Valid_rnd rnd } x,
  (x 0)%R
  (Rabs (round beta fexp rnd x - x) < ulp x)%R.

Theorem error_le_ulp :
   rnd { Zrnd : Valid_rnd rnd } x,
  (Rabs (round beta fexp rnd x - x) ulp x)%R.

Theorem error_le_half_ulp :
   choice x,
  (Rabs (round beta fexp (Znearest choice) x - x) /2 × ulp x)%R.

Theorem ulp_DN :
   x, (0 x)%R
  ulp (round beta fexp Zfloor x) = ulp x.

Theorem round_neq_0_negligible_exp :
  negligible_exp = None rnd { Zrnd : Valid_rnd rnd } x,
  (x 0)%R (round beta fexp rnd x 0)%R.

allows rnd x to be 0
Theorem error_lt_ulp_round :
   { Hm : Monotone_exp fexp } rnd { Zrnd : Valid_rnd rnd } x,
  (x 0)%R
  (Rabs (round beta fexp rnd x - x) < ulp (round beta fexp rnd x))%R.

allows both x and rnd x to be 0
Theorem error_le_half_ulp_round :
   { Hm : Monotone_exp fexp },
   choice x,
  (Rabs (round beta fexp (Znearest choice) x - x) /2 × ulp (round beta fexp (Znearest choice) x))%R.

Theorem pred_le :
   x y, F x F y (x y)%R
  (pred x pred y)%R.

Theorem succ_le :
   x y, F x F y (x y)%R
  (succ x succ y)%R.

Theorem pred_le_inv: x y, F x F y
    (pred x pred y)%R (x y)%R.

Theorem succ_le_inv: x y, F x F y
    (succ x succ y)%R (x y)%R.

Theorem pred_lt :
   x y, F x F y (x < y)%R
  (pred x < pred y)%R.

Theorem succ_lt :
   x y, F x F y (x < y)%R
  (succ x < succ y)%R.

Theorem round_DN_ge_UP_gt :
   x y, F y
  (y < round beta fexp Zceil x y round beta fexp Zfloor x)%R.

Theorem round_UP_le_DN_lt :
   x y, F y
  (round beta fexp Zfloor x < y round beta fexp Zceil x y)%R.

Theorem pred_UP_le_DN :
   x, (pred (round beta fexp Zceil x) round beta fexp Zfloor x)%R.

Theorem UP_le_succ_DN :
   x, (round beta fexp Zceil x succ (round beta fexp Zfloor x))%R.

Theorem pred_UP_eq_DN :
   x, ¬ F x
  (pred (round beta fexp Zceil x) = round beta fexp Zfloor x)%R.

Theorem succ_DN_eq_UP :
   x, ¬ F x
  (succ (round beta fexp Zfloor x) = round beta fexp Zceil x)%R.

Theorem round_DN_eq :
   x d, F d (d x < succ d)%R
  round beta fexp Zfloor x = d.

Theorem round_UP_eq :
   x u, F u (pred u < x u)%R
  round beta fexp Zceil x = u.

Lemma ulp_ulp_0 : {H : Exp_not_FTZ fexp},
  ulp (ulp 0) = ulp 0.

Lemma ulp_succ_pos : x, F x (0 < x)%R
   ulp (succ x) = ulp x succ x = bpow (mag beta x).

Lemma ulp_round_pos :
   { Not_FTZ_ : Exp_not_FTZ fexp},
    rnd { Zrnd : Valid_rnd rnd } x,
  (0 < x)%R ulp (round beta fexp rnd x) = ulp x
      round beta fexp rnd x = bpow (mag beta x).

Theorem ulp_round : { Not_FTZ_ : Exp_not_FTZ fexp},
    rnd { Zrnd : Valid_rnd rnd } x,
     ulp (round beta fexp rnd x) = ulp x
          Rabs (round beta fexp rnd x) = bpow (mag beta x).

Properties of rounding to nearest and ulp

Theorem round_N_le_midp: choice u v,
  F u (v < (u + succ u)/2)%R
       (round beta fexp (Znearest choice) v u)%R.

Theorem round_N_ge_midp: choice u v,
 F u ((u + pred u)/2 < v)%R
       (u round beta fexp (Znearest choice) v)%R.

Lemma round_N_eq_DN: choice x,
       let d:=round beta fexp Zfloor x in
       let u:=round beta fexp Zceil x in
      (x<(d+u)/2)%R
     round beta fexp (Znearest choice) x = d.

Lemma round_N_eq_DN_pt: choice x d u,
      Rnd_DN_pt F x d Rnd_UP_pt F x u
      (x<(d+u)/2)%R
     round beta fexp (Znearest choice) x = d.

Lemma round_N_eq_UP: choice x,
      let d:=round beta fexp Zfloor x in
      let u:=round beta fexp Zceil x in
     ((d+u)/2 < x)%R
     round beta fexp (Znearest choice) x = u.

Lemma round_N_eq_UP_pt: choice x d u,
      Rnd_DN_pt F x d Rnd_UP_pt F x u
      ((d+u)/2 < x)%R
     round beta fexp (Znearest choice) x = u.

End Fcore_ulp.