Skip to content

Commit 4d21e75

Browse files
lntuelntue
authored
[libc][math][c23] Add fmodl and fmodf128 math functions. (#84600)
- Allow `FMod` template to have different computational types and make it work for 80-bit long double. - Switch to use `uint64_t` as the intermediate computational types for `float`, significantly reduce the latency of `fmodf` when the exponent difference is large.
1 parent 4628e33 commit 4d21e75

File tree

23 files changed

+303
-105
lines changed

23 files changed

+303
-105
lines changed

‎libc/config/linux/aarch64/entrypoints.txt

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -335,6 +335,7 @@ set(TARGET_LIBM_ENTRYPOINTS
335335
libc.src.math.fminl
336336
libc.src.math.fmod
337337
libc.src.math.fmodf
338+
libc.src.math.fmodl
338339
libc.src.math.frexp
339340
libc.src.math.frexpf
340341
libc.src.math.frexpl
@@ -426,6 +427,7 @@ if(LIBC_TYPES_HAS_FLOAT128)
426427
libc.src.math.floorf128
427428
libc.src.math.fmaxf128
428429
libc.src.math.fminf128
430+
libc.src.math.fmodf128
429431
libc.src.math.frexpf128
430432
libc.src.math.ilogbf128
431433
libc.src.math.ldexpf128

‎libc/config/linux/riscv/entrypoints.txt

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -343,6 +343,7 @@ set(TARGET_LIBM_ENTRYPOINTS
343343
libc.src.math.fmaxl
344344
libc.src.math.fmod
345345
libc.src.math.fmodf
346+
libc.src.math.fmodl
346347
libc.src.math.frexp
347348
libc.src.math.frexpf
348349
libc.src.math.frexpl
@@ -434,6 +435,7 @@ if(LIBC_TYPES_HAS_FLOAT128)
434435
libc.src.math.floorf128
435436
libc.src.math.fmaxf128
436437
libc.src.math.fminf128
438+
libc.src.math.fmodf128
437439
libc.src.math.frexpf128
438440
libc.src.math.ilogbf128
439441
libc.src.math.ldexpf128

‎libc/config/linux/x86_64/entrypoints.txt

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -376,6 +376,7 @@ set(TARGET_LIBM_ENTRYPOINTS
376376
libc.src.math.fmaxl
377377
libc.src.math.fmod
378378
libc.src.math.fmodf
379+
libc.src.math.fmodl
379380
libc.src.math.frexp
380381
libc.src.math.frexpf
381382
libc.src.math.frexpl
@@ -469,6 +470,7 @@ if(LIBC_TYPES_HAS_FLOAT128)
469470
libc.src.math.floorf128
470471
libc.src.math.fmaxf128
471472
libc.src.math.fminf128
473+
libc.src.math.fmodf128
472474
libc.src.math.frexpf128
473475
libc.src.math.ilogbf128
474476
libc.src.math.ldexpf128

‎libc/config/windows/entrypoints.txt

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -155,6 +155,7 @@ set(TARGET_LIBM_ENTRYPOINTS
155155
libc.src.math.fmaxl
156156
libc.src.math.fmod
157157
libc.src.math.fmodf
158+
libc.src.math.fmodl
158159
libc.src.math.frexp
159160
libc.src.math.frexpf
160161
libc.src.math.frexpl

‎libc/docs/math/index.rst

Lines changed: 3 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -169,7 +169,9 @@ Basic Operations
169169
+--------------+---------+---------+---------+---------+---------+---------+---------+---------+---------+---------+---------+---------+
170170
| fmodf | |check| | |check| | |check| | |check| | |check| | | | |check| | |check| | |check| | | |
171171
+--------------+---------+---------+---------+---------+---------+---------+---------+---------+---------+---------+---------+---------+
172-
| fmodl | | | | | | | | | | | | |
172+
| fmodl | |check| | |check| | | |check| | |check| | | | |check| | | | | |
173+
+--------------+---------+---------+---------+---------+---------+---------+---------+---------+---------+---------+---------+---------+
174+
| fmodf128 | |check| | |check| | | |check| | | | | | | | | |
173175
+--------------+---------+---------+---------+---------+---------+---------+---------+---------+---------+---------+---------+---------+
174176
| frexp | |check| | |check| | |check| | |check| | |check| | | | |check| | |check| | |check| | | |
175177
+--------------+---------+---------+---------+---------+---------+---------+---------+---------+---------+---------+---------+---------+

‎libc/spec/stdc.td

Lines changed: 2 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -405,8 +405,9 @@ def StdC : StandardSpec<"stdc"> {
405405
FunctionSpec<"fmaf", RetValSpec<FloatType>, [ArgSpec<FloatType>, ArgSpec<FloatType>, ArgSpec<FloatType>]>,
406406

407407
FunctionSpec<"fmod", RetValSpec<DoubleType>, [ArgSpec<DoubleType>, ArgSpec<DoubleType>]>,
408-
409408
FunctionSpec<"fmodf", RetValSpec<FloatType>, [ArgSpec<FloatType>, ArgSpec<FloatType>]>,
409+
FunctionSpec<"fmodl", RetValSpec<LongDoubleType>, [ArgSpec<LongDoubleType>, ArgSpec<LongDoubleType>]>,
410+
GuardedFunctionSpec<"fmodf128", RetValSpec<Float128Type>, [ArgSpec<Float128Type>, ArgSpec<Float128Type>], "LIBC_TYPES_HAS_FLOAT128">,
410411

411412
FunctionSpec<"frexp", RetValSpec<DoubleType>, [ArgSpec<DoubleType>, ArgSpec<IntPtr>]>,
412413
FunctionSpec<"frexpf", RetValSpec<FloatType>, [ArgSpec<FloatType>, ArgSpec<IntPtr>]>,

‎libc/src/__support/FPUtil/FPBits.h

Lines changed: 10 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -640,6 +640,7 @@ struct FPRepImpl : public FPRepSem<fp_type, RetT> {
640640
using UP::EXP_MASK;
641641
using UP::FRACTION_MASK;
642642
using UP::SIG_LEN;
643+
using UP::SIG_MASK;
643644
using UP::SIGN_MASK;
644645
LIBC_INLINE_VAR static constexpr int MAX_BIASED_EXPONENT =
645646
(1 << UP::EXP_LEN) - 1;
@@ -729,6 +730,9 @@ struct FPRepImpl : public FPRepSem<fp_type, RetT> {
729730
bits = UP::merge(bits, mantVal, FRACTION_MASK);
730731
}
731732

733+
LIBC_INLINE constexpr void set_significand(StorageType sigVal) {
734+
bits = UP::merge(bits, sigVal, SIG_MASK);
735+
}
732736
// Unsafe function to create a floating point representation.
733737
// It simply packs the sign, biased exponent and mantissa values without
734738
// checking bound nor normalization.
@@ -755,20 +759,19 @@ struct FPRepImpl : public FPRepSem<fp_type, RetT> {
755759
// 4) "number" zero value is not processed correctly.
756760
// 5) Number is unsigned, so the result can be only positive.
757761
LIBC_INLINE static constexpr RetT make_value(StorageType number, int ep) {
758-
static_assert(fp_type != FPType::X86_Binary80,
759-
"This function is not tested for X86 Extended Precision");
760-
FPRepImpl result;
761-
// offset: +1 for sign, but -1 for implicit first bit
762-
int lz = cpp::countl_zero(number) - UP::EXP_LEN;
762+
FPRepImpl result(0);
763+
int lz =
764+
UP::FRACTION_LEN + 1 - (UP::STORAGE_LEN - cpp::countl_zero(number));
765+
763766
number <<= lz;
764767
ep -= lz;
765768

766769
if (LIBC_LIKELY(ep >= 0)) {
767770
// Implicit number bit will be removed by mask
768-
result.set_mantissa(number);
771+
result.set_significand(number);
769772
result.set_biased_exponent(ep + 1);
770773
} else {
771-
result.set_mantissa(number >> -ep);
774+
result.set_significand(number >> -ep);
772775
}
773776
return RetT(result.uintval());
774777
}

‎libc/src/__support/FPUtil/generic/FMod.h

Lines changed: 56 additions & 87 deletions
Original file line numberDiff line numberDiff line change
@@ -117,63 +117,9 @@ namespace generic {
117117
// be implemented in another handler.
118118
// Signaling NaN converted to quiet NaN with FE_INVALID exception.
119119
// https://www.open-std.org/JTC1/SC22/WG14/www/docs/n1011.htm
120-
template <typename T> struct FModExceptionalInputHandler {
121-
122-
static_assert(cpp::is_floating_point_v<T>,
123-
"FModCStandardWrapper instantiated with invalid type.");
124-
125-
LIBC_INLINE static bool pre_check(T x, T y, T &out) {
126-
using FPB = fputil::FPBits<T>;
127-
const T quiet_nan = FPB::quiet_nan().get_val();
128-
FPB sx(x), sy(y);
129-
if (LIBC_LIKELY(!sy.is_zero() && !sy.is_inf_or_nan() &&
130-
!sx.is_inf_or_nan())) {
131-
return false;
132-
}
133-
134-
if (sx.is_nan() || sy.is_nan()) {
135-
if ((sx.is_nan() && !sx.is_quiet_nan()) ||
136-
(sy.is_nan() && !sy.is_quiet_nan()))
137-
fputil::raise_except_if_required(FE_INVALID);
138-
out = quiet_nan;
139-
return true;
140-
}
141-
142-
if (sx.is_inf() || sy.is_zero()) {
143-
fputil::raise_except_if_required(FE_INVALID);
144-
fputil::set_errno_if_required(EDOM);
145-
out = quiet_nan;
146-
return true;
147-
}
148-
149-
if (sy.is_inf()) {
150-
out = x;
151-
return true;
152-
}
153-
154-
// case where x == 0
155-
out = x;
156-
return true;
157-
}
158-
};
159-
160-
template <typename T> struct FModFastMathWrapper {
161-
162-
static_assert(cpp::is_floating_point_v<T>,
163-
"FModFastMathWrapper instantiated with invalid type.");
164-
165-
static bool pre_check(T, T, T &) { return false; }
166-
};
167-
168-
template <typename T> class FModDivisionSimpleHelper {
169-
private:
170-
using StorageType = typename FPBits<T>::StorageType;
171-
172-
public:
173-
LIBC_INLINE constexpr static StorageType execute(int exp_diff,
174-
int sides_zeroes_count,
175-
StorageType m_x,
176-
StorageType m_y) {
120+
template <typename T> struct FModDivisionSimpleHelper {
121+
LIBC_INLINE constexpr static T execute(int exp_diff, int sides_zeroes_count,
122+
T m_x, T m_y) {
177123
while (exp_diff > sides_zeroes_count) {
178124
exp_diff -= sides_zeroes_count;
179125
m_x <<= sides_zeroes_count;
@@ -185,28 +131,21 @@ template <typename T> class FModDivisionSimpleHelper {
185131
}
186132
};
187133

188-
template <typename T> class FModDivisionInvMultHelper {
189-
private:
190-
using FPB = FPBits<T>;
191-
using StorageType = typename FPB::StorageType;
192-
193-
public:
194-
LIBC_INLINE constexpr static StorageType execute(int exp_diff,
195-
int sides_zeroes_count,
196-
StorageType m_x,
197-
StorageType m_y) {
134+
template <typename T> struct FModDivisionInvMultHelper {
135+
LIBC_INLINE constexpr static T execute(int exp_diff, int sides_zeroes_count,
136+
T m_x, T m_y) {
137+
constexpr int LENGTH = sizeof(T) * CHAR_BIT;
198138
if (exp_diff > sides_zeroes_count) {
199-
StorageType inv_hy = (cpp::numeric_limits<StorageType>::max() / m_y);
139+
T inv_hy = (cpp::numeric_limits<T>::max() / m_y);
200140
while (exp_diff > sides_zeroes_count) {
201141
exp_diff -= sides_zeroes_count;
202-
StorageType hd =
203-
(m_x * inv_hy) >> (FPB::TOTAL_LEN - sides_zeroes_count);
142+
T hd = (m_x * inv_hy) >> (LENGTH - sides_zeroes_count);
204143
m_x <<= sides_zeroes_count;
205144
m_x -= hd * m_y;
206145
while (LIBC_UNLIKELY(m_x > m_y))
207146
m_x -= m_y;
208147
}
209-
StorageType hd = (m_x * inv_hy) >> (FPB::TOTAL_LEN - exp_diff);
148+
T hd = (m_x * inv_hy) >> (LENGTH - exp_diff);
210149
m_x <<= exp_diff;
211150
m_x -= hd * m_y;
212151
while (LIBC_UNLIKELY(m_x > m_y))
@@ -219,22 +158,49 @@ template <typename T> class FModDivisionInvMultHelper {
219158
}
220159
};
221160

222-
template <typename T, class Wrapper = FModExceptionalInputHandler<T>,
223-
class DivisionHelper = FModDivisionSimpleHelper<T>>
161+
template <typename T, typename U = typename FPBits<T>::StorageType,
162+
typename DivisionHelper = FModDivisionSimpleHelper<U>>
224163
class FMod {
225-
static_assert(cpp::is_floating_point_v<T>,
164+
static_assert(cpp::is_floating_point_v<T> && cpp::is_unsigned_v<U> &&
165+
(sizeof(U) * CHAR_BIT > FPBits<T>::FRACTION_LEN),
226166
"FMod instantiated with invalid type.");
227167

228168
private:
229169
using FPB = FPBits<T>;
230170
using StorageType = typename FPB::StorageType;
231171

172+
LIBC_INLINE static bool pre_check(T x, T y, T &out) {
173+
using FPB = fputil::FPBits<T>;
174+
const T quiet_nan = FPB::quiet_nan().get_val();
175+
FPB sx(x), sy(y);
176+
if (LIBC_LIKELY(!sy.is_zero() && !sy.is_inf_or_nan() &&
177+
!sx.is_inf_or_nan()))
178+
return false;
179+
180+
if (sx.is_nan() || sy.is_nan()) {
181+
if (sx.is_signaling_nan() || sy.is_signaling_nan())
182+
fputil::raise_except_if_required(FE_INVALID);
183+
out = quiet_nan;
184+
return true;
185+
}
186+
187+
if (sx.is_inf() || sy.is_zero()) {
188+
fputil::raise_except_if_required(FE_INVALID);
189+
fputil::set_errno_if_required(EDOM);
190+
out = quiet_nan;
191+
return true;
192+
}
193+
194+
out = x;
195+
return true;
196+
}
197+
232198
LIBC_INLINE static constexpr FPB eval_internal(FPB sx, FPB sy) {
233199

234200
if (LIBC_LIKELY(sx.uintval() <= sy.uintval())) {
235201
if (sx.uintval() < sy.uintval())
236202
return sx; // |x|<|y| return x
237-
return FPB(FPB::zero()); // |x|=|y| return 0.0
203+
return FPB::zero(); // |x|=|y| return 0.0
238204
}
239205

240206
int e_x = sx.get_biased_exponent();
@@ -247,27 +213,29 @@ class FMod {
247213
StorageType m_y = sy.get_explicit_mantissa();
248214
StorageType d = (e_x == e_y) ? (m_x - m_y) : (m_x << (e_x - e_y)) % m_y;
249215
if (d == 0)
250-
return FPB(FPB::zero());
216+
return FPB::zero();
251217
// iy - 1 because of "zero power" for number with power 1
252218
return FPB::make_value(d, e_y - 1);
253219
}
254-
/* Both subnormal special case. */
220+
// Both subnormal special case.
255221
if (LIBC_UNLIKELY(e_x == 0 && e_y == 0)) {
256222
FPB d;
257223
d.set_mantissa(sx.uintval() % sy.uintval());
258224
return d;
259225
}
260226

261227
// Note that hx is not subnormal by conditions above.
262-
StorageType m_x = sx.get_explicit_mantissa();
228+
U m_x = static_cast<U>(sx.get_explicit_mantissa());
263229
e_x--;
264230

265-
StorageType m_y = sy.get_explicit_mantissa();
266-
int lead_zeros_m_y = FPB::EXP_LEN;
231+
U m_y = static_cast<U>(sy.get_explicit_mantissa());
232+
constexpr int DEFAULT_LEAD_ZEROS =
233+
sizeof(U) * CHAR_BIT - FPB::FRACTION_LEN - 1;
234+
int lead_zeros_m_y = DEFAULT_LEAD_ZEROS;
267235
if (LIBC_LIKELY(e_y > 0)) {
268236
e_y--;
269237
} else {
270-
m_y = sy.get_mantissa();
238+
m_y = static_cast<U>(sy.get_mantissa());
271239
lead_zeros_m_y = cpp::countl_zero(m_y);
272240
}
273241

@@ -286,26 +254,27 @@ class FMod {
286254

287255
{
288256
// Shift hx left until the end or n = 0
289-
int left_shift = exp_diff < int(FPB::EXP_LEN) ? exp_diff : FPB::EXP_LEN;
257+
int left_shift =
258+
exp_diff < DEFAULT_LEAD_ZEROS ? exp_diff : DEFAULT_LEAD_ZEROS;
290259
m_x <<= left_shift;
291260
exp_diff -= left_shift;
292261
}
293262

294263
m_x %= m_y;
295264
if (LIBC_UNLIKELY(m_x == 0))
296-
return FPB(FPB::zero());
265+
return FPB::zero();
297266

298267
if (exp_diff == 0)
299-
return FPB::make_value(m_x, e_y);
268+
return FPB::make_value(static_cast<StorageType>(m_x), e_y);
300269

301-
/* hx next can't be 0, because hx < hy, hy % 2 == 1 hx * 2^i % hy != 0 */
270+
// hx next can't be 0, because hx < hy, hy % 2 == 1 hx * 2^i % hy != 0
302271
m_x = DivisionHelper::execute(exp_diff, sides_zeroes_count, m_x, m_y);
303-
return FPB::make_value(m_x, e_y);
272+
return FPB::make_value(static_cast<StorageType>(m_x), e_y);
304273
}
305274

306275
public:
307276
LIBC_INLINE static T eval(T x, T y) {
308-
if (T out; Wrapper::pre_check(x, y, out))
277+
if (T out; LIBC_UNLIKELY(pre_check(x, y, out)))
309278
return out;
310279
FPB sx(x), sy(y);
311280
Sign sign = sx.sign();

‎libc/src/math/CMakeLists.txt

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -119,6 +119,8 @@ add_math_entrypoint_object(fminf128)
119119

120120
add_math_entrypoint_object(fmod)
121121
add_math_entrypoint_object(fmodf)
122+
add_math_entrypoint_object(fmodl)
123+
add_math_entrypoint_object(fmodf128)
122124

123125
add_math_entrypoint_object(frexp)
124126
add_math_entrypoint_object(frexpf)

‎libc/src/math/fmodf128.h

Lines changed: 20 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,20 @@
1+
//===-- Implementation header for fmodf128 ----------------------*- 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_MATH_FMODF128_H
10+
#define LLVM_LIBC_SRC_MATH_FMODF128_H
11+
12+
#include "src/__support/macros/properties/types.h"
13+
14+
namespace LIBC_NAMESPACE {
15+
16+
float128 fmodf128(float128 x, float128 y);
17+
18+
} // namespace LIBC_NAMESPACE
19+
20+
#endif // LLVM_LIBC_SRC_MATH_FMODF128_H

‎libc/src/math/fmodl.h

Lines changed: 18 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,18 @@
1+
//===-- Implementation header for fmodl -------------------------*- 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_MATH_FMODL_H
10+
#define LLVM_LIBC_SRC_MATH_FMODL_H
11+
12+
namespace LIBC_NAMESPACE {
13+
14+
long double fmodl(long double x, long double y);
15+
16+
} // namespace LIBC_NAMESPACE
17+
18+
#endif // LLVM_LIBC_SRC_MATH_FMODL_H

0 commit comments

Comments
 (0)