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/
