A Deterministic Ionization Algorithm for the OSIRIS Particle-in-Cell Framework

This paper presents the development, validation, and verification of a deterministic collisional ionization algorithm for the OSIRIS particle-in-cell framework that significantly improves ionization rate accuracy by up to two orders of magnitude while maintaining linear computational scaling.

The Gist

Imagine you are trying to simulate a storm inside a computer. In the world of physics, this "storm" is a plasma—a super-hot soup of atoms that have been ripped apart into electrons and ions. This happens in everything from the sun to fusion reactors to the lasers used in medical imaging.

To study these storms, scientists use a digital tool called a Particle-in-Cell (PIC) code. Think of this code as a massive, high-speed video game engine. Instead of rendering dragons and castles, it renders billions of tiny particles (electrons and ions) moving through a grid of space.

The problem? In the real world, when a fast electron zooms past an atom, it can knock an electron loose, turning that neutral atom into a new ion. This is called collisional ionization.

In the old video game engines (the previous versions of these codes), simulating this was like trying to guess if a coin flip would land on heads or tails by just rolling dice. It was a random process (Monte Carlo method). If you didn't roll enough dice (use enough particles), your simulation would be full of "noise"—random glitches that made the results look messy and inaccurate. To get a clear picture, you needed a massive amount of computer power.

The Paper's Big Idea: A Deterministic Algorithm

The authors of this paper, working with the OSIRIS simulation framework, built a new engine. Instead of rolling dice, they built a calculator.

Here is the simple breakdown of what they did, using some everyday analogies:

1. The "Rain Gauge" vs. The "Dice Roll"

• The Old Way (Dice Roll): Imagine you are trying to measure how much rain falls in a city. The old method was to stand on every street corner and flip a coin. If it's heads, it rains; if tails, it doesn't. To get an accurate average, you need thousands of people flipping coins. If you only have a few, your data is noisy and unreliable.
• The New Way (Rain Gauge): The new OSIRIS algorithm is like placing a giant, precise rain gauge on every street corner. Instead of guessing, it calculates exactly how much rain should fall based on the clouds above. It doesn't rely on luck; it relies on math. This is called a deterministic approach.

2. The "Grid" and the "Deposit"

In the simulation, space is divided into a grid (like a chessboard).

• The Old Way: Particles would bump into each other randomly. If two particles happened to be close, the code would check if they collided.
• The New Way: The code looks at every single particle and asks, "How much ionization are you causing right now?" It takes that answer and deposits it onto the grid squares, like pouring water into a bucket.
  • Because it calculates the total amount of water (ionization) for the whole bucket at once, it doesn't matter if you have 10 particles or 10,000. The result is smooth and accurate.
  • The Result: The authors found this new method is 100 times more accurate (two orders of magnitude) than the old random methods, especially when there aren't many particles to work with.

3. The "Toll Booth" for Energy

When an electron knocks another electron loose, it loses some of its own energy (like a car losing speed to push a stalled vehicle).

• The Old Way: The code would randomly decide how much energy was lost, which could lead to weird glitches where energy wasn't conserved.
• The New Way: The code acts like a precise toll booth. It calculates exactly how much "speed" (momentum) the fast electron loses and exactly how much "speed" the new, slow electron gains. It keeps a running ledger on the grid to ensure the math always adds up perfectly.

4. Two Modes: "Sticky" vs. "Moving" Ions

The paper describes two ways to handle the atoms being ionized:

• Immobile Ions (The "Sticky" Mode): Imagine the atoms are painted on the floor. They don't move, but they can change color (charge state) when hit. This is great for short, fast simulations where the atoms don't have time to move.
• Mobile Ions (The "Moving" Mode): Imagine the atoms are people walking around a room. When they get hit, they change their "charge" (maybe they put on a different colored hat) and keep walking. The new algorithm is smart enough to track these walking people, updating their "hat color" and giving them a little push when they get hit, all without losing track of the physics.

Why Does This Matter?

Think of this new algorithm as upgrading from a hand-drawn sketch to a high-definition 3D render.

• Accuracy: It removes the "static" or "snow" from the TV screen, giving scientists a crystal-clear view of how plasmas behave.
• Efficiency: Because it's so accurate, scientists don't need to use as many particles to get a good result. This saves massive amounts of computer time and money.
• Versatility: It works for everything from the cold, sparse gases in space (astrophysics) to the super-hot, dense fuel in fusion reactors.

In a Nutshell:
The authors took a messy, random guessing game and turned it into a precise, mathematical calculation. They built a system that knows exactly how much energy is transferred and how many new particles are created, ensuring that the digital simulations of our universe are as accurate as possible, even when the computer resources are limited.

Technical Summary

1. Problem Statement

Particle-in-Cell (PIC) simulations are essential for modeling plasma dynamics in astrophysics, fusion, and laser-plasma interactions. However, standard PIC methods often lack accurate treatments of collisional ionization (ionization caused by electron-impact collisions) and field ionization (ionization caused by strong electric fields).

• Limitations of Existing Methods: Most existing PIC codes (e.g., Smilei, Epoch) use Monte-Carlo (MC) approaches for collisional ionization. These methods pair macro-particles and use random number generation to determine if an ionization event occurs.
  • Noise: MC methods introduce statistical noise, requiring thousands of macro-particles per cell to achieve convergence.
  • Accuracy: Inaccuracies arise when colliding macro-particles have large discrepancies in weight or energy.
  • Time-Step Constraints: Some MC implementations struggle with large time steps, potentially allowing a single macro-particle to undergo multiple ionization events within one step, leading to unphysical density evolution.
