fvMeshGeometry.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) 2011-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 "fvMesh.H"
27 #include "Time.H"
28 #include "volFields.H"
29 #include "surfaceFields.H"
30 #include "slicedVolFields.H"
31 #include "slicedSurfaceFields.H"
32 #include "SubField.H"
33 #include "cyclicFvPatchFields.H"
34 #include "cyclicAMIFvPatchFields.H"
35 
36 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
37 
38 void Foam::fvMesh::makeSf() const
39 {
40  if (debug)
41  {
42  InfoInFunction << "Assembling face areas" << endl;
43  }
44 
45  // It is an error to attempt to recalculate
46  // if the pointer is already set
47  if (SfPtr_)
48  {
50  << "face areas already exist"
51  << abort(FatalError);
52  }
53 
54  SfPtr_ = new slicedSurfaceVectorField
55  (
56  IOobject
57  (
58  "Sf",
60  meshSubDir,
61  *this,
64  false
65  ),
66  *this,
67  dimArea,
68  faceAreas()
69  );
70 }
71 
72 
73 void Foam::fvMesh::makeMagSf() const
74 {
75  if (debug)
76  {
77  InfoInFunction << "Assembling mag face areas" << endl;
78  }
79 
80  // It is an error to attempt to recalculate
81  // if the pointer is already set
82  if (magSfPtr_)
83  {
85  << "mag face areas already exist"
86  << abort(FatalError);
87  }
88 
89  magSfPtr_ = new slicedSurfaceScalarField
90  (
91  IOobject
92  (
93  "magSf",
95  meshSubDir,
96  *this,
99  false
100  ),
101  *this,
102  dimArea,
103  magFaceAreas()
104  );
105 }
106 
107 
108 void Foam::fvMesh::makeC() const
109 {
110  if (debug)
111  {
112  InfoInFunction << "Assembling cell centres" << endl;
113  }
114 
115  // It is an error to attempt to recalculate
116  // if the pointer is already set
117  if (CPtr_)
118  {
120  << "cell centres already exist"
121  << abort(FatalError);
122  }
123 
124  // Construct as slices. Only preserve processor (not e.g. cyclic)
125 
126  CPtr_ = new slicedVolVectorField
127  (
128  IOobject
129  (
130  "C",
131  pointsInstance(),
132  meshSubDir,
133  *this,
136  false
137  ),
138  *this,
139  dimLength,
140  cellCentres(),
141  faceCentres(),
142  true, // preserveCouples
143  true // preserveProcOnly
144  );
145 }
146 
147 
148 void Foam::fvMesh::makeCf() const
149 {
150  if (debug)
151  {
152  InfoInFunction << "Assembling face centres" << endl;
153  }
154 
155  // It is an error to attempt to recalculate
156  // if the pointer is already set
157  if (CfPtr_)
158  {
160  << "face centres already exist"
161  << abort(FatalError);
162  }
163 
164  CfPtr_ = new slicedSurfaceVectorField
165  (
166  IOobject
167  (
168  "Cf",
169  pointsInstance(),
170  meshSubDir,
171  *this,
174  false
175  ),
176  *this,
177  dimLength,
178  faceCentres()
179  );
180 }
181 
182 
183 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
184 
186 {
187  if (!VPtr_)
188  {
189  if (debug)
190  {
192  << "Constructing from primitiveMesh::cellVolumes()" << endl;
193  }
194 
196  (
197  IOobject
198  (
199  "V",
200  time().timeName(),
201  *this,
204  false
205  ),
206  *this,
207  dimVolume,
208  cellVolumes()
209  );
210  }
211 
212  return *static_cast<slicedVolScalarField::Internal*>(VPtr_);
213 }
214 
215 
217 {
218  if (!V0Ptr_)
219  {
221  << "V0 is not available"
222  << abort(FatalError);
223  }
224 
225  return *V0Ptr_;
226 }
227 
228 
230 {
231  if (!V0Ptr_)
232  {
234  << "V0 is not available"
235  << abort(FatalError);
236  }
237 
238  return *V0Ptr_;
239 }
240 
241 
243 {
244  if (!V00Ptr_)
245  {
246  if (debug)
247  {
248  InfoInFunction << "Constructing from V0" << endl;
249  }
250 
252  (
253  IOobject
254  (
255  "V00",
256  time().timeName(),
257  *this,
260  false
261  ),
262  V0()
263  );
264  }
265 
266  return *V00Ptr_;
267 }
268 
269 
272 {
273  if (moving() && time().subCycling())
274  {
275  const TimeState& ts = time();
276  const TimeState& ts0 = time().prevTimeState();
277 
278  scalar tFrac =
279  (
280  ts.value() - (ts0.value() - ts0.deltaTValue())
281  )/ts0.deltaTValue();
282 
283  if (tFrac < (1 - small))
284  {
285  return V0() + tFrac*(V() - V0());
286  }
287  else
288  {
289  return V();
290  }
291  }
292  else
293  {
294  return V();
295  }
296 }
297 
298 
301 {
302  if (moving() && time().subCycling())
303  {
304  const TimeState& ts = time();
305  const TimeState& ts0 = time().prevTimeState();
306 
307  scalar t0Frac =
308  (
309  (ts.value() - ts.deltaTValue())
310  - (ts0.value() - ts0.deltaTValue())
311  )/ts0.deltaTValue();
312 
313  if (t0Frac > small)
314  {
315  return V0() + t0Frac*(V() - V0());
316  }
317  else
318  {
319  return V0();
320  }
321  }
322  else
323  {
324  return V0();
325  }
326 }
327 
328 
330 {
331  if (!SfPtr_)
332  {
333  makeSf();
334  }
335 
336  return *SfPtr_;
337 }
338 
339 
341 {
342  if (!magSfPtr_)
343  {
344  makeMagSf();
345  }
346 
347  return *magSfPtr_;
348 }
349 
350 
352 {
353  if (!CPtr_)
354  {
355  makeC();
356  }
357 
358  return *CPtr_;
359 }
360 
361 
363 {
364  if (!CfPtr_)
365  {
366  makeCf();
367  }
368 
369  return *CfPtr_;
370 }
371 
372 
374 {
375  if (debug)
376  {
377  InfoInFunction << "Calculating face deltas" << endl;
378  }
379 
381  (
383  (
384  "delta",
385  *this,
386  dimLength
387  )
388  );
389  surfaceVectorField& delta = tdelta.ref();
390 
391  const volVectorField& C = this->C();
392  const labelUList& owner = this->owner();
393  const labelUList& neighbour = this->neighbour();
394 
395  forAll(owner, facei)
396  {
397  delta[facei] = C[neighbour[facei]] - C[owner[facei]];
398  }
399 
400  surfaceVectorField::Boundary& deltabf =
401  delta.boundaryFieldRef();
402 
403  forAll(deltabf, patchi)
404  {
405  deltabf[patchi] = boundary()[patchi].delta();
406  }
407 
408  return tdelta;
409 }
410 
411 
413 {
414  if (!phiPtr_)
415  {
417  << "mesh flux field does not exist, is the mesh actually moving?"
418  << abort(FatalError);
419  }
420 
421  // Set zero current time
422  // mesh motion fluxes if the time has been incremented
423  if (!time().subCycling() && phiPtr_->timeIndex() != time().timeIndex())
424  {
425  (*phiPtr_) = dimensionedScalar(dimVolume/dimTime, 0);
426  }
427 
428  return *phiPtr_;
429 }
430 
431 
433 {
434  if (!phiPtr_)
435  {
437  << "mesh flux field does not exist, is the mesh actually moving?"
438  << abort(FatalError);
439  }
440 
441  return *phiPtr_;
442 }
443 
444 
445 // ************************************************************************* //
Foam::surfaceFields.
Graphite solid properties.
Definition: C.H:48
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
const surfaceVectorField & Sf() const
Return cell face area vectors.
const dimensionSet dimArea
bool moving() const
Is mesh moving.
Definition: polyMesh.H:512
The time value with time-stepping information, user-defined remapping, etc.
Definition: TimeState.H:48
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:323
const surfaceVectorField & Cf() const
Return face centres as surfaceVectorField.
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
static word meshSubDir
Return the mesh sub-directory name (usually "polyMesh")
Definition: polyMesh.H:312
IOobject(const word &name, const fileName &instance, const objectRegistry &registry, readOption r=NO_READ, writeOption w=NO_WRITE, bool registerObject=true)
Construct from name, instance, registry, io options.
Definition: IOobject.C:167
static tmp< GeometricField< vector, fvsPatchField, surfaceMesh > > New(const word &name, const Internal &, const PtrList< fvsPatchField< vector >> &)
Return a temporary field constructed from name,.
const surfaceScalarField & phi() const
Return cell face motion fluxes.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
const DimensionedField< scalar, volMesh > & V00() const
Return old-old-time cell volumes.
DimensionedField< scalar, volMesh > & setV0()
Return old-time cell volumes.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:239
const DimensionedField< scalar, volMesh > & V() const
Return cell volumes.
const TimeState & prevTimeState() const
Return previous TimeState if time is being sub-cycled.
Definition: Time.H:428
DimensionedField< Type, GeoMesh > Internal
Type of the internal field from which this GeometricField is derived.
const dimensionSet dimLength
const labelUList & neighbour() const
Internal face neighbour.
Definition: fvMesh.H:284
surfaceScalarField & setPhi()
Return cell face motion fluxes.
const scalarField & magFaceAreas() const
const dimensionSet dimTime
tmp< DimensionedField< scalar, volMesh > > Vsc() const
Return sub-cycle cell volumes.
const fileName & pointsInstance() const
Return the current instance directory for points.
Definition: polyMesh.C:802
const DimensionedField< scalar, volMesh > & V0() const
Return old-time cell volumes.
SlicedGeometricField< vector, fvsPatchField, slicedFvsPatchField, surfaceMesh > slicedSurfaceVectorField
scalar deltaTValue() const
Return time step value.
Definition: TimeStateI.H:41
const Type & value() const
Return const reference to value.
word timeName
Definition: getTimeIndex.H:3
SlicedGeometricField< vector, fvPatchField, slicedFvPatchField, volMesh > slicedVolVectorField
const vectorField & cellCentres() const
errorManip< error > abort(error &err)
Definition: errorManip.H:131
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
Definition: HashTable.H:60
tmp< DimensionedField< scalar, volMesh > > Vsc0() const
Return sub-cycl old-time cell volumes.
const surfaceScalarField & magSf() const
Return cell face area magnitudes.
tmp< surfaceVectorField > delta() const
Return face deltas as surfaceVectorField.
const labelUList & owner() const
Internal face owner.
Definition: fvMesh.H:278
label timeIndex() const
Return the time index of the field.
The internalField of a SlicedGeometricField.
const vectorField & faceCentres() const
label timeIndex() const
Return current time index.
Definition: TimeStateI.H:35
label patchi
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
const vectorField & faceAreas() const
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const dimensionSet dimVolume
const volVectorField & C() const
Return cell centres as volVectorField.
A class for managing temporary objects.
Definition: PtrList.H:53
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
SlicedGeometricField< scalar, fvsPatchField, slicedFvsPatchField, surfaceMesh > slicedSurfaceScalarField
const scalarField & cellVolumes() const
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:540
#define InfoInFunction
Report an information message using Foam::Info.