kOmegaSSTSato.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::kOmegaSSTSato
26 
27 Description
28  Implementation of the k-omega-SST turbulence model for dispersed bubbly
29  flows with Sato (1981) bubble induced turbulent viscosity model.
30 
31  Bubble induced turbulent viscosity model described in:
32  \verbatim
33  Sato, Y., Sadatomi, M.
34  "Momentum and heat transfer in two-phase bubble flow - I, Theory"
35  International Journal of Multiphase FLow 7, pp. 167-177, 1981.
36  \endverbatim
37 
38  Turbulence model described in:
39  \verbatim
40  Menter, F., Esch, T.
41  "Elements of Industrial Heat Transfer Prediction"
42  16th Brazilian Congress of Mechanical Engineering (COBEM),
43  Nov. 2001
44  \endverbatim
45 
46  with the addition of the optional F3 term for rough walls from
47  \verbatim
48  Hellsten, A.
49  "Some Improvements in Menter’s k-omega-SST turbulence model"
50  29th AIAA Fluid Dynamics Conference,
51  AIAA-98-2554,
52  June 1998.
53  \endverbatim
54 
55  Note that this implementation is written in terms of alpha diffusion
56  coefficients rather than the more traditional sigma (alpha = 1/sigma) so
57  that the blending can be applied to all coefficients in a consistent
58  manner. The paper suggests that sigma is blended but this would not be
59  consistent with the blending of the k-epsilon and k-omega models.
60 
61  Also note that the error in the last term of equation (2) relating to
62  sigma has been corrected.
63 
64  Wall-functions are applied in this implementation by using equations (14)
65  to specify the near-wall omega as appropriate.
66 
67  The blending functions (15) and (16) are not currently used because of the
68  uncertainty in their origin, range of applicability and that as y+ becomes
69  sufficiently small blending u_tau in this manner is clearly nonsense.
70 
71  The default model coefficients correspond to the following:
72  \verbatim
73  kOmegaSSTCoeffs
74  {
75  alphaK1 0.85034;
76  alphaK2 1.0;
77  alphaOmega1 0.5;
78  alphaOmega2 0.85616;
79  Prt 1.0; // only for compressible
80  beta1 0.075;
81  beta2 0.0828;
82  betaStar 0.09;
83  gamma1 0.5532;
84  gamma2 0.4403;
85  a1 0.31;
86  b1 1.0;
87  c1 10.0;
88  F3 no;
89  Cmub 0.6;
90  }
91  \endverbatim
92 
93 SourceFiles
94  kOmegaSSTSato.C
95 
96 \*---------------------------------------------------------------------------*/
97 
98 #ifndef kOmegaSSTSato_H
99 #define kOmegaSSTSato_H
100 
101 #include "kOmegaSST.H"
102 #include "phaseSystem.H"
103 #include "hashedWordList.H"
104 
105 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
106 
107 namespace Foam
108 {
109 namespace RASModels
110 {
111 
112 /*---------------------------------------------------------------------------*\
113  Class kOmegaSSTSato Declaration
114 \*---------------------------------------------------------------------------*/
115 
116 template<class BasicMomentumTransportModel>
117 class kOmegaSSTSato
118 :
119  public kOmegaSST<BasicMomentumTransportModel>
120 {
121  // Private Data
122 
123  //- The phase to which the turbulent model is applied to
124  const phaseModel& phase_;
125 
126  //- True if a list of dispersed phases is specified in the dictionary
127  bool hasDispersedPhaseNames_;
128 
129  //- Optional list of names of dispersed phases
130  hashedWordList dispersedPhaseNames_;
131 
132  // Private Member Functions
133 
134  //- Return the dispersed phase models
135  UPtrList<const phaseModel> getDispersedPhases() const;
136 
137 
138 protected:
139 
140  // Protected data
141 
142  // Model coefficients
145 
146 
147  // Protected Member Functions
148 
149  virtual void correctNut
150  (
151  const volScalarField& S2,
152  const volScalarField& F2
153  );
154 
155 
156 public:
158  typedef typename BasicMomentumTransportModel::alphaField alphaField;
159  typedef typename BasicMomentumTransportModel::rhoField rhoField;
160  typedef typename BasicMomentumTransportModel::transportModel transportModel;
161 
162 
163  //- Runtime type information
164  TypeName("kOmegaSSTSato");
165 
166 
167  // Constructors
168 
169  //- Construct from components
171  (
172  const alphaField& alpha,
173  const rhoField& rho,
174  const volVectorField& U,
175  const surfaceScalarField& alphaRhoPhi,
176  const surfaceScalarField& phi,
177  const transportModel& transport,
178  const word& type = typeName
179  );
180 
181  //- Disallow default bitwise copy construction
182  kOmegaSSTSato(const kOmegaSSTSato&) = delete;
183 
184 
185  //- Destructor
186  virtual ~kOmegaSSTSato()
187  {}
188 
189 
190  // Member Functions
191 
192  //- Read model coefficients if they have changed
193  virtual bool read();
194 
195  //- Solve the turbulence equations and correct the turbulence viscosity
196  virtual void correct();
197 
198 
199  // Member Operators
200 
201  //- Disallow default bitwise assignment
202  void operator=(const kOmegaSSTSato&) = delete;
203 };
204 
205 
206 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
207 
208 } // End namespace RASModels
209 } // End namespace Foam
210 
211 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
212 
213 #ifdef NoRepository
214  #include "kOmegaSSTSato.C"
215 #endif
216 
217 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
218 
219 #endif
220 
221 // ************************************************************************* //
void operator=(const kOmegaSSTSato &)=delete
Disallow default bitwise assignment.
Implementation of the k-omega-SST turbulence model for dispersed bubbly flows with Sato (1981) bubble...
kOmegaSSTSato(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.
Definition: kOmegaSSTSato.C:43
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
TypeName("kOmegaSSTSato")
Runtime type information.
BasicMomentumTransportModel::alphaField alphaField
A class for handling words, derived from string.
Definition: word.H:59
virtual ~kOmegaSSTSato()
Destructor.
phi
Definition: correctPhi.H:3
Specialisation for RAS of the generic kOmegaSSTBase base class. For more information, see Description of kOmegaSSTBase.H.
Definition: kOmegaSST.H:57
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: UPtrList.H:54
BasicMomentumTransportModel::rhoField rhoField
virtual bool read()
Read model coefficients if they have changed.
Definition: kOmegaSSTSato.C:93
U
Definition: pEqn.H:72
A wordList with hashed indices for faster lookup by name.
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
BasicMomentumTransportModel::transportModel transportModel
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
Definition: phaseModel.H:53
Namespace for OpenFOAM.