Cantera / cantera

Chemical kinetics, thermodynamics, and transport tool suite
https://cantera.org
Other
580 stars 341 forks source link

Equilibrium Constant of Irreversible Reactions is Non-Zero (and Reverse Reaction Rate) #1691

Closed jPurdue2 closed 3 weeks ago

jPurdue2 commented 2 months ago

Hello,

I am running windows 10, python 3.9, and observed this error with Cantera 2.6 and Cantera 3.0.

In Cantera I created a chemical mechanism with irreversible reactions using " => " instead of " <=> ". In Cantera's documentation it says that the equilibrium constant (Kc) of the irreversible reactions should be zero. However, when I displayed the Kc values they were non-zero. In addition, the reverse reaction rates were non-zero (which for an irreversible reaction they should be zero).

I observed this with the mechanism I created when using .cti and .yaml files for Cantera 2.6, and .yaml files for Cantera 3.0. I also observed this when running Cantera 2.6 though MATLAB (with both .cti and .yaml file types). See attached image for from Cantera 3.0.

In MATLAB I used "equil_Kc(gas)" to get the equilibrium constants. In Python I used gas.equilibrium_constants.

image

jPurdue2 commented 2 months ago

According to: https://cantera.org/documentation/docs-3.0/sphinx/html/matlab/kinetics.html the equilibrium constant for the function equil_Kc: "Returns a column vector of the equilibrium constants for all reactions. The vector has an entry for every reaction, whether reversible or not, but non-zero values occur only for the reversible reactions. If the output is not assigned to a variable, a bar graph is produced instead."

Running Cantera 2.6 in MATLAB does not result in non-zero equilibrium constants. I will try with the latest version of Cantera (3.0) and see if the result is different.

jPurdue2 commented 2 months ago

The issue is present in Cantera 3.0. Here is my gas: image

Here are the reactions: image

And the equilibrium constants from equil_Kc (after taking log10) image

As well as the output plot from equil_Kc image

jPurdue2 commented 2 months ago

After running some simulations it seems that the irreversible reaction is one directional (I simulated all O with only the forward reaction O2 + M => O + O +M) and there was no recombination of O unless I included the reverse reaction.

However, the non-zero rate constants and equilibrium constants are a little odd when using an irreversible reaction...

speth commented 2 months ago

While the Matlab equil_Kc function claims that it should return zero for reactions marked as irreversible, there is no similar specification for the Python Kinetics.equilibrium_constants property or the C++ BulkKinetics::getEquilibriumConstants method.

There is an ambiguity with "irreversible" reactions between what are otherwise two equivalent definitions for the equilibrium constant: it can either be defined based on the standard Gibbs free energy of the reaction, or it can be defined as the ratio of the forward and reverse rate constants. The first definition corresponds to how it is calculated in Cantera, and the documentation of BulkKinetics::getEquilibriumConstants. The second definition would give a result of infinity (not zero). Given the direction of calculations in Cantera, I think sticking to the first definition, as currently implemented, is preferable.

Regarding the reverse rate constants, the Python module correctly returns 0 for irreversible reactions. In the Matlab toolbox (both the now-removed legacy toolbox and the new/experimental toolbox), it returns non-zero values because it changes the default value of the argument to C++ BulkKinetics::getRevRateConstants that toggled how to handle irreversible reactions; the C++ default of doIrreversible = false seems reasonable to me.

In summary, I think the necessary fixes for this are just a couple of small changes to the Matlab toolbox: