PARADISE Particle Transport Code

• PARADISE is a stochastic focused-transport solver that propagates energetic particles through prescribed heliospheric and coronal MHD fields.
• It transforms solar wind, magnetic field, and shock data from models like EUHFORIA, Icarus, and COCONUT into time-dependent particle intensity profiles.
• The code employs an Itô stochastic differential equation framework, efficiently simulating particle movement without coupling back to the plasma dynamics.

PARADISE, short for PArticle Radiation Asset Directed at Interplanetary Space Exploration, is an energetic-particle transport code for heliospheric and, more recently, coronal applications. Across its published implementations, it is described as a stochastic focused-transport solver that propagates energetic particles as test particles through prescribed magnetohydrodynamic backgrounds, rather than as a self-consistent code in which particles modify the plasma. In the heliospheric studies, PARADISE takes solar-wind velocity, magnetic field, shock geometry, and related quantities from models such as EUHFORIA and Icarus and converts them into time-dependent particle intensities that can be compared directly with in-situ observations; in the coronal extension, it performs the same role on a time-dependent COCONUT background (Wijsen et al., 2023, Husidic et al., 2024, Husidic et al., 2024, Wijsen et al., 2022).

1. Concept and computational role

PARADISE is used as the energetic-particle transport component of a coupled modelling chain. Its principal function is to evolve energetic-particle distributions in structured, time-dependent plasma and magnetic-field environments supplied by external MHD models. In the 2021 October 9 solar energetic particle event study, it is described as the module that takes the heliospheric plasma and magnetic-field environment from EUHFORIA and turns it into predicted time profiles of SEP intensities at spacecraft, specifically producing the directional differential intensity on a 3D spherical mesh by solving the 5-dimensional focused transport equation through an equivalent set of Itô stochastic differential equations integrated forward in time (Wijsen et al., 2023).

The same underlying role appears in later framework papers. In Icarus+PARADISE, PARADISE propagates energetic particles as test particles through prescribed solar-wind configurations from 0.1 au onward and computes time-dependent intensity distributions; in COCONUT+PARADISE, it takes evolving coronal MHD fields and advances particles through them on a coronal grid (Husidic et al., 2024, Husidic et al., 2024). The 14 July 2012 energetic storm particle study likewise presents PARADISE as a stochastic focused-transport code coupled to EUHFORIA for a real-event simulation near Earth (Wijsen et al., 2022).

A recurring feature of these applications is that PARADISE is not presented as a first-principles microscopic acceleration model. Its main role is to transport particle populations injected at shocks or placed in coronal structures while accounting for field-aligned motion, pitch-angle scattering, cross-field diffusion, and adiabatic energy change. This suggests that PARADISE occupies the intermediate layer between global MHD structure and particle observations: it uses resolved heliospheric or coronal geometry to determine particle access, trapping, and intensity-time profiles, while leaving the background plasma dynamics to the MHD model and, in several studies, prescribing shock injection phenomenologically (Wijsen et al., 2023, Wijsen et al., 2022).

2. Focused-transport formulation

At its core, PARADISE solves the gyrotropic focused transport equation for the distribution function $f(\mathbf{x},p,\mu,t)$, where $\mathbf{x}$ is position, $p$ is momentum magnitude, $\mu=\cos\alpha$ is pitch-angle cosine, and $t$ is time (Husidic et al., 2024). In the form given for both the Icarus and COCONUT couplings, the equation is

$$\frac{\partial f}{\partial t} + \frac{d\mathbf{x}}{dt}\cdot\nabla f + \frac{d\mu}{dt}\frac{\partial f}{\partial \mu} + \frac{dp}{dt}\frac{\partial f}{\partial p} = \frac{\partial}{\partial \mu}\left(D_{\mu\mu}\frac{\partial f}{\partial \mu}\right) + \nabla\cdot\left(\boldsymbol{\kappa}_\perp\cdot\nabla f\right).$$

The phase-space characteristics used in PARADISE are written as

$$\frac{d\mathbf{x}}{dt} = \mathbf{V}_{\mathrm{sw}} + \mathbf{V}_{\mathrm{d}} + \mu v \mathbf{b},$$

$$\frac{d \mu}{dt} = \frac{1-\mu^2}{2} \left( v\,\nabla\cdot\mathbf{b} +\mu\,\nabla\cdot\mathbf{V}_{\mathrm{sw}} -3\mu\,\mathbf{b}\mathbf{b}:\nabla\mathbf{V}_{\mathrm{sw}} -\frac{2}{v}\mathbf{b}\cdot\frac{d\mathbf{V}_{\mathrm{sw}}}{dt} \right),$$

