ENH: Changes to input for the PDR models

applications/solvers/combustion/PDRFoam/PDRModels/dragModels/PDRDragModel/PDRDragModel.C

@@ -81,7 +81,7 @@ bool Foam::PDRDragModel::read(const dictionary& PDRProperties)
 {
     PDRDragModelCoeffs_ = PDRProperties.subDict(type() + "Coeffs");
 
-    PDRDragModelCoeffs_.lookup("PDRDragModel") >> on_;
+    PDRDragModelCoeffs_.lookup("drag") >> on_;
 
     return true;
 }

applications/solvers/combustion/PDRFoam/PDRModels/dragModels/basic/basic.C

@@ -91,24 +91,79 @@ Foam::PDRDragModels::basic::~basic()
 
 Foam::tmp<Foam::volSymmTensorField> Foam::PDRDragModels::basic::Dcu() const
 {
-    const volScalarField& betav =
+    tmp<volSymmTensorField> tDragDcu
-        U_.db().lookupObject<volScalarField>("betav");
+    (
+        new volSymmTensorField
+        (
+            IOobject
+            (
+                "tDragDcu",
+                U_.mesh().time().constant(),
+                U_.mesh(),
+                IOobject::NO_READ,
+                IOobject::NO_WRITE
+            ),
+            U_.mesh(),
+            dimensionedSymmTensor
+            (
+                "zero",
+                dimMass/dimTime/pow(dimLength, 3),
+                pTraits<symmTensor>::zero
+            ),
+            zeroGradientFvPatchSymmTensorField::typeName
+        )
+    );
 
-    return (0.5*rho_)*CR_*mag(U_) + (Csu*I)*betav*turbulence_.muEff()*sqr(Aw_);
+    volSymmTensorField& DragDcu = tDragDcu();
+
+    if (on_)
+    {
+        const volScalarField& betav =
+            U_.db().lookupObject<volScalarField>("betav");
+
+        DragDcu = (0.5*rho_)*CR_*mag(U_) + (Csu*I)*betav*turbulence_.muEff()*sqr(Aw_);
+    }
+
+    return tDragDcu;
 }
 
 
 Foam::tmp<Foam::volScalarField> Foam::PDRDragModels::basic::Gk() const
 {
-    const volScalarField& betav =
+    tmp<volScalarField> tGk
-        U_.db().lookupObject<volScalarField>("betav");
+    (
+        new volScalarField
+        (
+            IOobject
+            (
+                "tGk",
+                U_.mesh().time().constant(),
+                U_.mesh(),
+                IOobject::NO_READ,
+                IOobject::NO_WRITE
+            ),
+            U_.mesh(),
+            dimensionedScalar("zero", dimMass/dimLength/pow(dimTime, 3), 0.0),
+            zeroGradientFvPatchVectorField::typeName
+        )
+    );
 
-    const volSymmTensorField& CT =
+    volScalarField& Gk = tGk();
-        U_.db().lookupObject<volSymmTensorField>("CT");
 
-    return
+    if (on_)
-        (0.5*rho_)*mag(U_)*(U_ & CT & U_)
+    {
-      + Csk*betav*turbulence_.muEff()*sqr(Aw_)*magSqr(U_);
+        const volScalarField& betav =
+            U_.db().lookupObject<volScalarField>("betav");
+
+        const volSymmTensorField& CT =
+            U_.db().lookupObject<volSymmTensorField>("CT");
+
+        Gk =
+            (0.5*rho_)*mag(U_)*(U_ & CT & U_)
+          + Csk*betav*turbulence_.muEff()*sqr(Aw_)*magSqr(U_);
+    }
+
+    return tGk;
 }
 
 

applications/solvers/combustion/PDRFoam/XiModels/XiEqModels/Gulder/Gulder.C

@@ -50,7 +50,12 @@ Foam::XiEqModels::Gulder::Gulder
 :
     XiEqModel(XiEqProperties, thermo, turbulence, Su),
     XiEqCoef_(readScalar(XiEqModelCoeffs_.lookup("XiEqCoef"))),
-    SuMin_(0.01*Su.average())
+    SuMin_(0.01*Su.average()),
+    uPrimeCoef_(readScalar(XiEqModelCoeffs_.lookup("uPrimeCoef"))),
+    subGridSchelkin_
+    (
+        readBool(XiEqModelCoeffs_.lookup("subGridSchelkin"))
+    )
 {}
 
 
@@ -67,9 +72,9 @@ Foam::tmp<Foam::volScalarField> Foam::XiEqModels::Gulder::XiEq() const
     volScalarField up(sqrt((2.0/3.0)*turbulence_.k()));
     const volScalarField& epsilon = turbulence_.epsilon();
 
-    if (subGridSchelkin())
+    if (subGridSchelkin_)
     {
-        up.internalField() += calculateSchelkinEffect();
+        up.internalField() += calculateSchelkinEffect(uPrimeCoef_);
     }
 
     volScalarField tauEta(sqrt(mag(thermo_.muu()/(thermo_.rhou()*epsilon))));

applications/solvers/combustion/PDRFoam/XiModels/XiEqModels/Gulder/Gulder.H

@@ -61,6 +61,12 @@ class Gulder
         //- Minimum laminar burning velocity
         const dimensionedScalar SuMin_;
 
+        //- Schelkin effect Model constant
+        scalar uPrimeCoef_;
+
+        //- Use sub-grid Schelkin effect
+        bool subGridSchelkin_;
+
 
     // Private Member Functions
 

applications/solvers/combustion/PDRFoam/XiModels/XiEqModels/SCOPEXiEq/SCOPEXiEq.C

