May 2012 

NASA/TM-2012-217575 

NESC-RP-07-037 

 

Reacting Multi-Species Gas Capability for 

USM3D Flow Solver 

Neal T. Frink and David M. Schuster/NESC 

Langley Research Center, Hampton, Virginia 
NASA STI Program . . . in Profile

 Since its founding, NASA has been dedicated to 

the advancement of aeronautics and space science. 

The NASA scientific and technical information (STI) 

program plays a key part in helping NASA maintain 

this important role. 

 The NASA STI program operates under the 

auspices of the Agency Chief Information Officer. It 

collects, organizes, provides for archiving, and 

disseminates NASA’s STI. The NASA STI program 

provides access to the NASA Aeronautics and Space 

Database and its public interface, the NASA Technical 

Report Server, thus providing one of the largest 

collections of aeronautical and space science STI in 

the world. Results are published in both non-NASA 

channels and by NASA in the NASA STI Report 

Series, which includes the following report types: 

• TECHNICAL PUBLICATION. Reports of 

completed research or a major significant phase 

of research that present the results of NASA 

programs and include extensive data or 

theoretical analysis. Includes compilations of 

significant scientific and technical data and 

information deemed to be of continuing 

reference value. NASA counterpart of peerreviewed formal professional papers, but having 

less stringent limitations on manuscript length 

and extent of graphic presentations. 

• TECHNICAL MEMORANDUM. Scientific 

and technical findings that are preliminary or of 

specialized interest, e.g., quick release reports, 

working papers, and bibliographies that contain 

minimal annotation. Does not contain extensive 

analysis. 

• CONTRACTOR REPORT. Scientific and 

technical findings by NASA-sponsored 

contractors and grantees. 

• CONFERENCE PUBLICATION. Collected 

papers from scientific and technical 

conferences, symposia, seminars, or other 

meetings sponsored or co-sponsored by NASA. 

• SPECIAL PUBLICATION. Scientific, 

technical, or historical information from NASA 

programs, projects, and missions, often 

concerned with subjects having substantial 

public interest. 

• TECHNICAL TRANSLATION. Englishlanguage translations of foreign scientific and 

technical material pertinent to NASA’s mission. 

 Specialized services also include creating custom 

thesauri, building customized databases, and 

organizing and publishing research results. 

 For more information about the NASA STI 

program, see the following: 

• Access the NASA STI program home page at 

http://www.sti.nasa.gov

• E-mail your question via the Internet to 

help@sti.nasa.gov

• Fax your question to the NASA STI Help Desk 

at 443-757-5803 

• Phone the NASA STI Help Desk at 

443-757-5802 

• Write to: 

 NASA STI Help Desk 

 NASA Center for AeroSpace Information 

 7115 Standard Drive 

 Hanover, MD 21076-1320 
National Aeronautics and 

Space Administration 

Langley Research Center 

Hampton, Virginia 23681-2199 

 

May 2012 

NASA/TM-2012-217575 

NESC-RP-07-037 

 

Reacting Multi-Species Gas Capability for 

USM3D Flow Solver 

Neal T. Frink and David M. Schuster/NESC 

Langley Research Center, Hampton, Virginia 
Available from: 

NASA Center for AeroSpace Information 

7115 Standard Drive 

Hanover, MD 21076-1320 

443-757-5802 

The use of trademarks or names of manufacturers in the report is for accurate reporting and does not 

constitute an official endorsement, either expressed or implied, of such products or manufacturers by the 

National Aeronautics and Space Administration. 
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

1 of 78

NESC Request No.: 07-037-I

Reacting Multi-Species Gas Capability for USM3D Flow Solver

March 8, 2012
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

2 of 78

NESC Request No.: 07-037-I

Approval and Document Revision History 

NOTE: This document was approved at the March 8, 2012, NRB. This document was 

submitted to the NESC Director on April 24, 2012, for configuration control.

Approved: Original Signature on File 4/25/12

NESC Director Date

Version Description of Revision Office of Primary 

Responsibility

Effective Date

1.0 Initial Release Dr. Neal T. Frink, Senior 

Research Engineer, Langley 

Research Center

3/8/12
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

3 of 78

NESC Request No.: 07-037-I

Table of Contents

Technical Assessment Report

1.0 Notification and Authorization.........................................................................................5

2.0 Signature Page....................................................................................................................6

3.0 Team List............................................................................................................................7

4.0 Executive Summary ...........................................................................................................8

5.0 Assessment Plan .................................................................................................................9

6.0 Problem Description ..........................................................................................................9

7.0 Technical Approach.........................................................................................................10

7.1 Theoretical Formulation........................................................................................... 10

7.2 Numerical Implementation ...................................................................................... 22

8.0 Verification of Multi-Species Capability .......................................................................31

8.1 Comparison of Ideal and Real Gas Simulations...................................................... 32

8.2 Comparison of Flux Schemes with N2 Exhaust Gas................................................ 35

8.3 Results from Simulations using Hydrogen (H2) as the Exhaust Gas....................... 39

9.0 Findings, Observations, and NESC Recommendations...............................................42

9.1 Findings.................................................................................................................... 42

9.2 Observations ............................................................................................................ 43

9.3 NESC Recommendations......................................................................................... 43

10.0 Alternate Viewpoints.......................................................................................................44

11.0 Other Deliverables...........................................................................................................44

12.0 Lessons Learned...............................................................................................................44

13.0 Definition of Terms..........................................................................................................44

14.0 Acronyms List ..................................................................................................................45

14.1 Nomenclature........................................................................................................... 45

15.0 References.........................................................................................................................47

16.0 Appendices........................................................................................................................47
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

4 of 78

NESC Request No.: 07-037-I

List of Figures

Figure 7.1-1. Nomenclature of a Face and it's Neighbors .......................................................... 20

Figure 8.0-1. Triangulation on Boundary Plane of Tetrahedral Quadrant Grid for Nozzle Test 

Case....................................................................................................................... 32

Figure 8.1-1. Comparison of Results from Ideal and Real Air Simulations............................... 34

Figure 8.2-1. Results from the Simulations using N2 as the Exhaust Gas.................................. 39

Figure 8.3-1. Results from the Simulations using H2 as the Exhaust Gas.................................. 42
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

5 of 78

NESC Request No.: 07-037-I

Technical Assessment Report

1.0 Notification and Authorization 

The NASA Engineering and Safety Center (NESC) approved this assessment on 

June 14, 2007. Dr. Neal Frink, Senior Research Engineer at Langley Research Center (LaRC),

was selected to lead this assessment. The Phase I Assessment Plan was approved by the NESC 

Review Board (NRB) on July 12, 2007. The Phase II Assessment Plan was approved by the 

NRB on December 16, 2008.
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

6 of 78

NESC Request No.: 07-037-I

2.0 Signature Page 

Submitted by:

Team Signature Page on File – 5/4/12

Dr. Neal T. Frink Date Dr. David M. Schuster Date

Significant Contributors:

Dr. Gary Cheng Date Dr. Roy P. Koomullil Date

 

 

Signatories declare the findings, observations, and NESC recommendations compiled in the 

report are factually based from data extracted from Program/Project documents, contractor 

reports, and open literature, and/or generated from independently conducted tests, analyses, and 

inspections. 
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

7 of 78

NESC Request No.: 07-037-I

3.0 Team List

Name Discipline Organization/Location

Core Team

Neal Frink

Team Lead, Senior Research 

Engineer LaRC

Dave Schuster

NASA Technical Fellow for 

Aerosciences LaRC

Mohagna Pandya USM3D Expert LaRC

Gary Cheng Investigator

University of Alabama at 

Birmingham/Mechanical Engineering

Roy Koomullil Investigator

University of Alabama at 

Birmingham/Mechanical Engineering

Patricia Pahlavani MTSO Program Analyst LaRC

Administrative Support

Christina Williams Technical Writer ATK, LaRC
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

8 of 78

NESC Request No.: 07-037-I

4.0 Executive Summary

The USM3D Navier-Stokes flow solver contributed heavily to the NASA Constellation Project

(CxP) as a highly productive computational tool for generating the aerodynamic databases for 

the Ares I and V launch vehicles and Orion launch abort vehicle (LAV). At the peak of that 

effort, USM3D was consuming approximately 25 percent of the total available computer

resources on the 10,000-processor NASA Advanced Supercomputing Columbia machine at the 

NASA Ames Research Center for these projects.

Some of the largest database uncertainties were associated with inadequate aerodynamic 

modeling of jet impingement. For example, the aerodynamic modeling uncertainties from the 

abort motor plume impingement on the LAV flight vehicle continually contributed to failure of 

the flight simulations, which elevated the risk to flight safety for launch abort. Similar 

uncertainties were present in the jet impingement from the Ares roll control system (RCS). 

There is a need to improve this deficiency through implementation of more physically correct 

models of the hot-gas jet flows. Efforts to resolve these uncertainties within any of the flow 

solvers used in the CxP have not been successful to date.

USM3D is currently limited to ideal-gas flows, which are not adequate for modeling the 

chemistry or temperature effects of hot-gas jet flows. This task was initiated to create an 

efficient implementation of multi-species gas and equilibrium chemistry into the USM3D code to 

improve its predictive capabilities for hot jet impingement effects. A finite-rate chemistry model 

was also pursued to improve future predictions of re-entry heating rates of planetary entry 

vehicles. 

The effort was only partially successful as the investigators struggled with converting an existing

flow solver with tightly embedded ideal-gas assumptions to a real-gas chemistry code. At this 

time, the implementation of the real-gas and multi-species transport equations has been verified

on simple unit problems, but has not been validated against experimental data. Also, some

numerical solution stability issues are yet to be resolved. The modules for the equilibrium and 

reacting chemistry were developed, but were not implemented because of the prior difficulties.

Work is continuing through the support of the Space Launch System (SLS) Aero Task Lead, 

Dr. John Blevins at NASA Marshall Space Flight Center, using the services of subject matter 

experts, including: Dr. Michael Applebaum, Dr. William Eppard, and Mr. Les Hall at CRM 

Solutions, Inc., in Huntsville, Alabama. Assessments are underway to determine the best path 

forward. The near-term goal is to resolve the numerical stability issues with the multi-species 

capability as soon as possible and commence application to Multi-Purpose Crew Vehicle 

(MPCV) launch abort and SLS RCS jet cases. The long-term goal is to implement and test the 

equilibrium and reacting-gas capability.
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

9 of 78

NESC Request No.: 07-037-I

5.0 Assessment Plan

The goal of this assessment was to implement and validate a simulation capability to handle realgas effects in the USM3D code. This implementation was proposed as a three-step process: 

1. Modification of governing equations to handle variable ratio of specific heats. 

2. Addition of species transport equations to model the transport and mixing of multiple 

species in the computational domain. 

3. Addition of capability for different species to react with each other and produce new 

species. 

A detailed description of the USM3D code modifications, addition of new capabilities, and 

implementation model validation are described in the following sections.

6.0 Problem Description

The underlying physics associated with flow around a supersonic/hypersonic vehicle is complex, 

which involves an oblique/bow shock, air dissociation (or ionization), and large heat transfer. 

Accurate predictions of the aerothermodynamics of supersonic/hypersonic vehicles require 

efficient numerical models for resolving chemical reactions and modeling of transport of 

multiple species, and accurate estimation of heat transfer. The modeling of chemical reactions 

and species transport is also important in the numerical simulation of combustion. Challenges 

involved in the simulation of flows with chemical reactions are the stability and efficiency issues 

associated with the numerical methods. The stability of the numerical scheme for the species 

equation is hindered by disparity between the time scales for the mean flow and chemical 

reactions and the stiffness problem of the source terms. Various methods were developed to 

alleviate the stiffness problem and to efficiently solve a large set of species transport equations 

with momentum and energy equations. Since each method has its strength and weakness, the 

USM3D code needed to be analyzed with the assistance of the code developer to determine the 

numerical method most appropriate to efficiently and accurately solve the set of governing 

equations. This report describes the addition of the chemically-reacting flow capability into the 

USM3D code.
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

10 of 78

NESC Request No.: 07-037-I

7.0 Technical Approach

Section 14.1 provides the nomenclature definitions for equations used in the following sections.

7.1 Theoretical Formulation

Governing Equations

The governing equations for single-species real-gas flow can be expressed as:

Continuity equation:

0

ˆ ˆˆ

ˆ

ˆ

 











 





V

t







(1)

Momentum equation:

 

 





 



ˆ

ˆ ˆˆ ˆ ˆ ˆˆ

ˆ

ˆ

    











 





VV p

t

V

(2)

Total energy equation:

E V E p  T  V Q

t



  



ˆ ˆ ˆˆ ˆ ˆ ˆˆ ˆ)ˆ

( ˆ

ˆ ˆ ˆ ˆ

ˆ

 











     









  





   

(3)

where 

 

Pr

ˆ

ˆ

ˆ ˆ ˆ ˆ , and ˆ (thermalconductivity) ˆ

2 2 2

2

1





Cp

E  e u v w  .

The set of governing equations is non-dimensionalized with a set of reference quantities as:

(

ref ref ref ref ref ref ref  ,  , a , T , R ,  , and 

), such as:

  ref

3

ref ref ref ref

ref ref ref ref

ref ref v

ref

ref

v

v

ref ref ref

2

ref

2

2

ref

2

ref

2

ref

ref

2

ref ref ref ref ref ref

ref

/

ˆ

;

ˆ

; ;

1

1

ˆ

;

ˆ

;

ˆ

;

ˆ

;

2

ˆ

ˆ

ˆ

;

ˆ

ˆ

;

ˆ

;

ˆ

;

ˆ

;

ˆ

;

ˆ

















a

Q

Q

R

a R T

C R C

C

R

C

C

R

C

C

R

R

R

a

V

a

e

a

E

E

a

e

e

T

T

T

a

p

p

a

t

t

a

V

V

p

p

p

p

 



  







 





  



 

















       

       

 (4)

The reference length scale (

ref 

) was read into the USM3D code through input based on the unit 

of the mesh coordinates. If the air at the sea level (1 atm and 300K) was taken as reference state, 

then the reference temperature (

Tref

), pressure (

pref

), species mass fractions (

s,ref

), the 

reference gas constant (

Rref

), reference density (

 ref

), reference specific heats ratio (

ref 

), 

reference speed of sound (

aref

), and reference viscosity (

 ref

) were determined automatically. 

The values of all flow variables and mesh coordinates were non-dimensionalized based on this 

set of reference quantities once they were read into the code. The solution of the flowfield was
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

11 of 78

NESC Request No.: 07-037-I

calculated based on the non-dimensionalized governing equation. The Sutherland law for air 

was used as:

T S

T S

T

T























ˆ

ˆ

ˆ

ˆ

ˆ

ˆ

0

1.5

0 0







ref

0 ref

1.5

ref 0

0 0

1.5

ref 0

0

ref

0

ref 0

/

ˆ /

ˆ

ˆ

ˆ

ˆ ˆ ˆ

ˆ

ˆ ˆ

T S T

T S T

T

T

T S

T S

T

T











































   























(5)

If 

0 ref T

ˆ

 T

, then 

T01

The set of non-dimensionalized governing equations was obtained as:

    0





V

t







(6)

 

  





 



     





ref

ref

Re

M

VV p

t

V

(7)

     V  Q

M

T

M

E V E p

t



  



Re Re

( )

ref

ref

ref

ref         





   

(8)

where the non-dimensional  is defined as 

 





































ref Pr 1 









R

(9)

For the ease in discretization of these equations for unstructured meshes, the governing equations 

were written in the integral form as:

Qd F nd S F nd S

t

S

v

S

i



 

 

   





  

V

V

(10)

where



































E

w

v

u

Q













(11)
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

12 of 78

NESC Request No.: 07-037-I

  































 







 

t

z

y

x

i

E p pn

w pn

v pn

u pn

F n

 

 

 

 





(12)

 

 

 





































































  

 

 

 



































































 

z

T

n

y

T

n

x

T

u f v f wf n

n n n

n n n

n n n

M

f

f

f

f

f

M

F n

x y z

v v v

x z x y z y z z z

x y x y y y z y z

x x x y x y z x z

v

v

v

v

v

v



  

  

  

2 3 4

ref

ref

5

4

3

2

1

ref

ref

0

R e R e





(13)

where  is the contravariant velocity component of the fluid particle with respect to the grid and 

is defined as:

x y

wnznt un vn  

The thermal equation of state for perfect a gas is:

R

p

T



 ref



(14)

Relations between enthalpy per mole, specific heat per mole at constant pressure, and 

temperature are:

8

5

5 7

4 1

4 6

3 1

3 5

2 1

2 4

1

2 3

1

1

o

4

7

3

6

2

3 4 5

1

2

2

1

o

ˆ ˆ ˆ ˆ ˆ ˆ

l n ˆ ˆ

ˆ ˆ ˆ ˆ ˆ ˆ ˆ

h R a T a T a T a T a T a T a T a

C R a T a T a a T a T a T a T

u

p u

        

      



 

(15)

Species Continuity Equations

One of the prerequisites for the chemically-reacting flow modeling is the capability to solve 

species transport equations. This is a result of chemical reactions, where species mass fractions 

can be changed and transported through convection and diffusion processes. The local mass 

fraction of different species with the temperature and pressure govern the local reaction. The 

change of species mass fractions alters the mixture temperature and density at a given pressure 

and enthalpy. Therefore, an accurate modeling of the species transport equation is an important 

step. Appropriate equations of state (i.e., thermal and caloric) need to be employed to account 

for the real-gas effect involving multiple species. 
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

13 of 78

NESC Request No.: 07-037-I

The species transport equation can be written as:

    i i i i

i

V D

t

  







   





(16)

where 

V



 and

are the mean density and velocity vector, 

Diii, , and 

are diffusion coefficient, 

species mass fraction, and source term of i-th species, respectively. 

Using an eddy viscosity model, the diffusion coefficient can be written as:

i

c t

t

t

S

D 



 

,

(17)

Therefore, the species transport equation can be written as:

  i i

c t

t

c

i

i

S S

M

V

t

 

 

 









































    





Reref ,

 (18)

Integrating over the control volume and applying the divergence theorem, Equation 18 reduces 

to:

   

   

  

















   





ndA dV

S S

M

dV V ndA

t

i i

c t

t

c

i

i

 

 







 



Reref ,

 (19)

where i = 1, 2, ..., NS and NS is the number of species

For modeling chemical reactions, the species source term can be evaluated based on local species 

mass fractions and thermal state using either equilibrium or finite-rate chemistry. The species 

transport equation can be solved either independently (i.e., segregated approach) or by coupling 

with other governing equations (i.e., coupled approach). For the segregated approach, the 

species equation can be solved with the same time-step size as other governing equations, or with 

a smaller time-step size than other equations due to the disparate time-scale between flow and 

chemical reactions. The segregated approach requires sub-iterations to ensure convergence. The 

coupled approach does not need sub-iterations, but is limited to use the same time-step size for 

all equations, and thus requires large computer memory to store the Jacobian matrix, especially 

when there are a large number of species involved.
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

14 of 78

NESC Request No.: 07-037-I

For the source term of the species transport equation, it can be solved using an explicit or 

implicit approach. Each of these methods has advantages and disadvantages. In the case of 

explicit methods, the source term will be treated explicitly so that it is computationally efficient. 

However, the stability and accuracy are sacrificed due to the stiffness of the source term. In the 

case of implicit schemes, the sources terms will be coupled together and treated implicitly. This 

increases the slow solver stability, but the computational efficiency will be compromised 

because of the need for the solution of the denser matrix system. 

Equilibrium Chemistry Model

The general form of chemical reactions can be expressed as:

   

  

NS

i

ji i

b j

f j

NS

i

ji i

K

K

1 ,

,

1

 

(20)

where 

 ji  ji

 and 

are the stoichiometric coefficients of species “i” in forward and backward 

reactions “j”, 

Kf , jKb, jand

are the forward and backward reaction rates of the j

th reaction, Φi is 

the name of i-th species, and NS is the total number of species.

The equilibrium chemistry model assumes that the forward and backward reactions are in the 

equilibrium state. In this case, the species mass fractions are determined based on the local 

thermal state (i.e., pressure and mixture energy) and element balance (i.e., chemical elements 

cannot be created or destroyed). Based on these principles, a Gibbs energy minimization method 

can be devised to evaluate the equilibrium species mass fraction (

n1 i

). There are two ways to 

model the equilibrium chemistry with species transport. The first approach is to assume that 

reaction takes place instantly (i.e., infinitesimal time), where the species transport equation is 

solved with zero source term. After the species equation is solved, the equilibrium chemistry is 

used to determine the new species mass fraction based on the local mass fractions and thermal 

state. The second approach is to treat the convection, diffusion, and reaction processes that occur 

simultaneously, and evaluate the source term as:

t

n

i

n

i

i









( )

1

  

 (21)

The chemical equilibrium constant can be used to calculate the reverse rate constant if the 

forward rate constant is known.

Limitations of the equilibrium chemistry include: 

1. Equilibrium chemistry assumption is valid only at high temperature and high pressure 

when the reaction is extremely fast, 
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

15 of 78

NESC Request No.: 07-037-I

2. Fuel combustion (e.g., heavy hydrocarbon fuel) does not follow the equilibrium 

chemistry, especially at the fuel-rich condition, 

3. The occurrence of soot and/or nitrogen oxide is controlled by the kinetic chemistry, and 

is not in chemical equilibrium, 

4. Ignition and flame extinctions are not equilibrium processes, and are dependent on 

kinetic rates, and 

5. The equilibrium chemistry cannot be used to determine the flame front location.

Finite-rate Chemistry Model

For a finite-rate chemistry model, which assumes all chemical reactions occur at a finite rate, the 

source term in Equation 19 is calculated based on the reaction processes, their associated kinetic 

rates, local species mass fractions, pressure, temperature, and the presence of a catalyst or 

inhibitor. The source term due to chemical reaction can be expressed as:

























 























 









 

j l j l r

l

l

N S

l

b j

N R

j

r

l

l

N S

l

i i ji f j

M

K

M

M K

 

 

1

,

1

1

,

 (22)

where 

Mi

is the molecular weight of species i, and NR is the total number of reactions. The 

chemical reaction rates in Equation 22 controls the rate of species formation and destruction. For 

example, higher temperature and pressure will result in a faster reaction rate. The reaction rate 

formulated in the empirical Arrhenius form is:

K AT  E R T a u

B

f

 exp  (23)

where A, B, Ea are empirical constants. As described earlier, there are two ways of calculating 

the species source term. For the explicit method, since it is numerically unstable and inaccurate, 

some researchers employed the penalty function method to minimize the numerical error and 

improve the numerical stability. 

In the method, the time-step size of the species equation is determined by:

















 i a

i

i

t ( )

max

1

 



(24)
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

16 of 78

NESC Request No.: 07-037-I

where 

a

()

is the specified allowable maximum species mass fraction change. However, for a 

stiff reaction model, the calculated time-step size can become small, and thus requires a large 

number of sub-iterations. For the implicit method, since all the species involved in the source 

term are coupled and treated implicitly, the Taylor series expansion is employed to linearize the 

source term. Because the solution of a matrix is involved in the implicit scheme, it is less 

computationally efficient, but is more stable and accurate.

Linearization of the Inviscid Flux

Representing the conserved variable vector Q as:





































































5

4

3

2

1

Q

Q

Q

Q

Q

E

w

v

u

Q













(25)

The contravariant velocity can be represented as:

  x y z

Q n Q n Q n

Q

2 3 4

1

1

   

(26)

Using Equation 26, the inviscid flux vector can be written as:

 

 

 

   





















































      

    

    

    

   



































 

2 3 4

1

1 2 3 4

1

5

1 2 3 4

1

4

1 2 3 4

1

3

1 2 3 4

1

2

1 2 3 4

5

4

3

2

1

n Q n Q n Q

Q

p

n Q n Q n Q n Q

Q

Q

n Q n Q n Q n Q n p

Q

Q

n Q n Q n Q n Q n p

Q

Q

n Q n Q n Q n Q n p

Q

Q

n Q n Q n Q n Q

f

f

f

f

f

F n

t x y z x y z

t x y z z

t x y z y

t x y z x

t x y z

i

i

i

i

i

i





(27)

Equation 27 can be differentiated with respect to Equation 25 to get the inviscid flux Jacobian. 
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

17 of 78

NESC Request No.: 07-037-I

Each component of the inviscid flux vector can be represented as 



Ajk 

 f

j

i

Qk

, which will result 

in:



































5 1 5 2 5 3 5 4 5 5

4 1 4 2 4 3 4 4 4 5

3 1 3 2 3 3 3 4 3 5

2 1 2 2 2 3 2 4 2 5

1 1 1 2 1 3 1 4 1 5

A A A A A

A A A A A

A A A A A

A A A A A

A A A A A

A

(28)

For the linearization of the flux vector, it is assumed that the ratio of specific heats is locally 

constant. Based on this assumption, each term in Equation 28 can be expressed as:



A 

ntnxnynz0

nx  u nx

2   u U ny

u  nx

 1v nz

u  nx

 1w nx 1

ny  v nxv  ny

 1u ny2   v U nz

v  ny

 1w ny 1

nz  w nxw  nz

 1u ny

w  nz

 1v nz2   w U nz 1

2 

E



 

   

 

    nx

E



 

 

   

 

     1u ny

E



 

 

   

 

     1v nz

E



 

 

   

 

     1w   nt

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

(29)

where

x y

wnz  un  vn 

(30)

 

2 2 2

2

1

u  v  w









(31)

Linearization of Viscous Terms

The viscous terms involves the gradients of the primitive variables. Therefore, linearization of 

the viscous terms will be simpler if primitive variables are used as the independent variable. To 

linearize the viscous flux with respect to the conserved variables, the Jacobian, based on the 

primitive variables, can be multiplied by the linearization of the primitive variables with respect 

to the conserved variables. 
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

18 of 78

NESC Request No.: 07-037-I

The variable vector q selected for this linearization procedure is taken as:







































































































5

4 1

3 1

2 1

1

5

4

3

2

1

q

Q Q

Q Q

Q Q

Q

q

q

q

q

q

p

w

v

u

q





(32)

Hence, Jacobian resulting from the linearization of the primitive variable with respect to the 

conserved variable can be written as:



































































































































































5

5

4

5

3

5

2

5

1

5

5

4

4

4

3

4

2

4

1

4

5

3

4

3

3

3

2

3

1

3

5

2

4

2

3

2

2

2

1

2

5

1

4

1

3

1

2

1

1

1

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q





(33)

Using the relation between primitive and conserve variable, the above Jacobian can be reduced 

to:





































      

















( 1) ( 1) ( 1) 1

/ 0 0 1/ 0

/ 0 1/ 0 0

/ 1/ 0 0 0

1 0 0 0 0

   



 

 

 

u v w

p

w

v

u

Q

q





(34)

The details of this derivation are provided in Appendix A. The linearization of pressure for realgas is provided in Appendix B.
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

19 of 78

NESC Request No.: 07-037-I

From Equation 14, Equation 35 can be obtained as:

d p

R

d

T

dT









ref

  

(35)

Hence,



















R

T

q

T







ref

 0 0 0

(36)

Combining Equations 34 and 36 produces: 

 











        

















       















































      









































































R R

w

R

v

R

V u

R

T

R R

w

R

v

R

p u

R

T

u v w

p

w

v

u

R

T

Q

q

q

T

Q

T

Q

T

Q

T

Q

T

Q

T

Q

T



 



 



 



  









 



 



 



 

 





   



 

 

 







( 1) ( 1) ( 1) ( 1)

2

1

 

( 1) ( 1) ( 1) ( 1)

( 1) ( 1) ( 1) 1

/ 0 0 1/ 0

/ 0 1/ 0 0

/ 1/ 0 0 0

1 0 0 0 0

0 0 0

ref ref ref ref

2

ref

ref ref ref ref ref

ref

1 2 3 4 5





 

(37)

The spatial derivatives any primitive variable  at the interface between two adjacent cells can be 

calculated as:

n x n y n z

n

z

n x n y n z

n

y

n x n y n z

n

x

x y z

z

z z z

x y z

y

y y y

x y z

x

x x x

    

    







    

    







    

    







;

;

;



















 (38)
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

20 of 78

NESC Request No.: 07-037-I

and 

 RL

where  Rand  Lare the values of  at the center of right and left cells, as 

shown in Figure 7.1-1. 

Figure 7.1-1. Nomenclature of a Face and it's Neighbors

Hence, the spatial derivatives of a dependent variable at the cell face can be linearized by the 

conservative as:

 



















 

 



 

 





Q Q Q

R L

x x

x

  

 ( )  

(39)

Q


NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

21 of 78

NESC Request No.: 07-037-I

Therefore,





















 





R

R

x

R

x

Q Q

 

 

(40)





















 





L

L

x

L

x

Q Q

 

 

(41)

Equation 13 can be differentiated with respect to Equation 25 to get the viscous flux Jacobian. 

Each component of the viscous flux vector can be represented as 

k

v

j

jk Q

f

B







and this will result 

in:



































5 1 5 2 5 3 5 4 5 5

4 1 4 2 4 3 4 4 4 5

3 1 3 2 3 3 3 4 3 5

2 1 2 2 2 3 2 4 2 5

1 1 1 2 1 3 1 4 1 5

B B B B B

B B B B B

B B B B B

B B B B B

B B B B B

B

(42)

By substituting Equations 34, 37, and 39, each components of the viscous flux Jacobian can be 

expressed as:

B11  B12  B13  B14  B15  0

(43)

        







































c j

z x

z

c j

y x

y

c j

x x

x

c j

v

j

Q

n

Q

n

Q

n

Q

f

B

   2

2

(44)

       









































C j

z y

z

C j

y y

y

C j

x y

x

C j

v

j

Q

n

Q

n

Q

n

Q

f

B

  

3

3

(45)

       









































C j

z z

z

C j

y z

y

C j

x z

x

C j

v

j

Q

n

Q

n

Q

n

Q

f

B

   4

4

(46)

 

 

 

 

 

 

       

      































   



























































C j

z

z

C j

y

y

C j

x

f j f j f j x

C j

z

z

C j

y

y

C j

x

x

C j

v

C j

v

C j

v

f

C j

v

j

Q

T

n

Q

T

n

Q

T

u B v B w B n

Q

T

n

Q

T

n

Q

T

n

Q

w f

Q

v f

Q

u f

Q

f

B





2 3 4

5 2 3 4

5

(47)
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

22 of 78

NESC Request No.: 07-037-I

where C represents either cell on the left or the cell on the right-side of the face and j represents 

the column number.

7.2 Numerical Implementation

An overall flow diagram of the USM3D code main time loop with the current modifications to 

add real-gas capability is given below. The sections highlighted in bold letters represent new 

code additions. The USM3D code was modified to add variable ratio of specific heats. The 

modifications added to the different code modules are described in this section.

DO 25 ntsUN = 1, ntstep

Move grid

Qnode

UBC2

DO 10 nts = 1, niter

UCTime

Limiter

Turbulence model -> Update viscosity

Calculate species source terms based on equilibrium or finite chemistry

Solve species transport equations

Save ()

n+1

Viscous flux

Update (Implicit scheme)

Calculate RHS

Calculate LHS

Calculate Q

Update state variables

Update species mass fraction





n 1





n 1



n1

Calculate new species mass fraction based on equilibrium chemistry 

model if the effect is not included in the species source term
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

23 of 78

NESC Request No.: 07-037-I

Update temperature and pressure using new species mass fractions and 

state variables

Qnode

UBC2

Force

........

10 Continue

25 Continue

Roe’s Approximate Riemann Solver

Based on Roe’s approximate Riemann solver, the numerical flux that crosses a cell-face can be 

written as:

       FLR

 F QL F QR A QR QL

2

1

(48)

where L and R represent the cells on the left and right sides of the face respectively, Q is the 

conserved variable vector, and 

A

is the Roe averaged matrix. The flux vector F(Q) is defined 

as:



F(Q) 

U

uU  pnx

vU  pny

wU  pnz

U(E  p)

































(49)

where E is the total energy per unit mass and U is the contravariant velocity. The variables used 

in the Roe averaged matrix are evaluated as:

   L R

(50)

 

R L

uLuR R L

u

 

 

1 /

/







(51)

 

R L

L R R L

v v

v

 

 

1 /

/







(52)
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

24 of 78

NESC Request No.: 07-037-I

 

R L

wLwR R L

w

 

 

1 /

/







(53)

 

R L

h Lh R R L

h

 

 

1 /

/

0 0

0







(54)

The diffusion term (i.e., last term in Equation 48) can be split it into three components as:

  A QR

QL F1 F2 F3

 (55)

where, 

  



































































      

  

  

  







































 











 

   

u u v v w w U U

w n U

v n U

u n U

u v w

w

v

u

a

p

F U

z

y

x

0

2

1

1

2 2 2

1 2

 

(56)



















































   

  

h U a

w n a

v n a

u n a

a

p a U

F U a

z

y

x

0

4,5 2

1

2



(57)

where ∆ represents the change in the variable value between the cell on the right and cell on the 

left. 

The different components (i.e., mass, momentum, and energy) of the flux vector can be written:

 

f

1



UL

 UR

2



U 

2

 

 p

a 

2



  U  a 

 p  a U

4a 

2







 U  a 

 p  a U

4a 

2









(58)
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

25 of 78

NESC Request No.: 07-037-I



f

2



uUL

 uUR pL pR

nx

2



U 

2

u  

 p

a 

2

 

         u  nx

U

 

   

 

   

 U  a 

 p  a U

4a 

2

 

   

 

   u  nx

a  U  a 

 p  a U

4a 

2

 

   

 

   u  nx

a 

(59)



f

3



vUL

 vUR p L pR

ny

2



U 

2

v  

 p

a 

2

 

         v  ny

U

 

   

 

   

 U  a 

 p  a U

4a 

2

 

   

 

   v  ny

a  U  a 

 p  a U

4a 

2

 

   

 

   v  ny

a 

(60)



f

4



wUL

 wUR pL pR

nz

2



U 

2

w  

 p

a 

2

 

         w  nz

U

 

   

 

   

 U  a 

 p  a U

4a 

2

 

   

 

   w  nz

a  U  a 

 p  a U

4a 

2

 

   

 

   w  nz

a 

(61)



f

5



1

2

UE  pL

 UE  pR

U  

 p

a 

2

 

       1

2

u 

2

 v 

2

 w 

2

 

 

       u u  v v  w w U U

 

  

 

  

U  a 

 p  a U

2a 

2

 

  

 

  h  0

 U a  U  a 

 p  a U

2a 

2

 

      h  0

 U a 

 

 

 

 

 

 

 

 

  

 

 

 

 

 

 

 

 

  

(62)



f

5



Uh0L

 Uh0

R

2



U 

2

 

 p

a 

2

 

        V 2

2



U 

2

u u  v v  w w U U

 U  a 

 p  a U

4a 

2

 

   

 

   h  0

U a  U  a 

 p  a U

4a 

2

 

        h  0

U a 

(63)

For flows with multiple species, the ratio of specific heats used in the Roe averaged matrix is 

evaluated in different ways by different authors. Two of the commonly used methods are 

described in the following sections.
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

26 of 78

NESC Request No.: 07-037-I

Methods used by Luke [ref. 2]

The ratio of specific heats is evaluated as:

*

1

Cv

R

  

(64)

where,





NS

s

R sRs

1



(65)





NS

s

Cv sCvs

1

* *



, and (66)







R

L

T

T

vs

R L

vs C dT

T T

C

*

1

(67)





s



s,L

s,R R L 

1  R L

(68)



T 

TL

 TRR L 

1 R L

(69)



e s

es,L e  s,R R L 

1  R L

(70)

The speed of sound is calculated as:

    











      

N S

s

v s s

a h u v w C T e

1

2 2 2 *

0

2

2

1

 1 

(71)

Method used by Molvik and Merkle [ref. 3]

Molvik and Merkle used the following approach to calculate the speed of sound and the ratio of 

specific heats. The initial estimation for the speed of sound, ratio of specific heats, and the 

chemical source term are defined as:

2

L R

a a

a





(72)

2

L R

 







(73)
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

27 of 78

NESC Request No.: 07-037-I

2

sL sR

s

D D

D





(74)

These values are refined to meet Roe’s criteria as:

 

  

1 /

/

4

p p

D a p

D

s s

s

 

 



(75)

1

1 /

1









 





p p

(76)

where ∆ represents the change in the value between right and left states, and the other variables 

appearing in Equation 76 are calculated as:

     



       

N S

s

p p Ds sh p

1

   1 

(77)

   

2

1

2 2

a p

NS

s

  s 



 

(78)

The Roe averaged speed ( ) of sound is calculated as:

  

  

NS

s

a h sDs

1

2

 1 

(79)

Methodology Implemented in USM3D Code

Currently, in the USM3D code, the Roe averaged speed of sound is calculated based on constant 

specific heat assumption expressed as:

       











         

2 2 2

v

2

2

1

a  RT   1C T   1 e   1 E u v w

(80)

This has been modified to calculate the speed of sound based on variable ratio of specific heats.

The Roe averaged density, velocity components, and total energy were calculated using 

Equations 50 through 54. The other Roe averaged correlations used for different flow variables 

are listed as:





s



s,L

s,R R L 

1  R L

(81)

      

NS

s

s s

e E u v w e

1

2 2 2

2

1



(82)

a
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

28 of 78

NESC Request No.: 07-037-I

Since 

de  CvdT

, the Newton-Raphson method was used to solve for 

( )

1



n

T T

as:

n

n n

n n

C

e e

T T

v

1

1 

 





(83)

C R

R

e

R

e

T T

n

n n

n n

v

1

1





 





(84)

 

1 1 1

1

o o

v

     



 

N S

s u

s

p

s

u

p p p

R

C

n

R

C

R

C

R

C R

R

C

(85)

T

R

h

T n

R

h

T

R

h

R

h RT

R

e

N S

s u

s

s

u

     



 1

o o

(86)

 



 NS

s

s s

s s

s

M

M

n

1





(87)

    u s sCpR a T a T a a T a T a T a T

4

7

3

6

2

3 4 5

1

2

2

1

o

      

 

(88)

  s u s

h R aT a T aT a T a T a T a T a8

5

5 7

4 1

4 6

3 1

3 5

2 1

2 4

1

2 3

1

1

o

   ln      



(89)

Once the Roe-averaged temperature is obtained, 

a , can be calculated based on:

a   RT

(90)

where





NS

s

R sRs

1



(91)

s

u

s M

R

R 

(92)

( ) 1



C R

C R

p

p



(93)
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

29 of 78

NESC Request No.: 07-037-I

Modifications to Boundary Conditions

Boundary Condition 0: Freestream (Supersonic Inflow)

 Specify and fix density, velocity components, species mass fraction (i.e., mixture gas 

constant), and pressure.

 Calculate mixture temperature and Gamma (i.e., ratio of specific heats) based on 

pressure, mixture density, and mixture gas constant using the thermal equation of state.

 Calculate total energy from species mass fractions and mixture temperature, and velocity 

based on the caloric equation of state.

Boundary Condition 1: Tangent Flow (Reflection Plane)

 Set contravariant velocity component to zero. 

 Set pressure, mixture density, and species mass fractions at the boundary to be the same 

as the values of the cell values. 

 Calculate mixture temperature and Gamma based on pressure, mixture density, and 

mixture gas constant using the thermal equation of state.

 Calculate total energy from species mass fractions and mixture temperature, and velocity 

based on the caloric equation of state.

Boundary Condition 2: Full Extrapolation (Supersonic Outflow)

 Extrapolate velocity components, pressure and mixture density at the boundary from the 

cell values.

 Set species mass fractions at the boundary to be the same as the values of the cell inside.

 Calculate mixture temperature and Gamma at the boundary based on pressure, mixture 

density, and mixture gas constant using the thermal equation of state.

 Calculate total energy at the boundary from species mass fractions and mixture 

temperature, and velocity based on the caloric equation of state.

Boundary Condition 3: Characteristic Inflow/Outflow (Subsonic Outer Boundaries)

Appendix C provides a detailed description.

 Assumptions

o Species mass fractions at the boundary are the same as species mass fractions inside if 

the contravariant velocity is positive (i.e., outflow).

o Species mass fractions at the boundary are the same as freestream species mass 

fractions if the contravariant velocity is negative (i.e., inflow).

 Use Riemann invariants to calculate the boundary contravariant velocity and the speed of 

sound. Use the speed of sound to calculate the temperature. 

The Riemann invariants corresponding to the outgoing and incoming characteristics are given by:

1

2



  





a

U , and 

1

2



  





a

U
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

30 of 78

NESC Request No.: 07-037-I

where 

 

is corresponding to the eigenvalue, U+a, and 

 

is corresponding to the eigenvalue 

U-a. These Riemann invariants at the boundary can be written as:







  

  

 







 if 0

 if 0

freestream

inside

U a

U a

b

(94)







  

  

  





 if 0

 if 0

freestream

inside

U a

U a

b

(95)

The boundary contravariant velocity (Ub) and speed of sound (

ab

) can be calculated by solving 

Equations 94 and 95. The boundary species mass fractions (

s b

( )

) will be the same as those of 

the cell at the upwind side. Since:

  b s b b b s b

b b b

b

R f f T

R T

a ; ( ) ; and ,( )

,

ref

  





  

An iterative procedure was used to determine the boundary temperature and the ratio of specific 

heats based on the boundary speed of sound and species mass fractions. With the isentropic 

assumption between the upwind cell and the boundary, the boundary mixture density and 

pressure can be calculated based on the temperature (Tc) and mixture density (c) of the upwind 

cell as:

 

b b b b

c c

b

b

p R T

T

T b

b



 











 















1

; and

1

1

1

Total energy at the boundary was calculated from boundary temperature, species mass fraction, 

and velocity components based on the caloric equation of state.

Boundary Condition 4: No-Slip (Viscous Surfaces)

Appendix D provides a detailed description.

 Stationary wall: set velocity components at the wall boundary to be zero and velocity 

components at the ghost point to have the opposite direction to those at the neighboring 

interior cell.

 Moving wall: set velocity components at the wall boundary to be the prescribed speed of 

the moving wall and extrapolate velocity components at the ghost point from the 

neighboring interior cell and wall boundary.

 Extrapolate pressure at the wall boundary from neighboring interior cells, and extrapolate 

pressure at the ghost point from the neighboring interior cell and wall boundary.

 Set the species mass fractions to be the same as those of the neighboring interior cells.
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

31 of 78

NESC Request No.: 07-037-I

 Adiabatic: in the current implementation, temperature at the wall boundary is set to the 

same as its neighboring interior cell.

 Isothermal: temperature at the wall boundary remains unchanged.

 Calculate mixture density: total energy at the wall boundary is calculated from 

temperature, pressure, species mass fraction, and velocity components based on the 

thermal and caloric equation of state.

Boundary Condition 5: Tangent Flow (Inviscid Aerodynamic Surface)

 Same as boundary condition 1.

Boundary Condition 44: Special no-slip (Blunt base)

 Same as boundary condition 4.

 If wall function, then it is different than boundary condition 4. 

8.0 Verification of Multi-Species Capability

The implementation of the numerical models from Section 7.0 was verified using hot-gas flow 

conditions though a round convergent-divergent nozzle exiting into ambient air. The solutions

presented in the following sections will only serve as qualitative verification of the 

implementation of the real-gas energy and multi-species transport equations, and are not 

validated against experimental data. 

The configuration consists of a plenum chamber, which is represented as the small straight 

section on the left-side of the geometry shown in Figure 8.0-1, followed by the convergent 

divergent nozzle. A computational grid of 2,706,755 tetrahedral cells was constructed for one 

quadrant of the round nozzle configuration. Figure 8.0-1 depicts the grid clustering on one of the 

quadrant boundaries. The hot jet flow exits the nozzle at nominally Mach 2 into near-quiescent

ambient air. The freestream Mach number is set to 0.05 for numerical purposes. The plenum is 

set to total pressure, ptot, of 115 psi and total temperature, Ttot, of 2009 degR to produce a hot jet 

plume. The chemical composition of the plenum gas is varied to simulate different flow 

scenarios. The inflow boundary condition to the plenum is specified using Engine boundary 

condition in USM3D code. A description of the USM3D input parameters is presented in 

Appendix E.
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

32 of 78

NESC Request No.: 07-037-I

Figure 8.0-1. Triangulation on Boundary Plane of Tetrahedral Quadrant Grid for Nozzle Test Case

8.1 Comparison of Ideal and Real Gas Simulations

This section presents the comparison of the simulation results from the ideal- and real-gas 

versions of the USM3D code. The results confirm the expected variation of the ratio of specific 

heats with temperature from the modified energy equation. The multi-species transport 

equations are not solved in this example. In the following comparison, the real-gas is taken as 

air with a mixture of 73.5 percent nitrogen (N2) and 26.5 percent oxygen (O2). This mass 

fraction is assumed to be the same for the entire flowfield. In these simulations, the CourantFriedrichs-Lewy (CFL) number is ramped from 1 to 10 in 800 iterations. The ideal-gas 

simulation was run for 9,000 iterations. By using 9,000 iterations for real-gas simulation, the 

flow was not fully developed and contour deviations for temperature and Mach number near the 

axis of rotation between the ideal- and real-gas simulations were noted. Therefore, the real-gas 

version of the USM3D code was run for 30,000 iterations. The comparison of the pressure, 

temperature, Mach number, and ratio of specific heats is given in Figure 8.1-1. This figure 

shows that contours of pressure, temperature, and Mach number are similar for both ideal- and 

real-air simulations. The variation of ratio of specific heats for the real-gas simulation can be 

seen in Figure 8.1-1(d). As expected, the ratio of specific heats is lower in the settling chamber 

where the temperature is high, and defuses toward the ambient gamma ratio of 1.4 as the 

temperature of the jet plume decreases toward the ambient condition.
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

33 of 78

NESC Request No.: 07-037-I

(a) Pressure

(b) Temperature
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

34 of 78

NESC Request No.: 07-037-I

(c) Mach Number

(d) Ratio of Specific Heats

Figure 8.1-1. Comparison of Results from Ideal and Real Air Simulations
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

35 of 78

NESC Request No.: 07-037-I

8.2 Comparison of Flux Schemes with N2 Exhaust Gas

This section presents the comparison of the simulation results using Roe and van Leer flux 

schemes from the real-gas version of USM3D code. In the following comparison, the ambient 

gas is taken as a mixture of 73.5 percent N2 and 26.5 percent O2. The plenum/nozzle gas is 

prescribed as 100 percent N2. The species transport was solved to propagate the species into the 

domain, along with the energy equation that modifies the ratio of specific heat with temperature. 

Solutions were started with first-order spatial differencing and switched to second-order. A 

comparison of density, pressure, temperature, velocity magnitudes, Mach number, and species 

mass fractions are shown in Figure 8.2-1. The species mass fraction plots show the N2 mass 

fraction inside the settling chamber is 1.0, and O2 mass fraction is zero. It is observed that the 

distributions of flow variables computed by the Roe scheme and the van Leer scheme are in 

close agreement, thereby increasing confidence in the implementation of the species transport 

equations. An expected non-reacting transport of the N2 species from the engine exhaust into the 

ambient air was confirmed.

(a) Density
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

36 of 78

NESC Request No.: 07-037-I

(b) Pressure

(c) Temperature
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

37 of 78

NESC Request No.: 07-037-I

(d) Velocity Magnitude

(e) Mach Number
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

38 of 78

NESC Request No.: 07-037-I

(f) Ratio of Specific Heats

(g) Species Mass fraction: N2
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

39 of 78

NESC Request No.: 07-037-I

(h) Species Mass fraction: O2

Figure 8.2-1. Results from the Simulations using N2 as the Exhaust Gas

8.3 Results from Simulations using Hydrogen (H2) as the Exhaust Gas

This section presents the results from the simulation of the nozzle jet with H2 as the exhaust gas. 

With H2 exhaust gas, the equations are extremely stiff and required more effort to run the 

simulation. First, the solution was run with first-order differencing with a maximum CFL 

number of 5. The CFL ramping was set from 1 to 5 over the first 800 iterations. Second, the 

flow variables inside the plenum were initialized to the hot-gas H2 conditions to facilitate 

solution start up. While USM3D has a general capability for modifying prescribed cells, the gas 

chemistry extension was not added to that feature. For the current testing, the x-coordinate of the 

cell centroids was used to flag the cells in the plenum and initialize the variables. The number of 

iterations was set as 30,000, while maintaining the maximum value of CFL number as 5. The 

distribution of density, pressure, temperature, ratio of specific heats, and species mass fractions 

after 30,000 iterations are shown in Figure 8.3-1. This figure shows that the flow is not fully 

developed at the axis of symmetry due to the cell size. Additional iterations are needed to 

converge the solution and clean up the slight turn in the contour lines near the axis. 
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

40 of 78

NESC Request No.: 07-037-I

(a) Density

(b) Pressure

(c) Temperature

Air 

(N2,O2) 

H(N

2
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

41 of 78

NESC Request No.: 07-037-I

(d) Mach Number

(e) Ratio of Specific Heats

(f) H2 Mass Fraction
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

42 of 78

NESC Request No.: 07-037-I

(g) N2 Mass Fraction

(h) O2 Mass Fraction

Figure 8.3-1. Results from the Simulations using H2 as the Exhaust Gas

9.0 Findings, Observations, and NESC Recommendations

9.1 Findings

The following findings were identified:

F-1. The implemented real-gas energy equation solver enables the USM3D code to account 

for the effect of variable specific heats.

- Solutions require a smaller numerical stability parameter (CFL) to maintain numerical 

stability, which prolongs run times than its ideal-gas counterpart.
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

43 of 78

NESC Request No.: 07-037-I

F-2. The van Leer and Roe flux schemes were modified in the USM3D code for real-gas 

energy equation solver, with the van Leer flux scheme being the most stable.

F-3. The species transport equation enables USM3D code to model the multiple species effect. 

- A substantially reduced CFL number is needed to maintain numerical stability. The 

extent of CFL number reduction depends on the species system.

F-4. The numerical modules for solving chemical reactions were developed, but not 

implemented and tested because of the unresolved numerical stability issues for some 

multi-species mixing problems.

F-5. Continuation of this work is being supported by the SLS Program Aero Task Lead

through the services of Dr. Michael Applebaum, Dr. William Eppard, and Mr. Les Hall 

of CRM Solutions, Inc.

9.2 Observations

The following observations were identified by the Investigators:

O-1. The original USM3D flow solver was tightly embedded with the ideal-gas assumption 

making it difficult to locate every occurrence during the conversion to a real-gas model, 

thereby increasing the chance for inconsistencies.

O-2. Redundancy in the numerical flux scheme coding led to difficulties in implementing a 

consistent modification for real-gas equations. 

O-3. The USM3D code has several specialized boundary conditions that presented a challenge 

to adapt to real-gas formulation.

O-4. USM3D uses a one-dimensional array data structure throughout that requires one to have 

an intimate knowledge of the entire code for making modifications. A more modular 

approach with locally defined arrays would be desirable.

9.3 NESC Recommendations

The following NESC recommendations are directed to the NASA Technical Fellow for 

Aerosciences and the Aerothermodynamics Community:

R-1. Seek recommendation from CRM Solutions, Inc., on best approach for maturing a robust 

multi-species and gas chemistry capability in USM3D.

(F-1 through F-5, O-1 through O-4)

R-2. Monitor the SLS Program support of the USM3D code modification for application to the 

MPCV Program crew module launch abort analysis. (F-5)
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

44 of 78

NESC Request No.: 07-037-I

10.0 Alternate Viewpoints

There were no alternate viewpoints identified during the course of this assessment by the NESC 

team or the NRB quorum.

11.0 Other Deliverables

There were no other deliverables during the course of this assessment. 

12.0 Lessons Learned

No applicable lessons learned were identified for entry into the NASA Lessons Learned 

Information System.

13.0 Definition of Terms 

Corrective Actions Changes to design processes, work instructions, workmanship practices, 

training, inspections, tests, procedures, specifications, drawings, tools, 

equipment, facilities, resources, or material that result in preventing, 

minimizing, or limiting the potential for recurrence of a problem. 

Finding A conclusion based on facts established by the investigating authority. 

Lessons Learned Knowledge or understanding gained by experience. The experience may 

be positive, as in a successful test or mission, or negative, as in a mishap 

or failure. A lesson must be significant in that it has real or assumed 

impact on operations; valid in that it is factually and technically correct; 

and applicable in that it identifies a specific design, process, or decision 

that reduces or limits the potential for failures and mishaps, or reinforces a 

positive result. 

Observation A factor, event, or circumstance identified during the assessment that did 

not contribute to the problem, but if left uncorrected has the potential to 

cause a mishap, injury, or increase the severity should a mishap occur. 

Alternatively, an observation could be a positive acknowledgement of a 

Center/Program/Project/Organization’s operational structure, tools, and/or 

support provided.

Problem The subject of the independent technical assessment.

Proximate Cause The event(s) that occurred, including any condition(s) that existed 

immediately before the undesired outcome, directly resulted in its 

occurrence and, if eliminated or modified, would have prevented the 

undesired outcome.

Recommendation An action identified by the NESC to correct a root cause or deficiency 

identified during the investigation. The recommendations may be used by 
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

45 of 78

NESC Request No.: 07-037-I

the responsible Center/Program/Project/Organization in the preparation of 

a corrective action plan.

Root Cause One of multiple factors (events, conditions, or organizational factors) that 

contributed to or created the proximate cause and subsequent undesired 

outcome and, if eliminated or modified, would have prevented the 

undesired outcome. Typically, multiple root causes contribute to an 

undesired outcome.

14.0 Acronyms List

CFL Courant-Friedrichs-Lewy

CxP Constellation Project

H2 Hydrogen

LaRC Langley Research Center

LAV Launch Abort Vehicle

MPCV Multi-Purpose Crew Vehicle

N2 Nitrogen

NESC NASA Engineering and Safety Center

NRB NESC Review Board

O2 Oxygen

RCS Roll Control System

SLS Space Launch System

14.1 Nomenclature

A Inviscid flux Jacobian or empirical constant in the Arrhenius form of chemical 

reaction rate

Ajk Components of the inviscid flux Jacobian

B Viscous flux Jacobian or empirical constant in the Arrhenius form of chemical 

reaction rate

Bjk Components of the viscous flux Jacobian

Cp Specific heat at constant pressure

Cv Specific heat at constant volume

dV Elemental volume

dS Elemental surface area

E Total energy per unit mass

e Internal energy per unit mass

Ea Empirical constant in the Arrhenius form of chemical reaction rate



F i

Inviscid flux vector

i i i

f f f

1 2 5

, ,,

Components of inviscid flux vector

v

F



Viscous flux vector
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

46 of 78

NESC Request No.: 07-037-I

v v v

f f f

1 2 5

, ,,

Components of viscous flux vector

h Enthalpy per unit mass

h0 Total enthalpy per unit mass

Kb, j

Backward reaction rate of the j

th reaction

K f , j

Forward reaction rate of the j

th reaction

Mw Molecular mass

M Mach number

n



Outward pointing unit normal 

NR Number of reactions

NS Number of species

nt Contravariant velocity component of the grid speed and is defined as

t t x t y t

nzn  x n  y n  z

nx, ny, nz Components of the unit vector 

n



p Pressure 

Pr Prandtl number 

Q



Conserved variable vector

Q1, Q2, ..., Q5 Components of the conserved variable vector

R Mixture gas constant

Re Reynolds number

Ru Universal gas constant

S Sutherland constant

T Temperature

t Time

u, v, w Components of velocity in x, y, and z directions, respectively

xt, yt, zt Components of grids speed in x, y, and z directions, respectively

Greek Letters

α Species mass fraction defined as the ratio of species density to total density

β Contravariant velocity component of the fluid particle with respect to the grid and is 

defined as

    nt

xy

z

, ,

Geometric factors

 Coefficient of thermal conductivity and is defined as 

Pr





Cp



 Ratio of specific heats

Density

 Dynamic viscosity





Source term due to chemical reaction




NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

47 of 78

NESC Request No.: 07-037-I

 Contravariant velocity component of the fluid particle and is defined as 

x y

wnz  un  vn 



Stoichiometric coefficient

xx xy xz yy yz z z  , , , , ,

Component of the rate of shear stress tensor

Subscripts

b Backward reaction

f Forward reaction

ref Reference conditions

s Species index

x, y, z Components of a vector along x, y, and z directions, respectively

 Freestream conditions

Over

^ Dimensional quantities

 Vectors

Tensors

15.0 References

1) Frink, N. T., “Tetrahedral Unstructured Navier-Stokes Method for Turbulent Flows,” 

AIAA Journal, Vol. 36, No. 11, November 1998, pp. 1975-1982.

2) Luke, E. A., “A Rule-Based Specification System for Computational Fluid Dynamics,” 

