]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
parallel geometric multigrid: limit the level difference at vertices
authorheister <heister@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 9 Sep 2013 21:21:06 +0000 (21:21 +0000)
committerheister <heister@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 9 Sep 2013 21:21:06 +0000 (21:21 +0000)
git-svn-id: https://svn.dealii.org/trunk@30658 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/source/distributed/tria.cc

index 1af06931300af9d389fbc9e89f65286302e729bb..93145356aada3f2934f65560e65911a6da0e1eea 100644 (file)
@@ -2596,10 +2596,20 @@ namespace parallel
       // smoothing flag is used in the normal
       // refinement process.
       typename Triangulation<dim,spacedim>::MeshSmoothing
-      save_smooth = this->smooth_grid;
-      this->smooth_grid = dealii::Triangulation<dim,spacedim>::none;
-      bool mesh_changed = false;
+        save_smooth = this->smooth_grid;
+
+      // We will refine manually to match the p4est further down, which
+      // obeys a level difference of 2 at each vertex (see the balance call
+      // to p4est). We can disable this here so we store fewer artificial
+      // cells (in some cases). For geometric multigrid it turns out that
+      // we will miss level cells at shared vertices if we ignore this.
+      // See tests/mpi/mg_06.
+      if (settings & construct_multigrid_hierarchy)
+        this->smooth_grid = dealii::Triangulation<dim,spacedim>::none;
+      else
+        this->smooth_grid = dealii::Triangulation<dim,spacedim>::limit_level_difference_at_vertices;
 
+      bool mesh_changed = false;
 
       // remove all deal.II refinements. Note that we could skip this and
       // start from our current state, because the algorithm later coarsens as

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.