Stochastic Particle Method Research Articles

Stochastic variational principles for the collisional Vlasov-Maxwell and Vlasov-Poisson equations.

In this work, we recast the collisional Vlasov–Maxwell and Vlasov–Poisson equations as systems of coupled stochastic and partial differential equations, and we derive stochastic variational principles which underlie such reformulations. We also propose a stochastic particle method for the collisional Vlasov–Maxwell equations and provide a variational characterization of it, which can be used as a basis for a further development of stochastic structure-preserving particle-in-cell integrators.

An efficient algorithm of the unified stochastic particle Bhatnagar-Gross-Krook method for the simulation of multi-scale gas flows

The unified stochastic particle method based on the Bhatnagar-Gross-Krook model (USP-BGK) has been proposed recently to overcome the low accuracy and efficiency of the traditional stochastic particle methods, such as the direct simulation Monte Carlo (DSMC) method, for the simulation of multi-scale gas flows. However, running with extra virtual particles and space interpolation, the previous USP-BGK method cannot be directly transplanted into the existing DSMC codes. In this work, the implementation of USP-BGK is simplified using new temporal evolution and spatial reconstruction schemes. As a result, the present algorithm of the USP-BGK method is similar to the DSMC method and can be implemented efficiently based on any existing DSMC codes just by modifying the collision module.

Stochastic Reaction-Diffusion Modeling of Calcium Dynamics in 3D-Dendritic Spines of Purkinje Cells

Calcium (Ca2+) is a second messenger assumed to control changes in synaptic strength in the form of both long-term depression (LTD) and long-term potentiation (LTP) at Purkinje cell dendritic spine synapses via inositol trisphosphate (IP3) induced Ca2+ release. These Ca2+ transients happen in response to stimuli from parallel fibers (PF) from granule cells and climbing fibers (CF) from the inferior olivary nucleus. These events occur at low numbers of free Ca2+ requiring stochastic single particle methods when modeling them.

A hybrid particle approach based on the unified stochastic particle Bhatnagar-Gross-Krook and DSMC methods

For the simulation of multiscale gas flows, the numerical scheme should be valid and efficient in both rarefied and continuum regimes. For example, the Direct Simulation Monte Carlo (DSMC) method is not appropriate because of its huge computational cost in the continuum regime. Various kinds of hybrid methods with DSMC have been proposed to deal with this difficulty. One of them is the particle-particle hybrid method, which combines DSMC with another stochastic particle method that is based on an approximate kinetic model, such as the BGK-DSMC hybrid method. Since the same kind of computational particles is applied in the stochastic particle BGK (SPBGK) and DSMC methods, they can be implemented in a unified particle framework. Therefore, compared to the CFD-DSMC hybrid method, the BGK-DSMC hybrid method avoids difficulties caused by the amalgamation of two fundamentally different types of solvers. However, the traditional SPBGK method decouples the molecular motions and collisions in analogy to DSMC, and hence its transport properties deviate from physical values as the time step size increases. This defect significantly affects its computational accuracy and efficiency in the continuum regime. In the present paper, instead of the traditional SPBGK method, a unified stochastic particle BGK (USPBGK) method, which has second-order accuracy in time and space, is combined with DSMC. Comparing the computational performance of the USPBGK-DSMC and BGK-DSMC hybrid methods for numerical test cases of hypersonic flows past a cylinder and plume flows out of a planar micronozzle, one note that the proposed USPBGK-DSMC hybrid method can achieve higher efficiency for such multiscale gas flow simulations. In addition, an equilibrium-breakdown criterion based on the ratio of the time scales of the macroscopic flow field and molecule collision is proposed as a switching criterion for the new hybrid method. This criterion depends on the local Knudsen number and is easy to implement in a unified particle framework.

Ray effect in rarefied flow simulation

Ray effect usually appears in the computation of radiative transfer when discrete ordinates method (DOM) is used. The cause and remedy for the ray effect have been intensively investigated in the radiation community. For rarefied gas flow, the ray effect is also associated with the discrete velocity method (DVM). However, few studies have been carried out in the rarefied community. In this paper, we make a detailed investigation of the ray effect in the rarefied flow simulation. Starting from a few commonly used benchmark tests, the root of the ray effect will be analyzed theoretically and validated numerically. The mitigation strategy by adjusting the discretization in the particle velocity space is studied as well. The guidelines are proposed to optimize the velocity space discretization for effectively reducing the ray effect. However, the design of optimal velocity space discretization is problem dependent and it can be hardly obtained in a highly rarefied flow simulation with complicated geometry. Due to the intrinsic self-adaptation of particle velocity, the stochastic particle method could automatically provide an optimized velocity discretization and avoid the ray effect. In the term of ray effect mitigation, the use of stochastic particle method is recommended in the numerical simulation of highly rarefied flow.

