Boost logo

Boost-Commit :

From: john_at_[hidden]
Date: 2008-01-05 05:10:03


Author: johnmaddock
Date: 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
New Revision: 42474
URL: http://svn.boost.org/trac/boost/changeset/42474

Log:
Added initial support for interval arithmetic and mpfr extended-precision support.
Added:
   sandbox/interval_math_toolkit/boost/math/concepts/real_type_concept.hpp (contents, props changed)
   sandbox/interval_math_toolkit/libs/math/test/interval_concept_check.cpp (contents, props changed)
   sandbox/interval_math_toolkit/libs/math/test/mpfr_concept_check.cpp (contents, props changed)
   sandbox/interval_math_toolkit/libs/math/test/ntl_concept_check.cpp (contents, props changed)
Text files modified:
   sandbox/interval_math_toolkit/boost/math/bindings/mpfr.hpp | 1
   sandbox/interval_math_toolkit/boost/math/concepts/distributions.hpp | 2
   sandbox/interval_math_toolkit/boost/math/constants/constants.hpp | 16 +++
   sandbox/interval_math_toolkit/boost/math/distributions/binomial.hpp | 6
   sandbox/interval_math_toolkit/boost/math/distributions/cauchy.hpp | 2
   sandbox/interval_math_toolkit/boost/math/distributions/chi_squared.hpp | 6
   sandbox/interval_math_toolkit/boost/math/distributions/detail/inv_discrete_quantile.hpp | 9 +
   sandbox/interval_math_toolkit/boost/math/distributions/exponential.hpp | 4
   sandbox/interval_math_toolkit/boost/math/distributions/extreme_value.hpp | 4
   sandbox/interval_math_toolkit/boost/math/distributions/fisher_f.hpp | 16 ++--
   sandbox/interval_math_toolkit/boost/math/distributions/lognormal.hpp | 2
   sandbox/interval_math_toolkit/boost/math/distributions/negative_binomial.hpp | 10 +-
   sandbox/interval_math_toolkit/boost/math/distributions/normal.hpp | 6
   sandbox/interval_math_toolkit/boost/math/distributions/poisson.hpp | 8 +-
   sandbox/interval_math_toolkit/boost/math/distributions/rayleigh.hpp | 4
   sandbox/interval_math_toolkit/boost/math/distributions/triangular.hpp | 6
   sandbox/interval_math_toolkit/boost/math/distributions/weibull.hpp | 24 +++---
   sandbox/interval_math_toolkit/boost/math/policies/error_handling.hpp | 20 ++++-
   sandbox/interval_math_toolkit/boost/math/policies/policy.hpp | 13 +-
   sandbox/interval_math_toolkit/boost/math/special_functions/bessel.hpp | 12 +-
   sandbox/interval_math_toolkit/boost/math/special_functions/beta.hpp | 61 ++++++++-------
   sandbox/interval_math_toolkit/boost/math/special_functions/detail/bessel_ik.hpp | 10 +-
   sandbox/interval_math_toolkit/boost/math/special_functions/detail/bessel_jy.hpp | 27 +++---
   sandbox/interval_math_toolkit/boost/math/special_functions/detail/bessel_jy_asym.hpp | 12 ++-
   sandbox/interval_math_toolkit/boost/math/special_functions/detail/gamma_inva.hpp | 10 +-
   sandbox/interval_math_toolkit/boost/math/special_functions/detail/ibeta_inv_ab.hpp | 27 +++---
   sandbox/interval_math_toolkit/boost/math/special_functions/detail/ibeta_inverse.hpp | 19 ++--
   sandbox/interval_math_toolkit/boost/math/special_functions/detail/igamma_inverse.hpp | 11 +-
   sandbox/interval_math_toolkit/boost/math/special_functions/detail/lgamma_small.hpp | 4
   sandbox/interval_math_toolkit/boost/math/special_functions/detail/t_distribution_inv.hpp | 93 +++++++++++++++--------
   sandbox/interval_math_toolkit/boost/math/special_functions/detail/unchecked_factorial.hpp | 19 ++++
   sandbox/interval_math_toolkit/boost/math/special_functions/ellint_1.hpp | 6
   sandbox/interval_math_toolkit/boost/math/special_functions/ellint_2.hpp | 4
   sandbox/interval_math_toolkit/boost/math/special_functions/ellint_3.hpp | 8 +-
   sandbox/interval_math_toolkit/boost/math/special_functions/ellint_rj.hpp | 2
   sandbox/interval_math_toolkit/boost/math/special_functions/erf.hpp | 6
   sandbox/interval_math_toolkit/boost/math/special_functions/expint.hpp | 4
   sandbox/interval_math_toolkit/boost/math/special_functions/expm1.hpp | 2
   sandbox/interval_math_toolkit/boost/math/special_functions/factorials.hpp | 14 +-
   sandbox/interval_math_toolkit/boost/math/special_functions/gamma.hpp | 60 +++++++++------
   sandbox/interval_math_toolkit/boost/math/special_functions/legendre.hpp | 2
   sandbox/interval_math_toolkit/boost/math/special_functions/log1p.hpp | 4
   sandbox/interval_math_toolkit/boost/math/special_functions/sin_pi.hpp | 2
   sandbox/interval_math_toolkit/boost/math/special_functions/spherical_harmonic.hpp | 4
   sandbox/interval_math_toolkit/boost/math/special_functions/sqrt1pm1.hpp | 2
   sandbox/interval_math_toolkit/boost/math/special_functions/zeta.hpp | 4
   sandbox/interval_math_toolkit/boost/math/tools/config.hpp | 153 +++++++++++++++++++++++++++++++++++++++
   sandbox/interval_math_toolkit/boost/math/tools/fraction.hpp | 8 +-
   sandbox/interval_math_toolkit/boost/math/tools/minima.hpp | 4
   sandbox/interval_math_toolkit/boost/math/tools/rational.hpp | 12 +-
   sandbox/interval_math_toolkit/boost/math/tools/real_cast.hpp | 7 +
   sandbox/interval_math_toolkit/boost/math/tools/roots.hpp | 151 +++++++++++++++++++-------------------
   sandbox/interval_math_toolkit/boost/math/tools/series.hpp | 12 +-
   sandbox/interval_math_toolkit/boost/math/tools/test.hpp | 115 +++++++++++++++++++++++++++++
   sandbox/interval_math_toolkit/boost/math/tools/toms748_solve.hpp | 16 ++-
   sandbox/interval_math_toolkit/libs/math/test/handle_test_result.hpp | 47 ++++++++---
   sandbox/interval_math_toolkit/libs/math/test/log1p_expm1_test.cpp | 11 ++
   sandbox/interval_math_toolkit/libs/math/test/powm1_sqrtp1m1_test.cpp | 9 ++
   sandbox/interval_math_toolkit/libs/math/test/test_beta.cpp | 9 ++
   sandbox/interval_math_toolkit/libs/math/test/test_binomial_coeff.cpp | 9 ++
   sandbox/interval_math_toolkit/libs/math/test/test_carlson.cpp | 16 +++
   sandbox/interval_math_toolkit/libs/math/test/test_cbrt.cpp | 9 ++
   sandbox/interval_math_toolkit/libs/math/test/test_digamma.cpp | 9 ++
   sandbox/interval_math_toolkit/libs/math/test/test_ellint_1.cpp | 11 ++
   sandbox/interval_math_toolkit/libs/math/test/test_ellint_2.cpp | 10 ++
   sandbox/interval_math_toolkit/libs/math/test/test_ellint_3.cpp | 10 ++
   sandbox/interval_math_toolkit/libs/math/test/test_erf.cpp | 9 ++
   sandbox/interval_math_toolkit/libs/math/test/test_expint.cpp | 9 ++
   sandbox/interval_math_toolkit/libs/math/test/test_gamma.cpp | 9 ++
   sandbox/interval_math_toolkit/libs/math/test/test_hermite.cpp | 11 ++
   sandbox/interval_math_toolkit/libs/math/test/test_ibeta.cpp | 9 ++
   sandbox/interval_math_toolkit/libs/math/test/test_ibeta_inv.cpp | 12 +++
   72 files changed, 904 insertions(+), 371 deletions(-)

Modified: sandbox/interval_math_toolkit/boost/math/bindings/mpfr.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/bindings/mpfr.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/bindings/mpfr.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -705,7 +705,6 @@
     return value;
 }
 
-
 } // namespace detail
 
 }}

Modified: sandbox/interval_math_toolkit/boost/math/concepts/distributions.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/concepts/distributions.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/concepts/distributions.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -162,6 +162,7 @@
       v = quantile(complement(dist, d));
       v = hazard(dist, d);
       v = chf(dist, d);
+#ifndef TEST_MPFR
       long double ld = 1;
       v = cdf(dist, ld);
       v = cdf(complement(dist, ld));
@@ -170,6 +171,7 @@
       v = quantile(complement(dist, ld));
       v = hazard(dist, ld);
       v = chf(dist, ld);
+#endif
       int i = 1;
       v = cdf(dist, i);
       v = cdf(complement(dist, i));

Added: sandbox/interval_math_toolkit/boost/math/concepts/real_type_concept.hpp
==============================================================================
--- (empty file)
+++ sandbox/interval_math_toolkit/boost/math/concepts/real_type_concept.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -0,0 +1,111 @@
+// Copyright John Maddock 2007-8.
+// Use, modification and distribution are subject to the
+// Boost Software License, Version 1.0. (See accompanying file
+// LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
+
+#ifndef BOOST_MATH_REAL_TYPE_CONCEPT_HPP
+#define BOOST_MATH_REAL_TYPE_CONCEPT_HPP
+
+#ifdef BOOST_MSVC
+#pragma warning(push)
+#pragma warning(disable: 4100)
+#pragma warning(disable: 4510)
+#pragma warning(disable: 4610)
+#endif
+#include <boost/concept_check.hpp>
+#ifdef BOOST_MSVC
+#pragma warning(pop)
+#endif
+#include <boost/math/tools/config.hpp>
+#include <boost/math/tools/precision.hpp>
+
+
+namespace boost{ namespace math{ namespace concepts{
+
+template <class RealType>
+struct RealTypeConcept
+{
+ template <class Other>
+ void check_binary_ops(Other o)
+ {
+ using namespace boost::numeric::interval_math_compare;
+ RealType r(o);
+ r = o;
+ r -= o;
+ r += o;
+ r *= o;
+ r /= o;
+ r = r - o;
+ r = o - r;
+ r = r + o;
+ r = o + r;
+ r = o * r;
+ r = r * o;
+ r = r / o;
+ r = o / r;
+ bool b;
+ b = r == o;
+ b = o == r;
+ b = r != o;
+ b = o != r;
+ b = r <= o;
+ b = o <= r;
+ b = r >= o;
+ b = o >= r;
+ b = r < o;
+ b = o < r;
+ b = r > o;
+ b = o > r;
+ }
+
+ void constraints()
+ {
+ BOOST_MATH_STD_USING
+ using namespace boost::numeric::interval_math_compare;
+
+ RealType r;
+ check_binary_ops(r);
+ check_binary_ops(0.5f);
+ check_binary_ops(0.5);
+ //check_binary_ops(0.5L);
+ check_binary_ops(1);
+ //check_binary_ops(1u);
+ check_binary_ops(1L);
+ //check_binary_ops(1uL);
+#ifndef BOOST_HAS_LONG_LONG
+ check_binary_ops(1LL);
+#endif
+ RealType r2 = +r;
+ r2 = -r;
+
+ r2 = fabs(r);
+ r2 = abs(r);
+ r2 = ceil(r);
+ r2 = floor(r);
+ r2 = exp(r);
+ r2 = pow(r, r2);
+ r2 = sqrt(r);
+ r2 = log(r);
+ r2 = cos(r);
+ r2 = sin(r);
+ r2 = tan(r);
+ r2 = asin(r);
+ r2 = acos(r);
+ r2 = atan(r);
+ int i;
+ r2 = ldexp(r, i);
+ r2 = frexp(r, &i);
+ i = boost::math::tools::digits<RealType>();
+ r2 = boost::math::tools::max_value<RealType>();
+ r2 = boost::math::tools::min_value<RealType>();
+ r2 = boost::math::tools::log_max_value<RealType>();
+ r2 = boost::math::tools::log_min_value<RealType>();
+ r2 = boost::math::tools::epsilon<RealType>();
+ }
+}; // struct DistributionConcept
+
+
+}}} // namespaces
+
+#endif
+

Modified: sandbox/interval_math_toolkit/boost/math/constants/constants.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/constants/constants.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/constants/constants.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -35,12 +35,24 @@
     // (This is necessary because you can't use a numeric constant
     // since even a long double might not have enough digits).
 
+// TODO: this does not create propper intervals!
 
   #define BOOST_DEFINE_MATH_CONSTANT(name, x, y, exp)\
+ template <class T> T name(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE(T));\
+ namespace detail{\
+ template <class T> inline T name(boost::mpl::false_& BOOST_MATH_APPEND_EXPLICIT_TEMPLATE_TYPE(T))\
+ {\
+ static const T result = ::boost::lexical_cast<T>(BOOST_STRINGIZE(BOOST_JOIN(BOOST_JOIN(x, y), BOOST_JOIN(e, exp))));\
+ return result;\
+ }\
+ template <class T> inline T name(boost::mpl::true_& BOOST_MATH_APPEND_EXPLICIT_TEMPLATE_TYPE(T))\
+ {\
+ return boost::math::constants::name<typename T::base_type>();\
+ }\
+ }\
    template <class T> inline T name(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE(T))\
    {\
- static const T result = ::boost::lexical_cast<T>(BOOST_STRINGIZE(BOOST_JOIN(BOOST_JOIN(x, y), BOOST_JOIN(e, exp))));\
- return result;\
+ return detail::name<T>(boost::numeric::is_interval<T>());\
    }\
    template <> inline float name<float>(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(float))\
    { return BOOST_JOIN(BOOST_JOIN(x, BOOST_JOIN(e, exp)), F); }\

Modified: sandbox/interval_math_toolkit/boost/math/distributions/binomial.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/distributions/binomial.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/distributions/binomial.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -173,7 +173,7 @@
             // kurtosis:
             // T k = (1 - 6 * sf * (1 - sf) ) / (n * sf * (1 - sf));
             // Get the inverse of a std normal distribution:
- T x = boost::math::erfc_inv(p > q ? 2 * q : 2 * p, pol) * constants::root_two<T>();
+ T x = boost::math::erfc_inv(p > q ? static_cast<T>(2 * q) : static_cast<T>(2 * p), pol) * constants::root_two<T>();
             // Set the sign:
             if(p < 0.5)
                x = -x;
@@ -649,13 +649,13 @@
       template <class RealType, class Policy>
       inline RealType quantile(const binomial_distribution<RealType, Policy>& dist, const RealType& p)
       {
- return binomial_detail::quantile_imp(dist, p, 1-p);
+ return binomial_detail::quantile_imp(dist, p, static_cast<RealType>(1-p));
       } // quantile
 
       template <class RealType, class Policy>
       RealType quantile(const complemented2_type<binomial_distribution<RealType, Policy>, RealType>& c)
       {
- return binomial_detail::quantile_imp(c.dist, 1-c.param, c.param);
+ return binomial_detail::quantile_imp(c.dist, static_cast<RealType>(1-c.param), c.param);
       } // quantile
 
       template <class RealType, class Policy>

Modified: sandbox/interval_math_toolkit/boost/math/distributions/cauchy.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/distributions/cauchy.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/distributions/cauchy.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -140,7 +140,7 @@
       return location;
    }
    result = -scale / tan(constants::pi<RealType>() * P);
- return complement ? location - result : location + result;
+ return complement ? static_cast<RealType>(location - result) : static_cast<RealType>(location + result);
 } // quantile
 
 } // namespace detail

Modified: sandbox/interval_math_toolkit/boost/math/distributions/chi_squared.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/distributions/chi_squared.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/distributions/chi_squared.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -107,7 +107,7 @@
       }
    }
 
- return gamma_p_derivative(degrees_of_freedom / 2, chi_square / 2, Policy()) / 2;
+ return gamma_p_derivative(static_cast<RealType>(degrees_of_freedom / 2), static_cast<RealType>(chi_square / 2), Policy()) / 2;
 } // pdf
 
 template <class RealType, class Policy>
@@ -128,7 +128,7 @@
          function, "Chi Square parameter was %1%, but must be > 0 !", chi_square, Policy());
    }
 
- return boost::math::gamma_p(degrees_of_freedom / 2, chi_square / 2, Policy());
+ return boost::math::gamma_p(static_cast<RealType>(degrees_of_freedom / 2), static_cast<RealType>(chi_square / 2), Policy());
 } // cdf
 
 template <class RealType, class Policy>
@@ -165,7 +165,7 @@
          function, "Chi Square parameter was %1%, but must be > 0 !", chi_square, Policy());
    }
 
- return boost::math::gamma_q(degrees_of_freedom / 2, chi_square / 2, Policy());
+ return boost::math::gamma_q(static_cast<RealType>(degrees_of_freedom / 2), static_cast<RealType>(chi_square / 2), Policy());
 }
 
 template <class RealType, class Policy>

Modified: sandbox/interval_math_toolkit/boost/math/distributions/detail/inv_discrete_quantile.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/distributions/detail/inv_discrete_quantile.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/distributions/detail/inv_discrete_quantile.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -127,8 +127,9 @@
          }
          else
          {
+ using std::max;
             b = a;
- a = (std::max)(b - 1, value_type(0));
+ a = (max)(static_cast<value_type>(b - 1), value_type(0));
             if(a < min_bound)
                a = min_bound;
             fa = f(a);
@@ -158,7 +159,8 @@
       }
       else
       {
- b = (std::max)(a - adder, value_type(0));
+ using std::max;
+ b = max BOOST_PREVENT_MACRO_SUBSTITUTION(static_cast<value_type>(a - adder), value_type(0));
          if(b < min_bound)
             b = min_bound;
       }
@@ -182,7 +184,8 @@
          }
          else
          {
- b = (std::max)(a - adder, value_type(0));
+ using std::max;
+ b = max BOOST_PREVENT_MACRO_SUBSTITUTION(static_cast<value_type>(a - adder), value_type(0));
             if(b < min_bound)
                b = min_bound;
          }

Modified: sandbox/interval_math_toolkit/boost/math/distributions/exponential.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/distributions/exponential.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/distributions/exponential.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -124,7 +124,7 @@
       return result;
    if(0 == detail::verify_exp_x(function, x, &result, Policy()))
       return result;
- result = -boost::math::expm1(-x * lambda, Policy());
+ result = -boost::math::expm1(static_cast<RealType>(-x * lambda), Policy());
 
    return result;
 } // cdf
@@ -148,7 +148,7 @@
    if(p == 1)
       return policies::raise_overflow_error<RealType>(function, 0, Policy());
 
- result = -boost::math::log1p(-p, Policy()) / lambda;
+ result = -boost::math::log1p(static_cast<RealType>(-p), Policy()) / lambda;
    return result;
 } // quantile
 

Modified: sandbox/interval_math_toolkit/boost/math/distributions/extreme_value.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/distributions/extreme_value.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/distributions/extreme_value.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -153,7 +153,7 @@
    if(0 == detail::verify_scale_b("boost::math::cdf(const extreme_value_distribution<%1%>&, %1%)", b, &result, Policy()))
       return result;
 
- result = -boost::math::expm1(-exp((a-c.param)/b), Policy());
+ result = -boost::math::expm1(static_cast<RealType>(-exp((a-c.param)/b)), Policy());
 
    return result;
 }
@@ -179,7 +179,7 @@
    if(q == 1)
       return -policies::raise_overflow_error<RealType>(function, 0, Policy());
 
- result = a - log(-boost::math::log1p(-q, Policy())) * b;
+ result = a - log(-boost::math::log1p(static_cast<RealType>(-q), Policy())) * b;
 
    return result;
 }

Modified: sandbox/interval_math_toolkit/boost/math/distributions/fisher_f.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/distributions/fisher_f.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/distributions/fisher_f.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -115,13 +115,13 @@
    if(v1x > df2)
    {
       result = (df2 * df1) / ((df2 + v1x) * (df2 + v1x));
- result *= ibeta_derivative(df2 / 2, df1 / 2, df2 / (df2 + v1x), Policy());
+ result *= ibeta_derivative(static_cast<RealType>(df2 / 2), static_cast<RealType>(df1 / 2), static_cast<RealType>(df2 / (df2 + v1x)), Policy());
    }
    else
    {
       result = df2 + df1 * x;
       result = (result * df1 - x * df1 * df1) / (result * result);
- result *= ibeta_derivative(df1 / 2, df2 / 2, v1x / (df2 + v1x), Policy());
+ result *= ibeta_derivative(static_cast<RealType>(df1 / 2), static_cast<RealType>(df2 / 2), static_cast<RealType>(v1x / (df2 + v1x)), Policy());
    }
    return result;
 } // pdf
@@ -157,8 +157,8 @@
    // to switch things around so we're passing 1-z instead.
    //
    return v1x > df2
- ? boost::math::ibetac(df2 / 2, df1 / 2, df2 / (df2 + v1x), Policy())
- : boost::math::ibeta(df1 / 2, df2 / 2, v1x / (df2 + v1x), Policy());
+ ? boost::math::ibetac(static_cast<RealType>(df2 / 2), static_cast<RealType>(df1 / 2), static_cast<RealType>(df2 / (df2 + v1x)), Policy())
+ : boost::math::ibeta(static_cast<RealType>(df1 / 2), static_cast<RealType>(df2 / 2), static_cast<RealType>(v1x / (df2 + v1x)), Policy());
 } // cdf
 
 template <class RealType, class Policy>
