backgroundMeshDecomposition.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-2016 OpenFOAM Foundation
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26Class
27 Foam::backgroundMeshDecomposition
28
29Description
30 Store a background polyMesh to use for the decomposition of space and
31 queries for parallel conformalVoronoiMesh.
32
33 The requirements are:
34
35 - To have a decomposition of space which can quickly interrogate an
36 arbitrary location from any processor to reliably and unambiguously
37 determine which processor owns the space that the point is in, i.e. as
38 the vertices move, or need inserted as part of the surface conformation,
39 send them to the correct proc.
40
41 - To be able to be dynamically built, refined and redistributed to other
42 procs the partitioning as the meshing progresses to balance the load.
43
44 - To be able to query whether a sphere (the circumsphere of a Delaunay tet)
45 overlaps any part of the space defined by the structure, and whether a
46 ray (Voronoi edge) penetrates any part of the space defined by the
47 structure, this is what determines if points get referred to a processor.
48
49SourceFiles
50 backgroundMeshDecompositionI.H
51 backgroundMeshDecomposition.C
52
53\*---------------------------------------------------------------------------*/
54
55#ifndef backgroundMeshDecomposition_H
56#define backgroundMeshDecomposition_H
57
58#include "fvMesh.H"
59#include "hexRef8.H"
60#include "cellSet.H"
61#include "meshTools.H"
62#include "polyTopoChange.H"
63#include "mapPolyMesh.H"
64#include "decompositionMethod.H"
65#include "fvMeshDistribute.H"
66#include "removeCells.H"
68#include "globalIndex.H"
69#include "treeBoundBox.H"
70#include "primitivePatch.H"
71#include "face.H"
72#include "labelList.H"
73#include "pointField.H"
74#include "indexedOctree.H"
76#include "volumeType.H"
78
79// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
80
81namespace Foam
82{
86
87class Time;
88class Random;
90
91/*---------------------------------------------------------------------------*\
92 Class backgroundMeshDecomposition Declaration
93\*---------------------------------------------------------------------------*/
96{
97 // Private data
98
99 //- Method details dictionary
100 //dictionary coeffsDict_;
101
102 //- Reference to runtime
103 const Time& runTime_;
104
105 //- Reference to surface
106 const conformationSurfaces& geometryToConformTo_;
107
108 //- Random number generator
109 Random& rndGen_;
110
111 //- Mesh stored on for this processor, specifying the domain that it
112 //- is responsible for.
113 fvMesh mesh_;
114
115 //- Refinement object
116 hexRef8 meshCutter_;
117
118 //- Patch containing an independent representation of the surface of the
119 // mesh of this processor
120 autoPtr<bPatch> boundaryFacesPtr_;
121
122 //- Search tree for the boundaryFaces_ patch
124
125 //- The bounds of all background meshes on all processors
126 treeBoundBoxList allBackgroundMeshBounds_;
127
128 //- The overall bounds of all of the background meshes, used to test if
129 // a point that is not found on any processor is in the domain at all
130 treeBoundBox globalBackgroundBounds_;
131
132 //- merge distance required by fvMeshDistribute
133 scalar mergeDist_;
134
135 //- Scale of a cell span vs cell size used to decide to refine a cell
136 scalar spanScale_;
137
138 //- Smallest minimum cell size allowed, i.e. to avoid high initial
139 // refinement of areas of small size
140 scalar minCellSizeLimit_;
141
142 //- Minimum normal level of refinement
143 label minLevels_;
144
145 //- How fine should the initial sample of the volume a box be to
146 // investigate the local cell size
147 label volRes_;
148
149 //- Allowed factor above the average cell weight before a background
150 // cell needs to be split
151 scalar maxCellWeightCoeff_;
152
153
154 // Private Member Functions
155
156 void initialRefinement();
157
158 //- Print details of the decomposed mesh
159 void printMeshData(const polyMesh& mesh) const;
160
161 //- Estimate the number of vertices that will be in this cell, returns
162 // true if the cell is to be split because of the density ratio inside
163 // it
164 bool refineCell
165 (
166 label celli,
167 volumeType volType,
168 scalar& weightEstimate
169 ) const;
170
171 //- Select cells for refinement at the surface of the geometry to be
172 // meshed
173 labelList selectRefinementCells
174 (
175 List<volumeType>& volumeStatus,
176 volScalarField& cellWeights
177 ) const;
178
179 //- Build the surface patch and search tree
180 void buildPatchAndTree();
181
182 //- No copy construct
184 (
186 ) = delete;
187
188 //- No copy assignment
189 void operator=(const backgroundMeshDecomposition&) = delete;
190
191
192public:
193
194 //- Runtime type information
195 ClassName("backgroundMeshDecomposition");
196
197
198 // Constructors
199
200 //- Construct from components in foamyHexMesh operation
202 (
203 const Time& runTime,
204 Random& rndGen,
205 const conformationSurfaces& geometryToConformTo,
206 const dictionary& coeffsDict,
207 const fileName& decompDictFile = ""
208 );
209
210
211 //- Destructor
213
214
215 // Member Functions
216
217 //- Build a mapDistribute for the supplied destination processor data
218 static autoPtr<mapDistribute> buildMap(const List<label>& toProc);
219
220 //- Redistribute the background mesh based on a supplied weight field,
221 // returning a map to use to redistribute vertices.
223 (
224 volScalarField& cellWeights
225 );
226
227 //- Distribute supplied the points to the appropriate processor
228 template<class PointType>
230
231 //- Is the given position inside the domain of this decomposition
232 bool positionOnThisProcessor(const point& pt) const;
233
234 //- Are the given positions inside the domain of this decomposition
236
237 //- Does the given box overlap the faces of the boundary of this
238 // processor
239 bool overlapsThisProcessor(const treeBoundBox& box) const;
240
241 //- Does the given sphere overlap the faces of the boundary of this
242 // processor
244 (
245 const point& centre,
246 const scalar radiusSqr
247 ) const;
248
249 //- Find nearest intersection of line between start and end, (exposing
250 // underlying indexedOctree)
252 (
253 const point& start,
254 const point& end
255 ) const;
256
257 //- Find any intersection of line between start and end, (exposing
258 // underlying indexedOctree)
260 (
261 const point& start,
262 const point& end
263 ) const;
264
265 //- What processor is the given position on?
266 template<class PointType>
268
269 //- What is the nearest processor to the given position?
271
272 //- Which processors are intersected by the line segment, returns all
273 // processors whose boundary patch is intersected by the sphere. By
274 // default this does not return the processor that the query is
275 // launched from, it is assumed that the point is on that processor.
276 // The index data member of the pointIndexHit is replaced with the
277 // processor index.
279 (
280 const List<point>& starts,
281 const List<point>& ends,
282 bool includeOwnProcessor = false
283 ) const;
286 (
287 const point& centre,
288 const scalar& radiusSqr
289 ) const;
292 (
293 const point& centre,
294 const scalar radiusSqr
295 ) const;
296
297// //- Which processors overlap the given sphere, returns all processors
298// // whose boundary patch is touched by the sphere or whom the sphere
299// // is inside. By default this does not return the processor that the
300// // query is launched from, it is assumed that the point is on that
301// // processor.
302// labelListList overlapsProcessors
303// (
304// const List<point>& centres,
305// const List<scalar>& radiusSqrs,
306// const Delaunay& T,
307// bool includeOwnProcessor
308// ) const;
309
310 // Access
311
312 //- Return access to the underlying mesh
313 inline const fvMesh& mesh() const;
314
315 //- Return access to the underlying tree
316 inline const indexedOctree<treeDataBPatch>& tree() const;
317
318 //- Return the boundBox of this processor
319 inline const treeBoundBox& procBounds() const;
320
321 //- Return the cell level of the underlying mesh
322 inline const labelList& cellLevel() const;
323
324 //- Return the point level of the underlying mesh
325 inline const labelList& pointLevel() const;
326
327};
328
329
330// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
331
332} // End namespace Foam
333
334// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
335
337
338#ifdef NoRepository
340#endif
341
342// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
343
344#endif
345
346// ************************************************************************* //
CGAL data structures used for 3D Delaunay meshing.
This class describes the interaction of (usually) a face and a point. It carries the info of a succes...
Definition: PointIndexHit.H:66
A list of faces which address into the list of points.
Random number generator.
Definition: Random.H:60
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
Store a background polyMesh to use for the decomposition of space and queries for parallel conformalV...
ClassName("backgroundMeshDecomposition")
Runtime type information.
backgroundMeshDecomposition(const Time &runTime, Random &rndGen, const conformationSurfaces &geometryToConformTo, const dictionary &coeffsDict, const fileName &decompDictFile="")
Construct from components in foamyHexMesh operation.
labelList processorNearestPosition(const List< point > &pts) const
What is the nearest processor to the given position?
pointIndexHit findLineAny(const point &start, const point &end) const
Find any intersection of line between start and end, (exposing.
bool overlapsThisProcessor(const treeBoundBox &box) const
Does the given box overlap the faces of the boundary of this.
bool positionOnThisProcessor(const point &pt) const
Is the given position inside the domain of this decomposition.
const treeBoundBox & procBounds() const
Return the boundBox of this processor.
labelList processorPosition(const List< PointType > &pts) const
What processor is the given position on?
autoPtr< mapDistributePolyMesh > distribute(volScalarField &cellWeights)
Redistribute the background mesh based on a supplied weight field,.
boolList positionOnThisProcessor(const List< point > &pts) const
Are the given positions inside the domain of this decomposition.
static autoPtr< mapDistribute > buildMap(const List< label > &toProc)
Build a mapDistribute for the supplied destination processor data.
autoPtr< mapDistribute > distributePoints(List< PointType > &points) const
Distribute supplied the points to the appropriate processor.
List< List< pointIndexHit > > intersectsProcessors(const List< point > &starts, const List< point > &ends, bool includeOwnProcessor=false) const
Which processors are intersected by the line segment, returns all.
~backgroundMeshDecomposition()=default
Destructor.
const labelList & cellLevel() const
Return the cell level of the underlying mesh.
const fvMesh & mesh() const
Return access to the underlying mesh.
const indexedOctree< treeDataBPatch > & tree() const
Return access to the underlying tree.
const labelList & pointLevel() const
Return the point level of the underlying mesh.
pointIndexHit findLine(const point &start, const point &end) const
Find nearest intersection of line between start and end, (exposing.
bool overlapsOtherProcessors(const point &centre, const scalar &radiusSqr) const
labelList overlapProcessors(const point &centre, const scalar radiusSqr) const
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
A class for handling file names.
Definition: fileName.H:76
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:91
Refinement of (split) hexes using polyTopoChange.
Definition: hexRef8.H:68
Non-pointer based hierarchical recursive searching.
Definition: indexedOctree.H:74
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:81
Container with cells to refine. Refinement given as single direction.
Definition: refineCell.H:57
Standard boundBox with extra functionality for use in octree.
Definition: treeBoundBox.H:89
Encapsulation of data needed to search on PrimitivePatches.
An enumeration wrapper for classification of a location as being inside/outside of a volume.
Definition: volumeType.H:61
#define ClassName(TypeNameString)
Add typeName information from argument TypeNameString to a class.
Definition: className.H:67
engineTime & runTime
const pointField & points
Namespace for OpenFOAM.
treeDataPrimitivePatch< bPatch > treeDataBPatch
PrimitivePatch<::Foam::List< face >, const pointField > bPatch
Random rndGen
Definition: createFields.H:23