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-2021 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 "fvModels.H"
28 #include "fvConstraints.H"
29 #include "phaseSystem.H"
30 
31 
32 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36 namespace RASModels
37 {
38 
39 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
40 
41 template<class BasicMomentumTransportModel>
43 (
44  const alphaField& alpha,
45  const rhoField& rho,
46  const volVectorField& U,
47  const surfaceScalarField& alphaRhoPhi,
48  const surfaceScalarField& phi,
49  const transportModel& transport,
50  const word& type
51 )
52 :
54  (
55  alpha,
56  rho,
57  U,
58  alphaRhoPhi,
59  phi,
60  transport,
61  type
62  ),
63 
64  phase_(refCast<const phaseModel>(transport)),
65 
66  hasDispersedPhaseNames_(this->coeffDict_.found("dispersedPhases")),
67 
68  dispersedPhaseNames_
69  (
70  this->coeffDict_.lookupOrDefault("dispersedPhases", hashedWordList())
71  ),
72 
73  Cmub_
74  (
76  (
77  "Cmub",
78  this->coeffDict_,
79  0.6
80  )
81  )
82 {
83  if (type == typeName)
84  {
85  this->printCoeffs(type);
86  }
87 }
88 
89 
90 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
91 
92 template<class BasicMomentumTransportModel>
94 {
96  {
97  Cmub_.readIfPresent(this->coeffDict());
98 
99  return true;
100  }
101  else
102  {
103  return false;
104  }
105 }
106 
107 
108 template<class BasicMomentumTransportModel>
111 {
112  UPtrList<const phaseModel> dispersedPhases;
113 
114  const phaseSystem& fluid = phase_.fluid();
115 
116  if (hasDispersedPhaseNames_)
117  {
118  dispersedPhases.resize(dispersedPhaseNames_.size());
119 
120  forAll(dispersedPhaseNames_, dispersedPhasei)
121  {
122  dispersedPhases.set
123  (
124  dispersedPhasei,
125  &fluid.phases()[dispersedPhaseNames_[dispersedPhasei]]
126  );
127  }
128  }
129  else
130  {
131  dispersedPhases.resize(fluid.movingPhases().size() - 1);
132 
133  label dispersedPhasei = 0;
134 
135  forAll(fluid.movingPhases(), movingPhasei)
136  {
137  const phaseModel& otherPhase = fluid.movingPhases()[movingPhasei];
138 
139  if (&otherPhase != &phase_)
140  {
141  dispersedPhases.set
142  (
143  dispersedPhasei ++,
144  &otherPhase
145  );
146  }
147  }
148  }
149 
150  return dispersedPhases;
151 }
152 
153 
154 template<class BasicMomentumTransportModel>
156 (
157  const volScalarField& S2,
158  const volScalarField& F2
159 )
160 {
162  (
163  pow(this->betaStar_, 0.25)*this->y_*sqrt(this->k_)/this->nu()
164  );
165 
166  this->nut_ =
167  this->a1_*this->k_/max(this->a1_*this->omega_, this->b1_*F2*sqrt(S2));
168 
169  UPtrList<const phaseModel> dispersedPhases(getDispersedPhases());
170  forAllIter(UPtrList<const phaseModel>, dispersedPhases, phaseIter)
171  {
172  this->nut_ +=
173  sqr(1 - exp(-yPlus/16.0))
174  *Cmub_*phaseIter().d()*phaseIter()
175  *(mag(this->U_ - phaseIter().U()));
176  }
177 
178  this->nut_.correctBoundaryConditions();
179  fvConstraints::New(this->mesh_).constrain(this->nut_);
180 }
181 
182 
183 template<class BasicMomentumTransportModel>
185 {
187 }
188 
189 
190 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
191 
192 } // End namespace RASModels
193 } // End namespace Foam
194 
195 // ************************************************************************* //
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 > &)
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:43
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
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.
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:70
A class for handling words, derived from string.
Definition: word.H:59
phaseSystem & fluid
Definition: createFields.H:11
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
bool set(const label) const
Is element set.
Definition: UPtrListI.H:78
static autoPtr< dictionary > New(Istream &)
Construct top-level dictionary on freestore from Istream.
Definition: dictionaryIO.C:96
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:93
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 > &)
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