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