help-octave
[Top][All Lists]

## Re: precision of calculations

 From: Ted Harding Subject: Re: precision of calculations Date: Fri, 23 Jul 1999 16:54:28 +0100 (BST)

```On 23-Jul-99 A. Scottedward Hodel wrote:
> double precision arithmetic has about 16 significant digits
> (depending on hardware).  Hence, if (A-B) should be zero, but
> you get
>     norm(A-B)/norm(A) approximately = 1e-16
> you're o.k.  The relative error can creep up as problem dimension
> and problem condition (how sensitive the results are to small
> perturbations in data) increases.

Good reply. In those cases where these tiny discrepancies really matter,
I'd suggest doing the following kind of test at critical points (example
code only -- there'll be better ways of doing it in real life):

....

x = A-B;
if( abs(x) < 5*eps ), x=0; endif
....

In practice, I reckon that this will make the right correction more often
than the wrong one: in real life, (A-B) is much more likely to be
precisely zero, than to be that small but non-zero.

By the way, on my octave eps =  2.2204e-16.

I suggest a factor like 5 (but definitively) because smallish
factors arise very easily indeed; e.g.

2 - sqrt(2)*sqrt(2) = -4.4409e-16   [ = -2*eps ]

This will heavily complicate what might normally be straightforward code,
and may slow it down a lot; but when it matters, it matters!

Ted.

PS Interestingly,

2 - (2.^(1/3))*(2.^(1/3))*(2.^(1/3)) = 0

while

2 - (2.^(1/4))*(2.^(1/4))*(2.^(1/4))*(2.^(1/4)) = -2*eps again!

You never know, do you?

--------------------------------------------------------------------
Date: 23-Jul-99                                       Time: 16:54:28
------------------------------ XFMail ------------------------------

---------------------------------------------------------------------
Octave is freely available under the terms of the GNU GPL.  To ensure