Skip to content

Commit ba154a7

Browse files
ckormanyosckormanyos
committed
Preliminary sinh 64, 128
1 parent 49486d5 commit ba154a7

File tree

5 files changed

+206
-49
lines changed

5 files changed

+206
-49
lines changed

‎include/boost/decimal/detail/cmath/impl/cosh_impl.hpp‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -32,7 +32,7 @@ struct cosh_table_imp
3232
public:
3333
static constexpr d32_coeffs_t d32_coeffs =
3434
{{
35-
// Series[Cosh[x], {x, 0, 10}]
35+
// Series[Cosh[x], {x, 0, 12}]
3636
// (1), // * 1
3737
::boost::decimal::decimal32 { 5, -1 }, // * x^2
3838
::boost::decimal::decimal32 { UINT64_C(4166666666666666667), - 19 - 1 }, // * x^6

‎include/boost/decimal/detail/cmath/impl/sinh_impl.hpp‎

Lines changed: 126 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,126 @@
1+
// Copyright 2024 Matt Borland
2+
// Copyright 2024 Christopher Kormanyos
3+
// Distributed under the Boost Software License, Version 1.0.
4+
// https://www.boost.org/LICENSE_1_0.txt
5+
6+
#ifndef BOOST_DECIMAL_DETAIL_CMATH_IMPL_SINH_IMPL_HPP
7+
#define BOOST_DECIMAL_DETAIL_CMATH_IMPL_SINH_IMPL_HPP
8+
9+
#include <boost/decimal/detail/concepts.hpp>
10+
#include <boost/decimal/detail/cmath/impl/taylor_series_result.hpp>
11+
12+
#ifndef BOOST_DECIMAL_BUILD_MODULE
13+
#include <array>
14+
#include <cstddef>
15+
#include <cstdint>
16+
#endif
17+
18+
namespace boost {
19+
namespace decimal {
20+
namespace detail {
21+
22+
namespace sinh_detail {
23+
24+
template <bool b>
25+
struct sinh_table_imp
26+
{
27+
private:
28+
using d32_coeffs_t = std::array<decimal32, 6>;
29+
using d64_coeffs_t = std::array<decimal64, 9>;
30+
using d128_coeffs_t = std::array<decimal128, 17>;
31+
32+
public:
33+
static constexpr d32_coeffs_t d32_coeffs =
34+
{{
35+
// Series[Sinh[x], {x, 0, 13}]
36+
// (1), // * x
37+
::boost::decimal::decimal32 { UINT64_C(1666666666666666667), - 19 - 0 }, // * x^3
38+
::boost::decimal::decimal32 { UINT64_C(8333333333333333333), - 19 - 2 }, // * x^5
39+
::boost::decimal::decimal32 { UINT64_C(1984126984126984127), - 19 - 3 }, // * x^7
40+
::boost::decimal::decimal32 { UINT64_C(2755731922398589065), - 19 - 5 }, // * x^9
41+
::boost::decimal::decimal32 { UINT64_C(2505210838544171878), - 19 - 7 }, // * x^11
42+
::boost::decimal::decimal32 { UINT64_C(1605904383682161460), - 19 - 9 }, // * x^13
43+
}};
44+
45+
static constexpr d64_coeffs_t d64_coeffs =
46+
{{
47+
// Series[Sinh[x], {x, 0, 19}]
48+
// (1), // * x
49+
::boost::decimal::decimal64 { UINT64_C(1666666666666666667), - 19 - 0 }, // * x^3
50+
::boost::decimal::decimal64 { UINT64_C(8333333333333333333), - 19 - 2 }, // * x^5
51+
::boost::decimal::decimal64 { UINT64_C(1984126984126984127), - 19 - 3 }, // * x^7
52+
::boost::decimal::decimal64 { UINT64_C(2755731922398589065), - 19 - 5 }, // * x^9
53+
::boost::decimal::decimal64 { UINT64_C(2505210838544171878), - 19 - 7 }, // * x^11
54+
::boost::decimal::decimal64 { UINT64_C(1605904383682161460), - 19 - 9 }, // * x^13
55+
::boost::decimal::decimal64 { UINT64_C(7647163731819816476), - 19 - 12 }, // * x^15
56+
::boost::decimal::decimal64 { UINT64_C(2811457254345520763), - 19 - 14 }, // * x^17
57+
::boost::decimal::decimal64 { UINT64_C(8220635246624329717), - 19 - 17 } // * x^19
58+
}};
59+
60+
static constexpr d128_coeffs_t d128_coeffs =
61+
{{
62+
// Series[Cosh[x], {x, 0, 34}]
63+
// (1), // * 1
64+
::boost::decimal::decimal128 { 5, -1 }, // * x^2
65+
::boost::decimal::decimal128 { boost::decimal::detail::uint128 { UINT64_C(225875452601146), UINT64_C(13965751134118914724) }, -35 }, // * x^4
66+
::boost::decimal::decimal128 { boost::decimal::detail::uint128 { UINT64_C(75291817533715), UINT64_C(10804165069276155440) }, -36 }, // * x^6
67+
::boost::decimal::decimal128 { boost::decimal::detail::uint128 { UINT64_C(134449674167349), UINT64_C(4799281565792772746) }, -38 }, // * x^8
68+
::boost::decimal::decimal128 { boost::decimal::detail::uint128 { UINT64_C(149388526852610), UINT64_C(5332535073103080820) }, -40 }, // * x^10
69+
::boost::decimal::decimal128 { boost::decimal::detail::uint128 { UINT64_C(113173126403492), UINT64_C(11865690723015477068) }, -42 }, // * x^12
70+
::boost::decimal::decimal128 { boost::decimal::detail::uint128 { UINT64_C(62183036485435), UINT64_C(9560282387433155251) }, -44 }, // * x^14
71+
::boost::decimal::decimal128 { boost::decimal::detail::uint128 { UINT64_C(259095985355981), UINT64_C(6015479145837302244) }, -47 }, // * x^16
72+
::boost::decimal::decimal128 { boost::decimal::detail::uint128 { UINT64_C(84671890639209), UINT64_C(10767230553416093986) }, -49 }, // * x^18
73+
::boost::decimal::decimal128 { boost::decimal::detail::uint128 { UINT64_C(222820764840025), UINT64_C(4062785569898205740) }, -52 }, // * x^20
74+
::boost::decimal::decimal128 { boost::decimal::detail::uint128 { UINT64_C(482296027792262), UINT64_C(7037075391028107068) }, -55 }, // * x^22
75+
::boost::decimal::decimal128 { boost::decimal::detail::uint128 { UINT64_C(87372468802946), UINT64_C(1542176615384940434) }, -57 }, // * x^24
76+
::boost::decimal::decimal128 { boost::decimal::detail::uint128 { UINT64_C(134419182773763), UINT64_C(3791559721646796942) }, -60 }, // * x^26
77+
::boost::decimal::decimal128 { boost::decimal::detail::uint128 { UINT64_C(177803151817147), UINT64_C(1794430560736952558) }, -63 }, // * x^28
78+
::boost::decimal::decimal128 { boost::decimal::detail::uint128 { UINT64_C(204371438870284), UINT64_C(366311534299067156) }, -66 }, // * x^30
79+
::boost::decimal::decimal128 { boost::decimal::detail::uint128 { UINT64_C(206019595635366), UINT64_C(17625897212400736954) }, -69 }, // * x^32
80+
::boost::decimal::decimal128 { boost::decimal::detail::uint128 { UINT64_C(183618177928134), UINT64_C(9987905770721758456) }, -72 }, // * x^34
81+
}};
82+
};
83+
84+
#if !(defined(__cpp_inline_variables) && __cpp_inline_variables >= 201606L) && (!defined(_MSC_VER) || _MSC_VER != 1900)
85+
86+
template <bool b>
87+
constexpr typename sinh_table_imp<b>::d32_coeffs_t sinh_table_imp<b>::d32_coeffs;
88+
89+
template <bool b>
90+
constexpr typename sinh_table_imp<b>::d64_coeffs_t sinh_table_imp<b>::d64_coeffs;
91+
92+
template <bool b>
93+
constexpr typename sinh_table_imp<b>::d128_coeffs_t sinh_table_imp<b>::d128_coeffs;
94+
95+
#endif
96+
97+
} //namespace sinh_detail
98+
99+
using sinh_table = sinh_detail::sinh_table_imp<true>;
100+
101+
template <BOOST_DECIMAL_DECIMAL_FLOATING_TYPE T>
102+
constexpr auto sinh_series_expansion(T z) noexcept;
103+
104+
template <>
105+
constexpr auto sinh_series_expansion<decimal32>(decimal32 z2) noexcept
106+
{
107+
return taylor_series_result(z2, sinh_table::d32_coeffs);
108+
}
109+
110+
template <>
111+
constexpr auto sinh_series_expansion<decimal64>(decimal64 z2) noexcept
112+
{
113+
return taylor_series_result(z2, sinh_table::d64_coeffs);
114+
}
115+
116+
template <>
117+
constexpr auto sinh_series_expansion<decimal128>(decimal128 z2) noexcept
118+
{
119+
return taylor_series_result(z2, sinh_table::d128_coeffs);
120+
}
121+
122+
} //namespace detail
123+
} //namespace decimal
124+
} //namespace boost
125+
126+
#endif //BOOST_DECIMAL_DETAIL_CMATH_IMPL_SINH_IMPL_HPP

‎include/boost/decimal/detail/cmath/sinh.hpp‎

Lines changed: 5 additions & 43 deletions
Original file line numberDiff line numberDiff line change
@@ -1,15 +1,16 @@
1-
// Copyright 2023 Matt Borland
2-
// Copyright 2023 Christopher Kormanyos
1+
// Copyright 2023 - 2024 Matt Borland
2+
// Copyright 2023 - 2024 Christopher Kormanyos
33
// Distributed under the Boost Software License, Version 1.0.
44
// https://www.boost.org/LICENSE_1_0.txt
55

66
#ifndef BOOST_DECIMAL_DETAIL_CMATH_SINH_HPP
77
#define BOOST_DECIMAL_DETAIL_CMATH_SINH_HPP
88

99
#include <boost/decimal/fwd.hpp> // NOLINT(llvm-include-order)
10-
#include <boost/decimal/detail/type_traits.hpp>
10+
#include <boost/decimal/detail/cmath/impl/sinh_impl.hpp>
1111
#include <boost/decimal/detail/concepts.hpp>
1212
#include <boost/decimal/detail/config.hpp>
13+
#include <boost/decimal/detail/type_traits.hpp>
1314
#include <boost/decimal/numbers.hpp>
1415

1516
#ifndef BOOST_DECIMAL_BUILD_MODULE
@@ -57,48 +58,9 @@ constexpr auto sinh(T x) noexcept
5758

5859
if (x < one)
5960
{
60-
using coefficient_array_type = std::array<T, static_cast<std::size_t>(UINT8_C(9))>;
61-
62-
#if (defined(__clang__) && (__clang__ < 6))
63-
# pragma clang diagnostic push
64-
# pragma clang diagnostic ignored "-Wmissing-braces"
65-
#endif
66-
67-
constexpr auto coefficient_table =
68-
coefficient_array_type
69-
{
70-
// Series[Sinh[x], {x, 0, 19}]
71-
// (1), // * x
72-
T { UINT64_C(166666666666666667), -18 - 0 }, // * x^3
73-
T { UINT64_C(833333333333333333), -18 - 2 }, // * x^5
74-
T { UINT64_C(198412698412698413), -18 - 3 }, // * x^7
75-
T { UINT64_C(275573192239858907), -18 - 5 }, // * x^9
76-
T { UINT64_C(250521083854417188), -18 - 7 }, // * x^11
77-
T { UINT64_C(160590438368216146), -18 - 9 }, // * x^13
78-
T { UINT64_C(764716373181981648), -18 - 12 }, // * x^15
79-
T { UINT64_C(281145725434552076), -18 - 14 }, // * x^17
80-
T { UINT64_C(822063524662432972), -18 - 17 }, // * x^19
81-
};
82-
83-
#if (defined(__clang__) && (__clang__ < 6))
84-
# pragma clang diagnostic pop
85-
#endif
86-
87-
auto rit =
88-
coefficient_table.crbegin()
89-
+ static_cast<std::size_t>
90-
(
91-
(sizeof(T) == static_cast<std::size_t>(UINT8_C(4))) ? 4U : 0U
92-
);
93-
94-
result = *rit;
95-
9661
const auto xsq = x * x;
9762

98-
while(rit != coefficient_table.crend())
99-
{
100-
result = fma(result, xsq, *rit++);
101-
}
63+
result = detail::sinh_series_expansion(xsq);
10264

10365
result = x * fma(result, xsq, one);
10466
}

‎test/test_cosh.cpp‎

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -367,6 +367,8 @@ auto main() -> int
367367

368368
BOOST_TEST(result_pos64_is_ok);
369369

370+
BOOST_TEST(result_pos128_is_ok);
371+
370372
result_is_ok = (result_pos_is_ok && result_is_ok);
371373
result_is_ok = (result_neg_is_ok && result_is_ok);
372374

‎test/test_sinh.cpp‎

Lines changed: 72 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -1,5 +1,5 @@
1-
// Copyright 2023 Matt Borland
2-
// Copyright 2023 Christopher Kormanyos
1+
// Copyright 2023 - 2024 Matt Borland
2+
// Copyright 2023 - 2024 Christopher Kormanyos
33
// Distributed under the Boost Software License, Version 1.0.
44
// https://www.boost.org/LICENSE_1_0.txt
55

@@ -54,17 +54,31 @@ namespace local
5454

5555
auto result_is_ok = bool { };
5656

57+
NumericType delta { };
58+
5759
if(b == static_cast<NumericType>(0))
5860
{
59-
result_is_ok = (fabs(a - b) < tol);
61+
delta = fabs(a - b);
62+
63+
result_is_ok = (delta < tol);
6064
}
6165
else
6266
{
63-
const auto delta = fabs(1 - (a / b));
67+
delta = fabs(1 - (a / b));
6468

6569
result_is_ok = (delta < tol);
6670
}
6771

72+
// LCOV_EXCL_START
73+
if (!result_is_ok)
74+
{
75+
std::cerr << std::setprecision(std::numeric_limits<NumericType>::digits10) << "a: " << a
76+
<< "\nb: " << b
77+
<< "\ndelta: " << delta
78+
<< "\ntol: " << tol << std::endl;
79+
}
80+
// LCOV_EXCL_STOP
81+
6882
return result_is_ok;
6983
}
7084

@@ -205,6 +219,55 @@ namespace local
205219
return result_is_ok;
206220
}
207221

