# Free fall

This is just elementary calculus, but it is harder than I thought at first glance. A friend told me that there is no formula expressing the altitude of a falling body if you is the inverse square law for the force, rather than making the approximation of a constant gravitational field. As far as I am able to tell, he is right. One can however compute some things.

We study the motion of a mass $m$, moving along an axis, whose position is denoted by $x(t)$. At the time $t=0$, the mass is at rest in the position $x(0)=x_0$. There is a mass $M$, at the position $x=0$, that does not move and attracts $m$ by a force $f$ given by the inverse square law:
$\displaystyle f=-\frac{GMm}{x^2},$
$G$ being the gravitational constant.

Let us write the conservation of the mechanical energy for the mass $m$: at any time $t\ge 0$,
$\displaystyle \frac{1}{2}mx'(t)^2-\frac{GMm}{x(t)}=-\frac{GMm}{x_0}.$
We consider a falling body, therefore $x'(t)\le 0$, and thus
$x'(t)=-\sqrt{2GM}\sqrt{\frac{1}{x(t)}-\frac{1}{x_0}},$
This yields
$-\frac{x'(t)}{\sqrt{\frac{1}{x(t)}-\frac{1}{x_0}}}=\sqrt{2GM},$
and therefore
$\displaystyle -\frac{x'(t)x(t)}{\sqrt{x(t)(x_0-x(t))}}=k, (1)$
with
$\displaystyle k=\sqrt{\frac{2GM}{x_0}}$
To find $t$ as a function of $x(t)$, we compute
$\displaystyle \int\frac{x\,dx}{\sqrt{x(x_0-x)}}.$
We have
$\displaystyle \frac{x\,dx}{\sqrt{x(x_0-x)}}=\frac{-(-2x+x_0)}{2\sqrt{x(x_0-x)}}+\frac{x_0}{2\sqrt{x(x_0-x)}},$
and therefore
$\displaystyle \int\frac{x\,dx}{\sqrt{x(x_0-x)}}=-\sqrt{x(x_0-x)}+\frac{x_0}{2}\int \frac{dx}{\sqrt{x(x_0-x)}}.$
We make the change of variable $x=\frac{x_0}{1+u^2}$, which gives
$\displaystyle dx=\frac{2x_0u\,du}{(1+u^2)^2}$
and
$\displaystyle \sqrt{x(x_0-x)}-x_0\arctan\left(\sqrt{\frac{x_0}{x}-1}\right).$
We obtain
$\displaystyle \int\frac{x\,dx}{\sqrt{x(x_0-x)}}=-\sqrt{x(x_0-x)}-\frac{x_0}{2}\int_{u=\sqrt{\frac{x_{0}}{x}-1}} \frac{2\,du}{1+u^2}$
and finally
$\displaystyle \int\frac{x\,dx}{\sqrt{x(x_0-x)}}=-\sqrt{x(x_0-x)}-x_0\arctan\left(\sqrt{\frac{x_0}{x}-1}\right)+C.$
We now integrate  $(1)$ between $0$ and $t$:
$\displaystyle k\,t=\sqrt{x(x_0-x)}+x_0\arctan\left(\sqrt{\frac{x_0}{x}-1}\right). (2)$
This gives us the duration of the fall as a function of the position of the falling mass $m$.

We can find another formula for the antiderivative. Let us make the change of variable $x=x_0\frac{1+\cos(u)}{2}=x_0\cos^2\left(\frac{u}{2}\right)$. Then
$dx=-x_0\sin\left(\frac{u}{2}\right)\cos\left(\frac{u}{2}\right),$ and
$\displaystyle \sqrt{x(x_0-x)}=x_0\sqrt{\cos^2\left(\frac{u}{2}\right)\left(1-\cos^2\left(\frac{u}{2}\right)\right)}=x_0\sin\left(\frac{u}{2}\right)\cos\left(\frac{u}{2}\right),$
$u$ being taken between $0$ and $\pi$. We obtain
$\displaystyle\int\frac{x\,dx}{\sqrt{x(x_0-x)}}=-\sqrt{x(x_0-x)}-\frac{x_0}{2}\int_{u=2\arccos\left(\sqrt{\frac{x}{x_0}}\right)} du,$
and therefore
$\displaystyle k\,t=\sqrt{x(x_0-x)}+x_0\arccos\left(\sqrt{\frac{x}{x_0}}\right), (3)$
which is the answer that can be found in Wikipedia. Taking $x=0$, we obtain the free-fall time $t_f$, the time it takes for the falling mass to reach the center of attraction:
$\displaystyle t_f=\frac{\pi}{2}\sqrt{\frac{x_0^{3}}{2GM}}.$

How does we find these changes of variable (this is the only tricky part in these calculations)? Well, if we set $y=\sqrt{x(x_0-x)}$, the points $(x,y)$, with $x$ in $[0,x_0]$, are on the circle of radius $x_0/2$ centered at $(x_0/2,0)$. Now if $(x(u),y(u))$ is a parametric representation of this circle, the change of variable $x=x(u)$ allows us to write
$\displaystyle \int \frac{dx}{\sqrt{x(x_0-x)}}=\int_{u=x^{-1}(x)}\frac{x'(u)\,du}{y(u)}.$
In the first change of variable, we used a parametrization of the circle by rational functions:
$\displaystyle\left\{\begin{array}{lcl} x(u)&=&\frac{x_0}{1+u^2};\\ y(u)&=&\frac{x_0\,u}{1+u^2}. \end{array} \right.$
In the second, we used a parametrization by trigonometric functions:
$\displaystyle\left\{\begin{array}{lcl} x(u)&=&x_0\frac{1+\cos(u)}{2};\\ y(u)&=&\frac{x_0}{2}\sin(u). \end{array} \right.$
Both are natural choices (the second maybe more so than the first), and both produce a simple result.

If we set
$\displaystyle F(\xi)=\xi\sqrt{1-\xi^2}+\arccos(\xi),$
we obtain
$\displaystyle k\,t=x_0\,F\left(\sqrt{\frac{x}{x_0}}\right)$
and thus
$\displaystyle x(t)=x_0\,\left(F^{-1}\right)^2\left(\frac{k\,t}{x_0}\right). (4)$
Indeed, $F$ is continuous on $[0,1]$, differentiable on $[0,1[$, and
$\displaystyle F'(\xi)=\sqrt{1-\xi^2}-\frac{\xi^2}{\sqrt{1-\xi^2}}-\frac{1}{\sqrt{1-\xi^2}}=\frac{-2\,\xi^2}{\sqrt{1-\xi^2}}.$
Therefore $F$ is decreasing, and thus a bijection from $[0,1]$ to $[0,\pi/2]$. Its inverse $F^{-1}$ is continuous on $[0,\pi/2]$ and differentiable on $]0,\pi/2[$.

On the other hand, if we set
$\displaystyle G(\xi)=\frac{\xi}{1+\xi^2}+\arctan(\xi),$
we obtain
$k\,t=x_0\,G\left(\sqrt{\frac{x_0}{x}-1}\right)$
and thus
$\displaystyle x(t)=\frac{x_0}{1+\left(G^{-1}\right)^{2}\left(\frac{k\,t}{x_0}\right)}.$
The function $G$ is continuous and differentiable on $\mathbb{R}$ and
$G'(\xi)=\frac{2}{(1+\xi^2)^2}.$
Its is therefore a diffeomorphism from $\mathbb{R}$ to $]-\pi/2,\pi/2[$. Let us denote by $g$ the inverse $G^{-1}$. It satisfies the differential equation
$\displaystyle g'(\zeta)=\frac{1}{2}+g(\zeta)^2+\frac{1}{2}g(\zeta)^4,$
with $g(0)=0$. In particular, $g'(0)=1/2$. This, taken with Equation $(4)$, implies that
$\displaystyle x(t)=x_0\left(1-\frac{k^2t^2}{4x_0^2}+O\left(\frac{k^3t^3}{x_0^3}\right)\right)$
for $t$ close to $0$. We have approximately
$\displaystyle x(t)\simeq x_0-\frac{1}{2}\gamma_0t^2,$
whith
$\displaystyle \gamma_0=\frac{GM}{x_0^2}.$
The right-hand side is the result that we would obtain for a constant gravitational field $\gamma_0$.

We have
$\displaystyle \lim_{\zeta\to \pi/2}g(\zeta)=+\infty.$
Since
$\displaystyle G(\xi)=\frac{\xi}{1+\xi^2}+\arctan(\xi)=\frac{\pi}{2}-\frac{1/\xi}{1+(1/\xi)^2}-\arctan(1/\xi),$
we have
$\displaystyle \zeta=\frac{\pi}{2}-\frac{1/g(\zeta)}{1+(1/g(\zeta))^2}-\arctan(1/g(\zeta))\\ = \frac{\pi}{2}-\frac{1}{g(\zeta)}-\frac{1}{g(\zeta)^3}-\frac{1}{g(\zeta)}+\frac{1}{3\,g(\zeta)^3}+o\left(\frac{1}{g(\zeta)^3}\right)\sim -\frac{2}{3\,g(\zeta)^3},$
and therefore
$\displaystyle g(\zeta)\sim \left(\frac{3}{2}\left(\frac{\pi}{2}-\zeta\right)\right)^{-1/3}.$
Therefore, when $t$ tends to $t_f$,
$\displaystyle x(t)\sim x_0\left(\frac{3}{2}\left(\frac{\pi}{2}-\frac{k\,t}{x_0}\right)\right)^{2/3}=\left(\frac{9GM}{2}\right)^{1/3}(t_f-t)^{2/3}.$
Using the differential equation, we also find
$\displaystyle x'(t)\sim-\sqrt{2GM}\left(\frac{9GM}{2}\right)^{-1/6}(t_f-t)^{-1/3}\\ =-\frac{2}{3}\left(\frac{9}{2}\right)^{1/3}(GM)^{1/3}(t_f-t)^{-1/3},$
which is what we would obtain by formally differentiating the equivalent.

The free-fall time can be rewritten
$\displaystyle t_f=\frac{\pi}{2}\sqrt{\frac{x^{2}}{2GM}}\,\frac{x_0}{x}\sqrt{x_0}=\frac{\pi x_0}{2x}\sqrt{\frac{x_0}{2\gamma}}$
with $\gamma$ the accelaration of gravity at the altitude $x$. How many time does it take to fall from the altitude of the Moon? We take $x$ to be  the radius of the Earth, and $x_0$ is the distance form the Earth to the Moon. Wikipedia gives
$x_0=384399$ kilometers, $x=6378$ kilometers and $\gamma=9.780327$ meters per second squared. We find roughly $116$ hours, which is to say four days and twenty hours.

The Wikipedia article I linked above gives two references for the free fall:

1. From Moon-fall to motions under inverse square laws, by S. K. Foong, European Journal of Physics, 2008 29: 987–1003;
2. Radial motion of Two mutually attracting particles, by Carl E. Mungan, The Physics Teacher, 2009, 47: 502-507.

I skimmed the first article (unfortunately it is behind a paywall, but I could get access with my university subscription). The formula for $x \ll x_0$ and the numerical value for the free fall time from the orbit of the moon agree with what I found. There is also a series expansion, and probably other stuff. I will read in detail later.