Merge branch 'master' into cvm

Conflicts:
	src/meshTools/searchableSurface/closedTriSurfaceMesh.C
	src/meshTools/searchableSurface/closedTriSurfaceMesh.H

applications/solvers/combustion/dieselEngineFoam/Make/options

@@ -16,6 +16,7 @@ EXE_INC = \
     -I$(LIB_SRC)/ODE/lnInclude \
     -I$(LIB_SRC)/engine/lnInclude \
     -I$(LIB_SRC)/finiteVolume/lnInclude \
+    -I$(LIB_SRC)/combustionModels/lnInclude
 
 EXE_LIBS = \
     -lengine \
@@ -35,4 +36,5 @@ EXE_LIBS = \
     -llaminarFlameSpeedModels \
     -lchemistryModel \
     -lODE \
-    -ldistributionModels
+    -ldistributionModels \
+    -lcombustionModels

applications/solvers/combustion/dieselEngineFoam/YEqn.H

@@ -10,7 +10,8 @@ tmp<fv::convectionScheme<scalar> > mvConvection
 );
 
 {
-
+    combustion->correct();
+    dQ = combustion->dQ();
     label inertIndex = -1;
     volScalarField Yt(0.0*Y[0]);
 
@@ -20,17 +21,19 @@ tmp<fv::convectionScheme<scalar> > mvConvection
         {
             volScalarField& Yi = Y[i];
 
-            solve
+            fvScalarMatrix YiEqn
             (
                 fvm::ddt(rho, Yi)
               + mvConvection->fvmDiv(phi, Yi)
               - fvm::laplacian(turbulence->muEff(), Yi)
               ==
                 dieselSpray.evaporationSource(i)
-              + kappa*chemistry.RR(i),
+              + combustion->R(Yi)
-                mesh.solver("Yi")
             );
 
+            YiEqn.relax();
+            YiEqn.solve(mesh.solver("Yi"));
+
             Yi.max(0.0);
             Yt += Yi;
         }

applications/solvers/combustion/dieselEngineFoam/createFields.H

@@ -1,10 +1,14 @@
-Info<< nl << "Reading thermophysicalProperties" << endl;
+Info<< "Creating combustion model\n" << endl;
 
-autoPtr<psiChemistryModel> pChemistry
+autoPtr<combustionModels::psiChemistryCombustionModel> combustion
 (
-    psiChemistryModel::New(mesh)
+    combustionModels::psiChemistryCombustionModel::New
+    (
+        mesh
+    )
 );
-psiChemistryModel& chemistry = pChemistry();
+
+psiChemistryModel& chemistry = combustion->pChemistry();
 
 hsCombustionThermo& thermo = chemistry.thermo();
 
@@ -55,20 +59,6 @@ volScalarField& hs = thermo.hs();
 
 #include "compressibleCreatePhi.H"
 
-volScalarField kappa
-(
-    IOobject
-    (
-        "kappa",
-        runTime.timeName(),
-        mesh,
-        IOobject::NO_READ,
-        IOobject::AUTO_WRITE
-    ),
-    mesh,
-    dimensionedScalar("zero", dimless, 0.0)
-);
-
 Info << "Creating turbulence model.\n" << nl;
 autoPtr<compressible::turbulenceModel> turbulence
 (
@@ -81,6 +71,9 @@ autoPtr<compressible::turbulenceModel> turbulence
     )
 );
 
+// Set the turbulence into the combustion model
+combustion->setTurbulence(turbulence());
+
 Info<< "Creating field DpDt\n" << endl;
 volScalarField DpDt
 (
@@ -96,16 +89,16 @@ forAll(Y, i)
 }
 fields.add(hs);
 
-DimensionedField<scalar, volMesh> chemistrySh
+volScalarField dQ
 (
     IOobject
     (
-        "chemistry::Sh",
+        "dQ",
         runTime.timeName(),
         mesh,
         IOobject::NO_READ,
-        IOobject::NO_WRITE
+        IOobject::AUTO_WRITE
     ),
     mesh,
-    dimensionedScalar("chemistrySh", dimEnergy/dimTime/dimVolume, 0.0)
+    dimensionedScalar("dQ", dimEnergy/dimTime, 0.0)
 );

applications/solvers/combustion/dieselEngineFoam/dieselEngineFoam.C

@@ -32,8 +32,8 @@ Description
 #include "fvCFD.H"
 #include "engineTime.H"
 #include "engineMesh.H"
-#include "hCombustionThermo.H"
 #include "turbulenceModel.H"
+#include "psiChemistryCombustionModel.H"
 #include "spray.H"
 #include "psiChemistryModel.H"
 #include "chemistrySolver.H"
@@ -54,7 +54,6 @@ int main(int argc, char *argv[])
     #include "createEngineMesh.H"
     #include "createFields.H"
     #include "readGravitationalAcceleration.H"
-    #include "readCombustionProperties.H"
     #include "createSpray.H"
     #include "initContinuityErrs.H"
     #include "readEngineTimeControls.H"
@@ -82,29 +81,6 @@ int main(int argc, char *argv[])
 
         dieselSpray.evolve();
 
-        Info<< "Solving chemistry" << endl;
-
-        chemistry.solve
-        (
-            runTime.value() - runTime.deltaTValue(),
-            runTime.deltaTValue()
-        );
-
-        // turbulent time scale
-        {
-            volScalarField tk
-            (
-                Cmix*sqrt(turbulence->muEff()/rho/turbulence->epsilon())
-            );
-            volScalarField tc(chemistry.tc());
-
-            // Chalmers PaSR model
-            kappa = (runTime.deltaT() + tc)/(runTime.deltaT() + tc + tk);
-        }
-
-        chemistrySh = kappa*chemistry.Sh()();
-
-
         #include "rhoEqn.H"
 
         for (pimple.start(); pimple.loop(); pimple++)
@@ -130,10 +106,7 @@ int main(int argc, char *argv[])
 
         rho = thermo.rho();
 
-        if (runTime.write())
+        runTime.write();
-        {
-            chemistry.dQ()().write();
-        }
 
         Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
             << "  ClockTime = " << runTime.elapsedClockTime() << " s"

applications/solvers/combustion/dieselEngineFoam/hsEqn.H

@@ -1,14 +1,21 @@
 {
-    solve
+    fvScalarMatrix hsEqn
     (
         fvm::ddt(rho, hs)
       + mvConvection->fvmDiv(phi, hs)
       - fvm::laplacian(turbulence->alphaEff(), hs)
      ==
-       DpDt
+        DpDt
-     + dieselSpray.heatTransferSource()().dimensionedInternalField()
+      + combustion->Sh()
-     + chemistrySh
+      + dieselSpray.heatTransferSource()()
+
     );
 
+    hsEqn.relax();
+    hsEqn.solve();
+
     thermo.correct();
+
+    Info<< "min/max(T) = "
+        << min(T).value() << ", " << max(T).value() << endl;
 }

applications/solvers/combustion/dieselEngineFoam/readCombustionProperties.H

@@ -1,18 +0,0 @@
-Info<< "Reading combustion properties\n" << endl;
-
-IOdictionary combustionProperties
-(
-    IOobject
-    (
-        "combustionProperties",
-        runTime.constant(),
-        mesh,
-        IOobject::MUST_READ_IF_MODIFIED,
-        IOobject::NO_WRITE
-    )
-);
-
-dimensionedScalar Cmix
-(
-    combustionProperties.lookup("Cmix")
-);

applications/solvers/combustion/dieselFoam/Make/options

@@ -14,7 +14,8 @@ EXE_INC = \
     -I$(LIB_SRC)/thermophysicalModels/laminarFlameSpeed/lnInclude \
     -I$(LIB_SRC)/../applications/solvers/reactionThermo/XiFoam \
     -I$(LIB_SRC)/thermophysicalModels/chemistryModel/lnInclude \
-    -I$(LIB_SRC)/ODE/lnInclude
+    -I$(LIB_SRC)/ODE/lnInclude \
+    -I$(LIB_SRC)/combustionModels/lnInclude
 
 EXE_LIBS = \
     -lcompressibleTurbulenceModel \
@@ -33,4 +34,5 @@ EXE_LIBS = \
     -lchemistryModel \
     -lODE \
     -ldistributionModels \
-    -lfiniteVolume
+    -lfiniteVolume \
+    -lcombustionModels

applications/solvers/combustion/dieselFoam/dieselFoam.C

@@ -30,12 +30,11 @@ Description
 \*---------------------------------------------------------------------------*/
 
 #include "fvCFD.H"
-#include "hCombustionThermo.H"
+#include "psiChemistryCombustionModel.H"
 #include "turbulenceModel.H"
-#include "spray.H"
 #include "psiChemistryModel.H"
 #include "chemistrySolver.H"
-
+#include "spray.H"
 #include "multivariateScheme.H"
 #include "IFstream.H"
 #include "OFstream.H"
@@ -52,7 +51,6 @@ int main(int argc, char *argv[])
     #include "createMesh.H"
     #include "createFields.H"
     #include "readGravitationalAcceleration.H"
-    #include "readCombustionProperties.H"
     #include "createSpray.H"
     #include "initContinuityErrs.H"
     #include "readTimeControls.H"
@@ -79,26 +77,6 @@ int main(int argc, char *argv[])
 
         Info<< "Solving chemistry" << endl;
 
-        chemistry.solve
-        (
-            runTime.value() - runTime.deltaTValue(),
-            runTime.deltaTValue()
-        );
-
-        // turbulent time scale
-        {
-            volScalarField tk
-            (
-                Cmix*sqrt(turbulence->muEff()/rho/turbulence->epsilon())
-            );
-            volScalarField tc(chemistry.tc());
-
-            // Chalmers PaSR model
-            kappa = (runTime.deltaT() + tc)/(runTime.deltaT()+tc+tk);
-        }
-
-        chemistrySh = kappa*chemistry.Sh()();
-
         #include "rhoEqn.H"
 
         // --- Pressure-velocity PIMPLE corrector loop
@@ -124,10 +102,7 @@ int main(int argc, char *argv[])
 
         rho = thermo.rho();
 
-        if (runTime.write())
+        runTime.write();
-        {
-            chemistry.dQ()().write();
-        }
 
         Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
             << "  ClockTime = " << runTime.elapsedClockTime() << " s"

applications/solvers/combustion/fireFoam/YhsEqn.H

@@ -20,7 +20,6 @@ tmp<fv::convectionScheme<scalar> > mvConvection
         if (Y[i].name() != inertSpecie)
         {
             volScalarField& Yi = Y[i];
-            fvScalarMatrix R(combustion->R(Yi));
 
             fvScalarMatrix YiEqn
             (
@@ -30,7 +29,7 @@ tmp<fv::convectionScheme<scalar> > mvConvection
               ==
                 parcels.SYi(i, Yi)
               + surfaceFilm.Srho(i)
-              + R
+              + combustion->R(Yi)
             );
 
             YiEqn.relax();
@@ -55,7 +54,7 @@ tmp<fv::convectionScheme<scalar> > mvConvection
       - fvm::laplacian(turbulence->alphaEff(), hs)
      ==
         DpDt
-      + dQ
+      + combustion->Sh()
       + radiation->Shs(thermo)
       + parcels.Sh(hs)
       + surfaceFilm.Sh()
@@ -66,5 +65,6 @@ tmp<fv::convectionScheme<scalar> > mvConvection
 
     thermo.correct();
 
-    Info<< "min/max(T) = " << min(T).value() << ", " << max(T).value() << endl;
+    Info<< "min/max(T) = "
+        << min(T).value() << ", " << max(T).value() << endl;
 }

applications/solvers/combustion/fireFoam/createFields.H

@@ -1,10 +1,16 @@
+    Info<< "Creating combustion model\n" << endl;
+
+    autoPtr<combustionModels::psiCombustionModel> combustion
+    (
+        combustionModels::psiCombustionModel::New
+        (
+            mesh
+        )
+    );
+
     Info<< "Reading thermophysical properties\n" << endl;
 
-    autoPtr<hsCombustionThermo> pThermo
+    hsCombustionThermo& thermo = combustion->thermo();
-    (
-        hsCombustionThermo::New(mesh)
-    );
-    hsCombustionThermo& thermo = pThermo();
 
     SLGThermo slgThermo(mesh, thermo);
 
@@ -60,30 +66,8 @@
         )
     );
 
-    IOdictionary combustionProperties
+    // Set the turbulence into the combustion model
-    (
+    combustion->setTurbulence(turbulence());
-        IOobject
-        (
-            "combustionProperties",
-            runTime.constant(),
-            mesh,
-            IOobject::MUST_READ_IF_MODIFIED,
-            IOobject::NO_WRITE
-        )
-    );
-
-    Info<< "Creating combustion model\n" << endl;
-    autoPtr<combustionModel> combustion
-    (
-        combustionModel::combustionModel::New
-        (
-            combustionProperties,
-            thermo,
-            turbulence(),
-            phi,
-            rho
-        )
-    );
 
     volScalarField dQ
     (
@@ -96,7 +80,7 @@
             IOobject::AUTO_WRITE
         ),
         mesh,
-        dimensionedScalar("dQ", dimMass/pow3(dimTime)/dimLength, 0.0)
+        dimensionedScalar("dQ", dimEnergy/dimTime, 0.0)
     );
 
     Info<< "Creating field DpDt\n" << endl;

applications/solvers/combustion/fireFoam/fireFoam.C

@@ -38,9 +38,8 @@ Description
 #include "pyrolysisModel.H"
 #include "radiationModel.H"
 #include "SLGThermo.H"
