diff --git a/Documentation/doc/biblio/geom.bib b/Documentation/doc/biblio/geom.bib index 0d135a1805f..b8898f23596 100644 --- a/Documentation/doc/biblio/geom.bib +++ b/Documentation/doc/biblio/geom.bib @@ -151997,3 +151997,17 @@ pages = {179--189} year={1998}, organization={Blackwell Publishers} } + +@inproceedings{boissonnat2009Delaunay, + TITLE = {{Incremental construction of the Delaunay graph in medium dimension}}, + AUTHOR = {Boissonnat, Jean-Daniel and Devillers, Olivier and Hornus, Samuel}, + URL = {https://hal.inria.fr/inria-00412437}, + BOOKTITLE = {{Annual Symposium on Computational Geometry}}, + ADDRESS = {Aarhus, Denmark}, + PAGES = {208-216}, + YEAR = {2009}, + MONTH = Jun, + PDF = {https://hal.inria.fr/inria-00412437/file/socg09.pdf}, + HAL_ID = {inria-00412437}, + HAL_VERSION = {v1}, +} diff --git a/Kernel_d/doc/Kernel_d/CGAL/Epick_d.h b/Kernel_d/doc/Kernel_d/CGAL/Epick_d.h index 9be521a4af1..7d91973363e 100644 --- a/Kernel_d/doc/Kernel_d/CGAL/Epick_d.h +++ b/Kernel_d/doc/Kernel_d/CGAL/Epick_d.h @@ -37,6 +37,9 @@ icc 15 work. \cgalModels `Kernel_d` \cgalModels `DelaunayTriangulationTraits` +\cgalModels `RegularTriangulationTraits` +\cgalModels `SearchTraits` +\cgalModels `RangeSearchTraits` \sa `CGAL::Cartesian_d` \sa `CGAL::Homogeneous_d` @@ -74,9 +77,25 @@ Cartesian_const_iterator_d cartesian_begin()const; Cartesian_const_iterator_d cartesian_end()const; }; +/*! +represents a weighted point in the Euclidean space +\cgalModels `DefaultConstructible` +\cgalModels `Assignable` +*/ +class Weighted_point_d { +public: +/*! introduces a weighted point with point p and weight w. */ +Weighted_point_d(const Point_d& p, const double& w); +/*! extracts the point of a weighted point. */ +Point_d point() const; +/*! extracts the weight of a weighted point. */ +double weight() const; +}; + /*! \cgalModels `Kernel_d::Center_of_sphere_d` */ -struct Construct_circumcenter_d { +class Construct_circumcenter_d { +public: /*! returns the center of the sphere defined by `A=tuple[first,last)`. The sphere is centered in the affine hull of A and passes through all the points of A. The order of the points of A does not matter. \pre A is affinely independant. \tparam ForwardIterator has `Epick_d::Point_d` as value type. @@ -84,7 +103,8 @@ struct Construct_circumcenter_d { template Point_d operator()(ForwardIterator first, ForwardIterator last); }; -struct Compute_squared_radius_d { +class Compute_squared_radius_d { +public: /*! returns the radius of the sphere defined by `A=tuple[first,last)`. The sphere is centered in the affine hull of A and passes through all the points of A. The order of the points of A does not matter. \pre A is affinely independant. \tparam ForwardIterator has `Epick_d::Point_d` as value type. @@ -94,7 +114,8 @@ Point_d operator()(ForwardIterator first, ForwardIterator last); }; /*! \cgalModels `Kernel_d::Side_of_bounded_sphere_d` */ -struct Side_of_bounded_sphere_d { +class Side_of_bounded_sphere_d { +public: /*! returns the relative position of point p to the sphere defined by `A=tuple[first,last)`. The sphere is centered in the affine hull of A and passes through all the points of A. The order of the points of A does not matter. \pre A is affinely independant. \tparam ForwardIterator has `Epick_d::Point_d` as value type. diff --git a/Kernel_d/doc/Kernel_d/dependencies b/Kernel_d/doc/Kernel_d/dependencies index fdf381a7a17..eea7200c3d4 100644 --- a/Kernel_d/doc/Kernel_d/dependencies +++ b/Kernel_d/doc/Kernel_d/dependencies @@ -6,4 +6,4 @@ Circulator Stream_support Number_types Triangulation - +Spatial_searching diff --git a/NewKernel_d/include/CGAL/Epeck_d.h b/NewKernel_d/include/CGAL/Epeck_d.h new file mode 100644 index 00000000000..52bce84ca3f --- /dev/null +++ b/NewKernel_d/include/CGAL/Epeck_d.h @@ -0,0 +1,53 @@ +// Copyright (c) 2014 +// INRIA Saclay-Ile de France (France) +// +// This file is part of CGAL (www.cgal.org); you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public License as +// published by the Free Software Foundation; either version 3 of the License, +// or (at your option) any later version. +// +// Licensees holding a valid commercial license may use this file in +// accordance with the commercial license agreement provided with the software. +// +// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE +// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE. +// +// $URL$ +// $Id$ +// +// Author(s) : Marc Glisse + +#ifndef CGAL_EPECK_D_H +#define CGAL_EPECK_D_H +#include +#include +#include +#include + + +namespace CGAL { +#define CGAL_BASE \ + Cartesian_base_d::Type, Dim> +template +struct Epeck_d_help1 +: CGAL_BASE +{ + CGAL_CONSTEXPR Epeck_d_help1(){} + CGAL_CONSTEXPR Epeck_d_help1(int d):CGAL_BASE(d){} +}; +#undef CGAL_BASE +#define CGAL_BASE \ + Kernel_d_interface< \ + Cartesian_wrap< \ + Epeck_d_help1, \ + Epeck_d > > +template +struct Epeck_d +: CGAL_BASE +{ + CGAL_CONSTEXPR Epeck_d(){} + CGAL_CONSTEXPR Epeck_d(int d):CGAL_BASE(d){} +}; +#undef CGAL_BASE +} +#endif diff --git a/NewKernel_d/include/CGAL/Epick_d.h b/NewKernel_d/include/CGAL/Epick_d.h index cf780630c56..6443853994f 100644 --- a/NewKernel_d/include/CGAL/Epick_d.h +++ b/NewKernel_d/include/CGAL/Epick_d.h @@ -26,6 +26,7 @@ #include #include #include +#include namespace CGAL { diff --git a/NewKernel_d/include/CGAL/NewKernel_d/Cartesian_LA_base.h b/NewKernel_d/include/CGAL/NewKernel_d/Cartesian_LA_base.h index 00755eab878..c13a9801041 100644 --- a/NewKernel_d/include/CGAL/NewKernel_d/Cartesian_LA_base.h +++ b/NewKernel_d/include/CGAL/NewKernel_d/Cartesian_LA_base.h @@ -86,6 +86,7 @@ struct Cartesian_LA_base_d : public Dimension_base ::add::type ::add::type ::add::type + ::add::type Object_list; typedef typeset< Point_cartesian_const_iterator_tag>::type diff --git a/NewKernel_d/include/CGAL/NewKernel_d/Cartesian_filter_K.h b/NewKernel_d/include/CGAL/NewKernel_d/Cartesian_filter_K.h index 403a0e4d68c..179e97bf016 100644 --- a/NewKernel_d/include/CGAL/NewKernel_d/Cartesian_filter_K.h +++ b/NewKernel_d/include/CGAL/NewKernel_d/Cartesian_filter_K.h @@ -30,7 +30,7 @@ namespace CGAL { template < typename Base_, typename AK_, typename EK_ > struct Cartesian_filter_K : public Base_, - private Store_kernel, private Store_kernel2 + private Store_kernel, private Store_kernel2 { CGAL_CONSTEXPR Cartesian_filter_K(){} CGAL_CONSTEXPR Cartesian_filter_K(int d):Base_(d){} @@ -61,11 +61,11 @@ struct Cartesian_filter_K : public Base_, template struct Type : Get_type {}; template::type> struct Functor : - Inherit_functor {}; + Inherit_functor {}; template struct Functor { - typedef typename Get_functor::type AP; - typedef typename Get_functor::type EP; - typedef Filtered_predicate2 type; + typedef typename Get_functor::type AP; + typedef typename Get_functor::type EP; + typedef Filtered_predicate2 type; }; // TODO: // template struct Functor : diff --git a/NewKernel_d/include/CGAL/NewKernel_d/Coaffine.h b/NewKernel_d/include/CGAL/NewKernel_d/Coaffine.h index 1c2efe462d8..43015d24fef 100644 --- a/NewKernel_d/include/CGAL/NewKernel_d/Coaffine.h +++ b/NewKernel_d/include/CGAL/NewKernel_d/Coaffine.h @@ -89,6 +89,7 @@ template struct Construct_flat_orientation : private Store_kernel std::vector& rest=o.rest; rest.reserve(dim+1); for(int i=0; i),(Point_tag),()); +CGAL_KD_DEFAULT_FUNCTOR(Construct_ttag,(CartesianDKernelFunctors::Construct_weighted_point),(Weighted_point_tag,Point_tag),()); +CGAL_KD_DEFAULT_FUNCTOR(Point_drop_weight_tag,(CartesianDKernelFunctors::Point_drop_weight),(Weighted_point_tag,Point_tag),()); +CGAL_KD_DEFAULT_FUNCTOR(Point_weight_tag,(CartesianDKernelFunctors::Point_weight),(Weighted_point_tag,Point_tag),()); +CGAL_KD_DEFAULT_FUNCTOR(Power_side_of_power_sphere_tag,(CartesianDKernelFunctors::Power_side_of_power_sphere),(Weighted_point_tag),(Power_side_of_power_sphere_raw_tag,Point_drop_weight_tag,Point_weight_tag)); +CGAL_KD_DEFAULT_FUNCTOR(In_flat_power_side_of_power_sphere_tag,(CartesianDKernelFunctors::In_flat_power_side_of_power_sphere),(Weighted_point_tag),(In_flat_power_side_of_power_sphere_raw_tag,Point_drop_weight_tag,Point_weight_tag)); +CGAL_KD_DEFAULT_FUNCTOR(Power_distance_tag,(CartesianDKernelFunctors::Power_distance),(Weighted_point_tag,Point_tag),(Squared_distance_tag,Point_drop_weight_tag,Point_weight_tag)); +CGAL_KD_DEFAULT_FUNCTOR(Power_distance_to_point_tag,(CartesianDKernelFunctors::Power_distance_to_point),(Weighted_point_tag,Point_tag),(Squared_distance_tag,Point_drop_weight_tag,Point_weight_tag)); +CGAL_KD_DEFAULT_FUNCTOR(Power_center_tag,(CartesianDKernelFunctors::Power_center),(Weighted_point_tag,Point_tag),(Compute_point_cartesian_coordinate_tag,Construct_ttag,Construct_ttag,Point_dimension_tag,Squared_distance_to_origin_tag,Point_drop_weight_tag,Point_weight_tag,Power_distance_to_point_tag)); +} // namespace CGAL +#endif diff --git a/NewKernel_d/include/CGAL/NewKernel_d/Wrapper/Cartesian_wrap.h b/NewKernel_d/include/CGAL/NewKernel_d/Wrapper/Cartesian_wrap.h index e1f4610d9ea..44e9aa9621b 100644 --- a/NewKernel_d/include/CGAL/NewKernel_d/Wrapper/Cartesian_wrap.h +++ b/NewKernel_d/include/CGAL/NewKernel_d/Wrapper/Cartesian_wrap.h @@ -33,6 +33,7 @@ #include #include #include +#include #include @@ -111,6 +112,7 @@ CGAL_REGISTER_OBJECT_WRAPPER(Vector); CGAL_REGISTER_OBJECT_WRAPPER(Segment); CGAL_REGISTER_OBJECT_WRAPPER(Sphere); CGAL_REGISTER_OBJECT_WRAPPER(Hyperplane); +CGAL_REGISTER_OBJECT_WRAPPER(Weighted_point); #undef CGAL_REGISTER_OBJECT_WRAPPER // Note: this tends to be an all or nothing thing currently, wrapping diff --git a/NewKernel_d/include/CGAL/NewKernel_d/Wrapper/Weighted_point_d.h b/NewKernel_d/include/CGAL/NewKernel_d/Wrapper/Weighted_point_d.h new file mode 100644 index 00000000000..877eea21b0f --- /dev/null +++ b/NewKernel_d/include/CGAL/NewKernel_d/Wrapper/Weighted_point_d.h @@ -0,0 +1,129 @@ +// Copyright (c) 2014 +// INRIA Saclay-Ile de France (France) +// +// This file is part of CGAL (www.cgal.org); you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public License as +// published by the Free Software Foundation; either version 3 of the License, +// or (at your option) any later version. +// +// Licensees holding a valid commercial license may use this file in +// accordance with the commercial license agreement provided with the software. +// +// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE +// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE. +// +// $URL$ +// $Id$ +// +// Author(s) : Marc Glisse + +#ifndef CGAL_WRAPPER_WEIGHTED_POINT_D_H +#define CGAL_WRAPPER_WEIGHTED_POINT_D_H + +#include +#include +#include +#include +#include +#ifndef CGAL_CXX11 +#include +#endif +#include + +namespace CGAL { +namespace Wrap { + +template +class Weighted_point_d : public Get_type::type +{ + typedef typename Get_type::type FT_; + typedef typename R_::Kernel_base Kbase; + typedef typename Get_type::type Point_; + typedef typename Get_functor >::type CWPBase; + typedef typename Get_functor::type PDWBase; + typedef typename Get_functor::type PWBase; + + typedef Weighted_point_d Self; + BOOST_STATIC_ASSERT((boost::is_same::type>::value)); + +public: + + typedef Tag_true Is_wrapper; + typedef typename R_::Default_ambient_dimension Ambient_dimension; + typedef Dimension_tag<0> Feature_dimension; + + typedef typename Get_type::type Rep; + + const Rep& rep() const + { + return *this; + } + + Rep& rep() + { + return *this; + } + + typedef R_ R; + +#ifdef CGAL_CXX11 + template::type...>,std::tuple >::value>::type> explicit Weighted_point_d(U&&...u) + : Rep(CWPBase()(std::forward(u)...)){} + +// // called from Construct_point_d +// template explicit Point_d(Eval_functor&&,U&&...u) +// : Rep(Eval_functor(), std::forward(u)...){} + template explicit Weighted_point_d(Eval_functor&&,F&&f,U&&...u) + : Rep(std::forward(f)(std::forward(u)...)){} + +#if 0 + // the new standard may make this necessary + Point_d(Point_d const&)=default; + Point_d(Point_d &);//=default; + Point_d(Point_d &&)=default; +#endif + + // try not to use these + Weighted_point_d(Rep const& v) : Rep(v) {} + Weighted_point_d(Rep& v) : Rep(static_cast(v)) {} + Weighted_point_d(Rep&& v) : Rep(std::move(v)) {} + +#else + + Weighted_point_d() : Rep(CWPBase()()) {} + + Weighted_point_d(Rep const& v) : Rep(v) {} // try not to use it + +#define CGAL_CODE(Z,N,_) template \ + explicit Weighted_point_d(BOOST_PP_ENUM_BINARY_PARAMS(N,T,const&t)) \ + : Rep(CWPBase()( \ + BOOST_PP_ENUM_PARAMS(N,t))) {} \ + \ + template \ + Weighted_point_d(Eval_functor,F const& f,BOOST_PP_ENUM_BINARY_PARAMS(N,T,const&t)) \ + : Rep(f(BOOST_PP_ENUM_PARAMS(N,t))) {} + /* + template \ + Point_d(Eval_functor,BOOST_PP_ENUM_BINARY_PARAMS(N,T,const&t)) \ + : Rep(Eval_functor(), BOOST_PP_ENUM_PARAMS(N,t)) {} + */ + + BOOST_PP_REPEAT_FROM_TO(1,11,CGAL_CODE,_) +#undef CGAL_CODE + +#endif + + //TODO: use references? + Point_ point()const{ + return Point_(Eval_functor(),PDWBase(),rep()); + } + FT_ weight()const{ + return PWBase()(rep()); + } + +}; + +} //namespace Wrap +} //namespace CGAL + +#endif // CGAL_WRAPPER_SPHERE_D_H diff --git a/NewKernel_d/include/CGAL/NewKernel_d/function_objects_cartesian.h b/NewKernel_d/include/CGAL/NewKernel_d/function_objects_cartesian.h index a969916ba64..66ebfb9e6c2 100644 --- a/NewKernel_d/include/CGAL/NewKernel_d/function_objects_cartesian.h +++ b/NewKernel_d/include/CGAL/NewKernel_d/function_objects_cartesian.h @@ -554,6 +554,60 @@ template struct Orientation : private Store_kernel { } #endif +namespace CartesianDKernelFunctors { +template struct Power_side_of_power_sphere_raw : private Store_kernel { + CGAL_FUNCTOR_INIT_STORE(Power_side_of_power_sphere_raw) + typedef R_ R; + typedef typename Get_type::type RT; + typedef typename Get_type::type FT; + typedef typename Get_type::type Point; + typedef typename Get_type::type result_type; + typedef typename Increment_dimension::type D1; + typedef typename Increment_dimension::type D2; + typedef typename R::LA::template Rebind_dimension::Other LA; + typedef typename LA::Square_matrix Matrix; + + template + result_type operator()(IterP f, IterP const& e, IterW fw, Pt const& p0, Wt const& w0) const { + typedef typename Get_functor::type Sqdo; + typename Get_functor::type c(this->kernel()); + typename Get_functor::type pd(this->kernel()); + + int d=pd(p0); + Matrix m(d+1,d+1); + if(CGAL::Is_stored::value) { + Sqdo sqdo(this->kernel()); + FT const& h0 = sqdo(p0) - w0; + for(int i=0;f!=e;++f,++fw,++i) { + Point const& p=*f; + for(int j=0;j),(Point_tag),(Point_dimension_tag,Squared_distance_to_origin_tag,Compute_point_cartesian_coordinate_tag)); + +// TODO: make Side_of_oriented_sphere call Power_side_of_power_sphere_raw namespace CartesianDKernelFunctors { template struct Side_of_oriented_sphere : private Store_kernel { CGAL_FUNCTOR_INIT_STORE(Side_of_oriented_sphere) diff --git a/NewKernel_d/include/CGAL/NewKernel_d/functor_tags.h b/NewKernel_d/include/CGAL/NewKernel_d/functor_tags.h index a184f453ba3..b8e17886e17 100644 --- a/NewKernel_d/include/CGAL/NewKernel_d/functor_tags.h +++ b/NewKernel_d/include/CGAL/NewKernel_d/functor_tags.h @@ -172,6 +172,7 @@ namespace CGAL { CGAL_DECL_OBJ(Iso_box, Object); CGAL_DECL_OBJ(Bbox, Object); CGAL_DECL_OBJ(Aff_transformation, Object); + CGAL_DECL_OBJ(Weighted_point, Object); #undef CGAL_DECL_OBJ_ #undef CGAL_DECL_OBJ @@ -214,6 +215,9 @@ namespace CGAL { CGAL_DECL_COMPUTE(Scalar_product); CGAL_DECL_COMPUTE(Hyperplane_translation); CGAL_DECL_COMPUTE(Value_at); + CGAL_DECL_COMPUTE(Point_weight); + CGAL_DECL_COMPUTE(Power_distance); + CGAL_DECL_COMPUTE(Power_distance_to_point); #undef CGAL_DECL_COMPUTE #define CGAL_DECL_ITER_OBJ(X,Y,Z,C) struct X##_tag {}; \ @@ -266,6 +270,8 @@ namespace CGAL { CGAL_DECL_CONSTRUCT(Construct_min_vertex,Point); CGAL_DECL_CONSTRUCT(Construct_max_vertex,Point); CGAL_DECL_CONSTRUCT(Construct_circumcenter,Point); + CGAL_DECL_CONSTRUCT(Point_drop_weight,Point); + CGAL_DECL_CONSTRUCT(Power_center,Weighted_point); #undef CGAL_DECL_CONSTRUCT #if 0 #define CGAL_DECL_ITER_CONSTRUCT(X,Y) struct X##_tag {}; \ @@ -306,6 +312,10 @@ namespace CGAL { CGAL_DECL_PREDICATE(Affinely_independent); CGAL_DECL_PREDICATE(Contained_in_linear_hull); CGAL_DECL_PREDICATE(Contained_in_simplex); + CGAL_DECL_PREDICATE(Power_side_of_power_sphere_raw); + CGAL_DECL_PREDICATE(Power_side_of_power_sphere); + CGAL_DECL_PREDICATE(In_flat_power_side_of_power_sphere_raw); + CGAL_DECL_PREDICATE(In_flat_power_side_of_power_sphere); #undef CGAL_DECL_PREDICATE #define CGAL_DECL_MISC(X) struct X##_tag {}; \ diff --git a/NewKernel_d/include/CGAL/transforming_iterator.h b/NewKernel_d/include/CGAL/transforming_iterator.h index 186901bc086..15ea19a59ee 100644 --- a/NewKernel_d/include/CGAL/transforming_iterator.h +++ b/NewKernel_d/include/CGAL/transforming_iterator.h @@ -22,6 +22,10 @@ #include #include #include +#include +#include +#include +#include #include #include @@ -54,23 +58,31 @@ template struct Functor_as_base : public T { template class transforming_iterator_helper { + typedef std::iterator_traits Iter_traits; + typedef typename Iter_traits::reference Iter_ref; typedef typename Default::Get()(std::declval::reference>())) + decltype(std::declval()(std::declval())) #else - typename boost::result_of::value_type)>::type + typename boost::result_of::type // should be reference instead of value_type #endif - >::type reference; + >::type reference_; - typedef typename Default::Get::type>::type>::type value_type; + typedef typename Default::Get::type>::type>::type value_type; + + // Crappy heuristic. If we have *it that returns a Weighted_point and F that returns a reference to the Point contained in the Weighted_point it takes as argument, we do NOT want the transformed iterator to return a reference to the temporary *it. On the other hand, if *it returns an int n, and F returns a reference to array[n] it is not so good to lose the reference. This probably should be done elsewhere and should at least be made optional... + typedef typename boost::mpl::if_< + boost::mpl::or_, + boost::is_integral >, + reference_, value_type>::type reference; public: typedef boost::iterator_adaptor< Derived, Iter, value_type, - typename std::iterator_traits::iterator_category, + typename Iter_traits::iterator_category, reference > type; }; diff --git a/NewKernel_d/test/NewKernel_d/Epick_d.cpp b/NewKernel_d/test/NewKernel_d/Epick_d.cpp index 3bbe73dfb02..bed06b20a03 100644 --- a/NewKernel_d/test/NewKernel_d/Epick_d.cpp +++ b/NewKernel_d/test/NewKernel_d/Epick_d.cpp @@ -24,6 +24,8 @@ int main() #include #include #include +#include +#include //typedef CGAL::Cartesian_base_d > K0; //typedef CGAL::Cartesian_base_d > KA; @@ -84,6 +86,7 @@ void test2(){ typedef typename K1::Ray_d R; typedef typename K1::Iso_box_d IB; typedef typename K1::Flat_orientation_d FO; + typedef typename K1::Weighted_point_d WP; //typedef K1::Construct_point CP; typedef typename K1::Construct_point_d CP; @@ -135,6 +138,14 @@ void test2(){ typedef typename K1::Construct_min_vertex_d CmV; typedef typename K1::Construct_max_vertex_d CMV; typedef typename K1::Compute_squared_radius_d SR; + typedef typename K1::Translated_point_d TP; + typedef typename K1::Power_center_d PC; + typedef typename K1::Power_distance_d PoD; + typedef typename K1::Weighted_point_d WP; + typedef typename K1::Construct_weighted_point_d CWP; + //typedef typename K1::Point_drop_weight_d PDW; + typedef CP PDW; + typedef typename K1::Compute_weight_d PW; CGAL_USE_TYPE(AT); CGAL_USE_TYPE(D); @@ -196,6 +207,13 @@ void test2(){ CmV cmv Kinit(construct_min_vertex_d_object); CMV cMv Kinit(construct_max_vertex_d_object); SR sr Kinit(compute_squared_radius_d_object); + TP tp Kinit(translated_point_d_object); + PC pc Kinit(power_center_d_object); + CWP cwp Kinit(construct_weighted_point_d_object); + //PDW pdw Kinit(point_drop_weight_d_object); + PDW const& pdw = cp; + PW pw Kinit(compute_weight_d_object); + PoD pod Kinit(power_distance_d_object); CGAL_USE(bc); CGAL_USE(pol); @@ -203,11 +221,12 @@ void test2(){ CGAL_USE(cd); CGAL_USE(cli); CGAL_USE(cr); + using std::abs; P a=cp(3,4); assert(pd(a)==2); assert(pv(a)[1]==4); P b=vp(cv(5,6,7)); - assert(fabs(b[0]-5./7)<.0001); + assert(abs(b[0]-5./7)<.0001); assert(lc(b,a,1)); assert(!lc(a,b,0)); int rr[]={3,5,2}; @@ -221,8 +240,8 @@ void test2(){ assert(cl(a,c)==CGAL::SMALLER); assert(ll(b,c)); assert(cl(c,b)==CGAL::LARGER); - assert(fabs(m(a,c)[0]-3)<.0001); - assert(fabs(m(a,c)[1]-4.5)<.0001); + assert(abs(m(a,c)[0]-3)<.0001); + assert(abs(m(a,c)[1]-4.5)<.0001); P d=cp(r,r+3,CGAL::Homogeneous_tag()); S s=cs(c,d); std::cout << cc(a,1) << std::endl; @@ -277,9 +296,9 @@ void test2(){ assert(v.size()==1); assert(lr(tab3+0,tab3+2)==1); H h=ch(tab2+1,tab2+3,tab2[0]); - assert(fabs(va(h,x2)-1)<.0001); - assert(fabs(va(h,x3)-1)<.0001); - assert(fabs(va(h,x1)+1)<.0001); + assert(abs(va(h,x2)-1)<.0001); + assert(abs(va(h,x3)-1)<.0001); + assert(abs(va(h,x1)+1)<.0001); H h2=ch(tab2+1,tab2+3,x1,CGAL::ON_POSITIVE_SIDE); assert(hops(h2,x1)); assert(os(h2,x1)==CGAL::ON_POSITIVE_SIDE); @@ -340,20 +359,35 @@ void test2(){ Sp sp = csp(tabz+0,tabz+3); P cent0=cos(sp); P cent1=cos(tabz+0,tabz+3); - assert(fabs(cent0[0]-2)<.0001); - assert(fabs(cent0[1]+3)<.0001); - assert(fabs(cent1[0]-2)<.0001); - assert(fabs(cent1[1]+3)<.0001); - assert(fabs(sp.squared_radius()-25)<.0001); + assert(abs(cent0[0]-2)<.0001); + assert(abs(cent0[1]+3)<.0001); + assert(abs(cent1[0]-2)<.0001); + assert(abs(cent1[1]+3)<.0001); + assert(abs(sp.squared_radius()-25)<.0001); +#if 1 + // Fails for an exact kernel P psp0=ps(sp,0); P psp1=ps(sp,1); P psp2=ps(sp,2); assert(!ed(psp0,psp1)); assert(!ed(psp0,psp2)); assert(!ed(psp2,psp1)); - assert(fabs(sd(cent0,psp0)-25)<.0001); - assert(fabs(sd(cent0,psp1)-25)<.0001); - assert(fabs(sd(cent0,psp2)-25)<.0001); + assert(abs(sd(cent0,psp0)-25)<.0001); + assert(abs(sd(cent0,psp1)-25)<.0001); + assert(abs(sd(cent0,psp2)-25)<.0001); +#endif + P x2py1 = tp(x2,y1); + assert(x2py1[1]==-2); + WP tw[]={cwp(cp(5,0),1.5),cwp(cp(2,std::sqrt(3)),1),cwp(cp(2,-std::sqrt(3)),1)}; + WP xw=pc(tw+0,tw+3); + assert(abs(pod(xw,tw[0]))<.0001); + assert(abs(pod(xw,tw[1]))<.0001); + assert(abs(pod(xw,tw[2]))<.0001); + assert(pdw(xw)[0]<2.95); + assert(pdw(xw)[0]>2.5); + assert(pw(xw)<2.95); + assert(pw(xw)>2.5); + P tl=cp(2,5); P br=cp(4,-1); @@ -459,6 +493,7 @@ void test3(){ PD pd Kinit(point_dimension_d_object); AI ai Kinit(affinely_independent_d_object); SBDS sbds Kinit(side_of_bounded_sphere_d_object); + using std::abs; P a; // Triangulation needs this :-( a=cp(2,3,4); assert(pd(a)==3); @@ -482,7 +517,7 @@ void test3(){ std::cout << *i << ' '; std::cout << '\n'; P e=cp(-2,3,0); - assert(fabs(sd(e,a)-32)<.0001); + assert(abs(sd(e,a)-32)<.0001); P tab[]={a,b,c,d,e}; std::cout << po (&tab[0],tab+4) << ' '; std::cout << sos(&tab[0],tab+5) << ' '; @@ -535,6 +570,7 @@ void test3(){ P x4=cp(0,0,1); P x5=cp(0,0,0); P x6=cp(0,0,-1); + assert(!ed(x1,x2)); P tab2[]={x1,x2,x3,x4,x5}; assert(cis(tab2+0,tab2+4,x5)); assert(po(tab2+0,tab2+4)==CGAL::POSITIVE); @@ -592,6 +628,28 @@ void test3(){ assert(sbds(t1+0,t1+2,cp(2,2,3.415)) == CGAL::ON_UNBOUNDED_SIDE); assert(sbds(t1+0,t1+3,cp(2.1,3.5,1.9)) == CGAL::ON_BOUNDED_SIDE); assert(sbds(t1+0,t1+3,cp(10,10,10)) == CGAL::ON_UNBOUNDED_SIDE); + + typedef typename K1::Weighted_point_d WP; + typedef typename K1::Construct_weighted_point_d CWP; + //typedef typename K1::Point_drop_weight_d PDW; + typedef CP_ PDW; + typedef typename K1::Compute_weight_d PW; + typedef typename K1::Power_side_of_power_sphere_d PT; + typedef typename K1::In_flat_power_side_of_power_sphere_d IFPT; + CWP cwp Kinit(construct_weighted_point_d_object); + //PDW pdw Kinit(point_drop_weight_d_object); + PDW const& pdw = cp_; + PW pw Kinit(compute_weight_d_object); + PT pt Kinit(power_side_of_power_sphere_d_object); + IFPT ifpt Kinit(in_flat_power_side_of_power_sphere_d_object); + WP wp; + wp = cwp (x1, 2); + WP xw6 = cwp (x6, 0); + assert (pw(wp) == 2); + assert (ed(pdw(wp), x1)); + WP tabw[]={cwp(x1,0),cwp(x2,0),cwp(x3,0),cwp(x4,0),cwp(x5,0)}; + assert(pt(tabw+0,tabw+4,tabw[4])==CGAL::ON_POSITIVE_SIDE); + assert(ifpt(fo4,tabw+0,tabw+3,xw6)==CGAL::ON_POSITIVE_SIDE); } template struct CGAL::Epick_d >; template struct CGAL::Epick_d >; diff --git a/Number_types/include/CGAL/Gmpq.h b/Number_types/include/CGAL/Gmpq.h index 2fde983a7b3..c98532b3994 100644 --- a/Number_types/include/CGAL/Gmpq.h +++ b/Number_types/include/CGAL/Gmpq.h @@ -156,6 +156,17 @@ namespace Eigen { MulCost = 100 }; }; + + namespace internal { + template<> + struct significant_decimals_impl + { + static inline int run() + { + return 0; + } + }; + } } //since types are included by Gmp_coercion_traits.h: diff --git a/Number_types/include/CGAL/Interval_nt.h b/Number_types/include/CGAL/Interval_nt.h index e0605802ca5..e0969570391 100644 --- a/Number_types/include/CGAL/Interval_nt.h +++ b/Number_types/include/CGAL/Interval_nt.h @@ -1284,6 +1284,13 @@ namespace Eigen { MulCost = 10 }; }; + + namespace internal { + template struct significant_decimals_impl; + template + struct significant_decimals_impl > + : significant_decimals_impl::value_type> { }; + } } #endif // CGAL_INTERVAL_NT_H diff --git a/Spatial_searching/doc/Spatial_searching/Concepts/RangeSearchTraits.h b/Spatial_searching/doc/Spatial_searching/Concepts/RangeSearchTraits.h index 614451d084e..56281b811fc 100644 --- a/Spatial_searching/doc/Spatial_searching/Concepts/RangeSearchTraits.h +++ b/Spatial_searching/doc/Spatial_searching/Concepts/RangeSearchTraits.h @@ -10,6 +10,7 @@ range search queries in a model of `SpatialTree`. \cgalHasModel `CGAL::Cartesian_d` \cgalHasModel `CGAL::Homogeneous_d` +\cgalHasModel `CGAL::Epick_d` \cgalHasModel `CGAL::Search_traits_2` \cgalHasModel `CGAL::Search_traits_3` diff --git a/Spatial_searching/doc/Spatial_searching/Concepts/SearchTraits.h b/Spatial_searching/doc/Spatial_searching/Concepts/SearchTraits.h index 633e5b006dd..2abebf75c0b 100644 --- a/Spatial_searching/doc/Spatial_searching/Concepts/SearchTraits.h +++ b/Spatial_searching/doc/Spatial_searching/Concepts/SearchTraits.h @@ -6,7 +6,8 @@ The concept `SearchTraits` defines the requirements for the template parameter of the search classes. \cgalHasModel `CGAL::Cartesian_d` -\cgalHasModel `CGAL::Homogeneous_d` +\cgalHasModel `CGAL::Homogeneous_d` +\cgalHasModel `CGAL::Epick_d` \cgalHasModel `CGAL::Search_traits_2` \cgalHasModel `CGAL::Search_traits_3` \cgalHasModel `CGAL::Search_traits_d` diff --git a/Spatial_searching/include/CGAL/Search_traits_adapter.h b/Spatial_searching/include/CGAL/Search_traits_adapter.h index dc06e1c8080..428e8ff0d94 100644 --- a/Spatial_searching/include/CGAL/Search_traits_adapter.h +++ b/Spatial_searching/include/CGAL/Search_traits_adapter.h @@ -105,6 +105,17 @@ public: typename Base_traits::Cartesian_const_iterator_d operator()(const Point_with_info& p, int) const { return Base::operator() (get(ppmap,p),0); } + + // These 2 additional operators forward the call to Base_traits. + // This is needed because of an undocumented requirement of + // Orthogonal_k_neighbor_search and Orthogonal_incremental_neighbor_search: + // Traits::Construct_cartesian_const_iterator should be callable + // on the query point type + typename Base_traits::Cartesian_const_iterator_d operator()(const typename Base_traits::Point_d& p) const + { return Base::operator() (p); } + + typename Base_traits::Cartesian_const_iterator_d operator()(const typename Base_traits::Point_d& p, int) const + { return Base::operator() (p,0); } }; struct Construct_iso_box_d: public Base::Construct_iso_box_d{ diff --git a/Triangulation/TODO b/Triangulation/TODO index e7c7b34e5a5..fca72fad682 100644 --- a/Triangulation/TODO +++ b/Triangulation/TODO @@ -1,258 +1,3 @@ - --------------------------------------------------- -Problems to be solved from Sam's reading (september 2012) --------------------------------------------------- - -*) Substitute iterator : it seems to be comparing the point directly, instead -of comparing the iterator. That look like a big performance killer. - -*) I have marked with 'FIXME' all the lines where the code assumes, or seems to -assume that the vertex at infinity is at index 0 in the full cell. The code is -still working because the vertex at infinity is indeed at index 0, but this is -no more a requirement from the documentation, so that the code (Tr.h and -Delaunay_tr.h) should be changed. OR, we re-document 0 as the index of the -vertex at infinity. - - - -ALL Remarks below are done or not important - - --------------------------------------------------- -Problems to be solved from the reviews (beginning 2012) --------------------------------------------------- - -example delaunay does not execute properly - - SAM: seems to be a compiler bug --> low priority. - OD: I am afraid that it is still some strange bug in the code that show up only in - some compilation circumstances. - running both delaunay compiled debug/release, the diverges after the 6th insertion (in 2d round) - -in Triangulation_data_structure : -put a default value for dim in the constructor -(does not work, I do not understand why). OK done. - - SAM: I've modified the Concept constructor's documentation to reflect this. - If it is not satisfying, we might remove this doc. from the Concept and - move it to the class documentation. - -check that the perturbation scheme is independant of the order of insertion - - SAM: It is your and Monique's scheme. It should be independent, which - is a requirement anyway for the Delaunay::remove() function to work - properly. - OD: seems ok, lexicographic order. - -Add a template parameter Location_policy - -ambient dim vs max dim - OD: global replace done, remains to get the dim from the traits (and not from the point of the traits) - -check all is_valid function, precise in the doc what they are doing. -(do sth like 2/3 d) - - SAM: I suggest that the documentation of the function in the concepts only - states that "any validity check can be performed here (see model doc. for - details)", and that we re document the same function in the class documentation - with details on what exactly is performed by the implementation. - OD: seems reasonable - SAM: I've changed my mind a little. We need some easy mandatory checks in - concepts T...DSFullCell and T...DSVertex so that we can rely on these tests in - the implementation of T.._data_structure, instead of having to re-implement them. - So : the concepts lists simple validity checks taht must be present, and refers - to the documentation of the models for possible additional validity checks - that are implemented. - --> This scheme is done (code&doc) for TDS TDSFullCell and TDSVertex and Triangulation. - --> ALL DONE. - -small feature with iterator "all tuples" - -make the code and doc agree on Flat_* stuff (orientation in a flat) - -iterator on points in concept TriangulationVertex should be removed to keep requirements minimal - - REMOVED: But this makes the concept TriangulationFullCell empty and useless. - -make doc of TDS-FullCellStoragePolicy in user manual - - DONE. - -missing figures : Triangulation/fig/detail.png Triangulation/fig/illustration.png - - - - - - - - - --------------------------------------------------- -Old todo list (beginning 2011) --------------------------------------------------- -Les premiers points meritent d'etre examinés - -Je me replonge dans Triangulation en dim d - -ma liste de trucs à faire ou questions à résoudre: - - ---- Constructeur de triangulation - prévoir un constructeur qui prends d+1 points en position générale - ---- Doc 1.5 à faire - - ---- mirror_index - - policy a ajouter dans le user manual - - le code ne fait pas ce qui est dit dans la doc (ça ne retourne pas -1) - - PS: Samuel tu as benché la policy ? ça vaut vraiment la peine de s'embeter avec ça ? - - ---- Face - - pourquoi une face ne pourrait pas être de pleine dimension ? - c'est pas la peine d'avoir un truc générique pour les faces de toutes dim sauf un - OD: J'ai propose de viré cette restriction dans la doc, mais il faudrait vérifier le code. - - ---- TriangulationTraits et DelaunayTraits - lower dimensional predicates - le truc de l'orientation "consistante"est quand même un peu délicat - - d'une part, je proposerai un truc du genre: le noyau doit fournir un iterateur - de points en position générale - (un truc fixe, ça donne toujours les mêmes points e.g. l'origine et les points - avec une coordonnée à 1 et les autres à 0). - Ce générateur permet donc de compléter tout ensemble de points pour obtenir - des points en position générale. - - d'autre part ???? il faut un minimum en plus pour qu'on comprenne que - l'orientation du sous-espace est stocké qq part. - ---- TriangulationTraits has models - ... il faut voir ce que l'on met vraiment. - ---- Marque visité. - elle devrait etre dans le concept TDSFullCell - --------------------------------------------------- -ci dessous les points ont deja etes traites - ---- Index du sommet à l'infini. - Il est en ce moment prévu que le sommet à l'infini est l'index 0 dans les full - cell où il apparraît. Ça me semble très discutable, en particulier la tds ne sait - pas qui est à l'infini, et donc comment peut-elle garantir dans une manipulation - que toutes les cellules crées ou modifiées vont préserver cette position. - - je suis donc pour virer ce truc ce qui implique des modifs dans la doc - (1.3 - ? p54) et dans le code. - - en plus il y a des trucs pas cohérent, i.e. p. 50 dans la doc de locate - l'infini est en dernier, pas en premier. - - finite_vertex_iterator marche pas - SAM: J'ai corrigé le problème, mais c'est pas super propre. Le problème - vient de boost::filter_iterator qui suppose que le prédicat de filtrage - prend en argument |const value_type &| plutôt que |iterator|. Mais - quand on veut tester si un sommet est infini, on a envie de comparer - les iterator (ou handle)... ce qui n'est pas possible avec le - filter_iterator de boost. (Le filter_iterator de CGAL appele le - prédicat de filtrage sur l'itérateur, et non pas sur la valeur vers - lequel il pointe, comme le fait boost. Mais je ne sais pas si le filtre - d'iterateur de CGAL est destiné à durer) - - ---- Orientation des cellules. - les cellules sont orientées positivement par convention - - le dire dans la doc (deux full cell partagent une facet, elles doivent avoir - des orientation complémentaires - - dans Triangulation, il faut virer la méthode "orientation" - (ou au minimum, la mettre advanced) - ---- Triangulation is_infinite - j'avais viré les préconditions, et je maintiens. - en dim 0 on a une full cell finie et une infinie - en dim 1 on a 2 full cell infinie, 1 vertex/facet infini - et ok, en d=0, facet n'a pas de sens, mais l'utilisateur aura du mal a en avoir une - ---- Triangulation fonctions incident_* et star - est-il nécessaire de reprendre leur doc? (on pourrait pointer sur celle de TDS) - au moins pour incident_upper_faces - ---- dans Triangulation_data_structure.h insert_in_tagged_hole - on a l'air de supposer que full_cell(f) est marqué visité - (pourquoi ça serait pas l'autre représentation de la Facet f ?) - la doc de TDS::insert_in_hole est maintenant claire - - - - - - - -Older todo list (2010 ?) -__________________________________________________________________________RENAMING - -*) doc : --------- - -*) code not done : ------------------- - -bring part of Delaunay::remove into TDS (see comments in Delaunay_triangulation.h, and item FUTURE below) - -__________________________________________________________________________ALL - -*) FUTURE: better system for simplex's flags to know if we can use them or not - in Delaunay_triangulation::remove(Vertex_handle) - -_________________________________________________TRIANGULATION_DATA_STRUCTURE - -*) Un-recursify insert_in_tagged_hole : crashed in dimension 8 : stack overflow. - -*) TriangulationDataStructure: - - Should we put >> and << in the documentation of the class - Triangulation_data_structure ? - -*) Triangulation_data_structure.h: - - Ensure that it is topologically possible to collapse the Face in - collapse_face(Face) - -*) TriangulationDSVertex - - Should we put >> and << in the documentation of the class - Triangulation_ds_vertex instead of in that of the concept ? - -*) TDSFullCellStoragePolicy - what is it ? => see the doc. - it is undocumented. => it is documented (Triangulation_ds_full_cell) - if the aim is to choose between "full representation" and "1-skeleton" - I would prefer a different TDS class, at least the doc has to rewriten a sentence - like "class TDS explicitely stores its vertices and full cells" - does not apply to 1-skeleton representation - => no, this has nothing to do with 1-skeleton. - -*) default dim parameter in constructor (especially in the static case) -done - -________________________________________________________________TRIANGULATION - -*) default dim parameter in constructor (especially in the static case) -done - -*) Why vertex at infinity should have index 0 in all cells it appears ? - This is a convention that is enforced throughout the code. Makes it faster - to check for an infinite cell. - The "old" package "Convex_hull_d" does the same. -done - -_______________________________________________________DELAUNAY_TRIANGULATION - -________________________________________________________REGULAR_TRIANGULATION - -*) write Regular_triangulation.h (!) -*) write RegularTriangulationTraits.tex -*) write Regular_triangulation.tex - ------------------------------------------------------- RANDOM DESIGN IDEAS extracted from Convex_hull.h ------------------------------------------------------- @@ -263,4 +8,3 @@ ________________________________________________________REGULAR_TRIANGULATION In this second case, we must keeps the points that are inserted in the hull, as they may become part of the boundary later on, when some points are removed. - Constructor with range argument uses quickhull. -*/ diff --git a/Triangulation/applications/Triangulation/CMakeLists.txt b/Triangulation/applications/Triangulation/CMakeLists.txt new file mode 100644 index 00000000000..1f49e2c545e --- /dev/null +++ b/Triangulation/applications/Triangulation/CMakeLists.txt @@ -0,0 +1,69 @@ +# Created by the script cgal_create_cmake_script_with_options +# This is the CMake script for compiling a set of CGAL applications. + +project( Triangulation_apps ) + + +cmake_minimum_required(VERSION 2.6.2) +if("${CMAKE_MAJOR_VERSION}.${CMAKE_MINOR_VERSION}" VERSION_GREATER 2.6) + if("${CMAKE_MAJOR_VERSION}.${CMAKE_MINOR_VERSION}.${CMAKE_PATCH_VERSION}" VERSION_GREATER 2.8.3) + cmake_policy(VERSION 2.8.4) + else() + cmake_policy(VERSION 2.6) + endif() +endif() + +set( CMAKE_ALLOW_LOOSE_LOOP_CONSTRUCTS true ) + +if ( COMMAND cmake_policy ) + + cmake_policy( SET CMP0003 NEW ) + +endif() + +# CGAL and its components +find_package( CGAL QUIET COMPONENTS ) + +if ( NOT CGAL_FOUND ) + + message(STATUS "This project requires the CGAL library, and will not be compiled.") + return() + +endif() + +# include helper file +include( ${CGAL_USE_FILE} ) + + +# Boost and its components +find_package( Boost REQUIRED ) + +if ( NOT Boost_FOUND ) + + message(STATUS "This project requires the Boost library, and will not be compiled.") + + return() + +endif() + +find_package(Eigen3 3.1.0) +if (EIGEN3_FOUND) + include( ${EIGEN3_USE_FILE} ) +endif() + +# include for local directory +include_directories( BEFORE include ) + +# include for local package +include_directories( BEFORE ../../include ) + + +# Creating entries for all .cpp/.C files with "main" routine +# ########################################################## + +include( CGAL_CreateSingleSourceCGALProgram ) + +create_single_source_cgal_program( "points_to_RT_to_off.cpp" ) +create_single_source_cgal_program( "points_to_DT_to_off.cpp" ) + + diff --git a/Triangulation/applications/Triangulation/data/points_2 - extreme weights.cin b/Triangulation/applications/Triangulation/data/points_2 - extreme weights.cin new file mode 100644 index 00000000000..19b563f4835 --- /dev/null +++ b/Triangulation/applications/Triangulation/data/points_2 - extreme weights.cin @@ -0,0 +1,11 @@ +2 +0.0071 1.6899 0 +0.3272 1.3694 0.05 +1.3697 1.8296 0.1 +0.6722 0.3012 0.15 +1.1726 0.1899 0.2 +0.4374 2.8541 100.25 +2.5923 0.1904 0.3 +1.3083 2.5462 200.35 +1.4981 1.3929 0.4 +2.1304 2.055 0.45 \ No newline at end of file diff --git a/Triangulation/applications/Triangulation/data/points_2.cin b/Triangulation/applications/Triangulation/data/points_2.cin new file mode 100644 index 00000000000..5a087138dcd --- /dev/null +++ b/Triangulation/applications/Triangulation/data/points_2.cin @@ -0,0 +1,20 @@ +2 +0 0 6.28953 +-2.85086 -0.471442 6.12896 +1.90972 0.101219 0.988689 +0.637771 2.59367 5.80372 +2.22209 0.903198 2.19478 +-0.487202 -2.71506 4.90996 +1.1193 -1.91787 2.99626 +1.54714 0.109831 0 +0.44556 -2.73047 4.48142 +0.427936 1.28495 6.23624 +-2.67212 0.766674 5.29623 +1.5763 -1.59828 2.58905 +-0.476603 2.2546 6.04797 +1.57172 -0.514711 6.11405 +1.84528 2.10139 5.53936 +-2.99827 -0.101677 5.92246 +-0.482122 -2.39584 4.44264 +-2.25558 -1.492 6.23448 +0.128475 -1.75125 3.18916 \ No newline at end of file diff --git a/Triangulation/applications/Triangulation/data/points_3 - extreme weights.cin b/Triangulation/applications/Triangulation/data/points_3 - extreme weights.cin new file mode 100644 index 00000000000..d1a088200f2 --- /dev/null +++ b/Triangulation/applications/Triangulation/data/points_3 - extreme weights.cin @@ -0,0 +1,11 @@ +3 +0.0071 1.6899 2.521 0 +0.3272 1.3694 3.15 0.05 +1.3697 1.8296 2.654 0.1 +-10.6722 0.3012 0.1548 1000.15 +1.1726 0.1899 0.3658 0.2 +0.4374 20.8541 1.45894 2000.25 +2.5923 0.1904 0.6971 0.3 +10.3083 2.5462 1.3658 1000.35 +1.4981 1.3929 2.949 0.4 +2.1304 2.055 0.6597455 1.45 \ No newline at end of file diff --git a/Triangulation/applications/Triangulation/data/points_3 - no weights.cin b/Triangulation/applications/Triangulation/data/points_3 - no weights.cin new file mode 100644 index 00000000000..58e5426c6a9 --- /dev/null +++ b/Triangulation/applications/Triangulation/data/points_3 - no weights.cin @@ -0,0 +1,11 @@ +3 +0.0071 1.6899 2.521 0 +0.3272 1.3694 3.15 0 +1.3697 1.8296 2.654 0 +-10.6722 0.3012 0.1548 0 +1.1726 0.1899 0.3658 0 +0.4374 20.8541 1.45894 0 +2.5923 0.1904 0.6971 0 +10.3083 2.5462 1.3658 0 +1.4981 1.3929 2.949 0 +2.1304 2.055 0.6597455 0 \ No newline at end of file diff --git a/Triangulation/applications/Triangulation/data/points_3.cin b/Triangulation/applications/Triangulation/data/points_3.cin new file mode 100644 index 00000000000..d1a088200f2 --- /dev/null +++ b/Triangulation/applications/Triangulation/data/points_3.cin @@ -0,0 +1,11 @@ +3 +0.0071 1.6899 2.521 0 +0.3272 1.3694 3.15 0.05 +1.3697 1.8296 2.654 0.1 +-10.6722 0.3012 0.1548 1000.15 +1.1726 0.1899 0.3658 0.2 +0.4374 20.8541 1.45894 2000.25 +2.5923 0.1904 0.6971 0.3 +10.3083 2.5462 1.3658 1000.35 +1.4981 1.3929 2.949 0.4 +2.1304 2.055 0.6597455 1.45 \ No newline at end of file diff --git a/Triangulation/applications/Triangulation/points_to_DT_to_off.cpp b/Triangulation/applications/Triangulation/points_to_DT_to_off.cpp new file mode 100644 index 00000000000..6c6da039129 --- /dev/null +++ b/Triangulation/applications/Triangulation/points_to_DT_to_off.cpp @@ -0,0 +1,42 @@ +#include +#include +#include + +#include + +typedef CGAL::Epick_d K; +typedef CGAL::Delaunay_triangulation DT; + +void test(int dim) +{ + std::stringstream input_filename; + input_filename << "data/points_" << dim << ".cin"; + std::ifstream in(input_filename.str()); + + DT::Point p; + std::vector points; + + int dim_from_file; + in >> dim_from_file; + while(in >> p) + points.push_back(p); + + // Build the Regular Triangulation + DT dt(dim_from_file); + dt.insert(points.begin(), points.end()); + CGAL_assertion(dt.is_valid(true)); + + // Export + std::stringstream output_filename; + output_filename << "data/dt_dim" << dim << ".off"; + std::ofstream off_stream(output_filename.str()); + CGAL::export_triangulation_to_off(off_stream, dt); +} + +int main() +{ + //test(2); + //test(3); + test(10); + return 0; +} diff --git a/Triangulation/applications/Triangulation/points_to_RT_to_off.cpp b/Triangulation/applications/Triangulation/points_to_RT_to_off.cpp new file mode 100644 index 00000000000..c882330d392 --- /dev/null +++ b/Triangulation/applications/Triangulation/points_to_RT_to_off.cpp @@ -0,0 +1,41 @@ +#include +#include +#include + +#include + +typedef CGAL::Epick_d K; +typedef CGAL::Regular_triangulation RT; + +void test(int dim) +{ + std::stringstream input_filename; + input_filename << "data/points_" << dim << ".cin"; + std::ifstream in(input_filename.str()); + + RT::Weighted_point wp; + std::vector wpoints; + + int dim_from_file; + in >> dim_from_file; + while(in >> wp) + wpoints.push_back(wp); + + // Build the Regular Triangulation + RT rt(dim_from_file); + rt.insert(wpoints.begin(), wpoints.end()); + CGAL_assertion(rt.is_valid(true)); + + // Export + std::stringstream output_filename; + output_filename << "data/rt_dim" << dim << ".off"; + std::ofstream off_stream(output_filename.str()); + CGAL::export_triangulation_to_off(off_stream, rt); +} + +int main() +{ + test(2); + test(3); + return 0; +} diff --git a/Triangulation/benchmark/Triangulation/CMakeLists.txt b/Triangulation/benchmark/Triangulation/CMakeLists.txt index 9529aed0791..97c0fe9b2db 100644 --- a/Triangulation/benchmark/Triangulation/CMakeLists.txt +++ b/Triangulation/benchmark/Triangulation/CMakeLists.txt @@ -21,6 +21,7 @@ if ( CGAL_FOUND ) include_directories (BEFORE "../../include") include_directories (BEFORE "include") create_single_source_cgal_program( "delaunay.cpp" ) + create_single_source_cgal_program( "Td_vs_T2_and_T3.cpp" ) else() message(STATUS "NOTICE: Some of the executables in this directory need Eigen 3.1 (or greater) and will not be compiled.") diff --git a/Triangulation/benchmark/Triangulation/Td_vs_T2_and_T3.cpp b/Triangulation/benchmark/Triangulation/Td_vs_T2_and_T3.cpp new file mode 100644 index 00000000000..2c4ac1c98a9 --- /dev/null +++ b/Triangulation/benchmark/Triangulation/Td_vs_T2_and_T3.cpp @@ -0,0 +1,267 @@ +// To deactivate statics filters in the 2D/3D case +//#define CGAL_NO_STATIC_FILTERS + +#include +#include +#include + +#include +#include +#include +#include +#include + +#include +#include +#include +#include + +#include +#include +#include +#include +#include + +#include +#include +#include "console_color.h" + +template +struct Stats_getter; + +// T2 specialization +template +struct Stats_getter > +{ + typedef CGAL::Delaunay_triangulation_2 DT; + + Stats_getter(DT const& dt) : m_dt(dt) {} + + std::size_t number_of_vertices() { return m_dt.number_of_vertices(); } + std::size_t number_of_finite_cells() { return m_dt.number_of_faces(); } + + DT m_dt; +}; + +// RT2 specialization +template +struct Stats_getter > +{ + typedef CGAL::Regular_triangulation_2 DT; + + Stats_getter(DT const& dt) : m_dt(dt) {} + + std::size_t number_of_vertices() { return m_dt.number_of_vertices(); } + std::size_t number_of_finite_cells() { return m_dt.number_of_faces(); } + + DT m_dt; +}; + +// T3 specialization +template +struct Stats_getter > +{ + typedef CGAL::Delaunay_triangulation_3 DT; + + Stats_getter(DT const& dt) : m_dt(dt) {} + + std::size_t number_of_vertices() { return m_dt.number_of_vertices(); } + std::size_t number_of_finite_cells() { return m_dt.number_of_finite_cells(); } + + DT m_dt; +}; + +// RT3 specialization +template +struct Stats_getter > +{ + typedef CGAL::Regular_triangulation_3 DT; + + Stats_getter(DT const& dt) : m_dt(dt) {} + + std::size_t number_of_vertices() { return m_dt.number_of_vertices(); } + std::size_t number_of_finite_cells() { return m_dt.number_of_finite_cells(); } + + DT m_dt; +}; + + +template +void test( + int d, int N, Pt_d_range const& points_d, Pt_23_range const& points_23, + std::string const& DTd_static_or_dyn) +{ + // Td + { + DT_d dt(d); + CGAL::Timer timer; + timer.start(); + dt.insert(points_d.begin(), points_d.end()); + + std::cerr << " * Td: " << yellow << timer.time() << " s" + << white << std::endl; + std::cerr << " " << dt.number_of_vertices() << " vertices, " + << dt.number_of_finite_full_cells() << " finite cells." + << std::endl; + } + + // T2 or T3 + { + CGAL::Timer timer; + timer.start(); + + DT_23 dt; + dt.insert(points_23.begin(), points_23.end()); + + std::cerr << " * T" << d << ": " << yellow << timer.time() << " s" + << white << std::endl; + Stats_getter sg(dt); + std::cerr << " " << sg.number_of_vertices() << " vertices, " + << sg.number_of_finite_cells() << " finite cells." + << std::endl; + } +} + +template< int D, typename Dim_tag > +void go(const int N) +{ + CGAL_assertion(D == 2 || D == 3); + + // Generate points (in a common "array" format) + std::vector > coords; + coords.reserve(N); + for (int i = 0; i < N; ++i) + { + CGAL::cpp11::array pt; + for (int j = 0; j < D; ++j) + pt[j] = CGAL::default_random.get_double(-1., 1.); + coords.push_back(pt); + } + // Generate weights + std::vector weights; + weights.reserve(N); + for (int i = 0; i < N; ++i) + weights.push_back(CGAL::default_random.get_double(-10., 10.)); + + // DTd + typedef CGAL::Epick_d Kd; + typedef CGAL::Delaunay_triangulation DT_d; + typedef typename DT_d::Point Point_d; + + std::vector points_d; + points_d.reserve(N); + for (int i = 0; i < N; ++i) + points_d.push_back(Point_d(D, coords[i].begin(), coords[i].end())); + + // RTd + typedef CGAL::Regular_triangulation RT_d; + typedef typename RT_d::Bare_point Bare_point_d; + typedef typename RT_d::Point WPoint_d; + + std::vector wpoints_d; + wpoints_d.reserve(N); + for (int i = 0; i < N; ++i) + { + wpoints_d.push_back(WPoint_d( + Bare_point_d(D, coords[i].begin(), coords[i].end()), + weights[i])); + } + + // T2 or T3 + typedef CGAL::Exact_predicates_inexact_constructions_kernel K23; + if (D == 2) + { + // Delaunay + typedef CGAL::Delaunay_triangulation_2 DT_2; + typedef typename DT_2::Point Point; + + std::vector points; + points.reserve(N); + for (int i = 0; i < N; ++i) + points.push_back(Point(coords[i][0], coords[i][1])); + + std::cerr << std::endl << "DELAUNAY - dim " << D << " - " + << N << " points." << std::endl; + test(D, N, points_d, points, "static"); + + // Regular + typedef CGAL::Regular_triangulation_filtered_traits_2 Traits_2; + typedef CGAL::Regular_triangulation_2 RT_2; + typedef typename RT_2::Bare_point Bare_point; + typedef typename RT_2::Point WPoint; + + std::vector wpoints; + wpoints.reserve(N); + for (int i = 0; i < N; ++i) + { + wpoints.push_back(WPoint( + Bare_point(coords[i][0], coords[i][1]), + weights[i])); + } + + std::cerr << std::endl << "REGULAR - dim " << D << " - " + << N << " points." << std::endl; + test(D, N, wpoints_d, wpoints, "static"); + } + else if (D == 3) + { + typedef CGAL::Delaunay_triangulation_3 DT_3; + typedef typename DT_3::Point Point; + + std::vector points; + points.reserve(N); + for (int i = 0; i < N; ++i) + points.push_back(Point(coords[i][0], coords[i][1], coords[i][2])); + + std::cerr << std::endl << "DELAUNAY - dim " << D << " - " + << N << " points." << std::endl; + test(D, N, points_d, points, "static"); + + // Regular + typedef CGAL::Regular_triangulation_filtered_traits_3 Traits_3; + typedef CGAL::Regular_triangulation_3 RT_3; + typedef typename RT_3::Bare_point Bare_point; + typedef typename RT_3::Point WPoint; + + std::vector wpoints; + wpoints.reserve(N); + for (int i = 0; i < N; ++i) + { + wpoints.push_back(WPoint( + Bare_point(coords[i][0], coords[i][1], coords[i][2]), + weights[i])); + } + + std::cerr << std::endl << "REGULAR - dim " << D << " - " + << N << " points." << std::endl; + test(D, N, wpoints_d, wpoints, "static"); + } +} + +int main(int argc, char **argv) +{ + srand(static_cast(time(NULL))); +#ifdef _DEBUG + int N = 100; +#else + int N = 100000; +#endif + if (argc > 1) N = atoi(argv[1]); + + std::cerr << "-----------------------------------------" << std::endl; + std::cerr << "-- STATIC --" << std::endl; + std::cerr << "-----------------------------------------" << std::endl; + go<2, CGAL::Dimension_tag<2> >(N); + go<3, CGAL::Dimension_tag<3> >(N); + std::cerr << std::endl; + + std::cerr << "-----------------------------------------" << std::endl; + std::cerr << "-- DYNAMIC --" << std::endl; + std::cerr << "-----------------------------------------" << std::endl; + go<2, CGAL::Dynamic_dimension_tag>(N); + go<3, CGAL::Dynamic_dimension_tag>(N); + std::cerr << std::endl; + + return 0; +} diff --git a/Triangulation/benchmark/Triangulation/console_color.h b/Triangulation/benchmark/Triangulation/console_color.h new file mode 100644 index 00000000000..50554178257 --- /dev/null +++ b/Triangulation/benchmark/Triangulation/console_color.h @@ -0,0 +1,68 @@ +#ifndef CONSOLE_COLOR_H_ +#define CONSOLE_COLOR_H_ + +#include + +#if defined(WIN32) +#include +#endif + +inline std::ostream& blue(std::ostream &s) +{ +#if defined(WIN32) + HANDLE hStdout = GetStdHandle(STD_OUTPUT_HANDLE); + SetConsoleTextAttribute(hStdout, + FOREGROUND_BLUE|FOREGROUND_GREEN|FOREGROUND_INTENSITY); +#else + s << "\x1b[0;34m"; +#endif + return s; +} + +inline std::ostream& red(std::ostream &s) +{ +#if defined(WIN32) + HANDLE hStdout = GetStdHandle(STD_OUTPUT_HANDLE); + SetConsoleTextAttribute(hStdout, FOREGROUND_RED|FOREGROUND_INTENSITY); +#else + s << "\x1b[0;31m"; +#endif + return s; +} + +inline std::ostream& green(std::ostream &s) +{ +#if defined(WIN32) + HANDLE hStdout = GetStdHandle(STD_OUTPUT_HANDLE); + SetConsoleTextAttribute(hStdout, FOREGROUND_GREEN|FOREGROUND_INTENSITY); +#else + s << "\x1b[0;32m"; +#endif + return s; +} + +inline std::ostream& yellow(std::ostream &s) +{ +#if defined(WIN32) + HANDLE hStdout = GetStdHandle(STD_OUTPUT_HANDLE); + SetConsoleTextAttribute(hStdout, + FOREGROUND_GREEN|FOREGROUND_RED|FOREGROUND_INTENSITY); +#else + s << "\x1b[0;33m"; +#endif + return s; +} + +inline std::ostream& white(std::ostream &s) +{ +#if defined(WIN32) + HANDLE hStdout = GetStdHandle(STD_OUTPUT_HANDLE); + SetConsoleTextAttribute(hStdout, + FOREGROUND_RED|FOREGROUND_GREEN|FOREGROUND_BLUE); +#else + s << "\x1b[0;37m"; +#endif + return s; +} + +#endif diff --git a/Triangulation/benchmark/Triangulation/delaunay.cpp b/Triangulation/benchmark/Triangulation/delaunay.cpp index ea6c160c807..acf913d56f8 100644 --- a/Triangulation/benchmark/Triangulation/delaunay.cpp +++ b/Triangulation/benchmark/Triangulation/delaunay.cpp @@ -1,70 +1,128 @@ #include #include +#include #include #include #include +#include + #include #include #include #include #include +//#define USE_DYNAMIC_KERNEL +#define OUTPUT_STATS_IN_CSV +//#define EXPORT_POINTS_TO_A_FILE -template -void test(const int d, const std::string & type, const int N) +#ifdef OUTPUT_STATS_IN_CSV +static std::ofstream csv_file("stats.csv"); +#endif + +// Return the number of Bytes used +template +std::size_t compute_triangulation(std::size_t N) { - typedef typename DT::Vertex Vertex; - typedef typename DT::Vertex_handle Vertex_handle; - typedef typename DT::Full_cell Full_cell; - typedef typename DT::Full_cell_handle Full_cell_handle; - typedef typename DT::Facet Facet; - typedef typename DT::Point Point; - typedef typename DT::Geom_traits::RT RT; - typedef typename DT::Finite_full_cell_const_iterator Finite_full_cell_const_iterator; +#ifdef USE_DYNAMIC_KERNEL + typedef CGAL::Epick_d K; +#else + typedef CGAL::Epick_d > K; +#endif + typedef CGAL::Delaunay_triangulation DT; - typedef CGAL::Random_points_in_cube_d Random_points_iterator; - CGAL::Timer cost; // timer + typedef typename DT::Vertex Vertex; + typedef typename DT::Vertex_handle Vertex_handle; + typedef typename DT::Full_cell Full_cell; + typedef typename DT::Full_cell_handle Full_cell_handle; + typedef typename DT::Facet Facet; + typedef typename DT::Point Point; + typedef typename DT::Geom_traits::RT RT; + typedef typename DT::Finite_full_cell_const_iterator Finite_full_cell_const_iterator; - DT dt(d); - assert(dt.empty()); + typedef CGAL::Random_points_in_cube_d Random_points_iterator; + CGAL::Timer cost; // timer - std::vector points; - CGAL::Random rng; - Random_points_iterator rand_it(d, 2.0, rng); - CGAL::cpp11::copy_n(rand_it, N, std::back_inserter(points)); - cost.reset();cost.start(); - std::cout << " Delaunay triangulation of "< points; + CGAL::Random rng; + Random_points_iterator rand_it(D, 2.0, rng); + CGAL::cpp11::copy_n(rand_it, N, std::back_inserter(points)); + +#ifdef EXPORT_POINTS_TO_A_FILE + std::ofstream os("points.txt"); + for (auto const& p : points) + { + CGAL::Triangulation_IO::output_point(os, K(), p); + os << std::endl; + } +#endif + + std::size_t mem_before = CGAL::Memory_sizer().virtual_size(); + cost.reset(); + cost.start(); + + std::cout << "Delaunay triangulation of " << N << + " points in dim " << D << ":" << std::endl; + + DT dt(D); + dt.insert(points.begin(), points.end()); + + std::size_t mem = CGAL::Memory_sizer().virtual_size() - mem_before; + double timing = cost.time(); + std::cout << " Done in " << timing << " seconds." << std::endl; + std::cout << " Memory consumption: " << (mem >> 10) << " KB.\n"; + std::size_t nbfc= dt.number_of_finite_full_cells(); + std::size_t nbc= dt.number_of_full_cells(); + std::cout << " " << dt.number_of_vertices() << " vertices, " + << nbfc << " finite simplices and " + << (nbc-nbfc) << " convex hull Facets." + << std::endl; + +#ifdef OUTPUT_STATS_IN_CSV + csv_file + << D << ";" + << N << ";" + << timing << ";" + << mem << ";" + << nbfc << "\n" + << std::flush; +#endif + + + return mem; } -template< int D > -void go(const int N) +// Will compute triangulations of i*num_points_steps points, +// with i in [1, 2...], stopping after the last computation that takes +// more memory than mem_threshold_in_bytes +template +void go( + std::size_t num_points_increment, + std::size_t mem_threshold_in_MB = (3 << 10)) // 3 GB { - typedef CGAL::Epick_d > K; - //typedef CGAL::Epick_d K; - typedef CGAL::Delaunay_triangulation Triangulation; - test(D, "static", N); + std::size_t mem = 0; + for (std::size_t i = 1 ; mem < (mem_threshold_in_MB << 20) ; ++i) + { + mem = compute_triangulation(i*num_points_increment); + } } int main(int argc, char **argv) { - srand(static_cast(time(NULL))); - int N = 100; if( argc > 1 ) N = atoi(argv[1]); - go<2>(N); - go<3>(N); - go<4>(N); - go<5>(N); - go<6>(N); - go<7>(N); - go<8>(N); + srand(static_cast(time(NULL))); + //int N = 100; if( argc > 1 ) N = atoi(argv[1]); + go<2>(5000000); + //go<3>(1000000); + //go<4>(300000); + //go<5>(50000); + //go<6>(5000); + //go<7>(1000); + //go<8>(300); + //go<9>(100); + //go<10>(30); + //go<11>(20); + //go<12>(15); - - return 0; + return 0; } diff --git a/Triangulation/doc/Triangulation/CGAL/Delaunay_triangulation.h b/Triangulation/doc/Triangulation/CGAL/Delaunay_triangulation.h index 3524620b653..fe74b7336f9 100644 --- a/Triangulation/doc/Triangulation/CGAL/Delaunay_triangulation.h +++ b/Triangulation/doc/Triangulation/CGAL/Delaunay_triangulation.h @@ -7,7 +7,7 @@ namespace CGAL { This class is used to maintain the Delaunay triangulation of a set of points in \f$ \mathbb{R}^D \f$. It permits point insertion and -removal. The dimension \f$ D\f$ can be specified at compile-time or +removal. The maximal dimension \f$ D\f$ can be specified at compile-time or run-time. It should be kept reasonably small, see the performance section in the user manual for what reasonable means. @@ -20,35 +20,32 @@ A circumscribing ball of a simplex is a ball having all vertices of the simplex on its boundary. -\tparam DelaunayTriangulationTraits is the geometric traits class that provides the geometric types -and predicates needed by Delaunay triangulations. `DelaunayTriangulationTraits` must be a model of +\tparam DelaunayTriangulationTraits_ is the geometric traits class that provides the geometric types +and predicates needed by Delaunay triangulations. `DelaunayTriangulationTraits_` must be a model of the concept `DelaunayTriangulationTraits`. -\tparam TriangulationDataStructure must be a model of the concept +\tparam TriangulationDataStructure_ must be a model of the concept `TriangulationDataStructure`. This model is used to store -the faces of the triangulation. The parameter `TriangulationDataStructure` defaults to +the faces of the triangulation. The parameter `TriangulationDataStructure_` defaults to `Triangulation_data_structure` whose template parameters are instantiated as follows:
    -
  • `DelaunayTriangulationTraits::Dimension`
    • -
  • `Triangulation_vertex`
    • -
  • `Triangulation_full_cell`.
    • +
  • `DelaunayTriangulationTraits_::Dimension`
    • +
  • `Triangulation_vertex`
    • +
  • `Triangulation_full_cell`.
-The class template `Delaunay_triangulation` can +\tparam Delaunay_triangulation can be defined by specifying only the first parameter, or by using the tag `CGAL::Default` as the second parameter. -The class `Delaunay_triangulation` inherits all the types -defined in the base class `Triangulation`. Additionally, it -defines or overloads the following methods: - -\sa `Triangulation_data_structure` +\sa `Regular_triangulation` +\sa `Triangulation_data_structure` */ -template< typename DelaunayTriangulationTraits, typename TriangulationDataStructure > +template< typename DelaunayTriangulationTraits_, typename TriangulationDataStructure_ > class Delaunay_triangulation - : public Triangulation + : public Triangulation { public: @@ -62,7 +59,7 @@ at infinity). See the description of the inherited nested type the use of the parameter `dim`. The complex stores a copy of the geometric traits `gt`. */ -Delaunay_triangulation(const int dim, const Geom_traits gt = Geom_traits()); +Delaunay_triangulation(int dim, const Geom_traits > = Geom_traits()); /// @} @@ -118,7 +115,6 @@ Same as above but uses a vertex as starting place for the search. Vertex_handle insert(const Point & p, Vertex_handle hint); /*! -\cgalAdvancedBegin Inserts the point `p` in the Delaunay triangulation and ensures that the empty-ball property is preserved. Returns a handle to the @@ -128,44 +124,19 @@ value of `lt`:
`OUTSIDE_AFFINE_HULL`
Point `p` is inserted so as to increase the current dimension of the Delaunay -triangulation. The method `dt`.`insert_outside_affine_hull()` is called. +triangulation.
`ON_VERTEX`
The position of the vertex `v` described by `f` is set to `p`. `v` is returned.
Anything else
The point `p` is inserted. the full cell `c` is assumed to be in conflict with `p`. -(Roughly speaking, the method `dt`.`insert_in_conflicting_cell()` -is called.)
The parameters `lt`, `f`, `ft` and `c` must be consistent with the localization of point `p` in the Delaunay triangulation e.g. by a call to -`c = locate(p, lt, f, ft)`. -\cgalAdvancedEnd +`Triangulation::locate(const Point &, Locate_type &, Face &, Vertex_handle) const`. */ -Vertex_handle insert(const Point & p, const Locate_type lt, -const Face & f, const Facet & ft, const Full_cell_handle c); - -/*! -\cgalAdvancedBegin -Inserts the point `p` in the Delaunay triangulation. Returns a handle to the -(possibly newly created) vertex at that position. -\pre The point `p` -must lie outside the affine hull of the Delaunay triangulation. This implies that -`dt`.`current_dimension()` must be less that -`dt`.`maximal_dimension()`. -\cgalAdvancedEnd -*/ -Vertex_handle insert_outside_affine_hull(const Point & p); - -/*! -\cgalAdvancedBegin -Inserts the point `p` in the Delaunay triangulation. Returns a handle to the -(possibly newly created) vertex at that position. -\pre The point `p` must be in conflict with the full cell `c`. -\cgalAdvancedEnd -*/ -Vertex_handle insert_in_conflicting_cell(const Point & p, const -Full_cell_handle c); +Vertex_handle insert(const Point & p, Locate_type lt, +const Face & f, const Facet & ft, Full_cell_handle c); /// @} @@ -177,25 +148,20 @@ Returns `true` if and only if the point `p` is in (Delaunay) conflict with full cell `c` (i.e., the circumscribing ball of \f$ c\f$ contains \f$ p\f$ in its interior). */ -bool is_in_conflict(const Point & p, Full_cell_const_handle c) -const; +bool is_in_conflict(const Point & p, Full_cell_const_handle c) const; /*! -\cgalAdvancedBegin -Outputs handles to the full cells in confict with +Outputs handles to the full cells in conflict with point `p` into the `OutputIterator out`. The full cell `c` is used as a starting point for gathering the full cells in conflict with `p`. A facet `(cc,i)` on the boundary of the conflict zone with `cc` in conflict is returned. -\pre `c` is in conflict -with `p`. -`dt`.`current_dimension()`\f$ \geq2\f$. -\cgalAdvancedEnd +\pre `c` is in conflict with `p` and `dt`.`current_dimension()`\f$ \geq2\f$. */ template< typename OutputIterator > -Facet compute_conflict_zone(const Point & p, const Full_cell_handle c, -OutputIterator out) const; +Facet compute_conflict_zone(const Point & p, Full_cell_handle c, + OutputIterator out) const; /// @} diff --git a/Triangulation/doc/Triangulation/CGAL/Regular_triangulation.h b/Triangulation/doc/Triangulation/CGAL/Regular_triangulation.h new file mode 100644 index 00000000000..c22ed0d7e17 --- /dev/null +++ b/Triangulation/doc/Triangulation/CGAL/Regular_triangulation.h @@ -0,0 +1,171 @@ + +namespace CGAL { + +/*! +\ingroup PkgTriangulationsTriangulationClasses + +This class is used to maintain the +regular triangulation -- also known as weighted Delaunay triangulation -- +of a set of weighted points in \f$ \mathbb{R}^D \f$. +The maximal dimension \f$ D\f$ can be specified at compile-time or +run-time. It should be kept reasonably small -- see the performance +section in the user manual for what reasonable means. + +\warning The removal of points is not supported yet. + +\tparam RegularTriangulationTraits_ is the geometric traits class that provides the +geometric types and predicates needed by regular triangulations. +`RegularTriangulationTraits_` must be a model of the concept +`RegularTriangulationTraits`. + +\tparam TriangulationDataStructure_ must be a model of the concept +`TriangulationDataStructure`. This model is used to store +the faces of the triangulation. The parameter `TriangulationDataStructure_` +defaults to `Triangulation_data_structure` whose template parameters are +instantiated as follows: +
    +
  • `RegularTriangulationTraits_::Dimension`
    • +
  • `Triangulation_vertex >`
    • +
  • `Triangulation_full_cell >`.
    • +
+ +`Regular_triangulation` can +be defined by specifying only the first parameter, or by using the +tag `CGAL::Default` as the second parameter. + +\sa `Delaunay_triangulation` +\sa `Triangulation_data_structure` +\sa `Regular_triangulation_traits_adapter` + +*/ +template< typename RegularTriangulationTraits_, typename TriangulationDataStructure_ > +class Regular_triangulation + : public Triangulation, TriangulationDataStructure_> +{ +public: + +/// \name Types +/// @{ + +/*! +A point in Euclidean space with an associated weight. +*/ +typedef RegularTriangulationTraits_::Weighted_point_d Weighted_point; + +/// @} + +/// \name Creation +/// @{ + +/*! +Instantiates a regular triangulation with one vertex (the vertex +at infinity). See the description of the inherited nested type +`Triangulation::Maximal_dimension` for an explanation of +the use of the parameter `dim`. The triangulation stores a copy of the +geometric traits `gt`. +*/ +Regular_triangulation(int dim, const Geom_traits > = Geom_traits()); + +/// @} + +/// \name Point Insertion +/// @{ + +/*! +Inserts weighted point `p` in the triangulation and returns the corresponding +vertex. + +If this insertion creates a vertex, this vertex is returned. + +If `p` coincides with an existing vertex and has a greater weight, +then the existing weighted point becomes hidden and `p` replaces it as vertex +of the triangulation. + +If `p` coincides with an already existing vertex (both point and +weights being equal), then this vertex is returned and the triangulation +remains unchanged. + +Otherwise if `p` does not appear as a vertex of the triangulation, +then it is stored as a hidden point and this method returns the default +constructed handle. + +Prior to the actual insertion, `p` is located in the triangulation; +`hint` is used as a starting place for locating `p`. +*/ +Vertex_handle insert(const Weighted_point & p, Full_cell_handle hint + = Full_cell_handle()); + +/*! +Same as above but uses a vertex as starting place for the search. +*/ +Vertex_handle insert(const Weighted_point & p, Vertex_handle hint); + +/*! +Inserts weighted point `p` in the triangulation. +Similar to the above `insert()` function, but takes as additional +parameter the return values of a previous location query. See description of +`Triangulation::locate()`. +*/ +Vertex_handle insert(const Weighted_point & p, Locate_type lt, + const Face & f, const Facet & ft, Full_cell_handle c); + +/*! +Inserts the weighted points found in range `[s,e)` in the regular triangulation. +Returns the difference between the number of vertices after and before +the insertions (it may be negative due to hidden points). +Note that this function is not guaranteed to insert the points +following the order of `ForwardIterator` because `spatial_sort()` +is used to improve efficiency. + +\tparam ForwardIterator must be an input iterator with the value +type `Weighted_point`. +*/ +template< typename ForwardIterator > +std::ptrdiff_t insert(ForwardIterator s, ForwardIterator e); + +/// @} + +/// \name Queries +/// @{ + +/*! +Returns `true` if and only if the point `p` is in +conflict with full cell `c` (A weighted point `p` is said to be in conflict +with a cell `c` if it has a negative power distance to the power sphere of `c`.) +*/ +bool is_in_conflict(const Weighted_point & p, Full_cell_const_handle c) +const; + +/*! +Outputs handles to the full cells in conflict with +point `p` into the `OutputIterator out`. The full cell `c` is used +as a starting point for gathering the full cells in conflict with +`p`. +A facet `(cc,i)` on the boundary of the conflict zone with +`cc` in conflict is returned. +\pre `c` is in conflict with `p` and `rt`.`current_dimension()`\f$ \geq 1\f$. +*/ +template< typename OutputIterator > +Facet compute_conflict_zone(const Weighted_point & p, Full_cell_handle c, +OutputIterator out) const; + +/// @} + +/// \name Access Functions +/// @{ + +/*! +Returns the number of finite vertices that are not hidden. +*/ +size_type number_of_vertices() const; + +/*! +Returns the number of hidden vertices. +*/ +size_type number_of_hidden_vertices() const; + +/// @} + + +}; /* end regular_triangulation */ +} /* end namespace CGAL */ diff --git a/Triangulation/doc/Triangulation/CGAL/Regular_triangulation_traits_adapter.h b/Triangulation/doc/Triangulation/CGAL/Regular_triangulation_traits_adapter.h new file mode 100644 index 00000000000..c1654d8efb0 --- /dev/null +++ b/Triangulation/doc/Triangulation/CGAL/Regular_triangulation_traits_adapter.h @@ -0,0 +1,62 @@ + +namespace CGAL { + +/*! +\ingroup PkgTriangulationsTraitsClasses + +The class `Regular_triangulation_traits_adapter` is used internally by the +class `Regular_triangulation` to wrap its first template parameter +(`RegularTriangulationTraits_`) +so that the base class `Triangulation` manipulates weighted points instead +of bare points. + +\tparam K must be a model of the `RegularTriangulationTraits` concept. + +In addition to the types described below, the following predicates and functors +are adapted so that they can be called +with weighted points instead of bare points as parameters. +In practice, the functors from the base class `K` are called, +ignoring the weights. +- `Orientation_d` +- `Construct_flat_orientation_d` +- `In_flat_orientation_d` +- `Contained_in_affine_hull_d` +- `Compare_lexicographically_d` +- `Compute_coordinate_d` +- `Point_dimension_d` +- `Less_coordinate_d` + +*/ + +template +class Regular_triangulation_traits_adapter : public K { +public: + + /// \name Types + /// @{ + + /*! + The weighted point type. + */ + typedef typename K::Weighted_point_d Point_d; + + /*! + The (un-weighted) point type. + */ + typedef typename K::Point_d Bare_point_d; + + /// @} + + /// \name Creation + /// @{ + + /*! + The default constructor. + */ + Regular_triangulation_traits_adapter(); + + /// @} + +}; + +} /* end namespace CGAL */ diff --git a/Triangulation/doc/Triangulation/CGAL/Triangulation.h b/Triangulation/doc/Triangulation/CGAL/Triangulation.h index 6a78db853e3..bf01cfdce74 100644 --- a/Triangulation/doc/Triangulation/CGAL/Triangulation.h +++ b/Triangulation/doc/Triangulation/CGAL/Triangulation.h @@ -4,7 +4,7 @@ namespace CGAL { /*! \ingroup PkgTriangulationsTriangulationClasses -This class implements triangulations of point sets in dimensions \f$ d \f$. +This class implements triangulations of point sets in dimension \f$ d \f$. The triangulation covers the convex hull of the input points (the embedded vertices of the triangulation). @@ -17,25 +17,25 @@ incident to the infinite vertex and to an \f$ (i-1)\f$-simplex of the convex hull boundary. -\tparam TriangulationTraits is the geometric traits class that provides the geometric types -and predicates needed by triangulations. `TriangulationTraits` must be a model of the +\tparam TriangulationTraits_ is the geometric traits class that provides the geometric types +and predicates needed by triangulations. `TriangulationTraits_` must be a model of the concept `TriangulationTraits`. -\tparam TriangulationDataStructure must be a model of the concept +\tparam TriangulationDataStructure_ must be a model of the concept `TriangulationDataStructure`. This model is used to store -the faces of the triangulation. The parameter `TriangulationDataStructure` defaults to +the faces of the triangulation. The parameter `TriangulationDataStructure_` defaults to `Triangulation_data_structure` whose template parameters are instantiated as follows:
    -
  • `DelaunayTriangulationTraits::Dimension`
    • -
  • `Triangulation_vertex`
    • -
  • `Triangulation_full_cell`.
    • +
  • `TriangulationTraits_::Dimension`
    • +
  • `Triangulation_vertex`
    • +
  • `Triangulation_full_cell`.
The triangulation deduces its maximal dimension from the type -`TriangulationTraits::Dimension`. This dimension has to match +`TriangulationTraits_::Dimension`. This dimension has to match the dimension returned by -`TriangulationDataStructure::maximal_dimension()`. +`TriangulationDataStructure_::maximal_dimension()`. Input/Output -------------- @@ -47,62 +47,64 @@ full cell, plus the non-combinatorial information about each full cell, then the indices of the neighbors of each full cell, where the index corresponds to the preceding list of full cells. -\sa `Triangulation_data_structure` -\sa `Delaunay_triangulation` +\sa `Triangulation_data_structure` +\sa `Delaunay_triangulation` */ -template< typename TriangulationTraits, typename TriangulationDataStructure > +template< typename TriangulationTraits_, typename TriangulationDataStructure_> class Triangulation { public: /// \name Types /// @{ /*! -Type for the model of the `TriangulationTraits` concept. +Type for the model of the `TriangulationTraits_` concept. */ -typedef TriangulationTraits Geom_traits; +typedef TriangulationTraits_ Geom_traits; /*! -A point in Euclidean space. +A point in Euclidean space. Note that in the context of a +`Regular_triangulation` class (which derives from this class), +`TriangulationTraits_::Point_d` is a weighted point. */ -typedef TriangulationTraits::Point_d Point; +typedef TriangulationTraits_::Point_d Point; /*! This indicates whether the maximal dimension is static (i.e.\ if the type of `Maximal_dimension` is `CGAL::Dimension_tag`) or dynamic (i.e.\ if the type of `Maximal_dimension` is `CGAL::Dynamic_dimension_tag`). -In the latter case, the `dim` parameter passed to the class's constructor -is used. +In the latter case, the `dim` parameter passed to the constructor +of the class is used. */ -typedef TriangulationTraits::Dimension Maximal_dimension; +typedef TriangulationTraits_::Dimension Maximal_dimension; /*! The second template parameter: the triangulation data structure. */ -typedef TriangulationDataStructure Triangulation_ds; +typedef TriangulationDataStructure_ Triangulation_ds; /*! A model of the concept `TriangulationVertex`. */ -typedef TriangulationDataStructure::Vertex Vertex; +typedef TriangulationDataStructure_::Vertex Vertex; /*! A model of the concept `TriangulationFullCell`. */ -typedef TriangulationDataStructure::Full_cell Full_cell; +typedef TriangulationDataStructure_::Full_cell Full_cell; /*! The facet class */ -typedef TriangulationDataStructure::Facet Facet; +typedef TriangulationDataStructure_::Facet Facet; /*! A model of the concept `TriangulationDSFace`. */ -typedef TriangulationDataStructure::Face Face; +typedef TriangulationDataStructure_::Face Face; /// @} @@ -122,25 +124,25 @@ typedef TriangulationDataStructure::Face Face; /*! handle to a a vertex */ -typedef TriangulationDataStructure::Vertex_handle +typedef TriangulationDataStructure_::Vertex_handle Vertex_handle; /*! const handle to a a vertex */ -typedef TriangulationDataStructure::Vertex_const_handle +typedef TriangulationDataStructure_::Vertex_const_handle Vertex_const_handle; /*! iterator over all vertices (including the infinite one) */ -typedef TriangulationDataStructure::Vertex_iterator +typedef TriangulationDataStructure_::Vertex_iterator Vertex_iterator; /*! const iterator over all vertices (including the infinite one) */ -typedef TriangulationDataStructure::Vertex_const_iterator +typedef TriangulationDataStructure_::Vertex_const_iterator Vertex_const_iterator; /*! @@ -156,27 +158,27 @@ typedef unspecified_type Finite_vertex_const_iterator; /*! handle to a full cell */ -typedef TriangulationDataStructure::Full_cell_handle +typedef TriangulationDataStructure_::Full_cell_handle Full_cell_handle; /*! const handle to a full cell */ -typedef TriangulationDataStructure::Full_cell_const_handle +typedef TriangulationDataStructure_::Full_cell_const_handle Full_cell_const_handle; /*! iterator over all full cells (including the infinite ones) */ typedef -TriangulationDataStructure::Full_cell_iterator +TriangulationDataStructure_::Full_cell_iterator Full_cell_iterator; /*! const iterator over all full cells (including the infinite ones) */ typedef -TriangulationDataStructure::Full_cell_const_iterator +TriangulationDataStructure_::Full_cell_const_iterator Full_cell_const_iterator; /*! @@ -192,7 +194,7 @@ typedef unspecified_type Finite_full_cell_const_iterator; /*! iterator over all facets (including the infinite ones) */ -typedef TriangulationDataStructure::Facet_iterator +typedef TriangulationDataStructure_::Facet_iterator Facet_iterator; /*! @@ -201,19 +203,19 @@ iterator over finite facets typedef unspecified_type Finite_facet_iterator; /*! -Size type (an unsigned integral -type). +size type (an unsigned integral type) */ -typedef TriangulationDataStructure::size_type size_type; +typedef TriangulationDataStructure_::size_type size_type; /*! -Difference -type (a signed integral type). +difference +type (a signed integral type) */ -typedef TriangulationDataStructure::difference_type difference_type; +typedef TriangulationDataStructure_::difference_type difference_type; /*! -specifies which case occurs when locating a point in the triangulation. +\enum Locate_type +\brief Used by `Triangulation` to specify which case occurs when locating a point in the triangulation. */ enum Locate_type { ON_VERTEX=0, /*!< when the located point coincides with a vertex of the triangulation */ IN_FACE, /*!< when the point is in the interior of a face of dimension equal or less than `maximal_dimension()` - 2 */ @@ -234,8 +236,7 @@ description of the nested type `Maximal_dimension` above for an explanation of the use of the parameter `dim`. The triangulation stores a copy of the geometric traits `gt`. */ -Triangulation(const int dim, const Geom_traits & gt = -Geom_traits()); +Triangulation(int dim, const Geom_traits & gt = Geom_traits()); /*! The copy constructor. @@ -248,15 +249,15 @@ Triangulation(const Triangulation & t2); /// @{ /*! -\cgalAdvancedBegin Returns a const reference to the underlying triangulation data structure. -\cgalAdvancedEnd */ const Triangulation_ds & tds() const; /*! +\cgalAdvancedBegin Returns a non-const reference to the underlying triangulation data structure. +\cgalAdvancedEnd */ Triangulation_ds & tds(); @@ -322,13 +323,13 @@ size_type number_of_finite_full_cells() const; /*! Returns `true` if and only if the vertex `v` is the infinite vertex. */ -bool is_infinite(const Vertex_handle v) const; +bool is_infinite(Vertex_handle v) const; /*! Returns `true` if and only if `c` is incident to the infinite vertex. */ -bool is_infinite(const Full_cell_handle c) const; +bool is_infinite(Full_cell_handle c) const; /*! Returns `true` if and only if facet `ft` is incident to the infinite @@ -454,7 +455,7 @@ p_1, p_2, \ldots, p_d, \infty\}\f$ is returned such that the full cell \f$ (p_1, on the other side of facet \f$ (p_1, p_2, \ldots, p_d)\f$). */ Full_cell_handle locate(const Point & query, -Full_cell_handle hint = Full_cell_handle()) const; + Full_cell_const_handle hint = Full_cell_handle()) const; /*! Same as above but `hint` is a vertex and not a full cell. @@ -505,7 +506,7 @@ The parameter `hint` is ignored if it is a default constructed */ Full_cell_handle locate(const Point & query, Locate_type & loc_type, -Face & f, Vertex_handle hint) const; + Face & f, Vertex_handle hint) const; /// @} @@ -513,13 +514,11 @@ Face & f, Vertex_handle hint) const; /// @{ /*! -\cgalAdvancedBegin Contracts the `Face f` to a single vertex at position `p`. Returns a handle to that vertex. \pre The boundary of the union of full cells incident to `f` must be a triangulation of a sphere of dimension `tr`.`current_dimension()`). -\cgalAdvancedEnd */ Vertex_handle collapse_face(const Point & p, const Face & f); @@ -544,16 +543,15 @@ Inserts point `p` in the triangulation. Returns a Prior to the actual insertion, `p` is located in the triangulation; `hint` is used as a starting place for locating `p`. */ -Vertex_handle insert(const Point p, Full_cell_handle hint = +Vertex_handle insert(const Point &p, Full_cell_handle hint = Full_cell_handle()); /*! Same as above but uses a vertex `hint` as the starting place for the search. */ -Vertex_handle insert(const Point p, Vertex_handle hint); +Vertex_handle insert(const Point &p, Vertex_handle hint); /*! -\cgalAdvancedBegin Inserts point `p` into the triangulation and returns a handle to the `Vertex` at that position. The action taken depends on the value of `loc_type`: @@ -567,12 +565,10 @@ of `loc_type`, using the full cell `c`. This method is used internally by the other `insert()` methods. -\cgalAdvancedEnd */ -Vertex_handle insert(const Point p, Locate_type loc_type, Face & f, Facet & ft, Full_cell_handle c); +Vertex_handle insert(const Point &p, Locate_type loc_type, Face & f, Facet & ft, Full_cell_handle c); /*! -\cgalAdvancedBegin Removes the full cells in the range \f$ C=\f$`[s, e)`, inserts a vertex at position `p` and fills the hole by connecting each face of the boundary to `p`. @@ -582,62 +578,49 @@ defining full cell, `tr`.`full_cell(ft)` must lie inside \f$ C\f$. Handles to the newly created full cells are output in the `out` output iterator. \pre \f$C\f$ must be a topological ball, must contain `p` in its interior and must not contain any vertex of the triangulation. -\cgalAdvancedEnd */ template < typename ForwardIterator, typename OutputIterator > Vertex_handle insert_in_hole(const Point & p, ForwardIterator s, ForwardIterator e, const Facet & ft, OutputIterator out); /*! -\cgalAdvancedBegin Same as above, but the newly created full cells are not retrieved. -\cgalAdvancedEnd */ template < typename ForwardIterator > Vertex_handle insert_in_hole(const Point & p, ForwardIterator s, ForwardIterator e, const Facet & ft); /*! -\cgalAdvancedBegin Inserts point `p` in the triangulation. \pre `p` must lie in the relative interior of `f`. -\cgalAdvancedEnd */ Vertex_handle insert_in_face(const Point & p, const Face & f); /*! -\cgalAdvancedBegin Inserts point `p` in the triangulation. \pre `p` must lie in the relative interior of `ft`. -\cgalAdvancedEnd */ Vertex_handle insert_in_facet(const Point & p, const Facet & ft); /*! -\cgalAdvancedBegin Inserts point `p` in the triangulation. \pre `p` must lie in the interior of `c`. -\cgalAdvancedEnd */ Vertex_handle insert_in_full_cell(const Point & p, Full_cell_handle c); /*! -\cgalAdvancedBegin Inserts point `p` in the triangulation. \pre `p` must lie outside the convex hull of `tr`. The half-space defined by the infinite full cell `c` must contain `p`. -\cgalAdvancedEnd */ Vertex_handle insert_outside_convex_hull(const Point &, Full_cell_handle c); /*! -\cgalAdvancedBegin Inserts point `p` in the triangulation. \pre `p` must lie outside the affine hull of `tr`. -\cgalAdvancedEnd */ Vertex_handle insert_outside_affine_hull(const Point &); diff --git a/Triangulation/doc/Triangulation/CGAL/Triangulation_data_structure.h b/Triangulation/doc/Triangulation/CGAL/Triangulation_data_structure.h index 32807379952..bc2459e2dbf 100644 --- a/Triangulation/doc/Triangulation/CGAL/Triangulation_data_structure.h +++ b/Triangulation/doc/Triangulation/CGAL/Triangulation_data_structure.h @@ -9,29 +9,29 @@ of dimension \f$ d\leq D\f$ (`D` is the maximal dimension). \tparam Dimensionality can be either
    -
  • CGAL::`Dimension_tag` for some integer `D`. This +
  • `CGAL::Dimension_tag` for some integer `D`. This indicates that the triangulation data structure can store simplices (full cells) of dimension at most `D`. The maximal dimension `D` is known by the compiler, which triggers some optimizations. Or -
  • CGAL::`Dynamic_dimension_tag`. In this case, the maximum +
  • `CGAL::Dynamic_dimension_tag`. In this case, the maximum dimension of the simplices (full cells) is passed as an integer argument to an instance constructor (see `TriangulationDataStructure`).
-\tparam TriangulationDSVertex stands for a class to +\tparam `TriangulationDSVertex_` stands for a class to be used as the base `Vertex` type in the triangulation data structure. It must be a model of the concept `TriangulationDSVertex`. The class template `Triangulation_data_structure` can be defined by specifying only the first parameter. It also accepts the tag `CGAL::Default` as -second parameter. In both cases, `TriangulationDSVertex` defaults to +second parameter. In both cases, `TriangulationDSVertex_` defaults to `CGAL::Triangulation_ds_vertex<>`. -\tparam TriangulationDSFullCell stands for a class to +\tparam `TriangulationDSFullCell_` stands for a class to be used as the base `Full_cell` type in the triangulation data structure. It must be a model of the concept `TriangulationDSFullCell`. The class template `Triangulation_data_structure` accepts that no third parameter be specified. It also accepts the tag `CGAL::Default` as -third parameter. In both cases, `TriangulationDSFullCell` defaults to +third parameter. In both cases, `TriangulationDSFullCell_` defaults to `CGAL::Triangulation_ds_full_cell<>`. \cgalModels `TriangulationDataStructure`. In addition, the class @@ -41,7 +41,7 @@ methods. \sa `Triangulation_ds_vertex` \sa `Triangulation_ds_full_cell` */ -template< typename Dimensionality, typename TriangulationDSVertex, typename TriangulationDSFullCell > +template< typename Dimensionality, typename TriangulationDSVertex_, typename TriangulationDSFullCell_ > class Triangulation_data_structure { public: diff --git a/Triangulation/doc/Triangulation/CGAL/Triangulation_ds_full_cell.h b/Triangulation/doc/Triangulation/CGAL/Triangulation_ds_full_cell.h index c79e3e9996b..fccb180dfe4 100644 --- a/Triangulation/doc/Triangulation/CGAL/Triangulation_ds_full_cell.h +++ b/Triangulation/doc/Triangulation/CGAL/Triangulation_ds_full_cell.h @@ -17,7 +17,7 @@ This class can be used directly or can serve as a base to derive other classes with some additional attributes tuned for a specific application. -\tparam TriangulationDataStructure must be a model of the +\tparam `TriangulationDataStructure_` must be a model of the `TriangulationDataStructure` concept. \tparam TriangulationDSFullCellStoragePolicy indicates whether or not @@ -49,11 +49,11 @@ Rebind mechanism In case of derivation from that class, the nested class `Rebind_TDS` need to be provided in the derived class. -\sa `Triangulation_ds_vertex` -\sa `Triangulation_data_structure>` +\sa `Triangulation_ds_vertex` +\sa `Triangulation_data_structure` */ -template< typename TriangulationDataStructure, typename TriangulationDSFullCellStoragePolicy > +template< typename TriangulationDataStructure_, typename TriangulationDSFullCellStoragePolicy > class Triangulation_ds_full_cell { public: diff --git a/Triangulation/doc/Triangulation/CGAL/Triangulation_ds_vertex.h b/Triangulation/doc/Triangulation/CGAL/Triangulation_ds_vertex.h index fd9feffa3f3..6053146b13b 100644 --- a/Triangulation/doc/Triangulation/CGAL/Triangulation_ds_vertex.h +++ b/Triangulation/doc/Triangulation/CGAL/Triangulation_ds_vertex.h @@ -5,7 +5,7 @@ namespace CGAL { \ingroup PkgTriangulationsVertexCellClasses The class `Triangulation_ds_vertex` serves as the default vertex template parameter in the -class `Triangulation_data_structure`. +class `Triangulation_data_structure`. This class does not contain any geometric information but only combinatorial (adjacency) information. Thus, if the `Triangulation_data_structure` is @@ -18,7 +18,7 @@ with some additional attributes tuned for a specific application (a color for example). -\tparam TriangulationDataStructure must be a model of the +\tparam `TriangulationDataStructure_` must be a model of the `TriangulationDataStructure` concept. \cgalModels `TriangulationDSVertex` @@ -29,11 +29,11 @@ Rebind Mechanism In case of derivation from that class, the nested class `Rebind_TDS` need to be provided in the derived class. -\sa `Triangulation_ds_full_cell` -\sa `Triangulation_data_structure>` +\sa `Triangulation_ds_full_cell` +\sa `Triangulation_data_structure` */ -template< typename TriangulationDataStructure > +template< typename TriangulationDataStructure_ > class Triangulation_ds_vertex { public: diff --git a/Triangulation/doc/Triangulation/CGAL/Triangulation_face.h b/Triangulation/doc/Triangulation/CGAL/Triangulation_face.h index 46f6d438a29..798bec78323 100644 --- a/Triangulation/doc/Triangulation/CGAL/Triangulation_face.h +++ b/Triangulation/doc/Triangulation/CGAL/Triangulation_face.h @@ -5,7 +5,7 @@ namespace CGAL { A `Triangulation_face` is a model of the concept `TriangulationDSFace`. -\tparam TriangulationDataStructure must be a model of the concept +\tparam TriangulationDataStructure_ must be a model of the concept `TriangulationDataStructure`. Actually, `Triangulation_face` needs only that this concept defines the types `Full_cell_handle`, @@ -18,7 +18,7 @@ Actually, `Triangulation_face` needs only that this concept defines the types \sa `TriangulationDataStructure` */ -template< typename TriangulationDataStructure > +template< typename TriangulationDataStructure_ > class Triangulation_face { }; /* end Triangulation_face */ diff --git a/Triangulation/doc/Triangulation/CGAL/Triangulation_full_cell.h b/Triangulation/doc/Triangulation/CGAL/Triangulation_full_cell.h index e0abe85691c..a6d8a002a6c 100644 --- a/Triangulation/doc/Triangulation/CGAL/Triangulation_full_cell.h +++ b/Triangulation/doc/Triangulation/CGAL/Triangulation_full_cell.h @@ -6,15 +6,15 @@ namespace CGAL { The class `Triangulation_full_cell` is a model of the concept `TriangulationFullCell`. It is used by default for representing full cells in the class -`Triangulation`. +`Triangulation`. A `Triangulation_full_cell` stores handles to the vertices of the cell as well as handles to its adjacent cells. -\tparam TriangulationTraits must be a model of the concept `TriangulationTraits`. It +\tparam `TriangulationTraits_` must be a model of the concept `TriangulationTraits`. It provides geometric types and predicates for use in the -`Triangulation` class. +`Triangulation` class. \tparam Data is an optional type of data to be stored in the full cell class. The class template `Triangulation_full_cell` accepts that no second parameter be specified. In @@ -31,13 +31,13 @@ cases, `TriangulationDSFullCell_` defaults to `CGAL::Triangulation_ds_full_cell< `Triangulation_full_cell` provides the following types, constructors and methods: -\sa `Triangulation_vertex` -\sa `Triangulation_data_structure` -\sa `Triangulation` -\sa `Delaunay_triangulation` +\sa `Triangulation_vertex` +\sa `Triangulation_data_structure` +\sa `Triangulation` +\sa `Delaunay_triangulation` */ -template< typename TriangulationTraits, typename Data, typename TriangulationDSFullCell_ > +template< typename TriangulationTraits_, typename Data, typename TriangulationDSFullCell_ > class Triangulation_full_cell : public TriangulationDSFullCell_ { public: diff --git a/Triangulation/doc/Triangulation/CGAL/Triangulation_vertex.h b/Triangulation/doc/Triangulation/CGAL/Triangulation_vertex.h index d9a47a43428..ddb31b55090 100644 --- a/Triangulation/doc/Triangulation/CGAL/Triangulation_vertex.h +++ b/Triangulation/doc/Triangulation/CGAL/Triangulation_vertex.h @@ -6,14 +6,14 @@ namespace CGAL { The class `Triangulation_vertex` is a model of the concept `TriangulationVertex`. It is used by default for representing vertices in the class -`Triangulation`. +`Triangulation`. A `Triangulation_vertex` stores a point and an incident full cell. -\tparam TriangulationTraits must be a model of the concept `TriangulationTraits`. It +\tparam `TriangulationTraits_` must be a model of the concept `TriangulationTraits`. It provides geometric types and predicates for use in the -`Triangulation` class. It is of interest here for its +`Triangulation` class. It is of interest here for its declaration of the `Point` type. \tparam Data is an optional type of data to be stored in the vertex class. The @@ -22,27 +22,32 @@ this case, `Data` defaults to `CGAL::No_vertex_data`. `CGAL::No_vertex_data` can be explicitely specified to allow to access the third parameter. -\tparam TriangulationDSVertex must be a model of the concept `TriangulationDSVertex`. The +\tparam `TriangulationDSVertex_` must be a model of the concept `TriangulationDSVertex`. The class template `Triangulation_vertex` accepts that no third parameter be specified. It also accepts the tag `CGAL::Default` as third parameter. In both cases, -`TriangulationDSVertex` defaults to `CGAL::Triangulation_ds_vertex<>`. +`TriangulationDSVertex_` defaults to `CGAL::Triangulation_ds_vertex<>`. \cgalModels `TriangulationVertex` Additionally, the class `Triangulation_vertex` provides the following types, constructors and methods: -\sa `Triangulation_full_cell` -\sa `Triangulation_data_structure` -\sa `Triangulation` -\sa `Delaunay_triangulation` +\sa `Triangulation_full_cell` +\sa `Triangulation_data_structure` +\sa `Triangulation` +\sa `Delaunay_triangulation` */ -template< typename TriangulationTraits, typename Data, typename TriangulationDSVertex > -class Triangulation_vertex { +template< typename TriangulationTraits_, typename Data, typename TriangulationDSVertex_ > +class Triangulation_vertex { public: /// \name Types /// @{ +/*! +The point type. +*/ +typedef TriangulationTraits_::Point_d Point; + /*! The type of the additional data stored in the vertex. If you read a `Triangulation_vertex` from a stream (a file) or write a `Triangulation_vertex` to a stream, then streaming operators `<<` and `>>` must be provided for this diff --git a/Triangulation/doc/Triangulation/Concepts/DelaunayTriangulationTraits.h b/Triangulation/doc/Triangulation/Concepts/DelaunayTriangulationTraits.h index a65cf615751..f5d6345384c 100644 --- a/Triangulation/doc/Triangulation/Concepts/DelaunayTriangulationTraits.h +++ b/Triangulation/doc/Triangulation/Concepts/DelaunayTriangulationTraits.h @@ -5,7 +5,7 @@ This concept describes the geometric types and predicates required to build a Delaunay triangulation. It corresponds to the first template parameter of the class -`CGAL::Delaunay_triangulation`. +`CGAL::Delaunay_triangulation`. \cgalRefines `TriangulationTraits` @@ -32,7 +32,7 @@ defined by the points in range `[start,end)`. If the simplex is positively oriented, then the positive side of sphere corresponds geometrically to its bounded side. -\pre If `Dimension`=`CGAL::``Dimension_tag`, +\pre If `Dimension`=`CGAL::Dimension_tag`, then `std::distance(start,end)=D+1`. The points in range `[start,end)` must be affinely independent, i.e., the simplex must @@ -70,14 +70,16 @@ typedef unspecified_type In_flat_side_of_oriented_sphere_d; /// @{ /*! -The default constructor. +The default constructor (optional). +This is not required when an instance of the traits is provided +to the constructor of `CGAL::Delaunay_triangulation`. */ DelaunayTriangulationTraits(); /// @} /// \name Operations -/// The following methods permit access to the traits class's predicates: +/// The following methods permit access to the traits class's predicates and functors: /// @{ /*! diff --git a/Triangulation/doc/Triangulation/Concepts/RegularTriangulationTraits.h b/Triangulation/doc/Triangulation/Concepts/RegularTriangulationTraits.h new file mode 100644 index 00000000000..ca0e3a0d9cd --- /dev/null +++ b/Triangulation/doc/Triangulation/Concepts/RegularTriangulationTraits.h @@ -0,0 +1,132 @@ + +/*! +\ingroup PkgTriangulationsConcepts +\cgalConcept + +This concept describes the geometric types and predicates required to build +a regular triangulation. It corresponds to the first template parameter of the class +`CGAL::Regular_triangulation`. + +\cgalRefines `TriangulationTraits` + +\cgalHasModel `CGAL::Epick_d` + +\sa `TriangulationTraits` +*/ +class RegularTriangulationTraits { +public: + +/// \name Types +/// @{ + +/*! +A number type that is a model for `FieldNumberType`. +*/ +typedef unspecified_type FT; + +/*! +The weighted point type. +*/ +typedef unspecified_type Weighted_point_d; + +/*! +A function object that must provide the operator +`Point_d operator()(const Weighted_point_d & wp)`, returning +`wp` without its weight. +*/ +typedef unspecified_type Construct_point_d; + +/*! +A function object that must provide the operator +`FT operator()(const Weighted_point_d & wp)`, returning +the weight of `wp`. +*/ +typedef unspecified_type Compute_weight_d; + +/*! +A predicate object that must provide the templated operator +`template Oriented_side operator()(ForwardIterator start, ForwardIterator end, const Weighted_point_d & p)`. +Let \f$ S \f$ be the power sphere of the weighted points in range `[start,end)`. +The operator returns: +- `ON_ORIENTED_BOUNDARY` if `p` is orthogonal to +\f$ S \f$, + +- `ON_NEGATIVE_SIDE` if the power distance between `p` and \f$ S \f$ is +positive. + +- `ON_POSITIVE_SIDE` otherwise. + +\pre If `Dimension` is `CGAL::Dimension_tag`, +then `std::distance(start,end)=D+1`. +The weighted points in range +`[start,end)` must be affinely independent, i.e., the simplex must +not be flat. +*/ +typedef unspecified_type Power_side_of_power_sphere_d; + +/*! +A predicate object that must provide the templated operator +`template Oriented_side operator()(Flat_orientation_d orient, ForwardIterator start, ForwardIterator end, const Weighted_point_d & p)`. + +The points in range `[start,end)` and `p` are supposed to belong to the lower-dimensional flat +whose orientation is given by `orient`. + +Let \f$ S \f$ be the power sphere of the weighted points in range `[start,end)` +in this lower dimensional flat. +The operator returns: +- `ON_ORIENTED_BOUNDARY` if `p` is orthogonal to +\f$ S \f$, + +- `ON_NEGATIVE_SIDE` if the power distance between `p` and \f$ S \f$ is +positive. + +- `ON_POSITIVE_SIDE` otherwise. + +\pre `std::distance(start,end)=k+1` where \f$ k\f$ is the number of +points used to construct `orient` (dimension of the flat). +The points in range `[start,end)` must be affinely independent. +`p` must be in the flat generated by these points. +*/ +typedef unspecified_type In_flat_power_side_of_power_sphere_d; + +/// @} + +/// \name Creation +/// @{ + +/*! +The default constructor (optional). +This is not required when an instance of the traits is provided +to the constructor of `CGAL::Regular_triangulation`. +*/ +RegularTriangulationTraits(); + +/// @} + +/// \name Operations +/// The following methods permit access to the traits class's predicates and functors: +/// @{ + +/*! + +*/ +Construct_point_d construct_point_d_object() const; + +/*! + +*/ +Compute_weight_d compute_weight_d_object() const; + +/*! + +*/ +Power_side_of_power_sphere_d power_side_of_power_sphere_d_object() const; + +/*! + +*/ +In_flat_power_side_of_power_sphere_d in_flat_power_side_of_power_sphere_d_object() const; + +/// @} + +}; /* end RegularTriangulationTraits */ diff --git a/Triangulation/doc/Triangulation/Concepts/TriangulationDSFace.h b/Triangulation/doc/Triangulation/Concepts/TriangulationDSFace.h index dd0676f37db..c477323ee94 100644 --- a/Triangulation/doc/Triangulation/Concepts/TriangulationDSFace.h +++ b/Triangulation/doc/Triangulation/Concepts/TriangulationDSFace.h @@ -15,7 +15,7 @@ The dimension of a face is implicitely set when `TriangulationDSFace::set_index` is called two times to set the first two vertices (`i = 0` and `i = 1`), then the dimension is 1. -\cgalHasModel `CGAL::Triangulation_face` +\cgalHasModel `CGAL::Triangulation_face` \sa `TriangulationDSFullCell` \sa `TriangulationDSVertex` diff --git a/Triangulation/doc/Triangulation/Concepts/TriangulationDSFullCell.h b/Triangulation/doc/Triangulation/Concepts/TriangulationDSFullCell.h index f32f6e96c32..df05688010a 100644 --- a/Triangulation/doc/Triangulation/Concepts/TriangulationDSFullCell.h +++ b/Triangulation/doc/Triangulation/Concepts/TriangulationDSFullCell.h @@ -36,8 +36,8 @@ of `CGAL::Triangulation_data_structure::Vertex`. \cgalRefines `TriangulationDataStructure::FullCell` -\cgalHasModel ` CGAL::Triangulation_ds_full_cell` -\cgalHasModel `CGAL::Triangulation_full_cell` +\cgalHasModel `CGAL::Triangulation_ds_full_cell` +\cgalHasModel `CGAL::Triangulation_full_cell` \sa `TriangulationDSVertex` \sa `TriangulationDSFace` diff --git a/Triangulation/doc/Triangulation/Concepts/TriangulationDSVertex.h b/Triangulation/doc/Triangulation/Concepts/TriangulationDSVertex.h index 4cdd971a394..0366254ca89 100644 --- a/Triangulation/doc/Triangulation/Concepts/TriangulationDSVertex.h +++ b/Triangulation/doc/Triangulation/Concepts/TriangulationDSVertex.h @@ -36,8 +36,8 @@ of `CGAL::Triangulation_data_structure::Vertex`. \cgalRefines `TriangulationDataStructure::Vertex` -\cgalHasModel `CGAL::Triangulation_ds_vertex` -\cgalHasModel `CGAL::Triangulation_vertex` +\cgalHasModel `CGAL::Triangulation_ds_vertex` +\cgalHasModel `CGAL::Triangulation_vertex` \sa `TriangulationDSFullCell` \sa `TriangulationDSFace` diff --git a/Triangulation/doc/Triangulation/Concepts/TriangulationDataStructure.h b/Triangulation/doc/Triangulation/Concepts/TriangulationDataStructure.h index 7d20be8644d..3b1b98c991c 100644 --- a/Triangulation/doc/Triangulation/Concepts/TriangulationDataStructure.h +++ b/Triangulation/doc/Triangulation/Concepts/TriangulationDataStructure.h @@ -26,8 +26,8 @@ which is also the unique vertex and the unique full cell in the `TriangulationDataStructure`. In a geometric realization of the `TriangulationDataStructure` (e.g., in a -`Triangulation` or a -`Delaunay_triangulation`), this vertex +`Triangulation` or a +`Delaunay_triangulation`), this vertex corresponds to the vertex at infinity.
0
This corresponds to two vertices, each incident to one \f$ 0\f$-face; @@ -70,7 +70,7 @@ The classes `Vertex` and `Full_cell` have to provide the relevant I/O operators (possibly empty). -\cgalHasModel `CGAL::Triangulation_data_structure` +\cgalHasModel `CGAL::Triangulation_data_structure` \sa `TriangulationDataStructure::Vertex` \sa `TriangulationDataStructure::FullCell` @@ -257,11 +257,13 @@ The predicate must return `true` if the traversal of that `Facet` leads to a good full cell. All the good full cells are output into the last argument `out`. -\pre `start != Full_cell_handle()` and `start` is a good cell. +Returns a facet on the boundary of the set of cells. + +\pre `start != Full_cell_handle()` and `start` is a good cell. */ template< typename TraversalPredicate, typename OutputIterator > -void gather_full_cells(Full_cell_handle start, TraversalPredicate & tp, +Facet gather_full_cells(Full_cell_handle start, TraversalPredicate & tp, OutputIterator & out) const; /*! @@ -656,8 +658,8 @@ It sets requirements of combinatorial nature only, as geometry is not concerned here. In particular, we only require that the vertex holds a handle to a full cell incident to it in the triangulation. -\cgalHasModel `CGAL::Triangulation_ds_vertex` -\cgalHasModel `CGAL::Triangulation_vertex` +\cgalHasModel `CGAL::Triangulation_ds_vertex` +\cgalHasModel `CGAL::Triangulation_vertex` \sa `TriangulationDataStructure::FullCell` \sa `TriangulationDataStructure::Face` @@ -763,8 +765,8 @@ full cell as well as handles to the adjacent full cells. Two full cells are said to be adjacent when they share a facet. Adjacent full cells are called hereafter neighbors. -\cgalHasModel `CGAL::Triangulation_ds_full_cell` -\cgalHasModel `CGAL::Triangulation_full_cell` +\cgalHasModel `CGAL::Triangulation_ds_full_cell` +\cgalHasModel `CGAL::Triangulation_full_cell` \sa `TriangulationDataStructure::FullCell` \sa `TriangulationDataStructure::Face` diff --git a/Triangulation/doc/Triangulation/Concepts/TriangulationFullCell.h b/Triangulation/doc/Triangulation/Concepts/TriangulationFullCell.h index 4269b675b7d..8c0f1f76cf4 100644 --- a/Triangulation/doc/Triangulation/Concepts/TriangulationFullCell.h +++ b/Triangulation/doc/Triangulation/Concepts/TriangulationFullCell.h @@ -4,17 +4,17 @@ \cgalConcept The concept `TriangulationFullCell` describes the requirements on the type used by the -class `CGAL::Triangulation`, and its derived classes, to +class `CGAL::Triangulation`, and its derived classes, to represent a full cell. \cgalRefines `TriangulationDSFullCell` We only list below the additional specific requirements of `TriangulationFullCell`. -\cgalHasModel `CGAL::Triangulation_full_cell` +\cgalHasModel `CGAL::Triangulation_full_cell` -\sa `CGAL::Triangulation_full_cell` +\sa `CGAL::Triangulation_full_cell` \sa `TriangulationVertex` -\sa `CGAL::Triangulation` +\sa `CGAL::Triangulation` */ diff --git a/Triangulation/doc/Triangulation/Concepts/TriangulationTraits.h b/Triangulation/doc/Triangulation/Concepts/TriangulationTraits.h index d66dd83d7f1..9791e86d42f 100644 --- a/Triangulation/doc/Triangulation/Concepts/TriangulationTraits.h +++ b/Triangulation/doc/Triangulation/Concepts/TriangulationTraits.h @@ -5,7 +5,7 @@ This concept describes the geometric types and predicates required to build a triangulation. It corresponds to the first template parameter of the class -`CGAL::Triangulation`. +`CGAL::Triangulation`. \cgalRefines `SpatialSortingTraits_d` @@ -30,8 +30,8 @@ A type representing the dimension of the predicates (but not necessarily the one of `Point_d`). If \f$ n \f$ is the number of points required by the `Orientation_d` predicate, then `Dimension` \f$ = n - 1\f$. -It can be static (`Dimension`=`CGAL::``Dimension_tag`) or -dynamic (`Dimension`=`CGAL::``Dynamic_dimension_tag`). +It can be static (`Dimension`=`CGAL::Dimension_tag`) or +dynamic (`Dimension`=`CGAL::Dynamic_dimension_tag`). */ typedef unspecified_type Dimension; @@ -48,7 +48,7 @@ templated operator The operator returns the orientation of the simplex defined by the points in the range `[start, end)`; the value can be `CGAL::POSITIVE`, `CGAL::NEGATIVE` or `CGAL::COPLANAR`. -\pre If `Dimension`=`CGAL::``Dimension_tag`, then `std::distance(start,end)=D+1`. +\pre If `Dimension`=`CGAL::Dimension_tag`, then `std::distance(start,end)=D+1`. */ typedef unspecified_type Orientation_d; @@ -59,7 +59,7 @@ the templated operator The operator returns `true` if and only if point `p` is contained in the affine space spanned by the points in the range `[start, end)`. That affine space is also called the affine hull of the points in the range. -\pre If `Dimension`=`CGAL::``Dimension_tag`, +\pre If `Dimension`=`CGAL::Dimension_tag`, then `std::distance(start,end)=D+1`. The points in the range must be affinely independent. Note that in the CGAL kernels, this predicate @@ -97,7 +97,7 @@ the range `R=[start, end)` can be oriented in two different ways, the operator returns an object that allow to orient that flat so that `R=[start, end)` defines a positive simplex. -\pre If `Dimension`=`CGAL::``Dimension_tag`, +\pre If `Dimension`=`CGAL::Dimension_tag`, then `std::distance(start,end)=D+1`. The points in range `[start,end)` must be affinely independent. @@ -131,13 +131,15 @@ the same and `LARGER` otherwise. */ typedef unspecified_type Compare_lexicographically_d; -/// @} +/// @} /// \name Creation /// @{ /*! -The default constructor. +The default constructor (optional). +This is not required when an instance of the traits is provided +to the constructor of `CGAL::Triangulation`. */ TriangulationTraits(); diff --git a/Triangulation/doc/Triangulation/Concepts/TriangulationVertex.h b/Triangulation/doc/Triangulation/Concepts/TriangulationVertex.h index 2123bdefb02..bd564bebf08 100644 --- a/Triangulation/doc/Triangulation/Concepts/TriangulationVertex.h +++ b/Triangulation/doc/Triangulation/Concepts/TriangulationVertex.h @@ -4,7 +4,7 @@ \cgalConcept The concept `TriangulationVertex` describes the requirements on the type used by the -class `CGAL::Triangulation`, and its derived classes, to +class `CGAL::Triangulation`, and its derived classes, to represent a vertex. \cgalRefines `TriangulationDSVertex` @@ -12,7 +12,7 @@ We only list below the additional specific requirements of ::TriangulationVertex Compared to ::TriangulationDSVertex, the main difference is the addition of an association of the vertex with a geometric point. -\cgalHasModel `CGAL::Triangulation_vertex ` +\cgalHasModel `CGAL::Triangulation_vertex` Input/Output -------------- @@ -20,9 +20,9 @@ Input/Output These operators can be used directly and are called by the I/O operator of class `Triangulation`. -\sa `CGAL::Triangulation_vertex` +\sa `CGAL::Triangulation_vertex` \sa `TriangulationFullCell` -\sa `CGAL::Triangulation` +\sa `CGAL::Triangulation` */ @@ -36,7 +36,7 @@ public: The type of the point stored in the vertex. It must be the same as the point type `TriangulationTraits::Point` (or its refined concepts) when the `TriangulationVertex` is used in the class -`Triangulation` (or its derived classes). +`Triangulation` (or its derived classes). */ typedef unspecified_type Point; diff --git a/Triangulation/doc/Triangulation/PackageDescription.txt b/Triangulation/doc/Triangulation/PackageDescription.txt index 3b95eea1180..b3bf73dd86f 100644 --- a/Triangulation/doc/Triangulation/PackageDescription.txt +++ b/Triangulation/doc/Triangulation/PackageDescription.txt @@ -5,6 +5,9 @@ /// \defgroup PkgTriangulationsTriangulationClasses Triangulation Classes /// \ingroup PkgTriangulations +/// \defgroup PkgTriangulationsTraitsClasses Traits Classes +/// \ingroup PkgTriangulations + /// \defgroup PkgTriangulationsVertexCellClasses Vertex, Face and Cell Classes /// \ingroup PkgTriangulations @@ -13,7 +16,7 @@ \cgalPkgDescriptionBegin{dD Triangulations,PkgTriangulationsSummary} \cgalPkgPicture{Hypertriangle.png} \cgalPkgSummaryBegin -\cgalPkgAuthors{Samuel Hornus, Olivier Devillers and Clément Jamin} +\cgalPkgAuthors{Olivier Devillers, Samuel Hornus, and Clément Jamin} \cgalPkgDesc{This package provides classes for manipulating triangulations (pure simplicial complexes) in Euclidean spaces whose dimension can be specified at compile-time or at run-time. Specifically, it provides a @@ -84,6 +87,7 @@ is opposite to the vertex with the same index. - `TriangulationTraits` - `DelaunayTriangulationTraits` +- `RegularTriangulationTraits` - `TriangulationVertex` - `TriangulationFullCell` @@ -93,17 +97,22 @@ The latter two concepts are also abbreviated respectively as `TrVertex` and `TrF ## Triangulation Data Structure ## -- `CGAL::Triangulation_data_structure` -- `CGAL::Triangulation_ds_vertex` -- `CGAL::Triangulation_ds_full_cell` -- `CGAL::Triangulation_face` +- `CGAL::Triangulation_data_structure` +- `CGAL::Triangulation_ds_vertex` +- `CGAL::Triangulation_ds_full_cell` +- `CGAL::Triangulation_face` ## (Geometric) Triangulations ## -- `CGAL::Triangulation` -- `CGAL::Delaunay_triangulation` -- `CGAL::Triangulation_vertex` -- `CGAL::Triangulation_full_cell` +- `CGAL::Triangulation` +- `CGAL::Delaunay_triangulation` +- `CGAL::Regular_triangulation` +- `CGAL::Triangulation_vertex` +- `CGAL::Triangulation_full_cell` + +## Traits Classes ## + +- `CGAL::Regular_triangulation_traits_adapter` ## Enums ## diff --git a/Triangulation/doc/Triangulation/Triangulation.txt b/Triangulation/doc/Triangulation/Triangulation.txt index 2a814263f14..962ec2bfe97 100644 --- a/Triangulation/doc/Triangulation/Triangulation.txt +++ b/Triangulation/doc/Triangulation/Triangulation.txt @@ -6,17 +6,17 @@ namespace CGAL { \anchor Chapter_Triangulations \cgalAutoToc -\authors Samuel Hornus, Olivier Devillers and Clément Jamin. +\authors Olivier Devillers, Samuel Hornus, and Clément Jamin. This package proposes data structures and algorithms to compute -triangulations of points in any dimensions. +triangulations of points in any dimensions \cgalCite{boissonnat2009Delaunay}. The `Triangulation_data_structure` handles the combinatorial aspect of triangulations while the geometric classes -`Triangulation` and `Delaunay_triangulation` allows to -compute and maintain triangulations and Delaunay triangulations of -sets of points. +`Triangulation`, `Delaunay_triangulation` and `Regular_triangulation` allow to +compute and maintain triangulations, Delaunay triangulations, and +regular triangulations of sets of points. -# Introduction # +\section TriangulationSecIntro Introduction ## Some Definitions ## @@ -30,64 +30,28 @@ The sets in \f$ S\f$ (which are subsets of \f$ V\f$) are called singular of which is simplex). A simplex \f$ s\in S\f$ is maximal if it is not a proper subset of some other set in \f$ S\f$. -A simplex having \f$ d+1 \f$ vertices is said of dimension \f$ d \f$. +A simplex having \f$ k+1 \f$ vertices is said of dimension \f$ k \f$. +An \f$ k\f$-face denotes a \f$ k\f$-dimensional simplex, i.e., a simplex with \f$ k+1\f$ +vertices. The simplicial complex is pure if all the maximal simplices have the same dimension. + +A triangulation is a simplicial complex +that is pure, connected and without boundaries nor singularities. The +dimension of the triangulation is the dimension of its maximal +simplices. + In the sequel, we will call these maximal simplices full cells. A face of a simplex is a subset of this simplex. A proper face of a simplex is a strict subset of this simplex. +Two faces \f$ \sigma\f$ and \f$ \sigma'\f$ are incident if and only if +\f$ \sigma'\f$ is a proper face of \f$ \sigma\f$ or vice versa. A complex has no boundaries if any proper face of a simplex is also a proper face of another simplex. -If the vertices are embedded into Euclidean space \f$ \mathbb{R}^d\f$, -we deal with -finite simplicial complexes which have slightly different simplices -and additional requirements: -
    -
  • vertices corresponds to points in space. -
  • a simplex \f$ s\in S\f$ is the convex hull of its vertices. -
  • the vertices of a simplex \f$ s\in S\f$ are affinely independent. -
  • the intersection of any two simplices of \f$ S\f$ is a proper face of both -simplices (the empty set counts). -
-See the wikipedia -entry for more about simplicial complexes. - -## What's in this Package? ## - -This \cgal package provides three main classes -for creating and manipulating triangulations. - -The class `CGAL::Triangulation_data_structure` -models an abstract triangulation: vertices in this -class are not embedded in Euclidean space but are only of combinatorial -nature. It deals with simplicial complexes -which are pure, connected and without boundaries nor singularities. - -The class `CGAL::Triangulation` -describes an embedded triangulation that has as vertices a given set of points. -Methods are provided for the insertion of points in the triangulation, the -traversal of various elements of the triangulation, as well as the localization of a -query point inside the triangulation. -The triangulation covers the convex hull of the set of points. - -The class `CGAL::Delaunay_triangulation` -builds the Delaunay triangulation of a set of points. -In a Delaunay triangulation, each face has the so-called -Delaunay or empty-ball property: there exists a -circumscribing ball whose interior does not contain -any vertex of the triangulation. - -## Further Definitions ## - -An \f$ i\f$-face denotes an \f$ i\f$-dimensional simplex, or a simplex with \f$ i+1\f$ -vertices. When these vertices are embedded in Euclidean space, they must be -affinely independent. - -If the maximal dimension of a simplex in the triangulation is -\f$ d\f$, we use the following terminology:
    +If the triangulation is of dimension \f$ d \f$, we use the following terminology:
    • face: an \f$ i\f$-face for some \f$ i\in[0,d]\f$;
    • vertex: a \f$ 0\f$-face;
    • edge: a \f$ 1\f$-face; @@ -96,32 +60,69 @@ If the maximal dimension of a simplex in the triangulation is
    • full cell: a \f$ d\f$-face.
    -Two faces \f$ \sigma\f$ and \f$ \sigma'\f$ are incident if and only if -\f$ \sigma'\f$ is a proper sub-face of \f$ \sigma\f$ or vice versa. +If the vertices are embedded into Euclidean space \f$ \mathbb{R}^n\f$, +we deal with +finite simplicial complexes, which have slightly different simplices +and additional requirements: +
      +
    • vertices correspond to points in space. +
    • a simplex \f$ s\in S\f$ is the convex hull of its vertices. +
    • the vertices of a simplex \f$ s\in S\f$ are affinely independent. +
    • the intersection of any two simplices of \f$ S\f$ is a proper face of both +simplices (the empty set counts). +
    -# %Triangulation Data Structure # +## What is Provided in this Package? ## + +This \cgal package provides four main classes +for creating and manipulating triangulations. + +The class `CGAL::Triangulation_data_structure` +models an abstract triangulation: vertices in this +class are not embedded in Euclidean space but are only of combinatorial +nature. + +The class `CGAL::Triangulation` +describes an embedded triangulation that has as vertices a given set of points. +Methods are provided for the insertion of points in the triangulation, the +traversal of various elements of the triangulation, as well as the location of a +query point inside the triangulation. +The triangulation covers the convex hull of the set of points. + +The class `CGAL::Delaunay_triangulation` +builds the Delaunay triangulation of a set of points. +In a Delaunay triangulation, each face has the so-called +Delaunay or empty-ball property: there exists a +circumscribing ball whose interior does not contain +any vertex of the triangulation. + +The class `CGAL::Regular_triangulation` +builds the regular triangulation +-- also known as weighted Delaunay triangulation -- of a set of points. +A detailed definition of such a triangulation is available in section +\ref TriangulationSecRT. + +\section TriangulationSecTDS Triangulation Data Structure In this section, we describe the concept `TriangulationDataStructure` for which \cgal provides one model class: -`CGAL::Triangulation_data_structure`. +`CGAL::Triangulation_data_structure`. -A `TriangulationDataStructure` can represent an abstract pure complex -such that any facet is incident to exactly two full cells. +A triangulation data structure can represent an abstract triangulation. -A `TriangulationDataStructure` has a -maximal dimension which is a +The maximal dimension of a triangulation data structure is a positive integer equal to the maximum dimension a full cell can have. -This maximal dimension can be chosen by the user at the creation of a -`TriangulationDataStructure` and can then be queried using the method `tds.maximal_dimension()`. -A `TriangulationDataStructure` also knows the current dimension of its full cells, -which can be queried with `tds.current_dimension()`. In the sequel, let +This maximal dimension can be chosen by the user at the creation of the +triangulation data structure and can then be obtained using the method `tds.maximal_dimension()`. +A triangulation data structure also knows the current dimension of its full cells, +which can be obtained using `tds.current_dimension()`. In the sequel, let us denote the maximal dimension with \f$ D \f$ and the current dimension with \f$ d \f$. -The inequalities \f$ -2 \leq d \leq D\f$ and \f$ 0 \le D\f$ always hold. +The inequalities \f$ -2 \leq d \leq D\f$ and \f$ 0 < D\f$ always hold. The special meaning of negative values for \f$d\f$ is explained below. -### The Set of Faces ### +## The Set of Faces ## -The set of faces of a `TriangulationDataStructure` with +The set of faces of a triangulation data structure with current dimension \f$ d \f$ forms a triangulation of the topological sphere \f$ \mathbb{S}^d\f$. @@ -132,7 +133,7 @@ Possible values of \f$d\f$ (the current dimension of the triangulation) i
    \f$d=-2\f$
    This corresponds to an empty -`TriangulationDataStructure`. +triangulation data structure.
    \f$d=-1\f$
    This corresponds to an abstract simplicial complex reduced to a single vertex. @@ -149,16 +150,16 @@ the sphere \f$ \mathbb{S}^d\f$. ## The `Triangulation_data_structure` Class ## We give here some details about the class -`Triangulation_data_structure` +`Triangulation_data_structure` implementing the concept `TriangulationDataStructure`. ### Storage ### -A `TriangulationDataStructure` explicitly stores its vertices and full cells. +A triangulation data structure explicitly stores its vertices and full cells. Each vertex stores a reference to one of its incident full cells. -Each full cell stores references to its \f$ d+1\f$ vertices and +Each full cell \f$ \sigma \f$ stores references to its \f$ d+1\f$ vertices and neighbors. Its vertices and neighbors are indexed from \f$ 0\f$ to \f$ d \f$. The indices of its neighbors have the following meaning: the \f$ i\f$-th neighbor of \f$ \sigma\f$ is the unique neighbor of \f$ \sigma\f$ that does not contain the \f$ i\f$-th vertex of @@ -190,7 +191,7 @@ indices alongside the references to the vertices and neighbors in a full cell. This improves speed a little, but requires more memory. \cgalAdvanced \cgal provides the class template -`Triangulation_ds_full_cell` for representing full cells in a triangulation. Its second template parameter is used to specify wether or not the mirror indices should be kept in memory or computed @@ -200,41 +201,41 @@ documentation of that class template for specific details. ###Template Parameters### -The `Triangulation_data_structure` +The `Triangulation_data_structure` class is designed in such a way that its user can choose
    • the maximal dimension of the triangulation data structure by specifying the `Dimensionality` template parameter, -
    • the type used to represent vertices by specifying the `TriangulationDSVertex` +
    • the type used to represent vertices by specifying the `TriangulationDSVertex_` template parameter and
    • the type used to represent full cells by specifying the -`TriangulationDSFullCell` template parameter. +`TriangulationDSFullCell_` template parameter.
    The last two parameters have default values and are thus not necessary, unless -the user needs custom types (see the reference manual page for this class -template). The first template parameter, `Dimensionality`, must be -one of the following: +the user needs custom types (see `Triangulation_data_structure`). +The first template parameter, `Dimensionality`, must be one of the following:
      -
    • CGAL::`Dimension_tag` for some integer \f$ D \f$. This +
    • `CGAL::Dimension_tag` for some integer \f$ D \f$. This indicates that the triangulation can store full cells of dimension at most \f$ D \f$. The maximum dimension \f$ D \f$ is known by the compiler, which triggers some optimizations. -
    • CGAL::`Dynamic_dimension_tag`. In this case, the maximum +
    • `CGAL::Dynamic_dimension_tag`. In this case, the maximum dimension of the full cells must be passed as an integer argument to an instance constructor (see `TriangulationDataStructure`).
    -The `TriangulationDSVertex` and `TriangulationDSFullCell` parameters to the class template +The `TriangulationDSVertex_` and `TriangulationDSFullCell_` parameters to the class template must be models of the concepts `TriangulationDSVertex` and `TriangulationDSFullCell` respectively. \cgal provides models for these -concepts: `Triangulation_ds_vertex` and -`Triangulation_ds_full_cell`, which, as one -can see, take the `TriangulationDataStructure` as a template parameter in order to get access to -some nested types in `TriangulationDataStructure`. +concepts: `Triangulation_ds_vertex` and +`Triangulation_ds_full_cell`, which, as one +can see, take the triangulation data structure as a template parameter in order to get access to +some nested types in `TriangulationDataStructure_`. The default values are `CGAL::Triangulation_ds_vertex` and `CGAL::Triangulation_ds_full_cell` -where `TDS` is the current class `Triangulation_data_structure` +where `TDS` is the current class +`Triangulation_data_structure`. This creates a circular dependency, which we resolve in the same way as in the \cgal `Triangulation_2` and `Triangulation_3` packages (see Chapters \ref Chapter_2D_Triangulation_Data_Structure, \ref Chapter_2D_Triangulations, @@ -276,8 +277,7 @@ full cells adjacent to `c` are automatically subdivided to match the subdivision of the full cell `c`. The barycentric subdivision of `c` is obtained by enumerating all the faces of `c` in order of decreasing dimension, from the dimension of `c` to dimension 1, and inserting a new -vertex in each face. For the enumeration, we use a combination enumerator, -which is not documented, but provided in \cgal. +vertex in each face. \cgalFigureBegin{triangulationfigbarycentric,barycentric-subdivision.png} Barycentric subdivision in dimension \f$ d=2\f$. @@ -285,9 +285,9 @@ Barycentric subdivision in dimension \f$ d=2\f$. \cgalExample{barycentric_subdivision.cpp} -# Triangulations # +\section TriangulationSecTriangulations Triangulations -The class `CGAL::Triangulation` +The class `CGAL::Triangulation` maintains a triangulation embedded in Euclidean space. The triangulation covers the convex hull of the input points (the embedded vertices of the triangulation). @@ -300,39 +300,39 @@ Each infinite \f$ i\f$-simplex is incident to the infinite vertex and to an \f$ (i-1)\f$-simplex of the convex hull boundary. -See Chapters \ref Chapter_2D_Triangulations "2D Triangulations" and +See Chapters \ref Chapter_2D_Triangulations "2D Triangulations" or \ref Chapter_3D_Triangulations "3D Triangulations" for more details about infinite vertices and cells. Methods are provided for the insertion of points in the triangulation, the contraction of faces, the traversal of various elements of the triangulation -as well as the localization of a query point inside the triangulation. +as well as the location of a query point inside the triangulation. The ordering of the vertices of a full cell defines an orientation of that full cell. As long as no advanced class method is called, it is guaranteed -that all finite full cells have positive orientation. The infinite full -cells are oriented as if the infinite vertex was on the other side -of the hyperplane supported by the convex hull facets that the other points. +that all finite full cells have positive orientation. Each infinite full +cell is oriented as if its infinite vertex was on the side of +the hyperplane supported by its finite facet where there is no other point. ## Implementation ## -The class `CGAL::Triangulation` -stores a model of the concept `TriangulationDataStructure` which is +The class `CGAL::Triangulation` +stores a model of the concept `TriangulationDataStructure` that is instantiated with a vertex type that stores a point. -The template parameter `TriangulationTraits` must be a model of the concept -`TriangulationTraits` which provides the `Point` type as well +The template parameter `TriangulationTraits_` must be a model of the concept +`TriangulationTraits`, which provides the point type as well as various geometric predicates used by the `Triangulation` class. The `TriangulationTraits` concept includes a nested type -`TriangulationTraits::Dimension` which is the dimension of the predicates. -This dimension governs the number of points given as arguments to the -predicates. This type is either `CGAL::Dimension_tag` or -`CGAL::Dynamic_dimension_tag`. In any case, the dimension of the traits -must match the maximal dimension of the `TriangulationDataStructure`. +`TriangulationTraits::Dimension`. This dimension governs the number of points +given as arguments to the predicates. This type is either +`CGAL::Dimension_tag` or `CGAL::Dynamic_dimension_tag`. +In any case, the dimension of the traits +must match the maximal dimension of the triangulation data structure. -The template parameter `TriangulationDataStructure` must be a model of the concept +The template parameter `TriangulationDataStructure_` must be a model of the concept `TriangulationDataStructure` which provides the triangulation data structure as described in the previous section. @@ -382,11 +382,11 @@ One important difference between the two examples above is that the first uses visits only the infinite full cells but stores handles to them into the potentially big array infinite_full_cells. -# Delaunay Triangulations # +\section TriangulationSecDT Delaunay Triangulations -The class `CGAL::Delaunay_triangulation` derives from -`CGAL::Triangulation` -and represent Delaunay triangulations. +The class `CGAL::Delaunay_triangulation` derives from +`CGAL::Triangulation` +and represents Delaunay triangulations. A circumscribing ball of a simplex is a ball having all vertices of the simplex on its boundary. @@ -396,8 +396,8 @@ circumscribing ball whose interior does not contain any vertex of the triangulation. In case of degeneracies (co-spherical points) the triangulation is not -uniquely defined. Note however that the \cgal implementation computes a unique -triangulation even in these cases. +uniquely defined. Note however that the CGAL implementation computes a +unique triangulation even in these cases. When a new point `p` is inserted into a Delaunay triangulation, the full cells whose circumscribing ball contains `p` are said to @@ -409,20 +409,23 @@ in the conflict zone are removed, leaving a hole that contains `p`. That hole is ``star shaped'' around `p` and thus is re-triangulated using `p` as a center vertex. -Delaunay triangulations also support vertex removal. +Delaunay triangulations support insertion of points, removal of vertices, +and location of a query point inside the triangulation. +Note that inserting a large set of points at once is much faster +than inserting the same points one by one. ## Implementation ## -The class `CGAL::Delaunay_triangulation` derives from -`CGAL::Triangulation`. It thus stores a model of -the concept `TriangulationDataStructure` which is instantiated with a vertex +The class `CGAL::Delaunay_triangulation` derives from +`CGAL::Triangulation`. It thus stores a model of +the concept `TriangulationDataStructure`, which is instantiated with a vertex type that stores a geometric point and allows its retrieval. -The template parameter `DelaunayTriangulationTraits` must be a model of the concept +The template parameter `DelaunayTriangulationTraits_` must be a model of the concept `DelaunayTriangulationTraits` which provides the geometric `Point` type as well as various geometric predicates used by the `Delaunay_triangulation` class. The concept `DelaunayTriangulationTraits` refines the concept -`TriangulationTraits` by requiring a few other geometric predicates, necessary +`TriangulationTraits` by requiring a few additional geometric predicates, necessary for the computation of Delaunay triangulations. ## Examples ## @@ -438,55 +441,150 @@ retaining an efficient update of the Delaunay triangulation. \cgalExample{delaunay_triangulation.cpp} -# Complexity and Performances # +\section TriangulationSecRT Regular Triangulations -The current implementation locates points by walking in the -triangulation, and sorts the points with spatial sort to insert a -set of points. Thus the theoretical complexity are -\f$ O(n\log n)\f$ for inserting \f$ n\f$ random points and \f$ O(n^{\frac{1}{d}})\f$ -for inserting one point in a triangulation of \f$ n\f$ random points. -In the worst case, the expected complexity is -\f$ O(n^{\lceil\frac{d}{2}\rceil}+n\log n)\f$. +The class `CGAL::Regular_triangulation` derives from +`CGAL::Triangulation` +and represents regular triangulations. -We provide below (Figure \cgalFigureRef{Triangulationfigbenchmarks}) the +Regular triangulations are also known as weighted Delaunay triangulations. + +Let \f$ {S}^{(w)}\f$ be a set of weighted points in \f$ \mathbb{R}^D\f$. Let +\f$ {p}^{(w)}=(p,w_p), p\in\mathbb{R}^D, w_p\in\mathbb{R}\f$ and +\f$ {z}^{(w)}=(z,w_z), z\in\mathbb{R}^D, w_z\in\mathbb{R}\f$ +be two weighted points. +A weighted point \f$ {p}^{(w)}=(p,w_p)\f$ can also be seen as a sphere of +center \f$ p\f$ and radius \f$ \sqrt{w_p}\f$. +The power product (or power distance ) +between \f$ {p}^{(w)}\f$ and \f$ {z}^{(w)}\f$ is +defined as +\f[ \Pi({p}^{(w)},{z}^{(w)}) = {\|{p-z}\|^2-w_p-w_z} \f] +where \f$ \|{p-z}\|\f$ is the Euclidean distance between \f$ p\f$ and \f$ z\f$. +\f$ {p}^{(w)}\f$ and \f$ {z}^{(w)}\f$ +are said to be orthogonal if \f$ \Pi({p}^{(w)},{z}^{(w)}) += 0\f$. + +\f$D + 1\f$ weighted points have a unique common orthogonal weighted point +called the power sphere. A sphere \f$ {z}^{(w)}\f$ is said to be +regular if \f$ \forall {p}^{(w)}\in{S}^{(w)}, +\Pi({p}^{(w)},{z}^{(w)})\geq 0\f$. + +A triangulation of \f$ {S}^{(w)}\f$ is regular if the power spheres +of all simplices are regular. + +Note that as a result, some points can be hidden and do not result in vertices +in the triangulation. Those points are discarded and cannot be retrieved. + +A weighted point `p` is said to be in conflict +with a simplex `s` if it has a negative power distance to the power sphere of `s`. + +Regular triangulations support insertion of weighted points, +and location of a query point inside the triangulation. +Note that inserting a large set of points at once is much faster +than inserting the same points one by one. +\warning The removal of vertices is not supported yet. + + +## Implementation ## + +The class `CGAL::Regular_triangulation` derives from +`CGAL::Triangulation`. It thus stores a model of +the concept `TriangulationDataStructure_` which is instantiated with a vertex +type that stores a weighted point and allows its retrieval. + +The template parameter `RegularTriangulationTraits_` must be a model of the concept +`RegularTriangulationTraits`. It must provide the `%Weighted_point_d` +type as well as various geometric predicates used by the +`Regular_triangulation` class. +The concept `RegularTriangulationTraits` refines the concept +`TriangulationTraits`. + +## Example ## + +This simple example shows how to create a regular triangulation. + +\cgalExample{regular_triangulation.cpp} + +\section TriangulationSecPerf Complexity and Performances + +When inserting a batch of points into a Delaunay triangulation, +the current implementation starts by spatially sorting the points. +Then, for each point to insert, it locates it by walking in the triangulation, +using the previously inserted vertex as a "hint". Finally, the point is +inserted. +In the worst case scenario, without spatial sort, the expected complexity is +\f$ O(n^{\lceil\frac{d}{2}\rceil+1}) \f$. +When the algorithm is run on uniformly distributed points, the localization complexity is +\f$ O(n^{\frac{1}{d}}) \f$ and the size of the triangulation is \f$ O(n) \f$, which gives +a complexity of \f$ O(n^{1+\frac{1}{d}}) \f$ for the insertion. +With spatial sort and random points, one can expect a complexity of \f$ O(n\log n) \f$. +Please refer to \cgalCite{boissonnat2009Delaunay} for more details. + +We provide below (\cgalFigureRef{Triangulationfigbenchmarks100}, +\cgalFigureRef{Triangulationfigbenchmarks1000} and +\cgalFigureRef{triangulationfigbenchmarkchart}) the performance of the Delaunay triangulation on randomly distributed points. The machine used is a PC running Windows 7 64-bits with an Intel Xeon CPU clocked at 2.80 GHz with 32GB of RAM. -The program has been compiled with Microsoft Visual C++ 2012 in Release mode. - -\cgalFigureAnchor{Triangulationfigbenchmarks} +The program has been compiled with Microsoft Visual C++ 2013 in Release mode. +\cgalFigureAnchor{Triangulationfigbenchmarks100}
    - +
    + + + + + + + +
    Dimension23456789101112
    Time (s)0.0030.0070.030.140.562.711.3451856862390
    Memory (MB)< 1< 1< 113135318266221877156
    Number of maximal simplices1844871,5485,54819,59867,102230,375715,9842,570,6237,293,29321,235,615
    Number of convex hull facets14663081,1644,41016,97457,589238,406670,5452,574,3268,603,589
    +
    +
    +\cgalFigureCaptionBegin{Triangulationfigbenchmarks100} +Performance of the insertion of 100 points in a Delaunay triangulation. +\cgalFigureCaptionEnd + +\cgalFigureAnchor{Triangulationfigbenchmarks1000} +
    + - - + + + +
    Dimension2345678
    Inserting 100 points0.0030.0070.030.140.562.711.3
    Inserting 1000 points0.0150.0560.523.526.21851385
    Time (s)0.010.050.53.4241831365
    Memory (MB)< 1< 12.714814832827
    Number of maximal simplices1,9796,31525,845122,116596,9273,133,31816,403,337
    Number of convex hull facets191389636,18441,135241,5401,406,797
    -\cgalFigureCaptionBegin{Triangulationfigbenchmarks} -Running times in seconds for algorithms on Delaunay triangulations. +\cgalFigureCaptionBegin{Triangulationfigbenchmarks1000} +Performance of the insertion of 1000 points in a Delaunay triangulation. \cgalFigureCaptionEnd -# Design and Implementation History # +\cgalFigureBegin{triangulationfigbenchmarkchart,benchmark_DTd.png} +Running time wrt. number of maximal simplices, for dimensions for 2 to 12. +\cgalFigureEnd -This package is heavily inspired by the works of -Monique Teillaud and Sylvain Pion (`Triangulation_3`) -and Mariette Yvinec (`Triangulation_2`). -The first version was written by Samuel Hornus. The final version is a joint -work by Samuel Hornus, Olivier Devillers and Clément Jamin. - -Clément Jamin's work was supported by the -Advanced Grant of the European Research Council GUDHI -(Geometric Understanding in Higher Dimensions). +\section TriangulationSecDesign Design and Implementation History Starting with the version 2.3 of \cgal, a package written by Susan Hert and Michael Seel was the first able to deal with triangulation and convex hulls in arbitrary dimension. It is deprecated since the version 4.6 of \cgal and this package should be used instead. +This package is heavily inspired by the works of +Monique Teillaud and Sylvain Pion (`Triangulation_3`) +and Mariette Yvinec (`Triangulation_2`). +The first version was written by Samuel Hornus. The final version is a joint +work by Samuel Hornus, Olivier Devillers and Clément Jamin. In 2017, Clément +Jamin added the regular triangulations. + +Clément Jamin's work was supported by the +Advanced Grant of the European Research Council GUDHI +(Geometric Understanding in Higher Dimensions). + + */ } /* namespace CGAL */ diff --git a/Triangulation/doc/Triangulation/examples.txt b/Triangulation/doc/Triangulation/examples.txt index 7837f9302d2..36db1597b8f 100644 --- a/Triangulation/doc/Triangulation/examples.txt +++ b/Triangulation/doc/Triangulation/examples.txt @@ -1,6 +1,7 @@ /*! \example barycentric_subdivision.cpp \example delaunay_triangulation.cpp +\example regular_triangulation.cpp \example triangulation.cpp \example triangulation1.cpp \example triangulation2.cpp diff --git a/Triangulation/doc/Triangulation/fig/benchmark_DTd.png b/Triangulation/doc/Triangulation/fig/benchmark_DTd.png new file mode 100644 index 00000000000..76ccd5144a8 Binary files /dev/null and b/Triangulation/doc/Triangulation/fig/benchmark_DTd.png differ diff --git a/Triangulation/examples/Triangulation/CMakeLists.txt b/Triangulation/examples/Triangulation/CMakeLists.txt index 910fcacfb5e..a1999ba3863 100644 --- a/Triangulation/examples/Triangulation/CMakeLists.txt +++ b/Triangulation/examples/Triangulation/CMakeLists.txt @@ -22,6 +22,8 @@ if ( CGAL_FOUND ) create_single_source_cgal_program( "barycentric_subdivision.cpp" ) create_single_source_cgal_program( "delaunay_triangulation.cpp" ) + create_single_source_cgal_program( "convex_hull.cpp" ) + create_single_source_cgal_program( "regular_triangulation.cpp" ) create_single_source_cgal_program( "triangulation.cpp" ) create_single_source_cgal_program( "triangulation_data_structure_dynamic.cpp" ) create_single_source_cgal_program( "triangulation_data_structure_static.cpp" ) diff --git a/Triangulation/examples/Triangulation/barycentric_subdivision.cpp b/Triangulation/examples/Triangulation/barycentric_subdivision.cpp index 13fba765d3c..08a5ae80ba0 100644 --- a/Triangulation/examples/Triangulation/barycentric_subdivision.cpp +++ b/Triangulation/examples/Triangulation/barycentric_subdivision.cpp @@ -1,5 +1,5 @@ #include -#include +#include #include #include @@ -34,8 +34,8 @@ void barycentric_subdivide(TDS & tds, typename TDS::Full_cell_handle fc) face_vertices.resize(d+1); // The following class // enumerates all (d+1)-tuple of the set {0, 1, ..., dim} - CGAL::internal::Combination_enumerator combi(d+1, 0, dim); - while( ! combi.end() ) + CGAL::Combination_enumerator combi(d+1, 0, dim); + while ( !combi.finished() ) { for( int i = 0; i <= d; ++i ) face_vertices[i] = vertices[combi[i]]; diff --git a/Triangulation/examples/Triangulation/convex_hull.cpp b/Triangulation/examples/Triangulation/convex_hull.cpp new file mode 100644 index 00000000000..ba7a85a59d9 --- /dev/null +++ b/Triangulation/examples/Triangulation/convex_hull.cpp @@ -0,0 +1,70 @@ +#include +#include +#include +#include +#include +#include + +#include +#include +#include + +const int D = 10; +typedef CGAL::Epick_d< CGAL::Dimension_tag > K; +typedef CGAL::Delaunay_triangulation T; +// The triangulation uses the default instanciation of the +// TriangulationDataStructure template parameter + +int main(int argc, char **argv) +{ + int N = 100; // number of points + if (argc > 1) + N = atoi(argv[1]); + + CGAL::Timer cost; // timer + + // Generate N random points + typedef CGAL::Random_points_in_cube_d Random_points_iterator; + Random_points_iterator rand_it(D, 1.0, CGAL::get_default_random()); + std::vector points; + CGAL::cpp11::copy_n(rand_it, N, std::back_inserter(points)); + + T t(D); + CGAL_assertion(t.empty()); + + // insert the points in the triangulation, only if they are outside the + // convex hull + std::cout << " Convex hull of "<::iterator it_p = points.begin() ; + it_p != points.end() ; ++it_p) + { + T::Locate_type lt; + T::Face f(t.maximal_dimension()); + T::Facet ft; + T::Full_cell_handle fch = t.locate(*it_p, lt, f, ft, hint); + if (lt == T::OUTSIDE_CONVEX_HULL || lt == T::OUTSIDE_AFFINE_HULL) + { + hint = t.insert(*it_p, lt, f, ft, fch)->full_cell(); + ++c; + } + else + { + hint = fch; + } + } + + std::cout << " done in " << cost.time() << " seconds.\n"; + std::cout << c << " points where actually inserted.\n"; + CGAL_assertion( t.is_valid() ); + + return 0; +} diff --git a/Triangulation/examples/Triangulation/delaunay_triangulation.cpp b/Triangulation/examples/Triangulation/delaunay_triangulation.cpp index 0ab829b01f6..73fafb3d1c7 100644 --- a/Triangulation/examples/Triangulation/delaunay_triangulation.cpp +++ b/Triangulation/examples/Triangulation/delaunay_triangulation.cpp @@ -9,69 +9,42 @@ int main() #else #include -#include #include -#include -#include -#include -#include -#include -#include - -const int D=5; -typedef CGAL::Epick_d< CGAL::Dimension_tag > K; -typedef CGAL::Delaunay_triangulation T; -// The triangulation uses the default instanciation of the -// TriangulationDataStructure template parameter - -int main(int argc, char **argv) +int main() { - int N = 100; if( argc > 2 )N = atoi(argv[1]); // number of points - CGAL::Timer cost; // timer + double pointsIn[][7] = { + { 42.89, 0, 60.55, 30.72, 0, 0, 0 }, + { 45.65, 50.83, 50.37, 16.13, 0, 0, 0 }, + { 79.06, 57.84, 61.59, 2.52, 0, 0, 0 }, + { 44.47, 39.46, 39.53, 28.72, 0, 0, 0 }, + { 0, 100, 0, 0, 100, 0, 53.47 }, + { 66.95, 100, 33.6, 0, 0, 0, 0 }, + { 42.89, 0, 0, 30.72, 100, 0, 53.47 }, + { 100, 100, 100, 100, 100, 100, 100 } + }; + + typedef CGAL::Triangulation > > T; + T dt(7); - // Instanciate a random point generator - CGAL::Random rng(0); - typedef CGAL::Random_points_in_cube_d Random_points_iterator; - Random_points_iterator rand_it(D, 1.0, rng); - // Generate N random points std::vector points; - CGAL::cpp11::copy_n(rand_it, N, std::back_inserter(points)); - - T t(D); - CGAL_assertion(t.empty()); - - // insert the points in the triangulation - cost.reset();cost.start(); - std::cout << " Delaunay triangulation of "< Full_cells; - Full_cells zone, new_full_cells; - std::back_insert_iterator out(zone); - c = t.locate(*++rand_it, lt, f, ft, v); - // previously inserted vertex v is used as hint for point location (if defined) - T::Facet ftc = t.compute_conflict_zone(*rand_it, c, out); - std::cout< "<::iterator it = points.begin(); it != points.end(); ++it) { + if (T::Vertex_handle() != hint) { + hint = dt.insert(*it, hint); + } + else { + hint = dt.insert(*it); + } + printf("Processing: %d/%d\n", ++i, (int)points.size()); + } return 0; } diff --git a/Triangulation/examples/Triangulation/regular_triangulation.cpp b/Triangulation/examples/Triangulation/regular_triangulation.cpp new file mode 100644 index 00000000000..aa70e175d55 --- /dev/null +++ b/Triangulation/examples/Triangulation/regular_triangulation.cpp @@ -0,0 +1,41 @@ +#include +#include +#include +#include + +#include +#include +#include + +const int D = 5; // Dimension +const int N = 100; // Number of points +typedef CGAL::Epick_d< CGAL::Dimension_tag > K; +typedef CGAL::Regular_triangulation T; +typedef T::Bare_point Bare_point; +typedef T::Weighted_point Weighted_point; + +int main() +{ + // Instantiate a random point generator + CGAL::Random rng(0); + typedef CGAL::Random_points_in_cube_d Random_points_iterator; + Random_points_iterator rand_it(D, 1.0, rng); + + // Generate N random points + std::vector points; + for (int i = 0; i < N; ++i) + points.push_back(Weighted_point(*rand_it++, rng.get_double(0., 10.))); + + T t(D); + CGAL_assertion(t.empty()); + + // Insert the points in the triangulation + t.insert(points.begin(), points.end()); + CGAL_assertion( t.is_valid() ); + std::cout << "Regular triangulation successfully computed: " + << t.number_of_vertices() << " vertices, " + << t.number_of_finite_full_cells() << " finite cells." + << std::endl; + + return 0; +} diff --git a/Triangulation/include/CGAL/Delaunay_triangulation.h b/Triangulation/include/CGAL/Delaunay_triangulation.h index 3897b78e4ec..2bdc52a3f48 100644 --- a/Triangulation/include/CGAL/Delaunay_triangulation.h +++ b/Triangulation/include/CGAL/Delaunay_triangulation.h @@ -59,7 +59,7 @@ class Delaunay_triangulation public: // PUBLIC NESTED TYPES typedef DCTraits Geom_traits; - typedef typename Base::Triangulation_ds Triangulation_ds; + typedef typename Base::Triangulation_ds Triangulation_ds; typedef typename Base::Vertex Vertex; typedef typename Base::Full_cell Full_cell; @@ -75,21 +75,27 @@ public: // PUBLIC NESTED TYPES typedef typename Base::Vertex_const_handle Vertex_const_handle; typedef typename Base::Vertex_const_iterator Vertex_const_iterator; - typedef typename Base::Full_cell_handle Full_cell_handle; - typedef typename Base::Full_cell_iterator Full_cell_iterator; - typedef typename Base::Full_cell_const_handle Full_cell_const_handle; - typedef typename Base::Full_cell_const_iterator Full_cell_const_iterator; + typedef typename Base::Full_cell_handle Full_cell_handle; + typedef typename Base::Full_cell_iterator Full_cell_iterator; + typedef typename Base::Full_cell_const_handle Full_cell_const_handle; + typedef typename Base::Full_cell_const_iterator Full_cell_const_iterator; + typedef typename Base::Finite_full_cell_const_iterator + Finite_full_cell_const_iterator; typedef typename Base::size_type size_type; typedef typename Base::difference_type difference_type; typedef typename Base::Locate_type Locate_type; + //Tag to distinguish triangulations with weighted_points + typedef Tag_false Weighted_tag; + protected: // DATA MEMBERS public: - + + using typename Base::Rotor; using Base::maximal_dimension; using Base::are_incident_full_cells_valid; using Base::coaffine_orientation_predicate; @@ -99,11 +105,12 @@ public: //using Base::incident_full_cells; using Base::geom_traits; using Base::index_of_covertex; + //using Base::index_of_second_covertex; using Base::infinite_vertex; + using Base::rotate_rotor; using Base::insert_in_hole; using Base::insert_outside_convex_hull_1; using Base::is_infinite; - using Base::is_valid; using Base::locate; using Base::points_begin; using Base::set_neighbors; @@ -115,10 +122,12 @@ public: using Base::full_cell; using Base::full_cells_begin; using Base::full_cells_end; + using Base::finite_full_cells_begin; + using Base::finite_full_cells_end; using Base::vertices_begin; using Base::vertices_end; // using Base:: - + private: //*** Side_of_oriented_subsphere_d *** typedef typename Base::Flat_orientation_d Flat_orientation_d; @@ -147,49 +156,22 @@ private: }; public: -// - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - UTILITIES - - // A co-dimension 2 sub-simplex. called a Rotor because we can rotate - // the two "covertices" around the sub-simplex. Useful for traversing the - // boundary of a hole. NOT DOCUMENTED - typedef cpp11::tuple Rotor; - - /*Full_cell_handle full_cell(const Rotor & r) const // NOT DOCUMENTED - { - return cpp11::get<0>(r); - } - int index_of_covertex(const Rotor & r) const // NOT DOCUMENTED - { - return cpp11::get<1>(r); - } - int index_of_second_covertex(const Rotor & r) const // NOT DOCUMENTED - { - return cpp11::get<2>(r); - }*/ - Rotor rotate_rotor(Rotor & r) // NOT DOCUMENTED... - { - int opposite = cpp11::get<0>(r)->mirror_index(cpp11::get<1>(r)); - Full_cell_handle s = cpp11::get<0>(r)->neighbor(cpp11::get<1>(r)); - int new_second = s->index(cpp11::get<0>(r)->vertex(cpp11::get<2>(r))); - return Rotor(s, new_second, opposite); - } - // - - - - - - - - - - - - - - - - - - - - - - - - - - CREATION / CONSTRUCTORS - Delaunay_triangulation(int dim, const Geom_traits k = Geom_traits()) + Delaunay_triangulation(int dim, const Geom_traits &k = Geom_traits()) : Base(dim, k) { } // With this constructor, - // the user can specify a Flat_orientation_d object to be used for - // orienting simplices of a specific dimension + // the user can specify a Flat_orientation_d object to be used for + // orienting simplices of a specific dimension // (= preset_flat_orientation_.first) // It it used by the dark triangulations created by DT::remove Delaunay_triangulation( - int dim, + int dim, const std::pair &preset_flat_orientation, - const Geom_traits k = Geom_traits()) + const Geom_traits &k = Geom_traits()) : Base(dim, preset_flat_orientation, k) { } @@ -202,8 +184,8 @@ public: Side_of_oriented_subsphere_d side_of_oriented_subsphere_predicate() const { return Side_of_oriented_subsphere_d ( - flat_orientation_, - geom_traits().construct_flat_orientation_d_object(), + flat_orientation_, + geom_traits().construct_flat_orientation_d_object(), geom_traits().in_flat_side_of_oriented_sphere_d_object() ); } @@ -250,8 +232,8 @@ public: } return number_of_vertices() - n; } - Vertex_handle insert(const Point &, const Locate_type, const Face &, const Facet &, const Full_cell_handle); - Vertex_handle insert(const Point & p, const Full_cell_handle start = Full_cell_handle()) + Vertex_handle insert(const Point &, Locate_type, const Face &, const Facet &, Full_cell_handle); + Vertex_handle insert(const Point & p, Full_cell_handle start = Full_cell_handle()) { Locate_type lt; Face f(maximal_dimension()); @@ -259,13 +241,13 @@ public: Full_cell_handle s = locate(p, lt, f, ft, start); return insert(p, lt, f, ft, s); } - Vertex_handle insert(const Point & p, const Vertex_handle hint) + Vertex_handle insert(const Point & p, Vertex_handle hint) { CGAL_assertion( Vertex_handle() != hint ); return insert(p, hint->full_cell()); } Vertex_handle insert_outside_affine_hull(const Point &); - Vertex_handle insert_in_conflicting_cell(const Point &, const Full_cell_handle); + Vertex_handle insert_in_conflicting_cell(const Point &, Full_cell_handle); // - - - - - - - - - - - - - - - - - - - - - - - - - GATHERING CONFLICTING SIMPLICES @@ -275,7 +257,7 @@ public: Full_cell_const_handle, const OrientationPredicate &) const; template< typename OutputIterator > - Facet compute_conflict_zone(const Point &, const Full_cell_handle, OutputIterator) const; + Facet compute_conflict_zone(const Point &, Full_cell_handle, OutputIterator) const; template < typename OrientationPredicate, typename SideOfOrientedSpherePredicate > class Conflict_predicate @@ -315,7 +297,7 @@ public: Orientation o = ori_( boost::make_transform_iterator(s->vertices_begin(), spivi), - boost::make_transform_iterator(s->vertices_begin() + cur_dim_ + 1, + boost::make_transform_iterator(s->vertices_begin() + cur_dim_ + 1, spivi)); if( POSITIVE == o ) @@ -345,6 +327,10 @@ public: } }; +// - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - VALIDITY + + bool is_valid(bool verbose = false, int level = 0) const; + private: // Some internal types to shorten notation typedef typename Base::Coaffine_orientation_d Coaffine_orientation_d; @@ -357,30 +343,9 @@ private: Conflict_traversal_pred_in_subspace; typedef Conflict_traversal_predicate Conflict_traversal_pred_in_fullspace; - - // This is used in the |remove(v)| member function to manage sets of Full_cell_handles - template< typename FCH > - struct Full_cell_set : public std::vector - { - typedef std::vector Base_set; - using Base_set::begin; - using Base_set::end; - void make_searchable() - { // sort the full cell handles - std::sort(begin(), end()); - } - bool contains(const FCH & fch) const - { - return std::binary_search(begin(), end(), fch); - } - bool contains_1st_and_not_2nd(const FCH & fst, const FCH & snd) const - { - return ( ! contains(snd) ) && ( contains(fst) ); - } - }; }; -// = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = +// = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = // FUNCTIONS THAT ARE MEMBER METHODS: // - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - REMOVALS @@ -407,24 +372,14 @@ Delaunay_triangulation return Full_cell_handle(); } Full_cell_handle left = v->full_cell(); - if( is_infinite(left) && left->neighbor(0)->index(left) == 0 ) // we are on the infinite right. - left = left->neighbor(0); if( 0 == left->index(v) ) left = left->neighbor(1); CGAL_assertion( 1 == left->index(v) ); Full_cell_handle right = left->neighbor(0); - if( ! is_infinite(right) ) - { + tds().associate_vertex_with_full_cell(left, 1, right->vertex(1)); set_neighbors(left, 0, right->neighbor(0), right->mirror_index(0)); - } - else - { - tds().associate_vertex_with_full_cell(left, 1, left->vertex(0)); - tds().associate_vertex_with_full_cell(left, 0, infinite_vertex()); - set_neighbors(left, 0, left->neighbor(1), left->mirror_index(1)); - set_neighbors(left, 1, right->neighbor(1), right->mirror_index(1)); - } + tds().delete_vertex(v); tds().delete_full_cell(right); return left; @@ -432,7 +387,7 @@ Delaunay_triangulation // THE CASE cur_dim >= 2 // Gather the finite vertices sharing an edge with |v| - typedef Full_cell_set Simplices; + typedef typename Base::template Full_cell_set Simplices; Simplices simps; std::back_insert_iterator out(simps); tds().incident_full_cells(v, out); @@ -513,24 +468,16 @@ Delaunay_triangulation { int v_idx = (*it)->index(v); tds().associate_vertex_with_full_cell(*it, v_idx, infinite_vertex()); - if( v_idx != 0 ) - { - // we must put the infinite vertex at index 0. - // OK, now with the new convention that the infinite vertex - // does not have to be at index 0, this is not necessary, - // but still, I prefer to keep this piece of code here. [-- Samuel Hornus] - (*it)->swap_vertices(0, v_idx); - // Now, we preserve the positive orientation of the full_cell - (*it)->swap_vertices(current_dimension() - 1, current_dimension()); } - } // Make the handles to infinite full cells searchable infinite_simps.make_searchable(); // Then, modify the neighboring relation for( typename Simplices::iterator it = simps.begin(); it != simps.end(); ++it ) { - for( int i = 1; i <= current_dimension(); ++i ) + for( int i = 0; i <= current_dimension(); ++i ) { + if (is_infinite((*it)->vertex(i))) + continue; (*it)->vertex(i)->set_full_cell(*it); Full_cell_handle n = (*it)->neighbor(i); // Was |n| a finite full cell prior to removing |v| ? @@ -568,10 +515,10 @@ Delaunay_triangulation Dark_s_handle dark_ret_s = dark_s; Full_cell_handle ret_s; - typedef Full_cell_set Dark_full_cells; + typedef typename Base::template Full_cell_set Dark_full_cells; Dark_full_cells conflict_zone; std::back_insert_iterator dark_out(conflict_zone); - + dark_ft = dark_side.compute_conflict_zone(v->point(), dark_s, dark_out); // Make the dark simplices in the conflict zone searchable conflict_zone.make_searchable(); @@ -669,7 +616,7 @@ Delaunay_triangulation int li = light_s->index(dark_s->vertex(di)->data()); Rotor light_r(light_s, li, light_i); typename Dark_triangulation::Rotor dark_r(dark_s, di, dark_i); - + while (simps.contains(cpp11::get<0>(light_r)->neighbor(cpp11::get<1>(light_r)))) light_r = rotate_rotor(light_r); @@ -722,7 +669,7 @@ Delaunay_triangulation template< typename DCTraits, typename TDS > typename Delaunay_triangulation::Vertex_handle Delaunay_triangulation -::insert(const Point & p, const Locate_type lt, const Face & f, const Facet &, const Full_cell_handle s) +::insert(const Point & p, Locate_type lt, const Face & f, const Facet &, Full_cell_handle s) { switch( lt ) { @@ -753,6 +700,15 @@ Delaunay_triangulation } } +/* +[Undocumented function] + +Inserts the point `p` in the Delaunay triangulation. Returns a handle to the +(possibly newly created) vertex at that position. +\pre The point `p` +must lie outside the affine hull of the Delaunay triangulation. This implies that +`dt`.`current_dimension()` must be less than `dt`.`maximal_dimension()`. +*/ template< typename DCTraits, typename TDS > typename Delaunay_triangulation::Vertex_handle Delaunay_triangulation @@ -774,15 +730,53 @@ Delaunay_triangulation CGAL_assertion( ZERO != o ); if( NEGATIVE == o ) reorient_full_cells(); + + // We just inserted the second finite point and the right infinite + // cell is like : (inf_v, v), but we want it to be (v, inf_v) to be + // consistent with the rest of the cells + if (current_dimension() == 1) + { + // Is "inf_v_cell" the right infinite cell? + // Then inf_v_index should be 1 + if (inf_v_cell->neighbor(inf_v_index)->index(inf_v_cell) == 0 + && inf_v_index == 0) + { + inf_v_cell->swap_vertices( + current_dimension() - 1, current_dimension()); + } + // Otherwise, let's find the right infinite cell + else + { + inf_v_cell = inf_v_cell->neighbor((inf_v_index + 1) % 2); + inf_v_index = inf_v_cell->index(infinite_vertex()); + // Is "inf_v_cell" the right infinite cell? + // Then inf_v_index should be 1 + if (inf_v_cell->neighbor(inf_v_index)->index(inf_v_cell) == 0 + && inf_v_index == 0) + { + inf_v_cell->swap_vertices( + current_dimension() - 1, current_dimension()); + } + } + } } return v; } +/*! +[Undocumented function] + +Inserts the point `p` in the Delaunay triangulation. Returns a handle to the +(possibly newly created) vertex at that position. +\pre The point `p` must be in conflict with the full cell `c`. +*/ template< typename DCTraits, typename TDS > typename Delaunay_triangulation::Vertex_handle Delaunay_triangulation -::insert_in_conflicting_cell(const Point & p, const Full_cell_handle s) +::insert_in_conflicting_cell(const Point & p, Full_cell_handle s) { + CGAL_precondition(is_in_conflict(p, s)); + // for storing conflicting full_cells. typedef std::vector Full_cell_h_vector; CGAL_STATIC_THREAD_LOCAL_VARIABLE(Full_cell_h_vector,cs,0); @@ -876,7 +870,7 @@ template< typename DCTraits, typename TDS > template< typename OutputIterator > typename Delaunay_triangulation::Facet Delaunay_triangulation -::compute_conflict_zone(const Point & p, const Full_cell_handle s, OutputIterator out) const +::compute_conflict_zone(const Point & p, Full_cell_handle s, OutputIterator out) const { CGAL_precondition( 2 <= current_dimension() ); if( current_dimension() < maximal_dimension() ) @@ -895,6 +889,48 @@ Delaunay_triangulation } } +// - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - VALIDITY + +template< typename DCTraits, typename TDS > +bool +Delaunay_triangulation +::is_valid(bool verbose, int level) const +{ + if (!Base::is_valid(verbose, level)) + return false; + + int dim = current_dimension(); + if (dim == maximal_dimension()) + { + for (Finite_full_cell_const_iterator cit = this->finite_full_cells_begin() ; + cit != this->finite_full_cells_end() ; ++cit ) + { + Full_cell_const_handle ch = cit.base(); + for(int i = 0; i < dim+1 ; ++i ) + { + // If the i-th neighbor is not an infinite cell + Vertex_handle opposite_vh = + ch->neighbor(i)->vertex(ch->neighbor(i)->index(ch)); + if (!is_infinite(opposite_vh)) + { + Side_of_oriented_sphere_d side = + geom_traits().side_of_oriented_sphere_d_object(); + if (side(Point_const_iterator(ch->vertices_begin()), + Point_const_iterator(ch->vertices_end()), + opposite_vh->point()) == ON_BOUNDED_SIDE) + { + if (verbose) + CGAL_warning_msg(false, "Non-empty sphere"); + return false; + } + } + } + } + } + return true; +} + + } //namespace CGAL #endif // CGAL_DELAUNAY_COMPLEX_H diff --git a/Triangulation/include/CGAL/IO/Triangulation_off_ostream.h b/Triangulation/include/CGAL/IO/Triangulation_off_ostream.h new file mode 100644 index 00000000000..701f08202cd --- /dev/null +++ b/Triangulation/include/CGAL/IO/Triangulation_off_ostream.h @@ -0,0 +1,320 @@ +// Copyright (c) 2014 INRIA Sophia-Antipolis (France). +// All rights reserved. +// +// This file is part of CGAL (www.cgal.org); you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public License as +// published by the Free Software Foundation; either version 3 of the License, +// or (at your option) any later version. +// +// Licensees holding a valid commercial license may use this file in +// accordance with the commercial license agreement provided with the software. +// +// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE +// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE. +// +// $URL: $ +// $Id: $ +// +// Author(s) : Clement Jamin + + +#ifndef CGAL_TRIANGULATION_IO_H +#define CGAL_TRIANGULATION_IO_H + +#include +#include +#include +#include + +namespace CGAL { + +namespace Triangulation_IO +{ +// TODO: test if the stream is binary or text? +template +int +output_point(std::ostream & os, const Traits &traits, const P & p) +{ + typedef typename Traits::Compute_coordinate_d Ccd; + const Ccd ccd = traits.compute_coordinate_d_object(); + const int dim = traits.point_dimension_d_object()(p); + if (dim > 0) + { + os << ccd(p, 0); + for (int i = 1 ; i < dim ; ++i) + os << " " << CGAL::to_double(ccd(p, i)); + } + return dim; +} + +// TODO: test if the stream is binary or text? +template +int +output_weighted_point(std::ostream & os, const Traits &traits, const P & p, + bool output_weight = true) +{ + typedef typename Traits::Compute_coordinate_d Ccd; + typename Traits::Construct_point_d cp = + traits.construct_point_d_object(); + typename Traits::Compute_weight_d pt_weight = traits.compute_weight_d_object(); + const Ccd ccd = traits.compute_coordinate_d_object(); + const int dim = traits.point_dimension_d_object()(p); + if (dim > 0) + { + output_point(os, traits, p); + if (output_weight) + os << " " << pt_weight(p); + } + return dim; +} + +// TODO: test if the stream is binary or text? +template +void +output_full_cell(std::ostream & os, const Traits &traits, const FCH & fch, + bool output_weights = false) +{ + typename FCH::value_type::Vertex_handle_iterator vit = fch->vertices_begin(); + for( ; vit != fch->vertices_end(); ++vit ) + { + int dim; + if (output_weights) + dim = output_weighted_point(os, traits, (*vit)->point()); + else + dim = output_point(os, traits, (*vit)->point()); + if (dim > 0) + os << std::endl; + } +} + +// TODO: test if the stream is binary or text? +/*template +void +input_point(std::istream & is, const Traits &traits, P & p) +{ + typedef typename Traits::FT FT; + std::vector coords; + + std::string line; + for(;;) + { + if (!std::getline(is, line)) + return is; + if (line != "") + break; + } + std::stringstream line_sstr(line); + FT temp; + while (line_sstr >> temp) + coords.push_back(temp); + + p = traits.construct_point_d_object()(coords.begin(), coords.end()); +}*/ + +} // namespace Triangulation_IO + +/////////////////////////////////////////////////////////////// +// TODO: replace these operator>> by an "input_point" function +/////////////////////////////////////////////////////////////// + +// TODO: test if the stream is binary or text? +template +std::istream & +operator>>(std::istream &is, typename Wrap::Point_d & p) +{ + typedef typename Wrap::Point_d P; + typedef typename K::FT FT; + std::vector coords; + + std::string line; + for(;;) + { + if (!std::getline(is, line)) + return is; + if (line != "") + break; + } + std::stringstream line_sstr(line); + FT temp; + while (line_sstr >> temp) + coords.push_back(temp); + + p = P(coords.begin(), coords.end()); + return is; +} + +// TODO: test if the stream is binary or text? +template +std::istream & +operator>>(std::istream &is, typename Wrap::Weighted_point_d & wp) +{ + typedef typename Wrap::Point_d P; + typedef typename Wrap::Weighted_point_d WP; + typedef typename K::FT FT; + + std::string line; + for(;;) + { + if (!std::getline(is, line)) + return is; + if (line != "") + break; + } + std::stringstream line_sstr(line); + FT temp; + std::vector coords; + while (line_sstr >> temp) + coords.push_back(temp); + + typename std::vector::iterator last = coords.end() - 1; + P p = P(coords.begin(), last); + wp = WP(p, *last); + + return is; +} + +// TODO: test if the stream is binary or text? +template +std::istream & +operator>>(std::istream &is, typename Wrap::Vector_d & v) +{ + typedef typename Wrap::Vector_d V; + typedef typename K::FT FT; + std::vector coords; + + std::string line; + for (;;) + { + if (!std::getline(is, line)) + return is; + if (line != "") + break; + } + std::stringstream line_sstr(line); + FT temp; + while (line_sstr >> temp) + coords.push_back(temp); + + v = V(coords.begin(), coords.end()); + return is; +} + +template < class GT, class TDS > +std::ostream & +export_triangulation_to_off(std::ostream & os, + const Triangulation & tr, + bool in_3D_export_surface_only = false) +{ + typedef Triangulation Tr; + typedef typename Tr::Vertex_const_handle Vertex_handle; + typedef typename Tr::Finite_vertex_const_iterator Finite_vertex_iterator; + typedef typename Tr::Finite_full_cell_const_iterator Finite_full_cell_iterator; + typedef typename Tr::Full_cell_const_iterator Full_cell_iterator; + typedef typename Tr::Full_cell Full_cell; + typedef typename Full_cell::Vertex_handle_const_iterator Full_cell_vertex_iterator; + + if (tr.maximal_dimension() < 2 || tr.maximal_dimension() > 3) + { + std::cerr << "Warning: export_tds_to_off => dimension should be 2 or 3."; + os << "Warning: export_tds_to_off => dimension should be 2 or 3."; + return os; + } + + size_t n = tr.number_of_vertices(); + + std::stringstream output; + + // write the vertices + std::map index_of_vertex; + int i = 0; + for(Finite_vertex_iterator it = tr.finite_vertices_begin(); + it != tr.finite_vertices_end(); ++it, ++i) + { + Triangulation_IO::output_point(output, tr.geom_traits(), it->point()); + if (tr.maximal_dimension() == 2) + output << " 0"; + output << std::endl; + index_of_vertex[it.base()] = i; + } + CGAL_assertion( i == n ); + + size_t number_of_triangles = 0; + if (tr.maximal_dimension() == 2) + { + for (Finite_full_cell_iterator fch = tr.finite_full_cells_begin() ; + fch != tr.finite_full_cells_end() ; ++fch) + { + output << "3 "; + for (Full_cell_vertex_iterator vit = fch->vertices_begin() ; + vit != fch->vertices_end() ; ++vit) + { + output << index_of_vertex[*vit] << " "; + } + output << std::endl; + ++number_of_triangles; + } + } + else if (tr.maximal_dimension() == 3) + { + if (in_3D_export_surface_only) + { + // Parse boundary facets + for (Full_cell_iterator fch = tr.full_cells_begin() ; + fch != tr.full_cells_end() ; ++fch) + { + if (tr.is_infinite(fch)) + { + output << "3 "; + for (Full_cell_vertex_iterator vit = fch->vertices_begin() ; + vit != fch->vertices_end() ; ++vit) + { + if (!tr.is_infinite(*vit)) + output << index_of_vertex[*vit] << " "; + } + output << std::endl; + ++number_of_triangles; + } + } + } + else + { + // Parse finite cells + for (Finite_full_cell_iterator fch = tr.finite_full_cells_begin() ; + fch != tr.finite_full_cells_end() ; ++fch) + { + output << "3 " + << index_of_vertex[fch->vertex(0)] << " " + << index_of_vertex[fch->vertex(1)] << " " + << index_of_vertex[fch->vertex(2)] + << std::endl; + output << "3 " + << index_of_vertex[fch->vertex(0)] << " " + << index_of_vertex[fch->vertex(2)] << " " + << index_of_vertex[fch->vertex(3)] + << std::endl; + output << "3 " + << index_of_vertex[fch->vertex(1)] << " " + << index_of_vertex[fch->vertex(2)] << " " + << index_of_vertex[fch->vertex(3)] + << std::endl; + output << "3 " + << index_of_vertex[fch->vertex(0)] << " " + << index_of_vertex[fch->vertex(1)] << " " + << index_of_vertex[fch->vertex(3)] + << std::endl; + number_of_triangles += 4; + } + } + } + + os << "OFF \n" + << n << " " + << number_of_triangles << " 0\n" + << output.str(); + + return os; +} + +} //namespace CGAL + +#endif // CGAL_TRIANGULATION_IO_H diff --git a/Triangulation/include/CGAL/Regular_triangulation.h b/Triangulation/include/CGAL/Regular_triangulation.h new file mode 100644 index 00000000000..111c6ac9441 --- /dev/null +++ b/Triangulation/include/CGAL/Regular_triangulation.h @@ -0,0 +1,1169 @@ +// Copyright (c) 2014 INRIA Sophia-Antipolis (France). +// All rights reserved. +// +// This file is part of CGAL (www.cgal.org). +// You can redistribute it and/or modify it under the terms of the GNU +// General Public License as published by the Free Software Foundation, +// either version 3 of the License, or (at your option) any later version. +// +// Licensees holding a valid commercial license may use this file in +// accordance with the commercial license agreement provided with the software. +// +// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE +// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE. +// +// $URL$ +// $Id$ +// +// Author(s) : Clement Jamin + +#ifndef CGAL_REGULAR_TRIANGULATION_H +#define CGAL_REGULAR_TRIANGULATION_H + +#include +#include +#include +#include +#include + +#include + +namespace CGAL { + +template< typename Traits_, typename TDS_ = Default > +class Regular_triangulation +: public Triangulation< + Regular_triangulation_traits_adapter, + typename Default::Get< + TDS_, + Triangulation_data_structure< + typename Regular_triangulation_traits_adapter::Dimension, + Triangulation_vertex >, + Triangulation_full_cell > + > + >::type> +{ + typedef Regular_triangulation_traits_adapter RTTraits; + typedef typename RTTraits::Dimension Maximal_dimension_; + typedef typename Default::Get< + TDS_, + Triangulation_data_structure< + Maximal_dimension_, + Triangulation_vertex, + Triangulation_full_cell + > >::type TDS; + typedef Triangulation Base; + typedef Regular_triangulation Self; + + typedef typename RTTraits::Orientation_d Orientation_d; + typedef typename RTTraits::Power_side_of_power_sphere_d Power_side_of_power_sphere_d; + typedef typename RTTraits::In_flat_power_side_of_power_sphere_d + In_flat_power_side_of_power_sphere_d; + typedef typename RTTraits::Flat_orientation_d Flat_orientation_d; + typedef typename RTTraits::Construct_flat_orientation_d Construct_flat_orientation_d; + +public: // PUBLIC NESTED TYPES + + typedef RTTraits Geom_traits; + typedef typename Base::Triangulation_ds Triangulation_ds; + + typedef typename Base::Vertex Vertex; + typedef typename Base::Full_cell Full_cell; + typedef typename Base::Facet Facet; + typedef typename Base::Face Face; + + typedef Maximal_dimension_ Maximal_dimension; + typedef typename RTTraits::Bare_point_d Bare_point; + typedef typename RTTraits::Weighted_point_d Weighted_point; + + typedef typename Base::Point_const_iterator Point_const_iterator; + typedef typename Base::Vertex_handle Vertex_handle; + typedef typename Base::Vertex_iterator Vertex_iterator; + typedef typename Base::Vertex_const_handle Vertex_const_handle; + typedef typename Base::Vertex_const_iterator Vertex_const_iterator; + + typedef typename Base::Full_cell_handle Full_cell_handle; + typedef typename Base::Full_cell_iterator Full_cell_iterator; + typedef typename Base::Full_cell_const_handle Full_cell_const_handle; + typedef typename Base::Full_cell_const_iterator Full_cell_const_iterator; + typedef typename Base::Finite_full_cell_const_iterator + Finite_full_cell_const_iterator; + + typedef typename Base::size_type size_type; + typedef typename Base::difference_type difference_type; + + typedef typename Base::Locate_type Locate_type; + + //Tag to distinguish Delaunay from Regular triangulations + typedef Tag_true Weighted_tag; + +protected: // DATA MEMBERS + + +public: + + using typename Base::Rotor; + using Base::maximal_dimension; + using Base::are_incident_full_cells_valid; + using Base::coaffine_orientation_predicate; + using Base::reset_flat_orientation; + using Base::current_dimension; + using Base::geom_traits; + using Base::index_of_covertex; + //using Base::index_of_second_covertex; + using Base::rotate_rotor; + using Base::infinite_vertex; + using Base::insert_in_hole; + using Base::is_infinite; + using Base::locate; + using Base::points_begin; + using Base::set_neighbors; + using Base::new_full_cell; + using Base::number_of_vertices; + using Base::orientation; + using Base::tds; + using Base::reorient_full_cells; + using Base::full_cell; + using Base::full_cells_begin; + using Base::full_cells_end; + using Base::finite_full_cells_begin; + using Base::finite_full_cells_end; + using Base::vertices_begin; + using Base::vertices_end; + +private: + + // Wrapper + struct Power_side_of_power_sphere_for_non_maximal_dim_d + { + boost::optional* fop; + Construct_flat_orientation_d cfo; + In_flat_power_side_of_power_sphere_d ifpt; + + Power_side_of_power_sphere_for_non_maximal_dim_d( + boost::optional& x, + Construct_flat_orientation_d const&y, + In_flat_power_side_of_power_sphere_d const&z) + : fop(&x), cfo(y), ifpt(z) {} + + template + CGAL::Orientation operator()(Iter a, Iter b, const Weighted_point & p)const + { + if(!*fop) + *fop=cfo(a,b); + return ifpt(fop->get(),a,b,p); + } + }; + +public: + +// - - - - - - - - - - - - - - - - - - - - - - - - - - CREATION / CONSTRUCTORS + + Regular_triangulation(int dim, const Geom_traits &k = Geom_traits()) + : Base(dim, k) + { + } + + // With this constructor, + // the user can specify a Flat_orientation_d object to be used for + // orienting simplices of a specific dimension + // (= preset_flat_orientation_.first) + // It it used by the dark triangulations created by DT::remove + Regular_triangulation( + int dim, + const std::pair &preset_flat_orientation, + const Geom_traits &k = Geom_traits()) + : Base(dim, preset_flat_orientation, k) + { + } + + ~Regular_triangulation() {} + +// - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - ACCESS + + // Not Documented + Power_side_of_power_sphere_for_non_maximal_dim_d power_side_of_power_sphere_for_non_maximal_dim_predicate() const + { + return Power_side_of_power_sphere_for_non_maximal_dim_d ( + flat_orientation_, + geom_traits().construct_flat_orientation_d_object(), + geom_traits().in_flat_power_side_of_power_sphere_d_object() + ); + } + + + // - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - REMOVALS + + // Warning: these functions are not correct since they do not restore hidden + // vertices + + Full_cell_handle remove(Vertex_handle); + Full_cell_handle remove(const Weighted_point & p, Full_cell_handle hint = Full_cell_handle()) + { + Locate_type lt; + Face f(maximal_dimension()); + Facet ft; + Full_cell_handle s = locate(p, lt, f, ft, hint); + if( Base::ON_VERTEX == lt ) + { + return remove(s->vertex(f.index(0))); + } + return Full_cell_handle(); + } + + template< typename ForwardIterator > + void remove(ForwardIterator start, ForwardIterator end) + { + while( start != end ) + remove(*start++); + } + + // Not documented + void remove_decrease_dimension(Vertex_handle); + + // - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - INSERTIONS + + template< typename ForwardIterator > + std::ptrdiff_t insert(ForwardIterator start, ForwardIterator end) + { + size_type n = number_of_vertices(); + typedef std::vector WP_vec; + WP_vec points(start, end); + + spatial_sort(points.begin(), points.end(), geom_traits()); + + Full_cell_handle hint; + for(typename WP_vec::const_iterator p = points.begin(); p != points.end(); ++p ) + { + Locate_type lt; + Face f(maximal_dimension()); + Facet ft; + Full_cell_handle c = locate (*p, lt, f, ft, hint); + Vertex_handle v = insert (*p, lt, f, ft, c); + + hint = v == Vertex_handle() ? c : v->full_cell(); + } + return number_of_vertices() - n; + } + + Vertex_handle insert(const Weighted_point &, + Locate_type, + const Face &, + const Facet &, + Full_cell_handle); + + Vertex_handle insert(const Weighted_point & p, + Full_cell_handle start = Full_cell_handle()) + { + Locate_type lt; + Face f(maximal_dimension()); + Facet ft; + Full_cell_handle s = locate(p, lt, f, ft, start); + return insert(p, lt, f, ft, s); + } + + Vertex_handle insert(const Weighted_point & p, Vertex_handle hint) + { + CGAL_assertion( Vertex_handle() != hint ); + return insert(p, hint->full_cell()); + } + + Vertex_handle insert_outside_affine_hull(const Weighted_point &); + Vertex_handle insert_in_conflicting_cell( + const Weighted_point &, Full_cell_handle, + Vertex_handle only_if_this_vertex_is_in_the_cz = Vertex_handle()); + + Vertex_handle insert_if_in_star(const Weighted_point &, + Vertex_handle, + Locate_type, + const Face &, + const Facet &, + Full_cell_handle); + + Vertex_handle insert_if_in_star( + const Weighted_point & p, Vertex_handle star_center, + Full_cell_handle start = Full_cell_handle()) + { + Locate_type lt; + Face f(maximal_dimension()); + Facet ft; + Full_cell_handle s = locate(p, lt, f, ft, start); + return insert_if_in_star(p, star_center, lt, f, ft, s); + } + + Vertex_handle insert_if_in_star( + const Weighted_point & p, Vertex_handle star_center, + Vertex_handle hint) + { + CGAL_assertion( Vertex_handle() != hint ); + return insert_if_in_star(p, star_center, hint->full_cell()); + } + +// - - - - - - - - - - - - - - - - - - - - - - - - - GATHERING CONFLICTING SIMPLICES + + bool is_in_conflict(const Weighted_point &, Full_cell_const_handle) const; + + template< class OrientationPredicate > + Oriented_side perturbed_power_side_of_power_sphere(const Weighted_point &, + Full_cell_const_handle, const OrientationPredicate &) const; + + template< typename OutputIterator > + Facet compute_conflict_zone(const Weighted_point &, Full_cell_handle, OutputIterator) const; + + template < typename OrientationPredicate, typename PowerTestPredicate > + class Conflict_predicate + { + const Self & rt_; + const Weighted_point & p_; + OrientationPredicate ori_; + PowerTestPredicate power_side_of_power_sphere_; + int cur_dim_; + public: + Conflict_predicate( + const Self & rt, + const Weighted_point & p, + const OrientationPredicate & ori, + const PowerTestPredicate & power_side_of_power_sphere) + : rt_(rt), p_(p), ori_(ori), power_side_of_power_sphere_(power_side_of_power_sphere), cur_dim_(rt.current_dimension()) {} + + inline + bool operator()(Full_cell_const_handle s) const + { + bool ok; + if( ! rt_.is_infinite(s) ) + { + Oriented_side power_side_of_power_sphere = power_side_of_power_sphere_(rt_.points_begin(s), rt_.points_begin(s) + cur_dim_ + 1, p_); + if( ON_POSITIVE_SIDE == power_side_of_power_sphere ) + ok = true; + else if( ON_NEGATIVE_SIDE == power_side_of_power_sphere ) + ok = false; + else + ok = ON_POSITIVE_SIDE == rt_.perturbed_power_side_of_power_sphere(p_, s, ori_); + } + else + { + typedef typename Full_cell::Vertex_handle_const_iterator VHCI; + typedef Substitute_point_in_vertex_iterator F; + F spivi(rt_.infinite_vertex(), &p_); + + Orientation o = ori_( + boost::make_transform_iterator(s->vertices_begin(), spivi), + boost::make_transform_iterator(s->vertices_begin() + cur_dim_ + 1, + spivi)); + + if( POSITIVE == o ) + ok = true; + else if( o == NEGATIVE ) + ok = false; + else + ok = (*this)(s->neighbor( s->index( rt_.infinite_vertex() ) )); + } + return ok; + } + }; + + template < typename ConflictPredicate > + class Conflict_traversal_predicate + { + const Self & rt_; + const ConflictPredicate & pred_; + public: + Conflict_traversal_predicate(const Self & rt, const ConflictPredicate & pred) + : rt_(rt), pred_(pred) + {} + inline + bool operator()(const Facet & f) const + { + return pred_(rt_.full_cell(f)->neighbor(rt_.index_of_covertex(f))); + } + }; + +// - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - VALIDITY + + bool is_valid(bool verbose = false, int level = 0) const; + +// - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - MISC + + std::size_t number_of_hidden_vertices() const + { + return m_hidden_points.size(); + } + +private: + + template + bool + does_cell_range_contain_vertex(InputIterator cz_begin, InputIterator cz_end, + Vertex_handle vh) const + { + // Check all vertices + while(cz_begin != cz_end) + { + Full_cell_handle fch = *cz_begin; + for (int i = 0 ; i <= current_dimension() ; ++i) + { + if (fch->vertex(i) == vh) + return true; + } + ++cz_begin; + } + return false; + } + + template + void + process_conflict_zone(InputIterator cz_begin, InputIterator cz_end, + OutputIterator vertices_out) const + { + // Get all vertices + while(cz_begin != cz_end) + { + Full_cell_handle fch = *cz_begin; + for (int i = 0 ; i <= current_dimension() ; ++i) + { + Vertex_handle vh = fch->vertex(i); + if (vh->full_cell() != Full_cell_handle()) + { + (*vertices_out++) = vh; + vh->set_full_cell(Full_cell_handle()); + } + } + ++cz_begin; + } + } + + + template + void + process_cz_vertices_after_insertion(InputIterator vertices_begin, + InputIterator vertices_end) + { + // Get all vertices + while(vertices_begin != vertices_end) + { + Vertex_handle vh = *vertices_begin; + if (vh->full_cell() == Full_cell_handle()) + { + m_hidden_points.push_back(vh->point()); + tds().delete_vertex(vh); + } + ++vertices_begin; + } + } + +private: + // Some internal types to shorten notation + using typename Base::Coaffine_orientation_d; + using Base::flat_orientation_; + typedef Conflict_predicate + Conflict_pred_in_subspace; + typedef Conflict_predicate + Conflict_pred_in_fullspace; + typedef Conflict_traversal_predicate + Conflict_traversal_pred_in_subspace; + typedef Conflict_traversal_predicate + Conflict_traversal_pred_in_fullspace; + +// - - - - - - - - - - - - - - - - - - - - - - - - - - - - - MEMBER VARIABLES + std::vector m_hidden_points; + +}; // class Regular_triangulation + + +// = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = +// FUNCTIONS THAT ARE MEMBER METHODS: + +// - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - REMOVALS + + +// Warning: this function is not correct since it does not restore hidden +// vertices +template< typename Traits, typename TDS > +typename Regular_triangulation::Full_cell_handle +Regular_triangulation +::remove( Vertex_handle v ) +{ + CGAL_precondition( ! is_infinite(v) ); + CGAL_expensive_precondition( is_vertex(v) ); + + // THE CASE cur_dim == 0 + if( 0 == current_dimension() ) + { + remove_decrease_dimension(v); + return Full_cell_handle(); + } + else if( 1 == current_dimension() ) + { // THE CASE cur_dim == 1 + if( 2 == number_of_vertices() ) + { + remove_decrease_dimension(v); + return Full_cell_handle(); + } + Full_cell_handle left = v->full_cell(); + if( 0 == left->index(v) ) + left = left->neighbor(1); + CGAL_assertion( 1 == left->index(v) ); + Full_cell_handle right = left->neighbor(0); + tds().associate_vertex_with_full_cell(left, 1, right->vertex(1)); + set_neighbors(left, 0, right->neighbor(0), right->mirror_index(0)); + tds().delete_vertex(v); + tds().delete_full_cell(right); + return left; + } + + // THE CASE cur_dim >= 2 + // Gather the finite vertices sharing an edge with |v| + typedef typename Base::template Full_cell_set Simplices; + Simplices simps; + std::back_insert_iterator out(simps); + tds().incident_full_cells(v, out); + typedef std::set Vertex_set; + Vertex_set verts; + Vertex_handle vh; + for( typename Simplices::iterator it = simps.begin(); it != simps.end(); ++it ) + for( int i = 0; i <= current_dimension(); ++i ) + { + vh = (*it)->vertex(i); + if( is_infinite(vh) ) + continue; + if( vh == v ) + continue; + verts.insert(vh); + } + + // After gathering finite neighboring vertices, create their Dark Delaunay triangulation + typedef Triangulation_vertex Dark_vertex_base; + typedef Triangulation_full_cell< + Geom_traits, + internal::Triangulation::Dark_full_cell_data > Dark_full_cell_base; + typedef Triangulation_data_structure Dark_tds; + typedef Regular_triangulation Dark_triangulation; + typedef typename Dark_triangulation::Face Dark_face; + typedef typename Dark_triangulation::Facet Dark_facet; + typedef typename Dark_triangulation::Vertex_handle Dark_v_handle; + typedef typename Dark_triangulation::Full_cell_handle Dark_s_handle; + + // If flat_orientation_ is defined, we give it the Dark triangulation + // so that the orientation it uses for "current_dimension()"-simplices is + // coherent with the global triangulation + Dark_triangulation dark_side( + maximal_dimension(), + flat_orientation_ ? + std::pair(current_dimension(), flat_orientation_.get_ptr()) + : std::pair(std::numeric_limits::max(), NULL) ); + + Dark_s_handle dark_s; + Dark_v_handle dark_v; + typedef std::map Vertex_map; + Vertex_map light_to_dark; + typename Vertex_set::iterator vit = verts.begin(); + while( vit != verts.end() ) + { + dark_v = dark_side.insert((*vit)->point(), dark_s); + dark_s = dark_v->full_cell(); + dark_v->data() = *vit; + light_to_dark[*vit] = dark_v; + ++vit; + } + + if( dark_side.current_dimension() != current_dimension() ) + { + CGAL_assertion( dark_side.current_dimension() + 1 == current_dimension() ); + // Here, the finite neighbors of |v| span a affine subspace of + // dimension one less than the current dimension. Two cases are possible: + if( (size_type)(verts.size() + 1) == number_of_vertices() ) + { + remove_decrease_dimension(v); + return Full_cell_handle(); + } + else + { // |v| is strictly outside the convex hull of the rest of the points. This is an + // easy case: first, modify the finite full_cells, then, delete the infinite ones. + // We don't even need the Dark triangulation. + Simplices infinite_simps; + { + Simplices finite_simps; + for( typename Simplices::iterator it = simps.begin(); it != simps.end(); ++it ) + if( is_infinite(*it) ) + infinite_simps.push_back(*it); + else + finite_simps.push_back(*it); + simps.swap(finite_simps); + } // now, simps only contains finite simplices + // First, modify the finite full_cells: + for( typename Simplices::iterator it = simps.begin(); it != simps.end(); ++it ) + { + int v_idx = (*it)->index(v); + tds().associate_vertex_with_full_cell(*it, v_idx, infinite_vertex()); + } + // Make the handles to infinite full cells searchable + infinite_simps.make_searchable(); + // Then, modify the neighboring relation + for( typename Simplices::iterator it = simps.begin(); it != simps.end(); ++it ) + { + for( int i = 0 ; i <= current_dimension(); ++i ) + { + if (is_infinite((*it)->vertex(i))) + continue; + (*it)->vertex(i)->set_full_cell(*it); + Full_cell_handle n = (*it)->neighbor(i); + // Was |n| a finite full cell prior to removing |v| ? + if( ! infinite_simps.contains(n) ) + continue; + int n_idx = n->index(v); + set_neighbors(*it, i, n->neighbor(n_idx), n->neighbor(n_idx)->index(n)); + } + } + Full_cell_handle ret_s; + // Then, we delete the infinite full_cells + for( typename Simplices::iterator it = infinite_simps.begin(); it != infinite_simps.end(); ++it ) + tds().delete_full_cell(*it); + tds().delete_vertex(v); + return simps.front(); + } + } + else // From here on, dark_side.current_dimension() == current_dimension() + { + dark_side.infinite_vertex()->data() = infinite_vertex(); + light_to_dark[infinite_vertex()] = dark_side.infinite_vertex(); + } + + // Now, compute the conflict zone of v->point() in + // the dark side. This is precisely the set of full_cells + // that we have to glue back into the light side. + Dark_face dark_f(dark_side.maximal_dimension()); + Dark_facet dark_ft; + typename Dark_triangulation::Locate_type lt; + dark_s = dark_side.locate(v->point(), lt, dark_f, dark_ft); + CGAL_assertion( lt != Dark_triangulation::ON_VERTEX + && lt != Dark_triangulation::OUTSIDE_AFFINE_HULL ); + + // |ret_s| is the full_cell that we return + Dark_s_handle dark_ret_s = dark_s; + Full_cell_handle ret_s; + + typedef typename Base::template Full_cell_set Dark_full_cells; + Dark_full_cells conflict_zone; + std::back_insert_iterator dark_out(conflict_zone); + + dark_ft = dark_side.compute_conflict_zone(v->point(), dark_s, dark_out); + // Make the dark simplices in the conflict zone searchable + conflict_zone.make_searchable(); + + // THE FOLLOWING SHOULD MAYBE GO IN TDS. + // Here is the plan: + // 1. Pick any Facet from boundary of the light zone + // 2. Find corresponding Facet on boundary of dark zone + // 3. stitch. + + // 1. Build a facet on the boudary of the light zone: + Full_cell_handle light_s = *simps.begin(); + Facet light_ft(light_s, light_s->index(v)); + + // 2. Find corresponding Dark_facet on boundary of the dark zone + Dark_full_cells dark_incident_s; + for( int i = 0; i <= current_dimension(); ++i ) + { + if( index_of_covertex(light_ft) == i ) + continue; + Dark_v_handle dark_v = light_to_dark[full_cell(light_ft)->vertex(i)]; + dark_incident_s.clear(); + dark_out = std::back_inserter(dark_incident_s); + dark_side.tds().incident_full_cells(dark_v, dark_out); + for(typename Dark_full_cells::iterator it = dark_incident_s.begin(); + it != dark_incident_s.end(); + ++it) + { + (*it)->data().count_ += 1; + } + } + + for( typename Dark_full_cells::iterator it = dark_incident_s.begin(); it != dark_incident_s.end(); ++it ) + { + if( current_dimension() != (*it)->data().count_ ) + continue; + if( ! conflict_zone.contains(*it) ) + continue; + // We found a full_cell incident to the dark facet corresponding to the light facet |light_ft| + int ft_idx = 0; + while( light_s->has_vertex( (*it)->vertex(ft_idx)->data() ) ) + ++ft_idx; + dark_ft = Dark_facet(*it, ft_idx); + break; + } + // Pre-3. Now, we are ready to traverse both boundary and do the stiching. + + // But first, we create the new full_cells in the light triangulation, + // with as much adjacency information as possible. + + // Create new full_cells with vertices + for( typename Dark_full_cells::iterator it = conflict_zone.begin(); it != conflict_zone.end(); ++it ) + { + Full_cell_handle new_s = new_full_cell(); + (*it)->data().light_copy_ = new_s; + for( int i = 0; i <= current_dimension(); ++i ) + tds().associate_vertex_with_full_cell(new_s, i, (*it)->vertex(i)->data()); + if( dark_ret_s == *it ) + ret_s = new_s; + } + + // Setup adjacencies inside the hole + for( typename Dark_full_cells::iterator it = conflict_zone.begin(); it != conflict_zone.end(); ++it ) + { + Full_cell_handle new_s = (*it)->data().light_copy_; + for( int i = 0; i <= current_dimension(); ++i ) + if( conflict_zone.contains((*it)->neighbor(i)) ) + tds().set_neighbors(new_s, i, (*it)->neighbor(i)->data().light_copy_, (*it)->mirror_index(i)); + } + + // 3. Stitch + simps.make_searchable(); + typedef std::queue > Queue; + Queue q; + q.push(std::make_pair(light_ft, dark_ft)); + dark_s = dark_side.full_cell(dark_ft); + int dark_i = dark_side.index_of_covertex(dark_ft); + // mark dark_ft as visited: + // TODO try by marking with Dark_v_handle (vertex) + dark_s->neighbor(dark_i)->set_neighbor(dark_s->mirror_index(dark_i), Dark_s_handle()); + while( ! q.empty() ) + { + std::pair p = q.front(); + q.pop(); + light_ft = p.first; + dark_ft = p.second; + light_s = full_cell(light_ft); + int light_i = index_of_covertex(light_ft); + dark_s = dark_side.full_cell(dark_ft); + int dark_i = dark_side.index_of_covertex(dark_ft); + Full_cell_handle light_n = light_s->neighbor(light_i); + set_neighbors(dark_s->data().light_copy_, dark_i, light_n, light_s->mirror_index(light_i)); + for( int di = 0; di <= current_dimension(); ++di ) + { + if( di == dark_i ) + continue; + int li = light_s->index(dark_s->vertex(di)->data()); + Rotor light_r(light_s, li, light_i); + typename Dark_triangulation::Rotor dark_r(dark_s, di, dark_i); + + while( simps.contains(cpp11::get<0>(light_r)->neighbor(cpp11::get<1>(light_r))) ) + light_r = rotate_rotor(light_r); + + while( conflict_zone.contains(cpp11::get<0>(dark_r)->neighbor(cpp11::get<1>(dark_r))) ) + dark_r = dark_side.rotate_rotor(dark_r); + + Dark_s_handle dark_ns = cpp11::get<0>(dark_r); + int dark_ni = cpp11::get<1>(dark_r); + Full_cell_handle light_ns = cpp11::get<0>(light_r); + int light_ni = cpp11::get<1>(light_r); + // mark dark_r as visited: + // TODO try by marking with Dark_v_handle (vertex) + Dark_s_handle outside = dark_ns->neighbor(dark_ni); + Dark_v_handle mirror = dark_ns->mirror_vertex(dark_ni, current_dimension()); + int dn = outside->index(mirror); + if( Dark_s_handle() == outside->neighbor(dn) ) + continue; + outside->set_neighbor(dn, Dark_s_handle()); + q.push(std::make_pair(Facet(light_ns, light_ni), Dark_facet(dark_ns, dark_ni))); + } + } + tds().delete_full_cells(simps.begin(), simps.end()); + tds().delete_vertex(v); + return ret_s; +} + +template< typename Traits, typename TDS > +void +Regular_triangulation +::remove_decrease_dimension(Vertex_handle v) +{ + CGAL_precondition( current_dimension() >= 0 ); + tds().remove_decrease_dimension(v, infinite_vertex()); + // reset the predicates: + reset_flat_orientation(); + if( 1 <= current_dimension() ) + { + Full_cell_handle inf_v_cell = infinite_vertex()->full_cell(); + int inf_v_index = inf_v_cell->index(infinite_vertex()); + Full_cell_handle s = inf_v_cell->neighbor(inf_v_index); + Orientation o = orientation(s); + CGAL_assertion( ZERO != o ); + if( NEGATIVE == o ) + reorient_full_cells(); + } +} + +// - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - INSERTIONS + +template< typename Traits, typename TDS > +typename Regular_triangulation::Vertex_handle +Regular_triangulation +::insert(const Weighted_point & p, Locate_type lt, const Face & f, const Facet & ft, Full_cell_handle s) +{ + switch( lt ) + { + case Base::OUTSIDE_AFFINE_HULL: + return insert_outside_affine_hull(p); + break; + case Base::ON_VERTEX: + { + Vertex_handle v = s->vertex(f.index(0)); + typename RTTraits::Compute_weight_d pw = + geom_traits().compute_weight_d_object(); + + if (pw(p) == pw(v->point())) + return v; + // If dim == 0 and the new point has a bigger weight, + // we just replace the point, and the former point gets hidden + else if (current_dimension() == 0) + { + if (pw(p) > pw(v->point())) + { + m_hidden_points.push_back(v->point()); + v->set_point(p); + return v; + } + // Otherwise, the new point is hidden + else + { + m_hidden_points.push_back(p); + return Vertex_handle(); + } + } + // Otherwise, we apply the "normal" algorithm + + // !NO break here! + } + default: + return insert_in_conflicting_cell(p, s); + } +} + +/* +Inserts the point `p` in the regular triangulation. Returns a handle to the +newly created vertex at that position. +\pre The point `p` +must lie outside the affine hull of the regular triangulation. This implies that +`rt`.`current_dimension()` must be smaller than `rt`.`maximal_dimension()`. +*/ +template< typename Traits, typename TDS > +typename Regular_triangulation::Vertex_handle +Regular_triangulation +::insert_outside_affine_hull(const Weighted_point & p) +{ + // we don't use Base::insert_outside_affine_hull(...) because here, we + // also need to reset the side_of_oriented_subsphere functor. + CGAL_precondition( current_dimension() < maximal_dimension() ); + Vertex_handle v = tds().insert_increase_dimension(infinite_vertex()); + // reset the predicates: + reset_flat_orientation(); + v->set_point(p); + if( current_dimension() >= 1 ) + { + Full_cell_handle inf_v_cell = infinite_vertex()->full_cell(); + int inf_v_index = inf_v_cell->index(infinite_vertex()); + Full_cell_handle s = inf_v_cell->neighbor(inf_v_index); + Orientation o = orientation(s); + CGAL_assertion( ZERO != o ); + if( NEGATIVE == o ) + reorient_full_cells(); + + // We just inserted the second finite point and the right infinite + // cell is like : (inf_v, v), but we want it to be (v, inf_v) to be + // consistent with the rest of the cells + if (current_dimension() == 1) + { + // Is "inf_v_cell" the right infinite cell? Then inf_v_index should be 1 + if (inf_v_cell->neighbor(inf_v_index)->index(inf_v_cell) == 0 + && inf_v_index == 0) + { + inf_v_cell->swap_vertices(current_dimension() - 1, current_dimension()); + } + else + { + inf_v_cell = inf_v_cell->neighbor((inf_v_index + 1) % 2); + inf_v_index = inf_v_cell->index(infinite_vertex()); + // Is "inf_v_cell" the right infinite cell? Then inf_v_index should be 1 + if (inf_v_cell->neighbor(inf_v_index)->index(inf_v_cell) == 0 + && inf_v_index == 0) + { + inf_v_cell->swap_vertices(current_dimension() - 1, current_dimension()); + } + } + } + } + return v; +} + +template< typename Traits, typename TDS > +typename Regular_triangulation::Vertex_handle +Regular_triangulation +::insert_if_in_star(const Weighted_point & p, + Vertex_handle star_center, + Locate_type lt, + const Face & f, + const Facet & ft, + Full_cell_handle s) +{ + switch( lt ) + { + case Base::OUTSIDE_AFFINE_HULL: + return insert_outside_affine_hull(p); + break; + case Base::ON_VERTEX: + { + Vertex_handle v = s->vertex(f.index(0)); + typename RTTraits::Compute_weight_d pw = + geom_traits().compute_weight_d_object(); + if (pw(p) == pw(v->point())) + return v; + // If dim == 0 and the new point has a bigger weight, + // we replace the point + else if (current_dimension() == 0) + { + if (pw(p) > pw(v->point())) + v->set_point(p); + else + return v; + } + // Otherwise, we apply the "normal" algorithm + + // !NO break here! + } + default: + return insert_in_conflicting_cell(p, s, star_center); + } + + return Vertex_handle(); +} + +/* +[Undocumented function] + +Inserts the point `p` in the regular triangulation. `p` must be +in conflict with the second parameter `c`, which is used as a +starting point for `compute_conflict_zone`. +The function is faster than the standard `insert` function since +it does not need to call `locate`. + +If this insertion creates a vertex, this vertex is returned. + +If `p` coincides with an existing vertex and has a greater weight, +then the existing weighted point becomes hidden and `p` replaces it as vertex +of the triangulation. + +If `p` coincides with an already existing vertex (both point and +weights being equal), then this vertex is returned and the triangulation +remains unchanged. + +Otherwise if `p` does not appear as a vertex of the triangulation, +then it is stored as a hidden point and this method returns the default +constructed handle. + +\pre The point `p` must be in conflict with the full cell `c`. +*/ + +template< typename Traits, typename TDS > +typename Regular_triangulation::Vertex_handle +Regular_triangulation +::insert_in_conflicting_cell(const Weighted_point & p, + Full_cell_handle s, + Vertex_handle only_if_this_vertex_is_in_the_cz) +{ + typedef std::vector Full_cell_h_vector; + + bool in_conflict = is_in_conflict(p, s); + + // If p is not in conflict with s, then p is hidden + // => we don't insert it + if (!in_conflict) + { + m_hidden_points.push_back(p); + return Vertex_handle(); + } + else + { + Full_cell_h_vector cs; // for storing conflicting full_cells. + cs.reserve(64); + std::back_insert_iterator out(cs); + Facet ft = compute_conflict_zone(p, s, out); + + // Check if the CZ contains "only_if_this_vertex_is_in_the_cz" + if (only_if_this_vertex_is_in_the_cz != Vertex_handle() + && !does_cell_range_contain_vertex(cs.begin(), cs.end(), + only_if_this_vertex_is_in_the_cz)) + { + return Vertex_handle(); + } + + // Otherwise, proceed with the insertion + std::vector cz_vertices; + cz_vertices.reserve(64); + process_conflict_zone(cs.begin(), cs.end(), + std::back_inserter(cz_vertices)); + + Vertex_handle ret = insert_in_hole(p, cs.begin(), cs.end(), ft); + + process_cz_vertices_after_insertion(cz_vertices.begin(), cz_vertices.end()); + + return ret; + } +} + +// - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - GATHERING CONFLICTING SIMPLICES + +// NOT DOCUMENTED +template< typename Traits, typename TDS > +template< typename OrientationPred > +Oriented_side +Regular_triangulation +::perturbed_power_side_of_power_sphere(const Weighted_point & p, Full_cell_const_handle s, + const OrientationPred & ori) const +{ + CGAL_precondition_msg( ! is_infinite(s), "full cell must be finite"); + CGAL_expensive_precondition( POSITIVE == orientation(s) ); + typedef std::vector Points; + Points points(current_dimension() + 2); + int i(0); + for( ; i <= current_dimension(); ++i ) + points[i] = &(s->vertex(i)->point()); + points[i] = &p; + std::sort(points.begin(), points.end(), + internal::Triangulation::Compare_points_for_perturbation(*this)); + typename Points::const_reverse_iterator cut_pt = points.rbegin(); + Points test_points; + while( cut_pt != points.rend() ) + { + if( &p == *cut_pt ) + // because the full_cell "s" is assumed to be positively oriented + return ON_NEGATIVE_SIDE; // we consider |p| to lie outside the sphere + test_points.clear(); + Point_const_iterator spit = points_begin(s); + int adjust_sign = -1; + for( i = 0; i < current_dimension(); ++i ) + { + if( &(*spit) == *cut_pt ) + { + ++spit; + adjust_sign = (((current_dimension() + i) % 2) == 0) ? -1 : +1; + } + test_points.push_back(&(*spit)); + ++spit; + } + test_points.push_back(&p); + + typedef typename CGAL::Iterator_project< + typename Points::iterator, + internal::Triangulation::Point_from_pointer, + const Weighted_point &, const Weighted_point * + > Point_pointer_iterator; + + Orientation ori_value = ori( + Point_pointer_iterator(test_points.begin()), + Point_pointer_iterator(test_points.end())); + + if( ZERO != ori_value ) + return Oriented_side( - adjust_sign * ori_value ); + + ++cut_pt; + } + CGAL_assertion(false); // we should never reach here + return ON_NEGATIVE_SIDE; +} + +template< typename Traits, typename TDS > +bool +Regular_triangulation +::is_in_conflict(const Weighted_point & p, Full_cell_const_handle s) const +{ + CGAL_precondition( 1 <= current_dimension() ); + if( current_dimension() < maximal_dimension() ) + { + Conflict_pred_in_subspace c( + *this, p, + coaffine_orientation_predicate(), + power_side_of_power_sphere_for_non_maximal_dim_predicate()); + return c(s); + } + else + { + Orientation_d ori = geom_traits().orientation_d_object(); + Power_side_of_power_sphere_d side = geom_traits().power_side_of_power_sphere_d_object(); + Conflict_pred_in_fullspace c(*this, p, ori, side); + return c(s); + } +} + +template< typename Traits, typename TDS > +template< typename OutputIterator > +typename Regular_triangulation::Facet +Regular_triangulation +::compute_conflict_zone(const Weighted_point & p, Full_cell_handle s, OutputIterator out) const +{ + CGAL_precondition( 1 <= current_dimension() ); + if( current_dimension() < maximal_dimension() ) + { + Conflict_pred_in_subspace c( + *this, p, + coaffine_orientation_predicate(), + power_side_of_power_sphere_for_non_maximal_dim_predicate()); + Conflict_traversal_pred_in_subspace tp(*this, c); + return tds().gather_full_cells(s, tp, out); + } + else + { + Orientation_d ori = geom_traits().orientation_d_object(); + Power_side_of_power_sphere_d side = geom_traits().power_side_of_power_sphere_d_object(); + Conflict_pred_in_fullspace c(*this, p, ori, side); + Conflict_traversal_pred_in_fullspace tp(*this, c); + return tds().gather_full_cells(s, tp, out); + } +} + +// - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - VALIDITY + +template< typename Traits, typename TDS > +bool +Regular_triangulation +::is_valid(bool verbose, int level) const +{ + if (!Base::is_valid(verbose, level)) + return false; + + int dim = current_dimension(); + if (dim == maximal_dimension()) + { + for (Finite_full_cell_const_iterator cit = finite_full_cells_begin() ; + cit != finite_full_cells_end() ; ++cit ) + { + Full_cell_const_handle ch = cit.base(); + for(int i = 0; i < dim+1 ; ++i ) + { + // If the i-th neighbor is not an infinite cell + Vertex_handle opposite_vh = + ch->neighbor(i)->vertex(ch->neighbor(i)->index(ch)); + if (!is_infinite(opposite_vh)) + { + Power_side_of_power_sphere_d side = + geom_traits().power_side_of_power_sphere_d_object(); + if (side(Point_const_iterator(ch->vertices_begin()), + Point_const_iterator(ch->vertices_end()), + opposite_vh->point()) == ON_POSITIVE_SIDE) + { + if (verbose) + CGAL_warning_msg(false, "Non-empty sphere"); + return false; + } + } + } + } + } + return true; +} + +} //namespace CGAL + +#endif //CGAL_REGULAR_TRIANGULATION_H diff --git a/Triangulation/include/CGAL/Regular_triangulation_traits_adapter.h b/Triangulation/include/CGAL/Regular_triangulation_traits_adapter.h new file mode 100644 index 00000000000..78bb95a60aa --- /dev/null +++ b/Triangulation/include/CGAL/Regular_triangulation_traits_adapter.h @@ -0,0 +1,288 @@ +// Copyright (c) 2014 INRIA Sophia-Antipolis (France). +// All rights reserved. +// +// This file is part of CGAL (www.cgal.org). +// You can redistribute it and/or modify it under the terms of the GNU +// General Public License as published by the Free Software Foundation, +// either version 3 of the License, or (at your option) any later version. +// +// Licensees holding a valid commercial license may use this file in +// accordance with the commercial license agreement provided with the software. +// +// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE +// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE. +// +// $URL$ +// $Id$ +// +// Author(s) : Clement Jamin + +#ifndef CGAL_REGULAR_TRIANGULATION_TRAITS_ADAPTER_H +#define CGAL_REGULAR_TRIANGULATION_TRAITS_ADAPTER_H + +#include + +#include + +namespace CGAL { + +// Wrapper class to make a model of `RegularTriangulationTraits` easily usable +// by the `Regular_triangulation` class. By using this class: +// - Point_d (used by `Triangulation` and the TDS) becomes a weighted point +// - Predicates and functors such as Less_coordinate_d or Orientation_d +// can be called using weighted points instead of bare points (this is +// needed because `Weighted_point_d` is not convertible to `Point_d`) +// This way, `Triangulation` works perfectly well with weighted points. + +template +class Regular_triangulation_traits_adapter + : public K +{ +public: + typedef K Base; + + // Required by TriangulationTraits + typedef typename K::Dimension Dimension; + typedef typename K::FT FT; + typedef typename K::Flat_orientation_d Flat_orientation_d; + typedef typename K::Weighted_point_d Point_d; + + // Required by RegularTriangulationTraits + typedef typename K::Point_d Bare_point_d; + typedef typename K::Weighted_point_d Weighted_point_d; + typedef typename K::Construct_point_d Construct_point_d; + typedef typename K::Compute_weight_d Compute_weight_d; + typedef typename K::Power_side_of_power_sphere_d Power_side_of_power_sphere_d; + typedef typename K::In_flat_power_side_of_power_sphere_d + In_flat_power_side_of_power_sphere_d; + + //=========================================================================== + // Custom types + //=========================================================================== + + // Required by SpatialSortingTraits_d + class Less_coordinate_d + { + const K &m_kernel; + + public: + typedef bool result_type; + + Less_coordinate_d(const K &kernel) + : m_kernel(kernel) {} + + result_type operator()( + Weighted_point_d const& p, Weighted_point_d const& q, int i) const + { + Construct_point_d cp = m_kernel.construct_point_d_object(); + return m_kernel.less_coordinate_d_object() (cp(p), cp(q), i); + } + }; + + //=========================================================================== + + // Required by TriangulationTraits + class Orientation_d + { + const K &m_kernel; + + public: + typedef Orientation result_type; + + Orientation_d(const K &kernel) + : m_kernel(kernel) {} + + template + result_type operator()(ForwardIterator start, ForwardIterator end) const + { + Construct_point_d cp = m_kernel.construct_point_d_object(); + return m_kernel.orientation_d_object() ( + boost::make_transform_iterator(start, cp), + boost::make_transform_iterator(end, cp) + ); + } + }; + + //=========================================================================== + + // Required by TriangulationTraits + class Construct_flat_orientation_d + { + const K &m_kernel; + + public: + typedef Flat_orientation_d result_type; + + Construct_flat_orientation_d(const K &kernel) + : m_kernel(kernel) {} + + template + result_type operator()(ForwardIterator start, ForwardIterator end) const + { + Construct_point_d cp = m_kernel.construct_point_d_object(); + return m_kernel.construct_flat_orientation_d_object() ( + boost::make_transform_iterator(start, cp), + boost::make_transform_iterator(end, cp) + ); + } + }; + + + //=========================================================================== + + // Required by TriangulationTraits + class In_flat_orientation_d + { + const K &m_kernel; + + public: + typedef Orientation result_type; + + In_flat_orientation_d(const K &kernel) + : m_kernel(kernel) {} + + template + result_type operator()(Flat_orientation_d orient, + ForwardIterator start, ForwardIterator end) const + { + Construct_point_d cp = m_kernel.construct_point_d_object(); + return m_kernel.in_flat_orientation_d_object() ( + orient, + boost::make_transform_iterator(start, cp), + boost::make_transform_iterator(end, cp) + ); + } + }; + + //=========================================================================== + + // Required by TriangulationTraits + class Contained_in_affine_hull_d + { + const K &m_kernel; + + public: + typedef bool result_type; + + Contained_in_affine_hull_d(const K &kernel) + : m_kernel(kernel) {} + + template + result_type operator()(ForwardIterator start, ForwardIterator end, + const Weighted_point_d & p) const + { + Construct_point_d cp = m_kernel.construct_point_d_object(); + return m_kernel.contained_in_affine_hull_d_object() ( + boost::make_transform_iterator(start, cp), + boost::make_transform_iterator(end, cp), + cp(p) + ); + } + }; + + //=========================================================================== + + // Required by TriangulationTraits + class Compare_lexicographically_d + { + const K &m_kernel; + + public: + typedef Comparison_result result_type; + + Compare_lexicographically_d(const K &kernel) + : m_kernel(kernel) {} + + result_type operator()( + const Weighted_point_d & p, const Weighted_point_d & q) const + { + Construct_point_d cp = m_kernel.construct_point_d_object(); + return m_kernel.compare_lexicographically_d_object()(cp(p), cp(q)); + } + }; + + //=========================================================================== + + // Only for Triangulation_off_ostream.h (undocumented) + class Compute_coordinate_d + { + const K &m_kernel; + + public: + typedef FT result_type; + + Compute_coordinate_d(const K &kernel) + : m_kernel(kernel) {} + + result_type operator()( + const Weighted_point_d & p, const int i) const + { + Construct_point_d cp = m_kernel.construct_point_d_object(); + return m_kernel.compute_coordinate_d_object()(cp(p), i); + } + }; + + //=========================================================================== + + // To satisfy SpatialSortingTraits_d + // and also for Triangulation_off_ostream.h (undocumented) + class Point_dimension_d + { + const K &m_kernel; + + public: + typedef int result_type; + + Point_dimension_d(const K &kernel) + : m_kernel(kernel) {} + + result_type operator()( + const Weighted_point_d & p) const + { + Construct_point_d cp = m_kernel.construct_point_d_object(); + return m_kernel.point_dimension_d_object()(cp(p)); + } + }; + + //=========================================================================== + // Object creation + //=========================================================================== + + Less_coordinate_d less_coordinate_d_object() const + { + return Less_coordinate_d(*this); + } + Contained_in_affine_hull_d contained_in_affine_hull_d_object() const + { + return Contained_in_affine_hull_d(*this); + } + Orientation_d orientation_d_object() const + { + return Orientation_d(*this); + } + Construct_flat_orientation_d construct_flat_orientation_d_object() const + { + return Construct_flat_orientation_d(*this); + } + In_flat_orientation_d in_flat_orientation_d_object() const + { + return In_flat_orientation_d(*this); + } + Compare_lexicographically_d compare_lexicographically_d_object() const + { + return Compare_lexicographically_d(*this); + } + Compute_coordinate_d compute_coordinate_d_object() const + { + return Compute_coordinate_d(*this); + } + Point_dimension_d point_dimension_d_object() const + { + return Point_dimension_d(*this); + } +}; + + +} //namespace CGAL + +#endif // CGAL_REGULAR_TRIANGULATION_TRAITS_ADAPTER_H diff --git a/Triangulation/include/CGAL/Triangulation.h b/Triangulation/include/CGAL/Triangulation.h index a826187e662..969e39ba01a 100644 --- a/Triangulation/include/CGAL/Triangulation.h +++ b/Triangulation/include/CGAL/Triangulation.h @@ -32,6 +32,7 @@ #include #include #include +#include #include #include @@ -141,20 +142,20 @@ public: typedef Maximal_dimension_ Maximal_dimension; typedef typename Geom_traits::Point_d Point; - typedef typename TDS::Vertex_handle Vertex_handle; - typedef typename TDS::Vertex_iterator Vertex_iterator; - typedef typename TDS::Vertex_const_handle Vertex_const_handle; - typedef typename TDS::Vertex_const_iterator Vertex_const_iterator; + typedef typename TDS::Vertex_handle Vertex_handle; + typedef typename TDS::Vertex_iterator Vertex_iterator; + typedef typename TDS::Vertex_const_handle Vertex_const_handle; + typedef typename TDS::Vertex_const_iterator Vertex_const_iterator; - typedef typename TDS::Full_cell_handle Full_cell_handle; - typedef typename TDS::Full_cell_iterator Full_cell_iterator; - typedef typename TDS::Full_cell_const_handle Full_cell_const_handle; - typedef typename TDS::Full_cell_const_iterator Full_cell_const_iterator; + typedef typename TDS::Full_cell_handle Full_cell_handle; + typedef typename TDS::Full_cell_iterator Full_cell_iterator; + typedef typename TDS::Full_cell_const_handle Full_cell_const_handle; + typedef typename TDS::Full_cell_const_iterator Full_cell_const_iterator; - typedef typename TDS::Facet_iterator Facet_iterator; + typedef typename TDS::Facet_iterator Facet_iterator; - typedef typename TDS::size_type size_type; - typedef typename TDS::difference_type difference_type; + typedef typename TDS::size_type size_type; + typedef typename TDS::difference_type difference_type; /// The type of location a new point is found lying on enum Locate_type @@ -184,18 +185,18 @@ public: protected: // DATA MEMBERS - Triangulation_ds tds_; - const Geom_traits kernel_; - Vertex_handle infinity_; - mutable std::vector orientations_; + Triangulation_ds tds_; + const Geom_traits kernel_; + Vertex_handle infinity_; + mutable std::vector orientations_; mutable boost::optional flat_orientation_; // The user can specify a Flat_orientation_d object to be used for // orienting simplices of a specific dimension // (= preset_flat_orientation_.first) // preset_flat_orientation_.first = numeric_limits::max() otherwise) - std::pair preset_flat_orientation_; + std::pair preset_flat_orientation_; // for stochastic walk in the locate() function: - mutable Random rng_; + mutable Random rng_; #ifdef CGAL_TRIANGULATION_STATISTICS mutable unsigned long walk_size_; #endif @@ -229,10 +230,38 @@ public: { return tds().index_of_covertex(f); } + +// - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - UTILITIES + + // A co-dimension 2 sub-simplex. called a Rotor because we can rotate + // the two "covertices" around the sub-simplex. Useful for traversing the + // boundary of a hole. NOT DOCUMENTED + typedef cpp11::tuple Rotor; + // Commented out because it was causing "internal compiler error" in MSVC + /*Full_cell_handle full_cell(const Rotor & r) const // NOT DOCUMENTED + { + return cpp11::get<0>(r); + } + int index_of_covertex(const Rotor & r) const // NOT DOCUMENTED + { + return cpp11::get<1>(r); + } + int index_of_second_covertex(const Rotor & r) const // NOT DOCUMENTED + { + return cpp11::get<2>(r); + }*/ + Rotor rotate_rotor(Rotor & r) // NOT DOCUMENTED... + { + int opposite = cpp11::get<0>(r)->mirror_index(cpp11::get<1>(r)); + Full_cell_handle s = cpp11::get<0>(r)->neighbor(cpp11::get<1>(r)); + int new_second = s->index(cpp11::get<0>(r)->vertex(cpp11::get<2>(r))); + return Rotor(s, new_second, opposite); + } + // - - - - - - - - - - - - - - - - - - - - - - - - CREATION / CONSTRUCTORS - Triangulation(int dim, const Geom_traits k = Geom_traits()) + Triangulation(int dim, const Geom_traits &k = Geom_traits()) : tds_(dim) , kernel_(k) , infinity_() @@ -503,7 +532,7 @@ public: bool is_infinite(const Facet & ft) const { Full_cell_const_handle s = full_cell(ft); - CGAL_precondition(s != Full_cell_handle()); + CGAL_precondition(s != Full_cell_const_handle()); if( is_infinite(s) ) return (s->vertex(index_of_covertex(ft)) != infinite_vertex()); return false; @@ -512,7 +541,7 @@ public: bool is_infinite(const Face & f) const { Full_cell_const_handle s = f.full_cell(); - CGAL_precondition(s != Full_cell_handle()); + CGAL_precondition(s != Full_cell_const_handle()); if( is_infinite(s) ) { Vertex_handle v; @@ -542,7 +571,7 @@ public: } template< typename OutputIterator > - OutputIterator incident_faces(Vertex_const_handle v, int d, OutputIterator out) + OutputIterator incident_faces(Vertex_const_handle v, int d, OutputIterator out) const { return tds().incident_faces(v, d, out); } @@ -604,7 +633,12 @@ public: return tds().new_full_cell(); } - Vertex_handle new_vertex(const Point & p) + Vertex_handle new_vertex() + { + return tds().new_vertex(); + } + + Vertex_handle new_vertex(const Point & p) { return tds().new_vertex(p); } @@ -623,13 +657,13 @@ public: protected: template< typename OrientationPredicate > - Full_cell_handle do_locate( const Point &, Locate_type &, Face &, Facet &, - Full_cell_handle start, - const OrientationPredicate & o) const; + Full_cell_handle do_locate(const Point &, Locate_type &, Face &, Facet &, + Full_cell_handle start, + const OrientationPredicate & o) const; public: - Full_cell_handle locate( const Point &, Locate_type &, Face &, Facet &, + Full_cell_handle locate(const Point &, Locate_type &, Face &, Facet &, Full_cell_handle start = Full_cell_handle()) const; - Full_cell_handle locate( const Point &, Locate_type &, Face &, Facet &, + Full_cell_handle locate(const Point &, Locate_type &, Face &, Facet &, Vertex_handle) const; Full_cell_handle locate(const Point & p, Full_cell_handle s = Full_cell_handle()) const; Full_cell_handle locate(const Point & p, Vertex_handle v) const; @@ -654,7 +688,7 @@ public: } return number_of_vertices() - n; } - Vertex_handle insert(const Point &, const Locate_type, const Face &, const Facet &, const Full_cell_handle); + Vertex_handle insert(const Point &, Locate_type, const Face &, const Facet &, Full_cell_handle); Vertex_handle insert(const Point &, Full_cell_handle start = Full_cell_handle()); Vertex_handle insert(const Point &, Vertex_handle); template< typename ForwardIterator > @@ -709,6 +743,43 @@ public: // make sure all full_cells have positive orientation void reorient_full_cells(); +protected: + // This is used in the |remove(v)| member function to manage sets of Full_cell_handles + template< typename FCH > + struct Full_cell_set : public std::vector + { + typedef std::vector Base_set; + using Base_set::begin; + using Base_set::end; + void make_searchable() + { // sort the full cell handles + std::sort(begin(), end()); + } + bool contains(const FCH & fch) const + { + return std::binary_search(begin(), end(), fch); + } + bool contains_1st_and_not_2nd(const FCH & fst, const FCH & snd) const + { + return ( ! contains(snd) ) && ( contains(fst) ); + } + }; + + void display_all_full_cells__debugging() const + { + std::cerr << "ALL FULL CELLS:" << std::endl; + for (Full_cell_const_iterator cit = full_cells_begin() ; + cit != full_cells_end() ; ++cit ) + { + std::cerr << std::hex << &*cit << ": "; + for (int jj = 0 ; jj <= current_dimension() ; ++jj) + std::cerr << (is_infinite(cit->vertex(jj)) ? 0xFFFFFFFF : (unsigned int)&*cit->vertex(jj)) << " - "; + std::cerr << std::dec << std::endl; + } + std::cerr << std::endl; + } + + }; // Triangulation<...> // = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = = @@ -722,17 +793,15 @@ Triangulation { if( current_dimension() < 1 ) return; + Full_cell_iterator sit = full_cells_begin(); Full_cell_iterator send = full_cells_end(); - while( sit != send ) + for ( ; sit != send ; ++sit) { - if( is_infinite(sit) && (1 == current_dimension()) ) + if( ! (is_infinite(sit) && (1 == current_dimension())) ) { - ++sit; - continue; + sit->swap_vertices(current_dimension() - 1, current_dimension()); } - sit->swap_vertices(current_dimension() - 1, current_dimension()); - ++sit; } } @@ -757,7 +826,7 @@ Triangulation template < class TT, class TDS > typename Triangulation::Vertex_handle Triangulation -::insert(const Point & p, const Locate_type lt, const Face & f, const Facet & ft, const Full_cell_handle s) +::insert(const Point & p, Locate_type lt, const Face & f, const Facet & ft, Full_cell_handle s) { switch( lt ) { @@ -853,14 +922,8 @@ Triangulation // infinite one... CGAL_precondition( is_infinite(s) ); CGAL_precondition( 1 == current_dimension() ); - int inf_v_index = s->index(infinite_vertex()); - bool swap = (0 == s->neighbor(inf_v_index)->index(s)); Vertex_handle v = tds().insert_in_full_cell(s); v->set_point(p); - if( swap ) - { - s->swap_vertices(0, 1); - } return v; } @@ -917,6 +980,36 @@ Triangulation CGAL_assertion( COPLANAR != o ); if( NEGATIVE == o ) reorient_full_cells(); + + + // We just inserted the second finite point and the right infinite + // cell is like : (inf_v, v), but we want it to be (v, inf_v) to be + // consistent with the rest of the cells + if (current_dimension() == 1) + { + // Is "inf_v_cell" the right infinite cell? + // Then inf_v_index should be 1 + if (inf_v_cell->neighbor(inf_v_index)->index(inf_v_cell) == 0 + && inf_v_index == 0) + { + inf_v_cell->swap_vertices( + current_dimension() - 1, current_dimension()); + } + // Otherwise, let's find the right infinite cell + else + { + inf_v_cell = inf_v_cell->neighbor((inf_v_index + 1) % 2); + inf_v_index = inf_v_cell->index(infinite_vertex()); + // Is "inf_v_cell" the right infinite cell? + // Then inf_v_index should be 1 + if (inf_v_cell->neighbor(inf_v_index)->index(inf_v_cell) == 0 + && inf_v_index == 0) + { + inf_v_cell->swap_vertices( + current_dimension() - 1, current_dimension()); + } + } + } } return v; } @@ -928,12 +1021,12 @@ template < class TT, class TDS > template< typename OrientationPredicate > typename Triangulation::Full_cell_handle Triangulation -::do_locate( const Point & p, // query point +::do_locate(const Point & p, // query point Locate_type & loc_type,// type of result (full_cell, face, vertex) Face & face,// the face containing the query in its interior (when appropriate) Facet & facet,// the facet containing the query in its interior (when appropriate) - const Full_cell_handle start// starting full_cell for the walk - , OrientationPredicate const& orientation_pred + Full_cell_handle start, // starting full_cell for the walk + OrientationPredicate const& orientation_pred ) const { const int cur_dim = current_dimension(); diff --git a/Triangulation/include/CGAL/Triangulation_data_structure.h b/Triangulation/include/CGAL/Triangulation_data_structure.h index b0033a3fff3..097faaace53 100644 --- a/Triangulation/include/CGAL/Triangulation_data_structure.h +++ b/Triangulation/include/CGAL/Triangulation_data_structure.h @@ -157,13 +157,13 @@ private: template < class Dim_tag > struct get_maximal_dimension { - static int value(const int D) { return D; } + static int value(int D) { return D; } }; // specialization template < int D > struct get_maximal_dimension > { - static int value(const int) { return D; } + static int value(int) { return D; } }; public: @@ -211,7 +211,7 @@ public: protected: - bool check_range(const int i) const + bool check_range(int i) const { if( current_dimension() < 0 ) { @@ -245,19 +245,19 @@ public: Full_cell_container & full_cells() { return full_cells_; } const Full_cell_container & full_cells() const { return full_cells_; } - Vertex_handle vertex(const Full_cell_handle s, const int i) const /* Concept */ + Vertex_handle vertex(Full_cell_handle s, int i) const /* Concept */ { CGAL_precondition(s != Full_cell_handle() && check_range(i)); return s->vertex(i); } - Vertex_const_handle vertex(const Full_cell_const_handle s, const int i) const /* Concept */ + Vertex_const_handle vertex(Full_cell_const_handle s, int i) const /* Concept */ { CGAL_precondition(s != Full_cell_handle() && check_range(i)); return s->vertex(i); } - bool is_vertex(const Vertex_const_handle & v) const /* Concept */ + bool is_vertex(Vertex_const_handle v) const /* Concept */ { if( Vertex_const_handle() == v ) return false; @@ -267,7 +267,7 @@ public: return v == vit; } - bool is_full_cell(const Full_cell_const_handle & s) const /* Concept */ + bool is_full_cell(Full_cell_const_handle s) const /* Concept */ { if( Full_cell_const_handle() == s ) return false; @@ -277,43 +277,43 @@ public: return s == sit; } - Full_cell_handle full_cell(const Vertex_handle v) const /* Concept */ + Full_cell_handle full_cell(Vertex_handle v) const /* Concept */ { CGAL_precondition(v != Vertex_handle()); return v->full_cell(); } - Full_cell_const_handle full_cell(const Vertex_const_handle v) const /* Concept */ + Full_cell_const_handle full_cell(Vertex_const_handle v) const /* Concept */ { CGAL_precondition(Vertex_const_handle() != v); return v->full_cell(); } - Full_cell_handle neighbor(const Full_cell_handle s, const int i) const /* Concept */ + Full_cell_handle neighbor(Full_cell_handle s, int i) const /* Concept */ { CGAL_precondition(Full_cell_handle() != s && check_range(i)); return s->neighbor(i); } - Full_cell_const_handle neighbor(const Full_cell_const_handle s, const int i) const/* Concept */ + Full_cell_const_handle neighbor(Full_cell_const_handle s, int i) const/* Concept */ { CGAL_precondition(Full_cell_const_handle() != s && check_range(i)); return s->neighbor(i); } - int mirror_index(const Full_cell_handle s, const int i) const /* Concept */ + int mirror_index(Full_cell_handle s, int i) const /* Concept */ { CGAL_precondition(Full_cell_handle() != s && check_range(i)); return s->mirror_index(i); } - int mirror_index(const Full_cell_const_handle s, const int i) const + int mirror_index(Full_cell_const_handle s, int i) const { CGAL_precondition(Full_cell_const_handle() != s && check_range(i)); /* Concept */ return s->mirror_index(i); } - int mirror_vertex(const Full_cell_handle s, const int i) const /* Concept */ + int mirror_vertex(Full_cell_handle s, int i) const /* Concept */ { CGAL_precondition(Full_cell_handle() != s && check_range(i)); return s->mirror_vertex(i); @@ -368,7 +368,7 @@ public: // NICE UPDATE OPERATIONS protected: - void do_insert_increase_dimension(const Vertex_handle, const Vertex_handle); + void do_insert_increase_dimension(Vertex_handle, Vertex_handle); public: // - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - REMOVALS @@ -381,9 +381,9 @@ public: Vertex_handle insert_in_face(const Face &); /* Concept */ Vertex_handle insert_in_facet(const Facet &); /* Concept */ template< typename Forward_iterator > - Vertex_handle insert_in_hole(Forward_iterator, const Forward_iterator, Facet); /* Concept */ + Vertex_handle insert_in_hole(Forward_iterator, Forward_iterator, Facet); /* Concept */ template< typename Forward_iterator, typename OutputIterator > - Vertex_handle insert_in_hole(Forward_iterator, const Forward_iterator, Facet, OutputIterator); /* Concept */ + Vertex_handle insert_in_hole(Forward_iterator, Forward_iterator, Facet, OutputIterator); /* Concept */ template< typename OutputIterator > Full_cell_handle insert_in_tagged_hole(Vertex_handle, Facet, OutputIterator); @@ -420,7 +420,6 @@ private: void clear_visited_marks(Full_cell_handle) const; // - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - DANGEROUS UPDATE OPERATIONS - // - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - DANGEROUS UPDATE OPERATIONS private: @@ -449,13 +448,13 @@ public: dcur_ = -2; } - void set_current_dimension(const int d) /* Concept */ + void set_current_dimension(int d) /* Concept */ { CGAL_precondition(-2<=d && d<=maximal_dimension()); dcur_ = d; } - Full_cell_handle new_full_cell(const Full_cell_handle s) + Full_cell_handle new_full_cell(Full_cell_handle s) { return full_cells_.emplace(*s); } @@ -497,7 +496,7 @@ public: vertices_.erase(v); } - void associate_vertex_with_full_cell(Full_cell_handle s, const int i, Vertex_handle v) /* Concept */ + void associate_vertex_with_full_cell(Full_cell_handle s, int i, Vertex_handle v) /* Concept */ { CGAL_precondition(check_range(i)); CGAL_precondition(s != Full_cell_handle()); @@ -557,7 +556,7 @@ public: const Triangulation_data_structure & tds_; public: Incident_full_cell_traversal_predicate(const Triangulation_data_structure & tds, - const Face & f) + const Face & f) : f_(f), tds_(tds) { dim_ = f.face_dimension(); @@ -582,7 +581,7 @@ public: const Triangulation_data_structure & tds_; public: Star_traversal_predicate(const Triangulation_data_structure & tds, - const Face & f) + const Face & f) : f_(f), tds_(tds) { dim_ = f.face_dimension(); @@ -610,28 +609,28 @@ public: OutputIterator star(const Face &, OutputIterator) const; /* Concept */ #ifndef CGAL_CFG_NO_CPP0X_DEFAULT_TEMPLATE_ARGUMENTS_FOR_FUNCTION_TEMPLATES template< typename OutputIterator, typename Comparator = std::less > - OutputIterator incident_upper_faces(Vertex_const_handle v, const int dim, OutputIterator out, Comparator cmp = Comparator()) + OutputIterator incident_upper_faces(Vertex_const_handle v, int dim, OutputIterator out, Comparator cmp = Comparator()) { return incident_faces(v, dim, out, cmp, true); } template< typename OutputIterator, typename Comparator = std::less > - OutputIterator incident_faces(Vertex_const_handle, const int, OutputIterator, Comparator = Comparator(), bool = false); + OutputIterator incident_faces(Vertex_const_handle, int, OutputIterator, Comparator = Comparator(), bool = false) const; #else template< typename OutputIterator, typename Comparator > - OutputIterator incident_upper_faces(Vertex_const_handle v, const int dim, OutputIterator out, Comparator cmp = Comparator()) + OutputIterator incident_upper_faces(Vertex_const_handle v, int dim, OutputIterator out, Comparator cmp = Comparator()) { return incident_faces(v, dim, out, cmp, true); } template< typename OutputIterator > - OutputIterator incident_upper_faces(Vertex_const_handle v, const int dim, OutputIterator out) + OutputIterator incident_upper_faces(Vertex_const_handle v, int dim, OutputIterator out) { return incident_faces(v, dim, out, std::less(), true); } template< typename OutputIterator, typename Comparator > - OutputIterator incident_faces(Vertex_const_handle, const int, OutputIterator, Comparator = Comparator(), bool = false); + OutputIterator incident_faces(Vertex_const_handle, int, OutputIterator, Comparator = Comparator(), bool = false) const; template< typename OutputIterator > - OutputIterator incident_faces(Vertex_const_handle, const int, OutputIterator, - std::less = std::less(), bool = false); + OutputIterator incident_faces(Vertex_const_handle, int, OutputIterator, + std::less = std::less(), bool = false) const; #endif // - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - INPUT / OUTPUT @@ -726,8 +725,8 @@ template< class Dim, class Vb, class Fcb > template< typename OutputIterator > OutputIterator Triangulation_data_structure -::incident_faces(Vertex_const_handle v, const int dim, OutputIterator out, - std::less cmp, bool upper_faces) +::incident_faces(Vertex_const_handle v, int dim, OutputIterator out, + std::less cmp, bool upper_faces) const { return incident_faces >(v, dim, out, cmp, upper_faces); } @@ -737,7 +736,7 @@ template< class Dim, class Vb, class Fcb > template< typename OutputIterator, typename Comparator > OutputIterator Triangulation_data_structure -::incident_faces(Vertex_const_handle v, const int dim, OutputIterator out, Comparator cmp, bool upper_faces) +::incident_faces(Vertex_const_handle v, int dim, OutputIterator out, Comparator cmp, bool upper_faces) const { CGAL_precondition( 0 < dim ); if( dim >= current_dimension() ) @@ -791,13 +790,13 @@ Triangulation_data_structure // init state for enumerating all candidate faces: internal::Combination_enumerator f_idx(dim, v_idx + 1, current_dimension()); Face f(*s); - f.set_index(0, v_idx); + f.set_index(0, sorted_idx[v_idx]); while( ! f_idx.end() ) { - // check if face has already been found for( int i = 0; i < dim; ++i ) f.set_index(1 + i, sorted_idx[f_idx[i]]); - face_set.insert(f); + face_set.insert(f); // checks if face has already been found + // compute next sorted face (lexicographic enumeration) ++f_idx; } @@ -892,8 +891,7 @@ Triangulation_data_structure if( v_idx != current_dimension() ) { (*it)->swap_vertices(v_idx, current_dimension()); - if( ( ! (*it)->has_vertex(star) ) || (current_dimension() > 2) ) - (*it)->swap_vertices(current_dimension() - 2, current_dimension() - 1); + (*it)->swap_vertices(current_dimension() - 2, current_dimension() - 1); } (*it)->set_vertex(current_dimension(), Vertex_handle()); (*it)->set_neighbor(current_dimension(), Full_cell_handle()); @@ -971,10 +969,10 @@ Triangulation_data_structure CGAL_assertion_msg(is_boundary_facet(f), "starting facet should be on the hole boundary"); const int cur_dim = current_dimension(); - Full_cell_handle new_s; + Full_cell_handle new_s; - std::queue task_queue; - task_queue.push( + std::queue task_queue; + task_queue.push( IITH_task(f, mirror_index(full_cell(f), index_of_covertex(f))) ); while (!task_queue.empty()) @@ -1003,7 +1001,7 @@ Triangulation_data_structure associate_vertex_with_full_cell(new_s, facet_index, v); set_neighbors(new_s, facet_index, - neighbor(old_s, facet_index), + outside_neighbor, mirror_index(old_s, facet_index)); // add the new full_cell to the list of new full_cells @@ -1120,7 +1118,7 @@ Triangulation_data_structure template void Triangulation_data_structure -::do_insert_increase_dimension(const Vertex_handle x, const Vertex_handle star) +::do_insert_increase_dimension(Vertex_handle x, Vertex_handle star) { Full_cell_handle start = full_cells_begin(); Full_cell_handle swap_me; @@ -1142,11 +1140,6 @@ void Triangulation_data_structure for( int k = 1; k <= cur_dim; ++k ) associate_vertex_with_full_cell(S_new, k, vertex(S, k - 1)); } - else if( cur_dim == 2 ) - { // if cur. dim. is 2, we must take care of the 'rightmost' infinite vertex. - if( S->mirror_index(S->index(star)) == 0 ) - swap_me = S; - } } // now we setup the neighbors set_visited(start, false); @@ -1526,7 +1519,7 @@ operator>>(std::istream & is, Triangulation_data_structure & tr) // - the neighbors of each full_cell by their index in the preceding list { typedef Triangulation_data_structure TDS; - typedef typename TDS::Vertex_handle Vertex_handle; + typedef typename TDS::Vertex_handle Vertex_handle; // read current dimension and number of vertices std::size_t n; @@ -1576,8 +1569,8 @@ operator<<(std::ostream & os, const Triangulation_data_structure // - the neighbors of each full_cell by their index in the preceding list { typedef Triangulation_data_structure TDS; - typedef typename TDS::Vertex_const_handle Vertex_handle; - typedef typename TDS::Vertex_const_iterator Vertex_iterator; + typedef typename TDS::Vertex_const_handle Vertex_handle; + typedef typename TDS::Vertex_const_iterator Vertex_iterator; // outputs dimension and number of vertices std::size_t n = tr.number_of_vertices(); @@ -1598,6 +1591,8 @@ operator<<(std::ostream & os, const Triangulation_data_structure for( Vertex_iterator it = tr.vertices_begin(); it != tr.vertices_end(); ++it, ++i ) { os << *it; // write the vertex + if (is_ascii(os)) + os << std::endl; index_of_vertex[it] = i; } CGAL_assertion( (std::size_t) i == n ); diff --git a/Triangulation/include/CGAL/Triangulation_ds_vertex.h b/Triangulation/include/CGAL/Triangulation_ds_vertex.h index 3a18e1dffa7..5c680629851 100644 --- a/Triangulation/include/CGAL/Triangulation_ds_vertex.h +++ b/Triangulation/include/CGAL/Triangulation_ds_vertex.h @@ -64,7 +64,6 @@ public: /// Set 's' as an incident full_cell void set_full_cell(Full_cell_handle s) /* Concept */ { - CGAL_precondition( Full_cell_handle() != s ); full_cell_ = s; } diff --git a/Triangulation/include/CGAL/Triangulation_full_cell.h b/Triangulation/include/CGAL/Triangulation_full_cell.h index ce6d17c2fe0..1268d6fbce3 100644 --- a/Triangulation/include/CGAL/Triangulation_full_cell.h +++ b/Triangulation/include/CGAL/Triangulation_full_cell.h @@ -69,22 +69,6 @@ public: Triangulation_full_cell(const Self & s) : Base(s), data_(s.data_) {} - Point circumcenter() const - { - TriangulationTraits pct; - Vertex_handle_const_iterator vhit = vertices_begin(); - while( vhit != vertices_end() ) - { - if( *vhit == Vertex_const_handle() ) - { - CGAL_warning_msg(false, "too few points; can not compute circumcenter."); - return Point(); - } - ++vhit; - } - return pct.center_of_sphere_d_object()(points_begin(), points_end()); - } - const Data & data() const { return data_; diff --git a/Triangulation/include/CGAL/Triangulation_vertex.h b/Triangulation/include/CGAL/Triangulation_vertex.h index 74e93f28e6c..0a87d6ba3cd 100644 --- a/Triangulation/include/CGAL/Triangulation_vertex.h +++ b/Triangulation/include/CGAL/Triangulation_vertex.h @@ -25,7 +25,6 @@ #include #include -#include namespace CGAL { diff --git a/Triangulation/include/CGAL/internal/Triangulation/utilities.h b/Triangulation/include/CGAL/internal/Triangulation/utilities.h index 63e43ceef70..7c4d49e57fa 100644 --- a/Triangulation/include/CGAL/internal/Triangulation/utilities.h +++ b/Triangulation/include/CGAL/internal/Triangulation/utilities.h @@ -101,7 +101,7 @@ public: template< class T > struct Compare_points_for_perturbation { - typedef typename T::Point_d Point; + typedef typename T::Geom_traits::Point_d Point; const T & t_; @@ -122,8 +122,8 @@ public: template< class T > struct Point_from_pointer { - typedef const typename T::Point_d * argument_type; - typedef const typename T::Point_d result_type; + typedef const typename T::Geom_traits::Point_d * argument_type; + typedef const typename T::Geom_traits::Point_d result_type; result_type & operator()(argument_type & x) const { return (*x); diff --git a/Triangulation/test/Triangulation/CMakeLists.txt b/Triangulation/test/Triangulation/CMakeLists.txt index 0e1651af197..36ca91eea95 100644 --- a/Triangulation/test/Triangulation/CMakeLists.txt +++ b/Triangulation/test/Triangulation/CMakeLists.txt @@ -1,7 +1,6 @@ # Created by the script cgal_create_cmake_script # This is the CMake script for compiling a CGAL application. - project( Triangulation_Tests ) cmake_minimum_required(VERSION 2.8.11) @@ -21,18 +20,18 @@ if ( CGAL_FOUND ) include_directories (BEFORE "../../include") include_directories (BEFORE "include") + create_single_source_cgal_program( "test_triangulation.cpp" ) create_single_source_cgal_program( "test_delaunay.cpp" ) + create_single_source_cgal_program( "test_regular.cpp" ) create_single_source_cgal_program( "test_tds.cpp" ) create_single_source_cgal_program( "test_torture.cpp" ) - create_single_source_cgal_program( "test_triangulation.cpp" ) + create_single_source_cgal_program( "test_insert_if_in_star.cpp" ) else() message(STATUS "NOTICE: Some of the executables in this directory need Eigen 3.1 (or greater) and will not be compiled.") endif() else() - message(STATUS "This program requires the CGAL library, and will not be compiled.") - endif() diff --git a/Triangulation/test/Triangulation/test_delaunay.cpp b/Triangulation/test/Triangulation/test_delaunay.cpp index 6f9b8e729f2..bd2bb695ed3 100644 --- a/Triangulation/test/Triangulation/test_delaunay.cpp +++ b/Triangulation/test/Triangulation/test_delaunay.cpp @@ -34,102 +34,105 @@ void test(const int d, const string & type, const int N) typedef CGAL::Random_points_in_cube_d Random_points_iterator; - DC pc(d); + DC dt(d); cerr << "\nBuilding Delaunay triangulation of (" << type << d << ") dimension with " << N << " points"; - assert(pc.empty()); + assert(dt.empty()); vector points; - CGAL::Random rng; - Random_points_iterator rand_it(d, 2.0, rng); + //CGAL::Random rng; + //Random_points_iterator rand_it(d, 2.0, rng); //CGAL::cpp11::copy_n(rand_it, N, back_inserter(points)); - vector coords(d); + srand(10); for( int i = 0; i < N; ++i ) { + vector coords(d); for( int j = 0; j < d; ++j ) - coords[j] = rand() % 100000; + coords[j] = static_cast(rand() % 100000)/10000; points.push_back(Point(d, coords.begin(), coords.end())); } - pc.insert(points.begin(), points.end()); + dt.insert(points.begin(), points.end()); cerr << "\nChecking topology and geometry..."; - assert( pc.is_valid() ); + assert( dt.is_valid() ); cerr << "\nTraversing finite full_cells... "; size_t nbfs(0), nbis(0); - Finite_full_cell_const_iterator fsit = pc.finite_full_cells_begin(); - while( fsit != pc.finite_full_cells_end() ) + Finite_full_cell_const_iterator fsit = dt.finite_full_cells_begin(); + while( fsit != dt.finite_full_cells_end() ) ++fsit, ++nbfs; cerr << nbfs << " + "; vector infinite_full_cells; - pc.tds().incident_full_cells(pc.infinite_vertex(), back_inserter(infinite_full_cells)); + dt.tds().incident_full_cells(dt.infinite_vertex(), back_inserter(infinite_full_cells)); nbis = infinite_full_cells.size(); cerr << nbis << " = " << (nbis+nbfs) - << " = " << pc.number_of_full_cells(); - cerr << "\nThe triangulation has current dimension " << pc.current_dimension(); - CGAL_assertion( pc.number_of_full_cells() == nbis+nbfs); + << " = " << dt.number_of_full_cells(); + cerr << "\nThe triangulation has current dimension " << dt.current_dimension(); + CGAL_assertion( dt.number_of_full_cells() == nbis+nbfs); cerr << "\nTraversing finite vertices... "; size_t nbfv(0); - Finite_vertex_iterator fvit = pc.finite_vertices_begin(); - while( fvit != pc.finite_vertices_end() ) + Finite_vertex_iterator fvit = dt.finite_vertices_begin(); + while( fvit != dt.finite_vertices_end() ) ++fvit, ++nbfv; cerr << nbfv < 1 ) + if( dt.maximal_dimension() > 1 ) { typedef vector Faces; Faces edges; back_insert_iterator out(edges); - pc.tds().incident_faces(pc.infinite_vertex(), 1, out); + dt.tds().incident_faces(dt.infinite_vertex(), 1, out); cout << "\nThere are " << edges.size() << " vertices on the convex hull."; edges.clear(); } - else // pc.maximal_dimension() == 1 + else // dt.maximal_dimension() == 1 { typedef vector Cells; Cells cells; back_insert_iterator out(cells); - pc.tds().incident_full_cells(pc.infinite_vertex(), out); + dt.tds().incident_full_cells(dt.infinite_vertex(), out); cout << "\nThere are " << cells.size() << " vertices on the convex hull."; cells.clear(); } // Remove all ! - cerr << "\nBefore removal: " << pc.number_of_vertices() << " vertices. After: "; + cerr << "\nBefore removal: " << dt.number_of_vertices() << " vertices. After: "; random_shuffle(points.begin(), points.end()); - pc.remove(points.begin(), points.end()); - assert( pc.is_valid() ); - cerr << pc.number_of_vertices() << " vertices."; - // assert( pc.empty() ); NOT YET ! + dt.remove(points.begin(), points.end()); + assert( dt.is_valid() ); + cerr << dt.number_of_vertices() << " vertices."; + // assert( dt.empty() ); NOT YET ! // CLEAR - pc.clear(); - assert( -1 == pc.current_dimension() ); - assert( pc.empty() ); - assert( pc.is_valid() ); + dt.clear(); + assert( -1 == dt.current_dimension() ); + assert( dt.empty() ); + assert( dt.is_valid() ); } template< int D > void go(const int N) { - //typedef CGAL::Epick_d FK; typedef CGAL::Epick_d > FK; typedef CGAL::Delaunay_triangulation Triangulation; - //test(D, "dynamic", N); test(D, "static", N); + + typedef CGAL::Epick_d FK_dyn; + typedef CGAL::Delaunay_triangulation Triangulation_dyn; + test(D, "dynamic", N); } int main(int argc, char **argv) { srand(static_cast(time(NULL))); - int N = 100; + int N = 10; if( argc > 1 ) N = atoi(argv[1]); - go<5>(N); - go<4>(N); - go<3>(N); - go<2>(N); - go<1>(N); + //go<5>(N); + go<4>(N); + go<3>(N); + go<2>(N); + go<1>(N); cerr << endl; return 0; diff --git a/Triangulation/test/Triangulation/test_insert_if_in_star.cpp b/Triangulation/test/Triangulation/test_insert_if_in_star.cpp new file mode 100644 index 00000000000..4aeb5204566 --- /dev/null +++ b/Triangulation/test/Triangulation/test_insert_if_in_star.cpp @@ -0,0 +1,92 @@ +#include +#include +#include +#include +#include + +#include +#include +#include +#include +#include + +using namespace std; + +template +void test(const int d, const string & type, const int N) +{ + typedef typename RTri::Vertex_handle Vertex_handle; + typedef typename RTri::Point Point; + typedef typename RTri::Bare_point Bare_point; + + typedef CGAL::Random_points_in_cube_d Random_points_iterator; + + RTri rt(d); + RTri rt_star_only(d); + cerr << "\nBuilding Regular triangulation of (" << type << d + << ") dimension with " << N << " points\n"; + assert(rt.empty()); + assert(rt_star_only.empty()); + + srand(static_cast(time(NULL))); + + // Insert first point (0, 0...) + vector coords(d); + for( int j = 0; j < d; ++j ) + coords[j] = 0; + + Point p = Point( + Bare_point(d, coords.begin(), coords.end()), + static_cast(rand() % 10000)/100000); + + rt.insert(p); + Vertex_handle first_vertex = rt_star_only.insert(p); + + // Insert the other points + for( int i = 1 ; i < N ; ++i ) + { + for( int j = 0; j < d; ++j ) + coords[j] = 10.*(rand() % RAND_MAX)/RAND_MAX - 5.; + + p = Point( + Bare_point(d, coords.begin(), coords.end()), + static_cast(rand() % 10000)/1000000); + + rt.insert(p); + rt_star_only.insert_if_in_star(p, first_vertex); + } + + cerr << "\nChecking topology and geometry..." + << (rt.is_valid(true) ? "OK.\n" : "Error.\n"); + + cerr << "\nThe triangulation using 'insert' has current dimension " << rt.current_dimension() + << " and " << rt.number_of_full_cells() << " full cells\n"; + + cerr << "\nThe triangulation using 'insert_if_in_star' has current dimension " << rt.current_dimension() + << " and " << rt_star_only.number_of_full_cells() << " full cells\n"; + + // Export + if (d <= 3) + { + std::ofstream off_stream_all("data/test_insert_all.off"); + CGAL::export_triangulation_to_off(off_stream_all, rt); + std::ofstream off_stream_star_only("data/test_insert_if_in_star.off"); + CGAL::export_triangulation_to_off(off_stream_star_only, rt_star_only); + } +} + +template< int D > +void go(const int N) +{ + //typedef CGAL::Epick_d FK; + typedef CGAL::Epick_d > FK; + typedef CGAL::Regular_triangulation Triangulation; + //test(D, "dynamic", N); + test(D, "static", N); +} + +int main(int argc, char **argv) +{ + go<2>(100); + return 0; +} diff --git a/Triangulation/test/Triangulation/test_regular.cpp b/Triangulation/test/Triangulation/test_regular.cpp new file mode 100644 index 00000000000..2ffd84892b6 --- /dev/null +++ b/Triangulation/test/Triangulation/test_regular.cpp @@ -0,0 +1,289 @@ +#include +#include +#include +#include + +#include +#include +#include +#include +#include +#include + +using namespace std; + +template +void test(const int d, const string & type, const int N) +{ + typedef typename RTri::Full_cell_handle Full_cell_handle; + typedef typename RTri::Face Face; + typedef typename RTri::Point Point; + typedef typename RTri::Bare_point Bare_point; + typedef typename RTri::Finite_full_cell_const_iterator Finite_full_cell_const_iterator; + typedef typename RTri::Finite_vertex_iterator Finite_vertex_iterator; + + typedef CGAL::Random_points_in_cube_d Random_points_iterator; + + RTri rt(d); + cerr << "\nBuilding Regular triangulation of (" << type << d << ") dimension with " << N << " points"; + assert(rt.empty()); + + vector points; + + srand(10); + for( int i = 0; i < N; ++i ) + { + vector coords(d); + for( int j = 0; j < d; ++j ) + coords[j] = static_cast(rand() % 100000)/10000; + points.push_back(Point( + Bare_point(d, coords.begin(), coords.end()), + static_cast(rand() % 100000)/100000 + )); + } + rt.insert(points.begin(), points.end()); + cerr << "\nChecking topology and geometry..."; + assert( rt.is_valid(true) ); + + cerr << "\nTraversing finite full_cells... "; + size_t nbfs(0), nbis(0); + Finite_full_cell_const_iterator fsit = rt.finite_full_cells_begin(); + while( fsit != rt.finite_full_cells_end() ) + ++fsit, ++nbfs; + cerr << nbfs << " + "; + vector infinite_full_cells; + rt.tds().incident_full_cells(rt.infinite_vertex(), back_inserter(infinite_full_cells)); + nbis = infinite_full_cells.size(); + cerr << nbis << " = " << (nbis+nbfs) + << " = " << rt.number_of_full_cells(); + cerr << "\nThe triangulation has current dimension " << rt.current_dimension(); + CGAL_assertion( rt.number_of_full_cells() == nbis+nbfs); + + cerr << "\nTraversing finite vertices... "; + size_t nbfv(0); + Finite_vertex_iterator fvit = rt.finite_vertices_begin(); + while( fvit != rt.finite_vertices_end() ) + ++fvit, ++nbfv; + cerr << nbfv < 1 ) + { + typedef vector Faces; + Faces edges; + back_insert_iterator out(edges); + rt.tds().incident_faces(rt.infinite_vertex(), 1, out); + cout << "\nThere are " << edges.size() << " vertices on the convex hull."; + edges.clear(); + } + else // rt.maximal_dimension() == 1 + { + typedef vector Cells; + Cells cells; + back_insert_iterator out(cells); + rt.tds().incident_full_cells(rt.infinite_vertex(), out); + cout << "\nThere are " << cells.size() << " vertices on the convex hull."; + cells.clear(); + } + + // Remove all ! + cerr << "\nBefore removal: " << rt.number_of_vertices() << " vertices. After: "; + random_shuffle(points.begin(), points.end()); + rt.remove(points.begin(), points.end()); + assert( rt.is_valid() ); + //std::cerr << ((rt.is_valid(true)) ? "VALID!" : "NOT VALID :(") << std::endl; + cerr << rt.number_of_vertices() << " vertices."; + // assert( rt.empty() ); NOT YET ! + // CLEAR + rt.clear(); + assert( -1 == rt.current_dimension() ); + assert( rt.empty() ); + assert( rt.is_valid() ); + //std::cerr << ((rt.is_valid(true)) ? "VALID!" : "NOT VALID :(") << std::endl; +} + +template< int D > +void go(const int N) +{ + typedef CGAL::Epick_d > FK; + typedef CGAL::Regular_triangulation Triangulation; + test(D, "static", N); + + typedef CGAL::Epick_d FK_dyn; + typedef CGAL::Regular_triangulation Triangulation_dyn; + test(D, "dynamic", N); +} + +void test_inserting_points_at_the_same_position() +{ + const int DIM = 5; + typedef CGAL::Epick_d > FK; + typedef CGAL::Regular_triangulation RTri; + + typedef RTri::Vertex_handle Vertex_handle; + typedef RTri::Full_cell_handle Full_cell_handle; + typedef RTri::Face Face; + typedef RTri::Point Point; + typedef RTri::Bare_point Bare_point; + typedef RTri::Finite_full_cell_const_iterator Finite_full_cell_const_iterator; + typedef RTri::Finite_vertex_iterator Finite_vertex_iterator; + + RTri rt(DIM); + + cerr << "\nTesting insertion of points at the same position"; + assert(rt.empty()); + + std::vector pt; + pt.push_back(1.2); + pt.push_back(20.3); + pt.push_back(10.); + pt.push_back(8.); + pt.push_back(7.1); + + //======================= + + // First point + Vertex_handle vh = rt.insert(Point(Bare_point(pt.begin(), pt.end()), 0.)); + assert(rt.number_of_vertices() == 1); + assert(rt.number_of_hidden_vertices() == 0); + assert(vh != Vertex_handle()); + + // Same point + Vertex_handle vh2 = rt.insert(Point(Bare_point(pt.begin(), pt.end()), 0.)); + assert(rt.number_of_vertices() == 1); + assert(rt.number_of_hidden_vertices() == 0); + assert(vh == vh2); + + // Same position, bigger weight + vh = rt.insert(Point(Bare_point(pt.begin(), pt.end()), 1.3)); + assert(rt.number_of_vertices() == 1); + assert(rt.number_of_hidden_vertices() == 1); + assert(vh != Vertex_handle()); + + // Same point + vh2 = rt.insert(Point(Bare_point(pt.begin(), pt.end()), 1.3)); + assert(rt.number_of_vertices() == 1); + assert(rt.number_of_hidden_vertices() == 1); + assert(vh == vh2); + + // Same position, lower weight + vh = rt.insert(Point(Bare_point(pt.begin(), pt.end()), 1.1)); + assert(rt.number_of_vertices() == 1); + assert(rt.number_of_hidden_vertices() == 2); + assert(vh == Vertex_handle()); + + //======================= + + // New position + pt[3] = 0.78; + vh = rt.insert(Point(Bare_point(pt.begin(), pt.end()), 0.)); + assert(rt.number_of_vertices() == 2); + assert(rt.number_of_hidden_vertices() == 2); + assert(vh != Vertex_handle()); + + // Same point + vh2 = rt.insert(Point(Bare_point(pt.begin(), pt.end()), 0.)); + assert(rt.number_of_vertices() == 2); + assert(rt.number_of_hidden_vertices() == 2); + assert(vh == vh2); + + // Same position, bigger weight + vh = rt.insert(Point(Bare_point(pt.begin(), pt.end()), 1.3)); + assert(rt.number_of_vertices() == 2); + assert(rt.number_of_hidden_vertices() == 3); + assert(vh != Vertex_handle()); + + // Same point + vh2 = rt.insert(Point(Bare_point(pt.begin(), pt.end()), 1.3)); + assert(rt.number_of_vertices() == 2); + assert(rt.number_of_hidden_vertices() == 3); + assert(vh == vh2); + + // Same position, lower weight + vh = rt.insert(Point(Bare_point(pt.begin(), pt.end()), 1.1)); + assert(rt.number_of_vertices() == 2); + assert(rt.number_of_hidden_vertices() == 4); + assert(vh == Vertex_handle()); + + //======================= + + // New position + pt[4] = 1.78; + vh = rt.insert(Point(Bare_point(pt.begin(), pt.end()), 0.2)); + assert(rt.number_of_vertices() == 3); + assert(rt.number_of_hidden_vertices() == 4); + assert(vh != Vertex_handle()); + + //======================= + + // New position + pt[1] = 1.78; + vh = rt.insert(Point(Bare_point(pt.begin(), pt.end()), 0.8)); + assert(rt.number_of_vertices() == 4); + assert(rt.number_of_hidden_vertices() == 4); + assert(vh != Vertex_handle()); + + //======================= + + // New position + pt[2] = 1.78; + vh = rt.insert(Point(Bare_point(pt.begin(), pt.end()), 0.25)); + assert(rt.number_of_vertices() == 5); + assert(rt.number_of_hidden_vertices() == 4); + assert(vh != Vertex_handle()); + + //======================= + + // New position + pt[0] = 12.13; + vh = rt.insert(Point(Bare_point(pt.begin(), pt.end()), 1.25)); + assert(rt.number_of_vertices() == 6); + assert(rt.number_of_hidden_vertices() == 4); + assert(vh != Vertex_handle()); + + // Same position, bigger weight + vh = rt.insert(Point(Bare_point(pt.begin(), pt.end()), 1.3)); + assert(rt.number_of_vertices() == 6); + assert(rt.number_of_hidden_vertices() == 5); + assert(vh != Vertex_handle()); + + //======================= + + // New position + pt[0] = 9.13; + vh = rt.insert(Point(Bare_point(pt.begin(), pt.end()), 1.25)); + assert(rt.number_of_vertices() == 7); + assert(rt.number_of_hidden_vertices() == 5); + assert(vh != Vertex_handle()); + + //======================= + + cerr << "\nChecking topology and geometry..."; + assert(rt.is_valid(true)); + + rt.clear(); + assert(-1 == rt.current_dimension()); + assert(rt.empty()); + assert(rt.is_valid()); + //std::cerr << ((rt.is_valid(true)) ? "VALID!" : "NOT VALID :(") << std::endl; +} + + +int main(int argc, char **argv) +{ + srand(static_cast(time(NULL))); + int N = 10; + if( argc > 1 ) + N = atoi(argv[1]); + + test_inserting_points_at_the_same_position(); + + //go<5>(N); + go<4>(N); + go<3>(N); + go<2>(N); + go<1>(N); + + cerr << endl; + return 0; +} diff --git a/Triangulation/test/Triangulation/test_torture.cpp b/Triangulation/test/Triangulation/test_torture.cpp index 3ec215278b9..d75b560791a 100644 --- a/Triangulation/test/Triangulation/test_torture.cpp +++ b/Triangulation/test/Triangulation/test_torture.cpp @@ -82,11 +82,7 @@ void test(const int D, const int d, const int N, bool no_transform) for( int j = 0; j < d; ++j ) coords[i] = coords[i] + (*pit)[j] * aff[j][i]; } -#ifdef USE_NEW_KERNEL - points.push_back(Point(coords)); // this is for New_kernel_d -#else - points.push_back(Point(D, coords.begin(), coords.end())); // this is for Old_kernel_d -#endif + points.push_back(Point(D, coords.begin(), coords.end())); } assert( dc.is_valid() ); cout << " Inserting " << points.size() << " points."; diff --git a/Triangulation/test/Triangulation/test_triangulation.cpp b/Triangulation/test/Triangulation/test_triangulation.cpp index eecc6338171..4e225333054 100644 --- a/Triangulation/test/Triangulation/test_triangulation.cpp +++ b/Triangulation/test/Triangulation/test_triangulation.cpp @@ -53,7 +53,6 @@ void test(const int d, const string & type, int N) Finite_full_cell_const_iterator fsit = tri.finite_full_cells_begin(); while( fsit != tri.finite_full_cells_end() ) { - fsit->circumcenter(); ++fsit, ++nbfs; } cerr << nbfs << " + "; @@ -124,10 +123,11 @@ int main(int argc, char **argv) int N = 1000; if( argc > 1 ) N = atoi(argv[1]); - go<5>(N); - go<3>(N); - go<2>(N); - go<1>(N); + //go<5>(N); + go<4>(N); + go<3>(N); + go<2>(N); + go<1>(N); cerr << std::endl; return 0; diff --git a/Triangulation_2/applications/Triangulation_2/data/points.cin b/Triangulation_2/applications/Triangulation_2/data/points.cin new file mode 100644 index 00000000000..5e78e2349bc --- /dev/null +++ b/Triangulation_2/applications/Triangulation_2/data/points.cin @@ -0,0 +1,19 @@ +0 0 6.28953 +-2.85086 -0.471442 6.12896 +1.90972 0.101219 0.988689 +0.637771 2.59367 5.80372 +2.22209 0.903198 2.19478 +-0.487202 -2.71506 4.90996 +1.1193 -1.91787 2.99626 +1.54714 0.109831 0 +0.44556 -2.73047 4.48142 +0.427936 1.28495 6.23624 +-2.67212 0.766674 5.29623 +1.5763 -1.59828 2.58905 +-0.476603 2.2546 6.04797 +1.57172 -0.514711 6.11405 +1.84528 2.10139 5.53936 +-2.99827 -0.101677 5.92246 +-0.482122 -2.39584 4.44264 +-2.25558 -1.492 6.23448 +0.128475 -1.75125 3.18916 \ No newline at end of file diff --git a/Triangulation_2/applications/Triangulation_2/points_to_RT2_to_off.cpp b/Triangulation_2/applications/Triangulation_2/points_to_RT2_to_off.cpp new file mode 100644 index 00000000000..c2989d5a902 --- /dev/null +++ b/Triangulation_2/applications/Triangulation_2/points_to_RT2_to_off.cpp @@ -0,0 +1,28 @@ +#include +#include +#include +#include +#include + +#include + +typedef CGAL::Exact_predicates_inexact_constructions_kernel K; +typedef CGAL::Regular_triangulation_filtered_traits_2 Traits; +typedef CGAL::Regular_triangulation_2 Regular_triangulation; + +int main() +{ + std::ifstream in("data/points.cin"); + + Regular_triangulation::Weighted_point wp; + std::vector wpoints; + + while(in >> wp) + wpoints.push_back(wp); + + Regular_triangulation rt(wpoints.begin(), wpoints.end()); + CGAL_assertion(rt.is_valid(true)); + std::ofstream off_stream("data/rt2.off"); + CGAL::export_triangulation_2_to_off(off_stream, rt); + return 0; +} diff --git a/Triangulation_2/include/CGAL/IO/Triangulation_off_ostream_2.h b/Triangulation_2/include/CGAL/IO/Triangulation_off_ostream_2.h new file mode 100644 index 00000000000..a4cd70e8e0d --- /dev/null +++ b/Triangulation_2/include/CGAL/IO/Triangulation_off_ostream_2.h @@ -0,0 +1,79 @@ +// Copyright (c) 2014 INRIA Sophia-Antipolis (France). +// All rights reserved. +// +// This file is part of CGAL (www.cgal.org); you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public License as +// published by the Free Software Foundation; either version 3 of the License, +// or (at your option) any later version. +// +// Licensees holding a valid commercial license may use this file in +// accordance with the commercial license agreement provided with the software. +// +// This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE +// WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE. +// +// $URL: $ +// $Id: $ +// +// Author(s) : Clement Jamin + + +#ifndef CGAL_TRIANGULATION_OFF_OSTREAM_2_H +#define CGAL_TRIANGULATION_OFF_OSTREAM_2_H + +#include +#include +#include + +namespace CGAL { + +template < class GT, class TDS > +std::ostream & +export_triangulation_2_to_off(std::ostream & os, + const Triangulation_2 & tr) +{ + typedef Triangulation_2 Tr; + typedef typename Tr::Vertex_handle Vertex_handle; + typedef typename Tr::Vertex_iterator Vertex_iterator; + typedef typename Tr::Finite_vertices_iterator Finite_vertex_iterator; + typedef typename Tr::Finite_faces_iterator Finite_faces_iterator; + + size_t n = tr.number_of_vertices(); + + std::stringstream output; + + // write the vertices + std::map index_of_vertex; + int i = 0; + for(Finite_vertex_iterator it = tr.finite_vertices_begin(); + it != tr.finite_vertices_end(); ++it, ++i) + { + output << it->point().x() << " " << it->point().y() << " 0" << std::endl; + index_of_vertex[it.base()] = i; + } + CGAL_assertion( i == n ); + + size_t number_of_triangles = 0; + + for (Finite_faces_iterator fit = tr.finite_faces_begin() ; + fit != tr.finite_faces_end() ; ++fit) + { + output << "3 " + << index_of_vertex[fit->vertex(0)] << " " + << index_of_vertex[fit->vertex(1)] << " " + << index_of_vertex[fit->vertex(2)] + << std::endl; + ++number_of_triangles; + } + + os << "OFF \n" + << n << " " + << number_of_triangles << " 0\n" + << output.str(); + + return os; +} + +} //namespace CGAL + +#endif // CGAL_TRIANGULATION_OFF_OSTREAM_2_H \ No newline at end of file diff --git a/Triangulation_3/applications/Triangulation_3/data/points - extreme weights.cin b/Triangulation_3/applications/Triangulation_3/data/points - extreme weights.cin new file mode 100644 index 00000000000..10c0791ad4c --- /dev/null +++ b/Triangulation_3/applications/Triangulation_3/data/points - extreme weights.cin @@ -0,0 +1,10 @@ +0.0071 1.6899 2.521 0 +0.3272 1.3694 3.15 0.05 +1.3697 1.8296 2.654 0.1 +-10.6722 0.3012 0.1548 1000.15 +1.1726 0.1899 0.3658 0.2 +0.4374 20.8541 1.45894 2000.25 +2.5923 0.1904 0.6971 0.3 +10.3083 2.5462 1.3658 1000.35 +1.4981 1.3929 2.949 0.4 +2.1304 2.055 0.6597455 1.45 \ No newline at end of file diff --git a/Triangulation_3/applications/Triangulation_3/data/points - no weights.cin b/Triangulation_3/applications/Triangulation_3/data/points - no weights.cin new file mode 100644 index 00000000000..9dcea4973dc --- /dev/null +++ b/Triangulation_3/applications/Triangulation_3/data/points - no weights.cin @@ -0,0 +1,10 @@ +0.0071 1.6899 2.521 0 +0.3272 1.3694 3.15 0 +1.3697 1.8296 2.654 0 +-10.6722 0.3012 0.1548 0 +1.1726 0.1899 0.3658 0 +0.4374 20.8541 1.45894 0 +2.5923 0.1904 0.6971 0 +10.3083 2.5462 1.3658 0 +1.4981 1.3929 2.949 0 +2.1304 2.055 0.6597455 0 \ No newline at end of file diff --git a/Triangulation_3/applications/Triangulation_3/data/points.cin b/Triangulation_3/applications/Triangulation_3/data/points.cin new file mode 100644 index 00000000000..9dcea4973dc --- /dev/null +++ b/Triangulation_3/applications/Triangulation_3/data/points.cin @@ -0,0 +1,10 @@ +0.0071 1.6899 2.521 0 +0.3272 1.3694 3.15 0 +1.3697 1.8296 2.654 0 +-10.6722 0.3012 0.1548 0 +1.1726 0.1899 0.3658 0 +0.4374 20.8541 1.45894 0 +2.5923 0.1904 0.6971 0 +10.3083 2.5462 1.3658 0 +1.4981 1.3929 2.949 0 +2.1304 2.055 0.6597455 0 \ No newline at end of file diff --git a/Triangulation_3/applications/Triangulation_3/points_to_RT3_to_off.cpp b/Triangulation_3/applications/Triangulation_3/points_to_RT3_to_off.cpp new file mode 100644 index 00000000000..4dffb3ea265 --- /dev/null +++ b/Triangulation_3/applications/Triangulation_3/points_to_RT3_to_off.cpp @@ -0,0 +1,26 @@ +#include +#include +#include +#include + +#include + +typedef CGAL::Exact_predicates_inexact_constructions_kernel K; +typedef CGAL::Regular_triangulation_euclidean_traits_3 Traits; +typedef CGAL::Regular_triangulation_3 Regular_triangulation; + +int main() +{ + std::ifstream in("data/points.cin"); + + Regular_triangulation::Weighted_point wp; + std::vector wpoints; + + while(in >> wp) + wpoints.push_back(wp); + + Regular_triangulation rt(wpoints.begin(), wpoints.end()); + std::ofstream off_stream("data/rt3.off"); + CGAL::export_triangulation_3_to_off(off_stream, rt); + return 0; +} diff --git a/Triangulation_3/benchmark/Triangulation_3/Triangulation_benchmark_3.cpp b/Triangulation_3/benchmark/Triangulation_3/Triangulation_benchmark_3.cpp index f0db573b038..a7c0d66c373 100644 --- a/Triangulation_3/benchmark/Triangulation_3/Triangulation_benchmark_3.cpp +++ b/Triangulation_3/benchmark/Triangulation_3/Triangulation_benchmark_3.cpp @@ -303,9 +303,6 @@ void benchmark_remove() Time_accumulator tt(time); tr.remove(&vhs[0], &vhs[NUM_VERTICES_TO_REMOVE - 1]); ++iterations; - //std::cout<<"\b\b\b\b\b\b"< +#include +#include + +namespace CGAL { + +template < class GT, class TDS > +std::ostream & +export_triangulation_3_to_off(std::ostream & os, + const Triangulation_3 & tr, + bool export_surface_only = false) +{ + typedef Triangulation_3 Tr; + typedef typename Tr::Vertex_handle Vertex_handle; + typedef typename Tr::Vertex_iterator Vertex_iterator; + typedef typename Tr::Finite_vertices_iterator Finite_vertex_iterator; + typedef typename Tr::All_cells_iterator Cells_iterator; + typedef typename Tr::Finite_cells_iterator Finite_cells_iterator; + + size_t n = tr.number_of_vertices(); + + std::stringstream output; + + // write the vertices + std::map index_of_vertex; + int i = 0; + for(Finite_vertex_iterator it = tr.finite_vertices_begin(); + it != tr.finite_vertices_end(); ++it, ++i) + { + output << it->point().x() << " " + << it->point().y() << " " + << it->point().z() << std::endl; + index_of_vertex[it.base()] = i; + } + CGAL_assertion( i == n ); + + size_t number_of_triangles = 0; + + if (export_surface_only) + { + for (Cells_iterator cit = tr.cells_begin() ; + cit != tr.cells_end() ; ++cit) + { + if (tr.is_infinite(cit)) + { + output << "3 "; + for (int i = 0 ; i < 4 ; ++i) + { + if (!tr.is_infinite(cit->vertex(i))) + output << index_of_vertex[cit->vertex(i)] << " "; + } + output << std::endl; + ++number_of_triangles; + } + } + } + else + { + for (Finite_cells_iterator cit = tr.finite_cells_begin() ; + cit != tr.finite_cells_end() ; ++cit) + { + output << "3 " + << index_of_vertex[cit->vertex(0)] << " " + << index_of_vertex[cit->vertex(1)] << " " + << index_of_vertex[cit->vertex(2)] + << std::endl; + output << "3 " + << index_of_vertex[cit->vertex(0)] << " " + << index_of_vertex[cit->vertex(2)] << " " + << index_of_vertex[cit->vertex(3)] + << std::endl; + output << "3 " + << index_of_vertex[cit->vertex(1)] << " " + << index_of_vertex[cit->vertex(2)] << " " + << index_of_vertex[cit->vertex(3)] + << std::endl; + output << "3 " + << index_of_vertex[cit->vertex(0)] << " " + << index_of_vertex[cit->vertex(1)] << " " + << index_of_vertex[cit->vertex(3)] + << std::endl; + number_of_triangles += 4; + } + } + + os << "OFF \n" + << n << " " + << number_of_triangles << " 0\n" + << output.str(); + + return os; +} + +} //namespace CGAL + +#endif // CGAL_TRIANGULATION_OFF_OSTREAM_3_H \ No newline at end of file