Matlab Tutorial - University of Calgary in Alberta



Matlab Tutorial

1 - Introduction

This tutorial serves to provide a basic introduction to Matlab. Students taking this course should have some exposure to Matlab, or at least some other Computer Algebra System (CAS) such as Mathematica or Maple. Support or translation for other packages will not be provided.

When studying the examples in this tutorial, feel free to try different values or cases instead of entering code verbatim. Some concepts are continuously revisited. Other important commands are only briefly introduced. Hopefully you will develop enough familiarity with Matlab to be able to solve any problems that arise by yourself.

Many of the examples in this tutorial have been designed to be relevant to the course material. You may not have been exposed to some of the concepts used. If you are unsure of significance or use of some of the topics, just enter the code verbatim and check to see if the results match those that are provided. This is a Matlab tutorial for controls, not a controls tutorial for Matlab.

2 - Start-up and Help

The Matlab window is shown below. Important features include the main menu, the command window, the current directory and the workspace. Different windows can be enabled through View. Matlab commands are entered in the Command Window. The current directory is the directory in which your files will be saved and accessed. It can be set to any directory. The workspace shows any variables that are initialized.

[pic]

The first thing you should do when learning Matlab is become familiar with the 'help' function. Matlab has thousands of functions, each with its own syntax and options. There is no way to remember every function you might ever need. To find out how to use the on-line help, enter 'help help':

>> help help

HELP On-line help, display text at command line.

HELP, by itself, lists all primary help topics. Each primary topic

…..

Help can be used by itself to show the main topics. The contents of each of these topics can be seen with 'help TOPIC'. This will list all of the functions within that directory. 'help FUNCTION' will display information about the specific function. This information will include the full name of the function, an explanation of how it works, the basic syntax and the syntax to use various options. Some things you may want to try are 'help', 'help syntax', 'help ops', 'help elfun', 'help elmat' and 'help control'.

Another tool that is useful for finding commands is ‘lookfor’

>> lookfor sqrt

SQRT Square root.

SQRTM Matrix square root.

REALSQRT Real square root.

vsqrtm.m: % function out = vsqrtm(mat)

BLKSGNSQRT This block defines a function that returns the ...

BLKSQRT This block defines a function that returns the . . .

‘lookfor STRING’ searches for STRING in the first line of the help comments for each function and command in Matlab. Read ‘help lookfor’.

There are also many demos that illustrate features of Matlab and Simulink. The list of demos can be seen by entering 'help demos' or simply 'demos'. The demos usually run in a graphical user interface. Try 'help demos', 'help ctrldemos' and 'intro'.

3 - Command Line and Built-in Functions

Matlab allows you to both enter commands in a command line interface and write scripts to carry out operations. Basic arithmetic operations (+, -, *, /, ^EXP) can be performed such as:

>> 1+2+3+4+5+6+7+8

ans =

36

The answer has been stored in a variable ‘ans’. As with C, the ‘=’ character is for assignment, not for equal. If another operation is performed without being assigned to a specific variable, ans will be overwritten. Variables in Matlab are case sensitive. Using the character ‘;’ at the end of a command will suppress the output.

>> ANS = 1*2*3*4*5*6*7*8

ANS =

40320

>> who

Your variables are:

ANS ans

Useful commands with their basic syntax and use include:

|Name |Syntax |Use |

|who |who |lists all variables in the workspace |

|whos |whos |lists all variables, long form |

|clear |clear |clears variables from the workspace |

| |clear VAR |clears VAR from the workspace |

|dir |dir |lists files in the current directory |

| |dir DIRECTORY |lists files in DIRECTORY |

|cd |cd |Prints current directory |

| |cd .. |moves one directory up |

| |cd DIRECTORY |moves to DIRECTORY |

|save |save |saves all variables to MATLAB.MAT |

| |save FILE |saves all variables to FILE.MAT |

| |save FILE.EXT |saves all variables to FILE.EXT |

|load |load |loads variables from MATLAB.MAT |

| |load FILE |loads variables from FILE.MAT |

