Skip to content

Commit d9c771b

Browse files
ckormanyosckormanyos
committed
Preliminary finish and optimize log10
1 parent 37a578c commit d9c771b

File tree

3 files changed

+129
-26
lines changed

3 files changed

+129
-26
lines changed

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

Lines changed: 30 additions & 23 deletions
Original file line numberDiff line numberDiff line change
@@ -50,21 +50,28 @@ constexpr auto log10_impl(T x) noexcept
5050
{
5151
int exp10val { };
5252

53-
auto gn { frexp10(x, &exp10val) };
53+
using representation_type = std::conditional_t<std::is_same<T, decimal32>::value, std::uint32_t,
54+
std::conditional_t<std::is_same<T, decimal64>::value, std::uint64_t, detail::uint128>>;
5455

55-
const auto zeros_removal_result { detail::remove_trailing_zeros(gn) };
56+
const representation_type gn { frexp10(x, &exp10val) };
5657

57-
const bool is_pure { (zeros_removal_result.trimmed_number == one) };
58+
const remove_trailing_zeros_return<representation_type>
59+
zeros_removal
60+
{
61+
remove_trailing_zeros(gn)
62+
};
63+
64+
const bool is_pure
65+
{
66+
(zeros_removal.trimmed_number == static_cast<representation_type>(UINT8_C(1)))
67+
};
5868

5969
if(is_pure)
6070
{
6171
// Here, a pure power-of-10 argument gets a pure integral result.
62-
result =
63-
T
64-
{
65-
exp10val
66-
+ static_cast<int>(zeros_removal_result.number_of_removed_zeros)
67-
};
72+
const int p10 { exp10val + static_cast<int>(zeros_removal.number_of_removed_zeros) };
73+
74+
result = T { p10 };
6875
}
6976
else
7077
{
@@ -86,23 +93,25 @@ constexpr auto log10_impl(T x) noexcept
8693
}
8794
else if(x > one)
8895
{
89-
// The algorithm for base-10 logarithm is based on Chapter 5, pages 35-36
90-
// of Cody and Waite, Software Manual for the Elementary Functions,
91-
// Prentice Hall, 1980.
96+
// The algorithm for base-10 logarithm is based on Chapter 5,
97+
// pages 35-36 of Cody and Waite, "Software Manual for the
98+
// Elementary Functions", Prentice Hall, 1980.
9299

93-
T g { T { gn } / pow10(std::numeric_limits<T>::digits10) };
100+
// In this implementation, however, we use 2s (as for natural
101+
// logarithm in Cody and Waite) even though we are computing
102+
// the base-10 logarithm.
94103

95-
exp10val += std::numeric_limits<T>::digits10;
104+
T g { gn, -std::numeric_limits<T>::digits10 };
96105

97-
constexpr T inv_sqrt10 { UINT64_C(3162277660168379332), -19 };
106+
exp10val += std::numeric_limits<T>::digits10;
98107

99-
const bool reduce_sqrt10 { g < inv_sqrt10 };
108+
int reduce_sqrt2 { };
100109

101-
if (reduce_sqrt10)
110+
while (g < numbers::inv_sqrt2_v<T>)
102111
{
103-
constexpr T sqrt10 { UINT64_C(3162277660168379332), -18 };
112+
g += g;
104113

105-
g *= sqrt10;
114+
++reduce_sqrt2;
106115
}
107116

108117
const T s { (g - one) / (g + one) };
@@ -113,11 +122,9 @@ constexpr auto log10_impl(T x) noexcept
113122

114123
result /= numbers::ln10_v<T>;
115124

116-
if(reduce_sqrt10)
125+
for(int i = 0; i < reduce_sqrt2; ++i)
117126
{
118-
constexpr T half { 5, -1 };
119-
120-
result -= half;
127+
result -= numbers::log10_2_v<T>;
121128
}
122129

123130
result += static_cast<T>(exp10val);

‎include/boost/decimal/numbers.hpp‎

Lines changed: 9 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -38,6 +38,13 @@ BOOST_DECIMAL_EXPORT template <>
3838
BOOST_DECIMAL_CONSTEXPR_VARIABLE_SPECIALIZATION decimal128 log10e_v<decimal128> = decimal128{detail::uint128{UINT64_C(235431510388986),
3939
UINT64_C(2047877485384264674)}, -34};
4040

41+
BOOST_DECIMAL_EXPORT template <BOOST_DECIMAL_DECIMAL_FLOATING_TYPE Dec, std::enable_if_t<detail::is_decimal_floating_point_v<Dec>, bool> = true>
42+
BOOST_DECIMAL_CONSTEXPR_VARIABLE Dec log10_2_v = Dec{UINT64_C(3010299956639811952), -19};
43+
44+
BOOST_DECIMAL_EXPORT template <>
45+
BOOST_DECIMAL_CONSTEXPR_VARIABLE_SPECIALIZATION decimal128 log10_2_v<decimal128> = decimal128{detail::uint128{UINT64_C(163188687641095),
46+
UINT64_C(3612628795761985410)}, -34};
47+
4148
BOOST_DECIMAL_EXPORT template <BOOST_DECIMAL_DECIMAL_FLOATING_TYPE Dec, std::enable_if_t<detail::is_decimal_floating_point_v<Dec>, bool> = true>
4249
BOOST_DECIMAL_CONSTEXPR_VARIABLE Dec pi_v = Dec{UINT64_C(3141592653589793238), -18};
4350

@@ -127,8 +134,9 @@ BOOST_DECIMAL_CONSTEXPR_VARIABLE_SPECIALIZATION decimal128 phi_v<decimal128> = d
127134
UINT64_C(2061523135646567614)}, -33};
128135

