Library

AstroForceModels

A comprehensive Julia package for modeling astrodynamics forces affecting satellite orbital motion.

This package provides a unified framework for computing accelerations from various perturbative forces including atmospheric drag, solar radiation pressure, third-body gravity, relativistic effects, and gravitational harmonics. All force models implement a common acceleration interface and can be efficiently combined using the CentralBodyDynamicsModel system.

Key Features

• Unified Interface: All force models implement acceleration(state, params, time, model)
• Efficient Combination: Use CentralBodyDynamicsModel and build_dynamics_model for optimal performance
• SatelliteToolbox Integration: Built on the SatelliteToolbox.jl ecosystem
• Automatic Differentiation: All models support ForwardDiff.jl for gradient-based optimization
• High Performance: Type-stable implementations with compile-time optimizations

Force Models Available

• Gravity Models: Keplerian and spherical harmonics (via SatelliteToolboxGravityModels)
• Atmospheric Drag: Multiple atmospheric models (JB2008, JR1971, NRLMSISE00, Harris-Priester, Modified Harris-Priester, Exponential)
• Solar Radiation Pressure: With shadow modeling (conical, cylindrical)
• Earth Albedo Radiation Pressure: Reflected and thermal Earth radiation
• Third-Body Gravity: Sun, Moon, and planetary perturbations
• Solid Body Tides: Geopotential perturbations from tidal deformation (IERS 2010)
• Plasma Drag: Ionospheric ion drag with Chapman layer density model
• Thermal Emission: Spacecraft thermal re-radiation from anisotropic surface emission
• Geomagnetic Lorentz Force: Lorentz force on charged spacecraft from Earth's magnetic field
• Relativistic Effects: Schwarzschild, Lense-Thirring, and de Sitter effects
• Low-Thrust Propulsion: Constant, tangential, and user-defined thrust profiles

Example Usage

using AstroForceModels

# Create individual force models
gravity = GravityHarmonicsAstroModel(...)
drag = DragAstroModel(...)
srp = SRPAstroModel(...)

# Combine efficiently
dynamics = CentralBodyDynamicsModel(gravity, (drag, srp))

# Use in ODE integration
function satellite_dynamics!(du, u, p, t)
    du[1:3] = u[4:6]  # velocity
    du[4:6] = build_dynamics_model(u, p, t, dynamics)
end

See Also

• Usage Guide: Comprehensive examples
• API Reference: Complete function documentation
• SatelliteToolbox.jl: Ecosystem foundation

Abstract type for albedo radiation models used in albedo force calculations.

AbstractAstroForceModel

Abstract base type for all astrodynamics force models.

All force models must implement the acceleration(state, params, time, model) interface to compute the 3-dimensional acceleration vector acting on a spacecraft.

Subtypes

• AbstractPotentialBasedForce: Forces derivable from a potential (gravity)
• AbstractNonPotentialBasedForce: Non-conservative forces (drag, SRP, etc.)

AbstractDynamicsModel

Abstract base type for dynamics models that combine multiple force models.

Dynamics models provide efficient ways to compute total acceleration from multiple force sources. The primary implementation is CentralBodyDynamicsModel.

AbstractIonosphereModel

Abstract type for ionospheric density models used to compute ion mass density for plasma drag calculations.

Implementations

• ChapmanIonosphere: Chapman layer model for the F2 region
• ConstantIonosphere: Fixed ion density for testing and parametric studies
• NoIonosphere: Zero ion density (disables plasma drag)

AbstractNonPotentialBasedForce <: AbstractAstroForceModel

Abstract type for non-conservative force models that cannot be derived from a potential.

Examples include atmospheric drag, solar radiation pressure, and relativistic effects. These forces typically depend on velocity and other non-conservative factors.

AbstractPotentialBasedForce <: AbstractAstroForceModel

Abstract type for conservative force models derivable from a gravitational potential.

Examples include point-mass gravity, gravitational harmonics, and third-body gravity. These forces depend only on position and can be derived from a potential function.

Abstract satellite drag model to help compute drag forces.

AbstractSatellitePlasmaDragModel

Abstract type for satellite models used to compute ion ballistic coefficients for plasma drag.

Ion drag coefficients typically range from 2.0 to 4.0 depending on surface accommodation and electrostatic sheath effects, compared to 2.0-2.5 for neutral atmospheric drag.

Abstract satellite SRP model to help compute solar radiation pressure forces.

Abstract satellite thermal emission model.

AbstractSpacecraftChargeModel

Abstract type for spacecraft charge-to-mass ratio models.

Subtypes must implement charge_mass_ratio(u, p, t, model) returning q/m in C/kg.

AbstractThrustFrame

Abstract base type for thrust reference frame specifications.

The frame determines how the 3-component acceleration vector produced by an AbstractThrustModel is interpreted before being transformed into the inertial (ECI) frame.

Subtypes

• InertialFrame: Components are inertial (ECI) — no rotation needed
• RTNFrame: Radial / Transverse / Normal (orbit-fixed)
• VNBFrame: Velocity / Normal / Binormal (velocity-fixed)

AbstractThrustModel

Abstract base type for all low-thrust propulsion models.

Concrete subtypes define how the thrust acceleration vector is computed from the spacecraft state, simulation parameters, and time. Each subtype must implement:

thrust_acceleration(u, p, t, model) → SVector{3}

returning a 3-component thrust acceleration [km/s²]. The components are interpreted in the reference frame specified by the parent LowThrustAstroModel's frame field (see AbstractThrustFrame).

Albedo Astro Model struct Contains information to compute the acceleration from Earth albedo radiation pressure.

Surface positions are pre-computed at construction time on a spherical Earth using Lebedev quadrature points, eliminating expensive geodetic conversions from the integration loop.

Fields

• satellite_shape_model::AbstractSatelliteSRPModel: The satellite shape model for computing the ballistic coefficient.
• sun_data::ThirdBodyModel: The data to compute the Sun's position.
• body_albedo_model::AbstractAlbedoModel{<:Number, <:Number}: The Earth albedo radiation model.
• eop_data::EopIau1980: Earth orientation parameters.
• solar_irradiance::Number: Solar irradiance at 1 AU [W/m²].
• speed_of_light::Number: Speed of light [km/s].
• surface_positions_ecef::Vector{SVector{3,Float64}}: Pre-computed surface element positions in ECEF [km].
• weights::Vector{Float64}: Scaled Lebedev quadrature weights.

AtmosphericModelType

Abstract base type for atmospheric density model selectors used by compute_density.

Concrete subtypes: JB2008, JR1971, MSIS2000, ExpAtmo, HarrisPriester, HarrisPriesterModified, NoAtmosphere.

Cannonball Fixed Drag Model struct Contains information to compute the ballistic coefficient of a cannonball drag model with a fixed ballistic coefficient.

Fields

• radius::Number: The radius of the spacecraft.
• mass::Number: The mass of the spacecraft.
• drag_coeff::Number: The drag coefficient of the spacecraft.
• ballistic_coeff::Number: The fixed ballistic coefficient to use.

CannonballFixedDrag(radius::Number, mass::Number, drag_coeff::Number)

Constructor for a fixed cannonball ballistic coefficient drag model.

The ballistic coefficient is computed with the following equation:

            BC = CD * area/mass

where area is the 2D projection of a sphere

            area = π * r^2

Arguments

• radius::Number: The radius of the spacecraft.
• mass::Number: The mass of the spacecraft.
• drag_coeff::Number: The drag coefficient of the spacecraft.

Returns

• drag_model::CannonballFixedDrag: A fixed ballistic coefficient drag model.

CannonballFixedDrag(ballistic_coeff::Number)

Constructor for a fixed ballistic coefficient drag model.

Arguments

• ballistic_coeff::Number: The fixed ballistic coefficient to use.

Returns

• drag_model::CannonballFixedDrag: A fixed ballistic coefficient drag model.

CannonballFixedPlasmaDrag{RT,MT,DT,BT} <: AbstractSatellitePlasmaDragModel

Cannonball plasma drag model with a fixed ion ballistic coefficient.

The ion ballistic coefficient is computed as BCi = CD,i × A / m, where CD,i is the ion drag coefficient. For a sphere in the free molecular flow regime with full surface accommodation, CD,i ≈ 4.0 (Chapra 1961, Lafleur 2023).

Fields

• radius::RT: Spacecraft radius [m].
• mass::MT: Spacecraft mass [kg].
• drag_coeff::DT: Ion drag coefficient (dimensionless). Typical range 2.0-4.0.
• ballistic_coeff::BT: Precomputed C_D,i × A / m [m²/kg].

Example

# From physical properties (sphere r=1m, m=500kg, Cd_i=4.0)
sat = CannonballFixedPlasmaDrag(1.0, 500.0, 4.0)

# From precomputed ballistic coefficient
sat = CannonballFixedPlasmaDrag(0.025)

CannonballFixedPlasmaDrag(radius::Number, mass::Number, drag_coeff::Number)

Construct a cannonball plasma drag model from physical properties.

The ion ballistic coefficient is BCi = CD,i × π r² / m.

Arguments

• radius::Number: Spacecraft radius [m].
• mass::Number: Spacecraft mass [kg].
• drag_coeff::Number: Ion drag coefficient (dimensionless).

CannonballFixedPlasmaDrag(ballistic_coeff::Number)

Construct a fixed ion ballistic coefficient plasma drag model.

Arguments

• ballistic_coeff::Number: Precomputed C_D,i × A / m [m²/kg].

Cannonball Fixed SRP Model struct Contains information to compute the reflectivity ballistic coefficient of a cannonball SRP model with a fixed reflectivity coefficient.

Fields

• radius::Number: The radius of the spacecraft.
• mass::Number: The mass of the spacecraft.
• reflectivity_coeff::Number: The reflectivity coefficient of the spacecraft.
• reflectivity_ballistic_coeff::Number: The fixed ballistic coefficient to use.

CannonballFixedSRP(radius::Number, mass::Number, reflectivity_coeff::Number)

Constructor for a fixed cannonball reflectivity ballistic coefficient SRP model.

The reflectivity ballistic coefficient is computed with the following equation:

            RC = Cr * area/mass

where area is the 2D projection of a sphere

            area = π * r^2

Arguments

• radius::Number: The radius of the spacecraft.
• mass::Number: The mass of the spacecraft.
• reflectivity_coeff::Number: The reflectivity coefficient of the spacecraft.

