De La Salle University

Electronics and Communications Engineering Department

System Models
Laboratory Experiment 2
I.

Objectives :
[1.] To be able to use the ode23 and ode45 of Matlab in solving differential equation employing the
Runge-Kutta method.
[2.] To use MATLAB functions in finding the roots of the characteristic equation, poles and zeros of
transfer function, and transformation of poles and zeros to transfer function.
[3.] To represent an nth-order differential equation in state-space and transform the state-space to
poles and zeros and vice versa.
[4.] To use the control system toolbox in model conversion such as state-space to transfer function and
vice versa.

II. Introduction
Analysis of any system, and the process of design to achieve acceptable dynamic behavior,
requires a suitable mathematical model for the system. A mathematical model is an equation or set of
equations that describes the state of the system as a function of time.
2.1 Linear System
A linear differential equation is one which has constant coefficients and which contains no terms
involving products of variables of their derivatives. If a system is linear then the principle of superposition
applies, which means that when two forms of input function are applied simultaneously then the resulting
response is the sum of the responses to the inputs applied individually. As an example consider an RLC
network shown in figure 2-1.

Figure 2-1. RLC Circuit

Applying KVL, will result to

v i ( t ) i( t ) R L

di( t ) 1
i( t )dt 0
dt
C

Let the transfer function of the RLC circuit be defined as the output voltage over the input
voltage in the s-domain, the resulting transfer function can be obtained by first defining i(t) in terms of v o
(t).
Since the current flowing through the capacitor is given by

i( t ) C

dv o ( t )
dt

then the above differential equation can be rewritten as:

Feedlab

Laboratory Experiment 1

Page 2-1

dv o ( t )
d 2 v o (t)
v i ( t ) RC
LC
v o (t) 0
dt
dt 2
rearranging and multiplying both sides by 1/LC will result to,

d 2 v o (t)
dt 2

R dv o ( t )
1
1

v o (t)
v i (t)
L dt
LC
LC

Laplace transforming this equation:

H (s)

Vo (s)
1 / LC
2
Vi (s) s R / Ls 1 / LC

2.1.1 Polynomial Roots and Characteristic polynomial.

Let L = 0.2H, C = 0.01F, R = 2 then the transfer function is given by

H (s)

Vo (s)
500
2
Vi (s) s 10s 500

To encode this transfer function in MATLAB, we define two variables. The first variable refers to the
numerator polynomial while the second variable is for the denominator polynomial. The polynomial is
represented by a row vector containing the coefficients of a polynomial. For the given example, let num
and den be the variables for the numerator and denominator polynomial, respectively, then
>> num = [500];
>>den = [1 10 500];
To display the transfer function in MATLAB, you can use the command printsys. Printsys command will
writes out the transfer function equivalent to num(s)/den(s). To illustrate, entering
>>printsys(num,den)
will yield
num/den =
500
---------------s^2 + 10 s + 500
The roots of a polynomial can be obtained using the MATLAB function roots. For the denominator
polynomial, the roots are found by
>>r=roots(den)

The polynomial roots are obtained in column vector

r=
-5.0000+21.7945i

Feedlab

Laboratory Experiment 1

Page 2-2

-5.0000-21.7945i
Similarly, to construct a polynomial from its roots, use the MATLAB function poly.
>>d = poly(r);
Notice that the values of the vector d are equal to the values of vector den. If you use a matrix as the input
argument to poly, then the characteristic polynomial equation is returned.
To obtain the value of a polynomial at s = s 0, use the function polyval. For example, to obtain the value of
den(5)
>>den5 = polyval(den,5)
will result to
ans =
575
2.1.2 Poles and Zeros of a Transfer Function
In general an nth-order transfer function can be expressed:

H (s)

s M b M 1 s M 1 ... b 0

s N a N 1 s N 1 a N 2 s N 2 ... a 0

Factoring the numerator and denominator will result to a transfer function in terms of zeros, poles and
scalar gain.

H ( s) k

s z 0 s z 1 ... s z M 1
s p 0 s p 1 ... s p N 1

In MATLAB, the command tf2zp will convert the transfer function to zeros, poles and scalar gain. To
illustrate, consider a transfer function given below

H (s)

s 3 11s 2 30s
s 4 9s 3 45s 2 87s 50

>> num = [1 11 30 0];

>> den = [1 9 45 87 50];
>>[z, p ,k] = tf2zp(num,den)
The zeros, poles and gain are
z=
0
-6
-5
p=
-3.0000+ 4.0000i
-3.0000- 4.0000i
-2.0000
-1.0000

Feedlab

Laboratory Experiment 1

Page 2-3

k=
1
The command zp2tf is used to convert the zero-pole-gain system into a transfer function
described by the numerator polynomial and denominator polynomial.
2.1.3 State-space representation
Consider an nth-order linear system described by the differential equation:

dNc
dt

a N 1

d N 1 c
dt

N 1

... a 1

dc
a 0c r
dt

where c(t) is the output and r(t) is the input. The state variables, referred to as phase variables, is defined
as

x1
x2

d N 1 c
dt N 1
d N 2 c
dt N 2

xN c
expressing

dx k
dx N
from the given differential equation
x k 1 for k = 1, 2,...,N-1 and solving for
dt
dt

will result to the following equations

dx 1 d N c
N a N 2 x 2 a N 3 x 3 ... a 0 x N r
dt
dt
dx 2
d N 1 c
N 1 x 1
dt
dt

dx N
x N 1
dt

Expressing the above equation in matrix form

Feedlab

Laboratory Experiment 1

Page 2-4

 dx 1
dt a N 1 a N 1 a 1 a 0 x 1 1
dx
0
0
0 x 2 0
2 1
r
dt

0 0


0
1
0 x N 0
dx N 0
dt
The output can be defined as:

