10: Parachute Person

Parachute person is an example of a problem relevant to engineers and scientists that can be made fairly simple for beginning programmers. This problem can be analytically solved using differential equations, however, there are many problems that cannot be solved analytically but can be solved on a computer using the method (a simple numerical method which in this case makes an approximation of a differential equation) that will be presented here.

So why use this example if it can be solved analytically? Because it is simpler than an intractable problem and thus will be easier to grasp for a beginning programmer. Will real world problems be as simple? No. This example can be expanded to intractable real world problems, but a good understanding of numerical methods is necessary to truly do this sort of work (therefore it is wise for a student to take a numerical methods course or a signals and systems course that includes numerical methods).

First we discuss the physics and show how to design a computer program from a "word" problem. Then we show the working program and the results we can get.

A person has jumped out of a plane with a parachute, what is the velocity as the person reaches the ground?

As with all programming an initial program is made and then new versions enhance the initial program. Assumptions can systematically be removed or refined for new versions, but for now we will do the initial program. Because this is such a general problem the instructor has obviously :) given us liberty to make some reasonable assumptions. So let us list some assumptions.

• The ground is sufficiently far enough that the parachutist doesn't hit the ground before reaching terminal velocity
• The parachutist opens the parachute as soon as he jumps out of the plane (this of course would be one of the first assumptions to be revisited)
• The drag constant will include drag coefficient, parachute area size, design of parachute, and density of the air (more assumptions that would eventually be revisited)
• Weight of parachute is negligible
• Wind will be assumed to be zero in all directions (this assumption might be difficult to revisit, but give it a shot when you are more expert at programming)

So this is a physics problem at its heart. Since the parachutist is moving we recognize this as a dynamic problem.

• First step in a dynamic physics problem is to set up your force diagram.
• Before you draw your force diagram, draw your coordinate directions (which way does "x" go....which way does "y" go)
• Use the parachutist below to draw your forces. ¬ Parachute dude's force diagram (with coordinate system) Good, so we now see that the drag force ($$F_d$$) is negative (because of the coordinate system) and the gravity force ($$F_g$$) is positive. Now what do the forces equal? Assume a generic constant for the drag force.

¬ Forces defined

${\vec{F_g}} = \frac{d{\vec{p}}} {dt}$

Assuming m is constant we can rewrite this as...

${\vec{F_g}} = m \vec{a} = mg$

The simplest version of the drag force is this...

${\vec{F_d}} = c v^2$

Where c represents a constant that incorporates the drag coefficient, cross sectional area of the parachute, and the density of the atmosphere...

¬ Force equation

The generic static rule is that the sum of the forces in one direction equal zero (in addition to other equations introduced later)

$\sum{F_i} = 0$

However this problem is not static, it is dynamic, so...

$\sum{F_i} = ma$

Given our previously defined forces with their directions we have as a final equation

$F_g - F_d = mg - c v^2 = ma$

Using a differential for the acceleration we will develop a simple differential equation

$mg - c v^2 = m \frac{dv}{dt}$

While we could solve this differential equation using standard ideas introduced in differential equations (which most of you haven't taken), that is not the point of this exercise (whew!), so the next step is to produce an equation suitable for incorporating into a computer program.

¬ Force equation redefined...this is where it moves from physics to computer science

Starting with the differential equation previously derived

$m \frac{dv}{dt} = mg - c v^2$

and the definition of a differential

$\frac{dv}{dt} = \lim_{t \to 0} \frac{\Delta v}{\Delta t}$

we can write an approximate "algebraic" equation

$\frac{\Delta v}{\Delta t} = g - \frac{c}{m} v^2$

which we will rewrite as

${\Delta v} = (g - \frac{c}{m} v^2) {\Delta t}$

This is the form we will make into a usable computer form by using the definition of what the delta means in mathematics...which is

${\Delta v} = v_{i+1} - v_{i}$

Substituting into the previous equation we have

$v_{i+1} - v_{i} = (g - \frac{c}{m} (v_{i})^2) {\Delta t}$

Which we then rewrite in a form that is recognized as Euler's method in numerical methods...

$v_{i+1} = v_{i} + (g - \frac{c}{m} (v_{i})^2) {\Delta t}$

A vector (array) form is useful when writing a computer program as subscripts are not recognized by computer languages

$v(i+1) = v(i) + (g - c*(v(i))^2/m) * tsteps$

Now we are ready to program this...yes!!!

¬ Vector array, say what?

The vector array as noted previously is a variable that has many cells. Here is an example where t and v are vector arrays. This means that t and v are actually many variables as one. If we define v (in Fortran, say) as -- double, dimension(10) :: v -- then v is actually 10 distinct variables having the names v(1), v(2), v(3), v(4), v(5), v(6), v(7), v(8), v(9), and v(10). Each of these variables can be assigned a number as the example below shows.

 i is an index variable t(i) is time (seconds) v(i) is velocity 1 0.0 0.0 2 0.1 0.98 3 0.2 1.93253 4 0.3 2.80572 5 0.4 3.56058 6 0.5 4.17800 7 0.6 4.65876 8 0.7 5.01803 ... ... ...

The previous example fits with the final program below...

• Here is one version of parachute person
• The students should try to write this in a number of other languages (Fortran, C++, IDL, and Python)
• The student should try to improve upon this program (a number of real world concepts have been sacrificed here for simplicity)
%
% Function parachute person
% Syntax: [time,velocity] = paraperson(0,50,0,57.2,200,9.8,0.1)
% Author: Scott D. Johnson as an example to understanding MATLAB/octave
%
% Inputs:
%
% ti is initial time
% tf is final time
% vi is initial velocity
% c is coefficient of drag of the person and parachute combined (this is an equation in itself)
% m is mass of person and parachute
% g is gravity
% tsteps is the size of the time increments
%
% Outputs:
%
% t is a time vector
% v is a velocity vector (terminal velocity)
%
% Note:  This is an extremely simplified physics of skydiving problem...why?  So the student can program
%            the more complex version.
%
%-----------------------------------------------------------------------------------------------------------------------------------------
%
function [t,v]=paraperson(ti,tf,vi,c,m,g,tsteps)
numsteps = (tf-ti)/tsteps;
i=1;
v(i)=vi;
t(i)=ti;
for i = 1:numsteps-1
v(i+1)=v(i) + [g-c*v(i)^2/m]*tsteps;
t(i+1)=t(i) + tsteps;
end %for
return
end%
%
% This is the run of the function
%
[time,vel]=paraperson(0,50,0,57.2,200,9.8,0.1);
%
% Plotting the data
%
plot(time,vel,';Parachute Person;');
xlabel('Time in seconds');
ylabel('Velocity in m/s');
title('Velocity of person with a parachute jumping out of a plane')
axis([0,20,0,8])