Skip to content

Commit 1bf9657

Browse files
rok-cesnovarrok-cesnovar
authored
Merge pull request #2450 from bstatcomp/opencl_exponential_cdf
Opencl exponential and frechet (lc)cdf
2 parents 7c55ef7 + 55462da commit 1bf9657

19 files changed

Lines changed: 1368 additions & 6 deletions

‎stan/math/opencl/prim.hpp‎

Lines changed: 6 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -136,7 +136,13 @@
136136
#include <stan/math/opencl/prim/dot_self.hpp>
137137
#include <stan/math/opencl/prim/double_exponential_lpdf.hpp>
138138
#include <stan/math/opencl/prim/exp_mod_normal_lpdf.hpp>
139+
#include <stan/math/opencl/prim/exponential_cdf.hpp>
140+
#include <stan/math/opencl/prim/exponential_lccdf.hpp>
141+
#include <stan/math/opencl/prim/exponential_lcdf.hpp>
139142
#include <stan/math/opencl/prim/exponential_lpdf.hpp>
143+
#include <stan/math/opencl/prim/frechet_cdf.hpp>
144+
#include <stan/math/opencl/prim/frechet_lccdf.hpp>
145+
#include <stan/math/opencl/prim/frechet_lcdf.hpp>
140146
#include <stan/math/opencl/prim/frechet_lpdf.hpp>
141147
#include <stan/math/opencl/prim/gamma_lpdf.hpp>
142148
#include <stan/math/opencl/prim/gp_exp_quad_cov.hpp>

‎stan/math/opencl/prim/exponential_cdf.hpp‎

