How do I change the output of the graph from velocity to acceleration?

I need to graph the acceleration for a spring-mass oscillator, but I am a bit lost. This is the current code, where m = mass, c = damping constant, k = spring constant, and f = force.

m = 46; % mass

k = 3220; % spring constant

c = sqrt(184*k); % damping constant

f = 46*(9.8); %force

Time = 0:0.1:1; %1X2 matrix

IC = [1;0]; %2X1 matrix of initial conditions

myoptions = odeset('RelTol',1.e-8);

[x,y] = ode45(@Sprink, Time, IC , myoptions, m, c, k, f);

figure

plot(x,y)

title('Harmonic Oscillation')

xlabel('Time 0 to 1 with step size 0.1')

ylabel('Displacement (b) and Velocity (r)')

legend({'y = displacement','y = velocity'},'Location','northeast')

The ODE is broken down into a first order ODE by

function xprime = Sprink(t, x, m,c,k,f)

xprime = [x(2); (1/m)*(f-c*x(2)-k*x(1))];

Currently, the code only graphs the position and velocity of the system. How would I go about making it graph the acceleration as well?

Answer 1

[Dyuman Joshi]

m = 46; % mass
k = 3220; % spring constant
c = sqrt(184*k); % damping constant
f = 46*(9.8); %force
Time = 0:0.1:1; %1X2 matrix 
IC = [1;0]; %2X1 matrix of initial conditions
myoptions = odeset('RelTol',1.e-8); 
[x,y] = ode45(@Sprink, Time, IC , myoptions, m, c, k, f);

%Accelration
y(:,3)=(f-c*y(:,2)-k*y(:,1))/m;

figure
plot(x,y)
title('Harmonic Oscillation')
xlabel('Time 0 to 1 with step size 0.1')  
ylabel('Displacement (b) and Velocity (r)')
legend({'y = displacement','y = velocity','y = accelaration'},'Location','southeast')

function xprime = Sprink(t, x, m,c,k,f) 
xprime = [x(2); (1/m)*(f-c*x(2)-k*x(1))];
end

Answer 2

[Walter Roberson]

Use gradient(x, y(:, 2)) to estimate the acceleration