Stan Math Library  2.12.0
reverse mode automatic differentiation
gamma_p.hpp
Go to the documentation of this file.
1 #ifndef STAN_MATH_FWD_SCAL_FUN_GAMMA_P_HPP
2 #define STAN_MATH_FWD_SCAL_FUN_GAMMA_P_HPP
3 
4 #include <stan/math/fwd/core.hpp>
5 
7 
8 namespace stan {
9  namespace math {
10 
11  template <typename T>
12  inline
13  fvar<T>
14  gamma_p(const fvar<T>& x1, const fvar<T>& x2) {
15  using std::log;
16  using std::exp;
17  using std::pow;
18  using std::fabs;
19  using boost::math::tgamma;
21 
22  T u = gamma_p(x1.val_, x2.val_);
23 
24  T S = 0;
25  T s = 1;
26  T l = log(x2.val_);
27  T g = tgamma(x1.val_);
28  T dig = digamma(x1.val_);
29 
30  int k = 0;
31  T delta = s / (x1.val_ * x1.val_);
32 
33  while (fabs(delta) > 1e-6) {
34  S += delta;
35  ++k;
36  s *= -x2.val_ / k;
37  delta = s / ((k + x1.val_) * (k + x1.val_));
38  }
39 
40  T der1 = u * (dig - l) + exp(x1.val_ * l) * S / g;
41  T der2 = exp(-x2.val_) * pow(x2.val_, x1.val_ - 1.0) / g;
42 
43  return fvar<T>(u, x1.d_ * -der1 + x2.d_ * der2);
44  }
45 
46  template <typename T>
47  inline
48  fvar<T>
49  gamma_p(const fvar<T>& x1, const double x2) {
50  using std::log;
51  using std::exp;
52  using std::pow;
53  using std::fabs;
54  using boost::math::tgamma;
56 
57  T u = gamma_p(x1.val_, x2);
58 
59  T S = 0.0;
60  double s = 1.0;
61  double l = log(x2);
62  T g = tgamma(x1.val_);
63  T dig = digamma(x1.val_);
64 
65  int k = 0;
66  T delta = s / (x1.val_ * x1.val_);
67 
68  while (fabs(delta) > 1e-6) {
69  S += delta;
70  ++k;
71  s *= -x2 / k;
72  delta = s / ((k + x1.val_) * (k + x1.val_));
73  }
74 
75  T der1 = u * (dig - l) + exp(x1.val_ * l) * S / g;
76 
77  return fvar<T>(u, x1.d_ * -der1);
78  }
79 
80  template <typename T>
81  inline
82  fvar<T>
83  gamma_p(const double x1, const fvar<T>& x2) {
84  using std::exp;
85  using std::pow;
86 
87  T u = gamma_p(x1, x2.val_);
88 
89  double g = boost::math::tgamma(x1);
90 
91  T der2 = exp(-x2.val_) * pow(x2.val_, x1 - 1.0) / g;
92 
93  return fvar<T>(u, x2.d_ * der2);
94  }
95  }
96 }
97 #endif
fvar< T > fabs(const fvar< T > &x)
Definition: fabs.hpp:15
fvar< T > log(const fvar< T > &x)
Definition: log.hpp:14
fvar< T > exp(const fvar< T > &x)
Definition: exp.hpp:10
fvar< T > gamma_p(const fvar< T > &x1, const fvar< T > &x2)
Definition: gamma_p.hpp:14
double e()
Return the base of the natural logarithm.
Definition: constants.hpp:94
fvar< T > pow(const fvar< T > &x1, const fvar< T > &x2)
Definition: pow.hpp:17
fvar< T > tgamma(const fvar< T > &x)
Definition: tgamma.hpp:14
fvar< T > digamma(const fvar< T > &x)
Definition: digamma.hpp:15

     [ Stan Home Page ] © 2011–2016, Stan Development Team.