In the last lab, we learned how to start up Matlab and a little bit about using Matlab. I asked you to read a brief excerpt from the Matlab ``Getting Started'' information that is available on the Web from The Mathworks. This lab is intended to reinforce that information for using the Matlab language as a programming language. We will not be concerned with all of the details of this language, we will focus on those aspects of the language that make it ideal for numerical calculations.
But Matlab is a programming language and it is important to learn the basics of good programming so that you will be able to use Matlab for research and applications.
This lab will take two sessions. If you print it, you may find the Adobe pdf version more appropriate.
The best way to use Matlab is to use its scripting (programming) facility. With sequences of Matlab commands contained in files, it is easy to see what calculations were done to produce a certain result, and it is easy to show that the correct values were used in producing a graph. It is terribly embarrassing to produce a very nice plot that you show to your advisor only to discover later that you cannot reproduce it or anything like it for similar conditions or parameters. When the commands are in clear text files, with easily read, wellcommented code, you have a very good idea of how a particular result was obtained. And you will be able to reproduce it and similar calculations as often as you please.
The Matlab comment character is a percent sign (%). That is, lines starting with % are not read as Matlab commands and can contain any text. Similarly, any text on a line following a % can contain textual comments and not Matlab commands. It is important to include comments in script and function files to explain what is being accomplished.
Matlab commands are sometimes terminated with a semicolon (;) and sometimes not. The difference is that the result of a calculation is printed to the screen when there is no semicolon but no printing is done when there is a semicolon. It is a good idea to put semicolons at the ends of all calculational lines in a function file.
There are three kinds of files that Matlab can use:
A Matlab script mfile is a text file with the extension .m. Matlab script files should always start off with comments that identify the author, the date, and a brief description of the intent of the calculation that the file performs. Matlab script files are invoked by typing their names without the .m at the Matlab command line or by using their names inside another Matlab file. Invoking the script causes the commands in the script to be executed, in order.
Matlab function files are also text files with the extension .m, but the first noncomment line must start with the word function and be of the form
function y=sin(x)and the name of the file would be sin.m. The defining line, starting with the word function is called the ``signature'' of the function. If a function has no input parameters, they, and the parentheses, can be omitted. Similarly, a function need not have output parameters. There can be more than one output parameter, and the syntax for several output parameters is discussed later in this lab. The name of the function must be the same as the file name. It is best to have the first line of the function mfile be the signature line, starting with the word ``function.'' The lines following the signature should always contain comments with the following information.
The key differences between function and script files are that
When I am working on a task, I often start out using script files. As I discover just what tasks are repetitive or when I start to need the same calculation repeated for different parameters, or when I have many intermediate variables that might have the same names as variables in other parts of the calculation, I switch to function files. In these labs, I will specify function file or script file when it is important, and you are free to use what you like when I do not specify.
Because function files are intended to be used multiple times, it is a bad idea to have them print or plot things. Imagine what happens if you have a function that prints just one line of information that you think might be useful, and you put it into a loop that is executed a thousand times. Do you plan to read those lines?
Note: Matlab function names are casesensitive. This means that the function cos is different from Cos, coS, and COS. File names in Microsoft operating systems, however, are not strictly case sensitive. To avoid confusion for those students who might be using Matlab on a Microsoft system, in these labs we will use lowercase names for Matlab function and script files.
Matlab also supports data files. The Matlab save command will cause every variable in the workspace to be saved in a file called ``matlab.mat''. You can also name the file with the command save filename that will put everything into a file named ``filename.mat''. This command has many other options, and you can find more about it using the help facility. The inverse of the save command is load.
Matlab uses variable names to represent data. A variable name represents a matrix containing complex doubleprecision data. Of course, if you simply tell Matlab x=1, Matlab will understand that you mean a matrix and it is smart enough to print x out without its decimal and imaginary parts, but make no mistake: they are there. And x can just as easily turn into a matrix.
A variable can represent some important value in a program, or it can represent some sort of dummy or temporary value. Important quantities should be given names longer than a few letters, and the names should indicate the meaning of the quantity. For example, if you were using Matlab to generate a matrix containing a table of squares of numbers, you might name the table tableOfSquares. (The convention I am using here is that the first part of the variable name should be a noun and it should be lower case. Modifying words follow with upper case letters separating the words. This rule comes from the officially recommended naming of Java variables.)
Once you have used a variable name, it is bad practice to reuse it to mean something else. It is sometimes necessary to do so, however, and the statement
clear variableOne variableTwoshould be used to clear the two variables variableOne and variableTwo before they are reused. This same command is critical if you reuse a variable name but intend it to have smaller dimensions.
Matlab has a few reserved variable names. You should not use these variables in your mfiles. If you do use such variables as i or pi, they will lose their special meaning until you clear them. Reserved names include
a1='sqrt(4)' a2=sqrt(4)
We said that Matlab treats all its variables as though they were matrices. Important subclasses of matrices include row vectors (matrices with a single row and possibly several columns) and column vectors (matrices with a single column and possibly several rows). One important thing to remember is that you don't have to declare the size of your variable; Matlab decides how big the variable is when you try to put a value in it. The easiest way to define a row vector is to list its values inside of square brackets, and separated by spaces or commas:
rowVector = [ 0, 1, 3, 6, pi ]The easiest way to define a column vector is to list its values inside of square brackets, separated by semicolons or line breaks.
columnVector1 = [ 0; 1; 3; 6; pi ] columnVector2 = [ 0 1 9 36 100 ](It is not necessary to line the entries up as I have done, but it makes it look nicer.) Note that rowVector is not equal to columnVector1 even though each of their components is the same.
Matrices can be written using both commas and semicolons. The matrix
(1) 
A=[1,2,3;4,5,6;7,8,9]or, more clearly, as
A=[ 1 2 3 4 5 6 7 8 9];
Matlab has a special notation for generating a set of equally spaced values, which can be useful for plotting and other tasks. The format is:
start : increment : finishor
start : finishin which case the increment is understood to be 1. Both of these expressions result in row vectors. So we could define the even values from 10 to 20 by:
evens = 10 : 2 : 20
Sometimes, you'd prefer to specify the number of items in the list, rather than their spacing. In that case, you can use the linspace function, which has the form
linspace( firstValue, lastValue, numberOfValues )in which case we could generate six even numbers with the command:
evens = linspace ( 10, 20, 6 )or fifty evenlyspaced points in the interval [10,20] with
points = linspace ( 10, 20, 50 )As a general rule, use the colon notation when the increment is an integer or when you know what the increment is and use linspace when you know the number of values but not the increment.
Another nice thing about Matlab vector variables is that they are flexible. If you decide you want to add another entry to a vector, it's very easy to do so. To add the value 22 to the end of our evens vector:
evens = [ evens, 22 ]and you could just as easily have inserted a value 8 before the other entries, as well.
Even though the number of elements in a vector can change, Matlab always knows how many there are. You can request this value at any time by using the numel function. For instance,
numel ( evens )should yield the value 7 (the 6 original values of 10, 12, ... 20, plus the value 22 tacked on later). In the case of matrices with more than one nontrivial dimension, the numel function returns the total number of entries. the size function returns a vector containing two values: the number of rows and the number of columns. To get the number of rows of a variable v, use size(v,1) and to get the number of columns use size(v,2). For example, since evens is a row vector, size( evens, 1)=1 and size( evens, 2)=7, one row and seven columns.
To specify an individual entry of a vector, you need to use index notation, which uses round parentheses enclosing the index of an entry. The first element of an array has index 1 (as in Fortran, but not C or Java). Thus, if I want to alter the third element of evens, I could say
evens(3) = 7There is special syntax for the final index value. Since v is a vector with length 7, you can write v(end) to mean the same thing as v(7).
plot(meshPoints,sin(2*pi*meshPoints))Please save this plot as a jpeg (.jpg) file and send it along with your summary. You can use the ``FileSave as'' menu choice from the plot window. Use file type ``JPEG image.'' The command line expression to save a plot as JPEG is
print djpeg plotname.jpgwhere ``plotname.jpg'' is a name you choose for the file.
% Lab 2, exercise 2 % A sample script file. % Your name and the dateFollow the header comments with the commands containing exactly the commands you used in the earlier parts of this exercise. Test your script by using clear to clear your results and then execute the script from the command line by typing exer2.
Matlab provides a large assembly of tools for matrix and vector manipulation. We will investigate a few of these by trying them out.
rowVec1 = [ 1 4 9] colVec1 = [ 2 9 8 ] mat1 = [ 1 3 5 7 9 0 2 4 6 ]
colVec2 = (pi/4) * colVec1
colVec2 = cos( colVec2 )Note that the values of colVec2 have been overwritten. Are these the values you expect?
colVec3 = colVec1 + colVec2
illegal = colVec1 + rowVec1;Look carefully at the error message. You must recognize from the message what went wrong when you see it in the future.
norm(colVec3)
colvec4 = mat1 * colVec1
mat1Transpose = mat1' rowVec2 = colVec3'Warning: The single quote really means the complexconjugate transpose (or Hermitian adjoint). If you want a true transpose applied to a complex matrix you must use ``
.'
''.
mat2 = mat1 * mat1' % symmetric matrix rowVec3 = rowVec1 * mat1 dotProduct = colVec3' * colVec1 euclideanNorm = sqrt(colVec2' * colVec2)
determinant = det( mat1 ) tr = trace( mat1 )
min(rowVec1)
max(mat1)
max(abs(rowVec1))
A=magic(100); % please do not print all 10,000 entries.The matrix A has 100 row sums (one for each row), 100 column sums (one for each column) and two diagonal sums. These 202 sums should all be exactly the same, and you could verify that they are the same by printing them and ``seeing'' that they are the same. It is easy to miss small differences among so many numbers, though. Instead, verify that A is a magic square by constructing the 100 column sums (without printing them) and computing the maximum and minimum values of the column sums. Do the same for the 100 row sums, and compute the two diagonal sums. Check that these six values are the same. If the maximum and minimum values are the same, the flyswatter principle says that all values are the same.
integers = 0 : 10but now we'll get an error when we try to multiply the entries of integers by themselves.
squareIntegers = integers * integers
Realize that Matlab deals with vectors, and the default multiplication operation with vectors is rowbycolumn multiplication. What we want here is elementbyelement multiplication, so we need to place a period in front of the operator:
squareIntegers = integers .* integers
Now we can define cubeIntegers and fourthIntegers in a similar way.
cubeIntegers = squareIntegers .* integers fourthIntegers = squareIntegers .* squareIntegers
Finally, we would like to print them out as a table. integers, squareIntegers, etc. are row vectors, so make a matrix whose columns consist of these vectors and allow Matlab to print out the whole matrix at once.
tableOfPowers=[integers', squareIntegers', cubeIntegers', fourthIntegers']
sqIntegers = integers .^ 2and check that the two calculations agree with the command
norm(sqIntegerssquareIntegers)which should result in zero.
Remark: Addition, subtraction, and division or multiplication by a scalar never require the dot in front of the operator, although you will get the correct result if you use one.
tableOfCubes = tableOfPowers(:,[1,3]) tableOfEvenCubes = tableOfPowers(1:2:end,[1,3]) tableOfOddFourths = tableOfPowers(2:2:end,1:3:4)Note: [1:3:4] is the same as [1,4].
A = magic(10)What commands would be needed to generate the four matrices in the upper left quarter AUL, the upper right quarter AUR, the lower left quarter ALL, and the lower right quarter ALR.
B=[AUL AUR ALL ALR];and show that A and B are the same by showing that norm(AB) is zero.
It is critical to be able to ask questions and to perform repetitive calculations in mfiles. These topics are examples of ``flow control'' constructs in programming languages. Matlab provides two basic looping (repetition) constructs: for and while, and the if construct for asking questions. These statements each surround several Matlab statements with for, while or if at the top and end at the bottom.
Note: It is an excellent idea to indent the statements between
the for, while, or if lines and the end line.
This indentation strategy makes code immensely more readable.
Your mfiles will be expected to follow this convention.
Syntax  Example  

for loop 



while loop 



a simple if 



a compound if 


Note that elseif is one word! Using two words else if changes the statement into two nested if statements with possibly a very different meaning, and a different number of end statements.
The following Matlab code computes the trapezoid rule for the case that , and .
% Use the trapezoid rule to approximate the integral from 0 to 1 % of exp(x), using N (=100) intervals % Your name and the date N=40; h=1/(N1); x=h; % look at this trick approxIntegral=0.; for k=1:N % compute current x value x=x+h; % add the terms up if k==1  k==N approxIntegral=approxIntegral+(h/2)*exp(x); % ends of interval else approxIntegral=approxIntegral+h*exp(x); % middle of interval end end
approxIntegral=approxIntegral+(h/2)*exp(x); % ends of interval
approxIntegral=approxIntegral+h*exp(x); % middle of interval
If you have to type everything at the command line, you will not get very far. You need some sort of scripting capability to save the trouble of typing, to make editing easier, and to provide a record of what you have done. You also need the capability of making functions or your scripts will become too long to understand. In this section we will consider first a script file and later a function mfile. We will be using graphics in the script file, so you can pick up how graphics can be used in our work.
The Fourier series for the function on the interval
is
% compute NTERMS terms of the Fourier Series for y=x^3 % plot the result using NPOINTS points from 1 to 1. % Your name and the date NTERMS=10; NPOINTS=1000; x=linspace(1,1,NPOINTS); y=zeros(size(x)); for k=1:NTERMS term=2*(1)^k*(pi^2/k6/k^3)*sin(k*x); y=yterm end plot(x,y); hold on plot(x,x.^3,'g'); % 'g' is for a green line axis([1,1,2,2]); hold offIt is always good programming practice to define constants symbolically at the beginning of a program and then to use the symbols within the program. Sometimes these special constants are called ``magic numbers.'' By convention, symbolic constants are named using all upper case.
y=zeros(size(x));do? Hint: You can use Matlab help for zeros and size for more information.
Note: The command hold, without the ``on'' or ``off'', is a ``toggle.'' Each time it is used, it switches from ``on'' to ``off'' or from ``off'' to ``on.'' Using it that way is easier but you have to remember which state you are in.
In the following exercise you are going to modify the calculation so that it continues to add terms so long as the largest component of the next term remains larger in absolute value than, say, 0.05. Since the series is of alternating sign, this quantity is a legitimate estimate of the difference between the partial sum and the limit. The while statement is designed for this case. It would be used in the following way, replacing the for statement.
TOLERANCE=0.05; % the chosen tolerance value <<some lines of code from above>> k=0; term=TOLERANCE+1; % bigger than TOLERANCE while max(abs(term)) > TOLERANCE k = k + 1; <<some lines of code from above>> end disp( strcat('Number of iterations =',num2str(k)) ) <<some lines of code to plot results>>
term=TOLERANCE+1;
k = k + 1;
Hint: If you try this script and it does not quit (it stays ``busy'') you can interrupt the calculation by holding the Control key (``CTRL'') and pressing ``C''.
If you find that your code does not quit, you have a bug. For some reason, the value of the variable difference is not getting small. Look at your code carefully to find out why. If you cannot see it, use the debugger to watch execution and see what is happening to the value of difference.
The next task is to make a function mfile out of exer6.m in such a way that it takes as argument the desired tolerance and returns the number of iterations required to meet that tolerance.
function k = exer7( tolerance ) % k = exer7( tolerance ) % more comments % Your name and the date
exer7(0.05)will cause Matlab to print the result. Using the function on the right side of an equal sign with a variable on the left side of the equal sign causes the variable to take the value given by the function.
numItsRequired=exer7(0.05)causes both a printed value and assignment of a value to numItsRequired.
Unlike ordinary mathematical notation, Matlab allows a function to return two or more values instead of a single value. The syntax to accomplish this trick is similar to that of defining a vector, but the meaning is nothing like a vector. For a function named ``fnct'' that returns two variables depending on a single variable as input, the signature line would look like:
[y,z] = fnct( x )and to use the function, for the value at x = 3,
[y,z] = fnct( 3 )If you have the same function but wish only the first output variable, y, you would write
y = fnct( 3 )and if you wish only the second output variable, z, you would write
[~,z] = fnct( 3 )
y8=exer8(0.02); y9=exer9(0.02, @sin); % following difference should be small norm(y8y9)Note: The
@
symbol is used to identify sin as
a function instead of an ordinary variable. If you forget the @
,
you will receive a very strange error message.
You have now completed all the required work for this lab. There is an extra credit exercise below. If you wish to gain the extra credit, do that exercise now.
When you have finished your work, be sure to send it to me. It should include your summary file, the diary.txt file, each of the mfiles exer2.m, exer5.m, exer6.m, exer7.m, exer8.m and the plot files (.jpg) files you created. Also include your extra credit work if you chose to do it. Remember that you can zip all these files in your area from the Matlab command line with the command
!zip filename.zip *.m *.txt *.jpgwhere filename.zip is a file name of your choosing. Do not forget to close your diary file first by using the command ``diary off''. (If you are using your own computer, you can zip the files using any method you are familiar with.)
There is an additional command that is used for displaying results and that allows you a great deal of formatting flexibility. It is similar to the ``printf'' command in C and Java.
for i=0:12 disp(['The square root of ',num2str(i^2),' is ',num2str(i)]); end
Back to MATH2070 page.