From adb6bf3ad747623f2b899a61de59d83948752157 Mon Sep 17 00:00:00 2001 From: Marco Feder Date: Tue, 17 May 2022 10:05:18 +0200 Subject: [PATCH] variadic template for CGAL Mesh criteria --- include/deal.II/cgal/utilities.h | 60 +++++++++++++++++------ tests/cgal/cgal_mesh_criteria.cc | 72 ++++++++++++++++++++++++++++ tests/cgal/cgal_mesh_criteria.output | 2 + tests/cgal/cgal_surface_mesh_02.cc | 2 +- tests/cgal/cgal_triangulation_06.cc | 4 +- 5 files changed, 123 insertions(+), 17 deletions(-) create mode 100644 tests/cgal/cgal_mesh_criteria.cc create mode 100644 tests/cgal/cgal_mesh_criteria.output diff --git a/include/deal.II/cgal/utilities.h b/include/deal.II/cgal/utilities.h index 4c37536268..d0b0b638ee 100644 --- a/include/deal.II/cgal/utilities.h +++ b/include/deal.II/cgal/utilities.h @@ -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 + template void - cgal_surface_mesh_to_cgal_coarse_triangulation( + cgal_surface_mesh_to_cgal_triangulation( CGAL::Surface_mesh &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 + template void - cgal_surface_mesh_to_cgal_coarse_triangulation( + cgal_surface_mesh_to_cgal_triangulation( CGAL::Surface_mesh &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(domain, criteria, diff --git a/tests/cgal/cgal_mesh_criteria.cc b/tests/cgal/cgal_mesh_criteria.cc new file mode 100644 index 0000000000..f62e432fd8 --- /dev/null +++ b/tests/cgal/cgal_mesh_criteria.cc @@ -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 + +#include +#include + +#include +#include +#include + +#include "../tests.h" + +using K = CGAL::Exact_predicates_inexact_constructions_kernel; +using CGALPoint = CGAL::Point_3; +using namespace CGALWrappers; +using Mesh_domain = + CGAL::Polyhedral_mesh_domain_with_features_3>; +using Tr = typename CGAL:: + Mesh_triangulation_3::type; +using Mesh_criteria = CGAL::Mesh_criteria_3; +using C3t3 = CGAL::Mesh_complex_3_in_triangulation_3; + +void +test() +{ + CGAL::Surface_mesh 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 index 0000000000..6064c35fc2 --- /dev/null +++ b/tests/cgal/cgal_mesh_criteria.output @@ -0,0 +1,2 @@ + +DEAL::true diff --git a/tests/cgal/cgal_surface_mesh_02.cc b/tests/cgal/cgal_surface_mesh_02.cc index 9ed9935628..404c383f06 100644 --- a/tests/cgal/cgal_surface_mesh_02.cc +++ b/tests/cgal/cgal_surface_mesh_02.cc @@ -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); diff --git a/tests/cgal/cgal_triangulation_06.cc b/tests/cgal/cgal_triangulation_06.cc index aad82fc0c3..11a59d0583 100644 --- a/tests/cgal/cgal_triangulation_06.cc +++ b/tests/cgal/cgal_triangulation_06.cc @@ -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); -- 2.39.5