@@ -179,7 +179,7 @@
 
    RealType x, y;
 
- x = boost::math::ibeta_inv(df1 / 2, df2 / 2, p, &y, Policy());
+ x = boost::math::ibeta_inv(static_cast<RealType>(df1 / 2), static_cast<RealType>(df2 / 2), p, &y, Policy());
 
    return df2 * x / (df1 * y);
 } // quantile
@@ -216,8 +216,8 @@
    // to switch things around so we're passing 1-z instead.
    //
    return v1x > df2
- ? boost::math::ibeta(df2 / 2, df1 / 2, df2 / (df2 + v1x), Policy())
- : boost::math::ibetac(df1 / 2, df2 / 2, v1x / (df2 + v1x), Policy());
+ ? boost::math::ibeta(static_cast<RealType>(df2 / 2), static_cast<RealType>(df1 / 2), static_cast<RealType>(df2 / (df2 + v1x)), Policy())
+ : boost::math::ibetac(static_cast<RealType>(df1 / 2), static_cast<RealType>(df2 / 2), static_cast<RealType>(v1x / (df2 + v1x)), Policy());
 }
 
 template <class RealType, class Policy>
@@ -239,7 +239,7 @@
 
    RealType x, y;
 
- x = boost::math::ibetac_inv(df1 / 2, df2 / 2, p, &y, Policy());
+ x = boost::math::ibetac_inv(static_cast<RealType>(df1 / 2), static_cast<RealType>(df2 / 2), p, &y, Policy());
 
    return df2 * x / (df1 * y);
 }

Modified: sandbox/interval_math_toolkit/boost/math/distributions/lognormal.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/distributions/lognormal.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/distributions/lognormal.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -221,7 +221,7 @@
    if(0 == detail::check_scale("boost::math::variance(const lognormal_distribution<%1%>&)", sigma, &result, Policy()))
       return result;
 
- return boost::math::expm1(sigma * sigma, Policy()) * exp(2 * mu + sigma * sigma);
+ return boost::math::expm1(static_cast<RealType>(sigma * sigma), Policy()) * exp(2 * mu + sigma * sigma);
 }
 
 template <class RealType, class Policy>

Modified: sandbox/interval_math_toolkit/boost/math/distributions/negative_binomial.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/distributions/negative_binomial.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/distributions/negative_binomial.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -431,6 +431,7 @@
       // BUT Cephes/CodeCogs says: finds argument p (0 to 1) such that cdf(k, n, p) = y
       static const char* function = "boost::math::quantile(const negative_binomial_distribution<%1%>&, %1%)";
       BOOST_MATH_STD_USING // ADL of std functions.
+ using std::min;
 
       RealType p = dist.success_fraction();
       RealType r = dist.successes();
@@ -468,14 +469,14 @@
       RealType guess = 0;
       RealType factor = 5;
       if(r * r * r * P * p > 0.005)
- guess = detail::inverse_negative_binomial_cornish_fisher(r, p, 1-p, P, 1-P, Policy());
+ guess = detail::inverse_negative_binomial_cornish_fisher(r, p, static_cast<RealType>(1-p), P, static_cast<RealType>(1-P), Policy());
 
       if(guess < 10)
       {
          //
          // Cornish-Fisher Negative binomial approximation not accurate in this area:
          //
- guess = (std::min)(r * 2, RealType(10));
+ guess = min BOOST_PREVENT_MACRO_SUBSTITUTION(static_cast<RealType>(r * 2), RealType(10));
       }
       else
          factor = (1-P < sqrt(tools::epsilon<RealType>())) ? 2 : (guess < 20 ? 1.2f : 1.1f);
@@ -503,6 +504,7 @@
        // complement of the probability Q = 1 - P.
        static const char* function = "boost::math::quantile(const negative_binomial_distribution<%1%>&, %1%)";
        BOOST_MATH_STD_USING
+ using std::min;
 
        // Error checks:
        RealType Q = c.param;
@@ -545,14 +547,14 @@
        RealType guess = 0;
        RealType factor = 5;
        if(r * r * r * (1-Q) * p > 0.005)
- guess = detail::inverse_negative_binomial_cornish_fisher(r, p, 1-p, 1-Q, Q, Policy());
+ guess = detail::inverse_negative_binomial_cornish_fisher(r, p, static_cast<RealType>(1-p), static_cast<RealType>(1-Q), Q, Policy());
 
        if(guess < 10)
        {
           //
           // Cornish-Fisher Negative binomial approximation not accurate in this area:
           //
- guess = (std::min)(r * 2, RealType(10));
+ guess = min BOOST_PREVENT_MACRO_SUBSTITUTION(static_cast<RealType>(r * 2), RealType(10));
        }
        else
           factor = (Q < sqrt(tools::epsilon<RealType>())) ? 2 : (guess < 20 ? 1.2f : 1.1f);

Modified: sandbox/interval_math_toolkit/boost/math/distributions/normal.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/distributions/normal.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/distributions/normal.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -166,7 +166,7 @@
      return result;
    }
    RealType diff = (x - mean) / (sd * constants::root_two<RealType>());
- result = boost::math::erfc(-diff, Policy()) / 2;
+ result = boost::math::erfc(static_cast<RealType>(-diff), Policy()) / 2;
    return result;
 } // cdf
 
@@ -187,7 +187,7 @@
    if(false == detail::check_probability(function, p, &result, Policy()))
       return result;
 
- result= boost::math::erfc_inv(2 * p, Policy());
+ result= boost::math::erfc_inv(static_cast<RealType>(2 * p), Policy());
    result = -result;
    result *= sd * constants::root_two<RealType>();
    result += mean;
@@ -247,7 +247,7 @@
    RealType q = c.param;
    if(false == detail::check_probability(function, q, &result, Policy()))
       return result;
- result = boost::math::erfc_inv(2 * q, Policy());
+ result = boost::math::erfc_inv(static_cast<RealType>(2 * q), Policy());
    result *= sd * constants::root_two<RealType>();
    result += mean;
    return result;

Modified: sandbox/interval_math_toolkit/boost/math/distributions/poisson.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/distributions/poisson.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/distributions/poisson.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -347,7 +347,7 @@
         // (e ^ -mean * mean ^ k) / k!
         // == exp(log(e ^ -mean) + log (mean ^ k) - lgamma(k+1))
         // exp( -mean + log(mean) * k - lgamma(k+1))
- return exp(-mean + log(mean) * k - boost::math::lgamma(k+1, Policy()));
+ return exp(-mean + log(mean) * k - boost::math::lgamma(static_cast<RealType>(k+1), Policy()));
         // return gamma_p_derivative(k+1, mean); // equivalent & also passes tests.
       }
     } // pdf
@@ -444,7 +444,7 @@
       }
       if (k == 0)
       { // Avoid repeated checks on k and mean in gamma_p.
- return -boost::math::expm1(-mean, Policy());
+ return -boost::math::expm1(static_cast<RealType>(-mean), Policy());
       }
       // Unlike un-complemented cdf (sum from 0 to k),
       // can't use finite sum from k+1 to infinity for small integral k,
@@ -493,7 +493,7 @@
       if(z < 1)
          guess = z;
       else
- guess = boost::math::detail::inverse_poisson_cornish_fisher(z, p, 1-p, Policy());
+ guess = boost::math::detail::inverse_poisson_cornish_fisher(z, p, static_cast<RealType>(1-p), Policy());
       if(z > 5)
       {
          if(z > 1000)
@@ -561,7 +561,7 @@
       if(z < 1)
          guess = z;
       else
- guess = boost::math::detail::inverse_poisson_cornish_fisher(z, 1-q, q, Policy());
+ guess = boost::math::detail::inverse_poisson_cornish_fisher(z, static_cast<RealType>(1-q), q, Policy());
       if(z > 5)
       {
          if(z > 1000)

Modified: sandbox/interval_math_toolkit/boost/math/distributions/rayleigh.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/distributions/rayleigh.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/distributions/rayleigh.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -129,7 +129,7 @@
    {
       return result;
    }
- result = -boost::math::expm1(-x * x / ( 2 * sigma * sigma), Policy());
+ result = -boost::math::expm1(static_cast<RealType>(-x * x / ( 2 * sigma * sigma)), Policy());
    return result;
 } // cdf
 
@@ -154,7 +154,7 @@
    {
      return policies::raise_overflow_error<RealType>(function, 0, Policy());
    }
- result = sqrt(-2 * sigma * sigma * boost::math::log1p(-p, Policy()));
+ result = sqrt(-2 * sigma * sigma * boost::math::log1p(static_cast<RealType>(-p), Policy()));
    return result;
 } // quantile
 

Modified: sandbox/interval_math_toolkit/boost/math/distributions/triangular.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/distributions/triangular.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/distributions/triangular.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -212,11 +212,11 @@
     }
     if (x == lower)
     { // (mode - lower) == 0 which would lead to divide by zero!
- return (mode == lower) ? 2 / (upper - lower) : 0;
+ return (mode == lower) ? RealType(2) / (upper - lower) : RealType(0);
     }
     else if (x == upper)
     {
- return (mode == upper) ? 2 / (upper - lower) : 0;
+ return (mode == upper) ? RealType(2) / (upper - lower) : RealType(0);
     }
     else if (x <= mode)
     {
@@ -468,7 +468,7 @@
     RealType mode = dist.mode();
     RealType upper = dist.upper();
     RealType result; // of checks.
- if(false == detail::check_triangular(function,lower, mode, upper, &result, Policy()))
+ if(false == boost::math::detail::check_triangular(function,lower, mode, upper, &result, Policy()))
     {
       return result;
     }

Modified: sandbox/interval_math_toolkit/boost/math/distributions/weibull.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/distributions/weibull.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/distributions/weibull.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -158,7 +158,7 @@
    if(false == detail::check_weibull_x(function, x, &result, Policy()))
       return result;
 
- result = -boost::math::expm1(-pow(x / scale, shape), Policy());
+ result = -boost::math::expm1(static_cast<RealType>(-pow(x / scale, shape)), Policy());
 
    return result;
 }
@@ -182,7 +182,7 @@
    if(p == 1)
       return policies::raise_overflow_error<RealType>(function, 0, Policy());
 
- result = scale * pow(-boost::math::log1p(-p, Policy()), 1 / shape);
+ result = scale * pow(-boost::math::log1p(static_cast<RealType>(-p), Policy()), 1 / shape);
 
    return result;
 }
@@ -247,7 +247,7 @@
    if(false == detail::check_weibull(function, scale, shape, &result, Policy()))
       return result;
 
- result = scale * boost::math::tgamma(1 + 1 / shape, Policy());
+ result = scale * boost::math::tgamma(static_cast<RealType>(1 + 1 / shape), Policy());
    return result;
 }
 
@@ -264,9 +264,9 @@
    {
       return result;
    }
- result = boost::math::tgamma(1 + 1 / shape, Policy());
+ result = boost::math::tgamma(static_cast<RealType>(1 + 1 / shape), Policy());
    result *= -result;
- result += boost::math::tgamma(1 + 2 / shape, Policy());
+ result += boost::math::tgamma(static_cast<RealType>(1 + 2 / shape), Policy());
    result *= scale * scale;
    return result;
 }
@@ -327,9 +327,9 @@
    }
    RealType g1, g2, g3, d;
 
- g1 = boost::math::tgamma(1 + 1 / shape, Policy());
- g2 = boost::math::tgamma(1 + 2 / shape, Policy());
- g3 = boost::math::tgamma(1 + 3 / shape, Policy());
+ g1 = boost::math::tgamma(static_cast<RealType>(1 + 1 / shape), Policy());
+ g2 = boost::math::tgamma(static_cast<RealType>(1 + 2 / shape), Policy());
+ g3 = boost::math::tgamma(static_cast<RealType>(1 + 3 / shape), Policy());
    d = pow(g2 - g1 * g1, RealType(1.5));
 
    result = (2 * g1 * g1 * g1 - 3 * g1 * g2 + g3) / d;
@@ -352,10 +352,10 @@
 
    RealType g1, g2, g3, g4, d, g1_2, g1_4;
 
- g1 = boost::math::tgamma(1 + 1 / shape, Policy());
- g2 = boost::math::tgamma(1 + 2 / shape, Policy());
- g3 = boost::math::tgamma(1 + 3 / shape, Policy());
- g4 = boost::math::tgamma(1 + 4 / shape, Policy());
+ g1 = boost::math::tgamma(static_cast<RealType>(1 + 1 / shape), Policy());
+ g2 = boost::math::tgamma(static_cast<RealType>(1 + 2 / shape), Policy());
+ g3 = boost::math::tgamma(static_cast<RealType>(1 + 3 / shape), Policy());
+ g4 = boost::math::tgamma(static_cast<RealType>(1 + 4 / shape), Policy());
    g1_2 = g1 * g1;
    g1_4 = g1_2 * g1_2;
    d = g2 - g1_2;

Modified: sandbox/interval_math_toolkit/boost/math/policies/error_handling.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/policies/error_handling.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/policies/error_handling.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -18,6 +18,7 @@
 #include <boost/math/policies/policy.hpp>
 #include <boost/math/tools/precision.hpp>
 #include <boost/cstdint.hpp>
+#include <boost/tr1/tuple.hpp>
 #ifdef BOOST_MSVC
 # pragma warning(push) // Quiet warnings in boost/format.hpp
 # pragma warning(disable: 4996) // _SCL_SECURE_NO_DEPRECATE
@@ -55,6 +56,17 @@
 template <class T>
 T user_evaluation_error(const char* function, const char* message, const T& val);
 
