SF [2159499] Full bigfloat precision for Gamma after the second call



Dieter Kaiser wrote:
> As reported on Sourceforge.net I have studied the subtle bug that we do not get
> the full precision of the Gamma function in the first call. I can reproduce the
> bug with GCL 2.6.8 and CLISP 2.44.
> 
> I have studied the effect in more detail and detected that we have the
> problem only for a negative argument to bffac. That is true for
> gamma(0.5b0).
> 
> The reflection formula calculates:
> 
> bfloat(%pi*z/sin(%pi*z))/bffac(-z,fpprec)
> 

I looked into this a bit.  It seems to me that the problem is actually
in the value of pi.  I put a break point at *fpsin, and this is what I
see, when we change the precision from 64 to 128.  The arg to sin is:

((MTIMES SIMP)
       ((BIGFLOAT SIMP 428)

346583711765101857447301773017885462929554634421977071896309947576827663475703202879996800763017447262173901370175446478621769728
        0)
       $%PI))

The bigfloat there works out to be exactly 0.5b-1, as expected:

(displa '((BIGFLOAT SIMP 428)

346583711765101857447301773017885462929554634421977071896309947576827663475703202879996800763017447262173901370175446478621769728
        0))
5.0b-1

But look at the value of bigfloat%pi

((BIGFLOAT SIMP 428)
 544412421367563193456719682845547953364731116819812684167733550941339415157339997517605497896021017270535788777222242612081295775
 2)
0] (displa '((BIGFLOAT SIMP 428)
 544412421367563193456719682845547953364731116819812684167733550941339415157339997517605497896021017270535788777222242612081295775
 2))
3.1415926535897932384626433832795028841971693993751058209749445923078164062862\
089986280348256812002485034480173261900300090541701b0

Compared to the actual value of %pi:

(%i5) bfloat(%pi);
(%o5)
3.1415926535897932384626433832795028841971693993751058209749445923078164\
062862089986280348253421170679821480865132823066470938446b0

I think the problem is that when fpprec is set, fpprec1 is called.  But
this only updates the bfloat versions of 1, 0, +/- 1/2.  Perhaps we need
to update the bigfloat values of %pi at least, and probably also %e and
%gamma.

Ray