]> https://gitweb.dealii.org/ - dealii.git/commitdiff
MGTransferBlock changed to blocks instead of components
authorGuido Kanschat <dr.guido.kanschat@gmail.com>
Mon, 30 Jan 2006 14:58:09 +0000 (14:58 +0000)
committerGuido Kanschat <dr.guido.kanschat@gmail.com>
Mon, 30 Jan 2006 14:58:09 +0000 (14:58 +0000)
git-svn-id: https://svn.dealii.org/trunk@12211 0785d39b-7218-0410-832d-ea1e28bc413d

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

index dc7c48027a13542cb52508edf3de3b4304f4dfc9..0e9f26fa27314faf9977d9b05a7c853246542d33 100644 (file)
@@ -1089,6 +1089,24 @@ class DoFTools
                              std::vector<unsigned int>  target_component
                              = std::vector<unsigned int>());
 
+                                    /**
+                                     * Count the degrees of freedom
+                                     * in each block. This function
+                                     * is similar to
+                                     * count_dofs_per_component(),
+                                     * with the difference that the
+                                     * counting is done by
+                                     * blocks. See @ref GlossBlock
+                                     * "blocks" in the glossary for
+                                     * details.
+                                     */
+    template <int dim>
+    static void
+    count_dofs_per_block (const DoFHandler<dim>&     dof_handler,
+                         std::vector<unsigned int>& dofs_per_component,
+                         std::vector<unsigned int>  target_block
+                         = std::vector<unsigned int>());
+    
                                     /**
                                      * @deprecated See the previous
                                      * function with the same name
index b00f7c747c6af872ce64f51d35b24a2a592389db..b13e0534940437b39ef33d0a4d9fa7c8c96ead7e 100644 (file)
@@ -2318,8 +2318,8 @@ inline
 std::pair<unsigned int,unsigned int>
 FiniteElement<dim>::system_to_block_index (const unsigned int index) const
 {
-  Assert (index < this->n_blocks(),
-        ExcIndexRange(index, 0, this->n_blocks()));
+  Assert (index < this->dofs_per_cell,
+         ExcIndexRange(index, 0, this->dofs_per_cell));
                                   // The block is computed simply as
                                   // first block of this base plus
                                   // the index within the base blocks
index c944344a1fa747519a997d81a83ad824bed66fd4..f6172f7a1518e40b7ac184cde43e231707218f0c 100644 (file)
@@ -2540,6 +2540,70 @@ DoFTools::count_dofs_per_component (
 }
 
 
+template <int dim>
+void
+DoFTools::count_dofs_per_block (
+  const DoFHandler<dim>&     dof_handler,
+  std::vector<unsigned int>& dofs_per_block,
+  std::vector<unsigned int>  target_block)
+{
+  const FiniteElement<dim>& fe = dof_handler.get_fe();
+  const unsigned int n_blocks = fe.n_blocks();
+  dofs_per_block.resize (n_blocks);
+  std::fill (dofs_per_block.begin(), dofs_per_block.end(), 0U);
+  
+                                  // If the empty vector was given as
+                                  // default argument, set up this
+                                  // vector as identity.
+  if (target_block.size()==0)
+    {
+      target_block.resize(n_blocks);
+      for (unsigned int i=0;i<n_blocks;++i)
+       target_block[i] = i;
+    }
+  
+  Assert(target_block.size()==n_blocks,
+        ExcDimensionMismatch(target_block.size(),n_blocks));
+
+                                  // special case for only one
+                                  // block. treat this first
+                                  // since it does not require any
+                                  // computations
+  if (n_blocks == 1)
+    {
+      dofs_per_block[0] = dof_handler.n_dofs();
+      return;
+    } 
+                                  // otherwise determine the number
+                                  // of dofs in each block
+                                  // separately. do so in parallel
+  std::vector<std::vector<bool> >
+    dofs_in_block (n_blocks, std::vector<bool>(dof_handler.n_dofs(), false));
+  std::vector<std::vector<bool> >
+    block_select (n_blocks, std::vector<bool>(n_blocks, false));
+  Threads::ThreadGroup<> threads;
+  for (unsigned int i=0; i<n_blocks; ++i)
+    {
+      void (*fun_ptr) (const DoFHandler<dim>   &,
+                      const std::vector<bool> &,
+                      std::vector<bool>       &,
+                      bool)
+        = &DoFTools::template extract_dofs<dim>;
+      block_select[i][i] = true;
+      threads += Threads::spawn (fun_ptr)(dof_handler, block_select[i],
+                                          dofs_in_block[i], true);
+    };
+  threads.join_all ();
+
+                                  // next count what we got
+  for (unsigned int block=0;block<fe.n_blocks();++block)
+    dofs_per_block[target_block[block]]
+      += std::count(dofs_in_block[block].begin(),
+                   dofs_in_block[block].end(),
+                   true);
+}
+
+
 template <int dim>
 void
 DoFTools::count_dofs_per_component (
@@ -4100,6 +4164,12 @@ DoFTools::count_dofs_per_component<deal_II_dimension> (
   const DoFHandler<deal_II_dimension>&,
   std::vector<unsigned int>&, bool, std::vector<unsigned int>);
 
+template
+void
+DoFTools::count_dofs_per_block<deal_II_dimension> (
+  const DoFHandler<deal_II_dimension>&,
+  std::vector<unsigned int>&, std::vector<unsigned int>);
+
 template
 void
 DoFTools::count_dofs_per_component<deal_II_dimension> (

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.