turbulence models get thermal dissipation source term(s)

- thermalDissipation()
  corresponds to the energy lost due to viscous efffects and
  what exits the energy cascade via dissipation.

- thermalDissipationEff()
  corresponds to the energy lost due to effective viscous efffects.
  Everything that is lost from momentum. Thus essentially assumes
  turbulent equilibrium, but is what STAR-CD and Fluent seem to be using.
  Thus even if it's wrong, provide it anyhow.

- minor consistency update in comments

src/turbulenceModels/compressible/LES/LESModel/LESModel.C

@@ -76,6 +76,7 @@ LESModel::LESModel
         )
     ),
 
+    turbulence_(true),  // TODO: turbulence_(lookup("turbulence")),
     printCoeffs_(lookupOrDefault<Switch>("printCoeffs", false)),
     coeffDict_(subDictPtr(type + "Coeffs")),
 
@@ -180,6 +181,55 @@ autoPtr<LESModel> LESModel::New
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
+tmp<volScalarField> LESModel::thermalDissipation() const
+{
+    tmp<volTensorField> tgradU = fvc::grad(this->U());
+
+    return tmp<volScalarField>
+    (
+        new volScalarField
+        (
+            IOobject
+            (
+                "thermalDissipation",
+                runTime_.timeName(),
+                mesh_,
+                IOobject::NO_READ,
+                IOobject::NO_WRITE
+            ),
+            (
+                ( this->mu()*dev(twoSymm(tgradU())) ) && tgradU()
+            ) + this->rho() * this->epsilon()
+        )
+    );
+}
+
+
+tmp<volScalarField> LESModel::thermalDissipationEff() const
+{
+    tmp<volTensorField> tgradU = fvc::grad(this->U());
+
+    return tmp<volScalarField>
+    (
+        new volScalarField
+        (
+            IOobject
+            (
+                "thermalDissipationEff",
+                runTime_.timeName(),
+                mesh_,
+                IOobject::NO_READ,
+                IOobject::NO_WRITE
+            ),
+            (
+                this->muEff()*dev(twoSymm(tgradU()))
+              - ((2.0/3.0)*I) * this->rho() * this->k()
+            ) && tgradU()
+        )
+    );
+}
+
+
 void LESModel::correct(const tmp<volTensorField>&)
 {
     delta_().correct();

src/turbulenceModels/compressible/LES/LESModel/LESModel.H

@@ -33,13 +33,13 @@ Class
     Foam::compressible::LESModel
 
 Description
-    Class for all compressible flow LES SGS models.
+    Base class for all compressible flow LES SGS models.
 
-    This class defines the basic interface for a compressible flow SGS model,
+    This class defines the basic interface for a compressible flow SGS
-    and encapsulates data of value to all possible models. In particular
+    model, and encapsulates data of value to all possible models.
-    this includes references to all the dependent fields (rho, U, phi),
+    In particular this includes references to all the dependent fields
-    the physical viscosity mu, and the LESProperties dictionary,
+    (rho, U, phi), the physical viscosity mu, and the LESProperties
-    which contains the model selection and model coefficients.
+    dictionary, which contains the model selection and model coefficients.
 
 SourceFiles
     LESModel.C
@@ -80,6 +80,7 @@ protected:
 
     // Protected data
 
+        Switch turbulence_;
         Switch printCoeffs_;
         dictionary coeffDict_;
 
@@ -292,6 +293,15 @@ public:
             }
 
 
+        //- The source for the enthalpy equation resulting from
+        //  viscous and turbulent dissipation
+        virtual tmp<volScalarField> thermalDissipation() const;
+
+        //- The source for the enthalpy equation resulting from
+        //  the effective viscous dissipation
+        //  (ie, when turbulent production and dissipation are in equilibrium)
+        virtual tmp<volScalarField> thermalDissipationEff() const;
+
         //- Correct Eddy-Viscosity and related properties.
         //  This calls correct(const tmp<volTensorField>& gradU) by supplying
         //  gradU calculated locally.

src/turbulenceModels/compressible/RAS/RASModel/RASModel.C

