PIConGPU modeling of nanoplasma formation in helium nanodroplets irradiated by intense femtosecond laser pulses

Abstract: Helium nanodroplets provide a unique and versatile platform for investigating strong-field-driven nanoplasma dynamics. In this work, we present large-scale, GPU-accelerated particle-in-cell simulations using \textsc{PIConGPU} to study the interaction of pure helium nanodroplets containing up to $10<sup>{6}$ atoms with intense near-infrared femtosecond laser pulses, and compare the results with single-shot velocity-map electron imaging and ion measurements. The simulations describe the plasma evolution from the first ionization events to collective electron motion, nanoplasma formation, and early expansion. We show that the calculated electron and ion observables reproduce the main features of the measured spectra in systems with similar cluster sizes and laser intensities. Our results demonstrate that \textsc{PIConGPU} captures the essential physics of nanoplasma formation previously addressed mainly with molecular-dynamics or TDDFT approaches, while remaining computationally efficient and applicable to much larger systems. This establishes \textsc{PIConGPU} as a powerful and scalable tool for connecting nanoplasma theory with experimentally accessible observables.

• The paper demonstrates that GPU-accelerated PIConGPU simulations with FLYonPIC atomic kinetics can model nanoplasma formation in helium droplets containing up to 10^6 atoms, extending beyond typical molecular-dynamics limits.
• The simulations reveal a two-stage ignition process involving field-ionization seeding followed by electron-impact avalanche ionization, while coherent electron polarization fades after quasi-neutrality breaks down near 20 fs.
• The model reproduces key experimental features, including electron velocity-map image shapes, increasing electron energies, and ion momentum scales, but remains qualitative because it omits long-term expansion and absolute-yield validation.

Overview

This work presents large-scale, GPU-accelerated particle-in-cell (PIC) simulations, performed with PIConGPU and its FLYonPIC atomic physics library, of nanoplasma formation in pure helium nanodroplets containing up to $10^6$ atoms irradiated by intense 800 nm femtosecond laser pulses. The simulations are benchmarked against single-shot velocity-map imaging (VMI) electron measurements and ion momentum distributions recorded at the MAC endstation of ELI Beamlines. The central claim is that a fully electromagnetic PIC framework—previously applied mainly to relativistic laser–plasma problems—can reproduce the essential features of nanoplasma dynamics that have historically been the domain of molecular-dynamics (MD) and TDDFT approaches, while scaling to system sizes an order of magnitude or more beyond what fully microscopic MD can treat.

Motivation and positioning relative to prior approaches

Laser-driven clusters occupy an intermediate regime between isolated atoms and bulk matter: they combine near-solid density with nearly uniform illumination, since surface screening does not prevent field penetration as it does in extended condensed-phase targets. Helium nanodroplets are a particularly clean platform because of their homogeneous density, simple electronic structure, and controllable size distribution; at the same time, helium's high ionization potential makes avalanche ionization acutely sensitive to seed electrons, intensity, and droplet size—a property exploited extensively in dopant-induced ignition studies (2606.14300).

The methodological gap addressed here is one of scale. Fully microscopic MD simulations resolve short-range many-body interactions at the particle level but become impractical beyond clusters of a few thousand atoms. PIC methods instead evolve macro-particles coupled to grid-resolved self-consistent electromagnetic fields, which naturally accommodates large systems and electromagnetic propagation effects, at the cost of not resolving close charged-particle encounters atomistically. The paper positions PIConGPU—with Monte Carlo binary collisions and FLYonPIC atomic kinetics added on top of the standard PIC cycle—as filling this niche for $>10^6$-atom nanoplasmas.

Simulation methodology

The simulations use the explicit 3D3V electromagnetic PIC scheme in PIConGPU: Yee-solved Maxwell equations, Boris particle pusher, Esirkepov current deposition with exact charge conservation, fourth-order piecewise cubic spline (PQS) shape functions, trilinear field interpolation, and double precision throughout. Unresolved short-range collisions are modeled with a Monte Carlo binary collision plugin using a relativistic dynamic-log operator for all species pairs ($\mathrm{He}_e$–$\mathrm{He}_e$, $\mathrm{He}_e$–$\mathrm{He}_i$, $\mathrm{He}_i$–$\mathrm{He}_i$).

FLYonPIC tracks thirty bound atomic states of helium explicitly, solving non-equilibrium atomic rate equations with adaptive sub-stepping (stability factor $\alpha = 0.1$, up to 200 sub-steps per PIC step). Implemented processes include collisional de-/excitation, spontaneous deexcitation, electron-impact ionization, autonomous ionization, pressure ionization, and field ionization (BSI + ADK), all with Stewart-Pyatt ionization potential depression. The atomic kinetics depend on local plasma conditions without assuming thermal equilibrium, and the solver is energy- and charge-conserving.

The droplet is initialized with a SoftSphere/FCC profile at liquid-helium density ($2.1433 \times 10^{28}\,\mathrm{m^{-3}}$), radius 25 nm, 25 macro-particles per cell per species, on a uniform 0.5 nm grid with $>10^6$0 as. Two initial conditions are compared:

• Neutral ignition: initially neutral droplet driven at $>10^6$1 ($>10^6$2), where avalanche ionization develops self-consistently.
• Pre-seeded nanoplasma: singly pre-ionized droplet driven at $>10^6$3 ($>10^6$4), mimicking dopant-triggered ignition.

Both use 40 fs FWHM Gaussian pulses with plane-wave transverse profiles, so any transverse intensity variation arises only from plasma-induced field redistribution. Runs used 16 GPUs on JUWELS at roughly 36 hours per run—a modest computational cost that underscores the scalability argument.

Nanoplasma formation dynamics

The spatio-temporal evolution of the neutral case proceeds through distinct stages. Before the pulse arrives, the target is inert. At the leading edge ($>10^6$5 fs), first electrons appear via field ionization and a transient quasi-free population builds up; a large fraction remains trapped by the space-charge potential, forming a quasi-neutral plasma that screens the cluster interior. The trapped electrons are driven coherently along the polarization axis, producing pronounced surface-polarized modulations in $>10^6$6. During the main interaction ($>10^6$7 to 35 fs), wave-like field structures inside the droplet correlate with sub-cycle displacement of the electron cloud, consistent with rapid charge separation in a transiently overdense plasma—the authors note that a detailed mode analysis is left outside the scope of the work.

The global diagnostics support a two-stage ignition picture. The electron production rate shows a sharp early peak coinciding with depletion of the neutral population ($>10^6$8), followed by a broader contribution near the pulse maximum associated with growth of the doubly charged fraction via electron-impact ionization. The mean electron density inside the root-mean-square ion radius reaches strongly overdense conditions ($>10^6$9) around $\mathrm{He}_e$0 fs, consistent with strong screening and surface-localized field response. Ionic expansion, quantified by $\mathrm{He}_e$1, remains weak during ignition and becomes appreciable only after substantial net charging—on the simulated $\mathrm{He}_e$2 fs timescale the paper captures only the onset of Coulomb expansion, with the main ionic response expected later.

The angular analysis adds a quantitative signature of collective behavior. The hemisphere asymmetry parameter $\mathrm{He}_e$3 exhibits pronounced oscillations upon initial ionization that are not simply locked to the optical cycle, reflecting the interplay of the driving field and the restoring force of the positively charged core. Critically, these oscillations persist only while quasi-neutrality holds: around $\mathrm{He}_e$4 fs, $\mathrm{He}_e$5 relaxes rapidly toward zero even though the laser is still on, coinciding with electron escape and the breakdown of quasi-neutral polarization dynamics. This provides a clear simulation-level marker of the transition from coherent collective response to expansion-dominated behavior.

The energy histograms show complementary electron and ion behavior. The electron spectrum broadens rapidly during ignition and develops an extended high-energy tail; the median electron kinetic energy peaks near $\mathrm{He}_e$6 fs after the laser peak and then decreases as the plasma expands and hot electrons escape. Ions remain cold during the quasi-neutral stage (median kinetic energy near zero up to $\mathrm{He}_e$7 fs) and gain energy steadily thereafter as net positive charge accumulates. A small ion fraction reaches the keV range at late times, while median ion energies stay in the few-tens-of-eV range—an explicit statement that the high-energy tail contributes only a minor part of the total yield.

The pre-ionized case serves as a controlled entry point into the already-ignited regime: the early $\mathrm{He}_e$8 peak is absent, but subsequent overdense plasma formation, $\mathrm{He}_e$9 growth, and delayed expansion closely resemble the second stage of the neutral case. This supports interpreting neutral-droplet ignition as a two-stage process and demonstrates that the late-stage nanoplasma dynamics are largely insensitive to how seeding occurred—at least within the parameter pairs chosen here, which is itself an assumption worth noting, since the neutral and pre-ionized runs differ in intensity as well as initial charge state.

