3d plot of Tridimensional parabolic motion

Hi, I am trying to get a 3D plot of a parabolic motion. I am trying to simulate debris from an airplane falling from a certain altitude. Obviously the parabolic motion and the absence of other forces is a simplification.

These are the initial conditions.

h_0=14325 m

V_0= 200*.3048 m/s

gamma = 0.78 rad

Obviously I want my object to reach the ground. I created the following code, but the plot is a little strange. I don't understand where or what is wrong

g=9.8;
z0= 14325;
x0=0;
y0=0;
V_0= 200*0.3048;
gamma=0.78;
V_0x=V_0*cos(gamma);
V_0y=V_0*cos(gamma);
V_0z=V_0*sin(gamma);

t = linspace(0,2*V_0*sin(gamma)/g);
x=x0+V_0x*t;
y=y0+V_0y*t-.5*g*t.^2;
z=z0+V_0z-.5*g*t.^2;

indices=z>=0;
x = x(indices);
y = y(indices);
z = z(indices);

plot3(x,y,z)
hold on
plot3(x0,y0,z0,'ro','MarkerSize',5)
xlabel('Posizione x (m)')
ylabel('Posizione y (m)')
zlabel('Posizione z (m)')
grid on

Answer 1

[John D'Errico]

It should look strange. Why?

g=9.8;
z0= 14325;
x0=0;
y0=0;
V_0= 200*0.3048;
gamma=0.78;
V_0x=V_0*cos(gamma);
V_0y=V_0*cos(gamma);
V_0z=V_0*sin(gamma);
t = linspace(0,2*V_0*sin(gamma)/g);
x=x0+V_0x*t;
y=y0+V_0y*t-.5*g*t.^2;
z=z0+V_0z-.5*g*t.^2;
indices=z>=0;
x = x(indices);
y = y(indices);
z = z(indices);
plot3(x,y,z)
hold on
plot3(x0,y0,z0,'ro','MarkerSize',5)
xlabel('Posizione x (m)')
ylabel('Posizione y (m)')
zlabel('Posizione z (m)')
grid on
hold off

Your equations of motion are wrong. Let me explain...

First, you appear to have an angle gamma. Should we presume that gamma is the angle above the horizontal of the initial velocity?

rad2deg(gamma)

ans = 44.6907

So a roughly 45 degree angle above the horizontal.

And then you split that into a vertical component, AND a horizontal component. So we saw this:

V_0x=V_0*cos(gamma);
V_0y=V_0*cos(gamma);
V_0z=V_0*sin(gamma);

But what you did was set the SAME velocity in both the x and y directions. Essentially, that puts too much velocity into the horizontal component, That is wrong, because now the initial velocity is

norm([V_0x,V_0y,V_0z])

ans = 74.7947

But V_0 is

V_0

V_0 = 60.9600

Do you understand the problem? You needed to split the horizontal velocity into two components properly, NOT match the SAME velocities in both x and y.

Ok, but that was a relatively minor problem. Your initial velocity was wrong. What else is wrong?

x=x0+V_0x*t;
y=y0+V_0y*t-.5*g*t.^2;
z=z0+V_0z-.5*g*t.^2;

Do you see what you did?

Next, you have the z equations wrong. You dropped a t in there.

But worse, why is there an acceleration due to gravity in the y direction?

Finally, is there a reason why you only went out as far in time as you did?

2*V_0*sin(gamma)/g

That is, is that really the amount of time necessary for the projectile to hit the ground, having started up a long way above the ground?

I suppose you can make up equations of motion that will work however you want. But then don't expect them to look as if this projectile is moving in a planetary gravitational field. :)

What would I change to fix it?

First, I would decide how to split up the HORIZONTAL component of velocity into two pieces. So what directino is the projective moving between x and y?

g=9.8;
z0= 14325;
x0=0;
y0=0;
V_0= 200*0.3048;
phi = .5236; % the angle between x and y for the initial components of velocity
gamma=0.78;
V_0x=V_0*cos(gamma)*cos(phi);
V_0y=V_0*cos(gamma)*sin(phi);
V_0z=V_0*sin(gamma);

First, is this the correct overall velocity?

norm([V_0x,V_0y,V_0z])

ans = 60.9600

Now V_0 is consistent.

Next, when will the projectile hit the ground? How long? If the z-equation is correct, then this will do:

troots = roots([-.5*g, V_0z, z0])

troots = 2×1

58.6205 -49.8711

Taking the positive root, we can do this:

t = linspace(0,troots(troots > 0),1000); 
x=x0+V_0x*t;
y=y0+V_0y*t; % no acceleration in y due to gravity
z=z0+V_0z*t-.5*g*t.^2; % see the difference!

Did it work?

z(end)

ans = 0

% You don't need to worry about negative z anymore.
% indices=z>=0;
% x = x(indices);
% y = y(indices);
% z = z(indices);
plot3(x,y,z)
hold on
plot3(x0,y0,z0,'ro','MarkerSize',5)
xlabel('Posizione x (m)')
ylabel('Posizione y (m)')
zlabel('Posizione z (m)')
grid on

Now you could add in air resistance too. Make things interesting.

[Giovanni Curiazio]

thank you, yes please explain

[Giovanni Curiazio]

OK, I completely understood the errors. But how do I achieve a final altitude of 0? Thank you