]> https://gitweb.dealii.org/ - dealii.git/commitdiff
change BlockList to SparsityPattern
authorGuido Kanschat <dr.guido.kanschat@gmail.com>
Fri, 28 Oct 2011 21:05:42 +0000 (21:05 +0000)
committerGuido Kanschat <dr.guido.kanschat@gmail.com>
Fri, 28 Oct 2011 21:05:42 +0000 (21:05 +0000)
git-svn-id: https://svn.dealii.org/trunk@24684 0785d39b-7218-0410-832d-ea1e28bc413d

tests/deal.II/block_list_01.cc
tests/deal.II/block_list_02.cc

index ebd9e8a88f9af3917eba081e9c2cf4acb3a80669..36269e8125f42a881950d4c7d201b8cb08e4dcdf 100644 (file)
@@ -17,6 +17,7 @@
 #include <deal.II/multigrid/mg_dof_handler.h>
 #include <deal.II/grid/tria_iterator.h>
 #include <deal.II/multigrid/mg_dof_accessor.h>
+#include <dofs/dof_tools.h>
 #include <deal.II/fe/fe_q.h>
 #include <deal.II/fe/fe_dgq.h>
 #include <deal.II/fe/fe_raviart_thomas.h>
@@ -38,17 +39,18 @@ test_block_list(const Triangulation<dim>& tr, const FiniteElement<dim>& fe)
   MGDoFHandler<dim> dof;
   dof.initialize(tr, fe);
   
-  BlockList bl;
-
+  
   const unsigned int level = tr.n_levels()-1;
   
-  bl.initialize_mg(tr.n_cells(level), dof.begin(level), dof.end(level));
-
-  for (unsigned int i=0;i<bl.size();++i)
+  SparsityPattern bl(tr.n_cells(level), dof.n_dofs(level), fe.dofs_per_cell);
+  DoFTools::make_cell_patches(bl, dof, level);
+  bl.compress();
+  
+  for (unsigned int i=0;i<bl.n_rows();++i)
     {
       deallog << "Block " << std::setw(3) << i;
       std::vector<unsigned int> entries;
-      for (BlockList::const_iterator b = bl.begin(i);b != bl.end(i);++b)
+      for (SparsityPattern::row_iterator b = bl.row_begin(i);b != bl.row_end(i);++b)
        entries.push_back(*b);
 
       std::sort(entries.begin(), entries.end());
index 736cfcb4ab7cc68eb4b2693170ae0f5aef617387..c925707d94fa97102e637e5bd4e6aaa8a9085ea4 100644 (file)
@@ -17,6 +17,7 @@
 #include <deal.II/multigrid/mg_dof_handler.h>
 #include <deal.II/grid/tria_iterator.h>
 #include <deal.II/multigrid/mg_dof_accessor.h>
+#include <dofs/dof_tools.h>
 #include <deal.II/fe/fe_q.h>
 #include <deal.II/fe/fe_dgq.h>
 #include <deal.II/fe/fe_raviart_thomas.h>
 using namespace dealii;
 
 
-template <int dim>
 void
-test_block_list(const Triangulation<dim>& tr, const FiniteElement<dim>& fe)
+print_patches (const SparsityPattern& bl)
 {
-  deallog << fe.get_name() << std::endl;
-  
-  MGDoFHandler<dim> dof;
-  dof.initialize(tr, fe);
-  
-  BlockList bl;
-
-  const unsigned int level = tr.n_levels()-1;
-  const unsigned int n = bl.count_vertex_patches<2>(dof.begin(level), dof.end(level), true);
-  bl.initialize_vertex_patches_mg<dim>(n, dof.begin(level), dof.end(level));
-
-  for (unsigned int i=0;i<bl.size();++i)
+  for (unsigned int i=0;i<bl.n_rows();++i)
     {
       deallog << "Block " << std::setw(3) << i;
       std::vector<unsigned int> entries;
-      for (BlockList::const_iterator b = bl.begin(i);b != bl.end(i);++b)
+      for (SparsityPattern::row_iterator b = bl.row_begin(i);b != bl.row_end(i);++b)
        entries.push_back(*b);
-
+      
       std::sort(entries.begin(), entries.end());
-
+      
       for (unsigned int i=0;i<entries.size();++i)
        deallog << ' ' << std::setw(4) << entries[i];
       deallog << std::endl;
@@ -60,6 +49,65 @@ test_block_list(const Triangulation<dim>& tr, const FiniteElement<dim>& fe)
 }
 
 
+template <int dim>
+void
+test_block_list(const Triangulation<dim>& tr, const FiniteElement<dim>& fe)
+{
+  deallog << fe.get_name() << std::endl;
+  
+  MGDoFHandler<dim> dof;
+  dof.initialize(tr, fe); 
+  
+  const unsigned int level = tr.n_levels()-1;
+
+  {
+    deallog.push("ttt");
+    SparsityPattern bl;
+    DoFTools::make_vertex_patches(bl, dof, level, true, true, true);
+    bl.compress();
+    print_patches(bl);
+    deallog.pop();
+    deallog << std::endl;
+  }
+  {
+    deallog.push("ttf");
+    SparsityPattern bl;
+    DoFTools::make_vertex_patches(bl, dof, level, true, true, false);
+    bl.compress();
+    print_patches(bl);
+    deallog.pop(); 
+    deallog << std::endl;
+  }
+  {
+    deallog.push("tff");
+    SparsityPattern bl;
+    DoFTools::make_vertex_patches(bl, dof, level, true, false, false);
+    bl.compress();
+    print_patches(bl);
+    deallog.pop();
+    deallog << std::endl;
+  }
+  {
+    deallog.push("ftt");
+    SparsityPattern bl;
+    DoFTools::make_vertex_patches(bl, dof, level, false, true, true);
+    bl.compress();
+    print_patches(bl);
+    deallog.pop();
+    deallog << std::endl;
+  }
+  {
+    deallog.push("fff");
+    SparsityPattern bl;
+    DoFTools::make_vertex_patches(bl, dof, level, false, false, false);
+    bl.compress();
+    print_patches(bl);
+    deallog.pop();
+    deallog << std::endl;
+  }
+}
+
+
 int main()
 {
   std::string logname = JobIdentifier::base_name(__FILE__) + std::string("/output");

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.