]> https://gitweb.dealii.org/ - dealii.git/commitdiff
p:d:GridRefinement: accept other types of triangulations 14072/head
authorPeter Munch <peterrmuench@gmail.com>
Tue, 28 Jun 2022 19:46:26 +0000 (21:46 +0200)
committerPeter Munch <peterrmuench@gmail.com>
Thu, 30 Jun 2022 17:36:20 +0000 (19:36 +0200)
include/deal.II/distributed/grid_refinement.h
source/distributed/grid_refinement.cc
source/distributed/grid_refinement.inst.in

index 96e73eaff9e7242fa323889adfbd60fb9e436a5e..bf2adf4c5f65ac433e2bc87ed8819e47334e3b92 100644 (file)
@@ -98,16 +98,18 @@ namespace parallel
      * cycle of the adaptive finite element loop.
      *
      * In contrast to the functions in namespace dealii::GridRefinement,
-     * the functions in the current namespace are intended for distributed
-     * meshes, i.e., objects of type parallel::distributed::Triangulation.
+     * the functions in the current namespace are intended for parallel
+     * computations, i.e., computations, e.g., on objects of type
+     * parallel::distributed::Triangulation.
      *
      * @ingroup grid
      */
     namespace GridRefinement
     {
       /**
-       * Like dealii::GridRefinement::refine_and_coarsen_fixed_number, but for
-       * parallel distributed triangulations.
+       * Like dealii::GridRefinement::refine_and_coarsen_fixed_number, but
+       * designed for parallel computations, where each process has only
+       * information about locally owned cells.
        *
        * The vector of criteria needs to be a vector of refinement criteria
        * for all cells active on the current triangulation, i.e.,
@@ -152,8 +154,8 @@ namespace parallel
       template <int dim, typename Number, int spacedim>
       void
       refine_and_coarsen_fixed_number(
-        parallel::distributed::Triangulation<dim, spacedim> &tria,
-        const dealii::Vector<Number> &                       criteria,
+        Triangulation<dim, spacedim> & tria,
+        const dealii::Vector<Number> & criteria,
         const double                   top_fraction_of_cells,
         const double                   bottom_fraction_of_cells,
         const types::global_cell_index max_n_cells =
@@ -161,7 +163,8 @@ namespace parallel
 
       /**
        * Like dealii::GridRefinement::refine_and_coarsen_fixed_fraction, but
-       * for parallel distributed triangulations.
+       * designed for parallel computations, where each process only has
+       * information about locally owned cells.
        *
        * The vector of criteria needs to be a vector of refinement criteria
        * for all cells active on the current triangulation, i.e.,
@@ -205,10 +208,10 @@ namespace parallel
       template <int dim, typename Number, int spacedim>
       void
       refine_and_coarsen_fixed_fraction(
-        parallel::distributed::Triangulation<dim, spacedim> &tria,
-        const dealii::Vector<Number> &                       criteria,
-        const double                top_fraction_of_error,
-        const double                bottom_fraction_of_error,
+        Triangulation<dim, spacedim> &tria,
+        const dealii::Vector<Number> &criteria,
+        const double                  top_fraction_of_error,
+        const double                  bottom_fraction_of_error,
         const VectorTools::NormType norm_type = VectorTools::NormType::L1_norm);
     } // namespace GridRefinement
   }   // namespace distributed
index f5100915a23a158909c7327767f41443df13be5a..6a9a5ceff05e94d45ea3c73daee60376db6f7f80 100644 (file)
@@ -42,6 +42,18 @@ DEAL_II_NAMESPACE_OPEN
 
 namespace
 {
+  template <int dim, int spacedim>
+  unsigned int
+  n_locally_owned_active_cells(const Triangulation<dim, spacedim> &tria)
+  {
+    if (const auto parallel_tria =
+          dynamic_cast<const parallel::TriangulationBase<dim, spacedim> *>(
+            &tria))
+      return parallel_tria->n_locally_owned_active_cells();
+    else
+      return tria.n_active_cells();
+  }
+
   template <typename number>
   inline number
   max_element(const dealii::Vector<number> &criteria)
@@ -107,7 +119,7 @@ namespace
     Vector<Number> &locally_owned_indicators)
   {
     Assert(locally_owned_indicators.size() ==
-             tria.n_locally_owned_active_cells(),
+             n_locally_owned_active_cells(tria),
            ExcInternalError());
 
     unsigned int owned_index = 0;
@@ -118,7 +130,7 @@ namespace
           criteria(cell->active_cell_index());
         ++owned_index;
       }
-    Assert(owned_index == tria.n_locally_owned_active_cells(),
+    Assert(owned_index == n_locally_owned_active_cells(tria),
            ExcInternalError());
   }
 
@@ -207,8 +219,7 @@ namespace
   {
     // first extract from the vector of indicators the ones that correspond
     // to cells that we locally own
-    Vector<Number> locally_owned_indicators(
-      tria.n_locally_owned_active_cells());
+    Vector<Number> locally_owned_indicators(n_locally_owned_active_cells(tria));
     get_locally_owned_indicators(tria, criteria, locally_owned_indicators);
 
     MPI_Comm mpi_communicator = tria.get_communicator();
@@ -520,7 +531,7 @@ namespace parallel
         // first extract from the vector of indicators the ones that correspond
         // to cells that we locally own
         Vector<Number> locally_owned_indicators(
-          tria.n_locally_owned_active_cells());
+          n_locally_owned_active_cells(tria));
         get_locally_owned_indicators(tria, criteria, locally_owned_indicators);
 
         MPI_Comm mpi_communicator = tria.get_communicator();
index 48ab446542e1297b44d99f36f00d3cb6970b1c06..8f3a7e025dd8db231f6765e4a32b1efdbc3c0f4f 100644 (file)
@@ -65,7 +65,7 @@ for (S : REAL_SCALARS; deal_II_dimension : DIMENSIONS)
           refine_and_coarsen_fixed_number<deal_II_dimension,
                                           S,
                                           deal_II_dimension>(
-            parallel::distributed::Triangulation<deal_II_dimension> &,
+            Triangulation<deal_II_dimension> &,
             const dealii::Vector<S> &,
             const double,
             const double,
@@ -75,7 +75,7 @@ for (S : REAL_SCALARS; deal_II_dimension : DIMENSIONS)
           refine_and_coarsen_fixed_fraction<deal_II_dimension,
                                             S,
                                             deal_II_dimension>(
-            parallel::distributed::Triangulation<deal_II_dimension> &,
+            Triangulation<deal_II_dimension> &,
             const dealii::Vector<S> &,
             const double,
             const double,
@@ -97,8 +97,7 @@ for (S : REAL_SCALARS; deal_II_dimension : DIMENSIONS)
           refine_and_coarsen_fixed_number<deal_II_dimension - 1,
                                           S,
                                           deal_II_dimension>(
-            parallel::distributed::Triangulation<deal_II_dimension - 1,
-                                                 deal_II_dimension> &,
+            Triangulation<deal_II_dimension - 1, deal_II_dimension> &,
             const dealii::Vector<S> &,
             const double,
             const double,
@@ -108,8 +107,7 @@ for (S : REAL_SCALARS; deal_II_dimension : DIMENSIONS)
           refine_and_coarsen_fixed_fraction<deal_II_dimension - 1,
                                             S,
                                             deal_II_dimension>(
-            parallel::distributed::Triangulation<deal_II_dimension - 1,
-                                                 deal_II_dimension> &,
+            Triangulation<deal_II_dimension - 1, deal_II_dimension> &,
             const dealii::Vector<S> &,
             const double,
             const double,

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.