Ph.D. Dissertation, Mississippi State University, December 1999. 

3) Molvik, G. A., and Merkle, C. L., “A Set of Strongly Coupled, Upwind Algorithms for 

Computing Flows in Chemical Non-Equilibrium,” AIAA Paper 89-0199, presented at the 

27th Aerospace Sciences Meeting, Jan. 9-12, 1989, Reno, Nevada. 

16.0 Appendices

Appendix A. Derivation of Inviscid and Viscous Jacobians

Appendix B. Linearization of Pressure for Real Gases

Appendix C. Numerical Treatment of Characteristic Boundary Conditions in Subroutine UBC2

Appendix D. Numerical Treatment of Viscous Wall Boundary Conditions

Appendix E. Description of USM3D Input File
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

48 of 78

NESC Request No.: 07-037-I

Appendix A. Derivation of Inviscid and Viscous Jacobians

Representing the conserved variable vector Q as:





































































5

4

3

2

1

Q

Q

Q

Q

Q

E

w

v

u

Q













(A.1)

The contravariant velocity can be represented as:

  x y z

Q n Q n Q n

Q

2 3 4

1

1

   

(A.2)

Using Equation (5), the inviscid flux vector can be written as:

 

 

 

   





















































      

    

    

    

   



































 

2 3 4