@@ -53,6 +53,11 @@ Foam::XiEqModels::SCOPEXiEq::SCOPEXiEq
     XiEqExp_(readScalar(XiEqModelCoeffs_.lookup("XiEqExp"))),
     lCoef_(readScalar(XiEqModelCoeffs_.lookup("lCoef"))),
     SuMin_(0.01*Su.average()),
+    uPrimeCoef_(readScalar(XiEqModelCoeffs_.lookup("uPrimeCoef"))),
+    subGridSchelkin_
+    (
+        readBool(XiEqModelCoeffs_.lookup("subGridSchelkin"))
+    ),
     MaModel
     (
         IOdictionary
@@ -84,9 +89,9 @@ Foam::tmp<Foam::volScalarField> Foam::XiEqModels::SCOPEXiEq::XiEq() const
     const volScalarField& epsilon = turbulence_.epsilon();
 
     volScalarField up(sqrt((2.0/3.0)*k));
-    if (subGridSchelkin())
+    if (subGridSchelkin_)
     {
-        up.internalField() += calculateSchelkinEffect();
+        up.internalField() += calculateSchelkinEffect(uPrimeCoef_);
     }
 
     volScalarField l(lCoef_*sqrt(3.0/2.0)*up*k/epsilon);

applications/solvers/combustion/PDRFoam/XiModels/XiEqModels/SCOPEXiEq/SCOPEXiEq.H

@@ -70,6 +70,12 @@ class SCOPEXiEq
         //- Minimum Su
         dimensionedScalar SuMin_;
 
+        //- Schelkin effect Model constant
+        scalar uPrimeCoef_;
+
+        //- Use sub-grid Schelkin effect
+        bool subGridSchelkin_;
+
         //- The SCOPE laminar flame speed model used to obtain the
         //  Marstein number.  Note: the laminar flame speed need not be
         //  obtained form the same model.

applications/solvers/combustion/PDRFoam/XiModels/XiEqModels/XiEqModel/XiEqModel.C

@@ -119,7 +119,7 @@ void Foam::XiEqModel::writeFields() const
 }
 
 Foam::tmp<Foam::volScalarField>
-Foam::XiEqModel::calculateSchelkinEffect() const
+Foam::XiEqModel::calculateSchelkinEffect(const scalar uPrimeCoef) const
 {
     const fvMesh& mesh = Su_.mesh();
 
@@ -188,7 +188,7 @@ Foam::XiEqModel::calculateSchelkinEffect() const
 
     const scalarField cellWidth(pow(mesh.V(), 1.0/3.0));
 
-    const scalarField upLocal(uPrimeCoef_*sqrt((U & CT & U)*cellWidth));
+    const scalarField upLocal(uPrimeCoef*sqrt((U & CT & U)*cellWidth));
 
     const scalarField deltaUp(upLocal*(max(scalar(1.0), pow(nr, 0.5)) - 1.0));
 

applications/solvers/combustion/PDRFoam/XiModels/XiEqModels/XiEqModel/XiEqModel.H

@@ -80,12 +80,6 @@ protected:
         //
         volSymmTensorField nsv_;
 
-        //- Schelkin effect Model constant
-        scalar uPrimeCoef_;
-
-        //- Use sub-grid Schelkin effect
-        bool subGridSchelkin_;
-
 
 private:
 
@@ -162,14 +156,9 @@ public:
             return turbulence_.muEff();
         }
 
-        //- Return state of the sub-grid Schelkin effect
-        bool subGridSchelkin() const
-        {
-            return subGridSchelkin_;
-        }
 
         //- Return the sub-grid Schelkin effect
-        tmp<volScalarField> calculateSchelkinEffect() const;
+        tmp<volScalarField> calculateSchelkinEffect(const scalar) const;
 
         //- Update properties from given dictionary
         virtual bool read(const dictionary& XiEqProperties) = 0;

tutorials/combustion/PDRFoam/flamePropagationWithObstacles/constant/PDRProperties

@@ -50,26 +50,26 @@ instabilityCoeffs
             XiEqModelL
             {
                 XiEqModel       Gulder;
-                uPrimeCoef      1.0;
-                subGridSchelkin true;
 
                 GulderCoeffs
                 {
                     XiEqCoef   0.62;
+                    uPrimeCoef      1.0;
+                    subGridSchelkin true;
                 }
             }
 
             XiEqModelH
             {
                 XiEqModel       SCOPEXiEq;
-                uPrimeCoef      1.0;
-                subGridSchelkin true;
 
                 SCOPEXiEqCoeffs
                 {
                     XiEqCoef   1.6;
                     XiEqExp    0.33333;
                     lCoef      0.336;
+                    uPrimeCoef      1.0;
+                    subGridSchelkin true;
                 }
             }
         }

tutorials/combustion/PDRFoam/flamePropagationWithObstacles/system/fvSchemes

@@ -60,7 +60,7 @@ laplacianSchemes
     laplacian(muEff,U) Gauss linear limited 0.333;
     laplacian(DkEff,k) Gauss linear limited 0.333;
     laplacian(DepsilonEff,epsilon) Gauss linear limited 0.333;
-    laplacian((rho*inv((((1)*A(U))+((((0.5*rho)*CR)*mag(U))+((((Csu*(1))*betav)*muEff)*sqr(Aw)))))),p) Gauss linear limited 0.333;
+    laplacian((rho*inv((((1)*A(U))+tDragDcu))),p) Gauss linear limited 0.333;
     laplacian(Db,b) Gauss linear limited 0.333;
     laplacian(Db,ft) Gauss linear limited 0.333;
     laplacian(Db,h) Gauss linear limited 0.333;