[Boost-commit] svn:boost r83228 - in trunk: boost/math/constants boost/math/special_functions libs/math/test

Subject: [Boost-commit] svn:boost r83228 - in trunk: boost/math/constants boost/math/special_functions libs/math/test
From: john_at_[hidden]
Date: 2013-03-01 08:18:08

Author: johnmaddock
Date: 2013-03-01 08:18:07 EST (Fri, 01 Mar 2013)
New Revision: 83228
URL: http://svn.boost.org/trac/boost/changeset/83228

Log:
Fix bug in zeta that results in an incorrect domain_error for -eps < s < 0.
Add new constant needed by fix.
Text files modified: 
   trunk/boost/math/constants/calculate_constants.hpp |     8 ++++++++                                
   trunk/boost/math/constants/constants.hpp           |     1 +                                       
   trunk/boost/math/constants/generate.hpp            |     4 ++--                                    
   trunk/boost/math/special_functions/zeta.hpp        |     6 +++---                                  
   trunk/libs/math/test/test_constant_generate.cpp    |     1 +                                       
   trunk/libs/math/test/test_zeta.hpp                 |    32 ++++++++++++++++++++++++++++++++        
   6 files changed, 47 insertions(+), 5 deletions(-)

Modified: trunk/boost/math/constants/calculate_constants.hpp
==============================================================================
--- trunk/boost/math/constants/calculate_constants.hpp	(original)
+++ trunk/boost/math/constants/calculate_constants.hpp	2013-03-01 08:18:07 EST (Fri, 01 Mar 2013)
@@ -122,6 +122,14 @@
 
 template <class T>
 template<int N>
