From: Ralf Hartmann Date: Fri, 12 Jun 2009 06:07:59 +0000 (+0000) Subject: New hp::DoFHandler::set/get_active_fe_indices functions. X-Git-Tag: v8.0.0~7618 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=8c34d233584e5215558c2346ed06777e0df4e56c;p=dealii.git New hp::DoFHandler::set/get_active_fe_indices functions. git-svn-id: https://svn.dealii.org/trunk@18929 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/hp/dof_handler.h b/deal.II/deal.II/include/hp/dof_handler.h index 7636d27379..4a34f2a4b3 100644 --- a/deal.II/deal.II/include/hp/dof_handler.h +++ b/deal.II/deal.II/include/hp/dof_handler.h @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 2005, 2006, 2007, 2008 by the deal.II authors +// Copyright (C) 2005, 2006, 2007, 2008, 2009 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -172,8 +172,8 @@ namespace hp * Go through the triangulation and * distribute the degrees of freedoms * needed for the given finite element - * according to the given distribution - * method. + * according to the current distribution + * of active fe indices. * * A pointer of the transferred * finite element is @@ -189,6 +189,23 @@ namespace hp */ virtual void distribute_dofs (const hp::FECollection &fe); + /** + * Go through the triangulation and set + * the active FE indices of all active + * cells to the values given in @p + * active_fe_indices. + */ + void set_active_fe_indices (const std::vector& active_fe_indices); + + /** + * Go through the triangulation and + * store the active FE indices of all + * active cells to the vector @p + * active_fe_indices. This vector is + * resized, if necessary. + */ + void get_active_fe_indices (std::vector& active_fe_indices) const; + /** * Clear all data of this object and * especially delete the lock this object @@ -1084,7 +1101,7 @@ namespace hp * Store the number of dofs * created last time. */ - unsigned int used_dofs; + unsigned int used_dofs; /** * Array to store the indices diff --git a/deal.II/deal.II/source/hp/dof_handler.cc b/deal.II/deal.II/source/hp/dof_handler.cc index 3af05721cb..5d312bcfc7 100644 --- a/deal.II/deal.II/source/hp/dof_handler.cc +++ b/deal.II/deal.II/source/hp/dof_handler.cc @@ -3397,6 +3397,44 @@ namespace hp } } + + + template + void DoFHandler::set_active_fe_indices (const std::vector& active_fe_indices) + { + Assert(active_fe_indices.size()==get_tria().n_active_cells(), + ExcDimensionMismatch(active_fe_indices.size(), get_tria().n_active_cells())); + + create_active_fe_table (); + // we could set the values directly, since + // they are stored as protected data of + // this object, but for simplicity we use + // the cell-wise access. this way we also + // have to pass some debug-mode tests which + // we would have to duplicate ourselves + // otherwise + active_cell_iterator cell=begin_active(), + endc=end(); + for (unsigned int i=0; cell!=endc; ++cell, ++i) + cell->set_active_fe_index(active_fe_indices[i]); + } + + + + template + void DoFHandler::get_active_fe_indices (std::vector& active_fe_indices) const + { + active_fe_indices.resize(get_tria().n_active_cells()); + + // we could try to extract the values directly, since + // they are stored as protected data of + // this object, but for simplicity we use + // the cell-wise access. + active_cell_iterator cell=begin_active(), + endc=end(); + for (unsigned int i=0; cell!=endc; ++cell, ++i) + active_fe_indices[i]=cell->active_fe_index(); + } diff --git a/deal.II/doc/news/changes.h b/deal.II/doc/news/changes.h index 9c867c8093..e9ad4a1c06 100644 --- a/deal.II/doc/news/changes.h +++ b/deal.II/doc/news/changes.h @@ -135,6 +135,20 @@ inconvenience this causes.

deal.II

    +
  1. +

    + New: The new hp::DoFHandler::set_active_fe_indices function allows + to distribute all active FE indices at once based on a given + vector. This might be useful if this information is stored + somewhere and has to be reconstructed or else if two DoFHandler + objects with the same FE index distribution should be created. + There is now also a corresponding + hp::DoFHandler::get_active_fe_indices function. +
    + (Tobias Leicht, RH 2009/06/12) +

    +
  2. +
  3. Fix: The projection of quadrature points to subfaces in