Solution (source code)

= Solution

Let $f(r)=1-r_s/r$ and parametrize the <timelike geodesic> by <proper time> $\tau$. Time-translation and rotational symmetry conserve $E=f\dot t$ and $h=r^2\dot\phi$ in the equatorial plane. The normalized <four-velocity> obeys
$$
-1=-f\dot t^2+f^{-1}\dot r^2+r^2\dot\phi^2.
$$
Substitution and multiplication by $f$ give $\dot r^2+f(1+h^2/r^2)=E^2$, or
$$
\boxed{\tfrac12\dot r^2+V(r)=\tfrac12E^2,\qquad2V=(1-r_s/r)(1+h^2/r^2).}
$$
The dot denotes a proper-time derivative, and $E,h$ are <energy> and <angular momentum> per unit rest mass.