Commit 6b42dbae authored by Henry Weller's avatar Henry Weller
Browse files

TurbulenceModels: Created a general base-class and selection mechanism for laminar stress models

Renamed the original 'laminar' model to 'Stokes' to indicate it is a
linear stress model supporting both Newtonian and non-Newtonian
viscosity.

This general framework will support linear, non-linear, visco-elastic
etc. laminar transport models.

For backward compatibility the 'Stokes' laminar stress model can be
selected either the original 'laminar' 'simulationType'
specification in turbulenceProperties:

    simulationType laminar;

or using the new more general 'laminarModel' specification:

    simulationType laminar;

    laminar
    {
        laminarModel        Stokes;
    }

which allows other laminar stress models to be selected.
parent 49d3b6ea
......@@ -28,8 +28,8 @@ License
#include "addToRunTimeSelectionTable.H"
#include "makeTurbulenceModel.H"
#include "laminar.H"
#include "turbulentTransportModel.H"
#include "laminarModel.H"
#include "RASModel.H"
#include "LESModel.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
......@@ -52,6 +52,10 @@ makeBaseTurbulenceModel
singlePhaseTransportModel
);
#define makeLaminarModel(Type) \
makeTemplatedTurbulenceModel \
(singlePhaseTransportModelPhaseIncompressibleTurbulenceModel, laminar, Type)
#define makeRASModel(Type) \
makeTemplatedTurbulenceModel \
(singlePhaseTransportModelPhaseIncompressibleTurbulenceModel, RAS, Type)
......@@ -60,6 +64,9 @@ makeBaseTurbulenceModel
makeTemplatedTurbulenceModel \
(singlePhaseTransportModelPhaseIncompressibleTurbulenceModel, LES, Type)
#include "Stokes.H"
makeLaminarModel(Stokes);
#include "kEpsilon.H"
makeRASModel(kEpsilon);
......
......@@ -28,7 +28,7 @@ License
#include "addToRunTimeSelectionTable.H"
#include "makeTurbulenceModel.H"
#include "laminar.H"
#include "laminarModel.H"
#include "RASModel.H"
#include "LESModel.H"
......@@ -52,6 +52,14 @@ makeBaseTurbulenceModel
incompressibleTwoPhaseInteractingMixture
);
#define makeLaminarModel(Type) \
makeTemplatedTurbulenceModel \
( \
incompressibleTwoPhaseInteractingMixtureCompressibleTurbulenceModel, \
laminar, \
Type \
)
#define makeRASModel(Type) \
makeTemplatedTurbulenceModel \
( \
......@@ -68,6 +76,9 @@ makeBaseTurbulenceModel
Type \
)
#include "Stokes.H"
makeLaminarModel(Stokes);
#include "kEpsilon.H"
makeRASModel(kEpsilon);
......
......@@ -27,7 +27,7 @@ License
#include "addToRunTimeSelectionTable.H"
#include "makeTurbulenceModel.H"
#include "laminar.H"
#include "laminarModel.H"
#include "RASModel.H"
#include "LESModel.H"
......@@ -53,6 +53,10 @@ makeBaseTurbulenceModel
phaseModel
);
#define makeLaminarModel(Type) \
makeTemplatedLaminarModel \
(phaseModelPhaseCompressibleTurbulenceModel, laminar, Type)
#define makeRASModel(Type) \
makeTemplatedTurbulenceModel \
(phaseModelPhaseCompressibleTurbulenceModel, RAS, Type)
......@@ -61,6 +65,9 @@ makeBaseTurbulenceModel
makeTemplatedTurbulenceModel \
(phaseModelPhaseCompressibleTurbulenceModel, LES, Type)
#include "Stokes.H"
makeLaminarModel(Stokes);
#include "kEpsilon.H"
makeRASModel(kEpsilon);
......
......@@ -27,7 +27,7 @@ License
#include "addToRunTimeSelectionTable.H"
#include "makeTurbulenceModel.H"
#include "laminar.H"
#include "laminarModel.H"
#include "RASModel.H"
#include "LESModel.H"
......@@ -53,6 +53,10 @@ makeBaseTurbulenceModel
phaseModel
);
#define makeLaminarModel(Type) \
makeTemplatedLaminarModel \
(phaseModelPhaseCompressibleTurbulenceModel, laminar, Type)
#define makeRASModel(Type) \
makeTemplatedTurbulenceModel \
(phaseModelPhaseCompressibleTurbulenceModel, RAS, Type)
......@@ -61,6 +65,9 @@ makeBaseTurbulenceModel
makeTemplatedTurbulenceModel \
(phaseModelPhaseCompressibleTurbulenceModel, LES, Type)
#include "Stokes.H"
makeLaminarModel(Stokes);
#include "kEpsilon.H"
makeRASModel(kEpsilon);
......
......@@ -32,7 +32,7 @@ License
#include "ThermalDiffusivity.H"
#include "EddyDiffusivity.H"
#include "laminar.H"
#include "laminarModel.H"
#include "RASModel.H"
#include "LESModel.H"
......@@ -58,6 +58,10 @@ makeBaseTurbulenceModel
phaseModel
);
#define makeLaminarModel(Type) \
makeTemplatedLaminarModel \
(phaseModelPhaseCompressibleTurbulenceModel, laminar, Type)
#define makeRASModel(Type) \
makeTemplatedTurbulenceModel \
(phaseModelPhaseCompressibleTurbulenceModel, RAS, Type)
......@@ -66,6 +70,9 @@ makeBaseTurbulenceModel
makeTemplatedTurbulenceModel \
(phaseModelPhaseCompressibleTurbulenceModel, LES, Type)
#include "Stokes.H"
makeLaminarModel(Stokes);
#include "kEpsilon.H"
makeRASModel(kEpsilon);
......
......@@ -34,7 +34,8 @@ License
{ \
typedef TDModel<BaseModel<Transport>> \
Transport##BaseModel; \
typedef laminar<Transport##BaseModel> Laminar##Transport##BaseModel; \
typedef laminarModel<Transport##BaseModel> \
laminar##Transport##BaseModel; \
typedef RASModel<EddyDiffusivity<Transport##BaseModel>> \
RAS##Transport##BaseModel; \
typedef LESModel<EddyDiffusivity<Transport##BaseModel>> \
......@@ -65,14 +66,18 @@ License
Transport##BaseModel; \
\
\
typedef laminar<Transport##BaseModel> Laminar##Transport##BaseModel; \
typedef laminarModel<Transport##BaseModel> \
laminar##Transport##BaseModel; \
\
defineNamedTemplateTypeNameAndDebug(Laminar##Transport##BaseModel, 0); \
defineNamedTemplateTypeNameAndDebug(laminar##Transport##BaseModel, 0); \
\
defineTemplateRunTimeSelectionTable \
(laminar##Transport##BaseModel, dictionary); \
\
addToRunTimeSelectionTable \
( \
Transport##baseModel, \
Laminar##Transport##BaseModel, \
laminar##Transport##BaseModel, \
dictionary \
); \
\
......@@ -110,6 +115,27 @@ License
}
#define makeTemplatedLaminarModel(BaseModel, SType, Type) \
typedef Foam::SType##Models::Type<Foam::BaseModel> \
Type##SType##BaseModel; \
defineNamedTemplateTypeNameAndDebug(Type##SType##BaseModel, 0); \
\
namespace Foam \
{ \
namespace SType##Models \
{ \
typedef Type<BaseModel> Type##SType##BaseModel; \
\
addToRunTimeSelectionTable \
( \
SType##BaseModel, \
Type##SType##BaseModel, \
dictionary \
); \
} \
}
#define makeTemplatedTurbulenceModel(BaseModel, SType, Type) \
typedef Foam::SType##Models::Type<Foam::EddyDiffusivity<Foam::BaseModel>> \
Type##SType##BaseModel; \
......
......@@ -46,6 +46,7 @@ SourceFiles
#include "CompressibleTurbulenceModel.H"
#include "ThermalDiffusivity.H"
#include "EddyDiffusivity.H"
#include "laminarModel.H"
#include "RASModel.H"
#include "LESModel.H"
#include "fluidThermo.H"
......@@ -59,6 +60,7 @@ namespace Foam
typedef ThermalDiffusivity<CompressibleTurbulenceModel<fluidThermo>>
turbulenceModel;
typedef laminarModel<turbulenceModel> laminarModel;
typedef RASModel<EddyDiffusivity<turbulenceModel>> RASModel;
typedef LESModel<EddyDiffusivity<turbulenceModel>> LESModel;
......
......@@ -37,6 +37,15 @@ makeBaseTurbulenceModel
fluidThermo
);
// -------------------------------------------------------------------------- //
// Laminar models
// -------------------------------------------------------------------------- //
#include "Stokes.H"
makeLaminarModel(Stokes);
// -------------------------------------------------------------------------- //
// RAS models
// -------------------------------------------------------------------------- //
......
......@@ -32,7 +32,7 @@ License
#include "ThermalDiffusivity.H"
#include "EddyDiffusivity.H"
#include "laminar.H"
#include "laminarModel.H"
#include "RASModel.H"
#include "LESModel.H"
......@@ -48,6 +48,10 @@ makeTurbulenceModelTypes
fluidThermo
);
#define makeLaminarModel(Type) \
makeTemplatedLaminarModel \
(fluidThermoCompressibleTurbulenceModel, laminar, Type)
#define makeRASModel(Type) \
makeTemplatedTurbulenceModel \
(fluidThermoCompressibleTurbulenceModel, RAS, Type)
......
......@@ -44,6 +44,7 @@ SourceFiles
#define turbulentTransportModel_H
#include "IncompressibleTurbulenceModel.H"
#include "laminarModel.H"
#include "RASModel.H"
#include "LESModel.H"
#include "incompressible/transportModel/transportModel.H"
......@@ -56,6 +57,7 @@ namespace Foam
{
typedef IncompressibleTurbulenceModel<transportModel> turbulenceModel;
typedef laminarModel<turbulenceModel> laminarModel;
typedef RASModel<turbulenceModel> RASModel;
typedef LESModel<turbulenceModel> LESModel;
......
......@@ -36,6 +36,15 @@ makeBaseTurbulenceModel
transportModel
);
// -------------------------------------------------------------------------- //
// Laminar models
// -------------------------------------------------------------------------- //
#include "Stokes.H"
makeLaminarModel(Stokes);
// -------------------------------------------------------------------------- //
// RAS models
// -------------------------------------------------------------------------- //
......
......@@ -28,7 +28,7 @@ License
#include "addToRunTimeSelectionTable.H"
#include "makeTurbulenceModel.H"
#include "laminar.H"
#include "laminarModel.H"
#include "RASModel.H"
#include "LESModel.H"
......@@ -43,6 +43,10 @@ makeTurbulenceModelTypes
transportModel
);
#define makeLaminarModel(Type) \
makeTemplatedTurbulenceModel \
(transportModelIncompressibleTurbulenceModel, laminar, Type)
#define makeRASModel(Type) \
makeTemplatedTurbulenceModel \
(transportModelIncompressibleTurbulenceModel, RAS, Type)
......
......@@ -23,17 +23,24 @@ License
\*---------------------------------------------------------------------------*/
#include "laminar.H"
#include "Stokes.H"
#include "volFields.H"
#include "surfaceFields.H"
#include "fvcGrad.H"
#include "fvcDiv.H"
#include "fvmLaplacian.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
namespace laminarModels
{
// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
template<class BasicTurbulenceModel>
Foam::laminar<BasicTurbulenceModel>::laminar
Stokes<BasicTurbulenceModel>::Stokes
(
const alphaField& alpha,
const rhoField& rho,
......@@ -44,7 +51,7 @@ Foam::laminar<BasicTurbulenceModel>::laminar
const word& propertiesName
)
:
linearViscousStress<BasicTurbulenceModel>
linearViscousStress<laminarModel<BasicTurbulenceModel>>
(
typeName,
alpha,
......@@ -58,57 +65,26 @@ Foam::laminar<BasicTurbulenceModel>::laminar
{}
// * * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * //
template<class BasicTurbulenceModel>
Foam::autoPtr<Foam::laminar<BasicTurbulenceModel>>
Foam::laminar<BasicTurbulenceModel>::New
(
const alphaField& alpha,
const rhoField& rho,
const volVectorField& U,
const surfaceScalarField& alphaRhoPhi,
const surfaceScalarField& phi,
const transportModel& transport,
const word& propertiesName
)
{
return autoPtr<laminar>
(
new laminar
(
alpha,
rho,
U,
alphaRhoPhi,
phi,
transport,
propertiesName
)
);
}
// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
template<class BasicTurbulenceModel>
const Foam::dictionary&
Foam::laminar<BasicTurbulenceModel>::coeffDict() const
const dictionary&
Stokes<BasicTurbulenceModel>::coeffDict() const
{
return dictionary::null;
}
template<class BasicTurbulenceModel>
bool Foam::laminar<BasicTurbulenceModel>::read()
bool Stokes<BasicTurbulenceModel>::read()
{
return true;
}
template<class BasicTurbulenceModel>
Foam::tmp<Foam::volScalarField>
Foam::laminar<BasicTurbulenceModel>::nut() const
tmp<volScalarField>
Stokes<BasicTurbulenceModel>::nut() const
{
return tmp<volScalarField>
(
......@@ -131,8 +107,8 @@ Foam::laminar<BasicTurbulenceModel>::nut() const
template<class BasicTurbulenceModel>
Foam::tmp<Foam::scalarField>
Foam::laminar<BasicTurbulenceModel>::nut
tmp<scalarField>
Stokes<BasicTurbulenceModel>::nut
(
const label patchi
) const
......@@ -145,8 +121,8 @@ Foam::laminar<BasicTurbulenceModel>::nut
template<class BasicTurbulenceModel>
Foam::tmp<Foam::volScalarField>
Foam::laminar<BasicTurbulenceModel>::nuEff() const
tmp<volScalarField>
Stokes<BasicTurbulenceModel>::nuEff() const
{
return tmp<volScalarField>
(
......@@ -159,8 +135,8 @@ Foam::laminar<BasicTurbulenceModel>::nuEff() const
template<class BasicTurbulenceModel>
Foam::tmp<Foam::scalarField>
Foam::laminar<BasicTurbulenceModel>::nuEff
tmp<scalarField>
Stokes<BasicTurbulenceModel>::nuEff
(
const label patchi
) const
......@@ -170,8 +146,8 @@ Foam::laminar<BasicTurbulenceModel>::nuEff
template<class BasicTurbulenceModel>
Foam::tmp<Foam::volScalarField>
Foam::laminar<BasicTurbulenceModel>::k() const
tmp<volScalarField>
Stokes<BasicTurbulenceModel>::k() const
{
return tmp<volScalarField>
(
......@@ -194,8 +170,8 @@ Foam::laminar<BasicTurbulenceModel>::k() const
template<class BasicTurbulenceModel>
Foam::tmp<Foam::volScalarField>
Foam::laminar<BasicTurbulenceModel>::epsilon() const
tmp<volScalarField>
Stokes<BasicTurbulenceModel>::epsilon() const
{
return tmp<volScalarField>
(
......@@ -221,8 +197,8 @@ Foam::laminar<BasicTurbulenceModel>::epsilon() const
template<class BasicTurbulenceModel>
Foam::tmp<Foam::volSymmTensorField>
Foam::laminar<BasicTurbulenceModel>::R() const
tmp<volSymmTensorField>
Stokes<BasicTurbulenceModel>::R() const
{
return tmp<volSymmTensorField>
(
......@@ -248,10 +224,15 @@ Foam::laminar<BasicTurbulenceModel>::R() const
template<class BasicTurbulenceModel>
void Foam::laminar<BasicTurbulenceModel>::correct()
void Stokes<BasicTurbulenceModel>::correct()
{
BasicTurbulenceModel::correct();
laminarModel<BasicTurbulenceModel>::correct();
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
} // End namespace laminarModels
} // End namespace Foam
// ************************************************************************* //
......@@ -22,34 +22,37 @@ License
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
Class
Foam::laminar
Foam::Stokes
Description
Turbulence model for laminar flow.
Turbulence model for Stokes flow.
SourceFiles
laminar.C
Stokes.C
\*---------------------------------------------------------------------------*/
#ifndef laminar_H
#define laminar_H
#ifndef Stokes_H
#define Stokes_H
#include "laminarModel.H"
#include "linearViscousStress.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
namespace Foam
{
namespace laminarModels
{
/*---------------------------------------------------------------------------* \
Class laminar Declaration
Class Stokes Declaration
\*---------------------------------------------------------------------------*/
template<class BasicTurbulenceModel>
class laminar
class Stokes
:
public linearViscousStress<BasicTurbulenceModel>
public linearViscousStress<laminarModel<BasicTurbulenceModel>>
{
public:
......@@ -60,13 +63,13 @@ public:
//- Runtime type information
TypeName("laminar");
TypeName("Stokes");
// Constructors
//- Construct from components
laminar
Stokes
(
const alphaField& alpha,
const rhoField& rho,
......@@ -81,7 +84,7 @@ public:
// Selectors
//- Return a reference to the selected turbulence model
static autoPtr<laminar> New
static autoPtr<Stokes> New
(
const alphaField& alpha,
const rhoField& rho,
......@@ -94,7 +97,7 @@ public:
//- Destructor
virtual ~laminar()
virtual ~Stokes()
{}
......@@ -106,41 +109,42 @@ public:
//- Read turbulenceProperties dictionary
virtual bool read();
//- Return the turbulence viscosity, i.e. 0 for laminar flow
//- Return the turbulence viscosity, i.e. 0 for Stokes flow
virtual tmp<volScalarField> nut() const;
//- Return the turbulence viscosity on patch
virtual tmp<scalarField> nut(const label patchi) const;
//- Return the effective viscosity, i.e. the laminar viscosity