Returns

• srp_model::CannonballFixedSRP: A fixed reflectivity ballistic coefficient SRP model.

CannonballFixedSRP(reflectivity_ballistic_coeff::Number)

Constructor for a fixed reflectivity ballistic coefficient SRP model.

Arguments

• reflectivity_ballistic_coeff::Number: The fixed reflectivity ballistic coefficient to use.

Returns

• srp_model::CannonballFixedSRP: A fixed reflectivity ballistic coefficient SRP model.

CelestialBody{Name, T<:Number}

Celestial body representation with compile-time identity for type-stable dispatch.

Type Parameters

• Name::Symbol: Compile-time body identifier (e.g., :Sun, :Moon, :Earth)
• T<:Number: Numeric type for physical parameters

Fields

• central_body::Symbol: Name of the central body being orbited
• jpl_code::Int: NAIF ID Code
• μ::T: Gravitational parameter [km^3/s^2]
• Req::T: Equatorial radius [km]

CentralBodyDynamicsModel{N,GT,PT} <: AbstractDynamicsModel

A dynamics model structure that efficiently combines a central body gravity model with multiple perturbing force models to compute the total acceleration acting on a spacecraft.

The model treats gravity as the dominant central body force and adds perturbations from other effects such as atmospheric drag, solar radiation pressure, third-body gravity, and relativistic effects.

Type Parameters

• N::Int: Number of perturbing force models
• GT <: AbstractGravityAstroModel: Type of the gravity model
• PT <: Tuple: Type of the tuple containing perturbing models

Fields

• gravity_model::GT: The central body gravitational force model (e.g., KeplerianGravityAstroModel, GravityHarmonicsAstroModel)
• perturbing_models::PT: Tuple of perturbing force models (e.g., drag, SRP, third-body, relativistic effects)

Constructors

# With explicit gravity model and perturbations
CentralBodyDynamicsModel(gravity_model, (drag_model, srp_model, ...))

# With only perturbations (uses Keplerian gravity by default)
CentralBodyDynamicsModel((drag_model, srp_model, ...))

# With only gravity model (no perturbations)
CentralBodyDynamicsModel(gravity_model)

Example

using AstroForceModels

# Create force models
gravity = GravityHarmonicsAstroModel(...)
drag = DragAstroModel(...)
srp = SRPAstroModel(...)
third_body = ThirdBodyAstroModel(...)

# Combine into dynamics model
dynamics = CentralBodyDynamicsModel(gravity, (drag, srp, third_body))

# Use in ODE integration
function satellite_dynamics!(du, u, p, t)
    du[1:3] = u[4:6]  # velocity
    du[4:6] = build_dynamics_model(u, p, t, dynamics)
end

See Also

• build_dynamics_model: Function to compute total acceleration
• AbstractDynamicsModel: Parent abstract type
• acceleration: Individual force model acceleration interface

CentralBodyDynamicsModel(gravity_model)

Create a dynamics model with only a gravity model (no perturbations).

Arguments

• gravity_model::AbstractGravityAstroModel: Central body gravity model

CentralBodyDynamicsModel(perturbing_models)

Create a dynamics model with perturbing forces and default Keplerian gravity.

Arguments

• perturbing_models::NTuple{N,AbstractAstroForceModel}: Tuple of perturbing force models

CentralBodyDynamicsModel(gravity_model, perturbing_models)

Create a dynamics model with explicit gravity and perturbing force models.

Arguments

• gravity_model::AbstractGravityAstroModel: Central body gravity model
• perturbing_models::NTuple{N,AbstractAstroForceModel}: Tuple of perturbing force models

ChapmanIonosphere{NmT,HmT,HsT,MiT} <: AbstractIonosphereModel

Chapman layer model for the F2 region of the ionosphere.

The ion number density follows the Chapman production function profile:

nᵢ(h) = Nₘₐₓ ⋅ exp(½(1 - z - exp(-z)))

where z = (h - hₘₐₓ) / Hₛ, and the ion mass density is ρᵢ = nᵢ ⋅ mᵢ.

Default parameters represent typical mid-latitude, moderate solar activity conditions. The dominant ion in the F2 region (200-500 km) is O⁺.

Fields

• Nmax::NmT: Peak ion number density [m⁻³]. Typical range: 10¹⁰ (solar min) to 10¹² (solar max).
• hmax::HmT: Altitude of the F2 peak [km]. Typical range: 250-450 km.
• Hs::HsT: Scale height of the F2 layer [km]. Typical range: 50-80 km.
• mi::MiT: Mean ion mass [kg]. Default is O⁺ mass (2.6567628e-26 kg).

References

• Rishbeth & Garriott (1969), "Introduction to Ionospheric Physics"
• Lafleur (2023), "Charged aerodynamics", Acta Astronautica 212, 370-386

ConstantCartesianThrust{T1,T2,T3} <: AbstractThrustModel

Constant thrust acceleration expressed as a fixed 3-component vector.

The components are interpreted in the reference frame set on the parent LowThrustAstroModel:

• With InertialFrame: (ax, ay, az) are ECI components
• With RTNFrame: (ax, ay, az) are (a_R, a_T, a_N) components
• With VNBFrame: (ax, ay, az) are (a_V, a_N, a_B) components

Type Parameters

• T1 <: Number: Type of the first-component acceleration
• T2 <: Number: Type of the second-component acceleration
• T3 <: Number: Type of the third-component acceleration

Fields

• ax::T1: First-component thrust acceleration [km/s²]
• ay::T2: Second-component thrust acceleration [km/s²]
• az::T3: Third-component thrust acceleration [km/s²]

Constructors

ConstantCartesianThrust(ax, ay, az)

Create from acceleration components directly [km/s²].

ConstantCartesianThrust(direction, thrust, mass)

Create from a thrust direction vector (will be normalized), thrust magnitude [N], and spacecraft mass [kg]. The acceleration is computed as:

    a = (T / m) / 1000 * d̂   [km/s²]

ConstantCartesianThrust(direction::AbstractVector, thrust::Number, mass::Number)

Construct a constant Cartesian thrust model from a direction vector, thrust [N], and mass [kg].

ConstantIonosphere{DT} <: AbstractIonosphereModel

Ionosphere model with a fixed ion mass density, useful for testing and parametric studies.

Fields

• rho_i::DT: Fixed ion mass density [kg/m³].

Example

iono = ConstantIonosphere(rho_i=1e-17)

ConstantTangentialThrust{MT<:Number} <: AbstractThrustModel

Constant-magnitude thrust directed along (or opposite to) the velocity vector.

When the magnitude is positive, thrust is aligned with the velocity direction (orbit raising). When negative, thrust opposes the velocity direction (orbit lowering).

Type Parameters

• MT <: Number: Type of the acceleration magnitude

Fields

• magnitude::MT: Thrust acceleration magnitude [km/s²]. Positive = along velocity, negative = anti-velocity.

Constructors

ConstantTangentialThrust(magnitude)

Create from acceleration magnitude directly [km/s²].

ConstantTangentialThrust(thrust, mass)

Create from thrust [N] and spacecraft mass [kg]. The acceleration magnitude is computed as:

    |a| = T / (m × 1000)   [km/s²]

ConstantTangentialThrust(thrust::Number, mass::Number)

Construct a constant tangential thrust model from thrust [N] and mass [kg].

DipoleMagneticField <: GeomagneticFieldType

Simplified geomagnetic dipole model.

Assumes Earth's magnetic field is a perfect dipole. Less accurate than IGRF but sufficient for preliminary analysis where uncertainties are high.

DragAstroModel{ST,AT,EoT,RT,PT} <: AbstractNonPotentialBasedForce

Atmospheric drag force model for spacecraft orbital dynamics.

This model computes the acceleration due to atmospheric drag acting on a spacecraft in Earth's atmosphere. The drag force is proportional to the atmospheric density and the square of the relative velocity between the spacecraft and the atmosphere.

Type Parameters

• ST <: AbstractSatelliteDragModel: Type of the satellite drag model
• AT <: AtmosphericModelType: Type of the atmospheric model
• EoT <: Union{EopIau1980,EopIau2000A}: Type of Earth Orientation Parameters
• RT <: Union{Nothing,AbstractVector}: Type for RTS (optional)
• PT <: Union{Nothing,AbstractMatrix}: Type for P matrix (optional)

Fields

• satellite_drag_model::ST: The satellite drag model for computing ballistic coefficient (Cd*A/m)
• atmosphere_model::AT: The atmospheric model for computing density (e.g., JB2008, JR1971, MSIS2000, HarrisPriester, HarrisPriesterModified, ExpAtmo)
• eop_data::EoT: Earth Orientation Parameters for coordinate transformations and atmospheric modeling
• rts::RT: Optional RTS parameter for atmospheric model (defaults to nothing)
• P::PT: Optional P matrix parameter for atmospheric model (defaults to nothing)

Example

# Create satellite drag model
sat_drag = CannonballDragModel(
    area = 10.0,        # [m²]
    drag_coeff = 2.2,   # dimensionless
    mass = 1000.0       # [kg]
)

# Create drag force model
drag_model = DragAstroModel(
    satellite_drag_model = sat_drag,
    atmosphere_model = JB2008(),
    eop_data = eop_data
)

See Also

• acceleration: Compute drag acceleration
• ballistic_coefficient: Compute ballistic coefficient
• Atmospheric density computation via SatelliteToolboxAtmosphericModels.jl

ExpAtmo <: AtmosphericModelType

Simple exponential atmospheric density model. No space weather indices required.

FixedChargeMassRatio{QT<:Number} <: AbstractSpacecraftChargeModel

Fixed charge-to-mass ratio model. The user directly specifies q/m [C/kg].

Typical values range from 1e-6 to 1e-3 C/kg for natural spacecraft charging in Low Earth Orbit, and up to ~0.03 C/kg for proposed active Lorentz-augmented orbits.

Fields

• q_over_m::QT: Charge-to-mass ratio [C/kg].

FixedThermalEmission{CT<:Number} <: AbstractSatelliteThermalModel

Fixed thermal emission coefficient model. The user directly specifies the effective thermal emission coefficient C_thm [m²/kg].

The thermal emission acceleration magnitude is computed as:

|a| = ν * C_thm * Ψ * (AU / d_sun)²

