## Wednesday, 13 July 2011

### Simulation of a system described by its Dynamics

The present post and following few posts will talk about the simulation of any system described by its dynamics in MATLAB. We will see how we can simulate an open system and closed system with some small examples.

System Dynamics
Any system is described by its dynamics, Dynamics, in simple words, means how the outputs are dependent on current input and some internal variables called state along with how derivatives of  these states are behaving. Take this example of an electronic circuit with a Resistance R, an inductance L and a capacitance C connected as in Figure 1. Output is the voltage across L and input is the voltage source. This system has two internal parameters Vc and i which affect the output as

$V_o = V-(V_c + iR)$
Now this Vc and i are dynamic and their derivatives are given as
$\frac{dV_c}{dt}=\frac{i}{C}$
$\frac{di}{dt}=\frac{1}{L}(V-(V_c + iR))$

As you can see derivatives of states are dependent on state and inputs and so is output. This type of representation is also called state space model.
 Fig 1: LCR System with two states

For simulation purposes, we will take two models- one is simple integrator system and second is cruise system.

System1-Integrator model
This system has only one state y which is its output too. State dynamics are given as
$\frac{dy}{dt}=u(t)$
where u(t) is the input to this system.

To create this model, we will make a function sys1 with input u, y and t and which returns derivative term of y. This function describes the dynamics of system, or in other words it represents the system itself.
function dybydt = sys1(t,y,u)
dybydt = u ;
Building an Environment
Now build an environment which contains the model and its stimulus input. Let us take step signal as input for this example.
function dybydt = env1(t,y)
u=1*(t>0);
dybydt = sys1(t,y,u) ;
Simulating the system
The above build environment can be simulated using ode45. Create a blank script file and write the following in that.
%define timespan
tspan=[0 10];
%define initial value
y0=8;
[t y]= ode45(@env1,tspan,y0);
The output is defined as y and can be plotted by
plot(t,y);

as in Fig 2.
 Fig 2: Output of system 1

1. As an improvement for the env1 function I might suggest:

==================
function dybydt = env1(t,y)

step_finalvalue = 1;
step_time = 2;

u=step_finalvalue*(t>step_time);

%%%% incorporating the response into system
dybydt = sys1(t,y,u) ;
=============================

And a question, if I may. You use time span as tspan=[0 10]; But the question is - what is the granularity of the time span? I mean, the ODE is solved numerically - can we use fixed time step for this?

I'm wondering if we can do the same with discrete control systems, where sampling time is defined explicitly.

And again, huge thanks for the blog. Keep posting, it really helps!

2. Thanks for your suggestion,
In that way, you can also give different inputs to the system like sin wave or others :)

yes you can use fixed time step for this, but I wont recommend that,
See ode45 help, you can pass option structure with the options you want.

Discreet time systems are different than this, there we don't have any differential terms,
only difference terms are there

In difference terms, no time vector is involved, so accuracy of result does not depend on which time step you choose, while in differential terms (in continuous systems), time step is important, and variable time step is used because the variation in signal is not constant and may change much as time progresses.

You need not use ode45 for discreet time system. These are simpler system. Simple vector of previous states and LDE update of this states will do the job.

Do tell me if you need a blog on that

3. @Gupta, Abhishek Kumar said...

yes you can use fixed time step for this, but I wont recommend that
I know what you mean, but the reason I want to do this is to have a total control on the simulation. Let me explain.

I have simulated that example and compared it with Simulink - works fine. But I was interested in the values of controller's response (value v in this example). So I used dlmwrite funtion to output the value. My step time was 3 seconds and final value = 1. Now that was what I saw:
- ode45 chosen large time step, and proceeded from 0 to 3.5, than got LARGE values.
- ode45 then went back in time (!), chosen smaller step, went forwards. Still got large values.
- ode45 went backwards in time again, used even smaller values, and proceeded ahead with output, which was close to that one of Simulink.

On the final plot there is nothing like this - it looks neat and nice. But I don't like when too-smart-algorithms do something nasty behind my back :-)

Anyway, thanks again for the clarifications. In the meantime I have implemented my own Runge-Kutta solver, and I can visually study the influence of the step on the accuracy of solution.

Discreet time systems are different than this, there we don't have any differential terms, only difference terms are there
Yes, that's right. I asked just to check if I understand everything correctly.

You need not use ode45 for discreet time system. These are simpler system.
OK, I got it.

Do tell me if you need a blog on that
Well, if you plan to write about simulation of discrete systems, it would be great.

And thanks again for the help.

4. If you have your own solver instead of inbuilt like 0de45, there is one more way to solve this problem (different architecture basically).
Here, you are collecting all the states in to the environment and feeding it to the main solver(ode45). In the different approach you can solve the system locally and feed only outputs to the main environment (easier to implement and much clean code).

I will try to write one more article on comparison between these approaches.

5. Dr Abhishek Kumar Gupta,

As per the above :『For simulation purposes, we will take two models- one is simple integrator system and second is cruise system.』. May I know if there is case study "cruise system" ?

Thank you for your kind attention
Andrew