All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
velocityGroup.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) 2017-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 "velocityGroup.H"
27 #include "sizeGroup.H"
28 #include "populationBalanceModel.H"
31 
32 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
33 
34 namespace Foam
35 {
36 namespace diameterModels
37 {
38  defineTypeNameAndDebug(velocityGroup, 0);
39  addToRunTimeSelectionTable(diameterModel, velocityGroup, dictionary);
40 }
41 }
42 
43 
44 // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
45 
46 Foam::tmp<Foam::volScalarField> Foam::diameterModels::velocityGroup::dsm() const
47 {
48  tmp<volScalarField> tInvDsm
49  (
51  (
52  "invDsm",
53  phase().mesh(),
55  )
56  );
57 
58  volScalarField& invDsm = tInvDsm.ref();
59 
60  forAll(sizeGroups_, i)
61  {
62  const sizeGroup& fi = sizeGroups_[i];
63 
64  invDsm += fi.a()*fi/fi.x();
65  }
66 
67  return 6/tInvDsm;
68 }
69 
70 
71 Foam::tmp<Foam::volScalarField> Foam::diameterModels::velocityGroup::N() const
72 {
73  tmp<volScalarField> tN
74  (
76  (
77  "N",
78  phase().mesh(),
80  )
81  );
82 
83  volScalarField& N = tN.ref();
84 
85  forAll(sizeGroups_, i)
86  {
87  N += phase()*sizeGroups_[i]/sizeGroups_[i].x();
88  }
89 
90  return tN;
91 }
92 
93 
95 Foam::diameterModels::velocityGroup::fSum() const
96 {
97  tmp<volScalarField> tsumSizeGroups
98  (
100  (
101  "sumSizeGroups",
102  phase().mesh(),
104  )
105  );
106 
107  volScalarField& sumSizeGroups = tsumSizeGroups.ref();
108 
109  forAll(sizeGroups_, i)
110  {
111  sumSizeGroups += sizeGroups_[i];
112  }
113 
114  return tsumSizeGroups;
115 }
116 
117 
118 void Foam::diameterModels::velocityGroup::scale()
119 {
120  Info<< "Scaling sizeGroups for velocityGroup " << phase().name() << endl;
121 
122  forAll(sizeGroups_, i)
123  {
124  sizeGroups_[i].max(0);
125  };
126 
127  f_ = fSum();
128 
129  forAll(sizeGroups_, i)
130  {
131  sizeGroups_[i] /= f_;
132 
133  sizeGroups_[i].correctBoundaryConditions();
134  };
135 }
136 
137 
138 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
139 
141 (
142  const dictionary& diameterProperties,
143  const phaseModel& phase
144 )
145 :
146  diameterModel(diameterProperties, phase),
147  popBalName_(diameterProperties.lookup("populationBalance")),
148  f_
149  (
150  IOobject
151  (
152  IOobject::groupName
153  (
154  "f",
155  phase.name()
156  ),
157  phase.time().timeName(),
158  phase.mesh(),
159  IOobject::MUST_READ,
160  IOobject::AUTO_WRITE
161  ),
162  phase.mesh()
163  ),
164  sizeGroups_
165  (
166  diameterProperties.lookup("sizeGroups"),
167  sizeGroup::iNew(phase, *this)
168  ),
169  d_(IOobject::groupName("d", phase.name()), dsm())
170 {}
171 
172 
173 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
174 
176 {}
177 
178 
179 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
180 
182 {
183  return d_;
184 }
185 
186 
188 {
189  tmp<volScalarField> tA
190  (
192  (
193  "a",
194  phase().mesh(),
196  )
197  );
198 
199  volScalarField& a = tA.ref();
200 
201  forAll(sizeGroups_, i)
202  {
203  const sizeGroup& fi = sizeGroups_[i];
204 
205  a += fi.a()*fi/fi.x();
206  }
207 
208  return phase()*a;
209 }
210 
211 
213 {
214  forAll(sizeGroups_, i)
215  {
216  sizeGroups_[i].correct();
217  }
218 
219  if
220  (
221  phase().mesh().solverDict(popBalName_).lookupOrDefault<Switch>
222  (
223  "scale",
224  true
225  )
226  )
227  {
228  scale();
229  }
230 
231  f_ = fSum();
232 
233  f_.correctBoundaryConditions();
234 
235  Info<< phase().name() << " sizeGroups-sum volume fraction, min, max = "
236  << f_.weightedAverage(phase().mesh().V()).value()
237  << ' ' << min(f_).value()
238  << ' ' << max(f_).value()
239  << endl;
240 
241  d_ = dsm();
242 
243  Info<< this->phase().name() << " Sauter mean diameter, min, max = "
244  << d_.weightedAverage(d_.mesh().V()).value()
245  << ' ' << min(d_).value()
246  << ' ' << max(d_).value()
247  << endl;
248 }
249 
250 
252 (
253  const dictionary& phaseProperties
254 )
255 {
256  diameterModel::read(phaseProperties);
257 
258  return true;
259 }
260 
261 
262 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
virtual tmp< volScalarField > d() const
Get the diameter field.
dimensionedSphericalTensor inv(const dimensionedSphericalTensor &dt)
static tmp< GeometricField< scalar, fvPatchField, volMesh > > New(const word &name, const Internal &, const PtrList< fvPatchField< scalar >> &)
Return a temporary field constructed from name,.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
sizeGroup(const word &name, const dictionary &dict, const phaseModel &phase, const velocityGroup &velocityGroup, const fvMesh &mesh)
const word & name() const
Definition: phaseModel.H:110
const dimensionSet dimless
velocityGroup(const dictionary &diameterProperties, const phaseModel &phase)
Construct from dictionary and phase.
const phaseModel & phase() const
Return const-reference to the phase.
Definition: sizeGroupI.H:36
Macros for easy insertion into run-time selection tables.
virtual tmp< volScalarField > a() const
Get the surface area per unit volume field.
const dimensionSet dimLength
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:57
const tmp< volScalarField > a() const
Return representative surface area of the sizeGroup.
stressControl lookup("compactNormalStress") >> compactNormalStress
dynamicFvMesh & mesh
virtual bool read(const dictionary &phaseProperties)
Read phaseProperties dictionary.
virtual ~velocityGroup()
Destructor.
word timeName
Definition: getTimeIndex.H:3
void min(const dimensioned< Type > &)
static const zero Zero
Definition: zero.H:97
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
const Mesh & mesh() const
Return mesh.
defineTypeNameAndDebug(combustionModel, 0)
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
Internal & ref()
Return a reference to the dimensioned internal field.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
void max(const dimensioned< Type > &)
const dimensionSet dimVolume
messageStream Info
const autoPtr< label > & i() const
Return label of the sizeGroup within the population balance.
virtual bool read(const dictionary &diameterProperties)
Read diameterProperties dictionary.
A class for managing temporary objects.
Definition: PtrList.H:53
Namespace for OpenFOAM.
virtual void correct()
Correct the model.