package codex

  1. Overview
  2. Docs
Legend:
Page
Library
Module
Module type
Parameter
Class
Class type
Source

Source file fc_float.ml

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
(**************************************************************************)
(*  This file is part of the Codex semantics library.                     *)
(*                                                                        *)
(*  Copyright (C) 2013-2025                                               *)
(*    CEA (Commissariat à l'énergie atomique et aux énergies              *)
(*         alternatives)                                                  *)
(*                                                                        *)
(*  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, version 2.1.                                              *)
(*                                                                        *)
(*  It 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         *)
(*  GNU Lesser General Public License for more details.                   *)
(*                                                                        *)
(*  See the GNU Lesser General Public License version 2.1                 *)
(*  for more details (enclosed in the file LICENSE).                      *)
(*                                                                        *)
(**************************************************************************)

open Float_sig

let set_rounding_mode = function
  | Down -> Floating_point.set_round_downward ()
  | Up -> Floating_point.set_round_upward ()
  | Zero -> Floating_point.set_round_toward_zero ()
  | Near -> ()

let (>>%) round f =
  set_rounding_mode round;
  let result = f () in
  if round <> Near then Floating_point.set_round_nearest_even ();
  result

let is_single = function Single -> true | _ -> false

let round_to_precision prec f =
  if is_single prec then Floating_point.round_to_single_precision_float f else f

type t = float

(* let packed_descr = Structural_descr.p_float *)

let hash = Hashtbl.hash
let pretty = Floating_point.pretty

let is_exact = function
  | Single | Double -> true
  | Long_Double | Real -> false

let cmp_ieee = (compare: float -> float -> int)

(** NOTE: all floating-point comparisons using OCaml's standard operators
    do NOT distinguish between -0.0 and 0.0.
    Whenever floats are compared using them, it implies that negative zeroes
    are also considered, e.g. "if x < 0.0" is equivalent to "if x < -0.0",
    which is also equivalent to "F.compare x (-0.0) < 0".
    This 'compare' operator distinguishes -0. and 0. *)
(* replace "noalloc" with [@@noalloc] for OCaml version >= 4.03.0 *)
[@@@ warning "-3"]
external compare : float -> float -> int = "fc_ival_float_compare_total" "noalloc"
[@@@ warning "+3"]

let of_float round prec f = round >>% fun () -> round_to_precision prec f

let to_float f = f

let is_nan a = classify_float a = FP_nan
let is_infinite f = classify_float f = FP_infinite
let is_finite f = match classify_float f with
  | FP_nan | FP_infinite -> false
  | FP_normal | FP_subnormal | FP_zero -> true

(* replace "noalloc" with [@@noalloc] for OCaml version >= 4.03.0 *)
[@@@ warning "-3"]
external is_negative : float -> bool = "fc_ival_float_is_negative" "noalloc"
[@@@ warning "+3"]

let round_to_precision round prec t =
  if is_single prec
  then round >>% fun () -> Floating_point.round_to_single_precision_float t
  else t

(* Wrong for [next -0.] and [prev 0.], as this function uses the binary
   representation of floating-point values, which is not continuous at 0. *)
let next_previous int64fup int64fdown float =
  let r = Int64.bits_of_float float in
  let f = if r >= 0L then int64fup else int64fdown in
  let f = Int64.float_of_bits (f r) in
  match classify_float f with
  | FP_nan -> float (* can only be produced from an infinity or a nan *)
  | FP_infinite | FP_normal | FP_subnormal | FP_zero -> f

let next_float prec f =
  if not (is_exact prec) then f else
  if compare f (-0.) = 0 then 0. else
    let f = next_previous Int64.succ Int64.pred f in
    if is_single prec
    then Up >>% fun () -> Floating_point.round_to_single_precision_float f
    else f

let prev_float prec f =
  if not (is_exact prec) then f else
  if compare f 0. = 0 then -0. else
    let f = next_previous Int64.pred Int64.succ f in
    if is_single prec
    then Down >>% fun () -> Floating_point.round_to_single_precision_float f
    else f

let m_pi = 3.1415929794311523 (* single-precision *)
let m_pi_2 = 1.5707964897155761 (* single-precision *)
let max_single_precision_float = Floating_point.max_single_precision_float

let le f1 f2 = compare f1 f2 <= 0

let widen_up f =
  if le f (-0.) then -0.
  else if le f 0. then 0.
  else if le f 1. then 1.
  else if le f m_pi_2 then m_pi_2
  else if le f m_pi then m_pi
  else if le f 10. then 10.
  else if le f 1e10 then 1e10
  else if le f max_single_precision_float then max_single_precision_float
  else if le f 1e80 then 1e80
  else if le f max_float then max_float
  else infinity

let neg = (~-.)
let abs = abs_float

let floor = floor
let ceil = ceil
let trunc = Floating_point.trunc
let fround = Floating_point.fround

let binary op = fun round prec x y ->
  round >>% fun () ->
  if is_single prec
  then Floating_point.round_to_single_precision_float (op x y)
  else op x y

let add round = binary ( +. ) round
let sub round = binary ( -. ) round
let mul round = binary ( *. ) round
let div round = binary ( /. ) round
let fmod round = binary mod_float round

let generate single double =
  fun round prec ->
    round >>% fun () -> if is_single prec then single else double

let exp round = generate Floating_point.expf exp round
let log round = generate Floating_point.logf log round
let log10 round = generate Floating_point.log10f log10 round
let sqrt round = generate Floating_point.sqrtf sqrt round
let pow round = generate Floating_point.powf ( ** ) round

let cos round = generate Floating_point.cosf cos round
let sin round = generate Floating_point.sinf sin round
let atan2 round = generate Floating_point.atan2f atan2 round