2 dnl Copyright (C) 2011-2012 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 Determine FMAF_LIBM.
12 gl_MATHFUNC([fmaf], [float], [(float, float, float)])
13 if test $gl_cv_func_fmaf_no_libm = yes \
14 || test $gl_cv_func_fmaf_in_libm = yes; then
16 case "$gl_cv_func_fmaf_works" in
17 *no) REPLACE_FMAF=1 ;;
22 if test $HAVE_FMAF = 0 || test $REPLACE_FMAF = 1; then
23 dnl Find libraries needed to link lib/fmaf.c.
24 AC_REQUIRE([gl_FUNC_FREXPF])
25 AC_REQUIRE([gl_FUNC_LDEXPF])
26 AC_REQUIRE([gl_FUNC_FEGETROUND])
28 dnl Append $FREXPF_LIBM to FMAF_LIBM, avoiding gratuitous duplicates.
29 case " $FMAF_LIBM " in
30 *" $FREXPF_LIBM "*) ;;
31 *) FMAF_LIBM="$FMAF_LIBM $FREXPF_LIBM" ;;
33 dnl Append $LDEXPF_LIBM to FMAF_LIBM, avoiding gratuitous duplicates.
34 case " $FMAF_LIBM " in
35 *" $LDEXPF_LIBM "*) ;;
36 *) FMAF_LIBM="$FMAF_LIBM $LDEXPF_LIBM" ;;
38 dnl Append $FEGETROUND_LIBM to FMAF_LIBM, avoiding gratuitous duplicates.
39 case " $FMAF_LIBM " in
40 *" $FEGETROUND_LIBM "*) ;;
41 *) FMAF_LIBM="$FMAF_LIBM $FEGETROUND_LIBM" ;;
47 dnl Test whether fmaf() has any of the 7 known bugs of glibc 2.11.3 on x86_64.
48 AC_DEFUN([gl_FUNC_FMAF_WORKS],
50 AC_REQUIRE([AC_PROG_CC])
51 AC_REQUIRE([AC_CANONICAL_HOST]) dnl for cross-compiles
52 AC_REQUIRE([gl_FUNC_LDEXPF])
54 LIBS="$LIBS $FMAF_LIBM $LDEXPF_LIBM"
55 AC_CACHE_CHECK([whether fmaf works], [gl_cv_func_fmaf_works],
65 /* These tests fail with glibc 2.11.3 on x86_64. */
67 volatile float x = 1.5f; /* 3 * 2^-1 */
69 volatile float z = ldexpf (1.0f, FLT_MANT_DIG + 1); /* 2^25 */
70 /* x * y + z with infinite precision: 2^25 + 9 * 2^-2.
71 Lies between (2^23 + 0) * 2^2 and (2^23 + 1) * 2^2
72 and is closer to (2^23 + 1) * 2^2, therefore the rounding
73 must round up and produce (2^23 + 1) * 2^2. */
74 volatile float expected = z + 4.0f;
75 volatile float result = fmaf (x, y, z);
76 if (result != expected)
80 volatile float x = 1.25f; /* 2^0 + 2^-2 */
81 volatile float y = - x;
82 volatile float z = ldexpf (1.0f, FLT_MANT_DIG + 1); /* 2^25 */
83 /* x * y + z with infinite precision: 2^25 - 2^0 - 2^-1 - 2^-4.
84 Lies between (2^24 - 1) * 2^1 and 2^24 * 2^1
85 and is closer to (2^24 - 1) * 2^1, therefore the rounding
86 must round down and produce (2^24 - 1) * 2^1. */
87 volatile float expected = (ldexpf (1.0f, FLT_MANT_DIG) - 1.0f) * 2.0f;
88 volatile float result = fmaf (x, y, z);
89 if (result != expected)
93 volatile float x = 1.0f + ldexpf (1.0f, 1 - FLT_MANT_DIG); /* 2^0 + 2^-23 */
95 volatile float z = 4.0f; /* 2^2 */
96 /* x * y + z with infinite precision: 2^2 + 2^0 + 2^-22 + 2^-46.
97 Lies between (2^23 + 2^21) * 2^-21 and (2^23 + 2^21 + 1) * 2^-21
98 and is closer to (2^23 + 2^21 + 1) * 2^-21, therefore the rounding
99 must round up and produce (2^23 + 2^21 + 1) * 2^-21. */
100 volatile float expected = 4.0f + 1.0f + ldexpf (1.0f, 3 - FLT_MANT_DIG);
101 volatile float result = fmaf (x, y, z);
102 if (result != expected)
106 volatile float x = 1.0f + ldexpf (1.0f, 1 - FLT_MANT_DIG); /* 2^0 + 2^-23 */
107 volatile float y = - x;
108 volatile float z = 8.0f; /* 2^3 */
109 /* x * y + z with infinite precision: 2^2 + 2^1 + 2^0 - 2^-22 - 2^-46.
110 Lies between (2^23 + 2^22 + 2^21 - 1) * 2^-21 and
111 (2^23 + 2^22 + 2^21) * 2^-21 and is closer to
112 (2^23 + 2^22 + 2^21 - 1) * 2^-21, therefore the rounding
113 must round down and produce (2^23 + 2^22 + 2^21 - 1) * 2^-21. */
114 volatile float expected = 7.0f - ldexpf (1.0f, 3 - FLT_MANT_DIG);
115 volatile float result = fmaf (x, y, z);
116 if (result != expected)
120 volatile float x = 1.25f; /* 2^0 + 2^-2 */
121 volatile float y = - 0.75f; /* - 2^0 + 2^-2 */
122 volatile float z = ldexpf (1.0f, FLT_MANT_DIG); /* 2^24 */
123 /* x * y + z with infinite precision: 2^24 - 2^0 + 2^-4.
124 Lies between (2^24 - 2^0) and 2^24 and is closer to (2^24 - 2^0),
125 therefore the rounding must round down and produce (2^24 - 2^0). */
126 volatile float expected = ldexpf (1.0f, FLT_MANT_DIG) - 1.0f;
127 volatile float result = fmaf (x, y, z);
128 if (result != expected)
131 if ((FLT_MANT_DIG % 2) == 0)
133 volatile float x = 1.0f + ldexpf (1.0f, - FLT_MANT_DIG / 2); /* 2^0 + 2^-12 */
134 volatile float y = x;
135 volatile float z = ldexpf (1.0f, FLT_MIN_EXP - FLT_MANT_DIG); /* 2^-149 */
136 /* x * y + z with infinite precision: 2^0 + 2^-11 + 2^-24 + 2^-149.
137 Lies between (2^23 + 2^12 + 0) * 2^-23 and (2^23 + 2^12 + 1) * 2^-23
138 and is closer to (2^23 + 2^12 + 1) * 2^-23, therefore the rounding
139 must round up and produce (2^23 + 2^12 + 1) * 2^-23. */
140 volatile float expected =
141 1.0f + ldexpf (1.0f, 1 - FLT_MANT_DIG / 2) + ldexpf (1.0f, 1 - FLT_MANT_DIG);
142 volatile float result = fmaf (x, y, z);
143 if (result != expected)
147 float minus_inf = -1.0f / p0;
148 volatile float x = ldexpf (1.0f, FLT_MAX_EXP - 1);
149 volatile float y = ldexpf (1.0f, FLT_MAX_EXP - 1);
150 volatile float z = minus_inf;
151 volatile float result = fmaf (x, y, z);
152 if (!(result == minus_inf))
157 [gl_cv_func_fmaf_works=yes],
158 [gl_cv_func_fmaf_works=no],
159 [dnl Guess no, even on glibc systems.
160 gl_cv_func_fmaf_works="guessing no"
166 # Prerequisites of lib/fmaf.c.
167 AC_DEFUN([gl_PREREQ_FMAF], [:])