@@ -191,6 +191,55 @@ autoPtr<RASModel> RASModel::New
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
+tmp<volScalarField> RASModel::thermalDissipation() const
+{
+    tmp<volTensorField> tgradU = fvc::grad(this->U());
+
+    return tmp<volScalarField>
+    (
+        new volScalarField
+        (
+            IOobject
+            (
+                "thermalDissipation",
+                runTime_.timeName(),
+                mesh_,
+                IOobject::NO_READ,
+                IOobject::NO_WRITE
+            ),
+            (
+                ( this->mu()*dev(twoSymm(tgradU())) ) && tgradU()
+            ) + this->rho() * this->epsilon()
+        )
+    );
+}
+
+
+tmp<volScalarField> RASModel::thermalDissipationEff() const
+{
+    tmp<volTensorField> tgradU = fvc::grad(this->U());
+
+    return tmp<volScalarField>
+    (
+        new volScalarField
+        (
+            IOobject
+            (
+                "thermalDissipationEff",
+                runTime_.timeName(),
+                mesh_,
+                IOobject::NO_READ,
+                IOobject::NO_WRITE
+            ),
+            (
+                this->muEff()*dev(twoSymm(tgradU()))
+              - ((2.0/3.0)*I) * this->rho() * this->k()
+            ) && tgradU()
+        )
+    );
+}
+
+
 scalar RASModel::yPlusLam(const scalar kappa, const scalar E) const
 {
     scalar ypl = 11.0;

src/turbulenceModels/compressible/RAS/RASModel/RASModel.H

@@ -345,6 +345,15 @@ public:
         //- Return the source term for the momentum equation
         virtual tmp<fvVectorMatrix> divDevRhoReff(volVectorField& U) const = 0;
 
+        //- The source for the enthalpy equation resulting from
+        //  viscous and turbulent dissipation
+        virtual tmp<volScalarField> thermalDissipation() const;
+
+        //- The source for the enthalpy equation resulting from
+        //  the effective viscous dissipation
+        //  (ie, when turbulent production and dissipation are in equilibrium)
+        virtual tmp<volScalarField> thermalDissipationEff() const;
+
         //- Return yPlus for the given patch
         virtual tmp<scalarField> yPlus(const label patchI) const;
 

src/turbulenceModels/compressible/RAS/laminar/laminar.C

@@ -177,6 +177,50 @@ tmp<fvVectorMatrix> laminar::divDevRhoReff(volVectorField& U) const
 }
 
 
+tmp<volScalarField> laminar::thermalDissipation() const
+{
+    tmp<volTensorField> tgradU = fvc::grad(this->U());
+
+    return tmp<volScalarField>
+    (
+        new volScalarField
+        (
+            IOobject
+            (
+                "thermalDissipation",
+                runTime_.timeName(),
+                mesh_,
+                IOobject::NO_READ,
+                IOobject::NO_WRITE
+            ),
+            ( this->mu()*dev(twoSymm(tgradU())) ) && tgradU()
+        )
+    );
+}
+
+
+tmp<volScalarField> laminar::thermalDissipationEff() const
+{
+    tmp<volTensorField> tgradU = fvc::grad(this->U());
+
+    return tmp<volScalarField>
+    (
+        new volScalarField
+        (
+            IOobject
+            (
+                "thermalDissipationEff",
+                runTime_.timeName(),
+                mesh_,
+                IOobject::NO_READ,
+                IOobject::NO_WRITE
+            ),
+            ( this->mu()*dev(twoSymm(tgradU())) ) && tgradU()
+        )
+    );
+}
+
+
 bool laminar::read()
 {
     return RASModel::read();

src/turbulenceModels/compressible/RAS/laminar/laminar.H

@@ -30,7 +30,6 @@ Description
 
 SourceFiles
     laminar.C
-    laminarCorrect.C
 
 \*---------------------------------------------------------------------------*/
 
@@ -82,42 +81,51 @@ public:
     // Member Functions
 
         //- Return the turbulence viscosity, i.e. 0 for laminar flow
-        tmp<volScalarField> mut() const;
+        virtual tmp<volScalarField> mut() const;
 
         //- Return the effective viscosity, i.e. the laminar viscosity
-        tmp<volScalarField> muEff() const
+        virtual tmp<volScalarField> muEff() const
         {
             return tmp<volScalarField>(new volScalarField("muEff", mu()));
         }
 
         //- Return the effective turbulent thermal diffusivity,
         //  i.e. the laminar thermal diffusivity
-        tmp<volScalarField> alphaEff() const
+        virtual tmp<volScalarField> alphaEff() const
         {
             return tmp<volScalarField>(new volScalarField("alphaEff", alpha()));
         }
 
         //- Return the turbulence kinetic energy, i.e. 0 for laminar flow
-        tmp<volScalarField> k() const;
+        virtual tmp<volScalarField> k() const;
 
         //- Return the turbulence kinetic energy dissipation rate,
         //  i.e. 0 for laminar flow
-        tmp<volScalarField> epsilon() const;
+        virtual tmp<volScalarField> epsilon() const;
 
         //- Return the Reynolds stress tensor, i.e. 0 for laminar flow
-        tmp<volSymmTensorField> R() const;
+        virtual tmp<volSymmTensorField> R() const;
 
         //- Return the effective stress tensor, i.e. the laminar stress
-        tmp<volSymmTensorField> devRhoReff() const;
+        virtual tmp<volSymmTensorField> devRhoReff() const;
 
         //- Return the source term for the momentum equation
-        tmp<fvVectorMatrix> divDevRhoReff(volVectorField& U) const;
+        virtual tmp<fvVectorMatrix> divDevRhoReff(volVectorField& U) const;
+
+        //- The source for the enthalpy equation resulting from
+        //  viscous and turbulent dissipation
+        virtual tmp<volScalarField> thermalDissipation() const;
+
+        //- The source for the enthalpy equation resulting from
+        //  the effective viscous dissipation
+        //  (ie, when turbulent production and dissipation are in equilibrium)
+        virtual tmp<volScalarField> thermalDissipationEff() const;
 
         //- Correct the laminar viscosity
-        void correct();
+        virtual void correct();
 
         //- Read RASProperties dictionary
-        bool read();
+        virtual bool read();
 };
 
 

