(Partial) Answers to Review Questions of Chapters 0, 1, and 2

You may use a formula sheet: postscript version , and version in pdf format .

1. Consider the representation of floating-point numbers with base 2 and 3 digits in the fraction part. The values for the exponents are between -16 and +16.

   1. What is the machine precision is this number system?

   2. Represent the numbers 7 and 6 as floating point numbers and illustrate the calculation of 7+6 and 7-6.

   Answers:

   1. The machine precision is b^(-p) = 2^(-3).

   2. Adding 7 and 6 in this number system yields

                  7 = 111 = +.111 2^3
                  6 = 110 = +.110 2^3
                ----------------------
                 13 = 1101 = +.110 2^4 = 12 (when truncating)
                           = +.111 2^4 = 14 (when rounding)
      

      Notice that, if one works with a hidden bit (a normalized number always have leading bit equal to one), the result can be represented exactly.

      The subtraction 7 - 6 first yields +.001 2^3 which is normalized into +.100 e^1.

2. Let f(x) = x/(x-1). Consider the following values for x: a random number, a point close to the origin (i.e.: |x| < delta), and a point close to 1 (i.e.: |x - 1| < delta). Answer the following:

   1. What is the difference in relative and absolute error in the evaluation of f(x) for these three points?

   2. Assuming there is a small error of delta on the three points, what is effect of this error on the function value?

   3. Discuss the difference in numerical conditioning between the root finding and the evaluation problem.

   Answers:

   1. The absolute error is |y-f(x)|, where y is the computed value, and the relative error is |y-f(x)|/|f(x)|. We have a small absolute error for x around the root 0, but a large relative error when x is close to 0. Exactly the opposite happens when x is around 1: as f(x) goes to infinity, the absolute error will be large, but the relative error will be small as we divide by a huge number. For random numbers, the difference between absolute and relative error does not matter that much, except if we go to the extremes, close to 0 or 1, or if x becomes extremely large.

   2. The Taylor series expansion at x is

             f(x + delta) = f(x) + delta f'(x) + O(delta^2)
      

      which shows that the error on f(x) can be written as

             f(x + delta) - f(x) = delta f'(x) + O(delta^2)
      

      and depends mostly on f'(x) = -1/(x-1)^2. We see that if x is close to the root or a random point, f'(x) stays moderate. On the other hand, because f'(x) is so steep close to x = 1, we will see huge differences in the function values when there are even small errors on the input value x.

   3. Recall that when f'(x) is tiny around a root, the graph of f(x) is very flat and it is difficult to obtain an accurate numerical approximation of the root. Thus the condition of the root finding problem is bad when f'(x) is small, but good if f'(x) is large. Exactly the opposite situation happens for the evaluation problem, as illustrated by the function f(x).

