|
% Diary file for August 31, 1998
%-----------------------------------------------------------------
% This file was created using the "diary" command:
% diary aug31.out
% Here is the help information for the "diary" command:
help diary
DIARY Save text of MATLAB session.
DIARY file_name causes a copy of all subsequent terminal input
and most of the resulting output to be written on the named
file. DIARY OFF suspends it. DIARY ON turns it back on.
DIARY, by itself, toggles the diary state.
Use the functional form of DIARY, such as DIARY('file'),
when the file name is stored in a string.
%-----------------------------------------------------------------
% You can use Matlab to compute the function value in Project 1:
x = 5
x = 5
sqrt(1+x)
ans = 2.4495
% Matlab computes results using about 16 decimal digits (I have
% to say "about" because the computer uses binary arithmetic,
% and this leads to answers that have between 16 and 17 decimal
% digits). You can see all the digits of an answer by adjusting
% the "format":
format long
sqrt(1+x)
ans = 2.44948974278318
% The default format is "format short":
format short
sqrt(1+x)
ans = 2.4495
% You can learn more about formats:
help format
FORMAT Set output format.
All computations in MATLAB are done in double precision.
FORMAT may be used to switch between different output
display formats as follows:
FORMAT Default. Same as SHORT.
FORMAT SHORT Scaled fixed point format with 5 digits.
FORMAT LONG Scaled fixed point format with 15 digits.
FORMAT SHORT E Floating point format with 5 digits.
FORMAT LONG E Floating point format with 15 digits.
FORMAT SHORT G Best of fixed/floating format with 5 digits.
FORMAT LONG G Best of fixed/floating format with 15 digits.
FORMAT HEX Hexadecimal format.
FORMAT + The symbols +, - and blank are printed
for positive, negative and zero elements.
Imaginary parts are ignored.
FORMAT BANK Fixed format for dollars and cents.
FORMAT RAT Approximation by ratio of small integers.
Spacing:
FORMAT COMPACT Suppress extra line-feeds.
FORMAT LOOSE Puts the extra line-feeds back in.
%-----------------------------------------------------------------
% Most of today's lecture was related to Project 1, performing
% many of the same calculations using f(x) = exp(x), instead of
% f(x) = sqrt(1+x). Matlab has a built-in function to compute
% f(x) = exp(x) = e^x:
exp(1)
ans = 2.7183
% Matlab can be used to derive a Taylor series formula. (For
% further details, see page 18 of your text.) This is a
% "symbolic" calculation, meaning that the result is a formula,
% and not a number. To tell Matlab that you want to perform
% a symbolic calculation, you must declare a symbolic variable:
syms x
% And then you can use the "taylor" function to derive the
% Taylor series. [See the end of this file for information on
% what to do if you are using osf1.]
taylor(exp(x))
ans = 1+x+1/2*x^2+1/6*x^3+1/24*x^4+1/120*x^5
% You can make the result look better:
pretty(ans)
2 3 4 5
1 + x + 1/2 x + 1/6 x + 1/24 x + 1/120 x
% You can also get more terms of the series:
taylor(exp(x),10)
ans = 1+x+1/2*x^2+1/6*x^3+1/24*x^4+1/120*x^5+1/720*x^6+1/5040*x^7
+1/40320*x^8+1/362880*x^9
pretty(ans)
2 3 4 5 6 7
1 + x + 1/2 x + 1/6 x + 1/24 x + 1/120 x + 1/720 x + 1/5040 x
8 9
+ 1/40320 x + 1/362880 x
% The symbolic Taylor series is dramatically different
% than the numeric results we obtained earlier:
exp(1)
ans = 2.7183
exp(5)
ans = 148.4132
% For the rest of this lecture, I want to perform numeric
% calculations. Currently, x is a symbolic variable:
whos x
Name Size Bytes Class
x 1x1 126 sym object
Grand total is 2 elements using 126 bytes
% But I can "clear" this definition, and (from here on)
% use it to compute numerical results:
clear x
whos x
[Matlab prints no information on x, since it has
been "cleared"]
%-----------------------------------------------------------------
% My routine getexp.m (also on the web pages) implements
% the Taylor series formula for exp(x). The initial comments
% are printed if you ask for help:
help getexp
---------------------------------
Use Taylor series to evaluate
exp(x). Continue adding terms
until there is no change in the
answer.
---------------------------------
Usage: y = getexp(x)
---------------------------------
% You should always use the extension ".m" for your Matlab
% programs.
% In getexp.m, I compute each term of the series recursively:
% term = term*x/n
% (this is explained on the lecture overheads).
% An alternative would be to compute each term directly:
% term = x^n/gamma(n+1)
% where the built-in "gamma" function computes the factorial:
help gamma
GAMMA Gamma function.
Y = GAMMA(X) evaluates the gamma function for each element of X.
X must be real. The gamma function is defined as:
gamma(x) = integral from 0 to inf of t^(x-1) exp(-t) dt.
The gamma function interpolates the factorial function. For
integer n, gamma(n+1) = n! (n factorial) = prod(1:n).
See also GAMMALN, GAMMAINC.
Overloaded methods
help sym/gamma.m
% This is fine for small values of n (i.e., for the first
% few terms of the series):
x = 5; n = 3; term = x^n/gamma(n+1)
term = 20.8333
x = 5; n = 30; term = x^n/gamma(n+1)
term = 3.5111e-012
% This notation "3.5111e-012" is in scientific notation,
% and means 3.5111 times 10^(-12).
% For large values of n, the formula breaks down:
x = 5; n = 300; term = x^n/gamma(n+1)
term = 0
% The denominator is "infinite", because the computer
% can only store numbers up to (about) 10^(308):
gamma(n+1)
ans = Inf
% The recursive formula
% term = term*x/n
% avoids this problem.
%-----------------------------------------------------------------
% In Project 1, you have to use a more complicated condition
% for the "while" loop. I will talk about this on Wednesday,
% but you can learn more via:
help while
...
help ops
...
%-----------------------------------------------------------------
% Let's test my routine:
x = 1, y = exp(x), y1 = getexp(x), err = y-y1
x =
1
y =
2.7183
y1 =
2.7183
err =
-4.4409e-016
x = 5, y = exp(x), y1 = getexp(x), err = y-y1
x =
5
y =
148.4132
y1 =
148.4132
err =
5.6843e-014
% It appears that this calculation is less accurate
% when x=5 than when x=1, but if we examine the
% results more carefully, we find that both calculations
% are accurate to about 16 digits; see my lecture overheads.
% If we calculate the relative error instead, we can see
% this:
x = 5, y = exp(x), y1 = getexp(x), err = (y-y1)/y
x = 5
y = 148.4132
y1 = 148.4132
err = 3.8301e-016
x = 10, y = exp(x), y1 = getexp(x), err = (y-y1)/y
x = 10
y = 2.2026e+004
y1 = 2.2026e+004
err = 1.6516e-016
% So far, everything has been fine, but if we try
% negative values of x, the results can be terrible:
x = -10, y = exp(x), y1 = getexp(x), err = (y-y1)/y
x = -10
y = 4.5400e-005
y1 = 4.5400e-005
err = 3.0717e-009
x = -20, y = exp(x), y1 = getexp(x), err = (y-y1)/y
x = -20
y = 2.0612e-009
y1 = 5.6219e-009
err = -1.7275
x = -40, y = exp(x), y1 = getexp(x), err = (y-y1)/y
x = -40
y = 4.2484e-018
y1 = -3.1657
err = 7.4517e+017
% We can plot the function exp(x), and see if it is
% well behaved:
help fplot
FPLOT Plot function.
FPLOT(FUN,LIMS) plots the function specified by the string FUN
between the x-axis limits specified by LIMS = [XMIN XMAX]. Using
LIMS = [XMIN XMAX YMIN YMAX] also controls the y-axis limits. FUN
must be the name of an M-file function or a string with variable x
that may be passed to EVAL, such as 'sin(x)', 'diric(x,10)' or
'[sin(x),cos(x)]'.
...
% Here are a couple of plots, for different ranges of x values:
fplot('exp',[-10 10])
fplot('exp',[-10 1])
% The function is well behaved, but (as I showed above)
% my function getexp.m doesn't always work well:
x = -40, y = exp(x), y1 = getexp(x), err = (y-y1)/y
x = -40
y = 4.2484e-018
y1 = -3.1657
err = 7.4517e+017
% I'll explain why on Wednesday.
exit
38889 flops.
% The term "flops" stands for "floating point operations", and
% is a count on the number of calculations performed in this
% session.
%-----------------------------------------------------------------
% Using Maple on osf1
%-----------------------------------------------------------------
% If you are using Matlab on osf1, the symbolic functions are
% not available. Equivalent capabilities are available from Maple.
% If you are logged on to osf1, type
matlab
% to start Matlab, and
maple
% to start Maple.
% Once inside Maple, you can type
f:=exp(x);
taylor(f,x=0);
% to derive the Taylor series for f(x) = exp(x). To learn more
% about the taylor function, type
?taylor
% and to exit Maple, type
quit;
% Notice that every command ends with ";", and that ":=" is used
% for an assignment statement.
|