where ν is the shadow factor, Ψ is the solar radiation pressure at 1 AU, and d_sun is the Sun-spacecraft distance.

Fields

• thermal_emission_coeff::Number: The thermal emission coefficient [m²/kg].

FlatPlateThermalModel{AT,MT,EFT,EBT,ABT,XT,ET,CT} <: AbstractSatelliteThermalModel

Flat plate thermal emission model based on the solar panel thermal imbalance formulation of Duan & Hugentobler (2022) and Vigue et al. (1994).

Computes the thermal emission coefficient from physical surface properties using the equilibrium temperature assumption:

C_thm = (2/3) * (A/M) * [(ε_f - ξ·ε_b) / (ε_f + ξ·ε_b)] * (1 - η) * α

The resulting acceleration is directed along the Sun-spacecraft line.

Fields

• area::Number: Emitting surface area [m²]
• mass::Number: Spacecraft mass [kg]
• ε_front::Number: Front-side (Sun-facing) thermal emissivity
• ε_back::Number: Back-side thermal emissivity
• absorptivity::Number: Surface absorptivity for solar radiation
• ξ::Number: Temperature ratio factor between front and back sides (default: 1.0)
• η::Number: Electric conversion efficiency of solar panels (default: 0.0)
• thermal_emission_coeff::Number: Computed thermal emission coefficient [m²/kg]

FlatPlateThermalModel(; area, mass, ε_front, ε_back, absorptivity, ξ=1.0, η=0.0)

Construct a flat plate thermal model from physical surface properties.

The thermal emission coefficient is computed as:

C_thm = (2/3) * (A/M) * [(ε_f - ξ·ε_b) / (ε_f + ξ·ε_b)] * (1 - η) * α

Keyword Arguments

• area::Number: Emitting surface area [m²]
• mass::Number: Spacecraft mass [kg]
• ε_front::Number: Front-side (Sun-facing) thermal emissivity (0-1)
• ε_back::Number: Back-side thermal emissivity (0-1)
• absorptivity::Number: Surface absorptivity for solar radiation (0-1)
• ξ::Number: Temperature ratio factor between front and back sides (default: 1.0)
• η::Number: Electric conversion efficiency of solar panels (default: 0.0)

GeomagneticFieldType

Abstract type for geomagnetic field model selection.

Gravitational Harmonics Astro Model struct Contains information to compute the acceleration of a Gravitational Harmonics Model acting on a spacecraft.

Fields

• gravity_model::AbstractGravityModel: The gravitational potential model and coefficient data.
• eop_data::Union{EopIau1980,EopIau2000A}: The data compute the Earth's orientation.
• order::Int: The maximum order to compute the gravitational potential to, a value of -1 compute the maximum order of the supplied model. (Default=-1)
• degree::Int: The maximum degree to compute the gravitational potential to, a value of -1 compute the maximum degree of the supplied model. (Default=-1)

HarrisPriester <: AtmosphericModelType

Harris-Priester atmospheric density model. Valid from 100 to 1000 km altitude. Uses a static density table for mean solar activity conditions.

Fields

• n::Int: Cosine exponent for diurnal bulge modeling (2 ≤ n ≤ 6). Default: 4.

HarrisPriesterModified <: AtmosphericModelType

Modified Harris-Priester atmospheric density model (Hatten & Russell 2017). Provides C³ continuity, eliminates singularities, and uses cubic dependency on 81-day centered average F10.7 solar flux. Valid from 100 to 1000 km altitude.

Fields

• n::Number: Cosine exponent for diurnal bulge modeling. Default: 4.

IGRFField <: GeomagneticFieldType

International Geomagnetic Reference Field (IGRF) v14.

Uses spherical harmonic expansion up to degree 13 for high-fidelity geomagnetic field computation. Valid for dates between 1900 and 2035.

See SatelliteToolboxGeomagneticField.jl for details.

InertialFrame <: AbstractThrustFrame

Indicates that the thrust acceleration vector is expressed in the inertial (ECI/J2000) frame. No rotation is applied.

This is the default frame for LowThrustAstroModel.

JB2008 <: AtmosphericModelType

Jacchia-Bowman 2008 atmospheric density model. Valid from 100 to 1000 km altitude. Requires SpaceIndices.init() for automatic solar/geomagnetic index retrieval.

JR1971 <: AtmosphericModelType

Jacchia-Roberts 1971 atmospheric density model. Valid from 100 to 2500 km altitude. Requires SpaceIndices.init() for automatic solar/geomagnetic index retrieval.

Gravitational Keplerian Astro Model struct Contains information to compute the acceleration of a Gravitational Harmonics Model acting on a spacecraft.

Fields

• μ::Number: The gravitational potential constant of the central body.

LowThrustAstroModel{TM<:AbstractThrustModel,FR<:AbstractThrustFrame} <: AbstractNonPotentialBasedForce

Low-thrust propulsion force model for spacecraft orbital dynamics.

This model computes the acceleration due to low-thrust propulsion (e.g., electric/ion engines, solar sails, or any continuous thrust source). The thrust acceleration is determined by the underlying AbstractThrustModel, which defines the magnitude and direction of the thrust vector, and the AbstractThrustFrame, which specifies how the thrust components are oriented.

Type Parameters

• TM <: AbstractThrustModel: Type of the thrust model
• FR <: AbstractThrustFrame: Type of the reference frame

Fields

• thrust_model::TM: The thrust model defining the acceleration vector (e.g., ConstantCartesianThrust, ConstantTangentialThrust, StateThrustModel, PiecewiseConstantThrust)
• frame::FR: The reference frame in which the thrust model's output is expressed. Defaults to InertialFrame. Use RTNFrame or VNBFrame to specify thrust in orbital coordinates.

Constructors

LowThrustAstroModel(; thrust_model, frame=InertialFrame())

MSIS2000 <: AtmosphericModelType

NRLMSISE-00 atmospheric density model. Valid from 0 to 1000 km altitude. Requires SpaceIndices.init() for automatic solar/geomagnetic index retrieval.

MagneticFieldAstroModel{CT,GT,EoT,PT,DPT} <: AbstractNonPotentialBasedForce

Geomagnetic Lorentz force model for charged spacecraft.

A charged spacecraft moving through Earth's magnetic field experiences a perturbative Lorentz force. The magnetic field co-rotates with Earth, so the force depends on the spacecraft's velocity relative to the rotating field [1,2]:

𝐚 = (q/m) × (𝐯_rel × 𝐁)

where q/m is the charge-to-mass ratio, 𝐯_rel = 𝐯 - 𝛚×𝐫 is the velocity relative to the co-rotating magnetic field, and 𝐁 is the geomagnetic field vector.

The force is always perpendicular to both the relative velocity and the magnetic field, making it a non-dissipative perturbation that can modify orbital elements without adding or removing energy from the orbit (in the rotating frame) [3].

Type Parameters

• CT <: AbstractSpacecraftChargeModel: Charge-to-mass ratio model
• GT <: GeomagneticFieldType: Geomagnetic field model type
• EoT <: Union{EopIau1980,EopIau2000A}: Earth Orientation Parameters type
• PT <: Union{Nothing,AbstractMatrix}: Legendre polynomial buffer type
• DPT <: Union{Nothing,AbstractMatrix}: Legendre derivative buffer type

Fields

• spacecraft_charge_model::CT: Model providing the charge-to-mass ratio q/m [C/kg]
• geomagnetic_field_model::GT: Geomagnetic field model — IGRFField() or DipoleMagneticField()
• eop_data::EoT: Earth Orientation Parameters for J2000 ↔ ITRF transformations
• max_degree::Int: Maximum spherical harmonic degree for IGRF (ignored for dipole). Default: 13
• P::PT: Pre-allocated Legendre polynomial matrix for IGRF (reduces allocations). Default: nothing
• dP::DPT: Pre-allocated Legendre derivative matrix for IGRF (reduces allocations). Default: nothing

NoAtmosphere <: AtmosphericModelType

Sentinel type that returns zero density. Useful for disabling drag in a dynamics model without removing the DragAstroModel from the perturbation tuple.

NoIonosphere <: AbstractIonosphereModel

Null ionosphere model that returns zero density. Use to disable plasma drag while keeping the model in the perturbation tuple for interface consistency.

PiecewiseConstantThrust{N,TT<:Number,AT<:Number} <: AbstractThrustModel

A piecewise-constant thrust schedule defined by a sequence of time-tagged arcs.

Each arc specifies a constant 3-component acceleration vector that is active from its start time until the next arc begins. This is the natural representation for Sims-Flanagan trajectory segments, finite-burn maneuver schedules, and any thrust profile that changes discretely between arcs.

The 3-component vectors are interpreted in the frame specified by the parent LowThrustAstroModel.

Type Parameters

• N: Number of arcs (compile-time constant for type stability)
• TT <: Number: Element type of the time breakpoints
• AT <: Number: Element type of the acceleration components

Fields

• times::SVector{N,TT}: Sorted arc start times [s]. Must be strictly increasing.
• accelerations::SVector{N,SVector{3,AT}}: Acceleration vector for each arc [km/s²].

Before times[1] and after the last arc, the last arc's acceleration is held.

Constructors

PiecewiseConstantThrust(times, accelerations)

times is an AbstractVector of arc start times and accelerations is an AbstractVector of 3-component acceleration vectors (or SVector{3}).

PlasmaDragAstroModel{ST,IT,EoT} <: AbstractNonPotentialBasedForce

Ionospheric plasma drag force model for spacecraft in low-Earth orbit.

Plasma drag arises from direct collection of ionospheric ions (primarily O⁺ in the F2 region) and their momentum transfer to the spacecraft. The acceleration follows the same cannonball drag formulation as atmospheric drag, but uses the ion mass density from the ionosphere rather than the neutral atmospheric density:

𝐚 = -½ ⋅ BC_i ⋅ ρᵢ ⋅ |𝐯_app| ⋅ 𝐯_app

where BCi = C{D,i} A/m is the ion ballistic coefficient, ρᵢ is the ion mass density, and 𝐯_app is the spacecraft velocity relative to the co-rotating ionosphere [1, 3].

At LEO altitudes (300-600 km), plasma drag can contribute 5-35% of total aerodynamic force, with larger contributions at higher altitudes where neutral density decreases faster than ion density [1].

Type Parameters

