Evaluating $a(b + c)$ more accurately with FMA

Question

I'm using machine-precision floating-point arithmetic, and every so often it happens that I need to evaluate an expression of the form $a(b + c)$. I found that the accuracy can be improved using FMA (fused multiply-add), using one of these two forms:

1. $fma(a, b, a c)$
2. $fma(a, c, a b)$

Sometimes the first form is more accurate, and sometimes it's the second one.

Suppose $a$, $b$, $c$ belong to some floating-point format, I'd like conditions on $a$, $b$, $c$ to determine when the first form is the more accurate one, and when it's the second form.

Answer 1

This response doesn't specify conditions on $a$, $b$ or $c$ but does describe how you might go about computing the result you seek regardless of their values.

Let the exact sum of $b$ and $c$ be represented by $s$ and $t$, where $s$ and $t$ are floating-point numbers with $|t| <= ulp(s)/2$. (E.g., use the well-known 2Sum algorithm; see [1].) Thus $a(b + c)$ is exactly $a(s + t)$. Since $|s|$ is much larger than $|t|$, we can approximate the computation of $a(s + t)$ by fma(a, s, a * t).

Sample code might be

double f(double a, double b, double c) {
    double s = b + c;
    double t = (b - (s - c)) + (c - (s - (s - c)));
    return fma(a, s, a * t);
}

EDITED TO ADD:

Because the value of $a(s + t)$ is very close to $a \times s$ (recall: $|a * t| <= ulp(a * s)/2$), $a \times t$ is a relatively small correction term. Since fma(a, s, a * t) involves calculating $a \times s$ exactly without rounding and then adding $a * t$, it provides a very good approximation to $a(s + t)$.

If one were to calculate fma(a, t, a * s) instead, one would be calculating $a \times t$ exactly and then adding a * s (the rounded value of $a \times s$). Because of the relative magnitudes of $a * s$ and $a \times t$, this is essentially the same as calculating (a * s) + (a * t). In this case, no use is made of the additional precision provided by the fma operation.

I think the key is to visualize the alignment shift that takes place when, in fma(x, y, z), the addition of $z$ with the exact product of $x$ and $y$ is performed. In the case at hand, we want $x \times y$ to be the larger quantity so that $z$ is "shifted to the right" when the addition is performed.

END OF EDIT

[1] J.-M. Muller and L. Rideau, “Formalization of double-word arithmetic, and comments on ”Tight and rigorous error bounds for basic building blocks of double-word arithmetic”,” ACM Transactions on Mathematical Software, vol. 48, no. 1, pp. 1–24, Mar. 2022, doi: 10.1145/3484514. https://hal.archives-ouvertes.fr/hal-02972245