Icarus Model: Asteroid & Heliospheric MHD

• Icarus Model is a comprehensive framework that quantitatively characterizes NEA 1566 Icarus using multi-decadal radar and astrometric data, resolving its shape, spin, and Yarkovsky drift with high precision.
• The heliospheric MHD component simulates solar wind and CME propagation using adaptive mesh refinement and dynamic inner boundary conditions to capture fine-scale plasma dynamics.
• Operational advances in Icarus yield significant computational speed-ups and improved shock capture, enhancing real-time space-weather forecasting and particle acceleration modeling.

The term “Icarus Model” refers to rigorously quantitative and modular models, each significant in its respective domain—asteroid physical/rotational/orbital characterization, the 3D heliospheric magnetohydrodynamics (MHD) of solar wind and transient solar events, and operational space-weather forecasting. The following exposition pursues the full technical depth of the Icarus Model(s) as presently defined in prominent literature, focusing on the dynamical asteroid model (Asteroid 1566 Icarus) and the suite of heliospheric MHD models that share the same name, while referencing the distinct technical and modeling principles underlying each.

1. Icarus Model for Asteroid 1566 Icarus: Physical, Rotational, and Dynamical Description

The asteroid-focused Icarus Model synthesizes multidecadal radar and optical astrometry for the near-Earth asteroid 1566 Icarus ($a=1.08$ AU, $e=0.83$, $i=22.8^\circ$). This model resolves the asteroid’s debated shape, spin, surface scattering, and Yarkovsky-driven orbital drift with unprecedented precision (Greenberg et al., 2016).

1.1 Shape and Size Inversion

• The shape inversion used the “shape” software, fitting delay-Doppler radar images (0.1–0.2 μs range, 0.3–1 Hz Doppler) and CW spectra from the 2015 Arecibo/Goldstone campaign.
• The triaxial ellipsoid delivered the best fit:
  • $2a ≃ 1.61$ km, $2b ≃ 1.60$ km, $2c ≃ 1.17$ km
  • Equivalent volume diameter: $D_\text{eq} = 1.44$ km (±18%).

1.2 Spin–Pole Determination

• A 15°-spaced grid search in pole orientation uniquely isolated a retrograde solution: ecliptic longitude $\lambda=270^\circ\pm10^\circ$, latitude $\beta=-81^\circ\pm10^\circ$.
• This determination supersedes all prior lightcurve-based poles. It is required to reproduce both the daily bandwidth variations in radar data and the historical CW bandwidths across decades without invoking non-principal-axis rotation.

1.3 Surface Scattering and Radar Albedo

• The echo strength versus incidence $\theta$ fits a Hagfors surface:

$$\frac{d\sigma}{dA}(\theta) = H(\theta_0 - |\theta|)\ \frac{RC}{2} \left(\cos^4\theta + C\sin^2\theta\right)^{-3/2}$$

with $C \approx 13.8 \pm 1.4$, corresponding to decimeter-scale RMS slopes of $S_0 \sim 16^\circ$.

• S-band radar albedo is $0.02$—among the lowest observed—implying unusually high surface porosity ($\epsilon \approx 1.8$).

1.4 Integrated Orbit Fitting and Yarkovsky Drift

• The orbit-fit pipeline (IDOS, using JPL MONTE/DIVA and a Kalman filter), with all modern and historical astrometry, produces an eight-sigma measurement of Yarkovsky drift:

$$\dot{a} = (-4.62 \pm 0.48) \times 10^{-4}\ \mathrm{au/Myr}$$

• Model for transverse acceleration due to anisotropic thermal re-radiation:

$$\ddot{\mathbf{r}}(t) = \xi\,\frac{3}{8\pi}\,\frac{L_\odot}{c\,D\,\rho}\, \frac{X_{\hat{p}}(\phi)\,\mathbf{r}(t)}{|\mathbf{r}(t)|^3}$$

• Orbit-averaged semi-major-axis drift:

$$\frac{da}{dt} = \xi\,\frac{3 L_\odot}{4 \pi c \sqrt{G M_\odot}}\, \frac{\hat{\alpha}(e)}{D \rho \sqrt{a}}$$

where $\hat\alpha(e)\sim3.16$ for Icarus.

• The model reconciles previously discrepant Yarkovsky rates and confirms Icarus is a “clean retrograde Yarkovsky-drifter”—i.e., its sense of rotation and measured drift are in dynamical agreement.

Table 1. Summary: Physical Parameters

