1 #ifndef STAN_MATH_PRIM_SCAL_FUN_LOWER_REG_INC_GAMMA_HPP 2 #define STAN_MATH_PRIM_SCAL_FUN_LOWER_REG_INC_GAMMA_HPP 105 template <
typename T1,
typename T2>
107 const T1& a,
const T2& z,
double precision = 1
e-10,
int max_steps = 1e5) {
114 return std::numeric_limits<TP>::quiet_NaN();
122 if ((a < 0.8 && z > 15.0) || (a < 12.0 && z > 30.0)
136 T1 lgamma_a_plus_1 =
lgamma(a + 1);
137 T1 lgamma_a_plus_n_plus_1 = lgamma_a_plus_1;
140 term =
exp(a_plus_n * log_z - lgamma_a_plus_n_plus_1);
142 if (term <= precision)
145 domain_error(
"grad_reg_lower_inc_gamma",
"n (internal counter)",
146 max_steps,
"exceeded ",
147 " iterations, gamma_p(a,z) gradient (a) " 148 "did not converge.");
150 lgamma_a_plus_n_plus_1 +=
log1p(a_plus_n);
156 TP sum_b =
digamma(a + 1) *
exp(a * log_z - lgamma_a_plus_1);
157 lgamma_a_plus_n_plus_1 = lgamma_a_plus_1 +
log(a_plus_n);
159 term =
exp(a_plus_n * log_z - lgamma_a_plus_n_plus_1)
162 if (term <= precision)
163 return emz * (log_z * sum_a - sum_b);
165 domain_error(
"grad_reg_lower_inc_gamma",
"n (internal counter)",
166 max_steps,
"exceeded ",
167 " iterations, gamma_p(a,z) gradient (a) " 168 "did not converge.");
170 lgamma_a_plus_n_plus_1 +=
log1p(a_plus_n);
fvar< T > lgamma(const fvar< T > &x)
Return the natural logarithm of the gamma function applied to the specified argument.
fvar< T > sqrt(const fvar< T > &x)
double value_of_rec(const fvar< T > &v)
Return the value of the specified variable.
fvar< T > log(const fvar< T > &x)
bool is_any_nan(const T &x)
Returns true if the input is NaN and false otherwise.
void check_positive_finite(const char *function, const char *name, const T_y &y)
Check if y is positive and finite.
return_type< T1, T2 >::type grad_reg_inc_gamma(T1 a, T2 z, T1 g, T1 dig, double precision=1e-6, int max_steps=1e5)
Gradient of the regularized incomplete gamma functions igamma(a, z)
boost::math::tools::promote_args< double, typename scalar_type< T >::type, typename return_type< Types_pack... >::type >::type type
fvar< T > exp(const fvar< T > &x)
void domain_error(const char *function, const char *name, const T &y, const char *msg1, const char *msg2)
Throw a domain error with a consistently formatted message.
double e()
Return the base of the natural logarithm.
fvar< T > log1p(const fvar< T > &x)
fvar< T > pow(const fvar< T > &x1, const fvar< T > &x2)
fvar< T > tgamma(const fvar< T > &x)
Return the result of applying the gamma function to the specified argument.
return_type< T1, T2 >::type grad_reg_lower_inc_gamma(const T1 &a, const T2 &z, double precision=1e-10, int max_steps=1e5)
Computes the gradient of the lower regularized incomplete gamma function.
fvar< T > digamma(const fvar< T > &x)
Return the derivative of the log gamma function at the specified argument.