Public Member Functions | |
| SolutionTransfer (const DH &dof) | |
| ~SolutionTransfer () | |
| void | prepare_for_coarsening_and_refinement (const std::vector< const VECTOR * > &all_in) |
| void | prepare_for_coarsening_and_refinement (const VECTOR &in) |
| void | interpolate (VECTOR &out) |
| unsigned int | get_data_size () const |
| void | prepare_serialization (const VECTOR &in) |
| void | prepare_serialization (const std::vector< const VECTOR * > &all_in) |
| void | deserialize (VECTOR &in) |
| void | deserialize (std::vector< VECTOR * > &all_in) |
Private Member Functions | |
| void | pack_callback (const typename Triangulation< dim, dim >::cell_iterator &cell, const typename Triangulation< dim, dim >::CellStatus status, void *data) |
| void | unpack_callback (const typename Triangulation< dim, dim >::cell_iterator &cell, const typename Triangulation< dim, dim >::CellStatus status, const void *data, std::vector< VECTOR * > &all_out) |
Private Attributes | |
| SmartPointer< const DH, SolutionTransfer< dim, VECTOR, DH > > | dof_handler |
| std::vector< const VECTOR * > | input_vectors |
| unsigned int | offset |
Transfers a discrete FE function (like a solution vector) by interpolation while refining and/or coarsening a distributed grid and handles the necessary communication.
In a parallel computation PETSc or Trilinos vector may contain ghost elements or not. For reading in information with prepare_for_coarsening_and_refinement() or prepare_serialization() you need to supply vectors with ghost elements, so that all locally_active elements can be read. On the other hand, ghosted vectors are generally not writable, so for calls to interpolate() or deserialize() you need to supply distributed vectors without ghost elements.
Here VECTOR is your favorite vector type, e.g. PETScWrappers::MPI::Vector, TrilinosWrappers::MPI::Vector, or corresponding blockvectors.
SolutionTransfer<dim, VECTOR> 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 interpolated_solution;
//create VECTOR in the right size here
soltrans.interpolate(interpolated_solution);
This class can be used to serialize and later deserialize a distributed mesh with solution vectors to a file. If you use more than one DoFHandler and therefore more than one SolutionTransfer object, they need to be serialized and deserialized in the same order.
If vector has the locally relevant DoFs, serialization works as follows:
parallel::distributed::SolutionTransfer<dim,VECTOR> sol_trans(dof_handler); sol_trans.prepare_serialization (vector); triangulation.save(filename);
For deserialization the vector needs to be a distributed vector (without ghost elements):
//[create coarse mesh...] triangulation.load(filename); parallel::distributed::SolutionTransfer<dim,VECTOR> sol_trans(dof_handler); sol_trans.deserialize (distributed_vector);
Definition at line 100 of file solution_transfer.h.
| parallel::distributed::SolutionTransfer< dim, VECTOR, DH >::SolutionTransfer | ( | const DH & | dof | ) |
Constructor, takes the current DoFHandler as argument.
| parallel::distributed::SolutionTransfer< dim, VECTOR, DH >::~SolutionTransfer | ( | ) |
Destructor.
| void parallel::distributed::SolutionTransfer< dim, VECTOR, DH >::prepare_for_coarsening_and_refinement | ( | const std::vector< const VECTOR * > & | 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 parallel::distributed::SolutionTransfer< dim, VECTOR, DH >::prepare_for_coarsening_and_refinement | ( | const VECTOR & | in | ) |
Same as previous function but for only one discrete function to be interpolated.
| void parallel::distributed::SolutionTransfer< dim, VECTOR, DH >::interpolate | ( | VECTOR & | out | ) |
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 parallel::distributed::SolutionTransfer< dim, VECTOR, DH >::get_data_size | ( | ) | const |
Return the size in bytes that need to be stored per cell.
| void parallel::distributed::SolutionTransfer< dim, VECTOR, DH >::prepare_serialization | ( | const VECTOR & | in | ) |
Prepare the serialization of the given vector. The serialization is done by Triangulation::save(). The given vector needs all information on the locally active DoFs (it must be ghosted). See documentation of this class for more information.
| void parallel::distributed::SolutionTransfer< dim, VECTOR, DH >::prepare_serialization | ( | const std::vector< const VECTOR * > & | all_in | ) |
Same as the function above, only for a list of vectors.
| void parallel::distributed::SolutionTransfer< dim, VECTOR, DH >::deserialize | ( | VECTOR & | in | ) |
Execute the deserialization of the given vector. This needs to be done after calling Triangulation::load(). The given vector must be a fully distributed vector without ghost elements. See documentation of this class for more information.
| void parallel::distributed::SolutionTransfer< dim, VECTOR, DH >::deserialize | ( | std::vector< VECTOR * > & | all_in | ) |
Same as the function above, only for a list of vectors.
| void parallel::distributed::SolutionTransfer< dim, VECTOR, DH >::pack_callback | ( | const typename Triangulation< dim, dim >::cell_iterator & | cell, |
| const typename Triangulation< dim, dim >::CellStatus | status, | ||
| void * | data | ||
| ) | [private] |
A callback function used to pack the data on the current mesh into objects that can later be retrieved after refinement, coarsening and repartitioning.
| void parallel::distributed::SolutionTransfer< dim, VECTOR, DH >::unpack_callback | ( | const typename Triangulation< dim, dim >::cell_iterator & | cell, |
| const typename Triangulation< dim, dim >::CellStatus | status, | ||
| const void * | data, | ||
| std::vector< VECTOR * > & | all_out | ||
| ) | [private] |
A callback function used to unpack the data on the current mesh that has been packed up previously on the mesh before refinement, coarsening and repartitioning.
SmartPointer<const DH,SolutionTransfer<dim,VECTOR,DH> > parallel::distributed::SolutionTransfer< dim, VECTOR, DH >::dof_handler [private] |
Pointer to the degree of freedom handler to work with.
Definition at line 205 of file solution_transfer.h.
std::vector<const VECTOR*> parallel::distributed::SolutionTransfer< dim, VECTOR, DH >::input_vectors [private] |
A vector that stores pointers to all the vectors we are supposed to copy over from the old to the new mesh.
Definition at line 214 of file solution_transfer.h.
unsigned int parallel::distributed::SolutionTransfer< dim, VECTOR, DH >::offset [private] |
The offset that the Triangulation has assigned to this object starting at which we are allowed to write.
Definition at line 223 of file solution_transfer.h.
documentation generated on Tue May 22 2012 12:07:05 by
doxygen
1.7.3