Go to the documentation of this file.
39 template<
class Triangulation>
49 template<
class Triangulation>
61 Info<<
"Reading " << meshName <<
" from " <<
runTime.timeName() <<
endl;
118 List<Vb> pointsToInsert(pts.size());
120 forAll(pointsToInsert, pI)
134 pointsToInsert.begin(),
135 pointsToInsert.end(),
140 vertexCount_ = Triangulation::number_of_vertices();
147 template<
class Triangulation>
160 template<
class Triangulation>
163 Info<<
"Clearing triangulation" <<
endl;
169 Finite_vertices_iterator vit = Triangulation::finite_vertices_begin();
170 vit != Triangulation::finite_vertices_end();
187 vertices.last().fixed() = vit->fixed();
198 Info<<
"Inserted " << vertexCount() <<
" fixed points" <<
endl;
202 template<
class Triangulation>
209 return rangeInsertWithInfo
219 template<
class Triangulation>
227 return typename Gt::Less_x_3()(*(
p.first), *(q.first));
230 template<
class Triangulation>
238 return typename Gt::Less_y_3()(*(
p.first), *(q.first));
241 template<
class Triangulation>
249 return typename Gt::Less_z_3()(*(
p.first), *(q.first));
252 template<
class Triangulation>
260 template<
class Triangulation>
268 template<
class Triangulation>
277 template<
class Triangulation>
278 template<
class Po
intIterator>
294 > vectorPairPointIndex;
296 vectorPairPointIndex
points;
299 for (PointIterator it =
begin; it !=
end; ++it)
303 std::make_pair(&(it->point()),
count++)
313 Traits_for_spatial_sort()
318 Map<label> oldToNewIndex(
points.size());
322 typename vectorPairPointIndex::const_iterator
p =
points.begin();
327 const size_t checkInsertion = Triangulation::number_of_vertices();
329 hint = this->
insert(*(p->first), hint);
331 const Vb& vert = *(
begin +
p->second);
333 if (checkInsertion != Triangulation::number_of_vertices() - 1)
337 Vertex_handle nearV =
338 Triangulation::nearest_vertex(*(
p->first));
340 Pout<<
"Failed insertion : " << vert.
info()
341 <<
" nearest : " << nearV->info();
346 const label oldIndex = vert.
index();
347 hint->index() = getNewVertexIndex();
351 oldToNewIndex.insert(oldIndex, hint->index());
354 hint->type() = vert.
type();
361 return oldToNewIndex;
srcOptions insert("case", fileName(rootDirSource/caseDirSource))
Foam::tensor & alignment()
constexpr auto begin(C &c) -> decltype(c.begin())
Return iterator to the beginning of the container c.
CGAL::indexedVertex< K > Vb
static word meshSubDir
Return the mesh sub-directory name (usually "polyMesh")
IOField< label > labelIOField
labelField with IO.
An indexed form of CGAL::Triangulation_vertex_base_3<K> used to keep track of the Delaunay vertices i...
A HashTable to objects of type <T> with a label key.
Ostream & endl(Ostream &os)
Add newline and flush stream.
Map< label > rangeInsertWithInfo(PointIterator begin, PointIterator end, bool printErrors=false, bool reIndex=true)
Function inserting points into a triangulation and setting the.
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
#define forAll(list, i)
Loop across all elements in list.
A HashTable to objects of type <T> with a labelPair key. The hashing is based on labelPair (FixedList...
vectorIOField pointIOField
pointIOField is a vectorIOField.
messageStream Info
Information stream (uses stdout - output is on the master only)
Map< label > insertPoints(const List< Vb > &vertices, const bool reIndex)
Insert the list of vertices (calls rangeInsertWithInfo)
Foam::InfoProxy< indexedVertex< Gt, Vb > > info() const
Info proxy, to print information to a stream.
Foam::scalar & targetCellSize()
pointField vertices(const blockVertexList &bvl)
PointFrompoint toPoint(const Foam::point &p)
constexpr auto end(C &c) -> decltype(c.end())
Return iterator to the end of the container c.
~DelaunayMesh()
Destructor.
unsigned int count(const UList< bool > &bools, const bool val=true)
Count number of 'true' entries.
The vertex and cell classes must have an index defined.
void reset()
Clear the entire triangulation.