root/Math/Math.h

0001 // @(#)root/mathcore:$Id$
0002 // Author: L. Moneta Tue Nov 14 15:44:38 2006
0003
0004 /**********************************************************************
0005  *                                                                    *
0006  * Copyright (c) 2006  LCG ROOT Math Team, CERN/PH-SFT                *
0007  *                                                                    *
0008  *                                                                    *
0009  **********************************************************************/
0010
0011 // mathematical constants like Pi
0012
0013 #ifndef ROOT_Math_Math
0014 #define ROOT_Math_Math
0015
0016 #ifdef _MSC_VER
0017 #define _USE_MATH_DEFINES
0018 #endif
0019
0020 #include <cmath>
0021
0022 #if defined(__sun) || defined(_MSC_VER)
0023 //Microsoft and solaris definition of cmath does not include math.h which has the definitions of numerical constants
0024 #include <math.h> // for M_PI
0025 // TODO replace with std::numbers::pi once minimum version is C++20, and remove this code block
0026 #endif
0027
0028
0029 #ifdef HAVE_NO_EXPM1
0030 // needed to implement expm1
0031 #include <limits>
0032 #endif
0033
0034
0035 #ifndef M_PI
0036
0037 #define M_PI       3.14159265358979323846264338328      // Pi
0038 #endif
0039
0040 #ifndef M_PI_2
0041 #define M_PI_2     1.57079632679489661923132169164      // Pi/2
0042 #endif
0043
0044 #ifndef M_PI_4
0045 #define M_PI_4     0.78539816339744830961566084582      // Pi/4
0046 #endif
0047
0048 namespace ROOT {
0049
0050 /**
0051 \namespace Math
0052 Namespace for new Math classes and functions.
0053 See the \ref Math "Math Libraries" page for a detailed description.
0054 */
0055
0056 namespace Math {
0057
0058 /**
0059     Mathematical constants
0060 */
0061 inline double Pi()
0062 {
0063    return M_PI;
0064    }
0065
0066    /**
0067        declarations for functions which are not implemented by some compilers
0068    */
0069
0070    /// log(1+x) with error cancelation when x is small
0071    inline double log1p(double x)
0072    {
0073 #ifndef HAVE_NO_LOG1P
0074    return ::log1p(x);
0075 #else
0076    // if log1p is not in c math library
0077   volatile double y;
0078   y = 1 + x;
0079   return std::log(y) - ((y-1)-x)/y ;  /* cancels errors with IEEE arithmetic */
0080 #endif
0081 }
0082 /// exp(x) -1 with error cancellation when x is small
0083 inline double expm1( double x) {
0084 #ifndef HAVE_NO_EXPM1
0085    return ::expm1(x);
0086 #else
0087    // compute using taylor expansion until difference is less than epsilon
0088    // use for values smaller than 0.5 (for larger (exp(x)-1 is fine
0089    if (std::abs(x) < 0.5)
0090    {
0091        // taylor series S = x + (1/2!) x^2 + (1/3!) x^3 + ...
0092
0093       double i = 1.0;
0094       double sum = x;
0095       double term = x / 1.0;
0096       do {
0097          i++ ;
0098          term *= x/i;
0099          sum += term;
0100       }
0101       while (std::abs(term) > std::abs(sum) * std::numeric_limits<double>::epsilon() ) ;
0102
0103       return sum ;
0104    }
0105    else
0106    {
0107       return std::exp(x) - 1;
0108    }
0109 #endif
0110 }
0111
0112 } // end namespace Math
0113
0114 } // end namespace ROOT
0115
0116
0117
0118
0119
0120 #endif /* ROOT_Math_Math */