From: David Cournapeau on
On Sun, Dec 6, 2009 at 1:46 AM, Mark Dickinson <dickinsm(a)gmail.com> wrote:
> On Dec 5, 3:37 pm, Anton81 <gerenu...(a)googlemail.com> wrote:
>> I'd like to do calculations with floats and at some point equality of
>> two number will be checked.
>> What is the best way to make sure that equality of floats will be
>> detected, where I assume that mismatches beyond a certain point are
>> due to truncation errors?
>
> Well, it depends a lot on the details of the application, but
> a good general scheme is to allow both a fixed relative error
> and an absolute error, and to assert that your two values are
> 'nearly equal' if they're within *either* the relative error *or*
> the absolute error.  Something like, for example:
>
> def almostEqual(expected, actual, rel_err=1e-7, abs_err = 1e-20):
>    absolute_error = abs(actual-expected)
>    return absolute_error <= max(abs_err, rel_err * abs(expected))

If you can depend on IEEE 754 semantics, one relatively robust method
is to use the number of representable floats between two numbers. The
main advantage compared to the proposed methods is that it somewhat
automatically takes into account the amplitude of input numbers:

abs(x - y) <= N * spacing(max(abs(x), abs(y)))

Where spacing(a) is the smallest number such as a + spacing(a) != a.
Whether a and b are small or big, the same value of N can be used, and
it tells you how close two numbers are in terms of internal
representation.

Upcoming numpy 1.4.0 has an implementation for spacing - implementing
your own for double is not difficult, though,

cheers,

David
From: Dave Angel on


Carl Banks wrote:
> On Dec 6, 11:34 am, Anton81 <gerenu...(a)googlemail.com> wrote:
>
>> I do some linear algebra and whenever the prefactor of a vector turns
>> out to be zero, I want to remove it.
>>
>> I'd like to keep the system comfortable. So basically I should write a
>> new class for numbers that has it's own __eq__ operator?
>> Is there an existing module for that?
>>
>
> I highly recommend against it; among other things it invalidates the
> transitive property of equality:
>
> "If a =b and b == c, then a == c."
>
> It will also make the number non-hashable, and have several other
> negative consequences. Plus, it's not something that's never
> foolproof. What numbers are close enought to be condidered "equal"
> depends on the calculations.
>
> (I remember once struggling in a homework assignment over seemingly
> large discrepancies in a calculation I was doing, until i realized
> that the actual numbers were on the scale of 10**11, and the
> difference was around 10**1, so it really didn't matter.)
>
>
>
> Carl Banks
>
>
A few decades ago I implemented the math package (microcode) under the
machine language for a proprietary processor (this is when a processor
took 5 boards of circuitry to implement). I started with floating point
add and subtract, and continued all the way through the trig, log, and
even random functions. Anyway, a customer called asking whether a
particular problem he had was caused by his logic, or by errors in our
math. He was calculating the difference in height between an
always-level table and a perfectly flat table (between an arc of a great
circle around the earth, and a flat table that doesn't follow the
curvature.) In a couple of hundred feet of table, the difference was
measured in millionths of an inch, as I recall. Anyway it turned out
his calculation was effectively subtracting
(8000 miles plus a little bit) - (8000 miles)
and if he calculated it three different ways, he got three different
results, one was off in about the 3rd place, while the other was only
half the value. I was able to show him another way (through geometrical
transformations) to solve the problem that got the exact answer, or at
least to more digits than he could possibly measure. I think I recall
that the new solution also cancelled out the need for trig. Sometimes
the math package shouldn't hide the problem, but give it to you straight.

DaveA
From: dbd on
On Dec 6, 1:48 pm, sturlamolden <sturlamol...(a)yahoo.no> wrote:
> On 6 Des, 21:52, r0g <aioe....(a)technicalbloke.com> wrote:
>
> > > .> Right.  Using abs(x-y) < eps is the way to go.
> > > .>
> > > .> Raymond
>
> > > This only works when abs(x) and abs(y) are larger that eps, but not
> > > too much larger.
>
> > Okay, I'm confused now... I thought them being larger was entirely the
> > point.
>
> Yes. dbd got it wrong. If both a smaller than eps, the absolute
> difference is smaller than eps, so they are considered equal.

Small x,y failure case:
eps and even eps squared are representable as floats. If you have
samples of a sine wave with peak amplitude of one half eps, the "abs(x-
y) < eps" test would report all values on the sine wave as equal to
zero. This would not be correct.
Large x,y failure case:
If you have two calculation paths that symbolically should produce the
same value of size one over eps, valid floating point implementations
may differ by an lsb or more. An single lsb error would be 1, much
greater than the test allows as 'nearly equal' for floating point
comparison.

1.0 + eps is the smallest value greater than 1.0, distinguishable from
1.0. Long chains of floating point calculations that would
symbolically be expected to produce a value of 1.0 many be expected to
produce errors of an eps or more due to the inexactness of floating
point representation. These errors should be allowed in floating point
equality comparison. The value of the minimum representable error will
scale as the floating point number varies. A constant comparison value
is not appropriate.

Mark was right, DaveA's discussion explains a strategy to use.

Dale B. Dalrymple
From: Steven D'Aprano on
On Sun, 06 Dec 2009 14:54:37 -0800, Carl Banks wrote:

> (I remember once struggling in a homework assignment over seemingly
> large discrepancies in a calculation I was doing, until i realized that
> the actual numbers were on the scale of 10**11, and the difference was
> around 10**1, so it really didn't matter.)

Well that depends on the accuracy of the calculations, surely? If the
calculations were accurate to one part in 10**20, then an error around
10**1 is about ten trillion times larger than acceptable.

*wink*



--
Steven
From: sturlamolden on
On 7 Des, 06:43, dbd <d...(a)ieee.org> wrote:

> If you have
> samples of a sine wave with peak amplitude of one half eps, the "abs(x-
> y) < eps" test would report all values on the sine wave as equal to
> zero. This would not be correct.

You don't understand this at all do you?

If you have a sine wave with an amplitude less than the truncation
error, it will always be approximately equal to zero.

Numerical maths is about approximations, not symbolic equalities.


> 1.0 + eps is the smallest value greater than 1.0, distinguishable from
> 1.0.

Which is the reason 0.5*eps*sin(x) is never distinguishable from 0.


> A constant comparison value is not appropriate.

That require domain specific knowledge. Sometimes we look at a
significant number of digits; sometimes we look at a fixed number of
decimals; sometimes we look at abs(y/x). But there will always be a
truncation error of some sort, and differences less than that is never
significant.






>
> Mark was right, DaveA's discussion explains a strategy to use.
>
> Dale B. Dalrymple