]> https://gitweb.dealii.org/ - dealii.git/commitdiff
check counting by component
authorGuido Kanschat <dr.guido.kanschat@gmail.com>
Fri, 5 May 2006 07:44:36 +0000 (07:44 +0000)
committerGuido Kanschat <dr.guido.kanschat@gmail.com>
Fri, 5 May 2006 07:44:36 +0000 (07:44 +0000)
git-svn-id: https://svn.dealii.org/trunk@13055 0785d39b-7218-0410-832d-ea1e28bc413d

tests/multigrid/Makefile
tests/multigrid/count_01.cc [new file with mode: 0644]
tests/multigrid/count_01/cmp/generic [new file with mode: 0644]

index 0ac46ac8d5d14a554660629f3114d43267f756df..06559e9d0aa8d83cca0cd05ba3bc2f333d00a1f8 100644 (file)
@@ -1,6 +1,6 @@
 ############################################################
 # $Id$
-# Copyright (C) 2000, 2001, 2002, 2003, 2005 by the deal.II authors
+# Copyright (C) 2000, 2001, 2002, 2003, 2005, 2006 by the deal.II authors
 ############################################################
 
 ############################################################
@@ -22,7 +22,7 @@ default: run-tests
 ############################################################
 
 tests_x = cycles transfer transfer_select transfer_block \
-       smoother_block renumbering_*
+       smoother_block count_* renumbering_*
 
 # from above list of regular expressions, generate the real set of
 # tests
