.gitignore

@@ -82,20 +82,3 @@ doc/Doxygen/DTAGS
 
 # Ignore the test directory
 /tutorialsTest
-
-.*.swp
-*.gz
-log.*
-postProcessing
-*~
-processor*
-[0-9]*
-[0-9]*.[0-9]*
-!0.org
-constant/extendedFeatureEdgeMesh
-*.eMesh
-constant/polyMesh/*
-!constant/polyMesh/blockMeshDict
-constant/cellToRegion
-*.OpenFOAM
-*.foam

examples/FlameD_Equilibrium/constant/polyMesh/blockMeshDict

@@ -19,28 +19,30 @@ convertToMeters 0.001;
 vertices
 (
 
-(0 0 0)
-(0 576 0)
-(-0.15702979454  0 3.5965736)
-(-0.15702979454  576 3.5965736)
-(-0.167934643323 0 3.8463357)
-(-0.167934643323 576 3.8463357)
-(-0.396936424116 0 9.0913388)
-(-0.396936424116 576 9.0913388)
-(-0.412203210667 0 9.4410057)
-(-0.412203210667 576 9.4410057)
-(-6.5429081      0 149.8572332)
-(-6.5429081      576 149.8572332)
-(0.15702979454  0 3.5965736)
-(0.15702979454  576 3.5965736)
-(0.167934643323 0 3.8463357)
-(0.167934643323 576 3.8463357)
-(0.396936424116 0 9.0913388)
-(0.396936424116 576 9.0913388)
-(0.412203210667 0 9.4410057)
-(0.412203210667 576 9.4410057)
-(6.5429081      0 149.8572332)
-(6.5429081      576 149.8572332)
+	(0 0 0)
+	(0 576 0)
+
+	(-0.1570298 0 3.5965736)
+        (-0.1570298 576 3.5965736)
+	(-0.1679346 0 3.8463357)
+        (-0.1679346 576 3.8463357)
+	(-0.3969364 0 9.0913388)
+        (-0.3969364 576 9.0913388)
+	(-0.4122032 0 9.4410057)
+	(-0.4122032 576 9.4410057)
+	(-6.5429081 0 149.8572332)
+        (-6.5429081 576 149.8572332)
+
+        (0.1570298 0 3.5965736)
+        (0.1570298 576 3.5965736)
+        (0.1679346 0 3.8463357)
+        (0.1679346 576 3.8463357)
+        (0.3969364 0 9.0913388)
+        (0.3969364 576 9.0913388)
+        (0.4122032 0 9.4410057)
+        (0.4122032 576 9.4410057)
+        (6.5429081 0 149.8572332)
+        (6.5429081 576 149.8572332)
 );
 
 blocks

examples/FlameD_Equilibrium/constant/polyMesh/boundary

@@ -0,0 +1,76 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  1.7.x                                 |
+|   \\  /    A nd           | Web:      www.OpenFOAM.com                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       polyBoundaryMesh;
+    location    "constant/polyMesh";
+    object      boundary;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+9
+(
+    axis
+    {
+        type            symmetryPlane;
+        nFaces          0;
+        startFace       31843;
+    }
+    mainjet
+    {
+        type            patch;
+        nFaces          10;
+        startFace       31843;
+    }
+    pilot
+    {
+        type            patch;
+        nFaces          15;
+        startFace       31853;
+    }
+    coflow
+    {
+        type            patch;
+        nFaces          80;
+        startFace       31868;
+    }
+    outlet
+    {
+        type            patch;
+        nFaces          107;
+        startFace       31948;
+    }
+    burnerwall
+    {
+        type            wall;
+        nFaces          2;
+        startFace       32055;
+    }
+    slipwall
+    {
+        type            wall;
+        nFaces          150;
+        startFace       32057;
+    }
+    front
+    {
+        type            wedge;
+        nFaces          16050;
+        startFace       32207;
+    }
+    back
+    {
+        type            wedge;
+        nFaces          16050;
+        startFace       48257;
+    }
+)
+
+// ************************************************************************* //

examples/FlameD_LagrangianCMC/constant/polyMesh/blockMeshDict

@@ -19,28 +19,30 @@ convertToMeters 0.001;
 vertices
 (
 
-(0 0 0)
-(0 576 0)
-(-0.15702979454  0 3.5965736)
-(-0.15702979454  576 3.5965736)
-(-0.167934643323 0 3.8463357)
-(-0.167934643323 576 3.8463357)
-(-0.396936424116 0 9.0913388)
-(-0.396936424116 576 9.0913388)
-(-0.412203210667 0 9.4410057)
-(-0.412203210667 576 9.4410057)
-(-6.5429081      0 149.8572332)
-(-6.5429081      576 149.8572332)
-(0.15702979454  0 3.5965736)
-(0.15702979454  576 3.5965736)
-(0.167934643323 0 3.8463357)
-(0.167934643323 576 3.8463357)
-(0.396936424116 0 9.0913388)
-(0.396936424116 576 9.0913388)
-(0.412203210667 0 9.4410057)
-(0.412203210667 576 9.4410057)
-(6.5429081      0 149.8572332)
-(6.5429081      576 149.8572332)
+	(0 0 0)
+	(0 576 0)
+
+	(-0.1570298 0 3.5965736)
+        (-0.1570298 576 3.5965736)
+	(-0.1679346 0 3.8463357)
+        (-0.1679346 576 3.8463357)
+	(-0.3969364 0 9.0913388)
+        (-0.3969364 576 9.0913388)
+	(-0.4122032 0 9.4410057)
+	(-0.4122032 576 9.4410057)
+	(-6.5429081 0 149.8572332)
+        (-6.5429081 576 149.8572332)
+
+        (0.1570298 0 3.5965736)
+        (0.1570298 576 3.5965736)
+        (0.1679346 0 3.8463357)
+        (0.1679346 576 3.8463357)
+        (0.3969364 0 9.0913388)
+        (0.3969364 576 9.0913388)
+        (0.4122032 0 9.4410057)
+        (0.4122032 576 9.4410057)
+        (6.5429081 0 149.8572332)
+        (6.5429081 576 149.8572332)
 );
 
 blocks

examples/FlameD_LagrangianCMC/constant/polyMesh/boundary

@@ -0,0 +1,76 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  1.7.x                                 |
+|   \\  /    A nd           | Web:      www.OpenFOAM.com                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       polyBoundaryMesh;
+    location    "constant/polyMesh";
+    object      boundary;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+9
+(
+    axis
+    {
+        type            symmetryPlane;
+        nFaces          0;
+        startFace       31843;
+    }
+    mainjet
+    {
+        type            patch;
+        nFaces          10;
+        startFace       31843;
+    }
+    pilot
+    {
+        type            patch;
+        nFaces          15;
+        startFace       31853;
+    }
+    coflow
+    {
+        type            patch;
+        nFaces          80;
+        startFace       31868;
+    }
+    outlet
+    {
+        type            patch;
+        nFaces          107;
+        startFace       31948;
+    }
+    burnerwall
+    {
+        type            wall;
+        nFaces          2;
+        startFace       32055;
+    }
+    slipwall
+    {
+        type            wall;
+        nFaces          150;
+        startFace       32057;
+    }
+    front
+    {
+        type            wedge;
+        nFaces          16050;
+        startFace       32207;
+    }
+    back
+    {
+        type            wedge;
+        nFaces          16050;
+        startFace       48257;
+    }
+)
+
+// ************************************************************************* //

examples/FlameD_LagrangianCMC/constant/thermophysicalProperties

@@ -21,7 +21,7 @@ thermoType
     mixture         reactingMixture;
     transport       sutherland;
     thermo          janaf;
-    energy          absoluteEnthalpy;
+    energy          sensibleEnthalpy;
     equationOfState perfectGas;
     specie          specie;
 }

examples/FlameD_LagrangianCMC/system/fvSchemes

@@ -33,7 +33,7 @@ divSchemes
     div(phi,U)      Gauss linearUpwindV cellLimited Gauss linear 1;
     div(phi,Yi_h)   Gauss multivariateSelection
     {
-	ha linearUpwind cellLimited Gauss linear 1;
+	h linearUpwind cellLimited Gauss linear 1;
 	H2 linearUpwind cellLimited Gauss linear 1;
 	H linearUpwind cellLimited Gauss linear 1;
 	O linearUpwind cellLimited Gauss linear 1;

examples/FlameD_LagrangianCMC/system/fvSolution

@@ -46,14 +46,14 @@ solvers
         relTol           0.0;
     }
 
-    "(U|Yi|ha|k|epsilon|mf|mfVar|Fk)"
+    "(U|Yi|h|k|epsilon|mf|mfVar|Fk)"
     {
         solver          PBiCG;
         preconditioner  DILU;
         tolerance       1e-08;
         relTol          0;
     }
-    "(U|Yi|ha|k|epsilon|mf|mfVar|Fk)Final"
+    "(U|Yi|h|k|epsilon|mf|mfVar|Fk)Final"
     {
         solver          PBiCG;
         preconditioner  DILU;

examples/FlameD_SLFM/constant/polyMesh/blockMeshDict

@@ -19,28 +19,30 @@ convertToMeters 0.001;
 vertices
 (
 
-(0 0 0)
-(0 576 0)
-(-0.15702979454  0 3.5965736)
-(-0.15702979454  576 3.5965736)
-(-0.167934643323 0 3.8463357)
-(-0.167934643323 576 3.8463357)
-(-0.396936424116 0 9.0913388)
-(-0.396936424116 576 9.0913388)
-(-0.412203210667 0 9.4410057)
-(-0.412203210667 576 9.4410057)
-(-6.5429081      0 149.8572332)
-(-6.5429081      576 149.8572332)
-(0.15702979454  0 3.5965736)
-(0.15702979454  576 3.5965736)
-(0.167934643323 0 3.8463357)
-(0.167934643323 576 3.8463357)
-(0.396936424116 0 9.0913388)
-(0.396936424116 576 9.0913388)
-(0.412203210667 0 9.4410057)
-(0.412203210667 576 9.4410057)
-(6.5429081      0 149.8572332)
-(6.5429081      576 149.8572332)
+	(0 0 0)
+	(0 576 0)
+
+	(-0.1570298 0 3.5965736)
+        (-0.1570298 576 3.5965736)
+	(-0.1679346 0 3.8463357)
+        (-0.1679346 576 3.8463357)
+	(-0.3969364 0 9.0913388)
+        (-0.3969364 576 9.0913388)
+	(-0.4122032 0 9.4410057)
+	(-0.4122032 576 9.4410057)
+	(-6.5429081 0 149.8572332)
+        (-6.5429081 576 149.8572332)
+
+        (0.1570298 0 3.5965736)
+        (0.1570298 576 3.5965736)
+        (0.1679346 0 3.8463357)
+        (0.1679346 576 3.8463357)
+        (0.3969364 0 9.0913388)
+        (0.3969364 576 9.0913388)
+        (0.4122032 0 9.4410057)
+        (0.4122032 576 9.4410057)
+        (6.5429081 0 149.8572332)
+        (6.5429081 576 149.8572332)
 );
 
 blocks

examples/FlameD_SLFM/constant/polyMesh/boundary

@@ -0,0 +1,76 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  1.7.x                                 |
+|   \\  /    A nd           | Web:      www.OpenFOAM.com                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       polyBoundaryMesh;
+    location    "constant/polyMesh";
+    object      boundary;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+9
+(
+    axis
+    {
+        type            symmetryPlane;
+        nFaces          0;
+        startFace       31843;
+    }
+    mainjet
+    {
+        type            patch;
+        nFaces          10;
+        startFace       31843;
+    }
+    pilot
+    {
+        type            patch;
+        nFaces          15;
+        startFace       31853;
+    }
+    coflow
+    {
+        type            patch;
+        nFaces          80;
+        startFace       31868;
+    }
+    outlet
+    {
+        type            patch;
+        nFaces          107;
+        startFace       31948;
+    }
+    burnerwall
+    {
+        type            wall;
+        nFaces          2;
+        startFace       32055;
+    }
+    slipwall
+    {
+        type            wall;
+        nFaces          150;
+        startFace       32057;
+    }
+    front
+    {
+        type            wedge;
+        nFaces          16050;
+        startFace       32207;
+    }
+    back
+    {
+        type            wedge;
+        nFaces          16050;
+        startFace       48257;
+    }
+)
+
+// ************************************************************************* //

libs/Allwmake

@@ -5,6 +5,6 @@ targetType=libso
 . $WM_PROJECT_DIR/wmake/scripts/AllwmakeParseArguments
 
 wmake $targetType chemistryModel_POSTECH
-wmake $targetType thermos
+wmake $targetType combustionModels_POSTECH
 
 #------------------------------------------------------------------------------

libs/combustionModels_POSTECH/FSD/FSD.C

@@ -0,0 +1,370 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-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 <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "FSD.H"
+#include "addToRunTimeSelectionTable.H"
+#include "LESModel.H"
+#include "fvcGrad.H"
+#include "fvcDiv.H"
+
+namespace Foam
+{
+namespace combustionModels
+{
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+template<class CombThermoType, class ThermoType>
+FSD<CombThermoType, ThermoType>::FSD
+(
+    const word& modelType,
+    const fvMesh& mesh,
+    const word& phaseName
+)
+:
+    singleStepCombustion<CombThermoType, ThermoType>
+    (
+        modelType,
+        mesh,
+        phaseName
+    ),
+    reactionRateFlameArea_
+    (
+        reactionRateFlameArea::New
+        (
+            this->coeffs(),
+            this->mesh(),
+            *this
+        )
+    ),
+    ft_
+    (
+        IOobject
+        (
+            IOobject::groupName("ft", phaseName),
+            this->mesh().time().timeName(),
+            this->mesh(),
+            IOobject::NO_READ,
+            IOobject::AUTO_WRITE
+        ),
+        this->mesh(),
+        dimensionedScalar("zero", dimless, 0.0)
+    ),
+    YFuelFuelStream_(dimensionedScalar("YFuelStream", dimless, 1.0)),
+    YO2OxiStream_(dimensionedScalar("YOxiStream", dimless, 0.23)),
+    Cv_(readScalar(this->coeffs().lookup("Cv"))),
+    C_(5.0),
+    ftMin_(0.0),
+    ftMax_(1.0),
+    ftDim_(300),
+    ftVarMin_(readScalar(this->coeffs().lookup("ftVarMin")))
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+template<class CombThermoType, class ThermoType>
+FSD<CombThermoType, ThermoType>::~FSD()
+{}
+
+
+// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
+
+template<class CombThermoType, class ThermoType>
+void FSD<CombThermoType, ThermoType>::calculateSourceNorm()
+{
+    this->singleMixturePtr_->fresCorrect();
+
+    const label fuelI = this->singleMixturePtr_->fuelIndex();
+
+    const volScalarField& YFuel = this->thermoPtr_->composition().Y()[fuelI];
+
+    const volScalarField& YO2 = this->thermoPtr_->composition().Y("O2");
+
+    const dimensionedScalar s = this->singleMixturePtr_->s();
+
+    ft_ =
+        (s*YFuel - (YO2 - YO2OxiStream_))/(s*YFuelFuelStream_ + YO2OxiStream_);
+
+
+    volVectorField nft(fvc::grad(ft_));
+
+    volScalarField mgft(mag(nft));
+
+    surfaceVectorField SfHat(this->mesh().Sf()/this->mesh().magSf());
+
+    volScalarField cAux(scalar(1) - ft_);
+
+    dimensionedScalar dMgft = 1.0e-3*
+        (ft_*cAux*mgft)().weightedAverage(this->mesh().V())
+       /((ft_*cAux)().weightedAverage(this->mesh().V()) + SMALL)
+      + dimensionedScalar("ddMgft", mgft.dimensions(), SMALL);
+
+    mgft += dMgft;
+
+    nft /= mgft;
+
+    const volVectorField& U = YO2.db().lookupObject<volVectorField>("U");
+
+    const volScalarField sigma
+    (
+        (nft & nft)*fvc::div(U) - (nft & fvc::grad(U) & nft)
+    );
+
+    reactionRateFlameArea_->correct(sigma);
+
+    const volScalarField& omegaFuel = reactionRateFlameArea_->omega();
+
+
+    const scalar ftStoich =
+        YO2OxiStream_.value()
+       /(
+            s.value()*YFuelFuelStream_.value() + YO2OxiStream_.value()
+        );
+
+    tmp<volScalarField> tPc
+    (
+        new volScalarField
+        (
+            IOobject
+            (
+                IOobject::groupName("Pc", this->phaseName_),
+                U.time().timeName(),
+                U.db(),
+                IOobject::NO_READ,
+                IOobject::NO_WRITE
+            ),
+            U.mesh(),
+            dimensionedScalar("Pc", dimless, 0)
+        )
+    );
+
+    volScalarField& pc = tPc.ref();
+
+    tmp<volScalarField> tomegaFuel
+    (
+        new volScalarField
+        (
+            IOobject
+            (
+                IOobject::groupName("omegaFuelBar", this->phaseName_),
+                U.time().timeName(),
+                U.db(),
+                IOobject::NO_READ,
+                IOobject::NO_WRITE
+            ),
+            U.mesh(),
+            dimensionedScalar
+            (
+                "omegaFuelBar",
+                omegaFuel.dimensions(),
+                0
+            )
+        )
+    );
+
+    volScalarField& omegaFuelBar = tomegaFuel.ref();
+
+    // Calculation of the mixture fraction variance (ftVar)
+    const compressible::LESModel& lesModel =
+        YO2.db().lookupObject<compressible::LESModel>
+        (
+            turbulenceModel::propertiesName
+        );
+
+    const volScalarField& delta = lesModel.delta();
+    const volScalarField ftVar(Cv_*sqr(delta)*sqr(mgft));
+
+    // Thickened flame (average flame thickness for counterflow configuration
+    // is 1.5 mm)
+
+    volScalarField  deltaF
+    (
+        lesModel.delta()/dimensionedScalar("flame", dimLength, 1.5e-3)
+    );
+
+    // Linear correlation between delta and flame thickness
+    volScalarField omegaF(max(deltaF*(4.0/3.0) + (2.0/3.0), scalar(1)));
+
+    scalar deltaFt = 1.0/ftDim_;
+
+    forAll(ft_, celli)
+    {
+        if (ft_[celli] > ftMin_ && ft_[celli] < ftMax_)
+        {
+            scalar ftCell = ft_[celli];
+
+            if (ftVar[celli] > ftVarMin_) //sub-grid beta pdf of ft_
+            {
+                scalar ftVarc = ftVar[celli];
+                scalar a =
+                    max(ftCell*(ftCell*(1.0 - ftCell)/ftVarc - 1.0), 0.0);
+                scalar b = max(a/ftCell - a, 0.0);
+
+                for (int i=1; i<ftDim_; i++)
+                {
+                    scalar ft = i*deltaFt;
+                    pc[celli] += pow(ft, a-1.0)*pow(1.0 - ft, b - 1.0)*deltaFt;
+                }
+
+                for (int i=1; i<ftDim_; i++)
+                {
+                    scalar ft = i*deltaFt;
+                    omegaFuelBar[celli] +=
+                        omegaFuel[celli]/omegaF[celli]
+                       *exp
+                        (
+                           -sqr(ft - ftStoich)
+                           /(2.0*sqr(0.01*omegaF[celli]))
+                        )
+                       *pow(ft, a - 1.0)
+                       *pow(1.0 - ft, b - 1.0)
+                       *deltaFt;
+                }
+                omegaFuelBar[celli] /= max(pc[celli], 1e-4);
+            }
+            else
+            {
+                omegaFuelBar[celli] =
+                   omegaFuel[celli]/omegaF[celli]
+                  *exp(-sqr(ftCell - ftStoich)/(2.0*sqr(0.01*omegaF[celli])));
+            }
+        }
+        else
+        {
+            omegaFuelBar[celli] = 0.0;
+        }
+    }
+
+
+    // Combustion progress variable, c
+
+    List<label> productsIndex(2, label(-1));
+    {
+        label i = 0;
+        forAll(this->singleMixturePtr_->specieProd(), specieI)
+        {
+            if (this->singleMixturePtr_->specieProd()[specieI] < 0)
+            {
+                productsIndex[i] = specieI;
+                i++;
+            }
+        }
+    }
+
+
+    // Flamelet probability of the progress c based on IFC (reuse pc)
+    scalar YprodTotal = 0;
+    forAll(productsIndex, j)
+    {
+        YprodTotal += this->singleMixturePtr_->Yprod0()[productsIndex[j]];
+    }
+
+    forAll(ft_, celli)
+    {
+        if (ft_[celli] < ftStoich)
+        {
+            pc[celli] = ft_[celli]*(YprodTotal/ftStoich);
+        }
+        else
+        {
+            pc[celli] = (1.0 - ft_[celli])*(YprodTotal/(1.0 - ftStoich));
+        }
+    }
+
+    tmp<volScalarField> tproducts
+    (
+        new volScalarField
+        (
+            IOobject
+            (
+                IOobject::groupName("products", this->phaseName_),
+                U.time().timeName(),
+                U.db(),
+                IOobject::NO_READ,
+                IOobject::NO_WRITE
+            ),
+            U.mesh(),
+            dimensionedScalar("products", dimless, 0)
+        )
+    );
+
+    volScalarField& products = tproducts.ref();
+
+    forAll(productsIndex, j)
+    {
+        label specieI = productsIndex[j];
+        const volScalarField& Yp = this->thermoPtr_->composition().Y()[specieI];
+        products += Yp;
+    }
+
+    volScalarField c
+    (
+        max(scalar(1) - products/max(pc, scalar(1e-5)), scalar(0))
+    );
+
+    pc = min(C_*c, scalar(1));
+
+    const volScalarField fres(this->singleMixturePtr_->fres(fuelI));
+
+    this->wFuel_ == mgft*pc*omegaFuelBar;
+}
+
+
+template<class CombThermoType, class ThermoType>
+void FSD<CombThermoType, ThermoType>::correct()
+{
+    this->wFuel_ ==
+        dimensionedScalar("zero", dimMass/pow3(dimLength)/dimTime, 0.0);
+
+    if (this->active())
+    {
+        calculateSourceNorm();
+    }
+}
+
+
+template<class CombThermoType, class ThermoType>
+bool FSD<CombThermoType, ThermoType>::read()
+{
+    if (singleStepCombustion<CombThermoType, ThermoType>::read())
+    {
+        this->coeffs().lookup("Cv") >> Cv_ ;
+        this->coeffs().lookup("ftVarMin") >> ftVarMin_;
+        reactionRateFlameArea_->read(this->coeffs());
+        return true;
+    }
+    else
+    {
+        return false;
+    }
+}
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace combustionModels
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

libs/combustionModels_POSTECH/FSD/FSD.H

@@ -0,0 +1,174 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-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 <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::combustionModels::FSD
+
+Description
+
+    Flame Surface Dennsity (FDS) combustion model.
+
+    The fuel source term is given by mgft*pc*omegaFuelBar.
+
+    where:
+          mgft: filtered flame area.
+          pc:   probability of the combustion progress.
+          omegaFuelBar: filtered consumption speed per unit of flame area.
+
+    pc is considered from the IFC solution.
+    omegaFuelBar is calculated solving a relaxation equation which tends to
+    omegaEq. This omegaEq is obtained from the flamelet solution for
+    different strain rates and fit using a expential distribution.
+
+    The spacial distribution of the consumption speed (omega) is obtained also
+    from a strained flamelet solution and it is assumed to have a guassian
+    distribution.
+
+    If the grid resolution is not enough to resolve the flame, the consumption
+    speed distribution is linearly thickened conserving the overall heat
+    release.
+
+    If the turbulent fluctuation of the mixture fraction at the sub-grid level
+    is large (>1e-04) then a beta pdf is used for filtering.
+
+    At the moment the flame area combustion model is only fit to work in a LES
+    frame work. In RAS the subgrid fluctuation has to be solved by an extra
+    transport equation.
+
+SourceFiles
+    FSD.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef FSD_H
+#define FSD_H
+
+#include "singleStepCombustion.H"
+#include "reactionRateFlameArea.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace combustionModels
+{
+
+/*---------------------------------------------------------------------------*\
+                          Class FSD Declaration
+\*---------------------------------------------------------------------------*/
+
+template<class CombThermoType, class ThermoType>
+class FSD
+:
+    public singleStepCombustion <CombThermoType, ThermoType>
+{
+    // Private data
+
+        //- Auto pointer to consumption speed per unit of flame area model
+        autoPtr<reactionRateFlameArea> reactionRateFlameArea_;
+
+        //- Mixture fraction
+        volScalarField ft_;
+
+        //- Fuel mass concentration on the fuel stream
+        dimensionedScalar YFuelFuelStream_;
+
+        //- Oxygen mass concentration on the oxydizer stream
+        dimensionedScalar YO2OxiStream_;
+
+        //- Similarity constant for the sub-grid ft fluctuations
+        scalar Cv_;
+
+        //- Model constant
+        scalar C_;
+
+        //- Lower flammability limit
+        scalar ftMin_;
+
+        //- Upper flammability limit
+        scalar ftMax_;
+
+        //- Dimension of the ft space. Used to integrate the beta-pdf
+        scalar ftDim_;
+
+        //- Minimum mixture freaction variance to calculate pdf
+        scalar ftVarMin_;
+
+
+    // Private Member Functions
+
+        //- Calculate the normalised fuel source term
+        void calculateSourceNorm();
+
+        //- Disallow copy construct
+        FSD(const FSD&);
+
+        //- Disallow default bitwise assignment
+        void operator=(const FSD&);
+
+
+public:
+
+    //- Runtime type information
+    TypeName("FSD");
+
+
+    // Constructors
+
+        //- Construct from components
+        FSD(const word& modelType, const fvMesh& mesh, const word& phaseName);
+
+
+    // Destructor
+        virtual ~FSD();
+
+
+    // Evolution
+
+        //- Correct combustion rate
+        virtual void correct();
+
+
+    // IO
+
+        //- Update properties
+        virtual bool read();
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace combustionModels
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+    #include "FSD.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

libs/combustionModels_POSTECH/FSD/FSDs.C

@@ -0,0 +1,104 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-2015 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 <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "makeCombustionTypes.H"
+
+#include "thermoPhysicsTypes.H"
+#include "psiCombustionModel.H"
+#include "psiThermoCombustion.H"
+#include "rhoCombustionModel.H"
+#include "rhoThermoCombustion.H"
+#include "FSD.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+// Combustion models based on sensibleEnthalpy
+makeCombustionTypesThermo
+(
+    FSD,
+    psiThermoCombustion,
+    gasHThermoPhysics,
+    psiCombustionModel
+);
+
+makeCombustionTypesThermo
+(
+    FSD,
+    psiThermoCombustion,
+    constGasHThermoPhysics,
+    psiCombustionModel
+);
+
+makeCombustionTypesThermo
+(
+    FSD,
+    rhoThermoCombustion,
+    gasHThermoPhysics,
+    rhoCombustionModel
+);
+
+makeCombustionTypesThermo
+(
+    FSD,
+    rhoThermoCombustion,
+    constGasHThermoPhysics,
+    rhoCombustionModel
+);
+
+// Combustion models based on sensibleInternalEnergy
+makeCombustionTypesThermo
+(
+    FSD,
+    psiThermoCombustion,
+    gasEThermoPhysics,
+    psiCombustionModel
+);
+
+makeCombustionTypesThermo
+(
+    FSD,
+    psiThermoCombustion,
+    constGasEThermoPhysics,
+    psiCombustionModel
+);
+
+makeCombustionTypesThermo
+(
+    FSD,
+    rhoThermoCombustion,
+    gasEThermoPhysics,
+    rhoCombustionModel
+);
+
+makeCombustionTypesThermo
+(
+    FSD,
+    rhoThermoCombustion,
+    constGasEThermoPhysics,
+    rhoCombustionModel
+);
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

libs/combustionModels_POSTECH/FSD/reactionRateFlameAreaModels/consumptionSpeed/consumptionSpeed.C

@@ -0,0 +1,142 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-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 <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "consumptionSpeed.H"
+
+/* * * * * * * * * * * * * * * private static data * * * * * * * * * * * * * */
+
+namespace Foam
+{
+    defineTypeNameAndDebug(consumptionSpeed, 0);
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::consumptionSpeed::consumptionSpeed
+(
+    const dictionary& dict
+)
+:   omega0_(readScalar(dict.lookup("omega0"))),
+    eta_(readScalar(dict.lookup("eta"))),
+    sigmaExt_(readScalar(dict.lookup("sigmaExt"))),
+    omegaMin_(readScalar(dict.lookup("omegaMin")))
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+Foam::consumptionSpeed::~consumptionSpeed()
+{}
+
+
+// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
+
+Foam::scalar Foam::consumptionSpeed::omega0Sigma
+(
+    scalar sigma,
+    scalar a
+) const
+{
+    scalar omega0 = 0.0;
+
+    if (sigma < sigmaExt_)
+    {
+        omega0 = max
+        (
+            a*omega0_*(1.0 - exp(eta_*max(sigma, 0.0))),
+            omegaMin_
+        ) ;
+    }
+
+    return omega0;
+}
+
+
+Foam::tmp<Foam::volScalarField> Foam::consumptionSpeed::omega0Sigma
+(
+    const volScalarField& sigma
+)
+{
+    tmp<volScalarField> tomega0
+    (
+        new volScalarField
+        (
+            IOobject
+            (
+                "omega0",
+                sigma.time().timeName(),
+                sigma.db(),
+                IOobject::NO_READ,
+                IOobject::NO_WRITE
+            ),
+            sigma.mesh(),
+            dimensionedScalar
+            (
+                "omega0",
+                dimensionSet(1, -2, -1, 0, 0, 0, 0),
+                0
+            )
+        )
+    );
+
+    volScalarField& omega0 = tomega0.ref();
+
+    volScalarField::Internal& iomega0 = omega0;
+
+    forAll(iomega0, celli)
+    {
+        iomega0[celli] = omega0Sigma(sigma[celli], 1.0);
+    }
+
+    volScalarField::Boundary& bomega0 = omega0.boundaryFieldRef();
+
+    forAll(bomega0, patchi)
+    {
+        forAll(bomega0[patchi], facei)
+        {
+            bomega0[patchi][facei] =
+                omega0Sigma
+                (
+                    sigma.boundaryField()[patchi][facei],
+                    1.0
+                );
+        }
+    }
+
+    return tomega0;
+}
+
+
+void  Foam::consumptionSpeed::read(const dictionary& dict)
+{
+    dict.lookup("omega0") >> omega0_ ;
+    dict.lookup("eta") >> eta_ ;
+    dict.lookup("sigmaExt") >> sigmaExt_;
+    dict.lookup("omegaMin") >> omegaMin_;
+}
+
+
+// ************************************************************************* //

libs/combustionModels_POSTECH/FSD/reactionRateFlameAreaModels/consumptionSpeed/consumptionSpeed.H

@@ -0,0 +1,139 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-2015 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 <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::reactionRateFlameAreaModels::consumptionSpeed
+
+Description
+    Correlation function for laminar consumption speed obtained from flamelet
+    solution at increasing strain rates.
+
+SourceFiles
+    consumptionSpeed.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef consumptionSpeed_H
+#define consumptionSpeed_H
+
+#include "IOdictionary.H"
+#include "volFields.H"
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                           Class consumptionSpeed Declaration
+\*---------------------------------------------------------------------------*/
+
+class consumptionSpeed
+{
+    // Private Data
+
+        //- Maximum consumption speed
+        scalar omega0_;
+
+        //- Exponential factor
+        scalar eta_;
+
+        //- Extinction strain
+        scalar sigmaExt_;
+
+        //- Minimum consumption speed
+        scalar omegaMin_;
+
+
+    // Private member functions
+
+        //- Return consumption rate
+        scalar omega0Sigma(scalar sigma, scalar a) const;
+
+        //- Disallow copy construct
+        consumptionSpeed(const consumptionSpeed&);
+
+         //- Disallow default bitwise assignment
+        void operator=(const consumptionSpeed&);
+
+
+public:
+
+    //- Runtime type information
+    TypeName("consumptionSpeed");
+
+
+    // Constructors
+
+        //- Construct from dictionary
+        consumptionSpeed(const dictionary& dict);
+
+
+    //- Destructor
+        virtual ~consumptionSpeed();
+
+
+    // Member functions
+
+        //- Return speed consumption rate temp
+        tmp<volScalarField> omega0Sigma(const volScalarField& sigma);
+
+
+        // Access functions
+
+            scalar omega0() const
+            {
+                return omega0_;
+            }
+
+            scalar eta() const
+            {
+                return eta_;
+            }
+
+            scalar sigmaExt() const
+            {
+                return sigmaExt_;
+            }
+
+            scalar omegaMin() const
+            {
+                return omegaMin_;
+            }
+
+
+    // IO
+
+        //- Update properties
+        void read(const dictionary& dict);
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

libs/combustionModels_POSTECH/FSD/reactionRateFlameAreaModels/reactionRateFlameArea/reactionRateFlameArea.C

@@ -0,0 +1,108 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-2015 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 <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "reactionRateFlameArea.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+    defineTypeNameAndDebug(reactionRateFlameArea, 0);
+    defineRunTimeSelectionTable(reactionRateFlameArea, dictionary);
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::reactionRateFlameArea::reactionRateFlameArea
+(
+    const dictionary& dict,
+    const fvMesh& mesh,
+    const combustionModel& combModel
+)
+:
+    coeffDict_(dictionary::null),
+    mesh_(mesh),
+    combModel_(combModel),
+    fuel_(dict.lookup("fuel")),
+    omega_
+    (
+        IOobject
+        (
+            "omega",
+            mesh_.time().timeName(),
+            mesh_,
+            IOobject::MUST_READ,
+            IOobject::AUTO_WRITE
+        ),
+        mesh_
+    )
+{}
+
+
+Foam::reactionRateFlameArea::reactionRateFlameArea
+(
+    const word& modelType,
+    const dictionary& dict,
+    const fvMesh& mesh,
+    const combustionModel& combModel
+)
+:
+    coeffDict_(dict.subDict(modelType + "Coeffs")),
+    mesh_(mesh),
+    combModel_(combModel),
+    fuel_(dict.lookup("fuel")),
+    omega_
+    (
+        IOobject
+        (
+            "omega",
+            mesh_.time().timeName(),
+            mesh_,
+            IOobject::MUST_READ,
+            IOobject::AUTO_WRITE
+        ),
+        mesh_
+    )
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+Foam::reactionRateFlameArea::~reactionRateFlameArea()
+{}
+
+
+// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
+
+bool Foam::reactionRateFlameArea::read(const dictionary& dict)
+{
+    dict.lookup("fuel") >> fuel_;
+
+    return true;
+}
+
+
+// ************************************************************************* //

libs/combustionModels_POSTECH/FSD/reactionRateFlameAreaModels/reactionRateFlameArea/reactionRateFlameArea.H

@@ -0,0 +1,175 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-2013 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 <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::reactionRateFlameArea
+
+Description
+    Abstract class for reaction rate per flame area unit
+
+SourceFiles
+    reactionRateFlameArea.C
+    reactionRateFlameAreaNew.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef reactionRateFlameArea_H
+#define reactionRateFlameArea_H
+
+#include "runTimeSelectionTables.H"
+#include "dictionary.H"
+#include "autoPtr.H"
+#include "volFields.H"
+#include "combustionModel.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+class fvMesh;
+
+/*---------------------------------------------------------------------------*\
+                           Class reactionRateFlameArea Declaration
+\*---------------------------------------------------------------------------*/
+
+class reactionRateFlameArea
+{
+
+protected:
+
+    // Protected data
+
+        //- Dictionary
+        dictionary coeffDict_;
+
+        //- Mesh reference
+        const fvMesh& mesh_;
+
+        //- Combstion model owner
+        const combustionModel& combModel_;
+
+        //- Fuel name
+        word fuel_;
+
+        //- Fuel consumption rate per unit of flame area
+        volScalarField omega_;
+
+
+private:
+
+    // Private member functions
+
+         //- Disallow copy construct
+        reactionRateFlameArea(const reactionRateFlameArea&);
+
+         //- Disallow default bitwise assignment
+        void operator=(const reactionRateFlameArea&);
+
+
+public:
+
+    //- Runtime type information
+    TypeName("reactionRateFlameArea");
+
+
+    // Declare run-time constructor selection table
+
+        declareRunTimeSelectionTable
+        (
+            autoPtr,
+            reactionRateFlameArea,
+            dictionary,
+            (
+                const word modelType,
+                const dictionary& dict,
+                const fvMesh& mesh,
+                const combustionModel& combModel
+            ),
+            (modelType, dict, mesh, combModel)
+        );
+
+
+    // Constructors
+
+        //- Construct from dictionary and psiReactionThermo
+        reactionRateFlameArea
+        (
+            const dictionary& dict,
+            const fvMesh& mesh,
+            const combustionModel& combModel
+        );
+
+        //- Construct from components
+        reactionRateFlameArea
+        (
+            const word& modelType,
+            const dictionary& dict,
+            const fvMesh& mesh,
+            const combustionModel& combModel
+        );
+
+
+    // Selector
+
+        static autoPtr<reactionRateFlameArea> New
+        (
+            const dictionary& dict,
+            const fvMesh& mesh,
+            const combustionModel& combModel
+        );
+
+
+    // Destructor
+
+        virtual ~reactionRateFlameArea();
+
+
+    // Member functions
+
+        //- Access functions
+
+            //- Return omega
+            const volScalarField& omega() const
+            {
+                return omega_;
+            }
+
+
+        //- Correct omega
+        virtual void correct(const volScalarField& sigma) = 0;
+
+        //- Update from dictionary
+        virtual bool read(const dictionary& dictProperties);
+
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

libs/combustionModels_POSTECH/FSD/reactionRateFlameAreaModels/reactionRateFlameArea/reactionRateFlameAreaNew.C

@@ -0,0 +1,69 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-2015 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 <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "reactionRateFlameArea.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+Foam::autoPtr<Foam::reactionRateFlameArea> Foam::reactionRateFlameArea::New
+(
+    const dictionary& dict,
+    const fvMesh& mesh,
+    const combustionModel& combModel
+)
+{
+    word reactionRateFlameAreaType
+    (
+        dict.lookup("reactionRateFlameArea")
+    );
+
+    Info<< "Selecting reaction rate flame area correlation "
+        << reactionRateFlameAreaType << endl;
+
+    dictionaryConstructorTable::iterator cstrIter =
+        dictionaryConstructorTablePtr_->find(reactionRateFlameAreaType);
+
+    if (cstrIter == dictionaryConstructorTablePtr_->end())
+    {
+        FatalIOErrorInFunction
+        (
+            dict
+        )   << "Unknown reactionRateFlameArea type "
+            << reactionRateFlameAreaType << endl << endl
+            << "Valid reaction rate flame area types are :" << endl
+            << dictionaryConstructorTablePtr_->toc()
+            << exit(FatalIOError);
+    }
+
+    const label tempOpen = reactionRateFlameAreaType.find('<');
+
+    const word className = reactionRateFlameAreaType(0, tempOpen);
+
+    return autoPtr<reactionRateFlameArea>
+        (cstrIter()(className, dict, mesh, combModel));
+}
+
+
+// ************************************************************************* //

libs/combustionModels_POSTECH/FSD/reactionRateFlameAreaModels/relaxation/relaxation.C

@@ -0,0 +1,165 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-2015 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 <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "relaxation.H"
+#include "addToRunTimeSelectionTable.H"
+#include "fvm.H"
+#include "LESModel.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace reactionRateFlameAreaModels
+{
+    defineTypeNameAndDebug(relaxation, 0);
+    addToRunTimeSelectionTable
+    (
+        reactionRateFlameArea,
+        relaxation,
+        dictionary
+    );
+}
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::reactionRateFlameAreaModels::relaxation::relaxation
+(
+    const word modelType,
+    const dictionary& dict,
+    const fvMesh& mesh,
+    const combustionModel& combModel
+)
+:
+    reactionRateFlameArea(modelType, dict, mesh, combModel),
+    correlation_(dict.subDict(typeName + "Coeffs").subDict(fuel_)),
+    C_(readScalar(dict.subDict(typeName + "Coeffs").lookup("C"))),
+    alpha_(readScalar(dict.subDict(typeName + "Coeffs").lookup("alpha")))
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+Foam::reactionRateFlameAreaModels::relaxation::~relaxation()
+{}
+
+
+// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
+
+void Foam::reactionRateFlameAreaModels::relaxation::correct
+(
+    const volScalarField& sigma
+)
+{
+    dimensionedScalar omega0
+    (
+        "omega0",
+        dimensionSet(1, -2, -1, 0, 0, 0, 0),
+        correlation_.omega0()
+    );
+
+    dimensionedScalar sigmaExt
+    (
+        "sigmaExt",
+        dimensionSet(0, 0, -1, 0, 0, 0, 0),
+        correlation_.sigmaExt()
+    );
+
+    dimensionedScalar omegaMin
+    (
+        "omegaMin",
+        omega0.dimensions(),
+        1e-4
+    );
+
+    dimensionedScalar kMin
+    (
+        "kMin",
+        sqr(dimVelocity),
+        SMALL
+    );
+
+    const compressibleTurbulenceModel& turbulence = combModel_.turbulence();
+
+    // Total strain
+    const volScalarField sigmaTotal
+    (
+        sigma + alpha_*turbulence.epsilon()/(turbulence.k() + kMin)
+    );
+
+    const volScalarField omegaInf(correlation_.omega0Sigma(sigmaTotal));
+
+    dimensionedScalar sigma0("sigma0", sigma.dimensions(), 0.0);
+
+    const volScalarField tau(C_*mag(sigmaTotal));
+
+    volScalarField Rc
+    (
+        (tau*omegaInf*(omega0 - omegaInf) + sqr(omegaMin)*sigmaExt)
+       /(sqr(omega0 - omegaInf) + sqr(omegaMin))
+    );
+
+    const volScalarField rho(combModel_.rho());
+    const surfaceScalarField phi(combModel_.phi());
+
+    solve
+    (
+         fvm::ddt(rho, omega_)
+       + fvm::div(phi, omega_, "div(phi,omega)")
+      ==
+         rho*Rc*omega0
+       - fvm::SuSp(rho*(tau + Rc), omega_)
+    );
+
+    omega_.min(omega0);
+    omega_.max(0.0);
+}
+
+
+bool  Foam::reactionRateFlameAreaModels::relaxation::read
+(
+    const dictionary& dict
+)
+{
+    if (reactionRateFlameArea::read(dict))
+    {
+        coeffDict_ = dict.subDict(typeName + "Coeffs");
+        coeffDict_.lookup("C") >> C_;
+        coeffDict_.lookup("alpha") >> alpha_;
+        correlation_.read
+        (
+            coeffDict_.subDict(fuel_)
+        );
+        return true;
+    }
+    else
+    {
+        return false;
+    }
+}
+
+// ************************************************************************* //

libs/combustionModels_POSTECH/FSD/reactionRateFlameAreaModels/relaxation/relaxation.H

@@ -0,0 +1,122 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-2015 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 <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::reactionRateFlameAreaModels::relaxation
+
+Description
+    Consumption rate per unit of flame area obtained from a relaxation equation
+
+SourceFiles
+    relaxation.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef relaxation_H
+#define relaxation_H
+
+#include "reactionRateFlameArea.H"
+#include "consumptionSpeed.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace reactionRateFlameAreaModels
+{
+
+/*---------------------------------------------------------------------------*\
+                           Class relaxation Declaration
+\*---------------------------------------------------------------------------*/
+
+class relaxation
+:
+    public reactionRateFlameArea
+{
+    // Private Data
+
+        //- Correlation
+        consumptionSpeed correlation_;
+
+        //- Proportionality constant for time scale in the relaxation Eq.
+        scalar C_;
+
+        //- Proportionality constant for sub-grid strain
+        scalar alpha_;
+
+
+     // Private Member Functions
+
+        //- Disallow copy construct
+        relaxation(const relaxation&);
+
+        //- Disallow default bitwise assignment
+        void operator=(const relaxation&);
+
+
+public:
+
+    //- Runtime type information
+    TypeName("relaxation");
+
+
+    // Constructors
+
+        //- Construct from dictionary and psiReactionThermo
+        relaxation
+        (
+            const word modelType,
+            const dictionary& dictCoeffs,
+            const fvMesh& mesh,
+            const combustionModel& combModel
+        );
+
+
+    // Destructor
+
+        virtual ~relaxation();
+
+
+    // Member functions
+
+        //- Correct omega
+        virtual void correct(const volScalarField& sigma);
+
+
+    // IO
+
+        //- Update properties from given dictionary
+        virtual bool read(const dictionary& dictProperties);
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End reactionRateFlameAreaModels
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

libs/combustionModels_POSTECH/Make/files

@@ -0,0 +1,29 @@
+combustionModel/combustionModel.C
+
+psiCombustionModel/psiCombustionModel/psiCombustionModel.C
+psiCombustionModel/psiCombustionModel/psiCombustionModelNew.C
+psiCombustionModel/psiThermoCombustion/psiThermoCombustion.C
+psiCombustionModel/psiChemistryCombustion/psiChemistryCombustion.C
+
+rhoCombustionModel/rhoCombustionModel/rhoCombustionModel.C
+rhoCombustionModel/rhoCombustionModel/rhoCombustionModelNew.C
+rhoCombustionModel/rhoThermoCombustion/rhoThermoCombustion.C
+rhoCombustionModel/rhoChemistryCombustion/rhoChemistryCombustion.C
+
+diffusion/diffusions.C
+infinitelyFastChemistry/infinitelyFastChemistrys.C
+
+PaSR/PaSRs.C
+
+laminar/laminars.C
+
+FSD/reactionRateFlameAreaModels/consumptionSpeed/consumptionSpeed.C
+FSD/reactionRateFlameAreaModels/reactionRateFlameArea/reactionRateFlameArea.C
+FSD/reactionRateFlameAreaModels/reactionRateFlameArea/reactionRateFlameAreaNew.C
+FSD/reactionRateFlameAreaModels/relaxation/relaxation.C
+
+FSD/FSDs.C
+
+noCombustion/noCombustions.C
+
+LIB = $(FOAM_USER_LIBBIN)/libcombustionModels_POSTECH

libs/combustionModels_POSTECH/Make/options

@@ -0,0 +1,18 @@
+EXE_INC = \
+    -I$(LIB_SRC)/transportModels/compressible/lnInclude \
+    -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
+    -I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \
+    -I$(LIB_SRC)/thermophysicalModels/reactionThermo/lnInclude \
+    -I$(LIB_SRC)/thermophysicalModels/chemistryModel/lnInclude \
+    -I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
+    -I$(LIB_SRC)/TurbulenceModels/compressible/lnInclude \
+    -I$(LIB_SRC)/finiteVolume/lnInclude \
+    -I$(LIB_SRC)/meshTools/lnInclude
+
+LIB_LIBS = \
+    -lcompressibleTransportModels \
+    -lturbulenceModels \
+    -lcompressibleTurbulenceModels \
+    -lchemistryModel \
+    -lfiniteVolume \
+    -lmeshTools

libs/combustionModels_POSTECH/PaSR/PaSR.C

@@ -0,0 +1,162 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-2015 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 <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "PaSR.H"
+#include "fvmSup.H"
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+template<class Type>
+Foam::combustionModels::PaSR<Type>::PaSR
+(
+    const word& modelType,
+    const fvMesh& mesh,
+    const word& phaseName
+)
+:
+    laminar<Type>(modelType, mesh, phaseName),
+    Cmix_(readScalar(this->coeffs().lookup("Cmix"))),
+    turbulentReaction_(this->coeffs().lookup("turbulentReaction")),
+    kappa_
+    (
+        IOobject
+        (
+            IOobject::groupName("PaSR:kappa", phaseName),
+            mesh.time().timeName(),
+            mesh,
+            IOobject::NO_READ,
+            IOobject::AUTO_WRITE
+        ),
+        mesh,
+        dimensionedScalar("kappa", dimless, 0.0)
+    )
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+template<class Type>
+Foam::combustionModels::PaSR<Type>::~PaSR()
+{}
+
+
+// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
+
+template<class Type>
+void Foam::combustionModels::PaSR<Type>::correct()
+{
+    if (this->active())
+    {
+        laminar<Type>::correct();
+
+        if (turbulentReaction_)
+        {
+            tmp<volScalarField> tepsilon(this->turbulence().epsilon());
+            const volScalarField& epsilon = tepsilon();
+            tmp<volScalarField> tmuEff(this->turbulence().muEff());
+            const volScalarField& muEff = tmuEff();
+            tmp<volScalarField> ttc(this->tc());
+            const volScalarField& tc = ttc();
+            tmp<volScalarField> trho(this->rho());
+            const volScalarField& rho = trho();
+
+            forAll(epsilon, i)
+            {
+                scalar tk =
+                    Cmix_*sqrt(max(muEff[i]/rho[i]/(epsilon[i] + SMALL), 0));
+
+                if (tk > SMALL)
+                {
+                    kappa_[i] = tc[i]/(tc[i] + tk);
+                }
+                else
+                {
+                    kappa_[i] = 1.0;
+                }
+            }
+        }
+        else
+        {
+            kappa_ = 1.0;
+        }
+    }
+}
+
+
+template<class Type>
+Foam::tmp<Foam::fvScalarMatrix>
+Foam::combustionModels::PaSR<Type>::R(volScalarField& Y) const
+{
+    return kappa_*laminar<Type>::R(Y);
+}
+
+
+template<class Type>
+Foam::tmp<Foam::volScalarField>
+Foam::combustionModels::PaSR<Type>::dQ() const
+{
+    return tmp<volScalarField>
+    (
+        new volScalarField
+        (
+            IOobject::groupName("PaSR:dQ", this->phaseName_),
+            kappa_*laminar<Type>::dQ()
+        )
+    );
+}
+
+
+template<class Type>
+Foam::tmp<Foam::volScalarField>
+Foam::combustionModels::PaSR<Type>::Sh() const
+{
+    return tmp<volScalarField>
+    (
+        new volScalarField
+        (
+            IOobject::groupName("PaSR:Sh", this->phaseName_),
+            kappa_*laminar<Type>::Sh()
+        )
+    );
+}
+
+
+template<class Type>
+bool Foam::combustionModels::PaSR<Type>::read()
+{
+    if (laminar<Type>::read())
+    {
+        this->coeffs().lookup("Cmix") >> Cmix_;
+        this->coeffs().lookup("turbulentReaction") >> turbulentReaction_;
+        return true;
+    }
+    else
+    {
+        return false;
+    }
+}
+
+
+// ************************************************************************* //

libs/combustionModels_POSTECH/PaSR/PaSR.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2016 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -22,137 +22,111 @@ License
     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
 
 Class
-    Foam::BetaGrid
+    Foam::combustionModels::PaSR
 
 Description
+    Partially stirred reactor combustion model.  The model calculates a finite
+    rate, based on both turbulence and chemistry time scales.  Depending on
+    mesh resolution, the Cmix parameter can be used to scale the turbulence
+    mixing time scale.
 
 SourceFiles
-    BetaGridI.H
-    BetaGrid.C
-    BetaGridIO.C
+    PaSR.C
 
 \*---------------------------------------------------------------------------*/
 
-#ifndef BetaGrid_H
-#define BetaGrid_H
+#ifndef PaSR_H
+#define PaSR_H
 
-#include "scalarField.H"
-#include "labelList.H"
+#include "../laminar/laminar.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 namespace Foam
 {
-
-// Forward declaration of classes
-class Istream;
-class Ostream;
-
-// Forward declaration of friend functions and operators
-namespace CMC
+namespace combustionModels
 {
-    class BetaGrid;
-}
-class dictionary;
-
-Istream& operator>>(Istream&, CMC::BetaGrid&);
-Ostream& operator<<(Ostream&, const CMC::BetaGrid&);
-
 
 /*---------------------------------------------------------------------------*\
-                         Class BetaGrid Declaration
+                            Class PaSR Declaration
 \*---------------------------------------------------------------------------*/
 
-namespace CMC
-{
-
-class BetaGrid
+template<class Type>
+class PaSR
+:
+    public laminar<Type>
 {
     // Private data
 
-        //- Piecewise uniform interval points
-        scalarField etaCut_;
+        //- Mixing constant
+        scalar Cmix_;
 
-        //- Sizes of piecewise uniform intervals
-        labelList nSpacings_;
+        //- Turbulent reaction switch
+        Switch turbulentReaction_;
 
-        //- Mixture fraction eta grid points
-        scalarField etaSpace_;
+        //- Mixing parameter
+        volScalarField kappa_;
 
 
     // Private Member Functions
 
-        //- Disallow default bitwise copy construct
-        BetaGrid(const BetaGrid&);
+        //- Disallow copy construct
+        PaSR(const PaSR&);
 
         //- Disallow default bitwise assignment
-        void operator=(const BetaGrid&);
+        void operator=(const PaSR&);
 
 
 public:
 
-    // Static data members
-
-        //- Static data staticData
-        // static const scalar staticData;
+    //- Runtime type information
+    TypeName("PaSR");
 
 
     // Constructors
 
         //- Construct from components
-        BetaGrid(const dictionary& dict);
-
-        //- Construct from Istream
-        BetaGrid(Istream&);
-
-        //- Construct as copy
-        // BetaGrid(const BetaGrid&);
+        PaSR(const word& modelType, const fvMesh& mesh, const word& phaseName);
 
 
     //- Destructor
-    ~BetaGrid();
+    virtual ~PaSR();
 
 
     // Member Functions
 
-        // Access
-        const labelList& N() const {return nSpacings_;}
+        // Evolution
 
-        const scalarField& etaCut() const {return etaCut_;}
+            //- Correct combustion rate
+            virtual void correct();
 
-        const scalarField& etaSpace() const {return etaSpace_;}
+            //- Fuel consumption rate matrix.
+            virtual tmp<fvScalarMatrix> R(volScalarField& Y) const;
 
-        // Check
+            //- Heat release rate calculated from fuel consumption rate matrix
+            virtual tmp<volScalarField> dQ() const;
 
-        // Edit
-
-        // Write
+            //-  Return source for enthalpy equation [kg/m/s3]
+            virtual tmp<volScalarField> Sh() const;
 
 
-    // Member Operators
+    // IO
 
-        // void operator=(const BetaGrid&);
-
-
-    // Friend Functions
-
-    // Friend Operators
-
-    // IOstream Operators
-
-        friend Istream& ::Foam::operator>>(Istream&, BetaGrid&);
-        friend Ostream& ::Foam::operator<<(Ostream&, const BetaGrid&);
+            //- Update properties from given dictionary
+            virtual bool read();
 };
 
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
-} // End namespace CMC
+} // End namespace combustionModels
 } // End namespace Foam
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
-#include "BetaGridI.H"
+#ifdef NoRepository
+    #include "PaSR.C"
+#endif
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 

libs/combustionModels_POSTECH/PaSR/PaSRs.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2015 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -23,34 +23,16 @@ License
 
 \*---------------------------------------------------------------------------*/
 
-#include "BetaIntegrator.H"
-#include "IOstreams.H"
+#include "makeCombustionTypes.H"
 
-// * * * * * * * * * * * * * * * IOstream Operators  * * * * * * * * * * * * //
+#include "psiChemistryCombustion.H"
+#include "rhoChemistryCombustion.H"
+#include "PaSR.H"
 
-Foam::Istream& Foam::operator>>(Istream& is, CMC::BetaIntegrator&)
-{
-    // Check state of Istream
-    is.check
-    (
-        "Foam::Istream& Foam::operator>>(Foam::Istream&, Foam::CMC::BetaIntegrator&)"
-    );
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
-    return is;
-}
-
-
-Foam::Ostream& Foam::operator<<(Ostream& os, const CMC::BetaIntegrator&)
-{
-    // Check state of Ostream
-    os.check
-    (
-        "Foam::Ostream& Foam::operator<<(Foam::Ostream&, "
-        "const Foam::CMC::BetaIntegrator&)"
-    );
-
-    return os;
-}
+makeCombustionTypes(PaSR, psiChemistryCombustion, psiCombustionModel);
+makeCombustionTypes(PaSR, rhoChemistryCombustion, rhoCombustionModel);
 
 
 // ************************************************************************* //

libs/combustionModels_POSTECH/combustionModel/combustionModel.C

@@ -0,0 +1,114 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-2015 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 <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "combustionModel.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+    defineTypeNameAndDebug(combustionModel, 0);
+}
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::combustionModel::combustionModel
+(
+    const word& modelType,
+    const fvMesh& mesh,
+    const word& phaseName
+)
+:
+    IOdictionary
+    (
+        IOobject
+        (
+            IOobject::groupName("combustionProperties", phaseName),
+            mesh.time().constant(),
+            mesh,
+            IOobject::MUST_READ_IF_MODIFIED,
+            IOobject::NO_WRITE
+        )
+    ),
+    turbulencePtr_(),
+    mesh_(mesh),
+    active_(lookupOrDefault<Switch>("active", true)),
+    coeffs_(subDict(modelType + "Coeffs")),
+    modelType_(modelType),
+    phaseName_(phaseName)
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+Foam::combustionModel::~combustionModel()
+{
+    if (turbulencePtr_)
+    {
+        turbulencePtr_ = 0;
+    }
+}
+
+
+// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
+
+
+bool Foam::combustionModel::read()
+{
+    if (regIOobject::read())
+    {
+        this->lookup("active") >> active_;
+        coeffs_ = subDict(modelType_ + "Coeffs");
+        return true;
+    }
+    else
+    {
+        return false;
+    }
+}
+
+
+Foam::tmp<Foam::volScalarField> Foam::combustionModel::Sh() const
+{
+    return tmp<Foam::volScalarField>
+    (
+        new volScalarField
+        (
+            IOobject
+            (
+                IOobject::groupName("Sh", phaseName_),
+                mesh_.time().timeName(),
+                mesh_,
+                IOobject::NO_READ,
+                IOobject::NO_WRITE
+            ),
+            mesh_,
+            dimensionedScalar("zero", dimEnergy/dimVolume/dimTime, 0.0)
+        )
+    );
+}
+
+
+// ************************************************************************* //

libs/combustionModels_POSTECH/combustionModel/combustionModel.H

@@ -0,0 +1,170 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-2015 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 <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::combustionModel
+
+Description
+    Base class for combustion models
+
+SourceFiles
+    combustionModel.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef combustionModel_H
+#define combustionModel_H
+
+#include "IOdictionary.H"
+#include "turbulentFluidThermoModel.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                     Class combustionModel Declaration
+\*---------------------------------------------------------------------------*/
+
+class combustionModel
+:
+    public IOdictionary
+{
+    // Private Member Functions
+
+        //- Disallow copy construct
+        combustionModel(const combustionModel&);
+
+        //- Disallow default bitwise assignment
+        void operator=(const combustionModel&);
+
+
+protected:
+
+    // Protected data
+
+        //- Reference to the turbulence model
+        compressibleTurbulenceModel* turbulencePtr_;
+
+        //- Reference to the mesh database
+        const fvMesh& mesh_;
+
+        //- Active
+        Switch active_;
+
+        //- Dictionary of the model
+        dictionary coeffs_;
+
+        //- Model type
+        const word modelType_;
+
+        //- Phase name
+        const word phaseName_;
+
+
+public:
+
+    //- Runtime type information
+    TypeName("combustionModel");
+
+
+    // Constructors
+
+        //- Construct from components
+        combustionModel
+        (
+            const word& modelType,
+            const fvMesh& mesh,
+            const word& phaseName=word::null
+        );
+
+
+    //- Destructor
+    virtual ~combustionModel();
+
+
+    // Member Functions
+
+        // Access
+
+            //- Return const access to the mesh database
+            inline const fvMesh& mesh() const;
+
+            //- Return const access to phi
+            inline const surfaceScalarField& phi() const;
+
+            //- Return const access to rho
+            virtual tmp<volScalarField> rho() const = 0;
+
+            //- Return access to turbulence
+            inline const compressibleTurbulenceModel& turbulence() const;
+
+            //- Set turbulence
+            inline void setTurbulence
+            (
+                compressibleTurbulenceModel& turbModel
+            );
+
+            //- Is combustion active?
+            inline const Switch& active() const;
+
+            //- Return const dictionary of the model
+            inline const dictionary& coeffs() const;
+
+
+    // Evolution
+
+        //- Correct combustion rate
+        virtual void correct() = 0;
+
+        //- Fuel consumption rate matrix, i.e. source term for fuel equation
+        virtual tmp<fvScalarMatrix> R(volScalarField& Y) const = 0;
+
+        //- Heat release rate calculated from fuel consumption rate matrix
+        virtual tmp<volScalarField> dQ() const = 0;
+
+        //-  Return source for enthalpy equation [kg/m/s3]
+        virtual tmp<volScalarField> Sh() const;
+
+
+    // IO
+
+        //- Update properties from given dictionary
+        virtual bool read();
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#include "combustionModelI.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

libs/combustionModels_POSTECH/combustionModel/combustionModelI.H

@@ -0,0 +1,93 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-2015 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 <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+inline const Foam::fvMesh& Foam::combustionModel::mesh() const
+{
+    return mesh_;
+}
+
+
+inline const Foam::surfaceScalarField& Foam::combustionModel::phi() const
+{
+    if (turbulencePtr_)
+    {
+        return turbulencePtr_->phi();
+    }
+    else
+    {
+        FatalErrorInFunction
+            << "turbulencePtr_ is empty. Please use "
+            << "combustionModel::setTurbulence "
+            << "(compressibleTurbulenceModel& )"
+            << abort(FatalError);
+
+        return turbulencePtr_->phi();
+    }
+}
+
+
+inline const Foam::compressibleTurbulenceModel&
+Foam::combustionModel::turbulence() const
+{
+    if (turbulencePtr_)
+    {
+        return *turbulencePtr_;
+    }
+    else
+    {
+        FatalErrorInFunction
+            << "turbulencePtr_ is empty. Please use "
+            << "combustionModel::setTurbulence "
+            << "(compressibleTurbulenceModel& )"
+            << abort(FatalError);
+
+        return *turbulencePtr_;
+    }
+}
+
+
+inline const Foam::Switch& Foam::combustionModel::active() const
+{
+    return active_;
+}
+
+
+inline void Foam::combustionModel::setTurbulence
+(
+    compressibleTurbulenceModel& turbModel
+)
+{
+    turbulencePtr_ = &turbModel;
+}
+
+
+inline const Foam::dictionary& Foam::combustionModel::coeffs() const
+{
+    return coeffs_;
+}
+
+// ************************************************************************* //

libs/combustionModels_POSTECH/combustionModel/makeCombustionTypes.H

@@ -0,0 +1,93 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-2015 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 <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef makeCombustionTypes_H
+#define makeCombustionTypes_H
+
+#include "addToRunTimeSelectionTable.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#define makeCombustionTypesThermo(CombModel, CombType, Thermo, Table)          \
+                                                                               \
+    typedef Foam::combustionModels::CombModel                                  \
+        <Foam::combustionModels::CombType, Foam::Thermo>                       \
+        CombModel##CombType##Thermo;                                           \
+                                                                               \
+    defineTemplateTypeNameAndDebugWithName                                     \
+    (                                                                          \
+        CombModel##CombType##Thermo,                                           \
+        #CombModel"<"#CombType","#Thermo">",                                   \
+        0                                                                      \
+    );                                                                         \
+                                                                               \
+    namespace Foam                                                             \
+    {                                                                          \
+        namespace combustionModels                                             \
+        {                                                                      \
+            typedef CombModel<CombType, Thermo> CombModel##CombType##Thermo;   \
+            addToRunTimeSelectionTable                                         \
+            (                                                                  \
+                Table,                                                         \
+                CombModel##CombType##Thermo,                                   \
+                dictionary                                                     \
+            );                                                                 \
+        }                                                                      \
+    }
+
+
+#define makeCombustionTypes(CombModel, CombType, Table)                        \
+                                                                               \
+    typedef Foam::combustionModels::CombModel                                  \
+        <Foam::combustionModels::CombType>                                     \
+        CombModel##CombType;                                                   \
+                                                                               \
+    defineTemplateTypeNameAndDebugWithName                                     \
+    (                                                                          \
+        CombModel##CombType,                                                   \
+        #CombModel"<"#CombType">",                                             \
+        0                                                                      \
+    );                                                                         \
+                                                                               \
+    namespace Foam                                                             \
+    {                                                                          \
+        namespace combustionModels                                             \
+        {                                                                      \
+            typedef CombModel<CombType> CombModel##CombType;                   \
+                                                                               \
+            addToRunTimeSelectionTable                                         \
+            (                                                                  \
+                Table,                                                         \
+                CombModel##CombType,                                           \
+                dictionary                                                     \
+            );                                                                 \
+        }                                                                      \
+    }
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

libs/combustionModels_POSTECH/diffusion/diffusion.C

@@ -0,0 +1,114 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2012-2015 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 <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "diffusion.H"
+#include "fvcGrad.H"
+
+namespace Foam
+{
+namespace combustionModels
+{
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+template<class CombThermoType, class ThermoType>
+diffusion<CombThermoType, ThermoType>::diffusion
+(
+    const word& modelType,
+    const fvMesh& mesh,
+    const word& phaseName
+)
+:
+    singleStepCombustion<CombThermoType, ThermoType>
+    (
+        modelType,
+        mesh,
+        phaseName
+    ),
+    C_(readScalar(this->coeffs().lookup("C"))),
+    oxidantName_(this->coeffs().template lookupOrDefault<word>("oxidant", "O2"))
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+template<class CombThermoType, class ThermoType>
+diffusion<CombThermoType, ThermoType>::~diffusion()
+{}
+
+
+// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
+
+template<class CombThermoType, class ThermoType>
+void diffusion<CombThermoType, ThermoType>::correct()
+{
+    this->wFuel_ ==
+        dimensionedScalar("zero", dimMass/pow3(dimLength)/dimTime, 0.0);
+
+    if (this->active())
+    {
+        this->singleMixturePtr_->fresCorrect();
+
+        const label fuelI = this->singleMixturePtr_->fuelIndex();
+
+        const volScalarField& YFuel =
+            this->thermoPtr_->composition().Y()[fuelI];
+
+        if (this->thermoPtr_->composition().contains(oxidantName_))
+        {
+            const volScalarField& YO2 =
+                this->thermoPtr_->composition().Y(oxidantName_);
+
+            this->wFuel_ ==
+                C_*this->turbulence().muEff()
+               *mag(fvc::grad(YFuel) & fvc::grad(YO2))
+               *pos(YFuel)*pos(YO2);
+        }
+    }
+}
+
+
+template<class CombThermoType, class ThermoType>
+bool diffusion<CombThermoType, ThermoType>::read()
+{
+    if (singleStepCombustion<CombThermoType, ThermoType>::read())
+    {
+        this->coeffs().lookup("C") >> C_ ;
+        this->coeffs().readIfPresent("oxidant", oxidantName_);
+        return true;
+    }
+    else
+    {
+        return false;
+    }
+}
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace combustionModels
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

libs/combustionModels_POSTECH/diffusion/diffusion.H

@@ -0,0 +1,128 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2012-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 <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::combustionModels::diffusion
+
+Description
+    Simple diffusion-based combustion model based on the principle mixed is
+    burnt. Additional parameter C is used to distribute the heat release rate
+    in time.
+
+SourceFiles
+    diffusion.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef diffusion_H
+#define diffusion_H
+
+#include "singleStepCombustion.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace combustionModels
+{
+
+/*---------------------------------------------------------------------------*\
+                         Class diffusion Declaration
+\*---------------------------------------------------------------------------*/
+
+template<class CombThermoType, class ThermoType>
+class diffusion
+:
+    public singleStepCombustion<CombThermoType, ThermoType>
+{
+    // Private data
+
+        //- Model constant
+        scalar C_;
+
+        //- Name of oxidant - default is "O2"
+        word oxidantName_;
+
+
+    // Private Member Functions
+
+        //- Disallow copy construct
+        diffusion(const diffusion&);
+
+        //- Disallow default bitwise assignment
+        void operator=(const diffusion&);
+
+
+public:
+
+    //- Runtime type information
+    TypeName("diffusion");
+
+
+    // Constructors
+
+        //- Construct from components
+        diffusion
+        (
+            const word& modelType,
+            const fvMesh& mesh,
+            const word& phaseName
+        );
+
+
+    //- Destructor
+    virtual ~diffusion();
+
+
+    // Member Functions
+
+        // Evolution
+
+            //- Correct combustion rate
+            virtual void correct();
+
+
+        // IO
+
+            //- Update properties
+            virtual bool read();
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace combustionModels
+} // End namespace Foam
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+    #include "diffusion.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

libs/combustionModels_POSTECH/diffusion/diffusions.C

@@ -0,0 +1,103 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2012-2015 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 <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "makeCombustionTypes.H"
+
+#include "thermoPhysicsTypes.H"
+#include "psiThermoCombustion.H"
+#include "rhoThermoCombustion.H"
+#include "diffusion.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+// Combustion models based on sensibleEnthalpy
+makeCombustionTypesThermo
+(
+    diffusion,
+    psiThermoCombustion,
+    gasHThermoPhysics,
+    psiCombustionModel
+);
+
+makeCombustionTypesThermo
+(
+    diffusion,
+    psiThermoCombustion,
+    constGasHThermoPhysics,
+    psiCombustionModel
+);
+
+makeCombustionTypesThermo
+(
+    diffusion,
+    rhoThermoCombustion,
+    gasHThermoPhysics,
+    rhoCombustionModel
+);
+
+makeCombustionTypesThermo
+(
+    diffusion,
+    rhoThermoCombustion,
+    constGasHThermoPhysics,
+    rhoCombustionModel
+);
+
+ // Combustion models based on sensibleInternalEnergy
+
+makeCombustionTypesThermo
+(
+    diffusion,
+    psiThermoCombustion,
+    gasEThermoPhysics,
+    psiCombustionModel
+);
+
+makeCombustionTypesThermo
+(
+    diffusion,
+    psiThermoCombustion,
+    constGasEThermoPhysics,
+    psiCombustionModel
+);
+
+makeCombustionTypesThermo
+(
+    diffusion,
+    rhoThermoCombustion,
+    gasEThermoPhysics,
+    rhoCombustionModel
+);
+
+makeCombustionTypesThermo
+(
+    diffusion,
+    rhoThermoCombustion,
+    constGasEThermoPhysics,
+    rhoCombustionModel
+);
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

libs/combustionModels_POSTECH/infinitelyFastChemistry/infinitelyFastChemistry.C

@@ -0,0 +1,111 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-2015 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 <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "infinitelyFastChemistry.H"
+
+namespace Foam
+{
+namespace combustionModels
+{
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+template<class CombThermoType, class ThermoType>
+infinitelyFastChemistry<CombThermoType, ThermoType>::infinitelyFastChemistry
+(
+    const word& modelType,
+    const fvMesh& mesh,
+    const word& phaseName
+)
+:
+    singleStepCombustion<CombThermoType, ThermoType>
+    (
+        modelType,
+        mesh,
+        phaseName
+    ),
+    C_(readScalar(this->coeffs().lookup("C")))
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+template<class CombThermoType, class ThermoType>
+infinitelyFastChemistry<CombThermoType, ThermoType>::~infinitelyFastChemistry()
+{}
+
+
+// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
+
+template<class CombThermoType, class ThermoType>
+void infinitelyFastChemistry<CombThermoType, ThermoType>::correct()
+{
+    this->wFuel_ ==
+        dimensionedScalar("zero", dimMass/pow3(dimLength)/dimTime, 0.0);
+
+    if (this->active())
+    {
+        this->singleMixturePtr_->fresCorrect();
+
+        const label fuelI = this->singleMixturePtr_->fuelIndex();
+
+        const volScalarField& YFuel =
+            this->thermoPtr_->composition().Y()[fuelI];
+
+        const dimensionedScalar s = this->singleMixturePtr_->s();
+
+        if (this->thermoPtr_->composition().contains("O2"))
+        {
+            const volScalarField& YO2 = this->thermoPtr_->composition().Y("O2");
+
+            this->wFuel_ ==
+                this->rho()/(this->mesh().time().deltaT()*C_)
+               *min(YFuel, YO2/s.value());
+        }
+    }
+}
+
+
+template<class CombThermoType, class ThermoType>
+bool infinitelyFastChemistry<CombThermoType, ThermoType>::read()
+{
+    if (singleStepCombustion<CombThermoType, ThermoType>::read())
+    {
+        this->coeffs().lookup("C") >> C_ ;
+        return true;
+    }
+    else
+    {
+        return false;
+    }
+}
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace combustionModels
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

libs/combustionModels_POSTECH/infinitelyFastChemistry/infinitelyFastChemistry.H

@@ -0,0 +1,125 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-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 <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::combustionModels::infinitelyFastChemistry
+
+Description
+    Simple infinitely fast chemistry combustion model based on the principle
+    mixed is burnt. Additional parameter C is used to distribute the heat
+    release rate.in time
+
+SourceFiles
+    infinitelyFastChemistry.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef infinitelyFastChemistry_H
+#define infinitelyFastChemistry_H
+
+#include "singleStepCombustion.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace combustionModels
+{
+
+/*---------------------------------------------------------------------------*\
+                  Class infinitelyFastChemistry Declaration
+\*---------------------------------------------------------------------------*/
+
+template<class CombThermoType, class ThermoType>
+class infinitelyFastChemistry
+:
+    public singleStepCombustion<CombThermoType, ThermoType>
+{
+    // Private data
+
+        //- Model constant
+        scalar C_;
+
+
+    // Private Member Functions
+
+        //- Disallow copy construct
+        infinitelyFastChemistry(const infinitelyFastChemistry&);
+
+        //- Disallow default bitwise assignment
+        void operator=(const infinitelyFastChemistry&);
+
+
+public:
+
+    //- Runtime type information
+    TypeName("infinitelyFastChemistry");
+
+
+    // Constructors
+
+        //- Construct from components
+        infinitelyFastChemistry
+        (
+            const word& modelType,
+            const fvMesh& mesh,
+            const word& phaseName
+        );
+
+
+    //- Destructor
+    virtual ~infinitelyFastChemistry();
+
+
+    // Member Functions
+
+        // Evolution
+
+            //- Correct combustion rate
+            virtual void correct();
+
+
+        // IO
+
+            //- Update properties
+            virtual bool read();
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace combustionModels
+} // End namespace Foam
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+    #include "infinitelyFastChemistry.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

libs/combustionModels_POSTECH/infinitelyFastChemistry/infinitelyFastChemistrys.C

@@ -0,0 +1,104 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-2015 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 <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "makeCombustionTypes.H"
+
+#include "thermoPhysicsTypes.H"
+#include "psiThermoCombustion.H"
+#include "rhoThermoCombustion.H"
+#include "infinitelyFastChemistry.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+// Combustion models based on sensibleEnthalpy
+
+makeCombustionTypesThermo
+(
+    infinitelyFastChemistry,
+    psiThermoCombustion,
+    gasHThermoPhysics,
+    psiCombustionModel
+);
+
+makeCombustionTypesThermo
+(
+    infinitelyFastChemistry,
+    psiThermoCombustion,
+    constGasHThermoPhysics,
+    psiCombustionModel
+);
+
+makeCombustionTypesThermo
+(
+    infinitelyFastChemistry,
+    rhoThermoCombustion,
+    gasHThermoPhysics,
+    rhoCombustionModel
+);
+
+makeCombustionTypesThermo
+(
+    infinitelyFastChemistry,
+    rhoThermoCombustion,
+    constGasHThermoPhysics,
+    rhoCombustionModel
+);
+
+// Combustion models based on sensibleInternalEnergy
+
+makeCombustionTypesThermo
+(
+    infinitelyFastChemistry,
+    psiThermoCombustion,
+    gasEThermoPhysics,
+    psiCombustionModel
+);
+
+makeCombustionTypesThermo
+(
+    infinitelyFastChemistry,
+    psiThermoCombustion,
+    constGasEThermoPhysics,
+    psiCombustionModel
+);
+
+makeCombustionTypesThermo
+(
+    infinitelyFastChemistry,
+    rhoThermoCombustion,
+    gasEThermoPhysics,
+    rhoCombustionModel
+);
+
+makeCombustionTypesThermo
+(
+    infinitelyFastChemistry,
+    rhoThermoCombustion,
+    constGasEThermoPhysics,
+    rhoCombustionModel
+);
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

libs/combustionModels_POSTECH/laminar/laminar.C

@@ -0,0 +1,214 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013-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 <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "laminar.H"
+#include "fvmSup.H"
+#include "localEulerDdtScheme.H"
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+template<class Type>
+Foam::combustionModels::laminar<Type>::laminar
+(
+    const word& modelType,
+    const fvMesh& mesh,
+    const word& phaseName
+)
+:
+    Type(modelType, mesh, phaseName),
+    integrateReactionRate_
+    (
+        this->coeffs().lookupOrDefault("integrateReactionRate", true)
+    )
+{
+    if (integrateReactionRate_)
+    {
+        Info<< "    using integrated reaction rate" << endl;
+    }
+    else
+    {
+        Info<< "    using instantaneous reaction rate" << endl;
+    }
+}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+template<class Type>
+Foam::combustionModels::laminar<Type>::~laminar()
+{}
+
+
+// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
+
+template<class Type>
+Foam::tmp<Foam::volScalarField>
+Foam::combustionModels::laminar<Type>::tc() const
+{
+    return this->chemistryPtr_->tc();
+}
+
+
+template<class Type>
+void Foam::combustionModels::laminar<Type>::correct()
+{
+    if (this->active())
+    {
+        if (integrateReactionRate_)
+        {
+            if (fv::localEulerDdt::enabled(this->mesh()))
+            {
+                const scalarField& rDeltaT =
+                    fv::localEulerDdt::localRDeltaT(this->mesh());
+
+                if (this->coeffs().found("maxIntegrationTime"))
+                {
+                    scalar maxIntegrationTime
+                    (
+                        readScalar(this->coeffs().lookup("maxIntegrationTime"))
+                    );
+
+                    this->chemistryPtr_->solve
+                    (
+                        min(1.0/rDeltaT, maxIntegrationTime)()
+                    );
+                }
+                else
+                {
+                    this->chemistryPtr_->solve((1.0/rDeltaT)());
+                }
+            }
+            else
+            {
+                this->chemistryPtr_->solve(this->mesh().time().deltaTValue());
+            }
+        }
+        else
+        {
+            this->chemistryPtr_->calculate();
+        }
+    }
+}
+
+
+template<class Type>
+Foam::tmp<Foam::fvScalarMatrix>
+Foam::combustionModels::laminar<Type>::R(volScalarField& Y) const
+{
+    tmp<fvScalarMatrix> tSu(new fvScalarMatrix(Y, dimMass/dimTime));
+
+    fvScalarMatrix& Su = tSu.ref();
+
+    if (this->active())
+    {
+        const label specieI =
+            this->thermo().composition().species()[Y.member()];
+
+        Su += this->chemistryPtr_->RR(specieI);
+    }
+
+    return tSu;
+}
+
+
+template<class Type>
+Foam::tmp<Foam::volScalarField>
+Foam::combustionModels::laminar<Type>::dQ() const
+{
+    tmp<volScalarField> tdQ
+    (
+        new volScalarField
+        (
+            IOobject
+            (
+                IOobject::groupName(typeName + ":dQ", this->phaseName_),
+                this->mesh().time().timeName(),
+                this->mesh(),
+                IOobject::NO_READ,
+                IOobject::NO_WRITE,
+                false
+            ),
+            this->mesh(),
+            dimensionedScalar("dQ", dimEnergy/dimTime, 0.0)
+        )
+    );
+
+    if (this->active())
+    {
+        tdQ.ref() = this->chemistryPtr_->dQ();
+    }
+
+    return tdQ;
+}
+
+
+template<class Type>
+Foam::tmp<Foam::volScalarField>
+Foam::combustionModels::laminar<Type>::Sh() const
+{
+    tmp<volScalarField> tSh
+    (
+        new volScalarField
+        (
+            IOobject
+            (
+                IOobject::groupName(typeName + ":Sh", this->phaseName_),
+                this->mesh().time().timeName(),
+                this->mesh(),
+                IOobject::NO_READ,
+                IOobject::NO_WRITE,
+                false
+            ),
+            this->mesh(),
+            dimensionedScalar("zero", dimEnergy/dimTime/dimVolume, 0.0)
+        )
+    );
+
+    if (this->active())
+    {
+        tSh.ref() = this->chemistryPtr_->Sh();
+    }
+
+    return tSh;
+}
+
+
+template<class Type>
+bool Foam::combustionModels::laminar<Type>::read()
+{
+    if (Type::read())
+    {
+        this->coeffs().lookup("integrateReactionRate")
+            >> integrateReactionRate_;
+        return true;
+    }
+    else
+    {
+        return false;
+    }
+}
+
+
+// ************************************************************************* //

libs/combustionModels_POSTECH/laminar/laminar.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2013-2016 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -22,126 +22,114 @@ License
     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
 
 Class
-    Foam::LinearInterpolator1D
+    Foam::combustionModels::laminar
 
 Description
+    Laminar combustion model.
 
 SourceFiles
-    LinearInterpolator1DI.H
-    LinearInterpolator1D.C
-    LinearInterpolator1DIO.C
+    laminar.C
 
 \*---------------------------------------------------------------------------*/
 
-#ifndef LinearInterpolator1D_H
-#define LinearInterpolator1D_H
-
-#include "scalarField.H"
+#ifndef laminar_H
+#define laminar_H
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 namespace Foam
 {
-
-// Forward declaration of classes
-class Istream;
-class Ostream;
-
-// Forward declaration of friend functions and operators
-class LinearInterpolator1D;
-Istream& operator>>(Istream&, LinearInterpolator1D&);
-Ostream& operator<<(Ostream&, const LinearInterpolator1D&);
-
+namespace combustionModels
+{
 
 /*---------------------------------------------------------------------------*\
-                         Class LinearInterpolator1D Declaration
+                            Class laminar Declaration
 \*---------------------------------------------------------------------------*/
 
-class LinearInterpolator1D
+template<class Type>
+class laminar
+:
+    public Type
 {
     // Private data
 
-        //- Source grids
-        scalarField srcGrid_;
+        //- Integrate reaction rate over the time-step
+        //  using the selected ODE solver
+        bool integrateReactionRate_;
 
-        //- Destination grids
-        scalarField dstGrid_;
+protected:
 
-        //- Source grid to refer for each dst grid
-        labelList idxHigh_;
+    // Protected Member Functions
 
-        //- Coefficients for dst grid points
-        scalarField coef_;
+        //- Return the chemical time scale
+        tmp<volScalarField> tc() const;
 
+private:
 
     // Private Member Functions
 
-        //- Check grid sorted
-        bool _checkSorted(const scalarField& grid);
-
-        //- Disallow default bitwise copy construct
-        LinearInterpolator1D(const LinearInterpolator1D&);
+        //- Disallow copy construct
+        laminar(const laminar&);
 
         //- Disallow default bitwise assignment
-        void operator=(const LinearInterpolator1D&);
+        void operator=(const laminar&);
 
 
 public:
 
+    //- Runtime type information
+    TypeName("laminar");
+
 
     // Constructors
 
-        //- Construct null
-        LinearInterpolator1D();
-
         //- Construct from components
-        LinearInterpolator1D(const scalarField& srcGrid, const scalarField& dstGrid);
-
-        //- Construct from Istream
-        LinearInterpolator1D(Istream&);
-
+        laminar
+        (
+            const word& modelType,
+            const fvMesh& mesh,
+            const word& phaseName
+        );
 
 
     //- Destructor
-    ~LinearInterpolator1D();
+    virtual ~laminar();
 
 
     // Member Functions
 
-        // Access
+        // Evolution
 
-        // Check
+            //- Correct combustion rate
+            virtual void correct();
 
-        // Edit
+            //- Fuel consumption rate matrix.
+            virtual tmp<fvScalarMatrix> R(volScalarField& Y) const;
 
-        // Write
+            //- Heat release rate calculated from fuel consumption rate matrix
+            virtual tmp<volScalarField> dQ() const;
 
-        //- Evaluate interpolation
-        void interpolate(const scalarField& src, scalarField& dst);
+            //-  Return source for enthalpy equation [kg/m/s3]
+            virtual tmp<volScalarField> Sh() const;
 
 
-    // Member Operators
+    // IO
 
-
-
-    // Friend Functions
-
-    // Friend Operators
-
-    // IOstream Operators
-
-        friend Istream& operator>>(Istream&, LinearInterpolator1D&);
-        friend Ostream& operator<<(Ostream&, const LinearInterpolator1D&);
+            //- Update properties from given dictionary
+            virtual bool read();
 };
 
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
+} // End namespace combustionModels
 } // End namespace Foam
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
-#include "LinearInterpolator1DI.H"
+#ifdef NoRepository
+    #include "laminar.C"
+#endif
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 

libs/combustionModels_POSTECH/laminar/laminars.C

@@ -0,0 +1,38 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013-2015 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 <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "makeCombustionTypes.H"
+
+#include "psiChemistryCombustion.H"
+#include "rhoChemistryCombustion.H"
+#include "laminar.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+makeCombustionTypes(laminar, psiChemistryCombustion, psiCombustionModel);
+makeCombustionTypes(laminar, rhoChemistryCombustion, rhoCombustionModel);
+
+
+// ************************************************************************* //

libs/combustionModels_POSTECH/noCombustion/noCombustion.C

@@ -0,0 +1,141 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-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 <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "noCombustion.H"
+#include "fvmSup.H"
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+template<class CombThermoType>
+Foam::combustionModels::noCombustion<CombThermoType>::noCombustion
+(
+    const word& modelType,
+    const fvMesh& mesh,
+    const word& phaseName
+)
+:
+    CombThermoType(modelType, mesh, phaseName)
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+template<class CombThermoType>
+Foam::combustionModels::noCombustion<CombThermoType>::~noCombustion()
+{}
+
+
+// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
+
+template<class CombThermoType>
+void Foam::combustionModels::noCombustion<CombThermoType>::correct()
+{
+//  Do Nothing
+}
+
+
+template<class CombThermoType>
+Foam::tmp<Foam::fvScalarMatrix>
+Foam::combustionModels::noCombustion<CombThermoType>::R
+(
+    volScalarField& Y
+) const
+{
+    tmp<fvScalarMatrix> tSu
+    (
+        new fvScalarMatrix(Y, dimMass/dimTime)
+    );
+
+    return tSu;
+}
+
+
+template<class CombThermoType>
+Foam::tmp<Foam::volScalarField>
+Foam::combustionModels::noCombustion<CombThermoType>::dQ() const
+{
+    tmp<volScalarField> tdQ
+    (
+        new volScalarField
+        (
+            IOobject
+            (
+                IOobject::groupName("dQ", this->phaseName_),
+                this->mesh().time().timeName(),
+                this->mesh(),
+                IOobject::NO_READ,
+                IOobject::NO_WRITE,
+                false
+            ),
+            this->mesh(),
+            dimensionedScalar("dQ", dimEnergy/dimTime, 0.0)
+        )
+    );
+
+    return tdQ;
+}
+
+
+template<class CombThermoType>
+Foam::tmp<Foam::volScalarField>
+Foam::combustionModels::noCombustion<CombThermoType>::Sh() const
+{
+    tmp<volScalarField> tSh
+    (
+        new volScalarField
+        (
+            IOobject
+            (
+                IOobject::groupName("Sh", this->phaseName_),
+                this->mesh().time().timeName(),
+                this->mesh(),
+                IOobject::NO_READ,
+                IOobject::NO_WRITE,
+                false
+            ),
+            this->mesh(),
+            dimensionedScalar("zero", dimEnergy/dimTime/dimVolume, 0.0)
+        )
+    );
+
+    return tSh;
+}
+
+
+template<class CombThermoType>
+bool Foam::combustionModels::noCombustion<CombThermoType>::read()
+{
+    if (CombThermoType::read())
+    {
+        return true;
+    }
+    else
+    {
+        return false;
+    }
+}
+
+
+// ************************************************************************* //

libs/combustionModels_POSTECH/noCombustion/noCombustion.H

@@ -0,0 +1,121 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-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 <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::combustionModels::noCombustion
+
+Description
+    Dummy combustion model for 'no combustion'
+
+SourceFiles
+    noCombustion.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef noCombustion_H
+#define noCombustion_H
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace combustionModels
+{
+
+/*---------------------------------------------------------------------------*\
+                        Class noCombustion Declaration
+\*---------------------------------------------------------------------------*/
+
+template<class CombThermoType>
+class noCombustion
+:
+    public CombThermoType
+{
+
+    //- Disallow copy construct
+    noCombustion(const noCombustion&);
+
+    //- Disallow default bitwise assignment
+    void operator=(const noCombustion&);
+
+
+public:
+
+    //- Runtime type information
+    TypeName("noCombustion");
+
+
+    // Constructors
+
+        //- Construct from components
+        noCombustion
+        (
+            const word& modelType,
+            const fvMesh& mesh,
+            const word& phaseName
+        );
+
+
+    //- Destructor
+    virtual ~noCombustion();
+
+
+    // Member Functions
+
+        // Evolution
+
+            //- Correct combustion rate
+            virtual void correct();
+
+            //- Fuel consumption rate matrix
+            virtual tmp<fvScalarMatrix> R(volScalarField& Y) const;
+
+            //- Heat release rate calculated from fuel consumption rate matrix
+            virtual tmp<volScalarField> dQ() const;
+
+            //-  Return source for enthalpy equation [kg/m/s3]
+            virtual tmp<volScalarField> Sh() const;
+
+    // IO
+
+            //- Update properties from given dictionary
+            virtual bool read();
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace combustionModels
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+    #include "noCombustion.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

libs/combustionModels_POSTECH/noCombustion/noCombustions.C

@@ -0,0 +1,52 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2012-2015 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 <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "makeCombustionTypes.H"
+
+#include "psiCombustionModel.H"
+#include "rhoCombustionModel.H"
+#include "psiThermoCombustion.H"
+#include "rhoThermoCombustion.H"
+
+#include "noCombustion.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+makeCombustionTypes
+(
+    noCombustion,
+    psiThermoCombustion,
+    psiCombustionModel
+);
+
+makeCombustionTypes
+(
+    noCombustion,
+    rhoThermoCombustion,
+    rhoCombustionModel
+);
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

libs/combustionModels_POSTECH/psiCombustionModel/psiChemistryCombustion/psiChemistryCombustion.C

@@ -0,0 +1,77 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2012-2015 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 <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "psiChemistryCombustion.H"
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::combustionModels::psiChemistryCombustion::psiChemistryCombustion
+(
+    const word& modelType,
+    const fvMesh& mesh,
+    const word& phaseName
+)
+:
+    psiCombustionModel(modelType, mesh, phaseName),
+    chemistryPtr_(psiChemistryModel::New(mesh, phaseName))
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+Foam::combustionModels::psiChemistryCombustion::~psiChemistryCombustion()
+{}
+
+
+// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
+
+Foam::autoPtr<Foam::psiChemistryModel>
+Foam::combustionModels::psiChemistryCombustion::chem()
+{
+    return chemistryPtr_;
+}
+
+Foam::psiReactionThermo&
+Foam::combustionModels::psiChemistryCombustion::thermo()
+{
+    return chemistryPtr_->thermo();
+}
+
+
+const Foam::psiReactionThermo&
+Foam::combustionModels::psiChemistryCombustion::thermo() const
+{
+    return chemistryPtr_->thermo();
+}
+
+
+Foam::tmp<Foam::volScalarField>
+Foam::combustionModels::psiChemistryCombustion::rho() const
+{
+    return chemistryPtr_->thermo().rho();
+}
+
+
+// ************************************************************************* //

libs/combustionModels_POSTECH/psiCombustionModel/psiChemistryCombustion/psiChemistryCombustion.H

@@ -0,0 +1,113 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2012-2015 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 <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::psiChemistryCombustion
+
+Description
+    Compressibility-based chemistry model wrapper for combustion models
+
+SourceFiles
+    psiChemistryCombustion.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef psiChemistryCombustion_H
+#define psiChemistryCombustion_H
+
+#include "autoPtr.H"
+#include "psiCombustionModel.H"
+#include "psiChemistryModel.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace combustionModels
+{
+
+/*---------------------------------------------------------------------------*\
+                    class psiChemistryCombustion Declaration
+\*---------------------------------------------------------------------------*/
+
+class psiChemistryCombustion
+:
+    public psiCombustionModel
+{
+    // Private Member Functions
+
+        //- Construct as copy (not implemented)
+        psiChemistryCombustion(const psiChemistryCombustion&);
+
+        //- Disallow default bitwise assignment
+        void operator=(const psiChemistryCombustion&);
+
+
+protected:
+
+    // Protected data
+
+        //- Pointer to chemistry model
+        autoPtr<psiChemistryModel> chemistryPtr_;
+
+
+public:
+
+    // Constructors
+
+        //- Construct from components and thermo
+        psiChemistryCombustion
+        (
+            const word& modelType,
+            const fvMesh& mesh,
+            const word& phaseName
+        );
+
+
+    //- Destructor
+    virtual ~psiChemistryCombustion();
+
+
+    // Member Functions
+        autoPtr<psiChemistryModel> chem(); 
+        //- Return access to the thermo package
+        virtual psiReactionThermo& thermo();
+
+        //- Return const access to the thermo package
+        virtual const psiReactionThermo& thermo() const;
+
+        //- Return const access to the density field
+        virtual tmp<volScalarField> rho() const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace combustionModels
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

libs/combustionModels_POSTECH/psiCombustionModel/psiCombustionModel/psiCombustionModel.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2015 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -23,36 +23,51 @@ License
 
 \*---------------------------------------------------------------------------*/
 
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+#include "psiCombustionModel.H"
 
-// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+/* * * * * * * * * * * * * * * private static data * * * * * * * * * * * * * */
 
+namespace Foam
+{
+namespace combustionModels
+{
+    defineTypeNameAndDebug(psiCombustionModel, 0);
+    defineRunTimeSelectionTable(psiCombustionModel, dictionary);
+}
+}
 
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
+Foam::combustionModels::psiCombustionModel::psiCombustionModel
+(
+    const word& modelType,
+    const fvMesh& mesh,
+    const word& phaseName
+)
+:
+    combustionModel(modelType, mesh, phaseName)
+{}
+
 
 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
 
-
-// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+Foam::combustionModels::psiCombustionModel::~psiCombustionModel()
+{}
 
 
-// * * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * //
+// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
 
-
-// * * * * * * * * * * * * * * * Friend Functions  * * * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * * * Friend Operators  * * * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * * IOstream Operators  * * * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * * * Ostream Operator  * * * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+bool Foam::combustionModels::psiCombustionModel::read()
+{
+    if (combustionModel::read())
+    {
+        return true;
+    }
+    else
+    {
+        return false;
+    }
+}
 
 
 // ************************************************************************* //

libs/combustionModels_POSTECH/psiCombustionModel/psiCombustionModel/psiCombustionModel.H

@@ -0,0 +1,150 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-2015 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 <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::psiCombustionModel
+
+Description
+    Combustion models for compressibility-based thermodynamics
+
+SourceFiles
+    psiCombustionModelI.H
+    psiCombustionModel.C
+    psiCombustionModelNew.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef psiCombustionModel_H
+#define psiCombustionModel_H
+
+#include "combustionModel.H"
+#include "autoPtr.H"
+#include "runTimeSelectionTables.H"
+#include "psiReactionThermo.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace combustionModels
+{
+
+/*---------------------------------------------------------------------------*\
+                     class psiCombustionModel Declaration
+\*---------------------------------------------------------------------------*/
+
+class psiCombustionModel
+:
+    public combustionModel
+{
+
+    // Private Member Functions
+
+        //- Construct as copy (not implemented)
+        psiCombustionModel(const psiCombustionModel&);
+
+        //- Disallow default bitwise assignment
+        void operator=(const psiCombustionModel&);
+
+
+public:
+
+    //- Runtime type information
+    TypeName("psiCombustionModel");
+
+
+    //- Declare run-time constructor selection tables
+    declareRunTimeSelectionTable
+    (
+        autoPtr,
+        psiCombustionModel,
+        dictionary,
+        (
+            const word& modelType,
+            const fvMesh& mesh,
+            const word& phaseName
+        ),
+        (modelType, mesh, phaseName)
+    );
+
+
+    // Constructors
+
+        //- Construct from components
+        psiCombustionModel
+        (
+            const word& modelType,
+            const fvMesh& mesh,
+            const word& phaseName
+        );
+
+
+
+    //- Selector
+    static autoPtr<psiCombustionModel> New
+    (
+        const fvMesh& mesh,
+        const word& phaseName=word::null
+    );
+
+
+    //- Destructor
+    virtual ~psiCombustionModel();
+
+
+    // Member Functions
+
+        //- Return access to the thermo package
+        virtual psiReactionThermo& thermo() = 0;
+
+        //- Return const access to the thermo package
+        virtual const psiReactionThermo& thermo() const = 0;
+
+        //- Return reaction rate as volScalarField ,karam
+        //- defined at laminar.H and .C
+        virtual tmp<volScalarField> RR(label i) const
+        {
+            return tmp<volScalarField> (NULL);
+        }
+ 
+        //- Return tmp of rho
+        virtual tmp<volScalarField> rho() const = 0;
+
+
+    // IO
+
+        //- Update properties from given dictionary
+        virtual bool read();
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace combustionModels
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

libs/combustionModels_POSTECH/psiCombustionModel/psiCombustionModel/psiCombustionModelNew.C

@@ -0,0 +1,76 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-2015 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 <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "psiCombustionModel.H"
+
+// * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * * //
+
+Foam::autoPtr<Foam::combustionModels::psiCombustionModel>
+Foam::combustionModels::psiCombustionModel::New
+(
+    const fvMesh& mesh,
+    const word& phaseName
+)
+{
+    const word combModelName
+    (
+        IOdictionary
+        (
+            IOobject
+            (
+                IOobject::groupName("combustionProperties", phaseName),
+                mesh.time().constant(),
+                mesh,
+                IOobject::MUST_READ,
+                IOobject::NO_WRITE,
+                false
+            )
+        ).lookup("combustionModel")
+    );
+
+    Info<< "Selecting combustion model " << combModelName << endl;
+
+    dictionaryConstructorTable::iterator cstrIter =
+        dictionaryConstructorTablePtr_->find(combModelName);
+
+    if (cstrIter == dictionaryConstructorTablePtr_->end())
+    {
+        FatalErrorInFunction
+            << "Unknown psiCombustionModel type "
+            << combModelName << endl << endl
+            << "Valid  combustionModels are : " << endl
+            << dictionaryConstructorTablePtr_->toc()
+            << exit(FatalError);
+    }
+
+    const label tempOpen = combModelName.find('<');
+
+    const word className = combModelName(0, tempOpen);
+
+    return autoPtr<psiCombustionModel>(cstrIter()(className, mesh, phaseName));
+}
+
+
+// ************************************************************************* //

libs/combustionModels_POSTECH/psiCombustionModel/psiThermoCombustion/psiThermoCombustion.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2015 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -23,36 +23,49 @@ License
 
 \*---------------------------------------------------------------------------*/
 
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
-
+#include "psiThermoCombustion.H"
 
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
+Foam::combustionModels::psiThermoCombustion::psiThermoCombustion
+(
+    const word& modelType,
+    const fvMesh& mesh,
+    const word& phaseName
+)
+:
+    psiCombustionModel(modelType, mesh, phaseName),
+    thermoPtr_(psiReactionThermo::New(mesh, phaseName))
+{}
+
 
 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
 
-
-// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+Foam::combustionModels::psiThermoCombustion::~psiThermoCombustion()
+{}
 
 
-// * * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * //
+// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
+
+Foam::psiReactionThermo&
+Foam::combustionModels::psiThermoCombustion::thermo()
+{
+    return thermoPtr_();
+}
 
 
-// * * * * * * * * * * * * * * * Friend Functions  * * * * * * * * * * * * * //
+const Foam::psiReactionThermo&
+Foam::combustionModels::psiThermoCombustion::thermo() const
+{
+    return thermoPtr_();
+}
 
 
-// * * * * * * * * * * * * * * * Friend Operators  * * * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * * IOstream Operators  * * * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * * * Ostream Operator  * * * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+Foam::tmp<Foam::volScalarField>
+Foam::combustionModels::psiThermoCombustion::rho() const
+{
+    return thermoPtr_->rho();
+}
 
 
 // ************************************************************************* //

libs/combustionModels_POSTECH/psiCombustionModel/psiThermoCombustion/psiThermoCombustion.H

@@ -0,0 +1,113 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2012-2015 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 <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::psiThermoCombustion
+
+Description
+    Compressibility-based thermo model wrapper for combustion models
+
+SourceFiles
+    psiThermoCombustion.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef psiThermoCombustion_H
+#define psiThermoCombustion_H
+
+#include "autoPtr.H"
+#include "psiCombustionModel.H"
+#include "psiReactionThermo.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace combustionModels
+{
+
+/*---------------------------------------------------------------------------*\
+                    class psiThermoCombustion Declaration
+\*---------------------------------------------------------------------------*/
+
+class psiThermoCombustion
+:
+    public psiCombustionModel
+{
+    // Private Member Functions
+
+        //- Construct as copy (not implemented)
+        psiThermoCombustion(const psiThermoCombustion&);
+
+        //- Disallow default bitwise assignment
+        void operator=(const psiThermoCombustion&);
+
+
+protected:
+
+    // Protected data
+
+        //- Pointer to chemistry model
+        autoPtr<psiReactionThermo> thermoPtr_;
+
+
+public:
+
+    // Constructors
+
+        //- Construct from components
+        psiThermoCombustion
+        (
+            const word& modelType,
+            const fvMesh& mesh,
+            const word& phaseName
+        );
+
+
+    //- Destructor
+    virtual ~psiThermoCombustion();
+
+
+    // Member Functions
+
+        //- Return access to the thermo package
+        virtual psiReactionThermo& thermo();
+
+        //- Return const access to the thermo package
+        virtual const psiReactionThermo& thermo() const;
+
+        //- Return const access to the density field
+        virtual tmp<volScalarField> rho() const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace combustionModels
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

libs/combustionModels_POSTECH/rhoCombustionModel/rhoChemistryCombustion/rhoChemistryCombustion.C

@@ -0,0 +1,76 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2012-2015 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 <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "rhoChemistryCombustion.H"
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::combustionModels::rhoChemistryCombustion::rhoChemistryCombustion
+(
+    const word& modelType,
+    const fvMesh& mesh,
+    const word& phaseName
+)
+:
+    rhoCombustionModel(modelType, mesh, phaseName),
+    chemistryPtr_(rhoChemistryModel::New(mesh, phaseName))
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+Foam::combustionModels::rhoChemistryCombustion::~rhoChemistryCombustion()
+{}
+
+
+// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
+Foam::autoPtr<Foam::rhoChemistryModel>
+Foam::combustionModels::rhoChemistryCombustion::chem()
+{
+    return chemistryPtr_;
+}
+
+Foam::rhoReactionThermo&
+Foam::combustionModels::rhoChemistryCombustion::thermo()
+{
+    return chemistryPtr_->thermo();
+}
+
+
+const Foam::rhoReactionThermo&
+Foam::combustionModels::rhoChemistryCombustion::thermo() const
+{
+    return chemistryPtr_->thermo();
+}
+
+
+Foam::tmp<Foam::volScalarField>
+Foam::combustionModels::rhoChemistryCombustion::rho() const
+{
+    return chemistryPtr_->thermo().rho();
+}
+
+
+// ************************************************************************* //

libs/combustionModels_POSTECH/rhoCombustionModel/rhoChemistryCombustion/rhoChemistryCombustion.H

@@ -0,0 +1,113 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2012-2015 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 <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::rhoChemistryCombustion
+
+Description
+    Density-based chemistry model wrapper for combustion models
+
+SourceFiles
+    rhoChemistryCombustion.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef rhoChemistryCombustion_H
+#define rhoChemistryCombustion_H
+
+#include "autoPtr.H"
+#include "rhoCombustionModel.H"
+#include "rhoChemistryModel.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace combustionModels
+{
+
+/*---------------------------------------------------------------------------*\
+                    class rhoChemistryCombustion Declaration
+\*---------------------------------------------------------------------------*/
+
+class rhoChemistryCombustion
+:
+    public rhoCombustionModel
+{
+    // Private Member Functions
+
+        //- Construct as copy (not implemented)
+        rhoChemistryCombustion(const rhoChemistryCombustion&);
+
+        //- Disallow default bitwise assignment
+        void operator=(const rhoChemistryCombustion&);
+
+
+protected:
+
+    // Protected data
+
+        //- Pointer to chemistry model
+        autoPtr<rhoChemistryModel> chemistryPtr_;
+
+
+public:
+
+    // Constructors
+
+        //- Construct from components and thermo
+        rhoChemistryCombustion
+        (
+            const word& modelType,
+            const fvMesh& mesh,
+            const word& phaseName
+        );
+
+
+    //- Destructor
+    virtual ~rhoChemistryCombustion();
+
+
+    // Member Functions
+        autoPtr<rhoChemistryModel> chem(); 
+        //- Return access to the thermo package
+        virtual rhoReactionThermo& thermo();
+
+        //- Return const access to the thermo package
+        virtual const rhoReactionThermo& thermo() const;
+
+        //- Return const access to the density field
+        virtual tmp<volScalarField> rho() const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace combustionModels
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

libs/combustionModels_POSTECH/rhoCombustionModel/rhoCombustionModel/rhoCombustionModel.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2015 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -23,36 +23,51 @@ License
 
 \*---------------------------------------------------------------------------*/
 
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+#include "rhoCombustionModel.H"
 
-// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+/* * * * * * * * * * * * * * * private static data * * * * * * * * * * * * * */
 
+namespace Foam
+{
+namespace combustionModels
+{
+    defineTypeNameAndDebug(rhoCombustionModel, 0);
+    defineRunTimeSelectionTable(rhoCombustionModel, dictionary);
+}
+}
 
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
+Foam::combustionModels::rhoCombustionModel::rhoCombustionModel
+(
+    const word& modelType,
+    const fvMesh& mesh,
+    const word& phaseName
+)
+:
+    combustionModel(modelType, mesh, phaseName)
+{}
+
 
 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
 
-
-// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+Foam::combustionModels::rhoCombustionModel::~rhoCombustionModel()
+{}
 
 
-// * * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * //
+// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
 
-
-// * * * * * * * * * * * * * * * Friend Functions  * * * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * * * Friend Operators  * * * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * * IOstream Operators  * * * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * * * Ostream Operator  * * * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+bool Foam::combustionModels::rhoCombustionModel::read()
+{
+    if (combustionModel::read())
+    {
+        return true;
+    }
+    else
+    {
+        return false;
+    }
+}
 
 
 // ************************************************************************* //

libs/combustionModels_POSTECH/rhoCombustionModel/rhoCombustionModel/rhoCombustionModel.H

@@ -0,0 +1,150 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-2015 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 <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::rhoCombustionModel
+
+Description
+    Combustion models for rho-based thermodynamics
+
+SourceFiles
+    rhoCombustionModelI.H
+    rhoCombustionModel.C
+    rhoCombustionModelNew.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef rhoCombustionModel_H
+#define rhoCombustionModel_H
+
+#include "combustionModel.H"
+#include "autoPtr.H"
+#include "runTimeSelectionTables.H"
+#include "rhoReactionThermo.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace combustionModels
+{
+
+/*---------------------------------------------------------------------------*\
+                     class rhoCombustionModel Declaration
+\*---------------------------------------------------------------------------*/
+
+class rhoCombustionModel
+:
+    public combustionModel
+{
+
+    // Private Member Functions
+
+        //- Construct as copy (not implemented)
+        rhoCombustionModel(const rhoCombustionModel&);
+
+        //- Disallow default bitwise assignment
+        void operator=(const rhoCombustionModel&);
+
+
+public:
+
+    //- Runtime type information
+    TypeName("rhoCombustionModel");
+
+
+    //- Declare run-time constructor selection tables
+    declareRunTimeSelectionTable
+    (
+        autoPtr,
+        rhoCombustionModel,
+        dictionary,
+        (
+            const word& modelType,
+            const fvMesh& mesh,
+            const word& phaseName
+        ),
+        (modelType, mesh, phaseName)
+    );
+
+
+    // Constructors
+
+
+        //- Construct from components
+        rhoCombustionModel
+        (
+            const word& modelType,
+            const fvMesh& mesh,
+            const word& phaseName
+        );
+
+
+
+    //- Selector
+    static autoPtr<rhoCombustionModel> New
+    (
+        const fvMesh& mesh,
+        const word& phaseName=word::null
+    );
+
+
+    //- Destructor
+    virtual ~rhoCombustionModel();
+
+
+    // Member Functions
+
+        // Access functions
+
+            //- Access combustion dict
+             inline const dictionary& coeff() const;
+
+
+        //- Return access to the thermo package
+        virtual rhoReactionThermo& thermo() = 0;
+
+        //- Return const access to the thermo package
+        virtual const rhoReactionThermo& thermo() const = 0;
+
+        //- Return tmp of rho
+        virtual tmp<volScalarField> rho() const = 0;
+
+
+     // IO
+
+        //- Update properties from given dictionary
+        virtual bool read();
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace combustionModels
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

libs/combustionModels_POSTECH/rhoCombustionModel/rhoCombustionModel/rhoCombustionModelNew.C

@@ -0,0 +1,76 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-2015 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 <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "rhoCombustionModel.H"
+
+// * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * * //
+
+Foam::autoPtr<Foam::combustionModels::rhoCombustionModel>
+Foam::combustionModels::rhoCombustionModel::New
+(
+    const fvMesh& mesh,
+    const word& phaseName
+)
+{
+    const word combTypeName
+    (
+        IOdictionary
+        (
+            IOobject
+            (
+                IOobject::groupName("combustionProperties", phaseName),
+                mesh.time().constant(),
+                mesh,
+                IOobject::MUST_READ,
+                IOobject::NO_WRITE,
+                false
+            )
+        ).lookup("combustionModel")
+    );
+
+    Info<< "Selecting combustion model " << combTypeName << endl;
+
+    dictionaryConstructorTable::iterator cstrIter =
+        dictionaryConstructorTablePtr_->find(combTypeName);
+
+    if (cstrIter == dictionaryConstructorTablePtr_->end())
+    {
+        FatalErrorInFunction
+            << "Unknown rhoCombustionModel type "
+            << combTypeName << endl << endl
+            << "Valid  combustionModels are : " << endl
+            << dictionaryConstructorTablePtr_->toc()
+            << exit(FatalError);
+    }
+
+    const label tempOpen = combTypeName.find('<');
+
+    const word className = combTypeName(0, tempOpen);
+
+    return autoPtr<rhoCombustionModel> (cstrIter()(className, mesh, phaseName));
+}
+
+
+// ************************************************************************* //

libs/combustionModels_POSTECH/rhoCombustionModel/rhoThermoCombustion/rhoThermoCombustion.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2015 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -23,36 +23,49 @@ License
 
 \*---------------------------------------------------------------------------*/
 
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
-
+#include "rhoThermoCombustion.H"
 
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
+Foam::combustionModels::rhoThermoCombustion::rhoThermoCombustion
+(
+    const word& modelType,
+    const fvMesh& mesh,
+    const word& phaseName
+)
+:
+    rhoCombustionModel(modelType, mesh, phaseName),
+    thermoPtr_(rhoReactionThermo::New(mesh, phaseName))
+{}
+
 
 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
 
-
-// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+Foam::combustionModels::rhoThermoCombustion::~rhoThermoCombustion()
+{}
 
 
-// * * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * //
+// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
+
+Foam::rhoReactionThermo&
+Foam::combustionModels::rhoThermoCombustion::thermo()
+{
+    return thermoPtr_();
+}
 
 
-// * * * * * * * * * * * * * * * Friend Functions  * * * * * * * * * * * * * //
+const Foam::rhoReactionThermo&
+Foam::combustionModels::rhoThermoCombustion::thermo() const
+{
+    return thermoPtr_();
+}
 
 
-// * * * * * * * * * * * * * * * Friend Operators  * * * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * * IOstream Operators  * * * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * * * Ostream Operator  * * * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+Foam::tmp<Foam::volScalarField>
+Foam::combustionModels::rhoThermoCombustion::rho() const
+{
+    return thermoPtr_().rho();
+}
 
 
 // ************************************************************************* //

libs/combustionModels_POSTECH/rhoCombustionModel/rhoThermoCombustion/rhoThermoCombustion.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2015 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -21,48 +21,89 @@ License
     You should have received a copy of the GNU General Public License
     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
 
-Function
-    Foam::CMC::AMC
+Class
+    Foam::rhoThermoCombustion
 
 Description
+    Density-based thermo model wrapper for combustion models
 
 SourceFiles
-    AMC.C
+    rhoThermoCombustion.C
 
 \*---------------------------------------------------------------------------*/
 
-#ifndef AMC_H
-#define AMC_H
+#ifndef rhoThermoCombustion_H
+#define rhoThermoCombustion_H
 
-#include "tmp.H"
-#include "scalarField.H"
+#include "autoPtr.H"
+#include "rhoCombustionModel.H"
+#include "rhoChemistryModel.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 namespace Foam
 {
-namespace CMC
+namespace combustionModels
 {
 
 /*---------------------------------------------------------------------------*\
-                         Function AMC Declaration
+                    class rhoThermoCombustion Declaration
 \*---------------------------------------------------------------------------*/
 
-    const scalar pi = 3.141592;
+class rhoThermoCombustion
+:
+    public rhoCombustionModel
+{
+    // Private Member Functions
 
-    const scalar spi = Foam::sqrt(pi);
+        //- Construct as copy (not implemented)
+        rhoThermoCombustion(const rhoThermoCombustion&);
 
-    //Amplitude Mapping Closure (from KIVA)
-    tmp<scalarField> AMC(const scalarField& eta);
+        //- Disallow default bitwise assignment
+        void operator=(const rhoThermoCombustion&);
 
 
-    //Amplitude Mapping Closure (from KIVA)
-    //Define exp(-2*(erf^-1(2*eta - 1))^2)
-    inline scalar AMC(const scalar eta);
+protected:
+
+    // Protected data
+
+        //- Pointer to thermo model
+        autoPtr<rhoReactionThermo> thermoPtr_;
+
+
+public:
+
+    // Constructors
+
+        //- Construct from components
+        rhoThermoCombustion
+        (
+            const word& modelType,
+            const fvMesh& mesh,
+            const word& phaseName
+        );
+
+
+    //- Destructor
+    virtual ~rhoThermoCombustion();
+
+
+    // Member Functions
+
+        //- Return access to the thermo package
+        virtual rhoReactionThermo& thermo();
+
+        //- Return const access to the thermo package
+        virtual const rhoReactionThermo& thermo() const;
+
+        //- Return const access to the density field
+        virtual tmp<volScalarField> rho() const;
+};
+
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
-} // End namespace CMC
+} // End namespace combustionModels
 } // End namespace Foam
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

libs/combustionModels_POSTECH/singleStepCombustion/singleStepCombustion.C

@@ -0,0 +1,189 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-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 <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "singleStepCombustion.H"
+#include "fvmSup.H"
+
+namespace Foam
+{
+namespace combustionModels
+{
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+template<class CombThermoType, class ThermoType>
+singleStepCombustion<CombThermoType, ThermoType>::singleStepCombustion
+(
+    const word& modelType,
+    const fvMesh& mesh,
+    const word& phaseName
+)
+:
+    CombThermoType(modelType, mesh, phaseName),
+    singleMixturePtr_(NULL),
+    wFuel_
+    (
+        IOobject
+        (
+            IOobject::groupName("wFuel", phaseName),
+            this->mesh().time().timeName(),
+            this->mesh(),
+            IOobject::NO_READ,
+            IOobject::NO_WRITE
+        ),
+        this->mesh(),
+        dimensionedScalar("zero", dimMass/dimVolume/dimTime, 0.0)
+    ),
+    semiImplicit_(readBool(this->coeffs_.lookup("semiImplicit")))
+{
+    if (isA<singleStepReactingMixture<ThermoType>>(this->thermo()))
+    {
+        singleMixturePtr_ =
+            &dynamic_cast<singleStepReactingMixture<ThermoType>&>
+            (
+                this->thermo()
+            );
+    }
+    else
+    {
+        FatalErrorInFunction
+            << "Inconsistent thermo package for " << this->type() << " model:\n"
+            << "    " << this->thermo().type() << nl << nl
+            << "Please select a thermo package based on "
+            << "singleStepReactingMixture" << exit(FatalError);
+    }
+
+    if (semiImplicit_)
+    {
+        Info<< "Combustion mode: semi-implicit" << endl;
+    }
+    else
+    {
+        Info<< "Combustion mode: explicit" << endl;
+    }
+}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+template<class CombThermoType, class ThermoType>
+singleStepCombustion<CombThermoType, ThermoType>::~singleStepCombustion()
+{}
+
+
+// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
+
+template<class CombThermoType, class ThermoType>
+tmp<fvScalarMatrix> singleStepCombustion<CombThermoType, ThermoType>::R
+(
+    volScalarField& Y
+) const
+{
+    const label specieI =
+        this->thermoPtr_->composition().species()[Y.member()];
+
+    volScalarField wSpecie
+    (
+        wFuel_*singleMixturePtr_->specieStoichCoeffs()[specieI]
+    );
+
+    if (semiImplicit_)
+    {
+        const label fNorm = singleMixturePtr_->specieProd()[specieI];
+        const volScalarField fres(singleMixturePtr_->fres(specieI));
+        wSpecie /= max(fNorm*(Y - fres), scalar(1e-2));
+
+        return -fNorm*wSpecie*fres + fNorm*fvm::Sp(wSpecie, Y);
+    }
+    else
+    {
+        return wSpecie + fvm::Sp(0.0*wSpecie, Y);
+    }
+}
+
+
+template<class CombThermoType, class ThermoType>
+tmp<volScalarField>
+singleStepCombustion<CombThermoType, ThermoType>::Sh() const
+{
+    const label fuelI = singleMixturePtr_->fuelIndex();
+    volScalarField& YFuel =
+        const_cast<volScalarField&>(this->thermoPtr_->composition().Y(fuelI));
+
+    return -singleMixturePtr_->qFuel()*(R(YFuel) & YFuel);
+}
+
+
+template<class CombThermoType, class ThermoType>
+tmp<volScalarField>
+singleStepCombustion<CombThermoType, ThermoType>::dQ() const
+{
+    tmp<volScalarField> tdQ
+    (
+        new volScalarField
+        (
+            IOobject
+            (
+                IOobject::groupName("dQ", this->phaseName_),
+                this->mesh_.time().timeName(),
+                this->mesh_,
+                IOobject::NO_READ,
+                IOobject::NO_WRITE,
+                false
+            ),
+            this->mesh_,
+            dimensionedScalar("dQ", dimEnergy/dimTime, 0.0)
+        )
+    );
+
+    if (this->active())
+    {
+        volScalarField& dQ = tdQ.ref();
+        dQ.ref() = this->mesh().V()*Sh()();
+    }
+    return tdQ;
+}
+
+
+template<class CombThermoType, class ThermoType>
+bool singleStepCombustion<CombThermoType, ThermoType>::read()
+{
+    if (CombThermoType::read())
+    {
+        return true;
+    }
+    else
+    {
+        return false;
+    }
+}
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace combustionModels
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

libs/combustionModels_POSTECH/singleStepCombustion/singleStepCombustion.H

@@ -0,0 +1,132 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-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 <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::combustionModels::singleStepCombustion
+
+Description
+    Base class for combustion models using singleStepReactingMixture.
+
+SourceFiles
+    singleStepCombustion.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef singleStepCombustion_H
+#define singleStepCombustion_H
+
+#include "singleStepReactingMixture.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace combustionModels
+{
+
+/*---------------------------------------------------------------------------*\
+                    Class singleStepCombustion Declaration
+\*---------------------------------------------------------------------------*/
+
+template<class CombThermoType, class ThermoType>
+class singleStepCombustion
+:
+    public CombThermoType
+{
+    // Private Member Functions
+
+        //- Disallow copy construct
+        singleStepCombustion(const singleStepCombustion&);
+
+        //- Disallow default bitwise assignment
+        void operator=(const singleStepCombustion&);
+
+
+protected:
+
+    // Protected data
+
+        //- Pointer to singleStepReactingMixture mixture
+        singleStepReactingMixture<ThermoType>* singleMixturePtr_;
+
+        //- Fuel consumption rate
+        volScalarField wFuel_;
+
+        //- Semi-implicit (true) or explicit (false) treatment
+        bool semiImplicit_;
+
+
+public:
+
+    // Constructors
+
+        //- Construct from components
+        singleStepCombustion
+        (
+            const word& modelType,
+            const fvMesh& mesh,
+            const word& phaseName
+        );
+
+
+    //- Destructor
+    virtual ~singleStepCombustion();
+
+
+    // Member Functions
+
+        // Evolution
+
+            //- Fuel consumption rate matrix
+            virtual tmp<fvScalarMatrix> R(volScalarField& Y) const;
+
+            //- Heat release rate calculated from fuel consumption rate matrix
+            virtual tmp<volScalarField> dQ() const;
+
+            //- Sensible enthalpy source term
+            virtual tmp<volScalarField> Sh() const;
+
+
+        // IO
+
+            //- Update properties from given dictionary
+            virtual bool read();
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace combustionModels
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+    #include "singleStepCombustion.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

libs/thermos/Make/files

@@ -1,3 +0,0 @@
-thermoAdd.C
-
-LIB = $(FOAM_USER_LIBBIN)/libthermoAdd

libs/thermos/Make/options

@@ -1,29 +0,0 @@
-OFP_LIB_SRC = ../
-
-EXE_INC = \
-    -I$(LIB_SRC)/transportModels/compressible/lnInclude \
-    -I$(OFP_LIB_SRC)//chemistryModel_POSTECH/lnInclude \
-    -I$(LIB_SRC)/thermophysicalModels/reactionThermo/lnInclude \
-    -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
-    -I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \
-    -I$(LIB_SRC)/thermophysicalModels/functions/Polynomial \
-    -I$(LIB_SRC)/thermophysicalModels/thermophysicalFunctions/lnInclude \
-    -I$(LIB_SRC)/turbulenceModels/compressible/lnInclude \
-    -I$(LIB_SRC)/ODE/lnInclude \
-    -I$(LIB_SRC)/thermophysicalModels/solidSpecie/lnInclude \
-    -I$(LIB_SRC)/finiteVolume/lnInclude
-
-LIB_LIBS = \
-    -lcompressibleTransportModels \
-    -lfluidThermophysicalModels \
-    -lreactionThermophysicalModels \
-    -lspecie \
-    -lthermophysicalFunctions \
-    -lODE \
-    -lfiniteVolume \
-    -lcompressibleTransportModels \
-    -lfluidThermophysicalModels \
-    -lspecie \
-    -lsolidSpecie \
-    -lfiniteVolume \
-    -lOpenFOAM

libs/thermos/thermoAdd.C

@@ -1,384 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 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 <http://www.gnu.org/licenses/>.
-
-Typedefs
-    Foam::thermoPhysicsTypes
-
-Description
-    Type definitions for thermo-physics models
-
-\*---------------------------------------------------------------------------*/
-
-#include "specie.H"
-#include "perfectGas.H"
-#include "janafThermo.H"
-
-#include "absoluteEnthalpy.H"
-#include "absoluteInternalEnergy.H"
-#include "thermo.H"
-#include "sutherlandTransport.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-    // thermo physics types based on absoluteEnthalpy
-
-    typedef
-    sutherlandTransport
-    <
-        species::thermo
-        <
-            janafThermo
-            <
-                perfectGas<specie>
-            >,
-            absoluteEnthalpy
-        >
-    > gasHaThermoPhysics;
-
-
-    // thermo physics types based on absoluteInternalEnergy
-
-    typedef
-    sutherlandTransport
-    <
-        species::thermo
-        <
-            janafThermo
-            <
-                perfectGas<specie>
-            >,
-            absoluteInternalEnergy
-        >
-    > gasEaThermoPhysics;
-
-} // End namespace Foam
-
-/*---------------------------------------------------------------------------*\
-Typedefs
-    Foam::reactionTypes
-
-Description
-    Type definitions for reactions
-\*---------------------------------------------------------------------------*/
-
-#include "Reaction.H"
-
-namespace Foam
-{
-    // absolute enthalpy based reactions
-    typedef Reaction<gasHaThermoPhysics> gasHaReaction;
-
-    // absolute internal ennergy based reactions
-    typedef Reaction<gasEaThermoPhysics> gasEaReaction;
-
-}
-
-// ************************************************************************* //
-
-#include "makeReactionThermo.H"
-
-#include "chemistryReader.H"
-#include "foamChemistryReader.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-
-// Gas chemistry readers based on absoluteEnthalpy
-
-makeChemistryReader(gasHaThermoPhysics);
-
-makeChemistryReaderType(foamChemistryReader, gasHaThermoPhysics);
-
-
-
-// Gas chemistry readers based on absoluteInternalEnergy
-
-makeChemistryReader(gasEaThermoPhysics);
-
-makeChemistryReaderType(foamChemistryReader, gasEaThermoPhysics);
-
-
-} // End namespace Foam
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#include "makeReaction.H"
-
-#include "ArrheniusReactionRate.H"
-#include "infiniteReactionRate.H"
-#include "LandauTellerReactionRate.H"
-#include "thirdBodyArrheniusReactionRate.H"
-
-#include "ChemicallyActivatedReactionRate.H"
-#include "JanevReactionRate.H"
-#include "powerSeriesReactionRate.H"
-
-#include "FallOffReactionRate.H"
-#include "LindemannFallOffFunction.H"
-#include "SRIFallOffFunction.H"
-#include "TroeFallOffFunction.H"
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#define makeReactions(Thermo, Reaction)                                        \
-                                                                               \
-    defineTemplateTypeNameAndDebug(Reaction, 0);                               \
-    defineTemplateRunTimeSelectionTable(Reaction, Istream);                    \
-    defineTemplateRunTimeSelectionTable(Reaction, dictionary);                 \
-                                                                               \
-    makeIRNReactions(Thermo, ArrheniusReactionRate)                            \
-    makeIRNReactions(Thermo, infiniteReactionRate)                             \
-    makeIRNReactions(Thermo, LandauTellerReactionRate)                         \
-    makeIRNReactions(Thermo, thirdBodyArrheniusReactionRate)                   \
-                                                                               \
-    makeIRReactions(Thermo, JanevReactionRate)                                 \
-    makeIRReactions(Thermo, powerSeriesReactionRate)                           \
-                                                                               \
-    makePressureDependentReactions                                             \
-    (                                                                          \
-       Thermo,                                                                 \
-       ArrheniusReactionRate,                                                  \
-       LindemannFallOffFunction                                                \
-    )                                                                          \
-                                                                               \
-    makePressureDependentReactions                                             \
-    (                                                                          \
-       Thermo,                                                                 \
-       ArrheniusReactionRate,                                                  \
-       TroeFallOffFunction                                                     \
-    )                                                                          \
-                                                                               \
-    makePressureDependentReactions                                             \
-    (                                                                          \
-       Thermo,                                                                 \
-       ArrheniusReactionRate,                                                  \
-       SRIFallOffFunction                                                      \
-    )
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-    // absolute enthalpy based reactions
-    makeReactions(gasHaThermoPhysics, gasHaReaction)
-    makeReactions(gasEaThermoPhysics, gasEaReaction)
-}
-
-// ************************************************************************* //
-
-#include "makeReactionThermo.H"
-
-#include "psiReactionThermo.H"
-#include "hePsiThermo.H"
-
-#include "specie.H"
-#include "perfectGas.H"
-#include "janafThermo.H"
-#include "absoluteEnthalpy.H"
-#include "thermo.H"
-#include "sutherlandTransport.H"
-
-#include "homogeneousMixture.H"
-#include "inhomogeneousMixture.H"
-#include "veryInhomogeneousMixture.H"
-#include "multiComponentMixture.H"
-#include "reactingMixture.H"
-#include "singleStepReactingMixture.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-
-// sutherlandTransport, janafThermo
-
-makeReactionThermo
-(
-    psiThermo,
-    psiReactionThermo,
-    hePsiThermo,
-    homogeneousMixture,
-    sutherlandTransport,
-    absoluteEnthalpy,
-    janafThermo,
-    perfectGas,
-    specie
-);
-
-makeReactionThermo
-(
-    psiThermo,
-    psiReactionThermo,
-    hePsiThermo,
-    inhomogeneousMixture,
-    sutherlandTransport,
-    absoluteEnthalpy,
-    janafThermo,
-    perfectGas,
-    specie
-);
-
-makeReactionThermo
-(
-    psiThermo,
-    psiReactionThermo,
-    hePsiThermo,
-    veryInhomogeneousMixture,
-    sutherlandTransport,
-    absoluteEnthalpy,
-    janafThermo,
-    perfectGas,
-    specie
-);
-
-
-// Multi-component thermo for absolute enthalpy
-
-makeReactionMixtureThermo
-(
-    psiThermo,
-    psiReactionThermo,
-    hePsiThermo,
-    multiComponentMixture,
-    gasHaThermoPhysics
-);
-
-
-// Multi-component thermo for internal energy
-
-makeReactionMixtureThermo
-(
-    psiThermo,
-    psiReactionThermo,
-    hePsiThermo,
-    multiComponentMixture,
-    gasEaThermoPhysics
-);
-
-
-// Multi-component reaction thermo for absolute enthalpy
-
-makeReactionMixtureThermo
-(
-    psiThermo,
-    psiReactionThermo,
-    hePsiThermo,
-    reactingMixture,
-    gasHaThermoPhysics
-);
-
-makeReactionMixtureThermo
-(
-    psiThermo,
-    psiReactionThermo,
-    hePsiThermo,
-    singleStepReactingMixture,
-    gasHaThermoPhysics
-);
-
-
-// Multi-component reaction thermo for internal energy
-
-makeReactionMixtureThermo
-(
-    psiThermo,
-    psiReactionThermo,
-    hePsiThermo,
-    reactingMixture,
-    gasEaThermoPhysics
-);
-
-makeReactionMixtureThermo
-(
-    psiThermo,
-    psiReactionThermo,
-    hePsiThermo,
-    singleStepReactingMixture,
-    gasEaThermoPhysics
-);
-
-} // End namespace Foam
-
-/*---------------------------------------------------------------------------*\
-InClass
-    Foam::psiChemistryModel
-
-Description
-    Creates chemistry model instances templated on the type of thermodynamics
-\*---------------------------------------------------------------------------*/
-
-#include "makeChemistryModel.H"
-
-#include "psiChemistryModel.H"
-#include "chemistryModel.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-    // Chemistry moldels based on absoluteEnthalpy
-    makeChemistryModel
-    (
-        chemistryModel,
-        psiChemistryModel,
-        gasHaThermoPhysics
-    );
-
-    // Chemistry moldels based on absoluteInternalEnergy
-    makeChemistryModel
-    (
-        chemistryModel,
-        psiChemistryModel,
-        gasEaThermoPhysics
-    );
-}
-
-// ************************************************************************* //
-
-#include "makeChemistrySolverTypes.H"
-
-#include "psiChemistryModel.H"
-#include "rhoChemistryModel.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-    // Chemistry solvers based on absoluteEnthalpy
-    makeChemistrySolverTypes(psiChemistryModel, gasHaThermoPhysics);
-    // makeChemistrySolverTypes(rhoChemistryModel, gasHaThermoPhysics);
-
-    // Chemistry solvers based on absoluteInternalEnergy
-    makeChemistrySolverTypes(psiChemistryModel, gasEaThermoPhysics);
-    // makeChemistrySolverTypes(rhoChemistryModel, gasEaThermoPhysics);
-}
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-

solvers_post/AMC/AMC.C

@@ -1,68 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  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 <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-#include "AMC.H"
-
-#include "error.H"
-
-// * * * * * * * * * * * * * Functions * * * * * * * * * * * * //
-
-//Amplitude Mapping Closure (from KIVA)
-//Define exp(-2*(erf^-1(2*eta - 1))^2)
-inline Foam::scalar Foam::CMC::AMC(const scalar eta)
-{
-    const scalar a0 = (2.0*eta - 1.0);
-    scalar a1 = 0.5;
-    scalar slope;
-    scalar da = GREAT;
-
-    while(mag(da) > SMALL)
-    {
-        slope = (2.0/spi)*Foam::exp(-1.0*Foam::pow(a1,2.0));
-        da = (a0 - Foam::erf(a1))/slope;
-        a1 = a1 + da;
-    }
-
-    return Foam::exp(-2.0*Foam::pow(a1,2.0));
-}
-
-//Amplitude Mapping Closure (from KIVA)
-Foam::tmp<Foam::scalarField> Foam::CMC::AMC(const scalarField& eta)
-{
-    tmp<scalarField> tRes(new scalarField(eta.size(), 0.0));
-
-    scalarField& Res = tRes.ref();
-
-    for(label i=0 ; i<eta.size() ; i++)
-    {
-        Res[i] = AMC(eta[i]);
-    }
-
-    return tRes;
-}
-
-
-
-// ************************************************************************* //

solvers_post/BetaFunction/BetaFunction.C

@@ -1,120 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  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 <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-#include "BetaFunction.H"
-
-#include "error.H"
-
-// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
-
-// const scalar Foam::CMC::BetaFunction::staticData();
-
-
-// * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
-
-void Foam::CMC::BetaFunction::_setParameters()
-{
-        scalar gamma = mf_*(1.0 - mf_)/(mfVar_ + SMALL) - 1.0;
-
-        if (gamma < SMALL || (mfVar_)/(mf_*(1 - mf_) + SMALL) < 0.001)
-        {
-            fdelta_ = true;
-        }
-        else
-        {
-            a_ = max(mf_*gamma, 0.0);
-            b_ = max((1.0 - mf_)*gamma, 0.0);
-
-            if(a_ < 1.0)
-            {
-                delta0_ = true;
-            }
-            if(b_ < 1.0)
-            {
-                delta1_ = true;
-            }
-
-            _limitAB();
-        }
-}
-
-
-void Foam::CMC::BetaFunction::_limitAB()
-{
-    scalar fmax = 1.0/(1.0 + (b_ - 1.0)/(a_ - 1.0));
-
-    if(a_ > 500.0)
-    {
-        a_ = 500.0;
-        b_ = (a_ - 1.0 - fmax*(a_ - 2.0))/fmax;
-    }
-    else if(b_ > 500.0)
-    {
-        b_ = 500.0;
-        a_ = (1.0 + fmax*(b_ - 2.0))/(1.0 - fmax);
-    }
-}
-
-// * * * * * * * * * * * * Protected Member Functions  * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
-
-Foam::CMC::BetaFunction::BetaFunction(const scalar mf, const scalar mfVar)
-:
-    mf_(mf),
-    mfVar_(mfVar),
-    a_(0.0),
-    b_(0.0),
-    fdelta_(false),
-    delta0_(false),
-    delta1_(false)
-{
-    _setParameters();
-}
-
-
-// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
-
-Foam::CMC::BetaFunction::~BetaFunction()
-{}
-
-
-// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * * Friend Functions  * * * * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * * //
-
-
-// ************************************************************************* //

solvers_post/BetaFunction/BetaFunction.H

@@ -1,198 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  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 <http://www.gnu.org/licenses/>.
-
-Class
-    Foam::BetaFunction
-
-Description
-
-SourceFiles
-    BetaFunctionI.H
-    BetaFunction.C
-    BetaFunctionIO.C
-
-\*---------------------------------------------------------------------------*/
-
-#ifndef BetaFunction_H
-#define BetaFunction_H
-
-#include "scalar.H"
-#include "Switch.H"
-#include "scalarField.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-
-// Forward declaration of classes
-class Istream;
-class Ostream;
-
-// Forward declaration of friend functions and operators
-namespace CMC
-{
-    class BetaFunction;
-}
-Istream& operator>>(Istream&, CMC::BetaFunction&);
-Ostream& operator<<(Ostream&, const CMC::BetaFunction&);
-
-
-/*---------------------------------------------------------------------------*\
-                         Class BetaFunction Declaration
-\*---------------------------------------------------------------------------*/
-
-namespace CMC
-{
-
-class BetaFunction
-{
-    // Private data
-
-        //- Mean mixture fraction
-        scalar mf_;
-
-        //- Mixture fraction variance
-        scalar mfVar_;
-
-        //- Beta pdf parameter alpha
-        scalar a_;
-
-        //- Beta pdf parameter beta
-        scalar b_;
-
-        //- Flag beta pdf special case: forced-delta function 
-        Switch fdelta_;
-
-        //- delta function at pure oxidizer stream
-        Switch delta0_;
-
-        //- delta function at pure fuel stream
-        Switch delta1_;
-
-
-    // Private Member Functions
-
-        //- Calculate alpha and beta and set special case flags
-        void _setParameters();
-
-        //- Limit too large value of beta pdf parameters
-        void _limitAB();
-
-        //- Disallow default bitwise copy construct
-        // BetaFunction(const BetaFunction&);
-
-        //- Disallow default bitwise assignment
-        // void operator=(const BetaFunction&);
-
-
-public:
-
-    // Static data members
-
-        //- Static data staticData
-        // static const scalar staticData;
-
-
-    // Constructors
-
-        //- Construct from components
-        BetaFunction(const scalar mf, const scalar mfVar);
-
-        //- Construct from Istream
-        BetaFunction(Istream&);
-
-        //- Construct as copy
-        // BetaFunction(const BetaFunction&);
-
-
-    //- Destructor
-    ~BetaFunction();
-
-
-    // Member Functions
-
-        scalarField etaFunc(const UList<scalar>& eta) const
-        {
-            return pow(eta, a_-1.0)*pow(1.0-eta, b_-1.0);
-        }
-
-        scalar etaFunc(const scalar eta) const
-        {
-            return Foam::pow(eta, a_-1.0)*Foam::pow(1.0-eta, b_-1.0);
-        }
-
-
-        // Access
-
-            Switch fdelta() const {return fdelta_;}
-
-            Switch delta0() const {return delta0_;}
-
-            Switch delta1() const {return delta1_;}
-
-            scalar mf() const {return mf_;}
-
-            scalar mfVar() const {return mfVar_;}
-
-            scalar alpha() const {return a_;}
-
-            scalar beta() const {return b_;}
-
-        // Check
-
-        // Edit
-
-        // Write
-
-
-    // Member Operators
-
-        // void operator=(const BetaFunction&);
-
-
-    // Friend Functions
-
-    // Friend Operators
-
-    // IOstream Operators
-
-        friend Istream& ::Foam::operator>>(Istream&, BetaFunction&);
-        friend Ostream& ::Foam::operator<<(Ostream&, const BetaFunction&);
-};
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace CMC
-} // End namespace Foam
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#include "BetaFunctionI.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#endif
-
-// ************************************************************************* //

solvers_post/BetaFunction/BetaFunctionIO.C

@@ -1,65 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  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 <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-#include "BetaFunction.H"
-#include "IOstreams.H"
-
-// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
-
-Foam::CMC::BetaFunction::BetaFunction(Istream& is)
-{
-    // Check state of Istream
-    is.check("Foam::CMC::BetaFunction::BetaFunction(Foam::Istream&)");
-}
-
-
-// * * * * * * * * * * * * * * * IOstream Operators  * * * * * * * * * * * * //
-
-Foam::Istream& Foam::operator>>(Istream& is, CMC::BetaFunction&)
-{
-    // Check state of Istream
-    is.check
-    (
-        "Foam::Istream& Foam::operator>>(Foam::Istream&, Foam::CMC::BetaFunction&)"
-    );
-
-    return is;
-}
-
-
-Foam::Ostream& Foam::operator<<(Ostream& os, const CMC::BetaFunction&)
-{
-    // Check state of Ostream
-    os.check
-    (
-        "Foam::Ostream& Foam::operator<<(Foam::Ostream&, "
-        "const Foam::CMC::BetaFunction&)"
-    );
-
-    return os;
-}
-
-
-// ************************************************************************* //

solvers_post/BetaGrid/BetaGrid.C

@@ -1,126 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  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 <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-#include "BetaGrid.H"
-
-#include "error.H"
-#include "dictionary.H"
-#include "SubList.H"
-
-// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
-
-// const scalar Foam::CMC::BetaGrid::staticData();
-
-
-// * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * Protected Member Functions  * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
-
-Foam::CMC::BetaGrid::BetaGrid(const dictionary& dict)
-:
-    etaCut_(dict.lookup("detailedEta")),
-    nSpacings_(dict.lookup("detailedN")),
-    etaSpace_(sum(nSpacings_)+3)
-{
-    // Check for input list size
-    if (etaCut_.size() != nSpacings_.size() + 1)
-    {
-        FatalErrorInFunction
-            << "Number of grid interval points does not match number of grid sizes"
-            << abort(FatalError);
-    }
-
-    // Check for evenness of number of spacings in each interval
-    forAll (nSpacings_, i)
-    {
-        if ((nSpacings_[i] % 2) == 1)
-        {
-            FatalErrorInFunction
-                << "Number of spacings in intervals should be even"
-                << abort(FatalError);
-        }
-    }
-
-    forAll(nSpacings_, i)
-    {
-        const label ns = nSpacings_[i];
-        const labelList::subList prev(nSpacings_, i);
-        const label baseI = sum(prev) + 1;
-
-        const scalar low = etaCut_[i];
-        const scalar upp = etaCut_[i+1];
-        const scalar delta = (upp - low) / scalar(ns);
-
-        for (label j = 0; j < ns; j++)
-        {
-            etaSpace_[baseI+j] = low + delta * j;
-        }
-    }
-    etaSpace_[0] = 0.0;
-    etaSpace_[etaSpace_.size()-2] = etaCut_.last();
-    etaSpace_[etaSpace_.size()-1] = 1.0;
-}
-
-
-Foam::CMC::BetaGrid::BetaGrid(const BetaGrid&)
-{}
-
-
-// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
-
-Foam::CMC::BetaGrid::~BetaGrid()
-{}
-
-
-// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * * //
-
-void Foam::CMC::BetaGrid::operator=(const BetaGrid& rhs)
-{
-    // Check for assignment to self
-    if (this == &rhs)
-    {
-        FatalErrorInFunction
-            << "Attempted assignment to self"
-            << abort(FatalError);
-    }
-}
-
-// * * * * * * * * * * * * * * Friend Functions  * * * * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * * //
-
-
-// ************************************************************************* //

solvers_post/BetaGrid/BetaGridIO.C

@@ -1,65 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  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 <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-#include "BetaGrid.H"
-#include "IOstreams.H"
-
-// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
-
-Foam::CMC::BetaGrid::BetaGrid(Istream& is)
-{
-    // Check state of Istream
-    is.check("Foam::CMC::BetaGrid::BetaGrid(Foam::Istream&)");
-}
-
-
-// * * * * * * * * * * * * * * * IOstream Operators  * * * * * * * * * * * * //
-
-Foam::Istream& Foam::operator>>(Istream& is, CMC::BetaGrid&)
-{
-    // Check state of Istream
-    is.check
-    (
-        "Foam::Istream& Foam::operator>>(Foam::Istream&, Foam::CMC::BetaGrid&)"
-    );
-
-    return is;
-}
-
-
-Foam::Ostream& Foam::operator<<(Ostream& os, const CMC::BetaGrid&)
-{
-    // Check state of Ostream
-    os.check
-    (
-        "Foam::Ostream& Foam::operator<<(Foam::Ostream&, "
-        "const Foam::CMC::BetaGrid&)"
-    );
-
-    return os;
-}
-
-
-// ************************************************************************* //

solvers_post/BetaIntegrator/BetaIntegrator.C

@@ -1,161 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  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 <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-#include "BetaIntegrator.H"
-
-
-#include "error.H"
-#include "SubList.H"
-#include "SubField.H"
-#include "interpolateXY.H"
-
-// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
-
-// const scalar Foam::CMC::BetaIntegrator::staticData();
-
-
-// * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
-
-//extended Simpson's rule (Numerical recipes, 2nd Ed. p.128)
-//for equally spaced and even N intervals (or odd N+1 points) 
-Foam::scalar Foam::CMC::BetaIntegrator::simps
-(const scalar xl, const scalar xh, const label n, const UList<scalar>& fx)
-    const
-{
-    scalar evensum(0.0), oddsum(0.0), sum(0.0);
-
-    scalar h = (xh - xl)/scalar(n);
-
-    for(label i=0 ; i<fx.size() ; i++)
-    {
-        if(i%2 == 0)
-        {
-            evensum += fx[i];
-        }
-        else
-        {
-            oddsum += fx[i];
-        }
-    }
-
-    sum = -1.0*fx.first() + 2.0*evensum + 4.0*oddsum - 1.0*fx.last();
-
-    return sum*(h/3.0);
-}
-
-
-Foam::scalar Foam::CMC::BetaIntegrator::sumSimps(const UList<scalar>& f) const
-{
-    scalar total = 0.0;
-
-    const labelList &N(bg_.N());
-    const scalarField &etaCut(bg_.etaCut());
-
-    forAll(N, i)
-    {
-        const labelList::subList prev(N, i);
-        const label baseI = sum(prev) + 1;
-
-        const scalarField::subField subInterval(f, N[i]+1, baseI);
-
-        total += simps(etaCut[i], etaCut[i+1], N[i], subInterval);
-    }
-
-    return total;
-}
-
-// * * * * * * * * * * * * Protected Member Functions  * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
-
-Foam::CMC::BetaIntegrator::BetaIntegrator(const BetaFunction& bf, const BetaGrid &bg)
-:
-    etaSpace_(bg.etaSpace()),
-    pdfNum_(bg.etaSpace().size(), 0.0),
-    pdfDen_(0.0),
-    bf_(bf),
-    bg_(bg)
-{
-    const scalarField::subField domain(etaSpace_, etaSpace_.size()-2, 1);
-    scalarField::subField range(pdfNum_, etaSpace_.size()-2, 1);
-    range = bf_.etaFunc(domain);
-    pdfDen_ = sumSimps(pdfNum_)
-        + Foam::pow(etaSpace_[1], bf_.alpha())/(bf_.alpha() + SMALL)
-        + Foam::pow(etaSpace_[1], bf_.beta())/(bf_.beta() + SMALL);
-}
-
-
-// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
-
-Foam::CMC::BetaIntegrator::~BetaIntegrator()
-{}
-
-
-// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
-
-// beta-pdf weighted integration for given mf, mfVar and f
-Foam::scalar Foam::CMC::BetaIntegrator::betaIntegrate(const scalarField& f) const
-{
-    scalar result = 0.0;
-
-    if(bf_.fdelta())
-    {
-        result = interpolateXY(bf_.mf(), etaSpace_, f);
-    }
-    else
-    {
-        scalar num = sumSimps (f * pdfNum_)
-            + f.first()*Foam::pow(etaSpace_[1], bf_.alpha())/(bf_.alpha() + SMALL)
-            + f.last()*Foam::pow(etaSpace_[1], bf_.beta())/(bf_.beta() + SMALL);
-        result = num/pdfDen_;
-    }
-
-    return result;
-}
-
-// * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * * //
-
-void Foam::CMC::BetaIntegrator::operator=(const BetaIntegrator& rhs)
-{
-    // Check for assignment to self
-    if (this == &rhs)
-    {
-        FatalErrorInFunction
-            << "Attempted assignment to self"
-            << abort(FatalError);
-    }
-}
-
-// * * * * * * * * * * * * * * Friend Functions  * * * * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * * //
-
-
-// ************************************************************************* //

solvers_post/BetaIntegrator/BetaIntegrator.H

@@ -1,166 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  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 <http://www.gnu.org/licenses/>.
-
-Class
-    Foam::CMC::BetaIntegrator
-
-Description
-
-SourceFiles
-    BetaIntegratorI.H
-    BetaIntegrator.C
-    BetaIntegratorIO.C
-
-\*---------------------------------------------------------------------------*/
-
-#ifndef BetaIntegrator_H
-#define BetaIntegrator_H
-
-#include "scalarField.H"
-#include "BetaFunction.H"
-#include "BetaGrid.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-
-// Forward declaration of classes
-class Istream;
-class Ostream;
-
-// Forward declaration of friend functions and operators
-namespace CMC
-{
-    class BetaIntegrator;
-}
-Istream& operator>>(Istream&, CMC::BetaIntegrator&);
-Ostream& operator<<(Ostream&, const CMC::BetaIntegrator&);
-
-
-/*---------------------------------------------------------------------------*\
-                         Class BetaIntegrator Declaration
-\*---------------------------------------------------------------------------*/
-
-namespace CMC
-{
-
-class BetaIntegrator
-{
-    // Private data
-
-        //- Quadrature points
-        const scalarField &etaSpace_;
-
-        //- Beta pdf numerators
-        scalarField pdfNum_;
-
-        //- Beta pdf denominator
-        scalar pdfDen_;
-
-        //- Description of data_
-        BetaFunction bf_;
-
-        //- Description of data_
-        const BetaGrid &bg_;
-
-
-    // Private Member Functions
-
-        //- Calculate integral using Simpson's rule
-        scalar simps(const scalar xl, const scalar xh, const label n, const UList<scalar>& fx) const;
-
-        //- Calculate integral using Simpson's rule
-        scalar sumSimps(const UList<scalar>& f) const;
-
-        //- Disallow default bitwise copy construct
-        BetaIntegrator(const BetaIntegrator&);
-
-        //- Disallow default bitwise assignment
-        void operator=(const BetaIntegrator&);
-
-
-public:
-
-    // Static data members
-
-        //- Static data staticData
-        // static const scalar staticData;
-
-
-    // Constructors
-
-        //- Construct from components
-        BetaIntegrator(const BetaFunction& bf, const BetaGrid &bg);
-
-        //- Construct as copy
-        // BetaIntegrator(const BetaIntegrator&);
-
-
-    //- Destructor
-    ~BetaIntegrator();
-
-
-    // Member Functions
-
-        scalar betaIntegrate(const scalarField& f) const;
-
-        // Access
-
-        // Check
-
-        // Edit
-
-        // Write
-
-
-    // Member Operators
-
-        // void operator=(const BetaIntegrator&);
-
-
-    // Friend Functions
-
-    // Friend Operators
-
-    // IOstream Operators
-
-        friend Istream& ::Foam::operator>>(Istream&, BetaIntegrator&);
-        friend Ostream& ::Foam::operator<<(Ostream&, const BetaIntegrator&);
-};
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace CMC
-} // End namespace Foam
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#include "BetaIntegratorI.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#endif
-
-// ************************************************************************* //

solvers_post/FlameStructure/FlameStructure.C

@@ -1,265 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  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 <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-#include "FlameStructure.H"
-
-#include "IFstream.H"
-
-// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
-
-// const Foam::CMC::dataType Foam::CMC::FlameStructure::staticData();
-
-
-// * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * Protected Member Functions  * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
-
-Foam::CMC::FlameStructure::FlameStructure(word fileName)
-:
-    names_(),
-    Y_(),
-    W_(),
-    T_(),
-    Q_(),
-    h_(),
-    rho_(),
-    Rgas_()
-{
-    IFstream slfmFile (fileName);
-
-    word garbage;
-    word yname;
-
-    // number of eta grids and species
-    label nEta(0);
-    label nY(0);
-
-    slfmFile.getLine(garbage);      //INPUT FILE FOR THE SR-CMC ROUTINES
-    slfmFile.getLine(garbage);      //No.GRID POINTS No. SPECIES
-    slfmFile >> nEta >> nY;
-
-    Y_.setSize(nY);
-    W_.setSize(nY);
-
-    eta_.setSize(nEta);
-
-    T_.setSize(nEta);
-    Q_.setSize(nEta);
-    h_.setSize(nEta);
-    rho_.setSize(nEta);
-    Rgas_.setSize(nEta);
-
-    slfmFile.getLine(garbage);      //blank line (trailing white spaces)
-    slfmFile.getLine(garbage);      //PRESSURE(ATM)
-    slfmFile.getLine(garbage);      //1.0
-    slfmFile.getLine(garbage);      //MIXTURE FRACTION GRID
-
-    scalarField YTemp(nEta,0);
-    scalarField WTemp(nEta,0);
-
-
-    for(label j=0 ; j<nEta ; j++)
-    {
-        slfmFile >> eta_[j];
-    }
-
-    slfmFile.getLine(garbage);      //blank line (trailing white spaces)
-    slfmFile.getLine(garbage);      //INITIAL CONDITIONS
-
-
-    for(label i=0; i<nY ; i++)
-    {
-        //species name
-        slfmFile >> yname;
-        names_.append(yname);
-
-        //species mass fractions
-        for(label j=0 ; j<nEta ; j++)
-        {
-            slfmFile >> YTemp[j];
-        }
-        Y_.set(i, new scalarField(YTemp));
-
-        //species production rates
-        for(label j=0 ; j<nEta ; j++)
-        {
-            slfmFile >> WTemp[j];
-        }
-        W_.set(i, new scalarField(WTemp));
-
-        slfmFile.getLine(garbage);      //blank line (trailing white spaces)
-    }
-
-
-    //read temperature
-    slfmFile >> yname;
-
-    for(label j=0 ; j<nEta ; j++)
-    {
-        slfmFile >> T_[j];
-        T_[j] = max(0.0, T_[j]);
-    }
-
-    for(label j=0 ; j<nEta ; j++)
-    {
-        slfmFile >> Q_[j];
-    }
-
-
-    /* 
-    // calcuate rather than read
-
-    //read density
-    slfmFile >> yname;
-
-    for(label j=0 ; j<nEta ; j++)
-    {
-        slfmFile >> rho_[j];
-    }
-
-
-    //read enthalpy
-    slfmFile >> yname;
-
-    for(label j=0 ; j<nEta ; j++)
-    {
-        slfmFile >> Rgas_[j];
-    }
-    */
-
-
-
-
-
-
-
-
-
-
-
-
-
-    /*
-    Interpolation(etaValue, eta, 1, nY+1, yi_temp, yi);
-
-    for(label j = 0 ; j<=etamax ; j++)
-    {
-        for(label i = 0 ; i<Ysize ; i++)
-        {
-            QiCMC[j*Ysize+i] = yi[j*(nY+1)+i];
-        }
-        TCMC[j] = yi[j*(nY+1)+nY];
-    }
-
-    for(label k=1 ; k<nf ; k++)
-    {
-        for(label j = 0 ; j<=etamax ; j++)
-        {
-            for(label i = 0 ; i<Ysize ; i++)
-            {
-                QiCMC[(k*(etamax+1)+j)*Ysize+i] = QiCMC[j*Ysize+i];
-            }
-            TCMC[j+k*(etamax+1)] = TCMC[j];
-        }
-    }
-
-
-
-    for(label j=0 ; j <= etamax ; j++) //get conditional enthalpy 
-    {
-        for(label k = 0 ; k<nf ; k++)
-        {
-            scalar Tin = TCMC[ j+k*(etamax+1) ];
-            scalarField QiCMCin(Ysize, 0);
-
-            for(label i=0 ; i<Ysize ; i++)
-            {
-                QiCMCin[i] = QiCMC[(k*(etamax+1)+j)*Ysize + i];
-            }
-
-            rhoCMC[j+k*(etamax+1)] = chemistry->calculateRHOCMC(QiCMCin , Tin, Pin);
-
-            QhCMC[j+k*(etamax+1)] = chemistry->calculateHCMC(QiCMCin , Tin, rhoCMC[ j+k*(etamax+1)], Pin);
-
-        }
-    }
-    */
-}
-
-
-Foam::CMC::FlameStructure::FlameStructure(const FlameStructure&)
-// :
-    // baseClassName(),
-    // data_()
-{}
-
-
-// * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * * //
-
-/*
-Foam::autoPtr<Foam::CMC::FlameStructure>
-Foam::CMC::FlameStructure::New()
-{
-    return autoPtr<FlameStructure>(new FlameStructure);
-}
-*/
-
-
-// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
-
-Foam::CMC::FlameStructure::~FlameStructure()
-{}
-
-
-// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * * //
-
-void Foam::CMC::FlameStructure::operator=(const FlameStructure& rhs)
-{
-    // Check for assignment to self
-    if (this == &rhs)
-    {
-        FatalErrorInFunction
-            << "Attempted assignment to self"
-            << abort(FatalError);
-    }
-}
-
-// * * * * * * * * * * * * * * Friend Functions  * * * * * * * * * * * * * * //
-
-
-// * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * * //
-
-
-// ************************************************************************* //

solvers_post/FlameStructure/FlameStructure.H

@@ -1,242 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  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 <http://www.gnu.org/licenses/>.
-
-Class
-    Foam::CMC::FlameStructure
-
-Description
-
-SourceFiles
-    FlameStructureI.H
-    FlameStructure.C
-    FlameStructureIO.C
-
-\*---------------------------------------------------------------------------*/
-
-#ifndef FlameStructure_H
-#define FlameStructure_H
-
-#include "scalarField.H"
-#include "hashedWordList.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-
-// Forward declaration of classes
-class Istream;
-class Ostream;
-
-namespace CMC
-{
-class FlameStructure;
-}
-
-// Forward declaration of friend functions and operators
-Istream& operator>>(Istream&, CMC::FlameStructure&);
-Ostream& operator<<(Ostream&, const CMC::FlameStructure&);
-
-namespace CMC
-{
-
-/*---------------------------------------------------------------------------*\
-                         Class FlameStructure Declaration
-\*---------------------------------------------------------------------------*/
-
-class FlameStructure
-{
-    // Private data
-
-        //- Mixture fraction grid
-        scalarField eta_;
-
-        //- Species name table
-        hashedWordList names_;
-
-        //- Species mass fraction
-        PtrList<scalarField> Y_;
-        // scalarFieldArray2d Y_;
-
-        //- Species production rate
-        PtrList<scalarField> W_;
-        //scalarFieldArray2d W_;
-
-        //- Temperature
-        scalarField T_;
-
-        //- Heat source
-        scalarField Q_;
-
-        //- Enthalpy
-        scalarField h_;
-
-        //- Density
-        scalarField rho_;
-
-        //- Specific gas constant for mean W
-        scalarField Rgas_;
-
-
-    // Private Member Functions
-
-        //- Disallow Construct null
-        FlameStructure();
-
-        //- Disallow default bitwise copy construct
-        FlameStructure(const FlameStructure&);
-
-        //- Disallow default bitwise assignment
-        void operator=(const FlameStructure&);
-
-        //- limit value of alpha and beta
-        void limitAB();
-
-
-public:
-
-    // Static data members
-
-        //- Static data staticData
-        // static const dataType staticData;
-
-
-    // Constructors
-
-        //- Construct from components
-        FlameStructure(word fileName);
-
-        //- Construct from Istream
-        FlameStructure(Istream&);
-
-        //- Construct as copy
-        // FlameStructure(const FlameStructure&);
-
-
-    // Selectors
-
-        //- Select null constructed
-        // static autoPtr<FlameStructure> New();
-
-
-    //- Destructor
-    ~FlameStructure();
-
-
-    // Member Functions
-
-        // Access
-
-            //- Select null constructed
-            scalarField& eta() {return eta_;}
-
-            //- Select null constructed
-            hashedWordList& names() {return names_;}
-
-            //- Select null constructed
-            PtrList<scalarField>& Y() {return Y_;}
-
-            //- Select null constructed
-            PtrList<scalarField>& W() {return W_;}
-
-            //- Select null constructed
-            scalarField& T() {return T_;}
-
-            //- Select null constructed
-            scalarField& Q() {return Q_;}
-
-            //- Select null constructed
-            scalarField& rho() {return rho_;}
-
-            //- Select null constructed
-            scalarField& h() {return h_;}
-
-            //- Select null constructed
-            scalarField& Rgas() {return Rgas_;}
-
-            //- Select null constructed
-            const scalarField& eta() const {return eta_;}
-
-            //- Select null constructed
-            const hashedWordList& names() const {return names_;}
-
-            //- Select null constructed
-            const PtrList<scalarField>& Y() const {return Y_;}
-
-            //- Select null constructed
-            const PtrList<scalarField>& W() const {return W_;}
-
-            //- Select null constructed
-            const scalarField& T() const {return T_;}
-
-            //- Select null constructed
-            const scalarField& Q() const {return Q_;}
-
-            //- Select null constructed
-            const scalarField& rho() const {return rho_;}
-
-            //- Select null constructed
-            const scalarField& h() const {return h_;}
-
-            //- Select null constructed
-            const scalarField& Rgas() const {return Rgas_;}
-    /*
-    */
-
-        // Check
-
-        // Edit
-
-        // Write
-
-
-    // Member Operators
-
-        // void operator=(const FlameStructure&);
-
-
-    // Friend Functions
-
-    // Friend Operators
-
-    // IOstream Operators
-
-        friend Istream& ::Foam::operator>>(Istream&, FlameStructure&);
-        friend Ostream& ::Foam::operator<<(Ostream&, const FlameStructure&);
-};
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace CMC
-} // End namespace Foam
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#include "FlameStructureI.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#endif
-
-// ************************************************************************* //

solvers_post/FlameStructure/FlameStructureIO.C

@@ -1,79 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  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 <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-#include "FlameStructure.H"
-#include "IOstreams.H"
-#include "error.H"
-
-// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
-
-Foam::CMC::FlameStructure::FlameStructure(Istream& is)
-    /*
-:
-    base1(is),
-    base2(is),
-    member1(is),
-    member2(is)
-    */
-{
-    // Check state of Istream
-    is.check("Foam::CMC::FlameStructure::FlameStructure(Foam::Istream&)");
-}
-
-
-// * * * * * * * * * * * * * * * IOstream Operators  * * * * * * * * * * * * //
-
-Foam::Istream& Foam::operator>>(Istream& is, CMC::FlameStructure&)
-{
-    // Check state of Istream
-    is.check
-    (
-        "Foam::Istream& Foam::operator>>(Foam::Istream&, Foam::CMC::FlameStructure&)"
-    );
-
-    return is;
-}
-
-
-Foam::Ostream& Foam::operator<<(Ostream& os, const CMC::FlameStructure& fs)
-{
-    // Check state of Ostream
-    os.check
-    (
-        "Foam::Ostream& Foam::operator<<(Foam::Ostream&, "
-        "const Foam::CMC::FlameStructure&)"
-    );
-
-    os << fs.names()
-        << fs.eta()
-        << fs.Y()
-        << fs.T()
-        << endl;
-
-    return os;
-}
-
-
-// ************************************************************************* //

solvers_post/LagrangianCMCFoam/CMCintegration.H

@@ -1,37 +1,50 @@
-
+                
         forAll(rho, celli)
         {
-            if(F_total[celli] > 0)
-            {
-                scalarField Ytemp(Ysize, 0.0);
-
-                for(label i = 0 ; i < Ysize ; i++) //mass fraction integration over an eta-space for species i
+                if(F_total[celli] > 0)
                 {
-                    for(label k=0 ; k<=group ; k++)
-                    {
-                        for(label j=0 ; j<=etamax ; j++)
+			scalarField Ytemp(Ysize, 0.0);
+			
+                        for(label i=0 ; i < Ysize ; i++) //reaction rate integration over an eta-space for species i
                         {
-                            pdf[j] = Peta[j][celli];
-                            f[j] = QiCMC[(k*(etamax+1)+j)*Ysize + i];
-                        }
-                        Ytemp[i] += Ffrac[k][celli]*integration(deltaftn[celli], MFcut, Neta, pdf, f);
-                    }
-                }
+        			for(label k=0 ; k<=group ; k++)
+        			{
+					for(label j=0 ; j<= etamax ; j++)
+					{
+						pdf[j] = Peta[j][celli];
+						f[j] = QiCMC[(k*(etamax+1)+j)*Ysize + i];
+					}
+					Ytemp[i] += Ffrac[k][celli]*integration(deltaftn[celli], MFcut, Neta, pdf, f);
+                		}
+        		}
 
-                Switch ignoreOFIncrease = true;
-                Switch oxygenIncrease =  Ytemp[o2index] > Y[o2index][celli];
-                Switch fuelIncrease = Ytemp[fuelindex] > Y[fuelindex][celli];
-                if( ignoreOFIncrease && (oxygenIncrease || fuelIncrease))
-                {
-                    // Do not update Yi
+                        if( Ytemp[o2index] > Y[o2index][celli] || Ytemp[fuelindex] > Y[fuelindex][celli])
+                        {
+                                SUMSh[celli] = 0;
+                                //keep Yi value
+                        }
+                        else
+                        {
+                                scalarField ReactionRate(Ysize, 0.0);
+                                scalar tShtemp(0);
+
+                                for(label i = 0 ; i<Ysize ; i++)
+                                {
+                                        ReactionRate[i] = rho[celli] * (Ytemp[i] - Y[i][celli]) / runTime.deltaT().value();
+                                        tShtemp -= chemistry->calculateShCMC( ReactionRate , i);
+                                }
+                                SUMSh[celli] = tShtemp;
+
+                                for(label i = 0 ; i <Ysize ; i++)
+                                {
+                                        Y[i][celli] = Ytemp[i];
+                                }
+
+                        }
                 }
                 else
                 {
-                    for(label i = 0 ; i < Ysize ; i++)
-                    {
-                        Y[i][celli] = Ytemp[i];
-                    }
+                        SUMSh[celli] = 0;
                 }
-            }
         }
-
+ 

solvers_post/LagrangianCMCFoam/EEqn.H

@@ -6,7 +6,7 @@
         fvm::ddt(rho, he) + mvConvection->fvmDiv(phi, he)
       + fvc::ddt(rho, K) + fvc::div(phi, K)
       + (
-            he.name() == "ea"
+            he.name() == "e"
           ? fvc::div
             (
                 fvc::absolute(phi/fvc::interpolate(rho), U),
@@ -17,7 +17,8 @@
         )
       - fvm::laplacian(1.47*turbulence->mut(), he)
      ==
-        fvOptions(rho, he)
+        SUMSh
+       +fvOptions(rho, he)
     );
 
     EEqn.relax();

solvers_post/LagrangianCMCFoam/Flamegroup.H

@@ -19,20 +19,3 @@
         F_total += F[k];
     }
     
-    forAll(rho , celli)
-    {
-        if(F_total[celli] > 0)
-        {
-            for(label k=0 ; k<=group ; k++)
-            {
-                Ffrac[k][celli] = F[k][celli]/F_total[celli];
-            }
-        }
-        else if(F_total[celli] <= 0)
-        {
-            for(label k=0 ; k<=group ; k++)
-            {
-                Ffrac[k][celli] = 0.0;
-            }
-        }
-    }

solvers_post/LagrangianCMCFoam/InjNewFG.H

@@ -1,7 +1,7 @@
 if(fvc::domainIntegrate(F[group]*rho).value() >= fuelperFG) 
 {
 	NewFGindex = NewFGindex + 1;
-	if(group < (nf-1))
+	if(group < (nf-1.5))
 	{
 		group++;
         	forAll(mf.boundaryFieldRef(), patchi)
@@ -17,7 +17,7 @@ if(fvc::domainIntegrate(F[group]*rho).value() >= fuelperFG)
         	}
 		eventRecordFile<<runTime.value()<<tab<<"New flame group is injected"<<tab<<group<<endl;		
 	}
-	else if(group >= (nf-1))
+	else if(group > (nf-1.5))
 	{
 		//retain group, BC
 		forAll(rho, celli)

solvers_post/LagrangianCMCFoam/LagrangianCMCFoam.C

@@ -39,7 +39,7 @@ Contact
 
 #include "fvCFD.H"
 #include "turbulentFluidThermoModel.H"
-#include "psiChemistryModel.H"
+#include "psiChemistryCombustion.H"
 #include "multivariateScheme.H"
 #include "pisoControl.H"
 #include "fvOptions.H"
@@ -48,17 +48,6 @@ Contact
 #include "IFstream.H"
 #include "OFstream.H"
 #include "Math.H" //Mathematical functions for CMC (AMC, gammaln, TDMA)
-
-#include "LinearMesh.H"
-
-#include "FlameStructure.H"
-
-#include "AMC.H"
-
-#include "BetaFunction.H"
-#include "BetaGrid.H"
-#include "BetaIntegrator.H"
-
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 int main(int argc, char *argv[])
@@ -84,8 +73,6 @@ int main(int argc, char *argv[])
     #include "startSummary.H" //Make logSummary file for CMC
     #include "readCMCProperties.H" //Read and set fields for CMC calculation
 
-    #include "createLinearMesh.H"
-
     // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
     Info<< "\nStarting time loop\n" << endl;
@@ -100,8 +87,6 @@ int main(int argc, char *argv[])
 
         Info<< "Time = " << runTime.timeName() << nl << endl;
 
-        #include "QiEqn.H"
-
         #include "rhoEqn.H"
         #include "UEqn.H"
         if(init_start_CMC == true)
@@ -111,6 +96,7 @@ int main(int argc, char *argv[])
                 init_start_CMC = false;
         }
         #include "CMCequation.H" //Solve CMC equations
+        #include "YEqn.H"
         #include "CMCintegration.H" //Conditional field integration and get reaction rate
         #include "EEqn.H"
 

solvers_post/LagrangianCMCFoam/LinearMesh.H

@@ -1,337 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  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 <http://www.gnu.org/licenses/>.
-
-Class
-    Foam::LinearMesh
-
-Description
-
-SourceFiles
-    LinearMeshI.H
-    LinearMesh.C
-    LinearMeshIO.C
-
-\*---------------------------------------------------------------------------*/
-
-#ifndef LinearMesh_H
-#define LinearMesh_H
-
-#include "fvMesh.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-
-// Forward declaration of classes
-class Istream;
-class Ostream;
-
-// Forward declaration of friend functions and operators
-class LinearMesh;
-Istream& operator>>(Istream&, LinearMesh&);
-Ostream& operator<<(Ostream&, const LinearMesh&);
-
-
-/*---------------------------------------------------------------------------*\
-                         Class LinearMesh Declaration
-\*---------------------------------------------------------------------------*/
-
-class LinearMesh
-{
-    // Private data
-
-        //- Lower end value of 1-D domain
-        scalar x0_;
-
-        //- Upper end value of 1-D domain
-        scalar x1_;
-
-        //- Total number of cells
-        label nCells_;
-
-        //- Total number of cells
-        pointField points_;
-
-        //- Total number of cells
-        faceList faces_;
-
-        //- Total number of cells
-        labelList owner_;
-
-        //- Total number of cells
-        labelList neighbour_;
-
-        //- Total number of cells
-        label nPatches_;
-
-        //- Total number of cells
-        wordList patchNames_;
-
-        //- Total number of cells
-        labelList patchStartIndices_;
-
-        //- Total number of cells
-        labelList patchSizes_;
-
-
-    // Private Member Functions
-
-        //- Disallow Construct null
-        LinearMesh();
-
-        //- Disallow default bitwise copy construct
-        LinearMesh(const LinearMesh&);
-
-        //- Disallow default bitwise assignment
-        void operator=(const LinearMesh&);
-
-        //- Generate points and faces
-        void init(const scalarField &etaValue)
-        {
-            for(label i = 0; i < (nCells_+1); i++)
-            {
-                points_[4*i+0] = vector(etaValue[i], 0.0, 0.0);
-                points_[4*i+1] = vector(etaValue[i], 0.1, 0.0);
-                points_[4*i+2] = vector(etaValue[i], 0.1, 0.1);
-                points_[4*i+3] = vector(etaValue[i], 0.0, 0.1);
-            }
-
-            const label nBoundaryFaces = 2;
-            const label nInternalFaces = nCells_-1;
-            const label nEmptyFaces = 4*nCells_;
-
-            for(label i = 0; i < nInternalFaces; i++)
-            {
-                IStringStream faceFormatStream ("4(4 5 6 7)");
-                labelField faceFormat (faceFormatStream);
-                faces_[i] = face(faceFormat + 4*i);
-            }
-
-            {
-                IStringStream faceFormatStream ("4(0 3 2 1)");
-                faces_[nInternalFaces] = face(faceFormatStream);
-            }
-
-            {
-                IStringStream faceFormatStream ("4(4 5 6 7)");
-                labelField faceFormat (faceFormatStream);
-                faces_[nInternalFaces+1] = face(faceFormat + 4*nInternalFaces);
-            }
-
-            for(label i = 0; i < nCells_; i++)
-            {
-                IStringStream faceFormatStream ("4(4(0 1 5 4)4(1 2 6 5)4(2 3 7 6)4(3 0 4 7))");
-                List<labelField> faceFormats (faceFormatStream);
-
-                label fi = nInternalFaces + nBoundaryFaces + 4*i;
-                faces_[fi+0] = face(faceFormats[0] + 4*i);
-                faces_[fi+1] = face(faceFormats[1] + 4*i);
-                faces_[fi+2] = face(faceFormats[2] + 4*i);
-                faces_[fi+3] = face(faceFormats[3] + 4*i);
-            }
-
-
-            for(label i = 0; i < nInternalFaces; i++)
-            {
-                owner_[i] = i;
-                neighbour_[i] = i+1;
-            }
-
-            {
-                owner_[nInternalFaces] = 0;
-            }
-
-            {
-                owner_[nInternalFaces+1] = nCells_ - 1;
-            }
-
-            for(label i = 0; i < nCells_; i++)
-            {
-                label fi = nInternalFaces + nBoundaryFaces + 4*i;
-                owner_[fi+0] = i;
-                owner_[fi+1] = i;
-                owner_[fi+2] = i;
-                owner_[fi+3] = i;
-            }
-
-
-            nPatches_ = 3;
-
-            IStringStream sInd
-                ( 
-                 word("3(")
-                 + name(nInternalFaces) + " "
-                 + name(nInternalFaces+1) + " "
-                 + name(nInternalFaces+2) + " "
-                 + word(")")
-                );
-            patchStartIndices_ = labelList(sInd);
-
-            IStringStream sSize(word("3(1 1 ")+name(4*nCells_)+word(")"));
-            patchSizes_ = labelList(sSize);
-
-            IStringStream sName("3(lowerEnd upperEnd Sides)");
-            patchNames_ = wordList(sName);
-        }
-
-public:
-
-    // Static data members
-
-        //- Static data staticData
-        // static const dataType staticData;
-
-
-    // Constructors
-
-        //- Construct from components
-        LinearMesh(const scalar x0, const scalar x1, const label nCells)
-            :
-                x0_(x0),
-                x1_(x1),
-                nCells_(nCells),
-                points_(4*(nCells+1)),
-                faces_(5*nCells+1),
-                owner_(faces_.size()),
-                neighbour_(nCells-1)
-        {
-
-            scalarField etaValue(nCells+1, 0.0);
-            scalar delta = (x1 - x0) / scalar(nCells);
-            forAll(etaValue, i)
-            {
-                etaValue[i] = x0 + delta*i;
-            }
-
-
-            init(etaValue);
-
-        }
-
-        //- Construct from components
-        LinearMesh(const scalarField &x)
-            :
-                x0_(x[0]),
-                x1_(x.last()),
-                nCells_(x.size()-1),
-                points_(4*(nCells_+1)),
-                faces_(5*nCells_+1),
-                owner_(faces_.size()),
-                neighbour_(nCells_-1)
-        {
-
-            init(x);
-
-        }
-
-        //- Construct from Istream
-        LinearMesh(Istream&);
-
-        //- Construct as copy
-        // LinearMesh(const LinearMesh&);
-
-
-    // Selectors
-
-        //- Select null constructed
-        // static autoPtr<LinearMesh> New();
-
-
-    //- Destructor
-    ~LinearMesh() {}
-
-
-    // Member Functions
-
-        // Access
-            //- Lower end value of 1-D domain
-            scalar x0() const {return x0_;};
-
-            //- Upper end value of 1-D domain
-            scalar x1() const {return x1_;};
-
-            //- Total number of cells
-            label nCells() const {return nCells_;};
-
-            //- Total number of cells
-            const pointField &points() const {return points_;};
-
-            //- Total number of cells
-            const faceList &faces() const {return faces_;};
-
-            //- Total number of cells
-            const labelList &owner() const {return owner_;};
-
-            //- Total number of cells
-            const labelList &neighbour() const {return neighbour_;};
-
-            //- Total number of cells
-            label nPatches () const {return nPatches_;};
-
-            //- Total number of cells
-            word patchName (const label i) const {return patchNames_[i];};
-
-            //- Total number of cells
-            label patchStartIndex (const label i) const {return patchStartIndices_[i];};
-
-            //- Total number of cells
-            label patchSize (const label i) const {return patchSizes_[i];};
-
-
-        // Check
-
-        // Edit
-
-        // Write
-
-
-    // Member Operators
-
-        // void operator=(const LinearMesh&);
-
-
-    // Friend Functions
-
-    // Friend Operators
-
-    // IOstream Operators
-
-        friend Istream& operator>>(Istream&, LinearMesh&);
-        friend Ostream& operator<<(Ostream&, const LinearMesh&);
-};
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace Foam
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-// #include "LinearMeshI.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#endif
-
-// ************************************************************************* //

solvers_post/LagrangianCMCFoam/Make/files

@@ -1,20 +1,3 @@
-../FlameStructure/FlameStructure.C
-../FlameStructure/FlameStructureIO.C
-
-../AMC/AMC.C
-
-../LinearInterpolator1D/LinearInterpolator1D.C
-../LinearInterpolator1D/LinearInterpolator1DIO.C
-
-../BetaFunction/BetaFunction.C
-../BetaFunction/BetaFunctionIO.C
-
-../BetaGrid/BetaGrid.C
-../BetaGrid/BetaGridIO.C
-
-../BetaIntegrator/BetaIntegrator.C
-../BetaIntegrator/BetaIntegratorIO.C
-
 LagrangianCMCFoam.C
 
 EXE = $(FOAM_USER_APPBIN)/LagrangianCMCFoam

solvers_post/LagrangianCMCFoam/Make/options

@@ -1,12 +1,6 @@
 DEV_PATH=../../libs
 
 EXE_INC = \
-    -I../FlameStructure \
-    -I../LinearInterpolator1D \
-    -I../AMC \
-    -I../BetaFunction \
-    -I../BetaGrid \
-    -I../BetaIntegrator \
     -I$(LIB_SRC)/finiteVolume/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude \
     -I$(LIB_SRC)/sampling/lnInclude \
@@ -17,12 +11,13 @@ EXE_INC = \
     -I$(LIB_SRC)/transportModels/compressible/lnInclude \
     -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
     -I$(LIB_SRC)/ODE/lnInclude \
-    -I$(DEV_PATH)/chemistryModel_POSTECH/lnInclude
+    -I$(DEV_PATH)/chemistryModel_POSTECH/lnInclude \
+    -I$(DEV_PATH)/combustionModels_POSTECH/lnInclude
 
 EXE_LIBS = \
     -L$(FOAM_USER_LIBBIN) \
         -lchemistryModel_POSTECH \
-        -lthermoAdd \
+        -lcombustionModels_POSTECH \
     -lfiniteVolume \
     -lfvOptions \
     -lmeshTools \

solvers_post/LagrangianCMCFoam/QCMC_init.H

@@ -1,93 +1,97 @@
-IFstream slfmFile (SLFMinit);
-CMC::FlameStructure test (SLFMinit);
+		IFstream slfmFile
+		(
+		        SLFMinit
+		);
 
-word garbage, yname;
-scalar neta(0), ny(0);
+		word garbage, yname;
+		scalar neta(0), ny(0);
 
-slfmFile.getLine(garbage);//INPUT
-slfmFile.getLine(garbage);//NO.
-slfmFile>>neta>>ny;//number of eta and species
+	        slfmFile.getLine(garbage);//INPUT
+	        slfmFile.getLine(garbage);//NO.
+	        slfmFile>>neta>>ny;//number of eta and species
 
-slfmFile.getLine(garbage);//blank line
-slfmFile.getLine(garbage);//PRESSURE
-slfmFile.getLine(garbage);//1.0
-slfmFile.getLine(garbage);//MIXTURE
+	        slfmFile.getLine(garbage);//blank line
+	        slfmFile.getLine(garbage);//PRESSURE
+	        slfmFile.getLine(garbage);//1.0
+	        slfmFile.getLine(garbage);//MIXTURE
 
-scalarField eta(neta, 0);
-scalarField yi_temp(neta*(ny+1), 0), yi((etamax+1)*(ny+1),0), wi(neta*(ny+1), 0);
+	        scalarField eta(neta, 0);
+	        scalarField yi_temp(neta*(ny+1), 0), yi((etamax+1)*(ny+1),0), wi(neta*(ny+1), 0);
 
-for(label j=0 ; j<neta ; j++)
-{
-    slfmFile>>eta[j];
-}
+	        for(label j=0 ; j<neta ; j++)
+	        {
+	                slfmFile>>eta[j];
+	        }
 
-slfmFile.getLine(garbage); //blank line
-slfmFile.getLine(garbage); //INITIAL
+	        slfmFile.getLine(garbage); //blank line
+	        slfmFile.getLine(garbage); //INITIAL
 
-for(label i=0; i<ny ; i++) //species
-{
-    slfmFile>>yname;
-    for(label j=0 ; j<neta ; j++)
-    {
-        slfmFile>>yi_temp[j*(ny+1)+i];
-    }
+	        for(label i=0; i<ny ; i++) //species
+	        {
+	                slfmFile>>yname;
+			Info<<yname<<endl;
+			for(label j=0 ; j<neta ; j++)
+			{
+				slfmFile>>yi_temp[j*(ny+1)+i];
+			}
 
-    for(label j=0 ; j<neta ; j++)
-    {
-        slfmFile>>wi[j*(ny+1)+i]; //reaction rates
-    }
-    slfmFile.getLine(garbage);//blank line
-}
-//Temperature
+			for(label j=0 ; j<neta ; j++)
+			{
+				slfmFile>>wi[j*(ny+1)+i]; //reaction rates
+			}
+			slfmFile.getLine(garbage);//blank line
+		}
+		//Temperature
 
-slfmFile>>yname;
-for(label j=0 ; j<neta ; j++)
-{
-    slfmFile>>yi_temp[j*(ny+1)+ny];
+		slfmFile>>yname;
+		Info<<yname<<endl;
+		for(label j=0 ; j<etamax+1 ; j++)
+		{
+			slfmFile>>yi_temp[j*(ny+1)+ny];
 
-}
+		}
 
-Interpolation(etaValue, eta, 1, ny+1, yi_temp, yi);
+		Interpolation(etaValue, eta, 1, ny+1, yi_temp, yi);
 
-for(label j = 0 ; j<=etamax ; j++)
-{
-    for(label i = 0 ; i<Ysize ; i++)
-    {
-        QiCMC[j*Ysize+i] = yi[j*(ny+1)+i];
-    }
-    TCMC[j] = yi[j*(ny+1)+ny];
-}
+		for(label j = 0 ; j<=etamax ; j++)
+		{
+			for(label i = 0 ; i<Ysize ; i++)
+			{
+				QiCMC[j*Ysize+i] = yi[j*(ny+1)+i];
+			}
+			TCMC[j] = yi[j*(ny+1)+ny];
+		}
 
-for(label k=1 ; k<nf ; k++)
-{
-    for(label j = 0 ; j<=etamax ; j++)
-    {
-        for(label i = 0 ; i<Ysize ; i++)
-        {
-            QiCMC[(k*(etamax+1)+j)*Ysize+i] = QiCMC[j*Ysize+i];
-        }
-        TCMC[j+k*(etamax+1)] = TCMC[j];
-    }
-}
+		for(label k=1 ; k<nf ; k++)
+		{
+			for(label j = 0 ; j<=etamax ; j++)
+			{
+				for(label i = 0 ; i<Ysize ; i++)
+				{
+					QiCMC[(k*(etamax+1)+j)*Ysize+i] = QiCMC[j*Ysize+i];
+				}
+				TCMC[j+k*(etamax+1)] = TCMC[j];
+			}
+		}
 
 
 
 for(label j=0 ; j <= etamax ; j++) //get conditional enthalpy 
 {
-    for(label k = 0 ; k<nf ; k++)
-    {
-        scalar Tin = TCMC[ j+k*(etamax+1) ];
-        scalarField QiCMCin(Ysize, 0);
-
-        for(label i=0 ; i<Ysize ; i++)
+        for(label k = 0 ; k<nf ; k++)
         {
-            QiCMCin[i] = QiCMC[(k*(etamax+1)+j)*Ysize + i];
+                scalar Tin = TCMC[ j+k*(etamax+1) ];
+                scalarField QiCMCin(Ysize, 0);
+
+                for(label i=0 ; i<Ysize ; i++)
+                {
+                        QiCMCin[i] = QiCMC[(k*(etamax+1)+j)*Ysize + i];
+                }
+
+                rhoCMC[j+k*(etamax+1)] = chemistry->calculateRHOCMC(QiCMCin , Tin, Pin);
+
+                QhCMC[j+k*(etamax+1)] = chemistry->calculateHCMC(QiCMCin , Tin, rhoCMC[ j+k*(etamax+1)], Pin);
+
         }
-
-        rhoCMC[j+k*(etamax+1)] = chemistry->calculateRHOCMC(QiCMCin , Tin, Pin);
-
-        QhCMC[j+k*(etamax+1)] = chemistry->calculateHCMC(QiCMCin , Tin, rhoCMC[ j+k*(etamax+1)], Pin);
-
-    }
 }
 

solvers_post/LagrangianCMCFoam/QiEqn.H

@@ -1,29 +0,0 @@
-{
-    fvScalarMatrix testSource(testEtaMesh, testEtaMesh.dimensions()*dimVol/dimTime);
-
-    testSource.source() = 1.0;
-    
-    fvScalarMatrix QiEqn
-    (
-        fvm::ddt(testEtaMesh)
-//      + mvConvection->fvmDiv(phi, he)
-//      + fvc::ddt(rho, K) + fvc::div(phi, K)
-//      + (
-//            he.name() == "ea"
-//          ? fvc::div
-//            (
-//                fvc::absolute(phi/fvc::interpolate(rho), U),
-//                p,
-//                "div(phiv,p)"
-//            )
-//          : -dpdt
-//        )
-      - fvm::laplacian(DQi, testEtaMesh)
-      ==
-      testSource
-    );
-
-    QiEqn.relax();
-
-    QiEqn.solve();
-}

solvers_post/LagrangianCMCFoam/createFields.H

@@ -1,10 +1,16 @@
 Info<< "Creating reaction model\n" << endl;
 
-autoPtr<psiChemistryModel> chemistry(psiChemistryModel::New(mesh));
+combustionModels::psiCombustionModel* pRxnCMC = (combustionModels::psiCombustionModel::New(mesh)).ptr();
 
+autoPtr<combustionModels::psiChemistryCombustion> reaction
+(
+    dynamic_cast <combustionModels::psiChemistryCombustion*> (pRxnCMC)
+);
+
+autoPtr<psiChemistryModel> chemistry = reaction->chem();
 psiReactionThermo& thermo = chemistry->thermo();
 
-thermo.validate(args.executable(), "ha");
+thermo.validate(args.executable(), "h", "e");
 
 basicMultiComponentMixture& composition = thermo.composition();
 PtrList<volScalarField>& Y = composition.Y();
@@ -77,6 +83,9 @@ autoPtr<compressible::turbulenceModel> turbulence
     )
 );
 
+// Set the turbulence into the reaction model
+reaction->setTurbulence(turbulence());
+
 
 Info<< "Creating field dpdt\n" << endl;
 volScalarField dpdt

solvers_post/LagrangianCMCFoam/createLinearMesh.H

@@ -1,145 +0,0 @@
-Info<< "Constructing single cell mesh" << nl << endl;
-
-LinearMesh lMesh(0, 1, 100);
-
-pointField CMCPoints(lMesh.points());
-faceList CMCFaces(lMesh.faces());
-labelList CMCOwners(lMesh.owner());
-labelList CMCNeighbors(lMesh.neighbour());
-
-fvMesh mfMesh
-(
-    IOobject
-    (
-        "mfSpace",
-        runTime.timeName(),
-        runTime,
-        IOobject::READ_IF_PRESENT
-    ),
-    xferMove<Field<vector>>(CMCPoints),
-    CMCFaces.xfer(),
-    CMCOwners.xfer(),
-    CMCNeighbors.xfer()
-);
-
-if (mfMesh.boundaryMesh().empty())
-{
-    List<polyPatch*> mfPatches(3);
-
-    mfPatches[0] = new polyPatch
-    (
-        lMesh.patchName(0),
-        lMesh.patchSize(0),
-        lMesh.patchStartIndex(0),
-        0,
-        mfMesh.boundaryMesh(),
-        polyPatch::typeName
-    );
-
-    mfPatches[1] = new polyPatch
-    (
-        lMesh.patchName(1),
-        lMesh.patchSize(1),
-        lMesh.patchStartIndex(1),
-        1,
-        mfMesh.boundaryMesh(),
-        polyPatch::typeName
-    );
-
-    mfPatches[2] = new emptyPolyPatch
-    (
-        lMesh.patchName(2),
-        lMesh.patchSize(2),
-        lMesh.patchStartIndex(2),
-        2,
-        mfMesh.boundaryMesh(),
-        emptyPolyPatch::typeName
-    );
-
-    mfMesh.addFvPatches(mfPatches);
-}
-
-IStringStream tempStream("3(fixedValue fixedValue empty)");
-wordList testEtaPatchTypeNames(tempStream);
-
-volScalarField testEtaMesh
-(
-    IOobject
-    (
-        "testEtaMesh",
-        runTime.timeName(),
-        mfMesh,
-        IOobject::NO_READ,
-        IOobject::AUTO_WRITE
-    ),
-    mfMesh,
-    dimensionedScalar("testEtaMesh", dimless, 1.0),
-    testEtaPatchTypeNames
-);
-
-const volVectorField &cEtaMesh(mfMesh.C());
-forAll(cEtaMesh, celli)
-{
-    const scalar x = cEtaMesh[celli].x();
-    if (x > 0.5)
-    {
-        testEtaMesh[celli] = 1.0;
-    }
-    else
-    {
-        testEtaMesh[celli] = 0.0;
-    }
-    
-}
-
-volScalarField DQi
-(
-    IOobject
-    (
-        "DQi",
-        runTime.timeName(),
-        mfMesh,
-        IOobject::NO_READ,
-        IOobject::NO_WRITE
-    ),
-    mfMesh,
-    dimensionedScalar("DQi", dimArea/dimTime, 0.01),
-    testEtaPatchTypeNames
-);
-
-PtrList<PtrList<volScalarField> > FStruct (nf);
-for(label k = 0 ; k < nf ; k++)
-{
-    FStruct.set (k, new PtrList<volScalarField> (Ysize+3));
-
-    PtrList<volScalarField> &Fk(FStruct[k]);
-    Info << Fk.size() << endl;
-    for (label i = 0; i < Ysize + 3; i++)
-    {
-        Fk.set (i,
-            new volScalarField 
-            (
-                IOobject
-                (
-                    "test",
-                    runTime.timeName(),
-                    mfMesh,
-                    IOobject::NO_READ,
-                    IOobject::AUTO_WRITE
-                ),
-                mfMesh,
-                dimensionedScalar("test", dimless, 0.0),
-                testEtaPatchTypeNames
-            )
-        );
-
-    }
-}
-
-for(label k = 0 ; k < nf ; k++)
-{
-    for (label i = 0; i < Ysize + 3; i++)
-    {
-        FStruct[k][i].rename("test."+name(k)+"."+name(i));
-    }
-}

solvers_post/LagrangianCMCFoam/readCMCProperties.H

@@ -81,7 +81,7 @@ word fuelName(CMCdict.lookup("fuel_name"));
 word CMCresults_Time(CMCdict.lookup("CMCresults_Time"));
 scalar previousnf(readScalar(CMCdict.lookup("previousnf")));
 scalar fuelperFG(readScalar(CMCdict.lookup("fuelperFG")));
-label nf(readLabel(CMCdict.lookup("nf")));
+scalar nf(readScalar(CMCdict.lookup("nf")));
 
 scalar etamax( readScalar(CMCdict.lookup("etamax")) );
 scalarField etaValue(etamax+1, 0.0);
@@ -95,9 +95,6 @@ scalarField pdf(etamax+1,0.0);
 scalarField f(etamax+1, 0.0);
 Info<<"Set eta space grid"<<tab<<uni_eta<<endl;
 
-
-CMC::BetaGrid bg(CMCdict.subDict("betaGrid"));
-
 //Uniform eta space grid
 if(uni_eta == true)
 {
@@ -152,7 +149,7 @@ else if(uni_eta == false)
 }
 
 scalarField ChemDeltaT( (etamax+1)*nf, 0.0000001); //initial chemical time step;
-label group(0); //flame group number to add evaporation source (default = 0), Flamegroup.H line# 29
+scalar group(0); //flame group number to add evaporation source (default = 0), Flamegroup.H line# 29
 
 
 //Peta field

solvers_post/LagrangianCMCFoam/setSDRCMC.H

@@ -1,34 +1,50 @@
-// Calculate conditional SDR
+//get conditional SDR
+        scalar C2(0); 
+        scalar Airside(0);
+        scalar Fuelside(0);
 
-        scalarField amc(CMC::AMC(bg.etaSpace()));
-
-        forAll(rho, celli)
+        forAll(rho , celli) 
         {
-            scalar mfi = mf[celli];
-            scalar mfVari = mfVar[celli];
-
-            CMC::BetaFunction bf (mfi, mfVari);
-            CMC::BetaIntegrator bi (bf, bg);
-
-            //if there exist delta function at eta = 0 or 1
-            if( bf.delta0() || bf.delta1() )
-            {
-                for(label j = 0 ; j<= etamax ; j++)
+                Airside = Peta[0][celli]*deta[0]; 
+                Fuelside = Peta[etamax][celli]*deta[etamax]; 
+                if( (1-Airside) < 0.000001 || (1-Fuelside) < 0.000001 )//if there exist delta function at eta = 0 or 1
                 {
-                    CSDRCELL[j][celli] = 0;
+                        for(label j = 0 ; j<= etamax ; j++)
+                        {
+                                CSDRCELL[j][celli] = 0;
+                        }
                 }
-            }
-            else
-            {
-                //Eqn (13) from SH Kim C&F paper 120:75-90(20)
-                scalar C1 = bi.betaIntegrate(amc);
-                scalar C2 = SDR[celli] / C1;
-
-                for(label j = 0 ; j <= etamax ; j++)
+                else
                 {
-                    CSDRCELL[j][celli] = C2 * CMC::AMC(etaValue[j]);
+                        scalar C1(0);
+			for(label j=0 ; j<=etamax ;j++)
+			{
+				f[j] = AMC(etaValue[j]);
+				pdf[j] = Peta[j][celli];
+			}
+			C1 = integration(deltaftn[celli], MFcut, Neta, pdf, f);
+                        C2 = SDR[celli] / C1; //Eqn (13) from SH Kim C&F paper 120:75-90(20)
+
+                        for(label j = 0 ; j <= etamax ; j++)
+                        {
+                                scalar eta = etaValue[j];
+                                CSDRCELL[j][celli] = C2 * AMC(eta);
+                        }
                 }
-            }
+		if(F_total[celli] > 0)
+               	{
+			for(label k=0 ; k<=group ; k++)
+			{
+                		Ffrac[k][celli] = F[k][celli]/F_total[celli];
+			}
+              	}
+		else if(F_total[celli] <= 0)
+		{
+                        for(label k=0 ; k<=group ; k++)
+                        {
+                                Ffrac[k][celli] = 0.0;
+                        }			
+		}
         }
 
 //get density weighted averaged value of conditional SDR for each eta and flame-group 

solvers_post/LagrangianCMCSprayFoam/LagrangianCMCSprayFoam.C

@@ -34,7 +34,7 @@ Description
 #include "turbulentFluidThermoModel.H"
 #include "basicSprayCloud.H"
 //#include "psiCombustionModel.H"
-#include "psiChemistryModel.H"
+#include "psiChemistryCombustion.H" //karam
 #include "radiationModel.H"
 #include "SLGThermo.H"
 //#include "pimpleControl.H"

solvers_post/LagrangianCMCSprayFoam/Make/options

@@ -25,11 +25,13 @@ EXE_INC = \
     -I$(LIB_SRC)/ODE/lnInclude \
     -I$(LIB_SRC)/regionModels/regionModel/lnInclude \
     -I$(LIB_SRC)/regionModels/surfaceFilmModels/lnInclude \
-    -I$(DEV_PATH)/chemistryModel_POSTECH/lnInclude
+    -I$(DEV_PATH)/chemistryModel_POSTECH/lnInclude \
+    -I$(DEV_PATH)/combustionModels_POSTECH/lnInclude
 
 EXE_LIBS = \
     -L$(FOAM_USER_LIBBIN) \
         -lchemistryModel_POSTECH \
+        -lcombustionModels_POSTECH \
     -lturbulenceModels \
     -lcompressibleTurbulenceModels \
     -llagrangian \