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

van Nek wrote:
> Am 16 Oct 2008 um 10:35 hat Raymond Toy geschrieben:
>
[snip]
> (defun fprt18231_3 nil 
>   (let ((a 1823176476672000))   
>      (setq a (ash a (ash fpprec 1)))
>      (fpquotient (intofp (isqrt a)) (intofp (ash 1 fpprec))) ))
I think if we really wanted to optimize this, we shouldn't do the
fpquotient at all since we're dividing by a power of two, which is  the
base of the floating point numbers.

I'll have to check, but we might be able to take the output from intofp,
extract the exponent (the second part of the list) and just subtract
fpprec from it to get the correct answer.   That ought to speed things
up somewhat.  Same idea for version 4.
>
> (defun fprt18231_4 nil 
>   (let ((a 1823176476672000))   
>      (setq a (ash a (- (ash fpprec 1) 52))) ;; we already have at
> least 52 bits precision
>      (fpquotient (intofp (isqrt a)) (intofp (ash 1 (- fpprec 26)))) ))
I'm not sure this is quite right.  The value of a might have > 52 bits,
but isqrt will only give 26 bits, since it returns the integer part.  So
if fpprec = 52, we'd only get 26 bits from isqrt.  I think you want to
subtract 26, not 52.
>
> result: version 3 and 4 are about the same concerning time, version 2
> is significantly slower, so I tested version 1 und 3 more in detail.
>
> I always used precisions $fpprec, where fpprec is just 32*2^n + some
> bits, which is the worst case for version 1.
>
> On Windows / GCL:
>
> precs: [616,  1233, 2466,  4932,  9864,  19728,  39456]$
> time1: [0.17, 0.49, 1.22,  3.18,  8.6,   23.92,  70.59]$
> time3: [0.48, 1.17, 3.52, 11.54, 37.43, 116.63, 368.76]$
>
> (in all cases 1000 runs, so time for a single run in millisec)
>
> So version 1 is of type O(x^1.45) and version 3 of type O(x^1.6).
>
> With GCL there is a significant advantage for version 1.
>
> On Ubuntu / clisp (same computer):
>
> precs: [616,      1233,     2466,     4932,     9864,      19728,    
> 39456]$
> time1: [0.356022, 0.748047, 2.016126, 5.820365, 18.08513,  58.195637,
> 195.188199]$
> time3: [0.192013, 0.308018, 0.816051, 2.632164,  9.720607, 37.12232, 
> 148.445278]$
>
> Here version 1 is of type O(x^1.53) and version 3 of type O(x^1.66).
>
> If I assume that one might not always hit the worst case 32*2^n+1 for
> version 1, on clisp both versions are about the same. 1 is better in
> high precs and 3 in small precs.
>
Thanks for do the experiments.  I'll try to optimize version 3 a bit and
run them with cmucl and ecl, just to see if there are significant
differences between implementations.  Cmucl is supposed to have a pretty
fast isqrt, using some algorithm sent to comp.lang.lisp long ago.

Ray