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-2015 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 
28 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
29 
30 namespace Foam
31 {
32 namespace LESModels
33 {
34 
35 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
36 
37 template<class BasicTurbulenceModel>
38 SmagorinskyZhang<BasicTurbulenceModel>::SmagorinskyZhang
39 (
40  const alphaField& alpha,
41  const rhoField& rho,
42  const volVectorField& U,
43  const surfaceScalarField& alphaRhoPhi,
44  const surfaceScalarField& phi,
45  const transportModel& transport,
46  const word& propertiesName,
47  const word& type
48 )
49 :
51  (
52  alpha,
53  rho,
54  U,
55  alphaRhoPhi,
56  phi,
57  transport,
58  propertiesName,
59  type
60  ),
61 
62  gasTurbulencePtr_(NULL),
63 
64  Cmub_
65  (
67  (
68  "Cmub",
69  this->coeffDict_,
70  0.6
71  )
72  )
73 {
74  if (type == typeName)
75  {
76  // Cannot correct nut yet: construction of the phases is not complete
77  // correctNut();
78 
79  this->printCoeffs(type);
80  }
81 }
82 
83 
84 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
85 
86 template<class BasicTurbulenceModel>
88 {
90  {
91  Cmub_.readIfPresent(this->coeffDict());
92 
93  return true;
94  }
95  else
96  {
97  return false;
98  }
99 }
100 
101 
102 template<class BasicTurbulenceModel>
104 <
105  typename BasicTurbulenceModel::transportModel
106 >&
108 {
109  if (!gasTurbulencePtr_)
110  {
111  const volVectorField& U = this->U_;
112 
113  const transportModel& liquid = this->transport();
114  const twoPhaseSystem& fluid =
115  refCast<const twoPhaseSystem>(liquid.fluid());
116  const transportModel& gas = fluid.otherPhase(liquid);
117 
118  gasTurbulencePtr_ =
119  &U.db()
121  (
123  (
125  gas.name()
126  )
127  );
128  }
129 
130  return *gasTurbulencePtr_;
131 }
132 
133 
134 template<class BasicTurbulenceModel>
136 {
138  this->gasTurbulence();
139 
140  volScalarField k(this->k(fvc::grad(this->U_)));
141 
142  this->nut_ =
143  this->Ck_*sqrt(k)*this->delta()
144  + Cmub_*gasTurbulence.transport().d()*gasTurbulence.alpha()
145  *(mag(this->U_ - gasTurbulence.U()));
146 
147  this->nut_.correctBoundaryConditions();
148 }
149 
150 
151 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
152 
153 } // End namespace LESModels
154 } // End namespace Foam
155 
156 // ************************************************************************* //
dimensionedScalar sqrt(const dimensionedScalar &ds)
multiphaseSystem & fluid
Definition: createFields.H:10
dimensioned< scalar > mag(const dimensioned< Type > &)
A class for handling words, derived from string.
Definition: word.H:59
The Smagorinsky SGS model including bubble-generated turbulence.
Namespace for OpenFOAM.
static word groupName(Name name, const word &group)
Generic dimensioned Type class.
const phaseModel & otherPhase(const phaseModel &phase) const
Constant access the phase not given as an argument.
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
Class which solves the volume fraction equations for two phases.
BasicTurbulenceModel::transportModel transportModel
label k
Boltzmann constant.
BasicTurbulenceModel::alphaField alphaField
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
Templated abstract base class for multiphase compressible turbulence models.
scalar delta
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:589
surfaceScalarField & phi
The Smagorinsky SGS model.
Definition: Smagorinsky.H:89
virtual bool read()
Read model coefficients if they have changed.
static const word propertiesName
Default name of the turbulence properties dictionary.
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:52
BasicTurbulenceModel::rhoField rhoField
const objectRegistry & db() const
Return the local objectRegistry.
Definition: IOobject.C:245
virtual void correctNut()
Update the SGS eddy viscosity.
U
Definition: pEqn.H:82