Impact of geomagnetic storm intensity on thermospheric density and orbital decay of LEO satellites: a G1–G5 per-event analysis

Technical Summary: Impact of Geomagnetic Storm Intensity on Thermospheric Density and Orbital Decay of LEO Satellites

Problem Statement
The rapid expansion of Low-Earth Orbit (LEO) satellite constellations has heightened the operational need for accurate thermospheric density characterization. Atmospheric drag is the dominant non-gravitational perturbation for satellites below 600 km, and geomagnetic storms significantly enhance thermospheric density, thereby accelerating orbital decay. While historical events, such as the loss of 38 Starlink satellites in February 2022 and the premature re-entry of a Starlink satellite in October 2024, demonstrate the severe operational risks, existing statistical methods often rely on superposed epoch analysis (SEA). SEA averages responses across many events, potentially smoothing out the sharp peak responses critical for immediate orbit prediction and collision avoidance. Furthermore, empirical thermospheric models often struggle to accurately reproduce density enhancements during intense storm events.

Methodology
This study employs a "per-event peak" analysis approach to isolate the instantaneous peak response of 22 geomagnetic storms ranging from G1 to G5 (maximum $K_p$ from 5 to 9) spanning solar cycles 23, 24, and 25 (2003–2024). The dataset includes the May 2024 Gannon G5 storm.

• Data Sources: In-situ neutral mass density and precise orbit data were derived from accelerometer measurements on the GRACE (2002–2017) and GRACE-FO (2018–present) twin satellites. Geomagnetic indices ($K_p$, $Dst$, $SYM-H$, $AU$, $AL$) and solar wind parameters were obtained from NASA/GSFC OMNIWeb.
• Event Selection: 22 events were selected based on three criteria: (1) maximum $K_p \ge 5$ (G1 or higher); (2) continuous accelerometer data availability from 48 hours pre-Storm Sudden Commencement (SSC) to 48 hours post-storm peak; and (3) no overlapping preceding storm main phases. The selection includes five events each for G1–G4 categories and the two available G5 events (October 2003 and May 2024).
• Processing:
  1. Orbit Selection: For each event, the "storm orbit" was defined as the orbit closest to the storm peak (minimum $Dst$), and the "quiet reference orbit" was the orbit closest to the SSC occurring after at least 24 consecutive hours of $K_p \le 3$.
  2. Decay Rate Calculation: Raw 10-second altitude data were averaged per orbit, then smoothed using a 24-hour moving average to remove semi-diurnal modulation. Linear fits were applied to the quiet and storm periods to determine the altitude decay rates ($R_{quiet}$ and $R_{storm}$).
  3. Enhancement Metric: The storm-time orbital decay enhancement was calculated as the percentage increase: $\frac{R_{storm} - R_{quiet}}{R_{quiet}} \times 100$.
  4. Statistical Modeling: An empirical multi-linear regression (MLR) model was developed to predict storm-time decay rates using three accessible predictors: quiet-time decay rate, satellite operation altitude, and maximum $K_p$.

Key Results

• Monotonic Scaling with Intensity: The mean orbital decay rate enhancement increases monotonically with storm intensity. The observed enhancements are: 77.2% (G1), 98.5% (G2), 138.1% (G3), 149.6% (G4), and 166.5% (G5).
• Altitude as a Critical Modulator: Satellite operation altitude acts as an independent and significant factor. Low-altitude events (e.g., ~350 km) exhibited strong drag responses even under relatively weak geomagnetic forcing, often exceeding the mean enhancement of higher-intensity storms at higher altitudes.
• Event-Specific Variability: Significant scatter exists within intensity categories. Notably, the two G5 events (October 2003 and May 2024) showed a ~2.6-fold difference in decay enhancement despite identical $K_p$ values and similar altitudes. The authors attribute the higher enhancement in the 2024 Gannon storm to a lower background solar activity level (lower quiet-time decay rate denominator) and a prolonged duration driven by successive CMEs (>30 hours), compared to the shorter, more intense 2003 Halloween storm.
• MLR Model Performance: The empirical model, $D_s = 1.55D_q - 0.29(h_0 - 450) + 7.39Kp_{max} - 35.60$, achieved a coefficient of determination ($R^2$) of 0.91 and a root mean square error (RMSE) of 12.3 m/day across all 22 events. The model successfully captures the dependence on storm intensity and altitude without requiring satellite-specific physical parameters (mass, cross-section, drag coefficient).

Significance and Claims
The paper positions this work as a complementary benchmark to existing SEA-based analyses. By focusing on per-event peak responses, it preserves the sharp signals that averaging procedures obscure, offering a more direct validation reference for empirical thermospheric models. The study provides a dataset and a simple, operationally accessible MLR tool that allows satellite operators to estimate storm-time decay rates using only standard tracking data (altitude, quiet-time decay) and space weather indices ($K_p$). The authors claim these results support the development of improved orbit-prediction and collision-avoidance strategies for the expanding LEO population, particularly in the context of the current Solar Cycle 25.

The authors maintain a modest scope regarding the model's generalizability, noting that the MLR was derived from two specific satellite platforms (GRACE/GRACE-FO) with similar physical characteristics. They acknowledge that future work is required to incorporate satellite mass, cross-sectional area, and velocity to extend the tool's applicability to the broader and more diverse LEO satellite population.