-#include "hsCombustionThermo.H"
 #include "solidChemistryModel.H"
-#include "combustionModel.H"
+#include "psiCombustionModel.H"
 #include "pimpleControl.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -51,7 +50,6 @@ int main(int argc, char *argv[])
 
     #include "createTime.H"
     #include "createMesh.H"
-    #include "readChemistryProperties.H"
     #include "readGravitationalAcceleration.H"
     #include "createFields.H"
     #include "createClouds.H"

applications/solvers/combustion/fireFoam/readChemistryProperties.H

@@ -1,23 +0,0 @@
-Info<< "Reading chemistry properties\n" << endl;
-
-IOdictionary chemistryProperties
-(
-    IOobject
-    (
-        "chemistryProperties",
-        runTime.constant(),
-        mesh,
-        IOobject::MUST_READ_IF_MODIFIED,
-        IOobject::NO_WRITE,
-        false
-    )
-);
-
-Switch turbulentReaction(chemistryProperties.lookup("turbulentReaction"));
-
-dimensionedScalar Cmix("Cmix", dimless, 1.0);
-
-if (turbulentReaction)
-{
-    chemistryProperties.lookup("Cmix") >> Cmix;
-}

applications/solvers/combustion/reactingFoam/Make/options

@@ -5,7 +5,8 @@ EXE_INC = \
     -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
     -I$(LIB_SRC)/thermophysicalModels/chemistryModel/lnInclude \
     -I$(LIB_SRC)/ODE/lnInclude \
-    -I$(LIB_SRC)/finiteVolume/lnInclude
+    -I$(LIB_SRC)/finiteVolume/lnInclude \
+    -I$(LIB_SRC)/combustionModels/lnInclude
 
 EXE_LIBS = \
     -lcompressibleTurbulenceModel \
@@ -16,4 +17,5 @@ EXE_LIBS = \
     -lbasicThermophysicalModels \
     -lchemistryModel \
     -lODE \
-    -lfiniteVolume
+    -lfiniteVolume \
+    -lcombustionModels

applications/solvers/combustion/reactingFoam/YEqn.H

@@ -10,6 +10,8 @@ tmp<fv::convectionScheme<scalar> > mvConvection
 );
 
 {
+    combustion->correct();
+    dQ = combustion->dQ();
     label inertIndex = -1;
     volScalarField Yt(0.0*Y[0]);
 
@@ -19,16 +21,18 @@ tmp<fv::convectionScheme<scalar> > mvConvection
         {
             volScalarField& Yi = Y[i];
 
-            solve
+            fvScalarMatrix YiEqn
             (
                 fvm::ddt(rho, Yi)
               + mvConvection->fvmDiv(phi, Yi)
               - fvm::laplacian(turbulence->muEff(), Yi)
              ==
-                kappa*chemistry.RR(i),
+                combustion->R(Yi)
-                mesh.solver("Yi")
             );
 
+            YiEqn.relax();
+            YiEqn.solve(mesh.solver("Yi"));
+
             Yi.max(0.0);
             Yt += Yi;
         }

applications/solvers/combustion/reactingFoam/chemistry.H

@@ -1,44 +0,0 @@
-if (chemistry.chemistry())
-{
-    Info<< "Solving chemistry" << endl;
-
-    chemistry.solve
-    (
-        runTime.value() - runTime.deltaTValue(),
-        runTime.deltaTValue()
-    );
-
-    // turbulent time scale
-    if (turbulentReaction)
-    {
-        tmp<volScalarField> tepsilon(turbulence->epsilon());
-        const volScalarField& epsilon = tepsilon();
-        tmp<volScalarField> tmuEff(turbulence->muEff());
-        const volScalarField& muEff = tmuEff();
-        tmp<volScalarField> ttc(chemistry.tc());
-        const volScalarField& tc = ttc();
-
-        forAll(epsilon, i)
-        {
-            if (epsilon[i] > 0)
-            {
-                // Chalmers PaSR model
-                scalar tk = Cmix.value()*Foam::sqrt(muEff[i]/rho[i]/epsilon[i]);
-                kappa[i] =
-                    (runTime.deltaTValue() + tc[i])
-                   /(runTime.deltaTValue() + tc[i] + tk);
-            }
-            else
-            {
-                // Return to laminar combustion
-                kappa[i] = 1.0;
-            }
-        }
-    }
-    else
-    {
-        kappa = 1.0;
-    }
-
-    chemistrySh = kappa*chemistry.Sh()();
-}

applications/solvers/combustion/reactingFoam/createFields.H

@@ -1,9 +1,14 @@
-Info<< nl << "Reading thermophysicalProperties" << endl;
+Info<< "Creating combustion model\n" << endl;
-autoPtr<psiChemistryModel> pChemistry
+
+autoPtr<combustionModels::psiChemistryCombustionModel> combustion
 (
-    psiChemistryModel::New(mesh)
+    combustionModels::psiChemistryCombustionModel::New
+    (
+        mesh
+    )
 );
-psiChemistryModel& chemistry = pChemistry();
+
+psiChemistryModel& chemistry = combustion->pChemistry();
 
 hsCombustionThermo& thermo = chemistry.thermo();
 
@@ -45,20 +50,6 @@ const volScalarField& T = thermo.T();
 
 #include "compressibleCreatePhi.H"
 
-volScalarField kappa
-(
-    IOobject
-    (
-        "kappa",
-        runTime.timeName(),
-        mesh,
-        IOobject::NO_READ,
-        IOobject::AUTO_WRITE
-    ),
-    mesh,
-    dimensionedScalar("zero", dimless, 0.0)
-);
-
 Info << "Creating turbulence model.\n" << nl;
 autoPtr<compressible::turbulenceModel> turbulence
 (
@@ -71,6 +62,9 @@ autoPtr<compressible::turbulenceModel> turbulence
     )
 );
 
+// Set the turbulence into the combustion model
+combustion->setTurbulence(turbulence());
+
 Info<< "Creating field DpDt\n" << endl;
 volScalarField DpDt
 (
@@ -85,16 +79,16 @@ forAll(Y, i)
 }
 fields.add(hs);
 
-DimensionedField<scalar, volMesh> chemistrySh
+volScalarField dQ
 (
     IOobject
     (
-        "chemistry::Sh",
+        "dQ",
         runTime.timeName(),
         mesh,
         IOobject::NO_READ,
-        IOobject::NO_WRITE
+        IOobject::AUTO_WRITE
     ),
     mesh,
-    dimensionedScalar("chemistrySh", dimEnergy/dimTime/dimVolume, 0.0)
+    dimensionedScalar("dQ", dimEnergy/dimTime, 0.0)
 );

applications/solvers/combustion/reactingFoam/hsEqn.H

@@ -7,7 +7,7 @@
 //      - fvm::laplacian(turbulence->muEff(), hs)  // unit lewis no.
      ==
         DpDt
-      + chemistrySh
+      + combustion->Sh()
     );
 
     hsEqn.relax();

applications/solvers/combustion/reactingFoam/reactingFoam.C

@@ -30,10 +30,8 @@ Description
 \*---------------------------------------------------------------------------*/
 
 #include "fvCFD.H"
-#include "hCombustionThermo.H"
 #include "turbulenceModel.H"
-#include "psiChemistryModel.H"
+#include "psiChemistryCombustionModel.H"
-#include "chemistrySolver.H"
 #include "multivariateScheme.H"
 #include "pimpleControl.H"
 
@@ -44,7 +42,6 @@ int main(int argc, char *argv[])
     #include "setRootCase.H"
     #include "createTime.H"
     #include "createMesh.H"
-    #include "readChemistryProperties.H"
     #include "readGravitationalAcceleration.H"
     #include "createFields.H"
     #include "initContinuityErrs.H"
@@ -67,7 +64,6 @@ int main(int argc, char *argv[])
         runTime++;
         Info<< "Time = " << runTime.timeName() << nl << endl;
 
-        #include "chemistry.H"
         #include "rhoEqn.H"
 
         for (pimple.start(); pimple.loop(); pimple++)
@@ -88,11 +84,6 @@ int main(int argc, char *argv[])
             }
         }
 
-        if (runTime.write())
-        {
-            chemistry.dQ()().write();
-        }
-
         runTime.write();
 
         Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"

applications/solvers/combustion/reactingFoam/readChemistryProperties.H

@@ -1,23 +0,0 @@
-Info<< "Reading chemistry properties\n" << endl;
-
-IOdictionary chemistryProperties
-(
-    IOobject
-    (
-        "chemistryProperties",
-        runTime.constant(),
-        mesh,
-        IOobject::MUST_READ_IF_MODIFIED,
-        IOobject::NO_WRITE,
-        false
-    )
-);
-
-Switch turbulentReaction(chemistryProperties.lookup("turbulentReaction"));
-
-dimensionedScalar Cmix("Cmix", dimless, 1.0);
-
-if (turbulentReaction)
-{
-    chemistryProperties.lookup("Cmix") >> Cmix;
-}

applications/solvers/combustion/rhoReactingFoam/Make/options

@@ -6,7 +6,8 @@ EXE_INC = \
     -I$(LIB_SRC)/thermophysicalModels/chemistryModel/lnInclude \
     -I$(LIB_SRC)/ODE/lnInclude \
     -I$(LIB_SRC)/finiteVolume/lnInclude \
-    -I$(FOAM_SOLVERS)/combustion/reactingFoam
+    -I$(FOAM_SOLVERS)/combustion/reactingFoam \
+    -I$(LIB_SRC)/combustionModels/lnInclude
 
 
 EXE_LIBS = \
@@ -18,4 +19,5 @@ EXE_LIBS = \
     -lbasicThermophysicalModels \
     -lchemistryModel \
     -lODE \
-    -lfiniteVolume
+    -lfiniteVolume \
+    -lcombustionModels

applications/solvers/combustion/rhoReactingFoam/YEqn.H

@@ -10,6 +10,8 @@ tmp<fv::convectionScheme<scalar> > mvConvection
 );
 
 {
+    combustion->correct();
+    dQ = combustion->dQ();
     label inertIndex = -1;
     volScalarField Yt(0.0*Y[0]);
 
@@ -19,14 +21,13 @@ tmp<fv::convectionScheme<scalar> > mvConvection
         {
             volScalarField& Yi = Y[i];
 
-            solve
+            fvScalarMatrix YiEqn
             (
                 fvm::ddt(rho, Yi)
               + mvConvection->fvmDiv(phi, Yi)
               - fvm::laplacian(turbulence->muEff(), Yi)
              ==
-                kappa*chemistry.RR(i),
+                combustion->R(Yi)
-                mesh.solver("Yi")
             );
 
             Yi.max(0.0);

applications/solvers/combustion/rhoReactingFoam/createFields.H

@@ -1,9 +1,14 @@
-Info<< nl << "Reading thermophysicalProperties" << endl;
+Info<< "Creating combustion model\n" << endl;
-autoPtr<rhoChemistryModel> pChemistry
+
+autoPtr<combustionModels::rhoChemistryCombustionModel> combustion
 (
-    rhoChemistryModel::New(mesh)
+    combustionModels::rhoChemistryCombustionModel::New
+    (
+        mesh
+    )
 );
-rhoChemistryModel& chemistry = pChemistry();
+
+rhoChemistryModel& chemistry = combustion->pChemistry();
 
 hsReactionThermo& thermo = chemistry.thermo();
 
@@ -46,19 +51,6 @@ const volScalarField& T = thermo.T();
 
 #include "compressibleCreatePhi.H"
 
-volScalarField kappa
-(
-    IOobject
-    (
-        "kappa",
-        runTime.timeName(),
-        mesh,
-        IOobject::NO_READ,
-        IOobject::AUTO_WRITE
-    ),
-    mesh,
-    dimensionedScalar("zero", dimless, 0.0)
-);
 
 Info << "Creating turbulence model.\n" << nl;
 autoPtr<compressible::turbulenceModel> turbulence
@@ -72,6 +64,9 @@ autoPtr<compressible::turbulenceModel> turbulence
     )
 );
 
+// Set the turbulence into the combustion model
+combustion->setTurbulence(turbulence());
+
 Info<< "Creating field DpDt\n" << endl;
 volScalarField DpDt
 (
@@ -86,16 +81,16 @@ forAll(Y, i)
 }
 fields.add(hs);
 
-DimensionedField<scalar, volMesh> chemistrySh
+volScalarField dQ
 (
     IOobject
     (
-        "chemistry::Sh",
+        "dQ",
         runTime.timeName(),
         mesh,
         IOobject::NO_READ,
-        IOobject::NO_WRITE
+        IOobject::AUTO_WRITE
     ),
     mesh,
-    dimensionedScalar("chemistrySh", dimEnergy/dimTime/dimVolume, 0.0)
+    dimensionedScalar("dQ", dimEnergy/dimTime, 0.0)
 );

applications/solvers/combustion/rhoReactingFoam/hsEqn.H

@@ -4,9 +4,10 @@
         fvm::ddt(rho, hs)
       + mvConvection->fvmDiv(phi, hs)
       - fvm::laplacian(turbulence->alphaEff(), hs)
+//    - fvm::laplacian(turbulence->muEff(), hs)  // unit lewis no.
      ==
         DpDt
-      + chemistrySh
+      + combustion->Sh()
     );
 
     hsEqn.relax();
@@ -14,6 +15,6 @@
 
     thermo.correct();
 
-    Info<< "T gas min/max   = " << min(T).value() << ", "
+   Info<< "min/max(T) = "
-        << max(T).value() << endl;
+        << min(T).value() << ", " << max(T).value() << endl;
 }

