gdtk-uq/gdtk

Question on finite rate chem

Opened this issue · 7 comments

All,

I am wondering, at one point in my simulation, the run fails with a constantly decreasing dt (reaching a dt = 1e-25 or lower which is certainly a sign of instability in my solution - or in my choice of settings...) outputting dt_therm = -1 and locations of dt_chem being = -1. Is this indicative of the thermo dt being too low for the chemistry steps to properly update?

What steps would you recommend I follow to trace this error, and would it be possible to place a limiter on one of these two variables, so that the lower limit of the chemistry/thermo is set?

Are there other config settings that I should look into?

Thank you.

Hi,
There are a few dts mentioned here, so I'll go on the assumption that 'dt' relates to the fluid simulation timestep, and the other dts are obvious by their names.

If the fluid timestep is attempting to dive to something tiny, this is usually an indication that something is awry. I would try to identify what parts of the domain are causing that, or sometimes it is a poorly selected boundary condition. Some steps to debug the simulation, either in isolation or simultaneously: remove thermochemical effects; use low order reconstruction. This might at least give a good starting solution for a simulation with more complicated numerics and physics.

dt_chem = -1 and dt_therm = -1 does NOT signal an error. It signals a reset to do an internal estimate of timestep size for those updates. The chemistry update and energy exchange update is typically sub-cycled (has a smaller dt than the flow), but not always. When it sub-cycling mode, it tries to increase the sub-cycled dt by 10% when things are going well. This is a computational efficiency measure. We basically let the timestep grow until we've gone too far. Now when we've gone too far, we try to reduce that timestep by something like 50% and then let the process of 10% growth start repeating again. Occasionally, things are pushed too far. We signal a reset by setting the dt_chem or dt_therm to -1.0.

It is normal to see some instances of dt_chem = -1 and dt_therm = -1 peppered across the domain (you can plot these in Paraview, for example). If they are all clustered in one region, that's probably an indication that the interaction between gas dynamics and thermochemistry is not being well simulated in that region.

For tracing this problem, you need to determine if gas dynamics is triggering thermochemistry problems or vice versa. Start with no thermochemical effects. You can then look at the config variable: reaction_time_delay (see: https://gdtk.uqcloud.net/docs/eilmer/eilmer-reference-manual/#_finite_rate_thermo_chemistry )

A limiter on dt_chem and dt_thermo is not the issue. As mentioned, -1 is a valid signal. If those values are going very small, there are usually other issues at play.

What type of chemistry is this? Is the flow field going very cold? If so, there is a control to turn off chemistry below a certain static temperature (see: config.T_frozen at same link as above).

Some information on my simulation; the total run time is very short ~15e-5 sec, with the code exiting at around 2.5e-5. This is the exiting time with reactions enabled from the beginning of the simulation. I will try now to run from the beginning with the reaction delay enabled, perhaps until 3.0e-5 to see if it does/does not fail at 2.5e-5s, to give us an answer on whether its the reactions or the gasdynamics causing the exit. I am currently using the Gupta 5sp air model.

The flow field is somewhat cold, perhaps around 150 K is the lowest, with ~3000 K being the highest. Perhaps it would be useful to ensure chemical reactions are not occurring below the temperature at which the lowest energy bond is expected to break... Thank you for the suggestions. I will try this as well and return to this comment chain.

Setting a T_frozen value did not change the outcome this time around,... Having the same error in this case, could just be a need for a constant time step - do you think trying a fixed dt that is low enough would aid in our debugging? The case runs well with a lower refinement, so perhaps this is a CFL related exit; so I am thinking about setting fixed_time_step = true and setting that dt_init to the same value at which the lower refinement case became stable. I'll update here if that is a workaround; please let me know your thoughts.

Please let me know your thoughts and thank you for your comments and help!

I'd probably not use the fixed time step. If the case runs on a coarser mesh, you might want to just try a lower CFL on the finer one instead.

Also a quick question: Are you getting a solution file written when the solver fails? You should be able to look at that to get some idea of what the problem is.

You can also try changing the errTol in the chemistry file, one of mine looks like this:

Config{
   tempLimits = {lower=200.0, upper=30000.0},
   odeStep = {method='rkf', errTol=1.000000e-09},
   tightTempCoupling = false,
}

What type of chemistry is this? endothermic or exothermic?

When I saw Nick's suggestion, I was also reminded of the tightTempCoupling option for problematic chemistry schemes.

In general, I think you need to find out where and when the problem is occurring. Then look for flow field/grid problems there.

You mentioned low resolution works ok, so this could indicate a grid issue. Have you tried warm-starting the high resolution calculation with the low resolution calculation? This can help as flow features are mostly in place.

It is an endothermic reaction with N2 and O2 in the inflow. Currently the tightTempCoupling option is set to true, are you suggesting having that set to true or to false? At the time of the error, I do not see any outstanding features, but that could mean I need a smaller dt_plot increment to be able to capture it. My settings are:

species = {[0]='N2', [1]='O2', [2]='N', [3]='O', [4]='NO'}
config = {
tempLimits = {lower=300.000000, upper=30000.000000},
odeStep = {method='rkf', errTol=1.000000e-05},
tightTempCoupling = true,
maxSubcycles = 10000,
maxAttempts = 4
}

I will reduce the errorTol as well and keep this thread updated.

Endothermic chemistry tends to alter the flow field in a local sense less dramatically than exothermic (combustion). So, I do not think tightTempCoupling will help.

In fact, I don't think playing with chemistry settings will help.

I think a step back and look at what features of the fine grid are different from the coarse grid and where the failure is occurring will be your best bet. This sounds like a flow-field/grid/chemistry interaction and you'll need a more holistic approach to sort it out than just the chemistry settings.