x1
x
2
y 0 0 0 1
x3

x4
where y is the output of the system.
For linear system, the general state equation can be described as follows

X' Ax Bu

y Cx Du

where A is an N x N system matrix, B is an N x 1 input matrix, C is a 1 x N matrix (output matrix), and D

is a scalar called the direct transmission term.
To illustrate, consider the differential equation shown below:

d 3c
dt

d 2c
dt

dc
4c 10r
dt

To encode this in MATLAB, we have the following set on instructions

>> A = [-2 -3 4;1 0 0;0 1 0];
>> B = [1 0 0];
>> C = [0 0 1];
>> D = 0;
>> printsys(A,B,C,D)
a=
x1
x2
x3

Feedlab

x1
x2
-2.00000 -3.00000
1.00000
0
0
1.00000

x3
-4.00000
0
0

Laboratory Experiment 1

Page 2-5

b=
x1
x2
x3

u1
1.00000
0
0

c=
y1

x1
0

x2
0

x3
1.00000

d=
y1

u1
0

2.1.4 Changing from transfer function to state-space and state-space to transfer function
If the given differential equation contains the derivative of the input, r, the output y (state-space)
is not one of the state elements unlike in the previous discussion. In other words, it is not straightforward
to write the state-variable description of the system. In this regard, it is advisable to determine the transfer
function first then convert it into state-space using the MATLAB function tf2ss. To convert the state-space
to transfer function use the function ss2tf.
To express the system in terms of poles, zeros, and gain from state-space representation and vice
versa, the MATLAB function ss2zp and zp2ss, respectively, will do the conversion.
2.2 Non-linear System
A great majority of the system are linear within some range of the variables. However, all
systems ultimately become nonlinear as the ranges are increased without limit. For the nonlinear system,
the superposition principle of superposition does not apply. The MATLAB function ode23 and ode45
simply the task of solving a set of nonlinear differential equation. Another way of solving non-linear
system is by approximating it by linear equation through linearisation or small perturbation analysis
An example of a non-linear system is a simple pendulum shown in figure 2-2. Let the weight W
= mg kg is hung from a support by a weightless rod of length L meters. Assuming that the viscous
damping coefficient is B kg/m/sec and if in radians is the angle of deflection of the rod, the velocity of
the weight at the end will be L and the tangential force acting to increase the angle can be written as:

Figure 2-2

mL

Feedlab

d 2
dt

BL

d
W sin 0
dt

Laboratory Experiment 1

Page 2-6

Let x 1 and x 2

d
then
dt

dx 1
x2
dt
dx 2
B
W
x2
sin x 1
dt
m
mL
For W = 2; L = 0.6, B = 0.02; g = 9.81, m = W/g, the solution for x 1 and x2 can be obtained using the
following MATLAB code.
First create a function pendulum. The function pendulum will calculate the state derivative. The following
is the code of the function pendulum
function xdot = pendulum(t,x)
%returns the state derivative
W = 2; L =0.6; B = 0.02; g = 9.81; m = W/g
xdot = [x(2) ; -B/m*x(2) - W/(m*L)*sin(x(1)) ];
Second, use the ode23 to simulate the system over an interval 0 to 5 sec.
>>
>>
>>
>>
>>
>>
>>
>>
>>
>>
>>

t0 = 0; tfinal = 5
%time interval
x0 =[1 , 0]
%initial conditions
t01 = 0.0001; trace = 0;
[t , x] = ode23(pendulum,t0,tfinal,x0,t01,trace);
subplot(211); plot(t,x);
title(Time response of rigid pendulum);
xlabel(time - sec);
text(3.2,3.5,velocity), text(3.2,1.2,Angle - rad)
th = x(:,1); w=x(:,2); subplot(212),plot(th,w);
title(Phase plane plot of pendulum);
xlabel(Position in rad);ylabel(Angular velcity)

III. Laboratory Exercise

1.

Consider the function

H ( s)

N ( s)
D( s)

where
N (s) s 3 8s 2 21s 20
D(s) s 4 14s 3 46s 2 64s 40

1.1 Find N(-12), N(-10), N(-8), D(-12), N(-10), and D(-8)

1.2 Find H(-12), H(-10), H(-8)
1.3 Now consider s = a + bi, where a and b are integers, -15 a 0 and -2 b 2. Generate a 16 x
5 matrix S of all possible values of s.
1.4 Find H(s) for the values of s considered in 1.3.
1.5 For what values of s considered in 1.4 does H(s) = 0 ? For what values does H(s) blow up ?
2.

The roots of a polynomial are -1, -2, -3 j4. Determine the polynomial equation.

Feedlab

Laboratory Experiment 1

Page 2-7

3.

Find the poles and zeros of the following transfer function:

H (s)

s 3 11s 30s
s 4 9s 3 45s 2 87s 50

4.

A system has zeros and poles at -6, -5, 0, poles at -3 j4, -2 , -1 and a gain of 1. Determine the
system transfer function

5.

A system is described by the following state-space equations

dx 1
dt
dx 0 1 0
2 0 1 1

dt

3
dx 3
dt

x 1 10
x 0 u
2
x 3 0

x1
y 1 0 0 x 2
x 3
Find the transfer function.
6.

The function ode23 and ode45 are functions for numerical solution of ordinary differential equations.
Consider a second order differential equation known as the Van der Pol equation:

d2x
dt

dxdt x 0

x2 1

Rewriting the above differential equation as a system of first order differential equation will result to:

dx 1
x1 1 x 22 x 2
dt
dx 2
x1
dt
Create a function M-file containing these differential equations. Call it vdpol.m Solve the differential
equation numerically using the ode23 or ode45. Plot the time domain response and the phase plane
plot.

Feedlab

Laboratory Experiment 1

Page 2-8