[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
[lmi-commits] [lmi] master 290dc89 10/13: Remove a source of inaccuracy
From: |
Greg Chicares |
Subject: |
[lmi-commits] [lmi] master 290dc89 10/13: Remove a source of inaccuracy [285] |
Date: |
Fri, 9 Apr 2021 18:42:38 -0400 (EDT) |
branch: master
commit 290dc8918d38927beb4316f47dc134a662d7bd04
Author: Gregory W. Chicares <gchicares@sbcglobal.net>
Commit: Gregory W. Chicares <gchicares@sbcglobal.net>
Remove a source of inaccuracy [285]
Formerly, scale_fwd_ was something like 1000 or 0.001, and scale_back_
was its reciprocal, e.g., 1/1000 or 1/.001; but calculating 1/.001, the
reciprocal of the reciprocal of 1000, invites error. Now, they're e.g
1000 and 1/1000, or 1/1000 and 1000. This may not matter as long as the
factors use extended precision, but it's a needless impediment to the
possible future use of double precision.
Expunged an unnecessary cover function for nonstd::power(). It was not
used in enough places to warrant a separate function; worse, it masked
unintended outcomes like calculating 1/.001 above.
Abandoned the idea of favoring const data members. It is more commonly
advised that they be avoided. The idea that they would permit better
optimization seems implausible, and is hardly relevant because round_to
objects are used far more frequently than they're constructed.
---
round_to.hpp | 86 ++++++++++++++-----------------------------------------
round_to_test.cpp | 8 +++++-
2 files changed, 28 insertions(+), 66 deletions(-)
diff --git a/round_to.hpp b/round_to.hpp
index 2339732..a5aed11 100644
--- a/round_to.hpp
+++ b/round_to.hpp
@@ -53,50 +53,6 @@
// this alias-declaration to double.
using max_prec_real = long double;
-namespace detail
-{
-/// Raise 'r' to the integer power 'n'.
-///
-/// Motivation: To raise an integer-valued real to a positive integer
-/// power without any roundoff error as long as the result is exactly
-/// representable. See:
-/// https://lists.nongnu.org/archive/html/lmi/2016-12/msg00049.html
-///
-/// For negative 'n', the most accurate result possible is obtained by
-/// calculating power(r, -n), and returning its reciprocal calculated
-/// with the maximum available precision.
-///
-/// Because this template function is called only by the round_to
-/// constructor, efficiency here is not crucial in the contemplated
-/// typical case where a round_to object is created once and used to
-/// round many numbers, whereas it is crucial to avoid roundoff error.
-/// However, that does not justify gratuitous inefficiency, and the
-/// use of power() here means that the number of multiplications is
-/// O(log n), so this should be as fast as a library function that
-/// has been optimized for accuracy.
-///
-/// Fails to check for overflow or undeflow, but the round_to ctor
-/// does compare 'n' to the minimum and maximum decimal exponents,
-/// which suffices there because its 'r' is always ten.
-
-template<typename RealType>
-RealType int_pow(RealType r, int n)
-{
- if(0 == n)
- {
- return RealType(1.0);
- }
- if(n < 0)
- {
- return max_prec_real(1.0) / nonstd::power(r, -n);
- }
- else
- {
- return nonstd::power(r, n);
- }
-}
-} // namespace detail
-
inline rounding_style& default_rounding_style()
{
static rounding_style default_style = r_to_nearest;
@@ -251,12 +207,6 @@ class round_to
rounding_fn_t rounding_function_ {detail::erroneous_rounding_function};
};
-// Naran used const data members, reasoning that a highly optimizing
-// compiler could then calculate std::pow(10.0, n) at compile time.
-// Not all compilers do this. None available to the author do.
-// The data members were made non-const after profiling showed no
-// penalty on four available compilers around the year 2000.
-
// Division by an exact integer value should have slightly better
// accuracy in some cases. But profiling shows that multiplication by
// the reciprocal stored in scale_back_ makes a realistic application
@@ -268,27 +218,33 @@ template<typename RealType>
round_to<RealType>::round_to(int a_decimals, rounding_style a_style)
:decimals_ {a_decimals}
,style_ {a_style}
- ,scale_fwd_ {detail::int_pow(max_prec_real(10.0), decimals_)}
- ,scale_back_ {max_prec_real(1.0) / scale_fwd_}
,decimals_cents_ {decimals_ - currency::cents_digits}
- ,scale_fwd_cents_ {detail::int_pow(max_prec_real(10.0), decimals_cents_)}
- ,scale_back_cents_ {max_prec_real(1.0) / scale_fwd_cents_}
,rounding_function_ {select_rounding_function(style_)}
{
-/*
-// TODO ?? This might improve accuracy slightly, but would prevent
-// the data members from being const.
- if(0 <= a_decimals)
+ constexpr max_prec_real one( 1.0);
+ constexpr max_prec_real ten(10.0);
+
+ if(0 <= decimals_)
+ {
+ scale_fwd_ = nonstd::power(ten, decimals_);
+ scale_back_ = one / scale_fwd_;
+ }
+ else
+ {
+ scale_back_ = nonstd::power(ten, -decimals_);
+ scale_fwd_ = one / scale_back_;
+ }
+
+ if(0 <= decimals_cents_)
{
- scale_fwd_ = detail::int_pow(max_prec_real(10.0), a_decimals);
- scale_back_ = max_prec_real(1.0) / scale_fwd_;
+ scale_fwd_cents_ = nonstd::power(ten, decimals_cents_);
+ scale_back_cents_ = one / scale_fwd_cents_;
}
else
{
- scale_back_ = detail::int_pow(max_prec_real(10.0), -a_decimals);
- scale_fwd_ = max_prec_real(1.0) / scale_back_;
+ scale_back_cents_ = nonstd::power(ten, -decimals_cents_);
+ scale_fwd_cents_ = one / scale_back_cents_;
}
-*/
// This throws only if all use of the function object is invalid.
// Even if it doesn't throw, there are numbers that it cannot round
@@ -298,8 +254,8 @@ round_to<RealType>::round_to(int a_decimals, rounding_style
a_style)
// std::numeric_limits<RealType>::max_exponent10
// decimals.
if
- (a_decimals < std::numeric_limits<RealType>::min_exponent10
- || std::numeric_limits<RealType>::max_exponent10 <
a_decimals
+ (decimals_ < std::numeric_limits<RealType>::min_exponent10
+ || std::numeric_limits<RealType>::max_exponent10 < decimals_
)
{
throw std::domain_error("Invalid number of decimals.");
diff --git a/round_to_test.cpp b/round_to_test.cpp
index 3b331ca..6d37987 100644
--- a/round_to_test.cpp
+++ b/round_to_test.cpp
@@ -26,6 +26,7 @@
#include "currency.hpp" // currency::cents_digits
#include "fenv_lmi.hpp"
#include "miscellany.hpp" // floating_rep()
+#include "stl_extensions.hpp" // nonstd::power()
#include "test_tools.hpp"
#include <algorithm> // max()
@@ -299,7 +300,12 @@ void round_to_test::test_various_float_types
,long double expected
)
{
- long double factor = detail::int_pow(10.0L, -decimals);
+ int const inverse_decimals = -decimals;
+ long double factor =
+ (0 <= inverse_decimals)
+ ? nonstd::power(10.0L, inverse_decimals)
+ : 1.0L / nonstd::power(10.0L, -inverse_decimals)
+ ;
long double u = unrounded * factor;
long double e = expected * factor;
LMI_TEST((test_one_case(static_cast<float >(u), static_cast<float >(e),
decimals, style)));
- [lmi-commits] [lmi] master updated (8aa09b9 -> a32cee0), Greg Chicares, 2021/04/09
- [lmi-commits] [lmi] master 79741b1 01/13: Assert a precondition more consistently, Greg Chicares, 2021/04/09
- [lmi-commits] [lmi] master 358c5c9 06/13: Modernize, Greg Chicares, 2021/04/09
- [lmi-commits] [lmi] master edb098f 02/13: Strengthen unit tests, Greg Chicares, 2021/04/09
- [lmi-commits] [lmi] master 2a3d961 03/13: Measure cost of a needless transcendental calculation, Greg Chicares, 2021/04/09
- [lmi-commits] [lmi] master 6bc828a 05/13: Realign, Greg Chicares, 2021/04/09
- [lmi-commits] [lmi] master 290dc89 10/13: Remove a source of inaccuracy [285],
Greg Chicares <=
- [lmi-commits] [lmi] master a32cee0 13/13: Hoist a division, Greg Chicares, 2021/04/09
- [lmi-commits] [lmi] master 1e49ee3 04/13: Revise 'round_to' documentation [286], Greg Chicares, 2021/04/09
- [lmi-commits] [lmi] master 1bead9e 08/13: Purge unwanted, commented-out code, Greg Chicares, 2021/04/09
- [lmi-commits] [lmi] master ce9ed5c 07/13: Add an inchoate test of power-of-ten scaling, Greg Chicares, 2021/04/09
- [lmi-commits] [lmi] master 858d037 09/13: Include appropriate headers, and say why they're included, Greg Chicares, 2021/04/09
- [lmi-commits] [lmi] master ca63009 11/13: Expunge an obsolete comment, Greg Chicares, 2021/04/09
- [lmi-commits] [lmi] master c41faec 12/13: Resolve a marked defect [284], Greg Chicares, 2021/04/09