fvMeshGeometry.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2011-2017 OpenFOAM Foundation
9  Copyright (C) 2020-2021 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
12  This file is part of OpenFOAM.
13 
14  OpenFOAM is free software: you can redistribute it and/or modify it
15  under the terms of the GNU General Public License as published by
16  the Free Software Foundation, either version 3 of the License, or
17  (at your option) any later version.
18 
19  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22  for more details.
23 
24  You should have received a copy of the GNU General Public License
25  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26 
27 \*---------------------------------------------------------------------------*/
28 
29 #include "fvMesh.H"
30 #include "Time.H"
31 #include "volFields.H"
32 #include "surfaceFields.H"
33 #include "slicedVolFields.H"
34 #include "slicedSurfaceFields.H"
35 #include "SubField.H"
36 #include "cyclicFvPatchFields.H"
37 #include "cyclicAMIFvPatchFields.H"
38 
39 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
40 
42 {
43  DebugInFunction << "Assembling face areas" << endl;
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 
55  (
56  IOobject
57  (
58  "S",
60  meshSubDir,
61  *this,
64  false
65  ),
66  *this,
67  dimArea,
68  faceAreas()
69  );
70 
71  SfPtr_->setOriented();
72 }
73 
74 
76 {
77  DebugInFunction << "Assembling mag face areas" << endl;
78 
79  // It is an error to attempt to recalculate
80  // if the pointer is already set
81  if (magSfPtr_)
82  {
84  << "mag face areas already exist"
85  << abort(FatalError);
86  }
87 
88  // Note: Added stabilisation for faces with exactly zero area.
89  // These should be caught on mesh checking but at least this stops
90  // the code from producing Nans.
91  magSfPtr_ = new surfaceScalarField
92  (
93  IOobject
94  (
95  "magSf",
96  pointsInstance(),
97  meshSubDir,
98  *this,
101  false
102  ),
103  mag(Sf()) + dimensionedScalar("vs", dimArea, VSMALL)
104  );
105 }
106 
107 
109 {
110  DebugInFunction << "Assembling cell centres" << endl;
111 
112  // It is an error to attempt to recalculate
113  // if the pointer is already set
114  if (CPtr_)
115  {
117  << "cell centres already exist"
118  << abort(FatalError);
119  }
120 
121  // Construct as slices. Only preserve processor (not e.g. cyclic)
122 
123  CPtr_ = new slicedVolVectorField
124  (
125  IOobject
126  (
127  "C",
128  pointsInstance(),
129  meshSubDir,
130  *this,
133  false
134  ),
135  *this,
136  dimLength,
137  cellCentres(),
138  faceCentres(),
139  true, //preserveCouples
140  true //preserveProcOnly
141  );
142 }
143 
144 
146 {
147  DebugInFunction << "Assembling face centres" << endl;
148 
149  // It is an error to attempt to recalculate
150  // if the pointer is already set
151  if (CfPtr_)
152  {
154  << "face centres already exist"
155  << abort(FatalError);
156  }
157 
158  CfPtr_ = new slicedSurfaceVectorField
159  (
160  IOobject
161  (
162  "Cf",
163  pointsInstance(),
164  meshSubDir,
165  *this,
168  false
169  ),
170  *this,
171  dimLength,
172  faceCentres()
173  );
174 }
175 
176 
177 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
178 
180 {
181  if (!VPtr_)
182  {
184  << "Constructing from primitiveMesh::cellVolumes()" << endl;
185 
187  (
188  IOobject
189  (
190  "V",
191  time().timeName(),
192  *this,
195  false
196  ),
197  *this,
198  dimVolume,
199  cellVolumes()
200  );
201  }
202 
203  return *static_cast<slicedVolScalarField::Internal*>(VPtr_);
204 }
205 
206 
208 {
209  if (!V0Ptr_)
210  {
212  << "V0 is not available"
213  << abort(FatalError);
214  }
215 
216  return *V0Ptr_;
217 }
218 
219 
221 {
222  if (!V0Ptr_)
223  {
225  << "V0 is not available"
226  << abort(FatalError);
227  }
228 
229  return *V0Ptr_;
230 }
231 
232 
234 {
235  if (!V00Ptr_)
236  {
237  DebugInFunction << "Constructing from V0" << endl;
238 
240  (
241  IOobject
242  (
243  "V00",
244  time().timeName(),
245  *this,
248  true
249  ),
250  V0()
251  );
252 
253 
254  // If V00 is used then V0 should be stored for restart
255  V0Ptr_->writeOpt(IOobject::AUTO_WRITE);
256  }
257 
258  return *V00Ptr_;
259 }
260 
261 
263 {
264  if (!steady() && moving() && time().subCycling())
265  {
266  const TimeState& ts = time();
267  const TimeState& ts0 = time().prevTimeState();
268 
269  scalar tFrac =
270  (
271  ts.value() - (ts0.value() - ts0.deltaTValue())
272  )/ts0.deltaTValue();
273 
274  if (tFrac < (1 - SMALL))
275  {
276  return V0() + tFrac*(V() - V0());
277  }
278  else
279  {
280  return V();
281  }
282  }
283  else
284  {
285  return V();
286  }
287 }
288 
289 
291 {
292  if (!steady() && moving() && time().subCycling())
293  {
294  const TimeState& ts = time();
295  const TimeState& ts0 = time().prevTimeState();
296 
297  scalar t0Frac =
298  (
299  (ts.value() - ts.deltaTValue())
300  - (ts0.value() - ts0.deltaTValue())
301  )/ts0.deltaTValue();
302 
303  if (t0Frac > SMALL)
304  {
305  return V0() + t0Frac*(V() - V0());
306  }
307  else
308  {
309  return V0();
310  }
311  }
312  else
313  {
314  return V0();
315  }
316 }
317 
318 
320 {
321  if (!SfPtr_)
322  {
323  makeSf();
324  }
325 
326  return *SfPtr_;
327 }
328 
329 
331 {
332  if (!magSfPtr_)
333  {
334  makeMagSf();
335  }
336 
337  return *magSfPtr_;
338 }
339 
340 
342 {
343  if (!CPtr_)
344  {
345  makeC();
346  }
347 
348  return *CPtr_;
349 }
350 
351 
353 {
354  if (!CfPtr_)
355  {
356  makeCf();
357  }
358 
359  return *CfPtr_;
360 }
361 
362 
364 {
365  DebugInFunction << "Calculating face deltas" << endl;
366 
368  (
370  (
371  IOobject
372  (
373  "delta",
374  pointsInstance(),
375  meshSubDir,
376  *this,
379  false
380  ),
381  *this,
382  dimLength
383  )
384  );
385  surfaceVectorField& delta = tdelta.ref();
386  delta.setOriented();
387 
388  const volVectorField& C = this->C();
389  const labelUList& owner = this->owner();
390  const labelUList& neighbour = this->neighbour();
391 
392  forAll(owner, facei)
393  {
394  delta[facei] = C[neighbour[facei]] - C[owner[facei]];
395  }
396 
397  surfaceVectorField::Boundary& deltabf = delta.boundaryFieldRef();
398 
399  forAll(deltabf, patchi)
400  {
401  deltabf[patchi] = boundary()[patchi].delta();
402  }
403 
404  return tdelta;
405 }
406 
407 
409 {
410  if (!phiPtr_)
411  {
413  << "mesh flux field does not exist, is the mesh actually moving?"
414  << abort(FatalError);
415  }
416 
417  // Set zero current time
418  // mesh motion fluxes if the time has been incremented
419  if (!time().subCycling() && phiPtr_->timeIndex() != time().timeIndex())
420  {
421  (*phiPtr_) = dimensionedScalar(dimVolume/dimTime, Zero);
422  }
423 
424  phiPtr_->setOriented();
425 
426  return *phiPtr_;
427 }
428 
429 
431 {
432  if (!phiPtr_)
433  {
435  << "mesh flux field does not exist, is the mesh actually moving?"
436  << abort(FatalError);
437  }
438 
439  return *phiPtr_;
440 }
441 
442 
443 // ************************************************************************* //
Foam::IOobject::NO_WRITE
Definition: IOobject.H:195
volFields.H
Foam::IOobject
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:169
Foam::slicedSurfaceVectorField
SlicedGeometricField< vector, fvsPatchField, slicedFvsPatchField, surfaceMesh > slicedSurfaceVectorField
Definition: slicedSurfaceFieldsFwd.H:66
Foam::IOobject::AUTO_WRITE
Definition: IOobject.H:194
SubField.H
Foam::fvMesh::makeCf
void makeCf() const
Definition: fvMeshGeometry.C:145
Foam::slicedVolVectorField
SlicedGeometricField< vector, fvPatchField, slicedFvPatchField, volMesh > slicedVolVectorField
Definition: slicedVolFieldsFwd.H:66
Foam::dimLength
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:52
Foam::tmp
A class for managing temporary objects.
Definition: PtrList.H:61
Foam::Zero
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
Foam::fvMesh::Vsc0
tmp< DimensionedField< scalar, volMesh > > Vsc0() const
Return sub-cycle old-time cell volumes.
Definition: fvMeshGeometry.C:290
slicedVolFields.H
Foam::polyMesh::meshSubDir
static word meshSubDir
Return the mesh sub-directory name (usually "polyMesh")
Definition: polyMesh.H:321
Foam::fvMesh::setPhi
surfaceScalarField & setPhi()
Return cell face motion fluxes.
Definition: fvMeshGeometry.C:430
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:369
surfaceFields.H
Foam::surfaceFields.
Foam::dimensioned::value
const Type & value() const
Return const reference to value.
Definition: dimensionedType.C:434
Foam::fvMesh::V00
const DimensionedField< scalar, volMesh > & V00() const
Return old-old-time cell volumes.
Definition: fvMeshGeometry.C:233
Foam::fvMesh::makeC
void makeC() const
Definition: fvMeshGeometry.C:108
C
volScalarField & C
Definition: readThermalProperties.H:102
Foam::TimeState::deltaTValue
scalar deltaTValue() const noexcept
Return time step value.
Definition: TimeStateI.H:43
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
Foam::fvMesh::SfPtr_
slicedSurfaceVectorField * SfPtr_
Face area vectors.
Definition: fvMesh.H:121
Foam::fvMesh::V0
const DimensionedField< scalar, volMesh > & V0() const
Return old-time cell volumes.
Definition: fvMeshGeometry.C:207
Foam::polyMesh::pointsInstance
const fileName & pointsInstance() const
Return the current instance directory for points.
Definition: polyMesh.C:846
Foam::dimTime
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:53
Foam::fvMesh::makeSf
void makeSf() const
Definition: fvMeshGeometry.C:41
Foam::TimeState
The time value with time-stepping information, user-defined remapping, etc.
Definition: TimeState.H:51
Foam::tmp::ref
T & ref() const
Definition: tmpI.H:227
Foam::dimArea
const dimensionSet dimArea(sqr(dimLength))
Definition: dimensionSets.H:59
Foam::fvMesh::magSf
const surfaceScalarField & magSf() const
Return cell face area magnitudes.
Definition: fvMeshGeometry.C:330
DebugInFunction
#define DebugInFunction
Report an information message using Foam::Info.
Definition: messageStream.H:388
Foam::fvMesh::Vsc
tmp< DimensionedField< scalar, volMesh > > Vsc() const
Return sub-cycle cell volumes.
Definition: fvMeshGeometry.C:262
delta
scalar delta
Definition: LISASMDCalcMethod2.H:8
Foam::dimensionedScalar
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Definition: dimensionedScalarFwd.H:42
Foam::fvMesh::C
const volVectorField & C() const
Return cell centres as volVectorField.
Definition: fvMeshGeometry.C:341
timeName
word timeName
Definition: getTimeIndex.H:3
Foam::FatalError
error FatalError
fvMesh.H
Foam::abort
errorManip< error > abort(error &err)
Definition: errorManip.H:144
Foam::fvMesh::phi
const surfaceScalarField & phi() const
Return cell face motion fluxes.
Definition: fvMeshGeometry.C:408
Time.H
FatalErrorInFunction
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
Foam::surfaceScalarField
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
Definition: surfaceFieldsFwd.H:54
Foam::mag
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Foam::UList< label >
Foam::SlicedGeometricField::Internal
The internalField of a SlicedGeometricField.
Definition: SlicedGeometricField.H:196
slicedSurfaceFields.H
Foam::fvMesh::setV0
DimensionedField< scalar, volMesh > & setV0()
Return old-time cell volumes.
Definition: fvMeshGeometry.C:220
Foam::C
Graphite solid properties.
Definition: C.H:50
Foam::dimVolume
const dimensionSet dimVolume(pow3(dimLength))
Definition: dimensionSets.H:60
Foam::GeometricField< vector, fvsPatchField, surfaceMesh >
Foam::IOobject::NO_READ
Definition: IOobject.H:188
cyclicFvPatchFields.H
Foam::fvMesh::Cf
const surfaceVectorField & Cf() const
Return face centres as surfaceVectorField.
Definition: fvMeshGeometry.C:352
boundary
faceListList boundary
Definition: createBlockMesh.H:4
Foam::fvMesh::delta
tmp< surfaceVectorField > delta() const
Return face deltas as surfaceVectorField.
Definition: fvMeshGeometry.C:363
Foam::DimensionedField
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Definition: DimensionedField.H:54
cyclicAMIFvPatchFields.H
Foam::fvMesh::makeMagSf
void makeMagSf() const
Definition: fvMeshGeometry.C:75
Foam::fvMesh::V
const DimensionedField< scalar, volMesh > & V() const
Return cell volumes.
Definition: fvMeshGeometry.C:179
Foam::primitiveMesh::faceAreas
const vectorField & faceAreas() const
Definition: primitiveMeshFaceCentresAndAreas.C:89
Foam::fvMesh::Sf
const surfaceVectorField & Sf() const
Return cell face area vectors.
Definition: fvMeshGeometry.C:319