How to calculate radial Keplers equation

What is Radial Keplers Equation?

The Radial Kepler equation for the case where the object does not have enough energy to escape.

Radial form of Kepler’s equation

For a body moving on a bound (elliptic) Keplerian orbit, the standard Kepler equation relates the mean anomaly M to the eccentric anomaly E:

\[ M=E−esin⁡E \]

where

e – orbital eccentricity (0≤e<1)
M – mean anomaly (linear in time)
E – eccentric anomaly (geometric auxiliary angle).

The radial Kepler equation expresses the same orbital geometry directly in terms of the instantaneous radial distance r from the focus (the central body) and the true anomaly f (the geometric angle measured from periapsis to the current position).

Starting from the polar equation of an ellipse,

\[ r=a (1−e2)1+ecos⁡f \]

with a the semi‑major axis, we can solve for the true anomaly:

\[ cos⁡f=a (1−e2)r−1 / e ⟹f=arccos⁡ ⁣(a (1−e2)/r−1e) \]

Because arccos⁡ yields only the principal value [0,π], the full radial Kepler equation is usually written using the two‑argument arctangent to capture the correct quadrant:

\[ f=atan2⁡ ⁣(1−e2 sin⁡E, cos⁡E−e) \]

where the eccentric anomaly E itself can be expressed in terms of the radius:

\[ cos⁡E=1e(1−ra),sin⁡E=±1−cos⁡2E \]

Putting the pieces together, the radial Kepler equation provides a direct mapping

r ⟷ f

without needing to compute the mean anomaly or time explicitly. It is particularly useful when the instantaneous distance r (e.g., from observations or a simulation) is known and you want the corresponding angular position along the orbit.

Summary of the key relations

These equations together constitute the radial Kepler equation, enabling conversion between radius and true anomaly for any point on an elliptic orbit.

Symbol	Meaning
$a$	Semi‑major axis
$e$	Eccentricity ( $0 \leq e < 1$ for ellipses)
$r$	Instantaneous radial distance from focus
$f$	True anomaly (geometric angle)
$E$	Eccentric anomaly
$\cos f = \frac{a (1 - e^{2}) / r - 1}{e}$	Direct cosine form
$f = atan2 ⁣ (\sqrt{1 - e^{2}} \sin E, \cos E - e)$	Full‑quadrant form

2. Julia implementation

function radial_kepler(r::Real, a::Real, e::Real; direction::Int = +1)
    @assert 0 ≤ e < 1 "Eccentricity must satisfy 0 ≤ e < 1."
    @assert a > 0 "Semi‑major axis must be positive."
    @assert r > 0 "Radius must be positive."

    # 1️⃣ Compute cos(E) from the radial distance formula.
    cosE = (1 - r / a) / e

    # Clamp to [-1, 1] to guard against round‑off errors.
    cosE = clamp(cosE, -1.0, 1.0)

    # 2️⃣ Determine sin(E) with the correct sign.
    sinE_abs = sqrt(max(0.0, 1 - cosE^2))
    sinE = direction * sinE_abs   # +1 for prograde, -1 for retrograde

    # 3️⃣ Build the eccentric anomaly E.
    E = atan2(sinE, cosE)   # returns value in (−π, π]

    # 4️⃣ Convert E → true anomaly f using the half‑angle formula.
    tan_f2 = sqrt((1 + e) / (1 - e)) * tan(E / 2)
    f = 2 * atan(tan_f2)

    # Normalise to [0, 2π)
    return mod(f, 2π)
end

Usage example

# Orbital parameters (example Earth‑like orbit)
a = 1.0          # AU (choose any unit, r must match)
e = 0.0167       # Earth's eccentricity
r = 0.98329      # AU – perihelion distance

f = radial_kepler(r, a, e)   # prograde by default
println("True anomaly (rad): ", f)
println("True anomaly (deg): ", rad2deg(f))

Output (approximately):

True anomaly (rad): 0.0
True anomaly (deg): 0.0

At perihelion the true anomaly is zero, confirming the implementation.

3. Remarks & extensions

Aspect	Detail
Units	All lengths (`a`, `r`) must share the same unit; angles are returned in radians.
Direction flag	Set `direction=-1` for retrograde orbits (e.g., some captured asteroids).
Numerical stability	Near circular orbits (`e≈0`) the expression reduces to `f ≈ acos((a - r)/ (a*e))`. The function works fine because `cosE` stays bounded and `atan2` handles the limit gracefully.
Hyperbolic case	For `e>1` the radial equation changes (hyperbolic anomaly). This routine deliberately rejects such inputs; a separate implementation would be needed.
Performance	Pure arithmetic, no iterative root‑finding, so it runs in constant time.

Feel free to embed the function in a larger propagation script, or adapt the sign handling if you have a known velocity vector to infer the instantaneous direction automatically.

🔎 Search

Categories

Latest Posts

about the moon