+template <class T>
+inline T make_printable_type(const T& v)
+{ return v; }
+
+template <class T, class Policy>
+inline std::tr1::tuple<T, T>
+make_printable_type(const boost::numeric::interval<T, Policy>& inter)
+{
+ return std::tr1::tuple<T, T>(inter.lower(), inter.upper());
+}
+
 namespace detail
 {
 
@@ -89,7 +101,7 @@
   msg += message;
 
   int prec = 2 + (boost::math::policies::digits<T, boost::math::policies::policy<> >() * 30103UL) / 100000UL;
- msg = (boost::format(msg) % boost::io::group(std::setprecision(prec), val)).str();
+ msg = (boost::format(msg) % boost::io::group(std::setprecision(prec), make_printable_type(val))).str();
 
   E e(msg);
   boost::throw_exception(e);
@@ -429,7 +441,7 @@
 inline bool check_overflow(T val, R* result, const char* function, const Policy& pol)
 {
    BOOST_MATH_STD_USING
- if(fabs(val) > tools::max_value<R>())
+ if(boost::math::tools::maybe_greater(fabs(val), static_cast<T>(tools::max_value<R>())))
    {
       *result = static_cast<R>(boost::math::policies::detail::raise_overflow_error<R>(function, 0, pol));
       return true;
@@ -439,7 +451,7 @@
 template <class R, class T, class Policy>
 inline bool check_underflow(T val, R* result, const char* function, const Policy& pol)
 {
- if((val != 0) && (static_cast<R>(val) == 0))
+ if(!boost::math::tools::maybe_equal(val, T(0)) && boost::math::tools::maybe_equal(static_cast<R>(val), R(0)))
    {
       *result = static_cast<R>(boost::math::policies::detail::raise_underflow_error<R>(function, 0, pol));
       return true;
@@ -450,7 +462,7 @@
 inline bool check_denorm(T val, R* result, const char* function, const Policy& pol)
 {
    BOOST_MATH_STD_USING
- if((fabs(val) < static_cast<T>(tools::min_value<R>())) && (static_cast<R>(val) != 0))
+ if(boost::math::tools::maybe_less(fabs(val), static_cast<T>(tools::min_value<R>())) && !boost::math::tools::maybe_equal(static_cast<R>(val), R(0)))
    {
       *result = static_cast<R>(boost::math::policies::detail::raise_denorm_error<R>(function, 0, static_cast<R>(val), pol));
       return true;

Modified: sandbox/interval_math_toolkit/boost/math/policies/policy.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/policies/policy.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/policies/policy.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -730,25 +730,26 @@
 template <class Real, class Policy>
 struct precision
 {
+ typedef typename boost::numeric::interval_base_type<Real>::type base_type;
 #ifndef __BORLANDC__
    typedef typename Policy::precision_type precision_type;
    typedef typename mpl::if_c<
- ((::std::numeric_limits<Real>::is_specialized == 0) || (::std::numeric_limits<Real>::digits == 0)),
+ ((::std::numeric_limits<base_type>::is_specialized == 0) || (::std::numeric_limits<base_type>::digits == 0)),
       // Possibly unknown precision:
       precision_type,
       typename mpl::if_c<
- ((::std::numeric_limits<Real>::digits <= precision_type::value)
+ ((::std::numeric_limits<base_type>::digits <= precision_type::value)
          || (Policy::precision_type::value <= 0)),
          // Default case, full precision for RealType:
- digits2< ::std::numeric_limits<Real>::digits>,
+ digits2< ::std::numeric_limits<base_type>::digits>,
          // User customised precision:
          precision_type
>::type
>::type type;
 #else
    typedef typename Policy::precision_type precision_type;
- typedef mpl::int_< ::std::numeric_limits<Real>::digits> digits_t;
- typedef mpl::bool_< ::std::numeric_limits<Real>::is_specialized> spec_t;
+ typedef mpl::int_< ::std::numeric_limits<base_type>::digits> digits_t;
+ typedef mpl::bool_< ::std::numeric_limits<base_type>::is_specialized> spec_t;
    typedef typename mpl::if_<
       mpl::or_<mpl::equal_to<spec_t, mpl::false_>, mpl::equal_to<digits_t, mpl::int_<0> > >,
       // Possibly unknown precision:
@@ -756,7 +757,7 @@
       typename mpl::if_<
          mpl::or_<mpl::less_equal<digits_t, precision_type>, mpl::less_equal<precision_type, mpl::int_<0> > >,
          // Default case, full precision for RealType:
- digits2< ::std::numeric_limits<Real>::digits>,
+ digits2< ::std::numeric_limits<base_type>::digits>,
          // User customised precision:
          precision_type
>::type

Modified: sandbox/interval_math_toolkit/boost/math/special_functions/bessel.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/special_functions/bessel.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/special_functions/bessel.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -41,7 +41,7 @@
    {
       BOOST_MATH_STD_USING
       mult = x / 2;
- term = pow(mult, v) / boost::math::tgamma(v+1, Policy());
+ term = pow(mult, v) / boost::math::tgamma(static_cast<T>(v+1), Policy());
       mult *= -mult;
    }
    T operator()()
@@ -68,7 +68,7 @@
    {
       BOOST_MATH_STD_USING
       mult = x / 2;
- term = pow(mult, T(v)) / boost::math::tgamma(v+1+T(0.5f), Policy());
+ term = pow(mult, T(v)) / boost::math::tgamma(static_cast<T>(v+1+T(0.5f)), Policy());
       mult *= -mult;
    }
    T operator()()
@@ -126,7 +126,7 @@
       // better have integer v:
       if(floor(v) == v)
       {
- T r = cyl_bessel_j_imp(v, -x, t, pol);
+ T r = cyl_bessel_j_imp(v, static_cast<T>(-x), t, pol);
          if(tools::real_cast<int>(v) & 1)
             r = -r;
          return r;
@@ -207,7 +207,7 @@
    // Default case is just a naive evaluation of the definition:
    //
    return sqrt(constants::pi<T>() / (2 * x))
- * cyl_bessel_j_imp(T(n)+T(0.5f), x, bessel_no_int_tag(), pol);
+ * cyl_bessel_j_imp(static_cast<T>(T(n)+T(0.5f)), x, bessel_no_int_tag(), pol);
 }
 
 template <class T, class Policy>
@@ -225,7 +225,7 @@
       // better have integer v:
       if(floor(v) == v)
       {
- T r = cyl_bessel_i_imp(v, -x, pol);
+ T r = cyl_bessel_i_imp(v, static_cast<T>(-x), pol);
          if(tools::real_cast<int>(v) & 1)
             r = -r;
          return r;
@@ -378,7 +378,7 @@
    if(x < 2 * tools::min_value<T>())
       return -policies::raise_overflow_error<T>(function, 0, pol);
 
- T result = cyl_neumann_imp(T(v)+0.5f, x, bessel_no_int_tag(), pol);
+ T result = cyl_neumann_imp(static_cast<T>(T(v)+0.5f), x, bessel_no_int_tag(), pol);
    T tx = sqrt(constants::pi<T>() / (2 * x));
 
    if((tx > 1) && (tools::max_value<T>() / tx < result))

Modified: sandbox/interval_math_toolkit/boost/math/special_functions/beta.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/special_functions/beta.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/special_functions/beta.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -89,7 +89,7 @@
    {
       // Special case where the base of the power term is close to 1
       // compute (1+x)^y instead:
- result *= exp(ambh * boost::math::log1p(-b / cgh, pol));
+ result *= exp(ambh * boost::math::log1p(static_cast<T>(-b / cgh), pol));
    }
    else
    {
@@ -154,9 +154,10 @@
       std::swap(a, b);
 
    // set integration limits:
- T la = (std::max)(T(10), a);
- T lb = (std::max)(T(10), b);
- T lc = (std::max)(T(10), a+b);
+ using std::max;
+ T la = max BOOST_PREVENT_MACRO_SUBSTITUTION(T(10), a);
+ T lb = max BOOST_PREVENT_MACRO_SUBSTITUTION(T(10), b);
+ T lc = max BOOST_PREVENT_MACRO_SUBSTITUTION(T(10), static_cast<T>(a+b));
 
    // calculate the fraction parts:
    T sa = detail::lower_gamma_series(a, la, pol) / a;
@@ -200,6 +201,8 @@
                         const Policy& pol)
 {
    BOOST_MATH_STD_USING
+ using std::min;
+ using std::max;
 
    if(!normalised)
    {
@@ -221,11 +224,11 @@
    // l1 and l2 are the base of the exponents minus one:
    T l1 = (x * b - y * agh) / agh;
    T l2 = (y * a - x * bgh) / bgh;
- if(((std::min)(fabs(l1), fabs(l2)) < 0.2))
+ if((min BOOST_PREVENT_MACRO_SUBSTITUTION(fabs(l1), fabs(l2)) < 0.2))
    {
       // when the base of the exponent is very near 1 we get really
       // gross errors unless extra care is taken:
- if((l1 * l2 > 0) || ((std::min)(a, b) < 1))
+ if((l1 * l2 > 0) || (min BOOST_PREVENT_MACRO_SUBSTITUTION(a, b) < 1))
       {
          //
          // This first branch handles the simple cases where either:
@@ -249,7 +252,7 @@
          else
             result *= pow((y * cgh) / bgh, b);
       }
- else if((std::max)(fabs(l1), fabs(l2)) < 0.5)
+ else if(max BOOST_PREVENT_MACRO_SUBSTITUTION(fabs(l1), fabs(l2)) < 0.5)
       {
          //
          // Both exponents are near one and both the exponents are
@@ -274,14 +277,14 @@
          T ratio = b / a;
          if((small_a && (ratio * l2 < 0.1)) || (!small_a && (l1 / ratio > 0.1)))
          {
- T l3 = boost::math::expm1(ratio * boost::math::log1p(l2, pol), pol);
+ T l3 = boost::math::expm1(static_cast<T>(ratio * boost::math::log1p(l2, pol)), pol);
             l3 = l1 + l3 + l3 * l1;
             l3 = a * boost::math::log1p(l3, pol);
             result *= exp(l3);
          }
          else
          {
- T l3 = boost::math::expm1(boost::math::log1p(l1, pol) / ratio, pol);
+ T l3 = boost::math::expm1(static_cast<T>(boost::math::log1p(l1, pol) / ratio), pol);
             l3 = l2 + l3 + l3 * l2;
             l3 = b * boost::math::log1p(l3, pol);
             result *= exp(l3);
@@ -460,7 +463,7 @@
       T cgh = c + L::g() - T(0.5);
       result = L::lanczos_sum_expG_scaled(c) / (L::lanczos_sum_expG_scaled(a) * L::lanczos_sum_expG_scaled(b));
       if(a * b < bgh * 10)
- result *= exp((b - 0.5f) * boost::math::log1p(a / bgh, pol));
+ result *= exp((b - 0.5f) * boost::math::log1p(static_cast<T>(a / bgh), pol));
       else
          result *= pow(cgh / bgh, b - 0.5f);
       result *= pow(x * cgh / agh, a);
@@ -717,7 +720,7 @@
    T t = a + bm1 / 2;
    T lx, u;
    if(y < 0.35)
- lx = boost::math::log1p(-y, pol);
+ lx = boost::math::log1p(static_cast<T>(-y), pol);
    else
       lx = log(x);
    u = -t * lx;
@@ -821,7 +824,7 @@
    BOOST_MATH_STD_USING // ADL of std names
    T result = pow(x, n);
    T term = result;
- for(unsigned i = tools::real_cast<unsigned>(n - 1); i > k; --i)
+ for(unsigned i = tools::real_cast<unsigned>(static_cast<T>(n - 1)); i > k; --i)
    {
       term *= ((i + 1) * y) / ((n - i) * x) ;
       result += term;
@@ -843,6 +846,8 @@
    static const char* function = "boost::math::ibeta<%1%>(%1%, %1%, %1%)";
    typedef typename lanczos::lanczos<T, Policy>::type lanczos_type;
    BOOST_MATH_STD_USING // for ADL of std math functions.
+ using std::max;
+ using std::min;
 
    bool invert = inv;
    T fract;
@@ -873,7 +878,7 @@
    {
       if(p_derivative)
       {
- *p_derivative = (a == 1) ? 1 : (a < 1) ? tools::max_value<T>() / 2 : tools::min_value<T>() * 2;
+ *p_derivative = (a == 1) ? 1 : (a < 1) ? static_cast<T>(tools::max_value<T>() / 2) : static_cast<T>(tools::min_value<T>() * 2);
       }
       return (invert ? (normalised ? 1 : boost::math::beta(a, b, pol)) : 0);
    }
@@ -881,12 +886,12 @@
    {
       if(p_derivative)
       {
- *p_derivative = (b == 1) ? 1 : (b < 1) ? tools::max_value<T>() / 2 : tools::min_value<T>() * 2;
+ *p_derivative = (b == 1) ? 1 : (b < 1) ? static_cast<T>(tools::max_value<T>() / 2) : static_cast<T>(tools::min_value<T>() * 2);
       }
       return (invert == 0 ? (normalised ? 1 : boost::math::beta(a, b, pol)) : 0);
    }
 
- if((std::min)(a, b) <= 1)
+ if(min BOOST_PREVENT_MACRO_SUBSTITUTION(a, b) <= 1)
    {
       if(x > 0.5)
       {
@@ -894,10 +899,10 @@
          std::swap(x, y);
          invert = !invert;
       }
- if((std::max)(a, b) <= 1)
+ if(max BOOST_PREVENT_MACRO_SUBSTITUTION(a, b) <= 1)
       {
          // Both a,b < 1:
- if((a >= (std::min)(T(0.2), b)) || (pow(x, a) <= 0.9))
+ if((a >= min BOOST_PREVENT_MACRO_SUBSTITUTION(T(0.2), b)) || (pow(x, a) <= 0.9))
          {
             if(!invert)
                fract = ibeta_series(a, b, x, T(0), lanczos_type(), normalised, p_derivative, y, pol);
@@ -930,7 +935,7 @@
                T prefix;
                if(!normalised)
                {
- prefix = rising_factorial_ratio(a+b, a, 20);
+ prefix = rising_factorial_ratio(static_cast<T>(a+b), a, 20);
                }
                else
                {
@@ -938,12 +943,12 @@
                }
                fract = ibeta_a_step(a, b, x, y, 20, pol, normalised, p_derivative);
                if(!invert)
- fract = beta_small_b_large_a_series(a + 20, b, x, y, fract, prefix, pol, normalised);
+ fract = beta_small_b_large_a_series(static_cast<T>(a + 20), b, x, y, fract, prefix, pol, normalised);
                else
                {
                   fract -= (normalised ? 1 : boost::math::beta(a, b, pol));
                   invert = false;
- fract = -beta_small_b_large_a_series(a + 20, b, x, y, fract, prefix, pol, normalised);
+ fract = -beta_small_b_large_a_series(static_cast<T>(a + 20), b, x, y, fract, prefix, pol, normalised);
                }
             }
          }
@@ -996,7 +1001,7 @@
                T prefix;
                if(!normalised)
                {
- prefix = rising_factorial_ratio(a+b, a, 20);
+ prefix = rising_factorial_ratio(static_cast<T>(a+b), a, 20);
                }
                else
                {
@@ -1004,12 +1009,12 @@
                }
                fract = ibeta_a_step(a, b, x, y, 20, pol, normalised, p_derivative);
                if(!invert)
- fract = beta_small_b_large_a_series(a + 20, b, x, y, fract, prefix, pol, normalised);
+ fract = beta_small_b_large_a_series(static_cast<T>(a + 20), b, x, y, fract, prefix, pol, normalised);
                else
                {
                   fract -= (normalised ? 1 : boost::math::beta(a, b, pol));
                   invert = false;
- fract = -beta_small_b_large_a_series(a + 20, b, x, y, fract, prefix, pol, normalised);
+ fract = -beta_small_b_large_a_series(static_cast<T>(a + 20), b, x, y, fract, prefix, pol, normalised);
                }
             }
          }
@@ -1059,14 +1064,14 @@
          else if(a > 15)
          {
             // sidestep so we can use the series representation:
- int n = static_cast<int>(boost::math::tools::real_cast<long double>(floor(b)));
+ int n = tools::real_cast<int>(static_cast<T>(floor(b)));
             if(n == b)
                --n;
             T bbar = b - n;
             T prefix;
             if(!normalised)
             {
- prefix = rising_factorial_ratio(a+bbar, bbar, n);
+ prefix = rising_factorial_ratio(static_cast<T>(a+bbar), bbar, n);
             }
             else
             {
@@ -1081,7 +1086,7 @@
             // the formula here for the non-normalised case is tricky to figure
             // out (for me!!), and requires two pochhammer calculations rather
             // than one, so leave it for now....
- int n = static_cast<int>(boost::math::tools::real_cast<long double>(floor(b)));
+ int n = tools::real_cast<int>(static_cast<T>(floor(b)));
             T bbar = b - n;
             if(bbar <= 0)
             {
@@ -1093,7 +1098,7 @@
             if(invert)
                fract -= (normalised ? 1 : boost::math::beta(a, b, pol));
             //fract = ibeta_series(a+20, bbar, x, fract, l, normalised, p_derivative, y);
- fract = beta_small_b_large_a_series(a+20, bbar, x, y, fract, T(1), pol, normalised);
+ fract = beta_small_b_large_a_series(static_cast<T>(a+20), bbar, x, y, fract, T(1), pol, normalised);
             if(invert)
             {
                fract = -fract;
@@ -1166,7 +1171,7 @@
    // Now the regular cases:
    //
    typedef typename lanczos::lanczos<T, Policy>::type lanczos_type;
- T f1 = ibeta_power_terms(a, b, x, 1 - x, lanczos_type(), true, pol);
+ T f1 = ibeta_power_terms(a, b, x, static_cast<T>(1 - x), lanczos_type(), true, pol);
    T y = (1 - x) * x;
 
    if(f1 == 0)

Modified: sandbox/interval_math_toolkit/boost/math/special_functions/detail/bessel_ik.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/special_functions/detail/bessel_ik.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/special_functions/detail/bessel_ik.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -42,17 +42,17 @@
     BOOST_ASSERT(abs(v) <= 0.5f);
 
     T gp = boost::math::tgamma1pm1(v, pol);
- T gm = boost::math::tgamma1pm1(-v, pol);
+ T gm = boost::math::tgamma1pm1(static_cast<T>(-v), pol);
 
     a = log(x / 2);
     b = exp(v * a);
     sigma = -a * v;
     c = abs(v) < tools::epsilon<T>() ?
- 1 : boost::math::sin_pi(v) / (v * pi<T>());
+ static_cast<T>(1) : static_cast<T>(boost::math::sin_pi(v) / (v * pi<T>()));
     d = abs(sigma) < tools::epsilon<T>() ?
- 1 : sinh(sigma) / sigma;
+ static_cast<T>(1) : static_cast<T>(sinh(sigma) / sigma);
     gamma1 = abs(v) < tools::epsilon<T>() ?
- -euler<T>() : (0.5f / v) * (gp - gm) * c;
+ static_cast<T>(-euler<T>()) : static_cast<T>((0.5f / v) * (gp - gm) * c);
     gamma2 = (2 + gp + gm) * c / 2;
 
     // initial values
@@ -234,7 +234,7 @@
         v = -v; // v is non-negative from here
         kind |= need_k;
     }
- n = tools::real_cast<unsigned>(v + 0.5f);
+ n = tools::real_cast<unsigned>(static_cast<T>(v + 0.5f));
     u = v - n; // -1/2 <= u < 1/2
     BOOST_MATH_INSTRUMENT_VARIABLE(n);
     BOOST_MATH_INSTRUMENT_VARIABLE(u);

Modified: sandbox/interval_math_toolkit/boost/math/special_functions/detail/bessel_jy.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/special_functions/detail/bessel_jy.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/special_functions/detail/bessel_jy.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -46,21 +46,21 @@
     BOOST_ASSERT(fabs(v) <= 0.5f); // precondition for using this routine
 
     T gp = boost::math::tgamma1pm1(v, pol);
- T gm = boost::math::tgamma1pm1(-v, pol);
+ T gm = boost::math::tgamma1pm1(static_cast<T>(-v), pol);
     T spv = boost::math::sin_pi(v, pol);
- T spv2 = boost::math::sin_pi(v/2, pol);
+ T spv2 = boost::math::sin_pi(static_cast<T>(v/2), pol);
     T xp = pow(x/2, v);
 
     a = log(x / 2);
     sigma = -a * v;
     d = abs(sigma) < tools::epsilon<T>() ?
         T(1) : sinh(sigma) / sigma;
- e = abs(v) < tools::epsilon<T>() ? v*pi<T>()*pi<T>() / 2
- : 2 * spv2 * spv2 / v;
+ e = abs(v) < tools::epsilon<T>() ? static_cast<T>(v*pi<T>()*pi<T>() / 2)
+ : static_cast<T>(2 * spv2 * spv2 / v);
 
- T g1 = (v == 0) ? -euler<T>() : (gp - gm) / ((1 + gp) * (1 + gm) * 2 * v);
+ T g1 = (v == 0) ? static_cast<T>(-euler<T>()) : static_cast<T>((gp - gm) / ((1 + gp) * (1 + gm) * 2 * v));
     T g2 = (2 + gp + gm) / ((1 + gp) * (1 + gm) * 2);
- T vspv = (fabs(v) < tools::epsilon<T>()) ? 1/constants::pi<T>() : v / spv;
+ T vspv = (fabs(v) < tools::epsilon<T>()) ? static_cast<T>(1/constants::pi<T>()) : static_cast<T>(v / spv);
     f = (g1 * cosh(sigma) - g2 * a * d) * 2 * vspv;
 
     p = vspv / (xp * (1 + gm));
@@ -118,7 +118,7 @@
     tolerance = 2 * tools::epsilon<T>();
     tiny = sqrt(tools::min_value<T>());
     C = f = tiny; // b0 = 0, replace with tiny
- D = 0.0L;
+ D = 0;
     for (k = 1; k < policies::get_max_series_iterations<Policy>() * 100; k++)
     {
         a = -1;
@@ -131,7 +131,7 @@
         delta = C * D;
         f *= delta;
         if (D < 0) { s = -s; }
- if (abs(delta - 1.0L) < tolerance)
+ if (abs(delta - 1) < tolerance)
         { break; }
     }
     policies::check_series_iterations("boost::math::bessel_jy<%1%>(%1%,%1%) in CF1_jy", k / 100, pol);
@@ -158,7 +158,7 @@
     typedef typename complex_trait<T>::type complex_type;
 
     complex_type C, D, f, a, b, delta, one(1);
- T tiny, zero(0.0L);
+ T tiny, zero(0);
     unsigned long k;
 
     // |x| >= |v|, CF2_jy converges rapidly
@@ -169,7 +169,7 @@
     // Lentz, Applied Optics, vol 15, 668 (1976)
     T tolerance = 2 * tools::epsilon<T>();
     tiny = sqrt(tools::min_value<T>());
- C = f = complex_type(-0.5f/x, 1.0L);
+ C = f = complex_type(static_cast<T>(-0.5f/x), 1);
     D = 0;
     for (k = 1; k < policies::get_max_series_iterations<Policy>(); k++)
     {
@@ -225,7 +225,7 @@
         v = -v; // v is non-negative from here
         kind = need_j|need_y; // need both for reflection formula
     }
- n = real_cast<unsigned>(v + 0.5L);
+ n = real_cast<unsigned>(static_cast<T>(v + 0.5f));
     u = v - n; // -1/2 <= u < 1/2
 
     if (x == 0)
@@ -279,7 +279,8 @@
            lim = asymptotic_bessel_y_limit<T>(tag_type());
            break;
         default:
- lim = (std::max)(
+ using std::max;
+ lim = max BOOST_PREVENT_MACRO_SUBSTITUTION(
               asymptotic_bessel_j_limit<T>(v, tag_type()),
               asymptotic_bessel_y_limit<T>(tag_type()));
            break;
@@ -289,7 +290,7 @@
            if(kind&need_y)
            {
               Yu = asymptotic_bessel_y_large_x_2(u, x);
- Yu1 = asymptotic_bessel_y_large_x_2(u + 1, x);
+ Yu1 = asymptotic_bessel_y_large_x_2(static_cast<T>(u + 1), x);
            }
            else
               Yu = std::numeric_limits<T>::quiet_NaN(); // any value will do, we're not using it.

Modified: sandbox/interval_math_toolkit/boost/math/special_functions/detail/bessel_jy_asym.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/special_functions/detail/bessel_jy_asym.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/special_functions/detail/bessel_jy_asym.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -208,28 +208,32 @@
 {
    // default case:
    BOOST_MATH_STD_USING
- T v2 = (std::max)(T(3), v * v);
+ using std::max;
+ T v2 = max BOOST_PREVENT_MACRO_SUBSTITUTION(T(3), static_cast<T>(v * v));
    return v2 / pow(100 * tools::epsilon<T>() / T(2e-5f), T(0.17f));
 }
 template <class T>
 inline T asymptotic_bessel_j_limit(const T& v, const mpl::int_<53>&)
 {
    // double case:
- T v2 = (std::max)(T(3), v * v);
+ using std::max;
+ T v2 = max BOOST_PREVENT_MACRO_SUBSTITUTION(T(3), v * v);
    return v2 * 33 /*73*/;
 }
 template <class T>
 inline T asymptotic_bessel_j_limit(const T& v, const mpl::int_<64>&)
 {
    // 80-bit extended-double case:
- T v2 = (std::max)(T(3), v * v);
+ using std::max;
+ T v2 = max BOOST_PREVENT_MACRO_SUBSTITUTION(T(3), v * v);
    return v2 * 121 /*266*/;
 }
 template <class T>
 inline T asymptotic_bessel_j_limit(const T& v, const mpl::int_<113>&)
 {
    // 128-bit long double case:
- T v2 = (std::max)(T(3), v * v);
+ using std::max;
+ T v2 = max BOOST_PREVENT_MACRO_SUBSTITUTION(T(3), v * v);
    return v2 * 39154 /*85700*/;
 }
 

Modified: sandbox/interval_math_toolkit/boost/math/special_functions/detail/gamma_inva.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/special_functions/detail/gamma_inva.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/special_functions/detail/gamma_inva.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -48,7 +48,7 @@
    // kurtosis:
    // T k = 1/lambda;
    // Get the inverse of a std normal distribution:
- T x = boost::math::erfc_inv(p > q ? 2 * q : 2 * p, pol) * constants::root_two<T>();
+ T x = boost::math::erfc_inv(p > q ? static_cast<T>(2 * q) : static_cast<T>(2 * p), pol) * constants::root_two<T>();
    // Set the sign:
    if(p < 0.5)
       x = -x;
@@ -153,25 +153,25 @@
 template <class T, class Policy>
 inline T gamma_p_inva(T x, T p, const Policy& pol)
 {
- return detail::gamma_inva_imp(x, p, 1 - p, pol);
+ return detail::gamma_inva_imp(x, p, static_cast<T>(1 - p), pol);
 }
 
 template <class T, class Policy>
 inline T gamma_q_inva(T x, T q, const Policy& pol)
 {
- return detail::gamma_inva_imp(x, 1 - q, q, pol);
+ return detail::gamma_inva_imp(x, static_cast<T>(1 - q), q, pol);
 }
 
 template <class T>
 inline T gamma_p_inva(T x, T p)
 {
- return detail::gamma_inva_imp(x, p, 1 - p, policies::policy<>());
+ return detail::gamma_inva_imp(x, p, static_cast<T>(1 - p), policies::policy<>());
 }
 
 template <class T>
 inline T gamma_q_inva(T x, T q)
 {
- return detail::gamma_inva_imp(x, 1 - q, q, policies::policy<>());
+ return detail::gamma_inva_imp(x, static_cast<T>(1 - q), q, policies::policy<>());
 }
 
 } // namespace math

Modified: sandbox/interval_math_toolkit/boost/math/special_functions/detail/ibeta_inv_ab.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/special_functions/detail/ibeta_inv_ab.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/special_functions/detail/ibeta_inv_ab.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -51,7 +51,7 @@
    // kurtosis:
    T k = (6 - sf * (5+sfc)) / (n * (sfc));
    // Get the inverse of a std normal distribution:
- T x = boost::math::erfc_inv(p > q ? 2 * q : 2 * p, pol) * constants::root_two<T>();
+ T x = boost::math::erfc_inv(p > q ? static_cast<T>(2 * q) : static_cast<T>(2 * p), pol) * constants::root_two<T>();
    // Set the sign:
    if(p < 0.5)
       x = -x;
@@ -74,6 +74,7 @@
 T ibeta_inv_ab_imp(const T& b, const T& z, const T& p, const T& q, bool swap_ab, const Policy& pol)
 {
    BOOST_MATH_STD_USING // for ADL of std lib math functions
+ using std::min;
    //
    // Special cases first:
    //
@@ -120,11 +121,11 @@
       //
       if((p < q) != swap_ab)
       {
- guess = (std::min)(b * 2, T(1));
+ guess = min BOOST_PREVENT_MACRO_SUBSTITUTION(static_cast<T>(b * 2), T(1));
       }
       else
       {
- guess = (std::min)(b / 2, T(1));
+ guess = min BOOST_PREVENT_MACRO_SUBSTITUTION(static_cast<T>(b / 2), T(1));
       }
    }
    if(n * n * n * u * sf > 0.005)
@@ -137,11 +138,11 @@
       //
       if((p < q) != swap_ab)
       {
- guess = (std::min)(b * 2, T(10));
+ guess = min BOOST_PREVENT_MACRO_SUBSTITUTION(static_cast<T>(b * 2), T(10));
       }
       else
       {
- guess = (std::min)(b / 2, T(10));
+ guess = min BOOST_PREVENT_MACRO_SUBSTITUTION(static_cast<T>(b / 2), T(10));
       }
    }
    else
@@ -162,49 +163,49 @@
 template <class T, class Policy>
 inline T ibeta_inva(T b, T x, T p, const Policy& pol)
 {
- return detail::ibeta_inv_ab_imp(b, x, p, 1 - p, false, pol);
+ return detail::ibeta_inv_ab_imp(b, x, p, static_cast<T>(1 - p), false, pol);
 }
 
 template <class T, class Policy>
 inline T ibetac_inva(T b, T x, T q, const Policy& pol)
 {
- return detail::ibeta_inv_ab_imp(b, x, 1 - q, q, false, pol);
+ return detail::ibeta_inv_ab_imp(b, x, static_cast<T>(1 - q), q, false, pol);
 }
 
 template <class T, class Policy>
 inline T ibeta_invb(T b, T x, T p, const Policy& pol)
 {
- return detail::ibeta_inv_ab_imp(b, x, p, 1 - p, true, pol);
+ return detail::ibeta_inv_ab_imp(b, x, p, static_cast<T>(1 - p), true, pol);
 }
 
 template <class T, class Policy>
 inline T ibetac_invb(T b, T x, T q, const Policy& pol)
 {
- return detail::ibeta_inv_ab_imp(b, x, 1 - q, q, true, pol);
+ return detail::ibeta_inv_ab_imp(b, x, static_cast<T>(1 - q), q, true, pol);
 }
 
 template <class T>
 inline T ibeta_inva(T b, T x, T p)
 {
- return detail::ibeta_inv_ab_imp(b, x, p, 1 - p, false, policies::policy<>());
+ return detail::ibeta_inv_ab_imp(b, x, p, static_cast<T>(1 - p), false, policies::policy<>());
 }
 
 template <class T>
 inline T ibetac_inva(T b, T x, T q)
 {
- return detail::ibeta_inv_ab_imp(b, x, 1 - q, q, false, policies::policy<>());
+ return detail::ibeta_inv_ab_imp(b, x, static_cast<T>(1 - q), q, false, policies::policy<>());
 }
 
 template <class T>
 inline T ibeta_invb(T b, T x, T p)
 {
- return detail::ibeta_inv_ab_imp(b, x, p, 1 - p, true, policies::policy<>());
+ return detail::ibeta_inv_ab_imp(b, x, p, static_cast<T>(1 - p), true, policies::policy<>());
 }
 
 template <class T>
 inline T ibetac_invb(T b, T x, T q)
 {
- return detail::ibeta_inv_ab_imp(b, x, 1 - q, q, true, policies::policy<>());
+ return detail::ibeta_inv_ab_imp(b, x, static_cast<T>(1 - q), q, true, policies::policy<>());
 }
 
 } // namespace math

Modified: sandbox/interval_math_toolkit/boost/math/special_functions/detail/ibeta_inverse.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/special_functions/detail/ibeta_inverse.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/special_functions/detail/ibeta_inverse.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -66,7 +66,7 @@
    // get the first approximation for eta from the inverse
    // error function (Eq: 2.9 and 2.10).
    //
- T eta0 = boost::math::erfc_inv(2 * z, pol);
+ T eta0 = boost::math::erfc_inv(static_cast<T>(2 * z), pol);
    eta0 /= -sqrt(a / 2);
 
    T terms[4] = { eta0 };
@@ -106,7 +106,7 @@
    //
    // Bring them together to get a final estimate for eta:
    //
- T eta = tools::evaluate_polynomial(terms, 1/a, 4);
+ T eta = tools::evaluate_polynomial(terms, static_cast<T>(1/a), 4);
    //
    // now we need to convert eta to x, by solving the appropriate
    // quadratic equation:
@@ -143,7 +143,7 @@
    // Get first estimate for eta, see Eq 3.9 and 3.10,
    // but note there is a typo in Eq 3.10:
    //
- T eta0 = boost::math::erfc_inv(2 * z, pol);
+ T eta0 = boost::math::erfc_inv(static_cast<T>(2 * z), pol);
    eta0 /= -sqrt(r / 2);
 
    T s = sin(theta);
@@ -207,7 +207,7 @@
    // Bring the correction terms together to evaluate eta,
    // this is the last equation on page 151:
    //
- T eta = tools::evaluate_polynomial(terms, 1/r, 4);
+ T eta = tools::evaluate_polynomial(terms, static_cast<T>(1/r), 4);
    //
    // Now that we have eta we need to back solve for x,
    // we seek the value of x that gives eta in Eq 3.2.
@@ -453,6 +453,8 @@
 T ibeta_inv_imp(T a, T b, T p, T q, const Policy& pol, T* py)
 {
    BOOST_MATH_STD_USING // For ADL of math functions.
+ using std::max;
+ using std::min;
 
    //
    // The flag invert is set to true if we swap a for b and p for q,
@@ -521,8 +523,8 @@
          std::swap(p, q);
          invert = !invert;
       }
- T minv = (std::min)(a, b);
- T maxv = (std::max)(a, b);
+ T minv = min BOOST_PREVENT_MACRO_SUBSTITUTION(a, b);
+ T maxv = max BOOST_PREVENT_MACRO_SUBSTITUTION(a, b);
       if((sqrt(minv) > (maxv - minv)) && (minv > 5))
       {
          //
@@ -660,7 +662,7 @@
       //
       T lx = log(p * a * boost::math::beta(a, b, pol)) / a;
       x = exp(lx);
- y = x < 0.9 ? 1 - x : -boost::math::expm1(lx, pol);
+ y = x < 0.9 ? static_cast<T>(1 - x) : static_cast<T>(-boost::math::expm1(lx, pol));
 
       if((b < a) && (x < 0.2))
       {
@@ -734,7 +736,6 @@
          x = 1 - y;
       }
    }
-
    //
    // Now we have a guess for x (and for y) we can set things up for
    // iteration. If x > 0.5 it pays to swap things round:
@@ -804,7 +805,7 @@
    //BOOST_ASSERT(x != upper);
    //BOOST_ASSERT((x != lower) || (x == boost::math::tools::min_value<T>()) || (x == boost::math::tools::epsilon<T>()));
    //
- // Tidy up, if we "lower" was too high then zero is the best answer we have:
+ // Tidy up, if "lower" was too high then zero is the best answer we have:
    //
    if(x == lower)
       x = 0;

Modified: sandbox/interval_math_toolkit/boost/math/special_functions/detail/igamma_inverse.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/special_functions/detail/igamma_inverse.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/special_functions/detail/igamma_inverse.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -88,7 +88,7 @@
    // See equation 34.
    //
    BOOST_MATH_STD_USING
- T u = log(p) + boost::math::lgamma(a + 1, pol);
+ T u = log(p) + boost::math::lgamma(static_cast<T>(a + 1), pol);
    return exp((u + x - log(didonato_SN(a, x, N, tolerance))) / a);
 }
 
@@ -105,6 +105,7 @@
    // December 1986, Pages 377-393.
    //
    BOOST_MATH_STD_USING
+ using std::max;
 
    T result;
 
@@ -209,7 +210,7 @@
          }
          else
          {
- T D = (std::max)(T(2), a * (a - 1));
+ T D = max BOOST_PREVENT_MACRO_SUBSTITUTION(T(2), static_cast<T>(a * (a - 1)));
             T lg = boost::math::lgamma(a, pol);
             T lb = log(q) + lg;
             if(lb < -D * 2.3)
@@ -261,7 +262,7 @@
          {
             // DiDonato and Morris Eq 36:
             T zb = didonato_FN(p, a, z, 100, T(1e-4), pol);
- T u = log(p) + boost::math::lgamma(a + 1, pol);
+ T u = log(p) + boost::math::lgamma(static_cast<T>(a + 1), pol);
             result = zb * (1 - (a * log(zb) - zb - u + log(didonato_SN(a, z, 100, T(1e-4)))) / (a - zb));
          }
       }
@@ -357,7 +358,7 @@
       return tools::max_value<T>();
    if(p == 0)
       return 0;
- T guess = detail::find_inverse_gamma(a, p, 1 - p, pol);
+ T guess = detail::find_inverse_gamma(a, p, static_cast<T>(1 - p), pol);
    T lower = tools::min_value<T>();
    if(guess <= lower)
       guess = tools::min_value<T>();
@@ -399,7 +400,7 @@
       return tools::max_value<T>();
    if(q == 1)
       return 0;
- T guess = detail::find_inverse_gamma(a, 1 - q, q, pol);
+ T guess = detail::find_inverse_gamma(a, static_cast<T>(1 - q), q, pol);
    T lower = tools::min_value<T>();
    if(guess <= lower)
       guess = tools::min_value<T>();

Modified: sandbox/interval_math_toolkit/boost/math/special_functions/detail/lgamma_small.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/special_functions/detail/lgamma_small.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/special_functions/detail/lgamma_small.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -492,7 +492,7 @@
       // special case near 2:
       T dz = zm2;
       result = dz * log((z + L::g() - T(0.5)) / boost::math::constants::e<T>());
- result += boost::math::log1p(dz / (L::g() + T(1.5)), pol) * T(1.5);
+ result += boost::math::log1p(static_cast<T>(dz / (L::g() + T(1.5))), pol) * T(1.5);
       result += boost::math::log1p(L::lanczos_sum_near_2(dz), pol);
    }
    else
@@ -500,7 +500,7 @@
       // special case near 1:
       T dz = zm1;
       result = dz * log((z + L::g() - T(0.5)) / boost::math::constants::e<T>());
- result += boost::math::log1p(dz / (L::g() + T(0.5)), pol) / 2;
+ result += boost::math::log1p(static_cast<T>(dz / (L::g() + T(0.5))), pol) / 2;
       result += boost::math::log1p(L::lanczos_sum_near_1(dz), pol);
    }
    return result;

Modified: sandbox/interval_math_toolkit/boost/math/special_functions/detail/t_distribution_inv.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/special_functions/detail/t_distribution_inv.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/special_functions/detail/t_distribution_inv.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -30,7 +30,7 @@
    T a, b, c, d, q, x, y;
 
    if (ndf > 1e20f)
- return -boost::math::erfc_inv(2 * u, pol) * constants::root_two<T>();
+ return -boost::math::erfc_inv(static_cast<T>(2 * u), pol) * constants::root_two<T>();
 
    a = 1 / (ndf - 0.5f);
    b = 48 / (a * a);
@@ -43,14 +43,14 @@
       //
       // Asymptotic inverse expansion about normal:
       //
- x = -boost::math::erfc_inv(2 * u, pol) * constants::root_two<T>();
+ x = -boost::math::erfc_inv(static_cast<T>(2 * u), pol) * constants::root_two<T>();
       y = x * x;
 
       if (ndf < 5)
          c += 0.3f * (ndf - 4.5f) * (x + 0.6f);
       c += (((0.05f * d * x - 5) * x - 7) * x - 2) * x + b;
       y = (((((0.4f * y + 6.3f) * y + 36) * y + 94.5f) / c - y - 3) / b + 1) * x;
- y = boost::math::expm1(a * y * y, pol);
+ y = boost::math::expm1(static_cast<T>(a * y * y), pol);
    }
    else
    {
@@ -141,35 +141,60 @@
    // Figure out what the coefficients are, note these depend
    // only on the degrees of freedom (Eq 57 of Shaw):
    //
- c[2] = T(1) / 6 + T(1) / (6 * df);
+ c[2] = 0.16666666666666666667 + 0.16666666666666666667 / df;
    T in = 1 / df;
- c[3] = (((T(1) / 120) * in) + (T(1) / 15)) * in + (T(7) / 120);
- c[4] = ((((T(1) / 5040) * in + (T(1) / 560)) * in + (T(3) / 112)) * in + T(127) / 5040);
- c[5] = ((((T(1) / 362880) * in + (T(17) / 45360)) * in + (T(67) / 60480)) * in + (T(479) / 45360)) * in + (T(4369) / 362880);
- c[6] = ((((((T(1) / 39916800) * in + (T(2503) / 39916800)) * in + (T(11867) / 19958400)) * in + (T(1285) / 798336)) * in + (T(153161) / 39916800)) * in + (T(34807) / 5702400));
- c[7] = (((((((T(1) / 6227020800LL) * in + (T(37) / 2402400)) * in +
- (T(339929) / 2075673600LL)) * in + (T(67217) / 97297200)) * in +
- (T(870341) / 691891200LL)) * in + (T(70691) / 64864800LL)) * in +
- (T(20036983LL) / 6227020800LL));
- c[8] = (((((((T(1) / 1307674368000LL) * in + (T(1042243LL) / 261534873600LL)) * in +
- (T(21470159) / 435891456000LL)) * in + (T(326228899LL) / 1307674368000LL)) * in +
- (T(843620579) / 1307674368000LL)) * in + (T(332346031LL) / 435891456000LL)) * in +
- (T(43847599) / 1307674368000LL)) * in + (T(2280356863LL) / 1307674368000LL);
- c[9] = (((((((((T(1) / 355687428096000LL)) * in + (T(24262727LL) / 22230464256000LL)) * in +
- (T(123706507LL) / 8083805184000LL)) * in + (T(404003599LL) / 4446092851200LL)) * in +
- (T(51811946317LL) / 177843714048000LL)) * in + (T(91423417LL) / 177843714048LL)) * in +
- (T(32285445833LL) / 88921857024000LL)) * in + (T(531839683LL) / 1710035712000LL)) * in +
- (T(49020204823LL) / 50812489728000LL);
- c[10] = (((((((((T(1) / 121645100408832000LL) * in +
- (T(4222378423LL) / 13516122267648000LL)) * in +
- (T(49573465457LL) / 10137091700736000LL)) * in +
- (T(176126809LL) / 5304600576000LL)) * in +
- (T(44978231873LL) / 355687428096000LL)) * in +
- (T(5816850595639LL) / 20274183401472000LL)) * in +
- (T(73989712601LL) / 206879422464000LL)) * in +
- (T(26591354017LL) / 259925428224000LL)) * in +
- (T(14979648446341LL) / 40548366802944000LL)) * in +
- (T(65967241200001LL) / 121645100408832000LL);
+ c[3] = (0.0083333333333333333333 * in
+ + 0.066666666666666666667) * in
+ + 0.058333333333333333333;
+ c[4] = ((0.00019841269841269841270 * in
+ + 0.0017857142857142857143) * in
+ + 0.026785714285714285714) * in
+ + 0.025198412698412698413;
+ c[5] = (((2.7557319223985890653e10-6 * in
+ + 0.00037477954144620811287) * in
+ - 0.0011078042328042328042) * in
+ + 0.010559964726631393298) * in
+ + 0.012039792768959435626;
+ c[6] = ((((2.5052108385441718775e-8 * in
+ - 0.000062705427288760622094) * in
+ + 0.00059458674042007375341) * in
+ - 0.0016095979637646304313) * in
+ + 0.0061039211560044893378) * in
+ + 0.0038370059724226390893;
+ c[7] = (((((1.6059043836821614599e-10 * in
+ + 0.000015401265401265401265) * in
+ - 0.00016376804137220803887) * in
+ + 0.00069084207973096861986) * in
+ - 0.0012579159844784844785) * in
+ + 0.0010898206731540064873) * in
+ + 0.0032177478835464946576;
+ c[8] = ((((((7.6471637318198164759e-13 * in
+ - 3.9851014346715404916e-6) * in
+ + 0.000049255746366361445727) * in
+ - 0.00024947258047043099953) * in
+ + 0.00064513046951456342991) * in
+ - 0.00076245135440323932387) * in
+ + 0.000033530976880017885309) * in
+ + 0.0017438262298340009980;
+ c[9] = (((((((2.8114572543455207632e-15 * in
+ + 1.0914179173496789432e-6) * in
+ - 0.000015303004486655377567) * in
+ + 0.000090867107935219902229) * in
+ - 0.00029133414466938067350) * in
+ + 0.00051406605788341121363) * in
+ - 0.00036307660358786885787) * in
+ - 0.00031101086326318780412) * in
+ + 0.00096472747321388644237;
+ c[10] = ((((((((8.2206352466243297170e-18 * in
+ - 3.1239569599829868045e-7) * in
+ + 4.8903045291975346210e-6) * in
+ - 0.000033202652391372058698) * in
+ + 0.00012645437628698076975) * in
+ - 0.00028690924218514613987) * in
+ + 0.00035764655430568632777) * in
+ - 0.00010230378073700412687) * in
+ - 0.00036942667800009661203) * in
+ + 0.00054229262813129686486;
    //
    // The result is then a polynomial in v (see Eq 56 of Shaw):
    //
@@ -258,7 +283,7 @@
             //
             T a = 4 * (u - u * u);//1 - 4 * (u - 0.5f) * (u - 0.5f);
             T b = boost::math::cbrt(a);
- static const T c = 0.85498797333834849467655443627193L;
+ static const T c = 0.85498797333834849467655443627193;
             T p = 6 * (1 + c * (1 / b - 1));
             T p0;
             do{
@@ -368,7 +393,7 @@
          // where we use Shaw's tail series.
          // The crossover point is roughly exponential in -df:
          //
- T crossover = ldexp(1.0f, tools::real_cast<int>(df / -0.654f));
+ T crossover = ldexp(1.0f, tools::real_cast<int>(static_cast<T>(df / -0.654f)));
          if(u > crossover)
          {
             result = boost::math::detail::inverse_students_t_hill(df, u, pol);
@@ -385,7 +410,7 @@
 template <class T, class Policy>
 inline T find_ibeta_inv_from_t_dist(T a, T p, T q, T* py, const Policy& pol)
 {
- T u = (p > q) ? 0.5f - q / 2 : p / 2;
+ T u = (p > q) ? static_cast<T>(0.5f - q / 2) : static_cast<T>(p / 2);
    T v = 1 - u; // u < 0.5 so no cancellation error
    T df = a * 2;
    T t = boost::math::detail::inverse_students_t(df, u, v, pol);

Modified: sandbox/interval_math_toolkit/boost/math/special_functions/detail/unchecked_factorial.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/special_functions/detail/unchecked_factorial.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/special_functions/detail/unchecked_factorial.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -279,8 +279,16 @@
       value = ::boost::math::max_factorial<long double>::value);
 };
 
+namespace detail{
+
 template <class T>
-inline T unchecked_factorial(unsigned i BOOST_MATH_APPEND_EXPLICIT_TEMPLATE_TYPE_SPEC(T))
+inline T unchecked_factorial_imp(const mpl::true_&, unsigned i BOOST_MATH_APPEND_EXPLICIT_TEMPLATE_TYPE_SPEC(T))
+{
+ return boost::math::unchecked_factorial<typename T::base_type>(i);
+}
+
+template <class T>
+inline T unchecked_factorial_imp(const mpl::false_&, unsigned i BOOST_MATH_APPEND_EXPLICIT_TEMPLATE_TYPE_SPEC(T))
 {
    static const boost::array<T, 101> factorials = {{
       boost::lexical_cast<T>("1"),
@@ -389,6 +397,15 @@
    return factorials[i];
 }
 
+} // namespace detail
+
+template <class T>
+inline T unchecked_factorial(unsigned i BOOST_MATH_APPEND_EXPLICIT_TEMPLATE_TYPE_SPEC(T))
+{
+ typedef typename boost::numeric::is_interval<T>::type tag_type;
+ return detail::unchecked_factorial_imp<T>(tag_type(), i);
+}
+
 template <class T>
 struct max_factorial
 {

Modified: sandbox/interval_math_toolkit/boost/math/special_functions/ellint_1.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/special_functions/ellint_1.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/special_functions/ellint_1.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -88,9 +88,9 @@
        // so rewritten to use fmod instead:
        //
        BOOST_MATH_INSTRUMENT_CODE("pi/2 = " << constants::pi<T>() / 2);
- T rphi = boost::math::tools::fmod_workaround(phi, constants::pi<T>() / 2);
+ T rphi = boost::math::tools::fmod_workaround(phi, static_cast<T>(constants::pi<T>() / 2));
        BOOST_MATH_INSTRUMENT_VARIABLE(rphi);
- T m = 2 * (phi - rphi) / constants::pi<T>();
+ T m = floor(2 * (phi - rphi) / constants::pi<T>() + 0.5f);
        BOOST_MATH_INSTRUMENT_VARIABLE(m);
        int s = 1;
        if(boost::math::tools::fmod_workaround(m, T(2)) > 0.5)
@@ -104,7 +104,7 @@
        T cosp = cos(rphi);
        BOOST_MATH_INSTRUMENT_VARIABLE(sinp);
        BOOST_MATH_INSTRUMENT_VARIABLE(cosp);
- result = s * sinp * ellint_rf_imp(cosp * cosp, 1 - k * k * sinp * sinp, T(1), pol);
+ result = s * sinp * ellint_rf_imp(static_cast<T>(cosp * cosp), static_cast<T>(1 - k * k * sinp * sinp), T(1), pol);
        BOOST_MATH_INSTRUMENT_VARIABLE(result);
        if(m != 0)
        {

Modified: sandbox/interval_math_toolkit/boost/math/special_functions/ellint_2.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/special_functions/ellint_2.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/special_functions/ellint_2.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -74,8 +74,8 @@
        // but that fails if T has more digits than a long long,
        // so rewritten to use fmod instead:
        //
- T rphi = boost::math::tools::fmod_workaround(phi, constants::pi<T>() / 2);
- T m = 2 * (phi - rphi) / constants::pi<T>();
+ T rphi = boost::math::tools::fmod_workaround(phi, static_cast<T>(constants::pi<T>() / 2));
+ T m = floor(2 * (phi - rphi) / constants::pi<T>() + 0.5f);
        int s = 1;
        if(boost::math::tools::fmod_workaround(m, T(2)) > 0.5)
        {

Modified: sandbox/interval_math_toolkit/boost/math/special_functions/ellint_3.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/special_functions/ellint_3.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/special_functions/ellint_3.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -61,7 +61,7 @@
 
     T sphi = sin(fabs(phi));
 
- if(v > 1 / (sphi * sphi))
+ if(tools::maybe_greater(v, static_cast<T>(1 / (sphi * sphi))))
     {
         // Complex result is a domain error:
        return policies::raise_domain_error<T>(function,
@@ -102,7 +102,7 @@
           // v > 1:
           T vcr = sqrt(-vc);
           T arg = vcr * tan(phi);
- return (boost::math::log1p(arg, pol) - boost::math::log1p(-arg, pol)) / (2 * vcr);
+ return (boost::math::log1p(arg, pol) - boost::math::log1p(static_cast<T>(-arg), pol)) / (2 * vcr);
        }
     }
 
@@ -185,8 +185,8 @@
     }
     else
     {
- T rphi = boost::math::tools::fmod_workaround(fabs(phi), constants::pi<T>() / 2);
- T m = 2 * (fabs(phi) - rphi) / constants::pi<T>();
+ T rphi = boost::math::tools::fmod_workaround(fabs(phi), static_cast<T>(constants::pi<T>() / 2));
+ T m = floor(2 * (fabs(phi) - rphi) / constants::pi<T>() + 0.5f);
        int sign = 1;
        if(boost::math::tools::fmod_workaround(m, T(2)) > 0.5)
        {

Modified: sandbox/interval_math_toolkit/boost/math/special_functions/ellint_rj.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/special_functions/ellint_rj.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/special_functions/ellint_rj.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -94,7 +94,7 @@
        value = boost::math::ellint_rj(x, y, z, p, pol);
        value *= pmy;
        value -= 3 * boost::math::ellint_rf(x, y, z, pol);
- value += 3 * sqrt((x * y * z) / (x * z + p * q)) * boost::math::ellint_rc(x * z + p * q, p * q, pol);
+ value += 3 * sqrt((x * y * z) / (x * z + p * q)) * boost::math::ellint_rc(static_cast<T>(x * z + p * q), static_cast<T>(p * q), pol);
        value /= (y + q);
        return value;
     }

Modified: sandbox/interval_math_toolkit/boost/math/special_functions/erf.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/special_functions/erf.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/special_functions/erf.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -122,9 +122,9 @@
    if(z < 0)
    {
       if(!invert)
- return -erf_imp(-z, invert, pol, t);
+ return -erf_imp(static_cast<T>(-z), invert, pol, t);
       else
- return 1 + erf_imp(-z, false, pol, t);
+ return 1 + erf_imp(static_cast<T>(-z), false, pol, t);
    }
 
    T result;
@@ -142,7 +142,7 @@
       if(x < 0.6)
       {
          // Compute P:
- result = z * exp(-x);
+ result = z * exp(static_cast<T>(-x));
          result /= sqrt(boost::math::constants::pi<T>());
          if(result != 0)
             result *= 2 * detail::lower_gamma_series(T(0.5f), x, pol);

Modified: sandbox/interval_math_toolkit/boost/math/special_functions/expint.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/special_functions/expint.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/special_functions/expint.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -436,7 +436,7 @@
    if(z < 0)
       return policies::raise_domain_error<T>(function, "Function requires z >= 0 but got %1%.", z, pol);
    if(z == 0)
- return n == 1 ? policies::raise_overflow_error<T>(function, 0, pol) : 1 / (static_cast<T>(n - 1));
+ return n == 1 ? policies::raise_overflow_error<T>(function, 0, pol) : static_cast<T>(1 / (static_cast<T>(n - 1)));
 
    T result;
 
@@ -504,7 +504,7 @@
 {
    static const char* function = "boost::math::expint<%1%>(%1%)";
    if(z < 0)
- return -expint_imp(1, -z, pol, tag);
+ return -expint_imp(1, static_cast<T>(-z), pol, tag);
    if(z == 0)
       return -policies::raise_overflow_error<T>(function, 0, pol);
    return expint_i_as_series(z, pol);

Modified: sandbox/interval_math_toolkit/boost/math/special_functions/expm1.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/special_functions/expm1.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/special_functions/expm1.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -75,7 +75,7 @@
    BOOST_MATH_STD_USING
 
    T a = fabs(x);
- if(a > T(0.5L))
+ if(a > T(0.5f))
       return exp(x) - T(1);
    if(a < tools::epsilon<T>())
       return x;

Modified: sandbox/interval_math_toolkit/boost/math/special_functions/factorials.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/special_functions/factorials.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/special_functions/factorials.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -79,7 +79,7 @@
       // Fallthrough: i is too large to use table lookup, try the
       // gamma function instead.
       //
- T result = boost::math::tgamma(static_cast<T>(i) / 2 + 1, pol) / sqrt(constants::pi<T>());
+ T result = boost::math::tgamma(static_cast<T>(static_cast<T>(i) / 2 + 1), pol) / sqrt(constants::pi<T>());
       if(ldexp(tools::max_value<T>(), -static_cast<int>(i+1) / 2) > result)
          return ceil(result * ldexp(T(1), (i+1) / 2) - 0.5f);
    }
@@ -125,7 +125,7 @@
          n = -n;
          inv = true;
       }
- T result = ((n&1) ? -1 : 1) * falling_factorial(-x, n, pol);
+ T result = ((n&1) ? -1 : 1) * falling_factorial(static_cast<T>(-x), n, pol);
       if(inv)
          result = 1 / result;
       return result;
@@ -152,7 +152,7 @@
       // For x < 0 we really have a rising factorial
       // modulo a possible change of sign:
       //
- return (n&1 ? -1 : 1) * rising_factorial(-x, n, pol);
+ return (n&1 ? -1 : 1) * rising_factorial(static_cast<T>(-x), n, pol);
    }
    if(n == 0)
       return 1;
@@ -163,15 +163,15 @@
       // handle it, split the product up into three parts:
       //
       T xp1 = x + 1;
- unsigned n2 = tools::real_cast<unsigned>(floor(xp1));
+ unsigned n2 = tools::real_cast<unsigned>(static_cast<T>(floor(xp1)));
       if(n2 == xp1)
          return 0;
- T result = boost::math::tgamma_delta_ratio(xp1, -static_cast<T>(n2), pol);
+ T result = boost::math::tgamma_delta_ratio(xp1, static_cast<T>(-static_cast<T>(n2)), pol);
       x -= n2;
       result *= x;
       ++n2;
       if(n2 < n)
- result *= falling_factorial(x - 1, n - n2, pol);
+ result *= falling_factorial(static_cast<T>(x - 1), n - n2, pol);
       return result;
    }
    //
@@ -181,7 +181,7 @@
    // because tgamma_delta_ratio is alreay optimised
    // for that use case:
    //
- return boost::math::tgamma_delta_ratio(x + 1, -static_cast<T>(n), pol);
+ return boost::math::tgamma_delta_ratio(static_cast<T>(x + 1), static_cast<T>(-static_cast<T>(n)), pol);
 }
 
 } // namespace detail

Modified: sandbox/interval_math_toolkit/boost/math/special_functions/gamma.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/special_functions/gamma.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/special_functions/gamma.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -146,7 +146,7 @@
       return policies::raise_pole_error<T>(function, "Evaluation of tgamma at a negative integer %1%.", z, pol);
    if(z <= -20)
    {
- result = gamma_imp(-z, pol, l) * sinpx(z);
+ result = gamma_imp(static_cast<T>(-z), pol, l) * sinpx(z);
       if((fabs(result) < 1) && (tools::max_value<T>() * fabs(result) < boost::math::constants::pi<T>()))
          return policies::raise_overflow_error<T>(function, "Result of tgamma is too large to represent.", pol);
       result = -boost::math::constants::pi<T>() / result;
@@ -233,13 +233,17 @@
    {
       typedef typename policies::precision<T, Policy>::type precision_type;
       typedef typename mpl::if_<
- mpl::less_equal<precision_type, mpl::int_<64> >,
- mpl::int_<64>,
- typename mpl::if_<
- mpl::less_equal<precision_type, mpl::int_<113> >,
- mpl::int_<113>, mpl::int_<0> >::type
- >::type tag_type;
- result = lgamma_small_imp(z, z - 1, z - 2, tag_type(), pol, l);
+ mpl::less_equal<precision_type, mpl::int_<0> >,
+ mpl::int_<0>,
+ mpl::if_<
+ mpl::less_equal<precision_type, mpl::int_<64> >,
+ mpl::int_<64>,
+ typename mpl::if_<
+ mpl::less_equal<precision_type, mpl::int_<113> >,
+ mpl::int_<113>, mpl::int_<0> >::type
+ >::type
+ >::type tag_type;
+ result = lgamma_small_imp(z, static_cast<T>(z - 1), static_cast<T>(z - 2), tag_type(), pol, l);
    }
    else if((z >= 3) && (z < 100))
    {
@@ -346,7 +350,7 @@
       return policies::raise_pole_error<T>(function, "Evaluation of tgamma at a negative integer %1%.", z, pol);
    if(z <= -20)
    {
- T result = gamma_imp(-z, pol, l) * sinpx(z);
+ T result = gamma_imp(static_cast<T>(-z), pol, l) * sinpx(z);
       if((fabs(result) < 1) && (tools::max_value<T>() * fabs(result) < boost::math::constants::pi<T>()))
          return policies::raise_overflow_error<T>(function, "Result of tgamma is too large to represent.", pol);
       result = -boost::math::constants::pi<T>() / result;
@@ -414,7 +418,8 @@
    }
    else if((z != 1) && (z != 2))
    {
- T limit = (std::max)(z+1, T(10));
+ using std::max;
+ T limit = max BOOST_PREVENT_MACRO_SUBSTITUTION(static_cast<T>(z+1), T(10));
       T prefix = z * log(limit) - limit;
       T sum = detail::lower_gamma_series(z, limit, pol) / z;
       sum += detail::upper_gamma_fraction(z, limit, ::boost::math::policies::digits<T, Policy>());
@@ -456,14 +461,14 @@
       if(dz < -0.5)
       {
          // Best method is simply to subtract 1 from tgamma:
- result = boost::math::tgamma(1+dz, pol) - 1;
+ result = boost::math::tgamma(static_cast<T>(1+dz), pol) - 1;
          BOOST_MATH_INSTRUMENT_CODE(result);
       }
       else
       {
          // Use expm1 on lgamma:
- result = boost::math::expm1(-boost::math::log1p(dz, pol)
- + lgamma_small_imp(dz+2, dz + 1, dz, tag_type(), pol, l));
+ result = boost::math::expm1(static_cast<T>(-boost::math::log1p(dz, pol)
+ + lgamma_small_imp(static_cast<T>(dz+2), static_cast<T>(dz + 1), dz, tag_type(), pol, l)));
          BOOST_MATH_INSTRUMENT_CODE(result);
       }
    }
@@ -472,13 +477,13 @@
       if(dz < 2)
       {
          // Use expm1 on lgamma:
- result = boost::math::expm1(lgamma_small_imp(dz+1, dz, dz-1, tag_type(), pol, l), pol);
+ result = boost::math::expm1(lgamma_small_imp(static_cast<T>(dz+1), dz, static_cast<T>(dz-1), tag_type(), pol, l), pol);
          BOOST_MATH_INSTRUMENT_CODE(result);
       }
       else
       {
          // Best method is simply to subtract 1 from tgamma:
- result = boost::math::tgamma(1+dz, pol) - 1;
+ result = boost::math::tgamma(static_cast<T>(1+dz), pol) - 1;
          BOOST_MATH_INSTRUMENT_CODE(result);
       }
    }
@@ -496,7 +501,7 @@
    // algebra isn't easy for the general case....
    // Start by subracting 1 from tgamma:
    //
- T result = gamma_imp(1 + dz, pol, l) - 1;
+ T result = gamma_imp(static_cast<T>(1 + dz), pol, l) - 1;
    BOOST_MATH_INSTRUMENT_CODE(result);
    //
    // Test the level of cancellation error observed: we loose one bit
@@ -597,6 +602,9 @@
 T regularised_gamma_prefix(T a, T z, const Policy& pol, const L& l)
 {
    BOOST_MATH_STD_USING
+ using std::max;
+ using std::min;
+
    T agh = a + static_cast<T>(L::g()) - T(0.5);
    T prefix;
    T d = ((z - a) - static_cast<T>(L::g()) + T(0.5)) / agh;
@@ -639,16 +647,16 @@
       //
       T alz = a * log(z / agh);
       T amz = a - z;
- if(((std::min)(alz, amz) <= tools::log_min_value<T>()) || ((std::max)(alz, amz) >= tools::log_max_value<T>()))
+ if((min BOOST_PREVENT_MACRO_SUBSTITUTION(alz, amz) <= tools::log_min_value<T>()) || (max BOOST_PREVENT_MACRO_SUBSTITUTION(alz, amz) >= tools::log_max_value<T>()))
       {
          T amza = amz / a;
- if(((std::min)(alz, amz)/2 > tools::log_min_value<T>()) && ((std::max)(alz, amz)/2 < tools::log_max_value<T>()))
+ if((min BOOST_PREVENT_MACRO_SUBSTITUTION(alz, amz)/2 > tools::log_min_value<T>()) && (max BOOST_PREVENT_MACRO_SUBSTITUTION(alz, amz)/2 < tools::log_max_value<T>()))
          {
             // compute square root of the result and then square it:
             T sq = pow(z / agh, a / 2) * exp(amz / 2);
             prefix = sq * sq;
          }
- else if(((std::min)(alz, amz)/4 > tools::log_min_value<T>()) && ((std::max)(alz, amz)/4 < tools::log_max_value<T>()) && (z > a))
+ else if((min BOOST_PREVENT_MACRO_SUBSTITUTION(alz, amz)/4 > tools::log_min_value<T>()) && (max BOOST_PREVENT_MACRO_SUBSTITUTION(alz, amz)/4 < tools::log_max_value<T>()) && (z > a))
          {
             // compute the 4th root of the result then square it twice:
             T sq = pow(z / agh, a / 4) * exp(amz / 4);
@@ -679,8 +687,10 @@
 T regularised_gamma_prefix(T a, T z, const Policy& pol, const lanczos::undefined_lanczos&)
 {
    BOOST_MATH_STD_USING
+ using std::max;
+ using std::min;
 
- T limit = (std::max)(T(10), a);
+ T limit = max BOOST_PREVENT_MACRO_SUBSTITUTION(T(10), a);
    T sum = detail::lower_gamma_series(a, limit, pol) / a;
    sum += detail::upper_gamma_fraction(a, limit, ::boost::math::policies::digits<T, Policy>());
 
@@ -701,7 +711,7 @@
    T amz = a - z;
    T alzoa = a * log(zoa);
    T prefix;
- if(((std::min)(alzoa, amz) <= tools::log_min_value<T>()) || ((std::max)(alzoa, amz) >= tools::log_max_value<T>()))
+ if((min BOOST_PREVENT_MACRO_SUBSTITUTION(alzoa, amz) <= tools::log_min_value<T>()) || (max BOOST_PREVENT_MACRO_SUBSTITUTION(alzoa, amz) >= tools::log_max_value<T>()))
    {
       T amza = amz / a;
       if((amza <= tools::log_min_value<T>()) || (amza >= tools::log_max_value<T>()))
@@ -777,7 +787,7 @@
    // Calculates normalised Q when a is a half-integer:
    //
    BOOST_MATH_STD_USING
- T e = boost::math::erfc(sqrt(x), pol);
+ T e = boost::math::erfc(static_cast<T>(sqrt(x)), pol);
    if((e != 0) && (a > 1))
    {
       T term = exp(-x) / sqrt(constants::pi<T>() * x);
@@ -1025,7 +1035,7 @@
    T result;
    if(fabs(delta) < 10)
    {
- result = exp((constants::half<T>() - z) * boost::math::log1p(delta / zgh, pol));
+ result = exp((constants::half<T>() - z) * boost::math::log1p(static_cast<T>(delta / zgh), pol));
    }
    else
    {
@@ -1095,7 +1105,7 @@
          //
          if((z <= max_factorial<T>::value) && (z + delta <= max_factorial<T>::value))
          {
- return unchecked_factorial<T>(tools::real_cast<unsigned>(z) - 1) / unchecked_factorial<T>(tools::real_cast<unsigned>(z + delta) - 1);
+ return unchecked_factorial<T>(tools::real_cast<unsigned>(z) - 1) / unchecked_factorial<T>(tools::real_cast<unsigned>(static_cast<T>(z + delta)) - 1);
          }
       }
       if(fabs(delta) < 20)
@@ -1426,7 +1436,7 @@
       policies::discrete_quantile<>,
       policies::assert_undefined<> >::type forwarding_policy;
 
- return policies::checked_narrowing_cast<result_type, forwarding_policy>(detail::tgamma_delta_ratio_imp(static_cast<value_type>(a), static_cast<value_type>(b) - static_cast<value_type>(a), forwarding_policy()), "boost::math::tgamma_delta_ratio<%1%>(%1%, %1%)");
+ return policies::checked_narrowing_cast<result_type, forwarding_policy>(detail::tgamma_delta_ratio_imp(static_cast<value_type>(a), static_cast<value_type>(static_cast<value_type>(b) - static_cast<value_type>(a)), forwarding_policy()), "boost::math::tgamma_delta_ratio<%1%>(%1%, %1%)");
 }
 template <class T1, class T2>
 inline typename tools::promote_args<T1, T2>::type

Modified: sandbox/interval_math_toolkit/boost/math/special_functions/legendre.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/special_functions/legendre.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/special_functions/legendre.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -45,7 +45,7 @@
    if(second)
    {
       // A solution of the second kind (Q):
- p0 = (boost::math::log1p(x, pol) - boost::math::log1p(-x, pol)) / 2;
+ p0 = (boost::math::log1p(x, pol) - boost::math::log1p(static_cast<T>(-x), pol)) / 2;
       p1 = x * p0 - 1;
    }
    else

Modified: sandbox/interval_math_toolkit/boost/math/special_functions/log1p.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/special_functions/log1p.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/special_functions/log1p.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -88,7 +88,7 @@
          function, 0, pol);
 
    result_type a = abs(result_type(x));
- if(a > result_type(0.5L))
+ if(a > result_type(0.5f))
       return log(1 + result_type(x));
    // Note that without numeric_limits specialisation support,
    // epsilon just returns zero, and our "optimisation" will always fail:
@@ -248,7 +248,7 @@
          function, 0, pol);
 
    result_type a = abs(result_type(x));
- if(a > result_type(0.95L))
+ if(a > result_type(0.95f))
       return log(1 + result_type(x)) - result_type(x);
    // Note that without numeric_limits specialisation support,
    // epsilon just returns zero, and our "optimisation" will always fail:

Modified: sandbox/interval_math_toolkit/boost/math/special_functions/sin_pi.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/special_functions/sin_pi.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/special_functions/sin_pi.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -22,7 +22,7 @@
 {
    BOOST_MATH_STD_USING // ADL of std names
    if(x < 0)
- return -sin_pi(-x);
+ return -sin_pi(static_cast<T>(-x));
    // sin of pi*x:
    bool invert;
    if(x < 0.5)

Modified: sandbox/interval_math_toolkit/boost/math/special_functions/spherical_harmonic.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/special_functions/spherical_harmonic.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/special_functions/spherical_harmonic.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -61,7 +61,7 @@
    if(m&1)
    {
       // Check phase if theta is outside [0, PI]:
- T mod = boost::math::tools::fmod_workaround(theta, 2 * constants::pi<T>());
+ T mod = boost::math::tools::fmod_workaround(theta, static_cast<T>(2 * constants::pi<T>()));
       if(mod < 0)
          mod += 2 * constants::pi<T>();
       if(mod > constants::pi<T>())
@@ -88,7 +88,7 @@
    if(m&1)
    {
       // Check phase if theta is outside [0, PI]:
- T mod = boost::math::tools::fmod_workaround(theta, 2 * constants::pi<T>());
+ T mod = boost::math::tools::fmod_workaround(theta, static_cast<T>(2 * constants::pi<T>()));
       if(mod < 0)
          mod += 2 * constants::pi<T>();
       if(mod > constants::pi<T>())

Modified: sandbox/interval_math_toolkit/boost/math/special_functions/sqrt1pm1.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/special_functions/sqrt1pm1.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/special_functions/sqrt1pm1.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -28,7 +28,7 @@
 
    if(fabs(result_type(val)) > 0.75)
       return sqrt(1 + result_type(val)) - 1;
- return boost::math::expm1(boost::math::log1p(val, pol) / 2, pol);
+ return boost::math::expm1(static_cast<T>(boost::math::log1p(val, pol) / 2), pol);
 }
 
 template <class T>

Modified: sandbox/interval_math_toolkit/boost/math/special_functions/zeta.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/special_functions/zeta.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/special_functions/zeta.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -836,7 +836,7 @@
          result = 0;
       else
       {
- result = boost::math::sin_pi(0.5f * sc, pol)
+ result = boost::math::sin_pi(static_cast<T>(0.5f * sc), pol)
             * 2 * pow(2 * constants::pi<T>(), -s)
             * boost::math::tgamma(s, pol)
             * zeta_imp(s, sc, pol, tag);
@@ -884,7 +884,7 @@
 
    return policies::checked_narrowing_cast<result_type, forwarding_policy>(detail::zeta_imp(
       static_cast<value_type>(s),
- 1 - static_cast<value_type>(s),
+ static_cast<value_type>(1 - static_cast<value_type>(s)),
       forwarding_policy(),
       tag_type()), "boost::math::zeta<%1%>(%1%)");
 }

Modified: sandbox/interval_math_toolkit/boost/math/tools/config.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/tools/config.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/tools/config.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -12,7 +12,9 @@
 
 #include <boost/cstdint.hpp> // for boost::uintmax_t
 #include <boost/config.hpp>
+#include <boost/mpl/bool.hpp>
 #include <boost/detail/workaround.hpp>
+#include <boost/tr1/tuple.hpp>
 #include <algorithm> // for min and max
 #include <cmath>
 #include <climits>
@@ -159,6 +161,63 @@
 # define BOOST_MATH_CONTROL_FP
 #endif
 //
+// Forward declare interval type so we can handle it as a special
+// case where required.
+//
+namespace boost{ namespace numeric{
+
+template <class T, class Policy>
+class interval;
+
+template <class T>
+struct is_interval : public boost::mpl::false_{};
+
+template <class T, class Policy>
+struct is_interval<interval<T, Policy> > : public boost::mpl::true_{};
+
+template <class T>
+struct interval_base_type
+{
+ typedef T type;
+};
+
+template <class T, class Policy>
+struct interval_base_type<interval<T, Policy> >
+{
+ typedef T type;
+};
+
+// Namespace forward declarations:
+namespace interval_math_compare{}
+namespace interval_lib{ namespace compare{ namespace possible{} } }
+
+template <class T, class Policy>
+inline std::tr1::tuple<T,T> make_printable(const boost::numeric::interval<T, Policy>& v)
+{
+ return std::make_pair(v.lower(), v.upper());
+}
+
+template <class T>
+inline const T& make_printable(const T& v)
+{
+ return v;
+}
+
+#ifdef BOOST_MATH_INSTRUMENT
+//
+// We need a stream operator for interval types, this may well conflict
+// with a user provided definition:
+//
+template <class T, class Policy>
+std::ostream& operator << (std::ostream& os, const interval<T, Policy>& val)
+{
+ os << "{ " << val.lower() << ", " << val.upper() << " }";
+ return os;
+}
+#endif
+
+}} // namespaces
+//
 // Helper macro for using statements:
 //
 #define BOOST_MATH_STD_USING \
@@ -184,7 +243,8 @@
    using std::ceil;\
    using std::floor;\
    using std::log10;\
- using std::sqrt;
+ using std::sqrt;\
+ using namespace boost::numeric::interval_math_compare;
 
 
 namespace boost{ namespace math{
@@ -194,14 +254,101 @@
 template <class T>
 inline T max BOOST_PREVENT_MACRO_SUBSTITUTION(T a, T b, T c)
 {
- return (std::max)((std::max)(a, b), c);
+ using std::max;
+ return max BOOST_PREVENT_MACRO_SUBSTITUTION(max BOOST_PREVENT_MACRO_SUBSTITUTION(a, b), c);
 }
 
 template <class T>
 inline T max BOOST_PREVENT_MACRO_SUBSTITUTION(T a, T b, T c, T d)
 {
- return (std::max)((std::max)(a, b), (std::max)(c, d));
+ using std::max;
+ return max BOOST_PREVENT_MACRO_SUBSTITUTION(max BOOST_PREVENT_MACRO_SUBSTITUTION(a, b), max BOOST_PREVENT_MACRO_SUBSTITUTION(c, d));
+}
+//
+// Interval aware "soft" comparisons:
+//
+template <class T>
+inline bool maybe_equal(const T& a, const T& b)
+{
+ return a == b;
+}
+template <class T>
+inline bool maybe_less(const T& a, const T& b)
+{
+ return a < b;
+}
+template <class T>
+inline bool maybe_less_equal(const T& a, const T& b)
+{
+ return a <= b;
+}
+template <class T>
+inline bool maybe_greater(const T& a, const T& b)
+{
+ return a > b;
+}
+template <class T>
+inline bool maybe_greater_equal(const T& a, const T& b)
+{
+ return a >= b;
+}
+
+template <class T, class Policy>
+inline bool maybe_equal(const boost::numeric::interval<T, Policy>& a, const boost::numeric::interval<T, Policy>& b)
+{
+ using namespace boost::numeric::interval_lib::compare::possible;
+ return a == b;
+}
+template <class T, class Policy>
+inline bool maybe_less(const boost::numeric::interval<T, Policy>& a, const boost::numeric::interval<T, Policy>& b)
+{
+ using namespace boost::numeric::interval_lib::compare::possible;
+ return a < b;
+}
+template <class T, class Policy>
+inline bool maybe_less_equal(const boost::numeric::interval<T, Policy>& a, const boost::numeric::interval<T, Policy>& b)
+{
+ using namespace boost::numeric::interval_lib::compare::possible;
+ return a <= b;
 }
+template <class T, class Policy>
+inline bool maybe_greater(const boost::numeric::interval<T, Policy>& a, const boost::numeric::interval<T, Policy>& b)
+{
+ using namespace boost::numeric::interval_lib::compare::possible;
+ return a > b;
+}
+template <class T, class Policy>
+inline bool maybe_greater_equal(const boost::numeric::interval<T, Policy>& a, const boost::numeric::interval<T, Policy>& b)
+{
+ using namespace boost::numeric::interval_lib::compare::possible;
+ return a >= b;
+}
+
+//
+// Normalisation functions:
+//
+template <class T>
+inline void normalize_median(T&){}
+template <class T, class Policy>
+inline void normalize_median(boost::numeric::interval<T, Policy>& v)
+{
+ v = median(v);
+}
+template <class T>
+inline void normalize_up(T&){}
+template <class T, class Policy>
+inline void normalize_up(boost::numeric::interval<T, Policy>& v)
+{
+ v = v.upper();
+}
+template <class T>
+inline void normalize_down(T&){}
+template <class T, class Policy>
+inline void normalize_down(boost::numeric::interval<T, Policy>& v)
+{
+ v = v.lower();
+}
+
 } // namespace tools
 }} // namespace boost namespace math
 

Modified: sandbox/interval_math_toolkit/boost/math/tools/fraction.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/tools/fraction.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/tools/fraction.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -116,7 +116,7 @@
       D = 1/D;
       delta = C*D;
       f = f * delta;
- }while(fabs(delta - 1) > factor);
+ }while(!maybe_less_equal(fabs(delta - 1), factor));
 
    return f;
 }
@@ -155,7 +155,7 @@
       D = 1/D;
       delta = C*D;
       f = f * delta;
- }while((fabs(delta - 1) > factor) && --counter);
+ }while(!maybe_less_equal(fabs(delta - 1), factor) && --counter);
 
    max_terms = max_terms - counter;
 
@@ -209,7 +209,7 @@
       D = 1/D;
       delta = C*D;
       f = f * delta;
- }while(fabs(delta - 1) > factor);
+ }while(!maybe_less_equal(fabs(delta - 1), factor));
 
    return a0/f;
 }
@@ -249,7 +249,7 @@
       D = 1/D;
       delta = C*D;
       f = f * delta;
- }while((fabs(delta - 1) > factor) && --counter);
+ }while(!maybe_less_equal(fabs(delta - 1), factor) && --counter);
 
    max_terms = max_terms - counter;
 

Modified: sandbox/interval_math_toolkit/boost/math/tools/minima.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/tools/minima.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/tools/minima.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -23,7 +23,9 @@
 std::pair<T, T> brent_find_minima(F f, T min, T max, int bits, boost::uintmax_t& max_iter)
 {
    BOOST_MATH_STD_USING
- bits = (std::min)(policies::digits<T, policies::policy<> >() / 2, bits);
+ using std::min;
+
+ bits = min BOOST_PREVENT_MACRO_SUBSTITUTION(policies::digits<T, policies::policy<> >() / 2, bits);
    T tolerance = static_cast<T>(ldexp(1.0, 1-bits));
    T x; // minima so far
    T w; // second best point

Modified: sandbox/interval_math_toolkit/boost/math/tools/rational.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/tools/rational.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/tools/rational.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -98,19 +98,19 @@
 template <class T, class U>
 inline U evaluate_even_polynomial(const T* poly, U z, std::size_t count)
 {
- return evaluate_polynomial(poly, z*z, count);
+ return evaluate_polynomial(poly, static_cast<U>(z*z), count);
 }
 
 template <std::size_t N, class T, class V>
 inline V evaluate_even_polynomial(const T(&a)[N], const V& z)
 {
- return evaluate_polynomial(a, z*z);
+ return evaluate_polynomial(a, static_cast<V>(z*z));
 }
 
 template <std::size_t N, class T, class V>
 inline V evaluate_even_polynomial(const boost::array<T,N>& a, const V& z)
 {
- return evaluate_polynomial(a, z*z);
+ return evaluate_polynomial(a, static_cast<V>(z*z));
 }
 //
 // Odd polynomials come next:
@@ -118,21 +118,21 @@
 template <class T, class U>
 inline U evaluate_odd_polynomial(const T* poly, U z, std::size_t count)
 {
- return poly[0] + z * evaluate_polynomial(poly+1, z*z, count-1);
+ return poly[0] + z * evaluate_polynomial(poly+1, static_cast<U>(z*z), count-1);
 }
 
 template <std::size_t N, class T, class V>
 inline V evaluate_odd_polynomial(const T(&a)[N], const V& z)
 {
    typedef mpl::int_<N-1> tag_type;
- return a[0] + z * detail::evaluate_polynomial_c_imp(static_cast<const T*>(a) + 1, z*z, static_cast<tag_type const*>(0));
+ return a[0] + z * detail::evaluate_polynomial_c_imp(static_cast<const T*>(a) + 1, static_cast<V>(z*z), static_cast<tag_type const*>(0));
 }
 
 template <std::size_t N, class T, class V>
 inline V evaluate_odd_polynomial(const boost::array<T,N>& a, const V& z)
 {
    typedef mpl::int_<N-1> tag_type;
- return a[0] + z * detail::evaluate_polynomial_c_imp(static_cast<const T*>(a.data()) + 1, z*z, static_cast<tag_type const*>(0));
+ return a[0] + z * detail::evaluate_polynomial_c_imp(static_cast<const T*>(a.data()) + 1, static_cast<V>(z*z), static_cast<tag_type const*>(0));
 }
 
 template <class T, class U, class V>

Modified: sandbox/interval_math_toolkit/boost/math/tools/real_cast.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/tools/real_cast.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/tools/real_cast.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -10,6 +10,8 @@
 #pragma once
 #endif
 
+#include <boost/math/tools/config.hpp>
+
 namespace boost{ namespace math
 {
   namespace tools
@@ -19,6 +21,11 @@
     {
        return static_cast<To>(t);
     }
+ template <class To, class T, class Policy>
+ inline To real_cast(boost::numeric::interval<T, Policy> t)
+ {
+ return real_cast<To>(norm(t));
+ }
   } // namespace tools
 } // namespace math
 } // namespace boost

Modified: sandbox/interval_math_toolkit/boost/math/tools/roots.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/tools/roots.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/tools/roots.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -44,20 +44,20 @@
                             T& delta,
                             T& result,
                             T& guess,
- const T& min,
- const T& max)
+ const T& vmin,
+ const T& vmax)
 {
    if(last_f0 == 0)
    {
       // this must be the first iteration, pretend that we had a
- // previous one at either min or max:
- if(result == min)
+ // previous one at either vmin or vmax:
+ if(result == vmin)
       {
- guess = max;
+ guess = vmax;
       }
       else
       {
- guess = min;
+ guess = vmin;
       }
       last_f0 = std::tr1::get<0>(f(guess));
       delta = guess - result;
@@ -67,11 +67,11 @@
       // we've crossed over so move in opposite direction to last step:
       if(delta < 0)
       {
- delta = (result - min) / 2;
+ delta = (result - vmin) / 2;
       }
       else
       {
- delta = (result - max) / 2;
+ delta = (result - vmax) / 2;
       }
    }
    else
@@ -79,11 +79,11 @@
       // move in same direction as last step:
       if(delta < 0)
       {
- delta = (result - max) / 2;
+ delta = (result - vmax) / 2;
       }
       else
       {
- delta = (result - min) / 2;
+ delta = (result - vmin) / 2;
       }
    }
 }
@@ -91,28 +91,28 @@
 } // namespace
 
 template <class F, class T, class Tol, class Policy>
-std::pair<T, T> bisect(F f, T min, T max, Tol tol, boost::uintmax_t& max_iter, const Policy& pol)
+std::pair<T, T> bisect(F f, T vmin, T vmax, Tol tol, boost::uintmax_t& max_iter, const Policy& pol)
 {
- T fmin = f(min);
- T fmax = f(max);
+ T fmin = f(vmin);
+ T fmax = f(vmax);
    if(fmin == 0)
- return std::make_pair(min, min);
+ return std::make_pair(vmin, vmin);
    if(fmax == 0)
- return std::make_pair(max, max);
+ return std::make_pair(vmax, vmax);
 
    //
    // Error checking:
    //
    static const char* function = "boost::math::tools::bisect<%1%>";
- if(min >= max)
+ if(vmin >= vmax)
    {
       policies::raise_evaluation_error(function,
- "Arguments in wrong order in boost::math::tools::bisect (first arg=%1%)", min, pol);
+ "Arguments in wrong order in boost::math::tools::bisect (first arg=%1%)", vmin, pol);
    }
    if(fmin * fmax >= 0)
    {
       policies::raise_evaluation_error(function,
- "No change of sign in boost::math::tools::bisect, either there is no root to find, or there are multiple roots in the interval (f(min) = %1%).", fmin, pol);
+ "No change of sign in boost::math::tools::bisect, either there is no root to find, or there are multiple roots in the interval (f(vmin) = %1%).", fmin, pol);
    }
 
    //
@@ -124,25 +124,25 @@
    else
       count -= 3;
 
- while(count && (0 == tol(min, max)))
+ while(count && (0 == tol(vmin, vmax)))
    {
- T mid = (min + max) / 2;
+ T mid = (vmin + vmax) / 2;
       T fmid = f(mid);
- if((mid == max) || (mid == min))
+ if((mid == vmax) || (mid == vmin))
          break;
       if(fmid == 0)
       {
- min = max = mid;
+ vmin = vmax = mid;
          break;
       }
       else if(sign(fmid) * sign(fmin) < 0)
       {
- max = mid;
+ vmax = mid;
          fmax = fmid;
       }
       else
       {
- min = mid;
+ vmin = mid;
          fmin = fmid;
       }
       --count;
@@ -161,24 +161,24 @@
    }
 #endif
 
- return std::make_pair(min, max);
+ return std::make_pair(vmin, vmax);
 }
 
 template <class F, class T, class Tol>
-inline std::pair<T, T> bisect(F f, T min, T max, Tol tol, boost::uintmax_t& max_iter)
+inline std::pair<T, T> bisect(F f, T vmin, T vmax, Tol tol, boost::uintmax_t& max_iter)
 {
- return bisect(f, min, max, tol, max_iter, policies::policy<>());
+ return bisect(f, vmin, vmax, tol, max_iter, policies::policy<>());
 }
 
 template <class F, class T, class Tol>
-inline std::pair<T, T> bisect(F f, T min, T max, Tol tol)
+inline std::pair<T, T> bisect(F f, T vmin, T vmax, Tol tol)
 {
    boost::uintmax_t m = (std::numeric_limits<boost::uintmax_t>::max)();
- return bisect(f, min, max, tol, m, policies::policy<>());
+ return bisect(f, vmin, vmax, tol, m, policies::policy<>());
 }
 
 template <class F, class T>
-T newton_raphson_iterate(F f, T guess, T min, T max, int digits, boost::uintmax_t& max_iter)
+T newton_raphson_iterate(F f, T guess, T vmin, T vmax, int digits, boost::uintmax_t& max_iter)
 {
    BOOST_MATH_STD_USING
 
@@ -205,7 +205,7 @@
 #ifdef BOOST_MATH_INSTRUMENT
          std::cout << "Newton iteration, zero derivative found" << std::endl;
 #endif
- detail::handle_zero_derivative(f, last_f0, f0, delta, result, guess, min, max);
+ detail::handle_zero_derivative(f, last_f0, f0, delta, result, guess, vmin, vmax);
       }
       else
       {
@@ -217,29 +217,29 @@
       if(fabs(delta * 2) > fabs(delta2))
       {
          // last two steps haven't converged, try bisection:
- delta = (delta > 0) ? (result - min) / 2 : (result - max) / 2;
+ delta = (delta > 0) ? (result - vmin) / 2 : (result - vmax) / 2;
       }
       guess = result;
       result -= delta;
- if(result <= min)
+ if(result <= vmin)
       {
- delta = 0.5F * (guess - min);
+ delta = 0.5F * (guess - vmin);
          result = guess - delta;
- if((result == min) || (result == max))
+ if((result == vmin) || (result == vmax))
             break;
       }
- else if(result >= max)
+ else if(result >= vmax)
       {
- delta = 0.5F * (guess - max);
+ delta = 0.5F * (guess - vmax);
          result = guess - delta;
- if((result == min) || (result == max))
+ if((result == vmin) || (result == vmax))
             break;
       }
       // update brackets:
       if(delta > 0)
- max = guess;
+ vmax = guess;
       else
- min = guess;
+ vmin = guess;
    }while(--count && (fabs(result * factor) < fabs(delta)));
 
    max_iter -= count;
@@ -259,22 +259,23 @@
 }
 
 template <class F, class T>
-inline T newton_raphson_iterate(F f, T guess, T min, T max, int digits)
+inline T newton_raphson_iterate(F f, T guess, T vmin, T vmax, int digits)
 {
    boost::uintmax_t m = (std::numeric_limits<boost::uintmax_t>::max)();
- return newton_raphson_iterate(f, guess, min, max, digits, m);
+ return newton_raphson_iterate(f, guess, vmin, vmax, digits, m);
 }
 
 template <class F, class T>
-T halley_iterate(F f, T guess, T min, T max, int digits, boost::uintmax_t& max_iter)
+T halley_iterate(F f, T guess, T vmin, T vmax, int digits, boost::uintmax_t& max_iter)
 {
    BOOST_MATH_STD_USING
+ using std::max;
 
    T f0(0), f1, f2;
    T result = guess;
 
    T factor = static_cast<T>(ldexp(1.0, 1 - digits));
- T delta = (std::max)(10000000 * guess, T(10000000)); // arbitarily large delta
+ T delta = max BOOST_PREVENT_MACRO_SUBSTITUTION(static_cast<T>(10000000 * guess), T(10000000)); // arbitarily large delta
    T last_f0 = 0;
    T delta1 = delta;
    T delta2 = delta;
@@ -300,7 +301,7 @@
 #ifdef BOOST_MATH_INSTRUMENT
          std::cout << "Halley iteration, zero derivative found" << std::endl;
 #endif
- detail::handle_zero_derivative(f, last_f0, f0, delta, result, guess, min, max);
+ detail::handle_zero_derivative(f, last_f0, f0, delta, result, guess, vmin, vmax);
       }
       else
       {
@@ -331,7 +332,7 @@
       if((convergence > 0.8) && (convergence < 2))
       {
          // last two steps haven't converged, try bisection:
- delta = (delta > 0) ? (result - min) / 2 : (result - max) / 2;
+ delta = (delta > 0) ? (result - vmin) / 2 : (result - vmax) / 2;
          // reset delta2 so that this branch will *not* be taken on the
          // next iteration:
          delta2 = delta * 3;
@@ -340,53 +341,53 @@
       result -= delta;
 
       // check for out of bounds step:
- if(result < min)
+ if(result < vmin)
       {
- T diff = ((fabs(min) < 1) && (fabs(result) > 1) && (tools::max_value<T>() / fabs(result) < fabs(min))) ? 1000 : result / min;
+ T diff = ((fabs(vmin) < 1) && (fabs(result) > 1) && (tools::max_value<T>() / fabs(result) < fabs(vmin))) ? static_cast<T>(1000) : static_cast<T>(result / vmin);
          if(fabs(diff) < 1)
             diff = 1 / diff;
          if(!out_of_bounds_sentry && (diff > 0) && (diff < 3))
          {
             // Only a small out of bounds step, lets assume that the result
- // is probably approximately at min:
- delta = 0.99f * (guess - min);
+ // is probably approximately at vmin:
+ delta = 0.99f * (guess - vmin);
             result = guess - delta;
             out_of_bounds_sentry = true; // only take this branch once!
          }
          else
          {
- delta = (guess - min) / 2;
+ delta = (guess - vmin) / 2;
             result = guess - delta;
- if((result == min) || (result == max))
+ if((result == vmin) || (result == vmax))
                break;
          }
       }
- else if(result > max)
+ else if(result > vmax)
       {
- T diff = ((fabs(max) < 1) && (fabs(result) > 1) && (tools::max_value<T>() / fabs(result) < fabs(max))) ? 1000 : result / max;
+ T diff = ((fabs(vmax) < 1) && (fabs(result) > 1) && (tools::max_value<T>() / fabs(result) < fabs(vmax))) ? static_cast<T>(1000) : static_cast<T>(result / vmax);
          if(fabs(diff) < 1)
             diff = 1 / diff;
          if(!out_of_bounds_sentry && (diff > 0) && (diff < 3))
          {
             // Only a small out of bounds step, lets assume that the result
- // is probably approximately at min:
- delta = 0.99f * (guess - max);
+ // is probably approximately at vmin:
+ delta = 0.99f * (guess - vmax);
             result = guess - delta;
             out_of_bounds_sentry = true; // only take this branch once!
          }
          else
          {
- delta = (guess - max) / 2;
+ delta = (guess - vmax) / 2;
             result = guess - delta;
- if((result == min) || (result == max))
+ if((result == vmin) || (result == vmax))
                break;
          }
       }
       // update brackets:
       if(delta > 0)
- max = guess;
+ vmax = guess;
       else
- min = guess;
+ vmin = guess;
    }while(--count && (fabs(result * factor) < fabs(delta)));
 
    max_iter -= count;
@@ -406,14 +407,14 @@
 }
 
 template <class F, class T>
-inline T halley_iterate(F f, T guess, T min, T max, int digits)
+inline T halley_iterate(F f, T guess, T vmin, T vmax, int digits)
 {
    boost::uintmax_t m = (std::numeric_limits<boost::uintmax_t>::max)();
- return halley_iterate(f, guess, min, max, digits, m);
+ return halley_iterate(f, guess, vmin, vmax, digits, m);
 }
 
 template <class F, class T>
-T schroeder_iterate(F f, T guess, T min, T max, int digits, boost::uintmax_t& max_iter)
+T schroeder_iterate(F f, T guess, T vmin, T vmax, int digits, boost::uintmax_t& max_iter)
 {
    BOOST_MATH_STD_USING
 
@@ -444,7 +445,7 @@
 #ifdef BOOST_MATH_INSTRUMENT
          std::cout << "Halley iteration, zero derivative found" << std::endl;
 #endif
- detail::handle_zero_derivative(f, last_f0, f0, delta, result, guess, min, max);
+ detail::handle_zero_derivative(f, last_f0, f0, delta, result, guess, vmin, vmax);
       }
       else
       {
@@ -462,32 +463,32 @@
       if(fabs(delta * 2) > fabs(delta2))
       {
          // last two steps haven't converged, try bisection:
- delta = (delta > 0) ? (result - min) / 2 : (result - max) / 2;
+ delta = (delta > 0) ? (result - vmin) / 2 : (result - vmax) / 2;
       }
       guess = result;
       result -= delta;
 #ifdef BOOST_MATH_INSTRUMENT
       std::cout << "Halley iteration, delta = " << delta << std::endl;
 #endif
- if(result <= min)
+ if(result <= vmin)
       {
- delta = 0.5F * (guess - min);
+ delta = 0.5F * (guess - vmin);
          result = guess - delta;
- if((result == min) || (result == max))
+ if((result == vmin) || (result == vmax))
             break;
       }
- else if(result >= max)
+ else if(result >= vmax)
       {
- delta = 0.5F * (guess - max);
+ delta = 0.5F * (guess - vmax);
          result = guess - delta;
- if((result == min) || (result == max))
+ if((result == vmin) || (result == vmax))
             break;
       }
       // update brackets:
       if(delta > 0)
- max = guess;
+ vmax = guess;
       else
- min = guess;
+ vmin = guess;
    }while(--count && (fabs(result * factor) < fabs(delta)));
 
    max_iter -= count;
@@ -507,10 +508,10 @@
 }
 
 template <class F, class T>
-inline T schroeder_iterate(F f, T guess, T min, T max, int digits)
+inline T schroeder_iterate(F f, T guess, T vmin, T vmax, int digits)
 {
    boost::uintmax_t m = (std::numeric_limits<boost::uintmax_t>::max)();
- return schroeder_iterate(f, guess, min, max, digits, m);
+ return schroeder_iterate(f, guess, vmin, vmax, digits, m);
 }
 
 } // namespace tools