1

1 2 3 4

1

5

1 2 3 4

1

4

1 2 3 4

1

3

1 2 3 4

1

2

1 2 3 4

5

4

3

2

1

n Q n Q n Q

Q

p

n Q n Q n Q n Q

Q

Q

n Q n Q n Q n Q n p

Q

Q

n Q n Q n Q n Q n p

Q

Q

n Q n Q n Q n Q n p

Q

Q

n Q n Q n Q n Q

f

f

f

f

f

F n

t x y y x y y

t x y y z

t x y y y

t x y y x

t x y y

i

i

i

i

i

i





 (A.3)

Equation (A.3) can be differentiated with respect to A.1 to get the inviscid flux Jacobian. Each 

component of the inviscid flux vector can be represented as 



Ajk 

 f

j

i

Qk

and this will result in:



































5 1 5 2 5 3 5 4 5 5

4 1 4 2 4 3 4 4 4 5

3 1 3 2 3 3 3 4 3 5

2 1 2 2 2 3 2 4 2 5

1 1 1 2 1 3 1 4 1 5

A A A A A

A A A A A

A A A A A

A A A A A

A A A A A

A

(A.4)
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

49 of 78

NESC Request No.: 07-037-I



A 

ntnxnynz0

nx  u nx

2   u U ny

u  nx

 1v nz

u  nx

 1w nx 1

ny  v nxv  ny

 1u ny2   v U nz

v  ny

 1w ny 1

