[Bug-gsl] incorrect values by bessel_steed_array

[Hu Zhan]

Hi There,

When I use gsl_sf_bessel_jl_steed_array(LMAX, x, wrk) with a large LMAX,
e.g. LMAX = 80, the function returns wrong values at small l and x. My
system is Debian Linux with Opteron CPU.

Attached is a program demonstrating this problem. It writes five files
for l = 0 to 5 in the format of 

x  value_by_gsl_sf_bessel_jl_steed_array  value_by_gsl_sf_bessel_jl.

One can easily see the difference between the two functions. 

Hope to see a patch soon. Thanks.

Hu


#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <gsl/gsl_sf_bessel.h>

#define LMAX 80
#define LMAXP1 (LMAX + 1)

int main(void)
{
  char fName[100];
  int l;
  double x, dx = 0.1, wrk[LMAXP1];
  FILE *fp;

  for (l = 0; l < 6; l++) {
    sprintf(fName, "jl%d.txt", l);
    fp = fopen(fName, "w");
    for (x = 0.; x < 10.; x += dx) {
      gsl_sf_bessel_jl_steed_array(LMAX, x, wrk);
      fprintf(fp, "%g %g %g\n", x, wrk[l], gsl_sf_bessel_jl(l, x));
    }
    fclose(fp);
  }

  return 0;
}