ceph/src/boost/libs/math/include/boost/math/distributions/rayleigh.hpp

   1 //  Copyright Paul A. Bristow 2007.
   2 //  Use, modification and distribution are subject to the
   3 //  Boost Software License, Version 1.0. (See accompanying file
   4 //  LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
   5
   6 #ifndef BOOST_STATS_rayleigh_HPP
   7 #define BOOST_STATS_rayleigh_HPP
   8
   9 #include <boost/math/distributions/fwd.hpp>
  10 #include <boost/math/constants/constants.hpp>
  11 #include <boost/math/special_functions/log1p.hpp>
  12 #include <boost/math/special_functions/expm1.hpp>
  13 #include <boost/math/distributions/complement.hpp>
  14 #include <boost/math/distributions/detail/common_error_handling.hpp>
  15 #include <boost/config/no_tr1/cmath.hpp>
  16
  17 #ifdef BOOST_MSVC
  18 # pragma warning(push)
  19 # pragma warning(disable: 4702) // unreachable code (return after domain_error throw).
  20 #endif
  21
  22 #include <utility>
  23
  24 namespace boost{ namespace math{
  25
  26 namespace detail
  27 { // Error checks:
  28   template <class RealType, class Policy>
  29   inline bool verify_sigma(const char* function, RealType sigma, RealType* presult, const Policy& pol)
  30   {
  31      if((sigma <= 0) || (!(boost::math::isfinite)(sigma)))
  32      {
  33         *presult = policies::raise_domain_error<RealType>(
  34            function,
  35            "The scale parameter \"sigma\" must be > 0 and finite, but was: %1%.", sigma, pol);
  36         return false;
  37      }
  38      return true;
  39   } // bool verify_sigma
  40
  41   template <class RealType, class Policy>
  42   inline bool verify_rayleigh_x(const char* function, RealType x, RealType* presult, const Policy& pol)
  43   {
  44      if((x < 0) || (boost::math::isnan)(x))
  45      {
  46         *presult = policies::raise_domain_error<RealType>(
  47            function,
  48            "The random variable must be >= 0, but was: %1%.", x, pol);
  49         return false;
  50      }
  51      return true;
  52   } // bool verify_rayleigh_x
  53 } // namespace detail
  54
  55 template <class RealType = double, class Policy = policies::policy<> >
  56 class rayleigh_distribution
  57 {
  58 public:
  59    typedef RealType value_type;
  60    typedef Policy policy_type;
  61
  62    rayleigh_distribution(RealType l_sigma = 1)
  63       : m_sigma(l_sigma)
  64    {
  65       RealType err;
  66       detail::verify_sigma("boost::math::rayleigh_distribution<%1%>::rayleigh_distribution", l_sigma, &err, Policy());
  67    } // rayleigh_distribution
  68
  69    RealType sigma()const
  70    { // Accessor.
  71      return m_sigma;
  72    }
  73
  74 private:
  75    RealType m_sigma;
  76 }; // class rayleigh_distribution
  77
  78 typedef rayleigh_distribution<double> rayleigh;
  79
  80 template <class RealType, class Policy>
  81 inline const std::pair<RealType, RealType> range(const rayleigh_distribution<RealType, Policy>& /*dist*/)
  82 { // Range of permissible values for random variable x.
  83    using boost::math::tools::max_value;
  84    return std::pair<RealType, RealType>(static_cast<RealType>(0), std::numeric_limits<RealType>::has_infinity ? std::numeric_limits<RealType>::infinity() : max_value<RealType>());
  85 }
  86
  87 template <class RealType, class Policy>
  88 inline const std::pair<RealType, RealType> support(const rayleigh_distribution<RealType, Policy>& /*dist*/)
  89 { // Range of supported values for random variable x.
  90    // This is range where cdf rises from 0 to 1, and outside it, the pdf is zero.
  91    using boost::math::tools::max_value;
  92    return std::pair<RealType, RealType>(static_cast<RealType>(0),  max_value<RealType>());
  93 }
  94
  95 template <class RealType, class Policy>
  96 inline RealType pdf(const rayleigh_distribution<RealType, Policy>& dist, const RealType& x)
  97 {
  98    BOOST_MATH_STD_USING // for ADL of std function exp.
  99
 100    RealType sigma = dist.sigma();
 101    RealType result = 0;
 102    static const char* function = "boost::math::pdf(const rayleigh_distribution<%1%>&, %1%)";
 103    if(false == detail::verify_sigma(function, sigma, &result, Policy()))
 104    {
 105       return result;
 106    }
 107    if(false == detail::verify_rayleigh_x(function, x, &result, Policy()))
 108    {
 109       return result;
 110    }
 111    if((boost::math::isinf)(x))
 112    {
 113       return 0;
 114    }
 115    RealType sigmasqr = sigma * sigma;
 116    result = x * (exp(-(x * x) / ( 2 * sigmasqr))) / sigmasqr;
 117    return result;
 118 } // pdf
 119
 120 template <class RealType, class Policy>
 121 inline RealType cdf(const rayleigh_distribution<RealType, Policy>& dist, const RealType& x)
 122 {
 123    BOOST_MATH_STD_USING // for ADL of std functions
 124
 125    RealType result = 0;
 126    RealType sigma = dist.sigma();
 127    static const char* function = "boost::math::cdf(const rayleigh_distribution<%1%>&, %1%)";
 128    if(false == detail::verify_sigma(function, sigma, &result, Policy()))
 129    {
 130       return result;
 131    }
 132    if(false == detail::verify_rayleigh_x(function, x, &result, Policy()))
 133    {
 134       return result;
 135    }
 136    result = -boost::math::expm1(-x * x / ( 2 * sigma * sigma), Policy());
 137    return result;
 138 } // cdf
 139
 140 template <class RealType, class Policy>
 141 inline RealType quantile(const rayleigh_distribution<RealType, Policy>& dist, const RealType& p)
 142 {
 143    BOOST_MATH_STD_USING // for ADL of std functions
 144
 145    RealType result = 0;
 146    RealType sigma = dist.sigma();
 147    static const char* function = "boost::math::quantile(const rayleigh_distribution<%1%>&, %1%)";
 148    if(false == detail::verify_sigma(function, sigma, &result, Policy()))
 149       return result;
 150    if(false == detail::check_probability(function, p, &result, Policy()))
 151       return result;
 152
 153    if(p == 0)
 154    {
 155       return 0;
 156    }
 157    if(p == 1)
 158    {
 159      return policies::raise_overflow_error<RealType>(function, 0, Policy());
 160    }
 161    result = sqrt(-2 * sigma * sigma * boost::math::log1p(-p, Policy()));
 162    return result;
 163 } // quantile
 164
 165 template <class RealType, class Policy>
 166 inline RealType cdf(const complemented2_type<rayleigh_distribution<RealType, Policy>, RealType>& c)
 167 {
 168    BOOST_MATH_STD_USING // for ADL of std functions
 169
 170    RealType result = 0;
 171    RealType sigma = c.dist.sigma();
 172    static const char* function = "boost::math::cdf(const rayleigh_distribution<%1%>&, %1%)";
 173    if(false == detail::verify_sigma(function, sigma, &result, Policy()))
 174    {
 175       return result;
 176    }
 177    RealType x = c.param;
 178    if(false == detail::verify_rayleigh_x(function, x, &result, Policy()))
 179    {
 180       return result;
 181    }
 182    RealType ea = x * x / (2 * sigma * sigma);
 183    // Fix for VC11/12 x64 bug in exp(float):
 184    if (ea >= tools::max_value<RealType>())
 185       return 0;
 186    result =  exp(-ea);
 187    return result;
 188 } // cdf complement
 189
 190 template <class RealType, class Policy>
 191 inline RealType quantile(const complemented2_type<rayleigh_distribution<RealType, Policy>, RealType>& c)
 192 {
 193    BOOST_MATH_STD_USING // for ADL of std functions, log & sqrt.
 194
 195    RealType result = 0;
 196    RealType sigma = c.dist.sigma();
 197    static const char* function = "boost::math::quantile(const rayleigh_distribution<%1%>&, %1%)";
 198    if(false == detail::verify_sigma(function, sigma, &result, Policy()))
 199    {
 200       return result;
 201    }
 202    RealType q = c.param;
 203    if(false == detail::check_probability(function, q, &result, Policy()))
 204    {
 205       return result;
 206    }
 207    if(q == 1)
 208    {
 209       return 0;
 210    }
 211    if(q == 0)
 212    {
 213      return policies::raise_overflow_error<RealType>(function, 0, Policy());
 214    }
 215    result = sqrt(-2 * sigma * sigma * log(q));
 216    return result;
 217 } // quantile complement
 218
 219 template <class RealType, class Policy>
 220 inline RealType mean(const rayleigh_distribution<RealType, Policy>& dist)
 221 {
 222    RealType result = 0;
 223    RealType sigma = dist.sigma();
 224    static const char* function = "boost::math::mean(const rayleigh_distribution<%1%>&, %1%)";
 225    if(false == detail::verify_sigma(function, sigma, &result, Policy()))
 226    {
 227       return result;
 228    }
 229    using boost::math::constants::root_half_pi;
 230    return sigma * root_half_pi<RealType>();
 231 } // mean
 232
 233 template <class RealType, class Policy>
 234 inline RealType variance(const rayleigh_distribution<RealType, Policy>& dist)
 235 {
 236    RealType result = 0;
 237    RealType sigma = dist.sigma();
 238    static const char* function = "boost::math::variance(const rayleigh_distribution<%1%>&, %1%)";
 239    if(false == detail::verify_sigma(function, sigma, &result, Policy()))
 240    {
 241       return result;
 242    }
 243    using boost::math::constants::four_minus_pi;
 244    return four_minus_pi<RealType>() * sigma * sigma / 2;
 245 } // variance
 246
 247 template <class RealType, class Policy>
 248 inline RealType mode(const rayleigh_distribution<RealType, Policy>& dist)
 249 {
 250    return dist.sigma();
 251 }
 252
 253 template <class RealType, class Policy>
 254 inline RealType median(const rayleigh_distribution<RealType, Policy>& dist)
 255 {
 256    using boost::math::constants::root_ln_four;
 257    return root_ln_four<RealType>() * dist.sigma();
 258 }
 259
 260 template <class RealType, class Policy>
 261 inline RealType skewness(const rayleigh_distribution<RealType, Policy>& /*dist*/)
 262 {
 263   // using namespace boost::math::constants;
 264   return static_cast<RealType>(0.63111065781893713819189935154422777984404221106391L);
 265   // Computed using NTL at 150 bit, about 50 decimal digits.
 266   // return 2 * root_pi<RealType>() * pi_minus_three<RealType>() / pow23_four_minus_pi<RealType>();
 267 }
 268
 269 template <class RealType, class Policy>
 270 inline RealType kurtosis(const rayleigh_distribution<RealType, Policy>& /*dist*/)
 271 {
 272   // using namespace boost::math::constants;
 273   return static_cast<RealType>(3.2450893006876380628486604106197544154170667057995L);
 274   // Computed using NTL at 150 bit, about 50 decimal digits.
 275   // return 3 - (6 * pi<RealType>() * pi<RealType>() - 24 * pi<RealType>() + 16) /
 276   // (four_minus_pi<RealType>() * four_minus_pi<RealType>());
 277 }
 278
 279 template <class RealType, class Policy>
 280 inline RealType kurtosis_excess(const rayleigh_distribution<RealType, Policy>& /*dist*/)
 281 {
 282   //using namespace boost::math::constants;
 283   // Computed using NTL at 150 bit, about 50 decimal digits.
 284   return static_cast<RealType>(0.2450893006876380628486604106197544154170667057995L);
 285   // return -(6 * pi<RealType>() * pi<RealType>() - 24 * pi<RealType>() + 16) /
 286   //   (four_minus_pi<RealType>() * four_minus_pi<RealType>());
 287 } // kurtosis
 288
 289 } // namespace math
 290 } // namespace boost
 291
 292 #ifdef BOOST_MSVC
 293 # pragma warning(pop)
 294 #endif
 295
 296 // This include must be at the end, *after* the accessors
 297 // for this distribution have been defined, in order to
 298 // keep compilers that support two-phase lookup happy.
 299 #include <boost/math/distributions/detail/derived_accessors.hpp>
 300
 301 #endif // BOOST_STATS_rayleigh_HPP