Parameter	Value	Uncertainty
$D_\text{eq}$	1.44 km	±18%
$a$	1.08 AU	—
$e$	0.83	—
Spin pole ($\lambda$)	$270^\circ$ (retrograde)	±10°
Spin pole ($\beta$)	$-81^\circ$	±10°
Radar albedo	2%	—
Yarkovsky rate	$-4.62 \times 10^{-4}\,\mathrm{au/Myr}$	±0.48e-4 AU/Myr

2. The Icarus Heliospheric MHD Model: Numerical Architecture and Governing Equations

The heliospheric Icarus Model signifies a 3D, block-structured, adaptive mesh refinement (AMR) MHD solver for solar wind and CME/CIR propagation, implemented atop MPI-AMRVAC (Baratashvili et al., 2022, Baratashvili et al., 4 Jan 2024, Baratashvili et al., 25 Jul 2024, Baratashvili et al., 27 Jan 2025, Husidic et al., 18 Apr 2024, Baratashvili et al., 28 May 2024).

2.1 Mathematical Formulation

• Core system: ideal MHD equations in a solar co-rotating spherical domain, with explicit Coriolis and centrifugal source terms.
• Standard form:

$$\frac{\partial \rho}{\partial t} + \nabla \cdot (\rho \mathbf{v}) = 0$$

$$\frac{\partial (\rho \mathbf{v})}{\partial t} + \nabla \cdot [\rho \mathbf{v} \mathbf{v} + (p + \frac{B^2}{2}) \mathbf{I} - \mathbf{B}\mathbf{B}] = \rho \mathbf{g} - 2 \rho \Omega \times \mathbf{v} - \rho \Omega \times (\Omega \times \mathbf{r})$$

$$\frac{\partial \mathbf{B}}{\partial t} + \nabla \times (\mathbf{v} \times \mathbf{B}) = 0$$

$$\frac{\partial E}{\partial t} + \nabla \cdot [(E + p + \frac{B^2}{2})\mathbf{v} - (\mathbf{v} \cdot \mathbf{B})\mathbf{B}] = \rho \mathbf{g} \cdot \mathbf{v} - 2\rho \mathbf{v} \cdot (\Omega \times \mathbf{v})$$

• $E$: total (thermal + kinetic + magnetic) energy
• Divergence control: parabolic (Dedner) or hyperbolic “eight-wave” cleaning
• Adiabatic index: $\gamma = 5/3$ (standard), with physically motivated variants for coronal and wind regions.

2.2 Domain, Boundary and Initial Conditions

• Domain: $r \in [0.1, 2]$ AU, stretched radially to maintain aspect ratio; latitude [-60°, +60°]; longitude covers $360^\circ$.
• Inner boundary ($r=0.1$ AU): time-dependent, 2D maps of $v_r, n, T, B_r$ from WSA or full-MHD coronal models (COCONUT).
• Outer boundary ($r=2$ AU or beyond): zero-gradient outflow (Neumann).

3. Numerical Techniques: Grid Stretching, Adaptive Mesh Refinement, and Solver Selection

3.1 Grid Stretching

• Radial grid stretching achieves cubic aspect ratios out to 2 AU; mapping: $$r(\xi) = r_\text{in} \left( \tfrac{r_\text{out}}{r_\text{in}} \right)^{\xi}$$
• Base grid typically $60 \times 32 \times 96$ (r, θ, ϕ) with up to 5 AMR levels.

3.2 Adaptive Mesh Refinement (AMR)

• Shock detection: refine on $\nabla \cdot \mathbf{v} < 0$ (compression).
• CME/structure tracing: refine on passive tracers or strong density gradients.
• Combined criterion: logical OR of above, spatially limited to regions of interest (e.g., ±30° Sun–Earth cone).
• Each refinement level doubles resolution; AMR Level 4 yields $\Delta r \simeq 1.34\,R_\odot$ at 1 AU.
• AMR reduces wall-clock time by factors of 2.7–31 and memory footprint by 50–75% relative to uniform grids, with arrival-time accuracy ($\Delta t \lesssim 0.1$ hr) approaching full-resolution limits (Baratashvili et al., 4 Jan 2024, Baratashvili et al., 2022).

3.3 Numerical Solver

• Riemann solver: TVD Lax–Friedrichs (TVDLF), sometimes HLL; option for higher-order variants in development.
• Time-stepping: explicit strong-stability-preserving Runge–Kutta; CFL given by the fast magnetosonic speed.
• Slope limiter: Woodward, van Leer, or Koren.

4. Model Implementations for CME Initiation, Propagation, and Validation

4.1 CME Models

