Reference documentation for deal.II version 9.5.0
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
Public Member Functions | Static Public Member Functions | Private Member Functions | Private Attributes | List of all members
Utilities::MPI::Partitioner Class Reference

#include <deal.II/base/partitioner.h>

Inheritance diagram for Utilities::MPI::Partitioner:
[legend]

Public Member Functions

 Partitioner ()
 
 Partitioner (const unsigned int size)
 
 Partitioner (const types::global_dof_index local_size, const types::global_dof_index ghost_size, const MPI_Comm communicator)
 
 Partitioner (const IndexSet &locally_owned_indices, const IndexSet &ghost_indices_in, const MPI_Comm communicator_in)
 
 Partitioner (const IndexSet &locally_owned_indices, const MPI_Comm communicator_in)
 
virtual void reinit (const IndexSet &vector_space_vector_index_set, const IndexSet &read_write_vector_index_set, const MPI_Comm communicator) override
 
void set_owned_indices (const IndexSet &locally_owned_indices)
 
void set_ghost_indices (const IndexSet &ghost_indices, const IndexSet &larger_ghost_index_set=IndexSet())
 
types::global_dof_index size () const
 
unsigned int local_size () const
 
unsigned int locally_owned_size () const
 
const IndexSetlocally_owned_range () const
 
std::pair< types::global_dof_index, types::global_dof_indexlocal_range () const
 
bool in_local_range (const types::global_dof_index global_index) const
 
unsigned int global_to_local (const types::global_dof_index global_index) const
 
types::global_dof_index local_to_global (const unsigned int local_index) const
 
bool is_ghost_entry (const types::global_dof_index global_index) const
 
const IndexSetghost_indices () const
 
unsigned int n_ghost_indices () const
 
const std::vector< std::pair< unsigned int, unsigned int > > & ghost_indices_within_larger_ghost_set () const
 
const std::vector< std::pair< unsigned int, unsigned int > > & ghost_targets () const
 
const std::vector< std::pair< unsigned int, unsigned int > > & import_indices () const
 
unsigned int n_import_indices () const
 
const std::vector< std::pair< unsigned int, unsigned int > > & import_targets () const
 
bool is_compatible (const Partitioner &part) const
 
bool is_globally_compatible (const Partitioner &part) const
 
unsigned int this_mpi_process () const
 
unsigned int n_mpi_processes () const
 
virtual MPI_Comm get_mpi_communicator () const override
 
bool ghost_indices_initialized () const
 
template<typename Number , typename MemorySpaceType = MemorySpace::Host>
void export_to_ghosted_array_start (const unsigned int communication_channel, const ArrayView< const Number, MemorySpaceType > &locally_owned_array, const ArrayView< Number, MemorySpaceType > &temporary_storage, const ArrayView< Number, MemorySpaceType > &ghost_array, std::vector< MPI_Request > &requests) const
 
template<typename Number , typename MemorySpaceType = MemorySpace::Host>
void export_to_ghosted_array_finish (const ArrayView< Number, MemorySpaceType > &ghost_array, std::vector< MPI_Request > &requests) const
 
template<typename Number , typename MemorySpaceType = MemorySpace::Host>
void import_from_ghosted_array_start (const VectorOperation::values vector_operation, const unsigned int communication_channel, const ArrayView< Number, MemorySpaceType > &ghost_array, const ArrayView< Number, MemorySpaceType > &temporary_storage, std::vector< MPI_Request > &requests) const
 
template<typename Number , typename MemorySpaceType = MemorySpace::Host>
void import_from_ghosted_array_finish (const VectorOperation::values vector_operation, const ArrayView< const Number, MemorySpaceType > &temporary_storage, const ArrayView< Number, MemorySpaceType > &locally_owned_storage, const ArrayView< Number, MemorySpaceType > &ghost_array, std::vector< MPI_Request > &requests) const
 
std::size_t memory_consumption () const
 

Static Public Member Functions

static ::ExceptionBaseExcIndexNotPresent (types::global_dof_index arg1, unsigned int arg2)
 
static ::ExceptionBaseExcGhostIndexArrayHasWrongSize (unsigned int arg1, unsigned int arg2, unsigned int arg3)
 

Private Member Functions

void initialize_import_indices_plain_dev () const
 

Private Attributes

types::global_dof_index global_size
 
IndexSet locally_owned_range_data
 
std::pair< types::global_dof_index, types::global_dof_indexlocal_range_data
 
IndexSet ghost_indices_data
 
unsigned int n_ghost_indices_data
 
std::vector< std::pair< unsigned int, unsigned int > > ghost_targets_data
 
std::vector< std::pair< unsigned int, unsigned int > > import_indices_data
 
std::vector< Kokkos::View< unsigned int *, MemorySpace::Default::kokkos_space > > import_indices_plain_dev
 
unsigned int n_import_indices_data
 
std::vector< std::pair< unsigned int, unsigned int > > import_targets_data
 
std::vector< unsigned intimport_indices_chunks_by_rank_data
 
unsigned int n_ghost_indices_in_larger_set
 
std::vector< unsigned intghost_indices_subset_chunks_by_rank_data
 
std::vector< std::pair< unsigned int, unsigned int > > ghost_indices_subset_data
 
unsigned int my_pid
 
unsigned int n_procs
 
MPI_Comm communicator
 
bool have_ghost_indices
 

Detailed Description

This class defines a model for the partitioning of a vector (or, in fact, any linear data structure) among processors using MPI.

The partitioner stores the global vector size and the locally owned range as a half-open interval [lower, upper) on each process. In other words, it assumes that every process stores a contiguous subset of the array mentioned in the documentation of the base class Utilities::MPI::CommunicationPatternBase. (If you want to store non-contiguous parts of these arrays on each process, take a look at Utilities::MPI::NoncontiguousPartitioner.)

Furthermore, it includes a structure for the point-to-point communication patterns. It allows the inclusion of ghost indices (i.e. indices that a current processor needs to have access to, but are owned by another process) through an IndexSet. In addition, it also stores the other processors' ghost indices belonging to the current processor (see import_targets()), which are the indices where other processors might require information from. In a sense, these import indices form the dual of the ghost indices. This information is gathered once when constructing the partitioner, which obviates subsequent global communication steps when exchanging data.

The figure below gives an example of index space \([0,74)\) being split into four parts that are each owned by one MPI process:

The first row (above the thick black line) shows which process owns which elements. Below it, the next four lines indicate which elements of the overall array each processor wants to know about – this is generally a superset of the locally owned elements, with the difference being what are called "ghost elements".

To understand the remaining pieces of the figure (and this class), remember that in MPI, you can't just ask another process for data. (That's not quite true: There are mechanisms in newer MPI standards for this, but as a general rule it's true.) Rather, if you need information, you need to send another process a message, the other process needs to expect the message and respond as appropriate with a message of its own. In practice, it is therefore easier and faster if each process will already know what it will be asked and, at the appropriate time, just send that data. The remaining lines of information set up this kind of scheme.

To this end, note that process 0 will want to know about five elements it does not own itself: those with indices 20, 21 (owned by process 1); and 40, 41, 43 (owned by process 2). Similar information can be obtained from the following lines. To satisfy this need for knowledge, it would therefore be quite useful if process 1 stored that, at the appropriate time, it will have to send elements 20, 21 to process 0. Similarly, if you go through lines 2-4 below the thick black line, you will see that process 0 should know that it will need to send elements 1, 2, 13, 18, 19 to process 1; 18, 19 to process 2; and elements 1, 2, 13 to process 3. These are called "import indices" because other processes will want to import them. Instead of storing these indices as a set, it is often useful to use half-open index sets instead, and so the import indices listed above form the following collection of sets: [1,3), [13,14), [18,20), [18,20), [1,3), [13,14). This is how the import indices are shown in the figure above. We now only have to know which of these half-open sets are to be sent to whom. This is done in the line above it where we list the "import targets" (i.e., the target processes for an import operations): Process 1 will receive 5 elements (which are comprised of the first three half-open target index sets), process 2 will receive 2 indices (the fourth half-open interval), and process 3 will receive 3 indices (the remaining two half-open intervals). This information is encoded as the pairs {1,5}, {2,2}, {3,3} as the import targets. Similar considerations can be made for what processes 1, 2, and 3 will have to send out.

Finally, when receiving information, it is useful to know how many indices each process will receive from whom since then one can already pre-allocate buffers of the right size. This is listed in the last line under "ghost targets": Process 0 will receive two elements from process 1 (namely those with indices 20, 21), and three from process 2 (namely those with indices 40, 41, 43). This is encoded as pairs {1,2} and {2,3}. Again, similar considerations can be made for what processes 1, 2, and 3 should expect, and what is then shown in their respective columns.

The main purpose of this class is to set up these data structures knowing only which process owns which elements, and for which additional ghost elements each process needs knowledge.

Local and global numbering

The partitioner includes a mechanism for converting global to local and local to global indices. Internally, this class stores vector elements using the convention as follows: The local range is associated with local indices [0, locally_owned_size()), and ghost indices are stored consecutively in [locally_owned_size(), locally_owned_size() + n_ghost_indices()). The ghost indices are sorted according to their global index.

Parallel data exchange

This class also handles the ghost data exchange for partitioned arrays of objects – i.e., if a separate class stores the locally owned elements on each process, then this class facilitates the importation of those elements that are locally needed as ghosts but stored elsewhere. An example of where this class is used is the LinearAlgebra::distributed::Vector class.

The data exchange happens through four functions:

The MPI communication routines are point-to-point communication patterns.

Sending only selected ghost data

This partitioner class operates on a fixed set of ghost indices and must always be compatible with the ghost indices inside the array whose partitioning it represents. In some cases, one only wants to send around some of the ghost indices present in a vector, but without creating a copy of the vector with a suitable index set - think e.g. of local time stepping where different regions of a vector might be exchanged at different stages of a time step slice. This class supports that case by the following model: A vector is first created with the full ghosted index set. Then, a second Partitioner instance is created that sets ghost indices with a tighter index set as ghosts, but specifying the larger index set as the second argument to the set_ghost_indices() call. When data is exchanged, the export_to_ghosted_array_start() and import_from_ghosted_array_start() detect this case and only send the selected indices, taken from the full array of ghost entries.

Definition at line 197 of file partitioner.h.

Constructor & Destructor Documentation

◆ Partitioner() [1/5]

Utilities::MPI::Partitioner::Partitioner ( )

Default constructor.

Definition at line 30 of file partitioner.cc.

◆ Partitioner() [2/5]

Utilities::MPI::Partitioner::Partitioner ( const unsigned int  size)

Constructor with size argument. Creates an MPI_COMM_SELF structure where there is no real parallel layout.

Definition at line 45 of file partitioner.cc.

◆ Partitioner() [3/5]

Utilities::MPI::Partitioner::Partitioner ( const types::global_dof_index  local_size,
const types::global_dof_index  ghost_size,
const MPI_Comm  communicator 
)

Constructor that takes the number of locally-owned degrees of freedom local_size and the number of ghost degrees of freedom ghost_size.

The local index range is translated to global indices in an ascending and one-to-one fashion, i.e., the indices of process \(p\) sit exactly between the indices of the processes \(p-1\) and \(p+1\), respectively.

Note
Setting the ghost_size variable to an appropriate value provides memory space for the ghost data in a vector's memory allocation as and allows access to it via local_element(). However, the associated global indices must be handled externally in this case.

Definition at line 65 of file partitioner.cc.

◆ Partitioner() [4/5]

Utilities::MPI::Partitioner::Partitioner ( const IndexSet locally_owned_indices,
const IndexSet ghost_indices_in,
const MPI_Comm  communicator_in 
)

