Fri 14 Mar 2014 06:20:06 PM UTC, comment #4:
Caveat: I am not experienced enough to be a maintainer for this code and have just looked into this particular problem.
It needs a more professional eye.
This was analyzed with ddd/gdb on Linux Mint.
I have looked into this and the problem seems to be in:
hyperg_U.c
In this section:
gsl_sf_result sum;
int stat_sum = hyperg_U_finite_sum(N, a, b, x, xeps, &sum);
int stat_inf;
/* Evaluate infinite sum. */
if(fabs(xeps-1.0) > 0.5 ) {
stat_inf = hyperg_U_infinite_sum_stable(N, a, bint, b, beps, x, xeps, sum, result);
} else if (1+a-b < 0 && 1+a-b==floor(1+a-b) && beps != 0) {
stat_inf = hyperg_U_infinite_sum_simple(N, a, bint, b, beps, x, xeps, sum, result);
} else {
stat_inf = hyperg_U_infinite_sum_improved(N, a, bint, b, beps, x, xeps, sum, result);
}
What happens is that when a-b+1=-n then a finite condition 1-b+1= -n is satisfied and the finite sum succeeds; as it supposed to.
But the program continues to the hyperg_U_infinite_sum_simple where that condition is a no..no and bombs in gamma.c deep inside so it never gets back to check if the finite sum worked.
What should happen is the success of the finite sum should skip the latter tests
Recommendation: If the finite sum succeeded skip the infinite sum tests (and execution).
|
Tue 11 Mar 2014 07:54:28 PM UTC, comment #2:
Sorry for so many emails... I also tried the GSL SHELL software and it
seemse that sf.hypergU(1, -500, 0.1) can be calculated correctly. However,
I found that any negative non-integers of b for U(a, b, x) cannot be
calculated. These negative non-integers should be suitable arguments and
the manual states that b is of type double.
|
Tue 11 Mar 2014 07:54:03 PM UTC, comment #1:
Sorry for the typo but for example 1) listed in the previous email:
gsl_sf_hyperg_U(1, -500, 0.1) gives the correct value: ~1.99e-3.
The problem happens when calculating gsl_sf_hyperg_U(1, -500, 40): gsl:
gsl: hyperg.c:165: ERROR: overflow. The value should be ~0.00184 so I
suppose this is a bug.
|
Tue 11 Mar 2014 07:52:19 PM UTC, original submission:
Hi,
I have been using the gsl_sf_hyperg_U function for a while and I found
several bugs when calculating normal values. For example,
1) gsl_sf_hyperg_U(1, -500, 0.1) gives incorrect value: 9.3246e+03, which
actually should be around 0.19.
2) Many domain error bugs. For example, gsl_sf_hyperg_U(1, -6.67, 1) should
have value ~0.128 but will invoke domain error for gsl.
|