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

Make a few distributions CUDA-compatible.

Browse Source
This commit is contained in:
jzmaddock
2016-01-28 09:15:25 +00:00
parent ad22a6a00a
commit b1131fdc3f
16 changed files with 1385 additions and 76 deletions
Download Patch File Download Diff File Expand all files Collapse all files

68
include/boost/math/distributions/arcsine.hpp
View File

@@ -54,7 +54,7 @@ namespace boost
// Common error checking routines for arcsine distribution functions:
// Duplicating for x_min and x_max provides specific error messages.
template <class RealType, class Policy>
inline bool check_x_min(const char* function, const RealType& x, RealType* result, const Policy& pol)
inline BOOST_GPU_ENABLED bool check_x_min(const char* function, const RealType& x, RealType* result, const Policy& pol)
{
if (!(boost::math::isfinite)(x))
{
@@ -67,7 +67,7 @@ namespace boost
} // bool check_x_min
template <class RealType, class Policy>
inline bool check_x_max(const char* function, const RealType& x, RealType* result, const Policy& pol)
inline BOOST_GPU_ENABLED bool check_x_max(const char* function, const RealType& x, RealType* result, const Policy& pol)
{
if (!(boost::math::isfinite)(x))
{
@@ -81,11 +81,16 @@ namespace boost
template <class RealType, class Policy>
inline bool check_x_minmax(const char* function, const RealType& x_min, const RealType& x_max, RealType* result, const Policy& pol)
inline BOOST_GPU_ENABLED bool check_x_minmax(const char* function, const RealType& x_min, const RealType& x_max, RealType* result, const Policy& pol)
{ // Check x_min < x_max
if (x_min >= x_max)
{
std::string msg = "x_max argument is %1%, but must be > x_min = " + lexical_cast<std::string>(x_min) + "!";
#ifdef __CUDA_ARCH__
*result = policies::raise_domain_error<RealType>(
function,
"x_max value %1% is out of range", x_max, pol);
#else
std::string msg = "x_max argument is %1%, but must be > x_min = " + lexical_cast<std::string>(x_min) + "!";
*result = policies::raise_domain_error<RealType>(
function,
msg.c_str(), x_max, pol);
@@ -94,13 +99,14 @@ namespace boost
// But would require replication of all helpers functions in /policies/error_handling.hpp for two values,
// as well as two value versions of raise_error, raise_domain_error and do_format ...
// so use slightly hacky lexical_cast to string instead.
#endif
return false;
}
return true;
} // bool check_x_minmax
template <class RealType, class Policy>
inline bool check_prob(const char* function, const RealType& p, RealType* result, const Policy& pol)
inline BOOST_GPU_ENABLED bool check_prob(const char* function, const RealType& p, RealType* result, const Policy& pol)
{
if ((p < 0) || (p > 1) || !(boost::math::isfinite)(p))
{
@@ -113,7 +119,7 @@ namespace boost
} // bool check_prob
template <class RealType, class Policy>
inline bool check_x(const char* function, const RealType& x_min, const RealType& x_max, const RealType& x, RealType* result, const Policy& pol)
inline BOOST_GPU_ENABLED bool check_x(const char* function, const RealType& x_min, const RealType& x_max, const RealType& x, RealType* result, const Policy& pol)
{ // Check x finite and x_min < x < x_max.
if (!(boost::math::isfinite)(x))
{
@@ -137,7 +143,7 @@ namespace boost
} // bool check_x
template <class RealType, class Policy>
inline bool check_dist(const char* function, const RealType& x_min, const RealType& x_max, RealType* result, const Policy& pol)
inline BOOST_GPU_ENABLED bool check_dist(const char* function, const RealType& x_min, const RealType& x_max, RealType* result, const Policy& pol)
{ // Check both x_min and x_max finite, and x_min < x_max.
return check_x_min(function, x_min, result, pol)
&& check_x_max(function, x_max, result, pol)
@@ -145,14 +151,14 @@ namespace boost
} // bool check_dist
template <class RealType, class Policy>
inline bool check_dist_and_x(const char* function, const RealType& x_min, const RealType& x_max, RealType x, RealType* result, const Policy& pol)
inline BOOST_GPU_ENABLED bool check_dist_and_x(const char* function, const RealType& x_min, const RealType& x_max, RealType x, RealType* result, const Policy& pol)
{
return check_dist(function, x_min, x_max, result, pol)
&& arcsine_detail::check_x(function, x_min, x_max, x, result, pol);
} // bool check_dist_and_x
template <class RealType, class Policy>
inline bool check_dist_and_prob(const char* function, const RealType& x_min, const RealType& x_max, RealType p, RealType* result, const Policy& pol)
inline BOOST_GPU_ENABLED bool check_dist_and_prob(const char* function, const RealType& x_min, const RealType& x_max, RealType p, RealType* result, const Policy& pol)
{
return check_dist(function, x_min, x_max, result, pol)
&& check_prob(function, p, result, pol);
@@ -167,7 +173,7 @@ namespace boost
typedef RealType value_type;
typedef Policy policy_type;
arcsine_distribution(RealType x_min = 0, RealType x_max = 1) : m_x_min(x_min), m_x_max(x_max)
BOOST_GPU_ENABLED arcsine_distribution(RealType x_min = 0, RealType x_max = 1) : m_x_min(x_min), m_x_max(x_max)
{ // Default beta (alpha = beta = 0.5) is standard arcsine with x_min = 0, x_max = 1.
// Generalized to allow x_min and x_max to be specified.
RealType result;
@@ -178,11 +184,11 @@ namespace boost
&result, Policy());
} // arcsine_distribution constructor.
// Accessor functions:
RealType x_min() const
BOOST_GPU_ENABLED RealType x_min() const
{
return m_x_min;
}
RealType x_max() const
BOOST_GPU_ENABLED RealType x_max() const
{
return m_x_max;
}
@@ -197,21 +203,21 @@ namespace boost
template <class RealType, class Policy>
inline const std::pair<RealType, RealType> range(const arcsine_distribution<RealType, Policy>& dist)
inline BOOST_GPU_ENABLED const std::pair<RealType, RealType> range(const arcsine_distribution<RealType, Policy>& dist)
{ // Range of permissible values for random variable x.
using boost::math::tools::max_value;
return std::pair<RealType, RealType>(static_cast<RealType>(dist.x_min()), static_cast<RealType>(dist.x_max()));
}
template <class RealType, class Policy>
inline const std::pair<RealType, RealType> support(const arcsine_distribution<RealType, Policy>& dist)
inline BOOST_GPU_ENABLED const std::pair<RealType, RealType> support(const arcsine_distribution<RealType, Policy>& dist)
{ // Range of supported values for random variable x.
// This is range where cdf rises from 0 to 1, and outside it, the pdf is zero.
return std::pair<RealType, RealType>(static_cast<RealType>(dist.x_min()), static_cast<RealType>(dist.x_max()));
}
template <class RealType, class Policy>
inline RealType mean(const arcsine_distribution<RealType, Policy>& dist)
inline BOOST_GPU_ENABLED RealType mean(const arcsine_distribution<RealType, Policy>& dist)
{ // Mean of arcsine distribution .
RealType result;
RealType x_min = dist.x_min();
@@ -230,7 +236,7 @@ namespace boost
} // mean
template <class RealType, class Policy>
inline RealType variance(const arcsine_distribution<RealType, Policy>& dist)
inline BOOST_GPU_ENABLED RealType variance(const arcsine_distribution<RealType, Policy>& dist)
{ // Variance of standard arcsine distribution = (1-0)/8 = 0.125.
RealType result;
RealType x_min = dist.x_min();
@@ -248,7 +254,7 @@ namespace boost
} // variance
template <class RealType, class Policy>
inline RealType mode(const arcsine_distribution<RealType, Policy>& /* dist */)
inline BOOST_GPU_ENABLED RealType mode(const arcsine_distribution<RealType, Policy>& /* dist */)
{ //There are always [*two] values for the mode, at ['x_min] and at ['x_max], default 0 and 1,
// so instead we raise the exception domain_error.
return policies::raise_domain_error<RealType>(
@@ -259,7 +265,7 @@ namespace boost
} // mode
template <class RealType, class Policy>
inline RealType median(const arcsine_distribution<RealType, Policy>& dist)
inline BOOST_GPU_ENABLED RealType median(const arcsine_distribution<RealType, Policy>& dist)
{ // Median of arcsine distribution (a + b) / 2 == mean.
RealType x_min = dist.x_min();
RealType x_max = dist.x_max();
@@ -277,7 +283,7 @@ namespace boost
}
template <class RealType, class Policy>
inline RealType skewness(const arcsine_distribution<RealType, Policy>& dist)
inline BOOST_GPU_ENABLED RealType skewness(const arcsine_distribution<RealType, Policy>& dist)
{
RealType result;
RealType x_min = dist.x_min();
@@ -296,7 +302,7 @@ namespace boost
} // skewness
template <class RealType, class Policy>
inline RealType kurtosis_excess(const arcsine_distribution<RealType, Policy>& dist)
inline BOOST_GPU_ENABLED RealType kurtosis_excess(const arcsine_distribution<RealType, Policy>& dist)
{
RealType result;
RealType x_min = dist.x_min();
@@ -316,7 +322,7 @@ namespace boost
} // kurtosis_excess
template <class RealType, class Policy>
inline RealType kurtosis(const arcsine_distribution<RealType, Policy>& dist)
inline BOOST_GPU_ENABLED RealType kurtosis(const arcsine_distribution<RealType, Policy>& dist)
{
RealType result;
RealType x_min = dist.x_min();
@@ -336,12 +342,12 @@ namespace boost
} // kurtosis
template <class RealType, class Policy>
inline RealType pdf(const arcsine_distribution<RealType, Policy>& dist, const RealType& xx)
inline BOOST_GPU_ENABLED RealType pdf(const arcsine_distribution<RealType, Policy>& dist, const RealType& xx)
{ // Probability Density/Mass Function arcsine.
BOOST_FPU_EXCEPTION_GUARD
BOOST_MATH_STD_USING // For ADL of std functions.
static const char* function = "boost::math::pdf(arcsine_distribution<%1%> const&, %1%)";
BOOST_MATH_GPU_STATIC const char* function = "boost::math::pdf(arcsine_distribution<%1%> const&, %1%)";
RealType lo = dist.x_min();
RealType hi = dist.x_max();
@@ -362,11 +368,11 @@ namespace boost
} // pdf
template <class RealType, class Policy>
inline RealType cdf(const arcsine_distribution<RealType, Policy>& dist, const RealType& x)
inline BOOST_GPU_ENABLED RealType cdf(const arcsine_distribution<RealType, Policy>& dist, const RealType& x)
{ // Cumulative Distribution Function arcsine.
BOOST_MATH_STD_USING // For ADL of std functions.
static const char* function = "boost::math::cdf(arcsine_distribution<%1%> const&, %1%)";
BOOST_MATH_GPU_STATIC const char* function = "boost::math::cdf(arcsine_distribution<%1%> const&, %1%)";
RealType x_min = dist.x_min();
RealType x_max = dist.x_max();
@@ -395,10 +401,10 @@ namespace boost
} // arcsine cdf
template <class RealType, class Policy>
inline RealType cdf(const complemented2_type<arcsine_distribution<RealType, Policy>, RealType>& c)
inline BOOST_GPU_ENABLED RealType cdf(const complemented2_type<arcsine_distribution<RealType, Policy>, RealType>& c)
{ // Complemented Cumulative Distribution Function arcsine.
BOOST_MATH_STD_USING // For ADL of std functions.
static const char* function = "boost::math::cdf(arcsine_distribution<%1%> const&, %1%)";
BOOST_MATH_GPU_STATIC const char* function = "boost::math::cdf(arcsine_distribution<%1%> const&, %1%)";
RealType x = c.param;
arcsine_distribution<RealType, Policy> const& dist = c.dist;
@@ -431,7 +437,7 @@ namespace boost
} // arcine ccdf
template <class RealType, class Policy>
inline RealType quantile(const arcsine_distribution<RealType, Policy>& dist, const RealType& p)
inline BOOST_GPU_ENABLED RealType quantile(const arcsine_distribution<RealType, Policy>& dist, const RealType& p)
{
// Quantile or Percent Point arcsine function or
// Inverse Cumulative probability distribution function CDF.
@@ -445,7 +451,7 @@ namespace boost
using boost::math::constants::half_pi;
static const char* function = "boost::math::quantile(arcsine_distribution<%1%> const&, %1%)";
BOOST_MATH_GPU_STATIC const char* function = "boost::math::quantile(arcsine_distribution<%1%> const&, %1%)";
RealType result = 0; // of argument checks:
RealType x_min = dist.x_min();
@@ -475,7 +481,7 @@ namespace boost
} // quantile
template <class RealType, class Policy>
inline RealType quantile(const complemented2_type<arcsine_distribution<RealType, Policy>, RealType>& c)
inline BOOST_GPU_ENABLED RealType quantile(const complemented2_type<arcsine_distribution<RealType, Policy>, RealType>& c)
{
// Complement Quantile or Percent Point arcsine function.
// Return the number of expected x for a given
@@ -483,7 +489,7 @@ namespace boost
BOOST_MATH_STD_USING // For ADL of std functions.
using boost::math::constants::half_pi;
static const char* function = "boost::math::quantile(arcsine_distribution<%1%> const&, %1%)";
BOOST_MATH_GPU_STATIC const char* function = "boost::math::quantile(arcsine_distribution<%1%> const&, %1%)";
// Error checks:
RealType q = c.param;

53
include/boost/math/distributions/cauchy.hpp
View File

@@ -31,7 +31,7 @@ namespace detail
{
template <class RealType, class Policy>
RealType cdf_imp(const cauchy_distribution<RealType, Policy>& dist, const RealType& x, bool complement)
BOOST_GPU_ENABLED RealType cdf_imp(const cauchy_distribution<RealType, Policy>& dist, const RealType& x, bool complement)
{
//
// This calculates the cdf of the Cauchy distribution and/or its complement.
@@ -55,7 +55,7 @@ RealType cdf_imp(const cauchy_distribution<RealType, Policy>& dist, const RealTy
// to get the result.
//
BOOST_MATH_STD_USING // for ADL of std functions
static const char* function = "boost::math::cdf(cauchy<%1%>&, %1%)";
BOOST_MATH_GPU_STATIC const char* function = "boost::math::cdf(cauchy<%1%>&, %1%)";
RealType result = 0;
RealType location = dist.location();
RealType scale = dist.scale();
@@ -67,6 +67,12 @@ RealType cdf_imp(const cauchy_distribution<RealType, Policy>& dist, const RealTy
{
return result;
}
#ifdef __CUDA_ARCH__
if(x > tools::max_value<RealType>())
return static_cast<RealType>((complement) ? 0 : 1);
if(x < -tools::max_value<RealType>())
return static_cast<RealType>((complement) ? 1 : 0);
#else
if(std::numeric_limits<RealType>::has_infinity && x == std::numeric_limits<RealType>::infinity())
{ // cdf +infinity is unity.
return static_cast<RealType>((complement) ? 0 : 1);
@@ -75,6 +81,7 @@ RealType cdf_imp(const cauchy_distribution<RealType, Policy>& dist, const RealTy
{ // cdf -infinity is zero.
return static_cast<RealType>((complement) ? 1 : 0);
}
#endif
if(false == detail::check_x(function, x, &result, Policy()))
{ // Catches x == NaN
return result;
@@ -89,7 +96,7 @@ RealType cdf_imp(const cauchy_distribution<RealType, Policy>& dist, const RealTy
} // cdf
template <class RealType, class Policy>
RealType quantile_imp(
BOOST_GPU_ENABLED RealType quantile_imp(
const cauchy_distribution<RealType, Policy>& dist,
const RealType& p,
bool complement)
@@ -102,7 +109,7 @@ RealType quantile_imp(
// mid-point of the distribution. This is either added or subtracted
// from the location parameter depending on whether `complement` is true.
//
static const char* function = "boost::math::quantile(cauchy<%1%>&, %1%)";
BOOST_MATH_GPU_STATIC const char* function = "boost::math::quantile(cauchy<%1%>&, %1%)";
BOOST_MATH_STD_USING // for ADL of std functions
RealType result = 0;
@@ -152,20 +159,20 @@ public:
typedef RealType value_type;
typedef Policy policy_type;
cauchy_distribution(RealType l_location = 0, RealType l_scale = 1)
BOOST_GPU_ENABLED cauchy_distribution(RealType l_location = 0, RealType l_scale = 1)
: m_a(l_location), m_hg(l_scale)
{
static const char* function = "boost::math::cauchy_distribution<%1%>::cauchy_distribution";
BOOST_MATH_GPU_STATIC const char* function = "boost::math::cauchy_distribution<%1%>::cauchy_distribution";
RealType result;
detail::check_location(function, l_location, &result, Policy());
detail::check_scale(function, l_scale, &result, Policy());
} // cauchy_distribution
RealType location()const
BOOST_GPU_ENABLED RealType location()const
{
return m_a;
}
RealType scale()const
BOOST_GPU_ENABLED RealType scale()const
{
return m_hg;
}
@@ -178,7 +185,7 @@ private:
typedef cauchy_distribution<double> cauchy;
template <class RealType, class Policy>
inline const std::pair<RealType, RealType> range(const cauchy_distribution<RealType, Policy>&)
inline BOOST_GPU_ENABLED const std::pair<RealType, RealType> range(const cauchy_distribution<RealType, Policy>&)
{ // Range of permissible values for random variable x.
if (std::numeric_limits<RealType>::has_infinity)
{
@@ -192,7 +199,7 @@ inline const std::pair<RealType, RealType> range(const cauchy_distribution<RealT
}
template <class RealType, class Policy>
inline const std::pair<RealType, RealType> support(const cauchy_distribution<RealType, Policy>& )
inline BOOST_GPU_ENABLED const std::pair<RealType, RealType> support(const cauchy_distribution<RealType, Policy>& )
{ // Range of supported values for random variable x.
// This is range where cdf rises from 0 to 1, and outside it, the pdf is zero.
if (std::numeric_limits<RealType>::has_infinity)
@@ -207,11 +214,11 @@ inline const std::pair<RealType, RealType> support(const cauchy_distribution<Rea
}
template <class RealType, class Policy>
inline RealType pdf(const cauchy_distribution<RealType, Policy>& dist, const RealType& x)
inline BOOST_GPU_ENABLED RealType pdf(const cauchy_distribution<RealType, Policy>& dist, const RealType& x)
{
BOOST_MATH_STD_USING // for ADL of std functions
static const char* function = "boost::math::pdf(cauchy<%1%>&, %1%)";
BOOST_MATH_GPU_STATIC const char* function = "boost::math::pdf(cauchy<%1%>&, %1%)";
RealType result = 0;
RealType location = dist.location();
RealType scale = dist.scale();
@@ -244,31 +251,31 @@ inline RealType pdf(const cauchy_distribution<RealType, Policy>& dist, const Rea
} // pdf
template <class RealType, class Policy>
inline RealType cdf(const cauchy_distribution<RealType, Policy>& dist, const RealType& x)
inline BOOST_GPU_ENABLED RealType cdf(const cauchy_distribution<RealType, Policy>& dist, const RealType& x)
{
return detail::cdf_imp(dist, x, false);
} // cdf
template <class RealType, class Policy>
inline RealType quantile(const cauchy_distribution<RealType, Policy>& dist, const RealType& p)
inline BOOST_GPU_ENABLED RealType quantile(const cauchy_distribution<RealType, Policy>& dist, const RealType& p)
{
return detail::quantile_imp(dist, p, false);
} // quantile
template <class RealType, class Policy>
inline RealType cdf(const complemented2_type<cauchy_distribution<RealType, Policy>, RealType>& c)
inline BOOST_GPU_ENABLED RealType cdf(const complemented2_type<cauchy_distribution<RealType, Policy>, RealType>& c)
{
return detail::cdf_imp(c.dist, c.param, true);
} // cdf complement
template <class RealType, class Policy>
inline RealType quantile(const complemented2_type<cauchy_distribution<RealType, Policy>, RealType>& c)
inline BOOST_GPU_ENABLED RealType quantile(const complemented2_type<cauchy_distribution<RealType, Policy>, RealType>& c)
{
return detail::quantile_imp(c.dist, c.param, true);
} // quantile complement
template <class RealType, class Policy>
inline RealType mean(const cauchy_distribution<RealType, Policy>&)
inline BOOST_GPU_ENABLED RealType mean(const cauchy_distribution<RealType, Policy>&)
{ // There is no mean:
typedef typename Policy::assert_undefined_type assert_type;
BOOST_STATIC_ASSERT(assert_type::value == 0);
@@ -281,7 +288,7 @@ inline RealType mean(const cauchy_distribution<RealType, Policy>&)
}
template <class RealType, class Policy>
inline RealType variance(const cauchy_distribution<RealType, Policy>& /*dist*/)
inline BOOST_GPU_ENABLED RealType variance(const cauchy_distribution<RealType, Policy>& /*dist*/)
{
// There is no variance:
typedef typename Policy::assert_undefined_type assert_type;
@@ -295,18 +302,18 @@ inline RealType variance(const cauchy_distribution<RealType, Policy>& /*dist*/)
}
template <class RealType, class Policy>
inline RealType mode(const cauchy_distribution<RealType, Policy>& dist)
inline BOOST_GPU_ENABLED RealType mode(const cauchy_distribution<RealType, Policy>& dist)
{
return dist.location();
}
template <class RealType, class Policy>
inline RealType median(const cauchy_distribution<RealType, Policy>& dist)
inline BOOST_GPU_ENABLED RealType median(const cauchy_distribution<RealType, Policy>& dist)
{
return dist.location();
}
template <class RealType, class Policy>
inline RealType skewness(const cauchy_distribution<RealType, Policy>& /*dist*/)
inline BOOST_GPU_ENABLED RealType skewness(const cauchy_distribution<RealType, Policy>& /*dist*/)
{
// There is no skewness:
typedef typename Policy::assert_undefined_type assert_type;
@@ -320,7 +327,7 @@ inline RealType skewness(const cauchy_distribution<RealType, Policy>& /*dist*/)
}
template <class RealType, class Policy>
inline RealType kurtosis(const cauchy_distribution<RealType, Policy>& /*dist*/)
inline BOOST_GPU_ENABLED RealType kurtosis(const cauchy_distribution<RealType, Policy>& /*dist*/)
{
// There is no kurtosis:
typedef typename Policy::assert_undefined_type assert_type;
@@ -334,7 +341,7 @@ inline RealType kurtosis(const cauchy_distribution<RealType, Policy>& /*dist*/)
}
template <class RealType, class Policy>
inline RealType kurtosis_excess(const cauchy_distribution<RealType, Policy>& /*dist*/)
inline BOOST_GPU_ENABLED RealType kurtosis_excess(const cauchy_distribution<RealType, Policy>& /*dist*/)
{
// There is no kurtosis excess:
typedef typename Policy::assert_undefined_type assert_type;

20
include/boost/math/distributions/detail/common_error_handling.hpp
View File

@@ -18,7 +18,7 @@ namespace boost{ namespace math{ namespace detail
{
template <class RealType, class Policy>
inline bool check_probability(const char* function, RealType const& prob, RealType* result, const Policy& pol)
inline BOOST_GPU_ENABLED bool check_probability(const char* function, RealType const& prob, RealType* result, const Policy& pol)
{
if((prob < 0) || (prob > 1) || !(boost::math::isfinite)(prob))
{
@@ -31,7 +31,7 @@ inline bool check_probability(const char* function, RealType const& prob, RealTy
}
template <class RealType, class Policy>
inline bool check_df(const char* function, RealType const& df, RealType* result, const Policy& pol)
inline BOOST_GPU_ENABLED bool check_df(const char* function, RealType const& df, RealType* result, const Policy& pol)
{ // df > 0 but NOT +infinity allowed.
if((df <= 0) || !(boost::math::isfinite)(df))
{
@@ -44,7 +44,7 @@ inline bool check_df(const char* function, RealType const& df, RealType* result,
}
template <class RealType, class Policy>
inline bool check_df_gt0_to_inf(const char* function, RealType const& df, RealType* result, const Policy& pol)
inline BOOST_GPU_ENABLED bool check_df_gt0_to_inf(const char* function, RealType const& df, RealType* result, const Policy& pol)
{ // df > 0 or +infinity are allowed.
if( (df <= 0) || (boost::math::isnan)(df) )
{ // is bad df <= 0 or NaN or -infinity.
@@ -58,7 +58,7 @@ inline bool check_df_gt0_to_inf(const char* function, RealType const& df, RealTy
template <class RealType, class Policy>
inline bool check_scale(
inline BOOST_GPU_ENABLED bool check_scale(
const char* function,
RealType scale,
RealType* result,
@@ -75,7 +75,7 @@ inline bool check_scale(
}
template <class RealType, class Policy>
inline bool check_location(
inline BOOST_GPU_ENABLED bool check_location(
const char* function,
RealType location,
RealType* result,
@@ -92,7 +92,7 @@ inline bool check_location(
}
template <class RealType, class Policy>
inline bool check_x(
inline BOOST_GPU_ENABLED bool check_x(
const char* function,
RealType x,
RealType* result,
@@ -113,7 +113,7 @@ inline bool check_x(
} // bool check_x
template <class RealType, class Policy>
inline bool check_x_gt0(
inline BOOST_GPU_ENABLED bool check_x_gt0(
const char* function,
RealType x,
RealType* result,
@@ -134,7 +134,7 @@ inline bool check_x_gt0(
} // bool check_x_gt0
template <class RealType, class Policy>
inline bool check_positive_x(
inline BOOST_GPU_ENABLED bool check_positive_x(
const char* function,
RealType x,
RealType* result,
@@ -154,7 +154,7 @@ inline bool check_positive_x(
}
template <class RealType, class Policy>
inline bool check_non_centrality(
inline BOOST_GPU_ENABLED bool check_non_centrality(
const char* function,
RealType ncp,
RealType* result,
@@ -171,7 +171,7 @@ inline bool check_non_centrality(
}
template <class RealType, class Policy>
inline bool check_finite(
inline BOOST_GPU_ENABLED bool check_finite(
const char* function,
RealType x,
RealType* result,

24
include/boost/math/distributions/detail/derived_accessors.hpp
View File

@@ -39,24 +39,24 @@
namespace boost{ namespace math{
template <class Distribution>
typename Distribution::value_type variance(const Distribution& dist);
BOOST_GPU_ENABLED typename Distribution::value_type variance(const Distribution& dist);
template <class Distribution>
inline typename Distribution::value_type standard_deviation(const Distribution& dist)
inline BOOST_GPU_ENABLED typename Distribution::value_type standard_deviation(const Distribution& dist)
{
BOOST_MATH_STD_USING // ADL of sqrt.
return sqrt(variance(dist));
}
template <class Distribution>
inline typename Distribution::value_type variance(const Distribution& dist)
inline BOOST_GPU_ENABLED typename Distribution::value_type variance(const Distribution& dist)
{
typename Distribution::value_type result = standard_deviation(dist);
return result * result;
}
template <class Distribution, class RealType>
inline typename Distribution::value_type hazard(const Distribution& dist, const RealType& x)
inline BOOST_GPU_ENABLED typename Distribution::value_type hazard(const Distribution& dist, const RealType& x)
{ // hazard function
// http://www.itl.nist.gov/div898/handbook/eda/section3/eda362.htm#HAZ
typedef typename Distribution::value_type value_type;
@@ -75,7 +75,7 @@ inline typename Distribution::value_type hazard(const Distribution& dist, const
}
template <class Distribution, class RealType>
inline typename Distribution::value_type chf(const Distribution& dist, const RealType& x)
inline BOOST_GPU_ENABLED typename Distribution::value_type chf(const Distribution& dist, const RealType& x)
{ // cumulative hazard function.
// http://www.itl.nist.gov/div898/handbook/eda/section3/eda362.htm#HAZ
BOOST_MATH_STD_USING
@@ -83,7 +83,7 @@ inline typename Distribution::value_type chf(const Distribution& dist, const Rea
}
template <class Distribution>
inline typename Distribution::value_type coefficient_of_variation(const Distribution& dist)
inline BOOST_GPU_ENABLED typename Distribution::value_type coefficient_of_variation(const Distribution& dist)
{
typedef typename Distribution::value_type value_type;
typedef typename Distribution::policy_type policy_type;
@@ -104,19 +104,19 @@ inline typename Distribution::value_type coefficient_of_variation(const Distribu
// implementation with all arguments of the same type:
//
template <class Distribution, class RealType>
inline typename Distribution::value_type pdf(const Distribution& dist, const RealType& x)
inline BOOST_GPU_ENABLED typename Distribution::value_type pdf(const Distribution& dist, const RealType& x)
{
typedef typename Distribution::value_type value_type;
return pdf(dist, static_cast<value_type>(x));
}
template <class Distribution, class RealType>
inline typename Distribution::value_type cdf(const Distribution& dist, const RealType& x)
inline BOOST_GPU_ENABLED typename Distribution::value_type cdf(const Distribution& dist, const RealType& x)
{
typedef typename Distribution::value_type value_type;
return cdf(dist, static_cast<value_type>(x));
}
template <class Distribution, class RealType>
inline typename Distribution::value_type quantile(const Distribution& dist, const RealType& x)
inline BOOST_GPU_ENABLED typename Distribution::value_type quantile(const Distribution& dist, const RealType& x)
{
typedef typename Distribution::value_type value_type;
return quantile(dist, static_cast<value_type>(x));
@@ -130,21 +130,21 @@ inline typename Distribution::value_type chf(const Distribution& dist, const Rea
}
*/
template <class Distribution, class RealType>
inline typename Distribution::value_type cdf(const complemented2_type<Distribution, RealType>& c)
inline BOOST_GPU_ENABLED typename Distribution::value_type cdf(const complemented2_type<Distribution, RealType>& c)
{
typedef typename Distribution::value_type value_type;
return cdf(complement(c.dist, static_cast<value_type>(c.param)));
}
template <class Distribution, class RealType>
inline typename Distribution::value_type quantile(const complemented2_type<Distribution, RealType>& c)
inline BOOST_GPU_ENABLED typename Distribution::value_type quantile(const complemented2_type<Distribution, RealType>& c)
{
typedef typename Distribution::value_type value_type;
return quantile(complement(c.dist, static_cast<value_type>(c.param)));
}
template <class Dist>
inline typename Dist::value_type median(const Dist& d)
inline BOOST_GPU_ENABLED typename Dist::value_type median(const Dist& d)
{ // median - default definition for those distributions for which a
// simple closed form is not known,
// and for which a domain_error and/or NaN generating function is NOT defined.

108
test/cuda/arcsine_cdf_double.cu Normal file
View File

@@ -0,0 +1,108 @@
// Copyright John Maddock 2016.
// 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)
#define BOOST_MATH_OVERFLOW_ERROR_POLICY ignore_error
#include <iostream>
#include <iomanip>
#include <boost/math/distributions/arcsine.hpp>
#include <boost/math/special_functions/relative_difference.hpp>
#include <boost/random/mersenne_twister.hpp>
#include <boost/random/uniform_real_distribution.hpp>
#include "cuda_managed_ptr.hpp"
#include "stopwatch.hpp"
// For the CUDA runtime routines (prefixed with "cuda_")
#include <cuda_runtime.h>
typedef double float_type;
/**
* CUDA Kernel Device code
*
*/
__global__ void cuda_test(const float_type *in1, float_type *out, int numElements)
{
using std::cos;
int i = blockDim.x * blockIdx.x + threadIdx.x;
if (i < numElements)
{
out[i] = cdf(boost::math::arcsine_distribution<float_type>(), in1[i]);
}
}
/**
* Host main routine
*/
int main(void)
{
try{
// Error code to check return values for CUDA calls
cudaError_t err = cudaSuccess;
// Print the vector length to be used, and compute its size
int numElements = 50000;
std::cout << "[Vector operation on " << numElements << " elements]" << std::endl;
// Allocate the managed input vector A
cuda_managed_ptr<float_type> input_vector1(numElements);
// Allocate the managed output vector C
cuda_managed_ptr<float_type> output_vector(numElements);
boost::random::mt19937 gen;
boost::random::uniform_real_distribution<float_type> dist;
// Initialize the input vectors
for (int i = 0; i < numElements; ++i)
{
input_vector1[i] = dist(gen);
}
// Launch the Vector Add CUDA Kernel
int threadsPerBlock = 512;
int blocksPerGrid =(numElements + threadsPerBlock - 1) / threadsPerBlock;
std::cout << "CUDA kernel launch with " << blocksPerGrid << " blocks of " << threadsPerBlock << " threads" << std::endl;
watch w;
cuda_test<<<blocksPerGrid, threadsPerBlock>>>(input_vector1.get(), output_vector.get(), numElements);
std::cout << "CUDA kernal done in " << w.elapsed() << "s" << std::endl;
err = cudaGetLastError();
if (err != cudaSuccess)
{
std::cerr << "Failed to launch vectorAdd kernel (error code " << cudaGetErrorString(err) << ")!" << std::endl;
return EXIT_FAILURE;
}
// Verify that the result vector is correct
std::vector<float_type> results;
results.reserve(numElements);
w.reset();
for(int i = 0; i < numElements; ++i)
results.push_back(cdf(boost::math::arcsine_distribution<float_type>(), input_vector1[i]));
double t = w.elapsed();
// check the results
for(int i = 0; i < numElements; ++i)
{
if (boost::math::epsilon_difference(output_vector[i], results[i]) > 100.0)
{
std::cerr << "Result verification failed at element " << i << "!" << std::endl;
std::cerr << "Error rate was: " << boost::math::epsilon_difference(output_vector[i], results[i]) << "eps" << std::endl;
return EXIT_FAILURE;
}
}
std::cout << "Test PASSED with calculation time: " << t << "s" << std::endl;
std::cout << "Done\n";
}
catch(const std::exception& e)
{
std::cerr << "Stopped with exception: " << e.what() << std::endl;
}
return 0;
}

108
test/cuda/arcsine_cdf_float.cu Normal file
View File

@@ -0,0 +1,108 @@
// Copyright John Maddock 2016.
// 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)
#define BOOST_MATH_OVERFLOW_ERROR_POLICY ignore_error
#include <iostream>
#include <iomanip>
#include <boost/math/distributions/arcsine.hpp>
#include <boost/math/special_functions/relative_difference.hpp>
#include <boost/random/mersenne_twister.hpp>
#include <boost/random/uniform_real_distribution.hpp>
#include "cuda_managed_ptr.hpp"
#include "stopwatch.hpp"
// For the CUDA runtime routines (prefixed with "cuda_")
#include <cuda_runtime.h>
typedef float float_type;
/**
* CUDA Kernel Device code
*
*/
__global__ void cuda_test(const float_type *in1, float_type *out, int numElements)
{
using std::cos;
int i = blockDim.x * blockIdx.x + threadIdx.x;
if (i < numElements)
{
out[i] = cdf(boost::math::arcsine_distribution<float_type>(), in1[i]);
}
}
/**
* Host main routine
*/
int main(void)
{
try{
// Error code to check return values for CUDA calls
cudaError_t err = cudaSuccess;
// Print the vector length to be used, and compute its size
int numElements = 50000;
std::cout << "[Vector operation on " << numElements << " elements]" << std::endl;
// Allocate the managed input vector A
cuda_managed_ptr<float_type> input_vector1(numElements);
// Allocate the managed output vector C
cuda_managed_ptr<float_type> output_vector(numElements);
boost::random::mt19937 gen;
boost::random::uniform_real_distribution<float_type> dist;
// Initialize the input vectors
for (int i = 0; i < numElements; ++i)
{
input_vector1[i] = dist(gen);
}
// Launch the Vector Add CUDA Kernel
int threadsPerBlock = 512;
int blocksPerGrid =(numElements + threadsPerBlock - 1) / threadsPerBlock;
std::cout << "CUDA kernel launch with " << blocksPerGrid << " blocks of " << threadsPerBlock << " threads" << std::endl;
watch w;
cuda_test<<<blocksPerGrid, threadsPerBlock>>>(input_vector1.get(), output_vector.get(), numElements);
std::cout << "CUDA kernal done in " << w.elapsed() << "s" << std::endl;
err = cudaGetLastError();
if (err != cudaSuccess)
{
std::cerr << "Failed to launch vectorAdd kernel (error code " << cudaGetErrorString(err) << ")!" << std::endl;
return EXIT_FAILURE;
}
// Verify that the result vector is correct
std::vector<float_type> results;
results.reserve(numElements);
w.reset();
for(int i = 0; i < numElements; ++i)
results.push_back(cdf(boost::math::arcsine_distribution<float_type>(), input_vector1[i]));
double t = w.elapsed();
// check the results
for(int i = 0; i < numElements; ++i)
{
if (boost::math::epsilon_difference(output_vector[i], results[i]) > 100.0)
{
std::cerr << "Result verification failed at element " << i << "!" << std::endl;
std::cerr << "Error rate was: " << boost::math::epsilon_difference(output_vector[i], results[i]) << "eps" << std::endl;
return EXIT_FAILURE;
}
}
std::cout << "Test PASSED with calculation time: " << t << "s" << std::endl;
std::cout << "Done\n";
}
catch(const std::exception& e)
{
std::cerr << "Stopped with exception: " << e.what() << std::endl;
}
return 0;
}

108
test/cuda/arcsine_pdf_double.cu Normal file
View File

@@ -0,0 +1,108 @@
// Copyright John Maddock 2016.
// 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)
#define BOOST_MATH_OVERFLOW_ERROR_POLICY ignore_error
#include <iostream>
#include <iomanip>
#include <boost/math/distributions/arcsine.hpp>
#include <boost/math/special_functions/relative_difference.hpp>
#include <boost/random/mersenne_twister.hpp>
#include <boost/random/uniform_real_distribution.hpp>
#include "cuda_managed_ptr.hpp"
#include "stopwatch.hpp"
// For the CUDA runtime routines (prefixed with "cuda_")
#include <cuda_runtime.h>
typedef double float_type;
/**
* CUDA Kernel Device code
*
*/
__global__ void cuda_test(const float_type *in1, float_type *out, int numElements)
{
using std::cos;
int i = blockDim.x * blockIdx.x + threadIdx.x;
if (i < numElements)
{
out[i] = pdf(boost::math::arcsine_distribution<float_type>(), in1[i]);
}
}
/**
* Host main routine
*/
int main(void)
{
try{
// Error code to check return values for CUDA calls
cudaError_t err = cudaSuccess;
// Print the vector length to be used, and compute its size
int numElements = 50000;
std::cout << "[Vector operation on " << numElements << " elements]" << std::endl;
// Allocate the managed input vector A
cuda_managed_ptr<float_type> input_vector1(numElements);
// Allocate the managed output vector C
cuda_managed_ptr<float_type> output_vector(numElements);
boost::random::mt19937 gen;
boost::random::uniform_real_distribution<float_type> dist;
// Initialize the input vectors
for (int i = 0; i < numElements; ++i)
{
input_vector1[i] = dist(gen);
}
// Launch the Vector Add CUDA Kernel
int threadsPerBlock = 512;
int blocksPerGrid =(numElements + threadsPerBlock - 1) / threadsPerBlock;
std::cout << "CUDA kernel launch with " << blocksPerGrid << " blocks of " << threadsPerBlock << " threads" << std::endl;
watch w;
cuda_test<<<blocksPerGrid, threadsPerBlock>>>(input_vector1.get(), output_vector.get(), numElements);
std::cout << "CUDA kernal done in " << w.elapsed() << "s" << std::endl;
err = cudaGetLastError();
if (err != cudaSuccess)
{
std::cerr << "Failed to launch vectorAdd kernel (error code " << cudaGetErrorString(err) << ")!" << std::endl;
return EXIT_FAILURE;
}
// Verify that the result vector is correct
std::vector<float_type> results;
results.reserve(numElements);
w.reset();
for(int i = 0; i < numElements; ++i)
results.push_back(pdf(boost::math::arcsine_distribution<float_type>(), input_vector1[i]));
double t = w.elapsed();
// check the results
for(int i = 0; i < numElements; ++i)
{
if (boost::math::epsilon_difference(output_vector[i], results[i]) > 100.0)
{
std::cerr << "Result verification failed at element " << i << "!" << std::endl;
std::cerr << "Error rate was: " << boost::math::epsilon_difference(output_vector[i], results[i]) << "eps" << std::endl;
return EXIT_FAILURE;
}
}
std::cout << "Test PASSED with calculation time: " << t << "s" << std::endl;
std::cout << "Done\n";
}
catch(const std::exception& e)
{
std::cerr << "Stopped with exception: " << e.what() << std::endl;
}
return 0;
}

108
test/cuda/arcsine_pdf_float.cu Normal file
View File

@@ -0,0 +1,108 @@
// Copyright John Maddock 2016.
// 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)
#define BOOST_MATH_OVERFLOW_ERROR_POLICY ignore_error
#include <iostream>
#include <iomanip>
#include <boost/math/distributions/arcsine.hpp>
#include <boost/math/special_functions/relative_difference.hpp>
#include <boost/random/mersenne_twister.hpp>
#include <boost/random/uniform_real_distribution.hpp>
#include "cuda_managed_ptr.hpp"
#include "stopwatch.hpp"
// For the CUDA runtime routines (prefixed with "cuda_")
#include <cuda_runtime.h>
typedef float float_type;
/**
* CUDA Kernel Device code
*
*/
__global__ void cuda_test(const float_type *in1, float_type *out, int numElements)
{
using std::cos;
int i = blockDim.x * blockIdx.x + threadIdx.x;
if (i < numElements)
{
out[i] = pdf(boost::math::arcsine_distribution<float_type>(), in1[i]);
}
}
/**
* Host main routine
*/
int main(void)
{
try{
// Error code to check return values for CUDA calls
cudaError_t err = cudaSuccess;
// Print the vector length to be used, and compute its size
int numElements = 50000;
std::cout << "[Vector operation on " << numElements << " elements]" << std::endl;
// Allocate the managed input vector A
cuda_managed_ptr<float_type> input_vector1(numElements);
// Allocate the managed output vector C
cuda_managed_ptr<float_type> output_vector(numElements);
boost::random::mt19937 gen;
boost::random::uniform_real_distribution<float_type> dist;
// Initialize the input vectors
for (int i = 0; i < numElements; ++i)
{
input_vector1[i] = dist(gen);
}
// Launch the Vector Add CUDA Kernel
int threadsPerBlock = 512;
int blocksPerGrid =(numElements + threadsPerBlock - 1) / threadsPerBlock;
std::cout << "CUDA kernel launch with " << blocksPerGrid << " blocks of " << threadsPerBlock << " threads" << std::endl;
watch w;
cuda_test<<<blocksPerGrid, threadsPerBlock>>>(input_vector1.get(), output_vector.get(), numElements);
std::cout << "CUDA kernal done in " << w.elapsed() << "s" << std::endl;
err = cudaGetLastError();
if (err != cudaSuccess)
{
std::cerr << "Failed to launch vectorAdd kernel (error code " << cudaGetErrorString(err) << ")!" << std::endl;
return EXIT_FAILURE;
}
// Verify that the result vector is correct
std::vector<float_type> results;
results.reserve(numElements);
w.reset();
for(int i = 0; i < numElements; ++i)
results.push_back(pdf(boost::math::arcsine_distribution<float_type>(), input_vector1[i]));
double t = w.elapsed();
// check the results
for(int i = 0; i < numElements; ++i)
{
if (boost::math::epsilon_difference(output_vector[i], results[i]) > 100.0)
{
std::cerr << "Result verification failed at element " << i << "!" << std::endl;
std::cerr << "Error rate was: " << boost::math::epsilon_difference(output_vector[i], results[i]) << "eps" << std::endl;
return EXIT_FAILURE;
}
}
std::cout << "Test PASSED with calculation time: " << t << "s" << std::endl;
std::cout << "Done\n";
}
catch(const std::exception& e)
{
std::cerr << "Stopped with exception: " << e.what() << std::endl;
}
return 0;
}

108
test/cuda/arcsine_quan_double.cu Normal file
View File

@@ -0,0 +1,108 @@
// Copyright John Maddock 2016.
// 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)
#define BOOST_MATH_OVERFLOW_ERROR_POLICY ignore_error
#include <iostream>
#include <iomanip>
#include <boost/math/distributions/arcsine.hpp>
#include <boost/math/special_functions/relative_difference.hpp>
#include <boost/random/mersenne_twister.hpp>
#include <boost/random/uniform_real_distribution.hpp>
#include "cuda_managed_ptr.hpp"
#include "stopwatch.hpp"
// For the CUDA runtime routines (prefixed with "cuda_")
#include <cuda_runtime.h>
typedef double float_type;
/**
* CUDA Kernel Device code
*
*/
__global__ void cuda_test(const float_type *in1, float_type *out, int numElements)
{
using std::cos;
int i = blockDim.x * blockIdx.x + threadIdx.x;
if (i < numElements)
{
out[i] = quantile(boost::math::arcsine_distribution<float_type>(), in1[i]);
}
}
/**
* Host main routine
*/
int main(void)
{
try{
// Error code to check return values for CUDA calls
cudaError_t err = cudaSuccess;
// Print the vector length to be used, and compute its size
int numElements = 50000;
std::cout << "[Vector operation on " << numElements << " elements]" << std::endl;
// Allocate the managed input vector A
cuda_managed_ptr<float_type> input_vector1(numElements);
// Allocate the managed output vector C
cuda_managed_ptr<float_type> output_vector(numElements);
boost::random::mt19937 gen;
boost::random::uniform_real_distribution<float_type> dist;
// Initialize the input vectors
for (int i = 0; i < numElements; ++i)
{
input_vector1[i] = dist(gen);
}
// Launch the Vector Add CUDA Kernel
int threadsPerBlock = 512;
int blocksPerGrid =(numElements + threadsPerBlock - 1) / threadsPerBlock;
std::cout << "CUDA kernel launch with " << blocksPerGrid << " blocks of " << threadsPerBlock << " threads" << std::endl;
watch w;
cuda_test<<<blocksPerGrid, threadsPerBlock>>>(input_vector1.get(), output_vector.get(), numElements);
std::cout << "CUDA kernal done in " << w.elapsed() << "s" << std::endl;
err = cudaGetLastError();
if (err != cudaSuccess)
{
std::cerr << "Failed to launch vectorAdd kernel (error code " << cudaGetErrorString(err) << ")!" << std::endl;
return EXIT_FAILURE;
}
// Verify that the result vector is correct
std::vector<float_type> results;
results.reserve(numElements);
w.reset();
for(int i = 0; i < numElements; ++i)
results.push_back(quantile(boost::math::arcsine_distribution<float_type>(), input_vector1[i]));
double t = w.elapsed();
// check the results
for(int i = 0; i < numElements; ++i)
{
if (boost::math::epsilon_difference(output_vector[i], results[i]) > 100.0)
{
std::cerr << "Result verification failed at element " << i << "!" << std::endl;
std::cerr << "Error rate was: " << boost::math::epsilon_difference(output_vector[i], results[i]) << "eps" << std::endl;
return EXIT_FAILURE;
}
}
std::cout << "Test PASSED with calculation time: " << t << "s" << std::endl;
std::cout << "Done\n";
}
catch(const std::exception& e)
{
std::cerr << "Stopped with exception: " << e.what() << std::endl;
}
return 0;
}

108
test/cuda/arcsine_quan_float.cu Normal file
View File

@@ -0,0 +1,108 @@
// Copyright John Maddock 2016.
// 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)
#define BOOST_MATH_OVERFLOW_ERROR_POLICY ignore_error
#include <iostream>
#include <iomanip>
#include <boost/math/distributions/arcsine.hpp>
#include <boost/math/special_functions/relative_difference.hpp>
#include <boost/random/mersenne_twister.hpp>
#include <boost/random/uniform_real_distribution.hpp>
#include "cuda_managed_ptr.hpp"
#include "stopwatch.hpp"
// For the CUDA runtime routines (prefixed with "cuda_")
#include <cuda_runtime.h>
typedef float float_type;
/**
* CUDA Kernel Device code
*
*/
__global__ void cuda_test(const float_type *in1, float_type *out, int numElements)
{
using std::cos;
int i = blockDim.x * blockIdx.x + threadIdx.x;
if (i < numElements)
{
out[i] = quantile(boost::math::arcsine_distribution<float_type>(), in1[i]);
}
}
/**
* Host main routine
*/
int main(void)
{
try{
// Error code to check return values for CUDA calls
cudaError_t err = cudaSuccess;
// Print the vector length to be used, and compute its size
int numElements = 50000;
std::cout << "[Vector operation on " << numElements << " elements]" << std::endl;
// Allocate the managed input vector A
cuda_managed_ptr<float_type> input_vector1(numElements);
// Allocate the managed output vector C
cuda_managed_ptr<float_type> output_vector(numElements);
boost::random::mt19937 gen;
boost::random::uniform_real_distribution<float_type> dist;
// Initialize the input vectors
for (int i = 0; i < numElements; ++i)
{
input_vector1[i] = dist(gen);
}
// Launch the Vector Add CUDA Kernel
int threadsPerBlock = 512;
int blocksPerGrid =(numElements + threadsPerBlock - 1) / threadsPerBlock;
std::cout << "CUDA kernel launch with " << blocksPerGrid << " blocks of " << threadsPerBlock << " threads" << std::endl;
watch w;
cuda_test<<<blocksPerGrid, threadsPerBlock>>>(input_vector1.get(), output_vector.get(), numElements);
std::cout << "CUDA kernal done in " << w.elapsed() << "s" << std::endl;
err = cudaGetLastError();
if (err != cudaSuccess)
{
std::cerr << "Failed to launch vectorAdd kernel (error code " << cudaGetErrorString(err) << ")!" << std::endl;
return EXIT_FAILURE;
}
// Verify that the result vector is correct
std::vector<float_type> results;
results.reserve(numElements);
w.reset();
for(int i = 0; i < numElements; ++i)
results.push_back(quantile(boost::math::arcsine_distribution<float_type>(), input_vector1[i]));
double t = w.elapsed();
// check the results
for(int i = 0; i < numElements; ++i)
{
if (boost::math::epsilon_difference(output_vector[i], results[i]) > 100.0)
{
std::cerr << "Result verification failed at element " << i << "!" << std::endl;
std::cerr << "Error rate was: " << boost::math::epsilon_difference(output_vector[i], results[i]) << "eps" << std::endl;
return EXIT_FAILURE;
}
}
std::cout << "Test PASSED with calculation time: " << t << "s" << std::endl;
std::cout << "Done\n";
}
catch(const std::exception& e)
{
std::cerr << "Stopped with exception: " << e.what() << std::endl;
}
return 0;
}

108
test/cuda/cauchy_cdf_double.cu Normal file
View File

@@ -0,0 +1,108 @@
// Copyright John Maddock 2016.
// 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)
#define BOOST_MATH_OVERFLOW_ERROR_POLICY ignore_error
#include <iostream>
#include <iomanip>
#include <boost/math/distributions/cauchy.hpp>
#include <boost/math/special_functions/relative_difference.hpp>
#include <boost/random/mersenne_twister.hpp>
#include <boost/random/uniform_real_distribution.hpp>
#include "cuda_managed_ptr.hpp"
#include "stopwatch.hpp"
// For the CUDA runtime routines (prefixed with "cuda_")
#include <cuda_runtime.h>
typedef double float_type;
/**
* CUDA Kernel Device code
*
*/
__global__ void cuda_test(const float_type *in1, float_type *out, int numElements)
{
using std::cos;
int i = blockDim.x * blockIdx.x + threadIdx.x;
if (i < numElements)
{
out[i] = cdf(boost::math::cauchy_distribution<float_type>(), in1[i]);
}
}
/**
* Host main routine
*/
int main(void)
{
try{
// Error code to check return values for CUDA calls
cudaError_t err = cudaSuccess;
// Print the vector length to be used, and compute its size
int numElements = 50000;
std::cout << "[Vector operation on " << numElements << " elements]" << std::endl;
// Allocate the managed input vector A
cuda_managed_ptr<float_type> input_vector1(numElements);
// Allocate the managed output vector C
cuda_managed_ptr<float_type> output_vector(numElements);
boost::random::mt19937 gen;
boost::random::uniform_real_distribution<float_type> dist(-10000, 10000);
// Initialize the input vectors
for (int i = 0; i < numElements; ++i)
{
input_vector1[i] = dist(gen);
}
// Launch the Vector Add CUDA Kernel
int threadsPerBlock = 512;
int blocksPerGrid =(numElements + threadsPerBlock - 1) / threadsPerBlock;
std::cout << "CUDA kernel launch with " << blocksPerGrid << " blocks of " << threadsPerBlock << " threads" << std::endl;
watch w;
cuda_test<<<blocksPerGrid, threadsPerBlock>>>(input_vector1.get(), output_vector.get(), numElements);
std::cout << "CUDA kernal done in " << w.elapsed() << "s" << std::endl;
err = cudaGetLastError();
if (err != cudaSuccess)
{
std::cerr << "Failed to launch vectorAdd kernel (error code " << cudaGetErrorString(err) << ")!" << std::endl;
return EXIT_FAILURE;
}
// Verify that the result vector is correct
std::vector<float_type> results;
results.reserve(numElements);
w.reset();
for(int i = 0; i < numElements; ++i)
results.push_back(cdf(boost::math::cauchy_distribution<float_type>(), input_vector1[i]));
double t = w.elapsed();
// check the results
for(int i = 0; i < numElements; ++i)
{
if (boost::math::epsilon_difference(output_vector[i], results[i]) > 100.0)
{
std::cerr << "Result verification failed at element " << i << "!" << std::endl;
std::cerr << "Error rate was: " << boost::math::epsilon_difference(output_vector[i], results[i]) << "eps" << std::endl;
return EXIT_FAILURE;
}
}
std::cout << "Test PASSED with calculation time: " << t << "s" << std::endl;
std::cout << "Done\n";
}
catch(const std::exception& e)
{
std::cerr << "Stopped with exception: " << e.what() << std::endl;
}
return 0;
}

108
test/cuda/cauchy_cdf_float.cu Normal file
View File

@@ -0,0 +1,108 @@
// Copyright John Maddock 2016.
// 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)
#define BOOST_MATH_OVERFLOW_ERROR_POLICY ignore_error
#include <iostream>
#include <iomanip>
#include <boost/math/distributions/cauchy.hpp>
#include <boost/math/special_functions/relative_difference.hpp>
#include <boost/random/mersenne_twister.hpp>
#include <boost/random/uniform_real_distribution.hpp>
#include "cuda_managed_ptr.hpp"
#include "stopwatch.hpp"
// For the CUDA runtime routines (prefixed with "cuda_")
#include <cuda_runtime.h>
typedef float float_type;
/**
* CUDA Kernel Device code
*
*/
__global__ void cuda_test(const float_type *in1, float_type *out, int numElements)
{
using std::cos;
int i = blockDim.x * blockIdx.x + threadIdx.x;
if (i < numElements)
{
out[i] = cdf(boost::math::cauchy_distribution<float_type>(), in1[i]);
}
}
/**
* Host main routine
*/
int main(void)
{
try{
// Error code to check return values for CUDA calls
cudaError_t err = cudaSuccess;
// Print the vector length to be used, and compute its size
int numElements = 50000;
std::cout << "[Vector operation on " << numElements << " elements]" << std::endl;
// Allocate the managed input vector A
cuda_managed_ptr<float_type> input_vector1(numElements);
// Allocate the managed output vector C
cuda_managed_ptr<float_type> output_vector(numElements);
boost::random::mt19937 gen;
boost::random::uniform_real_distribution<float_type> dist(-10000, 10000);
// Initialize the input vectors
for (int i = 0; i < numElements; ++i)
{
input_vector1[i] = dist(gen);
}
// Launch the Vector Add CUDA Kernel
int threadsPerBlock = 512;
int blocksPerGrid =(numElements + threadsPerBlock - 1) / threadsPerBlock;
std::cout << "CUDA kernel launch with " << blocksPerGrid << " blocks of " << threadsPerBlock << " threads" << std::endl;
watch w;
cuda_test<<<blocksPerGrid, threadsPerBlock>>>(input_vector1.get(), output_vector.get(), numElements);
std::cout << "CUDA kernal done in " << w.elapsed() << "s" << std::endl;
err = cudaGetLastError();
if (err != cudaSuccess)
{
std::cerr << "Failed to launch vectorAdd kernel (error code " << cudaGetErrorString(err) << ")!" << std::endl;
return EXIT_FAILURE;
}
// Verify that the result vector is correct
std::vector<float_type> results;
results.reserve(numElements);
w.reset();
for(int i = 0; i < numElements; ++i)
results.push_back(cdf(boost::math::cauchy_distribution<float_type>(), input_vector1[i]));
double t = w.elapsed();
// check the results
for(int i = 0; i < numElements; ++i)
{
if (boost::math::epsilon_difference(output_vector[i], results[i]) > 100.0)
{
std::cerr << "Result verification failed at element " << i << "!" << std::endl;
std::cerr << "Error rate was: " << boost::math::epsilon_difference(output_vector[i], results[i]) << "eps" << std::endl;
return EXIT_FAILURE;
}
}
std::cout << "Test PASSED with calculation time: " << t << "s" << std::endl;
std::cout << "Done\n";
}
catch(const std::exception& e)
{
std::cerr << "Stopped with exception: " << e.what() << std::endl;
}
return 0;
}

108
test/cuda/cauchy_pdf_double.cu Normal file
View File

@@ -0,0 +1,108 @@
// Copyright John Maddock 2016.
// 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)
#define BOOST_MATH_OVERFLOW_ERROR_POLICY ignore_error
#include <iostream>
#include <iomanip>
#include <boost/math/distributions/cauchy.hpp>
#include <boost/math/special_functions/relative_difference.hpp>
#include <boost/random/mersenne_twister.hpp>
#include <boost/random/uniform_real_distribution.hpp>
#include "cuda_managed_ptr.hpp"
#include "stopwatch.hpp"
// For the CUDA runtime routines (prefixed with "cuda_")
#include <cuda_runtime.h>
typedef double float_type;
/**
* CUDA Kernel Device code
*
*/
__global__ void cuda_test(const float_type *in1, float_type *out, int numElements)
{
using std::cos;
int i = blockDim.x * blockIdx.x + threadIdx.x;
if (i < numElements)
{
out[i] = pdf(boost::math::cauchy_distribution<float_type>(), in1[i]);
}
}
/**
* Host main routine
*/
int main(void)
{
try{
// Error code to check return values for CUDA calls
cudaError_t err = cudaSuccess;
// Print the vector length to be used, and compute its size
int numElements = 50000;
std::cout << "[Vector operation on " << numElements << " elements]" << std::endl;
// Allocate the managed input vector A
cuda_managed_ptr<float_type> input_vector1(numElements);
// Allocate the managed output vector C
cuda_managed_ptr<float_type> output_vector(numElements);
boost::random::mt19937 gen;
boost::random::uniform_real_distribution<float_type> dist(-10000, 10000);
// Initialize the input vectors
for (int i = 0; i < numElements; ++i)
{
input_vector1[i] = dist(gen);
}
// Launch the Vector Add CUDA Kernel
int threadsPerBlock = 512;
int blocksPerGrid =(numElements + threadsPerBlock - 1) / threadsPerBlock;
std::cout << "CUDA kernel launch with " << blocksPerGrid << " blocks of " << threadsPerBlock << " threads" << std::endl;
watch w;
cuda_test<<<blocksPerGrid, threadsPerBlock>>>(input_vector1.get(), output_vector.get(), numElements);
std::cout << "CUDA kernal done in " << w.elapsed() << "s" << std::endl;
err = cudaGetLastError();
if (err != cudaSuccess)
{
std::cerr << "Failed to launch vectorAdd kernel (error code " << cudaGetErrorString(err) << ")!" << std::endl;
return EXIT_FAILURE;
}
// Verify that the result vector is correct
std::vector<float_type> results;
results.reserve(numElements);
w.reset();
for(int i = 0; i < numElements; ++i)
results.push_back(pdf(boost::math::cauchy_distribution<float_type>(), input_vector1[i]));
double t = w.elapsed();
// check the results
for(int i = 0; i < numElements; ++i)
{
if (boost::math::epsilon_difference(output_vector[i], results[i]) > 100.0)
{
std::cerr << "Result verification failed at element " << i << "!" << std::endl;
std::cerr << "Error rate was: " << boost::math::epsilon_difference(output_vector[i], results[i]) << "eps" << std::endl;
return EXIT_FAILURE;
}
}
std::cout << "Test PASSED with calculation time: " << t << "s" << std::endl;
std::cout << "Done\n";
}
catch(const std::exception& e)
{
std::cerr << "Stopped with exception: " << e.what() << std::endl;
}
return 0;
}

108
test/cuda/cauchy_pdf_float.cu Normal file
View File

@@ -0,0 +1,108 @@
// Copyright John Maddock 2016.
// 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)
#define BOOST_MATH_OVERFLOW_ERROR_POLICY ignore_error
#include <iostream>
#include <iomanip>
#include <boost/math/distributions/cauchy.hpp>
#include <boost/math/special_functions/relative_difference.hpp>
#include <boost/random/mersenne_twister.hpp>
#include <boost/random/uniform_real_distribution.hpp>
#include "cuda_managed_ptr.hpp"
#include "stopwatch.hpp"
// For the CUDA runtime routines (prefixed with "cuda_")
#include <cuda_runtime.h>
typedef float float_type;
/**
* CUDA Kernel Device code
*
*/
__global__ void cuda_test(const float_type *in1, float_type *out, int numElements)
{
using std::cos;
int i = blockDim.x * blockIdx.x + threadIdx.x;
if (i < numElements)
{
out[i] = pdf(boost::math::cauchy_distribution<float_type>(), in1[i]);
}
}
/**
* Host main routine
*/
int main(void)
{
try{
// Error code to check return values for CUDA calls
cudaError_t err = cudaSuccess;
// Print the vector length to be used, and compute its size
int numElements = 50000;
std::cout << "[Vector operation on " << numElements << " elements]" << std::endl;
// Allocate the managed input vector A
cuda_managed_ptr<float_type> input_vector1(numElements);
// Allocate the managed output vector C
cuda_managed_ptr<float_type> output_vector(numElements);
boost::random::mt19937 gen;
boost::random::uniform_real_distribution<float_type> dist(-10000, 10000);
// Initialize the input vectors
for (int i = 0; i < numElements; ++i)
{
input_vector1[i] = dist(gen);
}
// Launch the Vector Add CUDA Kernel
int threadsPerBlock = 512;
int blocksPerGrid =(numElements + threadsPerBlock - 1) / threadsPerBlock;
std::cout << "CUDA kernel launch with " << blocksPerGrid << " blocks of " << threadsPerBlock << " threads" << std::endl;
watch w;
cuda_test<<<blocksPerGrid, threadsPerBlock>>>(input_vector1.get(), output_vector.get(), numElements);
std::cout << "CUDA kernal done in " << w.elapsed() << "s" << std::endl;
err = cudaGetLastError();
if (err != cudaSuccess)
{
std::cerr << "Failed to launch vectorAdd kernel (error code " << cudaGetErrorString(err) << ")!" << std::endl;
return EXIT_FAILURE;
}
// Verify that the result vector is correct
std::vector<float_type> results;
results.reserve(numElements);
w.reset();
for(int i = 0; i < numElements; ++i)
results.push_back(pdf(boost::math::cauchy_distribution<float_type>(), input_vector1[i]));
double t = w.elapsed();
// check the results
for(int i = 0; i < numElements; ++i)
{
if (boost::math::epsilon_difference(output_vector[i], results[i]) > 100.0)
{
std::cerr << "Result verification failed at element " << i << "!" << std::endl;
std::cerr << "Error rate was: " << boost::math::epsilon_difference(output_vector[i], results[i]) << "eps" << std::endl;
return EXIT_FAILURE;
}
}
std::cout << "Test PASSED with calculation time: " << t << "s" << std::endl;
std::cout << "Done\n";
}
catch(const std::exception& e)
{
std::cerr << "Stopped with exception: " << e.what() << std::endl;
}
return 0;
}

108
test/cuda/cauchy_quan_double.cu Normal file
View File

@@ -0,0 +1,108 @@
// Copyright John Maddock 2016.
// 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)
#define BOOST_MATH_OVERFLOW_ERROR_POLICY ignore_error
#include <iostream>
#include <iomanip>
#include <boost/math/distributions/cauchy.hpp>
#include <boost/math/special_functions/relative_difference.hpp>
#include <boost/random/mersenne_twister.hpp>
#include <boost/random/uniform_real_distribution.hpp>
#include "cuda_managed_ptr.hpp"
#include "stopwatch.hpp"
// For the CUDA runtime routines (prefixed with "cuda_")
#include <cuda_runtime.h>
typedef double float_type;
/**
* CUDA Kernel Device code
*
*/
__global__ void cuda_test(const float_type *in1, float_type *out, int numElements)
{
using std::cos;
int i = blockDim.x * blockIdx.x + threadIdx.x;
if (i < numElements)
{
out[i] = quantile(boost::math::cauchy_distribution<float_type>(), in1[i]);
}
}
/**
* Host main routine
*/
int main(void)
{
try{
// Error code to check return values for CUDA calls
cudaError_t err = cudaSuccess;
// Print the vector length to be used, and compute its size
int numElements = 50000;
std::cout << "[Vector operation on " << numElements << " elements]" << std::endl;
// Allocate the managed input vector A
cuda_managed_ptr<float_type> input_vector1(numElements);
// Allocate the managed output vector C
cuda_managed_ptr<float_type> output_vector(numElements);
boost::random::mt19937 gen;
boost::random::uniform_real_distribution<float_type> dist;
// Initialize the input vectors
for (int i = 0; i < numElements; ++i)
{
input_vector1[i] = dist(gen);
}
// Launch the Vector Add CUDA Kernel
int threadsPerBlock = 512;
int blocksPerGrid =(numElements + threadsPerBlock - 1) / threadsPerBlock;
std::cout << "CUDA kernel launch with " << blocksPerGrid << " blocks of " << threadsPerBlock << " threads" << std::endl;
watch w;
cuda_test<<<blocksPerGrid, threadsPerBlock>>>(input_vector1.get(), output_vector.get(), numElements);
std::cout << "CUDA kernal done in " << w.elapsed() << "s" << std::endl;
err = cudaGetLastError();
if (err != cudaSuccess)
{
std::cerr << "Failed to launch vectorAdd kernel (error code " << cudaGetErrorString(err) << ")!" << std::endl;
return EXIT_FAILURE;
}
// Verify that the result vector is correct
std::vector<float_type> results;
results.reserve(numElements);
w.reset();
for(int i = 0; i < numElements; ++i)
results.push_back(quantile(boost::math::cauchy_distribution<float_type>(), input_vector1[i]));
double t = w.elapsed();
// check the results
for(int i = 0; i < numElements; ++i)
{
if (boost::math::epsilon_difference(output_vector[i], results[i]) > 100.0)
{
std::cerr << "Result verification failed at element " << i << "!" << std::endl;
std::cerr << "Error rate was: " << boost::math::epsilon_difference(output_vector[i], results[i]) << "eps" << std::endl;
return EXIT_FAILURE;
}
}
std::cout << "Test PASSED with calculation time: " << t << "s" << std::endl;
std::cout << "Done\n";
}
catch(const std::exception& e)
{
std::cerr << "Stopped with exception: " << e.what() << std::endl;
}
return 0;
}

108
test/cuda/cauchy_quan_float.cu Normal file
View File

@@ -0,0 +1,108 @@
// Copyright John Maddock 2016.
// 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)
#define BOOST_MATH_OVERFLOW_ERROR_POLICY ignore_error
#include <iostream>
#include <iomanip>
#include <boost/math/distributions/cauchy.hpp>
#include <boost/math/special_functions/relative_difference.hpp>
#include <boost/random/mersenne_twister.hpp>
#include <boost/random/uniform_real_distribution.hpp>
#include "cuda_managed_ptr.hpp"
#include "stopwatch.hpp"
// For the CUDA runtime routines (prefixed with "cuda_")
#include <cuda_runtime.h>
typedef float float_type;
/**
* CUDA Kernel Device code
*
*/
__global__ void cuda_test(const float_type *in1, float_type *out, int numElements)
{
using std::cos;
int i = blockDim.x * blockIdx.x + threadIdx.x;
if (i < numElements)
{
out[i] = quantile(boost::math::cauchy_distribution<float_type>(), in1[i]);
}
}
/**
* Host main routine
*/
int main(void)
{
try{
// Error code to check return values for CUDA calls
cudaError_t err = cudaSuccess;
// Print the vector length to be used, and compute its size
int numElements = 50000;
std::cout << "[Vector operation on " << numElements << " elements]" << std::endl;
// Allocate the managed input vector A
cuda_managed_ptr<float_type> input_vector1(numElements);
// Allocate the managed output vector C
cuda_managed_ptr<float_type> output_vector(numElements);
boost::random::mt19937 gen;
boost::random::uniform_real_distribution<float_type> dist;
// Initialize the input vectors
for (int i = 0; i < numElements; ++i)
{
input_vector1[i] = dist(gen);
}
// Launch the Vector Add CUDA Kernel
int threadsPerBlock = 512;
int blocksPerGrid =(numElements + threadsPerBlock - 1) / threadsPerBlock;
std::cout << "CUDA kernel launch with " << blocksPerGrid << " blocks of " << threadsPerBlock << " threads" << std::endl;
watch w;
cuda_test<<<blocksPerGrid, threadsPerBlock>>>(input_vector1.get(), output_vector.get(), numElements);
std::cout << "CUDA kernal done in " << w.elapsed() << "s" << std::endl;
err = cudaGetLastError();
if (err != cudaSuccess)
{
std::cerr << "Failed to launch vectorAdd kernel (error code " << cudaGetErrorString(err) << ")!" << std::endl;
return EXIT_FAILURE;
}
// Verify that the result vector is correct
std::vector<float_type> results;
results.reserve(numElements);
w.reset();
for(int i = 0; i < numElements; ++i)
results.push_back(quantile(boost::math::cauchy_distribution<float_type>(), input_vector1[i]));
double t = w.elapsed();
// check the results
for(int i = 0; i < numElements; ++i)
{
if (boost::math::epsilon_difference(output_vector[i], results[i]) > 100.0)
{
std::cerr << "Result verification failed at element " << i << "!" << std::endl;
std::cerr << "Error rate was: " << boost::math::epsilon_difference(output_vector[i], results[i]) << "eps" << std::endl;
return EXIT_FAILURE;
}
}
std::cout << "Test PASSED with calculation time: " << t << "s" << std::endl;
std::cout << "Done\n";
}
catch(const std::exception& e)
{
std::cerr << "Stopped with exception: " << e.what() << std::endl;
}
return 0;
}