applications/solvers/combustion/rhoReactingFoam/readChemistryProperties.H

@@ -1,23 +0,0 @@
-Info<< "Reading chemistry properties\n" << endl;
-
-IOdictionary chemistryProperties
-(
-    IOobject
-    (
-        "chemistryProperties",
-        runTime.constant(),
-        mesh,
-        IOobject::MUST_READ_IF_MODIFIED,
-        IOobject::NO_WRITE,
-        false
-    )
-);
-
-Switch turbulentReaction(chemistryProperties.lookup("turbulentReaction"));
-
-dimensionedScalar Cmix("Cmix", dimless, 1.0);
-
-if (turbulentReaction)
-{
-    chemistryProperties.lookup("Cmix") >> Cmix;
-}

applications/solvers/combustion/rhoReactingFoam/rhoReactingFoam.C

@@ -31,10 +31,8 @@ Description
 \*---------------------------------------------------------------------------*/
 
 #include "fvCFD.H"
-#include "hReactionThermo.H"
+#include "rhoChemistryCombustionModel.H"
 #include "turbulenceModel.H"
-#include "rhoChemistryModel.H"
-#include "chemistrySolver.H"
 #include "multivariateScheme.H"
 #include "pimpleControl.H"
 
@@ -45,7 +43,6 @@ int main(int argc, char *argv[])
     #include "setRootCase.H"
     #include "createTime.H"
     #include "createMesh.H"
-    #include "readChemistryProperties.H"
     #include "readGravitationalAcceleration.H"
     #include "createFields.H"
     #include "initContinuityErrs.H"
@@ -68,7 +65,6 @@ int main(int argc, char *argv[])
         runTime++;
         Info<< "Time = " << runTime.timeName() << nl << endl;
 
-        #include "chemistry.H"
         #include "rhoEqn.H"
 
         // --- Pressure-velocity PIMPLE corrector loop
@@ -92,10 +88,7 @@ int main(int argc, char *argv[])
 
         rho = thermo.rho();
 
-        if (runTime.write())
+        runTime.write();
-        {
-            chemistry.dQ()().write();
-        }
 
         Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
             << "  ClockTime = " << runTime.elapsedClockTime() << " s"

applications/solvers/heatTransfer/buoyantBaffleSimpleFoam/Make/options

@@ -1,5 +1,7 @@
 EXE_INC = \
     -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
+    -I$(LIB_SRC)/thermophysicalModels/basicSolidThermo/lnInclude \
+    -I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude \
     -I$(LIB_SRC)/turbulenceModels \
     -I$(LIB_SRC)/turbulenceModels/compressible/RAS/lnInclude \

applications/solvers/lagrangian/LTSReactingParcelFoam/LTSReactingParcelFoam.C

@@ -36,11 +36,9 @@ Description
 \*---------------------------------------------------------------------------*/
 
 #include "fvCFD.H"
-#include "hReactionThermo.H"
 #include "turbulenceModel.H"
 #include "basicReactingMultiphaseCloud.H"
-#include "rhoChemistryModel.H"
+#include "rhoChemistryCombustionModel.H"
-#include "chemistrySolver.H"
 #include "radiationModel.H"
 #include "porousZones.H"
 #include "timeActivatedExplicitSource.H"
@@ -75,7 +73,6 @@ int main(int argc, char *argv[])
 
     while (runTime.run())
     {
-        #include "readChemistryProperties.H"
         #include "readAdditionalSolutionControls.H"
         #include "readTimeControls.H"
 
@@ -85,7 +82,6 @@ int main(int argc, char *argv[])
 
         parcels.evolve();
 
-        #include "chemistry.H"
         #include "timeScales.H"
 
         #include "rhoEqn.H"
@@ -111,10 +107,7 @@ int main(int argc, char *argv[])
             }
         }
 
-        if (runTime.write())
+        runTime.write();
-        {
-            chemistry.dQ()().write();
-        }
 
         Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
             << "  ClockTime = " << runTime.elapsedClockTime() << " s"

applications/solvers/lagrangian/LTSReactingParcelFoam/Make/options

@@ -20,7 +20,8 @@ EXE_INC = \
     -I$(LIB_SRC)/ODE/lnInclude \
     -I$(LIB_SRC)/regionModels/regionModel/lnInclude \
     -I$(LIB_SRC)/regionModels/surfaceFilmModels/lnInclude \
-    -I$(LIB_SRC)/sampling/lnInclude
+    -I$(LIB_SRC)/sampling/lnInclude \
+    -I$(LIB_SRC)/combustionModels/lnInclude
 
 EXE_LIBS = \
     -lfiniteVolume \
@@ -44,4 +45,5 @@ EXE_LIBS = \
     -lODE \
     -lregionModels \
     -lsurfaceFilmModels \
-    -lsampling
+    -lsampling \
+    -lcombustionModels

applications/solvers/lagrangian/LTSReactingParcelFoam/YEqn.H

@@ -9,6 +9,9 @@ tmp<fv::convectionScheme<scalar> > mvConvection
     )
 );
 
+combustion->correct();
+dQ = combustion->dQ();
+
 if (solveSpecies)
 {
     label inertIndex = -1;
@@ -19,14 +22,15 @@ if (solveSpecies)
         if (Y[i].name() != inertSpecie)
         {
             volScalarField& Yi = Y[i];
+
             solve
             (
                 fvm::ddt(rho, Yi)
               + mvConvection->fvmDiv(phi, Yi)
               - fvm::laplacian(turbulence->muEff(), Yi)
-              ==
+             ==
                 parcels.SYi(i, Yi)
-              + kappa*chemistry.RR(i)().dimensionedInternalField()
+              + combustion->R(Yi)
               + massSource.Su(i),
                 mesh.solver("Yi")
             );

applications/solvers/lagrangian/LTSReactingParcelFoam/createFields.H

@@ -1,10 +1,14 @@
-    Info<< "Reading thermophysical properties\n" << endl;
+    Info<< "Creating combustion model\n" << endl;
 
-    autoPtr<rhoChemistryModel> pChemistry
+    autoPtr<combustionModels::rhoChemistryCombustionModel> combustion
     (
-        rhoChemistryModel::New(mesh)
+        combustionModels::rhoChemistryCombustionModel::New
+        (
+            mesh
+        )
     );
-    rhoChemistryModel& chemistry = pChemistry();
+
+    rhoChemistryModel& chemistry = combustion->pChemistry();
 
     hsReactionThermo& thermo = chemistry.thermo();
 
@@ -57,20 +61,6 @@
 
     #include "compressibleCreatePhi.H"
 
-    DimensionedField<scalar, volMesh> kappa
-    (
-        IOobject
-        (
-            "kappa",
-            runTime.timeName(),
-            mesh,
-            IOobject::NO_READ,
-            IOobject::AUTO_WRITE
-        ),
-        mesh,
-        dimensionedScalar("zero", dimless, 0.0)
-    );
-
     dimensionedScalar rhoMax
     (
         mesh.solutionDict().subDict("PIMPLE").lookup("rhoMax")
@@ -93,6 +83,9 @@
         )
     );
 
+    // Set the turbulence into the combustion model
+    combustion->setTurbulence(turbulence());
+
     Info<< "Creating multi-variate interpolation scheme\n" << endl;
     multivariateSurfaceInterpolationScheme<scalar>::fieldTable fields;
 
@@ -102,20 +95,21 @@
     }
     fields.add(hs);
 
-    DimensionedField<scalar, volMesh> chemistrySh
+    volScalarField dQ
     (
         IOobject
         (
-            "chemistry::Sh",
+            "dQ",
             runTime.timeName(),
             mesh,
             IOobject::NO_READ,
-            IOobject::NO_WRITE
+            IOobject::AUTO_WRITE
         ),
         mesh,
-        dimensionedScalar("chemistrySh", dimEnergy/dimTime/dimVolume, 0.0)
+        dimensionedScalar("dQ", dimEnergy/dimTime, 0.0)
     );
 
+
     volScalarField rDeltaT
     (
         IOobject

applications/solvers/lagrangian/LTSReactingParcelFoam/hsEqn.H

@@ -9,7 +9,7 @@
       + parcels.Sh(hs)
       + radiation->Shs(thermo)
       + energySource.Su()
-      + chemistrySh
+      + combustion->Sh()
     );
 
     hsEqn.solve();

applications/solvers/lagrangian/LTSReactingParcelFoam/readChemistryProperties.H

@@ -1,23 +0,0 @@
-// Info<< "Reading chemistry properties\n" << endl;
-
-IOdictionary chemistryProperties
-(
-    IOobject
-    (
-        "chemistryProperties",
-        runTime.constant(),
-        mesh,
-        IOobject::MUST_READ_IF_MODIFIED,
-        IOobject::NO_WRITE,
-        false
-    )
-);
-
-Switch turbulentReaction(chemistryProperties.lookup("turbulentReaction"));
-
-dimensionedScalar Cmix("Cmix", dimless, 1.0);
-
-if (turbulentReaction)
-{
-    chemistryProperties.lookup("Cmix") >> Cmix;
-}

applications/solvers/lagrangian/LTSReactingParcelFoam/timeScales.H

@@ -77,7 +77,7 @@ Info<< "Time scales min/max:" << endl;
                 DpDt
               + parcels.hsTrans()/(mesh.V()*runTime.deltaT())
               + energySource.Su()
-              + chemistrySh
+              + combustion->Sh()()
             )
            /rho
         );

applications/solvers/lagrangian/coalChemistryFoam/Make/options

@@ -21,6 +21,7 @@ EXE_INC = \
     -I$(LIB_SRC)/regionModels/surfaceFilmModels/lnInclude \
     -I$(LIB_SRC)/ODE/lnInclude \
     -I$(LIB_SRC)/sampling/lnInclude \
+    -I$(LIB_SRC)/combustionModels/lnInclude \
     -I$(FOAM_SOLVERS)/combustion/reactingFoam
 
 
@@ -47,4 +48,5 @@ EXE_LIBS = \
     -lregionModels \
     -lsurfaceFilmModels \
     -lODE \
-    -lsampling
+    -lsampling \
+    -lcombustionModels

applications/solvers/lagrangian/coalChemistryFoam/YEqn.H

@@ -11,6 +11,8 @@ tmp<fv::convectionScheme<scalar> > mvConvection
 
 
 {
+    combustion->correct();
+    dQ = combustion->dQ();
     label inertIndex = -1;
     volScalarField Yt(0.0*Y[0]);
 
@@ -19,16 +21,20 @@ tmp<fv::convectionScheme<scalar> > mvConvection
         if (Y[i].name() != inertSpecie)
         {
             volScalarField& Yi = Y[i];
-            solve
+
+            fvScalarMatrix YiEqn
             (
                 fvm::ddt(rho, Yi)
               + mvConvection->fvmDiv(phi, Yi)
               - fvm::laplacian(turbulence->muEff(), Yi)
               ==
                 coalParcels.SYi(i, Yi)
-              + kappa*chemistry.RR(i)().dimensionedInternalField()
+              + combustion->R(Yi)
             );
 
+            YiEqn.relax();
+            YiEqn.solve(mesh.solver("Yi"));
+
             Yi.max(0.0);
             Yt += Yi;
         }

applications/solvers/lagrangian/coalChemistryFoam/coalChemistryFoam.C

@@ -36,12 +36,10 @@ Description
 \*---------------------------------------------------------------------------*/
 
 #include "fvCFD.H"
-#include "hCombustionThermo.H"
 #include "turbulenceModel.H"
 #include "basicThermoCloud.H"
 #include "coalCloud.H"
-#include "psiChemistryModel.H"
+#include "psiChemistryCombustionModel.H"
-#include "chemistrySolver.H"
 #include "timeActivatedExplicitSource.H"
 #include "radiationModel.H"
 #include "SLGThermo.H"
@@ -55,7 +53,6 @@ int main(int argc, char *argv[])
 
     #include "createTime.H"
     #include "createMesh.H"
-    #include "readChemistryProperties.H"
     #include "readGravitationalAcceleration.H"
     #include "createFields.H"
     #include "createClouds.H"
@@ -88,7 +85,6 @@ int main(int argc, char *argv[])
 
         limestoneParcels.evolve();
 
-        #include "chemistry.H"
         #include "rhoEqn.H"
 
         // --- Pressure-velocity PIMPLE corrector loop
@@ -112,10 +108,7 @@ int main(int argc, char *argv[])
 
         rho = thermo.rho();
 
-        if (runTime.write())
+        runTime.write();
-        {
-            chemistry.dQ()().write();
-        }
 
         Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
             << "  ClockTime = " << runTime.elapsedClockTime() << " s"

applications/solvers/lagrangian/coalChemistryFoam/createFields.H

@@ -1,10 +1,14 @@
-    Info<< "Reading thermophysical properties\n" << endl;
+    Info<< "Creating combustion model\n" << endl;
 
-    autoPtr<psiChemistryModel> pChemistry
+    autoPtr<combustionModels::psiChemistryCombustionModel> combustion
     (
-        psiChemistryModel::New(mesh)
+        combustionModels::psiChemistryCombustionModel::New
+        (
+            mesh
+        )
     );
-    psiChemistryModel& chemistry = pChemistry();
+
+    psiChemistryModel& chemistry = combustion->pChemistry();
 
     hsCombustionThermo& thermo = chemistry.thermo();
 
@@ -96,20 +100,6 @@
 
     #include "compressibleCreatePhi.H"
 
-    DimensionedField<scalar, volMesh> kappa
-    (
-        IOobject
-        (
-            "kappa",
-            runTime.timeName(),
-            mesh,
-            IOobject::NO_READ,
-            IOobject::AUTO_WRITE
-        ),
-        mesh,
-        dimensionedScalar("zero", dimless, 0.0)
-    );
-
     Info<< "Creating turbulence model\n" << endl;
     autoPtr<compressible::turbulenceModel> turbulence
     (
@@ -122,6 +112,9 @@
         )
     );
 
