As the interpolation while coarsening is much more complicated to organize (see further documentation below) than interpolation while pure refinement, SolutionTransfer
offers two possible usages.
SolutionTransfer
as follows * SolutionTransfer<dim, double> soltrans(*dof_handler); * // flag some cells for refinement, e.g. * GridRefinement::refine_and_coarsen_fixed_fraction( * *tria, error_indicators, 0.3, 0); * // prepare the triangulation * // for refinement, * tria->prepare_coarsening_and_refinement(); * // tell the SolutionTransfer object * // that we intend to do pure refinement, * soltrans.prepare_for_pure_refinement(); * // actually execute the refinement, * tria->execute_coarsening_and_refinement(); * // and redistribute dofs. * dof_handler->distribute_dofs (fe); *
Then to proceed do
* // take a copy of the solution vector * Vector<double> solution_old(solution); * // resize solution vector to the correct * // size, as the @p refine_interpolate * // function requires the vectors to be * // of right sizes * solution.reinit(dof_handler->n_dofs()); * // and finally interpolate * soltrans.refine_interpolate(solution_old, solution); *
Although the refine_interpolate
functions are allowed to be called multiple times, e.g. for interpolating several solution vectors, there is following possibility of interpolating several functions simultaneously.
* vector<Vector<double> > solutions_old(n_vectors, Vector<double> (n)); * ... * vector<Vector<double> > solutions(n_vectors, Vector<double> (n)); * soltrans.refine_interpolate(solutions_old, solutions); *
SolutionTransfer
as follows * SolutionTransfer<dim, double> soltrans(*dof_handler); * // flag some cells for refinement * // and coarsening, e.g. * GridRefinement::refine_and_coarsen_fixed_fraction( * *tria, error_indicators, 0.3, 0.05); * // prepare the triangulation, * tria->prepare_coarsening_and_refinement(); * // prepare the SolutionTransfer object * // for coarsening and refinement and give * // the solution vector that we intend to * // interpolate later, * soltrans.prepare_for_coarsening_and_refinement(solution); * // actually execute the refinement, * tria->execute_coarsening_and_refinement (); * // redistribute dofs, * dof_handler->distribute_dofs (fe); * // and interpolate the solution * Vector<double> interpolate_solution(dof_handler->n_dofs()); * soltrans.interpolate(solution, interpolated_solution); *
Multiple calls to the function interpolate (const Vector<number> &in, Vector<number> &out)
are NOT allowed. Interpolating several functions can be performed in one step by using void interpolate (const vector<Vector<number> >&all_in, vector<Vector<number> >&all_out) const
, and using the respective prepare_for_coarsening_and_refinement
function taking several vectors as input before actually refining and coarsening the triangulation (see there).
For deleting all stored data in SolutionTransfer
and reinitializing it use the clear()
function.
The template argument number
denotes the data type of the vectors you want to transfer.
In the function refine_interpolate(in,out)
and on each cell where the pointer is set (i.e. the cells that were active in the original grid) we can now access the local values of the solution vector in
on that cell by using the stored DoF indices. These local values are interpolated and set into the vector out
that is at the end the discrete function in
interpolated on the refined mesh.
The refine_interpolate(in,out)
function can be called multiple times for arbitrary many discrete functions (solution vectors) on the original grid.
prepare_coarsening_and_refinement
the coarsen flags of either all or none of the children of a (father-)cell are set. While coarsening (Triangulationexecute_coarsening_and_refinement
) the cells that are not needed any more will be deleted from the Triangulation.
For the interpolation from the (to be coarsenend) children to their father the children cells are needed. Hence this interpolation and the storing of the interpolated values of each of the discrete functions that we want to interpolate needs to take place before these children cells are coarsened (and deleted!!). Again a pointers for the relevant cells is set to point to these values (see below). Additionally the DoF indices of the cells that will not be coarsened need to be stored according to the solution transfer while pure refinement (cf there). All this is performed by prepare_for_coarsening_and_refinement(all_in)
where the vector<Vector<number> >vector all_in
includes all discrete functions to be interpolated onto the new grid.
As we need two different kinds of pointers (vector<unsigned int> *
for the Dof indices and vector<Vector<number> > *
for the interpolated DoF values) we use the Pointerstruct
that includes both of these pointers and the pointer for each cell points to these Pointerstructs
. On each cell only one of the two different pointers is used at one time hence we could use a void * pointer
as vector<unsigned int> *
at one time and as vector<Vector<number> > *
at the other but using this Pointerstruct
in between makes the use of these pointers more safe and gives better possibility to expand their usage.
In interpolate(all_in, all_out)
the refined cells are treated according to the solution transfer while pure refinement. Additionally, on each cell that is coarsened (hence previously was a father cell), the values of the discrete functions in all_out
are set to the stored local interpolated values that are accessible due to the 'vector<Vector<number> > *' pointer in Pointerstruct
that is pointed to by the pointer of that cell. It is clear that interpolate(all_in, all_out)
only can be called with the vector<Vector<number> > all_in
that previously was the parameter of the prepare_for_coarsening_and_refinement(all_in)
function. Hence interpolate(all_in, all_out)
can (in contrast to refine_interpolate(in, out)
) only be called once.
enum SolutionTransfer::PreparationState [private] |
Declaration of PreparationState
that denotes the three possible states of the SolutionTransfer:
being prepared for 'pure refinement', prepared for 'coarsening and refinement' or not prepared.
SolutionTransfer< dim, VectorType, DH >::SolutionTransfer | ( | const DH & | dof | ) |
Constructor, takes the current DoFHandler as argument.
SolutionTransfer< dim, VectorType, DH >::~SolutionTransfer | ( | ) |
Destructor
void SolutionTransfer< dim, VectorType, DH >::clear | ( | ) |
Reinit this class to the state that it has directly after calling the Constructor
void SolutionTransfer< dim, VectorType, DH >::prepare_for_pure_refinement | ( | ) |
Prepares the SolutionTransfer
for pure refinement. It stores the dof indices of each cell. After calling this function only calling the refine_interpolate
functions is allowed.
void SolutionTransfer< dim, VectorType, DH >::prepare_for_coarsening_and_refinement | ( | const std::vector< VectorType > & | all_in | ) |
Prepares the SolutionTransfer
for coarsening and refinement. It stores the dof indices of each cell and stores the dof values of the vectors in all_in
in each cell that'll be coarsened. all_in
includes all vectors that are to be interpolated onto the new (refined and/or coarsenend) grid.
void SolutionTransfer< dim, VectorType, DH >::prepare_for_coarsening_and_refinement | ( | const VectorType & | in | ) |
Same as previous function but for only one discrete function to be interpolated.
void SolutionTransfer< dim, VectorType, DH >::refine_interpolate | ( | const VectorType & | in, | |
VectorType & | out | |||
) | const |
This function interpolates the discrete function in
, which is a vector on the grid before the refinement, to the function out
which then is a vector on the refined grid. It assumes the vectors having the right sizes (i.e. in.size()==n_dofs_old
, out.size()==n_dofs_refined
)
Calling this function is allowed only if prepare_for_pure_refinement
is called and the refinement is executed before. Multiple calling of this function is allowed. e.g. for interpolating several functions.
void SolutionTransfer< dim, VectorType, DH >::interpolate | ( | const std::vector< VectorType > & | all_in, | |
std::vector< VectorType > & | all_out | |||
) | const |
This function interpolates the discrete functions that are stored in all_out
onto the refined and/or coarsenend grid. It assumes the vectors in all_in
denote the same vectors as in all_in
as parameter of prepare_for_refinement_and_coarsening(all_in)
. However, there is no way of verifying this internally, so be careful here.
Calling this function is allowed only if first Triangulationprepare_coarsening_and_refinement
, second SolutionTransfer::prepare_for_coarsening_and_refinement
, an then third Triangulationexecute_coarsening_and_refinement
are called before. Multiple calling of this function is NOT allowed. Interpolating several functions can be performed in one step.
The number of output vectors is assumed to be the same as the number of input vectors. Also, the sizes of the output vectors are assumed to be of the right size (n_dofs_refined
). Otherwise an assertion will be thrown.
void SolutionTransfer< dim, VectorType, DH >::interpolate | ( | const VectorType & | in, | |
VectorType & | out | |||
) | const |
Same as the previous function. It interpolates only one function. It assumes the vectors having the right sizes (i.e. in.size()==n_dofs_old
, out.size()==n_dofs_refined
)
Multiple calling of this function is NOT allowed. Interpolating several functions can be performed in one step by using interpolate (all_in, all_out)
unsigned int SolutionTransfer< dim, VectorType, DH >::memory_consumption | ( | ) | const |
Determine an estimate for the memory consumption (in bytes) of this object.
SmartPointer<const DH> SolutionTransfer< dim, VectorType, DH >::dof_handler [private] |
Pointer to the degree of freedom handler to work with.
unsigned int SolutionTransfer< dim, VectorType, DH >::n_dofs_old [private] |
Stores the number of DoFs before the refinement and/or coarsening.
PreparationState SolutionTransfer< dim, VectorType, DH >::prepared_for [private] |
Definition of the respective variable.
std::vector<std::vector<unsigned int> > SolutionTransfer< dim, VectorType, DH >::indices_on_cell [private] |
Is used for prepare_for_refining
(of course also for repare_for_refining_and_coarsening
) and stores all dof indices of the cells that'll be refined
std::map<std::pair<unsigned int, unsigned int>, Pointerstruct> SolutionTransfer< dim, VectorType, DH >::cell_map [private] |
Map mapping from level and index of cell to the Pointerstructs
(cf. there). This map makes it possible to keep all the information needed to transfer the solution inside this object rather than using user pointers of the Triangulation for this purpose.
std::vector<std::vector<Vector<typename VectorType::value_type> > > SolutionTransfer< dim, VectorType, DH >::dof_values_on_cell [private] |
Is used for prepare_for_refining_and_coarsening
The interpolated dof values of all cells that'll be coarsened will be stored in this vector.