DelaunayMesh.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) 2012-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::DelaunayMesh
28
29Description
30 The vertex and cell classes must have an index defined
31
32SourceFiles
33 DelaunayMeshI.H
34 DelaunayMesh.C
35 DelaunayMeshIO.C
36
37\*---------------------------------------------------------------------------*/
38
39#ifndef DelaunayMesh_H
40#define DelaunayMesh_H
41
42#include "labelPairHashes.H"
43#include "boundBox.H"
44#include "indexedVertex.H"
46#include "Time.H"
47#include "autoPtr.H"
48
49// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
50
51namespace Foam
52{
53
54class fvMesh;
55
56/*---------------------------------------------------------------------------*\
57 Class DelaunayMesh Declaration
58\*---------------------------------------------------------------------------*/
59
60template<class Triangulation>
61class DelaunayMesh
62:
63 public Triangulation
64{
65public:
69 typedef typename Triangulation::Edge Edge;
70 typedef typename Triangulation::Point Point;
71 typedef typename Triangulation::Facet Facet;
72
79
80private:
81
82 // Private data
83
84 //- Keep track of the number of vertices that have been added.
85 // This allows a unique index to be assigned to each vertex.
86 mutable label vertexCount_;
87
88 //- Keep track of the number of cells that have been added.
89 // This allows a unique index to be assigned to each cell.
90 mutable label cellCount_;
91
92 //- Reference to Time
93 const Time& runTime_;
94
95 //- Spatial sort traits to use with a pair of point pointers and an int.
96 // Taken from a post on the CGAL lists: 2010-01/msg00004.html by
97 // Sebastien Loriot (Geometry Factory).
98 struct Traits_for_spatial_sort
99 :
100 public Triangulation::Geom_traits
101 {
102 typedef typename Triangulation::Geom_traits Gt;
103
104 typedef std::pair<const typename Triangulation::Point*, label>
105 Point_3;
107 struct Less_x_3
109 bool operator()(const Point_3& p, const Point_3& q) const;
110 };
112 struct Less_y_3
114 bool operator()(const Point_3& p, const Point_3& q) const;
115 };
117 struct Less_z_3
119 bool operator()(const Point_3& p, const Point_3& q) const;
120 };
121
122 Less_x_3 less_x_3_object() const;
123 Less_y_3 less_y_3_object() const;
124 Less_z_3 less_z_3_object() const;
125 };
126
127
128 // Private Member Functions
129
130 void sortFaces
131 (
132 faceList& faces,
133 labelList& owner,
134 labelList& neighbour
135 ) const;
136
137 void addPatches
138 (
139 const label nInternalFaces,
140 faceList& faces,
141 labelList& owner,
143 const List<DynamicList<face>>& patchFaces,
144 const List<DynamicList<label>>& patchOwners
145 ) const;
146
147 //- No copy construct
149
150 //- No copy assignment
151 void operator=(const DelaunayMesh<Triangulation>&) = delete;
152
153
154public:
155
156 // Constructors
157
158 //- Construct from components
159 explicit DelaunayMesh(const Time& runTime);
162 (
163 const Time& runTime,
164 const word& meshName
165 );
166
167
168 //- Destructor
170
171
172 // Member Functions
173
174 // Access
175
176 //- Return a reference to the Time object
177 inline const Time& time() const;
178
179
180 // Check
181
182 //- Write the cpuTime to screen
183 inline void timeCheck
184 (
185 const string& description,
186 const bool check = true
187 ) const;
188
189
190 // Indexing functions
191
192 //- Create a new unique cell index and return
193 inline label getNewCellIndex() const;
194
195 //- Create a new unique vertex index and return
196 inline label getNewVertexIndex() const;
197
198 //- Return the cell count (the next unique cell index)
199 inline label cellCount() const;
200
201 //- Return the vertex count (the next unique vertex index)
202 inline label vertexCount() const;
203
204 //- Set the cell count to zero
205 inline void resetCellCount();
206
207 //- Set the vertex count to zero
208 inline void resetVertexCount();
209
210
211 // Triangulation manipulation functions
212
213 //- Clear the entire triangulation
214 void reset();
215
216 //- Insert the list of vertices (calls rangeInsertWithInfo)
218 (
219 const List<Vb>& vertices,
220 const bool reIndex
221 );
222
223 //- Function inserting points into a triangulation and setting the
224 // index and type data of the point in the correct order. This is
225 // faster than inserting points individually.
226 //
227 // Adapted from a post on the CGAL lists: 2010-01/msg00004.html by
228 // Sebastien Loriot (Geometry Factory).
229 template<class PointIterator>
231 (
232 PointIterator begin,
233 PointIterator end,
234 bool printErrors = false,
235 bool reIndex = true
236 );
237
238
239 // Write
240
241 //- Write mesh statistics to stream
242 void printInfo(Ostream& os) const;
243
244 //- Write vertex statistics in the form of a table to stream
245 void printVertexInfo(Ostream& os) const;
246
247 //- Create an fvMesh from the triangulation.
248 // The mesh is not parallel consistent - only used for viewing
250 (
251 const fileName& name,
252 labelPairLookup& vertexMap,
253 labelList& cellMap,
254 const bool writeDelaunayData = true
255 ) const;
256};
257
258
259// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
260
261} // End namespace Foam
262
263// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
264
265#include "DelaunayMeshI.H"
266
267// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
268
269#ifdef NoRepository
270 #include "DelaunayMesh.C"
271#endif
272
273// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
274
275#endif
276
277// ************************************************************************* //
CGAL data structures used for 3D Delaunay meshing.
The vertex and cell classes must have an index defined.
Definition: DelaunayMesh.H:63
Triangulation::Vertex_handle Vertex_handle
Definition: DelaunayMesh.H:67
void resetVertexCount()
Set the vertex count to zero.
Triangulation::Facet Facet
Definition: DelaunayMesh.H:70
const Time & time() const
Return a reference to the Time object.
Definition: DelaunayMeshI.H:31
Map< label > insertPoints(const List< Vb > &vertices, const bool reIndex)
Insert the list of vertices (calls rangeInsertWithInfo)
Triangulation::Finite_vertices_iterator Finite_vertices_iterator
Definition: DelaunayMesh.H:73
Map< label > rangeInsertWithInfo(PointIterator begin, PointIterator end, bool printErrors=false, bool reIndex=true)
Function inserting points into a triangulation and setting the.
label vertexCount() const
Return the vertex count (the next unique vertex index)
label getNewCellIndex() const
Create a new unique cell index and return.
Definition: DelaunayMeshI.H:65
void printVertexInfo(Ostream &os) const
Write vertex statistics in the form of a table to stream.
Triangulation::Finite_facets_iterator Finite_facets_iterator
Definition: DelaunayMesh.H:77
DelaunayMesh(const Time &runTime, const word &meshName)
Triangulation::Cell_handle Cell_handle
Definition: DelaunayMesh.H:66
~DelaunayMesh()
Destructor.
autoPtr< polyMesh > createMesh(const fileName &name, labelPairLookup &vertexMap, labelList &cellMap, const bool writeDelaunayData=true) const
Create an fvMesh from the triangulation.
label getNewVertexIndex() const
Create a new unique vertex index and return.
Definition: DelaunayMeshI.H:80
Triangulation::Finite_cells_iterator Finite_cells_iterator
Definition: DelaunayMesh.H:75
Triangulation::Edge Edge
Definition: DelaunayMesh.H:68
DelaunayMesh(const Time &runTime)
Construct from components.
void reset()
Clear the entire triangulation.
label cellCount() const
Return the cell count (the next unique cell index)
Definition: DelaunayMeshI.H:95
void resetCellCount()
Set the cell count to zero.
void printInfo(Ostream &os) const
Write mesh statistics to stream.
Triangulation::Point Point
Definition: DelaunayMesh.H:69
void timeCheck(const string &description, const bool check=true) const
Write the cpuTime to screen.
Definition: DelaunayMeshI.H:39
A 1D vector of objects of type <T> that resizes itself as necessary to accept the new objects.
Definition: DynamicList.H:72
A HashTable to objects of type <T> with a label key.
Definition: Map.H:60
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:62
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
Definition: PtrList.H:73
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
CellSizeDelaunay::Cell_handle Cell_handle
CellSizeDelaunay::Vertex_handle Vertex_handle
CellSizeDelaunay::Point Point
A class for handling file names.
Definition: fileName.H:76
A class for handling words, derived from Foam::string.
Definition: word.H:68
volScalarField & p
engineTime & runTime
OBJstream os(runTime.globalPath()/outputName)
A HashTable to objects of type <T> with a labelPair key. The hashing is based on labelPair (FixedList...
Namespace for OpenFOAM.
static void check(const int retVal, const char *what)
pointField vertices(const blockVertexList &bvl)
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
PtrList< dictionary > patchDicts
Definition: readKivaGrid.H:532
bool operator()(const Point_3 &p, const Point_3 &q) const
bool operator()(const Point_3 &p, const Point_3 &q) const
bool operator()(const Point_3 &p, const Point_3 &q) const