Lines changed: 102 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,102 @@
1+
#ifndef STAN_MATH_OPENCL_PRIM_EXPONENTIAL_CDF_HPP
2+
#define STAN_MATH_OPENCL_PRIM_EXPONENTIAL_CDF_HPP
3+
#ifdef STAN_OPENCL
4+
5+
#include <stan/math/prim/meta.hpp>
6+
#include <stan/math/prim/err.hpp>
7+
#include <stan/math/prim/fun/constants.hpp>
8+
#include <stan/math/prim/fun/elt_divide.hpp>
9+
#include <stan/math/prim/fun/elt_multiply.hpp>
10+
#include <stan/math/opencl/kernel_generator.hpp>
11+
#include <stan/math/prim/functor/operands_and_partials.hpp>
12+
13+
namespace stan {
14+
namespace math {
15+
16+
/** \ingroup opencl
17+
* Calculates the exponential cumulative distribution function for
18+
* the given y and beta.
19+
*
20+
* Inverse scale parameter must be greater than 0.
21+
* y must be greater than or equal to 0.
22+
*
23+
* @tparam T_y_cl type of scalar
24+
* @tparam T_inv_scale_cl type of inverse scale
25+
* @param y a scalar variable
26+
* @param beta inverse scale parameter
27+
* @return The product of densities
28+
*/
29+
template <typename T_y_cl, typename T_inv_scale_cl,
30+
require_all_prim_or_rev_kernel_expression_t<
31+
T_y_cl, T_inv_scale_cl>* = nullptr,
32+
require_any_not_stan_scalar_t<T_y_cl, T_inv_scale_cl>* = nullptr>
33+
return_type_t<T_y_cl, T_inv_scale_cl> exponential_cdf(
34+
const T_y_cl& y, const T_inv_scale_cl& beta) {
35+
static const char* function = "exponential_cdf(OpenCL)";
36+
using T_partials_return = partials_return_t<T_y_cl, T_inv_scale_cl>;
37+
using std::isfinite;
38+
using std::isnan;
39+
40+
check_consistent_sizes(function, "Random variable", y,
41+
"Inverse scale parameter", beta);
42+
const size_t N = max_size(y, beta);
43+
if (N == 0) {
44+
return 1.0;
45+
}
46+
47+
const auto& y_col = as_column_vector_or_scalar(y);
48+
const auto& beta_col = as_column_vector_or_scalar(beta);
49+
50+
const auto& y_val = value_of(y_col);
51+
const auto& beta_val = value_of(beta_col);
52+
53+
auto check_y_nonnegative
54+
= check_cl(function, "Random variable", y_val, "nonnegative");
55+
auto y_nonnegative_expr = y_val >= 0.0;
56+
auto check_beta_positive_finite = check_cl(
57+
function, "Inverse scale parameter", beta_val, "positive finite");
58+
auto beta_positive_finite_expr = 0.0 < beta_val && isfinite(beta_val);
59+
60+
auto exp_val = exp(elt_multiply(-beta_val, y_val));
61+
auto one_m_exp = 1.0 - exp_val;
62+
auto cdf_expr = colwise_prod(one_m_exp);
63+
auto rep_deriv1 = elt_divide(exp_val, one_m_exp);
64+
65+
matrix_cl<double> cdf_cl;
66+
matrix_cl<double> y_deriv_cl;
67+
matrix_cl<double> beta_deriv_cl;
68+
69+
results(check_y_nonnegative, check_beta_positive_finite, cdf_cl,
70+
beta_deriv_cl)
71+
= expressions(
72+
y_nonnegative_expr, beta_positive_finite_expr, cdf_expr,
73+
calc_if<!is_constant_all<T_y_cl, T_inv_scale_cl>::value>(rep_deriv1));
74+
75+
T_partials_return cdf = (from_matrix_cl(cdf_cl)).prod();
76+
77+
auto rep_deriv2 = beta_deriv_cl * cdf;
78+
auto y_deriv = elt_multiply(
79+
static_select<is_constant<T_y_cl>::value>(0, beta_val), rep_deriv2);
80+
auto beta_deriv = elt_multiply(
81+
static_select<is_constant<T_inv_scale_cl>::value>(0, y_val), rep_deriv2);
82+
83+
results(y_deriv_cl, beta_deriv_cl)
84+
= expressions(calc_if<!is_constant<T_y_cl>::value>(y_deriv),
85+
calc_if<!is_constant<T_inv_scale_cl>::value>(beta_deriv));
86+
87+
operands_and_partials<decltype(y_col), decltype(beta_col)> ops_partials(
88+
y_col, beta_col);
89+
90+
if (!is_constant<T_y_cl>::value) {
91+
ops_partials.edge1_.partials_ = std::move(y_deriv_cl);
92+
}
93+
if (!is_constant<T_inv_scale_cl>::value) {
94+
ops_partials.edge2_.partials_ = std::move(beta_deriv_cl);
95+
}
96+
return ops_partials.build(cdf);
97+
}
98+
99+
} // namespace math
100+
} // namespace stan
101+
#endif
102+
#endif

‎stan/math/opencl/prim/exponential_lccdf.hpp‎

