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

   1 /* Round to int long double floating-point values.
   2    IBM extended format long double version.
   3    Copyright (C) 2006 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, write to the Free
  18    Software Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA
  19    02111-1307 USA.  */
  20
  21 /* This has been coded in assembler because GCC makes such a mess of it
  22    when it's coded in C.  */
  23
  24 #include <math.h>
  25 #include <math_ldbl_opt.h>
  26 #include <float.h>
  27 #include <ieee754.h>
  28
  29
  30 #ifdef __STDC__
  31 long double
  32 __rintl (long double x)
  33 #else
  34 long double
  35 __rintl (x)
  36      long double x;
  37 #endif
  38 {
  39   static const long double TWO52 = 4503599627370496.0L;
  40   union ibm_extended_long_double u;
  41   u.d = x;
  42
  43   if (fabs (u.dd[0]) < TWO52)
  44     {
  45       double high = u.dd[0];
  46       if (high > 0.0)
  47         {
  48           high += TWO52;
  49           high -= TWO52;
  50           if (high == -0.0) high = 0.0;
  51         }
  52       else if (high < 0.0)
  53         {
  54           high -= TWO52;
  55           high += TWO52;
  56           if (high == 0.0) high = -0.0;
  57         }
  58       u.dd[0] = high;
  59       u.dd[1] = 0.0;
  60     }
  61   else if (fabs (u.dd[1]) < TWO52 && u.dd[1] != 0.0)
  62     {
  63       double high, low, tau;
  64       /* In this case we have to round the low double and handle any
  65          adjustment to the high double that may be caused by rounding
  66          (up).  This is complicated by the fact that the high double
  67          may already be rounded and the low double may have the
  68          opposite sign to compensate.  */
  69       if (u.dd[0] > 0.0)
  70         {
  71           if (u.dd[1] > 0.0)
  72             {
  73               /* If the high/low doubles are the same sign then simply
  74                  round the low double.  */
  75               high = u.dd[0];
  76               low = u.dd[1];
  77             }
  78           else if (u.dd[1] < 0.0)
  79             {
  80               /* Else the high double is pre rounded and we need to
  81                  adjust for that.  */
  82
  83               tau = nextafter (u.dd[0], 0.0);
  84               tau = (u.dd[0] - tau) * 2.0;
  85               high = u.dd[0] - tau;
  86               low = u.dd[1] + tau;
  87             }
  88           low += TWO52;
  89           low -= TWO52;
  90         }
  91       else if (u.dd[0] < 0.0)
  92         {
  93           if (u.dd[1] < 0.0)
  94             {
  95               /* If the high/low doubles are the same sign then simply
  96                  round the low double.  */
  97               high = u.dd[0];
  98               low = u.dd[1];
  99             }
 100           else if (u.dd[1] > 0.0)
 101             {
 102               /* Else the high double is pre rounded and we need to
 103                  adjust for that.  */
 104               tau = nextafter (u.dd[0], 0.0);
 105               tau = (u.dd[0] - tau) * 2.0;
 106               high = u.dd[0] - tau;
 107               low = u.dd[1] + tau;
 108             }
 109           low = TWO52 - low;
 110           low = -(low - TWO52);
 111         }
 112       u.dd[0] = high + low;
 113       u.dd[1] = high - u.dd[0] + low;
 114     }
 115
 116   return u.d;
 117 }
 118
 119 long_double_symbol (libm, __rintl, rintl);