From: Marco Feder <marco.feder@sissa.it>
Date: Mon, 9 May 2022 22:04:34 +0000 (+0200)
Subject: Convert CGAL::Surface_mesh to coarse mesh
X-Git-Tag: v9.4.0-rc1~219^2
X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=fce63d5ca6090a62874133a75b7ebf23d2149956;p=dealii.git

Convert CGAL::Surface_mesh to coarse mesh
---

diff --git a/cmake/config/template-arguments.in b/cmake/config/template-arguments.in
index 4a9cc4aab1..a4231884ba 100644
--- a/cmake/config/template-arguments.in
+++ b/cmake/config/template-arguments.in
@@ -279,6 +279,5 @@ OUTPUT_FLAG_TYPES := { DXFlags; UcdFlags; GnuplotFlags; PovrayFlags; EpsFlags;
                        Deal_II_IntermediateFlags }
 
 // CGAL Kernels
-CGAL_KERNELS := {CGAL::Simple_cartesian<double>; CGAL::Exact_predicates_exact_constructions_kernel;
-                CGAL::Exact_predicates_exact_constructions_kernel_with_sqrt; 
+CGAL_KERNELS := {CGAL::Simple_cartesian<double>; CGAL::Exact_predicates_exact_constructions_kernel; 
                 CGAL::Exact_predicates_inexact_constructions_kernel }
diff --git a/include/deal.II/cgal/utilities.h b/include/deal.II/cgal/utilities.h
index 930ed2247a..b26cb32454 100644
--- a/include/deal.II/cgal/utilities.h
+++ b/include/deal.II/cgal/utilities.h
@@ -20,10 +20,23 @@
 
 #ifdef DEAL_II_WITH_CGAL
 #  include <CGAL/Cartesian.h>
+#  include <CGAL/Complex_2_in_triangulation_3.h>
 #  include <CGAL/Exact_predicates_exact_constructions_kernel.h>
-#  include <CGAL/Exact_predicates_exact_constructions_kernel_with_sqrt.h>
 #  include <CGAL/Exact_predicates_inexact_constructions_kernel.h>
+#  include <CGAL/Kernel_traits.h>
+#  include <CGAL/Mesh_complex_3_in_triangulation_3.h>
+#  include <CGAL/Mesh_criteria_3.h>
+#  include <CGAL/Mesh_triangulation_3.h>
+#  include <CGAL/Polygon_mesh_processing/triangulate_faces.h>
+#  include <CGAL/Polyhedral_mesh_domain_with_features_3.h>
 #  include <CGAL/Simple_cartesian.h>
+#  include <CGAL/Surface_mesh.h>
+#  include <CGAL/Triangulation_3.h>
+#  include <CGAL/make_mesh_3.h>
+#  include <CGAL/make_surface_mesh.h>
+
+#  include <type_traits>
+
 
 
 DEAL_II_NAMESPACE_OPEN
@@ -41,6 +54,12 @@ DEAL_II_NAMESPACE_OPEN
  */
 namespace CGALWrappers
 {
+#  ifdef CGAL_CONCURRENT_MESH_3
+  using ConcurrencyTag = CGAL::Parallel_tag;
+#  else
+  using ConcurrencyTag = CGAL::Sequential_tag;
+#  endif
+
   /**
    * Convert from deal.II Point to any compatible CGAL point.
    *
@@ -64,6 +83,26 @@ namespace CGALWrappers
   template <int dim, typename CGALPointType>
   inline dealii::Point<dim>
   cgal_point_to_dealii_point(const CGALPointType &p);
+
+  /**
+   * Given a closed CGAL::Surface_mesh, this function fills the
+   * region bounded by the surface with tets, keeping them as coarse as
+   * possible.
+   *
+   * 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.
+   */
+  template <typename C3t3>
+  void
+  cgal_surface_mesh_to_cgal_coarse_triangulation(
+    CGAL::Surface_mesh<typename C3t3::Point::Point> &surface_mesh,
+    C3t3 &                                           triangulation);
 } // namespace CGALWrappers
 
 #  ifndef DOXYGEN
@@ -106,6 +145,37 @@ namespace CGALWrappers
     else
       Assert(false, dealii::ExcNotImplemented());
   }