222+
auto test_sinh_64(const int tol_factor) -> bool
223+
{
224+
using decimal_type = boost::decimal::decimal64;
225+
226+
using val_ctrl_array_type = std::array<double, 19U>;
227+
228+
const val_ctrl_array_type ctrl_values =
229+
{{
230+
// Table[N[Sinh[n/10 + n/100], 17], {n, 1, 19, 1}]
231+
0.11022196758117152, 0.22177896631245117, 0.33602219751592705,
232+
0.45433539871409734, 0.57815160374345427, 0.70897049995516614,
233+
0.84837659273684347, 0.99805839736781424, 1.1598288906636083,
234+
1.3356474701241768, 1.5276436865595682, 1.7381430376475061,
235+
1.9696951348397458, 2.2251045847805740, 2.5074649592795473,
236+
2.8201962652897691, 3.1670863687357898, 3.5523368739250597,
237+
3.9806140142438027
238+
}};
239+
240+
std::array<decimal_type, std::tuple_size<val_ctrl_array_type>::value> sinh_values { };
241+
242+
int nx { 1 };
243+
244+
bool result_is_ok { true };
245+
246+
const decimal_type my_tol { std::numeric_limits<decimal_type>::epsilon() * static_cast<decimal_type>(tol_factor) };
247+
248+
for(auto i = static_cast<std::size_t>(UINT8_C(0)); i < std::tuple_size<val_ctrl_array_type>::value; ++i)
249+
{
250+
// Table[N[Sinh[n/10 + n/100], 17], {n, 1, 19, 1}]
251+
252+
const decimal_type
253+
x_arg
254+
{
255+
decimal_type { nx, -1 }
256+
+ decimal_type { nx, -2 }
257+
};
258+
259+
sinh_values[i] = sinh(x_arg);
260+
261+
++nx;
262+
263+
const auto result_sinh_is_ok = is_close_fraction(sinh_values[i], decimal_type(ctrl_values[i]), my_tol);
264+
265+
result_is_ok = (result_sinh_is_ok && result_is_ok);
266+
}
267+
268+
return result_is_ok;
269+
}
270+
208271
} // namespace local
209272

210273
auto main() -> int
@@ -222,6 +285,8 @@ auto main() -> int
222285

223286
const auto result_edge_is_ok = local::test_sinh_edge();
224287

288+
const auto result_pos64_is_ok = local::test_sinh_64(64);
289+
225290
BOOST_TEST(result_pos_is_ok);
226291
BOOST_TEST(result_neg_is_ok);
227292

@@ -231,7 +296,7 @@ auto main() -> int
231296
BOOST_TEST(result_pos_wide_is_ok);
232297
BOOST_TEST(result_neg_wide_is_ok);
233298

234-
BOOST_TEST(result_edge_is_ok);
299+
BOOST_TEST(result_pos64_is_ok);
235300

236301
result_is_ok = (result_pos_is_ok && result_is_ok);
237302
result_is_ok = (result_neg_is_ok && result_is_ok);
@@ -244,6 +309,8 @@ auto main() -> int
244309

245310
result_is_ok = (result_edge_is_ok && result_is_ok);
246311

312+
result_is_ok = (result_pos64_is_ok && result_is_ok);
313+
247314
result_is_ok = ((boost::report_errors() == 0) && result_is_ok);
248315

249316
return (result_is_ok ? 0 : -1);

0 commit comments

Comments
 (0)