From 8c4accf14ce070e106460eb9be7ec78defd4f5ad Mon Sep 17 00:00:00 2001 From: Praveen C Date: Sun, 19 Nov 2017 10:14:46 +0530 Subject: [PATCH] Add convert_generalized_support_point_values_to_dof_values This function is missing for FE_FaceQ and FE_TraceQ which is now implemented by copying from FE_DGQ. See https://github.com/dealii/dealii/issues/5449 --- include/deal.II/fe/fe_face.h | 12 +++++ include/deal.II/fe/fe_trace.h | 12 +++++ source/fe/fe_face.cc | 19 +++++++ source/fe/fe_trace.cc | 20 ++++++++ tests/fe/fe_faceq_interpolate.cc | 55 +++++++++++++++++++++ tests/fe/fe_faceq_interpolate.debug.output | 3 ++ tests/fe/fe_traceq_interpolate.cc | 55 +++++++++++++++++++++ tests/fe/fe_traceq_interpolate.debug.output | 3 ++ 8 files changed, 179 insertions(+) create mode 100644 tests/fe/fe_faceq_interpolate.cc create mode 100644 tests/fe/fe_faceq_interpolate.debug.output create mode 100644 tests/fe/fe_traceq_interpolate.cc create mode 100644 tests/fe/fe_traceq_interpolate.debug.output diff --git a/include/deal.II/fe/fe_face.h b/include/deal.II/fe/fe_face.h index 98d45b1649..b795fc1fcb 100644 --- a/include/deal.II/fe/fe_face.h +++ b/include/deal.II/fe/fe_face.h @@ -70,6 +70,18 @@ public: */ 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 > &support_point_values, + std::vector &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 diff --git a/include/deal.II/fe/fe_trace.h b/include/deal.II/fe/fe_trace.h index f4cc223760..fbadc0ffa4 100644 --- a/include/deal.II/fe/fe_trace.h +++ b/include/deal.II/fe/fe_trace.h @@ -62,6 +62,18 @@ public: std::unique_ptr > 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 > &support_point_values, + std::vector &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. diff --git a/source/fe/fe_face.cc b/source/fe/fe_face.cc index 448d2b3b6c..3c1cb96a29 100644 --- a/source/fe/fe_face.cc +++ b/source/fe/fe_face.cc @@ -315,7 +315,26 @@ FE_FaceQ::get_constant_modes () const (constant_modes, std::vector(1, 0)); } +template +void +FE_FaceQ:: +convert_generalized_support_point_values_to_dof_values (const std::vector > &support_point_values, + std::vector &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; idofs_per_cell; ++i) + { + AssertDimension (support_point_values[i].size(), 1); + nodal_values[i] = support_point_values[i](0); + } +} // ----------------------------- FE_FaceQ<1,spacedim> ------------------------ diff --git a/source/fe/fe_trace.cc b/source/fe/fe_trace.cc index dee1e7f869..4bc1bd6594 100644 --- a/source/fe/fe_trace.cc +++ b/source/fe/fe_trace.cc @@ -122,6 +122,26 @@ FE_TraceQ::get_constant_modes () const (constant_modes, std::vector(1, 0)); } +template +void +FE_TraceQ:: +convert_generalized_support_point_values_to_dof_values (const std::vector > &support_point_values, + std::vector &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; idofs_per_cell; ++i) + { + AssertDimension (support_point_values[i].size(), 1); + + nodal_values[i] = support_point_values[i](0); + } +} template diff --git a/tests/fe/fe_faceq_interpolate.cc b/tests/fe/fe_faceq_interpolate.cc new file mode 100644 index 0000000000..2f795a634a --- /dev/null +++ b/tests/fe/fe_faceq_interpolate.cc @@ -0,0 +1,55 @@ +// --------------------------------------------------------------------- +// +// 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 +#include +#include +#include +#include + +using namespace dealii; + +template +void test() +{ + Triangulation triangulation; + FE_FaceQ fe(2); + DoFHandler dof_handler(triangulation); + + GridGenerator::hyper_cube (triangulation, 0, 1); + triangulation.refine_global(6); + + dof_handler.distribute_dofs (fe); + Vector solution(dof_handler.n_dofs()); + + VectorTools::interpolate(dof_handler, + ZeroFunction(), + solution); + deallog << "Success, dim = " << dim << std::endl; +} + +int main() +{ + initlog(); + + test<2>(); + test<3>(); + return 0; +} diff --git a/tests/fe/fe_faceq_interpolate.debug.output b/tests/fe/fe_faceq_interpolate.debug.output new file mode 100644 index 0000000000..472a129ab0 --- /dev/null +++ b/tests/fe/fe_faceq_interpolate.debug.output @@ -0,0 +1,3 @@ + +DEAL::Success, dim = 2 +DEAL::Success, dim = 3 diff --git a/tests/fe/fe_traceq_interpolate.cc b/tests/fe/fe_traceq_interpolate.cc new file mode 100644 index 0000000000..bfca5bfaa7 --- /dev/null +++ b/tests/fe/fe_traceq_interpolate.cc @@ -0,0 +1,55 @@ +// --------------------------------------------------------------------- +// +// 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 +#include +#include +#include +#include + +using namespace dealii; + +template +void test() +{ + Triangulation triangulation; + FE_TraceQ fe(2); + DoFHandler dof_handler(triangulation); + + GridGenerator::hyper_cube (triangulation, 0, 1); + triangulation.refine_global(6); + + dof_handler.distribute_dofs (fe); + Vector solution(dof_handler.n_dofs()); + + VectorTools::interpolate(dof_handler, + ZeroFunction(), + solution); + deallog << "Success, dim = " << dim << std::endl; +} + +int main() +{ + initlog(); + + test<2>(); + test<3>(); + return 0; +} diff --git a/tests/fe/fe_traceq_interpolate.debug.output b/tests/fe/fe_traceq_interpolate.debug.output new file mode 100644 index 0000000000..472a129ab0 --- /dev/null +++ b/tests/fe/fe_traceq_interpolate.debug.output @@ -0,0 +1,3 @@ + +DEAL::Success, dim = 2 +DEAL::Success, dim = 3 -- 2.39.5