SmagorinskyZhang.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2013-2016 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 "SmagorinskyZhang.H"
27 #include "fvOptions.H"
28 
29 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
30 
31 namespace Foam
32 {
33 namespace LESModels
34 {
35 
36 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
37 
38 template<class BasicTurbulenceModel>
39 SmagorinskyZhang<BasicTurbulenceModel>::SmagorinskyZhang
40 (
41  const alphaField& alpha,
42  const rhoField& rho,
43  const volVectorField& U,
44  const surfaceScalarField& alphaRhoPhi,
45  const surfaceScalarField& phi,
46  const transportModel& transport,
47  const word& propertiesName,
48  const word& type
49 )
50 :
52  (
53  alpha,
54  rho,
55  U,
56  alphaRhoPhi,
57  phi,
58  transport,
59  propertiesName,
60  type
61  ),
62 
63  gasTurbulencePtr_(NULL),
64 
65  Cmub_
66  (
68  (
69  "Cmub",
70  this->coeffDict_,
71  0.6
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  Cmub_.readIfPresent(this->coeffDict());
90 
91  return true;
92  }
93  else
94  {
95  return false;
96  }
97 }
98 
99 
100 template<class BasicTurbulenceModel>
102 <
103  typename BasicTurbulenceModel::transportModel
104 >&
106 {
107  if (!gasTurbulencePtr_)
108  {
109  const volVectorField& U = this->U_;
110 
111  const transportModel& liquid = this->transport();
112  const twoPhaseSystem& fluid =
113  refCast<const twoPhaseSystem>(liquid.fluid());
114  const transportModel& gas = fluid.otherPhase(liquid);
115 
116  gasTurbulencePtr_ =
117  &U.db()
119  (
121  (
123  gas.name()
124  )
125  );
126  }
127 
128  return *gasTurbulencePtr_;
129 }
130 
131 
132 template<class BasicTurbulenceModel>
134 {
136  this->gasTurbulence();
137 
138  volScalarField k(this->k(fvc::grad(this->U_)));
139 
140  this->nut_ =
141  this->Ck_*sqrt(k)*this->delta()
142  + Cmub_*gasTurbulence.transport().d()*gasTurbulence.alpha()
143  *(mag(this->U_ - gasTurbulence.U()));
144 
145  this->nut_.correctBoundaryConditions();
146  fv::options::New(this->mesh_).correct(this->nut_);
147 
148  BasicTurbulenceModel::correctNut();
149 }
150 
151 
152 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
153 
154 } // End namespace LESModels
155 } // End namespace Foam
156 
157 // ************************************************************************* //
scalar delta
surfaceScalarField & phi
BasicTurbulenceModel::transportModel transportModel
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:52
const transportModel & transport() const
Access function to incompressible transport model.
U
Definition: pEqn.H:83
The Smagorinsky SGS model including bubble-generated turbulence.
multiphaseSystem & fluid
Definition: createFields.H:10
const objectRegistry & db() const
Return the local objectRegistry.
Definition: IOobject.C:221
void correct(GeometricField< Type, fvPatchField, volMesh > &field)
Apply correction to field.
const volVectorField & U() const
Access function to velocity field.
dimensionedScalar sqrt(const dimensionedScalar &ds)
label k
Boltzmann constant.
Generic dimensioned Type class.
The Smagorinsky SGS model.
Definition: Smagorinsky.H:89
Class which solves the volume fraction equations for two phases.
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
BasicTurbulenceModel::rhoField rhoField
static const word propertiesName
Default name of the turbulence properties dictionary.
A class for handling words, derived from string.
Definition: word.H:59
static word groupName(Name name, const word &group)
const phaseModel & otherPhase(const phaseModel &phase) const
Constant access the phase not given as an argument.
Templated abstract base class for multiphase compressible turbulence models.
virtual bool read()
Read model coefficients if they have changed.
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:461
dimensioned< scalar > mag(const dimensioned< Type > &)
static options & New(const fvMesh &mesh)
Construct fvOptions and register to datbase if not present.
Definition: fvOptions.C:103
virtual void correctNut()
Update the SGS eddy viscosity.
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
BasicTurbulenceModel::alphaField alphaField
Namespace for OpenFOAM.