Mercurial > hg > octave-lojdl > gnulib-hg

changeset 16697:6b65a767165c

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
sqrtl: Bypass broken implementation in OpenBSD 5.1/SPARC. * lib/math.in.h (sqrtl): Replace it if REPLACE_SQRTL is 1. * m4/sqrtl.m4 (gl_FUNC_SQRTL_WORKS): New macro. (gl_FUNC_SQRTL): Invoke it. Set REPLACE_SQRTL to 1 if sqrtl() produces too big rounding errors. * m4/math_h.m4 (gl_MATH_H_DEFAULTS): Initialize REPLACE_SQRTL. * modules/math (Makefile.am): Substitute REPLACE_SQRTL. * modules/sqrtl (configure.ac): Consider REPLACE_SQRTL. (Depends-on): Update conditions. * tests/test-sqrtl.c (my_ldexpl): New function. (main): Add test of a particular value. * doc/posix-functions/sqrtl.texi: Mention the OpenBSD 5.1/SPARC bug.
author Bruno Haible <bruno@clisp.org>
date Wed, 14 Mar 2012 01:51:10 +0100
parents a61aad62576e
children bf76aa416bc7
files ChangeLog doc/posix-functions/sqrtl.texi lib/math.in.h m4/math_h.m4 m4/sqrtl.m4 modules/math modules/sqrtl tests/test-sqrtl.c
diffstat 8 files changed, 128 insertions(+), 11 deletions(-) [+]
line wrap: on
line diff
--- a/ChangeLog
+++ b/ChangeLog
@@ -1,3 +1,18 @@
+2012-03-13  Bruno Haible  <bruno@clisp.org>
+
+	sqrtl: Bypass broken implementation in OpenBSD 5.1/SPARC.
+	* lib/math.in.h (sqrtl): Replace it if REPLACE_SQRTL is 1.
+	* m4/sqrtl.m4 (gl_FUNC_SQRTL_WORKS): New macro.
+	(gl_FUNC_SQRTL): Invoke it. Set REPLACE_SQRTL to 1 if sqrtl() produces
+	too big rounding errors.
+	* m4/math_h.m4 (gl_MATH_H_DEFAULTS): Initialize REPLACE_SQRTL.
+	* modules/math (Makefile.am): Substitute REPLACE_SQRTL.
+	* modules/sqrtl (configure.ac): Consider REPLACE_SQRTL.
+	(Depends-on): Update conditions.
+	* tests/test-sqrtl.c (my_ldexpl): New function.
+	(main): Add test of a particular value.
+	* doc/posix-functions/sqrtl.texi: Mention the OpenBSD 5.1/SPARC bug.
+
 2012-03-13  Pádraig Brady  <P@draigBrady.com>
 
 	doc: Update timer_* platform portability notes.
--- a/doc/posix-functions/sqrtl.texi
+++ b/doc/posix-functions/sqrtl.texi
@@ -17,6 +17,9 @@
 @item
 This function is not declared on some platforms:
 MacOS X 10.3.
+@item
+This function produces very imprecise results on some platforms:
+OpenBSD 5.1/SPARC.
 @end itemize
 
 Portability problems not fixed by Gnulib:
--- a/lib/math.in.h
+++ b/lib/math.in.h
@@ -1698,11 +1698,20 @@
 #endif
 
 #if @GNULIB_SQRTL@
-# if !@HAVE_SQRTL@ || !@HAVE_DECL_SQRTL@
-#  undef sqrtl
+# if @REPLACE_SQRTL@
+#  if !(defined __cplusplus && defined GNULIB_NAMESPACE)
+#   undef sqrtl
+#   define sqrtl rpl_sqrtl
+#  endif
+_GL_FUNCDECL_RPL (sqrtl, long double, (long double x));
+_GL_CXXALIAS_RPL (sqrtl, long double, (long double x));
+# else
+#  if !@HAVE_SQRTL@ || !@HAVE_DECL_SQRTL@
+#   undef sqrtl
 _GL_FUNCDECL_SYS (sqrtl, long double, (long double x));
+#  endif
+_GL_CXXALIAS_SYS (sqrtl, long double, (long double x));
 # endif
-_GL_CXXALIAS_SYS (sqrtl, long double, (long double x));
 _GL_CXXALIASWARN (sqrtl);
 #elif defined GNULIB_POSIXCHECK
 # undef sqrtl
