Stochastic particle methods

Numerical simulations are required for a wide range of flow regimes for various applications such as the entry of spacecraft into the atmosphere and electric propulsion systems. For example, the accurate design of thermal protection systems requires a correct estimation of the flow characteristics such as the heat flux that a vehicle must withstand during re-entry. The well-known Direct Simulation Monte Carlo method (DSMC) has been used for such numerical simulations for many years. While the undisturbed flow during atmospheric entry maneuvers may well be in the continuum range, non-equilibrium effects can occur in the wake of the spacecraft. Since the computational effort for the DSMC method increases significantly in the transition and continuum regimes, particle-based continuum methods are a promising solution for the simulation of multi-scale problems in different flow regimes. Recently, many methods have been developed, such as the Bhatnagar-Gross-Krook (BGK) models or the Fokker-Planck (FP) approach (see MEDUSA).

In these applications, the Boltzmann equation is the basis for describing the gas flow on a microscopic level:

This equation contains four terms that relate to different aspects of the gas flow. The leftmost term (I) describes the changes in particle position and velocity over time, term (II), the diffusion term, represents changes in the flow due to the particle motion itself. External forces, as described by term (III), are typically electromagnetic forces in the case of plasmas and can be introduced by coupling with the Particle-In-Cell (PIC) module of PICLas. The last term (IV) on the right-hand side, also called collision integral, takes into account the interactions between the particles.

Since an analytical solution of the Boltzmann equation is not possible and requires several assumptions, numerical methods are used. Particle methods such as the DSMC method, but also the Bhatnagar-Gross-Krook (BGK) or the Fokker-Planck (FP) methods are often used to approximate the Boltzmann collision integral in the simulation of dilute and non-equilibrium gas flows. Modules for all three methods are available in PICLas, whereby a simple coupling of the BGK or FP method with the DSMC method is also possible. Different coupling criteria can be used between the methods.

In such a particle simulation, the physical gas flow is represented by the movement and collision of the simulation particles. All three methods are stochastic in nature and describe the gas distribution function as a linear combination of delta functions. In order to reduce the calculation time, a simulation particle usually represents a larger number of real gas particles.

The particle motion and the interactions due to particle collisions are treated one after the other, with the motion step following Newton's laws. The treatment of particle collisions is different for the three methods mentioned. In a DSMC simulation, the collision pairs for the subsequent interaction are determined first. The collision probability is calculated for each particle pair and compared with a random number in order to decide whether an interaction takes place. If a collision occurs, internal energy can be exchanged or chemical reactions can take place between the two particles. The modeling of these processes in PICLas is based on phenomenological reaction rates.

After the collision step and the recalculation of the particle properties, the macroscopic flow parameters are determined by sampling and averaging the microscopic gas properties within a simulation cell. Various boundary conditions can be defined for the simulation area in order to simulate a gas inflow and outflow as well as the scattering process and catalytic reactions on a solid surface.

Due to the stochastic nature of particle methods, a certain number of simulation particles is required to minimize noise. To ensure physically meaningful results, certain quality factors such as the mean free path of the gas or the collision frequency must also be resolved by the simulation parameters. This is the reason for the increasing calculation time in transition or continuum regimes. Here, BGK and FP methods are less computationally expensive due to the lower resolution restrictions.

In the BGK method, the particles relax in the direction of a target distribution function instead of performing collisions as in the DSMC method in order to approximate the Boltzmann collision term.

In PICLas, the particles are sampled with a certain probability from a target distribution function that depends on the transport coefficients of the gas, the particle density, the temperature and the time step. There are various BGK models with different target distribution functions, e.g. the Ellipsoidal Statistical BGK (ESBGK) or Shakhov BGK (SBGK).

The FP approach simplifies the Boltzmann collision term to a Fokker-Planck equation, which depends on a drift (I) and a diffusion term (II). In contrast to the BGK method, the velocity of each particle in the simulation is relaxed with a stochastic differential equation in order to solve the corresponding Fokker-Planck equation.