]> https://gitweb.dealii.org/ - dealii.git/commitdiff
changes 2
authortcclevenger <tcleven@clemson.edu>
Fri, 17 Mar 2017 13:32:39 +0000 (09:32 -0400)
committertcclevenger <tcleven@clemson.edu>
Fri, 17 Mar 2017 13:32:39 +0000 (09:32 -0400)
include/deal.II/grid/grid_tools.h
source/distributed/shared_tria.cc
source/grid/grid_tools.cc
source/grid/grid_tools.inst.in
tests/sharedtria/tria_multigrid_02.cc
tests/sharedtria/tria_multigrid_02.mpirun=3.output
tests/sharedtria/tria_zorder_02.cc
tests/sharedtria/tria_zorder_02.mpirun=3.output

index 554fb1347dab24684a31bdb25626a74a2fcf0ac9..952d74a3bfee5847b373ebab5d98c44a606b1e9f 100644 (file)
@@ -715,16 +715,18 @@ namespace GridTools
 
 
   /**
-   * Extract and return the layer of cells around a subdomain on a given level
-   * in the @p mesh (i.e. those that share a common set of
-   * vertices with the level subdomain but are not a part of it).
+   * Extract and return the cell layer around a subdomain (set of
+   * cells) on a specified level of the @p mesh (i.e. those cells on
+   * that level that share a common set of vertices with the subdomain
+   * but are not a part of it). Here, the "subdomain" consists of exactly
+   * all of those cells for which the @p predicate returns @p true.
    */
   template <class MeshType>
   std::vector<typename MeshType::cell_iterator>
   compute_cell_halo_layer_on_level
-  (const MeshType     &mesh,
-   const unsigned int my_level_subdomain_id,
-   const unsigned int level);
+  (const MeshType                                                             &mesh,
+   const std_cxx11::function<bool (const typename MeshType::cell_iterator &)> &predicate,
+   const unsigned int                                                         level);
 
 
   /**
index fba438e439bbe42d97a3644568a3d0dbd80e9677..31903f812b81389c23fbe3265f302fb8cdddf3a6 100644 (file)
@@ -50,50 +50,6 @@ namespace parallel
 
 
 
-    namespace
-    {
-      /**
-       *  Helper function for partition() which determines halo
-       * layer cells for a given level
-       */
-      template <int dim, int spacedim>
-      std::vector<typename parallel::shared::Triangulation<dim,spacedim>::cell_iterator>
-      compute_cell_halo_layer_on_level
-      (parallel::shared::Triangulation<dim,spacedim>                                                                   &tria,
-       const std_cxx11::function<bool (const typename parallel::shared::Triangulation<dim,spacedim>::cell_iterator &)> &predicate,
-       const unsigned int                                                                                              level)
-      {
-        std::vector<typename parallel::shared::Triangulation<dim,spacedim>::cell_iterator> level_halo_layer;
-        std::vector<bool> locally_active_vertices_on_level_subdomain (tria.n_vertices(), false);
-
-        // Find the cells for which the predicate is true
-        // These are the cells around which we wish to construct
-        // the halo layer
-        for (typename parallel::shared::Triangulation<dim,spacedim>::cell_iterator
-             cell = tria.begin(level);
-             cell != tria.end(level); ++cell)
-          if (predicate(cell)) // True predicate -> part of subdomain
-            for (unsigned int v=0; v<GeometryInfo<dim>::vertices_per_cell; ++v)
-              locally_active_vertices_on_level_subdomain[cell->vertex_index(v)] = true;
-
-        // Find the cells that do not conform to the predicate
-        // but share a vertex with the selected level subdomain
-        // These comprise the halo layer
-        for (typename parallel::shared::Triangulation<dim,spacedim>::cell_iterator
-             cell = tria.begin(level);
-             cell != tria.end(level); ++cell)
-          if (!predicate(cell)) // False predicate -> possible halo layer cell
-            for (unsigned int v=0; v<GeometryInfo<dim>::vertices_per_cell; ++v)
-              if (locally_active_vertices_on_level_subdomain[cell->vertex_index(v)] == true)
-                {
-                  level_halo_layer.push_back(cell);
-                  break;
-                }
-
-        return level_halo_layer;
-      }
-    }
-
     template <int dim, int spacedim>
     void Triangulation<dim,spacedim>::partition()
     {
@@ -114,7 +70,8 @@ namespace parallel
           AssertThrow(false, ExcInternalError())
         }
 