Modified: sandbox/interval_math_toolkit/boost/math/tools/series.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/tools/series.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/tools/series.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -33,7 +33,7 @@
       next_term = func();
       result += next_term;
    }
- while(fabs(result) < fabs(factor * next_term));
+ while(!maybe_greater_equal(fabs(result), fabs(factor * next_term)));
    return result;
 }
 
@@ -53,7 +53,7 @@
       next_term = func();
       result += next_term;
    }
- while((fabs(result) < fabs(factor * next_term)) && --counter);
+ while(!maybe_greater_equal(fabs(result), fabs(factor * next_term)) && --counter);
 
    // set max_terms to the actual number of terms of the series evaluated:
    max_terms = max_terms - counter;
@@ -75,7 +75,7 @@
       next_term = func();
       result += next_term;
    }
- while(fabs(result) < fabs(factor * next_term));
+ while(!maybe_greater_equal(fabs(result), fabs(factor * next_term)));
 
    return result;
 }
@@ -96,7 +96,7 @@
       next_term = func();
       result += next_term;
    }
- while((fabs(result) < fabs(factor * next_term)) && --counter);
+ while(!maybe_greater_equal(fabs(result), fabs(factor * next_term)) && --counter);
 
    // set max_terms to the actual number of terms of the series evaluated:
    max_terms = max_terms - counter;