Simplified unified wave-particle method with quantified model-competition mechanism for numerical calculation of multiscale flows.

A quantified model-competition (QMC) mechanism for multiscale flows is extracted from the integral (analytical) solution of the Boltzmann-BGK model equation. In the QMC mechanism, the weight of the rarefied model and the weight of the continuum (aerodynamic and hydrodynamic) model are quantified. Then, a simplified unified wave-particle method (SUWP) is constructed based on the QMC mechanism. In the SUWP, the stochastic particle method and the continuum Navier-Stokes method are combined together. Their weights are determined by the QMC mechanism quantitatively in every discrete cell of the computational domain. The validity and accuracy of the present numerical method are examined using a series of test cases including the high nonequilibrium shock wave structure case, the unsteady Sod shock-tube case with a wide range of Kn number, the hypersonic flow around the circular cylinder from the free-molecular regime to the near continuum regime, and the viscous boundary layer case. In the construction process of the present method, an antidissipation effect in the continuum mechanism is also discussed.

A unified stochastic particle Bhatnagar-Gross-Krook method for multiscale gas flows

The stochastic particle method based on Bhatnagar-Gross-Krook (BGK) or ellipsoidal statistical BGK (ESBGK) model approximates the pairwise collisions in the Boltzmann equation using a relaxation process. Therefore, it is more efficient to simulate gas flows at small Knudsen numbers than the counterparts based on the original Boltzmann equation, such as the Direct Simulation Monte Carlo (DSMC) method. However, the traditional stochastic particle BGK method decouples the molecular motions and collisions in analogy to the DSMC method, and hence its transport properties deviate from physical values as the time step size increases. This defect significantly affects its computational accuracy and efficiency for the simulation of multiscale flows, especially when the transport processes in the continuum regime is important. In the present paper, we propose a unified stochastic particle ESBGK (USP-ESBGK) method by combining the molecular convection and collision effects. In the continuum regime, the proposed method can be applied using large temporal-spatial discretization and approaches to the Navier-Stokes solutions with second-order accuracy. Furthermore, it is capable to simulate both the small scale non-equilibrium flows and large scale continuum flows within a unified framework efficiently and accurately. The applications of USP-ESBGK method to a variety of benchmark problems, including Couette flow, thermal Couette flow, Poiseuille flow, Sod tube flow, cavity flow, Taylor-Green vortex, and flow through a slit, demonstrated that it is a promising tool to simulate multiscale gas flows ranging from rarefied to continuum regime.

Unified gas-kinetic wave-particle methods. II. Multiscale simulation on unstructured mesh

In this paper, we present a unified gas-kinetic wave-particle (UGKWP) method on unstructured mesh for the multiscale simulation of continuum and rarefied flow. Inheriting from the multiscale transport in the unified gas-kinetic scheme (UGKS), the integral solution of the kinetic model equation is employed in the construction of the UGKWP method to model the flow physics on the scales of cell size and time step. A novel wave-particle adaptive formulation is introduced in the UGKWP method to describe the flow dynamics in each control volume. The local gas evolution is constructed through the dynamical interaction of the deterministic hydrodynamic wave and the stochastic kinetic particle. To model the gas dynamics on the scales of cell size and time step, the decomposition, interaction, and evolution of the hydrodynamic wave and the kinetic particle depend on the ratio of time step to local collision time. In the rarefied flow regime, the UGKWP method recovers the nonequilibrium flow physics by discrete particles and performs as a stochastic particle method. In the continuum flow regime, the UGKWP method captures the flow behavior solely by macroscopic variable evolution and becomes a gas-kinetic hydrodynamic flow solver, the same as the gas-kinetic scheme, for viscous and heat-conducting Navier–Stokes solutions. In the transition regime, both kinetic particle and hydrodynamic wave contribute adaptively in the UGKWP to capture the peculiar nonequilibrium flow physics in a most efficient way. In different flow regimes, the Sod shock tube, lid-driven cavity flow, laminar boundary layer, and high-speed flow around a circular cylinder are computed to validate the UGKWP method on unstructured mesh. The UGKWP method obtains the same UGKS solutions in all Knudsen regimes. However, with an automatic wave-particle decomposition, the UGKWP method becomes very efficient. For example, at Mach number 30 and Knudsen number 0.1, the UGKWP has several-order-of-magnitude reductions in computational cost and memory requirement in comparison with UGKS. Overall, the UGKWP can capture the gas dynamics in all flow regimes efficiently and accurately from the free molecular transport to the Navier-Stokes flow evolution.

