There are two problems:
First, ‘x0’ is initially identically zero, and that creates an Inf result for the first value of ‘A’ and NaN for the rest, none of which will be plotted, so I added eps to ‘x0’.
Second, the system is ‘stiff’ so you need a ‘stiff’ solver, such as ode15s or ode23s (that I chose here, for no particular reason other than it being a stiff solver, either one would work here).
Then using both of those changes, it works —
all_units = 'all_units.xlsx';
maks = xlsread(all_units);
% Vehicle parameters
a = 1.024;
b = 1.602;
m = 1574.02;
cr = -25488;
cf = -184780;
I = 1790;
% Extract relevant columns from the data
v = maks(:, 24);
r = maks(:, 20);
u = maks(:,23);
t1 = maks(:, 25);
del = (maks(:, 18) + maks(:, 19)) / 2;
% Time span for the simulation
tspan = [min(t1), max(t1)];
% Initial condition
x0 = [0; 0]; % Initial values for v and r
x0 = x0 + eps;
% Solve the differential equation using ode45
[t, x] = ode23s(@(t, x) myVehicleODE(t, x, t1, del), tspan, x0);
% Plot the results
figure;
subplot(2, 1, 1);
plot(t, x(:, 1), 'b', 'LineWidth', 2);
title('Longitudinal Velocity (v)');
xlabel('Time');
ylabel('Velocity');
subplot(2, 1, 2);
plot(t, x(:, 2), 'r', 'LineWidth', 2);
title('Yaw Rate (r)');
xlabel('Time');
ylabel('Yaw Rate');
% Define the function myVehicleODE at the end of the script
function dxdt = myVehicleODE(t, x, t1, del)
% Parameters
a = 1.024;
b = 1.602;
m = 1574.02;
cr = -25488;
cf = -184780;
I = 1790;
% Interpolate the input function 'del' with respect to time 't'
del_interpolated = interp1(t1, del, t, 'linear', 'extrap');
% System matrices
a11 = (-cf - cr) / (m * x(1));
a12 = (((-cf * a) + (cr * b)) / (m * x(1))) - x(1);
a21 = ((-cf * a) + (cr * b)) / (I * x(1));
a22 = ((-cf * a^2) - (cr * b^2)) / (I * x(1));
b1 = cf / m;
b2 = (cf * a) / I;
% State-space matrices A and B
A = [a11, a12;
a21, a22];
B = [b1;
b2];
% System of differential equations
dxdt = A * x + B * del_interpolated;
end
.