@@ -135,7 +135,7 @@
       carry -= y;
       result = t;
    }
- while(fabs(result) < fabs(factor * next_term));
+ while(!maybe_greater_equal(fabs(result), fabs(factor * next_term)));
    return result;
 }
 
@@ -160,7 +160,7 @@
       carry -= y;
       result = t;
    }
- while((fabs(result) < fabs(factor * next_term)) && --counter);
+ while(!maybe_greater_equal(fabs(result), fabs(factor * next_term)) && --counter);
 
    // set max_terms to the actual number of terms of the series evaluated:
    max_terms = max_terms - counter;

Modified: sandbox/interval_math_toolkit/boost/math/tools/test.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/tools/test.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/tools/test.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -14,6 +14,7 @@
 #include <boost/math/tools/stats.hpp>
 #include <boost/math/special_functions/fpclassify.hpp>
 #include <boost/test/test_tools.hpp>
+#include <boost/mpl/bool.hpp>
 #include <stdexcept>
 
 namespace boost{ namespace math{ namespace tools{
@@ -158,7 +159,7 @@
    {
       if(i)
          std::cout << ", ";
- std::cout << row[i];
+ std::cout << boost::numeric::make_printable(row[i]);
    }
    std::cout << std::endl;
 }
@@ -171,7 +172,7 @@
 // expressions.
 //
 template <class A, class F1, class F2>
-test_result<typename calculate_result_type<A>::value_type> test(const A& a, F1 test_func, F2 expect_func)
+test_result<typename calculate_result_type<A>::value_type> test(const A& a, F1 test_func, F2 expect_func, const mpl::false_&)
 {
    typedef typename A::value_type row_type;
    typedef typename row_type::value_type value_type;
@@ -238,6 +239,116 @@
    return result;
 }
 
+template <class A, class F1, class F2>
+test_result<typename calculate_result_type<A>::value_type> test(const A& a, F1 test_func, F2 expect_func, const mpl::true_&)
+{
+ typedef typename A::value_type row_type;
+ typedef typename row_type::value_type value_type;
+
+ test_result<value_type> result;
+
+ for(unsigned i = 0; i < a.size(); ++i)
+ {
+ const row_type& row = a[i];
+ value_type point;
+ try
+ {
+ point = test_func(row);
+ }
+ catch(const std::underflow_error&)
+ {
+ point = 0;
+ }
+ catch(const std::overflow_error&)
+ {
+ point = std::numeric_limits<value_type>::has_infinity ?
+ std::numeric_limits<value_type>::infinity()
+ : tools::max_value<value_type>();
+ }
+ catch(const std::exception& e)
+ {
+ std::cerr << e.what() << std::endl;
+ print_row(row);
+ BOOST_ERROR("Unexpected exception.");
+ // so we don't get further errors:
+ point = expect_func(row);
+ }
+ value_type expected = expect_func(row);
+ value_type err = (median(point) != 0)
+ ? width(point) / median(point)
+ : width(point);
+
+#ifdef BOOST_INSTRUMENT
+ if(err != 0)
+ {
+ std::cout << std::setprecision(35);
+ std::cout << row[0] << " " << err;
+ std::cout << " (" << err.lower() / tools::epsilon<value_type>() << "eps)";
+ std::cout << std::endl;
+ }
+#endif
+ if(!(boost::math::isfinite)(point) && (boost::math::isfinite)(expected))
+ {
+ std::cout << std::setprecision(35);
+ std::cout << "CAUTION: Found non-finite result, when a finite value was expected at entry " << i << "\n";
+ std::cout << "Found: " << boost::numeric::make_printable(point) << " Expected " << boost::numeric::make_printable(expected) << " Error: " << boost::numeric::make_printable(err) << std::endl;
+ print_row(row);
+ BOOST_ERROR("Unexpected non-finite result");
+ }
+ if(err > 0.5)
+ {
+ std::cout << std::setprecision(35);
+ std::cout << "CAUTION: Gross error found at entry " << i << ".\n";
+ std::cout << "Found: " << boost::numeric::make_printable(point) << " Expected " << boost::numeric::make_printable(expected) << " Error: " << boost::numeric::make_printable(err) << std::endl;
+ print_row(row);
+ BOOST_ERROR("Gross error");
+ }
+ //
+ // Check whether the interval includes the result:
+ //
+ if(!overlap(point, expected))
+ {
+ //
+ // Since the transcendental functions do *not* guarentee
+ // to bracket the true value, we only raise an error if
+ // we're more than a few epsilon from the true value:
+ //
+ bool have_error = false;
+ int max_eps = 3;
+ typename value_type::base_type err = 0;
+ if(point.upper() < expected.lower())
+ err = (relative_error(point.upper(), expected.lower()) / epsilon<typename value_type::base_type>());
+ else if(point.lower() > expected.upper())
+ err = (relative_error(point.lower(), expected.upper()) / epsilon<typename value_type::base_type>());
+ have_error = max_eps < err;
+
+ if(have_error)
+ {
+ std::cout << std::setprecision(35);
+ std::cout << "Interval and expected result do not overlap.\n";
+ std::cout << "Found: " << boost::numeric::make_printable(point) << " Expected " << boost::numeric::make_printable(expected) << " Error: " << boost::numeric::make_printable(err) << std::endl;
+ print_row(row);
+ BOOST_ERROR("Non Overlapping Interval Error");
+ }
+ }
+
+ result.add(err);
+ if(maybe_greater_equal((result.max)(), err))
+ result.set_worst(i);
+ }
+ return result;
+}
+
+template <class A, class F1, class F2>
+inline test_result<typename calculate_result_type<A>::value_type> test(const A& a, F1 test_func, F2 expect_func)
+{
+ typedef typename A::value_type row_type;
+ typedef typename row_type::value_type value_type;
+ typedef typename boost::numeric::is_interval<value_type>::type tag_type;
+
+ return test(a, test_func, expect_func, tag_type());
+}
+
 } // namespace tools
 } // namespace math
 } // namespace boost

