New Distribution Testing Infastructure

include/boost/math/distributions/arcsine.hpp

@@ -229,11 +229,11 @@ namespace boost
         // Special cases:
         if (p == 0)
         {
-          return 0;
+          return x_min;
         }
         if (p == 1)
         {
-          return 1;
+          return x_max;
         }
         using boost::math::constants::half_pi;
         RealType sin2hpip = sin(half_pi<RealType>() * p);
@@ -247,11 +247,11 @@ namespace boost
         // Special cases:
         if (q == 1)
         {
-          return 0;
+          return x_min;
         }
         if (q == 0)
         {
-          return 1;
+          return x_max;
         }
         // Naive RealType p = 1 - q; result = sin(half_pi<RealType>() * p); loses accuracy, so use a cos alternative instead.
         //result = cos(half_pi<RealType>() * q); // for arcsine(0,1)

include/boost/math/distributions/exponential.hpp

@@ -109,10 +109,15 @@ template <class RealType, class Policy>
 BOOST_MATH_GPU_ENABLED inline boost::math::pair<RealType, RealType> support(const exponential_distribution<RealType, Policy>& /*dist*/)
 { // Range of supported values for random variable x.
    // This is range where cdf rises from 0 to 1, and outside it, the pdf is zero.
-   using boost::math::tools::max_value;
-   using boost::math::tools::min_value;
-   return boost::math::pair<RealType, RealType>(min_value<RealType>(),  max_value<RealType>());
-   // min_value<RealType>() to avoid a discontinuity at x = 0.
+   BOOST_MATH_IF_CONSTEXPR (boost::math::numeric_limits<RealType>::has_infinity)
+   { 
+      return boost::math::pair<RealType, RealType>(static_cast<RealType>(0), boost::math::numeric_limits<RealType>::infinity()); // - to + infinity.
+   }
+   else
+   { // Can only use max_value.
+      using boost::math::tools::max_value;
+      return boost::math::pair<RealType, RealType>(static_cast<RealType>(0), max_value<RealType>()); // - to + max.
+   }
 }
 
 template <class RealType, class Policy>
@@ -128,9 +133,6 @@ BOOST_MATH_GPU_ENABLED inline RealType pdf(const exponential_distribution<RealTy
       return result;
    if(0 == detail::verify_exp_x(function, x, &result, Policy()))
       return result;
-   // Workaround for VC11/12 bug:
-   if ((boost::math::isinf)(x))
-      return 0;
    result = lambda * exp(-lambda * x);
    return result;
 } // pdf
