diff --git a/src/TurbulenceModels/incompressible/turbulentTransportModels/turbulentTransportModels.C b/src/TurbulenceModels/incompressible/turbulentTransportModels/turbulentTransportModels.C index 508288db76f1182612e664a4eb15eb721e01015b..02669630353c819283940e8e3e685dcd02001729 100644 --- a/src/TurbulenceModels/incompressible/turbulentTransportModels/turbulentTransportModels.C +++ b/src/TurbulenceModels/incompressible/turbulentTransportModels/turbulentTransportModels.C @@ -82,5 +82,8 @@ makeLESModel(SpalartAllmarasDDES); #include "SpalartAllmarasIDDES.H" makeLESModel(SpalartAllmarasIDDES); +#include "DeardorffDiffStress.H" +makeLESModel(DeardorffDiffStress); + // ************************************************************************* // diff --git a/src/TurbulenceModels/turbulenceModels/LES/DeardorffDiffStress/DeardorffDiffStress.C b/src/TurbulenceModels/turbulenceModels/LES/DeardorffDiffStress/DeardorffDiffStress.C new file mode 100644 index 0000000000000000000000000000000000000000..9a69a065667dd2c095772b10c3e75b5ef6ddf3e8 --- /dev/null +++ b/src/TurbulenceModels/turbulenceModels/LES/DeardorffDiffStress/DeardorffDiffStress.C @@ -0,0 +1,215 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2011-2015 OpenFOAM Foundation + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software: you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>. + +\*---------------------------------------------------------------------------*/ + +#include "DeardorffDiffStress.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ +namespace LESModels +{ + +// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // + +template<class BasicTurbulenceModel> +void DeardorffDiffStress<BasicTurbulenceModel>::correctNut() +{ + this->nut_ = Ck_*sqrt(this->k())*this->delta(); + this->nut_.correctBoundaryConditions(); +} + + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +template<class BasicTurbulenceModel> +DeardorffDiffStress<BasicTurbulenceModel>::DeardorffDiffStress +( + const alphaField& alpha, + const rhoField& rho, + const volVectorField& U, + const surfaceScalarField& alphaRhoPhi, + const surfaceScalarField& phi, + const transportModel& transport, + const word& propertiesName, + const word& type +) +: + ReynoldsStress<LESModel<BasicTurbulenceModel> > + ( + type, + alpha, + rho, + U, + alphaRhoPhi, + phi, + transport, + propertiesName + ), + + Ck_ + ( + dimensioned<scalar>::lookupOrAddToDict + ( + "Ck", + this->coeffDict_, + 0.094 + ) + ), + Cm_ + ( + dimensioned<scalar>::lookupOrAddToDict + ( + "Cm", + this->coeffDict_, + 4.13 + ) + ), + Ce_ + ( + dimensioned<scalar>::lookupOrAddToDict + ( + "Ce", + this->coeffDict_, + 1.048 + ) + ), + Cs_ + ( + dimensioned<scalar>::lookupOrAddToDict + ( + "Cs", + this->coeffDict_, + 0.22 + ) + ) +{ + if (type == typeName) + { + this->boundNormalStress(this->R_); + correctNut(); + this->printCoeffs(type); + } +} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +template<class BasicTurbulenceModel> +bool DeardorffDiffStress<BasicTurbulenceModel>::read() +{ + if (ReynoldsStress<LESModel<BasicTurbulenceModel> >::read()) + { + Ck_.readIfPresent(this->coeffDict()); + Cm_.readIfPresent(this->coeffDict()); + Ce_.readIfPresent(this->coeffDict()); + Cs_.readIfPresent(this->coeffDict()); + + return true; + } + else + { + return false; + } +} + + +template<class BasicTurbulenceModel> +tmp<volScalarField> DeardorffDiffStress<BasicTurbulenceModel>::epsilon() const +{ + volScalarField k(this->k()); + + return tmp<volScalarField> + ( + new volScalarField + ( + IOobject + ( + IOobject::groupName("epsilon", this->U_.group()), + this->runTime_.timeName(), + this->mesh_, + IOobject::NO_READ, + IOobject::NO_WRITE + ), + this->Ce_*k*sqrt(k)/this->delta() + ) + ); +} + + +template<class BasicTurbulenceModel> +void DeardorffDiffStress<BasicTurbulenceModel>::correct() +{ + if (!this->turbulence_) + { + return; + } + + // Local references + const alphaField& alpha = this->alpha_; + const rhoField& rho = this->rho_; + const surfaceScalarField& alphaRhoPhi = this->alphaRhoPhi_; + const volVectorField& U = this->U_; + volSymmTensorField& R = this->R_; + + ReynoldsStress<LESModel<BasicTurbulenceModel> >::correct(); + + tmp<volTensorField> tgradU(fvc::grad(U)); + const volTensorField& gradU = tgradU(); + + volSymmTensorField D(symm(gradU)); + + volSymmTensorField P(-twoSymm(R & gradU)); + + volScalarField k(this->k()); + + tmp<fvSymmTensorMatrix> BEqn + ( + fvm::ddt(alpha, rho, R) + + fvm::div(alphaRhoPhi, R) + - fvm::laplacian(Cs_*(k/this->epsilon())*R, R) + + fvm::Sp(Cm_*alpha*rho*sqrt(k)/this->delta(), R) + == + alpha*rho*P + + (4.0/5.0)*alpha*rho*k*D // Deardorff + //- 0.6*alpha*rho*dev(P) // LRR + - ((2.0/3.0)*(1.0 - Cm_/this->Ce_)*I)*(alpha*rho*this->epsilon()) + ); + + BEqn().relax(); + BEqn().solve(); + + this->boundNormalStress(this->R_); + correctNut(); + this->correctWallShearStress(this->R_); +} + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace LESModels +} // End namespace Foam + +// ************************************************************************* // diff --git a/src/TurbulenceModels/turbulenceModels/LES/DeardorffDiffStress/DeardorffDiffStress.H b/src/TurbulenceModels/turbulenceModels/LES/DeardorffDiffStress/DeardorffDiffStress.H new file mode 100644 index 0000000000000000000000000000000000000000..451e47aca3aad207e27dd68ee2e206a6d1cba62a --- /dev/null +++ b/src/TurbulenceModels/turbulenceModels/LES/DeardorffDiffStress/DeardorffDiffStress.H @@ -0,0 +1,166 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 2011-2012 OpenFOAM Foundation + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software: you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>. + +Class + Foam::LESModels::DeardorffDiffStress + +Group + grpLESTurbulence + +Description + Differential SGS Stress Equation Model for incompressible and + compressible flows + + Reference: + \verbatim + Deardorff, J. W. (1973). + The use of subgrid transport equations in a three-dimensional model + of atmospheric turbulence. + Journal of Fluids Engineering, 95(3), 429-438. + \endverbatim + + This SGS model uses a full balance equation for the SGS stress tensor to + simulate the behaviour of B. + +SourceFiles + DeardorffDiffStress.C + +\*---------------------------------------------------------------------------*/ + +#ifndef DeardorffDiffStress_H +#define DeardorffDiffStress_H + +#include "LESModel.H" +#include "ReynoldsStress.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ +namespace LESModels +{ + +/*---------------------------------------------------------------------------*\ + Class DeardorffDiffStress Declaration +\*---------------------------------------------------------------------------*/ + +template<class BasicTurbulenceModel> +class DeardorffDiffStress +: + public ReynoldsStress<LESModel<BasicTurbulenceModel> > +{ + // Private Member Functions + + // Disallow default bitwise copy construct and assignment + DeardorffDiffStress(const DeardorffDiffStress&); + DeardorffDiffStress& operator=(const DeardorffDiffStress&); + + +protected: + + // Protected data + + // Model constants + + dimensionedScalar Ck_; + dimensionedScalar Cm_; + dimensionedScalar Ce_; + dimensionedScalar Cs_; + + + // Protected Member Functions + + //- Update the eddy-viscosity + virtual void correctNut(); + + +public: + + typedef typename BasicTurbulenceModel::alphaField alphaField; + typedef typename BasicTurbulenceModel::rhoField rhoField; + typedef typename BasicTurbulenceModel::transportModel transportModel; + + + //- Runtime type information + TypeName("DeardorffDiffStress"); + + + // Constructors + + //- Constructor from components + DeardorffDiffStress + ( + const alphaField& alpha, + const rhoField& rho, + const volVectorField& U, + const surfaceScalarField& alphaRhoPhi, + const surfaceScalarField& phi, + const transportModel& transport, + const word& propertiesName = turbulenceModel::propertiesName, + const word& type = typeName + ); + + + //- Destructor + virtual ~DeardorffDiffStress() + {} + + + // Member Functions + + //- Read model coefficients if they have changed + virtual bool read(); + + //- Return the turbulence kinetic energy dissipation rate + virtual tmp<volScalarField> epsilon() const; + + //- Return the effective diffusivity for B + tmp<volScalarField> DBEff() const + { + return tmp<volScalarField> + ( + new volScalarField("DBEff", this->nut_ + this->nu()) + ); + } + + //- Correct sub-grid stress, eddy-Viscosity and related properties + virtual void correct(); +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace LESModels +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#ifdef NoRepository +# include "DeardorffDiffStress.C" +#endif + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/src/TurbulenceModels/turbulenceModels/LES/WALE/WALE.H b/src/TurbulenceModels/turbulenceModels/LES/WALE/WALE.H index 7cae23f7da7e8f3bd0a70a7d1462e1bf09e9834e..53ae8fe5db1cc56ee03bb930349f8b130b7e5069 100644 --- a/src/TurbulenceModels/turbulenceModels/LES/WALE/WALE.H +++ b/src/TurbulenceModels/turbulenceModels/LES/WALE/WALE.H @@ -100,7 +100,7 @@ protected: // calculated from the given velocity gradient tmp<volScalarField> k(const volTensorField& gradU) const; - //- Update the SGS eddy viscosity + //- Update the SGS eddy-viscosity virtual void correctNut(); diff --git a/src/TurbulenceModels/turbulenceModels/LES/kEqn/kEqn.H b/src/TurbulenceModels/turbulenceModels/LES/kEqn/kEqn.H index b525666a798f6cd3aa24c3aacee9dd992fd0ccb2..c5052249eccf9b33c0ed89e37c02b385ee8c8c9c 100644 --- a/src/TurbulenceModels/turbulenceModels/LES/kEqn/kEqn.H +++ b/src/TurbulenceModels/turbulenceModels/LES/kEqn/kEqn.H @@ -50,7 +50,7 @@ Description The default model coefficients correspond to the following: \verbatim - NicenoKEqnCoeffs + kEqnCoeffs { Ck 0.094; Ce 1.048; @@ -165,7 +165,7 @@ public: ); } - //- Correct Eddy-Viscosity and related properties + //- Correct eddy-Viscosity and related properties virtual void correct(); }; diff --git a/src/TurbulenceModels/turbulenceModels/ReynoldsStress/ReynoldsStress.C b/src/TurbulenceModels/turbulenceModels/ReynoldsStress/ReynoldsStress.C index 345e69d2a3c09f9c6e631e0f90050c2cddce6513..4f4cba8d2719d04a71601ab57ec4de661372996e 100644 --- a/src/TurbulenceModels/turbulenceModels/ReynoldsStress/ReynoldsStress.C +++ b/src/TurbulenceModels/turbulenceModels/ReynoldsStress/ReynoldsStress.C @@ -26,6 +26,82 @@ License #include "ReynoldsStress.H" #include "fvc.H" #include "fvm.H" +#include "wallFvPatch.H" + +// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * // + +template<class BasicTurbulenceModel> +void Foam::ReynoldsStress<BasicTurbulenceModel>::boundNormalStress +( + volSymmTensorField& R +) const +{ + scalar kMin = this->kMin_.value(); + + R.max + ( + dimensionedSymmTensor + ( + "zero", + R.dimensions(), + symmTensor + ( + kMin, -GREAT, -GREAT, + kMin, -GREAT, + kMin + ) + ) + ); +} + + +template<class BasicTurbulenceModel> +void Foam::ReynoldsStress<BasicTurbulenceModel>::correctWallShearStress +( + volSymmTensorField& R +) const +{ + const fvPatchList& patches = this->mesh_.boundary(); + + forAll(patches, patchi) + { + const fvPatch& curPatch = patches[patchi]; + + if (isA<wallFvPatch>(curPatch)) + { + symmTensorField& Rw = R.boundaryField()[patchi]; + + const scalarField& nutw = this->nut_.boundaryField()[patchi]; + + const vectorField snGradU + ( + this->U_.boundaryField()[patchi].snGrad() + ); + + const vectorField& faceAreas + = this->mesh_.Sf().boundaryField()[patchi]; + + const scalarField& magFaceAreas + = this->mesh_.magSf().boundaryField()[patchi]; + + forAll(curPatch, facei) + { + // Calculate near-wall velocity gradient + tensor gradUw + = (faceAreas[facei]/magFaceAreas[facei])*snGradU[facei]; + + // Calculate near-wall shear-stress tensor + tensor tauw = -nutw[facei]*2*dev(symm(gradUw)); + + // Reset the shear components of the stress tensor + Rw[facei].xy() = tauw.xy(); + Rw[facei].xz() = tauw.xz(); + Rw[facei].yz() = tauw.yz(); + } + } + } +} + // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // @@ -113,7 +189,7 @@ template<class BasicTurbulenceModel> Foam::tmp<Foam::volScalarField> Foam::ReynoldsStress<BasicTurbulenceModel>::k() const { - tmp<Foam::volScalarField> tk(tr(R_)); + tmp<Foam::volScalarField> tk(0.5*tr(R_)); tk().rename("k"); return tk; } @@ -168,7 +244,7 @@ Foam::ReynoldsStress<BasicTurbulenceModel>::divDevRhoReff "laplacian(nuEff,U)" ) - fvm::laplacian(this->alpha_*this->rho_*this->nuEff(), U) - - fvc::div(this->alpha_*this->rho_*nu()*dev2(T(fvc::grad(U)))) + - fvc::div(this->alpha_*this->rho_*this->nu()*dev2(T(fvc::grad(U)))) ); } else @@ -183,7 +259,7 @@ Foam::ReynoldsStress<BasicTurbulenceModel>::divDevRhoReff "laplacian(nuEff,U)" ) - fvm::laplacian(this->alpha_*this->rho_*this->nuEff(), U) - - fvc::div(this->alpha_*this->rho_*nu()*dev2(T(fvc::grad(U)))) + - fvc::div(this->alpha_*this->rho_*this->nu()*dev2(T(fvc::grad(U)))) ); } diff --git a/src/TurbulenceModels/turbulenceModels/ReynoldsStress/ReynoldsStress.H b/src/TurbulenceModels/turbulenceModels/ReynoldsStress/ReynoldsStress.H index 9c27cc4c81d73dd8c3bb24342c0a7474de0d3288..b27468b5a335cadcf34b6c59d2533c8394759795 100644 --- a/src/TurbulenceModels/turbulenceModels/ReynoldsStress/ReynoldsStress.H +++ b/src/TurbulenceModels/turbulenceModels/ReynoldsStress/ReynoldsStress.H @@ -69,6 +69,10 @@ protected: // Protected Member Functions + void boundNormalStress(volSymmTensorField& R) const; + void correctWallShearStress(volSymmTensorField& R) const; + + //- Update the eddy-viscosity virtual void correctNut() = 0;