]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Enforce limit_level_difference_at_vertices if construct_multigrid_hierarchy is set 3068/head
authorDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Tue, 6 Sep 2016 15:40:13 +0000 (17:40 +0200)
committerDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Tue, 6 Sep 2016 23:34:13 +0000 (01:34 +0200)
include/deal.II/distributed/tria.h
include/deal.II/grid/tria.h
source/distributed/tria.cc

index d09927274fecc970a3b636aa7a6ab7ea6aaece02..960333b90455d17cf8929bf1e0de37d16fba3406 100644 (file)
@@ -414,6 +414,9 @@ namespace parallel
        * the kinds of smoothing operations that can be applied.
        *
        * @param settings See the description of the Settings enumerator.
+       * Providing <code>construct_multigrid_hierarchy</code> enforces
+       * <code>Triangulation::limit_level_difference_at_vertices</code>
+       * for smooth_grid.
        *
        * @note This class does not currently support the
        * <code>check_for_distorted_cells</code> argument provided by the base
index db5d54f6e15f0dfaf3a1c3c3dc7f4f0d358979cd..3a8f3252beb0b087c35ee15b317899f45cd2e130 100644 (file)
@@ -1555,7 +1555,7 @@ public:
   /**
    * Return the mesh smoothing requirements that are obeyed.
    */
-  virtual const Triangulation<dim, spacedim>::MeshSmoothing &get_mesh_smoothing() const;
+  virtual const MeshSmoothing &get_mesh_smoothing() const;
 
   /**
    * If @p dim==spacedim, assign a boundary object to a certain part of the
index 356ee3c52cf670c26c60760d6b8cc3ee0ba747e2..7eadc3450179ead2b188d842cf7f294fe750ba42 100644 (file)
@@ -2151,9 +2151,14 @@ namespace parallel
                    const typename dealii::Triangulation<dim,spacedim>::MeshSmoothing smooth_grid,
                    const Settings settings_)
       :
-      // do not check for distorted cells
+      // Do not check for distorted cells.
+      // For multigrid, we need limit_level_difference_at_vertices
+      // to make sure the transfer operators only need to consider two levels.
       dealii::parallel::Triangulation<dim,spacedim>
       (mpi_communicator,
+       (settings_ & construct_multigrid_hierarchy) ?
+       static_cast<typename dealii::Triangulation<dim,spacedim>::MeshSmoothing>
+       (smooth_grid | Triangulation<dim,spacedim>::limit_level_difference_at_vertices) :
        smooth_grid,
        false),
       settings(settings_),
@@ -3380,9 +3385,11 @@ namespace parallel
       // 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
+      // 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.
+      // See tests/mpi/mg_06. In particular, the flag is still necessary
+      // even though we force it for the original smooth_grid in the constructor.
       if (settings & construct_multigrid_hierarchy)
         this->smooth_grid = dealii::Triangulation<dim,spacedim>::limit_level_difference_at_vertices;
       else

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.