src/turbulenceModels/compressible/turbulenceModel/laminar/laminar.C

@@ -52,10 +52,10 @@ laminar::laminar
     const volScalarField& rho,
     const volVectorField& U,
     const surfaceScalarField& phi,
-    const basicThermo& thermoPhysicalModel
+    const basicThermo& thermophysicalModel
 )
 :
-    turbulenceModel(rho, U, phi, thermoPhysicalModel)
+    turbulenceModel(rho, U, phi, thermophysicalModel)
 {}
 
 
@@ -66,10 +66,10 @@ autoPtr<laminar> laminar::New
     const volScalarField& rho,
     const volVectorField& U,
     const surfaceScalarField& phi,
-    const basicThermo& thermoPhysicalModel
+    const basicThermo& thermophysicalModel
 )
 {
-    return autoPtr<laminar>(new laminar(rho, U, phi, thermoPhysicalModel));
+    return autoPtr<laminar>(new laminar(rho, U, phi, thermophysicalModel));
 }
 
 
@@ -96,18 +96,6 @@ tmp<volScalarField> laminar::mut() const
 }
 
 
-tmp<volScalarField> laminar::muEff() const
-{
-    return tmp<volScalarField>(new volScalarField("muEff", mu()));
-}
-
-
-tmp<volScalarField> laminar::alphaEff() const
-{
-    return tmp<volScalarField>(new volScalarField("alphaEff", alpha()));
-}
-
-
 tmp<volScalarField> laminar::k() const
 {
     return tmp<volScalarField>
@@ -207,6 +195,50 @@ tmp<fvVectorMatrix> laminar::divDevRhoReff(volVectorField& U) const
 }
 
 
+tmp<volScalarField> laminar::thermalDissipation() const
+{
+    tmp<volTensorField> tgradU = fvc::grad(this->U());
+
+    return tmp<volScalarField>
+    (
+        new volScalarField
+        (
+            IOobject
+            (
+                "thermalDissipation",
+                runTime_.timeName(),
+                mesh_,
+                IOobject::NO_READ,
+                IOobject::NO_WRITE
+            ),
+            ( this->mu()*dev(twoSymm(tgradU())) ) && tgradU()
+        )
+    );
+}
+
+
+tmp<volScalarField> laminar::thermalDissipationEff() const
+{
+    tmp<volTensorField> tgradU = fvc::grad(this->U());
+
+    return tmp<volScalarField>
+    (
+        new volScalarField
+        (
+            IOobject
+            (
+                "thermalDissipationEff",
+                runTime_.timeName(),
+                mesh_,
+                IOobject::NO_READ,
+                IOobject::NO_WRITE
+            ),
+            ( this->mu()*dev(twoSymm(tgradU())) ) && tgradU()
+        )
+    );
+}
+
+
 bool laminar::read()
 {
     return true;

src/turbulenceModels/compressible/turbulenceModel/laminar/laminar.H

@@ -67,7 +67,7 @@ public:
             const volScalarField& rho,
             const volVectorField& U,
             const surfaceScalarField& phi,
-            const basicThermo& thermoPhysicalModel
+            const basicThermo& thermophysicalModel
         );
 
 