Modified: sandbox/interval_math_toolkit/boost/math/tools/toms748_solve.hpp
==============================================================================
--- sandbox/interval_math_toolkit/boost/math/tools/toms748_solve.hpp (original)
+++ sandbox/interval_math_toolkit/boost/math/tools/toms748_solve.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -26,12 +26,14 @@
    eps_tolerance(unsigned bits)
    {
       BOOST_MATH_STD_USING
- eps = (std::max)(T(ldexp(1.0F, 1-bits)), 2 * tools::epsilon<T>());
+ using std::max;
+ eps = max BOOST_PREVENT_MACRO_SUBSTITUTION(T(ldexp(1.0F, 1-bits)), T(2 * 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;
+ using std::min;
+ return (fabs(a - b) / min BOOST_PREVENT_MACRO_SUBSTITUTION(fabs(a), fabs(b))) <= eps;
    }
 private:
    T eps;
@@ -193,9 +195,9 @@
    //
    // Start by obtaining the coefficients of the quadratic polynomial:
    //
- T B = safe_div(fb - fa, b - a, tools::max_value<T>());
- T A = safe_div(fd - fb, d - b, tools::max_value<T>());
- A = safe_div(A - B, d - a, T(0));
+ T B = safe_div(static_cast<T>(fb - fa), static_cast<T>(b - a), tools::max_value<T>());
+ T A = safe_div(static_cast<T>(fd - fb), static_cast<T>(d - b), tools::max_value<T>());
+ A = safe_div(static_cast<T>(A - B), static_cast<T>(d - a), T(0));
 
    if(a == 0)
    {
@@ -220,7 +222,7 @@
    for(unsigned i = 1; i <= count; ++i)
    {
       //c -= safe_div(B * c, (B + A * (2 * c - a - b)), 1 + c - a);
- c -= safe_div(fa+(B+A*(c-b))*(c-a), (B + A * (2 * c - a - b)), 1 + c - a);
+ c -= safe_div(static_cast<T>(fa+(B+A*(c-b))*(c-a)), static_cast<T>((B + A * (2 * c - a - b))), static_cast<T>(1 + c - a));
    }
    if((c <= a) || (c >= b))
    {
@@ -436,7 +438,7 @@
       //
       e = d;
       fe = fd;
- detail::bracket(f, a, b, a + (b - a) / 2, fa, fb, d, fd);
+ detail::bracket(f, a, b, static_cast<T>(a + (b - a) / 2), fa, fb, d, fd);
       --count;
       BOOST_MATH_INSTRUMENT_CODE("Not converging: Taking a bisection!!!!");
       BOOST_MATH_INSTRUMENT_CODE(" a = " << a << " b = " << b);

Modified: sandbox/interval_math_toolkit/libs/math/test/handle_test_result.hpp
==============================================================================
--- sandbox/interval_math_toolkit/libs/math/test/handle_test_result.hpp (original)
+++ sandbox/interval_math_toolkit/libs/math/test/handle_test_result.hpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -111,7 +111,7 @@
                        const char* test_name,
                        const char* group_name)
 {
- using namespace std; // To aid selection of the right pow.
+ BOOST_MATH_STD_USING
    T eps = boost::math::tools::epsilon<T>();
    std::cout << std::setprecision(4);
 
@@ -120,8 +120,8 @@
    //
    // Begin by printing the main tag line with the results:
    //
- std::cout << test_name << "<" << type_name << "> Max = " << max_error_found
- << " RMS Mean=" << mean_error_found;
+ std::cout << test_name << "<" << type_name << "> Max = " << boost::numeric::make_printable(max_error_found)
+ << " RMS Mean=" << boost::numeric::make_printable(mean_error_found);
    //
    // If the max error is non-zero, give the row of the table that
    // produced the worst error:
@@ -137,7 +137,7 @@
 #if defined(__SGI_STL_PORT)
          std::cout << boost::math::tools::real_cast<double>(worst[i]);
 #else
- std::cout << worst[i];
+ std::cout << boost::numeric::make_printable(worst[i]);
 #endif
       }
       std::cout << " }";
@@ -146,16 +146,19 @@
    //
    // Now verify that the results are within our expected bounds:
    //
- std::pair<boost::uintmax_t, boost::uintmax_t> const& bounds = get_max_errors(type_name, test_name, group_name);
- if(bounds.first < max_error_found)
+ if(!boost::numeric::is_interval<T>::value)
    {
- std::cerr << "Peak error greater than expected value of " << bounds.first << std::endl;
- BOOST_CHECK(bounds.first >= max_error_found);
- }
- if(bounds.second < mean_error_found)
- {
- std::cerr << "Mean error greater than expected value of " << bounds.second << std::endl;
- BOOST_CHECK(bounds.second >= mean_error_found);
+ std::pair<boost::uintmax_t, boost::uintmax_t> const& bounds = get_max_errors(type_name, test_name, group_name);
+ if(bounds.first < max_error_found)
+ {
+ std::cerr << "Peak error greater than expected value of " << bounds.first << std::endl;
+ BOOST_CHECK(bounds.first >= max_error_found);
+ }
+ if(bounds.second < mean_error_found)
+ {
+ std::cerr << "Mean error greater than expected value of " << bounds.second << std::endl;
+ BOOST_CHECK(bounds.second >= mean_error_found);
+ }
    }
    std::cout << std::endl;
 }
@@ -194,5 +197,23 @@
    std::cout << std::endl;
 }
 
+#ifdef TEST_INTERVAL
+
+#include <boost/math/bindings/interval.hpp>
+
+typedef boost::numeric::interval<double,
+ boost::numeric::interval_lib::policies<
+ boost::numeric::interval_lib::save_state<
+ boost::numeric::interval_lib::rounded_transc_nearest<
+ double,
+ boost::numeric::interval_lib::rounded_arith_opp<double>
+ >
+ >,
+ boost::numeric::interval_lib::checking_base<double>
+ >
+> interval_type;
+
+#endif
+
 #endif // BOOST_MATH_HANDLE_TEST_RESULT
 

Added: sandbox/interval_math_toolkit/libs/math/test/interval_concept_check.cpp
==============================================================================
--- (empty file)
+++ sandbox/interval_math_toolkit/libs/math/test/interval_concept_check.cpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -0,0 +1,37 @@
+// Copyright John Maddock 2007-8.
+// Use, modification and distribution are subject to the
+// Boost Software License, Version 1.0. (See accompanying file
+// LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
+
+//
+// This tests two things: that boost::interval meets our
+// conceptual requirements, and that we can instantiate
+// all our distributions and special functions on this type.
+//
+#define BOOST_MATH_INSTRUMENT
+
+#define BOOST_MATH_ASSERT_UNDEFINED_POLICY false
+
+#include <boost/math/bindings/interval.hpp>
+#include <boost/math/concepts/real_type_concept.hpp>
+#include "compile_test/instantiate.hpp"
+
+
+typedef boost::numeric::interval<double,
+ boost::numeric::interval_lib::policies<
+ boost::numeric::interval_lib::save_state<
+ boost::numeric::interval_lib::rounded_transc_opp<double> >,
+ boost::numeric::interval_lib::checking_base<double> > > interval_type;
+
+
+void foo()
+{
+ instantiate(interval_type());
+}
+
+int main()
+{
+ BOOST_CONCEPT_ASSERT((boost::math::concepts::RealTypeConcept<double>));
+ BOOST_CONCEPT_ASSERT((boost::math::concepts::RealTypeConcept<interval_type>));
+}
+

Modified: sandbox/interval_math_toolkit/libs/math/test/log1p_expm1_test.cpp
==============================================================================
--- sandbox/interval_math_toolkit/libs/math/test/log1p_expm1_test.cpp (original)
+++ sandbox/interval_math_toolkit/libs/math/test/log1p_expm1_test.cpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -5,7 +5,12 @@
 
 #define BOOST_MATH_OVERFLOW_ERROR_POLICY ignore_error
 
+#ifdef TEST_INTERVAL
+#include <boost/math/bindings/interval.hpp>
+#else
 #include <boost/math/concepts/real_concept.hpp>
+#endif
+
 #include <boost/test/included/test_exec_monitor.hpp>
 #include <boost/test/floating_point_comparison.hpp>
 #include <boost/math/special_functions/log1p.hpp>
@@ -104,6 +109,7 @@
 
    do_test(log1p_expm1_data, type_name, "expm1 and log1p");
 
+#ifndef TEST_INTERVAL
    //
    // C99 Appendix F special cases:
    static const T zero = 0;
@@ -117,11 +123,13 @@
       BOOST_CHECK_EQUAL(boost::math::expm1(-std::numeric_limits<T>::infinity()), m_one);
       BOOST_CHECK_EQUAL(boost::math::expm1(std::numeric_limits<T>::infinity()), std::numeric_limits<T>::infinity());
    }
