songsenand
/
cgal
mirror of https://github.com/CGAL/cgal
1
0
Fork 0
Code Issues Packages Projects Releases Wiki Activity

Use Marc's new Kernel_d

This commit is contained in:
Clement Jamin 2014-03-25 11:51:23 +01:00
parent 94ae9ec16a
commit 885741d466
6 changed files with 130 additions and 100 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

30
Triangulation/benchmark/Triangulation/delaunay.cpp
View File

@ -1,6 +1,4 @@
#include <CGAL/Cartesian_d.h>
//#include <CGAL/Simple_cartesian_d.h>
//#include <CGAL/Filtered_kernel_d.h>
#include <CGAL/Epick_d.h>
#include <CGAL/Delaunay_triangulation.h>
#include <CGAL/point_generators_d.h>
#include <CGAL/Timer.h>
@ -38,35 +36,33 @@ void test(const int d, const std::string & type, const int N)
std::cout << " Delaunay triangulation of "<<N<<" points in dim "<<d<< std::flush;
dt.insert(points.begin(), points.end());
std::cout << " done in "<<cost.time()<<" seconds." << std::endl;
int nbfc= dt.number_of_finite_full_cells();
int nbc= dt.number_of_full_cells();
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;
}
template< int D, typename RT >
template< int D >
void go(const int N)
{
//typedef CGAL::Simple_cartesian_d<RT, D> K;
typedef CGAL::Cartesian_d<RT> K;
//typedef CGAL::Filtered_kernel_d<K> FK;
typedef CGAL::Epick_d<CGAL::Dimension_tag<D> > K;
typedef CGAL::Delaunay_triangulation<K> Triangulation;
test<Triangulation>(D, "static", N);
}
int main(int argc, char **argv)
{
srand48(time(NULL));
srand(static_cast<unsigned int>(time(NULL)));
int N = 100; if( argc > 1 ) N = atoi(argv[1]);
go<2, double>(N);
go<3, double>(N);
go<4, double>(N);
go<5, double>(N);
go<6, double>(N);
go<7, double>(N);
go<8, double>(N);
go<2>(N);
go<3>(N);
go<4>(N);
go<5>(N);
go<6>(N);
go<7>(N);
go<8>(N);
return 0;

60
Triangulation/include/CGAL/Delaunay_triangulation.h
View File

@ -30,13 +30,12 @@ template< typename DCTraits, typename _TDS = Default >
class Delaunay_triangulation
: public Triangulation<DCTraits,
typename Default::Get<_TDS, Triangulation_data_structure<
typename Ambient_dimension<typename DCTraits::Point_d>::type,
typename DCTraits::Dimension,
Triangulation_vertex<DCTraits>,
Triangulation_full_cell<DCTraits> >
>::type >
{
typedef typename Ambient_dimension<typename DCTraits::Point_d>::type
Maximal_dimension_;
typedef typename DCTraits::Dimension Maximal_dimension_;
typedef typename Default::Get<_TDS, Triangulation_data_structure<
Maximal_dimension_,
Triangulation_vertex<DCTraits>,
@ -45,14 +44,36 @@ class Delaunay_triangulation
typedef Triangulation<DCTraits, TDS> Base;
typedef Delaunay_triangulation<DCTraits, _TDS> Self;
typedef typename DCTraits::Side_of_oriented_subsphere_d
Side_of_oriented_subsphere_d;
typedef typename DCTraits::Side_of_oriented_sphere_d
Side_of_oriented_sphere_d;
typedef typename DCTraits::Coaffine_orientation_d
Coaffine_orientation_d;
typedef typename DCTraits::Orientation_d Orientation_d;
//*** Side_of_oriented_subsphere_d ***
using Base::Flat_orientation_d;
using Base::Construct_flat_orientation_d;
typedef typename DCTraits::In_flat_side_of_oriented_sphere_d In_flat_side_of_oriented_sphere_d;
// Wrapper
struct Side_of_oriented_subsphere_d
{
boost::optional<Flat_orientation_d>* fop;
Construct_flat_orientation_d cfo;
In_flat_side_of_oriented_sphere_d ifsoos;
Side_of_oriented_subsphere_d(
boost::optional<Flat_orientation_d>& x,
Construct_flat_orientation_d const&y,
In_flat_side_of_oriented_sphere_d const&z)
: fop(&x), cfo(y), ifsoos(z) {}
template<class Iter>
CGAL::Orientation operator()(Iter a, Iter b, const Point & p)const
{
if(!*fop)
*fop=cfo(a,b);
return ifsoos(fop->get(),a,b,p);
}
};
public: // PUBLIC NESTED TYPES
typedef DCTraits Geom_traits;
@ -84,13 +105,13 @@ public: // PUBLIC NESTED TYPES
protected: // DATA MEMBERS
Side_of_oriented_subsphere_d side_of_oss_;
public:
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::star;
//using Base::incident_full_cells;
@ -145,7 +166,7 @@ public:
// - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - CREATION
Delaunay_triangulation(const int dim, const Geom_traits k = Geom_traits())
: Base(dim, k), side_of_oss_()
: Base(dim, k)
{
}
@ -154,16 +175,15 @@ public:
// - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - ACCESS
// Not Documented
Side_of_oriented_subsphere_d & side_of_oriented_subsphere_predicate()
Side_of_oriented_subsphere_d side_of_oriented_subsphere_predicate() const
{
return side_of_oss_;
return Side_of_oriented_subsphere_d (
flat_orientation_,
geom_traits().construct_flat_orientation_d_object(),
geom_traits().in_flat_side_of_oriented_sphere_d_object()
);
}
// Not Documented
const Side_of_oriented_subsphere_d & side_of_oriented_subsphere_predicate() const
{
return side_of_oss_;
}
// - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - REMOVALS
@ -307,6 +327,8 @@ public:
private:
// Some internal types to shorten notation
using Base::Coaffine_orientation_d;
using Base::flat_orientation_;
typedef Conflict_predicate<Coaffine_orientation_d, Side_of_oriented_subsphere_d>
Conflict_pred_in_subspace;
typedef Conflict_predicate<Orientation_d, Side_of_oriented_sphere_d>
@ -653,8 +675,7 @@ Delaunay_triangulation<DCTraits, TDS>
CGAL_precondition( current_dimension() >= 0 );
tds().remove_decrease_dimension(v, infinite_vertex());
// reset the predicates:
coaffine_orientation_predicate() = geom_traits().coaffine_orientation_d_object();
side_of_oriented_subsphere_predicate() = geom_traits().side_of_oriented_subsphere_d_object();
reset_flat_orientation();
if( 1 <= current_dimension() )
{
// FIXME: infinite vertex is NOT at index 0 a priori.
@ -712,8 +733,7 @@ Delaunay_triangulation<DCTraits, TDS>
CGAL_precondition( current_dimension() < maximal_dimension() );
Vertex_handle v = tds().insert_increase_dimension(infinite_vertex());
// reset the predicates:
coaffine_orientation_predicate() = geom_traits().coaffine_orientation_d_object();
side_of_oriented_subsphere_predicate() = geom_traits().side_of_oriented_subsphere_d_object();
reset_flat_orientation();
v->set_point(p);
if( current_dimension() >= 1 )
{

72
Triangulation/include/CGAL/Triangulation.h
View File

@ -37,17 +37,50 @@ namespace CGAL {
template < class TriangulationTraits, class TDS_ = Default >
class Triangulation
{
typedef typename Ambient_dimension<typename TriangulationTraits::Point_d>::type
Maximal_dimension_;
typedef typename TriangulationTraits::Dimension Maximal_dimension_;
typedef typename Default::Get<TDS_, Triangulation_data_structure
< Maximal_dimension_,
Triangulation_vertex<TriangulationTraits>,
Triangulation_full_cell<TriangulationTraits> >
>::type TDS;
typedef Triangulation<TriangulationTraits, TDS_> Self;
protected:
typedef typename TriangulationTraits::Flat_orientation_d Flat_orientation_d;
typedef typename TriangulationTraits::Construct_flat_orientation_d Construct_flat_orientation_d;
typedef typename TriangulationTraits::In_flat_orientation_d In_flat_orientation_d;
// Wrapper
struct Coaffine_orientation_d
{
boost::optional<Flat_orientation_d>* fop;
Construct_flat_orientation_d cfo;
In_flat_orientation_d ifo;
Coaffine_orientation_d(
boost::optional<Flat_orientation_d>& x,
Construct_flat_orientation_d const&y,
In_flat_orientation_d const&z)
: fop(&x), cfo(y), ifo(z) {}
template<class Iter>
CGAL::Orientation operator()(Iter a, Iter b) const
{
// CJTODO: temporary code
if (!*fop)
*fop = cfo(a,b);
return ifo(fop->get(),a,b);
// CJTODO: replace with this code:
/*if (*fop)
return ifo(fop->get(),a,b);
*fop = cfo(a,b);
CGAL_assertion(ifo(fop->get(),a,b) == CGAL::POSITIVE);
return CGAL::POSITIVE;*/
}
};
void reset_flat_orientation(){flat_orientation_=boost::none;}
typedef typename TriangulationTraits::Coaffine_orientation_d
Coaffine_orientation_d;
typedef typename TriangulationTraits::Orientation_d
Orientation_d;
@ -111,7 +144,7 @@ protected: // DATA MEMBERS
const Geom_traits kernel_;
Vertex_handle infinity_;
mutable std::vector<Oriented_side> orientations_;
Coaffine_orientation_d coaffine_orientation_;
mutable boost::optional<Flat_orientation_d> flat_orientation_;
// for stochastic walk in the locate() function:
mutable Random rng_;
#ifdef CGAL_TRIANGULATION_STATISTICS
@ -185,7 +218,7 @@ public:
// A full_cell has at most 1 + maximal_dimension() facets:
orientations_.resize(1 + maximal_dimension());
// Our coaffine orientation predicates HAS state member variables
coaffine_orientation_predicate() = geom_traits().coaffine_orientation_d_object();
reset_flat_orientation();
}
~Triangulation() {}
@ -215,14 +248,13 @@ public:
tds().set_visited(s, b);
} */
Coaffine_orientation_d & coaffine_orientation_predicate()
Coaffine_orientation_d coaffine_orientation_predicate() const
{
return coaffine_orientation_;
}
const Coaffine_orientation_d & coaffine_orientation_predicate() const
{
return coaffine_orientation_;
return Coaffine_orientation_d (
flat_orientation_,
geom_traits().construct_flat_orientation_d_object(),
geom_traits().in_flat_orientation_d_object()
);
}
const Triangulation_ds & tds() const
@ -482,7 +514,7 @@ public:
// A full_cell has at most 1 + maximal_dimension() facets:
orientations_.resize(1 + maximal_dimension());
// Our coaffine orientation predicates HAS state member variables
coaffine_orientation_predicate() = geom_traits().coaffine_orientation_d_object();
reset_flat_orientation();
#ifdef CGAL_TRIANGULATION_STATISTICS
walk_size_ = 0;
#endif
@ -519,7 +551,7 @@ protected:
template< typename OrientationPredicate >
Full_cell_handle do_locate( const Point &, Locate_type &, Face &, Facet &,
Full_cell_handle start = Full_cell_handle(),
OrientationPredicate & o = geom_traits().orientation_d_object()) const;
const OrientationPredicate & o = geom_traits().orientation_d_object()) const;
public:
Full_cell_handle locate( const Point &, Locate_type &, Face &, Facet &,
Full_cell_handle start = Full_cell_handle()) const;
@ -579,11 +611,11 @@ public:
{
Triangulation & t_;
const Point & p_;
OrientationPredicate & ori_;
OrientationPredicate const& ori_;
int cur_dim_;
public:
Outside_convex_hull_traversal_predicate(Triangulation & t, const Point & p,
OrientationPredicate & ori)
OrientationPredicate const& ori)
: t_(t), p_(p), ori_(ori), cur_dim_(t.current_dimension()) {}
// FUTURE change parameter to const reference
bool operator()(Facet f) const
@ -798,7 +830,7 @@ Triangulation<TT, TDS>
CGAL_precondition( current_dimension() < maximal_dimension() );
Vertex_handle v = tds().insert_increase_dimension(infinite_vertex());
// reset the orientation predicate:
coaffine_orientation_predicate() = geom_traits().coaffine_orientation_d_object();
reset_flat_orientation();
v->set_point(p);
if( current_dimension() >= 1 )
{
@ -824,7 +856,7 @@ Triangulation<TT, TDS>
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 & orientation_pred
, OrientationPredicate const& orientation_pred
) const
{
const int cur_dim = current_dimension();
@ -919,7 +951,7 @@ Triangulation<TT, TDS>
points_begin(s) + cur_dim + 1);
// restore the correct point for vertex |i| of the full_cell
s->vertex(i)->set_point(backup); /**/
s->vertex(i)->set_point(backup); */
if( orientations_[i] != NEGATIVE )
{

21
Triangulation/test/Triangulation/delaunay.cpp
View File

@ -1,7 +1,5 @@
#include <CGAL/Epick_d.h>
#include <CGAL/point_generators_d.h>
//#include <CGAL/Simple_cartesian_d.h>
#include <CGAL/Cartesian_d.h>
#include <CGAL/Filtered_kernel_d.h>
#include <CGAL/Delaunay_triangulation.h>
#include <CGAL/spatial_sort.h>
#include <CGAL/algorithm.h>
@ -43,11 +41,11 @@ void test(const int d, const string & type, const int N)
Random_points_iterator rand_it(d, 2.0, rng);
//CGAL::copy_n(rand_it, N, back_inserter(points));
vector<int> coords(d);
vector<int> coords(d);
for( int i = 0; i < N; ++i )
{
for( int j = 0; j < d; ++j )
coords[j] = lrand48() % 100000;
coords[j] = rand() % 100000;
points.push_back(Point(d, coords.begin(), coords.end()));
}
pc.insert(points.begin(), points.end());
@ -112,19 +110,20 @@ void test(const int d, const string & type, const int N)
template< int D >
void go(const int N)
{
typedef double RT;
//typedef double RT;
//typedef CGAL::Gmpq RT;
typedef CGAL::Cartesian_d<RT> K;
//typedef CGAL::Cartesian_d<RT> K;
//typedef CGAL::Simple_cartesian_d<RT, D> K;
typedef CGAL::Filtered_kernel_d<K> FK;
//typedef CGAL::Filtered_kernel_d<K> FK;
typedef CGAL::Epick_d<CGAL::Dimension_tag<D> > FK;
typedef CGAL::Delaunay_triangulation<FK> Triangulation;
test<Triangulation>(D, "dynamic", N);
//test<Triangulation>(D, "static", N);
//test<Triangulation>(D, "dynamic", N);
test<Triangulation>(D, "static", N);
}
int main(int argc, char **argv)
{
srand48(time(NULL));
srand(static_cast<unsigned int>(time(NULL)));
int N = 100;
if( argc > 1 )
N = atoi(argv[1]);

26
Triangulation/test/Triangulation/torture.cpp
View File

@ -1,12 +1,6 @@
#include <CGAL/Epick_d.h>
#include <CGAL/internal/Combination_enumerator.h>
#include <CGAL/point_generators_d.h>
//#define USE_NEW_KERNEL
#ifndef USE_NEW_KERNEL
#include <CGAL/Cartesian_d.h> // this is for Old_kernel_d
#else
#include <CGAL/Simple_cartesian_d.h> // this is for New_kernel_d
#endif
#include <CGAL/Filtered_kernel_d.h>
#include <CGAL/Delaunay_triangulation.h>
#include <CGAL/algorithm.h>
#include <tilted_grid.h>
@ -50,7 +44,7 @@ void test(const int D, const int d, const int N, bool no_transform)
CGAL::internal::Combination_enumerator combi(k, 0, D-1);
if( no_transform )
{ // choose a random set of axes:
int nb = lrand48() % 1000;
int nb = rand() % 1000;
for( int i = 0; i < nb; ++i )
{
++combi;
@ -70,7 +64,7 @@ void test(const int D, const int d, const int N, bool no_transform)
{
int c(0);
while( 0 == c )
c = (lrand48() % 11) - 5;
c = (rand() % 11) - 5;
aff[j].push_back(c);
}
}
@ -116,16 +110,8 @@ void test(const int D, const int d, const int N, bool no_transform)
template< int D >
void go(const int N, const int nb_trials)
{
typedef double RT;
//typedef CGAL::Quotient<CGAL::MP_Float> RT;
//typedef CGAL::Gmpq RT;
#ifndef USE_NEW_KERNEL
typedef CGAL::Cartesian_d<RT> K; // this is for Old_kernel_d
#else
typedef CGAL::Simple_cartesian_d<RT, D> K; // this is for New_kernel_d
#endif
typedef CGAL::Filtered_kernel_d<K> FK;
typedef CGAL::Delaunay_triangulation<FK> Triangulation;
typedef CGAL::Epick_d<CGAL::Dimension_tag<D> > K;
typedef CGAL::Delaunay_triangulation<K> Triangulation;
for( int d = 0; d <= D; ++d )
{
cout << "\n\n** Delaunay of " << d
@ -150,7 +136,7 @@ int main(int argc, char **argv)
if( argc > 3 )
rand_init = atoi(argv[3]);
std::cout<<argv[0]<<" "<<N<<" "<<nb_trials<<" "<<rand_init<<std::endl;
srand48(rand_init);
srand(rand_init);
go<1>(N, nb_trials);
go<2>(N, nb_trials);
go<3>(N, nb_trials);

21
Triangulation/test/Triangulation/triangulation.cpp
View File

@ -1,6 +1,5 @@
#include <CGAL/Epick_d.h>
#include <CGAL/point_generators_d.h>
#include <CGAL/Cartesian_d.h>
#include <CGAL/Filtered_kernel_d.h>
#include <CGAL/Triangulation.h>
#include <CGAL/algorithm.h>
#include <tilted_grid.h>
@ -34,10 +33,10 @@ void test(const int d, const string & type, int N)
assert(tri.empty());
vector<RT> coords(d);
vector<Point> points;
vector<Point> points;
CGAL::Random rng;
Random_points_iterator rand_it(d, 1.0, rng);
CGAL::copy_n(rand_it, N, std::back_inserter(points));
Random_points_iterator rand_it(d, 1.0, rng);
CGAL::copy_n(rand_it, N, std::back_inserter(points));
cerr << '\n' << points.size() << " points in the grid.";
@ -107,17 +106,15 @@ void test(const int d, const string & type, int N)
template< int D >
void go(int N)
{
typedef double RT;
typedef CGAL::Cartesian_d<RT> K;
typedef CGAL::Filtered_kernel_d<K> FK;
typedef CGAL::Triangulation<FK> Triangulation;
//test<Triangulation>(D, "static", N);
test<Triangulation>(D, "dynamic", N);
typedef CGAL::Epick_d<CGAL::Dimension_tag<D> > K;
typedef CGAL::Triangulation<K> Triangulation;
test<Triangulation>(D, "static", N);
//test<Triangulation>(D, "dynamic", N);
}
int main(int argc, char **argv)
{
srand48(time(NULL));
srand(static_cast<unsigned int>(time(NULL)));
int N = 1000;
if( argc > 1 )
N = atoi(argv[1]);