### Transcription

Dynamics and VibrationsMATLAB tutorialSchool of EngineeringBrown UniversityThis tutorial is intended to provide a crash-course on using a small subset of the features of MATLAB. Ifyou complete the whole of this tutorial, you will be able to use MATLAB to integrate equations of motionfor dynamical systems, plot the results, and use MATLAB optimizers and solvers to make designdecisions.You can work step-by-step through this tutorial, or if you prefer, you can brush up on topics from the listbelow.If you are working through the tutorial for the first time, you should complete sections 1-13.6. You cando the other sections later, when they are needed. If you have taken EN30, you should be familiar withthe material in Sections 1-10: you can review these if you like, or skip directly to Section 11.1. What is MATLAB2. Starting MATLAB3. Basic MATLAB windows4. Simple calculations using MATLAB5. MATLAB help6. Errors associated with floating point arithmetic (and an example of a basic loop)7. Vectors in MATLAB7.1 Creating Vectors with a Loop7.2 Dot notation for operations on vectors7.3 Operations on vectors using a loop8. Matrices in MATLAB9. Plotting and graphics in MATLAB10. Working with M-files11. MATLAB Functions12. Organizing complex calculations as functions in an M-file13. Solving ordinary differential equations (ODEs) using MATLAB13.1 What is a differential equation?13.2 Solving a basic differential equation13.3 Solving a basic differential equation in an M-file13.4 Solving a differential equation with adjustable parameters13.5 Common errors13.6 Solving simultaneous differential equations13.7 Controlling the accuracy of solutions to differential equations13.8 Looking for special events in a solution13.9 How the ODE solver works13.10 Other MATLAB differential equation solvers14. Using MATLAB solvers and optimizers to make design decisions14.1 Using fzero to solve equations14.2 Simple unconstrained optimization problem14.3 Optimizing with constraints15. Reading and writing data to/from files16. Movies and animation

Select a convenient directory where you will be able to save your files.4. Simple calculations using MATLABYou can use MATLAB as a calculator. Try this for yourself, by typing the following into the commandwindow. Press ‘enter’ at the end of each line. x 4 y x 2 z factorial(y) w log(z)*1.e-05 format long w format long eng w format short w sin(pi)MATLAB will display the solution to each step of the calculation just below the command. Do you noticeanything strange about the solution given to the last step?Once you have assigned a value to a variable, MATLAB remembers it forever. To remove a value from avariable you can use the ‘clear’ statement - try clear a aIf you type ‘clear’ and omit the variable, then everything gets cleared. Don’t do that now – but it isuseful when you want to start a fresh calculation.MATLAB can handle complex numbers. Try the following z x i*y real(z) imag(z) conj(z) angle(z) abs(z)You can even do things like log(z) sqrt(-1) log(-1)Notice that: Unlike MAPLE, Java, or C, you don’t need to type a semicolon at the end of the line (To properlyexpress your feelings about this, type load handel and then sound(y,Fs) in the commandwindow). If you do put a semicolon, the operation will be completed but MATLAB will not print the result.This can be useful when you want to do a sequence of calculations.

Special numbers, like pi’ and ‘i’ don’t need to be capitalized. But beware – you often use i as acounter in loops – and then the complex number i gets re-assigned as a number. You can also dodumb things like pi 3.2 (You may know that in 1897 a bill was submitted to the Indianalegislature to declare pi 3.2 but fortunately the bill did not pass). You can reset these specialvariables to their proper definitions by using clear i or clear piThe Command History window keeps track of everything you have typed. You can double leftclick on a line in the Command history window to repeat it, or right click it to see a list of otheroptions.The output in the command window looks like crap. MATLAB is not really supposed to be usedlike this. We will discuss a better approach later.If you screw up early on in a sequence of calculations, there is no quick way to fix your error,other than to type in the sequence of commands again. You can use the ‘up arrow’ key to scrollback through a sequence of commands. Again, there is a better way to use MATLAB that getsaround this problem.If you are really embarrassed by what you typed, you can right click the command window anddelete everything (but this will not reset variables). You can also delete lines from the Commandhistory, by right clicking the line and selecting Delete Selection. Or you can delete the entireCommand History.You can get help on MATLAB functions by highlighting the function, then right clicking the lineand selecting Help on Selection. Try this for the sqrt(-1) line.5. MATLAB helpHelp is available through the online manual – Click on the question-mark in the strip near the top right ofthe window).By default the help window opens inside theMATLAB GUI, but you can drag it out so itoccupies a new window on your desktop.If you already know the name of the MATLABfunction you want to use the help manual is quitegood – you can just enter the name of the functionin the search, and a page with a good number ofexamples usually comes up.It is morechallenging to find out how to do something, butmost of the functions you need can be found byclicking on the MATLAB link on the main pageand then following the menus that come up onsubsequent pages.

