Math'375 Fall 2005 4' Errors in Numerical Computation

1
Math.375 Fall 20054. Errors in Numerical
Computation

• Vageli Coutsias

2
Types of errors in numerical computation

• Roundoff errors
• Pi 3.14159
• Pi 3.1415926535897932384626
• Truncation errors
• Cos x 1 x2/2
• Cos x 1 x2/2 x4/4!
• Errors usually accumulate randomly
• (random walk)
• But they can also be systematic, and the reasons
• may be subtle!

3
x linspace(1-210-8,1210-8,401) f
x.7-7x.621x.5-35x.435x.3-21x.27x-1
plot(x,f)
CANCELLATION ERRORS
4
g -1x.(7x.(-21x.(35x.(-35x.(21x.(-7
x)))))) plot(x,g)
5
h (x-1).7 plot(x,h)
6
plot(x,h,x,g,x,f)
7
z(1) 0 z(2) 2 for k 229 z(k1)
2(k-1/2)(1-(1-4(1-k)z(k)2)(1/2))(1/2) end
semilogy(130,abs(z-pi)/pi)
8
Plot derivative error vs. h to
exhibit optimal step clear fname 'sin'
fname1 'cos' delta eps hmin 106eps
hmax 1010eps xmax 1 xpi/4 M21
estimate for second derivative hopt
2sqrt(delta/M2) h for minimum error yopt
(M2/2)hopt2delta/hopt d (feval(fname,xhopt
)-feval(fname,x))/hopt
9
h linspace(hmin,hmax,1000) y1 (M2/2)h
truncation error y2 2delta./h
roundoff error y y1 y2 total
error derivh (feval(fname,xh) -
feval(fname,x))./h abserr abs(feval(fname1,x)-d
erivh) loglog(h,y,h,y1,h,y2,h,abserr,hopt,yopt,'
bo') xlabel(sprintf('hopt 3.2d,
ds(3.2d)/dx 3.2d
/-3.2d',hopt,fname ,x,d,err)) ylabel('abserr,
maxerr') title(sprintf('derivative error as
function of h for
s(x)',fname)) axis(hmin hmax 10(-12) 10(-4))
10
(No Transcript)
11
Vocabulary

• Consistency correctness -
• Convergence (rate)
• Efficiency (operation counts)
• Numerical Stability (error amplification)
• Accuracy (relative vs. absolute error)
• Roundoff vs. Truncation

12
Due to limitations in computer arithmetic, need
to practice defensive coding.

• Mathematics?numerics implementation
• The sequence of computations carried out for the
  solution of a mathematical problem can be so
  complex as to require a mathematical analysis of
  its own correctness.
• How can we determine efficiency?

Text
13
Recurrence relations
14
Instability is inescapable But one can learn to
ride it!
15
(No Transcript)
16
(No Transcript)
17
APPENDIX

• NUMERICAL MONSTERS
• a collection of curious phenomena
• exhibited by computations at the
• limit of machine precision

18
The first monster
for k150, term xterm/k s stermerr(k)
abs(f(k) - s)end relerr err/exp(x)
semilogy(1nTerms,relerr) semilogy(1nTerms,err)
ylabel('Relative Error in Partial Sum.')
xlabel('Order of Partial Sum.')
title(sprintf('x 5.2f',x)) figure
semilogy(1nTerms,err)end term 1 s 1 f
exp(x)ones(nTerms,1) for k150, term
xterm/k s sterm err(k)
abs(f(k) - s) end relerr err/exp(x)
semilogy(1nTerms,relerr) semilogy(1nTerms,err)
ylabel('Relative Error in Partial Sum.')
xlabel('Order of Partial Sum.')
title(sprintf('x 5.2f',x)) figure
semilogy(1nTerms,err)end

• T(x) ex log(1e-1)
• For xgt30 a dramatic deviation from the
  theoretical value 1 is found.
• For xgt36.7 the plotted value drops (incorrectly)
  to zero
• Here
• x linspace(0,40,1000)
• y exp(x).log(1exp(-x))

19
(No Transcript)
20
(No Transcript)
21
(No Transcript)
22
Pet Monster I
23
(No Transcript)
24
(No Transcript)
25
(No Transcript)
26
III.MISCELLANEOUS OPERATIONS (1) Setting
ranges for axes, axis(-40 0 0 0.1) (2)
Superimposing plots plot(x,y,x,exp(x),'o')
(3) using "hold on/hold off" (4) Subplots
subplot(m,n,k)
27
Summary

• Roundoff and other errors
• Multiple plots

28
References

• C. Essex, M.Davidson, C. Schulzky,
• Numerical Monsters, preprint.
• Higham Higham, Matlab Guide, SIAM
• B5 Trailer http//www.scifi.com/b5rangers/