@@ -79,7 +79,7 @@ public:
             const volScalarField& rho,
             const volVectorField& U,
             const surfaceScalarField& phi,
-            const basicThermo& thermoPhysicalModel
+            const basicThermo& thermophysicalModel
         );
 
 
@@ -94,10 +94,17 @@ public:
         virtual tmp<volScalarField> mut() const;
 
         //- Return the effective viscosity, i.e. the laminar viscosity
-        virtual tmp<volScalarField> muEff() const;
+        virtual tmp<volScalarField> muEff() const
+        {
+            return tmp<volScalarField>(new volScalarField("muEff", mu()));
+        }
 
-        //- Return the effective turbulent thermal diffusivity
+        //- Return the effective turbulent thermal diffusivity,
-        virtual tmp<volScalarField> alphaEff() const;
+        //  i.e. the laminar thermal diffusivity
+        virtual tmp<volScalarField> alphaEff() const
+        {
+            return tmp<volScalarField>(new volScalarField("alphaEff", alpha()));
+        }
 
         //- Return the turbulence kinetic energy, i.e. 0 for laminar flow
         virtual tmp<volScalarField> k() const;
@@ -109,12 +116,21 @@ public:
         //- Return the Reynolds stress tensor, i.e. 0 for laminar flow
         virtual tmp<volSymmTensorField> R() const;
 
-        //- Return the effective stress tensor including the laminar stress
+        //- Return the effective stress tensor, i.e. the laminar stress
         virtual tmp<volSymmTensorField> devRhoReff() const;
 
         //- Return the source term for the momentum equation
         virtual tmp<fvVectorMatrix> divDevRhoReff(volVectorField& U) const;
 
+        //- The source for the enthalpy equation resulting from
+        //  viscous and turbulent dissipation
+        virtual tmp<volScalarField> thermalDissipation() const;
+
+        //- The source for the enthalpy equation resulting from
+        //  the effective viscous dissipation
+        //  (ie, when turbulent production and dissipation are in equilibrium)
+        virtual tmp<volScalarField> thermalDissipationEff() const;
+
         //- Correct the laminar viscosity
         virtual void correct();
 

src/turbulenceModels/compressible/turbulenceModel/turbulenceModel.C

@@ -117,12 +117,6 @@ autoPtr<turbulenceModel> turbulenceModel::New
 }
 
 
-// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
-
-turbulenceModel::~turbulenceModel()
-{}
-
-
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
 void turbulenceModel::correct()

src/turbulenceModels/compressible/turbulenceModel/turbulenceModel.H

