# How to compute rounding error bound?

I am using the IEEE754 half-precision floating point format, which has 11 significand bits.

My input is drawn randomly with values between 1.0 and 2.0.

I would like to approximate the maximum rounding error that could happen in the following computation.

(a1 + a2 + a3 + a4) * 0.25

I would only like to use units of least precision to calculate this rounding error.

Here's what I tried.

Let a1+a2 introduce rounding error e1.

Let a1+a2+e1+a3 introduce rounding error e2.

Let a1+a2+e1+a3+e2+a4 introduce rounding error e3.

Let (a1+a2+e1+a3+e2+a4+e3)*0.25 introduce rounding error e4.

The whole computation is

(a1+a2+e1+a3+e2+a4+e3)*0.25+e4
= ((a1+a2+a3+a4) + (e1+e2+e3))*0.25 + e4
= (a1+a2+a3+a4)*0.25 + (e1+e2+e3)*0.25 + e4


The rounding error is (e1+e2+e3)*0.25 + e4.

Largest possible e1 is about 0.5*ulp(2.0+2.0).

Largest possible e2 is about 0.5*ulp(2.0+2.0+2.0)

Largest possible e3 is about 0.5*ulp(2.0+2.0+2.0+2.0)

Largest possible e4 is about 0.5*ulp((2.0+2.0+2.0+2.0)*0.25)

So the maximum rounding error possible is then plugging these numbers in the total rounding error.

(e1+e2+e3)*0.25 + e4
= (0.5*ulp(4) + 0.5*ulp(6) + 0.5*ulp(8)) * 0.25 + 0.5*ulp(2)


I have several questions.

1. Am I doing it the right way?
2. Is there a way to over approximate the rounding error only by using one ulp and the number of operations?
3. If this is not the recommended way, how should I go about approximating the maximum rounding error?

My application does not require extreme accuracy.

EDIT: forgot to multiply ulp by 0.5.

• Are you considering your inputs to be exact and only want the rounding error induced by the computation, or do the inputs already have rounding error? Is your multiply by $0.25$ done as a floating multiply as opposed to a divide by an integer $4$? – Ross Millikan Jun 16 '18 at 1:23
• @RossMillikan I'm considering my inputs to be exact and only want the rounding error induced by the computation. My multiply is done as a floating point multiply, but if you have more things that I should concern with a divide by 4 and it's not too much effort, I'd love to learn more how that makes things different. Thanks. – Russell Jun 16 '18 at 2:21

Your general approach makes sense. I am not sure how you are defining $\operatorname{ulp}(n)$. Presumably $\operatorname{ulp}$ is units of the last place and $n$ is the magnitude of the quantity but then there are breaks as $n$ crosses a power of $2$. I will assume $\operatorname{ulp}(n)$ is a unit of the least place when $n$ is a power of $2$ and the value is between $n$ and $2n$. When you add $a1+a2$ the total is strictly less than $4$, so the exponent part will be $2^1$ if you always regard the mantissa as a number between $1$ and $2$.
When you add $a1+a2$ you get an exact $12$ bit result with the exponent $2^0$ and a mantissa between $2$ and $4$. You then shift the exponent by $1$ to normalize it and lose the LSB, which is at most $\operatorname{ulp}(1)$, which is your rounding error. When you add in $a3$ you throw away the LSB at the start to shift the exponent of $a3$ to match $a1+a2$. That is an additional $\operatorname{ulp}(1)$. What happens in the addition depends on whether $a1+a2+a3 \gt 4$. If it is less the addition is exact and we could be out $2\operatorname{ulp}(1)$. If it is greater we lose another LSB and could be out an additional $\operatorname{ulp}(2)$. The second is the worst case, so we will consider it and say we are out $4\operatorname{ulp}(1)$ so far. Now we add in $a4$, throwing away the two LSB to normalize, which can be an error of $3\operatorname{ulp}(1)$. We can't carry in this scenario, so there is no further error. We now have an error of up to $7\operatorname{ulp}(1)$.
When we multiply by $0.25$ it is better to divide by integer $4$. There is then no error and our worst case rounding error is $\frac 74\operatorname{ulp}(1)$ If you multiply by the float $0.25$ you need to worry whether you got $1\cdot 2^{-2}$ or $1.1111111111_2\cdot 2^{-3}$ The second will introduce a bit of rounding error.
You can do slightly better. It would be better to do $(a1+a2)+(a3+a4)$ and then divide by integer $4$. I assume they taught you in grade school that addition is associative, but computer addition is not. Here the $a3+a4$ addition happens with one smaller exponent, so the roundoff is reduced.
All of this is too much detail for a low precision application and depends on the fact that you know the magnitude of your addends is very close. For low precision applications the easy approach is to assume every operation has a relative error of $\epsilon$ and then to assume it is a random walk, so $k$ computations give a relative error of $\epsilon \sqrt k$. If that gets too large, use higher precision and call it good.
• It is like a random walk on a 1D line with fixed step size. The "average" distance from the origin after $n$steps is $\sqrt n$ – Ross Millikan Jun 16 '18 at 4:37