@@ -225,9 +227,6 @@ BOOST_MATH_GPU_ENABLED inline RealType cdf(const complemented2_type<exponential_
       return result;
    if(0 == detail::verify_exp_x(function, c.param, &result, Policy()))
       return result;
-   // Workaround for VC11/12 bug:
-   if (c.param >= tools::max_value<RealType>())
-      return 0;
    result = exp(-c.param * lambda);
 
    return result;
@@ -246,9 +245,6 @@ BOOST_MATH_GPU_ENABLED inline RealType logcdf(const complemented2_type<exponenti
       return result;
    if(0 == detail::verify_exp_x(function, c.param, &result, Policy()))
       return result;
-   // Workaround for VC11/12 bug:
-   if (c.param >= tools::max_value<RealType>())
-      return 0;
    result = -c.param * lambda;
 
    return result;

include/boost/math/distributions/fisher_f.hpp

@@ -140,10 +140,10 @@ BOOST_MATH_GPU_ENABLED inline RealType cdf(const fisher_f_distribution<RealType,
    RealType df2 = dist.degrees_of_freedom2();
    // Error check:
    RealType error_result = 0;
-   if(false == detail::check_df(
+   if(false == (detail::check_df(
          function, df1, &error_result, Policy())
          && detail::check_df(
-         function, df2, &error_result, Policy()))
+         function, df2, &error_result, Policy())))
       return error_result;
 
    if((x < 0) || !(boost::math::isfinite)(x))
@@ -201,10 +201,10 @@ BOOST_MATH_GPU_ENABLED inline RealType cdf(const complemented2_type<fisher_f_dis
    RealType x = c.param;
    // Error check:
    RealType error_result = 0;
-   if(false == detail::check_df(
+   if(false == (detail::check_df(
          function, df1, &error_result, Policy())
          && detail::check_df(
-         function, df2, &error_result, Policy()))
+         function, df2, &error_result, Policy())))
       return error_result;
 
    if((x < 0) || !(boost::math::isfinite)(x))
@@ -260,10 +260,10 @@ BOOST_MATH_GPU_ENABLED inline RealType mean(const fisher_f_distribution<RealType
    RealType df2 = dist.degrees_of_freedom2();
    // Error check:
    RealType error_result = 0;
-   if(false == detail::check_df(
+   if(false == (detail::check_df(
             function, df1, &error_result, Policy())
          && detail::check_df(
-            function, df2, &error_result, Policy()))
+            function, df2, &error_result, Policy())))
       return error_result;
    if(df2 <= 2)
    {
@@ -281,10 +281,10 @@ BOOST_MATH_GPU_ENABLED inline RealType variance(const fisher_f_distribution<Real
    RealType df2 = dist.degrees_of_freedom2();
    // Error check:
    RealType error_result = 0;
-   if(false == detail::check_df(
+   if(false == (detail::check_df(
             function, df1, &error_result, Policy())
          && detail::check_df(
-            function, df2, &error_result, Policy()))
+            function, df2, &error_result, Policy())))
       return error_result;
    if(df2 <= 4)
    {
@@ -302,10 +302,10 @@ BOOST_MATH_GPU_ENABLED inline RealType mode(const fisher_f_distribution<RealType
    RealType df2 = dist.degrees_of_freedom2();
    // Error check:
    RealType error_result = 0;
-   if(false == detail::check_df(
+   if(false == (detail::check_df(
             function, df1, &error_result, Policy())
          && detail::check_df(
-            function, df2, &error_result, Policy()))
+            function, df2, &error_result, Policy())))
       return error_result;
    if(df1 <= 2)
    {
@@ -333,10 +333,10 @@ BOOST_MATH_GPU_ENABLED inline RealType skewness(const fisher_f_distribution<Real
    RealType df2 = dist.degrees_of_freedom2();
    // Error check:
    RealType error_result = 0;
-   if(false == detail::check_df(
+   if(false == (detail::check_df(
             function, df1, &error_result, Policy())
          && detail::check_df(
-            function, df2, &error_result, Policy()))
+            function, df2, &error_result, Policy())))
       return error_result;
    if(df2 <= 6)
    {
@@ -364,10 +364,10 @@ BOOST_MATH_GPU_ENABLED inline RealType kurtosis_excess(const fisher_f_distributi
    RealType df2 = dist.degrees_of_freedom2();
    // Error check:
    RealType error_result = 0;
-   if(false == detail::check_df(
+   if(false == (detail::check_df(
             function, df1, &error_result, Policy())
          && detail::check_df(
-            function, df2, &error_result, Policy()))
+            function, df2, &error_result, Policy())))
       return error_result;
    if(df2 <= 8)
    {

include/boost/math/distributions/gamma.hpp

@@ -57,6 +57,22 @@ BOOST_MATH_GPU_ENABLED inline bool check_gamma_x(
    return true;
 }
 
+template <class RealType, class Policy>
+BOOST_MATH_GPU_ENABLED inline bool check_gamma_x_positive(
+      const char* function,
+      RealType const& x,
+      RealType* result, const Policy& pol)
+{
+   if((x < 0) || (boost::math::isnan)(x))
+   {
+      *result = policies::raise_domain_error<RealType>(
+         function,
+         "Random variate is %1% but must be >= 0 !", x, pol);
+      return false;
+   }
+   return true;
+}
+
 template <class RealType, class Policy>
 BOOST_MATH_GPU_ENABLED inline bool check_gamma(
       const char* function,
@@ -112,17 +128,30 @@ gamma_distribution(RealType,RealType)->gamma_distribution<typename boost::math::
 template <class RealType, class Policy>
 BOOST_MATH_GPU_ENABLED inline boost::math::pair<RealType, RealType> range(const gamma_distribution<RealType, Policy>& /* dist */)
 { // Range of permissible values for random variable x.
-   using boost::math::tools::max_value;
-   return boost::math::pair<RealType, RealType>(static_cast<RealType>(0), max_value<RealType>());
+  BOOST_MATH_IF_CONSTEXPR (boost::math::numeric_limits<RealType>::has_infinity)
+  { 
+      return boost::math::pair<RealType, RealType>(static_cast<RealType>(0), boost::math::numeric_limits<RealType>::infinity()); // 0 to + infinity.
+  }
+  else
+  {
+      using boost::math::tools::max_value;
+      return boost::math::pair<RealType, RealType>(static_cast<RealType>(0), max_value<RealType>()); // 0 to + max
+  }
 }
 
 template <class RealType, class Policy>
 BOOST_MATH_GPU_ENABLED inline boost::math::pair<RealType, RealType> support(const gamma_distribution<RealType, Policy>& /* dist */)
 { // Range of supported values for random variable x.
    // This is range where cdf rises from 0 to 1, and outside it, the pdf is zero.
-   using boost::math::tools::max_value;
-   using boost::math::tools::min_value;
-   return boost::math::pair<RealType, RealType>(min_value<RealType>(),  max_value<RealType>());
+  BOOST_MATH_IF_CONSTEXPR (boost::math::numeric_limits<RealType>::has_infinity)
+  { 
+      return boost::math::pair<RealType, RealType>(static_cast<RealType>(0), boost::math::numeric_limits<RealType>::infinity()); // 0 to + infinity.
+  }
+  else
+  {
+      using boost::math::tools::max_value;
+      return boost::math::pair<RealType, RealType>(static_cast<RealType>(0), max_value<RealType>()); // 0 to + max
+   }
 }
 
 template <class RealType, class Policy>
@@ -138,13 +167,18 @@ BOOST_MATH_GPU_ENABLED inline RealType pdf(const gamma_distribution<RealType, Po
    RealType result = 0;
    if(false == detail::check_gamma(function, scale, shape, &result, Policy()))
       return result;
-   if(false == detail::check_gamma_x(function, x, &result, Policy()))
+   if(false == detail::check_gamma_x_positive(function, x, &result, Policy()))
       return result;
 
    if(x == 0)
    {
       return 0;
    }
+   if ((boost::math::isinf)(x))
+   {
+      return 0;
+   }
+
    result = gamma_p_derivative(shape, x / scale, Policy()) / scale;
    return result;
 } // pdf
@@ -168,9 +202,8 @@ BOOST_MATH_GPU_ENABLED inline RealType logpdf(const gamma_distribution<RealType,
 
    if(x == 0)
    {
-      return boost::math::numeric_limits<RealType>::quiet_NaN();
+      return policies::raise_domain_error<RealType>(function, "Random variate is %1% but must be > 0 !", x, Policy());;
    }
-
    result = -k*log(theta) + (k-1)*log(x) - lgamma(k) - (x/theta);
 
    return result;
@@ -189,8 +222,17 @@ BOOST_MATH_GPU_ENABLED inline RealType cdf(const gamma_distribution<RealType, Po
    RealType result = 0;
    if(false == detail::check_gamma(function, scale, shape, &result, Policy()))
       return result;
-   if(false == detail::check_gamma_x(function, x, &result, Policy()))
+   if(false == detail::check_gamma_x_positive(function, x, &result, Policy()))
       return result;
+
+   if(x == 0)
+   {
+      return 0;
+   }
+   if ((boost::math::isinf)(x))
+   {
+      return 1;
+   }
 
    result = boost::math::gamma_p(shape, x / scale, Policy());
    return result;
@@ -233,9 +275,18 @@ BOOST_MATH_GPU_ENABLED inline RealType cdf(const complemented2_type<gamma_distri
    RealType result = 0;
    if(false == detail::check_gamma(function, scale, shape, &result, Policy()))
       return result;
-   if(false == detail::check_gamma_x(function, c.param, &result, Policy()))
+   if(false == detail::check_gamma_x_positive(function, c.param, &result, Policy()))
       return result;
 
+   if(c.param == 0)
+   {
+      return 1;
+   }
+   if ((boost::math::isinf)(c.param))
+   {
+      return 0;
+   }
+
    result = gamma_q(shape, c.param / scale, Policy());
 
    return result;

include/boost/math/distributions/geometric.hpp

@@ -145,9 +145,9 @@ namespace boost
         RealType result = 0;  // of error checks.
         RealType successes = 1;
         RealType failures = trials - successes;
-        if(false == detail::check_probability(function, alpha, &result, Policy())
+        if(false == (detail::check_probability(function, alpha, &result, Policy())
           && geometric_detail::check_dist_and_k(
-          function, RealType(0), failures, &result, Policy()))
+          function, RealType(0), failures, &result, Policy())))
         {
           return result;
         }
@@ -170,9 +170,9 @@ namespace boost
         RealType result = 0;  // of error checks.
         RealType successes = 1;
         RealType failures = trials - successes;
-        if(false == geometric_detail::check_dist_and_k(
+        if(false == (geometric_detail::check_dist_and_k(
           function, RealType(0), failures, &result, Policy())
-          && detail::check_probability(function, alpha, &result, Policy()))
+          && detail::check_probability(function, alpha, &result, Policy())))
         {
           return result;
         }
@@ -202,9 +202,9 @@ namespace boost
         constexpr auto function = "boost::math::geometric<%1%>::find_minimum_number_of_trials";
         // Error checks:
         RealType result = 0;
-        if(false == geometric_detail::check_dist_and_k(
+        if(false == (geometric_detail::check_dist_and_k(
           function, p, k, &result, Policy())
-          && detail::check_probability(function, alpha, &result, Policy()))
+          && detail::check_probability(function, alpha, &result, Policy())))
         {
           return result;
         }
@@ -220,9 +220,9 @@ namespace boost
         constexpr auto function = "boost::math::geometric<%1%>::find_maximum_number_of_trials";
         // Error checks:
         RealType result = 0;
-        if(false == geometric_detail::check_dist_and_k(
+        if(false == (geometric_detail::check_dist_and_k(
           function, p, k, &result, Policy())
-          &&  detail::check_probability(function, alpha, &result, Policy()))
+          &&  detail::check_probability(function, alpha, &result, Policy())))
         {
           return result;
         }
@@ -396,7 +396,7 @@ namespace boost
         k,
         &result, Policy()))
       {
-        return -boost::math::numeric_limits<RealType>::infinity();
+        return result;
       }
       if(k == 0)
       {
@@ -452,7 +452,7 @@ namespace boost
         k,
         &result, Policy()))
       {
-        return -boost::math::numeric_limits<RealType>::infinity();
+        return result;
       }
 
       return boost::math::log1p(-p, Policy()) * (k+1);
@@ -497,10 +497,6 @@ namespace boost
       { // p <= pdf(dist, 0) == cdf(dist, 0)
         return 0;
       }
-      if (x == 1)
-      {
-        return 0;
-      }
 
       // log(1-x) /log(1-success_fraction) -1; but use log1p in case success_fraction is small
       result = boost::math::log1p(-x, Policy()) / boost::math::log1p(-success_fraction, Policy()) - 1;