+
+
+
+  template <typename C3t3>
+  void
+  cgal_surface_mesh_to_cgal_coarse_triangulation(
+    CGAL::Surface_mesh<typename C3t3::Point::Point> &surface_mesh,
+    C3t3 &                                           triangulation)
+  {
+    using CGALPointType = typename C3t3::Point::Point;
+    Assert(CGAL::is_closed(surface_mesh),
+           ExcMessage("The surface mesh must be closed."));
+
+    using K           = typename CGAL::Kernel_traits<CGALPointType>::Kernel;
+    using Mesh_domain = CGAL::Polyhedral_mesh_domain_with_features_3<
+      K,
+      CGAL::Surface_mesh<CGALPointType>>;
+    using Tr = typename CGAL::
+      Mesh_triangulation_3<Mesh_domain, CGAL::Default, ConcurrencyTag>::type;
+    using Mesh_criteria = CGAL::Mesh_criteria_3<Tr>;
+
+    CGAL::Polygon_mesh_processing::triangulate_faces(surface_mesh);
+    Mesh_domain domain(surface_mesh);
+    domain.detect_features();
+    Mesh_criteria criteria;
+    // Mesh generation
+    triangulation = CGAL::make_mesh_3<C3t3>(domain,
+                                            criteria,
+                                            CGAL::parameters::no_perturb(),
+                                            CGAL::parameters::no_exude());
+  }
 } // namespace CGALWrappers
 #  endif
 