+#endif
 }
 
 
 int test_main(int, char* [])
 {
+#ifndef TEST_INTERVAL
    expected_results();
    BOOST_MATH_CONTROL_FP;
    test(float(0), "float");
@@ -142,5 +150,8 @@
       "not available at all, or because they are too inaccurate for these tests "
       "to pass.</note>" << std::cout;
 #endif
+#else // TEST_INTERVAL
+ test((interval_type)(0), "boost::interval<double>");
+#endif // TEST_INTERVAL
    return 0;
 }

Added: sandbox/interval_math_toolkit/libs/math/test/mpfr_concept_check.cpp
==============================================================================
--- (empty file)
+++ sandbox/interval_math_toolkit/libs/math/test/mpfr_concept_check.cpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -0,0 +1,30 @@
+// Copyright John Maddock 2007-8.
+// Use, modification and distribution are subject to the
+// Boost Software License, Version 1.0. (See accompanying file
+// LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
+
+//
+// This tests two things: that mpfr_class meets our
+// conceptual requirements, and that we can instantiate
+// all our distributions and special functions on this type.
+//
+#define BOOST_MATH_ASSERT_UNDEFINED_POLICY false
+#define TEST_MPFR
+
+#pragma warning(disable:4800)
+
+#include <boost/math/bindings/mpfr.hpp>
+#include <boost/math/concepts/real_type_concept.hpp>
+#include "compile_test/instantiate.hpp"
+
+void foo()
+{
+ //instantiate(mpfr_class());
+}
+
+int main()
+{
+ BOOST_CONCEPT_ASSERT((mpfr_class));
+ BOOST_CONCEPT_ASSERT((mpfr_class));
+}
+

Added: sandbox/interval_math_toolkit/libs/math/test/ntl_concept_check.cpp
==============================================================================
--- (empty file)
+++ sandbox/interval_math_toolkit/libs/math/test/ntl_concept_check.cpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -0,0 +1,27 @@
+// Copyright John Maddock 2007-8.
+// Use, modification and distribution are subject to the
+// Boost Software License, Version 1.0. (See accompanying file
+// LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
+
+//
+// This tests two things: that boost::math::ntl::RR meets our
+// conceptual requirements, and that we can instantiate
+// all our distributions and special functions on this type.
+//
+#define BOOST_MATH_ASSERT_UNDEFINED_POLICY false
+
+#include <boost/math/bindings/rr.hpp>
+#include <boost/math/concepts/real_type_concept.hpp>
+#include "compile_test/instantiate.hpp"
+
+void foo()
+{
+ instantiate(boost::math::ntl::RR());
+}
+
+int main()
+{
+ BOOST_CONCEPT_ASSERT((boost::math::ntl::RR));
+ BOOST_CONCEPT_ASSERT((boost::math::ntl::RR));
+}
+

