From 67e70b0c69535daed453652cf4c385e9a3104f5e Mon Sep 17 00:00:00 2001 From: Peter Munch Date: Mon, 18 Jan 2021 07:15:22 +0100 Subject: [PATCH] Implement ReferenceCell::vertex_faces() --- include/deal.II/grid/reference_cell.h | 75 +++++++++++++++++++++ source/grid/grid_tools.cc | 38 +++++------ source/grid/grid_tools_dof_handlers.cc | 38 +++++------ tests/simplex/reference_cell_kind_01.cc | 59 ++++++++++++++++ tests/simplex/reference_cell_kind_01.output | 40 +++++++++++ 5 files changed, 208 insertions(+), 42 deletions(-) create mode 100644 tests/simplex/reference_cell_kind_01.cc create mode 100644 tests/simplex/reference_cell_kind_01.output diff --git a/include/deal.II/grid/reference_cell.h b/include/deal.II/grid/reference_cell.h index ca56d20f89..dda8190245 100644 --- a/include/deal.II/grid/reference_cell.h +++ b/include/deal.II/grid/reference_cell.h @@ -18,6 +18,7 @@ #include +#include #include @@ -108,6 +109,12 @@ namespace ReferenceCell void serialize(Archive &archive, const unsigned int /*version*/); + /** + * Return a vector of faces a @p vertex belongs to. + */ + ArrayView + faces_for_given_vertex(const unsigned int vertex) const; + private: /** * The variable that stores what this object actually corresponds to. @@ -182,6 +189,74 @@ namespace ReferenceCell + inline ArrayView + Type::faces_for_given_vertex(const unsigned int vertex) const + { + if (kind == Type::Line) + { + AssertIndexRange(vertex, GeometryInfo<1>::vertices_per_cell); + return {&GeometryInfo<2>::vertex_to_face[vertex][0], 1}; + } + else if (kind == Type::Quad) + { + AssertIndexRange(vertex, GeometryInfo<2>::vertices_per_cell); + return {&GeometryInfo<2>::vertex_to_face[vertex][0], 2}; + } + else if (kind == Type::Hex) + { + AssertIndexRange(vertex, GeometryInfo<3>::vertices_per_cell); + return {&GeometryInfo<3>::vertex_to_face[vertex][0], 3}; + } + else if (kind == Type::Tri) + { + AssertIndexRange(vertex, 3); + static const std::array, 3> table = { + {{{0, 2}}, {{0, 1}}, {{1, 2}}}}; + + return table[vertex]; + } + else if (kind == Type::Tet) + { + AssertIndexRange(vertex, 4); + static const std::array, 4> table = { + {{{0, 1, 2}}, {{0, 1, 3}}, {{0, 2, 3}}, {{1, 2, 3}}}}; + + return table[vertex]; + } + else if (kind == Type::Wedge) + { + AssertIndexRange(vertex, 6); + static const std::array, 6> table = { + {{{0, 2, 4}}, + {{0, 2, 3}}, + {{0, 3, 4}}, + {{1, 2, 4}}, + {{1, 2, 3}}, + {{1, 3, 4}}}}; + + return table[vertex]; + } + else if (kind == Type::Pyramid) + { + AssertIndexRange(vertex, 5); + static const unsigned int X = numbers::invalid_unsigned_int; + static const std::array, 5> table = { + {{{0, 1, 3, X}}, + {{0, 2, 3, X}}, + {{0, 1, 4, X}}, + {{0, 2, 4, X}}, + {{1, 2, 3, 4}}}}; + + return {&table[vertex][0], vertex == 4 ? 4u : 3u}; + } + + Assert(false, ExcNotImplemented()); + + return {}; + } + + + /** * Return the dimension of the given reference-cell type @p type. */ diff --git a/source/grid/grid_tools.cc b/source/grid/grid_tools.cc index a4fdbb3c0c..1d163e6589 100644 --- a/source/grid/grid_tools.cc +++ b/source/grid/grid_tools.cc @@ -1545,27 +1545,23 @@ namespace GridTools // (possibly) coarser neighbor. if this is the case, // then we need to also add this neighbor if (dim >= 2) - for (unsigned int vface = 0; vface < dim; vface++) - { - const unsigned int face = - GeometryInfo::vertex_to_face[v][vface]; // TODO - - if (!cell->at_boundary(face) && - cell->neighbor(face)->is_active()) - { - // there is a (possibly) coarser cell behind a - // face to which the vertex belongs. the - // vertex we are looking at is then either a - // vertex of that coarser neighbor, or it is a - // hanging node on one of the faces of that - // cell. in either case, it is adjacent to the - // vertex, so add it to the list as well (if - // the cell was already in the list then the - // std::set makes sure that we get it only - // once) - adjacent_cells.insert(cell->neighbor(face)); - } - } + for (const auto face : + cell->reference_cell_type().faces_for_given_vertex(v)) + if (!cell->at_boundary(face) && + cell->neighbor(face)->is_active()) + { + // there is a (possibly) coarser cell behind a + // face to which the vertex belongs. the + // vertex we are looking at is then either a + // vertex of that coarser neighbor, or it is a + // hanging node on one of the faces of that + // cell. in either case, it is adjacent to the + // vertex, so add it to the list as well (if + // the cell was already in the list then the + // std::set makes sure that we get it only + // once) + adjacent_cells.insert(cell->neighbor(face)); + } // in any case, we have found a cell, so go to the next cell goto next_cell; diff --git a/source/grid/grid_tools_dof_handlers.cc b/source/grid/grid_tools_dof_handlers.cc index 339d64e1f0..8444430e82 100644 --- a/source/grid/grid_tools_dof_handlers.cc +++ b/source/grid/grid_tools_dof_handlers.cc @@ -266,27 +266,23 @@ namespace GridTools // (possibly) coarser neighbor. if this is the case, // then we need to also add this neighbor if (dim >= 2) - for (unsigned int vface = 0; vface < dim; vface++) - { - const unsigned int face = - GeometryInfo::vertex_to_face[v][vface]; - - if (!cell->at_boundary(face) && - cell->neighbor(face)->is_active()) - { - // there is a (possibly) coarser cell behind a - // face to which the vertex belongs. the - // vertex we are looking at is then either a - // vertex of that coarser neighbor, or it is a - // hanging node on one of the faces of that - // cell. in either case, it is adjacent to the - // vertex, so add it to the list as well (if - // the cell was already in the list then the - // std::set makes sure that we get it only - // once) - adjacent_cells.insert(cell->neighbor(face)); - } - } + for (const auto face : + cell->reference_cell_type().faces_for_given_vertex(v)) + if (!cell->at_boundary(face) && + cell->neighbor(face)->is_active()) + { + // there is a (possibly) coarser cell behind a + // face to which the vertex belongs. the + // vertex we are looking at is then either a + // vertex of that coarser neighbor, or it is a + // hanging node on one of the faces of that + // cell. in either case, it is adjacent to the + // vertex, so add it to the list as well (if + // the cell was already in the list then the + // std::set makes sure that we get it only + // once) + adjacent_cells.insert(cell->neighbor(face)); + } // in any case, we have found a cell, so go to the next cell goto next_cell; diff --git a/tests/simplex/reference_cell_kind_01.cc b/tests/simplex/reference_cell_kind_01.cc new file mode 100644 index 0000000000..a11f25d51b --- /dev/null +++ b/tests/simplex/reference_cell_kind_01.cc @@ -0,0 +1,59 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 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. +// +// --------------------------------------------------------------------- + + +// Test ReferenceCell::Kind::faces_for_given_vertex(). + + +#include + +#include + +#include "../tests.h" + +using namespace dealii; + +template +void +test(const ReferenceCell::Type &reference_cell) +{ + const auto kind = ReferenceCell::Type(reference_cell); + const auto &info = ReferenceCell::internal::Info::get_cell(reference_cell); + + for (const auto v : info.vertex_indices()) + { + deallog << v << ": "; + for (const auto i : kind.faces_for_given_vertex(v)) + deallog << i << " "; + deallog << std::endl; + } + deallog << std::endl; +} + + + +int +main() +{ + initlog(); + + test<2>(ReferenceCell::Type::Line); + test<2>(ReferenceCell::Type::Tri); + test<2>(ReferenceCell::Type::Quad); + test<3>(ReferenceCell::Type::Tet); + test<3>(ReferenceCell::Type::Pyramid); + test<3>(ReferenceCell::Type::Wedge); + test<3>(ReferenceCell::Type::Hex); +} diff --git a/tests/simplex/reference_cell_kind_01.output b/tests/simplex/reference_cell_kind_01.output new file mode 100644 index 0000000000..86988d0fb0 --- /dev/null +++ b/tests/simplex/reference_cell_kind_01.output @@ -0,0 +1,40 @@ + +DEAL::0: 0 +DEAL::1: 1 +DEAL:: +DEAL::0: 0 2 +DEAL::1: 0 1 +DEAL::2: 1 2 +DEAL:: +DEAL::0: 0 2 +DEAL::1: 1 2 +DEAL::2: 0 3 +DEAL::3: 1 3 +DEAL:: +DEAL::0: 0 1 2 +DEAL::1: 0 1 3 +DEAL::2: 0 2 3 +DEAL::3: 1 2 3 +DEAL:: +DEAL::0: 0 1 3 +DEAL::1: 0 2 3 +DEAL::2: 0 1 4 +DEAL::3: 0 2 4 +DEAL::4: 1 2 3 4 +DEAL:: +DEAL::0: 0 2 4 +DEAL::1: 0 2 3 +DEAL::2: 0 3 4 +DEAL::3: 1 2 4 +DEAL::4: 1 2 3 +DEAL::5: 1 3 4 +DEAL:: +DEAL::0: 0 2 4 +DEAL::1: 1 2 4 +DEAL::2: 0 3 4 +DEAL::3: 1 3 4 +DEAL::4: 0 2 5 +DEAL::5: 1 2 5 +DEAL::6: 0 3 5 +DEAL::7: 1 3 5 +DEAL:: -- 2.39.5