Re: [Maxima-discuss] e, pi, gamma, log(2)

[Raymond T. <toy...@gm...>]

>>>>> "Volker" == Volker van Nek <vol...@gm...> writes:

    Volker> Am 22.09.2014 18:27, schrieb Raymond Toy:
    >> 
    Volker> fpprec:30$
    Volker> x:0.5b0$
    Volker> asin(x);            --> 5.23598775598298873077107230547b−1
    Volker> bs_asin(x);         --> 5.23598775598298873077107230547b−1
    Volker> sin(asin(x));       --> 5.00000000000000000000000000001b−1
    Volker> bs_sin(bs_asin(x)); --> 5.00000000000000000023968760547b−1
    >> 
    >> Interesting. I'm sure this can be fixed, but it seems like an odd
    >> artifact.
    >> 

    Volker> I fixed it. Now it seems that the bs_ functions are more accurate.

What was the fix? (Did you ever post the code? I can't seem to find
it.)

    >> Can you give an example of the problem? I have not looked at fpround
    >> in a very long time and don't really remember how it works anymore.
    >> 

    Volker> fpprec:30$ x:0.5b0$

    Volker> asin(sin(x));        5.0b−1
    Volker> bs_asin(bs_sin(x));  2.5b−1  BUG


    Volker> The essential lines of code:

    Volker> (defun bs-asincos (x fn)
    Volker> (let ((phi
    Volker> (let* ((fpprec (+ fpprec 12))

    Volker> I leave out some lines here.
    Volker> phi is now computed with 12 extra bits and
    Volker> the following line is the correct result but with 12 bits extra.

    Volker> (print (bs-carg1 re-z im-z fpprec)) )))

    Volker> After returning from the block of higher precision I want to round back
    Volker> to the original precision:

    Volker> (list (print (fpround (car phi))) (print (cadr phi))) )) ;
    Volker> BUG

I would look at what bigfloatp does in the case where fpprec is lower
than the actual precision of the bigfloat. It adds the value of *m
into the exponent.

    Volker> print 1 : (20769187434139310514121985316880257 -1)  OK
    Volker> print 2 : 2535301200456458802993406410752  OK
    Volker> print 3 : -1  WRONG

    Volker> The exponent must be 0 instead of -1. The binary output shows the
    Volker> reason. fpround rounds from below.

Hmm, that seems not right. I would expect fpround to round, and
hopefully round to even.  Richard might remember how this works.

Ray