*/
virtual std::string get_name () const;
+ /**
+ * Implementation of the corresponding function in the FiniteElement
+ * class. Since the current element is interpolatory, the nodal
+ * values are exactly the support point values. Furthermore, since
+ * the current element is scalar, the support point values need to
+ * be vectors of length 1.
+ */
+ virtual
+ void
+ convert_generalized_support_point_values_to_dof_values (const std::vector<Vector<double> > &support_point_values,
+ std::vector<double> &nodal_values) const;
+
/**
* Return the matrix interpolating from a face of of one element to the face
* of the neighboring element. The size of the matrix is then
std::unique_ptr<FiniteElement<dim,spacedim> >
clone() const;
+ /**
+ * Implementation of the corresponding function in the FiniteElement
+ * class. Since the current element is interpolatory, the nodal
+ * values are exactly the support point values. Furthermore, since
+ * the current element is scalar, the support point values need to
+ * be vectors of length 1.
+ */
+ virtual
+ void
+ convert_generalized_support_point_values_to_dof_values (const std::vector<Vector<double> > &support_point_values,
+ std::vector<double> &nodal_values) const;
+
/**
* This function returns @p true, if the shape function @p shape_index has
* non-zero function values somewhere on the face @p face_index.
(constant_modes, std::vector<unsigned int>(1, 0));
}
+template <int dim, int spacedim>
+void
+FE_FaceQ<dim,spacedim>::
+convert_generalized_support_point_values_to_dof_values (const std::vector<Vector<double> > &support_point_values,
+ std::vector<double> &nodal_values) const
+{
+ AssertDimension (support_point_values.size(),
+ this->get_unit_support_points().size());
+ AssertDimension (support_point_values.size(),
+ nodal_values.size());
+ AssertDimension (this->dofs_per_cell,
+ nodal_values.size());
+
+ for (unsigned int i=0; i<this->dofs_per_cell; ++i)
+ {
+ AssertDimension (support_point_values[i].size(), 1);
+ nodal_values[i] = support_point_values[i](0);
+ }
+}
// ----------------------------- FE_FaceQ<1,spacedim> ------------------------
(constant_modes, std::vector<unsigned int>(1, 0));
}
+template <int dim, int spacedim>
+void
+FE_TraceQ<dim,spacedim>::
+convert_generalized_support_point_values_to_dof_values (const std::vector<Vector<double> > &support_point_values,
+ std::vector<double> &nodal_values) const
+{
+ AssertDimension (support_point_values.size(),
+ this->get_unit_support_points().size());
+ AssertDimension (support_point_values.size(),
+ nodal_values.size());
+ AssertDimension (this->dofs_per_cell,
+ nodal_values.size());
+
+ for (unsigned int i=0; i<this->dofs_per_cell; ++i)
+ {
+ AssertDimension (support_point_values[i].size(), 1);
+
+ nodal_values[i] = support_point_values[i](0);
+ }
+}
template <int dim, int spacedim>
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2015 - 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.
+//
+// ---------------------------------------------------------------------
+
+// Call VectorTools::interpolate for a function in FE_TraceQ. The purpose is
+// to check that all functions needed for interpolate to work exist.
+
+#include "../tests.h"
+
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/fe/fe_face.h>
+#include <deal.II/dofs/dof_handler.h>
+#include <deal.II/numerics/vector_tools.h>
+
+using namespace dealii;
+
+template <int dim>
+void test()
+{
+ Triangulation<dim> triangulation;
+ FE_FaceQ<dim> fe(2);
+ DoFHandler<dim> dof_handler(triangulation);
+
+ GridGenerator::hyper_cube (triangulation, 0, 1);
+ triangulation.refine_global(6);
+
+ dof_handler.distribute_dofs (fe);
+ Vector<double> solution(dof_handler.n_dofs());
+
+ VectorTools::interpolate(dof_handler,
+ ZeroFunction<dim>(),
+ solution);
+ deallog << "Success, dim = " << dim << std::endl;
+}
+
+int main()
+{
+ initlog();
+
+ test<2>();
+ test<3>();
+ return 0;
+}
--- /dev/null
+
+DEAL::Success, dim = 2
+DEAL::Success, dim = 3
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2015 - 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.
+//
+// ---------------------------------------------------------------------
+
+// Call VectorTools::interpolate for a function in FE_TraceQ. The purpose is
+// to check that all functions needed for interpolate to work exist.
+
+#include "../tests.h"
+
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/fe/fe_trace.h>
+#include <deal.II/dofs/dof_handler.h>
+#include <deal.II/numerics/vector_tools.h>
+
+using namespace dealii;
+
+template <int dim>
+void test()
+{
+ Triangulation<dim> triangulation;
+ FE_TraceQ<dim> fe(2);
+ DoFHandler<dim> dof_handler(triangulation);
+
+ GridGenerator::hyper_cube (triangulation, 0, 1);
+ triangulation.refine_global(6);
+
+ dof_handler.distribute_dofs (fe);
+ Vector<double> solution(dof_handler.n_dofs());
+
+ VectorTools::interpolate(dof_handler,
+ ZeroFunction<dim>(),
+ solution);
+ deallog << "Success, dim = " << dim << std::endl;
+}
+
+int main()
+{
+ initlog();
+
+ test<2>();
+ test<3>();
+ return 0;
+}
--- /dev/null
+
+DEAL::Success, dim = 2
+DEAL::Success, dim = 3