@@ -65,7 +65,7 @@ namespace compressible
 {
 
 /*---------------------------------------------------------------------------*\
-                           Class turbulenceModel Declaration
+                       Class turbulenceModel Declaration
 \*---------------------------------------------------------------------------*/
 
 class turbulenceModel
@@ -144,7 +144,8 @@ public:
 
 
     //- Destructor
-    virtual ~turbulenceModel();
+    virtual ~turbulenceModel()
+    {}
 
 
     // Member Functions
@@ -209,6 +210,15 @@ public:
         //- Return the source term for the momentum equation
         virtual tmp<fvVectorMatrix> divDevRhoReff(volVectorField& U) const = 0;
 
+        //- The source for the enthalpy equation resulting from
+        //  viscous and turbulent dissipation
+        virtual tmp<volScalarField> thermalDissipation() const = 0;
+
+        //- The source for the enthalpy equation resulting from
+        //  the effective viscous dissipation
+        //  (ie, when turbulent production and dissipation are in equilibrium)
+        virtual tmp<volScalarField> thermalDissipationEff() const = 0;
+
         //- Solve the turbulence equations and correct the turbulence viscosity
         virtual void correct() = 0;
 

src/turbulenceModels/incompressible/LES/LESModel/LESModel.C

@@ -125,8 +125,8 @@ autoPtr<LESModel> LESModel::New
     {
         FatalErrorIn
         (
-            "LESModel::select(const volVectorField&, const "
+            "LESModel::New(const volVectorField& U, const "
-            "surfaceScalarField&, transportModel&)"
+            "surfaceScalarField& phi, transportModel&)"
         )   << "Unknown LESModel type " << modelName
             << endl << endl
             << "Valid LESModel types are :" << endl
@@ -138,14 +138,57 @@ autoPtr<LESModel> LESModel::New
 }
 
 
-// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
-
-LESModel::~LESModel()
-{}
-
-
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
+tmp<volScalarField> LESModel::thermalDissipation() const
+{
+    tmp<volTensorField> tgradU = fvc::grad(this->U());
+
+    return tmp<volScalarField>
+    (
+        new volScalarField
+        (
+            IOobject
+            (
+                "thermalDissipation",
+                runTime_.timeName(),
+                mesh_,
+                IOobject::NO_READ,
+                IOobject::NO_WRITE
+            ),
+            (
+                ( this->nu()*dev(twoSymm(tgradU())) ) && tgradU()
+            ) + this->epsilon()
+        )
+    );
+}
+
+
+tmp<volScalarField> LESModel::thermalDissipationEff() const
+{
+    tmp<volTensorField> tgradU = fvc::grad(this->U());
+
+    return tmp<volScalarField>
+    (
+        new volScalarField
+        (
+            IOobject
+            (
+                "thermalDissipationEff",
+                runTime_.timeName(),
+                mesh_,
+                IOobject::NO_READ,
+                IOobject::NO_WRITE
+            ),
+            (
+                this->nuEff()*dev(twoSymm(tgradU()))
+              - ((2.0/3.0)*I) * this->k()
+            ) && tgradU()
+        )
+    );
+}
+
+
 void LESModel::correct(const tmp<volTensorField>&)
 {
     turbulenceModel::correct();

src/turbulenceModels/incompressible/LES/LESModel/LESModel.H

@@ -35,10 +35,10 @@ Description
     Base class for all incompressible flow LES SGS models.
 
     This class defines the basic interface for an incompressible flow SGS
-    model, and encapsulates data of value to all possible models. In
+    model, and encapsulates data of value to all possible models.
-    particular this includes references to all the dependent fields (U,
+    In particular this includes references to all the dependent fields
-    phi), the physical viscosity nu, and the LESProperties
+    (U, phi), the physical viscosity nu, and the LESProperties
-    dictionary which contains the model selection and model coefficients.
+    dictionary, which contains the model selection and model coefficients.
 
 SourceFiles
     LESModel.C
@@ -152,7 +152,8 @@ public:
 
 
     //- Destructor
-    virtual ~LESModel();
+    virtual ~LESModel()
+    {}
 
 
     // Member Functions
@@ -241,14 +242,23 @@ public:
             }
 
 
-        //- Correct Eddy-Viscosity and related properties
+        //- The source for the enthalpy equation resulting from
-        virtual void correct(const tmp<volTensorField>& gradU);
+        //  viscous and turbulent dissipation
+        virtual tmp<volScalarField> thermalDissipation() const;
+
+        //- The source for the enthalpy equation resulting from
+        //  the effective viscous dissipation
+        //  (ie, when turbulent production and dissipation are in equilibrium)
+        virtual tmp<volScalarField> thermalDissipationEff() const;
 
         //- Correct Eddy-Viscosity and related properties.
         //  This calls correct(const tmp<volTensorField>& gradU) by supplying
         //  gradU calculated locally.
         void correct();
 
+        //- Correct Eddy-Viscosity and related properties
+        virtual void correct(const tmp<volTensorField>& gradU);
+
         //- Read LESProperties dictionary
         virtual bool read() = 0;
 };

src/turbulenceModels/incompressible/LES/dynOneEqEddy/dynOneEqEddy.C

@@ -133,12 +133,6 @@ dynOneEqEddy::dynOneEqEddy
 }
 
 
-// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
-
-dynOneEqEddy::~dynOneEqEddy()
-{}
-
-
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
 void dynOneEqEddy::correct(const tmp<volTensorField>& gradU)

src/turbulenceModels/incompressible/LES/dynOneEqEddy/dynOneEqEddy.H

@@ -116,7 +116,8 @@ public:
 
 
     //- Destructor
