Merge remote-tracking branch 'cgal/6.0.x-branch' into 'cgal/master'

This commit is contained in:
Sébastien Loriot 2025-08-05 16:28:05 +02:00
commit bc0000919f
2 changed files with 100 additions and 76 deletions

View File

@ -40,17 +40,17 @@
namespace CGAL { namespace CGAL {
// parameterization of the hierarchy
//const float Triangulation_hierarchy_2__ratio = 30.0;
const int Triangulation_hierarchy_2__ratio = 30;
const int Triangulation_hierarchy_2__minsize = 20;
const int Triangulation_hierarchy_2__maxlevel = 5;
// maximal number of points is 30^5 = 24 millions !
template <class Tr_> template <class Tr_>
class Triangulation_hierarchy_2 class Triangulation_hierarchy_2
: public Tr_ : public Tr_
{ {
// parameterization of the hierarchy
// maximal number of points is 30^5 = 24 millions !
enum Ratio : int { ratio = 30 };
enum { minsize = 20};
enum { maxlevel = 5};
public: public:
typedef Tr_ Tr_Base; typedef Tr_ Tr_Base;
typedef typename Tr_Base::Geom_traits Geom_traits; typedef typename Tr_Base::Geom_traits Geom_traits;
@ -84,13 +84,13 @@ public:
private: private:
void init_hierarchy() { void init_hierarchy() {
hierarchy[0] = this; hierarchy[0] = this;
for(int i=1; i<Triangulation_hierarchy_2__maxlevel; ++i) for(int i=1; i<maxlevel; ++i)
hierarchy[i] = &hierarchy_triangulations[i-1]; hierarchy[i] = &hierarchy_triangulations[i-1];
} }
// here is the stack of triangulations which form the hierarchy // here is the stack of triangulations which form the hierarchy
std::array<Tr_Base,Triangulation_hierarchy_2__maxlevel-1> hierarchy_triangulations; std::array<Tr_Base,maxlevel-1> hierarchy_triangulations;
std::array<Tr_Base*,Triangulation_hierarchy_2__maxlevel> hierarchy; std::array<Tr_Base*,maxlevel> hierarchy;
boost::rand48 random; boost::rand48 random;
public: public:
@ -163,7 +163,7 @@ public:
// hints[i] is the face of the previously inserted point in level i. // hints[i] is the face of the previously inserted point in level i.
// Thanks to spatial sort, they are better hints than what the hierarchy // Thanks to spatial sort, they are better hints than what the hierarchy
// would give us. // would give us.
Face_handle hints[Triangulation_hierarchy_2__maxlevel]; Face_handle hints[maxlevel];
for (typename std::vector<Point>::const_iterator p = points.begin(), end = points.end(); for (typename std::vector<Point>::const_iterator p = points.begin(), end = points.end();
p != end; ++p) p != end; ++p)
{ {
@ -260,7 +260,7 @@ private:
int& li, int& li,
Face_handle loc, Face_handle loc,
Face_handle Face_handle
pos[Triangulation_hierarchy_2__maxlevel]) const; pos[maxlevel]) const;
int random_level(); int random_level();
@ -291,7 +291,7 @@ Triangulation_hierarchy_2<Tr_>::
Triangulation_hierarchy_2(const Geom_traits& traits) Triangulation_hierarchy_2(const Geom_traits& traits)
: Tr_Base(traits) : Tr_Base(traits)
, hierarchy_triangulations( , hierarchy_triangulations(
make_filled_array<Triangulation_hierarchy_2__maxlevel-1, make_filled_array<maxlevel-1,
Tr_Base>(traits)) Tr_Base>(traits))
{ {
init_hierarchy(); init_hierarchy();
@ -325,7 +325,7 @@ Triangulation_hierarchy_2<Tr_>::
copy_triangulation(const Triangulation_hierarchy_2<Tr_> &tr) copy_triangulation(const Triangulation_hierarchy_2<Tr_> &tr)
{ {
{ {
for(int i=0;i<Triangulation_hierarchy_2__maxlevel;++i) for(int i=0;i<maxlevel;++i)
hierarchy[i]->copy_triangulation(*tr.hierarchy[i]); hierarchy[i]->copy_triangulation(*tr.hierarchy[i]);
} }
@ -343,7 +343,7 @@ copy_triangulation(const Triangulation_hierarchy_2<Tr_> &tr)
add_hidden_vertices_into_map(Weighted_tag(), V); add_hidden_vertices_into_map(Weighted_tag(), V);
{ {
for(int i=1;i<Triangulation_hierarchy_2__maxlevel;++i) { for(int i=1;i<maxlevel;++i) {
for( Finite_vertices_iterator it=hierarchy[i]->finite_vertices_begin(); for( Finite_vertices_iterator it=hierarchy[i]->finite_vertices_begin();
it != hierarchy[i]->finite_vertices_end(); ++it) { it != hierarchy[i]->finite_vertices_end(); ++it) {
// down pointer goes in original instead in copied triangulation // down pointer goes in original instead in copied triangulation
@ -396,7 +396,7 @@ void
Triangulation_hierarchy_2<Tr_>:: Triangulation_hierarchy_2<Tr_>::
clear() clear()
{ {
for(int i=0;i<Triangulation_hierarchy_2__maxlevel;++i) for(int i=0;i<maxlevel;++i)
if(hierarchy[i]) hierarchy[i]->clear(); if(hierarchy[i]) hierarchy[i]->clear();
} }
@ -410,7 +410,7 @@ is_valid(bool verbose, int level) const
int i; int i;
Finite_vertices_iterator it; Finite_vertices_iterator it;
//verify correctness of triangulation at all levels //verify correctness of triangulation at all levels
for(i=0;i<Triangulation_hierarchy_2__maxlevel;++i) { for(i=0;i<maxlevel;++i) {
if(verbose) // pirnt number of vertices at each level if(verbose) // pirnt number of vertices at each level
std::cout << "number_of_vertices " std::cout << "number_of_vertices "
<< hierarchy[i]->number_of_vertices() << std::endl; << hierarchy[i]->number_of_vertices() << std::endl;
@ -421,13 +421,13 @@ is_valid(bool verbose, int level) const
it != hierarchy[0]->finite_vertices_end(); ++it) it != hierarchy[0]->finite_vertices_end(); ++it)
result = result && ( it->down() == Vertex_handle()); result = result && ( it->down() == Vertex_handle());
//verify that other levels have down pointer and reciprocal link is fine //verify that other levels have down pointer and reciprocal link is fine
for(i=1;i<Triangulation_hierarchy_2__maxlevel;++i) for(i=1;i<maxlevel;++i)
for( it = hierarchy[i]->finite_vertices_begin(); for( it = hierarchy[i]->finite_vertices_begin();
it != hierarchy[i]->finite_vertices_end(); ++it) it != hierarchy[i]->finite_vertices_end(); ++it)
result = result && result = result &&
( &*(it->down()->up()) == &*(it) ); ( &*(it->down()->up()) == &*(it) );
//verify that levels have up pointer and reciprocal link is fine //verify that levels have up pointer and reciprocal link is fine
for(i=0;i<Triangulation_hierarchy_2__maxlevel-1;++i) for(i=0;i<maxlevel-1;++i)
for( it = hierarchy[i]->finite_vertices_begin(); for( it = hierarchy[i]->finite_vertices_begin();
it != hierarchy[i]->finite_vertices_end(); ++it) it != hierarchy[i]->finite_vertices_end(); ++it)
result = result && ( it->up() == Vertex_handle() || result = result && ( it->up() == Vertex_handle() ||
@ -445,7 +445,7 @@ insert(const Point &p, Face_handle loc)
Locate_type lt; Locate_type lt;
int i; int i;
// locate using hierarchy // locate using hierarchy
Face_handle positions[Triangulation_hierarchy_2__maxlevel]; Face_handle positions[maxlevel];
locate_in_all(p,lt,i,loc,positions); locate_in_all(p,lt,i,loc,positions);
Vertex_handle vertex=hierarchy[0]->Tr_Base::insert(p,lt,positions[0],i); Vertex_handle vertex=hierarchy[0]->Tr_Base::insert(p,lt,positions[0],i);
Vertex_handle previous=vertex; Vertex_handle previous=vertex;
@ -480,7 +480,7 @@ insert(const Point& p,
// locate using hierarchy // locate using hierarchy
Locate_type ltt; Locate_type ltt;
int lii; int lii;
Face_handle positions[Triangulation_hierarchy_2__maxlevel]; Face_handle positions[maxlevel];
locate_in_all(p,ltt,lii,loc,positions); locate_in_all(p,ltt,lii,loc,positions);
//insert in higher levels //insert in higher levels
int level = 1; int level = 1;
@ -514,7 +514,7 @@ remove(Vertex_handle v )
while(1){ while(1){
hierarchy[l++]->remove(v); hierarchy[l++]->remove(v);
if (u == Vertex_handle()) break; if (u == Vertex_handle()) break;
if (l >= Triangulation_hierarchy_2__maxlevel) break; if (l >= maxlevel) break;
v=u; u=v->up(); v=u; u=v->up();
} }
} }
@ -531,7 +531,7 @@ remove_and_give_new_faces(Vertex_handle v, OutputItFaces fit)
if(l==0) hierarchy[l++]->remove_and_give_new_faces(v, fit); if(l==0) hierarchy[l++]->remove_and_give_new_faces(v, fit);
else hierarchy[l++]->remove(v); else hierarchy[l++]->remove(v);
if (u == Vertex_handle()) break; if (u == Vertex_handle()) break;
if (l >= Triangulation_hierarchy_2__maxlevel) break; if (l >= maxlevel) break;
v=u; u=v->up(); v=u; u=v->up();
} }
} }
@ -565,16 +565,19 @@ template <class Tr_>
typename Triangulation_hierarchy_2<Tr_>::Vertex_handle typename Triangulation_hierarchy_2<Tr_>::Vertex_handle
Triangulation_hierarchy_2<Tr_>:: Triangulation_hierarchy_2<Tr_>::
move_if_no_collision(Vertex_handle v, const Point &p) { move_if_no_collision(Vertex_handle v, const Point &p) {
Vertex_handle u=v->up(), norm = v;
int l = 0 ; CGAL_precondition(!this->is_infinite(v));
while(1) { if(v->point() == p) return v;
Vertex_handle w = hierarchy[l++]->move_if_no_collision(v, p); Vertex_handle ans = hierarchy[0]->move_if_no_collision(v, p);
if(w != v) return w; if(ans != v) return ans; // ans is an existing vertex at p and v was not changed
if (u == Vertex_handle()) break; for (int l = 1; l < maxlevel; ++l) {
if (l >= Triangulation_hierarchy_2__maxlevel) break; Vertex_handle u = v->up();
v=u; u=v->up(); if (u == Vertex_handle())
break;
hierarchy[l]->move_if_no_collision(u, p);
v = u;
} }
return norm; return ans;
} }
template <class Tr_> template <class Tr_>
@ -597,22 +600,18 @@ Triangulation_hierarchy_2<Tr_>::
move_if_no_collision_and_give_new_faces(Vertex_handle v, const Point &p, move_if_no_collision_and_give_new_faces(Vertex_handle v, const Point &p,
OutputItFaces oif) OutputItFaces oif)
{ {
Vertex_handle u=v->up(), norm = v; CGAL_precondition(!this->is_infinite(v));
int l = 0 ; if(v->point() == p) return v;
while(1){ Vertex_handle ans = hierarchy[0]->move_if_no_collision_and_give_new_faces(v, p, oif);
Vertex_handle w; if(ans != v) return ans; // ans is an existing vertex at p and v was not changed
if(l == 0) for (int l = 1; l < maxlevel; ++l) {
w = Vertex_handle u = v->up();
hierarchy[l++]->move_if_no_collision_and_give_new_faces(v, p, oif); if (u == Vertex_handle())
else w = hierarchy[l++]->move_if_no_collision(v, p); break;
hierarchy[l]->move_if_no_collision(u, p);
if(w != v) return w; v = u;
if (u == Vertex_handle()) break;
if (l >= Triangulation_hierarchy_2__maxlevel) break;
v=u; u=v->up();
} }
return norm; return ans;
} }
template <class Tr_> template <class Tr_>
@ -627,7 +626,7 @@ Triangulation_hierarchy_2<Tr_>::insert_and_give_new_faces(const Point &p,
Locate_type lt; Locate_type lt;
int i; int i;
// locate using hierarchy // locate using hierarchy
Face_handle positions[Triangulation_hierarchy_2__maxlevel]; Face_handle positions[maxlevel];
locate_in_all(p,lt,i,loc,positions); locate_in_all(p,lt,i,loc,positions);
Vertex_handle vertex= Vertex_handle vertex=
hierarchy[0]->Tr_Base::insert_and_give_new_faces(p,lt,positions[0],i,oif); hierarchy[0]->Tr_Base::insert_and_give_new_faces(p,lt,positions[0],i,oif);
@ -666,7 +665,7 @@ insert_and_give_new_faces(const Point &p,
// locate using hierarchy // locate using hierarchy
Locate_type ltt; Locate_type ltt;
int lii; int lii;
Face_handle positions[Triangulation_hierarchy_2__maxlevel]; Face_handle positions[maxlevel];
locate_in_all(p,ltt,lii,loc,positions); locate_in_all(p,ltt,lii,loc,positions);
//insert in higher levels //insert in higher levels
int level = 1; int level = 1;
@ -686,7 +685,7 @@ typename Triangulation_hierarchy_2<Tr_>::Face_handle
Triangulation_hierarchy_2<Tr_>:: Triangulation_hierarchy_2<Tr_>::
locate(const Point& p, Locate_type& lt, int& li, Face_handle loc) const locate(const Point& p, Locate_type& lt, int& li, Face_handle loc) const
{ {
Face_handle positions[Triangulation_hierarchy_2__maxlevel]; Face_handle positions[maxlevel];
locate_in_all(p,lt,li,loc,positions); locate_in_all(p,lt,li,loc,positions);
return positions[0]; return positions[0];
} }
@ -708,11 +707,11 @@ locate_in_all(const Point& p,
Locate_type& lt, Locate_type& lt,
int& li, int& li,
Face_handle loc, Face_handle loc,
Face_handle pos[Triangulation_hierarchy_2__maxlevel]) const Face_handle pos[maxlevel]) const
{ {
Face_handle position; Face_handle position;
Vertex_handle nearest; Vertex_handle nearest;
int level = Triangulation_hierarchy_2__maxlevel; int level = maxlevel;
typename Geom_traits::Compare_distance_2 typename Geom_traits::Compare_distance_2
closer = geom_traits().compare_distance_2_object(); closer = geom_traits().compare_distance_2_object();
@ -721,7 +720,7 @@ locate_in_all(const Point& p,
// find the highest level with enough vertices that is at the same time 2D // find the highest level with enough vertices that is at the same time 2D
while ( (hierarchy[--level]->number_of_vertices() while ( (hierarchy[--level]->number_of_vertices()
< static_cast<size_type> (Triangulation_hierarchy_2__minsize )) < static_cast<size_type> (minsize ))
|| (hierarchy[level]->dimension()<2) ){ || (hierarchy[level]->dimension()<2) ){
if ( ! level) break; // do not go below 0 if ( ! level) break; // do not go below 0
} }
@ -729,7 +728,7 @@ locate_in_all(const Point& p,
level--; level--;
} }
for (int i=level+1; i<Triangulation_hierarchy_2__maxlevel;++i) pos[i]=nullptr; for (int i=level+1; i<maxlevel;++i) pos[i]=nullptr;
while(level > 0) { while(level > 0) {
pos[level]=position=hierarchy[level]->locate(p, position); pos[level]=position=hierarchy[level]->locate(p, position);
// locate at that level from "position" // locate at that level from "position"
@ -770,10 +769,10 @@ int
Triangulation_hierarchy_2<Tr_>:: Triangulation_hierarchy_2<Tr_>::
random_level() random_level()
{ {
boost::geometric_distribution<> proba(1.0/Triangulation_hierarchy_2__ratio); boost::geometric_distribution<> proba(1.0/double(ratio));
boost::variate_generator<boost::rand48&, boost::geometric_distribution<> > die(random, proba); boost::variate_generator<boost::rand48&, boost::geometric_distribution<> > die(random, proba);
return (std::min)(die(), Triangulation_hierarchy_2__maxlevel)-1; return (std::min)(die(), (int)maxlevel)-1;
} }

View File

@ -62,7 +62,7 @@ class Triangulation_hierarchy_3
{ {
// parameterization of the hierarchy // parameterization of the hierarchy
// maximal number of points is 30^5 = 24 millions ! // maximal number of points is 30^5 = 24 millions !
enum { ratio = 30 }; enum Ratio : int { ratio = 30 };
enum { minsize = 20}; enum { minsize = 20};
enum { maxlevel = 5}; enum { maxlevel = 5};
@ -504,26 +504,54 @@ is_valid(bool verbose, int level) const
bool result = true; bool result = true;
// verify correctness of triangulation at all levels // verify correctness of triangulation at all levels
for(int i=0; i<maxlevel; ++i) for(int i=0; i<maxlevel; ++i){
result = result && hierarchy[i]->is_valid(verbose, level); result = result && hierarchy[i]->is_valid(verbose, level);
if(verbose && (! result)){
// verify that lower level has no down pointers std::cerr << "triangulation at level " << i << " invalid" << std::endl;
}
}
// verify that lowest level has no down pointers
for( Finite_vertices_iterator it = hierarchy[0]->finite_vertices_begin(), for( Finite_vertices_iterator it = hierarchy[0]->finite_vertices_begin(),
end = hierarchy[0]->finite_vertices_end(); it != end; ++it) end = hierarchy[0]->finite_vertices_end(); it != end; ++it){
result = result && (it->down() == Vertex_handle()); result = result && (it->down() == Vertex_handle());
if(verbose && (! result)){
std::cerr << "lowest level has a down pointer" << std::endl;
}
}
// verify that other levels has down pointer and reciprocal link is fine // verify that other levels have down pointer and reciprocal link is fine
for(int j=1; j<maxlevel; ++j) for(int j=1; j<maxlevel; ++j)
for( Finite_vertices_iterator it = hierarchy[j]->finite_vertices_begin(), for( Finite_vertices_iterator it = hierarchy[j]->finite_vertices_begin(),
end = hierarchy[j]->finite_vertices_end(); it != end; ++it) end = hierarchy[j]->finite_vertices_end(); it != end; ++it)
result = result && &*(it) == &*(it->down()->up()); {
result = result && (it->down() != Vertex_handle());
if(verbose && (! result)){
std::cerr << "missing down pointer" << std::endl;
}
if(it->down() == Vertex_handle()){
return false;
}
result = result && Vertex_handle(it) == Vertex_handle(it->down()->up());
if(verbose && (! result)){
std::cerr << "wrong reciprocal link with down()" << std::endl;
}
result = result && Vertex_handle(it)->point() == Vertex_handle(it->down())->point();
if(verbose && (! result)){
std::cerr << "inconsistent vertex positions" << std::endl;
}
}
// verify that other levels has down pointer and reciprocal link is fine // verify that all levels have up pointer and reciprocal link is fine
for(int k=0; k<maxlevel-1; ++k) for(int k=0; k<maxlevel-1; ++k)
for( Finite_vertices_iterator it = hierarchy[k]->finite_vertices_begin(), for( Finite_vertices_iterator it = hierarchy[k]->finite_vertices_begin(),
end = hierarchy[k]->finite_vertices_end(); it != end; ++it) end = hierarchy[k]->finite_vertices_end(); it != end; ++it)
{
result = result && ( it->up() == Vertex_handle() || result = result && ( it->up() == Vertex_handle() ||
&*it == &*(it->up())->down() ); Vertex_handle(it) == Vertex_handle(it->up())->down() );
if(verbose && (! result)){
std::cerr << "wrong reciprocal link with up()" << std::endl;
}
}
return result; return result;
} }
@ -718,14 +746,13 @@ move_if_no_collision(Vertex_handle v, const Point & p)
{ {
CGAL_precondition(!this->is_infinite(v)); CGAL_precondition(!this->is_infinite(v));
if(v->point() == p) return v; if(v->point() == p) return v;
Vertex_handle ans; Vertex_handle ans = hierarchy[0]->move_if_no_collision(v, p);
for (int l = 0; l < maxlevel; ++l) { if(ans != v) return ans; // ans is an existing vertex at p and v was not changed
for (int l = 1; l < maxlevel; ++l) {
Vertex_handle u = v->up(); Vertex_handle u = v->up();
if(l) hierarchy[l]->move_if_no_collision(v, p);
else ans = hierarchy[l]->move_if_no_collision(v, p);
if(ans != v) return ans;
if (u == Vertex_handle()) if (u == Vertex_handle())
break; break;
hierarchy[l]->move_if_no_collision(u, p);
v = u; v = u;
} }
return ans; return ans;
@ -753,17 +780,15 @@ Triangulation_hierarchy_3<Tr>::
move_if_no_collision_and_give_new_cells( move_if_no_collision_and_give_new_cells(
Vertex_handle v, const Point & p, OutputItCells fit) Vertex_handle v, const Point & p, OutputItCells fit)
{ {
CGAL_precondition(!is_infinite(v)); CGAL_precondition(!this->is_infinite(v));
if(v->point() == p) return v; if(v->point() == p) return v;
Vertex_handle ans; Vertex_handle ans = hierarchy[0]->move_if_no_collision_and_give_new_cells(v, p, fit);
for (int l = 0; l < maxlevel; ++l) { if(ans != v) return ans; // ans is an existing vertex at p and v was not changed
for (int l = 1; l < maxlevel; ++l) {
Vertex_handle u = v->up(); Vertex_handle u = v->up();
if(l) hierarchy[l]->move_if_no_collision(v, p);
else ans =
hierarchy[l]->move_if_no_collision_and_give_new_cells(v, p, fit);
if(ans != v) return ans;
if (u == Vertex_handle()) if (u == Vertex_handle())
break; break;
hierarchy[l]->move_if_no_collision(u, p);
v = u; v = u;
} }
return ans; return ans;