boost/math
2
0
Fork 0
mirror of https://github.com/boostorg/math.git synced 2026-01-28 19:32:08 +00:00
Code Issues Projects Releases Wiki Activity

Modified T7 verification routine to avoid underflows. Fixed speed issues

Browse Source 
with the tests.


[SVN r77263]
This commit is contained in:
Benjamin Sobotta
2012-03-08 13:35:53 +00:00
parent d21534ff06
commit 5dd0c19d85
2 changed files with 15 additions and 10 deletions
Download Patch File Download Diff File Expand all files Collapse all files

21
test/owens_t_T7.hpp
View File

@@ -95,22 +95,22 @@ namespace boost
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 @@ namespace boost
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)

4
test/test_owens_t.cpp
View File

@@ -159,8 +159,8 @@ void check_against_T7(RealType)
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);