## Development of a Hybrid Modular Particle-Continuum (MPC) Method for Nonequilibrium, Hypersonic Flows

### Motivation

Computation of the aerothermodynamics of hypersonic re-entry vehicles along their entire trajectory involves near-continuum conditions at low altitudes and rarefied, or non-equilibrium, conditions at high altitudes. Well-established simulation methods already exist for each of these flow regimes. For example, the continuum Navier-Stokes (NS) equations may be solved using algorithms from Computational Fluid Dynamics (CFD) and the particle-based, direct simulation Monte Carlo method (DSMC) may be used for the non-equilibrium, or rarefied, flows. A hybrid method that blends the CFD and DSMC techniques is an attractive approach for simulation of flows involving a mixture of both continuum and non-equilibrium flow regimes. This is the situation at intermediate altitudes where within a mostly continuum flow, there may be local regions of non-equilibrium flow generated by both the rapid expansion behid a re-entering capsule as well as by strong gradients in shock waves and boundary layers.

For dilute gases, such as the earth's atmosphere, the most popular numerical method for simulating non-equilibrium flow is the DSMC method first developed by Graeme Bird. Here, the trajectories of a large number of simulated particles are followed simultaneously through a grid of computational cells. For each iteration, particles are first moved along their trajectories without colliding, after which particles residing within the same cell are randomly selected for a collision process. A major limitation of DSMC is that in order for this collision process to be physically accurate, the cell size must be on the order of the mean free path, λ, while containing at least 20 simulated particles per cell. As a result, 2D and 3D DSMC simulations can require prohibitively high numbers of computational cells and therefore simulated particles, especially in regions where λ is very small. This can be seen in variation in mean free path around a planetary probe (seen below) where the characteristic flow length scale can be very small in the dense fore body region of the probe. It is precisely in these regions, where λ is small, that the continuum approximation is valid and the NS equations can be solved without the same restriction on cell size. This is the motivation behind developing a hybrid particle-continuum numerical scheme for hypersonic non-equilibrium flows.

**Basic Hybrid Methodology
**

A hybrid CFD-DSMC approach, developed at the
University of Michigan, is called the Modular Particle Continuum (MPC)
method and merges existing CFD and DSMC codes. Both codes remain
virtually unmodified and hybrid data structures and functions are
constructed to track macroscopic quantities, passs information between
solvers, and perform other required hybrid processes. The accuracy of
the MPC method relies on the proper positioning of the DSMC-CFD
interfaces. The interface must lie in near equilibrium regions where
solution of the NS equations will introduce minimal or no error.
Typically, particle and continuum regions are determined by evaluating
a continuum breakdown parameter in the flow field. The MPC method uses
the gradient-length Knudsen number, Kn_{GL-Q} =λ |∇Q/Q|,
where Q represents local flow quantities of interest such as density,
temperature, or velocity magnitude, λ is the local mean-free-path. It
has been shown for hypersonic flows that, in regions of the flow field
where Kn_{GL-MAX}
<0.05, the discrepancy between a NS and DSMC solution is less than
5%. Thus, these regions could be solved using a NS solver with little
error.

The MPC method begins with a NS solution of
the entire flow field and then uses gradient-length Knudsen number to
decompose the domain into CFD and DSMC regions. The method uses
state-based coupling to transfer information between particle and
continuum regions. After evaluation of the breakdown parameter, the
particle region is extended by a few extra cells into the continuum
domain to create an overlap region as seen in the figure below. Next,
one row of NS and two rows of DSMC boundary cells are initialized. The
DSMC and NS domains are then coupled by transferring information across
the interfaces. With state-based coupling, this involves updating the
boundary conditions of each solver. In this way, information transfer
into both the particle and continuum regions is handled through
existing boundary procedures already used by both solvers. The DSMC
boundary cells are continually filled with particles consistent with
the flow properties in the corresponding NS cell using the
Chapman-Enskog distribution based on the local macroscopic state and
gradients, known from the NS solver. As particles in the DSMC domain
interact, their distributions evolve in time towards the distribution
predicted by full DSMC. The MPC method also tracks the
macroscopic variation in each DSMC cell. In order to provide these
averaged properties with minimal statistical scatter, a mixture of
spatial and temporal averaging is used. Specifically, MPC uses the
subrelaxation technique proposed by Sun et al. These averaged DSMC
properties are then used to update boundary conditions for the NS
solver.

The DSMC-CFD interface locations for Mach 20
flow over a sting-mounted planetary probe with a global Knudsen number
of 0.01 are shown below along with the variation in flow length
scales. Due to the overlap regions, the interface locations have
moved from the initial position that was calculated using a full
Navier-Stokes solution to the final interface location based on the
physically correct, DSMC solution. The MPC method is able to
simulate steady-state, hypersonic flow that is in thermal (rotational
and vibrational) nonequilibrium. It is parallelized for
distributed memory systems using dynamic domain decomposition with
scaled efficiency remaining above 80% for typical near-continuum 2D and
axi-symmetric flows.

**Recent Results**

The following
figures compare translational and rotational temperature contours
predicted by full DSMC, full CFD, and the MPC method for Mach 20 flow
over a sting-mounted, planetary probe with a global Knudsen number of
0.01. In general the MPC method has improved agreement with full DSMC
results from the initial CFD solution across the entire flow field.
Even in regions where the CFD module is used, the MPC method has
obtained improved boundary conditions from DSMC, which has allowed the
CFD module to shift its result to excellent agreement with full DSMC.

The MPC method is also in near exact agreement
with full DSMC for surface properties. The following figure shows
surface heat transfer to the planetary probe predicted by full DSMC,
full CFD, and the MPC method compared with available experimental
measurements. Along the fore body where the flow is highly
collisional, all three methods are in good agreement with each other
and the experimental measurements. Despite this highly collisional
flow, CFD still slightly over predicts DSMC and MPC results. This is
due to the inability to correctly model the Knudsen layer within the
CFD solver. As the flow expands around the corner, full CFD over
predicts both DSMC and experimental measurements by over an order of
magnitude. In contrast, the MPC method remains in good agreement with
both the experimental measurements and the full DSMC results.
Similarly, the MPC method remains in excellent agreement with full DSMC
and experimental measurements along the sting mount, while full CFD
over predicts full DSMC along the entire sting mount.

In addition to reproducing full DSMC results, an objective of the MPC method is to reproduce DSMC at a reduced compuational cost compared to full DSMC. All flows considered have experienced a speedup compared to full DSMC. Speedup factors have ranged from 1.4 times faster for mostly rarefied flow to up to 30 times faster than DSMC for 2D and axi-symmetric, near-continuum flows, while still maintaining excellent agreement with full DSMC results. It is expected that current work of extending the method to simulate full 3 dimensional flow will experience even larger differentials for compuational cost between the MPC method and full DSMC.

### Recent Publications

### Acknowledgements

Financial support for this work is provided by the Constellation University Institute Project (CUIP), under NASA grant NCC3-989.