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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2017 OpenFOAM Foundation
9 Copyright (C) 2015-2022 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27Class
28 Foam::meshRefinement
29
30Description
31 Helper class which maintains intersections of (changing) mesh with
32 (static) surfaces.
33
34 Maintains
35 - per face any intersections of the cc-cc segment with any of the surfaces
36
37SourceFiles
38 meshRefinement.C
39 meshRefinementBaffles.C
40 meshRefinementGapRefine.C
41 meshRefinementMerge.C
42 meshRefinementProblemCells.C
43 meshRefinementRefine.C
44
45\*---------------------------------------------------------------------------*/
46
47#ifndef Foam_meshRefinement_H
48#define Foam_meshRefinement_H
49
50#include "hexRef8.H"
51#include "mapPolyMesh.H"
52#include "autoPtr.H"
53#include "labelPairHashes.H"
55#include "pointFieldsFwd.H"
56#include "Tuple2.H"
57#include "pointIndexHit.H"
58#include "wordPairHashes.H"
59#include "surfaceZonesInfo.H"
60#include "volumeType.H"
61#include "DynamicField.H"
62#include "coordSetWriter.H"
63
64// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
65
66namespace Foam
67{
68
69// Forward Declarations
70class fvMesh;
71class mapDistributePolyMesh;
72class mapDistribute;
73class decompositionMethod;
74class refinementSurfaces;
75class refinementFeatures;
76class shellSurfaces;
77class removeCells;
78class fvMeshDistribute;
79class removePoints;
80class localPointRegion;
81class snapParameters;
82
83/*---------------------------------------------------------------------------*\
84 Class meshRefinement Declaration
85\*---------------------------------------------------------------------------*/
88{
89public:
90
91 // Public Data Types
92
93 //- Enumeration for what to debug. Used as a bit-pattern.
94 enum debugType
95 {
96 MESH = (1 << 0),
97 OBJINTERSECTIONS = (1 << 1),
98 FEATURESEEDS = (1 << 2),
99 ATTRACTION = (1 << 3),
100 LAYERINFO = (1 << 4)
101 };
103 static const Enum<debugType> debugTypeNames;
104
106 //enum outputType
107 //{
108 // OUTPUTLAYERINFO = (1 << 0)
109 //};
110 //
111 //static const Enum<outputType> outputTypeNames;
112
113 //- Enumeration for what to write. Used as a bit-pattern.
114 enum writeType
116 WRITEMESH = (1 << 0),
118 WRITELEVELS = (1 << 2),
119 WRITELAYERSETS = (1 << 3),
121 };
123 static const Enum<writeType> writeTypeNames;
124
125 //- Enumeration for how userdata is to be mapped upon refinement.
126 enum mapType
129 KEEPALL = 2,
131 };
132
133 //- Enumeration for what to do with co-planar patch faces on a single
134 // cell
135 enum FaceMergeType
137 NONE, // no merging
138 GEOMETRIC, // use feature angle
139 IGNOREPATCH // use feature angle, allow merging of different
140 // patches
141 };
142
143
144private:
145
146 // Static Data Members
147
148 //- Control of writing level
149 static writeType writeLevel_;
150
152 //static outputType outputLevel_;
153
154
155 // Private Data
156
157 //- Reference to mesh
158 fvMesh& mesh_;
159
160 //- Tolerance used for sorting coordinates (used in 'less' routine)
161 const scalar mergeDistance_;
162
163 //- Overwrite the mesh?
164 const bool overwrite_;
165
166 //- Instance of mesh upon construction. Used when in overwrite_ mode.
167 const word oldInstance_;
168
169 //- All surface-intersection interaction
170 const refinementSurfaces& surfaces_;
171
172 //- All feature-edge interaction
173 const refinementFeatures& features_;
174
175 //- All shell-refinement interaction
176 const shellSurfaces& shells_;
177
178 //- All limit-refinement interaction
179 const shellSurfaces& limitShells_;
180
181 //- Are we operating in test mode?
182 const bool dryRun_;
183
184 //- Refinement engine
185 hexRef8 meshCutter_;
186
187 //- Per cc-cc vector the index of the surface hit
188 labelIOList surfaceIndex_;
189
190
191 // For baffle merging
192
193 //- Original patch for baffle faces that used to be on
194 // coupled patches
195 Map<label> faceToCoupledPatch_;
196
197
198 //- User supplied face based data.
199 List<Tuple2<mapType, labelList>> userFaceData_;
200
201 //- Meshed patches - are treated differently. Stored as wordList since
202 // order changes.
203 wordList meshedPatches_;
204
205 //- FaceZone to master patch name
206 HashTable<word> faceZoneToMasterPatch_;
207
208 //- FaceZone to slave patch name
209 HashTable<word> faceZoneToSlavePatch_;
210
211 //- FaceZone to method to handle faces
213
214
215 // Private Member Functions
216
217 //- Add patchfield of given type to all fields on mesh
218 template<class GeoField>
219 static void addPatchFields(fvMesh&, const word& patchFieldType);
220
221 //- Reorder patchfields of all fields on mesh
222 template<class GeoField>
223 static void reorderPatchFields(fvMesh&, const labelList& oldToNew);
224
225 //- Find out which faces have changed given cells (old mesh labels)
226 // that were marked for refinement.
227 static labelList getChangedFaces
228 (
229 const mapPolyMesh&,
230 const labelList& oldCellsToRefine
231 );
232
233 // Debug: count number of master-faces (and do some checking
234 // for consistency)
235 label globalFaceCount(const labelList& elems) const;
236
237 //- Calculate coupled boundary end vector and refinement level
238 void calcNeighbourData(labelList& neiLevel, pointField& neiCc) const;
239
240 //- Calculate rays from cell-centre to cell-centre and corresponding
241 // min cell refinement level
242 void calcCellCellRays
243 (
244 const pointField& neiCc,
245 const labelList& neiLevel,
246 const labelList& testFaces,
247 pointField& start,
248 pointField& end,
249 labelList& minLevel
250 ) const;
251
252 //- Get cells which are inside any closed surface. Note that
253 // all closed surfaces
254 // will have already been oriented to have keepPoint outside.
255 labelList getInsideCells(const word&) const;
256
257 //- Do all to remove inside cells
258 autoPtr<mapPolyMesh> removeInsideCells
259 (
260 const string& msg,
261 const label exposedPatchi
262 );
263
264
265 // Refinement candidate selection
266
267 //- Mark cell for refinement (if not already marked). Return false
268 // if refinelimit hit. Keeps running count (in nRefine) of cells
269 // marked for refinement
270 static bool markForRefine
271 (
272 const label markValue,
273 const label nAllowRefine,
274 label& cellValue,
275 label& nRefine
276 );
277
278 //- Mark every cell with level of feature passing through it
279 // (or -1 if not passed through). Uses tracking.
280 void markFeatureCellLevel
281 (
282 const pointField& keepPoints,
283 labelList& maxFeatureLevel
284 ) const;
285
286 //- Calculate list of cells to refine based on intersection of
287 // features.
288 label markFeatureRefinement
289 (
290 const pointField& keepPoints,
291 const label nAllowRefine,
292
294 label& nRefine
295 ) const;
296
297 //- Mark cells for distance-to-feature based refinement.
298 label markInternalDistanceToFeatureRefinement
299 (
300 const label nAllowRefine,
302 label& nRefine
303 ) const;
304
305 //- Mark cells for refinement-shells based refinement.
306 label markInternalRefinement
307 (
308 const label nAllowRefine,
310 label& nRefine
311 ) const;
312
313 //- Unmark cells for refinement based on limit-shells. Return number
314 // of limited cells.
315 label unmarkInternalRefinement
316 (
318 label& nRefine
319 ) const;
320
321 //- Collect faces that are intersected and whose neighbours aren't
322 // yet marked for refinement.
323 labelList getRefineCandidateFaces
324 (
325 const labelList& refineCell
326 ) const;
327
328 //- Mark cells for surface intersection based refinement.
329 label markSurfaceRefinement
330 (
331 const label nAllowRefine,
332 const labelList& neiLevel,
333 const pointField& neiCc,
335 label& nRefine
336 ) const;
337
338 //- Collect cells intersected by the surface that are candidates
339 // for gap checking. Used inside markSurfaceGapRefinement
340 void collectGapCandidates
341 (
342 const shellSurfaces& shells,
343 const labelList& testFaces,
344 const labelList& neiLevel,
345 const pointField& neiCc,
346 labelList& cellToCompact,
347 labelList& bFaceToCompact,
348 List<FixedList<label, 3>>& shellGapInfo,
349 List<volumeType>& shellGapMode
350 ) const;
351 void collectGapCells
352 (
353 const scalar planarCos,
354
355 const List<FixedList<label, 3>>& extendedGapLevel,
356 const List<volumeType>& extendedGapMode,
357 const labelList& testFaces,
358 const pointField& start,
359 const pointField& end,
360
361 const labelList& cellToCompact,
362 const labelList& bFaceToCompact,
363 const List<FixedList<label, 3>>& shellGapInfo,
364 const List<volumeType>& shellGapMode,
365
366 const label nAllowRefine,
367 const labelList& neiLevel,
368 const pointField& neiCc,
369
371 label& nRefine
372 ) const;
373
374
375 //- Mark cells intersected by the surface if they are inside
376 // close gaps
377 label markSurfaceGapRefinement
378 (
379 const scalar planarCos,
380 const label nAllowRefine,
381 const labelList& neiLevel,
382 const pointField& neiCc,
383
385 label& nRefine
386 ) const;
387
388 //- Generate single ray from nearPoint in direction of nearNormal
389 label generateRays
390 (
391 const point& nearPoint,
392 const vector& nearNormal,
393 const FixedList<label, 3>& gapInfo,
394 const volumeType& mode,
395
396 const label cLevel,
397
398 DynamicField<point>& start,
400 ) const;
401
402 //- Generate pairs of rays through cell centre
403 // Each ray pair has start, end, and expected gap size
404 label generateRays
405 (
406 const bool useSurfaceNormal,
407
408 const point& nearPoint,
409 const vector& nearNormal,
410 const FixedList<label, 3>& gapInfo,
411 const volumeType& mode,
412
413 const point& cc,
414 const label cLevel,
415
416 DynamicField<point>& start,
418 DynamicField<scalar>& gapSize,
419
420 DynamicField<point>& start2,
422 DynamicField<scalar>& gapSize2
423 ) const;
424
425 //- Select candidate cells (cells inside a shell with gapLevel
426 // specified)
427 void selectGapCandidates
428 (
429 const labelList& refineCell,
430 const label nRefine,
431
432 labelList& cellMap,
433 labelList& gapShell,
434 List<FixedList<label, 3>>& shellGapInfo,
435 List<volumeType>& shellGapMode
436 ) const;
437
438 //- Merge gap information coming from shell and from surface
439 // (surface wins)
440 void mergeGapInfo
441 (
442 const FixedList<label, 3>& shellGapInfo,
443 const volumeType shellGapMode,
444 const FixedList<label, 3>& surfGapInfo,
445 const volumeType surfGapMode,
446 FixedList<label, 3>& gapInfo,
447 volumeType& gapMode
448 ) const;
449
450 //- Mark cells for non-surface intersection based gap refinement
451 label markInternalGapRefinement
452 (
453 const scalar planarCos,
454 const bool spreadGapSize,
455 const label nAllowRefine,
457 label& nRefine,
458 labelList& numGapCells,
459 scalarField& gapSize
460 ) const;
461
462 //- Refine cells containing small gaps
463 label markSmallFeatureRefinement
464 (
465 const scalar planarCos,
466 const label nAllowRefine,
467 const labelList& neiLevel,
468 const pointField& neiCc,
469
471 label& nRefine
472 ) const;
473
474 //- Helper: count number of normals1 that are in normals2
475 label countMatches
476 (
477 const List<point>& normals1,
478 const List<point>& normals2,
479 const scalar tol = 1e-6
480 ) const;
481
482 //- Mark cells for surface curvature based refinement. Marks if
483 // local curvature > curvature or if on different regions
484 // (markDifferingRegions)
485 label markSurfaceCurvatureRefinement
486 (
487 const scalar curvature,
488 const label nAllowRefine,
489 const labelList& neiLevel,
490 const pointField& neiCc,
492 label& nRefine
493 ) const;
494
495 //- Mark cell if local a gap topology or
496 bool checkProximity
497 (
498 const scalar planarCos,
499 const label nAllowRefine,
500
501 const label surfaceLevel,
502 const vector& surfaceLocation,
503 const vector& surfaceNormal,
504
505 const label celli,
506
507 label& cellMaxLevel,
508 vector& cellMaxLocation,
509 vector& cellMaxNormal,
510
512 label& nRefine
513 ) const;
514
515 //- Mark cells for surface proximity based refinement.
516 label markProximityRefinement
517 (
518 const scalar curvature,
519
520 // Per region the min and max cell level
521 const labelList& surfaceMinLevel,
522 const labelList& surfaceMaxLevel,
523
524 const label nAllowRefine,
525 const labelList& neiLevel,
526 const pointField& neiCc,
527
529 label& nRefine
530 ) const;
531
532 //- Mark cells for surface proximity based refinement.
533 label markProximityRefinementWave
534 (
535 const scalar planarCos,
536 const labelList& blockedSurfaces,
537 const label nAllowRefine,
538 const labelList& neiLevel,
539 const pointField& neiCc,
540
542 label& nRefine
543 ) const;
544
545
546 // Baffle handling
547
548 //- Get faces to repatch. Returns map from face to patch.
549 Map<labelPair> getZoneBafflePatches
550 (
551 const bool allowBoundary,
552 const labelList& globalToMasterPatch,
553 const labelList& globalToSlavePatch
554 ) const;
555
556 //- Calculate intersections. Return per face -1 or the global
557 // surface region
558 void getIntersections
559 (
560 const labelList& surfacesToTest,
561 const pointField& neiCc,
562 const labelList& testFaces,
563
564 labelList& globalRegion1,
565 labelList& globalRegion2
566 ) const;
567
568 //- Calculate intersections on zoned faces. Return per face -1
569 // or the global region of the surface and the orientation
570 // w.r.t. surface
571 void getIntersections
572 (
573 const labelList& surfacesToTest,
574 const pointField& neiCc,
575 const labelList& testFaces,
576
577 labelList& namedSurfaceRegion,
578 bitSet& posOrientation
579 ) const;
580
581 //- Determine patches for baffles
582 void getBafflePatches
583 (
584 const label nErodeCellZones,
585 const labelList& globalToMasterPatch,
586 const pointField& locationsInMesh,
587 const wordList& regionsInMesh,
588 const pointField& locationsOutsideMesh,
589 const bool exitIfLeakPath,
590 const refPtr<coordSetWriter>& leakPathFormatter,
591
592 const labelList& neiLevel,
593 const pointField& neiCc,
594 labelList& ownPatch,
595 labelList& neiPatch
596 ) const;
597
598 autoPtr<mapPolyMesh> splitMesh
599 (
600 const label nBufferLayers,
601 const labelList& globalToMasterPatch,
602 const labelList& globalToSlavePatch,
603 labelList& cellRegion,
604 labelList& ownPatch,
605 labelList& neiPatch
606 );
607
608 //- Repatches external face or creates baffle for internal face
609 // with user specified patches (might be different for both sides).
610 // Returns label of added face.
611 label createBaffle
612 (
613 const label facei,
614 const label ownPatch,
615 const label neiPatch,
616 polyTopoChange& meshMod
617 ) const;
618
619 //- Write leak path
620 static fileName writeLeakPath
621 (
622 const polyMesh& mesh,
623 const pointField& locationsInMesh,
624 const pointField& locationsOutsideMesh,
625 const boolList& blockedFace,
626 coordSetWriter& leakPathWriter
627 );
628
629
630 // Problem cell handling
631
632 //- Helper function to mark face as being on 'boundary'. Used by
633 // markFacesOnProblemCells
634 void markBoundaryFace
635 (
636 const label facei,
637 boolList& isBoundaryFace,
638 boolList& isBoundaryEdge,
639 boolList& isBoundaryPoint
640 ) const;
641
642 void findNearest
643 (
644 const labelList& meshFaces,
645 List<pointIndexHit>& nearestInfo,
646 labelList& nearestSurface,
647 labelList& nearestRegion,
648 vectorField& nearestNormal
649 ) const;
650
651 Map<label> findEdgeConnectedProblemCells
652 (
653 const scalarField& perpendicularAngle,
654 const labelList&
655 ) const;
656
657 bool isCollapsedFace
658 (
659 const pointField&,
660 const pointField& neiCc,
661 const scalar minFaceArea,
662 const scalar maxNonOrtho,
663 const label facei
664 ) const;
665
666 bool isCollapsedCell
667 (
668 const pointField&,
669 const scalar volFraction,
670 const label celli
671 ) const;
672
673 //- Returns list with for every internal face -1 or the patch
674 // they should be baffled into. If removeEdgeConnectedCells is set
675 // removes cells based on perpendicularAngle.
676 void markFacesOnProblemCells
677 (
678 const dictionary& motionDict,
679 const bool removeEdgeConnectedCells,
680 const scalarField& perpendicularAngle,
681 const labelList& globalToMasterPatch,
682
683 labelList& facePatch,
685 ) const;
686
687 //- Calculates for every face the nearest 'start' face. Any
688 // unreached face does not get set (faceToStart[facei] = -1)
689 void nearestFace
690 (
691 const labelUList& startFaces,
692 const bitSet& isBlockedFace,
693
695 labelList& faceToStart,
696 const label nIter = labelMax
697 ) const;
698
699 //- Calculates for every face the label of the nearest
700 // patch its zone. Any unreached face (disconnected mesh?) becomes
701 // adaptPatchIDs[0]
702 void nearestPatch
703 (
704 const labelList& adaptPatchIDs,
705 labelList& nearestPatch,
706 labelList& nearestZone
707 ) const;
708
709 //- Returns list with for every face the label of the nearest
710 // patch. Any unreached face (disconnected mesh?) becomes
711 // adaptPatchIDs[0]
712 labelList nearestPatch(const labelList& adaptPatchIDs) const;
713
714 //- Returns list with for every face the label of the nearest
715 // (global) region. Any unreached face (disconnected mesh?) becomes
716 // defaultRegion
717 labelList nearestIntersection
718 (
719 const labelList& surfacesToTest,
720 const label defaultRegion
721 ) const;
722
723 //- Returns list with for every internal face -1 or the patch
724 // they should be baffled into.
725 void markFacesOnProblemCellsGeometric
726 (
727 const snapParameters& snapParams,
728 const dictionary& motionDict,
729 const labelList& globalToMasterPatch,
730 const labelList& globalToSlavePatch,
731
732 labelList& facePatch,
734 ) const;
735
736
737 // Baffle merging
738
739 //- Extract those baffles (duplicate) faces that are on the edge
740 // of a baffle region. These are candidates for merging.
741 List<labelPair> freeStandingBaffles
742 (
743 const List<labelPair>&,
744 const scalar freeStandingAngle
745 ) const;
746
747
748 // Zone handling
749
750 //- Finds zone per cell for cells inside closed named surfaces.
751 // (uses geometric test for insideness)
752 // Adapts namedSurfaceRegion so all faces on boundary of cellZone
753 // have corresponding faceZone.
754 void findCellZoneGeometric
755 (
756 const pointField& neiCc,
757 const labelList& closedNamedSurfaces,
758 labelList& namedSurfaceRegion,
759 const labelList& surfaceToCellZone,
760 labelList& cellToZone
761 ) const;
762
763 //- Finds zone per cell for cells inside region for which name
764 // is specified.
765 void findCellZoneInsideWalk
766 (
767 const pointField& locationsInMesh,
768 const labelList& zonesInMesh,
769 const labelList& blockedFace, // per face -1 or some index >= 0
770 labelList& cellToZone
771 ) const;
772
773 //- Finds zone per cell for cells inside region for which name
774 // is specified.
775 void findCellZoneInsideWalk
776 (
777 const pointField& locationsInMesh,
778 const wordList& zoneNamesInMesh,
779 const labelList& faceToZone, // per face -1 or some index >= 0
780 labelList& cellToZone
781 ) const;
782
783 //- Determines cell zone from cell region information.
784 bool calcRegionToZone
785 (
786 const label backgroundZoneID,
787 const label surfZoneI,
788 const label ownRegion,
789 const label neiRegion,
790
791 labelList& regionToCellZone
792 ) const;
793
794 //- Finds zone per cell. Uses topological walk with all faces
795 // marked in unnamedSurfaceRegion (intersections with unnamed
796 // surfaces) and namedSurfaceRegion (intersections with named
797 // surfaces) regarded as blocked.
798 void findCellZoneTopo
799 (
800 const label backgroundZoneID,
801 const pointField& locationsInMesh,
802 const labelList& unnamedSurfaceRegion,
803 const labelList& namedSurfaceRegion,
804 const labelList& surfaceToCellZone,
805 labelList& cellToZone
806 ) const;
807
808 //- Opposite of findCellTopo: finds assigned cell connected to
809 // an unassigned one and puts it in the background zone.
810 void erodeCellZone
811 (
812 const label nErodeCellZones,
813 const label backgroundZoneID,
814 const labelList& unnamedSurfaceRegion,
815 const labelList& namedSurfaceRegion,
816 labelList& cellToZone
817 ) const;
818
819 //- Variation of findCellZoneTopo: walks from cellZones but ignores
820 // face intersections (unnamedSurfaceRegion). WIP
821 void growCellZone
822 (
823 const label nGrowCellZones,
824 const label backgroundZoneID,
825 labelList& unnamedSurfaceRegion1,
826 labelList& unnamedSurfaceRegion2,
827 labelList& namedSurfaceRegion,
828 labelList& cellToZone
829 ) const;
830
831 //- Make namedSurfaceRegion consistent with cellToZone
832 // - clear out any blocked faces inbetween same cell zone.
833 void makeConsistentFaceIndex
834 (
835 const labelList& zoneToNamedSurface,
836 const labelList& cellToZone,
837 labelList& namedSurfaceRegion
838 ) const;
839
840 //- Calculate cellZone allocation
841 void zonify
842 (
843 const bool allowFreeStandingZoneFaces,
844 const label nErodeCellZones,
845 const label backgroundZoneID,
846 const pointField& locationsInMesh,
847 const wordList& zonesInMesh,
848 const pointField& locationsOutsideMesh,
849 const bool exitIfLeakPath,
850 const refPtr<coordSetWriter>& leakPathFormatter,
851
852 labelList& cellToZone,
853 labelList& unnamedRegion1,
854 labelList& unnamedRegion2,
855 labelList& namedSurfaceRegion,
856 bitSet& posOrientation
857 ) const;
858
859 //- Put cells into cellZone, faces into faceZone
860 void zonify
861 (
862 const bitSet& isMasterFace,
863 const labelList& cellToZone,
864 const labelList& neiCellZone,
865 const labelList& faceToZone,
866 const bitSet& meshFlipMap,
867 polyTopoChange& meshMod
868 ) const;
869
870 //- Allocate faceZoneName
871 void allocateInterRegionFaceZone
872 (
873 const label ownZone,
874 const label neiZone,
875 wordPairHashTable& zonesToFaceZone,
876 LabelPairMap<word>& zoneIDsToFaceZone
877 ) const;
878
879 //- Remove any loose standing cells
880 void handleSnapProblems
881 (
882 const snapParameters& snapParams,
883 const bool useTopologicalSnapDetection,
884 const bool removeEdgeConnectedCells,
885 const scalarField& perpendicularAngle,
886 const dictionary& motionDict,
887 Time& runTime,
888 const labelList& globalToMasterPatch,
889 const labelList& globalToSlavePatch
890 );
891
892
893 // Some patch utilities
894
895 //- Get all faces in faceToZone that have no cellZone on
896 // either side.
897 labelList freeStandingBaffleFaces
898 (
899 const labelList& faceToZone,
900 const labelList& cellToZone,
901 const labelList& neiCellZone
902 ) const;
903
904 //- Determine per patch edge the number of master faces. Used
905 // to detect non-manifold situations.
906 void calcPatchNumMasterFaces
907 (
908 const bitSet& isMasterFace,
909 const indirectPrimitivePatch& patch,
910 labelList& nMasterFaces
911 ) const;
912
913 //- Determine per patch face the (singly-) connected zone it
914 // is in. Return overall number of zones.
915 label markPatchZones
916 (
917 const indirectPrimitivePatch& patch,
918 const labelList& nMasterFaces,
919 labelList& faceToZone
920 ) const;
921
922 //- Make faces consistent.
923 void consistentOrientation
924 (
925 const bitSet& isMasterFace,
926 const indirectPrimitivePatch& patch,
927 const labelList& nMasterFaces,
928 const labelList& faceToZone,
929 const Map<label>& zoneToOrientation,
930 bitSet& meshFlipMap
931 ) const;
932
933
934 //- No copy construct
935 meshRefinement(const meshRefinement&) = delete;
936
937 //- No copy assignment
938 void operator=(const meshRefinement&) = delete;
939
940public:
941
942 //- Runtime type information
943 ClassName("meshRefinement");
944
945
946 // Constructors
947
948 //- Construct from components
950 (
951 fvMesh& mesh,
952 const scalar mergeDistance,
953 const bool overwrite,
954 const refinementSurfaces&,
955 const refinementFeatures&,
956 const shellSurfaces&, // omnidirectional refinement
957 const shellSurfaces&, // limit refinement
958 const labelUList& checkFaces, // initial faces to check
959 const bool dryRun
960 );
961
962
963 // Member Functions
964
965 // Access
966
967 //- Reference to mesh
968 const fvMesh& mesh() const
969 {
970 return mesh_;
972 fvMesh& mesh()
973 {
974 return mesh_;
975 }
977 scalar mergeDistance() const
978 {
979 return mergeDistance_;
980 }
981
982 //- Overwrite the mesh?
983 bool overwrite() const
984 {
985 return overwrite_;
986 }
987
988 //- (points)instance of mesh upon construction
989 const word& oldInstance() const
990 {
991 return oldInstance_;
992 }
993
994 //- Reference to surface search engines
995 const refinementSurfaces& surfaces() const
996 {
997 return surfaces_;
998 }
999
1000 //- Reference to feature edge mesh
1001 const refinementFeatures& features() const
1002 {
1003 return features_;
1004 }
1005
1006 //- Reference to refinement shells (regions)
1007 const shellSurfaces& shells() const
1008 {
1009 return shells_;
1010 }
1011
1012 //- Reference to limit shells (regions)
1013 const shellSurfaces& limitShells() const
1014 {
1015 return limitShells_;
1016 }
1017
1018 //- Reference to meshcutting engine
1019 const hexRef8& meshCutter() const
1020 {
1021 return meshCutter_;
1022 }
1023
1024 //- Per start-end edge the index of the surface hit
1025 const labelList& surfaceIndex() const;
1026
1028
1029 //- For faces originating from processor faces store the original
1030 // patch
1031 const Map<label>& faceToCoupledPatch() const
1032 {
1033 return faceToCoupledPatch_;
1034 }
1035
1036 //- Additional face data that is maintained across
1037 // topo changes. Every entry is a list over all faces.
1038 // Bit of a hack. Additional flag to say whether to maintain master
1039 // only (false) or increase set to account for face-from-face.
1041 {
1042 return userFaceData_;
1043 }
1046 {
1047 return userFaceData_;
1048 }
1049
1050
1051 // Other
1052
1053 //- Count number of intersections (local)
1054 label countHits() const;
1055
1056 //- Redecompose according to cell count
1057 // keepZoneFaces : find all faceZones from zoned surfaces and keep
1058 // owner and neighbour together
1059 // keepBaffles : find all baffles and keep them together
1061 (
1062 const bool keepZoneFaces,
1063 const bool keepBaffles,
1064 const scalarField& cellWeights,
1065 decompositionMethod& decomposer,
1066 fvMeshDistribute& distributor
1067 );
1068
1069 //- Get faces with intersection.
1071
1072 //- Get points on surfaces with intersection and boundary faces.
1074
1075 //- Create patch from set of patches
1077 (
1078 const polyMesh&,
1079 const labelList&
1080 );
1081
1082 //- Helper function to make a pointVectorField with correct
1083 // bcs for mesh movement:
1084 // - adaptPatchIDs : fixedValue
1085 // - processor : calculated (so free to move)
1086 // - cyclic/wedge/symmetry : slip
1087 // - other : slip
1089 (
1090 const pointMesh& pMesh,
1091 const labelList& adaptPatchIDs
1092 );
1093
1094 //- Helper function: check that face zones are synced
1095 static void checkCoupledFaceZones(const polyMesh&);
1096
1097 //- Helper: calculate edge weights (1/length)
1098 static void calculateEdgeWeights
1099 (
1100 const polyMesh& mesh,
1101 const bitSet& isMasterEdge,
1102 const labelList& meshPoints,
1103 const edgeList& edges,
1104 scalarField& edgeWeights,
1105 scalarField& invSumWeight
1106 );
1107
1108 //- Helper: weighted sum (over all subset of mesh points) by
1109 // summing contribution from (master) edges
1110 template<class Type>
1111 static void weightedSum
1112 (
1113 const polyMesh& mesh,
1114 const bitSet& isMasterEdge,
1115 const labelList& meshPoints,
1116 const edgeList& edges,
1117 const scalarField& edgeWeights,
1118 const Field<Type>& data,
1120 );
1121
1122
1123 // Refinement
1124
1125 //- Is local topology a small gap?
1126 bool isGap
1127 (
1128 const scalar,
1129 const vector&,
1130 const vector&,
1131 const vector&,
1132 const vector&
1133 ) const;
1134
1135 //- Is local topology a small gap normal to the test vector
1136 bool isNormalGap
1137 (
1138 const scalar planarCos,
1139 const label level0,
1140 const vector& point0,
1141 const vector& normal0,
1142 const label level1,
1143 const vector& point1,
1144 const vector& normal1
1145 ) const;
1146
1147 //- Calculate list of cells to refine.
1149 (
1150 const pointField& keepPoints,
1151 const scalar curvature,
1152 const scalar planarAngle,
1153
1154 const bool featureRefinement,
1155 const bool featureDistanceRefinement,
1156 const bool internalRefinement,
1157 const bool surfaceRefinement,
1158 const bool curvatureRefinement,
1159 const bool smallFeatureRefinement,
1160 const bool gapRefinement,
1161 const bool bigGapRefinement,
1162 const bool spreadGapSize,
1163 const label maxGlobalCells,
1164 const label maxLocalCells
1165 ) const;
1166
1167
1168 // Blocking cells
1169
1170 //- Mark faces on interface between set and rest
1171 // (and same cell level)
1172 void markOutsideFaces
1173 (
1174 const labelList& cellLevel,
1175 const labelList& neiLevel,
1176 const labelList& refineCell,
1177 bitSet& isOutsideFace
1178 ) const;
1179
1180 //- Count number of faces on cell that are in set
1181 label countFaceDirs
1182 (
1183 const bitSet& isOutsideFace,
1184 const label celli
1185 ) const;
1186
1187 //- Add one layer of cells to set
1188 void growSet
1189 (
1190 const labelList& neiLevel,
1191 const bitSet& isOutsideFace,
1193 label& nRefine
1194 ) const;
1195
1196 //- Detect gapRefinement cells and remove them
1198 (
1199 const scalar planarAngle,
1200 const labelList& minSurfaceLevel,
1201 const labelList& globalToMasterPatch,
1202 const label growIter
1203 );
1204
1205 // Blocking leaks (by blocking cells)
1206
1207 //- Faces currently on boundary or intersected by surface
1209 (
1210 const labelList& surfaces,
1211 boolList& isBlockedFace
1212 ) const;
1213
1214 //- Return list of cells to block by walking from the seedCells
1215 // until reaching a leak face
1217 (
1218 const boolList& isBlockedFace,
1219 const labelList& leakFaces,
1220 const labelList& seedCells
1221 ) const;
1222
1223 //- Remove minimum amount of cells to break any leak from
1224 // inside to outside
1226 (
1227 const labelList& globalToMasterPatch,
1228 const labelList& globalToSlavePatch,
1229 const pointField& locationsInMesh,
1230 const pointField& locationsOutsideMesh,
1231 const labelList& selectedSurfaces
1232 );
1233
1234 //- Baffle faces to break any leak from inside to outside
1236 (
1237 const labelList& globalToMasterPatch,
1238 const labelList& globalToSlavePatch,
1239 const pointField& locationsInMesh,
1240 const wordList& zonesInMesh,
1241 const pointField& locationsOutsideMesh,
1242 const labelList& selectedSurfaces
1243 );
1244
1245
1246 //- Refine some cells
1247 autoPtr<mapPolyMesh> refine(const labelList& cellsToRefine);
1248
1249 //- Refine some cells and rebalance
1251 (
1252 const string& msg,
1253 decompositionMethod& decomposer,
1254 fvMeshDistribute& distributor,
1255 const labelList& cellsToRefine,
1256 const scalar maxLoadUnbalance
1257 );
1258
1259 //- Balance before refining some cells
1261 (
1262 const string& msg,
1263 decompositionMethod& decomposer,
1264 fvMeshDistribute& distributor,
1265 const labelList& cellsToRefine,
1266 const scalar maxLoadUnbalance
1267 );
1268
1269 //- Calculate list of cells to directionally refine
1271 (
1272 const label maxGlobalCells,
1273 const label maxLocalCells,
1274 const labelList& currentLevel,
1275 const direction dir
1276 ) const;
1277
1278 //- Directionally refine in direction cmpt
1280 (
1281 const string& msg,
1282 const direction cmpt,
1283 const labelList& cellsToRefine
1284 );
1285
1286
1287 // Baffle handling
1288
1289 //- Split off unreachable areas of mesh.
1291 (
1292 const bool handleSnapProblems,
1293
1294 // How to remove problem snaps
1295 const snapParameters& snapParams,
1296 const bool useTopologicalSnapDetection,
1297 const bool removeEdgeConnectedCells,
1298 const scalarField& perpendicularAngle,
1299 const label nErodeCellZones,
1300 const dictionary& motionDict,
1301 Time& runTime,
1302 const labelList& globalToMasterPatch,
1303 const labelList& globalToSlavePatch,
1304 const pointField& locationsInMesh,
1305 const wordList& regionsInMesh,
1306 const pointField& locationsOutsideMesh,
1307 const bool exitIfLeakPath,
1308 const refPtr<coordSetWriter>& leakPathFormatter
1309 );
1310
1311 //- Merge free-standing baffles
1313 (
1314 const snapParameters& snapParams,
1315 const bool useTopologicalSnapDetection,
1316 const bool removeEdgeConnectedCells,
1317 const scalarField& perpendicularAngle,
1318 const scalar planarAngle,
1319 const dictionary& motionDict,
1320 Time& runTime,
1321 const labelList& globalToMasterPatch,
1322 const labelList& globalToSlavePatch,
1323 const pointField& locationsInMesh,
1324 const pointField& locationsOutsideMesh
1325 );
1326
1327 //- Split off (with optional buffer layers) unreachable areas
1328 // of mesh. Does not introduce baffles.
1329 autoPtr<mapPolyMesh> splitMesh
1330 (
1331 const label nBufferLayers,
1332 const label nErodeCellZones,
1333 const labelList& globalToMasterPatch,
1334 const labelList& globalToSlavePatch,
1335
1336 const pointField& locationsInMesh,
1337 const wordList& regionsInMesh,
1338 const pointField& locationsOutsideMesh,
1339 const bool exitIfLeakPath,
1340 const refPtr<coordSetWriter>& leakPathFormatter
1341 );
1342
1343 //- Remove cells from limitRegions if level -1
1345 (
1346 const label nBufferLayers,
1347 const label nErodeCellZones,
1348 const labelList& globalToMasterPatch,
1349 const labelList& globalToSlavePatch,
1350 const pointField& locationsInMesh,
1351 const wordList& regionsInMesh,
1352 const pointField& locationsOutsideMesh
1353 );
1354
1355 //- Find boundary points that connect to more than one cell
1356 // region and split them.
1358
1359 //- Find boundary points that connect to more than one cell
1360 // region and split them.
1362
1363 //- Find boundary points that are on faceZones of type boundary
1364 // and duplicate them
1366
1367 //- Merge duplicate points
1369 (
1370 const labelList& pointToDuplicate
1371 );
1372
1373 //- Create baffle for every internal face where ownPatch != -1.
1374 // External faces get repatched according to ownPatch (neiPatch
1375 // should be -1 for these)
1377 (
1378 const labelList& ownPatch,
1379 const labelList& neiPatch
1380 );
1381
1382 //- Get zones of given type
1384 (
1386 ) const;
1387
1388 //- Subset baffles according to zones
1390 (
1391 const polyMesh& mesh,
1392 const labelList& zoneIDs,
1393 const List<labelPair>& baffles
1394 );
1395
1396 //- Get per-face information (faceZone, master/slave patch)
1397 void getZoneFaces
1398 (
1399 const labelList& zoneIDs,
1401 labelList& ownPatch,
1402 labelList& neiPatch,
1403 labelList& nBaffles
1404 ) const;
1405
1406 //- Create baffles for faces on faceZones. Return created baffles
1407 // (= pairs of faces) and corresponding faceZone
1409 (
1410 const labelList& zoneIDs,
1411 List<labelPair>& baffles,
1412 labelList& originatingFaceZone
1413 );
1414
1415 //- Merge baffles. Gets pairs of faces and boundary faces to move
1416 // onto (coupled) patches
1418 (
1419 const List<labelPair>&,
1420 const Map<label>& faceToPatch
1421 );
1422
1423 //- Merge all baffles on faceZones
1425 (
1426 const bool doInternalZones,
1427 const bool doBaffleZones
1428 );
1429
1430 //- Put faces/cells into zones according to surface specification.
1431 // Returns null if no zone surfaces present. Regions containing
1432 // locationsInMesh/regionsInMesh will be put in corresponding
1433 // cellZone. keepPoints is for backwards compatibility and sets
1434 // all yet unassigned cells to be non-zoned (zone = -1)
1436 (
1437 const bool allowFreeStandingZoneFaces,
1438 const label nErodeCellZones,
1439 const pointField& locationsInMesh,
1440 const wordList& regionsInMesh,
1441 const pointField& locationsOutsideMesh,
1442 const bool exitIfLeakPath,
1443 const refPtr<coordSetWriter>& leakPathFormatter,
1444 wordPairHashTable& zonesToFaceZone
1445 );
1446
1447
1448 // Other topo changes
1449
1450 //- Helper:append patch to end of mesh.
1451 static label appendPatch
1452 (
1453 fvMesh&,
1454 const label insertPatchi,
1455 const word&,
1456 const dictionary&
1457 );
1458
1459 //- Helper:add patch to mesh. Update all registered fields.
1460 // Used by addMeshedPatch to add patches originating from surfaces.
1461 static label addPatch(fvMesh&, const word& name, const dictionary&);
1462
1463 //- Add patch originating from meshing. Update meshedPatches_.
1464 label addMeshedPatch(const word& name, const dictionary&);
1465
1466 //- Get patchIDs for patches added in addMeshedPatch.
1467 labelList meshedPatches() const;
1468
1469 //- Add/lookup faceZone and update information. Return index of
1470 // faceZone
1471 label addFaceZone
1472 (
1473 const word& fzName,
1474 const word& masterPatch,
1475 const word& slavePatch,
1476 const surfaceZonesInfo::faceZoneType& fzType
1477 );
1478
1479 //- Lookup faceZone information. Return false if no information
1480 // for faceZone
1481 bool getFaceZoneInfo
1482 (
1483 const word& fzName,
1484 label& masterPatchID,
1485 label& slavePatchID,
1487 ) const;
1488
1489 //- Add pointZone if does not exist. Return index of zone
1490 label addPointZone(const word& name);
1491
1492 //- Count number of faces per patch edge. Parallel consistent.
1494
1495 //- Select coupled faces that are not collocated
1497
1498 //- Find any intersection of surface. Store in surfaceIndex_.
1499 void updateIntersections(const labelList& changedFaces);
1500
1501 //- Calculate nearest intersection for selected mesh faces
1502 void nearestIntersection
1503 (
1504 const labelList& surfacesToTest,
1505 const labelList& testFaces,
1506
1507 labelList& surface1,
1508 List<pointIndexHit>& hit1,
1509 labelList& region1,
1510 labelList& surface2,
1511 List<pointIndexHit>& hit2,
1512 labelList& region2
1513 ) const;
1514
1515 //- Remove cells. Put exposedFaces into exposedPatchIDs.
1517 (
1518 const labelList& cellsToRemove,
1519 const labelList& exposedFaces,
1520 const labelList& exposedPatchIDs,
1521 removeCells& cellRemover
1522 );
1523
1524 //- Find cell point is in. Uses optional perturbation to re-test.
1525 // Returns -1 on processors that do not have the cell.
1526 static label findCell
1527 (
1528 const polyMesh&,
1529 const vector& perturbVec,
1530 const point& p
1531 );
1532
1533 //- Find region point is in. Uses optional perturbation to re-test.
1534 static label findRegion
1535 (
1536 const polyMesh&,
1537 const labelList& cellRegion,
1538 const vector& perturbVec,
1539 const point& p
1540 );
1541
1542 //- Find regions points are in.
1543 // \return number of cells to be removed
1544 static label findRegions
1545 (
1546 const polyMesh&,
1547 const vector& perturbVec,
1548 const pointField& locationsInMesh,
1549 const pointField& locationsOutsideMesh,
1550 const label nRegions,
1551 labelList& cellRegion,
1552 const boolList& blockedFace,
1553 // Leak-path
1554 const bool exitIfLeakPath,
1555 const refPtr<coordSetWriter>& leakPathFormatter
1556 );
1557
1558 //- Split mesh. Keep part containing point. Return empty map if
1559 // no cells removed.
1561 (
1562 const labelList& globalToMasterPatch,
1563 const labelList& globalToSlavePatch,
1564 const pointField& locationsInMesh,
1565 const pointField& locationsOutsideMesh,
1566 // Leak-path
1567 const bool exitIfLeakPath,
1568 const refPtr<coordSetWriter>& leakPathFormatter
1569 );
1570
1571 //- Split faces into two
1572 void doSplitFaces
1573 (
1574 const labelList& splitFaces,
1575 const labelPairList& splits,
1576 polyTopoChange& meshMod
1577 ) const;
1578
1579 //- Split faces along diagonal. Maintain mesh quality. Return
1580 // total number of faces split.
1581 label splitFacesUndo
1582 (
1583 const labelList& splitFaces,
1584 const labelPairList& splits,
1585 const dictionary& motionDict,
1586
1587 labelList& duplicateFace,
1588 List<labelPair>& baffles
1589 );
1590
1591 //- Update local numbering for mesh redistribution
1592 void distribute(const mapDistributePolyMesh&);
1593
1594 //- Update for external change to mesh. changedFaces are in new mesh
1595 // face labels.
1596 void updateMesh
1597 (
1598 const mapPolyMesh&,
1599 const labelList& changedFaces
1600 );
1601
1602 //- Helper: reorder list according to map.
1603 template<class T>
1604 static void updateList
1605 (
1606 const labelList& newToOld,
1607 const T& nullValue,
1608 List<T>& elems
1609 );
1610
1611
1612 // Restoring : is where other processes delete and reinsert data.
1613
1614 //- Signal points/face/cells for which to store data
1615 void storeData
1616 (
1617 const labelList& pointsToStore,
1618 const labelList& facesToStore,
1619 const labelList& cellsToStore
1620 );
1621
1622 //- Update local numbering + undo
1623 // Data to restore given as new pointlabel + stored pointlabel
1624 // (i.e. what was in pointsToStore)
1625 void updateMesh
1626 (
1627 const mapPolyMesh&,
1628 const labelList& changedFaces,
1629 const Map<label>& pointsToRestore,
1630 const Map<label>& facesToRestore,
1631 const Map<label>& cellsToRestore
1632 );
1633
1634 // Merging coplanar faces and edges
1635
1636 //- Merge coplanar faces if sets are of size mergeSize
1637 // (usually 4)
1638 label mergePatchFaces
1639 (
1640 const scalar minCos,
1641 const scalar concaveCos,
1642 const label mergeSize,
1643 const labelList& patchIDs,
1644 const meshRefinement::FaceMergeType mergeType
1645 );
1646
1647 //- Merge coplanar faces. preserveFaces is != -1 for faces
1648 // to be preserved
1650 (
1651 const scalar minCos,
1652 const scalar concaveCos,
1653 const labelList& patchIDs,
1654 const dictionary& motionDict,
1655 const labelList& preserveFaces,
1656 const meshRefinement::FaceMergeType mergeType
1657 );
1658
1660 (
1661 removePoints& pointRemover,
1662 const boolList& pointCanBeDeleted
1663 );
1664
1666 (
1667 removePoints& pointRemover,
1668 const labelList& facesToRestore
1669 );
1670
1672 (
1673 const labelList& candidateFaces,
1674 const labelHashSet& set
1675 ) const;
1676
1677 // Pick up faces of cells of faces in set.
1679 (
1680 const labelUList& set
1681 ) const;
1682
1683 // Pick up faces of cells of faces in set.
1685 (
1686 const labelHashSet& set
1687 ) const;
1688
1689 //- Merge edges, maintain mesh quality. Return global number
1690 // of edges merged
1691 label mergeEdgesUndo
1692 (
1693 const scalar minCos,
1694 const dictionary& motionDict
1695 );
1696
1697
1698 // Debug/IO
1699
1700 //- Debugging: check that all faces still obey start()>end()
1701 void checkData();
1702
1703 static void testSyncPointList
1704 (
1705 const string& msg,
1706 const polyMesh& mesh,
1707 const List<scalar>& fld
1708 );
1709
1710 static void testSyncPointList
1711 (
1712 const string& msg,
1713 const polyMesh& mesh,
1714 const List<point>& fld
1715 );
1716
1717 //- Compare two lists over all boundary faces
1718 template<class T>
1720 (
1721 const scalar mergeDistance,
1722 const string&,
1723 const UList<T>&,
1724 const UList<T>&
1725 ) const;
1726
1727 //- Print list according to (collected and) sorted coordinate
1728 template<class T>
1729 static void collectAndPrint
1730 (
1731 const UList<point>& points,
1732 const UList<T>& data
1733 );
1734
1735 //- Determine master point for subset of points. If coupled
1736 // chooses only one
1737 static bitSet getMasterPoints
1738 (
1739 const polyMesh& mesh,
1740 const labelList& meshPoints
1741 );
1742
1743 //- Determine master edge for subset of edges. If coupled
1744 // chooses only one
1745 static bitSet getMasterEdges
1746 (
1747 const polyMesh& mesh,
1748 const labelList& meshEdges
1749 );
1750
1751 //- Print some mesh stats.
1752 void printMeshInfo(const bool, const string&) const;
1753
1754 //- Replacement for Time::timeName() that returns oldInstance
1755 //- (if overwrite_)
1756 word timeName() const;
1757
1758 //- Set instance of all local IOobjects
1759 void setInstance(const fileName&);
1760
1761 //- Write mesh and all data
1762 bool write() const;
1763
1764 //- Write refinement level as volScalarFields for postprocessing
1765 void dumpRefinementLevel() const;
1766
1767 //- Debug: Write intersection information to OBJ format
1768 void dumpIntersections(const fileName& prefix) const;
1769
1770 //- Do any one of above IO functions
1771 void write
1772 (
1773 const debugType debugFlags,
1774 const writeType writeFlags,
1775 const fileName&
1776 ) const;
1777
1778 //- Helper: remove all relevant files from mesh instance
1779 static void removeFiles(const polyMesh&);
1780
1781 //- Helper: calculate average
1782 template<class T>
1783 static T gAverage
1784 (
1785 const bitSet& isMasterElem,
1786 const UList<T>& values
1787 );
1788
1789 //- Get/set write level
1790 static writeType writeLevel();
1791 static void writeLevel(const writeType);
1792
1794 //static outputType outputLevel();
1795 //static void outputLevel(const outputType);
1796
1797
1798 //- Helper: convert wordList into bit pattern using provided Enum
1799 template<class EnumContainer>
1800 static int readFlags
1801 (
1802 const EnumContainer& namedEnum,
1803 const wordList& words
1804 );
1805
1806 //- Wrapper around dictionary::get which does not exit
1807 template<class Type>
1808 static Type get
1809 (
1810 const dictionary& dict,
1811 const word& keyword,
1812 const bool noExit,
1813 enum keyType::option matchOpt = keyType::REGEX,
1814 const Type& deflt = Zero
1815 );
1816
1817 //- Wrapper around dictionary::subDict which does not exit
1818 static const dictionary& subDict
1819 (
1820 const dictionary& dict,
1821 const word& keyword,
1822 const bool noExit,
1823 enum keyType::option matchOpt = keyType::REGEX
1824 );
1825
1826 //- Wrapper around dictionary::lookup which does not exit
1827 static ITstream& lookup
1828 (
1829 const dictionary& dict,
1830 const word& keyword,
1831 const bool noExit,
1832 enum keyType::option matchOpt = keyType::REGEX
1833 );
1834};
1835
1836
1837// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
1838
1839} // End namespace Foam
1840
1841// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
1842
1843#ifdef NoRepository
1844 #include "meshRefinementTemplates.C"
1845#endif
1846
1847// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
1848
1849#endif
1850
1851// ************************************************************************* //
Info<< nl<< "Wrote faMesh in vtk format: "<< writer.output().name()<< nl;}{ vtk::lineWriter writer(aMesh.points(), aMesh.edges(), fileName(aMesh.mesh().time().globalPath()/"finiteArea-edges"));writer.writeGeometry();writer.beginCellData(4);writer.writeProcIDs();{ Field< scalar > fld(faMeshTools::flattenEdgeField(aMesh.magLe(), true))
Dynamically sized Field.
Definition: DynamicField.H:65
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
Definition: Enum.H:61
A 1D vector of objects of type <T> with a fixed length <N>.
Definition: FixedList.H:81
A HashTable similar to std::unordered_map.
Definition: HashTable.H:123
An input stream of tokens.
Definition: ITstream.H:56
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: List.H:77
A HashTable to objects of type <T> with a label key.
Definition: Map.H:60
A list of faces which address into the list of points.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:80
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
Definition: autoPtr.H:66
A bitSet stores bits (elements with only two states) in packed internal format and supports a variety...
Definition: bitSet.H:66
Base class for writing coordSet(s) and tracks with fields.
Database for solution data, solver performance and other reduced data.
Definition: data.H:58
Abstract base class for domain decomposition.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
A subset of mesh faces organised as a primitive patch.
Definition: faceZone.H:67
A class for handling file names.
Definition: fileName.H:76
Sends/receives parts of mesh+fvfields to neighbouring processors. Used in load balancing.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:91
Refinement of (split) hexes using polyTopoChange.
Definition: hexRef8.H:68
option
Enumeration for the data type and search/match modes (bitmask)
Definition: keyType.H:79
@ REGEX
Regular expression.
Definition: keyType.H:82
Takes mesh with 'baffles' (= boundary faces sharing points). Determines for selected points on bounda...
Class containing mesh-to-mesh mapping information after a mesh distribution where we send parts of me...
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Definition: mapPolyMesh.H:162
Helper class which maintains intersections of (changing) mesh with (static) surfaces.
autoPtr< mapDistributePolyMesh > balanceAndRefine(const string &msg, decompositionMethod &decomposer, fvMeshDistribute &distributor, const labelList &cellsToRefine, const scalar maxLoadUnbalance)
Balance before refining some cells.
static void calculateEdgeWeights(const polyMesh &mesh, const bitSet &isMasterEdge, const labelList &meshPoints, const edgeList &edges, scalarField &edgeWeights, scalarField &invSumWeight)
Helper: calculate edge weights (1/length)
label addFaceZone(const word &fzName, const word &masterPatch, const word &slavePatch, const surfaceZonesInfo::faceZoneType &fzType)
Add/lookup faceZone and update information. Return index of.
label splitFacesUndo(const labelList &splitFaces, const labelPairList &splits, const dictionary &motionDict, labelList &duplicateFace, List< labelPair > &baffles)
Split faces along diagonal. Maintain mesh quality. Return.
labelList getZones(const List< surfaceZonesInfo::faceZoneType > &fzTypes) const
Get zones of given type.
bool getFaceZoneInfo(const word &fzName, label &masterPatchID, label &slavePatchID, surfaceZonesInfo::faceZoneType &fzType) const
Lookup faceZone information. Return false if no information.
static T gAverage(const bitSet &isMasterElem, const UList< T > &values)
Helper: calculate average.
static void collectAndPrint(const UList< point > &points, const UList< T > &data)
Print list according to (collected and) sorted coordinate.
labelList intersectedFaces() const
Get faces with intersection.
labelList intersectedPoints() const
Get points on surfaces with intersection and boundary faces.
autoPtr< mapPolyMesh > removeLimitShells(const label nBufferLayers, const label nErodeCellZones, const labelList &globalToMasterPatch, const labelList &globalToSlavePatch, const pointField &locationsInMesh, const wordList &regionsInMesh, const pointField &locationsOutsideMesh)
Remove cells from limitRegions if level -1.
static void removeFiles(const polyMesh &)
Helper: remove all relevant files from mesh instance.
label mergeEdgesUndo(const scalar minCos, const dictionary &motionDict)
Merge edges, maintain mesh quality. Return global number.
static bitSet getMasterPoints(const polyMesh &mesh, const labelList &meshPoints)
Determine master point for subset of points. If coupled.
static void testSyncPointList(const string &msg, const polyMesh &mesh, const List< scalar > &fld)
void checkData()
Debugging: check that all faces still obey start()>end()
autoPtr< mapPolyMesh > blockLeakFaces(const labelList &globalToMasterPatch, const labelList &globalToSlavePatch, const pointField &locationsInMesh, const wordList &zonesInMesh, const pointField &locationsOutsideMesh, const labelList &selectedSurfaces)
Baffle faces to break any leak from inside to outside.
void updateIntersections(const labelList &changedFaces)
Find any intersection of surface. Store in surfaceIndex_.
static List< labelPair > subsetBaffles(const polyMesh &mesh, const labelList &zoneIDs, const List< labelPair > &baffles)
Subset baffles according to zones.
const List< Tuple2< mapType, labelList > > & userFaceData() const
Additional face data that is maintained across.
static Type get(const dictionary &dict, const word &keyword, const bool noExit, enum keyType::option matchOpt=keyType::REGEX, const Type &deflt=Zero)
Wrapper around dictionary::get which does not exit.
void updateMesh(const mapPolyMesh &, const labelList &changedFaces)
Update for external change to mesh. changedFaces are in new mesh.
writeType
Enumeration for what to write. Used as a bit-pattern.
const hexRef8 & meshCutter() const
Reference to meshcutting engine.
autoPtr< mapPolyMesh > mergePoints(const labelList &pointToDuplicate)
Merge duplicate points.
const labelList & surfaceIndex() const
Per start-end edge the index of the surface hit.
autoPtr< mapPolyMesh > doRestorePoints(removePoints &pointRemover, const labelList &facesToRestore)
autoPtr< mapPolyMesh > dupNonManifoldPoints()
Find boundary points that connect to more than one cell.
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 smallFeatureRefinement, const bool gapRefinement, const bool bigGapRefinement, const bool spreadGapSize, const label maxGlobalCells, const label maxLocalCells) const
Calculate list of cells to refine.
const shellSurfaces & limitShells() const
Reference to limit shells (regions)
label countFaceDirs(const bitSet &isOutsideFace, const label celli) const
Count number of faces on cell that are in set.
labelList growFaceCellFace(const labelUList &set) const
void printMeshInfo(const bool, const string &) const
Print some mesh stats.
void dumpRefinementLevel() const
Write refinement level as volScalarFields for postprocessing.
autoPtr< mapPolyMesh > removeGapCells(const scalar planarAngle, const labelList &minSurfaceLevel, const labelList &globalToMasterPatch, const label growIter)
Detect gapRefinement cells and remove them.
autoPtr< mapDistributePolyMesh > balance(const bool keepZoneFaces, const bool keepBaffles, const scalarField &cellWeights, decompositionMethod &decomposer, fvMeshDistribute &distributor)
Redecompose according to cell count.
autoPtr< mapPolyMesh > dupNonManifoldBoundaryPoints()
Find boundary points that are on faceZones of type boundary.
static const Enum< writeType > writeTypeNames
static int readFlags(const EnumContainer &namedEnum, const wordList &words)
Helper: convert wordList into bit pattern using provided Enum.
void distribute(const mapDistributePolyMesh &)
Update local numbering for mesh redistribution.
void getZoneFaces(const labelList &zoneIDs, labelList &faceZoneID, labelList &ownPatch, labelList &neiPatch, labelList &nBaffles) const
Get per-face information (faceZone, master/slave patch)
const Map< label > & faceToCoupledPatch() const
For faces originating from processor faces store the original.
autoPtr< mapPolyMesh > mergeBaffles(const List< labelPair > &, const Map< label > &faceToPatch)
Merge baffles. Gets pairs of faces and boundary faces to move.
static label addPatch(fvMesh &, const word &name, const dictionary &)
Helper:add patch to mesh. Update all registered fields.
void baffleAndSplitMesh(const bool handleSnapProblems, const snapParameters &snapParams, const bool useTopologicalSnapDetection, const bool removeEdgeConnectedCells, const scalarField &perpendicularAngle, const label nErodeCellZones, const dictionary &motionDict, Time &runTime, const labelList &globalToMasterPatch, const labelList &globalToSlavePatch, const pointField &locationsInMesh, const wordList &regionsInMesh, const pointField &locationsOutsideMesh, const bool exitIfLeakPath, const refPtr< coordSetWriter > &leakPathFormatter)
Split off unreachable areas of mesh.
label countHits() const
Count number of intersections (local)
void markOutsideFaces(const labelList &cellLevel, const labelList &neiLevel, const labelList &refineCell, bitSet &isOutsideFace) const
Mark faces on interface between set and rest.
ClassName("meshRefinement")
Runtime type information.
autoPtr< mapPolyMesh > doRemoveCells(const labelList &cellsToRemove, const labelList &exposedFaces, const labelList &exposedPatchIDs, removeCells &cellRemover)
Remove cells. Put exposedFaces into exposedPatchIDs.
mapType
Enumeration for how userdata is to be mapped upon refinement.
@ KEEPALL
have slaves (upon refinement) from master
@ REMOVE
set value to -1 any face that was refined
@ MASTERONLY
maintain master only
static const dictionary & subDict(const dictionary &dict, const word &keyword, const bool noExit, enum keyType::option matchOpt=keyType::REGEX)
Wrapper around dictionary::subDict which does not exit.
static label findCell(const polyMesh &, const vector &perturbVec, const point &p)
Find cell point is in. Uses optional perturbation to re-test.
const refinementSurfaces & surfaces() const
Reference to surface search engines.
label mergePatchFaces(const scalar minCos, const scalar concaveCos, const label mergeSize, const labelList &patchIDs, const meshRefinement::FaceMergeType mergeType)
Merge coplanar faces if sets are of size mergeSize.
debugType
Enumeration for what to debug. Used as a bit-pattern.
void doSplitFaces(const labelList &splitFaces, const labelPairList &splits, polyTopoChange &meshMod) const
Split faces into two.
static bitSet getMasterEdges(const polyMesh &mesh, const labelList &meshEdges)
Determine master edge for subset of edges. If coupled.
const shellSurfaces & shells() const
Reference to refinement shells (regions)
labelList meshedPatches() const
Get patchIDs for patches added in addMeshedPatch.
static void checkCoupledFaceZones(const polyMesh &)
Helper function: check that face zones are synced.
labelList countEdgeFaces(const uindirectPrimitivePatch &pp) const
Count number of faces per patch edge. Parallel consistent.
bool isNormalGap(const scalar planarCos, const label level0, const vector &point0, const vector &normal0, const label level1, const vector &point1, const vector &normal1) const
Is local topology a small gap normal to the test vector.
labelList collectFaces(const labelList &candidateFaces, const labelHashSet &set) const
static void weightedSum(const polyMesh &mesh, const bitSet &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.
autoPtr< mapPolyMesh > splitMeshRegions(const labelList &globalToMasterPatch, const labelList &globalToSlavePatch, const pointField &locationsInMesh, const pointField &locationsOutsideMesh, const bool exitIfLeakPath, const refPtr< coordSetWriter > &leakPathFormatter)
Split mesh. Keep part containing point. Return empty map if.
label addMeshedPatch(const word &name, const dictionary &)
Add patch originating from meshing. Update meshedPatches_.
List< Tuple2< mapType, labelList > > & userFaceData()
const fvMesh & mesh() const
Reference to mesh.
bool isGap(const scalar, const vector &, const vector &, const vector &, const vector &) const
Is local topology a small gap?
const refinementFeatures & features() const
Reference to feature edge mesh.
void testSyncBoundaryFaceList(const scalar mergeDistance, const string &, const UList< T > &, const UList< T > &) const
Compare two lists over all boundary faces.
void selectSeparatedCoupledFaces(boolList &) const
Select coupled faces that are not collocated.
static tmp< pointVectorField > makeDisplacementField(const pointMesh &pMesh, const labelList &adaptPatchIDs)
Helper function to make a pointVectorField with correct.
const word & oldInstance() const
(points)instance of mesh upon construction
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.
label addPointZone(const word &name)
Add pointZone if does not exist. Return index of zone.
void selectIntersectedFaces(const labelList &surfaces, boolList &isBlockedFace) const
Faces currently on boundary or intersected by surface.
static label findRegions(const polyMesh &, const vector &perturbVec, const pointField &locationsInMesh, const pointField &locationsOutsideMesh, const label nRegions, labelList &cellRegion, const boolList &blockedFace, const bool exitIfLeakPath, const refPtr< coordSetWriter > &leakPathFormatter)
Find regions points are in.
autoPtr< mapPolyMesh > doRemovePoints(removePoints &pointRemover, const boolList &pointCanBeDeleted)
scalar mergeDistance() const
FaceMergeType
Enumeration for what to do with co-planar patch faces on a single.
autoPtr< mapDistributePolyMesh > refineAndBalance(const string &msg, decompositionMethod &decomposer, fvMeshDistribute &distributor, const labelList &cellsToRefine, const scalar maxLoadUnbalance)
Refine some cells and rebalance.
autoPtr< mapPolyMesh > directionalRefine(const string &msg, const direction cmpt, const labelList &cellsToRefine)
Directionally refine in direction cmpt.
autoPtr< mapPolyMesh > refine(const labelList &cellsToRefine)
Refine some cells.
static autoPtr< indirectPrimitivePatch > makePatch(const polyMesh &, const labelList &)
Create patch from set of patches.
void storeData(const labelList &pointsToStore, const labelList &facesToStore, const labelList &cellsToStore)
Signal points/face/cells for which to store data.
static void updateList(const labelList &newToOld, const T &nullValue, List< T > &elems)
Helper: reorder list according to map.
bool overwrite() const
Overwrite the mesh?
bool write() const
Write mesh and all data.
autoPtr< mapPolyMesh > mergeZoneBaffles(const bool doInternalZones, const bool doBaffleZones)
Merge all baffles on faceZones.
labelList directionalRefineCandidates(const label maxGlobalCells, const label maxLocalCells, const labelList &currentLevel, const direction dir) const
Calculate list of cells to directionally refine.
void growSet(const labelList &neiLevel, const bitSet &isOutsideFace, labelList &refineCell, label &nRefine) const
Add one layer of cells to set.
label mergePatchFacesUndo(const scalar minCos, const scalar concaveCos, const labelList &patchIDs, const dictionary &motionDict, const labelList &preserveFaces, const meshRefinement::FaceMergeType mergeType)
Merge coplanar faces. preserveFaces is != -1 for faces.
autoPtr< mapPolyMesh > removeLeakCells(const labelList &globalToMasterPatch, const labelList &globalToSlavePatch, const pointField &locationsInMesh, const pointField &locationsOutsideMesh, const labelList &selectedSurfaces)
Remove minimum amount of cells to break any leak from.
void setInstance(const fileName &)
Set instance of all local IOobjects.
autoPtr< mapPolyMesh > createBaffles(const labelList &ownPatch, const labelList &neiPatch)
Create baffle for every internal face where ownPatch != -1.
static const Enum< debugType > debugTypeNames
static label appendPatch(fvMesh &, const label insertPatchi, const word &, const dictionary &)
Helper:append patch to end of mesh.
static writeType writeLevel()
Get/set write level.
autoPtr< mapPolyMesh > createZoneBaffles(const labelList &zoneIDs, List< labelPair > &baffles, labelList &originatingFaceZone)
Create baffles for faces on faceZones. Return created baffles.
void dumpIntersections(const fileName &prefix) const
Debug: Write intersection information to OBJ format.
labelList detectLeakCells(const boolList &isBlockedFace, const labelList &leakFaces, const labelList &seedCells) const
Return list of cells to block by walking from the seedCells.
void mergeFreeStandingBaffles(const snapParameters &snapParams, const bool useTopologicalSnapDetection, const bool removeEdgeConnectedCells, const scalarField &perpendicularAngle, const scalar planarAngle, const dictionary &motionDict, Time &runTime, const labelList &globalToMasterPatch, const labelList &globalToSlavePatch, const pointField &locationsInMesh, const pointField &locationsOutsideMesh)
Merge free-standing baffles.
Mesh representing a set of points created from polyMesh.
Definition: pointMesh.H:55
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:81
Direct mesh changes based on v1.3 polyTopoChange syntax.
Lookup type of boundary radiation properties.
Definition: lookup.H:66
A class for managing references or pointers (no reference counting)
Definition: refPtr.H:58
Container with cells to refine. Refinement given as single direction.
Definition: refineCell.H:57
Encapsulates queries for features.
Container for data on surfaces used for surface-driven refinement. Contains all the data about the le...
Given list of cells to remove, insert all the topology changes.
Definition: removeCells.H:64
Removes selected points from mesh and updates faces using these points.
Definition: removePoints.H:61
Encapsulates queries for volume refinement ('refine all cells within shell').
Definition: shellSurfaces.H:58
Simple container to keep together snap specific information.
Contains information about location on a triSurface.
faceZoneType
What to do with faceZone faces.
A class for managing temporary objects.
Definition: tmp.H:65
An enumeration wrapper for classification of a location as being inside/outside of a volume.
Definition: volumeType.H:61
A class for handling words, derived from Foam::string.
Definition: word.H:68
#define ClassName(TypeNameString)
Add typeName information from argument TypeNameString to a class.
Definition: className.H:67
volScalarField & p
const volScalarField & T
engineTime & runTime
const pointField & points
const labelIOList & zoneIDs
Definition: correctPhi.H:59
A HashTable to objects of type <T> with a labelPair key. The hashing is based on labelPair (FixedList...
Namespace for OpenFOAM.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
mode_t mode(const fileName &name, const bool followLink=true)
Return the file mode, normally following symbolic links.
Definition: MSwindows.C:572
constexpr label labelMax
Definition: label.H:61
uint8_t direction
Definition: direction.H:56
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
Forwards and collection of common point field types.
dictionary dict
volScalarField & e
Definition: createFields.H:11