Technische Universit at M unchen Institut f ur Informatik Prof. Dr. H.-J. Bungartz Dipl.-Tech. Math. S.

Schraufstetter

WT 2007/08

February 21, 2008

Numerical Programming I (for CSE)

Final Exam

Name: Matr.Nr.: Program:

.................................................................... .................................................................... ....................................................................

General Instructions
Material: You may only use one hand-written sheet of paper (size A4, on both pages). All other material including electronic devices of any kind are forbidden. Use the exam paper that was handed out to solve the exercises. It is recommendable to start a new page for every problem. Do not use a pencil, or red or green ink. General hint: Often, exercises b), c), etc. can be solved without the results from the previous exercise a): if you are stuck with exercise a, then dont immediately skip exercises b, c, etc. Show the working in your answers and give reasons for your statements or conclusions. Maximum score: The maximum score is 40 points plus a bonus of 5 points. Grades will be computed relatively to a maximum score of 40 points. 17 points are required to pass the exam. Working time: 90 minutes.

1) Floating Point Numbers and Rounding

a) Let f : R R be a mapping dened by f (x) = 5ex 5. Compute the relative rounding error rd(f (x)) f (x) , f (x)

( 4+2 = 6 points)

where rd(f (x)) takes the rounding errors into account when evaluating f (x). Negation does not produce any rounding errors. Is the computer-based evaluation of f stable? b) Consider the following two Matlab programs: function x = function A(p,q) x = -p/2 + sqrt(p^2-4*q)/2; end function x = function B(p,q) z = -p/2 - sqrt(p^2-4*q)/2; x = q/z; end

Both programs return a solution of the quadratic equation x2 + px + q = 0. Evaluating the formulas analytically, you would get the same result x in both cases. Evaluating the formulas with Matlab, you get: input
p 101 102 103 104 105 106 107 108 109 1010 1011 1012 q 1 1 1 1 1 1 1 1 1 1 1 1

output 1
x 0.09901951359278 0.00999900019995 0.00099999900000 0.00009999999900 0.00001000000000 0.00000100000000 0.00000010000000 0.00000001000000 0.00000000100000 0.00000000010000 0.00000000001000 0.00000000000100 x*x+p*x+q 0 0 0 0 0 0 0 0 0 0 0 0

output 2
x 0.09901951359278 0.00999900019995 0.00099999899999 0.00009999999929 0.00001000000339 0.00000100000761 0.00000009965152 0.00000000745058 0 0 0 0 x*x+p*x+q 0 -0.00000000000028 -0.00000000001507 0.00000000292777 0.00000033863576 0.00000761449437 -0.00348484516143 -0.25494194030762 -1.00000000000000 -1.00000000000000 -1.00000000000000 -1.00000000000000

In the second column of both outputs, you can see whether the solution x fullls x2 +px+q 0. Assign the outputs 1 and 2 to the programs function A and function B (give reasons for your decision!) and state which program you would choose for solving the equation x2 + px 1 = 0 for some p 0.

2) Interpolation
Consider the mapping f : [0, 2] [0, 1], f (x) = sin

( 3+3+1+1 = 8 points)
x . 2

a) To approximate the mapping f , nd the interpolant p(x) that interpolates f at the three support points P0 , P1 , and P2 with support abscissas x0 = 0, x1 = 1, x2 = 2.

Compute the divided dierences and give a closed representation of the interpolant p(x) = a0 + a1 (x x0 ) + a2 (x x0 )(x x1 ). with the help of Newtons interpolation formula. b) To get a higher accuracy, another support point P3 with support abscissa x3 = Find the new interpolant that interpolates the four points P0 , P1 , P2 , and P3 .
1 3

is added.

c) Let p1 (u) be the interpolant with general support abscissas {u0 , u1 , u2 , u3 } and p2 (u) be the interpolant with general support abscissas {u1 , u2 , u3 , u4 }. Find the interpolant p3 (u) with support abscissas {u0 , u1 , u2 , u3 , u4 } based on p1 (u) and p2 (u)! d) Which kind of interpolant would you choose, if 50 support points had to be interpolated? Give reasons for your answer!

3) Numerical Quadrature
Consider the integral
b

( 4+1+3+3 = 11 points)
f (x)dx,
a

where f is a function f : R R. a) Formulate a Matlab program function y = simpson sum(f,a,b,n) that computes an approximation of the value of the integral based on the Simpson sum with n subintervals, i.e. with n + 1 function evaluations of f . b) What does the Simpson sum with n = 4 subintervals return for f (x) = x2 , a = 0, and b = 4? c) Approximating the integral
4

I :=
0

x2 dx

with the trapezoidal sum leads to I T1 = 24 in case of n1 = 2 subintervals, I T2 = 22 in case of n2 = 4 subintervals. Compute an extrapolation step to get a better approximation! How many further extrapolation steps (e.g. using the trapezoidal sum with n3 = 8, n4 = 16,. . . subintervals) are reasonable? Why? d) How many partitions (subintervals) n are at least necessary to compute the integral
1

cos(2x)dx
0

with the trapezoidal sum with an error of not greater than You may use the estimates | sin(x)| 1, | cos(x)| 1.

1 ? 75

4) Solving Systems of Linear Equations

Consider the following Matlab code: function A = func(A) for k=1:n for j=1:k-1 A(k,k)=A(k,k)-A(k,j)^2; end A(k,k) = sqrt(A(k,k)); for i=k+1:n for j=1:k-1 A(i,k) = A(i,k)-A(i,j)*A(k,j); end A(i,k)=A(i,k)/A(k,k); end end end

( 2+3+2+2 = 9 points)

a) Which algorithm is implemented by the program func(A) and what does it return as result? What is the complexity of the algorithm? b) Find an LU factorization of 4 2 3 A = 2 3 0 . 4 4 9 c) Give two reasons for pivoting (in general)! d) Why are direct solving methods not suitable for large sparse systems of linear equations?

5) Solving Ordinary Dierential Equations

Consider the initial value problem y (t) = (y (t))2 y (1) = 1. It is known that y (t) > 0 for all t 1.

( 2+3+2+4 = 11 points)
t 1,

a) Solve the initial value problem analytically via separation of the variables. b) Let yk be a numerical approximation of y (tk ). Apply the implicit Euler method (backward Euler method) to nd an approximation yk+1 y (tk+1 ) in tk+1 = tk + h. 1 to nd an approximation Perform one step of the implicit Euler method with stepwidth h = 4 of the solution of the initial value problem at t = 1.25 (note that y (t) > 0). c) Instead of solving the implicit Euler formula directly, the value yk+1 could also be computed iteratively by solving the non-linear equation numerically. The three plots of gure 1 show the error reduction (depending on the number of iteration steps) for the following three methods exemplarily when solving the non-linear equation of b): (A) the bisection method (with initial interval [0, 10]), (B) Newtons method (with initial value 10) (C) the secant method (with initial interval [0, 10]). Assign every method one of the plots! d) Apply Newtons method to solve the non-linear equation of b) for yk+1 numerically: For this purpose, specify an appropriate function g and show that the application of Newtons method to the function g leads to the iteration rule
(j +1) yk+1

:=

h (yk+1 )2 + yk 2hyk+1 + 1
(j )

(j )

,
(0)

where j denotes the index of Newtons iteration. What might be a good initial value yk+1 ? Why is the choice of the initial point important?

plot 1

plot 2

plot 3 Figure 1: Solving non-linear equations: error plotted against the number of iteration steps (semilogarithmic plots)