--- /dev/null
+// ---------------------------------------------------------------------
+//
+// 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.
+//
+// ---------------------------------------------------------------------
+
+#ifndef dealii_cgal_additional_data_h
+#define dealii_cgal_additional_data_h
+
+#include <deal.II/base/config.h>
+
+#ifdef DEAL_II_WITH_CGAL
+
+# include <CGAL/boost/graph/Named_function_parameters.h>
+
+DEAL_II_NAMESPACE_OPEN
+namespace CGALWrappers
+{
+ /**
+ * Struct that must 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 Boost 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 constant
+ * providing 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 constant
+ * describing a uniform upper bound for the radii
+ * of the surface Delaunay balls.
+ * - `CGAL::parameters::facet_distance`: a constant
+ * describing 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 constant
+ * describing a uniform upper bound for the
+ * circumradii of the mesh tetrahedra.
+ *
+ * @note This struct must be instantiated with `dim=3`.
+ */
+ template <int dim>
+ struct AdditionalData
+ {
+ double facet_angle;
+ double facet_size;
+ double facet_distance;
+ double cell_radius_edge_ratio;
+ 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))
+ {
+ AssertThrow(
+ dim == 3,
+ ExcMessage(
+ "These struct can be instantiated with 3D Triangulations only."));
+ facet_angle = facet_a;
+ facet_size = facet_s;
+ facet_distance = facet_d;
+ cell_radius_edge_ratio = cell_radius_edge_r;
+ cell_size = cell_s;
+ }
+ };
+
+
+
+ /**
+ * Specialization of the above struct when the object to be constructed is a
+ * 2D triangulation embedded in the 3D space, i.e. a Triangulation<2,3>.
+ * Only three parameters are accepted:
+ * - `angular_bound` is a lower bound in degrees for the angles of mesh
+ * facets.
+ * - `radius_bound` is an 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.
+ * - `distance_bound` is an upper bound for the distance between the
+ * circumcenter of a mesh facet and the center of a surface Delaunay ball of
+ * this facet.
+ */
+ template <>
+ struct AdditionalData<2>
+ {
+ double angular_bound, radius_bound, distance_bound;
+ AdditionalData(double angular_b = 0.,
+ double radius_b = 0.,
+ double distance_b = 0.)
+ {
+ angular_bound = angular_b;
+ radius_bound = radius_b;
+ distance_bound = distance_b;
+ }
+ };
+
+} // namespace CGALWrappers
+DEAL_II_NAMESPACE_CLOSE
+
+#else
+
+DEAL_II_NAMESPACE_OPEN
+namespace CGALWrappers
+{
+ template <int dim>
+ struct AdditionalData
+ {};
+} // namespace CGALWrappers
+
+DEAL_II_NAMESPACE_CLOSE
+#endif
+
+#endif
#include <deal.II/base/config.h>
-#ifdef DEAL_II_WITH_CGAL
+#include <deal.II/base/exceptions.h>
+#include <deal.II/base/function.h>
+
+#include <deal.II/grid/tria.h>
-# include <deal.II/grid/tria.h>
+#include <deal.II/cgal/utilities.h>
+
+#ifdef DEAL_II_WITH_CGAL
# include <boost/hana.hpp>
+# include <CGAL/Complex_2_in_triangulation_3.h>
+# include <CGAL/IO/facets_in_complex_2_to_triangle_mesh.h>
+# include <CGAL/Implicit_surface_3.h>
+# include <CGAL/Labeled_mesh_domain_3.h>
+# include <CGAL/Mesh_complex_3_in_triangulation_3.h>
+# include <CGAL/Mesh_criteria_3.h>
+# include <CGAL/Mesh_triangulation_3.h>
# include <CGAL/Polyhedron_3.h>
# include <CGAL/Surface_mesh.h>
+# include <CGAL/Surface_mesh_default_triangulation_3.h>
# include <CGAL/Triangulation_2.h>
# include <CGAL/Triangulation_3.h>
-# include <deal.II/cgal/utilities.h>
+# include <CGAL/make_mesh_3.h>
+# include <CGAL/make_surface_mesh.h>
+# include <deal.II/cgal/surface_mesh.h>
DEAL_II_NAMESPACE_OPEN
// A face is identified by a cell and by the index within the cell
// of the opposite vertex. Loop over vertices, and retain only those
// that belong to this face.
- unsigned int j = 0;
- for (unsigned int i = 0; i < 4; ++i)
+ int j = 0;
+ for (int i = 0; i < 4; ++i)
if (i != cgal_vertex_face_index)
dealii_face.vertices[j++] =
cgal_to_dealii_vertex_map[cgal_cell->vertex(i)];
}
triangulation.create_triangulation(vertices, cells, subcell_data);
}
-# endif
} // namespace CGALWrappers
+# endif // doxygen
DEAL_II_NAMESPACE_CLOSE
#include <deal.II/base/config.h>
+#include <deal.II/cgal/additional_data.h>
+
#ifdef DEAL_II_WITH_CGAL
# include <deal.II/base/quadrature_lib.h>
- /**
- * Struct that must 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 Boost 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 constant
- * providing 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 constant
- * describing a uniform upper bound for the radii
- * of the surface Delaunay balls.
- * - `CGAL::parameters::facet_distance`: a constant
- * describing 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 constant
- * describing a uniform upper bound for the
- * circumradii of the mesh tetrahedra.
- *
- * @note This struct must be instantiated with `dim=3`.
- */
- template <int dim>
- struct AdditionalData
- {
- double facet_angle;
- double facet_size;
- double facet_distance;
- double cell_radius_edge_ratio;
- 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))
- {
- AssertThrow(
- dim == 3,
- ExcMessage(
- "These struct can be instantiated with 3D Triangulations only."));
- facet_angle = facet_a;
- facet_size = facet_s;
- facet_distance = facet_d;
- cell_radius_edge_ratio = cell_radius_edge_r;
- cell_size = cell_s;
- }
- };
-
-
-
- /**
- * Specialization of the above struct when the object to be constructed is a
- * 2D triangulation embedded in the 3D space, i.e. a Triangulation<2,3>.
- * Only three parameters are accepted:
- * - `angular_bound` is a lower bound in degrees for the angles of mesh
- * facets.
- * - `radius_bound` is an 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.
- * - `distance_bound` is an upper bound for the distance between the
- * circumcenter of a mesh facet and the center of a surface Delaunay ball of
- * this facet.
- */
- template <>
- struct AdditionalData<2>
- {
- double angular_bound, radius_bound, distance_bound;
- AdditionalData(double angular_b = 0.,
- double radius_b = 0.,
- double distance_b = 0.)
- {
- angular_bound = angular_b;
- radius_bound = radius_b;
- distance_bound = distance_b;
- }
- };
-
-
-
/**
* Given a closed CGAL::Surface_mesh, this function fills the
* region bounded by the surface with tets.
#include <deal.II/grid/reference_cell.h>
#include <deal.II/grid/tria.h>
-#ifdef DEAL_II_WITH_CGAL
-// All functions needed by the CGAL mesh generation utilities
-# include <CGAL/Complex_2_in_triangulation_3.h>
-# include <CGAL/IO/facets_in_complex_2_to_triangle_mesh.h>
-# include <CGAL/Implicit_surface_3.h>
-# include <CGAL/Labeled_mesh_domain_3.h>
-# include <CGAL/Mesh_complex_3_in_triangulation_3.h>
-# include <CGAL/Mesh_criteria_3.h>
-# include <CGAL/Mesh_triangulation_3.h>
-# include <CGAL/Surface_mesh.h>
-# include <CGAL/Surface_mesh_default_triangulation_3.h>
-# include <CGAL/make_mesh_3.h>
-# include <CGAL/make_surface_mesh.h>
-# include <deal.II/cgal/triangulation.h>
-# include <deal.II/cgal/utilities.h>
-#endif
+#include <deal.II/cgal/additional_data.h>
#include <array>
#include <map>
const std::string & grid_generator_function_name,
const std::string & grid_generator_function_arguments);
-#ifdef DEAL_II_WITH_CGAL
/**
* Generate a Triangulation from the zero level set of an implicit function,
* using the CGAL library.
CGALWrappers::AdditionalData<dim>{},
const Point<3> &interior_point = Point<3>(),
const double & outer_ball_radius = 1.0);
-#endif
///@}
/**
const double,
const bool);
-# ifdef DEAL_II_WITH_CGAL
- template <int dim>
- void
- implicit_function(Triangulation<dim, 3> &tria,
- const Function<3> & dealii_implicit_function,
- const CGALWrappers::AdditionalData<dim> &data,
- const Point<3> & interior_point,
- const double & outer_ball_radius)
- {
- Assert(dealii_implicit_function.n_components == 1,
- ExcMessage(
- "The implicit function must have exactly one component."));
- Assert(dealii_implicit_function.value(interior_point) < 0,
- ExcMessage(
- "The implicit function must be negative at the interior point."));
- Assert(outer_ball_radius > 0,
- ExcMessage("The outer ball radius must be positive."));
- Assert(tria.n_active_cells() == 0,
- ExcMessage("The triangulation must be empty."));
-
- if constexpr (dim == 3)
- {
- using K = CGAL::Exact_predicates_inexact_constructions_kernel;
- using NumberType = K::FT;
- using Point_3 = K::Point_3;
- using Sphere_3 = K::Sphere_3;
-
- using Mesh_domain = CGAL::Labeled_mesh_domain_3<K>;
- using Tr =
- CGAL::Mesh_triangulation_3<Mesh_domain,
- CGAL::Default,
- CGALWrappers::ConcurrencyTag>::type;
- using C3t3 = CGAL::Mesh_complex_3_in_triangulation_3<Tr, int, int>;
- using Mesh_criteria = CGAL::Mesh_criteria_3<Tr>;
-
-
- auto cgal_implicit_function = [&](const Point_3 &p) {
- return NumberType(
- dealii_implicit_function.value(Point<3>(p.x(), p.y(), p.z())));
- };
-
- Mesh_domain domain = Mesh_domain::create_implicit_mesh_domain(
- cgal_implicit_function,
- K::Sphere_3(
- Point_3(interior_point[0], interior_point[1], interior_point[2]),
- outer_ball_radius * outer_ball_radius));
-
- Mesh_criteria criteria(CGAL::parameters::facet_size = data.facet_size,
- CGAL::parameters::facet_angle = data.facet_angle,
- CGAL::parameters::facet_distance =
- data.facet_distance,
- CGAL::parameters::cell_radius_edge_ratio =
- data.cell_radius_edge_ratio,
- CGAL::parameters::cell_size = data.cell_size);
-
- auto cgal_triangulation = CGAL::make_mesh_3<C3t3>(domain, criteria);
- CGALWrappers::cgal_triangulation_to_dealii_triangulation(
- cgal_triangulation, tria);
- }
- else if constexpr (dim == 2)
- {
- // default triangulation for Surface_mesher
- using Tr = CGAL::Surface_mesh_default_triangulation_3;
- using C2t3 = CGAL::Complex_2_in_triangulation_3<Tr>;
- using GT = Tr::Geom_traits;
- using Sphere_3 = GT::Sphere_3;
- using Point_3 = GT::Point_3;
- using FT = GT::FT;
- typedef FT (*Function)(Point_3);
- using Surface_3 = CGAL::Implicit_surface_3<GT, Function>;
- using Surface_mesh = CGAL::Surface_mesh<Point_3>;
-
-
- auto cgal_implicit_function = [&](const Point_3 &p) {
- return FT(
- dealii_implicit_function.value(Point<3>(p.x(), p.y(), p.z())));
- };
-
- Surface_3 surface(cgal_implicit_function,
- Sphere_3(Point_3(interior_point[0],
- interior_point[1],
- interior_point[2]),
- outer_ball_radius * outer_ball_radius));
-
- Tr tr;
- C2t3 c2t3(tr);
- Surface_mesh mesh;
-
- CGAL::Surface_mesh_default_criteria_3<Tr> criteria(data.angular_bound,
- data.radius_bound,
- data.distance_bound);
- CGAL::make_surface_mesh(c2t3,
- surface,
- criteria,
- CGAL::Non_manifold_tag());
- CGAL::facets_in_complex_2_to_triangle_mesh(c2t3, mesh);
- CGALWrappers::cgal_surface_mesh_to_dealii_triangulation(mesh, tria);
- }
- else
- {
- Assert(false, ExcImpossibleInDim(dim));
- }
- }
-# endif
+
#endif
} // namespace GridGenerator
#include <deal.II/physics/transformations.h>
+#ifdef DEAL_II_WITH_CGAL
+// Functions needed by the CGAL mesh generation utilities are inside
+# include <deal.II/cgal/triangulation.h>
+#endif
+
#include <array>
#include <cmath>
#include <limits>
}
}
+
+
+ template <int dim>
+ void
+ implicit_function(Triangulation<dim, 3> &tria,
+ const Function<3> & dealii_implicit_function,
+ const CGALWrappers::AdditionalData<dim> &data,
+ const Point<3> & interior_point,
+ const double & outer_ball_radius)
+ {
+# ifdef DEAL_II_WITH_CGAL
+ Assert(dealii_implicit_function.n_components == 1,
+ ExcMessage(
+ "The implicit function must have exactly one component."));
+ Assert(dealii_implicit_function.value(interior_point) < 0,
+ ExcMessage(
+ "The implicit function must be negative at the interior point."));
+ Assert(outer_ball_radius > 0,
+ ExcMessage("The outer ball radius must be positive."));
+ Assert(tria.n_active_cells() == 0,
+ ExcMessage("The triangulation must be empty."));
+
+ if constexpr (dim == 3)
+ {
+ using K = CGAL::Exact_predicates_inexact_constructions_kernel;
+ using NumberType = K::FT;
+ using Point_3 = K::Point_3;
+ using Sphere_3 = K::Sphere_3;
+
+ using Mesh_domain = CGAL::Labeled_mesh_domain_3<K>;
+ using Tr =
+ CGAL::Mesh_triangulation_3<Mesh_domain,
+ CGAL::Default,
+ CGALWrappers::ConcurrencyTag>::type;
+ using C3t3 = CGAL::Mesh_complex_3_in_triangulation_3<Tr, int, int>;
+ using Mesh_criteria = CGAL::Mesh_criteria_3<Tr>;
+
+
+ auto cgal_implicit_function = [&](const Point_3 &p) {
+ return NumberType(
+ dealii_implicit_function.value(Point<3>(p.x(), p.y(), p.z())));
+ };
+
+ Mesh_domain domain = Mesh_domain::create_implicit_mesh_domain(
+ cgal_implicit_function,
+ K::Sphere_3(
+ Point_3(interior_point[0], interior_point[1], interior_point[2]),
+ outer_ball_radius * outer_ball_radius));
+
+ Mesh_criteria criteria(CGAL::parameters::facet_size = data.facet_size,
+ CGAL::parameters::facet_angle = data.facet_angle,
+ CGAL::parameters::facet_distance =
+ data.facet_distance,
+ CGAL::parameters::cell_radius_edge_ratio =
+ data.cell_radius_edge_ratio,
+ CGAL::parameters::cell_size = data.cell_size);
+
+ auto cgal_triangulation = CGAL::make_mesh_3<C3t3>(domain, criteria);
+ CGALWrappers::cgal_triangulation_to_dealii_triangulation(
+ cgal_triangulation, tria);
+ }
+ else if constexpr (dim == 2)
+ {
+ // default triangulation for Surface_mesher
+ using Tr = CGAL::Surface_mesh_default_triangulation_3;
+ using C2t3 = CGAL::Complex_2_in_triangulation_3<Tr>;
+ using GT = Tr::Geom_traits;
+ using Sphere_3 = GT::Sphere_3;
+ using Point_3 = GT::Point_3;
+ using FT = GT::FT;
+ typedef FT (*Function)(Point_3);
+ using Surface_3 = CGAL::Implicit_surface_3<GT, Function>;
+ using Surface_mesh = CGAL::Surface_mesh<Point_3>;
+
+
+ auto cgal_implicit_function = [&](const Point_3 &p) {
+ return FT(
+ dealii_implicit_function.value(Point<3>(p.x(), p.y(), p.z())));
+ };
+
+ Surface_3 surface(cgal_implicit_function,
+ Sphere_3(Point_3(interior_point[0],
+ interior_point[1],
+ interior_point[2]),
+ outer_ball_radius * outer_ball_radius));
+
+ Tr tr;
+ C2t3 c2t3(tr);
+ Surface_mesh mesh;
+
+ CGAL::Surface_mesh_default_criteria_3<Tr> criteria(data.angular_bound,
+ data.radius_bound,
+ data.distance_bound);
+ CGAL::make_surface_mesh(c2t3,
+ surface,
+ criteria,
+ CGAL::Non_manifold_tag());
+ CGAL::facets_in_complex_2_to_triangle_mesh(c2t3, mesh);
+ CGALWrappers::cgal_surface_mesh_to_dealii_triangulation(mesh, tria);
+ }
+ else
+ {
+ Assert(false, ExcImpossibleInDim(dim));
+ }
+
+# else
+
+ (void)tria;
+ (void)dealii_implicit_function;
+ (void)data;
+ (void)interior_point;
+ (void)outer_ball_radius;
+ AssertThrow(false, ExcMessage("This function needs CGAL to be installed."));
+
+# endif
+ }
} // namespace GridGenerator
// explicit instantiations
hyper_ball_balanced<deal_II_dimension>(Triangulation<deal_II_dimension> &,
const Point<deal_II_dimension> &,
const double);
+
+ template void
+ implicit_function<deal_II_dimension>(
+ Triangulation<deal_II_dimension, 3> &,
+ const Function<3> &,
+ const CGALWrappers::AdditionalData<deal_II_dimension> &,
+ const Point<3> &,
+ const double &);
#endif
\}
}
projects / dealii.git / commitdiff