]> https://gitweb.dealii.org/ - dealii.git/commitdiff
VectorNorm for fixed_fraction refinement. 11023/head
authorMarc Fehling <mafehling.git@gmail.com>
Mon, 5 Oct 2020 02:48:43 +0000 (20:48 -0600)
committerMarc Fehling <mafehling.git@gmail.com>
Tue, 27 Oct 2020 05:46:45 +0000 (23:46 -0600)
Added fixed fraction tests to compare serial/parallel implementations.

17 files changed:
doc/news/changes/minor/20201006Fehling [new file with mode: 0644]
include/deal.II/distributed/grid_refinement.h
include/deal.II/grid/grid_refinement.h
source/distributed/grid_refinement.cc
source/distributed/grid_refinement.inst.in
source/grid/grid_refinement.cc
source/grid/grid_refinement.inst.in
tests/grid/refine_and_coarsen_fixed_fraction_03.cc [new file with mode: 0644]
tests/grid/refine_and_coarsen_fixed_fraction_03.output [new file with mode: 0644]
tests/grid/refine_and_coarsen_fixed_fraction_04.cc [new file with mode: 0644]
tests/grid/refine_and_coarsen_fixed_fraction_04.output [new file with mode: 0644]
tests/mpi/refine_and_coarsen_fixed_fraction_08.cc [new file with mode: 0644]
tests/mpi/refine_and_coarsen_fixed_fraction_08.mpirun=1.with_p4est=true.output [new file with mode: 0644]
tests/mpi/refine_and_coarsen_fixed_fraction_08.mpirun=4.with_p4est=true.output [new file with mode: 0644]
tests/mpi/refine_and_coarsen_fixed_fraction_09.cc [new file with mode: 0644]
tests/mpi/refine_and_coarsen_fixed_fraction_09.mpirun=1.with_p4est=true.output [new file with mode: 0644]
tests/mpi/refine_and_coarsen_fixed_fraction_09.mpirun=4.with_p4est=true.output [new file with mode: 0644]

diff --git a/doc/news/changes/minor/20201006Fehling b/doc/news/changes/minor/20201006Fehling
new file mode 100644 (file)
index 0000000..39a950b
--- /dev/null
@@ -0,0 +1,6 @@
+New: GridRefinement::refine_and_coarsen_fixed_fraction() and
+parallel::distributed::GridRefinement::refine_and_coarsen_fixed_fraction()
+now allow to specify a VectorTools::NormType, which determines how
+combined errors on subsets of cells will be calculated.
+<br>
+(Marc Fehling, 2020/10/06)
index 7d8a1dc0997b3764ef7db85cea8c5b9976c95423..4ff1e2c3480d0402033cc9225a2047d7d43bbb59 100644 (file)
@@ -23,6 +23,8 @@
 
 #include <deal.II/distributed/tria.h>
 
+#include <deal.II/numerics/vector_tools_common.h>
+
 #include <limits>
 #include <vector>
 
@@ -124,6 +126,28 @@ namespace parallel
        * refined at all.
        *
        * The same is true for the fraction of cells that is coarsened.
+       *
+       * @param[in,out] tria The triangulation whose cells this function is
+       * supposed to mark for coarsening and refinement.
+       *
+       * @param[in] criteria The refinement criterion for each mesh cell active
+       * on the current triangulation. Entries may not be negative.
+       *
+       * @param[in] top_fraction_of_cells The fraction of cells to be refined.
+       * If this number is zero, no cells will be refined. If it equals one,
+       * the result will be flagging for global refinement.
+       *
+       * @param[in] bottom_fraction_of_cells The fraction of cells to be
+       * coarsened. If this number is zero, no cells will be coarsened.
+       *
+       * @param[in] max_n_cells This argument can be used to specify a maximal
+       * number of cells. If this number is going to be exceeded upon
+       * refinement, then refinement and coarsening fractions are going to be
+       * adjusted in an attempt to reach the maximum number of cells. Be aware
+       * though that through proliferation of refinement due to
+       * Triangulation::MeshSmoothing, this number is only an indicator. The
+       * default value of this argument is to impose no limit on the number of
+       * cells.
        */
       template <int dim, typename Number, int spacedim>
       void
@@ -157,14 +181,35 @@ namespace parallel
        * processors, no cells are refined at all.
        *
        * The same is true for the fraction of cells that is coarsened.
+       *
+       * @param[in,out] tria The triangulation whose cells this function is
+       * supposed to mark for coarsening and refinement.
+       *
+       * @param[in] criteria The refinement criterion computed on each mesh cell
+       * active on the current triangulation. Entries may not be negative.
+       *
+       * @param[in] top_fraction_of_error The fraction of the total estimate
+       * which should be refined. If this number is zero, no cells will be
+       * refined. If it equals one, the result will be flagging for global
+       * refinement.
+       *
+       * @param[in] bottom_fraction_of_error The fraction of the estimate
+       * coarsened. If this number is zero, no cells will be coarsened.
+       *
+       * @param[in] norm_type To determine thresholds, combined errors on
+       * subsets of cells are calculated as norms of the criteria on these
+       * cells. Different types of norms can be used for this purpose, from
+       * which VectorTools::NormType::L1_norm and
+       * VectorTools::NormType::L2_norm are currently supported.
        */
       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);
