--- /dev/null
+Improved: provide abstract interface to p::d::Triangulation::load()/save().
+<br>
+(Pasquale Claudio Africa, 2021/01/18)
virtual bool
is_multilevel_hierarchy_constructed() const override;
+ /**
+ * Save the triangulation 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. Additional cell-based data can be saved using
+ * register_data_attach().
+ */
+ virtual void
+ save(const std::string &filename) const override;
+
+ /**
+ * Load the triangulation saved with save() back in. The mesh
+ * must be empty before calling this function.
+ *
+ * You need to load with the same number of MPI processes that
+ * you saved with. Cell-based data that was saved with
+ * register_data_attach() can be read in with notify_ready_to_unpack()
+ * after calling load().
+ */
+ virtual void
+ load(const std::string &filename,
+ const bool autopartition = true) override;
+
private:
virtual unsigned int
coarse_cell_id_to_coarse_cell_index(
coarse_cell_index_to_coarse_cell_id(
const unsigned int coarse_cell_index) const override;
+ /**
+ * Go through all locally owned cells and store the relations between
+ * cells and their status in the private member local_cell_relations.
+ *
+ * The stored vector will be ordered by the occurrence of cells.
+ */
+ virtual void
+ update_cell_relations() override;
+
/**
* store the Settings.
*/
bool
is_multilevel_hierarchy_constructed() const override;
+ /**
+ * Transfer data across forests.
+ *
+ * Besides the actual @p parallel_forest, which has been already refined
+ * and repartitioned, this function also needs information about its
+ * previous state, i.e. the locally owned intervals in p4est's
+ * sc_array of each processor. This information needs to be memcopyied
+ * out of the old p4est object and has to be provided via the parameter
+ * @p previous_global_first_quadrant.
+ *
+ * Data has to be previously packed with
+ * DistributedTriangulationBase::DataTransfer::pack_data().
+ */
+ void
+ execute_transfer(
+ const typename dealii::internal::p4est::types<dim>::forest
+ *parallel_forest,
+ const typename dealii::internal::p4est::types<dim>::gloidx
+ *previous_global_first_quadrant);
+
/**
* Implementation of the same function as in the base class.
*
* 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. Additional
- * cell-based data can be saved using register_data_attach().
+ * cell-based data can be saved using
+ * DistributedTriangulationBase::DataTransfer::register_data_attach().
*/
- void
- save(const std::string &filename) const;
+ virtual void
+ save(const std::string &filename) const override;
/**
* Load the refinement information saved with save() back in. The mesh
* you saved with. Rather, if a mesh is loaded with a different number
* of MPI processes than used at the time of saving, the mesh is
* repartitioned appropriately. Cell-based data that was saved with
- * register_data_attach() can be read in with notify_ready_to_unpack()
+ * DistributedTriangulationBase::DataTransfer::register_data_attach() can
+ * be read in with
+ * DistributedTriangulationBase::DataTransfer::notify_ready_to_unpack()
* after calling load().
*
* If you use p4est version > 0.3.4.2 the @p autopartition flag tells
* set to false, the triangulation is only repartitioned if needed (i.e.
* if a different number of MPI processes is encountered).
*/
- void
- load(const std::string &filename, const bool autopartition = true);
-
- /**
- * Register a function that can be used to attach data of fixed size
- * to cells. This is useful for two purposes: (i) Upon refinement and
- * coarsening of a triangulation (in
- * parallel::distributed::Triangulation::execute_coarsening_and_refinement()),
- * one needs to be able to store one or more data vectors per cell that
- * characterizes the solution values on the cell so that this data can
- * then be transferred to the new owning processor of the cell (or
- * its parent/children) when the mesh is re-partitioned; (ii) when
- * serializing a computation to a file, it is necessary to attach
- * data to cells so that it can be saved (in
- * parallel::distributed::Triangulation::save()) along with the cell's
- * other information and, if necessary, later be reloaded from disk
- * with a different subdivision of cells among the processors.
- *
- * The way this function works is that it allows any number of interest
- * parties to register their intent to attach data to cells. One example
- * of classes that do this is parallel::distributed::SolutionTransfer
- * where each parallel::distributed::SolutionTransfer object that works
- * on the current Triangulation object then needs to register its intent.
- * Each of these parties registers a callback function (the first
- * argument here, @p pack_callback) that will be called whenever the
- * triangulation's execute_coarsening_and_refinement() or save()
- * functions are called.
- *
- * The current function then returns an integer handle that corresponds
- * to the number of data set that the callback provided here will attach.
- * While this number could be given a precise meaning, this is
- * not important: You will never actually have to do anything with
- * this number except return it to the notify_ready_to_unpack() function.
- * In other words, each interested party (i.e., the caller of the current
- * function) needs to store their respective returned handle for later use
- * when unpacking data in the callback provided to
- * notify_ready_to_unpack().
- *
- * Whenever @p pack_callback is then called by
- * execute_coarsening_and_refinement() or load() on a given cell, it
- * receives a number of arguments. In particular, the first
- * argument passed to the callback indicates the cell for which
- * it is supposed to attach data. This is always an active cell.
- *
- * The second, CellStatus, argument provided to the callback function
- * will tell you if the given cell will be coarsened, refined, or will
- * persist as is. (This status may be different than the refinement
- * or coarsening flags set on that cell, to accommodate things such as
- * the "one hanging node per edge" rule.). These flags need to be
- * read in context with the p4est quadrant they belong to, as their
- * relations are gathered in local_quadrant_cell_relations.
- *
- * Specifically, the values for this argument mean the following:
- *
- * - `CELL_PERSIST`: The cell won't be refined/coarsened, but might be
- * moved to a different processor. If this is the case, the callback
- * will want to pack up the data on this cell into an array and store
- * it at the provided address for later unpacking wherever this cell
- * may land.
- * - `CELL_REFINE`: This cell will be refined into 4 or 8 cells (in 2d
- * and 3d, respectively). However, because these children don't exist
- * yet, you cannot access them at the time when the callback is
- * called. Thus, in local_quadrant_cell_relations, the corresponding
- * p4est quadrants of the children cells are linked to the deal.II
- * cell which is going to be refined. To be specific, only the very
- * first child is marked with `CELL_REFINE`, whereas the others will be
- * marked with `CELL_INVALID`, which indicates that these cells will be
- * ignored by default during the packing or unpacking process. This
- * ensures that data is only transferred once onto or from the parent
- * cell. If the callback is called with `CELL_REFINE`, the callback
- * will want to pack up the data on this cell into an array and store
- * it at the provided address for later unpacking in a way so that
- * it can then be transferred to the children of the cell that will
- * then be available. In other words, if the data the callback
- * will want to pack up corresponds to a finite element field, then
- * the prolongation from parent to (new) children will have to happen
- * during unpacking.
- * - `CELL_COARSEN`: The children of this cell will be coarsened into the
- * given cell. These children still exist, so if this is the value
- * given to the callback as second argument, the callback will want
- * to transfer data from the children to the current parent cell and
- * pack it up so that it can later be unpacked again on a cell that
- * then no longer has any children (and may also be located on a
- * different processor). In other words, if the data the callback
- * will want to pack up corresponds to a finite element field, then
- * it will need to do the restriction from children to parent at
- * this point.
- * - `CELL_INVALID`: See `CELL_REFINE`.
- *
- * @note If this function is used for serialization of data
- * using save() and load(), then the cell status argument with which
- * the callback is called will always be `CELL_PERSIST`.
- *
- * The callback function is expected to return a memory chunk of the
- * format `std::vector<char>`, representing the packed data on a
- * certain cell.
- *
- * The second parameter @p returns_variable_size_data indicates whether
- * the returned size of the memory region from the callback function
- * varies by cell (<tt>=true</tt>) or stays constant on each one
- * throughout the whole domain (<tt>=false</tt>).
- *
- * @note The purpose of this function is to register intent to
- * attach data for a single, subsequent call to
- * execute_coarsening_and_refinement() and notify_ready_to_unpack(),
- * save(), load(). Consequently, notify_ready_to_unpack(), save(),
- * and load() all forget the registered callbacks once these
- * callbacks have been called, and you will have to re-register
- * them with a triangulation if you want them to be active for
- * another call to these functions.
- */
- unsigned int
- register_data_attach(
- const std::function<std::vector<char>(const cell_iterator &,
- const CellStatus)> &pack_callback,
- const bool returns_variable_size_data);
-
- /**
- * This function is the opposite of register_data_attach(). It is called
- * <i>after</i> the execute_coarsening_and_refinement() or save()/load()
- * functions are done when classes and functions that have previously
- * attached data to a triangulation for either transfer to other
- * processors, across mesh refinement, or serialization of data to
- * a file are ready to receive that data back. The important part about
- * this process is that the triangulation cannot do this right away from
- * the end of execute_coarsening_and_refinement() or load() via a
- * previously attached callback function (as the register_data_attach()
- * function does) because the classes that eventually want the data
- * back may need to do some setup between the point in time where the
- * mesh has been recreated and when the data can actually be received.
- * An example is the parallel::distributed::SolutionTransfer class
- * that can really only receive the data once not only the mesh is
- * completely available again on the current processor, but only
- * after a DoFHandler has been reinitialized and distributed
- * degrees of freedom. In other words, there is typically a significant
- * amount of set up that needs to happen in user space before the classes
- * that can receive data attached to cell are ready to actually do so.
- * When they are, they use the current function to tell the triangulation
- * object that now is the time when they are ready by calling the
- * current function.
- *
- * The supplied callback function is then called for each newly locally
- * owned cell. The first argument to the callback is an iterator that
- * designates the cell; the second argument indicates the status of the
- * cell in question; and the third argument localizes a memory area by
- * two iterators that contains the data that was previously saved from
- * the callback provided to register_data_attach().
- *
- * The CellStatus will indicate if the cell was refined, coarsened, or
- * persisted unchanged. The @p cell_iterator argument to the callback
- * will then either be an active,
- * locally owned cell (if the cell was not refined), or the immediate
- * parent if it was refined during execute_coarsening_and_refinement().
- * Therefore, contrary to during register_data_attach(), you can now
- * access the children if the status is `CELL_REFINE` but no longer for
- * callbacks with status `CELL_COARSEN`.
- *
- * The first argument to this function, `handle`, corresponds to
- * the return value of register_data_attach(). (The precise
- * meaning of what the numeric value of this handle is supposed
- * to represent is neither important, nor should you try to use
- * it for anything other than transmit information between a
- * call to register_data_attach() to the corresponding call to
- * notify_ready_to_unpack().)
- */
- void
- notify_ready_to_unpack(
- const unsigned int handle,
- const std::function<void(
- const cell_iterator &,
- const CellStatus,
- const boost::iterator_range<std::vector<char>::const_iterator> &)>
- &unpack_callback);
+ virtual void
+ load(const std::string &filename,
+ const bool autopartition = true) override;
/**
* Return a permutation vector for the order the coarse cells are handed
typename dealii::internal::p4est::types<dim>::ghost *parallel_ghost;
/**
- * A structure that stores information about the data that has been, or
- * will be, attached to cells via the register_data_attach() function
- * and later retrieved via notify_ready_to_unpack().
- */
- struct CellAttachedData
- {
- /**
- * number of functions that get attached to the Triangulation through
- * register_data_attach() for example SolutionTransfer.
- */
- unsigned int n_attached_data_sets;
-
- /**
- * number of functions that need to unpack their data after a call from
- * load()
- */
- unsigned int n_attached_deserialize;
-
- using pack_callback_t = std::function<std::vector<char>(
- typename Triangulation<dim, spacedim>::cell_iterator,
- CellStatus)>;
-
- /**
- * These callback functions will be stored in the order in which they
- * have been registered with the register_data_attach() function.
- */
- std::vector<pack_callback_t> pack_callbacks_fixed;
- std::vector<pack_callback_t> pack_callbacks_variable;
- };
-
- CellAttachedData cell_attached_data;
-
- /**
- * This auxiliary data structure stores the relation between a p4est
- * quadrant, a deal.II cell and its current CellStatus. For an extensive
- * description of the latter, see the documentation for the member
- * function register_data_attach().
- */
- using quadrant_cell_relation_t = typename std::tuple<
- typename dealii::internal::p4est::types<dim>::quadrant *,
- CellStatus,
- cell_iterator>;
-
- /**
- * Vector of tuples, which each contain a p4est quadrant, a deal.II cell
- * and their relation after refinement. To update its contents, use the
- * compute_quadrant_cell_relations member function.
- *
- * The size of this vector is assumed to be equal to the number of locally
- * owned quadrants in the parallel_forest object.
- */
- std::vector<quadrant_cell_relation_t> local_quadrant_cell_relations;
-
- /**
- * Go through all p4est trees and store the relations between locally
- * owned quadrants and cells in the private member
- * local_quadrant_cell_relations.
+ * Go through all p4est trees and store the relations between the status
+ * of locally owned quadrants and cells in the private member
+ * local_cell_relations.
*
* The stored vector will be ordered by the occurrence of quadrants in
* the corresponding local sc_array of the parallel_forest. p4est requires
* this specific ordering for its transfer functions.
*/
- void
- update_quadrant_cell_relations();
-
- /**
- * This class in the private scope of parallel::distributed::Triangulation
- * is dedicated to the data transfer across repartitioned meshes
- * and to the file system.
- *
- * It is designed to store all data buffers intended for transfer
- * separated from the parallel_forest and to interface with p4est
- * where it is absolutely necessary.
- */
- class DataTransfer
- {
- public:
- DataTransfer(const MPI_Comm &mpi_communicator);
-
- /**
- * Prepare data transfer by calling the pack callback functions on each
- * cell
- * in @p quad_cell_relations.
- *
- * All registered callback functions in @p pack_callbacks_fixed will write
- * into the fixed size buffer, whereas each entry of @p pack_callbacks_variable
- * will write its data into the variable size buffer.
- */
- void
- pack_data(
- const std::vector<quadrant_cell_relation_t> &quad_cell_relations,
- const std::vector<typename CellAttachedData::pack_callback_t>
- &pack_callbacks_fixed,
- const std::vector<typename CellAttachedData::pack_callback_t>
- &pack_callbacks_variable);
-
- /**
- * Transfer data across forests.
- *
- * Besides the actual @p parallel_forest, which has been already refined
- * and repartitioned, this function also needs information about its
- * previous state, i.e. the locally owned intervals in p4est's
- * sc_array of each processor. This information needs to be memcopyied
- * out of the old p4est object and has to be provided via the parameter
- * @p previous_global_first_quadrant.
- *
- * Data has to be previously packed with pack_data().
- */
- void
- execute_transfer(
- const typename dealii::internal::p4est::types<dim>::forest
- *parallel_forest,
- const typename dealii::internal::p4est::types<dim>::gloidx
- *previous_global_first_quadrant);
-
- /**
- * Unpack the CellStatus information on each entry of
- * @p quad_cell_relations.
- *
- * Data has to be previously transferred with execute_transfer()
- * or read from the file system via load().
- */
- void
- unpack_cell_status(
- std::vector<quadrant_cell_relation_t> &quad_cell_relations) const;
-
- /**
- * Unpack previously transferred data on each cell registered in
- * @p quad_cell_relations with the provided @p unpack_callback function.
- *
- * The parameter @p handle corresponds to the position where the
- * @p unpack_callback function is allowed to read from the memory. Its
- * value needs to be in accordance with the corresponding pack_callback
- * function that has been registered previously.
- *
- * Data has to be previously transferred with execute_transfer()
- * or read from the file system via load().
- */
- void
- unpack_data(
- const std::vector<quadrant_cell_relation_t> &quad_cell_relations,
- const unsigned int handle,
- const std::function<void(
- const typename dealii::Triangulation<dim, spacedim>::cell_iterator
- &,
- const typename dealii::Triangulation<dim, spacedim>::CellStatus &,
- const boost::iterator_range<std::vector<char>::const_iterator> &)>
- &unpack_callback) const;
-
- /**
- * Transfer data to file system.
- *
- * The data will be written in a separate file, whose name
- * consists of the stem @p filename and an attached identifier
- * <tt>-fixed.data</tt> for fixed size data and <tt>-variable.data</tt>
- * for variable size data.
- *
- * All processors write into these files simultaneously via MPIIO.
- * Each processor's position to write to will be determined
- * from the provided @p parallel_forest.
- *
- * Data has to be previously packed with pack_data().
- */
- void
- save(const typename dealii::internal::p4est::types<dim>::forest
- * parallel_forest,
- const std::string &filename) const;
-
- /**
- * Transfer data from file system.
- *
- * The data will be read from separate file, whose name
- * consists of the stem @p filename and an attached identifier
- * <tt>-fixed.data</tt> for fixed size data and <tt>-variable.data</tt>
- * for variable size data.
- * The @p n_attached_deserialize_fixed and @p n_attached_deserialize_variable
- * parameters are required to gather the memory offsets for each
- * callback.
- *
- * All processors read from these files simultaneously via MPIIO.
- * Each processor's position to read from will be determined
- * from the provided @p parallel_forest.
- *
- * After loading, unpack_data() needs to be called to finally
- * distribute data across the associated triangulation.
- */
- void
- load(const typename dealii::internal::p4est::types<dim>::forest
- * parallel_forest,
- const std::string &filename,
- const unsigned int n_attached_deserialize_fixed,
- const unsigned int n_attached_deserialize_variable);
-
- /**
- * Clears all containers and associated data, and resets member
- * values to their default state.
- *
- * Frees memory completely.
- */
- void
- clear();
-
- private:
- MPI_Comm mpi_communicator;
-
- /**
- * Flag that denotes if variable size data has been packed.
- */
- bool variable_size_data_stored;
-
- /**
- * Cumulative size in bytes that those functions that have called
- * register_data_attach() want to attach to each cell. This number
- * only pertains to fixed-sized buffers where the data attached to
- * each cell has exactly the same size.
- *
- * The last entry of this container corresponds to the data size
- * packed per cell in the fixed size buffer (which can be accessed
- * calling <tt>sizes_fixed_cumulative.back()</tt>).
- */
- std::vector<unsigned int> sizes_fixed_cumulative;
-
- /**
- * Consecutive buffers designed for the fixed size transfer
- * functions of p4est.
- */
- std::vector<char> src_data_fixed;
- std::vector<char> dest_data_fixed;
-
- /**
- * Consecutive buffers designed for the variable size transfer
- * functions of p4est.
- */
- std::vector<int> src_sizes_variable;
- std::vector<int> dest_sizes_variable;
- std::vector<char> src_data_variable;
- std::vector<char> dest_data_variable;
- };
-
- DataTransfer data_transfer;
+ virtual void
+ update_cell_relations() override;
/**
* Two arrays that store which p4est tree corresponds to which coarse
* This function is not implemented, but needs to be present for the
* compiler.
*/
- void
- load(const std::string &filename, const bool autopartition = true);
+ virtual void
+ load(const std::string &filename,
+ const bool autopartition = true) override;
/**
* This function is not implemented, but needs to be present for the
* compiler.
*/
- void
- save(const std::string &filename) const;
-
- bool
- is_multilevel_hierarchy_constructed() const override;
+ virtual void
+ save(const std::string &filename) const override;
/**
* This function is not implemented, but needs to be present for the
* compiler.
*/
- unsigned int
- register_data_attach(
- const std::function<std::vector<char>(
- const typename dealii::Triangulation<1, spacedim>::cell_iterator &,
- const typename dealii::Triangulation<1, spacedim>::CellStatus)>
- & pack_callback,
- const bool returns_variable_size_data);
+ virtual bool
+ is_multilevel_hierarchy_constructed() const override;
/**
* This function is not implemented, but needs to be present for the
* compiler.
*/
- void
- notify_ready_to_unpack(
- const unsigned int handle,
- const std::function<void(
- const typename dealii::Triangulation<1, spacedim>::cell_iterator &,
- const typename dealii::Triangulation<1, spacedim>::CellStatus,
- const boost::iterator_range<std::vector<char>::const_iterator> &)>
- &unpack_callback);
+ virtual void
+ update_cell_relations() override;
/**
* Dummy arrays. This class isn't usable but the compiler wants to see
const typename dealii::Triangulation<dim, spacedim>::MeshSmoothing
smooth_grid = (dealii::Triangulation<dim, spacedim>::none),
const bool check_for_distorted_cells = false);
+
+ /**
+ * Reset this triangulation into a virgin state by deleting all data.
+ *
+ * Note that this operation is only allowed if no subscriptions to this
+ * object exist any more, such as DoFHandler objects using it.
+ */
+ virtual void
+ clear() override;
+
+ using cell_iterator =
+ typename dealii::Triangulation<dim, spacedim>::cell_iterator;
+
+ using CellStatus =
+ typename dealii::Triangulation<dim, spacedim>::CellStatus;
+
+ /**
+ * Save the triangulation 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. Additional cell-based data can be saved using
+ * register_data_attach().
+ */
+ virtual void
+ save(const std::string &filename) const = 0;
+
+ /**
+ * Load the triangulation saved with save() back in. Cell-based data that
+ * was saved with register_data_attach() can be read in with
+ * notify_ready_to_unpack() after calling load().
+ */
+ virtual void
+ load(const std::string &filename, const bool autopartition = true) = 0;
+
+ /**
+ * Register a function that can be used to attach data of fixed size
+ * to cells. This is useful for two purposes: (i) Upon refinement and
+ * coarsening of a triangulation (@a e.g. in
+ * parallel::distributed::Triangulation::execute_coarsening_and_refinement()),
+ * one needs to be able to store one or more data vectors per cell that
+ * characterizes the solution values on the cell so that this data can
+ * then be transferred to the new owning processor of the cell (or
+ * its parent/children) when the mesh is re-partitioned; (ii) when
+ * serializing a computation to a file, it is necessary to attach
+ * data to cells so that it can be saved (@a e.g. in
+ * parallel::distributed::Triangulation::save()) along with the cell's
+ * other information and, if necessary, later be reloaded from disk
+ * with a different subdivision of cells among the processors.
+ *
+ * The way this function works is that it allows any number of interest
+ * parties to register their intent to attach data to cells. One example
+ * of classes that do this is parallel::distributed::SolutionTransfer
+ * where each parallel::distributed::SolutionTransfer object that works
+ * on the current Triangulation object then needs to register its intent.
+ * Each of these parties registers a callback function (the first
+ * argument here, @p pack_callback) that will be called whenever the
+ * triangulation's execute_coarsening_and_refinement() or save()
+ * functions are called.
+ *
+ * The current function then returns an integer handle that corresponds
+ * to the number of data set that the callback provided here will attach.
+ * While this number could be given a precise meaning, this is
+ * not important: You will never actually have to do anything with
+ * this number except return it to the notify_ready_to_unpack() function.
+ * In other words, each interested party (i.e., the caller of the current
+ * function) needs to store their respective returned handle for later use
+ * when unpacking data in the callback provided to
+ * notify_ready_to_unpack().
+ *
+ * Whenever @p pack_callback is then called by
+ * execute_coarsening_and_refinement() or load() on a given cell, it
+ * receives a number of arguments. In particular, the first
+ * argument passed to the callback indicates the cell for which
+ * it is supposed to attach data. This is always an active cell.
+ *
+ * The second, CellStatus, argument provided to the callback function
+ * will tell you if the given cell will be coarsened, refined, or will
+ * persist as is. (This status may be different than the refinement
+ * or coarsening flags set on that cell, to accommodate things such as
+ * the "one hanging node per edge" rule.). These flags need to be
+ * read in context with the p4est quadrant they belong to, as their
+ * relations are gathered in local_quadrant_cell_relations.
+ *
+ * Specifically, the values for this argument mean the following:
+ *
+ * - `CELL_PERSIST`: The cell won't be refined/coarsened, but might be
+ * moved to a different processor. If this is the case, the callback
+ * will want to pack up the data on this cell into an array and store
+ * it at the provided address for later unpacking wherever this cell
+ * may land.
+ * - `CELL_REFINE`: This cell will be refined into 4 or 8 cells (in 2d
+ * and 3d, respectively). However, because these children don't exist
+ * yet, you cannot access them at the time when the callback is
+ * called. Thus, in local_quadrant_cell_relations, the corresponding
+ * p4est quadrants of the children cells are linked to the deal.II
+ * cell which is going to be refined. To be specific, only the very
+ * first child is marked with `CELL_REFINE`, whereas the others will be
+ * marked with `CELL_INVALID`, which indicates that these cells will be
+ * ignored by default during the packing or unpacking process. This
+ * ensures that data is only transferred once onto or from the parent
+ * cell. If the callback is called with `CELL_REFINE`, the callback
+ * will want to pack up the data on this cell into an array and store
+ * it at the provided address for later unpacking in a way so that
+ * it can then be transferred to the children of the cell that will
+ * then be available. In other words, if the data the callback
+ * will want to pack up corresponds to a finite element field, then
+ * the prolongation from parent to (new) children will have to happen
+ * during unpacking.
+ * - `CELL_COARSEN`: The children of this cell will be coarsened into the
+ * given cell. These children still exist, so if this is the value
+ * given to the callback as second argument, the callback will want
+ * to transfer data from the children to the current parent cell and
+ * pack it up so that it can later be unpacked again on a cell that
+ * then no longer has any children (and may also be located on a
+ * different processor). In other words, if the data the callback
+ * will want to pack up corresponds to a finite element field, then
+ * it will need to do the restriction from children to parent at
+ * this point.
+ * - `CELL_INVALID`: See `CELL_REFINE`.
+ *
+ * @note If this function is used for serialization of data
+ * using save() and load(), then the cell status argument with which
+ * the callback is called will always be `CELL_PERSIST`.
+ *
+ * The callback function is expected to return a memory chunk of the
+ * format `std::vector<char>`, representing the packed data on a
+ * certain cell.
+ *
+ * The second parameter @p returns_variable_size_data indicates whether
+ * the returned size of the memory region from the callback function
+ * varies by cell (<tt>=true</tt>) or stays constant on each one
+ * throughout the whole domain (<tt>=false</tt>).
+ *
+ * @note The purpose of this function is to register intent to
+ * attach data for a single, subsequent call to
+ * execute_coarsening_and_refinement() and notify_ready_to_unpack(),
+ * save(), load(). Consequently, notify_ready_to_unpack(), save(),
+ * and load() all forget the registered callbacks once these
+ * callbacks have been called, and you will have to re-register
+ * them with a triangulation if you want them to be active for
+ * another call to these functions.
+ */
+ unsigned int
+ register_data_attach(
+ const std::function<std::vector<char>(const cell_iterator &,
+ const CellStatus)> &pack_callback,
+ const bool returns_variable_size_data);
+
+ /**
+ * This function is the opposite of register_data_attach(). It is called
+ * <i>after</i> the execute_coarsening_and_refinement() or save()/load()
+ * functions are done when classes and functions that have previously
+ * attached data to a triangulation for either transfer to other
+ * processors, across mesh refinement, or serialization of data to
+ * a file are ready to receive that data back. The important part about
+ * this process is that the triangulation cannot do this right away from
+ * the end of execute_coarsening_and_refinement() or load() via a
+ * previously attached callback function (as the register_data_attach()
+ * function does) because the classes that eventually want the data
+ * back may need to do some setup between the point in time where the
+ * mesh has been recreated and when the data can actually be received.
+ * An example is the parallel::distributed::SolutionTransfer class
+ * that can really only receive the data once not only the mesh is
+ * completely available again on the current processor, but only
+ * after a DoFHandler has been reinitialized and distributed
+ * degrees of freedom. In other words, there is typically a significant
+ * amount of set up that needs to happen in user space before the classes
+ * that can receive data attached to cell are ready to actually do so.
+ * When they are, they use the current function to tell the triangulation
+ * object that now is the time when they are ready by calling the
+ * current function.
+ *
+ * The supplied callback function is then called for each newly locally
+ * owned cell. The first argument to the callback is an iterator that
+ * designates the cell; the second argument indicates the status of the
+ * cell in question; and the third argument localizes a memory area by
+ * two iterators that contains the data that was previously saved from
+ * the callback provided to register_data_attach().
+ *
+ * The CellStatus will indicate if the cell was refined, coarsened, or
+ * persisted unchanged. The @p cell_iterator argument to the callback
+ * will then either be an active,
+ * locally owned cell (if the cell was not refined), or the immediate
+ * parent if it was refined during execute_coarsening_and_refinement().
+ * Therefore, contrary to during register_data_attach(), you can now
+ * access the children if the status is `CELL_REFINE` but no longer for
+ * callbacks with status `CELL_COARSEN`.
+ *
+ * The first argument to this function, `handle`, corresponds to
+ * the return value of register_data_attach(). (The precise
+ * meaning of what the numeric value of this handle is supposed
+ * to represent is neither important, nor should you try to use
+ * it for anything other than transmit information between a
+ * call to register_data_attach() to the corresponding call to
+ * notify_ready_to_unpack().)
+ */
+ void
+ notify_ready_to_unpack(
+ const unsigned int handle,
+ const std::function<
+ void(const cell_iterator &,
+ const CellStatus,
+ const boost::iterator_range<std::vector<char>::const_iterator> &)>
+ &unpack_callback);
+
+ protected:
+ /**
+ * Save additional cell-attached data into the given file. The first
+ * arguments are used to determine the offsets where to write buffers to.
+ *
+ * Called by @ref save.
+ */
+ void
+ save_attached_data(const unsigned int global_first_cell,
+ const unsigned int global_num_cells,
+ const std::string &filename) const;
+
+ /**
+ * Load additional cell-attached data from the given file, if any was saved.
+ * The first arguments are used to determine the offsets where to read
+ * buffers from.
+ *
+ * Called by @ref load.
+ */
+ void
+ load_attached_data(const unsigned int global_first_cell,
+ const unsigned int global_num_cells,
+ const unsigned int local_num_cells,
+ const std::string &filename,
+ const unsigned int n_attached_deserialize_fixed,
+ const unsigned int n_attached_deserialize_variable);
+
+ /**
+ * Go through all cells and store the relations between locally
+ * owned quadrants and cells in the private member
+ * local_cell_relations.
+ *
+ * The stored vector will be ordered by the occurrence of quadrants.
+ */
+ virtual void
+ update_cell_relations() = 0;
+
+ /**
+ * This auxiliary data structure stores the relation between
+ * a deal.II cell and its current CellStatus. For an extensive
+ * description of the latter, see the documentation for the member
+ * function register_data_attach().
+ */
+ using cell_relation_t = typename std::tuple<CellStatus, cell_iterator>;
+
+ /**
+ * Vector of tuples, which each contain a deal.II cell
+ * and their relation after refinement. To update its contents, use the
+ * compute_cell_relations member function.
+ *
+ * The size of this vector is assumed to be equal to the number of locally
+ * owned quadrants in the parallel_forest object.
+ */
+ std::vector<cell_relation_t> local_cell_relations;
+
+ /**
+ * A structure that stores information about the data that has been, or
+ * will be, attached to cells via the register_data_attach() function
+ * and later retrieved via notify_ready_to_unpack().
+ */
+ struct CellAttachedData
+ {
+ /**
+ * number of functions that get attached to the Triangulation through
+ * register_data_attach() for example SolutionTransfer.
+ */
+ unsigned int n_attached_data_sets;
+
+ /**
+ * number of functions that need to unpack their data after a call from
+ * load()
+ */
+ unsigned int n_attached_deserialize;
+
+ using pack_callback_t = std::function<std::vector<char>(
+ typename dealii::Triangulation<dim, spacedim>::cell_iterator,
+ typename dealii::Triangulation<dim, spacedim>::CellStatus)>;
+
+ /**
+ * These callback functions will be stored in the order in which they
+ * have been registered with the register_data_attach() function.
+ */
+ std::vector<pack_callback_t> pack_callbacks_fixed;
+ std::vector<pack_callback_t> pack_callbacks_variable;
+ };
+
+ CellAttachedData cell_attached_data;
+
+ /**
+ * This class in the private scope of parallel::DistributedTriangulationBase
+ * is dedicated to the data transfer across repartitioned meshes
+ * and to/from the file system.
+ *
+ * It is designed to store all data buffers intended for transfer.
+ */
+ class DataTransfer
+ {
+ public:
+ DataTransfer(const MPI_Comm &mpi_communicator);
+
+ /**
+ * Prepare data transfer by calling the pack callback functions on each
+ * cell
+ * in @p cell_relations.
+ *
+ * All registered callback functions in @p pack_callbacks_fixed will write
+ * into the fixed size buffer, whereas each entry of @p pack_callbacks_variable
+ * will write its data into the variable size buffer.
+ */
+ void
+ pack_data(const std::vector<cell_relation_t> &cell_relations,
+ const std::vector<typename CellAttachedData::pack_callback_t>
+ &pack_callbacks_fixed,
+ const std::vector<typename CellAttachedData::pack_callback_t>
+ &pack_callbacks_variable);
+
+
+
+ /**
+ * Unpack the CellStatus information on each entry of
+ * @p cell_relations.
+ *
+ * Data has to be previously transferred with execute_transfer()
+ * or read from the file system via load().
+ */
+ void
+ unpack_cell_status(std::vector<cell_relation_t> &cell_relations) const;
+
+ /**
+ * Unpack previously transferred data on each cell registered in
+ * @p cell_relations with the provided @p unpack_callback function.
+ *
+ * The parameter @p handle corresponds to the position where the
+ * @p unpack_callback function is allowed to read from the memory. Its
+ * value needs to be in accordance with the corresponding pack_callback
+ * function that has been registered previously.
+ *
+ * Data has to be previously transferred with execute_transfer()
+ * or read from the file system via load().
+ */
+ void
+ unpack_data(
+ const std::vector<cell_relation_t> &cell_relations,
+ const unsigned int handle,
+ const std::function<void(
+ const typename dealii::Triangulation<dim, spacedim>::cell_iterator &,
+ const typename dealii::Triangulation<dim, spacedim>::CellStatus &,
+ const boost::iterator_range<std::vector<char>::const_iterator> &)>
+ &unpack_callback) const;
+
+ /**
+ * Transfer data to file system.
+ *
+ * The data will be written in a separate file, whose name
+ * consists of the stem @p filename and an attached identifier
+ * <tt>_fixed.data</tt> for fixed size data and <tt>_variable.data</tt>
+ * for variable size data.
+ *
+ * All processors write into these files simultaneously via MPIIO.
+ * Each processor's position to write to will be determined
+ * from the provided input parameters.
+ *
+ * Data has to be previously packed with pack_data().
+ */
+ void
+ save(const unsigned int global_first_cell,
+ const unsigned int global_num_cells,
+ const std::string &filename) const;
+
+ /**
+ * Transfer data from file system.
+ *
+ * The data will be read from separate file, whose name
+ * consists of the stem @p filename and an attached identifier
+ * <tt>_fixed.data</tt> for fixed size data and <tt>_variable.data</tt>
+ * for variable size data.
+ * The @p n_attached_deserialize_fixed and @p n_attached_deserialize_variable
+ * parameters are required to gather the memory offsets for each
+ * callback.
+ *
+ * All processors read from these files simultaneously via MPIIO.
+ * Each processor's position to read from will be determined
+ * from the provided input arguments.
+ *
+ * After loading, unpack_data() needs to be called to finally
+ * distribute data across the associated triangulation.
+ */
+ void
+ load(const unsigned int global_first_cell,
+ const unsigned int global_num_cells,
+ const unsigned int local_num_cells,
+ const std::string &filename,
+ const unsigned int n_attached_deserialize_fixed,
+ const unsigned int n_attached_deserialize_variable);
+
+ /**
+ * Clears all containers and associated data, and resets member
+ * values to their default state.
+ *
+ * Frees memory completely.
+ */
+ void
+ clear();
+
+ /**
+ * Flag that denotes if variable size data has been packed.
+ */
+ bool variable_size_data_stored;
+
+ /**
+ * Cumulative size in bytes that those functions that have called
+ * register_data_attach() want to attach to each cell. This number
+ * only pertains to fixed-sized buffers where the data attached to
+ * each cell has exactly the same size.
+ *
+ * The last entry of this container corresponds to the data size
+ * packed per cell in the fixed size buffer (which can be accessed
+ * calling <tt>sizes_fixed_cumulative.back()</tt>).
+ */
+ std::vector<unsigned int> sizes_fixed_cumulative;
+
+ /**
+ * Consecutive buffers designed for the fixed size transfer
+ * functions of p4est.
+ */
+ std::vector<char> src_data_fixed;
+ std::vector<char> dest_data_fixed;
+
+ /**
+ * Consecutive buffers designed for the variable size transfer
+ * functions of p4est.
+ */
+ std::vector<int> src_sizes_variable;
+ std::vector<int> dest_sizes_variable;
+ std::vector<char> src_data_variable;
+ std::vector<char> dest_data_variable;
+
+ private:
+ MPI_Comm mpi_communicator;
+ };
+
+ DataTransfer data_transfer;
};
} // namespace parallel
#include <deal.II/grid/grid_tools.h>
+#include <fstream>
#include <memory>
DEAL_II_NAMESPACE_OPEN
}
this->update_number_cache();
+ // this->update_cell_relations();
}
}
+ template <int dim, int spacedim>
+ void
+ Triangulation<dim, spacedim>::update_cell_relations()
+ {
+ AssertThrow(false, ExcNotImplemented());
+ }
+
+ template <int dim, int spacedim>
+ void
+ Triangulation<dim, spacedim>::save(const std::string &filename) const
+ {
+ (void)filename;
+
+ AssertThrow(false, ExcNotImplemented());
+ }
+
+
+
+ template <int dim, int spacedim>
+ void
+ Triangulation<dim, spacedim>::load(const std::string &filename,
+ const bool autopartition)
+ {
+ (void)filename;
+ (void)autopartition;
+
+ AssertThrow(false, ExcNotImplemented());
+ }
+
} // namespace fullydistributed
} // namespace parallel
return *this_object->current_pointer++;
}
-
-
template <int dim, int spacedim>
- using quadrant_cell_relation_t = typename std::tuple<
- typename dealii::internal::p4est::types<dim>::quadrant *,
+ using cell_relation_t = typename std::tuple<
typename dealii::Triangulation<dim, spacedim>::CellStatus,
typename dealii::Triangulation<dim, spacedim>::cell_iterator>;
-
-
/**
* Adds a tuple of a p4est quadrant, @p status and @p dealii_cell
- * to the vector containing all relations @p quad_cell_rel.
+ * to the vector containing all relations @p cell_rel.
* The tuple will be inserted in the position corresponding to the one
* of the p4est quadrant in the underlying p4est sc_array. The position
* will be determined from @p idx, which is the position of the quadrant
*/
template <int dim, int spacedim>
inline void
- add_single_quadrant_cell_relation(
- std::vector<quadrant_cell_relation_t<dim, spacedim>> & quad_cell_rel,
+ add_single_cell_relation(
+ std::vector<cell_relation_t<dim, spacedim>> & cell_rel,
const typename dealii::internal::p4est::types<dim>::tree & tree,
const unsigned int idx,
const typename Triangulation<dim, spacedim>::cell_iterator &dealii_cell,
{
const unsigned int local_quadrant_index = tree.quadrants_offset + idx;
- const auto q =
- static_cast<typename dealii::internal::p4est::types<dim>::quadrant *>(
- sc_array_index(const_cast<sc_array_t *>(&tree.quadrants), idx));
-
// check if we will be writing into valid memory
- Assert(local_quadrant_index < quad_cell_rel.size(), ExcInternalError());
+ Assert(local_quadrant_index < cell_rel.size(), ExcInternalError());
// store relation
- quad_cell_rel[local_quadrant_index] =
- std::make_tuple(q, status, dealii_cell);
+ cell_rel[local_quadrant_index] = std::make_tuple(status, dealii_cell);
}
/**
* This is the recursive part of the member function
- * update_quadrant_cell_relations().
+ * update_cell_relations().
*
* Find the relation between the @p p4est_cell and the @p dealii_cell in the
* corresponding @p tree. Depending on the CellStatus relation between the two,
- * a new entry will either be inserted in @p quad_cell_rel or the recursion
+ * a new entry will either be inserted in @p cell_rel or the recursion
* will be continued.
*/
template <int dim, int spacedim>
void
- update_quadrant_cell_relations_recursively(
- std::vector<quadrant_cell_relation_t<dim, spacedim>> & quad_cell_rel,
+ update_cell_relations_recursively(
+ std::vector<cell_relation_t<dim, spacedim>> & cell_rel,
const typename dealii::internal::p4est::types<dim>::tree & tree,
const typename Triangulation<dim, spacedim>::cell_iterator & dealii_cell,
const typename dealii::internal::p4est::types<dim>::quadrant &p4est_cell)
for (unsigned int c = 0; c < GeometryInfo<dim>::max_children_per_cell;
++c)
{
- update_quadrant_cell_relations_recursively<dim, spacedim>(
- quad_cell_rel, tree, dealii_cell->child(c), p4est_child[c]);
+ update_cell_relations_recursively<dim, spacedim>(
+ cell_rel, tree, dealii_cell->child(c), p4est_child[c]);
}
}
else if (!p4est_has_children && !dealii_cell->has_children())
{
// this active cell didn't change
// save tuple into corresponding position
- add_single_quadrant_cell_relation<dim, spacedim>(
- quad_cell_rel,
+ add_single_cell_relation<dim, spacedim>(
+ cell_rel,
tree,
idx,
dealii_cell,
cell_status = (i == 0) ? Triangulation<dim, spacedim>::CELL_REFINE :
Triangulation<dim, spacedim>::CELL_INVALID;
- add_single_quadrant_cell_relation<dim, spacedim>(
- quad_cell_rel, tree, child_idx, dealii_cell, cell_status);
+ add_single_cell_relation<dim, spacedim>(
+ cell_rel, tree, child_idx, dealii_cell, cell_status);
}
}
else // based on the conditions above, we know that p4est_cell has no
{
// its children got coarsened into this cell in p4est,
// but the dealii_cell still has its children
- add_single_quadrant_cell_relation<dim, spacedim>(
- quad_cell_rel,
+ add_single_cell_relation<dim, spacedim>(
+ cell_rel,
tree,
idx,
dealii_cell,
{
namespace distributed
{
- /* ------------------ class DataTransfer<dim,spacedim> ----------------- */
-
-
- template <int dim, int spacedim>
- Triangulation<dim, spacedim>::DataTransfer::DataTransfer(
- const MPI_Comm &mpi_communicator)
- : mpi_communicator(mpi_communicator)
- , variable_size_data_stored(false)
- {}
-
-
-
- template <int dim, int spacedim>
- void
- Triangulation<dim, spacedim>::DataTransfer::pack_data(
- const std::vector<quadrant_cell_relation_t> &quad_cell_relations,
- const std::vector<typename CellAttachedData::pack_callback_t>
- &pack_callbacks_fixed,
- const std::vector<typename CellAttachedData::pack_callback_t>
- &pack_callbacks_variable)
- {
- Assert(src_data_fixed.size() == 0,
- ExcMessage("Previously packed data has not been released yet!"));
- Assert(src_sizes_variable.size() == 0, ExcInternalError());
-
- const unsigned int n_callbacks_fixed = pack_callbacks_fixed.size();
- const unsigned int n_callbacks_variable = pack_callbacks_variable.size();
-
- // Store information that we packed variable size data in
- // a member variable for later.
- variable_size_data_stored = (n_callbacks_variable > 0);
-
- // If variable transfer is scheduled, we will store the data size that
- // each variable size callback function writes in this auxiliary
- // container. The information will be stored by each cell in this vector
- // temporarily.
- std::vector<unsigned int> cell_sizes_variable_cumulative(
- n_callbacks_variable);
-
- // Prepare the buffer structure, in which each callback function will
- // store its data for each active cell.
- // The outmost shell in this container construct corresponds to the
- // data packed per cell. The next layer resembles the data that
- // each callback function packs on the corresponding cell. These
- // buffers are chains of chars stored in an std::vector<char>.
- // A visualisation of the data structure:
- /* clang-format off */
- // | cell_1 | | cell_2 | ...
- // || callback_1 || callback_2 |...| || callback_1 || callback_2 |...| ...
- // |||char|char|...|||char|char|...|...| |||char|char|...|||char|char|...|...| ...
- /* clang-format on */
- std::vector<std::vector<std::vector<char>>> packed_fixed_size_data(
- quad_cell_relations.size());
- std::vector<std::vector<std::vector<char>>> packed_variable_size_data(
- variable_size_data_stored ? quad_cell_relations.size() : 0);
-
- //
- // --------- Pack data for fixed and variable size transfer ---------
- //
- // Iterate over all cells, call all callback functions on each cell,
- // and store their data in the corresponding buffer scope.
- {
- auto quad_cell_rel_it = quad_cell_relations.cbegin();
- auto data_cell_fixed_it = packed_fixed_size_data.begin();
- auto data_cell_variable_it = packed_variable_size_data.begin();
- for (; quad_cell_rel_it != quad_cell_relations.cend();
- ++quad_cell_rel_it, ++data_cell_fixed_it)
- {
- const auto &cell_status = std::get<1>(*quad_cell_rel_it);
- const auto &dealii_cell = std::get<2>(*quad_cell_rel_it);
-
- // Assertions about the tree structure.
- switch (cell_status)
- {
- case parallel::distributed::Triangulation<dim, spacedim>::
- CELL_PERSIST:
- case parallel::distributed::Triangulation<dim, spacedim>::
- CELL_REFINE:
- // double check the condition that we will only ever attach
- // data to active cells when we get here
- Assert(dealii_cell->is_active(), ExcInternalError());
- break;
-
- case parallel::distributed::Triangulation<dim, spacedim>::
- CELL_COARSEN:
- // double check the condition that we will only ever attach
- // data to cells with children when we get here. however, we
- // can only tolerate one level of coarsening at a time, so
- // check that the children are all active
- Assert(dealii_cell->is_active() == false, ExcInternalError());
- for (unsigned int c = 0;
- c < GeometryInfo<dim>::max_children_per_cell;
- ++c)
- Assert(dealii_cell->child(c)->is_active(),
- ExcInternalError());
- break;
-
- case parallel::distributed::Triangulation<dim, spacedim>::
- CELL_INVALID:
- // do nothing on invalid cells
- break;
-
- default:
- Assert(false, ExcInternalError());
- break;
- }
-
- // Reserve memory corresponding to the number of callback
- // functions that will be called.
- // If variable size transfer is scheduled, we need to leave
- // room for an array that holds information about how many
- // bytes each of the variable size callback functions will
- // write.
- // On cells flagged with CELL_INVALID, only its CellStatus
- // will be stored.
- const unsigned int n_fixed_size_data_sets_on_cell =
- 1 +
- ((cell_status ==
- parallel::distributed::Triangulation<dim,
- spacedim>::CELL_INVALID) ?
- 0 :
- ((variable_size_data_stored ? 1 : 0) + n_callbacks_fixed));
- data_cell_fixed_it->resize(n_fixed_size_data_sets_on_cell);
-
- // We continue with packing all data on this specific cell.
- auto data_fixed_it = data_cell_fixed_it->begin();
-
- // First, we pack the CellStatus information.
- // to get consistent data sizes on each cell for the fixed size
- // transfer, we won't allow compression
- *data_fixed_it =
- Utilities::pack(cell_status, /*allow_compression=*/false);
- ++data_fixed_it;
-
- // Proceed with all registered callback functions.
- // Skip cells with the CELL_INVALID flag.
- if (cell_status !=
- parallel::distributed::Triangulation<dim,
- spacedim>::CELL_INVALID)
- {
- // Pack fixed size data.
- for (auto callback_it = pack_callbacks_fixed.cbegin();
- callback_it != pack_callbacks_fixed.cend();
- ++callback_it, ++data_fixed_it)
- {
- *data_fixed_it = (*callback_it)(dealii_cell, cell_status);
- }
-
- // Pack variable size data.
- // If we store variable size data, we need to transfer
- // the sizes of each corresponding callback function
- // via fixed size transfer as well.
- if (variable_size_data_stored)
- {
- const unsigned int n_variable_size_data_sets_on_cell =
- ((cell_status ==
- parallel::distributed::Triangulation<dim, spacedim>::
- CELL_INVALID) ?
- 0 :
- n_callbacks_variable);
- data_cell_variable_it->resize(
- n_variable_size_data_sets_on_cell);
-
- auto callback_it = pack_callbacks_variable.cbegin();
- auto data_variable_it = data_cell_variable_it->begin();
- auto sizes_variable_it =
- cell_sizes_variable_cumulative.begin();
- for (; callback_it != pack_callbacks_variable.cend();
- ++callback_it, ++data_variable_it, ++sizes_variable_it)
- {
- *data_variable_it =
- (*callback_it)(dealii_cell, cell_status);
-
- // Store data sizes for each callback function first.
- // Make it cumulative below.
- *sizes_variable_it = data_variable_it->size();
- }
-
- // Turn size vector into its cumulative representation.
- std::partial_sum(cell_sizes_variable_cumulative.begin(),
- cell_sizes_variable_cumulative.end(),
- cell_sizes_variable_cumulative.begin());
-
- // Serialize cumulative variable size vector value-by-value.
- // This way we can circumvent the overhead of storing the
- // container object as a whole, since we know its size by
- // the number of registered callback functions.
- data_fixed_it->resize(n_callbacks_variable *
- sizeof(unsigned int));
- for (unsigned int i = 0; i < n_callbacks_variable; ++i)
- std::memcpy(&(data_fixed_it->at(i *
- sizeof(unsigned int))),
- &(cell_sizes_variable_cumulative.at(i)),
- sizeof(unsigned int));
-
- ++data_fixed_it;
- }
-
- // Double check that we packed everything we wanted
- // in the fixed size buffers.
- Assert(data_fixed_it == data_cell_fixed_it->end(),
- ExcInternalError());
- }
-
- // Increment the variable size data iterator
- // only if we actually pack this kind of data
- // to avoid getting out of bounds.
- if (variable_size_data_stored)
- ++data_cell_variable_it;
- } // loop over quad_cell_relations
- }
-
- //
- // ----------- Gather data sizes for fixed size transfer ------------
- //
- // Generate a vector which stores the sizes of each callback function,
- // including the packed CellStatus transfer.
- // Find the very first cell that we wrote to with all callback
- // functions (i.e. a cell that was not flagged with CELL_INVALID)
- // and store the sizes of each buffer.
- //
- // To deal with the case that at least one of the processors does not own
- // any cell at all, we will exchange the information about the data sizes
- // among them later. The code in between is still well-defined, since the
- // following loops will be skipped.
- std::vector<unsigned int> local_sizes_fixed(
- 1 + n_callbacks_fixed + (variable_size_data_stored ? 1 : 0));
- for (const auto &data_cell : packed_fixed_size_data)
- {
- if (data_cell.size() == local_sizes_fixed.size())
- {
- auto sizes_fixed_it = local_sizes_fixed.begin();
- auto data_fixed_it = data_cell.cbegin();
- for (; data_fixed_it != data_cell.cend();
- ++data_fixed_it, ++sizes_fixed_it)
- {
- *sizes_fixed_it = data_fixed_it->size();
- }
-
- break;
- }
- }
-
- // Check if all cells have valid sizes.
- for (auto data_cell_fixed_it = packed_fixed_size_data.cbegin();
- data_cell_fixed_it != packed_fixed_size_data.cend();
- ++data_cell_fixed_it)
- {
- Assert((data_cell_fixed_it->size() == 1) ||
- (data_cell_fixed_it->size() == local_sizes_fixed.size()),
- ExcInternalError());
- }
-
- // Share information about the packed data sizes
- // of all callback functions across all processors, in case one
- // of them does not own any cells at all.
- std::vector<unsigned int> global_sizes_fixed(local_sizes_fixed.size());
- Utilities::MPI::max(local_sizes_fixed,
- this->mpi_communicator,
- global_sizes_fixed);
-
- // Construct cumulative sizes, since this is the only information
- // we need from now on.
- sizes_fixed_cumulative.resize(global_sizes_fixed.size());
- std::partial_sum(global_sizes_fixed.begin(),
- global_sizes_fixed.end(),
- sizes_fixed_cumulative.begin());
-
- //
- // ---------- Gather data sizes for variable size transfer ----------
- //
- if (variable_size_data_stored)
- {
- src_sizes_variable.reserve(packed_variable_size_data.size());
- for (const auto &data_cell : packed_variable_size_data)
- {
- int variable_data_size_on_cell = 0;
-
- for (const auto &data : data_cell)
- variable_data_size_on_cell += data.size();
-
- src_sizes_variable.push_back(variable_data_size_on_cell);
- }
- }
-
- //
- // ------------------------ Build buffers ---------------------------
- //
- const unsigned int expected_size_fixed =
- quad_cell_relations.size() * sizes_fixed_cumulative.back();
- const unsigned int expected_size_variable =
- std::accumulate(src_sizes_variable.begin(),
- src_sizes_variable.end(),
- std::vector<int>::size_type(0));
-
- // Move every piece of packed fixed size data into the consecutive buffer.
- src_data_fixed.reserve(expected_size_fixed);
- for (const auto &data_cell_fixed : packed_fixed_size_data)
- {
- // Move every fraction of packed data into the buffer
- // reserved for this particular cell.
- for (const auto &data_fixed : data_cell_fixed)
- std::move(data_fixed.begin(),
- data_fixed.end(),
- std::back_inserter(src_data_fixed));
-
- // If we only packed the CellStatus information
- // (i.e. encountered a cell flagged CELL_INVALID),
- // fill the remaining space with invalid entries.
- // We can skip this if there is nothing else to pack.
- if ((data_cell_fixed.size() == 1) &&
- (sizes_fixed_cumulative.size() > 1))
- {
- const std::size_t bytes_skipped =
- sizes_fixed_cumulative.back() - sizes_fixed_cumulative.front();
-
- src_data_fixed.insert(src_data_fixed.end(),
- bytes_skipped,
- static_cast<char>(-1)); // invalid_char
- }
- }
-
- // Move every piece of packed variable size data into the consecutive
- // buffer.
- if (variable_size_data_stored)
- {
- src_data_variable.reserve(expected_size_variable);
- for (const auto &data_cell : packed_variable_size_data)
- {
- // Move every fraction of packed data into the buffer
- // reserved for this particular cell.
- for (const auto &data : data_cell)
- std::move(data.begin(),
- data.end(),
- std::back_inserter(src_data_variable));
- }
- }
-
- // Double check that we packed everything correctly.
- Assert(src_data_fixed.size() == expected_size_fixed, ExcInternalError());
- Assert(src_data_variable.size() == expected_size_variable,
- ExcInternalError());
- }
-
-
-
- template <int dim, int spacedim>
- void
- Triangulation<dim, spacedim>::DataTransfer::execute_transfer(
- const typename dealii::internal::p4est::types<dim>::forest
- *parallel_forest,
- const typename dealii::internal::p4est::types<dim>::gloidx
- *previous_global_first_quadrant)
- {
- Assert(sizes_fixed_cumulative.size() > 0,
- ExcMessage("No data has been packed!"));
-
- // Resize memory according to the data that we will receive.
- dest_data_fixed.resize(parallel_forest->local_num_quadrants *
- sizes_fixed_cumulative.back());
-
- // Execute non-blocking fixed size transfer.
- typename dealii::internal::p4est::types<dim>::transfer_context
- *tf_context;
- tf_context =
- dealii::internal::p4est::functions<dim>::transfer_fixed_begin(
- parallel_forest->global_first_quadrant,
- previous_global_first_quadrant,
- parallel_forest->mpicomm,
- 0,
- dest_data_fixed.data(),
- src_data_fixed.data(),
- sizes_fixed_cumulative.back());
-
- if (variable_size_data_stored)
- {
- // Resize memory according to the data that we will receive.
- dest_sizes_variable.resize(parallel_forest->local_num_quadrants);
-
- // Execute fixed size transfer of data sizes for variable size
- // transfer.
- dealii::internal::p4est::functions<dim>::transfer_fixed(
- parallel_forest->global_first_quadrant,
- previous_global_first_quadrant,
- parallel_forest->mpicomm,
- 1,
- dest_sizes_variable.data(),
- src_sizes_variable.data(),
- sizeof(int));
- }
-
- dealii::internal::p4est::functions<dim>::transfer_fixed_end(tf_context);
-
- // Release memory of previously packed data.
- src_data_fixed.clear();
- src_data_fixed.shrink_to_fit();
-
- if (variable_size_data_stored)
- {
- // Resize memory according to the data that we will receive.
- dest_data_variable.resize(
- std::accumulate(dest_sizes_variable.begin(),
- dest_sizes_variable.end(),
- std::vector<int>::size_type(0)));
-
-# if DEAL_II_P4EST_VERSION_GTE(2, 0, 65, 0)
-# else
- // ----- WORKAROUND -----
- // An assertion in p4est prevents us from sending/receiving no data
- // at all, which is mandatory if one of our processes does not own
- // any quadrant. This bypasses the assertion from being triggered.
- // - see: https://github.com/cburstedde/p4est/issues/48
- if (src_sizes_variable.size() == 0)
- src_sizes_variable.resize(1);
- if (dest_sizes_variable.size() == 0)
- dest_sizes_variable.resize(1);
-# endif
-
- // Execute variable size transfer.
- dealii::internal::p4est::functions<dim>::transfer_custom(
- parallel_forest->global_first_quadrant,
- previous_global_first_quadrant,
- parallel_forest->mpicomm,
- 1,
- dest_data_variable.data(),
- dest_sizes_variable.data(),
- src_data_variable.data(),
- src_sizes_variable.data());
-
- // Release memory of previously packed data.
- src_sizes_variable.clear();
- src_sizes_variable.shrink_to_fit();
- src_data_variable.clear();
- src_data_variable.shrink_to_fit();
- }
- }
-
-
-
- template <int dim, int spacedim>
- void
- Triangulation<dim, spacedim>::DataTransfer::unpack_cell_status(
- std::vector<quadrant_cell_relation_t> &quad_cell_relations) const
- {
- Assert(sizes_fixed_cumulative.size() > 0,
- ExcMessage("No data has been packed!"));
- if (quad_cell_relations.size() > 0)
- {
- Assert(dest_data_fixed.size() > 0,
- ExcMessage("No data has been received!"));
- }
-
- // Size of CellStatus object that will be unpacked on each cell.
- const unsigned int size = sizes_fixed_cumulative.front();
-
- // Iterate over all cells and overwrite the CellStatus
- // information from the transferred data.
- // Proceed buffer iterator position to next cell after
- // each iteration.
- auto quad_cell_rel_it = quad_cell_relations.begin();
- auto dest_fixed_it = dest_data_fixed.cbegin();
- for (; quad_cell_rel_it != quad_cell_relations.end();
- ++quad_cell_rel_it, dest_fixed_it += sizes_fixed_cumulative.back())
- {
- std::get<1>(*quad_cell_rel_it) = // cell_status
- Utilities::unpack<typename parallel::distributed::
- Triangulation<dim, spacedim>::CellStatus>(
- dest_fixed_it,
- dest_fixed_it + size,
- /*allow_compression=*/false);
- }
- }
-
-
-
- template <int dim, int spacedim>
- void
- Triangulation<dim, spacedim>::DataTransfer::unpack_data(
- const std::vector<quadrant_cell_relation_t> &quad_cell_relations,
- const unsigned int handle,
- const std::function<void(
- const typename dealii::Triangulation<dim, spacedim>::cell_iterator &,
- const typename dealii::Triangulation<dim, spacedim>::CellStatus &,
- const boost::iterator_range<std::vector<char>::const_iterator> &)>
- &unpack_callback) const
- {
- // We decode the handle returned by register_data_attach() back into
- // a format we can use. All even handles belong to those callback
- // functions which write/read variable size data, all odd handles interact
- // with fixed size buffers.
- const bool callback_variable_transfer = (handle % 2 == 0);
- const unsigned int callback_index = handle / 2;
-
- // Cells will always receive fixed size data (i.e., CellStatus
- // information), but not necessarily variable size data (e.g., with a
- // ParticleHandler a cell might not contain any particle at all).
- // Thus it is sufficient to check if fixed size data has been received.
- Assert(sizes_fixed_cumulative.size() > 0,
- ExcMessage("No data has been packed!"));
- if (quad_cell_relations.size() > 0)
- {
- Assert(dest_data_fixed.size() > 0,
- ExcMessage("No data has been received!"));
- }
-
- std::vector<char>::const_iterator dest_data_it;
- std::vector<char>::const_iterator dest_sizes_cell_it;
-
- // Depending on whether our callback function unpacks fixed or
- // variable size data, we have to pursue different approaches
- // to localize the correct fraction of the buffer from which
- // we are allowed to read.
- unsigned int offset = numbers::invalid_unsigned_int;
- unsigned int size = numbers::invalid_unsigned_int;
- unsigned int data_increment = numbers::invalid_unsigned_int;
-
- if (callback_variable_transfer)
- {
- // For the variable size data, we need to extract the
- // data size from the fixed size buffer on each cell.
- //
- // We packed this information last, so the last packed
- // object in the fixed size buffer corresponds to the
- // variable data sizes.
- //
- // The last entry of sizes_fixed_cumulative corresponds
- // to the size of all fixed size data packed on the cell.
- // To get the offset for the last packed object, we need
- // to get the next-to-last entry.
- const unsigned int offset_variable_data_sizes =
- sizes_fixed_cumulative[sizes_fixed_cumulative.size() - 2];
-
- // This iterator points to the data size that the
- // callback_function packed for each specific cell.
- // Adjust buffer iterator to the offset of the callback
- // function so that we only have to advance its position
- // to the next cell after each iteration.
- dest_sizes_cell_it = dest_data_fixed.cbegin() +
- offset_variable_data_sizes +
- callback_index * sizeof(unsigned int);
-
- // Let the data iterator point to the correct buffer.
- dest_data_it = dest_data_variable.cbegin();
- }
- else
- {
- // For the fixed size data, we can get the information about
- // the buffer location on each cell directly from the
- // sizes_fixed_cumulative vector.
- offset = sizes_fixed_cumulative[callback_index];
- size = sizes_fixed_cumulative[callback_index + 1] - offset;
- data_increment = sizes_fixed_cumulative.back();
-
- // Let the data iterator point to the correct buffer.
- // Adjust buffer iterator to the offset of the callback
- // function so that we only have to advance its position
- // to the next cell after each iteration.
- dest_data_it = dest_data_fixed.cbegin() + offset;
- }
-
- // Iterate over all cells and unpack the transferred data.
- auto quad_cell_rel_it = quad_cell_relations.begin();
- auto dest_sizes_it = dest_sizes_variable.cbegin();
- for (; quad_cell_rel_it != quad_cell_relations.end();
- ++quad_cell_rel_it, dest_data_it += data_increment)
- {
- const auto &cell_status = std::get<1>(*quad_cell_rel_it);
- const auto &dealii_cell = std::get<2>(*quad_cell_rel_it);
-
- if (callback_variable_transfer)
- {
- // Update the increment according to the whole data size
- // of the current cell.
- data_increment = *dest_sizes_it;
-
- if (cell_status !=
- parallel::distributed::Triangulation<dim,
- spacedim>::CELL_INVALID)
- {
- // Extract the corresponding values for offset and size from
- // the cumulative sizes array stored in the fixed size buffer.
- if (callback_index == 0)
- offset = 0;
- else
- std::memcpy(&offset,
- &(*(dest_sizes_cell_it - sizeof(unsigned int))),
- sizeof(unsigned int));
-
- std::memcpy(&size,
- &(*dest_sizes_cell_it),
- sizeof(unsigned int));
-
- size -= offset;
-
- // Move the data iterator to the corresponding position
- // of the callback function and adjust the increment
- // accordingly.
- dest_data_it += offset;
- data_increment -= offset;
- }
-
- // Advance data size iterators to the next cell.
- dest_sizes_cell_it += sizes_fixed_cumulative.back();
- ++dest_sizes_it;
- }
-
- switch (cell_status)
- {
- case parallel::distributed::Triangulation<dim,
- spacedim>::CELL_PERSIST:
- case parallel::distributed::Triangulation<dim,
- spacedim>::CELL_COARSEN:
- unpack_callback(dealii_cell,
- cell_status,
- boost::make_iterator_range(dest_data_it,
- dest_data_it +
- size));
- break;
-
- case parallel::distributed::Triangulation<dim,
- spacedim>::CELL_REFINE:
- unpack_callback(dealii_cell->parent(),
- cell_status,
- boost::make_iterator_range(dest_data_it,
- dest_data_it +
- size));
- break;
-
- case parallel::distributed::Triangulation<dim,
- spacedim>::CELL_INVALID:
- // Skip this cell.
- break;
-
- default:
- Assert(false, ExcInternalError());
- break;
- }
- }
- }
-
-
-
- template <int dim, int spacedim>
- void
- Triangulation<dim, spacedim>::DataTransfer::save(
- const typename dealii::internal::p4est::types<dim>::forest
- * parallel_forest,
- const std::string &filename) const
- {
- // Large fractions of this function have been copied from
- // DataOutInterface::write_vtu_in_parallel.
- // TODO: Write general MPIIO interface.
-
- Assert(sizes_fixed_cumulative.size() > 0,
- ExcMessage("No data has been packed!"));
-
- const int myrank = Utilities::MPI::this_mpi_process(mpi_communicator);
-
- //
- // ---------- Fixed size data ----------
- //
- {
- const std::string fname_fixed = std::string(filename) + "_fixed.data";
-
- MPI_Info info;
- int ierr = MPI_Info_create(&info);
- AssertThrowMPI(ierr);
-
- MPI_File fh;
- ierr = MPI_File_open(mpi_communicator,
- DEAL_II_MPI_CONST_CAST(fname_fixed.c_str()),
- MPI_MODE_CREATE | MPI_MODE_WRONLY,
- info,
- &fh);
- AssertThrowMPI(ierr);
-
- ierr = MPI_File_set_size(fh, 0); // delete the file contents
- AssertThrowMPI(ierr);
- // this barrier is necessary, because otherwise others might already
- // write while one core is still setting the size to zero.
- ierr = MPI_Barrier(mpi_communicator);
- AssertThrowMPI(ierr);
- ierr = MPI_Info_free(&info);
- AssertThrowMPI(ierr);
- // ------------------
-
- // Check if number of processors is lined up with p4est partitioning.
- Assert(myrank < parallel_forest->mpisize, ExcInternalError());
-
- // Write cumulative sizes to file.
- // Since each processor owns the same information about the data sizes,
- // it is sufficient to let only the first processor perform this task.
- if (myrank == 0)
- {
- const unsigned int *data = sizes_fixed_cumulative.data();
-
- ierr = MPI_File_write_at(fh,
- 0,
- DEAL_II_MPI_CONST_CAST(data),
- sizes_fixed_cumulative.size(),
- MPI_UNSIGNED,
- MPI_STATUS_IGNORE);
- AssertThrowMPI(ierr);
- }
-
- // Write packed data to file simultaneously.
- const unsigned int offset_fixed =
- sizes_fixed_cumulative.size() * sizeof(unsigned int);
-
- const char *data = src_data_fixed.data();
-
- ierr = MPI_File_write_at(
- fh,
- offset_fixed +
- parallel_forest->global_first_quadrant[myrank] *
- sizes_fixed_cumulative.back(), // global position in file
- DEAL_II_MPI_CONST_CAST(data),
- src_data_fixed.size(), // local buffer
- MPI_CHAR,
- MPI_STATUS_IGNORE);
- AssertThrowMPI(ierr);
-
- ierr = MPI_File_close(&fh);
- AssertThrowMPI(ierr);
- }
-
- //
- // ---------- Variable size data ----------
- //
- if (variable_size_data_stored)
- {
- const std::string fname_variable =
- std::string(filename) + "_variable.data";
-
- MPI_Info info;
- int ierr = MPI_Info_create(&info);
- AssertThrowMPI(ierr);
-
- MPI_File fh;
- ierr = MPI_File_open(mpi_communicator,
- DEAL_II_MPI_CONST_CAST(fname_variable.c_str()),
- MPI_MODE_CREATE | MPI_MODE_WRONLY,
- info,
- &fh);
- AssertThrowMPI(ierr);
-
- ierr = MPI_File_set_size(fh, 0); // delete the file contents
- AssertThrowMPI(ierr);
- // this barrier is necessary, because otherwise others might already
- // write while one core is still setting the size to zero.
- ierr = MPI_Barrier(mpi_communicator);
- AssertThrowMPI(ierr);
- ierr = MPI_Info_free(&info);
- AssertThrowMPI(ierr);
-
- // Write sizes of each cell into file simultaneously.
- {
- const int *data = src_sizes_variable.data();
- ierr =
- MPI_File_write_at(fh,
- parallel_forest->global_first_quadrant[myrank] *
- sizeof(int), // global position in file
- DEAL_II_MPI_CONST_CAST(data),
- src_sizes_variable.size(), // local buffer
- MPI_INT,
- MPI_STATUS_IGNORE);
- AssertThrowMPI(ierr);
- }
-
-
- const unsigned int offset_variable =
- parallel_forest->global_num_quadrants * sizeof(int);
-
- // Gather size of data in bytes we want to store from this processor.
- const unsigned int size_on_proc = src_data_variable.size();
-
- // Compute prefix sum
- unsigned int prefix_sum = 0;
- ierr = MPI_Exscan(DEAL_II_MPI_CONST_CAST(&size_on_proc),
- &prefix_sum,
- 1,
- MPI_UNSIGNED,
- MPI_SUM,
- mpi_communicator);
- AssertThrowMPI(ierr);
-
- const char *data = src_data_variable.data();
-
- // Write data consecutively into file.
- ierr = MPI_File_write_at(fh,
- offset_variable +
- prefix_sum, // global position in file
- DEAL_II_MPI_CONST_CAST(data),
- src_data_variable.size(), // local buffer
- MPI_CHAR,
- MPI_STATUS_IGNORE);
- AssertThrowMPI(ierr);
-
- ierr = MPI_File_close(&fh);
- AssertThrowMPI(ierr);
- }
- }
-
-
-
- template <int dim, int spacedim>
- void
- Triangulation<dim, spacedim>::DataTransfer::load(
- const typename dealii::internal::p4est::types<dim>::forest
- * parallel_forest,
- const std::string &filename,
- const unsigned int n_attached_deserialize_fixed,
- const unsigned int n_attached_deserialize_variable)
- {
- // Large fractions of this function have been copied from
- // DataOutInterface::write_vtu_in_parallel.
- // TODO: Write general MPIIO interface.
-
- Assert(dest_data_fixed.size() == 0,
- ExcMessage("Previously loaded data has not been released yet!"));
-
- variable_size_data_stored = (n_attached_deserialize_variable > 0);
-
- const int myrank = Utilities::MPI::this_mpi_process(mpi_communicator);
-
- //
- // ---------- Fixed size data ----------
- //
- {
- const std::string fname_fixed = std::string(filename) + "_fixed.data";
-
- MPI_Info info;
- int ierr = MPI_Info_create(&info);
- AssertThrowMPI(ierr);
-
- MPI_File fh;
- ierr = MPI_File_open(mpi_communicator,
- DEAL_II_MPI_CONST_CAST(fname_fixed.c_str()),
- MPI_MODE_RDONLY,
- info,
- &fh);
- AssertThrowMPI(ierr);
-
- ierr = MPI_Info_free(&info);
- AssertThrowMPI(ierr);
-
- // Check if number of processors is lined up with p4est partitioning.
- Assert(myrank < parallel_forest->mpisize, ExcInternalError());
-
- // Read cumulative sizes from file.
- // Since all processors need the same information about the data sizes,
- // let each of them retrieve it by reading from the same location in
- // the file.
- sizes_fixed_cumulative.resize(1 + n_attached_deserialize_fixed +
- (variable_size_data_stored ? 1 : 0));
- ierr = MPI_File_read_at(fh,
- 0,
- sizes_fixed_cumulative.data(),
- sizes_fixed_cumulative.size(),
- MPI_UNSIGNED,
- MPI_STATUS_IGNORE);
- AssertThrowMPI(ierr);
-
- // Allocate sufficient memory.
- dest_data_fixed.resize(parallel_forest->local_num_quadrants *
- sizes_fixed_cumulative.back());
-
- // Read packed data from file simultaneously.
- const unsigned int offset =
- sizes_fixed_cumulative.size() * sizeof(unsigned int);
-
- ierr = MPI_File_read_at(
- fh,
- offset + parallel_forest->global_first_quadrant[myrank] *
- sizes_fixed_cumulative.back(), // global position in file
- dest_data_fixed.data(),
- dest_data_fixed.size(), // local buffer
- MPI_CHAR,
- MPI_STATUS_IGNORE);
- AssertThrowMPI(ierr);
-
- ierr = MPI_File_close(&fh);
- AssertThrowMPI(ierr);
- }
-
- //
- // ---------- Variable size data ----------
- //
- if (variable_size_data_stored)
- {
- const std::string fname_variable =
- std::string(filename) + "_variable.data";
-
- MPI_Info info;
- int ierr = MPI_Info_create(&info);
- AssertThrowMPI(ierr);
-
- MPI_File fh;
- ierr = MPI_File_open(mpi_communicator,
- DEAL_II_MPI_CONST_CAST(fname_variable.c_str()),
- MPI_MODE_RDONLY,
- info,
- &fh);
- AssertThrowMPI(ierr);
-
- ierr = MPI_Info_free(&info);
- AssertThrowMPI(ierr);
-
- // Read sizes of all locally owned cells.
- dest_sizes_variable.resize(parallel_forest->local_num_quadrants);
- ierr =
- MPI_File_read_at(fh,
- parallel_forest->global_first_quadrant[myrank] *
- sizeof(int),
- dest_sizes_variable.data(),
- dest_sizes_variable.size(),
- MPI_INT,
- MPI_STATUS_IGNORE);
- AssertThrowMPI(ierr);
-
- const unsigned int offset =
- parallel_forest->global_num_quadrants * sizeof(int);
-
- const unsigned int size_on_proc =
- std::accumulate(dest_sizes_variable.begin(),
- dest_sizes_variable.end(),
- 0);
-
- // share information among all processors by prefix sum
- unsigned int prefix_sum = 0;
- ierr = MPI_Exscan(DEAL_II_MPI_CONST_CAST(&size_on_proc),
- &prefix_sum,
- 1,
- MPI_UNSIGNED,
- MPI_SUM,
- mpi_communicator);
- AssertThrowMPI(ierr);
-
- dest_data_variable.resize(size_on_proc);
- ierr = MPI_File_read_at(fh,
- offset + prefix_sum,
- dest_data_variable.data(),
- dest_data_variable.size(),
- MPI_CHAR,
- MPI_STATUS_IGNORE);
- AssertThrowMPI(ierr);
-
- ierr = MPI_File_close(&fh);
- AssertThrowMPI(ierr);
- }
- }
-
-
-
- template <int dim, int spacedim>
- void
- Triangulation<dim, spacedim>::DataTransfer::clear()
- {
- variable_size_data_stored = false;
-
- // free information about data sizes
- sizes_fixed_cumulative.clear();
- sizes_fixed_cumulative.shrink_to_fit();
-
- // free fixed size transfer data
- src_data_fixed.clear();
- src_data_fixed.shrink_to_fit();
-
- dest_data_fixed.clear();
- dest_data_fixed.shrink_to_fit();
-
- // free variable size transfer data
- src_sizes_variable.clear();
- src_sizes_variable.shrink_to_fit();
-
- src_data_variable.clear();
- src_data_variable.shrink_to_fit();
-
- dest_sizes_variable.clear();
- dest_sizes_variable.shrink_to_fit();
-
- dest_data_variable.clear();
- dest_data_variable.shrink_to_fit();
- }
-
-
-
- /* ----------------- class Triangulation<dim,spacedim> ----------------- */
-
-
+ /*----------------- class Triangulation<dim,spacedim> ---------------\*/
template <int dim, int spacedim>
Triangulation<dim, spacedim>::Triangulation(
const MPI_Comm &mpi_communicator,
, triangulation_has_content(false)
, connectivity(nullptr)
, parallel_forest(nullptr)
- , cell_attached_data({0, 0, {}, {}})
- , data_transfer(mpi_communicator)
{
parallel_ghost = nullptr;
}
{
triangulation_has_content = false;
- cell_attached_data = {0, 0, {}, {}};
- data_transfer.clear();
-
if (parallel_ghost != nullptr)
{
dealii::internal::p4est::functions<dim>::ghost_destroy(
coarse_cell_to_p4est_tree_permutation.resize(0);
p4est_tree_to_coarse_cell_permutation.resize(0);
- dealii::Triangulation<dim, spacedim>::clear();
+ dealii::parallel::DistributedTriangulationBase<dim, spacedim>::clear();
this->update_number_cache();
}
+ template <int dim, int spacedim>
+ void
+ Triangulation<dim, spacedim>::execute_transfer(
+ const typename dealii::internal::p4est::types<dim>::forest
+ *parallel_forest,
+ const typename dealii::internal::p4est::types<dim>::gloidx
+ *previous_global_first_quadrant)
+ {
+ Assert(this->data_transfer.sizes_fixed_cumulative.size() > 0,
+ ExcMessage("No data has been packed!"));
+
+ // Resize memory according to the data that we will receive.
+ this->data_transfer.dest_data_fixed.resize(
+ parallel_forest->local_num_quadrants *
+ this->data_transfer.sizes_fixed_cumulative.back());
+
+ // Execute non-blocking fixed size transfer.
+ typename dealii::internal::p4est::types<dim>::transfer_context
+ *tf_context;
+ tf_context =
+ dealii::internal::p4est::functions<dim>::transfer_fixed_begin(
+ parallel_forest->global_first_quadrant,
+ previous_global_first_quadrant,
+ parallel_forest->mpicomm,
+ 0,
+ this->data_transfer.dest_data_fixed.data(),
+ this->data_transfer.src_data_fixed.data(),
+ this->data_transfer.sizes_fixed_cumulative.back());
+
+ if (this->data_transfer.variable_size_data_stored)
+ {
+ // Resize memory according to the data that we will receive.
+ this->data_transfer.dest_sizes_variable.resize(
+ parallel_forest->local_num_quadrants);
+
+ // Execute fixed size transfer of data sizes for variable size
+ // transfer.
+ dealii::internal::p4est::functions<dim>::transfer_fixed(
+ parallel_forest->global_first_quadrant,
+ previous_global_first_quadrant,
+ parallel_forest->mpicomm,
+ 1,
+ this->data_transfer.dest_sizes_variable.data(),
+ this->data_transfer.src_sizes_variable.data(),
+ sizeof(unsigned int));
+ }
+
+ dealii::internal::p4est::functions<dim>::transfer_fixed_end(tf_context);
+
+ // Release memory of previously packed data.
+ this->data_transfer.src_data_fixed.clear();
+ this->data_transfer.src_data_fixed.shrink_to_fit();
+
+ if (this->data_transfer.variable_size_data_stored)
+ {
+ // Resize memory according to the data that we will receive.
+ this->data_transfer.dest_data_variable.resize(
+ std::accumulate(this->data_transfer.dest_sizes_variable.begin(),
+ this->data_transfer.dest_sizes_variable.end(),
+ std::vector<int>::size_type(0)));
+
+# if DEAL_II_P4EST_VERSION_GTE(2, 0, 65, 0)
+# else
+ // ----- WORKAROUND -----
+ // An assertion in p4est prevents us from sending/receiving no data
+ // at all, which is mandatory if one of our processes does not own
+ // any quadrant. This bypasses the assertion from being triggered.
+ // - see: https://github.com/cburstedde/p4est/issues/48
+ if (this->data_transfer.src_sizes_variable.size() == 0)
+ this->data_transfer.src_sizes_variable.resize(1);
+ if (this->data_transfer.dest_sizes_variable.size() == 0)
+ this->data_transfer.dest_sizes_variable.resize(1);
+# endif
+
+ // Execute variable size transfer.
+ dealii::internal::p4est::functions<dim>::transfer_custom(
+ parallel_forest->global_first_quadrant,
+ previous_global_first_quadrant,
+ parallel_forest->mpicomm,
+ 1,
+ this->data_transfer.dest_data_variable.data(),
+ this->data_transfer.dest_sizes_variable.data(),
+ this->data_transfer.src_data_variable.data(),
+ this->data_transfer.src_sizes_variable.data());
+
+ // Release memory of previously packed data.
+ this->data_transfer.src_sizes_variable.clear();
+ this->data_transfer.src_sizes_variable.shrink_to_fit();
+ this->data_transfer.src_data_variable.clear();
+ this->data_transfer.src_data_variable.shrink_to_fit();
+ }
+ }
+
+
template <int dim, int spacedim>
bool
Triangulation<dim, spacedim>::has_hanging_nodes() const
Triangulation<dim, spacedim>::save(const std::string &filename) const
{
Assert(
- cell_attached_data.n_attached_deserialize == 0,
+ this->cell_attached_data.n_attached_deserialize == 0,
ExcMessage(
"not all SolutionTransfer's got deserialized after the last load()"));
Assert(this->n_cells() > 0,
ExcMessage("Can not save() an empty Triangulation."));
+ const int myrank =
+ Utilities::MPI::this_mpi_process(this->mpi_communicator);
+
// signal that serialization is going to happen
this->signals.pre_distributed_save();
<< std::endl
<< 4 << " "
<< Utilities::MPI::n_mpi_processes(this->mpi_communicator) << " "
- << cell_attached_data.pack_callbacks_fixed.size() << " "
- << cell_attached_data.pack_callbacks_variable.size() << " "
+ << this->cell_attached_data.pack_callbacks_fixed.size() << " "
+ << this->cell_attached_data.pack_callbacks_variable.size() << " "
<< this->n_cells(0) << std::endl;
}
// each cell should have been flagged `CELL_PERSIST`
- for (const auto &quad_cell_rel : local_quadrant_cell_relations)
+ for (const auto &cell_rel : this->local_cell_relations)
{
- (void)quad_cell_rel;
+ (void)cell_rel;
Assert(
- (std::get<1>(quad_cell_rel) == // cell_status
+ (std::get<0>(cell_rel) == // cell_status
parallel::distributed::Triangulation<dim, spacedim>::CELL_PERSIST),
ExcInternalError());
}
- if (cell_attached_data.n_attached_data_sets > 0)
- {
- // cast away constness
- auto tria = const_cast<
- dealii::parallel::distributed::Triangulation<dim, spacedim> *>(
- this);
-
- // pack attached data first
- tria->data_transfer.pack_data(
- local_quadrant_cell_relations,
- cell_attached_data.pack_callbacks_fixed,
- cell_attached_data.pack_callbacks_variable);
-
- // then store buffers in file
- tria->data_transfer.save(parallel_forest, filename);
-
- // and release the memory afterwards
- tria->data_transfer.clear();
- }
+ // Save cell attached data.
+ this->save_attached_data(parallel_forest->global_first_quadrant[myrank],
+ parallel_forest->global_num_quadrants,
+ filename);
dealii::internal::p4est::functions<dim>::save(filename.c_str(),
parallel_forest,
false);
- // clear all of the callback data, as explained in the documentation of
- // register_data_attach()
- {
- dealii::parallel::distributed::Triangulation<dim, spacedim> *tria =
- const_cast<
- dealii::parallel::distributed::Triangulation<dim, spacedim> *>(
- this);
-
- tria->cell_attached_data.n_attached_data_sets = 0;
- tria->cell_attached_data.pack_callbacks_fixed.clear();
- tria->cell_attached_data.pack_callbacks_variable.clear();
- }
-
// signal that serialization has finished
this->signals.post_distributed_save();
}
ExcMessage(
"Triangulation may only contain coarse cells when calling load()."));
+ const int myrank =
+ Utilities::MPI::this_mpi_process(this->mpi_communicator);
+
// signal that de-serialization is going to happen
this->signals.pre_distributed_load();
// clear all of the callback data, as explained in the documentation of
// register_data_attach()
- cell_attached_data.n_attached_data_sets = 0;
- cell_attached_data.n_attached_deserialize =
+ this->cell_attached_data.n_attached_data_sets = 0;
+ this->cell_attached_data.n_attached_deserialize =
attached_count_fixed + attached_count_variable;
parallel_forest = dealii::internal::p4est::functions<dim>::load_ext(
Assert(false, ExcInternalError());
}
- // load saved data, if any was stored
- if (cell_attached_data.n_attached_deserialize > 0)
- {
- data_transfer.load(parallel_forest,
- filename,
- attached_count_fixed,
- attached_count_variable);
-
- data_transfer.unpack_cell_status(local_quadrant_cell_relations);
+ // Load attached cell data, if any was stored.
+ this->load_attached_data(parallel_forest->global_first_quadrant[myrank],
+ parallel_forest->global_num_quadrants,
+ parallel_forest->local_num_quadrants,
+ filename,
+ attached_count_fixed,
+ attached_count_variable);
- // the CellStatus of all stored cells should always be CELL_PERSIST.
- for (const auto &quad_cell_rel : local_quadrant_cell_relations)
- {
- (void)quad_cell_rel;
- Assert(
- (std::get<1>(quad_cell_rel) == // cell_status
- parallel::distributed::Triangulation<dim,
- spacedim>::CELL_PERSIST),
- ExcInternalError());
- }
- }
+ // signal that de-serialization is finished
+ this->signals.post_distributed_load();
this->update_periodic_face_map();
this->update_number_cache();
-
- // signal that de-serialization is finished
- this->signals.post_distributed_load();
}
this->smooth_grid = save_smooth;
// finally, after syncing the parallel_forest with the triangulation,
- // also update the quadrant_cell_relations, which will be used for
+ // also update the cell_relations, which will be used for
// repartitioning, further refinement/coarsening, and unpacking
// of stored or transferred data.
- update_quadrant_cell_relations();
+ update_cell_relations();
}
// since refinement and/or coarsening on the parallel forest
// has happened, we need to update the quadrant cell relations
- update_quadrant_cell_relations();
+ update_cell_relations();
// before repartitioning the mesh, store the current distribution
// of the p4est quadrants and let others attach mesh related info
previous_global_first_quadrant;
// pack data only if anything has been attached
- if (cell_attached_data.n_attached_data_sets > 0)
+ if (this->cell_attached_data.n_attached_data_sets > 0)
{
- data_transfer.pack_data(local_quadrant_cell_relations,
- cell_attached_data.pack_callbacks_fixed,
- cell_attached_data.pack_callbacks_variable);
+ this->data_transfer.pack_data(
+ this->local_cell_relations,
+ this->cell_attached_data.pack_callbacks_fixed,
+ this->cell_attached_data.pack_callbacks_variable);
// before repartitioning the p4est object, save a copy of the
// positions of the global first quadrants for data transfer later
// transfer data
// only if anything has been attached
- if (cell_attached_data.n_attached_data_sets > 0)
+ if (this->cell_attached_data.n_attached_data_sets > 0)
{
// execute transfer after triangulation got updated
- data_transfer.execute_transfer(parallel_forest,
- previous_global_first_quadrant.data());
+ this->execute_transfer(parallel_forest,
+ previous_global_first_quadrant.data());
// also update the CellStatus information on the new mesh
- data_transfer.unpack_cell_status(local_quadrant_cell_relations);
+ this->data_transfer.unpack_cell_status(this->local_cell_relations);
}
# ifdef DEBUG
previous_global_first_quadrant;
// pack data only if anything has been attached
- if (cell_attached_data.n_attached_data_sets > 0)
+ if (this->cell_attached_data.n_attached_data_sets > 0)
{
- data_transfer.pack_data(local_quadrant_cell_relations,
- cell_attached_data.pack_callbacks_fixed,
- cell_attached_data.pack_callbacks_variable);
+ this->data_transfer.pack_data(
+ this->local_cell_relations,
+ this->cell_attached_data.pack_callbacks_fixed,
+ this->cell_attached_data.pack_callbacks_variable);
// before repartitioning the p4est object, save a copy of the
// positions of quadrant for data transfer later
// transfer data
// only if anything has been attached
- if (cell_attached_data.n_attached_data_sets > 0)
+ if (this->cell_attached_data.n_attached_data_sets > 0)
{
// execute transfer after triangulation got updated
- data_transfer.execute_transfer(parallel_forest,
- previous_global_first_quadrant.data());
+ this->execute_transfer(parallel_forest,
+ previous_global_first_quadrant.data());
}
this->update_periodic_face_map();
- template <int dim, int spacedim>
- unsigned int
- Triangulation<dim, spacedim>::register_data_attach(
- const std::function<std::vector<char>(const cell_iterator &,
- const CellStatus)> &pack_callback,
- const bool returns_variable_size_data)
- {
- unsigned int handle = numbers::invalid_unsigned_int;
-
- // Add new callback function to the corresponding register.
- // Encode handles according to returns_variable_size_data.
- if (returns_variable_size_data)
- {
- handle = 2 * cell_attached_data.pack_callbacks_variable.size();
- cell_attached_data.pack_callbacks_variable.push_back(pack_callback);
- }
- else
- {
- handle = 2 * cell_attached_data.pack_callbacks_fixed.size() + 1;
- cell_attached_data.pack_callbacks_fixed.push_back(pack_callback);
- }
-
- // Increase overall counter.
- ++cell_attached_data.n_attached_data_sets;
-
- return handle;
- }
-
-
-
- template <int dim, int spacedim>
- void
- Triangulation<dim, spacedim>::notify_ready_to_unpack(
- const unsigned int handle,
- const std::function<
- void(const cell_iterator &,
- const CellStatus,
- const boost::iterator_range<std::vector<char>::const_iterator> &)>
- &unpack_callback)
- {
- Assert(cell_attached_data.n_attached_data_sets > 0,
- ExcMessage("The notify_ready_to_unpack() has already been called "
- "once for each registered callback."));
-
- // check if local_quadrant_cell_relations have been previously gathered
- // correctly
- Assert(local_quadrant_cell_relations.size() ==
- static_cast<unsigned int>(parallel_forest->local_num_quadrants),
- ExcInternalError());
-
-# ifdef DEBUG
- // check validity of handle and deregister pack_callback function.
- // first reset with invalid entries to preserve ambiguity of
- // handles, then free memory when all were unpacked (see below).
- const unsigned int callback_index = handle / 2;
- if (handle % 2 == 0)
- {
- Assert(callback_index <
- cell_attached_data.pack_callbacks_variable.size(),
- ExcMessage("Invalid handle."));
-
- Assert(cell_attached_data.pack_callbacks_variable[callback_index] !=
- nullptr,
- ExcInternalError());
- cell_attached_data.pack_callbacks_variable[callback_index] = nullptr;
- }
- else
- {
- Assert(callback_index <
- cell_attached_data.pack_callbacks_fixed.size(),
- ExcMessage("Invalid handle."));
-
- Assert(cell_attached_data.pack_callbacks_fixed[callback_index] !=
- nullptr,
- ExcInternalError());
- cell_attached_data.pack_callbacks_fixed[callback_index] = nullptr;
- }
-# endif
-
- // perform unpacking
- data_transfer.unpack_data(local_quadrant_cell_relations,
- handle,
- unpack_callback);
-
- // decrease counters
- --cell_attached_data.n_attached_data_sets;
- if (cell_attached_data.n_attached_deserialize > 0)
- --cell_attached_data.n_attached_deserialize;
-
- // important: only remove data if we are not in the deserialization
- // process. There, each SolutionTransfer registers and unpacks before
- // the next one does this, so n_attached_data_sets is only 1 here. This
- // would destroy the saved data before the second SolutionTransfer can
- // get it. This created a bug that is documented in
- // tests/mpi/p4est_save_03 with more than one SolutionTransfer.
- if (cell_attached_data.n_attached_data_sets == 0 &&
- cell_attached_data.n_attached_deserialize == 0)
- {
- // everybody got their data, time for cleanup!
- cell_attached_data.pack_callbacks_fixed.clear();
- cell_attached_data.pack_callbacks_variable.clear();
- data_transfer.clear();
-
- // reset all cell_status entries after coarsening/refinement
- for (auto &quad_cell_rel : local_quadrant_cell_relations)
- std::get<1>(quad_cell_rel) =
- parallel::distributed::Triangulation<dim, spacedim>::CELL_PERSIST;
- }
- }
-
-
-
template <int dim, int spacedim>
const std::vector<types::global_dof_index> &
Triangulation<dim, spacedim>::get_p4est_tree_to_coarse_cell_permutation()
MemoryConsumption::memory_consumption(connectivity) +
MemoryConsumption::memory_consumption(parallel_forest) +
MemoryConsumption::memory_consumption(
- cell_attached_data.n_attached_data_sets) +
+ this->cell_attached_data.n_attached_data_sets) +
// MemoryConsumption::memory_consumption(cell_attached_data.pack_callbacks_fixed)
// +
// MemoryConsumption::memory_consumption(cell_attached_data.pack_callbacks_variable)
other_tria_x->coarse_cell_to_p4est_tree_permutation;
p4est_tree_to_coarse_cell_permutation =
other_tria_x->p4est_tree_to_coarse_cell_permutation;
- cell_attached_data = other_tria_x->cell_attached_data;
- data_transfer = other_tria_x->data_transfer;
+ this->cell_attached_data = other_tria_x->cell_attached_data;
+ this->data_transfer = other_tria_x->data_transfer;
settings = other_tria_x->settings;
}
template <int dim, int spacedim>
void
- Triangulation<dim, spacedim>::update_quadrant_cell_relations()
+ Triangulation<dim, spacedim>::update_cell_relations()
{
- // reorganize memory for local_quadrant_cell_relations
- local_quadrant_cell_relations.resize(
- parallel_forest->local_num_quadrants);
- local_quadrant_cell_relations.shrink_to_fit();
+ // reorganize memory for local_cell_relations
+ this->local_cell_relations.resize(parallel_forest->local_num_quadrants);
+ this->local_cell_relations.shrink_to_fit();
// recurse over p4est
for (typename Triangulation<dim, spacedim>::cell_iterator cell =
typename dealii::internal::p4est::types<dim>::tree *tree =
init_tree(cell->index());
- update_quadrant_cell_relations_recursively<dim, spacedim>(
- local_quadrant_cell_relations, *tree, cell, p4est_coarse_cell);
+ update_cell_relations_recursively<dim, spacedim>(
+ this->local_cell_relations, *tree, cell, p4est_coarse_cell);
}
}
std::vector<unsigned int>
Triangulation<dim, spacedim>::get_cell_weights() const
{
- // check if local_quadrant_cell_relations have been previously gathered
+ // check if local_cell_relations have been previously gathered
// correctly
- Assert(local_quadrant_cell_relations.size() ==
+ Assert(this->local_cell_relations.size() ==
static_cast<unsigned int>(parallel_forest->local_num_quadrants),
ExcInternalError());
// Note that we need to follow the p4est ordering
// instead of the deal.II ordering to get the cell_weights
// in the same order p4est will encounter them during repartitioning.
- for (const auto &quad_cell_rel : local_quadrant_cell_relations)
+ for (const auto &cell_rel : this->local_cell_relations)
{
- const auto &cell_status = std::get<1>(quad_cell_rel);
- const auto &cell_it = std::get<2>(quad_cell_rel);
+ const auto &cell_status = std::get<0>(cell_rel);
+ const auto &cell_it = std::get<1>(cell_rel);
switch (cell_status)
{
- template <int spacedim>
- unsigned int
- Triangulation<1, spacedim>::register_data_attach(
- const std::function<std::vector<char>(
- const typename dealii::Triangulation<1, spacedim>::cell_iterator &,
- const typename dealii::Triangulation<1, spacedim>::CellStatus)>
- & /*pack_callback*/,
- const bool /*returns_variable_size_data*/)
- {
- Assert(false, ExcNotImplemented());
- return 0;
- }
-
-
-
- template <int spacedim>
- void
- Triangulation<1, spacedim>::notify_ready_to_unpack(
- const unsigned int /*handle*/,
- const std::function<
- void(const typename dealii::Triangulation<1, spacedim>::cell_iterator &,
- const typename dealii::Triangulation<1, spacedim>::CellStatus,
- const boost::iterator_range<std::vector<char>::const_iterator> &)>
- & /*unpack_callback*/)
- {
- Assert(false, ExcNotImplemented());
- }
-
-
-
template <int spacedim>
const std::vector<types::global_dof_index> &
Triangulation<1, spacedim>::get_p4est_tree_to_coarse_cell_permutation()
}
+
template <int spacedim>
void
Triangulation<1, spacedim>::load(const std::string &, const bool)
return false;
}
+
+
+ template <int spacedim>
+ void
+ Triangulation<1, spacedim>::update_cell_relations()
+ {
+ Assert(false, ExcNotImplemented());
+ }
+
} // namespace distributed
} // namespace parallel
void
TriangulationBase<dim, spacedim>::update_number_cache()
{
- Assert(false, ExcNotImplemented());
+ Assert(false, ExcNeedsMPI());
}
#endif
return number_cache.level_cell_index_partitioners[level];
}
+
+
template <int dim, int spacedim>
DistributedTriangulationBase<dim, spacedim>::DistributedTriangulationBase(
const MPI_Comm &mpi_communicator,
mpi_communicator,
smooth_grid,
check_for_distorted_cells)
+ , cell_attached_data({0, 0, {}, {}})
+ , data_transfer(mpi_communicator)
{}
+
+
+ template <int dim, int spacedim>
+ void
+ DistributedTriangulationBase<dim, spacedim>::clear()
+ {
+ cell_attached_data = {0, 0, {}, {}};
+ data_transfer.clear();
+
+ dealii::Triangulation<dim, spacedim>::clear();
+ }
+
+
+
+ template <int dim, int spacedim>
+ void
+ DistributedTriangulationBase<dim, spacedim>::save_attached_data(
+ const unsigned int global_first_cell,
+ const unsigned int global_num_cells,
+ const std::string &filename) const
+ {
+ // cast away constness
+ auto tria = const_cast<
+ dealii::parallel::DistributedTriangulationBase<dim, spacedim> *>(this);
+
+ if (this->cell_attached_data.n_attached_data_sets > 0)
+ {
+ // pack attached data first
+ tria->data_transfer.pack_data(
+ tria->local_cell_relations,
+ tria->cell_attached_data.pack_callbacks_fixed,
+ tria->cell_attached_data.pack_callbacks_variable);
+
+ // then store buffers in file
+ tria->data_transfer.save(global_first_cell, global_num_cells, filename);
+
+ // and release the memory afterwards
+ tria->data_transfer.clear();
+ }
+
+ // clear all of the callback data, as explained in the documentation of
+ // register_data_attach()
+ {
+ tria->cell_attached_data.n_attached_data_sets = 0;
+ tria->cell_attached_data.pack_callbacks_fixed.clear();
+ tria->cell_attached_data.pack_callbacks_variable.clear();
+ }
+ }
+
+
+ template <int dim, int spacedim>
+ void
+ DistributedTriangulationBase<dim, spacedim>::load_attached_data(
+ const unsigned int global_first_cell,
+ const unsigned int global_num_cells,
+ const unsigned int local_num_cells,
+ const std::string &filename,
+ const unsigned int n_attached_deserialize_fixed,
+ const unsigned int n_attached_deserialize_variable)
+ {
+ // load saved data, if any was stored
+ if (this->cell_attached_data.n_attached_deserialize > 0)
+ {
+ this->data_transfer.load(global_first_cell,
+ global_num_cells,
+ local_num_cells,
+ filename,
+ n_attached_deserialize_fixed,
+ n_attached_deserialize_variable);
+
+ this->data_transfer.unpack_cell_status(this->local_cell_relations);
+
+ // the CellStatus of all stored cells should always be CELL_PERSIST.
+ for (const auto &cell_rel : this->local_cell_relations)
+ {
+ (void)cell_rel;
+ Assert(
+ (std::get<0>(cell_rel) == // cell_status
+ parallel::DistributedTriangulationBase<dim,
+ spacedim>::CELL_PERSIST),
+ ExcInternalError());
+ }
+ }
+ }
+
+ template <int dim, int spacedim>
+ unsigned int
+ DistributedTriangulationBase<dim, spacedim>::register_data_attach(
+ const std::function<std::vector<char>(const cell_iterator &,
+ const CellStatus)> &pack_callback,
+ const bool returns_variable_size_data)
+ {
+ unsigned int handle = numbers::invalid_unsigned_int;
+
+ // Add new callback function to the corresponding register.
+ // Encode handles according to returns_variable_size_data.
+ if (returns_variable_size_data)
+ {
+ handle = 2 * cell_attached_data.pack_callbacks_variable.size();
+ cell_attached_data.pack_callbacks_variable.push_back(pack_callback);
+ }
+ else
+ {
+ handle = 2 * cell_attached_data.pack_callbacks_fixed.size() + 1;
+ cell_attached_data.pack_callbacks_fixed.push_back(pack_callback);
+ }
+
+ // Increase overall counter.
+ ++cell_attached_data.n_attached_data_sets;
+
+ return handle;
+ }
+
+
+ template <int dim, int spacedim>
+ void
+ DistributedTriangulationBase<dim, spacedim>::notify_ready_to_unpack(
+ const unsigned int handle,
+ const std::function<
+ void(const cell_iterator &,
+ const CellStatus,
+ const boost::iterator_range<std::vector<char>::const_iterator> &)>
+ &unpack_callback)
+ {
+ // perform unpacking
+ data_transfer.unpack_data(local_cell_relations, handle, unpack_callback);
+
+ // decrease counters
+ --cell_attached_data.n_attached_data_sets;
+ if (cell_attached_data.n_attached_deserialize > 0)
+ --cell_attached_data.n_attached_deserialize;
+
+ // important: only remove data if we are not in the deserialization
+ // process. There, each SolutionTransfer registers and unpacks before
+ // the next one does this, so n_attached_data_sets is only 1 here. This
+ // would destroy the saved data before the second SolutionTransfer can
+ // get it. This created a bug that is documented in
+ // tests/mpi/p4est_save_03 with more than one SolutionTransfer.
+ if (cell_attached_data.n_attached_data_sets == 0 &&
+ cell_attached_data.n_attached_deserialize == 0)
+ {
+ // everybody got their data, time for cleanup!
+ cell_attached_data.pack_callbacks_fixed.clear();
+ cell_attached_data.pack_callbacks_variable.clear();
+ }
+
+ if (this->cell_attached_data.n_attached_data_sets == 0 &&
+ this->cell_attached_data.n_attached_deserialize == 0)
+ {
+ // reset all cell_status entries after coarsening/refinement
+ for (auto &cell_rel : local_cell_relations)
+ std::get<0>(cell_rel) =
+ parallel::DistributedTriangulationBase<dim, spacedim>::CELL_PERSIST;
+ }
+ }
+
+
+
+ /* ------------------ class DataTransfer<dim,spacedim> ----------------- */
+
+
+ template <int dim, int spacedim>
+ DistributedTriangulationBase<dim, spacedim>::DataTransfer::DataTransfer(
+ const MPI_Comm &mpi_communicator)
+ : variable_size_data_stored(false)
+ , mpi_communicator(mpi_communicator)
+ {}
+
+
+
+ template <int dim, int spacedim>
+ void
+ DistributedTriangulationBase<dim, spacedim>::DataTransfer::pack_data(
+ const std::vector<cell_relation_t> &cell_relations,
+ const std::vector<typename CellAttachedData::pack_callback_t>
+ &pack_callbacks_fixed,
+ const std::vector<typename CellAttachedData::pack_callback_t>
+ &pack_callbacks_variable)
+ {
+ Assert(src_data_fixed.size() == 0,
+ ExcMessage("Previously packed data has not been released yet!"));
+ Assert(src_sizes_variable.size() == 0, ExcInternalError());
+
+ const unsigned int n_callbacks_fixed = pack_callbacks_fixed.size();
+ const unsigned int n_callbacks_variable = pack_callbacks_variable.size();
+
+ // Store information that we packed variable size data in
+ // a member variable for later.
+ variable_size_data_stored = (n_callbacks_variable > 0);
+
+ // If variable transfer is scheduled, we will store the data size that
+ // each variable size callback function writes in this auxiliary
+ // container. The information will be stored by each cell in this vector
+ // temporarily.
+ std::vector<unsigned int> cell_sizes_variable_cumulative(
+ n_callbacks_variable);
+
+ // Prepare the buffer structure, in which each callback function will
+ // store its data for each active cell.
+ // The outmost shell in this container construct corresponds to the
+ // data packed per cell. The next layer resembles the data that
+ // each callback function packs on the corresponding cell. These
+ // buffers are chains of chars stored in an std::vector<char>.
+ // A visualisation of the data structure:
+ /* clang-format off */
+ // | cell_1 | | cell_2 | ...
+ // || callback_1 || callback_2 |...| || callback_1 || callback_2 |...| ...
+ // |||char|char|...|||char|char|...|...| |||char|char|...|||char|char|...|...| ...
+ /* clang-format on */
+ std::vector<std::vector<std::vector<char>>> packed_fixed_size_data(
+ cell_relations.size());
+ std::vector<std::vector<std::vector<char>>> packed_variable_size_data(
+ variable_size_data_stored ? cell_relations.size() : 0);
+
+ //
+ // --------- Pack data for fixed and variable size transfer ---------
+ //
+ // Iterate over all cells, call all callback functions on each cell,
+ // and store their data in the corresponding buffer scope.
+ {
+ auto cell_rel_it = cell_relations.cbegin();
+ auto data_cell_fixed_it = packed_fixed_size_data.begin();
+ auto data_cell_variable_it = packed_variable_size_data.begin();
+ for (; cell_rel_it != cell_relations.cend(); ++cell_rel_it)
+ {
+ const auto &cell_status = std::get<0>(*cell_rel_it);
+ const auto &dealii_cell = std::get<1>(*cell_rel_it);
+
+ // Assertions about the tree structure.
+ switch (cell_status)
+ {
+ case parallel::DistributedTriangulationBase<dim, spacedim>::
+ CELL_PERSIST:
+ case parallel::DistributedTriangulationBase<dim, spacedim>::
+ CELL_REFINE:
+ // double check the condition that we will only ever attach
+ // data to active cells when we get here
+ Assert(dealii_cell->is_active(), ExcInternalError());
+ break;
+
+ case parallel::DistributedTriangulationBase<dim, spacedim>::
+ CELL_COARSEN:
+ // double check the condition that we will only ever attach
+ // data to cells with children when we get here. however, we
+ // can only tolerate one level of coarsening at a time, so
+ // check that the children are all active
+ Assert(dealii_cell->is_active() == false, ExcInternalError());
+ for (unsigned int c = 0;
+ c < GeometryInfo<dim>::max_children_per_cell;
+ ++c)
+ Assert(dealii_cell->child(c)->is_active(),
+ ExcInternalError());
+ break;
+
+ case parallel::DistributedTriangulationBase<dim, spacedim>::
+ CELL_INVALID:
+ // do nothing on invalid cells
+ break;
+
+ default:
+ Assert(false, ExcInternalError());
+ break;
+ }
+
+ // Reserve memory corresponding to the number of callback
+ // functions that will be called.
+ // If variable size transfer is scheduled, we need to leave
+ // room for an array that holds information about how many
+ // bytes each of the variable size callback functions will
+ // write.
+ // On cells flagged with CELL_INVALID, only its CellStatus
+ // will be stored.
+ const unsigned int n_fixed_size_data_sets_on_cell =
+ 1 +
+ ((cell_status ==
+ parallel::DistributedTriangulationBase<dim,
+ spacedim>::CELL_INVALID) ?
+ 0 :
+ ((variable_size_data_stored ? 1 : 0) + n_callbacks_fixed));
+ data_cell_fixed_it->resize(n_fixed_size_data_sets_on_cell);
+
+ // We continue with packing all data on this specific cell.
+ auto data_fixed_it = data_cell_fixed_it->begin();
+
+ // First, we pack the CellStatus information.
+ // to get consistent data sizes on each cell for the fixed size
+ // transfer, we won't allow compression
+ *data_fixed_it =
+ Utilities::pack(cell_status, /*allow_compression=*/false);
+ ++data_fixed_it;
+
+ // Proceed with all registered callback functions.
+ // Skip cells with the CELL_INVALID flag.
+ if (cell_status !=
+ parallel::DistributedTriangulationBase<dim,
+ spacedim>::CELL_INVALID)
+ {
+ // Pack fixed size data.
+ for (auto callback_it = pack_callbacks_fixed.cbegin();
+ callback_it != pack_callbacks_fixed.cend();
+ ++callback_it, ++data_fixed_it)
+ {
+ *data_fixed_it = (*callback_it)(dealii_cell, cell_status);
+ }
+
+ // Pack variable size data.
+ // If we store variable size data, we need to transfer
+ // the sizes of each corresponding callback function
+ // via fixed size transfer as well.
+ if (variable_size_data_stored)
+ {
+ const unsigned int n_variable_size_data_sets_on_cell =
+ ((cell_status ==
+ parallel::DistributedTriangulationBase<dim, spacedim>::
+ CELL_INVALID) ?
+ 0 :
+ n_callbacks_variable);
+ data_cell_variable_it->resize(
+ n_variable_size_data_sets_on_cell);
+
+ auto callback_it = pack_callbacks_variable.cbegin();
+ auto data_variable_it = data_cell_variable_it->begin();
+ auto sizes_variable_it =
+ cell_sizes_variable_cumulative.begin();
+ for (; callback_it != pack_callbacks_variable.cend();
+ ++callback_it, ++data_variable_it, ++sizes_variable_it)
+ {
+ *data_variable_it =
+ (*callback_it)(dealii_cell, cell_status);
+
+ // Store data sizes for each callback function first.
+ // Make it cumulative below.
+ *sizes_variable_it = data_variable_it->size();
+ }
+
+ // Turn size vector into its cumulative representation.
+ std::partial_sum(cell_sizes_variable_cumulative.begin(),
+ cell_sizes_variable_cumulative.end(),
+ cell_sizes_variable_cumulative.begin());
+
+ // Serialize cumulative variable size vector value-by-value.
+ // This way we can circumvent the overhead of storing the
+ // container object as a whole, since we know its size by
+ // the number of registered callback functions.
+ data_fixed_it->resize(n_callbacks_variable *
+ sizeof(unsigned int));
+ for (unsigned int i = 0; i < n_callbacks_variable; ++i)
+ std::memcpy(&(data_fixed_it->at(i * sizeof(unsigned int))),
+ &(cell_sizes_variable_cumulative.at(i)),
+ sizeof(unsigned int));
+
+ ++data_fixed_it;
+ }
+
+ // Double check that we packed everything we wanted
+ // in the fixed size buffers.
+ Assert(data_fixed_it == data_cell_fixed_it->end(),
+ ExcInternalError());
+ }
+
+ ++data_cell_fixed_it;
+
+ // Increment the variable size data iterator
+ // only if we actually pack this kind of data
+ // to avoid getting out of bounds.
+ if (variable_size_data_stored)
+ ++data_cell_variable_it;
+ } // loop over cell_relations
+ }
+
+ //
+ // ----------- Gather data sizes for fixed size transfer ------------
+ //
+ // Generate a vector which stores the sizes of each callback function,
+ // including the packed CellStatus transfer.
+ // Find the very first cell that we wrote to with all callback
+ // functions (i.e. a cell that was not flagged with CELL_INVALID)
+ // and store the sizes of each buffer.
+ //
+ // To deal with the case that at least one of the processors does not own
+ // any cell at all, we will exchange the information about the data sizes
+ // among them later. The code in between is still well-defined, since the
+ // following loops will be skipped.
+ std::vector<unsigned int> local_sizes_fixed(
+ 1 + n_callbacks_fixed + (variable_size_data_stored ? 1 : 0));
+ for (const auto &data_cell : packed_fixed_size_data)
+ {
+ if (data_cell.size() == local_sizes_fixed.size())
+ {
+ auto sizes_fixed_it = local_sizes_fixed.begin();
+ auto data_fixed_it = data_cell.cbegin();
+ for (; data_fixed_it != data_cell.cend();
+ ++data_fixed_it, ++sizes_fixed_it)
+ {
+ *sizes_fixed_it = data_fixed_it->size();
+ }
+
+ break;
+ }
+ }
+
+ // Check if all cells have valid sizes.
+ for (auto data_cell_fixed_it = packed_fixed_size_data.cbegin();
+ data_cell_fixed_it != packed_fixed_size_data.cend();
+ ++data_cell_fixed_it)
+ {
+ Assert((data_cell_fixed_it->size() == 1) ||
+ (data_cell_fixed_it->size() == local_sizes_fixed.size()),
+ ExcInternalError());
+ }
+
+ // Share information about the packed data sizes
+ // of all callback functions across all processors, in case one
+ // of them does not own any cells at all.
+ std::vector<unsigned int> global_sizes_fixed(local_sizes_fixed.size());
+ Utilities::MPI::max(local_sizes_fixed,
+ this->mpi_communicator,
+ global_sizes_fixed);
+
+ // Construct cumulative sizes, since this is the only information
+ // we need from now on.
+ sizes_fixed_cumulative.resize(global_sizes_fixed.size());
+ std::partial_sum(global_sizes_fixed.begin(),
+ global_sizes_fixed.end(),
+ sizes_fixed_cumulative.begin());
+
+ //
+ // ---------- Gather data sizes for variable size transfer ----------
+ //
+ if (variable_size_data_stored)
+ {
+ src_sizes_variable.reserve(packed_variable_size_data.size());
+ for (const auto &data_cell : packed_variable_size_data)
+ {
+ int variable_data_size_on_cell = 0;
+
+ for (const auto &data : data_cell)
+ variable_data_size_on_cell += data.size();
+
+ src_sizes_variable.push_back(variable_data_size_on_cell);
+ }
+ }
+
+ //
+ // ------------------------ Build buffers ---------------------------
+ //
+ const unsigned int expected_size_fixed =
+ cell_relations.size() * sizes_fixed_cumulative.back();
+ const unsigned int expected_size_variable =
+ std::accumulate(src_sizes_variable.begin(),
+ src_sizes_variable.end(),
+ std::vector<int>::size_type(0));
+
+ // Move every piece of packed fixed size data into the consecutive buffer.
+ src_data_fixed.reserve(expected_size_fixed);
+ for (const auto &data_cell_fixed : packed_fixed_size_data)
+ {
+ // Move every fraction of packed data into the buffer
+ // reserved for this particular cell.
+ for (const auto &data_fixed : data_cell_fixed)
+ std::move(data_fixed.begin(),
+ data_fixed.end(),
+ std::back_inserter(src_data_fixed));
+
+ // If we only packed the CellStatus information
+ // (i.e. encountered a cell flagged CELL_INVALID),
+ // fill the remaining space with invalid entries.
+ // We can skip this if there is nothing else to pack.
+ if ((data_cell_fixed.size() == 1) &&
+ (sizes_fixed_cumulative.size() > 1))
+ {
+ const std::size_t bytes_skipped =
+ sizes_fixed_cumulative.back() - sizes_fixed_cumulative.front();
+
+ src_data_fixed.insert(src_data_fixed.end(),
+ bytes_skipped,
+ static_cast<char>(-1)); // invalid_char
+ }
+ }
+
+ // Move every piece of packed variable size data into the consecutive
+ // buffer.
+ if (variable_size_data_stored)
+ {
+ src_data_variable.reserve(expected_size_variable);
+ for (const auto &data_cell : packed_variable_size_data)
+ {
+ // Move every fraction of packed data into the buffer
+ // reserved for this particular cell.
+ for (const auto &data : data_cell)
+ std::move(data.begin(),
+ data.end(),
+ std::back_inserter(src_data_variable));
+ }
+ }
+
+ // Double check that we packed everything correctly.
+ Assert(src_data_fixed.size() == expected_size_fixed, ExcInternalError());
+ Assert(src_data_variable.size() == expected_size_variable,
+ ExcInternalError());
+ }
+
+
+
+ template <int dim, int spacedim>
+ void
+ DistributedTriangulationBase<dim, spacedim>::DataTransfer::unpack_cell_status(
+ std::vector<cell_relation_t> &cell_relations) const
+ {
+ Assert(sizes_fixed_cumulative.size() > 0,
+ ExcMessage("No data has been packed!"));
+ if (cell_relations.size() > 0)
+ {
+ Assert(dest_data_fixed.size() > 0,
+ ExcMessage("No data has been received!"));
+ }
+
+ // Size of CellStatus object that will be unpacked on each cell.
+ const unsigned int size = sizes_fixed_cumulative.front();
+
+ // Iterate over all cells and overwrite the CellStatus
+ // information from the transferred data.
+ // Proceed buffer iterator position to next cell after
+ // each iteration.
+ auto cell_rel_it = cell_relations.begin();
+ auto dest_fixed_it = dest_data_fixed.cbegin();
+ for (; cell_rel_it != cell_relations.end();
+ ++cell_rel_it, dest_fixed_it += sizes_fixed_cumulative.back())
+ {
+ std::get<0>(*cell_rel_it) = // cell_status
+ Utilities::unpack<typename parallel::DistributedTriangulationBase<
+ dim,
+ spacedim>::CellStatus>(dest_fixed_it,
+ dest_fixed_it + size,
+ /*allow_compression=*/false);
+ }
+ }
+
+
+
+ template <int dim, int spacedim>
+ void
+ DistributedTriangulationBase<dim, spacedim>::DataTransfer::unpack_data(
+ const std::vector<cell_relation_t> &cell_relations,
+ const unsigned int handle,
+ const std::function<
+ void(const typename dealii::Triangulation<dim, spacedim>::cell_iterator &,
+ const typename dealii::Triangulation<dim, spacedim>::CellStatus &,
+ const boost::iterator_range<std::vector<char>::const_iterator> &)>
+ &unpack_callback) const
+ {
+ // We decode the handle returned by register_data_attach() back into
+ // a format we can use. All even handles belong to those callback
+ // functions which write/read variable size data, all odd handles interact
+ // with fixed size buffers.
+ const bool callback_variable_transfer = (handle % 2 == 0);
+ const unsigned int callback_index = handle / 2;
+
+ // Cells will always receive fixed size data (i.e., CellStatus
+ // information), but not necessarily variable size data (e.g., with a
+ // ParticleHandler a cell might not contain any particle at all).
+ // Thus it is sufficient to check if fixed size data has been received.
+ Assert(sizes_fixed_cumulative.size() > 0,
+ ExcMessage("No data has been packed!"));
+ if (cell_relations.size() > 0)
+ {
+ Assert(dest_data_fixed.size() > 0,
+ ExcMessage("No data has been received!"));
+ }
+
+ std::vector<char>::const_iterator dest_data_it;
+ std::vector<char>::const_iterator dest_sizes_cell_it;
+
+ // Depending on whether our callback function unpacks fixed or
+ // variable size data, we have to pursue different approaches
+ // to localize the correct fraction of the buffer from which
+ // we are allowed to read.
+ unsigned int offset = numbers::invalid_unsigned_int;
+ unsigned int size = numbers::invalid_unsigned_int;
+ unsigned int data_increment = numbers::invalid_unsigned_int;
+
+ if (callback_variable_transfer)
+ {
+ // For the variable size data, we need to extract the
+ // data size from the fixed size buffer on each cell.
+ //
+ // We packed this information last, so the last packed
+ // object in the fixed size buffer corresponds to the
+ // variable data sizes.
+ //
+ // The last entry of sizes_fixed_cumulative corresponds
+ // to the size of all fixed size data packed on the cell.
+ // To get the offset for the last packed object, we need
+ // to get the next-to-last entry.
+ const unsigned int offset_variable_data_sizes =
+ sizes_fixed_cumulative[sizes_fixed_cumulative.size() - 2];
+
+ // This iterator points to the data size that the
+ // callback_function packed for each specific cell.
+ // Adjust buffer iterator to the offset of the callback
+ // function so that we only have to advance its position
+ // to the next cell after each iteration.
+ dest_sizes_cell_it = dest_data_fixed.cbegin() +
+ offset_variable_data_sizes +
+ callback_index * sizeof(unsigned int);
+
+ // Let the data iterator point to the correct buffer.
+ dest_data_it = dest_data_variable.cbegin();
+ }
+ else
+ {
+ // For the fixed size data, we can get the information about
+ // the buffer location on each cell directly from the
+ // sizes_fixed_cumulative vector.
+ offset = sizes_fixed_cumulative[callback_index];
+ size = sizes_fixed_cumulative[callback_index + 1] - offset;
+ data_increment = sizes_fixed_cumulative.back();
+
+ // Let the data iterator point to the correct buffer.
+ // Adjust buffer iterator to the offset of the callback
+ // function so that we only have to advance its position
+ // to the next cell after each iteration.
+ dest_data_it = dest_data_fixed.cbegin() + offset;
+ }
+
+ // Iterate over all cells and unpack the transferred data.
+ auto cell_rel_it = cell_relations.begin();
+ auto dest_sizes_it = dest_sizes_variable.cbegin();
+ for (; cell_rel_it != cell_relations.end(); ++cell_rel_it)
+ {
+ const auto &cell_status = std::get<0>(*cell_rel_it);
+ const auto &dealii_cell = std::get<1>(*cell_rel_it);
+
+ if (callback_variable_transfer)
+ {
+ // Update the increment according to the whole data size
+ // of the current cell.
+ data_increment = *dest_sizes_it;
+
+ if (cell_status !=
+ parallel::DistributedTriangulationBase<dim,
+ spacedim>::CELL_INVALID)
+ {
+ // Extract the corresponding values for offset and size from
+ // the cumulative sizes array stored in the fixed size buffer.
+ if (callback_index == 0)
+ offset = 0;
+ else
+ std::memcpy(&offset,
+ &(*(dest_sizes_cell_it - sizeof(unsigned int))),
+ sizeof(unsigned int));
+
+ std::memcpy(&size,
+ &(*dest_sizes_cell_it),
+ sizeof(unsigned int));
+
+ size -= offset;
+
+ // Move the data iterator to the corresponding position
+ // of the callback function and adjust the increment
+ // accordingly.
+ dest_data_it += offset;
+ data_increment -= offset;
+ }
+
+ // Advance data size iterators to the next cell.
+ dest_sizes_cell_it += sizes_fixed_cumulative.back();
+ ++dest_sizes_it;
+ }
+
+ switch (cell_status)
+ {
+ case parallel::DistributedTriangulationBase<dim,
+ spacedim>::CELL_PERSIST:
+ case parallel::DistributedTriangulationBase<dim,
+ spacedim>::CELL_COARSEN:
+ unpack_callback(dealii_cell,
+ cell_status,
+ boost::make_iterator_range(dest_data_it,
+ dest_data_it + size));
+ break;
+
+ case parallel::DistributedTriangulationBase<dim,
+ spacedim>::CELL_REFINE:
+ unpack_callback(dealii_cell->parent(),
+ cell_status,
+ boost::make_iterator_range(dest_data_it,
+ dest_data_it + size));
+ break;
+
+ case parallel::DistributedTriangulationBase<dim,
+ spacedim>::CELL_INVALID:
+ // Skip this cell.
+ break;
+
+ default:
+ Assert(false, ExcInternalError());
+ break;
+ }
+
+ dest_data_it += data_increment;
+ }
+ }
+
+
+
+ template <int dim, int spacedim>
+ void
+ DistributedTriangulationBase<dim, spacedim>::DataTransfer::save(
+ const unsigned int global_first_cell,
+ const unsigned int global_num_cells,
+ const std::string &filename) const
+ {
+#ifdef DEAL_II_WITH_MPI
+ // Large fractions of this function have been copied from
+ // DataOutInterface::write_vtu_in_parallel.
+ // TODO: Write general MPIIO interface.
+
+ Assert(sizes_fixed_cumulative.size() > 0,
+ ExcMessage("No data has been packed!"));
+
+ const int myrank = Utilities::MPI::this_mpi_process(mpi_communicator);
+
+ //
+ // ---------- Fixed size data ----------
+ //
+ {
+ const std::string fname_fixed = std::string(filename) + "_fixed.data";
+
+ MPI_Info info;
+ int ierr = MPI_Info_create(&info);
+ AssertThrowMPI(ierr);
+
+ MPI_File fh;
+ ierr = MPI_File_open(mpi_communicator,
+ DEAL_II_MPI_CONST_CAST(fname_fixed.c_str()),
+ MPI_MODE_CREATE | MPI_MODE_WRONLY,
+ info,
+ &fh);
+ AssertThrowMPI(ierr);
+
+ ierr = MPI_File_set_size(fh, 0); // delete the file contents
+ AssertThrowMPI(ierr);
+ // this barrier is necessary, because otherwise others might already
+ // write while one core is still setting the size to zero.
+ ierr = MPI_Barrier(mpi_communicator);
+ AssertThrowMPI(ierr);
+ ierr = MPI_Info_free(&info);
+ AssertThrowMPI(ierr);
+ // ------------------
+
+ // Write cumulative sizes to file.
+ // Since each processor owns the same information about the data sizes,
+ // it is sufficient to let only the first processor perform this task.
+ if (myrank == 0)
+ {
+ const unsigned int *data = sizes_fixed_cumulative.data();
+
+ ierr = MPI_File_write_at(fh,
+ 0,
+ DEAL_II_MPI_CONST_CAST(data),
+ sizes_fixed_cumulative.size(),
+ MPI_UNSIGNED,
+ MPI_STATUS_IGNORE);
+ AssertThrowMPI(ierr);
+ }
+
+ // Write packed data to file simultaneously.
+ const unsigned int offset_fixed =
+ sizes_fixed_cumulative.size() * sizeof(unsigned int);
+
+ const char *data = src_data_fixed.data();
+
+ ierr = MPI_File_write_at(
+ fh,
+ offset_fixed +
+ global_first_cell *
+ sizes_fixed_cumulative.back(), // global position in file
+ DEAL_II_MPI_CONST_CAST(data),
+ src_data_fixed.size(), // local buffer
+ MPI_CHAR,
+ MPI_STATUS_IGNORE);
+ AssertThrowMPI(ierr);
+
+ ierr = MPI_File_close(&fh);
+ AssertThrowMPI(ierr);
+ }
+
+ //
+ // ---------- Variable size data ----------
+ //
+ if (variable_size_data_stored)
+ {
+ const std::string fname_variable =
+ std::string(filename) + "_variable.data";
+
+ MPI_Info info;
+ int ierr = MPI_Info_create(&info);
+ AssertThrowMPI(ierr);
+
+ MPI_File fh;
+ ierr = MPI_File_open(mpi_communicator,
+ DEAL_II_MPI_CONST_CAST(fname_variable.c_str()),
+ MPI_MODE_CREATE | MPI_MODE_WRONLY,
+ info,
+ &fh);
+ AssertThrowMPI(ierr);
+
+ ierr = MPI_File_set_size(fh, 0); // delete the file contents
+ AssertThrowMPI(ierr);
+ // this barrier is necessary, because otherwise others might already
+ // write while one core is still setting the size to zero.
+ ierr = MPI_Barrier(mpi_communicator);
+ AssertThrowMPI(ierr);
+ ierr = MPI_Info_free(&info);
+ AssertThrowMPI(ierr);
+
+ // Write sizes of each cell into file simultaneously.
+ {
+ const int *data = src_sizes_variable.data();
+ ierr =
+ MPI_File_write_at(fh,
+ global_first_cell *
+ sizeof(unsigned int), // global position in file
+ DEAL_II_MPI_CONST_CAST(data),
+ src_sizes_variable.size(), // local buffer
+ MPI_INT,
+ MPI_STATUS_IGNORE);
+ AssertThrowMPI(ierr);
+ }
+
+
+ const unsigned int offset_variable =
+ global_num_cells * sizeof(unsigned int);
+
+ // Gather size of data in bytes we want to store from this processor.
+ const unsigned int size_on_proc = src_data_variable.size();
+
+ // Compute prefix sum
+ unsigned int prefix_sum = 0;
+ ierr = MPI_Exscan(DEAL_II_MPI_CONST_CAST(&size_on_proc),
+ &prefix_sum,
+ 1,
+ MPI_UNSIGNED,
+ MPI_SUM,
+ mpi_communicator);
+ AssertThrowMPI(ierr);
+
+ const char *data = src_data_variable.data();
+
+ // Write data consecutively into file.
+ ierr = MPI_File_write_at(fh,
+ offset_variable +
+ prefix_sum, // global position in file
+ DEAL_II_MPI_CONST_CAST(data),
+ src_data_variable.size(), // local buffer
+ MPI_CHAR,
+ MPI_STATUS_IGNORE);
+ AssertThrowMPI(ierr);
+
+ ierr = MPI_File_close(&fh);
+ AssertThrowMPI(ierr);
+ }
+#else
+ (void)global_first_cell;
+ (void)global_num_cells;
+ (void)filename;
+
+ AssertThrow(false, ExcNeedsMPI());
+#endif
+ }
+
+
+
+ template <int dim, int spacedim>
+ void
+ DistributedTriangulationBase<dim, spacedim>::DataTransfer::load(
+ const unsigned int global_first_cell,
+ const unsigned int global_num_cells,
+ const unsigned int local_num_cells,
+ const std::string &filename,
+ const unsigned int n_attached_deserialize_fixed,
+ const unsigned int n_attached_deserialize_variable)
+ {
+#ifdef DEAL_II_WITH_MPI
+ // Large fractions of this function have been copied from
+ // DataOutInterface::write_vtu_in_parallel.
+ // TODO: Write general MPIIO interface.
+
+ Assert(dest_data_fixed.size() == 0,
+ ExcMessage("Previously loaded data has not been released yet!"));
+
+ variable_size_data_stored = (n_attached_deserialize_variable > 0);
+
+ //
+ // ---------- Fixed size data ----------
+ //
+ {
+ const std::string fname_fixed = std::string(filename) + "_fixed.data";
+
+ MPI_Info info;
+ int ierr = MPI_Info_create(&info);
+ AssertThrowMPI(ierr);
+
+ MPI_File fh;
+ ierr = MPI_File_open(mpi_communicator,
+ DEAL_II_MPI_CONST_CAST(fname_fixed.c_str()),
+ MPI_MODE_RDONLY,
+ info,
+ &fh);
+ AssertThrowMPI(ierr);
+
+ ierr = MPI_Info_free(&info);
+ AssertThrowMPI(ierr);
+
+ // Read cumulative sizes from file.
+ // Since all processors need the same information about the data sizes,
+ // let each of them retrieve it by reading from the same location in
+ // the file.
+ sizes_fixed_cumulative.resize(1 + n_attached_deserialize_fixed +
+ (variable_size_data_stored ? 1 : 0));
+ ierr = MPI_File_read_at(fh,
+ 0,
+ sizes_fixed_cumulative.data(),
+ sizes_fixed_cumulative.size(),
+ MPI_UNSIGNED,
+ MPI_STATUS_IGNORE);
+ AssertThrowMPI(ierr);
+
+ // Allocate sufficient memory.
+ dest_data_fixed.resize(local_num_cells * sizes_fixed_cumulative.back());
+
+ // Read packed data from file simultaneously.
+ const unsigned int offset =
+ sizes_fixed_cumulative.size() * sizeof(unsigned int);
+
+ ierr = MPI_File_read_at(
+ fh,
+ offset + global_first_cell *
+ sizes_fixed_cumulative.back(), // global position in file
+ dest_data_fixed.data(),
+ dest_data_fixed.size(), // local buffer
+ MPI_CHAR,
+ MPI_STATUS_IGNORE);
+ AssertThrowMPI(ierr);
+
+ ierr = MPI_File_close(&fh);
+ AssertThrowMPI(ierr);
+ }
+
+ //
+ // ---------- Variable size data ----------
+ //
+ if (variable_size_data_stored)
+ {
+ const std::string fname_variable =
+ std::string(filename) + "_variable.data";
+
+ MPI_Info info;
+ int ierr = MPI_Info_create(&info);
+ AssertThrowMPI(ierr);
+
+ MPI_File fh;
+ ierr = MPI_File_open(mpi_communicator,
+ DEAL_II_MPI_CONST_CAST(fname_variable.c_str()),
+ MPI_MODE_RDONLY,
+ info,
+ &fh);
+ AssertThrowMPI(ierr);
+
+ ierr = MPI_Info_free(&info);
+ AssertThrowMPI(ierr);
+
+ // Read sizes of all locally owned cells.
+ dest_sizes_variable.resize(local_num_cells);
+ ierr = MPI_File_read_at(fh,
+ global_first_cell * sizeof(unsigned int),
+ dest_sizes_variable.data(),
+ dest_sizes_variable.size(),
+ MPI_INT,
+ MPI_STATUS_IGNORE);
+ AssertThrowMPI(ierr);
+
+ const unsigned int offset = global_num_cells * sizeof(unsigned int);
+
+ const unsigned int size_on_proc =
+ std::accumulate(dest_sizes_variable.begin(),
+ dest_sizes_variable.end(),
+ 0);
+
+ // share information among all processors by prefix sum
+ unsigned int prefix_sum = 0;
+ ierr = MPI_Exscan(DEAL_II_MPI_CONST_CAST(&size_on_proc),
+ &prefix_sum,
+ 1,
+ MPI_UNSIGNED,
+ MPI_SUM,
+ mpi_communicator);
+ AssertThrowMPI(ierr);
+
+ dest_data_variable.resize(size_on_proc);
+ ierr = MPI_File_read_at(fh,
+ offset + prefix_sum,
+ dest_data_variable.data(),
+ dest_data_variable.size(),
+ MPI_CHAR,
+ MPI_STATUS_IGNORE);
+ AssertThrowMPI(ierr);
+
+ ierr = MPI_File_close(&fh);
+ AssertThrowMPI(ierr);
+ }
+#else
+ (void)global_first_cell;
+ (void)global_num_cells;
+ (void)local_num_cells;
+ (void)filename;
+ (void)n_attached_deserialize_fixed;
+ (void)n_attached_deserialize_variable;
+
+ AssertThrow(false, ExcNeedsMPI());
+#endif
+ }
+
+
+
+ template <int dim, int spacedim>
+ void
+ DistributedTriangulationBase<dim, spacedim>::DataTransfer::clear()
+ {
+ variable_size_data_stored = false;
+
+ // free information about data sizes
+ sizes_fixed_cumulative.clear();
+ sizes_fixed_cumulative.shrink_to_fit();
+
+ // free fixed size transfer data
+ src_data_fixed.clear();
+ src_data_fixed.shrink_to_fit();
+
+ dest_data_fixed.clear();
+ dest_data_fixed.shrink_to_fit();
+
+ // free variable size transfer data
+ src_sizes_variable.clear();
+ src_sizes_variable.shrink_to_fit();
+
+ src_data_variable.clear();
+ src_data_variable.shrink_to_fit();
+
+ dest_sizes_variable.clear();
+ dest_sizes_variable.shrink_to_fit();
+
+ dest_data_variable.clear();
+ dest_data_variable.shrink_to_fit();
+ }
+
} // end namespace parallel
projects / dealii.git / commitdiff