]> https://gitweb.dealii.org/ - dealii.git/commitdiff
test dg numbering
authorGuido Kanschat <dr.guido.kanschat@gmail.com>
Thu, 3 Jul 2003 11:49:00 +0000 (11:49 +0000)
committerGuido Kanschat <dr.guido.kanschat@gmail.com>
Thu, 3 Jul 2003 11:49:00 +0000 (11:49 +0000)
git-svn-id: https://svn.dealii.org/trunk@7835 0785d39b-7218-0410-832d-ea1e28bc413d

tests/deal.II/dof_renumbering.cc

index 72c55b5dcfa3de6aceed92c386a4e63b17d4e081..4968214a0ab2cc81190fc455ace75be08e934528 100644 (file)
@@ -29,6 +29,7 @@
 #include <multigrid/mg_dof_handler.h>
 #include <fe/fe_q.h>
 #include <fe/fe_dgq.h>
+#include <fe/fe_dgp.h>
 #include <fe/fe_system.h>
 
 #include <fstream>
@@ -70,44 +71,40 @@ print_dofs (const MGDoFHandler<dim> &dof, unsigned int level)
 }
 
 
-
 template <int dim>
 void
-check ()
+check_renumbering(MGDoFHandler<dim>& mgdof, bool discontinuous)
 {
-  Functions::CosineFunction<dim> cosine;
-  
-  Triangulation<dim> tr;  
-  if (dim==2)
-    GridGenerator::hyper_ball(tr, Point<dim>(), 1);
-  else
-    GridGenerator::hyper_cube(tr, -1,1);
-  tr.refine_global (1);
-  tr.begin_active()->set_refine_flag ();
-  tr.execute_coarsening_and_refinement ();
-  if (dim==1)
-    tr.refine_global(2);
-  
-  FESystem<dim> element (FE_Q<dim>(2), 2, FE_DGQ<dim>(1), 1);
-  MGDoFHandler<dim> mgdof(tr);
+  const FiniteElement<dim>& element = mgdof.get_fe();
   DoFHandler<dim>& dof = mgdof;
-  dof.distribute_dofs(element);
-
+  
                                   // Prepare a reordering of
                                   // components for later use
   std::vector<unsigned int> order(element.n_components());
   for (unsigned int i=0; i<order.size(); ++i) order[i] = order.size()-i-1;
 
+  Point<dim> direction;
+  for (unsigned int i=0;i<dim;++i)
+    direction(i) = pow(10.,i);
+  
                                   // Check global ordering
   print_dofs (dof);
   
-  DoFRenumbering::Cuthill_McKee (dof);
-  print_dofs (dof);
-  DoFRenumbering::Cuthill_McKee (dof, true);
-  print_dofs (dof);
-  DoFRenumbering::Cuthill_McKee (dof, true, true);
-  print_dofs (dof);
-
+  if (discontinuous)
+    {
+      DoFRenumbering::downstream_dg(dof, direction);
+    }
+  else
+    {
+      DoFRenumbering::Cuthill_McKee (dof);
+      print_dofs (dof);
+      DoFRenumbering::Cuthill_McKee (dof, true);
+      print_dofs (dof);
+      DoFRenumbering::Cuthill_McKee (dof, true, true);
+      print_dofs (dof);
+    }
+  
+  
   DoFRenumbering::component_wise (dof, order);
   print_dofs (dof);
 
@@ -117,7 +114,7 @@ check ()
   print_dofs (dof);
 
                                   // Check level ordering
-  for (unsigned int level=0;level<tr.n_levels();++level)
+  for (unsigned int level=0;level<dof.get_tria().n_levels();++level)
     {
       print_dofs (mgdof, level);
 
@@ -127,9 +124,44 @@ check ()
 //        DoFRenumbering::Cuthill_McKee (mgdof, level, true);
 //        print_dofs (mgdof, level);
 
+      if (discontinuous)
+       {
+         DoFRenumbering::downstream_dg(mgdof, level, direction);
+       }
+  
       DoFRenumbering::component_wise (mgdof, order);
       print_dofs (mgdof, level);
     }
+
+}
+
+
+template <int dim>
+void
+check ()
+{
+  Triangulation<dim> tr;  
+  if (dim==2)
+    GridGenerator::hyper_ball(tr, Point<dim>(), 1);
+  else
+    GridGenerator::hyper_cube(tr, -1,1);
+  tr.refine_global (1);
+  tr.begin_active()->set_refine_flag ();
+  tr.execute_coarsening_and_refinement ();
+  if (dim==1)
+    tr.refine_global(2);
+
+  MGDoFHandler<dim> mgdof(tr);
+  
+  FESystem<dim> e1 (FE_Q<dim>(2), 2, FE_DGQ<dim>(1), 1);
+  mgdof.distribute_dofs(e1);
+  check_renumbering(mgdof, false);
+  mgdof.clear();
+  
+  FESystem<dim> e2 (FE_DGP<dim>(2), 2, FE_DGQ<dim>(1), 1);
+  mgdof.distribute_dofs(e2);
+  check_renumbering(mgdof, true);
+  mgdof.clear();
 }
 
 

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.