--- /dev/null
+Changed: The hp-variants of KellyErrorEstimator::estimate now take
+hp::MappingCollection objects. The versions that require Mapping objects
+are now deprecated.
+<br>
+(Marc Fehling, 2025/06/20)
{
template <int>
class QCollection;
-}
+
+ template <int, int>
+ class MappingCollection;
+} // namespace hp
#endif
/**
* Equivalent to the set of functions above, except that this one takes a
- * quadrature collection for hp-finite element dof handlers.
+ * mapping and quadrature collection for hp-finite element dof handlers.
*/
template <typename Number>
static void
estimate(
+ const hp::MappingCollection<dim, spacedim> &mapping,
+ const DoFHandler<dim, spacedim> &dof,
+ const hp::QCollection<dim - 1> &quadrature,
+ const std::map<types::boundary_id, const Function<spacedim, Number> *>
+ &neumann_bc,
+ const ReadVector<Number> &solution,
+ Vector<float> &error,
+ const ComponentMask &component_mask = {},
+ const Function<spacedim> *coefficients = nullptr,
+ const unsigned int n_threads = numbers::invalid_unsigned_int,
+ const types::subdomain_id subdomain_id = numbers::invalid_subdomain_id,
+ const types::material_id material_id = numbers::invalid_material_id,
+ const Strategy strategy = cell_diameter_over_24);
+
+
+ /**
+ * Equivalent to the set of functions above, except that this one takes a
+ * quadrature collection for hp-finite element dof handlers.
+ *
+ * @deprecated Use the version of this function which takes a
+ * hp::MappingCollection instead.
+ */
+ template <typename Number>
+ DEAL_II_DEPRECATED_EARLY_WITH_COMMENT(
+ "Use the version of this function which takes a hp::MappingCollection instead.")
+ static void estimate(
const Mapping<dim, spacedim> &mapping,
const DoFHandler<dim, spacedim> &dof,
const hp::QCollection<dim - 1> &quadrature,
/**
* Equivalent to the set of functions above, except that this one takes a
- * quadrature collection for hp-finite element dof handlers.
+ * mapping and quadrature collection for hp-finite element dof handlers.
*/
template <typename Number>
static void
estimate(
+ const hp::MappingCollection<dim, spacedim> &mapping,
+ const DoFHandler<dim, spacedim> &dof,
+ const hp::QCollection<dim - 1> &quadrature,
+ const std::map<types::boundary_id, const Function<spacedim, Number> *>
+ &neumann_bc,
+ const ArrayView<const ReadVector<Number> *> &solutions,
+ ArrayView<Vector<float> *> &errors,
+ const ComponentMask &component_mask = {},
+ const Function<spacedim> *coefficients = nullptr,
+ const unsigned int n_threads = numbers::invalid_unsigned_int,
+ const types::subdomain_id subdomain_id = numbers::invalid_subdomain_id,
+ const types::material_id material_id = numbers::invalid_material_id,
+ const Strategy strategy = cell_diameter_over_24);
+
+
+ /**
+ * Equivalent to the set of functions above, except that this one takes a
+ * quadrature collection for hp-finite element dof handlers.
+ *
+ * @deprecated Use the version of this function which takes a
+ * hp::MappingCollection instead.
+ */
+ template <typename Number>
+ DEAL_II_DEPRECATED_EARLY_WITH_COMMENT(
+ "Use the version of this function which takes a hp::MappingCollection instead.")
+ static void estimate(
const Mapping<dim, spacedim> &mapping,
const DoFHandler<dim, spacedim> &dof,
const hp::QCollection<dim - 1> &quadrature,
/**
* Equivalent to the set of functions above, except that this one takes a
- * quadrature collection for hp-finite element dof handlers.
+ * mapping and quadrature collection for hp-finite element dof handlers.
*/
template <typename Number>
static void
estimate(
+ const hp::MappingCollection<1, spacedim> &mapping,
+ const DoFHandler<1, spacedim> &dof,
+ const hp::QCollection<0> &quadrature,
+ const std::map<types::boundary_id, const Function<spacedim, Number> *>
+ &neumann_bc,
+ const ReadVector<Number> &solution,
+ Vector<float> &error,
+ const ComponentMask &component_mask = {},
+ const Function<spacedim> *coefficients = nullptr,
+ const unsigned int n_threads = numbers::invalid_unsigned_int,
+ const types::subdomain_id subdomain_id = numbers::invalid_subdomain_id,
+ const types::material_id material_id = numbers::invalid_material_id,
+ const Strategy strategy = cell_diameter_over_24);
+
+
+
+ /**
+ * Equivalent to the set of functions above, except that this one takes a
+ * quadrature collection for hp-finite element dof handlers.
+ *
+ * @deprecated Use the version of this function which takes a
+ * hp::MappingCollection instead.
+ */
+ template <typename Number>
+ DEAL_II_DEPRECATED_EARLY_WITH_COMMENT(
+ "Use the version of this function which takes a hp::MappingCollection instead.")
+ static void estimate(
const Mapping<1, spacedim> &mapping,
const DoFHandler<1, spacedim> &dof,
const hp::QCollection<0> &quadrature,
/**
* Equivalent to the set of functions above, except that this one takes a
- * quadrature collection for hp-finite element dof handlers.
+ * mapping and quadrature collection for hp-finite element dof handlers.
*/
template <typename Number>
static void
estimate(
+ const hp::MappingCollection<1, spacedim> &mapping,
+ const DoFHandler<1, spacedim> &dof,
+ const hp::QCollection<0> &quadrature,
+ const std::map<types::boundary_id, const Function<spacedim, Number> *>
+ &neumann_bc,
+ const ArrayView<const ReadVector<Number> *> &solutions,
+ ArrayView<Vector<float> *> &errors,
+ const ComponentMask &component_mask = {},
+ const Function<spacedim> *coefficients = nullptr,
+ const unsigned int n_threads = numbers::invalid_unsigned_int,
+ const types::subdomain_id subdomain_id = numbers::invalid_subdomain_id,
+ const types::material_id material_id = numbers::invalid_material_id,
+ const Strategy strategy = cell_diameter_over_24);
+
+
+
+ /**
+ * Equivalent to the set of functions above, except that this one takes a
+ * quadrature collection for hp-finite element dof handlers.
+ *
+ * @deprecated Use the version of this function which takes a
+ * hp::MappingCollection instead.
+ */
+ template <typename Number>
+ DEAL_II_DEPRECATED_EARLY_WITH_COMMENT(
+ "Use the version of this function which takes a hp::MappingCollection instead.")
+ static void estimate(
const Mapping<1, spacedim> &mapping,
const DoFHandler<1, spacedim> &dof,
const hp::QCollection<0> &quadrature,
}
+template <int dim, int spacedim>
+template <typename Number>
+void
+KellyErrorEstimator<dim, spacedim>::estimate(
+ const hp::MappingCollection<dim, spacedim> &mapping,
+ const DoFHandler<dim, spacedim> &dof_handler,
+ const hp::QCollection<dim - 1> &quadrature,
+ const std::map<types::boundary_id, const Function<spacedim, Number> *>
+ &neumann_bc,
+ const ReadVector<Number> &solution,
+ Vector<float> &error,
+ const ComponentMask &component_mask,
+ const Function<spacedim> *coefficients,
+ const unsigned int n_threads,
+ const types::subdomain_id subdomain_id,
+ const types::material_id material_id,
+ const Strategy strategy)
+{
+ // just pass on to the other function
+ std::vector<const ReadVector<Number> *> solutions(1, &solution);
+ std::vector<Vector<float> *> errors(1, &error);
+ ArrayView<Vector<float> *> error_view = make_array_view(errors);
+ estimate(mapping,
+ dof_handler,
+ quadrature,
+ neumann_bc,
+ make_array_view(solutions),
+ error_view,
+ component_mask,
+ coefficients,
+ n_threads,
+ subdomain_id,
+ material_id,
+ strategy);
+}
+
+
template <int dim, int spacedim>
template <typename Number>
void
const types::material_id material_id,
const Strategy strategy)
{
+ // DEPRECATED
// just pass on to the other function
std::vector<const ReadVector<Number> *> solutions(1, &solution);
std::vector<Vector<float> *> errors(1, &error);
ArrayView<Vector<float> *> error_view = make_array_view(errors);
- estimate(mapping,
+ const hp::MappingCollection<dim, spacedim> mapping_collection(mapping);
+ estimate(mapping_collection,
dof_handler,
quadrature,
neumann_bc,
const types::material_id material_id,
const Strategy strategy)
{
- estimate(get_default_linear_mapping(dof_handler.get_triangulation()),
+ const hp::MappingCollection<dim, spacedim> mapping(
+ get_default_linear_mapping(dof_handler.get_triangulation()));
+ estimate(mapping,
dof_handler,
quadrature,
neumann_bc,
template <typename Number>
void
KellyErrorEstimator<dim, spacedim>::estimate(
- const Mapping<dim, spacedim> &mapping,
- const DoFHandler<dim, spacedim> &dof_handler,
- const hp::QCollection<dim - 1> &face_quadratures,
+ const hp::MappingCollection<dim, spacedim> &mapping,
+ const DoFHandler<dim, spacedim> &dof_handler,
+ const hp::QCollection<dim - 1> &face_quadratures,
const std::map<types::boundary_id, const Function<spacedim, Number> *>
&neumann_bc,
const ArrayView<const ReadVector<Number> *> &solutions,
// all the data needed in the error estimator by each of the threads is
// gathered in the following structures
- const hp::MappingCollection<dim, spacedim> mapping_collection(mapping);
const internal::ParallelData<dim, spacedim, Number> parallel_data(
dof_handler.get_fe_collection(),
face_quadratures,
- mapping_collection,
+ mapping,
(!neumann_bc.empty() || (coefficients != nullptr)),
solutions.size(),
subdomain_id,
+template <int dim, int spacedim>
+template <typename Number>
+void
+KellyErrorEstimator<dim, spacedim>::estimate(
+ const Mapping<dim, spacedim> &mapping,
+ const DoFHandler<dim, spacedim> &dof_handler,
+ const hp::QCollection<dim - 1> &quadrature,
+ const std::map<types::boundary_id, const Function<spacedim, Number> *>
+ &neumann_bc,
+ const ArrayView<const ReadVector<Number> *> &solutions,
+ ArrayView<Vector<float> *> &errors,
+ const ComponentMask &component_mask,
+ const Function<spacedim> *coefficients,
+ const unsigned int n_threads,
+ const types::subdomain_id subdomain_id,
+ const types::material_id material_id,
+ const Strategy strategy)
+{
+ // DEPRECATED
+ // just pass on to the other function
+ const hp::MappingCollection<dim, spacedim> mapping_collection(mapping);
+ estimate(mapping_collection,
+ dof_handler,
+ quadrature,
+ neumann_bc,
+ solutions,
+ errors,
+ component_mask,
+ coefficients,
+ n_threads,
+ subdomain_id,
+ material_id,
+ strategy);
+}
+
+
+
template <int dim, int spacedim>
template <typename Number>
void
const types::material_id material_id,
const Strategy strategy)
{
- // forward to the function with the QCollection
- estimate(mapping,
+ const hp::MappingCollection<dim, spacedim> mapping_collection(mapping);
+ const hp::QCollection<dim - 1> quadrature_collection(quadrature);
+ estimate(mapping_collection,
dof_handler,
- hp::QCollection<dim - 1>(quadrature),
+ quadrature_collection,
neumann_bc,
solutions,
errors,
const types::material_id material_id,
const Strategy strategy)
{
- estimate(get_default_linear_mapping(dof_handler.get_triangulation()),
+ const hp::MappingCollection<dim, spacedim> mapping(
+ get_default_linear_mapping(dof_handler.get_triangulation()));
+ estimate(mapping,
dof_handler,
quadrature,
neumann_bc,
const KellyErrorEstimator<deal_II_dimension,
deal_II_space_dimension>::Strategy);
+ template void
+ KellyErrorEstimator<deal_II_dimension, deal_II_space_dimension>::estimate<
+ S>(
+ const hp::MappingCollection<deal_II_dimension, deal_II_space_dimension> &,
+ const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+ const hp::QCollection<deal_II_dimension - 1> &,
+ const std::map<types::boundary_id,
+ const Function<deal_II_space_dimension, S> *> &,
+ const ReadVector<S> &,
+ Vector<float> &,
+ const ComponentMask &,
+ const Function<deal_II_space_dimension> *,
+ const unsigned int,
+ const types::subdomain_id,
+ const types::material_id,
+ const KellyErrorEstimator<deal_II_dimension,
+ deal_II_space_dimension>::Strategy);
+
+ // DEPRECATED
template void
KellyErrorEstimator<deal_II_dimension, deal_II_space_dimension>::estimate<
S>(const Mapping<deal_II_dimension, deal_II_space_dimension> &,
const KellyErrorEstimator<deal_II_dimension,
deal_II_space_dimension>::Strategy);
+ template void
+ KellyErrorEstimator<deal_II_dimension, deal_II_space_dimension>::estimate<
+ S>(
+ const hp::MappingCollection<deal_II_dimension, deal_II_space_dimension> &,
+ const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+ const hp::QCollection<deal_II_dimension - 1> &,
+ const std::map<types::boundary_id,
+ const Function<deal_II_space_dimension, S> *> &,
+ const ArrayView<const ReadVector<S> *> &,
+ ArrayView<Vector<float> *> &,
+ const ComponentMask &,
+ const Function<deal_II_space_dimension> *,
+ const unsigned int,
+ const types::subdomain_id,
+ const types::material_id,
+ const KellyErrorEstimator<deal_II_dimension,
+ deal_II_space_dimension>::Strategy);
+
+ // DEPRECATED
template void
KellyErrorEstimator<deal_II_dimension, deal_II_space_dimension>::estimate<
S>(const Mapping<deal_II_dimension, deal_II_space_dimension> &,
template <typename Number>
void
KellyErrorEstimator<1, spacedim>::estimate(
- const Mapping<1, spacedim> &mapping,
- const DoFHandler<1, spacedim> &dof_handler,
- const hp::QCollection<0> &quadrature,
+ const hp::MappingCollection<1, spacedim> &mapping,
+ const DoFHandler<1, spacedim> &dof_handler,
+ const hp::QCollection<0> &quadrature,
const std::map<types::boundary_id, const Function<spacedim, Number> *>
&neumann_bc,
const ReadVector<Number> &solution,
+template <int spacedim>
+template <typename Number>
+void
+KellyErrorEstimator<1, spacedim>::estimate(
+ const Mapping<1, spacedim> &mapping,
+ const DoFHandler<1, spacedim> &dof_handler,
+ const hp::QCollection<0> &quadrature,
+ const std::map<types::boundary_id, const Function<spacedim, Number> *>
+ &neumann_bc,
+ const ReadVector<Number> &solution,
+ Vector<float> &error,
+ const ComponentMask &component_mask,
+ const Function<spacedim> *coefficients,
+ const unsigned int n_threads,
+ const types::subdomain_id subdomain_id,
+ const types::material_id material_id,
+ const Strategy strategy)
+{
+ // DEPRECATED
+ // just pass on to the other function
+ std::vector<const ReadVector<Number> *> solutions(1, &solution);
+ std::vector<Vector<float> *> errors(1, &error);
+ ArrayView<Vector<float> *> error_view = make_array_view(errors);
+ const hp::MappingCollection<1, spacedim> mapping_collection(mapping);
+ estimate(mapping_collection,
+ dof_handler,
+ quadrature,
+ neumann_bc,
+ make_array_view(solutions),
+ error_view,
+ component_mask,
+ coefficients,
+ n_threads,
+ subdomain_id,
+ material_id,
+ strategy);
+}
+
+
+
template <int spacedim>
template <typename Number>
void
const Strategy strategy)
{
const auto reference_cell = ReferenceCells::Line;
- estimate(reference_cell.template get_default_linear_mapping<1, spacedim>(),
+ const hp::MappingCollection<1, spacedim> mapping(
+ reference_cell.template get_default_linear_mapping<1, spacedim>());
+ estimate(mapping,
dof_handler,
quadrature,
neumann_bc,
const Strategy strategy)
{
const auto reference_cell = ReferenceCells::Line;
- estimate(reference_cell.template get_default_linear_mapping<1, spacedim>(),
+ const hp::MappingCollection<1, spacedim> mapping(
+ reference_cell.template get_default_linear_mapping<1, spacedim>());
+ estimate(mapping,
dof_handler,
quadrature,
neumann_bc,
template <typename Number>
void
KellyErrorEstimator<1, spacedim>::estimate(
- const Mapping<1, spacedim> &mapping,
- const DoFHandler<1, spacedim> &dof_handler,
+ const hp::MappingCollection<1, spacedim> &mapping,
+ const DoFHandler<1, spacedim> &dof_handler,
const hp::QCollection<0> &,
const std::map<types::boundary_id, const Function<spacedim, Number> *>
&neumann_bc,
const hp::FECollection<1, spacedim> &fe = dof_handler.get_fe_collection();
- hp::MappingCollection<1, spacedim> mapping_collection;
- mapping_collection.push_back(mapping);
- hp::FEValues<1, spacedim> fe_values(mapping_collection,
+ hp::FEValues<1, spacedim> fe_values(mapping,
fe,
q_collection,
update_gradients);
hp::FEFaceValues<1, spacedim> fe_face_values(
- /*mapping_collection,*/ fe, q_face_collection, update_normal_vectors);
+ /*mapping,*/ fe, q_face_collection, update_normal_vectors);
// loop over all cells and do something on the cells which we're told to
// work on. note that the error indicator is only a sum over the two
+template <int spacedim>
+template <typename Number>
+void
+KellyErrorEstimator<1, spacedim>::estimate(
+ const Mapping<1, spacedim> &mapping,
+ const DoFHandler<1, spacedim> &dof_handler,
+ const hp::QCollection<0> &quadrature,
+ const std::map<types::boundary_id, const Function<spacedim, Number> *>
+ &neumann_bc,
+ const ArrayView<const ReadVector<Number> *> &solutions,
+ ArrayView<Vector<float> *> &errors,
+ const ComponentMask &component_mask,
+ const Function<spacedim> *coefficients,
+ const unsigned int n_threads,
+ const types::subdomain_id subdomain_id,
+ const types::material_id material_id,
+ const Strategy strategy)
+{
+ // DEPRECATED
+ // just pass on to the other function
+ const hp::MappingCollection<1, spacedim> mapping_collection(mapping);
+ estimate(mapping_collection,
+ dof_handler,
+ quadrature,
+ neumann_bc,
+ solutions,
+ errors,
+ component_mask,
+ coefficients,
+ n_threads,
+ subdomain_id,
+ material_id,
+ strategy);
+}
+
+
+
template <int spacedim>
template <typename Number>
void
const types::material_id material_id,
const Strategy strategy)
{
- const hp::QCollection<0> quadrature_collection(quadrature);
- estimate(mapping,
+ const hp::MappingCollection<1, spacedim> mapping_collection(mapping);
+ const hp::QCollection<0> quadrature_collection(quadrature);
+ estimate(mapping_collection,
dof_handler,
quadrature_collection,
neumann_bc,
const types::material_id,
const Strategy);
+ template void
+ KellyErrorEstimator<deal_II_dimension, deal_II_space_dimension>::estimate<
+ S>(
+ const hp::MappingCollection<deal_II_dimension, deal_II_space_dimension> &,
+ const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+ const hp::QCollection<deal_II_dimension - 1> &,
+ const std::map<types::boundary_id,
+ const Function<deal_II_space_dimension, S> *> &,
+ const ReadVector<S> &,
+ Vector<float> &,
+ const ComponentMask &,
+ const Function<deal_II_space_dimension> *,
+ const unsigned int,
+ const types::subdomain_id,
+ const types::material_id,
+ const Strategy);
+
+ // DEPRECATED
template void
KellyErrorEstimator<deal_II_dimension, deal_II_space_dimension>::estimate<
S>(const Mapping<deal_II_dimension, deal_II_space_dimension> &,
const types::material_id,
const Strategy);
+ template void
+ KellyErrorEstimator<deal_II_dimension, deal_II_space_dimension>::estimate<
+ S>(
+ const hp::MappingCollection<deal_II_dimension, deal_II_space_dimension> &,
+ const DoFHandler<deal_II_dimension, deal_II_space_dimension> &,
+ const hp::QCollection<deal_II_dimension - 1> &,
+ const std::map<types::boundary_id,
+ const Function<deal_II_space_dimension, S> *> &,
+ const ArrayView<const ReadVector<S> *> &,
+ ArrayView<Vector<float> *> &,
+ const ComponentMask &,
+ const Function<deal_II_space_dimension> *,
+ const unsigned int,
+ const types::subdomain_id,
+ const types::material_id,
+ const Strategy);
+
+ // DEPRECATED
template void
KellyErrorEstimator<deal_II_dimension, deal_II_space_dimension>::estimate<
S>(const Mapping<deal_II_dimension, deal_II_space_dimension> &,
--- /dev/null
+// ------------------------------------------------------------------------
+//
+// SPDX-License-Identifier: LGPL-2.1-or-later
+// Copyright (C) 2025 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// Part of the source code is dual licensed under Apache-2.0 WITH
+// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
+// governing the source code and code contributions can be found in
+// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
+//
+// ------------------------------------------------------------------------
+
+
+// Test that we pick the correct mapping in KellyErrorEstimator
+// for hp applications.
+
+
+#include <deal.II/base/function_lib.h>
+#include <deal.II/base/quadrature_lib.h>
+
+#include <deal.II/dofs/dof_handler.h>
+
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/fe/mapping_q.h>
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/tria.h>
+
+#include <deal.II/hp/fe_collection.h>
+#include <deal.II/hp/mapping_collection.h>
+#include <deal.II/hp/q_collection.h>
+
+#include <deal.II/lac/vector.h>
+
+#include <deal.II/numerics/error_estimator.h>
+#include <deal.II/numerics/vector_tools.h>
+
+#include "../tests.h"
+
+
+template <int dim>
+void
+test()
+{
+ // We investigate the effect of different mappings in an hp-context.
+ // Thus we only feed the MappingCollection with different mappings,
+ // and keep the FEs and quadratures the same.
+ hp::FECollection<dim> fes(FE_Q<dim>(1), FE_Q<dim>(1));
+ hp::QCollection<dim - 1> quads_face(QGauss<dim - 1>(2), QGauss<dim - 1>(2));
+ hp::MappingCollection<dim> mappings(MappingQ<dim>(1), MappingQ<dim>(2));
+
+ // Set up hyper_ball and assign different active FE indices
+ // to inner and outer cells.
+ Triangulation<dim> tria;
+ DoFHandler<dim> dofh(tria);
+
+ GridGenerator::hyper_ball(tria);
+ for (const auto &cell : dofh.active_cell_iterators())
+ if (cell->center() != Point<dim>())
+ cell->set_active_fe_index(1);
+ tria.refine_global(1);
+
+ dofh.distribute_dofs(fes);
+
+ // Interpolate a transcendental function and estimate the error via Kelly
+ // with different mappings.
+ Functions::CosineFunction<dim> function;
+ Vector<double> interpolation(dofh.n_dofs());
+
+ Vector<float> error(tria.n_active_cells());
+
+ // MappingQ1
+ hp::MappingCollection<dim> mapping_q1(mappings[0]);
+ VectorTools::interpolate(mapping_q1, dofh, function, interpolation);
+ KellyErrorEstimator<dim>::estimate(
+ mapping_q1, dofh, quads_face, {}, interpolation, error);
+ deallog << "MappingQ1 : global error estimate=" << error.l2_norm()
+ << std::endl;
+
+ // MappingQ2
+ hp::MappingCollection<dim> mapping_q2(mappings[1]);
+ VectorTools::interpolate(mapping_q2, dofh, function, interpolation);
+ KellyErrorEstimator<dim>::estimate(
+ mapping_q2, dofh, quads_face, {}, interpolation, error);
+ deallog << "MappingQ2 : global error estimate=" << error.l2_norm()
+ << std::endl;
+
+ // MappingQ1&2
+ VectorTools::interpolate(mappings, dofh, function, interpolation);
+ KellyErrorEstimator<dim>::estimate(
+ mappings, dofh, quads_face, {}, interpolation, error);
+ deallog << "MappingQ1&2: global error estimate=" << error.l2_norm()
+ << std::endl;
+}
+
+
+int
+main()
+{
+ initlog();
+
+ deallog.push("2d");
+ test<2>();
+ deallog.pop();
+ deallog.push("3d");
+ test<3>();
+ deallog.pop();
+}
--- /dev/null
+
+DEAL:2d::MappingQ1 : global error estimate=0.585534
+DEAL:2d::MappingQ2 : global error estimate=0.465513
+DEAL:2d::MappingQ1&2: global error estimate=0.465513
+DEAL:3d::MappingQ1 : global error estimate=0.922872
+DEAL:3d::MappingQ2 : global error estimate=0.603987
+DEAL:3d::MappingQ1&2: global error estimate=0.603987