+        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
 } // namespace parallel
index 8306ea82d29c3cc775eff944429be532af1ae3a8..0a5741e34103ebce89871ab8fff0f6ea9b90dd63 100644 (file)
@@ -21,6 +21,8 @@
 
 #include <deal.II/base/exceptions.h>
 
+#include <deal.II/numerics/vector_tools_common.h>
+
 #include <limits>
 
 DEAL_II_NAMESPACE_OPEN
@@ -219,6 +221,12 @@ namespace GridRefinement
    * through proliferation of refinement due to Triangulation::MeshSmoothing,
    * this number is only an indicator. The default value of this argument is
    * to impose no limit on the number of cells.
+   *
+   * @param[in] norm_type To determine thresholds, combined errors on
+   * subsets of cells are calculated as norms of the criteria on these
+   * cells. Different types of norms can be used for this purpose, from
+   * which VectorTools::NormType::L1_norm and
+   * VectorTools::NormType::L2_norm are currently supported.
    */
   template <int dim, typename Number, int spacedim>
   void
@@ -227,7 +235,8 @@ namespace GridRefinement
     const Vector<Number> &        criteria,
     const double                  top_fraction,
     const double                  bottom_fraction,
-    const unsigned int max_n_cells = std::numeric_limits<unsigned int>::max());
+    const unsigned int max_n_cells = std::numeric_limits<unsigned int>::max(),
+    const VectorTools::NormType norm_type = VectorTools::NormType::L1_norm);
 
 
 
index 9c438ee9fd3a390365ca6eeb3e21731051f358c1..b3042b6c07549972afde1c13f7d2427cca9f781e 100644 (file)
@@ -186,6 +186,67 @@ namespace
           cell->clear_coarsen_flag();
         }
   }
+
+
+
+  /**
+   * Fixed fraction algorithm without a specified vector norm.
+   *
+   * Entries of the criteria vector and fractions are taken as is, so this
+   * function basically evaluates norms on the vector or its subsets as
+   * l1-norms.
+   */
+  template <int dim, int spacedim, typename Number>
+  void
+  refine_and_coarsen_fixed_fraction_via_l1_norm(
+    parallel::distributed::Triangulation<dim, spacedim> &tria,
+    const dealii::Vector<Number> &                       criteria,
+    const double                                         top_fraction_of_error,
+    const double bottom_fraction_of_error)
+  {
+    // 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());
+    get_locally_owned_indicators(tria, criteria, locally_owned_indicators);
+
+    MPI_Comm mpi_communicator = tria.get_communicator();
+
+    // figure out the global max and min of the indicators. we don't need it
+    // here, but it's a collective communication call
+    const std::pair<double, double> global_min_and_max =
+      dealii::internal::parallel::distributed::GridRefinement::
+        compute_global_min_and_max_at_root(locally_owned_indicators,
+                                           mpi_communicator);
+
+    const double total_error =
+      compute_global_sum(locally_owned_indicators, mpi_communicator);
+
+    double top_target_error    = top_fraction_of_error * total_error,
+           bottom_target_error = (1. - bottom_fraction_of_error) * total_error;
+
+    double top_threshold, bottom_threshold;
+    top_threshold = dealii::internal::parallel::distributed::GridRefinement::
+      RefineAndCoarsenFixedFraction::compute_threshold(locally_owned_indicators,
+                                                       global_min_and_max,
+                                                       top_target_error,
+                                                       mpi_communicator);
+
+    // compute bottom threshold only if necessary. otherwise use the lowest
+    // threshold possible
+    if (bottom_fraction_of_error > 0)
+      bottom_threshold = dealii::internal::parallel::distributed::
+        GridRefinement::RefineAndCoarsenFixedFraction::compute_threshold(
+          locally_owned_indicators,
+          global_min_and_max,
+          bottom_target_error,
+          mpi_communicator);
+    else
+      bottom_threshold = std::numeric_limits<Number>::lowest();
+
+    // now refine the mesh
+    mark_cells(tria, criteria, top_threshold, bottom_threshold);
+  }
 } // namespace
 
 