Constructor with index set arguments. This constructor creates a distributed layout based on a given communicators, an IndexSet describing the locally owned range and another one for describing ghost indices that are owned by other processors, but that we need to have read or write access to.

Definition at line 101 of file partitioner.cc.

◆ Partitioner() [5/5]

Utilities::MPI::Partitioner::Partitioner ( const IndexSet locally_owned_indices,
const MPI_Comm  communicator_in 
)

Constructor with one index set argument. This constructor creates a distributed layout based on a given communicator, and an IndexSet describing the locally owned range. It allows to set the ghost indices at a later time. Apart from this, it is similar to the other constructor with two index sets.

Definition at line 120 of file partitioner.cc.

Member Function Documentation

◆ reinit()

void Utilities::MPI::Partitioner::reinit ( const IndexSet vector_space_vector_index_set,
const IndexSet read_write_vector_index_set,
const MPI_Comm  communicator 
)
overridevirtual

Reinitialize the communication pattern. The first argument vector_space_vector_index_set is the index set associated to a VectorSpaceVector object. The second argument read_write_vector_index_set is the index set associated to a ReadWriteVector object.

Implements Utilities::MPI::CommunicationPatternBase.

Definition at line 138 of file partitioner.cc.

◆ set_owned_indices()

void Utilities::MPI::Partitioner::set_owned_indices ( const IndexSet locally_owned_indices)

Set the locally owned indices. Used in the constructor.

Definition at line 151 of file partitioner.cc.

◆ set_ghost_indices()

void Utilities::MPI::Partitioner::set_ghost_indices ( const IndexSet ghost_indices,
const IndexSet larger_ghost_index_set = IndexSet() 
)

Set the ghost indices after the constructor has been called.

The optional parameter larger_ghost_index_set allows defining an indirect addressing into a larger set of ghost indices. This setup is useful if a distributed vector is based on that larger ghost index set but only a tighter subset should be communicated according to ghost_indices.

Definition at line 192 of file partitioner.cc.

◆ size()

types::global_dof_index Utilities::MPI::Partitioner::size ( ) const

Return the global size.

◆ local_size()

unsigned int Utilities::MPI::Partitioner::local_size ( ) const

Return the number of locally owned indices, i.e., local_range().second minus local_range().first. The returned numbers need to add up to the total number of indices when summed over all processes

Deprecated:
Use the more clearly named function locally_owned_size() instead.

◆ locally_owned_size()

unsigned int Utilities::MPI::Partitioner::locally_owned_size ( ) const

Return the number of locally owned indices, i.e., local_range().second minus local_range().first. The returned numbers need to add up to the total number of indices when summed over all processes

◆ locally_owned_range()

const IndexSet & Utilities::MPI::Partitioner::locally_owned_range ( ) const

Return an IndexSet representation of the local range. This class only supports contiguous local ranges, so the IndexSet actually only consists of one single range of data, and is equivalent to the result of local_range().

◆ local_range()

std::pair< types::global_dof_index, types::global_dof_index > Utilities::MPI::Partitioner::local_range ( ) const

Return the local range. The returned pair consists of the index of the first element and the index of the element one past the last locally owned one.

◆ in_local_range()

bool Utilities::MPI::Partitioner::in_local_range ( const types::global_dof_index  global_index) const

Return true if the given global index is in the local range of this processor.

◆ global_to_local()

unsigned int Utilities::MPI::Partitioner::global_to_local ( const types::global_dof_index  global_index) const

Return the local index corresponding to the given global index. If the given global index is neither locally owned nor a ghost, an exception is thrown.

Note that the returned local index for locally owned indices will be between 0 and locally_owned_size() - 1, and the local index for ghosts is between locally_owned_size() and locally_owned_size() + n_ghost_indices() - 1.

◆ local_to_global()

types::global_dof_index Utilities::MPI::Partitioner::local_to_global ( const unsigned int  local_index) const

Return the global index corresponding to the given local index.

