1 #ifndef STAN_MATH_REV_SCAL_FUN_IBETA_HPP
2 #define STAN_MATH_REV_SCAL_FUN_IBETA_HPP
4 #include <boost/math/special_functions/digamma.hpp>
5 #include <boost/math/special_functions/gamma.hpp>
18 double ibeta_hypergeometric_helper(
double a,
double b,
double z,
19 double precision = 1
e-8,
20 double max_steps = 1000) {
37 class ibeta_vvv_vari :
public op_vvv_vari {
39 ibeta_vvv_vari(vari* avi, vari* bvi, vari* xvi) :
40 op_vvv_vari(
stan::math::
ibeta(avi->val_, bvi->val_, xvi->val_),
44 double a = avi_->val_;
45 double b = bvi_->val_;
46 double c = cvi_->val_;
55 using stan::math::ibeta_hypergeometric_helper;
60 * ibeta_hypergeometric_helper(a, 1-b, c);
63 * ibeta_hypergeometric_helper(b, 1-a, 1-c)
67 boost::math::ibeta_derivative(a, b, c);
70 class ibeta_vvd_vari :
public op_vvd_vari {
72 ibeta_vvd_vari(vari* avi, vari* bvi,
double x) :
73 op_vvd_vari(
stan::math::
ibeta(avi->val_, bvi->val_, x), avi, bvi, x) {
76 double a = avi_->val_;
77 double b = bvi_->val_;
87 using stan::math::ibeta_hypergeometric_helper;
92 * ibeta_hypergeometric_helper(a, 1-b, c);
95 * ibeta_hypergeometric_helper(b, 1-a, 1-c)
100 class ibeta_vdv_vari :
public op_vdv_vari {
102 ibeta_vdv_vari(vari* avi,
double b, vari* xvi) :
103 op_vdv_vari(
stan::math::
ibeta(avi->val_, b, xvi->val_), avi, b, xvi) {
106 double a = avi_->val_;
108 double c = cvi_->val_;
117 using stan::math::ibeta_hypergeometric_helper;
122 * ibeta_hypergeometric_helper(a, 1-b, c);
124 boost::math::ibeta_derivative(a, b, c);
127 class ibeta_vdd_vari :
public op_vdd_vari {
129 ibeta_vdd_vari(vari* avi,
double b,
double x) :
130 op_vdd_vari(
stan::math::
ibeta(avi->val_, b, x), avi, b, x) {
133 double a = avi_->val_;
144 using stan::math::ibeta_hypergeometric_helper;
149 * ibeta_hypergeometric_helper(a, 1-b, c);
152 class ibeta_dvv_vari :
public op_dvv_vari {
154 ibeta_dvv_vari(
double a, vari* bvi, vari* xvi) :
155 op_dvv_vari(
stan::math::
ibeta(a, bvi->val_, xvi->val_), a, bvi, xvi) {
159 double b = bvi_->val_;
160 double c = cvi_->val_;
169 using stan::math::ibeta_hypergeometric_helper;
172 * ibeta_hypergeometric_helper(b, 1-a, 1-c)
176 boost::math::ibeta_derivative(a, b, c);
179 class ibeta_dvd_vari :
public op_dvd_vari {
181 ibeta_dvd_vari(
double a, vari* bvi,
double x) :
182 op_dvd_vari(
stan::math::
ibeta(a, bvi->val_, x), a, bvi, x) {
186 double b = bvi_->val_;
196 using stan::math::ibeta_hypergeometric_helper;
199 * ibeta_hypergeometric_helper(b, 1-a, 1-c)
204 class ibeta_ddv_vari :
public op_ddv_vari {
206 ibeta_ddv_vari(
double a,
double b, vari* xvi) :
207 op_ddv_vari(
stan::math::
ibeta(a, b, xvi->val_), a, b, xvi) {
212 double c = cvi_->val_;
215 boost::math::ibeta_derivative(a, b, c);
fvar< T > abs(const fvar< T > &x)
double ibeta(const double a, const double b, const double x)
The normalized incomplete beta function of a, b, and x.
fvar< T > log(const fvar< T > &x)
Independent (input) and dependent (output) variables for gradients.
int isnan(const stan::math::var &a)
Checks if the given number is NaN.
fvar< T > sin(const fvar< T > &x)
vari * vi_
Pointer to the implementation of this variable.
double e()
Return the base of the natural logarithm.
var ibeta(const var &a, const var &b, const var &x)
The normalized incomplete beta function of a, b, and x.
double pi()
Return the value of pi.
double adj_
The adjoint of this variable, which is the partial derivative of this variable with respect to the ro...
fvar< T > pow(const fvar< T > &x1, const fvar< T > &x2)
fvar< T > tgamma(const fvar< T > &x)
fvar< T > digamma(const fvar< T > &x)