nz  w nxw  nz

 1u ny

w  nz

 1v nz2   w U nz 1

2 

E



 

   

 

    nx

E



 

 

   

 

     1u ny

E



 

 

   

 

     1v nz

E



 

 

   

 

     1w   nt

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

 

and



A1 1 

f

1

i

Q1

nt; A1 2 

f

1

i

Q2

 nx; A1 3 

f

1

i

Q3

 ny; A1 4 

f

1

i

Q4

 nz; A1 5 

f

1

i

Q5

0

( 1)

( 1)

( 1)

(2 )

5

2

2 5

1

4

1

2

4

2

2 4

1

3

1

2

3

2

2 3

1

2

2

2

2 2

1

2

1

2 4

1

2 3

1

2

2

2

1 1

2

2 1

 







  







  







  











   















    



















 





x

i

z x

i

y x

i

x

i

t x y z x x

i

n

Q

f

A

Q

Q

n

Q

Q

n

Q

f

A

Q

Q

n

Q

Q

n

Q

f

A

Q

Q

n

Q

f

A

p

n

Q

Q

n p

Q

Q Q

n

Q

Q Q

n

Q

Q

n Q n

Q Q

f

A

( 1)

( 1)

(2 )

( 1)

5

3

3 5

1

4

1

3

4

3

3 4

1

3

3

3

3 3

1

2

1

3

2

3

3 2

1

3

1

3 4

1

2

3

1

2 3

3

1 1

3

3 1

 







  







  







  











   















    

















 







y

i

z y

i

y

i

x y

i

t x y z y y

i

n

Q

f

A

Q

Q

n

Q

Q

n

Q

f

A

Q

Q

n

Q

f

A

Q

Q

n

Q

Q

n

Q

f

A

p

n

Q

Q

n p

Q

Q Q

n

Q

Q

n

Q

Q Q

n Q n

Q Q

f

A
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

50 of 78

NESC Request No.: 07-037-I

( 1)

(2 )

( 1)

( 1)

5

4

4 5

1

4

4

4

4 4

1

3

1

4

3

4

4 3

1

2

1

4

2

4

4 2

1

4

1

2

4

1

3 4

1

2 4

4

1 1

4

4 1

 







  







  







  











   















    















 









z

i

z

i

y z

i

x z

i

t x y z z z

i

n

Q

f

A

Q

Q

n

Q

f

A

Q

Q

n

Q

Q

n

Q

f

A

Q

Q

n

Q

Q

n

Q

f

A

p

n

Q

Q

n p

Q

Q

n

Q

Q Q

n

Q

Q Q

n Q n

Q Q

f

A



A5 1 

 f

5

i

Q1





Q1

ntQ5 nx

Q2Q5

Q1

 ny

Q3Q5

Q1

 nz

Q4Q5

Q1

 nxp

Q2

Q1

 nyp

Q3

Q1

 nzp

Q4

Q1









  

Q5

Q1



 1

Q1

Q2

2

 Q3

2

 Q4

2

 









A5 2 

 f

5

i

Q2

 nx

Q5 p

Q1

 ( 1)

Q2

Q1

A5 3 

 f

5

i

Q3

 ny

Q5 p

Q1

 ( 1)

Q3

Q1

A5 4 

 f

5

i

Q4

 nz

Q5 p

Q1

 ( 1)

Q4

Q1

A5 5 

 f

5

i

Q5

   nt

where 

 

2

4

2

3

2

5 2

1 1

1

Q Q Q Q

Q Q

p p

   



 



 



is derived as shown in Appendix B. For constant 

specific heat assumption:

   

2 2 2 2

4

2

3

2

2

1

2

1

2

1

Q Q Q u v w

Q

p

 



  







  


NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

51 of 78

NESC Request No.: 07-037-I

Also, representing the dependent variable vector q as:







































































































5

4 1

3 1

2 1

1

5

4

3

2

1

q

Q Q

Q Q

Q Q

Q

q

q

q

q

q

p

w

v

u

q





Hence,



































































































































































5

5

4

5

3

5

2

5

1

5

5

4

4

4

3

4

2

4

1

4

5

3

4

3

3

3

2

3

1

3

5

2

4

2

3

2

2

2

1

2

5

1

4

1

3

1

2

1

1

1

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q

Q

q









































      

















( 1) ( 1) ( 1) 1

/ 0 0 1/ 0

/ 0 1/ 0 0

/ 1/ 0 0 0

1 0 0 0 0

   



 

 

 

u v w

p

w

v

u

Q

q





where:

1; 0

5

1

4

1

3

1

2

1

1

1

1

1





































Q

q

Q

q

Q

q

Q

q

Q

Q

Q

q

; 0

1 1

;

5

2

4

2

3

2

2 1

2 1

2

2

2

1

2

1

2 1

1

2



















 









 

  











Q

q

Q

q

Q

q

Q Q

Q Q

Q

u q

Q

Q

Q

Q Q

Q

q

 

; 0

1 1

;

5

3

4

3

2

3

3 1

3 1

3

3

2

1

3

1

3 1

1

3



















 









 

  











Q

q

Q

q

Q

q

Q Q

Q Q

Q

v q

Q

Q

Q

Q Q

Q

q

 
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

52 of 78

NESC Request No.: 07-037-I

; 0

1 1

;

5

4

3

4

2

4

4 1

4 1

4

4

2

1

4

1

4 1

1

4



















 









 

  











Q

q

Q

q

Q

q

Q Q

Q Q

Q

w q

Q

Q

Q

Q Q

Q

q

 

And q5/Q1, q5/Q2, q5/Q3, q5/Q4, q5/Q5 can be obtained as shown in Appendix B and 

shown as the following:

   

   

   

1

1 1

1 1

1 1

1

1

5

5

1

4

4

5

1

3

3

5

1

2

2

5

2

4

2

3

2

5 2

1 1

5

 





     





     





     





















   



















 

 

 







Q

q

Q

Q

w

Q

q

Q

Q

v

Q

q

Q

Q

u

Q

q

Q Q Q Q

p

Q

P

Q

q

Hence,





































      

















( 1) ( 1) ( 1) 1

/ 0 0 1/ 0

/ 0 1/ 0 0

/ 1/ 0 0 0

1 0 0 0 0

   



 

 

 

u v w

p

w

v

u

Q

q





(A.5)

From Equation 14, the team can obtain:

d p

R

d

T

dT









ref

  

Hence,



















R

T

q

T







ref

 0 0 0

(A.6)
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

53 of 78

NESC Request No.: 07-037-I

Combine Equations (A.5) and (A.6), the team can obtain:

 











        

















       















































      









































































R R

w

R

v

R

V u

R

T

R R

w

R

v

R

p u

R

T

u v w

p

w

v

u

R

T

Q

q

q

T

Q

T

Q

T

Q

T

Q

T

Q

T

Q

T



 



 



 



  









 



 



 



 

 





   



 

 

 







( 1) ( 1) ( 1) ( 1)

2

1

 

( 1) ( 1) ( 1) ( 1)

( 1) ( 1) ( 1) 1

/ 0 0 1/ 0

/ 0 1/ 0 0

/ 1/ 0 0 0

1 0 0 0 0

0 0 0

ref ref ref ref

2

ref

ref ref ref ref ref

ref

1 2 3 4 5





 

(A.7)

The spatial derivatives a dependent variable  at the interface between two adjacent cells can be 

calculated as:

n x n y n z

n

z

n x n y n z

n

y

n x n y n z

n

x

x y z

z

z z z

x y z

y

y y y

x y z

x

x x x

    

    







    

    







    

    







;

;

;



















and 

 RL

where  Rand  Lare the values of  at the center of right and left. Hence, the 

spatial derivatives of a dependent variable at the cell face can be linearized by the conservative 

Q



as:

 



















 

 



 

 





Q Q Q

R L

x x

x

  

 ( )  

(A.8)

Therefore,





















 





R

R

x

R

x

Q Q

 

 
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

54 of 78

NESC Request No.: 07-037-I





















 





L

L

x

L

x

Q Q

 

 

Equation (13) can be differentiated with respect to (A.1) to get the viscous flux Jacobian. Each 

component of the viscous flux vector can be represented as 

k

v

j

jk Q

f

B







and this will result in:



































5 1 5 2 5 3 5 4 5 5

4 1 4 2 4 3 4 4 4 5

3 1 3 2 3 3 3 4 3 5

2 1 2 2 2 3 2 4 2 5

1 1 1 2 1 3 1 4 1 5

B B B B B

B B B B B

B B B B B

B B B B B

B B B B B

B

(A.9)

Substituting Equations (A.8), (A.10), and (A.11), each components of the viscous flux Jacobian 

can be expressed as:

B11  B12  B13  B14  B15  0

        







































c j

z x

z

c j

y x

y

c j

x x

x

c j

v

j

Q

n

Q

n

Q

n

Q

f

B

   2

2

       









































C j

z y

z

C j

y y

y

C j

x y

x

C j

v

j

Q

n

Q

n

Q

n

Q

f

B

  

3

3

       









































C j

z z

z

C j

y z

y

C j

x z

x

C j

v

j

Q

n

Q

n

Q

n

Q

f

B

   4

4

 

 

 

 

 

 

       

      































   



























































C j

z

z

C j

y

y

C j

x

f j f j f j x

C j

z

z

C j

y

y

C j

x

x

C j

v

C j

v

C j

v

f

C j

v

j

Q

T

n

Q

T

n

Q

T

u B v B w B n

Q

T

n

Q

T

n

Q

T

n

Q

w f

Q

v f

Q

u f

Q

f

B





2 3 4

5 2 3 4

5

where C represents either cell on the left or cell on the right side of the face and j represents the 

column number.
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

55 of 78

NESC Request No.: 07-037-I

Appendix B. Linearization of Pressure for Real Gases

For thermally perfect, single-species flow:

ref 

RT

p 

Hence,

dT Rd

RT

d p

ref ref 







 

(B.1)

Since the specific internal energy (e) is defined as:

ref ref v

v

1

C

d e d T d T C

d e  

 

(B.2)

Also, the specific total energy (E) is defined as:

 

2 2 2

2

E  e

1

u v w

Hence,

dE  deuduvdvwdw, or de  dEuduvdvwdw

(B.3)

For a general flow variable 

 

 



 

  

d d

d d d d

 

       

Apply the above equation for flow variables, E, u, v, and w to Eq. (B.3), the team can obtain:

          



d E Ed u v w d u d u vd v wd w

de       

2 2 2

(B.4)

Substitute Equation (B.4) into Equation (B.2):

          

v

2 2 2

ref

1

C

d E Ed u v w d ud u vd v wd w

dT



     



      



(B.5)

Substitute Equation (B.5) into Equation (B.1), the team can obtain:
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

56 of 78

NESC Request No.: 07-037-I

    

 

5

5

4

4

3

3

2

2

1

1

2 2 2

2 2 2

ref

v

v

2 2 2

v

( 1) ( ) ( 1) ( )

( 1) ( ) ( 1) ( )

( 1)

1

( ) ( ) ( ) ( )

( ) ( ) ( ) ( )

dQ

Q

p

dQ

Q

p

dQ

Q

p

dQ

Q

p

dQ

Q

p

wd w d E

E u v w d ud u vd v

p

E u v w d ud u vd v wd w d E

C T

C

R

d E Ed u v w d ud u vd v wd w

C

R

d p RTd































   

    















   



 













    















    

        

   

    

 



    



      

(B.6)

where





































































E

w

v

u

Q

Q

Q

Q

Q

Q











5

4

3

2

1



and

   

   

   

   

1

1 1

1 1

1 1

1

( 1)

1

5

1

4

4

1

3

3

1

2

2

2

4

2

3

2

5 2

1 1

2 2 2

1

 





     





     





     









    



  















   



 







 

 

 





 



Q

p

Q

Q

w

Q

p

Q

Q

v

Q

p

Q

Q

u

Q

p

p

Q Q Q Q

Q Q

p

E u v w

p

Q

p

(B.7)
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

57 of 78

NESC Request No.: 07-037-I

Appendix C. Numerical Treatment of Characteristic Boundary 

Conditions in Subroutine UBC2

Current formulation (constant  = ):

Figure C-1. Orientation of Boundary Face for Characteristic Variable based Boundary Condition

Based on the orientation of a boundary face and the nomenclature used in Figure C-1, the 

Riemann invariants can be written as: 

1

2

1

2

1

2

1

2



 



  



 



  





 

 

b

b

L

L

b

b

R

R

a

U

a

U

a

U

a

U

where U and a are the contravariant velocity and speed of sound, and

  



UR Uc U  s ; ULUgU



Subscript ∞ denote the freestream condition. Hence, the condition at the boundary can be 

obtained as:

   

    

   

4

 

1

2

1

;



Ubab

If 

UbUm0

where Um is the contravariant velocity of the moving boundary, then the 

information is traveling from the right-hand-side. Hence,
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

58 of 78

NESC Request No.: 07-037-I

 

 

  b R z b R

b R y b R

b R x b R

w w n U U

v v n U U

u u n U U

  

  

  

To obtain the thermal properties at the boundary, isentropic process between R and b is 

assumed, and thus

b b b

R R

b

b

R

R

b

b

p T

T

T T T



 



 



  

1

and

1

1

1 1 1

 











  



  

where

p p p s s

p

T a T R c R c

R

R

b b R

 

           



;  ; ;

2

If 

UbUm0

, then the information is traveling from the left-hand-side. Hence,

 

 

  





   

   

   

w w n U U ; w w

v v n U U ; v v

u u n U U ; u u

b L z b L L

b L y b L L

b L x b L L

To obtain the thermal properties at the boundary, isentropic process between L and b is 

assumed, and thus

b b b

L L

b

b

L

L

b

b

p T

T

T T T



 



 



  

1

and

1

1

1 1 1 













  



  

where

 



   



 L L

L

L

L

p p

p

T ; ;

Modified implementation:

1

2

1

2

1

2

1

2



 



  



 



  





b

b

b

L

L

L

b

b

b

R

R

R

a

U

a

U

a

U

a

U

 

 

(C.1)

and

        



 







URUcU s; ULUgU ; s Ls gs; RLR



(C.2)
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

59 of 78

NESC Request No.: 07-037-I

where αi is the mass fraction of i-th species and R is the mixture gas constant. Subscript ∞

denotes the freestream condition. Hence, the condition at the boundary can be obtained as:

   

   

 



    

4

1

;

2

1 b

Ubab



(C.3)

If 

UbUm0

where Um is the contravariant velocity of the moving boundary, then the 

information is traveling from the right-hand-side. Hence,

 

 

 

        s b s R s c b R s R

b R z b R

b R y b R

b R x b R

R R f

w w n U U

v v n U U

u u n U U

     ;   ( )

  

  

  

(C.4)

Once the speed of sound and the mixture gas constant at the boundary are obtained, an 

iterative procedure can be used to determine the temperature and the ratio of specific heats at 

the boundary based on the following correlations. 

 

 

 

  1

;

ref 





 

b

p

b

p

b bp

b

p

b

b b b

b

C R

C R

C R

R T C

a 





(C.5)

 

  











































 















N S

s

s s

s s

s

b

N S

s u

s

p

s

b

u

p

b

p

M

M

n

R

C

n

R

C

R

C

1

1

o o

;





(C.6)

    u s sCpR a T a T a a T a T a T a T

4

7

3

6

2

3 4 5

1

2

2

1

o

      

 

(C.7)

With the isentropic assumption between the upwind cell and the boundary, the mixture 

density and pressure at the boundary can be calculated based on the temperature (Tc) and 

mixture density (c) of the upwind cell as:

 

b b b b

c c

b

b

p R T

T

T b

b



 







ref

1

1

1

1



; and 

















(C.8)

Total energy at the boundary is then calculated from temperature, species mass fraction, and 

velocity components at the boundary based on the caloric equation of state.
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

60 of 78

NESC Request No.: 07-037-I

If 

UbUm0

, then the information is traveling from the left-hand-side. Hence,

 

 

 

      







   

   

   

   

R R R

w w n U U w w

v v n U U ; v v

u u n U U ; u u

s b s L s b R

b L z b L L

b L y b L L

b L x b L L

;

;

  

(C.9)

The procedure stated above used to solve Equation (C.5) through (C.7) can be used to obtain 

the thermal properties at the boundary. Finally, the total enthalpy at the boundary can be 

calculated as:

     

    w E  EE C R dT u v b b b g

T

b

b p

b

and 2 2 2

2

1

0

(C.10)
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

61 of 78

NESC Request No.: 07-037-I

Appendix D. Numerical Treatment of Viscous Wall Boundary 

Conditions

Figure D-1. Orientation of a Wall Boundary Face

Based on the notations used in Figure D-1, the boundary conditions for a solid wall boundary can 

be summarized as:

b c g b c

g b c g b c g b c

b w b w b w

p p p s p p p

n

p

u u u v v v w w w

u u v v w w

        





     

  

0 ; 2

2 ; 2 ; 2

; ;



Current formulation:

1) adiabatic wall:

g

g

g

b

b

b

b c g c

T

p

T

p

T T T T

n

T

ref ref ;

0 ;

     

   





2) isothermal wall:

g b c

b

b

b w b

T

p

T T

  

 

 

  

2

ref
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

62 of 78

NESC Request No.: 07-037-I

Modified implementation:

s b s c s g s c

s

n

0 ( ) ( ) ; ( ) ( )



   





1) adiabatic wall:

¶T

¶n

= 0 Þ Tb = Tc; Tg = Tc

Þ rb = g ref

pb

Tb

; eb = CvdT

0

Tb

ò

rg = g ref

pg

Tg

; eg = CvdT

0

Tg

ò

and Eb = eb +

1

2

ub

2

+ vb

2

+wb

2

( ) ; Eg = eg +

1

2

ug

2

+ vg

2

+wg

2

( )

2) isothermal wall:

 

 

2 2 2

2

1

0

ref v

ref

2 2 2

2

1

0

v

2 and

and

g g g

T

g

g

g

g b c g

b

b

b b b b

T

b w b

E C d T u v w

p

T

T

p

T T E C d T u v w

g

b

       

      







   

 
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

63 of 78

NESC Request No.: 07-037-I

Appendix E. Description of USM3D Input File

USM3D (v6.0) - Sample Input File

 Mach alpha beta ReUe,mil Tinf,dR itwall Tw/Tinf ipwall

 0.84000 3.06000 0.00000 11.7 460.00000 0 -1.0 0

 sref cref bref xmc ymc zmc

 0.52550 0.66700 1.00000 0.00000 0.00000 0.00000

 ioverset impl cfl1 iramp cfl2 cflmin GS_tol crelax

 0 1 -150.000 50 150. 5.0 -20 1.00

 itimeacc delta_t ntstep res_step imvgrd isolavg timebgn

 -2 0.3 500 -3.0 1 0 1500.0

 irest mstage iresmth dqmax p_break pmin limiter

 0 3 1 0.5 0.050 0.001 0.

 nupdate nwrest nwflo nwflobgn ipltqn idiagnos nodeypl

 50 75 100 1500 2 0 0

 iorder lapl-avg high-bc ifds ivisc itrip ev_lim

 2 1 1 1 3 0 0.0

 ncyc nengines nsinkbc nrotor compF&M p_bc1002 cldes

 10 1 1 1 0 0.0 0.0

 jet1 fuel gammaj pjet p0jet Rratio T0jet

 103 0.02000 1.40000 0.71430 1.42857 0.77400 3.17390

 jet1 vec(1) vec(2) vec(3) engsim nstation idirec

 103 1.00000 0.00000 0.00000 0 10 1

 jet1 fuel gammaj pjet p0jet Rratio T0jet

 102 0.02000 1.40000 0.71430 1.42857 0.77400 3.17390

 jet1 vec(1) vec(2) vec(3) engsim nstation idirec

 102 1.00000 0.00000 0.00000 1 10 1

 jet1 inner_R outer_R cen_X cen_Y cen_Z

 102 1.0114 1.39846 23.646 8.713 -3.322

(r-inner_R)/(outer_R-inner_R)

0.168,0.350,0.500,0.600,0.700,0.800,0.900,0.950,0.975,1.00

p0jet

1.07143 1.07143 1.07143 1.07143 1.07143 1.07143 1.07143 1.07143 1.07143 1.07143

T0jet

1.32610 1.32610 1.32610 1.32610 1.32610 1.32610 1.32610 1.32610 1.32610 1.32610

V_circular/Vmag

0.0168 0.035 0.050 0.060 0.070 0.080 0.090 0.095 0.0975 0.1

V_radial/Vmag

0. 0. 0. 0. 0. 0. 0. 0. 0. 0.

 sinkid outer_R cen-X cen-Y cen-Z massflx omega

 111 1.0 0.0 0.0 0.0 3.0 -0.5

 rotid inner_R outer_R cen-X cen-Y cen-Z

 502 0.00000 24.0000 0.0 0.0 0.0

 Advratio thrst_c torq_c irotat rotsim

 0.20 0.00861 0.00 1 0

 Xrcg Yrcg Zrcg reflen iroll ipitch iyaw

 0.0 0.0 0.0 1.0 0 1 0

 ang_mean ang_ampl red_freq ideformg imvMOMCEN

 0.00 2.41 0.0808 0 0

 ikeord icons nstagek t_dtfact t_intsity mut/mul ratiokp dkemax

 1 1 5 1.0 1.e-3 1.0 0.00 0.25

 inl ilhg iwallf icompCorr itempCorr itk isk idt_proc

 7 -14 0 0 0 2 1 0
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

64 of 78

NESC Request No.: 07-037-I

 f1kemax <=== THIS IS A TEMPORARY VARIABLE FOR CAPPING MAX f1 in LB MODEL

 1.0

$SPECIES

 H2O O2 CO CO2 C7H8 H2 O H OH N2

$REACTIONS

 13

REACTION:H2O, O2, CO, CO2, C7H8, H2, O, H, OH, N2,

 1, 4.4963E9, 1.00, 2.679E4, 0, 2,

 0., -3.5, 7., 0., -1., 4., 0., 0., 0., 0.,

 0., -1.0, 0., 0., -0.5, 0., 0., 0., 0., 0.,

 2, 1.7000E13, 0.00,24070., 0, 0,

 0., -1., 0., 0., 0., -1., 0., 0., 2., 0.,

 3, 2.1900E13, 0.00, 2590., 0, 0,

 1., 0., 0., 0., 0., -1., 0., 1., -1., 0.,

 4, 6.0230E12, 0.00, 550., 0, 0,

 1., 0., 0., 0., 0., 0., 1., 0., -2., 0.,

 5, 1.8000E10, 1.00, 4480., 0, 0,

 0., 0., 0., 0., 0., -1., -1., 1., 1., 0.,

 6, 1.2200E17, -0.91, 8369., 0, 0,

 0., -1., 0., 0., 0., 0., 1., -1., 1., 0.,

 7, 4.0000E12, 0.00,10030., 0, 0,

 0., 0., -1., 1., 0., 0., 0., 1., -1., 0.,

 8, 3.0000E12, 0.00,2.50E4, 0, 0,

 0., -1., -1., 1., 0., 0., 1., 0., 0., 0.,

 9, 1.0000E16, 0.00, 0.,999, 0,

 0., 0., 0., 0., 0., 0., -1., -1., 1., 0.,

10, 2.5500E18, -1.00,59390.,999, 0,

 0., 1., 0., 0., 0., 0., -2., 0., 0., 0.,

11, 5.0000E15, 0.00, 0.,999, 0,

 0., 0., 0., 0., 0., 1., 0., -2., 0., 0.,

12, 8.4000E21, -2.00, 0.,999, 0,

 1., 0., 0., 0., 0., 0., 0., -1., -1., 0.,

13, 6.0000E13, 0.00, 0.,999, 0,

 0., 0., -1., 1., 0., 0., -1., 0., 0., 0.,

c==============================================================================

c LINE 1

c Title - free form upto 80 characters

c==============================================================================

c LINE 2

c xmach - freestream Mach number

c for xmach < 0, velocity will be initialized by 

c a fixed mach value of 0.4. Some usefulness for supersonic cases,

c however *not* recommended.

c alpha - angle of attack, deg. (in X-Z plane)

c beta - side-slip angle, deg. (in X-Y plane)

c reue - freestream Reynolds number per unit length (millions)

c tinf - freestream temperature (degrees Rankine)

c itwall - wall temperature boundary condition flag

c = 0 adiabatic wall temperature

c = 1 specified wall temperature

c twtinf - wall temperature (temperature at wall divided by temperature

c of freestream)

c if twtinf<=0, twtinf taken as freestream stagnation temperature)
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

65 of 78

NESC Request No.: 07-037-I

c ipwall - type of extrapolation of pressure to wall

c = 0 zeroth order

c = 1 first order (linear extrapolation)

c

c-----------------------------------------------------------------------

c LINE 2

c igmodel - flag to specify gas model

c = 0 Air (calorically perfect)

c = 1 Any gas with specified sutherland const and gamma

c (calorically perfect)

c = 2,-2 Any gas with specified sutherland const, gamma, and

c gas name (to set temperature variation of gamma)

c = 3,-3 Two-gas model with specified gamma, sutherland const,

c stagnation enthalpies, and gas names

c -2 or -3 is Recommended

c gamma1,gamma2 - specific heat ratios (gamma)

c htot1,htot2 - total enthalpies (non-dimensional)

c numspec - total number of species constituting 1- or 2-gas system

c NOTE : For igmodel values of 2 or 3, user must have "thermo.inp"

c file of McBride and Gordon.

c : For igmodel=2 or 3, you can set gamma1<=0; code will 

c calculate the value internally based on freestream temp.

c : For igmodel=3, you can set Htot1<= 0; code will calculate 

c the value internally. 

c : For igmodel=3, and exit BC -102, you can set Htot2<= 0; 

c code will calculate the value internally. 

c

c-----------------------------------------------------------------------

c LINE 3

c gasname - Name of a gas (must be part of McBride-Gordon table)

c ispec_grp - specie group index, must be 1 or 2 (2-gas model)

c mole_frac - Mole fraction of the gas in a given specie group 

c-----------------------------------------------------------------------

c LINE 5

c sref - reference area (half-area for half configuration)

c cref - reference length

c bref - reference span

c xmc - moment center in x-direction

c ymc - moment center in y-direction

c zmc - moment center in z-direction

c

c-----------------------------------------------------------------------

c LINE 6

c ioverset - overset grid flag

c = 0 single grid

c = 1 overset grid

c impl - implicit time stepping flag

c = 0 explicit Runge-Kutta time stepping 

c = 1 implicit Gauss-Seidel time stepping 

c = 1 UPDATQ, timestep scaling; UPDATQ enforces positivec definite updates for density & internal energy

c = 10 UPDATQ, timestep scaling off

c = 2 UPDATQ1, timestep scaling; UPDATQ1 *does not*

c change solver-provided solution updates

c = 20 UPDATQ1, timestep scaling off
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

66 of 78

NESC Request No.: 07-037-I

c = 3 UPDATQ2, timestep scaling

c = 30 UPDATQ2, timestep scaling off

c code will internally reset impl to 1 after raising

c appropriate flags.

c cfl1 - time step or initial clf number

c < 0 local time stepping or ramping

c iramp - number of cycles for ramping cfl number

c cfl2 - final value of cfl number for ramping

c cflmin - minimum allowable CFL for variable time step strategy

c GS_tol = Tolerance criterion for automatic termination of Gauss-Seidel

c iteration process termination for Flow and S-A equations

c recommended value (1-dq_prev/dq_new) is 0.001-0.005.

c = Negative value will be interpreted as G-S stages

c by the code for solving linear problem for mean flow.

c Minimum threshold value for stages is set to 5.

c crelax = Under/over-relaxation coefficient for Newton's scheme.

c = 0.7 as recommended by Dr. Paul Pao. 

c Without this relaxation, Newton's scheme produced jagged 

c hysteresis of Lift and moment coefficients for 

c NACA 0012 viscous pitching case. 

c

c *CAUTION* itimeacc=0, 1 or 2 options also will use relaxation 

c coefficient. For itimeacc=2, relaxation is not typically

c needed. Therefore it is important to set it to 1.0 when

c not using Newton's scheme (itimeacc=-2).

c-----------------------------------------------------------------------

c LINE 7

c itimeacc = 0 Not time accurate; assumes steady state final solution

c = 1 Time-accurate solution with subiterations using 

c Crank-Nicholson

c = 2 Time-accurate solution with subiterations using 

c 3-point backward differencing and pseudo-time 

c =-2 Time-accurate solution with subiterations using

c 3-point backward differencing without pseudo-time (Newton) 

c delta_t - time step increment (nondimensional) 

c HOW DOES ONE DETERMINE VALUE FOR TIMESTEP?

c -> Think of it as how many timesteps one needs for a fluid

c particle to travel some characteristic length, Lchar, 

c e.g. wing chord, at freestream velocity. A rule of thumb 

c from Jim Forsythe, Cobalt Solutions, is N=50 steps. 

c Scott Morton of the Air Force Academy suggests between 

c 100 and 1000 steps.

c -> For USM3D nondimensionalization, the relationship for

c timestep and number of steps, N, to traverse Lchar:

c delta_t = Lchar/xmach/N

c ntstep - number of timestep increments

c (Note that total "iterations" will be approx. ntstep*ncyc)

c res_step - level of subiteration convergence before taking next step

c Typical value might be -2 or -3.

c imvgrd - invokes prescribed grid motion 

c (additional input required, see LINE 5a below)

c isolavg - solution averaging flag

c = 0, turn-off solution averaging (for steady-state, it will

c be turned off by the code)

c = 1, perform time averaging for mean flow solution qunatities
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

67 of 78

NESC Request No.: 07-037-I

c (extra memory equivalent to 5*nnode). Also, mean values

c of aerodynamic quantities will be printed. 

c = 2, In addition to 1 above, perform solution averaging for 

