View Issue Details
| ID | Project | Category | View Status | Date Submitted | Last Update |
|---|---|---|---|---|---|
| 0000327 | OpenFOAM | Bug | public | 2011-10-27 02:04 | 2015-01-29 09:51 |
| Reporter | russelldub | Assigned To | henry | ||
| Priority | normal | Severity | major | Reproducibility | always |
| Status | resolved | Resolution | fixed | ||
| Platform | ALL | ||||
| Summary | 0000327: Mixtures and Reactions using janafThermo sometimes produce incorrect thermodynamic properties | ||||
| Description | When the janafThermo class is used to construct mixtures and Reactions the JANAF coefficients of individual species are combined algebraically (cf. janafThermoI.H, multiComponentMixture.C, Reaction.C). If the combination has the same Tcommon_ as all constituent species then this will produce the same results for thermodynamic properties as calculating each individually and then combining. This does not happen in practice which means JANAF polynomial fits are used outside of their specified ranges, leading to errors in calculated properties. For mixtures: Tcommon_ is taken as a weighted average of that for the components. While this is likely the best possible way to enable the algebraic combination of janaf coefficients, it still can lead to unacceptable errors (details below). For reactions: Tcommon_ is calculated as a weighted difference of reactants and products. Depending on the order of the reaction this can lead to Tcommon_ being zero, or some very large number, meaning that the high temperature coefficients will always be used or vice versa. This affects reactions which have their reverse rates calculated by equilibrium (i.e. "reversibleReaction"s). | ||||
| Steps To Reproduce | Download accompanying gzipped tarball. Compile example programs chemFoamEx2 and chemFoamEx3. Run chemFoamEx2 on case nc7h16_chemFoamEx2.* Note that errors in temperature of greater than 1 K are possible. Run chemFoamEx3 on case h2_chemFoamEx3.** Note that Tcommon_ for most reactions is 0. (Third number on line after those starting with "H2") This means that high temperature coefficients will be used to calculate the reverse reaction rate for all these reactions regardless of temperature. * The only material difference with tutorial nc7h16 was an increase in the concentration of nc7h16. ** This case is effectively identical to tutorial h2. | ||||
| Additional Information | Proposed fix: Modification of janafThermo (or creation of new thermo (e.g. janafMixThermo)) that saves relative quantities of all constituents in addition to their JANAF coefficients and temperature limits. | ||||
| Tags | No tags attached. | ||||
|
|
|
|
|
What range of Tcommon do you have for your set of species? Do you have or can you generate JANAF coefficients for all you species which use the same Tcommon? Is the averaging of Tcommon acceptable for the mixtures? The main issue appears to be with the way in which Tcommon is calculated for a reaction, this should be easy to fix and we will look into it. JANAF thermodynamics is already expensive, how much overhead would your proposal add for a system containing a large number of species? |
|
|
For reactions with an equal number of molecules on each side the average JANAF coeffcients and Tcommon cannot be calculated: H2 1e-15 0.5 200 5000 0 -7.02007e+14 1.34225e+12 -5.89e+08 107584 -7.46174 -7.47943e+18 2.1621e+15 -1.04855e+15 2.46495e+12 -3.86423e+09 4.67613e+06 -2076.93 -7.33074e+18 4.06535e+15 Note the 1e-15 in the number ef moles! However these coefficients are not used, only the information concerning the rates of reaction and equilibrium coefficients are used. The temperature difference obtained using JANAF functions with different Tcommon can be resolved by refitting the functions so that all species use the same Tcommon, thus removing the ambigutiy and preserving the relative efficiency of the current method. |
|
|
henry, your analysis of the problem concerning reactions is incorrect, and I'm unsatisfied with your solution concerning mixtures. Reactions: The coefficients are most certainly used. They are used in calculating the equilibrium constants (cf. specieThermoI.H and janafThermoI.H). I've attached a new example application and case which help to show this more explicitly. Ten (10) % error in an equilibrium constant (and hence a rate coefficient) is not in my opinion an acceptable error. Note that the example checks only one reaction of one chemical mechanism. The discrepancies could be much larger considering the fits are polynomials and are used way outside of the region to which they were fitted. Mixtures: Your resolution of this issue without changes leaves the OpenFOAM user three options: 1. Accept (or more likely be blissfully ignorant of) the innaccuracy caused by the frequent use of unequal Tcommon in thermodynamic databases used in the modeling community (including those distributed in OpenFOAM tutorials). An error one cannot fully evaluate without significant effort for each situation one may be trying to model. 2. Perform the onerous task of refitting each database one might use to have Tcommon be equal for all species that could be combined in a mixture in a specific thermodynamic database. Keeping in mind that if Tcommon was chosen different from standard the likely reason was to maintain acceptable error between the fitted polynomials and the measured properties. 3. Modify, test, and maintain, in parallel to OpenFOAM development, one's own thermodynamic class, which corrects the issue and doesn't risk significant errors in thermodynamic properties of mixtures. At the very least the code should warn when a thermodynamic database is used that is at risk of this problem. While calculating the thermodynamics in the correct way may be significantly more expensive than the current implementation (at worst nSpecies times), it should still not be a significant contribution to the cost of the overall simulation considering that the integration of chemical reactors and the solution of the mass, momentum, species, and energy conservation equations are much more numerically intensive. Thank you for your consideration of these arguments. |
|
|
|
|
|
Thank your for the further clarification of the problem. Adding a warning concerning the averaging of non-equal Tcommon is trivial and I am happy to do that. Changing the implementation so that the thermodynamics is nSpecies more expensive is actually a big change as for many cases we run the cost of the thermodynamics is already significant and we are looking at ways to reduce this cost. I will also reconsider the method of averaging Tcommon and report back. |
|
|
I have corrected the averaging of Tcommon for reactions and put in a warning for when Tcommon requires averaging. See commit 3354f73de21c75339ed7bf489cbd170956f3dcc4 in OpenFOAM-2.0.x Your test cases now appear to run correctly. Could you please try these changes yourself and let us know. |
|
|
Averaging Tcommon is not quite right with this fix: When using the -= operator in Reaction.C for the combination of thermodynamic data the number of moles of the reaction (which in the end is the change in moles for the reaction) is not the correct quantity to weight by for averaging. As an example, see the Tcommon for reaction 15 (CH3 + OH = CH3OH) of the nc7h16 tutorial which should have averaged Tcommon = (2*1000 + 1357)/3 = 1119 but instead has Tcommon = 1357 (plus digits of order 1e-15). It's necessary to somehow keep track of the total number of moles involved in a reaction while applying the correct algebra to the coefficients. The Tcommon of janafThermo mixtures created by use of subtraction operators are not going to behave in a consistent manner either. Consider adding 1 mole of NC7H16 to 1 mole of O2 now with Tcommon = (1391 + 1000)/2 = 1195.5. Now subtract 1 mole of O2 so that Tcommon = (2*1195.5 + 1000)/3 = 1130.333. If this is repeated, the problem only gets worse. The warning is not exactly what I was imagining and seems a bit overbearing now. I was actually envisioning something in the chemistry reader that would give a one time warning. As it is, running chemFoam on the nc7h16 tutorial and directing output to a file results in a 54 MB file which is 99.9 % warnings. How difficult would it be to keep the current janafThermo class and implement a new one that does the safer/more expensive/more accurate calculation, which could then be used at the modeler's discretion. |
|
|
The more I think about this, the more I believe that unequal Tcommon should be a FatalError with the current implementation of janafThermo. There's just no way to know (without testing across all possible mixtures and temperatures) how bad an answer you're getting. As mentioned before, you're using a high order polynomial fit outside of where it was specified. Depending on the particular combination of polynomials and how wildly they behave outside of the fitted range, this could lead to a very large deviation. What types of simulations are you running for which janafThermo is a significant part of the cost? |
|
|
> The more I think about this, the more I believe that unequal Tcommon > should be a FatalError with the current implementation of janafThermo. I have also been rethinking the averaging strategy and I agree with you. If averaging is appropriate it should be declared to be so but the user by providing the input thermo data with pre-averaged and unified Tcommon. I will make this change. > What types of simulations are you running for which janafThermo is a > significant part of the cost? Any of the combustion calculations we do with JANAF thermodynamics using dieselEngineFoam, dieselFoam, engineFoam, fireFoam, PDRFoam, reactingFoam, rhoReactingFoam or XiFoam. |
|
|
With the change (commit 329e824b) I've confirmed that the following tutorials fail with FatalErrorIn janafThermo: chemFoam/gri chemFoam/nc7h16 chemFoam/ic8h18 dieselFoam/aachenBomb Should I file that as a separate bug? Is an enhancement request for a utility to automatically refit a thermo database appropriate (which could be used to fix the above tutorials)? How about one requesting a thermo class that doesn't require averaging? Should this discussion be held in some other forum? |
|
|
Yes indeed, I have not yet updated the tutorials and will do so and once done I will close this bug report. Your other requests are for additional features rather than bug fixes. If you would like to develop these enhancements yourself and provide them we can consider including and maintaining them in future releases of OpenFOAM. Alternatively we undertake contract developments and could take on this work for you. |
|
|
I've done some more investigation on this topic and gotten some interesting timing results. Please find attached a new example solver which compares timing for calculation of thermodynamic properties (cp, h, and s) using the mixture paradigm versus calculated for each species and combined afterwards. For janafThermo on my test machine the mixture calculation is 1.7-1.9 times SLOWER than the individual calculation for the chemFoam tutorial cases (h2, gri, nc7h16, ic8h18). Note that I arbitrarily set Tcommon = 1000 for the thermo databases in these tutorials to get all of them to run in the latest commits on github. If the thermo classes were re-worked to use combinations of individual species properties instead of the mixture idea then the whole calculation would be sped up and there would be no need to refit thermo databases to have equal Tcommon. Win and win. From looking at it, I'm not sure how easy it is to implement (it seems it may need more than a new species thermo class but some changes in multiComponentMixture and perhaps higher), but I believe it to be a worthwhile effort. Thanks again for your consideration of this information. (Please excuse the inclusion of the Debug object files in the example solver. Timings were conducted with the optimized binary.) |
|
|
|
|
|
Already resolved |
| Date Modified | Username | Field | Change |
|---|---|---|---|
| 2011-10-27 02:04 | russelldub | New Issue | |
| 2011-10-27 02:04 | russelldub | File Added: janafThermoBug.tgz | |
| 2011-10-27 09:41 | henry | Note Added: 0000757 | |
| 2011-10-27 17:03 | henry | Note Added: 0000764 | |
| 2011-10-27 17:03 | henry | Status | new => resolved |
| 2011-10-27 17:03 | henry | Resolution | open => fixed |
| 2011-10-27 17:03 | henry | Assigned To | => henry |
| 2011-10-28 03:11 | russelldub | Note Added: 0000765 | |
| 2011-10-28 03:11 | russelldub | Status | resolved => feedback |
| 2011-10-28 03:11 | russelldub | Resolution | fixed => reopened |
| 2011-10-28 03:12 | russelldub | File Added: janafThermoBug2.tgz | |
| 2011-10-28 03:13 | russelldub | Status | feedback => assigned |
| 2011-10-28 03:14 | russelldub | Note Edited: 0000765 | |
| 2011-10-28 11:45 | henry | Note Added: 0000767 | |
| 2011-10-28 15:14 | henry | Note Added: 0000771 | |
| 2011-10-29 00:40 | russelldub | Note Added: 0000773 | |
| 2011-10-29 02:16 | russelldub | Note Added: 0000774 | |
| 2011-10-29 10:12 | henry | Note Added: 0000775 | |
| 2011-10-30 20:11 | russelldub | Note Added: 0000776 | |
| 2011-10-30 20:55 | henry | Note Added: 0000777 | |
| 2011-11-05 01:23 | russelldub | Note Added: 0000788 | |
| 2011-11-05 01:24 | russelldub | File Added: janafTiming.tgz | |
| 2011-11-05 01:26 | russelldub | Note Edited: 0000788 | |
| 2015-01-28 19:43 | henry | Note Added: 0003610 | |
| 2015-01-28 19:43 | henry | Status | assigned => resolved |
| 2015-01-28 19:43 | henry | Resolution | reopened => fixed |