• The Gap: There is a need for a deterministic, grid-based algorithm that improves accuracy, reduces noise at low particle counts, and scales efficiently for modeling complex plasma regimes ranging from cold astrophysical settings to hot, dense fusion systems.

2. Methodology

The authors developed a new deterministic ionization algorithm integrated into the OSIRIS PIC framework. Unlike MC methods that rely on particle-particle interactions, this algorithm deposits ionization rates onto a spatial grid and advances ion densities deterministically.

A. Theoretical Framework

• Collisional Ionization: Uses Modified Relativistic Binary Encounter Bethe (MRBEB) cross-sections. These cross-sections are calculated based on incident electron kinetic energy, binding energies, and shell occupation numbers. The model accounts for:
  • Energy loss of the incident electron.
  • Energy transfer to the newly ionized (secondary) electron.
  • Momentum balance (assuming collinear momentum transfer for simplicity).
• Field Ionization: Implements a piecewise model covering three regimes:
  1. Tunneling Regime: Uses the ADK model (instantaneous rate).
  2. Intermediate Regime: Uses a barrier-suppression model (RBM).
  3. Barrier-Suppression Regime: Uses a classical/quantum analytic model for extremely high fields.
  • Transition thresholds ($E_1, E_2$) are pre-calculated to ensure continuity.

B. Algorithmic Implementation

The algorithm operates within the PIC loop in three main steps:

1. Rate Calculation:
   • Grid-Based Approach: Instead of pairing particles, the algorithm iterates over all macro-particles and deposits their contribution to the ionization rate ($R_{CI}$) onto a grid.
   • Deterministic Advancement: The ion number densities ($n_i$) are advanced using a semi-implicit, center-differenced scheme (Runge-Kutta style) to solve the coupled differential equations:
     $$\frac{dn_i}{dt} = n_{i-1}R_{i-1} - n_i R_i$$
   • Momentum Tracking: The algorithm calculates the momentum lost by incident electrons and the momentum transferred to secondary electrons, depositing these values onto momentum grids.
2. Ion Density Advancement:
   • Immobile Ions: Ionization states are tracked via grid variables ($Z_i$) representing the fraction of density in each charge state.
   • Mobile Ions: A novel infrastructure tracks ionizable moving ions using macro-particles with a continuous ionization variable ($\zeta$). When $\zeta$ drops below a random threshold, the particle is "promoted" to the next charge state.
3. Particle Injection:
   • New macro-particles (electrons or promoted ions) are injected when sufficient charge density has accumulated.
   • Momentum Assignment: Newly created particles receive momentum based on the accumulated values in the momentum grids, ensuring energy and momentum conservation.

3. Key Contributions

• Deterministic vs. Stochastic: The shift from Monte-Carlo pairing to a grid-based deterministic deposition eliminates statistical noise in ionization rates, allowing for accurate simulations with significantly fewer macro-particles per cell.
• Dual Infrastructure: The code supports both immobile ions (fixed background) and mobile ions (dynamic tracking of charge states), the latter being a significant advancement for modeling ion motion in plasma.
• Improved Momentum Conservation: The algorithm explicitly tracks momentum loss and transfer, using a Runge-Kutta scheme to prevent numerical instability (e.g., removing more momentum than a particle possesses) and ensuring physical consistency.
• Comprehensive Field Ionization: The implementation includes a seamless transition between tunneling, intermediate, and barrier-suppression ionization regimes, crucial for high-intensity laser-plasma interactions.

4. Results and Verification

The algorithm was verified and validated against theoretical models and other PIC codes (Smilei and Epoch).

• Collisional Ionization Rates:
  • OSIRIS results matched theoretical numerical integrations of the rate equations with minimal noise, even with low particle counts (e.g., 10 particles/cell).
  • In contrast, Monte-Carlo codes (Smilei/Epoch) required ~10,000 particles/cell to achieve similar convergence.
• Momentum Loss (Stopping Power):
  • The energy loss of incident electrons matched the theoretical stopping power formulation (Rohrlich & Carlson) for Hydrogen and Aluminum across various energies (1 keV to 1 MeV).
• Momentum Transfer:
  • The average momentum assigned to newly ionized electrons matched theoretical predictions.
  • Convergence: The deterministic method showed a sharp drop in error as particle count increased, whereas MC methods followed the expected $1/\sqrt{N}$ statistical scaling.
• Benchmarking:
  • When compared to Smilei and Epoch, OSIRIS showed excellent agreement once the "multiple ionizations per timestep" artifact in the MC codes was disabled.
  • OSIRIS remained accurate at larger time steps ($\Delta t = 5/\omega_p$) where MC codes began to deviate due to numerical constraints.
• Performance:
  • Execution time scales linearly with the number of particles per cell, similar to MC codes.
  • Currently, OSIRIS is slightly slower due to lack of vectorization (SIMD) in the ionization module, but the authors note this can be optimized.

5. Significance

This work provides a robust, high-fidelity tool for simulating plasma dynamics where collisional ionization is critical.

• Efficiency: By reducing the noise floor, the algorithm allows researchers to simulate complex plasma scenarios (e.g., ionization fronts, filamentation) with fewer macro-particles, reducing computational cost.
• Accuracy: The deterministic approach eliminates the "noise" that can seed unphysical instabilities in ionization fronts, a common issue in MC-based PIC simulations.
• Extensibility: The framework is designed to be modular, allowing for the future inclusion of recombination processes (radiative, three-body) and non-collinear momentum scattering, making it a versatile platform for next-generation plasma physics research.

In summary, the paper presents a significant advancement in PIC methodology by replacing stochastic collisional ionization with a deterministic, grid-based approach, achieving superior accuracy and convergence properties while maintaining compatibility with existing PIC infrastructure.