@@ -499,13 +560,15 @@ 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)
+        const double                top_fraction_of_error,
+        const double                bottom_fraction_of_error,
+        const VectorTools::NormType norm_type)
       {
         Assert(criteria.size() == tria.n_active_cells(),
                ExcDimensionMismatch(criteria.size(), tria.n_active_cells()));
@@ -519,45 +582,44 @@ namespace parallel
         Assert(criteria.is_non_negative(),
                dealii::GridRefinement::ExcNegativeCriteria());
 
-        // 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());
-        get_locally_owned_indicators(tria, criteria, locally_owned_indicators);
-
-        MPI_Comm mpi_communicator = tria.get_communicator();
-
-        // figure out the global max and min of the indicators. we don't need it
-        // here, but it's a collective communication call
-        const std::pair<double, double> global_min_and_max =
-          dealii::internal::parallel::distributed::GridRefinement::
-            compute_global_min_and_max_at_root(locally_owned_indicators,
-                                               mpi_communicator);
-
-        const double total_error =
-          compute_global_sum(locally_owned_indicators, mpi_communicator);
-        double top_threshold, bottom_threshold;
-        top_threshold = dealii::internal::parallel::distributed::
-          GridRefinement::RefineAndCoarsenFixedFraction::compute_threshold(
-            locally_owned_indicators,
-            global_min_and_max,
-            top_fraction_of_error * total_error,
-            mpi_communicator);
-
-        // compute bottom threshold only if necessary. otherwise use the lowest
-        // threshold possible
-        if (bottom_fraction_of_error > 0)
-          bottom_threshold = dealii::internal::parallel::distributed::
-            GridRefinement::RefineAndCoarsenFixedFraction::compute_threshold(
-              locally_owned_indicators,
-              global_min_and_max,
-              (1. - bottom_fraction_of_error) * total_error,
-              mpi_communicator);
-        else
-          bottom_threshold = std::numeric_limits<Number>::lowest();
+        switch (norm_type)
+          {
+            case VectorTools::NormType::L1_norm:
+              // evaluate norms on subsets and compare them as
+              //   c_0 + c_1 + ... < fraction * l1-norm(c)
+              refine_and_coarsen_fixed_fraction_via_l1_norm(
+                tria,
+                criteria,
+                top_fraction_of_error,
+                bottom_fraction_of_error);
+              break;
+
+            case VectorTools::NormType::L2_norm:
+              {
+                // we do not want to evaluate norms on subsets as:
+                //   sqrt(c_0^2 + c_1^2 + ...) < fraction * l2-norm(c)
+                // instead take the square of both sides of the equation
+                // and evaluate:
+                //   c_0^2 + c_1^2 + ... < fraction^2 * l1-norm(c.c)
+                // we adjust all parameters accordingly
+                Vector<Number> criteria_squared(criteria.size());
+                std::transform(criteria.begin(),
+                               criteria.end(),
+                               criteria_squared.begin(),
+                               [](Number c) { return c * c; });
+
+                refine_and_coarsen_fixed_fraction_via_l1_norm(
+                  tria,
+                  criteria_squared,
+                  top_fraction_of_error * top_fraction_of_error,
+                  bottom_fraction_of_error * bottom_fraction_of_error);
+              }
+              break;
 
-        // now refine the mesh
-        mark_cells(tria, criteria, top_threshold, bottom_threshold);
+            default:
+              Assert(false, ExcNotImplemented());
+              break;
+          }
       }
     } // namespace GridRefinement
   }   // namespace distributed
index 00204672eee97c331e5fce04563c77b90eb54c70..c653059f77518cbc0b597db7add9b1322e83a4b9 100644 (file)
@@ -78,7 +78,8 @@ for (S : REAL_SCALARS; deal_II_dimension : DIMENSIONS)
             parallel::distributed::Triangulation<deal_II_dimension> &,
             const dealii::Vector<S> &,
             const double,
-            const double);
+            const double,
+            const VectorTools::NormType);
         \}
       \}
     \}
@@ -109,7 +110,8 @@ for (S : REAL_SCALARS; deal_II_dimension : DIMENSIONS)
                                                  deal_II_dimension> &,
             const dealii::Vector<S> &,
             const double,
-            const double);
+            const double,
+            const VectorTools::NormType);
         \}
       \}
     \}
index 8c3afb22172a3480791e84fba2af090944b3e425..fe0a1c3e0a53a64ca19c4230c13a9da356f8c6e3 100644 (file)
 
 DEAL_II_NAMESPACE_OPEN
 
