In this article, we’ll introduce Direct Simulation Monte Carlo (DSMC), which is a method for simulating interparticle collisions. It is frequently used to simulate neutral gas flows and to add collisions to plasma simulations. This method has been popularized by G.A. Bird, and to this day, his book Molecular Gas Dynamics and the Direct Simulation of Gas Flows serves as the bible for the DSMC method. Prof. Bird also runs a website containing various DSMC resources. Please make sure to check out these references for the details of DSMC, the write up here is definitely not exhaustive!
DSMC is the analogue to the particle in cell method for neutral gases. It is a method for computing the motion of neutral flows directly using simulation particles instead of solving fluid conservation equations, as is done in the CFD approach. The main advantage over CFD is that DSMC does not assume a priori any distribution function. The CFD method is generally applicable only to thermalized flows – ones in which gas molecule velocities follow the Maxwellian distribution function. The DSMC method works by moving the simulation particles through a finite time step through a discretized computational domain. After each push, particles are paired up with other particles located in the same computational mesh and collision probability is computed for these pairs. Those pairs that are colliding then have their velocities modified in a way that conserves energy and momentum.
Since by itself, the PIC method does not account for particle collisions, the particle in cell method must be coupled with a collisional algorithm to account for momentum transfer and charge exchange collisions. Previously we have introduced one way of adding collisions: the Monte Carlo Collisions (MCC) method. The problem with MCC is that it is not momentum conserving, and is thus physically suitable only when collisions have only a negligible effect on the target population. This significantly restricts what we can do with it. But no worries, DSMC is here to the rescue!
Quick DSMC Theory Crash Course
Let’s consider a single simulation cell containing a number of simulation macroparticles. Bird gives the probability of a collision for a macroparticle in this cell as
where is the macroparticle weight, is the total collision cross-section, is the relative velocity, is the simulation time step, and is the cell volume. You can think of this expression as
or the ratio of the volume swept by the cross-section to the total volume of the cell scaled by the macroparticle weight.
No Time Counter Method (NTC)
One order to check for collisions would be to iterate over all particles, and compute probability with all remaining particles. This would result in pairs. For a large number of particles, this method clearly becomes very computationally inefficient. Bird’s No Time Counter (NTC) method was designed to help with this issue. It allows us to estimate ahead of time the maximum number of pairs that need to be checked. The maximum collision probability is
where is a parameter chosen ahead of time using approximate predictions of cross-section and velocity. As Bird points out in his book, the actual value is not all that important since it ends up getting cancelled out. The number of pairs to check is then
Here the second is the average particle count. The average is used to reduce the statistical time-step to time-step oscillations. For each pair, we then compute the probability as follows
this value is compared to a random number and the collision occurs if .
Variable Hard Sphere
The collision cross-section for a hard-sphere, billiard ball collision is . This model is not physically realistic for real molecules. Bird proposed a new model which takes into account the decrease in cross-section as relative velocity between molecules increases. The diameter of the sphere is allowed to vary according to
All that’s needed now is the code for actually performing the collision. Bird goes into the details of deriving the equations in Chapter 2 of his book. In short, momentum is conserved in the collision such that
Bird also shows that energy conservation dictates that the magnitude of relative velocity remains constant during the collision. We can also write the post-collision velocity as
Since we know the magnitude of the relative velocity from the pre-collision velocities, we just need to select a new direction for the post-collision velocity vector. In general, the post-collision scattering is anisotropic – the direction will depend on the pre-impact velocities and differential cross-sections are then used to select the in-plane angle. However, here for simplicity we assume that the scattering is isotropic. We then let
where the two R1 and R2 are two random numbers. The relative velocity is then given by
Zero D Example