From 0d6a9e8367fc04f243da5d6006f7efea0cbdeac9 Mon Sep 17 00:00:00 2001 From: Jane Tournois Date: Thu, 16 May 2019 16:46:43 +0200 Subject: [PATCH] add a benchmark for simplex traverser and cell traverser --- .../benchmark/Triangulation_3/CMakeLists.txt | 29 ++++ .../segment_traverser_benchmark.cpp | 131 ++++++++++++++++++ .../CGAL/Triangulation_segment_traverser_3.h | 6 + 3 files changed, 166 insertions(+) create mode 100644 Triangulation_3/benchmark/Triangulation_3/CMakeLists.txt create mode 100644 Triangulation_3/benchmark/Triangulation_3/segment_traverser_benchmark.cpp diff --git a/Triangulation_3/benchmark/Triangulation_3/CMakeLists.txt b/Triangulation_3/benchmark/Triangulation_3/CMakeLists.txt new file mode 100644 index 00000000000..c7c212da684 --- /dev/null +++ b/Triangulation_3/benchmark/Triangulation_3/CMakeLists.txt @@ -0,0 +1,29 @@ + +project( Triangulation_3_Benchmark ) + +cmake_minimum_required(VERSION 2.8.11) + +find_package(CGAL QUIET) + +if ( CGAL_FOUND ) + + include( ${CGAL_USE_FILE} ) + +# find_package( TBB QUIET ) + + include( CGAL_CreateSingleSourceCGALProgram ) + + include_directories (BEFORE "../../include") + + create_single_source_cgal_program( "incident_edges.cpp" ) + create_single_source_cgal_program( "simple.cpp" ) + create_single_source_cgal_program( "simple_2.cpp" ) + create_single_source_cgal_program( "Triangulation_benchmark_3.cpp" ) + create_single_source_cgal_program( "segment_traverser_benchmark.cpp" ) + +else() + + message(STATUS "This program requires the CGAL library, and will not be compiled.") + +endif() + diff --git a/Triangulation_3/benchmark/Triangulation_3/segment_traverser_benchmark.cpp b/Triangulation_3/benchmark/Triangulation_3/segment_traverser_benchmark.cpp new file mode 100644 index 00000000000..580d9d7e649 --- /dev/null +++ b/Triangulation_3/benchmark/Triangulation_3/segment_traverser_benchmark.cpp @@ -0,0 +1,131 @@ +// Copyright (c) 2019 GeometryFactory (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$ +// +// + +#include +#include + +#include +#include + +#include +#include + +#include +#include + +// Define the kernel. +typedef CGAL::Exact_predicates_exact_constructions_kernel Epeck; +typedef CGAL::Exact_predicates_inexact_constructions_kernel Epick; + + +template +void bench_segment_traverser(const int nb_queries, + const int nbv, + const double rad, + CGAL::Random& rng) +{ + typedef CGAL::Delaunay_triangulation_3 DT; + typedef CGAL::Triangulation_segment_simplex_iterator_3
Simplex_traverser; + typedef CGAL::Triangulation_segment_cell_iterator_3
Cell_traverser; + typedef typename DT::Point_3 Point_3; + typedef typename DT::Cell Cell; + + std::cout << "\nBench :\t " << nb_queries << " queries," << std::endl + << "\t in triangulation of size " << nbv << std::endl + << "\t lying in the sphere (-" << rad << "; " << rad << ")" << std::endl + << "\t and Kernel : " << typeid(Kernel()).name() << std::endl; + + // Collect random points to build a triangulation + std::vector points(nbv); + int i = 0; + while (i < nbv) + { + points[i++] = Point_3(rng.get_double(-rad, rad), + rng.get_double(-rad, rad), + rng.get_double(-rad, rad)); + } + + // Construct the Delaunay triangulation. + DT dt(points.begin(), points.end()); + assert(dt.is_valid()); + + // Collect random segments for queries in [-2*rad, 2*rad] + std::vector segments(2 * nb_queries); + i = 0; + while (i < 2 * nb_queries) + { + segments[i++] = Point_3(rng.get_double(-2 * rad, 2 * rad), + rng.get_double(-2 * rad, 2 * rad), + rng.get_double(-2 * rad, 2 * rad)); + } + + CGAL::Timer timer_st, timer_ct; + timer_st.reset(); + timer_ct.reset(); + + for (i = 0; i < nb_queries; ++i) + { + //Simplex traverser + timer_st.start(); + Simplex_traverser st(dt, segments[2*i], segments[2*i + 1]); + + // Count the number of finite cells traversed. + unsigned int inf = 0, fin = 0; + for (; st != st.end(); ++st) + { + Cell c = st.cell(); +// if (dt.is_infinite(c)) ++inf; +// else ++fin; + } + timer_st.stop(); + + //Cell traverser + timer_ct.start(); + Cell_traverser ct(dt, segments[2*i], segments[2*i + 1]); + + // Count the number of finite cells traversed. + inf = 0, fin = 0; + for (; ct != ct.end(); ++ct) + { + Cell c = ct.cell(); +// if (dt.is_infinite(c)) ++inf; +// else ++fin; + } + timer_ct.stop(); + } + + std::cout << "Simplex traverser took " << timer_st.time() << " seconds." << std::endl; + std::cout << "Cell traverser took " << timer_ct.time() << " seconds." << std::endl; +} + +int main(int argc, char* argv[]) +{ + CGAL::Random rng; + std::cout << "CGAL::Random seed is " << rng.get_seed() << std::endl; + + //nb of segments tested + const int nb_queries = (argc > 1) ? atoi(argv[1]) : 100; + //nb of vertices in triangulation + const int nbv = (argc > 2) ? atoi(argv[2]) : 1000; + //radius of sphere including the triangulation + const double rad = (argc > 3) ? atof(argv[3]) : 10.; + +// bench_segment_traverser(nb_queries, nbv, rad, rng); + bench_segment_traverser(nb_queries, nbv, rad, rng); +} diff --git a/Triangulation_3/include/CGAL/Triangulation_segment_traverser_3.h b/Triangulation_3/include/CGAL/Triangulation_segment_traverser_3.h index 0ec766475de..396df844073 100644 --- a/Triangulation_3/include/CGAL/Triangulation_segment_traverser_3.h +++ b/Triangulation_3/include/CGAL/Triangulation_segment_traverser_3.h @@ -505,6 +505,7 @@ public: // \{ typedef typename SCI::Vertex_handle Vertex_handle;//< defines the type of a handle for a vertex in the triangulation typedef typename SCI::Cell_handle Cell_handle; //< defines the type of a handle for a cell in the triangulation. + typedef typename SCI::Cell Cell; //< defines the type of a handle for a cell in the triangulation. typedef typename SCI::Triangulation::Edge Edge; //< defines the type of an edge in the triangulation. typedef typename SCI::Triangulation::Facet Facet; //< defines the type of a facet in the triangulation. typedef typename SCI::Locate_type Locate_type; //< defines the simplex type returned from location. @@ -813,6 +814,11 @@ public: bool is_facet() const { return _curr_simplex.dimension() == 2; } bool is_cell() const { return _curr_simplex.dimension() == 3; } + const Cell cell() const + { + return _cell_iterator.cell(); + } + const Simplex_type& get_simplex() const { return _curr_simplex; } Vertex_handle get_vertex() const {