+namespace
+{
+  /**
+   * Fixed fraction algorithm without a specified vector norm.
+   *
+   * Entries of the criteria vector and fractions are taken as is, so this
+   * function basically evaluates norms on the vector or its subsets as
+   * l1-norms.
+   */
+  template <int dim, int spacedim, typename Number>
+  void
+  refine_and_coarsen_fixed_fraction_via_l1_norm(
+    Triangulation<dim, spacedim> &tria,
+    const Vector<Number> &        criteria,
+    const double                  top_fraction,
+    const double                  bottom_fraction,
+    const unsigned int            max_n_cells)
+  {
+    // sort the criteria in descending order in an auxiliary vector, which we
+    // have to sum up and compare with @p{fraction_of_error*total_error}
+    Vector<Number> criteria_sorted = criteria;
+    std::sort(criteria_sorted.begin(),
+              criteria_sorted.end(),
+              std::greater<double>());
+
+    const double total_error = criteria_sorted.l1_norm();
+
+    // compute thresholds
+    typename Vector<Number>::const_iterator pp = criteria_sorted.begin();
+    for (double sum = 0;
+         (sum < top_fraction * total_error) && (pp != criteria_sorted.end());
+         ++pp)
+      sum += *pp;
+    double top_threshold =
+      (pp != criteria_sorted.begin() ? (*pp + *(pp - 1)) / 2 : *pp);
+
+    typename Vector<Number>::const_iterator qq = criteria_sorted.end() - 1;
+    for (double sum = 0; (sum < bottom_fraction * total_error) &&
+                         (qq != criteria_sorted.begin() - 1);
+         --qq)
+      sum += *qq;
+    double bottom_threshold =
+      ((qq != criteria_sorted.end() - 1) ? (*qq + *(qq + 1)) / 2 : 0.);
+
+    // we now have an idea how many cells we
+    // are going to refine and coarsen. we use
+    // this information to see whether we are
+    // over the limit and if so use a function
+    // that knows how to deal with this
+    // situation
+
+    // note, that at this point, we have no
+    // information about anisotropically refined
+    // cells, thus use the situation of purely
+    // isotropic refinement as guess for a mixed
+    // refinemnt as well.
+    const unsigned int refine_cells  = pp - criteria_sorted.begin(),
+                       coarsen_cells = criteria_sorted.end() - 1 - qq;
+
+    if (static_cast<unsigned int>(
+          tria.n_active_cells() +
+          refine_cells * (GeometryInfo<dim>::max_children_per_cell - 1) -
+          (coarsen_cells * (GeometryInfo<dim>::max_children_per_cell - 1) /
+           GeometryInfo<dim>::max_children_per_cell)) > max_n_cells)
+      {
+        GridRefinement::refine_and_coarsen_fixed_number(tria,
+                                                        criteria,
+                                                        1. * refine_cells /
+                                                          criteria.size(),
+                                                        1. * coarsen_cells /
+                                                          criteria.size(),
+                                                        max_n_cells);
+        return;
+      }
+
+    // in some rare cases it may happen that
+    // both thresholds are the same (e.g. if
+    // there are many cells with the same
+    // error indicator). That would mean that
+    // all cells will be flagged for
+    // refinement or coarsening, but some will
+    // be flagged for both, namely those for
+    // which the indicator equals the
+    // thresholds. This is forbidden, however.
+    //
+    // In some rare cases with very few cells
+    // we also could get integer round off
+    // errors and get problems with
+    // the top and bottom fractions.
+    //
+    // In these case we arbitrarily reduce the
+    // bottom threshold by one permille below
+    // the top threshold
+    //
+    // Finally, in some cases
+    // (especially involving symmetric
+    // solutions) there are many cells
+    // with the same error indicator
+    // values. if there are many with
+    // indicator equal to the top
+    // threshold, no refinement will
+    // take place below; to avoid this
+    // case, we also lower the top
+    // threshold if it equals the
+    // largest indicator and the
+    // top_fraction!=1
+    const double max_criterion = *(criteria_sorted.begin()),
+                 min_criterion = *(criteria_sorted.end() - 1);
+
+    if ((top_threshold == max_criterion) && (top_fraction != 1))
+      top_threshold *= 0.999;
+
+    if (bottom_threshold >= top_threshold)
+      bottom_threshold = 0.999 * top_threshold;
+
+    // actually flag cells
+    if (top_threshold < max_criterion)
+      GridRefinement::refine(tria, criteria, top_threshold, refine_cells);
+
+    if (bottom_threshold > min_criterion)
+      GridRefinement::coarsen(tria, criteria, bottom_threshold);
+  }
+} // namespace
+
+
 
 template <int dim, typename Number, int spacedim>
 void
