fvMesh.H
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-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 Class
25  Foam::fvMesh
26 
27 Description
28  Mesh data needed to do the Finite Volume discretisation.
29 
30  NOTE ON USAGE:
31  fvMesh contains all the topological and geometric information
32  related to the mesh. It is also responsible for keeping the data
33  up-to-date. This is done by deleting the cell volume, face area,
34  cell/face centre, addressing and other derived information as
35  required and recalculating it as necessary. The fvMesh therefore
36  reserves the right to delete the derived information upon every
37  topological (mesh refinement/morphing) or geometric change (mesh
38  motion). It is therefore unsafe to keep local references to the
39  derived data outside of the time loop.
40 
41 SourceFiles
42  fvMesh.C
43  fvMeshGeometry.C
44 
45 \*---------------------------------------------------------------------------*/
46 
47 #ifndef fvMesh_H
48 #define fvMesh_H
49 
50 #include "polyMesh.H"
51 #include "lduMesh.H"
52 #include "primitiveMesh.H"
53 #include "fvBoundaryMesh.H"
54 #include "surfaceInterpolation.H"
55 #include "fvSchemes.H"
56 #include "fvSolution.H"
57 #include "data.H"
58 #include "DimensionedField.H"
59 #include "volFieldsFwd.H"
60 #include "surfaceFieldsFwd.H"
61 #include "pointFieldsFwd.H"
62 #include "slicedVolFieldsFwd.H"
63 #include "slicedSurfaceFieldsFwd.H"
64 #include "className.H"
65 
66 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
67 
68 namespace Foam
69 {
70 
71 class fvMeshLduAddressing;
72 class volMesh;
73 class mapDistributePolyMesh;
74 
75 /*---------------------------------------------------------------------------*\
76  Class fvMesh Declaration
77 \*---------------------------------------------------------------------------*/
78 
79 class fvMesh
80 :
81  public polyMesh,
82  public lduMesh,
83  public surfaceInterpolation,
84  public fvSchemes,
85  public fvSolution,
86  public data
87 {
88  // Private Data
89 
90  //- Boundary mesh
91  fvBoundaryMesh boundary_;
92 
93 
94  // Demand-driven data
95 
96  mutable fvMeshLduAddressing* lduPtr_;
97 
98  //- Current time index for cell volumes
99  // Note. The whole mechanism will be replaced once the
100  // dimensionedField is created and the dimensionedField
101  // will take care of the old-time levels.
102  mutable label curTimeIndex_;
103 
104  //- Cell volumes old time level
105  mutable void* VPtr_;
106 
107  //- Cell volumes old time level
108  mutable DimensionedField<scalar, volMesh>* V0Ptr_;
109 
110  //- Cell volumes old-old time level
111  mutable DimensionedField<scalar, volMesh>* V00Ptr_;
112 
113  //- Face area vectors
114  mutable slicedSurfaceVectorField* SfPtr_;
115 
116  //- Mag face area vectors
117  mutable slicedSurfaceScalarField* magSfPtr_;
118 
119  //- Cell centres
120  mutable slicedVolVectorField* CPtr_;
121 
122  //- Face centres
123  mutable slicedSurfaceVectorField* CfPtr_;
124 
125  //- Face motion fluxes
126  mutable surfaceScalarField* phiPtr_;
127 
128 
129  // Private Member Functions
130 
131  // Storage management
132 
133  //- Clear geometry but not the old-time cell volumes
134  void clearGeomNotOldVol();
135 
136  //- Clear geometry like clearGeomNotOldVol but recreate any
137  // geometric demand-driven data that was set
138  void updateGeomNotOldVol();
139 
140  //- Clear geometry
141  void clearGeom();
142 
143  //- Clear addressing
144  void clearAddressing(const bool isMeshUpdate = false);
145 
146  //- Preserve old volume(s)
147  void storeOldVol(const scalarField&);
148 
149 
150  // Make geometric data
151 
152  void makeSf() const;
153  void makeMagSf() const;
154 
155  void makeC() const;
156  void makeCf() const;
157 
158 
159 public:
160 
161  // Public Typedefs
163  typedef fvMesh Mesh;
165 
166 
167  // Declare name of the class and its debug switch
168  ClassName("fvMesh");
169 
170 
171  // Constructors
172 
173  //- Construct from IOobject
174  explicit fvMesh(const IOobject& io);
175 
176  //- Construct from cellShapes with boundary.
177  fvMesh
178  (
179  const IOobject& io,
180  pointField&& points,
181  const cellShapeList& shapes,
182  const faceListList& boundaryFaces,
183  const wordList& boundaryPatchNames,
184  const PtrList<dictionary>& boundaryDicts,
185  const word& defaultBoundaryPatchName,
186  const word& defaultBoundaryPatchType,
187  const bool syncPar = true
188  );
189 
190  //- Construct from components without boundary.
191  // Boundary is added using addFvPatches() member function
192  fvMesh
193  (
194  const IOobject& io,
195  pointField&& points,
196  faceList&& faces,
197  labelList&& allOwner,
198  labelList&& allNeighbour,
199  const bool syncPar = true
200  );
201 
202  //- Construct without boundary from cells rather than owner/neighbour.
203  // Boundary is added using addPatches() member function
204  fvMesh
205  (
206  const IOobject& io,
207  pointField&& points,
208  faceList&& faces,
209  cellList&& cells,
210  const bool syncPar = true
211  );
212 
213  //- Disallow default bitwise copy construction
214  fvMesh(const fvMesh&);
215 
216 
217  //- Destructor
218  virtual ~fvMesh();
219 
220 
221  // Member Functions
222 
223  // Helpers
224 
225  //- Add boundary patches. Constructor helper
226  void addFvPatches
227  (
228  const List<polyPatch*>&,
229  const bool validBoundary = true
230  );
231 
232  //- Update the mesh based on the mesh files saved in time
233  // directories
234  virtual readUpdateState readUpdate();
235 
236 
237  // Access
238 
239  //- Return the top-level database
240  const Time& time() const
241  {
242  return polyMesh::time();
243  }
244 
245  //- Return the object registry - resolve conflict polyMesh/lduMesh
246  virtual const objectRegistry& thisDb() const
247  {
248  return polyMesh::thisDb();
249  }
250 
251  //- Return reference to name
252  // Note: name() is currently ambiguous due to derivation from
253  // surfaceInterpolation
254  const word& name() const
255  {
256  return polyMesh::name();
257  }
258 
259  //- Return reference to boundary mesh
260  const fvBoundaryMesh& boundary() const;
261 
262  //- Return ldu addressing
263  virtual const lduAddressing& lduAddr() const;
264 
265  //- Return a list of pointers for each patch
266  // with only those pointing to interfaces being set
267  virtual lduInterfacePtrsList interfaces() const
268  {
269  return boundary().interfaces();
270  }
271 
272  //- Return communicator used for parallel communication
273  virtual label comm() const
274  {
275  return polyMesh::comm();
276  }
277 
278  //- Internal face owner
279  const labelUList& owner() const
280  {
281  return lduAddr().lowerAddr();
282  }
283 
284  //- Internal face neighbour
285  const labelUList& neighbour() const
286  {
287  return lduAddr().upperAddr();
288  }
289 
290  //- Return cell volumes
291  const DimensionedField<scalar, volMesh>& V() const;
292 
293  //- Return old-time cell volumes
294  const DimensionedField<scalar, volMesh>& V0() const;
295 
296  //- Return old-old-time cell volumes
297  const DimensionedField<scalar, volMesh>& V00() const;
298 
299  //- Return sub-cycle cell volumes
301 
302  //- Return sub-cycl old-time cell volumes
304 
305  //- Return cell face area vectors
306  const surfaceVectorField& Sf() const;
307 
308  //- Return cell face area magnitudes
309  const surfaceScalarField& magSf() const;
310 
311  //- Return cell face motion fluxes
312  const surfaceScalarField& phi() const;
313 
314  //- Return cell centres as volVectorField
315  const volVectorField& C() const;
316 
317  //- Return face centres as surfaceVectorField
318  const surfaceVectorField& Cf() const;
319 
320  //- Return face deltas as surfaceVectorField
322 
323  //- Return a labelType of valid component indicators
324  // 1 : valid (solved)
325  // -1 : invalid (not solved)
326  template<class Type>
327  typename pTraits<Type>::labelType validComponents() const;
328 
329 
330  // Edit
331 
332  //- Clear all geometry and addressing
333  void clearOut();
334 
335  //- Update mesh corresponding to the given map
336  virtual void updateMesh(const mapPolyMesh& mpm);
337 
338  //- Update mesh corresponding to the given distribution map
339  // This is a prototype implementation without field mapping
340  // or handling of old-time volumes for mesh-morphing
341  virtual void updateMesh(const mapDistributePolyMesh& mdpm);
342 
343  //- Move points, returns volumes swept by faces in motion
344  virtual tmp<scalarField> movePoints(const pointField&);
345 
346  //- Map all fields in time using given map.
347  virtual void mapFields(const mapPolyMesh& mpm);
348 
349  //- Add/insert single patch. If validBoundary the new situation
350  // is consistent across processors.
351  virtual void addPatch
352  (
353  const label insertPatchi,
354  const polyPatch& patch,
355  const dictionary& patchFieldDict,
356  const word& defaultPatchFieldType,
357  const bool validBoundary
358  );
359 
360  //- Reorder and trim existing patches. If validBoundary the new
361  // situation is consistent across processors
362  virtual void reorderPatches
363  (
364  const labelUList& newToOld,
365  const bool validBoundary
366  );
367 
368  //- Remove boundary patches. Warning: fvPatchFields hold ref to
369  // these fvPatches.
370  void removeFvBoundary();
371 
372  //- Return cell face motion fluxes
374 
375  //- Return old-time cell volumes
377 
378 
379  // Write
380 
381  //- Write the underlying polyMesh and other data
382  virtual bool writeObject
383  (
387  const bool write = true
388  ) const;
389 
390  //- Write mesh using IO settings from time
391  virtual bool write(const bool write = true) const;
392 
393 
394  // Member Operators
395 
396  //- Disallow default bitwise assignment
397  void operator=(const fvMesh&);
398 
399  bool operator!=(const fvMesh&) const;
400  bool operator==(const fvMesh&) const;
401 };
402 
403 
404 template<>
406 fvMesh::validComponents<sphericalTensor>() const;
407 
408 
409 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
410 
411 } // End namespace Foam
412 
413 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
414 
415 #ifdef NoRepository
416  #include "fvMeshTemplates.C"
417 #endif
418 
419 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
420 
421 #endif
422 
423 // ************************************************************************* //
void clearOut()
Clear all geometry and addressing.
Definition: fvMesh.C:240
const surfaceVectorField & Sf() const
Return cell face area vectors.
void clearAddressing()
Clear topological data.
void removeFvBoundary()
Remove boundary patches. Warning: fvPatchFields hold ref to.
Definition: fvMesh.C:473
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 word & name() const
Return name.
Definition: IOobject.H:303
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
void addFvPatches(const List< polyPatch *> &, const bool validBoundary=true)
Add boundary patches. Constructor helper.
Definition: fvMesh.C:455
const surfaceVectorField & Cf() const
Return face centres as surfaceVectorField.
virtual void addPatch(const label insertPatchi, const polyPatch &patch, const dictionary &patchFieldDict, const word &defaultPatchFieldType, const bool validBoundary)
Add/insert single patch. If validBoundary the new situation.
Definition: fvMesh.C:913
virtual void reorderPatches(const labelUList &newToOld, const bool validBoundary)
Reorder and trim existing patches. If validBoundary the new.
Definition: fvMesh.C:1045
const surfaceScalarField & phi() const
Return cell face motion fluxes.
Traits class for primitives.
Definition: pTraits.H:50
Abstract base class for meshes which provide LDU addressing for the construction of lduMatrix and LDU...
Definition: lduMesh.H:59
bool movePoints()
Do what is necessary if the mesh has moved.
const DimensionedField< scalar, volMesh > & V00() const
Return old-old-time cell volumes.
Class containing mesh-to-mesh mapping information after a mesh distribution where we send parts of me...
fvBoundaryMesh BoundaryMesh
Definition: fvMesh.H:163
const cellList & cells() const
Cell to surface interpolation scheme. Included in fvMesh.
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.
virtual bool write(const bool write=true) const
Write mesh using IO settings from time.
Definition: fvMesh.C:1093
virtual const labelUList & lowerAddr() const =0
Return lower addressing.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:68
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Definition: mapPolyMesh.H:158
virtual bool writeObject(IOstream::streamFormat fmt, IOstream::versionNumber ver, IOstream::compressionType cmp, const bool write=true) const
Write the underlying polyMesh and other data.
Definition: fvMesh.C:1070
virtual const objectRegistry & thisDb() const
Return the object registry - resolve conflict polyMesh/lduMesh.
Definition: fvMesh.H:245
const labelUList & neighbour() const
Internal face neighbour.
Definition: fvMesh.H:284
fvMesh Mesh
Definition: fvMesh.H:162
virtual const pointField & points() const
Return raw points.
Definition: polyMesh.C:1131
surfaceScalarField & setPhi()
Return cell face motion fluxes.
virtual ~fvMesh()
Destructor.
Definition: fvMesh.C:446
bool operator!=(const fvMesh &) const
Definition: fvMesh.C:1109
virtual label comm() const
Return communicator used for parallel communication.
Definition: fvMesh.H:272
tmp< DimensionedField< scalar, volMesh > > Vsc() const
Return sub-cycle cell volumes.
virtual void updateMesh(const mapPolyMesh &mpm)
Update mesh corresponding to the given map.
Definition: fvMesh.C:795
virtual readUpdateState readUpdate()
Update the mesh based on the mesh files saved in time.
Definition: fvMesh.C:489
const DimensionedField< scalar, volMesh > & V0() const
Return old-time cell volumes.
A class for handling words, derived from string.
Definition: word.H:59
virtual const labelUList & upperAddr() const =0
Return upper addressing.
pTraits< Type >::labelType validComponents() const
Return a labelType of valid component indicators.
virtual const lduAddressing & lduAddr() const
Return ldu addressing.
Definition: fvMesh.C:546
streamFormat
Enumeration for the format of data in the stream.
Definition: IOstream.H:86
virtual const faceList & faces() const
Return raw faces.
Definition: polyMesh.C:1156
label comm() const
Return communicator used for parallel communication.
Definition: polyMesh.C:1412
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.
compressionType
Enumeration for the format of data in the stream.
Definition: IOstream.H:193
const surfaceScalarField & magSf() const
Return cell face area magnitudes.
Foam::fvMeshLduAddressing.
tmp< surfaceVectorField > delta() const
Return face deltas as surfaceVectorField.
const Time & time() const
Return time.
Database for solution and other reduced data.
Definition: data.H:51
lduInterfacePtrsList interfaces() const
Return a list of pointers for each patch.
const labelUList & owner() const
Internal face owner.
Definition: fvMesh.H:278
fvMesh(const IOobject &io)
Construct from IOobject.
Definition: fvMesh.C:256
const word & name() const
Return reference to name.
Definition: fvMesh.H:253
Selector class for finite volume solution solution. fvMesh is derived from fvSolution so that all fie...
Definition: fvSolution.H:47
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:70
bool operator==(const fvMesh &) const
Definition: fvMesh.C:1115
Foam::fvBoundaryMesh.
Selector class for finite volume differencing schemes. fvMesh is derived from fvSchemes so that all f...
Definition: fvSchemes.H:50
Macro definitions for declaring ClassName(), NamespaceName(), etc.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Version number type.
Definition: IOstream.H:96
const objectRegistry & thisDb() const
Return the object registry.
Definition: polyMesh.H:497
The class contains the addressing required by the lduMatrix: upper, lower and losort.
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
const volVectorField & C() const
Return cell centres as volVectorField.
readUpdateState
Enumeration defining the state of the mesh after a read update.
Definition: polyMesh.H:88
A class for managing temporary objects.
Definition: PtrList.H:53
Registry of regIOobjects.
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
Specialisation of GeometricField which holds slices of given complete fields in a form that they act ...
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
void operator=(const fvMesh &)
Disallow default bitwise assignment.
virtual lduInterfacePtrsList interfaces() const
Return a list of pointers for each patch.
Definition: fvMesh.H:266
Namespace for OpenFOAM.
virtual void mapFields(const mapPolyMesh &mpm)
Map all fields in time using given map.
Definition: fvMesh.C:557
const fvBoundaryMesh & boundary() const
Return reference to boundary mesh.
Definition: fvMesh.C:540
ClassName("fvMesh")