import new upstream nautilus stable release 14.2.8

[ceph.git] / ceph / src / boost / boost / multiprecision / cpp_bin_float / transcendental.hpp

diff --git a/ceph/src/boost/boost/multiprecision/cpp_bin_float/transcendental.hpp b/ceph/src/boost/boost/multiprecision/cpp_bin_float/transcendental.hpp
index 5c969716a544ed24f1c0b1213b14ce88a72c2b1a..7de9ef9f028c6bea78d1a2d2bad8b648c7638a4e 100644 (file)
--- a/ceph/src/boost/boost/multiprecision/cpp_bin_float/transcendental.hpp
+++ b/ceph/src/boost/boost/multiprecision/cpp_bin_float/transcendental.hpp
@@ -1,15 +1,15 @@
 ///////////////////////////////////////////////////////////////
 //  Copyright 2013 John Maddock. Distributed under the Boost
 //  Software License, Version 1.0. (See accompanying file
-//  LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_
+//  LICENSE_1_0.txt or copy at https://www.boost.org/LICENSE_1_0.txt
 
 #ifndef BOOST_MULTIPRECISION_CPP_BIN_FLOAT_TRANSCENDENTAL_HPP
 #define BOOST_MULTIPRECISION_CPP_BIN_FLOAT_TRANSCENDENTAL_HPP
 
-namespace boost{ namespace multiprecision{ namespace backends{
+namespace boost { namespace multiprecision { namespace backends {
 
 template <unsigned Digits, digit_base_type DigitBase, class Allocator, class Exponent, Exponent MinE, Exponent MaxE>
-void eval_exp_taylor(cpp_bin_float<Digits, DigitBase, Allocator, Exponent, MinE, MaxE> &res, const cpp_bin_float<Digits, DigitBase, Allocator, Exponent, MinE, MaxE> &arg)
+void eval_exp_taylor(cpp_bin_float<Digits, DigitBase, Allocator, Exponent, MinE, MaxE>& res, const cpp_bin_float<Digits, DigitBase, Allocator, Exponent, MinE, MaxE>& arg)
 {
    static const int bits = cpp_bin_float<Digits, DigitBase, Allocator, Exponent, MinE, MaxE>::bit_count;
    //
@@ -20,19 +20,19 @@ void eval_exp_taylor(cpp_bin_float<Digits, DigitBase, Allocator, Exponent, MinE,
    denom = limb_type(1);
    eval_add(res, num);
 
-   for(unsigned k = 2; ; ++k)
+   for (unsigned k = 2;; ++k)
    {
       eval_multiply(denom, k);
       eval_multiply(num, arg);
       eval_divide(t, num, denom);
       eval_add(res, t);
-      if(eval_is_zero(t) || (res.exponent() - bits > t.exponent()))
+      if (eval_is_zero(t) || (res.exponent() - bits > t.exponent()))
          break;
    }
 }
 
 template <unsigned Digits, digit_base_type DigitBase, class Allocator, class Exponent, Exponent MinE, Exponent MaxE>
-void eval_exp(cpp_bin_float<Digits, DigitBase, Allocator, Exponent, MinE, MaxE> &res, const cpp_bin_float<Digits, DigitBase, Allocator, Exponent, MinE, MaxE> &arg)
+void eval_exp(cpp_bin_float<Digits, DigitBase, Allocator, Exponent, MinE, MaxE>& res, const cpp_bin_float<Digits, DigitBase, Allocator, Exponent, MinE, MaxE>& arg)
 {
    //
    // This is based on MPFR's method, let:
@@ -60,35 +60,36 @@ void eval_exp(cpp_bin_float<Digits, DigitBase, Allocator, Exponent, MinE, MaxE>
    // Then add the final 1 at the end, given that e0 is small, this effectively wipes
    // out the error in the last step.
    //
-   using default_ops::eval_multiply;
-   using default_ops::eval_subtract;
    using default_ops::eval_add;
    using default_ops::eval_convert_to;
+   using default_ops::eval_increment;
+   using default_ops::eval_multiply;
+   using default_ops::eval_subtract;
 
-   int type = eval_fpclassify(arg);
+   int  type  = eval_fpclassify(arg);
    bool isneg = eval_get_sign(arg) < 0;
-   if(type == (int)FP_NAN)
+   if (type == (int)FP_NAN)
    {
-      res = arg;
+      res   = arg;
       errno = EDOM;
       return;
    }
-   else if(type == (int)FP_INFINITE)
+   else if (type == (int)FP_INFINITE)
    {
       res = arg;
-      if(isneg)
+      if (isneg)
          res = limb_type(0u);
-      else 
+      else
          res = arg;
       return;
    }
-   else if(type == (int)FP_ZERO)
+   else if (type == (int)FP_ZERO)
    {
       res = limb_type(1);
       return;
    }
    cpp_bin_float<Digits, DigitBase, Allocator, Exponent, MinE, MaxE> t, n;
-   if(isneg)
+   if (isneg)
    {
       t = arg;
       t.negate();
@@ -104,7 +105,14 @@ void eval_exp(cpp_bin_float<Digits, DigitBase, Allocator, Exponent, MinE, MaxE>
    eval_multiply(t, n, default_ops::get_constant_ln2<cpp_bin_float<Digits, DigitBase, Allocator, Exponent, MinE, MaxE> >());
    eval_subtract(t, arg);
    t.negate();
-   if(eval_get_sign(t) < 0)
+   if (t.compare(default_ops::get_constant_ln2<cpp_bin_float<Digits, DigitBase, Allocator, Exponent, MinE, MaxE> >()) > 0)
+   {
+      // There are some rare cases where the multiply rounds down leaving a remainder > ln2
+      // See https://github.com/boostorg/multiprecision/issues/120
+      eval_increment(n);
+      t = limb_type(0);
+   }
+   if (eval_get_sign(t) < 0)
    {
       // There are some very rare cases where arg/ln2 is an integer, and the subsequent multiply
       // rounds up, in that situation t ends up negative at this point which breaks our invariants below:
@@ -131,7 +139,7 @@ void eval_exp(cpp_bin_float<Digits, DigitBase, Allocator, Exponent, MinE, MaxE>
    //
    // Square 1 + res k times:
    //
-   for(Exponent s = 0; s < k; ++s)
+   for (Exponent s = 0; s < k; ++s)
    {
       t.swap(res);
       eval_multiply(res, t, t);
@@ -142,7 +150,6 @@ void eval_exp(cpp_bin_float<Digits, DigitBase, Allocator, Exponent, MinE, MaxE>
    eval_ldexp(res, res, nn);
 }
 
-}}} // namespaces
+}}} // namespace boost::multiprecision::backends
 
 #endif
-