Date: Wed, 1 Aug 2012 18:12:59 +1000 (EST) From: Bruce Evans <brde@optusnet.com.au> To: Stephen Montgomery-Smith <stephen@missouri.edu> Cc: freebsd-numerics@FreeBSD.org, Bruce Evans <brde@optusnet.com.au> Subject: Re: Bad error computing 1/z Message-ID: <20120801180703.B1153@besplex.bde.org> In-Reply-To: <5018441C.5040909@missouri.edu> References: <5015DA8A.10200@missouri.edu> <20120730025118.GB7958@server.rulingia.com> <5015F9E9.9070408@missouri.edu> <20120730131929.Y876@besplex.bde.org> <5018441C.5040909@missouri.edu>
next in thread | previous in thread | raw e-mail | index | archive | help
On Tue, 31 Jul 2012, Stephen Montgomery-Smith wrote: > On 07/29/12 22:35, Bruce Evans wrote: > >> 1/z calculated by gcc is much more accurate than the singular function >> z*z :-). The real part of the latter is x*x-y*y. This may have an >> error of almost 2**MANT_DIG ulps due to cancelation. > > Could one avoid the cancellation errors in computing x*x-y*y by using the > formula (x-y)*(x+y) instead? Yes, but there seems to be no similar trick for general complex multiplication or higher powers of z. > Here is my proposed argument: Case 1: x and y have the same sign. If they > are nearly equal, then the computation (x-y) should be exact. The error in > computing (x+y) shouldn't be worse than 0.5 ULP. So the whole calculation > could probably be performed with an error less than 1 ULP, maybe even 0.5 ULP > if one took a lot of care. Case 2: x and y have different signs - just > switch the roles of (x+y) and (x-y) in the preceding argument. 1.5 or even 2.0 ulps I think (since there are 3 operations). > Do you know for a fact that gcc uses x*x-y*y? To use (x-y)*(x+y) would be no > worse performance wise - perhaps better if addition is cheaper than > multiplication. It would be the simplest programming switch in the world, as > long as one knew where to find the code fragment that needed to be changed. I checked, and it doesn't seem to. Bruce
Want to link to this message? Use this URL: <https://mail-archive.FreeBSD.org/cgi/mid.cgi?20120801180703.B1153>