c turbulent-stress-related qunatities (extra memory 

c equivalent to 11*nodes)

c timebgn - beginning timestep for initiating time-averaging of solution

c (based on the "ntstep" timestep numbering, not "ncyc")

c

c NOTES: 1) "ncyc" becomes the number of subiterations+1 used to converge

c solution to Q_(n+1)

c 2) Need to monitor hist.subit for subiteration convergence

c

c-----------------------------------------------------------------------

c LINE 8

c irest = 0 no restart

c = 1 restart

c mstage - number of Runge-Kutta stages (for explicit scheme)

c For implicit scheme mstage < 0 specifies interval at which 

c all the cells are reset as second-order cells

c iresmth - implicit residual smoothing (for explicit scheme)

c = 1 smoothing on every m-stage cycle

c = -1 smoothing only on odd m-stage cycles

c For implicit scheme iresmth specifies technique to mark 

c first order cells

c = 0 => do not set cells to first order locally

c = 1 => use p<=pminimum or rho<=pminimum switch to set

c first order cells

c = 2 => use p<=pminimum or rho<=pminimum or pt>=1.2*pt0

c switch to set first order cells iresmth=2 will

c only be used for non-engine, non-rotor cases

c dqmax - maximum allowable (-dp/p, -drho/rho) in solution

c (used to scale CFL number)

c = 0.5 suggested

c = greater than zero --> No Limiter Applied

c pbrk = threshold pressure for lowering local CFL toward CFLMIN

c pminimum= Minimum allowable level of non-dimensional pressure and

c threshold for setting cells to 1st order

c <0 also resets 1st-order flagged cells to 2nd-order

c bcoef - Limiter coefficent

c = 0. no limiting

c = 1. MINMOD limiter

c = 2. SUPERBEE limiter with variable compression factor

c = -2. will involke automatic selection of limiter 

c

c-----------------------------------------------------------------------

c LINE 9

c nupdate - number of iterations between updates of the time step 

c (not used for time accurate computations)

c nwrest - number of iterations between updates of the binary

c restart file

c - for UNSTEADY runs (itimeacc.NE.0), this is restart file

c updates are based on the number of "ntstep", not "ncyc"

c nwflo - number of iterations between writing of instantaneous 

