]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Fixing default values in AdditionalData
authorMarco Feder <marco.feder@sissa.it>
Thu, 9 Jun 2022 22:58:22 +0000 (00:58 +0200)
committerMarco Feder <marco.feder@sissa.it>
Fri, 10 Jun 2022 08:34:20 +0000 (10:34 +0200)
include/deal.II/cgal/additional_data.h
tests/cgal/implicit_triangulation_01.cc
tests/cgal/implicit_triangulation_02.cc

index 2fcdadfede34aaa249c50b0dd8e8e152f2405eba..8bc74614fe84d69455b98865ec0166fc2259ddc9 100644 (file)
 
 #ifdef DEAL_II_WITH_CGAL
 
-#  include <CGAL/boost/graph/Named_function_parameters.h>
+#  include <CGAL/Mesh_facet_topology.h>
+
 
 DEAL_II_NAMESPACE_OPEN
 namespace CGALWrappers
 {
+  enum class FacetTopology
+  {
+    /**
+     * Each vertex of the facet have to be on the surface, on a curve, or on a
+     * corner.
+     */
+    facet_vertices_on_surface = CGAL::FACET_VERTICES_ON_SURFACE,
+    /**
+     * The three vertices of a facet belonging to a surface patch s have to be
+     * on the same surface patch s, on a curve or on a corner.
+     */
+    facet_vertices_on_same_surface_patch =
+      CGAL::FACET_VERTICES_ON_SAME_SURFACE_PATCH,
+    /**
+     * The three vertices of a facet belonging to a surface patch s have to be
+     * on the same surface patch s, or on a curve incident to the surface patch
+     * s or on a corner incident to the surface patch s.
+     */
+    facet_vertices_on_same_surface_patch_with_adjacency_check =
+      CGAL::FACET_VERTICES_ON_SAME_SURFACE_PATCH_WITH_ADJACENCY_CHECK,
+
+
+  };
+
   /**
    * Struct that must be used to pass additional
    * arguments to the CGAL::Mesh_criteria_3 class (see
@@ -50,7 +75,7 @@ namespace CGALWrappers
    * 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
+   * `CGAL::FACET_VERTICES_ON_SURFACE`. See the enum @p FacetToplogy 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.
@@ -63,26 +88,46 @@ namespace CGALWrappers
   template <int dim>
   struct AdditionalData
   {
+    // a constant providing a uniform upper bound for the lengths ofcurve edges.
+    // This parameter has to be set to a positive value when1-dimensional
+    // features protection is used.
+    double edge_size;
+    //   lower bound for the angles (in degrees) of the surface mesh facets.
     double facet_angle;
+    // constant describing a uniform upper bound for the radii of the surface
+    // Delaunay balls.
     double facet_size;
+    //  a constant describing a uniform upper bound for the distance between the
+    //  facet circumcenter and the center of its surface Delaunay ball.
     double facet_distance;
+    // set of topological constraints which have to be verified by each surface
+    // facet.
+    FacetTopology facet_topology;
+    // upper bound for the radius-edge ratio of the mesh tetrahedra.
     double cell_radius_edge_ratio;
+    // a constant describing a uniform upper bound for the circumradii of the
+    // mesh tetrahedra.
     double cell_size;
 
     AdditionalData(
-      double facet_a            = CGAL::parameters::is_default_parameter(true),
-      double facet_s            = CGAL::parameters::is_default_parameter(true),
-      double facet_d            = CGAL::parameters::is_default_parameter(true),
-      double cell_radius_edge_r = CGAL::parameters::is_default_parameter(true),
-      double cell_s             = CGAL::parameters::is_default_parameter(true))
+      double        edge_s  = std::numeric_limits<double>::max(),
+      double        facet_a = 0.,
+      double        facet_s = 0.,
+      double        facet_d = 0.,
+      FacetTopology facet_t =
+        dealii::CGALWrappers::FacetTopology::facet_vertices_on_surface,
+      double cell_radius_edge_r = 0.,
+      double cell_s             = 0.)
     {
       AssertThrow(
         dim == 3,
         ExcMessage(
           "These struct can be instantiated with 3D Triangulations only."));
+      edge_size              = edge_s;
       facet_angle            = facet_a;
       facet_size             = facet_s;
       facet_distance         = facet_d;
+      facet_topology         = facet_t;
       cell_radius_edge_ratio = cell_radius_edge_r;
       cell_size              = cell_s;
     }
@@ -106,7 +151,15 @@ namespace CGALWrappers
   template <>
   struct AdditionalData<2>
   {
-    double angular_bound, radius_bound, distance_bound;
+    // lower bound in degrees for the angles of mesh facets.
+    double angular_bound;
+    // upper bound on the radii of surface Delaunay balls. A surface Delaunay
+    // ball is a ball circumscribing a mesh facet and centered on the surface.
+    double radius_bound;
+    // upper bound for the distance between the circumcenter of a mesh facet and
+    // the center of a surface Delaunay ball of this facet.
+    double distance_bound;
+
     AdditionalData(double angular_b  = 0.,
                    double radius_b   = 0.,
                    double distance_b = 0.)
index 054fc51ea20dc8076d6c8f541eb85b1eb1e6e0ed..a5ee98a785fc61a77dc41b589bdc564a2eb55102 100644 (file)
@@ -25,7 +25,6 @@
 
 #include "../tests.h"
 
-using namespace CGAL::parameters;
 
 class ImplicitFunction : public Function<3>
 {
@@ -54,7 +53,7 @@ main()
     std::ofstream of("tria.vtk");
     go.write_vtk(tria, of);
   }
-  // remove(/"tria.vtk");
+  remove("tria.vtk");
   //  If we got here, everything was ok, including writing the grid.
   deallog << "OK" << std::endl;
 }
index 5e8456b7ffc51d0b4ca7aadc8b97ebf2f30e22ae..474cb359b4ed1ff3a1c34285681bbac3ed5ef1ed 100644 (file)
@@ -25,7 +25,6 @@
 
 #include "../tests.h"
 
-using namespace CGAL::parameters;
 
 class ImplicitFunction : public Function<3>
 {
@@ -53,7 +52,7 @@ main()
     tria, implicit_function, data, Point<3>(1, 0, 0), 10.);
   {
     GridOut       go;
-    std::ofstream of("tria_ancora.vtk");
+    std::ofstream of("tria.vtk");
     go.write_vtk(tria, of);
   }
   remove("tria.vtk");

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.