X-Git-Url: http://erislabs.net/gitweb/?a=blobdiff_plain;f=lib%2Fceil.c;h=fbcca33a4a4dfcbdae30e57ae4d759c0552d278d;hb=fa1db0dd22768f09a507674a30beb5b8a87bb35f;hp=2538bf64d412720275960a58995f09d5b3bdc0fb;hpb=57fdfd3f8ec62b105c53bcdf6f127c35c7fe7391;p=gnulib.git diff --git a/lib/ceil.c b/lib/ceil.c index 2538bf64d..fbcca33a4 100644 --- a/lib/ceil.c +++ b/lib/ceil.c @@ -1,5 +1,5 @@ /* Round towards positive infinity. - Copyright (C) 2007 Free Software Foundation, Inc. + Copyright (C) 2007, 2010-2013 Free Software Foundation, Inc. This program is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by @@ -16,32 +16,52 @@ /* Written by Bruno Haible , 2007. */ -#include +#if ! defined USE_LONG_DOUBLE +# include +#endif /* Specification. */ #include #include +#undef MIN + #ifdef USE_LONG_DOUBLE # define FUNC ceill # define DOUBLE long double # define MANT_DIG LDBL_MANT_DIG +# define MIN LDBL_MIN # define L_(literal) literal##L #elif ! defined USE_FLOAT # define FUNC ceil # define DOUBLE double # define MANT_DIG DBL_MANT_DIG +# define MIN DBL_MIN # define L_(literal) literal #else /* defined USE_FLOAT */ # define FUNC ceilf # define DOUBLE float # define MANT_DIG FLT_MANT_DIG +# define MIN FLT_MIN # define L_(literal) literal##f #endif +/* -0.0. See minus-zero.h. */ +#if defined __hpux || defined __sgi || defined __ICC +# define MINUS_ZERO (-MIN * MIN) +#else +# define MINUS_ZERO L_(-0.0) +#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 + /* 2^(MANT_DIG-1). */ -static const double TWO_MANT_DIG = +static const DOUBLE TWO_MANT_DIG = /* Assume MANT_DIG <= 5 * 31. Use the identity n = floor(n/5) + floor((n+1)/5) + ... + floor((n+4)/5). */ @@ -63,20 +83,35 @@ FUNC (DOUBLE x) volatile DOUBLE y = x; volatile DOUBLE z = y; - /* Round to the next integer (nearest or up or down, doesn't matter). */ if (z > L_(0.0)) { - z += TWO_MANT_DIG; - z -= TWO_MANT_DIG; + /* Avoid rounding errors for values near 2^k, where k >= MANT_DIG-1. */ + if (z < TWO_MANT_DIG) + { + /* Round to the next integer (nearest or up or down, doesn't matter). */ + z += TWO_MANT_DIG; + z -= TWO_MANT_DIG; + /* Enforce rounding up. */ + if (z < y) + z += L_(1.0); + } } else if (z < L_(0.0)) { - z -= TWO_MANT_DIG; - z += TWO_MANT_DIG; + /* For -1 < x < 0, return -0.0 regardless of the current rounding + mode. */ + if (z > L_(-1.0)) + z = MINUS_ZERO; + /* Avoid rounding errors for values near -2^k, where k >= MANT_DIG-1. */ + else if (z > - TWO_MANT_DIG) + { + /* Round to the next integer (nearest or up or down, doesn't matter). */ + z -= TWO_MANT_DIG; + z += TWO_MANT_DIG; + /* Enforce rounding up. */ + if (z < y) + z += L_(1.0); + } } - /* Enforce rounding up. */ - if (z < y) - z += L_(1.0); - return z; }