Unexpected zeroing out of Q_REAC

[mcgratta]

The case Species/hrrpuv_reac_soot.fds has periods during which the Q_REAC terms are zero.

[rmcdermo]

This should be fixed with PR #13837 . Here is the result I get...

hrrpuv_reac_soot.pdf

[rmcdermo]

@cxp484 @drjfloyd This case will not work with ODE_SOLVER='CVODE'. Looks like there are a couple of issues. One is whether we can mix a fast reaction with an Arrhenius reaction and use CVODE for the finite-rate part. The next issue is that we set Q_REAC_SUB = 0._EB in CVODE section, so HRRPUV REAC will definitely break.

I understand that this was not a priority when getting CVODE off the ground. The question is whether it is something we want to pursue. If we do not do this for CVODE then I don't see how we can get rid of RK2_RICHARDSON. Well, perhaps we get rid of it and just live with EXPLICIT EULER for cases like this. Something to think about.

[drjfloyd]

I think for CVODE plus fast we would just have to pass in the fast reaction as another Arrhenius reaction. Guess try and come up with some A and E that is fast but not make CVODE unstable fast.

CVODE integrates all the reactions over the timestep and we just get the end state for the species from which we get total HRRPUV by the net enthalpy change. With reversible reactions, competing reactions making the same product, etc., I don't know if we would be able to easily untangle a per reaction HRRPUV for some really complex scheme. Also each reaction would need its own Q_REAC array. So Smooke with 16 species would need 50 Q_REAC slots (assuming you want to know both forward and reversible).

[rmcdermo]

I agree, was thinking the same thing later for the "fast" reaction mixed with finite-rate.

If CVODE is working now, I think we can get rid of RK2_RICHARDSON. The soot case works the same with EXPLICIT EULER. We can throw an error if someone tries to output HRRPUV REAC with CVODE.

Before purging RK2_RICHARDSON, I think setting the hrrpuv_reac_soot case up with FE and CVODE (sans hrrpuv reac) would be a good test.

[marcosvanella]

I understand CVODE needs some work to run with REACs where N_S is defined.

[rmcdermo]

OK. Also, CVODE does not do the mixing step. I guess let's get those working before we consider removing RK2.

[drjfloyd]

@marcosvanella that might be just in read processing. I don't know if there is anything in chem that presumes N_S are nice integer values.

[marcosvanella]

Ok, maybe @cxp484 can chip in on this. I noted Species/jones_lindstedt.fds was giving substantially different equilibrium results using RK2 and CVODE as ODE_SOLVER, and discussing with Chandan this came up as a potential reason.

[cxp484]

@marcosvanella, @drjfloyd, and @rmcdermo, I have not yet tested CVODE with two-step finite-rate chemistry. My understanding is that in this approach, the exponents in the Arrhenius reaction rate can be fractional, unlike detailed chemistry where they are typically stoichiometric coefficients. A quick code check suggests that these exponents are stored in the same N_S variables (as Jason pointed out). I need to test this with CVODE to see where things are failing.

Overall, I think next steps to make CVODE work like RK2 would be:

1. Enable support for 2-3 step reaction mechanisms with fractional exponents.
2. Incorporate mixing during CVODE integration substeps. This likely involves adding additional terms to the derivatives and Jacobian to account for mixing. I need to sit and do the math on this. Once this will be done, we can probably do LES with detailed chemical mechanism.

I am going to work on this. It is in my next to do list.

[drjfloyd]

If we had a reaction equation like d[B]/dt = A [C]^c [D]^d e^(-E/RT) where c and d were fractional, I think the logic in chem should work as long as the correct values are present for c and d. In read though, to make life simpler when first starting down this path, detailed chemistry is a bit seperate from the EULER and RK2 for some of the input processsing. Need to make sure we are correctly assigning N_S if it is given in the input file when CVODE is called. There might be the issue that we are replacing the REAC input N_S values with stoichiometric values when CVODE is being used.

[rmcdermo]

I am going to close this issue since the hrrpuv_reac_soot seems to have worked. We can continue the CVODE conversation offline or start a new discussion.