From 0fe915e3bfc6632baa222e54057a74fe721331c9 Mon Sep 17 00:00:00 2001 From: Jonathan Robey Date: Sun, 24 Jul 2016 10:41:40 -0700 Subject: [PATCH] Add test of functionality --- .../deal.II/numerics/vector_tools.templates.h | 28 +-- tests/dofs/interpolate_q_system_mixed_01.cc | 189 ++++++++++++++++++ ...interpolate_q_system_mixed_01.debug.output | 120 +++++++++++ 3 files changed, 319 insertions(+), 18 deletions(-) create mode 100644 tests/dofs/interpolate_q_system_mixed_01.cc create mode 100644 tests/dofs/interpolate_q_system_mixed_01.debug.output diff --git a/include/deal.II/numerics/vector_tools.templates.h b/include/deal.II/numerics/vector_tools.templates.h index ab6c5ff009..cc8987ebc6 100644 --- a/include/deal.II/numerics/vector_tools.templates.h +++ b/include/deal.II/numerics/vector_tools.templates.h @@ -122,29 +122,21 @@ namespace VectorTools { if (component_mask[component_index] == true) { - std::pair base_index = fe[fe_index].component_to_base_index(component_index); - - AssertThrow ((fe[fe_index].base_element(base_index.first).has_support_points()) || - (fe[fe_index].base_element(base_index.first).dofs_per_cell == 0), - ExcNonInterpolatingFE()); + Assert ((fe[fe_index].base_element(fe[fe_index].component_to_base_index(component_index).first).has_support_points()) || + (fe[fe_index].base_element(fe[fe_index].component_to_base_index(component_index).first).dofs_per_cell == 0), + ExcNonInterpolatingFE()); } } } - // Find the support points on a cell that - // are mentioned multiple times, and ony add - // each once. - // Each multiple point gets to know the dof - // index of its representative point by the - // dof_to_rep_dof_table. + // Find the support points on a cell that are mentioned multiple times, and + // only add each once. Each multiple point gets to know the dof index of + // its representative point by the dof_to_rep_dof_table. // the following vector collects all unit support points p[i], - // 0<=i > > rep_unit_support_points (fe.size()); // the following table converts a dof i // to the rep index. @@ -160,7 +152,7 @@ namespace VectorTools = fe[fe_index].system_to_component_index(i).first; if (component_mask[component] == true) { - Point dof_support_point = fe[fe_index].unit_support_point(i); + const Point dof_support_point = fe[fe_index].unit_support_point(i); bool representative=true; // the following loop is looped diff --git a/tests/dofs/interpolate_q_system_mixed_01.cc b/tests/dofs/interpolate_q_system_mixed_01.cc new file mode 100644 index 0000000000..0b5570c08f --- /dev/null +++ b/tests/dofs/interpolate_q_system_mixed_01.cc @@ -0,0 +1,189 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2006 - 2015 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 that VectorTools::interpolate with a component mask works for +// FE_System(FE_Q(p)) elements correctly on a uniformly refined mesh for +// functions of degree q when a non-interpolating (FE_DGP(p)) element is also +// included in the system. + +#include "../tests.h" +#include +#include +#include +#include + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include +#include + + +template +class F : public Function +{ +public: + F (const unsigned int q, const double adj) : Function(4), q(q), adj(adj) {} + + virtual void vector_value (const Point &p, + Vector &v) const + { + for (unsigned int c=0; c +void test () +{ + Triangulation triangulation; + GridGenerator::hyper_cube (triangulation); + triangulation.refine_global (3); + + for (unsigned int p=1; p<6-dim; ++p) + { + FE_Q fe_1(p); + FE_Q fe_2(p+1); + FE_DGP fe_3(p); + FESystem fe(fe_1, 2, + fe_2, 1, + fe_3, 1); + DoFHandler dof_handler(triangulation); + dof_handler.distribute_dofs (fe); + + // Use constant offset to distinguish between masks + const double adj1 = 0.3; + ComponentSelectFunction select_mask1(0,4); + ComponentMask mask1(4, false); + mask1.set(0, true); + + const double adj2 = 1.7; + ComponentSelectFunction select_mask2(std::make_pair(1,3),4); + ComponentMask mask2(4, false); + mask2.set(1, true); + mask2.set(2, true); + + ComponentMask mask_fail(4, false); + mask_fail.set(3, true); + + Vector interpolant (dof_handler.n_dofs()); + Vector error (triangulation.n_active_cells()); + for (unsigned int q=0; q<=p+2; ++q) + { + // interpolate the function with mask 1 + VectorTools::interpolate (dof_handler, + F (q, adj1), + interpolant, + mask1); + + // interpolate the function with mask 2 + VectorTools::interpolate (dof_handler, + F (q, adj2), + interpolant, + mask2); + + // then compute the interpolation error for mask 1 + VectorTools::integrate_difference (dof_handler, + interpolant, + F (q, adj1), + error, + QGauss(q+2), + VectorTools::L2_norm, + &select_mask1); + if (q<=p) + Assert (error.l2_norm() < 1e-12*interpolant.l2_norm(), + ExcInternalError()); + + deallog << fe.get_name() << ", Mask 1, P_" << q + << ", rel. error=" << error.l2_norm() / interpolant.l2_norm() + << std::endl; + + // then compute the interpolation error for mask 2 + VectorTools::integrate_difference (dof_handler, + interpolant, + F (q, adj2), + error, + QGauss(q+2), + VectorTools::L2_norm, + &select_mask2); + if (q<=p) + Assert (error.l2_norm() < 1e-12*interpolant.l2_norm(), + ExcInternalError()); + + deallog << fe.get_name() << ", Mask 2, P_" << q + << ", rel. error=" << error.l2_norm() / interpolant.l2_norm() + << std::endl; + } + + // Test for correct failure + deallog << fe.get_name() + << ", Fails for including non-interpolating FE "; + try + { + VectorTools::interpolate (dof_handler, + F (0, 0.0), + interpolant, + mask_fail); + } + catch (const ExceptionBase &e) + { + deallog << "OK" << std::endl; + deallog << "\tFails with " << e.get_exc_name() << std::endl; + } + deallog << std::endl; + } +} + + + +int main () +{ + deal_II_exceptions::disable_abort_on_exception(); + + std::ofstream logfile("output"); + deallog << std::setprecision (3); + + deallog.attach(logfile); + deallog.threshold_double(1.e-10); + + test<1>(); + test<2>(); + test<3>(); +} + diff --git a/tests/dofs/interpolate_q_system_mixed_01.debug.output b/tests/dofs/interpolate_q_system_mixed_01.debug.output new file mode 100644 index 0000000000..a8843b62ea --- /dev/null +++ b/tests/dofs/interpolate_q_system_mixed_01.debug.output @@ -0,0 +1,120 @@ + +DEAL::FESystem<1>[FE_Q<1>(1)^2-FE_Q<1>(2)-FE_DGP<1>(1)], Mask 1, P_0, rel. error=0 +DEAL::FESystem<1>[FE_Q<1>(1)^2-FE_Q<1>(2)-FE_DGP<1>(1)], Mask 2, P_0, rel. error=0 +DEAL::FESystem<1>[FE_Q<1>(1)^2-FE_Q<1>(2)-FE_DGP<1>(1)], Mask 1, P_1, rel. error=0 +DEAL::FESystem<1>[FE_Q<1>(1)^2-FE_Q<1>(2)-FE_DGP<1>(1)], Mask 2, P_1, rel. error=0 +DEAL::FESystem<1>[FE_Q<1>(1)^2-FE_Q<1>(2)-FE_DGP<1>(1)], Mask 1, P_2, rel. error=0.000124 +DEAL::FESystem<1>[FE_Q<1>(1)^2-FE_Q<1>(2)-FE_DGP<1>(1)], Mask 2, P_2, rel. error=0.000124 +DEAL::FESystem<1>[FE_Q<1>(1)^2-FE_Q<1>(2)-FE_DGP<1>(1)], Mask 1, P_3, rel. error=0.000296 +DEAL::FESystem<1>[FE_Q<1>(1)^2-FE_Q<1>(2)-FE_DGP<1>(1)], Mask 2, P_3, rel. error=0.000296 +DEAL::FESystem<1>[FE_Q<1>(1)^2-FE_Q<1>(2)-FE_DGP<1>(1)], Fails for including non-interpolating FE OK +DEAL:: Fails with ExcNonInterpolatingFE() +DEAL:: +DEAL::FESystem<1>[FE_Q<1>(2)^2-FE_Q<1>(3)-FE_DGP<1>(2)], Mask 1, P_0, rel. error=0 +DEAL::FESystem<1>[FE_Q<1>(2)^2-FE_Q<1>(3)-FE_DGP<1>(2)], Mask 2, P_0, rel. error=0 +DEAL::FESystem<1>[FE_Q<1>(2)^2-FE_Q<1>(3)-FE_DGP<1>(2)], Mask 1, P_1, rel. error=0 +DEAL::FESystem<1>[FE_Q<1>(2)^2-FE_Q<1>(3)-FE_DGP<1>(2)], Mask 2, P_1, rel. error=0 +DEAL::FESystem<1>[FE_Q<1>(2)^2-FE_Q<1>(3)-FE_DGP<1>(2)], Mask 1, P_2, rel. error=0 +DEAL::FESystem<1>[FE_Q<1>(2)^2-FE_Q<1>(3)-FE_DGP<1>(2)], Mask 2, P_2, rel. error=0 +DEAL::FESystem<1>[FE_Q<1>(2)^2-FE_Q<1>(3)-FE_DGP<1>(2)], Mask 1, P_3, rel. error=2.31e-06 +DEAL::FESystem<1>[FE_Q<1>(2)^2-FE_Q<1>(3)-FE_DGP<1>(2)], Mask 2, P_3, rel. error=2.31e-06 +DEAL::FESystem<1>[FE_Q<1>(2)^2-FE_Q<1>(3)-FE_DGP<1>(2)], Mask 1, P_4, rel. error=6.92e-06 +DEAL::FESystem<1>[FE_Q<1>(2)^2-FE_Q<1>(3)-FE_DGP<1>(2)], Mask 2, P_4, rel. error=6.92e-06 +DEAL::FESystem<1>[FE_Q<1>(2)^2-FE_Q<1>(3)-FE_DGP<1>(2)], Fails for including non-interpolating FE OK +DEAL:: Fails with ExcNonInterpolatingFE() +DEAL:: +DEAL::FESystem<1>[FE_Q<1>(3)^2-FE_Q<1>(4)-FE_DGP<1>(3)], Mask 1, P_0, rel. error=0 +DEAL::FESystem<1>[FE_Q<1>(3)^2-FE_Q<1>(4)-FE_DGP<1>(3)], Mask 2, P_0, rel. error=0 +DEAL::FESystem<1>[FE_Q<1>(3)^2-FE_Q<1>(4)-FE_DGP<1>(3)], Mask 1, P_1, rel. error=0 +DEAL::FESystem<1>[FE_Q<1>(3)^2-FE_Q<1>(4)-FE_DGP<1>(3)], Mask 2, P_1, rel. error=0 +DEAL::FESystem<1>[FE_Q<1>(3)^2-FE_Q<1>(4)-FE_DGP<1>(3)], Mask 1, P_2, rel. error=0 +DEAL::FESystem<1>[FE_Q<1>(3)^2-FE_Q<1>(4)-FE_DGP<1>(3)], Mask 2, P_2, rel. error=0 +DEAL::FESystem<1>[FE_Q<1>(3)^2-FE_Q<1>(4)-FE_DGP<1>(3)], Mask 1, P_3, rel. error=0 +DEAL::FESystem<1>[FE_Q<1>(3)^2-FE_Q<1>(4)-FE_DGP<1>(3)], Mask 2, P_3, rel. error=0 +DEAL::FESystem<1>[FE_Q<1>(3)^2-FE_Q<1>(4)-FE_DGP<1>(3)], Mask 1, P_4, rel. error=5.66e-08 +DEAL::FESystem<1>[FE_Q<1>(3)^2-FE_Q<1>(4)-FE_DGP<1>(3)], Mask 2, P_4, rel. error=5.66e-08 +DEAL::FESystem<1>[FE_Q<1>(3)^2-FE_Q<1>(4)-FE_DGP<1>(3)], Mask 1, P_5, rel. error=2.03e-07 +DEAL::FESystem<1>[FE_Q<1>(3)^2-FE_Q<1>(4)-FE_DGP<1>(3)], Mask 2, P_5, rel. error=2.03e-07 +DEAL::FESystem<1>[FE_Q<1>(3)^2-FE_Q<1>(4)-FE_DGP<1>(3)], Fails for including non-interpolating FE OK +DEAL:: Fails with ExcNonInterpolatingFE() +DEAL:: +DEAL::FESystem<1>[FE_Q<1>(4)^2-FE_Q<1>(5)-FE_DGP<1>(4)], Mask 1, P_0, rel. error=0 +DEAL::FESystem<1>[FE_Q<1>(4)^2-FE_Q<1>(5)-FE_DGP<1>(4)], Mask 2, P_0, rel. error=0 +DEAL::FESystem<1>[FE_Q<1>(4)^2-FE_Q<1>(5)-FE_DGP<1>(4)], Mask 1, P_1, rel. error=0 +DEAL::FESystem<1>[FE_Q<1>(4)^2-FE_Q<1>(5)-FE_DGP<1>(4)], Mask 2, P_1, rel. error=0 +DEAL::FESystem<1>[FE_Q<1>(4)^2-FE_Q<1>(5)-FE_DGP<1>(4)], Mask 1, P_2, rel. error=0 +DEAL::FESystem<1>[FE_Q<1>(4)^2-FE_Q<1>(5)-FE_DGP<1>(4)], Mask 2, P_2, rel. error=0 +DEAL::FESystem<1>[FE_Q<1>(4)^2-FE_Q<1>(5)-FE_DGP<1>(4)], Mask 1, P_3, rel. error=0 +DEAL::FESystem<1>[FE_Q<1>(4)^2-FE_Q<1>(5)-FE_DGP<1>(4)], Mask 2, P_3, rel. error=0 +DEAL::FESystem<1>[FE_Q<1>(4)^2-FE_Q<1>(5)-FE_DGP<1>(4)], Mask 1, P_4, rel. error=0 +DEAL::FESystem<1>[FE_Q<1>(4)^2-FE_Q<1>(5)-FE_DGP<1>(4)], Mask 2, P_4, rel. error=0 +DEAL::FESystem<1>[FE_Q<1>(4)^2-FE_Q<1>(5)-FE_DGP<1>(4)], Mask 1, P_5, rel. error=1.50e-09 +DEAL::FESystem<1>[FE_Q<1>(4)^2-FE_Q<1>(5)-FE_DGP<1>(4)], Mask 2, P_5, rel. error=1.50e-09 +DEAL::FESystem<1>[FE_Q<1>(4)^2-FE_Q<1>(5)-FE_DGP<1>(4)], Mask 1, P_6, rel. error=6.27e-09 +DEAL::FESystem<1>[FE_Q<1>(4)^2-FE_Q<1>(5)-FE_DGP<1>(4)], Mask 2, P_6, rel. error=6.27e-09 +DEAL::FESystem<1>[FE_Q<1>(4)^2-FE_Q<1>(5)-FE_DGP<1>(4)], Fails for including non-interpolating FE OK +DEAL:: Fails with ExcNonInterpolatingFE() +DEAL:: +DEAL::FESystem<2>[FE_Q<2>(1)^2-FE_Q<2>(2)-FE_DGP<2>(1)], Mask 1, P_0, rel. error=0 +DEAL::FESystem<2>[FE_Q<2>(1)^2-FE_Q<2>(2)-FE_DGP<2>(1)], Mask 2, P_0, rel. error=0 +DEAL::FESystem<2>[FE_Q<2>(1)^2-FE_Q<2>(2)-FE_DGP<2>(1)], Mask 1, P_1, rel. error=0 +DEAL::FESystem<2>[FE_Q<2>(1)^2-FE_Q<2>(2)-FE_DGP<2>(1)], Mask 2, P_1, rel. error=0 +DEAL::FESystem<2>[FE_Q<2>(1)^2-FE_Q<2>(2)-FE_DGP<2>(1)], Mask 1, P_2, rel. error=4.17e-05 +DEAL::FESystem<2>[FE_Q<2>(1)^2-FE_Q<2>(2)-FE_DGP<2>(1)], Mask 2, P_2, rel. error=4.17e-05 +DEAL::FESystem<2>[FE_Q<2>(1)^2-FE_Q<2>(2)-FE_DGP<2>(1)], Mask 1, P_3, rel. error=9.68e-05 +DEAL::FESystem<2>[FE_Q<2>(1)^2-FE_Q<2>(2)-FE_DGP<2>(1)], Mask 2, P_3, rel. error=9.68e-05 +DEAL::FESystem<2>[FE_Q<2>(1)^2-FE_Q<2>(2)-FE_DGP<2>(1)], Fails for including non-interpolating FE OK +DEAL:: Fails with ExcNonInterpolatingFE() +DEAL:: +DEAL::FESystem<2>[FE_Q<2>(2)^2-FE_Q<2>(3)-FE_DGP<2>(2)], Mask 1, P_0, rel. error=0 +DEAL::FESystem<2>[FE_Q<2>(2)^2-FE_Q<2>(3)-FE_DGP<2>(2)], Mask 2, P_0, rel. error=0 +DEAL::FESystem<2>[FE_Q<2>(2)^2-FE_Q<2>(3)-FE_DGP<2>(2)], Mask 1, P_1, rel. error=0 +DEAL::FESystem<2>[FE_Q<2>(2)^2-FE_Q<2>(3)-FE_DGP<2>(2)], Mask 2, P_1, rel. error=0 +DEAL::FESystem<2>[FE_Q<2>(2)^2-FE_Q<2>(3)-FE_DGP<2>(2)], Mask 1, P_2, rel. error=0 +DEAL::FESystem<2>[FE_Q<2>(2)^2-FE_Q<2>(3)-FE_DGP<2>(2)], Mask 2, P_2, rel. error=0 +DEAL::FESystem<2>[FE_Q<2>(2)^2-FE_Q<2>(3)-FE_DGP<2>(2)], Mask 1, P_3, rel. error=4.87e-07 +DEAL::FESystem<2>[FE_Q<2>(2)^2-FE_Q<2>(3)-FE_DGP<2>(2)], Mask 2, P_3, rel. error=4.87e-07 +DEAL::FESystem<2>[FE_Q<2>(2)^2-FE_Q<2>(3)-FE_DGP<2>(2)], Mask 1, P_4, rel. error=1.46e-06 +DEAL::FESystem<2>[FE_Q<2>(2)^2-FE_Q<2>(3)-FE_DGP<2>(2)], Mask 2, P_4, rel. error=1.46e-06 +DEAL::FESystem<2>[FE_Q<2>(2)^2-FE_Q<2>(3)-FE_DGP<2>(2)], Fails for including non-interpolating FE OK +DEAL:: Fails with ExcNonInterpolatingFE() +DEAL:: +DEAL::FESystem<2>[FE_Q<2>(3)^2-FE_Q<2>(4)-FE_DGP<2>(3)], Mask 1, P_0, rel. error=0 +DEAL::FESystem<2>[FE_Q<2>(3)^2-FE_Q<2>(4)-FE_DGP<2>(3)], Mask 2, P_0, rel. error=0 +DEAL::FESystem<2>[FE_Q<2>(3)^2-FE_Q<2>(4)-FE_DGP<2>(3)], Mask 1, P_1, rel. error=0 +DEAL::FESystem<2>[FE_Q<2>(3)^2-FE_Q<2>(4)-FE_DGP<2>(3)], Mask 2, P_1, rel. error=0 +DEAL::FESystem<2>[FE_Q<2>(3)^2-FE_Q<2>(4)-FE_DGP<2>(3)], Mask 1, P_2, rel. error=0 +DEAL::FESystem<2>[FE_Q<2>(3)^2-FE_Q<2>(4)-FE_DGP<2>(3)], Mask 2, P_2, rel. error=0 +DEAL::FESystem<2>[FE_Q<2>(3)^2-FE_Q<2>(4)-FE_DGP<2>(3)], Mask 1, P_3, rel. error=0 +DEAL::FESystem<2>[FE_Q<2>(3)^2-FE_Q<2>(4)-FE_DGP<2>(3)], Mask 2, P_3, rel. error=0 +DEAL::FESystem<2>[FE_Q<2>(3)^2-FE_Q<2>(4)-FE_DGP<2>(3)], Mask 1, P_4, rel. error=1.02e-08 +DEAL::FESystem<2>[FE_Q<2>(3)^2-FE_Q<2>(4)-FE_DGP<2>(3)], Mask 2, P_4, rel. error=1.02e-08 +DEAL::FESystem<2>[FE_Q<2>(3)^2-FE_Q<2>(4)-FE_DGP<2>(3)], Mask 1, P_5, rel. error=3.68e-08 +DEAL::FESystem<2>[FE_Q<2>(3)^2-FE_Q<2>(4)-FE_DGP<2>(3)], Mask 2, P_5, rel. error=3.68e-08 +DEAL::FESystem<2>[FE_Q<2>(3)^2-FE_Q<2>(4)-FE_DGP<2>(3)], Fails for including non-interpolating FE OK +DEAL:: Fails with ExcNonInterpolatingFE() +DEAL:: +DEAL::FESystem<3>[FE_Q<3>(1)^2-FE_Q<3>(2)-FE_DGP<3>(1)], Mask 1, P_0, rel. error=0 +DEAL::FESystem<3>[FE_Q<3>(1)^2-FE_Q<3>(2)-FE_DGP<3>(1)], Mask 2, P_0, rel. error=0 +DEAL::FESystem<3>[FE_Q<3>(1)^2-FE_Q<3>(2)-FE_DGP<3>(1)], Mask 1, P_1, rel. error=0 +DEAL::FESystem<3>[FE_Q<3>(1)^2-FE_Q<3>(2)-FE_DGP<3>(1)], Mask 2, P_1, rel. error=0 +DEAL::FESystem<3>[FE_Q<3>(1)^2-FE_Q<3>(2)-FE_DGP<3>(1)], Mask 1, P_2, rel. error=1.26e-05 +DEAL::FESystem<3>[FE_Q<3>(1)^2-FE_Q<3>(2)-FE_DGP<3>(1)], Mask 2, P_2, rel. error=1.26e-05 +DEAL::FESystem<3>[FE_Q<3>(1)^2-FE_Q<3>(2)-FE_DGP<3>(1)], Mask 1, P_3, rel. error=2.89e-05 +DEAL::FESystem<3>[FE_Q<3>(1)^2-FE_Q<3>(2)-FE_DGP<3>(1)], Mask 2, P_3, rel. error=2.89e-05 +DEAL::FESystem<3>[FE_Q<3>(1)^2-FE_Q<3>(2)-FE_DGP<3>(1)], Fails for including non-interpolating FE OK +DEAL:: Fails with ExcNonInterpolatingFE() +DEAL:: +DEAL::FESystem<3>[FE_Q<3>(2)^2-FE_Q<3>(3)-FE_DGP<3>(2)], Mask 1, P_0, rel. error=0 +DEAL::FESystem<3>[FE_Q<3>(2)^2-FE_Q<3>(3)-FE_DGP<3>(2)], Mask 2, P_0, rel. error=0 +DEAL::FESystem<3>[FE_Q<3>(2)^2-FE_Q<3>(3)-FE_DGP<3>(2)], Mask 1, P_1, rel. error=0 +DEAL::FESystem<3>[FE_Q<3>(2)^2-FE_Q<3>(3)-FE_DGP<3>(2)], Mask 2, P_1, rel. error=0 +DEAL::FESystem<3>[FE_Q<3>(2)^2-FE_Q<3>(3)-FE_DGP<3>(2)], Mask 1, P_2, rel. error=0 +DEAL::FESystem<3>[FE_Q<3>(2)^2-FE_Q<3>(3)-FE_DGP<3>(2)], Mask 2, P_2, rel. error=0 +DEAL::FESystem<3>[FE_Q<3>(2)^2-FE_Q<3>(3)-FE_DGP<3>(2)], Mask 1, P_3, rel. error=1.03e-07 +DEAL::FESystem<3>[FE_Q<3>(2)^2-FE_Q<3>(3)-FE_DGP<3>(2)], Mask 2, P_3, rel. error=1.03e-07 +DEAL::FESystem<3>[FE_Q<3>(2)^2-FE_Q<3>(3)-FE_DGP<3>(2)], Mask 1, P_4, rel. error=3.08e-07 +DEAL::FESystem<3>[FE_Q<3>(2)^2-FE_Q<3>(3)-FE_DGP<3>(2)], Mask 2, P_4, rel. error=3.08e-07 +DEAL::FESystem<3>[FE_Q<3>(2)^2-FE_Q<3>(3)-FE_DGP<3>(2)], Fails for including non-interpolating FE OK +DEAL:: Fails with ExcNonInterpolatingFE() +DEAL:: -- 2.39.5