-    virtual ~dynOneEqEddy();
+    virtual ~dynOneEqEddy()
+    {}
 
 
     // Member Functions

src/turbulenceModels/incompressible/LES/dynSmagorinsky/dynSmagorinsky.C

@@ -124,12 +124,6 @@ dynSmagorinsky::dynSmagorinsky
 }
 
 
-// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
-
-dynSmagorinsky::~dynSmagorinsky()
-{}
-
-
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
 void dynSmagorinsky::correct(const tmp<volTensorField>& gradU)

src/turbulenceModels/incompressible/LES/dynSmagorinsky/dynSmagorinsky.H

@@ -125,7 +125,8 @@ public:
 
 
     //- Destructor
-    virtual ~dynSmagorinsky();
+    virtual ~dynSmagorinsky()
+    {}
 
 
     // Member Functions

src/turbulenceModels/incompressible/LES/locDynOneEqEddy/locDynOneEqEddy.C

@@ -126,12 +126,6 @@ locDynOneEqEddy::locDynOneEqEddy
 }
 
 
-// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
-
-locDynOneEqEddy::~locDynOneEqEddy()
-{}
-
-
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
 void locDynOneEqEddy::correct(const tmp<volTensorField>& gradU)

src/turbulenceModels/incompressible/LES/locDynOneEqEddy/locDynOneEqEddy.H

@@ -138,7 +138,8 @@ public:
 
 
     //- Destructor
-    virtual ~locDynOneEqEddy();
+    virtual ~locDynOneEqEddy()
+    {}
 
 
     // Member Functions

src/turbulenceModels/incompressible/LES/scaleSimilarity/scaleSimilarity.C

@@ -56,12 +56,6 @@ scaleSimilarity::scaleSimilarity
 }
 
 
-// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
-
-scaleSimilarity::~scaleSimilarity()
-{}
-
-
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
 tmp<volScalarField> scaleSimilarity::k() const

src/turbulenceModels/incompressible/LES/scaleSimilarity/scaleSimilarity.H

@@ -90,7 +90,8 @@ public:
 
 
     //- Destructor
-    virtual ~scaleSimilarity();
+    virtual ~scaleSimilarity()
+    {}
 
 
     // Member Functions

src/turbulenceModels/incompressible/RAS/RASModel/RASModel.C

@@ -174,14 +174,57 @@ autoPtr<RASModel> RASModel::New
 }
 
 
-// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
-
-RASModel::~RASModel()
-{}
-
-
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
+tmp<volScalarField> RASModel::thermalDissipation() const
+{
+    tmp<volTensorField> tgradU = fvc::grad(this->U());
+
+    return tmp<volScalarField>
+    (
+        new volScalarField
+        (
+            IOobject
+            (
+                "thermalDissipation",
+                runTime_.timeName(),
+                mesh_,
+                IOobject::NO_READ,
+                IOobject::NO_WRITE
+            ),
+            (
+                ( this->nu()*dev(twoSymm(tgradU())) ) && tgradU()
+            ) + this->epsilon()
+        )
+    );
+}
+
+
+tmp<volScalarField> RASModel::thermalDissipationEff() const
+{
+    tmp<volTensorField> tgradU = fvc::grad(this->U());
+
+    return tmp<volScalarField>
+    (
+        new volScalarField
+        (
+            IOobject
+            (
+                "thermalDissipationEff",
+                runTime_.timeName(),
+                mesh_,
+                IOobject::NO_READ,
+                IOobject::NO_WRITE
+            ),
+            (
+                this->nuEff()*dev(twoSymm(tgradU()))
+              - ((2.0/3.0)*I) * this->k()
+            ) && tgradU()
+        )
+    );
+}
+
+
 scalar RASModel::yPlusLam(const scalar kappa, const scalar E) const
 {
     scalar ypl = 11.0;

src/turbulenceModels/incompressible/RAS/RASModel/RASModel.H

@@ -182,7 +182,8 @@ public:
 
 
     //- Destructor
-    virtual ~RASModel();
+    virtual ~RASModel()
+    {}
 
 
     // Member Functions
@@ -310,9 +311,6 @@ public:
             );
         }
 
-        //- Return yPlus for the given patch
-        virtual tmp<scalarField> yPlus(const label patchI) const;
-
         //- Return the turbulence kinetic energy
         virtual tmp<volScalarField> k() const = 0;
 
