Updated boost to v1.55.0

[rsem.git] / boost / math / tools / toms748_solve.hpp

diff --git a/boost/math/tools/toms748_solve.hpp b/boost/math/tools/toms748_solve.hpp
index e7ada17435fbc8f00f53d6c260142f1ea1496ab4..1310ff29a69f3cd3b1e3688ea759ab34f08ce6e9 100644 (file)
--- a/boost/math/tools/toms748_solve.hpp
+++ b/boost/math/tools/toms748_solve.hpp
@@ -26,12 +26,12 @@ public:
    eps_tolerance(unsigned bits)
    {
       BOOST_MATH_STD_USING
-      eps = (std::max)(T(ldexp(1.0F, 1-bits)), T(2 * tools::epsilon<T>()));
+      eps = (std::max)(T(ldexp(1.0F, 1-bits)), T(4 * tools::epsilon<T>()));
    }
    bool operator()(const T& a, const T& b)
    {
       BOOST_MATH_STD_USING
-      return (fabs(a - b) / (std::min)(fabs(a), fabs(b))) <= eps;
+      return fabs(a - b) <= (eps * (std::min)(fabs(a), fabs(b)));
    }
 private:
    T eps;
@@ -197,7 +197,7 @@ T quadratic_interpolate(const T& a, const T& b, T const& d,
    T A = safe_div(T(fd - fb), T(d - b), tools::max_value<T>());
    A = safe_div(T(A - B), T(d - a), T(0));
 
-   if(a == 0)
+   if(A == 0)
    {
       // failure to determine coefficients, try a secant step:
       return secant_interpolate(a, b, fa, fb);
@@ -422,9 +422,10 @@ std::pair<T, T> toms748_solve(F f, const T& ax, const T& bx, const T& fax, const
       e = d;
       fe = fd;
       detail::bracket(f, a, b, c, fa, fb, d, fd);
+      BOOST_MATH_INSTRUMENT_CODE(" a = " << a << " b = " << b);
+      BOOST_MATH_INSTRUMENT_CODE(" tol = " << T((fabs(a) - fabs(b)) / fabs(a)));
       if((0 == --count) || (fa == 0) || tol(a, b))
          break;
-      BOOST_MATH_INSTRUMENT_CODE(" a = " << a << " b = " << b);
       //
       // And finally... check to see if an additional bisection step is 
       // to be taken, we do this if we're not converging fast enough: