add function to get max possible mg coarse level

processor

diff --git a/include/deal.II/multigrid/mg_tools.h b/include/deal.II/multigrid/mg_tools.h
index fdcf167e2d0de5a2f44cd10a3a077354b45d98b3..fddc736ed8ec7d068468b064866b28eaa394d1b8 100644 (file)
--- a/include/deal.II/multigrid/mg_tools.h
+++ b/include/deal.II/multigrid/mg_tools.h
@@ -247,6 +247,18 @@ namespace MGTools
   void
   extract_non_interface_dofs (const DoFHandler<dim,spacedim> &mg_dof_handler,
                               std::vector<std::set<types::global_dof_index> > &non_interface_dofs);
+
+  /**
+   * Return the highest possible level that can be used as the coarsest level in
+   * a Multigrid computation, that is, the highest level in the hierarchy whose mesh
+   * covers the entire domain. This corresponds to the minimum level of a cell on
+   * the active mesh. Since each processor only has a local view of the mesh, each
+   * processor must call this function. Note that this is a global minimum over the
+   * entire mesh and therefore each processor will return the same value.
+   */
+  template <int dim, int spacedim>
+  unsigned int
+  max_level_for_coarse_mesh (const Triangulation<dim,spacedim> &tria);
 }
 
 /* @} */

diff --git a/source/multigrid/mg_tools.cc b/source/multigrid/mg_tools.cc
index b0ec9651fea6fe3a9f47bf98e72ef0f61ab88028..d3c44bc18f48a0b23af528884ebdfe34f5c0f935 100644 (file)
--- a/source/multigrid/mg_tools.cc
+++ b/source/multigrid/mg_tools.cc
@@ -1547,6 +1547,40 @@ namespace MGTools
       }
 
   }
+
+
+
+  template <int dim, int spacedim>
+  unsigned int
+  max_level_for_coarse_mesh (const Triangulation<dim,spacedim> &tria)
+  {
+    // Find minimum level for an active cell in
+    // this locally owned subdomain
+    // Note: with the way active cells are traversed,
+    // the first locally owned cell we find will have
+    // the lowest level in the particular subdomain.
+    unsigned int min_level = tria.n_global_levels();
+    typename Triangulation<dim,spacedim>::active_cell_iterator cell = tria.begin_active(),
+                                                               endc = tria.end();
+    for (; cell!=endc; ++cell)
+      if (cell->is_locally_owned())
+        {
+          min_level = cell->level();
+          break;
+        }
+
+    unsigned int global_min = min_level;
+    // If necessary, communicate to find minimum
+    // level for an active cell over all subdomains
+    if (const parallel::Triangulation<dim,spacedim> *tr
+        = dynamic_cast<const parallel::Triangulation<dim,spacedim> *>
+        (&tria))
+      global_min = Utilities::MPI::min(min_level, tr->get_communicator());
+
+    AssertIndexRange(global_min, tria.n_global_levels());
+
+    return global_min;
+  }
 }
 
 

diff --git a/source/multigrid/mg_tools.inst.in b/source/multigrid/mg_tools.inst.in
index 45e5c10cad58e7a0904022faf7e0a7a2f02cb576..b42950abf5b3fa7472686d7d546d34d077cf0764 100644 (file)
--- a/source/multigrid/mg_tools.inst.in
+++ b/source/multigrid/mg_tools.inst.in
@@ -126,6 +126,7 @@ for (deal_II_dimension : DIMENSIONS)
     extract_non_interface_dofs (const DoFHandler<deal_II_dimension> & mg_dof_handler,
                                 std::vector<std::set<types::global_dof_index> > &non_interface_dofs);
 
+
 #if deal_II_dimension < 3
     template void count_dofs_per_block (
         const DoFHandler<deal_II_dimension,deal_II_dimension+1>&,
@@ -140,3 +141,15 @@ for (deal_II_dimension : DIMENSIONS)
     \}
 }
 
