sysdeps/ieee754/ldbl-128ibm/s_llrintl.c

   1 /* Round to long long int long double floating-point values.
   2    IBM extended format long double version.
   3    Copyright (C) 2006-2014 Free Software Foundation, Inc.
   4    This file is part of the GNU C Library.
   5
   6    The GNU C Library is free software; you can redistribute it and/or
   7    modify it under the terms of the GNU Lesser General Public
   8    License as published by the Free Software Foundation; either
   9    version 2.1 of the License, or (at your option) any later version.
  10
  11    The GNU C Library is distributed in the hope that it will be useful,
  12    but WITHOUT ANY WARRANTY; without even the implied warranty of
  13    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  14    Lesser General Public License for more details.
  15
  16    You should have received a copy of the GNU Lesser General Public
  17    License along with the GNU C Library; if not, see
  18    <http://www.gnu.org/licenses/>.  */
  19
  20 #include <math.h>
  21 #include <fenv_libc.h>
  22 #include <math_ldbl_opt.h>
  23 #include <float.h>
  24 #include <ieee754.h>
  25
  26
  27 long long
  28 __llrintl (long double x)
  29 {
  30   double xh, xl;
  31   long long res, hi, lo;
  32   int save_round;
  33
  34   ldbl_unpack (x, &xh, &xl);
  35
  36   /* Limit the range of values handled by the conversion to long long.
  37      We do this because we aren't sure whether that conversion properly
  38      raises FE_INVALID.  */
  39   if (__builtin_expect
  40       ((__builtin_fabs (xh) <= -(double) (-__LONG_LONG_MAX__ - 1)), 1)
  41 #if !defined (FE_INVALID)
  42       || 1
  43 #endif
  44     )
  45     {
  46       save_round = __fegetround ();
  47
  48       if (__builtin_expect ((xh == -(double) (-__LONG_LONG_MAX__ - 1)), 0))
  49         {
  50           /* When XH is 9223372036854775808.0, converting to long long will
  51              overflow, resulting in an invalid operation.  However, XL might
  52              be negative and of sufficient magnitude that the overall long
  53              double is in fact in range.  Avoid raising an exception.  In any
  54              case we need to convert this value specially, because
  55              the converted value is not exactly represented as a double
  56              thus subtracting HI from XH suffers rounding error.  */
  57           hi = __LONG_LONG_MAX__;
  58           xh = 1.0;
  59         }
  60       else
  61         {
  62           hi = (long long) xh;
  63           xh -= hi;
  64         }
  65       ldbl_canonicalize (&xh, &xl);
  66
  67       lo = (long long) xh;
  68
  69       /* Peg at max/min values, assuming that the above conversions do so.
  70          Strictly speaking, we can return anything for values that overflow,
  71          but this is more useful.  */
  72       res = hi + lo;
  73
  74       /* This is just sign(hi) == sign(lo) && sign(res) != sign(hi).  */
  75       if (__builtin_expect (((~(hi ^ lo) & (res ^ hi)) < 0), 0))
  76         goto overflow;
  77
  78       xh -= lo;
  79       ldbl_canonicalize (&xh, &xl);
  80
  81       hi = res;
  82       switch (save_round)
  83         {
  84         case FE_TONEAREST:
  85           if (fabs (xh) < 0.5
  86               || (fabs (xh) == 0.5
  87                   && ((xh > 0.0 && xl < 0.0)
  88                       || (xh < 0.0 && xl > 0.0)
  89                       || (xl == 0.0 && (res & 1) == 0))))
  90             return res;
  91
  92           if (xh < 0.0)
  93             res -= 1;
  94           else
  95             res += 1;
  96           break;
  97
  98         case FE_TOWARDZERO:
  99           if (res > 0 && (xh < 0.0 || (xh == 0.0 && xl < 0.0)))
 100             res -= 1;
 101           else if (res < 0 && (xh > 0.0 || (xh == 0.0 && xl > 0.0)))
 102             res += 1;
 103           return res;
 104           break;
 105
 106         case FE_UPWARD:
 107           if (xh > 0.0 || (xh == 0.0 && xl > 0.0))
 108             res += 1;
 109           break;
 110
 111         case FE_DOWNWARD:
 112           if (xh < 0.0 || (xh == 0.0 && xl < 0.0))
 113             res -= 1;
 114           break;
 115         }
 116
 117       if (__builtin_expect (((~(hi ^ (res - hi)) & (res ^ hi)) < 0), 0))
 118         goto overflow;
 119
 120       return res;
 121     }
 122   else
 123     {
 124       if (xh > 0.0)
 125         hi = __LONG_LONG_MAX__;
 126       else if (xh < 0.0)
 127         hi = -__LONG_LONG_MAX__ - 1;
 128       else
 129         /* Nan */
 130         hi = 0;
 131     }
 132
 133 overflow:
 134 #ifdef FE_INVALID
 135   feraiseexcept (FE_INVALID);
 136 #endif
 137   return hi;
 138 }
 139
 140 long_double_symbol (libm, __llrintl, llrintl);