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

Add logcdf specialization to rayleigh distribution

Browse Source
This commit is contained in:
Matt Borland
2023-02-11 18:07:27 -08:00
parent 614445505a
commit 7b678a9001
2 changed files with 69 additions and 0 deletions
Download Patch File Download Diff File Expand all files Collapse all files

47
include/boost/math/distributions/rayleigh.hpp
View File

@@ -1,4 +1,5 @@
// Copyright Paul A. Bristow 2007.
// Copyright Matt Borland 2023.
// Use, modification and distribution are subject to the
// Boost Software License, Version 1.0. (See accompanying file
// LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
@@ -19,6 +20,7 @@
#endif
#include <utility>
#include <limits>
#include <cmath>
namespace boost{ namespace math{
@@ -168,6 +170,26 @@ inline RealType cdf(const rayleigh_distribution<RealType, Policy>& dist, const R
return result;
} // cdf
template <class RealType, class Policy>
inline RealType logcdf(const rayleigh_distribution<RealType, Policy>& dist, const RealType& x)
{
BOOST_MATH_STD_USING // for ADL of std functions
RealType result = 0;
RealType sigma = dist.sigma();
static const char* function = "boost::math::logcdf(const rayleigh_distribution<%1%>&, %1%)";
if(false == detail::verify_sigma(function, sigma, &result, Policy()))
{
return -std::numeric_limits<RealType>::infinity();
}
if(false == detail::verify_rayleigh_x(function, x, &result, Policy()))
{
return -std::numeric_limits<RealType>::infinity();
}
result = log1p(-exp(-x * x / ( 2 * sigma * sigma)), Policy());
return result;
} // logcdf
template <class RealType, class Policy>
inline RealType quantile(const rayleigh_distribution<RealType, Policy>& dist, const RealType& p)
{
@@ -218,6 +240,31 @@ inline RealType cdf(const complemented2_type<rayleigh_distribution<RealType, Pol
return result;
} // cdf complement
template <class RealType, class Policy>
inline RealType logcdf(const complemented2_type<rayleigh_distribution<RealType, Policy>, RealType>& c)
{
BOOST_MATH_STD_USING // for ADL of std functions
RealType result = 0;
RealType sigma = c.dist.sigma();
static const char* function = "boost::math::logcdf(const rayleigh_distribution<%1%>&, %1%)";
if(false == detail::verify_sigma(function, sigma, &result, Policy()))
{
return -std::numeric_limits<RealType>::infinity();
}
RealType x = c.param;
if(false == detail::verify_rayleigh_x(function, x, &result, Policy()))
{
return -std::numeric_limits<RealType>::infinity();
}
RealType ea = x * x / (2 * sigma * sigma);
// Fix for VC11/12 x64 bug in exp(float):
if (ea >= tools::max_value<RealType>())
return 0;
result = -ea;
return result;
} // logcdf complement
template <class RealType, class Policy>
inline RealType quantile(const complemented2_type<rayleigh_distribution<RealType, Policy>, RealType>& c)
{

22
test/test_rayleigh.cpp
View File

@@ -1,4 +1,5 @@
// Copyright John Maddock 2006.
// Copyright Matt Borland 2023.
// Use, modification and distribution are subject to the
// Boost Software License, Version 1.0.
@@ -28,10 +29,19 @@
using std::setprecision;
#include <cmath>
using std::log;
#include <type_traits>
template <class RealType>
void test_spot(RealType s, RealType x, RealType p, RealType q, RealType tolerance)
{
RealType logtolerance = tolerance;
BOOST_IF_CONSTEXPR (std::is_same<RealType, long double>::value ||
std::is_same<RealType, boost::math::concepts::real_concept>::value)
{
logtolerance *= 100;
}
BOOST_CHECK_CLOSE(
::boost::math::cdf(
rayleigh_distribution<RealType>(s),
@@ -44,6 +54,18 @@ void test_spot(RealType s, RealType x, RealType p, RealType q, RealType toleranc
x)),
q,
tolerance); // %
BOOST_CHECK_CLOSE(
::boost::math::logcdf(
rayleigh_distribution<RealType>(s),
x),
log(p),
tolerance); // %
BOOST_CHECK_CLOSE(
::boost::math::logcdf(
complement(rayleigh_distribution<RealType>(s),
x)),
log(q),
tolerance); // %
// Special extra tests for p and q near to unity.
if(p < 0.999)
{