• ST <: AbstractSatellitePlasmaDragModel: Satellite plasma drag model type
• IT <: AbstractIonosphereModel: Ionospheric density model type
• EoT <: Union{EopIau1980,EopIau2000A}: Earth Orientation Parameters type

Fields

• satellite_plasma_drag_model::ST: Model for computing ion ballistic coefficient
• ionosphere_model::IT: Model for computing ion mass density
• eop_data::EoT: Earth Orientation Parameters for coordinate transformations

RTNFrame <: AbstractThrustFrame

Indicates that the thrust acceleration vector is expressed in the RTN (Radial–Transverse–Normal) frame, also known as RIC (Radial–In-track–Cross-track).

The RTN basis vectors are defined as:

• R̂ (Radial): Along the position vector, away from the central body: r / |r|
• N̂ (Normal): Along the orbital angular momentum: (r × v) / |r × v|
• T̂ (Transverse): Completes the right-handed triad: N̂ × R̂

The acceleration components [a_R, a_T, a_N] are converted to inertial via:

a_inertial = a_R R̂ + a_T T̂ + a_N N̂

Relativity Astro Model struct Contains information to compute the acceleration of relativity acting on a spacecraft.

Fields

• central_body::ThirdBodyModel: The data to compute the central body's gravitational parameter.
• sun_body::ThirdBodyModel: The data to compute the Sun's position, velocity, and gravitational parameter.
• eop_data::Union{EopIau1980,EopIau2000A}: Earth Orientation Parameter data.
• c::Number: The speed of light [km/s].
• γ::Number: Post-Newtonian Parameterization parameter. γ=1 in General Relativity.
• β::Number: Post-Newtonian Parameterization parameter. β=1 in General Relativity.
• schwarzschild_effect::Bool: Include the Schwarzschild relativity effect.
• lense_thirring_effect::Bool: Include the Lense Thirring relativity effect.
• de_Sitter_effect::Bool: Include the De Sitter relativity effect.

RelativityModel(eop_data; kwargs...)

Convenience constructor that shares a single EOP dataset across all sub-models, avoiding redundant fetch_iers_eop() calls. Any of central_body, sun_body, or eop_data can be overridden via keyword arguments.

SRPAstroModel{ST,SDT,EoT,SMT,RST,ROT,PT,AUT} <: AbstractNonPotentialBasedForce

Solar Radiation Pressure (SRP) force model for spacecraft orbital dynamics.

This model computes the acceleration due to solar radiation pressure acting on a spacecraft. The force is proportional to the solar flux, inversely proportional to the square of the Sun-spacecraft distance, and includes shadow effects from the central body (typically Earth).

Type Parameters

• ST <: AbstractSatelliteSRPModel: Type of the satellite SRP model
• SDT <: ThirdBodyModel: Type of the Sun data model
• EoT <: Union{EopIau1980,EopIau2000A}: Type of Earth Orientation Parameters
• SMT <: ShadowModelType: Type of the shadow model

Fields

• satellite_srp_model::ST: The satellite SRP model for computing the reflectivity coefficient and cross-sectional area
• sun_data::SDT: The data model to compute the Sun's position and velocity
• eop_data::EoT: Earth Orientation Parameters for coordinate transformations
• shadow_model::SMT: Shadow model type (Conical, Cylindrical, etc.) - defaults to Conical()
• R_Sun::RST: Radius of the Sun [km] - defaults to R_SUN constant
• R_Occulting::ROT: Radius of the occulting body [km] - defaults to R_EARTH constant
• Ψ::PT: Solar flux constant at 1 AU [W/m²] - defaults to SOLAR_FLUX constant
• AU::AUT: Astronomical Unit [km] - defaults to ASTRONOMICAL_UNIT / 1E3

Example

# Create satellite SRP model
sat_srp = CannonballFixedSRP(
    radius = 1.0,           # [m]
    mass = 1000.0,          # [kg] 
    reflectivity_coeff = 1.3
)

# Create Sun position model
sun_model = ThirdBodyModel(body = SunBody(), eop_data = eop_data)

# Create SRP force model
srp_model = SRPAstroModel(
    satellite_srp_model = sat_srp,
    sun_data = sun_model,
    eop_data = eop_data,
    shadow_model = Conical()
)

See Also

• acceleration: Compute SRP acceleration
• reflectivity_ballistic_coefficient: Compute reflectivity coefficient
• Shadow modeling with Conical, Cylindrical, and other shadow types

SmoothedConical{CST,CTT} <: ShadowModelType

A smoothed (differentiable) shadow model based on Aziz et al. [2] that uses a logistic sigmoid function to provide a continuous transition between sunlight and shadow.

This model is fully compatible with automatic differentiation since it avoids branching (if/else) on shadow state, making it ideal for gradient-based optimization.

Fields

• cs::CST: Sharpness coefficient controlling the slope of the transition (default: 289.78)
• ct::CTT: Transition coefficient controlling where the curve is centered (default: 1.0)

The default values of cs = 289.78 and ct = 1.0 were found by Aziz et al. to minimize error relative to the discontinuous eclipse model for geocentric orbits.

SolidBodyTidesModel{N,TBT,K2T,K3T,KP2T,RET,BT} <: AbstractNonPotentialBasedForce

Solid body tides acceleration model following IERS Conventions (2010), Section 6.2.

Models the perturbation to the geopotential caused by the tidal deformation of a central body due to the gravitational attraction of tide-raising bodies. Implements Step 1 of the IERS procedure using frequency-independent Love numbers and the Legendre polynomial addition theorem for efficient computation directly in the inertial frame.

The formulation is general and works for any central body (Earth, Mars, etc.) with any number of tide-raising bodies. For Earth, the dominant tide raisers are the Moon and Sun; for Mars they would be the Sun, Phobos, and Deimos.

The tidal perturbation potential at satellite position $\mathbf{r}$ from a tide-raising body $j$ at degree $n$ is:

\[\Delta V_j^{(n)} = k_n \frac{GM_j R_e^{2n+1}}{r_j^{n+1} r^{n+1}} P_n(\cos\gamma_j)\]

where $k_n$ is the Love number, $R_e$ is the central body's equatorial radius, $r_j$ is the distance to body $j$, and $\gamma_j$ is the geocentric angle between the satellite and body $j$.

Fields

• tide_raising_bodies::NTuple{N,ThirdBodyModel}: Tuple of ThirdBodyModels for each tide-raising body. Gravitational parameters are obtained from each body's CelestialBody.
• k2::Number: Degree-2 Love number (default: 0.30190, IERS 2010 Table 6.3, anelastic).
• k3::Number: Degree-3 Love number (default: 0.093, IERS 2010 Table 6.3).
• k2_plus::Number: Degree-2 → degree-4 coupling Love number $k_2^{(+)}$ (default: -0.00089, IERS 2010 Table 6.3, anelastic m=0). Set to 0 to disable.
• R_e::Number: Central body equatorial radius in km (default: R_EARTH).
• include_degree_3::Bool: Include degree-3 tidal contribution (default: true).

SolidBodyTidesModel(eop_data; kwargs...)

Convenience constructor for Earth that creates Sun and Moon tide-raising body models sharing a single EOP dataset.

Arguments

• eop_data::Union{EopIau1980,EopIau2000A}: Earth Orientation Parameter data.

Keyword Arguments

All fields of SolidBodyTidesModel may be overridden.

StateChargeModel{FT<:Function} <: AbstractSpacecraftChargeModel

State-dependent charge-to-mass ratio model. The user provides a function f(u, p, t) -> q/m that returns the charge-to-mass ratio [C/kg] as a function of the current state, parameters, and time.

Fields

• charge_function::FT: Function (u, p, t) -> q/m [C/kg].

StateDragModel{IT<:Integer} <: AbstractSatelliteDragModel

Drag model that reads the ballistic coefficient from the state vector. Intended for estimation problems where the ballistic coefficient is a solve-for parameter.

Fields

• state_index::IT: Index into the state vector u where the ballistic coefficient is stored.

```

StatePlasmaDragModel <: AbstractSatellitePlasmaDragModel

Plasma drag model that reads the ion ballistic coefficient from the state vector. Intended for estimation problems where the ballistic coefficient is a solve-for parameter.

StateSRPModel{IT<:Integer} <: AbstractSatelliteSRPModel

SRP model that reads the reflectivity ballistic coefficient from the state vector. Intended for estimation problems where the coefficient is a solve-for parameter.

Fields

• state_index::IT: Index into the state vector u where the reflectivity ballistic coefficient is stored.

```

StateThrustModel{F} <: AbstractThrustModel

A user-defined thrust model where the thrust acceleration vector is computed by a provided function.

This is the most flexible thrust model, allowing arbitrary state- and time-dependent thrust profiles such as feedback control laws, scheduled thrust arcs, or optimization-based thrust histories.

The 3-component return value is interpreted in the frame specified by the parent LowThrustAstroModel.

Type Parameters

• F: Type of the callable (function, functor, or closure)

Fields

• f::F: Callable with signature f(u, p, t) → SVector{3} returning thrust acceleration [km/s²]

ThermalEmissionAstroModel{TT,SDT,EoT,SMT,RST,ROT,PT,AUT} <: AbstractNonPotentialBasedForce

Spacecraft thermal emission (thermal re-radiation) force model.

When a spacecraft absorbs solar radiation, the surface heats up and re-emits thermal infrared radiation. If the emission is anisotropic (e.g., due to different emissivities on front and back surfaces of solar panels), the net photon recoil produces a small but persistent perturbative acceleration. This effect is most significant during eclipse season transitions for GNSS satellites, where it can cause orbit errors of 1-3 cm [1].

In the steady-state (equilibrium) approximation, the thermal emission acceleration has the same functional form as solar radiation pressure:

𝐚 = ν * C_thm * Ψ * (AU / d_sun)² * d̂_sun

where ν is the shadow factor, Cthm is the thermal emission coefficient [m²/kg], Ψ is the solar flux at 1 AU [N/m²], dsun is the Sun-spacecraft distance, and d̂_sun is the unit vector from the Sun to the spacecraft.

For a flat plate with front emissivity εf and back emissivity εb, the thermal emission coefficient in equilibrium is [1,2]:

C_thm = (2/3) * (A/M) * [(ε_f - ξ·ε_b) / (ε_f + ξ·ε_b)] * (1 - η) * α

Type Parameters

