computing expontential as repeated multiplication for small integer exponent
Subject: computing expontential as repeated multiplication for small integer exponent
From: Daniel Lakeland
Date: Mon, 30 Jul 2007 09:37:23 -0700
On Mon, Jul 30, 2007 at 09:35:41AM -0600, Robert Dodier wrote:
> On 7/30/07, Raymond Toy <raymond.toy at ericsson.com> wrote:
>
> > Your example of 1.2345678901234567d0^3 just illustrates that expt and
> > multiplication are different. As best as I can tell, (expt x 3) from
> > cmucl (and sbcl?) is more accurate than repeated multiplication.
>
> At this point, I think it is more important to get various Lisp
> implementations to agree, rather than getting maximally-accurate
> results for some implementations. I say this after having crunched
> more than a few numbers in various contexts.
I think maximally accurate results are more important than consistency
between implementations *FOR THE USER*. The usual user will be using a
single implementation, and will want accurate computations.
> Well, this fuzz factor business is problematic. There is already
> a fuzz factor; it is 8 * (64 bit float epsilon). Test #69 comes out
> with a greater error because it is the result of an iterative algorithm
> (namely find_root); a small discrepancy in EXPTB causes find_root
> to take a different number of steps for different Lisps.
> In general it is difficult to gauge an acceptable discrepancy.
> For one of the test cases for lsquares which I was recently working
> on, different Lisps agree to only 3 decimal places.
> When floating point operations aren't guaranteed to come
> out the same, we can't find a suitable fuzz factor without
> considering each case, which is tedious.
This is undoubtedly due to the non-stability of naive least
squares. a^2+b^2+c^2... quickly gets large and then individual terms
are small compared to the accumulation, leading to poor roundoff
behavior. Due to the nature of a sum of squares, the minimum will
always be at the bottom of a locally paraboloid surface. Finding that
minimum will be ill-conditioned.
I think that floating point test cases should specify an interval into
which the floating point result should fall, and this interval should
be calculated from analytical considerations when at all
possible. After all, we have a CAS to help us get exact results, and
also to help us compute error terms. So for example, x^3 should fall
in an interval of width 3 * x^2 * float_eps * 2 or a little wider.
--
Daniel Lakeland
dlakelan at street-artists.org
http://www.street-artists.org/~dlakelan