+    // Set the turbulence into the combustion model
+    combustion->setTurbulence(turbulence());
+
     Info<< "Creating field DpDt\n" << endl;
     volScalarField DpDt
     (
@@ -138,16 +131,16 @@
         "hs"
     );
 
-    DimensionedField<scalar, volMesh> chemistrySh
+    volScalarField dQ
     (
         IOobject
         (
-            "chemistry::Sh",
+            "dQ",
             runTime.timeName(),
             mesh,
             IOobject::NO_READ,
-            IOobject::NO_WRITE
+            IOobject::AUTO_WRITE
         ),
         mesh,
-        dimensionedScalar("chemistrySh", dimEnergy/dimTime/dimVolume, 0.0)
+        dimensionedScalar("dQ", dimEnergy/dimTime, 0.0)
     );

applications/solvers/lagrangian/coalChemistryFoam/hsEqn.H

@@ -6,15 +6,14 @@
       - fvm::laplacian(turbulence->alphaEff(), hs)
      ==
         DpDt
+     +  combustion->Sh()
      +  coalParcels.Sh(hs)
      +  limestoneParcels.Sh(hs)
      +  enthalpySource.Su()
      +  radiation->Shs(thermo)
-     +  chemistrySh
     );
 
     hsEqn.relax();
-
     hsEqn.solve();
 
     thermo.correct();

applications/solvers/lagrangian/coalChemistryFoam/readChemistryProperties.H

@@ -1,23 +0,0 @@
-Info<< "Reading chemistry properties\n" << endl;
-
-IOdictionary chemistryProperties
-(
-    IOobject
-    (
-        "chemistryProperties",
-        runTime.constant(),
-        mesh,
-        IOobject::MUST_READ_IF_MODIFIED,
-        IOobject::NO_WRITE,
-        false
-    )
-);
-
-Switch turbulentReaction(chemistryProperties.lookup("turbulentReaction"));
-
-dimensionedScalar Cmix("Cmix", dimless, 1.0);
-
-if (turbulentReaction)
-{
-    chemistryProperties.lookup("Cmix") >> Cmix;
-}

applications/solvers/lagrangian/porousExplicitSourceReactingParcelFoam/Make/options

@@ -21,6 +21,7 @@ EXE_INC = \
     -I$(LIB_SRC)/regionModels/regionModel/lnInclude \
     -I$(LIB_SRC)/regionModels/surfaceFilmModels/lnInclude \
     -I$(LIB_SRC)/sampling/lnInclude \
+    -I$(LIB_SRC)/combustionModels/lnInclude \
     -I$(FOAM_SOLVERS)/combustion/reactingFoam
 
 
@@ -46,4 +47,5 @@ EXE_LIBS = \
     -lODE \
     -lregionModels \
     -lsurfaceFilmModels \
-    -lsampling
+    -lsampling \
+    -lcombustionModels

applications/solvers/lagrangian/porousExplicitSourceReactingParcelFoam/YEqn.H

@@ -10,6 +10,8 @@ tmp<fv::convectionScheme<scalar> > mvConvection
     )
 );
 
+combustion->correct();
+dQ = combustion->dQ();
 
 if (solveSpecies)
 {
@@ -21,6 +23,7 @@ if (solveSpecies)
         if (Y[i].name() != inertSpecie)
         {
             volScalarField& Yi = Y[i];
+
             solve
             (
                 fvm::ddt(rho, Yi)
@@ -28,7 +31,7 @@ if (solveSpecies)
               - fvm::laplacian(turbulence->muEff(), Yi)
               ==
                 parcels.SYi(i, Yi)
-              + kappa*chemistry.RR(i)().dimensionedInternalField()
+              + combustion->R(Yi)
               + massSource.Su(i),
                 mesh.solver("Yi")
             );

applications/solvers/lagrangian/porousExplicitSourceReactingParcelFoam/createFields.H

@@ -1,10 +1,14 @@
-    Info<< "Reading thermophysical properties\n" << endl;
+    Info<< "Creating combustion model\n" << endl;
 
-    autoPtr<rhoChemistryModel> pChemistry
+    autoPtr<combustionModels::rhoChemistryCombustionModel> combustion
     (
-        rhoChemistryModel::New(mesh)
+        combustionModels::rhoChemistryCombustionModel::New
+        (
+            mesh
+        )
     );
-    rhoChemistryModel& chemistry = pChemistry();
+
+    rhoChemistryModel& chemistry = combustion->pChemistry();
 
     hsReactionThermo& thermo = chemistry.thermo();
 
@@ -57,20 +61,6 @@
 
     #include "compressibleCreatePhi.H"
 
-    DimensionedField<scalar, volMesh> kappa
-    (
-        IOobject
-        (
-            "kappa",
-            runTime.timeName(),
-            mesh,
-            IOobject::NO_READ,
-            IOobject::AUTO_WRITE
-        ),
-        mesh,
-        dimensionedScalar("zero", dimless, 0.0)
-    );
-
     Info<< "Creating turbulence model\n" << endl;
     autoPtr<compressible::turbulenceModel> turbulence
     (
@@ -83,6 +73,9 @@
         )
     );
 
+     // Set the turbulence into the combustion model
+    combustion->setTurbulence(turbulence());
+
     Info<< "Creating multi-variate interpolation scheme\n" << endl;
     multivariateSurfaceInterpolationScheme<scalar>::fieldTable fields;
 
@@ -92,16 +85,16 @@
     }
     fields.add(hs);
 
-    DimensionedField<scalar, volMesh> chemistrySh
+    volScalarField dQ
     (
         IOobject
         (
-            "chemistry::Sh",
+            "dQ",
             runTime.timeName(),
             mesh,
             IOobject::NO_READ,
-            IOobject::NO_WRITE
+            IOobject::AUTO_WRITE
         ),
         mesh,
-        dimensionedScalar("chemistrySh", dimEnergy/dimTime/dimVolume, 0.0)
+        dimensionedScalar("dQ", dimEnergy/dimTime, 0.0)
     );

applications/solvers/lagrangian/porousExplicitSourceReactingParcelFoam/hsEqn.H

@@ -39,7 +39,7 @@
           + parcels.Sh(hs)
           + radiation->Shs(thermo)
           + energySource.Su()
-          + chemistrySh
+          + combustion->Sh()
         );
 
         thermo.correct();

applications/solvers/lagrangian/porousExplicitSourceReactingParcelFoam/porousExplicitSourceReactingParcelFoam.C

@@ -40,11 +40,9 @@ Description
 \*---------------------------------------------------------------------------*/
 
 #include "fvCFD.H"
-#include "hReactionThermo.H"
 #include "turbulenceModel.H"
 #include "basicReactingMultiphaseCloud.H"
-#include "rhoChemistryModel.H"
+#include "rhoChemistryCombustionModel.H"
-#include "chemistrySolver.H"
 #include "radiationModel.H"
 #include "porousZones.H"
 #include "timeActivatedExplicitSource.H"
@@ -59,7 +57,6 @@ int main(int argc, char *argv[])
 
     #include "createTime.H"
     #include "createMesh.H"
-    #include "readChemistryProperties.H"
     #include "readGravitationalAcceleration.H"
     #include "createFields.H"
     #include "createRadiationModel.H"
@@ -90,7 +87,6 @@ int main(int argc, char *argv[])
 
         parcels.evolve();
 
-        #include "chemistry.H"
         #include "rhoEqn.H"
 
         // --- Pressure-velocity PIMPLE corrector loop
@@ -114,10 +110,7 @@ int main(int argc, char *argv[])
 
         rho = thermo.rho();
 
-        if (runTime.write())
+        runTime.write();
-        {
-            chemistry.dQ()().write();
-        }
 
         Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
             << "  ClockTime = " << runTime.elapsedClockTime() << " s"

applications/solvers/lagrangian/porousExplicitSourceReactingParcelFoam/readChemistryProperties.H

@@ -1,23 +0,0 @@
-Info<< "Reading chemistry properties\n" << endl;
-
-IOdictionary chemistryProperties
-(
-    IOobject
-    (
-        "chemistryProperties",
-        runTime.constant(),
-        mesh,
-        IOobject::MUST_READ_IF_MODIFIED,
-        IOobject::NO_WRITE,
-        false
-    )
-);
-
-Switch turbulentReaction(chemistryProperties.lookup("turbulentReaction"));
-
-dimensionedScalar Cmix("Cmix", dimless, 1.0);
-
-if (turbulentReaction)
-{
-    chemistryProperties.lookup("Cmix") >> Cmix;
-}

applications/solvers/lagrangian/reactingParcelFilmFoam/Make/options

@@ -20,6 +20,7 @@ EXE_INC = \
     -I$(LIB_SRC)/lagrangian/intermediate/lnInclude \
     -I$(LIB_SRC)/ODE/lnInclude \
     -I$(LIB_SRC)/sampling/lnInclude \
+    -I$(LIB_SRC)/combustionModels/lnInclude \
     -I$(FOAM_SOLVERS)/combustion/reactingFoam
 
 
@@ -41,4 +42,5 @@ EXE_LIBS = \
     -lsurfaceFilmModels \
     -llagrangianIntermediate \
     -lODE \
-    -lsampling
+    -lsampling \
+    -lcombustionModels

applications/solvers/lagrangian/reactingParcelFilmFoam/YEqn.H

@@ -11,6 +11,8 @@ tmp<fv::convectionScheme<scalar> > mvConvection
 
 
 {
+    combustion->correct();
+    dQ = combustion->dQ();
     label inertIndex = -1;
     volScalarField Yt(0.0*Y[0]);
 
@@ -19,6 +21,7 @@ tmp<fv::convectionScheme<scalar> > mvConvection
         if (Y[i].name() != inertSpecie)
         {
             volScalarField& Yi = Y[i];
+
             solve
             (
                 fvm::ddt(rho, Yi)
@@ -27,7 +30,7 @@ tmp<fv::convectionScheme<scalar> > mvConvection
               ==
                 parcels.SYi(i, Yi)
               + surfaceFilm.Srho(i)
-              + kappa*chemistry.RR(i)().dimensionedInternalField(),
+              + combustion->R(Yi),
                 mesh.solver("Yi")
             );
 

applications/solvers/lagrangian/reactingParcelFilmFoam/createFields.H

@@ -1,10 +1,14 @@
-    Info<< "Reading thermophysical properties\n" << endl;
+    Info<< "Creating combustion model\n" << endl;
 
-    autoPtr<psiChemistryModel> pChemistry
+    autoPtr<combustionModels::psiChemistryCombustionModel> combustion
     (
-        psiChemistryModel::New(mesh)
+        combustionModels::psiChemistryCombustionModel::New
+        (
+            mesh
+        )
     );
-    psiChemistryModel& chemistry = pChemistry();
+
+    psiChemistryModel& chemistry = combustion->pChemistry();
 
     hsCombustionThermo& thermo = chemistry.thermo();
 
@@ -50,21 +54,6 @@
 
     #include "compressibleCreatePhi.H"
 
-    Info<< "Creating field kappa\n" << endl;
-    DimensionedField<scalar, volMesh> kappa
-    (
-        IOobject
-        (
-            "kappa",
-            runTime.timeName(),
-            mesh,
-            IOobject::NO_READ,
-            IOobject::AUTO_WRITE
-        ),
-        mesh,
-        dimensionedScalar("zero", dimless, 0.0)
-    );
-
     Info<< "Creating turbulence model\n" << endl;
     autoPtr<compressible::turbulenceModel> turbulence
     (
@@ -77,6 +66,9 @@
         )
     );
 
+    // Set the turbulence into the combustion model
+    combustion->setTurbulence(turbulence());
+
     Info<< "Creating field DpDt\n" << endl;
     volScalarField DpDt
     (
@@ -131,16 +123,16 @@
         additionalControlsDict.lookup("solvePrimaryRegion")
     );
 
-    DimensionedField<scalar, volMesh> chemistrySh
+    volScalarField dQ
     (
         IOobject
         (
-            "chemistry::Sh",
+            "dQ",
             runTime.timeName(),
             mesh,
             IOobject::NO_READ,
-            IOobject::NO_WRITE
+            IOobject::AUTO_WRITE
         ),
         mesh,
-        dimensionedScalar("chemistrySh", dimEnergy/dimTime/dimVolume, 0.0)
+        dimensionedScalar("dQ", dimEnergy/dimTime, 0.0)
     );

applications/solvers/lagrangian/reactingParcelFilmFoam/hsEqn.H

@@ -9,7 +9,7 @@
       + parcels.Sh(hs)
       + surfaceFilm.Sh()
       + radiation->Shs(thermo)
-      + chemistrySh
+      + combustion->Sh()
     );
 
     hsEqn.relax();

applications/solvers/lagrangian/reactingParcelFilmFoam/reactingParcelFilmFoam.C

@@ -31,12 +31,10 @@ Description
 \*---------------------------------------------------------------------------*/
 
 #include "fvCFD.H"
-#include "hCombustionThermo.H"
 #include "turbulenceModel.H"
 #include "basicReactingCloud.H"
 #include "surfaceFilmModel.H"
-#include "psiChemistryModel.H"
+#include "psiChemistryCombustionModel.H"
-#include "chemistrySolver.H"
 #include "radiationModel.H"
 #include "SLGThermo.H"
 #include "pimpleControl.H"
