From e2d83467ff20a8fea555e7a0650dbb486e769117 Mon Sep 17 00:00:00 2001 From: Matthias Maier Date: Fri, 10 Jun 2022 21:49:47 -0500 Subject: [PATCH] GridGenerator: put CGAL implementation into own compilation unit --- source/grid/CMakeLists.txt | 18 +-- source/grid/grid_generator_cgal.cc | 155 ++++++++++++++++++++++++ source/grid/grid_generator_cgal.inst.in | 30 +++++ 3 files changed, 195 insertions(+), 8 deletions(-) create mode 100644 source/grid/grid_generator_cgal.cc create mode 100644 source/grid/grid_generator_cgal.inst.in diff --git a/source/grid/CMakeLists.txt b/source/grid/CMakeLists.txt index 88a236d06e..0d13509ef3 100644 --- a/source/grid/CMakeLists.txt +++ b/source/grid/CMakeLists.txt @@ -21,12 +21,12 @@ INCLUDE_DIRECTORIES(BEFORE ${CMAKE_CURRENT_BINARY_DIR}) # 100k characters produced by one of the CGAL headers. Unfortunately, # guarding the include by DEAL_II_DISABLE_EXTRA_DIAGNOSTICS is not enough # due to a longstanding bug in gcc. Thus, simply wet -# -Wno-misleading-indentation on the command line for the grid_generator.cc -# compilation unit. +# -Wno-misleading-indentation on the command line for the +# grid_generator_cgal.cc compilation unit. # IF(DEAL_II_WITH_CGAL) ENABLE_IF_SUPPORTED(_flag "-Wno-misleading-indentation") - SET_PROPERTY(SOURCE "grid_generator.cc" + SET_PROPERTY(SOURCE "grid_generator_cgal.cc" APPEND_STRING PROPERTY COMPILE_FLAGS " ${_flag}" ) ENDIF() @@ -49,15 +49,16 @@ SET(_unity_include_src ) SET(_separate_src + grid_generator.cc + grid_generator_cgal.cc grid_generator_from_name.cc + grid_generator_pipe_junction.cc grid_in.cc grid_out.cc - grid_generator.cc - grid_generator_pipe_junction.cc - grid_tools.cc grid_tools_cache.cc - grid_tools_nontemplates.cc + grid_tools.cc grid_tools_dof_handlers.cc + grid_tools_nontemplates.cc tria.cc ) @@ -72,8 +73,9 @@ SETUP_SOURCE_LIST("${_unity_include_src}" SET(_inst cell_id.inst.in - grid_generator.inst.in + grid_generator_cgal.inst.in grid_generator_from_name.inst.in + grid_generator.inst.in grid_generator_pipe_junction.inst.in grid_in.inst.in grid_out.inst.in diff --git a/source/grid/grid_generator_cgal.cc b/source/grid/grid_generator_cgal.cc new file mode 100644 index 0000000000..ed9e7481ae --- /dev/null +++ b/source/grid/grid_generator_cgal.cc @@ -0,0 +1,155 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 1999 - 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. +// +// --------------------------------------------------------------------- + +#include + +#include + +#ifdef DEAL_II_WITH_CGAL +// Functions needed by the CGAL mesh generation utilities are inside +# include +#endif + + +DEAL_II_NAMESPACE_OPEN + +// work around the problem that doxygen for some reason lists all template +// specializations in this file +#ifndef DOXYGEN + +namespace GridGenerator +{ + template + void + implicit_function(Triangulation &tria, + const Function<3> & dealii_implicit_function, + const CGALWrappers::AdditionalData &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; + using Tr = + CGAL::Mesh_triangulation_3::type; + using C3t3 = CGAL::Mesh_complex_3_in_triangulation_3; + using Mesh_criteria = CGAL::Mesh_criteria_3; + + + 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(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; + 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; + using Surface_mesh = CGAL::Surface_mesh; + + + 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 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 +# include "grid_generator_cgal.inst" + +#endif // DOXYGEN + +DEAL_II_NAMESPACE_CLOSE diff --git a/source/grid/grid_generator_cgal.inst.in b/source/grid/grid_generator_cgal.inst.in new file mode 100644 index 0000000000..0db59c879d --- /dev/null +++ b/source/grid/grid_generator_cgal.inst.in @@ -0,0 +1,30 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 1999 - 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. +// +// --------------------------------------------------------------------- + +for (deal_II_dimension : DIMENSIONS) + { + namespace GridGenerator + \{ +#if deal_II_dimension >= 2 + template void + implicit_function( + Triangulation &, + const Function<3> &, + const CGALWrappers::AdditionalData &, + const Point<3> &, + const double &); +#endif + \} + } -- 2.39.5