summaryrefslogtreecommitdiff
path: root/boost/math/special_functions/detail/erf_inv.hpp
diff options
context:
space:
mode:
Diffstat (limited to 'boost/math/special_functions/detail/erf_inv.hpp')
-rw-r--r--boost/math/special_functions/detail/erf_inv.hpp51
1 files changed, 38 insertions, 13 deletions
diff --git a/boost/math/special_functions/detail/erf_inv.hpp b/boost/math/special_functions/detail/erf_inv.hpp
index d51db9d52f..77aa72fc26 100644
--- a/boost/math/special_functions/detail/erf_inv.hpp
+++ b/boost/math/special_functions/detail/erf_inv.hpp
@@ -50,7 +50,7 @@ T erf_inv_imp(const T& p, const T& q, const Policy&, const boost::mpl::int_<64>*
BOOST_MATH_BIG_CONSTANT(T, 64, -0.00538772965071242932965)
};
static const T Q[] = {
- 1,
+ BOOST_MATH_BIG_CONSTANT(T, 64, 1.0),
BOOST_MATH_BIG_CONSTANT(T, 64, -0.970005043303290640362),
BOOST_MATH_BIG_CONSTANT(T, 64, -1.56574558234175846809),
BOOST_MATH_BIG_CONSTANT(T, 64, 1.56221558398423026363),
@@ -92,7 +92,7 @@ T erf_inv_imp(const T& p, const T& q, const Policy&, const boost::mpl::int_<64>*
BOOST_MATH_BIG_CONSTANT(T, 64, -3.67192254707729348546)
};
static const T Q[] = {
- BOOST_MATH_BIG_CONSTANT(T, 64, 1),
+ BOOST_MATH_BIG_CONSTANT(T, 64, 1.0),
BOOST_MATH_BIG_CONSTANT(T, 64, 6.24264124854247537712),
BOOST_MATH_BIG_CONSTANT(T, 64, 3.9713437953343869095),
BOOST_MATH_BIG_CONSTANT(T, 64, -28.6608180499800029974),
@@ -147,7 +147,7 @@ T erf_inv_imp(const T& p, const T& q, const Policy&, const boost::mpl::int_<64>*
BOOST_MATH_BIG_CONSTANT(T, 64, -0.681149956853776992068e-9)
};
static const T Q[] = {
- 1,
+ BOOST_MATH_BIG_CONSTANT(T, 64, 1.0),
BOOST_MATH_BIG_CONSTANT(T, 64, 3.46625407242567245975),
BOOST_MATH_BIG_CONSTANT(T, 64, 5.38168345707006855425),
BOOST_MATH_BIG_CONSTANT(T, 64, 4.77846592945843778382),
@@ -176,7 +176,7 @@ T erf_inv_imp(const T& p, const T& q, const Policy&, const boost::mpl::int_<64>*
BOOST_MATH_BIG_CONSTANT(T, 64, 0.266339227425782031962e-11)
};
static const T Q[] = {
- BOOST_MATH_BIG_CONSTANT(T, 64, 1),
+ BOOST_MATH_BIG_CONSTANT(T, 64, 1.0),
BOOST_MATH_BIG_CONSTANT(T, 64, 1.3653349817554063097),
BOOST_MATH_BIG_CONSTANT(T, 64, 0.762059164553623404043),
BOOST_MATH_BIG_CONSTANT(T, 64, 0.220091105764131249824),
@@ -204,7 +204,7 @@ T erf_inv_imp(const T& p, const T& q, const Policy&, const boost::mpl::int_<64>*
BOOST_MATH_BIG_CONSTANT(T, 64, 0.99055709973310326855e-16)
};
static const T Q[] = {
- BOOST_MATH_BIG_CONSTANT(T, 64, 1),
+ BOOST_MATH_BIG_CONSTANT(T, 64, 1.0),
BOOST_MATH_BIG_CONSTANT(T, 64, 0.591429344886417493481),
BOOST_MATH_BIG_CONSTANT(T, 64, 0.138151865749083321638),
BOOST_MATH_BIG_CONSTANT(T, 64, 0.0160746087093676504695),
@@ -231,7 +231,7 @@ T erf_inv_imp(const T& p, const T& q, const Policy&, const boost::mpl::int_<64>*
BOOST_MATH_BIG_CONSTANT(T, 64, -0.116765012397184275695e-17)
};
static const T Q[] = {
- BOOST_MATH_BIG_CONSTANT(T, 64, 1),
+ BOOST_MATH_BIG_CONSTANT(T, 64, 1.0),
BOOST_MATH_BIG_CONSTANT(T, 64, 0.207123112214422517181),
BOOST_MATH_BIG_CONSTANT(T, 64, 0.0169410838120975906478),
BOOST_MATH_BIG_CONSTANT(T, 64, 0.000690538265622684595676),
@@ -258,7 +258,7 @@ T erf_inv_imp(const T& p, const T& q, const Policy&, const boost::mpl::int_<64>*
BOOST_MATH_BIG_CONSTANT(T, 64, -0.348890393399948882918e-21)
};
static const T Q[] = {
- BOOST_MATH_BIG_CONSTANT(T, 64, 1),
+ BOOST_MATH_BIG_CONSTANT(T, 64, 1.0),
BOOST_MATH_BIG_CONSTANT(T, 64, 0.0845746234001899436914),
BOOST_MATH_BIG_CONSTANT(T, 64, 0.00282092984726264681981),
BOOST_MATH_BIG_CONSTANT(T, 64, 0.468292921940894236786e-4),
@@ -282,7 +282,7 @@ struct erf_roots
BOOST_MATH_STD_USING
T derivative = sign * (2 / sqrt(constants::pi<T>())) * exp(-(guess * guess));
T derivative2 = -2 * guess * derivative;
- return boost::math::make_tuple(((sign > 0) ? boost::math::erf(guess, Policy()) : boost::math::erfc(guess, Policy())) - target, derivative, derivative2);
+ return boost::math::make_tuple(((sign > 0) ? static_cast<T>(boost::math::erf(guess, Policy()) - target) : static_cast<T>(boost::math::erfc(guess, Policy())) - target), derivative, derivative2);
}
erf_roots(T z, int s) : target(z), sign(s) {}
private:
@@ -331,18 +331,34 @@ struct erf_inv_initializer
{
do_init();
}
+ static bool is_value_non_zero(T);
static void do_init()
{
boost::math::erf_inv(static_cast<T>(0.25), Policy());
boost::math::erf_inv(static_cast<T>(0.55), Policy());
boost::math::erf_inv(static_cast<T>(0.95), Policy());
boost::math::erfc_inv(static_cast<T>(1e-15), Policy());
- if(static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 1e-130)) != 0)
+ // These following initializations must not be called if
+ // type T can not hold the relevant values without
+ // underflow to zero. We check this at runtime because
+ // some tools such as valgrind silently change the precision
+ // of T at runtime, and numeric_limits basically lies!
+ if(is_value_non_zero(static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 1e-130))))
boost::math::erfc_inv(static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 1e-130)), Policy());
- if(static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 1e-800)) != 0)
+
+ // Some compilers choke on constants that would underflow, even in code that isn't instantiated
+ // so try and filter these cases out in the preprocessor:
+#if LDBL_MAX_10_EXP >= 800
+ if(is_value_non_zero(static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 1e-800))))
boost::math::erfc_inv(static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 1e-800)), Policy());
- if(static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 1e-900)) != 0)
+ if(is_value_non_zero(static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 1e-900))))
boost::math::erfc_inv(static_cast<T>(BOOST_MATH_BIG_CONSTANT(T, 64, 1e-900)), Policy());
+#else
+ if(is_value_non_zero(static_cast<T>(BOOST_MATH_HUGE_CONSTANT(T, 64, 1e-800))))
+ boost::math::erfc_inv(static_cast<T>(BOOST_MATH_HUGE_CONSTANT(T, 64, 1e-800)), Policy());
+ if(is_value_non_zero(static_cast<T>(BOOST_MATH_HUGE_CONSTANT(T, 64, 1e-900))))
+ boost::math::erfc_inv(static_cast<T>(BOOST_MATH_HUGE_CONSTANT(T, 64, 1e-900)), Policy());
+#endif
}
void force_instantiate()const{}
};
@@ -356,6 +372,15 @@ struct erf_inv_initializer
template <class T, class Policy>
const typename erf_inv_initializer<T, Policy>::init erf_inv_initializer<T, Policy>::initializer;
+template <class T, class Policy>
+bool erf_inv_initializer<T, Policy>::init::is_value_non_zero(T v)
+{
+ // This needs to be non-inline to detect whether v is non zero at runtime
+ // rather than at compile time, only relevant when running under valgrind
+ // which changes long double's to double's on the fly.
+ return v != 0;
+}
+
} // namespace detail
template <class T, class Policy>
@@ -368,7 +393,7 @@ typename tools::promote_args<T>::type erfc_inv(T z, const Policy& pol)
//
static const char* function = "boost::math::erfc_inv<%1%>(%1%, %1%)";
if((z < 0) || (z > 2))
- policies::raise_domain_error<result_type>(function, "Argument outside range [0,2] in inverse erfc function (got p=%1%).", z, pol);
+ return policies::raise_domain_error<result_type>(function, "Argument outside range [0,2] in inverse erfc function (got p=%1%).", z, pol);
if(z == 0)
return policies::raise_overflow_error<result_type>(function, 0, pol);
if(z == 2)
@@ -432,7 +457,7 @@ typename tools::promote_args<T>::type erf_inv(T z, const Policy& pol)
//
static const char* function = "boost::math::erf_inv<%1%>(%1%, %1%)";
if((z < -1) || (z > 1))
- policies::raise_domain_error<result_type>(function, "Argument outside range [-1, 1] in inverse erf function (got p=%1%).", z, pol);
+ return policies::raise_domain_error<result_type>(function, "Argument outside range [-1, 1] in inverse erf function (got p=%1%).", z, pol);
if(z == 1)
return policies::raise_overflow_error<result_type>(function, 0, pol);
if(z == -1)