[Boost-commit] svn:boost r77263 - sandbox/math/libs/math/test

Subject: [Boost-commit] svn:boost r77263 - sandbox/math/libs/math/test
From: ben_at_[hidden]
Date: 2012-03-08 08:35:54

Author: bensobotta
Date: 2012-03-08 08:35:53 EST (Thu, 08 Mar 2012)
New Revision: 77263
URL: http://svn.boost.org/trac/boost/changeset/77263

Log:
Modified T7 verification routine to avoid underflows. Fixed speed issues
with the tests.

Text files modified: 
   sandbox/math/libs/math/test/owens_t_T7.hpp   |    21 +++++++++++++--------                   
   sandbox/math/libs/math/test/test_owens_t.cpp |     4 ++--                                    
   2 files changed, 15 insertions(+), 10 deletions(-)

Modified: sandbox/math/libs/math/test/owens_t_T7.hpp
==============================================================================
--- sandbox/math/libs/math/test/owens_t_T7.hpp	(original)
+++ sandbox/math/libs/math/test/owens_t_T7.hpp	2012-03-08 08:35:53 EST (Thu, 08 Mar 2012)
@@ -95,22 +95,22 @@
         const RealType hs = h*h;
         const RealType as = a*a;
 
-	RealType u = one_div_two_pi<RealType>() * a * exp(-half<RealType>()*hs*(static_cast<RealType>(1)+as));
+	//RealType u = one_div_two_pi<RealType>() * a * exp(-half<RealType>()*hs*(static_cast<RealType>(1)+as));
         RealType v = c2[0];
 
-	RealType val = u*v;
+	RealType val = v;
         RealType last_val = val+1; // last_val must not be the same as val
 
-	unsigned k = 0;
+	int k = 0;
 
         std::vector<RealType> memory;
-	memory.push_back(u*v);
+	memory.push_back(v);
 
-	while((val != last_val) || (k<c2.size())) // use all c2
+	while((val != last_val) || (k<51)) // use all c2
           {
             last_val = val;
             k++;
-	    u *= as;
+	    const RealType u = std::pow(as,k);
             if(k < c2.size())
               {
                 v = (hs*v + c2[k])/(static_cast<RealType>(2*k+1));
@@ -133,12 +133,17 @@
         std::sort(ndx_4_sorted_data.begin(), ndx_4_sorted_data.end(),
                   boost::bind(owens_t_sort_proxy<typename std::vector<RealType>::size_type, RealType>, _1, _2, &memory[0]));
 
-	val = static_cast<RealType>(0);
-	for(unsigned i = 0; i != memory.size(); i++)
+	val = memory[ndx_4_sorted_data[0]];
+	for(unsigned i = 1; i != memory.size(); i++)
           {
             val+=memory[ndx_4_sorted_data[i]];
           }
 
+        // split the exponential to avoid values that go below the minimum floating pt value
+        val *= exp(-half<RealType>()*hs*as);
+        val *= exp(-half<RealType>()*hs);
+        val *= one_div_two_pi<RealType>() * a;
+        
         return val;
       } // RealType compute_owens_t_T7(const RealType h, const RealType a)
       

Modified: sandbox/math/libs/math/test/test_owens_t.cpp
==============================================================================
--- sandbox/math/libs/math/test/test_owens_t.cpp	(original)
+++ sandbox/math/libs/math/test/test_owens_t.cpp	2012-03-08 08:35:53 EST (Thu, 08 Mar 2012)
@@ -159,8 +159,8 @@
   using namespace std; // ADL of std names.
 
   // apply log scale because points near zero are more interesting
-  for(RealType a = static_cast<RealType>(-10.0l); a < static_cast<RealType>(3l); a+= static_cast<RealType>(0.1l))
-    for(RealType h = static_cast<RealType>(-10.0l); h < static_cast<RealType>(4l); h+= static_cast<RealType>(0.1l))
+  for(RealType a = static_cast<RealType>(-10.0l); a < static_cast<RealType>(3l); a+= static_cast<RealType>(0.2l))
+    for(RealType h = static_cast<RealType>(-10.0l); h < static_cast<RealType>(4l); h+= static_cast<RealType>(0.2l))
     {
       const RealType expa = exp(a);
       const RealType exph = exp(h);