]> https://gitweb.dealii.org/ - dealii.git/commitdiff
counting by blocks
authorGuido Kanschat <dr.guido.kanschat@gmail.com>
Sun, 29 Jan 2006 20:21:16 +0000 (20:21 +0000)
committerGuido Kanschat <dr.guido.kanschat@gmail.com>
Sun, 29 Jan 2006 20:21:16 +0000 (20:21 +0000)
git-svn-id: https://svn.dealii.org/trunk@12205 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/dofs/dof_tools.h
deal.II/deal.II/source/dofs/dof_tools.cc

index b4a3c7f86ad03ba16471ef24c6cdaa119db52423..dc7c48027a13542cb52508edf3de3b4304f4dfc9 100644 (file)
@@ -751,27 +751,29 @@ class DoFTools
                                      * Extract the indices of the
                                      * degrees of freedom belonging
                                      * to certain vector components
-                                     * of a vector-valued finite
-                                     * element. The bit vector
-                                     * @p select defines, which
-                                     * components of an
-                                     * FESystem are to be
-                                     * extracted from the
-                                     * DoFHandler @p dof. The
-                                     * entries in @p selected_dofs
-                                     * corresponding to degrees of
-                                     * freedom belonging to these
-                                     * components are then flagged
-                                     * @p true, while all others are
-                                     * set to @p false.
+                                     * or blocks (if the last
+                                     * argument is <tt>true</tt>) of
+                                     * a vector-valued finite
+                                     * element. The bit vector @p
+                                     * select defines, which
+                                     * components or blocks of an
+                                     * FESystem are to be extracted
+                                     * from the DoFHandler @p
+                                     * dof. The entries in @p
+                                     * selected_dofs corresponding to
+                                     * degrees of freedom belonging
+                                     * to these components are then
+                                     * flagged @p true, while all
+                                     * others are set to @p false.
                                      *
-                                     * The size of
-                                     * @p component_select shall
+                                     * The size of @p select must
                                      * equal the number of components
-                                     * in the finite element used by
-                                     * @p dof. The size of
-                                     * @p selected_dofs shall equal
-                                     * <tt>dof_handler.n_dofs()</tt>. Previous
+                                     * or blocks in the FiniteElement
+                                     * used by @p dof, depending on
+                                     * the argument
+                                     * <tt>blocks</tt>. The size of
+                                     * @p selected_dofs must equal
+                                     * DoFHandler::n_dofs(). Previous
                                      * contents of this array are
                                      * overwritten.
                                      *
@@ -781,8 +783,8 @@ class DoFTools
                                      * of its shape functions are
                                      * non-zero in more than one
                                      * vector component (which holds,
-                                     * for example, for Nedelec or
-                                     * Raviart-Thomas elements), then
+                                     * for example, for FE_Nedelec or
+                                     * FE_RaviartThomas elements), then
                                      * shape functions cannot be
                                      * associated with a single
                                      * vector component. In this
@@ -799,12 +801,13 @@ class DoFTools
     template <int dim>
     static void
     extract_dofs (const DoFHandler<dim>   &dof_handler,
-                 const std::vector<bool> &component_select,
-                 std::vector<bool>       &selected_dofs);
+                 const std::vector<bool> &select,
+                 std::vector<bool>       &selected_dofs,
+                 bool blocks = false);
 
                                     /**
                                      * Do the same thing as
-                                     * @p extract_dofs for one level
+                                     * extract_dofs() for one level
                                      * of a multi-grid DoF numbering.
                                      */
     template <int dim>
@@ -812,7 +815,8 @@ class DoFTools
     extract_level_dofs (const unsigned int       level,
                        const MGDoFHandler<dim> &dof,
                        const std::vector<bool> &select,
-                       std::vector<bool>       &selected_dofs);
+                       std::vector<bool>       &selected_dofs,
+                       bool blocks = false);
 
                                     /**
                                      * Extract all degrees of freedom
index 0c345d6bd0778e15597c44605ac54afd02ed8c3e..c944344a1fa747519a997d81a83ad824bed66fd4 100644 (file)
@@ -1732,12 +1732,25 @@ void DoFTools::distribute_cell_to_dof_vector (
 
 template <int dim>
 void
-DoFTools::extract_dofs (const DoFHandler<dim>   &dof,
-                       const std::vector<bool> &component_select,
-                       std::vector<bool>       &selected_dofs)
+DoFTools::extract_dofs (
+  const DoFHandler<dim>   &dof,
+  const std::vector<bool> &component_select,
+  std::vector<bool>       &selected_dofs,
+  bool count_by_blocks)
 {
-  Assert(component_select.size() == n_components(dof),
-        ExcDimensionMismatch(component_select.size(), n_components(dof)));
+  const FiniteElement<dim> &fe = dof.get_fe();
+  
+  if (count_by_blocks)
+    {
+      Assert(component_select.size() == fe.n_blocks(),
+            ExcDimensionMismatch(component_select.size(), fe.n_blocks()));
+    }
+  else
+    {
+      Assert(component_select.size() == n_components(dof),
+            ExcDimensionMismatch(component_select.size(), n_components(dof)));
+    }
+  
   Assert(selected_dofs.size() == dof.n_dofs(),
         ExcDimensionMismatch(selected_dofs.size(), dof.n_dofs()));
 
@@ -1759,8 +1772,6 @@ DoFTools::extract_dofs (const DoFHandler<dim>   &dof,
     };
   
 
-  const FiniteElement<dim> &fe = dof.get_fe();
-  
                                   // preset all values by false
   std::fill_n (selected_dofs.begin(), dof.n_dofs(), false);
 
@@ -1770,50 +1781,54 @@ DoFTools::extract_dofs (const DoFHandler<dim>   &dof,
                                    // something interesting or not
   std::vector<bool> local_selected_dofs (fe.dofs_per_cell, false);
   for (unsigned int i=0; i<fe.dofs_per_cell; ++i)
-    if (fe.is_primitive(i))
+    if (count_by_blocks)
       local_selected_dofs[i]
-        = component_select[fe.system_to_component_index(i).first];
+        = component_select[fe.system_to_block_index(i).first];
     else
-                                       // if this shape function is
-                                       // not primitive, then we have
-                                       // to work harder. we have to
-                                       // find out whether _any_ of
-                                       // the vector components of
-                                       // this element is selected or
-                                       // not
-                                       //
-                                       // to do so, get the first and
-                                       // last vector components of
-                                       // the base element to which
-                                       // the local dof with index i
-                                       // belongs
-      {
-        unsigned int first_comp = 0;
-        const unsigned int this_base = fe.system_to_base_index(i).first.first;
-        const unsigned int this_multiplicity
-          = fe.system_to_base_index(i).first.second;
-        
-        for (unsigned int b=0; b<this_base; ++b)
-          first_comp += fe.base_element(b).n_components() *
-                        fe.element_multiplicity(b);
-        for (unsigned int m=0; m<this_multiplicity; ++m)
-          first_comp += fe.base_element(this_base).n_components();
-        const unsigned int end_comp = first_comp +
-                                      fe.base_element(this_base).n_components();
-
-        Assert (first_comp < fe.n_components(), ExcInternalError());
-        Assert (end_comp <= fe.n_components(),  ExcInternalError());
-
-                                         // now check whether any of
-                                         // the components in between
-                                         // is set
-        for (unsigned int c=first_comp; c<end_comp; ++c)
-          if (component_select[c] == true)
-            {
-              local_selected_dofs[i] = true;
-              break;
-            };
-      };  
+      if (fe.is_primitive(i))
+       local_selected_dofs[i]
+         = component_select[fe.system_to_component_index(i).first];
+      else
+                                        // if this shape function is
+                                        // not primitive, then we have
+                                        // to work harder. we have to
+                                        // find out whether _any_ of
+                                        // the vector components of
+                                        // this element is selected or
+                                        // not
+                                        //
+                                        // to do so, get the first and
+                                        // last vector components of
+                                        // the base element to which
+                                        // the local dof with index i
+                                        // belongs
+       {
+         unsigned int first_comp = 0;
+         const unsigned int this_base = fe.system_to_base_index(i).first.first;
+         const unsigned int this_multiplicity
+           = fe.system_to_base_index(i).first.second;
+         
+         for (unsigned int b=0; b<this_base; ++b)
+           first_comp += fe.base_element(b).n_components() *
+                         fe.element_multiplicity(b);
+         for (unsigned int m=0; m<this_multiplicity; ++m)
+           first_comp += fe.base_element(this_base).n_components();
+         const unsigned int end_comp = first_comp +
+                                       fe.base_element(this_base).n_components();
+         
+         Assert (first_comp < fe.n_components(), ExcInternalError());
+         Assert (end_comp <= fe.n_components(),  ExcInternalError());
+         
+                                          // now check whether any of
+                                          // the components in between
+                                          // is set
+         for (unsigned int c=first_comp; c<end_comp; ++c)
+           if (component_select[c] == true)
+             {
+               local_selected_dofs[i] = true;
+               break;
+             }
+       }
   
                                    // then loop over all cells and do
                                    // the work
@@ -1830,14 +1845,26 @@ DoFTools::extract_dofs (const DoFHandler<dim>   &dof,
 
 template<int dim>
 void
-DoFTools::extract_level_dofs(const unsigned int       level,
-                            const MGDoFHandler<dim> &dof,
-                            const std::vector<bool> &component_select,
-                            std::vector<bool>       &selected_dofs)
+DoFTools::extract_level_dofs(
+  const unsigned int       level,
+  const MGDoFHandler<dim> &dof,
+  const std::vector<bool> &component_select,
+  std::vector<bool>       &selected_dofs,
+  bool count_by_blocks)
 {
   const FiniteElement<dim>& fe = dof.get_fe();
-  Assert(component_select.size() == fe.n_components(),
-        ExcDimensionMismatch(component_select.size(), fe.n_components()));
+  
+  if (count_by_blocks)
+    {
+      Assert(component_select.size() == fe.n_blocks(),
+            ExcDimensionMismatch(component_select.size(), fe.n_blocks()));
+    }
+  else
+    {
+      Assert(component_select.size() == fe.n_components(),
+            ExcDimensionMismatch(component_select.size(), fe.n_components()));
+    }
+  
   Assert(selected_dofs.size() == dof.n_dofs(level),
         ExcDimensionMismatch(selected_dofs.size(), dof.n_dofs(level)));
 
@@ -1867,51 +1894,55 @@ DoFTools::extract_level_dofs(const unsigned int       level,
                                    // something interesting or not
   std::vector<bool> local_selected_dofs (fe.dofs_per_cell, false);
   for (unsigned int i=0; i<fe.dofs_per_cell; ++i)
-    if (fe.is_primitive(i))
+    if (count_by_blocks)
       local_selected_dofs[i]
-        = component_select[fe.system_to_component_index(i).first];
+        = component_select[fe.system_to_block_index(i).first];
     else
-                                       // if this shape function is
-                                       // not primitive, then we have
-                                       // to work harder. we have to
-                                       // find out whether _any_ of
-                                       // the vector components of
-                                       // this element is selected or
-                                       // not
-                                       //
-                                       // to do so, get the first and
-                                       // last vector components of
-                                       // the base element to which
-                                       // the local dof with index i
-                                       // belongs
-      {
-        unsigned int first_comp = 0;
-        const unsigned int this_base = fe.system_to_base_index(i).first.first;
-        const unsigned int this_multiplicity
-          = fe.system_to_base_index(i).first.second;
-        
-        for (unsigned int b=0; b<this_base; ++b)
-          first_comp += fe.base_element(b).n_components() *
-                        fe.element_multiplicity(b);
-        for (unsigned int m=0; m<this_multiplicity; ++m)
-          first_comp += fe.base_element(this_base).n_components();
-        const unsigned int end_comp = first_comp +
-                                      fe.base_element(this_base).n_components();
-
-        Assert (first_comp < fe.n_components(), ExcInternalError());
-        Assert (end_comp <= fe.n_components(),  ExcInternalError());
-
-                                         // now check whether any of
-                                         // the components in between
-                                         // is set
-        for (unsigned int c=first_comp; c<end_comp; ++c)
-          if (component_select[c] == true)
-            {
-              local_selected_dofs[i] = true;
-              break;
-            };
-      };
-
+      if (fe.is_primitive(i))
+       local_selected_dofs[i]
+         = component_select[fe.system_to_component_index(i).first];
+      else
+                                        // if this shape function is
+                                        // not primitive, then we have
+                                        // to work harder. we have to
+                                        // find out whether _any_ of
+                                        // the vector components of
+                                        // this element is selected or
+                                        // not
+                                        //
+                                        // to do so, get the first and
+                                        // last vector components of
+                                        // the base element to which
+                                        // the local dof with index i
+                                        // belongs
+       {
+         unsigned int first_comp = 0;
+         const unsigned int this_base = fe.system_to_base_index(i).first.first;
+         const unsigned int this_multiplicity
+           = fe.system_to_base_index(i).first.second;
+         
+         for (unsigned int b=0; b<this_base; ++b)
+           first_comp += fe.base_element(b).n_components() *
+                         fe.element_multiplicity(b);
+         for (unsigned int m=0; m<this_multiplicity; ++m)
+           first_comp += fe.base_element(this_base).n_components();
+         const unsigned int end_comp = first_comp +
+                                       fe.base_element(this_base).n_components();
+         
+         Assert (first_comp < fe.n_components(), ExcInternalError());
+         Assert (end_comp <= fe.n_components(),  ExcInternalError());
+         
+                                          // now check whether any of
+                                          // the components in between
+                                          // is set
+         for (unsigned int c=first_comp; c<end_comp; ++c)
+           if (component_select[c] == true)
+             {
+               local_selected_dofs[i] = true;
+               break;
+             }
+       }
+  
                                    // then loop over all cells and do
                                    // work
   std::vector<unsigned int> indices(fe.dofs_per_cell);  
@@ -2465,11 +2496,12 @@ DoFTools::count_dofs_per_component (
     {
       void (*fun_ptr) (const DoFHandler<dim>   &,
                       const std::vector<bool> &,
-                      std::vector<bool>       &)
+                      std::vector<bool>       &,
+                      bool)
         = &DoFTools::template extract_dofs<dim>;
       component_select[i][i] = true;
       threads += Threads::spawn (fun_ptr)(dof_handler, component_select[i],
-                                          dofs_in_component[i]);
+                                          dofs_in_component[i], false);
     };
   threads.join_all ();
 
@@ -4007,15 +4039,12 @@ DoFTools::distribute_cell_to_dof_vector<deal_II_dimension>
 
 
 template void DoFTools::extract_dofs<deal_II_dimension>
-(const DoFHandler<deal_II_dimension>& dof,
- const std::vector<bool>& component_select,
- std::vector<bool>& selected_dofs);
+(const DoFHandler<deal_II_dimension>&,
+ const std::vector<bool>&, std::vector<bool>&, bool);
 
 template void DoFTools::extract_level_dofs<deal_II_dimension>
-(const unsigned int level,
- const MGDoFHandler<deal_II_dimension>& dof,
- const std::vector<bool>& component_select,
- std::vector<bool>& selected_dofs);
+(const unsigned int level, const MGDoFHandler<deal_II_dimension>&,
+ const std::vector<bool>&, std::vector<bool>&, bool);
 
 #if deal_II_dimension != 1
 template

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.