From: heister Date: Tue, 12 Feb 2013 20:12:13 +0000 (+0000) Subject: rename dof_indices() function to get_active_or_mg_dof_indices() and document a couple... X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=ed6ca71a7f1dd1bb6391bed6c51d475555738563;p=dealii-svn.git rename dof_indices() function to get_active_or_mg_dof_indices() and document a couple of functions git-svn-id: https://svn.dealii.org/trunk@28338 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/include/deal.II/dofs/dof_accessor.h b/deal.II/include/deal.II/dofs/dof_accessor.h index beebb35966..bac83d687b 100644 --- a/deal.II/include/deal.II/dofs/dof_accessor.h +++ b/deal.II/include/deal.II/dofs/dof_accessor.h @@ -136,7 +136,7 @@ namespace internal * the dimension of the Triangulation this object refers to as well as * the dimension of the space into which it is embedded. Finally, the * template argument level_dof_access governs the - * behavior of the function dof_indices(). See the section on Generic loops below. + * behavior of the function get_active_or_mg_dof_indices(). See the section on Generic loops below. * *

Typedefs

* @@ -152,7 +152,7 @@ namespace internal * dofs of the active cells of the Triangulation or on the level dodfs * of a single level or the whole grid hierarchy. In order to use * polymorphism in such loops, they access degrees of freedom through - * the function dof_indices(), which changes behavior according to the + * the function get_active_or_mg_dof_indices(), which changes behavior according to the * third template argument. If the argument is false, then the active * dofs of active cells are accessed. If it is true, the level dofs * are used. DoFHandler has functions, for instance begin() and @@ -160,7 +160,7 @@ namespace internal * they can be cast into each other, in case this is needed, since * they access the same data. * - * It is highly recommended to use the function dof_indices() in + * It is highly recommended to use the function get_active_or_mg_dof_indices() in * generic loops in lieu of get_dof_indices() or get_mg_dof_indices(). * *

Inheritance

