1 #include "../inc/util.h"
4 namespace vertex_lcfi {
namespace util{
6 double gamma(
double a,
double x)
8 if (a <= 0 || x <= 0)
return 0;
10 if (x < (a+1))
return gamSer(a,x);
11 else return gamCf(a,x);
18 double fpmin = 1.e-30;
20 if (a <= 0 || x <= 0)
return 0;
28 for (
int i=1; i<=itmax; i++) {
29 an = double(-i)*(double(i)-a);
32 if (fabs(d) < fpmin) d = fpmin;
34 if (fabs(c) < fpmin) c = fpmin;
38 if (fabs(del-1) < eps)
break;
41 double v = exp(-x+a*log(x)-gln)*h;
51 if (a <= 0 || x <= 0)
return 0;
57 for (
int n=1; n<=itmax; n++) {
61 if ( fabs(del) < fabs(sum*eps))
break;
64 double v = sum*exp(-x+a*log(x)-gln);
73 double c[7] = { 2.5066282746310005, 76.18009172947146, -86.50532032941677
74 ,24.01409824083091, -1.231739572450155, 0.1208650973866179e-2
75 ,-0.5395239384953e-5};
80 tmp = (x+0.5)*log(tmp)-tmp;
81 double ser = 1.000000000190015;
82 for (
int i=1; i<7; i++) {
86 double v = tmp+log(c[0]*ser/x);
double gamCf(double a, double x)
Computation of the incomplete gamma function P(a,x) via its continued fraction representation.
double gamSer(double a, double x)
Incomplete gamma function P(a,x) via its series representation.
double lnGamma(double z)
Computation of ln[gamma(z)] for all z>0.
double gamma(double a, double x)
Computation of the upper incomplete gamma function P(a,x)