The following equations are used by CFPOST for deriving a requested variable that is not present in the solution file. If a variable cannot be derived by one of the following equations or by the simple multiplication or division by ρ, then the requested variable must exist within the solution file. Note that multiple derivations may be used to derive a variable. For example, when deriving the total temperature T_{0}, the temperature T and Mach number M may themselves be derived.
The derivation of variables that include the gas constant R or the ratio of specific heats γ will take real gas effects into account if sufficient information is available. For frozen and finite-rate chemistry flows, species properties are read from a chemistry data (.chm) file specified with the chemistry command. Otherwise, real gas effects are modeled by multiplying R_{∞} by the compressibility factor Z, and using the effective specific heat ratio β in place of γ.
Note that some variables require the specification of a surface, or normal direction into the flow field, and others require knowledge of "up" and "down" axes orientation. Variables denoted with (s) indicate that the normal direction must be specified properly using the subset or surface command. Variables denoted with (o) are affected by the axes orientation set by the orientation command.
Several variables may be computed in multiple ways, depending on the information available in the solution file. The equations are shown below in the order that CFPOST attempts to use them. If any of the needed quantities are missing from the solution file and cannot be derived, CFPOST will try using the subsequent equation, etc. Note though, that for frozen and finite-rate chemistry flows, if the necessary information isn't available CFPOST will abort with an error message, rather than use the perfect gas equations.
rho - Density
ρ | = | p / (R_{∞}TZ) | |
= | p / (R_{∞}T) |
Z contains the local real gas effects, R/R_{∞}.
rho0 - Stagnation density
For frozen and finite-rate chemistry,
ρ_{0} | = | p_{0} / (R_{∞}T_{0}Z) |
For a perfect gas,
ρ_{0} | = | ρ [ 1 + (γ_{∞}−1) M^{ 2}/2 ] ^{1 / (&gamma∞−1)} |
p - Pressure
p | = | ρR_{∞}TZ | |
= | ρR_{∞}T | ||
= | (β−1) (ρe_{0}− ρV^{ 2}/2− ρh_{f}) | ||
= | (γ_{∞}−1) (ρe_{0}− ρV^{ 2}/2− ρh_{f}) |
Z contains the local real gas effects, R/R_{∞}.
p0 - Stagnation pressure
For frozen and finite-rate chemistry,
p_{0} | = | p exp{(1/R_{mix}) ∑^{ns}_{i=1} C_{i} (R_{u} / MW_{i}) [Sfun_{i,T0} − Sfun_{i,T}]} |
where Sfun is a function used to compute the entropy for a species as a function of temperature, from the curve fit data in the .chm file.
For a perfect gas,
p_{0} | = | p [ 1 + (γ_{∞}−1) M^{ 2}/2 ] ^{&gamma∞ / (&gamma∞−1)} |
For a value in a rotating reference frame (i.e., p0r), the total temperature T_{0} and Mach number M in the above equations are replaced by the corresponding values in the rotating reference frame.
Cp - Pressure coefficient
C_{p} | = | (p−p_{∞}) / q_{∞} |
where q_{∞} = ρ_{∞}(M_{∞}a_{∞})^{2}/2
deltap - Delta pressure
Δp | = | p − p_{∞} |
pp - Pitot pressure
For M ≤ 1,
p_{p} | = | p [ 1 + (γ_{∞}−1) M^{ 2}/2 ] ^{&gamma∞ / (&gamma∞−1)} |
For M > 1,
p_{p} | = | p [ (γ_{∞}+1) M^{ 2}/2 ] ^{&gamma∞ / (&gamma∞−1)} / [ 2γ_{∞} M^{ 2}/(γ_{∞}+1) − (γ_{∞}−1)/(γ_{∞}+1) ] ^{1 / (γ∞−1)} |
Cpt - Total pressure coefficient
C_{pt} | = | (p_{0} − p_{∞}) / q_{∞} |
where q_{∞} = ρ_{∞}(M_{∞}a_{∞})^{2}/2
q - Dynamic pressure
q | = | ρV^{ 2}/2 |
T - Temperature
T | = | p / (ρR_{∞}Z) | |
= | p / (ρR_{∞}) |
Z contains the local real gas effects, R/R_{∞}.
T0 - Stagnation temperature
For frozen and finite-rate chemistry, iterate on T_{0} until
h_{mix,T0} | = | h_{0,mix} |
∑^{ns}_{i=1} C_{i} (R_{u} / MW_{i}) Hfun_{i,T0} | = | ∑^{ns}_{i=1} C_{i} (R_{u} / MW_{i}) Hfun_{i,T} + V^{ 2}/2 |
where Hfun is a function used to compute the enthalpy for a species as a function of temperature, from the curve fit data in the .chm file.
For a perfect gas,
T_{0} | = | [ 1 + (γ_{∞}−1) M^{ 2}/2 ] |
For a value in a rotating reference frame (i.e., T0r), the values of T_{0}, h_{0,mix}, V^{ 2}, and M in the above equations are replaced by the corresponding values in the rotating reference frame.
h - Sensible enthalpy
h | = | βe_{i} | |
= | γ_{∞}e_{i} |
h0 - Sensible stagnation enthalpy
h_{0} | = | βe_{i} + e_{k} | |
= | γ_{∞}e_{i} + e_{k} |
ei - Sensible internal energy per unit volume
e_{i} | = | p / [ρ (β−1)] | |
= | p / [ρ (γ_{∞}−1)] |
ek - Kinetic energy per unit volume
e_{k} | = | V^{ 2}/2 |
e0 - Absolute stagnation energy per unit volume
e_{0} | = | (ρe_{0}) / ρ | |
= | [ ρ(e_{0})_{r} + ρ(ω × r)^{2}/2 − ρ(g ⋅ r)/2 ] / ρ | ||
= | e_{i} + e_{k} + h_{f} |
Ma - Equivalent isentropic Mach number
M_{a} | = | { [2/(β−1)] [ [1 + (β−1)M^{ 2}_{∞} / 2] [1 + βC_{p}M^{ 2}_{∞} / 2]^{(1−β)/β} − 1 ] }^{1/2} | |
= | { [2/(γ_{∞}−1)] [ [1 + (γ_{∞}−1)M^{ 2}_{∞} / 2] [1 + γ_{∞}C_{p}M^{ 2}_{∞} / 2]^{(1−γ∞)/γ∞} − 1 ] }^{1/2} |
where C_{p} is the pressure coefficient.
u - x velocity
In a non-rotating reference frame,
u | = | (ρu)/ρ | |
= | u_{r} + (ω × r) ⋅ e_{x} | ||
= | 0 |
In a rotating reference frame,
u_{r} | = | (ρu_{r})/ρ | |
= | 0 |
v - y velocity
In a non-rotating reference frame,
v | = | (ρv)/ρ | |
= | v_{r} + (ω × r) ⋅ e_{y} | ||
= | 0 |
In a rotating reference frame,
v_{r} | = | (ρv_{r})/ρ | |
= | 0 |
w - z velocity
In a non-rotating reference frame,
w | = | (ρw)/ρ | |
= | w_{r} + (ω × r) ⋅ e_{z} | ||
= | 0 |
In a rotating reference frame,
w_{r} | = | (ρw_{r})/ρ | |
= | 0 |
V - Velocity magnitude
V | = | (V^{ 2})^{1/2} |
For a value in a rotating reference frame (i.e., Vr), V^{ 2} in the above equation is replaced by the corresponding value in the rotating reference frame.
The value of V^{ 2} in a non-rotating reference frame is computed using one of the following equations. (As noted earlier, the equation used depends on the information available in the solution file. If any of the needed quantities in the first equation listed below are missing from the solution file and cannot be derived, CFPOST will try using the subsequent equation, etc.)
V^{ 2} | = | |V^{2}| | |
= | u^{2} + v^{2} + w^{2} | ||
= | (ρ |V|)^{2} / ρ^{2} = [(ρu)^{2} + (ρv)^{2} + (ρw)^{2}] / ρ^{2} | ||
= | (|V_{r}| + |ω × r|)^{2} | ||
= | (ρ |V_{r}| + ρ |ω × r|)^{2} / ρ^{2} |
In a rotating reference frame, V_{r}^{2} is computed using one of the following equations.
V_{r}^{ 2} | = | |V_{r}|^{2} | |
= | u_{r}^{2} + v_{r}^{2} + w_{r}^{2} | ||
= | (ρ |V_{r}|)^{2} / ρ^{2} = [(ρu_{r})^{2} + (ρv_{r})^{2} + (ρw_{r})^{2}] / ρ^{2} |
Vxy, Yyz, Vxz - Crossflow velocity magnitudes
V_{xy} | = | (u^{2} + v^{2})^{1/2} | |
V_{yz} | = | (v^{2} + w^{2})^{1/2} | |
V_{xz} | = | (u^{2} + w^{2})^{1/2} |
M - Mach number
M | = | V / a |
For a value in a rotating reference frame (i.e., Mr), V in the above equation is replaced by the corresponding value in the rotating reference frame.
a - Speed of sound
a | = | (γ_{∞}R_{∞}TZ)^{1/2} | |
= | (γ_{∞}R_{∞}T)^{1/2} |
s - Entropy
For frozen and finite-rate chemistry,
s − s_{∞} | = | ∑^{ns}_{i=1} C_{i} (R_{u} / MW_{i}) [Sfun_{i,T} − Sfun_{i,T∞}] − R_{mix} ln (p/p_{∞}) |
where Sfun is a function used to compute the entropy for a species as a function of temperature, from the curve fit data in the .chm file.
For a perfect gas,
s − s_{∞} | = | (c_{v})_{∞} ln [(p/p_{∞}) / (ρ/ρ_{∞})^{γ∞}] |
where (c_{v})_{∞} = R_{∞}/(γ_{∞}−1). Here s_{∞} may be assumed to be zero.
S - Strain rate magnitude
S | = | (2 S_{ij} S_{ij})^{1/2} |
Sxx, Sxy, Sxz, Syx, Syy, Syz, Syz, Szx, Szy, Szz - Strain rate components
S_{ij} | = | (∂u_{i}/∂x_{j} + ∂u_{j}/∂x_{i} )/2 | |
W - Rotation rate magnitude or vorticity magnitude
Ω | = | (2 W_{ij} W_{ij})^{1/2} | = | (ω_{x}^{2} + ω_{y}^{2} + ω_{z}^{2})^{1/2} |
The magnitude of the rotation rate tensor is equivalent to the vorticity magnitude.
Wxx, Wxy, Wxz, Wyx, Wyy, Wyz, Wyz, Wzx, Wzy, Wzz - Rotation rate components
W_{ij} | = | (∂u_{i}/∂x_{j} − ∂u_{j}/∂x_{i} )/2 | |
wx, wy, wz - Vorticity components
ω_{x} | = | ∂w/∂y − ∂v/∂z | |
ω_{y} | = | ∂u/∂z − ∂w/∂x | |
ω_{z} | = | ∂v/∂x − ∂u/∂y |
helicity - Helicity
helicity | = | Ω ⋅ V = ω_{x}u + ω_{y}v + ω_{z}w |
swirl - Swirl
swirl | = | Ω ⋅ V / (ρV^{ 2}) = (ω_{x}u + ω_{y}v + ω_{z}w) / (ρV^{ 2}) |
dila - Dilatation
dila | = | ∇ ⋅ V = ∂u/∂x + ∂v/∂y + ∂w/∂z |
localpha - Local α (o)
Let v_{ds} be the velocity in the "downstream" direction and v_{up} be the velocity in the "upward" direction as defined by the aerodynamic axes set in the grid file or specified with the orientation command. Local α is then defined as:
α | = | arctan (v_{up} / v_{ds}) |
locbeta - Local β (o)
Let v_{ss} be the velocity in the "sidestream" direction as defined by the aerodynamic axes set in the grid file or specified with the orientation command and let V be the magnitude of the velocity vector. Local β is then defined as:
β | = | arcsin (v_{ss} / V) |
mu - Total viscosity coefficient
μ | = | μ_{l} + μ_{t} |
mul - Laminar viscosity coefficient
For frozen and finite-rate chemistry,
μ_{l} | = | ∑^{ns}_{i=1} [(X_{i} μ_{i}) / (∑^{ns}_{j=1} X_{j} φ_{ij})] |
where
μ_{i} | = | μ_{0} (T/T_{μ})^{3/2} [(T_{μ} + S_{μ}) / (T + S_{μ})] | |
φ_{ij} | = | [ 1 + (μ_{i} / μ_{j})^{1/2} (MW_{j} / MW_{i})^{1/4} ]^{2} / ( 8 + 8MW_{i} / MW_{j} )^{1/2} |
The molecular weight of each species is obtained from the thermodynamic data in the chemistry (.chm) file, and the values of μ_{0}, T_{μ}, and S_{μ} for each species are obtained from the transport data in the .chm file.
For a perfect gas,
μ_{i} | = | μ_{0} (T/T_{0})^{3/2} [(T_{0} + S_{1}) / (T + S_{1})] |
where [Pratt and Whitney Aircraft Company "Aeronautical Vestpocket Handbook," 19th. ed., Part Number 79500]:
μ_{0} | = | 0.37445×10^{−6} slug/ft-sec =
1.79287×10^{−5} kg/m-sec | |
T_{0} | = | 518.67 °R = 288.15 K | |
S_{0} | = | 216.0 °R = 120 K |
kappa - Thermal conductivity coefficient
For frozen and finite-rate chemistry,
κ | = | ∑^{ns}_{i=1} [(X_{i} κ_{i}) / (∑^{ns}_{j=1} X_{j} φ_{ij})] |
where
κ | = | κ_{0} (T/T_{κ})^{3/2} [(T_{κ} + S_{κ}) / (T + S_{κ})] | |
μ_{i} | = | μ_{0} (T/T_{μ})^{3/2} [(T_{μ} + S_{μ}) / (T + S_{μ})] | |
φ_{ij} | = | [ 1 + (μ_{i} / μ_{j})^{1/2} (MW_{j} / MW_{i})^{1/4} ]^{2} / ( 8 + 8MW_{i} / MW_{j} )^{1/2} |
The molecular weight of each species is obtained from the thermodynamic data in the chemistry (.chm) file, and the values of κ_{0}, T_{κ}, S_{κ}, μ_{0}, T_{μ}, and S_{μ} for each species are obtained from the transport data in the .chm file.
For a perfect gas,
κ | = | μ_{l} c_{p} / P_{r} |
delta - Boundary layer thickness (s)
The edge of the boundary layer is determined according to the criteria established by the boundary layer edge command.
delta1 or delta* - Boundary layer displacement thickness (s)
δ_{1} | = | ∫_{0}^{δ} [ 1 − (ρV) / (ρ_{δ}V_{δ}) ] ds = ∫_{0}^{δ} ds − 1 / (ρ_{δ}V_{δ}) ∫_{0}^{δ} ρV ds |
where V is the total velocity magnitude.
Starting at the wall and advancing one point at a time along the grid line intersecting the wall, the integrals in the second equation are accumulated. At each grid point, the values of the subscripted variables are assumed to be the value at the current grid point. The value of the equation is then computed and compared to the value at the previous iteration. If the change meets the criteria specified by the boundary layer edge command, the iteration stops and the thickness will be set to the accumulated arc length along the grid line (i.e., the value of the first integral). (Algorithm compliments of R. H. Bush.)
delta2 or THETA - Boundary layer momentum thickness (s)
δ_{2} | = | ∫_{0}^{δ} (ρV) / (ρ_{δ}V_{δ}) ( 1 − V / V_{δ} ) ds = 1 / (ρ_{δ}V_{δ}) ∫_{0}^{δ} ρV ds − 1 / (ρ_{δ}V_{δ}^{2}) ∫_{0}^{δ} ρV^{ 2} ds |
The method described for the boundary layer displacement thickness delta1 is used to evaluate the integral.
delta3 - Boundary layer energy thickness (s)
δ_{3} | = | ∫_{0}^{δ} (ρV) / (ρ_{δ}V_{δ}) ( 1 − V^{ 2} / V_{δ}^{ 2} ) ds = 1 / (ρ_{δ}V_{δ}) ∫_{0}^{δ} ρV ds − 1 / (ρ_{δ}V_{δ}^{3}) ∫_{0}^{δ} ρV^{ 3} ds |
The method described for the boundary layer displacement thickness delta1 is used to evaluate the integral.
deltau - Boundary layer velocity thickness (s)
δ_{u} | = | ∫_{0}^{δ} [ 1 − V / V_{δ} ] ds = ∫_{0}^{δ} ds − 1 / V_{δ} ∫_{0}^{δ} V ds |
The method described for the boundary layer displacement thickness delta1 is used to evaluate the integral. Note that the velocity thickness is the same as the incompressible displacement thickness.
THETAinc - Boundary layer incompressible momentum thickness (s)
θ_{inc} | = | ∫_{0}^{δ} V / V_{δ} [ 1 − V / V_{δ} ] ds = ∫_{0}^{δ} V / V_{δ} ds − 1 / V_{δ}^{2} ∫_{0}^{δ} V^{2} ds |
The method described for the boundary layer displacement thickness delta1 is used to evaluate the integral. The incompressible momentum thickness is often used to compute an incompressible shape factor for evaluating flow control devices.
Q - Heat transfer rate (s)
Q | = | κ ∂T/∂n = κ ∇T ⋅ n = κ (∂T/∂x ∂x/∂n + ∂T/∂y ∂y/∂n + ∂T/∂z ∂z/∂n) |
The equation is evaluated using the same process described in the derivation of tau.
tau, taux, tauy, tauz - Total (laminar+turbulent) shear stress (s)
τ | = | μ ∂V_{tan} / ∂n |
where μ is the total (laminar+turbulent) viscosity
μ | = | μ_{l} + μ_{t} |
and V_{tan} is the tangential velocity vector determined as follows:
V_{tan} | = | V − (V ⋅ n) n |
The magnitude of the shear stress is computed as
τ | = | μ [(∂u_{tan} / ∂n)^{2} + (∂v_{tan} / ∂n)^{2} + (∂w_{tan} / ∂n)^{2}]^{1/2} |
where the partial derivatives are evaluated using the chain rule:
∂f / ∂n | = | (∂f / ∂x) (∂x / ∂n) + (∂f / ∂y) (∂y / ∂n) + (∂f / ∂z) (∂z / ∂n) |
The partial derivatives are computed using the equations presented in the section on "Derivatives of Variables With Respect to Cartesian Coordinates".
The shear stress acts in the direction of the velocity, so the components are:
τ_{x} | = | τ u_{tan} / V_{t} | |
τ_{y} | = | τ v_{tan} / V_{t} | |
τ_{z} | = | τ w_{tan} / V_{t} |
If wall function boundary conditions were used in the flow solver, CFPOST will also use wall functions when computing the shear stress.
Cf, Cfx, Cfy, Cfz - Skin friction coefficient (s)
The total skin friction coefficient is the shear stress magnitude normalized by the freestream dynamic pressure.
C_{f} | = | τ / (ρ_{∞}V_{∞}^{2}/2) |
Similarly, the Cartesian components of the skin friction coefficient are defined as
(C_{f} )_{x} | = | τ_{x} / (ρ_{∞}V_{∞}^{2}/2) | |
(C_{f} )_{y} | = | τ_{y} / (ρ_{∞}V_{∞}^{2}/2) | |
(C_{f} )_{z} | = | τ_{z} / (ρ_{∞}V_{∞}^{2}/2) |
Note that the freestream values used are the reference values specified with the Wind-US FREESTREAM keyword and may differ substantially from those at the local boundary layer edge.
uu, uv, uw, vu, vv, vw, wu, wv, ww - Turbulent stress components (s)
u_{i}"u_{j}" | = | {ρu_{i}"u_{j}"} / ρ | = | −[2μ_{t}( S_{ij} − S_{kk}δ_{ij}/3 ) − 2/3*ρkδ_{ij}] |
Turbulent stresses can be computed for any of the structured grid turbulence models. The formula above represents the Boussinesq approximation used by most models. For the Spalart-Allmaras turbulence model, the turbulent kinetic energy is estimated by means of a generalized form of Bradshaw's relation:
ρk | = | μ_{t} (2S_{ij}S_{ij})^{1/2} / 0.31 |
y+ - Non-dimensional boundary layer wall coordinate (s)
y^{+} | = | y u_{τ} / ν_{wall} |
where
u_{τ} | = | (τ_{wall} / ρ_{wall})^{1/2} |
ν_{wall} | = | μ_{wall} / ρ_{wall} |
and here y is the distance from the wall, measured along the grid line intersecting the wall.
u+ - Non-dimensional boundary layer velocity (s)
u^{+} | = | V / u_{τ} |
where u_{τ} = (τ_{wall} / ρ_{wall})^{1/2}.
k+ - Non-dimensional turbulent kinetic energy (s)
k^{+} | = | k / u_{τ}^{2} |
where u_{τ} = (τ_{wall} / ρ_{wall})^{1/2}.
epsilon+ - Non-dimensional turbulent dissipation rate (s)
ε^{+} | = | ε μ_{l} / (ρ u_{τ}^{4}) |
where u_{τ} = (τ_{wall} / ρ_{wall})^{1/2}.
uu+, uv+, uw+, vu+, vv+, vw+, wu+, wv+, ww+ - Non-dimensional turbulent stress components (s)
(u_{i}"u_{j}")^{+} | = | {ρu_{i}"u_{j}"} / (ρu_{τ}^{2}) |
where u_{τ} = (τ_{wall} / ρ_{wall})^{1/2}.
Mt - Turbulent Mach number
Mt | = | (2 k / a^{2})^{1/2} |
Rt - Turbulent Reynolds number
[Ladd, J. A., and Kral, L. D. (1992) "Development and Application of a Zonal k-ε Turbulence Model for Complex 3-D Flowfields," AIAA/SAE/ASME/ASEE 28th. Joint Propulsion Conference and Exhibit, AIAA Paper 92-3176, p. 3]
Rt | = | ρk^{2} / (μ_{l} ε) |
Ry - Turbulent Reynolds number based on y (s)
[Ladd, J. A., and Kral, L. D. (1992) "Development and Application of a Zonal k-ε Turbulence Model for Complex 3-D Flowfields," AIAA/SAE/ASME/ASEE 28th. Joint Propulsion Conference and Exhibit, AIAA Paper 92-3176, p. 3]
Ry | = | yρk^{1/2}/μ_{l} |
where here y is the distance from the wall, measured along the grid line intersecting the wall.
fmuJL, fmuCH, fmuSP, fmuLB - Turbulence model damping functions (s)
[Ladd, J. A., and Kral, L. D. (1992) "Development and Application of a Zonal k-ε Turbulence Model for Complex 3-D Flowfields," AIAA/SAE/ASME/ASEE 28th. Joint Propulsion Conference and Exhibit, AIAA Paper 92-3176, p. 3]
(f_{μ})_{JL} | = | e^{−2.5/(1 + Rt / 50)} | |
(f_{μ})_{CH} | = | 1 − e^{−0.0115y+} | |
(f_{μ})_{SP} | = | (1 + 3.45 / Rt^{1/2}) tanh (y^{+}/70) | |
(f_{μ})_{LB} | = | (1 + 20.5 / Rt) (1 − e^{−0.0165Ry})^{2} |
C-species - Mass fraction of species
Using O2 as an example for species,
C-O2 | = | ρ_{O2} / ρ |
X-species - Mole fraction of species
Using O2 as an example for species,
X-O2 | = | (C-O2/MW_{O2}) ∑^{ns}_{i=1} MW_{i} |
Veff - Effective collision frequency
V_{eff} | = | C (p/T) T_{e}^{1/2} |
where C is a function of the gas, p is the pressure in torr (i.e., millimeters of mercury), T is the temperature in K, and T_{e} is the electron temperature. In CFPOST, C is assumed to be 7.3×10^{9}, the value for air, and T_{e} is assumed equal to T.
MW - Molecular weight
For frozen and finite-rate chemistry,
MW | = | 1 / (∑^{ns}_{i=1} C_{i} / MW_{i}) |
where C_{i} is the mass fraction of species i. The molecular weights for all the species may be listed using the show chemistry command.
For a perfect gas,
MW | = | R_{u} / R_{∞} |
R - Gas constant
For frozen and finite-rate chemistry,
R | = | R_{u} / MW |
For a perfect gas,
R | = | R_{∞} |
cp - Specific heat at constant pressure
For frozen and finite-rate chemistry,
c_{p} | = | ∑^{ns}_{i=1} C_{i} (R_{u} / MW_{i}) Cpfun_{i,T} |
where Cpfun is a function used to compute the specific heat at constant pressure for a species as a function of temperature, from the curve fit data in the .chm file.
For a perfect gas,
c_{p} | = | γ_{∞}R_{∞} / (γ_{∞} − 1) |
cv - Specific heat at constant volume
c_{v} | = | c_{p} − R |
gamma - Ratio of specific heats
For frozen and finite-rate chemistry,
γ | = | c_{p} / c_{v} |
For a perfect gas,
γ | = | γ_{∞} |
hf - Heat of formation
For frozen and finite-rate chemistry,
h_{f} | = | ∑^{ns}_{i=1} C_{i} [h_{f(298)}]_{i} |
For a perfect gas,
h_{f} | = | 0 |
The functions Cpfun, Hfun, and Sfun appearing in some of the above equations for frozen and finite-rate chemistry flows are used to compute thermodynamic properties for a species as a function of temperature. Properties are computed using curve fit equations, with the coefficients supplied in a chemistry data (.chm) file. The .chm file to be used must be specified using the CFPOST chemistry command.
Three different (but related) formats are available for the type of curve fits being specified in the .chm file - SPARKCRV, WINDNASA, and NASA3287.
Cpfun - Molar specific heat at constant pressure
Function Cpfun computes the specific heat per mole at constant pressure for a species. For the SPARKCRV and WINDNASA formats, the curve fit formula is
C_{p} / R_{u} | = | a_{1} + a_{2}T + a_{3}T^{2} + a_{4}T^{3} + a_{5}T^{4} |
For the NASA3287 format, the formula is
C_{p} / R_{u} | = | a_{1}T^{−2} + a_{2}T^{−1} + a_{3} + a_{4}T + a_{5}T^{2} + a_{6}T^{3} + a_{7}T^{4} + a_{8}T^{5} |
Note that the value returned by Cpfun is dimensionless.
Hfun - Molar absolute enthalpy
Function Hfun computes the enthalpy per mole for a species. For the SPARKCRV and WINDNASA formats, the curve fit formula is
H / R_{u} | = | a_{1}T + a_{2}T^{2}/2 + a_{3}T^{3}/3 + a_{4}T^{4}/4 + a_{5}T^{5}/5 + b_{1} |
For the NASA3287 format, the formula is
H / R_{u} | = | −a_{1}T^{−1} + a_{2} ln T + a_{3}T + a_{4}T^{2}/2 + a_{5}T^{3}/3 + a_{6}T^{4}/4 + a_{7}T^{5}/5 + a_{8}T^{6}/6 + b_{1} |
For the SPARKCRV format, the enthalpy reference state is at T = 0 K. I.e., the curve fit for H/R_{u} computes
H / R_{u} | = | (H_{f(0)} + ΔH_{(0→T)}) / R_{u} |
where H_{f(0)} is the heat of formation at 0 K, and ΔH_{(0→T)} is the change in enthalpy between 0 K and T.
However, for the WINDNASA and NASA3287 formats the enthalpy reference state is at T = 298.15 K. I.e., for these formats the curve fit computes
H' / R_{u} | = | (H_{f(298)} + ΔH_{(298→T)}) / R_{u} |
where H_{f(298)} is the heat of formation at 298.15 K, and ΔH_{(298→T)} is the change in enthalpy between 298.15 K and T.
Thus, for the WINDNASA and NASA3287 formats, function Hfun adjusts the value to shift the reference state to 0 K. I.e., the actual returned value is
H / R_{u} | = | H' / R_{u} + H_{(shift)} / R_{u} | |
= | (H_{f(298)} + ΔH_{(298→T)}) / R_{u} + (H_{f(0)} − H_{f(298)} + ΔH_{(0→298)}) / R_{u} | ||
= | (H_{f(0)} + ΔH_{(0→T)}) / R_{u} |
Note that the value returned by Hfun has units of K.
Sfun - Molar entropy
Function Sfun computes the entropy per mole for a species. For the SPARKCRV and WINDNASA formats, the curve fit formula is
S / R_{u} | = | a_{1} ln T + a_{2}T + a_{3}T^{2}/2 + a_{4}T^{3}/3 + a_{5}T^{4}/4 + b_{2} |
For the NASA3287 format, the formula is
S / R_{u} | = | −a_{1}T^{−2}/2 − a_{2}T^{−1} + a_{3} ln T + a_{4}T + a_{5}T^{2}/2 + a_{6}T^{3}/3 + a_{7}T^{4}/4 + a_{8}T^{5}/5 + b_{2} |
It should be noted that the entropy is defined by the equation
dS | = | C_{p} dT/T − R_{u} dp/p |
Thus
S / R_{u} | = | ∫ (1/R_{u}) dS | |
= | ∫ (C_{p} / R_{u}) dT/T − ∫ dp/p |
and Sfun computes the temperature portion of the above equation, ∫ (C_{p} / R_{u}) dT/T.
Note also that the value returned by Sfun is dimensionless.
The area to be integrated is considered to be a collection of individual polygons. A surface may be defined directly by a subset or surface command. If the dimensions of the surface is i points by j points, there will be (i−1) × (j−1) polygons that define the surface. The perimeter of the surface is described by the line segments that join adjacent grid points at the extreme bounds of the subset. A surface may also be defined as the intersection of a cutting plane defined by a cut command with a computational volume defined by a subset command. The polygons created by passing the cutting plane through the cells that make up the volume will be those processed by the integration.
For surfaces defined directly, each polygon will have exactly four verticies defined by the grid points that make up the polygon. The area will be determined by computing the magnitude of the cross products of the diagonals. For surfaces defined by a cutting plane, each polygon will have from three to six verticies defined by the intersection of the cutting plane with the cell edges. The area will be determined by decomposing the polygon into triangles and adding their individual areas. The unit normal vector for each polygon will be defined to be the unit vector in the direction of the vector sum of the normal vectors at the vertices.
For any cell, the value of a particular property for that cell (pressure, temperature, density, etc.) will be defined to be the average of the values at the verticies for that cell; i.e., for a cell with N_{v} vertices,
f_{cell} | = | (1/N_{v}) ∑ f_{i} |
where the summation is from i = 1 to N_{v}.
Averages over a group of N_{c} cells are computed as shown below. In these formulas the summations are from i = 1 to N_{c}, and the variables being summed are cell averages as defined above.
Area Average
f_{area} | = | (∑ f_{i} A_{i}) / (∑ A_{i}) |
Mass Average
f_{mass} | = | (∑ ρ_{i} f_{i} A_{i}) / (∑ ρ_{i} A_{i}) |
Mass Flux Average
f_{massflux} | = | [ ∑ ρ_{i} f_{i} (V_{i} ⋅ A_{i}) ] / [ ∑ ρ_{i} (V_{i} ⋅ A_{i}) ] |
Arithmetic Mean
f | = | (1/N_{c}) ∑ f_{i} |
σ | = | { [ 1/(N_{c}−1) ] ∑ (f_{i} − f )^{2} }^{1/2} | |
= | { [ 1/(N_{c}−1) ] [ (∑ f_{i}^{2}) − N_{c}f^{ 2} ] }^{1/2} |
Mass Flux
m | = | ρ (V⋅n) A |
Momentum Flux
F_{momentum} | = | ρV (V⋅n) A |
Pressure Flux
F_{pressure} | = | pA n |
Total Thrust
F_{total} | = | F_{momentum} + F_{pressure} |
Gross Thrust
F_{G} | = | (p − p_{∞}) A n |
F | = | F_{pressure} + F_{viscous} |
Computing partial derivatives of variables with respect to Cartesian coordinates requires the calculation of the metric coefficients. Using subscripts to denote partial differentiation, the Jacobian is defined as:
J = x_{ξ} (y_{η}z_{ζ} − y_{ζ}z_{η}) + x_{η} (y_{ζ}z_{ξ} − y_{ξ}z_{ζ}) + x_{ζ} (y_{ξ}z_{η} − y_{η}z_{ξ})
The metric coefficients are then computed as follows:
ξ_{x} =
J^{ −1}
(y_{η}z_{ζ} −
y_{ζ}z_{η})
ξ_{y} =
J^{ −1}
(x_{ζ}z_{η} −
x_{η}z_{ζ})
ξ_{z} =
J^{ −1}
(x_{η}y_{ζ} −
x_{ζ}y_{η})
η_{x} =
J^{ −1}
(y_{ζ}z_{ξ} −
y_{ξ}z_{ζ})
η_{y} =
J^{ −1}
(x_{ξ}z_{ζ} −
x_{ζ}z_{ξ})
η_{z} =
J^{ −1}
(x_{ζ}y_{ξ} −
x_{ξ}y_{ζ})
ζ_{x} =
J^{ −1}
(y_{ξ}z_{η} −
y_{η}z_{ξ})
ζ_{y} =
J^{ −1}
(x_{η}z_{ξ} −
x_{ξ}z_{η})
ζ_{z} =
J^{ −1}
(x_{ξ}y_{η} −
x_{η}y_{ξ})
Given the above, it is then possible to compute the partial derivative of a variable with respect to a Cartesian coordinate. Using the chain rule:
f_{x} = f_{ξ} ξ_{x} + f_{η} η_{x} + f_{ζ} ζ_{x}
It is often necessary to know the vector normal to computational surface. First we define vectors tangent to a line in which only one of the computational coordinates vary:
V_{ξ} =
x_{ξ} i +
y_{ξ} j +
z_{ξ} k
V_{η} =
x_{η} i +
y_{η} j +
z_{η} k
V_{ζ} =
x_{ζ} i +
y_{ζ} j +
z_{ζ} k
The vector normal to a constant computational surface is then the cross product of the two vectors whose computational coordinates vary on the surface:
N_{ξ} =
V_{η} × V_{ζ}
N_{η} =
V_{ζ} × V_{ξ}
N_{ζ} =
V_{ξ} × V_{η}
Given a normal vector, one can then simply form a unit normal vector and directional derivatives.
n = N / |N| = (∂x/∂n) i + (∂y/∂n) j + (∂z/∂n) k
Derivatives of values with respect to computational coordinates may be either first or second order accurate depending on the need. Where normal vectors are not of a concern, second order central differences are used except where a boundary is involved and then first order forward or backward differences are used:
∂f/∂ξ ≈ δf =
(f_{i+1} −
f_{i−1}) / 2
Central difference
∂f/∂ξ ≈ Δf =
f_{i+1} − f_{i}
Forward difference
∂f/∂ξ ≈ ∇f =
f_{i} − f_{i−1}
Backward difference
For normal derivatives, first or second order may be selected by the first order or second order command or parameter on the integrate force command. Because of the internal structure of CFPOST, these derivatives are always forward derivatives:
∂f/∂ξ ≈
f_{i+1} − f_{i}
First order
∂f/∂ξ ≈
(−f_{i+2} +
4f_{i+1} −
3f_{i}) / 4
Second order
A low-pressure region is defined as a range of rake probes on a particular ring whose total pressures fall below the average total pressure for that ring.
Circumferential Intensity
Intensity = (ΔPC / P)_{Ring i} = (PAV_{Ring i} − PAVLOW_{Ring i}) / PAV_{Ring i}
where
PAV_{Ring i} = Ring-averaged total pressure
PAVLOW_{Ring i} = Ring-averaged total pressure over
all low-pressure regions
Circumferential Extent
Extent = Cumulative angular extent of all low-pressure regions, in degrees
Multiple-per-revolution (MPR)
MPR_{Ring i} = ∑^{Q}_{k=1} [ (ΔPC / P)_{i} k θ_{i}k^{−} ] / max [ (ΔPC / P)_{i} k θ_{i}k^{−} ]
where
Q = Number of low-pressure regions
&theta_{i}k^{−} =
Extent of low-pressure region k on the i'th ring
Radial Intensity
Intensity = (ΔPR / P)_{Ring i} = (PFAV − PAV_{Ring i}) / PFAV
where
PFAV = Area-weighted face-average total pressure
Circumferential Intensity
Intensity = (PAV_{Ring i} − PMIN_{Ring i}) / PFAV
where
PAV_{Ring i} = Area-averaged total pressure of ring i
PMIN_{Ring i} = Minimum total pressure of ring i
PAV = Area-averaged total pressure over the complete face
Circumferential Extent
The circumferential extent is the same as the SAE definition.
Radial Intensity
The radial intensity is the same as the SAE definition.
Shape Factor
Shape Factor = (GE ring circumferential intensity) / (SAE ring circumferential intensity)
Last updated 2 June 2016