Stokes.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2013-2020 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 Class
25  Foam::laminarModels::Stokes
26 
27 Description
28  Turbulence model for Stokes flow.
29 
30 SourceFiles
31  Stokes.C
32 
33 \*---------------------------------------------------------------------------*/
34 
35 #ifndef Stokes_H
36 #define Stokes_H
37 
38 #include "laminarModel.H"
39 #include "linearViscousStress.H"
40 
41 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
42 
43 namespace Foam
44 {
45 namespace laminarModels
46 {
47 
48 /*---------------------------------------------------------------------------* \
49  Class Stokes Declaration
50 \*---------------------------------------------------------------------------*/
51 
52 template<class BasicMomentumTransportModel>
53 class Stokes
54 :
55  public linearViscousStress<laminarModel<BasicMomentumTransportModel>>
56 {
57 
58 public:
59 
60  typedef typename BasicMomentumTransportModel::alphaField alphaField;
61  typedef typename BasicMomentumTransportModel::rhoField rhoField;
62  typedef typename BasicMomentumTransportModel::transportModel transportModel;
63 
64 
65  //- Runtime type information
66  TypeName("Stokes");
67 
68 
69  // Constructors
70 
71  //- Construct from components
72  Stokes
73  (
74  const alphaField& alpha,
75  const rhoField& rho,
76  const volVectorField& U,
77  const surfaceScalarField& alphaRhoPhi,
78  const surfaceScalarField& phi,
79  const transportModel& transport
80  );
81 
82 
83  // Selectors
84 
85  //- Return a reference to the selected turbulence model
86  static autoPtr<Stokes> New
87  (
88  const alphaField& alpha,
89  const rhoField& rho,
90  const volVectorField& U,
91  const surfaceScalarField& alphaRhoPhi,
92  const surfaceScalarField& phi,
93  const transportModel& transport
94  );
95 
96 
97  //- Destructor
98  virtual ~Stokes()
99  {}
100 
101 
102  // Member Functions
103 
104  //- Const access to the coefficients dictionary
105  virtual const dictionary& coeffDict() const;
106 
107  //- Read momentumTransport dictionary
108  virtual bool read();
109 
110  //- Return the turbulence viscosity, i.e. 0 for Stokes flow
111  virtual tmp<volScalarField> nut() const;
112 
113  //- Return the turbulence viscosity on patch
114  virtual tmp<scalarField> nut(const label patchi) const;
115 
116  //- Return the effective viscosity, i.e. the Stokes viscosity
117  virtual tmp<volScalarField> nuEff() const;
118 
119  //- Return the effective viscosity on patch
120  virtual tmp<scalarField> nuEff(const label patchi) const;
121 
122  //- Return the turbulence kinetic energy, i.e. 0 for Stokes flow
123  virtual tmp<volScalarField> k() const;
124 
125  //- Return the turbulence kinetic energy dissipation rate,
126  // i.e. 0 for Stokes flow
127  virtual tmp<volScalarField> epsilon() const;
128 
129  //- Return the stress tensor [m^2/s^2], i.e. 0 for Stokes flow
130  virtual tmp<volSymmTensorField> sigma() const;
131 
132  //- Correct the Stokes viscosity
133  virtual void correct();
134 };
135 
136 
137 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
138 
139 } // End namespace laminarModels
140 } // End namespace Foam
141 
142 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
143 
144 #ifdef NoRepository
145  #include "Stokes.C"
146 #endif
147 
148 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
149 
150 #endif
151 
152 // ************************************************************************* //
virtual tmp< volSymmTensorField > sigma() const
Return the stress tensor [m^2/s^2], i.e. 0 for Stokes flow.
Definition: Stokes.C:161
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
virtual tmp< volScalarField > nuEff() const
Return the effective viscosity, i.e. the Stokes viscosity.
Definition: Stokes.C:112
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
virtual bool read()
Read momentumTransport dictionary.
Definition: Stokes.C:77
Linear viscous stress turbulence model base class.
virtual const dictionary & coeffDict() const
Const access to the coefficients dictionary.
Definition: Stokes.C:70
virtual tmp< volScalarField > epsilon() const
Return the turbulence kinetic energy dissipation rate,.
Definition: Stokes.C:148
BasicMomentumTransportModel::rhoField rhoField
Definition: Stokes.H:60
phi
Definition: pEqn.H:104
Stokes(const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport)
Construct from components.
Definition: Stokes.C:44
TypeName("Stokes")
Runtime type information.
BasicMomentumTransportModel::transportModel transportModel
Definition: Stokes.H:61
virtual ~Stokes()
Destructor.
Definition: Stokes.H:97
BasicMomentumTransportModel::alphaField alphaField
Definition: Stokes.H:59
label patchi
U
Definition: pEqn.H:72
static autoPtr< Stokes > New(const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport)
Return a reference to the selected turbulence model.
virtual tmp< volScalarField > k() const
Return the turbulence kinetic energy, i.e. 0 for Stokes flow.
Definition: Stokes.C:135
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
Turbulence model for Stokes flow.
Definition: Stokes.H:52
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
virtual void correct()
Correct the Stokes viscosity.
Definition: Stokes.C:173
A class for managing temporary objects.
Definition: PtrList.H:53
virtual tmp< volScalarField > nut() const
Return the turbulence viscosity, i.e. 0 for Stokes flow.
Definition: Stokes.C:85
Namespace for OpenFOAM.