meshRefinement.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-2022 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::meshRefinement
26 
27 Description
28  Helper class which maintains intersections of (changing) mesh with
29  (static) surfaces.
30 
31  Maintains
32  - per face any intersections of the cc-cc segment with any of the surfaces
33 
34 SourceFiles
35  meshRefinement.C
36  meshRefinementBaffles.C
37  meshRefinementMerge.C
38  meshRefinementProblemCells.C
39  meshRefinementRefine.C
40 
41 \*---------------------------------------------------------------------------*/
42 
43 #ifndef meshRefinement_H
44 #define meshRefinement_H
45 
46 #include "hexRef8.H"
47 #include "polyTopoChangeMap.H"
48 #include "autoPtr.H"
49 #include "labelPair.H"
50 #include "indirectPrimitivePatch.H"
51 #include "pointFieldsFwd.H"
52 #include "Tuple2.H"
53 #include "pointIndexHit.H"
54 #include "refinementParameters.H"
55 
56 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
57 
58 namespace Foam
59 {
60 
61 // Class forward declarations
62 class fvMesh;
63 class polyDistributionMap;
64 class decompositionMethod;
65 class refinementSurfaces;
66 class refinementFeatures;
67 class refinementRegions;
68 class removeCells;
69 class fvMeshDistribute;
70 class searchableSurface;
71 class regionSplit;
72 class globalIndex;
73 class removePoints;
74 class localPointRegion;
75 
76 class snapParameters;
77 
78 /*---------------------------------------------------------------------------*\
79  Class meshRefinement Declaration
80 \*---------------------------------------------------------------------------*/
81 
82 class meshRefinement
83 {
84 
85 public:
86 
87  // Public data types
88 
89  //- Enumeration for what to debug
90  enum IOdebugType
91  {
97  };
98 
100  enum debugType
101  {
102  MESH = 1<<IOMESH,
106  LAYERINFO = 1<<IOLAYERINFO
107  };
108 
109  //- Enumeration for what to output
110  enum IOoutputType
111  {
113  };
116  enum outputType
117  {
119  };
120 
121  //- Enumeration for what to write
122  enum IOwriteType
123  {
129  };
132  enum writeType
133  {
138  WRITELAYERFIELDS = 1<<IOWRITELAYERFIELDS
139  };
140 
141  //- Enumeration for how the userdata is to be mapped upon refinement.
142  enum mapType
143  {
145  KEEPALL = 2,
146  REMOVE = 4
147  };
148 
149 
150 private:
151 
152  // Static Data Members
153 
154  //- Control of writing level
155  static writeType writeLevel_;
156 
157  //- Control of output/log level
158  static outputType outputLevel_;
159 
160 
161  // Private Data
162 
163  //- Reference to mesh
164  fvMesh& mesh_;
165 
166  //- Tolerance used for sorting coordinates (used in 'less' routine)
167  const scalar mergeDistance_;
168 
169  //- Overwrite the mesh?
170  const bool overwrite_;
171 
172  //- Instance of mesh upon construction. Used when in overwrite_ mode.
173  const word oldInstance_;
174 
175  //- All surface-intersection interaction
176  const refinementSurfaces& surfaces_;
177 
178  //- All feature-edge interaction
179  const refinementFeatures& features_;
180 
181  //- All shell-refinement interaction
182  const refinementRegions& shells_;
183 
184  //- Refinement engine
185  hexRef8 meshCutter_;
186 
187  //- Per cc-cc vector the index of the surface hit
188  labelIOList surfaceIndex_;
189 
190  //- User supplied face based data.
191  List<Tuple2<mapType, labelList>> userFaceData_;
192 
193  //- Meshed patches - are treated differently. Stored as wordList since
194  // order changes.
195  wordList meshedPatches_;
196 
197 
198  // Private Member Functions
199 
200  //- Find out which faces have changed given cells (old mesh labels)
201  // that were marked for refinement.
202  static labelList getChangedFaces
203  (
204  const polyTopoChangeMap&,
205  const labelList& oldCellsToRefine
206  );
207 
208  //- Calculate coupled boundary end vector and refinement level
209  void calcNeighbourData
210  (
211  labelList& neiLevel,
212  pointField& neiCc
213  ) const;
214 
215  //- Find any intersection of surface. Store in surfaceIndex_.
216  void updateIntersections(const labelList& changedFaces);
217 
218  //- Remove cells. Put exposedFaces into exposedPatchIDs.
219  autoPtr<polyTopoChangeMap> doRemoveCells
220  (
221  const labelList& cellsToRemove,
222  const labelList& exposedFaces,
223  const labelList& exposedPatchIDs,
224  removeCells& cellRemover
225  );
226 
227  // Get cells which are inside any closed surface. Note that
228  // all closed surfaces
229  // will have already been oriented to have insidePoint outside.
230  labelList getInsideCells(const word&) const;
231 
232  // Do all to remove inside cells
233  autoPtr<polyTopoChangeMap> removeInsideCells
234  (
235  const string& msg,
236  const label exposedPatchi
237  );
238 
239 
240  // Refinement candidate selection
241 
242  //- Mark cell for refinement (if not already marked). Return false
243  // if refinelimit hit. Keeps running count (in nRefine) of cells
244  // marked for refinement
245  static bool markForRefine
246  (
247  const label markValue,
248  const label nAllowRefine,
249  label& cellValue,
250  label& nRefine
251  );
252 
253  //- Mark every cell with level of feature passing through it
254  // (or -1 if not passed through). Uses tracking.
255  void markFeatureCellLevel
256  (
257  const List<point>& insidePoints,
258  labelList& maxFeatureLevel
259  ) const;
260 
261  //- Calculate list of cells to refine based on intersection of
262  // features.
263  label markFeatureRefinement
264  (
265  const List<point>& insidePoints,
266  const label nAllowRefine,
267 
269  label& nRefine
270  ) const;
271 
272  //- Mark cells for distance-to-feature based refinement.
273  label markInternalDistanceToFeatureRefinement
274  (
275  const label nAllowRefine,
277  label& nRefine
278  ) const;
279 
280  //- Mark cells for refinement-shells based refinement.
281  label markInternalRefinement
282  (
283  const label nAllowRefine,
285  label& nRefine
286  ) const;
287 
288  //- Collect faces that are intersected and whose neighbours aren't
289  // yet marked for refinement.
290  labelList getRefineCandidateFaces
291  (
292  const labelList& refineCell
293  ) const;
294 
295  //- Mark cells for surface intersection based refinement.
296  label markSurfaceRefinement
297  (
298  const label nAllowRefine,
299  const labelList& neiLevel,
300  const pointField& neiCc,
302  label& nRefine
303  ) const;
304 
305  //- Helper: count number of normals1 that are in normals2
306  label countMatches
307  (
308  const List<point>& normals1,
309  const List<point>& normals2,
310  const scalar tol = 1e-6
311  ) const;
312 
313  //- Mark cells for surface curvature based refinement. Marks if
314  // local curvature > curvature or if on different regions
315  // (markDifferingRegions)
316  label markSurfaceCurvatureRefinement
317  (
318  const scalar curvature,
319  const label nAllowRefine,
320  const labelList& neiLevel,
321  const pointField& neiCc,
323  label& nRefine
324  ) const;
325 
326  //- Mark cell if local a gap topology or
327  bool checkProximity
328  (
329  const scalar planarCos,
330  const label nAllowRefine,
331 
332  const label surfaceLevel,
333  const vector& surfaceLocation,
334  const vector& surfaceNormal,
335 
336  const label celli,
337 
338  label& cellMaxLevel,
339  vector& cellMaxLocation,
340  vector& cellMaxNormal,
341 
343  label& nRefine
344  ) const;
345 
346  //- Mark cells for surface proximity based refinement.
347  label markProximityRefinement
348  (
349  const scalar curvature,
350  const label nAllowRefine,
351  const labelList& neiLevel,
352  const pointField& neiCc,
353 
355  label& nRefine
356  ) const;
357 
358  // Baffle handling
359 
360  //- Get faces to repatch. Returns map from face to patch.
361  Map<labelPair> getZoneBafflePatches
362  (
363  const bool allowBoundary,
364  const labelList& globalToMasterPatch,
365  const labelList& globalToSlavePatch
366  ) const;
367 
368  //- Determine patches for baffles on all intersected unnamed faces
369  void getBafflePatches
370  (
371  const labelList& globalToMasterPatch,
372  const labelList& neiLevel,
373  const pointField& neiCc,
374  labelList& ownPatch,
375  labelList& nbrPatch
376  ) const;
377 
378  //- Repatches external face or creates baffle for internal face
379  // with user specified patches (might be different for both sides).
380  // Returns label of added face.
381  label createBaffle
382  (
383  const label facei,
384  const label ownPatch,
385  const label nbrPatch,
386  polyTopoChange& meshMod
387  ) const;
388 
389  // Problem cell handling
390 
391  //- Helper function to mark face as being on 'boundary'. Used by
392  // markFacesOnProblemCells
393  void markBoundaryFace
394  (
395  const label facei,
396  boolList& isBoundaryFace,
397  boolList& isBoundaryEdge,
398  boolList& isBoundaryPoint
399  ) const;
400 
401  void findNearest
402  (
403  const labelList& meshFaces,
404  List<pointIndexHit>& nearestInfo,
405  labelList& nearestSurface,
406  labelList& nearestRegion,
407  vectorField& nearestNormal
408  ) const;
409 
410  Map<label> findEdgeConnectedProblemCells
411  (
412  const scalarField& perpendicularAngle,
413  const labelList&
414  ) const;
415 
416  bool isCollapsedFace
417  (
418  const pointField&,
419  const pointField& neiCc,
420  const scalar minFaceArea,
421  const scalar maxNonOrtho,
422  const label facei
423  ) const;
424 
425  bool isCollapsedCell
426  (
427  const pointField&,
428  const scalar volFraction,
429  const label celli
430  ) const;
431 
432  //- Returns list with for every internal face -1 or the patch
433  // they should be baffled into. If removeEdgeConnectedCells is set
434  // removes cells based on perpendicularAngle.
435  //
436  // Gets run after createBaffles so all the unzoned surface
437  // intersections have already been turned into baffles. (Note:
438  // zoned surfaces are not baffled at this stage) Used to remove
439  // cells by baffling all their faces and have the splitMeshRegions
440  // chuck away non used regions.
441  labelList markFacesOnProblemCells
442  (
443  const dictionary& motionDict,
444  const bool removeEdgeConnectedCells,
445  const scalarField& perpendicularAngle,
446  const labelList& globalToMasterPatch
447  ) const;
448 
449  //- Returns list with for every face the label of the nearest
450  // patch. Any unreached face (disconnected mesh?) becomes
451  // adaptPatchIDs[0]
452  labelList nearestPatch(const labelList& adaptPatchIDs) const;
453 
454  //- Returns list with for every internal face -1 or the patch
455  // they should be baffled into.
456  labelList markFacesOnProblemCellsGeometric
457  (
458  const snapParameters& snapParams,
459  const dictionary& motionDict
460  ) const;
461 
462 
463  // Baffle merging
464 
465  //- Extract those baffles (duplicate) faces that are on the edge
466  // of a baffle region. These are candidates for merging.
467  List<labelPair> freeStandingBaffles
468  (
469  const List<labelPair>&,
470  const scalar freeStandingAngle
471  ) const;
472 
473 
474  // Zone handling
475 
476  //- Finds zone per cell for cells inside closed named surfaces.
477  // (uses geometric test for insideness)
478  // Adapts namedSurfaceIndex so all faces on boundary of cellZone
479  // have corresponding faceZone.
480  void findCellZoneGeometric
481  (
482  const pointField& neiCc,
483  const labelList& closedNamedSurfaces,
484  labelList& namedSurfaceIndex,
485  const labelList& surfaceToCellZone,
486  labelList& cellToZone
487  ) const;
488 
489  //- Finds zone per cell for cells inside named surfaces that have
490  // an inside point specified.
491  void findCellZoneInsideWalk
492  (
493  const labelList& locationSurfaces,
494  const labelList& namedSurfaceIndex,
495  const labelList& surfaceToCellZone,
496  labelList& cellToZone
497  ) const;
498 
499  //- Determines cell zone from cell region information.
500  bool calcRegionToZone
501  (
502  const label surfZoneI,
503  const label ownRegion,
504  const label neiRegion,
505 
506  labelList& regionToCellZone
507  ) const;
508 
509  //- Finds zone per cell. Uses topological walk with all faces
510  // marked in namedSurfaceIndex regarded as blocked.
511  void findCellZoneTopo
512  (
513  const List<point>& insidePoints,
514  const labelList& namedSurfaceIndex,
515  const labelList& surfaceToCellZone,
516  labelList& cellToZone
517  ) const;
518 
519  //- Make namedSurfaceIndex consistent with cellToZone
520  // - clear out any blocked faces in between same cell zone.
521  void makeConsistentFaceIndex
522  (
523  const labelList& cellToZone,
524  labelList& namedSurfaceIndex
525  ) const;
526 
527  //- Remove any loose standing cells
528  void handleSnapProblems
529  (
530  const snapParameters& snapParams,
531  const bool useTopologicalSnapDetection,
532  const bool removeEdgeConnectedCells,
533  const scalarField& perpendicularAngle,
534  const dictionary& motionDict,
535  Time& runTime,
536  const labelList& globalToMasterPatch,
537  const labelList& globalToSlavePatch
538  );
539 
540 
541  // Some patch utilities
542 
543  //- Get all faces in faceToZone that have no cellZone on
544  // either side.
545  labelList freeStandingBaffleFaces
546  (
547  const labelList& faceToZone,
548  const labelList& cellToZone,
549  const labelList& neiCellZone
550  ) const;
551 
552  //- Determine per patch edge the number of master faces. Used
553  // to detect non-manifold situations.
554  void calcPatchNumMasterFaces
555  (
556  const PackedBoolList& isMasterFace,
557  const indirectPrimitivePatch& patch,
558  labelList& nMasterFaces
559  ) const;
560 
561  //- Determine per patch face the (singly-) connected zone it
562  // is in. Return overall number of zones.
563  label markPatchZones
564  (
565  const indirectPrimitivePatch& patch,
566  const labelList& nMasterFaces,
567  labelList& faceToZone
568  ) const;
569 
570  //- Make faces consistent.
571  void consistentOrientation
572  (
573  const PackedBoolList& isMasterFace,
574  const indirectPrimitivePatch& patch,
575  const labelList& nMasterFaces,
576  const labelList& faceToZone,
577  const Map<label>& zoneToOrientation,
578  boolList& meshFlipMap
579  ) const;
580 
581 
582 public:
583 
584  //- Runtime type information
585  ClassName("meshRefinement");
586 
587 
588  // Constructors
589 
590  //- Construct from components
592  (
593  fvMesh& mesh,
594  const dictionary& refineDict,
595  const scalar mergeDistance,
596  const bool overwrite,
598  const refinementFeatures&,
599  const refinementRegions&
600  );
601 
602  //- Disallow default bitwise copy construction
603  meshRefinement(const meshRefinement&) = delete;
604 
605 
606  // Member Functions
607 
608  // Access
609 
610  //- Reference to mesh
611  const fvMesh& mesh() const
612  {
613  return mesh_;
614  }
616  fvMesh& mesh()
617  {
618  return mesh_;
619  }
621  scalar mergeDistance() const
622  {
623  return mergeDistance_;
624  }
625 
626  //- Overwrite the mesh?
627  bool overwrite() const
628  {
629  return overwrite_;
630  }
631 
632  //- (points)instance of mesh upon construction
633  const word& oldInstance() const
634  {
635  return oldInstance_;
636  }
637 
638  //- Reference to surface search engines
639  const refinementSurfaces& surfaces() const
640  {
641  return surfaces_;
642  }
643 
644  //- Reference to feature edge mesh
645  const refinementFeatures& features() const
646  {
647  return features_;
648  }
649 
650  //- Reference to refinement shells (regions)
651  const refinementRegions& shells() const
652  {
653  return shells_;
654  }
655 
656  //- Reference to meshcutting engine
657  const hexRef8& meshCutter() const
658  {
659  return meshCutter_;
660  }
661 
662  //- Per start-end edge the index of the surface hit
663  const labelList& surfaceIndex() const
664  {
665  return surfaceIndex_;
666  }
669  {
670  return surfaceIndex_;
671  }
672 
673  //- Additional face data that is maintained across
674  // topo changes. Every entry is a list over all faces.
675  // Bit of a hack. Additional flag to say whether to maintain master
676  // only (false) or increase set to account for face-from-face.
678  {
679  return userFaceData_;
680  }
683  {
684  return userFaceData_;
685  }
686 
687 
688  // Other
689 
690  //- Count number of intersections (local)
691  label countHits() const;
692 
693  //- Redecompose according to cell count
694  // keepZoneFaces : find all faceZones from zoned surfaces and keep
695  // owner and neighbour together
696  // keepBaffles : find all baffles and keep them together
698  (
699  const bool keepZoneFaces,
700  const bool keepBaffles,
701  const scalarField& cellWeights,
702  decompositionMethod& decomposer,
703  fvMeshDistribute& distributor
704  );
705 
706  //- Get faces with intersection.
707  labelList intersectedFaces() const;
708 
709  //- Get points on surfaces with intersection and boundary faces.
711 
712  //- Create patch from set of patches
714  (
715  const polyMesh&,
716  const labelList&
717  );
718 
719  //- Helper function to make a pointVectorField with correct
720  // bcs for mesh movement:
721  // - adaptPatchIDs : fixedValue
722  // - processor : calculated (so free to move)
723  // - cyclic/wedge/symmetry : slip
724  // - other : slip
726  (
727  const pointMesh& pMesh,
728  const labelList& adaptPatchIDs
729  );
730 
731  //- Helper function: check that face zones are synced
732  static void checkCoupledFaceZones(const polyMesh&);
733 
734  //- Helper: calculate edge weights (1/length)
735  static void calculateEdgeWeights
736  (
737  const polyMesh& mesh,
738  const PackedBoolList& isMasterEdge,
739  const labelList& meshPoints,
740  const edgeList& edges,
741  scalarField& edgeWeights,
742  scalarField& invSumWeight
743  );
744 
745  //- Helper: weighted sum (over all subset of mesh points) by
746  // summing contribution from (master) edges
747  template<class Type>
748  static void weightedSum
749  (
750  const polyMesh& mesh,
751  const PackedBoolList& isMasterEdge,
752  const labelList& meshPoints,
753  const edgeList& edges,
754  const scalarField& edgeWeights,
755  const Field<Type>& data,
757  );
758 
759 
760  // Refinement
761 
762  //- Is local topology a small gap?
763  bool isGap
764  (
765  const scalar,
766  const vector&,
767  const vector&,
768  const vector&,
769  const vector&
770  ) const;
771 
772  //- Is local topology a small gap normal to the test vector
773  bool isNormalGap
774  (
775  const scalar,
776  const vector&,
777  const vector&,
778  const vector&,
779  const vector&
780  ) const;
781 
782  //- Calculate list of cells to refine.
784  (
785  const List<point>& insidePoints,
786  const scalar curvature,
787  const scalar planarAngle,
788 
789  const bool featureRefinement,
790  const bool featureDistanceRefinement,
791  const bool internalRefinement,
792  const bool surfaceRefinement,
793  const bool curvatureRefinement,
794  const bool gapRefinement,
795  const label maxGlobalCells,
796  const label maxLocalCells
797  ) const;
798 
799  //- Refine some cells
800  autoPtr<polyTopoChangeMap> refine(const labelList& cellsToRefine);
801 
802  //- Refine some cells and rebalance
804  (
805  const string& msg,
806  decompositionMethod& decomposer,
807  fvMeshDistribute& distributor,
808  const labelList& cellsToRefine,
809  const scalar maxLoadUnbalance
810  );
811 
812  //- Balance before refining some cells
814  (
815  const string& msg,
816  decompositionMethod& decomposer,
817  fvMeshDistribute& distributor,
818  const labelList& cellsToRefine,
819  const scalar maxLoadUnbalance
820  );
821 
822 
823  // Baffle handling
824 
825  //- Split off unreachable areas of mesh.
826  void baffleAndSplitMesh
827  (
828  const bool handleSnapProblems,
829 
830  // How to remove problem snaps
831  const snapParameters& snapParams,
832  const bool useTopologicalSnapDetection,
833  const bool removeEdgeConnectedCells,
834  const scalarField& perpendicularAngle,
835 
836  // How to handle free-standing baffles
837  const bool mergeFreeStanding,
838  const scalar freeStandingAngle,
839 
840  const dictionary& motionDict,
841  Time& runTime,
842  const labelList& globalToMasterPatch,
843  const labelList& globalToSlavePatch,
844  const refinementParameters::cellSelectionPoints& selectionPoints
845  );
846 
847  //- Split off (with optional buffer layers) unreachable areas
848  // of mesh. Does not introduce baffles.
850  (
851  const label nBufferLayers,
852  const labelList& globalToMasterPatch,
853  const labelList& globalToSlavePatch,
854  const refinementParameters::cellSelectionPoints& selectionPoints
855  );
856 
857  //- Find boundary points that connect to more than one cell
858  // region and split them.
860  (
861  const localPointRegion&
862  );
863 
864  //- Find boundary points that connect to more than one cell
865  // region and split them.
867 
868  //- Create baffle for every internal face where ownPatch != -1.
869  // External faces get repatched according to ownPatch (nbrPatch
870  // should be -1 for these)
872  (
873  const labelList& ownPatch,
874  const labelList& nbrPatch
875  );
876 
877  //- Debug helper: check faceZones are not on processor patches
878  void checkZoneFaces() const;
879 
880  //- Create baffles for faces straddling zoned surfaces. Return
881  // baffles.
883  (
884  const labelList& globalToMasterPatch,
885  const labelList& globalToSlavePatch,
887  );
888 
889  //- Merge baffles. Gets pairs of faces.
891 
892  //- Put faces/cells into zones according to surface specification.
893  // Returns null if no zone surfaces present. Region containing
894  // the insidePoint will not be put into a cellZone.
896  (
897  const List<point>& insidePoints,
898  const bool allowFreeStandingZoneFaces
899  );
900 
901 
902  // Other topo changes
903 
904  //- Add patch originating from meshing. Update meshedPatches_.
905  label addMeshedPatch(const word& name, const dictionary&);
906 
907  //- Get patchIDs for patches added in addMeshedPatch.
908  labelList meshedPatches() const;
909 
910  //- Select coupled faces that are not collocated
912 
913  //- Find region point is in. Uses optional perturbation to re-test.
914  static label findRegion
915  (
916  const polyMesh&,
917  const labelList& cellRegion,
918  const vector& perturbVec,
919  const point& location
920  );
921 
922  //- Find regions points are in and update cellRegion
923  // Uses optional perturbation to re-test.
924  static void findRegions
925  (
926  const polyMesh&,
927  labelList& cellRegion,
928  const vector& perturbVec,
929  const refinementParameters::cellSelectionPoints& selectionPoints
930  );
931 
932  //- Split mesh according to selectionPoints
934  (
935  const labelList& globalToMasterPatch,
936  const labelList& globalToSlavePatch,
937  const refinementParameters::cellSelectionPoints& selectionPoints
938  );
939 
940  //- Split faces into two
942  (
943  const labelList& splitFaces,
944  const labelPairList& splits
945  );
946 
947  //- Update local numbering for mesh redistribution
948  void distribute(const polyDistributionMap&);
949 
950  //- Update for external change to mesh. changedFaces are in new mesh
951  // face labels.
952  void topoChange
953  (
954  const polyTopoChangeMap&,
955  const labelList& changedFaces
956  );
957 
958  //- Helper: reorder list according to map.
959  template<class T>
960  static void updateList
961  (
962  const labelList& newToOld,
963  const T& nullValue,
964  List<T>& elems
965  );
966 
967 
968  // Restoring : is where other processes delete and reinsert data.
969 
970  //- Signal points/face/cells for which to store data
971  void storeData
972  (
973  const labelList& pointsToStore,
974  const labelList& facesToStore,
975  const labelList& cellsToStore
976  );
977 
978  //- Update local numbering + undo
979  // Data to restore given as new pointlabel + stored pointlabel
980  // (i.e. what was in pointsToStore)
981  void topoChange
982  (
983  const polyTopoChangeMap&,
984  const labelList& changedFaces,
985  const Map<label>& pointsToRestore,
986  const Map<label>& facesToRestore,
987  const Map<label>& cellsToRestore
988  );
989 
990  // Merging coplanar faces and edges
991 
992  //- Merge coplanar faces. preserveFaces is != -1 for faces
993  // to be preserved
995  (
996  const scalar minCos,
997  const scalar concaveCos,
998  const labelHashSet& patchIDs,
999  const dictionary& motionDict,
1000  const labelList& preserveFaces
1001  );
1002 
1004  (
1005  removePoints& pointRemover,
1006  const boolList& pointCanBeDeleted
1007  );
1008 
1010  (
1011  removePoints& pointRemover,
1012  const labelList& facesToRestore
1013  );
1014 
1016  (
1017  const labelList& candidateFaces,
1018  const labelHashSet& set
1019  ) const;
1020 
1021  // Pick up faces of cells of faces in set.
1023  (
1024  const labelHashSet& set
1025  ) const;
1026 
1027  //- Merge edges, maintain mesh quality. Return global number
1028  // of edges merged
1030  (
1031  const scalar minCos,
1032  const dictionary& motionDict
1033  );
1034 
1035 
1036  // Debug/IO
1037 
1038  //- Debugging: check that all faces still obey start()>end()
1039  void checkData();
1040 
1041  static void testSyncPointList
1042  (
1043  const string& msg,
1044  const polyMesh& mesh,
1045  const List<scalar>& fld
1046  );
1047 
1048  static void testSyncPointList
1049  (
1050  const string& msg,
1051  const polyMesh& mesh,
1052  const List<point>& fld
1053  );
1054 
1055  //- Compare two lists over all boundary faces
1056  template<class T>
1058  (
1059  const scalar mergeDistance,
1060  const string&,
1061  const UList<T>&,
1062  const UList<T>&
1063  ) const;
1064 
1065  //- Print list according to (collected and) sorted coordinate
1066  template<class T>
1067  static void collectAndPrint
1068  (
1069  const UList<point>& points,
1070  const UList<T>& data
1071  );
1072 
1073  //- Determine master point for subset of points. If coupled
1074  // chooses only one
1076  (
1077  const polyMesh& mesh,
1078  const labelList& meshPoints
1079  );
1080 
1081  //- Determine master edge for subset of edges. If coupled
1082  // chooses only one
1084  (
1085  const polyMesh& mesh,
1086  const labelList& meshEdges
1087  );
1088 
1089  //- Print some mesh stats.
1090  void printMeshInfo(const bool, const string&) const;
1091 
1092  //- Replacement for Time::timeName() : return oldInstance (if
1093  // overwrite_)
1094  word timeName() const;
1095 
1096  //- Set instance of all local IOobjects
1097  void setInstance(const fileName&);
1098 
1099  //- Write mesh and all data
1100  bool write() const;
1101 
1102  //- Write refinement level as volScalarFields for postprocessing
1103  void dumpRefinementLevel() const;
1104 
1105  //- Debug: Write intersection information to OBJ format
1106  void dumpIntersections(const fileName& prefix) const;
1107 
1108  //- Do any one of above IO functions
1109  void write
1110  (
1111  const debugType debugFlags,
1112  const writeType writeFlags,
1113  const fileName&
1114  ) const;
1115 
1116  //- Helper: calculate average
1117  template<class T>
1118  static T gAverage
1119  (
1120  const PackedBoolList& isMasterElem,
1121  const UList<T>& values
1122  );
1123 
1124  //- Get/set write level
1125  static writeType writeLevel();
1126  static void writeLevel(const writeType);
1127 
1128  //- Get/set output level
1129  static outputType outputLevel();
1130  static void outputLevel(const outputType);
1131 
1132 
1133  //- Helper: convert wordList into bit pattern using provided
1134  // NamedEnum
1135  template<class Enum>
1136  static int readFlags(const Enum& namedEnum, const wordList&);
1137 
1138 
1139  // Member Operators
1140 
1141  //- Disallow default bitwise assignment
1142  void operator=(const meshRefinement&) = delete;
1143 };
1144 
1145 
1146 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
1147 
1148 } // End namespace Foam
1149 
1150 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
1151 
1152 #ifdef NoRepository
1153  #include "meshRefinementTemplates.C"
1154 #endif
1155 
1156 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
1157 
1158 #endif
1159 
1160 // ************************************************************************* //
const labelList & surfaceIndex() const
Per start-end edge the index of the surface hit.
void printMeshInfo(const bool, const string &) const
Print some mesh stats.
const refinementSurfaces & surfaces() const
Reference to surface search engines.
Class containing mesh-to-mesh mapping information after a mesh distribution where we send parts of me...
void checkZoneFaces() const
Debug helper: check faceZones are not on processor patches.
IOoutputType
Enumeration for what to output.
A class for handling file names.
Definition: fileName.H:79
const List< Tuple2< mapType, labelList > > & userFaceData() const
Additional face data that is maintained across.
autoPtr< polyTopoChangeMap > createBaffles(const labelList &ownPatch, const labelList &nbrPatch)
Create baffle for every internal face where ownPatch != -1.
const word & oldInstance() const
(points)instance of mesh upon construction
label countHits() const
Count number of intersections (local)
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
autoPtr< polyTopoChangeMap > zonify(const List< point > &insidePoints, const bool allowFreeStandingZoneFaces)
Put faces/cells into zones according to surface specification.
IOdebugType
Enumeration for what to debug.
bool overwrite() const
Overwrite the mesh?
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
label mergePatchFacesUndo(const scalar minCos, const scalar concaveCos, const labelHashSet &patchIDs, const dictionary &motionDict, const labelList &preserveFaces)
Merge coplanar faces. preserveFaces is != -1 for faces.
void topoChange(const polyTopoChangeMap &, const labelList &changedFaces)
Update for external change to mesh. changedFaces are in new mesh.
static const NamedEnum< IOoutputType, 1 > IOoutputTypeNames
void operator=(const meshRefinement &)=delete
Disallow default bitwise assignment.
Given list of cells to remove insert all the topology changes.
Definition: removeCells.H:59
autoPtr< polyTopoChangeMap > splitFaces(const labelList &splitFaces, const labelPairList &splits)
Split faces into two.
autoPtr< polyTopoChangeMap > dupNonManifoldPoints()
Find boundary points that connect to more than one cell.
labelList intersectedFaces() const
Get faces with intersection.
static writeType writeLevel()
Get/set write level.
autoPtr< polyTopoChangeMap > createZoneBaffles(const labelList &globalToMasterPatch, const labelList &globalToSlavePatch, List< labelPair > &)
Create baffles for faces straddling zoned surfaces. Return.
Encapsulates queries for volume refinement (&#39;refine all cells within shell&#39;).
set value to -1 any face that was refined
labelList intersectedPoints() const
Get points on surfaces with intersection and boundary faces.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:69
autoPtr< polyTopoChangeMap > splitMesh(const label nBufferLayers, const labelList &globalToMasterPatch, const labelList &globalToSlavePatch, const refinementParameters::cellSelectionPoints &selectionPoints)
Split off (with optional buffer layers) unreachable areas.
autoPtr< polyDistributionMap > balanceAndRefine(const string &msg, decompositionMethod &decomposer, fvMeshDistribute &distributor, const labelList &cellsToRefine, const scalar maxLoadUnbalance)
Balance before refining some cells.
labelList collectFaces(const labelList &candidateFaces, const labelHashSet &set) const
Takes mesh with &#39;baffles&#39; (= boundary faces sharing points). Determines for selected points on bounda...
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
const fvMesh & mesh() const
Reference to mesh.
IOwriteType
Enumeration for what to write.
Container for data on surfaces used for surface-driven refinement. Contains all the data about the le...
Mesh representing a set of points created from polyMesh.
Definition: pointMesh.H:48
label mergeEdgesUndo(const scalar minCos, const dictionary &motionDict)
Merge edges, maintain mesh quality. Return global number.
meshRefinement(fvMesh &mesh, const dictionary &refineDict, const scalar mergeDistance, const bool overwrite, refinementSurfaces &, const refinementFeatures &, const refinementRegions &)
Construct from components.
labelList refineCandidates(const List< point > &insidePoints, const scalar curvature, const scalar planarAngle, const bool featureRefinement, const bool featureDistanceRefinement, const bool internalRefinement, const bool surfaceRefinement, const bool curvatureRefinement, const bool gapRefinement, const label maxGlobalCells, const label maxLocalCells) const
Calculate list of cells to refine.
Encapsulates queries for features.
static void updateList(const labelList &newToOld, const T &nullValue, List< T > &elems)
Helper: reorder list according to map.
A list of faces which address into the list of points.
Sends/receives parts of mesh+fvfields to neighbouring processors. Used in load balancing.
autoPtr< polyTopoChangeMap > refine(const labelList &cellsToRefine)
Refine some cells.
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().x()<< ' ';}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().y()<< ' ';}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){ gmvFile<< iter().position().z()<< ' ';}gmvFile<< nl;forAll(lagrangianScalarNames, i){ const word &name=lagrangianScalarNames[i];IOField< scalar > fld(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
ClassName("meshRefinement")
Runtime type information.
const pointField & points
static void weightedSum(const polyMesh &mesh, const PackedBoolList &isMasterEdge, const labelList &meshPoints, const edgeList &edges, const scalarField &edgeWeights, const Field< Type > &data, Field< Type > &sum)
Helper: weighted sum (over all subset of mesh points) by.
static void calculateEdgeWeights(const polyMesh &mesh, const PackedBoolList &isMasterEdge, const labelList &meshPoints, const edgeList &edges, scalarField &edgeWeights, scalarField &invSumWeight)
Helper: calculate edge weights (1/length)
Refinement of (split) hexes using polyTopoChange.
Definition: hexRef8.H:64
A class for handling words, derived from string.
Definition: word.H:59
static PackedBoolList getMasterEdges(const polyMesh &mesh, const labelList &meshEdges)
Determine master edge for subset of edges. If coupled.
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Simple container to keep together snap specific information.
autoPtr< polyDistributionMap > refineAndBalance(const string &msg, decompositionMethod &decomposer, fvMeshDistribute &distributor, const labelList &cellsToRefine, const scalar maxLoadUnbalance)
Refine some cells and rebalance.
static const NamedEnum< IOdebugType, 5 > IOdebugTypeNames
static void checkCoupledFaceZones(const polyMesh &)
Helper function: check that face zones are synced.
void testSyncBoundaryFaceList(const scalar mergeDistance, const string &, const UList< T > &, const UList< T > &) const
Compare two lists over all boundary faces.
Abstract base class for decomposition.
Container with cells to refine. Refinement given as single direction.
Definition: refineCell.H:56
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 refinementFeatures & features() const
Reference to feature edge mesh.
word timeName() const
Replacement for Time::timeName() : return oldInstance (if.
autoPtr< polyTopoChangeMap > doRemovePoints(removePoints &pointRemover, const boolList &pointCanBeDeleted)
mapType
Enumeration for how the userdata is to be mapped upon refinement.
const hexRef8 & meshCutter() const
Reference to meshcutting engine.
void checkData()
Debugging: check that all faces still obey start()>end()
insidePoints((1 2 3))
static int readFlags(const Enum &namedEnum, const wordList &)
Helper: convert wordList into bit pattern using provided.
Database for solution and other reduced data.
Definition: data.H:51
label addMeshedPatch(const word &name, const dictionary &)
Add patch originating from meshing. Update meshedPatches_.
labelList meshedPatches() const
Get patchIDs for patches added in addMeshedPatch.
autoPtr< polyDistributionMap > balance(const bool keepZoneFaces, const bool keepBaffles, const scalarField &cellWeights, decompositionMethod &decomposer, fvMeshDistribute &distributor)
Redecompose according to cell count.
void baffleAndSplitMesh(const bool handleSnapProblems, const snapParameters &snapParams, const bool useTopologicalSnapDetection, const bool removeEdgeConnectedCells, const scalarField &perpendicularAngle, const bool mergeFreeStanding, const scalar freeStandingAngle, const dictionary &motionDict, Time &runTime, const labelList &globalToMasterPatch, const labelList &globalToSlavePatch, const refinementParameters::cellSelectionPoints &selectionPoints)
Split off unreachable areas of mesh.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
static autoPtr< indirectPrimitivePatch > makePatch(const polyMesh &, const labelList &)
Create patch from set of patches.
static outputType outputLevel()
Get/set output level.
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
scalar mergeDistance() const
Helper class which maintains intersections of (changing) mesh with (static) surfaces.
void setInstance(const fileName &)
Set instance of all local IOobjects.
static void collectAndPrint(const UList< point > &points, const UList< T > &data)
Print list according to (collected and) sorted coordinate.
void distribute(const polyDistributionMap &)
Update local numbering for mesh redistribution.
A bit-packed bool list.
static T gAverage(const PackedBoolList &isMasterElem, const UList< T > &values)
Helper: calculate average.
autoPtr< polyTopoChangeMap > splitMeshRegions(const labelList &globalToMasterPatch, const labelList &globalToSlavePatch, const refinementParameters::cellSelectionPoints &selectionPoints)
Split mesh according to selectionPoints.
void storeData(const labelList &pointsToStore, const labelList &facesToStore, const labelList &cellsToStore)
Signal points/face/cells for which to store data.
static void testSyncPointList(const string &msg, const polyMesh &mesh, const List< scalar > &fld)
have slaves (upon refinement) from master
labelList growFaceCellFace(const labelHashSet &set) const
const refinementRegions & shells() const
Reference to refinement shells (regions)
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:95
Direct mesh changes based on v1.3 polyTopoChange syntax.
void selectSeparatedCoupledFaces(boolList &) const
Select coupled faces that are not collocated.
static void findRegions(const polyMesh &, labelList &cellRegion, const vector &perturbVec, const refinementParameters::cellSelectionPoints &selectionPoints)
Find regions points are in and update cellRegion.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
const doubleScalar e
Elementary charge.
Definition: doubleScalar.H:105
Removes selected points from mesh and updates faces using these points.
Definition: removePoints.H:58
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:76
bool isGap(const scalar, const vector &, const vector &, const vector &, const vector &) const
Is local topology a small gap?
Contains information about location on a triSurface.
bool write() const
Write mesh and all data.
static tmp< pointVectorField > makeDisplacementField(const pointMesh &pMesh, const labelList &adaptPatchIDs)
Helper function to make a pointVectorField with correct.
A class for managing temporary objects.
Definition: PtrList.H:53
static PackedBoolList getMasterPoints(const polyMesh &mesh, const labelList &meshPoints)
Determine master point for subset of points. If coupled.
autoPtr< polyTopoChangeMap > mergeBaffles(const List< labelPair > &)
Merge baffles. Gets pairs of faces.
void dumpIntersections(const fileName &prefix) const
Debug: Write intersection information to OBJ format.
static label findRegion(const polyMesh &, const labelList &cellRegion, const vector &perturbVec, const point &location)
Find region point is in. Uses optional perturbation to re-test.
void dumpRefinementLevel() const
Write refinement level as volScalarFields for postprocessing.
bool isNormalGap(const scalar, const vector &, const vector &, const vector &, const vector &) const
Is local topology a small gap normal to the test vector.
Namespace for OpenFOAM.
static const NamedEnum< IOwriteType, 5 > IOwriteTypeNames
A HashTable to objects of type <T> with a label key.
Definition: Map.H:49
Class to hold the points to select cells inside and outside.
autoPtr< polyTopoChangeMap > doRestorePoints(removePoints &pointRemover, const labelList &facesToRestore)