| |load FILE.EXT |loads variables from FILE.EXT |

|delete |delete FILE.EXT |deletes FILE.EXT from current directory |

|clc |clc |clears the command window |

4 - Vectors and Matrices

Brackets, ‘[]’, are used to form vectors and matrices. A row vector can be entered in Matlab as:

>> lint = [1 2 3 4 5 6 7 8]

lint =

1 2 3 4 5 6 7 8

>> sum(lint)

ans =

36

To place elements in another row to make either column vectors or matrices one uses the semi-colon ‘;’. A colon ‘:’ can be used to select elements of a vector or matrix. The transpose can be obtained either with the 'transpose(VAR)' command, or just an apostrophe:

>> spleen = [1 2 3; 4 5 6]

spleen =

1 2 3

4 5 6

>> spleen'

1 4

2 5

3 6

>> spleen(1,:)

ans =

1 2 3

>> spleen(:,2)

2

5

Parentheses, ‘()’, are used for function arguments, order of operation and to select vector or matrix components. The ‘size(VAR)’ command returns the number of rows and columns of an M x N matrix:

>> spleen_size = size(spleen)

spleen_size =

2 3

>> spleen_size(2)

ans =

3

>> spleen_size(3)

??? Index exceeds matrix dimensions

To generate a vector with equally spaced elements between two values one can use:

>> t = 0:3:15

t =

0 3 6 9 12 15