• TT <: AbstractSatelliteThermalModel: Type of the satellite thermal model
• SDT <: ThirdBodyModel: Type of the Sun data model
• EoT <: Union{EopIau1980,EopIau2000A}: Type of Earth Orientation Parameters
• SMT <: ShadowModelType: Type of the shadow model

Fields

• satellite_thermal_model::TT: Satellite thermal model providing the emission coefficient
• sun_data::SDT: Model to compute the Sun's position and velocity
• eop_data::EoT: Earth Orientation Parameters for coordinate transformations
• shadow_model::SMT: Shadow model type — defaults to Conical()
• R_Sun::RST: Radius of the Sun [km] — defaults to R_SUN
• R_Occulting::ROT: Radius of the occulting body [km] — defaults to R_EARTH
• Ψ::PT: Solar flux constant at 1 AU [N/m²] — defaults to SOLAR_FLUX
• AU::AUT: Astronomical Unit [km] — defaults to ASTRONOMICAL_UNIT / 1E3

Third Body Model Astro Model struct Contains information to compute the acceleration of a third body force acting on a spacecraft.

Fields

• body::CelestialBody: Celestial body acting on the craft.
• ephem_type::AbstractEphemerisType: Ephemeris type used to compute body's position. Options are currently Vallado().

Convenience to compute the ephemeris position of a CelestialBody in a ThirdBodyModel Wraps get_position().

Arguments

• time::Number: Current time of the simulation in seconds.

Returns

• body_position: SVector{3}: The 3-dimensional third body position in the J2000 frame.

Uniform Albedo Model struct Simple uniform albedo model with constant reflection and emission coefficients.

Fields

• visible_albedo::Number: Visible light albedo coefficient (0-1)
• infrared_emissivity::Number: Infrared emissivity coefficient (0-1)

Default Values

The default values (visiblealbedo=0.3, infraredemissivity=0.7) represent commonly accepted Earth-averaged values used in astrodynamics literature. These are based on Earth radiation budget studies and are referenced in works such as Knocke et al. (1988), Borderies & Longaretti (1990), and other albedo radiation pressure studies. The 0.3 value represents approximately 30% of incoming solar radiation being reflected, while 0.7 represents Earth's average thermal emission characteristics.

VNBFrame <: AbstractThrustFrame

Indicates that the thrust acceleration vector is expressed in the VNB (Velocity–Normal–Binormal) frame.

The VNB basis vectors are defined as:

• V̂ (Velocity): Along the velocity vector: v / |v|
• N̂ (Normal): In the orbital plane, perpendicular to V̂: B̂ × V̂
• B̂ (Binormal): Along the orbital angular momentum: (r × v) / |r × v|

The acceleration components [a_V, a_N, a_B] are converted to inertial via:

a_inertial = a_V V̂ + a_N N̂ + a_B B̂

Degree-2 tidal acceleration from a single perturbing body.

Gradient of the degree-2 tidal perturbation potential: $\Delta V^{(2)} = k_2 GM_j R_e^5 / (r_j^3 r^3) \cdot P_2(\cos\gamma)$

Degree-3 tidal acceleration from a single perturbing body.

Gradient of the degree-3 tidal perturbation potential: $\Delta V^{(3)} = k_3 GM_j R_e^7 / (r_j^4 r^4) \cdot P_3(\cos\gamma)$

acceleration(u::AbstractVector, p::ComponentVector, t::Number, albedo_model::AlbedoAstroModel)

Computes the albedo acceleration acting on a spacecraft given an albedo model and current state and parameters of an object.

Arguments

• u::AbstractVector: Current State of the simulation.
• p::ComponentVector: Current parameters of the simulation.
• t::Number: Current time of the simulation.
• albedo_model::AlbedoAstroModel: Albedo model struct containing the relevant information to compute the acceleration.

Returns

• acceleration: SVector{3}: The 3-dimensional albedo acceleration acting on the spacecraft.

acceleration(u::AbstractVector, p::ComponentVector, t::Number, drag_model::DragAstroModel)

Computes the drag acceleration acting on a spacecraft given a drag model and current state and parameters of an object.

Arguments

• u::AbstractVector: Current State of the simulation.
• p::ComponentVector: Current parameters of the simulation.
• t::Number: Current time of the simulation.
• drag_model::DragAstroModel: Drag model struct containing the relevant information to compute the acceleration.

Returns

• acceleration: SVector{3}: The 3-dimensional drag acceleration acting on the spacecraft.

acceleration(u::AbstractVector, p::ComponentVector, t::Number, grav_model::GravityHarmonicsAstroModel)

Computes the gravitational acceleration acting on a spacecraft given a gravity model and current state and parameters of an object.

Arguments

• u::AbstractVector: Current State of the simulation.
• p::ComponentVector: Current parameters of the simulation.
• t::Number: Current time of the simulation.
• gravity_model::GravityHarmonicsAstroModel: Gravity model struct containing the relevant information to compute the acceleration.

Returns

• acceleration: SVector{3}: The 3-dimensional gravity acceleration acting on the spacecraft.

acceleration(u::AbstractVector, p::ComponentVector, t::Number, grav_model::KeplerianGravityAstroModel)

Computes the gravitational acceleration acting on a spacecraft given a gravity model and current state and parameters of an object.

Arguments

• u::AbstractVector: Current State of the simulation.
• p::ComponentVector: Current parameters of the simulation.
• t::Number: Current time of the simulation.
• gravity_model::KeplerianGravityAstroModel: Gravity model struct containing the relevant information to compute the acceleration.

Returns

• acceleration: SVector{3}: The 3-dimensional gravity acceleration acting on the spacecraft.

acceleration(u::AbstractVector, p::ComponentVector, t::Number, lt_model::LowThrustAstroModel)

Computes the low-thrust acceleration acting on a spacecraft in the inertial frame.

The thrust model first produces a 3-component acceleration in the model's reference frame, which is then transformed to the inertial frame using transform_to_inertial.

Arguments

• u::AbstractVector: Current state of the simulation [x, y, z, vx, vy, vz] in km and km/s.
• p::ComponentVector: Current parameters of the simulation.
• t::Number: Current time of the simulation [s].
• lt_model::LowThrustAstroModel: Low-thrust model containing the thrust configuration.

Returns

• SVector{3}: The 3-dimensional thrust acceleration in the inertial frame [km/s²].

acceleration(u::AbstractVector, p::ComponentVector, t::Number, model::MagneticFieldAstroModel)

Compute the geomagnetic Lorentz force acceleration on a charged spacecraft.

Arguments

• u::AbstractVector: Current state [rx, ry, rz, vx, vy, vz] in km and km/s (J2000).
• p::ComponentVector: Parameters including JD (Julian Date at epoch).
• t::Number: Elapsed time since epoch [s].
• model::MagneticFieldAstroModel: Magnetic field force model configuration.

Returns

• SVector{3}: Lorentz force acceleration [km/s²] in J2000 frame.

acceleration(u, p, t, model::PlasmaDragAstroModel) -> SVector{3}

Compute the plasma drag acceleration on a spacecraft.

Arguments

• u::AbstractVector: State vector [r; v] in J2000 ECI [km; km/s].
• p::ComponentVector: Parameters; must include p.JD (Julian Date at epoch).
• t::Number: Elapsed time since epoch [s].
• model::PlasmaDragAstroModel: Plasma drag force model.

Returns

• SVector{3}: Plasma drag acceleration in J2000 ECI [km/s²].

acceleration(u::AbstractVector, p::ComponentVector, t::Number, srp_model::SRPAstroModel)

Computes the srp acceleration acting on a spacecraft given a srp model and current state and parameters of an object.

Arguments

• u::AbstractVector: Current State of the simulation.
• p::ComponentVector: Current parameters of the simulation.
• t::Number: Current time of the simulation.
• srp_model::SRPAstroModel: SRP model struct containing the relevant information to compute the acceleration.

Returns

• acceleration: SVector{3}: The 3-dimensional srp acceleration acting on the spacecraft.

acceleration(u, p, t, model::SolidBodyTidesModel) -> SVector{3}

Compute the solid body tides acceleration perturbation on a satellite.

Arguments

• u::AbstractVector: Spacecraft state [rx, ry, rz, vx, vy, vz] in km and km/s (J2000 ECI).
• p::ComponentVector: Simulation parameters (must contain JD).
• t::Number: Elapsed time since epoch [s].
• model::SolidBodyTidesModel: Solid body tides model configuration.

Returns

• SVector{3}: Acceleration perturbation [km/s²] in J2000 ECI.

acceleration(u::AbstractVector, p::ComponentVector, t::Number, model::ThermalEmissionAstroModel)

Compute the acceleration from spacecraft thermal emission given a thermal model and the current state and parameters of the spacecraft.

Arguments

• u::AbstractVector: Current state of the simulation [km, km/s].
• p::ComponentVector: Current parameters of the simulation.
• t::Number: Current time of the simulation [s].
• model::ThermalEmissionAstroModel: Thermal emission model.

Returns

• SVector{3}: The 3-dimensional thermal emission acceleration [km/s²].

acceleration(u::AbstractVector, p::ComponentVector, t::Number, third_body_model::ThirdBodyModel)

Computes the third-body gravitational acceleration acting on a spacecraft given a third body model and current state and parameters of an object.

Arguments

• u::AbstractVector: Current State of the simulation.
• p::ComponentVector: Current parameters of the simulation.
• t::Number: Current time of the simulation.
• third_body_model::ThirdBodyModel: Third body model struct containing the relevant information to compute the acceleration.

Returns

• acceleration: SVector{3}: The 3-dimensional third-body acceleration acting on the spacecraft.

acceleration(u::AbstractVector, p::ComponentVector, t::Number, relativity_model::RelativityModel)

Computes the relativistic acceleration acting on a spacecraft given a relativity model and current state and parameters of an object.

Arguments

• u::AbstractVector: Current State of the simulation.
• p::ComponentVector: Current parameters of the simulation.
• t::Number: Current time of the simulation.
• relativity_model::RelativityModel: Relativity model struct containing the relevant information to compute the acceleration.

Returns

• acceleration: SVector{3}: The 3-dimensional relativity acceleration acting on the spacecraft.

