FreeFOAM The Cross-Platform CFD Toolkit
meshRefinement.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) 1991-2010 OpenCFD Ltd.
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 <dynamicMesh/hexRef8.H>
47 #include <OpenFOAM/mapPolyMesh.H>
48 #include <OpenFOAM/autoPtr.H>
49 #include <OpenFOAM/labelPair.H>
52 #include <OpenFOAM/Tuple2.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 shellSurfaces;
66 class removeCells;
67 class featureEdgeMesh;
68 class fvMeshDistribute;
69 class searchableSurface;
70 class regionSplit;
71 class globalIndex;
72 
73 /*---------------------------------------------------------------------------*\
74  Class meshRefinement Declaration
75 \*---------------------------------------------------------------------------*/
76 
78 {
79 
80 public:
81 
82  // Public data types
83 
84  //- Enumeration for debug dumping
85  enum writeFlag
86  {
87  MESH = 1,
90  };
91 
92 
93  //- Enumeration for how the userdata is to be mapped upon refinement.
94  enum mapType
95  {
96  MASTERONLY = 1,
97  KEEPALL = 2,
98  REMOVE = 4
99  };
100 
101 
102 private:
103 
104  // Private data
105 
106  //- Reference to mesh
107  fvMesh& mesh_;
108 
109  //- tolerance used for sorting coordinates (used in 'less' routine)
110  const scalar mergeDistance_;
111 
112  //- overwrite the mesh?
113  const bool overwrite_;
114 
115  //- Instance of mesh upon construction. Used when in overwrite_ mode.
116  const word oldInstance_;
117 
118  //- All surface-intersection interaction
119  const refinementSurfaces& surfaces_;
120 
121  //- All shell-refinement interaction
122  const shellSurfaces& shells_;
123 
124  //- refinement engine
125  hexRef8 meshCutter_;
126 
127  //- per cc-cc vector the index of the surface hit
128  labelIOList surfaceIndex_;
129 
130  //- user supplied face based data.
131  List<Tuple2<mapType, labelList> > userFaceData_;
132 
133  //- Meshed patches - are treated differently. Stored as wordList since
134  // order changes.
135  wordList meshedPatches_;
136 
137 
138  // Private Member Functions
139 
140  //- Reorder list according to map.
141  template<class T>
142  static void updateList
143  (
144  const labelList& newToOld,
145  const T& nullValue,
146  List<T>& elems
147  );
148 
149  //- Add patchfield of given type to all fields on mesh
150  template<class GeoField>
151  static void addPatchFields(fvMesh&, const word& patchFieldType);
152 
153  //- Reorder patchfields of all fields on mesh
154  template<class GeoField>
155  static void reorderPatchFields(fvMesh&, const labelList& oldToNew);
156 
157  //- Find out which faces have changed given cells (old mesh labels)
158  // that were marked for refinement.
159  static labelList getChangedFaces
160  (
161  const mapPolyMesh&,
162  const labelList& oldCellsToRefine
163  );
164 
165  //- Calculate coupled boundary end vector and refinement level
166  void calcNeighbourData
167  (
168  labelList& neiLevel,
169  pointField& neiCc
170  ) const;
171 
172  //- Find any intersection of surface. Store in surfaceIndex_.
173  void updateIntersections(const labelList& changedFaces);
174 
175  //- Remove cells. Put exposedFaces into exposedPatchIDs.
176  autoPtr<mapPolyMesh> doRemoveCells
177  (
178  const labelList& cellsToRemove,
179  const labelList& exposedFaces,
180  const labelList& exposedPatchIDs,
181  removeCells& cellRemover
182  );
183 
184  // Get cells which are inside any closed surface. Note that
185  // all closed surfaces
186  // will have already been oriented to have keepPoint outside.
187  labelList getInsideCells(const word&) const;
188 
189  // Do all to remove inside cells
190  autoPtr<mapPolyMesh> removeInsideCells
191  (
192  const string& msg,
193  const label exposedPatchI
194  );
195 
196  // For decomposeCombineRegions
197 
198  //- Used in decomposeCombineRegions. Given global region per cell
199  // determines master processor/cell for regions straddling
200  // procboundaries.
201  void getCoupledRegionMaster
202  (
203  const globalIndex& globalCells,
204  const boolList& blockedFace,
205  const regionSplit& globalRegion,
206  Map<label>& regionToMaster
207  ) const;
208 
209  //- Determine regions that are local to me or coupled ones that
210  // are owned by me. Determine representative location.
211  void calcLocalRegions
212  (
213  const globalIndex& globalCells,
214  const labelList& globalRegion,
215  const Map<label>& coupledRegionToMaster,
216  const scalarField& cellWeights,
217 
218  Map<label>& globalToLocalRegion,
219  pointField& localPoints,
220  scalarField& localWeights
221  ) const;
222 
223  //- Convert region into global index.
224  static label getShiftedRegion
225  (
226  const globalIndex& indexer,
227  const Map<label>& globalToLocalRegion,
228  const Map<label>& coupledRegionToShifted,
229  const label globalRegion
230  );
231 
232  //- helper: add element if not in list. Linear search.
233  static void addUnique(const label, labelList&);
234 
235  //- Calculate region connectivity. Major communication.
236  void calcRegionRegions
237  (
238  const labelList& globalRegion,
239  const Map<label>& globalToLocalRegion,
240  const Map<label>& coupledRegionToMaster,
241  labelListList& regionRegions
242  ) const;
243 
244 
245  // Refinement candidate selection
246 
247  //- Mark cell for refinement (if not already marked). Return false
248  // if refinelimit hit. Keeps running count (in nRefine) of cells
249  // marked for refinement
250  static bool markForRefine
251  (
252  const label markValue,
253  const label nAllowRefine,
254  label& cellValue,
255  label& nRefine
256  );
257 
258  //- Calculate list of cells to refine based on intersection of
259  // features.
260  label markFeatureRefinement
261  (
262  const point& keepPoint,
263  const PtrList<featureEdgeMesh>& featureMeshes,
264  const labelList& featureLevels,
265  const label nAllowRefine,
266 
268  label& nRefine
269  ) const;
270 
271  //- Mark cells for refinement-shells based refinement.
272  label markInternalRefinement
273  (
274  const label nAllowRefine,
276  label& nRefine
277  ) const;
278 
279  //- Collect faces that are intersected and whose neighbours aren't
280  // yet marked for refinement.
281  labelList getRefineCandidateFaces
282  (
283  const labelList& refineCell
284  ) const;
285 
286  //- Mark cells for surface intersection based refinement.
287  label markSurfaceRefinement
288  (
289  const label nAllowRefine,
290  const labelList& neiLevel,
291  const pointField& neiCc,
293  label& nRefine
294  ) const;
295 
296  //- Mark cell if local curvature > curvature or
297  // markDifferingRegions = true and intersections with different
298  // regions.
299  bool checkCurvature
300  (
301  const scalar curvature,
302  const label nAllowRefine,
303  const label surfaceLevel,
304  const vector& surfaceNormal,
305  const label cellI,
306  label& cellMaxLevel,
307  vector& cellMaxNormal,
309  label& nRefine
310  ) const;
311 
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  // Baffle handling
327 
328  //- Determine patches for baffles
329  void getBafflePatches
330  (
331  const labelList& globalToPatch,
332  const labelList& neiLevel,
333  const pointField& neiCc,
334  labelList& ownPatch,
335  labelList& neiPatch
336  ) const;
337 
338  //- Determine patch for baffle using some heuristic (and not
339  // surface)
340  label getBafflePatch
341  (
342  const labelList& facePatch,
343  const label faceI
344  ) const;
345 
346  //- Repatches external face or creates baffle for internal face
347  // with user specified patches (might be different for both sides).
348  // Returns label of added face.
349  label createBaffle
350  (
351  const label faceI,
352  const label ownPatch,
353  const label neiPatch,
354  polyTopoChange& meshMod
355  ) const;
356 
357  // Problem cell handling
358 
359  //- Helper function to mark face as being on 'boundary'. Used by
360  // markFacesOnProblemCells
361  void markBoundaryFace
362  (
363  const label faceI,
364  boolList& isBoundaryFace,
365  boolList& isBoundaryEdge,
366  boolList& isBoundaryPoint
367  ) const;
368 
369  void findNearest
370  (
371  const labelList& meshFaces,
372  List<pointIndexHit>& nearestInfo,
373  labelList& nearestSurface,
374  labelList& nearestRegion,
375  vectorField& nearestNormal
376  ) const;
377 
378  Map<label> findEdgeConnectedProblemCells
379  (
380  const scalarField& perpendicularAngle,
381  const labelList&
382  ) const;
383 
384  bool isCollapsedFace
385  (
386  const pointField&,
387  const pointField& neiCc,
388  const scalar minFaceArea,
389  const scalar maxNonOrtho,
390  const label faceI
391  ) const;
392 
393  bool isCollapsedCell
394  (
395  const pointField&,
396  const scalar volFraction,
397  const label cellI
398  ) const;
399 
400  //- Returns list with for every internal face -1 or the patch
401  // they should be baffled into. If removeEdgeConnectedCells is set
402  // removes cells based on perpendicularAngle.
403  labelList markFacesOnProblemCells
404  (
405  const dictionary& motionDict,
406  const bool removeEdgeConnectedCells,
407  const scalarField& perpendicularAngle,
408  const labelList& globalToPatch
409  ) const;
410 
412  //labelList markFacesOnProblemCellsGeometric
413  //(
414  // const dictionary& motionDict
415  //) const;
416 
417 
418  // Baffle merging
419 
420  //- Extract those baffles (duplicate) faces that are on the edge
421  // of a baffle region. These are candidates for merging.
422  List<labelPair> filterDuplicateFaces(const List<labelPair>&) const;
423 
424 
425  // Zone handling
426 
427  //- Finds zone per cell for cells inside closed named surfaces.
428  // (uses geometric test for insideness)
429  // Adapts namedSurfaceIndex so all faces on boundary of cellZone
430  // have corresponding faceZone.
431  void findCellZoneGeometric
432  (
433  const labelList& closedNamedSurfaces,
434  labelList& namedSurfaceIndex,
435  const labelList& surfaceToCellZone,
436  labelList& cellToZone
437  ) const;
438 
439  //- Determines cell zone from cell region information.
440  bool calcRegionToZone
441  (
442  const label surfZoneI,
443  const label ownRegion,
444  const label neiRegion,
445 
446  labelList& regionToCellZone
447  ) const;
448 
449  //- Finds zone per cell. Uses topological walk with all faces
450  // marked in namedSurfaceIndex regarded as blocked.
451  void findCellZoneTopo
452  (
453  const point& keepPoint,
454  const labelList& namedSurfaceIndex,
455  const labelList& surfaceToCellZone,
456  labelList& cellToZone
457  ) const;
458 
459  void makeConsistentFaceIndex
460  (
461  const labelList& cellToZone,
462  labelList& namedSurfaceIndex
463  ) const;
464 
465 
466  //- Disallow default bitwise copy construct
468 
469  //- Disallow default bitwise assignment
470  void operator=(const meshRefinement&);
471 
472 public:
473 
474  //- Runtime type information
475  ClassName("meshRefinement");
476 
477 
478  // Constructors
479 
480  //- Construct from components
482  (
483  fvMesh& mesh,
484  const scalar mergeDistance,
485  const bool overwrite,
486  const refinementSurfaces&,
487  const shellSurfaces&
488  );
489 
490 
491  // Member Functions
492 
493  // Access
494 
495  //- reference to mesh
496  const fvMesh& mesh() const
497  {
498  return mesh_;
499  }
501  {
502  return mesh_;
503  }
504 
505  scalar mergeDistance() const
506  {
507  return mergeDistance_;
508  }
509 
510  //- Overwrite the mesh?
511  bool overwrite() const
512  {
513  return overwrite_;
514  }
515 
516  //- (points)instance of mesh upon construction
517  const word& oldInstance() const
518  {
519  return oldInstance_;
520  }
521 
522  //- reference to surface search engines
524  {
525  return surfaces_;
526  }
527 
528  //- reference to refinement shells (regions)
529  const shellSurfaces& shells() const
530  {
531  return shells_;
532  }
533 
534  //- reference to meshcutting engine
535  const hexRef8& meshCutter() const
536  {
537  return meshCutter_;
538  }
539 
540  //- per start-end edge the index of the surface hit
541  const labelList& surfaceIndex() const
542  {
543  return surfaceIndex_;
544  }
545 
547  {
548  return surfaceIndex_;
549  }
550 
551  //- Additional face data that is maintained across
552  // topo changes. Every entry is a list over all faces.
553  // Bit of a hack. Additional flag to say whether to maintain master
554  // only (false) or increase set to account for face-from-face.
556  {
557  return userFaceData_;
558  }
559 
561  {
562  return userFaceData_;
563  }
564 
565 
566  // Other
567 
568  //- Count number of intersections (local)
569  label countHits() const;
570 
571  //- Helper function to get decomposition such that all connected
572  // regions get moved onto one processor. Used to prevent baffles
573  // straddling processor boundaries. explicitConnections is to
574  // keep pairs of non-coupled boundary faces together
575  // (e.g. to keep baffles together)
577  (
578  const scalarField& cellWeights,
579  const boolList& blockedFace,
580  const List<labelPair>& explicitConnections,
582  ) const;
583 
584  //- Redecompose according to cell count
585  // keepZoneFaces : find all faceZones from zoned surfaces and keep
586  // owner and neighbour together
587  // keepBaffles : find all baffles and keep them together
589  (
590  const bool keepZoneFaces,
591  const bool keepBaffles,
592  const scalarField& cellWeights,
593  decompositionMethod& decomposer,
594  fvMeshDistribute& distributor
595  );
596 
597  //- Get faces with intersection.
598  labelList intersectedFaces() const;
599 
600  //- Get points on surfaces with intersection and boundary faces.
602 
603  //- Create patch from set of patches
605  (
606  const polyMesh&,
607  const labelList&
608  );
609 
610  //- Helper function to make a pointVectorField with correct
611  // bcs for mesh movement:
612  // - adaptPatchIDs : fixedValue
613  // - processor : calculated (so free to move)
614  // - cyclic/wedge/symmetry : slip
615  // - other : slip
617  (
618  const pointMesh& pMesh,
619  const labelList& adaptPatchIDs
620  );
621 
622  //- Helper function: check that face zones are synced
623  static void checkCoupledFaceZones(const polyMesh&);
624 
625 
626  // Refinement
627 
628  //- Calculate list of cells to refine.
630  (
631  const point& keepPoint,
632  const scalar curvature,
633 
634  const PtrList<featureEdgeMesh>& featureMeshes,
635  const labelList& featureLevels,
636 
637  const bool featureRefinement,
638  const bool internalRefinement,
639  const bool surfaceRefinement,
640  const bool curvatureRefinement,
641  const label maxGlobalCells,
642  const label maxLocalCells
643  ) const;
644 
645  //- Refine some cells
646  autoPtr<mapPolyMesh> refine(const labelList& cellsToRefine);
647 
648  //- Refine some cells and rebalance
650  (
651  const string& msg,
652  decompositionMethod& decomposer,
653  fvMeshDistribute& distributor,
654  const labelList& cellsToRefine,
655  const scalar maxLoadUnbalance
656  );
657 
658  //- Balance before refining some cells
660  (
661  const string& msg,
662  decompositionMethod& decomposer,
663  fvMeshDistribute& distributor,
664  const labelList& cellsToRefine,
665  const scalar maxLoadUnbalance
666  );
667 
668 
669  // Baffle handling
670 
671  //- Split off unreachable areas of mesh.
672  void baffleAndSplitMesh
673  (
674  const bool handleSnapProblems,
675  const bool removeEdgeConnectedCells,
676  const scalarField& perpendicularAngle,
677  const bool mergeFreeStanding,
678  const dictionary& motionDict,
679  Time& runTime,
680  const labelList& globalToPatch,
681  const point& keepPoint
682  );
683 
684  //- Split off (with optional buffer layers) unreachable areas
685  // of mesh. Does not introduce baffles.
687  (
688  const label nBufferLayers,
689  const labelList& globalToPatch,
690  const point& keepPoint
691  );
692 
693  //- Find boundary points that connect to more than one cell
694  // region and split them.
696 
697  //- Create baffle for every internal face where ownPatch != -1.
698  // External faces get repatched according to ownPatch (neiPatch
699  // should be -1 for these)
701  (
702  const labelList& ownPatch,
703  const labelList& neiPatch
704  );
705 
706  //- Return a list of coupled face pairs, i.e. faces that
707  // use the same vertices.
708  List<labelPair> getDuplicateFaces(const labelList& testFaces) const;
709 
710  //- Merge baffles. Gets pairs of faces.
712 
713  //- Put faces/cells into zones according to surface specification.
714  // Returns null if no zone surfaces present. Region containing
715  // the keepPoint will not be put into a cellZone.
717  (
718  const point& keepPoint,
719  const bool allowFreeStandingZoneFaces
720  );
721 
722 
723  // Other topo changes
724 
725  //- Helper:add patch to mesh. Update all registered fields.
726  // Use addMeshedPatch to add patches originating from surfaces.
727  static label addPatch(fvMesh&, const word& name, const word& type);
728 
729  //- Add patch originating from meshing. Update meshedPatches_.
730  label addMeshedPatch(const word& name, const word& type);
731 
732  //- Get patchIDs for patches added in addMeshedPatch.
733  labelList meshedPatches() const;
734 
735  //- Split mesh. Keep part containing point.
736  autoPtr<mapPolyMesh> splitMeshRegions(const point& keepPoint);
737 
738  //- Update local numbering for mesh redistribution
739  void distribute(const mapDistributePolyMesh&);
740 
741  //- Update for external change to mesh. changedFaces are in new mesh
742  // face labels.
743  void updateMesh
744  (
745  const mapPolyMesh&,
746  const labelList& changedFaces
747  );
748 
749 
750  // Restoring : is where other processes delete and reinsert data.
751 
752  //- Signal points/face/cells for which to store data
753  void storeData
754  (
755  const labelList& pointsToStore,
756  const labelList& facesToStore,
757  const labelList& cellsToStore
758  );
759 
760  //- Update local numbering + undo
761  // Data to restore given as new pointlabel + stored pointlabel
762  // (i.e. what was in pointsToStore)
763  void updateMesh
764  (
765  const mapPolyMesh&,
766  const labelList& changedFaces,
767  const Map<label>& pointsToRestore,
768  const Map<label>& facesToRestore,
769  const Map<label>& cellsToRestore
770  );
771 
772 
773  //- Merge faces on the same patch (usually from exposing refinement)
774  // Returns global number of faces merged.
775  label mergePatchFaces
776  (
777  const scalar minCos,
778  const scalar concaveCos,
779  const labelList& patchIDs
780  );
781 
782  //- Remove points not used by any face or points used
783  // by only two faces where the edges are in line
784  autoPtr<mapPolyMesh> mergeEdges(const scalar minCos);
785 
786 
787  // Debug/IO
788 
789  //- Debugging: check that all faces still obey start()>end()
790  void checkData();
791 
792  //- Compare two lists over all boundary faces
793  template<class T>
795  (
796  const scalar mergeDistance,
797  const string&,
798  const UList<T>&,
799  const UList<T>&
800  ) const;
801 
802  //- Print some mesh stats.
803  void printMeshInfo(const bool, const string&) const;
804 
805  //- Replacement for Time::timeName() : return oldInstance (if
806  // overwrite_)
807  word timeName() const;
808 
809  //- Set instance of all local IOobjects
810  void setInstance(const fileName&);
811 
812  //- Write mesh and all data
813  bool write() const;
814 
815  //- Write refinement level as volScalarFields for postprocessing
816  void dumpRefinementLevel() const;
817 
818  //- Debug: Write intersection information to OBJ format
819  void dumpIntersections(const fileName& prefix) const;
820 
821  //- Do any one of above IO functions. flag is combination of
822  // writeFlag values.
823  void write(const label flag, const fileName&) const;
824 };
825 
826 
827 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
828 
829 } // End namespace Foam
830 
831 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
832 
833 #ifdef NoRepository
835 #endif
836 
837 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
838 
839 #endif
840 
841 // ************************ vim: set sw=4 sts=4 et: ************************ //