--- /dev/null
+New: The method parallel::distributed::Triagnulation::communicate_locally_moved_vertices()
+has been refactored and moved to parallel::DistributedTriangulationBase so that it can now also be
+used for parallel::fullydistributed::Triangulation.
+<br>
+(Daniel Arndt, Ivan Fumagalli, Peter Munch, 2021/01/15)
+
+/**
+ * A class that always returns a given value.
+ * This is needed as a workaround for lambdas used as default parameters
+ * some compilers struggle to deal with.
+ */
+template <typename ArgType, typename ValueType>
+struct always_return
+{
+ ValueType value;
+ ValueType
+ operator()(const ArgType &)
+ {
+ return value;
+ }
+};
+
+
+
/**
* A class to perform comparisons of arbitrary pointers for equality. In some
* circumstances, one would like to make sure that two arguments to a function
void
repartition();
- /**
- * When vertices have been moved locally, for example using code like
- * @code
- * cell->vertex(0) = new_location;
- * @endcode
- * then this function can be used to update the location of vertices
- * between MPI processes.
- *
- * All the vertices that have been moved and might be in the ghost layer
- * of a process have to be reported in the @p vertex_locally_moved
- * argument. This ensures that that part of the information that has to
- * be send between processes is actually sent. Additionally, it is quite
- * important that vertices on the boundary between processes are
- * reported on exactly one process (e.g. the one with the highest id).
- * Otherwise we could expect undesirable results if multiple processes
- * move a vertex differently. A typical strategy is to let processor $i$
- * move those vertices that are adjacent to cells whose owners include
- * processor $i$ but no other processor $j$ with $j<i$; in other words,
- * for vertices at the boundary of a subdomain, the processor with the
- * lowest subdomain id "owns" a vertex.
- *
- * @note It only makes sense to move vertices that are either located on
- * locally owned cells or on cells in the ghost layer. This is because
- * you can be sure that these vertices indeed exist on the finest mesh
- * aggregated over all processors, whereas vertices on artificial cells
- * but not at least in the ghost layer may or may not exist on the
- * globally finest mesh. Consequently, the @p vertex_locally_moved
- * argument may not contain vertices that aren't at least on ghost
- * cells.
- *
- * @note This function moves vertices in such a way that on every
- * processor, the vertices of every locally owned and ghost cell is
- * consistent with the corresponding location of these cells on other
- * processors. On the other hand, the locations of artificial cells will
- * in general be wrong since artificial cells may or may not exist on
- * other processors and consequently it is not possible to determine
- * their location in any way. This is not usually a problem since one
- * never does anything on artificial cells. However, it may lead to
- * problems if the mesh with moved vertices is refined in a later step.
- * If that's what you want to do, the right way to do it is to save the
- * offset applied to every vertex, call this function, and before
- * refining or coarsening the mesh apply the opposite offset and call
- * this function again.
- *
- * @param vertex_locally_moved A bitmap indicating which vertices have
- * been moved. The size of this array must be equal to
- * Triangulation::n_vertices() and must be a subset of those vertices
- * flagged by GridTools::get_locally_owned_vertices().
- *
- * @see This function is used, for example, in
- * GridTools::distort_random().
- */
- void
- communicate_locally_moved_vertices(
- const std::vector<bool> &vertex_locally_moved);
-
/**
* Return true if the triangulation has hanging nodes.
const std::vector<types::global_dof_index> &
get_p4est_tree_to_coarse_cell_permutation() const;
- /**
- * When vertices have been moved locally, for example using code like
- * @code
- * cell->vertex(0) = new_location;
- * @endcode
- * then this function can be used to update the location of vertices
- * between MPI processes.
- *
- * All the vertices that have been moved and might be in the ghost layer
- * of a process have to be reported in the @p vertex_locally_moved
- * argument. This ensures that that part of the information that has to
- * be send between processes is actually sent. Additionally, it is quite
- * important that vertices on the boundary between processes are
- * reported on exactly one process (e.g. the one with the highest id).
- * Otherwise we could expect undesirable results if multiple processes
- * move a vertex differently. A typical strategy is to let processor $i$
- * move those vertices that are adjacent to cells whose owners include
- * processor $i$ but no other processor $j$ with $j<i$; in other words,
- * for vertices at the boundary of a subdomain, the processor with the
- * lowest subdomain id "owns" a vertex.
- *
- * @note It only makes sense to move vertices that are either located on
- * locally owned cells or on cells in the ghost layer. This is because
- * you can be sure that these vertices indeed exist on the finest mesh
- * aggregated over all processors, whereas vertices on artificial cells
- * but not at least in the ghost layer may or may not exist on the
- * globally finest mesh. Consequently, the @p vertex_locally_moved
- * argument may not contain vertices that aren't at least on ghost
- * cells.
- *
- * @see This function is used, for example, in
- * GridTools::distort_random().
- */
- void
- communicate_locally_moved_vertices(
- const std::vector<bool> &vertex_locally_moved);
-
/**
* This function is not implemented, but needs to be present for the
* compiler.
virtual std::vector<types::manifold_id>
get_manifold_ids() const override;
+ /**
+ * When vertices have been moved locally, for example using code like
+ * @code
+ * cell->vertex(0) = new_location;
+ * @endcode
+ * then this function can be used to update the location of vertices
+ * between MPI processes.
+ *
+ * All the vertices that have been moved and might be in the ghost layer
+ * of a process have to be reported in the @p vertex_locally_moved
+ * argument. This ensures that that part of the information that has to
+ * be send between processes is actually sent. Additionally, it is quite
+ * important that vertices on the boundary between processes are
+ * reported on exactly one process (e.g. the one with the highest id).
+ * Otherwise we could expect undesirable results if multiple processes
+ * move a vertex differently. A typical strategy is to let processor $i$
+ * move those vertices that are adjacent to cells whose owners include
+ * processor $i$ but no other processor $j$ with $j<i$; in other words,
+ * for vertices at the boundary of a subdomain, the processor with the
+ * lowest subdomain id "owns" a vertex.
+ *
+ * @note It only makes sense to move vertices that are either located on
+ * locally owned cells or on cells in the ghost layer. This is because
+ * you can be sure that these vertices indeed exist on the finest mesh
+ * aggregated over all processors, whereas vertices on artificial cells
+ * but not at least in the ghost layer may or may not exist on the
+ * globally finest mesh. Consequently, the @p vertex_locally_moved
+ * argument may not contain vertices that aren't at least on ghost
+ * cells.
+ *
+ * @note This function moves vertices in such a way that on every
+ * processor, the vertices of every locally owned and ghost cell is
+ * consistent with the corresponding location of these cells on other
+ * processors. On the other hand, the locations of artificial cells will
+ * in general be wrong since artificial cells may or may not exist on
+ * other processors and consequently it is not possible to determine
+ * their location in any way. This is not usually a problem since one
+ * never does anything on artificial cells. However, it may lead to
+ * problems if the mesh with moved vertices is refined in a later step.
+ * If that's what you want to do, the right way to do it is to save the
+ * offset applied to every vertex, call this function, and before
+ * refining or coarsening the mesh apply the opposite offset and call
+ * this function again.
+ *
+ * @param vertex_locally_moved A bitmap indicating which vertices have
+ * been moved. The size of this array must be equal to
+ * Triangulation::n_vertices() and must be a subset of those vertices
+ * flagged by GridTools::get_locally_owned_vertices().
+ *
+ * @see This function is used, for example, in
+ * GridTools::distort_random().
+ */
+ void
+ communicate_locally_moved_vertices(
+ const std::vector<bool> &vertex_locally_moved);
+
protected:
/**
* MPI communicator to be used for the triangulation. We create a unique
const DataType &)> & unpack,
const std::function<bool(const typename MeshType::active_cell_iterator &)>
&cell_filter =
- [](const typename MeshType::active_cell_iterator &) { return true; });
+ always_return<typename MeshType::active_cell_iterator, bool>{true});
/**
* Exchange arbitrary data of type @p DataType provided by the function
const typename MeshType::level_cell_iterator &)> &pack,
const std::function<void(const typename MeshType::level_cell_iterator &,
const DataType &)> & unpack,
- const std::function<bool(const typename MeshType::level_cell_iterator &)>
- &cell_filter =
- [](const typename MeshType::level_cell_iterator &) { return true; });
+ const std::function<bool(const typename MeshType::level_cell_iterator &)> &
+ cell_filter = always_return<typename MeshType::level_cell_iterator, bool>{
+ true});
/* Exchange with all processors of the MPI communicator @p mpi_communicator the vector of bounding
* boxes @p local_bboxes.
template <typename DataType,
typename MeshType,
typename MeshCellIteratorType>
- void
+ inline void
exchange_cell_data(
const MeshType &mesh,
const std::function<
} // namespace internal
template <typename DataType, typename MeshType>
- void
+ inline void
exchange_cell_data_to_ghosts(
const MeshType & mesh,
const std::function<std_cxx17::optional<DataType>(
template <typename DataType, typename MeshType>
- void
+ inline void
exchange_cell_data_to_level_ghosts(
const MeshType & mesh,
const std::function<std_cxx17::optional<DataType>(
- // This anonymous namespace contains utility for
- // the function Triangulation::communicate_locally_moved_vertices
- namespace CommunicateLocallyMovedVertices
- {
- namespace
- {
- /**
- * A list of tree+quadrant and their vertex indices.
- * The bool vector describes which vertices are of interest
- * and should be set on the receiving processes.
- */
- template <int dim, int spacedim>
- struct CellInfo
- {
- // store all the tree_indices we send/receive consecutively (n_cells
- // entries)
- std::vector<unsigned int> tree_index;
- // store all the quadrants we send/receive consecutively (n_cells
- // entries)
- std::vector<typename dealii::internal::p4est::types<dim>::quadrant>
- quadrants;
- // store for each cell the number of vertices we send/receive
- // and then the vertex indices (for each cell: n_vertices+1 entries)
- std::vector<unsigned int> vertex_indices;
- // store for each cell the vertices we send/receive
- // (for each cell n_vertices entries)
- std::vector<dealii::Point<spacedim>> vertices;
- // for receiving and unpacking data we need to store pointers to the
- // first vertex and vertex_index on each cell additionally
- // both vectors have as many entries as there are cells
- std::vector<unsigned int *> first_vertex_indices;
- std::vector<dealii::Point<spacedim> *> first_vertices;
-
- unsigned int
- bytes_for_buffer() const
- {
- return sizeof(unsigned int) +
- tree_index.size() * sizeof(unsigned int) +
- quadrants.size() *
- sizeof(
- typename dealii::internal::p4est::types<dim>::quadrant) +
- vertex_indices.size() * sizeof(unsigned int) +
- vertices.size() * sizeof(dealii::Point<spacedim>);
- }
-
- void
- pack_data(std::vector<char> &buffer) const
- {
- buffer.resize(bytes_for_buffer());
-
- char *ptr = buffer.data();
-
- const unsigned int num_cells = tree_index.size();
- std::memcpy(ptr, &num_cells, sizeof(unsigned int));
- ptr += sizeof(unsigned int);
-
- std::memcpy(ptr,
- tree_index.data(),
- num_cells * sizeof(unsigned int));
- ptr += num_cells * sizeof(unsigned int);
-
- std::memcpy(
- ptr,
- quadrants.data(),
- num_cells *
- sizeof(typename dealii::internal::p4est::types<dim>::quadrant));
- ptr +=
- num_cells *
- sizeof(typename dealii::internal::p4est::types<dim>::quadrant);
-
- std::memcpy(ptr,
- vertex_indices.data(),
- vertex_indices.size() * sizeof(unsigned int));
- ptr += vertex_indices.size() * sizeof(unsigned int);
- std::memcpy(ptr,
- vertices.data(),
- vertices.size() * sizeof(dealii::Point<spacedim>));
- ptr += vertices.size() * sizeof(dealii::Point<spacedim>);
-
- Assert(ptr == buffer.data() + buffer.size(), ExcInternalError());
- }
-
- void
- unpack_data(const std::vector<char> &buffer)
- {
- const char * ptr = buffer.data();
- unsigned int cells;
- memcpy(&cells, ptr, sizeof(unsigned int));
- ptr += sizeof(unsigned int);
-
- tree_index.resize(cells);
- memcpy(tree_index.data(), ptr, sizeof(unsigned int) * cells);
- ptr += sizeof(unsigned int) * cells;
-
- quadrants.resize(cells);
- memcpy(quadrants.data(),
- ptr,
- sizeof(
- typename dealii::internal::p4est::types<dim>::quadrant) *
- cells);
- ptr +=
- sizeof(typename dealii::internal::p4est::types<dim>::quadrant) *
- cells;
-
- vertex_indices.clear();
- first_vertex_indices.resize(cells);
- std::vector<unsigned int> n_vertices_on_cell(cells);
- std::vector<unsigned int> first_indices(cells);
- for (unsigned int c = 0; c < cells; ++c)
- {
- // The first 'vertex index' is the number of vertices.
- // Additionally, we need to store the pointer to this
- // vertex index with respect to the std::vector
- const unsigned int *const vertex_index =
- reinterpret_cast<const unsigned int *>(ptr);
- first_indices[c] = vertex_indices.size();
- vertex_indices.push_back(*vertex_index);
- n_vertices_on_cell[c] = *vertex_index;
- ptr += sizeof(unsigned int);
- // Now copy all the 'real' vertex_indices
- vertex_indices.resize(vertex_indices.size() +
- n_vertices_on_cell[c]);
- memcpy(&vertex_indices[vertex_indices.size() -
- n_vertices_on_cell[c]],
- ptr,
- n_vertices_on_cell[c] * sizeof(unsigned int));
- ptr += n_vertices_on_cell[c] * sizeof(unsigned int);
- }
- for (unsigned int c = 0; c < cells; ++c)
- first_vertex_indices[c] = &vertex_indices[first_indices[c]];
-
- vertices.clear();
- first_vertices.resize(cells);
- for (unsigned int c = 0; c < cells; ++c)
- {
- first_indices[c] = vertices.size();
- vertices.resize(vertices.size() + n_vertices_on_cell[c]);
- memcpy(&vertices[vertices.size() - n_vertices_on_cell[c]],
- ptr,
- n_vertices_on_cell[c] * sizeof(dealii::Point<spacedim>));
- ptr += n_vertices_on_cell[c] * sizeof(dealii::Point<spacedim>);
- }
- for (unsigned int c = 0; c < cells; ++c)
- first_vertices[c] = &vertices[first_indices[c]];
-
- Assert(ptr == buffer.data() + buffer.size(), ExcInternalError());
- }
- };
-
-
-
- // This function is responsible for gathering the information
- // we want to send to each process.
- // For each dealii cell on the coarsest level the corresponding
- // p4est_cell has to be provided when calling this function.
- // By recursing through all children we consider each active cell.
- // vertices_with_ghost_neighbors tells us which vertices
- // are in the ghost layer and for which processes they might
- // be interesting.
- // Whether a vertex has actually been updated locally is
- // stored in vertex_locally_moved. Only those are considered
- // for sending.
- // The gathered information is saved into needs_to_get_cell.
- template <int dim, int spacedim>
- void
- fill_vertices_recursively(
- const typename parallel::distributed::Triangulation<dim, spacedim>
- & tria,
- const unsigned int tree_index,
- const typename Triangulation<dim, spacedim>::cell_iterator
- &dealii_cell,
- const typename dealii::internal::p4est::types<dim>::quadrant
- &p4est_cell,
- const std::map<unsigned int, std::set<dealii::types::subdomain_id>>
- & vertices_with_ghost_neighbors,
- const std::vector<bool> &vertex_locally_moved,
- std::map<dealii::types::subdomain_id, CellInfo<dim, spacedim>>
- &needs_to_get_cell)
- {
- // see if we have to
- // recurse...
- if (dealii_cell->has_children())
- {
- typename dealii::internal::p4est::types<dim>::quadrant
- p4est_child[GeometryInfo<dim>::max_children_per_cell];
- dealii::internal::p4est::init_quadrant_children<dim>(p4est_cell,
- p4est_child);
-
-
- for (unsigned int c = 0;
- c < GeometryInfo<dim>::max_children_per_cell;
- ++c)
- fill_vertices_recursively<dim, spacedim>(
- tria,
- tree_index,
- dealii_cell->child(c),
- p4est_child[c],
- vertices_with_ghost_neighbors,
- vertex_locally_moved,
- needs_to_get_cell);
- return;
- }
-
- // We're at a leaf cell. If the cell is locally owned, we may
- // have to send its vertices to other processors if any of
- // its vertices is adjacent to a ghost cell and has been moved
- //
- // If one of the vertices of the cell is interesting,
- // send all moved vertices of the cell to all processors
- // adjacent to all cells adjacent to this vertex
- if (dealii_cell->is_locally_owned())
- {
- std::set<dealii::types::subdomain_id> send_to;
- for (const unsigned int v : GeometryInfo<dim>::vertex_indices())
- {
- const std::map<unsigned int,
- std::set<dealii::types::subdomain_id>>::
- const_iterator neighbor_subdomains_of_vertex =
- vertices_with_ghost_neighbors.find(
- dealii_cell->vertex_index(v));
-
- if (neighbor_subdomains_of_vertex !=
- vertices_with_ghost_neighbors.end())
- {
- Assert(neighbor_subdomains_of_vertex->second.size() != 0,
- ExcInternalError());
- send_to.insert(
- neighbor_subdomains_of_vertex->second.begin(),
- neighbor_subdomains_of_vertex->second.end());
- }
- }
-
- if (send_to.size() > 0)
- {
- std::vector<unsigned int> vertex_indices;
- std::vector<dealii::Point<spacedim>> local_vertices;
- for (const unsigned int v :
- GeometryInfo<dim>::vertex_indices())
- if (vertex_locally_moved[dealii_cell->vertex_index(v)])
- {
- vertex_indices.push_back(v);
- local_vertices.push_back(dealii_cell->vertex(v));
- }
-
- if (vertex_indices.size() > 0)
- for (const auto subdomain : send_to)
- {
- // get an iterator to what needs to be sent to that
- // subdomain (if already exists), or create such an
- // object
- const typename std::map<
- dealii::types::subdomain_id,
- CellInfo<dim, spacedim>>::iterator p =
- needs_to_get_cell
- .insert(std::make_pair(subdomain,
- CellInfo<dim, spacedim>()))
- .first;
-
- p->second.tree_index.push_back(tree_index);
- p->second.quadrants.push_back(p4est_cell);
-
- p->second.vertex_indices.push_back(
- vertex_indices.size());
- p->second.vertex_indices.insert(
- p->second.vertex_indices.end(),
- vertex_indices.begin(),
- vertex_indices.end());
-
- p->second.vertices.insert(p->second.vertices.end(),
- local_vertices.begin(),
- local_vertices.end());
- }
- }
- }
- }
-
-
-
- // After the cell data has been received this function is responsible
- // for moving the vertices in the corresponding ghost layer locally.
- // As in fill_vertices_recursively for each dealii cell on the
- // coarsest level the corresponding p4est_cell has to be provided
- // when calling this function. By recursing through through all
- // children we consider each active cell.
- // Additionally, we need to give a pointer to the first vertex indices
- // and vertices. Since the first information saved in vertex_indices
- // is the number of vertices this all the information we need.
- template <int dim, int spacedim>
- void
- set_vertices_recursively(
- const parallel::distributed::Triangulation<dim, spacedim> &tria,
- const typename dealii::internal::p4est::types<dim>::quadrant
- &p4est_cell,
- const typename Triangulation<dim, spacedim>::cell_iterator
- &dealii_cell,
- const typename dealii::internal::p4est::types<dim>::quadrant
- & quadrant,
- const dealii::Point<spacedim> *const vertices,
- const unsigned int *const vertex_indices)
- {
- if (dealii::internal::p4est::quadrant_is_equal<dim>(p4est_cell,
- quadrant))
- {
- Assert(!dealii_cell->is_artificial(), ExcInternalError());
- Assert(dealii_cell->is_active(), ExcInternalError());
- Assert(!dealii_cell->is_locally_owned(), ExcInternalError());
-
- const unsigned int n_vertices = vertex_indices[0];
-
- // update dof indices of cell
- for (unsigned int i = 0; i < n_vertices; ++i)
- dealii_cell->vertex(vertex_indices[i + 1]) = vertices[i];
-
- return;
- }
-
- if (dealii_cell->is_active())
- return;
-
- if (!dealii::internal::p4est::quadrant_is_ancestor<dim>(p4est_cell,
- quadrant))
- return;
-
- typename dealii::internal::p4est::types<dim>::quadrant
- p4est_child[GeometryInfo<dim>::max_children_per_cell];
- dealii::internal::p4est::init_quadrant_children<dim>(p4est_cell,
- p4est_child);
-
- for (unsigned int c = 0; c < GeometryInfo<dim>::max_children_per_cell;
- ++c)
- set_vertices_recursively<dim, spacedim>(tria,
- p4est_child[c],
- dealii_cell->child(c),
- quadrant,
- vertices,
- vertex_indices);
- }
- } // namespace
- } // namespace CommunicateLocallyMovedVertices
-
-
-
template <int dim, int spacedim>
void
Triangulation<dim, spacedim>::clear()
- template <int dim, int spacedim>
- void
- Triangulation<dim, spacedim>::communicate_locally_moved_vertices(
- const std::vector<bool> &vertex_locally_moved)
- {
- Assert(vertex_locally_moved.size() == this->n_vertices(),
- ExcDimensionMismatch(vertex_locally_moved.size(),
- this->n_vertices()));
-# ifdef DEBUG
- {
- const std::vector<bool> locally_owned_vertices =
- dealii::GridTools::get_locally_owned_vertices(*this);
- for (unsigned int i = 0; i < locally_owned_vertices.size(); ++i)
- Assert((vertex_locally_moved[i] == false) ||
- (locally_owned_vertices[i] == true),
- ExcMessage("The vertex_locally_moved argument must not "
- "contain vertices that are not locally owned"));
- }
-# endif
-
- // First find out which process should receive which vertices.
- // These are specifically the ones that are located on cells at the
- // boundary of the subdomain this process owns and the receiving
- // process taking periodic faces into account.
- // Here, it is sufficient to collect all vertices that are located
- // at that boundary.
- const std::map<unsigned int, std::set<dealii::types::subdomain_id>>
- vertices_with_ghost_neighbors =
- GridTools::compute_vertices_with_ghost_neighbors(*this);
-
- // now collect cells and their vertices
- // for the interested neighbors
- using cellmap_t =
- std::map<dealii::types::subdomain_id,
- CommunicateLocallyMovedVertices::CellInfo<dim, spacedim>>;
- cellmap_t needs_to_get_cells;
-
- for (typename Triangulation<dim, spacedim>::cell_iterator cell =
- this->begin(0);
- cell != this->end(0);
- ++cell)
- {
- typename dealii::internal::p4est::types<dim>::quadrant
- p4est_coarse_cell;
- dealii::internal::p4est::init_coarse_quadrant<dim>(p4est_coarse_cell);
-
- CommunicateLocallyMovedVertices::fill_vertices_recursively<dim,
- spacedim>(
- *this,
- this->get_coarse_cell_to_p4est_tree_permutation()[cell->index()],
- cell,
- p4est_coarse_cell,
- vertices_with_ghost_neighbors,
- vertex_locally_moved,
- needs_to_get_cells);
- }
-
- // Send information
-
- // We need to protect this communication below using a mutex:
- static Utilities::MPI::CollectiveMutex mutex;
- Utilities::MPI::CollectiveMutex::ScopedLock lock(
- mutex, this->get_communicator());
-
- const int mpi_tag = Utilities::MPI::internal::Tags::
- triangulation_communicate_locally_moved_vertices;
-
- std::vector<std::vector<char>> sendbuffers(needs_to_get_cells.size());
- std::vector<std::vector<char>>::iterator buffer = sendbuffers.begin();
- std::vector<MPI_Request> requests(needs_to_get_cells.size());
- std::vector<unsigned int> destinations;
-
- unsigned int idx = 0;
-
- for (typename cellmap_t::iterator it = needs_to_get_cells.begin();
- it != needs_to_get_cells.end();
- ++it, ++buffer, ++idx)
- {
- const unsigned int num_cells = it->second.tree_index.size();
- (void)num_cells;
- destinations.push_back(it->first);
-
- Assert(num_cells == it->second.quadrants.size(), ExcInternalError());
- Assert(num_cells > 0, ExcInternalError());
-
- // pack all the data into
- // the buffer for this
- // recipient and send
- // it. keep data around
- // till we can make sure
- // that the packet has been
- // received
- it->second.pack_data(*buffer);
- const int ierr = MPI_Isend(buffer->data(),
- buffer->size(),
- MPI_BYTE,
- it->first,
- mpi_tag,
- this->get_communicator(),
- &requests[idx]);
- AssertThrowMPI(ierr);
- }
-
- Assert(destinations.size() == needs_to_get_cells.size(),
- ExcInternalError());
-
- // collect the neighbors
- // that are going to send stuff to us
- const unsigned int n_senders =
- Utilities::MPI::compute_n_point_to_point_communications(
- this->get_communicator(), destinations);
-
- // receive ghostcelldata
- std::vector<char> receive;
- CommunicateLocallyMovedVertices::CellInfo<dim, spacedim> cellinfo;
- for (unsigned int i = 0; i < n_senders; ++i)
- {
- MPI_Status status;
- int ierr = MPI_Probe(MPI_ANY_SOURCE,
- mpi_tag,
- this->get_communicator(),
- &status);
- AssertThrowMPI(ierr);
-
- int len;
- ierr = MPI_Get_count(&status, MPI_BYTE, &len);
- AssertThrowMPI(ierr);
- receive.resize(len);
-
- char *ptr = receive.data();
- ierr = MPI_Recv(ptr,
- len,
- MPI_BYTE,
- status.MPI_SOURCE,
- status.MPI_TAG,
- this->get_communicator(),
- &status);
- AssertThrowMPI(ierr);
-
- cellinfo.unpack_data(receive);
- const unsigned int cells = cellinfo.tree_index.size();
- for (unsigned int c = 0; c < cells; ++c)
- {
- typename dealii::parallel::distributed::
- Triangulation<dim, spacedim>::cell_iterator cell(
- this,
- 0,
- this->get_p4est_tree_to_coarse_cell_permutation()
- [cellinfo.tree_index[c]]);
-
- typename dealii::internal::p4est::types<dim>::quadrant
- p4est_coarse_cell;
- dealii::internal::p4est::init_coarse_quadrant<dim>(
- p4est_coarse_cell);
-
- CommunicateLocallyMovedVertices::set_vertices_recursively<
- dim,
- spacedim>(*this,
- p4est_coarse_cell,
- cell,
- cellinfo.quadrants[c],
- cellinfo.first_vertices[c],
- cellinfo.first_vertex_indices[c]);
- }
- }
-
- // complete all sends, so that we can
- // safely destroy the buffers.
- if (requests.size() > 0)
- {
- const int ierr =
- MPI_Waitall(requests.size(), requests.data(), MPI_STATUSES_IGNORE);
- AssertThrowMPI(ierr);
- }
-
- // check all msgs got sent and received
- Assert(Utilities::MPI::sum(needs_to_get_cells.size(),
- this->get_communicator()) ==
- Utilities::MPI::sum(n_senders, this->get_communicator()),
- ExcInternalError());
- }
-
-
-
template <int dim, int spacedim>
const std::vector<types::global_dof_index> &
Triangulation<dim, spacedim>::get_p4est_tree_to_coarse_cell_permutation()
- template <int spacedim>
- void
- Triangulation<1, spacedim>::communicate_locally_moved_vertices(
- const std::vector<bool> & /*vertex_locally_moved*/)
- {
- 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 dim, int spacedim>
+ void
+ TriangulationBase<dim, spacedim>::communicate_locally_moved_vertices(
+ const std::vector<bool> &vertex_locally_moved)
+ {
+ AssertDimension(vertex_locally_moved.size(), this->n_vertices());
+#ifdef DEBUG
+ {
+ const std::vector<bool> locally_owned_vertices =
+ dealii::GridTools::get_locally_owned_vertices(*this);
+ for (unsigned int i = 0; i < locally_owned_vertices.size(); ++i)
+ Assert((vertex_locally_moved[i] == false) ||
+ (locally_owned_vertices[i] == true),
+ ExcMessage("The vertex_locally_moved argument must not "
+ "contain vertices that are not locally owned"));
+ }
+#endif
+
+ Point<spacedim> invalid_point;
+ for (int d = 0; d < spacedim; ++d)
+ invalid_point[d] = std::numeric_limits<double>::quiet_NaN();
+
+ const auto pack = [&](const auto &cell) {
+ std::vector<Point<spacedim>> vertices(cell->n_vertices());
+
+ for (const auto v : cell->vertex_indices())
+ if (vertex_locally_moved[cell->vertex_index(v)])
+ vertices[v] = cell->vertex(v);
+ else
+ vertices[v] = invalid_point;
+
+ return vertices;
+ };
+
+ const auto unpack = [&](const auto &cell, const auto &vertices) {
+ for (const auto v : cell->vertex_indices())
+ if (std::isnan(vertices[v][0]) == false)
+ cell->vertex(v) = vertices[v];
+ };
+
+ if (this->is_multilevel_hierarchy_constructed())
+ GridTools::exchange_cell_data_to_level_ghosts<
+ std::vector<Point<spacedim>>>(*this, pack, unpack);
+ else
+ GridTools::exchange_cell_data_to_ghosts<std::vector<Point<spacedim>>>(
+ *this, pack, unpack);
+ }
+
+
+
template <int dim, int spacedim>
const Utilities::MPI::Partitioner &
TriangulationBase<dim, spacedim>::global_active_cell_index_partitioner() const
// If the triangulation is distributed, we need to
// exchange the moved vertices across mpi processes
- if (parallel::distributed::Triangulation<dim, spacedim>
- *distributed_triangulation =
- dynamic_cast<parallel::distributed::Triangulation<dim, spacedim> *>(
- &triangulation))
+ if (auto distributed_triangulation =
+ dynamic_cast<parallel::DistributedTriangulationBase<dim, spacedim> *>(
+ &triangulation))
{
const std::vector<bool> locally_owned_vertices =
get_locally_owned_vertices(triangulation);
}
}
-#ifdef DEAL_II_WITH_P4EST
distributed_triangulation->communicate_locally_moved_vertices(
locally_owned_vertices);
-#else
- (void)distributed_triangulation;
- Assert(false, ExcInternalError());
-#endif
}
else
// if this is a sequential triangulation, we could in principle
vertex_to_cell_map(triangulation);
// Create a local index for the locally "owned" vertices
- types::global_vertex_index next_index = 0;
- unsigned int max_cellid_size = 0;
+ types::global_vertex_index next_index = 0;
+ unsigned int max_cellid_size = 0;
std::set<std::pair<types::subdomain_id, types::global_vertex_index>>
- vertices_added;
+ vertices_added;
std::map<types::subdomain_id, std::set<unsigned int>> vertices_to_recv;
std::map<types::subdomain_id,
std::vector<std::tuple<types::global_vertex_index,
types::global_vertex_index,
std::string>>>
- vertices_to_send;
+ vertices_to_send;
active_cell_iterator cell = triangulation.begin_active(),
endc = triangulation.end();
std::set<active_cell_iterator> missing_vert_cells;
- std::set<unsigned int> used_vertex_index;
+ std::set<unsigned int> used_vertex_index;
for (; cell != endc; ++cell)
{
if (cell->is_locally_owned())
types::subdomain_id lowest_subdomain_id = cell->subdomain_id();
typename std::set<active_cell_iterator>::iterator
adjacent_cell = vertex_to_cell[cell->vertex_index(i)].begin(),
- end_adj_cell = vertex_to_cell[cell->vertex_index(i)].end();
+ end_adj_cell = vertex_to_cell[cell->vertex_index(i)].end();
for (; adjacent_cell != end_adj_cell; ++adjacent_cell)
lowest_subdomain_id =
std::min(lowest_subdomain_id,
// Make indices global by getting the number of vertices owned by each
// processors and shifting the indices accordingly
types::global_vertex_index shift = 0;
- int ierr = MPI_Exscan(&next_index,
+ int ierr = MPI_Exscan(&next_index,
&shift,
1,
DEAL_II_VERTEX_INDEX_MPI_TYPE,
AssertThrowMPI(ierr);
std::map<unsigned int, types::global_vertex_index>::iterator
- global_index_it = local_to_global_vertex_index.begin(),
+ global_index_it = local_to_global_vertex_index.begin(),
global_index_end = local_to_global_vertex_index.end();
for (; global_index_it != global_index_end; ++global_index_it)
global_index_it->second += shift;
std::vector<std::tuple<types::global_vertex_index,
types::global_vertex_index,
std::string>>>::iterator
- vert_to_send_it = vertices_to_send.begin(),
+ vert_to_send_it = vertices_to_send.begin(),
vert_to_send_end = vertices_to_send.end();
for (unsigned int i = 0; vert_to_send_it != vert_to_send_end;
++vert_to_send_it, ++i)
{
- int destination = vert_to_send_it->first;
- const unsigned int n_vertices = vert_to_send_it->second.size();
- const int buffer_size = 2 * n_vertices;
+ int destination = vert_to_send_it->first;
+ const unsigned int n_vertices = vert_to_send_it->second.size();
+ const int buffer_size = 2 * n_vertices;
vertices_send_buffers[i].resize(buffer_size);
// fill the buffer
std::vector<std::vector<types::global_vertex_index>> vertices_recv_buffers(
vertices_to_recv.size());
typename std::map<types::subdomain_id, std::set<unsigned int>>::iterator
- vert_to_recv_it = vertices_to_recv.begin(),
+ vert_to_recv_it = vertices_to_recv.begin(),
vert_to_recv_end = vertices_to_recv.end();
for (unsigned int i = 0; vert_to_recv_it != vert_to_recv_end;
++vert_to_recv_it, ++i)
{
- int source = vert_to_recv_it->first;
- const unsigned int n_vertices = vert_to_recv_it->second.size();
- const int buffer_size = 2 * n_vertices;
+ int source = vert_to_recv_it->first;
+ const unsigned int n_vertices = vert_to_recv_it->second.size();
+ const int buffer_size = 2 * n_vertices;
vertices_recv_buffers[i].resize(buffer_size);
// Receive the message
for (unsigned int i = 0; vert_to_send_it != vert_to_send_end;
++vert_to_send_it, ++i)
{
- int destination = vert_to_send_it->first;
- const unsigned int n_vertices = vert_to_send_it->second.size();
- const int buffer_size = max_cellid_size * n_vertices;
+ int destination = vert_to_send_it->first;
+ const unsigned int n_vertices = vert_to_send_it->second.size();
+ const int buffer_size = max_cellid_size * n_vertices;
cellids_send_buffers[i].resize(buffer_size);
// fill the buffer
for (unsigned int i = 0; vert_to_recv_it != vert_to_recv_end;
++vert_to_recv_it, ++i)
{
- int source = vert_to_recv_it->first;
- const unsigned int n_vertices = vert_to_recv_it->second.size();
- const int buffer_size = max_cellid_size * n_vertices;
+ int source = vert_to_recv_it->first;
+ const unsigned int n_vertices = vert_to_recv_it->second.size();
+ const int buffer_size = max_cellid_size * n_vertices;
cellids_recv_buffers[i].resize(buffer_size);
// Receive the message
&cellids_recv_buffers[i][max_cellid_size * j] + max_cellid_size);
bool found = false;
typename std::set<active_cell_iterator>::iterator
- cell_set_it = missing_vert_cells.begin(),
+ cell_set_it = missing_vert_cells.begin(),
end_cell_set = missing_vert_cells.end();
for (; (found == false) && (cell_set_it != end_cell_set);
++cell_set_it)
// are owned by other processors -- either because the vertex is
// on an artificial cell, or because it is on a ghost cell with
// a smaller subdomain
- if (const parallel::distributed::Triangulation<dim, spacedim> *tr =
- dynamic_cast<const parallel::distributed::Triangulation<dim, spacedim>
- *>(&triangulation))
+ if (const auto *tr = dynamic_cast<
+ const parallel::DistributedTriangulationBase<dim, spacedim> *>(
+ &triangulation))
for (const auto &cell : triangulation.active_cell_iterators())
if (cell->is_artificial() ||
(cell->is_ghost() &&
guessed_points = GridTools::guess_point_owner(global_bboxes, local_points);
// Preparing to call compute_point_locations on points which may be local
- const auto &guess_loc_idx = std::get<0>(guessed_points)[my_rank];
+ const auto & guess_loc_idx = std::get<0>(guessed_points)[my_rank];
const unsigned int n_local_guess = guess_loc_idx.size();
// Vector containing points which are probably local
for (const auto &probable_owners : points_to_probable_owners)
{
- const auto &point_idx = probable_owners.first;
+ const auto &point_idx = probable_owners.first;
const auto &probable_owner_ranks = probable_owners.second;
if (!std::binary_search(found_point_indices.begin(),
found_point_indices.end(),
for (const auto &received_ghost_point : received_ghost_points)
{
// Transforming CellsIds into Tria iterators
- const auto &cell_ids = std::get<0>(received_ghost_point.second);
- const unsigned int n_cells = cell_ids.size();
+ const auto & cell_ids = std::get<0>(received_ghost_point.second);
+ const unsigned int n_cells = cell_ids.size();
std::vector<typename Triangulation<dim, spacedim>::active_cell_iterator>
cell_iter(n_cells);
for (unsigned int c = 0; c < n_cells; ++c)
// Copying data from the unordered map to the tuple
// and returning output
const unsigned int size_output = found_points.size();
- auto &out_cells = std::get<0>(output_tuple);
- auto &out_qpoints = std::get<1>(output_tuple);
- auto &out_maps = std::get<2>(output_tuple);
- auto &out_points = std::get<3>(output_tuple);
- auto &out_ranks = std::get<4>(output_tuple);
+ auto & out_cells = std::get<0>(output_tuple);
+ auto & out_qpoints = std::get<1>(output_tuple);
+ auto & out_maps = std::get<2>(output_tuple);
+ auto & out_points = std::get<3>(output_tuple);
+ auto & out_ranks = std::get<4>(output_tuple);
out_cells.resize(size_output);
out_qpoints.resize(size_output);
unsigned int c = 0;
for (const auto &cell_and_data : found_points)
{
- out_cells[c] = cell_and_data.first;
+ out_cells[c] = cell_and_data.first;
out_qpoints[c] = std::get<0>(cell_and_data.second);
- out_maps[c] = std::get<1>(cell_and_data.second);
- out_points[c] = std::get<2>(cell_and_data.second);
- out_ranks[c] = std::get<3>(cell_and_data.second);
+ out_maps[c] = std::get<1>(cell_and_data.second);
+ out_points[c] = std::get<2>(cell_and_data.second);
+ out_ranks[c] = std::get<3>(cell_and_data.second);
++c;
}
// Step 4: create the array of bboxes for output
std::vector<std::vector<BoundingBox<spacedim>>> global_bboxes(n_procs);
- unsigned int begin_idx = 0;
+ unsigned int begin_idx = 0;
for (unsigned int i = 0; i < n_procs; ++i)
{
// Number of local bounding boxes
projects / dealii.git / commitdiff