-      // custom partition require manual partitioning of level cells
+      // do not partition multigrid levels if user is
+      // defining a custom partition
       if ((settings & construct_multigrid_hierarchy) && !(settings & partition_custom_signal))
         dealii::GridTools::partition_multigrid_levels(*this);
 
@@ -155,7 +112,7 @@ namespace parallel
               for (unsigned int lvl=0; lvl<this->n_levels(); ++lvl)
                 {
                   const std::vector<typename parallel::shared::Triangulation<dim,spacedim>::cell_iterator>
-                  level_halo_layer_vector = compute_cell_halo_layer_on_level (*this, predicate, lvl);
+                  level_halo_layer_vector = GridTools::compute_cell_halo_layer_on_level (*this, predicate, lvl);
                   std::set<typename parallel::shared::Triangulation<dim,spacedim>::cell_iterator>
                   level_halo_layer(level_halo_layer_vector.begin(), level_halo_layer_vector.end());
 
@@ -164,10 +121,11 @@ namespace parallel
                   endc = this->end(lvl);
                   for (; cell != endc; cell++)
                     {
-                      // for active cells we must keep level subdomain id of all neighbors,
-                      // not just neighbors that exist on the same level.
+                      // for active cells we must keep level subdomain id of all neighbors
+                      // to our subdomain, not just cells that share a vertex on the same level.
                       // if the cells subdomain id was not set to artitficial above, we will
-                      // also keep its level subdomain id.
+                      // also keep its level subdomain id since it is either owned by this processor
+                      // or in the ghost layer of the active mesh.
                       if (!cell->has_children() && cell->subdomain_id() != numbers::artificial_subdomain_id)
                         continue;
                       if (!cell->is_locally_owned_on_level() &&
index f4989746a87a085604223b41c96932bbee850cd2..373b279613d7691e36ed132ca2334cdb2af90a35 100644 (file)
@@ -1607,6 +1607,46 @@ next_cell:
 
 
 
+  template <class MeshType>
+  std::vector<typename MeshType::cell_iterator>
+  compute_cell_halo_layer_on_level
+  (const MeshType                                                             &mesh,
+   const std_cxx11::function<bool (const typename MeshType::cell_iterator &)> &predicate,
+   const unsigned int                                                         level)
+  {
+    std::vector<typename MeshType::cell_iterator> level_halo_layer;
+    std::vector<bool> locally_active_vertices_on_level_subdomain (mesh.get_triangulation().n_vertices(),
+        false);
+
+    // Find the cells for which the predicate is true
+    // These are the cells around which we wish to construct
+    // the halo layer
+    for (typename MeshType::cell_iterator
+         cell = mesh.begin(level);
+         cell != mesh.end(level); ++cell)
+      if (predicate(cell)) // True predicate --> Part of subdomain
+        for (unsigned int v=0; v<GeometryInfo<MeshType::dimension>::vertices_per_cell; ++v)
+          locally_active_vertices_on_level_subdomain[cell->vertex_index(v)] = true;
+
+    // Find the cells that do not conform to the predicate
+    // but share a vertex with the selected subdomain on that level
+    // These comprise the halo layer
+    for (typename MeshType::cell_iterator
+         cell = mesh.begin(level);
+         cell != mesh.end(level); ++cell)
+      if (!predicate(cell)) // False predicate --> Potential halo cell
+        for (unsigned int v=0; v<GeometryInfo<MeshType::dimension>::vertices_per_cell; ++v)
+          if (locally_active_vertices_on_level_subdomain[cell->vertex_index(v)] == true)
+            {
+              level_halo_layer.push_back(cell);
+              break;
+            }
+
+    return level_halo_layer;
+  }
+
+
+
   template <class MeshType>
   std::vector<typename MeshType::active_cell_iterator>
   compute_ghost_cell_halo_layer (const MeshType &mesh)
index 8bc7c7cf3561aae3c4aa81f177624f97118aa004..b8fbeba4295e7511d7bf928d9931d2dec1c256a2 100644 (file)
@@ -45,6 +45,12 @@ for (X : TRIANGULATION_AND_DOFHANDLERS; deal_II_dimension : DIMENSIONS ; deal_II
     compute_active_cell_halo_layer (const X &,
                                     const std_cxx11::function<bool (const dealii::internal::ActiveCellIterator<deal_II_dimension, deal_II_space_dimension, X>::type&)> &);
 
+    template
+    std::vector<X::cell_iterator>
+    compute_cell_halo_layer_on_level (const X &,
+                                      const std_cxx11::function<bool (const X::cell_iterator&)> &,
+                                      const unsigned int);
+
     template
     std::vector<dealii::internal::ActiveCellIterator<deal_II_dimension, deal_II_space_dimension, X>::type>
     compute_ghost_cell_halo_layer (const X &);
index 0e92e9d2aa4b320c640510923e134a099cbe694d..de3a0ae24296dd90f24f7e1979df8ac4ea33d05d 100644 (file)
@@ -53,10 +53,13 @@ void test()
     if (cell->center().norm() > 0.3 && cell->center().norm() < 0.42)
       cell->set_refine_flag();
   shared_tria.execute_coarsening_and_refinement();
-  for (typename Triangulation<dim>::active_cell_iterator cell=shared_tria.begin_active(); cell != shared_tria.end(); ++cell)
-    if (cell->at_boundary() && (cell->center()[0] < 0 || cell->center()[1] < 0))
-      cell->set_refine_flag();
-  shared_tria.execute_coarsening_and_refinement();
+  if (dim != 1)
+    {
+      for (typename Triangulation<dim>::active_cell_iterator cell=shared_tria.begin_active(); cell != shared_tria.end(); ++cell)
+        if (cell->at_boundary() && (cell->center()[0] < 0 || cell->center()[1] < 0))
+          cell->set_refine_flag();
+      shared_tria.execute_coarsening_and_refinement();
+    }
 
   deallog << "(CellId,level_subdomain_id) for each active cell:" << std::endl;
   for (unsigned int lvl=0; lvl<shared_tria.n_levels(); ++lvl)
@@ -75,6 +78,9 @@ int main(int argc, char *argv[])
   Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
   MPILogInitAll all;
 
+  deallog.push("1d");
+  test<1>();
+  deallog.pop();
   deallog.push("2d");
   test<2>();
   deallog.pop();
index aa5f75a3ab69dea7ad232717804f461f351da721..06ad6d31684129e34cb2feb439cb1cd868dd31e8 100644 (file)
@@ -1,4 +1,19 @@
 
+DEAL:0:1d::(CellId,level_subdomain_id) for each active cell:
+DEAL:0:1d::(0_0:,0)
+DEAL:0:1d::(1_0:,1)
+DEAL:0:1d::(0_1:0,0)
+DEAL:0:1d::(0_1:1,0)
+DEAL:0:1d::(1_1:0,1)
+DEAL:0:1d::(0_2:00,0)
+DEAL:0:1d::(0_2:01,0)
+DEAL:0:1d::(0_2:10,0)
+DEAL:0:1d::(0_2:11,1)
+DEAL:0:1d::(0_3:100,0)
+DEAL:0:1d::(0_3:101,0)
+DEAL:0:1d::(0_3:110,1)
+DEAL:0:1d::(0_4:1010,0)
+DEAL:0:1d::(0_4:1011,0)
 DEAL:0:2d::(CellId,level_subdomain_id) for each active cell:
 DEAL:0:2d::(0_0:,0)
 DEAL:0:2d::(1_0:,0)
@@ -1756,6 +1771,24 @@ DEAL:0:3d::(7_4:0100,2)
 DEAL:0:3d::(7_4:0101,2)
 DEAL:0::OK
 
+DEAL:1:1d::(CellId,level_subdomain_id) for each active cell:
+DEAL:1:1d::(0_0:,0)
+DEAL:1:1d::(1_0:,1)
+DEAL:1:1d::(0_1:1,0)
+DEAL:1:1d::(1_1:0,1)
+DEAL:1:1d::(1_1:1,2)
+DEAL:1:1d::(0_2:10,0)
+DEAL:1:1d::(0_2:11,1)
+DEAL:1:1d::(1_2:00,1)
+DEAL:1:1d::(1_2:01,2)
+DEAL:1:1d::(0_3:101,0)
+DEAL:1:1d::(0_3:110,1)
+DEAL:1:1d::(0_3:111,1)
+DEAL:1:1d::(1_3:000,1)
+DEAL:1:1d::(1_3:001,1)
+DEAL:1:1d::(1_3:010,2)
+DEAL:1:1d::(0_4:1011,0)
+DEAL:1:1d::(1_4:0100,2)
 DEAL:1:2d::(CellId,level_subdomain_id) for each active cell:
 DEAL:1:2d::(0_0:,0)
 DEAL:1:2d::(1_0:,0)
@@ -3829,6 +3862,18 @@ DEAL:1:3d::(7_4:0404,2)
 DEAL:1::OK
 
 
+DEAL:2:1d::(CellId,level_subdomain_id) for each active cell:
+DEAL:2:1d::(1_1:0,1)
+DEAL:2:1d::(1_1:1,2)
+DEAL:2:1d::(1_2:00,1)
+DEAL:2:1d::(1_2:01,2)
+DEAL:2:1d::(1_2:10,2)
+DEAL:2:1d::(1_2:11,2)
+DEAL:2:1d::(1_3:001,1)
+DEAL:2:1d::(1_3:010,2)
+DEAL:2:1d::(1_3:011,2)
+DEAL:2:1d::(1_4:0100,2)
+DEAL:2:1d::(1_4:0101,2)
 DEAL:2:2d::(CellId,level_subdomain_id) for each active cell:
 DEAL:2:2d::(0_0:,0)
 DEAL:2:2d::(1_0:,0)
index 7af1191c7ada7be64c4fff24d85b9225585fc781..919d6645bbaf205de11cc20541738debc1c2f192 100644 (file)
@@ -52,10 +52,13 @@ void test()
     if (cell->center().norm() > 0.3 && cell->center().norm() < 0.42)
       cell->set_refine_flag();
   shared_tria.execute_coarsening_and_refinement();
-  for (typename Triangulation<dim>::active_cell_iterator cell=shared_tria.begin_active(); cell != shared_tria.end(); ++cell)
-    if (cell->at_boundary() && (cell->center()[0] < 0 || cell->center()[1] < 0))
-      cell->set_refine_flag();
-  shared_tria.execute_coarsening_and_refinement();
+  if (dim != 1)
+    {
+      for (typename Triangulation<dim>::active_cell_iterator cell=shared_tria.begin_active(); cell != shared_tria.end(); ++cell)
+        if (cell->at_boundary() && (cell->center()[0] < 0 || cell->center()[1] < 0))
+          cell->set_refine_flag();
+      shared_tria.execute_coarsening_and_refinement();
+    }
 
   deallog << "(CellId,subdomain_id) for each active cell:" << std::endl;
   typename Triangulation<dim>::active_cell_iterator
@@ -71,6 +74,9 @@ int main(int argc, char *argv[])
   Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
   MPILogInitAll all;
 
+  deallog.push("1d");
+  test<1>();
+  deallog.pop();
   deallog.push("2d");
   test<2>();
   deallog.pop();
index bf66b82b53963cfd758ce8a5a8595d5cf862788c..4e9ad0d4657e4256c5df654c98cab1bb617fc044 100644 (file)
@@ -1,4 +1,10 @@
 
+DEAL:0:1d::(CellId,subdomain_id) for each active cell:
+DEAL:0:1d::(0_2:00,0)
+DEAL:0:1d::(0_2:01,0)
+DEAL:0:1d::(0_3:100,0)
+DEAL:0:1d::(0_4:1010,0)
+DEAL:0:1d::(0_4:1011,0)
 DEAL:0:2d::(CellId,subdomain_id) for each active cell:
 DEAL:0:2d::(0_2:03,0)
 DEAL:0:2d::(0_2:12,0)
@@ -1161,6 +1167,11 @@ DEAL:0:3d::(1_4:6766,0)
 DEAL:0:3d::(1_4:6767,0)
 DEAL:0::OK
 
+DEAL:1:1d::(CellId,subdomain_id) for each active cell:
+DEAL:1:1d::(0_3:110,1)
+DEAL:1:1d::(0_3:111,1)
+DEAL:1:1d::(1_3:000,1)
+DEAL:1:1d::(1_3:001,1)
 DEAL:1:2d::(CellId,subdomain_id) for each active cell:
 DEAL:1:2d::(1_2:12,1)
 DEAL:1:2d::(1_2:30,1)
@@ -2318,6 +2329,12 @@ DEAL:1:3d::(4_4:3737,1)
 DEAL:1::OK
 
 
+DEAL:2:1d::(CellId,subdomain_id) for each active cell:
+DEAL:2:1d::(1_2:10,2)
+DEAL:2:1d::(1_2:11,2)
+DEAL:2:1d::(1_3:011,2)
+DEAL:2:1d::(1_4:0100,2)
+DEAL:2:1d::(1_4:0101,2)
 DEAL:2:2d::(CellId,subdomain_id) for each active cell:
 DEAL:2:2d::(2_2:21,2)
 DEAL:2:2d::(2_2:30,2)

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.