From 5d0c434efa1700e2a9624523a56d532109615d7a Mon Sep 17 00:00:00 2001 From: David Wells <wellsd2@rpi.edu> Date: Sun, 19 Feb 2017 21:15:32 -0500 Subject: [PATCH] Add missing DoFAccessor<0, ..> member functions. --- doc/news/changes/minor/20170219DavidWells | 5 + include/deal.II/dofs/dof_accessor.templates.h | 117 +++++++++++++++ tests/hp/accessor_0.cc | 134 ++++++++++++++++++ tests/hp/accessor_0.output | 40 ++++++ 4 files changed, 296 insertions(+) create mode 100644 doc/news/changes/minor/20170219DavidWells create mode 100644 tests/hp/accessor_0.cc create mode 100644 tests/hp/accessor_0.output diff --git a/doc/news/changes/minor/20170219DavidWells b/doc/news/changes/minor/20170219DavidWells new file mode 100644 index 0000000000..3ae812b83a --- /dev/null +++ b/doc/news/changes/minor/20170219DavidWells @@ -0,0 +1,5 @@ +Fixed: The DoFAccessor specialization for <code>structdim == 0</code> and +<code>dimension == 1</code> was missing several definitions, which were +required for use of hp::DoFHandler in 1D. These have now been added. +<br> +(David Wells, 2017/02/19) diff --git a/include/deal.II/dofs/dof_accessor.templates.h b/include/deal.II/dofs/dof_accessor.templates.h index c89b068d65..25d39b02aa 100644 --- a/include/deal.II/dofs/dof_accessor.templates.h +++ b/include/deal.II/dofs/dof_accessor.templates.h @@ -2285,6 +2285,33 @@ DoFAccessor<0,DoFHandlerType<1,spacedim>, level_dof_access>::set_dof_handler +template <template <int, int> class DoFHandlerType, int spacedim, bool level_dof_access> +inline +void +DoFAccessor<0,DoFHandlerType<1,spacedim>, level_dof_access>::set_dof_index +(const unsigned int /*i*/, + const types::global_dof_index /*index*/, + const unsigned int /*fe_index*/) const +{ + Assert (false, ExcNotImplemented()); +} + + + +template <template <int, int> class DoFHandlerType, int spacedim, bool level_dof_access> +inline +void +DoFAccessor<0,DoFHandlerType<1,spacedim>, level_dof_access>::set_vertex_dof_index +(const unsigned int /*vertex*/, + const unsigned int /*i*/, + const types::global_dof_index /*index*/, + const unsigned int /*fe_index*/) const +{ + Assert (false, ExcNotImplemented()); +} + + + template <template <int, int> class DoFHandlerType, int spacedim, bool level_dof_access> inline const DoFHandlerType<1,spacedim> & @@ -2352,6 +2379,70 @@ vertex_dof_index (const unsigned int vertex, +template <template <int, int> class DoFHandlerType, int spacedim, bool level_dof_access> +inline +types::global_dof_index +DoFAccessor<0,DoFHandlerType<1,spacedim>, level_dof_access>:: +dof_index (const unsigned int i, + const unsigned int fe_index) const +{ + return dealii::internal::DoFAccessor::Implementation::get_vertex_dof_index + (*this->dof_handler, + this->vertex_index(0), + fe_index, + i); +} + + + +template <template <int, int> class DoFHandlerType, int spacedim, bool level_dof_access> +inline +unsigned int +DoFAccessor<0,DoFHandlerType<1,spacedim>, level_dof_access>:: +n_active_fe_indices() const +{ + Assert(false, ExcNotImplemented()); + return numbers::invalid_unsigned_int; +} + + + +template <template <int, int> class DoFHandlerType, int spacedim, bool level_dof_access> +inline +unsigned int +DoFAccessor<0,DoFHandlerType<1,spacedim>, level_dof_access>:: +nth_active_fe_index(const unsigned int /*n*/) const +{ + Assert(false, ExcNotImplemented()); + return numbers::invalid_unsigned_int; +} + + + +template <template <int, int> class DoFHandlerType, int spacedim, bool level_dof_access> +inline +bool +DoFAccessor<0,DoFHandlerType<1,spacedim>, level_dof_access>:: +fe_index_is_active (const unsigned int fe_index) const +{ + Assert(false, ExcNotImplemented()); + return false; +} + + + +template <template <int, int> class DoFHandlerType, int spacedim, bool level_dof_access> +inline +const FiniteElement<DoFHandlerType<1,spacedim>::dimension,DoFHandlerType<1,spacedim>::space_dimension> & +DoFAccessor<0,DoFHandlerType<1,spacedim>, level_dof_access>:: +get_fe (const unsigned int fe_index) const +{ + Assert (this->dof_handler != 0, ExcInvalidObject()); + return dof_handler->get_fe()[fe_index]; +} + + + template <template <int, int> class DoFHandlerType, int spacedim, bool level_dof_access> inline void @@ -2387,6 +2478,32 @@ DoFAccessor<0,DoFHandlerType<1,spacedim>, level_dof_access>::child (const unsign +template <template <int, int> class DoFHandlerType, int spacedim, bool level_dof_access> +inline +typename dealii::internal::DoFHandler::Iterators<DoFHandlerType<1,spacedim>, level_dof_access>::line_iterator +DoFAccessor<0,DoFHandlerType<1,spacedim>, level_dof_access>::line +(const unsigned int /*c*/) const +{ + Assert(false, ExcNotImplemented()); + return typename dealii::internal:: + DoFHandler::Iterators<DoFHandlerType<1,spacedim>, level_dof_access>::line_iterator(); +} + + + +template <template <int, int> class DoFHandlerType, int spacedim, bool level_dof_access> +inline +typename dealii::internal::DoFHandler::Iterators<DoFHandlerType<1,spacedim>, level_dof_access>::quad_iterator +DoFAccessor<0,DoFHandlerType<1,spacedim>, level_dof_access>::quad +(const unsigned int /*c*/) const +{ + Assert(false, ExcNotImplemented()); + return typename dealii::internal:: + DoFHandler::Iterators<DoFHandlerType<1,spacedim>, level_dof_access>::quad_iterator(); +} + + + template <template <int, int> class DoFHandlerType, int spacedim, bool level_dof_access> template <int dim2, class DoFHandlerType2, bool level_dof_access2> inline diff --git a/tests/hp/accessor_0.cc b/tests/hp/accessor_0.cc new file mode 100644 index 0000000000..623cae7d9b --- /dev/null +++ b/tests/hp/accessor_0.cc @@ -0,0 +1,134 @@ +// --------------------------------------------------------------------- +// +// 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. +// +// --------------------------------------------------------------------- + +// Ensure that the newly added instantiations for DoFAccessor<0, ...> +// (dof_index and get_fe) work correctly for an hp::DoFHandler. + +#include <deal.II/fe/fe_q.h> + +#include <deal.II/grid/grid_generator.h> +#include <deal.II/grid/tria.h> + +#include <deal.II/hp/dof_handler.h> +#include <deal.II/hp/fe_collection.h> + +#include "../tests.h" + +int main() +{ + std::ofstream logfile("output"); + deallog.attach(logfile); + deallog << std::boolalpha; + + using namespace dealii; + + Triangulation<1> triangulation; + GridGenerator::hyper_cube(triangulation); + triangulation.refine_global(2); + + hp::FECollection<1> fe_collection; + hp::DoFHandler<1> dof_handler(triangulation); + fe_collection.push_back(FE_Q<1>(2)); + fe_collection.push_back(FE_Q<1>(4)); + fe_collection.push_back(FE_Q<1>(6)); + + const unsigned int n_fe_indices = 3; + { + typename hp::DoFHandler<1>::active_cell_iterator + cell = dof_handler.begin_active(); + dof_handler.begin_active()->set_active_fe_index(1); + ++cell; // go to cell 1 + ++cell; // go to cell 2 + ++cell; // go to cell 3 + cell->set_active_fe_index(2); + } + dof_handler.distribute_dofs(fe_collection); + + // dof_index + // get_fe + + std::vector<types::global_dof_index> dof_indices; + + typename hp::DoFHandler<1>::active_cell_iterator + cell = dof_handler.begin_active(), + endc = dof_handler.end(); + for (; cell != endc; ++cell) + { + deallog << "===================================" + << std::endl; + deallog << "cell center: " + << cell->center() + << std::endl; + + dof_indices.resize(fe_collection[cell->active_fe_index()].dofs_per_cell); + cell->get_dof_indices(dof_indices); + deallog << "cell dofs: "; + for (unsigned int dof_n = 0; dof_n < dof_indices.size(); ++dof_n) + { + deallog << dof_indices[dof_n]; + if (dof_n != dof_indices.size() - 1) + { + deallog << ", "; + } + } + deallog << std::endl; + + // see if we have a neighbor on the right. If so, the common vertex + // should be associated with two FE indices. + const typename hp::DoFHandler<1>::active_cell_iterator neighbor = cell->neighbor(1); + if (neighbor != dof_handler.end()) + { + const unsigned int current_index = cell->active_fe_index(); + const unsigned int neighbor_index = neighbor->active_fe_index(); + deallog << "dof index (current cell, current index): " + << cell->face(1)->dof_index(0, current_index) + << " (neighbor cell, current index): " + << neighbor->face(0)->dof_index(0, current_index) + << std::endl + << "dof index (current cell, neighbor index): " + << cell->face(1)->dof_index(0, neighbor_index) + << " (neighbor cell, neighbor index): " + << neighbor->face(0)->dof_index(0, neighbor_index) + << std::endl; + } + + for (unsigned int fe_index = 0; fe_index < n_fe_indices; ++fe_index) + { + const bool index_is_active = cell->fe_index_is_active(fe_index); + deallog << "cell uses fe index " << fe_index << ": " + << index_is_active + << std::endl; + + for (unsigned int face_n = 0; face_n < GeometryInfo<1>::faces_per_cell; + ++face_n) + { + AssertThrow (&cell->face(face_n)->get_fe(fe_index) == &fe_collection[fe_index], + ExcMessage("The result of get_fe should always return" + " a known finite element.")); + + if (index_is_active) + { + deallog << "vertex dof index: " + << cell->face(face_n)->dof_index(0, fe_index) + << std::endl; + } + } + } + } + + deallog << "OK" << std::endl; + + return 0; +} diff --git a/tests/hp/accessor_0.output b/tests/hp/accessor_0.output new file mode 100644 index 0000000000..fdc3a9b9ed --- /dev/null +++ b/tests/hp/accessor_0.output @@ -0,0 +1,40 @@ + +DEAL::=================================== +DEAL::cell center: 0.125000 +DEAL::cell dofs: 0, 4, 1, 2, 3 +DEAL::dof index (current cell, current index): 4 (neighbor cell, current index): 4 +DEAL::dof index (current cell, neighbor index): 4 (neighbor cell, neighbor index): 4 +DEAL::cell uses fe index 0: false +DEAL::cell uses fe index 1: true +DEAL::vertex dof index: 0 +DEAL::vertex dof index: 4 +DEAL::cell uses fe index 2: false +DEAL::=================================== +DEAL::cell center: 0.375000 +DEAL::cell dofs: 4, 5, 6 +DEAL::dof index (current cell, current index): 5 (neighbor cell, current index): 5 +DEAL::dof index (current cell, neighbor index): 5 (neighbor cell, neighbor index): 5 +DEAL::cell uses fe index 0: true +DEAL::vertex dof index: 4 +DEAL::vertex dof index: 5 +DEAL::cell uses fe index 1: false +DEAL::cell uses fe index 2: false +DEAL::=================================== +DEAL::cell center: 0.625000 +DEAL::cell dofs: 5, 7, 8 +DEAL::dof index (current cell, current index): 7 (neighbor cell, current index): 7 +DEAL::dof index (current cell, neighbor index): 7 (neighbor cell, neighbor index): 7 +DEAL::cell uses fe index 0: true +DEAL::vertex dof index: 5 +DEAL::vertex dof index: 7 +DEAL::cell uses fe index 1: false +DEAL::cell uses fe index 2: false +DEAL::=================================== +DEAL::cell center: 0.875000 +DEAL::cell dofs: 7, 9, 10, 11, 12, 13, 14 +DEAL::cell uses fe index 0: false +DEAL::cell uses fe index 1: false +DEAL::cell uses fe index 2: true +DEAL::vertex dof index: 7 +DEAL::vertex dof index: 9 +DEAL::OK -- 2.39.5