VectorTools::interpolate fails for FE_FaceQ

[cpraveen]

With deal.II version c99312d, VectorTools::interpolate fails for FE_FaceQ.

Here is an example code

#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;

int main()
{
   Triangulation<2>    triangulation;
   FE_FaceQ<2>         fe(2);
   DoFHandler<2>       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<2>(),
                            solution);
}

and the output

--------------------------------------------------------
An error occurred in line <1157> of file </home/praveen/Applications/deal.II/git/source/fe/fe.cc> in function
    void dealii::FiniteElement<<anonymous>, <anonymous> >::convert_generalized_support_point_values_to_dof_values(const std::vector<dealii::Vector<double> >&, std::vector<double>&) const [with int dim = 2; int spacedim = 2]
The violated condition was: 
    false
Additional information: 
    You are trying to use functionality in deal.II that is currently not implemented. In many cases, this indicates that there simply didn't appear much of a need for it, or that the author of the original code did not have the time to implement a particular case. If you hit this exception, it is therefore worth the time to look into the code to find out whether you may be able to implement the missing functionality. If you do, please consider providing a patch to the deal.II development sources (see the deal.II website on how to contribute).

Stacktrace:
-----------
#0  /opt/deal.II/lib/libdeal_II.g.so.9.0.0-pre: dealii::FiniteElement<2, 2>::convert_generalized_support_point_values_to_dof_values(std::vector<dealii::Vector<double>, std::allocator<dealii::Vector<double> > > const&, std::vector<double, std::allocator<double> >&) const
#1  /opt/deal.II/lib/libdeal_II.g.so.9.0.0-pre: 
#2  /opt/deal.II/lib/libdeal_II.g.so.9.0.0-pre: void dealii::FETools::convert_generalized_support_point_values_to_dof_values<2, 2, double>(dealii::FiniteElement<2, 2> const&, std::vector<dealii::Vector<double>, std::allocator<dealii::Vector<double> > > const&, std::vector<double, std::allocator<double> >&)
#3  /opt/deal.II/lib/libdeal_II.g.so.9.0.0-pre: 
#4  /opt/deal.II/lib/libdeal_II.g.so.9.0.0-pre: void dealii::VectorTools::interpolate<2, 2, dealii::Vector<double>, dealii::DoFHandler>(dealii::Mapping<2, 2> const&, dealii::DoFHandler<2, 2> const&, dealii::Function<2, dealii::Vector<double>::value_type> const&, dealii::Vector<double>&, dealii::ComponentMask const&)
#5  /opt/deal.II/lib/libdeal_II.g.so.9.0.0-pre: void dealii::VectorTools::interpolate<2, 2, dealii::Vector<double>, dealii::DoFHandler>(dealii::DoFHandler<2, 2> const&, dealii::Function<2, dealii::Vector<double>::value_type> const&, dealii::Vector<double>&, dealii::ComponentMask const&)
#6  ./main: main
--------------------------------------------------------

[bangerth]

I can confirm this.

FE_FaceQ is an interpolating element, i.e., the function FE_FaceQ::convert_generalized_support_point_values_to_dof_values that is not implemented should -- as far as I can see -- be trivial and simply copy its input array to its output array. Can you try what happens if you just steal the corresponding function from one of the other interpolating elements -- say FE_Q or the base class of that class that implements this function?

[bangerth]

As a matter of fact, I think the following (adjusted for the different class name) should work:

template <int dim, int spacedim>
void
FE_DGQ<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);
    }
}

[cpraveen]

I added the above function to FE_FaceQ and now interpolate works. Similar function is needed in FE_TraceQ also.

[masterleinad]

Great! Would you be able to turn this into a pull request with a corresponding test?

[bangerth]

Indeed, that would be great!

[masterleinad]

Fixed by #5503.