Problem 3. No instantaneous transitions...
The acceleration of over 200 ft/sec^2 (over 6 gs) could cause a dangerous (if not fatal) force on the sky diver. Of course, this doesn't really occur since in reality there is a short (but measurable) time interval over which the transitions from vertical to spread eagle flight and especially from spread eagle to open chute parts of flight. We account for this below by adding two more time intervals into the problem: one for each transition.
We'll make the easiest assumption: over the transition period, the k's vary linearly.
Time Interval 1. From plane to the transition to spread eagle.
Input :=
v0 = 0;
temp = Flatten[DSolve[{de1, v[0] == v0}, v[t], t]]
Output =
278.8
{v[t] -> 278.8 - -----------}
0.114778 t
E
Input :=
vel1[t_] = v[t] /. temp
Output =
278.8
278.8 - -----------
0.114778 t
E
Integrate to get the position function, pos1(t).
Input :=
p01 = .;
pos1[t_] = Integrate[vel1[t], t] + p01
Output =
2429.04
----------- + p01 + 278.8 t
0.114778 t
E
Since the velocity is always positive, pos1[t] is increasing, so any guess for Newton's method will work: we guess 0. (We assume that the jump height, 0, corresponds to 0 feet from the plane and 5000 feet above ground.)
Input :=
FindRoot[pos1[0], {p01, 0}]
Output =
{p01 -> -2429.04}
Input :=
p01 = p01 /. %
Output =
-2429.04
Next, we find the time when pos1[t] = 2000.
Input :=
Plot[pos1[t], {t, 0, 10}]
Output =
-Graphics-
Input :=
FindRoot[pos1[t] == 2000, {t, 10}]
Output =
{t -> 14.1736}
Input :=
t1 = t /. %
Output =
14.1736
Time Interval 2: Assume that this is t2=1 second long.
Then, k(t) = ((k2 - k1)/t2) t + k1
Input :=
t2 = 1;
ktrans1[t_] = (k2 - k1)/t2 t + k1
Output =
0.555954 + 0.303248 t
Input :=
de2 = m v'[t] == m g - ktrans1[t] v[t];
Input :=
v1 = vel1[t1];
temp = Flatten[DSolve[{de2, v[0] == v1}, v[t], t]]
Output =
2
-0.114778 t - 0.031303 t
{v[t] -> 169.281 E +
2
-0.114778 t - 0.031303 t
81.4018 E Sqrt[Pi]
Erfi[0.176927 (1.83333 + t)]}
Input :=
vel2[t_] = v[t] /. temp
Output =
2
-0.114778 t - 0.031303 t
169.281 E +
2
-0.114778 t - 0.031303 t
81.4018 E Sqrt[Pi]
Erfi[0.176927 (1.83333 + t)]
Integrate to get the position function, pos2(t).
Input :=
pos2[t_] := NIntegrate[vel2[s], {s, 0, t}] + 2000;
Input :=
Plot[pos2[t], {t, 0, t2}]
Output =
-Graphics-
Time Interval 3: Falling "spread eagle" .
Determine the new position at the end of the first transition period
Input :=
p03 = pos2[t2]
Output =
2224.75
Input :=
v3 = vel2[t2]
de3 = m v'[t] == m g - k2 v[t];
temp = Flatten[DSolve[{de3, v[0] == v3}, v[t], t]]
Output =
146.273 + 43.3864 Sqrt[Pi]
Output =
42.7739
{v[t] -> 180.4 + -----------}
0.177384 t
E
Input :=
vel3[t_] = v[t] /. temp
Output =
42.7739
180.4 + -----------
0.177384 t
E
Input :=
pos3[t_] := Integrate[vel3[s], {s, 0, t}] + p03;
Input :=
FindRoot[pos3[t] == 4500, {t, 12}]
Output =
{t -> 11.4509}
Input :=
t3 = t /. %
Output =
11.4509
Time Interval 4: 2nd transition period.
Assume that this is t4 = 2 seconds long
(between spread eagle and open chute).
Then, k(t) = ((k3 - k2)/t4) t + k2
Input :=
t4 = 2;
ktrans2[t_] = (k3 - k2)/t4 t + k2
Output =
0.859202 + 3.20548 t
Input :=
de4 = m v'[t] == m g - ktrans2[t] v[t];
Input :=
v3 = vel3[t3];
temp = Flatten[DSolve[{de4, v[0] == v3}, v[t], t]]
Output =
{v[t] -> 177.5682906621567
2.
-0.1773835920177384 t - 0.33088862357155 t
E +
49.300803920419 Power[E,
-0.02377305872402679 - 0.1773835920177384 t -
2.
0.33088862357155 t ]
Erfi[0.575229192210853 (0.2680412371134021 + t)]
}
Input :=
vel4[t_] = v[t] /. temp
Output =
177.5682906621567 E
2.
-0.1773835920177384 t - 0.33088862357155 t
+
49.300803920419 Power[E,
-0.02377305872402679 - 0.1773835920177384 t -
2.
0.33088862357155 t ]
Erfi[0.575229192210853 (0.2680412371134021 + t)]
Integrate to get the position function, pos4(t).
Input :=
pos4[t_] := NIntegrate[vel4[s], {s, 0, t}] + 4500;
Input :=
Plot[pos4[t], {t, 0, t4}]
Output =
-Graphics-
Time Interval 5: Falling with chute open.
Determine the new position at the end of the second transition period:
Input :=
p05 = pos4[t4]
Output =
4763.36
Input :=
v5 = vel4[t4]
de5 = m v'[t] == m g - k3 v[t];
temp = Flatten[DSolve[{de5, v[0] == v5}, v[t], t]]
Output =
-15
59.974578681858 + 0. 10 I
Output =
38.65457868186
{v[t] -> 21.32 + ------------------}
1.5009380863039 t
E
Hmmm. Why does the DSolver give 0. 10^(-15) I? An interesting aside for some!
Input :=
vel5[t_] = v[t] /. temp
Output =
38.65457868186
21.32 + ------------------
1.5009380863039 t
E
Input :=
pos5[t_] := Integrate[vel5[s], {s, 0, t}] + p05;
Input :=
FindRoot[pos5[t] == 5000, {t, 12}]
Output =
{t -> 9.89164}
Input :=
t5 = t /. %
Output =
9.89164
So, the Sky Diver will hit the ground in
Input :=
t1 + t2 + t3 + t4 + t5
Output =
38.5162
seconds with an average speed of
Input :=
5000 ft/(% sec) 3600 sec/(1 hr) 1 miles/(5280 ft)
Output =
88.5106 miles
-------------
hr
Below are the graphs of position, velocity and acceleration of the diver during her fall.
Input :=
pos[t_] := If[t < t1, pos1[t],
If[t < t1+t2, pos2[t-t1],
If[t < t1+t2+t3, pos3[t-t1-t2],
If[t < t1+t2+t3+t4, pos4[t-t1-t2-t3],
pos5[t-t1-t2-t3-t4]
]
]
]
];
vel[t_] := If[t < t1, vel1[t],
If[t < t1+t2, vel2[t-t1],
If[t < t1+t2+t3, vel3[t-t1-t2],
If[t < t1+t2+t3+t4, vel4[t-t1-t2-t3],
vel5[t-t1-t2-t3-t4]
]
]
]
];
Input :=
acc1[t_] = D[vel1[t], t];
acc2[t_] = D[vel2[t], t];
acc3[t_] = D[vel3[t], t];
acc4[t_] = D[vel4[t], t];
acc5[t_] = D[vel5[t], t];
acc[t_] := If[t < t1, acc1[t],
If[t < t1+t2, acc2[t-t1],
If[t < t1+t2+t3, acc3[t-t1-t2],
If[t < t1+t2+t3+t4, acc4[t-t1-t2-t3],
acc5[t-t1-t2-t3-t4]
]
]
]
];
Input :=
Plot[pos[t], {t, 0, t1 + t2 + t3 + t4 + t5},
PlotRange -> All,
AxesLabel -> {"t [sec]", "distance fallen [ft]"}]
Output =
-Graphics-
Input :=
Plot[vel[t], {t, 0, t1 + t2 + t3 + t4 + t5},
PlotRange -> All,
AxesLabel -> {"t [sec]", "speed [ft/sec]"}]
Output =
-Graphics-
Input :=
Plot[acc[t], {t, 0, t1 + t2 + t3 + t4 + t5},
PlotRange -> All,
AxesLabel -> {"t [sec]", "acceleration [ft/sec^2]"}]
Output =
-Graphics-
Note that now the maximum acceleration due to gravity is approximately 83 ft/sec^2 which is less than 3 g's. Of course, the solution is dependent on the length of the intervals. It's up to the solver to determine a reasonable length for each transition period.
It's also interesting to note that even though the diver had a faster flight according the second model, it was a far more comfortable one.