One problem is with the ‘dxdt’ assignment. MATLAB interprets the spaces as delimiters, so the second row has 3 columns, as MATLAB sees it. Put parentheses around the second row terms (to preserve readability) and after providing ‘F’ as an argument rather than internally as a vector, and a few other tweaks, it works —
t = 0:20; % Time duration
x0 = [1; 6]; % Initial conditions
F = 1:20; % For simplicity I have taken force as an array
for k = 1:numel(F)
[t, x] = ode45(@(t,x)eqm(t,x,F(k)),t,x0); % ode45 command
tv{k} = t;
xv{k} = x;
end
figure
for k = 1:numel(F)
subplot(5,4,k)
plot(tv{k}, xv{k})
title(sprintf('k = %2d',k))
grid
end
% Function
function dxdt = eqm(t, x, F)
m = 1; % mass
c = 1; % damping
k = 1; % stiffness
dxdt = [x(2); (F/m - c/m*x(2) - k*x(1))];
end
.
