>Andrew Bennett
>Sat, 11 Mar 2006 17:24:51 -0800
>
>On Sat, 11 Mar 2006 13:44:52 +0100, Pianoman wrote:
>
>>Hi, I need to perform very fast squareroot (sqrt) operations on double
>> type nukmbers.
>>I tried this
>>function mysqrt(a:double):double;
>> [...]
>>yn:=(y*y+a)/(2*y);
>>until abs(yn-y) < 10e-16;
>>mysqrt:=yn;
>>end;
>>but it is quite slow (FPC uses Heron iterations) which is 10 times
>>faster
>>than this code but I need even faster sqrt routin.
>>Any ideas for optimizing the function written above would be
>>appreciated.
>>>Pianoman
>>Once upon a time, I did this in single precision:
>
>1) I made a crude 1st guess by manufacturing, in
>Assembler, a linear fit mantissa - 2 different fits
>depending if exponent was odd or even - and the
>appropriate exponent. The worst case error of this
>is easily calculated. Process takes 2 floating ops and
>some integer stuff taking out and putting back the
>exponent.
>
>2) I then explicitly coded the iterations required to
>converge in this worst case: just 3 floating ops per
>(or is it 2 plus integer exponent manipulation?)
>and it didn't take very many. No test. Saving the test
>saved more time than was wasted doing unnecessary
>iterations, at least on the machine I did it on. This
>may not be so in double precision ...
For some reason I recently have been interested in fast sqrt computation
as well, and found the following (C) pseudo-code, which looks quite fast.
It approximates 1/sqrt(x), which can be done without a divide, just
multiplies. Maybe it is of help, don't know the Heron iteration thing
you mention...
(See especially section "B. sqrt(x) by Reciproot Iteration",
http://www.netlib.org/fdlibm/e_sqrt.c)
Regards,
Thomas
_______________________________________________
fpc-pascal maillist - fpc-pascal@lists.freepascal.org
http://lists.freepascal.org/mailman/listinfo/fpc-pascal
