Computer Arithmetic
CS/SE 4X03
Ned Nedialkov
McMaster University
January 19, 2021
Canclellations
Consider x − y.
Assume no roundoff in the subtraction, i.e. fl (x − y) = fl (x) − fl (y). The relative error is
fl(x−y)−(x−y) fl(x)−fl(y)−(x−y)
x−y=x−y (fl(x) − x) − (fl(y) − y)
= x − y ≤ |fl(x)−x|+|fl(y)−y|
|x−y|
If x ≈ y this ratio can be large.
Copyright © 2021 N. Nedialkov. All rights reserved.
2/9
Example 1. Consider a decimal FP system with 5 digits. Let x = 9.23450001 and y = 9.23455001.
Assuming rounding to the nearest, what is the relative error in (a) fl(x+y), (b) fl(x−y)?
x and y are stored as x = fl (x) = 9.2345 and y = fl (y) = 9.2346. (a) x + y = 1.84691 × 10 and fl (x + y) = fl (x + y) = 1.8469 × 10
fl(x+y)−(x+y)= 1.8469×10−1.846905002×10 x + y 1.846905002 × 10
≈ 2.7 × 10−6
Copyright © 2021 N. Nedialkov. All rights reserved. 3/9
Example 1. cont.
(b) x − y = −1.0000 × 10−4
fl (x − y) = fl (x − y) = −1.0000 × 10−4
fl(x−y)−(x−y)= −1.0000×10−4 −(−5.0000×10−5) x−y −5.0000×10−5
−10−(−5) = − 5
=1
Copyright © 2021 N. Nedialkov. All rights reserved. 4/9
Example 2. How to evaluate √x + 1 − √x to avoid cancellations? For large x, √x + 1 ≈ √x.
√ √ √x+1+√x
( x+1− x)√ √ =√
x+1+ x √1√
1 x+1+ x
√
Evaluate
x+1+ x
Let x = 100 000. In a 5-digit decima arithmetic,
x+1=1.0000×105 =100001roundsto1.0000×105. Then √x+1−√x gives 0, but
1 1 −3 √x+1+√x = 1.0000×105 +1.0000×105 = 6.3246×10
Copyright © 2021 N. Nedialkov. All rights reserved. 5/9
Example 3. Which one is better to compute
x2 −y2 or (x−y)(x+y) ?
The following
clear all;
h = 1e-10;
x = sin(1+h); y = sin(1-h);
a = (x-y)*(x+y); b = x^2-y^2;
% evaluate in higher precision
x = vpa(x, 32);
acc = (x-y)*(x+y);
fprintf(’a = %.8e\n’, a);
fprintf(’b = %.8e\n’, b);
fprintf(’accurate = %.8e\n’, acc);
fprintf(’rel error a=%.1e\n’, abs((acc-a)/acc));
fprintf(’rel error b=%.1e\n’, abs((acc-b)/acc));
outputs
a = 1.81859466e-10
b = 1.81859416e-10
accurate = 1.81859466e-10
rel error a=1.0e-17
rel error b=2.7e-07
Copyright © 2021 N. Nedialkov. All rights reserved. 6/9
Example 4. Consider approximating e−x for x > 0 by
x2x3 xk e−x ≈1−x+ − +···(−1)k
2!3! k! From e−x = 1/ex, it is better to approximate
for some k
x x2x3 xk e ≈1+x+2!+3!+···+k!
and then compute 1/ex
Copyright © 2021 N. Nedialkov. All rights reserved. 7/9
Solving ax2 + bx + c
Example 5. Compute the roots of ax2 + bx + c = 0
√
x1,2 = −b± b2 −4ac
2a If b2 ≫ 4|ac|, there may be cancellations
Consider 4-digit decimal arithmetic and take a = 1.01, b = 98.73, c = 4.03
b2 = 9748, 4ac = 16.28, b2 −4ac = 9732
d = b2 − 4ac = 98.65
−b+d = −98.73+98.65 = −0.08, −b−d = −98.73−98.71 = −197.4
x1 = (−b + d)/(2a) = −3.960 × 10−2 x2 = (−b − d)/(2a) = −97.71
Exact roots rounded to 4 digits x1 = −4.084 × 10−2 x2 = −97.71
Copyright © 2021 N. Nedialkov. All rights reserved. 8/9
Example 5. cont. x1x2 = c/a Compute using
√
d= b2−4ac
if b ≥ 0
x1 = (−b − d)/(2a) x2 = c/(ax1)
else
x1 = (−b + d)/(2a) x2 = c/(ax1)
x1 = −97.71, x2 = −4.084 × 10−2
Copyright © 2021 N. Nedialkov. All rights reserved. 9/9