$$\frac{dp}{dt} = \left[ \frac{1-3\mu^2}{2}\left(\mathbf{b}\mathbf{b}:\nabla\mathbf{V}_{\mathrm{sw}}\right) -\frac{1-\mu^2}{2}\nabla\cdot\mathbf{V}_{\mathrm{sw}} -\frac{\mu}{v}\mathbf{b}\cdot\frac{d\mathbf{V}_{\mathrm{sw}}}{dt} \right]p.$$

These terms encode convection with the solar wind, magnetic focusing, pitch-angle evolution, adiabatic deceleration or acceleration, and perpendicular transport (Husidic et al., 2024, Husidic et al., 2024). The 2022 event study explicitly notes that the FTE contains the physics needed for acceleration in converging, accelerating, and shear flows, and that the momentum term also captures energy changes due to particle drifts along the motional electric field; however, that study does not include a momentum-space diffusion term $D_{pp}$, so stochastic second-order acceleration downstream of the shock is not modelled (Wijsen et al., 2022).

PARADISE solves this equation stochastically by converting the deterministic-plus-diffusive problem into an equivalent set of stochastic differential equations and integrating them with Itô calculus. The resulting Monte Carlo or stochastic test-particle approach is central to the code’s design: it avoids solving the full phase-space transport equation on a fixed grid while retaining realistic 3D background structure (Wijsen et al., 2023, Husidic et al., 2024). In the Icarus validation paper, the equations are written in mixed coordinates, with spatial variables in the inertial frame and momentum and pitch angle in the solar-wind comoving frame; terms of order $\mathbf{x}$0 are neglected because $\mathbf{x}$1, and the authors explicitly set $\mathbf{x}$2, so magnetic drifts are not included in that particular study (Husidic et al., 2024).

3. Coupling to MHD backgrounds

PARADISE has been coupled to several global plasma models, each of which supplies the background through which particles are transported.

Framework	Background model	Coupling characteristics
EUHFORIA+PARADISE	3D heliospheric ideal-MHD solar wind and CME model	Uses EUHFORIA background fields and shock tracer; inner boundary at 0.1 au
Icarus+PARADISE	3D MHD code based on MPI-AMRVAC	Supports AMR and grid stretching; uses trilinear interpolation and on-demand block queries
COCONUT+PARADISE	Global coronal MHD model on an unstructured prism grid	Time interpolation between snapshots of cadence 289 s; adapted particle tracing on unstructured mesh

In the EUHFORIA-based studies, the workflow is explicit. EUHFORIA is run first to model the background solar wind, structures such as an HSS or SIR, and the evolving CME-driven shock. PARADISE is then run on top of this background, taking the solar-wind velocity field $\mathbf{x}$3, the magnetic field, and the evolving shock surface and properties from EUHFORIA output. The shock is identified with EUHFORIA’s shock tracer, and PARADISE injects particles continuously along the forward shock, but only from the CME-driven shock beyond the inner boundary at 0.1 au (Wijsen et al., 2023). In the 14 July 2012 ESP study, the EUHFORIA shock-finding and shock-tracing tool additionally determines the magnetically connected shock point, shock angle $\mathbf{x}$4, shock speed, and compression ratio, which are then supplied to PARADISE (Wijsen et al., 2022).

The Icarus coupling preserves the same conceptual division but changes the numerics. Icarus supplies density, velocity, and magnetic field on a 3D grid, and PARADISE must locate each test particle inside that grid, read the local MHD variables from surrounding cells, and interpolate them at each time step. The paper states that trilinear interpolation is used, as in the earlier EUHFORIA+PARADISE coupling. Because Icarus supports adaptive mesh refinement and grid stretching, the coupling must also handle nested block boundaries; the authors therefore query neighboring blocks on demand rather than storing ghost cells for all blocks, which they judge infeasible because it would create very large data files (Husidic et al., 2024).

The COCONUT extension shifts PARADISE into the low corona. COCONUT provides a sequence of time-dependent MHD snapshots with a cadence of 289 s, and PARADISE linearly interpolates the background fields in time between successive snapshots. Because COCONUT uses an unstructured prism grid, PARADISE was adapted to trace particles and compute gradients on that mesh. The authors also developed an interpolation to a structured EUHFORIA-like grid as a faster alternative, although the main results use the original unstructured output. For production runs, the minimum PARADISE timestep is $\mathbf{x}$5, and the stochastic scheme is parallelized by distributing pseudo-particles over cores (Husidic et al., 2024).

