From 1be223c9a48fa02cd2b72fa05af78e71d970138a Mon Sep 17 00:00:00 2001 From: Peter Munch Date: Sun, 8 May 2022 19:45:22 +0200 Subject: [PATCH] Add GenericDoFsPerObject::generate() --- include/deal.II/fe/fe_base.h | 87 ++++++++++++++++++++++ tests/fe/generic_dofs_per_object_01.cc | 78 +++++++++++++++++++ tests/fe/generic_dofs_per_object_01.output | 79 ++++++++++++++++++++ 3 files changed, 244 insertions(+) create mode 100644 tests/fe/generic_dofs_per_object_01.cc create mode 100644 tests/fe/generic_dofs_per_object_01.output diff --git a/include/deal.II/fe/fe_base.h b/include/deal.II/fe/fe_base.h index f817e16345..adc0aed987 100644 --- a/include/deal.II/fe/fe_base.h +++ b/include/deal.II/fe/fe_base.h @@ -28,6 +28,12 @@ DEAL_II_NAMESPACE_OPEN +// Forward declarations: +#ifndef DOXYGEN +template +class FiniteElementData; +#endif + /** * A namespace solely for the purpose of defining the Domination enum as well * as associated operators. @@ -192,6 +198,13 @@ namespace internal * First index of an object within a face. */ std::vector> first_object_index_on_face; + + /** + * Function that fills the fields based on a provided finite element. + */ + template + static GenericDoFsPerObject + generate(const FiniteElementData &fe); }; } // namespace internal @@ -959,6 +972,80 @@ FiniteElementData::get_first_face_quad_index( face_no]; } +template +internal::GenericDoFsPerObject +internal::GenericDoFsPerObject::generate(const FiniteElementData &fe) +{ + const auto reference_cell = fe.reference_cell(); + + internal::GenericDoFsPerObject result; + + result.dofs_per_object_exclusive.resize(dim + 1); + result.dofs_per_object_inclusive.resize(dim + 1); + result.object_index.resize(dim + 1); + + unsigned int counter = 0; + + for (unsigned int v : reference_cell.vertex_indices()) + { + const auto c = fe.template n_dofs_per_object<0>(v); + + result.dofs_per_object_exclusive[0].emplace_back(c); + result.dofs_per_object_inclusive[0].emplace_back(c); + result.object_index[0].emplace_back(counter); + + counter += c; + } + + if (dim >= 2) + for (unsigned int l : reference_cell.line_indices()) + { + const auto c = fe.template n_dofs_per_object<1>(l); + + result.dofs_per_object_exclusive[1].emplace_back(c); + result.dofs_per_object_inclusive[1].emplace_back( + c + 2 * fe.template n_dofs_per_object<0>()); + result.object_index[1].emplace_back(counter); + + counter += c; + } + + if (dim == 3) + for (unsigned int f : reference_cell.face_indices()) + { + const auto c = fe.template n_dofs_per_object<2>(f); + + result.dofs_per_object_exclusive[2].emplace_back(c); + result.dofs_per_object_inclusive[2].emplace_back(fe.n_dofs_per_face(f)); + result.object_index[2].emplace_back(counter); + + counter += c; + } + + { + result.dofs_per_object_exclusive[dim].emplace_back( + fe.template n_dofs_per_object()); + result.dofs_per_object_inclusive[dim].emplace_back(fe.n_dofs_per_cell()); + result.object_index[dim].emplace_back(counter); + } + + result.first_object_index_on_face.resize(dim); + for (unsigned int face_no : reference_cell.face_indices()) + { + result.first_object_index_on_face[0].emplace_back(0); + + if (dim >= 2) + result.first_object_index_on_face[1].emplace_back( + fe.get_first_face_line_index(face_no)); + + if (dim == 3) + result.first_object_index_on_face[2].emplace_back( + fe.get_first_face_quad_index(face_no)); + } + + return result; +} + #endif // DOXYGEN diff --git a/tests/fe/generic_dofs_per_object_01.cc b/tests/fe/generic_dofs_per_object_01.cc new file mode 100644 index 0000000000..854b58c18d --- /dev/null +++ b/tests/fe/generic_dofs_per_object_01.cc @@ -0,0 +1,78 @@ +// --------------------------------------------------------------------- +// +// 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. +// +// --------------------------------------------------------------------- + +// Test internal::GenericDoFsPerObject::generate(). + +#include +#include +#include +#include + +#include "../tests.h" + +using namespace dealii; + +template +void +test(const FiniteElement &fe) +{ + const auto dpo = internal::GenericDoFsPerObject::generate(fe); + + deallog << fe.get_name() << std::endl; + + for (const auto &row : dpo.dofs_per_object_exclusive) + for (const auto &i : row) + deallog << i << " "; + deallog << std::endl; + + for (const auto &row : dpo.dofs_per_object_inclusive) + for (const auto &i : row) + deallog << i << " "; + deallog << std::endl; + + for (const auto &row : dpo.object_index) + for (const auto &i : row) + deallog << i << " "; + deallog << std::endl; + + for (const auto &row : dpo.first_object_index_on_face) + for (const auto &i : row) + deallog << i << " "; + deallog << std::endl; + deallog << std::endl; +} + +int +main() +{ + initlog(); + + test(FE_Q<1>(1)); + test(FE_Q<1>(2)); + test(FE_Q<2>(1)); + test(FE_Q<2>(2)); + test(FE_Q<3>(1)); + test(FE_Q<3>(2)); + + test(FE_SimplexP<2>(1)); + test(FE_SimplexP<2>(2)); + test(FE_SimplexP<3>(1)); + test(FE_SimplexP<3>(2)); + + test(FE_PyramidP<3>(1)); + + test(FE_WedgeP<3>(1)); + test(FE_WedgeP<3>(2)); +} diff --git a/tests/fe/generic_dofs_per_object_01.output b/tests/fe/generic_dofs_per_object_01.output new file mode 100644 index 0000000000..9ad3e7cf86 --- /dev/null +++ b/tests/fe/generic_dofs_per_object_01.output @@ -0,0 +1,79 @@ + +DEAL::FE_Q<1>(1) +DEAL::1 1 0 +DEAL::1 1 2 +DEAL::0 1 2 +DEAL::0 0 +DEAL:: +DEAL::FE_Q<1>(2) +DEAL::1 1 1 +DEAL::1 1 3 +DEAL::0 1 2 +DEAL::0 0 +DEAL:: +DEAL::FE_Q<2>(1) +DEAL::1 1 1 1 0 0 0 0 0 +DEAL::1 1 1 1 2 2 2 2 4 +DEAL::0 1 2 3 4 4 4 4 4 +DEAL::0 0 0 0 2 2 2 2 +DEAL:: +DEAL::FE_Q<2>(2) +DEAL::1 1 1 1 1 1 1 1 1 +DEAL::1 1 1 1 3 3 3 3 9 +DEAL::0 1 2 3 4 5 6 7 8 +DEAL::0 0 0 0 2 2 2 2 +DEAL:: +DEAL::FE_Q<3>(1) +DEAL::1 1 1 1 1 1 1 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +DEAL::1 1 1 1 1 1 1 1 2 2 2 2 2 2 2 2 2 2 2 2 4 4 4 4 4 4 8 +DEAL::0 1 2 3 4 5 6 7 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 8 +DEAL::0 0 0 0 0 0 4 4 4 4 4 4 4 4 4 4 4 4 +DEAL:: +DEAL::FE_Q<3>(2) +DEAL::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 +DEAL::1 1 1 1 1 1 1 1 3 3 3 3 3 3 3 3 3 3 3 3 9 9 9 9 9 9 27 +DEAL::0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 +DEAL::0 0 0 0 0 0 4 4 4 4 4 4 8 8 8 8 8 8 +DEAL:: +DEAL::FE_SimplexP<2>(1) +DEAL::1 1 1 0 0 0 0 +DEAL::1 1 1 2 2 2 3 +DEAL::0 1 2 3 3 3 3 +DEAL::0 0 0 2 2 2 +DEAL:: +DEAL::FE_SimplexP<2>(2) +DEAL::1 1 1 1 1 1 0 +DEAL::1 1 1 3 3 3 6 +DEAL::0 1 2 3 4 5 6 +DEAL::0 0 0 2 2 2 +DEAL:: +DEAL::FE_SimplexP<3>(1) +DEAL::1 1 1 1 0 0 0 0 0 0 0 0 0 0 0 +DEAL::1 1 1 1 2 2 2 2 2 2 3 3 3 3 4 +DEAL::0 1 2 3 4 4 4 4 4 4 4 4 4 4 4 +DEAL::0 0 0 0 3 3 3 3 3 3 3 3 +DEAL:: +DEAL::FE_SimplexP<3>(2) +DEAL::1 1 1 1 1 1 1 1 1 1 0 0 0 0 0 +DEAL::1 1 1 1 3 3 3 3 3 3 6 6 6 6 10 +DEAL::0 1 2 3 4 5 6 7 8 9 10 10 10 10 10 +DEAL::0 0 0 0 3 3 3 3 6 6 6 6 +DEAL:: +DEAL::FE_PyramidP<3>(1) +DEAL::1 1 1 1 1 0 0 0 0 0 0 0 0 0 0 0 0 0 +DEAL::1 1 1 1 1 2 2 2 2 2 2 2 4 3 3 3 3 5 +DEAL::0 1 2 3 4 5 5 5 5 5 5 5 5 5 5 5 5 5 +DEAL::0 0 0 0 0 4 3 3 3 3 4 3 3 3 3 +DEAL:: +DEAL::FE_WedgeP<3>(1) +DEAL::1 1 1 1 1 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +DEAL::1 1 1 1 1 1 2 2 2 2 2 2 2 2 2 3 3 4 4 4 6 +DEAL::0 1 2 3 4 5 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 +DEAL::0 0 0 0 0 3 3 4 4 4 3 3 4 4 4 +DEAL:: +DEAL::FE_WedgeP<3>(2) +DEAL::1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 0 0 1 1 1 0 +DEAL::1 1 1 1 1 1 3 3 3 3 3 3 3 3 3 6 6 9 9 9 18 +DEAL::0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 15 15 16 17 18 +DEAL::0 0 0 0 0 3 3 4 4 4 6 6 8 8 8 +DEAL:: -- 2.39.5