ACMS 40390: Fall 2012 - Numerical Analysis

Course Information

Instructor: Zhiliang Xu, zxu2@nd.edu Hayes-Healy 226

Office phone number: (574)631-3423

Class time: M, W, F, 8:30am - 9:20am.

Classroom: Edward J. DeBartolo Hall 117

Teaching assistant: Dong Lu, Dong.Lv.2@nd.edu, Hayes-Healy 215

Office hours:

Xu: M 5:00pm - 6:30pm or by appointment and drop-in (when available) at HH226
Lu: R 3:00pm - 4:00pm at HH215 Textbook: Richard L. Burden and J. Douglas Faires, Numerical Analysis, 9th Edition.

Prerequisites:

MATH 20750 or MATH 20860 or MATH 30650 or ACMS 20750 or PHYS 20452

The course requires a moderate amount of programming. C, C++ or Fortran programming languages are preferred. However, you may also use software programs including Matlab, Mathematica.

Syllabus is here.

Tips about programming:

Important Dates
First test review Monday, 09/24 in class
First in class test Wednesday, 09/26 in class
Second test review Monday, 10/22 in class
Second in class test Wednesday, 10/24 in class
Final test review Wednesday, 12/05 in class
Final test Wednesday, 12/12 117 DeBartolo Hall
Office Hrs for Final Test 1:00pm - 7:00pm, Monday (12/10); 1:00pm - 7:00pm (12/11) Hayes-Healy Hayes-Healy 226

