From 3dc8cd6b71b83b81e0d279c3a32d3ee35bd64c26 Mon Sep 17 00:00:00 2001 From: Elias Dejene Date: Thu, 17 Dec 2020 21:55:39 +0100 Subject: [PATCH] Implement inheritance of manifold id_s in 2D. --- source/grid/grid_generator.cc | 306 +++++++++++++++--- tests/simplex/conv_hex_to_simplex.cc | 7 + ..._to_simplex.with_simplex_support=on.output | 38 ++- 3 files changed, 299 insertions(+), 52 deletions(-) diff --git a/source/grid/grid_generator.cc b/source/grid/grid_generator.cc index 852f918279..da8586353d 100644 --- a/source/grid/grid_generator.cc +++ b/source/grid/grid_generator.cc @@ -7101,8 +7101,16 @@ namespace GridGenerator Assert(dim > 1, ExcNotImplemented()); - // static tables with the definitions of cells and boundary faces for 2D - // and 3D + /* 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 + * 3D, also inner-edges and boundary-edges need to be defined. + */ + + /* Cell definition 2D: + * A quadrilateral element is converted to 8 simplices elements. Each + * triangle is defined by 3 vertices. + */ static const std::array, 8> table_2D_cell = { {{{0, 6, 4}}, {{8, 4, 6}}, @@ -7113,24 +7121,39 @@ namespace GridGenerator {{8, 5, 7}}, {{3, 7, 5}}}}; - static const std::array, 24> table_3D_cell = { - {{{0, 1, 12, 10}}, {{2, 3, 11, 12}}, {{7, 6, 11, 13}}, - {{5, 4, 13, 10}}, {{0, 2, 8, 12}}, {{4, 6, 13, 8}}, - {{5, 13, 7, 9}}, {{1, 9, 3, 12}}, {{0, 8, 4, 10}}, - {{1, 5, 9, 10}}, {{3, 7, 11, 9}}, {{2, 6, 8, 11}}, - {{12, 13, 10, 9}}, {{12, 13, 9, 11}}, {{12, 13, 11, 8}}, - {{12, 13, 8, 10}}, {{13, 8, 10, 4}}, {{13, 10, 9, 5}}, - {{13, 9, 11, 7}}, {{13, 11, 8, 6}}, {{10, 12, 9, 1}}, - {{9, 12, 11, 3}}, {{11, 12, 8, 2}}, {{8, 12, 10, 0}}}}; - + /* Cell definition 3D: + * A hexahedron element is converted to 24 tetrahedron elements. Each + * tetrahedron is defined by 4 vertices. + */ + static const std::array, 24> + vertex_ids_for_cells_3d = { + {{{0, 1, 12, 10}}, {{2, 3, 11, 12}}, {{7, 6, 11, 13}}, + {{5, 4, 13, 10}}, {{0, 2, 8, 12}}, {{4, 6, 13, 8}}, + {{5, 13, 7, 9}}, {{1, 9, 3, 12}}, {{0, 8, 4, 10}}, + {{1, 5, 9, 10}}, {{3, 7, 11, 9}}, {{2, 6, 8, 11}}, + {{12, 13, 10, 9}}, {{12, 13, 9, 11}}, {{12, 13, 11, 8}}, + {{12, 13, 8, 10}}, {{13, 8, 10, 4}}, {{13, 10, 9, 5}}, + {{13, 9, 11, 7}}, {{13, 11, 8, 6}}, {{10, 12, 9, 1}}, + {{9, 12, 11, 3}}, {{11, 12, 8, 2}}, {{8, 12, 10, 0}}}}; + + /* Boundary-faces 2D: + * After converting, each of the 4 quadrilateral faces is defined by faces + * of 2 different triangles, i.e., lines. Note that lines are defined by 2 + * vertices. + */ static const std::array, 2>, 4> - table_2D_boundary_faces = {{{{{{0, 4}}, {{4, 2}}}}, - {{{{1, 5}}, {{5, 3}}}}, - {{{{0, 6}}, {{6, 1}}}}, - {{{{2, 7}}, {{7, 3}}}}}}; - + vertex_ids_for_boundary_faces_2d = {{{{{{0, 4}}, {{4, 2}}}}, + {{{{1, 5}}, {{5, 3}}}}, + {{{{0, 6}}, {{6, 1}}}}, + {{{{2, 7}}, {{7, 3}}}}}}; + + /* Boundary-faces 3D: + * After converting, each of the 6 hexahedron faces corresponds to faces of + * 4 different tetrahedron faces, i.e., triangles. Note that a triangle is + * defined by 3 vertices. + */ static const std::array, 4>, 6> - table_3D_boundary_faces = { + vertex_ids_for_boundary_faces_3d = { {{{{{0, 4, 8}}, {{4, 8, 6}}, {{8, 6, 2}}, {{0, 2, 8}}}}, {{{{1, 3, 9}}, {{3, 9, 7}}, {{9, 7, 5}}, {{1, 9, 5}}}}, {{{{0, 1, 10}}, {{1, 10, 5}}, {{10, 5, 4}}, {{0, 10, 4}}}}, @@ -7138,6 +7161,120 @@ namespace GridGenerator {{{{0, 1, 12}}, {{1, 12, 3}}, {{12, 3, 2}}, {{0, 12, 2}}}}, {{{{4, 5, 13}}, {{5, 13, 7}}, {{13, 7, 6}}, {{4, 13, 6}}}}}}; + /* Inner-faces 2D: + * The converted triangulation based on simplices has 8 faces that do not + * form the boundary, i.e. inner-faces, each defined by 2 vertices. + */ + static const std::array, 8> + vertex_ids_for_inner_faces_2d = {{{{6, 4}}, + {{6, 8}}, + {{6, 5}}, + {{4, 8}}, + {{8, 5}}, + {{7, 4}}, + {{7, 8}}, + {{7, 5}}}}; + + /* Inner-faces 3D: + * The converted triangulation based on simplices has 72 faces that do not + * form the boundary, i.e. inner-faces, each defined by 3 vertices. + */ + static const std::array, 72> + vertex_ids_for_inner_faces_3d = { + {{{0, 12, 10}}, {{12, 1, 10}}, {{12, 1, 9}}, {{12, 3, 9}}, + {{12, 2, 11}}, {{12, 3, 11}}, {{12, 0, 8}}, {{12, 2, 8}}, + {{9, 13, 5}}, {{13, 7, 9}}, {{11, 7, 13}}, {{11, 6, 13}}, + {{4, 8, 13}}, {{6, 8, 13}}, {{4, 13, 10}}, {{13, 5, 10}}, + {{10, 9, 5}}, {{10, 9, 1}}, {{11, 9, 7}}, {{11, 9, 3}}, + {{8, 11, 2}}, {{8, 11, 6}}, {{8, 10, 0}}, {{8, 10, 4}}, + {{12, 3, 9}}, {{12, 9, 11}}, {{12, 3, 11}}, {{3, 9, 11}}, + {{2, 12, 8}}, {{2, 12, 11}}, {{2, 11, 8}}, {{8, 12, 11}}, + {{0, 12, 10}}, {{0, 12, 8}}, {{0, 8, 10}}, {{8, 10, 12}}, + {{12, 1, 10}}, {{12, 1, 9}}, {{1, 10, 9}}, {{10, 9, 12}}, + {{10, 8, 4}}, {{10, 8, 13}}, {{4, 13, 8}}, {{4, 13, 10}}, + {{10, 9, 13}}, {{10, 9, 5}}, {{13, 5, 10}}, {{13, 5, 9}}, + {{13, 7, 9}}, {{13, 7, 11}}, {{9, 11, 13}}, {{9, 11, 7}}, + {{8, 11, 13}}, {{8, 11, 6}}, {{6, 13, 8}}, {{6, 13, 11}}, + {{12, 13, 10}}, {{12, 13, 8}}, {{8, 10, 13}}, {{8, 10, 12}}, + {{12, 13, 10}}, {{12, 13, 9}}, {{10, 9, 13}}, {{10, 9, 12}}, + {{12, 13, 9}}, {{12, 13, 11}}, {{9, 11, 13}}, {{9, 11, 12}}, + {{12, 13, 11}}, {{12, 13, 8}}, {{8, 11, 13}}, {{8, 11, 12}}}}; + + /* Inner-edges 3D: + * The converted triangulation based on simplices has 60 edges that do not + * coincide with the boundary, i.e. inner-edges, each defined by 2 vertices. + */ + static const std::array, 60> + vertex_ids_for_inner_edges_3d = { + {{{12, 10}}, {{12, 9}}, {{12, 11}}, {{12, 8}}, {{9, 13}}, + {{11, 13}}, {{8, 13}}, {{10, 13}}, {{10, 9}}, {{9, 11}}, + {{11, 8}}, {{8, 10}}, {{12, 9}}, {{12, 11}}, {{11, 9}}, + {{12, 8}}, {{12, 11}}, {{11, 8}}, {{12, 8}}, {{12, 10}}, + {{10, 8}}, {{12, 10}}, {{12, 9}}, {{9, 10}}, {{13, 10}}, + {{13, 8}}, {{8, 10}}, {{13, 10}}, {{13, 9}}, {{9, 10}}, + {{13, 11}}, {{13, 9}}, {{11, 9}}, {{13, 11}}, {{13, 8}}, + {{11, 8}}, {{12, 13}}, {{8, 10}}, {{8, 13}}, {{10, 13}}, + {{8, 12}}, {{10, 12}}, {{12, 13}}, {{10, 9}}, {{10, 13}}, + {{9, 13}}, {{10, 12}}, {{9, 12}}, {{12, 13}}, {{9, 11}}, + {{9, 13}}, {{11, 13}}, {{9, 12}}, {{11, 12}}, {{12, 13}}, + {{11, 8}}, {{11, 13}}, {{8, 13}}, {{11, 12}}, {{8, 12}}}}; + + /* Boundary-edges 3D: + * For each of the 6 boundary-faces of the hexahedron, there are 8 edges (of + * different tetrahedrons) that coincide with the boundary, i.e. + * boundary-edges. Each boundary-edge is defined by 2 vertices. + */ + static const std::array, 8>, 6> + vertex_ids_for_boundary_edges_3d = {{{{{{4, 6}}, + {{4, 8}}, + {{6, 8}}, + {{4, 0}}, + {{6, 2}}, + {{0, 8}}, + {{2, 8}}, + {{0, 2}}}}, + {{{{5, 7}}, + {{5, 9}}, + {{7, 9}}, + {{5, 1}}, + {{7, 3}}, + {{1, 9}}, + {{3, 9}}, + {{1, 3}}}}, + {{{{4, 5}}, + {{4, 10}}, + {{5, 10}}, + {{4, 0}}, + {{5, 1}}, + {{0, 10}}, + {{1, 10}}, + {{0, 1}}}}, + {{{{6, 7}}, + {{6, 11}}, + {{7, 11}}, + {{6, 2}}, + {{7, 3}}, + {{2, 11}}, + {{3, 11}}, + {{2, 3}}}}, + {{{{2, 3}}, + {{2, 12}}, + {{3, 12}}, + {{2, 0}}, + {{3, 1}}, + {{0, 12}}, + {{1, 12}}, + {{0, 1}}}}, + {{{{6, 7}}, + {{6, 13}}, + {{7, 13}}, + {{6, 4}}, + {{7, 5}}, + {{4, 13}}, + {{5, 14}}, + {{4, 5}}}}}}; + + std::vector> vertices; std::vector> cells; SubCellData subcell_data; @@ -7205,22 +7342,44 @@ namespace GridGenerator // helper function for creating cells and subcells const auto add_cell = [&](const unsigned int struct_dim, const auto & index_vertices, - const unsigned int material_or_boundary_id) { + const unsigned int material_or_boundary_id, + const unsigned int manifold_id = 0) { + // sub-cell data only has to be stored if the information differs + // from the default + if (struct_dim < dim && + (material_or_boundary_id == numbers::internal_face_boundary_id && + manifold_id == numbers::flat_manifold_id)) + return; + if (struct_dim == dim) // cells { + if (dim == 2) + { + AssertDimension(index_vertices.size(), 3); + } + else if (dim == 3) + { + AssertDimension(index_vertices.size(), 4); + } + CellData cell_data(index_vertices.size()); for (unsigned int i = 0; i < index_vertices.size(); ++i) { cell_data.vertices[i] = local_vertex_indices[index_vertices[i]]; - cell_data.material_id = material_or_boundary_id; + cell_data.material_id = + material_or_boundary_id; // inherit material id + cell_data.manifold_id = + manifold_id; // inherit cell-manifold id } cells.push_back(cell_data); } - else if (dim == 2 && struct_dim == 1) // an edge of a triangle + else if (dim == 2 && struct_dim == 1) // an edge of a simplex { + Assert(index_vertices.size() == 2, ExcInternalError()); CellData<1> boundary_line(2); boundary_line.boundary_id = material_or_boundary_id; + boundary_line.manifold_id = manifold_id; for (unsigned int i = 0; i < index_vertices.size(); ++i) { boundary_line.vertices[i] = @@ -7228,10 +7387,12 @@ namespace GridGenerator } subcell_data.boundary_lines.push_back(boundary_line); } - else if (dim == 3 && struct_dim == 2) // a face of a thetrahedron + else if (dim == 3 && struct_dim == 2) // a face of a tetrahedron { + Assert(index_vertices.size() == 3, ExcInternalError()); CellData<2> boundary_quad(3); - boundary_quad.boundary_id = material_or_boundary_id; + boundary_quad.material_id = material_or_boundary_id; + boundary_quad.manifold_id = manifold_id; for (unsigned int i = 0; i < index_vertices.size(); ++i) { boundary_quad.vertices[i] = @@ -7239,40 +7400,101 @@ namespace GridGenerator } subcell_data.boundary_quads.push_back(boundary_quad); } + else if (dim == 3 && struct_dim == 1) // an edge of a tetrahedron + { + Assert(index_vertices.size() == 2, ExcInternalError()); + CellData<1> boundary_line(2); + boundary_line.manifold_id = manifold_id; + for (unsigned int i = 0; i < index_vertices.size(); ++i) + { + boundary_line.vertices[i] = + local_vertex_indices[index_vertices[i]]; + } + subcell_data.boundary_lines.push_back(boundary_line); + } else { Assert(false, ExcNotImplemented()); } }; - const auto current_material_id = cell.material_id(); + const auto material_id_cell = cell.material_id(); // create cells one by one - if (dim == 2) - for (const auto &cell : table_2D_cell) - add_cell(dim, cell, current_material_id); + { + // get cell-manifold id from current quad cell + const auto manifold_id_cell = cell.manifold_id(); + // inherit cell manifold + for (const auto &cell_vertices : table_2D_cell) + add_cell(dim, cell_vertices, material_id_cell, manifold_id_cell); + + // inherit inner manifold (faces) + for (const auto &face_vertices : vertex_ids_for_inner_faces_2d) + // set inner tri-faces according to cell-manifold of quad + // element, set invalid b_id, since we do not want to modify + // b_id inside + add_cell(1, + face_vertices, + numbers::internal_face_boundary_id, + manifold_id_cell); + } else if (dim == 3) - for (const auto &cell : table_3D_cell) - add_cell(dim, cell, current_material_id); + { + // get cell-manifold id from current quad cell + const auto manifold_id_cell = cell.manifold_id(); + // inherit cell manifold + for (const auto &cell_vertices : vertex_ids_for_cells_3d) + add_cell(dim, cell_vertices, material_id_cell, manifold_id_cell); + + // set manifold of inner FACES of tets according to + // hex-cell-manifold + for (const auto &face_vertices : vertex_ids_for_inner_faces_3d) + add_cell(2, + face_vertices, + numbers::internal_face_boundary_id, + manifold_id_cell); + + // set manifold of inner EDGES of tets according to + // hex-cell-manifold + for (const auto &edge_vertices : vertex_ids_for_inner_edges_3d) + add_cell(1, + edge_vertices, + numbers::internal_face_boundary_id, + manifold_id_cell); + } else Assert(false, ExcNotImplemented()); // Set up sub-cell data. for (const auto f : cell.face_indices()) - if (cell.face(f)->boundary_id() != - numbers::internal_face_boundary_id) // face at boundary - { - const auto current_boundary_id = cell.face(f)->boundary_id(); - if (dim == 2) // 2D boundary_lines - for (const auto &face : table_2D_boundary_faces[f]) - add_cell(1, face, current_boundary_id); - else if (dim == 3) // 3D boundary_quads - for (const auto &face : table_3D_boundary_faces[f]) - add_cell(2, face, current_boundary_id); - else - Assert(false, ExcNotImplemented()); - } + { + const auto bid = cell.face(f)->boundary_id(); + const auto mid = cell.face(f)->manifold_id(); + + // process boundary-faces: set boundary and manifold ids + if (dim == 2) // 2D boundary-faces + for (const auto &face_vertices : + vertex_ids_for_boundary_faces_2d[f]) + add_cell(1, face_vertices, bid, mid); + + else if (dim == 3) // 3D boundary-faces + { + // set manifold id of tet-boundary-faces according to + // hex-boundary-faces + for (const auto &face_vertices : + vertex_ids_for_boundary_faces_3d[f]) + add_cell(2, face_vertices, bid, mid); + // set manifold id of tet-boundary-edges according to + // hex-boundary-faces + for (const auto &edge_vertices : + vertex_ids_for_boundary_edges_3d[f]) + add_cell(1, edge_vertices, bid, mid); + } + + else + Assert(false, ExcNotImplemented()); + } } out_tria.create_triangulation(vertices, cells, subcell_data); diff --git a/tests/simplex/conv_hex_to_simplex.cc b/tests/simplex/conv_hex_to_simplex.cc index 94e58315af..10dab721d2 100644 --- a/tests/simplex/conv_hex_to_simplex.cc +++ b/tests/simplex/conv_hex_to_simplex.cc @@ -80,6 +80,13 @@ check_file() // for dim = spaceim 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)); + + // out_tria.refine_global(2); + // write 2 outputs (total mesh and only surface mesh) const auto grid_out = [](const auto &tria, const bool surface_mesh_only = false) { diff --git a/tests/simplex/conv_hex_to_simplex.with_simplex_support=on.output b/tests/simplex/conv_hex_to_simplex.with_simplex_support=on.output index 2b4d9978ac..aa358d6d45 100644 --- a/tests/simplex/conv_hex_to_simplex.with_simplex_support=on.output +++ b/tests/simplex/conv_hex_to_simplex.with_simplex_support=on.output @@ -75,7 +75,7 @@ LOOKUP_TABLE default SCALARS ManifoldID int 1 LOOKUP_TABLE default -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 -1 +-1 -1 -1 0 0 -1 0 -1 0 -1 # vtk DataFile Version 3.0 Triangulation generated with deal.II @@ -128,7 +128,7 @@ LOOKUP_TABLE default SCALARS ManifoldID int 1 LOOKUP_TABLE default --1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 -1 +-1 -1 -1 -1 -1 0 0 -1 0 -1 0 -1 DEAL:2D: conversion triangulation with quad elements to tri elements: ::OK! # vtk DataFile Version 3.0 Triangulation generated with deal.II @@ -568,7 +568,7 @@ POINTS 33 double 0.293743 0.00000 0.680892 0.364976 0.364976 0.856495 -CELLS 140 656 +CELLS 158 710 4 0 1 12 10 4 2 3 11 12 4 7 6 11 13 @@ -709,24 +709,42 @@ CELLS 140 656 3 24 29 30 3 4 30 29 3 16 32 17 +2 14 15 +2 14 20 +2 15 20 +2 15 17 +2 17 20 +2 17 25 +2 15 25 +2 16 14 +2 16 20 +2 29 32 +2 17 16 +2 17 32 +2 25 23 +2 24 17 +2 24 32 +2 24 25 +2 24 23 +2 29 24 -CELL_TYPES 140 +CELL_TYPES 158 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 10 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 5 +3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 - -CELL_DATA 140 +CELL_DATA 158 SCALARS MaterialID int 1 LOOKUP_TABLE default 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 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 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 1 2 2 1 3 2 8 7 5 2 5 8 1 9 10 1 7 9 3 5 4 5 4 6 10 3 4 10 11 4 3 8 7 6 8 11 6 7 6 10 11 9 9 11 - +0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 SCALARS ManifoldID int 1 LOOKUP_TABLE default -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 -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 -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 -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 -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 -1 0 0 -1 -1 0 -1 0 0 -1 -1 -1 0 -1 0 0 -1 0 -1 0 -1 -1 0 +0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 # vtk DataFile Version 3.0 Triangulation generated with deal.II ASCII @@ -828,5 +846,5 @@ LOOKUP_TABLE default SCALARS ManifoldID int 1 LOOKUP_TABLE default --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 -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 -1 -1 -1 -1 -1 -1 0 -1 0 0 -1 -1 0 -1 0 0 -1 -1 -1 0 -1 0 0 -1 0 -1 0 -1 -1 0 DEAL:3D: conversion triangulation with tet elements to hex elements: ::OK! -- 2.39.5