The code does not define ‘u’ anywhere except in the ‘springmass’ function, and since functions have their own workspaces, and it is not returned by ‘springmass’ as an output, lsim never sees it. (The lsim function uses the ‘t’ vector from the ode45 result.) The length of the ‘u’ vector (in a single-input system) or matrix (in a multiple-input system) must match the number of elements in ‘t’.
It might be better to define ‘u’ in the calling script as:
u = A*cos(omega*t);
using the ‘t’ output from ode45 to calculate it before the lsim call.
The tf2ss function is a Signal Processing Toolbox function, and will not return a system object that the Control System Toolbox can use. Using the ss funciton is an option, however once the system is created, just use it as it exists. There is no particualr reason to convert it to state-space.
Try this —
m = 250; m1 = m; m2 = m; % masses (all equal)
k = 50; k1 = k; k2 = k; k3 = k; % spring constants
c = 10; % damping
A = 0.00315; omega = 0.75; % forcing function
%
% tspan=[0 500]; % time range for simulation
tspan = linspace(0, 500, 1500); % time range for simulation
y0 = [0; 0; 0; 0]; % initial conditions
% Call ode45 routine
[t,y] = ode45(@springmass, tspan, y0, [], k1, k2, k3, m1, m2, c, A, omega);
%disp(y)
% Plot the input and outputs over entire period
u = A*cos(omega*t); % <— ADDED HERE
figure(1); clf
plot(t, u, t, y(:,1), t, y(:,2));
%
%
N1 = [12500 500 5500 20 150];
D1 = [0 0 0 0 50];
%Gm1 = tf(D1,N1);
N2 = [1250 50 1000 20 1500];
D2 = [0 0 250 0 100];
%Gm2 = tf(D2,N2);
%bode(Gm1,Gm2);
%------------------------------------------------------
sys1 = tf(D1,N1)
sys2 = tf(D2,N2);
%u = A*cos(omega*t);
y1 = lsim(sys1, u, t);
figure(4)
plot(t,y1)
grid
function dydt = springmass(t, y, k1, k2, k3, m1, m2, c, A, omega)
% compute the input to drive the system
u = A*cos(omega*t);
% compute the time derivative of the state vector
dydt = [
y(3);
y(4);
-(k1+k2)/m1*y(1) + k2/m1*y(2);
k2/m2*y(1) - (k2+k3)/m2*y(2) - c/m2*y(4) + k3/m2*u
];
end
.
