Error in dynamic system integration

Hello everyone, I am writing a control synthesis program for a dynamic system. First I discretize it, then I find a control and apply it. When I integrate the system, a strange error occurs, which says that I have symbolic variables somewhere in the integration. But I checked, they are not there. Here is my code:

clear
set(0,'DefaultTextInterpreter', 'latex');
m = 0.127;
M = 1.206;
I = 0.001;
L = 0.178;
Bc = 5.4;
Bp = 0.002;
g = 9.8;
h=0.1;

A_0 = [m+M -m*L; -m*L I+m*L^2];
A_1 = [Bc 0; 0 Bp];
A_2 = [0 0; 0 -m*g*L];
A = [zeros(2) eye(2); -inv(A_0)*A_2 -inv(A_0)*A_1];
B = [0;0; inv(A_0)*[1;0]];

A_d=expm(A*h);
syms s;
A_exps=eye(4);
for k=1:40
    A_exps=A_exps+(A*s)^k;
end
B_d=int(A_exps*B, s, 0, h);

evs = eig(A_d)

C=[B_d A_d*B_d A_d^2*B_d A_d^3*B_d];

evs(2) = 0.19343;
Xi_A1 = (A_d-evs(1)*eye(4))*(A_d-evs(2)*eye(4))*(A_d-evs(3)*eye(4))*(A_d-evs(4)*eye(4));
Theta1 = -[0 0 0 1]*inv(C)*Xi_A1;

X0 = [0.0;-0.6; 0.0; 0.0];
t_k2=0.0;
count = 10;
t_all=[];x_all=[];
for k=1:count
    t_k1=t_k2
    t_k2 = t_k1+h
    [t_tmp,X_tmp]=ode45(@(t,X)pendulum(t,X,Theta1),[t_k1 t_k2],X0) %This is line 71
    t_all=[t_all;t_tmp];
    x_all=[x_all;X_tmp];
    X0=X_tmp(end,:)';
end

figure; 
subplot(1,1,1); 
plot(t_all,x_all(:,1));

function dXdt = pendulum(~,X,Theta)
m = 0.127; M = 1.206; I = 0.1e-2; L = .178; Bc = 5.4; Bp = 0.2e-2;g=9.8;
x = X(1);
phi = X(2);
dx = X(3);
dphi = X(4);

A0=[m+M -m*L*cos(phi);-m*L*cos(phi) I+m*L^2];
b=[Theta*X-Bc*dx-m*L*dphi^2*sin(phi); -Bp*dphi+m*g*L*sin(phi)];

dXdt = [dx;dphi;inv(A0)*b];
end             

And here's the error:

Error using odearguments (line 113)
Inputs must be floats, namely single or double.

Error in ode45 (line 115)
  odearguments(FcnHandlesUsed, solver_name, ode, tspan, y0, options, varargin);

Error in LW21 (line 71)
    [t_tmp,X_tmp]=ode45(@(t,X)pendulum(t,X,Theta1),[t_k1 t_k2],X0)
 

Can you please tell me what could be causing this and how can I fix it?

Answer 1

[Walter Roberson]

set(0,'DefaultTextInterpreter', 'latex');
m = 0.127;
M = 1.206;
I = 0.001;
L = 0.178;
Bc = 5.4;
Bp = 0.002;
g = 9.8;
h=0.1;

A_0 = [m+M -m*L; -m*L I+m*L^2];
A_1 = [Bc 0; 0 Bp];
A_2 = [0 0; 0 -m*g*L];
A = [zeros(2) eye(2); -inv(A_0)*A_2 -inv(A_0)*A_1];
B = [0;0; inv(A_0)*[1;0]];

A_d=expm(A*h);
syms s;
A_exps=eye(4);
for k=1:40
    A_exps=A_exps+(A*s)^k;
end
B_d=int(A_exps*B, s, 0, h);

evs = eig(A_d)

evs = 4×1

1.0000 1.9343 0.4718 0.6769

C=[B_d A_d*B_d A_d^2*B_d A_d^3*B_d];

evs(2) = 0.19343;
Xi_A1 = (A_d-evs(1)*eye(4))*(A_d-evs(2)*eye(4))*(A_d-evs(3)*eye(4))*(A_d-evs(4)*eye(4));
Theta1 = -[0 0 0 1]*inv(C)*Xi_A1;

X0 = [0.0;-0.6; 0.0; 0.0];
t_k2=0.0;
count = 10;
t_all=[];x_all=[];
for k=1:count
    t_k1=t_k2;
    t_k2 = t_k1+h;
    [t_tmp,X_tmp]=ode45(@(t,X)pendulum(t,X,Theta1),[t_k1 t_k2],X0); %This is line 71
    t_all=[t_all;t_tmp];
    x_all=[x_all;X_tmp];
    X0=X_tmp(end,:)';
end

figure; 
subplot(1,1,1); 
plot(t_all,x_all(:,1));

function dXdt = pendulum(~,X,Theta)
m = 0.127; M = 1.206; I = 0.1e-2; L = .178; Bc = 5.4; Bp = 0.2e-2;g=9.8;
x = X(1);
phi = X(2);
dx = X(3);
dphi = X(4);

A0=[m+M -m*L*cos(phi);-m*L*cos(phi) I+m*L^2];
b=[Theta*X-Bc*dx-m*L*dphi^2*sin(phi); -Bp*dphi+m*g*L*sin(phi)];

dXdt = double([dx;dphi;inv(A0)*b]);
end             

[Egor Semenuta]

And due to the fact that a double is used there, the accuracy will not be lost?

[Walter Roberson]

I = 0.1e-2

That is only one significant digit. g=9.8 is only two significant digits. Why are you concerned about whether double precision (15 to 16 significant digits) is accurate enough when your equations have high error built in to them?