Week 1:

  • Wednesday (08/22): Sections 1.1 HW:
    Computer project: Setup your account for computer projects. Compile and run the C++ code on the machine "darrow" and submit your work following the instruction at Instruction and help on Unix systems . The code for practising is ex0-1.cpp .

    Section 1.1 Review of Calculus

    Section 1.2 Arithematic

    Matlab code for finding extreme values of functions .
    Computer Lab: When: Monday 08/27/2012 from 4:00pm to 5:00pm. Where: Fitzpatrick 149 (the Engineering Library)
    Increase AFS space quota: email oithelp@nd.edu to ask for 200 MB AFS disk space if you are NOT engineering students.
  • Friday (08/24): Sections 1.1 HW: Section 1.1: 1(b,c), 2(b), 4(c, you can also use Matlab to solve for the equation f'(x) = 0. To do this, take the above Matlab code and modify.), 9, 16 (Use Matlab or Maple), 21.
    Short review of binary numbers.

    A Short Tutorial on The Binary System

    In the decimal system, we use the digits 0-9 to represent numbers, and things are organized into columns. Take the decimal number 285 for example: 

        H | T | O
    2 | 8 | 5
    such that "H" is the hundreds column, "T" is the tens column, and "O" is the ones column. So the decimal number "285" is 2-hundreds plus 8-tens plus 5-ones.

    The ones column stands 10^0, the tens column stands 10^1, the hundreds column stands 10^2 and so on, so 

          10^2|10^1|10^0
        2 |  8 |  5
    the number 285 is really {(2*10^2)+(8*10^1)+(5*10^0)}.

    The binary system works under the exact same way as the decimal system, except it operates in the base 2 rather than the base 10. In other words, instead of columns being

           10^2|10^1|10^0
    they are 
            2^2|2^1|2^0

    Instead of using the digits 0-9, we only use 0-1 digits.

    Examples: What will the binary number 1011 be in the decimal system?

    Answer: By using the columns, we have 

            2^3|2^2|2^1|2^0
         1 | 0 | 1 | 1
    Therefore 
      1011=(1*2^3)+(0*2^2)+(1*2^1)+(1*2^0)
      = (1*8) + (0*4) + (1*2) + (1*1)
  = 11 (in decimal system)

    Week 2:

  • Monday (08/27): Sections 1.2 HW: Section 1.2: 1(a,b,c), 3(b,d), 15(a,b), 16(use the number given in 15(a,b)).
    Computer Lab: When: Tuesday 08/28/2012 from 4:00pm to 5:00pm. Where: Fitzpatrick 149 (the Engineering Library)
  • Wednesday (08/29): Sections 1.2 HW: Section 1.2: 4(b,c, d), 22, 24

    Section 1.3 Algorithms and Convergence

    Summary of Pseudo-code Language Constructions:

    An algorithm is an ordered sequence of unambiguous and well-defined instructions that performs some tasks.

    Three Categories of Algorithmic Operations

    1. sequential operations(Sequence) - instructions are executed in order
    2. conditional  ("question asking") operations - a control structure that asks a true/false question and then selects the next instruction based on the answer
    3. iterative operations (loops) - a control structure that repeats the execution of a block of instructions

    Computation/Assignment 

            set the value of "variable" to :"arithmetic expression" or
            "variable" equals "expression" or
            "variable" = "expression"

    Input/Output

            get "variable", "variable", ...
            display "variable", "variable", ...

    Conditional

                if  "condition" then
                            (subordinate) statement 1 
                             etc ...
                else
                            (subordinate) statement 2
                            etc ...

    Iterative

                while "condition" 
                            (subordinate) statement 1
                            (subordinate) statement 2 ...

                for "iteration bounds" 
                            (subordinate) statement 1
                            (subordinate) statement 2 ...

  • Friday (08/31): Sections 1.2, 1.3 HW: Section 1.2: 12, 13(a,d), 19(a), 21; Section 1.3: 1(b).

    Section 2.1 Bisection method

    Week 3:

  • Monday (09/03): Sections 1.3, 2.1 HW: Section 1.3: 6(a,c), 7(a, d), 9, 10, 14
    Matlab code for bisection method .
    The sample code for the bisection method is ALG021_bisection.cpp .

    Section 2.2 Fixed Point Iteration

    Main Matlab code for fixed-point iteration , Matlab code of Routine for fixed-point .
  • Wednesday (09/05): Section 2.1 HW: Section 2.1: 9, 11, 14, 19(Computer project)
    Computer project:
    Use Alg. 2.1 (the bisection method) to solve Exercise 2.1.19.
    Submit the code in acms40390/hw1 with ex1-1.cpp.
  • Friday (09/07): Section 2.2 HW: Section 2.2: 2a, 3, 6(Computer project), 11b
    Computer project:
    Use Alg. 2.2 (fixed point iteration) to solve Exercise 2.2.6.
    The sample code for the fixed-point iteration method is ALG022_fixed_pt_iter.cpp
    Submit the code in acms40390/hw1 with ex1-2.cpp.

    Week 4:

  • Monday (09/10): Sections 2.2, 2.3
    HW:
    Section 2.2: 8 (To find the fixed-point numerically, use the stopping criterion which is difference between two successive approximations is less than 10^{-4}.), 16;
    Section 2.3: 2 (write down each step), 6(a,b,c)(You can modify the sample code to solve these problems, and write down the approximations from the last two iterations for each of the problems).

  • Section 2.3 Newton's Method

    Algorithm of the Newton's method.
    The sample code for the Newton's method is ALG023_Newton.cpp
  • Algorithm of the Secant method.
    The sample code for the Secant method is ALG024_Secant.cpp
  • Wednesday (09/12): Sections 2.3, 2.4
    HW:
    Section 2.3: 3(a), 4(a)
    Section 2.4: 6, 7, 8(a)
    Computer project:
    Use Newton's method and Secant method respectively to solve Exercise 2.3.26. Which one converges faster?
    Submit the code in acms40390/hw2 with ex2-1.cpp (for Newton's method) and ex2-2.cpp (for Secant method) respectively.

  • Section 2.4 Error Analysis for Iterative Methods

  • Friday (09/14): Sections 2.4, 2.5, 2.6
    HW:
    Section 2.4: 4 (Do exercises 2(a,b). Use modified Newton's method to compute approximations for p0, p1, and p2 for each of these problems); 10.
    Section 2.5: 4, 5.

  • Section 2.5 Accelerating Convergence

  • Section 2.6 Zeros of polynomials and Horner's method

  • Algorithm for the Horner's method is here.
    The sample code for the Horner's method is ALG027_Horner_2.cpp

    Week 5:

  • Monday (09/17): Sections 2.6, 3.1
    HW:
    Section 2.6: 2(a), Let x0=1, compute f(x0) by Horner's method and find the expression of Horner's decomposition, which is the last eqn of Theorem 2.19.
    2(h), Let x0=2, find the expression of Horner's decomposition, which is the last eqn of Theorem 2.19.
    Section 3.1: 2(a,c).

  • Section 3.1 Interpolation and the Lagrange Polynomial


    The Matlab code for the Lagrange Interpolating Polynomial is Lagrange Interpolating Polynomial and Lagrange basis Polynomial .
  • To use this sample code to solve Example 3.1.1, do the following:
    Download these two files and in Matlab add paths where we save these two files;
    In the command window, excute the following statement to see the curve of the 2nd degree Lagrange polynomial for interpolating points taken from 1/x:
    xi = [2.0, 2.5, 4.0];
    fi = 1./xi;
    lp = lagrange_2(xi, fi);
    xx = 0.5:0.01:5;
    plot(xx, polyval(lp, xx), xi, fi, 'or', xx, 1./xx, '--g');
    To evaluate the value of the interpolating polynomial at a number x, use the command: polyval(lp, x).
  • Remark: the return values from calling " lagrange_2(xi, fi)" are coefficents of the interpolating polynomial from the highest degree to zeroth degree in the form of P_n(x) = a_n*x^n + a_{n-1} *x^{n-1} + ... + a_0. (a_n, ..., a_0) are return values.
  • Wednesday (09/19): Section 3.1
    HW: Section 3.1: 4(analyze exercises 2(a,c)), 6(a, justify how you choose points), 14, 19 (a. Use the Matlab code).

  • Section 3.3 Divided Differences

    Matlab function for computing Newton divided differences table
  • Friday (09/21): Section 3.3
    HW: Section 3.3: 2, 16, 17, and problem: what is the error term of Newton's divided difference interpolating polynomial Pn(x) for interpolating function f(x) at (n+1) distinct points? Justify your answer.

  • Section 3.4 Hermite Polynomial; Section 3.5 Cubic Splines

    Oscillation at the edges of an interval when using high degree polynomial interpolation ( known as Runge's phenomenon ).
    Consider the function: f(x) = 1/(1+25*x*x) on [-1, 1].
    Run the following Matlab code:
    xi = -1.0:0.1:1.0;
    fi = 1./(1.0+25.0*xi.*xi);
    lp = lagrange_2(xi, fi);
    xx = -1.0:0.05:1.0;
    plot(xx, polyval(lp, xx), xi, fi, 'or', xx, 1./(1.0+25*xx.*xx), '--g');
  • Interpolation oscillates toward the end of the interval.

    Week 6:

  • Monday (09/24): Exam review
  • Wednesday (09/26): 1st Exam
  • Friday (09/28) Sections 3.4, 3.5
    HW: Section 3.4: 2(a,c, use Divided differences to construct polynomial), 4(a,c).
  • The cubic spline Matlab code is here . The input data for Ex3.5.3 is here .
  • To run this cubic spline code, open it in Matlab and click on "run" button. Then you will see options to input data set. If you want to input data from a file. The format of the data set in the file should follow the one in the sample data file.

    Week 7:

  • Monday (10/01): Sections 3.5, 4.1
    HW: Section 3.5: 2, 11.

  • Section 4.1 Numerical Differentiations

  • Wednesday (10/03): Section 4.1
    HW: Section 4.1: 2(a), 6(a,c), 8(a,c), 13, 18(a), 26(use the most accurate formula for each of time points).

  • Section 4.3 Elements of Numerical Integration

  • Friday (10/05): Section 4.3
    HW: Section 4.3: 2(a,b,c), 4, 6, 8, 10, 12 (for exercises 4, 6, 8, 10, 12, do problems 2(a,b,c)), 14, 16, 17.

  • Section 4.4 Composite Numerical Integration

    Week 8:

  • Monday (10/08): Sections 4.4, 4.6
    HW: Section 4.4: 1(a,b,c), 3, 5 (for exercises 3, 5, do problems 1(a,b,c)), 9, 11, 22 (use composite Simpson's rule).

  • Section 4.6 Adaptive Quadrature Methods


    The sample code for the adaptive quadrature is ALG043.cpp .

  • Section 4.7 Gaussian Quadrature

  • Wednesday (10/10):

  • Section 5.1 Theory of IVP; Section 5.2 Euler's Method


    HW: Section 4.6: 1(a,c), 2(do 1(a)), 9. Section 4.7: 1(a,c), 2, 3 (do problems 1(a,c) for exercises 2 and 3).
  • Friday (10/12):
    HW: Section 4.7: 5. Section 5.1: 1(a,b,c).

    Week 9:

  • Fall break

    Week 10:

  • Pizza party (10/25)

  • Monday (10/22): Review for the 2nd exam
  • Wednesday (10/24): 2nd exam

  • Friday (10/26):
  • Euler's Method Matlab code .

  • Taylor method of order two to solve Example 1(section 5.3) (Matlab code)


    HW: Section 5.2: 2(a,b), 4 (a,b), 9 (note: part a could be done by using the Matlab code), 11(a) (note: turn in the matlab plots for different h values), 12.

  • Section 5.3 Higher-Order Taylor Method

    Week 11:

  • Monday (10/29):
    HW: Section 5.3: 2(a,b), 4(a,b), 9(a, b, c. Note: modify the sample Matlab code to solve a and c; turn in the Matlab plots of the solutions.).

  • Section 5.4 Runge-Kutta Methods

    4th order Runge-Kutta method to solve Example 3(section 5.4) (Matlab code)

    Modified Euler's method (C++ code)

  • Wednesday (10/31):
    HW: Section 5.4: 2(a,c), 4(a, implement a Matlab code to solve. Turn in the solution plot which consists of both approximate and exact solutions. Hint: take the sample 4th order Runge-Kutta Matlab code and modify it), 6(do problems 2(a,c)), 14 (do problems 2(a,c)).

    Computer project:
    1. Do Exercise 5.4.28(a). Submit the code in acms40390/hw5 with ex5-1.cpp (Hint: Use the modified Euler's method C code as the base and see how different stages are implemented in the 4th order Matlab code)
    2. Do Exercise 5.4.8(do Prob 5.4.4a). Submit the code in acms40390/hw5 with ex5-2.cpp (Hint: Modify the sample modified Euler method C code)
  • Friday (11/02):
    HW: Section 5.6: 2(a,b) (Use two-step and four-step Adams-Bashforth explicit methods respectively. You do not need to compute solutions over theentire time interval. The two-step explicit method computes w2, w3, and w4; The four-step explicit method computes w4, w5 and w6. Use 4th order RK method to compute the starting values for each method respectively)

  • Section 5.6 Multistep Methods & Section 5.10 Stability (part 1)

    The Adam-Bashforth four step explicit method C++ code .
    The 4th order predictor-Correction method C++ code .

    Week 12:

  • Monday (11/05): Section 5.10.
    HW: Section 5.10: 3.

    Computer project:
    Use Alg. 5.4 (predictor-corrector method) to compute approximate solution to Extercise 5.6.3(b). Compare with exact solution. Submit the code in acms40390/hw6 with ex6-1.cpp.
  • Wednesday (11/07): Sections 5.10, 5.11.
    HW: Section 5.10: 7, 8. Section 5.11: 2(a,b, use the Matlab code to do the computation), 10.

  • Section 5.10 Stability (part 2) & Section 5.11

  • RK4 to solve system of stiff ODEs (Matlab code)

  • Section 6.1 & Section 6.2

  • Friday (11/09): Sections 5.11.
    HW: Section 5.11: 12, 15.

    Week 13:

  • Monday (11/12): Sections 6.1 and 6.2.
    HW: Section 6.1: 3(a), 6(a,c), 9. Section 6.2: 4(do problems 2(a,c)), 14(do problem 10(a)), 15(do problem 9(b)).
    Algorithm for Gaussian elimination with backward substitution method.
    Gaussian elimination with backward substitution for solving linear systems of eqns code .
    Gaussian elimination with partial pivoting C++ code .
  • Wednesday (11/14): Section 6.2, 6.3, 6.4, 6.5
    HW: Section 6.2: 20 (do problems 10(a,b)), 31. Section 6.3: 6(a,c). Section 6.4: 6.

    Computer project:
    Use Alg. 6.3 (Gaussian elimination with scaled partial pivoting) to solve exercise 10.c of Section 6.2 (Hint: Modify the code for Gaussian elimination with partial pivoting). Submit the code in acms40390/hw7 with ex7-1.cpp.
  • Friday (11/16): Section 6.5
    HW: Section 6.5: 4(b,c), 9(b,c).

  • Section 6.3 - Section 6.5

    Week 14:

  • Monday (11/19): Section 6.6
    HW: Section 6.6: 2(a,b,c),4(d, Solve by computer code), 6(do problem 4(d), solve by computer code), 12(c), 17.

    Computer project:
    Use Alg. 6.4 (LU factorization) to solve problem 6.5.7d. ( Hint: Modify the LU factorization code by adding forward and backward substitution steps). Submit the code in acms40390/hw7 with ex7-2.cpp.

  • Section 6.6 Special Types of Matrices

    Algorithms for LU, LL^t, LDL^t and Crout method for tridiagonal matrices


    LU factorization C++ code
    LDL^t factorization C++ code
    LL^t Cholesky factorization C++ code
    Crout factorization for tridiagonal linear system C++ code

  • Sections 7.1 & 7.2

  • Wednesday (11/28): Section 7.3
    HW: Section 7.3: 2(b,c), 4 (do exercises 2(b,c)).
    Computer project:
    Use Alg. 7.1 (Jacobi iterative method) to solve Exercise 7.3.2(d) with TOL = 10^{-3} in l_2 norm. Submit the code in acms40390/hw8 with ex8-1.cpp.
    Use Alg. 7.2 (Gauss-Seidel iterative method) to solve Exercise 7.3.2(d) with TOL = 10^{-3} in l_2 norm. Submit the code in acms40390/hw8 with ex8-2.cpp.

    Section 7.3 Jacobi and Gauss-Siedel Methods

    Jacobi iteration method C++ code .
    Gauss Seidel iteration C++ code .
  • Friday (11/30): Section 7.4
    HW: Section 7.4: 1(b,c), 10(b, Modify the sample code to solve. Use the stopping criterion given in Step 4 of SOR Alg. 7.3).

  • Section 7.4 SOR Methods & 7.5 Error Bounds


    SOR C++ code .

    Week 15:

  • Monday (12/03): Sections 7.5 and 7.6
    HW: Section 7.5: 2(a,c), 4(a,c), 7.

  • Section 7.6 Conjugate Gradient Method


    Preconditioned Conjugate gradient method C++ code .