diff --git a/tests/multigrid/count_01.cc b/tests/multigrid/count_01.cc
new file mode 100644 (file)
index 0000000..d2f8ebf
--- /dev/null
@@ -0,0 +1,117 @@
+//----------------------------------------------------------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2006 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.
+//
+//----------------------------------------------------------------------------
+
+// check MGTools::count_dofs_per_component
+
+#include <base/logstream.h>
+#include <lac/vector.h>
+#include <lac/block_vector.h>
+#include <grid/tria.h>
+#include <grid/tria_iterator.h>
+#include <grid/tria_accessor.h>
+#include <grid/grid_generator.h>
+#include <fe/fe_dgq.h>
+#include <fe/fe_q.h>
+#include <fe/fe_raviart_thomas.h>
+#include <fe/fe_system.h>
+#include <multigrid/mg_dof_handler.h>
+#include <multigrid/mg_tools.h>
+
+#include <fstream>
+#include <iostream>
+#include <iomanip>
+#include <algorithm>
+
+using namespace std;
+
+void log_vector (const std::vector<std::vector<unsigned int> >& count)
+{
+  for (unsigned int l=0;l<count.size();++l)
+    {
+      deallog << "Level " << l;
+      for (unsigned int c=0;c<count[l].size();++c)
+       deallog << '\t' << count[l][c];
+      deallog << std::endl;
+    }
+}
+
+template <int dim>
+void check_fe(FiniteElement<dim>& fe)
+{
+  deallog << fe.get_name() << std::endl;
+  
+  Triangulation<dim> tr;
+  GridGenerator::hyper_cube(tr);
+  tr.refine_global(1);
+  tr.begin_active()->set_refine_flag();
+  tr.execute_coarsening_and_refinement ();
+  tr.refine_global(1);
+  
+  MGDoFHandler<dim> mgdof(tr);
+  mgdof.distribute_dofs(fe);
+  
+  std::vector<std::vector<unsigned int> > count(tr.n_levels());
+  MGTools::count_dofs_per_component(mgdof, count, false);
+  log_vector(count);
+  MGTools::count_dofs_per_component(mgdof, count, true);
+  log_vector(count);
+
+  std::vector<unsigned int> target(fe.n_components());
+  for (unsigned int i=0;i<target.size();++i)
+    target[i] = i/3;
+  deallog << std::endl << "Target";
+  for (unsigned int i=0;i<target.size();++i)
+    deallog << '\t' << target[i];
+  deallog << std::endl;
+  
+  MGTools::count_dofs_per_component(mgdof, count, false, target);
+  log_vector(count);
+  MGTools::count_dofs_per_component(mgdof, count, true, target);
+  log_vector(count);
+  
+}
+
+
+template <int dim>
+void check()
+{
+  FE_Q<dim> q1(1);
+  FE_Q<dim> q2(2);
+  FE_DGQ<dim> dq1(1);
+  
+  FESystem<dim> s1(q1, 2, q2,1);
+
+  check_fe(s1);
+  if (dim>1)
+    {
+      FE_RaviartThomas<dim> rt(1);
+      FESystem<dim> s2(rt, 2, dq1,1);
+      FESystem<dim> s3(rt,1, s1, 2);
+      
+      check_fe(s2);
+      check_fe(s3);
+    }
+}
+
+int main()
+{
+  std::ofstream logfile("count_01/output");
+  logfile.precision(3);
+  deallog.attach(logfile);
+  deallog.depth_console(10);
+  deallog.threshold_double(1.e-10);
+
+  check<1> ();
+  check<2> ();
+  check<3> ();
+}
diff --git a/tests/multigrid/count_01/cmp/generic b/tests/multigrid/count_01/cmp/generic
new file mode 100644 (file)
index 0000000..4326ee4
--- /dev/null
@@ -0,0 +1,134 @@
+
+DEAL::FESystem<1>[FE_Q<1>(1)^2-FE_Q<1>(2)]
+DEAL::Level 0  2       2       3
+DEAL::Level 1  3       3       5
+DEAL::Level 2  5       5       9
+DEAL::Level 3  5       5       9
+DEAL::Level 0  2       2       3
+DEAL::Level 1  3       3       5
+DEAL::Level 2  5       5       9
+DEAL::Level 3  5       5       9
+DEAL::
+DEAL::Target   0       0       0
+DEAL::Level 0  7       0       0
+DEAL::Level 1  11      0       0
+DEAL::Level 2  19      0       0
+DEAL::Level 3  19      0       0
+DEAL::Level 0  7       0       0
+DEAL::Level 1  11      0       0
+DEAL::Level 2  19      0       0
+DEAL::Level 3  19      0       0
+DEAL::FESystem<2>[FE_Q<2>(1)^2-FE_Q<2>(2)]
+DEAL::Level 0  4       4       9
+DEAL::Level 1  9       9       25
+DEAL::Level 2  25      25      81
+DEAL::Level 3  25      25      81
+DEAL::Level 0  4       4       9
+DEAL::Level 1  9       9       25
+DEAL::Level 2  25      25      81
+DEAL::Level 3  25      25      81
+DEAL::
+DEAL::Target   0       0       0
+DEAL::Level 0  17      0       0
+DEAL::Level 1  43      0       0
+DEAL::Level 2  131     0       0
+DEAL::Level 3  131     0       0
+DEAL::Level 0  17      0       0
+DEAL::Level 1  43      0       0
+DEAL::Level 2  131     0       0
+DEAL::Level 3  131     0       0
+DEAL::FESystem<2>[FE_RaviartThomas<2>(1)^2-FE_DGQ<2>(1)]
+DEAL::Level 0  12      12      12      12      4
+DEAL::Level 1  40      40      40      40      16
+DEAL::Level 2  144     144     144     144     64
+DEAL::Level 3  144     144     144     144     64
+DEAL::Level 0  12      0       12      0       4
+DEAL::Level 1  40      0       40      0       16
+DEAL::Level 2  144     0       144     0       64
+DEAL::Level 3  144     0       144     0       64
+DEAL::
+DEAL::Target   0       0       0       1       1
+DEAL::Level 0  36      16      0       0       0
+DEAL::Level 1  120     56      0       0       0
+DEAL::Level 2  432     208     0       0       0
+DEAL::Level 3  432     208     0       0       0
+DEAL::Level 0  24      4       0       0       0
+DEAL::Level 1  80      16      0       0       0
+DEAL::Level 2  288     64      0       0       0
+DEAL::Level 3  288     64      0       0       0
+DEAL::FESystem<2>[FE_RaviartThomas<2>(1)-FESystem<2>[FE_Q<2>(1)^2-FE_Q<2>(2)]^2]
+DEAL::Level 0  12      12      4       4       9       4       4       9
+DEAL::Level 1  40      40      9       9       25      9       9       25
+DEAL::Level 2  144     144     25      25      81      25      25      81
+DEAL::Level 3  144     144     25      25      81      25      25      81
+DEAL::Level 0  12      0       4       4       9       4       4       9
+DEAL::Level 1  40      0       9       9       25      9       9       25
+DEAL::Level 2  144     0       25      25      81      25      25      81
+DEAL::Level 3  144     0       25      25      81      25      25      81
+DEAL::
+DEAL::Target   0       0       0       1       1       1       2       2
+DEAL::Level 0  28      17      13      0       0       0       0       0
+DEAL::Level 1  89      43      34      0       0       0       0       0
+DEAL::Level 2  313     131     106     0       0       0       0       0
+DEAL::Level 3  313     131     106     0       0       0       0       0
+DEAL::Level 0  16      17      13      0       0       0       0       0
+DEAL::Level 1  49      43      34      0       0       0       0       0
+DEAL::Level 2  169     131     106     0       0       0       0       0
+DEAL::Level 3  169     131     106     0       0       0       0       0
+DEAL::FESystem<3>[FE_Q<3>(1)^2-FE_Q<3>(2)]
+DEAL::Level 0  8       8       27
+DEAL::Level 1  27      27      125
+DEAL::Level 2  125     125     729
+DEAL::Level 3  125     125     729
+DEAL::Level 0  8       8       27
+DEAL::Level 1  27      27      125
+DEAL::Level 2  125     125     729
+DEAL::Level 3  125     125     729
+DEAL::
+DEAL::Target   0       0       0
+DEAL::Level 0  43      0       0
+DEAL::Level 1  179     0       0
+DEAL::Level 2  979     0       0
+DEAL::Level 3  979     0       0
+DEAL::Level 0  43      0       0
+DEAL::Level 1  179     0       0
+DEAL::Level 2  979     0       0
+DEAL::Level 3  979     0       0
+DEAL::FESystem<3>[FE_RaviartThomas<3>(1)^2-FE_DGQ<3>(1)]
+DEAL::Level 0  36      36      36      36      36      36      8
+DEAL::Level 1  240     240     240     240     240     240     64
+DEAL::Level 2  1728    1728    1728    1728    1728    1728    512
+DEAL::Level 3  1728    1728    1728    1728    1728    1728    512
+DEAL::Level 0  36      0       0       36      0       0       8
+DEAL::Level 1  240     0       0       240     0       0       64
+DEAL::Level 2  1728    0       0       1728    0       0       512
+DEAL::Level 3  1728    0       0       1728    0       0       512
+DEAL::
+DEAL::Target   0       0       0       1       1       1       2
+DEAL::Level 0  108     108     8       0       0       0       0
+DEAL::Level 1  720     720     64      0       0       0       0
+DEAL::Level 2  5184    5184    512     0       0       0       0
+DEAL::Level 3  5184    5184    512     0       0       0       0
+DEAL::Level 0  36      36      8       0       0       0       0
+DEAL::Level 1  240     240     64      0       0       0       0
+DEAL::Level 2  1728    1728    512     0       0       0       0
+DEAL::Level 3  1728    1728    512     0       0       0       0
+DEAL::FESystem<3>[FE_RaviartThomas<3>(1)-FESystem<3>[FE_Q<3>(1)^2-FE_Q<3>(2)]^2]
+DEAL::Level 0  36      36      36      8       8       27      8       8       27
+DEAL::Level 1  240     240     240     27      27      125     27      27      125
+DEAL::Level 2  1728    1728    1728    125     125     729     125     125     729
+DEAL::Level 3  1728    1728    1728    125     125     729     125     125     729
+DEAL::Level 0  36      0       0       8       8       27      8       8       27
+DEAL::Level 1  240     0       0       27      27      125     27      27      125
+DEAL::Level 2  1728    0       0       125     125     729     125     125     729
+DEAL::Level 3  1728    0       0       125     125     729     125     125     729
+DEAL::
+DEAL::Target   0       0       0       1       1       1       2       2       2
+DEAL::Level 0  108     43      43      0       0       0       0       0       0
+DEAL::Level 1  720     179     179     0       0       0       0       0       0
+DEAL::Level 2  5184    979     979     0       0       0       0       0       0
+DEAL::Level 3  5184    979     979     0       0       0       0       0       0
+DEAL::Level 0  36      43      43      0       0       0       0       0       0
+DEAL::Level 1  240     179     179     0       0       0       0       0       0
+DEAL::Level 2  1728    979     979     0       0       0       0       0       0
+DEAL::Level 3  1728    979     979     0       0       0       0       0       0

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.