buoyantKEpsilon.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) 2014-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::RASModels::buoyantKEpsilon
26 
27 Description
28  Additional buoyancy generation/dissipation term applied to the
29  k and epsilon equations of the standard k-epsilon model.
30 
31  Reference:
32  \verbatim
33  Henkes, R.A.W.M., Van Der Vlugt, F.F. & Hoogendoorn, C.J. (1991).
34  Natural Convection Flow in a Square Cavity Calculated with
35  Low-Reynolds-Number Turbulence Models.
36  Int. J. Heat Mass Transfer, 34, 1543-1557.
37  \endverbatim
38 
39  This implementation is based on the density rather than temperature gradient
40  extending the applicability to systems in which the density gradient may be
41  generated by variation of composition rather than temperature. Further, the
42  1/Prt coefficient is replaced by Cg to provide more general control of
43  model.
44 
45  The default model coefficients are
46  \verbatim
47  buoyantKEpsilonCoeffs
48  {
49  Cg 1.0;
50  }
51  \endverbatim
52 
53 See also
54  Foam::RASModels::kEpsilon
55 
56 SourceFiles
57  buoyantKEpsilon.C
58 
59 \*---------------------------------------------------------------------------*/
60 
61 #ifndef buoyantKEpsilon_H
62 #define buoyantKEpsilon_H
63 
64 #include "kEpsilon.H"
65 
66 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
67 
68 namespace Foam
69 {
70 namespace RASModels
71 {
72 
73 /*---------------------------------------------------------------------------*\
74  Class buoyantKEpsilon Declaration
75 \*---------------------------------------------------------------------------*/
76 
77 template<class BasicMomentumTransportModel>
78 class buoyantKEpsilon
79 :
80  public kEpsilon<BasicMomentumTransportModel>
81 {
82 protected:
83 
84  // Protected data
85 
86  // Model coefficients
87 
89 
90 
91  // Protected Member Functions
92 
93  tmp<volScalarField> Gcoef() const;
94 
95  virtual tmp<fvScalarMatrix> kSource() const;
96  virtual tmp<fvScalarMatrix> epsilonSource() const;
97 
98 
99 public:
101  typedef typename BasicMomentumTransportModel::alphaField alphaField;
102  typedef typename BasicMomentumTransportModel::rhoField rhoField;
103  typedef typename BasicMomentumTransportModel::transportModel transportModel;
104 
105 
106  //- Runtime type information
107  TypeName("buoyantKEpsilon");
108 
109 
110  // Constructors
111 
112  //- Construct from components
114  (
115  const alphaField& alpha,
116  const rhoField& rho,
117  const volVectorField& U,
118  const surfaceScalarField& alphaRhoPhi,
119  const surfaceScalarField& phi,
120  const transportModel& transport,
121  const word& type = typeName
122  );
123 
124  //- Disallow default bitwise copy construction
125  buoyantKEpsilon(const buoyantKEpsilon&) = delete;
126 
127 
128  //- Destructor
129  virtual ~buoyantKEpsilon()
130  {}
131 
132 
133  // Member Functions
134 
135  //- Re-read model coefficients if they have changed
136  virtual bool read();
137 
138 
139  // Member Operators
140 
141  //- Disallow default bitwise assignment
142  void operator=(const buoyantKEpsilon&) = delete;
143 };
144 
145 
146 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
147 
148 } // End namespace RASModels
149 } // End namespace Foam
150 
151 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
152 
153 #ifdef NoRepository
154  #include "buoyantKEpsilon.C"
155 #endif
156 
157 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
158 
159 #endif
160 
161 // ************************************************************************* //
tmp< volScalarField > Gcoef() const
void operator=(const buoyantKEpsilon &)=delete
Disallow default bitwise assignment.
virtual ~buoyantKEpsilon()
Destructor.
Additional buoyancy generation/dissipation term applied to the k and epsilon equations of the standar...
virtual tmp< fvScalarMatrix > epsilonSource() const
virtual tmp< fvScalarMatrix > kSource() const
phi
Definition: pEqn.H:104
virtual bool read()
Re-read model coefficients if they have changed.
TypeName("buoyantKEpsilon")
Runtime type information.
Standard k-epsilon turbulence model for incompressible and compressible flows including rapid distort...
Definition: kEpsilon.H:83
A class for handling words, derived from string.
Definition: word.H:59
buoyantKEpsilon(const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport, const word &type=typeName)
Construct from components.
U
Definition: pEqn.H:72
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
BasicMomentumTransportModel::rhoField rhoField
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
A class for managing temporary objects.
Definition: PtrList.H:53
BasicMomentumTransportModel::transportModel transportModel
Namespace for OpenFOAM.
BasicMomentumTransportModel::alphaField alphaField