compute_union = 1 << 3,
};
+
+
+ /**
+ * 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;
+ }
+ };
+
+
+
/**
* Convert from a deal.II Point to any compatible CGAL point.
*
#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
* using the CGAL library.
*
* This function is only implemented for `dim` equal to two or three, and
- * requires that deal.II is configured using `DEAL_II_WITH_CGAL`. When
- * `dim` is equal to three, the @p implicit_function is supposed to be
- * negative in the interior of the domain, and positive outside, and the
- * triangulation that is generated covers the volume bounded by the zero level
- * set of the implicit function and the boundary of a ball of radius,
- * @p outer_ball_radius where the @p implicit_function is negative.
+ * requires that deal.II is configured using `DEAL_II_WITH_CGAL`. When `dim`
+ * is equal to three, the @p implicit_function is supposed to be negative in
+ * the interior of the domain, positive outside, and to be entirely enclosed
+ * in a ball of radius @p outer_ball_radius centered at the point
+ * @p interior_point. The triangulation that is generated covers the volume
+ * bounded by the zero level set of the implicit function where the
+ * @p implicit_function is negative.
*
* When `dim` is equal to two, the generated surface triangulation is the zero
- * level set of the @p implicit_function, with or without boundary, depending
- * on whether the zero level set intersects the boundary of the outer ball.
+ * level set of the @p implicit_function, oriented such that the surface
+ * triangulation has normals pointing towards the region where
+ * @p implicit_function is positive.
*
- * The orientation is such that the surface triangulation has normals pointing
- * towards the region where @p implicit_function is positive.
- *
- * The optional variable arguments @p cgal_args can be used to pass additional
+ * The struct @p data 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 fine control on the size, quality, and distribution of
- * the cells of the final triangulation. CGAL uses named parameters for these
- * arguments, i.e., they must be specified with the syntax
- * `CGAL::parameters::parameter_name=parameter_value`. Accepted parameters
- * are:
- *
- * - CGAL::parameters::facet_angl,
- * - CGAL::parameters::facet_size,
- * - CGAL::parameters::facet_distance,
- * - CGAL::parameters::cell_radius_edge_ratio,
- * - CGAL::parameters::cell_size;
+ * https://doc.cgal.org/latest/Mesh_3/index.html for more information.)
*
* An example usage of this function is given by
*
* @code
- * Triangulation<3> tria;
- * FunctionParser<3> implicit_function("(1-sqrt(x^2+y^2))^2+z^2-.25");
- * GridGenerator::implicit_function(
- * tria,
- * implicit_function,
- * Point<3>(1, 0, 0),
- * 10.0,
- * CGAL::parameters::cell_size = 0.2);
+ * Triangulation<dim, 3> tria;
+ * FunctionParser<3> my_function("(1-sqrt(x^2+y^2))^2+z^2-.25");
+ * GridGenerator::implicit_function( tria, my_function,
+ * Point<3>(.5, 0, 0), 1.0, cell_size = 0.2);
* @endcode
*
- * The above snippet of code generates the following image:
+ * The above snippet of code generates the following grid for `dim` equal to
+ * two and three respectively
+ *
+ * @image html grid_generator_implicit_function_2d.png
*
* @image html grid_generator_implicit_function_3d.png
*
+ * @ingroup simplex
*
* @param[out] tria The output triangulation
* @param[in] implicit_function The implicit function
+ * @param[in] data Additional parameters to pass to the CGAL::make_mesh_3
+ * function and to the CGAL::make_surface_mesh functions
* @param[in] interior_point A point in the interior of the domain, for which
* @p implicit_function is negative
* @param[in] outer_ball_radius The radius of the ball that will contain the
* generated Triangulation object
- * @param[in] cgal_args Additional parameters to pass to the CGAL::make_mesh_3
- * function and to the CGAL::make_surface_mesh functions
*/
- template <int dim, typename... Args>
+ template <int dim>
void
- implicit_function(Triangulation<dim, 3> &tria,
- const Function<3> & implicit_function,
- const Point<3> & interior_point = Point<3>(),
- const double & outer_ball_radius = 1.0,
- Args... cgal_args);
+ implicit_function(Triangulation<dim, 3> & tria,
+ const Function<3> & implicit_function,
+ const CGALWrappers::AdditionalData<dim> &data =
+ CGALWrappers::AdditionalData<dim>{},
+ const Point<3> &interior_point = Point<3>(),
+ const double & outer_ball_radius = 1.0);
#endif
///@}
const bool);
# ifdef DEAL_II_WITH_CGAL
- template <int dim, typename... Args>
+ template <int dim>
void
implicit_function(Triangulation<dim, 3> &tria,
- const Function<3> & implicit_function,
- const Point<3> & interior_point,
- const double & outer_ball_radius,
- Args... cgal_args)
+ const Function<3> & dealii_implicit_function,
+ const CGALWrappers::AdditionalData<dim> &data,
+ const Point<3> & interior_point,
+ const double & outer_ball_radius)
{
- Assert(implicit_function.n_components == 1,
+ Assert(dealii_implicit_function.n_components == 1,
ExcMessage(
"The implicit function must have exactly one component."));
- Assert(implicit_function.value(interior_point) < 0,
+ Assert(dealii_implicit_function.value(interior_point) < 0,
ExcMessage(
"The implicit function must be negative at the interior point."));
Assert(outer_ball_radius > 0,
Assert(tria.n_active_cells() == 0,
ExcMessage("The triangulation must be empty."));
- 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>;
-
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,
auto cgal_implicit_function = [&](const Point_3 &p) {
return NumberType(
- implicit_function.value(Point<3>(p.x(), p.y(), p.z())));
+ dealii_implicit_function.value(Point<3>(p.x(), p.y(), p.z())));
};
Mesh_domain domain = Mesh_domain::create_implicit_mesh_domain(
Point_3(interior_point[0], interior_point[1], interior_point[2]),
outer_ball_radius * outer_ball_radius));
- Mesh_criteria criteria(cgal_args...);
+ 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));