kOmegaSSTSato.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) 2014-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 "kOmegaSSTSato.H"
28 #include "twoPhaseSystem.H"
29 #include "dragModel.H"
30 
31 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 namespace RASModels
36 {
37 
38 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
39 
40 template<class BasicTurbulenceModel>
41 kOmegaSSTSato<BasicTurbulenceModel>::kOmegaSSTSato
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  gasTurbulencePtr_(NULL),
66 
67  Cmub_
68  (
70  (
71  "Cmub",
72  this->coeffDict_,
73  0.6
74  )
75  )
76 {
77  if (type == typeName)
78  {
79  // Cannot correct nut yet: construction of the phases is not complete
80  // correctNut();
81  this->printCoeffs(type);
82  }
83 }
84 
85 
86 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
87 
88 template<class BasicTurbulenceModel>
90 {
92  {
93  Cmub_.readIfPresent(this->coeffDict());
94 
95  return true;
96  }
97  else
98  {
99  return false;
100  }
101 }
102 
103 template<class BasicTurbulenceModel>
105 <
106  typename BasicTurbulenceModel::transportModel
107 >&
109 {
110  if (!gasTurbulencePtr_)
111  {
112  const volVectorField& U = this->U_;
113 
114  const transportModel& liquid = this->transport();
115  const twoPhaseSystem& fluid =
116  refCast<const twoPhaseSystem>(liquid.fluid());
117  const transportModel& gas = fluid.otherPhase(liquid);
118 
119  gasTurbulencePtr_ =
120  &U.db()
122  (
124  (
126  gas.name()
127  )
128  );
129  }
130 
131  return *gasTurbulencePtr_;
132 }
133 
134 
135 template<class BasicTurbulenceModel>
137 {
139  this->gasTurbulence();
140 
142  (
143  pow(this->betaStar_, 0.25)*this->y_*sqrt(this->k_)/this->nu()
144  );
145 
146  this->nut_ =
147  this->a1_*this->k_
148  /max
149  (
150  this->a1_*this->omega_,
151  this->F23()*sqrt(2.0)*mag(symm(fvc::grad(this->U_)))
152  )
153  + sqr(1 - exp(-yPlus/16.0))
154  *Cmub_*gasTurbulence.transport().d()*gasTurbulence.alpha()
155  *(mag(this->U_ - gasTurbulence.U()));
156 
157  this->nut_.correctBoundaryConditions();
158 }
159 
160 
161 template<class BasicTurbulenceModel>
163 {
165 }
166 
167 
168 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
169 
170 } // End namespace RASModels
171 } // End namespace Foam
172 
173 // ************************************************************************* //
dimensionedScalar sqrt(const dimensionedScalar &ds)
BasicTurbulenceModel::rhoField rhoField
BasicTurbulenceModel::alphaField alphaField
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
Definition: kOmegaSST.C:390
BasicTurbulenceModel::transportModel transportModel
multiphaseSystem & fluid
Definition: createFields.H:10
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
Implementation of the k-omega-SST turbulence model for incompressible and compressible flows...
Definition: kOmegaSST.H:120
dimensioned< scalar > mag(const dimensioned< Type > &)
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
A class for handling words, derived from string.
Definition: word.H:59
dimensionedScalar exp(const dimensionedScalar &ds)
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.
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
Implementation of the k-omega-SST turbulence model for dispersed bubbly flows with Sato (1981) bubble...
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.
Macros for easy insertion into run-time selection tables.
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:589
surfaceScalarField & phi
scalar yPlus
static const word propertiesName
Default name of the turbulence properties dictionary.
virtual bool read()
Read model coefficients if they have changed.
Definition: kOmegaSSTSato.C:89
volScalarField & nu
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:52
dimensionedSymmTensor symm(const dimensionedSymmTensor &dt)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const objectRegistry & db() const
Return the local objectRegistry.
Definition: IOobject.C:245
U
Definition: pEqn.H:82