6. Errors associated with floating point arithmetic (and a basic loop)If you have not already done so, use MATLAB to calculate sin(pi)The answer, of course, should be zero, but MATLAB returns a small, but finite, number. This is becauseMATLAB (and any other program) stores floating point numbers as sequences of binary digits with afinite length. Obviously, it is impossible to store the exact value of π in this way.More surprisingly, perhaps, it is not possible even to store a simple decimal number like 0.1 as a finitenumber of binary digits. Try typing the following simple MATLAB program into the command window a 0; for n 1:10 a a 0.1; end a a – 1Here, the line that reads “for n 1:10 a a 0.1; end” is called a “loop.” This is a very common operationin most computer programs. It can be interpreted as the command: “for each of the discrete values of theinteger variable n between 1 and 10 (inclusive), calculate the variable “a” by adding 0.1 to the previousvalue of “a” The loop starts with the value n 1 and ends with n 10. Loops can be used to repeatcalculations many times – we will see lots more examples in later parts of the tutorial.Thus, the for end loop therefore adds 0.1 to the variable a ten times.approximately 1. But when you compute a-1, you don’t end up with zero.It gives an answer that isOf course -1.1102e-016 is not a big error compared to 1, and this kind of accuracy is good enough forgovernment work. But if someone subtracted 1.1102e-016 from your bank account every time afinancial transaction occurred around the world, you would burn up your bank account pretty fast.Perhaps even faster than you do by paying your tuition bills.You can minimize errors caused by floating point arithmetic by careful programming, but you can’teliminate them altogether. As a user of MATLAB they are mostly out of your control, but you need toknow that they exist, and try to check the accuracy of your computations as carefully as possible.7. Vectors in MATLABMATLAB can do all vector operations completely painlessly. Try the following commands a [6,3,4] a(1) a(2) a(3) b [3,1,-6] c a b c dot(a,b) c cross(a,b)Calculate the magnitude of c (you should be able to do this with a dot product. MATLAB also has abuilt-in function called norm’ that calculates the magnitude of a vector)A vector in MATLAB need not be three dimensional. For example, try

a [9,8,7,6,5,4,3,2,1] b [1,2,3,4,5,6,7,8,9]You can add, subtract, and evaluate the dot product of vectors that are not 3D, but you can’t take a crossproduct. Try the following a b dot(a,b) cross(a,b)In MATLAB, vectors can be stored as either a row of numbers, or a column of numbers. So you couldalso enter the vector a as a [9;8;7;6;5;4;3;2;1]to produce a column vector.You can turn a row vector into a column vector, and vice-versa by b transpose(b)7.1 Creating vectors with a loopYou can create a vector containing regularly spaced data points very quickly with a loop. Try for i 1:11v(i) 0.1*(i-1);end vThe for end loop repeats the calculation with each value of i from 1 to 11. Here, the “counter”variable i now is used to refer to the ith entry in the vector v, and also is used in the formula itself.As another example, suppose you want to create a vector v of 101 equally spaced points, starting at 3 andending at 2*pi, you would use for i 1:101v(i) 3 (2*pi-3)*(i-1)/100;end v7.2 Dot notation for operations on vectorsYou can often manipulate the numbers in a vector using simple MATLAB commands. For example, ifyou type sin(v)MATLAB will compute the sin of every number stored in the vector v and return the result as anothervector. This is useful for plots – see section 11.You have to be careful to distinguish between operations on a vector (or matrix, see later) and operationson the components of the vector. For example, try typing v 2This will cause MATLAB to bomb, because the square of a row vector is not defined. But you can type v. 2(there is a period after the v, and no space). This squares every element within v.You can also do things like v. /(1 v)

