Skip to content

Commit cd04653

Browse files
wldfngrsovermighty
wldfngrs
and
overmighty
authored
[libc][math][c23] Add sinf16 C23 math function (#116674)
Co-authored-by: OverMighty <[email protected]>
1 parent fc11b67 commit cd04653

File tree

13 files changed

+299
-7
lines changed

13 files changed

+299
-7
lines changed

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

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -700,6 +700,7 @@ if(LIBC_TYPES_HAS_FLOAT16)
700700
libc.src.math.scalbnf16
701701
libc.src.math.setpayloadf16
702702
libc.src.math.setpayloadsigf16
703+
libc.src.math.sinf16
703704
libc.src.math.sinhf16
704705
libc.src.math.sinpif16
705706
libc.src.math.sqrtf16

‎libc/docs/math/index.rst‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -336,7 +336,7 @@ Higher Math Functions
336336
+-----------+------------------+-----------------+------------------------+----------------------+------------------------+------------------------+----------------------------+
337337
| rsqrt | | | | | | 7.12.7.9 | F.10.4.9 |
338338
+-----------+------------------+-----------------+------------------------+----------------------+------------------------+------------------------+----------------------------+
339-
| sin | |check| | |check| | | | | 7.12.4.6 | F.10.1.6 |
339+
| sin | |check| | |check| | | |check| | | 7.12.4.6 | F.10.1.6 |
340340
+-----------+------------------+-----------------+------------------------+----------------------+------------------------+------------------------+----------------------------+
341341
| sincos | |check| | |check| | | | | | |
342342
+-----------+------------------+-----------------+------------------------+----------------------+------------------------+------------------------+----------------------------+

‎libc/newhdrgen/yaml/math.yaml‎

Lines changed: 7 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -2339,6 +2339,13 @@ functions:
23392339
return_type: float
23402340
arguments:
23412341
- type: float
2342+
- name: sinf16
2343+
standards:
2344+
- stdc
2345+
return_type: _Float16
2346+
arguments:
2347+
- type: _Float16
2348+
guard: LIBC_TYPES_HAS_FLOAT16
23422349
- name: sinhf
23432350
standards:
23442351
- stdc

‎libc/src/math/CMakeLists.txt‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -484,6 +484,7 @@ add_math_entrypoint_object(sincosf)
484484

485485
add_math_entrypoint_object(sin)
486486
add_math_entrypoint_object(sinf)
487+
add_math_entrypoint_object(sinf16)
487488
add_math_entrypoint_object(sinpif)
488489
add_math_entrypoint_object(sinpif16)
489490

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

Lines changed: 21 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -498,6 +498,27 @@ add_entrypoint_object(
498498
${libc_opt_high_flag}
499499
)
500500

501+
add_entrypoint_object(
502+
sinf16
503+
SRCS
504+
sinf16.cpp
505+
HDRS
506+
../sinf16.h
507+
DEPENDS
508+
.sincosf16_utils
509+
libc.hdr.errno_macros
510+
libc.hdr.fenv_macros
511+
libc.src.__support.FPUtil.cast
512+
libc.src.__support.FPUtil.fenv_impl
513+
libc.src.__support.FPUtil.fp_bits
514+
libc.src.__support.FPUtil.except_value_utils
515+
libc.src.__support.FPUtil.multiply_add
516+
libc.src.__support.macros.optimization
517+
libc.src.__support.macros.properties.types
518+
COMPILE_OPTIONS
519+
-O3
520+
)
521+
501522
add_entrypoint_object(
502523
sincos
503524
SRCS

‎libc/src/math/generic/sincosf16_utils.h‎

Lines changed: 42 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -11,6 +11,7 @@
1111

1212
#include "src/__support/FPUtil/FPBits.h"
1313
#include "src/__support/FPUtil/PolyEval.h"
14+
#include "src/__support/FPUtil/cast.h"
1415
#include "src/__support/FPUtil/nearest_integer.h"
1516
#include "src/__support/common.h"
1617
#include "src/__support/macros/config.h"
@@ -46,10 +47,31 @@ LIBC_INLINE int32_t range_reduction_sincospif16(float x, float &y) {
4647
return static_cast<int32_t>(kf);
4748
}
4849

49-
LIBC_INLINE void sincospif16_eval(float xf, float &sin_k, float &cos_k,
50-
float &sin_y, float &cosm1_y) {
51-
float y;
52-
int32_t k = range_reduction_sincospif16(xf, y);
50+
// Recall, range reduction:
51+
// k = round(x * 32/pi)
52+
// y = x * 32/pi - k
53+
//
54+
// The constant 0x1.45f306dc9c883p3 is 32/pi rounded to double-precision.
55+
// 32/pi is generated by Sollya with the following commands:
56+
// > display = hexadecimal;
57+
// > round(32/pi, D, RN);
58+
//
59+
// The precision choice of 'double' is to minimize rounding errors
60+
// in this initial scaling step, preserving enough bits so errors accumulated
61+
// while computing the subtraction: y = x * 32/pi - round(x * 32/pi)
62+
// are beyond the least-significant bit of single-precision used during
63+
// further intermediate computation.
64+
LIBC_INLINE int32_t range_reduction_sincosf16(float x, float &y) {
65+
double prod = x * 0x1.45f306dc9c883p3;
66+
double kf = fputil::nearest_integer(prod);
67+
y = static_cast<float>(prod - kf);
68+
69+
return static_cast<int32_t>(kf);
70+
}
71+
72+
static LIBC_INLINE void sincosf16_poly_eval(int32_t k, float y, float &sin_k,
73+
float &cos_k, float &sin_y,
74+
float &cosm1_y) {
5375

5476
sin_k = SIN_K_PI_OVER_32[k & 63];
5577
cos_k = SIN_K_PI_OVER_32[(k + 16) & 63];
@@ -72,6 +94,22 @@ LIBC_INLINE void sincospif16_eval(float xf, float &sin_k, float &cos_k,
7294
0x1.a6f7a2p-29f);
7395
}
7496

97+
LIBC_INLINE void sincosf16_eval(float xf, float &sin_k, float &cos_k,
98+
float &sin_y, float &cosm1_y) {
99+
float y;
100+
int32_t k = range_reduction_sincosf16(xf, y);
101+
102+
sincosf16_poly_eval(k, y, sin_k, cos_k, sin_y, cosm1_y);
103+
}
104+
105+
LIBC_INLINE void sincospif16_eval(float xf, float &sin_k, float &cos_k,
106+
float &sin_y, float &cosm1_y) {
107+
float y;
108+
int32_t k = range_reduction_sincospif16(xf, y);
109+
110+
sincosf16_poly_eval(k, y, sin_k, cos_k, sin_y, cosm1_y);
111+
}
112+
75113
} // namespace LIBC_NAMESPACE_DECL
76114

77115
#endif // LLVM_LIBC_SRC_MATH_GENERIC_SINCOSF16_UTILS_H

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

Lines changed: 108 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,108 @@
1+
//===-- Half-precision sin(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/sinf16.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 = 4;
23+
24+
constexpr fputil::ExceptValues<float16, N_EXCEPTS> SINF16_EXCEPTS{{
25+
// (input, RZ output, RU offset, RD offset, RN offset)
26+
{0x2b45, 0x2b43, 1, 0, 1},
27+
{0x585c, 0x3ba3, 1, 0, 1},
28+
{0x5cb0, 0xbbff, 0, 1, 0},
29+
{0x51f5, 0xb80f, 0, 1, 0},
30+
}};
31+
32+
LLVM_LIBC_FUNCTION(float16, sinf16, (float16 x)) {
33+
using FPBits = fputil::FPBits<float16>;
34+
FPBits xbits(x);
35+
36+
uint16_t x_u = xbits.uintval();
37+
uint16_t x_abs = x_u & 0x7fff;
38+
float xf = x;
39+
40+
// Range reduction:
41+
// For |x| > pi/32, we perform range reduction as follows:
42+
// Find k and y such that:
43+
// x = (k + y) * pi/32
44+
// k is an integer, |y| < 0.5
45+
//
46+
// This is done by performing:
47+
// k = round(x * 32/pi)
48+
// y = x * 32/pi - k
49+
//
50+
// Once k and y are computed, we then deduce the answer by the sine of sum
51+
// formula:
52+
// sin(x) = sin((k + y) * pi/32)
53+
// = sin(k * pi/32) * cos(y * pi/32) +
54+
// sin(y * pi/32) * cos(k * pi/32)
55+
56+
// Handle exceptional values
57+
if (LIBC_UNLIKELY(x_abs == 0x585c || x_abs == 0x5cb0 || x_abs == 0x51f5 ||
58+
x_abs == 0x2b45)) {
59+
bool x_sign = x_u >> 15;
60+
if (auto r = SINF16_EXCEPTS.lookup_odd(x_abs, x_sign);
61+
LIBC_UNLIKELY(r.has_value()))
62+
return r.value();
63+
}
64+
65+
int rounding = fputil::quick_get_round();
66+
67+
// Exhaustive tests show that for |x| <= 0x1.f4p-11, 1ULP rounding errors
68+
// occur. To fix this, the following apply:
69+
if (LIBC_UNLIKELY(x_abs <= 0x13d0)) {
70+
// sin(+/-0) = +/-0
71+
if (LIBC_UNLIKELY(x_abs == 0U))
72+
return x;
73+
74+
// When x > 0, and rounding upward, sin(x) == x.
75+
// When x < 0, and rounding downward, sin(x) == x.
76+
if ((rounding == FE_UPWARD && xbits.is_pos()) ||
77+
(rounding == FE_DOWNWARD && xbits.is_neg()))
78+
return x;
79+
80+
// When x < 0, and rounding upward, sin(x) == (x - 1ULP)
81+
if (rounding == FE_UPWARD && xbits.is_neg()) {
82+
x_u--;
83+
return FPBits(x_u).get_val();
84+
}
85+
}
86+
87+
if (xbits.is_inf_or_nan()) {
88+
if (xbits.is_inf()) {
89+
fputil::set_errno_if_required(EDOM);
90+
fputil::raise_except_if_required(FE_INVALID);
91+
}
92+
93+
return x + FPBits::quiet_nan().get_val();
94+
}
95+
96+
float sin_k, cos_k, sin_y, cosm1_y;
97+
sincosf16_eval(xf, sin_k, cos_k, sin_y, cosm1_y);
98+
99+
if (LIBC_UNLIKELY(sin_y == 0 && sin_k == 0))
100+
return FPBits::zero(xbits.sign()).get_val();
101+
102+
// Since, cosm1_y = cos_y - 1, therfore:
103+
// sin(x) = cos_k * sin_y + sin_k + (cosm1_y * sin_k)
104+
return fputil::cast<float16>(fputil::multiply_add(
105+
sin_y, cos_k, fputil::multiply_add(cosm1_y, sin_k, sin_k)));
106+
}
107+
108+
} // namespace LIBC_NAMESPACE_DECL

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

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -21,7 +21,7 @@ namespace LIBC_NAMESPACE_DECL {
2121

2222
constexpr size_t N_EXCEPTS = 21;
2323

24-
constexpr fputil::ExceptValues<float16, N_EXCEPTS> TANF16_EXCEPTS{{
24+
constexpr fputil::ExceptValues<float16, N_EXCEPTS> TANPIF16_EXCEPTS{{
2525
// (input, RZ output, RU offset, RD offset, RN offset)
2626
{0x07f2, 0x0e3d, 1, 0, 0}, {0x086a, 0x0eee, 1, 0, 1},
2727
{0x08db, 0x0fa0, 1, 0, 0}, {0x094c, 0x1029, 1, 0, 0},
@@ -49,7 +49,7 @@ LLVM_LIBC_FUNCTION(float16, tanpif16, (float16 x)) {
4949
return x;
5050

5151
bool x_sign = x_u >> 15;
52-
if (auto r = TANF16_EXCEPTS.lookup_odd(x_abs, x_sign);
52+
if (auto r = TANPIF16_EXCEPTS.lookup_odd(x_abs, x_sign);
5353
LIBC_UNLIKELY(r.has_value()))
5454
return r.value();
5555
}

‎libc/src/math/sinf16.h‎

Lines changed: 21 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,21 @@
1+
//===-- Implementation header for sinf16 ------------------------*- 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_SINF16_H
10+
#define LLVM_LIBC_SRC_MATH_SINF16_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 sinf16(float16 x);
18+
19+
} // namespace LIBC_NAMESPACE_DECL
20+
21+
#endif // LLVM_LIBC_SRC_MATH_SINF16_H

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

Lines changed: 11 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -85,6 +85,17 @@ add_fp_unittest(
8585
libc.src.__support.FPUtil.fp_bits
8686
)
8787

88+
add_fp_unittest(
89+
sinf16_test
90+
NEED_MPFR
91+
SUITE
92+
libc-math-unittests
93+
SRCS
94+
sinf16_test.cpp
95+
DEPENDS
96+
libc.src.math.sinf16
97+
)
98+
8899
add_fp_unittest(
89100
sinpif_test
90101
NEED_MPFR

0 commit comments

Comments
 (0)