Simulink integrator error, please help

What should I put into the integrator to make ti work?

Controller

%sinistra uscita----destra ingresso
function [u,V]=controller(x)

%unpack inputs
x1=x(1);
x2=x(2);
x3=x(3);
x4=x(4);

%constant parameters
k=100;
g=9.81;
m=1;

%Lyapunov function
V=0.5*(x1^2+x2^2);

%controller law
u=(k*x2^2-10*x1^2+x1*x2+(m*x1^2*(-x2-1)-x2)*x2)/(g*(x1^2+x2^2+2)*x2);

end

BackStepping1

%sinistra uscita----destra ingresso
function [u1,V1]=backstepping1(phi,dphi,x,V)

%unpack inputs
x1=x(1);
x2=x(2);
x3=x(3);
x4=x(4);

%constant parameters
k=100;
g=9.81;
m=1;

%first Lyapunov function
V1=V+0.5*(x3-phi)'*(x3-phi);

%first controller law
u1=1/(1+x1^2)*(dphi+g*x2*(x1^2+x2^2+2)-x1^2*(pi-sin(x1)*(x2)-k*(x3-phi)));

end

Backstepping2

function [u2,V2]=backstepping2(u,du,x,V1)

%phi-->u
%dphi--->du

%unpack inputs
x1=x(1);
x2=x(2);
x3=x(3);
x4=x(4);

%constant parameters
k=100;
g=9.81;
m=1;

%second Lyapunov function
V2=V1+0.5*(x4-u)'*(x4-u);

%second controller law
u2=1/(1+x1^2)*(du-x3*(1+x1^2)-k*(x4-u));

end

Dynamic

function dx=dynamic(x,u)

%unpack inputs
x1=x(1);
x2=x(2);
x3=x(3);
x4=x(4);

%constant parameters
k=100;
g=9.81;
m=1;

%dinamica del sistema

dx1=-10*x1+x2;
dx2=-g*x3*(x1^2+x2^2+2)+m*x1^2*(-x2-1)-x2;
dx3= x1^2*(pi-sin(x1)*x2)+(1+x1^2)*x4;
dx4=(1+x1^2)*u;

dx=[dx1; dx2; dx3; dx4];

end

Here's a photo of the simulink scheme

Error:An error occurred while running the simulation and the simulation was terminated

Caused by:

Derivative of state '2' in block 'altra_prova_backstep/Integrator' at time 0.007 is not finite. The simulation will be stopped. There may be a singularity in the solution. If not, try reducing the step size (either by reducing the fixed step size or by tightening the error tolerances)

[Sam Chak]

Hi @Federico Ennio,

Based on the knowledge of system dynamics, there doesn't appear to be anything wrong with the Integrator block in your Simulink model. The reason is that the system itself (function name: dynamic) is initially open-loop unstable. As a result, the responses of the states grow rapidly (diverging) until the solver in Simulink fails to meet the integration tolerances, leading to the termination of the simulation.

To address this issue, it's important to ensure that the signal u is properly designed to stabilize the system. Without validating your stability proof, it becomes difficult for us to understand what's going on and verify if your signal u is correctly designed in the code.

[t, x] = ode45(@dynamic, [0 10], [1 0 0 0]);

Warning: Failure at t=1.504024e+00. Unable to meet integration tolerances without reducing the step size below the smallest value allowed (3.552714e-15) at time t.

plot(t, x), grid on

function dx=dynamic(t, x)
    %unpack inputs
    x1=x(1);
    x2=x(2);
    x3=x(3);
    x4=x(4);
    %constant parameters
    k=100;
    g=9.81;
    m=1;
    
    %% set u = 0 to test open-loop stability
    u   = 0;
    % ----------------------
    
    %dinamica del sistema
    dx1=-10*x1+x2;
    dx2=-g*x3*(x1^2+x2^2+2)+m*x1^2*(-x2-1)-x2;
    dx3= x1^2*(pi-sin(x1)*x2)+(1+x1^2)*x4;
    dx4=(1+x1^2)*u;
    dx=[dx1; dx2; dx3; dx4];
end

Answer 1

[Sam Chak]

Hi @Federico Ennio,

Considering the system dynamics,

the controller can be designed as follows

,

where

,

,

.

However, it's important to note that this controller does not guarantee global stability. The initial values of and have the potential to destabilize the system.

[t, x] = ode45(@dynamic, [0 10], [1 0 0 0]);

TL = tiledlayout(2,2,'TileSpacing','Compact');

nexttile
plot(t, x(:,1)), grid on
title('x_{1}')

nexttile
plot(t, x(:,2)), grid on
title('x_{2}')

nexttile
plot(t, x(:,3)), grid on
title('x_{3}')

nexttile
plot(t, x(:,4)), grid on
title('x_{4}')

title( TL, 'Responses of the States')
xlabel(TL, 'Time (sec)')
ylabel(TL, 'Magnitude')

function dx=dynamic(t, x)
    % unpack inputs
    x1  = x(1);
    x2  = x(2);
    x3  = x(3);
    x4  = x(4);
    
    % constant parameters
    k   = 100;
    g   = 9.81;
    m   = 1;
    
    % Controller
    wn  = 0.5;
    x2d = 0;
    k2  = 4*wn^3;
    x3d = - (m*(x1^2)*(x2 + 1) + x2 - k2*(x2 - x2d))/(g*(x1^2 + x2^2 + 2));
    k3  = 6*wn^2;
    x4d = 1/(1 + x1^2)*(- k3*(x3 - x3d) + (x1^2)*(x2*sin(x1) - pi));
    k4  = 4*wn;
    u   = - k4*(x4 - x4d)/(1 + x1^2);
    
    % dinamica del sistema
    dx1 = - 10*x1 + x2;
    dx2 = - g*x3*(x1^2 + x2^2 + 2) + m*x1^2*(- x2 - 1) - x2;
    dx3 =  x1^2*(pi - sin(x1)*x2) + (1 + x1^2)*x4;
    dx4 =  (1 + x1^2)*u;
    dx  =  [dx1; dx2; dx3; dx4];
end