\documentclass[aps,12pt]{revtex4}
\usepackage{color,graphics}
\usepackage{verbatim}
\usepackage{amsmath}

\newcommand{\dt}{\Delta t}
\newcommand{\tm}{\Delta t/2}
\newcommand{\vt}{v_t}

\begin{document}

\begin{center}
The Effect of Air Resistance\\
Harvey Gould\\
Physics 127\\
2/06/07
\end{center}

\section{Introduction}
I simulated the motion of
a falling body near the earth's surface in the presence of gravity
and air resistance. I used the Euler, Euler-Cromer, and
Euler-Richardson algorithms and considered both a linear and
quadratic dependence of the drag force on the velocity. Of
particular interest is the rise and fall of a pebble and the
relative times of ascent and descent.

\section{Method}
As shown in the text, Newton's laws of motion can be
written as two coupled first-order differential equations. The
Euler algorithm for the numerical solution of the latter equations
can be written in the form
\begin{subequations}
\label{eq:euler}
\begin{align}
v_{n+1} & = v_n + a_n \dt 	\\
y_{n+1} & = y_n + v_n \dt .
\end{align}
\end{subequations}
The Euler-Cromer algorithm modifies Eq.~(\ref{eq:euler}) slightly
and can be expressed as
\begin{subequations}
\begin{align}
v_{n+1} & = v_n + a_n \dt 	\\
y_{n+1}& = y_n + v_{n+1} \dt .
\end{align}
\end{subequations}
A more sophisticated version of the Euler algorithm evaluates the
acceleration and the velocity at the time $\tm = t + \dt/2$. The
Euler-Richardson algorithm can be written as follows:
\begin{subequations}
\begin{align}
a_m & = F(y_m,v_m,\tm)/m \\
v_m & = v_n + a_n \dt/2 	\\
y_m & = y_n + v_n \dt/2,
\end{align}
\end{subequations}
and
\begin{subequations}
\begin{align}
v_{n+1} & = v_n + a_m \, \dt	\\
y_{n+1} & = y_n + v_m \, \dt .
\end{align}
\end{subequations}

The total force on a falling body in the presence of air
resistance can be written in the convenient form
\begin{subequations}
\label{eq:drag}
\begin{align}
F & = -mg(1 - v^2/\vt^2)	\qquad \mbox{(quadratic
drag force)} \\
\noalign{\noindent or}
F & = -mg(1 - v/\vt)	\quad \qquad \mbox{(linear drag force)}
\end{align}
\end{subequations}
The signs in Eq.~(\ref{eq:drag}) are for a falling body with
the vertical coordinate increasing upward.

\section{Verification}
I made various modifications of so that I could
determine whether the Euler, Euler-Cromer, or Euler-Richardson
gave the most accurate results for the position of a freely
falling body. 
My numerical results using the three algorithms
are summarized in Table~\ref{tab:compare}. Note that the three 
algorithms give the exact result for the velocity because the
acceleration is constant for free fall. We see that the
Euler-Richardson algorithm gives the most accurate result for a
given value of
$\dt$. For this reason I will adopt the Euler-Richardson
algorithm with $\dt = 0.001$ for the rest of this report. I also
note that the Euler and Euler-Cromer algorithms overestimate and
underestimate the exact answer, respectively, and that their
average value gives the exact answer for $y$.

\begin{table}[h]
\begin{center}
\begin{tabular}{|l|l|l|l|}
\hline
$\dt$	& $y$ (Euler)	& $y$ (Euler-Cromer)	& $y$ (Euler-Richardson)
\\
\hline
0.1	  & 5.5900	& 4.6100	& 4.8550 \\
0.01	 & 5.1490	& 5.0510	& 5.0755 \\
0.001	& 5.1049	& 5.0951	& 5.0976 \\
\hline
\end{tabular}
\caption{\label{tab:compare}Comparison of the results of the Euler,
Euler-Cromer, and Euler-Richardson algorithms with the analytical
result, $y(t) = y_0 - v_0 t - gt^2/2 = 5.1$ at $t = 1$ with $y_0 =
10$ and $v_0 = 0$. Note that the Euler-Richardson algorithm gives
the most accurate result for a given value of $\dt$.}
\end{center}
\end{table}

