Merge pull request #12790 from drwells/simplex-convert-refined

Permit conversion of multi-level triangulations to simplices.

source/grid/grid_generator.cc

@@ -7631,11 +7631,17 @@ namespace GridGenerator
   convert_hypercube_to_simplex_mesh(const Triangulation<dim, spacedim> &in_tria,
                                     Triangulation<dim, spacedim> &out_tria)
   {
-    Assert(in_tria.n_global_levels() == 1,
-           ExcMessage("Number of global levels has to be 1."));
-
     Assert(dim > 1, ExcNotImplemented());
 
+    Triangulation<dim, spacedim> temp_tria;
+    if (in_tria.n_global_levels() > 1)
+      {
+        AssertThrow(!in_tria.has_hanging_nodes(), ExcNotImplemented());
+        flatten_triangulation(in_tria, temp_tria);
+      }
+    const Triangulation<dim, spacedim> &ref_tria =
+      in_tria.n_global_levels() > 1 ? temp_tria : in_tria;
+
     /* static tables with the definitions of cells, faces and edges by its
      * vertices for 2D and 3D. For the inheritance of the manifold_id,
      * definitions of inner-faces and boundary-faces are required. In case of
@@ -7811,9 +7817,9 @@ namespace GridGenerator
     // store for each vertex and face the assigned index so that we only
     // assign them a value once
     std::vector<unsigned int> old_to_new_vertex_indices(
-      in_tria.n_vertices(), numbers::invalid_unsigned_int);
+      ref_tria.n_vertices(), numbers::invalid_unsigned_int);
     std::vector<unsigned int> face_to_new_vertex_indices(
-      in_tria.n_faces(), numbers::invalid_unsigned_int);
+      ref_tria.n_faces(), numbers::invalid_unsigned_int);
 
     // We first have to create all of the new vertices. To do this, we loop over
     // all cells and on each cell
@@ -7822,7 +7828,7 @@ namespace GridGenerator
     // (ii) create new midpoint vertex locations for each face (and record their
     // new indices in the 'face_to_new_vertex_indices' vector),
     // (iii) create new midpoint vertex locations for each cell (dim = 2 only)
-    for (const auto &cell : in_tria)
+    for (const auto &cell : ref_tria)
       {
         // temporary array storing the global indices of each cell entity in the
         // sequence: vertices, edges/faces, cell

tests/simplex/conv_hex_to_simplex_02.cc

@@ -0,0 +1,133 @@
+/* ---------------------------------------------------------------------
+ *
+ * Copyright (C) 2020 - 2021 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.
+ *
+ * ---------------------------------------------------------------------
+ *
+ * Like conv_hex_to_simplex, but also refines the grid once.
+ */
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_out.h>
+#include <deal.II/grid/grid_tools.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/tria_iterator.h>
+
+#include <fstream>
+#include <iostream>
+
+#include "../tests.h"
+
+using namespace dealii;
+
+template <int dim, int spacedim>
+void
+create_triangulation(Triangulation<dim, spacedim> &triangulation)
+{
+  GridGenerator::subdivided_hyper_cube(triangulation, 4);
+}
+
+template <int dim>
+void
+create_triangulation(Triangulation<dim, dim> &triangulation)
+{
+  GridGenerator::quarter_hyper_ball(triangulation);
+}
+
+template <int dim, int spacedim>
+void
+check_file() // for dim = spaceim
+{
+  Triangulation<dim, spacedim> in_tria, out_tria;
+  create_triangulation(in_tria);
+
+  // make each cell a different material id
+  unsigned int m_id = 0;
+  for (const auto &cell : in_tria)
+    {
+      cell.set_material_id(m_id++);
+    }
+
+  // set different boundary ids and output
+  unsigned int b_id = 0;
+  for (const auto &cell : in_tria)
+    {
+      for (const auto f : cell.face_indices())
+        {
+          if (cell.face(f)->at_boundary())
+            {
+              cell.face(f)->set_boundary_id(b_id);
+              b_id++;
+            }
+        }
+    }
+  in_tria.refine_global(1);
+
+  GridGenerator::convert_hypercube_to_simplex_mesh(in_tria, out_tria);
+
+  // copy manifolds to test global refining
+  for (const auto i : in_tria.get_manifold_ids())
+    if (i != numbers::flat_manifold_id)
+      out_tria.set_manifold(i, in_tria.get_manifold(i));
+
+  // write 2 outputs (total mesh and only surface mesh)
+  const auto grid_out = [](const auto &tria,
+                           const bool  surface_mesh_only = false) {
+    GridOutFlags::Vtk flags;
+
+    if (surface_mesh_only)
+      {
+        flags.output_cells         = false;
+        flags.output_faces         = true;
+        flags.output_edges         = false;
+        flags.output_only_relevant = false;
+      }
+
+    GridOut grid_out;
+    grid_out.set_flags(flags);
+
+    grid_out.write_vtk(tria, deallog.get_file_stream());
+  };
+
+  grid_out(out_tria);       // total mesh
+  grid_out(out_tria, true); // only surface mesh
+
+  deallog << "OK!" << std::endl;
+}
+
+
+int
+main()
+{
+  initlog();
+  // TRIANGULAR ELEMENTS
+  // dim = spacedim = 2
+  deallog.push(
+    "2D: conversion triangulation with quad elements to tri elements: ");
+  check_file<2, 2>();
+  deallog.pop();
+
+  // TETRAHEDRAL ELEMENTS
+  // dim = 2, spacedim = 2
+  deallog.push(
+    "2D: conversion triangulation with quad elements to tri elements: ");
+  check_file<2, 3>();
+  deallog.pop();
+
+  // TETRAHEDRAL ELEMENTS
+  // dim = spacedim = 3
+  deallog.push(
+    "3D: conversion triangulation with tet elements to hex elements: ");
+  check_file<3, 3>();
+  deallog.pop();
+}