]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
put these auxiliary functions into a reasonable place
authorguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 29 Jun 2006 12:11:27 +0000 (12:11 +0000)
committerguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 29 Jun 2006 12:11:27 +0000 (12:11 +0000)
instead of duplicating them

git-svn-id: https://svn.dealii.org/trunk@13311 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/dofs/dof_tools.h
deal.II/deal.II/include/numerics/vectors.templates.h
deal.II/deal.II/source/dofs/dof_tools.cc
deal.II/deal.II/source/numerics/vectors.cc

index c264a55a952771abacc7fb69338f0443dbbc85e4..6de8092e5cb4e8e7cc522465052fd7299871f6ac 100644 (file)
@@ -196,6 +196,110 @@ class DoFTools
                                             */
          nonzero
     };
+    
+                                    /**
+                                     * @name Auxiliary functions
+                                     * @{
+                                     */
+                                    /**
+                                     * Maximal number of degrees of
+                                     * freedom on a cell. This is
+                                     * just
+                                     * FiniteElementData::dofs_per_cell,
+                                     * but allows for a common
+                                     * interface with hp::DoFHandler.
+                                     */
+    template <int dim>
+    static unsigned int
+    max_dofs_per_cell (const DoFHandler<dim> &dh);
+    
+    
+                                    /**
+                                     * Maximal number of degrees of
+                                     * freedom on a face. This is
+                                     * just
+                                     * FiniteElementData::dofs_per_face,
+                                     * but allows for a common
+                                     * interface with hp::DoFHandler.
+                                     */
+    template <int dim>
+    static unsigned int
+    max_dofs_per_face (const DoFHandler<dim> &dh);
+
+                                    /**
+                                     * Maximal number of degrees of
+                                     * freedom on a vertex. This is
+                                     * just
+                                     * FiniteElementData::dofs_per_vertex,
+                                     * but allows for a common
+                                     * interface with hp::DoFHandler.
+                                     */
+    template <int dim>
+    static unsigned int
+    max_dofs_per_vertex (const DoFHandler<dim> &dh);
+    
+                                    /**
+                                     * Number of components in an
+                                     * hp-conforming way.
+                                     */
+    template <int dim>
+    static unsigned int
+    n_components (const DoFHandler<dim> &dh);
+    
+                                    /**
+                                     * Find out if a FiniteElement is
+                                     * primitive in an hp-conforming
+                                     * way.
+                                     */
+    template <int dim>
+    static unsigned int
+    fe_is_primitive (const DoFHandler<dim> &dh);
+    
+                                    /**
+                                     * Maximal number of degrees of
+                                     * freedom on a cell in an hp hierarchy.
+                                     */
+    template <int dim>
+    static unsigned int
+    max_dofs_per_cell (const hp::DoFHandler<dim> &dh);
+    
+                                    /**
+                                     * Maximal number of degrees of
+                                     * freedom on a face in an hp hierarchy.
+                                     */
+    template <int dim>
+    static unsigned int
+    max_dofs_per_face (const hp::DoFHandler<dim> &dh);
+    
+                                    /**
+                                     *Maximal number of degrees of
+                                     * freedom on a vertex in an hp hierarchy.
+                                     */
+    template <int dim>
+    static unsigned int
+    max_dofs_per_vertex (const hp::DoFHandler<dim> &dh);
+    
+                                    /**
+                                     * Number of components in an
+                                     * hp-conforming way.
+                                     */
+    template <int dim>
+    static unsigned int
+    n_components (const hp::DoFHandler<dim> &dh);
+    
+                                    /**
+                                     * Find out if an hp::FECollection is
+                                     * primitive in an hp-conforming
+                                     * way.
+                                     */
+    template <int dim>
+    static unsigned int
+    fe_is_primitive (const hp::DoFHandler<dim> &dh);
+    
+                                    /**
+                                     * @}
+                                     */
+    
                                     /**
                                      * @name Sparsity Pattern Generation
                                      * @{
@@ -1644,6 +1748,87 @@ DoFTools::Coupling operator | (const DoFTools::Coupling c1,
 
 // ---------------------- inline and template functions --------------------
 
+template <int dim>
+inline unsigned int
+DoFTools::max_dofs_per_cell (const DoFHandler<dim> &dh)
+{
+  return dh.get_fe().dofs_per_cell;
+}
+
+
+template <int dim>
+inline unsigned int
+DoFTools::max_dofs_per_face (const DoFHandler<dim> &dh) 
+{
+  return dh.get_fe().dofs_per_face;
+}
+
+
+template <int dim>
+inline unsigned int
+DoFTools::max_dofs_per_vertex (const DoFHandler<dim> &dh) 
+{
+  return dh.get_fe().dofs_per_vertex;
+}
+
+
+template <int dim>
+inline unsigned int
+DoFTools::n_components (const DoFHandler<dim> &dh) 
+{
+  return dh.get_fe().n_components();
+}
+
+
+
+template <int dim>
+inline unsigned int
+DoFTools::fe_is_primitive (const DoFHandler<dim> &dh) 
+{
+  return dh.get_fe().is_primitive();
+}
+
+
+template <int dim>
+inline unsigned int
+DoFTools::max_dofs_per_cell (const hp::DoFHandler<dim> &dh) 
+{
+  return dh.get_fe().max_dofs_per_cell ();
+}
+
+
+template <int dim>
+inline unsigned int
+DoFTools::max_dofs_per_face (const hp::DoFHandler<dim> &dh) 
+{
+  return dh.get_fe().max_dofs_per_face ();
+}
+
+
+template <int dim>
+inline unsigned int
+DoFTools::max_dofs_per_vertex (const hp::DoFHandler<dim> &dh) 
+{
+  return dh.get_fe().max_dofs_per_vertex ();
+}
+
+
+template <int dim>
+inline unsigned int
+DoFTools::n_components (const hp::DoFHandler<dim> &dh) 
+{
+  return dh.get_fe()[0].n_components();
+}
+
+
+template <int dim>
+inline unsigned int
+DoFTools::fe_is_primitive (const hp::DoFHandler<dim> &dh) 
+{
+  return dh.get_fe()[0].is_primitive();
+}
+
+
 template <class DH, class SparsityPattern>
 inline
 void
index 67068eb5b894157dcc35769399089bfe20b2a422..33fdc7bca7f89e6499dfaf89597f8ff2080208b3 100644 (file)
@@ -896,7 +896,7 @@ interpolate_boundary_values (const Mapping<dim>            &mapping,
   Assert (function_map.find(255) == function_map.end(),
          ExcInvalidBoundaryIndicator());
 
-  const unsigned int        n_components = get_n_components(dof);
+  const unsigned int        n_components = DoFTools::n_components(dof);
   const bool                fe_is_system = (n_components != 1);
 
   for (typename FunctionMap<dim>::type::const_iterator i=function_map.begin();
@@ -915,11 +915,11 @@ interpolate_boundary_values (const Mapping<dim>            &mapping,
 
                                   // field to store the indices
   std::vector<unsigned int> face_dofs;
-  face_dofs.reserve (max_dofs_per_face(dof));
+  face_dofs.reserve (DoFTools::max_dofs_per_face(dof));
   std::fill (face_dofs.begin (), face_dofs.end (), DoFHandler<dim>::invalid_dof_index);
 
   std::vector<Point<dim> >  dof_locations;
-  dof_locations.reserve (max_dofs_per_face(dof));
+  dof_locations.reserve (DoFTools::max_dofs_per_face(dof));
   std::fill (dof_locations.begin(), dof_locations.end (), Point<dim>());
 
                                   // array to store the values of
@@ -930,8 +930,8 @@ interpolate_boundary_values (const Mapping<dim>            &mapping,
                                   // respectively
   std::vector<double>          dof_values_scalar;
   std::vector<Vector<double> > dof_values_system;
-  dof_values_scalar.reserve (max_dofs_per_face (dof));
-  dof_values_system.reserve (max_dofs_per_face (dof));
+  dof_values_scalar.reserve (DoFTools::max_dofs_per_face (dof));
+  dof_values_system.reserve (DoFTools::max_dofs_per_face (dof));
 
   typename DH<dim>::active_cell_iterator cell = dof.begin_active(),
                                         endc = dof.end();
index 00a6d85ddf3af673dc9eaa234839ab2e33f7758c..7f2f251e3910e19fac9871e936fef10ada719b92 100644 (file)
 #include <numeric>
 
 
-namespace 
-{
-  // Functions for the DoFHandler
-  template <int dim>
-  unsigned int
-  max_dofs_per_cell (const DoFHandler<dim> &dh)
-  {
-    return dh.get_fe().dofs_per_cell;
-  }
-
-
-
-  template <int dim>
-  unsigned int
-  max_dofs_per_face (const DoFHandler<dim> &dh) 
-  {
-    return dh.get_fe().dofs_per_face;
-  }
-
-
-
-  template <int dim>
-  unsigned int
-  max_dofs_per_vertex (const DoFHandler<dim> &dh) 
-  {
-    return dh.get_fe().dofs_per_vertex;
-  }
-
-
-
-  template <int dim>
-  unsigned int
-  n_components (const DoFHandler<dim> &dh) 
-  {
-    return dh.get_fe().n_components();
-  }
-
-
-
-  template <int dim>
-  unsigned int
-  fe_is_primitive (const DoFHandler<dim> &dh) 
-  {
-    return dh.get_fe().is_primitive();
-  }
-
-  // Functions for the hp::DoFHandler
-  template <int dim>
-  unsigned int
-  max_dofs_per_cell (const hp::DoFHandler<dim> &dh) 
-  {
-    return dh.get_fe().max_dofs_per_cell ();
-  }
-
-
-
-  template <int dim>
-  unsigned int
-  max_dofs_per_face (const hp::DoFHandler<dim> &dh) 
-  {
-    return dh.get_fe().max_dofs_per_face ();
-  }
-
-
-
-  template <int dim>
-  unsigned int
-  max_dofs_per_vertex (const hp::DoFHandler<dim> &dh) 
-  {
-    return dh.get_fe().max_dofs_per_vertex ();
-  }
-
-
-
-  template <int dim>
-  unsigned int
-  n_components (const hp::DoFHandler<dim> &dh) 
-  {
-//TODO:[?] Verify that this is really correct
-    return dh.get_fe()[0].n_components();
-  }
-
-
-
-  template <int dim>
-  unsigned int
-  fe_is_primitive (const hp::DoFHandler<dim> &dh) 
-  {
-//TODO:[?] Verify that this is really correct
-    return dh.get_fe()[0].is_primitive();
-  }
-}
-
-
-
 #if  deal_II_dimension == 1
 
 // Specialization for 1D
index d5d5e6f3a3a22d029d5cf64886600e3286bd1427..0f3a57ab376b8c761abd96d07b632fb9d2de9641 100644 (file)
 #include <dofs/dof_handler.h>
 #include <dofs/hp_dof_handler.h>
 
-namespace 
-{
-  // Functions for the DoFHandler
-  template <int dim>
-  unsigned int
-  max_dofs_per_cell (const DoFHandler<dim> &dh)
-  {
-    return dh.get_fe().dofs_per_cell;
-  }
-
-
-  template <int dim>
-  unsigned int
-  max_dofs_per_face (const DoFHandler<dim> &dh) 
-  {
-    return dh.get_fe().dofs_per_face;
-  }
-
-
-  template <int dim>
-  unsigned int
-  get_n_components (const DoFHandler<dim> &dh) 
-  {
-    return dh.get_fe().n_components();
-  }
-
-
-  // Functions for the hp::DoFHandler
-  template <int dim>
-  unsigned int
-  max_dofs_per_cell (const hp::DoFHandler<dim> &dh) 
-  {
-    return dh.get_fe().max_dofs_per_cell ();
-  }
-
-
-  template <int dim>
-  unsigned int
-  max_dofs_per_face (const hp::DoFHandler<dim> &dh) 
-  {
-    return dh.get_fe().max_dofs_per_face ();
-  }
-
-
-
-  template <int dim>
-  unsigned int
-  get_n_components (const hp::DoFHandler<dim> &dh) 
-  {
-    return dh.get_fe().n_components();
-  }
-}
-
-
 #include<numerics/vectors.templates.h>
 
 // explicit instantiations

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.