kOmegaSSTSato.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-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 \*---------------------------------------------------------------------------*/
25 
26 #include "kOmegaSSTSato.H"
27 #include "fvOptions.H"
28 #include "phaseSystem.H"
29 
30 
31 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 namespace RASModels
36 {
37 
38 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
39 
40 template<class BasicMomentumTransportModel>
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& type
50 )
51 :
53  (
54  alpha,
55  rho,
56  U,
57  alphaRhoPhi,
58  phi,
59  transport,
60  type
61  ),
62 
63  phase_(transport),
64 
65  hasDispersedPhaseNames_(this->coeffDict_.found("dispersedPhases")),
66 
67  dispersedPhaseNames_
68  (
69  this->coeffDict_.lookupOrDefault("dispersedPhases", hashedWordList())
70  ),
71 
72  Cmub_
73  (
75  (
76  "Cmub",
77  this->coeffDict_,
78  0.6
79  )
80  )
81 {
82  if (type == typeName)
83  {
84  this->printCoeffs(type);
85  }
86 }
87 
88 
89 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
90 
91 template<class BasicMomentumTransportModel>
93 {
95  {
96  Cmub_.readIfPresent(this->coeffDict());
97 
98  return true;
99  }
100  else
101  {
102  return false;
103  }
104 }
105 
106 
107 template<class BasicMomentumTransportModel>
110 {
111  UPtrList<const phaseModel> dispersedPhases;
112 
113  const phaseSystem& fluid = phase_.fluid();
114 
115  if (hasDispersedPhaseNames_)
116  {
117  dispersedPhases.resize(dispersedPhaseNames_.size());
118 
119  forAll(dispersedPhaseNames_, dispersedPhasei)
120  {
121  dispersedPhases.set
122  (
123  dispersedPhasei,
124  &fluid.phases()[dispersedPhaseNames_[dispersedPhasei]]
125  );
126  }
127  }
128  else
129  {
130  dispersedPhases.resize(fluid.movingPhases().size() - 1);
131 
132  label dispersedPhasei = 0;
133 
134  forAll(fluid.movingPhases(), movingPhasei)
135  {
136  const phaseModel& otherPhase = fluid.movingPhases()[movingPhasei];
137 
138  if (&otherPhase != &phase_)
139  {
140  dispersedPhases.set
141  (
142  dispersedPhasei ++,
143  &otherPhase
144  );
145  }
146  }
147  }
148 
149  return dispersedPhases;
150 }
151 
152 
153 template<class BasicMomentumTransportModel>
155 (
156  const volScalarField& S2,
157  const volScalarField& F2
158 )
159 {
161  (
162  pow(this->betaStar_, 0.25)*this->y_*sqrt(this->k_)/this->nu()
163  );
164 
165  this->nut_ =
166  this->a1_*this->k_/max(this->a1_*this->omega_, this->b1_*F2*sqrt(S2));
167 
168  UPtrList<const phaseModel> dispersedPhases(getDispersedPhases());
169  forAllIter(UPtrList<const phaseModel>, dispersedPhases, phaseIter)
170  {
171  this->nut_ +=
172  sqr(1 - exp(-yPlus/16.0))
173  *Cmub_*phaseIter().d()*phaseIter()
174  *(mag(this->U_ - phaseIter().U()));
175  }
176 
177  this->nut_.correctBoundaryConditions();
178  fv::options::New(this->mesh_).correct(this->nut_);
179 }
180 
181 
182 template<class BasicMomentumTransportModel>
184 {
186 }
187 
188 
189 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
190 
191 } // End namespace RASModels
192 } // End namespace Foam
193 
194 // ************************************************************************* //
BasicMomentumTransportModel::alphaField alphaField
Definition: kOmegaSST.H:68
void resize(const label)
Reset size of UPtrList. This can only be used to set the size.
Definition: UPtrListI.H:71
const phaseModelPartialList & movingPhases() const
Return the models for phases that are moving.
Definition: phaseSystemI.H:49
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
Implementation of the k-omega-SST turbulence model for dispersed bubbly flows with Sato (1981) bubble...
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
void correct(GeometricField< Type, fvPatchField, volMesh > &field)
Apply correction to field.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
#define forAllIter(Container, container, iter)
Iterate across all elements in the container object of type.
Definition: UList.H:459
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:42
dimensionedScalar sqrt(const dimensionedScalar &ds)
BasicMomentumTransportModel::transportModel transportModel
Definition: kOmegaSST.H:70
#define F2(B, C, D)
Definition: SHA1.C:174
Generic dimensioned Type class.
phi
Definition: pEqn.H:104
const phaseModelList & phases() const
Return the phase models.
Definition: phaseSystemI.H:35
dimensionedScalar exp(const dimensionedScalar &ds)
Class to represent a system of phases and model interfacial transfers between them.
Definition: phaseSystem.H:67
A class for handling words, derived from string.
Definition: word.H:59
phaseSystem & fluid
Definition: createFields.H:11
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
bool set(const label) const
Is element set.
Definition: UPtrListI.H:78
label size() const
Return the number of elements in the UPtrList.
Definition: UPtrListI.H:29
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
virtual bool read()
Read model coefficients if they have changed.
Definition: kOmegaSSTSato.C:92
U
Definition: pEqn.H:72
A wordList with hashed indices for faster lookup by name.
scalar yPlus
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
dimensioned< scalar > mag(const dimensioned< Type > &)
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
static options & New(const fvMesh &mesh)
Construct fvOptions and register to datbase if not present.
Definition: fvOptions.C:101
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
Definition: phaseModel.H:53
Namespace for OpenFOAM.
BasicMomentumTransportModel::rhoField rhoField
Definition: kOmegaSST.H:69