brintos

brintos / llvm-project-archived public Read only

0
0
Text · 12.6 KiB · 9316661 Raw
412 lines · c
1//===-- Nearest integer floating-point operations ---------------*- C++ -*-===//2//3// Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions.4// See https://llvm.org/LICENSE.txt for license information.5// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception6//7//===----------------------------------------------------------------------===//8 9#ifndef LLVM_LIBC_SRC___SUPPORT_FPUTIL_NEARESTINTEGEROPERATIONS_H10#define LLVM_LIBC_SRC___SUPPORT_FPUTIL_NEARESTINTEGEROPERATIONS_H11 12#include "FEnvImpl.h"13#include "FPBits.h"14#include "rounding_mode.h"15 16#include "hdr/math_macros.h"17#include "src/__support/CPP/type_traits.h"18#include "src/__support/common.h"19#include "src/__support/macros/config.h"20 21namespace LIBC_NAMESPACE_DECL {22namespace fputil {23 24template <typename T, cpp::enable_if_t<cpp::is_floating_point_v<T>, int> = 0>25LIBC_INLINE T trunc(T x) {26  using StorageType = typename FPBits<T>::StorageType;27  FPBits<T> bits(x);28 29  // If x is infinity or NaN, return it.30  // If it is zero also we should return it as is, but the logic31  // later in this function takes care of it. But not doing a zero32  // check, we improve the run time of non-zero values.33  if (bits.is_inf_or_nan())34    return x;35 36  int exponent = bits.get_exponent();37 38  // If the exponent is greater than the most negative mantissa39  // exponent, then x is already an integer.40  if (exponent >= static_cast<int>(FPBits<T>::FRACTION_LEN))41    return x;42 43  // If the exponent is such that abs(x) is less than 1, then return 0.44  if (exponent <= -1)45    return FPBits<T>::zero(bits.sign()).get_val();46 47  int trim_size = FPBits<T>::FRACTION_LEN - exponent;48  StorageType trunc_mantissa =49      static_cast<StorageType>((bits.get_mantissa() >> trim_size) << trim_size);50  bits.set_mantissa(trunc_mantissa);51  return bits.get_val();52}53 54template <typename T, cpp::enable_if_t<cpp::is_floating_point_v<T>, int> = 0>55LIBC_INLINE T ceil(T x) {56  using StorageType = typename FPBits<T>::StorageType;57  FPBits<T> bits(x);58 59  // If x is infinity NaN or zero, return it.60  if (bits.is_inf_or_nan() || bits.is_zero())61    return x;62 63  bool is_neg = bits.is_neg();64  int exponent = bits.get_exponent();65 66  // If the exponent is greater than the most negative mantissa67  // exponent, then x is already an integer.68  if (exponent >= static_cast<int>(FPBits<T>::FRACTION_LEN))69    return x;70 71  if (exponent <= -1) {72    if (is_neg)73      return T(-0.0);74    else75      return T(1.0);76  }77 78  uint32_t trim_size = FPBits<T>::FRACTION_LEN - exponent;79  StorageType x_u = bits.uintval();80  StorageType trunc_u =81      static_cast<StorageType>((x_u >> trim_size) << trim_size);82 83  // If x is already an integer, return it.84  if (trunc_u == x_u)85    return x;86 87  bits.set_uintval(trunc_u);88  T trunc_value = bits.get_val();89 90  // If x is negative, the ceil operation is equivalent to the trunc operation.91  if (is_neg)92    return trunc_value;93 94  return trunc_value + T(1.0);95}96 97template <typename T, cpp::enable_if_t<cpp::is_floating_point_v<T>, int> = 0>98LIBC_INLINE T floor(T x) {99  FPBits<T> bits(x);100  if (bits.is_neg()) {101    return -ceil(-x);102  } else {103    return trunc(x);104  }105}106 107template <typename T, cpp::enable_if_t<cpp::is_floating_point_v<T>, int> = 0>108LIBC_INLINE T round(T x) {109  using StorageType = typename FPBits<T>::StorageType;110  FPBits<T> bits(x);111 112  // If x is infinity NaN or zero, return it.113  if (bits.is_inf_or_nan() || bits.is_zero())114    return x;115 116  int exponent = bits.get_exponent();117 118  // If the exponent is greater than the most negative mantissa119  // exponent, then x is already an integer.120  if (exponent >= static_cast<int>(FPBits<T>::FRACTION_LEN))121    return x;122 123  if (exponent == -1) {124    // Absolute value of x is greater than equal to 0.5 but less than 1.125    return FPBits<T>::one(bits.sign()).get_val();126  }127 128  if (exponent <= -2) {129    // Absolute value of x is less than 0.5.130    return FPBits<T>::zero(bits.sign()).get_val();131  }132 133  uint32_t trim_size = FPBits<T>::FRACTION_LEN - exponent;134  bool half_bit_set =135      bool(bits.get_mantissa() & (StorageType(1) << (trim_size - 1)));136  StorageType x_u = bits.uintval();137  StorageType trunc_u =138      static_cast<StorageType>((x_u >> trim_size) << trim_size);139 140  // If x is already an integer, return it.141  if (trunc_u == x_u)142    return x;143 144  bits.set_uintval(trunc_u);145  T trunc_value = bits.get_val();146 147  if (!half_bit_set) {148    // Franctional part is less than 0.5 so round value is the149    // same as the trunc value.150    return trunc_value;151  } else {152    return bits.is_neg() ? trunc_value - T(1.0) : trunc_value + T(1.0);153  }154}155 156template <typename T>157LIBC_INLINE constexpr cpp::enable_if_t<cpp::is_floating_point_v<T>, T>158round_using_specific_rounding_mode(T x, int rnd) {159  using StorageType = typename FPBits<T>::StorageType;160  FPBits<T> bits(x);161 162  // If x is infinity NaN or zero, return it.163  if (bits.is_inf_or_nan() || bits.is_zero())164    return x;165 166  bool is_neg = bits.is_neg();167  int exponent = bits.get_exponent();168 169  // If the exponent is greater than the most negative mantissa170  // exponent, then x is already an integer.171  if (exponent >= static_cast<int>(FPBits<T>::FRACTION_LEN))172    return x;173 174  if (exponent <= -1) {175    switch (rnd) {176    case FP_INT_DOWNWARD:177      return is_neg ? T(-1.0) : T(0.0);178    case FP_INT_UPWARD:179      return is_neg ? T(-0.0) : T(1.0);180    case FP_INT_TOWARDZERO:181      return is_neg ? T(-0.0) : T(0.0);182    case FP_INT_TONEARESTFROMZERO:183      if (exponent < -1)184        return is_neg ? T(-0.0) : T(0.0); // abs(x) < 0.5185      return is_neg ? T(-1.0) : T(1.0);   // abs(x) >= 0.5186    case FP_INT_TONEAREST:187    default:188      if (exponent <= -2 || bits.get_mantissa() == 0)189        return is_neg ? T(-0.0) : T(0.0); // abs(x) <= 0.5190      else191        return is_neg ? T(-1.0) : T(1.0); // abs(x) > 0.5192    }193  }194 195  uint32_t trim_size = FPBits<T>::FRACTION_LEN - exponent;196  StorageType x_u = bits.uintval();197  StorageType trunc_u =198      static_cast<StorageType>((x_u >> trim_size) << trim_size);199 200  // If x is already an integer, return it.201  if (trunc_u == x_u)202    return x;203 204  FPBits<T> new_bits(trunc_u);205  T trunc_value = new_bits.get_val();206 207  StorageType trim_value =208      bits.get_mantissa() &209      static_cast<StorageType>(((StorageType(1) << trim_size) - 1));210  StorageType half_value =211      static_cast<StorageType>((StorageType(1) << (trim_size - 1)));212  // If exponent is 0, trimSize will be equal to the mantissa width, and213  // truncIsOdd` will not be correct. So, we handle it as a special case214  // below.215  StorageType trunc_is_odd =216      new_bits.get_mantissa() & (StorageType(1) << trim_size);217 218  switch (rnd) {219  case FP_INT_DOWNWARD:220    return is_neg ? trunc_value - T(1.0) : trunc_value;221  case FP_INT_UPWARD:222    return is_neg ? trunc_value : trunc_value + T(1.0);223  case FP_INT_TOWARDZERO:224    return trunc_value;225  case FP_INT_TONEARESTFROMZERO:226    if (trim_value >= half_value)227      return is_neg ? trunc_value - T(1.0) : trunc_value + T(1.0);228    return trunc_value;229  case FP_INT_TONEAREST:230  default:231    if (trim_value > half_value) {232      return is_neg ? trunc_value - T(1.0) : trunc_value + T(1.0);233    } else if (trim_value == half_value) {234      if (exponent == 0)235        return is_neg ? T(-2.0) : T(2.0);236      if (trunc_is_odd)237        return is_neg ? trunc_value - T(1.0) : trunc_value + T(1.0);238      else239        return trunc_value;240    } else {241      return trunc_value;242    }243  }244}245 246template <typename T>247LIBC_INLINE cpp::enable_if_t<cpp::is_floating_point_v<T>, T>248round_using_current_rounding_mode(T x) {249  int rounding_mode = quick_get_round();250 251  switch (rounding_mode) {252  case FE_DOWNWARD:253    return round_using_specific_rounding_mode(x, FP_INT_DOWNWARD);254  case FE_UPWARD:255    return round_using_specific_rounding_mode(x, FP_INT_UPWARD);256  case FE_TOWARDZERO:257    return round_using_specific_rounding_mode(x, FP_INT_TOWARDZERO);258  case FE_TONEAREST:259    return round_using_specific_rounding_mode(x, FP_INT_TONEAREST);260  default:261    __builtin_unreachable();262  }263}264 265template <bool IsSigned, typename T>266LIBC_INLINE constexpr cpp::enable_if_t<cpp::is_floating_point_v<T>, T>267fromfp(T x, int rnd, unsigned int width) {268  using StorageType = typename FPBits<T>::StorageType;269 270  constexpr StorageType EXPLICIT_BIT =271      FPBits<T>::SIG_MASK - FPBits<T>::FRACTION_MASK;272 273  if (width == 0U) {274    raise_except_if_required(FE_INVALID);275    return FPBits<T>::quiet_nan().get_val();276  }277 278  FPBits<T> bits(x);279 280  if (bits.is_inf_or_nan()) {281    raise_except_if_required(FE_INVALID);282    return FPBits<T>::quiet_nan().get_val();283  }284 285  T rounded_value = round_using_specific_rounding_mode(x, rnd);286 287  if constexpr (IsSigned) {288    // T can't hold a finite number >= 2.0 * 2^EXP_BIAS.289    if (width - 1 > FPBits<T>::EXP_BIAS)290      return rounded_value;291 292    StorageType range_exp =293        static_cast<StorageType>(width - 1 + FPBits<T>::EXP_BIAS);294    // rounded_value < -2^(width - 1)295    T range_min =296        FPBits<T>::create_value(Sign::NEG, range_exp, EXPLICIT_BIT).get_val();297    if (rounded_value < range_min) {298      raise_except_if_required(FE_INVALID);299      return FPBits<T>::quiet_nan().get_val();300    }301    // rounded_value > 2^(width - 1) - 1302    T range_max =303        FPBits<T>::create_value(Sign::POS, range_exp, EXPLICIT_BIT).get_val() -304        T(1.0);305    if (rounded_value > range_max) {306      raise_except_if_required(FE_INVALID);307      return FPBits<T>::quiet_nan().get_val();308    }309 310    return rounded_value;311  }312 313  if (rounded_value < T(0.0)) {314    raise_except_if_required(FE_INVALID);315    return FPBits<T>::quiet_nan().get_val();316  }317 318  // T can't hold a finite number >= 2.0 * 2^EXP_BIAS.319  if (width > FPBits<T>::EXP_BIAS)320    return rounded_value;321 322  StorageType range_exp = static_cast<StorageType>(width + FPBits<T>::EXP_BIAS);323  // rounded_value > 2^width - 1324  T range_max =325      FPBits<T>::create_value(Sign::POS, range_exp, EXPLICIT_BIT).get_val() -326      T(1.0);327  if (rounded_value > range_max) {328    raise_except_if_required(FE_INVALID);329    return FPBits<T>::quiet_nan().get_val();330  }331 332  return rounded_value;333}334 335template <bool IsSigned, typename T>336LIBC_INLINE constexpr cpp::enable_if_t<cpp::is_floating_point_v<T>, T>337fromfpx(T x, int rnd, unsigned int width) {338  T rounded_value = fromfp<IsSigned>(x, rnd, width);339  FPBits<T> bits(rounded_value);340 341  if (!bits.is_nan() && rounded_value != x)342    raise_except_if_required(FE_INEXACT);343 344  return rounded_value;345}346 347namespace internal {348 349template <typename FloatType, typename IntType,350          cpp::enable_if_t<cpp::is_floating_point_v<FloatType> &&351                               cpp::is_integral_v<IntType>,352                           int> = 0>353LIBC_INLINE IntType rounded_float_to_signed_integer(FloatType x) {354  constexpr IntType INTEGER_MIN = (IntType(1) << (sizeof(IntType) * 8 - 1));355  constexpr IntType INTEGER_MAX = -(INTEGER_MIN + 1);356  FPBits<FloatType> bits(x);357  auto set_domain_error_and_raise_invalid = []() {358    set_errno_if_required(EDOM);359    raise_except_if_required(FE_INVALID);360  };361 362  if (bits.is_inf_or_nan()) {363    set_domain_error_and_raise_invalid();364    return bits.is_neg() ? INTEGER_MIN : INTEGER_MAX;365  }366 367  int exponent = bits.get_exponent();368  constexpr int EXPONENT_LIMIT = sizeof(IntType) * 8 - 1;369  if (exponent > EXPONENT_LIMIT) {370    set_domain_error_and_raise_invalid();371    return bits.is_neg() ? INTEGER_MIN : INTEGER_MAX;372  } else if (exponent == EXPONENT_LIMIT) {373    if (bits.is_pos() || bits.get_mantissa() != 0) {374      set_domain_error_and_raise_invalid();375      return bits.is_neg() ? INTEGER_MIN : INTEGER_MAX;376    }377    // If the control reaches here, then it means that the rounded378    // value is the most negative number for the signed integer type IntType.379  }380 381  // For all other cases, if `x` can fit in the integer type `IntType`,382  // we just return `x`. static_cast will convert the floating383  // point value to the exact integer value.384  return static_cast<IntType>(x);385}386 387} // namespace internal388 389template <typename FloatType, typename IntType,390          cpp::enable_if_t<cpp::is_floating_point_v<FloatType> &&391                               cpp::is_integral_v<IntType>,392                           int> = 0>393LIBC_INLINE IntType round_to_signed_integer(FloatType x) {394  return internal::rounded_float_to_signed_integer<FloatType, IntType>(395      round(x));396}397 398template <typename FloatType, typename IntType,399          cpp::enable_if_t<cpp::is_floating_point_v<FloatType> &&400                               cpp::is_integral_v<IntType>,401                           int> = 0>402LIBC_INLINE IntType403round_to_signed_integer_using_current_rounding_mode(FloatType x) {404  return internal::rounded_float_to_signed_integer<FloatType, IntType>(405      round_using_current_rounding_mode(x));406}407 408} // namespace fputil409} // namespace LIBC_NAMESPACE_DECL410 411#endif // LLVM_LIBC_SRC___SUPPORT_FPUTIL_NEARESTINTEGEROPERATIONS_H412