diff --git a/tests/cgal/cgal_surface_mesh_02.cc b/tests/cgal/cgal_surface_mesh_02.cc
new file mode 100644
index 0000000000..f9a1d228d1
--- /dev/null
+++ b/tests/cgal/cgal_surface_mesh_02.cc
@@ -0,0 +1,74 @@
+// ---------------------------------------------------------------------
+//
+// 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.
+//
+// ---------------------------------------------------------------------
+
+// Read a Surface_mesh from an .off file, then create a coarse mesh filled with
+// tets
+
+#include <deal.II/base/config.h>
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/tria.h>
+
+#include <CGAL/IO/File_medit.h>
+#include <CGAL/IO/io.h>
+#include <deal.II/cgal/utilities.h>
+
+#include "../tests.h"
+
+// Create a Surface_mesh from an .off file, then fill it with tets and print
+// vertices.
+
+using K         = CGAL::Exact_predicates_inexact_constructions_kernel;
+using CGALPoint = CGAL::Point_3<K>;
+using namespace CGALWrappers;
+using Mesh_domain =
+  CGAL::Polyhedral_mesh_domain_with_features_3<K,
+                                               CGAL::Surface_mesh<CGALPoint>>;
+using Tr = typename CGAL::
+  Mesh_triangulation_3<Mesh_domain, CGAL::Default, ConcurrencyTag>::type;
+using Mesh_criteria = CGAL::Mesh_criteria_3<Tr>;
+using C3t3          = CGAL::Mesh_complex_3_in_triangulation_3<Tr,
+                                                     Mesh_domain::Corner_index,
+                                                     Mesh_domain::Curve_index>;
+
+void
+test()
+{
+  const std::vector<std::string> fnames{"input_grids/cube.off",
+                                        "input_grids/tetrahedron.off"};
+  CGAL::Surface_mesh<CGALPoint>  sm;
+  C3t3                           tria;
+  for (const auto &name : fnames)
+    {
+      std::ifstream input(name);
+      input >> sm;
+      cgal_surface_mesh_to_cgal_coarse_triangulation(sm, tria);
+      {
+        std::ofstream off_file_medit("coarse.mesh");
+        tria.output_to_medit(off_file_medit, false);
+      }
+      cat_file("coarse.mesh");
+      sm.clear(); // reset surface
+      tria.clear();
+      deallog << std::endl << std::endl;
+    }
+}
+
+int
+main()
+{
+  initlog();
+  test();
+}
diff --git a/tests/cgal/cgal_surface_mesh_02.output b/tests/cgal/cgal_surface_mesh_02.output
new file mode 100644
index 0000000000..8330e7b04e
--- /dev/null
+++ b/tests/cgal/cgal_surface_mesh_02.output
@@ -0,0 +1,407 @@
+
+MeshVersionFormatted 1
+Dimension 3
+Vertices
+32
+-1 -1 -1 -1
+-1 -1 1 -1
+-1 1 -1 -1
+-1 1 1 -1
+1 -1 -1 -1
+1 -1 1 -1
+1 1 -1 -1
+1 1 1 -1
+-1 -0.33333333333333337 -1 -1
+-1 0.33333333333333326 -1 -1
+-0.33333333333333337 -1 -1 -1
+0.33333333333333326 -1 -1 -1
+-1 -1 -0.33333333333333337 -1
+-1 -1 0.33333333333333326 -1
+-0.33333333333333337 1 -1 -1
+0.33333333333333326 1 -1 -1
+-1 1 -0.33333333333333337 -1
+-1 1 0.33333333333333326 -1
+1 0.33333333333333337 -1 -1
+1 -0.33333333333333326 -1 -1
+1 1 -0.33333333333333337 -1
+1 1 0.33333333333333326 -1
+1 -1 -0.33333333333333337 -1
+1 -1 0.33333333333333326 -1
+-1 -0.33333333333333337 1 -1
+-1 0.33333333333333326 1 -1
+-0.33333333333333337 -1 1 -1
+0.33333333333333326 -1 1 -1
+-0.33333333333333337 1 1 -1
+0.33333333333333326 1 1 -1
+1 0.33333333333333337 1 -1
+1 -0.33333333333333326 1 -1
+Triangles
+120
+30 22 8 1
+8 22 30 1
+25 14 26 1
+26 14 25 1
+12 20 16 1
+16 20 12 1
+12 5 20 1
+20 5 12 1
+31 8 22 1
+22 8 31 1
+31 30 8 1
+8 30 31 1
+1 13 11 1
+11 13 1 1
+10 11 12 1
+12 11 10 1
+5 23 20 1
+20 23 5 1
+13 1 9 1
+9 1 13 1
+26 30 32 1
+32 30 26 1
+25 26 28 1
+28 26 25 1
+16 19 7 1
+7 19 16 1
+10 12 16 1
+16 12 10 1
+21 22 16 1
+16 22 21 1
+3 10 15 1
+15 10 3 1
+24 32 20 1
+20 32 24 1
+14 12 11 1
+11 12 14 1
+14 9 10 1
+10 9 14 1
+31 32 30 1
+30 32 31 1
+22 20 32 1
+32 20 22 1
+10 3 17 1
+17 3 10 1
+21 7 19 1
+19 7 21 1
+21 16 7 1
+7 16 21 1
+15 10 16 1
+16 10 15 1
+10 17 18 1
+18 17 10 1
+24 28 6 1
+6 28 24 1
+14 28 12 1
+12 28 14 1
+22 30 16 1
+16 30 22 1
+10 9 11 1
+11 9 10 1
+2 14 25 1
+25 14 2 1
+6 32 24 1
+24 32 6 1
+32 6 28 1
+28 6 32 1
+25 27 2 1
+2 27 25 1
+29 18 30 1
+30 18 29 1
+28 26 32 1
+32 26 28 1
+12 28 24 1
+24 28 12 1
+20 22 21 1
+21 22 20 1
+18 16 30 1
+30 16 18 1
+26 10 18 1
+18 10 26 1
+32 31 22 1
+22 31 32 1
+15 17 3 1
+3 17 15 1
+24 23 12 1
+12 23 24 1
+14 2 27 1
+27 2 14 1
+29 4 18 1
+18 4 29 1
+26 18 4 1
+4 18 26 1
+5 12 23 1
+23 12 5 1
+23 24 20 1
+20 24 23 1
+30 26 29 1
+29 26 30 1
+11 13 14 1
+14 13 11 1
+14 13 9 1
+9 13 14 1
+28 14 27 1
+27 14 28 1
+1 11 9 1
+9 11 1 1
+29 26 4 1
+4 26 29 1
+15 16 18 1
+18 16 15 1
+20 19 16 1
+16 19 20 1
+28 27 25 1
+25 27 28 1
+14 10 26 1
+26 10 14 1
+19 20 21 1
+21 20 19 1
+18 17 15 1
+15 17 18 1
+Tetrahedra
+36
+30 22 8 31 1
+26 18 30 32 1
+24 16 20 12 1
+24 10 16 12 1
+32 16 20 24 1
+16 32 20 22 1
+16 19 7 21 1
+10 32 28 24 1
+32 18 10 28 1
+32 10 16 24 1
+32 16 30 22 1
+32 6 28 24 1
+20 16 22 21 1
+18 16 30 32 1
+10 3 17 15 1
+28 26 14 25 1
+23 24 20 12 1
+32 30 31 22 1
+9 11 14 13 1
+11 1 13 9 1
+15 10 16 18 1
+26 18 4 29 1
+32 26 18 28 1
+10 18 26 28 1
+28 10 24 12 1
+18 10 16 32 1
+23 5 12 20 1
+12 14 10 11 1
+14 2 27 25 1
+10 26 14 28 1
+18 29 26 30 1
+28 14 27 25 1
+16 20 19 21 1
+14 10 28 12 1
+11 14 10 9 1
+17 10 15 18 1
+End
+
+DEAL::
+DEAL::
+MeshVersionFormatted 1
+Dimension 3
+Vertices
+34
+0 0 0 -1
+0 0 1 -1
+0 1 0 -1
+1 0 0 -1
+0.25 0.75 0 -1
+0.33333333333333337 0.66666666666666663 0 -1
+0.75 0.25 0 -1
+0 0.66666666666666674 0 -1
+0 0.33333333333333337 0 -1
+0 0.75 0.25 -1
+0 0.33333333333333348 0.66666666666666652 -1
+0 0.25 0.75 -1
+0.66666666666666674 0 0 -1
+0.33333333333333337 0 0 -1
+0.75 0 0.25 -1
+0.33333333333333348 0 0.66666666666666652 -1
+0.25 0 0.75 -1
+0 0 0.33333333333333331 -1
+0 0 0.66666666666666663 -1
+0.41666666666666685 0 0.58333333333333315 -1
+0.50000000000000022 0 0.49999999999999983 -1
+0.58333333333333348 0 0.41666666666666652 -1
+0.66666666666666685 0 0.33333333333333315 -1
+0 0.41666666666666685 0.58333333333333315 -1
+0 0.50000000000000022 0.49999999999999983 -1
+0 0.58333333333333348 0.41666666666666652 -1
+0 0.66666666666666685 0.33333333333333315 -1
+0.41666666666666663 0.58333333333333337 0 -1
+0.50000000000000011 0.49999999999999989 0 -1
+0.58333333333333337 0.41666666666666663 0 -1
+0.66666666666666674 0.33333333333333326 0 -1
+0.54454601052364293 0 0 -1
+0 0 0.54454601052364282 -1
+0 0.54454601052364293 0 -1
+Triangles
+128
+10 27 34 1
+34 27 10 1
+34 8 10 1
+10 8 34 1
+25 26 28 1
+28 26 25 1
+26 25 9 1
+9 25 26 1
+22 30 31 1
+31 30 22 1
+21 30 22 1
+22 30 21 1
+3 8 5 1
+5 8 3 1
+4 15 7 1
+7 15 4 1
+1 18 14 1
+14 18 1 1
+18 1 9 1
+9 1 18 1
+15 4 13 1
+13 4 15 1
+7 13 4 1
+4 13 7 1
+32 15 13 1
+13 15 32 1
+1 14 9 1
+9 14 1 1
+18 16 20 1
+20 16 18 1
+13 7 32 1
+32 7 13 1
+19 17 33 1
+33 17 19 1
+33 12 19 1
+19 12 33 1
+3 10 8 1
+8 10 3 1
+17 16 33 1
+33 16 17 1
+17 19 2 1
+2 19 17 1
+7 15 23 1
+23 15 7 1
+29 9 14 1
+14 9 29 1
+32 31 14 1
+14 31 32 1
+14 23 32 1
+32 23 14 1
+2 12 17 1
+17 12 2 1
+23 22 31 1
+31 22 23 1
+12 2 19 1
+19 2 12 1
+21 29 30 1
+30 29 21 1
+21 25 29 1
+29 25 21 1
+11 17 12 1
+12 17 11 1
+27 10 5 1
+5 10 27 1
+26 27 6 1
+6 27 26 1
+27 26 9 1
+9 26 27 1
+24 16 11 1
+11 16 24 1
+31 30 14 1
+14 30 31 1
+14 21 22 1
+22 21 14 1
+14 22 23 1
+23 22 14 1
+33 18 11 1
+11 18 33 1
+16 18 33 1
+33 18 16 1
+24 20 16 1
+16 20 24 1
+5 10 3 1
+3 10 5 1
+21 14 18 1
+18 14 21 1
+18 20 21 1
+21 20 18 1
+6 27 5 1
+5 27 6 1
+34 27 9 1
+9 27 34 1
+26 6 28 1
+28 6 26 1
+9 6 34 1
+34 6 9 1
+34 5 8 1
+8 5 34 1
+11 12 33 1
+33 12 11 1
+7 31 32 1
+32 31 7 1
+17 11 16 1
+16 11 17 1
+24 11 18 1
+18 11 24 1
+32 23 15 1
+15 23 32 1
+6 9 28 1
+28 9 6 1
+25 28 29 1
+29 28 25 1
+20 24 25 1
+25 24 20 1
+25 21 20 1
+20 21 25 1
+25 24 18 1
+18 24 25 1
+29 28 9 1
+9 28 29 1
+18 9 25 1
+25 9 18 1
+30 29 14 1
+14 29 30 1
+7 23 31 1
+31 23 7 1
+34 6 5 1
+5 6 34 1
+Tetrahedra
+32
+26 25 9 28 1
+7 13 4 15 1
+1 18 14 9 1
+15 7 13 32 1
+12 17 19 33 1
+23 14 31 32 1
+17 19 2 12 1
+27 10 5 34 1
+27 26 9 6 1
+16 11 18 33 1
+18 16 20 24 1
+3 10 8 5 1
+30 22 14 31 1
+9 18 14 29 1
+14 22 23 31 1
+9 27 6 34 1
+5 10 8 34 1
+21 14 18 29 1
+12 11 17 33 1
+11 16 18 24 1
+23 7 15 32 1
+21 18 25 29 1
+6 26 9 28 1
+20 18 24 25 1
+18 20 21 25 1
+18 9 25 29 1
+16 17 11 33 1
+28 25 9 29 1
+29 21 14 30 1
+22 14 21 30 1
+7 23 31 32 1
+6 27 5 34 1
+End
+
+DEAL::
+DEAL::
diff --git a/tests/cgal/input_grids/cube.off b/tests/cgal/input_grids/cube.off
new file mode 100644
index 0000000000..e200aa44bd
--- /dev/null
+++ b/tests/cgal/input_grids/cube.off
@@ -0,0 +1,22 @@
+OFF
+8 12 0
+-1 -1 -1
+-1 1 -1
+1 1 -1
+1 -1 -1
+-1 -1 1
+-1 1 1
+1 1 1
+1 -1 1
+3  0 1 3
+3  3 1 2
+3  0 4 1
+3  1 4 5
+3  3 2 7
+3  7 2 6
+3  4 0 3
+3  7 4 3
+3  6 4 7
+3  6 5 4
+3  1 5 6
+3  2 1 6
\ No newline at end of file
diff --git a/tests/cgal/input_grids/tetrahedron.off b/tests/cgal/input_grids/tetrahedron.off
new file mode 100644
index 0000000000..9e02601667
--- /dev/null
+++ b/tests/cgal/input_grids/tetrahedron.off
@@ -0,0 +1,11 @@
+OFF
+4 4 0
+
+0 1 0
+1 0 0
+0 0 0
+0 0 1
+3  0 1 2
+3  2 3 0
+3  1 3 2
+3  0 3 1