129136
BOOST_DECIMAL_EXPORT BOOST_DECIMAL_CONSTEXPR_VARIABLE auto e {e_v<decimal64>};
130-
BOOST_DECIMAL_EXPORT BOOST_DECIMAL_CONSTEXPR_VARIABLE auto log2e {log2e_v<decimal64>};
137+
BOOST_DECIMAL_EXPORT BOOST_DECIMAL_CONSTEXPR_VARIABLE auto log10_2 {log10_2_v<decimal64>};
131138
BOOST_DECIMAL_EXPORT BOOST_DECIMAL_CONSTEXPR_VARIABLE auto log10e {log10e_v<decimal64>};
139+
BOOST_DECIMAL_EXPORT BOOST_DECIMAL_CONSTEXPR_VARIABLE auto log2e {log2e_v<decimal64>};
132140
BOOST_DECIMAL_EXPORT BOOST_DECIMAL_CONSTEXPR_VARIABLE auto pi {pi_v<decimal64>};
133141
BOOST_DECIMAL_EXPORT BOOST_DECIMAL_CONSTEXPR_VARIABLE auto inv_pi {inv_pi_v<decimal64>};
134142
BOOST_DECIMAL_EXPORT BOOST_DECIMAL_CONSTEXPR_VARIABLE auto inv_sqrtpi {inv_sqrtpi_v<decimal64>};

‎test/test_log10.cpp‎

Lines changed: 90 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -98,7 +98,7 @@ namespace local
9898
auto dis_r =
9999
std::uniform_real_distribution<float_type>
100100
{
101-
static_cast<float_type>(1.2L),
101+
static_cast<float_type>(1.4L),
102102
static_cast<float_type>(8.9L)
103103
};
104104

@@ -294,6 +294,75 @@ namespace local
294294
return result_is_ok;
295295
}
296296

