Species Net Production Rates of Reactor

[Michael Hegetschweiler]

Hello Cantera Users!

I am trying to understand what the function "getNetProductionRates" (of Kinetics.h) exactly does.

I have a reactor network with just one IdealGasConstPressureReactor without any inlets, outlets or walls. The reactor is filled
with an IdealGasMix.
Then I evolve the reactor network for a certain time dt. After that I would like to get all species net production rates in kmol/m3/s (basically the species concentration change over dt).
First I thought the function "getNetProductionRates" will give me the exactly this information but I realized when I let evolve the
reactor for a long time such that nothing happens anymore then the output of getNetProductionRates is basically zero for all the species.
Therefore I conclude that this function doesn't calculate the integrated species net production rates from the initial to the final reactor state. Am I correct? what else it does calculate?

I figure that what I actually want I can calculate "by hand" by subtracting the initial specie concentrations by the final specie concentration and dived by dt?

Thanks for some clarification!

Best
Michael 

[Bryan W. Weber]

Michael,

That function gives the instantaneous net production rates of all the species. Thus, at equilibrium, all the net production rates are zero (by definition). If you want the integrated net production rate, I'd suggest you store the instantaneous production rate in an array of some sort and use a numerical integration technique to find the integrated production rate (for instance, the trapezoid rule is very simple).

Best,

Bryan

[Michael Hegetschweiler]

Thanks Bryan!

I wrote a small test program which integrates which calculates the the overall net production rates of the species in two different ways. First as you suggested by integrating the instantaneous rate using the function getNetProductionRates and second, by taking the differences between the initial species concentration and the final concentration and dividing by the overall elapsed reactor time.

I see some differences between those two results even for major species as CH4, CO2 or H20 (I use a initially stoichiometric CH4 - air mixture). I suspect that these could be numerical integration errors.
Could you have a look at the program and tell me if things are done right?

Appreciate the help very much!

Michael

[Bryan W. Weber]

Michael,

I don't think your formula for the trapezoidal integration is correct (presuming that's what you're implementing). If the array rates is storing the integrated rates, then I think the formula should be

rates[k] += (tmp_rates0[k] + tmp_rates1[k])*(t1 - t0)/2

Then, at the end of the integration, you don't need to divide everything by the total time. I don't think this will be very close the the "beginning minus end divided by Delta t" approach, because the former is an integration while the latter is a differentiation... What are you actually trying to calculate here?

Best,

Bryan

[Michael Hegetschweiler]

Good Morning Bryan

You are right my implementation of the trapezoidal formula was wrong and it should be as you stated.

What I try to reach is, I need the source term for the species conservation equations (in kmol/m3/s or kg/m3/s). I use OpenFoam as the flow solver and combine it with Cantera. Therefore I need to extract species source terms, enthalpy source term and species thermo and transport properties.
For extracting the species source term I initially thought I could use the function getNetProductionRates of Cantera and wanted to double-check
and found the "inconsistency" I described in my first post.

Best
Michael

[Bryan W. Weber]

Michael,

If you need the source term for the CFD solver, don't you want the instantaneous production rates of the species (in which case getNetProductionRates does exactly what you want)? Why do you need the integrated production rate?

Best,

Bryan

[Michael Hegetschweiler]

So, I create a reactor for each grid cell. Then I let evolve the reactor for prescribed time interval. For that I use the advance(time) of the ReactorNetwork class. Then I call getNetProductionRates and this should give me the the species source term. Is that what you would suggest?

What I suspected is that with the advance() function the solver uses several internal substeps and when I call getNetProductionRates I get just the instantaneous rates of the last substep...

Bet
Michael

[Michael Hegetschweiler]

What I have to mention is that in the original OpenFoam chemical solver the species source term was calculated as species concentrations at the end of the time step - species concentrations at the beginning of the time step divided by the the time step.
That is why I wanted to compare this against the output of getNetProductionRates and found it is not consistent.

Best
Michael

[Ray Speth]

Michael,

the advance() function the solver uses several internal substeps and when I call getNetProductionRates I get just the instantaneous rates

This is precisely what is happening, and I think this is the most straightforward meaning of the functions that Cantera provides. A Phase / Solution / IdealGasMix object represents a single thermodynamic state and has no notion of time, so the only meaning that “net production rates” can have is the instantaneous one, which is reinforced by the term rates and the units associated with this value (kmol/m^3/s). A difference between two discrete states would only have dimensions of kmol/m^3.

For your CFD application, if the intent is to use the ReactorNet model to integrate the chemical source term as part of an operator split integration, then the finite difference approach that your code was using is the correct one. If you are implementing your own integration of the chemical source term, then you would want to use the netProductionRates function (which is what the ReactorNet model uses when it is integrating the equations).

Regards,
Ray

[Michael Hegetschweiler]

Thank you for the clarifications, Ray.

In my case I use the OpenFoam solver for all conservation equations. For the kinetics part I use Cantera, since I believe that Cantera provides a better and more robust kinetic solver. Therefore I need the chemical reaction rates for species and enthalpy to feed into the species and enthalpy conservation equations. The goal is to perform transient flame calculations and therefore I have to resolve the chemical time scales well enough. Because of that I guess I can use the instantaneous net production rate directly and I don't have to solve a reactorNet.

Best,
Michael