From: Daniel Arndt Date: Tue, 6 Sep 2016 15:40:13 +0000 (+0200) Subject: Enforce limit_level_difference_at_vertices if construct_multigrid_hierarchy is set X-Git-Tag: v8.5.0-rc1~697^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=849245cac8b230ec8193f46c162c059939ba1b6b;p=dealii.git Enforce limit_level_difference_at_vertices if construct_multigrid_hierarchy is set --- diff --git a/include/deal.II/distributed/tria.h b/include/deal.II/distributed/tria.h index d09927274f..960333b904 100644 --- a/include/deal.II/distributed/tria.h +++ b/include/deal.II/distributed/tria.h @@ -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 construct_multigrid_hierarchy enforces + * Triangulation::limit_level_difference_at_vertices + * for smooth_grid. * * @note This class does not currently support the * check_for_distorted_cells argument provided by the base diff --git a/include/deal.II/grid/tria.h b/include/deal.II/grid/tria.h index db5d54f6e1..3a8f3252be 100644 --- a/include/deal.II/grid/tria.h +++ b/include/deal.II/grid/tria.h @@ -1555,7 +1555,7 @@ public: /** * Return the mesh smoothing requirements that are obeyed. */ - virtual const Triangulation::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 diff --git a/source/distributed/tria.cc b/source/distributed/tria.cc index 356ee3c52c..7eadc34501 100644 --- a/source/distributed/tria.cc +++ b/source/distributed/tria.cc @@ -2151,9 +2151,14 @@ namespace parallel const typename dealii::Triangulation::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 (mpi_communicator, + (settings_ & construct_multigrid_hierarchy) ? + static_cast::MeshSmoothing> + (smooth_grid | Triangulation::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::limit_level_difference_at_vertices; else