matlab – Spring/Damper Calculation & Plotting
matlab – Spring/Damper Calculation & Plotting
Answer updated to match Simulink model implementation
To use ode45, you first need to write a function that computes the derivative of you input vector (i.e. your differential equation), and store that function in a separate file with the function name as the filename. Please note that the ode solvers can only solve first-order differential equations, so you first need to do a bit of work to convert your second-order differential equation to a first-order one. For more details, see the documentation on ode45.
Based on what you have done in your Simulink model, D2y is known for all values of t (its the step input), so we need to integrate it with respect to time to get Dy and y. So our state vector is X = [x; Dx; y; Dy] and our function looks like (stored in diff_eqn.m):
function dX = diff_eqn(t,X)
m1=500;
c=1200;
k1=25000;
dX(1) = X(2); % Dx
dX(2) = -(1/m1)*(c*(X(2)-X(4)/m1) + k1*(X(1)-X(3)/m1));; % D2x
dX(3) = X(4); % Dy
if t<2
dX(4) = 0; % D2y
else
dX(4) = 0.5;
end
as dX = [Dx; D2x; Dy; D2y].
In your script or your MATLAB command window, you can then call the ode solver (initial conditions all being equal to zero for Dx, x, Dy and y, as per your Simulink model):
[t,X] = ode45(@diff_eqn,[0 20],[0; 0; 0; 0]);
Adjust the ode solver options (e.g. max step size, etc…) to get results with more data points. To get the same plot as in your Simulink model, you can then process the results from the ode solver:
D2x = diff(X(:,2))./diff(t);
D2x = [0; D2x];
D2y = zeros(size(D2x));
D2y(t>=2) = 0.5;
plot(t,[D2y 500*D2x])
grid on
xlabel(Time [s])
legend(D2y,m1*D2x,Location,NorthEast)
Which gives the following plot, matching the results from your Simulink model:
