1 /*---------------------------------------------------------------------------*\
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
5 \\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
7 -------------------------------------------------------------------------------
9 This file is part of OpenFOAM.
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 the
13 Free Software Foundation; either version 2 of the License, or (at your
14 option) any later version.
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
21 You should have received a copy of the GNU General Public License
22 along with OpenFOAM; if not, write to the Free Software Foundation,
23 Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
29 Helper class which maintains intersections of (changing) mesh with
33 - per face any intersections of the cc-cc segment with any of the surfaces
37 meshRefinementBaffles.C
39 meshRefinementProblemCells.C
40 meshRefinementRefine.C
42 \*---------------------------------------------------------------------------*/
44 #ifndef meshRefinement_H
45 #define meshRefinement_H
48 #include "mapPolyMesh.H"
50 #include "labelPair.H"
51 #include "indirectPrimitivePatch.H"
52 #include "pointFieldsFwd.H"
54 #include "pointIndexHit.H"
56 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
61 // Class forward declarations
63 class mapDistributePolyMesh;
64 class decompositionMethod;
65 class refinementSurfaces;
68 class featureEdgeMesh;
69 class fvMeshDistribute;
70 class searchableSurface;
74 /*---------------------------------------------------------------------------*\
75 Class meshRefinement Declaration
76 \*---------------------------------------------------------------------------*/
85 //- Enumeration for debug dumping
94 //- Enumeration for how the userdata is to be mapped upon refinement.
97 MASTERONLY = 1, /*!< maintain master only */
98 KEEPALL = 2, /*!< have slaves (upon refinement) from master */
99 REMOVE = 4 /*!< set value to -1 any face that was refined */
107 //- Reference to mesh
110 //- tolerance used for sorting coordinates (used in 'less' routine)
111 const scalar mergeDistance_;
113 //- overwrite the mesh?
114 const bool overwrite_;
116 //- Instance of mesh upon construction. Used when in overwrite_ mode.
117 const word oldInstance_;
119 //- All surface-intersection interaction
120 const refinementSurfaces& surfaces_;
122 //- All shell-refinement interaction
123 const shellSurfaces& shells_;
125 //- refinement engine
128 //- per cc-cc vector the index of the surface hit
129 labelIOList surfaceIndex_;
131 //- user supplied face based data.
132 List<Tuple2<mapType, labelList> > userFaceData_;
134 //- Meshed patches - are treated differently. Stored as wordList since
136 wordList meshedPatches_;
139 // Private Member Functions
141 //- Reorder list according to map.
143 static void updateList
145 const labelList& newToOld,
150 //- Add patchfield of given type to all fields on mesh
151 template<class GeoField>
152 static void addPatchFields(fvMesh&, const word& patchFieldType);
154 //- Reorder patchfields of all fields on mesh
155 template<class GeoField>
156 static void reorderPatchFields(fvMesh&, const labelList& oldToNew);
158 //- Find out which faces have changed given cells (old mesh labels)
159 // that were marked for refinement.
160 static labelList getChangedFaces
163 const labelList& oldCellsToRefine
166 //- Calculate coupled boundary end vector and refinement level
167 void calcNeighbourData
173 //- Find any intersection of surface. Store in surfaceIndex_.
174 void updateIntersections(const labelList& changedFaces);
176 //- Remove cells. Put exposedFaces into exposedPatchIDs.
177 autoPtr<mapPolyMesh> doRemoveCells
179 const labelList& cellsToRemove,
180 const labelList& exposedFaces,
181 const labelList& exposedPatchIDs,
182 removeCells& cellRemover
185 // Get cells which are inside any closed surface. Note that
186 // all closed surfaces
187 // will have already been oriented to have keepPoint outside.
188 labelList getInsideCells(const word&) const;
190 // Do all to remove inside cells
191 autoPtr<mapPolyMesh> removeInsideCells
194 const label exposedPatchI
197 // For decomposeCombineRegions
199 //- Used in decomposeCombineRegions. Given global region per cell
200 // determines master processor/cell for regions straddling
202 void getCoupledRegionMaster
204 const globalIndex& globalCells,
205 const boolList& blockedFace,
206 const regionSplit& globalRegion,
207 Map<label>& regionToMaster
210 //- Determine regions that are local to me or coupled ones that
211 // are owned by me. Determine representative location.
212 void calcLocalRegions
214 const globalIndex& globalCells,
215 const labelList& globalRegion,
216 const Map<label>& coupledRegionToMaster,
218 Map<label>& globalToLocalRegion,
219 pointField& localPoints
222 //- Convert region into global index.
223 static label getShiftedRegion
225 const globalIndex& indexer,
226 const Map<label>& globalToLocalRegion,
227 const Map<label>& coupledRegionToShifted,
228 const label globalRegion
231 //- helper: add element if not in list. Linear search.
232 static void addUnique(const label, labelList&);
234 //- Calculate region connectivity. Major communication.
235 void calcRegionRegions
237 const labelList& globalRegion,
238 const Map<label>& globalToLocalRegion,
239 const Map<label>& coupledRegionToMaster,
240 labelListList& regionRegions
244 // Refinement candidate selection
246 //- Mark cell for refinement (if not already marked). Return false
247 // if refinelimit hit. Keeps running count (in nRefine) of cells
248 // marked for refinement
249 static bool markForRefine
251 const label markValue,
252 const label nAllowRefine,
257 //- Calculate list of cells to refine based on intersection of
259 label markFeatureRefinement
261 const point& keepPoint,
262 const PtrList<featureEdgeMesh>& featureMeshes,
263 const labelList& featureLevels,
264 const label nAllowRefine,
266 labelList& refineCell,
270 //- Mark cells for refinement-shells based refinement.
271 label markInternalRefinement
273 const label nAllowRefine,
274 labelList& refineCell,
278 //- Collect faces that are intersected and whose neighbours aren't
279 // yet marked for refinement.
280 labelList getRefineCandidateFaces
282 const labelList& refineCell
285 //- Mark cells for surface intersection based refinement.
286 label markSurfaceRefinement
288 const label nAllowRefine,
289 const labelList& neiLevel,
290 const pointField& neiCc,
291 labelList& refineCell,
295 //- Mark cell if local curvature > curvature or
296 // markDifferingRegions = true and intersections with different
300 const scalar curvature,
301 const label nAllowRefine,
302 const label surfaceLevel,
303 const vector& surfaceNormal,
306 vector& cellMaxNormal,
307 labelList& refineCell,
312 //- Mark cells for surface curvature based refinement. Marks if
313 // local curvature > curvature or if on different regions
314 // (markDifferingRegions)
315 label markSurfaceCurvatureRefinement
317 const scalar curvature,
318 const label nAllowRefine,
319 const labelList& neiLevel,
320 const pointField& neiCc,
321 labelList& refineCell,
327 //- Determine patches for baffles
328 void getBafflePatches
330 const labelList& globalToPatch,
331 const labelList& neiLevel,
332 const pointField& neiCc,
337 //- Determine patch for baffle using some heuristic (and not
341 const labelList& facePatch,
345 //- Repatches external face or creates baffle for internal face
346 // with user specified patches (might be different for both sides).
347 // Returns label of added face.
351 const label ownPatch,
352 const label neiPatch,
353 polyTopoChange& meshMod
356 // Problem cell handling
358 //- Helper function to mark face as being on 'boundary'. Used by
359 // markFacesOnProblemCells
360 void markBoundaryFace
363 boolList& isBoundaryFace,
364 boolList& isBoundaryEdge,
365 boolList& isBoundaryPoint
370 const labelList& meshFaces,
371 List<pointIndexHit>& nearestInfo,
372 labelList& nearestSurface,
373 labelList& nearestRegion,
374 vectorField& nearestNormal
377 Map<label> findEdgeConnectedProblemCells
379 const scalarField& perpendicularAngle,
386 const pointField& neiCc,
387 const scalar minFaceArea,
388 const scalar maxNonOrtho,
395 const scalar volFraction,
399 //- Returns list with for every internal face -1 or the patch
400 // they should be baffled into. If removeEdgeConnectedCells is set
401 // removes cells based on perpendicularAngle.
402 labelList markFacesOnProblemCells
404 const dictionary& motionDict,
405 const bool removeEdgeConnectedCells,
406 const scalarField& perpendicularAngle,
407 const labelList& globalToPatch
410 ////- Initial test of marking faces using geometric information.
411 //labelList markFacesOnProblemCellsGeometric
413 // const dictionary& motionDict
419 //- Extract those baffles (duplicate) faces that are on the edge
420 // of a baffle region. These are candidates for merging.
421 List<labelPair> filterDuplicateFaces(const List<labelPair>&) const;
426 //- Finds zone per cell for cells inside closed named surfaces.
427 // (uses geometric test for insideness)
428 // Adapts namedSurfaceIndex so all faces on boundary of cellZone
429 // have corresponding faceZone.
430 void findCellZoneGeometric
432 const labelList& closedNamedSurfaces,
433 labelList& namedSurfaceIndex,
434 const labelList& surfaceToCellZone,
435 labelList& cellToZone
438 //- Determines cell zone from cell region information.
439 bool calcRegionToZone
441 const label surfZoneI,
442 const label ownRegion,
443 const label neiRegion,
445 labelList& regionToCellZone
448 //- Finds zone per cell. Uses topological walk with all faces
449 // marked in namedSurfaceIndex regarded as blocked.
450 void findCellZoneTopo
452 const point& keepPoint,
453 const labelList& namedSurfaceIndex,
454 const labelList& surfaceToCellZone,
455 labelList& cellToZone
458 void makeConsistentFaceIndex
460 const labelList& cellToZone,
461 labelList& namedSurfaceIndex
465 //- Disallow default bitwise copy construct
466 meshRefinement(const meshRefinement&);
468 //- Disallow default bitwise assignment
469 void operator=(const meshRefinement&);
473 //- Runtime type information
474 ClassName("meshRefinement");
479 //- Construct from components
483 const scalar mergeDistance,
484 const bool overwrite,
485 const refinementSurfaces&,
494 //- reference to mesh
495 const fvMesh& mesh() const
504 scalar mergeDistance() const
506 return mergeDistance_;
509 //- Overwrite the mesh?
510 bool overwrite() const
515 //- (points)instance of mesh upon construction
516 const word& oldInstance() const
521 //- reference to surface search engines
522 const refinementSurfaces& surfaces() const
527 //- reference to refinement shells (regions)
528 const shellSurfaces& shells() const
533 //- reference to meshcutting engine
534 const hexRef8& meshCutter() const
539 //- per start-end edge the index of the surface hit
540 const labelList& surfaceIndex() const
542 return surfaceIndex_;
545 labelList& surfaceIndex()
547 return surfaceIndex_;
550 //- Additional face data that is maintained across
551 // topo changes. Every entry is a list over all faces.
552 // Bit of a hack. Additional flag to say whether to maintain master
553 // only (false) or increase set to account for face-from-face.
554 const List<Tuple2<mapType, labelList> >& userFaceData() const
556 return userFaceData_;
559 List<Tuple2<mapType, labelList> >& userFaceData()
561 return userFaceData_;
567 //- Count number of intersections (local)
568 label countHits() const;
570 //- Helper function to get decomposition such that all connected
571 // regions get moved onto one processor. Used to prevent baffles
572 // straddling processor boundaries. explicitConnections is to
573 // keep pairs of non-coupled boundary faces together
574 // (e.g. to keep baffles together)
575 labelList decomposeCombineRegions
577 const boolList& blockedFace,
578 const List<labelPair>& explicitConnections,
582 //- Redecompose according to cell count
583 // keepZoneFaces : find all faceZones from zoned surfaces and keep
584 // owner and neighbour together
585 // keepBaffles : find all baffles and keep them together
586 autoPtr<mapDistributePolyMesh> balance
588 const bool keepZoneFaces,
589 const bool keepBaffles,
590 decompositionMethod& decomposer,
591 fvMeshDistribute& distributor
594 //- Get faces with intersection.
595 labelList intersectedFaces() const;
597 //- Get points on surfaces with intersection and boundary faces.
598 labelList intersectedPoints() const;
600 //- Create patch from set of patches
601 static autoPtr<indirectPrimitivePatch> makePatch
607 //- Helper function to make a pointVectorField with correct
608 // bcs for mesh movement:
609 // - adaptPatchIDs : fixedValue
610 // - processor : calculated (so free to move)
611 // - cyclic/wedge/symmetry : slip
613 static tmp<pointVectorField> makeDisplacementField
615 const pointMesh& pMesh,
616 const labelList& adaptPatchIDs
619 //- Helper function: check that face zones are synced
620 static void checkCoupledFaceZones(const polyMesh&);
625 //- Calculate list of cells to refine.
626 labelList refineCandidates
628 const point& keepPoint,
629 const scalar curvature,
631 const PtrList<featureEdgeMesh>& featureMeshes,
632 const labelList& featureLevels,
634 const bool featureRefinement,
635 const bool internalRefinement,
636 const bool surfaceRefinement,
637 const bool curvatureRefinement,
638 const label maxGlobalCells,
639 const label maxLocalCells
642 //- Refine some cells
643 autoPtr<mapPolyMesh> refine(const labelList& cellsToRefine);
645 //- Refine some cells and rebalance
646 autoPtr<mapDistributePolyMesh> refineAndBalance
649 decompositionMethod& decomposer,
650 fvMeshDistribute& distributor,
651 const labelList& cellsToRefine
657 //- Split off unreachable areas of mesh.
658 void baffleAndSplitMesh
660 const bool handleSnapProblems,
661 const bool removeEdgeConnectedCells,
662 const scalarField& perpendicularAngle,
663 const bool mergeFreeStanding,
664 const dictionary& motionDict,
666 const labelList& globalToPatch,
667 const point& keepPoint
670 //- Split off (with optional buffer layers) unreachable areas
671 // of mesh. Does not introduce baffles.
672 autoPtr<mapPolyMesh> splitMesh
674 const label nBufferLayers,
675 const labelList& globalToPatch,
676 const point& keepPoint
679 //- Find boundary points that connect to more than one cell
680 // region and split them.
681 autoPtr<mapPolyMesh> dupNonManifoldPoints();
683 //- Create baffle for every internal face where ownPatch != -1.
684 // External faces get repatched according to ownPatch (neiPatch
685 // should be -1 for these)
686 autoPtr<mapPolyMesh> createBaffles
688 const labelList& ownPatch,
689 const labelList& neiPatch
692 //- Return a list of coupled face pairs, i.e. faces that
693 // use the same vertices.
694 List<labelPair> getDuplicateFaces(const labelList& testFaces) const;
696 //- Merge baffles. Gets pairs of faces.
697 autoPtr<mapPolyMesh> mergeBaffles(const List<labelPair>&);
699 //- Put faces/cells into zones according to surface specification.
700 // Returns null if no zone surfaces present. Region containing
701 // the keepPoint will not be put into a cellZone.
702 autoPtr<mapPolyMesh> zonify(const point& keepPoint);
705 // Other topo changes
707 //- Helper:add patch to mesh. Update all registered fields.
708 // Use addMeshedPatch to add patches originating from surfaces.
709 static label addPatch(fvMesh&, const word& name, const word& type);
711 //- Add patch originating from meshing. Update meshedPatches_.
712 label addMeshedPatch(const word& name, const word& type);
714 //- Get patchIDs for patches added in addMeshedPatch.
715 labelList meshedPatches() const;
717 //- Split mesh. Keep part containing point.
718 autoPtr<mapPolyMesh> splitMeshRegions(const point& keepPoint);
720 //- Update local numbering for mesh redistribution
721 void distribute(const mapDistributePolyMesh&);
723 //- Update for external change to mesh. changedFaces are in new mesh
728 const labelList& changedFaces
732 // Restoring : is where other processes delete and reinsert data.
734 //- Signal points/face/cells for which to store data
737 const labelList& pointsToStore,
738 const labelList& facesToStore,
739 const labelList& cellsToStore
742 //- Update local numbering + undo
743 // Data to restore given as new pointlabel + stored pointlabel
744 // (i.e. what was in pointsToStore)
748 const labelList& changedFaces,
749 const Map<label>& pointsToRestore,
750 const Map<label>& facesToRestore,
751 const Map<label>& cellsToRestore
755 //- Merge faces on the same patch (usually from exposing refinement)
756 // Returns global number of faces merged.
757 label mergePatchFaces
760 const scalar concaveCos,
761 const labelList& patchIDs
764 //- Remove points not used by any face or points used
765 // by only two faces where the edges are in line
766 autoPtr<mapPolyMesh> mergeEdges(const scalar minCos);
771 //- Debugging: check that all faces still obey start()>end()
774 //- Compare two lists over all boundary faces
776 void testSyncBoundaryFaceList
778 const scalar mergeDistance,
784 //- Print some mesh stats.
785 void printMeshInfo(const bool, const string&) const;
787 //- Replacement for Time::timeName() : return oldInstance (if
789 word timeName() const;
791 //- Set instance of all local IOobjects
792 void setInstance(const fileName&);
794 //- Write mesh and all data
797 //- Write refinement level as volScalarFields for postprocessing
798 void dumpRefinementLevel() const;
800 //- Debug: Write intersection information to OBJ format
801 void dumpIntersections(const fileName& prefix) const;
803 //- Do any one of above IO functions. flag is combination of
805 void write(const label flag, const fileName&) const;
809 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
811 } // End namespace Foam
813 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
816 # include "meshRefinementTemplates.C"
819 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
823 // ************************************************************************* //