From 33330b801f50044116211635e3667ee484877168 Mon Sep 17 00:00:00 2001 From: Matthias Maier Date: Tue, 26 Sep 2017 11:40:25 -0500 Subject: [PATCH] add a test --- .../deal.II/numerics/vector_tools.templates.h | 20 ++++++--- .../generalized_interpolation_04.output | 4 +- .../numerics/generalized_interpolation_05.cc | 45 +++++++++++++++++++ .../generalized_interpolation_05.output | 3 ++ 4 files changed, 63 insertions(+), 9 deletions(-) create mode 100644 tests/numerics/generalized_interpolation_05.cc create mode 100644 tests/numerics/generalized_interpolation_05.output diff --git a/include/deal.II/numerics/vector_tools.templates.h b/include/deal.II/numerics/vector_tools.templates.h index 04966b4bff..c1c98a6208 100644 --- a/include/deal.II/numerics/vector_tools.templates.h +++ b/include/deal.II/numerics/vector_tools.templates.h @@ -238,15 +238,21 @@ namespace VectorTools // A small helper function to transform a component range starting // at offset from the real to the unit cell according to the - // supplied conformity. + // supplied conformity. The function_values vector is transformed + // in place. // // FIXME: This should be refactored into the mapping (i.e. - // implement the inverse function of Mapping::transform). Further, the finite element should make - // the information about the correct mapping directly accessible - - // fe.conforming_space is not the right call (thing about BDM). - const auto transform = [&](const typename FiniteElementData::Conformity conformity, - unsigned int offset) + // implement the inverse function of Mapping::transform). + // Further, the finite element should make the information about + // the correct mapping directly accessible (i.e. which MappingType + // should be used). Using fe.conforming_space might be a bit of a + // problem because we only support doing nothing, Hcurl, and Hdiv + // conforming mappings. + // + + const auto transform = [&function_values, &fe_values_jacobians, &cell]( + const typename FiniteElementData::Conformity conformity, + const unsigned int offset) { switch (conformity) { diff --git a/tests/numerics/generalized_interpolation_04.output b/tests/numerics/generalized_interpolation_04.output index 59027680b3..ce616a49ad 100644 --- a/tests/numerics/generalized_interpolation_04.output +++ b/tests/numerics/generalized_interpolation_04.output @@ -1,3 +1,3 @@ -DEAL::dim 2 FESystem<2>[FE_RaviartThomas<2>(1)^2-FE_Q<2>(1)-FESystem<2>[FE_Nedelec<2>(2)^2]-FESystem<2>[FE_Q<2>(1)^2]] -DEAL::Check projection property: 20.8500 +DEAL::dim 2 FESystem<2>[FE_RaviartThomas<2>(1)-FE_Q<2>(1)-FE_Nedelec<2>(2)^2-FESystem<2>[FE_Q<2>(1)^2]] +DEAL::Check projection property: 0.00000 diff --git a/tests/numerics/generalized_interpolation_05.cc b/tests/numerics/generalized_interpolation_05.cc new file mode 100644 index 0000000000..7d82fe3918 --- /dev/null +++ b/tests/numerics/generalized_interpolation_05.cc @@ -0,0 +1,45 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2017 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 at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + +// Check projection property of VectorTools::interpolate for a complex, +// staggered system of Hdiv / Hcurl / L2 conforming spaces. + +#include +#include +#include +#include + +#include "../tests.h" +#include "generalized_interpolation.h" + +int main () +{ + initlog(); + + static const int dim = 2; + + FESystem fe(FE_RaviartThomas(1), + 2, + FE_Q(1), + 1, + FESystem(FE_Nedelec(2), 2), + 1, + FESystem(FE_Q(1), dim), + 1); + + const unsigned int n_comp = fe.n_components(); + + test<2>(fe, F<2>(n_comp, 1), 1, false); +} diff --git a/tests/numerics/generalized_interpolation_05.output b/tests/numerics/generalized_interpolation_05.output new file mode 100644 index 0000000000..e6dc071228 --- /dev/null +++ b/tests/numerics/generalized_interpolation_05.output @@ -0,0 +1,3 @@ + +DEAL::dim 2 FESystem<2>[FE_RaviartThomas<2>(1)^2-FE_Q<2>(1)-FESystem<2>[FE_Nedelec<2>(2)^2]-FESystem<2>[FE_Q<2>(1)^2]] +DEAL::Check projection property: 0.00000 -- 2.39.5