]> https://gitweb.dealii.org/ - dealii.git/commitdiff
For the moment, a last wave of restructuring and consolidation of code.
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Mon, 1 May 2006 22:30:11 +0000 (22:30 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Mon, 1 May 2006 22:30:11 +0000 (22:30 +0000)
git-svn-id: https://svn.dealii.org/trunk@12946 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/dofs/dof_accessor.templates.h
deal.II/deal.II/include/dofs/dof_levels.h
deal.II/deal.II/include/dofs/hp_dof_levels.h
deal.II/deal.II/source/dofs/dof_accessor.cc
deal.II/deal.II/source/dofs/dof_levels.cc

index a3de91ff911da1e7ca6bb004f5fbfae7afe847a5..8946d4ca8ace2539673de1f2553b480bb0ef16e8 100644 (file)
@@ -190,6 +190,38 @@ DoFObjectAccessor (const Triangulation<DH::dimension> *tria,
 
 
 
+template <class DH>
+inline
+unsigned int
+DoFObjectAccessor<1,DH>::dof_index (const unsigned int i,
+                                   const unsigned int fe_index) const
+{
+  return this->dof_handler->levels[this->present_level]
+    ->get_line_dof_index (*this->dof_handler,
+                          this->present_index,
+                          fe_index,
+                          i);
+}
+
+
+
+template <class DH>
+inline
+void
+DoFObjectAccessor<1,DH>::set_dof_index (const unsigned int i,
+                                       const unsigned int index,
+                                       const unsigned int fe_index) const
+{
+  this->dof_handler->levels[this->present_level]
+    ->set_line_dof_index (*this->dof_handler,
+                          this->present_index,
+                          fe_index,
+                          i,
+                          index);
+}
+
+
+
 template <class DH>
 inline
 void
@@ -266,6 +298,38 @@ DoFObjectAccessor<2,DH>::line (const unsigned int i) const
 
 
 
+template <class DH>
+inline
+unsigned int
+DoFObjectAccessor<2,DH>::dof_index (const unsigned int i,
+                                   const unsigned int fe_index) const
+{
+  return this->dof_handler->levels[this->present_level]
+    ->get_quad_dof_index (*this->dof_handler,
+                          this->present_index,
+                          fe_index,
+                          i);
+}
+
+
+
+template <class DH>
+inline
+void
+DoFObjectAccessor<2,DH>::set_dof_index (const unsigned int i,
+                                       const unsigned int index,
+                                       const unsigned int fe_index) const
+{
+  this->dof_handler->levels[this->present_level]
+    ->set_quad_dof_index (*this->dof_handler,
+                          this->present_index,
+                          fe_index,
+                          i,
+                          index);
+}
+
+
+
 template <class DH>
 inline
 void
@@ -362,6 +426,38 @@ DoFObjectAccessor<3,DH>::quad (const unsigned int i) const
 
 
 
+template <class DH>
+inline
+unsigned int
+DoFObjectAccessor<3,DH>::dof_index (const unsigned int i,
+                                   const unsigned int fe_index) const
+{
+  return this->dof_handler->levels[this->present_level]
+    ->get_hex_dof_index (*this->dof_handler,
+                        this->present_index,
+                        fe_index,
+                        i);
+}
+
+
+
+template <class DH>
+inline
+void
+DoFObjectAccessor<3,DH>::set_dof_index (const unsigned int i,
+                                       const unsigned int index,
+                                       const unsigned int fe_index) const
+{
+  this->dof_handler->levels[this->present_level]
+    ->set_hex_dof_index (*this->dof_handler,
+                        this->present_index,
+                        fe_index,
+                        i,
+                        index);
+}
+
+
+
 template <class DH>
 inline
 void
@@ -416,363 +512,6 @@ DoFObjectAccessor<3,DH>::get_dof_indices (std::vector<unsigned int> &dof_indices
 
 
 
-/*--------------- Functions: DoFObjectAccessor<1,DoFHandler> -----------*/
-
-
-template <>
-inline
-unsigned int
-DoFObjectAccessor<1,DoFHandler<1> >::dof_index (const unsigned int i,
-                                               const unsigned int fe_index) const
-{
-  Assert (fe_index == DoFHandler<1>::default_fe_index,
-         ExcMessage ("Only the default FE index is allowed for non-hp DoFHandler objects"));
-  Assert (static_cast<unsigned int>(this->present_level) < this->dof_handler->levels.size(),
-          ExcMessage ("DoFHandler not initialized"));
-
-  Assert (this->dof_handler != 0, BaseClass::ExcInvalidObject());
-                                  // make sure a FE has been selected
-                                  // and enough room was reserved
-  Assert (&this->dof_handler->get_fe() != 0, BaseClass::ExcInvalidObject());
-  Assert (i<this->dof_handler->get_fe().dofs_per_line,
-         ExcIndexRange (i, 0, this->dof_handler->get_fe().dofs_per_line));
-
-  return this->dof_handler->levels[this->present_level]
-    ->line_dofs[this->present_index*this->dof_handler->get_fe().dofs_per_line+i];
-}
-
-
-
-template <>
-inline
-unsigned int
-DoFObjectAccessor<1,DoFHandler<2> >::dof_index (const unsigned int i,
-                                               const unsigned int fe_index) const
-{
-  Assert (fe_index == DoFHandler<2>::default_fe_index,
-         ExcMessage ("Only the default FE index is allowed for non-hp DoFHandler objects"));
-  Assert (static_cast<unsigned int>(this->present_level) < this->dof_handler->levels.size(),
-          ExcMessage ("DoFHandler not initialized"));
-
-  Assert (this->dof_handler != 0, BaseClass::ExcInvalidObject());
-                                  // make sure a FE has been selected
-                                  // and enough room was reserved
-  Assert (&this->dof_handler->get_fe() != 0, BaseClass::ExcInvalidObject());
-  Assert (i<this->dof_handler->get_fe().dofs_per_line,
-         ExcIndexRange (i, 0, this->dof_handler->get_fe().dofs_per_line));
-
-  return this->dof_handler->levels[this->present_level]
-    ->line_dofs[this->present_index*this->dof_handler->get_fe().dofs_per_line+i];
-}
-
-
-
-template <>
-inline
-unsigned int
-DoFObjectAccessor<1,DoFHandler<3> >::dof_index (const unsigned int i,
-                                               const unsigned int fe_index) const
-{
-  Assert (fe_index == DoFHandler<3>::default_fe_index,
-         ExcMessage ("Only the default FE index is allowed for non-hp DoFHandler objects"));
-  Assert (static_cast<unsigned int>(this->present_level) < this->dof_handler->levels.size(),
-          ExcMessage ("DoFHandler not initialized"));
-
-  Assert (this->dof_handler != 0, BaseClass::ExcInvalidObject());
-                                  // make sure a FE has been selected
-                                  // and enough room was reserved
-  Assert (&this->dof_handler->get_fe() != 0, BaseClass::ExcInvalidObject());
-  Assert (i<this->dof_handler->get_fe().dofs_per_line,
-         ExcIndexRange (i, 0, this->dof_handler->get_fe().dofs_per_line));
-
-  return this->dof_handler->levels[this->present_level]
-    ->line_dofs[this->present_index*this->dof_handler->get_fe().dofs_per_line+i];
-}
-
-
-
-/*--------------- Functions: DoFObjectAccessor<2,DoFHandler> -----------*/
-
-template <>
-inline
-unsigned int
-DoFObjectAccessor<2,DoFHandler<2> >::dof_index (const unsigned int i,
-                                               const unsigned int fe_index) const
-{
-  Assert (fe_index == DoFHandler<2>::default_fe_index,
-         ExcMessage ("Only the default FE index is allowed for non-hp DoFHandler objects"));
-  Assert (this->dof_handler != 0,
-         BaseClass::ExcInvalidObject());
-                                  // make sure a FE has been selected
-                                  // and enough room was reserved
-  Assert (&this->dof_handler->get_fe() != 0,
-         BaseClass::ExcInvalidObject());
-  Assert (i<this->dof_handler->get_fe().dofs_per_quad,
-         ExcIndexRange (i, 0, this->dof_handler->get_fe().dofs_per_quad));
-  Assert (static_cast<unsigned int>(this->present_level) < this->dof_handler->levels.size(),
-          ExcMessage ("DoFHandler not initialized"));
-
-  return this->dof_handler->levels[this->present_level]
-    ->quad_dofs[this->present_index*this->dof_handler->get_fe().dofs_per_quad+i];
-}
-
-
-
-template <>
-inline
-unsigned int
-DoFObjectAccessor<2,DoFHandler<3> >::dof_index (const unsigned int i,
-                                               const unsigned int fe_index) const
-{
-  Assert (fe_index == DoFHandler<3>::default_fe_index,
-         ExcMessage ("Only the default FE index is allowed for non-hp DoFHandler objects"));
-  Assert (this->dof_handler != 0,
-         BaseClass::ExcInvalidObject());
-                                  // make sure a FE has been selected
-                                  // and enough room was reserved
-  Assert (&this->dof_handler->get_fe() != 0,
-         BaseClass::ExcInvalidObject());
-  Assert (i<this->dof_handler->get_fe().dofs_per_quad,
-         ExcIndexRange (i, 0, this->dof_handler->get_fe().dofs_per_quad));
-  Assert (static_cast<unsigned int>(this->present_level) < this->dof_handler->levels.size(),
-          ExcMessage ("DoFHandler not initialized"));
-
-  return this->dof_handler->levels[this->present_level]
-    ->quad_dofs[this->present_index*this->dof_handler->get_fe().dofs_per_quad+i];
-}
-
-
-
-/*--------------- Functions: DoFObjectAccessor<3,DoFHandler> -----------*/
-
-
-
-
-template <>
-inline
-unsigned int
-DoFObjectAccessor<3,DoFHandler<3> >::dof_index (const unsigned int i,
-                                               const unsigned int fe_index) const
-{
-  Assert (fe_index == DoFHandler<3>::default_fe_index,
-         ExcMessage ("Only the default FE index is allowed for non-hp DoFHandler objects"));
-  Assert (this->dof_handler != 0,
-         BaseClass::ExcInvalidObject());
-                                  // make sure a FE has been selected
-                                  // and enough room was reserved
-  Assert (&this->dof_handler->get_fe() != 0,
-         BaseClass::ExcInvalidObject());
-  Assert (i<this->dof_handler->get_fe().dofs_per_hex,
-         ExcIndexRange (i, 0, this->dof_handler->get_fe().dofs_per_hex));
-  Assert (static_cast<unsigned int>(this->present_level) < this->dof_handler->levels.size(),
-          ExcMessage ("DoFHandler not initialized"));
-
-  return this->dof_handler->levels[this->present_level]
-    ->hex_dofs[this->present_index*this->dof_handler->get_fe().dofs_per_hex+i];
-}
-
-
-
-/* -------------- hp vertex dofs stuff ------------------------------------ */
-
-
-/*--------------- Functions: DoFObjectAccessor<1,hp::DoFHandler> -----------*/
-
-template <>
-inline
-unsigned int
-DoFObjectAccessor<1,hp::DoFHandler<1> >::dof_index (const unsigned int i,
-                                                   const unsigned int fe_index) const
-{
-  return this->dof_handler->levels[this->present_level]
-    ->get_line_dof_index (*this->dof_handler,
-                          this->present_index,
-                          fe_index,
-                          i);
-}
-
-
-
-template <>
-inline
-unsigned int
-DoFObjectAccessor<1,hp::DoFHandler<2> >::dof_index (const unsigned int i,
-                                                   const unsigned int fe_index) const
-{
-  return this->dof_handler->levels[this->present_level]
-    ->get_line_dof_index (*this->dof_handler,
-                          this->present_index,
-                          fe_index,
-                          i);
-}
-
-
-
-template <>
-inline
-unsigned int
-DoFObjectAccessor<1,hp::DoFHandler<3> >::dof_index (const unsigned int i,
-                                                   const unsigned int fe_index) const
-{
-  return this->dof_handler->levels[this->present_level]
-    ->get_line_dof_index (*this->dof_handler,
-                          this->present_index,
-                          fe_index,
-                          i);
-}
-
-
-
-template <>
-inline
-void
-DoFObjectAccessor<1,hp::DoFHandler<1> >::set_dof_index (const unsigned int i,
-                                                       const unsigned int index,
-                                                       const unsigned int fe_index) const
-{
-  this->dof_handler->levels[this->present_level]
-    ->set_line_dof_index (*this->dof_handler,
-                          this->present_index,
-                          fe_index,
-                          i,
-                          index);
-}
-
-
-
-template <>
-inline
-void
-DoFObjectAccessor<1, hp::DoFHandler<2> >::set_dof_index (const unsigned int i,
-                                                        const unsigned int index,
-                                                        const unsigned int fe_index) const
-{
-  this->dof_handler->levels[this->present_level]
-    ->set_line_dof_index (*this->dof_handler,
-                          this->present_index,
-                          fe_index,
-                          i,
-                          index);
-}
-
-
-
-template <>
-inline
-void
-DoFObjectAccessor<1, hp::DoFHandler<3> >::set_dof_index (const unsigned int i,
-                                                        const unsigned int index,
-                                                        const unsigned int fe_index) const
-{
-  this->dof_handler->levels[this->present_level]
-    ->set_line_dof_index (*this->dof_handler,
-                          this->present_index,
-                          fe_index,
-                          i,
-                          index);
-}
-
-
-
-
-/*------------- Functions: DoFObjectAccessor<2,hp::DoFHandler> -------------*/
-
-template <>
-inline
-unsigned int DoFObjectAccessor<2,hp::DoFHandler<2> >::dof_index (const unsigned int i,
-                                                                const unsigned int fe_index) const
-{
-  return this->dof_handler->levels[this->present_level]
-    ->get_quad_dof_index (*this->dof_handler,
-                          this->present_index,
-                          fe_index,
-                          i);
-}
-
-
-
-template <>
-inline
-unsigned int DoFObjectAccessor<2,hp::DoFHandler<3> >::dof_index (const unsigned int i,
-                                                                const unsigned int fe_index) const
-{
-  return this->dof_handler->levels[this->present_level]
-    ->get_quad_dof_index (*this->dof_handler,
-                          this->present_index,
-                          fe_index,
-                          i);
-}
-
-
-
-template <>
-inline
-void
-DoFObjectAccessor<2, hp::DoFHandler<2> >::set_dof_index (const unsigned int i,
-                                                        const unsigned int index,
-                                                        const unsigned int fe_index) const
-{
-  this->dof_handler->levels[this->present_level]
-    ->set_quad_dof_index (*this->dof_handler,
-                          this->present_index,
-                          fe_index,
-                          i,
-                          index);
-}
-
-
-
-template <>
-inline
-void
-DoFObjectAccessor<2, hp::DoFHandler<3> >::set_dof_index (const unsigned int i,
-                                                        const unsigned int index,
-                                                        const unsigned int fe_index) const
-{
-  this->dof_handler->levels[this->present_level]
-    ->set_quad_dof_index (*this->dof_handler,
-                          this->present_index,
-                          fe_index,
-                          i,
-                          index);
-}
-
-
-
-
-/*------------- Functions: DoFObjectAccessor<3,hp::DoFHandler> -------------*/
-
-template <>
-inline
-unsigned int
-DoFObjectAccessor<3,hp::DoFHandler<3> >::dof_index (const unsigned int i,
-                                                   const unsigned int fe_index) const
-{
-  return this->dof_handler->levels[this->present_level]
-    ->get_hex_dof_index (*this->dof_handler,
-                         this->present_index,
-                         fe_index,
-                         i);
-}
-
-
-
-template <>
-inline
-void
-DoFObjectAccessor<3, hp::DoFHandler<3> >::set_dof_index (const unsigned int i,
-                                                        const unsigned int index,
-                                                        const unsigned int fe_index) const
-{
-  this->dof_handler->levels[this->present_level]
-    ->set_hex_dof_index (*this->dof_handler,
-                         this->present_index,
-                         fe_index,
-                         i,
-                         index);
-}
-
-
 
 /*------------------------- Functions: DoFCellAccessor -----------------------*/
 
index a361c9a0fde4d4a5c5affc8d330de97aadec435c..2e1a3f1ca42d20e53bccb48a143fd4fd49df69bb 100644 (file)
@@ -17,6 +17,8 @@
 #include <base/config.h>
 #include <vector>
 
+template <int> class DoFHandler;
+
 
 namespace internal
 {
@@ -141,6 +143,21 @@ namespace internal
                                           */
         std::vector<unsigned int> line_dofs;
 
+        template <int dim>
+        void
+        set_line_dof_index (const ::DoFHandler<dim> &dof_handler,
+                            const unsigned int       line_index,
+                            const unsigned int       fe_index,
+                            const unsigned int       local_index,
+                            const unsigned int       global_index);
+
+        template <int dim>
+        unsigned int
+        get_line_dof_index (const ::DoFHandler<dim> &dof_handler,
+                            const unsigned int       line_index,
+                            const unsigned int       fe_index,
+                            const unsigned int       local_index) const;
+
                                          /**
                                           * Determine an estimate for the
                                           * memory consumption (in bytes)
@@ -170,6 +187,21 @@ namespace internal
                                           */
         std::vector<unsigned int> quad_dofs;
 
+        template <int dim>
+        void
+        set_quad_dof_index (const ::DoFHandler<dim> &dof_handler,
+                            const unsigned int       quad_index,
+                            const unsigned int       fe_index,
+                            const unsigned int       local_index,
+                            const unsigned int       global_index);
+
+        template <int dim>
+        unsigned int
+        get_quad_dof_index (const ::DoFHandler<dim> &dof_handler,
+                            const unsigned int       quad_index,
+                            const unsigned int       fe_index,
+                            const unsigned int       local_index) const;
+
                                          /**
                                           * Determine an estimate for the
                                           * memory consumption (in bytes)
@@ -199,6 +231,21 @@ namespace internal
                                           */
         std::vector<unsigned int> hex_dofs;
 
+        template <int dim>
+        void
+        set_hex_dof_index (const ::DoFHandler<dim> &dof_handler,
+                           const unsigned int       hex_index,
+                           const unsigned int       fe_index,
+                           const unsigned int       local_index,
+                           const unsigned int       global_index);
+
+        template <int dim>
+        unsigned int
+        get_hex_dof_index (const ::DoFHandler<dim> &dof_handler,
+                           const unsigned int       hex_index,
+                           const unsigned int       fe_index,
+                           const unsigned int       local_index) const;
+
                                          /**
                                           * Determine an estimate for the
                                           * memory consumption (in bytes)
index 42fc106fad52fa67cc46bdeecc1f897fa34e7ed0..e7cd7768714a5391f586940d9dbbea4dcdc92d80 100644 (file)
@@ -295,6 +295,7 @@ namespace internal
 
 
     
+    
   } // namespace hp
   
 } // namespace internal
index e18099bb2c84a1c42a5012e427f3f8108e5866f9..7bee86dcbb71535320ed5811f6b76c2b0f618683 100644 (file)
 #include <vector>
 
 
-/*------------------------- Functions: DoFObjectAccessor<1,dim> -----------------------*/
-
-
-template <class DH>
-void DoFObjectAccessor<1, DH>::set_dof_index (const unsigned int i,
-                                              const unsigned int index,
-                                             const unsigned int fe_index) const
-{
-  Assert (fe_index == DoFHandler<1>::default_fe_index,
-         ExcMessage ("Only the default FE index is allowed for non-hp DoFHandler objects"));
-  Assert (this->dof_handler != 0,
-         typename BaseClass::ExcInvalidObject());
-                                  // make sure a FE has been selected
-                                  // and enough room was reserved
-  Assert (&this->dof_handler->get_fe() != 0,
-         typename BaseClass::ExcInvalidObject());
-  Assert (i<this->dof_handler->get_fe().dofs_per_line,
-         ExcIndexRange (i, 0, this->dof_handler->get_fe().dofs_per_line));
-
-  this->dof_handler->levels[this->present_level]
-    ->line_dofs[this->present_index*this->dof_handler->get_fe().dofs_per_line+i]
-    = index;
-}
-
-
-
-
-/*------------------------- Functions: DoFObjectAccessor<2,dim> -----------------------*/
-
-
-template <class DH>
-void DoFObjectAccessor<2, DH>::set_dof_index (const unsigned int i,
-                                              const unsigned int index,
-                                             const unsigned int fe_index) const
-{
-  Assert (fe_index == DoFHandler<1>::default_fe_index,
-         ExcMessage ("Only the default FE index is allowed for non-hp DoFHandler objects"));
-  Assert (this->dof_handler != 0,
-         typename BaseClass::ExcInvalidObject());
-                                  // make sure a FE has been selected
-                                  // and enough room was reserved
-  Assert (&this->dof_handler->get_fe() != 0,
-         typename BaseClass::ExcInvalidObject());
-  Assert (i<this->dof_handler->get_fe().dofs_per_quad,
-         ExcIndexRange (i, 0, this->dof_handler->get_fe().dofs_per_quad));
-
-  this->dof_handler->levels[this->present_level]
-    ->quad_dofs[this->present_index*this->dof_handler->get_fe().dofs_per_quad+i] = index;
-}
-
-
-
-
-
-/*------------------------- Functions: DoFObjectAccessor<3,dim> -----------------------*/
-
-
-template <class DH>
-void DoFObjectAccessor<3, DH>::set_dof_index (const unsigned int i,
-                                              const unsigned int index,
-                                             const unsigned int fe_index) const
-{
-  Assert (fe_index == DoFHandler<1>::default_fe_index,
-         ExcMessage ("Only the default FE index is allowed for non-hp DoFHandler objects"));
-  Assert (this->dof_handler != 0,
-         typename BaseClass::ExcInvalidObject());
-                                  // make sure a FE has been selected
-                                  // and enough room was reserved
-  Assert (&this->dof_handler->get_fe() != 0,
-         typename BaseClass::ExcInvalidObject());
-  Assert (i<this->dof_handler->get_fe().dofs_per_hex,
-         ExcIndexRange (i, 0, this->dof_handler->get_fe().dofs_per_hex));
-
-  this->dof_handler->levels[this->present_level]
-    ->hex_dofs[this->present_index*this->dof_handler->get_fe().dofs_per_hex+i] = index;
-}
-
-
 /*------------------------- Functions: DoFCellAccessor -----------------------*/
 
 
index 2f97a2661b7cb822ec52ca54ecc5b50f1d6aac79..3d2cef25383c67aad493e91d73a55eeb058668ef 100644 (file)
 //---------------------------------------------------------------------------
 
 
+#include <base/exceptions.h>
 #include <base/memory_consumption.h>
 #include <dofs/dof_levels.h>
+#include <dofs/dof_handler.h>
+#include <fe/fe.h>
 
 
 namespace internal
@@ -52,6 +55,203 @@ namespace internal
       return (DoFLevel<2>::memory_consumption () +
               MemoryConsumption::memory_consumption (hex_dofs));
     }
+
+
+
+    template <int dim>
+    unsigned int
+    DoFLevel<1>::
+    get_line_dof_index (const ::DoFHandler<dim> &dof_handler,
+                        const unsigned int       line_index,
+                        const unsigned int       fe_index,
+                        const unsigned int       local_index) const
+    {
+      Assert (fe_index == ::DoFHandler<dim>::default_fe_index,
+             ExcMessage ("Only the default FE index is allowed for non-hp DoFHandler objects"));
+      Assert (local_index<dof_handler.get_fe().dofs_per_line,
+             ExcIndexRange (local_index, 0, dof_handler.get_fe().dofs_per_line));
+      Assert (line_index * dof_handler.get_fe().dofs_per_line+local_index
+             <
+             line_dofs.size(),
+             ExcInternalError());
+      
+      return line_dofs[line_index * dof_handler.get_fe().dofs_per_line+local_index];
+    }
+
+
+
+    template <int dim>
+    void
+    DoFLevel<1>::
+    set_line_dof_index (const ::DoFHandler<dim> &dof_handler,
+                        const unsigned int       line_index,
+                        const unsigned int       fe_index,
+                        const unsigned int       local_index,
+                        const unsigned int       global_index)
+    {
+      Assert (fe_index == ::DoFHandler<dim>::default_fe_index,
+             ExcMessage ("Only the default FE index is allowed for non-hp DoFHandler objects"));
+      Assert (local_index<dof_handler.get_fe().dofs_per_line,
+             ExcIndexRange (local_index, 0, dof_handler.get_fe().dofs_per_line));
+      Assert (line_index * dof_handler.get_fe().dofs_per_line+local_index
+             <
+             line_dofs.size(),
+             ExcInternalError());
+      
+      line_dofs[line_index * dof_handler.get_fe().dofs_per_line+local_index] = global_index;
+    }
+
+
+
+    template <int dim>
+    unsigned int
+    DoFLevel<2>::
+    get_quad_dof_index (const ::DoFHandler<dim> &dof_handler,
+                        const unsigned int       quad_index,
+                        const unsigned int       fe_index,
+                        const unsigned int       local_index) const
+    {
+      Assert (fe_index == ::DoFHandler<dim>::default_fe_index,
+             ExcMessage ("Only the default FE index is allowed for non-hp DoFHandler objects"));
+      Assert (local_index<dof_handler.get_fe().dofs_per_quad,
+             ExcIndexRange (local_index, 0, dof_handler.get_fe().dofs_per_quad));
+      Assert (quad_index * dof_handler.get_fe().dofs_per_quad+local_index
+             <
+             quad_dofs.size(),
+             ExcInternalError());
+      
+      return quad_dofs[quad_index * dof_handler.get_fe().dofs_per_quad+local_index];
+    }
+
+
+
+    template <int dim>
+    void
+    DoFLevel<2>::
+    set_quad_dof_index (const ::DoFHandler<dim> &dof_handler,
+                        const unsigned int       quad_index,
+                        const unsigned int       fe_index,
+                        const unsigned int       local_index,
+                        const unsigned int       global_index)
+    {
+      Assert (fe_index == ::DoFHandler<dim>::default_fe_index,
+             ExcMessage ("Only the default FE index is allowed for non-hp DoFHandler objects"));
+      Assert (local_index<dof_handler.get_fe().dofs_per_quad,
+             ExcIndexRange (local_index, 0, dof_handler.get_fe().dofs_per_quad));
+      Assert (quad_index * dof_handler.get_fe().dofs_per_quad+local_index
+             <
+             quad_dofs.size(),
+             ExcInternalError());
+      
+      quad_dofs[quad_index * dof_handler.get_fe().dofs_per_quad+local_index] = global_index;
+    }
+
+
+
+    template <int dim>
+    unsigned int
+    DoFLevel<3>::
+    get_hex_dof_index (const ::DoFHandler<dim> &dof_handler,
+                       const unsigned int       hex_index,
+                       const unsigned int       fe_index,
+                       const unsigned int       local_index) const
+    {
+      Assert (fe_index == ::DoFHandler<dim>::default_fe_index,
+             ExcMessage ("Only the default FE index is allowed for non-hp DoFHandler objects"));
+      Assert (local_index<dof_handler.get_fe().dofs_per_hex,
+             ExcIndexRange (local_index, 0, dof_handler.get_fe().dofs_per_hex));
+      Assert (hex_index * dof_handler.get_fe().dofs_per_hex+local_index
+             <
+             hex_dofs.size(),
+             ExcInternalError());
+      
+      return hex_dofs[hex_index * dof_handler.get_fe().dofs_per_hex+local_index];
+    }
+
+
+
+    template <int dim>
+    void
+    DoFLevel<3>::
+    set_hex_dof_index (const ::DoFHandler<dim> &dof_handler,
+                       const unsigned int       hex_index,
+                       const unsigned int       fe_index,
+                       const unsigned int       local_index,
+                       const unsigned int       global_index)
+    {
+      Assert (fe_index == ::DoFHandler<dim>::default_fe_index,
+             ExcMessage ("Only the default FE index is allowed for non-hp DoFHandler objects"));
+      Assert (local_index<dof_handler.get_fe().dofs_per_hex,
+             ExcIndexRange (local_index, 0, dof_handler.get_fe().dofs_per_hex));
+      Assert (hex_index * dof_handler.get_fe().dofs_per_hex+local_index
+             <
+             hex_dofs.size(),
+             ExcInternalError());
+      
+      hex_dofs[hex_index * dof_handler.get_fe().dofs_per_hex+local_index] = global_index;
+    }
+    
+    
+
+// explicit instantiations
+    template
+    unsigned int
+    DoFLevel<1>::
+    get_line_dof_index (const ::DoFHandler<deal_II_dimension> &dof_handler,
+                        const unsigned int       line_index,
+                        const unsigned int       fe_index,
+                        const unsigned int       local_index) const;
+    
+    template
+    void
+    DoFLevel<1>::
+    set_line_dof_index (const ::DoFHandler<deal_II_dimension> &dof_handler,
+                        const unsigned int       line_index,
+                        const unsigned int       fe_index,
+                        const unsigned int       local_index,
+                        const unsigned int       global_index);
+
+#if deal_II_dimension >= 2
+
+    template
+    unsigned int
+    DoFLevel<2>::
+    get_quad_dof_index (const ::DoFHandler<deal_II_dimension> &dof_handler,
+                        const unsigned int       quad_index,
+                        const unsigned int       fe_index,
+                        const unsigned int       local_index) const;
+    
+    template
+    void
+    DoFLevel<2>::
+    set_quad_dof_index (const ::DoFHandler<deal_II_dimension> &dof_handler,
+                        const unsigned int       quad_index,
+                        const unsigned int       fe_index,
+                        const unsigned int       local_index,
+                        const unsigned int       global_index);
+
+#endif
+
+#if deal_II_dimension >= 3
+    
+    template
+    unsigned int
+    DoFLevel<3>::
+    get_hex_dof_index (const ::DoFHandler<deal_II_dimension> &dof_handler,
+                       const unsigned int       hex_index,
+                       const unsigned int       fe_index,
+                       const unsigned int       local_index) const;
+    
+    template
+    void
+    DoFLevel<3>::
+    set_hex_dof_index (const ::DoFHandler<deal_II_dimension> &dof_handler,
+                       const unsigned int       hex_index,
+                       const unsigned int       fe_index,
+                       const unsigned int       local_index,
+                       const unsigned int       global_index);
+
+#endif
     
   }
   

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.