+for (deal_II_dimension : DIMENSIONS;
+        deal_II_space_dimension: SPACE_DIMENSIONS)
+{
+    namespace MGTools
+    \{
+#if deal_II_dimension <= deal_II_space_dimension
+    template
+    unsigned int
+    max_level_for_coarse_mesh (const Triangulation<deal_II_dimension, deal_II_space_dimension> & tria);
+#endif
+    \}
+}

diff --git a/tests/multigrid/max_level_for_coarse_mesh.cc b/tests/multigrid/max_level_for_coarse_mesh.cc
new file mode 100644 (file)
index 0000000..dead182
--- /dev/null
+++ b/tests/multigrid/max_level_for_coarse_mesh.cc
@@ -0,0 +1,87 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2006 - 2017 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+
+// check function MGTools::max_level_for_coarse_mesh()
+
+#include "../tests.h"
+#include <deal.II/base/utilities.h>
+
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/tria_iterator.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_out.h>
+
+#include <deal.II/numerics/vector_tools.h>
+
+#include <deal.II/base/index_set.h>
+#include <deal.II/distributed/tria.h>
+#include <deal.II/distributed/grid_refinement.h>
+
+#include <deal.II/multigrid/mg_constrained_dofs.h>
+#include <deal.II/multigrid/multigrid.h>
+
+template <int dim>
+void test ()
+{
+  parallel::distributed::Triangulation<dim> tria(MPI_COMM_WORLD,Triangulation<dim>::
+                                                 limit_level_difference_at_vertices,
+                                                 parallel::distributed::Triangulation<dim>::construct_multigrid_hierarchy);
+  GridGenerator::hyper_cube(tria,0,1);
+  tria.refine_global(1);
+  for (unsigned int cycle=0; cycle<2; ++cycle)
+    {
+      for (typename parallel::distributed::Triangulation<dim>::active_cell_iterator cell = tria.begin_active();
+           cell != tria.end(); ++cell)
+        for (unsigned int v=0; v<GeometryInfo<dim>::vertices_per_cell; ++v)
+          {
+            if (dim == 2)
+              if (cell->vertex(v)[0] < 0.25 && cell->vertex(v)[1] < 0.25)
+                {
+                  cell->set_refine_flag();
+                  break;
+                }
+            if (dim == 3)
+              if (cell->vertex(v)[0] < 0.25 && cell->vertex(v)[1] < 0.25 && cell->vertex(v)[2] < 0.25)
+                {
+                  cell->set_refine_flag();
+                  break;
+                }
+          }
+      tria.execute_coarsening_and_refinement();
+    }
+
+  const unsigned int max_possible_level = MGTools::max_level_for_coarse_mesh(tria);
+  deallog << "Max possible level for coarse mesh: "
+          << max_possible_level << std::endl;
+}
+
+
+
+int main (int argc, char *argv[])
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
+  MPILogInitAll all;
+
+  deallog.push("2d");
+  test<2>();
+  deallog.pop();
+  deallog.push("3d");
+  test<3>();
+  deallog.pop();
+
+  deallog << "OK" << std::endl;
+}

diff --git a/tests/multigrid/max_level_for_coarse_mesh.with_mpi=true.with_p4est=true.mpirun=3.output b/tests/multigrid/max_level_for_coarse_mesh.with_mpi=true.with_p4est=true.mpirun=3.output
new file mode 100644 (file)
index 0000000..0537bfd
--- /dev/null
+++ b/tests/multigrid/max_level_for_coarse_mesh.with_mpi=true.with_p4est=true.mpirun=3.output
@@ -0,0 +1,14 @@
+
+DEAL:0:2d::Max possible level for coarse mesh: 1
+DEAL:0:3d::Max possible level for coarse mesh: 1
+DEAL:0::OK
+
+DEAL:1:2d::Max possible level for coarse mesh: 1
+DEAL:1:3d::Max possible level for coarse mesh: 1
+DEAL:1::OK
+
+
+DEAL:2:2d::Max possible level for coarse mesh: 1
+DEAL:2:3d::Max possible level for coarse mesh: 1
+DEAL:2::OK
+