Mercurial > hg > octave-lyh > gnulib-hg
changeset 16565:0badbd4d62e4
fmod, fmodl: Fix computation for large quotients x / y.
* lib/fmod.c: Completely rewritten.
* lib/fmodl.c (fmodl): Use implementation of fmod.c with
USE_LONG_DOUBLE.
* modules/fmod (Depends-on): Add isfinite, signbit, fabs, frexp, ldexp,
isnand. Remove fma.
* modules/fmodl (Depends-on): Add float, isfinite, signbit, fabsl,
frexpl, ldexpl, isnanl. Remove fma.
* m4/fmod.m4 (gl_FUNC_FMOD): Update computation of FMOD_LIBM.
* m4/fmodl.m4 (gl_FUNC_FMODL): Update computation of FMODL_LIBM.
author | Bruno Haible <bruno@clisp.org> |
---|---|
date | Sun, 04 Mar 2012 20:56:32 +0100 |
parents | a5e3d6ba8dee |
children | bd4785c10d56 |
files | ChangeLog lib/fmod.c lib/fmodl.c m4/fmod.m4 m4/fmodl.m4 modules/fmod modules/fmodl |
diffstat | 7 files changed, 410 insertions(+), 111 deletions(-) [+] |
line wrap: on
line diff
--- a/ChangeLog +++ b/ChangeLog @@ -1,3 +1,16 @@ +2012-03-04 Bruno Haible <bruno@clisp.org> + + fmod, fmodl: Fix computation for large quotients x / y. + * lib/fmod.c: Completely rewritten. + * lib/fmodl.c (fmodl): Use implementation of fmod.c with + USE_LONG_DOUBLE. + * modules/fmod (Depends-on): Add isfinite, signbit, fabs, frexp, ldexp, + isnand. Remove fma. + * modules/fmodl (Depends-on): Add float, isfinite, signbit, fabsl, + frexpl, ldexpl, isnanl. Remove fma. + * m4/fmod.m4 (gl_FUNC_FMOD): Update computation of FMOD_LIBM. + * m4/fmodl.m4 (gl_FUNC_FMODL): Update computation of FMODL_LIBM. + 2012-03-03 Bruno Haible <bruno@clisp.org> fmod* tests: More tests.
--- a/lib/fmod.c +++ b/lib/fmod.c @@ -14,59 +14,347 @@ You should have received a copy of the GNU General Public License along with this program. If not, see <http://www.gnu.org/licenses/>. */ -#include <config.h> +#if ! defined USE_LONG_DOUBLE +# include <config.h> +#endif /* Specification. */ #include <math.h> -double -fmod (double x, double y) +#include <float.h> +#include <stdlib.h> + +#ifdef USE_LONG_DOUBLE +# define FMOD fmodl +# define DOUBLE long double +# define MANT_DIG LDBL_MANT_DIG +# define L_(literal) literal##L +# define FREXP frexpl +# define LDEXP ldexpl +# define FABS fabsl +# define TRUNC truncl +# define ISNAN isnanl +#else +# define FMOD fmod +# define DOUBLE double +# define MANT_DIG DBL_MANT_DIG +# define L_(literal) literal +# define FREXP frexp +# define LDEXP ldexp +# define FABS fabs +# define TRUNC trunc +# define ISNAN isnand +#endif + +/* MSVC with option -fp:strict refuses to compile constant initializers that + contain floating-point operations. Pacify this compiler. */ +#ifdef _MSC_VER +# pragma fenv_access (off) +#endif + +#undef NAN +#if defined _MSC_VER +static DOUBLE zero; +# define NAN (zero / zero) +#else +# define NAN (L_(0.0) / L_(0.0)) +#endif + +/* To avoid rounding errors during the computation of x - q * y, + there are three possibilities: + - Use fma (- q, y, x). + - Have q be a single bit at a time, and compute x - q * y + through a subtraction. + - Have q be at most MANT_DIG/2 bits at a time, and compute + x - q * y by splitting y into two halves: + y = y1 * 2^(MANT_DIG/2) + y0 + x - q * y = (x - q * y1 * 2^MANT_DIG/2) - q * y0. + The latter is probably the most efficient. */ + +/* Number of bits in a limb. */ +#define LIMB_BITS (MANT_DIG / 2) + +/* 2^LIMB_BITS. */ +static const DOUBLE TWO_LIMB_BITS = + /* Assume LIMB_BITS <= 3 * 31. + Use the identity + n = floor(n/3) + floor((n+1)/3) + floor((n+2)/3). */ + (DOUBLE) (1U << (LIMB_BITS / 3)) + * (DOUBLE) (1U << ((LIMB_BITS + 1) / 3)) + * (DOUBLE) (1U << ((LIMB_BITS + 2) / 3)); + +/* 2^-LIMB_BITS. */ +static const DOUBLE TWO_LIMB_BITS_INVERSE = + /* Assume LIMB_BITS <= 3 * 31. + Use the identity + n = floor(n/3) + floor((n+1)/3) + floor((n+2)/3). */ + L_(1.0) / ((DOUBLE) (1U << (LIMB_BITS / 3)) + * (DOUBLE) (1U << ((LIMB_BITS + 1) / 3)) + * (DOUBLE) (1U << ((LIMB_BITS + 2) / 3))); + +DOUBLE +FMOD (DOUBLE x, DOUBLE y) { - if (isinf (y)) - return x; + if (isfinite (x) && isfinite (y) && y != L_(0.0)) + { + if (x == L_(0.0)) + /* Return x, regardless of the sign of y. */ + return x; + + { + int negate = ((!signbit (x)) ^ (!signbit (y))); + + /* Take the absolute value of x and y. */ + x = FABS (x); + y = FABS (y); + + /* Trivial case that requires no computation. */ + if (x < y) + return (negate ? - x : x); + + { + int yexp; + DOUBLE ym; + DOUBLE y1; + DOUBLE y0; + int k; + DOUBLE x2; + DOUBLE x1; + DOUBLE x0; + + /* Write y = 2^yexp * (y1 * 2^-LIMB_BITS + y0 * 2^-(2*LIMB_BITS)) + where y1 is an integer, 2^(LIMB_BITS-1) <= y1 < 2^LIMB_BITS, + y1 has at most LIMB_BITS bits, + 0 <= y0 < 2^LIMB_BITS, + y0 has at most (MANT_DIG + 1) / 2 bits. */ + ym = FREXP (y, &yexp); + ym = ym * TWO_LIMB_BITS; + y1 = TRUNC (ym); + y0 = (ym - y1) * TWO_LIMB_BITS; + + /* Write + x = 2^(yexp+(k-3)*LIMB_BITS) + * (x2 * 2^(2*LIMB_BITS) + x1 * 2^LIMB_BITS + x0) + where x2, x1, x0 are each integers >= 0, < 2^LIMB_BITS. */ + { + int xexp; + DOUBLE xm = FREXP (x, &xexp); + /* Since we know x >= y, we know xexp >= yexp. */ + xexp -= yexp; + /* Compute k = ceil(xexp / LIMB_BITS). */ + k = (xexp + LIMB_BITS - 1) / LIMB_BITS; + /* Note: (k - 1) * LIMB_BITS + 1 <= xexp <= k * LIMB_BITS. */ + /* Note: 0.5 <= xm < 1.0. */ + xm = LDEXP (xm, xexp - (k - 1) * LIMB_BITS); + /* Note: Now xm < 2^(xexp - (k - 1) * LIMB_BITS) <= 2^LIMB_BITS + and xm >= 0.5 * 2^(xexp - (k - 1) * LIMB_BITS) >= 1.0 + and xm has at most MANT_DIG <= 2*LIMB_BITS+1 bits. */ + x2 = TRUNC (xm); + x1 = (xm - x2) * TWO_LIMB_BITS; + /* Split off x0 from x1 later. */ + } + + /* Test whether [x2,x1,0] >= 2^LIMB_BITS * [y1,y0]. */ + if (x2 > y1 || (x2 == y1 && x1 >= y0)) + { + /* Subtract 2^LIMB_BITS * [y1,y0] from [x2,x1,0]. */ + x2 -= y1; + x1 -= y0; + if (x1 < L_(0.0)) + { + if (!(x2 >= L_(1.0))) + abort (); + x2 -= L_(1.0); + x1 += TWO_LIMB_BITS; + } + } + + /* Split off x0 from x1. */ + { + DOUBLE x1int = TRUNC (x1); + x0 = TRUNC ((x1 - x1int) * TWO_LIMB_BITS); + x1 = x1int; + } + + for (; k > 0; k--) + { + /* Multiprecision division of the limb sequence [x2,x1,x0] + by [y1,y0]. */ + /* Here [x2,x1,x0] < 2^LIMB_BITS * [y1,y0]. */ + /* The first guess takes into account only [x2,x1] and [y1]. + + By Knuth's theorem, we know that + q* = min (floor ([x2,x1] / [y1]), 2^LIMB_BITS - 1) + and + q = floor ([x2,x1,x0] / [y1,y0]) + are not far away: + q* - 2 <= q <= q* + 1. + + Proof: + a) q* * y1 <= floor ([x2,x1] / [y1]) * y1 <= [x2,x1]. + Hence + [x2,x1,x0] - q* * [y1,y0] + = 2^LIMB_BITS * ([x2,x1] - q* * [y1]) + x0 - q* * y0 + >= x0 - q* * y0 + >= - q* * y0 + > - 2^(2*LIMB_BITS) + >= - 2 * [y1,y0] + So + [x2,x1,x0] > (q* - 2) * [y1,y0]. + b) If q* = floor ([x2,x1] / [y1]), then + [x2,x1] < (q* + 1) * y1 + Hence + [x2,x1,x0] - q* * [y1,y0] + = 2^LIMB_BITS * ([x2,x1] - q* * [y1]) + x0 - q* * y0 + <= 2^LIMB_BITS * (y1 - 1) + x0 - q* * y0 + <= 2^LIMB_BITS * (2^LIMB_BITS-2) + (2^LIMB_BITS-1) - 0 + < 2^(2*LIMB_BITS) + <= 2 * [y1,y0] + So + [x2,x1,x0] < (q* + 2) * [y1,y0]. + and so + q < q* + 2 + which implies + q <= q* + 1. + In the other case, q* = 2^LIMB_BITS - 1. Then trivially + q < 2^LIMB_BITS = q* + 1. + + We know that floor ([x2,x1] / [y1]) >= 2^LIMB_BITS if and + only if x2 >= y1. */ + DOUBLE q = + (x2 >= y1 + ? TWO_LIMB_BITS - L_(1.0) + : TRUNC ((x2 * TWO_LIMB_BITS + x1) / y1)); + if (q > L_(0.0)) + { + /* Compute + [x2,x1,x0] - q* * [y1,y0] + = 2^LIMB_BITS * ([x2,x1] - q* * [y1]) + x0 - q* * y0. */ + DOUBLE q_y1 = q * y1; /* exact, at most 2*LIMB_BITS bits */ + DOUBLE q_y1_1 = TRUNC (q_y1 * TWO_LIMB_BITS_INVERSE); + DOUBLE q_y1_0 = q_y1 - q_y1_1 * TWO_LIMB_BITS; + DOUBLE q_y0 = q * y0; /* exact, at most MANT_DIG bits */ + DOUBLE q_y0_1 = TRUNC (q_y0 * TWO_LIMB_BITS_INVERSE); + DOUBLE q_y0_0 = q_y0 - q_y0_1 * TWO_LIMB_BITS; + x2 -= q_y1_1; + x1 -= q_y1_0; + x1 -= q_y0_1; + x0 -= q_y0_0; + /* Move negative carry from x0 to x1 and from x1 to x2. */ + if (x0 < L_(0.0)) + { + x0 += TWO_LIMB_BITS; + x1 -= L_(1.0); + } + if (x1 < L_(0.0)) + { + x1 += TWO_LIMB_BITS; + x2 -= L_(1.0); + if (x1 < L_(0.0)) /* not sure this can happen */ + { + x1 += TWO_LIMB_BITS; + x2 -= L_(1.0); + } + } + if (x2 < L_(0.0)) + { + /* Reduce q by 1. */ + x1 += y1; + x0 += y0; + /* Move overflow from x0 to x1 and from x1 to x0. */ + if (x0 >= TWO_LIMB_BITS) + { + x0 -= TWO_LIMB_BITS; + x1 += L_(1.0); + } + if (x1 >= TWO_LIMB_BITS) + { + x1 -= TWO_LIMB_BITS; + x2 += L_(1.0); + } + if (x2 < L_(0.0)) + { + /* Reduce q by 1 again. */ + x1 += y1; + x0 += y0; + /* Move overflow from x0 to x1 and from x1 to x0. */ + if (x0 >= TWO_LIMB_BITS) + { + x0 -= TWO_LIMB_BITS; + x1 += L_(1.0); + } + if (x1 >= TWO_LIMB_BITS) + { + x1 -= TWO_LIMB_BITS; + x2 += L_(1.0); + } + if (x2 < L_(0.0)) + /* Shouldn't happen, because we proved that + q >= q* - 2. */ + abort (); + } + } + } + if (x2 > L_(0.0) + || x1 > y1 + || (x1 == y1 && x0 >= y0)) + { + /* Increase q by 1. */ + x1 -= y1; + x0 -= y0; + /* Move negative carry from x0 to x1 and from x1 to x2. */ + if (x0 < L_(0.0)) + { + x0 += TWO_LIMB_BITS; + x1 -= L_(1.0); + } + if (x1 < L_(0.0)) + { + x1 += TWO_LIMB_BITS; + x2 -= L_(1.0); + } + if (x2 < L_(0.0)) + abort (); + if (x2 > L_(0.0) + || x1 > y1 + || (x1 == y1 && x0 >= y0)) + /* Shouldn't happen, because we proved that + q <= q* + 1. */ + abort (); + } + /* Here [x2,x1,x0] < [y1,y0]. */ + /* Next round. */ + x2 = x1; +#if (MANT_DIG + 1) / 2 > LIMB_BITS /* y0 can have a fractional bit */ + x1 = TRUNC (x0); + x0 = (x0 - x1) * TWO_LIMB_BITS; +#else + x1 = x0; + x0 = L_(0.0); +#endif + /* Here [x2,x1,x0] < 2^LIMB_BITS * [y1,y0]. */ + } + /* Here k = 0. + The result is + 2^(yexp-3*LIMB_BITS) + * (x2 * 2^(2*LIMB_BITS) + x1 * 2^LIMB_BITS + x0). */ + { + DOUBLE r = + LDEXP ((x2 * TWO_LIMB_BITS + x1) * TWO_LIMB_BITS + x0, + yexp - 3 * LIMB_BITS); + return (negate ? - r : r); + } + } + } + } else { - double q = - trunc (x / y); - double r = fma (q, y, x); /* = x + q * y, computed in one step */ - /* Correct possible rounding errors in the quotient x / y. */ - if (y >= 0) - { - if (x >= 0) - { - /* Expect 0 <= r < y. */ - if (r < 0) - q += 1, r = fma (q, y, x); - else if (r >= y) - q -= 1, r = fma (q, y, x); - } - else - { - /* Expect - y < r <= 0. */ - if (r > 0) - q -= 1, r = fma (q, y, x); - else if (r <= - y) - q += 1, r = fma (q, y, x); - } - } + if (ISNAN (x) || ISNAN (y)) + return x + y; /* NaN */ + else if (isinf (y)) + return x; else - { - if (x >= 0) - { - /* Expect 0 <= r < - y. */ - if (r < 0) - q -= 1, r = fma (q, y, x); - else if (r >= - y) - q += 1, r = fma (q, y, x); - } - else - { - /* Expect y < r <= 0. */ - if (r > 0) - q += 1, r = fma (q, y, x); - else if (r <= y) - q -= 1, r = fma (q, y, x); - } - } - return r; + /* x infinite or y zero */ + return NAN; } }
--- a/lib/fmodl.c +++ b/lib/fmodl.c @@ -29,56 +29,7 @@ #else -long double -fmodl (long double x, long double y) -{ - if (isinf (y)) - return x; - else - { - long double q = - truncl (x / y); - long double r = fmal (q, y, x); /* = x + q * y, computed in one step */ - /* Correct possible rounding errors in the quotient x / y. */ - if (y >= 0) - { - if (x >= 0) - { - /* Expect 0 <= r < y. */ - if (r < 0) - q += 1, r = fmal (q, y, x); - else if (r >= y) - q -= 1, r = fmal (q, y, x); - } - else - { - /* Expect - y < r <= 0. */ - if (r > 0) - q -= 1, r = fmal (q, y, x); - else if (r <= - y) - q += 1, r = fmal (q, y, x); - } - } - else - { - if (x >= 0) - { - /* Expect 0 <= r < - y. */ - if (r < 0) - q -= 1, r = fmal (q, y, x); - else if (r >= - y) - q += 1, r = fmal (q, y, x); - } - else - { - /* Expect y < r <= 0. */ - if (r > 0) - q += 1, r = fmal (q, y, x); - else if (r <= y) - q -= 1, r = fmal (q, y, x); - } - } - return r; - } -} +# define USE_LONG_DOUBLE +# include "fmod.c" #endif
--- a/m4/fmod.m4 +++ b/m4/fmod.m4 @@ -1,4 +1,4 @@ -# fmod.m4 serial 2 +# fmod.m4 serial 3 dnl Copyright (C) 2011-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, @@ -67,18 +67,36 @@ ]) if test $REPLACE_FMOD = 1; then dnl Find libraries needed to link lib/fmod.c. + AC_REQUIRE([gl_FUNC_FABS]) + AC_REQUIRE([gl_FUNC_FREXP]) AC_REQUIRE([gl_FUNC_TRUNC]) - AC_REQUIRE([gl_FUNC_FMA]) + AC_REQUIRE([gl_FUNC_LDEXP]) + AC_REQUIRE([gl_FUNC_ISNAND]) FMOD_LIBM= + dnl Append $FABS_LIBM to FMOD_LIBM, avoiding gratuitous duplicates. + case " $FMOD_LIBM " in + *" $FABS_LIBM "*) ;; + *) FMOD_LIBM="$FMOD_LIBM $FABS_LIBM" ;; + esac + dnl Append $FREXP_LIBM to FMOD_LIBM, avoiding gratuitous duplicates. + case " $FMOD_LIBM " in + *" $FREXP_LIBM "*) ;; + *) FMOD_LIBM="$FMOD_LIBM $FREXP_LIBM" ;; + esac dnl Append $TRUNC_LIBM to FMOD_LIBM, avoiding gratuitous duplicates. case " $FMOD_LIBM " in *" $TRUNC_LIBM "*) ;; *) FMOD_LIBM="$FMOD_LIBM $TRUNC_LIBM" ;; esac - dnl Append $FMA_LIBM to FMOD_LIBM, avoiding gratuitous duplicates. + dnl Append $LDEXP_LIBM to FMOD_LIBM, avoiding gratuitous duplicates. case " $FMOD_LIBM " in - *" $FMA_LIBM "*) ;; - *) FMOD_LIBM="$FMOD_LIBM $FMA_LIBM" ;; + *" $LDEXP_LIBM "*) ;; + *) FMOD_LIBM="$FMOD_LIBM $LDEXP_LIBM" ;; + esac + dnl Append $ISNAND_LIBM to FMOD_LIBM, avoiding gratuitous duplicates. + case " $FMOD_LIBM " in + *" $ISNAND_LIBM "*) ;; + *) FMOD_LIBM="$FMOD_LIBM $ISNAND_LIBM" ;; esac fi ])
--- a/m4/fmodl.m4 +++ b/m4/fmodl.m4 @@ -1,4 +1,4 @@ -# fmodl.m4 serial 2 +# fmodl.m4 serial 3 dnl Copyright (C) 2011-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, @@ -82,18 +82,36 @@ if test $HAVE_SAME_LONG_DOUBLE_AS_DOUBLE = 1; then FMODL_LIBM="$FMOD_LIBM" else + AC_REQUIRE([gl_FUNC_FABSL]) + AC_REQUIRE([gl_FUNC_FREXPL]) AC_REQUIRE([gl_FUNC_TRUNCL]) - AC_REQUIRE([gl_FUNC_FMAL]) + AC_REQUIRE([gl_FUNC_LDEXPL]) + AC_REQUIRE([gl_FUNC_ISNANL]) FMODL_LIBM= + dnl Append $FABSL_LIBM to FMODL_LIBM, avoiding gratuitous duplicates. + case " $FMODL_LIBM " in + *" $FABSL_LIBM "*) ;; + *) FMODL_LIBM="$FMODL_LIBM $FABSL_LIBM" ;; + esac + dnl Append $FREXPL_LIBM to FMODL_LIBM, avoiding gratuitous duplicates. + case " $FMODL_LIBM " in + *" $FREXPL_LIBM "*) ;; + *) FMODL_LIBM="$FMODL_LIBM $FREXPL_LIBM" ;; + esac dnl Append $TRUNCL_LIBM to FMODL_LIBM, avoiding gratuitous duplicates. case " $FMODL_LIBM " in *" $TRUNCL_LIBM "*) ;; *) FMODL_LIBM="$FMODL_LIBM $TRUNCL_LIBM" ;; esac - dnl Append $FMAL_LIBM to FMODL_LIBM, avoiding gratuitous duplicates. + dnl Append $LDEXPL_LIBM to FMODL_LIBM, avoiding gratuitous duplicates. case " $FMODL_LIBM " in - *" $FMAL_LIBM "*) ;; - *) FMODL_LIBM="$FMODL_LIBM $FMAL_LIBM" ;; + *" $LDEXPL_LIBM "*) ;; + *) FMODL_LIBM="$FMODL_LIBM $LDEXPL_LIBM" ;; + esac + dnl Append $ISNANL_LIBM to FMODL_LIBM, avoiding gratuitous duplicates. + case " $FMODL_LIBM " in + *" $ISNANL_LIBM "*) ;; + *) FMODL_LIBM="$FMODL_LIBM $ISNANL_LIBM" ;; esac fi fi
--- a/modules/fmod +++ b/modules/fmod @@ -8,9 +8,14 @@ Depends-on: math -isinf [test $REPLACE_FMOD = 1] +isfinite [test $REPLACE_FMOD = 1] +signbit [test $REPLACE_FMOD = 1] +fabs [test $REPLACE_FMOD = 1] +frexp [test $REPLACE_FMOD = 1] trunc [test $REPLACE_FMOD = 1] -fma [test $REPLACE_FMOD = 1] +ldexp [test $REPLACE_FMOD = 1] +isnand [test $REPLACE_FMOD = 1] +isinf [test $REPLACE_FMOD = 1] configure.ac: gl_FUNC_FMOD
--- a/modules/fmodl +++ b/modules/fmodl @@ -9,9 +9,15 @@ Depends-on: math fmod [{ test $HAVE_FMODL = 0 || test $REPLACE_FMODL = 1; } && test $HAVE_SAME_LONG_DOUBLE_AS_DOUBLE = 1] -isinf [{ test $HAVE_FMODL = 0 || test $REPLACE_FMODL = 1; } && test $HAVE_SAME_LONG_DOUBLE_AS_DOUBLE = 0] +float [{ test $HAVE_FMODL = 0 || test $REPLACE_FMODL = 1; } && test $HAVE_SAME_LONG_DOUBLE_AS_DOUBLE = 0] +isfinite [{ test $HAVE_FMODL = 0 || test $REPLACE_FMODL = 1; } && test $HAVE_SAME_LONG_DOUBLE_AS_DOUBLE = 0] +signbit [{ test $HAVE_FMODL = 0 || test $REPLACE_FMODL = 1; } && test $HAVE_SAME_LONG_DOUBLE_AS_DOUBLE = 0] +fabsl [{ test $HAVE_FMODL = 0 || test $REPLACE_FMODL = 1; } && test $HAVE_SAME_LONG_DOUBLE_AS_DOUBLE = 0] +frexpl [{ test $HAVE_FMODL = 0 || test $REPLACE_FMODL = 1; } && test $HAVE_SAME_LONG_DOUBLE_AS_DOUBLE = 0] truncl [{ test $HAVE_FMODL = 0 || test $REPLACE_FMODL = 1; } && test $HAVE_SAME_LONG_DOUBLE_AS_DOUBLE = 0] -fmal [{ test $HAVE_FMODL = 0 || test $REPLACE_FMODL = 1; } && test $HAVE_SAME_LONG_DOUBLE_AS_DOUBLE = 0] +ldexpl [{ test $HAVE_FMODL = 0 || test $REPLACE_FMODL = 1; } && test $HAVE_SAME_LONG_DOUBLE_AS_DOUBLE = 0] +isnanl [{ test $HAVE_FMODL = 0 || test $REPLACE_FMODL = 1; } && test $HAVE_SAME_LONG_DOUBLE_AS_DOUBLE = 0] +isinf [{ test $HAVE_FMODL = 0 || test $REPLACE_FMODL = 1; } && test $HAVE_SAME_LONG_DOUBLE_AS_DOUBLE = 0] configure.ac: gl_FUNC_FMODL