Modified: sandbox/interval_math_toolkit/libs/math/test/powm1_sqrtp1m1_test.cpp
==============================================================================
--- sandbox/interval_math_toolkit/libs/math/test/powm1_sqrtp1m1_test.cpp (original)
+++ sandbox/interval_math_toolkit/libs/math/test/powm1_sqrtp1m1_test.cpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -3,7 +3,12 @@
 // Boost Software License, Version 1.0. (See accompanying file
 // LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
 
+#ifdef TEST_INTERVAL
+#include <boost/math/bindings/interval.hpp>
+#else
 #include <boost/math/concepts/real_concept.hpp>
+#endif
+
 #include <boost/test/included/test_exec_monitor.hpp>
 #include <boost/test/floating_point_comparison.hpp>
 #include <boost/math/special_functions/sqrt1pm1.hpp>
@@ -1637,6 +1642,7 @@
 
 int test_main(int, char* [])
 {
+#ifndef TEST_INTERVAL
    expected_results();
    BOOST_MATH_CONTROL_FP;
    test_powm1_sqrtp1m1(1.0F, "float");
@@ -1652,6 +1658,9 @@
       "not available at all, or because they are too inaccurate for these tests "
       "to pass.</note>" << std::cout;
 #endif
+#else // TEST_INTERVAL
+ test_powm1_sqrtp1m1(interval_type(), "boost::numeric::interval<double>");
+#endif // TEST_INTERVAL
    return 0;
 }
 

Modified: sandbox/interval_math_toolkit/libs/math/test/test_beta.cpp
==============================================================================
--- sandbox/interval_math_toolkit/libs/math/test/test_beta.cpp (original)
+++ sandbox/interval_math_toolkit/libs/math/test/test_beta.cpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -5,7 +5,12 @@
 // Boost Software License, Version 1.0. (See accompanying file
 // LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
 
+#ifdef TEST_INTERVAL
+#include <boost/math/bindings/interval.hpp>
+#else
 #include <boost/math/concepts/real_concept.hpp>
+#endif
+
 #include <boost/test/included/test_exec_monitor.hpp>
 #include <boost/test/floating_point_comparison.hpp>
 #include <boost/math/special_functions/beta.hpp>
@@ -186,6 +191,7 @@
 
 int test_main(int, char* [])
 {
+#ifndef TEST_INTERVAL
    expected_results();
    BOOST_MATH_CONTROL_FP;
    test_spots(0.0F);
@@ -217,6 +223,9 @@
       "not available at all, or because they are too inaccurate for these tests "
       "to pass.</note>" << std::cout;
 #endif
+#else // TEST_INTERVAL
+ test_beta(interval_type(), "boost::numeric::interval<double>");
+#endif // TEST_INTERVAL
    return 0;
 }
 

Modified: sandbox/interval_math_toolkit/libs/math/test/test_binomial_coeff.cpp
==============================================================================
--- sandbox/interval_math_toolkit/libs/math/test/test_binomial_coeff.cpp (original)
+++ sandbox/interval_math_toolkit/libs/math/test/test_binomial_coeff.cpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -3,7 +3,12 @@
 // Boost Software License, Version 1.0. (See accompanying file
 // LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
 
+#ifdef TEST_INTERVAL
+#include <boost/math/bindings/interval.hpp>
+#else
 #include <boost/math/concepts/real_concept.hpp>
+#endif
+
 #include <boost/test/included/test_exec_monitor.hpp>
 #include <boost/test/floating_point_comparison.hpp>
 #include <boost/math/special_functions/binomial.hpp>
@@ -163,6 +168,7 @@
 
 int test_main(int, char* [])
 {
+#ifndef TEST_INTERVAL
    expected_results();
    BOOST_MATH_CONTROL_FP;
 
@@ -186,6 +192,9 @@
       "not available at all, or because they are too inaccurate for these tests "
       "to pass.</note>" << std::cout;
 #endif
+#else // TEST_INTERVAL
+ test_binomial(interval_type(), "boost:math::interval<double>");
+#endif // TEST_INTERVAL
    return 0;
 }
 

Modified: sandbox/interval_math_toolkit/libs/math/test/test_carlson.cpp
==============================================================================
--- sandbox/interval_math_toolkit/libs/math/test/test_carlson.cpp (original)
+++ sandbox/interval_math_toolkit/libs/math/test/test_carlson.cpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -5,7 +5,12 @@
 // Boost Software License, Version 1.0. (See accompanying file
 // LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
 
+#ifdef TEST_INTERVAL
+#include <boost/math/bindings/interval.hpp>
+#else
 #include <boost/math/concepts/real_concept.hpp>
+#endif
+
 #include <boost/test/included/test_exec_monitor.hpp>
 #include <boost/test/floating_point_comparison.hpp>
 #include <boost/math/special_functions/ellint_rf.hpp>
@@ -195,7 +200,7 @@
 
    std::cout << "Testing: " << test << std::endl;
 
- value_type (*fp)(value_type, value_type, value_type) = boost::math::ellint_rd;
+ value_type (*fp)(value_type, value_type, value_type) = boost::math::ellint_rd;
     boost::math::tools::test_result<value_type> result;
  
     result = boost::math::tools::test(
@@ -214,6 +219,8 @@
 {
    using namespace boost::math;
    using namespace std;
+
+#ifndef TEST_INTERVAL
    // Spot values from Numerical Computation of Real or Complex
    // Elliptic Integrals, B. C. Carlson: http://arxiv.org/abs/math.CA/9409227
    // RF:
@@ -281,7 +288,7 @@
       s2 = ellint_rj(x, y, z, z);
       BOOST_CHECK_CLOSE_FRACTION(s1, s2, eps40);
    }
-
+#endif // TEST_INTERVAL
    //
    // Now random spot values:
    //
@@ -304,6 +311,7 @@
 
 int test_main(int, char* [])
 {
+#ifndef TEST_INTERVAL
     expected_results();
     BOOST_MATH_CONTROL_FP;
 
@@ -322,7 +330,9 @@
       "not available at all, or because they are too inaccurate for these tests "
       "to pass.</note>" << std::cout;
 #endif
-
+#else // TEST_INTERVAL
+ test_spots(interval_type(0), "boost::numeric::interval<double>");
+#endif // TEST_INTERVAL
     return 0;
 }
 

Modified: sandbox/interval_math_toolkit/libs/math/test/test_cbrt.cpp
==============================================================================
--- sandbox/interval_math_toolkit/libs/math/test/test_cbrt.cpp (original)
+++ sandbox/interval_math_toolkit/libs/math/test/test_cbrt.cpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -3,7 +3,12 @@
 // Boost Software License, Version 1.0. (See accompanying file
 // LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
 
+#ifdef TEST_INTERVAL
+#include <boost/math/bindings/interval.hpp>
+#else
 #include <boost/math/concepts/real_concept.hpp>
+#endif
+
 #include <boost/test/included/test_exec_monitor.hpp>
 #include <boost/test/floating_point_comparison.hpp>
 #include <boost/math/tools/stats.hpp>
@@ -115,6 +120,7 @@
 
 int test_main(int, char* [])
 {
+#ifndef TEST_INTERVAL
    expected_results();
    BOOST_MATH_CONTROL_FP;
    test_cbrt(0.1F, "float");
@@ -125,6 +131,9 @@
    test_cbrt(boost::math::concepts::real_concept(0.1), "real_concept");
 #endif
 #endif
+#else // TEST_INTERVAL
+ test_cbrt(interval_type(0.1), "boost::numeric::interval<double>");
+#endif // TEST_INTERVAL
    return 0;
 }
 

Modified: sandbox/interval_math_toolkit/libs/math/test/test_digamma.cpp
==============================================================================
--- sandbox/interval_math_toolkit/libs/math/test/test_digamma.cpp (original)
+++ sandbox/interval_math_toolkit/libs/math/test/test_digamma.cpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -3,7 +3,12 @@
 // Boost Software License, Version 1.0. (See accompanying file
 // LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
 
+#ifdef TEST_INTERVAL
+#include <boost/math/bindings/interval.hpp>
+#else
 #include <boost/math/concepts/real_concept.hpp>
+#endif
+
 #include <boost/test/included/test_exec_monitor.hpp>
 #include <boost/test/floating_point_comparison.hpp>
 #include <boost/math/special_functions/digamma.hpp>
@@ -149,6 +154,7 @@
 
 int test_main(int, char* [])
 {
+#ifndef TEST_INTERVAL
    BOOST_MATH_CONTROL_FP;
    test_spots(0.0F, "float");
    test_spots(0.0, "double");
@@ -172,6 +178,9 @@
       "not available at all, or because they are too inaccurate for these tests "
       "to pass.</note>" << std::cout;
 #endif
+#else // TEST_INTERVAL
+ test_digamma(interval_type(0.1), "boost::numeric::interval<double>");
+#endif
    return 0;
 }
 

Modified: sandbox/interval_math_toolkit/libs/math/test/test_ellint_1.cpp
==============================================================================
--- sandbox/interval_math_toolkit/libs/math/test/test_ellint_1.cpp (original)
+++ sandbox/interval_math_toolkit/libs/math/test/test_ellint_1.cpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -11,7 +11,12 @@
 // Constants are too big for float case, but this doesn't matter for test.
 #endif
 
+#ifdef TEST_INTERVAL
+#include <boost/math/bindings/interval.hpp>
+#else
 #include <boost/math/concepts/real_concept.hpp>
+#endif
+
 #include <boost/test/included/test_exec_monitor.hpp>
 #include <boost/test/floating_point_comparison.hpp>
 #include <boost/math/special_functions/ellint_1.hpp>
@@ -192,6 +197,8 @@
 int test_main(int, char* [])
 {
     expected_results();
+
+#ifndef TEST_INTERVAL
     BOOST_MATH_CONTROL_FP;
 
     test_spots(0.0F, "float");
@@ -207,6 +214,8 @@
       "not available at all, or because they are too inaccurate for these tests "
       "to pass.</note>" << std::cout;
 #endif
-
+#else // TEST_INTERVAL
+ test_spots(interval_type(0), "boost::numeric::interval<T>");
+#endif // TEST_INTERVAL
     return 0;
 }

Modified: sandbox/interval_math_toolkit/libs/math/test/test_ellint_2.cpp
==============================================================================
--- sandbox/interval_math_toolkit/libs/math/test/test_ellint_2.cpp (original)
+++ sandbox/interval_math_toolkit/libs/math/test/test_ellint_2.cpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -11,7 +11,12 @@
 // Constants are too big for float case, but this doesn't matter for test.
 #endif
 
+#ifdef TEST_INTERVAL
+#include <boost/math/bindings/interval.hpp>
+#else
 #include <boost/math/concepts/real_concept.hpp>
+#endif
+
 #include <boost/test/included/test_exec_monitor.hpp>
 #include <boost/test/floating_point_comparison.hpp>
 #include <boost/math/special_functions/ellint_2.hpp>
@@ -183,6 +188,7 @@
 int test_main(int, char* [])
 {
     expected_results();
+#ifndef TEST_INTERVAL
     BOOST_MATH_CONTROL_FP;
 #ifndef BOOST_MATH_BUGGY_LARGE_FLOAT_CONSTANTS
     test_spots(0.0F, "float");
@@ -199,6 +205,8 @@
       "not available at all, or because they are too inaccurate for these tests "
       "to pass.</note>" << std::cout;
 #endif
-
+#else // TEST_INTERVAL
+ test_spots(interval_type(0), "boost::numeric::interval<double>");
+#endif // TEST_INTERVAL
     return 0;
 }

Modified: sandbox/interval_math_toolkit/libs/math/test/test_ellint_3.cpp
==============================================================================
--- sandbox/interval_math_toolkit/libs/math/test/test_ellint_3.cpp (original)
+++ sandbox/interval_math_toolkit/libs/math/test/test_ellint_3.cpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -11,7 +11,12 @@
 // Constants are too big for float case, but this doesn't matter for test.
 #endif
 
+#ifdef TEST_INTERVAL
+#include <boost/math/bindings/interval.hpp>
+#else
 #include <boost/math/concepts/real_concept.hpp>
+#endif
+
 #include <boost/test/included/test_exec_monitor.hpp>
 #include <boost/test/floating_point_comparison.hpp>
 #include <boost/math/special_functions/ellint_3.hpp>
@@ -220,6 +225,7 @@
 int test_main(int, char* [])
 {
     expected_results();
+#ifndef TEST_INTERVAL
     BOOST_MATH_CONTROL_FP;
     test_spots(0.0F, "float");
     test_spots(0.0, "double");
@@ -234,6 +240,8 @@
       "not available at all, or because they are too inaccurate for these tests "
       "to pass.</note>" << std::cout;
 #endif
-
+#else // TEST_INTERVAL
+ test_spots(interval_type(0.0), "boost::numeric::interval<double>");
+#endif // TEST_INTERVAL
     return 0;
 }

Modified: sandbox/interval_math_toolkit/libs/math/test/test_erf.cpp
==============================================================================
--- sandbox/interval_math_toolkit/libs/math/test/test_erf.cpp (original)
+++ sandbox/interval_math_toolkit/libs/math/test/test_erf.cpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -4,7 +4,12 @@
 // Boost Software License, Version 1.0. (See accompanying file
 // LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
 
+#ifdef TEST_INTERVAL
+#include <boost/math/bindings/interval.hpp>
+#else
 #include <boost/math/concepts/real_concept.hpp>
+#endif
+
 #include <boost/test/included/test_exec_monitor.hpp>
 #include <boost/test/floating_point_comparison.hpp>
 #include <boost/math/special_functions/erf.hpp>
@@ -317,6 +322,7 @@
 
 int test_main(int, char* [])
 {
+#ifndef TEST_INTERVAL
    BOOST_MATH_CONTROL_FP;
    test_spots(0.0F, "float");
    test_spots(0.0, "double");
@@ -344,6 +350,9 @@
       "not available at all, or because they are too inaccurate for these tests "
       "to pass.</note>" << std::cout;
 #endif
+#else // TEST_INTERVAL
+ test_erf(interval_type(0), "boost::numeric::interval<double>");
+#endif // TEST_INTERVAL
    return 0;
 }
 

Modified: sandbox/interval_math_toolkit/libs/math/test/test_expint.cpp
==============================================================================
--- sandbox/interval_math_toolkit/libs/math/test/test_expint.cpp (original)
+++ sandbox/interval_math_toolkit/libs/math/test/test_expint.cpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -3,7 +3,12 @@
 // Boost Software License, Version 1.0. (See accompanying file
 // LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
 
+#ifdef TEST_INTERVAL
+#include <boost/math/bindings/interval.hpp>
+#else
 #include <boost/math/concepts/real_concept.hpp>
+#endif
+
 #include <boost/math/special_functions/expint.hpp>
 #include <boost/test/included/test_exec_monitor.hpp>
 #include <boost/test/floating_point_comparison.hpp>
@@ -293,6 +298,7 @@
 int test_main(int, char* [])
 {
    expected_results();
+#ifndef TEST_INTERVAL
    BOOST_MATH_CONTROL_FP;
 
    boost::math::expint(114.7);
@@ -320,6 +326,9 @@
       "not available at all, or because they are too inaccurate for these tests "
       "to pass.</note>" << std::cout;
 #endif
+#else // TEST_INTERVAL
+ test_expint(interval_type(), "boost::numeric::interval<double>");
+#endif // TEST_INTERVAL
    return 0;
 }
 

Modified: sandbox/interval_math_toolkit/libs/math/test/test_gamma.cpp
==============================================================================
--- sandbox/interval_math_toolkit/libs/math/test/test_gamma.cpp (original)
+++ sandbox/interval_math_toolkit/libs/math/test/test_gamma.cpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -5,7 +5,12 @@
 
 #define BOOST_MATH_OVERFLOW_ERROR_POLICY ignore_error
 
+#ifdef TEST_INTERVAL
+#include <boost/math/bindings/interval.hpp>
+#else
 #include <boost/math/concepts/real_concept.hpp>
+#endif
+
 #include <boost/math/special_functions/gamma.hpp>
 #include <boost/test/included/test_exec_monitor.hpp>
 #include <boost/test/floating_point_comparison.hpp>
@@ -450,6 +455,7 @@
 int test_main(int, char* [])
 {
    expected_results();
+#ifndef TEST_INTERVAL
    BOOST_MATH_CONTROL_FP;
 
 #ifndef BOOST_MATH_BUGGY_LARGE_FLOAT_CONSTANTS
@@ -478,6 +484,9 @@
       "not available at all, or because they are too inaccurate for these tests "
       "to pass.</note>" << std::cout;
 #endif
+#else // TEST_INTERVAL
+ test_gamma(interval_type(), "boost::numeric::interval<double>");
+#endif // TEST_INTERVAL
    return 0;
 }
 

Modified: sandbox/interval_math_toolkit/libs/math/test/test_hermite.cpp
==============================================================================
--- sandbox/interval_math_toolkit/libs/math/test/test_hermite.cpp (original)
+++ sandbox/interval_math_toolkit/libs/math/test/test_hermite.cpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -12,7 +12,12 @@
 // Constants are too big for float case, but this doesn't matter for test.
 #endif
 
+#ifdef TEST_INTERVAL
+#include <boost/math/bindings/interval.hpp>
+#else
 #include <boost/math/concepts/real_concept.hpp>
+#endif
+
 #include <boost/test/included/test_exec_monitor.hpp>
 #include <boost/test/floating_point_comparison.hpp>
 #include <boost/math/special_functions/hermite.hpp>
@@ -168,8 +173,7 @@
 {
    BOOST_MATH_CONTROL_FP;
 
- boost::math::hermite(51, 915.0);
-
+#ifndef TEST_INTERVAL
 #ifndef BOOST_MATH_BUGGY_LARGE_FLOAT_CONSTANTS
    test_spots(0.0F, "float");
 #endif
@@ -196,6 +200,9 @@
       "not available at all, or because they are too inaccurate for these tests "
       "to pass.</note>" << std::cout;
 #endif
+#else // TEST_INTERVAL
+ test_hermite(interval_type(), "boost::numeric::interval<double>");
+#endif // TEST_INTERVAL
    return 0;
 }
 

Modified: sandbox/interval_math_toolkit/libs/math/test/test_ibeta.cpp
==============================================================================
--- sandbox/interval_math_toolkit/libs/math/test/test_ibeta.cpp (original)
+++ sandbox/interval_math_toolkit/libs/math/test/test_ibeta.cpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -3,7 +3,12 @@
 // Boost Software License, Version 1.0. (See accompanying file
 // LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
 
+#ifdef TEST_INTERVAL
+#include <boost/math/bindings/interval.hpp>
+#else
 #include <boost/math/concepts/real_concept.hpp>
+#endif
+
 #include <boost/test/included/test_exec_monitor.hpp>
 #include <boost/test/floating_point_comparison.hpp>
 #include <boost/math/special_functions/beta.hpp>
@@ -515,6 +520,7 @@
 int test_main(int, char* [])
 {
    expected_results();
+#ifndef TEST_INTERVAL
    BOOST_MATH_CONTROL_FP;
 #ifdef TEST_GSL
    gsl_set_error_handler_off();
@@ -557,6 +563,9 @@
       "not available at all, or because they are too inaccurate for these tests "
       "to pass.</note>" << std::cout;
 #endif
+#else // TEST_INTERVAL
+ test_beta(interval_type(), "boost::numeric::interval<double>");
+#endif // TEST_INTERVAL
    return 0;
 }
 

Modified: sandbox/interval_math_toolkit/libs/math/test/test_ibeta_inv.cpp
==============================================================================
--- sandbox/interval_math_toolkit/libs/math/test/test_ibeta_inv.cpp (original)
+++ sandbox/interval_math_toolkit/libs/math/test/test_ibeta_inv.cpp 2008-01-05 05:09:55 EST (Sat, 05 Jan 2008)
@@ -3,7 +3,12 @@
 // Boost Software License, Version 1.0. (See accompanying file
 // LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
 
+#ifdef TEST_INTERVAL
+#include <boost/math/bindings/interval.hpp>
+#else
 #include <boost/math/concepts/real_concept.hpp>
+#endif
+
 #include <boost/test/included/test_exec_monitor.hpp>
 #include <boost/test/floating_point_comparison.hpp>
 #include <boost/math/special_functions/beta.hpp>
@@ -270,6 +275,7 @@
    // The contents are as follows, each row of data contains
    // five items, input value a, input value b, integration limits x, beta(a, b, x) and ibeta(a, b, x):
    //
+#ifndef TEST_INTERVAL
 # include "ibeta_small_data.ipp"
 
    test_inverses(ibeta_small_data);
@@ -282,6 +288,8 @@
 
    test_inverses(ibeta_large_data);
 
+#endif // TEST_INTERVAL
+
 # include "ibeta_inv_data.ipp"
 
    test_inverses2(ibeta_inv_data, name, "Inverse incomplete beta");
@@ -324,6 +332,7 @@
 {
    BOOST_MATH_CONTROL_FP;
    expected_results();
+#ifndef TEST_INTERVAL
 #ifdef TEST_GSL
    gsl_set_error_handler_off();
 #endif
@@ -363,6 +372,9 @@
       "not available at all, or because they are too inaccurate for these tests "
       "to pass.</note>" << std::cout;
 #endif
+#else // TEST_INTERVAL
+ test_beta(interval_type(), "boost::numeric::interval<double>");
+#endif // TEST_INTERVAL
    return 0;
 }
 


Boost-Commit list run by bdawes at acm.org, david.abrahams at rcn.com, gregod at cs.rpi.edu, cpdaniel at pacbell.net, john at johnmaddock.co.uk