@@ -267,10 +392,10 @@ GridRefinement::refine_and_coarsen_fixed_fraction(
   const Vector<Number> &        criteria,
   const double                  top_fraction,
   const double                  bottom_fraction,
-  const unsigned int            max_n_cells)
+  const unsigned int            max_n_cells,
+  const VectorTools::NormType   norm_type)
 {
-  // correct number of cells is
-  // checked in @p{refine}
+  // correct number of cells is checked in @p{refine}
   Assert((top_fraction >= 0) && (top_fraction <= 1),
          ExcInvalidParameterValue());
   Assert((bottom_fraction >= 0) && (bottom_fraction <= 1),
@@ -280,108 +405,43 @@ GridRefinement::refine_and_coarsen_fixed_fraction(
          ExcInvalidParameterValue());
   Assert(criteria.is_non_negative(), ExcNegativeCriteria());
 
-  // let tmp be the cellwise square of the
-  // error, which is what we have to sum
-  // up and compare with
-  // @p{fraction_of_error*total_error}.
-  Vector<Number> tmp;
-  tmp                      = criteria;
-  const double total_error = tmp.l1_norm();
-
-  // sort the largest criteria to the
-  // beginning of the vector
-  std::sort(tmp.begin(), tmp.end(), std::greater<double>());
-
-  // compute thresholds
-  typename Vector<Number>::const_iterator pp = tmp.begin();
-  for (double sum = 0; (sum < top_fraction * total_error) && (pp != tmp.end());
-       ++pp)
-    sum += *pp;
-  double top_threshold = (pp != tmp.begin() ? (*pp + *(pp - 1)) / 2 : *pp);
-  typename Vector<Number>::const_iterator qq = (tmp.end() - 1);
-  for (double sum = 0;
-       (sum < bottom_fraction * total_error) && (qq != tmp.begin() - 1);
-       --qq)
-    sum += *qq;
-  double bottom_threshold = (qq != (tmp.end() - 1) ? (*qq + *(qq + 1)) / 2 : 0);
-
-  // we now have an idea how many cells we
-  // are going to refine and coarsen. we use
-  // this information to see whether we are
-  // over the limit and if so use a function
-  // that knows how to deal with this
-  // situation
-
-  // note, that at this point, we have no
-  // information about anisotropically refined
-  // cells, thus use the situation of purely
-  // isotropic refinement as guess for a mixed
-  // refinemnt as well.
-  {
-    const unsigned int refine_cells  = pp - tmp.begin(),
-                       coarsen_cells = tmp.end() - 1 - qq;
-
-    if (static_cast<unsigned int>(
-          tria.n_active_cells() +
-          refine_cells * (GeometryInfo<dim>::max_children_per_cell - 1) -
-          (coarsen_cells * (GeometryInfo<dim>::max_children_per_cell - 1) /
-           GeometryInfo<dim>::max_children_per_cell)) > max_n_cells)
-      {
-        refine_and_coarsen_fixed_number(tria,
-                                        criteria,
-                                        1. * refine_cells / criteria.size(),
-                                        1. * coarsen_cells / criteria.size(),
-                                        max_n_cells);
-        return;
-      }
-  }
-
+  switch (norm_type)
+    {
+      case VectorTools::NormType::L1_norm:
+        // evaluate norms on subsets and compare them as
+        //   c_0 + c_1 + ... < fraction * l1-norm(c)
+        refine_and_coarsen_fixed_fraction_via_l1_norm(
+          tria, criteria, top_fraction, bottom_fraction, max_n_cells);
+        break;
+
+      case VectorTools::NormType::L2_norm:
+        {
+          // we do not want to evaluate norms on subsets as:
+          //   sqrt(c_0^2 + c_1^2 + ...) < fraction * l2-norm(c)
+          // instead take the square of both sides of the equation
+          // and evaluate:
+          //   c_0^2 + c_1^2 + ... < fraction^2 * l1-norm(c.c)
+          // we adjust all parameters accordingly
+          Vector<Number> criteria_squared(criteria.size());
+          std::transform(criteria.begin(),
+                         criteria.end(),
+                         criteria_squared.begin(),
+                         [](Number c) { return c * c; });
+
+          refine_and_coarsen_fixed_fraction_via_l1_norm(tria,
+                                                        criteria_squared,
+                                                        top_fraction *
+                                                          top_fraction,
+                                                        bottom_fraction *
+                                                          bottom_fraction,
+                                                        max_n_cells);
+        }
+        break;
 
-  // in some rare cases it may happen that
-  // both thresholds are the same (e.g. if
-  // there are many cells with the same
-  // error indicator). That would mean that
-  // all cells will be flagged for
-  // refinement or coarsening, but some will
-  // be flagged for both, namely those for
-  // which the indicator equals the
-  // thresholds. This is forbidden, however.
-  //
-  // In some rare cases with very few cells
-  // we also could get integer round off
-  // errors and get problems with
-  // the top and bottom fractions.
-  //
-  // In these case we arbitrarily reduce the
-  // bottom threshold by one permille below
-  // the top threshold
-  //
-  // Finally, in some cases
-  // (especially involving symmetric
-  // solutions) there are many cells
-  // with the same error indicator
-  // values. if there are many with
-  // indicator equal to the top
-  // threshold, no refinement will
-  // take place below; to avoid this
-  // case, we also lower the top
-  // threshold if it equals the
-  // largest indicator and the
-  // top_fraction!=1
-  const auto minmax_element =
-    std::minmax_element(criteria.begin(), criteria.end());
-  if ((top_threshold == *minmax_element.second) && (top_fraction != 1))
-    top_threshold *= 0.999;
-
-  if (bottom_threshold >= top_threshold)
-    bottom_threshold = 0.999 * top_threshold;
-
-  // actually flag cells
-  if (top_threshold < *minmax_element.second)
-    refine(tria, criteria, top_threshold, pp - tmp.begin());
-
-  if (bottom_threshold > *minmax_element.first)
-    coarsen(tria, criteria, bottom_threshold);
+      default:
+        Assert(false, ExcNotImplemented());
+        break;
+    }
 }
 
 
index 1a6ff7f790898a74ff0fcd598913a207a2211152..aad99dcb06de22ec670aad7b075a6f9877618f29 100644 (file)
@@ -45,7 +45,8 @@ for (S : REAL_SCALARS; deal_II_dimension : DIMENSIONS)
                          const dealii::Vector<S> &,
                          const double,
                          const double,
-                         const unsigned int);
+                         const unsigned int,
+                         const VectorTools::NormType);
 
     template void GridRefinement::
       refine_and_coarsen_optimize<deal_II_dimension, S, deal_II_dimension>(
@@ -85,7 +86,8 @@ for (S : REAL_SCALARS; deal_II_dimension : DIMENSIONS)
         const dealii::Vector<S> &,
         const double,
         const double,
-        const unsigned int);
+        const unsigned int,
+        const VectorTools::NormType);
 
     template void GridRefinement::
       refine_and_coarsen_optimize<deal_II_dimension, S, deal_II_dimension + 1>(
diff --git a/tests/grid/refine_and_coarsen_fixed_fraction_03.cc b/tests/grid/refine_and_coarsen_fixed_fraction_03.cc
new file mode 100644 (file)
index 0000000..103d6ec
--- /dev/null
@@ -0,0 +1,107 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2020 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.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+// Verify fixed fraction algorithm with l1-norm and l2-norm
+// Equidistant indicators
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_refinement.h>
+#include <deal.II/grid/tria.h>
+
+#include <deal.II/lac/vector.h>
+
+#include <limits>
+
+#include "../tests.h"
+
+
+template <int dim>
+void
+count_flags(const Triangulation<dim> &tria)
+{
+  unsigned int n_refine = 0, n_coarsen = 0;
+  for (const auto &cell : tria.active_cell_iterators())
+    {
+      if (cell->refine_flag_set())
+        ++n_refine;
+      else if (cell->coarsen_flag_set())
+        ++n_coarsen;
+    }
+
+  deallog << "n_refine_flags: " << n_refine
+          << ", n_coarsen_flags: " << n_coarsen;
+}
+
+
+
+template <int dim>
+void
+test()
+{
+  Triangulation<dim> tria;
+  GridGenerator::hyper_cube(tria);
+  tria.refine_global(2);
+
+  Vector<float> indicator(tria.n_active_cells());
+  // assign each cell a globally unique cellid
+  for (const auto &cell : tria.active_cell_iterators())
+    {
+      const std::string  cellid = cell->id().to_string();
+      const unsigned int fine_cellid =
+        std::stoul(cellid.substr(cellid.find(':') + 1, std::string::npos));
+
+      indicator[cell->active_cell_index()] = fine_cellid + 1;
+    }
+
+  deallog << "l1-norm: ";
+  GridRefinement::refine_and_coarsen_fixed_fraction(
+    tria,
+    indicator,
+    0.3,
+    0.3,
+    std::numeric_limits<unsigned int>::max(),
+    VectorTools::NormType::L1_norm);
+  count_flags(tria);
+  deallog << std::endl;
+
+  // reset refinement flags
+  for (const auto &cell : tria.active_cell_iterators())
+    {
+      cell->clear_refine_flag();
+      cell->clear_coarsen_flag();
+    }
+
+  deallog << "l2-norm: ";
+  GridRefinement::refine_and_coarsen_fixed_fraction(
+    tria,
+    indicator,
+    0.3,
+    0.3,
+    std::numeric_limits<unsigned int>::max(),
+    VectorTools::NormType::L2_norm);
+  count_flags(tria);
+  deallog << std::endl;
+}
+
+
+
+int
+main(int argc, const char *argv[])
+{
+  initlog();
+
+  test<2>();
+}
diff --git a/tests/grid/refine_and_coarsen_fixed_fraction_03.output b/tests/grid/refine_and_coarsen_fixed_fraction_03.output
new file mode 100644 (file)
index 0000000..4a513e9
--- /dev/null
@@ -0,0 +1,3 @@
+
+DEAL::l1-norm: n_refine_flags: 3, n_coarsen_flags: 10
+DEAL::l2-norm: n_refine_flags: 1, n_coarsen_flags: 8
diff --git a/tests/grid/refine_and_coarsen_fixed_fraction_04.cc b/tests/grid/refine_and_coarsen_fixed_fraction_04.cc
new file mode 100644 (file)
index 0000000..ca8655c
--- /dev/null
@@ -0,0 +1,108 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2020 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.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+// Verify fixed fraction algorithm with l1-norm and l2-norm
+// Random indicators
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_refinement.h>
+#include <deal.II/grid/tria.h>
+
+#include <deal.II/lac/vector.h>
+
+#include <limits>
+
+#include "../tests.h"
+
+
+template <int dim>
+void
+count_flags(const Triangulation<dim> &tria)
+{
+  unsigned int n_refine = 0, n_coarsen = 0;
+  for (const auto &cell : tria.active_cell_iterators())
+    {
+      if (cell->refine_flag_set())
+        ++n_refine;
+      else if (cell->coarsen_flag_set())
+        ++n_coarsen;
+    }
+
+  deallog << "n_refine_flags: " << n_refine
+          << ", n_coarsen_flags: " << n_coarsen;
+}
+
+
+
+template <int dim>
+void
+test()
+{
+  Triangulation<dim> tria;
+  GridGenerator::hyper_cube(tria);
+  tria.refine_global(2);
+
+  Vector<float> indicator(tria.n_active_cells());
+  // assign each cell a globally unique cellid
+  for (const auto &cell : tria.active_cell_iterators())
+    {
+      const std::string  cellid = cell->id().to_string();
+      const unsigned int fine_cellid =
+        std::stoul(cellid.substr(cellid.find(':') + 1, std::string::npos));
+
+      Testing::srand(fine_cellid);
+      indicator[cell->active_cell_index()] = random_value<float>();
+    }
+
+  deallog << "l1-norm: ";
+  GridRefinement::refine_and_coarsen_fixed_fraction(
+    tria,
+    indicator,
+    0.3,
+    0.3,
+    std::numeric_limits<unsigned int>::max(),
+    VectorTools::NormType::L1_norm);
+  count_flags(tria);
+  deallog << std::endl;
+
+  // reset refinement flags
+  for (const auto &cell : tria.active_cell_iterators())
+    {
+      cell->clear_refine_flag();
+      cell->clear_coarsen_flag();
+    }
+
+  deallog << "l2-norm: ";
+  GridRefinement::refine_and_coarsen_fixed_fraction(
+    tria,
+    indicator,
+    0.3,
+    0.3,
+    std::numeric_limits<unsigned int>::max(),
+    VectorTools::NormType::L2_norm);
+  count_flags(tria);
+  deallog << std::endl;
+}
+
+
+
+int
+main(int argc, const char *argv[])
+{
+  initlog();
+
+  test<2>();
+}
diff --git a/tests/grid/refine_and_coarsen_fixed_fraction_04.output b/tests/grid/refine_and_coarsen_fixed_fraction_04.output
new file mode 100644 (file)
index 0000000..860181d
--- /dev/null
@@ -0,0 +1,3 @@
+
+DEAL::l1-norm: n_refine_flags: 4, n_coarsen_flags: 8
+DEAL::l2-norm: n_refine_flags: 1, n_coarsen_flags: 5
diff --git a/tests/mpi/refine_and_coarsen_fixed_fraction_08.cc b/tests/mpi/refine_and_coarsen_fixed_fraction_08.cc
new file mode 100644 (file)
index 0000000..7de5a22
--- /dev/null
@@ -0,0 +1,115 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2020 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.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+// Verify fixed fraction algorithm with l1-norm and l2-norm
+// Equidistant indicators
+
+#include <deal.II/base/utilities.h>
+
+#include <deal.II/distributed/grid_refinement.h>
+#include <deal.II/distributed/tria.h>
+
+#include <deal.II/grid/grid_generator.h>
+
+#include <deal.II/lac/vector.h>
+
+#include "../tests.h"
+
+
+template <int dim>
+void
+count_flags(const parallel::distributed::Triangulation<dim> &tria)
+{
+  unsigned int n_refine = 0, n_coarsen = 0;
+  for (const auto &cell : tria.active_cell_iterators())
+    if (cell->is_locally_owned())
+      {
+        if (cell->refine_flag_set())
+          ++n_refine;
+        else if (cell->coarsen_flag_set())
+          ++n_coarsen;
+      }
+
+  const unsigned int n_refine_global =
+                       Utilities::MPI::sum(n_refine, MPI_COMM_WORLD),
+                     n_coarsen_global =
+                       Utilities::MPI::sum(n_coarsen, MPI_COMM_WORLD);
+
+  deallog << "n_refine_flags: " << n_refine_global
+          << ", n_coarsen_flags: " << n_coarsen_global;
+}
+
+
+
+template <int dim>
+void
+test()
+{
+  parallel::distributed::Triangulation<dim> tria(MPI_COMM_WORLD);
+  GridGenerator::hyper_cube(tria);
+  tria.refine_global(2);
+
+  Vector<float> indicator(tria.n_active_cells());
+  // assign each cell a globally unique cellid
+  for (const auto &cell : tria.active_cell_iterators())
+    if (cell->is_locally_owned())
+      {
+        const std::string  cellid = cell->id().to_string();
+        const unsigned int fine_cellid =
+          std::stoul(cellid.substr(cellid.find(':') + 1, std::string::npos));
+
+        indicator[cell->active_cell_index()] = fine_cellid + 1;
+      }
+
+  deallog << "l1-norm: ";
+  parallel::distributed::GridRefinement::refine_and_coarsen_fixed_fraction(
+    tria, indicator, 0.3, 0.3, VectorTools::NormType::L1_norm);
+  count_flags(tria);
+  deallog << std::endl;
+
+  // reset refinement flags
+  for (const auto &cell : tria.active_cell_iterators())
+    if (cell->is_locally_owned())
+      {
+        cell->clear_refine_flag();
+        cell->clear_coarsen_flag();
+      }
+
+  deallog << "l2-norm: ";
+  parallel::distributed::GridRefinement::refine_and_coarsen_fixed_fraction(
+    tria, indicator, 0.3, 0.3, VectorTools::NormType::L2_norm);
+  count_flags(tria);
+  deallog << std::endl;
+}
+
+
+
+int
+main(int argc, char *argv[])
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
+
+  const unsigned int myid = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD);
+
+  if (myid == 0)
+    {
+      initlog();
+
+      test<2>();
+    }
+  else
+    test<2>();
+}
diff --git a/tests/mpi/refine_and_coarsen_fixed_fraction_08.mpirun=1.with_p4est=true.output b/tests/mpi/refine_and_coarsen_fixed_fraction_08.mpirun=1.with_p4est=true.output
new file mode 100644 (file)
index 0000000..c8ef4f8
--- /dev/null
@@ -0,0 +1,3 @@
+
+DEAL::l1-norm: n_refine_flags: 2, n_coarsen_flags: 9
+DEAL::l2-norm: n_refine_flags: 1, n_coarsen_flags: 7
diff --git a/tests/mpi/refine_and_coarsen_fixed_fraction_08.mpirun=4.with_p4est=true.output b/tests/mpi/refine_and_coarsen_fixed_fraction_08.mpirun=4.with_p4est=true.output
new file mode 100644 (file)
index 0000000..c8ef4f8
--- /dev/null
@@ -0,0 +1,3 @@
+
+DEAL::l1-norm: n_refine_flags: 2, n_coarsen_flags: 9
+DEAL::l2-norm: n_refine_flags: 1, n_coarsen_flags: 7
diff --git a/tests/mpi/refine_and_coarsen_fixed_fraction_09.cc b/tests/mpi/refine_and_coarsen_fixed_fraction_09.cc
new file mode 100644 (file)
index 0000000..1fedd2e
--- /dev/null
@@ -0,0 +1,116 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2020 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.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+// Verify fixed fraction algorithm with l1-norm and l2-norm
+// Random indicators
+
+#include <deal.II/base/utilities.h>
+
+#include <deal.II/distributed/grid_refinement.h>
+#include <deal.II/distributed/tria.h>
+
+#include <deal.II/grid/grid_generator.h>
+
+#include <deal.II/lac/vector.h>
+
+#include "../tests.h"
+
+
+template <int dim>
+void
+count_flags(const parallel::distributed::Triangulation<dim> &tria)
+{
+  unsigned int n_refine = 0, n_coarsen = 0;
+  for (const auto &cell : tria.active_cell_iterators())
+    if (cell->is_locally_owned())
+      {
+        if (cell->refine_flag_set())
+          ++n_refine;
+        else if (cell->coarsen_flag_set())
+          ++n_coarsen;
+      }
+
+  const unsigned int n_refine_global =
+                       Utilities::MPI::sum(n_refine, MPI_COMM_WORLD),
+                     n_coarsen_global =
+                       Utilities::MPI::sum(n_coarsen, MPI_COMM_WORLD);
+
+  deallog << "n_refine_flags: " << n_refine_global
+          << ", n_coarsen_flags: " << n_coarsen_global;
+}
+
+
+
+template <int dim>
+void
+test()
+{
+  parallel::distributed::Triangulation<dim> tria(MPI_COMM_WORLD);
+  GridGenerator::hyper_cube(tria);
+  tria.refine_global(2);
+
+  Vector<float> indicator(tria.n_active_cells());
+  // assign each cell a globally unique cellid
+  for (const auto &cell : tria.active_cell_iterators())
+    if (cell->is_locally_owned())
+      {
+        const std::string  cellid = cell->id().to_string();
+        const unsigned int fine_cellid =
+          std::stoul(cellid.substr(cellid.find(':') + 1, std::string::npos));
+
+        Testing::srand(fine_cellid);
+        indicator[cell->active_cell_index()] = random_value<float>();
+      }
+
+  deallog << "l1-norm: ";
+  parallel::distributed::GridRefinement::refine_and_coarsen_fixed_fraction(
+    tria, indicator, 0.3, 0.3, VectorTools::NormType::L1_norm);
+  count_flags(tria);
+  deallog << std::endl;
+
+  // reset refinement flags
+  for (const auto &cell : tria.active_cell_iterators())
+    if (cell->is_locally_owned())
+      {
+        cell->clear_refine_flag();
+        cell->clear_coarsen_flag();
+      }
+
+  deallog << "l2-norm: ";
+  parallel::distributed::GridRefinement::refine_and_coarsen_fixed_fraction(
+    tria, indicator, 0.3, 0.3, VectorTools::NormType::L2_norm);
+  count_flags(tria);
+  deallog << std::endl;
+}
+
+
+
+int
+main(int argc, char *argv[])
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
+
+  const unsigned int myid = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD);
+
+  if (myid == 0)
+    {
+      initlog();
+
+      test<2>();
+    }
+  else
+    test<2>();
+}
diff --git a/tests/mpi/refine_and_coarsen_fixed_fraction_09.mpirun=1.with_p4est=true.output b/tests/mpi/refine_and_coarsen_fixed_fraction_09.mpirun=1.with_p4est=true.output
new file mode 100644 (file)
index 0000000..860181d
--- /dev/null
@@ -0,0 +1,3 @@
+
+DEAL::l1-norm: n_refine_flags: 4, n_coarsen_flags: 8
+DEAL::l2-norm: n_refine_flags: 1, n_coarsen_flags: 5
diff --git a/tests/mpi/refine_and_coarsen_fixed_fraction_09.mpirun=4.with_p4est=true.output b/tests/mpi/refine_and_coarsen_fixed_fraction_09.mpirun=4.with_p4est=true.output
new file mode 100644 (file)
index 0000000..860181d
--- /dev/null
@@ -0,0 +1,3 @@
+
+DEAL::l1-norm: n_refine_flags: 4, n_coarsen_flags: 8
+DEAL::l2-norm: n_refine_flags: 1, n_coarsen_flags: 5

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.