/ src / ref / ceil.c
ceil.c
  1  /*
  2   * Copyright (C) 2008-2020 Advanced Micro Devices, Inc. All rights reserved.
  3   *
  4   * Redistribution and use in source and binary forms, with or without modification,
  5   * are permitted provided that the following conditions are met:
  6   * 1. Redistributions of source code must retain the above copyright notice,
  7   *    this list of conditions and the following disclaimer.
  8   * 2. Redistributions in binary form must reproduce the above copyright notice,
  9   *    this list of conditions and the following disclaimer in the documentation
 10   *    and/or other materials provided with the distribution.
 11   * 3. Neither the name of the copyright holder nor the names of its contributors
 12   *    may be used to endorse or promote products derived from this software without
 13   *    specific prior written permission.
 14   *
 15   * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
 16   * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
 17   * WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED.
 18   * IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT,
 19   * INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
 20   * BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA,
 21   * OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY,
 22   * WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
 23   * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
 24   * POSSIBILITY OF SUCH DAMAGE.
 25   *
 26   */
 27  
 28  #include "fn_macros.h"
 29  #include "libm_util_amd.h"
 30  #include "libm_special.h"
 31  
 32  double FN_PROTOTYPE_REF(ceil)(double x)
 33  {
 34    double r;
 35    long long rexp, xneg;
 36    unsigned long long ux, ax, ur, mask;
 37  
 38    GET_BITS_DP64(x, ux);
 39    /*ax is |x|*/
 40    ax = ux & (~SIGNBIT_DP64);
 41    /*xneg stores the sign of the input x*/
 42    xneg = (ux != ax);
 43    /*The range is divided into
 44      > 2^53. ceil will either the number itself or Nan
 45              always returns a QNan. Raises exception if input is a SNan
 46      < 1.0   If 0.0 then return with the appropriate sign
 47              If input is less than -0.0 and greater than -1.0 then return -0.0
 48              If input is greater than 0.0 and less than 1.0 then return 1.0
 49      1.0 < |x| < 2^53
 50              appropriately check the exponent and set the return Value by shifting
 51              */
 52    if (ax >= 0x4340000000000000) /* abs(x) > 2^53*/
 53      {
 54        /* abs(x) is either NaN, infinity, or >= 2^53 */
 55        if (ax > 0x7ff0000000000000)
 56          /* x is NaN */
 57          #ifdef WINDOWS
 58              return __amd_handle_error("ceil", __amd_ceil, ux|0x0008000000000000, _DOMAIN, 0, EDOM, x, 0.0, 1);
 59          #else
 60             {
 61             if(!(ax & 0x0008000000000000))// x is snan
 62                 return __amd_handle_error("ceil", __amd_ceil, ux|0x0008000000000000, _DOMAIN, AMD_F_INVALID, EDOM, x, 0.0, 1);
 63             else
 64                 return x;
 65  		    }
 66          #endif
 67  
 68        else
 69          return x;
 70      }
 71    else if (ax < 0x3ff0000000000000) /* abs(x) < 1.0 */
 72      {
 73        if (ax == 0x0000000000000000)
 74          /* x is +zero or -zero; return the same zero */
 75            return x;
 76        else if (xneg) /* x < 0.0; return -0.0 */
 77          {
 78            PUT_BITS_DP64(0x8000000000000000, r);
 79            return r;
 80          }
 81        else
 82          return 1.0;
 83      }
 84    else
 85      {
 86        /*Get the exponent for the floating point number. Should be between 0 and 53*/
 87        rexp = ((ux & EXPBITS_DP64) >> EXPSHIFTBITS_DP64) - EXPBIAS_DP64;
 88        /* Mask out the bits of r that we don't want */
 89        mask = 1;
 90        mask = (mask << (EXPSHIFTBITS_DP64 - rexp)) - 1;
 91        /*Keeps the exponent part and the required mantissa.*/
 92        ur = (ux & ~mask);
 93        PUT_BITS_DP64(ur, r);
 94        if (xneg || (ur == ux))
 95          return r;
 96        else
 97          /* We threw some bits away and x was positive */
 98          return r + 1.0;
 99      }
100  
101  }
102  
103