This notebook is due on Friday, December 9th, 2016 at 11:59 p.m.. Please make sure to get started early, and come by the instructors' office hours if you have any questions. Office hours and locations can be found in the course syllabus. IMPORTANT: While it's fine if you talk to other people in class about this homework - and in fact we encourage it! - you are responsible for creating the solutions for this homework on your own, and each student must submit their own homework assignment.
FOR THIS HOMEWORK: In addition to the correctness of your answers, we will be grading you on:
To that end:
Put your name here!
Emboldened by his success in parachuting from a helicopter that's hovering in place, your professor has decided to start jumping from actual moving airplanes. In this section, we're going to model the behavior of a skydiver that is jumping from a moving plane and opening their parachute at some point above the ground. Your professor has requested your help in trying to figure out where he is going to land, knowing some basic information about how high above the ground he is and how fast he's going when he jumps out of the plane.
An object's position ($\vec{x}$) and velocity ($\vec{v}$) evolve according to a straightforward set of mathematical rules:
$\vec{x} = \vec{x}_0 + \int_0^t \vec{v} dt$
$\vec{v} = \vec{v}_0 + \int_0^t \vec{a} dt$
$\vec{a} = \frac{d\vec{v}}{dt} = \frac{1}{m}\sum \vec{F}$
Note that $\vec{a}$ is the acceleration, m is the mass of the skydiver, and $\sum \vec{F}$ are the sum of all of the forces acting on the sky diver. $\vec{x}_0$ and $\vec{v}_0$ are the position and velocity at $t=0$, respectively. The little arrows on top of letters indicate vector quantities - i.e., quantities that have direction as well as magnitude. With vector quantities, you can separate the components and solve them independently - for example, you only need to know the x-velocity to calculate the x-position (though ultimately you may need some quantities that need information from both - you need the speed, which includes x- and y-velocity, to calculate the acceleration due to air resistance). In the last term of the third equation, the sum is over all forces acting on the object (i.e., the sky diver).
The forces acting on the object are gravity, which points downward with a magnitude of $g=9.81$m s$^{-2}$, and air resistance, which can be written as:
$\vec{F}_{air} = -\frac{1}{2} C \rho A v^2 \hat{v}$
Where C is a constant corresponding to the shape of the object ($C=0.47$ for a sphere, $C=2.0$ for a parachute), $\rho$ is the density of the atmosphere, 1.2 kg/m$^3$, A is the cross-sectional area of the object, v is the magnitude of the velocity, and $\hat{v}$ is the unit vector of the velocity. Note the minus sign - this means that the force of air resistance always opposes the motion of the object, regardless of its direction.
Your task is to create a code that implements the equations presented above to model the position and velocity of the sky diver as a function of time, assuming the sky diver jumps out of the plane at some height above the ground Y and traveling at some horizontal speed V, and that the parachute opens at some time T. You have been provided a function that calculates the acceleration on the sky diver due to both air resistance and gravity; you need to write a function that takes in the height at which sky diver jumps out of the plane, his initial horizontal speed, and the time that the parachute opens, and returns the sky diver's position and velocity over time (so position, velocity, and time). Assume that the sky diver jumps out of the plane 2,000 meters above the ground while traveling at 50 m/s in the horizontal direction, and pulls the rip cord when he is 500 meters above the ground.
After you create the model, plot the sky diver's horizontal and vertical position and velocity as a function of time. Assuming that he jumps out of the plane at $t=0$ and $x=0$, how far does he land from the place he jumped out of the plane (in terms of horizontal distance), and at what time does he reach the ground?
PUT YOUR ANSWSERS TO THE QUESTIONS HERE! (And put your plots below, after your model)
In [ ]:
%matplotlib inline
import matplotlib.pyplot as plt
import numpy as np
import math
def calc_accel(velocity, parachute_open=False, use_resistance=True):
'''
This function calculates and returns the acceleration vector
for a sky diver, given an input velocity. You can indicate whether the
parachute is open or not (which affects the air resistance). Optionally,
you can turn on and off air resistance for comparison purposes with
resistance=True or False.
Inputs are velocity, a Boolean argument indicating if the parachute is open
(default False, so the parachute is assumed to be closed), and a Boolean argument
indicating if air resistance is to be used.
Output is a 2-element numpy array that indicates velocity.
Note 1: inputs and outputs are in SI units - meters per second for
velocity and meters per second^2 for acceleration.
Note 2: both input velocity and returned acceleration are 2-element
numpy arrays, with the horizontal (x) direction in the first element
and vertical (y) direction in the second element of the array.
'''
C_diver = 0.47 # skydiver, parachute closed
C_chute = 1.0 # skydiver, parachute open
rho_air = 1.225 # kg/m^3
mass_obj = 100.0 # mass of skydiver in kg, including all gear
# speed = magnitude of velocity vector
speed = (velocity[0]**2 + velocity[1]**2)**0.5
accel_vector = np.zeros(2)
# if we have air resistance turned on, compute the force due to air
# resistance and add it to the acceleration.
if use_resistance==True:
'''
Depending on whether the parachute is open or not,
that radically changes our air resistance! Check to
see if parachute's open and set parameters accordingly.
'''
if parachute_open==True:
rad_obj = 3.0
C_obj = C_chute
else:
rad_obj = 0.75
C_obj = C_diver
'''
Calculates acceleration due to air resistance
|F_air| = 0.5 * C * rho * A * v^2
C = constant based on shape
rho = density of air
A = cross-sectional area (depends on whether parachute is open or closed)
v = speed (magnitude of velocity)
Note that we calculate the magnitude in the first step and turn it into a vector in the second step.
'''
accel_air_resist = -0.5 * (C_obj/mass_obj) * rho_air * math.pi * rad_obj**2 * speed**2
# turn this into a vector quantity - velocity[0]/speed gives x-component of velocity unit vector
# and velocity[1]/speed gives y-component of velocity unit vector
accel_vector[0] += accel_air_resist*velocity[0]/speed
accel_vector[1] += accel_air_resist*velocity[1]/speed
# quick check to make sure that the acceleration magnitude isn't too ridiculous
# (this is for numerical stability, and is a slight hack. Forgive me!)
accel_mag = (accel_vector[0]**2 + accel_vector[1]**2)**0.5
max_accel = 40.0
if accel_mag > max_accel:
accel_vector[0] *= max_accel / accel_mag
accel_vector[1] *= max_accel / accel_mag
# we always have gravity, which always points downward.
accel_vector[1] -= 9.81
return accel_vector
In [ ]:
# Put your model here. Add additional cells if necessary!
In [ ]:
# Put the code that generates your plots here. Add additional cells if necessary!
In this part of the homework, we're going to extend our fire model to make it more realistic. This model, like all models, makes a variety of simplifying assumptions. One that is particularly important is that there are many factors that affect the spread of fire in a forest, such as wind, the type of wood, the amount of dead underbrush, and even how close the branches of adjoining trees are to each other.
Given that there are so many possible complications, we're going to make a small modification of our model to take all of these uncertainties into account. We know that it's not certain that fire will spread, but we don't know exactly what that probability is. So, we're going to make this a new parameter for our system, called prob_spread. By experimenting with values of prob_spread (which is a fraction that should vary between 0 and 1), we can see what effect the probability of the fire spreading from cell to cell has on the overall spread of the fire.
The plan: Modify the fire model you developed in class as follows. In the original model, if a cell has trees and an adjacent cell is on fire, that cell automatically catches on fire in the next turn. Now, instead of it automatically catchin on fire, you generate a random number and compare it to prob_spread. If the number is smaller than prob_spread, the trees catch on fire.
Question: In the original model, you see a rapid transition in the amount of the forest that catches on fire when more than approximately 59% of cells are occupied by trees (this is referred to as the "tipping point" or "critical threshold"). As you try different values of prob_spread (say, values of 0.1 to 1.0 by steps of 0.1), how does the behavior differ from the standard model? At what value of prob_spread is it essentially indistinguishable from the standard model? Answer the question below, and support your answer with plots!
Put your answer here!
In [ ]:
# Put your model here. Add additional cells if necessary.
In [ ]:
# Put the code that generates your plots here. Add additional cells if necessary.
In [ ]:
from IPython.display import HTML
HTML(
"""
<iframe
src="https://goo.gl/forms/OU4CVKhXdrdTGVVx1?embedded=true"
width="80%"
height="1200px"
frameborder="0"
marginheight="0"
marginwidth="0">
Loading...
</iframe>
"""
)