Plot the acceleration with ODE45

Hello, I have the problem to plot accelaration i matlab with ODE45 function.

Here is the question:

I wrote the function like this:

function dxdt = vdp1(t,x)

dxdt = [x(2); (9.8 - 0.00324*sign(x(2))*x(2)^2-(1/7)*(x(1)-150)*heaviside(x(1)-150))];
end

For part 1 and 2,to get results I wrote this code:

%Part 1

[t,x] = ode45(@vdp1,[0 11.47], [0;0])
plot(t,x(:,1),'-o',t,x(:,2),'-o')
title('Solution with ODE45');
xlabel('Time t');
ylabel('Solution x');
 legend('x_1','x_2')

% Part 2
[t,x] = ode45(@vdp1,[0 5.988], [0;0])
plot(t,x(:,1),'-o',t,x(:,2),'-o')
title('Solution with ODE45');
xlabel('Time t');
ylabel('Solution x');
% legend('x_1','x_2')

[t,x] = ode45(@vdp1,[0 11.47], [0;0])
plot(t,x(:,1),'-o',t,x(:,2),'-o')
title('Solution with ODE45');
xlabel('Time t');
ylabel('Solution x');
legend('x_1','x_2')

Now, I don't know how to find accelaration with with this way. Please can anyone help me. Thanks in advance.

Answer 1

[James Tursa]

You find acceleration by taking the x solution from ode45( ) and feeding that back into your derivative function vdp1( ). Since your derivative function is not vectorized, you will have to do this in a loop.

[Elnur Shahbalayev]

Thank you so much for your answer. But I am still a bit confused. If possible, please can you show with code?

[James Tursa]

E.g., something like this:

[t,x] = ode45(@vdp1,[0 11.47], [0;0]); % solve the ode
n = size(x,1); % figure out how many rows are in the solution
xdot = zeros(size(x)); % allocate the derivative matrix
for k=1:n % loop through each solution vector to fill in the derivative matrix
    xdot(k,:) = vdp1(t(k),x(k,:)); % calculate the derivative at time t(k)
end
acceleration = xdot(:,2); % pick off the acceleration values

[Elnur Shahbalayev]

Thank you so much for your help!