Note that the local index for locally owned indices must be between 0 and locally_owned_size() - 1, and the local index for ghosts must be between locally_owned_size() and locally_owned_size() + n_ghost_indices() - 1.

◆ is_ghost_entry()

bool Utilities::MPI::Partitioner::is_ghost_entry ( const types::global_dof_index  global_index) const

Return whether the given global index is a ghost index on the present processor. Returns false for indices that are owned locally and for indices not present at all.

◆ ghost_indices()

const IndexSet & Utilities::MPI::Partitioner::ghost_indices ( ) const

Return an IndexSet representation of all ghost indices.

◆ n_ghost_indices()

unsigned int Utilities::MPI::Partitioner::n_ghost_indices ( ) const

Return the number of ghost indices. Same as ghost_indices().n_elements(), but cached for simpler access.

◆ ghost_indices_within_larger_ghost_set()

const std::vector< std::pair< unsigned int, unsigned int > > & Utilities::MPI::Partitioner::ghost_indices_within_larger_ghost_set ( ) const

In case the partitioner was built to define ghost indices as a subset of indices in a larger set of ghosts, this function returns the numbering in terms of ranges within that set. Similar structure as in an IndexSet, but tailored to be iterated over.

In case the partitioner did not take a second set of ghost indices into account, this subset is simply defined as the half-open interval [0, n_ghost_indices()).

◆ ghost_targets()

const std::vector< std::pair< unsigned int, unsigned int > > & Utilities::MPI::Partitioner::ghost_targets ( ) const

Return a list of processors (first entry) and the number of ghost degrees of freedom owned by that processor (second entry). The sum of the latter over all processors equals n_ghost_indices().

◆ import_indices()

const std::vector< std::pair< unsigned int, unsigned int > > & Utilities::MPI::Partitioner::import_indices ( ) const

Return a vector of ranges of local indices that we are importing during compress(), i.e., others' ghosts that belong to the local range. Similar structure as in an IndexSet, but tailored to be iterated over, and some indices may be duplicated. The returned pairs consists of the index of the first element and the index of the element one past the last one in a range.

◆ n_import_indices()

unsigned int Utilities::MPI::Partitioner::n_import_indices ( ) const

Number of import indices, i.e., indices that are ghosts on other processors and we will receive data from.

◆ import_targets()

const std::vector< std::pair< unsigned int, unsigned int > > & Utilities::MPI::Partitioner::import_targets ( ) const

Return a list of processors (first entry) and the number of degrees of freedom imported from it during compress() operation (second entry) for all the processors that data is obtained from, i.e., locally owned indices that are ghosts on other processors.

Note
The returned vector only contains those processor id's for which the second entry is non-zero.

◆ is_compatible()

bool Utilities::MPI::Partitioner::is_compatible ( const Partitioner part) const

Check whether the given partitioner is compatible with the current partitioner. Two partitioners are compatible if they have the same local sizes and the same ghost indices. They do not necessarily need to correspond to the same data that is stored based on these partioner objects. This is a local operation only, i.e., if only some processors decide that the partitioning is not compatible, only these processors will return false, whereas the other processors will return true.

Definition at line 476 of file partitioner.cc.

◆ is_globally_compatible()

bool Utilities::MPI::Partitioner::is_globally_compatible ( const Partitioner part) const

Check whether the given partitioner is compatible with the current partitioner. Two partitioners are compatible if they have the same local size and the same ghost indices. They do not necessarily need to correspond to the same data that is stored based on these partioner objects. As opposed to is_compatible(), this method checks for compatibility among all processors and the method only returns true if the partitioner is the same on all processors. In other words, it does a global "and" operation over the results returned by is_compatible() on all involved processes.

This method performs global communication, so make sure to use it only in a context where all processors call it the same number of times.

Definition at line 503 of file partitioner.cc.

◆ this_mpi_process()

unsigned int Utilities::MPI::Partitioner::this_mpi_process ( ) const

Return the MPI ID of the calling processor. Cached to have simple access.

◆ n_mpi_processes()

unsigned int Utilities::MPI::Partitioner::n_mpi_processes ( ) const

Return the total number of MPI processor participating in the given partitioner. Cached to have simple access.

◆ get_mpi_communicator()

virtual MPI_Comm Utilities::MPI::Partitioner::get_mpi_communicator ( ) const
overridevirtual

Return the underlying MPI communicator.

Implements Utilities::MPI::CommunicationPatternBase.

◆ ghost_indices_initialized()

bool Utilities::MPI::Partitioner::ghost_indices_initialized ( ) const

Return whether ghost indices have been explicitly added as a ghost_indices argument. Only true if a reinit() call or constructor provided that argument.

◆ export_to_ghosted_array_start()

template<typename Number , typename MemorySpaceType = MemorySpace::Host>
void Utilities::MPI::Partitioner::export_to_ghosted_array_start ( const unsigned int  communication_channel,
const ArrayView< const Number, MemorySpaceType > &  locally_owned_array,
const ArrayView< Number, MemorySpaceType > &  temporary_storage,
const ArrayView< Number, MemorySpaceType > &  ghost_array,
std::vector< MPI_Request > &  requests 
) const

Start the exportation of the data in a locally owned array to the range described by the ghost indices of this class.

Parameters
communication_channelSets an offset to the MPI_Isend and MPI_Irecv calls that avoids interference with other ongoing export_to_ghosted_array_start() calls on different entries. Typically handled within the blocks of a block vector. Any value less than 200 is a valid value.
locally_owned_arrayThe array of data from which the data is extracted and sent to the ghost entries on a remote processor.
temporary_storageA temporary storage array of length n_import_indices() that is used to hold the packed data from the locally_owned_array to be sent. Note that this array must not be touched until the respective export_to_ghosted_array_finish() call has been made because the model uses non-blocking communication.
ghost_arrayThe array that will receive the exported data, i.e., the entries that a remote processor sent to the calling process. Its size must either be n_ghost_indices() or equal the number of ghost indices in the larger index set that was given as second argument to set_ghost_indices(). In case only selected indices are sent, no guarantee is made regarding the entries that do not get set. Some of them might be used to organize the transfer and later reset to zero, so make sure you do not use them in computations.
requestsThe list of MPI requests for the ongoing non-blocking communication that will be finalized in the export_to_ghosted_array_finish() call.

This functionality is used in LinearAlgebra::distributed::Vector::update_ghost_values().

◆ export_to_ghosted_array_finish()

template<typename Number , typename MemorySpaceType = MemorySpace::Host>
void Utilities::MPI::Partitioner::export_to_ghosted_array_finish ( const ArrayView< Number, MemorySpaceType > &  ghost_array,
std::vector< MPI_Request > &  requests 
) const

Finish the exportation of the data in a locally owned array to the range described by the ghost indices of this class.

Parameters
ghost_arrayThe array that will receive the exported data started in the export_to_ghosted_array_start(). This must be the same array as passed to that function, otherwise the behavior is undefined.
requestsThe list of MPI requests for the ongoing non-blocking communication that were started in the export_to_ghosted_array_start() call. This must be the same array as passed to that function, otherwise MPI will likely throw an error.

This functionality is used in LinearAlgebra::distributed::Vector::update_ghost_values().

◆ import_from_ghosted_array_start()

template<typename Number , typename MemorySpaceType = MemorySpace::Host>
void Utilities::MPI::Partitioner::import_from_ghosted_array_start ( const VectorOperation::values  vector_operation,
const unsigned int  communication_channel,
const ArrayView< Number, MemorySpaceType > &  ghost_array,
const ArrayView< Number, MemorySpaceType > &  temporary_storage,
std::vector< MPI_Request > &  requests 
) const

Start importing the data on an array indexed by the ghost indices of this class that is later accumulated into a locally owned array with import_from_ghosted_array_finish().

