c2d

Convert model from continuous to discrete time

Contents

Syntax

SYSD = c2d(SYSC,TS) computes a discrete-time model SYSD with time step TS using one step forward euler approximation from an continuous-time model SYSC.

Description

For the discretization the one step forward euler approximation is used: xd = (xp - x)/TS. The algorithm directly implements additional auxilary equations for all state predictions and replaces the state derivatives by auxiliary variables by manipulation of the parameter tensor (the strucutre and parameter matrix of the CP decomposition). Afterwards it will be checked if those additional auxilary variables and equations can be reduced algebraicly, if so, they will be eliminated by algbraic elimination by direct manipulation of the parameter tensor (the strucutre and parameter matrix of the CP decomposition).

Examples

Discretization and simulation of the Dadras Attractor using the dmss-class and c2d function

The Dadras attractor, introduced in [1], can be described by the following set of ODEs

eq = ["xp1 = x2 - a*x1 + b*x2*x3";...
    "xp2 = c*x2 - x1*x3 + x3";...
    "xp3 = d*x1*x2 - h*x3"];

which can be converted into a continuous time dmss object

ts = 0;
sys = sym2dmss(eq, ts);
disp(sys.symbolicEquations)
Reduced by 2 column(s) and 0 equation(s) with trivial Reduction due to duplications.
xp1 - x2 + a*x1 - b*x2*x3 == 0
  xp2 - x3 - c*x2 + x1*x3 == 0
     xp3 + h*x3 - d*x1*x2 == 0
 

For simulation the symbolic parameters should be replace by numeric values

symParameters = sym(["a", "b", "c", "d", "h"]);
numParameters = [3, 2.7, 4.7, 2, 9];
sys = sys.replaceSymbolicParameters(symParameters, numParameters);
disp(sys.symbolicEquations)
3*x1 - x2 + xp1 - (27*x2*x3)/10 == 0
  xp2 - x3 - (47*x2)/10 + x1*x3 == 0
           9*x3 + xp3 - 2*x1*x2 == 0
 

The model can be converted into discrete time model by

tsD = 0.01;
sysD = c2d(sys, tsD);
disp(sysD.symbolicEquations)
System reduced algebraicly by 3 algebraic variable(s) and 3 equation(s).
         xp1 - x2/100 - (97*x1)/100 - (27*x2*x3)/1000 == 0
(1000*xp2)/1047 - (10*x3)/1047 - x2 + (10*x1*x3)/1047 == 0
                       xp3 - (91*x3)/100 - (x1*x2)/50 == 0
 

With simulation the phase potrait can be computed for the continuous- and discrete-time models. Also the simulation times can be compared.

x0 = [5, 0, -4];
tend = 1e2;

tic()
simoutC = dmsim(sys, x0, [0, tend]);
timeC = toc();

tic()
simoutD = dmsim(sysD, x0, [0:tsD:tend]);
timeD = toc();

figure()
plot3(simoutC.x(:,1),simoutC.x(:,2),simoutC.x(:,3), "r-", LineWidth=0.01);
hold on
plot3(simoutD.x(:,1),simoutD.x(:,2),simoutD.x(:,3), "k--", LineWidth=0.01);
view([45, 30])
grid on
legend([sprintf("continuous model: %.3f sec.", timeC), sprintf("discrete model: %.3f sec.", timeD)], "Location","southeast")
hold off
Seperation in 3 subset(s) with 3 subproblem(s), from which 3 are explicit solvable.

In the phase potrait it can be seen that the discretization is stable and behaves similar (within a similar magnitude) as the continuous time model. Due to the chaotic behaviour of the system it is extremely unlikely to generate exactly similar trends.

References

[1] S. Dadras and H. R. Momeni, “A novel three-dimensional autonomous chaotic system generating two, three and four-scroll attractors,” Physics Letters A, vol. 373, no. 40, pp. 3637–3642, Sep. 2009. doi:10.1016/j.physleta.2009.07.088

See also

d2c, d2d, dmss, dmsim, sym2dmss

Author(s): Torben Warnecke