Another approach would be to use ‘linspace(START,END,#_OF_POINTS)’. The command ‘logspace(START,END)’ will generate a vector with logarithmically equally spaced points between 10^START and 10^END. The default number of points for ‘linspace’ and ‘logspace’ are 100 and 50 respectively:

>> t2 = linspace(0, 15, 6)

t2 =

0 3 6 9 12 15

>> t3 = logspace(0,2,5)

t3 =

1.0000 3.1623 10.0000 31.6228 100.0000

In general, as long as the dimensions of the variables agree, vectors and matrices can be manipulated as easily as scalars:

>> kidney = transpose(spleen - 5)

kidney =

-4 -1

-3 0

-2 1

>> kidney*spleen

ans =

-8 -13 -18

-3 -6 -9

2 1 0

>> spleen*kidney

ans =

-16 2

-43 2

The 'find' command can be used to search vectors and matrices. It returns the indices of non-zero elements. To find elements of spleen greater than or equal to five:

>> [r c] = find(spleen >= 5)

r =

2

2

c =

2

3

Useful commands with their basic syntax and use include:

|Name |Syntax |Use |

|zeros |zeros(N) |creates an N x N zero matrix |

|ones |ones(N) |creates an N x N matrix of ones |

|eye |eye(N) |creates an N x N identity matrix |

|rand |rand(N) |creates an N x N matrix with random elements between 0 and 1 |

|sum |sum(VAR) |sum of elements in each row of VAR |

|prod |prod(VAR) |product of elements in each row of VAR |

|diag |diag(VAR) |the main diagonal of VAR |

|eig |eig(VAR) |eigenvalues of VAR |

|inv |inv(VAR) |inverse of VAR |

|rank |rank(VAR) |rank of VAR |

|max |max(VAR) |maximum value in VAR |

|mean |mean(VAR) |mean of VAR |

|abs |abs(VAR) |absolute value of elements of VAR |

|sort |sort(VAR) |sorts elements of VAR in ascending order |

Things to try:

- Create a 2 x 2 matrix m2= [1 2; 3 4]. Use the 'inv' command to get the inverse. Also, get the inverse using 'eye(2)/m2'. What is the rank of m2?

- Do the same as before with m3 = [1 2 3; 4 5 6;7 8 9]. What happens? What is the rank of m3? What will be the rank of the 15 x 15 matrix m15? These matrices can be made with the command: 'mN=transpose(reshape(1:N^2, N, N))'

5 - Polynomials

Many Matlab functions that will be used in this course take polynomials as arguments. Polynomials in Matlab are represented as vectors. The second degree polynomial [pic] is entered as:

>> bird = [1 2 3];

Note that 'size(POLY) -1' is equal to the degree of the polynomial. Any missing coefficients must be entered are zeros. Mathematically, vector multiplication is not a standard operation. Matlab allows component wise multiplication of equal length vectors with the times command or '.*'. This is not how polynomials are multiplied. Polynomial multiplication is performed with 'conv'. Polynomial addition can be performed using regular addition as long as you add sufficient leading zeros to the polynomial of smaller degree.

>> silo = [-4 5];

>> squish = conv(bird, silo)

squish =

-4 -3 -2 15

>> roots(squish)

ans =

-1.0000 + 1.4142i

-1.0000 - 1.4142i

1.2500

Other useful commands to explore can be found by looking at 'help polyfun'.

6 - Plotting

Matlab has many options for displaying data in two and three dimensional graphics. The function that will be used most often is 'plot'. The command 'plot(X,Y)' will plot the elements of Y versus the corresponding element of X. X and Y must be of the same length. 'plot(X)' will plot X versus the element indices. If X is a matrix, each row is plotted.

>> t = 0:pi/10:2*pi; x = sin(t); y = cos(t); plot(t,x)

[pic]

Matlab plots are saved as .fig files, but they can be exported in a variety of formats. The following command plots x vs t, y vs t and y vs x:

>> plot(t,x,t,y,x,y)

[pic]

There are many options that can be used to add or customize titles, labels, text, axes, legends and more. Look at 'help graphics' and 'help graph2d'. The following plots show a random matrix, its inverse, the m15 matrix (from section 4) and its inverse, versus each of their diagonals.

>> squiggles = rand(15); x = diag(squiggles);

>> sqinv = inv(squiggles); xinv = diag(sqinv);

>>

>> m15 = transpose(reshape(1:15^2,15,15)); m15diag = diag(m15);

>> m15pinv = pinv(m15); m15pinvdiag = diag(m15pinv);

>>

>> subplot(2,2,1), plot(x, squiggles), title('squiggles vs diagonal')

>> subplot(2,2,2), plot(xinv, sqinv), title('sqinv vs diagonal'), grid on

>> subplot(2,2,3), plot(m15diag, m15), title('m15 vs diagonal'), box off

>> subplot(2,2,4), plot(m15pinvdiag, m15pinv), title('m15pinv vs diagonal')

[pic]

If this image is displayed or printed in colour, you will notice that Matlab will automatically assign different colours each row of the matrix that it plots. It only has seven default colours and starts repeating after that. The plot in the lower righthand corner is now the current plot. Any new plotting commands entered will update that corner alone. The whole figure can be cleared with 'clf'.

Things to try:

-Plot the polynomial squish from the previous section. Are the roots where you expect them to be? The function 'polyval' will be useful for this exercise.

-There is a way to set the aspect ratio of a plot to 1. Find out how to do this.

7 - Relational Operators, Logic and Control Flow

Matlab is very similar to C in its implementation of relational and logical operators. They can be used either as a symbol or a command.

|Name |Matlab Symbol |Command |

|Equal |A == B |eq(A,B) |

|Not equal |A ~= B |ne(A,B) |

|Less than |A < B |lt(A,B) |

|Greater than |A > B |gt(A,B) |

|Less than or equal |A = B |ge(A,B) |

|And |A & B |and(A,B) |

|Or |A | B |or(A,B) |

|Not |~A |not(A) |

For loops, if statements, while loops and switches should be familiar to anyone with programming experience. A for loop and if statement are demonstrated in the following example. The m15 matrix can be made in the following manner:

>> for i=1:15

for j=1:15

m15nested(i,j) = 15*(i-1) + j;

end

end;

If the variable m15 from previous work is still in the workspace:

>> if m15nested == m15

'True'

else

'You messed up'

end

ans =

True

Note that, when within the for loop or if statement, evaluation is not carried out right away. 'end' is needed to terminate any of the control statements. Also, character strings are entered in Matlab as text between apostrophes. The command line is not the best place to implement for loops or any other complicated series of commands.

8 - m-Files and Functions

The examples presented so far have been entered in the command line. Scripts can be created to perform series of commands. Files with a .m extension will be evaluated in Matlab. These files can be either written as a text file with any text editor, or through Matlab with MATLAB Editor. At the command line, enter:

>> edit invspeed

MATLAB Editor will open. Write the following:

%invspeed.m

%invspeed.m calculates the time to perform two types of matrix inversions;

%it then plots the time versus matrix dimension.

for i=1:300

id=eye(i);

mat=rand(i);

tic;

id/mat;

x(i,1)=toc;

tic;

inv(mat);

x(i,2)=toc;

end

plot(1000*x)

title('Inversion Time vs Matrix Dimension')

legend('1/A', 'inv(A)')

xlabel('Matrix Dimension')

ylabel('Time (milliseconds)')

text(20, 200,['1/A total = ', num2str(sum(x(:,1))), ' seconds'])

text(20, 180,['inv(A) total = ', num2str(sum(x(:,2))), ' seconds'])

[pic]

The text in MATLAB Editor will be colour coded. Comments are green, commands for control flow (for, if, end …) as well as ‘function’ are blue and strings are red. New material used in this example includes:

- The use of '%' to enter comments. Everything on a line following a percentage sign is considered a comment. This is also useful to temporarily disable commands.

- 'tic' starts a timer. 'toc' ends the timer and returns the time elapsed. This is useful for determining the computational efficiency of your code. Computationally intensive calculations may be detrimental in real-time control applications.

- Some string manipulation was used to create the text with the totals. An array of strings was created, with the middle string converted from the numerical sum. Look at 'help strings'.

In addition to scripts for performing a series of commands, m-files are ideal for writing your own functions. A function can be written in the command line, but m-files are much easier for debugging. The basic syntax is demonstrated in the following example:

function lot(N, R)

%lot picks lottery numbers. N numbers from 1 to R are chosen.

if (N 0)

tot = randperm(R);

x = tot(1:N);

disp(['Your lottery numbers are: ', num2str(x)])

end

The variables ‘tot’ and ‘x’ are local, within the function. This function displays text but does not return a value. In order to return the variable ‘x’ to the workspace, we need to write ‘function x = lot(N, R) …’. The next example shows how to return several variables, as well as perform error checking on the function arguments.

function [msg, x] = nalot(S, A)

%nalot picks lottery numbers.

%nalot(S) uses your name to pick lottery numbers.

%nalot(S, A) uses your name and age to pick lottery numbers.

tot = randperm(49);

tot = cat(2,tot,tot(1:5));

t=1;

d=0;

if nargin==0

d=1;

elseif (nargin==1 & ischar(S))

n=double(S);

t=ceil(49*randperm(sum(n))/sum(n));

d=2;

elseif (nargin==2 & ischar(S) & isnumeric(A))

n=A*double(S);

t=ceil(49*randperm(sum(n))/sum(n));

d=3;

end

x=tot(t(1):t(1)+5);

xn=num2str(x);

switch d

case 1

msg=['Your lottery numbers are: ', xn];

case 2

msg=['Your lottery numbers, ', S, ', are: ', xn];

case 3

msg=['Your lottery numbers, ', S, ', age ', num2str(A) ', are: ', xn];

otherwise

msg=['You cannot win the lottery if your name is not a string and/or your age is not a

number'];

x=0;

end

The command ‘[txt, nums]=nalot(‘Edward Pepsi’, 25)’ will initialize two variables containing the text message, msg, and the lottery numbers, x. If the results are to be assigned to a single variable, only the text message will be assigned. Make sure that you understand how ‘randperm’, ‘nargin’, ‘switch’ and the ‘is’ functions work.

Other examples throughout this tutorial will continue to provide exposure to writing functions.

Things to try:

-Write a function that inverts a matrix using Gaussian Elimination. Compare it's performance to '1/A' and 'inv(A)'. The built in functions are written in a lower level language than Matlab.

9 - Control Toolbox

Mathworks develops specialized toolboxes that build upon Matlab. One toolbox that is very useful for control engineering (and is available both in the computer labs and on the computers in MEB 524) is the Control System Toolbox. ‘help control’ lists the functions within this toolbox. We will briefly introduce some of the functions that will be relevant to this course.

Most Matlab functions have no problem operating on scalars, vectors or matrices. The functions in the Control System Toolbox operate on transfer function or state-space objects. Where as transfer functions represent differential equations in the S-domain, a state-space representation converts systems of higher order differential equations into systems of first order differential equations. A numerator and a denominator define a transfer function. Four matrices describing the dynamics, input, output and any direct influence of the input on the output, characterize a state-space model (see section 2.3 and chapter 7 of Franklin for more information about state-space representations). Any transfer function can be converted to a state-space model, but not all state-space models can be converted into transfer functions. The following example illustrates initialization of a transfer function for the system described by the second order differential equation [pic] with an output [pic].

>> num=[0.01 1]; den=[1 7 10];

>> sys1=tf(num, den)

Transfer function:

0.01 s + 1

--------------

s^2 + 7 s + 10

The second order differential equation can be transformed into a system of two first order differential equations.

[pic]

>> a2=[-7 -10; 1 0]; b2=[1; 0]; c2=[0.01 1]; d2=[0];

>> sys2=ss(a2,b2,c2,d2)

a =

x1 x2

x1 -7 -10

x2 1 0

b =

u1

x1 1

x2 0

c =

x1 x2

y1 0.01 1

d =

u1

y1 0

Continuous-time model.

>> [num2 den2]=ss2tf(a2,b2,c2,d2); tf(num2,den2)

Transfer function:

0.01 s + 1

--------------

s^2 + 7 s + 10

The ‘feedback’ command can be used to convert and open loop system to a closed loop one. Impulse and step responses can be viewed using ‘impulse’ and ‘step’. Bode plots and root locus diagrams are displayed with ‘bode’ and ‘rlocus’. The response of a linear time invariant model to an arbitrary input can be simulated using the ‘lsim’ command. ‘gensig’ can be used to generate sine, square or pulse signals. ‘ltiview’ displays various responses and diagrams.

The SISO design tool is opened with ‘sisotool’. This is a GUI that has several useful features for the design of compensators for SISO systems. Different plants, sensor, filters and compensators can be inserted into the system. Root loci, bode plots and various loop responses can be viewed. ‘rltool’ opens sisotool set to view the root locus diagram.

10 - Numerical and Control Examples

The following examples demonstrate how Matlab can be prove useful for problems involving numerical or control system analysis. You may find that the comments do not always provide enough information for you to completely understand how the code works. Use help to find out how any new functions work. Lab work and numerical assignments will use and build upon this material.

Bode Plot:

A bode plot displays the magnitude and phase diagrams of the frequency response of a system. The following code shows one implementation for creating bode plots and compares the results to the built-in function bode. Enter the code into Matlab and note (at least) one significant difference between the two plots for arbitrary transfer functions.

function [m,p] = bode2(num,den)

%bode2(num,den) calculates the magnitude and phase of a system

%for sinusoidal inputs with frequencies between 0.001 and 1000 rad/s.

%The magnitude and phase diagrams are displayed along with the

%results of the built-in function bode.

%Note the use of .^ and ./ instead of ^ and / this ensures that

%the operation works on vector elements, not the entire vector.

sys=tf(num,den)

w=logspace(-3,3,1000);

n=zeros(1,length(w)); d=zeros(1,length(w));

for k=length(num):-1:1

n=n + num(k)*(j*w).^(length(num)-k);

end

for k=length(den):-1:1

d = d + den(k)*(j*w).^(length(den)-k);

end

g=n./d;

for i=1:length(w)

m(i) = norm(g(i)); %magnitude

p(i) = atan(imag(g(i))./real(g(i))); %phase

end

m = 20*log10(m); %conversion to decibels

p = 180/pi*p; %conversion to degrees

figure(1)

subplot(2,1,1)

semilogx(w,m)

subplot(2,1,2)

semilogx(w,p)

figure(2)

bode(sys, w)

The following plot is for the transfer function: [pic]

[pic]

Zero-Order Hold:

Problems can occur when converting an analog signal to a digital representation. One approach for analog to digital (A/D) conversion is zero-order hold sampling. The magnitude of the signal is measured at fixed intervals and is held constant between the sampling periods. Note that the original signal is itself sampled with a frequency of 0.001 and the plot joins points with straight lines. This example demonstrates manipulation of vectors.

function zoh(fs)

%zoh(fs)

%A sine wave of unit frequency is sampled with

%a frequency of fs. The original and reconstructed

%signals are plotted.

t = 0:0.001:8*pi;

sig = sin(t);

ts = 0:fs:8*pi;

ts = sort(cat(2,ts,ts));

sampledsig = sin(ts);

first = 2:length(ts);

last = 1:length(ts)-1;

ts = ts(first);

sampledsig=sampledsig(last);

plot(t,sig,ts,sampledsig);

This is the plot generated with the command ‘zoh(1)’.

[pic]

Numerical Integration and Differentiation:

There are many algorithms that have been developed for numerical integration and differentiation. This example demonstrates an implementation using the following discrete approximations:

[pic]

The continuous variable t is replaced by Tk, where T is the sampling period and k is the index of the current sample. The results displayed are for ‘ts=0.1;sig=sin(0:ts:10);nint(sig,ts)’.

function [int, diff] = nint(signal,ts)

%nint(signal)

%This function numerically integrates and differentiates

%a signal sampled with a period of ts.

int = zeros(1,length(signal));

diff = zeros(1,length(signal));

int(1) = 0;

diff(1) = 0;

t = 1:length(signal);

t = (t-1)*ts;

for i = 2:length(signal)

int(i) = int(i-1) + signal(i)*ts;

diff(i) = (signal(i) - signal(i-1))/ts;

end

plot(t,signal,t,int,t,diff)

[pic]

Spring Damper:

Use Matlab to simulate the behaviour of the system shown below.

[pic]

A mass is suspended from a datum by a spring and a damper. The equation of motion for this system is:

[pic]

The following code solves this equation using the Backward Euler method shown in the previous example:

function [x,v] = spring_damp(mass, spring_k, damp_c)

g = -9.81;

km = spring_k/mass;

cm = damp_c/mass;

x(1) = 0; %initial position

v(1) = 0; %initial velocity

t_max = 30;

dt = 0.01; %step size

t = 0:dt:t_max;

for j = 2:length(t);

x(j) = x(j-1) + v(j-1)*dt; %update position

v(j) = v(j-1) + (g - cm*v(j-1) - km*x(j-1))*dt; %update velocity

end

figure(1)

subplot(2,1,1)

plot(t,x)

title('Position')

subplot(2,1,2)

plot(t,v)

title('Velocity')

This simulation can also be performed using functions from the Control Toolbox. The equation of motion can be transformed into the S-domain using the Laplace transform and rewritten as a transfer function:

[pic]

Using ‘tf’ and ‘lsim’, the system can be simulated with the following code:

function [x,v] = spring_damp_tf(mass, spring_k, damp_c)

g = -9.81;

km = spring_k/mass;

cm = damp_c/mass;

t_max = 30;

dt = 0.01; %step size

t = 0:dt:t_max;

u = g*ones(1,length(t)); %step input

nump = [1];

numv = [1 0];

den = [1 cm km];

sysp = tf(nump,den);

sysv = tf(numv,den);

x = lsim(sysp,u,t);

v = lsim(sysv,u,t);

figure(1)

subplot(2,1,1)

plot(t,x)

title('Position')

subplot(2,1,2)

plot(t,v)

title('Velocity')

figure(2)

rlocus(sysp,sysv) %plot sysp blue, sysv green

Save these two files as spring_damp.m and spring_damp_tf.m.

................
................

In order to avoid copyright disputes, this page is only a partial summary.

Google Online Preview   Download