]> https://gitweb.dealii.org/ - dealii.git/commitdiff
add maximal cell number limit to parallel version of refine_and_coarsen_fixed_number
authorLei Qiao <qiaol618@gmail.com>
Sat, 30 May 2015 16:39:50 +0000 (11:39 -0500)
committerLei Qiao <qiaol618@gmail.com>
Wed, 3 Jun 2015 13:11:23 +0000 (08:11 -0500)
include/deal.II/distributed/grid_refinement.h
source/distributed/grid_refinement.cc
source/distributed/grid_refinement.inst.in

index 42e82af1a12b4424da9e83cebe05d965891b53f8..862ba4af84b7982ae972a411440db6a58e8ef02a 100644 (file)
@@ -70,8 +70,9 @@ namespace parallel
       refine_and_coarsen_fixed_number (
         parallel::distributed::Triangulation<dim,spacedim> &tria,
         const Vector                &criteria,
-        const double                top_fraction_of_cells,
-        const double                bottom_fraction_of_cells);
+        const double                 top_fraction_of_cells,
+        const double                 bottom_fraction_of_cells,
+        const unsigned int           max_n_cells = std::numeric_limits<unsigned int>::max());
 
       /**
        * Like dealii::GridRefinement::refine_and_coarsen_fixed_fraction, but
index 731f505a38175e9b6d5f3cf69257a6825621c1e2..1ee88da7df067cd78ec406d54d96e3f471f97d8a 100644 (file)
@@ -543,8 +543,9 @@ namespace parallel
       refine_and_coarsen_fixed_number (
         parallel::distributed::Triangulation<dim,spacedim> &tria,
         const Vector                &criteria,
-        const double                top_fraction_of_cells,
-        const double                bottom_fraction_of_cells)
+        const double                 top_fraction_of_cells,
+        const double                 bottom_fraction_of_cells,
+        const unsigned int           max_n_cells)
       {
         Assert ((top_fraction_of_cells>=0) && (top_fraction_of_cells<=1),
                 dealii::GridRefinement::ExcInvalidParameterValue());
@@ -555,6 +556,13 @@ namespace parallel
         Assert (criteria.is_non_negative (),
                 dealii::GridRefinement::ExcNegativeCriteria());
 
+        const std::pair<double, double> adjusted_fractions =
+          dealii::GridRefinement::adjust_refine_and_coarsen_number_fraction<dim> (
+            tria.n_global_active_cells(),
+            max_n_cells,
+            top_fraction_of_cells,
+            bottom_fraction_of_cells);
+
         // first extract from the
         // vector of indicators the
         // ones that correspond to
@@ -590,7 +598,7 @@ namespace parallel
                 master_compute_threshold (locally_owned_indicators,
                                           global_min_and_max,
                                           static_cast<unsigned int>
-                                          (top_fraction_of_cells *
+                                          (adjusted_fractions.first *
                                            tria.n_global_active_cells()),
                                           mpi_communicator);
 
@@ -600,14 +608,14 @@ namespace parallel
             // use a threshold lower
             // than the smallest
             // value we have locally
-            if (bottom_fraction_of_cells > 0)
+            if (adjusted_fractions.second > 0)
               bottom_threshold
                 =
                   RefineAndCoarsenFixedNumber::
                   master_compute_threshold (locally_owned_indicators,
                                             global_min_and_max,
                                             static_cast<unsigned int>
-                                            ((1-bottom_fraction_of_cells) *
+                                            ((1-adjusted_fractions.second) *
                                              tria.n_global_active_cells()),
                                             mpi_communicator);
             else
@@ -629,7 +637,7 @@ namespace parallel
             // compute bottom
             // threshold only if
             // necessary
-            if (bottom_fraction_of_cells > 0)
+            if (adjusted_fractions.second > 0)
               bottom_threshold
                 =
                   RefineAndCoarsenFixedNumber::
index c9527a2b490be57c3ff90f52797ed8c3ac00ef32..f2fddc5b131b464c96170abbbb537c184ce96a9a 100644 (file)
@@ -31,7 +31,8 @@ namespace parallel
       (parallel::distributed::Triangulation<deal_II_dimension> &,
        const dealii::Vector<S> &,
        const double,
-       const double);
+       const double,
+       const unsigned int);
 
       template
       void
@@ -60,7 +61,8 @@ namespace parallel
       (parallel::distributed::Triangulation<deal_II_dimension-1,deal_II_dimension> &,
        const dealii::Vector<S> &,
        const double,
-       const double);
+       const double,
+       const unsigned int);
 
       template
       void

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.