+2007-03-21 Bruno Haible <bruno@clisp.org>
+
+ * modules/frexp: New file.
+ * lib/frexp.c: New file.
+ * m4/frexp.m4: New file.
+ * lib/math_.h (frexp): New declaration.
+ * m4/math_h.m4 (gl_MATH_H_DEFAULTS): Also initialize GNULIB_FREXP and
+ REPLACE_FREXP.
+ * modules/math (math.h): Also substitute GNULIB_FREXP, REPLACE_FREXP.
+
2007-03-21 Bruno Haible <bruno@clisp.org>
* modules/isnanl-tests: New file.
--- /dev/null
+/* Split a double into fraction and mantissa.
+ Copyright (C) 2007 Free Software Foundation, Inc.
+
+ This program is free software; you can redistribute it and/or modify
+ it under the terms of the GNU General Public License as published by
+ the Free Software Foundation; either version 2, or (at your option)
+ any later version.
+
+ This program 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 General Public License for more details.
+
+ You should have received a copy of the GNU General Public License along
+ with this program; if not, write to the Free Software Foundation,
+ Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA. */
+
+#include <config.h>
+
+#if !(defined USE_LONG_DOUBLE && !HAVE_LONG_DOUBLE)
+
+/* Specification. */
+# include <math.h>
+
+# include <float.h>
+# ifdef USE_LONG_DOUBLE
+# include "isnanl-nolibm.h"
+# else
+# include "isnan.h"
+# endif
+
+/* This file assumes FLT_RADIX = 2. If FLT_RADIX is a power of 2 greater
+ than 2, or not even a power of 2, some rounding errors can occur, so that
+ then the returned mantissa is only guaranteed to be <= 1.0, not < 1.0. */
+
+# ifdef USE_LONG_DOUBLE
+# define FUNC frexpl
+# define DOUBLE long double
+# define ISNAN isnanl
+# define L_(literal) literal##L
+# else
+# define FUNC frexp
+# define DOUBLE double
+# define ISNAN isnan
+# define L_(literal) literal
+# endif
+
+DOUBLE
+FUNC (DOUBLE x, int *exp)
+{
+ int sign;
+ int exponent;
+
+ /* Test for NaN, infinity, and zero. */
+ if (ISNAN (x) || x + x == x)
+ {
+ *exp = 0;
+ return x;
+ }
+
+ sign = 0;
+ if (x < 0)
+ {
+ x = - x;
+ sign = -1;
+ }
+
+ if (0)
+ {
+ /* Implementation contributed by Paolo Bonzini.
+ Disabled because it's under GPL and doesn't pass the tests. */
+
+ /* Since the exponent is an 'int', it fits in 64 bits. Therefore the
+ loops are executed no more than 64 times. */
+ DOUBLE exponents[64];
+ DOUBLE *next;
+ int bit;
+
+ exponent = 0;
+ if (x >= L_(1.0))
+ {
+ for (next = exponents, exponents[0] = L_(2.0), bit = 1;
+ *next <= x + x;
+ bit <<= 1, next[1] = next[0] * next[0], next++);
+
+ for (; next >= exponents; bit >>= 1, next--)
+ if (x + x >= *next)
+ {
+ x /= *next;
+ exponent |= bit;
+ }
+ }
+ else if (x < L_(0.5))
+ {
+ for (next = exponents, exponents[0] = L_(0.5), bit = 1;
+ *next > x;
+ bit <<= 1, next[1] = next[0] * next[0], next++);
+
+ for (; next >= exponents; bit >>= 1, next--)
+ if (x < *next)
+ {
+ x /= *next;
+ exponent |= bit;
+ }
+ exponent = - exponent;
+ }
+ }
+ else
+ {
+ /* Implementation contributed by Bruno Haible. */
+
+ /* Since the exponent is an 'int', it fits in 64 bits. Therefore the
+ loops are executed no more than 64 times. */
+ DOUBLE pow2[64]; /* pow2[i] = 2^2^i */
+ DOUBLE powh[64]; /* powh[i] = 2^-2^i */
+ int i;
+
+ exponent = 0;
+ if (x >= L_(1.0))
+ {
+ /* A positive exponent. */
+ DOUBLE pow2_i; /* = pow2[i] */
+ DOUBLE powh_i; /* = powh[i] */
+
+ /* Invariants: pow2_i = 2^2^i, powh_i = 2^-2^i,
+ x * 2^exponent = argument, x >= 1.0. */
+ for (i = 0, pow2_i = L_(2.0), powh_i = L_(0.5);
+ ;
+ i++, pow2_i = pow2_i * pow2_i, powh_i = powh_i * powh_i)
+ {
+ if (x >= pow2_i)
+ {
+ exponent += (1 << i);
+ x *= powh_i;
+ }
+ else
+ break;
+
+ pow2[i] = pow2_i;
+ powh[i] = powh_i;
+ }
+ /* Avoid making x too small, as it could become a denormalized
+ number and thus lose precision. */
+ while (i > 0 && x < pow2[i - 1])
+ {
+ i--;
+ powh_i = powh[i];
+ }
+ exponent += (1 << i);
+ x *= powh_i;
+ /* Here 2^-2^i <= x < 1.0. */
+ }
+ else
+ {
+ /* A negative or zero exponent. */
+ DOUBLE pow2_i; /* = pow2[i] */
+ DOUBLE powh_i; /* = powh[i] */
+
+ /* Invariants: pow2_i = 2^2^i, powh_i = 2^-2^i,
+ x * 2^exponent = argument, x < 1.0. */
+ for (i = 0, pow2_i = L_(2.0), powh_i = L_(0.5);
+ ;
+ i++, pow2_i = pow2_i * pow2_i, powh_i = powh_i * powh_i)
+ {
+ if (x < powh_i)
+ {
+ exponent -= (1 << i);
+ x *= pow2_i;
+ }
+ else
+ break;
+
+ pow2[i] = pow2_i;
+ powh[i] = powh_i;
+ }
+ /* Here 2^-2^i <= x < 1.0. */
+ }
+
+ /* Invariants: x * 2^exponent = argument, and 2^-2^i <= x < 1.0. */
+ while (i > 0)
+ {
+ i--;
+ if (x < powh[i])
+ {
+ exponent -= (1 << i);
+ x *= pow2[i];
+ }
+ }
+ /* Here 0.5 <= x < 1.0. */
+ }
+
+ *exp = exponent;
+ return (sign < 0 ? - x : x);
+}
+
+#else
+
+/* This declaration is solely to ensure that after preprocessing
+ this file is never empty. */
+typedef int dummy;
+
+#endif
#endif
+/* Write x as
+ x = mantissa * 2^exp
+ where
+ If x finite and nonzero: 0.5 <= |mantissa| < 1.0.
+ If x is zero: mantissa = x, exp = 0.
+ If x is infinite or NaN: mantissa = x, exp unspecified.
+ Store exp and return mantissa. */
+#if @GNULIB_FREXP@
+# if @REPLACE_FREXP@
+# define frexp rpl_frexp
+extern double frexp (double x, int *exp);
+# endif
+#elif defined GNULIB_POSIXCHECK
+# undef frexp
+# define frexp(x,e) \
+ (GL_LINK_WARNING ("frexp is unportable - " \
+ "use gnulib module frexp for portability"), \
+ frexp (x, e))
+#endif
+
+
#if @GNULIB_MATHL@ || !@HAVE_DECL_ACOSL@
extern long double acosl (long double x);
#endif
--- /dev/null
+# frexp.m4 serial 1
+dnl Copyright (C) 2007 Free Software Foundation, Inc.
+dnl This file is free software; the Free Software Foundation
+dnl gives unlimited permission to copy and/or distribute it,
+dnl with or without modifications, as long as this notice is preserved.
+
+AC_DEFUN([gl_FUNC_FREXP],
+[
+ AC_REQUIRE([gl_MATH_H_DEFAULTS])
+ FREXP_LIBM=
+ AC_CACHE_CHECK([whether frexp() can be used without linking with libm],
+ [gl_cv_func_frexp_no_libm],
+ [
+ AC_TRY_LINK([#include <math.h>
+ long double x;],
+ [int e; return frexp (x, &e) > 0;],
+ [gl_cv_func_frexp_no_libm=yes],
+ [gl_cv_func_frexp_no_libm=no])
+ ])
+ if test $gl_cv_func_frexp_no_libm = no; then
+ AC_CACHE_CHECK([whether frexp() can be used with libm],
+ [gl_cv_func_frexp_in_libm],
+ [
+ save_LIBS="$LIBS"
+ LIBS="$LIBS -lm"
+ AC_TRY_LINK([#include <math.h>
+ long double x;],
+ [int e; return frexp (x, &e) > 0;],
+ [gl_cv_func_frexp_in_libm=yes],
+ [gl_cv_func_frexp_in_libm=no])
+ LIBS="$save_LIBS"
+ ])
+ if test $gl_cv_func_frexp_in_libm = yes; then
+ FREXP_LIBM=-lm
+ fi
+ fi
+ if test $gl_cv_func_frexp_no_libm = yes \
+ || test $gl_cv_func_frexp_in_libm = yes; then
+ save_LIBS="$LIBS"
+ LIBS="$LIBS $FREXP_LIBM"
+ gl_FUNC_FREXP_WORKS
+ LIBS="$save_LIBS"
+ case "$gl_cv_func_frexp_works" in
+ *yes) gl_func_frexp=yes ;;
+ *) gl_func_frexp=no; REPLACE_FREXP=1; FREXP_LIBM= ;;
+ esac
+ else
+ gl_func_frexp=no
+ fi
+ if test $gl_func_frexp = yes; then
+ AC_DEFINE([HAVE_FREXP], 1,
+ [Define if the frexp() function is available and works.])
+ else
+ AC_LIBOBJ([frexp])
+ fi
+])
+
+dnl Test whether frexp() works also on denormalized numbers.
+dnl This test fails e.g. on NetBSD.
+AC_DEFUN([gl_FUNC_FREXP_WORKS],
+[
+ AC_REQUIRE([AC_PROG_CC])
+ AC_REQUIRE([AC_CANONICAL_HOST]) dnl for cross-compiles
+ AC_CACHE_CHECK([whether frexp works], [gl_cv_func_frexp_works],
+ [
+ AC_TRY_RUN([
+#include <float.h>
+#include <math.h>
+int main()
+{
+ int i;
+ volatile double x;
+ for (i = 1, x = 1.0; i >= DBL_MIN_EXP; i--, x *= 0.5)
+ ;
+ if (x > 0.0)
+ {
+ int exp;
+ double y = frexp (x, &exp);
+ /* On machines with IEEE754 arithmetic: x = 1.11254e-308, exp = -1022.
+ On NetBSD: y = 0.75. Correct: y = 0.5. */
+ if (y != 0.5)
+ return 1;
+ }
+ return 0;
+}], [gl_cv_func_frexp_works=yes], [gl_cv_func_frexp_works=no],
+ [case "$host_os" in
+ netbsd*) gl_cv_func_frexp_works="guessing no";;
+ *) gl_cv_func_frexp_works="guessing yes";;
+ esac
+ ])
+ ])
+])
-# math_h.m4 serial 2
+# math_h.m4 serial 3
dnl Copyright (C) 2007 Free Software Foundation, Inc.
dnl This file is free software; the Free Software Foundation
dnl gives unlimited permission to copy and/or distribute it,
AC_DEFUN([gl_MATH_H_DEFAULTS],
[
+ GNULIB_FREXP=0; AC_SUBST([GNULIB_FREXP])
GNULIB_MATHL=0; AC_SUBST([GNULIB_MATHL])
dnl Assume proper GNU behavior unless another module says otherwise.
HAVE_DECL_ACOSL=1; AC_SUBST([HAVE_DECL_ACOSL])
HAVE_DECL_SINL=1; AC_SUBST([HAVE_DECL_SINL])
HAVE_DECL_SQRTL=1; AC_SUBST([HAVE_DECL_SQRTL])
HAVE_DECL_TANL=1; AC_SUBST([HAVE_DECL_TANL])
+ REPLACE_FREXP=0; AC_SUBST([REPLACE_FREXP])
])
--- /dev/null
+Description:
+frexp() function: split a double into its constituents.
+
+Files:
+lib/frexp.c
+m4/frexp.m4
+
+Depends-on:
+math
+isnan-nolibm
+
+configure.ac:
+gl_FUNC_FREXP
+gl_MATH_MODULE_INDICATOR([frexp])
+
+Makefile.am:
+
+Include:
+<math.h>
+
+License:
+LGPL
+
+Maintainer:
+Bruno Haible, Paolo Bonzini
+
rm -f $@-t $@
{ echo '/* DO NOT EDIT! GENERATED AUTOMATICALLY! */' && \
sed -e 's|@''ABSOLUTE_MATH_H''@|$(ABSOLUTE_MATH_H)|g' \
+ -e 's|@''GNULIB_FREXP''@|$(GNULIB_FREXP)|g' \
-e 's|@''GNULIB_MATHL''@|$(GNULIB_MATHL)|g' \
-e 's|@''HAVE_DECL_ACOSL''@|$(HAVE_DECL_ACOSL)|g' \
-e 's|@''HAVE_DECL_ASINL''@|$(HAVE_DECL_ASINL)|g' \
-e 's|@''HAVE_DECL_SINL''@|$(HAVE_DECL_SINL)|g' \
-e 's|@''HAVE_DECL_SQRTL''@|$(HAVE_DECL_SQRTL)|g' \
-e 's|@''HAVE_DECL_TANL''@|$(HAVE_DECL_TANL)|g' \
+ -e 's|@''REPLACE_FREXP''@|$(REPLACE_FREXP)|g' \
-e '/definition of GL_LINK_WARNING/r $(LINK_WARNING_H)' \
< $(srcdir)/math_.h; \
} > $@-t
projects / pspp / commitdiff