From 2e453c2c7563544a09c923dfeb3724afe3fc2368 Mon Sep 17 00:00:00 2001 From: Simon Sticko Date: Tue, 1 Jun 2021 09:01:01 +0200 Subject: [PATCH] Implement face interpolation between FE_Bernstein and elements with support points. By using the implementation of get_face_interpolation_matrix in the base class (FE_Q_Base) when the incoming element isn't a Bernstein element. --- include/deal.II/fe/fe_bernstein.h | 8 +-- source/fe/fe_bernstein.cc | 12 ++-- .../face_interpolation_fe_Bernstein_2_fe_q.cc | 65 +++++++++++++++++++ ...e_interpolation_fe_Bernstein_2_fe_q.output | 10 +++ 4 files changed, 84 insertions(+), 11 deletions(-) create mode 100644 tests/fe/face_interpolation_fe_Bernstein_2_fe_q.cc create mode 100644 tests/fe/face_interpolation_fe_Bernstein_2_fe_q.output diff --git a/include/deal.II/fe/fe_bernstein.h b/include/deal.II/fe/fe_bernstein.h index 362d0d7ec4..4d471eef78 100644 --- a/include/deal.II/fe/fe_bernstein.h +++ b/include/deal.II/fe/fe_bernstein.h @@ -110,8 +110,8 @@ public: * the neighboring element. The size of the matrix is then * source.dofs_per_face times this->dofs_per_face. The * FE_Bernstein element family only provides interpolation matrices for - * elements of the same type and FE_Nothing. For all other elements, an - * exception of type + * elements of the same type, for elements that have support points, and + * FE_Nothing. For all other elements, an exception of type * FiniteElement::ExcInterpolationNotImplemented is thrown. */ virtual void @@ -124,8 +124,8 @@ public: * the neighboring element. The size of the matrix is then * source.dofs_per_face times this->dofs_per_face. The * FE_Bernstein element family only provides interpolation matrices for - * elements of the same type and FE_Nothing. For all other elements, an - * exception of type + * elements of the same type, for elements that have support points, and + * FE_Nothing. For all other elements, an exception of type * FiniteElement::ExcInterpolationNotImplemented is thrown. */ virtual void diff --git a/source/fe/fe_bernstein.cc b/source/fe/fe_bernstein.cc index cff7465e48..b25f172e44 100644 --- a/source/fe/fe_bernstein.cc +++ b/source/fe/fe_bernstein.cc @@ -195,15 +195,13 @@ FE_Bernstein::get_subface_interpolation_matrix( Assert(std::fabs(sum - 1) < eps, ExcInternalError()); } } - else if (dynamic_cast *>(&x_source_fe) != nullptr) + else { - // nothing to do here, the FE_Nothing has no degrees of freedom anyway + // When the incoming element is not FE_Bernstein we can just delegate to + // the base class to create the interpolation matrix. + FE_Q_Base::get_subface_interpolation_matrix( + x_source_fe, subface, interpolation_matrix, face_no); } - else - AssertThrow( - false, - (typename FiniteElement::ExcInterpolationNotImplemented())); } diff --git a/tests/fe/face_interpolation_fe_Bernstein_2_fe_q.cc b/tests/fe/face_interpolation_fe_Bernstein_2_fe_q.cc new file mode 100644 index 0000000000..8dc698a658 --- /dev/null +++ b/tests/fe/face_interpolation_fe_Bernstein_2_fe_q.cc @@ -0,0 +1,65 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2021 - 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 that we can use FE_Bernstein::get_face_interpolation_matrix to create +// an interpolation matrix between FE_Bernstein and an element that has support +// points. + +#include +#include + +#include + +#include "../tests.h" + + +// Set up an FE_Q element and an FE_Bernstein element with the incoming order. +// Call get_face_interpolation_matrix to create the interpolation matrix +// between these and write it to deallog. +template +void +create_and_print_face_interpolation_matrix(const unsigned int order) +{ + deallog << "dim = " << dim << ", order = " << order << std::endl; + + const FE_Q fe_q(order); + const FE_Bernstein fe_bernstein(order); + + FullMatrix interpolation_matrix(fe_q.n_dofs_per_face(), + fe_bernstein.n_dofs_per_face()); + const unsigned int face_index = 0; + + + fe_bernstein.get_face_interpolation_matrix(fe_q, + interpolation_matrix, + face_index); + + interpolation_matrix.print(deallog); + deallog << std::endl; +} + + + +int +main() +{ + initlog(); + + const int dim = 2; + + const std::vector orders = {1, 2}; + for (const unsigned int order : orders) + create_and_print_face_interpolation_matrix(order); +} diff --git a/tests/fe/face_interpolation_fe_Bernstein_2_fe_q.output b/tests/fe/face_interpolation_fe_Bernstein_2_fe_q.output new file mode 100644 index 0000000000..9096def154 --- /dev/null +++ b/tests/fe/face_interpolation_fe_Bernstein_2_fe_q.output @@ -0,0 +1,10 @@ + +DEAL::dim = 2, order = 1 +DEAL::1.0 0.0 +DEAL::0.0 1.0 +DEAL:: +DEAL::dim = 2, order = 2 +DEAL::1.0 0.0 0.0 +DEAL::0.0 1.0 0.0 +DEAL::0.25 0.25 0.50 +DEAL:: -- 2.39.5