Floating point subtraction rounding error (NOT display error)

[Steven D'Aprano]

On Thu, 13 Dec 2007 14:30:18 -0800, Keflavich wrote:

> Hey, I have a bit of code that died on a domain error when doing an
> arcsin, and apparently it's because floating point subtraction is having
> problems.

I'm not convinced that your diagnosis is correct. Unless you're using 
some weird, uncommon hardware, it's unlikely that a bug in the floating 
point subtraction routines has escaped detection. (Unlikely, but not 
impossible.) Can you tell us what values give you incorrect results?



> I know about the impossibility of storing floating point
> numbers precisely, but I was under the impression that the standard used
> for that last digit would prevent subtraction errors from compounding.

What gave you that impression? Are you referring to guard digits?


> Is there a simple solution to this problem, or do I need to run some
> sort of check at every subtraction to make sure that my float does not
> deviate?  I'm not sure I know even how to do that.

I still don't quite understand your problem. If you think your float is 
deviating from the correct value, that implies you know what the correct 
value should be. How do you know what the correct value should be?

I should also mention that of course your answer will deviate, due to the 
finite precision of floats.


> A sample of the failure:
> ipdb>1.0000000000000001 == 1
> True
> ipdb>R
> 0.69999999999999996
> ipdb>R==.7
> True
> ipdb>y2
> 3.2999999999999998
> ipdb>y2 == 3.3
> True
> ipdb>cirY-y2
> 0.70000000000000018

What's cirY? How do you know this is the incorrect value?

> ipdb>cirY-y2 == .7
> False

Obviously not. As you've already shown, the correct value is 
0.70000000000000018, not 0.69999999999999996 (the closest floating point 
value to 0.7).


> I was unable to find solutions when searching the web because all of the
> hits I got were discussing display issues, which I'm not concerned with.

Have you read this?
http://docs.sun.com/source/806-3568/ncg_goldberg.html




-- 
Steven.