G12ISC Introduction to Scientific Computation — Coursework 4 (5%) — Submission deadline: 4pm, Tuesday, 18 Feb 2020 This coursework contributes 5% towards the overall grade for the module. Rules: • Each student is to submit their own coursework. • You are allowed to work together and discuss in small groups (2 to 3 people), but you must write your own coursework and program all code by yourself. Coursework Aim: To develop Matlab code related to polynomial interpolation. Getting Started: • Download the contents of the “Coursework 4 Pack” folder from Moodle. • Do the coursework by simply adding Matlab code to the downloaded file cw4.m, and by completing the incomplete function m-files: LagrPolyn.m, PolynInterp.m, PolyInterpolError.m, GridCheb.m • Your responses to the numbered questions below should be written in cw4.m, except the responses to questions 4, 6, and 8, which involve the completion of the other m-files. Advice: • A total of 20 points can be obtained in this coursework. The weighting is enclosed in brackets for each question, e.g. [2 / 20]. • Write your code as neatly and concisely as possible so that it is easy to follow. • Add some comments to your scripts that indicate what a piece of code does. • Remember to format output using for example disp, num2str or fprintf. • Some marking criteria: Correctness of numerical output, clarity of code and comments, com- pleteness of figures (title, axis labels, legend, etc) and tables, clear and concise answers to open questions, etc. See Moodle: Assessed Coursework 4 for all grading criteria per question. Submission Procedure: The coursework is to be submitted using the “Assessed Course- work 4” assignment on Moodle. A complete submission consists of the following files: • All m-files: LagrPolyn.m, PolynInterp.m, PolyInterpolError.m, GridCheb.m • 1 PDF file: cw4.pdf This PDF file is created using Matlab “Publish” on cw4.m (see for ex- ample Section 1.2.9 of the Matlab Guide (by Tom Wicks) for instructions on publishing scripts). I Introductory Questions Let f(x) = ex − 1, and let x0 = −1, x1 = 0 and x2 = 2 be three points. 1 [2 / 20]Find the (quadratic) Lagrange polynomials L0(x), L1(x), L2(x) based on the points x0, x1 and x2. (You can do this by hand calculation, so no need to use Matlab.) • Use Matlab to plot these 3 Lagrange polynomials in one figure. Clearly in- dicate which graph corresponds to which Lagrange polynomial (using for exam- School of Mathematical Sciences — Version: 29th January, 2020 — Lecturer: Dr Andrea Cangiani Page 1 of 4 ple the Matlab command “legend” and different plotting styles per graph; see “help plot”). 2 [2 / 20]Use Lagrange’s form to find the polynomial interpolant of degree at most 2 that interpolates f at x0, x1 and x2. (You can do this by hand calculation, so no need to use Matlab.) Next, write down the Taylor (Maclaurin) Polynomial of the function f(x) = ex − 1 of degree 2 centred at x = 0. • Use Matlab to plot f , the Lagrange’s form of the polynomial interpolant of degree at most 2, and the Taylor polynomial of degree 2 in one figure. Hint: Use a sufficient number of evaluation points x in the interval [−1, 2] (e.g. 151 uniformly-distributed points) to evaluate and plot each function. I General Lagrange Polynomials We now proceed to the general case of polynomial interpolation using n+1 points. Let the interval [a, b] contain n+1 uniformly-distributed (equally-spaced) points: zj := a+ (b− a) jn for j = 0, 1, . . . , n. (1) Let L0(x), L1(x), . . . , Ln(x) denote the corresponding Lagrange polynomials. 3 [1 / 20]By using theMATLAB command linspace or otherwise, generate the vector of points z = (z0, . . . , zn) storing the n + 1 uniformly distributed points in [a, b] defined in (1). You may set n = 3 and [a, b] = [−1, 2] and check that the result is as expected. 4 [3 / 20]Modify the function m-file “LagrPolyn.m”, which implements Lk(x), and has the following prototype: function y = LagrPolyn(k,x,z) This function should construct the Lagrange polynomial Lk based on the n+1 points z = (z0, . . . , zn), and evaluate it in the point x, thereby returning the value y = Lk(x). Make sure your function also works for a vector of points x, in which case it should return the vector of values y = Lk(x). Hint: Recall that the formula for Lk(x) contains a product of n terms. There- fore, it is recommended to use a “for” loop to gradually construct Lk. The file LagrPolyn.m from Moodle should help you get started. 5 [1 / 20]• Using your m-file with n = 3, a = −1 and b = 2, plot L0(x), L1(x), L2(x), L3(x) in one figure. Clearly indicate which graph corre- sponds to which Lagrange polynomial. I The General Polynomial Interpolant Let Pn(x) denote the interpolating polynomial of degree at most n that interpolates a function f at n+1 points {zj}. 6 [2 / 20]Modify the function m-file PolynInterp.m that implements Pn(x). The function Coursework 3 Page 2 of 4 must call your other m-file “LagrPolyn.m”, and has the following prototype: function y = PolynInterp(f,x,z) where f is a function handle and z is a vector with n + 1 entries (containing points {zj} used to construct the Lagrange polynomials). This function should evaluate Pn in the point x, thereby returning the value y = Pn(x). Make sure this function also works for a vector of points x, in which case it should return the vector of values y = Pn(x). 7 [2 / 20]Using your m-file PolynInterp.m, plot in one figure the function f(x) = ex − 1 and the different interpolants Pn(x) for n = 1, 2, 4, and 8 generated using n + 1 equally spaced interpolation points in the interval [−1, 2]. 8 [1 / 20]Modify the function m-file called PolyInterpolError.m that computes a numerical approximation to the∞-norm of the error: ‖f − Pn‖∞,[a,b] = max a≤x≤ b ∣∣f(x)− Pn(x)∣∣ . This function has the following prototype: function e = PolyInterpolError(a,b,f,z) where f is a function handle, z is a vector with n+1 entries (containing the points {zj} used to construct the Lagrange polynomials). The output of the function is e = ‖f − Pn‖∞,[a,b]. Hint: Use a sufficient number of points to sample the error f(x)−Pn(x). The Matlab functions “abs” and “max” may be helpful. 9 (?) [2 / 20]• Using your m-file PolyInterpolError.m, plot ‖f − Pn‖∞,[−1,2] versus a sufficiently large range of n, using a logarithmic y-scale. Hint: “help semilogy”. • [Bonus question, not counting for marks.] Let Tn(x) be the Taylor (Maclaurin) polynomial of f centred at x = 0 of order n. Compute the error ‖f − Tn‖∞,[−1,2] for the same range of n and add to the above plot (again in logarithmic y-scale) the computed errors versus n. What do you observe? • Does Pn converge to f , in other words, does ‖f − Pn‖∞,[−1,2] → 0 as n → ∞? If yes, explain how you got to this conclusion. If not, it may be that ‖f − Pn‖∞,[−c,c] → 0 for a sub-interval [−c, c] ⊂ [−1, 2]. In that case, determine (approximately) the largest value of c. 10 [2 / 20]Consider two new functions f = sin(2pix) and g(x) = sin(2pix) + 5 104 sin(10pix) both defined for x ∈ [−1, 1]. Plot the Lagrange polynomial of order n = 22 of the functions f and g. What do you observe? Hint: Think of the function g as a tiny perturbation of the function f . 11 (?) [2 / 20]Instead of using uniformly-distributed points, let us now consider the following Coursework 3 Page 3 of 4 special interpolation points (so-called Chebyshev points): zj = 1 2 (a+ b) + 1 2 (b− a) cos ( 2j + 1 2n+ 2 pi ) for j = 0, 1, . . . , n . (2) Modify the function m-file “GridCheb.m”, which produces Chebyshev points {zj}nj=0 in [a, b] defined via (2). The function has the following prototype: function z = GridCheb(n,a,b) • Repeat Question 10, keeping n = 22, but instead of using uniformly-distributed points, use Chebyshev points computed with GridCheb.m from Question 11. What do you observe? Coursework 3 Page 4 of 4