Parameters
vector_operationDefines how the data sent to the owner should be combined with the existing entries, e.g., added into.
communication_channelSets an offset to the MPI_Isend and MPI_Irecv calls that avoids interference with other ongoing import_from_ghosted_array_start() calls on different entries. Typically handled within the blocks of a block vector. Any value less than 200 is a valid value.
ghost_arrayThe array of ghost data that is sent to a remote owner of the respective index in a vector. Its size must either be n_ghost_indices() or equal the number of ghost indices in the larger index set that was given as second argument to set_ghost_indices(). This or the subsequent import_from_ghosted_array_finish() function, the order is implementation-dependent, will set all data entries behind ghost_array to zero.
temporary_storageA temporary storage array of length n_import_indices() that is used to hold the packed data from MPI communication that will later be written into the locally owned array. Note that this array must not be touched until the respective import_from_ghosted_array_finish() call has been made because the model uses non-blocking communication.
requestsThe list of MPI requests for the ongoing non-blocking communication that will be finalized in the export_to_ghosted_array_finish() call.

This functionality is used in LinearAlgebra::distributed::Vector::compress().

◆ import_from_ghosted_array_finish()

template<typename Number , typename MemorySpaceType = MemorySpace::Host>
void Utilities::MPI::Partitioner::import_from_ghosted_array_finish ( const VectorOperation::values  vector_operation,
const ArrayView< const Number, MemorySpaceType > &  temporary_storage,
const ArrayView< Number, MemorySpaceType > &  locally_owned_storage,
const ArrayView< Number, MemorySpaceType > &  ghost_array,
std::vector< MPI_Request > &  requests 
) const

Finish importing the data from an array indexed by the ghost indices of this class into a specified locally owned array, combining the results according to the given input vector_operation.

Parameters
vector_operationDefines how the data sent to the owner should be combined with the existing entries, e.g., added into.
temporary_storageThe same array given to the import_from_ghosted_array_start() call that contains the packed data from MPI communication. In thus function, it is combined at the corresponding entries described by the ghost relations according to vector_operation.
ghost_arrayThe array of ghost data that is sent to a remote owner of the respective index in a vector. Its size must either be n_ghost_indices() or equal the number of ghost indices in the larger index set that was given as second argument to set_ghost_indices(). This function will set all data entries behind ghost_array to zero for the implementation-dependent cases when it was not already done in the import_from_ghosted_array_start() call.
locally_owned_storageThe array of data where the resulting data sent by remote processes to the calling process will be accumulated into.
requestsThe list of MPI requests for the ongoing non-blocking communication that have been initiated in the import_to_ghosted_array_finish() call. This must be the same array as passed to that function, otherwise MPI will likely throw an error.

This functionality is used in LinearAlgebra::distributed::Vector::compress().

◆ memory_consumption()

std::size_t Utilities::MPI::Partitioner::memory_consumption ( ) const

Compute the memory consumption of this structure.

Definition at line 512 of file partitioner.cc.

◆ initialize_import_indices_plain_dev()

void Utilities::MPI::Partitioner::initialize_import_indices_plain_dev ( ) const
private

Initialize import_indices_plain_dev from import_indices_data. This function is only used when using device-aware MPI.

Definition at line 544 of file partitioner.cc.

Member Data Documentation

◆ global_size

types::global_dof_index Utilities::MPI::Partitioner::global_size
private

The global size of the vector over all processors

Definition at line 688 of file partitioner.h.

◆ locally_owned_range_data

IndexSet Utilities::MPI::Partitioner::locally_owned_range_data
private

The range of the vector that is stored locally.

Definition at line 693 of file partitioner.h.

◆ local_range_data

std::pair<types::global_dof_index, types::global_dof_index> Utilities::MPI::Partitioner::local_range_data
private

The range of the vector that is stored locally. Extracted from locally_owned_range for performance reasons.

Definition at line 700 of file partitioner.h.

◆ ghost_indices_data

IndexSet Utilities::MPI::Partitioner::ghost_indices_data
private

The set of indices to which we need to have read access but that are not locally owned

