From: bangerth Date: Thu, 27 Jul 2006 16:21:52 +0000 (+0000) Subject: Provide nth_active_fe_index. X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=8607e6e957096ed14b2401d9118cfc2f238bb610;p=dealii-svn.git Provide nth_active_fe_index. git-svn-id: https://svn.dealii.org/trunk@13448 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/dofs/hp_dof_handler.h b/deal.II/deal.II/include/dofs/hp_dof_handler.h index b3a4c00f93..cf76ae6d81 100644 --- a/deal.II/deal.II/include/dofs/hp_dof_handler.h +++ b/deal.II/deal.II/include/dofs/hp_dof_handler.h @@ -1085,6 +1085,16 @@ namespace hp unsigned int n_active_fe_indices (const unsigned int obj_level, const unsigned int obj_index) const; + /** + * Return the fe index of the + * n-th active finite element on + * this object. + */ + template + unsigned int nth_active_fe_index (const unsigned int obj_level, + const unsigned int obj_index, + const unsigned int n) const; + /** * Return whether fe-index is * an active fe, calls the diff --git a/deal.II/deal.II/include/dofs/hp_dof_objects.h b/deal.II/deal.II/include/dofs/hp_dof_objects.h index a3d6198a76..d386fb09f5 100644 --- a/deal.II/deal.II/include/dofs/hp_dof_objects.h +++ b/deal.II/deal.II/include/dofs/hp_dof_objects.h @@ -201,6 +201,18 @@ namespace internal n_active_fe_indices (const ::hp::DoFHandler &dof_handler, const unsigned int obj_index) const; + /** + * Return the fe_index of the + * n-th active finite element + * on this object. + */ + template + unsigned int + nth_active_fe_index (const ::hp::DoFHandler &dof_handler, + const unsigned int obj_level, + const unsigned int obj_index, + const unsigned int n) const; + /** * Check whether a given * finite element index is diff --git a/deal.II/deal.II/source/dofs/hp_dof_handler.cc b/deal.II/deal.II/source/dofs/hp_dof_handler.cc index a322dbfe49..ffc74f793c 100644 --- a/deal.II/deal.II/source/dofs/hp_dof_handler.cc +++ b/deal.II/deal.II/source/dofs/hp_dof_handler.cc @@ -1426,6 +1426,21 @@ namespace hp return this->levels[obj_level]->lines.n_active_fe_indices (*this, obj_index); } + + + + template <> + template <> + unsigned int + DoFHandler<1>::nth_active_fe_index<1> (const unsigned int obj_level, + const unsigned int obj_index, + const unsigned int n) const + { + return this->levels[obj_level]->lines.nth_active_fe_index (*this, + obj_level, + obj_index, + n); + } #endif @@ -1492,6 +1507,21 @@ namespace hp } + template <> + template <> + unsigned int + DoFHandler<2>::nth_active_fe_index<1> (const unsigned int obj_level, + const unsigned int obj_index, + const unsigned int n) const + { + return this->faces->lines.nth_active_fe_index (*this, + obj_level, + obj_index, + n); + } + + + template <> template <> unsigned int @@ -1552,6 +1582,21 @@ namespace hp return this->levels[obj_level]->quads.n_active_fe_indices (*this, obj_index); } + + + + template <> + template <> + unsigned int + DoFHandler<2>::nth_active_fe_index<2> (const unsigned int obj_level, + const unsigned int obj_index, + const unsigned int n) const + { + return this->levels[obj_level]->quads.nth_active_fe_index (*this, + obj_level, + obj_index, + n); + } #endif #if deal_II_dimension == 3 @@ -1618,6 +1663,21 @@ namespace hp + template <> + template <> + unsigned int + DoFHandler<3>::nth_active_fe_index<1> (const unsigned int obj_level, + const unsigned int obj_index, + const unsigned int n) const + { + return this->faces->lines.nth_active_fe_index (*this, + obj_level, + obj_index, + n); + } + + + template <> template <> unsigned int @@ -1729,6 +1789,23 @@ namespace hp obj_index); } + + + template <> + template <> + unsigned int + DoFHandler<3>::nth_active_fe_index<2> (const unsigned int obj_level, + const unsigned int obj_index, + const unsigned int n) const + { + return this->faces->quads.nth_active_fe_index (*this, + obj_level, + obj_index, + n); + } + + + template <> template <> unsigned int @@ -1738,6 +1815,21 @@ namespace hp return this->levels[obj_level]->hexes.n_active_fe_indices (*this, obj_index); } + + + + template <> + template <> + unsigned int + DoFHandler<3>::nth_active_fe_index<3> (const unsigned int obj_level, + const unsigned int obj_index, + const unsigned int n) const + { + return this->levels[obj_level]->hexes.nth_active_fe_index (*this, + obj_level, + obj_index, + n); + } #endif diff --git a/deal.II/deal.II/source/dofs/hp_dof_objects.cc b/deal.II/deal.II/source/dofs/hp_dof_objects.cc index 84f8f0093a..c3bf9c0ac1 100644 --- a/deal.II/deal.II/source/dofs/hp_dof_objects.cc +++ b/deal.II/deal.II/source/dofs/hp_dof_objects.cc @@ -18,6 +18,8 @@ #include +//TODO: Several of the functions in this file are identical. some template trickery should make it possible to merge them + namespace internal { namespace hp @@ -160,6 +162,40 @@ namespace internal + template <> + template <> + unsigned int + DoFObjects<1>:: + nth_active_fe_index<1> (const ::hp::DoFHandler<1> &dof_handler, + const unsigned int obj_level, + const unsigned int obj_index, + const unsigned int n) const + { + Assert (&dof_handler != 0, + ExcMessage ("No DoFHandler is specified for this iterator")); + Assert (&dof_handler.get_fe() != 0, + ExcMessage ("No finite element collection is associated with " + "this DoFHandler")); + Assert (obj_index < dof_offsets.size(), + ExcIndexRange (obj_index, 0, dof_offsets.size())); + + // make sure we are on an + // object for which DoFs have + // been allocated at all + Assert (dof_offsets[obj_index] != deal_II_numbers::invalid_unsigned_int, + ExcMessage ("You are trying to access degree of freedom " + "information for an object on which no such " + "information is available")); + + // this is a cell, so there is + // only a single fe_index + Assert (n == 0, ExcIndexRange (n, 0, 1)); + + return dof_handler.levels[obj_level]->active_fe_indices[obj_index]; + } + + + template <> template <> bool @@ -385,6 +421,66 @@ namespace internal + template <> + template <> + unsigned int + DoFObjects<1>:: + nth_active_fe_index<2> (const ::hp::DoFHandler<2> &dof_handler, + const unsigned int , + const unsigned int obj_index, + const unsigned int n) const + { + Assert (&dof_handler != 0, + ExcMessage ("No DoFHandler is specified for this iterator")); + Assert (&dof_handler.get_fe() != 0, + ExcMessage ("No finite element collection is associated with " + "this DoFHandler")); + Assert (obj_index < dof_offsets.size(), + ExcIndexRange (obj_index, 0, dof_offsets.size())); + + // make sure we are on an + // object for which DoFs have + // been allocated at all + Assert (dof_offsets[obj_index] != deal_II_numbers::invalid_unsigned_int, + ExcMessage ("You are trying to access degree of freedom " + "information for an object on which no such " + "information is available")); + + Assert (n < n_active_fe_indices(dof_handler, obj_index), + ExcIndexRange (n, 0, + n_active_fe_indices(dof_handler, obj_index))); + + // we are in higher space + // dimensions, so there may + // be multiple finite + // elements associated with + // this object. hop along + // the list of index sets + // until we find the one + // with the correct + // fe_index, and then poke + // into that part. trigger + // an exception if we can't + // find a set for this + // particular fe_index + const unsigned int starting_offset = dof_offsets[obj_index]; + const unsigned int *pointer = &dofs[starting_offset]; + unsigned int counter = 0; + while (true) + { + Assert (*pointer != deal_II_numbers::invalid_unsigned_int, + ExcInternalError()); + + if (counter == n) + return *pointer; + + ++counter; + pointer += dof_handler.get_fe()[*pointer].dofs_per_line + 1; + } + } + + + template <> template <> bool @@ -627,6 +723,66 @@ namespace internal + template <> + template <> + unsigned int + DoFObjects<1>:: + nth_active_fe_index<3> (const ::hp::DoFHandler<3> &dof_handler, + const unsigned int , + const unsigned int obj_index, + const unsigned int n) const + { + Assert (&dof_handler != 0, + ExcMessage ("No DoFHandler is specified for this iterator")); + Assert (&dof_handler.get_fe() != 0, + ExcMessage ("No finite element collection is associated with " + "this DoFHandler")); + Assert (obj_index < dof_offsets.size(), + ExcIndexRange (obj_index, 0, dof_offsets.size())); + + // make sure we are on an + // object for which DoFs have + // been allocated at all + Assert (dof_offsets[obj_index] != deal_II_numbers::invalid_unsigned_int, + ExcMessage ("You are trying to access degree of freedom " + "information for an object on which no such " + "information is available")); + + Assert (n < n_active_fe_indices(dof_handler, obj_index), + ExcIndexRange (n, 0, + n_active_fe_indices(dof_handler, obj_index))); + + // we are in higher space + // dimensions, so there may + // be multiple finite + // elements associated with + // this object. hop along + // the list of index sets + // until we find the one + // with the correct + // fe_index, and then poke + // into that part. trigger + // an exception if we can't + // find a set for this + // particular fe_index + const unsigned int starting_offset = dof_offsets[obj_index]; + const unsigned int *pointer = &dofs[starting_offset]; + unsigned int counter = 0; + while (true) + { + Assert (*pointer != deal_II_numbers::invalid_unsigned_int, + ExcInternalError()); + + if (counter == n) + return *pointer; + + ++counter; + pointer += dof_handler.get_fe()[*pointer].dofs_per_line + 1; + } + } + + + template <> template <> bool @@ -815,6 +971,40 @@ namespace internal + template <> + template <> + unsigned int + DoFObjects<2>:: + nth_active_fe_index<2> (const ::hp::DoFHandler<2> &dof_handler, + const unsigned int obj_level, + const unsigned int obj_index, + const unsigned int n) const + { + Assert (&dof_handler != 0, + ExcMessage ("No DoFHandler is specified for this iterator")); + Assert (&dof_handler.get_fe() != 0, + ExcMessage ("No finite element collection is associated with " + "this DoFHandler")); + Assert (obj_index < dof_offsets.size(), + ExcIndexRange (obj_index, 0, dof_offsets.size())); + + // make sure we are on an + // object for which DoFs have + // been allocated at all + Assert (dof_offsets[obj_index] != deal_II_numbers::invalid_unsigned_int, + ExcMessage ("You are trying to access degree of freedom " + "information for an object on which no such " + "information is available")); + + // this is a cell, so there is + // only a single fe_index + Assert (n == 0, ExcIndexRange (n, 0, 1)); + + return dof_handler.levels[obj_level]->active_fe_indices[obj_index]; + } + + + template <> template <> bool @@ -1038,6 +1228,66 @@ namespace internal + template <> + template <> + unsigned int + DoFObjects<2>:: + nth_active_fe_index<3> (const ::hp::DoFHandler<3> &dof_handler, + const unsigned int , + const unsigned int obj_index, + const unsigned int n) const + { + Assert (&dof_handler != 0, + ExcMessage ("No DoFHandler is specified for this iterator")); + Assert (&dof_handler.get_fe() != 0, + ExcMessage ("No finite element collection is associated with " + "this DoFHandler")); + Assert (obj_index < dof_offsets.size(), + ExcIndexRange (obj_index, 0, dof_offsets.size())); + + // make sure we are on an + // object for which DoFs have + // been allocated at all + Assert (dof_offsets[obj_index] != deal_II_numbers::invalid_unsigned_int, + ExcMessage ("You are trying to access degree of freedom " + "information for an object on which no such " + "information is available")); + + Assert (n < n_active_fe_indices(dof_handler, obj_index), + ExcIndexRange (n, 0, + n_active_fe_indices(dof_handler, obj_index))); + + // we are in higher space + // dimensions, so there may + // be multiple finite + // elements associated with + // this object. hop along + // the list of index sets + // until we find the one + // with the correct + // fe_index, and then poke + // into that part. trigger + // an exception if we can't + // find a set for this + // particular fe_index + const unsigned int starting_offset = dof_offsets[obj_index]; + const unsigned int *pointer = &dofs[starting_offset]; + unsigned int counter = 0; + while (true) + { + Assert (*pointer != deal_II_numbers::invalid_unsigned_int, + ExcInternalError()); + + if (counter == n) + return *pointer; + + ++counter; + pointer += dof_handler.get_fe()[*pointer].dofs_per_line + 1; + } + } + + + template <> template <> bool @@ -1225,6 +1475,40 @@ namespace internal + template <> + template <> + unsigned int + DoFObjects<3>:: + nth_active_fe_index<3> (const ::hp::DoFHandler<3> &dof_handler, + const unsigned int obj_level, + const unsigned int obj_index, + const unsigned int n) const + { + Assert (&dof_handler != 0, + ExcMessage ("No DoFHandler is specified for this iterator")); + Assert (&dof_handler.get_fe() != 0, + ExcMessage ("No finite element collection is associated with " + "this DoFHandler")); + Assert (obj_index < dof_offsets.size(), + ExcIndexRange (obj_index, 0, dof_offsets.size())); + + // make sure we are on an + // object for which DoFs have + // been allocated at all + Assert (dof_offsets[obj_index] != deal_II_numbers::invalid_unsigned_int, + ExcMessage ("You are trying to access degree of freedom " + "information for an object on which no such " + "information is available")); + + // this is a cell, so there is + // only a single fe_index + Assert (n == 0, ExcIndexRange (n, 0, 1)); + + return dof_handler.levels[obj_level]->active_fe_indices[obj_index]; + } + + + template <> template <> bool