Lab 1A: Introduction to Iteration.
January 16, 2022
In this lab you will create a script to compute a Taylor polynomial1 and use the plotting utilities in Matlab
to output the result. You will also learn how to implement Horner’s method to evaluate these polynomials.
The Taylor polynomial of a function f(x) about a point x0 is given by
k! (x−x0) . (1) ∑n xk
To write a Matlab function to evaluate this Taylor polynomial about a given point x0 we start by rewriting the formula in a form that is a bit closer to something the computer would understand. Pseudo-code is basically a list of step-by-step instructions in quasi-English (not in Matlab ) to describe an algorithm. Initially it can be quite vague, and then becomes more precise, as it gets closer to computer code. Here we are given the formula for Taylor series, so our pseudo-code is an intermediate step to help us implement the algorithm. Take a few minutes to try to hand-write rough pseudo-code to evaluate (2). See the end of this document for an example (but don’t peek until you’ve tried first).
INITIALIZE Step 1 Set Step 2
Step 3 Output
Implementation By implementation, we mean turning your algorithm into code that will run on a computer. For this problem we will guide you through the process of implementing the solution in Matlab so that you have an idea of what is expected in the future.
1. Start Matlab and create a new m-file. Put in appropriate comments at the top such as your name, date, and what the file will contain (i.e. evaluate Taylor series for ex about the point x0 = 0).
2. We are now ready to convert the main part of the algorithm into Matlab statements. To do this we need to know how to implement a loop statement and how to write a function. To help you get started, the following lines in your m-file will implement the main block of code:
1You can review the theory of Taylor series, truncation, and the error caused by truncation, in your old Calculus textbook: e. g. Stewart Chapter 11.10.
For example, if f(x) = ex, the Taylor polynomial about x0 = 0 of degree n is given by
∑n f ( k ) ( x 0 ) k
Lab 1A: Introduction to Iteration. 2
% =========================================================================
% AUTHOR ….. [Your Name]
% UPDATED …. [Today’s date]
% Evaluate the truncated Taylor series for exp(x) about the point x0 = 0
% x …. Vector of values to evaluate the Taylor polynomial at
% n …. Integer of last term to evaluate in Taylor polynomial
% T : Evaluated Taylor polynomial at points given by x degree n
% ========================================================================= function T = expTaylorPoly(x, n)
% Initialize sum as 0
% Loop over terms in series
for k = 0:n
T = T + x.^k ./ factorial(k);
Try this out, at the command line type expTaylorPoly(1, 10), to see that you get a number close to e1 ≈ 2.718… This will evaluate the first 10 terms of the Taylor series at x = 1. It is important to understand every part of the code so that you can modify it later. A few points to note are:
1. The = is reserved for assigning a new value to a variable so the statement T=T+…; replaces the previous value of the variable T by whatever the old value was added to the rest of the statement. Such a statement only makes sense if the variable T is preset, or initialized, before the statement is encountered. Thus, we need the initialization statement T=0; at the beginning. After we have typed in all the statements that we want executed in the loop (in this case there is only one, but we will want more later), we close the for statement with an end.
2. We use a semicolon ; at the end of each statement. This is to suppress the output of these commands.
3. You should indent your code as you go along. If you are not sure whether you are doing this in a reasonable way, you can get the Matlab editor to auto indent by highlighting a section of your code with the mouse and hitting ctrl-i. This can also be very useful for debugging your program. That is, finding an error. If the auto-indent indents your code in a different way from what you expect, it may indicate the location of an error in your code.
4. Make sure you save your m-file with the name expTaylorPoly.m, so that it corresponds with the name specified in the function definition in the first line of your file.
One of the powerful features of Matlab is its ability to operate on entire vectors of values of x at once. This is the reason that the power in the code has a .ˆ (period caret) rather than just ˆ (caret). This tells Matlab that this operation should be done element-wise rather than being interpreted as a matrix power. Similarly this is the reason the code has a dot before the divides sign / (i.e. ./).
The statement
linspace(-1, 1, 2019)
will generate a vector of 2019 values evenly spaced between −1 and 1 (Try it in the command window). If we assign this vector to a variable x then we will be able to pass it as an argument to our expTaylorPoly function by typing:
Lab 1A: Introduction to Iteration. 3
x = linspace(-1, 1, 2019); T = expTaylorPoly(x, 10);
The 10 indicates that we are taking 10 terms of the Taylor polynomial.
Most Matlab operations are the same for vectors and for scalar numbers. However, some operations can
be ambiguous for vector expressions. For example,
is now asking Matlab to take a 2019 × 1 matrix (i.e. a vector) and multiply it by a 2019 × 1 matrix, an operation that is not defined in linear algebra. In this case, we actually want the operation to be done on each element of the vector x. To do this we use a . in front of the operator:
Now that we have assigned values for x and evaluated our function expTaylorPoly to evaluate the Taylor
series at these values. Then we should plot the results. This can be done by typing:
plot(x, T)
The plot command can also take multiple inputs, so the command plot(x, T, x, x)
will plot the Taylor series, along with the function y = x on the same plot. Try it!
Horner’s Method
A function f(x) about a point x0 can be expressed and evaluated using the following formula:
( ( ( ( ))))
Hn(x)=f(x0)+(x−x0) f′(x)+(x−x0) f′′(x)+(x−x0) ···+(x−x0) f(n)(x0)··· 2 n!
We can express the Taylor polynomial of f(x) = ex about x0 = 0 of degree n using Horner’s method as ( ( 1 ( 1 ( ( 1 )))))
Hn(x)=1+x 1+x 2!+x 3!+x ···+x n! . (4) One more elegant, but useful, refinement makes things even nicer:
x( x( x( x( x )))) Hn(x)=1+1 1+2 1+3 1+4 ···+n(1) .
(5) (check for yourself that this factoring gets all the factorials correct) which we can then write recursively.
Starting with HN = 1,
Hk−1 = 1 + x · Hk , (6) k
where k = N,N − 1,…,1. This may seem confusing at first, but what this means is that we will be working backwards from k = N to k = 1. You should try to create a rough pseudo-code for this now (see the end of this document for an example … but don’t peek now … be creative!:
INITIALIZE Step 1 Set Step 2
Step 3 Output
Lab 1A: Introduction to Iteration. 4
Exercise 1
You should try this. Write your own function expHorner based on the Horner’s method pseudo-code given to you in this lab. Remember you are not submitting this lab, but try to complete these exercises.
Your function should have the following header
% =========================================================================
% AUTHOR ….. [Your Name]
% UPDATED …. [Today’s date]
% Evaluate exp(x) about the point x0 = 0 using Horner’s method
% x …. Vector of values to evaluate the Taylor polynomial at
% n …. Integer of last term to evaluate in Taylor polynomial
% H : Evaluated Taylor polynomial at points given by x degree n
% ========================================================================= function H = expHorner(x, n)
Create a script file (m-file) called plotOfExp.m that produces a single plot on the interval [−4,4] using 50 equally-spaced nodes with the following curves:
• H2(x) • H3(x) • H5(x)
Use different line styles or colours or both for each curve. Your plot should also include a legend (see the legend function in Matlab) that clearly indicates what each curve is. Your plot must include an appropriate title and axis labels. Save the resulting plot as expHornerPlot.pdf
Remark We’re just easing you in to the numerics, here. It’s easy to break this code by looking at using high degree polynomials to try to evaluate exp(x) for x large and negative. We’ll look at this later, but for instance your code should accurately evaluate (say) exp(30) with a degree 90 polynomial, but utterly fail to evaluate exp(−30) accurately with the same polynomial. That’s weird, and we’ll explain it later.
Matlab Tip It is often useful to start script files with the following 3 lines:
• clf will clear any previous figures you have created
• close all will close any plot windows you have open
• clear all will clear all previous variables you have defined
You do not need these 3 lines at the beginning of function files, only script files. Make sure you know the difference. If you include this for all your script files you can be confident that when your code is marked the same results will be produced.
close all;
clear all;
Lab 1A: Introduction to Iteration. 5
Appendix: Examples of pseudo-code
INITIALIZE n,x
Step 1 Step 2
Set T = 0;
Loop k from 0 to n do
Set T = T + x.^k ./ factorial(k) Endloop
Output(‘exp(x) is approximately ‘, T)
% Horner for exponential function
INITIALIZE n,x
Step 1 Step 2
Set H = 1;
Loop k from n to 1 do
Set H = 1 + x. * H ./ k Endloop
Output(‘exp(x) is approximately ‘, H)