<h3>General</h3>
<ol>
-<li> None so far.
+<li> New: parallel::distributed::Triangulation::save()/load() to store
+the refinement information to disk. Also supports saving solution vectors
+using the SolutionTransfer class. (Timo Heister, 2011/10/12)
</ol>
* coarsening a distributed grid and
* handles the necessary communication.
*
- * <h3>Usage</h3>
+ * @note: It is important to note, that if you use more than one
+ * SolutionTransfer object at the same time, that the calls to
+ * prepare_*() and interpolate()/deserialize() need to be in the same order.
+ *
+ * <h3>Note on ghost elements</h3>
+ * 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.
+ *
+ * <h3>Transfering a solution</h3>
+ * Here VECTOR is your favorite vector type, e.g. PETScWrappers::MPI::Vector,
+ * TrilinosWrappers::MPI::Vector, or corresponding blockvectors.
* @verbatim
- * SolutionTransfer<dim, Vector<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> interpolated_solution(dof_handler->n_dofs());
- * soltrans.interpolate(interpolated_solution);
- * @endverbatim
+ 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);
+@endverbatim
+ *
+ * <h3>Use for Serialization</h3>
+ *
+ * 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:
+ *@verbatim
+
+parallel::distributed::SolutionTransfer<dim,VECTOR> sol_trans(dof_handler);
+sol_trans.prepare_serialization (vector);
+
+triangulation.save(filename);
+@endverbatim
+* For deserialization the vector needs to be a distributed vector (without ghost elements):
+* @verbatim
+//[create coarse mesh...]
+triangulation.load(filename);
+
+parallel::distributed::SolutionTransfer<dim,VECTOR> sol_trans(dof_handler);
+sol_trans.deserialize (distributed_vector);
+@endverbatim
+ *
* @ingroup distributed
- * @author Timo Heister, 2009
+ * @author Timo Heister, 2009-2011
*/
template<int dim, typename VECTOR, class DH=DoFHandler<dim> >
class SolutionTransfer
/**
- * return the size in bytes that need
+ * Return the size in bytes that need
* to be stored per cell.
*/
unsigned int get_data_size() const;
+ /**
+ * 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 prepare_serialization(const VECTOR &in);
+
+
+ /**
+ * Same as the function above, only
+ * for a list of vectors.
+ */
+ void prepare_serialization(const std::vector<const VECTOR*> &all_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 deserialize(VECTOR &in);
+
+
+ /**
+ * Same as the function above, only
+ * for a list of vectors.
+ */
+ void deserialize(std::vector<VECTOR*> &all_in);
private:
/**
const void* data,
std::vector<VECTOR*> &all_out);
+
+ /**
+ *
+ */
+ void register_data_attach(const std::size_t size);
+
};
*/
unsigned int get_checksum () const;
+ /**
+ * Save the refinement information from the coarse mesh into the given
+ * file. This file needs to be reachable from all nodes in the computation
+ * on a shared network file system. See the SolutionTransfer class
+ * on how to store solution vectors into this file.
+ */
+ void save(const char* filename) const;
+
+ /**
+ * Load the refinement information saved with save() back in. The mesh
+ * must contain the same coarse mesh that was used in save(). You need
+ * to load with the same number of CPUs that you saved with.
+ */
+ void load(const char* filename);
+
/**
* Used to inform in the callbacks of
* register_data_attach() and
SolutionTransfer<dim, VECTOR, DH>::
prepare_for_coarsening_and_refinement (const std::vector<const VECTOR*> &all_in)
{
- Assert(all_in.size() > 0, ExcMessage("Please transfer atleast one vector!"));
+
input_vectors = all_in;
- SolutionTransfer<dim, VECTOR, DH> *ptr = this;
+ register_data_attach( get_data_size() * input_vectors.size() );
+ }
+
+
+ template<int dim, typename VECTOR, class DH>
+ void
+ SolutionTransfer<dim, VECTOR, DH>::
+ register_data_attach(const std::size_t size)
+ {
+ Assert(size > 0, ExcMessage("Please transfer at least one vector!"));
+
//TODO: casting away constness is bad
parallel::distributed::Triangulation<dim> * tria
= (dynamic_cast<parallel::distributed::Triangulation<dim>*>
Assert (tria != 0, ExcInternalError());
offset
- = tria->register_data_attach(static_cast<size_t>
- (get_data_size() * input_vectors.size()),
+ = tria->register_data_attach(size,
std_cxx1x::bind(&SolutionTransfer<dim, VECTOR, DH>::pack_callback,
- ptr,
+ this,
std_cxx1x::_1,
std_cxx1x::_2,
std_cxx1x::_3));
- }
-
+
+ }
+
template<int dim, typename VECTOR, class DH>
void
+ template<int dim, typename VECTOR, class DH>
+ void
+ SolutionTransfer<dim, VECTOR, DH>::
+ prepare_serialization(const VECTOR &in)
+ {
+ std::vector<const VECTOR*> all_in(1, &in);
+ prepare_serialization(all_in);
+ }
+
+
+
+ template<int dim, typename VECTOR, class DH>
+ void
+ SolutionTransfer<dim, VECTOR, DH>::
+ prepare_serialization(const std::vector<const VECTOR*> &all_in)
+ {
+ prepare_for_coarsening_and_refinement (all_in);
+ }
+
+
+
+ template<int dim, typename VECTOR, class DH>
+ void
+ SolutionTransfer<dim, VECTOR, DH>::
+ deserialize(VECTOR &in)
+ {
+ std::vector<VECTOR*> all_in(1, &in);
+ deserialize(all_in);
+ }
+
+
+
+ template<int dim, typename VECTOR, class DH>
+ void
+ SolutionTransfer<dim, VECTOR, DH>::
+ deserialize(std::vector<VECTOR*> &all_in)
+ {
+ register_data_attach( get_data_size() * all_in.size() );
+
+ // this makes interpolate() happy
+ input_vectors.resize(all_in.size());
+
+ interpolate(all_in);
+ }
+
+
template<int dim, typename VECTOR, class DH>
void
SolutionTransfer<dim, VECTOR, DH>::
#include <algorithm>
#include <numeric>
#include <iostream>
+#include <fstream>
DEAL_II_NAMESPACE_OPEN
}
+ template <int dim, int spacedim>
+ void
+ Triangulation<dim,spacedim>::
+ save(const char* filename) const
+ {
+ int real_data_size = 0;
+ if (attached_data_size>0)
+ real_data_size = attached_data_size+sizeof(CellStatus);
+
+ if (my_subdomain==0)
+ {
+ std::string fname=std::string(filename)+".info";
+ std::ofstream f(fname.c_str());
+ f << Utilities::System::get_n_mpi_processes (mpi_communicator) << " "
+ << real_data_size << " "
+ << attached_data_pack_callbacks.size() << std::endl;
+ }
+
+ if (attached_data_size>0)
+ {
+ const_cast<dealii::parallel::distributed::Triangulation<dim, spacedim>*>(this)
+ ->attach_mesh_data();
+ }
+
+ dealii::internal::p4est::functions<dim>::save(filename, parallel_forest, attached_data_size>0);
+
+ dealii::parallel::distributed::Triangulation<dim, spacedim>* tria
+ = const_cast<dealii::parallel::distributed::Triangulation<dim, spacedim>*>(this);
+
+ tria->n_attached_datas = 0;
+ tria->attached_data_size = 0;
+ tria->attached_data_pack_callbacks.clear();
+
+ // and release the data
+ void * userptr = parallel_forest->user_pointer;
+ dealii::internal::p4est::functions<dim>::reset_data (parallel_forest, 0, NULL, NULL);
+ parallel_forest->user_pointer = userptr;
+ }
+
+
+ template <int dim, int spacedim>
+ void
+ Triangulation<dim,spacedim>::
+ load(const char* filename)
+ {
+ dealii::internal::p4est::functions<dim>::destroy (parallel_forest);
+ parallel_forest = 0;
+ dealii::internal::p4est::functions<dim>::connectivity_destroy (connectivity);
+ connectivity=0;
+
+ unsigned int numcpus, attached_size, attached_count;
+ {
+ std::string fname=std::string(filename)+".info";
+ std::ifstream f(fname.c_str());
+ f >> numcpus >> attached_size >> attached_count;
+ if (numcpus != Utilities::System::get_n_mpi_processes (mpi_communicator))
+ throw ExcInternalError();
+ }
+
+ attached_data_size = 0;
+ n_attached_datas = 0;
+
+ parallel_forest = dealii::internal::p4est::functions<dim>::load (
+ filename, mpi_communicator,
+ attached_size, attached_size>0,
+ this,
+ &connectivity);
+
+ try
+ {
+ copy_local_forest_to_triangulation ();
+ }
+ catch (const typename Triangulation<dim>::DistortedCellList &)
+ {
+ // the underlying
+ // triangulation should not
+ // be checking for
+ // distorted cells
+ AssertThrow (false, ExcInternalError());
+ }
+
+ update_number_cache ();
+ }
+
+
template <int dim, int spacedim>
unsigned int
Triangulation<dim,spacedim>::get_checksum () const