# Why is subtraction a stable operation?

In Numerical Linear Algebra by Trefethen & Bau, it is claimed that subtraction is backward stable. Here is the proof:

Let $$f(x, y) = x-y$$ and let $$\tilde f(x,y)$$ be the answer you get when doing $$x-y$$ on a computer. We have \begin{align*}\tilde f(x, y) &= (x(1+\epsilon_1) - y(1 + \epsilon_2))(1 + \epsilon_3)\\ \end{align*} where $$|\epsilon_1|, |\epsilon_2|, |\epsilon_3| < \epsilon_m = \epsilon_{machine}$$. Now distributing, we have \begin{align*}(1+\epsilon_2)(1+\epsilon_3) &= 1 + \epsilon_2 + \epsilon_3 + \epsilon_2 \epsilon_3\\ &=1 + \epsilon_4 \end{align*} where $$|\epsilon_4| \le 2 \epsilon_m + \epsilon_m^2 \le 3 \epsilon_m$$, if $$\epsilon_m < 1.$$ Similarly, $$(1+\epsilon_2)(1+\epsilon_3) \le 1 + \epsilon_5$$ where $$\epsilon_5 \le 3 \epsilon_m$$ if $$\epsilon_m \le 1.$$

Thus \begin{align*} \tilde f(x, y) = x(1 + \epsilon_4) - y(1 + \epsilon_5) = f(\tilde x, \tilde y). \end{align*}

Now here's they key part. Using the 2-norm, the relative error is $$\frac {||(x, y) - (\tilde x, \tilde y)||}{||(x, y)||} = \frac {||(\epsilon_4, \epsilon_5)||}{||(x, y)||} \le \frac {3 \sqrt2}{||(x-y)||} \epsilon_m = c(x, y) \epsilon_m$$ The problem is that the coefficient of $$\epsilon_m$$ above depends on the input $$(x, y)$$. There doesn't exist one $$x$$ that works for all $$(x, y)$$. But this is what is needed for backward stability. So why is subtraction backward stable?

• Is your last equation line correct in the first equality? I see $(x,y)-(\tilde x,\tilde y)=(-ϵ_4x,-ϵ_5y)$ so that the norm has an upper bound of $3ϵ_m\|(x,y)\|$ and the coefficient gets to be independent of $(x,y)$. Nov 3, 2021 at 8:23
• @LutzLehmann Thanks! Jun 30, 2022 at 18:37

(My apologies for the terrible typesetting)

Let me start with the definition of backward stability of an algorithm:

An algorithm $$A$$ to solve a problem $$P$$ is called backward stable, if for all $$x$$ in the input space there exists $$y$$ such that $$||y-x||/||x||\approx \epsilon_m$$ and $$A(x) = P(y)$$.

Trefethen & Bau summarizes this as "Backward stable algorithms give exact answers to nearly right questions" (or something like that, I don't have the book. I am looking at my lecture notes, and I am not a good note-taker).

Now, the problem of subtraction can be written as $$P(X)=x_1-x_2$$ where $$X=(x_1,x_2)$$. So you need to show that given the algorithm you described, i.e. $$A(X) = \tilde{f}(x_1,x_2)$$, for all $$X$$ in the input space there exists a $$Y$$ such that $$||X-Y||/||X||\approx \epsilon_m$$ and $$A(X)=P(Y)$$.

If we choose $$||\cdot||$$ as the $$\infty-$$norm, then the proof you provided above (except for the last step) gives us the answer:

$$y_1 = x_1(1+\epsilon_4)$$, which implies $$y_1-x_1 = x_1\epsilon_4$$ which in turn implies $$|y_1-x_1|/|x_1| = |\epsilon_4|$$ (similar calculation for $$y_2$$, $$x_2$$). From this observation, we can conclude $$||X-Y||/||X||=\frac{\max\{|y_1-x_1|,|y_2-x_2|\}}{\max\{|x_1|,|x_2|\}}\approx \epsilon_m$$

To see why this holds, let's consider all 4 cases:

1. $$\max\{|y_1-x_1|,|y_2-x_2|\} = |y_1-x_1|$$ and $$\max\{|x_1|,|x_2|\} = |x_1|$$: $$||X-Y||/||X||=\frac{\max\{|y_1-x_1|,|y_2-x_2|\}}{\max\{|x_1|,|x_2|\}} = \frac{|y_1-x_1|}{|x_1|}=|\epsilon_4|$$

2. $$\max\{|y_1-x_1|,|y_2-x_2|\} = |y_2-x_2|$$ and $$\max\{|x_1|,|x_2|\} = |x_2|$$: $$||X-Y||/||X||=\frac{\max\{|y_1-x_1|,|y_2-x_2|\}}{\max\{|x_1|,|x_2|\}} = \frac{|y_2-x_2|}{|x_2|}=|\epsilon_5|$$

3. $$\max\{|y_1-x_1|,|y_2-x_2|\} = |y_1-x_1|$$ and $$\max\{|x_1|,|x_2|\} = |x_2|$$: $$||X-Y||/||X||=\frac{\max\{|y_1-x_1|,|y_2-x_2|\}}{\max\{|x_1|,|x_2|\}} = \frac{|y_1-x_1|}{|x_2|}\leq \frac{|y_1-x_1|}{|x_1|} =|\epsilon_4|$$

4. $$\max\{|y_1-x_1|,|y_2-x_2|\} = |y_2-x_2|$$ and $$\max\{|x_1|,|x_2|\} = |x_1|$$: $$||X-Y||/||X||=\frac{\max\{|y_1-x_1|,|y_2-x_2|\}}{\max\{|x_1|,|x_2|\}} = \frac{|y_2-x_2|}{|x_1|}\leq \frac{|y_2-x_2|}{|x_2|} =|\epsilon_5|$$

• This is very clear and easy to follow. Note that the relative errors $\epsilon_4$ and $\epsilon_5$ can be negative, so technically an absolute value sign is needed. Nov 3, 2021 at 13:23
• Thank you for the very detailed response. Nov 3, 2021 at 13:47
• I think using the 2-norm would save us some work, since we don't have to split into cases. Nov 3, 2021 at 13:55
• @CarlChristian Good observation. I corrected it now. Thank you :) Nov 4, 2021 at 1:37
• @user56202 Even the 2-norm proof will have some sort of $\max$ in it. We, mathematicians, have a dangerous but simple trick to avoid checking cases: "w.l.o.g. assume that $\max\{a,b\} = a$" :D Nov 4, 2021 at 1:42