[Bug-gsl] lnpoch_pos (specfunc) does not handle large inputs properly

[Reed A. Cartwright]

lnpoch_pos(a,x) calculates Log[Pochhammer[a,x]], where a >0 and a+x > 0

However, for some situations it does not calculate the proper value due to
improper conditionals.

lnpoch_pos(2**60,20) should return 831.77; however, it returns 0.

This is because the first if statement is true, which (eventually) results
in lnpoch_pos calculating the result via lngamma(a+x) - lngamma(a), which
will be wrong because a+x == a due to machine precision.

lnpoch_pos does appear to have a code (the stirling approximation) for this
situation.  This code path is used if `absx < 0.1*a && a > 15.0`.  However
it is not reached because earlier, `absx*log(GSL_MAX_DBL(a,2.0)) > 0.1` is
tested and found to be true.

I'm not sure why the earlier test exists, but it does seem like it needs to
be modified so that the stirling approximation will be used for appropriate
inputs.

 --
Reed A. Cartwright, PhD
Barrett Honors Faculty
Assistant Professor of Genomics, Evolution, and Bioinformatics
School of Life Sciences
Human and Comparative Genomics Laboratory
The Biodesign Institute
Arizona State University
==================
Address: The Biodesign Institute, PO Box 875301, Tempe, AZ 85287-5301 USA
Packages: The Biodesign Institute, 1001 S. McAllister Ave, Tempe, AZ
85287-5301 USA
Office: Biodesign A-224A, 1-480-965-9949
Website: http://cartwrig.ht/