File: coq_float64.c

package info (click to toggle)
coq 8.20.1%2Bdfsg-1
  • links: PTS, VCS
  • area: main
  • in suites: forky, sid, trixie
  • size: 44,116 kB
  • sloc: ml: 234,160; sh: 4,301; python: 3,270; ansic: 2,644; makefile: 882; lisp: 172; javascript: 63; xml: 24; sed: 2
file content (74 lines) | stat: -rw-r--r-- 2,847 bytes parent folder | download | duplicates (4)
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
/************************************************************************/
/*         *   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 <math.h>
#include <stdint.h>

#define CAML_INTERNALS
#include <caml/alloc.h>

#include "coq_values.h"

union double_bits {
  double d;
  uint64_t u;
};

static double next_up(double x) {
  union double_bits bits;
  if (!(x < INFINITY)) return x; // x is +oo or NaN
  bits.d = x;
  int64_t i = bits.u;
  if (i >= 0) ++bits.u; // x >= +0.0, go away from zero
  else if (bits.u + bits.u == 0) bits.u = 1; // x is -0.0, should almost never happen
  else --bits.u; // x < 0.0, go toward zero
  return bits.d;
}

static double next_down(double x) {
  union double_bits bits;
  if (!(x > -INFINITY)) return x; // x is -oo or NaN
  bits.d = x;
  int64_t i = bits.u;
  if (i == 0) bits.u = INT64_MIN + 1; // x is +0.0
  else if (i < 0) ++bits.u; // x <= -0.0, go away from zero
  else --bits.u; // x > 0.0, go toward zero
  return bits.d;
}

#define DECLARE_FBINOP(name, e)                                         \
  double coq_##name(double x, double y) {                               \
    return e;                                                           \
  }                                                                     \
                                                                        \
  value coq_##name##_byte(value x, value y) {                           \
    return caml_copy_double(coq_##name(Double_val(x), Double_val(y)));  \
  }

#define DECLARE_FUNOP(name, e)                                          \
  double coq_##name(double x) {                                         \
    return e;                                                           \
  }                                                                     \
                                                                        \
  value coq_##name##_byte(value x) {                                    \
    return caml_copy_double(coq_##name(Double_val(x)));                 \
  }

DECLARE_FBINOP(fmul, x * y)
DECLARE_FBINOP(fadd, x + y)
DECLARE_FBINOP(fsub, x - y)
DECLARE_FBINOP(fdiv, x / y)
DECLARE_FUNOP(fsqrt, sqrt(x))
DECLARE_FUNOP(next_up, next_up(x))
DECLARE_FUNOP(next_down, next_down(x))

value coq_is_double(value x) {
  return Val_long(Is_double(x));
}