Reverse reaction rate in ReversibleReaction artificially restricted - might deliver wrong results
Summary
Reversible reactions in OpenFOAM deliver wrong results when the backward/ reverse reaction acts strongly.
Steps to reproduce
Run appended test case and display generated lineplot.png file.
Example case
Below test case features the reacting flow in a channel (0,1 x 0,1 x 1 m³) at T_\mathrm{in}=300\,\mathrm{°C}, U_\mathrm{in}=10\frac{\mathrm m}{\mathrm s}, (x_{\mathrm H_2}=x_{\mathrm H_2\mathrm O}=x_\mathrm{CO}=x_{\mathrm{CO}2})_\mathrm{in}=0.25 and adiabatic free-slip walls considering both water-gas shift as well as methane-steam reforming reaction:
What is the current bug behaviour?
The backward reaction rate is currently limited by limiting the concentration-based equilibrium constant Kc to an (arbitrary) value of 1e-6:
Foam::scalar Foam::ReversibleReaction
<
ReactionType,
ReactionThermo,
ReactionRate
>::kr
(
const scalar kfwd,
const scalar p,
const scalar T,
const scalarField& c
) const
{
return kfwd/max(this->Kc(p, T), 1e-6);
}
The unit of this value differs with depending the stoichiometric coefficients of the reaction educts and products, which further points out the arbitrariness of this value:
K_c=K_{(p)}\left(\frac{p^\ominus}{\bar RT}\right)^{\sum_k{\nu_k''-\nu_k'}}=e^{-\frac{\Delta_{\mathrm R}G^\ominus}{\bar RT}}\left(\frac{p^\ominus}{\bar RT}\right)^{\sum_i{\nu_i''}-\sum_k{\nu_k'}}
[K_c]=\left(\frac{\mathrm{kmol}}{\mathrm m^3}\right)^{\sum_i{\nu_i''}-\sum_k{\nu_k'}}
In contrast to this, the forward reaction coefficient is not restricted at all.
Especially for reactions, where the direction of reaction is not clear at a first glance or strongly depends on the temperature, this can lead to wrong results which are hard to identify.
What is the expected correct behavior?
The backward reaction rate should remain unrestricted just as the forward reaction rate.
Relevant logs and/or images
The left picture shows the results of the test case for the original reversible reaction implementation in OpenFOAM. The reaction does nearly not take place as the backward reaction of the methane steam reaction is restricted. This also leads to nearly no change in the temperature. The right picture shows the results of the same test case for the correct reversible reaction implementation. Both the development of the species and the temperature are in perfect agreement with results generated in Ansys CFX, Ansys Fluent, StarCCM+ and analytical calculation.
 |
 |
Environment information
- OpenFOAM version : v2212, 10 and all previous
- Operating system : Ubuntu 2204 LTS @ WSL2.0 and all other systems
Possible fixes
In $FOAM_SRC/thermophysicalModels/specie/reaction/Reactions/ReversibleReaction/ReversibleReaction.C, the 1e-6 in line 131 should simply be replaced by VSMALL:
Foam::scalar Foam::ReversibleReaction
<
ReactionType,
ReactionThermo,
ReactionRate
>::kr
(
const scalar kfwd,
const scalar p,
const scalar T,
const scalarField& c
) const
{
return kfwd/max(this->Kc(p, T), VSMALL);
}