New insns for the s390 backend (2)
Bradley Lucier
lucier@math.purdue.edu
Thu Aug 28 02:20:00 GMT 2003
On Wednesday, August 27, 2003, at 04:43 PM, Roger Sayle wrote:
> I did however manage to find the following URL of interest:
> http://www.active-web.cc/html/research/sine/sin-cos.txt
> All of the algorithms for precomputing FFT sin/cos tables given
> on that page appear to be unaffected by -ffast-math.
Roger, look at the accuracy requirements for sin/cos in Percival's
paper as I applied it to my code---each entry can be off by four
rounding errors, max, in tables with tens of millions of values. Your
reference is not as complete as
http://www.dattalo.com/technical/theory/sinewave.html
The error for method 4 here grows linearly; the one for method 5 grows
quadratically. If you want to use the best method I know of
sin(a + da) = cos(da) * sin(a) + sin(da) * cos(a)
= sin(a) + ((cos(da)-1) * sin(a) + sin(da) * cos(a))
= sin(a) + ((-2*sin^2(da/2)) * sin(a) + sin(da) * cos(a))
and similarly for cos, you can make a probabilistic argument that the
error grows like the square root of the number of iterations. So, *if*
you have 80 bit intermediate values you *probably* have < 4 rounding
errors at double precision with tens of millions of table entries (and
on PPC, SPARC, ..., you don't have that extra precision). But that
doesn't cut it for the number theory folks---they (I) need guaranteed <
4 rounding error accuracy, so I just do library calls. On x86 I could
use the above iteration in blocks of 2048 table entries *if* I could be
assured that gcc wouldn't spill any temporaries to doubles, but then
I'd still have jumps in values at the boundaries of the blocks, and I
just gave up with that.
> Then just to double check, I've just built and run "make check"
> on GNU gmp-4.1.2 as configured by default and then again with
> CC="gcc -ffast-math". No failures/differences. This is especially
> relevant to your concerns, as mpn/generic/mul_fft.c contains code to
> perform bignum multiplications using Karatsuba, 3-way Toom-Cook,
> and Fermat FFT.
I don't believe that gmp uses floating-point FFTs at all, see the code
in mpn/generic/mul_fft.c, they use FFTs in other fields. My impression
of the gmp folks is that they don't trust floating-point arithmetic,
numerical analysis (or analysts ;-), or mathematical libraries, which
is fine, but I'm not surprised to find out that -ffast-math doesn't
affect the test results.
I think this is getting way afield of the topic of gcc-patches; perhaps
we should take any more discussion off list.
Brad
More information about the Gcc-patches
mailing list