Reference documentation for deal.II version GIT relicensing-214-g6e74dec06b 2024-03-27 18:10:01+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
Functions
parallel::distributed::GridRefinement Namespace Reference

Functions

template<int dim, typename Number , int spacedim>
void refine_and_coarsen_fixed_number (::Triangulation< dim, spacedim > &tria, const ::Vector< Number > &criteria, const double top_fraction_of_cells, const double bottom_fraction_of_cells, const types::global_cell_index max_n_cells=std::numeric_limits< types::global_cell_index >::max())
 
template<int dim, typename Number , int spacedim>
void refine_and_coarsen_fixed_fraction (::Triangulation< dim, spacedim > &tria, const ::Vector< Number > &criteria, const double top_fraction_of_error, const double bottom_fraction_of_error, const VectorTools::NormType norm_type=VectorTools::L1_norm)
 

Detailed Description

This namespace provides a collection of functions that aid in refinement and coarsening of triangulations. Despite the name of the namespace, the functions do not actually refine the triangulation, but only mark cells for refinement or coarsening. In other words, they perform the "mark" part of the typical "solve-estimate-mark-refine" cycle of the adaptive finite element loop.

In contrast to the functions in namespace GridRefinement, the functions in the current namespace are intended for parallel computations, i.e., computations, e.g., on objects of type parallel::distributed::Triangulation.

Function Documentation

◆ refine_and_coarsen_fixed_number()

template<int dim, typename Number , int spacedim>
void parallel::distributed::GridRefinement::refine_and_coarsen_fixed_number ( ::Triangulation< dim, spacedim > &  tria,
const ::Vector< Number > &  criteria,
const double  top_fraction_of_cells,
const double  bottom_fraction_of_cells,
const types::global_cell_index  max_n_cells = std::numeric_limits<types::global_cell_index>::max() 
)

Like 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., 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 a fraction 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.

The same is true for the fraction of cells that is coarsened.

Parameters
[in,out]triaThe triangulation whose cells this function is supposed to mark for coarsening and refinement.
[in]criteriaThe refinement criterion for each mesh cell active on the current triangulation. Entries may not be negative.
[in]top_fraction_of_cellsThe 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.
[in]bottom_fraction_of_cellsThe fraction of cells to be coarsened. If this number is zero, no cells will be coarsened.
[in]max_n_cellsThis 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.

Definition at line 503 of file grid_refinement.cc.

◆ refine_and_coarsen_fixed_fraction()

template<int dim, typename Number , int spacedim>
void parallel::distributed::GridRefinement::refine_and_coarsen_fixed_fraction ( ::Triangulation< dim, spacedim > &  tria,
const ::Vector< Number > &  criteria,
const double  top_fraction_of_error,
const double  bottom_fraction_of_error,
const VectorTools::NormType  norm_type = VectorTools::L1_norm 
)

Like GridRefinement::refine_and_coarsen_fixed_fraction, but 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., 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_owned_active_cells() on each processor individually) is refined that together make up a total of top_fraction_of_error of the total error. In other words, it may be that on some processors, no cells are refined at all.

The same is true for the fraction of cells that is coarsened.

Parameters
[in,out]triaThe triangulation whose cells this function is supposed to mark for coarsening and refinement.
[in]criteriaThe refinement criterion computed on each mesh cell active on the current triangulation. Entries may not be negative.
[in]top_fraction_of_errorThe 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.
[in]bottom_fraction_of_errorThe fraction of the estimate coarsened. If this number is zero, no cells will be coarsened.
[in]norm_typeTo 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::L1_norm and VectorTools::L2_norm are currently supported.

Definition at line 576 of file grid_refinement.cc.