Skip to content

Commit 94ff9a7

Browse files
committed
[libc] Add float log_eval.
1 parent 81c5cf7 commit 94ff9a7

File tree

4 files changed

+120
-2
lines changed

4 files changed

+120
-2
lines changed

libc/src/math/generic/atanhf16.cpp

Lines changed: 6 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -22,13 +22,17 @@
2222

2323
namespace LIBC_NAMESPACE_DECL {
2424

25-
static constexpr size_t N_EXCEPTS = 2;
25+
static constexpr size_t N_EXCEPTS = 4;
2626
static constexpr fputil::ExceptValues<float16, N_EXCEPTS> ATANHF16_EXCEPTS{{
2727
// (input, RZ output, RU offset, RD offset, RN offset)
2828
// x = 0x1.a5cp-4, atanhf16(x) = 0x1.a74p-4 (RZ)
2929
{0x2E97, 0x2E9D, 1, 0, 0},
3030
// x = -0x1.a5cp-4, atanhf16(x) = -0x1.a74p-4 (RZ)
3131
{0xAE97, 0xAE9D, 0, 1, 0},
32+
// x = -0x1.99cp-4, atanhf16(x) = -0x1.9bp-4 (RZ)
33+
{0xAE67, 0xAE6C, 0, 1, 1},
34+
// x = -0x1.b8cp-3, atanhf16(x) = -0x1.bfcp-3 (RZ)
35+
{0xB2E3, 0xB2FF, 0, 1, 0},
3236
}};
3337

3438
LLVM_LIBC_FUNCTION(float16, atanhf16, (float16 x)) {
@@ -89,7 +93,7 @@ LLVM_LIBC_FUNCTION(float16, atanhf16, (float16 x)) {
8993
}
9094

9195
float xf = x;
92-
return fputil::cast<float16>(0.5 * log_eval((xf + 1.0f) / (xf - 1.0f)));
96+
return fputil::cast<float16>(0.5 * log_eval_f((xf + 1.0f) / (xf - 1.0f)));
9397
}
9498

9599
} // namespace LIBC_NAMESPACE_DECL

libc/src/math/generic/common_constants.cpp

Lines changed: 72 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -13,6 +13,42 @@
1313

1414
namespace LIBC_NAMESPACE_DECL {
1515

16+
// Lookup table for logf(f) = logf(1 + n*2^(-7)) where n = 0..127,
17+
// computed and stored as float precision constants.
18+
const float ONE_OVER_F_FLOAT[128] = {
19+
0x1p0f, 0x1.fc07fp-1f, 0x1.f81f82p-1f, 0x1.f4465ap-1f,
20+
0x1.f07c2p-1f, 0x1.ecc07cp-1f, 0x1.e9131ap-1f, 0x1.e573acp-1f,
21+
0x1.e1e1e2p-1f, 0x1.de5d6ep-1f, 0x1.dae608p-1f, 0x1.d77b66p-1f,
22+
0x1.d41d42p-1f, 0x1.d0cb58p-1f, 0x1.cd8568p-1f, 0x1.ca4b3p-1f,
23+
0x1.c71c72p-1f, 0x1.c3f8fp-1f, 0x1.c0e07p-1f, 0x1.bdd2b8p-1f,
24+
0x1.bacf92p-1f, 0x1.b7d6c4p-1f, 0x1.b4e81cp-1f, 0x1.b20364p-1f,
25+
0x1.af286cp-1f, 0x1.ac5702p-1f, 0x1.a98ef6p-1f, 0x1.a6d01ap-1f,
26+
0x1.a41a42p-1f, 0x1.a16d4p-1f, 0x1.9ec8eap-1f, 0x1.9c2d14p-1f,
27+
0x1.99999ap-1f, 0x1.970e5p-1f, 0x1.948b1p-1f, 0x1.920fb4p-1f,
28+
0x1.8f9c18p-1f, 0x1.8d3018p-1f, 0x1.8acb9p-1f, 0x1.886e6p-1f,
29+
0x1.861862p-1f, 0x1.83c978p-1f, 0x1.818182p-1f, 0x1.7f406p-1f,
30+
0x1.7d05f4p-1f, 0x1.7ad22p-1f, 0x1.78a4c8p-1f, 0x1.767dcep-1f,
31+
0x1.745d18p-1f, 0x1.724288p-1f, 0x1.702e06p-1f, 0x1.6e1f76p-1f,
32+
0x1.6c16c2p-1f, 0x1.6a13cep-1f, 0x1.681682p-1f, 0x1.661ec6p-1f,
33+
0x1.642c86p-1f, 0x1.623fa8p-1f, 0x1.605816p-1f, 0x1.5e75bcp-1f,
34+
0x1.5c9882p-1f, 0x1.5ac056p-1f, 0x1.58ed24p-1f, 0x1.571ed4p-1f,
35+
0x1.555556p-1f, 0x1.539094p-1f, 0x1.51d07ep-1f, 0x1.501502p-1f,
36+
0x1.4e5e0ap-1f, 0x1.4cab88p-1f, 0x1.4afd6ap-1f, 0x1.49539ep-1f,
37+
0x1.47ae14p-1f, 0x1.460cbcp-1f, 0x1.446f86p-1f, 0x1.42d662p-1f,
38+
0x1.414142p-1f, 0x1.3fb014p-1f, 0x1.3e22ccp-1f, 0x1.3c995ap-1f,
39+
0x1.3b13b2p-1f, 0x1.3991c2p-1f, 0x1.381382p-1f, 0x1.3698ep-1f,
40+
0x1.3521dp-1f, 0x1.33ae46p-1f, 0x1.323e34p-1f, 0x1.30d19p-1f,
41+
0x1.2f684cp-1f, 0x1.2e025cp-1f, 0x1.2c9fb4p-1f, 0x1.2b404ap-1f,
42+
0x1.29e412p-1f, 0x1.288b02p-1f, 0x1.27350cp-1f, 0x1.25e228p-1f,
43+
0x1.24924ap-1f, 0x1.234568p-1f, 0x1.21fb78p-1f, 0x1.20b47p-1f,
44+
0x1.1f7048p-1f, 0x1.1e2ef4p-1f, 0x1.1cf06ap-1f, 0x1.1bb4a4p-1f,
45+
0x1.1a7b96p-1f, 0x1.194538p-1f, 0x1.181182p-1f, 0x1.16e068p-1f,
46+
0x1.15b1e6p-1f, 0x1.1485fp-1f, 0x1.135c82p-1f, 0x1.12358ep-1f,
47+
0x1.111112p-1f, 0x1.0fef02p-1f, 0x1.0ecf56p-1f, 0x1.0db20ap-1f,
48+
0x1.0c9714p-1f, 0x1.0b7e6ep-1f, 0x1.0a681p-1f, 0x1.0953f4p-1f,
49+
0x1.08421p-1f, 0x1.07326p-1f, 0x1.0624dep-1f, 0x1.05198p-1f,
50+
0x1.041042p-1f, 0x1.03091cp-1f, 0x1.020408p-1f, 0x1.010102p-1f};
51+
1652
// Lookup table for (1/f) where f = 1 + n*2^(-7), n = 0..127.
1753
const double ONE_OVER_F[128] = {
1854
0x1.0000000000000p+0, 0x1.fc07f01fc07f0p-1, 0x1.f81f81f81f820p-1,
@@ -59,6 +95,42 @@ const double ONE_OVER_F[128] = {
5995
0x1.05197f7d73404p-1, 0x1.0410410410410p-1, 0x1.03091b51f5e1ap-1,
6096
0x1.0204081020408p-1, 0x1.0101010101010p-1};
6197

98+
// Lookup table for (1/f) where f = 1 + n*2^(-7), n = 0..127,
99+
// computed and stored as float precision constants.
100+
const float LOG_F_FLOAT[128] = {
101+
0.0f, 0x1.fe02a6p-8f, 0x1.fc0a8cp-7f, 0x1.7b91bp-6f,
102+
0x1.f829bp-6f, 0x1.39e87cp-5f, 0x1.77459p-5f, 0x1.b42dd8p-5f,
103+
0x1.f0a30cp-5f, 0x1.16536ep-4f, 0x1.341d7ap-4f, 0x1.51b074p-4f,
104+
0x1.6f0d28p-4f, 0x1.8c345ep-4f, 0x1.a926d4p-4f, 0x1.c5e548p-4f,
105+
0x1.e27076p-4f, 0x1.fec914p-4f, 0x1.0d77e8p-3f, 0x1.1b72aep-3f,
106+
0x1.29553p-3f, 0x1.371fc2p-3f, 0x1.44d2b6p-3f, 0x1.526e5ep-3f,
107+
0x1.5ff308p-3f, 0x1.6d60fep-3f, 0x1.7ab89p-3f, 0x1.87fa06p-3f,
108+
0x1.9525aap-3f, 0x1.a23bc2p-3f, 0x1.af3c94p-3f, 0x1.bc2868p-3f,
109+
0x1.c8ff7cp-3f, 0x1.d5c216p-3f, 0x1.e27076p-3f, 0x1.ef0adcp-3f,
110+
0x1.fb9186p-3f, 0x1.04025ap-2f, 0x1.0a324ep-2f, 0x1.1058cp-2f,
111+
0x1.1675cap-2f, 0x1.1c898cp-2f, 0x1.22942p-2f, 0x1.2895a2p-2f,
112+
0x1.2e8e2cp-2f, 0x1.347ddap-2f, 0x1.3a64c6p-2f, 0x1.404308p-2f,
113+
0x1.4618bcp-2f, 0x1.4be5fap-2f, 0x1.51aad8p-2f, 0x1.576772p-2f,
114+
0x1.5d1bdcp-2f, 0x1.62c83p-2f, 0x1.686c82p-2f, 0x1.6e08eap-2f,
115+
0x1.739d8p-2f, 0x1.792a56p-2f, 0x1.7eaf84p-2f, 0x1.842d1ep-2f,
116+
0x1.89a338p-2f, 0x1.8f11e8p-2f, 0x1.947942p-2f, 0x1.99d958p-2f,
117+
0x1.9f323ep-2f, 0x1.a4840ap-2f, 0x1.a9cecap-2f, 0x1.af1294p-2f,
118+
0x1.b44f78p-2f, 0x1.b9858ap-2f, 0x1.beb4dap-2f, 0x1.c3dd7ap-2f,
119+
0x1.c8ff7cp-2f, 0x1.ce1afp-2f, 0x1.d32fe8p-2f, 0x1.d83e72p-2f,
120+
0x1.dd46ap-2f, 0x1.e24882p-2f, 0x1.e74426p-2f, 0x1.ec399ep-2f,
121+
0x1.f128f6p-2f, 0x1.f6124p-2f, 0x1.faf588p-2f, 0x1.ffd2ep-2f,
122+
0x1.02552ap-1f, 0x1.04bdfap-1f, 0x1.0723e6p-1f, 0x1.0986f4p-1f,
123+
0x1.0be72ep-1f, 0x1.0e4498p-1f, 0x1.109f3ap-1f, 0x1.12f71ap-1f,
124+
0x1.154c3ep-1f, 0x1.179eacp-1f, 0x1.19ee6cp-1f, 0x1.1c3b82p-1f,
125+
0x1.1e85f6p-1f, 0x1.20cdcep-1f, 0x1.23130ep-1f, 0x1.2555bcp-1f,
126+
0x1.2795e2p-1f, 0x1.29d38p-1f, 0x1.2c0e9ep-1f, 0x1.2e4744p-1f,
127+
0x1.307d74p-1f, 0x1.32b134p-1f, 0x1.34e28ap-1f, 0x1.37117cp-1f,
128+
0x1.393e0ep-1f, 0x1.3b6844p-1f, 0x1.3d9026p-1f, 0x1.3fb5b8p-1f,
129+
0x1.41d8fep-1f, 0x1.43f9fep-1f, 0x1.4618bcp-1f, 0x1.48353ep-1f,
130+
0x1.4a4f86p-1f, 0x1.4c679ap-1f, 0x1.4e7d82p-1f, 0x1.50913cp-1f,
131+
0x1.52a2d2p-1f, 0x1.54b246p-1f, 0x1.56bf9ep-1f, 0x1.58cadcp-1f,
132+
0x1.5ad404p-1f, 0x1.5cdb1ep-1f, 0x1.5ee02ap-1f, 0x1.60e33p-1f};
133+
62134
// Lookup table for log(f) = log(1 + n*2^(-7)) where n = 0..127.
63135
const double LOG_F[128] = {
64136
0x0.0000000000000p+0, 0x1.fe02a6b106788p-8, 0x1.fc0a8b0fc03e3p-7,

libc/src/math/generic/common_constants.h

Lines changed: 8 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -15,9 +15,17 @@
1515

1616
namespace LIBC_NAMESPACE_DECL {
1717

18+
// Lookup table for (1/f) where f = 1 + n*2^(-7), n = 0..127,
19+
// computed and stored as float precision constants.
20+
extern const float ONE_OVER_F_FLOAT[128];
21+
1822
// Lookup table for (1/f) where f = 1 + n*2^(-7), n = 0..127.
1923
extern const double ONE_OVER_F[128];
2024

25+
// Lookup table for log(f) = log(1 + n*2^(-7)) where n = 0..127,
26+
// computed and stored as float precision constants.
27+
extern const float LOG_F_FLOAT[128];
28+
2129
// Lookup table for log(f) = log(1 + n*2^(-7)) where n = 0..127.
2230
extern const double LOG_F[128];
2331

libc/src/math/generic/explogxf.h

Lines changed: 34 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -297,6 +297,40 @@ LIBC_INLINE static double log2_eval(double x) {
297297
return result;
298298
}
299299

300+
// x should be positive, normal finite value
301+
LIBC_INLINE static float log_eval_f(float x) {
302+
// For x = 2^ex * (1 + mx), logf(x) = ex * logf(2) + logf(1 + mx).
303+
using FPB = fputil::FPBits<float>;
304+
FPB bs(x);
305+
306+
float ex = static_cast<float>(bs.get_exponent());
307+
// p1 is the leading 7 bits of mx, i.e.
308+
// p1 * 2^(-7) <= m_x < (p1 + 1) * 2^(-7).
309+
int p1 = static_cast<int>(bs.get_mantissa() >> (FPB::FRACTION_LEN - 7));
310+
311+
// Set bs to (1 + (mx - p1*2^(-7))
312+
bs.set_uintval(bs.uintval() & (FPB::FRACTION_MASK >> 7));
313+
bs.set_biased_exponent(FPB::EXP_BIAS);
314+
// dx = (mx - p1*2^(-7)) / (1 + p1*2^(-7)).
315+
float dx = (bs.get_val() - 1.0f) * ONE_OVER_F_FLOAT[p1];
316+
317+
// Minimax polynomial of log(1 + dx) generated by Sollya with:
318+
// > P = fpminimax(log(1 + x)/x, 6, [|D...|], [0, 2^-7]);
319+
const float COEFFS[6] = {-0x1.fffffep-2f, 0x1.555556p-2f, -0x1.fffefep-3f,
320+
0x1.99999ap-3f, -0x1.554318p-3f, 0x1.1dc5c4p-3f};
321+
322+
float dx2 = dx * dx;
323+
324+
float c1 = fputil::multiply_add(dx, COEFFS[1], COEFFS[0]);
325+
float c2 = fputil::multiply_add(dx, COEFFS[3], COEFFS[2]);
326+
float c3 = fputil::multiply_add(dx, COEFFS[5], COEFFS[4]);
327+
328+
float p = fputil::polyeval(dx2, dx, c1, c2, c3);
329+
330+
float result = fputil::multiply_add(ex, 0x1.62e42ep-1f, LOG_F_FLOAT[p1] + p);
331+
return result;
332+
}
333+
300334
// x should be positive, normal finite value
301335
LIBC_INLINE static double log_eval(double x) {
302336
// For x = 2^ex * (1 + mx)

0 commit comments

Comments
 (0)