m4/fmal.m4

   1 # fmal.m4 serial 1
   2 dnl Copyright (C) 2011 Free Software Foundation, Inc.
   3 dnl This file is free software; the Free Software Foundation
   4 dnl gives unlimited permission to copy and/or distribute it,
   5 dnl with or without modifications, as long as this notice is preserved.
   6
   7 AC_DEFUN([gl_FUNC_FMAL],
   8 [
   9   AC_REQUIRE([gl_MATH_H_DEFAULTS])
  10   AC_REQUIRE([gl_LONG_DOUBLE_VS_DOUBLE])
  11
  12   dnl Determine FMAL_LIBM.
  13   gl_MATHFUNC([fmal], [long double], [(long double, long double, long double)])
  14   if test $gl_cv_func_fmal_no_libm = yes \
  15      || test $gl_cv_func_fmal_in_libm = yes; then
  16     gl_FUNC_FMAL_WORKS
  17     case "$gl_cv_func_fmal_works" in
  18       *no) REPLACE_FMAL=1 ;;
  19     esac
  20   else
  21     HAVE_FMAL=0
  22   fi
  23   if test $HAVE_FMAL = 0 || test $REPLACE_FMAL = 1; then
  24     dnl Find libraries needed to link lib/fmal.c.
  25     if test $HAVE_SAME_LONG_DOUBLE_AS_DOUBLE = 1; then
  26       AC_REQUIRE([gl_FUNC_FMA])
  27       FMAL_LIBM="$FMA_LIBM"
  28     else
  29       AC_REQUIRE([gl_FUNC_FREXPL])
  30       AC_REQUIRE([gl_FUNC_LDEXPL])
  31       AC_REQUIRE([gl_FUNC_FEGETROUND])
  32       FMAL_LIBM=
  33       dnl Append $FREXPL_LIBM to FMAL_LIBM, avoiding gratuitous duplicates.
  34       case " $FMAL_LIBM " in
  35         *" $FREXPL_LIBM "*) ;;
  36         *) FMAL_LIBM="$FMAL_LIBM $FREXPL_LIBM" ;;
  37       esac
  38       dnl Append $LDEXPL_LIBM to FMAL_LIBM, avoiding gratuitous duplicates.
  39       case " $FMAL_LIBM " in
  40         *" $LDEXPL_LIBM "*) ;;
  41         *) FMAL_LIBM="$FMAL_LIBM $LDEXPL_LIBM" ;;
  42       esac
  43       dnl Append $FEGETROUND_LIBM to FMAL_LIBM, avoiding gratuitous duplicates.
  44       case " $FMAL_LIBM " in
  45         *" $FEGETROUND_LIBM "*) ;;
  46         *) FMAL_LIBM="$FMAL_LIBM $FEGETROUND_LIBM" ;;
  47       esac
  48     fi
  49   fi
  50   AC_SUBST([FMAL_LIBM])
  51 ])
  52
  53 dnl Test whether fmal() has any of the 15 known bugs of glibc 2.11.3 on x86_64.
  54 AC_DEFUN([gl_FUNC_FMAL_WORKS],
  55 [
  56   AC_REQUIRE([AC_PROG_CC])
  57   AC_REQUIRE([AC_CANONICAL_HOST]) dnl for cross-compiles
  58   AC_REQUIRE([gl_FUNC_LDEXPL])
  59   save_LIBS="$LIBS"
  60   LIBS="$LIBS $FMAL_LIBM $LDEXPL_LIBM"
  61   AC_CACHE_CHECK([whether fmal works], [gl_cv_func_fmal_works],
  62     [
  63       AC_RUN_IFELSE(
  64         [AC_LANG_SOURCE([[
  65 #include <float.h>
  66 #include <math.h>
  67 /* Override the values of <float.h>, like done in float.in.h.  */
  68 #if defined __i386__ && (defined __BEOS__ || defined __OpenBSD__)
  69 # undef LDBL_MANT_DIG
  70 # define LDBL_MANT_DIG   64
  71 # undef LDBL_MIN_EXP
  72 # define LDBL_MIN_EXP    (-16381)
  73 # undef LDBL_MAX_EXP
  74 # define LDBL_MAX_EXP    16384
  75 #endif
  76 #if defined __i386__ && defined __FreeBSD__
  77 # undef LDBL_MANT_DIG
  78 # define LDBL_MANT_DIG   64
  79 # undef LDBL_MIN_EXP
  80 # define LDBL_MIN_EXP    (-16381)
  81 # undef LDBL_MAX_EXP
  82 # define LDBL_MAX_EXP    16384
  83 #endif
  84 #if (defined _ARCH_PPC || defined _POWER) && defined _AIX && (LDBL_MANT_DIG == 106) && defined __GNUC__
  85 # undef LDBL_MIN_EXP
  86 # define LDBL_MIN_EXP DBL_MIN_EXP
  87 #endif
  88 #if defined __sgi && (LDBL_MANT_DIG >= 106)
  89 # undef LDBL_MANT_DIG
  90 # define LDBL_MANT_DIG 106
  91 # if defined __GNUC__
  92 #  undef LDBL_MIN_EXP
  93 #  define LDBL_MIN_EXP DBL_MIN_EXP
  94 # endif
  95 #endif
  96 long double p0 = 0.0L;
  97 int main()
  98 {
  99   int failed_tests = 0;
 100   /* This test fails on glibc 2.11 powerpc.  */
 101   {
 102     volatile long double x = 1.5L; /* 3 * 2^-1 */
 103     volatile long double y = x;
 104     volatile long double z = ldexpl (1.0L, LDBL_MANT_DIG + 1); /* 2^65 */
 105     /* x * y + z with infinite precision: 2^65 + 9 * 2^-2.
 106        Lies between (2^63 + 0) * 2^2 and (2^63 + 1) * 2^2
 107        and is closer to (2^63 + 1) * 2^2, therefore the rounding
 108        must round up and produce (2^63 + 1) * 2^2.  */
 109     volatile long double expected = z + 4.0L;
 110     volatile long double result = fmal (x, y, z);
 111     if (result != expected)
 112       failed_tests |= 1;
 113   }
 114   /* This test fails on glibc 2.11 powerpc.  */
 115   {
 116     volatile long double x = 1.25L; /* 2^0 + 2^-2 */
 117     volatile long double y = - x;
 118     volatile long double z = ldexpl (1.0L, LDBL_MANT_DIG + 1); /* 2^65 */
 119     /* x * y + z with infinite precision: 2^65 - 2^0 - 2^-1 - 2^-4.
 120        Lies between (2^64 - 1) * 2^1 and 2^64 * 2^1
 121        and is closer to (2^64 - 1) * 2^1, therefore the rounding
 122        must round down and produce (2^64 - 1) * 2^1.  */
 123     volatile long double expected = (ldexpl (1.0L, LDBL_MANT_DIG) - 1.0L) * 2.0L;
 124     volatile long double result = fmal (x, y, z);
 125     if (result != expected)
 126       failed_tests |= 1;
 127   }
 128   /* This test fails on glibc 2.11 x86,x86_64,powerpc, glibc 2.7 hppa,sparc,
 129      OSF/1 5.1, mingw.  */
 130   {
 131     volatile long double x = 1.0L + ldexpl (1.0L, 1 - LDBL_MANT_DIG); /* 2^0 + 2^-63 */
 132     volatile long double y = x;
 133     volatile long double z = 4.0L; /* 2^2 */
 134     /* x * y + z with infinite precision: 2^2 + 2^0 + 2^-62 + 2^-126.
 135        Lies between (2^63 + 2^61) * 2^-61 and (2^63 + 2^61 + 1) * 2^-61
 136        and is closer to (2^63 + 2^61 + 1) * 2^-61, therefore the rounding
 137        must round up and produce (2^63 + 2^61 + 1) * 2^-61.  */
 138     volatile long double expected = 4.0L + 1.0L + ldexpl (1.0L, 3 - LDBL_MANT_DIG);
 139     volatile long double result = fmal (x, y, z);
 140     if (result != expected)
 141       failed_tests |= 2;
 142   }
 143   /* This test fails on glibc 2.11 x86,x86_64,powerpc glibc 2.7 hppa,sparc,
 144      OSF/1 5.1, mingw.  */
 145   {
 146     volatile long double x = 1.0L + ldexpl (1.0L, 1 - LDBL_MANT_DIG); /* 2^0 + 2^-63 */
 147     volatile long double y = - x;
 148     volatile long double z = 8.0L; /* 2^3 */
 149     /* x * y + z with infinite precision: 2^2 + 2^1 + 2^0 - 2^-62 - 2^-126.
 150        Lies between (2^63 + 2^62 + 2^61 - 1) * 2^-61 and
 151        (2^63 + 2^62 + 2^61) * 2^-61 and is closer to
 152        (2^63 + 2^62 + 2^61 - 1) * 2^-61, therefore the rounding
 153        must round down and produce (2^63 + 2^62 + 2^61 - 1) * 2^-61.  */
 154     volatile long double expected = 7.0L - ldexpl (1.0L, 3 - LDBL_MANT_DIG);
 155     volatile long double result = fmal (x, y, z);
 156     if (result != expected)
 157       failed_tests |= 2;
 158   }
 159   /* This test fails on glibc 2.11 powerpc.  */
 160   {
 161     volatile long double x = 1.25L; /* 2^0 + 2^-2 */
 162     volatile long double y = - 0.75L; /* - 2^0 + 2^-2 */
 163     volatile long double z = ldexpl (1.0L, LDBL_MANT_DIG); /* 2^64 */
 164     /* x * y + z with infinite precision: 2^64 - 2^0 + 2^-4.
 165        Lies between (2^64 - 2^0) and 2^64 and is closer to (2^64 - 2^0),
 166        therefore the rounding must round down and produce (2^64 - 2^0).  */
 167     volatile long double expected = ldexpl (1.0L, LDBL_MANT_DIG) - 1.0L;
 168     volatile long double result = fmal (x, y, z);
 169     if (result != expected)
 170       failed_tests |= 1;
 171   }
 172   if ((LDBL_MANT_DIG % 2) == 1)
 173     {
 174       /* These tests fail on glibc 2.7 hppa,sparc, OSF/1 5.1.  */
 175       {
 176         volatile long double x = 1.0L + ldexpl (1.0L, - (LDBL_MANT_DIG + 1) / 2); /* 2^0 + 2^-27 */
 177         volatile long double y = 1.0L - ldexpl (1.0L, - (LDBL_MANT_DIG + 1) / 2); /* 2^0 - 2^-27 */
 178         volatile long double z = - ldexpl (1.0L, LDBL_MIN_EXP - LDBL_MANT_DIG); /* - 2^-1074 */
 179         /* x * y + z with infinite precision: 2^0 - 2^-54 - 2^-1074.
 180            Lies between (2^53 - 1) * 2^-53 and 2^53 * 2^-53 and is closer to
 181            (2^53 - 1) * 2^-53, therefore the rounding must round down and
 182            produce (2^53 - 1) * 2^-53.  */
 183         volatile long double expected = 1.0L - ldexpl (1.0L, - LDBL_MANT_DIG);
 184         volatile long double result = fmal (x, y, z);
 185         if (result != expected)
 186           failed_tests |= 4;
 187       }
 188       {
 189         volatile long double x = 1.0L + ldexpl (1.0L, - (LDBL_MANT_DIG + 1) / 2); /* 2^0 + 2^-57 */
 190         volatile long double y = x;
 191         volatile long double z = ldexpl (1.0L, - LDBL_MANT_DIG); /* 2^-113 */
 192         /* x * y + z with infinite precision: 2^0 + 2^-56 + 2^-113 + 2^-114.
 193            Lies between (2^112 + 2^56) * 2^-112 and (2^112 + 2^56 + 1) * 2^-112
 194            and is closer to (2^112 + 2^56 + 1) * 2^-112, therefore the rounding
 195            must round up and produce (2^112 + 2^56 + 1) * 2^-112.  */
 196         volatile long double expected =
 197           1.0L + ldexpl (1.0L, - (LDBL_MANT_DIG - 1) / 2) + ldexpl (1.0L, 1 - LDBL_MANT_DIG);
 198         volatile long double result = fmal (x, y, z);
 199         if (result != expected)
 200           failed_tests |= 4;
 201       }
 202     }
 203   else
 204     {
 205       /* These tests fail on glibc 2.11 x86,x86_64,powerpc, mingw.  */
 206       {
 207         volatile long double x = 1.0L + ldexpl (1.0L, - LDBL_MANT_DIG / 2); /* 2^0 + 2^-32 */
 208         volatile long double y = x;
 209         volatile long double z = ldexpl (1.0L, LDBL_MIN_EXP - LDBL_MANT_DIG); /* 2^-16445 */
 210         /* x * y + z with infinite precision: 2^0 + 2^-31 + 2^-64 + 2^-16445.
 211            Lies between (2^63 + 2^32 + 0) * 2^-63 and (2^63 + 2^32 + 1) * 2^-63
 212            and is closer to (2^63 + 2^32 + 1) * 2^-63, therefore the rounding
 213            must round up and produce (2^63 + 2^32 + 1) * 2^-63.  */
 214         volatile long double expected =
 215           1.0L + ldexpl (1.0L, 1 - LDBL_MANT_DIG / 2) + ldexpl (1.0L, 1 - LDBL_MANT_DIG);
 216         volatile long double result = fmal (x, y, z);
 217         if (result != expected)
 218           failed_tests |= 8;
 219       }
 220       {
 221         volatile long double x = 1.0L + ldexpl (1.0L, - LDBL_MANT_DIG / 2); /* 2^0 + 2^-32 */
 222         volatile long double y = x;
 223         volatile long double z = ldexpl (1.0L, - LDBL_MANT_DIG); /* 2^-64 */
 224         /* x * y + z with infinite precision: 2^0 + 2^-31 + 2^-63.
 225            Rounding must return this value unchanged.  */
 226         volatile long double expected = 1.0L + ldexpl (1.0L, 1 - LDBL_MANT_DIG / 2) + ldexpl (1.0L, 1 - LDBL_MANT_DIG);
 227         volatile long double result = fmal (x, y, z);
 228         if (result != expected)
 229           failed_tests |= 8;
 230       }
 231       {
 232         volatile long double x = 1.0L + ldexpl (1.0L, - LDBL_MANT_DIG / 2); /* 2^0 + 2^-32 */
 233         volatile long double y = x;
 234         volatile long double z = ldexpl (1.0L, 1 - LDBL_MANT_DIG); /* 2^-63 */
 235         /* x * y + z with infinite precision: 2^0 + 2^-31 + 2^-63 + 2^-64.
 236            Lies between (2^63 + 2^32 + 1) * 2^-63 and (2^63 + 2^32 + 2) * 2^-63
 237            and is at the same distance from each.  According to the round-to-even
 238            rule, the rounding must round up and produce (2^63 + 2^32 + 2) * 2^-63.  */
 239         volatile long double expected = 1.0L + ldexpl (1.0L, -31) + ldexpl (1.0L, -62);
 240         volatile long double result = fmal (x, y, z);
 241         if (result != expected)
 242           failed_tests |= 8;
 243       }
 244       {
 245         volatile long double x = 1.0L + ldexpl (1.0L, - LDBL_MANT_DIG / 2); /* 2^0 + 2^-32 */
 246         volatile long double y = x;
 247         volatile long double z = ldexpl (1.0L, LDBL_MANT_DIG / 2 + 1); /* 2^33 */
 248         /* x * y + z with infinite precision: 2^33 + 2^0 + 2^-31 + 2^-64.
 249            Lies between (2^63 + 2^30) * 2^-30 and (2^63 + 2^30 + 1) * 2^-30
 250            and is closer to (2^63 + 2^30 + 1) * 2^-30, therefore the rounding
 251            must round up and produce (2^63 + 2^30 + 1) * 2^-30.  */
 252         volatile long double expected = z + 1.0L + ldexp (1.0L, 2 - LDBL_MANT_DIG / 2);
 253         volatile long double result = fmal (x, y, z);
 254         if (result != expected)
 255           failed_tests |= 8;
 256       }
 257       {
 258         volatile long double x = 1.0L + ldexpl (1.0L, - LDBL_MANT_DIG / 2); /* 2^0 + 2^-32 */
 259         volatile long double y = x;
 260         volatile long double z = - ldexpl (1.0, 1 - LDBL_MANT_DIG); /* - 2^-63 */
 261         /* x * y + z with infinite precision: 2^0 + 2^-31 - 2^-64.
 262            Lies between (2^63 + 2^32 - 1) * 2^-63 and (2^63 + 2^32) * 2^-63
 263            and is at the same distance from each.  According to the round-to-even
 264            rule, the rounding must round up and produce (2^63 + 2^32) * 2^-63.  */
 265         volatile long double expected = 1.0L + ldexpl (1.0L, 1 - LDBL_MANT_DIG / 2);
 266         volatile long double result = fmal (x, y, z);
 267         if (result != expected)
 268           failed_tests |= 8;
 269       }
 270       {
 271         volatile long double x = 1.0L + ldexpl (1.0L, - LDBL_MANT_DIG / 2); /* 2^0 + 2^-32 */
 272         volatile long double y = x;
 273         volatile long double z = - 1.0L; /* - 2^0 */
 274         /* x * y + z with infinite precision: 2^-31 + 2^-64.
 275            Rounding must return this value unchanged.  */
 276         volatile long double expected = ldexpl (1.0L, 1 - LDBL_MANT_DIG / 2) + ldexpl (1.0L, - LDBL_MANT_DIG);
 277         volatile long double result = fmal (x, y, z);
 278         if (result != expected)
 279           failed_tests |= 8;
 280       }
 281       {
 282         volatile long double x = 1.0L + ldexpl (1.0L, - LDBL_MANT_DIG / 2); /* 2^0 + 2^-32 */
 283         volatile long double y = - x;
 284         volatile long double z = 2.0L; /* 2^1 */
 285         /* x * y + z with infinite precision: 2^0 - 2^31 - 2^-64.
 286            Rounding must return this value unchanged.  */
 287         volatile long double expected = 1.0L - ldexpl (1.0L, 1 - LDBL_MANT_DIG / 2) - ldexpl (1.0L, - LDBL_MANT_DIG);
 288         volatile long double result = fmal (x, y, z);
 289         if (result != expected)
 290           failed_tests |= 8;
 291       }
 292       {
 293         volatile long double x = 1.0L + ldexpl (1.0L, - LDBL_MANT_DIG / 2); /* 2^0 + 2^-32 */
 294         volatile long double y = - x;
 295         volatile long double z = ldexpl (1.0L, LDBL_MANT_DIG / 2 + 2); /* 2^34 */
 296         /* x * y + z with infinite precision: 2^34 - (2^0 + 2^-31 + 2^-64).
 297            Lies between (2^64 - 2^30 - 1) * 2^-30 and (2^64 - 2^30) * 2^-30
 298            and is closer to (2^64 - 2^30 - 1) * 2^-30, therefore the rounding
 299            must round down and produce (2^64 - 2^30 - 1) * 2^-30.  */
 300         volatile long double expected = z - 1.0L - ldexpl (1.0L, 2 - LDBL_MANT_DIG / 2);
 301         volatile long double result = fmal (x, y, z);
 302         if (result != expected)
 303           failed_tests |= 8;
 304       }
 305       {
 306         volatile long double x = 1.0L + ldexpl (1.0L, - LDBL_MANT_DIG / 2 - 1); /* 2^0 + 2^-33 */
 307         volatile long double y = 1.0L - ldexpl (1.0L, - LDBL_MANT_DIG / 2 - 1); /* 2^0 - 2^-33 */
 308         volatile long double z = - ldexpl (1.0L, - LDBL_MANT_DIG - 1); /* 2^-65 */
 309         /* x * y + z with infinite precision: 2^0 - 2^-65 - 2^-66.
 310            Lies between (2^64 - 1) * 2^-64 and 2^64 * 2^-64 and is closer to
 311            (2^64 - 1) * 2^-64, therefore the rounding must round down and
 312            produce (2^64 - 1) * 2^-64.  */
 313         volatile long double expected = 1.0L - ldexpl (1.0L, - LDBL_MANT_DIG);
 314         volatile long double result = fmal (x, y, z);
 315         if (result != expected)
 316           failed_tests |= 8;
 317       }
 318       {
 319         volatile long double x = 1.0L + ldexpl (1.0L, - LDBL_MANT_DIG / 2 - 1); /* 2^0 + 2^-33 */
 320         volatile long double y = 1.0L - ldexpl (1.0L, - LDBL_MANT_DIG / 2 - 1); /* 2^0 - 2^-33 */
 321         volatile long double z = - 1.0L; /* 2^0 */
 322         /* x * y + z with infinite precision: - 2^-66.
 323            Rounding must return this value unchanged.  */
 324         volatile long double expected = - ldexpl (1.0L, - LDBL_MANT_DIG - 2);
 325         volatile long double result = fmal (x, y, z);
 326         if (result != expected)
 327           failed_tests |= 8;
 328       }
 329     }
 330   /* This test fails on glibc 2.11 x86,x86_64,powerpc, glibc 2.7 hppa,sparc,
 331      FreeBSD 6.4 x86, mingw.  */
 332   {
 333     long double minus_inf = -1.0L / p0;
 334     volatile long double x = ldexpl (1.0L, LDBL_MAX_EXP - 1);
 335     volatile long double y = ldexpl (1.0L, LDBL_MAX_EXP - 1);
 336     volatile long double z = minus_inf;
 337     volatile long double result = fmal (x, y, z);
 338     if (!(result == minus_inf))
 339       failed_tests |= 16;
 340   }
 341   /* This test fails on glibc 2.11 x86,x86_64,powerpc glibc 2.7 hppa,sparc,
 342      MacOS X 10.5, FreeBSD 6.4 x86, OSF/1 5.1, mingw.  */
 343   {
 344     volatile long double x = ldexpl (1.0L, LDBL_MAX_EXP - 1);
 345     volatile long double y = 2.0L;
 346     volatile long double z =
 347       - ldexpl (ldexpl (1.0L, LDBL_MAX_EXP - 1) - ldexpl (1.0L, LDBL_MAX_EXP - LDBL_MANT_DIG - 1), 1);
 348     volatile long double expected = ldexpl (1.0L, LDBL_MAX_EXP - LDBL_MANT_DIG);
 349     volatile long double result = fmal (x, y, z);
 350     if (result != expected)
 351       failed_tests |= 32;
 352   }
 353   return failed_tests;
 354 }]])],
 355         [gl_cv_func_fmal_works=yes],
 356         [gl_cv_func_fmal_works=no],
 357         [dnl Guess no, even on glibc systems.
 358          gl_cv_func_fmal_works="guessing no"
 359         ])
 360     ])
 361   LIBS="$save_LIBS"
 362 ])
 363
 364 # Prerequisites of lib/fmal.c.
 365 AC_DEFUN([gl_PREREQ_FMAL], [:])