Finite Rate Chemistry Setup for Hypersonic Flow
 

I've tried to provide a thorough description of my problem without overloading it with detail. If anyone has specific questions I'd be glad to provide as much info as I can.

For the past several months I've been putting together a solver setup for Fluent (12.0.16) to handle external hypersonic flow. The application problems will be 2D airfoil sections, 3D wings, and possibly simplified 3D vehicles.

I've gotten reasonable results for most of the validation cases I've run, at least as far as flow features: shock types and angles look right, expansion fans, boundary layer formation, etc show up where they're supposed to. The general flowfield looks good. However, the magnitude of surface profiles are way off, such as temperature and pressure distributions being off by an order of magnitude or more.

The biggest suspect to me is the chemistry. I looked into including finite-rate reactions early on, but had a lot of trouble getting them to run properly, so I dropped them. The target simulations will be run at fairly low ambient pressures (<1,000 Pa), and toward the lower hypersonic range (M=5.0-10.0). I was hoping that the errors would be acceptable, but that doesn't seem to be the case.

So I've turned back to including the finite-rate chemistry. I can converge to a reasonable solution with species transport, everything except reactions, but the moment I turn on the reactions Fluent fails. I was hoping someone here would be able to point me in the right direction or find a flaw in my chemistry model.

The general solver and model settings I'm using are:

• Density-Based Solver, Steady-State
• Energy Equation On
• k-w SST Turbulence Model
• Discrete Ordinate Radiation Model
• Species Transport, Volumetric Reactions, Eddy-Dissipation Concept Model
• Aungier-Redlich-Kwong Real Gas Model for Density
• Pressure Inlet/Outlet BCs, Wall BCs
• Boundary Turbulence Specification: I = 0.05%, Turbulent Viscosity Ratio = 1.0
• Calculating in Absolute Pressure (0 Pa Operating Pressure)
• AUSM Gradient Method
• All equations second order
• Solution Steering - Hypersonic, No FMG Initialization

If I left anything out, it's probably left as default, but feel free to ask.

I've poured over the Fluent docs to make sure I'm not missing some incompatibility between these models, and worked through all the tutorials I can find that relate. I have had trouble finding a good tutorial for setting up a custom chemistry model, or a solid hypersonic tutorial, although I have seen them referenced in several places.

The chemistry model I'm using is a 5-species, 5-reaction model:

Species Transport, Volumetric Reactions, Eddy Dissipation Concept
Third Body Efficiencies (for dissociation reactions), Calculate Backward Rates

Species:
N2, O2, N, O, NO

Reactions and Arrhenius Parameters:
(A = Pre-Exponential Factor, beta = Temperature Exponent, Ea = Activation Energy)

Oxygen Dissociation:

• O2 + M' -> O + O + M'
• A = 2.70 e+19, beta = -1.00, Ea = 494,000 J/kgmol
• Rate Exponent: [O2] = 2.00

Nitrogen Dissociation:

• N2 + M'' -> N + N + M''
• A = 3.70 e+21, beta = -1.60, Ea = 941,000 J/kgmol
• Rate Exponent: [N2] = 2.00

NO Dissociation:

• NO + M''' -> N + O + M'''
• A = 2.90 e+15, beta = 0.00, Ea = 621,000 J/kgmol
• Rate Exponent: [NO] = 2.00

N2-O Exchange:

• N2 + O -> NO + N
• A = 1.82 e+14, beta = 0.00, Ea = 319,000 J/kgmol
• Rate Exponents: [N2] = 1.00, [O] = 1.00

NO-O Exchange:

• NO + O -> N + O2
• A = 3.80 e+9, beta = 0.00, Ea = 173,100 J/kgmol
• Rate Exponents: [NO] = 1.00, [O] = 1.00

Third Body Efficiencies:
M' vs:

• O2 = 1.00
• N2 = 0.10
• O = 2.80
• N = 0.10
• NO = 0.10

M'' vs:

• O2 = 0.10
• N2 = 1.00
• O = 0.10
• N = 2.80
• NO = 0.10