297+
auto test_log10_128(const int tol_factor) -> bool
298+
{
299+
using decimal_type = boost::decimal::decimal128;
300+
301+
using str_ctrl_array_type = std::array<const char*, 28U>;
302+
303+
const str_ctrl_array_type ctrl_strings =
304+
{{
305+
// Table[N[Log[10, 456 10^n], 36], {n, -3, 24, 1}]
306+
"-0.341035157335565015527421936814476293",
307+
"0.658964842664434984472578063185523707",
308+
"1.65896484266443498447257806318552371",
309+
"2.65896484266443498447257806318552371",
310+
"3.65896484266443498447257806318552371",
311+
"4.65896484266443498447257806318552371",
312+
"5.65896484266443498447257806318552371",
313+
"6.65896484266443498447257806318552371",
314+
"7.65896484266443498447257806318552371",
315+
"8.65896484266443498447257806318552371",
316+
"9.65896484266443498447257806318552371",
317+
"10.6589648426644349844725780631855237",
318+
"11.6589648426644349844725780631855237",
319+
"12.6589648426644349844725780631855237",
320+
"13.6589648426644349844725780631855237",
321+
"14.6589648426644349844725780631855237",
322+
"15.6589648426644349844725780631855237",
323+
"16.6589648426644349844725780631855237",
324+
"17.6589648426644349844725780631855237",
325+
"18.6589648426644349844725780631855237",
326+
"19.6589648426644349844725780631855237",
327+
"20.6589648426644349844725780631855237",
328+
"21.6589648426644349844725780631855237",
329+
"22.6589648426644349844725780631855237",
330+
"23.6589648426644349844725780631855237",
331+
"24.6589648426644349844725780631855237",
332+
"25.6589648426644349844725780631855237",
333+
"26.6589648426644349844725780631855237",
334+
}};
335+
336+
std::array<decimal_type, std::tuple_size<str_ctrl_array_type>::value> log_values { };
337+
std::array<decimal_type, std::tuple_size<str_ctrl_array_type>::value> ctrl_values { };
338+
339+
int nx { -3 };
340+
341+
bool result_is_ok { true };
342+
343+
const decimal_type my_tol { std::numeric_limits<decimal_type>::epsilon() * static_cast<decimal_type>(tol_factor) };
344+
345+
for(auto i = static_cast<std::size_t>(UINT8_C(0)); i < std::tuple_size<str_ctrl_array_type>::value; ++i)
346+
{
347+
const decimal_type x_arg { 456, nx };
348+
349+
++nx;
350+
351+
log_values[i] = log10(x_arg);
352+
353+
static_cast<void>
354+
(
355+
from_chars(ctrl_strings[i], ctrl_strings[i] + std::strlen(ctrl_strings[i]), ctrl_values[i])
356+
);
357+
358+
const auto result_log_is_ok = is_close_fraction(log_values[i], ctrl_values[i], my_tol);
359+
360+
result_is_ok = (result_log_is_ok && result_is_ok);
361+
}
362+
363+
return result_is_ok;
364+
}
365+
297366
} // namespace local
298367

299368
auto main() -> int
@@ -304,7 +373,18 @@ auto main() -> int
304373
using decimal_type = boost::decimal::decimal32;
305374
using float_type = float;
306375

307-
const auto test_log10_is_ok = local::test_log10<decimal_type, float_type>(128);
376+
const auto test_log10_is_ok = local::test_log10<decimal_type, float_type>(64);
377+
378+
BOOST_TEST(test_log10_is_ok);
379+
380+
result_is_ok = (test_log10_is_ok && result_is_ok);
381+
}
382+
383+
{
384+
using decimal_type = boost::decimal::decimal64;
385+
using float_type = double;
386+
387+
const auto test_log10_is_ok = local::test_log10<decimal_type, float_type>(256);
308388

309389
BOOST_TEST(test_log10_is_ok);
310390

@@ -333,6 +413,14 @@ auto main() -> int
333413
result_is_ok = (test_log10_edge_is_ok && result_is_ok);
334414
}
335415

416+
{
417+
const auto result_pos128_is_ok = local::test_log10_128(1'400'000);
418+
419+
BOOST_TEST(result_pos128_is_ok);
420+
421+
result_is_ok = (result_pos128_is_ok && result_is_ok);
422+
}
423+
336424
result_is_ok = ((boost::report_errors() == 0) && result_is_ok);
337425

338426
return (result_is_ok ? 0 : -1);

0 commit comments

Comments
 (0)