How to determine that 2.0 == (float) 2 in a front-end.

[Geoff Keating]

kaih@khms.westfalen.de (Kai Henningsen) writes:

> dewar@gnat.com (Robert Dewar)  wrote on 12.07.02 in <20020712163252.D2579F2CB6@nile.gnat.com>:
> 
> > <<Why is that?  Naively, I would think that x*x*x*x (evaluated left to
> > right) has 3 multiplications in a row, while (x**2)**2 has only two.  So,
> > from an error propagation point of view, the latter should be more accurate.
> > >>
> >
> > It is always hard to get an intuitive feel for these results (the above
> > result is fairly easy to prove). Here is an attempt. In the case of
> > x*x*x*x, you have 3 multiplications but in the 2nd and 3rd, one of the
> > operands at least is exact. In the case of (x**2)**2 the second
> > multiplication has two potentially inexact operands so the error is
> > multiplied.
> 
> Let's see. (And let's assume x>0 so we can ignore the sign, which makes no  
> difference to the error anyway.)
> 
> With eps being the largest possible error for a multiplication:
> 
> x*x     = x^2+/-eps
> x*x*x   = x^3+/-eps*x+/-eps = x^3+/-eps*(x+1)
> x*x*x*x = x^4+/-eps*(x+1)*x+/-eps = x^4+/-eps*(x^2+x+1)
> 
> (x*x)*(x*x) = x^4+/-2*eps+eps^2+/-eps = x^4+/-eps*(3+eps)
> 
> Assuming I did not make a mistake here ... So, what is larger - x^2+x+1 or  
> eps+3? It seems there is no one true answer - for small x, x*x*x*x is  
> better, whereas for large x, (x*x)*(x*x) is better. |x|=1 seems to be  
> close to the dividing point.

I think I'll have to object to Robert's claim that this is easy.

The obvious way to do it is to model the multiplication as

x*y -> x*y*(1 +/- eps)

but then of course you discover that it doesn't matter in what order
you multiply, the error is related only to the (logical) number of
multiplications, and there are three multiplications both ways.

Here's a better analysis, in cases.  I'll write *_ to mean the
IEEE multiplication.

Let 'eps' be such that all values are of the form (1 + n*eps) * 2^m,
for integers m and n.  Let x = (1 + n_x*eps) * 2^(m_x).

We can assume that m_x=0.

The maximum error in a single multiplication is (0.5*eps)*2^m (where
'm' is from the result).

Let a = x *_ x.
Let b = a *_ a.
Let c = a *_ x.
Let d = c *_ x.

We can ignore the error introduced by the final multiplication, since
it'll be the same for b and d; it won't be mentioned at all below.

Now, there are these cases:

1. n_x == 0.  In this case the result will be completely accurate,
   it will be 1 + 2^(4*m_x).

2. 0 < 1 + n_x * eps < 2^(1/3).  In this case, the the maximum errors will 
   be:

   a = x^2 +/- 0.5*eps
   b = x^4 +/- 0.5*eps*(2*x^2+0.5*eps)
   c = x^3 +/- 0.5*eps*(x+1)
   d = x^4 +/- 0.5*eps*(x^2+x)

   Since we know x^2 + eps > x, 'd' is more accurate, by up to
   0.5*eps*0.327.

3. 2^(1/3) <= 1 + n_x * eps < 2^(1/2)

   a = x^2 +/- 0.5*eps
   b = x^4 +/- 0.5*eps*(2*x^2+0.5*eps)
   c = x^3 +/- 0.5*eps*(x+2)
   d = x^4 +/- 0.5*eps*(x^2+2*x)

   Here, 'b' is more accurate because x^2 < 2*x, by up to 0.5*eps*0.964.

4. 2^(1/2) <= 1 + n_x * eps < 2^(2/3)

   a = x^2 +/- 0.5*eps*2
   b = x^4 +/- 0.5*eps*2*(2*x^2+0.5*eps*2)
   c = x^3 +/- 0.5*eps*2*(x+1)
   d = x^4 +/- 0.5*eps*2*(x^2+x)

   Here, 'd' is more accurate again, by up to 0.5*eps*2*0.932.

5. 2^(2/3) <= 1 + n_x * eps < 2

   a = x^2 +/- 0.5*eps*2
   b = x^4 +/- 0.5*eps*2*(2*x^2+0.5*eps*2)
   c = x^3 +/- 0.5*eps*2*(x+2)
   d = x^4 +/- 0.5*eps*2*(x^2+2*x)

   Here, 'b' is more accurate again, by up to 0.5*eps*2*0.655.

My conclusion from all this is that there's no clear winner.
Sometimes, the last multiplication for 'd' is better because one of
its arguments is more accurate; sometimes, 'b' is better because it
avoids repeated rounding.

And the moral of the story is, you can keep analysing as long as you
like and keep getting different results :-).

-- 
- Geoffrey Keating <geoffk@geoffk.org> <geoffk@redhat.com>