polyMesh.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-2019 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::polyMesh
26 
27 Description
28  Mesh consisting of general polyhedral cells.
29 
30 SourceFiles
31  polyMesh.C
32  polyMeshInitMesh.C
33  polyMeshClear.C
34  polyMeshFromShapeMesh.C
35  polyMeshIO.C
36  polyMeshUpdate.C
37  polyMeshCheck.C
38 
39 \*---------------------------------------------------------------------------*/
40 
41 #ifndef polyMesh_H
42 #define polyMesh_H
43 
44 #include "objectRegistry.H"
45 #include "primitiveMesh.H"
46 #include "pointField.H"
47 #include "faceList.H"
48 #include "cellList.H"
49 #include "cellShapeList.H"
50 #include "pointIOField.H"
51 #include "faceIOList.H"
52 #include "labelIOList.H"
53 #include "polyBoundaryMesh.H"
54 #include "boundBox.H"
55 #include "pointZoneMesh.H"
56 #include "faceZoneMesh.H"
57 #include "cellZoneMesh.H"
58 
59 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
60 
61 namespace Foam
62 {
63 
64 // Forward declaration of classes
65 class globalMeshData;
66 class mapPolyMesh;
67 class polyMeshTetDecomposition;
68 class treeDataCell;
69 template<class Type> class indexedOctree;
70 
71 /*---------------------------------------------------------------------------*\
72  Class polyMesh Declaration
73 \*---------------------------------------------------------------------------*/
74 
75 class polyMesh
76 :
77  public objectRegistry,
78  public primitiveMesh
79 {
80 
81 public:
82 
83  // Public data types
84 
85  //- Enumeration defining the state of the mesh after a read update.
86  // Used for post-processing applications, where the mesh
87  // needs to update based on the files written in time
88  // directories
89  enum readUpdateState
90  {
95  };
96 
97  //- Enumeration defining the decomposition of the cell for
98  // inside/outside test
100  {
101  FACE_PLANES, //- Faces considered as planes
103  FACE_CENTRE_TRIS, //- Faces decomposed into triangles
104  // using face-centre
106  FACE_DIAG_TRIS, //- Faces decomposed into triangles diagonally
108  CELL_TETS //- Cell decomposed into tets
109  };
110 
111 
112 private:
113 
114  // Permanent data
115 
116  // Primitive mesh data
117 
118  //- Points
119  pointIOField points_;
120 
121  //- Faces
122  faceCompactIOList faces_;
123 
124  //- Face owner
125  labelIOList owner_;
126 
127  //- Face neighbour
128  labelIOList neighbour_;
129 
130  //- Have the primitives been cleared
131  bool clearedPrimitives_;
132 
133 
134  //- Boundary mesh
135  mutable polyBoundaryMesh boundary_;
136 
137  //- Mesh bounding-box.
138  // Created from points on construction, updated when the mesh moves
139  boundBox bounds_;
140 
141  //- Communicator used for parallel communication
142  label comm_;
143 
144  //- Vector of non-constrained directions in mesh
145  // defined according to the presence of empty and wedge patches
146  mutable Vector<label> geometricD_;
147 
148  //- Vector of valid directions in mesh
149  // defined according to the presence of empty patches
150  mutable Vector<label> solutionD_;
151 
152  //- Base point for face decomposition into tets
153  mutable autoPtr<labelIOList> tetBasePtIsPtr_;
154 
155  //- Search tree to allow spatial cell searching
156  mutable autoPtr<indexedOctree<treeDataCell>> cellTreePtr_;
157 
158 
159  // Zoning information
160 
161  //- Point zones
162  pointZoneMesh pointZones_;
163 
164  //- Face zones
165  faceZoneMesh faceZones_;
166 
167  //- Cell zones
168  cellZoneMesh cellZones_;
169 
170 
171  //- Parallel info
172  mutable autoPtr<globalMeshData> globalMeshDataPtr_;
173 
174 
175  // Mesh motion related data
176 
177  //- Is the mesh moving
178  bool moving_;
179 
180  //- Is the mesh topology changing
181  bool topoChanging_;
182 
183  //- Current time index for mesh motion
184  mutable label curMotionTimeIndex_;
185 
186  //- Old points (for the last mesh motion)
187  mutable autoPtr<pointField> oldPointsPtr_;
188 
189  //- Old cell centres (for the last mesh motion)
190  mutable autoPtr<pointField> oldCellCentresPtr_;
191 
192  //- Whether or not to store the old cell centres
193  mutable bool storeOldCellCentres_;
194 
195 
196  // Private Member Functions
197 
198  //- Initialise the polyMesh from the primitive data
199  void initMesh();
200 
201  //- Initialise the polyMesh from the given set of cells
202  void initMesh(cellList& c);
203 
204  //- Calculate the valid directions in the mesh from the boundaries
205  void calcDirections() const;
206 
207  //- Calculate the cell shapes from the primitive
208  // polyhedral information
209  void calcCellShapes() const;
210 
211  //- Read and return the tetBasePtIs
212  autoPtr<labelIOList> readTetBasePtIs() const;
213 
214 
215  // Helper functions for constructor from cell shapes
216 
217  labelListList cellShapePointCells(const cellShapeList&) const;
218 
219  labelList facePatchFaceCells
220  (
221  const faceList& patchFaces,
222  const labelListList& pointCells,
223  const faceListList& cellsFaceShapes,
224  const label patchID
225  ) const;
226 
227  void setTopology
228  (
229  const cellShapeList& cellsAsShapes,
230  const faceListList& boundaryFaces,
231  const wordList& boundaryPatchNames,
232  labelList& patchSizes,
233  labelList& patchStarts,
234  label& defaultPatchStart,
235  label& nFaces,
236  cellList& cells
237  );
238 
239 
240  // Geometry checks
241 
242  //- Check non-orthogonality
243  bool checkFaceOrthogonality
244  (
245  const vectorField& fAreas,
246  const vectorField& cellCtrs,
247  const bool report,
248  const bool detailedReport,
249  labelHashSet* setPtr
250  ) const;
251 
252  //- Check face skewness
253  bool checkFaceSkewness
254  (
255  const pointField& points,
256  const vectorField& fCtrs,
257  const vectorField& fAreas,
258  const vectorField& cellCtrs,
259  const bool report,
260  const bool detailedReport,
261  labelHashSet* setPtr
262  ) const;
263 
264  bool checkEdgeAlignment
265  (
266  const pointField& p,
267  const bool report,
268  const Vector<label>& directions,
269  labelHashSet* setPtr
270  ) const;
271 
272  bool checkCellDeterminant
273  (
274  const vectorField& faceAreas,
275  const bool report,
276  labelHashSet* setPtr,
277  const Vector<label>& meshD
278  ) const;
279 
280  bool checkFaceWeight
281  (
282  const vectorField& fCtrs,
283  const vectorField& fAreas,
284  const vectorField& cellCtrs,
285  const bool report,
286  const scalar minWeight,
287  labelHashSet* setPtr
288  ) const;
289 
290  bool checkVolRatio
291  (
292  const scalarField& cellVols,
293  const bool report,
294  const scalar minRatio,
295  labelHashSet* setPtr
296  ) const;
297 
298 public:
299 
300  // Public Typedefs
302  typedef polyMesh Mesh;
304 
305 
306  //- Runtime type information
307  TypeName("polyMesh");
308 
309  //- Return the default region name
310  static word defaultRegion;
311 
312  //- Return the mesh sub-directory name (usually "polyMesh")
313  static word meshSubDir;
314 
315 
316  // Constructors
317 
318  //- Construct from IOobject
319  explicit polyMesh(const IOobject& io);
320 
321  //- Move construct from IOobject or from components.
322  // Boundary is added using addPatches() member function
323  polyMesh
324  (
325  const IOobject& io,
326  pointField&& points,
327  faceList&& faces,
328  labelList&& owner,
329  labelList&& neighbour,
330  const bool syncPar = true
331  );
332 
333  //- Move construct without boundary with cells rather than
334  // owner/neighbour.
335  // Boundary is added using addPatches() member function
336  polyMesh
337  (
338  const IOobject& io,
339  pointField&& points,
340  faceList&& faces,
341  cellList&& cells,
342  const bool syncPar = true
343  );
344 
345  //- Move construct from cell shapes
346  polyMesh
347  (
348  const IOobject& io,
349  pointField&& points,
350  const cellShapeList& shapes,
351  const faceListList& boundaryFaces,
352  const wordList& boundaryPatchNames,
353  const wordList& boundaryPatchTypes,
354  const word& defaultBoundaryPatchName,
355  const word& defaultBoundaryPatchType,
356  const wordList& boundaryPatchPhysicalTypes,
357  const bool syncPar = true
358  );
359 
360  //- Move construct from cell shapes with patch information in dictionary
361  // format.
362  polyMesh
363  (
364  const IOobject& io,
365  pointField&& points,
366  const cellShapeList& shapes,
367  const faceListList& boundaryFaces,
368  const wordList& boundaryPatchNames,
369  const PtrList<dictionary>& boundaryDicts,
370  const word& defaultBoundaryPatchName,
371  const word& defaultBoundaryPatchType,
372  const bool syncPar = true
373  );
374 
375  //- Disallow default bitwise copy construction
376  polyMesh(const polyMesh&);
377 
378 
379  //- Destructor
380  virtual ~polyMesh();
381 
382 
383  // Member Functions
384 
385  // Database
386 
387  //- Override the objectRegistry dbDir for a single-region case
388  virtual const fileName& dbDir() const;
389 
390  //- Return the local mesh directory (dbDir()/meshSubDir)
391  fileName meshDir() const;
392 
393  //- Return the current instance directory for points
394  // Used in the consruction of gemometric mesh data dependent
395  // on points
396  const fileName& pointsInstance() const;
397 
398  //- Return the current instance directory for faces
399  const fileName& facesInstance() const;
400 
401  //- Set the instance for mesh files
402  void setInstance(const fileName&);
403 
404 
405  // Access
406 
407  //- Return raw points
408  virtual const pointField& points() const;
409 
410  //- Return true if io is up-to-date with points
411  virtual bool upToDatePoints(const regIOobject& io) const;
412 
413  //- Set io to be up-to-date with points
414  virtual void setUpToDatePoints(regIOobject& io) const;
415 
416  //- Return raw faces
417  virtual const faceList& faces() const;
418 
419  //- Return face owner
420  virtual const labelList& faceOwner() const;
421 
422  //- Return face neighbour
423  virtual const labelList& faceNeighbour() const;
424 
425  //- Return old points for mesh motion
426  virtual const pointField& oldPoints() const;
427 
428  //- Return old points for mesh motion
429  virtual const pointField& oldCellCentres() const;
430 
431  //- Return IO object for points0
432  static IOobject points0IO(const polyMesh& mesh);
433 
434  //- Return boundary mesh
435  const polyBoundaryMesh& boundaryMesh() const
436  {
437  return boundary_;
438  }
439 
440  //- Return mesh bounding box
441  const boundBox& bounds() const
442  {
443  return bounds_;
444  }
445 
446  //- Return the vector of geometric directions in mesh.
447  // Defined according to the presence of empty and wedge patches.
448  // 1 indicates unconstrained direction and -1 a constrained
449  // direction.
450  const Vector<label>& geometricD() const;
451 
452  //- Return the number of valid geometric dimensions in the mesh
453  label nGeometricD() const;
454 
455  //- Return the vector of solved-for directions in mesh.
456  // Differs from geometricD in that it includes for wedge cases
457  // the circumferential direction in case of swirl.
458  // 1 indicates valid direction and -1 an invalid direction.
459  const Vector<label>& solutionD() const;
460 
461  //- Return the number of valid solved-for dimensions in the mesh
462  label nSolutionD() const;
463 
464  //- Return the tetBasePtIs
465  const labelIOList& tetBasePtIs() const;
466 
467  //- Return the cell search tree
468  const indexedOctree<treeDataCell>& cellTree() const;
469 
470  //- Return point zone mesh
471  const pointZoneMesh& pointZones() const
472  {
473  return pointZones_;
474  }
475 
476  //- Return face zone mesh
477  const faceZoneMesh& faceZones() const
478  {
479  return faceZones_;
480  }
481 
482  //- Return cell zone mesh
483  const cellZoneMesh& cellZones() const
484  {
485  return cellZones_;
486  }
487 
488  //- Return parallel info
489  const globalMeshData& globalData() const;
490 
491  //- Return communicator used for parallel communication
492  label comm() const;
493 
494  //- Return communicator used for parallel communication
495  label& comm();
496 
497  //- Return the object registry
498  const objectRegistry& thisDb() const
499  {
500  return *this;
501  }
502 
503 
504  // Mesh motion
505 
506  //- Is mesh dynamic
507  virtual bool dynamic() const
508  {
509  return false;
510  }
511 
512  //- Is mesh moving
513  bool moving() const
514  {
515  return moving_;
516  }
517 
518  //- Set the mesh to be moving
519  bool moving(const bool m)
520  {
521  bool m0 = moving_;
522  moving_ = m;
523  return m0;
524  }
525 
526  //- Is mesh topology changing
527  bool topoChanging() const
528  {
529  return topoChanging_;
530  }
531 
532  //- Set the mesh topology to be changing
533  bool topoChanging(const bool c)
534  {
535  bool c0 = topoChanging_;
536  topoChanging_ = c;
537  return c0;
538  }
539 
540  //- Is mesh changing (topology changing and/or moving)
541  bool changing() const
542  {
543  return moving()||topoChanging();
544  }
545 
546  //- Move points, returns volumes swept by faces in motion
547  virtual tmp<scalarField> movePoints(const pointField&);
548 
549  //- Reset motion
550  void resetMotion() const;
551 
552 
553  // Topological change
554 
555  //- Return non-const access to the pointZones
557  {
558  return pointZones_;
559  }
560 
561  //- Return non-const access to the faceZones
563  {
564  return faceZones_;
565  }
566 
567  //- Return non-const access to the cellZones
569  {
570  return cellZones_;
571  }
572 
573  //- Add boundary patches
574  void addPatches
575  (
576  const List<polyPatch*>&,
577  const bool validBoundary = true
578  );
579 
580  //- Add mesh zones
581  void addZones
582  (
583  const List<pointZone*>& pz,
584  const List<faceZone*>& fz,
585  const List<cellZone*>& cz
586  );
587 
588  //- Add/insert single patch. If validBoundary the new situation
589  // is consistent across processors.
590  virtual void addPatch
591  (
592  const label insertPatchi,
593  const polyPatch& patch,
594  const dictionary& patchFieldDict,
595  const word& defaultPatchFieldType,
596  const bool validBoundary
597  );
598 
599  //- Reorder and trim existing patches. If validBoundary the new
600  // situation is consistent across processors
601  virtual void reorderPatches
602  (
603  const labelUList& newToOld,
604  const bool validBoundary
605  );
606 
607  //- Update the mesh based on the mesh files saved in
608  // time directories
609  virtual readUpdateState readUpdate();
610 
611  //- Update the mesh corresponding to given map
612  virtual void updateMesh(const mapPolyMesh& mpm);
613 
614  //- Remove boundary patches
615  void removeBoundary();
616 
617  //- Reset mesh primitive data. Assumes all patch info correct
618  // (so does e.g. parallel communication). If not use
619  // validBoundary=false
620  void resetPrimitives
621  (
622  pointField&& points,
623  faceList&& faces,
624  labelList&& owner,
625  labelList&& neighbour,
626  const labelList& patchSizes,
627  const labelList& patchStarts,
628  const bool validBoundary = true
629  );
630 
631 
632  // Storage management
633 
634  //- Clear geometry
635  void clearGeom();
636 
637  //- Clear addressing
638  void clearAddressing(const bool isMeshUpdate = false);
639 
640  //- Clear all geometry and addressing unnecessary for CFD
641  void clearOut();
642 
643  //- Clear primitive data (points, faces and cells)
644  void clearPrimitives();
645 
646  //- Clear tet base points
647  void clearTetBasePtIs();
648 
649  //- Clear cell tree data
650  void clearCellTree();
651 
652  //- Remove all files from mesh instance
653  void removeFiles(const fileName& instanceDir) const;
654 
655  //- Remove all files from mesh instance()
656  void removeFiles() const;
657 
658 
659  // Geometric checks. Selectively override primitiveMesh functionality.
660 
661  //- Check non-orthogonality
662  virtual bool checkFaceOrthogonality
663  (
664  const bool report = false,
665  labelHashSet* setPtr = nullptr
666  ) const;
667 
668  //- Check face skewness
669  virtual bool checkFaceSkewness
670  (
671  const bool report = false,
672  labelHashSet* setPtr = nullptr
673  ) const;
674 
675  //- Check edge alignment for 1D/2D cases
676  virtual bool checkEdgeAlignment
677  (
678  const bool report,
679  const Vector<label>& directions,
680  labelHashSet* setPtr
681  ) const;
682 
683  virtual bool checkCellDeterminant
684  (
685  const bool report,
686  labelHashSet* setPtr
687  ) const;
688 
689  //- Check mesh motion for correctness given motion points
690  virtual bool checkMeshMotion
691  (
692  const pointField& newPoints,
693  const bool report = false,
694  const bool detailedReport = false
695  ) const;
696 
697  //- Check for face weights
698  virtual bool checkFaceWeight
699  (
700  const bool report,
701  const scalar minWeight = 0.05,
702  labelHashSet* setPtr = nullptr
703  ) const;
704 
705  //- Check for neighbouring cell volumes
706  virtual bool checkVolRatio
707  (
708  const bool report,
709  const scalar minRatio = 0.01,
710  labelHashSet* setPtr = nullptr
711  ) const;
712 
713 
714  // Position search functions
715 
716  //- Find the cell, tetFacei and tetPti for point p
717  void findCellFacePt
718  (
719  const point& p,
720  label& celli,
721  label& tetFacei,
722  label& tetPti
723  ) const;
724 
725  //- Find the tetFacei and tetPti for point p in celli.
726  // tetFacei and tetPtI are set to -1 if not found
727  void findTetFacePt
728  (
729  const label celli,
730  const point& p,
731  label& tetFacei,
732  label& tetPti
733  ) const;
734 
735  //- Test if point p is in the celli
736  bool pointInCell
737  (
738  const point& p,
739  label celli,
741  ) const;
742 
743  //- Find cell enclosing this location and return index
744  // If not found -1 is returned
746  (
747  const point& p,
749  ) const;
750 
751 
752  // Member Operators
753 
754  //- Disallow default bitwise assignment
755  void operator=(const polyMesh&) = delete;
756 };
757 
758 
759 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
760 
761 } // End namespace Foam
762 
763 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
764 
765 #endif
766 
767 // ************************************************************************* //
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:434
virtual void reorderPatches(const labelUList &newToOld, const bool validBoundary)
Reorder and trim existing patches. If validBoundary the new.
Definition: polyMesh.C:1011
bool changing() const
Is mesh changing (topology changing and/or moving)
Definition: polyMesh.H:540
polyMesh Mesh
Definition: polyMesh.H:301
virtual tmp< scalarField > movePoints(const pointField &)
Move points, returns volumes swept by faces in motion.
Definition: polyMesh.C:1287
void clearAddressing()
Clear topological data.
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
A class for handling file names.
Definition: fileName.H:79
const faceZoneMesh & faceZones() const
Return face zone mesh.
Definition: polyMesh.H:476
bool moving() const
Is mesh moving.
Definition: polyMesh.H:512
const fileName & facesInstance() const
Return the current instance directory for faces.
Definition: polyMesh.C:808
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
virtual const labelList & faceNeighbour() const
Return face neighbour.
Definition: polyMesh.C:1175
const labelIOList & tetBasePtIs() const
Return the tetBasePtIs.
Definition: polyMesh.C:848
virtual void setUpToDatePoints(regIOobject &io) const
Set io to be up-to-date with points.
Definition: polyMesh.C:1150
Cell-face mesh analysis engine.
Definition: primitiveMesh.H:74
void resetMotion() const
Reset motion.
Definition: polyMesh.C:1386
label nFaces() const
Various mesh related information for a parallel run. Upon construction, constructs all info using par...
void clearOut()
Clear all geometry and addressing unnecessary for CFD.
static word meshSubDir
Return the mesh sub-directory name (usually "polyMesh")
Definition: polyMesh.H:312
const Vector< label > & solutionD() const
Return the vector of solved-for directions in mesh.
Definition: polyMesh.C:831
Set of directions for each cell in the mesh. Either uniform and size=1 or one set of directions per c...
Definition: directions.H:64
static word defaultRegion
Return the default region name.
Definition: polyMesh.H:309
Foam::pointZoneMesh.
A bounding box defined in terms of the points at its extremities.
Definition: boundBox.H:58
const cellList & cells() const
virtual const fileName & dbDir() const
Override the objectRegistry dbDir for a single-region case.
Definition: polyMesh.C:783
bool pointInCell(const point &p, label celli, const cellDecomposition=CELL_TETS) const
Test if point p is in the celli.
Definition: polyMesh.C:1496
cellDecomposition
Enumeration defining the decomposition of the cell for.
Definition: polyMesh.H:98
label nGeometricD() const
Return the number of valid geometric dimensions in the mesh.
Definition: polyMesh.C:825
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Definition: mapPolyMesh.H:158
virtual const pointField & points() const
Return raw points.
Definition: polyMesh.C:1131
const Vector< label > & geometricD() const
Return the vector of geometric directions in mesh.
Definition: polyMesh.C:814
Foam::cellZoneMesh.
const fileName & pointsInstance() const
Return the current instance directory for points.
Definition: polyMesh.C:802
dynamicFvMesh & mesh
void removeFiles() const
Remove all files from mesh instance()
Definition: polyMesh.C:1448
void clearPrimitives()
Clear primitive data (points, faces and cells)
virtual void updateMesh(const mapPolyMesh &mpm)
Update the mesh corresponding to given map.
const pointZoneMesh & pointZones() const
Return point zone mesh.
Definition: polyMesh.H:470
A class for handling words, derived from string.
Definition: word.H:59
const cellZoneMesh & cellZones() const
Return cell zone mesh.
Definition: polyMesh.H:482
void findCellFacePt(const point &p, label &celli, label &tetFacei, label &tetPti) const
Find the cell, tetFacei and tetPti for point p.
Definition: polyMesh.C:1455
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: polyMesh.C:1049
virtual const pointField & oldPoints() const
Return old points for mesh motion.
Definition: polyMesh.C:1181
virtual const labelList & faceOwner() const
Return face owner.
Definition: polyMesh.C:1169
const globalMeshData & globalData() const
Return parallel info.
Definition: polyMesh.C:1394
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
void addZones(const List< pointZone *> &pz, const List< faceZone *> &fz, const List< cellZone *> &cz)
Add mesh zones.
Definition: polyMesh.C:953
void addPatches(const List< polyPatch *> &, const bool validBoundary=true)
Add boundary patches.
Definition: polyMesh.C:909
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
const labelListList & pointCells() const
Foam::polyBoundaryMesh.
void clearGeom()
Clear geometry.
Definition: polyMeshClear.C:53
void removeBoundary()
Remove boundary patches.
Definition: polyMeshClear.C:36
void resetPrimitives(pointField &&points, faceList &&faces, labelList &&owner, labelList &&neighbour, const labelList &patchSizes, const labelList &patchStarts, const bool validBoundary=true)
Reset mesh primitive data. Assumes all patch info correct.
Definition: polyMesh.C:672
const scalarField & cellVols
void clearTetBasePtIs()
Clear tet base points.
void setInstance(const fileName &)
Set the instance for mesh files.
Definition: polyMeshIO.C:32
static IOobject points0IO(const polyMesh &mesh)
Return IO object for points0.
Definition: polyMesh.C:1220
label nSolutionD() const
Return the number of valid solved-for dimensions in the mesh.
Definition: polyMesh.C:842
virtual bool upToDatePoints(const regIOobject &io) const
Return true if io is up-to-date with points.
Definition: polyMesh.C:1144
void operator=(const polyMesh &)=delete
Disallow default bitwise assignment.
TypeName("polyMesh")
Runtime type information.
const indexedOctree< treeDataCell > & cellTree() const
Return the cell search tree.
Definition: polyMesh.C:882
const boundBox & bounds() const
Return mesh bounding box.
Definition: polyMesh.H:440
virtual const pointField & oldCellCentres() const
Return old points for mesh motion.
Definition: polyMesh.C:1199
Non-pointer based hierarchical recursive searching.
Definition: treeDataEdge.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
const vectorField & faceAreas() const
const dimensionedScalar c
Speed of light in a vacuum.
bool topoChanging() const
Is mesh topology changing.
Definition: polyMesh.H:526
fileName meshDir() const
Return the local mesh directory (dbDir()/meshSubDir)
Definition: polyMesh.C:796
virtual bool dynamic() const
Is mesh dynamic.
Definition: polyMesh.H:506
void findTetFacePt(const label celli, const point &p, label &tetFacei, label &tetPti) const
Find the tetFacei and tetPti for point p in celli.
Definition: polyMesh.C:1480
label findCell(const point &p, const cellDecomposition=CELL_TETS) const
Find cell enclosing this location and return index.
Definition: polyMesh.C:1606
polyBoundaryMesh BoundaryMesh
Definition: polyMesh.H:302
virtual ~polyMesh()
Destructor.
Definition: polyMesh.C:774
regIOobject is an abstract class derived from IOobject to handle automatic object registration with t...
Definition: regIOobject.H:55
const objectRegistry & thisDb() const
Return the object registry.
Definition: polyMesh.H:497
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
readUpdateState
Enumeration defining the state of the mesh after a read update.
Definition: polyMesh.H:88
volScalarField & p
A class for managing temporary objects.
Definition: PtrList.H:53
virtual readUpdateState readUpdate()
Update the mesh based on the mesh files saved in.
Definition: polyMeshIO.C:71
Registry of regIOobjects.
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
Foam::faceZoneMesh.
A primitive field of type <T> with automated input and output.
Definition: IOField.H:50
void clearCellTree()
Clear cell tree data.
virtual bool checkMeshMotion(const pointField &newPoints, const bool report=false, const bool detailedReport=false) const
Check mesh motion for correctness given motion points.
Namespace for OpenFOAM.
polyMesh(const IOobject &io)
Construct from IOobject.
Definition: polyMesh.C:163