Short-circuit Logic

[Chris Angelico]

py> y = 1e17 + x # x is not zero, so y should be > 1e17
py> 1/(1e17 - y)
Traceback (most recent call last):
File "<stdin>", line 1, in <module>
ZeroDivisionError: float division by zero

You don't even need to go for 1e17. By definition:
True

Therefore the same can be done with 2 as you did with 1e17.
Traceback (most recent call last):
File "<pyshell#182>", line 1, in <module>
1/(2-y)
ZeroDivisionError: float division by zero

Of course, since we're working with a number greater than epsilon, we
need to go a little further, but we can still work with small numbers:
Traceback (most recent call last):
File "<pyshell#191>", line 1, in <module>
1/(4-y)
ZeroDivisionError: float division by zero


ChrisA

 

[Grant Edwards]

Can you show me a value of x where x == 0.0 returns False, but x actually
isn't zero?

I'm confused. Don't all non-zero values satisfy your conditions?
False

 

[Steven D'Aprano]

I'm confused. Don't all non-zero values satisfy your conditions?

Of course they do :-(

I meant "but x actually *is* zero". Sorry for the confusion. I blame the
terrists.

 

[Ahmed Abdulshafy]

Incorrect. Testing for equality is always precise, and exact. The problem

is not the *equality test*, but that you don't always have the number

that you think you have. The problem lies elsewhere, not equality!


Steven

Well, this is taken from my python shell>
True

Anyway, man, those were not my words anyway, most programming books I've read state so. Here's an excerpt from the Python book, I'm currently reading>

">>> 0.0, 5.4, -2.5, 8.9e-4
(0.0, 5.4000000000000004, -2.5, 0.00088999999999999995)


The inexactness is not a problem specific to Python—all programming languages have this problem with floating-point numbers."

 

[Chris Angelico]

Well, this is taken from my python shell>

(0.3345585735242628, 0.3345585735242628)

They're not representably different.

ChrisA

 

[rusi]

Well, this is taken from my python shell>


True

Anyway, man, those were not my words anyway, most programming books I've read state so. Here's an excerpt from the Python book, I'm currently reading>

">>> 0.0, 5.4, -2.5, 8.9e-4
(0.0, 5.4000000000000004, -2.5, 0.00088999999999999995)

The inexactness is not a problem specific to Python—all programming languages have this problem with floating-point numbers."

0.0 == 0.0 implies 5.4 == 5.4
is not a true statement is what (I think) Steven is saying.
0 (or if you prefer 0.0) is special and is treated specially.

Naturally if you reach (nearabout) 0.0 by some numerical process thats
another matter...

 

[Dave Angel]

It has nothing to do with 0 being special. A floating point number
will always equal itself (except for nan, which is even more special),
and in particular 5.4 == 5.4. But if you have two different
calculations that produce 0, or two different calculations that
produce 5.4, you might actually get two different numbers that
approximate 0 or 5.4 thanks to rounding error. If you then compare
those two ever-so-slightly different numbers, you will find them
unequal.

Rounding error is just one of the problems. Usually less obvious is
quantization error. If you represent a floating number in decimal, but
you're using a binary floating point representation, it just might change.

Another error is roundoff error. Even in a pure decimal system of (say)
40 digits, I could type in a 42 digit number and it would get quantized.
So just because two 42 digit numbers are different doesn't imply that
the 40 digit internal format would be.

 

[Steven D'Aprano]

It has nothing to do with 0 being special. A floating point number will
always equal itself (except for nan, which is even more special), and in
particular 5.4 == 5.4. But if you have two different calculations that
produce 0, or two different calculations that produce 5.4, you might
actually get two different numbers that approximate 0 or 5.4 thanks to
rounding error. If you then compare those two ever-so-slightly
different numbers, you will find them unequal.

EXACTLY!

The problem does not lie with the *equality operator*, it lies with the
calculations. And that is an intractable problem -- in general, floating
point is *hard*. So the problem occurs when we start with a perfectly
good statement of the facts:

"If you naively test the results of a calculation for equality without
understanding what you are doing, you will often get surprising results"

which then turns into a general heuristic that is often, but not always,
reasonable:

"In general, you should test for floating point *approximate* equality,
in some appropriate sense, rather than exact equality"

which then gets mangled to:

"Never test floating point numbers for equality"

and then implemented badly by people who have no clue what they are doing
and have misunderstood the nature of the problem, leading to either:

* de facto exact equality testing, only slower and with the *illusion* of
avoiding equality, e.g. "abs(x-y) < sys.float_info.epsilon" is just a
long and slow way of saying "x == y" when both numbers are sufficiently
large;

* incorrectly accepting non-equal numbers as "equal" just because they
happen to be "close".


The problem is that there is *no one right answer*, except "have everyone
become an expert in floating point, then judge every case on its merits",
which will never happen.

But if nothing else, I wish that we can get past the rank superstition
that you should "never" test floats for equality. That would be a step
forward.

 

[Chris Angelico]

* de facto exact equality testing, only slower and with the *illusion* of
avoiding equality, e.g. "abs(x-y) < sys.float_info.epsilon" is just a
long and slow way of saying "x == y" when both numbers are sufficiently
large;

The problem here, I think, is that "epsilon" has two meanings:

* sys.float_info.epsilon, which is an extremely specific value (the
smallest x such that 1.0+x != x)

* the mathematical concept, which is where the other got its name from.

Let's suppose someone is told to compare floating point numbers by
seeing if the absolute value of the difference is less than some
epsilon. They look up "absolute value" and find abs(); they look up
"epsilon" and think they've found it. Trouble is, they've found the
wrong epsilon... and really, there's an engineering issue here too.
Here's one of my favourite examples of equality comparisons:

http://xkcd.com/1047/

# Let's say we measured this accurately to one part in 40
x = one_light_year_in_meters

y = pow(99,8)
x == y # False
abs(x-y) < x/40 # True

Measurement accuracy is usually far FAR worse than floating-point
accuracy. It's pretty pointless to compare for some kind of "equality"
that ignores this. Say you measure the diameter and circumference of a
circle, accurate to one meter, and got values of 79 and 248; does this
mean that pi is less than 3.14? No - in fact:

pi = 248/79
# math.pi = 3.141592653589793
abs(pi-math.pi) < pi/79 # True

Worst error is 1 in 79, so all comparisons are done with epsilon
derived from that.

ChrisA

 

[Steven D'Aprano]

Well, this is taken from my python shell>

True

This is an excellent example of misunderstanding what you are seeing.
Both 0.33455857352426283 and 0.33455857352426282 represent the same
float, so it is hardly a surprise that they compare equal -- they compare
equal because they are equal.

py> a, b = 0.33455857352426283, 0.33455857352426282
py> a.as_integer_ratio()
(6026871468229899, 18014398509481984)
py> b.as_integer_ratio()
(6026871468229899, 18014398509481984)

You've made a common error: neglecting to take into account the finite
precision of floats. Floats are not mathematical "real numbers", with
infinite precision. The error is more obvious if we exaggerate it:

py> 0.3 == 0.300000000000000000000000000000000000000000000000000001
True

Most people who have seen an ordinary four-function calculator will
realise that the issue here is *not* that the equality operator == is
wrongly stating that two unequal numbers are equal, but that just because
you enter 0.300...00001 doesn't mean that all those decimal places are
actually used.

Anyway, man, those were not my words anyway, most programming books I've
read state so. Here's an excerpt from the Python book, I'm currently
reading>

">>> 0.0, 5.4, -2.5, 8.9e-4
(0.0, 5.4000000000000004, -2.5, 0.00088999999999999995)


The inexactness is not a problem specific to Python—all programming
languages have this problem with floating-point numbers."

I'm not denying that floats are tricky to use correctly, or that testing
for exact equality is *sometimes* the wrong thing to do:

# Wrong, don't do this!
x = 0.1
while x != 17.3:
print(x)
x += 0.1


I'm just saying that a simple minded comparison with
sys.float_info.epsilon is *also* often wrong.

 

[Steven D'Aprano]

Even in a pure decimal system of (say)
40 digits, I could type in a 42 digit number and it would get quantized.
So just because two 42 digit numbers are different doesn't imply that
the 40 digit internal format would be.

Correct, and we can demonstrate it using Python:

py> from decimal import *
py> getcontext().prec = 3
py> a = +Decimal('1.0000')
py> b = +Decimal('1.0009')
py> a == b
True


(By default, the Decimal constructor does not honour the current
precision. To force it to do so, use the unary + operator.)

 

[Chris Angelico]

# Wrong, don't do this!
x = 0.1
while x != 17.3:
print(x)
x += 0.1

Actually, I wouldn't do that with integers either. There are too many
ways that a subsequent edit could get it wrong and go infinite, so I'd
*always* use an inequality for that:

x = 1
while x < 173:
print(x)
x += 1

Well, in Python I'd use for/range, but the equivalent still applies. A
range() is still based on an inequality:

list(range(1,6)) [1, 2, 3, 4, 5]
list(range(1,6,3))

[1, 4]

Stops once it's no longer less than the end. That's safe, since Python
can't do integer wraparound.

ChrisA

 

[Steven D'Aprano]

Let's suppose someone is told to compare floating point numbers by
seeing if the absolute value of the difference is less than some
epsilon.

Which is usually the wrong way to do it! Normally one would prefer
*relative* error, not absolute:

# absolute error:
abs(a - b) < epsilon


# relative error:
abs(a - b)/a < epsilon


One problem with absolute error is that it can give an entirely spurious
image of "fuzziness", when in reality it is actually performing the same
exact equality as == only slower and more verbosely. If a and b are
sufficiently large, the smallest possible difference between a and b may
be greater than epsilon (for whichever epsilon you pick). When that
happens, you might as well just use == and be done with it.

But using relative error also raises questions:

- what if a is negative?

- why relative to a instead of relative to b?

- what if a is zero?

The first, at least, is easy to solve: take the absolute value of a. But
strangely, you rarely see programming books mention that, so I expect
that there is a lot of code in the real world that assumes a is positive
and does the wrong thing when it isn't.

Here's another way, mathematically equivalent (although not necessarily
equivalent using floating point computations!) which avoids the divide-by-
zero problem:

abs(a - b) < epsilon*a


Whichever method you choose, there are gotchas to watch out for.

http://xkcd.com/1047/

Nice!

 

[Jussi Piitulainen]

Which is usually the wrong way to do it! Normally one would prefer
*relative* error, not absolute:

# absolute error:
abs(a - b) < epsilon


# relative error:
abs(a - b)/a < epsilon

....

I wonder why floating-point errors are not routinely discussed in
terms of ulps (units in last position). There is a recipe for
calculating the difference of two floating point numbers in ulps, and
it's possible to find the previous or next floating point number, but
I don't know of any programming language having built-in support for
these.

Why isn't this considered the most natural measure of a floating point
result being close to a given value? The meaning is roughly this: how
many floating point numbers there are between these two.

"close enough" if abs(ulps(a, b)) < 3 else "not close enough"

"equal" if ulps(a, b) == 0 else "not equal"

There must be some subtle technical issues here, too, but it puzzles
me that this measure of closeness is not often even discussed when
absolute and relative error are discussed - and computed using the
same approximate arithmetic whose accuracy is being measured. Scary.

Got light?

 

[Chris Angelico]

Which is usually the wrong way to do it! Normally one would prefer
*relative* error, not absolute:

# absolute error:
abs(a - b) < epsilon


# relative error:
abs(a - b)/a < epsilon

I was picking an epsilon based on a, though, which comes to pretty
much the same thing as the relative error calculation you're using.

But using relative error also raises questions:

- what if a is negative?

- why relative to a instead of relative to b?

- what if a is zero?

The first, at least, is easy to solve: take the absolute value of a.

One technique I saw somewhere is to use the average of a and b. But
probably better is to take the lower absolute value (ie the larger
epsilon). However, there's still the question of what epsilon should
be - what percentage of a or b you take to mean equal - and that one
is best answered by looking at the original inputs.

Take these guys, for instance. Doing the same thing I was, only with
more accuracy.

ChrisA

 

[Steven D'Aprano]

I wonder why floating-point errors are not routinely discussed in terms
of ulps (units in last position). There is a recipe for calculating the
difference of two floating point numbers in ulps, and it's possible to
find the previous or next floating point number, but I don't know of any
programming language having built-in support for these.

That is an excellent question!

I think it is because the traditional recipes for "close enough" equality
either pre-date any standardization of floating point types, or because
they're written by people who are thinking about abstract floating point
numbers and not considering the implementation.

Prior to most compiler and hardware manufacturers standardizing on IEEE
754, there was no real way to treat float's implementation in a machine
independent way. Every machine laid their floats out differently, or used
different number of bits. Some even used decimal, and in the case of a
couple of Russian machines, trinary. (Although that's going a fair way
back.)

But we now have IEEE 754, and C has conquered the universe, so it's
reasonable for programming languages to offer an interface for accessing
floating point objects in terms of ULPs. Especially for a language like
Python, which only has a single float type.

I have a module that works with ULPs. I may clean it up and publish it.
Would there be interest in seeing it in the standard library?

Why isn't this considered the most natural measure of a floating point
result being close to a given value? The meaning is roughly this: how
many floating point numbers there are between these two.

There are some subtleties here also. Firstly, how many ULP should you
care about? Three, as you suggest below, is awfully small, and chances
are most practical, real-world calculations could not justify 3 ULP.
Numbers that we normally care about, like "0.01mm", probably can justify
thousands of ULP when it comes to C-doubles, which Python floats are.

Another subtlety: small-but-positive numbers are millions of ULP away
from small-but-negative numbers. Also, there are issues to do with +0.0
and -0.0, NANs and the INFs.

 

[Jussi Piitulainen]

....

But we now have IEEE 754, and C has conquered the universe, so it's
reasonable for programming languages to offer an interface for
accessing floating point objects in terms of ULPs. Especially for a
language like Python, which only has a single float type.

Yes, that's what I'm thinking, that there is now a ubiquitous floating
point format or two, so the properties of the format could be used.

I have a module that works with ULPs. I may clean it up and publish it.
Would there be interest in seeing it in the standard library?

Yes, please.

There are some subtleties here also. Firstly, how many ULP should
you care about? Three, as you suggest below, is awfully small, and
chances are most practical, real-world calculations could not
justify 3 ULP. Numbers that we normally care about, like "0.01mm",
probably can justify thousands of ULP when it comes to C-doubles,
which Python floats are.

I suppose this depends on the complexity of the process and the amount
of data that produced the numbers of interest. Many individual
floating point operations are required to be within an ulp or two of
the mathematically correct result, I think, and the rounding error
when parsing a written representation of a number should be similar.
Either these add up to produce large errors, or the computation is
approximate in other ways in addition to using floating point.

One could develop a kind of sense for such differences. Ulps could be
a tangible measure when comparing different algorithms. (That's what I
tried to do with them in the first place. And that's how I began to
notice their absence when floating point errors are discussed.)

Another subtlety: small-but-positive numbers are millions of ULP
away from small-but-negative numbers. Also, there are issues to do
with +0.0 and -0.0, NANs and the INFs.

The usual suspects ^_^ and no reason to dismiss the ulp when the
competing kinds of error have their corresponding subtleties. A matter
of education, I'd say.

Thank you much for an illuminating discussion.

 

[Roy Smith]

Actually, I wouldn't do that with integers either. There are too many
ways that a subsequent edit could get it wrong and go infinite, so I'd
*always* use an inequality for that:

x = 1
while x < 173:
print(x)
x += 1

There's a big difference between these two. In the first case, using
less-than instead of testing for equality, you are protecting against
known and expected floating point behavior.

In the second case, you're protecting against some vague, unknown,
speculative future programming botch. So, what *is* the right behavior
if somebody were to accidentally drop three zeros into the source code:

x = 1000
while x < 173:
print(x)
x += 1

should the loop just quietly not execute (which is what it will do
here)? Will that make your program correct again, or will it simply
turn this into a difficult to find bug? If you're really worried about
that, why not:

 

[Roy Smith]

I wonder why floating-point errors are not routinely discussed in
terms of ulps (units in last position).

Analysis of error is a complicated topic (and is much older than digital
computers). These sorts of things come up in the real world, too. For
example, let's say I have two stakes driven into the ground 1000 feet
apart. One of them is near me and is my measurement datum.

I want to drive a third stake which is 1001 feet away from the datum.
Do I measure 1 foot from the second stake, or do I take out my
super-long tape measure and measure 1001 feet from the datum?

 

[Chris Angelico]

if somebody were to accidentally drop three zeros into the source code:


should the loop just quietly not execute (which is what it will do
here)? Will that make your program correct again, or will it simply
turn this into a difficult to find bug? If you're really worried about
that, why not:

If you iterate from 1000 to 173, you get nowhere. This is the expected
behaviour; this is what a C-style for loop would be written as, it's
what range() does, it's the normal thing. Going from a particular
starting point to a particular ending point that's earlier than the
start results in no iterations. The alternative would be an infinite
number of iterations, which is far far worse.

ChrisA