@@ -178,7 +178,7 @@ namespace internal * * @ingroup dofs * @ingroup Accessors - * @author Wolfgang Bangerth, 1998, 2006, 2008, Timo Heister, Guido Kanschat, 2012 + * @author Wolfgang Bangerth, 1998, 2006, 2008, Timo Heister, Guido Kanschat, 2012, 2013 */ template class DoFAccessor : public dealii::internal::DoFAccessor::Inheritance::BaseClass @@ -306,7 +306,7 @@ public: /** * Tell the caller whether - * dof_indices() accesses active + * get_active_or_mg_dof_indices() accesses active * or level dofs. */ static bool is_level_cell(); @@ -441,7 +441,23 @@ public: void get_dof_indices (std::vector &dof_indices, const unsigned int fe_index = DH::default_fe_index) const; + /** + * Return the indices of the dofs of this + * object in the standard ordering: dofs + * on vertex 0, dofs on vertex 1, ... + * dofs on line 0, dofs on line 1, ..., + * then quads, then hexes. + * + * It is assumed that the vector already + * has the right size beforehand. The + * indices refer to the local numbering + * for the level this line lives on. + */ void get_mg_dof_indices (const int level, std::vector &dof_indices, const unsigned int fe_index = DH::default_fe_index) const; + + /** + * Sets the level DoF indices that are returned by get_mg_dof_indices. + */ void set_mg_dof_indices (const int level, const std::vector &dof_indices, const unsigned int fe_index = DH::default_fe_index); /** @@ -484,6 +500,10 @@ public: const unsigned int i, const unsigned int fe_index = DH::default_fe_index) const; + /** + * Returns the Global DoF index of the i degree associated with the @p + * vertexth vertex on the level @p level. Also see vertex_dof_index(). + */ unsigned int mg_vertex_dof_index (const int level, const unsigned int vertex, const unsigned int i, const unsigned int fe_index = DH::default_fe_index) const; /** @@ -543,6 +563,9 @@ public: unsigned int dof_index (const unsigned int i, const unsigned int fe_index = DH::default_fe_index) const; + /** + * Returns the dof_index on the given level. Also see dof_index. + */ unsigned int mg_dof_index (const int level, const unsigned int i) const; /** @@ -2076,11 +2099,9 @@ public: * * Examples for this use are in the implementation of DoFRenumbering. */ - void dof_indices (std::vector &dof_indices) const; + void get_active_or_mg_dof_indices (std::vector &dof_indices) const; /** - * @deprecated Use dof_indices() instead. - * * Return the indices of the dofs of this * quad in the standard ordering: dofs * on vertex 0, dofs on vertex 1, etc, @@ -2111,11 +2132,13 @@ public: * * This is a function which requires that the cell be active. * + * Also see get_active_or_mg_dof_indices(). + * * @deprecated Currently, this function can also be called for non-active cells, if all degrees of freedom of the FiniteElement are located in vertices. This functionality will vanish in a future release. */ void get_dof_indices (std::vector &dof_indices) const; /** - * @deprecated Use dof_indices() with level_cell_iterator returned from begin_mg(). + * @deprecated Use get_active_or_mg_dof_indices() with level_cell_iterator returned from begin_mg(). * * Retrieve the global indices of the degrees of freedom on this cell in the level vector associated to the level of the cell. */ @@ -2170,6 +2193,11 @@ public: * given DoF handler class. */ void set_dof_indices (const std::vector &dof_indices); + + /** + * Set the Level DoF indices of this + * cell to the given values. + */ void set_mg_dof_indices (const std::vector &dof_indices); /** diff --git a/deal.II/include/deal.II/dofs/dof_accessor.templates.h b/deal.II/include/deal.II/dofs/dof_accessor.templates.h index 10aef3789f..2f654a07a4 100644 --- a/deal.II/include/deal.II/dofs/dof_accessor.templates.h +++ b/deal.II/include/deal.II/dofs/dof_accessor.templates.h @@ -3635,7 +3635,7 @@ void DoFCellAccessor::set_mg_dof_indices (const std::vector inline -void DoFCellAccessor::dof_indices (std::vector &dof_indices) const +void DoFCellAccessor::get_active_or_mg_dof_indices (std::vector &dof_indices) const { if (lda) get_mg_dof_indices (dof_indices); diff --git a/deal.II/source/dofs/dof_renumbering.cc b/deal.II/source/dofs/dof_renumbering.cc index c5ac2bbc76..3f0d782f33 100644 --- a/deal.II/source/dofs/dof_renumbering.cc +++ b/deal.II/source/dofs/dof_renumbering.cc @@ -287,7 +287,7 @@ namespace DoFRenumbering dofs_on_this_cell.resize (dofs_per_cell); - cell->dof_indices (dofs_on_this_cell); + cell->get_active_or_mg_dof_indices (dofs_on_this_cell); for (unsigned int i=0; i dofs_on_this_cell[j]) @@ -668,7 +668,7 @@ namespace DoFRenumbering const unsigned int fe_index = cell->active_fe_index(); const unsigned int dofs_per_cell = fe_collection[fe_index].dofs_per_cell; local_dof_indices.resize (dofs_per_cell); - cell->dof_indices (local_dof_indices); + cell->get_active_or_mg_dof_indices (local_dof_indices); for (unsigned int i=0; iget_dof_handler().locally_owned_dofs().is_element(local_dof_indices[i])) component_to_dof_map[component_list[fe_index][i]]. @@ -974,7 +974,7 @@ namespace DoFRenumbering const unsigned int fe_index = cell->active_fe_index(); const unsigned int dofs_per_cell =fe_collection[fe_index].dofs_per_cell; local_dof_indices.resize (dofs_per_cell); - cell->dof_indices (local_dof_indices); + cell->get_active_or_mg_dof_indices (local_dof_indices); for (unsigned int i=0; iget_dof_handler().locally_owned_dofs().is_element(local_dof_indices[i])) block_to_dof_map[block_list[fe_index][i]]. @@ -1386,7 +1386,7 @@ namespace DoFRenumbering unsigned int n_cell_dofs = (*cell)->get_fe().n_dofs_per_cell(); cell_dofs.resize(n_cell_dofs); - (*cell)->dof_indices(cell_dofs); + (*cell)->get_active_or_mg_dof_indices(cell_dofs); // Sort here to make sure that // degrees of freedom inside a @@ -1459,7 +1459,7 @@ namespace DoFRenumbering { Assert ((*cell)->level() == (int) level, ExcInternalError()); - (*cell)->dof_indices(cell_dofs); + (*cell)->get_active_or_mg_dof_indices(cell_dofs); std::sort(cell_dofs.begin(), cell_dofs.end()); for (unsigned int i=0; i &fe_values = hp_fe_values.get_present_fe_values (); - begin->dof_indices(local_dof_indices); + begin->get_active_or_mg_dof_indices(local_dof_indices); const std::vector > &points = fe_values.get_quadrature_points (); for (unsigned int i=0; i::cell_iterator &begin_tria = begin; - begin->dof_indices(local_dof_indices); + begin->get_active_or_mg_dof_indices(local_dof_indices); fe_values.reinit (begin_tria); const std::vector > &points = fe_values.get_quadrature_points ();