conformalVoronoiMesh.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2012-2016 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::conformalVoronoiMesh
26 
27 Description
28 
29 SourceFiles
30  conformalVoronoiMeshI.H
31  conformalVoronoiMesh.C
32  conformalVoronoiMeshZones.C
33  conformalVoronoiMeshIO.C
34  conformalVoronoiMeshConformToSurface.C
35  conformalVoronoiMeshFeaturePoints.C
36  conformalVoronoiMeshCalcDualMesh.C
37  conformalVoronoiMeshTemplates.C
38 
39 \*---------------------------------------------------------------------------*/
40 
41 #ifndef conformalVoronoiMesh_H
42 #define conformalVoronoiMesh_H
43 
44 // Include uint.H before CGAL headers to define __STDC_LIMIT_MACROS
45 #include "uint.H"
47 #include "searchableSurfaces.H"
48 #include "conformationSurfaces.H"
49 #include "cellShapeControl.H"
50 #include "cvControls.H"
51 #include "DynamicList.H"
52 #include "PackedBoolList.H"
53 #include "Time.H"
54 #include "polyMesh.H"
55 #include "plane.H"
56 #include "SortableList.H"
57 #include "meshTools.H"
58 #include "dynamicIndexedOctree.H"
59 #include "dynamicTreeDataPoint.H"
60 #include "indexedOctree.H"
61 #include "treeDataPoint.H"
62 #include "unitConversion.H"
63 #include "transform.H"
64 #include "volFields.H"
65 #include "fvMesh.H"
66 #include "labelPair.H"
67 #include "HashSet.H"
68 #include "memInfo.H"
69 #include "point.H"
70 #include "cellSet.H"
71 #include "wallPolyPatch.H"
72 #include "processorPolyPatch.H"
74 #include "globalIndex.H"
75 #include "pointFeatureEdgesTypes.H"
76 #include "pointConversion.H"
77 #include "Pair.H"
79 #include "featurePointConformer.H"
80 #include "pointPairs.H"
81 
82 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
83 
84 namespace Foam
85 {
86 
87 // Forward declaration of classes
88 class initialPointsMethod;
89 class relaxationModel;
90 class faceAreaWeightModel;
91 class backgroundMeshDecomposition;
92 class OBJstream;
93 
94 /*---------------------------------------------------------------------------*\
95  Class conformalVoronoiMesh Declaration
96 \*---------------------------------------------------------------------------*/
97 
99 :
100  public DistributedDelaunayMesh<Delaunay>
101 {
102 public:
104  typedef Delaunay::Vertex_handle Vertex_handle;
105  typedef Delaunay::Cell_handle Cell_handle;
106  typedef Delaunay::Edge Edge;
107  typedef Delaunay::Facet Facet;
108  typedef Delaunay::Point Point;
115 
116  // Static data
118  enum dualMeshPointType
119  {
120  internal = 0,
121  surface = 1,
124  constrained = 4
125  };
128 
129 
130 private:
131 
132  // Static data
133 
134  static const scalar searchConeAngle;
135 
136  static const scalar searchAngleOppositeSurface;
137 
138 
139  // Private data
140 
141  //- The time registry of the application
142  const Time& runTime_;
143 
144  //- Random number generator
145  mutable Random rndGen_;
146 
147  //- Controls for the conformal Voronoi meshing process
148  cvControls foamyHexMeshControls_;
149 
150  //- All geometry of the meshing process, including surfaces to be
151  // conformed to and those to be used for refinement
152  searchableSurfaces allGeometry_;
153 
154  //- The surfaces to conform to
155  conformationSurfaces geometryToConformTo_;
156 
157  //- Background mesh decomposition, only available in parallel.
159 
160  //- The cell shape control object
161  cellShapeControl cellShapeControl_;
162 
163  //- Limiting bound box before infinity begins
164  treeBoundBox limitBounds_;
165 
166  //-
167  mutable pointPairs<Delaunay> ptPairs_;
168 
169  //-
170  featurePointConformer ftPtConformer_;
171 
172  //- Search tree for edge point locations
174  edgeLocationTreePtr_;
175 
176  mutable DynamicList<Foam::point> existingEdgeLocations_;
177 
178  //- Search tree for surface point locations
180  surfacePtLocationTreePtr_;
181 
182  mutable DynamicList<Foam::point> existingSurfacePtLocations_;
183 
184  //- Store the surface and feature edge conformation locations to be
185  // reinserted
186  List<Vb> surfaceConformationVertices_;
187 
188  //- Method for inserting initial points. Runtime selectable.
189  autoPtr<initialPointsMethod> initialPointsMethod_;
190 
191  //- Relaxation coefficient model. Runtime selectable.
192  autoPtr<relaxationModel> relaxationModel_;
193 
194  //- Face area weight function. Runtime selectable.
195  autoPtr<faceAreaWeightModel> faceAreaWeightModel_;
196 
197 
198  // Private Member Functions
199 
200  inline scalar defaultCellSize() const;
201 
202  //- Return the local target cell size at the given location. Takes
203  // boolean argument to allow speed-up of queries if the point is going
204  // to be on a surface.
205  inline scalar targetCellSize(const Foam::point& pt) const;
206 
207  //- Return the target cell size from that stored on a pair of
208  // Delaunay vertices, including the possibility that one of
209  // them is not an internalOrBoundaryPoint, and so will not
210  // have valid data.
211  inline scalar averageAnyCellSize
212  (
213  const Vertex_handle& vA,
214  const Vertex_handle& vB
215  ) const;
216 
217  //- The average target cell size of a Delaunay facet, i.e., of
218  // a dual edge
219  inline scalar averageAnyCellSize
220  (
221  const Delaunay::Finite_facets_iterator& fit
222  ) const;
223 
224  //- Insert Delaunay vertices using the CGAL range insertion method,
225  // optionally check processor occupancy and distribute to other
226  // processors
227  void insertInternalPoints
228  (
230  const bool distribute = false
231  );
232 
233  Map<label> insertPointPairs
234  (
235  List<Vb>& vertices,
236  bool distribute,
237  bool reIndex
238  );
239 
240  //- Create a point-pair at a ppDist distance either side of
241  // surface point surfPt, in the direction n
242  inline void createPointPair
243  (
244  const scalar ppDist,
245  const Foam::point& surfPt,
246  const vector& n,
247  const bool ptPair,
248  DynamicList<Vb>& pts
249  ) const;
250 
251  inline Foam::point perturbPoint(const Foam::point& pt) const;
252 
253  //- Create a point-pair at a ppDist distance either side of
254  // surface point surfPt, in the direction n
255  inline void createBafflePointPair
256  (
257  const scalar ppDist,
258  const Foam::point& surfPt,
259  const vector& n,
260  const bool ptPair,
261  DynamicList<Vb>& pts
262  ) const;
263 
264  //- Check internal point is completely inside the meshable region
265  inline bool internalPointIsInside(const Foam::point& pt) const;
266 
267  //- Insert pairs of points on the surface with the given normals, at the
268  // specified spacing
269  void insertSurfacePointPairs
270  (
271  const pointIndexHitAndFeatureList& surfaceHits,
272  const fileName fName,
273  DynamicList<Vb>& pts
274  );
275 
276  //- Insert groups of points to conform to an edge given a list of
277  // pointIndexHits specifying the location and edge index of the point
278  // to be conformed to on the corresponding entry in featureHit
279  void insertEdgePointGroups
280  (
281  const pointIndexHitAndFeatureList& edgeHits,
282  const fileName fName,
283  DynamicList<Vb>& pts
284  );
285 
286  void createEdgePointGroupByCirculating
287  (
288  const extendedFeatureEdgeMesh& feMesh,
289  const pointIndexHit& edHit,
290  DynamicList<Vb>& pts
291  ) const;
292 
293  bool meshableRegion
294  (
295  const plane::side side,
297  ) const;
298 
299  bool regionIsInside
300  (
302  const vector& normalA,
304  const vector& normalB,
305  const vector& masterPtVec
306  ) const;
307 
308  //- Create points to conform to an external edge
309  void createExternalEdgePointGroup
310  (
311  const extendedFeatureEdgeMesh& feMesh,
312  const pointIndexHit& edHit,
313  DynamicList<Vb>& pts
314  ) const;
315 
316  //- Create points to conform to an internal edge
317  void createInternalEdgePointGroup
318  (
319  const extendedFeatureEdgeMesh& feMesh,
320  const pointIndexHit& edHit,
321  DynamicList<Vb>& pts
322  ) const;
323 
324  //- Create points to conform to a flat edge
325  void createFlatEdgePointGroup
326  (
327  const extendedFeatureEdgeMesh& feMesh,
328  const pointIndexHit& edHit,
329  DynamicList<Vb>& pts
330  ) const;
331 
332  //- Create points to conform to an open edge
333  void createOpenEdgePointGroup
334  (
335  const extendedFeatureEdgeMesh& feMesh,
336  const pointIndexHit& edHit,
337  DynamicList<Vb>& pts
338  ) const;
339 
340  //- Create points to conform to multiply connected edge
341  void createMultipleEdgePointGroup
342  (
343  const extendedFeatureEdgeMesh& feMesh,
344  const pointIndexHit& edHit,
345  DynamicList<Vb>& pts
346  ) const;
347 
348  //- Determine and insert point groups at the feature points
349  void insertFeaturePoints(bool distribute = false);
350 
351  //- Check if a location is in exclusion range around a feature point
352  bool nearFeaturePt(const Foam::point& pt) const;
353 
354  //- Check if a surface point is in exclusion range around a feature edge
355  bool surfacePtNearFeatureEdge(const Foam::point& pt) const;
356 
357  //- Insert the initial points into the triangulation, based on the
358  // initialPointsMethod
359  void insertInitialPoints();
360 
361  //- In parallel redistribute the backgroundMeshDecomposition and
362  // vertices to balance the number of vertices on each processor.
363  // Returns true if the background mesh changes as this removes all
364  // referred vertices, so the parallel interface may need rebuilt.
365  template<class Triangulation>
366  bool distributeBackground(const Triangulation& mesh);
367 
368  // Test for full containment
369  void cellSizeMeshOverlapsBackground() const;
370 
371  //-
372  void distribute();
373 
374  void buildCellSizeAndAlignmentMesh();
375 
376  //- Set the size and alignment data for each vertex
377  void setVertexSizeAndAlignment();
378 
379  //- Builds a dual face by circulating around the supplied edge.
380  face buildDualFace
381  (
382  const Delaunay::Finite_edges_iterator& eit
383  ) const;
384 
385  boolList dualFaceBoundaryPoints
386  (
387  const Delaunay::Finite_edges_iterator& eit
388  ) const;
389 
390  //- Finds the maximum filterCount of the dual vertices
391  // (Delaunay cells) that form the dual face produced by the
392  // supplied edge
393  label maxFilterCount
394  (
395  const Delaunay::Finite_edges_iterator& eit
396  ) const;
397 
398  //- Determines the owner and neighbour labels for dual cells
399  // corresponding to the dual face formed by the supplied
400  // Delaunay vertices. If the dual face is a boundary face
401  // then neighbour = -1. Returns true if the dual face
402  // created by vA -> vB needs to be reversed to be correctly
403  // orientated.
404  bool ownerAndNeighbour
405  (
406  Vertex_handle vA,
407  Vertex_handle vB,
408  label& owner,
409  label& neighbour
410  ) const;
411 
412  //- Insert the necessary point pairs to conform to the surface, either
413  // from stored results, or trigger a re-conformation
414  void conformToSurface();
415 
416  //- Decision making function for when to rebuild the surface
417  // conformation
418  bool reconformToSurface() const;
419 
420  //- Determines geometrically whether a vertex is close to a surface
421  // This is an optimisation
422  label findVerticesNearBoundaries();
423 
424  //- Create and insert the necessary point pairs to conform to the
425  // surface, then store the result
426  void buildSurfaceConformation();
427 
428  label synchroniseEdgeTrees
429  (
430  const DynamicList<label>& edgeToTreeShape,
431  pointIndexHitAndFeatureList& featureEdgeHits
432  );
433 
434  label synchroniseSurfaceTrees
435  (
436  const DynamicList<label>& surfaceToTreeShape,
437  pointIndexHitAndFeatureList& surfaceHits
438  );
439 
440  bool surfaceLocationConformsToInside
441  (
442  const pointIndexHitAndFeature& info
443  ) const;
444 
445  //- Check to see if dual cell specified by given vertex iterator
446  // intersects the boundary and hence reqires a point-pair
447  bool dualCellSurfaceAnyIntersection
448  (
449  const Delaunay::Finite_vertices_iterator& vit
450  ) const;
451 
452  //- Return all intersections
453  bool dualCellSurfaceAllIntersections
454  (
455  const Delaunay::Finite_vertices_iterator& vit,
456  pointIndexHitAndFeatureDynList& info
457  ) const;
458 
459  //- Return false if the line is entirely outside the current processor
460  // domain, true is either point is inside, or the processor domain
461  // bounadry is intersected (i.e. the points are box outside but the
462  // line cuts. The points will be moved onto the box where they
463  // intersect.
464  bool clipLineToProc
465  (
466  const Foam::point& pt,
467  Foam::point& a,
468  Foam::point& b
469  ) const;
470 
471  //- Find the "worst" protrusion of a dual cell through the surface,
472  // subject to the maxSurfaceProtrusion tolerance
473  void dualCellLargestSurfaceProtrusion
474  (
475  const Delaunay::Finite_vertices_iterator& vit,
476  pointIndexHit& surfHit,
477  label& hitSurface
478  ) const;
479 
480  void dualCellLargestSurfaceIncursion
481  (
482  const Delaunay::Finite_vertices_iterator& vit,
483  pointIndexHit& surfHit,
484  label& hitSurface
485  ) const;
486 
487  //- Write out vertex-processor occupancy information for debugging
488  void reportProcessorOccupancy();
489 
490  //- Write out debugging information about the surface conformation
491  // quality
492 // void reportSurfaceConformationQuality();
493 
494  //- Limit the displacement of a point so that it doesn't penetrate the
495  // surface to be meshed or come too close to it
496  void limitDisplacement
497  (
498  const Delaunay::Finite_vertices_iterator& vit,
499  vector& displacement,
500  label callCount = 0
501  ) const;
502 
503  //- Find angle between the normals of two close surface points.
504  scalar angleBetweenSurfacePoints(Foam::point pA, Foam::point pB) const;
505 
506  //- Check if a surface point is near another.
507  bool nearSurfacePoint
508  (
509  pointIndexHitAndFeature& pHit
510  ) const;
511 
512  //- Append a point to the surface point tree and the existing list
513  bool appendToSurfacePtTree
514  (
515  const Foam::point& pt
516  ) const;
517 
518  //- Append a point to the edge location tree and the existing list
519  bool appendToEdgeLocationTree
520  (
521  const Foam::point& pt
522  ) const;
523 
524  //- Return a list of the nearest feature edge locations
525  List<pointIndexHit> nearestFeatureEdgeLocations
526  (
527  const Foam::point& pt
528  ) const;
529 
530  //- Check if a point is near any feature edge points.
531  bool pointIsNearFeatureEdgeLocation(const Foam::point& pt) const;
532 
533  bool pointIsNearFeatureEdgeLocation
534  (
535  const Foam::point& pt,
536  pointIndexHit& info
537  ) const;
538 
539  //- Check if a point is near any surface conformation points.
540  bool pointIsNearSurfaceLocation(const Foam::point& pt) const;
541 
542  bool pointIsNearSurfaceLocation
543  (
544  const Foam::point& pt,
545  pointIndexHit& info
546  ) const;
547 
548  //- Check if a location is in the exclusion range of an existing feature
549  //- Edge conformation location
550  bool nearFeatureEdgeLocation
551  (
552  const pointIndexHit& pHit,
553  pointIndexHit& nearestEdgeHit
554  ) const;
555 
556  //- Build or rebuild the edge location tree
557  void buildEdgeLocationTree
558  (
559  const DynamicList<Foam::point>& existingEdgeLocations
560  ) const;
561 
562  //- Build or rebuild the surface point location tree
563  void buildSurfacePtLocationTree
564  (
565  const DynamicList<Foam::point>& existingSurfacePtLocations
566  ) const;
567 
568  //- Process the surface conformation locations to decide which surface
569  // and edge conformation locations to add
570  void addSurfaceAndEdgeHits
571  (
572  const Foam::point& vit,
573  const pointIndexHitAndFeatureDynList& surfaceIntersections,
574  scalar surfacePtReplaceDistCoeffSqr,
575  scalar edgeSearchDistCoeffSqr,
576  pointIndexHitAndFeatureDynList& surfaceHits,
577  pointIndexHitAndFeatureDynList& featureEdgeHits,
578  DynamicList<label>& surfaceToTreeShape,
579  DynamicList<label>& edgeToTreeShape,
580  Map<scalar>& surfacePtToEdgePtDist,
581  bool firstPass
582  ) const;
583 
584  //- Store the surface conformation with the indices offset to be
585  // relative to zero
586  void storeSurfaceConformation();
587 
588  //- Reinsert the surface conformation re-offsetting indices to be
589  // relative to new number of internal vertices
590  void reinsertSurfaceConformation();
591 
592  void checkCells();
593 
594  void checkDuals();
595 
596  void checkVertices();
597 
598  void checkCoPlanarCells() const;
599 
600  //- Dual calculation
601  void calcDualMesh
602  (
604  labelList& boundaryPts,
605  faceList& faces,
606  labelList& owner,
607  labelList& neighbour,
610  pointField& cellCentres,
611  labelList& cellToDelaunayVertex,
612  labelListList& patchToDelaunayVertex,
613  PackedBoolList& boundaryFacesToRemove
614  );
615 
616  void calcNeighbourCellCentres
617  (
618  const polyMesh& mesh,
619  const pointField& cellCentres,
620  pointField& neiCc
621  ) const;
622 
623  void selectSeparatedCoupledFaces
624  (
625  const polyMesh& mesh,
626  boolList& selected
627  ) const;
628 
629  //- From meshRefinementBaffles.C. Use insidePoint for a surface to
630  // determine the cell zone.
631  void findCellZoneInsideWalk
632  (
633  const polyMesh& mesh,
634  const labelList& locationSurfaces,
635  const labelList& faceToSurface,
636  labelList& cellToSurface
637  ) const;
638 
639  //- Calculate the cell zones from cellCentres using all closed surfaces
640  labelList calcCellZones(const pointField& cellCentres) const;
641 
642  //- Calculate the face zones
643  void calcFaceZones
644  (
645  const polyMesh& mesh,
646  const pointField& cellCentres,
647  const labelList& cellToSurface,
648  labelList& faceToSurface,
649  boolList& flipMap
650  ) const;
651 
652  //- Add zones to the polyMesh
653  void addZones(polyMesh& mesh, const pointField& cellCentres) const;
654 
655  //- Tet mesh calculation
656  void calcTetMesh
657  (
659  labelList& pointToDelaunayVertex,
660  faceList& faces,
661  labelList& owner,
662  labelList& neighbour,
665  );
666 
667  //- Determines if the dual face constructed by the Delaunay
668  // edge is a boundary face
669  inline bool isBoundaryDualFace
670  (
671  const Delaunay::Finite_edges_iterator& eit
672  ) const;
673 
674  //- Which processors are attached to the dual edge represented by this
675  // Delaunay facet
676  inline List<label> processorsAttached
677  (
678  const Delaunay::Finite_facets_iterator& fit
679  ) const;
680 
681  //- Determines if the edge constructed from the face is on
682  // a processor patch
683  inline bool isParallelDualEdge
684  (
685  const Delaunay::Finite_facets_iterator& fit
686  ) const;
687 
688  //- Determines if the dual face constructed by the Delaunay
689  // edge is a processor boundary face
690  inline bool isProcBoundaryEdge
691  (
692  const Delaunay::Finite_edges_iterator& eit
693  ) const;
694 
695  //- Merge vertices that are identical
696  void mergeIdenticalDualVertices
697  (
698  const pointField& pts,
699  labelList& boundaryPts
700  );
701 
702  label mergeIdenticalDualVertices
703  (
704  const pointField& pts,
705  Map<label>& dualPtIndexMap
706  ) const;
707 
708  //- Identify the face labels of the deferred collapse faces
709  void deferredCollapseFaceSet
710  (
711  labelList& owner,
712  labelList& neighbour,
713  const HashSet<labelPair, labelPair::Hash<>>& deferredCollapseFaces
714  ) const;
715 
716  //- Check whether the cell sizes are fine enough. Creates a polyMesh.
717  void checkCellSizing();
718 
719  //- Find all cells with a patch face that is not near the surface. The
720  // allowed offset is the fraction of the target cell size.
721  labelHashSet findOffsetPatchFaces
722  (
723  const polyMesh& mesh,
724  const scalar allowedOffset
725  ) const;
726 
727  //- Create a polyMesh and check its quality, reports which
728  // elements damage the mesh quality, allowing backtracking.
729  labelHashSet checkPolyMeshQuality(const pointField& pts) const;
730 
731  label classifyBoundaryPoint(Cell_handle cit) const;
732 
733  //- Index all of the the Delaunay cells and calculate their
734  //- Dual points
735  void indexDualVertices
736  (
737  pointField& pts,
738  labelList& boundaryPts
739  );
740 
741  //- Re-index all of the Delaunay cells
742  void reindexDualVertices
743  (
744  const Map<label>& dualPtIndexMap,
745  labelList& boundaryPts
746  );
747 
748  label createPatchInfo
749  (
752  ) const;
753 
754  vector calcSharedPatchNormal(Cell_handle c1, Cell_handle c2) const;
755 
756  bool boundaryDualFace(Cell_handle c1, Cell_handle c2) const;
757 
758  //- Create all of the internal and boundary faces
759  void createFacesOwnerNeighbourAndPatches
760  (
761  const pointField& pts,
762  faceList& faces,
763  labelList& owner,
764  labelList& neighbour,
767  labelListList& patchPointPairSlaves,
768  PackedBoolList& boundaryFacesToRemove,
769  bool includeEmptyPatches = false
770  ) const;
771 
772  //- Sort the faces, owner and neighbour lists into
773  // upper-triangular order. For internal faces only, use
774  // before adding patch faces
775  void sortFaces
776  (
777  faceList& faces,
778  labelList& owner,
779  labelList& neighbour
780  ) const;
781 
782  //- Sort the processor patches so that the faces are in the same order
783  // on both processors
784  void sortProcPatches
785  (
786  List<DynamicList<face>>& patchFaces,
787  List<DynamicList<label>>& patchOwners,
788  List<DynamicList<label>>& patchPointPairSlaves,
789  labelPairPairDynListList& patchSortingIndices
790  ) const;
791 
792  //- Add the faces and owner information for the patches
793  void addPatches
794  (
795  const label nInternalFaces,
796  faceList& faces,
797  labelList& owner,
799  PackedBoolList& boundaryFacesToRemove,
800  const List<DynamicList<face>>& patchFaces,
801  const List<DynamicList<label>>& patchOwners,
802  const List<DynamicList<bool>>& indirectPatchFace
803  ) const;
804 
805  //- Remove points that are no longer used by any faces
806  void removeUnusedPoints
807  (
808  faceList& faces,
809  pointField& pts,
810  labelList& boundaryPts
811  ) const;
812 
813  //- Remove dual cells that are not used by any faces. Return compaction
814  // map.
815  labelList removeUnusedCells
816  (
817  labelList& owner,
818  labelList& neighbour
819  ) const;
820 
821  //- Create an empty fvMesh
822  autoPtr<fvMesh> createDummyMesh
823  (
824  const IOobject& io,
825  const wordList& patchNames,
827  ) const;
828 
829  //- Create a polyMesh from points.
830  autoPtr<polyMesh> createPolyMeshFromPoints(const pointField& pts) const;
831 
832  void checkProcessorPatchesMatch
833  (
835  ) const;
836 
837  void reorderPoints
838  (
840  labelList& boundaryPts,
841  faceList& faces,
842  const label nInternalFaces
843  ) const;
844 
845  //- Rotate the faces on processor patches if necessary
846  void reorderProcessorPatches
847  (
848  const word& meshName,
849  const fileName& instance,
850  const pointField& points,
851  faceList& faces,
852  const wordList& patchNames,
854  ) const;
855 
856  void writePointPairs(const fileName& fName) const;
857 
858  //- Disallow default bitwise copy construct
860 
861  //- Disallow default bitwise assignment
862  void operator=(const conformalVoronoiMesh&);
863 
864 
865 public:
866 
867  //- Runtime type information
868  ClassName("conformalVoronoiMesh");
869 
870 
871  // Constructors
872 
873  //- Construct from Time and foamyHexMeshDict
875  (
876  const Time& runTime,
877  const dictionary& foamyHexMeshDict
878  );
879 
880 
881  //- Destructor
883 
884 
885  // Member Functions
886 
887  void initialiseForMotion();
888 
890 
891  //- Move the vertices according to the controller, re-conforming to the
892  // surface as required
893  void move();
894 
895 // //- Which other processors does each sphere overlap
896 // labelListList overlapsProc
897 // (
898 // const List<Foam::point>& centres,
899 // const List<scalar>& radiusSqrs
900 // ) const;
901 
902 
903  // Access
904 
905  //- Return the Time object
906  inline const Time& time() const;
907 
908  //- Return the random number generator
909  inline Random& rndGen() const;
910 
911  //- Return the allGeometry object
912  inline const searchableSurfaces& allGeometry() const;
913 
914  //- Return the conformationSurfaces object
915  inline const conformationSurfaces& geometryToConformTo() const;
916 
917  //- Return the backgroundMeshDecomposition
918  inline const backgroundMeshDecomposition& decomposition() const;
919 
920  //- Return the cellShapeControl object
921  inline const cellShapeControl& cellShapeControls() const;
922 
923  //- Return the foamyHexMeshControls object
924  inline const cvControls& foamyHexMeshControls() const;
925 
926 
927  // Query
928 
929  //- Return the local point pair separation at the given location
930  inline scalar pointPairDistance(const Foam::point& pt) const;
931 
932  //- Return the local mixed feature point placement distance
933  inline scalar mixedFeaturePointDistance
934  (
935  const Foam::point& pt
936  ) const;
937 
938  //- Return the square of the local feature point exclusion distance
940  (
941  const Foam::point& pt
942  ) const;
943 
944  //- Return the square of the local feature edge exclusion distance
945  inline scalar featureEdgeExclusionDistanceSqr
946  (
947  const Foam::point& pt
948  ) const;
949 
950  //- Return the square of the local surface point exclusion distance
951  inline scalar surfacePtExclusionDistanceSqr
952  (
953  const Foam::point& pt
954  ) const;
955 
956  //- Return the square of the local surface search distance
957  inline scalar surfaceSearchDistanceSqr(const Foam::point& pt) const;
958 
959  //- Return the local maximum surface protrusion distance
960  inline scalar maxSurfaceProtrusion(const Foam::point& pt) const;
961 
962  //- Call the appropriate function to conform to an edge
964  (
965  const extendedFeatureEdgeMesh& feMesh,
966  const pointIndexHit& edHit,
967  DynamicList<Vb>& pts
968  ) const;
969 
970 
971  // Write
972 
973  //- Write the elapsedCpuTime and memory usage, with an optional
974  // description
975  static void timeCheck
976  (
977  const Time& runTime,
978  const string& description = string::null,
979  const bool check = true
980  );
981 
982  void timeCheck
983  (
984  const string& description = string::null
985  ) const;
986 
987  //- Prepare data and call writeMesh for polyMesh and
988  // tetDualMesh
989  void writeMesh(const fileName& instance);
990 
991  //- Write mesh to disk
992  void writeMesh
993  (
994  const word& meshName,
995  const fileName& instance,
997  labelList& boundaryPts,
998  faceList& faces,
999  labelList& owner,
1000  labelList& neighbour,
1001  const wordList& patchNames,
1003  const pointField& cellCentres,
1004  PackedBoolList& boundaryFacesToRemove
1005  ) const;
1006 
1007  //- Calculate and write a field of the target cell size,
1008  // target cell volume, actual cell volume and equivalent
1009  // actual cell size (cbrt(actual cell volume)).
1010  void writeCellSizes(const fvMesh& mesh) const;
1011 
1012  void writeCellAlignments(const fvMesh& mesh) const;
1013 
1014  //- Calculate and write the cell centres.
1015  void writeCellCentres(const fvMesh& mesh) const;
1016 
1017  //- Find the cellSet of the boundary cells which have points that
1018  // protrude out of the surface beyond a tolerance.
1020 };
1021 
1022 
1023 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
1024 
1025 } // End namespace Foam
1026 
1027 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
1028 
1029 #include "conformalVoronoiMeshI.H"
1030 
1031 #ifdef NoRepository
1033 #endif
1034 
1035 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
1036 
1037 #endif
1038 
1039 // ************************************************************************* //
Delaunay::Cell_handle Cell_handle
A HashTable with keys but without contents.
Definition: HashSet.H:59
~conformalVoronoiMesh()
Destructor.
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:69
Delaunay::Vertex_handle Vertex_handle
scalar featurePointExclusionDistanceSqr(const Foam::point &pt) const
Return the square of the local feature point exclusion distance.
const cvControls & foamyHexMeshControls() const
Return the foamyHexMeshControls object.
A face is a list of labels corresponding to mesh vertices.
Definition: face.H:75
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
void writeCellSizes(const fvMesh &mesh) const
Calculate and write a field of the target cell size,.
A 2-tuple for storing two objects of different types.
Definition: Tuple2.H:47
labelHashSet findRemainingProtrusionSet(const polyMesh &mesh) const
Find the cellSet of the boundary cells which have points that.
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: HashTable.H:59
scalar pointPairDistance(const Foam::point &pt) const
Return the local point pair separation at the given location.
Unit conversion functions.
const dimensionedScalar c2
Second radiation constant: default SI units: [m.K].
static void timeCheck(const Time &runTime, const string &description=string::null, const bool check=true)
Write the elapsedCpuTime and memory usage, with an optional.
void writeCellCentres(const fvMesh &mesh) const
Calculate and write the cell centres.
ClassName("conformalVoronoiMesh")
Runtime type information.
const cellShapeControl & cellShapeControls() const
Return the cellShapeControl object.
This class describes the interaction of (usually) a face and a point. It carries the info of a succes...
Definition: PointIndexHit.H:53
List< pointIndexHitAndFeature > pointIndexHitAndFeatureList
void writeCellAlignments(const fvMesh &mesh) const
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:68
void writeMesh(const fileName &instance)
Prepare data and call writeMesh for polyMesh and.
void move()
Move the vertices according to the controller, re-conforming to the.
An ordered pair of two objects of type <T> with first() and second() elements.
Definition: contiguous.H:49
dynamicFvMesh & mesh
static char meshName[]
Definition: globalFoam.H:7
Tuple2< pointIndexHit, label > pointIndexHitAndFeature
const pointField & points
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects...
Definition: DynamicList.H:56
System uinteger.
const dimensionedScalar b
Wien displacement law constant: default SI units: [m.K].
Definition: createFields.H:27
A class for handling words, derived from string.
Definition: word.H:59
3D tensor transformation operations.
scalar maxSurfaceProtrusion(const Foam::point &pt) const
Return the local maximum surface protrusion distance.
const backgroundMeshDecomposition & decomposition() const
Return the backgroundMeshDecomposition.
wordList patchNames(nPatches)
Container for searchableSurfaces.
sideVolumeType
Normals point to the outside.
Simple random number generator.
Definition: Random.H:49
const searchableSurfaces & allGeometry() const
Return the allGeometry object.
Controls for the conformalVoronoiMesh mesh generator.
Definition: cvControls.H:53
static const string null
An empty string.
Definition: string.H:86
DynamicList< pointIndexHitAndFeature > pointIndexHitAndFeatureDynList
Random & rndGen() const
Return the random number generator.
scalar surfacePtExclusionDistanceSqr(const Foam::point &pt) const
Return the square of the local surface point exclusion distance.
The Delaunay vertices required to conform to a feature point can be determined upon initialisation be...
A bit-packed bool list.
side
Side of the plane.
Definition: plane.H:65
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:62
scalar surfaceSearchDistanceSqr(const Foam::point &pt) const
Return the square of the local surface search distance.
const Time & time() const
Return the Time object.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
const dimensionedScalar c1
First radiation constant: default SI units: [W/m2].
scalar mixedFeaturePointDistance(const Foam::point &pt) const
Return the local mixed feature point placement distance.
Standard boundBox + extra functionality for use in octree.
Definition: treeBoundBox.H:87
static const NamedEnum< dualMeshPointType, 5 > dualMeshPointTypeNames_
label n
Store a background polyMesh to use for the decomposition of space and queries for parallel conformalV...
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:53
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
CGAL data structures used for 3D Delaunay meshing.
PtrList< dictionary > patchDicts
Definition: readKivaGrid.H:537
void createEdgePointGroup(const extendedFeatureEdgeMesh &feMesh, const pointIndexHit &edHit, DynamicList< Vb > &pts) const
Call the appropriate function to conform to an edge.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
List< DynamicList< Pair< labelPair > > > labelPairPairDynListList
scalar featureEdgeExclusionDistanceSqr(const Foam::point &pt) const
Return the square of the local feature edge exclusion distance.
const conformationSurfaces & geometryToConformTo() const
Return the conformationSurfaces object.
Namespace for OpenFOAM.