Skip to content

[libc][math] Refactor acoshf16 implementation to header-only in src/__support/math folder. #148568

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
merged 1 commit into from
Jul 21, 2025
Merged

[libc][math] Refactor acoshf16 implementation to header-only in src/__support/math folder. #148568

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
1 change: 1 addition & 0 deletions libc/shared/math.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -15,6 +15,7 @@
#include "math/acosf.h"
#include "math/acosf16.h"
#include "math/acoshf.h"
#include "math/acoshf16.h"
#include "math/erff.h"
#include "math/exp.h"
#include "math/exp10.h"
Expand Down
29 changes: 29 additions & 0 deletions libc/shared/math/acoshf16.h
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1,29 @@
//===-- Shared acoshf16 function --------------------------------*- C++ -*-===//
//
// Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions.
// See https://llvm.org/LICENSE.txt for license information.
// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
//
//===----------------------------------------------------------------------===//

#ifndef LLVM_LIBC_SHARED_MATH_ACOSHF16_H
#define LLVM_LIBC_SHARED_MATH_ACOSHF16_H

#include "include/llvm-libc-macros/float16-macros.h"
#include "shared/libc_common.h"

#ifdef LIBC_TYPES_HAS_FLOAT16

#include "src/__support/math/acoshf16.h"

namespace LIBC_NAMESPACE_DECL {
namespace shared {

using math::acoshf16;

} // namespace shared
} // namespace LIBC_NAMESPACE_DECL

#endif // LIBC_TYPES_HAS_FLOAT16

#endif // LLVM_LIBC_SHARED_MATH_ACOSHF16_H
16 changes: 16 additions & 0 deletions libc/src/__support/math/CMakeLists.txt
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -79,6 +79,22 @@ add_header_library(
libc.src.__support.macros.optimization
)

add_header_library(
acoshf16
HDRS
acoshf16.h
DEPENDS
.acoshf_utils
libc.src.__support.FPUtil.cast
libc.src.__support.FPUtil.except_value_utils
libc.src.__support.FPUtil.fenv_impl
libc.src.__support.FPUtil.fp_bits
libc.src.__support.FPUtil.multiply_add
libc.src.__support.FPUtil.polyeval
libc.src.__support.FPUtil.sqrt
libc.src.__support.macros.optimization
)

add_header_library(
asin_utils
HDRS
Expand Down
123 changes: 123 additions & 0 deletions libc/src/__support/math/acoshf16.h
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1,123 @@
//===-- Implementation header for acoshf16 ----------------------*- C++ -*-===//
//
// Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions.
// See https://llvm.org/LICENSE.txt for license information.
// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
//
//===----------------------------------------------------------------------===//

#ifndef LLVM_LIBC_SRC___SUPPORT_MATH_ACOSHF_H
#define LLVM_LIBC_SRC___SUPPORT_MATH_ACOSHF_H

#include "include/llvm-libc-macros/float16-macros.h"

#ifdef LIBC_TYPES_HAS_FLOAT16

#include "acoshf_utils.h"
#include "src/__support/FPUtil/FEnvImpl.h"
#include "src/__support/FPUtil/FPBits.h"
#include "src/__support/FPUtil/PolyEval.h"
#include "src/__support/FPUtil/cast.h"
#include "src/__support/FPUtil/except_value_utils.h"
#include "src/__support/FPUtil/multiply_add.h"
#include "src/__support/FPUtil/sqrt.h"
#include "src/__support/macros/config.h"
#include "src/__support/macros/optimization.h"

namespace LIBC_NAMESPACE_DECL {

namespace math {

static constexpr float16 acoshf16(float16 x) {

using namespace acoshf_internal;
constexpr size_t N_EXCEPTS = 2;
constexpr fputil::ExceptValues<float16, N_EXCEPTS> ACOSHF16_EXCEPTS{{
// (input, RZ output, RU offset, RD offset, RN offset)
// x = 0x1.6dcp+1, acoshf16(x) = 0x1.b6p+0 (RZ)
{0x41B7, 0x3ED8, 1, 0, 0},
// x = 0x1.39p+0, acoshf16(x) = 0x1.4f8p-1 (RZ)
{0x3CE4, 0x393E, 1, 0, 1},
}};

using FPBits = fputil::FPBits<float16>;
FPBits xbits(x);
uint16_t x_u = xbits.uintval();

// Check for NaN input first.
if (LIBC_UNLIKELY(xbits.is_inf_or_nan())) {
if (xbits.is_signaling_nan()) {
fputil::raise_except_if_required(FE_INVALID);
return FPBits::quiet_nan().get_val();
}
if (xbits.is_neg()) {
fputil::set_errno_if_required(EDOM);
fputil::raise_except_if_required(FE_INVALID);
return FPBits::quiet_nan().get_val();
}
return x;
}

// Domain error for inputs less than 1.0.
if (LIBC_UNLIKELY(x <= 1.0f)) {
if (x == 1.0f)
return FPBits::zero().get_val();
fputil::set_errno_if_required(EDOM);
fputil::raise_except_if_required(FE_INVALID);
return FPBits::quiet_nan().get_val();
}

if (auto r = ACOSHF16_EXCEPTS.lookup(xbits.uintval());
LIBC_UNLIKELY(r.has_value()))
return r.value();

float xf = x;
// High-precision polynomial approximation for inputs close to 1.0
// ([1, 1.25)).
//
// Brief derivation:
// 1. Expand acosh(1 + delta) using Taylor series around delta=0:
// acosh(1 + delta) ≈ sqrt(2 * delta) * [1 - delta/12 + 3*delta^2/160
// - 5*delta^3/896 + 35*delta^4/18432 + ...]
// 2. Truncate the series to fit accurately for delta in [0, 0.25].
// 3. Polynomial coefficients (from sollya) used here are:
// P(delta) ≈ 1 - 0x1.555556p-4 * delta + 0x1.333334p-6 * delta^2
// - 0x1.6db6dcp-8 * delta^3 + 0x1.f1c71cp-10 * delta^4
// 4. The Sollya commands used to generate these coefficients were:
// > display = hexadecimal;
// > round(1/12, SG, RN);
// > round(3/160, SG, RN);
// > round(5/896, SG, RN);
// > round(35/18432, SG, RN);
// With hexadecimal display mode enabled, the outputs were:
// 0x1.555556p-4
// 0x1.333334p-6
// 0x1.6db6dcp-8
// 0x1.f1c71cp-10
// 5. The maximum absolute error, estimated using:
// dirtyinfnorm(acosh(1 + x) - sqrt(2*x) * P(x), [0, 0.25])
// is:
// 0x1.d84281p-22
if (LIBC_UNLIKELY(x_u < 0x3D00U)) {
float delta = xf - 1.0f;
float sqrt_2_delta = fputil::sqrt<float>(2.0 * delta);
float pe = fputil::polyeval(delta, 0x1p+0f, -0x1.555556p-4f, 0x1.333334p-6f,
-0x1.6db6dcp-8f, 0x1.f1c71cp-10f);
float approx = sqrt_2_delta * pe;
return fputil::cast<float16>(approx);
}

// acosh(x) = log(x + sqrt(x^2 - 1))
float sqrt_term = fputil::sqrt<float>(fputil::multiply_add(xf, xf, -1.0f));
float result = static_cast<float>(log_eval(xf + sqrt_term));

return fputil::cast<float16>(result);
}

} // namespace math

} // namespace LIBC_NAMESPACE_DECL

#endif // LIBC_TYPES_HAS_FLOAT16

#endif // LLVM_LIBC_SRC___SUPPORT_MATH_ACOSHF_H
14 changes: 2 additions & 12 deletions libc/src/math/generic/CMakeLists.txt
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -3878,18 +3878,8 @@ add_entrypoint_object(
HDRS
../acoshf16.h
DEPENDS
.explogxf
libc.hdr.errno_macros
libc.hdr.fenv_macros
libc.src.__support.FPUtil.cast
libc.src.__support.FPUtil.except_value_utils
libc.src.__support.FPUtil.fenv_impl
libc.src.__support.FPUtil.fp_bits
libc.src.__support.FPUtil.multiply_add
libc.src.__support.FPUtil.polyeval
libc.src.__support.FPUtil.sqrt
libc.src.__support.macros.optimization
libc.src.__support.macros.properties.types
libc.src.__support.math.acoshf16
libc.src.errno.errno
)

add_entrypoint_object(
Expand Down
99 changes: 2 additions & 97 deletions libc/src/math/generic/acoshf16.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -7,105 +7,10 @@
//===----------------------------------------------------------------------===//

#include "src/math/acoshf16.h"
#include "explogxf.h"
#include "hdr/errno_macros.h"
#include "hdr/fenv_macros.h"
#include "src/__support/FPUtil/FEnvImpl.h"
#include "src/__support/FPUtil/FPBits.h"
#include "src/__support/FPUtil/PolyEval.h"
#include "src/__support/FPUtil/cast.h"
#include "src/__support/FPUtil/except_value_utils.h"
#include "src/__support/FPUtil/multiply_add.h"
#include "src/__support/FPUtil/sqrt.h"
#include "src/__support/common.h"
#include "src/__support/macros/config.h"
#include "src/__support/macros/optimization.h"
#include "src/__support/math/acoshf16.h"

namespace LIBC_NAMESPACE_DECL {

static constexpr size_t N_EXCEPTS = 2;
static constexpr fputil::ExceptValues<float16, N_EXCEPTS> ACOSHF16_EXCEPTS{{
// (input, RZ output, RU offset, RD offset, RN offset)
// x = 0x1.6dcp+1, acoshf16(x) = 0x1.b6p+0 (RZ)
{0x41B7, 0x3ED8, 1, 0, 0},
// x = 0x1.39p+0, acoshf16(x) = 0x1.4f8p-1 (RZ)
{0x3CE4, 0x393E, 1, 0, 1},
}};

LLVM_LIBC_FUNCTION(float16, acoshf16, (float16 x)) {
using namespace acoshf_internal;
using FPBits = fputil::FPBits<float16>;
FPBits xbits(x);
uint16_t x_u = xbits.uintval();

// Check for NaN input first.
if (LIBC_UNLIKELY(xbits.is_inf_or_nan())) {
if (xbits.is_signaling_nan()) {
fputil::raise_except_if_required(FE_INVALID);
return FPBits::quiet_nan().get_val();
}
if (xbits.is_neg()) {
fputil::set_errno_if_required(EDOM);
fputil::raise_except_if_required(FE_INVALID);
return FPBits::quiet_nan().get_val();
}
return x;
}

// Domain error for inputs less than 1.0.
if (LIBC_UNLIKELY(x <= 1.0f)) {
if (x == 1.0f)
return FPBits::zero().get_val();
fputil::set_errno_if_required(EDOM);
fputil::raise_except_if_required(FE_INVALID);
return FPBits::quiet_nan().get_val();
}

if (auto r = ACOSHF16_EXCEPTS.lookup(xbits.uintval());
LIBC_UNLIKELY(r.has_value()))
return r.value();

float xf = x;
// High-precision polynomial approximation for inputs close to 1.0
// ([1, 1.25)).
//
// Brief derivation:
// 1. Expand acosh(1 + delta) using Taylor series around delta=0:
// acosh(1 + delta) ≈ sqrt(2 * delta) * [1 - delta/12 + 3*delta^2/160
// - 5*delta^3/896 + 35*delta^4/18432 + ...]
// 2. Truncate the series to fit accurately for delta in [0, 0.25].
// 3. Polynomial coefficients (from sollya) used here are:
// P(delta) ≈ 1 - 0x1.555556p-4 * delta + 0x1.333334p-6 * delta^2
// - 0x1.6db6dcp-8 * delta^3 + 0x1.f1c71cp-10 * delta^4
// 4. The Sollya commands used to generate these coefficients were:
// > display = hexadecimal;
// > round(1/12, SG, RN);
// > round(3/160, SG, RN);
// > round(5/896, SG, RN);
// > round(35/18432, SG, RN);
// With hexadecimal display mode enabled, the outputs were:
// 0x1.555556p-4
// 0x1.333334p-6
// 0x1.6db6dcp-8
// 0x1.f1c71cp-10
// 5. The maximum absolute error, estimated using:
// dirtyinfnorm(acosh(1 + x) - sqrt(2*x) * P(x), [0, 0.25])
// is:
// 0x1.d84281p-22
if (LIBC_UNLIKELY(x_u < 0x3D00U)) {
float delta = xf - 1.0f;
float sqrt_2_delta = fputil::sqrt<float>(2.0 * delta);
float pe = fputil::polyeval(delta, 0x1p+0f, -0x1.555556p-4f, 0x1.333334p-6f,
-0x1.6db6dcp-8f, 0x1.f1c71cp-10f);
float approx = sqrt_2_delta * pe;
return fputil::cast<float16>(approx);
}

// acosh(x) = log(x + sqrt(x^2 - 1))
float sqrt_term = fputil::sqrt<float>(fputil::multiply_add(xf, xf, -1.0f));
float result = static_cast<float>(log_eval(xf + sqrt_term));

return fputil::cast<float16>(result);
}
LLVM_LIBC_FUNCTION(float16, acoshf16, (float16 x)) { return math::acoshf16(x); }

} // namespace LIBC_NAMESPACE_DECL
1 change: 1 addition & 0 deletions libc/test/shared/CMakeLists.txt
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -11,6 +11,7 @@ add_fp_unittest(
libc.src.__support.math.acosf
libc.src.__support.math.acosf16
libc.src.__support.math.acoshf
libc.src.__support.math.acoshf16
libc.src.__support.math.erff
libc.src.__support.math.exp
libc.src.__support.math.exp10
Expand Down
2 changes: 2 additions & 0 deletions libc/test/shared/shared_math_test.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -14,6 +14,8 @@
TEST(LlvmLibcSharedMathTest, AllFloat16) {
int exponent;

EXPECT_FP_EQ(0x0p+0f, LIBC_NAMESPACE::shared::acoshf16(1.0f));

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

EXPECT_FP_EQ(0x1p+0f16, LIBC_NAMESPACE::shared::expf16(0.0f16));
Expand Down
24 changes: 24 additions & 0 deletions utils/bazel/llvm-project-overlay/libc/BUILD.bazel
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -2141,6 +2141,22 @@ libc_support_library(
],
)

libc_support_library(
name = "__support_math_acoshf16",
hdrs = ["src/__support/math/acoshf16.h"],
deps = [
":__support_math_acoshf_utils",
":__support_fputil_cast",
":__support_fputil_except_value_utils",
":__support_fputil_fenv_impl",
":__support_fputil_fp_bits",
":__support_fputil_multiply_add",
":__support_fputil_polyeval",
":__support_fputil_sqrt",
":__support_macros_optimization",
],
)

libc_support_library(
name = "__support_math_asin_utils",
hdrs = ["src/__support/math/asin_utils.h"],
Expand Down Expand Up @@ -2682,6 +2698,14 @@ libc_math_function(
],
)

libc_math_function(
name = "acoshf16",
additional_deps = [
":__support_math_acoshf16",
":errno",
],
)

libc_math_function(
name = "asinf",
additional_deps = [
Expand Down
Loading