--- a/m4/math_h.m4
+++ b/m4/math_h.m4
@@ -1,4 +1,4 @@
-# math_h.m4 serial 103
+# math_h.m4 serial 104
 dnl Copyright (C) 2007-2012 Free Software Foundation, Inc.
 dnl This file is free software; the Free Software Foundation
 dnl gives unlimited permission to copy and/or distribute it,
@@ -295,6 +295,7 @@
   REPLACE_ROUNDL=0;            AC_SUBST([REPLACE_ROUNDL])
   REPLACE_SIGNBIT=0;           AC_SUBST([REPLACE_SIGNBIT])
   REPLACE_SIGNBIT_USING_GCC=0; AC_SUBST([REPLACE_SIGNBIT_USING_GCC])
+  REPLACE_SQRTL=0;             AC_SUBST([REPLACE_SQRTL])
   REPLACE_TRUNC=0;             AC_SUBST([REPLACE_TRUNC])
   REPLACE_TRUNCF=0;            AC_SUBST([REPLACE_TRUNCF])
   REPLACE_TRUNCL=0;            AC_SUBST([REPLACE_TRUNCL])
--- a/m4/sqrtl.m4
+++ b/m4/sqrtl.m4
@@ -1,4 +1,4 @@
-# sqrtl.m4 serial 7
+# sqrtl.m4 serial 8
 dnl Copyright (C) 2010-2012 Free Software Foundation, Inc.
 dnl This file is free software; the Free Software Foundation
 dnl gives unlimited permission to copy and/or distribute it,
@@ -58,9 +58,20 @@
     dnl Also check whether it's declared.
     dnl MacOS X 10.3 has sqrtl() in libc but doesn't declare it in <math.h>.
     AC_CHECK_DECL([sqrtl], , [HAVE_DECL_SQRTL=0], [[#include <math.h>]])
+
+    save_LIBS="$LIBS"
+    LIBS="$LIBS $SQRTL_LIBM"
+    gl_FUNC_SQRTL_WORKS
+    LIBS="$save_LIBS"
+    case "$gl_cv_func_sqrtl_works" in
+      *yes) ;;
+      *) REPLACE_SQRTL=1 ;;
+    esac
   else
     HAVE_DECL_SQRTL=0
     HAVE_SQRTL=0
+  fi
+  if test $HAVE_SQRTL = 0 || test $REPLACE_SQRTL = 1; then
     dnl Find libraries needed to link lib/sqrtl.c.
     if test $HAVE_SAME_LONG_DOUBLE_AS_DOUBLE = 1; then
       AC_REQUIRE([gl_FUNC_SQRT])
@@ -94,3 +105,56 @@
   fi
   AC_SUBST([SQRTL_LIBM])
 ])
+
+dnl Test whether sqrtl() works.
+dnl On OpenBSD 5.1/SPARC, sqrtl(8.1974099812331540680810141969554806865L) has
+dnl rounding errors that eat up the last 8 to 9 decimal digits.
+AC_DEFUN([gl_FUNC_SQRTL_WORKS],
+[
+  AC_REQUIRE([AC_PROG_CC])
+  AC_REQUIRE([AC_CANONICAL_HOST]) dnl for cross-compiles
+  AC_CACHE_CHECK([whether sqrtl works], [gl_cv_func_sqrtl_works],
+    [
+      AC_RUN_IFELSE(
+        [AC_LANG_SOURCE([[
+#include <float.h>
+#include <math.h>
+extern
+#ifdef __cplusplus
+"C"
+#endif
+long double sqrtl (long double);
+static long double
+my_ldexpl (long double x, int d)
+{
+  for (; d > 0; d--)
+    x *= 2.0L;
+  for (; d < 0; d++)
+    x *= 0.5L;
+  return x;
+}
+volatile long double x;
+volatile long double y;
+long double z;
+int main ()
+{
+  x = 8.1974099812331540680810141969554806865L;
+  y = sqrtl (x);
+  z = y * y - x;
+  z = my_ldexpl (z, LDBL_MANT_DIG);
+  if (z < 0)
+    z = - z;
+  if (z > 100.0L)
+    return 1;
+  return 0;
+}
+]])],
+        [gl_cv_func_sqrtl_works=yes],
+        [gl_cv_func_sqrtl_works=no],
+        [case "$host_os" in
+           osf*) gl_cv_func_sqrtl_works="guessing no";;
+           *)    gl_cv_func_sqrtl_works="guessing yes";;
+         esac
+        ])
+    ])
+])
--- a/modules/math
+++ b/modules/math
@@ -262,6 +262,7 @@
 	      -e 's|@''REPLACE_ROUNDL''@|$(REPLACE_ROUNDL)|g' \
 	      -e 's|@''REPLACE_SIGNBIT''@|$(REPLACE_SIGNBIT)|g' \
 	      -e 's|@''REPLACE_SIGNBIT_USING_GCC''@|$(REPLACE_SIGNBIT_USING_GCC)|g' \