Definition at line 706 of file partitioner.h.

◆ n_ghost_indices_data

unsigned int Utilities::MPI::Partitioner::n_ghost_indices_data
private

A variable caching the number of ghost indices. It would be expensive to use ghost_indices.n_elements() to compute this.

Definition at line 712 of file partitioner.h.

◆ ghost_targets_data

std::vector<std::pair<unsigned int, unsigned int> > Utilities::MPI::Partitioner::ghost_targets_data
private

An array that contains information which processors my ghost indices belong to and how many those indices are

Definition at line 718 of file partitioner.h.

◆ import_indices_data

std::vector<std::pair<unsigned int, unsigned int> > Utilities::MPI::Partitioner::import_indices_data
private

The set of (local) indices that we are importing during compress(), i.e., others' ghosts that belong to the local range. Similar structure as in an IndexSet, but tailored to be iterated over, and some indices may be duplicates.

Definition at line 726 of file partitioner.h.

◆ import_indices_plain_dev

std::vector< Kokkos::View<unsigned int *, MemorySpace::Default::kokkos_space> > Utilities::MPI::Partitioner::import_indices_plain_dev
mutableprivate

The set of (local) indices that we are importing during compress(), i.e., others' ghosts that belong to the local range. The data stored is the same as in import_indices_data but the data is expanded in plain arrays. This variable is only used when using device-aware MPI.

Definition at line 738 of file partitioner.h.

◆ n_import_indices_data

unsigned int Utilities::MPI::Partitioner::n_import_indices_data
private

A variable caching the number of ghost indices. It would be expensive to compute it by iterating over the import indices and accumulate them.

Definition at line 744 of file partitioner.h.

◆ import_targets_data

std::vector<std::pair<unsigned int, unsigned int> > Utilities::MPI::Partitioner::import_targets_data
private

The set of processors and length of data field which send us their ghost data

Definition at line 750 of file partitioner.h.

◆ import_indices_chunks_by_rank_data

std::vector<unsigned int> Utilities::MPI::Partitioner::import_indices_chunks_by_rank_data
private

An array that caches the number of chunks in the import indices per MPI rank. The length is import_indices_data.size()+1.

Definition at line 756 of file partitioner.h.

◆ n_ghost_indices_in_larger_set

unsigned int Utilities::MPI::Partitioner::n_ghost_indices_in_larger_set
private

A variable caching the number of ghost indices in a larger set of indices given by the optional argument to set_ghost_indices().

Definition at line 762 of file partitioner.h.

◆ ghost_indices_subset_chunks_by_rank_data

std::vector<unsigned int> Utilities::MPI::Partitioner::ghost_indices_subset_chunks_by_rank_data
private

An array that caches the number of chunks in the import indices per MPI rank. The length is ghost_indices_subset_data.size()+1.

Definition at line 768 of file partitioner.h.

◆ ghost_indices_subset_data

std::vector<std::pair<unsigned int, unsigned int> > Utilities::MPI::Partitioner::ghost_indices_subset_data
private

The set of indices that appear for an IndexSet that is a subset of a larger set. Similar structure as in an IndexSet within all ghost indices, but tailored to be iterated over.

Definition at line 776 of file partitioner.h.

◆ my_pid

unsigned int Utilities::MPI::Partitioner::my_pid
private

The ID of the current processor in the MPI network

Definition at line 781 of file partitioner.h.

◆ n_procs

unsigned int Utilities::MPI::Partitioner::n_procs
private

The total number of processors active in the problem

Definition at line 786 of file partitioner.h.

◆ communicator

MPI_Comm Utilities::MPI::Partitioner::communicator
private

The MPI communicator involved in the problem

Definition at line 791 of file partitioner.h.

◆ have_ghost_indices

bool Utilities::MPI::Partitioner::have_ghost_indices
private

A variable storing whether the ghost indices have been explicitly set.

Definition at line 796 of file partitioner.h.


The documentation for this class was generated from the following files: