Skip to content

Commit 520398e

Browse files
bassiounixbassiounix
authored
[libc][math] Refactor acoshf16 implementation to header-only in src/__support/math folder. (llvm#148568)
Part of llvm#147386 in preparation for: https://discourse.llvm.org/t/rfc-make-clang-builtin-math-functions-constexpr-with-llvm-libc-to-support-c-23-constexpr-math-functions/86450
1 parent dcffa3d commit 520398e

File tree

9 files changed

+200
-109
lines changed

9 files changed

+200
-109
lines changed

‎libc/shared/math.h

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -15,6 +15,7 @@
1515
#include "math/acosf.h"
1616
#include "math/acosf16.h"
1717
#include "math/acoshf.h"
18+
#include "math/acoshf16.h"
1819
#include "math/erff.h"
1920
#include "math/exp.h"
2021
#include "math/exp10.h"

‎libc/shared/math/acoshf16.h

Lines changed: 29 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,29 @@
1+
//===-- Shared acoshf16 function --------------------------------*- 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_SHARED_MATH_ACOSHF16_H
10+
#define LLVM_LIBC_SHARED_MATH_ACOSHF16_H
11+
12+
#include "include/llvm-libc-macros/float16-macros.h"
13+
#include "shared/libc_common.h"
14+
15+
#ifdef LIBC_TYPES_HAS_FLOAT16
16+
17+
#include "src/__support/math/acoshf16.h"
18+
19+
namespace LIBC_NAMESPACE_DECL {
20+
namespace shared {
21+
22+
using math::acoshf16;
23+
24+
} // namespace shared
25+
} // namespace LIBC_NAMESPACE_DECL
26+
27+
#endif // LIBC_TYPES_HAS_FLOAT16
28+
29+
#endif // LLVM_LIBC_SHARED_MATH_ACOSHF16_H

‎libc/src/__support/math/CMakeLists.txt

Lines changed: 16 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -79,6 +79,22 @@ add_header_library(
7979
libc.src.__support.macros.optimization
8080
)
8181

82+
add_header_library(
83+
acoshf16
84+
HDRS
85+
acoshf16.h
86+
DEPENDS
87+
.acoshf_utils
88+
libc.src.__support.FPUtil.cast
89+
libc.src.__support.FPUtil.except_value_utils
90+
libc.src.__support.FPUtil.fenv_impl
91+
libc.src.__support.FPUtil.fp_bits
92+
libc.src.__support.FPUtil.multiply_add
93+
libc.src.__support.FPUtil.polyeval
94+
libc.src.__support.FPUtil.sqrt
95+
libc.src.__support.macros.optimization
96+
)
97+
8298
add_header_library(
8399
asin_utils
84100
HDRS

‎libc/src/__support/math/acoshf16.h

Lines changed: 123 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,123 @@
1+
//===-- Implementation header for acoshf16 ----------------------*- 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_MATH_ACOSHF_H
10+
#define LLVM_LIBC_SRC___SUPPORT_MATH_ACOSHF_H
11+
12+
#include "include/llvm-libc-macros/float16-macros.h"
13+
14+
#ifdef LIBC_TYPES_HAS_FLOAT16
15+
16+
#include "acoshf_utils.h"
17+
#include "src/__support/FPUtil/FEnvImpl.h"
18+
#include "src/__support/FPUtil/FPBits.h"
19+
#include "src/__support/FPUtil/PolyEval.h"
20+
#include "src/__support/FPUtil/cast.h"
21+
#include "src/__support/FPUtil/except_value_utils.h"
22+
#include "src/__support/FPUtil/multiply_add.h"
23+
#include "src/__support/FPUtil/sqrt.h"
24+
#include "src/__support/macros/config.h"
25+
#include "src/__support/macros/optimization.h"
26+
27+
namespace LIBC_NAMESPACE_DECL {
28+
29+
namespace math {
30+
31+
static constexpr float16 acoshf16(float16 x) {
32+
33+
using namespace acoshf_internal;
34+
constexpr size_t N_EXCEPTS = 2;
35+
constexpr fputil::ExceptValues<float16, N_EXCEPTS> ACOSHF16_EXCEPTS{{
36+
// (input, RZ output, RU offset, RD offset, RN offset)
37+
// x = 0x1.6dcp+1, acoshf16(x) = 0x1.b6p+0 (RZ)
38+
{0x41B7, 0x3ED8, 1, 0, 0},
39+
// x = 0x1.39p+0, acoshf16(x) = 0x1.4f8p-1 (RZ)
40+
{0x3CE4, 0x393E, 1, 0, 1},
41+
}};
42+
43+
using FPBits = fputil::FPBits<float16>;
44+
FPBits xbits(x);
45+
uint16_t x_u = xbits.uintval();
46+
47+
// Check for NaN input first.
48+
if (LIBC_UNLIKELY(xbits.is_inf_or_nan())) {
49+
if (xbits.is_signaling_nan()) {
50+
fputil::raise_except_if_required(FE_INVALID);
51+
return FPBits::quiet_nan().get_val();
52+
}
53+
if (xbits.is_neg()) {
54+
fputil::set_errno_if_required(EDOM);
55+
fputil::raise_except_if_required(FE_INVALID);
56+
return FPBits::quiet_nan().get_val();
57+
}
58+
return x;
59+
}
60+
61+
// Domain error for inputs less than 1.0.
62+
if (LIBC_UNLIKELY(x <= 1.0f)) {
63+
if (x == 1.0f)
64+
return FPBits::zero().get_val();
65+
fputil::set_errno_if_required(EDOM);
66+
fputil::raise_except_if_required(FE_INVALID);
67+
return FPBits::quiet_nan().get_val();
68+
}
69+
70+
if (auto r = ACOSHF16_EXCEPTS.lookup(xbits.uintval());
71+
LIBC_UNLIKELY(r.has_value()))
72+
return r.value();
73+
74+
float xf = x;
75+
// High-precision polynomial approximation for inputs close to 1.0
76+
// ([1, 1.25)).
77+
//
78+
// Brief derivation:
79+
// 1. Expand acosh(1 + delta) using Taylor series around delta=0:
80+
// acosh(1 + delta) ≈ sqrt(2 * delta) * [1 - delta/12 + 3*delta^2/160
81+
// - 5*delta^3/896 + 35*delta^4/18432 + ...]
82+
// 2. Truncate the series to fit accurately for delta in [0, 0.25].
83+
// 3. Polynomial coefficients (from sollya) used here are:
84+
// P(delta) ≈ 1 - 0x1.555556p-4 * delta + 0x1.333334p-6 * delta^2
85+
// - 0x1.6db6dcp-8 * delta^3 + 0x1.f1c71cp-10 * delta^4
86+
// 4. The Sollya commands used to generate these coefficients were:
87+
// > display = hexadecimal;
88+
// > round(1/12, SG, RN);
89+
// > round(3/160, SG, RN);
90+
// > round(5/896, SG, RN);
91+
// > round(35/18432, SG, RN);
92+
// With hexadecimal display mode enabled, the outputs were:
93+
// 0x1.555556p-4
94+
// 0x1.333334p-6
95+
// 0x1.6db6dcp-8
96+
// 0x1.f1c71cp-10
97+
// 5. The maximum absolute error, estimated using:
98+
// dirtyinfnorm(acosh(1 + x) - sqrt(2*x) * P(x), [0, 0.25])
99+
// is:
100+
// 0x1.d84281p-22
101+
if (LIBC_UNLIKELY(x_u < 0x3D00U)) {
102+
float delta = xf - 1.0f;
103+
float sqrt_2_delta = fputil::sqrt<float>(2.0 * delta);
104+
float pe = fputil::polyeval(delta, 0x1p+0f, -0x1.555556p-4f, 0x1.333334p-6f,
105+
-0x1.6db6dcp-8f, 0x1.f1c71cp-10f);
106+
float approx = sqrt_2_delta * pe;
107+
return fputil::cast<float16>(approx);
108+
}
109+
110+
// acosh(x) = log(x + sqrt(x^2 - 1))
111+
float sqrt_term = fputil::sqrt<float>(fputil::multiply_add(xf, xf, -1.0f));
112+
float result = static_cast<float>(log_eval(xf + sqrt_term));
113+
114+
return fputil::cast<float16>(result);
115+
}
116+
117+
} // namespace math
118+
119+
} // namespace LIBC_NAMESPACE_DECL
120+
121+
#endif // LIBC_TYPES_HAS_FLOAT16
122+
123+
#endif // LLVM_LIBC_SRC___SUPPORT_MATH_ACOSHF_H

‎libc/src/math/generic/CMakeLists.txt

Lines changed: 2 additions & 12 deletions
Original file line numberDiff line numberDiff line change
@@ -3878,18 +3878,8 @@ add_entrypoint_object(
38783878
HDRS
38793879
../acoshf16.h
38803880
DEPENDS
3881-
.explogxf
3882-
libc.hdr.errno_macros
3883-
libc.hdr.fenv_macros
3884-
libc.src.__support.FPUtil.cast
3885-
libc.src.__support.FPUtil.except_value_utils
3886-
libc.src.__support.FPUtil.fenv_impl
3887-
libc.src.__support.FPUtil.fp_bits
3888-
libc.src.__support.FPUtil.multiply_add
3889-
libc.src.__support.FPUtil.polyeval
3890-
libc.src.__support.FPUtil.sqrt
3891-
libc.src.__support.macros.optimization
3892-
libc.src.__support.macros.properties.types
3881+
libc.src.__support.math.acoshf16
3882+
libc.src.errno.errno
38933883
)
38943884

38953885
add_entrypoint_object(

‎libc/src/math/generic/acoshf16.cpp

Lines changed: 2 additions & 97 deletions
Original file line numberDiff line numberDiff line change
@@ -7,105 +7,10 @@
77
//===----------------------------------------------------------------------===//
88

99
#include "src/math/acoshf16.h"
10-
#include "explogxf.h"
11-
#include "hdr/errno_macros.h"
12-
#include "hdr/fenv_macros.h"
13-
#include "src/__support/FPUtil/FEnvImpl.h"
14-
#include "src/__support/FPUtil/FPBits.h"
15-
#include "src/__support/FPUtil/PolyEval.h"
16-
#include "src/__support/FPUtil/cast.h"
17-
#include "src/__support/FPUtil/except_value_utils.h"
18-
#include "src/__support/FPUtil/multiply_add.h"
19-
#include "src/__support/FPUtil/sqrt.h"
20-
#include "src/__support/common.h"
21-
#include "src/__support/macros/config.h"
22-
#include "src/__support/macros/optimization.h"
10+
#include "src/__support/math/acoshf16.h"
2311

2412
namespace LIBC_NAMESPACE_DECL {
2513

26-
static constexpr size_t N_EXCEPTS = 2;
27-
static constexpr fputil::ExceptValues<float16, N_EXCEPTS> ACOSHF16_EXCEPTS{{
28-
// (input, RZ output, RU offset, RD offset, RN offset)
29-
// x = 0x1.6dcp+1, acoshf16(x) = 0x1.b6p+0 (RZ)
30-
{0x41B7, 0x3ED8, 1, 0, 0},
31-
// x = 0x1.39p+0, acoshf16(x) = 0x1.4f8p-1 (RZ)
32-
{0x3CE4, 0x393E, 1, 0, 1},
33-
}};
34-
35-
LLVM_LIBC_FUNCTION(float16, acoshf16, (float16 x)) {
36-
using namespace acoshf_internal;
37-
using FPBits = fputil::FPBits<float16>;
38-
FPBits xbits(x);
39-
uint16_t x_u = xbits.uintval();
40-
41-
// Check for NaN input first.
42-
if (LIBC_UNLIKELY(xbits.is_inf_or_nan())) {
43-
if (xbits.is_signaling_nan()) {
44-
fputil::raise_except_if_required(FE_INVALID);
45-
return FPBits::quiet_nan().get_val();
46-
}
47-
if (xbits.is_neg()) {
48-
fputil::set_errno_if_required(EDOM);
49-
fputil::raise_except_if_required(FE_INVALID);
50-
return FPBits::quiet_nan().get_val();
51-
}
52-
return x;
53-
}
54-
55-
// Domain error for inputs less than 1.0.
56-
if (LIBC_UNLIKELY(x <= 1.0f)) {
57-
if (x == 1.0f)
58-
return FPBits::zero().get_val();
59-
fputil::set_errno_if_required(EDOM);
60-
fputil::raise_except_if_required(FE_INVALID);
61-
return FPBits::quiet_nan().get_val();
62-
}
63-
64-
if (auto r = ACOSHF16_EXCEPTS.lookup(xbits.uintval());
65-
LIBC_UNLIKELY(r.has_value()))
66-
return r.value();
67-
68-
float xf = x;
69-
// High-precision polynomial approximation for inputs close to 1.0
70-
// ([1, 1.25)).
71-
//
72-
// Brief derivation:
73-
// 1. Expand acosh(1 + delta) using Taylor series around delta=0:
74-
// acosh(1 + delta) ≈ sqrt(2 * delta) * [1 - delta/12 + 3*delta^2/160
75-
// - 5*delta^3/896 + 35*delta^4/18432 + ...]
76-
// 2. Truncate the series to fit accurately for delta in [0, 0.25].
77-
// 3. Polynomial coefficients (from sollya) used here are:
78-
// P(delta) ≈ 1 - 0x1.555556p-4 * delta + 0x1.333334p-6 * delta^2
79-
// - 0x1.6db6dcp-8 * delta^3 + 0x1.f1c71cp-10 * delta^4
80-
// 4. The Sollya commands used to generate these coefficients were:
81-
// > display = hexadecimal;
82-
// > round(1/12, SG, RN);
83-
// > round(3/160, SG, RN);
84-
// > round(5/896, SG, RN);
85-
// > round(35/18432, SG, RN);
86-
// With hexadecimal display mode enabled, the outputs were:
87-
// 0x1.555556p-4
88-
// 0x1.333334p-6
89-
// 0x1.6db6dcp-8
90-
// 0x1.f1c71cp-10
91-
// 5. The maximum absolute error, estimated using:
92-
// dirtyinfnorm(acosh(1 + x) - sqrt(2*x) * P(x), [0, 0.25])
93-
// is:
94-
// 0x1.d84281p-22
95-
if (LIBC_UNLIKELY(x_u < 0x3D00U)) {
96-
float delta = xf - 1.0f;
97-
float sqrt_2_delta = fputil::sqrt<float>(2.0 * delta);
98-
float pe = fputil::polyeval(delta, 0x1p+0f, -0x1.555556p-4f, 0x1.333334p-6f,
99-
-0x1.6db6dcp-8f, 0x1.f1c71cp-10f);
100-
float approx = sqrt_2_delta * pe;
101-
return fputil::cast<float16>(approx);
102-
}
103-
104-
// acosh(x) = log(x + sqrt(x^2 - 1))
105-
float sqrt_term = fputil::sqrt<float>(fputil::multiply_add(xf, xf, -1.0f));
106-
float result = static_cast<float>(log_eval(xf + sqrt_term));
107-
108-
return fputil::cast<float16>(result);
109-
}
14+
LLVM_LIBC_FUNCTION(float16, acoshf16, (float16 x)) { return math::acoshf16(x); }
11015

11116
} // namespace LIBC_NAMESPACE_DECL

‎libc/test/shared/CMakeLists.txt

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -11,6 +11,7 @@ add_fp_unittest(
1111
libc.src.__support.math.acosf
1212
libc.src.__support.math.acosf16
1313
libc.src.__support.math.acoshf
14+
libc.src.__support.math.acoshf16
1415
libc.src.__support.math.erff
1516
libc.src.__support.math.exp
1617
libc.src.__support.math.exp10

‎libc/test/shared/shared_math_test.cpp

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -14,6 +14,8 @@
1414
TEST(LlvmLibcSharedMathTest, AllFloat16) {
1515
int exponent;
1616

17+
EXPECT_FP_EQ(0x0p+0f, LIBC_NAMESPACE::shared::acoshf16(1.0f));
18+
1719
EXPECT_FP_EQ(0x1p+0f16, LIBC_NAMESPACE::shared::exp10f16(0.0f16));
1820

1921
EXPECT_FP_EQ(0x1p+0f16, LIBC_NAMESPACE::shared::expf16(0.0f16));

‎utils/bazel/llvm-project-overlay/libc/BUILD.bazel

Lines changed: 24 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -2141,6 +2141,22 @@ libc_support_library(
21412141
],
21422142
)
21432143

2144+
libc_support_library(
2145+
name = "__support_math_acoshf16",
2146+
hdrs = ["src/__support/math/acoshf16.h"],
2147+
deps = [
2148+
":__support_math_acoshf_utils",
2149+
":__support_fputil_cast",
2150+
":__support_fputil_except_value_utils",
2151+
":__support_fputil_fenv_impl",
2152+
":__support_fputil_fp_bits",
2153+
":__support_fputil_multiply_add",
2154+
":__support_fputil_polyeval",
2155+
":__support_fputil_sqrt",
2156+
":__support_macros_optimization",
2157+
],
2158+
)
2159+
21442160
libc_support_library(
21452161
name = "__support_math_asin_utils",
21462162
hdrs = ["src/__support/math/asin_utils.h"],
@@ -2682,6 +2698,14 @@ libc_math_function(
26822698
],
26832699
)
26842700

2701+
libc_math_function(
2702+
name = "acoshf16",
2703+
additional_deps = [
2704+
":__support_math_acoshf16",
2705+
":errno",
2706+
],
2707+
)
2708+
26852709
libc_math_function(
26862710
name = "asinf",
26872711
additional_deps = [

0 commit comments

Comments
 (0)