Hypergeometric 2F1() specfunc: cases x=1

Hypergeometric 2F1() specfunc: cases x=1

[Richard Mathar]

It would be nice if the implementation of the Gaussian hypergeometric
functions would catch the special case of unit argument and compute it
via the standard formula with the Gamma functions (Abramowitz & Stegun 15.1.20).
Currently, all x>0.995 are rejects as being too close to the limit where
the branch cut prohibits the standard series to converge. So my proposal is
to implement in gsl_sf_hyperg_2F1_e() in specfunc/hyperg_2F1.c just
at the start of the function a trap for this case:

....
  const int c_neg_integer = ( c < 0.0  &&  fabs(c - rintc) < locEPS );

  result->val = 0.0;
  result->err = 0.0;

/* start proposed handling of x=1.0 RJM */
   if ( fabs(x - 1.0) < locEPS && c-a-b > 0. && c != 0. && ! c_neg_integer )
   {
        gsl_sf_result lnGamc, lnGamcab, lnGamca, lnGamcb  ;
        int stat = gsl_sf_lngamma_e(c, & lnGamc);
        if ( stat != GSL_SUCCESS )
        {
                DOMAIN_ERROR(result);
        }
        stat = gsl_sf_lngamma_e(c-a-b, & lnGamcab);
        if ( stat != GSL_SUCCESS )
        {
                DOMAIN_ERROR(result);
        }
        stat = gsl_sf_lngamma_e(c-a, & lnGamca);
        if ( stat != GSL_SUCCESS )
        {
                DOMAIN_ERROR(result);
        }
        stat = gsl_sf_lngamma_e(c-b, & lnGamcb);
        if ( stat != GSL_SUCCESS )
        {
                DOMAIN_ERROR(result);
        }
        return gsl_sf_exp_err_e(lnGamc.val+lnGamcab.val-lnGamca.val
                -lnGamcb.val,lnGamc.err+lnGamcab.err+lnGamca.err+lnGamcb.err,result) ;
   }
/*end proposed handling of x=1.0 RJM */

  if(x < -1.0 || 1.0 <= x) {
....

Richard Mathar, http://www.strw.leidenuniv.nl/~mathar

Re: Hypergeometric 2F1() specfunc: cases x=1

[Brian Gough]

At Tue, 11 Dec 2007 22:14:27 +0100,
Richard Mathar wrote:
> It would be nice if the implementation of the Gaussian hypergeometric
> functions would catch the special case of unit argument and compute it
> via the standard formula with the Gamma functions (Abramowitz & Stegun 15.1.20).

Sounds reasonable --- can you write some test cases with high
precision values for specfunc/test_hyperg.c to make sure it's working.

-- 
Brian Gough