MRFZoneList.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) 2012-2023 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 "MRFZoneList.H"
27 #include "volFields.H"
29 
30 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
31 
33 (
34  const fvMesh& mesh,
35  const dictionary& dict
36 )
37 :
38  PtrList<MRFZone>(),
39  mesh_(mesh)
40 {
41  reset(dict);
42 }
43 
44 
45 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
46 
48 {}
49 
50 
51 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
52 
54 {
55  label count = 0;
57  {
58  if (iter().isDict())
59  {
60  count++;
61  }
62  }
63 
64  this->setSize(count);
65  label i = 0;
67  {
68  if (iter().isDict())
69  {
70  const word& name = iter().keyword();
71  const dictionary& modelDict = iter().dict();
72 
73  Info<< " creating MRF zone: " << name << endl;
74 
75  this->set
76  (
77  i++,
78  new MRFZone(name, mesh_, modelDict)
79  );
80  }
81  }
82 }
83 
84 
86 {
87  bool allOk = true;
88  forAll(*this, i)
89  {
90  MRFZone& pm = this->operator[](i);
91  bool ok = pm.read(dict.subDict(pm.name()));
92  allOk = (allOk && ok);
93  }
94  return allOk;
95 }
96 
97 
99 (
100  const volVectorField& U
101 ) const
102 {
104  (
106  (
107  "MRFZoneList:DDt",
108  U.mesh(),
109  dimensionedVector(U.dimensions()/dimTime, Zero)
110  )
111  );
112  volVectorField& DDt = tDDt.ref();
113 
114  forAll(*this, i)
115  {
116  operator[](i).addCoriolis(U, DDt);
117  }
118 
119  return tDDt;
120 }
121 
122 
124 (
125  const volScalarField& rho,
126  const volVectorField& U
127 ) const
128 {
129  return rho*DDt(U);
130 }
131 
132 
135 {
136  tmp<volVectorField> tcentrifugalAcceleration
137  (
139  (
140  "MRFZoneList:centrifugalAcceleration",
141  mesh_,
143  )
144  );
145  volVectorField& centrifugalAcceleration = tcentrifugalAcceleration.ref();
146 
147  forAll(*this, i)
148  {
149  operator[](i).addCentrifugalAcceleration(centrifugalAcceleration);
150  }
151 
152  return tcentrifugalAcceleration;
153 }
154 
155 
156 
158 {
159  forAll(*this, i)
160  {
161  operator[](i).makeRelative(U);
162  }
163 }
164 
165 
167 {
168  forAll(*this, i)
169  {
170  operator[](i).makeRelative(phi);
171  }
172 }
173 
174 
176 (
177  const tmp<surfaceScalarField>& tphi
178 ) const
179 {
180  if (size())
181  {
183  (
184  New
185  (
186  tphi,
187  "relative(" + tphi().name() + ')',
188  tphi().dimensions(),
189  true
190  )
191  );
192 
193  makeRelative(rphi.ref());
194 
195  tphi.clear();
196 
197  return rphi;
198  }
199  else
200  {
201  return tmp<surfaceScalarField>(tphi, true);
202  }
203 }
204 
205 
208 (
210 ) const
211 {
212  if (size())
213  {
214  tmp<FieldField<fvsPatchField, scalar>> rphi(New(tphi, true));
215 
216  forAll(*this, i)
217  {
218  operator[](i).makeRelative(rphi.ref());
219  }
220 
221  tphi.clear();
222 
223  return rphi;
224  }
225  else
226  {
227  return tmp<FieldField<fvsPatchField, scalar>>(tphi, true);
228  }
229 }
230 
231 
234 (
235  const tmp<Field<scalar>>& tphi,
236  const label patchi
237 ) const
238 {
239  if (size())
240  {
241  tmp<Field<scalar>> rphi(New(tphi, true));
242 
243  forAll(*this, i)
244  {
245  operator[](i).makeRelative(rphi.ref(), patchi);
246  }
247 
248  tphi.clear();
249 
250  return rphi;
251  }
252  else
253  {
254  return tmp<Field<scalar>>(tphi, true);
255  }
256 }
257 
258 
260 (
261  const surfaceScalarField& rho,
262  surfaceScalarField& phi
263 ) const
264 {
265  forAll(*this, i)
266  {
267  operator[](i).makeRelative(rho, phi);
268  }
269 }
270 
271 
273 {
274  forAll(*this, i)
275  {
276  operator[](i).makeAbsolute(U);
277  }
278 }
279 
280 
282 {
283  forAll(*this, i)
284  {
285  operator[](i).makeAbsolute(phi);
286  }
287 }
288 
289 
291 (
292  const tmp<surfaceScalarField>& tphi
293 ) const
294 {
295  if (size())
296  {
298  (
299  New
300  (
301  tphi,
302  "absolute(" + tphi().name() + ')',
303  tphi().dimensions(),
304  true
305  )
306  );
307 
308  makeAbsolute(rphi.ref());
309 
310  tphi.clear();
311 
312  return rphi;
313  }
314  else
315  {
316  return tmp<surfaceScalarField>(tphi, true);
317  }
318 }
319 
320 
322 (
323  const surfaceScalarField& rho,
324  surfaceScalarField& phi
325 ) const
326 {
327  forAll(*this, i)
328  {
329  operator[](i).makeAbsolute(rho, phi);
330  }
331 }
332 
333 
335 (
336  const tmp<surfaceScalarField>& tphi,
337  const volScalarField& rho
338 ) const
339 {
340  if (size())
341  {
343  (
344  New
345  (
346  tphi,
347  "absolute(" + tphi().name() + ')',
348  tphi().dimensions(),
349  true
350  )
351  );
352 
354 
355  tphi.clear();
356 
357  return rphi;
358  }
359  else
360  {
361  return tmp<surfaceScalarField>(tphi, true);
362  }
363 }
364 
365 
367 {
368  if (mesh_.topoChanged())
369  {
370  forAll(*this, i)
371  {
372  operator[](i).update();
373  }
374  }
375 }
376 
377 
378 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object of type.
Definition: UList.H:477
Generic field type.
Definition: FieldField.H:77
Generic GeometricField class.
static tmp< GeometricField< Type, PatchField, GeoMesh > > New(const word &name, const Internal &, const PtrList< PatchField< Type >> &, const HashPtrTable< Source > &=HashPtrTable< Source >())
Return a temporary field constructed from name,.
tmp< surfaceScalarField > relative(const tmp< surfaceScalarField > &phi) const
Return the given absolute flux relative within the MRF region.
Definition: MRFZoneList.C:176
void reset(const dictionary &dict)
Reset the source list.
Definition: MRFZoneList.C:53
tmp< volVectorField > centrifugalAcceleration() const
Return the centrifugal acceleration.
Definition: MRFZoneList.C:134
tmp< surfaceScalarField > absolute(const tmp< surfaceScalarField > &phi) const
Return the given relative flux absolute within the MRF region.
Definition: MRFZoneList.C:291
bool read(const dictionary &dict)
Read dictionary.
Definition: MRFZoneList.C:85
void makeAbsolute(volVectorField &U) const
Make the given relative velocity absolute within the MRF region.
Definition: MRFZoneList.C:272
void makeRelative(volVectorField &U) const
Make the given absolute velocity relative within the MRF region.
Definition: MRFZoneList.C:157
MRFZoneList(const fvMesh &mesh, const dictionary &dict)
Definition: MRFZoneList.C:33
void update()
Update MRFZone faces if the mesh topology changes.
Definition: MRFZoneList.C:366
~MRFZoneList()
Destructor.
Definition: MRFZoneList.C:47
tmp< volVectorField > DDt(const volVectorField &U) const
Return the Coriolis acceleration.
Definition: MRFZoneList.C:99
MRF zone definition based on cell zone and parameters obtained from a control dictionary constructed ...
Definition: MRFZone.H:56
bool read(const dictionary &dict)
Read MRF dictionary.
Definition: MRFZone.C:362
const word & name() const
Return const access to the MRF region name.
Definition: MRFZoneI.H:26
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: PtrList.H:75
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
Definition: dictionary.C:843
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
A class for managing temporary objects.
Definition: tmp.H:55
void clear() const
If object pointer points to valid object:
Definition: tmpI.H:237
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
A class for handling words, derived from string.
Definition: word.H:62
label patchi
U
Definition: pEqn.H:72
tmp< VolField< Type > > DDt(const surfaceScalarField &phi, const VolField< Type > &psi)
Definition: fvcDDt.C:45
static tmp< SurfaceField< Type > > interpolate(const VolField< Type > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
void makeAbsolute(surfaceScalarField &phi, const volVectorField &U)
Make the given flux absolute.
Definition: fvcMeshPhi.C:128
void makeRelative(surfaceScalarField &phi, const volVectorField &U)
Make the given flux relative.
Definition: fvcMeshPhi.C:89
static const zero Zero
Definition: zero.H:97
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
const HashTable< dimensionSet > & dimensions()
Get the table of dimension sets.
Definition: dimensionSets.C:96
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
messageStream Info
const dimensionSet dimAcceleration
const dimensionSet dimTime
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
dimensioned< vector > dimensionedVector
Dimensioned vector obtained from generic dimensioned type.
label count(const ListType &l, typename ListType::const_reference x)
Count the number of occurrences of a value in a list.
points setSize(newPointi)
dictionary dict