Chapter 1 Introduction
1.3. Specific Objectives of the Thesis
The objectives of the thesis are briefly summarized as follows:
1. To verify the implementation of transient sub-cell technique in a parallelized DSMC code (PDSC) [Wu, et al].
2. To simulate driven cavity flows, including M=1.1-4 of the speed of top driven plate, and Kn=10-0.0033, based on the average number density and wall temperature.
3. To discuss the effects of Knudsen Number and Mach number of the driven plate on the flow fields.
The organization of the thesis is stated as follows. Chapter 1 describes the Introduction, Chapter 2 describes the Numerical Method, Chapter 3 describes the verification of implementation of transient sub-cells, and followed by the Results and Discussion, and finally Chapter 4 describes the conclusion and recommendation of future work..
Chapter 2
Numerical Method
2.1. The Boltzmann Equation
The Boltzmann equation is one of the most important transport equations in non-equilibrium statistical mechanics, which deals with systems far from thermodynamics equilibrium. There are some assumptions made in the derivation of the Boltzmann equation which defines limits of applicability. They are summarized as follows:
1. Molecular chaos is assumed which is valid when the intermolecular forces are short range. It allows the representation of the two particles distribution function as a product of the two single particle distribution functions.
2. Distribution functions do not change before particle collision. This implies that the encounter is of short time duration in comparison to the mean free collision time.
3. All collisions are binary collisions.
4. Particles are uninfluenced by intermolecular potentials external to an interaction.
According to these assumptions, the Boltzmann equation is derived and shown as Eq. (2.1)
Meaning of particle phase-space distribution function f is the number of particles with
center of mass located within a small volume d r near the point r , and velocity within a 3
range d u , at time t . 3 F is an external force per unit mass and t is the time and i u is the i molecular velocity. σ is the differential cross section and dΩ is an element of solid angle.
The prime denotes the post-collision quantities and the subscript 1 denotes the collision partner. Meaning of each term in Eq. (2.1) is described in the following;
1. The first term on the left hand side of the equation represents the time variation of the distribution function of the particles (unsteady term).
2. The second term gives the spatial variation of the distribution function (flux term).
3. The third term describes the effect of a force on the particles (force term).
4. The term at right hand side of the equation is called the collision integral (collision term). It is the source of most of the difficulties in obtaining solutions of the Boltzmann equation.
In general, it is very hard to solve the Boltzmann equation directly using numerical method because the difficulties of correctly modeling the integral collision term. Instead, the DSMC method was used to simulated problems involving rarefied gas dynamics, which is the simulation tool used in the current thesis.
2.2. General Description of the Standard DSMC
According to the expected rarefaction caused by the rarefied gas flows, the direct simulation Monte Carlo (DSMC) method, which is a particle-based method developed by
Bird during the 1960s, and was widely used as an efficient technique to simulate rarefied gas regime [Bird, 1976 and Bird, 1994]. In the DSMC method, a large number of particles are generated in the flow field to represent real physical molecules rather than a mathematical foundation and it has been proved that the DSMC method is equivalent to solving the Boltzmann equation [Nanbu, 1986 and Wagner, 1992]. The assumptions of molecular chaos and a dilute gas are shared by both the Boltzmann equation and the DSMC method [Bird, 1976 and Bird, 1994]. An important feature of DSMC is that the molecular motion and the intermolecular collisions are uncoupled over the time intervals that are much smaller than the mean collision time. Both the collision between molecules and the interaction between molecules and solid boundaries are computed on a probabilistic basis and, hence, this method makes extensive use of random numbers. In most practical applications, the number of simulated molecules is extremely small as compared to the number of real molecules. The general procedures of the DSMC method are described in the next section, and the consequences of the computational approximations can be found in Bird [Bird, 1976 and Bird, 1994].
In DSMC, there are three popular molecular collision models for real physical behavior and imitate the real particle collision. They include the Hard Sphere (HS), Variable Hard Sphere (VHS) and Variable Soft Sphere (VSS) molecular models [Bird, 1994]. No time counter (NTC) method is often used to determined the probability of collision within a
samping cell based on the acceptance-rejection method. This method can yield the exact collision rate in both simple gases and gas mixtures, and under either equilibrium or non-equilibrium conditions.
Fig. 2.2 is a general flowchart of the DSMC method. Important steps of the DSMC method include setting up the initial conditions, moving all the simulated particles, indexing (or sorting) all the particles, colliding between particles and sampling the molecules within cells to determine the macroscopic quantities. The details of each step are described in the following:
z Initialization
The first step to use the DSMC method is to set up the geometry and flow conditions. A physical space is discredited into a network of cells and the boundary conditions at domain boundaries have to be assigned according to the flow conditions. An important rule of thumb of selecting the cell size is it should be much smaller than the local mean free path. At the same time the distance of the molecular movement per time step should be smaller than the cell dimension, which is similar to the CFL condition in CFD. After the data of geometry and flow conditions have been read in the code, the number of simulation particles each cell is calculated according to the initial number density and the cell volume. The initial particle velocities are normally assigned to each particle based on the Maxwell-Boltzmann distribution according to the initial velocities and temperature. The position of each particle is
then randomly allocated within the cell.
z Molecular Movement
After initialization process, the molecule begins move one by one, and the molecule moves in a straight line over the time step. If the particle collides with a solid surface, then it has to reflected back based on the boundary conditions at the surface. If the particle passes through a inflow/outflow boundary, then it has to be removed from the simulation immediately. After the particle reaches its final destination, the cell with which the particle affiliated is recorded.
z Indexing
Relation between particles and cells are reordered according to the order of the
increasing particle number in each cell, which is used for collision step.
z Gas Phase Collision
The other most important phase of the DSMC method is the gas phase collision. The DSMC method uses the no time counter (NTC) method to determine the correct collision rate in the collision cell. Number of collision pairs within a cell of volume Vc over a time interval
∆t is calculated by the following equation:
where N and N are fluctuating and average number of simulated particles, respectively.
FN
is the particle weight, which is the number of real particles that a simulated particle
represents. σT and cr are the cross section and the relative speed, respectively. The collision probability for each pair of particle chosen can be expressed as
(σTcr)/(σTcr)max (2.3)
Collision is accepted if the above value for the pair is greater than a random fraction. Each cell is treated independently and the collision partners for interactions are chosen randomly, regardless of their positions within the cell. Collision procedures are summarized as follows:
1. Number of collision pairs is calculated according to the NTC method, Eq. (2.2), in each cell.
2. The first particle is chosen randomly from the list of particles within a collision cell.
3. The other collision partner is also chosen randomly within the same cell.
4. The collision is accepted if the computed probability, Eq. (2.3), is greater than a random number.
5. If the collision pair is accepted, then the post-collision velocities are calculated using the mechanics of elastic collision. If the collision pair is not to collide, then proceed to choose the next collision pair.
6. If the collision pair is polyatomic gas, the translational and internal energy can be redistributed by the Larsen-Borgnakke model [1975], which assumes the state is in equilibrium.
The collision process will be finished until all the collision pairs are handled for all cells and then progress to the next step.
z Sampling
After the particle movement and collision process complete, the particle has updated positions and velocities. Macroscopic flow properties in each cell are assumed to be constant over the cell volume and are sampled from the microscopic properties of all the particles within the cell. Macroscopic properties, including number density, velocities and temperatures, are calculated following the equations as shown below: where n, m are the number density and molecule mass, receptively. c, co, and c’ are the total velocity, mean velocity, and random velocity, respectively. In addition, Ttr, Trot, Tv and Ttot are translational, rotational, vibration and total temperature, respectively. εrot and εv are the rotational and vibration energy, respectively. ζrot and ζv are the number of degree of freedom of rotation and vibration, respectively. If the simulated particle is monatomic gas, the
translational temperature is regarded simply as the total temperature. Vibration effect can be neglect if the temperature of the flow is low enough.
The flow will be monitored if steady state is reached. If the flow is still in the transient phase, then sampling of the properties should be reset until the flow reaches steady state. As the flow reaches physically steady state, time averaging in each cell is used to obtain macroscopic properties in each cell. As a rule of thumb, sampling of particles starts when the number of molecules in the calculation domain becomes approximately constant.
2.3. General Description of the PDSC
In the past few years, we have developed a parallelized direct simulation Monte Carlo code, named PDSC, which is portable on distributed memory machines. Important features of the PDSC code can be found in the references [Wu et al.’s JCP paper, 2006; JCP paper submitted in June 2007; JFM paper under preparation, 2007] and only several of them are briefly described in the following.
1) 2D/2D-axisymmetric/3-D unstructured-grid topology: PDSC can accept either 2D/2D-axisymmetric (triangular, quadrilateral or hybrid triangular-quadrilateral) or 3D (tetrahedral, hexahedral or hybrid tetrahedral-hexahedral) mesh [Wu et al.’s JCP paper, 2006]. Computational cost of particle tracking for the unstructured mesh is generally higher than that for the structured mesh. However, the use of the unstructured mesh, which provides excellent flexibility of handling boundary conditions with complicated geometry and of parallel computing using dynamic domain decomposition based on load balancing, is highly justified.
2) Parallel computing using dynamic domain decomposition: Load balancing of PDSC is achieved by repeatedly repartitioning the computational domain using a multi-level graph-partitioning tool, PMETIS [Wu and Tseng, 2005], by taking advantage of the unstructured mesh topology employed in the code. A decision policy for repartition with a concept of Stop-At-Rise (SAR) [Wu and Tseng, 2005] or constant period of time (fixed number of time steps) can be used to decide when to repartition the domain. Capability of repartitioning of the domain at constant or variable time interval is also provided in PDSC. Resulting parallel performance is excellent if the problem size is comparably large. Details can be found in Wu and Tseng [Wu and Tseng, 2005].
3) Spatial variable time-step scheme: PDSC employs a spatial variable time-step scheme (or equivalently a variable cell-weighting scheme), based on particle flux (mass, momentum, energy) conservation when particles pass interface between cells.
This strategy can greatly reduce both the number of iterations towards the steady state, and the required number of simulated particles for an acceptable statistical uncertainty. Past experience shows this scheme is very effective when coupled with an adaptive mesh refinement technique [Wu et al.’s CPC paper, 2004].
4) Unsteady flow simulation: An unsteady sampling routine is implemented in PDSC, allowing the simulation of time-dependent flow problems in the near continuum range [JCP paper submitted in June 2007]. A post-processing procedure called DSMC Rapid Ensemble Averaging Method (DREAM) is developed to improve the statistical scatter in the results while minimizing both memory and simulation time.
In addition, a temporal variable time-step (TVTS) scheme is also developed to speed up the unsteady flow simulation using PDSC. More details can be found in [JCP paper submitted in June 2007].
5) Transient Sub-cells: Recently, transient sub-cells are implemented in PDSC directly on the unstructured grid, in which the nearest-neighbor collision can be enforced, whilst maintaining minimal computational overhead [JFM paper under preparation, 2007]. Details of the idea and implementation are described next.
2.4. Transient Sub-Cells for PDSC
[JFM paper under preparation, 2007]The implementation of sub-cells in DSMC allows the maintenance of good collision quality within the simulation, even for grids which are “under-resolved” (that is, if the sampling cells are bigger than the recommended setting of 1/3~ 1/2 times the local mean free path). Running simulations with under-resolved sampling cells which employ sub-cells results in a reduction in the computational and memory requirements of the simulation, albeit at the cost of a reduction in the possible sampling resolution of the macroscopic properties, but without sacrificing simulation accuracy.
In PDSC, unstructured grids are used, requiring an adaptation of the transient sub-cells scheme, which was originally promoted by Bird [DS2V code by Bird]. In PDSC, the sampling cells are divided into sub-cells during the collision routine. Because the sub-cells only exist in one sampling cell at a time, and only during the collision routine, they can be considered
“transient sub-cells” which will have negligible computer memory overhead. In every case, these sub-cells are quadrilateral in 2D or hexahedral in 3D which reduces the complexity of sub-dividing the sampling cell and greatly facilitates particle indexing. The size of the sub-cells is indirectly controlled by the user, who inputs the desired averaged number of
particles per sub-cell, P. The program then determines the dimensions of the sub-cell array based on the number of particles with the cell, Nparts. Briefly, the total number of sub-cells are computed by the rule for the 2-D case Fig. 2.4 shows the way in which both rectangular and triangular sampling cells are
divided into sub-cells. As can be seen, in the unstructured case, there may be sub-cells which are entirely outside the boundary of the sampling cell, however this has no affect on the collision routine. In both cases, the concept is easily extended to three-dimensional sampling cells.
During the collision routine, a particle is chosen at random from some point within the whole sampling cell. The sub-cell in which the particle lies is then determined and if another particle is in the same sub-cell then these particles are chosen for collision. If the first particle is alone within the sub-cell, then adjacent sub-cells are scanned for a collision partner.
These sub-cell routines ensure nearest neighbor collisions, even within under-resolved sampling cells, with minimal computational and memory overhead.
Bird has also shown that preventing particles from colliding again their last collision partner, reduces the error in some variables such as heat transfer and shear stress by up to 5%
[ref-Bird manual of DS2V code]. The basis of this is that collisions between particles which just collided with each other is unphysical, since the particle must be moving away from each
other after the first collision. A minor modification was made to PDSC to prevent particles colliding with their last collision partner. This involved the creation of an array in which the last collision partner for every particle is stored and if the two particles are subsequently chosen for collision without having collided with any other particle, the collision is rejected.
Chapter 3
Results and Discussion
3.1 Benchmark of Transient Sub-cell Method
Although this thesis is concerned with supersonic square driven cavity flows in rarefied regime, a subsonic flow case with M=0.5 and Kn=0.01 is used as the benchmark case. The rationale is that a subsonic flow represents a more stringent test problem than a supersonic flow for the DSMC method in terms of statistical uncertainties.
3.1.1 Problem Description and Test Conditions
Flow and simulation conditions are summarized in Table I. Fig. 3.1 sketches the 2D square driven cavity flow. PDSC simulations are conducted, respectively, using a coarse mesh (100×100 cells, 100 particles per cell) with transient sub-cells and a finer mesh (400×400 cells, 25 particles per cell) without transient sub-cells. Results using the finer mesh are considered to be a benchmark solution since the simulation conditions satisfy all requirements for a good DSMC computation. Since the overhead of indexing the particles using transient sub-cells is minimal (at most 15% as compared to that without using transient sub-cells if less than 100 particles per cell), resulting computational time using transient sub-cells is about only 1/4 of that without transient sub-cells, considering the similar accuracy of DSMC
simulations, which is shown next. If 100 particles per cell is used in the case without transient sub-cells, more benefit of computational time can gained by using transient sub-cells.
3.1.2 Simulation Results of Benchmark
Fig. 3.5 (a)-(d) (a)-(d) present the profiles of various properties ((a) u-velocity, (b) v-velocity, (c)number density and (d)temperature) along the central vertical line (x/L=0.5). In these figures, up and L represents the top plate speed and length (also width) of the driven cavity, respectively. Results show that simulation data with and without transient sub-cells are in good agreement considering the inherent statistical uncertainties of the DSMC method.
Statistical uncertainties with transient sub-cells are even lower than those without transient sub-cells, especially for the v velocity, which is generally low as compared to the u velocity along the x/L=0.5 line. It clearly demonstrates that with transient sub-cells accuracy of DSMC simulation can be well preserved using cell size even close to local mean path. Similar trends can be found in Fig. 3.6 (a)-(d), in which the profiles of various properties are presented along the central horizontal line (y/L=0.5). In the simulation, desired averaged number of particles per sub-cell, P, is set to be 1. Not unexpectedly, the merit of collision (=mean collision distance/local mean free path), which represents the quality of collisions in a DSMC simulation, is generally much smaller by using transient sub-cells than that without using transient sub-cells, as illustrated in Fig. 3.11 and Fig. 3.12. Bird [???] has argued that the
merit of collision in the range of 0.1~0.2 represents a good DSMC simulation. Result show that the simulation used coarse mesh with transient sub-cell approximate that used finer mesh without transient sub-cell. Therefore, the simulation with transient sub-cell keeping minimal computational overhead and memory requirement simultaneously, but still has simulation accuracy.
3.2 Simulation Results of Driven Cavity Flows
According to previous section, we have found the transient sub-cell method can greatly reduce the computational cost and still has good quality of solution. In this section, we have used this technique in PDSC to simulate steady supersonic square driven cavity flows in a systematic way. Thus, we can understand more about the effects of rarefaction and compressibility in such conditions.
3.2.1 Problem Description and Test Conditions
Fig. 3.1 shows the sketch of the 2D square (L/H=1) driven cavity flow with moving top plate. Flow and simulation conditions are summarized in Table II and Table III. Argon is used as the working gas, while the initial temperatures inside the cavity and wall temperatures
Fig. 3.1 shows the sketch of the 2D square (L/H=1) driven cavity flow with moving top plate. Flow and simulation conditions are summarized in Table II and Table III. Argon is used as the working gas, while the initial temperatures inside the cavity and wall temperatures