\section{Data}
My first goal is to try to determine the form of the velocity dependence of the drag resistance on a 
styrofoam ball. To do so I will first use the data in
Table~\ref{tab:values} (the same as Table~3.1 in the text) to try 
to find the velocity dependence of the acceleration directly. The
velocity and acceleration were found from the approximations
\begin{align}
v(t) & = [y(t+\dt) - y(t)]/\dt \label{eq:vt} \\
a(t) & = [y(t+\dt) - 2y(t) + y(t-\dt)]/(\dt)^2 . \label{eq:at}
\end{align}
From Table~\ref{tab:values} we see that it is impossible to
determine the velocity-dependence of the acceleration. The reason
is that an error is introduced in trying to use the approximate
formulas (\ref{eq:vt}) and (\ref{eq:at}) to estimate $v$ and $a$.
For this reason I now assume a functional form of the velocity
dependence and try to find the best fit to the measured values of
the position given in Table~\ref{tab:values}.

\begin{table}[h]
\begin{center}
\begin{tabular}{|l|c|c|c|}
\hline
time	& position	& velocity	& acceleration \\
\hline
0	  & 0.075	& --	   & -- \\
0.1	  & 0.26	 & 2.65	& 8 \\
0.2	  & 0.525	& 3.45	& 8 \\
0.3	  & 0.870	& 4.00	& 5.5 \\
0.4	  & 1.27	 & 4.60	& 6 \\
0.5	  & 1.73	 & 5.00	& 4 \\
0.6	  & 2.23	 & 5.40	& 4 \\
0.7	  & 2.77	 & 5.80	& 4 \\
0.8	  & 3.35	 & --	& -- \\
\hline
\end{tabular}
\caption{\label{tab:values}The derived values of the velocity and
the acceleration from the empirical values of the position $y$.
The values of $y$ are taken from Table 3.1 of the text. Note that
the acceleration does not have an obvious dependence on the
velocity.}
\end{center}
\end{table}

I now use drag to determine the best fit to the empirical data. My results are 
summarized in Fig.~\ref{fig:fall}.

\begin{figure}[h]
\begin{center}
\includegraphics{fall}
\caption{\label{fig:fall}Comparison of the measured distance (open
squares) fallen by a styrofoam ball with the calculated distance
assuming a quadratic drag (open circles) and a linear drag (open
triangles).  The terminal velocity was assumed to be 6.5\,m/s in
both cases. Note that the quadratic form fits much better.}
\end{center}
\end{figure}

I wrote the class {\tt Pebble} so that I could simulate the time for
a pebble to ascend and descend. The main programming change I made
was to use the {\tt Math.abs} method
to account for the fact that the drag force changes sign.

I assumed a quadratic drag force and that the
ball is thrown vertically upward with an initial velocity $v0 =
50$. The drag force can be written in terms of the terminal
velocity even if the ball is moving upwards. I assumed that $\vt =
30$, a value much bigger than the value for the styrofoam ball. I
found that the time up was approximately 3.15 and the time down was
3.93. I also found that the speed of the ball when it returned to
the earth was 25.73, a value much less than the initial speed. All
units are SI.

\section{Analysis}
I did not try to optimize the value of the terminal velocity and used the Euler-Richardson algorithm 
with $\dt = 0.001$ for all runs. Because it was clear from my
results (see Fig.~1) that the quadratic form of the drag
resistance fits the empirical data better than the linear form, no
more analysis was really necessary.

\section{Interpretation}
I initially thought that the time up would be shorter than the time
down because the pebble starts out with a big velocity. I was
somewhat surprised to find that my intuition was correct. The ball 
reached a height of 61.03\,m, much less than the height of
$v_0^2/2g = 127.6$\,m it would have reached in the absence of drag
resistance.

\section{Critique}
It was fun to model the behavior of a system that I almost
understand. I was surprised that I understood that the time up
would be less than the time down. I would have liked to simulate
two- dimensional trajectories so I could have used my newfound
knowledge of Java's graphics capabilities, but I didn't have
time.

\end{document}