]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Implement DoFTools::count_dofs_per_component for hp::DoFHandler as well.
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 14 Nov 2011 05:14:46 +0000 (05:14 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 14 Nov 2011 05:14:46 +0000 (05:14 +0000)
git-svn-id: https://svn.dealii.org/trunk@24747 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/doc/news/changes.h
deal.II/include/deal.II/dofs/dof_tools.h
deal.II/source/dofs/dof_tools.cc
deal.II/source/dofs/dof_tools.inst.in
tests/hp/count_dofs_per_block_01.cc [new file with mode: 0644]
tests/hp/count_dofs_per_block_01/cmp/generic [new file with mode: 0644]
tests/hp/count_dofs_per_block_02.cc [new file with mode: 0644]
tests/hp/count_dofs_per_block_02/cmp/generic [new file with mode: 0644]

index 5023a1fc9bb52e3a4cd9374a73610bdcd2296c59..cc9bcb4456dd7871aba19bdf307eb5c2a61cb2fb 100644 (file)
@@ -34,6 +34,11 @@ inconvenience this causes.
 <h3>General</h3>
 
 <ol>
+<li> New: The DoFTools::count_dofs_per_block function now also works
+for objects of type hp::DoFHandler.
+<br>
+(Jason Sheldon, 2011/11/13)
+
 <li> Fixed: The Intel compiler complains that it can't copy Trilinos vector
 reference objects, preventing the compiling of step-32. This is now fixed.
 <br>
index bc9870f6c1a597b6af9f2ac80a3c0ccc5a51d512..e5202f3d3627f91f5fb6591c5a3993ab4f41842a 100644 (file)
@@ -1491,7 +1491,7 @@ namespace DoFTools
                         const unsigned int level, 
                          const std::vector<bool>& selected_dofs = std::vector<bool>(), 
                          unsigned int offset = 0);
-  
+
                                    /**
                                     * Create a sparsity pattern which
                                     * in each row lists the degrees
@@ -1558,7 +1558,7 @@ namespace DoFTools
                            const bool boundary_patches = false,
                            const bool level_boundary_patches = false,
                            const bool single_cell_patches = false);
-  
+
                                    /**
                                     * Create a sparsity pattern which
                                     * in each row lists the degrees of
@@ -1587,7 +1587,7 @@ namespace DoFTools
                                     * conditions.
                                     *
                                     * Since the patches are defined
-                                    * through refinement, th 
+                                    * through refinement, th
                                     *
                                     * @arg <tt>block_list</tt>: the
                                     * SparsityPattern into which the
@@ -1654,7 +1654,7 @@ namespace DoFTools
                          const DH& dof_handler,
                          const unsigned int level,
                          const bool interior_dofs_only = false);
-  
+
                                   //@}
                                   /**
                                    * Extract a vector that represents the
@@ -1843,12 +1843,12 @@ namespace DoFTools
                                    * in the target_blocks argument
                                    * if given.
                                    */
-  template <int dim, int spacedim>
+  template <class DH>
   void
