Force Model¶
The satkit numerical propagator integrates the forces acting on a satellite to produce a change in velocity, then integrates the velocity to produce a change in position:
This page describes each force in \(\vec{a}\left (t, \vec{p}, \vec{v}\right )\) that satkit models. The integration mechanics live in the ODE Integrators page; state vectors, the state transition matrix, and covariance propagation live in State Vectors, STM & Covariance.
The force model follows the treatment in Montenbruck & Gill (2000), Satellite Orbits: Models, Methods, Applications, Chapter 3 — consult that book for more depth on any term. Every model below cites its primary source; the full list is on the References page.
Summary¶
| Force | Default | Setting | Order of magnitude (LEO) |
|---|---|---|---|
| Earth gravity (spherical harmonics) | on | gravity_degree, gravity_order, gravity_model |
\(10^0\) m/s² |
| Sun third-body | on | use_sun_gravity |
\(10^{-6}\) m/s² |
| Moon third-body | on | use_moon_gravity |
\(10^{-6}\) m/s² |
| Atmospheric drag (NRLMSISE-00) | when alt < 1,000 km | use_spaceweather, satproperties.cd_a_over_m |
\(10^{-7}\) to \(10^{-3}\) m/s² |
| Solar radiation pressure (cannonball) | when craoverm > 0 |
satproperties.craoverm |
\(10^{-8}\) to \(10^{-7}\) m/s² |
| Solar radiation pressure (ECOM, experimental — see Empirical SRP: ECOM) | when ecom is set |
satproperties.ecom |
\(10^{-7}\) m/s² (D0), \(10^{-9}\) (Y, B) |
| Solid Earth tides (IERS 2010 §6.2.1 Step 1) | on | tide_model |
\(10^{-7}\) m/s² |
| General relativity (IERS 2010 §10.3: Schwarzschild + geodesic + Lense–Thirring) | on | use_relativistic_correction |
\(10^{-9}\) m/s² (LEO) |
| Continuous thrust | when configured | satproperties.thrusts |
user-specified |
The Forces-vs-altitude plot at the bottom of this page shows how each contribution scales with orbital altitude.
Earth Gravity¶
Earth's gravity dominates by many orders of magnitude. The Earth is not a point mass; its non-spherical mass distribution is captured by an expansion in spherical harmonics with coefficients \(\bar{C}_{nm}, \bar{S}_{nm}\):
The \(\bar{C}_{20}\) term (commonly called J2) captures Earth's equatorial bulge and is responsible for orbital precession. The expansion and the recursive evaluation of the normalized Legendre functions follow Montenbruck & Gill (2000), §3.2, Eqs. 3.28–3.33.
Coefficient files come from ICGEM (Ince et al. 2019). Five models are selectable via gravity_model:
| Model | Description | Tide system | Provided |
|---|---|---|---|
egm96 |
Earth Gravitational Model 1996 — Lemoine et al. (1998) | tide-free | compiled in |
egm2008 |
Earth Gravitational Model 2008 (default) — Pavlis et al. (2012) | tide-free | compiled in |
jgm3 |
Joint Gravity Model 3 — Tapley et al. (1996) | zero-tide | compiled in |
jgm2 |
Joint Gravity Model 2 — Nerem et al. (1994) | tide-free | compiled in |
itugrace16 |
ITU_GRACE16, GRACE-only — Akyilmaz et al. (2016); CC BY 4.0 | zero-tide | downloaded on first use |
The compiled-in models are truncated to degree 70, the evaluator's cap, so nothing is lost. At degree ≤ 40 the post-GRACE models agree with each other to ~10⁻¹¹ in the coefficients; the choice among them changes a one-day LEO propagation by a few metres and a GPS one by centimetres. Truncation matters more than the choice of model: at 400 km the terms from degree 41 to 70 are worth tens of metres per day, at 800 km a few metres, above 1200 km nothing.
The gravity_degree and gravity_order parameters cap the expansion (default 4×4, maximum 70×70). For high-precision work, degree-8 to degree-20 is typical; gains beyond ~degree-20 are small for satellites above ~500 km (satkit's recommendation — see the order-of-magnitude discussion of gravity-model truncation in Montenbruck & Gill 2000, §3.2). Order may be set lower than degree to zero out the longitudinal (tesseral) terms.
Third-Body Gravity (Sun, Moon)¶
The Sun and Moon each act as point-mass attractors. Their pull on the Earth must be subtracted to express acceleration in the geocentric frame:
The Moon expression has the same form (Montenbruck & Gill 2000, §3.3.1, Eq. 3.37). Body positions come from the JPL DE-series ephemerides (default DE440; Park et al. 2021). Disable either with use_sun_gravity=False / use_moon_gravity=False.
Solid Earth Tides¶
The Sun and Moon deform the solid body of the Earth, which in turn perturbs the geopotential. satkit implements IERS 2010 §6.2.1 Step 1 (Petit & Luzum 2010, Eqs. 6.6–6.7, with the anelastic Love numbers of Table 6.3) — the frequency-independent Love-number response, which in satkit's own measurements carries about 99% of the total solid-tide signal at ~5% per-step CPU overhead. The correction modifies the degree-2, degree-3, and a small degree-4 set of Stokes coefficients \(\Delta\bar{C}_{nm}, \Delta\bar{S}_{nm}\) as a function of Sun and Moon ITRF positions and the IERS 2010 nominal Love numbers.
Frequency-dependent corrections (IERS 2010 §6.2.2 Step 2, Tables 6.5a–c, 71 tidal constituents) are reserved for TideModel::SolidFull — currently a placeholder that falls back to Step 1. The Step 2 contribution is sub-mm class at LEO.
Set with tide_model:
| Value | Effect |
|---|---|
tidemodel.solid_step1 (default) |
IERS §6.2.1 Step 1 — recommended |
tidemodel.solid_full |
Step 1 + Step 2 (Step 2 not yet implemented; behaves as Step 1) |
tidemodel.none |
Disable solid tides (use for reproducibility with pre-tide versions) |
Tide system of the gravity model
Step 1 includes the permanent tide (its time average, \(A_0 H_0 k_{20} \approx -4.2\times10^{-9}\) in \(\bar{C}_{20}\); Petit & Luzum 2010, §6.2.2). A tide-free model (egm96, egm2008, jgm2; \(\bar{C}_{20} \approx -4.841654\times10^{-4}\)) has that part removed from its coefficients, so Step 1 completes it. A zero-tide model (jgm3, itugrace16; \(\bar{C}_{20} \approx -4.841695\times10^{-4}\)) already carries it. satkit reads each model's tide system on load (from the ICGEM tide_system header, or from the \(\bar{C}_{20}\) value for the older files that lack one) and, for a zero-tide model, removes the permanent part from the Step 1 correction (IERS 2010 Eq. 6.13) instead of counting it twice. The published coefficients are never altered, so with tidemodel.none every model propagates exactly as published — and a zero-tide model with tides off still differs from a tide-free one by the permanent tide, about 8 m/day at 500 km with J2 alone. Any model can therefore be paired with any tide setting; before 0.23, jgm3 or itugrace16 with solid_step1 double-counted the permanent tide (3–12 m/day at LEO, 0.3 m/day at GPS).
Atmospheric Drag¶
At altitudes below 1,000 km the residual atmosphere imposes a drag force:
with \(C_d\) the drag coefficient (typically 1.5–3; Montenbruck & Gill 2000, §3.5), \(A/m\) the area-to-mass ratio, \(\rho\) the atmospheric density, and \(\vec{v}_r\) the satellite velocity relative to the co-rotating atmosphere (assumed at rest in the Earth-fixed frame).
Density comes from the NRLMSISE-00 thermosphere model (Picone et al. 2002; a pure-Rust port of Dominik Brodowski's C implementation), which reads the space-weather indices from satkit's table (the GFZ Potsdam observed record, then the NOAA/SWPC 45-day and NASA MSAFE monthly forecasts — see Data coverage) automatically, following the model's own interface: the observed F10.7 of the previous UTC day, the observed 81-day centred average F10.7A of the current day, and the 7-element 3-hourly geomagnetic history (the current day's daily Ap, the current 3-hourly ap and the three before it, and the 12–33 h and 36–57 h means — NRLMSISE-00 switch 9 = −1), which lets the density respond to a geomagnetic storm within hours rather than at the next UTC midnight. When that history cannot be assembled (a day missing from the table, or a monthly forecast row, whose single climatological Ap has no 3-hourly structure) the model falls back to the current day's daily Ap. Each index stands on its own: a day missing its F10.7 in the observed record takes the most recent measured flux up to three days back (else the 81-day average) and keeps its Ap; only an index the table cannot supply at all — before the table starts, or F10.7 before 1947 — takes the nominal value below, with a one-time warning. Disable space-weather lookup with use_spaceweather=False to fall back to fixed nominal indices (F10.7 = F10.7A = 150, Ap = 4).
Drag is skipped more than 1,000 km above the equatorial radius — the upper limit of NRLMSISE-00's validity — regardless of settings. (Before 0.24 the cutoff was 700 km.)
Solar Radiation Pressure¶
Solar photons absorbed or scattered by the satellite transfer momentum, producing a force away from the Sun:
where \(P_\text{sun} = 1367\,\text{W/m}^2 / c \approx 4.56 \times 10^{-6}\) N/m² is the radiation pressure at 1 AU (consts.solar_pressure_1au), \(d\) is the satellite–Sun distance (\(\hat{p}_\text{sun}\) is the satellite→Sun direction), \(C_R A/m\) is the satellite's radiation susceptibility (the user-supplied satproperties.craoverm), and \(\nu(\vec{p}, \vec{p}_\text{sun}) \in [0, 1]\) is a shadow function that vanishes when the satellite is in Earth's umbra — the conical umbra/penumbra model of Montenbruck & Gill (2000), §3.4.2; the cannonball force itself is their §3.4, Eq. 3.75.
The \((\text{AU}/d)^2\) factor makes the flux follow Earth's eccentric orbit: it is 3.4 % above the 1 AU value at perihelion (early January) and 3.3 % below at aphelion (early July). satkit versions before 0.23.1 held it at the 1 AU value year-round.
Solar constant: satkit uses 1367 W/m², not the measured 1361 W/m²
| Solar constant | \(P_\text{sun}\) at 1 AU | Source |
|---|---|---|
| 1367 W/m² (satkit) | \(4.560 \times 10^{-6}\) N/m² | 1980s–90s radiometry (Nimbus-7, ACRIM, VIRGO); WMO standard; GMAT and STK default |
| 1361 W/m² (measured) | \(4.540 \times 10^{-6}\) N/m² | SORCE/TIM (2003–); Kopp & Lean (2011); IAU 2015 Resolution B3 nominal value |
The Sun did not dim: the older radiometers read high because light scattered past their precision aperture into the detector cavity, which ground tests at the NIST/LASP TSI Radiometer Facility later confirmed. satkit keeps 1367 W/m² so SRP results compare directly with GMAT and STK. The 0.4 % difference is well inside typical \(C_R\) uncertainty and is absorbed whenever craoverm is estimated; to model 1361 W/m² exactly, multiply craoverm by 1361/1367.
satkit's default is a cannonball model — the satellite's surface is treated as if its normal points toward the Sun, and the force acts along the satellite→Sun line. A physical box-wing model (Rodriguez-Solano et al. 2012) is not provided.
For GNSS-class work satkit also offers the Empirical CODE Orbit Model (ECOM) — the empirical, Sun-oriented D/Y/B parameterization used by CODE and most IGS analysis centres — as an experimental addition to the cannonball term, enabled by supplying coefficients through satproperties(ecom=...). Its equations, coefficient table, sign and eclipse conventions, and measured performance are on the Empirical SRP: ECOM page.
General-Relativistic Correction¶
The full IERS 2010 §10.3 Eq. 10.12 correction (Petit & Luzum 2010) with PPN parameters \(\beta = \gamma = 1\), three terms:
where \(\vec{\Omega} = \tfrac{3}{2}\,\frac{GM_\odot}{c^2 R_\oplus^3}\,(\vec{R}_\oplus\times\dot{\vec{R}}_\oplus)\) is the geodesic (de Sitter) precession of the geocentric frame from the Earth's heliocentric state \(\vec{R}_\oplus, \dot{\vec{R}}_\oplus\) (\(|\Omega| \approx 1.92''\)/century), and \(\vec{J} = \tfrac{2}{5}R^2\omega\,\hat{z}_\text{ITRF}\) is the Earth's spin angular momentum per unit mass (homogeneous rigid-sphere approximation, \(\approx 1.19\times10^9\) m²/s; the real moment-of-inertia factor is ≈0.33 rather than 0.4, a ~20% overstatement of a term that is ≤\(10^{-10}\) m/s² and mostly periodic — GMAT uses the same approximation, GMAT Mathematical Specifications, §4.2.6).
The Schwarzschild term dominates below GEO (~\(10^{-9}\) m/s² at LEO; in satkit runs, omitting it costs ~1 m/day at GPS altitude and ~3 m/day at GEO). The geodesic term is a near-constant \(4\times10^{-11}\) m/s² Coriolis-like acceleration at LEO and becomes the largest relativistic term beyond ~100,000 km (~1 m over 7 days at 200,000 km, measured in the GMAT comparison). Lense–Thirring is ~\(10^{-10}\) m/s² at LEO, mostly periodic, and falls as \(1/r^3\). This is the same formulation as GMAT's RelativisticCorrection (GMAT Mathematical Specifications, §4.2.6).
Toggle with use_relativistic_correction (default True).
Continuous Thrust¶
Constant acceleration in a chosen frame (used to model low-thrust maneuvers or — in orbit-determination contexts — empirical accelerations: a fitted catch-all that absorbs un-modeled physics). See satkit.thrust and the GPS Example tutorial for usage.
Future Propagation¶
When propagating into the future beyond the date range of downloaded data files:
- Earth Orientation Parameters (\(\Delta UT1\), \(x_p\), \(y_p\)): the last available values are held constant, with a one-time warning. This is much more accurate than defaulting to zero, but still drifts by ~0.1″ / ~10 ms over a few months (metres at LEO); check
satkit.frametransform.eop_status(t)or setpropsettings.require_eop_coverage = Trueto make the propagator raise instead (it also raises for a span before the table's 1973-01-02 start, where zero EOP would otherwise be used), and refresh withsatkit.utils.update_datafiles(). A propagation with no EOP table loaded at all is refused. - Space Weather (F10.7 solar flux, Ap geomagnetic index): past the observed record the table continues with the NOAA/SWPC 45-day forecast and then NASA's MSAFE monthly forecast, which carries a climatological Ap — see Data coverage.
Forces vs Altitude¶
The plot below, modeled on Fig. 3.1 of Montenbruck & Gill (2000), shows each force's order of magnitude vs orbital altitude:
See Also¶
- Tutorial: GPS Example — fits a GPS orbit against ESA SP3 truth and walks through the empirical-acceleration concept.
- Theory: ODE Integrators for the integration mechanics; State Vectors, STM & Covariance for state representation and covariance propagation.
- Validation: GMAT Comparison — how this force model agrees with NASA GMAT across seven orbital regimes, and which terms account for the remaining differences.
- API:
satkit.propagate,satkit.propsettings,satkit.satproperties. - References: the sources for every model on this page are collected on the References page.