• Cone CME: hydrodynamic, unmagnetized, geometry after Scolini et al. (spherical cap, user-tuned half-width, latitude/longitude, launch speed).
• Linear force-free spheromak: $\nabla \times \mathbf{B} = \alpha \mathbf{B}$, Chandrasekhar’s solution, injected at 0.1 AU using the EUHFORIA spheromak workflow (Baratashvili et al., 4 Jan 2024, Baratashvili et al., 28 May 2024).
• Parameter tuning: performed against in situ spacecraft (e.g., ACE, MESSENGER, OMNI) using observed CME arrival, field strength, and shock characteristics.

4.2 Validation and Performance Metrics

• CME shock arrival, peak field, and density jumps at L1 are robustly reproduced—best with AMR Level 4 (error in arrival time $\Delta t < 0.5$ hr, field match within 2–10%).
• For multi-spacecraft studies (MESSENGER/Earth), AMR Level 3 suffices near Mercury ($\Delta r\sim 1\,R_\odot$), AMR Level 4 at Earth ($1.34\,R_\odot$).
• Wall-clock times: AMR Level 4 ($\sim$2.5 h on 36 cores) versus uniform grid ($\sim$7 h Icarus, $>$18 h EUHFORIA) for Sun–2 AU.
• Limitation: In multiple-CME interaction scenarios, CME deceleration is underestimated; a drag force or improved expansion model is recommended.

Table 2. Wall-Clock Time vs. Grid/AMR Level (36 cores, July 12 2012 Event)

Simulation	Wall time	Speed-up vs. uniform
EUHFORIA uniform	18 h 34 m	1×
Icarus uniform	6 h 52 m	2.7×
Icarus AMR 2	11 m	101×
Icarus AMR 4	2 h 34 m	7.2×
Icarus AMR 5	4 h 26 m	4.2×

5. Dynamic Inner Boundary and Sun-to-Earth Operational Forecasting

5.1 Time-Dependent Boundary Conditions (Icarus 3.0)

• Inner boundary at 0.1 AU is driven by high-cadence ($\geq1$ h) magnetogram-based coronal model output (WSA or COCONUT).
• These are packed into 3D VTK tables (θ, φ, t); at each MHD timestep, boundary states are interpolated in time (“bc_data” module).
• Dynamic-driving offers improved reproduction of timing and fine structure in wind/CME vs. static boundary (Baratashvili et al., 27 Jan 2025, Baratashvili et al., 25 Jul 2024).

5.2 Coronal–Heliospheric Coupling

• Automated pipeline fetches magnetograms, computes coronal model, stacks boundary conditions, and launches Icarus with dynamic boundary updates.
• Relaxation (“spin-up”) window (typically 8 days) followed by forecast interval.

5.3 Validation and Operational Implications

• Dynamic boundary more accurately represents solar cycle and stream transients, yielding improvement in density, velocity, and shock signatures at observer locations (e.g., OMNI, STEREO, MESSENGER).
• Wall-times for dynamic runs (∼1–2 h on 100+ cores) support real-time/ensemble operational use.

6. Limitations and Current Development Frontiers

• CME scheme: currently, spheromak and cone implementations; planned extensions include spheromak/flux-rope with evolving self-similar expansion and on-the-fly inner-heliospheric coupling.
• Physical realism is ultimately restricted by the accuracy of the upstream coronal boundary model (COCONUT, WSA); strong sensitivity is seen in verification studies.
• AMR criteria and grid tuning are active development areas, with the goal to further localize refinement and reduce computational cost without sacrificing forecast accuracy.

7. Integration with Particle Acceleration Modeling

• Icarus+PARADISE allows hybrid MHD–stochastic test-particle kinetic studies of solar energetic particles (SEPs) in realistic CME/CIR-driven wind (Husidic et al., 18 Apr 2024).
• AMR in Icarus leads to improved shock capture and higher particle acceleration efficiency in coupled kinetic runs, as sharper shocks concentrate acceleration zones, in contrast to more diffusive uniform grids.
• The data interface allows trilinear interpolation of wind/B fields to moving ensemble particles, with the AMR arrangement ensuring adequate coverage at CME/CIR fronts.

---

In summary, the Icarus Model denotes a self-consistent, high-performance, and well-calibrated framework for both small-body dynamical/physical characterization (specifically for NEA 1566 Icarus) and 3D/AMR-enabled heliospheric MHD modeling of the inner solar system plasma environment. The innovations in shape inversion, radar-driven material property inference, and precision Yarkovsky drift for asteroids directly parallel the model’s refinement and computational advances for global, real-time, and operational MHD space-weather forecasting. Each Icarus framework provides benchmarked, reproducible simulations and parameter inferences for subsequent physical paper or mission planning.