m4/fmaf.m4

   1 # fmaf.m4 serial 3
   2 dnl Copyright (C) 2011-2014 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_FMAF],
   8 [
   9   AC_REQUIRE([gl_MATH_H_DEFAULTS])
  10
  11   dnl Persuade glibc <math.h> to declare fmaf().
  12   AC_REQUIRE([gl_USE_SYSTEM_EXTENSIONS])
  13
  14   dnl Determine FMAF_LIBM.
  15   gl_MATHFUNC([fmaf], [float], [(float, float, float)],
  16     [extern
  17      #ifdef __cplusplus
  18      "C"
  19      #endif
  20      float fmaf (float, float, float);
  21     ])
  22   if test $gl_cv_func_fmaf_no_libm = yes \
  23      || test $gl_cv_func_fmaf_in_libm = yes; then
  24     dnl Also check whether it's declared.
  25     dnl IRIX 6.5 has fmaf() in libm but doesn't declare it in <math.h>,
  26     dnl and the function is likely buggy.
  27     AC_CHECK_DECL([fmaf], , [REPLACE_FMAF=1], [[#include <math.h>]])
  28     if test $REPLACE_FMAF = 0; then
  29       gl_FUNC_FMAF_WORKS
  30       case "$gl_cv_func_fmaf_works" in
  31         *no) REPLACE_FMAF=1 ;;
  32       esac
  33     fi
  34   else
  35     HAVE_FMAF=0
  36   fi
  37   if test $HAVE_FMAF = 0 || test $REPLACE_FMAF = 1; then
  38     dnl Find libraries needed to link lib/fmaf.c.
  39     AC_REQUIRE([gl_FUNC_FREXPF])
  40     AC_REQUIRE([gl_FUNC_LDEXPF])
  41     AC_REQUIRE([gl_FUNC_FEGETROUND])
  42     FMAF_LIBM=
  43     dnl Append $FREXPF_LIBM to FMAF_LIBM, avoiding gratuitous duplicates.
  44     case " $FMAF_LIBM " in
  45       *" $FREXPF_LIBM "*) ;;
  46       *) FMAF_LIBM="$FMAF_LIBM $FREXPF_LIBM" ;;
  47     esac
  48     dnl Append $LDEXPF_LIBM to FMAF_LIBM, avoiding gratuitous duplicates.
  49     case " $FMAF_LIBM " in
  50       *" $LDEXPF_LIBM "*) ;;
  51       *) FMAF_LIBM="$FMAF_LIBM $LDEXPF_LIBM" ;;
  52     esac
  53     dnl Append $FEGETROUND_LIBM to FMAF_LIBM, avoiding gratuitous duplicates.
  54     case " $FMAF_LIBM " in
  55       *" $FEGETROUND_LIBM "*) ;;
  56       *) FMAF_LIBM="$FMAF_LIBM $FEGETROUND_LIBM" ;;
  57     esac
  58   fi
  59   AC_SUBST([FMAF_LIBM])
  60 ])
  61
  62 dnl Test whether fmaf() has any of the 7 known bugs of glibc 2.11.3 on x86_64.
  63 AC_DEFUN([gl_FUNC_FMAF_WORKS],
  64 [
  65   AC_REQUIRE([AC_PROG_CC])
  66   AC_REQUIRE([AC_CANONICAL_HOST]) dnl for cross-compiles
  67   AC_REQUIRE([gl_FUNC_LDEXPF])
  68   save_LIBS="$LIBS"
  69   LIBS="$LIBS $FMAF_LIBM $LDEXPF_LIBM"
  70   AC_CACHE_CHECK([whether fmaf works], [gl_cv_func_fmaf_works],
  71     [
  72       AC_RUN_IFELSE(
  73         [AC_LANG_SOURCE([[
  74 #include <float.h>
  75 #include <math.h>
  76 float p0 = 0.0f;
  77 int main()
  78 {
  79   int failed_tests = 0;
  80   /* These tests fail with glibc 2.11.3 on x86_64.  */
  81   {
  82     volatile float x = 1.5f; /* 3 * 2^-1 */
  83     volatile float y = x;
  84     volatile float z = ldexpf (1.0f, FLT_MANT_DIG + 1); /* 2^25 */
  85     /* x * y + z with infinite precision: 2^25 + 9 * 2^-2.
  86        Lies between (2^23 + 0) * 2^2 and (2^23 + 1) * 2^2
  87        and is closer to (2^23 + 1) * 2^2, therefore the rounding
  88        must round up and produce (2^23 + 1) * 2^2.  */
  89     volatile float expected = z + 4.0f;
  90     volatile float result = fmaf (x, y, z);
  91     if (result != expected)
  92       failed_tests |= 1;
  93   }
  94   {
  95     volatile float x = 1.25f; /* 2^0 + 2^-2 */
  96     volatile float y = - x;
  97     volatile float z = ldexpf (1.0f, FLT_MANT_DIG + 1); /* 2^25 */
  98     /* x * y + z with infinite precision: 2^25 - 2^0 - 2^-1 - 2^-4.
  99        Lies between (2^24 - 1) * 2^1 and 2^24 * 2^1
 100        and is closer to (2^24 - 1) * 2^1, therefore the rounding
 101        must round down and produce (2^24 - 1) * 2^1.  */
 102     volatile float expected = (ldexpf (1.0f, FLT_MANT_DIG) - 1.0f) * 2.0f;
 103     volatile float result = fmaf (x, y, z);
 104     if (result != expected)
 105       failed_tests |= 2;
 106   }
 107   {
 108     volatile float x = 1.0f + ldexpf (1.0f, 1 - FLT_MANT_DIG); /* 2^0 + 2^-23 */
 109     volatile float y = x;
 110     volatile float z = 4.0f; /* 2^2 */
 111     /* x * y + z with infinite precision: 2^2 + 2^0 + 2^-22 + 2^-46.
 112        Lies between (2^23 + 2^21) * 2^-21 and (2^23 + 2^21 + 1) * 2^-21
 113        and is closer to (2^23 + 2^21 + 1) * 2^-21, therefore the rounding
 114        must round up and produce (2^23 + 2^21 + 1) * 2^-21.  */
 115     volatile float expected = 4.0f + 1.0f + ldexpf (1.0f, 3 - FLT_MANT_DIG);
 116     volatile float result = fmaf (x, y, z);
 117     if (result != expected)
 118       failed_tests |= 4;
 119   }
 120   {
 121     volatile float x = 1.0f + ldexpf (1.0f, 1 - FLT_MANT_DIG); /* 2^0 + 2^-23 */
 122     volatile float y = - x;
 123     volatile float z = 8.0f; /* 2^3 */
 124     /* x * y + z with infinite precision: 2^2 + 2^1 + 2^0 - 2^-22 - 2^-46.
 125        Lies between (2^23 + 2^22 + 2^21 - 1) * 2^-21 and
 126        (2^23 + 2^22 + 2^21) * 2^-21 and is closer to
 127        (2^23 + 2^22 + 2^21 - 1) * 2^-21, therefore the rounding
 128        must round down and produce (2^23 + 2^22 + 2^21 - 1) * 2^-21.  */
 129     volatile float expected = 7.0f - ldexpf (1.0f, 3 - FLT_MANT_DIG);
 130     volatile float result = fmaf (x, y, z);
 131     if (result != expected)
 132       failed_tests |= 8;
 133   }
 134   {
 135     volatile float x = 1.25f; /* 2^0 + 2^-2 */
 136     volatile float y = - 0.75f; /* - 2^0 + 2^-2 */
 137     volatile float z = ldexpf (1.0f, FLT_MANT_DIG); /* 2^24 */
 138     /* x * y + z with infinite precision: 2^24 - 2^0 + 2^-4.
 139        Lies between (2^24 - 2^0) and 2^24 and is closer to (2^24 - 2^0),
 140        therefore the rounding must round down and produce (2^24 - 2^0).  */
 141     volatile float expected = ldexpf (1.0f, FLT_MANT_DIG) - 1.0f;
 142     volatile float result = fmaf (x, y, z);
 143     if (result != expected)
 144       failed_tests |= 16;
 145   }
 146   if ((FLT_MANT_DIG % 2) == 0)
 147     {
 148       volatile float x = 1.0f + ldexpf (1.0f, - FLT_MANT_DIG / 2); /* 2^0 + 2^-12 */
 149       volatile float y = x;
 150       volatile float z = ldexpf (1.0f, FLT_MIN_EXP - FLT_MANT_DIG); /* 2^-149 */
 151       /* x * y + z with infinite precision: 2^0 + 2^-11 + 2^-24 + 2^-149.
 152          Lies between (2^23 + 2^12 + 0) * 2^-23 and (2^23 + 2^12 + 1) * 2^-23
 153          and is closer to (2^23 + 2^12 + 1) * 2^-23, therefore the rounding
 154          must round up and produce (2^23 + 2^12 + 1) * 2^-23.  */
 155       volatile float expected =
 156         1.0f + ldexpf (1.0f, 1 - FLT_MANT_DIG / 2) + ldexpf (1.0f, 1 - FLT_MANT_DIG);
 157       volatile float result = fmaf (x, y, z);
 158       if (result != expected)
 159         failed_tests |= 32;
 160     }
 161   {
 162     float minus_inf = -1.0f / p0;
 163     volatile float x = ldexpf (1.0f, FLT_MAX_EXP - 1);
 164     volatile float y = ldexpf (1.0f, FLT_MAX_EXP - 1);
 165     volatile float z = minus_inf;
 166     volatile float result = fmaf (x, y, z);
 167     if (!(result == minus_inf))
 168       failed_tests |= 64;
 169   }
 170   return failed_tests;
 171 }]])],
 172         [gl_cv_func_fmaf_works=yes],
 173         [gl_cv_func_fmaf_works=no],
 174         [dnl Guess no, even on glibc systems.
 175          gl_cv_func_fmaf_works="guessing no"
 176         ])
 177     ])
 178   LIBS="$save_LIBS"
 179 ])
 180
 181 # Prerequisites of lib/fmaf.c.
 182 AC_DEFUN([gl_PREREQ_FMAF], [:])