@@ -328,6 +326,18 @@ public:
         //- Return the source term for the momentum equation
         virtual tmp<fvVectorMatrix> divDevReff(volVectorField& U) const = 0;
 
+        //- The source for the enthalpy equation resulting from
+        //  viscous and turbulent dissipation
+        virtual tmp<volScalarField> thermalDissipation() const;
+
+        //- The source for the enthalpy equation resulting from
+        //  the effective viscous dissipation
+        //  (ie, when turbulent production and dissipation are in equilibrium)
+        virtual tmp<volScalarField> thermalDissipationEff() const;
+
+        //- Return yPlus for the given patch
+        virtual tmp<scalarField> yPlus(const label patchI) const;
+
         //- Solve the turbulence equations and correct the turbulence viscosity
         virtual void correct() = 0;
 

src/turbulenceModels/incompressible/RAS/laminar/laminar.C

@@ -77,12 +77,6 @@ tmp<volScalarField> laminar::nut() const
 }
 
 
-tmp<volScalarField> laminar::nuEff() const
-{
-    return tmp<volScalarField>(new volScalarField("nuEff", nu()));
-}
-
-
 tmp<volScalarField> laminar::k() const
 {
     return tmp<volScalarField>
@@ -182,6 +176,50 @@ tmp<fvVectorMatrix> laminar::divDevReff(volVectorField& U) const
 }
 
 
+tmp<volScalarField> laminar::thermalDissipation() const
+{
+    tmp<volTensorField> tgradU = fvc::grad(this->U());
+
+    return tmp<volScalarField>
+    (
+        new volScalarField
+        (
+            IOobject
+            (
+                "thermalDissipation",
+                runTime_.timeName(),
+                mesh_,
+                IOobject::NO_READ,
+                IOobject::NO_WRITE
+            ),
+            ( this->nu()*dev(twoSymm(tgradU())) ) && tgradU()
+        )
+    );
+}
+
+
+tmp<volScalarField> laminar::thermalDissipationEff() const
+{
+    tmp<volTensorField> tgradU = fvc::grad(this->U());
+
+    return tmp<volScalarField>
+    (
+        new volScalarField
+        (
+            IOobject
+            (
+                "thermalDissipationEff",
+                runTime_.timeName(),
+                mesh_,
+                IOobject::NO_READ,
+                IOobject::NO_WRITE
+            ),
+            ( this->nu()*dev(twoSymm(tgradU())) ) && tgradU()
+        )
+    );
+}
+
+
 bool laminar::read()
 {
     return RASModel::read();

src/turbulenceModels/incompressible/RAS/laminar/laminar.H

@@ -83,7 +83,10 @@ public:
         virtual tmp<volScalarField> nut() const;
 
         //- Return the effective viscosity, i.e. the laminar viscosity
-        virtual tmp<volScalarField> nuEff() const;
+        virtual tmp<volScalarField> nuEff() const
+        {
+            return tmp<volScalarField>(new volScalarField("nuEff", nu()));
+        }
 
         //- Return the turbulence kinetic energy, i.e. 0 for laminar flow
         virtual tmp<volScalarField> k() const;
@@ -95,12 +98,21 @@ public:
         //- Return the Reynolds stress tensor, i.e. 0 for laminar flow
         virtual tmp<volSymmTensorField> R() const;
 
-        //- Return the effective stress tensor including the laminar stress
+        //- Return the effective stress tensor, i.e. the laminar stress
         virtual tmp<volSymmTensorField> devReff() const;
 
         //- Return the source term for the momentum equation
         virtual tmp<fvVectorMatrix> divDevReff(volVectorField& U) const;
 
+        //- The source for the enthalpy equation resulting from
+        //  viscous and turbulent dissipation
+        virtual tmp<volScalarField> thermalDissipation() const;
+
+        //- The source for the enthalpy equation resulting from
+        //  the effective viscous dissipation
+        //  (ie, when turbulent production and dissipation are in equilibrium)
+        virtual tmp<volScalarField> thermalDissipationEff() const;
+
         //- Correct the laminar viscosity
         virtual void correct();
 

src/turbulenceModels/incompressible/turbulenceModel/laminar/laminar.C

@@ -199,6 +199,50 @@ tmp<fvVectorMatrix> laminar::divDevReff(volVectorField& U) const
 }
 
 