(see if you can work out what this does).7.3 Operations on vectors using a loopI personally avoid using the dot notation – mostly because it makes code hard to read. Instead, I generallydo operations on vector elements using loops. For example, instead of writing w v. 2, I would use for i 1:length(v)w(i) v(i) 2;end wHere, ‘for i 1:length(v)’ repeats the calculation for every element in the vector v. The functionlength(vector) determines how many components the vector v has (101 in this case).Using loops is not elegant programming, and slows down MATLAB. Purists (like CS40 TAs) object to it.But I don’t care. For any seriously computationally intensive calculations I would use a programminglanguage like C or Fortran95, not MATLAB. Programming languages like Java, C , or Python arebetter if you need to work with complicated data structures.8. Matrices in MATLABHopefully you know what a matrix is If not, it doesn’t matter - for now, it is enough to know that amatrix is a set of numbers, arranged in rows and columns, as shown belowColumn 2 1 5 3 9543206082 6 5 7 row 1row 4Column 3A matrix need not necessarily have the same numbers of rows as columns, but most of the matrices wewill encounter in this course do. A matrix of this kind is called square. (My generation used to call ourprofessors and parents square too, but with hindsight it is hard to see why. ‘Lumpy’ would have beenmore accurate).You can create a matrix in MATLAB by entering the numbers one row at a time, separated bysemicolons, as follows A [1,5,0,2; 5,4,6,6;3,3,0,5;9,2,8,7]You can extract the numbers from the matrix using the convention A(row #, col #). Try A(1,3) A(3,1)You can also assign values of individual array elements A(1,3) 1000

There are some short-cuts for creating special matrices. Try the following B ones(1,4) C pascal(6) D eye(4,4) E zeros(3,3)The ‘eye’ command creates the ‘identity matrix’ – this is the matrix version of the number 1. You canuse help pascalto find out what pascal does.MATLAB can help you do all sorts of things to matrices, if you are the sort of person that enjoys doingthings to matrices. For example1. You can flip rows and columns with B transpose(A)2. You can add matrices (provided they have the same number of rows and columns C A BTry also C – transpose(C)A matrix that is equal to its transpose is called symmetric3. You can multiply matrices – this is a rather complicated operation, which will be discussed inmore detail elsewhere. But in MATLAB you need only to type D A*B to find the product ofA and B. Also try the following E A*B-B*A F eye(4,4)*A - A4. You can do titillating things like calculate the determinant of a matrix; the inverse of a matrix, theeigenvalues and eigenvectors of a matrix. If you want to try these things det(A) inv(A) eig(A) [W,D] eig(A)You can find out more about these functions, and also get a full list of MATLAB matrixoperations in Linear Algebra page of the help menu.

MATLAB can also calculate the product of a matrix and a vector. This operation is used very frequentlyin engineering calculations. For example, you can multiply a 3D column vector by a matrix with 3 rowsand 3 columns, as follows v [4;3;6] A [3,1,9;2,10,4;6,8,2] w A*vThe result is a 3D column vector. Notice that you can’t multiply a 3D row vector by a 3x3 matrix. Trythis v [4,3,6] w A*vIf you accidentally enter a vector as a row, you can convert it into a column vector by using v transpose(v)MATLAB is also very good at solving systems of linear equations. For example, consider the equations3 x1 4 x2 7 x3 65 x1 2 x2 9 x3 1 x1 13 x2 3 x3 8This system of equations can be expressed in matrix form asAx b x1 3 4 7 6 A 5 2 9 x x2 b 1 x3 1 13 3 8 To solve these in MATLAB, you would simply type A [3,4,7;5,2,-9;-1,13,3] b [6;1;8] x A\b(note the forward-slash, not the back-slash or divide sign) You can check your answer by calculating A*xThe notation here is supposed to tell you that x is b ‘divided’ by A – although division’ by a matrix hasto be interpreted rather carefully. Try also x transpose(b)/AThe notation transpose(b)/A solves the equations xA b , where x and b are row vectors. Again, you cancheck this with x*A(The answer should equal b, (as a row vector) of course)MATLAB can quickly solve huge systems of equations, which makes it useful for many engineeringcalculations. The feature has to be used carefully because systems of equations may not have a solution,or may have many solutions – MATLAB has procedures for dealing with both these situations but if youdon’t know what it’s doing you can get yourself into trouble. For more info on linear equations check thesection of the manual below

9. Plotting and graphics in MATLAB This is an important section of this tutorialPlotting data in MATLAB is very simple. Try the following for i 1:101 x(i) 2*pi*(i-1)/100; end y sin(x) plot(x,y)MATLAB should produce something that lookslike this (the 2014b MATLAB version producescurves with a different color scheme and linethickness to 2012b)MATLAB lets you edit and annotate a graphdirectly from the window. For example, youcan go to Tools Edit Plot, then double-clickthe plot. A menu should open up that will allowyou to add x and y axis labels, change the rangeof the x-y axes; add a title to the plot, and so on.You can change axis label fonts, the linethickness and color, the background, and so on – usually by double-clicking what you want to change andthen using the pop-up editor. You can export figures in many different formats from the File Save Asmenu – and you can also print your plot directly. Play with the figure for yourself to see what you can do.

Helpful Hint: If you annotate a plot by hand, youcan make MATLAB generate the code to re-createyour annotated plot automatically by going to File Generate Code. This is useful when you write yourown scripts and want to label plots automatically.To plot multiple lines on the same plot you can use y sin(x) plot(x,y) hold on y sin(2*x) plot(x,y)Alternatively, you can use y(1,:) sin(x); y(2,:) sin(2*x); y(3,:) sin(3*x); plot(x,y)Here, y is a matrix. The notation y(1,:) fills the first row of y, y(2,:) fills the second, and so on. Thecolon : ensures that the number of columns is equal to the number of terms in the vector x. If you prefer,you could accomplish the same calculation in a loop: for i 1:length(x) y(1,i) sin(x(i)); y(2,i) sin(2*x(i)); y(3,i) sin(3*x(i)); end plot(x,y)Notice that when you make a new plot, it always wipes out the old one. If you want to create a new plotwithout over-writing old ones, you can use figure plot(x,y)The ‘figure’ command will open a new window and will assign a new number to it (in this case, figure 2).If you want to go back and re-plot an earlier figure you can use figure(1) plot(x,y)If you like, you can display multiple plots in the same figure, by typing newaxes axes; plot(x,y)The new plot appears over the top of the old one, but you can drag it away by clicking on the arrow tooland then clicking on any axis or border of new plot. You can also re-size the plots in the figure windowto display them side by side. The statement ‘newaxes axes’ gives a name (or ‘handle’) to the new axes,so you can select them again later. For example, if you were to create a third set of axes yetmoreaxes axes; plot(x,y)you can then go back and re-plot newaxes’ by typing axes(newaxes); plot(x,y)Doing parametric plots is easy. For example, try for i 1:101 t(i) 2*pi*(i-1)/100; end x sin(t); y cos(t); plot(x,y)

MATLAB has vast numbers of different 2D and 3D plots. For example, to draw a filled contour plot ofthe function z sin(2π x)sin(2π y ) for 0 x 1, 0 y 1 , you can use for i 1:51 x(i) (i-1)/50; y(i) x(i); end z transpose(sin(2*pi*y))*sin(2*pi*x); figure contourf(x,y,z)The first two lines of this sequence should be familiar: they create row vectors of equally spaced points.The third needs some explanation – this operation constructs a matrix z, whose rows and columns satisfyz(i,j) sin(2*pi*y(i)) * sin(2*pi*x(j)) for each value of x(i) and y(j).If you like, you can change the number of contour levels contourf(x,y,z,15)You can also plot this data as a 3D surface using surface(x,y,z)The result will look a bit strange, but you can click on the ‘rotation 3D’ button (the little box with acircular arrow around it ) near the top of the figure window, and then rotate the view in the figure withyour mouse to make it look more sensible.You can find out more about the different kinds of MATLAB plots in the section of the manual shownbelow

10. Working with M filesSo far, we’ve run MATLAB by typing into the command window. This is OK for simple calculations,but it is usually better to type the list of commands you plan to execute into a file (called an M-file), andthen tell MATLAB to run all the commands together. One benefit of doing this is that you can fix errorseasily by editing and re-evaluating the file. But more importantly, you can use M-files to write simpleprograms and functions using MATLAB.To create an M-file, simply press the ‘New Script’ button on the main MATLAB window.This will open a new window for the matlab script editor, as shown belowNow, type the following lines into the editor (see the screenshot):for i 1:101theta(i) -pi 2*pi*(i-1)/100;rho(i) 2*sin(5*theta(i));endfigurepolar(theta,rho)You can make MATLAB execute these statements by:1. Pressing the green arrow labeled ‘Run’ near the top of the editor window – this will first save thefile (MATLAB will prompt you for a file name – save the script in a file called myscript.m), andwill then execute the file.

2. You can save the file yourself (e.g. in a file called myscript.m). You can then run the script fromthe command window, by typing myscriptYou don’t have to use the MATLAB editor to create M-files – any text editor will work. As long as yousave the file in plain text (ascii) format in a file with a .m extension, MATLAB will be able to find andexecute the file. To do this, you must open the directory containing the file with MATLAB (e.g. byentering the path in the field at the top of the window). Then, typing filenamein the command window will execute the file. Usually it’s more convenient to use the MATLAB editor but if you happen to be stuck somewhere without access to MATLAB this is a useful alternative. (Thenagain, this means you can’t use lack of access to MATLAB to avoid doing homework, so maybe there is adown-side)11. MATLAB functionsThis is an important section of this tutorialYou can also use M-files to define new MATLAB functions – these are programs that can accept userdefined data and use it to produce a new result. For example, to create a function that computes themagnitude of a vector:1. Open a new M-file with the editor (press the ‘New’ button on the editor)2. Type the following into the M-filefunction y magnitude(v)% Function to compute the magnitude of a vectory sqrt(dot(v,v));endNote that MATLAB ignores any lines that start with a % - this is to allow you to type commentsinto your programs that will help you, or users of your programs, to understand them.3. Save the file (accept the default file name, which is magnitude.m)4. Type the following into the command window x [1,2,3]; magnitude(x)Note the syntax for defining a function – it must always have the formfunction solution variable function name(input variables ) script that computes the function, ending with the command:endsolution variable .You can visualize the way a function worksusing the picture.Think of the ‘inputvariables’ as grass (the stuff cows eat, that is,not the other kind of grass); the body of thescript as the cow’s many stomachs, and thesolution variable as what comes out of theother end of the cow.jSolutionvariable(s)Computationsinside functionInputvariable(s)A cow can eat many different vegetables, and you can put what comes out of the cow in many differentcontainers. A function is the same: you can pass different variables or numbers into a function, and you

can assign the output to any variable you like – try result1 magnitude([3,4])magnitude([2,4,6,8]) for example.or result2 If you are writing functions for other people to use, it is a good idea to provide some help for the function,and some error checking to stop users doing stupid things. MATLAB treats the comment lines thatappear just under the function name in a function as help, so if you type help magnitudeinto the command window, MATLAB will tell you what magnitude’ does. You could also edit theprogram a bit to make sure that v is a vector, as followsfunction y magnitude(v)% Function to compute the magnitude of a vector% Check that v is a vector[m,n] size(v);% [m,n] are the number of rows and columns of v% The next line checks if both m and n are 1, in which case v is% a matrix, or if both m 1 and n 1, in which case v scalar; if either,% an error is returnedif ( ( (m 1) & (n 1) ) (m 1 & n 1) )error('Input must be a vector')endy sqrt(dot(v,v));endThis program shows another important programming concept. The “if end” set of commands iscalled a “conditional” statement. You use it like this:if (a test)calculation to do if the test is trueendIn this case the calculation is only done if test is true – otherwise the calculation is skipped altogether.You can also doif (a test)calculation to do if the test is trueelse if (a second test)calculation to do if the first test is false and the second test is trueelsecalculation to do if both tests are falseendIn the example (m 1) & (n 1) means “m 1 and n 1” while (m 1 & n 1) means“or m 1 and n 1”. The symbols & and are shorthand for and’ and ‘or’.Functions can accept or return data as numbers, arrays, strings, or even more abstract data types. Forexample, the program below is used by many engineers as a substitute for social skillsfunction s pickup(person)% function to say hello to someone cutes ['Hello ' person ' you are cute']beep;end

(If you try to cut and paste this function into MATLAB the quotation marks will probably come outwrong in the Matlab file and the function won’t work. The quotes are all ‘close single quote’ on yourkeyboard – you will need to correct them by hand.)You can try this out with pickup(‘Janet’)(Janet happens to be my wife’s name. If you are also called Janet please don’t attach any significance tothis example)You can also pass functions into other functions. For example, create an M-file with this programfunction plotit(func to plot,xmin,xmax,npoints)% plot a graph of a function of a single variable f(x)% from xmin to xmax, with npoints data pointsfor n 1:npointsx(n) xmin (xmax-xmin)*(n-1)/(npoints-1);v(n) func to plot(x(n));endfigure;plot(x,v);endThen save the M-file, and try plotting a cosine, by using plotit(@cos,0,4*pi,100)Several new concepts have been introduced here – firstly, notice that variables (like func to plot in thisexample) don’t necessarily have to contain numbers, or even strings. They can be more complicatedthings called “objects.” We won’t discuss objects and object oriented programming here, but the exampleshows that a variable can contain a function. To see this, try the following v @exp v(1)This assigns the exponential function to a variable called v – which then behaves just like the exponentialfunction. The ‘@’ before the exp is called a ‘function handle’ – it tells MATLAB that the functionexp(x), should be assigned to v, instead of just a variable called exp.The ‘@’ in front of cos “plotit(@cos,0,4*pi,100)” assigns the cosine function to the variable calledfunc to plot inside the function.Although MATLAB is not really intended to be a programming language, you can use it to write somequite complicated code, and it is widely used by engineers for this purpose. CS4 will give you a goodintroduction to MATLAB programming. But if you want to write a real program you should use agenuine programming language, like C, C , Java, lisp, etc.IMPORTANT: You may have got used to writing scripts in EN30 without defining them as functions.For the calculations you will be doing in EN40, it will be essential to make all your scripts define afunction. This means you must start the script with ‘function name of function’ and terminate it with an‘end’ statement. An example is shown below.function My Design Project%UNTITLED3 Summary of this function goes here%Detailed explanation goes hereend

If you like, you can automatically create a new function by using the New Function menu, (You willusually be able to delete the input and output arguments for the function).12. Organizing complicated calculations as functions in an M fileIt can be very helpful to divide a long and complicated calculation into a sequence of smaller ones, whichare done by separate functions in an M file. As a very simple example, the M file shown below createsplots a pseudo-random function of time, then computes the root-mean-square of the signal. Copy andpaste the example into a new MATLAB M-file, and press the green arrow to see what it does.function randomsignal% Function to plot a random signal and compute its RMS.close all % This closes any open figure windowsnpoints 100; % No. points in the plotdtime 0.01; % Time interval between points% Compute vectors of time and the value of the function.% This example shows how a function can calculate several% things at the same time.[time,function value] create random function(npoints,dtime);% Compute the rms value of the functionrms compute rms(time,function value);% Plot the functionplot(time,function value);% Write the rms value as a label on the plotlabel strcat('rms signal el},'FontSize',16,.'BackgroundColor',[1 1 1],.'Position',[0.3499 0.6924 0.3944 0.1]);end%function [t,y] create random function(npoints,time interval)% Create a vector of equally spaced times t(i), and% a vector y(i), of random values with i 1 npointsfor i 1:npointst(i) time interval*(i-1);% The rand function computes a random number between 0 and 1y(i) rand-0.5;endendfunction r compute rms(t,y)

%%%%Function to compute

MATLAB tutorial . School of Engineering . Brown University . This tutorial is intended to provide a crash-course on using a small subset of the features of MATLAB. If you complete the whole of this tutorial, you will be able to use MATLAB