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-exception
6//
7//===----------------------------------------------------------------------===//
8
9#ifndef LLVM_LIBC_SRC___SUPPORT_FPUTIL_NEARESTINTEGEROPERATIONS_H
10#define LLVM_LIBC_SRC___SUPPORT_FPUTIL_NEARESTINTEGEROPERATIONS_H
11
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
20namespace LIBC_NAMESPACE {
21namespace fputil {
22
23template <typename T, cpp::enable_if_t<cpp::is_floating_point_v<T>, int> = 0>
24LIBC_INLINE T trunc(T x) {
25 FPBits<T> bits(x);
26
27 // If x is infinity or NaN, return it.
28 // If it is zero also we should return it as is, but the logic
29 // later in this function takes care of it. But not doing a zero
30 // check, we improve the run time of non-zero values.
31 if (bits.is_inf_or_nan())
32 return x;
33
34 int exponent = bits.get_exponent();
35
36 // If the exponent is greater than the most negative mantissa
37 // exponent, then x is already an integer.
38 if (exponent >= static_cast<int>(FPBits<T>::FRACTION_LEN))
39 return x;
40
41 // If the exponent is such that abs(x) is less than 1, then return 0.
42 if (exponent <= -1)
43 return FPBits<T>::zero(bits.sign()).get_val();
44
45 int trim_size = FPBits<T>::FRACTION_LEN - exponent;
46 bits.set_mantissa((bits.get_mantissa() >> trim_size) << trim_size);
47 return bits.get_val();
48}
49
50template <typename T, cpp::enable_if_t<cpp::is_floating_point_v<T>, int> = 0>
51LIBC_INLINE T ceil(T x) {
52 FPBits<T> bits(x);
53
54 // If x is infinity NaN or zero, return it.
55 if (bits.is_inf_or_nan() || bits.is_zero())
56 return x;
57
58 bool is_neg = bits.is_neg();
59 int exponent = bits.get_exponent();
60
61 // If the exponent is greater than the most negative mantissa
62 // exponent, then x is already an integer.
63 if (exponent >= static_cast<int>(FPBits<T>::FRACTION_LEN))
64 return x;
65
66 if (exponent <= -1) {
67 if (is_neg)
68 return T(-0.0);
69 else
70 return T(1.0);
71 }
72
73 uint32_t trim_size = FPBits<T>::FRACTION_LEN - exponent;
74 bits.set_mantissa((bits.get_mantissa() >> trim_size) << trim_size);
75 T trunc_value = bits.get_val();
76
77 // If x is already an integer, return it.
78 if (trunc_value == x)
79 return x;
80
81 // If x is negative, the ceil operation is equivalent to the trunc operation.
82 if (is_neg)
83 return trunc_value;
84
85 return trunc_value + T(1.0);
86}
87
88template <typename T, cpp::enable_if_t<cpp::is_floating_point_v<T>, int> = 0>
89LIBC_INLINE T floor(T x) {
90 FPBits<T> bits(x);
91 if (bits.is_neg()) {
92 return -ceil(-x);
93 } else {
94 return trunc(x);
95 }
96}
97
98template <typename T, cpp::enable_if_t<cpp::is_floating_point_v<T>, int> = 0>
99LIBC_INLINE T round(T x) {
100 using StorageType = typename FPBits<T>::StorageType;
101 FPBits<T> bits(x);
102
103 // If x is infinity NaN or zero, return it.
104 if (bits.is_inf_or_nan() || bits.is_zero())
105 return x;
106
107 int exponent = bits.get_exponent();
108
109 // If the exponent is greater than the most negative mantissa
110 // exponent, then x is already an integer.
111 if (exponent >= static_cast<int>(FPBits<T>::FRACTION_LEN))
112 return x;
113
114 if (exponent == -1) {
115 // Absolute value of x is greater than equal to 0.5 but less than 1.
116 return FPBits<T>::one(bits.sign()).get_val();
117 }
118
119 if (exponent <= -2) {
120 // Absolute value of x is less than 0.5.
121 return FPBits<T>::zero(bits.sign()).get_val();
122 }
123
124 uint32_t trim_size = FPBits<T>::FRACTION_LEN - exponent;
125 bool half_bit_set =
126 bool(bits.get_mantissa() & (StorageType(1) << (trim_size - 1)));
127 bits.set_mantissa((bits.get_mantissa() >> trim_size) << trim_size);
128 T trunc_value = bits.get_val();
129
130 // If x is already an integer, return it.
131 if (trunc_value == x)
132 return x;
133
134 if (!half_bit_set) {
135 // Franctional part is less than 0.5 so round value is the
136 // same as the trunc value.
137 return trunc_value;
138 } else {
139 return bits.is_neg() ? trunc_value - T(1.0) : trunc_value + T(1.0);
140 }
141}
142
143template <typename T>
144LIBC_INLINE constexpr cpp::enable_if_t<cpp::is_floating_point_v<T>, T>
145round_using_specific_rounding_mode(T x, int rnd) {
146 using StorageType = typename FPBits<T>::StorageType;
147 FPBits<T> bits(x);
148
149 // If x is infinity NaN or zero, return it.
150 if (bits.is_inf_or_nan() || bits.is_zero())
151 return x;
152
153 bool is_neg = bits.is_neg();
154 int exponent = bits.get_exponent();
155
156 // If the exponent is greater than the most negative mantissa
157 // exponent, then x is already an integer.
158 if (exponent >= static_cast<int>(FPBits<T>::FRACTION_LEN))
159 return x;
160
161 if (exponent <= -1) {
162 switch (rnd) {
163 case FP_INT_DOWNWARD:
164 return is_neg ? T(-1.0) : T(0.0);
165 case FP_INT_UPWARD:
166 return is_neg ? T(-0.0) : T(1.0);
167 case FP_INT_TOWARDZERO:
168 return is_neg ? T(-0.0) : T(0.0);
169 case FP_INT_TONEARESTFROMZERO:
170 if (exponent < -1)
171 return is_neg ? T(-0.0) : T(0.0); // abs(x) < 0.5
172 return is_neg ? T(-1.0) : T(1.0); // abs(x) >= 0.5
173 case FP_INT_TONEAREST:
174 default:
175 if (exponent <= -2 || bits.get_mantissa() == 0)
176 return is_neg ? T(-0.0) : T(0.0); // abs(x) <= 0.5
177 else
178 return is_neg ? T(-1.0) : T(1.0); // abs(x) > 0.5
179 }
180 }
181
182 uint32_t trim_size = FPBits<T>::FRACTION_LEN - exponent;
183 FPBits<T> new_bits = bits;
184 new_bits.set_mantissa((bits.get_mantissa() >> trim_size) << trim_size);
185 T trunc_value = new_bits.get_val();
186
187 // If x is already an integer, return it.
188 if (trunc_value == x)
189 return x;
190
191 StorageType trim_value =
192 bits.get_mantissa() & ((StorageType(1) << trim_size) - 1);
193 StorageType half_value = (StorageType(1) << (trim_size - 1));
194 // If exponent is 0, trimSize will be equal to the mantissa width, and
195 // truncIsOdd` will not be correct. So, we handle it as a special case
196 // below.
197 StorageType trunc_is_odd =
198 new_bits.get_mantissa() & (StorageType(1) << trim_size);
199
200 switch (rnd) {
201 case FP_INT_DOWNWARD:
202 return is_neg ? trunc_value - T(1.0) : trunc_value;
203 case FP_INT_UPWARD:
204 return is_neg ? trunc_value : trunc_value + T(1.0);
205 case FP_INT_TOWARDZERO:
206 return trunc_value;
207 case FP_INT_TONEARESTFROMZERO:
208 if (trim_value >= half_value)
209 return is_neg ? trunc_value - T(1.0) : trunc_value + T(1.0);
210 return trunc_value;
211 case FP_INT_TONEAREST:
212 default:
213 if (trim_value > half_value) {
214 return is_neg ? trunc_value - T(1.0) : trunc_value + T(1.0);
215 } else if (trim_value == half_value) {
216 if (exponent == 0)
217 return is_neg ? T(-2.0) : T(2.0);
218 if (trunc_is_odd)
219 return is_neg ? trunc_value - T(1.0) : trunc_value + T(1.0);
220 else
221 return trunc_value;
222 } else {
223 return trunc_value;
224 }
225 }
226}
227
228template <typename T>
229LIBC_INLINE cpp::enable_if_t<cpp::is_floating_point_v<T>, T>
230round_using_current_rounding_mode(T x) {
231 int rounding_mode = quick_get_round();
232
233 switch (rounding_mode) {
234 case FE_DOWNWARD:
235 return round_using_specific_rounding_mode(x, FP_INT_DOWNWARD);
236 case FE_UPWARD:
237 return round_using_specific_rounding_mode(x, FP_INT_UPWARD);
238 case FE_TOWARDZERO:
239 return round_using_specific_rounding_mode(x, FP_INT_TOWARDZERO);
240 case FE_TONEAREST:
241 return round_using_specific_rounding_mode(x, FP_INT_TONEAREST);
242 default:
243 __builtin_unreachable();
244 }
245}
246
247template <bool IsSigned, typename T>
248LIBC_INLINE constexpr cpp::enable_if_t<cpp::is_floating_point_v<T>, T>
249fromfp(T x, int rnd, unsigned int width) {
250 using StorageType = typename FPBits<T>::StorageType;
251
252 constexpr StorageType EXPLICIT_BIT =
253 FPBits<T>::SIG_MASK - FPBits<T>::FRACTION_MASK;
254
255 if (width == 0U) {
256 raise_except_if_required(FE_INVALID);
257 return FPBits<T>::quiet_nan().get_val();
258 }
259
260 FPBits<T> bits(x);
261
262 if (bits.is_inf_or_nan()) {
263 raise_except_if_required(FE_INVALID);
264 return FPBits<T>::quiet_nan().get_val();
265 }
266
267 T rounded_value = round_using_specific_rounding_mode(x, rnd);
268
269 if constexpr (IsSigned) {
270 // T can't hold a finite number >= 2.0 * 2^EXP_BIAS.
271 if (width - 1 > FPBits<T>::EXP_BIAS)
272 return rounded_value;
273
274 StorageType range_exp = width - 1U + FPBits<T>::EXP_BIAS;
275 // rounded_value < -2^(width - 1)
276 T range_min =
277 FPBits<T>::create_value(Sign::NEG, range_exp, EXPLICIT_BIT).get_val();
278 if (rounded_value < range_min) {
279 raise_except_if_required(FE_INVALID);
280 return FPBits<T>::quiet_nan().get_val();
281 }
282 // rounded_value > 2^(width - 1) - 1
283 T range_max =
284 FPBits<T>::create_value(Sign::POS, range_exp, EXPLICIT_BIT).get_val() -
285 T(1.0);
286 if (rounded_value > range_max) {
287 raise_except_if_required(FE_INVALID);
288 return FPBits<T>::quiet_nan().get_val();
289 }
290
291 return rounded_value;
292 }
293
294 if (rounded_value < T(0.0)) {
295 raise_except_if_required(FE_INVALID);
296 return FPBits<T>::quiet_nan().get_val();
297 }
298
299 // T can't hold a finite number >= 2.0 * 2^EXP_BIAS.
300 if (width > FPBits<T>::EXP_BIAS)
301 return rounded_value;
302
303 StorageType range_exp = width + FPBits<T>::EXP_BIAS;
304 // rounded_value > 2^width - 1
305 T range_max =
306 FPBits<T>::create_value(Sign::POS, range_exp, EXPLICIT_BIT).get_val() -
307 T(1.0);
308 if (rounded_value > range_max) {
309 raise_except_if_required(FE_INVALID);
310 return FPBits<T>::quiet_nan().get_val();
311 }
312
313 return rounded_value;
314}
315
316template <bool IsSigned, typename T>
317LIBC_INLINE constexpr cpp::enable_if_t<cpp::is_floating_point_v<T>, T>
318fromfpx(T x, int rnd, unsigned int width) {
319 T rounded_value = fromfp<IsSigned>(x, rnd, width);
320 FPBits<T> bits(rounded_value);
321
322 if (!bits.is_nan() && rounded_value != x)
323 raise_except_if_required(FE_INEXACT);
324
325 return rounded_value;
326}
327
328namespace internal {
329
330template <typename F, typename I,
331 cpp::enable_if_t<cpp::is_floating_point_v<F> && cpp::is_integral_v<I>,
332 int> = 0>
333LIBC_INLINE I rounded_float_to_signed_integer(F x) {
334 constexpr I INTEGER_MIN = (I(1) << (sizeof(I) * 8 - 1));
335 constexpr I INTEGER_MAX = -(INTEGER_MIN + 1);
336 FPBits<F> bits(x);
337 auto set_domain_error_and_raise_invalid = []() {
338 set_errno_if_required(EDOM);
339 raise_except_if_required(FE_INVALID);
340 };
341
342 if (bits.is_inf_or_nan()) {
343 set_domain_error_and_raise_invalid();
344 return bits.is_neg() ? INTEGER_MIN : INTEGER_MAX;
345 }
346
347 int exponent = bits.get_exponent();
348 constexpr int EXPONENT_LIMIT = sizeof(I) * 8 - 1;
349 if (exponent > EXPONENT_LIMIT) {
350 set_domain_error_and_raise_invalid();
351 return bits.is_neg() ? INTEGER_MIN : INTEGER_MAX;
352 } else if (exponent == EXPONENT_LIMIT) {
353 if (bits.is_pos() || bits.get_mantissa() != 0) {
354 set_domain_error_and_raise_invalid();
355 return bits.is_neg() ? INTEGER_MIN : INTEGER_MAX;
356 }
357 // If the control reaches here, then it means that the rounded
358 // value is the most negative number for the signed integer type I.
359 }
360
361 // For all other cases, if `x` can fit in the integer type `I`,
362 // we just return `x`. static_cast will convert the floating
363 // point value to the exact integer value.
364 return static_cast<I>(x);
365}
366
367} // namespace internal
368
369template <typename F, typename I,
370 cpp::enable_if_t<cpp::is_floating_point_v<F> && cpp::is_integral_v<I>,
371 int> = 0>
372LIBC_INLINE I round_to_signed_integer(F x) {
373 return internal::rounded_float_to_signed_integer<F, I>(round(x));
374}
375
376template <typename F, typename I,
377 cpp::enable_if_t<cpp::is_floating_point_v<F> && cpp::is_integral_v<I>,
378 int> = 0>
379LIBC_INLINE I round_to_signed_integer_using_current_rounding_mode(F x) {
380 return internal::rounded_float_to_signed_integer<F, I>(
381 round_using_current_rounding_mode(x));
382}
383
384} // namespace fputil
385} // namespace LIBC_NAMESPACE
386
387#endif // LLVM_LIBC_SRC___SUPPORT_FPUTIL_NEARESTINTEGEROPERATIONS_H
388

source code of libc/src/__support/FPUtil/NearestIntegerOperations.h