Describing the motion of a composite body using system of second order differential equations

I am trying to graph the solution to the following system of differential equations.

t''= t'^2(tan(t)) - x''(4sec(t))

t'^2 (sec(t)-tan(t)) = x''*sec(t)*(1/0.625-4)

Initial conditions: t(0) = pi/3; t 0) = 0; x'(0) = x(0) = 0

The function file contains the following code.

function xy=PackageMotion(t,x)
% the differential equation soltuion 
xy=zeros(4,1);


xy(1)=x(2);   xy(3)=x(4);                               % The position x(t)
xy(2)=x(2)*x(2)*tan(x(1))-xy(4)*4*sec(x(1));         % The velcoity along x axis xdot(t)
                                % The position y(t)
xy(4)=(x(2)*x(2)*(sec(x(1))-tan(x(1))))/(1/0.625-4); 
    

end

And the code to graph the solution is as follows.

tspan=[0 10];     % tspan=[t0  tf]

xy0=[pi/3 0 0 0];

[t,y_sol]=ode45(@PackageMotion,tspan,xy0);


figure(),
plot(y_sol(:,1),y_sol(:,3));grid on
xlabel('Time (sc)');
ylabel('x(t) vs y(t)');
xlim([0 length(t)])

The plot is a either a vertical line or not visible. Can someone explain what is wrong about the code? Also, is it possible to obtain a closed form solution of the system of equations so I can understand what is being plotted?

Answer 1

[James Tursa]

In this line:

xy(2)=x(2)*x(2)*tan(x(1))-xy(4)*4*sec(x(1)); 

you are using xy(4) before it is defined. That is, you have t'' defined in terms of another highest order term x'' that has no value assigned to it yet. You can't do that.

To use ode45( ) on this, you first need to solve your equations for the highest order derivatives. I.e., you need to have equations that look like this:

t'' = stuff that is known

x'' = stuff that is known

Then code up your derivative function from these equations.

It looks like you could solve that 2nd equation for x'' first, evaluate it, and then use the x'' result to plug into the t'' equation.