CHW3, Math 1322, Fall 97, Due Tuesday, 9/23
So far most of the computations you have done with Matlab could
have been done with a good calculator. The power of software like Matlab
becomes apparent when your calculations involve multiple Matlab commands
as, for example, they do when calculating an approximation to an integral
by Simpson's Rule. Suppose further that you want to perform Simpson's Rule
several times, for different values of n, and then you want to carry it out
on different functions and over different intervals. This homework and the
next will introduce the use of scripts with what is called flow control.
You should read the chapter "Flow Control" in the "Getting Started"
document in "helpdesk". Also read the subsections "M-Files" and
Concatenation" in the chapter "Working with Matrices".
Example. Consider the integral from 2 to 5 of the function y =
xexp(-x/5). Verify that an antiderivative of y is F = -5exp(-x/5)(x + 5),
so that the integral is equal to A = F(5) - F(2) = -5(10exp(-1)
-7exp(-2/5)). We want to approximate this intergral by Mn using the
Midpoint Rule for n = 5, 10, 20, 40. We then want to calculate the error
En = |A - Mn| and see how this changes with n by calculating En/E(n+1). By
the error bounds in (5) on page 460, we know that these ratios should all
be roughly equal to 4, since we are doubling n at each step. We want to
exhibit the calculations in a table similar to the tables on page 460.
Namely, we want three rows and four columns with the top row being the
values of n, the second row being the values of Mn and the last row being
the values of En. Another table with one row and three columns should
display the three ratios. Here is a script that does all this.
%Script for CHW3 to carry out Midpoint Rule for the integral of
%xexp(-x/5) from a to b for the cases of n = 5, 10, 20, 40,
%then calculate the errors and the ratios of successive errors.
%This script will be put on your diskette in the a directory as a
%script m-file with name chw3p1a.m
A = -5*(10*exp(-1) - 7*exp(-2/5)) %the actual value of the integral
N = []; M = []; E = []; R = []; T = []; %rows of the tables start empty
a = 2; b = 5; %the endpoints of the interval
%format long %Use this statement if you want more than 4 places beyond the decimal
%loop to calculate Midpoint Rule approximations for each n
for k = 1:4,
n = 5*2^(k-1); %the values of n
h = (b-a)/n; %deltax
%take the regular partition of [a,b] into n subintervals
x = [a + h/2: h: b - h/2]; % the intercalation = the midpoints
y = x.*exp(-x/5); %values of the function at the midpoints
N = [N n]; %first row of the table by concatenation
M = [M h*sum(y)]; %second row, the Midpoint Rule approximations
E = [E abs(A - M(k))]; %third row, the errors
end; %end of the loop calculating the Midpoint Rule approximations
T = [N;M;E] %the table
for k = 1:3, %loop to calculate error ratios
R = [R E(k)/E(k+1)];
end;
R %table of error ratios
One types this script into a script m-file by going to the Matlab
command window, clicking on the icon on the far left of the tool bar (piece
of paper with upper right corner folded down) for New m-file and then
typing in the script as it appears above. You would save it to your
diskette with a name following some systematic naming scheme. For example,
a:/chw3p1a.m is a name indicating computer homework 3, problem 1 and a for
this preliminary sample script as opposed to the one you will write. Now
the name, without the extension and without the path (provided a: is on
your path), written on the Matlab command line, or as a Matlab command in
your notebook, will execute this script by going through it line by line.
Let's try it here.
chw3p1a
A =
5.0672
T =
5.0000 10.0000 20.0000 40.0000
5.0733 5.0687 5.0676 5.0673
0.0060 0.0015 0.0004 0.0001
R =
3.9968 3.9992 3.9998
It worked!
Recall that the actual value of the integral, found by the Fundamental
Theorem of Calculus, is A. The number of subintervals in the regular
partition is the value of n in the first row of T. The approximations by
the Midpoint Rule for the values of n are in the second row of T . The
values of the errors are given in the third row of T. The vector R gives
the ratio of errors as n is doubled. Notice that the ratios are all very
close to 4 = 2^2, as one should expect from the error bound in (5) on page
460.
The Assignment
In this homework you are to write a script which carries out
Simpson's Rule for n = 6, 12, 24 and 48 to approximate A, the integral from
1 to 3 of x*ln(x). Detailed steps are below. Next week we will consider
important elaborations. For example, how can we adjust the script so that
the user can enter the endpoints a and b? Shouldn't we be able to enter
the function, too? Can we set up for the user to choose the number of
steps? Maybe we should have two scripts, one to perform Simpson's Rule for
a value of n, a function and interval all to be input by the user, and the
second script to make the desired table of approximations, errors and error
ratios. The separate Simpson's Rule script might be useful in other
contexts.
1. Verify that F(x) = x^2*(ln(x) - .5)/2 is an antiderivative of x*ln(x).
Use it to calculate the value of A by hand.
2. From the Matlab command
window open a new m-file on which you will type your script. Name it and
save it to your diskette. Make sure the a: drive is on Matlab's path.
Your script should have lots of comment lines in it to explain what each
step is doing. It should do the following:
3. Use a for-loop to: carry
out Simpson's Rule to approximate the integral A from 1 to 3 of xln(x), for
the cases where the number of subintervals is n = 6, 12, 24, 48;
4. make a
table whose first row consists of the values of n, whose second row
consists of the values of Sn and whose third row consists of the values of
the error |A - Sn|.
5. Use a for-loop to make a table of the values of
the ratios of the successive errors.
6. In your notebook enter the input
cell filename, where filename is the name you gave your script m-file.
After the output, comment on how the table of ratios compares with what one
expects because of formula (7) on page 464.
Hand in your notebook with a printout of the script attached.