# Associativity in floating point arithmetic failing by two values

Cross-posting from math.stackexchange, since there might be people here familiar with this topic.

Assume working in floating point arithmetic with finite precision, bounded exponent and rounding to nearest.

Let $x,y$ be positive. It is not hard to find examples such that

$$(x+y)-x > y$$

In order to construct a test for a software I wanted to find an example or prove that such example doesn't exist of following:

Let $s(x)$ be the successor of $x$.

Is it possible to have both \begin{align}(x+y)-x&>y\\(x+y)-s(x)&>y\end{align}

It is easy to write a program that searches for such an example, but also it is unfeasible to test all possibilities and show that the example doesn't exist in this manner. So far my code hasn't got any example.

Example: In case seeing an example of $(x+y)-x>y$ helps somehow, take \begin{align} x&=1.1234567891234568\\ y&=1e-5\text{ ( denoting }10^{-5}) \end{align} Then $(x+y)-x=1.0000000000065512e-05 > y$. Examples of the first inequality there are many.

These kinds of floating point counterexamples are notoriously tricky to discover by hand. Usually you do expect random brute-force search to work, especially on 32-bit floats, but sometimes it fails.

These days, you have SMT solvers, like Z3 available to you, so you can try formulating your formula in Z3 to see if it can find a counterexample. This is the "high-tech" method of solving your problem. This approach is nice because it doesn't require you to know very much about floating point arithmetic to be able to use it.

I tried this, using Haskell, with the SBV package using Z3 as the backend. It looks only at the range $x\in(1,2)$, where $s(x) = x + \epsilon$

import Data.SBV

ε = 1.1920928955078125e-07

f :: SFloat -> SFloat -> SBool
f x y = 1 .< x &&& x .< 2 &&& fpIsNormal y &&& (x + y) - (x + ε) .> y


and it immediately found this counterexample:

λ> sat f
Satisfiable. Model:
s0 =       1.5 :: Float
s1 = 3.9999998 :: Float


This satisfies $(x+y)-s(x) > y$:

julia> x, y = 1.5f0, 3.9999998f0
julia> ((x + y) - nextfloat(x)) - y
2.3841858f-7
julia> (x + y) - nextfloat(x) > y
true


Note Looking at your example, it's unclear to me if you have restrictions on the magnitudes on $x$ and $y$. For example, if you restrict it to $x\in(1,2)$, $y\in(0,1)$, then there is no way to satisfy $(x+y)-s(x)>y$. But restricting only $x\in(1,2)$, but not $y$, it seems you get very large $n$ with $(x+y)-(x+n\epsilon) > y$.

Looking at some small examples, it seems the worst case can be reached by taking (in binary, single precision) $$x = 1.5 = (1.1)_2, \qquad y = 2^k\times (0.111111111111111111111111)_2 = 2^k\times (1-2^{-24}),$$ which gives $n=2^{k-1}$, so $(x+y)-(x+2^{k-1}\epsilon)=y$. In double precision it works to take $x=1.5$ and $y=(1-2^{-53})\times 2^k$ (which is just all ones in binary), and the relationship only starts to break down from $k\geq51$.

• These floats are single precision? – EEE Oct 3 '17 at 23:15
• @EEE Yes, that makes it less work for the SMT solver. But if the inequality doesn't hold for single precision floats, it will very likely also not hold for double precision floats—these things only very rarely depend on the number of bits of precision. – Kirill Oct 3 '17 at 23:17
• I am not familiar with Haskel, SBV or Z3, but is the nextafter float being computed by adding $\epsilon$? If so, this is strange to me. As I understand the intervals between floating points are not constant, are they? Oh, maybe it is constant in the interval $(1,2)$. Nevermind, it is constant between consecutive powers of $2$. – EEE Oct 3 '17 at 23:20
• @EEE Yes, for a specific exponent $e$, $x=s\times 2^e$, the distance between consecutive numbers is constant, I was relying on that to make the code more simple. – Kirill Oct 3 '17 at 23:23
• By the way, the claim that there are no solutions for $x\in(1,2)$ and $y\in(0,1)$, was it obtained by Z3 or you made it through an argument by hand? – EEE Oct 4 '17 at 1:12