Matlab loops help (with physics)

A skydiver with a mass of 80 kg (including parachute) jumps out of an airplane at an altitude of 4000 m, the initial vertical velocity is 0 m/s. He free-falls for 40 seconds before opening a parachute. The parachute has a 10 m diameter when completely opened (you can assume it opens instantaneously). Create a program in MATLAB that calculates and plots the skydiver’s vertical speed and altitude as a function of time using the following equations and values. Remember that the skydiver cannot exceed the terminal velocity and the parachute has a failsafe to open at 1000 m.

Hints: There are different terminal velocities with and without the parachute open, assume that the parachute opens instantaneously.

Constant values g = 9.81 [m/s^2 ] (Gravitational Acceleration) p_air = 1.225 [kg/m^3 ] (Density of air)

Before Opening: A = 0.5 [m^2 ] Cd = 0.7

After Opening: A = pi/4 * D^2 [m2 ] Cd = 1.4

Calculation of Terminal Velocity

Vt=0.99*sqrt((2*mass*gravity)/(p_air*A*Cd))

Calculation of Downward Acceleration (derived from F=m*a)

ai=(1/2)*(p_air/mass)*(V^2_i-1)*(A)*(Cd)-g

Kinematics Equations Setup for a loop

Vi=Vi-1+ai*(delta t) t is time

Yi=Yi-1+Vi-1*(delta t)+(1/2)*ai*(delta t)

Setup the diameter of the parachute, initial altitude, time of opening the parachute, and the time step as user inputs.

In addition to the plots, output the answers to the following questions the command window: 1. Is terminal velocity reached before the skydiver opens the parachute? If so, at what time did the skydiver reach terminal velocity? 2. At what altitude did the parachute open? 3. Does the skydiver open the parachute or did the failsafe mechanism open it?

Use While-end loops to solve for the equations, must use these

[James Tursa]

The general procedure on this forum is for you to make a legitimate attempt at solving the problem. If you have problems, then post your code and ask specific questions about it and we can help.

[Gerald Melins]

all I figured out was declaring the variables and what goes inside the while loop.

for the loop while (time > 40 && altitude > 1000)

end

I'm sorry I don't know that much, im completely new to matlab and haven't done physics for 2 years

[James Tursa]

So start with something simple. Can you code up a simulation that starts at 4000 m and 0 m/s velocity and does a free fall with constant acceleration down of 9.81 m/s^2 for 10 s with no other forces involved? I.e, ignoring parachute and atmospheric drag. E.g., just code up the following two equations in m-code and run it for 10 s:

Vi=Vi-1+ai*(delta t) t is time

Yi=Yi-1+Vi-1*(delta t)+(1/2)*ai*(delta t)

In the above, the Vi-1 is interpreted as the V from the previous time point, and Yi-1 the same. I.e., the equations would look something like this in MATLAB

V(i) = V(i-1) + a(i) * (delta_t)

Y(i) = Y(i-1) + V(i-1) * (delta_t) + (1/2) * a(i)*(delta_t)

Once you get that working then you can start to add in the other effects such as parachute and drag. As a check to see if you have this part working, you can compare your result to the exact solution:

The change in Y at future time t (assuming you start at t=0) should be (1/2)*a*t^2, where a = constant 9.81 for our simple test case. Since it is falling down, you would actually use a = -9.81 in your simulation and Y should decrease by the amount shown.

[Gerald Melins]

posted code below

[Gerald Melins]

I also got this error Attempted to access V(0); index must be a positive integer or logical.

[Image Analyst]

Gerald, see the FAQ

Answer 1

[James Tursa]

OK, now we have something to work with. Thanks for making the coding attempt and posting it.

My suggestion is to first simplify this to a test case that you can compare with the exact solution for the constant acceleration case. To do that, get rid of all the complicated parachute & drag stuff and while loop conditions and do the simple constant acceleration free fall (copy your code to a new file name first) for a fixed number of steps. E.g., simply extracting parts of your above code for a simple free-fall simulation:

% Constants
g=-9.81; % m/s^2

% Initial conditions
a(1) = g; % m/s^2
V(1) = 0; % m/s
Y(1) = 4000; % m

% Simulation step size
delta_t = 0.25; % s

% Initial time
  t = 0; % s

% Simulation
for i=2:41 % start at i=2, not i=1
    a(i) = -g; % Constant down acceleration
    V(i) = V(i-1) + a(i) * delta_t; % Update velocity
    Y(i) = Y(i-1) + V(i-1) * delta_t + (1/2) * a(i) * delta_t; % Update position
t = t + delta_t; % Update time
end

Now just compare the change in position to the exact solution:

Y(end) - Y(1) is the total change in position