A blob method for diffusion

As a counterpoint to classical stochastic particle methods for diffusion, we develop a deterministic particle method for linear and nonlinear diffusion. At first glance, deterministic particle methods are incompatible with diffusive partial differential equations since initial data given by sums of Dirac masses would be smoothed instantaneously: particles do not remain particles. Inspired by classical vortex blob methods, we introduce a nonlocal regularization of our velocity field that ensures particles do remain particles and apply this to develop a numerical blob method for a range of diffusive partial differential equations of Wasserstein gradient flow type, including the heat equation, the porous medium equation, the Fokker–Planck equation, and the Keller–Segel equation and its variants. Our choice of regularization is guided by the Wasserstein gradient flow structure, and the corresponding energy has a novel form, combining aspects of the well-known interaction and potential energies. In the presence of a confining drift or interaction potential, we prove that minimizers of the regularized energy exist and, as the regularization is removed, converge to the minimizers of the unregularized energy. We then restrict our attention to nonlinear diffusion of porous medium type with at least quadratic exponent. Under sufficient regularity assumptions, we prove that gradient flows of the regularized porous medium energies converge to solutions of the porous medium equation. As a corollary, we obtain convergence of our numerical blob method. We conclude by considering a range of numerical examples to demonstrate our method’s rate of convergence to exact solutions and to illustrate key qualitative properties preserved by the method, including asymptotic behavior of the Fokker–Planck equation and critical mass of the two-dimensional Keller–Segel equation.

Fully conservative spectral Galerkin–Petrov method for the inhomogeneous Boltzmann equation

In this paper, we present an application of a Galerkin-Petrov method to the spatially one-dimensional Boltzmann equation. The three-dimensional velocity space is discretised by a spectral method. The space of the test functions is spanned by polynomials, which includes the collision invariants. This automatically insures the exact conservation of mass, momentum and energy. The resulting system of hyperbolic PDEs is solved with a finite volume method. We illustrate our method with two standard tests, namely the Fourier and the Sod shock tube problems. Our results are validated with the help of a stochastic particle method.

Coupling of state-resolved rovibrational coarse-grain model for nitrogen to stochastic particle method for simulating internal energy excitation and dissociation.

We propose to couple a state-resolved rovibrational coarse-grain model to a stochastic particle method for simulating internal energy excitation and dissociation of a molecular gas. A coarse-grained model for a rovibrational reaction mechanism of an ab initio database developed at the NASA Ames Research Center for the N2-N system is modified based on variably spaced energy bins. The thermodynamic properties of the modified coarse-grained model allow us to closely match those obtained with the full set of rovibrational levels over a wide temperature range, while using a number of bins significantly smaller than the complete mechanism. The chemical-kinetic behavior of equally and variably spaced bin formulations is compared by simulating internal energy excitation and dissociation of nitrogen in an adiabatic, isochoric reactor. We find that the variably spaced formulation is better suited for reproducing the dynamics of the full database at conditions of interest in the Earth atmospheric entry. Also in this paper, we discuss the details of our particle method implementation for the uniform rovibrational collisional bin model and describe changes to the Direct Simulation Monte Carlo (DSMC) collision algorithm, which become necessary to accommodate our state-resolved reaction mechanism for excitation and dissociation reactions. The DSMC code is then verified against equivalent master equation calculations. In these simulations, state-resolved cross sections are used in analytical form. These cross sections verify micro-reversibility relations for the rovibrational bins and allow for fast execution of the DSMC code. In our verification calculations, we obtain very close agreement for the concentrations profiles of N and N2, as well as the translational and rovibrational mode temperatures obtained independently through both methods. In addition to macroscopic moments, we compare discrete internal energy populations predicted at selected time steps via DSMC and the master equations. We observe good agreement between the two sets of results within the limits imposed by statistical scatter, which is inherent to particle-based DSMC solutions. As future work, the rovibrational coarse-grain model coupled to the particle method will allow us to study 3D reentry flow configurations.

A stochastic algorithm without time discretization error for the Wigner equation

Stochastic particle methods for the numerical treatment of the Wigner equation are considered. The approximation properties of these methods depend on several numerical parameters. Such parameters are the number of particles, a time step (if transport and other processes are treated separately) and the grid size (used for the discretization of the position and the wave-vector). A stochastic algorithm without time discretization error is introduced. Its derivation is based on the theory of piecewise deterministic Markov processes. Numerical experiments are performed in a one-dimensional test case. Approximation properties with respect to the grid size and the number of particles are studied. Convergence of a time-splitting scheme to the no-splitting algorithm is demonstrated. The no-splitting algorithm is shown to be more efficient in terms of computational effort.

