Skip to content
Snippets Groups Projects
Commit 32561d3c authored by cprutean's avatar cprutean
Browse files

Update file Checkpoint5.ipynb

parent f32b0044
No related branches found
No related tags found
No related merge requests found
%% Cell type:markdown id:abaa6e4a tags:
## Checkpoint 5
### Aim
To write two Python programs to numerically analyse the path of a projectile under gravity and drag force. There are two tasks, one to display the trajectory and the second to display the relation between final kinetic energy and launch angle.
### Numerical Integration
In this checkpoint we will use the simplest numerical integration technique called First Order Euler. This involves iteratively calculating the position and velocity at a set of times t$_i$=iΔt for i=0,1,2,… where we assume that the position and velocity at time t$_{i+1}$ depend only on the position and velocity at time t$_i$. In particular if at time t$_i$ the particle has
$$\text{Position} = \vec{x}_i$$
$$\text{Velocity} = \vec{v}_i$$
$$\text{Acceleration} = \vec{a}_i$$
then at time t$$_{i+1}$$=t$$_i$$+Δt the position and velocity of the particle can be approximated by
then at time $t_{i+1}=t_i+Δt$ the position and velocity of the particle can be approximated by
$$\vec{x}_{i+1} = \vec{x}_i+Δt\vec{v}_i$$
$$\vec{v}_{i+1} = \vec{v}_i+Δt\vec{a}_i$$
So if we start at time t$_0$=0 at position $\vec{x}_0$ and velocity $\vec{v}_0$, then provided we can calculate the force on the particle, and hence the acceleration $\vec{a}_i$, and use small enough time steps Δt, we can iteratively trace out the path by repeated use of the above relations.
So if we start at time $t_0=0$ at position $\vec{x}_0$ and velocity $\vec{v}_0$, then provided we can calculate the force on the particle, and hence the acceleration $\vec{a}_i$, and use small enough time steps Δt, we can iteratively trace out the path by repeated use of the above relations.
### Forces on Particle
![ch5.png](ch5.png)
If we consider a particle of mass m travelling with velocity $\vec{v}$ subject to gravity and drag from a fluid as shown below, then the drag force is given by
$$ \vec{F}_d=−\frac{1}{2}\rho AC_Dv^2\hat{v}$$
where $\rho$ is the density of the fluid, A is the cross-sectional area of the particle, C$_D$ is the drag coefficient and $\hat{v}=\vec{v}/|v|$, being the unit vector in the direction of the velocity.
Thus the acceleration of of the particle is given by $\vec{a}=\frac{1}{m}\vec{F}_d−g\hat{j}$ where g is acceleration due to gravity, being 9.81ms$^{−2}$. We can write this as:
$$ \vec{a} =−\beta v^2\hat{v}−g\hat{j}$$
where $\beta=\frac{\rho AC_D}{2m}$, where $\beta$ is the normalised drag coefficient for the system.
In this problem the position, velocity and acceleration are all two dimensional vectors so they have a horizontal and vertical component that you have to treat separately in the calculation.
So in terms of components in x (horizontal) and y (vertical) we can write:
$\vec{x}=(xx,xy); \vec{v}=(vx,vy); \vec{a}=(ax,ay)$
where the components of acceleration from the above equations are:
$a_x=−\beta vv_x$ and $a_y=−\beta vv_y−g$
where v is the magnitude of the velocity, so being $v=\sqrt{v^2_x+v^2_y}$
### Task
**Task 1:** Write a Python program to plot out the trajectory of a projectile starting at position $\vec{x}_0=(0,0)$ for specified initial velocity and normalised drag coefficient. You program should prompt for
- Magnitude of the initial velocity v$_0$ in ms$^{−1}$.
- Angle θ from the horizontal of initial velocity in degrees.
- β the normalised drag coefficient.
- Δt the step interval in seconds.
The path should be traced and displayed using Matplotlib from the inital starting position until the particle re-crosses the x-axis. Your graph should have a suitable title and axis labels. Your program should also print the range of the projectile to the terminal.
Test your program with:
- v$_0$=10 ms$^{−1}$, $\theta$=45$^{\circ}$, $\beta$=0.0, $\Delta$t=0.001
- v$_0$=10 ms$^{−1}$, $\theta$=60$^{\circ}$, $\beta$=0.3, $\Delta$t=0.001
**Task 2:** Modify the above program to produce a plot of the ratio of final kinetic energy to initial kinetic energy against lauch angle θ in the range 0$\to$90$^{\circ}$. Your program should prompt for :
- Initial velocity v$_0$.
- Normalised drag coefficient $\beta$.
- Timestep $\Delta$t.
and plot out, using Matplotlib, a graph of K$_f$/K$_i$ against $\theta$.
### Checkpoint
- To help you with this checkpoint we've already provided you with a program skeletong structure below.
- Test your program with the following conditions:
- v$_0$=10 ms$^{−1}$, $\beta$=0.02, $\Delta$t=0.001
- v$_0$=10 ms$^{−1}$, $\beta$=0.3, $\Delta$t=0.001
- Call a demonstrator, show them your code, and run your programs for the two different conditions.
### Assessment
To pass this checkpoint:
- At least one of the progammes works fully with the supplied test examples.
- A good attempt has been made for the second program.
- The output graphs have the correct title and labels and the code is well structured with sensible use of functions called from the main() function.
- There are sensible header comments and sufficient in-line comments to explain the operation.
### Next steps
Once you have completed this checkpoint, please continue with the [Week5](../WeeklyTasks/Week5.ipynb) tasklist.
%% Cell type:code id:2539ac4d tags:
``` python
### Task 1 - plotting the trajectory
### Define here your imports
g = 9.81 #m/s**2
def set_initial(v_initial, theta):
""" Set inital conditions
"""
def acceleration(vx, vy, beta):
""" Calculate the acceleration
"""
def step_forward(x, y, vx, vy, beta, delta_t):
""" Do a forward step
"""
def Checkpoint5_1(v_initial, theta, beta, delta_t):
""" Main program to read in from terminal, do iteration, and plot out graph
"""
print("Trajectory for (v0, theta, beta, delta_t) = ({0}, {1}, {2}, {3})".format(v_initial, theta, beta, delta_t))
...
t = 0.0
x, y, vx, vy = set_initial(v_initial, theta)
...
x, y, vx, vy = step_forward(x, y, vx, vy, beta, delta_t)
while y>0.0:
#Propagate motion
plt.plot(x_values,y_values,"r") # Plot in red
plt.title("Trajectory of particle") # Add title/lables
plt.xlabel("Distance (m)")
plt.ylabel("Height")
plt.grid()
plt.show()
#Run interactively
v_initial = float(input("v0 (initial velocity) : "))
theta = float(input("theta (angle from the horizontal) : "))
beta = float(input("beta (drag coefficient) : "))
delta_t = float(input("timestep (delta_t) : "))
Checkpoint5_1(v_initial, theta, beta, delta_t)
```
%% Cell type:code id:a60f660b tags:
``` python
# Task 2 - ratio of kinetic energies
def kinetic_energy(v):
### Calculate kinetic energy
def Checkpoint5_2(v_initial, beta, delta_t):
print("K_f/K_i ratios for (v0, beta, delta_t) = ({0}, {1}, {2})".format(v_initial, beta, delta_t))
...
#Control resolution for angle input (total # calculations = 90 * angle_res)
angle_res = 10
#Scan all launch angles theta from 0 to 90
for theta in range(90*angle_res+1):
...
x, y, vx, vy = set_initial(v_initial, theta/angle_res)
x, y, vx, vy = step_forward(x, y, vx, vy, beta, delta_t)
while y>0.0:
#Propagate motion
#Plot graph
plt.plot(thetas, ke_ratios,"r") # Plot in red
plt.title("Ratio of Final-to-initial Kinetic Energies") # Add title/lables
plt.xlabel(r'Launch angles $\theta$ (degrees)')
plt.ylabel(r'K$_f$/K$_i$')
plt.xlim(0,90)
plt.grid()
plt.show()
#Run interactively
v_initial = float(input("v0 (initial velocity) : "))
beta = float(input("beta (drag coefficient) : "))
delta_t = float(input("timestep (delta_t) : "))
Checkpoint5_2(v_initial, beta, delta_t)
```
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment