Why are the round-off errors when solving the linear system $Ax = b$ of order $\varepsilon_\text{mach} x_j$?

Question

I was reading a paper on arXiv where, in Section 2.4, the authors are discussing the error that arises in the solution of a linear system $$Ax = b,$$ or, to match up better with the paper, $$\Phi \alpha = u,$$ due to machine precision. The system is overdetermined so they solve it with least squares.

To make a long story short, they state that large coefficients will lead to severe round-off errors that hamper the convergence of a numerical method...

Now for the details. I will write the problem in a general form for clarity. Let $M>N$ and let $\varepsilon_\text{mach} \approx 10^{-16}$ represent machine precision.

They have a system matrix $\Phi=(\phi_{i,j}) \in \mathbb{R}^{M\times N}$, a coefficients vector $\alpha = (\alpha_1,\dots, \alpha_N) \in \mathbb{R}^N$ and an input vector $u^{(N)} = (u_1,\dots, u_M)\in \mathbb{R}^M$, giving the linear system: $$ \sum_{j=1}^N \alpha_j \phi_{i,j} = u_i. $$

The vector $u^{(N)}$ is a numerical approximation of a continuous function $u$, see equation (2) in the paper. They are interested in the convergence of $u^{(N)}$ to $u$ as $N$ increases.

In Section 2.4 they say that because the elements of the system matrix are of $O(1)$, each coefficient $\alpha_j$ will result in round-off errors of approximately $\varepsilon_\text{mach} \alpha_j$ in the numerical approximation $u^{(N)}$. Thus the total error in $u^{(N)}$ due to machine precision is of order $\varepsilon_\text{mach}||\alpha||_{l^2}$.

They say this machine precision error places a limit on the achievable minimimum error of the numerical method. That is, as $N$ is increased eventually the approximation error $t(\alpha) = ||u-u^{(N)}||$ will decrease to approximately $\varepsilon_\text{mach}||\alpha||_{l^2}$ and then stop decreasing.

Thus if the norm of the coefficients $||\alpha||_{l^2}$ is very large this halting of convergence will happen very quickly and really impact the performance of the numerical method.

I get the overall idea but what I don't understand is how they can say that:

...because the elements of the system matrix are of $O(1)$, each coefficient $\alpha_j$ will result in round-off errors of approximately $\varepsilon_\text{mach} \alpha_j$ in the numerical approximation $u^{(N)}$.

Can this statement that the round-off errors are of size $\varepsilon_\text{mach} \alpha_j$ be proven? Even a rough proof? And/or, is it possible to show a simple explicit example that demonstrates that the round-off errors are indeed of this size?

I have spent several hours searching for further information on this, but although I have found many papers/books/documents, that state that large coefficients lead to round-off errors, I have not seen anything that proves or at least gives a detailed explanation that the size of the error due to each coefficient is of $\varepsilon_\text{mach} \alpha_j$?!

Answer 1

The previous and next IEEE machine numbers to $\alpha_j$ are at a distance $\approx |\alpha_j| \varepsilon_{mach}$ from each other; hence $fl(\alpha_j)$ (the closest machine number to $\alpha$) is at a distance at most $|\alpha_j| \varepsilon_{mach}$ from $\alpha_j$, but that's just an upper bound: if $\alpha=0$ or $\alpha=1$, for instance, then it is represented exactly and the error is 0.

I think that what they mean in Section 2.4 is in the reverse direction of what is usually done in error analysis: a perturbation of $\alpha$ to $\alpha+f$, with $\|f\| \leq \varepsilon \|\alpha\|$, gives rise to a perturbed $u+g =\Phi (\alpha+f) = u + \Phi f$ with $\|g\| = \|\Phi f\| \leq \|\Phi\|\|f\| \leq O(1) \varepsilon \|\alpha\|$, hence one expects that perturbations in $u$ of size smaller than that will not change the value of $\alpha$ further, or result only in perturbations in the last significant digit. This is confirmed by the experimental results, and it shows that it makes little sense to compute $u$ with very high precision.

It is not hard to construct an example where all $\leq$ become equalities, for instance take $N=1$ and $\Phi=1$.

It could be rephrased in the more familiar terms of condition numbers: if $\Phi\alpha=u$, then $\sigma_{\min}(\Phi) \leq \|u\| / \|\alpha\|$, and since $\sigma_{\max}(\Phi)=O(1)$, we get $\kappa(\Phi) \geq O(1) \|\alpha\| / \|u\|$. The relative error amplification in solving the linear system is then given by $\kappa(\Phi)$, which (if you use the classical condition number bound and convert between relative and absolute errors properly), shows that absolute errors in $u$ of size $\varepsilon \|\alpha\|$ become absolute errors in $\alpha$ of size at least $\varepsilon \|\alpha\|$, that is, relative errors of the size of roundoff precision. But that's just overcomplicating things, if you look at it in retrospect.