A single-column particle-resolved model for simulating the vertical distribution of aerosol mixing state: WRF-PartMC-MOSAIC-SCM v1.0

Abstract. The PartMC-MOSAIC particle-resolved aerosol model was previously developed to predict the aerosol mixing state as it evolves in the atmosphere. However, the modeling framework was limited to a zero-dimensional box model approach without resolving spatial gradients in aerosol concentrations. This paper presents the development of stochastic particle methods to simulate turbulent diffusion and dry deposition of aerosol particles in a vertical column within the planetary boundary layer. The new model, WRF-PartMC-MOSAIC-SCM, resolves the vertical distribution of aerosol mixing state. We verified the new algorithms with analytical solutions for idealized test cases and illustrate the capabilities with results from a 2-day urban scenario that shows the evolution of black carbon mixing state in a vertical column.

A stochastic Lagrangian micromixing model for the dispersion of reactive scalars in turbulent flows: role of concentration fluctuations and improvements to the conserved scalar theory under non-homogeneous conditions

Several reaction schemes, based on the conserved scalar theory, are implemented within a stochastic Lagrangian micromixing model to simulate the dispersion of reactive scalars in turbulent flows. In particular, the formulation of the reaction-dominated limit (RDL) reaction scheme is here extended to improve the model performance under non-homogeneous conditions (NHRDL scheme). The validation of the stochastic model is obtained by comparison with the available measurements of reactive pollutant concentrations in a grid-generated turbulent flow. This test case describes the dispersion of two atmospheric reactant species (NO and O3) and their reaction product (NO2) in an unbounded turbulent flow. Model inter-comparisons are also assessed, by considering the results of state-of-the-art models for pollutant dispersion. The present validation shows that RDL reaction scheme provides a systematic overestimation (relative error of ca. 85% around the centreline) in computing the local reactant consumption/production rate, whereas the NHRDL scheme drastically reduces this gap (relative error lower than 5% around the centreline). In terms of NO2 production (or reactant consumption), neglecting concentration fluctuations determines overestimations of the product mean of around 100% and a NO2 local production of one order of magnitude higher than the reference simulation. In terms of standard deviations, the concentration fluctuations of both the passive and reactive scalars are generally of the same order of magnitude or up to 1 or 2 orders of magnitudes higher than the corresponding ensemble mean values, except for the background reactant close to the plume edges. The study highlights the importance of modelling pollutant reactions depending on the instantaneous instead of the mean concentrations of the reactants, thus quantifying the role of the turbulent fluctuations of concentration, in terms of scalar statistics (mean, standard deviation, intensity of fluctuations, skewness and kurtosis of concentration, segregation coefficient, simulated reaction rate). This stochastic particle method represents an efficient numerical technique to solve the convection–diffusion equation for reactive scalars and involves several application fields: micro-scale air quality (urban and street-canyon scales), accidental releases, impact of odours, water quality and fluid flow industrial processes (e.g. combustion).

Importance sampling variance reduction for the Fokker–Planck rarefied gas particle method

The Fokker–Planck approximation to the Boltzmann equation, solved numerically by stochastic particle schemes, is used to provide estimates for rarefied gas flows. This paper presents a variance reduction technique for a stochastic particle method that is able to greatly reduce the uncertainty of the estimated flow fields when the characteristic speed of the flow is small in comparison to the thermal velocity of the gas. The method relies on importance sampling, requiring minimal changes to the basic stochastic particle scheme. We test the importance sampling scheme on a homogeneous relaxation, planar Couette flow and a lid-driven-cavity flow, and find that our method is able to greatly reduce the noise of estimated quantities. Significantly, we find that as the characteristic speed of the flow decreases, the variance of the noisy estimators becomes independent of the characteristic speed.

A stochastic model and Monte Carlo algorithm for fluctuation-induced H2 formation on the surface of interstellar dust grains

A stochastic algorithm for simulation of fluctuation-induced kinetics of H2 formation on grain surfaces is suggested as a generalization of the technique developed in our recent studies [1] where this method was developed to describe the annihilation of spatially separate electrons and holes in a disordered semiconductor. The stochastic model is based on the spatially inhomogeneous, nonlinear integro-differential Smoluchowski equations with random source term. In this paper we derive the general system of Smoluchowski type equations for the formation of H2 from two hydrogen atoms on the surface of interstellar dust grains with physisorption and chemisorption sites. We focus in this study on the spatial distribution, and numerically investigate the segregation in the case of a source with a continuous generation in time and randomly distributed in space. The stochastic particle method presented is based on a probabilistic interpretation of the underlying process as a stochastic Markov process of interacting particle system in discrete but randomly progressed time instances. The segregation is analyzed through the correlation analysis of the vector random field of concentrations which appears to be isotropic in space and stationary in time.