Lines changed: 84 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,84 @@
1+
#ifndef STAN_MATH_OPENCL_PRIM_EXPONENTIAL_LCCDF_HPP
2+
#define STAN_MATH_OPENCL_PRIM_EXPONENTIAL_LCCDF_HPP
3+
#ifdef STAN_OPENCL
4+
5+
#include <stan/math/prim/meta.hpp>
6+
#include <stan/math/prim/err.hpp>
7+
#include <stan/math/prim/fun/constants.hpp>
8+
#include <stan/math/prim/fun/elt_divide.hpp>
9+
#include <stan/math/prim/fun/elt_multiply.hpp>
10+
#include <stan/math/opencl/kernel_generator.hpp>
11+
#include <stan/math/prim/functor/operands_and_partials.hpp>
12+
13+
namespace stan {
14+
namespace math {
15+
16+
/** \ingroup opencl
17+
* Calculates the log exponential cumulative distribution function for
18+
* the given y and beta.
19+
*
20+
* Inverse scale parameter must be greater than 0.
21+
* y must be greater than or equal to 0.
22+
*
23+
* @tparam T_y_cl type of scalar
24+
* @tparam T_inv_scale_cl type of inverse scale
25+
* @param y a scalar variable
26+
* @param beta inverse scale parameter
27+
* @return The log of the product of densities
28+
*/
29+
template <typename T_y_cl, typename T_inv_scale_cl,
30+
require_all_prim_or_rev_kernel_expression_t<
31+
T_y_cl, T_inv_scale_cl>* = nullptr,
32+
require_any_not_stan_scalar_t<T_y_cl, T_inv_scale_cl>* = nullptr>
33+
return_type_t<T_y_cl, T_inv_scale_cl> exponential_lccdf(
34+
const T_y_cl& y, const T_inv_scale_cl& beta) {
35+
static const char* function = "exponential_lccdf(OpenCL)";
36+
using T_partials_return = partials_return_t<T_y_cl, T_inv_scale_cl>;
37+
using std::isfinite;
38+
using std::isnan;
39+
40+
check_consistent_sizes(function, "Random variable", y,
41+
"Inverse scale parameter", beta);
42+
const size_t N = max_size(y, beta);
43+
if (N == 0) {
44+
return 0.0;
45+
}
46+
47+
const auto& y_col = as_column_vector_or_scalar(y);
48+
const auto& beta_col = as_column_vector_or_scalar(beta);
49+
50+
const auto& y_val = value_of(y_col);
51+
const auto& beta_val = value_of(beta_col);
52+
53+
auto check_y_nonnegative
54+
= check_cl(function, "Random variable", y_val, "nonnegative");
55+
auto y_nonnegative_expr = y_val >= 0.0;
56+
auto check_beta_positive_finite = check_cl(
57+
function, "Inverse scale parameter", beta_val, "positive finite");
58+
auto beta_positive_finite_expr = 0.0 < beta_val && isfinite(beta_val);
59+
60+
auto lccdf_expr = colwise_sum(elt_multiply(-beta_val, y_val));
61+
62+
matrix_cl<double> lccdf_cl;
63+
64+
results(check_y_nonnegative, check_beta_positive_finite, lccdf_cl)
65+
= expressions(y_nonnegative_expr, beta_positive_finite_expr, lccdf_expr);
66+
67+
T_partials_return lccdf = (from_matrix_cl(lccdf_cl)).sum();
68+
69+
operands_and_partials<decltype(y_col), decltype(beta_col)> ops_partials(
70+
y_col, beta_col);
71+
72+
if (!is_constant<T_y_cl>::value) {
73+
ops_partials.edge1_.partials_ = constant(0.0, N, 1) - beta_val;
74+
}
75+
if (!is_constant<T_inv_scale_cl>::value) {
76+
ops_partials.edge2_.partials_ = constant(0.0, N, 1) - y_val;
77+
}
78+
return ops_partials.build(lccdf);
79+
}
80+
81+
} // namespace math
82+
} // namespace stan
83+
#endif
84+
#endif

‎stan/math/opencl/prim/exponential_lcdf.hpp‎