4. Scattering, diffusion, and source prescriptions

A consistent theme in PARADISE applications is the modular treatment of scattering and diffusion. The Icarus paper states explicitly that diffusion prescriptions can be added, removed, or replaced, which is useful because different transport problems require different turbulence models (Husidic et al., 2024). In heliospheric applications, pitch-angle diffusion is commonly represented with a quasi-linear-theory-based expression. One form used in the Icarus framework is

$\mathbf{x}$6

with $\mathbf{x}$7 and a constant parallel mean free path $\mathbf{x}$8 au for the 1 MeV protons in that study. The same paper uses a constant perpendicular mean free path $\mathbf{x}$9 au, giving $p$0, so transport is strongly field-aligned with weak cross-field spreading (Husidic et al., 2024). The 2021 October 9 SEP-event paper adopts the same $p$1 and a rigidity-dependent heliospheric parallel mean free path,

$p$2

with $p$3 and $p$4, the Kolmogorov turbulence spectral index (Wijsen et al., 2023).

The coronal COCONUT study uses a related but more general pitch-angle diffusion coefficient,

$p$5

with $p$6 and $p$7. It defines

$p$8

and implements perpendicular transport through

$p$9

That paper compares a constant perpendicular mean free path with a Larmor radius-dependent prescription and caps $\mu=\cos\alpha$0 at $\mu=\cos\alpha$1 in very weak magnetic fields (Husidic et al., 2024).

Source prescriptions in PARADISE are generally phenomenological. In the 14 July 2012 ESP study, continuous injection of 50 keV protons at the CME-driven shock is scaled with local solar-wind compression according to

$\mu=\cos\alpha$2

The authors state explicitly that this does not solve the injection problem from first principles, but is a phenomenological prescription motivated by the idea that strong compression regions are more likely to contain preheated seed particles and to facilitate injection into diffusive shock acceleration (Wijsen et al., 2022). The 2023 SEP-event study uses a closely related compression-based source in which the injected differential intensity scales with $\mu=\cos\alpha$3, an energy factor $\mu=\cos\alpha$4, and a radial factor $\mu=\cos\alpha$5, with normalization fixed by the ACE 68–115 keV EPAM/LEMS120 intensity at shock arrival (Wijsen et al., 2023). The Icarus validation study uses a different idealized prescription,

$\mu=\cos\alpha$6

implemented by distributing 1 MeV protons uniformly across shock regions identified by $\mu=\cos\alpha$7 and assigning a statistical weight proportional to $\mu=\cos\alpha$8 (Husidic et al., 2024).

A noteworthy refinement appears in the 2012 ESP paper: the foreshock and sheath are represented by a spatially varying parallel mean free path that decreases toward the shock and is estimated from in-situ observations. The shock-adjacent value inferred for 88 keV protons is about $\mu=\cos\alpha$9 au, and the authors use this observation-tuned profile to strengthen scattering in the foreshock and sheath (Wijsen et al., 2022). By contrast, the 2023 study explicitly notes that true diffusive shock acceleration is not especially efficient in its transport setup because the mean free path is not reduced near the shock enough to trap particles strongly there (Wijsen et al., 2023).

5. Demonstrated applications

PARADISE has been used in real-event modelling, method validation, and coronal transport studies.

In the observation-based modelling of the 14 July 2012 energetic storm particle event, the coupled EUHFORIA+PARADISE chain reproduces several low-energy features observed by ACE near Earth. For energies below about $t$0 MeV, the simulation captures the ESP onset time, the pre-shock intensity rise, the spectral shape at the shock, a secondary peak just before flux-rope onset, and the sharp drop in intensities at entry into the magnetic cloud. The paper emphasizes that the downstream drop is reproduced because the CME is represented as a magnetised structure using EUHFORIA’s spheromak model, illustrating the importance of a magnetised CME description for flux-rope shielding (Wijsen et al., 2022).

