From bug-octave-request at bevo dot che dot wisc dot edu  Tue Jul 11 10:08:29 2000
Subject: Re: Octave produces erroneous values when adding and subtracting
From: "Ross A. Lippert" <ripper at cs dot sandia dot gov>
To: Thomas Walter <walter at pctc dot chemie dot uni-erlangen dot de>
CC: etienne at isr dot isr dot ist dot utl dot pt, bug-octave@bevo.che.wisc.edu,         etienne@isr.ist.utl.pt
Date: Tue, 11 Jul 2000 09:06:08 -0600

There is an anecdote about the missile defense system used in a recent US war which
needed to be rebooted bc it had been programmed to count time in intervals on 1 tenth
of a second.  The problem is, of course that
.1 + .1 + .... + .1 != 1.0 be .1 cannot represented exactly in IEEE floating point.
The result was clock skew which eventually caused the system to malfunction and cost
lives (I credit Nick Higham on U of Manchester with this story).

Checking octave expression evaluation with C/C++ compilers is bound to produce
differences in evaluation.  All other things equal, octave will left associate
operations of equal precedence (check (1+-1+1e-16) for example).  All other
things being equal C/C++ compilers give no guarrantee that the compiler associate
in a given way.  In fact, C/C++ compilers are free to ignore parentheses when
reassociating, possibly compiling ((1+-1)+1e-16) to (1+(-1+1e-16)) (I credit the
O'reilly high performance computing book with this tidbit).

I hope that no one's life is relying upon the correct evaluation of expressions
like A-(1+A-1).  Whatever it is your are trying to compute via this blackmagic should
be done better in your code.


-r



Thomas Walter wrote:
> 
> >>>>> "etienne" == etienne grossmann <etienne at anonimo dot isr dot ist dot utl dot pt> writes:
> 
> [snip]
> 
>     >>> A=1.2; A - (1 + A - 1)
>     etienne> ans =  -2.2204e-16
> 
>     etienne> #  I checked it with a simple C and a simple C++ program because I also
>     etienne> #  depend on this while doing fitting and finding local minima.  And
>     etienne> #  these testprograms are correct as long as A > 1e-9.
> 
>     etienne> #  Therefore I think it is a bug in octave.  I don't know where but plain
>     etienne> #  C/C++ has no problems.
> 
>     etienne>   Maybe your C compiler is "smart" and simplifies constants when it
>     etienne> finds them (i.e. it sees that 1+A-1 == A and your program never
>     etienne> actually adds and substracts 1). That's just a guess.
> 
> Hello,
> your guess is correct!
> I jumped into the common pit with too simple tests for the compiler. 8-(
> And it is another proof to request for
> if (abs(A-B) < (100 * eps))
>    printf("A and B are equal");
> endif
> 
> --
> Was gibt sieben mal sieben?  Ganz feinen Sand. 8-)
> 
> ----------------------------------------------
> Dipl. Phys. Thomas Walter
> Inst. f. Physiklische Chemie II
> Egerlandstr. 3                          Tel.: ++9131-85 27326 / 27330
> 91058 Erlangen, Germany                 email: walter at pctc dot chemie dot uni-erlangen dot de
> 
> -----------------------------------------------------------------------
> Octave is freely available under the terms of the GNU GPL.
> 
> Octave's home on the web:  http://www.che.wisc.edu/octave/octave.html
> How to fund new projects:  http://www.che.wisc.edu/octave/funding.html
> Subscription information:  http://www.che.wisc.edu/octave/archive.html
> -----------------------------------------------------------------------



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

Octave's home on the web:  http://www.che.wisc.edu/octave/octave.html
How to fund new projects:  http://www.che.wisc.edu/octave/funding.html
Subscription information:  http://www.che.wisc.edu/octave/archive.html
-----------------------------------------------------------------------