albedo_accel(
    u::AbstractVector, 
    sun_pos::AbstractVector, 
    RC::Number, 
    current_time::Number,
    body_albedo_model::AbstractAlbedoModel,
    eop_data::EopIau1980,
    surface_positions_ecef::AbstractVector,
    weights::AbstractVector;
    solar_irradiance::Number=SOLAR_IRRADIANCE,
    speed_of_light::Number=SPEED_OF_LIGHT,
    AU::Number=ASTRONOMICAL_UNIT / 1E3
)

Compute the acceleration from Earth albedo radiation pressure using Lebedev quadrature.

Earth albedo radiation pressure arises from two sources:

1. Solar radiation reflected by Earth's surface (shortwave, albedo component)
2. Thermal radiation emitted by Earth (longwave, infrared component)

The total acceleration is computed by integrating over the Earth's surface visible to the satellite (field of view), considering both reflected and emitted radiation. Surface element positions are pre-computed in ECEF at construction time and rotated to ECI once per evaluation step, using dot products for all angle computations.

Mathematical formulation based on Knocke et al. (1988): aalbedo = RC * ∫∫ [FSW + F_LW] * cos(α) * dΩ / (π * c * r²) / 1E3

Arguments

• u::AbstractVector: The current state of the spacecraft in the central body inertial frame [km, km/s].
• sun_pos::AbstractVector: The current position of the Sun [km].
• RC::Number: The reflectivity ballistic coefficient of the satellite [m²/kg].
• current_time::Number: The current Julian date.
• body_albedo_model::AbstractAlbedoModel: Earth albedo radiation model.
• eop_data::EopIau1980: Earth orientation parameters.
• surface_positions_ecef::AbstractVector: Pre-computed ECEF surface element positions [km].
• weights::AbstractVector: Pre-computed scaled integration weights.

Keyword Arguments

• solar_irradiance::Number: Solar irradiance at 1 AU [W/m²]. Default: SOLAR_IRRADIANCE.
• speed_of_light::Number: Speed of light [km/s]. Default: SPEED_OF_LIGHT.
• AU::Number: Astronomical Unit [km]. Default: ASTRONOMICAL_UNIT / 1E3.

Returns

• SVector{3}: Inertial acceleration from albedo radiation pressure [km/s²].

angle_between_vectors(
    v1::AbstractVector{T1}, v2::AbstractVector{T2}
) where {T1<:Number,T2<:Number}

Computes the angle between two vectors in a more numerically stable way than dot product.

Arguments

-v1::AbstractVector{<:Number}: The first vector of the computation -v2::AbstractVector{<:Number}: The second vector of the computation

Returns

-angle::Number: The angle between the two vectors

ballistic_coefficient( u::AbstractVector, p::ComponentVector, t::Number, model::CannonballFixedDrag)

Returns the ballistic coefficient for a drag model given the model and current state of the simulation.

Arguments

- `u::AbstractVector`: The current state of the simulation.
- `p::ComponentVector`: The parameters of the simulation.
- `t::Number`: The current time of the simulation.
- `model::CannonballFixedDrag`: Drag model for the spacecraft.

Returns

-`ballistic_coeff::Number`: The current ballistic coefficient of the spacecraft.

ballistic_coefficient(u, p, t, model::StateDragModel)

Returns the ballistic coefficient from the state vector at model.state_index.

build_dynamics_model(u::AbstractVector, p::ComponentVector, t::Number, models::CentralBodyDynamicsModel)

Compute the total acceleration acting on a spacecraft using a CentralBodyDynamicsModel.

This function efficiently combines the central body gravity acceleration with all perturbing accelerations to produce the total acceleration vector. It is optimized for use in ODE integration routines and provides better performance than manually summing individual force model accelerations.

Arguments

• u::AbstractVector: Current spacecraft state vector [rx, ry, rz, vx, vy, vz] in km and km/s
• p::ComponentVector: Simulation parameters including time references (JD), spacecraft properties, etc.
• t::Number: Current simulation time (typically seconds since epoch)
• models::CentralBodyDynamicsModel: Combined dynamics model containing gravity and perturbing forces

Returns

• SVector{3}: The 3-dimensional total acceleration vector [ax, ay, az] in km/s²

Performance Notes

The function is marked @inline for performance and uses compile-time optimizations through the tuple-based storage of perturbing models. This approach provides:

• Zero-cost abstraction over manual force summation
• Type stability for all force combinations
• Efficient memory access patterns

Example

# Create dynamics model
dynamics = CentralBodyDynamicsModel(gravity_model, (drag_model, srp_model))

# Compute acceleration at current state
state = [6678.137, 0.0, 0.0, 0.0, 7.66, 0.0]  # km, km/s
params = ComponentVector(JD = 2.460310e6)  # Julian Date
time = 0.0

total_accel = build_dynamics_model(state, params, time, dynamics)

See Also

• CentralBodyDynamicsModel: Dynamics model structure
• acceleration: Individual force model interface
• sum_accelerations: Internal acceleration summation function

charge_mass_ratio(u, p, t, model::FixedChargeMassRatio)

Return the fixed charge-to-mass ratio [C/kg].

charge_mass_ratio(u, p, t, model::StateChargeModel)

Evaluate the state-dependent charge-to-mass ratio function [C/kg].

compute_density(JD::Number, u::AbstractVector, eop_data::EopIau1980, AtmosphereType::AtmosphericModelType)

Computes the atmospheric density at a point given the date, position, eop_data, and atmosphere type.

Arguments

• JD::Number: The current time of the simulation in Julian days.
• u::AbstractVector: The current state of the simulation.
• eop_data::EopIau1980: The earth orientation parameters.
• AtmosphereType::AtmosphericModelType: The type of atmospheric model used to compute the density. Available options are Jacchia-Bowman 2008 (JB2008), Jacchia-Roberts 1971 (JR1971), NRL MSIS 2000 (MSIS2000), Exponential (ExpAtmo), Harris-Priester (HarrisPriester), Modified Harris-Priester (HarrisPriesterModified), and no atmosphere (NoAtmosphere).

Returns

• rho::Number: The density of the atmosphere at the provided time and point [kg/m^3].

compute_earth_radiation_fluxes(
    body_albedo_model::UniformAlbedoModel, 
    cos_solar_zenith::Number, 
    irradiance_at_earth::Number
)

Compute the total radiation flux (shortwave + longwave) from a surface element.

Arguments

• body_albedo_model::UniformAlbedoModel: Earth albedo model with uniform coefficients.
• cos_solar_zenith::Number: Cosine of solar zenith angle at surface element.
• irradiance_at_earth::Number: Solar irradiance at Earth's distance [W/m²].

Returns

• Number: Total outgoing flux [W/m²].

compute_ion_density(JD, u, eop_data, model::ChapmanIonosphere)

Compute the ion mass density at the spacecraft position using a Chapman layer profile.

The geodetic altitude is computed from the J2000 state vector via ECEF transformation, and the Chapman function is evaluated for the F2 layer. Density is zero above 1500 km where the ionosphere becomes negligible.

Arguments

• JD::Number: Current Julian Date.
• u::AbstractVector: Spacecraft state vector [r; v] in J2000 ECI [km; km/s].
• eop_data: Earth Orientation Parameters for coordinate transformations.
• model::ChapmanIonosphere: Chapman layer model parameters.

Returns

• rho_i::Number: Ion mass density [kg/m³].

compute_magnetic_field(r_ecef::SVector{3}, JD::Number, model::MagneticFieldAstroModel{CT,DipoleMagneticField}) where CT

Compute the geomagnetic field vector in ECEF [T] using the simplified dipole model.

compute_magnetic_field(r_ecef::SVector{3}, JD::Number, model::MagneticFieldAstroModel{CT,IGRFField}) where CT

Compute the geomagnetic field vector in ECEF [T] using the IGRF model.

Converts ECEF Cartesian coordinates to geocentric spherical, calls the IGRF model to get B in the geocentric NED frame, then rotates back to ECEF.

current_jd(p, t)

Compute the current Julian Date from parameters p and elapsed time t [s].

de_Sitter_acceleration(
    u::AbstractVector,
    r_sun::AbstractVector,
    v_sun::AbstractVector,
    μ_Sun::Number;
    c::Number=SPEED_OF_LIGHT,
    γ::Number=1.0,
)

Computes the relativity acceleration acting on a spacecraft given a relativity model and current state and parameters of an object.

Arguments

• u::AbstractVector: Current State of the simulation.
• r_sun::AbstractVector: The position of the sun in the Earth inertial frame.
• v_sun::AbstractVector: The velocity of the sun in the Earth inertial frame.
• μ_Sun::Number: Gravitation Parameter of the Sun. [km^3/s^2]
• c::Number: Speed of Light [km/s]
• γ::Number: Post-Newtonian Parameterization parameter. γ=1 in General Relativity.

Returns

• de_Sitter_acceleration: SVector{3}: The 3-dimensional de Sitter acceleration acting on the spacecraft.

drag_accel(u::AbstractVector, rho::Number, BC::Number, ω_vec::AbstractVector, t::Number, [DragModel]) -> SVector{3}{Number}

Compute the Acceleration Atmospheric Drag

The atmosphere is treated as a solid revolving with the Earth and the apparent velocity of the satellite is computed using the transport theorem

            𝐯_app = 𝐯 - 𝛚 x 𝐫

The acceleration from drag is then computed with a cannonball model as

            𝐚 = 1/2 * ρ * BC * |𝐯_app|₂^2 * v̂

Arguments

• u::AbstractVector: The current state of the spacecraft in the central body's inertial frame.
• rho::Number: Atmospheric density at (t, u) [kg/m^3].
• BC::Number: The ballistic coefficient of the satellite – (area/mass) * drag coefficient [kg/m^2].
• ω_vec::AbstractVector: The angular velocity vector of Earth. Typically approximated as [0.0; 0.0; ω_Earth]
• t::Number: Current time of the simulation.

Returns

• SVector{3}{Number}: Inertial acceleration from drag

err_iau2006(JD_TT::Number) -> Float64

Compute the Earth Rotation Rate (time derivative of the Greenwich Mean Sidereal Time) according to the IAU 2006 Conventions.

Args

• JD_TT: Terrestrial Time (TT) Julian date.

Returns

• Float64: The Earth Rotation Rate [rad/s].

Computes the position of the celestial body using Vallado's ephemeris

Arguments

• ephem_type::Vallado: Ephemeris type used to compute body's position.
• body::CelestialBody: Celestial body acting on the craft.
• time::Number: Current time of the simulation in Julian days.

