Equatorial null orbit in Schwarzschild spacetime (source code)

= Equatorial null orbit in Schwarzschild spacetime

For inverse radius $y=1/r$ and azimuth $\phi$, a Schwarzschild null geodesic obeys
$$
(y')^2+y^2-2My^3=\frac1{b^2},
\qquad
y''+y=3My^2,
$$
where $b=L/E$ is the <impact parameter>.