M''' vs:

• O2 = 0.05
• N2 = 0.05
• O = 1.00
• N = 1.00
• NO = 1.00

Boundary Species Distribution:

• NO - 0.0
• N - 0.0
• O - 0.0
• O2 - 0.209476
• N2 - Bulk Species

I can attach screenshots of exactly how I've set these up in Fluent if it will help. I pulled the parameters and model from a paper I can dig up if someone wants. I pulled the reaction orders from the NIST Kinetics Database.

My grids vary; some are structured, some are unstructured. All use near-wall placement of cells to capture boundary layer effects (Y+ < 10). Some examples of the cases I've been running are:

• Diamond airfoils with various wedge angles
• Cylinder cases
• Ramp in wind tunnel
• Step in wind tunnel
• etc...

I can provide more details on the specific validation cases I'm running if requested. For now I'm focusing on the ramp case, which can be found here: http://www.grc.nasa.gov/WWW/wind/val...p/hypramp.html

The solution residuals will converge to 1e-6 without much trouble. General flow field and, where applicable, drag monitors also indicate convergence. When reactions are turned on, I receive the following error:

dasac failure at temperature = ...

I've seen this error pop up many times on these forums, but it never seems to get a solid reply. Hopefully I've provided enough information here for someone to point me in the right direction. This site has been a fantastic resource for sorting out problems, and I appreciate any help you can provide.

I thought I would update this for anyone else having a similar problem.

The root of the failures seems to have been the ISAT algorithm. I eventually disabled it and relied instead on direct integration, after which my simulations stopped failing. I also did a few studies on the effects of various chemistry settings (for the Eddy Dissipation Concept model), and settled on the following for my cases:

Flow Iterations/Update: 5
Volume Fraction Constant: 0.01
Time Scale Constant: 0.9
Thermal Diffusion: No
Multicomponent Diffusion: No
Stiff Chemistry Solver (from console): No

I also found that enabling backward reactions in the chemistry setup had a large impact on the solution chemistry.

After some more digging through databases, I revised the chemistry model a bit:

O2 -> O- + O-
A: 1.14 e+16 (cm^3/mol/s
beta: -1.00
Ea: 494,000 (J/mol)
Order: 2

N2 -> N + N
A: 3.48 e+18
beta: -1.70
Ea: 940,000
Order: 2

NO -> N + O
A: 1.45 e+15
beta: 0.00
Ea: 620,000
Order: 2

N2 + O -> NO + N
A: 2.41 e+14
beta: -0.125
Ea: 323,000
Order: 2

O2 + N -> NO + O
A: 1.81 e+14
beta: 1.00
Ea: 27,190
Order: 2

Obviously this tweaking is all problem dependent, but the key seemed to be disabling ISAT. I thought I would share since this was an especially frustrating roadblock.

What is difference between running simulation without enabling ISAT and running ISAT with DIRECT INTEGRATION ENABLED ? What is the role of aggrasiveness factor?

It's been a while since I've looked at this in detail, so if anyone spots an error feel free to correct me.

When solving for finite-rate chemistry using the eddy-dissipation concept (EDC) model, the reactions have to be integrated. This is done using direct integration (DI), ie the equations for each cell are numerically integrated. This integration is computationally expensive.

ISAT (In-Situ Adaptive Tabulation) is an algorithm which reduces the number of integrations that have to be performed. It does this via linear interpolation, and only performs direct integration when the interpolated result has too much error.

So the reaction equations are being integrated in both cases, but enabling ISAT reduces the number of total integrations required, thus speeding up the simulation. By disabling ISAT, I mean that the simulation is performing direct integrations all the time, no interpolation.

The aggressiveness factor is part of the EDC model itself; it looks like a control for the numerics of the method (robustness vs convergence). My guess would be that it is something akin to a under/over-relaxation factor, but I couldn't find any specific details on how it's integrated into the model.

Again, I had a great deal of trouble getting these simulations to run and I don't think I ever got them working quite as well as I'd like. I was reporting my results here mostly to help anyone else who stumbled into the same problem.

You should also dig into the Fluent Theory Manual and User Guide, which both discuss using ISAT (or not).

I couldnt find any information on how the aggresiveness factor is factored into the solver also but from... https://www.sharcnet.ca/Software/Flu...ies_Model.html
it states it "is a numerical factor which controls the robustness and the convergence speed. This value ranges between 0 and 1, where 0 (the default) is the most robust, but results in the slowest convergence."

AdeeMij,
I am doing similar simulations to your setup. I have a few comments and questions that I would like to ask.

1. Are you getting risidual convergence or are the results just unacceptable? If not getting residual convergence, try running it using 1st order accuracy first, then 2nd order accuracy. Or mess with courant # and relaxation values. If the results are unacceptable, specifically, what values are unacceptable?

2. Have you tried a simulation with ideal gas instead of the real gas models? I'm not sure if that would change much, but the conditions that you state are valid with the ideal gas model. from....
https://www.sharcnet.ca/Software/Flu...ug/node334.htm
you can see the limits that the ideal gas model is acceptable for.

3. I would try running it without the radiation model. Again, I am not sure how this will change the results, but a simpler simulation in terms of numerical calculations might help you see where the errors are coming from. Hypersonic speeds becomes radiatively dominated approximately around 12 km/s (Mach 35+). Turning it off for M= 10 simulations will not reduce accuracy.

4. How come you change your volume fraction constant? When i look up this constant's meaning, everywhere says its 2.13777 which is its default value and shouldn't vary.

5. Did you adjust your solution control -> limits? maximum static T limit is 5000 K which is kind of low. Not sure if your simulations should exceed that Temp.

6. Did you adjust your thermal conductivity, viscosity, and mass diffusivity? The default values that are given I believe is air at 300K. They should not be constant, since they depend on the momentum, energy and mass of a particle as it travels through the gas. Anderson in "Hypersonic and High Temp Gas Dynamics" uses kinetic theory to solve for these 3 values. Mass diffusivity in FLUENT has a kinetic theory option, but the others don't. They have mixing-law options which might be a better approximation. I am still reading into which would be the best way to approach this and I will keep you updated if I find anything useful.

7. Just a minor comment, the species B.C.'s you state are for the mole fractions. The default in FLUENT is mass fractions, so make sure the "specify species in mole fractions" is checked.


Keep me updated!

Panboy

Also,

you state your y+ is < 10. How close are they to 1? I am aiming mine close to 1 as thats what literature has taught me to accurately capture boundary layer.

Panboy

AdeeMij,

Also, I would carefully check your Arrhenius parameters.
check this website for a compilation of experimental values.
http://www.nist.gov/data/nsrds/NSRDS-NBS67.pdf
The first values you set up were from the article "Different Levels of Air Dissociation Chemistry & Its Coupling with Flow Models" I assume, as I ran into that article and they exactly match your first inputs. The second inputs I am not sure how you are getting those values.
However, if you check the link above, you can get another set of data. They closely match the data from the article (your first inputs).

For anyone who is using reaction in their simulation, its a great resource for rate constants.


Panboy

Third Body efficiencies
 

@ AdeeMij

In the model you have suggested,
for the following reaction :

O2 -> O- + O-
A: 1.14 e+16 (cm^3/mol/s
beta: -1.00
Ea: 494,000 (J/mol)
Order: 2

What would be the third body efficiencies?

Thanks in advance :)

gaurav_prabhudesai....

If you PM me your email address I can send you an article that states the Dunn-Kang and Park Model coefficients with 3rd body efficiencies.

Panboy

hey,

If anyone here is still active on the Hypersonics research, I am simulating a 3-D re-entry capsule for boundary and shock capturing on CFX. I am using the total energy fluid model and a Mach no. of 32. I have defined a no-slip wall, static pressure and velocity at the far field, static pressure at outlet, and the material is air as an ideal gas. Can someone advise me as to how I should go about setting up the reactions, the wall definition(if anything other than no-slip is necessary) and the thermal flux conditions(the monte carlo method is causing divergence)?

I tried running the simulation(steady) with the mentioned conditions but the solution runs to about 300 iterations with constantly increasing Mach Number going to the order of 10^6 and then, as would be expected, fails. Could someone kindly help me out?

thanks

I would suggest looking into the contour plots of Mach and seeing where on the geometry/mesh this divergence is initially occurring. Let me know what you find.

Also, the Mach number is a very high level value. Can you look into what else is diverging that is causing the Mach number to diverge? for example, temperature? Also if possible, plot species reaction variables and see if those values make sense.

Panboy

Hello,

I'm trying to simulate dissociated air at 5380K flowing over a non-symmetric wedge using fluent. I used Park's 5-species 5-reaction model in my setup. For cp of each individual species, I used my own piecewise-polynomial coefficients. For the cp of the mixture, I used the piecewise-polynomial coefficients for air at 300K to 20,000K. When I ran this case, my total temperature at the inlet, with static tempt = 5380K, is about 2000K lower. I also used mixing law for cp of the mixture, but that resulted in gamma about 1.4, which is incorrect. Could anyone give me some thoughts?

Thanks

Thanhho,

Are you sure it has converged? Are you using time dependent solver or no? What is your time step? Have you tried to change your time steps and see if you get same results?

Panboy,

The solution converged and the solver is set to steady state. With piecewise-polynomial coefficients entered for the 5 species, I got " Incorrect Specific Heat Value" was a warning message in the console window when I used mixing-law for Cp of the mixture.

This is very interesting. Is there anything odd looking about the specific heat value when you look at a contour plot of it? What values do you get for the specific heat value near the wedge and freestream.

One thing I would make sure is to double check your specific heat inputs for each species. A high number of significant digits is needed and any rounding can cause a large change in the Cp.

Let me know what you find out about the Cp

I curve fit the data using matlab with 4 decimal places for the polynomial coefficients. the contour of cp is attached and the values are very low. I used mixing law for cp in the mixture, and the cp for each individual species is set to polynomial with coefficients from matlab. The gamma in the freestream is ~1.40 when in fact it should be lower. I've attached the cp plot. Please let me know if you have any other thoughts. I pretty much tried everything I could think of and still the gamma is at ~1.40. Thank you.

In order to truly understand how FLUENT is getting 1.4, I would do the hand calculation for the freestream myself. It seems to be an input error on your part since all FLUENT does in plug and chug the values you give it using the equations you tell it to. What cp do you get when you do the freestream cp?

If you know the cp of each species as well, you can check to see which species is the one that is incorrect.

If you don't mind, can you send me your case and data file? I would like to look at it. I'll PM you

Hello,

I am using the above mentioned reactions for my re-entry capsule studies in CFX. Can someone please throw some light on the importance of the third body efficiency and the backward reaction rates?
The method I am using to incorporate the reactions is as follows:
-Define a the reactions with the activation energy, Arrhenius parameter, etc.
-Define a new reacting mixture using the reactions.
-Use this 'reacting mixture' as the fluid model.

Thanks

Third body efficiencies are very important to some, and not too important to others. They might make a difference in your simulation or not at all.

I think it would be best if you described what you know about third body reactions and why you believe you need them and then we can go from there. I say this because the importance of them is dependent on your simulation and what you are looking for to get out of your simulation. I look forward to hearing from you

Thanks Panboy. I am trying to capture the boundary layer transition on a re-entry capsule fro which all effects such as the finite rate chemistry, shock standoff distance, etc, need to be accurately predicted. I am new to modelling chemistry in flows so am not sure as to what the third body efficiency means and am not using them for now.

The problem I am facing when I try to run the solver, seems to be a technical one. When I try to run after defining the reaction, the error I get is that the Arrhenius number does not have the correct dimensions. The problem is that the CFX-pre does not have an option to specify any dimensions for the Arrhenius number.