+tmp<volScalarField> laminar::thermalDissipation() const
+{
+    tmp<volTensorField> tgradU = fvc::grad(this->U());
+
+    return tmp<volScalarField>
+    (
+        new volScalarField
+        (
+            IOobject
+            (
+                "thermalDissipation",
+                runTime_.timeName(),
+                mesh_,
+                IOobject::NO_READ,
+                IOobject::NO_WRITE
+            ),
+            ( this->nu()*dev(twoSymm(tgradU())) ) && tgradU()
+        )
+    );
+}
+
+
+tmp<volScalarField> laminar::thermalDissipationEff() const
+{
+    tmp<volTensorField> tgradU = fvc::grad(this->U());
+
+    return tmp<volScalarField>
+    (
+        new volScalarField
+        (
+            IOobject
+            (
+                "thermalDissipationEff",
+                runTime_.timeName(),
+                mesh_,
+                IOobject::NO_READ,
+                IOobject::NO_WRITE
+            ),
+            ( this->nu()*dev(twoSymm(tgradU())) ) && tgradU()
+        )
+    );
+}
+
+
 bool laminar::read()
 {
     return true;

src/turbulenceModels/incompressible/turbulenceModel/laminar/laminar.H

@@ -104,12 +104,21 @@ public:
         //- Return the Reynolds stress tensor, i.e. 0 for laminar flow
         virtual tmp<volSymmTensorField> R() const;
 
-        //- Return the effective stress tensor including the laminar stress
+        //- Return the effective stress tensor, i.e. the laminar stress
         virtual tmp<volSymmTensorField> devReff() const;
 
         //- Return the source term for the momentum equation
         virtual tmp<fvVectorMatrix> divDevReff(volVectorField& U) const;
 
+        //- The source for the enthalpy equation resulting from
+        //  viscous and turbulent dissipation
+        virtual tmp<volScalarField> thermalDissipation() const;
+
+        //- The source for the enthalpy equation resulting from
+        //  the effective viscous dissipation
+        //  (ie, when turbulent production and dissipation are in equilibrium)
+        virtual tmp<volScalarField> thermalDissipationEff() const;
+
         //- Correct the laminar viscosity
         virtual void correct();
 

src/turbulenceModels/incompressible/turbulenceModel/turbulenceModel.C

@@ -110,12 +110,6 @@ autoPtr<turbulenceModel> turbulenceModel::New
 }
 
 
-// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
-
-turbulenceModel::~turbulenceModel()
-{}
-
-
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
 void turbulenceModel::correct()

src/turbulenceModels/incompressible/turbulenceModel/turbulenceModel.H

@@ -139,7 +139,8 @@ public:
 
 
     //- Destructor
-    virtual ~turbulenceModel();
+    virtual ~turbulenceModel()
+    {}
 
 
     // Member Functions
@@ -189,6 +190,15 @@ public:
         //- Return the source term for the momentum equation
         virtual tmp<fvVectorMatrix> divDevReff(volVectorField& U) const = 0;
 
+        //- The source for the enthalpy equation resulting from
+        //  viscous and turbulent dissipation
+        virtual tmp<volScalarField> thermalDissipation() const = 0;
+
+        //- The source for the enthalpy equation resulting from
+        //  the effective viscous dissipation
+        //  (ie, when turbulent production and dissipation are in equilibrium)
+        virtual tmp<volScalarField> thermalDissipationEff() const = 0;
+
         //- Solve the turbulence equations and correct the turbulence viscosity
         virtual void correct() = 0;
 

tutorials/combustion/XiFoam/les/pitzDaily/constant/LESProperties

@@ -19,6 +19,8 @@ LESModel        oneEqEddy;
 
 delta           cubeRootVol;
 
+turbulence      on;
+
 printCoeffs     on;
 
 laminarCoeffs

tutorials/combustion/XiFoam/les/pitzDaily3D/constant/LESProperties

@@ -19,6 +19,8 @@ LESModel        oneEqEddy;
 
 delta           cubeRootVol;
 
+turbulence      on;
+
 printCoeffs     on;
 
 laminarCoeffs

tutorials/compressible/rhoPisoFoam/les/pitzDaily/constant/LESProperties

@@ -17,6 +17,8 @@ FoamFile
 
 LESModel        oneEqEddy;
 
+turbulence      on;
+
 printCoeffs     on;
 
 delta           cubeRootVol;