2 dnl Copyright (C) 2011-2017 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.
7 AC_DEFUN([gl_FUNC_FMAF],
9 AC_REQUIRE([gl_MATH_H_DEFAULTS])
11 dnl Persuade glibc <math.h> to declare fmaf().
12 AC_REQUIRE([gl_USE_SYSTEM_EXTENSIONS])
14 dnl Determine FMAF_LIBM.
15 gl_MATHFUNC([fmaf], [float], [(float, float, float)],
20 float fmaf (float, float, float);
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
30 case "$gl_cv_func_fmaf_works" in
31 *no) REPLACE_FMAF=1 ;;
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])
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" ;;
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" ;;
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" ;;
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],
65 AC_REQUIRE([AC_PROG_CC])
66 AC_REQUIRE([AC_CANONICAL_HOST]) dnl for cross-compiles
67 AC_REQUIRE([gl_FUNC_LDEXPF])
69 LIBS="$LIBS $FMAF_LIBM $LDEXPF_LIBM"
70 AC_CACHE_CHECK([whether fmaf works], [gl_cv_func_fmaf_works],
80 /* These tests fail with glibc 2.11.3 on x86_64. */
82 volatile float x = 1.5f; /* 3 * 2^-1 */
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)
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)
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)
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)
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)
146 if ((FLT_MANT_DIG % 2) == 0)
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)
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))
172 [gl_cv_func_fmaf_works=yes],
173 [gl_cv_func_fmaf_works=no],
174 [dnl Guess yes on native Windows with MSVC.
175 dnl Otherwise guess no, even on glibc systems.
176 gl_cv_func_fmaf_works="guessing no"
179 AC_EGREP_CPP([Known], [
183 ], [gl_cv_func_fmaf_works="guessing yes"])
191 # Prerequisites of lib/fmaf.c.
192 AC_DEFUN([gl_PREREQ_FMAF], [:])