Tidy N-body orbital mechanics for R.
Version 1.0 — the function names, arguments, and defaults are now stable; anything that changes in a future version will be deprecated first. Feedback, bug reports, and contributions are welcome on GitHub.
Full documentation, examples, and guides at orbit-r.com
# Install from CRAN:
install.packages("orbitr")
# Or install the development version from GitHub:
# install.packages("devtools")
devtools::install_github("DRosenman/orbitr")
library(orbitr)
load_solar_system() |>
simulate_system(time_step = seconds_per_day, duration = seconds_per_year) |>
plot_orbits(three_d = FALSE)

load_solar_system() builds the Sun, all eight planets, the Moon, and
Pluto with real orbital data from JPL — eccentricities, inclinations,
and all. One line to build it, one line to simulate and plot.
Don’t need every body? Drop them with remove_body():
load_solar_system() |>
remove_body(c("Pluto", "Moon")) |>
simulate_system(time_step = seconds_per_day, duration = seconds_per_year) |>
plot_orbits(three_d = FALSE)
library(orbitr)
sim <- create_system() |>
add_sun() |>
add_body("Earth", mass = mass_earth, x = distance_earth_sun, vy = speed_earth) |>
simulate_system(time_step = seconds_per_day, duration = seconds_per_year)
sim |> plot_orbits()

Built-in constants like mass_sun, distance_earth_sun, and
speed_earth are real-world values in SI units — no Googling needed.
simulate_system() returns a tidy tibble, ready for dplyr, ggplot2,
plotly, or anything else.
animate_system(sim, fps = 15, duration = 5)

add_planet()Pick and choose real solar system bodies without looking up any numbers:
create_system() |>
add_sun() |>
add_planet("Earth", parent = "Sun") |>
add_planet("Mars", parent = "Sun") |>
simulate_system(time_step = seconds_per_day, duration = seconds_per_year * 2) |>
plot_orbits(three_d = FALSE)

Every planet’s mass, eccentricity, inclination, and orbital orientation are filled in from JPL data. Override any element to explore “what if” scenarios:
# What if Mars had a perfectly circular orbit?
create_system() |>
add_sun() |>
add_planet("Mars", parent = "Sun", e = 0) |>
simulate_system(time_step = seconds_per_day, duration = seconds_per_year * 2) |>
plot_orbits(three_d = FALSE)

create_system() |>
add_body("Earth", mass = mass_earth) |>
add_planet("Moon", parent = "Earth") |>
simulate_system(time_step = seconds_per_hour, duration = seconds_per_day * 28) |>
plot_orbits()

For full control, add_body_keplerian() lets you specify orbits using
classical elements — semi-major axis, eccentricity, inclination, and
orientation angles — instead of raw positions and velocities:
# A comet on a highly eccentric, tilted orbit
create_system() |>
add_sun() |>
add_planet("Earth", parent = "Sun") |>
add_body_keplerian(
"Comet", mass = 1e13, parent = "Sun",
a = 3 * distance_earth_sun, e = 0.85,
i = 60, lan = 45, arg_pe = 90, nu = 0
) |>
simulate_system(time_step = seconds_per_hour, duration = seconds_per_year * 5) |>
plot_orbits()

Kepler-16b orbits two stars — a real-life Tatooine.
G <- gravitational_constant
AU <- distance_earth_sun
m_A <- 0.68 * mass_sun
m_B <- 0.20 * mass_sun
a_bin <- 0.22 * AU
r_A <- a_bin * m_B / (m_A + m_B)
r_B <- a_bin * m_A / (m_A + m_B)
v_A <- sqrt(G * m_B^2 / ((m_A + m_B) * a_bin))
v_B <- sqrt(G * m_A^2 / ((m_A + m_B) * a_bin))
r_planet <- 0.7048 * AU
v_planet <- sqrt(G * (m_A + m_B) / r_planet)
create_system() |>
add_body("Star A", mass = m_A, x = r_A, vy = v_A) |>
add_body("Star B", mass = m_B, x = -r_B, vy = -v_B) |>
add_body("Kepler-16b", mass = 0.333 * mass_jupiter, x = r_planet, vy = v_planet) |>
simulate_system(time_step = seconds_per_hour, duration = seconds_per_day * 228.8 * 3) |>
plot_orbits()

load_solar_system() and add_planet()
use JPL orbital elements, no manual setup neededadd_body_keplerian() lets you define orbits
with semi-major axis, eccentricity, inclination, and orientation
anglesplot_orbits() auto-dispatches to interactive plotly
when any body has Z-axis motion; force either with
three_d = TRUE/FALSEanimate_system() renders GIFs with fading trails
via gganimateshift_reference_frame("Earth") re-centers
everything on any body, or on the system's center of mass with
"barycenter"add_body_keplerian() builds orbits
from elements (including hyperbolic flybys with e > 1);
get_orbital_elements() recovers the osculating elements of any body
at every time stepget_energy(), get_momentum(), and
get_angular_momentum() give the system totals at every time step;
conserved_quantities() joins them with their errors so you can see
whether a run is trustworthycontinue_simulation() picks up where a run left
off with a different time step, for comets and close encounterssave_system() / load_system() round-trip a
system to .rds; export_bodies() writes the body table to CSVadd_body_keplerian() and add_planet()The physics and numerics in orbitr were informed by Classical
Dynamics of Particles and Systems by Thornton & Marion, Computational
Physics by Mark Newman, and Data Structures & Algorithms in Python by
Canning, Broder & Lafore. The C++ engine was built with help from the
official Rcpp documentation on
CRAN. See the full Further
Reading article for
details on how each resource was used.
Shoutout to The College of New Jersey, where I first got hooked on this stuff. Go Lions.
MIT
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.