Fix tests and benchmark

Isosurfacing_3/benchmark/Isosurfacing_3/benchmark.cpp

@@ -1,9 +1,9 @@
 #include <CGAL/Cartesian.h>
 #include <CGAL/Cartesian_grid_3.h>
-#include <CGAL/Cartesian_grid_domain.h>
+#include <CGAL/Default_gradients.h>
 #include <CGAL/Dual_contouring_3.h>
 #include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
-#include <CGAL/Implicit_domain.h>
+#include <CGAL/Isosurfacing_domains.h>
 #include <CGAL/Marching_cubes_3.h>
 #include <CGAL/Polygon_mesh_processing/polygon_soup_to_polygon_mesh.h>
 #include <CGAL/Simple_cartesian.h>
@@ -39,13 +39,16 @@ struct SphereGradient {
 template <class GeomTraits>
 struct Implicit_sphere {
 
-    typedef CGAL::Isosurfacing::Implicit_domain<GeomTraits, SphereValue<GeomTraits>, SphereGradient<GeomTraits>> Domain;
+    typedef CGAL::Isosurfacing::Implicit_cartesian_grid_domain_with_gradient<GeomTraits, SphereValue<GeomTraits>,
+                                                                             SphereGradient<GeomTraits>>
+        Domain;
     typedef typename GeomTraits::Vector_3 Vector;
 
     Implicit_sphere(const std::size_t N) : res(2.0 / N, 2.0 / N, 2.0 / N) {}
 
     Domain domain() const {
-        return Domain({-1, -1, -1, 1, 1, 1}, res, val, grad);
+        return CGAL::Isosurfacing::create_implicit_cartesian_grid_domain<GeomTraits>({-1, -1, -1, 1, 1, 1}, res, val,
+                                                                                     grad);
     }
 
     typename GeomTraits::FT iso() const {
@@ -94,13 +97,16 @@ struct IWPGradient {
 template <class GeomTraits>
 struct Implicit_iwp {
 
-    typedef CGAL::Isosurfacing::Implicit_domain<GeomTraits, IWPValue<GeomTraits>, IWPGradient<GeomTraits>> Domain;
+    typedef CGAL::Isosurfacing::Implicit_cartesian_grid_domain_with_gradient<GeomTraits, IWPValue<GeomTraits>,
+                                                                             IWPGradient<GeomTraits>>
+        Domain;
     typedef typename GeomTraits::Vector_3 Vector;
 
     Implicit_iwp(const std::size_t N) : res(2.0 / N, 2.0 / N, 2.0 / N) {}
 
     Domain domain() const {
-        return Domain({-1, -1, -1, 1, 1, 1}, res, val, grad);
+        return CGAL::Isosurfacing::create_implicit_cartesian_grid_domain<GeomTraits>({-1, -1, -1, 1, 1, 1}, res, val,
+                                                                                     grad);
     }
 
     typename GeomTraits::FT iso() const {
@@ -116,7 +122,8 @@ private:
 template <class GeomTraits>
 struct Grid_sphere {
 
-    typedef CGAL::Isosurfacing::Cartesian_grid_domain<GeomTraits> Domain;
+    typedef CGAL::Isosurfacing::Explicit_cartesian_grid_gradient<GeomTraits> Gradient;
+    typedef CGAL::Isosurfacing::Explicit_cartesian_grid_domain_with_gradient<GeomTraits, Gradient> Domain;
     typedef CGAL::Cartesian_grid_3<GeomTraits> Grid;
     typedef typename GeomTraits::FT FT;
     typedef typename GeomTraits::Point_3 Point;
@@ -143,7 +150,7 @@ struct Grid_sphere {
     }
 
     Domain domain() const {
-        return Domain(grid);
+        return CGAL::Isosurfacing::create_explicit_cartesian_grid_domain<GeomTraits>(grid, Gradient(grid));
     }
 
     typename GeomTraits::FT iso() const {
@@ -157,11 +164,12 @@ private:
 template <class GeomTraits>
 struct Skull_image {
 
-    typedef CGAL::Isosurfacing::Cartesian_grid_domain<GeomTraits> Domain;
+    typedef CGAL::Isosurfacing::Explicit_cartesian_grid_gradient<GeomTraits> Gradient;
+    typedef CGAL::Isosurfacing::Explicit_cartesian_grid_domain_with_gradient<GeomTraits, Gradient> Domain;
     typedef CGAL::Cartesian_grid_3<GeomTraits> Grid;
 
     Skull_image(const std::size_t N) : grid(2, 2, 2, {-1, -1, -1, 1, 1, 1}) {
-        const std::string fname = "../data/skull_2.9.inr";
+        const std::string fname = CGAL::data_file_path("images/skull_2.9.inr");
         CGAL::Image_3 image;
         if (!image.read(fname)) {
             std::cerr << "Error: Cannot read file " << fname << std::endl;
@@ -171,7 +179,7 @@ struct Skull_image {
     }
 
     Domain domain() const {
-        return Domain(grid);
+        return CGAL::Isosurfacing::create_explicit_cartesian_grid_domain<GeomTraits>(grid, Gradient(grid));
     }
 
     typename GeomTraits::FT iso() const {

Isosurfacing_3/include/CGAL/Isosurfacing_3/internal/Grid_topology.h

@@ -125,10 +125,13 @@ public:
 #ifdef CGAL_LINKED_WITH_TBB
     template <typename Functor>
     void iterate_vertices(Functor& f, Parallel_tag) const {
-        auto iterator = [&f, size_i, size_j, size_k](const tbb::blocked_range<std::size_t>& r) {
+        const std::size_t sj = size_j;
+        const std::size_t sk = size_k;
+
+        auto iterator = [&f, sj, sk](const tbb::blocked_range<std::size_t>& r) {
             for (std::size_t i = r.begin(); i != r.end(); i++) {
-                for (std::size_t j = 0; j < size_j; j++) {
+                for (std::size_t j = 0; j < sj; j++) {
-                    for (std::size_t k = 0; k < size_k; k++) {
+                    for (std::size_t k = 0; k < sk; k++) {
                         f({i, j, k});
                     }
                 }
@@ -140,10 +143,13 @@ public:
 
     template <typename Functor>
     void iterate_edges(Functor& f, Parallel_tag) const {
-        auto iterator = [&f, size_i, size_j, size_k](const tbb::blocked_range<std::size_t>& r) {
+        const std::size_t sj = size_j;
+        const std::size_t sk = size_k;
+
+        auto iterator = [&f, sj, sk](const tbb::blocked_range<std::size_t>& r) {
             for (std::size_t i = r.begin(); i != r.end(); i++) {
-                for (std::size_t j = 0; j < size_j - 1; j++) {
+                for (std::size_t j = 0; j < sj - 1; j++) {
-                    for (std::size_t k = 0; k < size_k - 1; k++) {
+                    for (std::size_t k = 0; k < sk - 1; k++) {
                         f({i, j, k, 0});
                         f({i, j, k, 1});
                         f({i, j, k, 2});
@@ -157,10 +163,13 @@ public:
 
     template <typename Functor>
     void iterate_cells(Functor& f, Parallel_tag) const {
-        auto iterator = [&f, size_i, size_j, size_k](const tbb::blocked_range<std::size_t>& r) {
+        const std::size_t sj = size_j;
+        const std::size_t sk = size_k;
+
+        auto iterator = [&f, sj, sk](const tbb::blocked_range<std::size_t>& r) {
             for (std::size_t i = r.begin(); i != r.end(); i++) {
-                for (std::size_t j = 0; j < size_j - 1; j++) {
+                for (std::size_t j = 0; j < sj - 1; j++) {
-                    for (std::size_t k = 0; k < size_k - 1; k++) {
+                    for (std::size_t k = 0; k < sk - 1; k++) {
                         f({i, j, k});
                     }
                 }

Isosurfacing_3/test/Isosurfacing_3/CMakeLists.txt

@@ -9,13 +9,9 @@ find_package(CGAL REQUIRED)
 create_single_source_cgal_program( "test_dual_contouring.cpp" )
 create_single_source_cgal_program( "test_marching_cubes.cpp" )
 
-find_package(OpenMP)
-
 find_package(TBB)
 include(CGAL_TBB_support)
 if(TARGET CGAL::TBB_support)
   target_link_libraries(test_dual_contouring PUBLIC CGAL::TBB_support)
-  target_link_libraries(test_dual_contouring PRIVATE OpenMP::OpenMP_CXX)
-
   target_link_libraries(test_marching_cubes PUBLIC CGAL::TBB_support)
 endif()

Isosurfacing_3/test/Isosurfacing_3/test_dual_contouring.cpp

@@ -1,29 +1,24 @@
 #include <CGAL/Cartesian_grid_3.h>
-#include <CGAL/Cartesian_grid_domain.h>
-#include <CGAL/Cartesian_grid_domain_old.h>
 #include <CGAL/Dual_contouring_3.h>
-#include <CGAL/Implicit_domain.h>
+#include <CGAL/Isosurfacing_domains.h>
-#include <CGAL/Implicit_domain_old.h>
 #include <CGAL/Marching_cubes_3.h>
 #include <CGAL/Octree_wrapper.h>
 #include <CGAL/Polygon_mesh_processing/polygon_soup_to_polygon_mesh.h>
 #include <CGAL/Simple_cartesian.h>
 #include <CGAL/Surface_mesh.h>
-#include <CGAL/TC_marching_cubes_3.h>
 #include <CGAL/boost/graph/IO/OFF.h>
-#include <tbb/concurrent_vector.h>
 
 #include "Timer.h"
 
-typedef CGAL::Simple_cartesian<float> Kernel;
+typedef CGAL::Simple_cartesian<double> Kernel;
 typedef typename Kernel::Vector_3 Vector;
 typedef typename Kernel::Point_3 Point;
 
 typedef CGAL::Surface_mesh<Point> Mesh;
 typedef CGAL::Cartesian_grid_3<Kernel> Grid;
 
-typedef tbb::concurrent_vector<Point> Point_range;
+typedef std::vector<Point> Point_range;
-typedef tbb::concurrent_vector<std::vector<std::size_t>> Polygon_range;
+typedef std::vector<std::vector<std::size_t>> Polygon_range;
 
 int main() {
     const Vector spacing(0.002, 0.002, 0.02);
@@ -35,8 +30,8 @@ int main() {
 
     typedef CGAL::Isosurfacing::Finite_difference_gradient<Kernel, decltype(sphere_function)> Gradient;
 
-    CGAL::Isosurfacing::Implicit_domain<Kernel, decltype(sphere_function), Gradient> implicit_domain(
+    auto implicit_domain = CGAL::Isosurfacing::create_implicit_cartesian_grid_domain<Kernel>(
-        {-1, -1, -1, 1, 1, 1}, spacing, sphere_function, Gradient(sphere_function, 0.0001));  // TODO: this is ugly
+        bbox, spacing, sphere_function, Gradient(sphere_function, 0.0001));
 
     const std::size_t nx = static_cast<std::size_t>(2.0 / spacing.x());
     const std::size_t ny = static_cast<std::size_t>(2.0 / spacing.y());
@@ -65,7 +60,7 @@ int main() {
     //}
     // Grid grid(image);
 
-    CGAL::Isosurfacing::Cartesian_grid_domain<Kernel> grid_domain(grid);
+    auto grid_domain = CGAL::Isosurfacing::create_explicit_cartesian_grid_domain<Kernel>(grid);
 
     Point_range points;
     Polygon_range polygons;

Isosurfacing_3/test/Isosurfacing_3/test_marching_cubes.cpp

@@ -22,8 +22,7 @@ void test_implicit_sphere() {
     const Vector spacing(0.2, 0.2, 0.2);
     const CGAL::Bbox_3 bbox = {-1, -1, -1, 1, 1, 1};
 
-    CGAL::Isosurfacing::Implicit_domain<Kernel, Sphere_function> domain(bbox, spacing,
+    auto domain = CGAL::Isosurfacing::create_implicit_cartesian_grid_domain<Kernel>(bbox, spacing, Sphere_function());
-                                                                        Sphere_function());  // TODO: this is ugly
 
     Point_range points;
     Polygon_range polygons;
@@ -64,7 +63,7 @@ void test_grid_sphere(const std::size_t n) {
         }
     }
 
-    CGAL::Isosurfacing::Cartesian_grid_domain<Kernel> domain(grid);
+    auto domain = CGAL::Isosurfacing::create_explicit_cartesian_grid_domain<Kernel>(grid);
 
     Point_range points;
     Polygon_range polygons;

Isosurfacing_3/test/Isosurfacing_3/test_util.h

@@ -2,16 +2,14 @@
 #define CGAL_ISOSURFACING_TEST_UTIL_H
 
 #include <CGAL/Cartesian_grid_3.h>
-#include <CGAL/Cartesian_grid_domain.h>
 #include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
-#include <CGAL/Implicit_domain.h>
+#include <CGAL/Isosurfacing_domains.h>
 #include <CGAL/Polygon_mesh_processing/distance.h>
 #include <CGAL/Polygon_mesh_processing/manifoldness.h>
 #include <CGAL/Polygon_mesh_processing/polygon_soup_to_polygon_mesh.h>
 #include <CGAL/Polygon_mesh_processing/repair.h>
 #include <CGAL/Polygon_mesh_processing/repair_polygon_soup.h>
 #include <CGAL/Surface_mesh.h>
-#include <tbb/concurrent_vector.h>
 
 #include <cassert>
 #include <iostream>