c NEW flow-file (prj_name.flo.iter#). Typically used for

c time-accurate runs.
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

68 of 78

NESC Request No.: 07-037-I

c nwflobgn - number of iterations beyond which writing of instantaneous

c NEW flow-file (prj_name.flo.iter#) will be triggered. 

c The parameters nwflo and nwflobgn relate to outer loop 

c time steps (ntstep) for an unsteady flow analysis 

c and to inner loop time steps (ncyc) for steady flow analysis.

c ipltqn = 0 suppress output of post analysis of flow field file

c = 1 write unformatted q-node file for post analysis of flow field

c = 2 write formatted q-node file for post analysis of flow field

c = 3 write formatted q-node file for postprocessing by TECPLOT

c = 4 write binary TECPLOT file (NOT CODED YET)

c = 5 write unformatted q-node file for postprocessing by FIELDVIEW

c X =13 write unformatted q-node file for using Tecplot converter

c X =23 write formatted q-node file for using Tecplot converter

c =15 write unformatted q-node file for using FIELDVIEW converter

c =25 write formatted q-node file for using FIELDVIEW converter

c =-25 write solution restart file for using FIELDVIEW converter

c also write formatted q-node file and do not delete

c X =16 write unformatted q-node file for using Ensight converter

c X =26 write formatted q-node file for using Ensight converter

c X = **indicates that this option has not been tested and**

c **not recommended for use.**

c idiagnos= 0 no output

c 1 write history of minimum static pressure, maximum 

c Mach number, maximum eddy viscosity in info.diag file.

c 2 write a flow file with 5 auxiliary variables in place of Qs

c for unsteady flow simulations. These variables are,

c instantaneous and time-average values of vorticity,

c eddy viscosity, and local CFL number. 

c 3 write info.diag2 file in Tecplot format that contains

c several variables at cell-centers for the cells associated

c with user-specified boundary faces. These boundary faces 

c must be supplied in project.trip file. You will need to

c create a project.trip file for idiagnos != 3 as well, if

c itrip != 0 (below). NOTE: Diagnostics related to 

c idiagnos = 3 is available only for ivisc > 0.

c 4 write info.probe.N file in Tecplot format that contains

c history of Qs for each probe specified in project.probe

c file. N will be substituted by the probe number.

c 5 write "Sol.Volume.USM3D" file in PLOT3D format that contains

c Qs at desired locations in the grid as read from the file

c "Grid.Volume.TEACC" in a pre-processing step. In addition,

c "Sol.Exit.TEACC" file will be read at appropriate time to

c derive exit pressures for appropriate locations on the

c boundary. These locations are written out in the file

c "Grid.Exit.USM3D".

c nodeypl= # Boundary node number where output of u+/y+ is desired

c 0 no output

c

c-----------------------------------------------------------------------

c LINE 10

c iorder = order of flux computation

c = <0 (INTEGER only) Automatic selection of order based on

c prescribed magnitude of mean-flow residual error drop

c = 0 Automatic selection of order based on residual error drop

c by one order
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

69 of 78

NESC Request No.: 07-037-I

c = 1 first order

c = 2 Taylor series with full gradient computation

c laplavg = 0 Inverse-distance averaging to all nodes (Vers. 1.0)

c = 1 Pseudo-Laplacian averaging to all nodes

c = -1 (Ps. Lapl. on interior, Inv. Dist. on boundary nodes)

c ihibc = 0 First-order accurate boundary conditions (Vers. 1.0)

c = 1 Second-order accurate boundary conditions

c ifds - spatial differencing parameter for Euler fluxes

c = 0 flux-vector splitting

c = 1 flux-difference splitting (Roe scheme) (recommended)

c = 2 AUSM+ (need to evaluate for USM3D (as of 06/17/04). 

c = 3 HLLC (need to evaluate, only for stationary cases

c as of 07/12/05)

c ivisc - viscous/inviscid interaction flag 

c = 0 inviscid

c = 1 laminar (not active in version 2.0)

c = 2 Spalart-Almaras, full viscous

c = -2 Spalart-Almaras, full viscous with DES

c = 3 Spalart-Almaras with wall function

c = -3 Spalart-Almaras with wall function with DES

c = 6 K-Epsilon turbulence model

c = -6 K-Epsilon turbulence model with DES

c = -60 K-Epsilon turbulence model with Multiscale capability

c = 7 K-Epsilon turbulence model, Carlson's formulation

c = 8 SST turbulence model

c = -8 SST turbulence model with DES capability

c = -80 SST turbulence model with Multiscale capability

c itrip - flag to activate flow tripping

c = 0 no tripping

c = 1 activates tripping (useful mostly for k-e turbulence model

c when used for low Reynolds number flows)

c = -1 enforces laminar flow region (useful mostly for SA and

c SST model when a certain flow region has to be forced

c to a laminar flow state. These two models generally 

c provide a fully turbulent flow on the body). 

c NOTE: for itrip value other than 0, user will need to

c (1) define trip patches in the existing grid 

c using PREDISC program (contact: Richard Campbell,

c LaRC)

c (2) create a file project.trip and list all the 

c trip patches and appropriate trip values.

c Refer to the format of project.trip file on

c the USM3D website. You will need to create 

c project.trip file for itrip = 0 as well, if

c idiagnos = 3 (above).

c ev_lim = Eigenvalue limiter value

c

c NOTE: iorder automatic switch should be used only for steady-state cases. 

c If you want to force first order computations based on

c. previous second-order solution upon restart, you must 

c. set iorder = 1. Negative iorder value will not let you

c. restart the solution in first order mode.

c

c-----------------------------------------------------------------------

c LINE 11
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

70 of 78

NESC Request No.: 07-037-I

c

c ncyc = STEADY-STATE PROBLEMS: number of time step cycles

c = UNSTEADY PROBLEMS: number of subiterations to "ntstep"

c ngine = number of jet engines (maximum of 4) pmh

c (additional input required, see LINE 9a below)

c nsinkbc = number of special intakes with specified mass flux

c required bc type in .mapbc: 111,211,311,411.

c (additional input required, see LINE 9b below)

c nrotor = number of rotors/propellers

c (additional input required, see LINE 9c below)

c compF&M = 0 do not include "comp.FandM" file; no special integration

c = 1 integrate component F&M's over patches prescribed in

c file "comp.FandM' at the end of a run. write on to tet.out

c = -1 integrate component F&M's over patches prescribed in

c file "comp.FandM' continuosly during run. write the

c CL,CD,CDv,CM, and 6 components of body forces and 

c moment history on "project.FMhistCOMPN" file where, 

c N stands for component index. You will get as many files

c as the number of components in "comp.FandM" file. 

c Also, write the last iteration F&M information on to

c tet.out file 

c = 2 automatically treat every existing patch in the grid as a 

c component and write forces and moments on to tet.out

c at the end of a run. "comp.FandM" file not read/needed.

c = -2 automatically treat every existing patch in the grid as a 

c component and write history of forces and moments in 

c the seprate files. You will get as many files as the 

c number of active patches in the grid. Also, write 

c last iteration F&M information on to tet.out. "comp.FandM"

c file not read/needed.

c = 3 merger of features activated by compF&M=1 and compF&M=2

c = -3 merger of features activated by compF&M=-1 and compF&M=-2

c = -12 merger of features activated by compF&M=-1 and compF&M=2

c = -13 merger of features activated by compF&M=-1 and compF&M=3

c p_bc1002 = non-dimensional plenum pressure for the BC 1002

c <= 0.0 Use freestream pressure (1/gamma = 0.714)

c > 0.0 Use the specified value

c cldes = prescribed CL for matching lift coefficient via feedback loop

c (cldes=0. means no feedback and maintain prescribed alpha)

c

c-----------------------------------------------------------------------

c LINE 11a

c****************************************************************

c*** NOTE: These lines are included ONLY IF ngine.ne.0 ***

c****************************************************************

c

c ibctype = boundary condition flag

c Second jet or fan exhaust: bctype = n*100+3 (optional)

c Primary jet exhaust by: bctype = n*100+2

c where "n" is the engine number

cvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv

c RULE: If second jet or fan exhaust is utilized, it's

c characteristics must be prescribed before those 

c for the primary jet b.c.

c^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

71 of 78

NESC Request No.: 07-037-I

c fuel = fuel fraction - in % of total exhaust mass flow

c gammaj = ratio of specific heats for (hot) jet

c pjet = static nozzle pressure (usually set to p_infinity)

c p0jet = stagnation pressure of jet

c Rratio = spec. gas constant of air/spec. gas constant of jet

c T0jet = stagnation temperature of jet

c vec(i) = directional cosines of exhaust jet (i=1,2,3 for x,y,z)

c engsim = 0 uniform exhaust properties such as p0jet,t0jet

c = 1 radially varying exhaust properties 

c will need radial distribution of total pressure, 

c total tempreature, tangential and radial velocities. 

c nstation = number of radial locations at which exhaust properties

c specified (only useful for engsim=1). LIMIT: nstation <=20.

c idirec = direction of tangential velocity (only useful for engsim=1)

c = 1 clockwise swirl on engine exhaust looking downstream

c = -1 counter-clockwise swirl on engine exhaust looking downstream

c-----

c NPR = nozzle pressure ratio, p0jet/pjet

c TTR = total temperature ratio, T0jet/T0inf

c RMACHJ = exit mach number of jet

c

c-------------------------------------------------------------------------

c LINE 11b

c sinkid = sink bc identification

c sink 1(BC 111), sink 2(BC 211), sink 3(BC 311)...

c #of sink bcs specified earlier in the input as nsinkbc,

c nsinkbc <= 4

c outer_R = outer radius of sink

c cen-X

c cen-Y

c cen-Z = sink BC's center coordinates

c massflx = mass flow rate through the sink (unit ???)

c -> nondimensional_massflx = dimensional{mdot/(rhoinf*ainf)}

c omega = ratio of fan-tip velocity to freestream velocity

c positive for clockwise looking downstream

c negative for counter-clockwise looking downstream

c

c------------------------------------------------------------------------

c LINE 11c

c rotid = rotor identification

c rotor 1 - BC 501,502

c rotor 2 - BC 601,602

c rotor 3 - BC 701,702

c rotor 4 - BC 801,802

c inner_R = inner radius of the rotor

c outer_R = outer radius of the rotor

c cen-X,cen-Y,cen-Z = rotor center coordinates

c Advratio = advance ratio (free stream velocity/rotor tip speed)

c thrst_c = thrust coefficient of rotor

c torq_c = torque coefficient of rotor

c irotat = direction of rotation (viewing from front/top)

c = 1 clockwise rotation

c = -1 counterclockwise rotation

c rotsim = 0 rotor bc with uniform jump in pressure and/or swirl

c = 1 rotor bc with pressure and/or swirl jump as f(r)
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

72 of 78

NESC Request No.: 07-037-I

c-----------------------------------------------------------------------

c LINE 12a

c..CAUTION: These two lines of input read/needed ONLY IF imvgrd != 0 in LINE 4

c..CAUTION: Current implementation set up to allow only 1 DOF for rotation.

c

c xrcg,yrcg,zrcg - origin of the rotation axis

c reflen - reference length, used for calculating reduced

c frequency, based on grid coordinates

c iroll - grid rotation about x-axis (rotation axis || x-axis)

c = 1 for activating rolling motion

c = 0 for NOT activating rolling motion

c ipitch - grid rotation about y-axis (rotation axis || y-axis)

c = 1 for activating pitching motion

c = 0 for NOT activating pitching motion

c iyaw - grid rotation about z-axis (rotation axis || z-axis)

c = 1 for activating yawing motion

c = 0 for NOT activating yawing motion

c angmean - mean angle of rotation (in degrees)

c angampl - amplitude of rotation (in degrees)

c rfreq - reduced frequency of oscillations (rfreq=k/Lchar) where

c k=2*pi*f*Lchar/(2*V_inf), where f is frequency, Lchar is

c characteristic length, where all quantities are in

c dimensional form

c - The number of "steps" required to cycle through one

c 2*pi period of motion can be computed via.

c Nperiod = pi*Lchar/(k*xmach*delta_t)

c ideformg - grid deformation flag, important only for itimeacc=2.

c = 0 for solid-body rotation of grid,

c value of 0 will switch off Geometric Conservation Law.

c = 1 for grid deformation to accommodate body motion,

c imvMOMCEN = 0 do not rotate moment refernce center with the grid

c = 1 rotate moment refernce center with the grid 

c (recommended)

c

cNOTE

c..angle for pitch (positive angle for nose up)

c..angle for roll (positive angle for right wing up)

c..angle psi for (positive angle for left turn, wind blowing in pilot's

c..left ear)

c-----------------------------------------------------------------------

NOTE: VGRID generated proj.cogsg file will now be protected. (Aug 31, 2005)

c-----------------------------------------------------------------------

c LINE 13a

c..CAUTION: These three lines of input read/needed ONLY IF ivisc = 6,7,or,8

c

c ikeord - spatial order of accuracy for 2-equation 

c turbulence models

c = 1, first order

c = 2, second order

c Recommended value is 1

c icons - Flag to activate either incompressible or compressible

c formulations of for 2-equation turbulence models

c = -1, incompressible

c = 1, compressible, conservative formulation 
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

73 of 78

NESC Request No.: 07-037-I

c Recommended value is 1

c nstagek - Gauss-seidel subiterations for the implicit solution 

c of 2-equation turbulence models

c Recommended value is 10

c t_dtfact - parameter to decrease/increse CFL number for the

c solution of 2-equation turbulence models

c > 1.0 decreases CFL number

c < 1.0 increases CFL number

c Recommended value is 1.0

c t_intsity - freestream turbulence intensity

c Recommended value is 1.0e-3

c mut/mul - ratio of freestream turbulence viscosity to the

c freestream laminar viscosity

c *NOTE*: You can specify either turbulence intensity

c (t_intsity) < 0 or ratio of turbulence-to-laminar 

c viscosity (mut/mul) < 0 to activate USM3D's default 

c values of freestream turbulence intensity 

c [=sqrt(9.e-9)] and mut/mul[=0.009]

c ratiokp - Upper limit on turbulence kinetic energy in relation 

c to static pressure

c Recommended value: 0.0 for first attempt at solution

c If convergence problems are encountered specify

c a value in the range of 0.05-0.1

c dkemax - Maximum allowed update of turbulence model variables

c at any given time step (=Delta(var)/var)

c Recommended value: 0.25

c-----------------------------------------------------------------------

c LINE 13b

c inl - Integer to activate different formulations of two

c equation turbulence models

c 0 = Basic turbulence stress model, assumes that 

c stresses are linearly related to the strain

c 1 = Nonlinear model of Shih, Zhu, and Lumley

c Ref. AIAA-95-2246

c 7 = Nonlinear model of Girimaji 

c Ref. NASA CR-198243, 1995

c *NOTE*: Currently, only k-eps model has provision to

c activate different formulations. For SST model, only

c one value (inl=0) is allowed and user input will be 

c overridden if an inconsistent value is specified for 

c this model

c ilhg - Near-wall modification of k-eps turbulence model

c -14 = Formulation proposed by Launder and Sharma

c -3 = Formulation proposed by Launder and Sharma with

c wall function 

c +3 = No near-wall modification of k-eps turbulence 

c model but use wall function

c Recommended value: -14

c iwallf - Flag to specify wall function grid 

c 0 = grid suitable for integration all the way to wall

c 1 = wall function grid

c Recommended value: 0

c icompCorr - Formulations for correction to the turbulence models 

c related to compressibility at high speeds

c 0 = no correction applied
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

74 of 78

NESC Request No.: 07-037-I

c 1 = Sarkar's formulation, Ref: J. of Fluid Mechanics,

c Vol. 227, 1991

c 2 = Zeman's formulation, Ref: Physics of Fluids, 

c Vol. 2, 1990

c 5 = Wilcox's formulation, Ref: AIAA Paper 92-0439

c Recommended value: 1

c itempCorr - Flag to adjust (for linear k-e model only) eddy 

c viscosity for high-temperature flowfields,

c for example, hot jets.

c 0 = no adjustment made

c 1 = temperature correction model of Hamid

c Ref: J. of Fluids Engg., Vol. 126, No. 4, 2004. 

c itk - Strategies for applying wall boundary condition for

c the turbulence dissipation variable of k-e model

c 0 = isotropic component of dissipation set to zero

c (Launder and Sharma model)

c 1 = zero gradient of dissipation rate at the wall

c (Lam and Bremhorst model)

c 2 = dissipation proportional to turbulence kinetic

c energy grdient at the wall

c Recommended value: 2

c isk - Flag to apply correction from turbulence kintetic

c energy equation to the mean-flow energy equation

c 0 = no modification

c 1 = apply corrections

c Recommended value: 0

c idt_proc - Flag to activate a specific procedure to compute local

c time-steps for steady-state solutions

c 0 = use face-based subroutine (UCTIME_MOD)

c 1 = use node-based subroutine (UCTIME_KE)

c For SA and SST models, idt_proc is reset to 0 

c internally. 

c NOTE: This input variable is introduced in light of

c the observed sensitivity of a CEV solution (Run 314)

c to the differences in the procedure used to compute

c local time steps. This option may be eventually 

c dropped from the export version 6.0.

c-----------------------------------------------------------------------

c LINE 13c

c f1kemax - Maximum allowed value for function f1 in dissipation

c rate equation for Lam Bremhorst model (itk=1) 

c Recommended value : 1.0

c for itk = 0 and 2 this input variable is not used.

Supplemental Description of Input for Real Gas and Chemical Reactions

$SPECIES: flag to activate real-gas model for single or multiple chemical species

The above line is followed by a list of species names involved in this test case. The name of 

each species can be either separated by blank space(s) or comma. An example is shown below:

H2O O2 CO CO2 C7H8 H2 O H OH N2

$REACTIONS: flag to activate chemical reactions
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

75 of 78

NESC Request No.: 07-037-I

The above line is followed by the number of chemical reaction steps (NR). If NR is negative, 

then equilibrium chemistry will be used. If the specified value is greater than 0, then the finiterate chemistry will be used where the value indicates the number of reaction steps. The format 

of finite-rate chemistry model is listed below. As shown in Equation (E.1), a general system of 

chemical reactions can be written in terms of its stoichiometric coefficients:

(

 i j  i j

 and 

)and the i-th chemical species (

i j i j  and 

) for the j-th reaction:

   

  

NS

i

ji i

b j

f j

NS

i

ji i

K

K

1 ,

,

1

 

(E.1)

The source term of i-th species (

i



), due to chemical reaction, can be evaluated as:

























 























 









 

jl jl r

l

l

ns

l

b j

nr

j

r

l

l

ns

l

i i ji f j M

K

M

M K

 

 

1

,

1

1

,

 (E.2)

where 

Mi

is the molecular weight of species “i”, nr is the total number of reactions,

 ji

is the 

stoichiometric coefficient, 

Kf , j Kb, jand

are the forward and backward reaction rates of the j

th

reaction, and ns is the total number of species. 

i j i j

r

 and r

are the power dependence of the 

reactant and product, respectively. The reaction rate, in general, is formulated in the empirical 

Arrhenius form as:

B E R T

f j j

j a j u K A T e,

,





(E.3)

where Aj, Bj, Ea,j are empirical constants for j-th reaction. The input format is shown in the 

following table, where E/R is Ea,j/Ru, STCOEF is 

 ji

, and STCOEG is 

ji r .
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

76 of 78

NESC Request No.: 07-037-I

Format

(2 or 3 

lines*NR)

REACTION: Species names, i = 1  ns (only 1 line as a title)

j, A, B, E/R, ITHIRD, IGLOB,

(STOCEF(i, j), i = 1  ns),

(STOCEG(i, j), i = 1  ns)

*

* needed if IGLOB = 2

Definition j reaction step counter

A reaction rate leading constant

B reaction rate temperature exponent

E/R reaction rate activation energy constant

ITHIRD third-body reaction indicator

= 0: deactivated

= n: for using the n-th species as third-body

= 999: for global (every species) third-body

IGLOB global reaction model indicator

= 0: elementary reactions with the rate of backward 

reaction is calculated from equilibrium constant and 

forward reaction rate

= 1: one-way reaction (either forward or backward 

reaction controlled by the sign of STCOEF); need 

only one input line of STCOEF 

= 2: one-way reaction with power dependency; need input 

line for STCOEF and STCOEG 

STCOEF stoichiometric coefficients of elementary reactions (negative 

signs apply to reactants and positive signs are for the 

products)

STCOEG power dependency coefficients (negative signs apply to 

reactants and positive signs are for the products)

For a sample reaction model:

13

REACTION:H2O, O2, CO, CO2, C7H8, H2, O, H, OH, N2,

 1, 4.4963E9, 1.00, 2.679E4, 0, 2,

 0., -3.5, 7., 0., -1., 4., 0., 0., 0., 0.,

 0., -1.0, 0., 0., -0.5, 0., 0., 0., 0., 0.,

 2, 1.7000E13, 0.00,24070., 0, 0,
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

77 of 78

NESC Request No.: 07-037-I

 0., -1., 0., 0., 0., -1., 0., 0., 2., 0.,

 3, 2.1900E13, 0.00, 2590., 0, 0,

 1., 0., 0., 0., 0., -1., 0., 1., -1., 0.,

 4, 6.0230E12, 0.00, 550., 0, 0,

 1., 0., 0., 0., 0., 0., 1., 0., -2., 0.,

 5, 1.8000E10, 1.00, 4480., 0, 0,

 0., 0., 0., 0., 0., -1., -1., 1., 1., 0.,

 6, 1.2200E17, -0.91, 8369., 0, 0,

 0., -1., 0., 0., 0., 0., 1., -1., 1., 0.,

 7, 4.0000E12, 0.00,10030., 0, 0,

 0., 0., -1., 1., 0., 0., 0., 1., -1., 0.,

 8, 3.0000E12, 0.00,2.50E4, 0, 0,

 0., -1., -1., 1., 0., 0., 1., 0., 0., 0.,

 9, 1.0000E16, 0.00, 0.,999, 0,

 0., 0., 0., 0., 0., 0., -1., -1., 1., 0.,

10, 2.5500E18, -1.00,59390.,999, 0,

 0., 1., 0., 0., 0., 0., -2., 0., 0., 0.,

11, 5.0000E15, 0.00, 0.,999, 0,

 0., 0., 0., 0., 0., 1., 0., -2., 0., 0.,

12, 8.4000E21, -2.00, 0.,999, 0,

 1., 0., 0., 0., 0., 0., 0., -1., -1., 0.,

13, 6.0000E13, 0.00, 0.,999, 0,

 0., 0., -1., 1., 0., 0., -1., 0., 0., 0.,

Reaction Step A B Ea /Ru

Power Dependencies 

(rij)

C7H8 + 3.5O2  7CO + 4H2 4.4963E9 1 2.679E4 [C7H8]

0.5 [O2]1

H2 + O2  OH + OH 1.7E13 0 2.407E4 ---

OH + H2  H2O + H 2.19E13 0 2.59E3 ---

OH + OH  O + H2O 6.023E12 0 5.5E2 ---

O + H2  H + OH 1.8E10 1.0 4.48E3 ---

H + O2  O + OH 1.22E17 -0.91 8.369E3 ---

H + O+ M  OH + M 1.0E16 0 0 ---

O + O + M  O2 + M 2.55E18 -1.0 5.939E4 ---

H + H + M  H2 + M 5.0E15 0 0 ---
NASA Engineering and Safety Center 

Technical Assessment Report

Document #:

NESC-RP07-037

Version:

1.0

Title:

Reacting Multi-Species Gas Capability for USM3D Flow 

Solver

Page #:

78 of 78

NESC Request No.: 07-037-I

H + OH + M  H2O + M 8.4E21 -2.0 0 ---

CO + OH  H + CO2 4.0E12 0 4.03E3 ---

CO + O2  CO2 + O 3.0E12 0 2.5E4 ---

CO + O + M  CO2 + M 6.0E13 0 0 ---

Note: M represents the third body, which includes all species in the above reaction model. If the

power dependency is not specified, then the stoichiometric coefficient will be used.
REPORT DOCUMENTATION PAGE Form Approved

OMB No. 0704-0188

2. REPORT TYPE 

Technical Memorandum

 4. TITLE AND SUBTITLE

Reacting Multi-Species Gas Capability for USM3D Flow Solver

5a. CONTRACT NUMBER

 6. AUTHOR(S)

Frink, Neal T.; Schuster, David M.

 7. PERFORMING ORGANIZATION NAME(S) AND ADDRESS(ES)

NASA Langley Research Center

Hampton, VA 23681-2199

 9. SPONSORING/MONITORING AGENCY NAME(S) AND ADDRESS(ES)

National Aeronautics and Space Administration

Washington, DC 20546-0001

 8. PERFORMING ORGANIZATION

 REPORT NUMBER

L-20147 NESC-RP-07-037

10. SPONSOR/MONITOR'S ACRONYM(S)

NASA

13. SUPPLEMENTARY NOTES

12. DISTRIBUTION/AVAILABILITY STATEMENT

Unclassified - Unlimited

Subject Category 16 Space Transportation and Safety

Availability: NASA CASI (443) 757-5802

19a. NAME OF RESPONSIBLE PERSON

STI Help Desk (email: help@sti.nasa.gov)

14. ABSTRACT

The USM3D Navier-Stokes flow solver contributed heavily to the NASA Constellation Project (CxP) as a highly productive 

computational tool for generating the aerodynamic databases for the Ares I and V launch vehicles and Orion launch abort 

vehicle (LAV). USM3D is currently limited to ideal-gas flows, which are not adequate for modeling the chemistry or 

temperature effects of hot-gas jet flows. This task was initiated to create an efficient implementation of multi-species gas and 

equilibrium chemistry into the USM3D code to improve its predictive capabilities for hot jet impingement effects. The goal of 

this NASA Engineering and Safety Center (NESC) assessment was to implement and validate a simulation capability to handle 

real-gas effects in the USM3D code. This document contains the outcome of the NESC assessment.

15. SUBJECT TERMS

NASA Engineering and Safety Center; USM3D Navier-Stokes flow solver; Launch abort vehicle; Modeling

18. NUMBER

 OF 

 PAGES

83

19b. TELEPHONE NUMBER (Include area code)

(443) 757-5802

a. REPORT

U

c. THIS PAGE

U

b. ABSTRACT

U

17. LIMITATION OF 

 ABSTRACT

UU

Prescribed by ANSI Std. Z39.18

Standard Form 298 (Rev. 8-98)

3. DATES COVERED (From - To)

December 2008 - March 2012

5b. GRANT NUMBER

5c. PROGRAM ELEMENT NUMBER

5d. PROJECT NUMBER

5e. TASK NUMBER

5f. WORK UNIT NUMBER

869021.04.07.01.14

11. SPONSOR/MONITOR'S REPORT

 NUMBER(S)

NASA/TM-2012-217575

16. SECURITY CLASSIFICATION OF:

The public reporting burden for this collection of information is estimated to average 1 hour per response, including the time for reviewing instructions, searching existing data sources, 

gathering and maintaining the data needed, and completing and reviewing the collection of information. Send comments regarding this burden estimate or any other aspect of this 

collection of information, including suggestions for reducing this burden, to Department of Defense, Washington Headquarters Services, Directorate for Information Operations and 

Reports (0704-0188), 1215 Jefferson Davis Highway, Suite 1204, Arlington, VA 22202-4302. Respondents should be aware that notwithstanding any other provision of law, no person 

shall be subject to any penalty for failing to comply with a collection of information if it does not display a currently valid OMB control number.

PLEASE DO NOT RETURN YOUR FORM TO THE ABOVE ADDRESS.

1. REPORT DATE (DD-MM-YYYY)

01- 05 - 2012