kEqn.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) 2011-2021 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::LESModels::kEqn
26 
27 Description
28  One equation eddy-viscosity model
29 
30  Eddy viscosity SGS model using a modeled balance equation to simulate the
31  behaviour of k.
32 
33  Reference:
34  \verbatim
35  Yoshizawa, A. (1986).
36  Statistical theory for compressible turbulent shear flows,
37  with the application to subgrid modeling.
38  Physics of Fluids (1958-1988), 29(7), 2152-2164.
39  \endverbatim
40 
41  The default model coefficients are
42  \verbatim
43  kEqnCoeffs
44  {
45  Ck 0.094;
46  Ce 1.048;
47  }
48  \endverbatim
49 
50 SourceFiles
51  kEqn.C
52 
53 \*---------------------------------------------------------------------------*/
54 
55 #ifndef kEqn_H
56 #define kEqn_H
57 
58 #include "LESeddyViscosity.H"
59 
60 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
61 
62 namespace Foam
63 {
64 namespace LESModels
65 {
66 
67 /*---------------------------------------------------------------------------*\
68  Class kEqn Declaration
69 \*---------------------------------------------------------------------------*/
70 
71 template<class BasicMomentumTransportModel>
72 class kEqn
73 :
74  public LESeddyViscosity<BasicMomentumTransportModel>
75 {
76 protected:
77 
78  // Protected data
79 
80  // Fields
81 
83 
84 
85  // Model constants
86 
88 
89 
90  // Protected Member Functions
91 
92  virtual void correctNut();
93  virtual tmp<fvScalarMatrix> kSource() const;
94 
95 
96 public:
97 
98  typedef typename BasicMomentumTransportModel::alphaField alphaField;
99  typedef typename BasicMomentumTransportModel::rhoField rhoField;
100 
101 
102  //- Runtime type information
103  TypeName("kEqn");
104 
105 
106  // Constructors
107 
108  //- Constructor from components
109  kEqn
110  (
111  const alphaField& alpha,
112  const rhoField& rho,
113  const volVectorField& U,
114  const surfaceScalarField& alphaRhoPhi,
115  const surfaceScalarField& phi,
116  const viscosity& viscosity,
117  const word& type = typeName
118  );
119 
120  //- Disallow default bitwise copy construction
121  kEqn(const kEqn&) = delete;
122 
123 
124  //- Destructor
125  virtual ~kEqn()
126  {}
127 
128 
129  // Member Functions
130 
131  //- Read model coefficients if they have changed
132  virtual bool read();
133 
134  //- Return SGS kinetic energy
135  virtual tmp<volScalarField> k() const
136  {
137  return k_;
138  }
139 
140  //- Return sub-grid disipation rate
141  virtual tmp<volScalarField> epsilon() const;
142 
143  //- Return the effective diffusivity for k
144  tmp<volScalarField> DkEff() const
145  {
146  return volScalarField::New
147  (
148  "DkEff",
149  this->nut_ + this->nu()
150  );
151  }
152 
153  //- Correct eddy-Viscosity and related properties
154  virtual void correct();
155 
156 
157  // Member Operators
158 
159  //- Disallow default bitwise assignment
160  void operator=(const kEqn&) = delete;
161 };
162 
163 
164 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
165 
166 } // End namespace LESModels
167 } // End namespace Foam
168 
169 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
170 
171 #ifdef NoRepository
172  #include "kEqn.C"
173 #endif
174 
175 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
176 
177 #endif
178 
179 // ************************************************************************* //
Generic GeometricField class.
static tmp< GeometricField< Type, PatchField, GeoMesh > > New(const word &name, const Internal &, const PtrList< PatchField< Type >> &, const HashPtrTable< Source > &=HashPtrTable< Source >())
Return a temporary field constructed from name,.
Eddy viscosity LES SGS model base class.
One equation eddy-viscosity model.
Definition: kEqn.H:74
BasicMomentumTransportModel::alphaField alphaField
Definition: kEqn.H:97
tmp< volScalarField > DkEff() const
Return the effective diffusivity for k.
Definition: kEqn.H:143
volScalarField k_
Definition: kEqn.H:81
dimensionedScalar Ck_
Definition: kEqn.H:86
virtual tmp< volScalarField > k() const
Return SGS kinetic energy.
Definition: kEqn.H:134
virtual void correct()
Correct eddy-Viscosity and related properties.
Definition: kEqn.C:150
void operator=(const kEqn &)=delete
Disallow default bitwise assignment.
virtual tmp< volScalarField > epsilon() const
Return sub-grid disipation rate.
Definition: kEqn.C:139
virtual ~kEqn()
Destructor.
Definition: kEqn.H:124
virtual void correctNut()
Definition: kEqn.C:40
virtual tmp< fvScalarMatrix > kSource() const
Definition: kEqn.C:49
TypeName("kEqn")
Runtime type information.
kEqn(const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const viscosity &viscosity, const word &type=typeName)
Constructor from components.
Definition: kEqn.C:67
virtual bool read()
Read model coefficients if they have changed.
Definition: kEqn.C:123
BasicMomentumTransportModel::rhoField rhoField
Definition: kEqn.H:98
A class for managing temporary objects.
Definition: tmp.H:55
Abstract base class for all fluid physical properties.
Definition: viscosity.H:50
A class for handling words, derived from string.
Definition: word.H:62
U
Definition: pEqn.H:72
volScalarField alpha(IOobject("alpha", runTime.name(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
Namespace for OpenFOAM.
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488