Imported Upstream version 1.72.0
[platform/upstream/boost.git] / boost / multiprecision / detail / functions / trig.hpp
index c001f95..38649fb 100644 (file)
 // in ACM TOMS, {VOL 37, ISSUE 4, (February 2011)} (C) ACM, 2011. http://doi.acm.org/10.1145/1916461.1916469
 //
 // This file has no include guards or namespaces - it's expanded inline inside default_ops.hpp
-// 
+//
 
 #ifdef BOOST_MSVC
 #pragma warning(push)
-#pragma warning(disable:6326)  // comparison of two constants
+#pragma warning(disable : 6326) // comparison of two constants
 #endif
 
 template <class T>
 void hyp0F1(T& result, const T& b, const T& x)
 {
-   typedef typename boost::multiprecision::detail::canonical<boost::int32_t, T>::type si_type;
+   typedef typename boost::multiprecision::detail::canonical<boost::int32_t, T>::type  si_type;
    typedef typename boost::multiprecision::detail::canonical<boost::uint32_t, T>::type ui_type;
 
    // Compute the series representation of Hypergeometric0F1 taken from
@@ -28,8 +28,8 @@ void hyp0F1(T& result, const T& b, const T& x)
    // There are no checks on input range or parameter boundaries.
 
    T x_pow_n_div_n_fact(x);
-   T pochham_b         (b);
-   T bp                (b);
+   T pochham_b(b);
+   T bp(b);
 
    eval_divide(result, x_pow_n_div_n_fact, pochham_b);
    eval_add(result, ui_type(1));
@@ -40,15 +40,16 @@ void hyp0F1(T& result, const T& b, const T& x)
    tol = ui_type(1);
    eval_ldexp(tol, tol, 1 - boost::multiprecision::detail::digits2<number<T, et_on> >::value());
    eval_multiply(tol, result);
-   if(eval_get_sign(tol) < 0)
+   if (eval_get_sign(tol) < 0)
       tol.negate();
    T term;
 
-   const int series_limit = 
-      boost::multiprecision::detail::digits2<number<T, et_on> >::value() < 100
-      ? 100 : boost::multiprecision::detail::digits2<number<T, et_on> >::value();
+   const int series_limit =
+       boost::multiprecision::detail::digits2<number<T, et_on> >::value() < 100
+           ? 100
+           : boost::multiprecision::detail::digits2<number<T, et_on> >::value();
    // Series expansion of hyperg_0f1(; b; x).
-   for(n = 2; n < series_limit; ++n)
+   for (n = 2; n < series_limit; ++n)
    {
       eval_multiply(x_pow_n_div_n_fact, x);
       eval_divide(x_pow_n_div_n_fact, n);
@@ -59,22 +60,21 @@ void hyp0F1(T& result, const T& b, const T& x)
       eval_add(result, term);
 
       bool neg_term = eval_get_sign(term) < 0;
-      if(neg_term)
+      if (neg_term)
          term.negate();
-      if(term.compare(tol) <= 0)
+      if (term.compare(tol) <= 0)
          break;
    }
 
-   if(n >= series_limit)
+   if (n >= series_limit)
       BOOST_THROW_EXCEPTION(std::runtime_error("H0F1 Failed to Converge"));
 }
 
-
 template <class T>
 void eval_sin(T& result, const T& x)
 {
    BOOST_STATIC_ASSERT_MSG(number_category<T>::value == number_kind_floating_point, "The sin function is only valid for floating point types.");
-   if(&result == &x)
+   if (&result == &x)
    {
       T temp;
       eval_sin(temp, x);
@@ -82,18 +82,18 @@ void eval_sin(T& result, const T& x)
       return;
    }
 
-   typedef typename boost::multiprecision::detail::canonical<boost::int32_t, T>::type si_type;
+   typedef typename boost::multiprecision::detail::canonical<boost::int32_t, T>::type  si_type;
    typedef typename boost::multiprecision::detail::canonical<boost::uint32_t, T>::type ui_type;
-   typedef typename mpl::front<typename T::float_types>::type fp_type;
+   typedef typename mpl::front<typename T::float_types>::type                          fp_type;
 
-   switch(eval_fpclassify(x))
+   switch (eval_fpclassify(x))
    {
    case FP_INFINITE:
    case FP_NAN:
-      if(std::numeric_limits<number<T, et_on> >::has_quiet_NaN)
+      if (std::numeric_limits<number<T, et_on> >::has_quiet_NaN)
       {
          result = std::numeric_limits<number<T, et_on> >::quiet_NaN().backend();
-         errno = EDOM;
+         errno  = EDOM;
       }
       else
          BOOST_THROW_EXCEPTION(std::domain_error("Result is undefined or complex and there is no NaN for this number type."));
@@ -101,7 +101,7 @@ void eval_sin(T& result, const T& x)
    case FP_ZERO:
       result = x;
       return;
-   default: ;
+   default:;
    }
 
    // Local copy of the argument
@@ -112,7 +112,7 @@ void eval_sin(T& result, const T& x)
    // The argument xx will be reduced to 0 <= xx <= pi/2.
    bool b_negate_sin = false;
 
-   if(eval_get_sign(x) < 0)
+   if (eval_get_sign(x) < 0)
    {
       xx.negate();
       b_negate_sin = !b_negate_sin;
@@ -120,7 +120,7 @@ void eval_sin(T& result, const T& x)
 
    T n_pi, t;
    // Remove even multiples of pi.
-   if(xx.compare(get_constant_pi<T>()) > 0)
+   if (xx.compare(get_constant_pi<T>()) > 0)
    {
       eval_divide(n_pi, xx, get_constant_pi<T>());
       eval_trunc(n_pi, n_pi);
@@ -140,7 +140,7 @@ void eval_sin(T& result, const T& x)
       BOOST_MATH_INSTRUMENT_CODE(n_pi.str(0, std::ios_base::scientific));
 
       // Adjust signs if the multiple of pi is not even.
-      if(!b_n_pi_is_even)
+      if (!b_n_pi_is_even)
       {
          b_negate_sin = !b_negate_sin;
       }
@@ -148,7 +148,7 @@ void eval_sin(T& result, const T& x)
 
    // Reduce the argument to 0 <= xx <= pi/2.
    eval_ldexp(t, get_constant_pi<T>(), -1);
-   if(xx.compare(t) > 0)
+   if (xx.compare(t) > 0)
    {
       eval_subtract(xx, get_constant_pi<T>(), xx);
       BOOST_MATH_INSTRUMENT_CODE(xx.str(0, std::ios_base::scientific));
@@ -159,18 +159,19 @@ void eval_sin(T& result, const T& x)
    const bool b_pi_half = eval_get_sign(t) == 0;
 
    // Check if the reduced argument is very close to 0 or pi/2.
-   const bool    b_near_zero    = xx.compare(fp_type(1e-1)) < 0;
-   const bool    b_near_pi_half = t.compare(fp_type(1e-1)) < 0;;
+   const bool b_near_zero    = xx.compare(fp_type(1e-1)) < 0;
+   const bool b_near_pi_half = t.compare(fp_type(1e-1)) < 0;
+   ;
 
-   if(b_zero)
+   if (b_zero)
    {
       result = ui_type(0);
    }
-   else if(b_pi_half)
+   else if (b_pi_half)
    {
       result = ui_type(1);
    }
-   else if(b_near_zero)
+   else if (b_near_zero)
    {
       eval_multiply(t, xx, xx);
       eval_divide(t, si_type(-4));
@@ -180,7 +181,7 @@ void eval_sin(T& result, const T& x)
       BOOST_MATH_INSTRUMENT_CODE(result.str(0, std::ios_base::scientific));
       eval_multiply(result, xx);
    }
-   else if(b_near_pi_half)
+   else if (b_near_pi_half)
    {
       eval_multiply(t, t);
       eval_divide(t, si_type(-4));
@@ -196,7 +197,7 @@ void eval_sin(T& result, const T& x)
       // divide by three identity a certain number of times.
       // Here we use division by 3^9 --> (19683 = 3^9).
 
-      static const si_type n_scale = 9;
+      static const si_type n_scale           = 9;
       static const si_type n_three_pow_scale = static_cast<si_type>(19683L);
 
       eval_divide(xx, n_three_pow_scale);
@@ -211,7 +212,7 @@ void eval_sin(T& result, const T& x)
       eval_multiply(result, xx);
 
       // Convert back using multiple angle identity.
-      for(boost::int32_t k = static_cast<boost::int32_t>(0); k < n_scale; k++)
+      for (boost::int32_t k = static_cast<boost::int32_t>(0); k < n_scale; k++)
       {
          // Rescale the cosine value using the multiple angle identity.
          eval_multiply(t2, result, ui_type(3));
@@ -222,7 +223,7 @@ void eval_sin(T& result, const T& x)
       }
    }
 
-   if(b_negate_sin)
+   if (b_negate_sin)
       result.negate();
 }
 
@@ -230,7 +231,7 @@ template <class T>
 void eval_cos(T& result, const T& x)
 {
    BOOST_STATIC_ASSERT_MSG(number_category<T>::value == number_kind_floating_point, "The cos function is only valid for floating point types.");
-   if(&result == &x)
+   if (&result == &x)
    {
       T temp;
       eval_cos(temp, x);
@@ -238,18 +239,18 @@ void eval_cos(T& result, const T& x)
       return;
    }
 
-   typedef typename boost::multiprecision::detail::canonical<boost::int32_t, T>::type si_type;
+   typedef typename boost::multiprecision::detail::canonical<boost::int32_t, T>::type  si_type;
    typedef typename boost::multiprecision::detail::canonical<boost::uint32_t, T>::type ui_type;
-   typedef typename mpl::front<typename T::float_types>::type fp_type;
+   typedef typename mpl::front<typename T::float_types>::type                          fp_type;
 
-   switch(eval_fpclassify(x))
+   switch (eval_fpclassify(x))
    {
    case FP_INFINITE:
    case FP_NAN:
-      if(std::numeric_limits<number<T, et_on> >::has_quiet_NaN)
+      if (std::numeric_limits<number<T, et_on> >::has_quiet_NaN)
       {
          result = std::numeric_limits<number<T, et_on> >::quiet_NaN().backend();
-         errno = EDOM;
+         errno  = EDOM;
       }
       else
          BOOST_THROW_EXCEPTION(std::domain_error("Result is undefined or complex and there is no NaN for this number type."));
@@ -257,7 +258,7 @@ void eval_cos(T& result, const T& x)
    case FP_ZERO:
       result = ui_type(1);
       return;
-   default: ;
+   default:;
    }
 
    // Local copy of the argument
@@ -268,14 +269,14 @@ void eval_cos(T& result, const T& x)
    // The argument xx will be reduced to 0 <= xx <= pi/2.
    bool b_negate_cos = false;
 
-   if(eval_get_sign(x) < 0)
+   if (eval_get_sign(x) < 0)
    {
       xx.negate();
    }
 
    T n_pi, t;
    // Remove even multiples of pi.
-   if(xx.compare(get_constant_pi<T>()) > 0)
+   if (xx.compare(get_constant_pi<T>()) > 0)
    {
       eval_divide(t, xx, get_constant_pi<T>());
       eval_trunc(n_pi, t);
@@ -303,7 +304,7 @@ void eval_cos(T& result, const T& x)
       eval_fmod(t, n_pi, t);
       const bool b_n_pi_is_even = eval_get_sign(t) == 0;
 
-      if(!b_n_pi_is_even)
+      if (!b_n_pi_is_even)
       {
          b_negate_cos = !b_negate_cos;
       }
@@ -312,7 +313,7 @@ void eval_cos(T& result, const T& x)
    // Reduce the argument to 0 <= xx <= pi/2.
    eval_ldexp(t, get_constant_pi<T>(), -1);
    int com = xx.compare(t);
-   if(com > 0)
+   if (com > 0)
    {
       eval_subtract(xx, get_constant_pi<T>(), xx);
       b_negate_cos = !b_negate_cos;
@@ -323,17 +324,17 @@ void eval_cos(T& result, const T& x)
    const bool b_pi_half = com == 0;
 
    // Check if the reduced argument is very close to 0.
-   const bool    b_near_zero    = xx.compare(fp_type(1e-1)) < 0;
+   const bool b_near_zero = xx.compare(fp_type(1e-1)) < 0;
 
-   if(b_zero)
+   if (b_zero)
    {
       result = si_type(1);
    }
-   else if(b_pi_half)
+   else if (b_pi_half)
    {
       result = si_type(0);
    }
-   else if(b_near_zero)
+   else if (b_near_zero)
    {
       eval_multiply(t, xx, xx);
       eval_divide(t, si_type(-4));
@@ -346,7 +347,7 @@ void eval_cos(T& result, const T& x)
       eval_subtract(t, xx);
       eval_sin(result, t);
    }
-   if(b_negate_cos)
+   if (b_negate_cos)
       result.negate();
 }
 
@@ -354,7 +355,7 @@ template <class T>
 void eval_tan(T& result, const T& x)
 {
    BOOST_STATIC_ASSERT_MSG(number_category<T>::value == number_kind_floating_point, "The tan function is only valid for floating point types.");
-   if(&result == &x)
+   if (&result == &x)
    {
       T temp;
       eval_tan(temp, x);
@@ -370,19 +371,19 @@ void eval_tan(T& result, const T& x)
 template <class T>
 void hyp2F1(T& result, const T& a, const T& b, const T& c, const T& x)
 {
-  // Compute the series representation of hyperg_2f1 taken from
-  // Abramowitz and Stegun 15.1.1.
-  // There are no checks on input range or parameter boundaries.
+   // Compute the series representation of hyperg_2f1 taken from
+   // Abramowitz and Stegun 15.1.1.
+   // There are no checks on input range or parameter boundaries.
 
    typedef typename boost::multiprecision::detail::canonical<boost::uint32_t, T>::type ui_type;
 
    T x_pow_n_div_n_fact(x);
-   T pochham_a         (a);
-   T pochham_b         (b);
-   T pochham_c         (c);
-   T ap                (a);
-   T bp                (b);
-   T cp                (c);
+   T pochham_a(a);
+   T pochham_b(b);
+   T pochham_c(c);
+   T ap(a);
+   T bp(b);
+   T cp(c);
 
    eval_multiply(result, pochham_a, pochham_b);
    eval_divide(result, pochham_c);
@@ -392,17 +393,18 @@ void hyp2F1(T& result, const T& a, const T& b, const T& c, const T& x)
    T lim;
    eval_ldexp(lim, result, 1 - boost::multiprecision::detail::digits2<number<T, et_on> >::value());
 
-   if(eval_get_sign(lim) < 0)
+   if (eval_get_sign(lim) < 0)
       lim.negate();
 
    ui_type n;
-   T term;
+   T       term;
 
-   const unsigned series_limit = 
-      boost::multiprecision::detail::digits2<number<T, et_on> >::value() < 100
-      ? 100 : boost::multiprecision::detail::digits2<number<T, et_on> >::value();
+   const unsigned series_limit =
+       boost::multiprecision::detail::digits2<number<T, et_on> >::value() < 100
+           ? 100
+           : boost::multiprecision::detail::digits2<number<T, et_on> >::value();
    // Series expansion of hyperg_2f1(a, b; c; x).
-   for(n = 2; n < series_limit; ++n)
+   for (n = 2; n < series_limit; ++n)
    {
       eval_multiply(x_pow_n_div_n_fact, x);
       eval_divide(x_pow_n_div_n_fact, n);
@@ -419,12 +421,12 @@ void hyp2F1(T& result, const T& a, const T& b, const T& c, const T& x)
       eval_multiply(term, x_pow_n_div_n_fact);
       eval_add(result, term);
 
-      if(eval_get_sign(term) < 0)
+      if (eval_get_sign(term) < 0)
          term.negate();
-      if(lim.compare(term) >= 0)
+      if (lim.compare(term) >= 0)
          break;
    }
-   if(n > series_limit)
+   if (n > series_limit)
       BOOST_THROW_EXCEPTION(std::runtime_error("H2F1 failed to converge."));
 }
 
@@ -433,23 +435,23 @@ void eval_asin(T& result, const T& x)
 {
    BOOST_STATIC_ASSERT_MSG(number_category<T>::value == number_kind_floating_point, "The asin function is only valid for floating point types.");
    typedef typename boost::multiprecision::detail::canonical<boost::uint32_t, T>::type ui_type;
-   typedef typename mpl::front<typename T::float_types>::type fp_type;
+   typedef typename mpl::front<typename T::float_types>::type                          fp_type;
 
-   if(&result == &x)
+   if (&result == &x)
    {
       T t(x);
       eval_asin(result, t);
       return;
    }
 
-   switch(eval_fpclassify(x))
+   switch (eval_fpclassify(x))
    {
    case FP_NAN:
    case FP_INFINITE:
-      if(std::numeric_limits<number<T, et_on> >::has_quiet_NaN)
+      if (std::numeric_limits<number<T, et_on> >::has_quiet_NaN)
       {
          result = std::numeric_limits<number<T, et_on> >::quiet_NaN().backend();
-         errno = EDOM;
+         errno  = EDOM;
       }
       else
          BOOST_THROW_EXCEPTION(std::domain_error("Result is undefined or complex and there is no NaN for this number type."));
@@ -457,37 +459,37 @@ void eval_asin(T& result, const T& x)
    case FP_ZERO:
       result = x;
       return;
-   default: ;
+   default:;
    }
 
    const bool b_neg = eval_get_sign(x) < 0;
 
    T xx(x);
-   if(b_neg)
+   if (b_neg)
       xx.negate();
 
    int c = xx.compare(ui_type(1));
-   if(c > 0)
+   if (c > 0)
    {
-      if(std::numeric_limits<number<T, et_on> >::has_quiet_NaN)
+      if (std::numeric_limits<number<T, et_on> >::has_quiet_NaN)
       {
          result = std::numeric_limits<number<T, et_on> >::quiet_NaN().backend();
-         errno = EDOM;
+         errno  = EDOM;
       }
       else
          BOOST_THROW_EXCEPTION(std::domain_error("Result is undefined or complex and there is no NaN for this number type."));
       return;
    }
-   else if(c == 0)
+   else if (c == 0)
    {
       result = get_constant_pi<T>();
       eval_ldexp(result, result, -1);
-      if(b_neg)
+      if (b_neg)
          result.negate();
       return;
    }
 
-   if(xx.compare(fp_type(1e-4)) < 0)
+   if (xx.compare(fp_type(1e-4)) < 0)
    {
       // http://functions.wolfram.com/ElementaryFunctions/ArcSin/26/01/01/
       eval_multiply(xx, xx);
@@ -498,7 +500,7 @@ void eval_asin(T& result, const T& x)
       eval_multiply(result, x);
       return;
    }
-   else if(xx.compare(fp_type(1 - 1e-4f)) > 0)
+   else if (xx.compare(fp_type(1 - 1e-4f)) > 0)
    {
       T dx1;
       T t1, t2;
@@ -513,7 +515,7 @@ void eval_asin(T& result, const T& x)
       eval_ldexp(t1, get_constant_pi<T>(), -1);
       result.negate();
       eval_add(result, t1);
-      if(b_neg)
+      if (b_neg)
          result.negate();
       return;
    }
@@ -528,16 +530,16 @@ void eval_asin(T& result, const T& x)
 
    result = (guess_type)(std::asin(dd));
 
-   // Newton-Raphson iteration, we should double our precision with each iteration, 
+   // Newton-Raphson iteration, we should double our precision with each iteration,
    // in practice this seems to not quite work in all cases... so terminate when we
-   // have at least 2/3 of the digits correct on the assumption that the correction 
+   // have at least 2/3 of the digits correct on the assumption that the correction
    // we've just added will finish the job...
 
    boost::intmax_t current_precision = eval_ilogb(result);
-   boost::intmax_t target_precision = current_precision - 1 - (std::numeric_limits<number<T> >::digits * 2) / 3;
+   boost::intmax_t target_precision  = current_precision - 1 - (std::numeric_limits<number<T> >::digits * 2) / 3;
 
    // Newton-Raphson iteration
-   while(current_precision > target_precision)
+   while (current_precision > target_precision)
    {
       T sine, cosine;
       eval_sin(sine, result);
@@ -546,10 +548,10 @@ void eval_asin(T& result, const T& x)
       eval_divide(sine, cosine);
       eval_subtract(result, sine);
       current_precision = eval_ilogb(sine);
-      if(current_precision <= (std::numeric_limits<typename T::exponent_type>::min)() + 1)
+      if (current_precision <= (std::numeric_limits<typename T::exponent_type>::min)() + 1)
          break;
    }
-   if(b_neg)
+   if (b_neg)
       result.negate();
 }
 
@@ -559,14 +561,14 @@ inline void eval_acos(T& result, const T& x)
    BOOST_STATIC_ASSERT_MSG(number_category<T>::value == number_kind_floating_point, "The acos function is only valid for floating point types.");
    typedef typename boost::multiprecision::detail::canonical<boost::uint32_t, T>::type ui_type;
 
-   switch(eval_fpclassify(x))
+   switch (eval_fpclassify(x))
    {
    case FP_NAN:
    case FP_INFINITE:
-      if(std::numeric_limits<number<T, et_on> >::has_quiet_NaN)
+      if (std::numeric_limits<number<T, et_on> >::has_quiet_NaN)
       {
          result = std::numeric_limits<number<T, et_on> >::quiet_NaN().backend();
-         errno = EDOM;
+         errno  = EDOM;
       }
       else
          BOOST_THROW_EXCEPTION(std::domain_error("Result is undefined or complex and there is no NaN for this number type."));
@@ -580,20 +582,20 @@ inline void eval_acos(T& result, const T& x)
    eval_abs(result, x);
    int c = result.compare(ui_type(1));
 
-   if(c > 0)
+   if (c > 0)
    {
-      if(std::numeric_limits<number<T, et_on> >::has_quiet_NaN)
+      if (std::numeric_limits<number<T, et_on> >::has_quiet_NaN)
       {
          result = std::numeric_limits<number<T, et_on> >::quiet_NaN().backend();
-         errno = EDOM;
+         errno  = EDOM;
       }
       else
          BOOST_THROW_EXCEPTION(std::domain_error("Result is undefined or complex and there is no NaN for this number type."));
       return;
    }
-   else if(c == 0)
+   else if (c == 0)
    {
-      if(eval_get_sign(x) < 0)
+      if (eval_get_sign(x) < 0)
          result = get_constant_pi<T>();
       else
          result = ui_type(0);
@@ -611,21 +613,21 @@ template <class T>
 void eval_atan(T& result, const T& x)
 {
    BOOST_STATIC_ASSERT_MSG(number_category<T>::value == number_kind_floating_point, "The atan function is only valid for floating point types.");
-   typedef typename boost::multiprecision::detail::canonical<boost::int32_t, T>::type si_type;
+   typedef typename boost::multiprecision::detail::canonical<boost::int32_t, T>::type  si_type;
    typedef typename boost::multiprecision::detail::canonical<boost::uint32_t, T>::type ui_type;
-   typedef typename mpl::front<typename T::float_types>::type fp_type;
+   typedef typename mpl::front<typename T::float_types>::type                          fp_type;
 
-   switch(eval_fpclassify(x))
+   switch (eval_fpclassify(x))
    {
    case FP_NAN:
       result = x;
-      errno = EDOM;
+      errno  = EDOM;
       return;
    case FP_ZERO:
       result = x;
       return;
    case FP_INFINITE:
-      if(eval_get_sign(x) < 0)
+      if (eval_get_sign(x) < 0)
       {
          eval_ldexp(result, get_constant_pi<T>(), -1);
          result.negate();
@@ -633,16 +635,16 @@ void eval_atan(T& result, const T& x)
       else
          eval_ldexp(result, get_constant_pi<T>(), -1);
       return;
-   default: ;
+   default:;
    }
 
    const bool b_neg = eval_get_sign(x) < 0;
 
    T xx(x);
-   if(b_neg)
+   if (b_neg)
       xx.negate();
 
-   if(xx.compare(fp_type(0.1)) < 0)
+   if (xx.compare(fp_type(0.1)) < 0)
    {
       T t1, t2, t3;
       t1 = ui_type(1);
@@ -655,7 +657,7 @@ void eval_atan(T& result, const T& x)
       return;
    }
 
-   if(xx.compare(fp_type(10)) > 0)
+   if (xx.compare(fp_type(10)) > 0)
    {
       T t1, t2, t3;
       t1 = fp_type(0.5f);
@@ -665,31 +667,30 @@ void eval_atan(T& result, const T& x)
       eval_divide(xx, si_type(-1), xx);
       hyp2F1(result, t1, t2, t3, xx);
       eval_divide(result, x);
-      if(!b_neg)
+      if (!b_neg)
          result.negate();
       eval_ldexp(t1, get_constant_pi<T>(), -1);
       eval_add(result, t1);
-      if(b_neg)
+      if (b_neg)
          result.negate();
       return;
    }
 
-
    // Get initial estimate using standard math function atan.
    fp_type d;
    eval_convert_to(&d, xx);
    result = fp_type(std::atan(d));
 
-   // Newton-Raphson iteration, we should double our precision with each iteration, 
+   // Newton-Raphson iteration, we should double our precision with each iteration,
    // in practice this seems to not quite work in all cases... so terminate when we
-   // have at least 2/3 of the digits correct on the assumption that the correction 
+   // have at least 2/3 of the digits correct on the assumption that the correction
    // we've just added will finish the job...
 
    boost::intmax_t current_precision = eval_ilogb(result);
-   boost::intmax_t target_precision = current_precision - 1 - (std::numeric_limits<number<T> >::digits * 2) / 3;
+   boost::intmax_t target_precision  = current_precision - 1 - (std::numeric_limits<number<T> >::digits * 2) / 3;
 
    T s, c, t;
-   while(current_precision > target_precision)
+   while (current_precision > target_precision)
    {
       eval_sin(s, result);
       eval_cos(c, result);
@@ -698,10 +699,10 @@ void eval_atan(T& result, const T& x)
       eval_multiply(s, t, c);
       eval_add(result, s);
       current_precision = eval_ilogb(s);
-      if(current_precision <= (std::numeric_limits<typename T::exponent_type>::min)() + 1)
+      if (current_precision <= (std::numeric_limits<typename T::exponent_type>::min)() + 1)
          break;
    }
-   if(b_neg)
+   if (b_neg)
       result.negate();
 }
 
@@ -709,13 +710,13 @@ template <class T>
 void eval_atan2(T& result, const T& y, const T& x)
 {
    BOOST_STATIC_ASSERT_MSG(number_category<T>::value == number_kind_floating_point, "The atan2 function is only valid for floating point types.");
-   if(&result == &y)
+   if (&result == &y)
    {
       T temp(y);
       eval_atan2(result, temp, x);
       return;
    }
-   else if(&result == &x)
+   else if (&result == &x)
    {
       T temp(x);
       eval_atan2(result, y, temp);
@@ -724,82 +725,82 @@ void eval_atan2(T& result, const T& y, const T& x)
 
    typedef typename boost::multiprecision::detail::canonical<boost::uint32_t, T>::type ui_type;
 
-   switch(eval_fpclassify(y))
+   switch (eval_fpclassify(y))
    {
    case FP_NAN:
       result = y;
-      errno = EDOM;
+      errno  = EDOM;
       return;
    case FP_ZERO:
+   {
+      if (eval_signbit(x))
       {
-         if(eval_signbit(x))
-         {
-            result = get_constant_pi<T>();
-            if(eval_signbit(y))
-               result.negate();
-         }
-         else
-         {
-            result = y; // Note we allow atan2(0,0) to be +-zero, even though it's mathematically undefined
-         }
-         return;
+         result = get_constant_pi<T>();
+         if (eval_signbit(y))
+            result.negate();
       }
+      else
+      {
+         result = y; // Note we allow atan2(0,0) to be +-zero, even though it's mathematically undefined
+      }
+      return;
+   }
    case FP_INFINITE:
+   {
+      if (eval_fpclassify(x) == FP_INFINITE)
       {
-         if(eval_fpclassify(x) == FP_INFINITE)
+         if (eval_signbit(x))
          {
-            if(eval_signbit(x))
-            {
-               // 3Pi/4
-               eval_ldexp(result, get_constant_pi<T>(), -2);
-               eval_subtract(result, get_constant_pi<T>());
-               if(eval_get_sign(y) >= 0)
-                  result.negate();
-            }
-            else
-            {
-               // Pi/4
-               eval_ldexp(result, get_constant_pi<T>(), -2);
-               if(eval_get_sign(y) < 0)
-                  result.negate();
-            }
+            // 3Pi/4
+            eval_ldexp(result, get_constant_pi<T>(), -2);
+            eval_subtract(result, get_constant_pi<T>());
+            if (eval_get_sign(y) >= 0)
+               result.negate();
          }
          else
          {
-            eval_ldexp(result, get_constant_pi<T>(), -1);
-            if(eval_get_sign(y) < 0)
+            // Pi/4
+            eval_ldexp(result, get_constant_pi<T>(), -2);
+            if (eval_get_sign(y) < 0)
                result.negate();
          }
-         return;
       }
+      else
+      {
+         eval_ldexp(result, get_constant_pi<T>(), -1);
+         if (eval_get_sign(y) < 0)
+            result.negate();
+      }
+      return;
+   }
    }
 
-   switch(eval_fpclassify(x))
+   switch (eval_fpclassify(x))
    {
    case FP_NAN:
       result = x;
-      errno = EDOM;
+      errno  = EDOM;
       return;
    case FP_ZERO:
-      {
-         eval_ldexp(result, get_constant_pi<T>(), -1);
-         if(eval_get_sign(y) < 0)
-            result.negate();
-         return;
-      }
+   {
+      eval_ldexp(result, get_constant_pi<T>(), -1);
+      if (eval_get_sign(y) < 0)
+         result.negate();
+      return;
+   }
    case FP_INFINITE:
-      if(eval_get_sign(x) > 0)
+      if (eval_get_sign(x) > 0)
          result = ui_type(0);
       else
          result = get_constant_pi<T>();
-      if(eval_get_sign(y) < 0)
+      if (eval_get_sign(y) < 0)
          result.negate();
       return;
    }
 
    T xx;
    eval_divide(xx, y, x);
-   if(eval_get_sign(xx) < 0)
+   if (eval_get_sign(xx) < 0)
       xx.negate();
 
    eval_atan(result, xx);
@@ -808,33 +809,33 @@ void eval_atan2(T& result, const T& y, const T& x)
    const bool y_neg = eval_get_sign(y) < 0;
    const bool x_neg = eval_get_sign(x) < 0;
 
-   if(y_neg != x_neg)
+   if (y_neg != x_neg)
       result.negate();
 
-   if(x_neg)
+   if (x_neg)
    {
-      if(y_neg)
+      if (y_neg)
          eval_subtract(result, get_constant_pi<T>());
       else
          eval_add(result, get_constant_pi<T>());
    }
 }
-template<class T, class A> 
+template <class T, class A>
 inline typename enable_if<is_arithmetic<A>, void>::type eval_atan2(T& result, const T& x, const A& a)
 {
-   typedef typename boost::multiprecision::detail::canonical<A, T>::type canonical_type;
+   typedef typename boost::multiprecision::detail::canonical<A, T>::type          canonical_type;
    typedef typename mpl::if_<is_same<A, canonical_type>, T, canonical_type>::type cast_type;
-   cast_type c;
+   cast_type                                                                      c;
    c = a;
    eval_atan2(result, x, c);
 }
 
-template<class T, class A> 
+template <class T, class A>
 inline typename enable_if<is_arithmetic<A>, void>::type eval_atan2(T& result, const A& x, const T& a)
 {
-   typedef typename boost::multiprecision::detail::canonical<A, T>::type canonical_type;
+   typedef typename boost::multiprecision::detail::canonical<A, T>::type          canonical_type;
    typedef typename mpl::if_<is_same<A, canonical_type>, T, canonical_type>::type cast_type;
-   cast_type c;
+   cast_type                                                                      c;
    c = x;
    eval_atan2(result, c, a);
 }