Skip to content

Commit ecf4f95

Browse files
authored
[libc][math][c23] Add tanf16 function (#121018)
- Implementation of tan for 16-bit floating point inputs. - Exhaustive tests across the 16-bit input range
1 parent 7bf1cb7 commit ecf4f95

File tree

14 files changed

+275
-16
lines changed

14 files changed

+275
-16
lines changed

libc/config/linux/x86_64/entrypoints.txt

+1
Original file line numberDiff line numberDiff line change
@@ -721,6 +721,7 @@ if(LIBC_TYPES_HAS_FLOAT16)
721721
libc.src.math.sinhf16
722722
libc.src.math.sinpif16
723723
libc.src.math.sqrtf16
724+
libc.src.math.tanf16
724725
libc.src.math.tanhf16
725726
libc.src.math.tanpif16
726727
libc.src.math.totalorderf16

libc/docs/headers/math/index.rst

+1-1
Original file line numberDiff line numberDiff line change
@@ -346,7 +346,7 @@ Higher Math Functions
346346
+-----------+------------------+-----------------+------------------------+----------------------+------------------------+------------------------+----------------------------+
347347
| sqrt | |check| | |check| | |check| | |check| | |check| | 7.12.7.10 | F.10.4.10 |
348348
+-----------+------------------+-----------------+------------------------+----------------------+------------------------+------------------------+----------------------------+
349-
| tan | |check| | |check| | | | | 7.12.4.7 | F.10.1.7 |
349+
| tan | |check| | |check| | | |check| | | 7.12.4.7 | F.10.1.7 |
350350
+-----------+------------------+-----------------+------------------------+----------------------+------------------------+------------------------+----------------------------+
351351
| tanh | |check| | | | |check| | | 7.12.5.6 | F.10.2.6 |
352352
+-----------+------------------+-----------------+------------------------+----------------------+------------------------+------------------------+----------------------------+

libc/include/math.yaml

+7
Original file line numberDiff line numberDiff line change
@@ -2418,6 +2418,13 @@ functions:
24182418
return_type: float
24192419
arguments:
24202420
- type: float
2421+
- name: tanf16
2422+
standards:
2423+
- stdc
2424+
return_type: _Float16
2425+
arguments:
2426+
- type: _Float16
2427+
guard: LIBC_TYPES_HAS_FLOAT16
24212428
- name: tanhf
24222429
standards:
24232430
- stdc

libc/src/math/CMakeLists.txt

+1
Original file line numberDiff line numberDiff line change
@@ -501,6 +501,7 @@ add_math_entrypoint_object(sqrtf128)
501501

502502
add_math_entrypoint_object(tan)
503503
add_math_entrypoint_object(tanf)
504+
add_math_entrypoint_object(tanf16)
504505

505506
add_math_entrypoint_object(tanh)
506507
add_math_entrypoint_object(tanhf)

libc/src/math/generic/CMakeLists.txt

+19
Original file line numberDiff line numberDiff line change
@@ -661,6 +661,25 @@ add_entrypoint_object(
661661
${libc_opt_high_flag}
662662
)
663663

664+
add_entrypoint_object(
665+
tanf16
666+
SRCS
667+
tanf16.cpp
668+
HDRS
669+
../tanf16.h
670+
DEPENDS
671+
.sincosf16_utils
672+
libc.hdr.errno_macros
673+
libc.hdr.fenv_macros
674+
libc.src.__support.FPUtil.cast
675+
libc.src.__support.FPUtil.fenv_impl
676+
libc.src.__support.FPUtil.fp_bits
677+
libc.src.__support.FPUtil.except_value_utils
678+
libc.src.__support.FPUtil.multiply_add
679+
libc.src.__support.macros.optimization
680+
libc.src.__support.macros.properties.types
681+
)
682+
664683
add_entrypoint_object(
665684
tanpif16
666685
SRCS

libc/src/math/generic/sincosf16_utils.h

+12-13
Original file line numberDiff line numberDiff line change
@@ -47,24 +47,23 @@ LIBC_INLINE int32_t range_reduction_sincospif16(float x, float &y) {
4747

4848
// Recall, range reduction:
4949
// k = round(x * 32/pi)
50-
// y = x * 32/pi - k
5150
//
52-
// The constant 0x1.45f306dc9c883p3 is 32/pi rounded to double-precision.
53-
// 32/pi is generated by Sollya with the following commands:
54-
// > display = hexadecimal;
55-
// > round(32/pi, D, RN);
56-
//
57-
// The precision choice of 'double' is to minimize rounding errors
58-
// in this initial scaling step, preserving enough bits so errors accumulated
59-
// while computing the subtraction: y = x * 32/pi - round(x * 32/pi)
51+
// The precision choice of 'double' in the following function is to minimize
52+
// rounding errors in this initial scaling step,
53+
// preserving enough bits so errors accumulated while computing the subtraction:
54+
// y = x * 32/pi - round(x * 32/pi)
6055
// are beyond the least-significant bit of single-precision used during
6156
// further intermediate computation.
6257
LIBC_INLINE int32_t range_reduction_sincosf16(float x, float &y) {
63-
double prod = x * 0x1.45f306dc9c883p3;
64-
double kf = fputil::nearest_integer(prod);
65-
y = static_cast<float>(prod - kf);
58+
// Generated by Sollya with:
59+
// > D(32/pi);
60+
constexpr double THIRTYTWO_OVER_PI = 0x1.45f306dc9c883p3;
6661

67-
return static_cast<int32_t>(kf);
62+
double prod = x * THIRTYTWO_OVER_PI;
63+
double kd = fputil::nearest_integer(prod);
64+
y = static_cast<float>(prod - kd);
65+
66+
return static_cast<int32_t>(kd);
6867
}
6968

7069
static LIBC_INLINE void sincosf16_poly_eval(int32_t k, float y, float &sin_k,

libc/src/math/generic/tanf16.cpp

+115
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,115 @@
1+
//===-- Half-precision tan(x) function ------------------------------------===//
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+
#include "src/math/tanf16.h"
10+
#include "hdr/errno_macros.h"
11+
#include "hdr/fenv_macros.h"
12+
#include "sincosf16_utils.h"
13+
#include "src/__support/FPUtil/FEnvImpl.h"
14+
#include "src/__support/FPUtil/FPBits.h"
15+
#include "src/__support/FPUtil/cast.h"
16+
#include "src/__support/FPUtil/except_value_utils.h"
17+
#include "src/__support/FPUtil/multiply_add.h"
18+
#include "src/__support/macros/optimization.h"
19+
20+
namespace LIBC_NAMESPACE_DECL {
21+
22+
constexpr size_t N_EXCEPTS = 9;
23+
24+
constexpr fputil::ExceptValues<float16, N_EXCEPTS> TANF16_EXCEPTS{{
25+
// (input, RZ output, RU offset, RD offset, RN offset)
26+
{0x2894, 0x2894, 1, 0, 1},
27+
{0x3091, 0x3099, 1, 0, 0},
28+
{0x3098, 0x30a0, 1, 0, 0},
29+
{0x55ed, 0x3911, 1, 0, 0},
30+
{0x607b, 0xc638, 0, 1, 1},
31+
{0x674e, 0x3b7d, 1, 0, 0},
32+
{0x6807, 0x4014, 1, 0, 1},
33+
{0x6f4d, 0xbe19, 0, 1, 1},
34+
{0x7330, 0xcb62, 0, 1, 0},
35+
}};
36+
37+
LLVM_LIBC_FUNCTION(float16, tanf16, (float16 x)) {
38+
using FPBits = fputil::FPBits<float16>;
39+
FPBits xbits(x);
40+
41+
uint16_t x_u = xbits.uintval();
42+
uint16_t x_abs = x_u & 0x7fff;
43+
bool x_sign = x_u >> 15;
44+
float xf = x;
45+
46+
// Handle exceptional values
47+
if (auto r = TANF16_EXCEPTS.lookup_odd(x_abs, x_sign);
48+
LIBC_UNLIKELY(r.has_value()))
49+
return r.value();
50+
51+
// |x| <= 0x1.d1p-5
52+
if (LIBC_UNLIKELY(x_abs <= 0x2b44)) {
53+
// |x| <= 0x1.398p-11
54+
if (LIBC_UNLIKELY(x_abs <= 0x10e6)) {
55+
// tan(+/-0) = +/-0
56+
if (LIBC_UNLIKELY(x_abs == 0))
57+
return x;
58+
59+
int rounding = fputil::quick_get_round();
60+
61+
// Exhaustive tests show that, when:
62+
// x > 0, and rounding upward or
63+
// x < 0, and rounding downward then,
64+
// tan(x) = x * 2^-11 + x
65+
if ((xbits.is_pos() && rounding == FE_UPWARD) ||
66+
(xbits.is_neg() && rounding == FE_DOWNWARD))
67+
return fputil::cast<float16>(fputil::multiply_add(xf, 0x1.0p-11f, xf));
68+
return x;
69+
}
70+
71+
float xsq = xf * xf;
72+
73+
// Degree-6 minimax odd polynomial of tan(x) generated by Sollya with:
74+
// > P = fpminimax(tan(x)/x, [|0, 2, 4, 6|], [|1, SG...|], [0, pi/32]);
75+
float result = fputil::polyeval(xsq, 0x1p0f, 0x1.555556p-2f, 0x1.110ee4p-3f,
76+
0x1.be80f6p-5f);
77+
78+
return fputil::cast<float16>(xf * result);
79+
}
80+
81+
// tan(+/-inf) = NaN, and tan(NaN) = NaN
82+
if (LIBC_UNLIKELY(x_abs >= 0x7c00)) {
83+
// x = +/-inf
84+
if (x_abs == 0x7c00) {
85+
fputil::set_errno_if_required(EDOM);
86+
fputil::raise_except_if_required(FE_INVALID);
87+
}
88+
89+
return x + FPBits::quiet_nan().get_val();
90+
}
91+
92+
// Range reduction:
93+
// For |x| > pi/32, we perform range reduction as follows:
94+
// Find k and y such that:
95+
// x = (k + y) * pi/32;
96+
// k is an integer, |y| < 0.5
97+
//
98+
// This is done by performing:
99+
// k = round(x * 32/pi)
100+
// y = x * 32/pi - k
101+
//
102+
// Once k and y are computed, we then deduce the answer by the formula:
103+
// tan(x) = sin(x) / cos(x)
104+
// = (sin_y * cos_k + cos_y * sin_k) / (cos_y * cos_k - sin_y * sin_k)
105+
float sin_k, cos_k, sin_y, cosm1_y;
106+
sincosf16_eval(xf, sin_k, cos_k, sin_y, cosm1_y);
107+
108+
// Note that, cosm1_y = cos_y - 1:
109+
using fputil::multiply_add;
110+
return fputil::cast<float16>(
111+
multiply_add(sin_y, cos_k, multiply_add(cosm1_y, sin_k, sin_k)) /
112+
multiply_add(sin_y, -sin_k, multiply_add(cosm1_y, cos_k, cos_k)));
113+
}
114+
115+
} // namespace LIBC_NAMESPACE_DECL

libc/src/math/generic/tanpif16.cpp

+1-1
Original file line numberDiff line numberDiff line change
@@ -79,7 +79,7 @@ LLVM_LIBC_FUNCTION(float16, tanpif16, (float16 x)) {
7979
// k = round(x * 32)
8080
// y = x * 32 - k
8181
//
82-
// Once k and y are computed, we then deduce the answer by tthe formula:
82+
// Once k and y are computed, we then deduce the answer by the formula:
8383
// tan(x) = sin(x) / cos(x)
8484
// = (sin_y * cos_k + cos_y * sin_k) / (cos_y * cos_k - sin_y * sin_k)
8585
float xf = x;

libc/src/math/tanf16.h

+21
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,21 @@
1+
//===-- Implementation header for tanf16 ------------------------*- 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_TANF16_H
10+
#define LLVM_LIBC_SRC_MATH_TANF16_H
11+
12+
#include "src/__support/macros/config.h"
13+
#include "src/__support/macros/properties/types.h"
14+
15+
namespace LIBC_NAMESPACE_DECL {
16+
17+
float16 tanf16(float16 x);
18+
19+
} // namespace LIBC_NAMESPACE_DECL
20+
21+
#endif // LLVM_LIBC_SRC_MATH_TANF16_H

libc/test/src/math/CMakeLists.txt

+11
Original file line numberDiff line numberDiff line change
@@ -190,6 +190,17 @@ add_fp_unittest(
190190
libc.src.__support.FPUtil.fp_bits
191191
)
192192

193+
add_fp_unittest(
194+
tanf16_test
195+
NEED_MPFR
196+
SUITE
197+
libc-math-unittests
198+
SRCS
199+
tanf16_test.cpp
200+
DEPENDS
201+
libc.src.math.tanf16
202+
)
203+
193204
add_fp_unittest(
194205
tanpif16_test
195206
NEED_MPFR

libc/test/src/math/cosf16_test.cpp

+1-1
Original file line numberDiff line numberDiff line change
@@ -17,7 +17,7 @@ namespace mpfr = LIBC_NAMESPACE::testing::mpfr;
1717

1818
// Range: [0, Inf]
1919
static constexpr uint16_t POS_START = 0x0000U;
20-
static constexpr uint16_t POS_STOP = 0x7c00u;
20+
static constexpr uint16_t POS_STOP = 0x7c00U;
2121

2222
// Range: [-Inf, 0]
2323
static constexpr uint16_t NEG_START = 0x8000U;

libc/test/src/math/smoke/CMakeLists.txt

+11
Original file line numberDiff line numberDiff line change
@@ -121,6 +121,17 @@ add_fp_unittest(
121121
libc.src.__support.FPUtil.fp_bits
122122
)
123123

124+
add_fp_unittest(
125+
tanf16_test
126+
SUITE
127+
libc-math-smoke-tests
128+
SRCS
129+
tanf16_test.cpp
130+
DEPENDS
131+
libc.src.errno.errno
132+
libc.src.math.tanf16
133+
)
134+
124135
add_fp_unittest(
125136
tanpif16_test
126137
SUITE
+34
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,34 @@
1+
//===-- Unittests for tanf16 ----------------------------------------------===//
2+
//
3+
//
4+
// Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions.
5+
// See https://llvm.org/LICENSE.txt for license information.
6+
// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception.
7+
//
8+
//===----------------------------------------------------------------------===//
9+
10+
#include "src/errno/libc_errno.h"
11+
#include "src/math/tanf16.h"
12+
#include "test/UnitTest/FPMatcher.h"
13+
#include "test/UnitTest/Test.h"
14+
15+
using LlvmLibcTanf16Test = LIBC_NAMESPACE::testing::FPTest<float16>;
16+
17+
TEST_F(LlvmLibcTanf16Test, SpecialNumbers) {
18+
LIBC_NAMESPACE::libc_errno = 0;
19+
20+
EXPECT_FP_EQ(aNaN, LIBC_NAMESPACE::tanf16(aNaN));
21+
EXPECT_MATH_ERRNO(0);
22+
23+
EXPECT_FP_EQ(zero, LIBC_NAMESPACE::tanf16(zero));
24+
EXPECT_MATH_ERRNO(0);
25+
26+
EXPECT_FP_EQ(neg_zero, LIBC_NAMESPACE::tanf16(neg_zero));
27+
EXPECT_MATH_ERRNO(0);
28+
29+
EXPECT_FP_EQ(aNaN, LIBC_NAMESPACE::tanf16(inf));
30+
EXPECT_MATH_ERRNO(EDOM);
31+
32+
EXPECT_FP_EQ(aNaN, LIBC_NAMESPACE::tanf16(neg_inf));
33+
EXPECT_MATH_ERRNO(EDOM);
34+
}

libc/test/src/math/tanf16_test.cpp

+40
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,40 @@
1+
//===-- Exhaustive test for tanf16 ----------------------------------------===//
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+
#include "src/math/tanf16.h"
10+
#include "test/UnitTest/FPMatcher.h"
11+
#include "test/UnitTest/Test.h"
12+
#include "utils/MPFRWrapper/MPFRUtils.h"
13+
14+
using LlvmLibcTanf16Test = LIBC_NAMESPACE::testing::FPTest<float16>;
15+
16+
namespace mpfr = LIBC_NAMESPACE::testing::mpfr;
17+
18+
// Range: [0, Inf]
19+
static constexpr uint16_t POS_START = 0x0000U;
20+
static constexpr uint16_t POS_STOP = 0x7c00U;
21+
22+
// Range: [-Inf, 0]
23+
static constexpr uint16_t NEG_START = 0x8000U;
24+
static constexpr uint16_t NEG_STOP = 0xfc00U;
25+
26+
TEST_F(LlvmLibcTanf16Test, PositiveRange) {
27+
for (uint16_t v = POS_START; v <= POS_STOP; ++v) {
28+
float16 x = FPBits(v).get_val();
29+
EXPECT_MPFR_MATCH_ALL_ROUNDING(mpfr::Operation::Tan, x,
30+
LIBC_NAMESPACE::tanf16(x), 0.5);
31+
}
32+
}
33+
34+
TEST_F(LlvmLibcTanf16Test, NegativeRange) {
35+
for (uint16_t v = NEG_START; v <= NEG_STOP; ++v) {
36+
float16 x = FPBits(v).get_val();
37+
EXPECT_MPFR_MATCH_ALL_ROUNDING(mpfr::Operation::Tan, x,
38+
LIBC_NAMESPACE::tanf16(x), 0.5);
39+
}
40+
}

0 commit comments

Comments
 (0)