Returns

• body_position::SVector{3}: The 3-dimensional third body position in the J2000 frame [m].

Computes the velocity of the celestial body using Vallado's ephemeris

Arguments

• ephem_type::Vallado: Ephemeris type used to compute body's velocity.
• body::CelestialBody: Celestial body acting on the craft.
• time::Number: Current time of the simulation in Julian days.

Returns

• body_velocity::SVector{3}: The 3-dimensional third body velocity in the J2000 frame.

ion_ballistic_coefficient(u, p, t, model::CannonballFixedPlasmaDrag)

Returns the fixed ion ballistic coefficient C_D,i × A / m [m²/kg].

ion_ballistic_coefficient(u, p, t, model::StatePlasmaDragModel)

Returns the ion ballistic coefficient from the state vector at model.state_index.

lense_thirring_acceleration(
    u::AbstractVector,
    μ_body::Number,
    J::AbstractVector;
    c::Number=SPEED_OF_LIGHT,
    γ::Number=1.0,
)

Computes the lense thirring relativity acceleration acting on a spacecraft given a relativity model and current state and parameters of an object.

Arguments

• u::AbstractVector: Current State of the simulation.
• μ_body::Number: Gravitation Parameter of the central body.
• J::AbstractVector: Angular momentum vector per unit mass of the central body. [km^3/s^2]
• c::Number: Speed of Light [km/s]
• γ::Number: Post-Newtonian Parameterization parameter. γ=1 in General Relativity.

Returns

• lense_thirring_acceleration: SVector{3}: The 3-dimensional lense thirring acceleration acting on the spacecraft.

magnetic_field_accel(u::AbstractVector, q_over_m::Number, B_eci::SVector{3}) -> SVector{3}

Compute the geomagnetic Lorentz force acceleration given the state, charge-to-mass ratio, and magnetic field vector in the inertial frame.

The Lorentz force on a charged body moving through a magnetic field is [1,2]:

𝐅 = q(𝐯_rel × 𝐁)

where 𝐯_rel = 𝐯 - 𝛚×𝐫 is the velocity relative to the co-rotating magnetic field. The acceleration is:

𝐚 = (q/m)(𝐯_rel × 𝐁)

The magnetic field co-rotates with the Earth, so the relative velocity accounts for Earth's rotation via the transport theorem, identical to the apparent velocity used in atmospheric drag computations.

Arguments

• u::AbstractVector: Spacecraft state [rx, ry, rz, vx, vy, vz] in km and km/s (J2000).
• q_over_m::Number: Charge-to-mass ratio [C/kg].
• B_eci::SVector{3}: Geomagnetic field vector in J2000 frame [T].

Returns

• SVector{3}: Lorentz force acceleration [km/s²] in J2000 frame.

plasma_drag_accel(u, rho_i, BC_i, ω_vec) -> SVector{3}

Low-level computation of plasma drag acceleration.

The ionospheric plasma is treated as co-rotating with the Earth. The apparent velocity of the spacecraft relative to the plasma is:

𝐯_app = 𝐯 - 𝛚 × 𝐫

and the plasma drag acceleration is:

𝐚 = -½ ⋅ BC_i ⋅ ρᵢ ⋅ |𝐯_app| ⋅ 𝐯_app

The factor of 1E3 converts from m/s² to km/s² (density is in kg/m³, velocity in km/s).

Arguments

• u::AbstractVector: State [r; v] in J2000 ECI [km; km/s].
• rho_i::Number: Ion mass density [kg/m³].
• BC_i::Number: Ion ballistic coefficient C_{D,i} A/m [m²/kg].
• ω_vec::AbstractVector: Earth angular velocity vector [rad/s].

Returns

• SVector{3}: Plasma drag acceleration [km/s²].

potential(u::AbstractVector, p::ComponentVector, t::Number, grav_model::GravityHarmonicsAstroModel)

Computes the gravitational potential acting on a spacecraft given a gravity model and current state and parameters of an object.

Arguments

• u::AbstractVector: Current State of the simulation.
• p::ComponentVector: Current parameters of the simulation.
• t::Number: Current time of the simulation.
• gravity_model::GravityHarmonicsAstroModel: Gravity model struct containing the relevant information to compute the potential.

Returns

• potential: Number: The gravitational potential acting on the spacecraft.

potential(u::AbstractVector, p::ComponentVector, t::Number, grav_model::KeplerianGravityAstroModel)

Computes the gravitational potential acting on a spacecraft given a gravity model and current state and parameters of an object.

Arguments

• u::AbstractVector: Current State of the simulation.
• p::ComponentVector: Current parameters of the simulation.
• t::Number: Current time of the simulation.
• gravity_model::KeplerianGravityAstroModel: Gravity model struct containing the relevant information to compute the potential.

Returns

• potential: Number: The gravitational potential acting on the spacecraft.

potential_time_derivative(u::AbstractVector, p::ComponentVector, t::Number, grav_model::GravityHarmonicsAstroModel)

Computes the time derivative of the gravitational potential acting on a spacecraft given a gravity model and current state and parameters of an object. Based on the IAU 2006 precession-nutation model, implementation from [1].

[1] Amato, Davide. "THALASSA: Orbit propagator for near-Earth and cislunar space." Astrophysics Source Code Library (2019): ascl-1905.

Arguments

• u::AbstractVector: Current State of the simulation.
• p::ComponentVector: Current parameters of the simulation.
• t::Number: Current time of the simulation.
• gravity_model::GravityHarmonicsAstroModel: Gravity model struct containing the relevant information to compute the potential time derivative.

Returns

• potential_time_derivative: Number: The time derivative of the gravitational potential acting on the spacecraft.

potential_time_derivative(u::AbstractVector, p::ComponentVector, t::Number, grav_model::KeplerianGravityAstroModel)

Computes the time derivative of the gravitational potential acting on a spacecraft given a gravity model and current state and parameters of an object.

Arguments

• u::AbstractVector: Current State of the simulation.
• p::ComponentVector: Current parameters of the simulation.
• t::Number: Current time of the simulation.
• gravity_model::KeplerianGravityAstroModel: Gravity model struct containing the relevant information to compute the potential time derivative.

Returns

• potential_time_derivative: Number: The time derivative of the gravitational potential acting on the spacecraft.

reflectivityballisticcoefficient( u::AbstractVector, p::ComponentVector, t::Number, model::CannonballFixedSRP)

Returns the ballistic coefficient for a SRP model given the model and current state of the simulation.

Arguments

- `u::AbstractVector`: The current state of the simulation.
- `p::ComponentVector`: The parameters of the simulation.
- `t::Number`: The current time of the simulation.
- `model::CannonballFixedSRP`: SRP model for the spacecraft.

Returns

-`ballistic_coeff::Number`: The current ballistic coefficient of the spacecraft.

reflectivity_ballistic_coefficient(u, p, t, model::StateSRPModel)

Returns the reflectivity ballistic coefficient from the state vector at model.state_index.

relativity_accel(
    u::AbstractVector,
    r_sun::AbstractVector,
    v_sun::AbstractVector,
    μ_body::Number,
    μ_Sun::Number,
    J::AbstractVector;
    c::Number=SPEED_OF_LIGHT,
    γ::Number=1.0,
    β::Number=1.0,
    schwarzschild_effect::Bool=true,
    lense_thirring_effect::Bool=true,
    de_Sitter_effect::Bool=true,
)

Computes the relativity acceleration acting on a spacecraft given a relativity model and current state and parameters of an object.

Arguments

• u::AbstractVector: Current State of the simulation.
• r_sun::AbstractVector: The position of the sun in the Earth inertial frame.
• v_sun::AbstractVector: The velocity of the sun in the Earth inertial frame.
• μ_body::Number: Gravitation Parameter of the central body.
• μ_Sun::Number: Gravitation Parameter of the Sun. [km^3/s^2]
• J::AbstractVector: Angular momentum vector per unit mass of the central body. [km^3/s^2]
• c::Number: Speed of Light [km/s]
• γ::Number: Post-Newtonian Parameterization parameter. γ=1 in General Relativity.
• β::Number: Post-Newtonian Parameterization parameter. β=1 in General Relativity.
• schwarzschild_effect::Bool: Include the Schwarzschild relativity effect.
• lense_thirring_effect::Bool: Include the Lense Thirring relativity effect.
• de_Sitter_effect::Bool: Include the De Sitter relativity effect.

Returns

• acceleration: SVector{3}: The 3-dimensional relativity acceleration acting on the spacecraft.

schwarzschild_acceleration(
    u::AbstractVector, μ_body::Number; c::Number=SPEED_OF_LIGHT, γ::Number=1.0, β::Number=1.0
)

Computes the relativity acceleration acting on a spacecraft given a relativity model and current state and parameters of an object.

Arguments

• u::AbstractVector: Current State of the simulation.
• μ_body::Number: Gravitation Parameter of the central body.
• c::Number: Speed of Light [km/s]
• γ::Number: Post-Newtonian Parameterization parameter. γ=1 in General Relativity.
• β::Number: Post-Newtonian Parameterization parameter. β=1 in General Relativity.

Returns

• schwarzschild_acceleration: SVector{3}: The 3-dimensional schwarzschild acceleration acting on the spacecraft.

shadow_model(sat_pos::AbstractVector, sun_pos::AbstractVector, R_Sun::Number, R_Occulting::Number, t::Number, ShadowModel::ShadowModelType)

Computes the Lighting Factor of the Sun occur from the Umbra and Prenumbra of Earth's Shadow

Arguments

• sat_pos::AbstractVector: The current satellite position.
• sun_pos::AbstractVector: The current Sun position.
• R_Sun::Number: The radius of the Sun.
• R_Occulting::Number: The radius of the Occulting Body.
• ShadowModel::ShadowModelType: The Earth shadow model to use. Current Options – Cylindrical, Conical, SmoothedConical, NoShadow

Returns

• Number: Shadow factor between 0.0 (full shadow) and 1.0 (full sunlight)

shadow_model(sat_pos, sun_pos, ::SmoothedConical; R_Sun, R_Occulting)

Compute the sunlight fraction using the smoothed eclipse model from Aziz et al. [2].

The model computes apparent angular semi-diameters of the Sun and occulting body as seen from the spacecraft, then applies a logistic sigmoid to produce a smooth, differentiable transition between full sunlight (γ ≈ 1) and full shadow (γ ≈ 0).