-  count_dofs_per_block (const DoFHandler<dim,spacedim>&     dof_handler,
-                       std::vector<unsigned int>& dofs_per_block,
-                       std::vector<unsigned int>  target_blocks
-                       = std::vector<unsigned int>());
+  count_dofs_per_block (const DH &dof,
+                        std::vector<unsigned int> &dofs_per_block,
+                        std::vector<unsigned int>  target_block
+                          = std::vector<unsigned int>());
 
                                   /**
                                    * @deprecated See the previous
index 560e752a660ef0672d3181cc90ece1d1cc264e67..cf6181d11c5d54dd36c92b807afc91bab0a65150 100644 (file)
@@ -4671,78 +4671,84 @@ namespace DoFTools
 
 
 
-  template <int dim, int spacedim>
+  template <class DH>
   void
-  count_dofs_per_block (const DoFHandler<dim,spacedim>& dof_handler,
-                       std::vector<unsigned int> &dofs_per_block,
-                       std::vector<unsigned int>  target_block)
+  count_dofs_per_block (const DH &dof_handler,
+                        std::vector<unsigned int> &dofs_per_block,
+                        std::vector<unsigned int>  target_block)
   {
-    const FiniteElement<dim,spacedim>& fe = dof_handler.get_fe();
-
-    std::fill (dofs_per_block.begin(), dofs_per_block.end(), 0U);
+    const dealii::hp::FECollection<DH::dimension,DH::space_dimension>
+      fe_collection (dof_handler.get_fe());
+    Assert (fe_collection.size() < 256, ExcNotImplemented());
 
-                                    // If the empty vector was given as
-                                    // default argument, set up this
-                                    // vector as identity.
-    if (target_block.size()==0)
+    for (unsigned int this_fe=0; this_fe<fe_collection.size(); ++this_fe)
       {
-       target_block.resize(fe.n_blocks());
-       for (unsigned int i=0; i<fe.n_blocks(); ++i)
-         target_block[i] = i;
-      }
-    else
-      Assert (target_block.size()==fe.n_blocks(),
-             ExcDimensionMismatch(target_block.size(),
-                                  fe.n_blocks()));
+        const FiniteElement<DH::dimension,DH::space_dimension>& fe = fe_collection[this_fe];
+        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(fe.n_blocks());
+            for (unsigned int i=0; i<fe.n_blocks(); ++i)
+             target_block[i] = i;
+         }
+        else
+         Assert (target_block.size()==fe.n_blocks(),
+                 ExcDimensionMismatch(target_block.size(),
+                                      fe.n_blocks()));
 
 
-    const unsigned int max_block
-      = *std::max_element (target_block.begin(),
-                          target_block.end());
-    const unsigned int n_target_blocks = max_block + 1;
-    const unsigned int n_blocks = fe.n_blocks();
 
-    AssertDimension (dofs_per_block.size(), n_target_blocks);
+        const unsigned int max_block
+         = *std::max_element (target_block.begin(),
+                              target_block.end());
+        const unsigned int n_target_blocks = max_block + 1;
+        const unsigned int n_blocks = fe.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.
-    std::vector<unsigned char> dofs_by_block (dof_handler.n_locally_owned_dofs());
-    internal::extract_dofs_by_component (dof_handler, std::vector<bool>(),
-                                        true, dofs_by_block);
+        AssertDimension (dofs_per_block.size(), n_target_blocks);
 
-                                    // next count what we got
-    for (unsigned int block=0; block<fe.n_blocks(); ++block)
-      dofs_per_block[target_block[block]]
-       += std::count(dofs_by_block.begin(), dofs_by_block.end(),
-                     block);
+                                        // 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.
+        std::vector<unsigned char> dofs_by_block (dof_handler.n_locally_owned_dofs());
+        internal::extract_dofs_by_component (dof_handler, std::vector<bool>(),
+                                            true, dofs_by_block);
+
+                                        // next count what we got
+        for (unsigned int block=0; block<fe.n_blocks(); ++block)
+         dofs_per_block[target_block[block]]
+           += std::count(dofs_by_block.begin(), dofs_by_block.end(),
+                         block);
 
 #ifdef DEAL_II_USE_P4EST
 #if DEAL_II_COMPILER_SUPPORTS_MPI
-                                    // if we are working on a parallel
-                                    // mesh, we now need to collect
-                                    // this information from all
-                                    // processors
-    if (const parallel::distributed::Triangulation<dim> * tria
-       = (dynamic_cast<const parallel::distributed::Triangulation<dim>*>
-          (&dof_handler.get_tria())))
-      {
-       std::vector<unsigned int> local_dof_count = dofs_per_block;
-       MPI_Allreduce ( &local_dof_count[0], &dofs_per_block[0], n_target_blocks,
-                       MPI_UNSIGNED, MPI_SUM, tria->get_communicator());
-      }
+                                        // if we are working on a parallel
+                                        // mesh, we now need to collect
+                                        // this information from all
+                                        // processors
+        if (const parallel::distributed::Triangulation<DH::dimension,DH::space_dimension> * tria
+           = (dynamic_cast<const parallel::distributed::Triangulation<DH::dimension,DH::space_dimension>*>
+              (&dof_handler.get_tria())))
+         {
+            std::vector<unsigned int> local_dof_count = dofs_per_block;
+            MPI_Allreduce ( &local_dof_count[0], &dofs_per_block[0], n_target_blocks,
+                           MPI_UNSIGNED, MPI_SUM, tria->get_communicator());
+         }
 #endif
 #endif
+      }
   }
 
 
@@ -5871,7 +5877,7 @@ namespace DoFTools
     typename DH::cell_iterator cell;
     typename DH::cell_iterator endc = dof_handler.end(level);
     std::vector<unsigned int> indices;
-    
+
     unsigned int i=0;
     for (cell=dof_handler.begin(level); cell != endc; ++i, ++cell)
       {
@@ -5893,8 +5899,8 @@ namespace DoFTools
         }
       }
   }
-  
-  
+
+
   template <class DH>
   void make_single_patch(
     SparsityPattern& block_list,
@@ -5906,15 +5912,15 @@ namespace DoFTools
       block_list.reinit(1, dof_handler.n_dofs(level), dof_handler.n_dofs(level));
       typename DH::cell_iterator cell;
       typename DH::cell_iterator endc = dof_handler.end(level);
-      
+
       std::vector<unsigned int> indices;
       std::vector<bool> exclude;
-      
+
       for (cell=dof_handler.begin(level); cell != endc;++cell)
       {
        indices.resize(cell->get_fe().dofs_per_cell);
        cell->get_mg_dof_indices(indices);
-       
+
        if (interior_only)
          {
                                             // Exclude degrees of
@@ -5923,7 +5929,7 @@ namespace DoFTools
            exclude.resize(fe.dofs_per_cell);
            std::fill(exclude.begin(), exclude.end(), false);
            const unsigned int dpf = fe.dofs_per_face;
-           
+
            for (unsigned int face=0;face<GeometryInfo<DH::dimension>::faces_per_cell;++face)
              if (cell->at_boundary(face) || cell->neighbor(face)->level() != cell->level())
                for (unsigned int i=0;i<dpf;++i)
@@ -5940,7 +5946,7 @@ namespace DoFTools
       }
     }
 
-  
+
   template <class DH>
   void make_child_patches(
     SparsityPattern& block_list,
@@ -5955,10 +5961,10 @@ namespace DoFTools
     typename DH::cell_iterator cell;
     typename DH::cell_iterator pcell = dof_handler.begin(level-1);
     typename DH::cell_iterator endc = dof_handler.end(level-1);
-    
+
     std::vector<unsigned int> indices;
     std::vector<bool> exclude;
-    
+
     for (unsigned int block = 0;pcell != endc;++pcell)
       {
        if (!pcell->has_children()) continue;
@@ -5982,7 +5988,7 @@ namespace DoFTools
                                                 // which are on faces
                                                 // of the parent
                const unsigned int dpf = fe.dofs_per_face;
-               
+
                for (unsigned int d=0;d<DH::dimension;++d)
                  {
                    const unsigned int face = GeometryInfo<DH::dimension>::vertex_to_face[child][d];
@@ -6000,7 +6006,7 @@ namespace DoFTools
                      for (unsigned int i=0;i<dpf;++i)
                        exclude[fe.face_to_cell_index(i,face)] = false;
              }
-           
+
            for (unsigned int i=0;i<n_dofs;++i)
              if (!exclude[i])
                block_list.add(block, indices[i]);
@@ -6008,9 +6014,9 @@ namespace DoFTools
        ++block;
       }
   }
-  
-  
-    
+
+
+
   template <class DH>
   void make_vertex_patches(
     SparsityPattern& block_list,
@@ -6037,7 +6043,7 @@ namespace DoFTools
                                     // Estimate for the number of
                                     // dofs at this point
     std::vector<unsigned int> vertex_dof_count(dof_handler.get_tria().n_vertices(), 0);
-    
+
                                     // Identify all vertices active
                                     // on this level and remember
                                     // some data about them
@@ -6061,7 +6067,7 @@ namespace DoFTools
                                     // From now on, only vertices
                                     // with positive dof count are
                                     // "in".
-    
+
                                     // Remove vertices at boundaries
                                     // or in corners
     for (unsigned int vg=0;vg<vertex_dof_count.size();++vg)
@@ -6081,12 +6087,12 @@ namespace DoFTools
     for (unsigned int vg=0;vg<vertex_mapping.size();++vg)
       if (vertex_dof_count[vg] != 0)
        vertex_dof_count[vertex_mapping[vg]] = vertex_dof_count[vg];
-    
+
                                     // Now that we have all the data,
                                     // we reduce it to the part we
                                     // actually want
     vertex_dof_count.resize(i);
-    
+
                                     // At this point, the list of
                                     // patches is ready. Now we enter
                                     // the dofs into the sparsity
@@ -6097,23 +6103,23 @@ namespace DoFTools
     for (unsigned int i=0;i<vertex_dof_count.size();++i)
       deallog << ' ' << vertex_dof_count[i];
     deallog << std::endl;
-    
+
     std::vector<unsigned int> indices;
     std::vector<bool> exclude;
-    
+
     for (cell=dof_handler.begin(level); cell != endc; ++cell)
       {
        const FiniteElement<DH::dimension>& fe = cell->get_fe();
        indices.resize(fe.dofs_per_cell);
        cell->get_mg_dof_indices(indices);
-       
+
        for (unsigned int v=0;v<GeometryInfo<DH::dimension>::vertices_per_cell;++v)
          {
            const unsigned int vg = cell->vertex_index(v);
            const unsigned int block = vertex_mapping[vg];
            if (block == numbers::invalid_unsigned_int)
              continue;
-           
+
            if (interior_only)
              {
                                                 // Exclude degrees of
@@ -6122,7 +6128,7 @@ namespace DoFTools
                exclude.resize(fe.dofs_per_cell);
                std::fill(exclude.begin(), exclude.end(), false);
                const unsigned int dpf = fe.dofs_per_face;
-               
+
                for (unsigned int d=0;d<DH::dimension;++d)
                  {
                    const unsigned int a_face = GeometryInfo<DH::dimension>::vertex_to_face[v][d];
index 866832b649cea8ee89e904468fe6001c77177a1b..996ba002a9c558c0574327558755801d926da015 100644 (file)
@@ -543,10 +543,21 @@ DoFTools::extract_level_dofs<deal_II_dimension, deal_II_dimension+1>
 
 template
 void
-DoFTools::count_dofs_per_block<deal_II_dimension> (
+DoFTools::count_dofs_per_block<DoFHandler<deal_II_dimension> > (
   const DoFHandler<deal_II_dimension>&,
   std::vector<unsigned int>&, std::vector<unsigned int>);
 
+template
+void
+DoFTools::count_dofs_per_block<hp::DoFHandler<deal_II_dimension> > (
+  const hp::DoFHandler<deal_II_dimension>&,
+  std::vector<unsigned int>&, std::vector<unsigned int>);
+
+template
+void
+DoFTools::count_dofs_per_block<MGDoFHandler<deal_II_dimension> > (
+  const MGDoFHandler<deal_II_dimension>&,
+  std::vector<unsigned int>&, std::vector<unsigned int>);
 
 template
 void
@@ -619,10 +630,22 @@ DoFTools::map_dofs_to_support_points<deal_II_dimension,deal_II_dimension+1>
 
 template
 void
-DoFTools::count_dofs_per_block<deal_II_dimension,deal_II_dimension+1> (
+DoFTools::count_dofs_per_block<DoFHandler<deal_II_dimension,deal_II_dimension+1> > (
   const DoFHandler<deal_II_dimension,deal_II_dimension+1>&,
   std::vector<unsigned int>&, std::vector<unsigned int>);
 
+template
+void
+DoFTools::count_dofs_per_block<hp::DoFHandler<deal_II_dimension,deal_II_dimension+1> > (
+  const hp::DoFHandler<deal_II_dimension,deal_II_dimension+1>&,
+  std::vector<unsigned int>&, std::vector<unsigned int>);
+
+template
+void
+DoFTools::count_dofs_per_block<MGDoFHandler<deal_II_dimension,deal_II_dimension+1> > (
+  const MGDoFHandler<deal_II_dimension,deal_II_dimension+1>&,
+  std::vector<unsigned int>&, std::vector<unsigned int>);
+
 #endif
 
 template
diff --git a/tests/hp/count_dofs_per_block_01.cc b/tests/hp/count_dofs_per_block_01.cc
new file mode 100644 (file)
index 0000000..a921f8e
--- /dev/null
@@ -0,0 +1,154 @@
+//----------------------------  count_dofs_per_block_01.cc  ---------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2009, 2011 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//----------------------------  count_dofs_per_block_01.cc  ---------------------------
+
+// like bits/count_dofs_per_block_01 but for hp::DoFHandler
+
+
+#include "../tests.h"
+#include <deal.II/base/logstream.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/tria_iterator.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/dofs/dof_accessor.h>
+#include <deal.II/hp/dof_handler.h>
+#include <deal.II/hp/fe_collection.h>
+#include <deal.II/dofs/dof_tools.h>
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/fe/fe_dgq.h>
+#include <deal.II/fe/fe_dgp.h>
+#include <deal.II/fe/fe_nedelec.h>
+#include <deal.II/fe/fe_raviart_thomas.h>
+#include <deal.II/fe/fe_system.h>
+
+#include <fstream>
+#include <iomanip>
+#include <iomanip>
+#include <string>
+
+
+std::string output_file_name = "count_dofs_per_block_01/output";
+
+
+
+void print (const std::vector<unsigned int> &v)
+{
+  deallog << v.size();
+  for (unsigned int i=0; i<v.size(); ++i)
+    deallog << ' ' << v[i];
+  deallog << std::endl;
+}
+
+
+
+template <int dim>
+void
+check ()
+{
+                                   // create tria and dofhandler
+                                   // objects. set different boundary
+                                   // and sub-domain ids
+  Triangulation<dim> tria;
+  GridGenerator::hyper_cube(tria, 0., 1.);
+  tria.refine_global (1);
+  for (int i=0; i<2; ++i)
+    {
+      tria.begin_active()->set_refine_flag();
+      tria.execute_coarsening_and_refinement ();
+    }
+
+                                  // produce two copies of the Taylor-Hood
+                                  // element
+  hp::FECollection<dim> fe;
+  fe.push_back (FESystem<dim> (FE_Q<dim>(1), dim, FE_DGQ<dim>(0), 1));
+  fe.push_back (FESystem<dim> (FE_Q<dim>(1), dim, FE_DGQ<dim>(0), 1));
+
+  hp::DoFHandler<dim> dof_handler (tria);
+
+  typename hp::DoFHandler<dim>::active_cell_iterator
+    cell = dof_handler.begin_active(),
+    endc = dof_handler.end();
+  for (; cell!=endc; ++cell)
+    cell->set_active_fe_index (std::rand() % fe.size());
+
+  dof_handler.distribute_dofs (fe);
+
+                                  // no grouping
+  {
+    std::vector<unsigned int> dpc(dim+1);
+    DoFTools::count_dofs_per_component (dof_handler, dpc);
+    print (dpc);
+  }
+
+  {
+    std::vector<unsigned int> dpc(dim+1);
+    DoFTools::count_dofs_per_block (dof_handler, dpc);
+    print (dpc);
+  }
+
+
+                                  // grouping into less groups than
+                                  // components
+  {
+    std::vector<unsigned int> group(dim+1, 0);
+    group[dim] = 1;
+    std::vector<unsigned int> dpc(2);
+    DoFTools::count_dofs_per_component (dof_handler, dpc, false, group);
+    Assert (dpc.size() == 2, ExcInternalError());
+    print (dpc);
+  }
+
+  {
+    std::vector<unsigned int> group(dim+1, 0);
+    group[dim] = 1;
+    std::vector<unsigned int> dpc(2);
+    DoFTools::count_dofs_per_block (dof_handler, dpc, group);
+    Assert (dpc.size() == 2, ExcInternalError());
+    print (dpc);
+  }
+
+                                  // grouping into more groups than
+                                  // components
+  {
+    std::vector<unsigned int> group(dim+1, 2*dim);
+    group[dim] = 0;
+    std::vector<unsigned int> dpc(2*dim+1);
+    DoFTools::count_dofs_per_component (dof_handler, dpc, false, group);
+    Assert (dpc.size() == 2*dim+1, ExcInternalError());
+    print (dpc);
+  }
+
+  {
+    std::vector<unsigned int> group(dim+1, 2*dim);
+    group[dim] = 0;
+    std::vector<unsigned int> dpc(2*dim+1);
+    DoFTools::count_dofs_per_block (dof_handler, dpc, group);
+    Assert (dpc.size() == 2*dim+1, ExcInternalError());
+    print (dpc);
+  }
+}
+
+
+
+int main()
+{
+  std::ofstream logfile(output_file_name.c_str());
+  logfile << std::setprecision (2);
+  deallog << std::setprecision (2);
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-10);
+
+  check<1> ();
+  check<2> ();
+  check<3> ();
+}
diff --git a/tests/hp/count_dofs_per_block_01/cmp/generic b/tests/hp/count_dofs_per_block_01/cmp/generic
new file mode 100644 (file)
index 0000000..1728b32
--- /dev/null
@@ -0,0 +1,19 @@
+
+DEAL::2 5 4
+DEAL::2 5 4
+DEAL::2 5 4
+DEAL::2 5 4
+DEAL::3 4 0 5
+DEAL::3 4 0 5
+DEAL::3 18 18 10
+DEAL::3 18 18 10
+DEAL::2 36 10
+DEAL::2 36 10
+DEAL::5 10 0 0 0 36
+DEAL::5 10 0 0 0 36
+DEAL::4 60 60 60 22
+DEAL::4 60 60 60 22
+DEAL::2 180 22
+DEAL::2 180 22
+DEAL::7 22 0 0 0 0 0 180
+DEAL::7 22 0 0 0 0 0 180
diff --git a/tests/hp/count_dofs_per_block_02.cc b/tests/hp/count_dofs_per_block_02.cc
new file mode 100644 (file)
index 0000000..9a2bb0d
--- /dev/null
@@ -0,0 +1,154 @@
+//----------------------------  count_dofs_per_block_02.cc  ---------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2009, 2011 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//----------------------------  count_dofs_per_block_02.cc  ---------------------------
+
+// like _01 but with different elements
+
+
+#include "../tests.h"
+#include <deal.II/base/logstream.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/tria_iterator.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/dofs/dof_accessor.h>
+#include <deal.II/hp/dof_handler.h>
+#include <deal.II/hp/fe_collection.h>
+#include <deal.II/dofs/dof_tools.h>
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/fe/fe_dgq.h>
+#include <deal.II/fe/fe_dgp.h>
+#include <deal.II/fe/fe_nedelec.h>
+#include <deal.II/fe/fe_raviart_thomas.h>
+#include <deal.II/fe/fe_system.h>
+
+#include <fstream>
+#include <iomanip>
+#include <iomanip>
+#include <string>
+
+
+std::string output_file_name = "count_dofs_per_block_02/output";
+
+
+
+void print (const std::vector<unsigned int> &v)
+{
+  deallog << v.size();
+  for (unsigned int i=0; i<v.size(); ++i)
+    deallog << ' ' << v[i];
+  deallog << std::endl;
+}
+
+
+
+template <int dim>
+void
+check ()
+{
+                                   // create tria and dofhandler
+                                   // objects. set different boundary
+                                   // and sub-domain ids
+  Triangulation<dim> tria;
+  GridGenerator::hyper_cube(tria, 0., 1.);
+  tria.refine_global (1);
+  for (int i=0; i<2; ++i)
+    {
+      tria.begin_active()->set_refine_flag();
+      tria.execute_coarsening_and_refinement ();
+    }
+
+                                  // produce two copies of the Taylor-Hood
+                                  // element
+  hp::FECollection<dim> fe;
+  fe.push_back (FESystem<dim> (FE_Q<dim>(1), dim, FE_DGQ<dim>(0), 1));
+  fe.push_back (FESystem<dim> (FE_Q<dim>(2), dim, FE_DGQ<dim>(1), 1));
+
+  hp::DoFHandler<dim> dof_handler (tria);
+
+  typename hp::DoFHandler<dim>::active_cell_iterator
+    cell = dof_handler.begin_active(),
+    endc = dof_handler.end();
+  for (; cell!=endc; ++cell)
+    cell->set_active_fe_index (std::rand() % fe.size());
+
+  dof_handler.distribute_dofs (fe);
+
+                                  // no grouping
+  {
+    std::vector<unsigned int> dpc(dim+1);
+    DoFTools::count_dofs_per_component (dof_handler, dpc);
+    print (dpc);
+  }
+
+  {
+    std::vector<unsigned int> dpc(dim+1);
+    DoFTools::count_dofs_per_block (dof_handler, dpc);
+    print (dpc);
+  }
+
+
+                                  // grouping into less groups than
+                                  // components
+  {
+    std::vector<unsigned int> group(dim+1, 0);
+    group[dim] = 1;
+    std::vector<unsigned int> dpc(2);
+    DoFTools::count_dofs_per_component (dof_handler, dpc, false, group);
+    Assert (dpc.size() == 2, ExcInternalError());
+    print (dpc);
+  }
+
+  {
+    std::vector<unsigned int> group(dim+1, 0);
+    group[dim] = 1;
+    std::vector<unsigned int> dpc(2);
+    DoFTools::count_dofs_per_block (dof_handler, dpc, group);
+    Assert (dpc.size() == 2, ExcInternalError());
+    print (dpc);
+  }
+
+                                  // grouping into more groups than
+                                  // components
+  {
+    std::vector<unsigned int> group(dim+1, 2*dim);
+    group[dim] = 0;
+    std::vector<unsigned int> dpc(2*dim+1);
+    DoFTools::count_dofs_per_component (dof_handler, dpc, false, group);
+    Assert (dpc.size() == 2*dim+1, ExcInternalError());
+    print (dpc);
+  }
+
+  {
+    std::vector<unsigned int> group(dim+1, 2*dim);
+    group[dim] = 0;
+    std::vector<unsigned int> dpc(2*dim+1);
+    DoFTools::count_dofs_per_block (dof_handler, dpc, group);
+    Assert (dpc.size() == 2*dim+1, ExcInternalError());
+    print (dpc);
+  }
+}
+
+
+
+int main()
+{
+  std::ofstream logfile(output_file_name.c_str());
+  logfile << std::setprecision (2);
+  deallog << std::setprecision (2);
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-10);
+
+  check<1> ();
+  check<2> ();
+  check<3> ();
+}
diff --git a/tests/hp/count_dofs_per_block_02/cmp/generic b/tests/hp/count_dofs_per_block_02/cmp/generic
new file mode 100644 (file)
index 0000000..d5d6551
--- /dev/null
@@ -0,0 +1,19 @@
+
+DEAL::2 8 7
+DEAL::2 8 7
+DEAL::2 8 7
+DEAL::2 8 7
+DEAL::3 7 0 8
+DEAL::3 7 0 8
+DEAL::3 45 45 28
+DEAL::3 45 45 28
+DEAL::2 90 28
+DEAL::2 90 28
+DEAL::5 28 0 0 0 90
+DEAL::5 28 0 0 0 90
+DEAL::4 205 205 205 92
+DEAL::4 205 205 205 92
+DEAL::2 615 92
+DEAL::2 615 92
+DEAL::7 92 0 0 0 0 0 615
+DEAL::7 92 0 0 0 0 0 615

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.