Lines changed: 94 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,94 @@
1+
#ifndef STAN_MATH_OPENCL_PRIM_EXPONENTIAL_LCDF_HPP
2+
#define STAN_MATH_OPENCL_PRIM_EXPONENTIAL_LCDF_HPP
3+
#ifdef STAN_OPENCL
4+
5+
#include <stan/math/prim/meta.hpp>
6+
#include <stan/math/prim/err.hpp>
7+
#include <stan/math/prim/fun/constants.hpp>
8+
#include <stan/math/prim/fun/elt_divide.hpp>
9+
#include <stan/math/prim/fun/elt_multiply.hpp>
10+
#include <stan/math/opencl/kernel_generator.hpp>
11+
#include <stan/math/prim/functor/operands_and_partials.hpp>
12+
13+
namespace stan {
14+
namespace math {
15+
16+
/** \ingroup opencl
17+
* Calculates the log exponential cumulative distribution function for
18+
* the given y and beta.
19+
*
20+
* Inverse scale parameter must be greater than 0.
21+
* y must be greater than or equal to 0.
22+
*
23+
* @tparam T_y_cl type of scalar
24+
* @tparam T_inv_scale_cl type of inverse scale
25+
* @param y a scalar variable
26+
* @param beta inverse scale parameter
27+
* @return The log of the product of densities
28+
*/
29+
template <typename T_y_cl, typename T_inv_scale_cl,
30+
require_all_prim_or_rev_kernel_expression_t<
31+
T_y_cl, T_inv_scale_cl>* = nullptr,
32+
require_any_not_stan_scalar_t<T_y_cl, T_inv_scale_cl>* = nullptr>
33+
return_type_t<T_y_cl, T_inv_scale_cl> exponential_lcdf(
34+
const T_y_cl& y, const T_inv_scale_cl& beta) {
35+
static const char* function = "exponential_lcdf(OpenCL)";
36+
using T_partials_return = partials_return_t<T_y_cl, T_inv_scale_cl>;
37+
using std::isfinite;
38+
using std::isnan;
39+
40+
check_consistent_sizes(function, "Random variable", y,
41+
"Inverse scale parameter", beta);
42+
const size_t N = max_size(y, beta);
43+
if (N == 0) {
44+
return 0.0;
45+
}
46+
47+
const auto& y_col = as_column_vector_or_scalar(y);
48+
const auto& beta_col = as_column_vector_or_scalar(beta);
49+
50+
const auto& y_val = value_of(y_col);
51+
const auto& beta_val = value_of(beta_col);
52+
53+
auto check_y_nonnegative
54+
= check_cl(function, "Random variable", y_val, "nonnegative");
55+
auto y_nonnegative_expr = y_val >= 0.0;
56+
auto check_beta_positive_finite = check_cl(
57+
function, "Inverse scale parameter", beta_val, "positive finite");
58+
auto beta_positive_finite_expr = 0.0 < beta_val && isfinite(beta_val);
59+
60+
auto exp_val = exp(elt_multiply(-beta_val, y_val));
61+
auto lcdf_expr = colwise_sum(log1m(exp_val));
62+
63+
auto rep_deriv = elt_divide(exp_val, 1.0 - exp_val);
64+
auto y_deriv = elt_multiply(beta_val, rep_deriv);
65+
auto beta_deriv = elt_multiply(y_val, rep_deriv);
66+
67+
matrix_cl<double> lcdf_cl;
68+
matrix_cl<double> y_deriv_cl;
69+
matrix_cl<double> beta_deriv_cl;
70+
71+
results(check_y_nonnegative, check_beta_positive_finite, lcdf_cl, y_deriv_cl,
72+
beta_deriv_cl)
73+
= expressions(y_nonnegative_expr, beta_positive_finite_expr, lcdf_expr,
74+
calc_if<!is_constant<T_y_cl>::value>(y_deriv),
75+
calc_if<!is_constant<T_inv_scale_cl>::value>(beta_deriv));
76+
77+
T_partials_return lcdf = (from_matrix_cl(lcdf_cl)).sum();
78+
79+
operands_and_partials<decltype(y_col), decltype(beta_col)> ops_partials(
80+
y_col, beta_col);
81+
82+
if (!is_constant<T_y_cl>::value) {
83+
ops_partials.edge1_.partials_ = std::move(y_deriv);
84+
}
85+
if (!is_constant<T_inv_scale_cl>::value) {
86+
ops_partials.edge2_.partials_ = std::move(beta_deriv);
87+
}
88+
return ops_partials.build(lcdf);
89+
}
90+
91+
} // namespace math
92+
} // namespace stan
93+
#endif
94+
#endif

0 commit comments

Comments
 (0)