This is the mail archive of the
gsl-discuss@sources.redhat.com
mailing list for the GSL project.
Bug fix (proposed): Incomplete gamma function, large parameters equal to each other
- From: gsl at boonstra dot org
- To: gsl-discuss at sources dot redhat dot com
- Date: Wed, 9 Feb 2005 12:47:06 -0800 (PST)
- Subject: Bug fix (proposed): Incomplete gamma function, large parameters equal to each other
Hi all
Special thanks to Gerard Jungman for his work on these functions.
I'd like to propose a bugfix for gamma_inc_Q_asymp_unif() in
specfunc/gamma_inc.c near line 150. Note that
const double eps = (x-a)/a;
but soon after we set
eta = eps * sqrt(-2.0*ln_term.val/(eps*eps))
which of course does not work if eps is zero (or effectively zero). It
seems to me the fix is simple, just by inserting something like:
if (eps<GSL_SQRT_DBL_EPSILON) {
eta = sqrt(-2.0*ln_term.val);
}
but I must say I'm far from sure GSL_SQRT_DBL_EPSILON is the correct bound.
Brian