Equations (Aziz et al. Eq. 16-19)

• a_SR = asin(R_Sun / ||r_sun/sc||) — apparent angular radius of the Sun
• a_BR = asin(R_B / ||r_B/sc||) — apparent angular radius of the occulting body
• a_D = acos(r̂_B/sc · r̂_sun/sc) — angular separation
• γ = 1 / (1 + exp(-cs * [a_D - ct * (a_SR + a_BR)])) — sunlight fraction

solid_body_tides_accel(
    u, bodies;
    k2=0.30190, k3=0.093,
    R_e=R_EARTH, include_degree_3=true,
) -> SVector{3}

Compute the solid body tides acceleration on a satellite from tide-raising body data.

Uses the gradient of the tidal perturbation potential derived from IERS Conventions (2010) Equation 6.6 with the Legendre polynomial addition theorem.

Degree 2 acceleration from body $j$:

\[\mathbf{a}_j^{(2)} = \frac{3 k_2 GM_j R_e^5}{r_j^3 r^4} \left[\frac{1 - 5\xi^2}{2}\hat{\mathbf{r}} + \xi\hat{\mathbf{r}}_j\right]\]

Degree 3 acceleration from body $j$:

\[\mathbf{a}_j^{(3)} = \frac{k_3 GM_j R_e^7}{2 r_j^4 r^5} \left[(15\xi - 35\xi^3)\hat{\mathbf{r}} + (15\xi^2 - 3)\hat{\mathbf{r}}_j\right]\]

where $\xi = \hat{\mathbf{r}} \cdot \hat{\mathbf{r}}_j$.

Arguments

• u::AbstractVector: Spacecraft state vector [rx, ry, rz, vx, vy, vz] in km, km/s.
• bodies::Tuple: Tuple of (r_body, μ_body) pairs, where r_body is the body position [km] and μ_body is the gravitational parameter [km³/s²].

Keyword Arguments

• k2::Number: Degree-2 Love number (default: 0.30190).
• k3::Number: Degree-3 Love number (default: 0.093).
• R_e::Number: Central body equatorial radius in km (default: R_EARTH).
• include_degree_3::Bool: Include degree-3 contribution (default: true).

Returns

• SVector{3}: Tidal acceleration [km/s²].

srp_accel(u::AbstractVector, sun_pos::AbstractVector, RC::Number; ShadowModel, R_Sun, R_Occulting, Ψ, AU)

Compute the acceleration from Solar Radiation Pressure.

Radiation from the Sun reflects off the satellite's surface and transfers momentum perturbing the satellite's trajectory. This force can be computed using a cannonball model with the following equation:

            𝐚 = F * RC * Ψ * (AU/(R_sc_Sun))^2 * R̂_sc_Sun

Arguments

• u::AbstractVector: The current state of the spacecraft in the central body's inertial frame [km, km/s].
• sun_pos::AbstractVector: The current position of the Sun [km].
• RC::Number: The reflectivity ballistic coefficient of the satellite – (Area/mass) * Reflectivity Coefficient [m^2/kg].

Keyword Arguments

• ShadowModel::ShadowModelType: Shadow model to use. Options: Conical(), Cylindrical(), SmoothedConical(), NoShadow(). Default: Conical().
• R_Sun::Number: The radius of the Sun [km]. Default: R_SUN.
• R_Occulting::Number: The radius of the occulting body [km]. Default: R_EARTH.
• Ψ::Number: Solar radiation pressure at 1 AU [N/m^2]. Default: SOLAR_FLUX.
• AU::Number: Astronomical Unit [km]. Default: ASTRONOMICAL_UNIT / 1E3.

Returns

• SVector{3}: Inertial acceleration from SRP [km/s^2].

sum_accelerations(u::AbstractVector, p::ComponentVector, t::Number, models::Tuple)

Internal recursive function to efficiently sum accelerations from multiple force models.

This function uses compile-time recursion over the tuple of force models to compute the total perturbation acceleration. The recursive approach allows the compiler to unroll the loop and inline all acceleration computations, providing optimal performance.

Arguments

• u::AbstractVector: Current spacecraft state vector
• p::ComponentVector: Simulation parameters
• t::Number: Current simulation time
• models::Tuple: Tuple of force models to sum over

Returns

• SVector{3}: Sum of all acceleration vectors from the force models

Implementation Notes

• Uses tail recursion for compile-time unrolling
• Base case returns zero acceleration for empty tuple
• Each recursive call processes one force model and continues with the tail

thermal_emission_accel(u, sun_pos, C_thm; ShadowModel, R_Sun, R_Occulting, Ψ, AU)

Compute the acceleration from spacecraft thermal emission (thermal re-radiation).

Spacecraft surfaces absorb solar radiation, heat up, and re-emit thermal photons. When emission is anisotropic (different emissivities on opposite sides), the net photon recoil produces a perturbative force. In the steady-state equilibrium approximation, the acceleration has the form [1]:

𝐚 = ν * C_thm * Ψ * (AU / d_sun)² * d̂_sun

This is structurally similar to solar radiation pressure but captures the distinct physical mechanism of absorbed-then-re-emitted thermal radiation rather than reflected/scattered photons. The thermal emission coefficient C_thm encodes the surface emissivity asymmetry and absorptivity properties of the spacecraft.

Arguments

• u::AbstractVector: Current state of the spacecraft in the central body inertial frame [km, km/s].
• sun_pos::AbstractVector: Current position of the Sun [km].
• C_thm::Number: Thermal emission coefficient [m²/kg].

Keyword Arguments

• ShadowModel::ShadowModelType: Shadow model. Default: Conical().
• R_Sun::Number: Radius of the Sun [km]. Default: R_SUN.
• R_Occulting::Number: Radius of the occulting body [km]. Default: R_EARTH.
• Ψ::Number: Solar radiation pressure at 1 AU [N/m²]. Default: SOLAR_FLUX.
• AU::Number: Astronomical Unit [km]. Default: ASTRONOMICAL_UNIT / 1E3.

Returns

• SVector{3}: Inertial acceleration from thermal emission [km/s²].

thermal_emission_coefficient(u, p, t, model::FixedThermalEmission)

Return the thermal emission coefficient for a fixed thermal emission model.

thermal_emission_coefficient(u, p, t, model::FlatPlateThermalModel)

Return the pre-computed thermal emission coefficient for a flat plate model.

third_body_accel(u::AbstractVector, μ_body::Number, body_pos::AbstractVector, h::Number) -> SVector{3}{Number}

Compute the Acceleration from a 3rd Body Represented as a Point Mass

Since the central body is also being acted upon by the third body, the acceleration of body 𝐁 acting on spacecraft 𝐀 in the orbiting body's 𝐂 is part of the force not acting on the central body

            a = ∇UB(rA) - ∇UB(rC)

Arguments

• u::AbstractVector: The current state of the spacecraft in the central body's inertial frame.
• μ_body: Gravitation Parameter of the 3rd body.
• body_pos::AbstractVector: The current position of the 3rd body in the central body's inertial frame.

Returns

• SVector{3}{Number}: Inertial acceleration from the 3rd body

thrust_acceleration(u::AbstractVector, p::AbstractVector, t::Number, model::ConstantCartesianThrust)

Returns the constant Cartesian thrust acceleration vector [km/s²].

Arguments

• u::AbstractVector: Current state of the simulation.
• p::AbstractVector: Parameters of the simulation.
• t::Number: Current time of the simulation.
• model::ConstantCartesianThrust: Constant Cartesian thrust model.

Returns

• SVector{3}: Thrust acceleration [km/s²] in the frame specified by the parent model.

thrust_acceleration(u::AbstractVector, p::AbstractVector, t::Number, model::StateThrustModel)

Returns the thrust acceleration vector computed by the user-provided function [km/s²].

Arguments

• u::AbstractVector: Current state of the simulation.
• p::AbstractVector: Parameters of the simulation.
• t::Number: Current time of the simulation.
• model::StateThrustModel: User-defined thrust model.

Returns

• SVector{3}: Thrust acceleration [km/s²] in the frame specified by the parent model.

thrust_acceleration(u::AbstractVector, p::AbstractVector, t::Number, model::ConstantTangentialThrust)

Returns the tangential thrust acceleration vector directed along the velocity [km/s²].

The result is always in the inertial frame regardless of the parent model's frame setting.

Arguments

• u::AbstractVector: Current state of the simulation.
• p::AbstractVector: Parameters of the simulation.
• t::Number: Current time of the simulation.
• model::ConstantTangentialThrust: Constant tangential thrust model.

Returns

• SVector{3}: Thrust acceleration in the inertial frame [km/s²].

thrust_acceleration(u::AbstractVector, p::AbstractVector, t::Number, model::PiecewiseConstantThrust)

Returns the piecewise-constant thrust acceleration for the arc containing time t [km/s²].

Uses a reverse linear scan to find the active arc (last arc whose start time ≤ t). If t is before the first arc, the first arc's acceleration is returned.

Arguments

• u::AbstractVector: Current state of the simulation.
• p::AbstractVector: Parameters of the simulation.
• t::Number: Current time of the simulation [s].
• model::PiecewiseConstantThrust: Piecewise-constant thrust schedule.

Returns

• SVector{3}: Thrust acceleration [km/s²] in the frame specified by the parent model.

transform_to_inertial(a_local::SVector{3}, u::AbstractVector, frame::InertialFrame)

Identity transformation — the acceleration is already in the inertial frame.

transform_to_inertial(a_local::SVector{3}, u::AbstractVector, frame::RTNFrame)

Transform a thrust acceleration from the RTN frame to the inertial frame.

The RTN basis is constructed from the spacecraft state u = [r; v] as:

R̂ = r / |r|
N̂ = (r × v) / |r × v|
T̂ = N̂ × R̂

The inertial acceleration is: a_R R̂ + a_T T̂ + a_N N̂.

transform_to_inertial(a_local::SVector{3}, u::AbstractVector, frame::VNBFrame)

Transform a thrust acceleration from the VNB frame to the inertial frame.

The VNB basis is constructed from the spacecraft state u = [r; v] as:

V̂ = v / |v|
B̂ = (r × v) / |r × v|
N̂ = B̂ × V̂

The inertial acceleration is: a_V V̂ + a_N N̂ + a_B B̂.