The Wigner equation in the presence of electromagnetic potentials

An analysis of the possible formulations of the Wigner equation under a general gauge for the electric field is presented with an emphasis on the computational aspects of the problem. The numerical peculiarities of those formulations enable alternative computational strategies based on existing numerical methods applied in the Wigner formalism, such as finite difference or stochastic particle methods. The phase space formulation of the problem along with certain relations to classical mechanics offers an insight about the role of the gauge transforms in quantum mechanics.

Cell Size Error in Stochastic Particle Methods for Coagulation Equations with Advection

The error arising from the delocalization of the coagulation interaction in stochastic particle methods is studied. The model under consideration includes advection, coagulation, and inception. Stochastic particle methods depend on several discretization parameters, due to the finite number of particles, the decoupling of transport and interaction (time step), and the spatial delocalization of the interaction (cell size). The paper studies the dependence on the cell size of the steady state solution obtained in the infinite-particle-number limit. Sufficient conditions for second order approximation are provided. Examples are given, where only first order approximation is observed.

Probability density function approach for modelling multi-phase flow with ganglia in porous media

AbstractA probabilistic approach to model macroscopic behaviour of non-wetting-phase ganglia or blobs in multi-phase flow through porous media is proposed. The key idea is to consider a set of stochastic Markov processes that can mimic the microscopic multi-phase dynamics. These processes are characterized by equilibrium probability density functions (PDFs) and correlation times, which can be obtained from micro-scale simulation studies or experiments. A Lagrangian viewpoint is adopted, where stochastic particles represent infinitesimal fluid elements and evolve in the physical and probability space. Ganglion mobilization and trapping are modelled by a two-state jump process with transition probabilities given as functions of ganglion size. Coalescence and breakup of ganglia influence the ganglion size distribution, which is modelled by a Langevin type equation. The joint probability density function (JPDF) of the chosen stochastic variables is governed by a high-dimensional Chapman–Kolmogorov equation. This equation can be used to derive moment (e.g. saturation, mean mobility etc.) transport equations, which in general do not form a closed system. However, in some special cases, which arise in the limit of one time scale being smaller or larger than the others, a closed set of moment transport equations can be obtained. For slowly varying and quasi-uniform flows, the saturation transport equation appears in closed form with the mean mobility fully determined, if the equilibrium PDFs are known. Furthermore, it is shown how statistical parameters such as mobilization and trapping rates and equilibrium PDFs can be obtained from the birth–death type approach, in which ganglia breakup and coalescence are explicitly considered. A two-equation transport model (one equation for the total saturation and one for the trapped saturation) is obtained in the limit of very fast coalescence and breakup processes. This model is employed to mimic hysteresis in relative permeability–saturation curves; a well known phenomenon observed in the successive processes of imbibition and drainage. For the general case, the JPDF-equation is solved using the stochastic particle method, which was proposed in our previous paper (Tyagi et al. J. Comput. Phys. 227, 2008, 6696–6714). Several one- and two-dimensional numerical simulation results are presented to show the influence of correlation times on the averaged macroscopic flow behaviour.

Stochastic flux freezing and magnetic dynamo

Magnetic flux conservation in turbulent plasmas at high magnetic Reynolds numbers is argued neither to hold in the conventional sense nor to be entirely broken, but instead to be valid in a statistical sense associated to the "spontaneous stochasticity" of Lagrangian particle trajectories. The latter phenomenon is due to the explosive separation of particles undergoing turbulent Richardson diffusion, which leads to a breakdown of Laplacian determinism for classical dynamics. Empirical evidence is presented for spontaneous stochasticity, including numerical results. A Lagrangian path-integral approach is then exploited to establish stochastic flux freezing for resistive hydromagnetic equations and to argue, based on the properties of Richardson diffusion, that flux conservation must remain stochastic at infinite magnetic Reynolds number. An important application of these results is the kinematic, fluctuation dynamo in nonhelical, incompressible turbulence at magnetic Prandtl number (Pr(m)) equal to unity. Numerical results on the Lagrangian dynamo mechanisms by a stochastic particle method demonstrate a strong similarity between the Pr(m)=1 and 0 dynamos. Stochasticity of field-line motion is an essential ingredient of both. Finally, some consequences for nonlinear magnetohydrodynamic turbulence, dynamo, and reconnection are briefly considered.