@@ -49,7 +47,6 @@ int main(int argc, char *argv[])
 
     #include "createTime.H"
     #include "createMesh.H"
-    #include "readChemistryProperties.H"
     #include "readGravitationalAcceleration.H"
     #include "createFields.H"
     #include "createClouds.H"
@@ -83,7 +80,6 @@ int main(int argc, char *argv[])
 
         if (solvePrimaryRegion)
         {
-            #include "chemistry.H"
             #include "rhoEqn.H"
 
             // --- PIMPLE loop
@@ -107,10 +103,7 @@ int main(int argc, char *argv[])
 
             rho = thermo.rho();
 
-            if (runTime.write())
+            runTime.write();
-            {
-                chemistry.dQ()().write();
-            }
         }
         else
         {

applications/solvers/lagrangian/reactingParcelFilmFoam/readChemistryProperties.H

@@ -1,23 +0,0 @@
-Info<< "Reading chemistry properties\n" << endl;
-
-IOdictionary chemistryProperties
-(
-    IOobject
-    (
-        "chemistryProperties",
-        runTime.constant(),
-        mesh,
-        IOobject::MUST_READ_IF_MODIFIED,
-        IOobject::NO_WRITE,
-        false
-    )
-);
-
-Switch turbulentReaction(chemistryProperties.lookup("turbulentReaction"));
-
-dimensionedScalar Cmix("Cmix", dimless, 1.0);
-
-if (turbulentReaction)
-{
-    chemistryProperties.lookup("Cmix") >> Cmix;
-}

applications/solvers/lagrangian/reactingParcelFoam/Make/options

@@ -20,6 +20,7 @@ EXE_INC = \
     -I$(LIB_SRC)/regionModels/regionModel/lnInclude \
     -I$(LIB_SRC)/regionModels/surfaceFilmModels/lnInclude \
     -I$(LIB_SRC)/sampling/lnInclude \
+    -I$(LIB_SRC)/combustionModels/lnInclude \
     -I$(FOAM_SOLVERS)/combustion/reactingFoam
 
 
@@ -45,4 +46,5 @@ EXE_LIBS = \
     -lODE \
     -lregionModels \
     -lsurfaceFilmModels \
-    -lsampling
+    -lsampling \
+    -lcombustionModels

applications/solvers/lagrangian/reactingParcelFoam/YEqn.H

@@ -11,6 +11,8 @@ tmp<fv::convectionScheme<scalar> > mvConvection
 
 
 {
+    combustion->correct();
+    dQ = combustion->dQ();
     label inertIndex = -1;
     volScalarField Yt(0.0*Y[0]);
 
@@ -19,17 +21,20 @@ tmp<fv::convectionScheme<scalar> > mvConvection
         if (Y[i].name() != inertSpecie)
         {
             volScalarField& Yi = Y[i];
-            solve
+
+            fvScalarMatrix YiEqn
             (
                 fvm::ddt(rho, Yi)
               + mvConvection->fvmDiv(phi, Yi)
               - fvm::laplacian(turbulence->muEff(), Yi)
               ==
                 parcels.SYi(i, Yi)
-              + kappa*chemistry.RR(i)().dimensionedInternalField(),
+              + combustion->R(Yi)
-                mesh.solver("Yi")
             );
 
+            YiEqn.relax();
+            YiEqn.solve(mesh.solver("Yi"));
+
             Yi.max(0.0);
             Yt += Yi;
         }

applications/solvers/lagrangian/reactingParcelFoam/createFields.H

@@ -1,10 +1,14 @@
-    Info<< "Reading thermophysical properties\n" << endl;
+    Info<< "Creating combustion model\n" << endl;
 
-    autoPtr<psiChemistryModel> pChemistry
+    autoPtr<combustionModels::psiChemistryCombustionModel> combustion
     (
-        psiChemistryModel::New(mesh)
+        combustionModels::psiChemistryCombustionModel::New
+        (
+            mesh
+        )
     );
-    psiChemistryModel& chemistry = pChemistry();
+
+    psiChemistryModel& chemistry = combustion->pChemistry();
 
     hsCombustionThermo& thermo = chemistry.thermo();
 
@@ -57,20 +61,6 @@
 
     #include "compressibleCreatePhi.H"
 
-    DimensionedField<scalar, volMesh> kappa
-    (
-        IOobject
-        (
-            "kappa",
-            runTime.timeName(),
-            mesh,
-            IOobject::NO_READ,
-            IOobject::AUTO_WRITE
-        ),
-        mesh,
-        dimensionedScalar("zero", dimless, 0.0)
-    );
-
     Info<< "Creating turbulence model\n" << endl;
     autoPtr<compressible::turbulenceModel> turbulence
     (
@@ -83,6 +73,9 @@
         )
     );
 
+    // Set the turbulence into the combustion model
+    combustion->setTurbulence(turbulence());
+
     Info<< "Creating field DpDt\n" << endl;
     volScalarField DpDt
     (
@@ -98,16 +91,16 @@
     }
     fields.add(hs);
 
-    DimensionedField<scalar, volMesh> chemistrySh
+    volScalarField dQ
     (
         IOobject
         (
-            "chemistry::Sh",
+            "dQ",
             runTime.timeName(),
             mesh,
             IOobject::NO_READ,
-            IOobject::NO_WRITE
+            IOobject::AUTO_WRITE
         ),
         mesh,
-        dimensionedScalar("chemistry::Sh", dimEnergy/dimTime/dimVolume, 0.0)
+        dimensionedScalar("dQ", dimEnergy/dimTime, 0.0)
     );

applications/solvers/lagrangian/reactingParcelFoam/hsEqn.H

@@ -8,7 +8,7 @@
         DpDt
      +  parcels.Sh(hs)
      +  radiation->Shs(thermo)
-     +  chemistrySh
+     +  combustion->Sh()
     );
 
     hEqn.relax();

applications/solvers/lagrangian/reactingParcelFoam/reactingParcelFoam.C

@@ -31,11 +31,9 @@ Description
 \*---------------------------------------------------------------------------*/
 
 #include "fvCFD.H"
-#include "hCombustionThermo.H"
 #include "turbulenceModel.H"
 #include "basicReactingCloud.H"
-#include "psiChemistryModel.H"
+#include "psiChemistryCombustionModel.H"
-#include "chemistrySolver.H"
 #include "radiationModel.H"
 #include "SLGThermo.H"
 #include "pimpleControl.H"
@@ -48,7 +46,6 @@ int main(int argc, char *argv[])
 
     #include "createTime.H"
     #include "createMesh.H"
-    #include "readChemistryProperties.H"
     #include "readGravitationalAcceleration.H"
     #include "createFields.H"
     #include "createClouds.H"
@@ -76,7 +73,6 @@ int main(int argc, char *argv[])
 
         parcels.evolve();
 
-        #include "chemistry.H"
         #include "rhoEqn.H"
 
         // --- Pressure-velocity PIMPLE corrector loop
@@ -100,10 +96,7 @@ int main(int argc, char *argv[])
 
         rho = thermo.rho();
 
-        if (runTime.write())
+        runTime.write();
-        {
-            chemistry.dQ()().write();
-        }
 
         Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
             << "  ClockTime = " << runTime.elapsedClockTime() << " s"

applications/solvers/lagrangian/reactingParcelFoam/readChemistryProperties.H

@@ -1,23 +0,0 @@
-Info<< "Reading chemistry properties\n" << endl;
-
-IOdictionary chemistryProperties
-(
-    IOobject
-    (
-        "chemistryProperties",
-        runTime.constant(),
-        mesh,
-        IOobject::MUST_READ_IF_MODIFIED,
-        IOobject::NO_WRITE,
-        false
-    )
-);
-
-Switch turbulentReaction(chemistryProperties.lookup("turbulentReaction"));
-
-dimensionedScalar Cmix("Cmix", dimless, 1.0);
-
-if (turbulentReaction)
-{
-    chemistryProperties.lookup("Cmix") >> Cmix;
-}

applications/solvers/lagrangian/sprayFoam/Make/options

@@ -21,6 +21,7 @@ EXE_INC = \
     -I$(LIB_SRC)/regionModels/regionModel/lnInclude \
     -I$(LIB_SRC)/regionModels/surfaceFilmModels/lnInclude \
     -I$(LIB_SRC)/sampling/lnInclude \
+    -I$(LIB_SRC)/combustionModels/lnInclude \
     -I$(FOAM_SOLVERS)/lagrangian/reactingParcelFoam
 
 
@@ -47,4 +48,5 @@ EXE_LIBS = \
     -lODE \
     -lregionModels \
     -lsurfaceFilmModels \
-    -lsampling
+    -lsampling \
+    -lcombustionModels

applications/solvers/lagrangian/sprayFoam/chemistry.H

@@ -1,44 +0,0 @@
-if (chemistry.chemistry())
-{
-    Info<< "Solving chemistry" << endl;
-
-    chemistry.solve
-    (
-        runTime.value() - runTime.deltaTValue(),
-        runTime.deltaTValue()
-    );
-
-    // turbulent time scale
-    if (turbulentReaction)
-    {
-        tmp<volScalarField> tepsilon(turbulence->epsilon());
-        const volScalarField& epsilon = tepsilon();
-        tmp<volScalarField> tmuEff(turbulence->muEff());
-        const volScalarField& muEff = tmuEff();
-        tmp<volScalarField> ttc(chemistry.tc());
-        const volScalarField& tc = ttc();
-
-        forAll(epsilon, i)
-        {
-            if (epsilon[i] > 0)
-            {
-                // Chalmers PaSR model
-                scalar tk = Cmix.value()*Foam::sqrt(muEff[i]/rho[i]/epsilon[i]);
-                kappa[i] =
-                    (runTime.deltaTValue() + tc[i])
-                   /(runTime.deltaTValue() + tc[i] + tk);
-            }
-            else
-            {
-                // Return to laminar combustion
-                kappa[i] = 1.0;
-            }
-        }
-    }
-    else
-    {
-        kappa = 1.0;
-    }
-
-    chemistrySh = kappa*chemistry.Sh()();
-}

applications/solvers/lagrangian/sprayFoam/sprayFoam.C

@@ -31,11 +31,9 @@ Description
 \*---------------------------------------------------------------------------*/
 
 #include "fvCFD.H"
-#include "hCombustionThermo.H"
 #include "turbulenceModel.H"
 #include "basicSprayCloud.H"
-#include "psiChemistryModel.H"
+#include "psiChemistryCombustionModel.H"
-#include "chemistrySolver.H"
 #include "radiationModel.H"
 #include "SLGThermo.H"
 #include "pimpleControl.H"
@@ -48,7 +46,6 @@ int main(int argc, char *argv[])
 
     #include "createTime.H"
     #include "createMesh.H"
-    #include "readChemistryProperties.H"
     #include "readGravitationalAcceleration.H"
     #include "createFields.H"
     #include "createClouds.H"
@@ -76,7 +73,6 @@ int main(int argc, char *argv[])
 
         parcels.evolve();
 
-        #include "chemistry.H"
         #include "rhoEqn.H"
 
         // --- Pressure-velocity PIMPLE corrector loop

applications/solvers/multiphase/bubbleFoam/bubbleFoam.C

@@ -34,6 +34,7 @@ Description
 #include "nearWallDist.H"
 #include "wallFvPatch.H"
 #include "Switch.H"
+
 #include "pimpleControl.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -47,6 +48,9 @@ int main(int argc, char *argv[])
     #include "readGravitationalAcceleration.H"
     #include "createFields.H"
     #include "initContinuityErrs.H"
+    #include "readTimeControls.H"
+    #include "CourantNo.H"
+    #include "setInitialDeltaT.H"
 
     pimpleControl pimple(mesh);
 
@@ -54,16 +58,23 @@ int main(int argc, char *argv[])
 
     Info<< "\nStarting time loop\n" << endl;
 
-    while (runTime.loop())
+    while (runTime.run())
     {
-        Info<< "Time = " << runTime.timeName() << nl << endl;
-
         #include "readBubbleFoamControls.H"
         #include "CourantNo.H"
+        #include "setDeltaT.H"
+
+        runTime++;
+        Info<< "Time = " << runTime.timeName() << nl << endl;
 
         // --- Pressure-velocity PIMPLE corrector loop
         for (pimple.start(); pimple.loop(); pimple++)
         {
+            if (pimple.nOuterCorr() != 1)
+            {
+                p.storePrevIter();
+            }
+
             #include "alphaEqn.H"
             #include "liftDragCoeffs.H"
             #include "UEqns.H"
@@ -73,7 +84,7 @@ int main(int argc, char *argv[])
             {
                 #include "pEqn.H"
 
-                if (correctAlpha)
+                if (correctAlpha && !pimple.finalIter())
                 {
                     #include "alphaEqn.H"
                 }

applications/solvers/multiphase/interMixingFoam/incompressibleThreePhaseMixture/threePhaseMixture.C

@@ -21,9 +21,6 @@ License
     You should have received a copy of the GNU General Public License
     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
 
-Class
-    threePhaseMixture
-
 \*---------------------------------------------------------------------------*/
 
 #include "threePhaseMixture.H"

applications/solvers/multiphase/twoPhaseEulerFoam/twoPhaseEulerFoam.C

@@ -26,7 +26,7 @@ Application
 
 Description
     Solver for a system of 2 incompressible fluid phases with one phase
-    dispersed, e.g. gas bubbles in a liquid.
+    dispersed, e.g. gas bubbles in a liquid or solid particles in a gas.
 
 \*---------------------------------------------------------------------------*/
 
@@ -78,10 +78,13 @@ int main(int argc, char *argv[])
         // --- Pressure-velocity PIMPLE corrector loop
         for (pimple.start(); pimple.loop(); pimple++)
         {
+            if (pimple.nOuterCorr() != 1)
+            {
+                p.storePrevIter();
+            }
+
             #include "alphaEqn.H"
-
             #include "liftDragCoeffs.H"
-
             #include "UEqns.H"
 
             // --- PISO loop

applications/utilities/mesh/conversion/gmshToFoam/gmshToFoam.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2004-2010 OpenCFD Ltd.
+    \\  /    A nd           | Copyright (C) 2004-2011 OpenCFD Ltd.
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -757,9 +757,23 @@ int main(int argc, char *argv[])
         "retain raw orientation for prisms/hexs"
     );
 
+#   include "addRegionOption.H"
+
 #   include "setRootCase.H"
 #   include "createTime.H"
 
+    Foam::word regionName;
+
+    if (args.optionReadIfPresent("region", regionName))
+    {
+        Foam::Info
+            << "Creating polyMesh for region " << regionName << endl;
+    }
+    else
+    {
+        regionName = Foam::polyMesh::defaultRegion;
+    }
+
     const bool keepOrientation = args.optionFound("keepOrientation");
     IFstream inFile(args[1]);
 
@@ -894,7 +908,7 @@ int main(int argc, char *argv[])
     (
         IOobject
         (
-            polyMesh::defaultRegion,
+            regionName,
             runTime.constant(),
             runTime
         ),
@@ -984,6 +998,7 @@ int main(int argc, char *argv[])
 
 
     //Get polyMesh to write to constant
+
     runTime.setTime(instant(runTime.constant()), 0);
 
     repatcher.repatch();
@@ -1079,6 +1094,32 @@ int main(int argc, char *argv[])
         mesh.addZones(List<pointZone*>(0), fz, cz);
     }
 
+    // Remove empty defaultFaces
+    label defaultPatchID = mesh.boundaryMesh().findPatchID(defaultFacesName);
+    if (mesh.boundaryMesh()[defaultPatchID].size() == 0)
+    {
+        List<polyPatch*> newPatchPtrList((mesh.boundaryMesh().size() - 1));
+        label newPatchI = 0;
+        forAll(mesh.boundaryMesh(), patchI)
+        {
+            if (patchI != defaultPatchID)
+            {
+                const polyPatch& patch = mesh.boundaryMesh()[patchI];
+
+                newPatchPtrList[newPatchI] = patch.clone
+                (
+                    mesh.boundaryMesh(),
+                    newPatchI,
+                    patch.size(),
+                    patch.start()
+                ).ptr();
+
+                newPatchI++;
+            }
+        }
+        repatcher.changePatches(newPatchPtrList);
+    }
+
     mesh.write();
 
     Info<< "End\n" << endl;

applications/utilities/mesh/generation/extrude/extrudeToRegionMesh/extrudeToRegionMeshDict

@@ -20,6 +20,9 @@ region  liquidFilm;
 // FaceZones to extrude
 faceZones (f0);
 
+// FaceZone shadow
+//faceZonesShadow (fBaffleShadow);
+
 // Adapt the original mesh to have directMapped patches at where the
 // faceZones are?
 // If true:
@@ -32,6 +35,10 @@ adaptMesh true;
 // Extrude 1D-columns of cells?
 oneD false;
 
+// If oneD is true. Specify which boundary is wanted between the layers
+//oneDPolyPatchType emptyPolyPatch; //wedgePolyPatch
+
+
 //- Extrusion model to use. The only logical choice is linearNormal?
 
 //- Linear extrusion in normal direction

applications/utilities/mesh/manipulation/polyDualMesh/meshDualiser.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2004-2010 OpenCFD Ltd.
+    \\  /    A nd           | Copyright (C) 2004-2011 OpenCFD Ltd.
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -21,9 +21,6 @@ License
     You should have received a copy of the GNU General Public License
     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
 
-Class
-    meshDualiser
-
 \*---------------------------------------------------------------------------*/
 
 #include "meshDualiser.H"

applications/utilities/mesh/manipulation/setSet/writeFuns.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2004-2010 OpenCFD Ltd.
+    \\  /    A nd           | Copyright (C) 2004-2011 OpenCFD Ltd.
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -22,18 +22,18 @@ License
     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
 
 Class
-    Foam::writeFunctions
+    Foam::writeFuns
 
 Description
     Various functions for collecting and writing binary data.
 
 SourceFiles
-    writeFunctions.C
+    writeFuns.C
 
 \*---------------------------------------------------------------------------*/
 
-#ifndef writeFunctions_H
+#ifndef writeFuns_H
-#define writeFunctions_H
+#define writeFuns_H
 
 #include "labelList.H"
 #include "floatScalar.H"

applications/utilities/miscellaneous/foamDebugSwitches/foamDebugSwitches.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2004-2010 OpenCFD Ltd.
+    \\  /    A nd           | Copyright (C) 2004-2011 OpenCFD Ltd.
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -61,7 +61,12 @@ int main(int argc, char *argv[])
 
     if (args.optionFound("old") || args.optionFound("new"))
     {
-        dictionary controlDict(IFstream(findEtcFile("controlDict", true))());
+        fileNameList controlDictFiles = findEtcFiles("controlDict", true);
+        dictionary controlDict;
+        forAllReverse(controlDictFiles, cdfi)
+        {
+            controlDict.merge(dictionary(IFstream(controlDictFiles[cdfi])()));
+        }
 
         wordHashSet oldDebug
         (

applications/utilities/postProcessing/dataConversion/foamToVTK/writeFuns.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2004-2010 OpenCFD Ltd.
+    \\  /    A nd           | Copyright (C) 2004-2011 OpenCFD Ltd.
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -22,18 +22,18 @@ License
     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
 
 Class
-    Foam::writeFunctions
+    Foam::writeFuns
 
 Description
     Various functions for collecting and writing binary data.
 
 SourceFiles
-    writeFunctions.C
+    writeFuns.C
 
 \*---------------------------------------------------------------------------*/
 
-#ifndef writeFunctions_H
+#ifndef writeFuns_H
-#define writeFunctions_H
+#define writeFuns_H
 
 #include "floatScalar.H"
 #include "DynamicList.H"

applications/utilities/postProcessing/sampling/sample/sampleDict

@@ -49,6 +49,9 @@ formatOptions
 //      cell          : use cell-centre value only; constant over cells (default)
 //      cellPoint     : use cell-centre and vertex values
 //      cellPointFace : use cell-centre, vertex and face values.
+//      pointMVC      : use point values only (Mean Value Coordinates)
+//      cellPatchConstrained : use cell-centre except on boundary faces where
+//        it uses the boundary value. For use with e.g. patchCloudSet.
 // 1] vertex values determined from neighbouring cell-centre values
 // 2] face values determined using the current face interpolation scheme
 //    for the field (linear, gamma, etc.)
@@ -83,6 +86,7 @@ fields
 //      uniform, face, midPoint, midPointAndFace : start and end coordinate
 //      uniform: extra number of sampling points
 //      polyLine, cloud: list of coordinates
+//      patchCloud: list of coordinates and set of patches to look for nearest
 sets
 (
     lineX1
@@ -113,8 +117,21 @@ sets
         points  ((0.049 0.049 0.00501)(0.051 0.049 0.00501));
     }
 
+    somePatchPoints
+    {
+        // Sample nearest points on selected patches. Use with
+        // interpolations:
+        // - cell (cell value)
+        // - cellPatchConstrained (boundary value)
+        // - cellPoint (interpolated boundary value)
+        type    patchCloud;
+        axis    xyz;
+        points  ((0.049 0.099 0.005)(0.051 0.054 0.005));
+        patches (".*Wall.*");
+    }
 );
 
+
 // Surface sampling definition
 //
 // 1] patches are not triangulated by default
@@ -241,4 +258,5 @@ surfaces
     }
 );
 
+
 // *********************************************************************** //

applications/utilities/preProcessing/faceAgglomerate/faceAgglomerate.C

@@ -104,7 +104,7 @@ int main(int argc, char *argv[])
         label patchI = pp.index();
         finalAgglom[patchI].setSize(pp.size(), 0);
 
-        if (pp.size() > 0 && !pp.coupled())
+        if (!pp.coupled())
         {
             if (agglomDict.found(pp.name()))
             {

applications/utilities/preProcessing/viewFactorsGen/shootRays.H

@@ -2,8 +2,7 @@
 // Pre-size by assuming a certain percentage is visible.
 
 // Maximum lenght for dynamicList
-const label maxDynListLenght = 10000;
+const label maxDynListLength = 10000;
-//label lenghtCount = 0;
 
 for (label procI = 0; procI < Pstream::nProcs(); procI++)
 {
@@ -22,69 +21,51 @@ for (label procI = 0; procI < Pstream::nProcs(); procI++)
     const pointField& remoteArea = remoteCoarseSf[procI];
     const pointField& remoteFc = remoteCoarseCf[procI];
 
-    if (myFc.size()*remoteFc.size() > 0)
+    label i = 0;
+    label j = 0;
+    do
     {
-        forAll(myFc, i)
+        for (; i < myFc.size(); i++)
         {
             const point& fc = myFc[i];
             const vector& fA = myArea[i];
 
-            forAll(remoteFc, j)
+            for (; j < remoteFc.size(); j++)//
             {
                 if (procI != Pstream::myProcNo() || i != j)
                 {
                     const point& remFc = remoteFc[j];
                     const vector& remA = remoteArea[j];
-                    const vector& d = remFc-fc;
+                    const vector& d = remFc - fc;
 
                     if (((d & fA) < 0.) && ((d & remA) > 0))
                     {
-                        //lenghtCount ++;
                         start.append(fc + 0.0001*d);
                         startIndex.append(i);
                         end.append(fc + 0.9999*d);
                         label globalI = globalNumbering.toGlobal(procI, j);
                         endIndex.append(globalI);
-
+                        if (startIndex.size() > maxDynListLength)
-                        if (startIndex.size() > maxDynListLenght)
                         {
-                            List<pointIndexHit> hitInfo(startIndex.size());
+                            break;
-                            surfacesMesh.findLine
-                            (
-                                start,
-                                end,
-                                hitInfo
-                            );
-                            surfacesMesh.findLine(start, end, hitInfo);
-                            forAll (hitInfo, rayI)
-                            {
-                                if (!hitInfo[rayI].hit())
-                                {
-                                    rayStartFace.append(startIndex[rayI]);
-                                    rayEndFace.append(endIndex[rayI]);
-                                }
-                            }
-                            //lenghtCount = 0;
-                            start.clear();
-                            startIndex.clear();
-                            end.clear();
-                            endIndex.clear();
                         }
                     }
                 }
             }
+            if (startIndex.size() > maxDynListLength)
+            {
+                break;
+            }
+
+            if (j == remoteFc.size())
+            {
+                j = 0;
+            }
         }
-    }
+
-    if (!start.empty())
-    {
         List<pointIndexHit> hitInfo(startIndex.size());
-        surfacesMesh.findLine
-        (
-            start,
-            end,
-            hitInfo
-        );
         surfacesMesh.findLine(start, end, hitInfo);
+
         forAll (hitInfo, rayI)
         {
             if (!hitInfo[rayI].hit())
@@ -93,5 +74,12 @@ for (label procI = 0; procI < Pstream::nProcs(); procI++)
                 rayEndFace.append(endIndex[rayI]);
             }
         }
-    }
+
+        start.clear();
+        startIndex.clear();
+        end.clear();
+        endIndex.clear();
+
+    }while (returnReduce(i < myFc.size(), orOp<bool>()));
+
 }

bin/paraFoam

@@ -70,6 +70,8 @@ export LC_ALL=C
 # reader extension
 extension=OpenFOAM
 
+requirePV=1
+
 # parse options
 while [ "$#" -gt 0 ]
 do
@@ -97,10 +99,12 @@ do
         ;;
     -touch)
         optTouch=true
+        requirePV=0
         shift
         ;;
     -touchAll)
         optTouch=all
+        requirePV=0
         shift
         ;;
     --)
@@ -117,6 +121,24 @@ do
 done
 
 
+#
+# check that reader module has been built
+#
+if [ $requirePV -eq 1 -a ! -f $PV_PLUGIN_PATH/libPV3FoamReader_SM.so ]
+then
+    cat<< BUILDREADER
+
+FATAL ERROR: ParaView reader module libraries do not exist
+
+Please build the reader module before continuing:
+cd \$FOAM_UTILITIES/postProcessing/graphics/PV3Readers
+./Allwclean
+./Allwmake
+
+BUILDREADER
+    exit 1
+fi
+
 #
 # check for --data=... argument
 #

etc/bashrc

@@ -32,7 +32,7 @@
 #------------------------------------------------------------------------------
 
 export WM_PROJECT=OpenFOAM
-export WM_PROJECT_VERSION=dev
+export WM_PROJECT_VERSION=dev.cvm
 
 ################################################################################
 # USER EDITABLE PART: Changes made here may be lost with the next upgrade

etc/config/settings.csh

@@ -209,10 +209,10 @@ case ThirdParty:
         breaksw
     case Gcc46:
     case Gcc46++0x:
-        set gcc_version=gcc-4.6.0
+        set gcc_version=gcc-4.6.1
-        set gmp_version=gmp-5.0.1
+        set gmp_version=gmp-5.0.2
-        set mpfr_version=mpfr-2.4.2
+        set mpfr_version=mpfr-3.0.1
-        set mpc_version=mpc-0.8.1
+        set mpc_version=mpc-0.9
         breaksw
     case Gcc45:
     case Gcc45++0x:
@@ -236,8 +236,8 @@ case ThirdParty:
         # using clang - not gcc
         setenv WM_CC 'clang'
         setenv WM_CXX 'clang++'
-        set clang_version=llvm-2.9
+        #set clang_version=llvm-2.9
-        #set clang_version=llvm-svn
+        set clang_version=llvm-svn
         breaksw
     default:
         echo

etc/config/settings.sh

@@ -228,10 +228,11 @@ OpenFOAM | ThirdParty)
         mpfr_version=mpfr-2.4.2
         ;;
     Gcc46 | Gcc46++0x)
-        gcc_version=gcc-4.6.0
+        gcc_version=gcc-4.6.1
-        gmp_version=gmp-5.0.1
+        gmp_version=gmp-5.0.2
-        mpfr_version=mpfr-2.4.2
+        mpfr_version=mpfr-3.0.1
-        mpc_version=mpc-0.8.1
+        mpc_version=mpc-0.9
+        gmpPACKAGE=gmp-5.0.2
         ;;
     Gcc45 | Gcc45++0x)
         gcc_version=gcc-4.5.2
@@ -253,8 +254,8 @@ OpenFOAM | ThirdParty)
         # using clang - not gcc
         export WM_CC='clang'
         export WM_CXX='clang++'
-        clang_version=llvm-2.9
+        #clang_version=llvm-2.9
-        #clang_version=llvm-svn
+        clang_version=llvm-svn
         ;;
     *)
         echo

etc/controlDict

@@ -24,7 +24,6 @@ Documentation
         "$WM_PROJECT_USER_DIR/html"
         "~OpenFOAM/html"
         "$WM_PROJECT_DIR/doc/Doxygen/html"
-        "$WM_PROJECT_DIR/doc/doxygen/html"
     );
     doxySourceFileExts
     (
@@ -38,6 +37,7 @@ InfoSwitches
 {
     writePrecision  6;
     writeJobInfo    0;
+    writeDictionaries 0;
 
     // Allow case-supplied C++ code (#codeStream, codedFixedValue)
     allowSystemOperations   0;

etc/cshrc

@@ -31,7 +31,7 @@
 #------------------------------------------------------------------------------
 
 setenv WM_PROJECT OpenFOAM
-setenv WM_PROJECT_VERSION dev
+setenv WM_PROJECT_VERSION dev.cvm   #dev
 
 ################################################################################
 # USER EDITABLE PART: Changes made here may be lost with the next upgrade

src/OSspecific/POSIX/POSIX.C

@@ -264,10 +264,16 @@ bool Foam::chDir(const fileName& dir)
 }
 
 
-Foam::fileName Foam::findEtcFile(const fileName& name, bool mandatory)
+Foam::fileNameList Foam::findEtcFiles
+(
+    const fileName& name,
+    bool mandatory,
+    bool findFirst
+)
 {
-    //
+    fileNameList results;
-    // search for user files in
+
+    // Search for user files in
     // * ~/.OpenFOAM/VERSION
     // * ~/.OpenFOAM
     //
@@ -277,19 +283,25 @@ Foam::fileName Foam::findEtcFile(const fileName& name, bool mandatory)
         fileName fullName = searchDir/FOAMversion/name;
         if (isFile(fullName))
         {
-            return fullName;
+            results.append(fullName);
+            if (findFirst)
+            {
+                return results;
+            }
         }
 
         fullName = searchDir/name;
         if (isFile(fullName))
         {
-            return fullName;
+            results.append(fullName);
+            if (findFirst)
+            {
+                return results;
+            }
         }
     }
 
-
+    // Search for group (site) files in
-    //
-    // search for group (site) files in
     // * $WM_PROJECT_SITE/VERSION
     // * $WM_PROJECT_SITE
     //
@@ -301,19 +313,26 @@ Foam::fileName Foam::findEtcFile(const fileName& name, bool mandatory)
             fileName fullName = searchDir/FOAMversion/name;
             if (isFile(fullName))
             {
-                return fullName;
+                results.append(fullName);
+                if (findFirst)
+                {
+                    return results;
+                }
             }
 
             fullName = searchDir/name;
             if (isFile(fullName))
             {
-                return fullName;
+                results.append(fullName);
+                if (findFirst)
+                {
+                    return results;
+                }
             }
         }
     }
     else
     {
-        //
         // OR search for group (site) files in
         // * $WM_PROJECT_INST_DIR/site/VERSION
         // * $WM_PROJECT_INST_DIR/site
@@ -324,20 +343,26 @@ Foam::fileName Foam::findEtcFile(const fileName& name, bool mandatory)
             fileName fullName = searchDir/"site"/FOAMversion/name;
             if (isFile(fullName))
             {
-                return fullName;
+                results.append(fullName);
+                if (findFirst)
+                {
+                    return results;
+                }
             }
 
             fullName = searchDir/"site"/name;
             if (isFile(fullName))
             {
-                return fullName;
+                results.append(fullName);
+                if (findFirst)
+                {
+                    return results;
+                }
             }
         }
     }
 
-
+    // Search for other (shipped) files in
-    //
-    // search for other (shipped) files in
     // * $WM_PROJECT_DIR/etc
     //
     searchDir = getEnv("WM_PROJECT_DIR");
@@ -346,24 +371,45 @@ Foam::fileName Foam::findEtcFile(const fileName& name, bool mandatory)
         fileName fullName = searchDir/"etc"/name;
         if (isFile(fullName))
         {
-            return fullName;
+            results.append(fullName);
+            if (findFirst)
+            {
+                return results;
+            }
         }
     }
 
     // Not found
-    // abort if the file is mandatory, otherwise return null
+    if (results.empty())
-    if (mandatory)
     {
-        std::cerr
+        // Abort if the file is mandatory, otherwise return null
-            << "--> FOAM FATAL ERROR in Foam::findEtcFile() :"
+        if (mandatory)
-               " could not find mandatory file\n    '"
+        {
-            << name.c_str() << "'\n\n" << std::endl;
+            std::cerr
-        ::exit(1);
+                << "--> FOAM FATAL ERROR in Foam::findEtcFiles() :"
+                   " could not find mandatory file\n    '"
+                << name.c_str() << "'\n\n" << std::endl;
+            ::exit(1);
+        }
     }
 
-    // Return null-constructed fileName rather than fileName::null
+    // Return list of matching paths or empty list if none found
-    // to avoid cyclic dependencies in the construction of globals
+    return results;
-    return fileName();
+}
+
+
+Foam::fileName Foam::findEtcFile(const fileName& name, bool mandatory)
+{
+    fileNameList results(findEtcFiles(name, mandatory, true));
+
+    if (results.size())
+    {
+        return results[0];
+    }
+    else
+    {
+        return fileName();
+    }
 }
 
 

src/OSspecific/POSIX/fileMonitor.C

@@ -21,9 +21,6 @@ License
     You should have received a copy of the GNU General Public License
     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
 
-Class
-    fileMonitor
-
 \*----------------------------------------------------------------------------*/
 
 #include "fileMonitor.H"

src/OpenFOAM/db/IOobjects/IOdictionary/IOdictionary.C

@@ -21,12 +21,6 @@ License
     You should have received a copy of the GNU General Public License
     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
 
-Description
-    IOdictionary is derived from dictionary and IOobject to give the
-    dictionary automatic IO functionality via the objectRegistry.  To facilitate
-    IO, IOdictioanry is provided with a constructor from IOobject and writeData
-    and write functions.
-
 \*---------------------------------------------------------------------------*/
 
 #include "IOdictionary.H"
@@ -37,84 +31,10 @@ Description
 
 defineTypeNameAndDebug(Foam::IOdictionary, 0);
 
-
+bool Foam::IOdictionary::writeDictionaries
-// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+(
-
+    Foam::debug::infoSwitch("writeDictionaries", 0)
-// Parallel aware reading, using non-virtual type information (typeName instead
+);
-// of type()) because of use in constructor.
-void Foam::IOdictionary::readFile(const bool masterOnly)
-{
-    if (Pstream::master() || !masterOnly)
-    {
-        if (debug)
-        {
-            Pout<< "IOdictionary : Reading " << objectPath()
-                << " from file " << endl;
-        }
-        readStream(typeName) >> *this;
-        close();
-    }
-
-    if (masterOnly && Pstream::parRun())
-    {
-        // Scatter master data using communication scheme
-
-        const List<Pstream::commsStruct>& comms =
-        (
-            (Pstream::nProcs() < Pstream::nProcsSimpleSum)
-          ? Pstream::linearCommunication()
-          : Pstream::treeCommunication()
-        );
-
-        // Master reads headerclassname from file. Make sure this gets
-        // transfered as well as contents.
-        Pstream::scatter(comms, const_cast<word&>(headerClassName()));
-        Pstream::scatter(comms, note());
-
-        // Get my communication order
-        const Pstream::commsStruct& myComm = comms[Pstream::myProcNo()];
-
-        // Reveive from up
-        if (myComm.above() != -1)
-        {
-            if (debug)
-            {
-                Pout<< "IOdictionary : Reading " << objectPath()
-                    << " from processor " << myComm.above() << endl;
-            }
-
-            // Note: use ASCII for now - binary IO of dictionaries is
-            // not currently supported
-            IPstream fromAbove
-            (
-                Pstream::scheduled,
-                myComm.above(),
-                0,
-                IOstream::ASCII
-            );
-            IOdictionary::readData(fromAbove);
-        }
-
-        // Send to my downstairs neighbours
-        forAll(myComm.below(), belowI)
-        {
-            if (debug)
-            {
-                Pout<< "IOdictionary : Sending " << objectPath()
-                    << " to processor " << myComm.below()[belowI] << endl;
-            }
-            OPstream toBelow
-            (
-                Pstream::scheduled,
-                myComm.below()[belowI],
-                0,
-                Pstream::msgType(),
-                IOstream::ASCII
-            );
-            IOdictionary::writeData(toBelow);
-        }
-    }
-}
 
 
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //

src/OpenFOAM/db/IOobjects/IOdictionary/IOdictionary.H

@@ -25,10 +25,10 @@ Class
     Foam::IOdictionary
 
 Description
-    IOdictionary is derived from dictionary and IOobject to give the
+    IOdictionary is derived from dictionary and IOobject to give the dictionary
-    dictionary automatic IO functionality via the objectRegistry.  To facilitate
+    automatic IO functionality via the objectRegistry.  To facilitate IO,
-    IO, IOdictioanry is provided with a constructor from IOobject and writeData
+    IOdictioanry is provided with a constructor from IOobject and writeData and
-    and write functions.
+    write functions.
 
 SourceFiles
     IOdictionary.C
@@ -56,12 +56,17 @@ class IOdictionary
     public regIOobject,
     public dictionary
 {
+    // Private data
+
+        static bool writeDictionaries;
+
 
     // Private Member Functions
 
-        //- read dictionary from file
+        //- Read dictionary from file
         void readFile(const bool);
 
+
 public:
 
     TypeName("dictionary");

src/OpenFOAM/db/IOobjects/IOdictionary/IOdictionaryIO.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2004-2010 OpenCFD Ltd.
+    \\  /    A nd           | Copyright (C) 2004-2011 OpenCFD Ltd.
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -21,21 +21,116 @@ License
     You should have received a copy of the GNU General Public License
     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
 
-Description
-    IOdictionary is derived from dictionary and IOobject to give the
-    dictionary automatic IO functionality via the objectRegistry.  To facilitate
-    IO, IOdictioanry is provided with a constructor from IOobject and writeData
-    and write functions.
-
 \*---------------------------------------------------------------------------*/
 
 #include "IOdictionary.H"
+#include "Pstream.H"
 
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+
+// Parallel aware reading, using non-virtual type information (typeName instead
+// of type()) because of use in constructor.
+void Foam::IOdictionary::readFile(const bool masterOnly)
+{
+    if (Pstream::master() || !masterOnly)
+    {
+        if (debug)
+        {
+            Pout<< "IOdictionary : Reading " << objectPath()
+                << " from file " << endl;
+        }
+        readStream(typeName) >> *this;
+        close();
+
+        if (writeDictionaries && Pstream::master())
+        {
+            Sout<< nl
+                << "--- IOdictionary " << name()
+                << ' ' << objectPath() << ":" << nl;
+            writeHeader(Sout);
+            writeData(Sout);
+            Sout<< "--- End of IOdictionary " << name() << nl << endl;
+        }
+    }
+
+    if (masterOnly && Pstream::parRun())
+    {
+        // Scatter master data using communication scheme
+
+        const List<Pstream::commsStruct>& comms =
+        (
+            (Pstream::nProcs() < Pstream::nProcsSimpleSum)
+          ? Pstream::linearCommunication()
+          : Pstream::treeCommunication()
+        );
+
+        // Master reads headerclassname from file. Make sure this gets
+        // transfered as well as contents.
+        Pstream::scatter(comms, const_cast<word&>(headerClassName()));
+        Pstream::scatter(comms, note());
+
+        // Get my communication order
+        const Pstream::commsStruct& myComm = comms[Pstream::myProcNo()];
+
+        // Reveive from up
+        if (myComm.above() != -1)
+        {
+            if (debug)
+            {
+                Pout<< "IOdictionary : Reading " << objectPath()
+                    << " from processor " << myComm.above() << endl;
+            }
+
+            // Note: use ASCII for now - binary IO of dictionaries is
+            // not currently supported
+            IPstream fromAbove
+            (
+                Pstream::scheduled,
+                myComm.above(),
+                0,
+                IOstream::ASCII
+            );
+            IOdictionary::readData(fromAbove);
+        }
+
+        // Send to my downstairs neighbours
+        forAll(myComm.below(), belowI)
+        {
+            if (debug)
+            {
+                Pout<< "IOdictionary : Sending " << objectPath()
+                    << " to processor " << myComm.below()[belowI] << endl;
+            }
+            OPstream toBelow
+            (
+                Pstream::scheduled,
+                myComm.below()[belowI],
+                0,
+                Pstream::msgType(),
+                IOstream::ASCII
+            );
+            IOdictionary::writeData(toBelow);
+        }
+    }
+}
+
+
+// * * * * * * * * * * * * * * * Members Functions * * * * * * * * * * * * * //
 
 bool Foam::IOdictionary::readData(Istream& is)
 {
     is >> *this;
+
+    if (writeDictionaries && Pstream::master() && !is.bad())
+    {
+        Sout<< nl
+            << "--- IOdictionary " << name()
+            << ' ' << objectPath() << ":" << nl;
+        writeHeader(Sout);
+        writeData(Sout);
+        Sout<< "--- End of IOdictionary " << name() << nl << endl;
+    }
+
     return !is.bad();
 }
 
@@ -46,4 +141,5 @@ bool Foam::IOdictionary::writeData(Ostream& os) const
     return os.good();
 }
 
+
 // ************************************************************************* //

src/OpenFOAM/db/dictionary/entry/entry.C

@@ -39,12 +39,14 @@ int Foam::entry::disableFunctionEntries
 
 Foam::entry::entry(const keyType& keyword)
 :
+    IDLList<entry>::link(),
     keyword_(keyword)
 {}
 
 
 Foam::entry::entry(const entry& e)
 :
+    IDLList<entry>::link(),
     keyword_(e.keyword_)
 {}
 

src/OpenFOAM/db/dictionary/functionEntries/codeStream/codeStream.C

@@ -197,22 +197,25 @@ Foam::functionEntries::codeStream::getFunction
             }
         }
 
-        // all processes must wait for compile to finish - except if this
+        //- We don't know whether this code was from IOdictionary
-        // file is only read on the master
+        //  (possibly read on master only) or from e.g. Field so cannot
-        bool masterOnly =
+        //  decide here.
-            (
+        //// all processes must wait for compile to finish - except if this
-                regIOobject::fileModificationChecking
+        //// file is only read on the master
-             == regIOobject::timeStampMaster
+        //bool masterOnly =
-            )
+        //    (
-         || (
+        //        regIOobject::fileModificationChecking
-                regIOobject::fileModificationChecking
+        //     == regIOobject::timeStampMaster
-             == regIOobject::inotifyMaster
+        //    )
-            );
+        // || (
-
+        //        regIOobject::fileModificationChecking
-        if (!masterOnly)
+        //     == regIOobject::inotifyMaster
-        {
+        //    );
+        //
+        //if (!masterOnly)
+        //{
             reduce(create, orOp<bool>());
-        }
+        //}
 
         if (isA<IOdictionary>(topDict(parentDict)))
         {

src/OpenFOAM/db/error/StaticAssert.H

@@ -22,7 +22,7 @@ License
     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
 
 Class
-    Foam::StaticAssertFailed
+    Foam::StaticAssertionFailed
 
 Description
     Macros and classes to provide static (compile-time) assertions.

src/OpenFOAM/db/regIOobject/regIOobjectRead.C

@@ -98,7 +98,7 @@ Foam::Istream& Foam::regIOobject::readStream()
         }
     }
 
-    // Mark as uptodate if read succesfully
+    // Mark as uptodate if read successfully
     if (watchIndex_ != -1)
     {
         time().setUnmodified(watchIndex_);

src/OpenFOAM/db/typeInfo/className.H

@@ -119,7 +119,6 @@ public:                                                                       \
     defineTypeNameWithName(Type, Name)
 //- Define the typeName as \a Name for template sub-classes
 # define defineTemplate2TypeNameWithName(Type, Name)                          \
-    template<>                                                                \
     template<>                                                                \
     defineTypeNameWithName(Type, Name)
 #endif
@@ -161,7 +160,6 @@ public:                                                                       \
     defineDebugSwitchWithName(Type, Name, DebugSwitch)
 //- Define the debug information for templates sub-classes, lookup as \a Name
 # define defineTemplate2DebugSwitchWithName(Type, Name, DebugSwitch)          \
-    template<>                                                                \
     template<>                                                                \
     defineDebugSwitchWithName(Type, Name, DebugSwitch)
 #endif

src/OpenFOAM/fields/GeometricFields/geometricOneField/geometricOneField.H

@@ -22,7 +22,7 @@ License
     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
 
 Class
-    Foam::GeometricField
+    Foam::geometricOneField
 
 Description
     A class representing the concept of a GeometricField of 1 used to avoid

src/OpenFOAM/fields/GeometricFields/geometricZeroField/geometricZeroField.H

@@ -22,7 +22,7 @@ License
     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
 
 Class
-    Foam::GeometricField
+    Foam::geometricZeroField
 
 Description
     A class representing the concept of a GeometricField of 1 used to avoid

src/OpenFOAM/global/debug/debug.C

@@ -75,10 +75,15 @@ Foam::dictionary& Foam::debug::controlDict()
 {
     if (!controlDictPtr_)
     {
-        controlDictPtr_ = new dictionary
+        fileNameList controlDictFiles = findEtcFiles("controlDict", true);
-        (
+        controlDictPtr_ = new dictionary();
-            IFstream(findEtcFile("controlDict", true))()
+        forAllReverse(controlDictFiles, cdfi)
-        );
+        {
+            controlDictPtr_->merge
+            (
+                dictionary(IFstream(controlDictFiles[cdfi])())
+            );
+        }
     }
 
     return *controlDictPtr_;

src/OpenFOAM/include/OSspecific.H

@@ -93,7 +93,7 @@ fileName cwd();
 //  else return false
 bool chDir(const fileName& dir);
 
-//- Search for a file from user/group/shipped directories.
+//- Search for files from user/group/shipped directories.
 //  The search scheme allows for version-specific and
 //  version-independent files using the following hierarchy:
 //  - \b user settings:
@@ -108,8 +108,21 @@ bool chDir(const fileName& dir);
 //  - \b other (shipped) settings:
 //    - $WM_PROJECT_DIR/etc/
 //
-//  \return The full path name or fileName() if the name cannot be found
+//  \return The list of full paths of all the matching files or
-//  Optionally abort if the file cannot be found
+//  an empty list if the name cannot be found.
+//  Optionally abort if the file cannot be found.
+//  Optionally stop search after the first file has been found.
+fileNameList findEtcFiles
+(
+    const fileName&,
+    bool mandatory=false,
+    bool findFirst=false
+);
+
+//- Search for a file using findEtcFiles.
+//  \return The full path name of the first file found in the
+//  search hierarchy or an empty fileName if the name cannot be found.
+//  Optionally abort if the file cannot be found.
 fileName findEtcFile(const fileName&, bool mandatory=false);
 
 //- Make a directory and return an error if it could not be created

src/OpenFOAM/interpolations/interpolationTable/tableReaders/tableReaders.H

@@ -21,9 +21,6 @@ License
     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::tableReader
-
 SourceFiles
     tableReaders.C
 

src/OpenFOAM/meshes/polyMesh/syncTools/syncToolsTemplates.C

@@ -1352,10 +1352,10 @@ void Foam::syncTools::syncBoundaryFaceList
                 label sz = cycPatch.size();
 
                 // Transform (copy of) data on both sides
-                Field<T> ownVals(SubField<T>(faceValues, sz, ownStart));
+                Field<T> ownVals = SubField<T>(faceValues, sz, ownStart);
                 top(nbrPatch, ownVals);
 
-                Field<T> nbrVals(SubField<T>(faceValues, sz, nbrStart));
+                Field<T> nbrVals = SubField<T>(faceValues, sz, nbrStart);
                 top(cycPatch, nbrVals);
 
                 label i0 = ownStart;

src/OpenFOAM/meshes/polyMesh/zones/zone/zone.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2004-2010 OpenCFD Ltd.
+    \\  /    A nd           | Copyright (C) 2004-2011 OpenCFD Ltd.
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -22,7 +22,7 @@ License
     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
 
 Class
-    Foam::cellZone
+    Foam::zone
 
 Description
     Base class for zones

src/OpenFOAM/primitives/globalIndexAndTransform/globalIndexAndTransformI.H

@@ -512,7 +512,7 @@ Foam::labelList Foam::globalIndexAndTransform::transformIndicesForPatches
         }
     }
 
-    label nUsedTrans = sum(mag(permutation));
+    label nUsedTrans = round(sum(mag(permutation)));
 
     if (nUsedTrans == 0)
     {

src/combustionModels/FSD/FSD.C

@@ -0,0 +1,366 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-2011 OpenCFD Ltd.
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+
+    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"
+
+namespace Foam
+{
+namespace combustionModels
+{
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+template<class CombThermoType, class ThermoType>
+FSD<CombThermoType, ThermoType>::FSD
+(
+    const word& modelType, const fvMesh& mesh
+)
+:
+    singleStepCombustion<CombThermoType, ThermoType>(modelType, mesh),
+    reactionRateFlameArea_
+    (
+        reactionRateFlameArea::New
+        (
+            this->coeffs(),
+            this->mesh(),
+            *this
+        )
+    ),
+    ft_
+    (
+        IOobject
+        (
+            "ft",
+            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")))
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructors * * * * * * * * * * * * * * * //
+
+template<class CombThermoType, class ThermoType>
+FSD<CombThermoType, ThermoType>::~FSD()
+{}
+
+
+// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
+
+template<class CombThermoType, class ThermoType>
+void FSD<CombThermoType, ThermoType>::calculateSourceNorm()
+{
+    this->singleMixture_.fresCorrect();
+
+    const label fuelI = this->singleMixture_.fuelIndex();
+
+    const volScalarField& YFuel = this->thermo_->composition().Y()[fuelI];
+
+    const volScalarField& YO2 = this->thermo_->composition().Y("O2");
+
+    const dimensionedScalar s = this->singleMixture_.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
+            (
+                "Pc",
+                U.time().timeName(),
+                U.db(),
+                IOobject::NO_READ,
+                IOobject::NO_WRITE
+            ),
+            U.mesh(),
+            dimensionedScalar("Pc", dimless, 0)
+        )
+    );
+
+    volScalarField& pc = tPc();
+
+    tmp<volScalarField> tomegaFuel
+    (
+        new volScalarField
+        (
+            IOobject
+            (
+                "omegaFuelBar",
+                U.time().timeName(),
+                U.db(),
+                IOobject::NO_READ,
+                IOobject::NO_WRITE
+            ),
+            U.mesh(),
+            dimensionedScalar
+            (
+                "omegaFuelBar",
+                omegaFuel.dimensions(),
+                0
+            )
+        )
+    );
+
+    volScalarField& omegaFuelBar = tomegaFuel();
+
+    // Calculation of the mixture fraction variance (ftVar)
+    const compressible::LESModel& lesModel =
+        YO2.db().lookupObject<compressible::LESModel>("LESProperties");
+
+    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), 1.0));
+
+    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, -1);
+    {
+        label i = 0;
+        forAll (this->singleMixture_.specieProd(), specieI)
+        {
+            if (this->singleMixture_.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->singleMixture_.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
+            (
+                "products",
+                U.time().timeName(),
+                U.db(),
+                IOobject::NO_READ,
+                IOobject::NO_WRITE
+            ),
+            U.mesh(),
+            dimensionedScalar("products", dimless, 0)
+        )
+    );
+
+    volScalarField& products = tproducts();
+
+    forAll (productsIndex, j)
+    {
+        label specieI = productsIndex[j];
+        const volScalarField& Yp = this->thermo_->composition().Y()[specieI];
+        products += Yp;
+    }
+
+    volScalarField c(max(scalar(1.0) - products/max(pc, 1e-5), 0.0));
+
+    pc = min(C_*c, scalar(1.0));
+
+    const volScalarField fres(this->singleMixture_.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
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

src/combustionModels/FSD/FSD.H

@@ -0,0 +1,173 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-2011 OpenCFD Ltd.
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+
+    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 fluctiuation 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);
+
+
+    // Destructor
+        virtual ~FSD();
+
+
+    // Evolution
+
+        //- Correct combustion rate
+        virtual void correct();
+
+
+    // I-O
+
+        //- Update properties
+        virtual bool read();
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace combustionModels
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+#   include "FSD.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/combustionModels/FSD/FSDs.C

@@ -0,0 +1,49 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-2011 OpenCFD Ltd.
+     \\/     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 "FSD.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace combustionModels
+{
+    makeCombustionTypesThermo
+    (
+        FSD,
+        psiCombustionModel,
+        gasThermoPhysics
+    );
+
+}
+}
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

src/combustionModels/FSD/reactionRateFlameAreaModels/consumptionSpeed/consumptionSpeed.C

@@ -0,0 +1,137 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-2011 OpenCFD Ltd.
+     \\/     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 2 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, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+\*---------------------------------------------------------------------------*/
+
+#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();
+
+    forAll(omega0, celli)
+    {
+        omega0[celli] = omega0Sigma(sigma[celli], 1.0);
+    }
+
+    forAll(omega0.boundaryField(), patchi)
+    {
+        forAll(omega0.boundaryField()[patchi], facei)
+        {
+            omega0.boundaryField()[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_;
+}
+
+// ************************************************************************* //