From 1d572696807bc3b353b6a3ef36ead2d15b2aa6fd Mon Sep 17 00:00:00 2001 From: Henry Weller Date: Sun, 17 Jul 2016 15:13:54 +0100 Subject: [PATCH] TDACChemistryModel: New chemistry model providing Tabulation of Dynamic Adaptive Chemistry MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Provides efficient integration of complex laminar reaction chemistry, combining the advantages of automatic dynamic specie and reaction reduction with ISAT (in situ adaptive tabulation). The advantages grow as the complexity of the chemistry increases. References: Contino, F., Jeanmart, H., Lucchini, T., & D’Errico, G. (2011). Coupling of in situ adaptive tabulation and dynamic adaptive chemistry: An effective method for solving combustion in engine simulations. Proceedings of the Combustion Institute, 33(2), 3057-3064. Contino, F., Lucchini, T., D'Errico, G., Duynslaegher, C., Dias, V., & Jeanmart, H. (2012). Simulations of advanced combustion modes using detailed chemistry combined with tabulation and mechanism reduction techniques. SAE International Journal of Engines, 5(2012-01-0145), 185-196. Contino, F., Foucher, F., Dagaut, P., Lucchini, T., D’Errico, G., & Mounaïm-Rousselle, C. (2013). Experimental and numerical analysis of nitric oxide effect on the ignition of iso-octane in a single cylinder HCCI engine. Combustion and Flame, 160(8), 1476-1483. Contino, F., Masurier, J. B., Foucher, F., Lucchini, T., D’Errico, G., & Dagaut, P. (2014). CFD simulations using the TDAC method to model iso-octane combustion for a large range of ozone seeding and temperature conditions in a single cylinder HCCI engine. Fuel, 137, 179-184. Two tutorial cases are currently provided: + tutorials/combustion/chemFoam/ic8h18_TDAC + tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI_TDAC the first of which clearly demonstrates the advantage of dynamic adaptive chemistry providing ~10x speedup, the second demonstrates ISAT on the modest complex GRI mechanisms for methane combustion, providing a speedup of ~4x. More tutorials demonstrating TDAC on more complex mechanisms and cases will be provided soon in addition to documentation for the operation and settings of TDAC. Also further updates to the TDAC code to improve consistency and integration with the rest of OpenFOAM and further optimize operation can be expected. Original code providing all algorithms for chemistry reduction and tabulation contributed by Francesco Contino, Tommaso Lucchini, Gianluca D’Errico, Hervé Jeanmart, Nicolas Bourgeois and Stéphane Backaert. Implementation updated, optimized and integrated into OpenFOAM-dev by Henry G. Weller, CFD Direct Ltd with the help of Francesco Contino. --- .../solvers/combustion/fireFoam/YEEqn.H | 2 +- .../solvers/combustion/reactingFoam/YEqn.H | 2 +- .../lagrangian/coalChemistryFoam/YEqn.H | 2 +- .../lagrangian/reactingParcelFilmFoam/YEqn.H | 2 +- .../lagrangian/reactingParcelFoam/YEqn.H | 2 +- .../simpleReactingParcelFoam/YEqn.H | 2 +- .../MultiComponentPhaseModel.C | 8 +- bin/tools/CleanFunctions | 1 + .../chemistryModel/Make/files | 3 + .../TDACChemistryModel/TDACChemistryModel.C | 843 ++ .../TDACChemistryModel/TDACChemistryModel.H | 280 + .../TDACChemistryModel/TDACChemistryModelI.H | 152 + .../TDACChemistryModel/reduction/DAC/DAC.C | 746 ++ .../TDACChemistryModel/reduction/DAC/DAC.H | 172 + .../TDACChemistryModel/reduction/DRG/DRG.C | 347 + .../TDACChemistryModel/reduction/DRG/DRG.H | 108 + .../reduction/DRGEP/DRGEP.C | 746 ++ .../reduction/DRGEP/DRGEP.H | 179 + .../reduction/DRGEP/SortableListDRGEP.C | 153 + .../reduction/DRGEP/SortableListDRGEP.H | 142 + .../TDACChemistryModel/reduction/EFA/EFA.C | 568 ++ .../TDACChemistryModel/reduction/EFA/EFA.H | 105 + .../reduction/EFA/SortableListEFA.C | 155 + .../reduction/EFA/SortableListEFA.H | 141 + .../TDACChemistryModel/reduction/PFA/PFA.C | 442 + .../TDACChemistryModel/reduction/PFA/PFA.H | 107 + .../chemistryReductionMethod.C | 58 + .../chemistryReductionMethod.H | 177 + .../chemistryReductionMethodI.H | 84 + .../chemistryReductionMethodNew.C | 99 + .../reduction/makeChemistryReductionMethods.C | 98 + .../reduction/makeChemistryReductionMethods.H | 122 + .../noChemistryReduction.C | 65 + .../noChemistryReduction.H | 101 + .../TDACChemistryModel/tabulation/ISAT/ISAT.C | 622 ++ .../TDACChemistryModel/tabulation/ISAT/ISAT.H | 256 + .../tabulation/ISAT/binaryNode/binaryNode.C | 180 + .../tabulation/ISAT/binaryNode/binaryNode.H | 207 + .../tabulation/ISAT/binaryTree/binaryTree.C | 821 ++ .../tabulation/ISAT/binaryTree/binaryTree.H | 250 + .../ISAT/chemPointISAT/chemPointISAT.C | 985 ++ .../ISAT/chemPointISAT/chemPointISAT.H | 460 + .../chemistryTabulationMethod.C | 55 + .../chemistryTabulationMethod.H | 184 + .../chemistryTabulationMethodNew.C | 99 + .../makeChemistryTabulationMethods.C | 101 + .../makeChemistryTabulationMethods.H | 89 + .../noChemistryTabulation.C | 54 + .../noChemistryTabulation.H | 141 + .../basicChemistryModelTemplates.C | 24 +- .../chemistryModel/chemistryModel.C | 147 +- .../chemistryModel/chemistryModel.H | 25 +- .../psiChemistryModel/psiChemistryModels.C | 76 +- .../rhoChemistryModel/rhoChemistryModels.C | 75 +- .../makeChemistrySolverTypes.H | 23 +- .../basicMultiComponentMixture.C | 86 +- .../basicMultiComponentMixture.H | 23 +- .../basicMultiComponentMixtureI.H | 44 +- .../reactingMixture/reactingMixture.C | 4 + .../reactingMixture/reactingMixture.H | 16 +- .../combustion/chemFoam/ic8h18_TDAC/Allclean | 11 + .../combustion/chemFoam/ic8h18_TDAC/Allrun | 14 + .../chemFoam/ic8h18_TDAC/chemkin/chem.inp | 7902 +++++++++++++++++ .../chemFoam/ic8h18_TDAC/chemkin/senk.inp | 13 + .../chemFoam/ic8h18_TDAC/chemkin/senk.out | 3775 ++++++++ .../chemFoam/ic8h18_TDAC/chemkin/therm.dat | 5959 +++++++++++++ .../ic8h18_TDAC/chemkin/transportProperties | 27 + .../ic8h18_TDAC/constant/chemistryProperties | 75 + .../ic8h18_TDAC/constant/initialConditions | 34 + .../constant/thermophysicalProperties | 33 + .../chemFoam/ic8h18_TDAC/system/controlDict | 55 + .../chemFoam/ic8h18_TDAC/system/fvSchemes | 32 + .../chemFoam/ic8h18_TDAC/system/fvSolution | 30 + .../chemFoam/ic8h18_TDAC/validation/Allrun | 14 + .../ic8h18_TDAC/validation/createGraph | 25 + .../laminar/counterFlowFlame2D/0/T | 6 +- .../laminar/counterFlowFlame2DLTS/0/N2 | 2 +- .../laminar/counterFlowFlame2DLTS/0/T | 6 +- .../laminar/counterFlowFlame2D_GRI/0/CH4 | 48 + .../laminar/counterFlowFlame2D_GRI/0/CO2 | 48 + .../laminar/counterFlowFlame2D_GRI/0/H2O | 48 + .../laminar/counterFlowFlame2D_GRI/0/N2 | 48 + .../laminar/counterFlowFlame2D_GRI/0/O2 | 47 + .../laminar/counterFlowFlame2D_GRI/0/T | 47 + .../laminar/counterFlowFlame2D_GRI/0/U | 46 + .../laminar/counterFlowFlame2D_GRI/0/Ydefault | 48 + .../laminar/counterFlowFlame2D_GRI/0/alphat | 45 + .../laminar/counterFlowFlame2D_GRI/0/p | 44 + .../constant/chemistryProperties | 35 + .../constant/combustionProperties | 27 + .../counterFlowFlame2D_GRI/constant/reactions | 28 + .../constant/reactionsGRI | 5590 ++++++++++++ .../constant/thermo.compressibleGas | 152 + .../constant/thermo.compressibleGasGRI | 1391 +++ .../constant/thermophysicalProperties | 36 + .../constant/turbulenceProperties | 21 + .../system/blockMeshDict | 82 + .../counterFlowFlame2D_GRI/system/controlDict | 53 + .../system/decomposeParDict | 29 + .../counterFlowFlame2D_GRI/system/fvSchemes | 57 + .../counterFlowFlame2D_GRI/system/fvSolution | 69 + .../laminar/counterFlowFlame2D_GRI_TDAC/0/CH4 | 48 + .../laminar/counterFlowFlame2D_GRI_TDAC/0/CO2 | 48 + .../laminar/counterFlowFlame2D_GRI_TDAC/0/H2O | 48 + .../laminar/counterFlowFlame2D_GRI_TDAC/0/N2 | 48 + .../laminar/counterFlowFlame2D_GRI_TDAC/0/O2 | 47 + .../laminar/counterFlowFlame2D_GRI_TDAC/0/T | 47 + .../laminar/counterFlowFlame2D_GRI_TDAC/0/U | 46 + .../counterFlowFlame2D_GRI_TDAC/0/Ydefault | 48 + .../counterFlowFlame2D_GRI_TDAC/0/alphat | 45 + .../laminar/counterFlowFlame2D_GRI_TDAC/0/p | 44 + .../constant/chemistryProperties | 134 + .../constant/chemistryProperties.new | 127 + .../constant/combustionProperties | 27 + .../constant/reactions | 28 + .../constant/reactionsGRI | 5590 ++++++++++++ .../constant/thermo.compressibleGas | 152 + .../constant/thermo.compressibleGasGRI | 1391 +++ .../constant/thermophysicalProperties | 36 + .../constant/turbulenceProperties | 21 + .../system/blockMeshDict | 82 + .../system/controlDict | 53 + .../system/decomposeParDict | 29 + .../system/fvSchemes | 57 + .../system/fvSolution | 69 + 125 files changed, 45393 insertions(+), 185 deletions(-) create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/TDACChemistryModel.C create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/TDACChemistryModel.H create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/TDACChemistryModelI.H create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/reduction/DAC/DAC.C create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/reduction/DAC/DAC.H create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/reduction/DRG/DRG.C create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/reduction/DRG/DRG.H create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/reduction/DRGEP/DRGEP.C create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/reduction/DRGEP/DRGEP.H create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/reduction/DRGEP/SortableListDRGEP.C create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/reduction/DRGEP/SortableListDRGEP.H create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/reduction/EFA/EFA.C create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/reduction/EFA/EFA.H create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/reduction/EFA/SortableListEFA.C create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/reduction/EFA/SortableListEFA.H create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/reduction/PFA/PFA.C create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/reduction/PFA/PFA.H create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/reduction/chemistryReductionMethod/chemistryReductionMethod.C create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/reduction/chemistryReductionMethod/chemistryReductionMethod.H create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/reduction/chemistryReductionMethod/chemistryReductionMethodI.H create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/reduction/chemistryReductionMethod/chemistryReductionMethodNew.C create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/reduction/makeChemistryReductionMethods.C create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/reduction/makeChemistryReductionMethods.H create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/reduction/noChemistryReduction/noChemistryReduction.C create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/reduction/noChemistryReduction/noChemistryReduction.H create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/tabulation/ISAT/ISAT.C create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/tabulation/ISAT/ISAT.H create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/tabulation/ISAT/binaryNode/binaryNode.C create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/tabulation/ISAT/binaryNode/binaryNode.H create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/tabulation/ISAT/binaryTree/binaryTree.C create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/tabulation/ISAT/binaryTree/binaryTree.H create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/tabulation/ISAT/chemPointISAT/chemPointISAT.C create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/tabulation/ISAT/chemPointISAT/chemPointISAT.H create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/tabulation/chemistryTabulationMethod/chemistryTabulationMethod.C create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/tabulation/chemistryTabulationMethod/chemistryTabulationMethod.H create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/tabulation/chemistryTabulationMethod/chemistryTabulationMethodNew.C create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/tabulation/makeChemistryTabulationMethods.C create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/tabulation/makeChemistryTabulationMethods.H create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/tabulation/noChemistryTabulation/noChemistryTabulation.C create mode 100644 src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/tabulation/noChemistryTabulation/noChemistryTabulation.H create mode 100755 tutorials/combustion/chemFoam/ic8h18_TDAC/Allclean create mode 100755 tutorials/combustion/chemFoam/ic8h18_TDAC/Allrun create mode 100644 tutorials/combustion/chemFoam/ic8h18_TDAC/chemkin/chem.inp create mode 100644 tutorials/combustion/chemFoam/ic8h18_TDAC/chemkin/senk.inp create mode 100644 tutorials/combustion/chemFoam/ic8h18_TDAC/chemkin/senk.out create mode 100644 tutorials/combustion/chemFoam/ic8h18_TDAC/chemkin/therm.dat create mode 100644 tutorials/combustion/chemFoam/ic8h18_TDAC/chemkin/transportProperties create mode 100644 tutorials/combustion/chemFoam/ic8h18_TDAC/constant/chemistryProperties create mode 100644 tutorials/combustion/chemFoam/ic8h18_TDAC/constant/initialConditions create mode 100644 tutorials/combustion/chemFoam/ic8h18_TDAC/constant/thermophysicalProperties create mode 100644 tutorials/combustion/chemFoam/ic8h18_TDAC/system/controlDict create mode 100644 tutorials/combustion/chemFoam/ic8h18_TDAC/system/fvSchemes create mode 100644 tutorials/combustion/chemFoam/ic8h18_TDAC/system/fvSolution create mode 100755 tutorials/combustion/chemFoam/ic8h18_TDAC/validation/Allrun create mode 100755 tutorials/combustion/chemFoam/ic8h18_TDAC/validation/createGraph create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI/0/CH4 create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI/0/CO2 create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI/0/H2O create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI/0/N2 create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI/0/O2 create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI/0/T create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI/0/U create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI/0/Ydefault create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI/0/alphat create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI/0/p create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI/constant/chemistryProperties create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI/constant/combustionProperties create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI/constant/reactions create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI/constant/reactionsGRI create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI/constant/thermo.compressibleGas create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI/constant/thermo.compressibleGasGRI create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI/constant/thermophysicalProperties create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI/constant/turbulenceProperties create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI/system/blockMeshDict create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI/system/controlDict create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI/system/decomposeParDict create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI/system/fvSchemes create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI/system/fvSolution create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI_TDAC/0/CH4 create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI_TDAC/0/CO2 create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI_TDAC/0/H2O create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI_TDAC/0/N2 create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI_TDAC/0/O2 create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI_TDAC/0/T create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI_TDAC/0/U create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI_TDAC/0/Ydefault create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI_TDAC/0/alphat create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI_TDAC/0/p create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI_TDAC/constant/chemistryProperties create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI_TDAC/constant/chemistryProperties.new create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI_TDAC/constant/combustionProperties create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI_TDAC/constant/reactions create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI_TDAC/constant/reactionsGRI create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI_TDAC/constant/thermo.compressibleGas create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI_TDAC/constant/thermo.compressibleGasGRI create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI_TDAC/constant/thermophysicalProperties create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI_TDAC/constant/turbulenceProperties create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI_TDAC/system/blockMeshDict create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI_TDAC/system/controlDict create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI_TDAC/system/decomposeParDict create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI_TDAC/system/fvSchemes create mode 100644 tutorials/combustion/reactingFoam/laminar/counterFlowFlame2D_GRI_TDAC/system/fvSolution diff --git a/applications/solvers/combustion/fireFoam/YEEqn.H b/applications/solvers/combustion/fireFoam/YEEqn.H index b84d655a0d..1c9c2176c5 100644 --- a/applications/solvers/combustion/fireFoam/YEEqn.H +++ b/applications/solvers/combustion/fireFoam/YEEqn.H @@ -16,7 +16,7 @@ tmp> mvConvection forAll(Y, i) { - if (i != inertIndex) + if (i != inertIndex && composition.active(i)) { volScalarField& Yi = Y[i]; diff --git a/applications/solvers/combustion/reactingFoam/YEqn.H b/applications/solvers/combustion/reactingFoam/YEqn.H index 246f7119ab..7c3f0f117b 100644 --- a/applications/solvers/combustion/reactingFoam/YEqn.H +++ b/applications/solvers/combustion/reactingFoam/YEqn.H @@ -16,7 +16,7 @@ tmp> mvConvection forAll(Y, i) { - if (i != inertIndex) + if (i != inertIndex && composition.active(i)) { volScalarField& Yi = Y[i]; diff --git a/applications/solvers/lagrangian/coalChemistryFoam/YEqn.H b/applications/solvers/lagrangian/coalChemistryFoam/YEqn.H index 52e44e2c47..4edbce9b42 100644 --- a/applications/solvers/lagrangian/coalChemistryFoam/YEqn.H +++ b/applications/solvers/lagrangian/coalChemistryFoam/YEqn.H @@ -17,7 +17,7 @@ tmp> mvConvection forAll(Y, i) { - if (i != inertIndex) + if (i != inertIndex && composition.active(i)) { volScalarField& Yi = Y[i]; diff --git a/applications/solvers/lagrangian/reactingParcelFilmFoam/YEqn.H b/applications/solvers/lagrangian/reactingParcelFilmFoam/YEqn.H index 1fbb98445b..555b5bc0bc 100644 --- a/applications/solvers/lagrangian/reactingParcelFilmFoam/YEqn.H +++ b/applications/solvers/lagrangian/reactingParcelFilmFoam/YEqn.H @@ -17,7 +17,7 @@ tmp> mvConvection forAll(Y, i) { - if (i != inertIndex) + if (i != inertIndex && composition.active(i)) { volScalarField& Yi = Y[i]; diff --git a/applications/solvers/lagrangian/reactingParcelFoam/YEqn.H b/applications/solvers/lagrangian/reactingParcelFoam/YEqn.H index 7259a1ca2a..db6628936c 100644 --- a/applications/solvers/lagrangian/reactingParcelFoam/YEqn.H +++ b/applications/solvers/lagrangian/reactingParcelFoam/YEqn.H @@ -16,7 +16,7 @@ tmp> mvConvection forAll(Y, i) { - if (i != inertIndex) + if (i != inertIndex && composition.active(i)) { volScalarField& Yi = Y[i]; diff --git a/applications/solvers/lagrangian/reactingParcelFoam/simpleReactingParcelFoam/YEqn.H b/applications/solvers/lagrangian/reactingParcelFoam/simpleReactingParcelFoam/YEqn.H index 6498fd73f2..7c7cf9a4c0 100644 --- a/applications/solvers/lagrangian/reactingParcelFoam/simpleReactingParcelFoam/YEqn.H +++ b/applications/solvers/lagrangian/reactingParcelFoam/simpleReactingParcelFoam/YEqn.H @@ -16,7 +16,7 @@ tmp> mvConvection forAll(Y, i) { - if (i != inertIndex) + if (i != inertIndex && composition.active(i)) { volScalarField& Yi = Y[i]; diff --git a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseModel/MultiComponentPhaseModel/MultiComponentPhaseModel.C b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseModel/MultiComponentPhaseModel/MultiComponentPhaseModel.C index efa486a0e4..081e5c893a 100644 --- a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseModel/MultiComponentPhaseModel/MultiComponentPhaseModel.C +++ b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/phaseModel/MultiComponentPhaseModel/MultiComponentPhaseModel.C @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 2015 OpenFOAM Foundation + \\ / A nd | Copyright (C) 2015-2016 OpenFOAM Foundation \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -142,6 +142,12 @@ Foam::MultiComponentPhaseModel::YiEqn this->name() ) ) + || ( + !this->thermo_->composition().active + ( + this->thermo_->composition().species()[Yi.member()] + ) + ) ) { return tmp(); diff --git a/bin/tools/CleanFunctions b/bin/tools/CleanFunctions index dde53bafc3..3df910567d 100644 --- a/bin/tools/CleanFunctions +++ b/bin/tools/CleanFunctions @@ -57,6 +57,7 @@ cleanCase() rm -rf processor* > /dev/null 2>&1 rm -rf postProcessing > /dev/null 2>&1 + rm -rf TDAC > /dev/null 2>&1 rm -rf probes* > /dev/null 2>&1 rm -rf forces* > /dev/null 2>&1 rm -rf graphs* > /dev/null 2>&1 diff --git a/src/thermophysicalModels/chemistryModel/Make/files b/src/thermophysicalModels/chemistryModel/Make/files index ef8c6d84c8..d0782c33b9 100644 --- a/src/thermophysicalModels/chemistryModel/Make/files +++ b/src/thermophysicalModels/chemistryModel/Make/files @@ -6,6 +6,9 @@ chemistryModel/psiChemistryModel/psiChemistryModels.C chemistryModel/rhoChemistryModel/rhoChemistryModel.C chemistryModel/rhoChemistryModel/rhoChemistryModels.C +chemistryModel/TDACChemistryModel/reduction/makeChemistryReductionMethods.C +chemistryModel/TDACChemistryModel/tabulation/makeChemistryTabulationMethods.C + chemistrySolver/chemistrySolver/makeChemistrySolvers.C LIB = $(FOAM_LIBBIN)/libchemistryModel diff --git a/src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/TDACChemistryModel.C b/src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/TDACChemistryModel.C new file mode 100644 index 0000000000..de13a49d01 --- /dev/null +++ b/src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/TDACChemistryModel.C @@ -0,0 +1,843 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2016 OpenFOAM Foundation + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software: you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM. If not, see . + +\*---------------------------------------------------------------------------*/ + +#include "TDACChemistryModel.H" +#include "UniformField.H" +#include "clockTime.H" + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +template +Foam::TDACChemistryModel::TDACChemistryModel +( + const fvMesh& mesh, + const word& phaseName +) +: + chemistryModel(mesh, phaseName), + NsDAC_(this->nSpecie_), + completeC_(this->nSpecie_,0.0), + reactionsDisabled_(this->reactions_.size(), false), + completeToSimplifiedIndex_(this->nSpecie_,-1), + simplifiedToCompleteIndex_(this->nSpecie_), + specieComp_(this->nSpecie_) +{ + basicMultiComponentMixture& composition = this->thermo().composition(); + + // Store the species composition according to the species index + speciesTable speciesTab = composition.species(); + + const HashTable>& specComp = + dynamicCast&>(this->thermo()) + .specieComposition(); + + forAll(specieComp_,i) + { + specieComp_[i] = specComp[this->Y()[i].name()]; + } + + mechRed_ = chemistryReductionMethod::New + ( + *this, + *this + ); + + // When the mechanism reduction method is used, the 'active' flag for every + // species should be initialized (by default 'active' is true) + if (mechRed_->active()) + { + forAll(this->Y(), i) + { + IOobject header + ( + this->Y()[i].name(), + mesh.time().timeName(), + mesh, + IOobject::NO_READ + ); + + // Check if the species file is provided, if not set inactive + // and NO_WRITE + if (!header.headerOk()) + { + composition.setInactive(i); + this->Y()[i].writeOpt() = IOobject::NO_WRITE; + } + } + } + + tabulation_ = chemistryTabulationMethod::New + ( + *this, + *this + ); + + if (mechRed_->log()) + { + cpuReduceFile_ = logFile("cpu_reduce.out"); + nActiveSpeciesFile_ = logFile("nActiveSpecies.out"); + } + + if (tabulation_->log()) + { + cpuAddFile_ = logFile("cpu_add.out"); + cpuRetrieveFile_ = logFile("cpu_retrieve.out"); + } + + if (mechRed_->log() || tabulation_->log()) + { + cpuSolveFile_ = logFile("cpu_solve.out"); + } +} + + +// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // + +template +Foam::TDACChemistryModel::~TDACChemistryModel() +{} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +template +void Foam::TDACChemistryModel::omega +( + const scalarField& c, // Contains all species even when mechRed is active + const scalar T, + const scalar p, + scalarField& dcdt +) const +{ + const bool reduced = mechRed_->active(); + + scalar pf, cf, pr, cr; + label lRef, rRef; + + dcdt = Zero; + + forAll(this->reactions_, i) + { + if (!reactionsDisabled_[i]) + { + const Reaction& R = this->reactions_[i]; + + scalar omegai = omega + ( + R, c, T, p, pf, cf, lRef, pr, cr, rRef + ); + + forAll(R.lhs(), s) + { + label si = R.lhs()[s].index; + if (reduced) + { + si = completeToSimplifiedIndex_[si]; + } + + const scalar sl = R.lhs()[s].stoichCoeff; + dcdt[si] -= sl*omegai; + } + forAll(R.rhs(), s) + { + label si = R.rhs()[s].index; + if (reduced) + { + si = completeToSimplifiedIndex_[si]; + } + + const scalar sr = R.rhs()[s].stoichCoeff; + dcdt[si] += sr*omegai; + } + } + } +} + + +template +Foam::scalar Foam::TDACChemistryModel::omega +( + const Reaction& R, + const scalarField& c, // Contains all species even when mechRed is active + const scalar T, + const scalar p, + scalar& pf, + scalar& cf, + label& lRef, + scalar& pr, + scalar& cr, + label& rRef +) const +{ + const scalar kf = R.kf(p, T, c); + const scalar kr = R.kr(kf, p, T, c); + + const label Nl = R.lhs().size(); + const label Nr = R.rhs().size(); + + label slRef = 0; + lRef = R.lhs()[slRef].index; + + pf = kf; + for (label s=1; s SMALL) + { + pf *= pow(cf, exp - 1); + } + else + { + pf = 0; + } + } + else + { + pf *= pow(cf, exp - 1); + } + } + + label srRef = 0; + rRef = R.rhs()[srRef].index; + + // Find the matrix element and element position for the rhs + pr = kr; + for (label s=1; sSMALL) + { + pr *= pow(cr, exp - 1); + } + else + { + pr = 0; + } + } + else + { + pr *= pow(cr, exp - 1); + } + } + + return pf*cf - pr*cr; +} + + +template +void Foam::TDACChemistryModel::derivatives +( + const scalar time, + const scalarField& c, + scalarField& dcdt +) const +{ + const bool reduced = mechRed_->active(); + + const scalar T = c[this->nSpecie_]; + const scalar p = c[this->nSpecie_ + 1]; + + if (reduced) + { + // When using DAC, the ODE solver submit a reduced set of species the + // complete set is used and only the species in the simplified mechanism + // are updated + this->c_ = completeC_; + + // Update the concentration of the species in the simplified mechanism + // the other species remain the same and are used only for third-body + // efficiencies + for (label i=0; ic_[simplifiedToCompleteIndex_[i]] = max(0.0, c[i]); + } + } + else + { + for (label i=0; inSpecie(); i++) + { + this->c_[i] = max(0.0, c[i]); + } + } + + omega(this->c_, T, p, dcdt); + + // Constant pressure + // dT/dt = ... + scalar rho = 0; + for (label i=0; ic_.size(); i++) + { + const scalar W = this->specieThermo_[i].W(); + rho += W*this->c_[i]; + } + + scalar cp = 0; + for (label i=0; ic_.size(); i++) + { + // cp function returns [J/(kmol K)] + cp += this->c_[i]*this->specieThermo_[i].cp(p, T); + } + cp /= rho; + + // When mechanism reduction is active + // dT is computed on the reduced set since dcdt is null + // for species not involved in the simplified mechanism + scalar dT = 0; + for (label i=0; inSpecie_; i++) + { + label si; + if (reduced) + { + si = simplifiedToCompleteIndex_[i]; + } + else + { + si = i; + } + + // ha function returns [J/kmol] + const scalar hi = this->specieThermo_[si].ha(p, T); + dT += hi*dcdt[i]; + } + dT /= rho*cp; + + dcdt[this->nSpecie_] = -dT; + + // dp/dt = ... + dcdt[this->nSpecie_ + 1] = 0; +} + + +template +void Foam::TDACChemistryModel::jacobian +( + const scalar t, + const scalarField& c, + scalarSquareMatrix& dfdc +) const +{ + const bool reduced = mechRed_->active(); + + // If the mechanism reduction is active, the computed Jacobian + // is compact (size of the reduced set of species) + // but according to the informations of the complete set + // (i.e. for the third-body efficiencies) + + const scalar T = c[this->nSpecie_]; + const scalar p = c[this->nSpecie_ + 1]; + + if (reduced) + { + this->c_ = completeC_; + for (label i=0; ic_[simplifiedToCompleteIndex_[i]] = max(0.0, c[i]); + } + } + else + { + forAll(this->c_, i) + { + this->c_[i] = max(c[i], 0.0); + } + } + + dfdc = Zero; + + forAll(this->reactions_, ri) + { + if (!reactionsDisabled_[ri]) + { + const Reaction& R = this->reactions_[ri]; + + const scalar kf0 = R.kf(p, T, this->c_); + const scalar kr0 = R.kr(kf0, p, T, this->c_); + + forAll(R.lhs(), j) + { + label sj = R.lhs()[j].index; + if (reduced) + { + sj = completeToSimplifiedIndex_[sj]; + } + scalar kf = kf0; + forAll(R.lhs(), i) + { + const label si = R.lhs()[i].index; + const scalar el = R.lhs()[i].exponent; + if (i == j) + { + if (el < 1) + { + if (this->c_[si] > SMALL) + { + kf *= el*pow(this->c_[si] + VSMALL, el - 1); + } + else + { + kf = 0; + } + } + else + { + kf *= el*pow(this->c_[si], el - 1); + } + } + else + { + kf *= pow(this->c_[si], el); + } + } + + forAll(R.lhs(), i) + { + label si = R.lhs()[i].index; + if (reduced) + { + si = completeToSimplifiedIndex_[si]; + } + const scalar sl = R.lhs()[i].stoichCoeff; + dfdc(si, sj) -= sl*kf; + } + forAll(R.rhs(), i) + { + label si = R.rhs()[i].index; + if (reduced) + { + si = completeToSimplifiedIndex_[si]; + } + const scalar sr = R.rhs()[i].stoichCoeff; + dfdc(si, sj) += sr*kf; + } + } + + forAll(R.rhs(), j) + { + label sj = R.rhs()[j].index; + if (reduced) + { + sj = completeToSimplifiedIndex_[sj]; + } + scalar kr = kr0; + forAll(R.rhs(), i) + { + const label si = R.rhs()[i].index; + const scalar er = R.rhs()[i].exponent; + if (i == j) + { + if (er < 1) + { + if (this->c_[si] > SMALL) + { + kr *= er*pow(this->c_[si] + VSMALL, er - 1); + } + else + { + kr = 0; + } + } + else + { + kr *= er*pow(this->c_[si], er - 1); + } + } + else + { + kr *= pow(this->c_[si], er); + } + } + + forAll(R.lhs(), i) + { + label si = R.lhs()[i].index; + if (reduced) + { + si = completeToSimplifiedIndex_[si]; + } + const scalar sl = R.lhs()[i].stoichCoeff; + dfdc(si, sj) += sl*kr; + } + forAll(R.rhs(), i) + { + label si = R.rhs()[i].index; + if (reduced) + { + si = completeToSimplifiedIndex_[si]; + } + const scalar sr = R.rhs()[i].stoichCoeff; + dfdc(si, sj) -= sr*kr; + } + } + } + } + + // Calculate the dcdT elements numerically + const scalar delta = 1e-3; + + omega(this->c_, T + delta, p, this->dcdt_); + for (label i=0; inSpecie_; i++) + { + dfdc(i, this->nSpecie_) = this->dcdt_[i]; + } + + omega(this->c_, T - delta, p, this->dcdt_); + for (label i=0; inSpecie_; i++) + { + dfdc(i, this->nSpecie_) = + 0.5*(dfdc(i, this->nSpecie_) - this->dcdt_[i])/delta; + } + + dfdc(this->nSpecie_, this->nSpecie_) = 0; + dfdc(this->nSpecie_ + 1, this->nSpecie_) = 0; +} + + +template +void Foam::TDACChemistryModel::jacobian +( + const scalar t, + const scalarField& c, + scalarField& dcdt, + scalarSquareMatrix& dfdc +) const +{ + jacobian(t, c, dfdc); + + const scalar T = c[this->nSpecie_]; + const scalar p = c[this->nSpecie_ + 1]; + + // Note: Uses the c_ field initialized by the call to jacobian above + omega(this->c_, T, p, dcdt); +} + + +template +template +Foam::scalar Foam::TDACChemistryModel::solve +( + const DeltaTType& deltaT +) +{ + const bool reduced = mechRed_->active(); + + basicMultiComponentMixture& composition = this->thermo().composition(); + + // CPU time analysis + const clockTime clockTime_= clockTime(); + clockTime_.timeIncrement(); + scalar reduceMechCpuTime_ = 0; + scalar addNewLeafCpuTime_ = 0; + scalar solveChemistryCpuTime_ = 0; + scalar searchISATCpuTime_ = 0; + + // Average number of active species + scalar nActiveSpecies = 0; + scalar nAvg = 0; + + CompType::correct(); + + scalar deltaTMin = GREAT; + + if (!this->chemistry_) + { + return deltaTMin; + } + + const volScalarField rho + ( + IOobject + ( + "rho", + this->time().timeName(), + this->mesh(), + IOobject::NO_READ, + IOobject::NO_WRITE, + false + ), + this->thermo().rho() + ); + + const scalarField& T = this->thermo().T(); + const scalarField& p = this->thermo().p(); + + scalarField c(this->nSpecie_); + scalarField c0(this->nSpecie_); + + // Composition vector (Yi, T, p) + scalarField phiq(this->nEqns()); + + scalarField Rphiq(this->nEqns()); + + forAll(rho, celli) + { + const scalar rhoi = rho[celli]; + scalar pi = p[celli]; + scalar Ti = T[celli]; + + for (label i=0; inSpecie_; i++) + { + c[i] = rhoi*this->Y_[i][celli]/this->specieThermo_[i].W(); + c0[i] = c[i]; + phiq[i] = this->Y()[i][celli]; + } + phiq[this->nSpecie()]=Ti; + phiq[this->nSpecie()+1]=pi; + + // Initialise time progress + scalar timeLeft = deltaT[celli]; + + // Not sure if this is necessary + Rphiq = Zero; + + clockTime_.timeIncrement(); + + // When tabulation is active (short-circuit evaluation for retrieve) + // It first tries to retrieve the solution of the system with the + // information stored through the tabulation method + if (tabulation_->active() && tabulation_->retrieve(phiq, Rphiq)) + { + // Retrieved solution stored in Rphiq + for (label i=0; inSpecie(); i++) + { + c[i] = rhoi*Rphiq[i]/this->specieThermo_[i].W(); + } + + searchISATCpuTime_ += clockTime_.timeIncrement(); + } + // This position is reached when tabulation is not used OR + // if the solution is not retrieved. + // In the latter case, it adds the information to the tabulation + // (it will either expand the current data or add a new stored poin). + else + { + clockTime_.timeIncrement(); + if (reduced) + { + // Reduce mechanism change the number of species (only active) + mechRed_->reduceMechanism(c, Ti, pi); + nActiveSpecies += mechRed_->NsSimp(); + nAvg++; + } + + reduceMechCpuTime_ += clockTime_.timeIncrement(); + + // Calculate the chemical source terms + while (timeLeft > SMALL) + { + scalar dt = timeLeft; + if (reduced) + { + // completeC_ used in the overridden ODE methods + // to update only the active species + completeC_ = c; + + // Solve the reduced set of ODE + this->solve + ( + simplifiedC_, Ti, pi, dt, this->deltaTChem_[celli] + ); + + for (label i=0; isolve(c, Ti, pi, dt, this->deltaTChem_[celli]); + } + timeLeft -= dt; + } + + solveChemistryCpuTime_ += clockTime_.timeIncrement(); + + // If tabulation is used, we add the information computed here to + // the stored points (either expand or add) + if (tabulation_->active()) + { + forAll(c, i) + { + Rphiq[i] = c[i]/rhoi*this->specieThermo_[i].W(); + } + + Rphiq[Rphiq.size()-2] = Ti; + Rphiq[Rphiq.size()-1] = pi; + tabulation_->add(phiq, Rphiq, rhoi); + } + + addNewLeafCpuTime_ += clockTime_.timeIncrement(); + + // When operations are done and if mechanism reduction is active, + // the number of species (which also affects nEqns) is set back + // to the total number of species (stored in the mechRed object) + if (reduced) + { + this->nSpecie_ = mechRed_->nSpecie(); + } + deltaTMin = min(this->deltaTChem_[celli], deltaTMin); + } + + // Set the RR vector (used in the solver) + for (label i=0; inSpecie_; i++) + { + this->RR_[i][celli] = + (c[i] - c0[i])*this->specieThermo_[i].W()/deltaT[celli]; + } + } + + if (mechRed_->log() || tabulation_->log()) + { + cpuSolveFile_() + << this->time().timeOutputValue() + << " " << solveChemistryCpuTime_ << endl; + } + + if (mechRed_->log()) + { + cpuReduceFile_() + << this->time().timeOutputValue() + << " " << reduceMechCpuTime_ << endl; + } + + if (tabulation_->active()) + { + // Every time-step, look if the tabulation should be updated + tabulation_->update(); + + // Write the performance of the tabulation + tabulation_->writePerformance(); + + if (tabulation_->log()) + { + cpuRetrieveFile_() + << this->time().timeOutputValue() + << " " << searchISATCpuTime_ << endl; + + cpuAddFile_() + << this->time().timeOutputValue() + << " " << addNewLeafCpuTime_ << endl; + } + } + + if (reduced && nAvg && mechRed_->log()) + { + // Write average number of species + nActiveSpeciesFile_() + << this->time().timeOutputValue() + << " " << nActiveSpecies/nAvg << endl; + } + + if (Pstream::parRun()) + { + List active(composition.active()); + Pstream::listCombineGather(active, orEqOp()); + Pstream::listCombineScatter(active); + + forAll(active, i) + { + if (active[i]) + { + composition.setActive(i); + } + } + } + + forAll(this->Y(), i) + { + if (composition.active(i)) + { + this->Y()[i].writeOpt() = IOobject::AUTO_WRITE; + } + } + + return deltaTMin; +} + + +template +Foam::scalar Foam::TDACChemistryModel::solve +( + const scalar deltaT +) +{ + // Don't allow the time-step to change more than a factor of 2 + return min + ( + this->solve>(UniformField(deltaT)), + 2*deltaT + ); +} + + +template +Foam::scalar Foam::TDACChemistryModel::solve +( + const scalarField& deltaT +) +{ + return this->solve(deltaT); +} + + +// ************************************************************************* // diff --git a/src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/TDACChemistryModel.H b/src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/TDACChemistryModel.H new file mode 100644 index 0000000000..52a8f0860f --- /dev/null +++ b/src/thermophysicalModels/chemistryModel/chemistryModel/TDACChemistryModel/TDACChemistryModel.H @@ -0,0 +1,280 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2016 OpenFOAM Foundation + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software: you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM. If not, see . + +Class + Foam::TDACChemistryModel + +Description + Extends chemistryModel by adding the TDAC method. + + References: + \verbatim + Contino, F., Jeanmart, H., Lucchini, T., & D’Errico, G. (2011). + Coupling of in situ adaptive tabulation and dynamic adaptive chemistry: + An effective method for solving combustion in engine simulations. + Proceedings of the Combustion Institute, 33(2), 3057-3064. + + Contino, F., Lucchini, T., D'Errico, G., Duynslaegher, C., + Dias, V., & Jeanmart, H. (2012). + Simulations of advanced combustion modes using detailed chemistry + combined with tabulation and mechanism reduction techniques. + SAE International Journal of Engines, + 5(2012-01-0145), 185-196. + + Contino, F., Foucher, F., Dagaut, P., Lucchini, T., D’Errico, G., & + Mounaïm-Rousselle, C. (2013). + Experimental and numerical analysis of nitric oxide effect on the + ignition of iso-octane in a single cylinder HCCI engine. + Combustion and Flame, 160(8), 1476-1483. + + Contino, F., Masurier, J. B., Foucher, F., Lucchini, T., D’Errico, G., & + Dagaut, P. (2014). + CFD simulations using the TDAC method to model iso-octane combustion + for a large range of ozone seeding and temperature conditions + in a single cylinder HCCI engine. + Fuel, 137, 179-184. + \endverbatim + +SourceFiles + TDACChemistryModelI.H + TDACChemistryModel.C + +\*---------------------------------------------------------------------------*/ + +#ifndef TDACChemistryModel_H +#define TDACChemistryModel_H + +#include "chemistryModel.H" +#include "chemistryReductionMethod.H" +#include "chemistryTabulationMethod.H" +#include "OFstream.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +/*---------------------------------------------------------------------------*\ + Class TDACChemistryModel Declaration +\*---------------------------------------------------------------------------*/ + +template +class TDACChemistryModel +: + public chemistryModel +{ + // Private member data + + // Mechanism reduction + label NsDAC_; + scalarField completeC_; + scalarField simplifiedC_; + Field reactionsDisabled_; + Field