Analysis Propagation

Optimization missions are built with Mission and Phase. For propagation and initial-condition studies that do not need a solve, start with the octavian.propagate namespace:

from octavian import propagate

It collects the package's analysis propagators without hiding which physical model is selected.

Inertial Propagation With Optional Perturbations

propagate.inertial is the Earth-centered counterpart to propagate.cr3bp: it takes a Cartesian initial state and output times and returns [r, v, t] rows in SI units. It integrates Cartesian equations directly, including for circular and equatorial states. No solver or mission is needed.

import numpy as np
from octavian import EARTH, Perturbations, propagate, state

initial = state([7_000_000.0, 0.0, 0.0], [0.0, 7_500.0, 200.0])
times_s = np.linspace(0.0, 6_000.0, 361)
reference = propagate.inertial(initial, times_s)
j2_reference = propagate.inertial(
    initial, times_s, central_body=EARTH,
    perturbations=Perturbations(j2=True), max_step_s=10.0,
)

The default force is point-mass gravity. Optional forces share the existing relative propagator's implementation: spherical harmonics or J2, Moon/Sun gravity, cannonball drag, and cannonball SRP. Moon/Sun gravity or SRP requires initial_epoch and uses Earth-centered TOD ephemerides; supply the initial state in that frame. Drag/SRP also requires spacecraft with mass and positive corresponding cannonball areas. Spacecraft mass remains constant; this is unpowered propagation.

Times must be strictly monotonic with zero at the first or last sample. Backward histories are supported. Fixed-step RK4 shortens its internal steps to hit the requested times; max_step_s (default 10 seconds) controls accuracy independently of output spacing. Reduce the step to check convergence for long or demanding trajectories. ephemeris_step_s controls Sun/Moon interpolation.

Pass the history directly as reference_trajectories=[{"name": "Reference", "traj": reference}] to a trajectory visualizer. For an elliptic orbit, its nominal Keplerian period is a useful display span; perturbed orbits need not close after that time. See the reference plotting guide.

Available Models

Call	Model	Return
`propagate.inertial(...)`	Cartesian RK4 with optional harmonics or J2, third-body gravity, drag, and SRP	`[r, v, time]` array
`propagate.two_body(...)`	Elliptic point-mass inertial motion	`[r, v, time]` array
`propagate.cwh(...)`	Linear circular-chief CWH	`[RIC state, time]` array
`propagate.nonlinear_ric(...)`	Exact circular-chief RIC before linearization	`[RIC state, time]` array
`propagate.relative(...)`	Coupled chief/deputy absolute states	`RelativePropagationResult`
`propagate.relative_elements(...)`	D'Amico or classical relative elements	`RelativeElementPropagationResult`
`propagate.cr3bp(...)`	Dimensional or canonical synodic CR3BP	`[state, time]` array

Every ordinary state-history array has seven columns, with time last. The relative result objects keep multiple useful views together instead of discarding information.

Inertial And Reduced Relative Histories

import numpy as np
from octavian import CR3BPSystem, EARTH, Perturbations, propagate, state

radius_m = EARTH.mean_radius_m + 500_000.0
chief = state(
    [radius_m, 0.0, 0.0],
    [0.0, np.sqrt(EARTH.mu_m3ps2 / radius_m), 0.0],
)
times_s = np.linspace(0.0, 600.0, 13)

inertial = propagate.inertial(chief, times_s)

mean_motion = np.sqrt(EARTH.mu_m3ps2 / radius_m**3)
relative = propagate.cwh(
    [100.0, -500.0, 0.0, 0.0, 0.02, 0.0],
    times_s,
    mean_motion_radps=mean_motion,
)

CWH is deliberately named in the call. Use nonlinear_ric when the full circular-chief relative equations are needed, or relative to propagate an arbitrary chief and deputy as absolute Cartesian states.

Coupled Relative Propagation

relative_initial = state(
    [100.0, -500.0, 0.0],
    [0.0, 0.02, 0.0],
)
result = propagate.relative(
    chief,
    relative_initial,
    times_s,
    perturbations=Perturbations(j2=True),
)

ric = result.relative_trajectory_ric
chief_eci = result.chief_trajectory_eci
deputy_eci = result.deputy_trajectory_eci

The same call accepts Moon/Sun gravity, drag, and SRP. Supply chief_spacecraft and deputy_spacecraft for differential cannonball forces, and initial_epoch whenever an ephemeris is required.

Relative Elements Without Duplicate Propagation

relative_elements returns both representations from one coupled numerical history:

result = propagate.relative_elements(
    initial_roe,
    times_s,
    chief_initial_state_eci=chief,
    mu_m3ps2=EARTH.mu_m3ps2,
    perturbations=Perturbations(j2=True, drag=True),
    chief_spacecraft=chief_vehicle,
    deputy_spacecraft=deputy_vehicle,
)

osculating_elements = result.elements
ric_for_plotting = result.ric

result.representation records whether the native rows are D'Amico or classical relative elements, and result.times_s returns their shared time vector.

CR3BP

system = CR3BPSystem.earth_moon()
l4 = state(
    system.lagrange_points(dimensional=False)["L4"],
    [0.0, 0.0, 0.0],
)
history = propagate.cr3bp(
    l4,
    [0.0, 0.01],
    system=system,
    dimensional=False,
)

The dimensional flag applies to both state and time. See the Cislunar CR3BP tutorial for conversions, Jacobi diagnostics, and composable optimization.

Specialized Functions Remain Available

The umbrella namespace is additive. Existing imports such as octavian.relative.propagate_cwh, octavian.relative.propagate_relative_numerical, and octavian.cislunar.propagate_cr3bp remain supported for code that benefits from the domain-specific module layout.

The complete executable overview is examples/quick/08_propagation_models.py.

Native Two-Body Integration

For adaptive numerical point-mass propagation, TwoBodyECI defaults to ASSET's specialized C++ Kepler ODE and integrator. Its entire integration loop runs in C++, with no Python force callbacks. This uses the existing ASSET wheel and needs no additional extension build.

from octavian import EARTH
from octavian.dynamics import TwoBodyECI

ode = TwoBodyECI(mu_m3ps2=EARTH.mu_m3ps2)
integrator = ode.integrator(10.0)  # Initial step in seconds.
integrator.setAbsTol(1e-7)
history = integrator.integrate_dense(
    [7e6, 0.0, 0.0, 0.0, 7500.0, 0.0, 0.0], 6000.0, 401
)

Rows contain [r, v, time] in SI units. Pass backend="expression" to compare against the generic ASSET expression integrator. Both are compiled; the specialized backend also specializes the integration loop for Kepler dynamics. Existing phase creation and variable groups remain available. propagate.two_body already uses ASSET's analytical C++ Cartesian Kepler propagation when ASSET is installed; it is useful for direct evaluation at arbitrary times. The numerical integrator above is a separate option.

native/spherical_harmonics/runtime_comparison.py reports raw integration times and checks agreement between the two numerical backends.