In the 2021 October 9 gradual SEP event, PARADISE is applied to protons from 50 keV to 6 MeV, with emphasis on the below-5 MeV ions measured by ACE near Earth and BepiColombo at 0.33 au. The study finds that a modest high-speed stream affected the low-energy ion component by deforming the CME-driven shock and making its compression ratio and obliquity highly non-uniform. When the particle source strength is scaled to the local solar-wind compression at the shock, the model achieves an excellent match to the measured $t$1 MeV ions. At ACE, it reproduces the first peak, the near-constant low-energy intensities while the spacecraft remains in the HSS and magnetically connected to the compressed region, and the general decline in the sheath. At BepiColombo, it reproduces a two-peak structure before shock arrival; the modelled time series is shifted by 4 hours earlier to align with the observed shock arrival, indicating that the simulated shock arrival was late by 4 hours at Bepi (Wijsen et al., 2023).

The Icarus+PARADISE paper is primarily a validation and resolution study rather than a real-event reconstruction. Using a synthetic CIR-rich wind, it reproduces the expected double-peaked intensity profile at $t$2 au and $t$3 au, with the right peak associated with the forward shock, the left peak with the reverse shock, and the valley aligned with the stream interface. The comparison with EUHFORIA+PARADISE shows good agreement in the lowest-energy channels, while higher-energy intensities differ somewhat because Icarus is described as slightly more diffusive in that setup (Husidic et al., 2024).

The COCONUT+PARADISE study extends the code to coronal flux-rope confinement. It injects 100 keV protons, isotropic in pitch angle, near the central axis of one negative-polarity leg of the Titov–Démoulin flux rope at $t$4 h after flux-rope insertion, with 3.6 million pseudo-particles and absorbing boundaries. Without cross-field diffusion, particles remain trapped within the TDFR, stream upward along interior field lines, and may later reach the opposite footpoint and fall back, but they do not escape onto the open field lines ahead of the CME. With cross-field diffusion, even relatively small perpendicular mean free paths allow particles on the outer layers to escape, preferentially in the CME propagation direction where reconnection opens field lines near the nose (Husidic et al., 2024).

6. Interpretation, limitations, and development trajectory

The published studies present PARADISE as a capable but deliberately non-self-consistent transport framework. In the Icarus paper, energetic particles are treated under the test-particle assumption, so they do not feed back on the background MHD flow (Husidic et al., 2024). In the 2012 ESP-event study, the authors explicitly identify the lack of self-consistent wave–particle coupling as a limitation: the parallel mean free path is prescribed and observation-tuned rather than evolved in response to the energetic-particle distribution, and they suggest coupling to a wave–particle code such as SOLPACS as a future improvement (Wijsen et al., 2022).

Another recurrent limitation concerns shock injection and acceleration physics. In the 2023 gradual SEP-event paper, the source term depends on compression but does not explicitly include a formula depending on shock obliquity $t$5, even though many acceleration theories do; the authors intentionally avoid imposing an explicit obliquity dependence because the exact physics is uncertain and the observations are already reproduced well with a compression-based source (Wijsen et al., 2023). The 2012 study similarly treats injection phenomenologically rather than from first principles (Wijsen et al., 2022). These choices should not be confused with a claim that PARADISE derives injection efficiencies microscopically; the papers instead use compression-scaled source terms as practical surrogates.

Resolution effects in the background MHD solution are also important. The Icarus+PARADISE study shows that AMR strongly affects shock acceleration: as the AMR level increases, the shock becomes better resolved and high-energy intensities increase, especially in higher-energy channels and at the forward shock. The paper interprets this as a consequence of the ratio of particle mean free path to numerical shock width, noting that MHD shocks are typically numerically broader than real interplanetary shocks and that better shock resolution favors more efficient acceleration (Husidic et al., 2024). The 2012 ESP paper discusses the same issue from a different angle and therefore lets particles cross the modelled shock scatter free through the numerical shock layer, conserving magnetic moment unless mirrored, to prevent the broadened MHD shock from producing merely compressional acceleration (Wijsen et al., 2022).

The development trajectory of PARADISE is therefore twofold. One direction is the extension from EUHFORIA-based heliospheric event studies to the newer Icarus+PARADISE framework with AMR and grid stretching, which is intended to improve shock representation and resolution control (Husidic et al., 2024). The other is the extension inward from the heliosphere to the low corona in COCONUT+PARADISE, with the explicit long-term goal of a consistent coronal-to-heliospheric particle model (Husidic et al., 2024). A plausible implication is that PARADISE has evolved from a heliospheric SEP and ESP transport code into a more general transport layer for multi-domain space-weather modelling, while retaining its defining traits: stochastic solution of the focused transport equation, test-particle dynamics in prescribed MHD backgrounds, and modular scattering and diffusion prescriptions.