--- /dev/null
+// ---------------------------------------------------------------------
+//
+// 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.
+//
+// ---------------------------------------------------------------------
+
+// This header contains a non-trivial testfunction F and a small test
+// procedure shared among all generalized_interpolation_* tests.
+
+#include <deal.II/dofs/dof_handler.h>
+#include <deal.II/fe/fe.h>
+#include <deal.II/fe/mapping_q.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_tools.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/lac/vector.h>
+#include <deal.II/numerics/vector_tools.h>
+#include <deal.II/numerics/fe_field_function.h>
+
+template <int dim>
+class F : public Function<dim>
+{
+public:
+ F(const unsigned int n_comp, const unsigned int q)
+ : Function<dim>(n_comp), q(q)
+ {
+ }
+
+ virtual double value(const Point<dim> &p, const unsigned int) const
+ {
+ double v = 0;
+ for (unsigned int d = 0; d < dim; ++d)
+ for (unsigned int i = 0; i <= q; ++i)
+ v += (d + 1) * (i + 1) * std::pow(p[d], 1. * i);
+ return v;
+ }
+
+private:
+ const unsigned int q;
+};
+
+template <int dim, typename T>
+void test(const FiniteElement<dim> &fe,
+ const T &f,
+ const unsigned int order_mapping,
+ bool distort_mesh,
+ bool print_function_values = false)
+{
+ deallog << "dim " << dim << " " << fe.get_name() << std::endl;
+
+ Triangulation<dim> triangulation;
+ GridGenerator::hyper_cube(triangulation, -0.3, 0.7);
+ triangulation.refine_global(dim == 2 ? 2 : 1);
+ if (distort_mesh)
+ GridTools::distort_random(0.03, triangulation);
+
+ MappingQ<dim> mapping(order_mapping);
+
+ DoFHandler<dim> dof_handler(triangulation);
+ dof_handler.distribute_dofs(fe);
+
+ Vector<double> interpolant(dof_handler.n_dofs());
+ VectorTools::interpolate(mapping, dof_handler, f, interpolant);
+
+ // Print function value in origin:
+ if (print_function_values)
+ {
+ Functions::FEFieldFunction<dim> f2(dof_handler, interpolant, mapping);
+ deallog << "Function value at (0.0,0.0): ";
+ for (unsigned int i = 0; i < fe.n_components(); ++i)
+ deallog << f2.value(Point<dim>(), i) << " ";
+ deallog << std::endl;
+ }
+
+ // Check that VectorTools::interpolate is in fact a
+ // projection, i.e. applying the interpolation twice results in the same
+ // vector:
+
+ Functions::FEFieldFunction<dim> f2(dof_handler, interpolant, mapping);
+
+ Vector<double> interpolant2(dof_handler.n_dofs());
+ VectorTools::interpolate(
+ mapping, dof_handler, f2, interpolant2);
+
+ interpolant2 -= interpolant;
+ deallog << "Check projection property: " << interpolant2.linfty_norm()
+ << std::endl;
+}
// Check that VectorTools::interpolate correctly recovers a
// constant vector field for H1, Hdiv and Hcurl conforming elements.
-#include <deal.II/base/quadrature_lib.h>
-#include <deal.II/dofs/dof_handler.h>
#include <deal.II/fe/fe_nedelec.h>
#include <deal.II/fe/fe_q.h>
#include <deal.II/fe/fe_raviart_thomas.h>
-#include <deal.II/fe/mapping_q.h>
-#include <deal.II/grid/grid_generator.h>
-#include <deal.II/grid/grid_tools.h>
-#include <deal.II/grid/tria.h>
-#include <deal.II/lac/vector.h>
-#include <deal.II/numerics/vector_tools.h>
-#include <deal.II/numerics/fe_field_function.h>
+#include <deal.II/fe/fe_system.h>
#include "../tests.h"
-
-template <int dim>
-void test(const FiniteElement<dim> &fe,
- const unsigned int n_comp,
- const unsigned int order_mapping,
- bool distort_mesh)
-{
- deallog << "dim " << dim << " " << fe.get_name() << std::endl;
-
- Triangulation<dim> triangulation;
- GridGenerator::hyper_cube(triangulation, -0.3, 0.7);
- triangulation.refine_global(dim == 2 ? 2 : 1);
- if (distort_mesh)
- GridTools::distort_random(0.03, triangulation);
-
- ConstantFunction<dim> f (1.0, n_comp);
-
- MappingQ<dim> mapping(order_mapping);
-
- DoFHandler<dim> dof_handler(triangulation);
- dof_handler.distribute_dofs(fe);
-
- Vector<double> interpolant(dof_handler.n_dofs());
- VectorTools::interpolate(mapping, dof_handler, f, interpolant);
-
- // Check that the interpoland really returns 1.0...
-
- Functions::FEFieldFunction<dim> f2(dof_handler, interpolant, mapping);
- deallog << "Function value at (0.0,0.0): ";
- for (unsigned int i = 0; i < n_comp; ++i)
- deallog << f2.value(Point<dim>(), i) << " ";
- deallog << std::endl;
-
- // Check that VectorTools::interpolate is in fact a
- // projection, i.e. applying the interpolation twice results in the same
- // vector:
-
- Vector<double> interpolant2(dof_handler.n_dofs());
- VectorTools::interpolate(
- mapping, dof_handler, f2, interpolant2);
-
- interpolant2 -= interpolant;
- deallog << "Check projection property: " << interpolant2.linfty_norm()
- << std::endl;
-}
-
+#include "generalized_interpolation.h"
int main ()
{
- deallog.depth_console(3);
+ initlog();
- test<2>(FE_Q<2>(1), 1, 1, false);
- test<2>(FE_Q<2>(2), 1, 2, false);
- test<3>(FE_Q<3>(3), 1, 3, false);
+ test<2>(FE_Q<2>(1), ConstantFunction<2>(1.0, 1), 1, false, true);
+ test<2>(FE_Q<2>(2), ConstantFunction<2>(1.0, 1), 2, false, true);
+ test<3>(FE_Q<3>(3), ConstantFunction<3>(1.0, 1), 3, false, true);
- test<2>(FE_Q<2>(1), 1, 1, true);
- test<2>(FE_Q<2>(2), 1, 2, true);
- test<3>(FE_Q<3>(3), 1, 3, true);
+ test<2>(FE_Q<2>(1), ConstantFunction<2>(1.0, 1), 1, true, true);
+ test<2>(FE_Q<2>(2), ConstantFunction<2>(1.0, 1), 2, true, true);
+ test<3>(FE_Q<3>(3), ConstantFunction<3>(1.0, 1), 3, true, true);
- test<2>(FE_RaviartThomas<2>(0), 2, 1, false);
- test<2>(FE_RaviartThomas<2>(1), 2, 2, false);
- test<2>(FE_RaviartThomas<2>(2), 2, 3, false);
- test<3>(FE_RaviartThomas<3>(0), 3, 1, false);
- test<3>(FE_RaviartThomas<3>(1), 3, 2, false);
+ test<2>(FE_RaviartThomas<2>(0), ConstantFunction<2>(1.0, 2), 1, false, true);
+ test<2>(FE_RaviartThomas<2>(1), ConstantFunction<2>(1.0, 2), 2, false, true);
+ test<2>(FE_RaviartThomas<2>(2), ConstantFunction<2>(1.0, 2), 3, false, true);
+ test<3>(FE_RaviartThomas<3>(0), ConstantFunction<3>(1.0, 3), 1, false, true);
+ test<3>(FE_RaviartThomas<3>(1), ConstantFunction<3>(1.0, 3), 2, false, true);
- test<2>(FE_RaviartThomas<2>(0), 2, 1, true);
- test<2>(FE_RaviartThomas<2>(1), 2, 2, true);
- test<2>(FE_RaviartThomas<2>(2), 2, 3, true);
+ test<2>(FE_RaviartThomas<2>(0), ConstantFunction<2>(1.0, 2), 1, true, true);
+ test<2>(FE_RaviartThomas<2>(1), ConstantFunction<2>(1.0, 2), 2, true, true);
+ test<2>(FE_RaviartThomas<2>(2), ConstantFunction<2>(1.0, 2), 3, true, true);
// lowest order RT in 3D does not contain constant 1 function on a
// distorted mesh.
- test<3>(FE_RaviartThomas<3>(1), 3, 2, true);
-
- test<2>(FE_Nedelec<2>(0), 2, 1, false);
- test<2>(FE_Nedelec<2>(1), 2, 2, false);
- test<2>(FE_Nedelec<2>(2), 2, 3, false);
- test<2>(FE_Nedelec<2>(0), 2, 1, true);
- test<2>(FE_Nedelec<2>(1), 2, 2, true);
- test<2>(FE_Nedelec<2>(2), 2, 3, true);
+ test<3>(FE_RaviartThomas<3>(1), ConstantFunction<3>(1.0, 3), 2, true, true);
+
+ test<2>(FE_Nedelec<2>(0), ConstantFunction<2>(1.0, 2), 1, false, true);
+ test<2>(FE_Nedelec<2>(1), ConstantFunction<2>(1.0, 2), 2, false, true);
+ test<2>(FE_Nedelec<2>(2), ConstantFunction<2>(1.0, 2), 3, false, true);
+ test<2>(FE_Nedelec<2>(0), ConstantFunction<2>(1.0, 2), 1, true, true);
+ test<2>(FE_Nedelec<2>(1), ConstantFunction<2>(1.0, 2), 2, true, true);
+ test<2>(FE_Nedelec<2>(2), ConstantFunction<2>(1.0, 2), 3, true, true);
}
+
DEAL::dim 2 FE_Q<2>(1)
DEAL::Function value at (0.0,0.0): 1.00000
DEAL::Check projection property: 0.00000
// ---------------------------------------------------------------------
// Check projection property of VectorTools::interpolate for
-// Hdiv conforming spaces on something nontrivial.
+// Hdiv and Hcurl conforming spaces on something nontrivial.
-#include <deal.II/base/quadrature_lib.h>
-#include <deal.II/dofs/dof_handler.h>
#include <deal.II/fe/fe_nedelec.h>
-#include <deal.II/fe/fe_q.h>
#include <deal.II/fe/fe_raviart_thomas.h>
-#include <deal.II/fe/mapping_q.h>
-#include <deal.II/grid/grid_generator.h>
-#include <deal.II/grid/grid_tools.h>
-#include <deal.II/grid/tria.h>
-#include <deal.II/lac/vector.h>
-#include <deal.II/numerics/vector_tools.h>
-#include <deal.II/numerics/fe_field_function.h>
#include "../tests.h"
-
-template <int dim>
-class F : public Function<dim>
-{
-public:
- F(const unsigned int q) : Function<dim>(dim), q(q) {}
-
- virtual double value(const Point<dim> &p, const unsigned int) const
- {
- double v = 0;
- for (unsigned int d = 0; d < dim; ++d)
- for (unsigned int i = 0; i <= q; ++i)
- v += (d + 1) * (i + 1) * std::pow(p[d], 1. * i);
- return v;
- }
-
-private:
- const unsigned int q;
-};
-
-template <int dim, typename T>
-void test(const FiniteElement<dim> &fe,
- const T &f,
- const unsigned int n_comp,
- const unsigned int order_mapping,
- bool distort_mesh)
-{
- deallog << "dim " << dim << " " << fe.get_name() << std::endl;
-
- Triangulation<dim> triangulation;
- GridGenerator::hyper_cube(triangulation, -0.3, 0.7);
- triangulation.refine_global(dim == 2 ? 2 : 1);
- if (distort_mesh)
- GridTools::distort_random(0.03, triangulation);
-
- MappingQ<dim> mapping(order_mapping);
-
- DoFHandler<dim> dof_handler(triangulation);
- dof_handler.distribute_dofs(fe);
-
- Vector<double> interpolant(dof_handler.n_dofs());
- VectorTools::interpolate(mapping, dof_handler, f, interpolant);
-
- // Check that VectorTools::interpolate is in fact a
- // projection, i.e. applying the interpolation twice results in the same
- // vector:
-
- Functions::FEFieldFunction<dim> f2(dof_handler, interpolant, mapping);
-
- Vector<double> interpolant2(dof_handler.n_dofs());
- VectorTools::interpolate(
- mapping, dof_handler, f2, interpolant2);
-
- interpolant2 -= interpolant;
- deallog << "Check projection property: " << interpolant2.linfty_norm()
- << std::endl;
-}
-
+#include "generalized_interpolation.h"
int main ()
{
deallog.depth_console(3);
- test<2>(FE_RaviartThomas<2>(0), F<2>(1), 2, 1, false);
- test<2>(FE_RaviartThomas<2>(1), F<2>(0), 2, 2, false);
- test<2>(FE_RaviartThomas<2>(1), F<2>(2), 2, 2, false);
+ test<2>(FE_RaviartThomas<2>(0), F<2>(2, 1), 1, false);
+ test<2>(FE_RaviartThomas<2>(1), F<2>(2, 0), 2, false);
+ test<2>(FE_RaviartThomas<2>(1), F<2>(2, 2), 2, false);
- test<3>(FE_RaviartThomas<3>(0), F<3>(0), 3, 1, false);
- test<3>(FE_RaviartThomas<3>(1), F<3>(0), 3, 2, false);
- test<3>(FE_RaviartThomas<3>(1), F<3>(2), 3, 2, false);
+ test<3>(FE_RaviartThomas<3>(0), F<3>(3, 0), 1, false);
+ test<3>(FE_RaviartThomas<3>(1), F<3>(3, 0), 2, false);
+ test<3>(FE_RaviartThomas<3>(1), F<3>(3, 2), 2, false);
}
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// 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 <deal.II/fe/fe_nedelec.h>
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/fe/fe_raviart_thomas.h>
+#include <deal.II/fe/fe_system.h>
+
+#include "../tests.h"
+#include "generalized_interpolation.h"
+
+int main ()
+{
+ initlog();
+
+ static const int dim = 2;
+
+ FESystem<dim> fe(FE_RaviartThomas<dim>(1),
+ 1,
+ FE_Q<dim>(1),
+ 1,
+ FE_Nedelec<dim>(2),
+ 2,
+ FESystem<dim>(FE_Q<dim>(1), dim),
+ 1);
+
+ const unsigned int n_comp = fe.n_components();
+
+ test<2>(fe, F<2>(n_comp, 1), 1, false);
+}
--- /dev/null
+
+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