Download Lab Exercise: Polynomial Approximation and Quadrature Methods and more Lab Reports Mathematical Methods for Numerical Analysis and Optimization in PDF only on Docsity!
MTH 351 Winter 2009 – Lab 5
Due: Mar 13 by 5pm
- Use lsq.m from the website to find the degree m least squares polynomial approximation to a function via the solution of a linear system. Explore how the error in the linear system solution method GEpivot affects the error in the approximation. For example, test an easy function such as exp(x) on the interval [0, 1].
(a) Try m = 1, 2, 3, and comment on the trend of maximum error on the interval.
(b) Increment m, for what value of m does the maximum error stop decreasing?
(c) Increment m, for what value of m are there visible errors in the approximation?
(d) Increment m a couple more times, comment on the trend of maximum error on the interval.
(e) Can we conclude that low degree polynomial approximations are better than high degree approximations?
(f) Could Jacobi iteration be applied to the linear system? Why or why not?
- Use gaussweights.m from the website to compute the weights for Gaussian quadrature using n nodes on the interval [− 1 , 1]. The code uses the exact values of the nodes (from a table) to determine the linear system to be solved. It applies GEpivot to solve the system, and compares the computed weights to the exact values graphically.
(a) Try n = 2, 4, and 8, and compare the results with Table 5.7 in the book.
(b) Try n = 32 and 64. Explain what goes wrong. It may be helpful to compute the condition number by adding cond(A,inf) to the code. Include these plots in your work to turn in.
(c) Could Jacobi iteration be applied to the linear system? Why or why not?
- Use cubspline.m from the website to solve the tridiagonal system corresponding to natural cubic spline interpolation. The code uses the example of approximating the Legendre polynomial of degree m using n (uniformly spaced) nodes. The linear system is solved using GEpivot. The polynomial and its approximation are plotted, and the residual and computation times are output. In this example, use m = 15. Try n = 8 and 16 to get an idea of what the code is doing.
(a) Try n = 64, 128, and 256. Record the time required to solve the linear system. Notes: Computation time varies due to many factors. Consider running each example a couple times, or modify the code to use cputime if you prefer.
(b) Modify the code to call the file Jacobi.m in order to use Jacobi iteration for the linear solve. Record the time required, and the number of iterations, for n = 64, 128, and 256. Compare the runtimes to those of GEpivot. Notes: Use error tolerance delta=eps to get accuracy upto round-off. Use a sufficiently large maximum number of iterations.
(c) Modify the code to call the file GS.m in order to use Gauss-Seidel iteration for the linear solve. Record the number of iterations required for n = 64, 128, and 256 and compare to those of Jacobi. Notes: Again, use delta=eps. In general one would write Gauss-Seidel in such a way that the time for each iteration is the same as a Jacobi iteration, but this has not been done here, so ignore the timings for Gauss-Seidel.
(d) Modify the code to call the file CG.m in order to use Conjugate-Gradient iteration for the linear solve. Record the runtimes and the number of iterations required for n = 64, 128, and 256, and compare to those of Jacobi. Notes: Use delta=eps. Although we did not discuss Conjugate-Gradient in detail, you can look at the code to get an idea of how it works.
(e) Modify the code to call the file tridiag.m in order to use the tridiagonal version of LU factorization and backsubstitution. Record the runtimes for n = 64, 128, and 256 and compare to those of all of the previous methods. Just for fun, try n = 2048 and compare to the other methods on n = 64. Do not try GEpivot for n = 2048.
You do not need to include these plots in your work to turn in.
