[Boost-commit] svn:boost r65444 - in trunk: boost/math/distributions libs/math/test

Subject: [Boost-commit] svn:boost r65444 - in trunk: boost/math/distributions libs/math/test
From: john_at_[hidden]
Date: 2010-09-17 14:13:48

Author: johnmaddock
Date: 2010-09-17 14:13:45 EDT (Fri, 17 Sep 2010)
New Revision: 65444
URL: http://svn.boost.org/trac/boost/changeset/65444

Log:
Improve error handling and fix Linux errors.
Text files modified: 
   trunk/boost/math/distributions/inverse_gamma.hpp         |    30 ++++++++++++++++++++++++++----          
   trunk/libs/math/test/test_inverse_gamma_distribution.cpp |    30 +++++++++++++++++-------------          
   2 files changed, 43 insertions(+), 17 deletions(-)

Modified: trunk/boost/math/distributions/inverse_gamma.hpp
==============================================================================
--- trunk/boost/math/distributions/inverse_gamma.hpp	(original)
+++ trunk/boost/math/distributions/inverse_gamma.hpp	2010-09-17 14:13:45 EDT (Fri, 17 Sep 2010)
@@ -1,7 +1,7 @@
 // inverse_gamma.hpp
 
-//  Copyright John Maddock 2006 - 2010.
 //  Copyright Paul A. Bristow 2010.
+//  Copyright John Maddock 2010.
 //  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)
@@ -164,7 +164,23 @@
    { // x bad.
       return result;
    }
-   result = gamma_p_derivative(shape, scale / x, Policy()) * scale / (x * x);
+   result = gamma_p_derivative(shape, scale / x, Policy()) * scale;
+   if(0 != result)
+   {
+      if(x < 0)
+      {
+         // x * x may under or overflow, likewise our result,
+         // so be extra careful about the arithmetic:
+         RealType lim = tools::max_value<RealType>() * x;
+         if(lim < result)
+            return policies::raise_overflow_error<RealType, Policy>(function, "PDF is infinite.", Policy());
+         result /= x;
+         if(lim < result)
+            return policies::raise_overflow_error<RealType, Policy>(function, "PDF is infinite.", Policy());
+         result /= x;
+      }
+      result /= (x * x);
+   }
    // better than naive
    // result = (pow(scale, shape) * pow(x, (-shape -1)) * exp(-scale/x) ) / tgamma(shape);
    return result;
@@ -217,7 +233,10 @@
    {
       return policies::raise_overflow_error<RealType>(function, 0, Policy());
    }
-   result = scale / gamma_q_inv(shape, p, Policy());
+   result = gamma_q_inv(shape, p, Policy());
+   if((result < 1) && (result * tools::max_value<RealType>() < scale))
+      return policies::raise_overflow_error<RealType, Policy>(function, "Value of random variable in inverse gamma distribution quantile is infinite.", Policy());
+   result = scale / result;
    return result;
 }
 
@@ -263,7 +282,10 @@
    {
       return policies::raise_overflow_error<RealType>(function, 0, Policy());
    }
-   result = scale / gamma_p_inv(shape, q, Policy());
+   result = gamma_p_inv(shape, q, Policy());
+   if((result < 1) && (result * tools::max_value<RealType>() < scale))
+      return policies::raise_overflow_error<RealType, Policy>(function, "Value of random variable in inverse gamma distribution quantile is infinite.", Policy());
+   result = scale / result;
    return result;
 }
 

Modified: trunk/libs/math/test/test_inverse_gamma_distribution.cpp
==============================================================================
--- trunk/libs/math/test/test_inverse_gamma_distribution.cpp	(original)
+++ trunk/libs/math/test/test_inverse_gamma_distribution.cpp	2010-09-17 14:13:45 EDT (Fri, 17 Sep 2010)
@@ -228,20 +228,24 @@
 
   // Special and limit cases:
 
-  RealType mx = (std::numeric_limits<RealType>::max)();
-  RealType mi = (std::numeric_limits<RealType>::min)();
+  if(std::numeric_limits<RealType>::is_specialized)
+  {
+    RealType mx = (std::numeric_limits<RealType>::max)();
+    RealType mi = (std::numeric_limits<RealType>::min)();
+
+     BOOST_CHECK_EQUAL(
+     pdf(inverse_gamma_distribution<RealType>(1),
+       static_cast<RealType>(mx)), // max()
+       static_cast<RealType>(0)
+       );
+
+     BOOST_CHECK_EQUAL(
+     pdf(inverse_gamma_distribution<RealType>(1),
+       static_cast<RealType>(mi)), // min()
+       static_cast<RealType>(0)
+       );
 
-  BOOST_CHECK_EQUAL(
-  pdf(inverse_gamma_distribution<RealType>(1),
-    static_cast<RealType>(mx)), // max()
-    static_cast<RealType>(0)
-    );
-
-  BOOST_CHECK_EQUAL(
-  pdf(inverse_gamma_distribution<RealType>(1),
-    static_cast<RealType>(mi)), // min()
-    static_cast<RealType>(0)
-    );
+  }
 
   BOOST_CHECK_EQUAL(
     pdf(inverse_gamma_distribution<RealType>(1), static_cast<RealType>(0)), static_cast<RealType>(0));