+	      -e 's|@''REPLACE_SQRTL''@|$(REPLACE_SQRTL)|g' \
 	      -e 's|@''REPLACE_TRUNC''@|$(REPLACE_TRUNC)|g' \
 	      -e 's|@''REPLACE_TRUNCF''@|$(REPLACE_TRUNCF)|g' \
 	      -e 's|@''REPLACE_TRUNCL''@|$(REPLACE_TRUNCL)|g' \
--- a/modules/sqrtl
+++ b/modules/sqrtl
@@ -8,15 +8,15 @@
 Depends-on:
 math
 extensions
-sqrt            [test $HAVE_SQRTL = 0]
-float           [test $HAVE_SQRTL = 0 && test $HAVE_SAME_LONG_DOUBLE_AS_DOUBLE = 0]
-isnanl          [test $HAVE_SQRTL = 0 && test $HAVE_SAME_LONG_DOUBLE_AS_DOUBLE = 0]
-frexpl          [test $HAVE_SQRTL = 0 && test $HAVE_SAME_LONG_DOUBLE_AS_DOUBLE = 0]
-ldexpl          [test $HAVE_SQRTL = 0 && test $HAVE_SAME_LONG_DOUBLE_AS_DOUBLE = 0]
+sqrt            [{ test $HAVE_SQRTL = 0 || test $REPLACE_SQRTL = 1; }]
+float           [{ test $HAVE_SQRTL = 0 || test $REPLACE_SQRTL = 1; } && test $HAVE_SAME_LONG_DOUBLE_AS_DOUBLE = 0]
+isnanl          [{ test $HAVE_SQRTL = 0 || test $REPLACE_SQRTL = 1; } && test $HAVE_SAME_LONG_DOUBLE_AS_DOUBLE = 0]
+frexpl          [{ test $HAVE_SQRTL = 0 || test $REPLACE_SQRTL = 1; } && test $HAVE_SAME_LONG_DOUBLE_AS_DOUBLE = 0]
+ldexpl          [{ test $HAVE_SQRTL = 0 || test $REPLACE_SQRTL = 1; } && test $HAVE_SAME_LONG_DOUBLE_AS_DOUBLE = 0]
 
 configure.ac:
 gl_FUNC_SQRTL
-if test $HAVE_SQRTL = 0; then
+if test $HAVE_SQRTL = 0 || test $REPLACE_SQRTL = 1; then
   AC_LIBOBJ([sqrtl])
 fi
 gl_MATH_MODULE_INDICATOR([sqrtl])
--- a/tests/test-sqrtl.c
+++ b/tests/test-sqrtl.c
@@ -35,6 +35,16 @@
 #define RANDOM randoml
 #include "test-sqrt.h"
 
+static long double
+my_ldexpl (long double x, int d)
+{
+  for (; d > 0; d--)
+    x *= 2.0L;
+  for (; d < 0; d++)
+    x *= 0.5L;
+  return x;
+}
+
 int
 main ()
 {
@@ -47,6 +57,20 @@
   y = sqrtl (x);
   ASSERT (y >= 0.7745966692L && y <= 0.7745966693L);
 
+  /* Another particular value.  */
+  {
+    long double z;
+    long double err;
+
+    x = 8.1974099812331540680810141969554806865L;
+    y = sqrtl (x);
+    z = y * y - x;
+    err = my_ldexpl (z, LDBL_MANT_DIG);
+    if (err < 0)
+      err = - err;
+    ASSERT (err <= 100.0L);
+  }
+
   test_function ();
 
   return 0;