+inline T constant_log_root_two_pi<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
+{
+   BOOST_MATH_STD_USING
+   return log(root_two_pi<T, policies::policy<policies::digits2<N> > >());
+}
+
+template <class T>
+template<int N>
 inline T constant_root_ln_four<T>::compute(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE_SPEC(mpl::int_<N>))
 {
    BOOST_MATH_STD_USING

Modified: trunk/boost/math/constants/constants.hpp
==============================================================================
--- trunk/boost/math/constants/constants.hpp	(original)
+++ trunk/boost/math/constants/constants.hpp	2013-03-01 08:18:07 EST (Fri, 01 Mar 2013)
@@ -247,6 +247,7 @@
   BOOST_DEFINE_MATH_CONSTANT(root_pi, 1.772453850905516027298167483341145182e+00, "1.77245385090551602729816748334114518279754945612238712821380778985291128459103218137495065673854466541622682362e+00")
   BOOST_DEFINE_MATH_CONSTANT(root_half_pi, 1.253314137315500251207882642405522626e+00, "1.25331413731550025120788264240552262650349337030496915831496178817114682730392098747329791918902863305800498633e+00")
   BOOST_DEFINE_MATH_CONSTANT(root_two_pi, 2.506628274631000502415765284811045253e+00, "2.50662827463100050241576528481104525300698674060993831662992357634229365460784197494659583837805726611600997267e+00")
+  BOOST_DEFINE_MATH_CONSTANT(log_root_two_pi, 9.189385332046727417803297364056176398e-01, "9.18938533204672741780329736405617639861397473637783412817151540482765695927260397694743298635954197622005646625e-01")
   BOOST_DEFINE_MATH_CONSTANT(one_div_root_pi, 5.641895835477562869480794515607725858e-01, "5.64189583547756286948079451560772585844050629328998856844085721710642468441493414486743660202107363443028347906e-01")
   BOOST_DEFINE_MATH_CONSTANT(root_one_div_pi, 5.641895835477562869480794515607725858e-01, "5.64189583547756286948079451560772585844050629328998856844085721710642468441493414486743660202107363443028347906e-01")
   BOOST_DEFINE_MATH_CONSTANT(pi_minus_three, 1.415926535897932384626433832795028841e-01, "1.41592653589793238462643383279502884197169399375105820974944592307816406286208998628034825342117067982148086513e-01")

Modified: trunk/boost/math/constants/generate.hpp
==============================================================================
--- trunk/boost/math/constants/generate.hpp	(original)
+++ trunk/boost/math/constants/generate.hpp	2013-03-01 08:18:07 EST (Fri, 01 Mar 2013)
@@ -17,7 +17,7 @@
 #elif defined(USE_MPREAL)
 #include <boost/math/bindings/mpreal.hpp>
 #elif defined(USE_CPP_FLOAT)
-#include <boost/multiprecision/cpp_float.hpp>
+#include <boost/multiprecision/cpp_dec_float.hpp>
 #else
 #include <boost/math/bindings/rr.hpp>
 #endif
@@ -29,7 +29,7 @@
 #elif defined(USE_MPREAL)
 typedef mpfr::mpreal generator_type;
 #elif defined(USE_CPP_FLOAT)
-typedef boost::multiprecision::mp_number<boost::multiprecision::cpp_float<500> > generator_type;
+typedef boost::multiprecision::number<boost::multiprecision::cpp_dec_float<500> > generator_type;
 #else
 typedef ntl::RR generator_type;
 #endif

Modified: trunk/boost/math/special_functions/zeta.hpp
==============================================================================
--- trunk/boost/math/special_functions/zeta.hpp	(original)
+++ trunk/boost/math/special_functions/zeta.hpp	2013-03-01 08:18:07 EST (Fri, 01 Mar 2013)
@@ -866,15 +866,15 @@
 T zeta_imp(T s, T sc, const Policy& pol, const Tag& tag)
 {
    BOOST_MATH_STD_USING
-   if(s == 1)
+   if(sc == 0)
       return policies::raise_pole_error<T>(
          "boost::math::zeta<%1%>", 
          "Evaluation of zeta function at pole %1%", 
          s, pol);
    T result;
-   if(s == 0)
+   if(fabs(s) < tools::root_epsilon<T>())
    {
-      result = -0.5;
+      result = -0.5f - constants::log_root_two_pi<T, Policy>() * s;
    }
    else if(s < 0)
    {

Modified: trunk/libs/math/test/test_constant_generate.cpp
==============================================================================
--- trunk/libs/math/test/test_constant_generate.cpp	(original)
+++ trunk/libs/math/test/test_constant_generate.cpp	2013-03-01 08:18:07 EST (Fri, 01 Mar 2013)
@@ -50,6 +50,7 @@
    BOOST_CONSTANTS_GENERATE(root_pi);
    BOOST_CONSTANTS_GENERATE(root_half_pi);
    BOOST_CONSTANTS_GENERATE(root_two_pi);
+   BOOST_CONSTANTS_GENERATE(log_root_two_pi);
    BOOST_CONSTANTS_GENERATE(two_div_pi);
    BOOST_CONSTANTS_GENERATE(root_two_div_pi);
 

Modified: trunk/libs/math/test/test_zeta.hpp
==============================================================================
--- trunk/libs/math/test/test_zeta.hpp	(original)
+++ trunk/libs/math/test/test_zeta.hpp	2013-03-01 08:18:07 EST (Fri, 01 Mar 2013)
@@ -118,5 +118,37 @@
    BOOST_CHECK_CLOSE(::boost::math::zeta(static_cast<T>(-20)), static_cast<T>(0), tolerance * 100);
    BOOST_CHECK_CLOSE(::boost::math::zeta(static_cast<T>(-21)), static_cast<T>(-281.46014492753623188405797101449275362318840579710L), tolerance * 100);
    BOOST_CHECK_CLOSE(::boost::math::zeta(static_cast<T>(-30.125)), static_cast<T>(2.2762941726834511267740045451463455513839970804578e7L), tolerance * 100);
+   // Very small values:
+   BOOST_CHECK_CLOSE(::boost::math::zeta(ldexp(static_cast<T>(1), -20)), static_cast<T>(-0.500000876368989859479646132126454890645615288202492097957612L), tolerance);
+   BOOST_CHECK_CLOSE(::boost::math::zeta(ldexp(static_cast<T>(1), -21)), static_cast<T>(-0.500000438184266833093492063649184012943132422189989164545507L), tolerance);
+   BOOST_CHECK_CLOSE(::boost::math::zeta(ldexp(static_cast<T>(1), -22)), static_cast<T>(-0.500000219092076392425852854644256723571669269957526445270374L), tolerance);
+   BOOST_CHECK_CLOSE(::boost::math::zeta(ldexp(static_cast<T>(1), -23)), static_cast<T>(-0.500000109546023940187789325464529558825433290921168958481804L), tolerance);
+   BOOST_CHECK_CLOSE(::boost::math::zeta(ldexp(static_cast<T>(1), -24)), static_cast<T>(-0.500000054773008406088246161057525197302821575823476487961574L), tolerance);
+   BOOST_CHECK_CLOSE(::boost::math::zeta(ldexp(static_cast<T>(1), -25)), static_cast<T>(-0.500000027386503312042790426817221131071450407798601059264341L), tolerance);
+   BOOST_CHECK_CLOSE(::boost::math::zeta(ldexp(static_cast<T>(1), -26)), static_cast<T>(-0.500000013693251433271071983943082871935521396740331377486886L), tolerance);
+   BOOST_CHECK_CLOSE(::boost::math::zeta(ldexp(static_cast<T>(1), -27)), static_cast<T>(-0.500000006846625660947956426350389518286874288247134329498289L), tolerance);
+   BOOST_CHECK_CLOSE(::boost::math::zeta(ldexp(static_cast<T>(1), -28)), static_cast<T>(-0.500000003423312816552083476988056486473169377162409806781384L), tolerance);
+   BOOST_CHECK_CLOSE(::boost::math::zeta(ldexp(static_cast<T>(1), -29)), static_cast<T>(-0.500000001711656404795568073849512135664960180586820144333542L), tolerance);
+   BOOST_CHECK_CLOSE(::boost::math::zeta(ldexp(static_cast<T>(1), -30)), static_cast<T>(-0.500000000855828201527665623188910582717329375986726355164261L), tolerance);
+   BOOST_CHECK_CLOSE(::boost::math::zeta(ldexp(static_cast<T>(1), -31)), static_cast<T>(-0.500000000427914100546303208463654361814800355929815322493143L), tolerance);
+   BOOST_CHECK_CLOSE(::boost::math::zeta(ldexp(static_cast<T>(1), -32)), static_cast<T>(-0.500000000213957050218769203487022003676593508474107873788445L), tolerance);
+   BOOST_CHECK_CLOSE(::boost::math::zeta(ldexp(static_cast<T>(1), -33)), static_cast<T>(-0.500000000106978525095789001562046589421133388262409441738089L), tolerance);
+   BOOST_CHECK_CLOSE(::boost::math::zeta(ldexp(static_cast<T>(1), -34)), static_cast<T>(-0.500000000053489262544495600736249301842352101231724731340202L), tolerance);
+
+   BOOST_CHECK_CLOSE(::boost::math::zeta(-ldexp(static_cast<T>(1), -20)), static_cast<T>(-0.499999123632834911086872289657767335473025908373776645822722L), tolerance);
+   BOOST_CHECK_CLOSE(::boost::math::zeta(-ldexp(static_cast<T>(1), -21)), static_cast<T>(-0.499999561816189359548137231641582253243376087534976981434190L), tolerance);
+   BOOST_CHECK_CLOSE(::boost::math::zeta(-ldexp(static_cast<T>(1), -22)), static_cast<T>(-0.499999780908037655734554449793729262345041281451929584703788L), tolerance);
+   BOOST_CHECK_CLOSE(::boost::math::zeta(-ldexp(static_cast<T>(1), -23)), static_cast<T>(-0.499999890454004571852312499433422838864632848598847415933664L), tolerance);
+   BOOST_CHECK_CLOSE(::boost::math::zeta(-ldexp(static_cast<T>(1), -24)), static_cast<T>(-0.499999945226998721921779295091241395945379526155584220813497L), tolerance);
+   BOOST_CHECK_CLOSE(::boost::math::zeta(-ldexp(static_cast<T>(1), -25)), static_cast<T>(-0.499999972613498469959715937215237923104705216198368099221577L), tolerance);
+   BOOST_CHECK_CLOSE(::boost::math::zeta(-ldexp(static_cast<T>(1), -26)), static_cast<T>(-0.499999986306749012229554607064736104475024094525587925697276L), tolerance);
+   BOOST_CHECK_CLOSE(::boost::math::zeta(-ldexp(static_cast<T>(1), -27)), static_cast<T>(-0.499999993153374450427200221401546739119918746163907954406855L), tolerance);
+   BOOST_CHECK_CLOSE(::boost::math::zeta(-ldexp(static_cast<T>(1), -28)), static_cast<T>(-0.499999996576687211291705684949926422460038672790251466963619L), tolerance);
+   BOOST_CHECK_CLOSE(::boost::math::zeta(-ldexp(static_cast<T>(1), -29)), static_cast<T>(-0.499999998288343602165379216634983519354686193860717726606017L), tolerance);
+   BOOST_CHECK_CLOSE(::boost::math::zeta(-ldexp(static_cast<T>(1), -30)), static_cast<T>(-0.499999999144171800212571199432213326524228740247618955829902L), tolerance);
+   BOOST_CHECK_CLOSE(::boost::math::zeta(-ldexp(static_cast<T>(1), -31)), static_cast<T>(-0.499999999572085899888755997191626615213504580792674808876724L), tolerance);
+   BOOST_CHECK_CLOSE(::boost::math::zeta(-ldexp(static_cast<T>(1), -32)), static_cast<T>(-0.499999999786042949889995597926798240562852438685508646794693L), tolerance);
+   BOOST_CHECK_CLOSE(::boost::math::zeta(-ldexp(static_cast<T>(1), -33)), static_cast<T>(-0.499999999893021474931402198791408471637626205588681812641711L), tolerance);
+   BOOST_CHECK_CLOSE(::boost::math::zeta(-ldexp(static_cast<T>(1), -34)), static_cast<T>(-0.499999999946510737462302199352114463422268928922372277519378L), tolerance);
 }