]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Provide nth_active_fe_index.
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 27 Jul 2006 16:21:52 +0000 (16:21 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 27 Jul 2006 16:21:52 +0000 (16:21 +0000)
git-svn-id: https://svn.dealii.org/trunk@13448 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/dofs/hp_dof_handler.h
deal.II/deal.II/include/dofs/hp_dof_objects.h
deal.II/deal.II/source/dofs/hp_dof_handler.cc
deal.II/deal.II/source/dofs/hp_dof_objects.cc

index b3a4c00f93af804a7d33417c7e996782f4499641..cf76ae6d81f2e8a7aeb9c48cd6a5588f52145153 100644 (file)
@@ -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 <int structdim>
+      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
index a3d6198a76934711890af7f8a8469dfdf1be4658..d386fb09f50ce50008fffeb48d093bb70d6d997e 100644 (file)
@@ -201,6 +201,18 @@ namespace internal
         n_active_fe_indices (const ::hp::DoFHandler<spacedim> &dof_handler,
                              const unsigned int               obj_index) const;
 
+                                        /**
+                                         * Return the fe_index of the
+                                         * n-th active finite element
+                                         * on this object.
+                                         */
+        template <int spacedim>
+        unsigned int
+        nth_active_fe_index (const ::hp::DoFHandler<spacedim> &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
index a322dbfe49462aec05068d5ef3a4a6f3614dab5b..ffc74f793ca812ac34f29fbddb8aa6b836511bdc 100644 (file)
@@ -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
 
 
index 84f8f0093a027c089aeea4437bda984856fbbe7b..c3bf9c0ac1172d828d5ebe048c6a73167c419e6b 100644 (file)
@@ -18,6 +18,8 @@
 #include <dofs/hp_dof_handler.h>
 
 
+//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

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.