Comparison with experiment

Experiments were performed at the MAC endstation of ELI Beamlines using helium nanodroplets of average size $\mathrm{He}_e$0 atoms, generated by supersonic expansion at 50 bar through a 5 μm cryogenic nozzle, intersected by focused 800 nm pulses. Electrons were detected by VMI with an MCP/Timepix readout operated at low gain to maintain linear response despite large shot-to-shot fluctuations; characteristic electron energies were extracted from Gaussian fits to radial intensity distributions. Direct VMI reconstruction was not possible for ions due to their high kinetic energies, so ion momenta were reconstructed from Timepix impact positions and times of arrival under a free-flight model, assuming emission from the laser focus and approximate symmetry of the energy distribution so that the detector-limited angular subset represents the full distribution.

The comparison is explicitly framed as qualitative rather than one-to-one. Three observables are compared:

Observable	Simulation	Experiment
Electron VMI morphology	Compact central distribution plus broader hot component	Compact central distribution; hot component partially outside detector acceptance
Max electron KE vs. pulse energy	Characteristic energies extracted with same radial analysis (neutral point divided by three for visibility)	Approximately linear increase with pulse energy
Ion momentum distribution	Low-momentum at $\mathrm{He}_e$1, extending to higher momenta at 106 fs; summed curve shifted to experimental peak	Broad accumulated distribution over many single shots

The experimental electron energies increase approximately linearly with pulse energy, indicating efficient plasma heating at stronger driving fields, consistent with earlier single-shot dopant-induced studies. For ions, the summed early-plus-late-time simulated distribution—shifted to the experimental peak position—is offered as a phenomenological illustration that the measured signal reflects a superposition of nanoplasmas with different heating efficiencies and extraction stages, not a prediction for any single event. The authors also note that continued ion expansion beyond the latest simulated extraction time likely means the experimental distribution corresponds to a later dynamical stage than the earliest simulated snapshot.

Overall, the simulations reproduce the main experimentally accessible features—VMI morphology, characteristic energy ranges, and ion momentum scales—with residual differences attributed to single-shot averaging over varying droplet sizes, ignition conditions, and effective intensities, plus detector and extraction-field effects absent from the model. It should be emphasized that the agreement claimed is at the level of spectral shapes and energy scales; no quantitative absolute-yield or shot-by-shot validation is attempted, and the rescaling of the neutral-run data point by a factor of three signals that simulation and experiment are not yet on a common absolute footing.

Limitations and open questions

Several limitations are acknowledged directly. The PIC treatment does not resolve close charged-particle encounters at atomistic resolution, unlike MicPIC, so short-range correlation and quantum effects relevant to small clusters are outside its scope. The comparison with experiment is qualitative: experiments average over many single-shot events with fluctuating conditions, while each simulation represents an idealized individual case with fixed droplet size and intensity. The ion analysis relies on the assumption of an approximately symmetric energy distribution to generalize from the detector's limited angular acceptance. The simulated timescale ($\mathrm{He}_e$2 fs) captures only the onset of ionic expansion, leaving the full Coulomb explosion and long-timescale evolution unaddressed. The internal wave-like field structures observed during the main interaction are described only phenomenologically; a mode analysis is explicitly deferred. Finally, whether the similarity between neutral and pre-ionized pathways holds across a wider range of intensity–seed combinations remains untested.

Conclusion

This work demonstrates that PIConGPU, augmented with Monte Carlo binary collisions and FLYonPIC atomic kinetics, captures the essential physics of nanoplasma formation in helium nanodroplets—from field-ionization seeding through avalanche ionization, overdense plasma buildup, coherent electron polarization, and the onset of Coulomb-driven expansion—at experimentally relevant sizes of order $\mathrm{He}_e$3 atoms, at a computational cost of tens of GPU-hours per run. The simulations reproduce the principal features of measured electron VMI images, characteristic electron energies, and ion momentum distributions, while MD-class methods remain necessary for small clusters where short-range many-body and quantum effects dominate. The paper leaves open quantitative benchmarking against single-shot statistics, extension to full expansion dynamics, detailed identification of the observed internal field structures, and application to doped droplets and shock-wave phenomena within nanodroplets.