Elements of an Orbit

Kepler nailed it. Conventional names and ordering of the elements: a, e, i, Ω, ϖ, L.

The elements of the solar system, (.) may be used to compute Cartesian coordinates and velocities of the planets. (discrepancy, now ,) Here is how.

We use this formula for an ellipse:

r = A/(1 + e cos θ) In some orthogonal coordinate system in some plane that includes the attractor at the origin. (A = a(1−e²))

The origin is at a focus of this ellipse.
The above reference uses these terms:

a = semi-major axis (a>0): Half the long axis of the ellipse (a = A(1/(1+e) + 1/(1−e))/2 = A/(1−e²))
e = eccentricity (0≤e<1): Departure from circle; see e in equation.
i = inclination (0≤i<π/2): angle between nominal plane of planet’s orbit and plane of Earth’s orbit.
Ω = longitude of the ascending node (−π<Ω≤π): Where the planet passes to the North side of the Earth’s orbital plane. Measured as an angle from a somewhat arbitrary point in the Earth’s orbit.
ϖ = longitude of perihelion (−π<ϖ≤π): Where the planet is closest to the Sun. Measured from the same arbitrary point. This is the sum of two angles, the first in the Earth’s orbit, to the ascending node, and the second in the planet’s orbit. This value would be awkward for a retrograde orbit, but nonesuch are of interest.
L = mean longitude (0≤L<2π): The angle which the planet would be at now, if it had been moving at its constant average angular velocity since perihelion. (That’s weird!) Measured from the same arbitrary point. L changes by 2π each orbit and is linear in time.

See too JPL’s Glossary.

Limits on a and e mean the formulae fail. The limitation of i to π/2 is only because these conventions don’t work well for retrograde orbits. Other angular limits are merely conventional. Traditionally they are given in degrees. We will use radians.

We use the ‘astronomical unit’ (AU) as the length unit. For Earth, a=1. We take the Sun as a ‘unit attractor’ which means that when a=1 and e=0 the velocity = 1 and it is one year between t=0 and t=2π. This is confusing but I find the alternatives more confusing.

Here is the plan, and code, to compute x, y, z, and their time derivatives:

Use e: Consider a unit ellipse with a focus at the origin, eccentricity = e, and X-axis as major axis. I.e. convert r = A/(1 + e cos θ) to Cartesian coordinates.
Use L: Locate body (x, y) on orbit iterating Kepler’s equal area equation.
Use a: Transform coordinates by scale factor A/(1−e²).
Use ϖ−Ω: Rotate ellipse.
Use i: Rotate in 3D about line of nodes by inclination.
Use Ω: Rotate line of nodes by longitude of the ascending node.

The last four steps are matrix multiplies.

As e varies, the orbit r = 1/(1+e cos θ) has the same specific angular momentum. The period of the orbit is 2π(1−e²)^−3/2 which is the area of the ellipse which is just the entire area swept out in one period. N.B. When e > 0 → (1−e²)^−3/2 > 1.