libc/src/__support/FPUtil/DivisionAndRemainderOperations.h

   1 //===-- Floating point divsion and remainder 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-exception
   6 //
   7 //===----------------------------------------------------------------------===//
   8
   9 #ifndef LLVM_LIBC_SRC_SUPPORT_FPUTIL_DIVISION_AND_REMAINDER_OPERATIONS_H
  10 #define LLVM_LIBC_SRC_SUPPORT_FPUTIL_DIVISION_AND_REMAINDER_OPERATIONS_H
  11
  12 #include "FPBits.h"
  13 #include "ManipulationFunctions.h"
  14 #include "NormalFloat.h"
  15
  16 #include "src/__support/CPP/type_traits.h"
  17 #include "src/__support/common.h"
  18
  19 namespace __llvm_libc {
  20 namespace fputil {
  21
  22 static constexpr int QUOTIENT_LSB_BITS = 3;
  23
  24 // The implementation is a bit-by-bit algorithm which uses integer division
  25 // to evaluate the quotient and remainder.
  26 template <typename T, cpp::enable_if_t<cpp::is_floating_point_v<T>, int> = 0>
  27 LIBC_INLINE T remquo(T x, T y, int &q) {
  28   FPBits<T> xbits(x), ybits(y);
  29   if (xbits.is_nan())
  30     return x;
  31   if (ybits.is_nan())
  32     return y;
  33   if (xbits.is_inf() || ybits.is_zero())
  34     return FPBits<T>::build_quiet_nan(1);
  35
  36   if (xbits.is_zero()) {
  37     q = 0;
  38     return __llvm_libc::fputil::copysign(T(0.0), x);
  39   }
  40
  41   if (ybits.is_inf()) {
  42     q = 0;
  43     return x;
  44   }
  45
  46   bool result_sign = (xbits.get_sign() == ybits.get_sign() ? false : true);
  47
  48   // Once we know the sign of the result, we can just operate on the absolute
  49   // values. The correct sign can be applied to the result after the result
  50   // is evaluated.
  51   xbits.set_sign(0);
  52   ybits.set_sign(0);
  53
  54   NormalFloat<T> normalx(xbits), normaly(ybits);
  55   int exp = normalx.exponent - normaly.exponent;
  56   typename NormalFloat<T>::UIntType mx = normalx.mantissa,
  57                                     my = normaly.mantissa;
  58
  59   q = 0;
  60   while (exp >= 0) {
  61     unsigned shift_count = 0;
  62     typename NormalFloat<T>::UIntType n = mx;
  63     for (shift_count = 0; n < my; n <<= 1, ++shift_count)
  64       ;
  65
  66     if (static_cast<int>(shift_count) > exp)
  67       break;
  68
  69     exp -= shift_count;
  70     if (0 <= exp && exp < QUOTIENT_LSB_BITS)
  71       q |= (1 << exp);
  72
  73     mx = n - my;
  74     if (mx == 0) {
  75       q = result_sign ? -q : q;
  76       return __llvm_libc::fputil::copysign(T(0.0), x);
  77     }
  78   }
  79
  80   NormalFloat<T> remainder(exp + normaly.exponent, mx, 0);
  81
  82   // Since NormalFloat to native type conversion is a truncation operation
  83   // currently, the remainder value in the native type is correct as is.
  84   // However, if NormalFloat to native type conversion is updated in future,
  85   // then the conversion to native remainder value should be updated
  86   // appropriately and some directed tests added.
  87   T native_remainder(remainder);
  88   T absy = T(ybits);
  89   int cmp = remainder.mul2(1).cmp(normaly);
  90   if (cmp > 0) {
  91     q = q + 1;
  92     if (x >= T(0.0))
  93       native_remainder = native_remainder - absy;
  94     else
  95       native_remainder = absy - native_remainder;
  96   } else if (cmp == 0) {
  97     if (q & 1) {
  98       q += 1;
  99       if (x >= T(0.0))
 100         native_remainder = -native_remainder;
 101     } else {
 102       if (x < T(0.0))
 103         native_remainder = -native_remainder;
 104     }
 105   } else {
 106     if (x < T(0.0))
 107       native_remainder = -native_remainder;
 108   }
 109
 110   q = result_sign ? -q : q;
 111   if (native_remainder == T(0.0))
 112     return __llvm_libc::fputil::copysign(T(0.0), x);
 113   return native_remainder;
 114 }
 115
 116 } // namespace fputil
 117 } // namespace __llvm_libc
 118
 119 #endif // LLVM_LIBC_SRC_SUPPORT_FPUTIL_DIVISION_AND_REMAINDER_OPERATIONS_H