Skip to content
Snippets Groups Projects
Commit 9106a4f5 authored by andy's avatar andy
Browse files

ENH: Updated SpalartAllmaras family of DES models

parent 1d009cfb
Branches
Tags
No related merge requests found
......@@ -118,7 +118,7 @@ SpalartAllmaras::SpalartAllmaras
const word& modelName
)
:
LESModel(modelName, rho, U, phi, thermoPhysicalModel, turbulenceModelName),
DESModel(modelName, rho, U, phi, thermoPhysicalModel, turbulenceModelName),
sigmaNut_
(
......@@ -364,6 +364,32 @@ bool SpalartAllmaras::read()
}
tmp<volScalarField> SpalartAllmaras::LESRegion() const
{
volScalarField wd(wallDist(mesh_).y());
tmp<volScalarField> tLESRegion
(
new volScalarField
(
IOobject
(
"DES::LESRegion",
mesh_.time().timeName(),
mesh_,
IOobject::NO_READ,
IOobject::NO_WRITE
),
neg(min(CDES_*delta(), wd) - wd)
// mesh_,
// dimensionedScalar("zero", dimless, 0.0)
)
);
return tLESRegion;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace LESModels
......
......@@ -25,7 +25,7 @@ Class
Foam::compressible::LESModels::SpalartAllmaras
Group
grpCmpLESTurbulence
grpCmpDESTurbulence
Description
SpalartAllmaras for compressible flows
......@@ -38,7 +38,7 @@ SourceFiles
#ifndef compressibleSpalartAllmaras_H
#define compressibleSpalartAllmaras_H
#include "LESModel.H"
#include "DESModel.H"
#include "volFields.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
......@@ -56,7 +56,7 @@ namespace LESModels
class SpalartAllmaras
:
public LESModel
public DESModel
{
// Private data
......@@ -169,6 +169,9 @@ public:
//- Read LESProperties dictionary
virtual bool read();
//- Return the LES field indicator
virtual tmp<volScalarField> LESRegion() const;
};
......
......@@ -152,7 +152,7 @@ SpalartAllmaras::SpalartAllmaras
const word& modelName
)
:
LESModel(modelName, U, phi, transport, turbulenceModelName),
DESModel(modelName, U, phi, transport, turbulenceModelName),
sigmaNut_
(
......@@ -395,6 +395,30 @@ bool SpalartAllmaras::read()
}
tmp<volScalarField> SpalartAllmaras::LESRegion() const
{
tmp<volScalarField> tLESRegion
(
new volScalarField
(
IOobject
(
"DES::LESRegion",
mesh_.time().timeName(),
mesh_,
IOobject::NO_READ,
IOobject::NO_WRITE
),
neg(dTilda(S(fvc::grad(U_))) - y_)
// mesh_,
// dimensionedScalar("zero", dimless, 0.0)
)
);
return tLESRegion;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace LESModels
......
......@@ -25,7 +25,7 @@ Class
Foam::incompressible::LESModels::SpalartAllmaras
Group
grpIcoLESTurbulence
grpIcoDESTurbulence
Description
SpalartAllmaras DES (SA + LES) turbulence model for incompressible flows
......@@ -38,7 +38,7 @@ SourceFiles
#ifndef SpalartAllmaras_H
#define SpalartAllmaras_H
#include "LESModel.H"
#include "DESModel.H"
#include "volFields.H"
#include "wallDist.H"
......@@ -57,7 +57,7 @@ namespace LESModels
class SpalartAllmaras
:
public LESModel
public DESModel
{
// Private Member Functions
......@@ -194,6 +194,9 @@ public:
//- Read LESProperties dictionary
virtual bool read();
//- Return the LES field indicator
virtual tmp<volScalarField> LESRegion() const;
};
......
......@@ -25,7 +25,7 @@ Class
Foam::incompressible::LESModels::SpalartAllmarasDDES
Group
grpIcoLESTurbulence
grpIcoDESTurbulence
Description
SpalartAllmaras DDES LES turbulence model for incompressible flows
......@@ -55,7 +55,7 @@ namespace LESModels
{
/*---------------------------------------------------------------------------*\
Class SpalartAllmarasDDES Declaration
Class SpalartAllmarasDDES Declaration
\*---------------------------------------------------------------------------*/
class SpalartAllmarasDDES
......
......@@ -25,7 +25,7 @@ Class
Foam::incompressible::LESModels::SpalartAllmarasIDDES
Group
grpIcoLESTurbulence
grpIcoDESTurbulence
Description
SpalartAllmarasIDDES LES turbulence model for incompressible flows
......
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