From cf0817acd9ddb91e630a7c8d5d4f799c0cac8b80 Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Thu, 8 Sep 2016 03:15:22 -0600 Subject: [PATCH] Update documentation for grid refinement functions. --- include/deal.II/distributed/grid_refinement.h | 30 +++++----- include/deal.II/grid/grid_refinement.h | 58 +++++++++---------- 2 files changed, 45 insertions(+), 43 deletions(-) diff --git a/include/deal.II/distributed/grid_refinement.h b/include/deal.II/distributed/grid_refinement.h index 1d09221300..e090a18c48 100644 --- a/include/deal.II/distributed/grid_refinement.h +++ b/include/deal.II/distributed/grid_refinement.h @@ -1,6 +1,6 @@ // --------------------------------------------------------------------- // -// Copyright (C) 2009 - 2015 by the deal.II authors +// Copyright (C) 2009 - 2016 by the deal.II authors // // This file is part of the deal.II library. // @@ -43,18 +43,20 @@ namespace parallel { /** * Like dealii::GridRefinement::refine_and_coarsen_fixed_number, but for - * parallel distributed triangulation. + * parallel distributed triangulations. * * The vector of criteria needs to be a vector of refinement criteria - * for all cells active on the current triangulation, i.e. - * tria.n_active_cells() (and not - * tria.n_locally_owned_active_cells()). However, the + * for all cells active on the current triangulation, i.e., + * it needs to be of length tria.n_active_cells() (and not + * tria.n_locally_owned_active_cells()). In other words, + * the vector needs to include entries for ghost and artificial + * cells. However, the current * function will only look at the indicators that correspond to those * cells that are actually locally owned, and ignore the indicators for * all other cells. The function will then coordinate among all - * processors that store part of the triangulation so that at the end @p - * top_fraction_of_cells are refined, where the fraction is enforced as - * a fraction of Triangulation::n_global_active_cells, not + * processors that store part of the triangulation so that at the end + * a fraction @p top_fraction_of_cells of all Triangulation::n_global_active_cells() + * active cells are refined, rather than a fraction of the * Triangulation::n_locally_active_cells on each processor individually. * In other words, it may be that on some processors, no cells are * refined at all. @@ -72,18 +74,20 @@ namespace parallel /** * Like dealii::GridRefinement::refine_and_coarsen_fixed_fraction, but - * for parallel distributed triangulation. + * for parallel distributed triangulations. * * The vector of criteria needs to be a vector of refinement criteria - * for all cells active on the current triangulation, - * tria.n_active_cells() (and not - * tria.n_locally_owned_active_cells()). However, the + * for all cells active on the current triangulation, i.e., + * it needs to be of length tria.n_active_cells() (and not + * tria.n_locally_owned_active_cells()). In other words, + * the vector needs to include entries for ghost and artificial + * cells. However, the current * function will only look at the indicators that correspond to those * cells that are actually locally owned, and ignore the indicators for * all other cells. The function will then coordinate among all * processors that store part of the triangulation so that at the end * the smallest fraction of Triangulation::n_global_active_cells (not - * Triangulation::n_locally_active_cells on each processor individually) + * Triangulation::n_locally_owned_active_cells() on each processor individually) * is refined that together make up a total of @p top_fraction_of_error * of the total error. In other words, it may be that on some * processors, no cells are refined at all. diff --git a/include/deal.II/grid/grid_refinement.h b/include/deal.II/grid/grid_refinement.h index 1f60f38c82..a55e9fc77e 100644 --- a/include/deal.II/grid/grid_refinement.h +++ b/include/deal.II/grid/grid_refinement.h @@ -64,11 +64,11 @@ namespace GridRefinement * default value of this argument is to impose no limit on the number of * cells. * - * @param[in] top_fraction_of_cells The requested fraction of cells to be - * refined. + * @param[in] top_fraction_of_cells The requested fraction of active + * cells to be refined. * - * @param[in] bottom_fraction_of_cells The requested fraction of cells to be - * coarsened. + * @param[in] bottom_fraction_of_cells The requested fraction of + * active 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 @@ -83,8 +83,8 @@ namespace GridRefinement const double bottom_fraction_of_cells); /** - * This function provides a refinement strategy with predictable growth of - * the mesh. + * This function provides a refinement strategy with predictable growth in + * the size of the mesh by refining a given fraction of all cells. * * The function takes a vector of refinement @p criteria and two values * between zero and one denoting the fractions of cells to be refined and @@ -96,31 +96,29 @@ namespace GridRefinement * *
  • Sort the cells according to descending values of @p criteria. * - *
  • Set the refinement threshold to be the criterion belonging to the - * cell at position @p top_fraction_of_cells times - * Triangulation::n_active_cells(). + *
  • Mark the @p top_fraction_of_cells times + * Triangulation::n_active_cells() active cells with the largest + * refinement criteria for refinement. * - *
  • Set the coarsening threshold accordingly using the cell @p - * bottom_fraction_of_cells times Triangulation::n_active_cells() from the - * end of the sorted list. - * - *
  • Use these two thresholds in calls to refine() and coarsen(), - * respectively. + *
  • Mark the @p bottom_fraction_of_cells times + * Triangulation::n_active_cells() active cells with the smallest + * refinement criteria for coarsening. * * * * As an example, with no coarsening, setting @p top_fraction_of_cells to * 1/3 will result in approximately doubling the number of cells in two - * dimensions. The same effect in three dimensions is achieved by refining - * 1/7th of the cells. These values are good initial guesses, but should be - * adjusted depending on the singularity of approximated function. - * - * The sorting of criteria is not done actually, since we only need the - * threshold values in order to call refine() and coarsen(). The order of - * cells with higher and of those with lower criteria is irrelevant. Getting - * this value is accomplished by the @p nth_element function of the - * C++ standard library, which takes only linear time in the number - * of elements, rather than N log N for sorting all values. + * dimensions. That is because each of these 1/3 of cells will be replaced by + * its four children, resulting in $4\times \frac 13 N$ cells, whereas the + * remaining 2/3 of cells remains untouched -- thus yielding a total of + * $4\times \frac 13 N + \frac 23 N = 2N$ cells. + * The same effect in three dimensions is achieved by refining + * 1/7th of the cells. These values are therefore frequently used because + * they ensure that the cost of computations on subsequent meshes become + * expensive sufficiently quickly that the fraction of time spent on + * the coarse meshes is not too large. On the other hand, the fractions + * are small enough that mesh adaptation does not refine too many cells + * in each step. * * @note This function only sets the coarsening and refinement flags. The * mesh is not changed until you call @@ -160,11 +158,11 @@ namespace GridRefinement * This function provides a refinement strategy controlling the reduction of * the error estimate. * - * Also known as the bulk criterion, this function computes the - * thresholds for refinement and coarsening such that the @p criteria of - * cells getting flagged for refinement make up for a certain fraction of - * the total error. We explain its operation for refinement, coarsening - * works analogously. + * Also known as the bulk criterion or Dörfler marking, + * this function computes the thresholds for refinement and coarsening + * such that the @p criteria of cells getting flagged for refinement make + * up for a certain fraction of the total error. We explain its operation + * for refinement, coarsening works analogously. * * Let cK be the criterion of cell K. Then the * total error estimate is computed by the formula -- 2.39.5