diff options
| author | Pierre Roux | 2020-10-06 16:52:03 +0200 |
|---|---|---|
| committer | Pierre Roux | 2020-10-06 18:26:38 +0200 |
| commit | 6fe8c44ff828ef4ec89b49ada634ce87639f384f (patch) | |
| tree | 43642a3c3fb5bdb5817afea42cd608c527b7044c /kernel/float64_63.ml | |
| parent | 6d3a9220204de22e0b81dc961d2eb269128b5c2e (diff) | |
Use OCaml floating-point operations on 64 bits arch
C functions were used for floating-point arithmetic operations, by
fear of double rounding that could happen on old x87 on 32 bits
architecture. This commit uses OCaml floating-point operations on 64
bits architectures.
The following snippet is made 17% faster by this commit.
From Coq Require Import Int63 BinPos PrimFloat.
Definition foo n :=
let eps := sub (next_up one) one in
Pos.iter (fun x => add x eps) two n.
Time Eval native_compute in foo 1000000000.
Diffstat (limited to 'kernel/float64_63.ml')
| -rw-r--r-- | kernel/float64_63.ml | 35 |
1 files changed, 35 insertions, 0 deletions
diff --git a/kernel/float64_63.ml b/kernel/float64_63.ml new file mode 100644 index 0000000000..0025531cb1 --- /dev/null +++ b/kernel/float64_63.ml @@ -0,0 +1,35 @@ +(************************************************************************) +(* * The Coq Proof Assistant / The Coq Development Team *) +(* v * Copyright INRIA, CNRS and contributors *) +(* <O___,, * (see version control and CREDITS file for authors & dates) *) +(* \VV/ **************************************************************) +(* // * This file is distributed under the terms of the *) +(* * GNU Lesser General Public License Version 2.1 *) +(* * (see LICENSE file for the text of the license) *) +(************************************************************************) + +include Float64_common + +let mul (x : float) (y : float) : float = x *. y +[@@ocaml.inline always] + +let add (x : float) (y : float) : float = x +. y +[@@ocaml.inline always] + +let sub (x : float) (y : float) : float = x -. y +[@@ocaml.inline always] + +let div (x : float) (y : float) : float = x /. y +[@@ocaml.inline always] + +let sqrt (x : float) : float = sqrt x +[@@ocaml.inline always] + +(*** Test at runtime that no harmful double rounding seems to + be performed with an intermediate 80 bits representation (x87). *) +let () = + let b = ldexp 1. 53 in + let s = add 1. (ldexp 1. (-52)) in + if add b s <= b || add b 1. <> b || ldexp 1. (-1074) <= 0. then + failwith "Detected non IEEE-754 compliant architecture (or wrong \ + rounding mode). Use of Float is thus unsafe." |