Add air resistance to projectile motion

I am given an initial x and y position and initial velocity and I was asked to graph the trajectory in 1 second intervals.

This is what I have so far:

If $x_0 = 1, v_{0x} = 70, y_0 = 0, v_{0y} = 80, a_x = 0, a_y = -9.8$, and time will be 0,1,2,3... and so on.

Using these equations on every second you can find the plot will be a bell shaped with the highest point being ~ 325 m at about 600 seconds:

$$x = x0 + (v_{0x})t + 1/2((a_x)t^2)$$ $$y = y0 + (v_{0y})t + 1/2((a_y)t^2)$$

Usually in physics, we are taught in perfect condition with no air resistance. But what if there was air resistance? How would it affect this problem? How can i add it to my calculations and see what the difference is?

• With v^2 drag, there is no closed form solution. I posted (x,y) expressions for linear drag on some previous question. If you want a realistic model, there is not much alternative to numerical integration. Finding the top altitude will also be iterative. Sep 5 '14 at 21:31
• More on projectile motion with air drag. Sep 5 '14 at 22:04
• One should probably add, that if there is air resistance, then there is not only drag, but also lift, and sideways deflection if the projectile is rotating. Before you know it, you will be replicating a hundred years of work on realistic ballistic equations. Sep 6 '14 at 0:51
• Initial vertical velocity of $80 m/s$ and acceleration of gravity = $-9.8 m/s^2$ implies you will reach a maximum height after a little over 8 seconds. Not sure where your 600 seconds comes from? Fix that before worrying about drag... Sep 6 '14 at 2:36
• Air resistance is a function of cross section area and $\mathrm {velocity}^2$. So the you will either need to integrate the formula and recalculate it, or divide your time into smaller sections and calculate a new velocity for each point.
– LDC3
Sep 6 '14 at 3:34

With drag force $- \alpha \left|{\dot{\bf r}}\right| {\dot{\bf r}}$ and gravitational force $-mg {\bf {\hat{y}}}$, the equations of motion are (see my answer to this question) \begin{align} {\ddot{x}} &= - \beta {\dot{x}} \sqrt{{\dot x}^2+{\dot y}^2} \\ {\ddot{y}} &= - g - \beta {\dot{y}} \sqrt{{\dot x}^2+{\dot y}^2} \end{align} where $\beta = \alpha / m$, $\alpha = \rho C_d A / 2$, $\rho$ is the air density, $C_d$ is the drag coefficient, and $A$ is the cross-sectional area of the projectile.

I've never seen the above equations solved analytically. Here's some Python code that plots the solution with and without air resistance. To run, copy into a text file myFile.py and run

python myFile.py

from the command line.

from pylab import *
import math

# Physical constants
g = 9.8
m = 1.0
rho = 1.0
Cd = 1.0
A = math.pi * pow(0.01,2.0)
alpha = rho * Cd * A / 2.0
beta = alpha / m

# Initial conditions
X0 = 1.0
Y0 = 0.0
Vx0 = 70.0
Vy0 = 80.0

# Time steps
steps = 1000
t_HIT = 2.0*Vy0/g
dt = t_HIT / steps

# No drag
X_ND = list()
Y_ND = list()

for i in xrange(steps+1):
X_ND.append(X0 + Vx0 * dt * i)
Y_ND.append(Y0 + Vy0 * dt * i - 0.5 * g * pow(dt * i,2.0))

# With drag
X_WD = list()
Y_WD = list()
Vx_WD = list()
Vy_WD = list()

for i in xrange(steps+1):
X_ND.append(X0 + Vx0 * dt * i)
Y_ND.append(Y0 + Vy0 * dt * i - 0.5 * g * pow(dt * i,2.0))

# With drag
X_WD = list()
Y_WD = list()
Vx_WD = list()
Vy_WD = list()

X_WD.append(X0)
Y_WD.append(Y0)
Vx_WD.append(Vx0)
Vy_WD.append(Vy0)

stop = 0
for i in range(1,steps+1):
if stop != 1:
speed = pow(pow(Vx_WD[i-1],2.0)+pow(Vy_WD[i-1],2.0),0.5)

# First calculate velocity
Vx_WD.append(Vx_WD[i-1] * (1.0 - beta * speed * dt))
Vy_WD.append(Vy_WD[i-1] + ( - g - beta * Vy_WD[i-1] * speed) * dt)

# Now calculate position
X_WD.append(X_WD[i-1] + Vx_WD[i-1] * dt)
Y_WD.append(Y_WD[i-1] + Vy_WD[i-1] * dt)

# Stop if hits ground
if Y_WD[i] <= 0.0:
stop = 1

# Plot results
plot(X_ND, Y_ND)
plot(X_WD, Y_WD)
show()
• I think the OP does not know how to simulate the above equations (even with the most basic Euler method) and so maybe you can add a reference to this step in your answer to get him/her pointed in the right direction. Sep 10 '14 at 13:36
• I'll try to update my answer some time today with a python script to solve the equations of motion. Sep 10 '14 at 15:19