From 594dc87afdb7cf442510c5348ac7909b4a0883a2 Mon Sep 17 00:00:00 2001 From: Lei Qiao Date: Sat, 30 May 2015 10:58:33 -0500 Subject: [PATCH] separate refinement and coarsening fraction adjustment into an individual function --- include/deal.II/grid/grid_refinement.h | 35 ++++++++++ source/grid/grid_refinement.cc | 91 ++++++++++++++++---------- source/grid/grid_refinement.inst.in | 13 ++++ 3 files changed, 105 insertions(+), 34 deletions(-) diff --git a/include/deal.II/grid/grid_refinement.h b/include/deal.II/grid/grid_refinement.h index cb49c85184..f0127282ff 100644 --- a/include/deal.II/grid/grid_refinement.h +++ b/include/deal.II/grid/grid_refinement.h @@ -47,6 +47,41 @@ template class Vector; */ namespace GridRefinement { + /** + * Return a pair of double values of which the first is adjusted refinement + * fraction of cells and the second is adjusted coarsening fraction of cells. + * + * + * @arg @p current_n_cells is current cell number. + * + * @arg @p max_n_cells is the maximal number of cells. If current cell number + * @p current_n_cells is already exceeded maximal cell number @p max_n_cells, + * refinement fraction of cells will be set to zero and coarsening fraction + * of cells will be adjusted to reduce cell number to @ max_n_cells. + * If cell number is going to be exceeded only upon refinement, then refinement + * and coarsening fractions are going to be adjusted with a same ratio 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. + * + * @arg @p top_fraction is the requested fraction of cells to be refined. + * + * @arg @p bottom_fraction is the requested fraction of cells to be coarsened. + * + * @note Usually you do not need to call this function explicitly. Pass + * @p max_n_cells to function refine_and_coarsen_fixed_number() or function + * refine_and_coarsen_fixed_fraction() and they will call this function if + * necessary. + */ + template + std::pair + adjust_refine_and_coarsen_number_fraction (const unsigned int current_n_cells, + const unsigned int max_n_cells, + const double top_fraction_of_cells, + const double bottom_fraction_of_cells); + /** * This function provides a refinement strategy with predictable growth of * the mesh. diff --git a/source/grid/grid_refinement.cc b/source/grid/grid_refinement.cc index 3282956d50..8b05b2455f 100644 --- a/source/grid/grid_refinement.cc +++ b/source/grid/grid_refinement.cc @@ -277,30 +277,30 @@ void GridRefinement::coarsen (Triangulation &tria, cell->set_coarsen_flag(); } - - -template -void -GridRefinement::refine_and_coarsen_fixed_number (Triangulation &tria, - const Vector &criteria, - const double top_fraction, - const double bottom_fraction, - const unsigned int max_n_cells) +template +std::pair +GridRefinement::adjust_refine_and_coarsen_number_fraction (const unsigned int current_n_cells, + const unsigned int max_n_cells, + const double top_fraction, + const double bottom_fraction) { - // correct number of cells is - // checked in @p{refine} - Assert ((top_fraction>=0) && (top_fraction<=1), ExcInvalidParameterValue()); - Assert ((bottom_fraction>=0) && (bottom_fraction<=1), ExcInvalidParameterValue()); + Assert (top_fraction>=0, ExcInvalidParameterValue()); + Assert (top_fraction<=1, ExcInvalidParameterValue()); + Assert (bottom_fraction>=0, ExcInvalidParameterValue()); + Assert (bottom_fraction<=1, ExcInvalidParameterValue()); Assert (top_fraction+bottom_fraction <= 1, ExcInvalidParameterValue()); - Assert (criteria.is_non_negative (), ExcNegativeCriteria()); - int refine_cells = static_cast(top_fraction*criteria.size()); - int coarsen_cells = static_cast(bottom_fraction*criteria.size()); + double refine_cells = current_n_cells * top_fraction; + double coarsen_cells = current_n_cells * bottom_fraction; + + const double cell_increase_on_refine = GeometryInfo::max_children_per_cell - 1.0; + const double cell_decrease_on_coarsen = 1.0 - 1.0/GeometryInfo::max_children_per_cell; + std::pair adjusted_fractions(top_fraction, bottom_fraction); // first we have to see whether we // currently already exceed the target // number of cells - if (tria.n_active_cells() >= max_n_cells) + if (current_n_cells >= max_n_cells) { // if yes, then we need to stop // refining cells and instead try to @@ -312,10 +312,10 @@ GridRefinement::refine_and_coarsen_fixed_number (Triangulation &tr // anisotropically, assume isotropic // refinement here, though that may // result in a worse approximation - refine_cells = 0; - coarsen_cells = (tria.n_active_cells() - max_n_cells) * - GeometryInfo::max_children_per_cell / - (GeometryInfo::max_children_per_cell - 1); + adjusted_fractions.first = 0; + coarsen_cells = (current_n_cells - max_n_cells) / + cell_decrease_on_coarsen; + adjusted_fractions.second = std::min(coarsen_cells/current_n_cells, 1.0); } // otherwise, see if we would exceed the // maximum desired number of cells with the @@ -331,11 +331,9 @@ GridRefinement::refine_and_coarsen_fixed_number (Triangulation &tr // refined cells. we take this as an // approximation of a mixed refinement. else if (static_cast - (tria.n_active_cells() - + refine_cells * (GeometryInfo::max_children_per_cell - 1) - - (coarsen_cells * - (GeometryInfo::max_children_per_cell - 1) / - GeometryInfo::max_children_per_cell)) + (current_n_cells + + refine_cells * cell_increase_on_refine + - coarsen_cells / cell_decrease_on_coarsen) > max_n_cells) { @@ -346,19 +344,44 @@ GridRefinement::refine_and_coarsen_fixed_number (Triangulation &tr // fractions and the resulting number // of cells to be equal to // max_n_cells. this leads to the - // following equation for lambda + // following equation for alpha const double alpha = 1. * - (max_n_cells - tria.n_active_cells()) + (max_n_cells - current_n_cells) / - (refine_cells * (GeometryInfo::max_children_per_cell - 1) - - (coarsen_cells * - (GeometryInfo::max_children_per_cell - 1) / - GeometryInfo::max_children_per_cell)); - refine_cells = static_cast (refine_cells * alpha); - coarsen_cells = static_cast (coarsen_cells * alpha); + (refine_cells * cell_increase_on_refine + - coarsen_cells / cell_decrease_on_coarsen); + + adjusted_fractions.first = alpha * top_fraction; + adjusted_fractions.second = alpha * bottom_fraction; } + return (adjusted_fractions); +} + +template +void +GridRefinement::refine_and_coarsen_fixed_number (Triangulation &tria, + const Vector &criteria, + const double top_fraction, + const double bottom_fraction, + const unsigned int max_n_cells) +{ + // correct number of cells is + // checked in @p{refine} + Assert ((top_fraction>=0) && (top_fraction<=1), ExcInvalidParameterValue()); + Assert ((bottom_fraction>=0) && (bottom_fraction<=1), ExcInvalidParameterValue()); + Assert (top_fraction+bottom_fraction <= 1, ExcInvalidParameterValue()); + Assert (criteria.is_non_negative (), ExcNegativeCriteria()); + + const std::pair adjusted_fractions = + adjust_refine_and_coarsen_number_fraction (criteria.size(), + max_n_cells, + top_fraction, + bottom_fraction); + + const int refine_cells = static_cast(adjusted_fractions.first * criteria.size()); + const int coarsen_cells = static_cast(adjusted_fractions.second * criteria.size()); if (refine_cells || coarsen_cells) { diff --git a/source/grid/grid_refinement.inst.in b/source/grid/grid_refinement.inst.in index bd90fb65b2..ad4123bd09 100644 --- a/source/grid/grid_refinement.inst.in +++ b/source/grid/grid_refinement.inst.in @@ -110,3 +110,16 @@ for (S : REAL_SCALARS; deal_II_dimension : DIMENSIONS) const unsigned int); #endif } + +for (deal_II_dimension : DIMENSIONS) +{ + template + std::pair + GridRefinement:: + adjust_refine_and_coarsen_number_fraction + (const unsigned int, + const unsigned int, + const double, + const double); +} + -- 2.39.5