buoyantKEpsilon.C
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-2018 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 \*---------------------------------------------------------------------------*/
25 
26 #include "buoyantKEpsilon.H"
28 #include "fvcGrad.H"
30 
31 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 namespace RASModels
36 {
37 
38 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
39 
40 template<class BasicTurbulenceModel>
42 (
43  const alphaField& alpha,
44  const rhoField& rho,
45  const volVectorField& U,
46  const surfaceScalarField& alphaRhoPhi,
47  const surfaceScalarField& phi,
48  const transportModel& transport,
49  const word& propertiesName,
50  const word& type
51 )
52 :
54  (
55  alpha,
56  rho,
57  U,
58  alphaRhoPhi,
59  phi,
60  transport,
61  propertiesName,
62  type
63  ),
64 
65  Cg_
66  (
68  (
69  "Cg",
70  this->coeffDict_,
71  1.0
72  )
73  )
74 {
75  if (type == typeName)
76  {
77  this->printCoeffs(type);
78  }
79 }
80 
81 
82 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
83 
84 template<class BasicTurbulenceModel>
86 {
88  {
89  Cg_.readIfPresent(this->coeffDict());
90 
91  return true;
92  }
93  else
94  {
95  return false;
96  }
97 }
98 
99 
100 template<class BasicTurbulenceModel>
103 {
105  this->mesh_.objectRegistry::template
106  lookupObject<uniformDimensionedVectorField>("g");
107 
108  return
109  (Cg_*this->Cmu_)*this->alpha_*this->k_*(g & fvc::grad(this->rho_))
110  /(this->epsilon_ + this->epsilonMin_);
111 }
112 
113 
114 template<class BasicTurbulenceModel>
117 {
119  this->mesh_.objectRegistry::template
120  lookupObject<uniformDimensionedVectorField>("g");
121 
122  if (mag(g.value()) > small)
123  {
124  return -fvm::SuSp(Gcoef(), this->k_);
125  }
126  else
127  {
129  }
130 }
131 
132 
133 template<class BasicTurbulenceModel>
136 {
138  this->mesh_.objectRegistry::template
139  lookupObject<uniformDimensionedVectorField>("g");
140 
141  if (mag(g.value()) > small)
142  {
143  vector gHat(g.value()/mag(g.value()));
144 
145  volScalarField v(gHat & this->U_);
147  (
148  mag(this->U_ - gHat*v)
150  );
151 
152  return -fvm::SuSp(this->C1_*tanh(mag(v)/u)*Gcoef(), this->epsilon_);
153  }
154  else
155  {
157  }
158 }
159 
160 
161 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
162 
163 } // End namespace RASModels
164 } // End namespace Foam
165 
166 // ************************************************************************* //
tmp< volScalarField > Gcoef() const
virtual tmp< fvScalarMatrix > epsilonSource() const
Definition: kEpsilon.C:66
dimensionedScalar tanh(const dimensionedScalar &ds)
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:52
tmp< fvMatrix< Type > > SuSp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
surfaceScalarField & phi
virtual tmp< fvScalarMatrix > kSource() const
Definition: kEpsilon.C:51
BasicTurbulenceModel::rhoField rhoField
virtual tmp< fvScalarMatrix > epsilonSource() const
Generic dimensioned Type class.
virtual tmp< fvScalarMatrix > kSource() const
Macros for easy insertion into run-time selection tables.
virtual bool read()
Re-read model coefficients if they have changed.
Standard k-epsilon turbulence model for incompressible and compressible flows including rapid distort...
Definition: kEpsilon.H:83
Calculate the gradient of the given field.
A class for handling words, derived from string.
Definition: word.H:59
const Type & value() const
Return const reference to value.
BasicTurbulenceModel::alphaField alphaField
U
Definition: pEqn.H:72
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
dimensioned< scalar > mag(const dimensioned< Type > &)
buoyantKEpsilon(const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport, const word &propertiesName=turbulenceModel::propertiesName, const word &type=typeName)
Construct from components.
BasicTurbulenceModel::transportModel transportModel
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
const dimensionedVector & g
Namespace for OpenFOAM.
const dimensionSet dimVelocity