Hands On: State Space Descriptions

Table of Contents

Introduction A Linear Continuous Time System: a Mechanical System A Linear Discrete Time System: a Supply Chain System A Linear Discrete Time System: a linear Finite Difference Equation A Nonlinear Continuous-Time Dynamic System: a Double Pendulum The Model System's Parameters How to Solve ODEs in MATLAB Critical Remark: from an N-th Order ODEs Set to a Set of First-Order ODEs The MATLAB Command ode Important Remark on Notation How to Obtain a Set of First-Order ODEs The IVP Initial Conditions

Introduction

This live script illustrates some models of dynamical systems (linear or nonlinear, time-invariant or time-variant, continuous or discrete time).

Inspired by the examples illustrated in the live script StateSpaceDescriptionExamples.m, you will able to create a Matlab live script to simulate the behavior of a continuous time or discrete time dynamical system.

What you will learn:

Given a state space description of a dynamical system, how to build a MATLAB model for such a system;
Tune the model, setting the model parameters (or the simulation parameters) using a Matlab script;
How to run a MATLAB model and retrieve the results.

A Linear Continuous Time System: a Mechanical System

The dynamics of the system is given by the state-space model derived from conservation balances (of mass, energy, etc.). Applying Newton's second law, we obtain the following state-space description:

where the state variables are

This state-space description may be easily implemented in Matlab.

Proposal: implement a MATLAB LTI model and a Matlab live script to analyse the evolution of the above described mechanical system. Use the following parameters to configure the system model:

point masses [kg]: , and ;
viscous damping coefficients [N s/m]: ;
spring constants [N/m] : ;
initial conditions:
external forces (inputs):
simulation time interval [s]:

% Create a state-space description using a MATLAB LTI model, then assign the model parameters,

% the initial state, the simulation configuration.

% Using the same live script, run the model, and retrieve the results.

% Plot the state and output movements.

% help initial

N = 20000;

tSET = linspace(0,30,N);

test = (tSET >=10);

test2 = -2*test;

U2 = test2;

U1 = 2*ones(size(tSET));

whos U1 U2

  Name      Size                Bytes  Class     Attributes

  U1        1x20000            160000  double              
  U2        1x20000            160000  double              

U = [U1; U2];

whos U

  Name      Size                Bytes  Class     Attributes

  U         2x20000            320000  double              

figure;plot(tSET, test2);ylim([-3.5, +1.5]);xlim([0,30])

legend('$x_1$','$x_2$','$x_3$','$x_4$','Interpreter','latex')

Warning: Ignoring extra legend entries.

whos test test2

  Name       Size                Bytes  Class      Attributes

  test       1x20000             20000  logical              
  test2      1x20000            160000  double               

% [yFREE, ~, xFREE] = initial(LTIsys, x0,tSET);

% U1 = 2 * ones(size(tSET));

% [yLSIM, ~,xLSIM] = lsim(LTIsys, U, tSET);

% xTOT = xFREE+xLSIM;

% figure;

% plot(tSET,xTOT);

A Linear Discrete Time System: a Supply Chain System

The dynamics of the system is given by the state-space model derived from conservation balances (of mass, energy, etc.).

Spurchases the quantity of raw material at each month k;
A fraction of raw material is discarded; a fraction is shipped to producer P;
A fraction of product is sold by P to retailer R, a fraction is discarded;
Retailer Rreturns a fraction of defective products every month, and sells a fraction to customers.

Paying attention to the material flow in each work step, we obtain the following description:

where

month counter
raw material in S
products in P
products in R
products sold to customers

The following constraints apply to the model parameters

and

Proposal: create a MATLAB model and a live script to analyse the evolution of this discrete time dynamical system. Use the following parameters to configure the system model:

choose random values for the system parameters , ensuring that constraints on the parameters are respected;
assuming that , and can take on non-integer values, choose for the initial state three random values in the range ;
input of raw material:
simulation time interval [months]:

% Create a state-space description using a MATLAB model, then using 
% a live script assign the model parameters, the initial state,
% the simulation configuration. 
% Using the same live script, run the model, retrieve the result
% and plot the state and output movements.

A Linear Discrete Time System: a linear Finite Difference Equation

The state-space model is derived from a difference equation model (an Input/Output model):

Letting

one gets:

Proposal: create a MATLAB model and a live script to analyse the evolution of this discrete time dynamical system. Use the following parameters to configure the system model:

choose random values for the initial state
assuming as input the sequence:
simulation time interval:

% Create a state-space description using a MATLAB model, then using 
% a live script assign the model parameters, the initial state,
% the simulation configuration. 
% Using the same live script, run the model, retrieve the result
% and plot the state and output movements.
 

A Nonlinear Continuous-Time Dynamic System: a Double Pendulum

The double pendulum consists of two pendulums and two rigid rods. The pendulums are attached to each other, as in the following figure

Assumptions regarding the dynamic system's model:

the pendulum rods are mass-less and rigid;
the pendulum masses are point masses;
the only force acting on pendulums is the gravity force.

The Model

The ODEs describing the system's motion can be derived using the direct Newtonian method (refer to [3], for example):

where

and are the angles of the first and latter pendulum, respectively;
and are the (constant) lengths of the two rods;
and are the (constant) masses of the pendulums;
g is the gravitational constant.

Note: for simplicity, it is assumed that there is no input. Gravity acceleration is NOT considered an input of the system but simply a parameter of the ODEs.

The dynamic model of the double pendulum is a set of two second-order nonlinear Ordinary Differential Equations (ODEs).

As outputs of the system consider simply:

System's Parameters

clear variables
 
g = 9.80665;  % Acceleration due to gravity (m/s^2)
 
% first set of parameters
% L1 = 1.0;  % Length of the first pendulum arm (m)
% L2 = 1.0;  % Length of the second pendulum arm (m)
% m1 = 1.0;  % Mass of the first pendulum bob [= the ball, attached to the end of the rod]  (kg)
% m2 = 1.0;  % Mass of the second pendulum bob (kg)
 
% a different set
L1 = 5.0;  % Length of the first pendulum arm (m)
L2 = 2.5;  % Length of the second pendulum arm (m)
m1 = 1.0;  % Mass of the first pendulum bob [= the ball, attached to the end of the rod]  (kg)
m2 = 2.0;  % Mass of the second pendulum bob (kg)

How to Solve ODEs in MATLAB

Critical Remark: from an N-th Order ODEs Set to a Set of First-Order ODEs

The order of the ODE equations is crucial: to be able to solve numerically in MATLAB a system of (linear or nonlinear) ODEs of order n (with

), we must always rewrite the n-th order ODEs set as a set of first-order ODEs, properly coupled each to the others. Once we have obtained the set of first-order ODEs, we can

describe this set using the MATLAB command ode, in case the ODEs are nonlinear and then solve the ODEs;
exploit the Control System Toolbox functions (as, for example, lsim) in case the ODEs are describing a linear time-invariant system (refer to the LTI Systems devoted live script in this collection).

The MATLAB Command ode

For help on using the command ode, and examples, refer to the online documentation or execute the following piece of code in the MATLAB Command Window

doc ode

Briefly, the ode command allows us to solve Initial Value Problems (IVPs) of the form

where y is a (vector) function (the solution of the ODE we are looking for) and p is a vector of parameters.

Thus, since the ODEs

are already a set of first-order ODEs (nonlinear ODEs), we can describe the ODEs set structure using the command ode.

Important Remark on Notation

It's not mandatory to obey the notation

. Indeed, we will use a different notation for IPVs, when they describe the evolution of a dynamic system (linear or nonlinear), keeping the notation x for the state variables (the variables affected by the ODEs), and using the notation y for the system's outputs. Thus, our default notation for an IVP is

How to Obtain a Set of First-Order ODEs

Given the ODEs set

, let's transform it into a set of first-order ODEs with proper initial conditions:

The dynamic system is mechanical; thus, the positions and velocities (linear, angular) of rigid bodies’ centers of mass are candidates for state variables.
The rods are mass-less, so the state variables are the angular positions and velocities of the two masses.

Using this set of variables, let's transform the Eq.s

into a set of first-order ODEs, i.e. a state-space description

The IVP, we would like to solve in MATLAB, is described by the Eqs

together with the initial conditions

where

is the initial time instant.

The IVP Initial Conditions

t0    =  0; % starting time instant [s]
tStop = 20; % end of the simulation [s]
 
%
barTheta1 = pi/2; % Initial angle of the first pendulum (radians)
barOmega1 = 0;    % Initial angular velocity of the first pendulum (rad/s)
 
barTheta2 = pi/2; % Initial angle of the second pendulum (radians)
barOmega2 = 0;    % Initial angular velocity of the second pendulum (rad/s)
%}
 
%{
% a different IC configuration
barTheta1 = pi/2; % Initial angle of the first pendulum (radians)
barOmega1 = 0;    % Initial angular velocity of the first pendulum (rad/s)
 
barTheta2 = deg2rad(-10); % Initial angle of the second pendulum (radians) *** try -8 deg ***
barOmega2 = 0;    % Initial angular velocity of the second pendulum (rad/s)
%}

% Create a state-space description using a MATLAB model, then using 
% a live script assign the model parameters, the initial state,
% the simulation configuration. 
% Using the same live script, run the model, retrieve the result
% and plot the state and output movements.