(1/2) * g * t^2 is the exact solution

First, get this working and understand it.

Then you can work on adding in the other stuff. E.g., replace the for loop with a while loop with appropriate test condition(s). And you will be using a different equation for a(i) at each step depending on the current states of the parachute and velocity etc. In other words, you will likely have an "if" check somewhere near the top of the loop to determine how to calculate a(i). So you will need to figure out how to code up that "if" test. And how to account for the terminal velocity stuff. Etc.

[Gerald Melins]

I'm supposed to use a while loop instead of for. So I don't think I will have to use an if test if I have the while loop checking that hes above 1000 feet and after 40 seconds have passed, right?

[James Tursa]

Yes, in your final code you will have a while loop per instructions. But you don't understand how to code it yet which is why I suggest you start with something simple first. E.g., here is what you have above:

while (Ui > 0 && altitude >= 1000)

But Ui is not defined at all, and altitude is not changing during the loop. These are both problems that will prevent the code from running properly. Ui I suspect is a typo and is probably meant to be Y(i), and altitude is probably just Y(i). So these errors are probably not going to be difficult to fix. But rather than address that up front I thought it best to get you started on a simple, related problem first.

[Gerald Melins]

ok thanks, I ran a few tests and think I'm doing better now. I changed the while loop to while (t > 0 && Y(i) > 1000) so it checks the conditions the parachute opens. Then I did if and else statements to check whether the parachute opened because the time ran out or because he reached the failsafe height.

I made a comparison of the terminal velocity, the one from inside the loop and another once the parachute opens. So the greater value would be the terminal velocity. I was able to print out the final altitude and the else/if statement said the reason for the parachute opening

[Gerald Melins]

g=-9.81; % m/s^2
airdens=1.225; % g/m^3
mass=80; %kg
diameter=10; %m

%initial equations
a(1) = g; % m/s^2
V(1) = 0; % m/s
Y(1) = 4000; % m
delta_t = 1; % s

% Initial time
  t = 0; %s

while (t > 0 && Y(i) > 1000) %starting loop to check time and altitude, parachute is not open yet

a(i) = (1/2)*(airdens/mass)*V(i)^2*0.5*0.7 -g; % Constant downward acceleration
V(i) = V(i-1) + a(i) * delta_t; % checking velocity
Y(i) = Y(i-1) + V(i-1) * delta_t + (1/2) *a(i) * delta_t; % checking position
t = t + delta_t; % checking time
end

Vterminalopen=0.99*sqrt((2*mass*g)/(airdens*0.5*0.7)); %terminal velocity, parachute unopened
Aopen=(pi/4)*diameter^2; % m^2
Vterminalclose=0.99*sqrt((2*mass*g)/(airdens*Aopen*1.4)); %terminal velocity, parachute opened

if Vterminalopen < Vterminalclose
    fprintf('The diver reached terminal velocity before the parachute opened \n');
    fprintf('The terminal velocity is %.2f \n',Vterminalclosed);
else
    fprintf('The diver reached terminal velocity after the parachute opened \n');
    fprintf('The terminal velocity is %.2f \n',Vterminalopen);
end

if(Y<=1000) %condition the parachute would have opened for, failsafe height
fprintf('The parachute opened once the skydiver reached 1000 meters \n');
disp(Y(i));
disp(t(i));

elseif t <=0 %another condition the parachute opens for, 40 seconds have passed
fprintf('The parachute opened after the skydiver was freefalling for 40 seconds \n');
disp(Y(i));
disp(t(i));
end

not sure what I am doing wrong

[James Tursa]

For starters, there is this

    t = 0; %s
  while (t > 0

The loop will never run because t does not start with a positive value.

But why are you checking t at all in your while loop check? You shouldn't, since it doesn't matter what t is. Get rid of the t > 0 check.

Then there are these problems:

You don't initialize i = 1 before starting the loop.

You don't increment i inside the while loop.

You don't have all of the other stuff (did parachute open, what A to use for the a(i) calculation, etc) inside the while loop.

The pseudo-code for the loop should look something like this:

A = value for parachute closed
while( some test, e.g., is altitude positive? )
    if( parachute is not open )
        if( conditions are met for parachute opening )
            open the parachute
            A = value for an open parachute
        end
    end
    insert kinematics equations here
    insert terminal velocity checks here
end

Note that all of the kinematics, parachute stuff, and terminal velocity stuff are all inside the loop. You may want to have flags in your code to keep track of parachute status, etc.

[Gerald Melins]

I need the time check in the while loop, the parachute deploys after 40 seconds

[James Tursa]

The time check for the parachute goes in this line:

if( conditions are met for parachute opening )