summaryrefslogtreecommitdiff
path: root/sysdeps/ieee754
diff options
context:
space:
mode:
Diffstat (limited to 'sysdeps/ieee754')
-rw-r--r--sysdeps/ieee754/dbl-64/wordsize-64/s_ceil.c67
-rw-r--r--sysdeps/ieee754/dbl-64/wordsize-64/s_isinf.c30
-rw-r--r--sysdeps/ieee754/dbl-64/wordsize-64/s_isnan.c39
-rw-r--r--sysdeps/ieee754/dbl-64/wordsize-64/s_llround.c1
-rw-r--r--sysdeps/ieee754/dbl-64/wordsize-64/s_lround.c80
-rw-r--r--sysdeps/ieee754/dbl-64/wordsize-64/s_nearbyint.c86
-rw-r--r--sysdeps/ieee754/dbl-64/wordsize-64/s_rint.c79
-rw-r--r--sysdeps/ieee754/dbl-64/wordsize-64/s_round.c77
-rw-r--r--sysdeps/ieee754/dbl-64/wordsize-64/s_scalbn.c68
-rw-r--r--sysdeps/ieee754/dbl-64/wordsize-64/s_trunc.c56
10 files changed, 583 insertions, 0 deletions
diff --git a/sysdeps/ieee754/dbl-64/wordsize-64/s_ceil.c b/sysdeps/ieee754/dbl-64/wordsize-64/s_ceil.c
new file mode 100644
index 0000000000..9123fdc7bd
--- /dev/null
+++ b/sysdeps/ieee754/dbl-64/wordsize-64/s_ceil.c
@@ -0,0 +1,67 @@
+/* @(#)s_ceil.c 5.1 93/09/24 */
+/*
+ * ====================================================
+ * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
+ *
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+/*
+ * ceil(x)
+ * Return x rounded toward -inf to integral value
+ * Method:
+ * Bit twiddling.
+ * Exception:
+ * Inexact flag raised if x not equal to ceil(x).
+ */
+
+#include "math.h"
+#include "math_private.h"
+
+#ifdef __STDC__
+static const double huge = 1.0e300;
+#else
+static double huge = 1.0e300;
+#endif
+
+#ifdef __STDC__
+ double __ceil(double x)
+#else
+ double __ceil(x)
+ double x;
+#endif
+{
+ int64_t i0,i;
+ int32_t j0;
+ EXTRACT_WORDS64(i0,x);
+ j0 = ((i0>>52)&0x7ff)-0x3ff;
+ if(j0<=51) {
+ if(j0<0) { /* raise inexact if x != 0 */
+ if(huge+x>0.0) {/* return 0*sign(x) if |x|<1 */
+ if(i0<0) {i0=INT64_C(0x8000000000000000);}
+ else if(i0!=0) { i0=INT64_C(0x3ff0000000000000);}
+ }
+ } else {
+ i = INT64_C(0x000fffffffffffff)>>j0;
+ if((i0&i)==0) return x; /* x is integral */
+ if(huge+x>0.0) { /* raise inexact flag */
+ if(i0>0) i0 += UINT64_C(0x0010000000000000)>>j0;
+ i0 &= (~i);
+ }
+ }
+ } else {
+ if(j0==0x400) return x+x; /* inf or NaN */
+ else return x; /* x is integral */
+ }
+ INSERT_WORDS64(x,i0);
+ return x;
+}
+weak_alias (__ceil, ceil)
+#ifdef NO_LONG_DOUBLE
+strong_alias (__ceil, __ceill)
+weak_alias (__ceil, ceill)
+#endif
diff --git a/sysdeps/ieee754/dbl-64/wordsize-64/s_isinf.c b/sysdeps/ieee754/dbl-64/wordsize-64/s_isinf.c
new file mode 100644
index 0000000000..c48e979a84
--- /dev/null
+++ b/sysdeps/ieee754/dbl-64/wordsize-64/s_isinf.c
@@ -0,0 +1,30 @@
+/*
+ * Written by J.T. Conklin <jtc@netbsd.org>.
+ * Changed to return -1 for -Inf by Ulrich Drepper <drepper@cygnus.com>.
+ * Public domain.
+ */
+
+/*
+ * isinf(x) returns 1 is x is inf, -1 if x is -inf, else 0;
+ * no branching!
+ */
+
+#include "math.h"
+#include "math_private.h"
+
+int
+__isinf (double x)
+{
+ int64_t ix;
+ EXTRACT_WORDS64(ix,x);
+ int64_t t = ix & UINT64_C(0x7fffffffffffffff);
+ t ^= UINT64_C(0x7ff0000000000000);
+ t |= -t;
+ return ~(t >> 63) & (ix >> 62);
+}
+hidden_def (__isinf)
+weak_alias (__isinf, isinf)
+#ifdef NO_LONG_DOUBLE
+strong_alias (__isinf, __isinfl)
+weak_alias (__isinf, isinfl)
+#endif
diff --git a/sysdeps/ieee754/dbl-64/wordsize-64/s_isnan.c b/sysdeps/ieee754/dbl-64/wordsize-64/s_isnan.c
new file mode 100644
index 0000000000..3b08c54dd4
--- /dev/null
+++ b/sysdeps/ieee754/dbl-64/wordsize-64/s_isnan.c
@@ -0,0 +1,39 @@
+/* @(#)s_isnan.c 5.1 93/09/24 */
+/*
+ * ====================================================
+ * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
+ *
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+/*
+ * isnan(x) returns 1 is x is nan, else 0;
+ * no branching!
+ */
+
+#include "math.h"
+#include "math_private.h"
+
+#ifdef __STDC__
+ int __isnan(double x)
+#else
+ int __isnan(x)
+ double x;
+#endif
+{
+ int64_t hx;
+ EXTRACT_WORDS64(hx,x);
+ hx &= UINT64_C(0x7fffffffffffffff);
+ hx = UINT64_C(0x7ff0000000000000) - hx;
+ return (int)(((uint64_t)hx)>>63);
+}
+hidden_def (__isnan)
+weak_alias (__isnan, isnan)
+#ifdef NO_LONG_DOUBLE
+strong_alias (__isnan, __isnanl)
+weak_alias (__isnan, isnanl)
+#endif
diff --git a/sysdeps/ieee754/dbl-64/wordsize-64/s_llround.c b/sysdeps/ieee754/dbl-64/wordsize-64/s_llround.c
new file mode 100644
index 0000000000..f7c9ea568f
--- /dev/null
+++ b/sysdeps/ieee754/dbl-64/wordsize-64/s_llround.c
@@ -0,0 +1 @@
+/* The code is the same as lround. Use an alias, see l_round.c. */
diff --git a/sysdeps/ieee754/dbl-64/wordsize-64/s_lround.c b/sysdeps/ieee754/dbl-64/wordsize-64/s_lround.c
new file mode 100644
index 0000000000..30ea5db7cc
--- /dev/null
+++ b/sysdeps/ieee754/dbl-64/wordsize-64/s_lround.c
@@ -0,0 +1,80 @@
+/* Round double value to long int.
+ Copyright (C) 1997, 2004, 2009 Free Software Foundation, Inc.
+ This file is part of the GNU C Library.
+ Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997.
+
+ The GNU C 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 2.1 of the License, or (at your option) any later version.
+
+ The GNU C 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 GNU
+ Lesser General Public License for more details.
+
+ You should have received a copy of the GNU Lesser General Public
+ License along with the GNU C Library; if not, write to the Free
+ Software Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA
+ 02111-1307 USA. */
+
+#define llround __hidden_llround
+#define __llround __hidden___llround
+
+#include <math.h>
+
+#include "math_private.h"
+
+
+long int
+__lround (double x)
+{
+ int32_t j0;
+ int64_t i0;
+ long int result;
+ int sign;
+
+ EXTRACT_WORDS64 (i0, x);
+ j0 = ((i0 >> 52) & 0x7ff) - 0x3ff;
+ sign = i0 < 0 ? -1 : 1;
+ i0 &= UINT64_C(0xfffffffffffff);
+ i0 |= UINT64_C(0x10000000000000);
+
+ if (j0 < (int32_t) (8 * sizeof (long int)) - 1)
+ {
+ if (j0 < 0)
+ return j0 < -1 ? 0 : sign;
+ else if (j0 >= 52)
+ result = i0 << (j0 - 52);
+ else
+ {
+ i0 += UINT64_C(0x8000000000000) >> j0;
+
+ result = i0 >> (52 - j0);
+ }
+ }
+ else
+ {
+ /* The number is too large. It is left implementation defined
+ what happens. */
+ return (long int) x;
+ }
+
+ return sign * result;
+}
+
+weak_alias (__lround, lround)
+#ifdef NO_LONG_DOUBLE
+strong_alias (__lround, __lroundl)
+weak_alias (__lround, lroundl)
+#endif
+
+/* long long has the same width as long on 64-bit machines. */
+#undef llround
+#undef __llround
+strong_alias (__lround, __llround)
+weak_alias (__lround, llround)
+#ifdef NO_LONG_DOUBLE
+strong_alias (__lround, __llroundl)
+weak_alias (__lround, llroundl)
+#endif
diff --git a/sysdeps/ieee754/dbl-64/wordsize-64/s_nearbyint.c b/sysdeps/ieee754/dbl-64/wordsize-64/s_nearbyint.c
new file mode 100644
index 0000000000..cb49019ddb
--- /dev/null
+++ b/sysdeps/ieee754/dbl-64/wordsize-64/s_nearbyint.c
@@ -0,0 +1,86 @@
+/* Adapted for use as nearbyint by Ulrich Drepper <drepper@cygnus.com>. */
+/*
+ * ====================================================
+ * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
+ *
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+/*
+ * rint(x)
+ * Return x rounded to integral value according to the prevailing
+ * rounding mode.
+ * Method:
+ * Using floating addition.
+ * Exception:
+ * Inexact flag raised if x not equal to rint(x).
+ */
+
+#include <fenv.h>
+#include "math.h"
+#include "math_private.h"
+
+#ifdef __STDC__
+static const double
+#else
+static double
+#endif
+TWO52[2]={
+ 4.50359962737049600000e+15, /* 0x43300000, 0x00000000 */
+ -4.50359962737049600000e+15, /* 0xC3300000, 0x00000000 */
+};
+
+#ifdef __STDC__
+ double __nearbyint(double x)
+#else
+ double __nearbyint(x)
+ double x;
+#endif
+{
+ fenv_t env;
+ int64_t i0,sx;
+ int32_t j0;
+ EXTRACT_WORDS64(i0,x);
+ sx = (i0>>63)&1;
+ j0 = ((i0>>52)&0x7ff)-0x3ff;
+ if(j0<52) {
+ if(j0<0) {
+ if((i0&UINT64_C(0x7fffffffffffffff))==0) return x;
+ uint64_t i = i0 & UINT64_C(0xfffffffffffff);
+ i0 &= UINT64_C(0xfffe000000000000);
+ i0 |= (((i|-i) >> 12) & UINT64_C(0x8000000000000));
+ INSERT_WORDS64(x,i0);
+ feholdexcept (&env);
+ double w = TWO52[sx]+x;
+ double t = w-TWO52[sx];
+ fesetenv (&env);
+ EXTRACT_WORDS64(i0,t);
+ INSERT_WORDS64(t,(i0&UINT64_C(0x7fffffffffffffff))|(sx<<63));
+ return t;
+ } else {
+ uint64_t i = UINT64_C(0x000fffffffffffff)>>j0;
+ if((i0&i)==0) return x; /* x is integral */
+ i>>=1;
+ if((i0&i)!=0)
+ i0 = (i0&(~i))|(UINT64_C(0x4000000000000)>>j0);
+ }
+ } else {
+ if(j0==0x400) return x+x; /* inf or NaN */
+ else return x; /* x is integral */
+ }
+ INSERT_WORDS64(x,i0);
+ feholdexcept (&env);
+ double w = TWO52[sx]+x;
+ double t = w-TWO52[sx];
+ fesetenv (&env);
+ return t;
+}
+weak_alias (__nearbyint, nearbyint)
+#ifdef NO_LONG_DOUBLE
+strong_alias (__nearbyint, __nearbyintl)
+weak_alias (__nearbyint, nearbyintl)
+#endif
diff --git a/sysdeps/ieee754/dbl-64/wordsize-64/s_rint.c b/sysdeps/ieee754/dbl-64/wordsize-64/s_rint.c
new file mode 100644
index 0000000000..4a60aa3278
--- /dev/null
+++ b/sysdeps/ieee754/dbl-64/wordsize-64/s_rint.c
@@ -0,0 +1,79 @@
+/* @(#)s_rint.c 5.1 93/09/24 */
+/*
+ * ====================================================
+ * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
+ *
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+/*
+ * rint(x)
+ * Return x rounded to integral value according to the prevailing
+ * rounding mode.
+ * Method:
+ * Using floating addition.
+ * Exception:
+ * Inexact flag raised if x not equal to rint(x).
+ */
+
+#include "math.h"
+#include "math_private.h"
+
+#ifdef __STDC__
+static const double
+#else
+static double
+#endif
+TWO52[2]={
+ 4.50359962737049600000e+15, /* 0x43300000, 0x00000000 */
+ -4.50359962737049600000e+15, /* 0xC3300000, 0x00000000 */
+};
+
+#ifdef __STDC__
+ double __rint(double x)
+#else
+ double __rint(x)
+ double x;
+#endif
+{
+ int64_t i0,sx;
+ int32_t j0;
+ EXTRACT_WORDS64(i0,x);
+ sx = (i0>>63)&1;
+ j0 = ((i0>>52)&0x7ff)-0x3ff;
+ if(j0<52) {
+ if(j0<0) {
+ if((i0 & UINT64_C(0x7fffffffffffffff))==0) return x;
+ uint64_t i = i0 & UINT64_C(0xfffffffffffff);
+ i0 &= UINT64_C(0xfffe000000000000);
+ i0 |= (((i|-i) >> 12) & UINT64_C(0x8000000000000));
+ INSERT_WORDS64(x,i0);
+ double w = TWO52[sx]+x;
+ double t = w-TWO52[sx];
+ EXTRACT_WORDS64(i0,t);
+ INSERT_WORDS64(t,(i0&UINT64_C(0x7fffffffffffffff))|(sx<<63));
+ return t;
+ } else {
+ uint64_t i = UINT64_C(0x000fffffffffffff)>>j0;
+ if((i0&i)==0) return x; /* x is integral */
+ i>>=1;
+ if((i0&i)!=0)
+ i0 = (i0&(~i))|(UINT64_C(0x4000000000000)>>j0);
+ }
+ } else {
+ if(j0==0x400) return x+x; /* inf or NaN */
+ else return x; /* x is integral */
+ }
+ INSERT_WORDS64(x,i0);
+ double w = TWO52[sx]+x;
+ return w-TWO52[sx];
+}
+weak_alias (__rint, rint)
+#ifdef NO_LONG_DOUBLE
+strong_alias (__rint, __rintl)
+weak_alias (__rint, rintl)
+#endif
diff --git a/sysdeps/ieee754/dbl-64/wordsize-64/s_round.c b/sysdeps/ieee754/dbl-64/wordsize-64/s_round.c
new file mode 100644
index 0000000000..5bd857910c
--- /dev/null
+++ b/sysdeps/ieee754/dbl-64/wordsize-64/s_round.c
@@ -0,0 +1,77 @@
+/* Round double to integer away from zero.
+ Copyright (C) 1997, 2009 Free Software Foundation, Inc.
+ This file is part of the GNU C Library.
+ Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997.
+
+ The GNU C 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 2.1 of the License, or (at your option) any later version.
+
+ The GNU C 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 GNU
+ Lesser General Public License for more details.
+
+ You should have received a copy of the GNU Lesser General Public
+ License along with the GNU C Library; if not, write to the Free
+ Software Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA
+ 02111-1307 USA. */
+
+#include <math.h>
+
+#include "math_private.h"
+
+
+static const double huge = 1.0e300;
+
+
+double
+__round (double x)
+{
+ int64_t i0, j0;
+
+ EXTRACT_WORDS64 (i0, x);
+ j0 = ((i0 >> 52) & 0x7ff) - 0x3ff;
+ if (__builtin_expect (j0 < 52, 1))
+ {
+ if (j0 < 0)
+ {
+ if (huge + x > 0.0)
+ {
+ i0 &= UINT64_C(0x8000000000000000);
+ if (j0 == -1)
+ i0 |= UINT64_C(0x3ff0000000000000);
+ }
+ }
+ else
+ {
+ uint64_t i = UINT64_C(0x000fffffffffffff) >> j0;
+ if ((i0 & i) == 0)
+ /* X is integral. */
+ return x;
+ if (huge + x > 0.0)
+ {
+ /* Raise inexact if x != 0. */
+ i0 += UINT64_C(0x0008000000000000) >> j0;
+ i0 &= ~i;
+ }
+ }
+ }
+ else
+ {
+ if (j0 == 0x400)
+ /* Inf or NaN. */
+ return x + x;
+ else
+ return x;
+ }
+
+ INSERT_WORDS64 (x, i0);
+ return x;
+}
+weak_alias (__round, round)
+#ifdef NO_LONG_DOUBLE
+strong_alias (__round, __roundl)
+weak_alias (__round, roundl)
+#endif
diff --git a/sysdeps/ieee754/dbl-64/wordsize-64/s_scalbn.c b/sysdeps/ieee754/dbl-64/wordsize-64/s_scalbn.c
new file mode 100644
index 0000000000..25cf3b1699
--- /dev/null
+++ b/sysdeps/ieee754/dbl-64/wordsize-64/s_scalbn.c
@@ -0,0 +1,68 @@
+/* @(#)s_scalbn.c 5.1 93/09/24 */
+/*
+ * ====================================================
+ * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
+ *
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+/*
+ * scalbn (double x, int n)
+ * scalbn(x,n) returns x* 2**n computed by exponent
+ * manipulation rather than by actually performing an
+ * exponentiation or a multiplication.
+ */
+
+#include "math.h"
+#include "math_private.h"
+
+#ifdef __STDC__
+static const double
+#else
+static double
+#endif
+two54 = 1.80143985094819840000e+16, /* 0x43500000, 0x00000000 */
+twom54 = 5.55111512312578270212e-17, /* 0x3C900000, 0x00000000 */
+huge = 1.0e+300,
+tiny = 1.0e-300;
+
+#ifdef __STDC__
+ double __scalbn (double x, int n)
+#else
+ double __scalbn (x,n)
+ double x; int n;
+#endif
+{
+ int64_t ix;
+ int64_t k;
+ EXTRACT_WORDS64(ix,x);
+ k = (ix >> 52) & 0x7ff; /* extract exponent */
+ if (k==0) { /* 0 or subnormal x */
+ if ((ix & UINT64_C(0xfffffffffffff))==0) return x; /* +-0 */
+ x *= two54;
+ EXTRACT_WORDS64(ix,x);
+ k = ((ix >> 52) & 0x7ff) - 54;
+ }
+ if (k==0x7ff) return x+x; /* NaN or Inf */
+ k = k+n;
+ if (n> 50000 || k > 0x7fe)
+ return huge*__copysign(huge,x); /* overflow */
+ if (n< -50000) return tiny*__copysign(tiny,x); /*underflow*/
+ if (k > 0) /* normal result */
+ {INSERT_WORDS64(x,(ix&UINT64_C(0x800fffffffffffff))|(k<<52));
+ return x;}
+ if (k <= -54)
+ return tiny*__copysign(tiny,x); /*underflow*/
+ k += 54; /* subnormal result */
+ INSERT_WORDS64(x,(ix&INT64_C(0x800fffffffffffff))|(k<<52));
+ return x*twom54;
+}
+weak_alias (__scalbn, scalbn)
+#ifdef NO_LONG_DOUBLE
+strong_alias (__scalbn, __scalbnl)
+weak_alias (__scalbn, scalbnl)
+#endif
diff --git a/sysdeps/ieee754/dbl-64/wordsize-64/s_trunc.c b/sysdeps/ieee754/dbl-64/wordsize-64/s_trunc.c
new file mode 100644
index 0000000000..9add5ade71
--- /dev/null
+++ b/sysdeps/ieee754/dbl-64/wordsize-64/s_trunc.c
@@ -0,0 +1,56 @@
+/* Truncate argument to nearest integral value not larger than the argument.
+ Copyright (C) 1997, 1998, 2009 Free Software Foundation, Inc.
+ This file is part of the GNU C Library.
+ Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997.
+
+ The GNU C 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 2.1 of the License, or (at your option) any later version.
+
+ The GNU C 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 GNU
+ Lesser General Public License for more details.
+
+ You should have received a copy of the GNU Lesser General Public
+ License along with the GNU C Library; if not, write to the Free
+ Software Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA
+ 02111-1307 USA. */
+
+#include <math.h>
+
+#include "math_private.h"
+
+
+double
+__trunc (double x)
+{
+ int64_t i0, j0;
+ int64_t sx;
+
+ EXTRACT_WORDS64 (i0, x);
+ sx = i0 & UINT64_C(0x8000000000000000);
+ j0 = ((i0 >> 52) & 0x7ff) - 0x3ff;
+ if (j0 < 52)
+ {
+ if (j0 < 0)
+ /* The magnitude of the number is < 1 so the result is +-0. */
+ INSERT_WORDS64 (x, sx);
+ else
+ INSERT_WORDS64 (x, sx | (i0 & ~(UINT64_C(0x000fffffffffffff) >> j0)));
+ }
+ else
+ {
+ if (j0 == 0x400)
+ /* x is inf or NaN. */
+ return x + x;
+ }
+
+ return x;
+}
+weak_alias (__trunc, trunc)
+#ifdef NO_LONG_DOUBLE
+strong_alias (__trunc, __truncl)
+weak_alias (__trunc, truncl)
+#endif
generated by cgit v1.2.3 (git 2.39.1) at 2024-05-23 05:16:51 +0000