3. Find an approximation for the square root of 4.12.

   Answers:

   Note that we applying the methods to the equation f(x) = x^2 - 4.12 = 0.

   1. Apply three steps of the bisection method, starting with the interval [2,2.1]. Complete the following table, using four decimal places:

      ---------------------------------------------------------------------------
        k   a(k)    b(k)   (a(k)+b(k))/2   f(a(k))   f(b(k))   f((a(k)+b(k)/2)
      ===========================================================================
        0   2.000   2.100     2.050        -0.1200   0.2900      0.0825
        1   2.000   2.050     2.025        -0.1200   0.0825     -0.01938
        2   2.025   2.050     2.038         0.0825  -0.01938     0.01314
      ---------------------------------------------------------------------------
      

   2. Apply three steps of the secant method, starting with x_0 = 1 and x_1 = 2. Complete the following table, using four decimal places:

      -------------------------------------------------------------------
        k   x(k)   x(k+1)  x(k+2)   f(x(k))   f(x(k+1))   f((x(k+2))
      ===================================================================
        0   1.000   2.000   2.040   -3.120     -0.1200    0.0416
        1   2.000   2.040   2.030   -0.1200     0.0416    -3.059E-4
        2   2.040   2.030   2.030    0.0416    -3.059E-4    ---
      -------------------------------------------------------------------
      

      Note that the precision of four decimal places is insufficient to make a difference in the end.

   3. Use Newton's method, starting at 2, until you are sure you have eight decimal places correct. Make a table with four columns, containing the respective values for k, x(k), f(x(k)), and f(x(k))/f'(x(k)). Write x(k) with as many decimal places as your calculator shows. Give the values for f(x(k)) and f(x(k))/f'(x(k)) in scientific format, using four decimal places.

      ------------------------------------------------------- 
        k       x(k)             f(x(k))   f(x(k))/f'(x(k))
      =======================================================
        0     2                -1.200E-1      -3.000E-2
        1     2.03              9.000E-4       2.217E-4
        2     2.029778325       4.914E-8       1.210E-8
        3     2.029778313
      --------------------------------------------------------
      

      We see eight corresponding decimal places between x(2) and x(3); the last correction term in the table is of order 1.0E-8.

4. Apply three steps of the golden section search method to find then minimum of f(x) = x^3 - x in the interval [0,1]. Write the values for a, b, x1, x2, f(x1), and f(x2) in the table (4 decimal places):

   ---------------------------------------------------------
     k      a       b      x1      x2      f(x1)    f(x2)
   =========================================================
     0   0.0000  1.0000  0.3820  0.6180  -0.3263  -0.3820
     1   0.3820  1.0000  0.6180  0.7639  -0.3820  -0.3181
     2   0.3820  0.7639  0.5279  0.6180  -0.3808  -0.3820
   ---------------------------------------------------------
   

5. What is error propagation? Illustrate the propagation of errors with a method we have seen during the course.

   Answer:

   Error propagation is the accumulation of roundoff error during the execution of a numerical algorithm. We observed error propagation in deflation methods to find all roots of a polynomial (see the project on Bairstow's method). Other instances of error propagation are Gaussian elimination without pivoting and the solution of ODEs (see the last project).

6. The following fixed-point iterations all have the same fixed point, i.e.: 1.

   1. x(k+1) = 3 - 2 x(k)^3, k=0,1,...

   2. x(k+1) = 3/(2 x(k)^2 + x(k)), k=0,1,...

   3. x(k+1) = (4 x(k)^3 + 3)/(6 x(k)^2+1), k=0,1,...

   For each of the three fixed-point iterations, make a cobweb picture illustrating the convergence (or divergence), starting at x(0) = 0.9. Compute the convergence (or divergence) rate for each iteration. Which iteration is best?

   Answers:

   Before making the cobweb pictures, we first best examine the rate of convergence, because this rate determines the slope of the function g(x) at the fixed point, where g(x) defines the fixed-point iteration x(k+1) = g(x(k)).

   1. g(x) = 3 - 2 x^3; g'(x) = - 6 x^2; g'(1) = - 6. As |g'(1)| > 1, this method diverges.

   2. g(x) = 2/(2 x^2 + x); g'(x) = -3(2 x^2 + x)^(-2)(4 x + 1); g'(1) = -15/9. As |g'(1)|, also this method diverges.

   3. g(x) = (4 x^3 + 3)/(6 x^2 + 1); g'(x) = 12 x^2/(6 x^2 + 1) - 12 x (4 x^3 + 3)/(6 x^2 + 1)^2; g'(1) = 0. This method converges (at least) quadratically.

   We do not list the plots here.

7. Consider the linear system

                  -1.000 x_1 - 1.000 x_2 + 1.000 x_3 = -1.000 
                   1.000 x_1 - 1.000 x_2 + 0.000 x_3 = 0.000 
                   2.000 x_1 + 0.000 x_2 + 1.000 x_3 = 3.000
   

   1. Use Gaussian elimination to compute an LU decomposition in the following two ways:

      1. without partial pivoting;

      2. with partial pivoting;

   2. Solve the system two times, using the two LU decompositions obtained from above.

   3. Compute the determinant using the two LU decompositions obtained from above.

   4. Use the second LU decomposition you obtained to compute A^(-1) and cond(A) with ||.||_1.

(Partial) Answers:

1. Here is how we can do the partial pivoting with MATLAB:

   > a = [-1 -1 1; 1 -1 0; 2 0 1]
   
   a =
   
       -1    -1     1
        1    -1     0
        2     0     1
   
   >> b = [-1;0;3]
   
   b =
   
       -1
        0
        3
   
   >> a\b
   
   ans =
   
        1
        1
        1
   
   >> [l,u,p] = lu(a)
   
   l =
   
       1.0000         0         0
       0.5000    1.0000         0
      -0.5000    1.0000    1.0000
   
   
   u =
   
       2.0000         0    1.0000
            0   -1.0000   -0.5000
            0         0    2.0000
   
   
   p =
   
        0     0     1
        0     1     0
        1     0     0
   

2. If we use partial pivoting, then we compute the determinant from det(u)*det(p) = -4*(-1) = 4.

3. The formula for the condition number we use is ||A||_1 * ||A^(-1)||_1. To compute the inverse of A, we use the LU decomposition, solving three linear systems: A x_1 = e_1, A x_2 = e_2, and A x_3 = e3, where x_1,x_2, and x_3 are the columns of the inverse and e_1,e_2,e_3 are the columns of the identity matrix.

FINAL EXAM is in BH 0309 on Monday 29 April 2002, from 1 till 3PM.