Commit c10dec17 authored by Henry's avatar Henry
Browse files

turbulenceModels/RAS/SSG: Added Speziale, Sarkar and Gatski pressure-strain...

turbulenceModels/RAS/SSG: Added Speziale, Sarkar and Gatski pressure-strain based Reynolds-stress turbulence model
parent c61529d0
......@@ -81,6 +81,9 @@ makeRASModel(v2f);
#include "LRR.H"
makeRASModel(LRR);
#include "SSG.H"
makeRASModel(SSG);
#include "Smagorinsky.H"
makeLESModel(Smagorinsky);
......
......@@ -76,6 +76,9 @@ makeRASModel(v2f);
#include "LRR.H"
makeRASModel(LRR);
#include "SSG.H"
makeRASModel(SSG);
#include "Smagorinsky.H"
makeLESModel(Smagorinsky);
......
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 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 "SSG.H"
#include "wallFvPatch.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
namespace RASModels
{
// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
template<class BasicTurbulenceModel>
void SSG<BasicTurbulenceModel>::correctNut()
{
this->nut_ = this->Cmu_*sqr(k_)/epsilon_;
this->nut_.correctBoundaryConditions();
BasicTurbulenceModel::correctNut();
}
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
template<class BasicTurbulenceModel>
SSG<BasicTurbulenceModel>::SSG
(
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<RASModel<BasicTurbulenceModel> >
(
type,
alpha,
rho,
U,
alphaRhoPhi,
phi,
transport,
propertiesName
),
Cmu_
(
dimensioned<scalar>::lookupOrAddToDict
(
"Cmu",
this->coeffDict_,
0.09
)
),
C1_
(
dimensioned<scalar>::lookupOrAddToDict
(
"C1",
this->coeffDict_,
3.4
)
),
C1s_
(
dimensioned<scalar>::lookupOrAddToDict
(
"C1s",
this->coeffDict_,
1.8
)
),
C2_
(
dimensioned<scalar>::lookupOrAddToDict
(
"C2",
this->coeffDict_,
4.2
)
),
C3_
(
dimensioned<scalar>::lookupOrAddToDict
(
"C3",
this->coeffDict_,
0.8
)
),
C3s_
(
dimensioned<scalar>::lookupOrAddToDict
(
"C3s",
this->coeffDict_,
1.3
)
),
C4_
(
dimensioned<scalar>::lookupOrAddToDict
(
"C4",
this->coeffDict_,
1.25
)
),
C5_
(
dimensioned<scalar>::lookupOrAddToDict
(
"C5",
this->coeffDict_,
0.4
)
),
Ceps1_
(
dimensioned<scalar>::lookupOrAddToDict
(
"Ceps1",
this->coeffDict_,
1.44
)
),
Ceps2_
(
dimensioned<scalar>::lookupOrAddToDict
(
"Ceps2",
this->coeffDict_,
1.92
)
),
Cs_
(
dimensioned<scalar>::lookupOrAddToDict
(
"Cs",
this->coeffDict_,
0.25
)
),
Ceps_
(
dimensioned<scalar>::lookupOrAddToDict
(
"Ceps",
this->coeffDict_,
0.15
)
),
k_
(
IOobject
(
"k",
this->runTime_.timeName(),
this->mesh_,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
0.5*tr(this->R_)
),
epsilon_
(
IOobject
(
"epsilon",
this->runTime_.timeName(),
this->mesh_,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
this->mesh_
)
{
if (type == typeName)
{
this->boundNormalStress(this->R_);
bound(epsilon_, this->epsilonMin_);
k_ = 0.5*tr(this->R_);
correctNut();
this->printCoeffs(type);
}
}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
template<class BasicTurbulenceModel>
bool SSG<BasicTurbulenceModel>::read()
{
if (ReynoldsStress<RASModel<BasicTurbulenceModel> >::read())
{
Cmu_.readIfPresent(this->coeffDict());
C1_.readIfPresent(this->coeffDict());
C1s_.readIfPresent(this->coeffDict());
C2_.readIfPresent(this->coeffDict());
C3_.readIfPresent(this->coeffDict());
C3s_.readIfPresent(this->coeffDict());
C4_.readIfPresent(this->coeffDict());
C5_.readIfPresent(this->coeffDict());
Ceps1_.readIfPresent(this->coeffDict());
Ceps2_.readIfPresent(this->coeffDict());
Cs_.readIfPresent(this->coeffDict());
Ceps_.readIfPresent(this->coeffDict());
return true;
}
else
{
return false;
}
}
template<class BasicTurbulenceModel>
void SSG<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<RASModel<BasicTurbulenceModel> >::correct();
tmp<volTensorField> tgradU(fvc::grad(U));
const volTensorField& gradU = tgradU();
volSymmTensorField P(-twoSymm(R & gradU));
volScalarField G(this->GName(), 0.5*mag(tr(P)));
// Update epsilon and G at the wall
epsilon_.boundaryField().updateCoeffs();
// Dissipation equation
tmp<fvScalarMatrix> epsEqn
(
fvm::ddt(alpha, rho, epsilon_)
+ fvm::div(alphaRhoPhi, epsilon_)
- fvm::laplacian(Ceps_*alpha*rho*(k_/epsilon_)*R, epsilon_)
==
Ceps1_*alpha*rho*G*epsilon_/k_
- fvm::Sp(Ceps2_*alpha*rho*epsilon_/k_, epsilon_)
);
epsEqn().relax();
epsEqn().boundaryManipulate(epsilon_.boundaryField());
solve(epsEqn);
bound(epsilon_, this->epsilonMin_);
// Correct the trace of the tensorial production to be consistent
// with the near-wall generation from the wall-functions
const fvPatchList& patches = this->mesh_.boundary();
forAll(patches, patchi)
{
const fvPatch& curPatch = patches[patchi];
if (isA<wallFvPatch>(curPatch))
{
forAll(curPatch, facei)
{
label faceCelli = curPatch.faceCells()[facei];
P[faceCelli] *= min
(
G[faceCelli]/(0.5*mag(tr(P[faceCelli])) + SMALL),
1.0
);
}
}
}
volSymmTensorField b(dev(R)/(2*k_));
volSymmTensorField S(symm(gradU));
volTensorField Omega(skew(gradU));
// Reynolds stress equation
tmp<fvSymmTensorMatrix> REqn
(
fvm::ddt(alpha, rho, R)
+ fvm::div(alphaRhoPhi, R)
- fvm::laplacian(Cs_*alpha*rho*(k_/epsilon_)*R, R)
+ fvm::Sp(((C1_/2)*epsilon_ + (C1s_/2)*G)*alpha*rho/k_, R)
==
alpha*rho*P
- ((1.0/3.0)*I)*(((2.0 - C1_)*epsilon_ - C1s_*G)*alpha*rho)
+ (C2_*(alpha*rho*epsilon_))*dev(symm(b&b)) // symm should not be needed
+ alpha*rho*k_
*(
(C3_ - C3s_*mag(b))*S
+ C4_*dev(twoSymm(b&S))
+ C5_*twoSymm(b&Omega)
)
);
REqn().relax();
solve(REqn);
this->boundNormalStress(R);
k_ = 0.5*tr(R);
correctNut();
// Correct wall shear-stresses when applying wall-functions
this->correctWallShearStress(R);
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace RASModels
} // End namespace Foam
// ************************************************************************* //
/*---------------------------------------------------------------------------*\
========= |
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
\\ / O peration |
\\ / A nd | Copyright (C) 2015 OpenFOAM Foundation
\\/ M anipulation |
-------------------------------------------------------------------------------
License
This file is part of OpenFOAM.
OpenFOAM is free software: you can redistribute it and/or modify it
under the terms of the GNU General Public License as published by
the Free Software Foundation, either version 3 of the License, or
(at your option) any later version.
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
for more details.
You should have received a copy of the GNU General Public License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
Class
Foam::RASModels::SSG
Group
grpRASTurbulence
Description
Speziale, Sarkar and Gatski Reynolds-stress turbulence model for
incompressible and compressible flows.
Reference:
\verbatim
Speziale, C. G., Sarkar, S., & Gatski, T. B. (1991).
Modelling the pressure–strain correlation of turbulence:
an invariant dynamical systems approach.
Journal of Fluid Mechanics, 227, 245-272.
\endverbatim
Including the generalized gradient diffusion model of
Daly and Harlow:
\verbatim
Daly, B. J., & Harlow, F. H. (1970).
Transport equations in turbulence.
Physics of Fluids (1958-1988), 13(11), 2634-2649.
\endverbatim
The default model coefficients are:
\verbatim
SSGCoeffs
{
Cmu 0.09;
C1 3.4;
C1s 1.8;
C2 4.2;
C3 0.8;
C3s 1.3;
C4 1.25;
C5 0.4;
Ceps1 1.44;
Ceps2 1.92;
Cs 0.25;
Ceps 0.15;
couplingFactor 0.0;
}
\endverbatim
SourceFiles
SSG.C
\*---------------------------------------------------------------------------*/
#ifndef SSG_H
#define SSG_H
#include "RASModel.H"
#include "ReynoldsStress.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
namespace RASModels
{
/*---------------------------------------------------------------------------*\
Class SSG Declaration
\*---------------------------------------------------------------------------*/
template<class BasicTurbulenceModel>
class SSG
:
public ReynoldsStress<RASModel<BasicTurbulenceModel> >
{
// Private Member Functions
// Disallow default bitwise copy construct and assignment
SSG(const SSG&);
SSG& operator=(const SSG&);
protected:
// Protected data
// Model coefficients
dimensionedScalar Cmu_;
dimensionedScalar C1_;
dimensionedScalar C1s_;
dimensionedScalar C2_;
dimensionedScalar C3_;
dimensionedScalar C3s_;
dimensionedScalar C4_;
dimensionedScalar C5_;
dimensionedScalar Ceps1_;
dimensionedScalar Ceps2_;
dimensionedScalar Cs_;
dimensionedScalar Ceps_;
// Fields
volScalarField k_;
volScalarField epsilon_;
// 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("SSG");
// Constructors
//- Construct from components
SSG
(
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 ~SSG()
{}
// Member Functions
//- Read model coefficients if they have changed
virtual bool read();
//- Return the turbulence kinetic energy
virtual tmp<volScalarField> k() const
{
return k_;
}
//- Return the turbulence kinetic energy dissipation rate
virtual tmp<volScalarField> epsilon() const
{
return epsilon_;
}
//- Solve the turbulence equations and correct eddy-Viscosity and
// related properties
virtual void correct();
};
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace RASModels
} // End namespace Foam
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#ifdef NoRepository
# include "SSG.C"
#endif
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
#endif
// ************************************************************************* //
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment