]> https://gitweb.dealii.org/ - dealii.git/commitdiff
variadic template for CGAL Mesh criteria 13753/head
authorMarco Feder <marco.feder@sissa.it>
Tue, 17 May 2022 08:05:18 +0000 (10:05 +0200)
committerMarco Feder <marco.feder@sissa.it>
Sun, 22 May 2022 17:22:31 +0000 (19:22 +0200)
include/deal.II/cgal/utilities.h
tests/cgal/cgal_mesh_criteria.cc [new file with mode: 0644]
tests/cgal/cgal_mesh_criteria.output [new file with mode: 0644]
tests/cgal/cgal_surface_mesh_02.cc
tests/cgal/cgal_triangulation_06.cc

index 4c375362689c51e35e864f0ac6635471be9289f1..d0b0b638eecc33bb60d917b1950031c96ee6e592 100644 (file)
@@ -130,23 +130,54 @@ namespace CGALWrappers
 
   /**
    * Given a closed CGAL::Surface_mesh, this function fills the
-   * region bounded by the surface with tets, keeping them as coarse as
-   * possible.
+   * region bounded by the surface with tets.
+   *
+   * The optional variable arguments @p cgal_args can be used to pass additional
+   * arguments to the CGAL::Mesh_criteria_3 class (see
+   * https://doc.cgal.org/latest/Mesh_3/index.html) for more information.
+   *
+   * The arguments allow for fine control on the size, quality, and distribution
+   * of the cells of the final triangulation. CGAL uses named parameters for
+   * these arguments in dimension three, i.e., they must be specified with the
+   * syntax `CGAL::parameters::parameter_name=parameter_value`, irrespective of
+   * their order. Accepted parameters are:
+   *
+   * - CGAL::parameters::edge_size: a scalar field (resp. a constant) providing
+   * a space varying (resp. a uniform) upper bound for the lengths of curve
+   * edges. This parameter has to be set to a positive value when 1-dimensional
+   * features protection is used.
+   * - CGAL::parameters::facet_angle: a lower bound for the angles (in degrees)
+   * of the surface mesh facets.
+   * - CGAL::parameters::facet_size: a scalar field (resp. a constant)
+   * describing a space varying (resp. a uniform) upper-bound or for the radii
+   * of the surface Delaunay balls.
+   * - CGAL::parameters::facet_distance: a scalar field (resp. a constant)
+   * describing a space varying (resp. a uniform) upper bound for the distance
+   * between the facet circumcenter and the center of its surface Delaunay ball.
+   * - CGAL::parameters::facet_topology: the set of topological constraints
+   * which have to be verified by each surface facet. The default value is
+   * CGAL::FACET_VERTICES_ON_SURFACE. See CGAL::Mesh_facet_topology manual page
+   * to get all possible values.
+   * - CGAL::parameters::cell_radius_edge_ratio: an upper bound for the
+   * radius-edge ratio of the mesh tetrahedra.
+   * - CGAL::parameters::cell_size: a scalar field (resp. a constant) describing
+   * a space varying (resp. a uniform) upper-bound for the circumradii of the
+   * mesh tetrahedra.
+   *
+   * If no parameters are specified, all the values will be set to `ignored` by
+   * CGAL.
    *
-   * The number of the generated tetrahedrons depends on the refinement level of
-   * surface mesh. This function does not attempt to construct a fine
-   * triangulation, nor to smooth the final result. If you need finer control on
-   * the resulting triangulation, you should consider using directly
-   * CGAL::make_mesh_3().
    *
    * @param [in] surface_mesh The (closed) surface mesh bounding the volume that has to be filled.
    * @param [out] triangulation The output triangulation filled with tetrahedra.
+   * @param [in] cgal_args Additional parameters to pass to the CGAL::make_mesh_3 function.
    */
-  template <typename C3t3>
+  template <typename C3t3, typename... Args>
   void
-  cgal_surface_mesh_to_cgal_coarse_triangulation(
+  cgal_surface_mesh_to_cgal_triangulation(
     CGAL::Surface_mesh<typename C3t3::Point::Point> &surface_mesh,
-    C3t3 &                                           triangulation);
+    C3t3 &                                           triangulation,
+    Args... cgal_args);
 
   /**
    * Given two triangulated surface meshes that bound two volumes, execute a
@@ -243,11 +274,12 @@ namespace CGALWrappers
 
 
 
-  template <typename C3t3>
+  template <typename C3t3, typename... Args>
   void
-  cgal_surface_mesh_to_cgal_coarse_triangulation(
+  cgal_surface_mesh_to_cgal_triangulation(
     CGAL::Surface_mesh<typename C3t3::Point::Point> &surface_mesh,
-    C3t3 &                                           triangulation)
+    C3t3 &                                           triangulation,
+    Args... cgal_args)
   {
     using CGALPointType = typename C3t3::Point::Point;
     Assert(CGAL::is_closed(surface_mesh),
@@ -264,7 +296,7 @@ namespace CGALWrappers
     CGAL::Polygon_mesh_processing::triangulate_faces(surface_mesh);
     Mesh_domain domain(surface_mesh);
     domain.detect_features();
-    Mesh_criteria criteria;
+    Mesh_criteria criteria(cgal_args...);
     // Mesh generation
     triangulation = CGAL::make_mesh_3<C3t3>(domain,
                                             criteria,
diff --git a/tests/cgal/cgal_mesh_criteria.cc b/tests/cgal/cgal_mesh_criteria.cc
new file mode 100644 (file)
index 0000000..f62e432
--- /dev/null
@@ -0,0 +1,72 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2022 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, 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 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+// Check that different `cell_size` parameter gives different number of cells.
+
+#include <deal.II/base/config.h>
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/tria.h>
+
+#include <CGAL/IO/File_medit.h>
+#include <CGAL/IO/io.h>
+#include <deal.II/cgal/utilities.h>
+
+#include "../tests.h"
+
+using K         = CGAL::Exact_predicates_inexact_constructions_kernel;
+using CGALPoint = CGAL::Point_3<K>;
+using namespace CGALWrappers;
+using Mesh_domain =
+  CGAL::Polyhedral_mesh_domain_with_features_3<K,
+                                               CGAL::Surface_mesh<CGALPoint>>;
+using Tr = typename CGAL::
+  Mesh_triangulation_3<Mesh_domain, CGAL::Default, ConcurrencyTag>::type;
+using Mesh_criteria = CGAL::Mesh_criteria_3<Tr>;
+using C3t3          = CGAL::Mesh_complex_3_in_triangulation_3<Tr,
+                                                     Mesh_domain::Corner_index,
+                                                     Mesh_domain::Curve_index>;
+
+void
+test()
+{
+  CGAL::Surface_mesh<CGALPoint> sm;
+  C3t3                          tria;
+  std::ifstream                 input("input_grids/cube.off");
+  input >> sm;
+  cgal_surface_mesh_to_cgal_triangulation(sm,
+                                          tria,
+                                          CGAL::parameters::cell_size = 0.1);
+  const unsigned int n_intial_facets = tria.number_of_facets_in_complex();
+  tria.clear();
+  cgal_surface_mesh_to_cgal_triangulation(sm,
+                                          tria,
+                                          CGAL::parameters::cell_size = 0.5);
+
+  Assert(
+    n_intial_facets > tria.number_of_facets_in_complex(),
+    ExcMessage(
+      "The number of facets in the finer mesh must be greater than the number of facets in the coarse mesh."));
+  deallog << std::boolalpha
+          << (n_intial_facets > tria.number_of_facets_in_complex())
+          << std::endl;
+}
+
+int
+main()
+{
+  initlog();
+  test();
+}
diff --git a/tests/cgal/cgal_mesh_criteria.output b/tests/cgal/cgal_mesh_criteria.output
new file mode 100644 (file)
index 0000000..6064c35
--- /dev/null
@@ -0,0 +1,2 @@
+
+DEAL::true
index 9ed9935628006117f45ca4b7d9336598a849baba..404c383f06a804e5e26e331bfb1cdc9e16ac0952 100644 (file)
@@ -56,7 +56,7 @@ test()
     {
       std::ifstream input(name);
       input >> sm;
-      cgal_surface_mesh_to_cgal_coarse_triangulation(sm, tria);
+      cgal_surface_mesh_to_cgal_triangulation(sm, tria);
       {
         std::ofstream off_file_medit("coarse.mesh");
         tria.output_to_medit(off_file_medit, false);
index aad82fc0c3dcf95b9e054cde448b8b5c83626501..11a59d0583e819454b384b774db8e5237aeae924 100644 (file)
@@ -53,8 +53,8 @@ main()
   }
 
   C3t3 cgal_triangulation;
-  cgal_surface_mesh_to_cgal_coarse_triangulation(cgal_surface_mesh,
-                                                 cgal_triangulation);
+  cgal_surface_mesh_to_cgal_triangulation(cgal_surface_mesh,
+                                          cgal_triangulation);
 
   Triangulation<3> tria;
   cgal_triangulation_to_dealii_triangulation(cgal_triangulation, tria);

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.