Tpetra parallel linear algebra
Version of the Day
|
Class that encapulates linear problem (Ax = b). More...
#include <Tpetra_LinearProblem_decl.hpp>
Public Types | |
Typedefs | |
using | map_type = Map< LocalOrdinal, GlobalOrdinal, Node > |
using | row_matrix_type = RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > |
using | multivector_type = MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > |
using | vector_type = Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > |
using | linear_problem_type = LinearProblem< Scalar, LocalOrdinal, GlobalOrdinal, Node > |
Typedefs | |
using | packet_type = typename::Kokkos::ArithTraits< Scalar >::val_type |
The type of each datum being sent or received in an Import or Export. More... | |
using | local_ordinal_type = LocalOrdinal |
The type of local indices. More... | |
using | global_ordinal_type = GlobalOrdinal |
The type of global indices. More... | |
using | node_type = Node |
The Node type. If you don't know what this is, don't use it. More... | |
using | device_type = typename Node::device_type |
The Kokkos Device type. More... | |
using | execution_space = typename device_type::execution_space |
The Kokkos execution space. More... | |
Public Member Functions | |
Constructors/Destructor | |
LinearProblem () | |
Default Constructor. More... | |
LinearProblem (const Teuchos::RCP< row_matrix_type > &A, const Teuchos::RCP< multivector_type > &X, const Teuchos::RCP< multivector_type > &B) | |
Constructor with a matrix. More... | |
LinearProblem (const LinearProblem< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Problem) | |
Copy Constructor. More... | |
virtual | ~LinearProblem ()=default |
LinearProblem Destructor. More... | |
Integrity check method | |
void | checkInput () const |
Check input parameters for existence and size consistency. More... | |
Implementation of DistObject interface | |
virtual bool | checkSizes (const SrcDistObject &source) override |
Compare the source and target (this) objects for compatibility. More... | |
Set methods | |
void | setMatrix (Teuchos::RCP< row_matrix_type > A) |
Set Matrix A of linear problem AX = B using a RowMatrix. More... | |
void | setLHS (Teuchos::RCP< multivector_type > X) |
Set left-hand-side X of linear problem AX = B. More... | |
void | setRHS (Teuchos::RCP< multivector_type > B) |
Set right-hand-side B of linear problem AX = B. More... | |
Computational methods | |
void | leftScale (const Teuchos::RCP< const vector_type > &D, Teuchos::ETransp mode=Teuchos::NO_TRANS) |
Perform left scaling of a linear problem. More... | |
void | rightScale (const Teuchos::RCP< const vector_type > &D, Teuchos::ETransp mode=Teuchos::NO_TRANS) |
Perform right scaling of a linear problem. More... | |
Accessor methods | |
Teuchos::RCP< row_matrix_type > | getMatrix () const |
Get an RCP to the matrix A. More... | |
Teuchos::RCP< multivector_type > | getLHS () const |
Get an RCP to the left-hand-side X. More... | |
Teuchos::RCP< multivector_type > | getRHS () const |
Get an RCP to the right-hand-side B. More... | |
Public methods for redistributing data | |
void | doImport (const SrcDistObject &source, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, const CombineMode CM, const bool restrictedMode=false) |
Import data into this object using an Import object ("forward mode"). More... | |
void | doImport (const SrcDistObject &source, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, const CombineMode CM, const bool restrictedMode=false) |
Import data into this object using an Export object ("reverse mode"). More... | |
void | doExport (const SrcDistObject &source, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, const CombineMode CM, const bool restrictedMode=false) |
Export data into this object using an Export object ("forward mode"). More... | |
void | doExport (const SrcDistObject &source, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, const CombineMode CM, const bool restrictedMode=false) |
Export data into this object using an Import object ("reverse mode"). More... | |
void | beginImport (const SrcDistObject &source, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, const CombineMode CM, const bool restrictedMode=false) |
void | beginImport (const SrcDistObject &source, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, const CombineMode CM, const bool restrictedMode=false) |
void | beginExport (const SrcDistObject &source, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, const CombineMode CM, const bool restrictedMode=false) |
void | beginExport (const SrcDistObject &source, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, const CombineMode CM, const bool restrictedMode=false) |
void | endImport (const SrcDistObject &source, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, const CombineMode CM, const bool restrictedMode=false) |
void | endImport (const SrcDistObject &source, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, const CombineMode CM, const bool restrictedMode=false) |
void | endExport (const SrcDistObject &source, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, const CombineMode CM, const bool restrictedMode=false) |
void | endExport (const SrcDistObject &source, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, const CombineMode CM, const bool restrictedMode=false) |
bool | transferArrived () const |
Whether the data from an import/export operation has arrived, and is ready for the unpack and combine step. More... | |
Attribute accessor methods | |
bool | isDistributed () const |
Whether this is a globally distributed object. More... | |
virtual Teuchos::RCP< const map_type > | getMap () const |
The Map describing the parallel distribution of this object. More... | |
I/O methods | |
void | print (std::ostream &os) const |
Print this object to the given output stream. More... | |
Implementation of Teuchos::Describable | |
virtual std::string | description () const |
One-line descriptiion of this object. More... | |
virtual void | describe (Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const |
Print a descriptiion of this object to the given output stream. More... | |
Methods for use only by experts | |
virtual void | removeEmptyProcessesInPlace (const Teuchos::RCP< const map_type > &newMap) |
Remove processes which contain no entries in this object's Map. More... | |
Protected Member Functions | |
virtual size_t | constantNumberOfPackets () const |
Whether the implementation's instance promises always to have a constant number of packets per LID (local index), and if so, how many packets per LID there are. More... | |
virtual void | doTransfer (const SrcDistObject &src, const ::Tpetra::Details::Transfer< local_ordinal_type, global_ordinal_type, node_type > &transfer, const char modeString[], const ReverseOption revOp, const CombineMode CM, const bool restrictedMode) |
Redistribute data across (MPI) processes. More... | |
virtual bool | reallocArraysForNumPacketsPerLid (const size_t numExportLIDs, const size_t numImportLIDs) |
Reallocate numExportPacketsPerLID_ and/or numImportPacketsPerLID_, if necessary. More... | |
void | beginTransfer (const SrcDistObject &src, const ::Tpetra::Details::Transfer< local_ordinal_type, global_ordinal_type, node_type > &transfer, const char modeString[], const ReverseOption revOp, const CombineMode CM, const bool restrictedMode) |
Implementation detail of doTransfer. More... | |
Methods implemented by subclasses and used by doTransfer(). | |
The doTransfer() method uses the subclass' implementations of these methods to implement data transfer. Subclasses of DistObject must implement these methods. This is an instance of the Template Method Pattern. ("Template" here doesn't mean "C++ template"; it means "pattern with holes that are filled in by the subclass' method implementations.") | |
Teuchos::RCP< const map_type > | map_ |
The Map over which this object is distributed. More... | |
Kokkos::DualView< packet_type *, buffer_device_type > | imports_ |
Buffer into which packed data are imported (received from other processes). More... | |
Kokkos::DualView< size_t *, buffer_device_type > | numImportPacketsPerLID_ |
Number of packets to receive for each receive operation. More... | |
Kokkos::DualView< packet_type *, buffer_device_type > | exports_ |
Buffer from which packed data are exported (sent to other processes). More... | |
Kokkos::DualView< size_t *, buffer_device_type > | numExportPacketsPerLID_ |
Number of packets to send for each send operation. More... | |
virtual void | copyAndPermute (const SrcDistObject &source, const size_t numSameIDs, const Kokkos::DualView< const local_ordinal_type *, buffer_device_type > &permuteToLIDs, const Kokkos::DualView< const local_ordinal_type *, buffer_device_type > &permuteFromLIDs, const CombineMode CM) |
Perform copies and permutations that are local to the calling (MPI) process. More... | |
virtual void | copyAndPermute (const SrcDistObject &source, const size_t numSameIDs, const Kokkos::DualView< const local_ordinal_type *, buffer_device_type > &permuteToLIDs, const Kokkos::DualView< const local_ordinal_type *, buffer_device_type > &permuteFromLIDs, const CombineMode CM, const execution_space &space) |
Same as copyAndPermute, but do operations in space . More... | |
virtual void | packAndPrepare (const SrcDistObject &source, const Kokkos::DualView< const local_ordinal_type *, buffer_device_type > &exportLIDs, Kokkos::DualView< packet_type *, buffer_device_type > &exports, Kokkos::DualView< size_t *, buffer_device_type > numPacketsPerLID, size_t &constantNumPackets) |
Pack data and metadata for communication (sends). More... | |
virtual void | packAndPrepare (const SrcDistObject &source, const Kokkos::DualView< const local_ordinal_type *, buffer_device_type > &exportLIDs, Kokkos::DualView< packet_type *, buffer_device_type > &exports, Kokkos::DualView< size_t *, buffer_device_type > numPacketsPerLID, size_t &constantNumPackets, const execution_space &space) |
Same as packAndPrepare, but in an execution space instance. More... | |
virtual void | unpackAndCombine (const Kokkos::DualView< const local_ordinal_type *, buffer_device_type > &importLIDs, Kokkos::DualView< packet_type *, buffer_device_type > imports, Kokkos::DualView< size_t *, buffer_device_type > numPacketsPerLID, const size_t constantNumPackets, const CombineMode combineMode) |
Perform any unpacking and combining after communication. More... | |
virtual void | unpackAndCombine (const Kokkos::DualView< const local_ordinal_type *, buffer_device_type > &importLIDs, Kokkos::DualView< packet_type *, buffer_device_type > imports, Kokkos::DualView< size_t *, buffer_device_type > numPacketsPerLID, const size_t constantNumPackets, const CombineMode combineMode, const execution_space &space) |
std::unique_ptr< std::string > | createPrefix (const char className[], const char methodName[]) const |
virtual bool | reallocImportsIfNeeded (const size_t newSize, const bool verbose, const std::string *prefix, const bool remoteLIDsContiguous=false, const CombineMode CM=INSERT) |
Reallocate imports_ if needed. More... | |
Class that encapulates linear problem (Ax = b).
The LinearProblem class is a wrapper that encapsulates the general information needed for solving a linear system of equations. Currently it accepts a Tpetra matrix/operator, initial guess and RHS and returns the solution.
Scalar | The type of the numerical entries of the matrix. (You can use real-valued or complex-valued types here.) |
LocalOrdinal | The type of local indices. See the documentation of Map for requirements. |
GlobalOrdinal | The type of global indices. See the documentation of Map for requirements. |
Node | The Kokkos Node type. See the documentation of Map for requirements. |
Definition at line 47 of file Tpetra_LinearProblem_decl.hpp.
|
inherited |
The type of each datum being sent or received in an Import or Export.
Note that this type does not always correspond to the Scalar
template parameter of subclasses.
Definition at line 304 of file Tpetra_DistObject_decl.hpp.
|
inherited |
The type of local indices.
Definition at line 306 of file Tpetra_DistObject_decl.hpp.
|
inherited |
The type of global indices.
Definition at line 308 of file Tpetra_DistObject_decl.hpp.
|
inherited |
The Node type. If you don't know what this is, don't use it.
Definition at line 310 of file Tpetra_DistObject_decl.hpp.
|
inherited |
The Kokkos Device type.
Definition at line 313 of file Tpetra_DistObject_decl.hpp.
|
inherited |
The Kokkos execution space.
Definition at line 315 of file Tpetra_DistObject_decl.hpp.
Tpetra::LinearProblem< Scalar, LocalOrdinal, GlobalOrdinal, Node >::LinearProblem | ( | ) |
Default Constructor.
Creates an empty LinearProblem instance. The matrix A, left-hand-side X and right-hand-side B must be set use the setMatrix(), SetLHS() and SetRHS() methods respectively.
Definition at line 31 of file Tpetra_LinearProblem_def.hpp.
Tpetra::LinearProblem< Scalar, LocalOrdinal, GlobalOrdinal, Node >::LinearProblem | ( | const Teuchos::RCP< row_matrix_type > & | A, |
const Teuchos::RCP< multivector_type > & | X, | ||
const Teuchos::RCP< multivector_type > & | B | ||
) |
Constructor with a matrix.
Creates a LinearProblem instance with a matrix.
Definition at line 41 of file Tpetra_LinearProblem_def.hpp.
Tpetra::LinearProblem< Scalar, LocalOrdinal, GlobalOrdinal, Node >::LinearProblem | ( | const LinearProblem< Scalar, LocalOrdinal, GlobalOrdinal, Node > & | Problem | ) |
Copy Constructor.
Definition at line 53 of file Tpetra_LinearProblem_def.hpp.
|
virtualdefault |
LinearProblem Destructor.
void Tpetra::LinearProblem< Scalar, LocalOrdinal, GlobalOrdinal, Node >::checkInput | ( | ) | const |
Check input parameters for existence and size consistency.
Returns 0 if all input parameters are valid. Returns +1 if operator is not a matrix. This is not necessarily an error, but no scaling can be done if the user passes in an operator that is not an matrix.
Definition at line 102 of file Tpetra_LinearProblem_def.hpp.
|
overridevirtual |
Compare the source and target (this) objects for compatibility.
Implements Tpetra::DistObject< Scalar, LocalOrdinal, GlobalOrdinal, Node >.
Definition at line 126 of file Tpetra_LinearProblem_def.hpp.
|
inline |
Set Matrix A of linear problem AX = B using a RowMatrix.
Sets an RCP to a RowMatrix. No copy of the operator is made.
Definition at line 123 of file Tpetra_LinearProblem_decl.hpp.
|
inline |
Set left-hand-side X of linear problem AX = B.
Sets an RCP to a MultiVector. No copy of the object is made.
Definition at line 129 of file Tpetra_LinearProblem_decl.hpp.
|
inline |
Set right-hand-side B of linear problem AX = B.
Sets an RCP to a MultiVector. No copy of the object is made.
Definition at line 134 of file Tpetra_LinearProblem_decl.hpp.
void Tpetra::LinearProblem< Scalar, LocalOrdinal, GlobalOrdinal, Node >::leftScale | ( | const Teuchos::RCP< const vector_type > & | D, |
Teuchos::ETransp | mode = Teuchos::NO_TRANS |
||
) |
Perform left scaling of a linear problem.
Applies the scaling vector D to the left side of the matrix A() and to the right hand side B().
In | D - Vector containing scaling values. D[i] will be applied to the ith row of A() and B(). mode - Indicating if transposed. |
Definition at line 63 of file Tpetra_LinearProblem_def.hpp.
void Tpetra::LinearProblem< Scalar, LocalOrdinal, GlobalOrdinal, Node >::rightScale | ( | const Teuchos::RCP< const vector_type > & | D, |
Teuchos::ETransp | mode = Teuchos::NO_TRANS |
||
) |
Perform right scaling of a linear problem.
Applies the scaling vector D to the right side of the matrix A(). Apply the inverse of D to the initial guess.
In | D - Vector containing scaling values. D[i] will be applied to the ith row of A(). 1/D[i] will be applied to the ith row of B(). mode - Indicating if transposed. |
Definition at line 83 of file Tpetra_LinearProblem_def.hpp.
|
inline |
Get an RCP to the matrix A.
Definition at line 177 of file Tpetra_LinearProblem_decl.hpp.
|
inline |
Get an RCP to the left-hand-side X.
Definition at line 179 of file Tpetra_LinearProblem_decl.hpp.
|
inline |
Get an RCP to the right-hand-side B.
Definition at line 181 of file Tpetra_LinearProblem_decl.hpp.
|
inherited |
Import data into this object using an Import object ("forward mode").
The input DistObject is always the source of the data redistribution operation, and the *this
object is always the target.
If you don't know the difference between forward and reverse mode, then you probably want forward mode. Use this method with your precomputed Import object if you want to do an Import, else use doExport() with a precomputed Export object.
"Restricted Mode" does two things:
*this
, in a "locallyFitted" sense. This cannot be used if (2) is not true, OR there are permutes. The "source" maps still need to match.
source | [in] The "source" object for redistribution. |
importer | [in] Precomputed data redistribution plan. Its source Map must be the same as the input DistObject's Map, and its target Map must be the same as this->getMap() . |
CM | [in] How to combine incoming data with the same global index. |
|
inherited |
Import data into this object using an Export object ("reverse mode").
The input DistObject is always the source of the data redistribution operation, and the *this
object is always the target.
If you don't know the difference between forward and reverse mode, then you probably want forward mode. Use the version of doImport() that takes a precomputed Import object in that case.
"Restricted Mode" does two things:
*this
, in a "locallyFitted" sense. This cannot be used if (2) is not true, OR there are permutes. The "source" maps still need to match.
source | [in] The "source" object for redistribution. |
exporter | [in] Precomputed data redistribution plan. Its target Map must be the same as the input DistObject's Map, and its source Map must be the same as this->getMap() . (Note the difference from forward mode.) |
CM | [in] How to combine incoming data with the same global index. |
|
inherited |
Export data into this object using an Export object ("forward mode").
The input DistObject is always the source of the data redistribution operation, and the *this
object is always the target.
If you don't know the difference between forward and reverse mode, then you probably want forward mode. Use this method with your precomputed Export object if you want to do an Export, else use doImport() with a precomputed Import object.
"Restricted Mode" does two things:
*this
, in a "locallyFitted" sense. This cannot be used if (2) is not true, OR there are permutes. The "source" maps still need to match.
source | [in] The "source" object for redistribution. |
exporter | [in] Precomputed data redistribution plan. Its source Map must be the same as the input DistObject's Map, and its target Map must be the same as this->getMap() . |
CM | [in] How to combine incoming data with the same global index. |
|
inherited |
Export data into this object using an Import object ("reverse mode").
The input DistObject is always the source of the data redistribution operation, and the *this
object is always the target.
If you don't know the difference between forward and reverse mode, then you probably want forward mode. Use the version of doExport() that takes a precomputed Export object in that case.
"Restricted Mode" does two things:
*this
, in a "locallyFitted" sense. This cannot be used if (2) is not true, OR there are permutes. The "source" maps still need to match.
source | [in] The "source" object for redistribution. |
importer | [in] Precomputed data redistribution plan. Its target Map must be the same as the input DistObject's Map, and its source Map must be the same as this->getMap() . (Note the difference from forward mode.) |
CM | [in] How to combine incoming data with the same global index. |
|
inherited |
Whether the data from an import/export operation has arrived, and is ready for the unpack and combine step.
|
inherited |
Whether this is a globally distributed object.
For a definition of "globally distributed" (and its opposite, "locally replicated"), see the documentation of Map's isDistributed() method.
|
inlinevirtualinherited |
The Map describing the parallel distribution of this object.
Note that some Tpetra objects might be distributed using multiple Map objects. For example, CrsMatrix has both a row Map and a column Map. It is up to the subclass to decide which Map to use when invoking the DistObject constructor.
Definition at line 559 of file Tpetra_DistObject_decl.hpp.
|
inherited |
Print this object to the given output stream.
We generally assume that all MPI processes can print to the given stream.
|
virtualinherited |
One-line descriptiion of this object.
We declare this method virtual so that subclasses of DistObject may override it.
Reimplemented in Tpetra::MultiVector< SC, LO, GO, NO >, Tpetra::MultiVector< Scalar, LO, GO, Node >, Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node >, and Tpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node >.
|
virtualinherited |
Print a descriptiion of this object to the given output stream.
We declare this method virtual so that subclasses of Distobject may override it.
Reimplemented in Tpetra::MultiVector< SC, LO, GO, NO >, Tpetra::MultiVector< Scalar, LO, GO, Node >, Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node >, and Tpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node >.
|
virtualinherited |
Remove processes which contain no entries in this object's Map.
On input, this object is distributed over the Map returned by getMap() (the "original Map," with its communicator, the "original communicator"). The input newMap
of this method must be the same as the result of calling getMap()->removeEmptyProcesses()
. On processes in the original communicator which contain zero entries ("excluded processes," as opposed to "included processes"), the input newMap
must be Teuchos::null
(which is what getMap()->removeEmptyProcesses()
returns anyway).
On included processes, reassign this object's Map (that would be returned by getMap()) to the input newMap
, and do any work that needs to be done to restore correct semantics. On excluded processes, free any data that needs freeing, and do any other work that needs to be done to restore correct semantics.
This method has collective semantics over the original communicator. On exit, the only method of this object which is safe to call on excluded processes is the destructor. This implies that subclasses' destructors must not contain communication operations.
|
protectedvirtualinherited |
Whether the implementation's instance promises always to have a constant number of packets per LID (local index), and if so, how many packets per LID there are.
If this method returns zero, the instance says that it might possibly have a different number of packets for each LID (local index) to send or receive. If it returns nonzero, the instance promises that the number of packets is the same for all LIDs, and that the return value is this number of packets per LID.
The default implementation of this method returns zero. This does not affect the behavior of doTransfer() in any way. If a nondefault implementation returns nonzero, doTransfer() will use this information to avoid unnecessary allocation and / or resizing of arrays.
Reimplemented in Tpetra::MultiVector< SC, LO, GO, NO >, Tpetra::MultiVector< Scalar, LO, GO, Node >, and Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node >.
|
protectedvirtualinherited |
Redistribute data across (MPI) processes.
src | [in] The source object, to redistribute into the target object, which is *this object. |
transfer | [in] The Export or Import object representing the communication pattern. (Details::Transfer is the common base class of these two objects.) |
modeString | [in] Human-readable string, for verbose debugging output and error output, explaining what function called this method. Example: "doImport (forward)", "doExport (reverse)". |
revOp | [in] Whether to do a forward or reverse mode redistribution. |
CM | [in] The combine mode that describes how to combine values that map to the same global ID on the same process. |
|
protectedvirtualinherited |
Reallocate numExportPacketsPerLID_ and/or numImportPacketsPerLID_, if necessary.
numExportLIDs | [in] Number of entries in the exportLIDs input array argument of doTransfer(). |
numImportLIDs | [in] Number of entries in the remoteLIDs input array argument of doTransfer(). |
|
protectedinherited |
Implementation detail of doTransfer.
LID DualViews come from the Transfer object given to doTransfer. They are always sync'd on both host and device. Users must never attempt to modify or sync them.
|
protectedvirtualinherited |
Perform copies and permutations that are local to the calling (MPI) process.
Subclasses must reimplement this function. Its default implementation does nothing. Note that the <t>target object of the Export or Import, namely *this
, packs the source object's data.
permuteToLIDs.need_sync_host()
, permuteToLIDs.need_sync_device()
, permuteFromLIDs.need_sync_host()
, and permuteFromLIDs.need_sync_device()
are all false.source | [in] On entry, the source object of the Export or Import operation. |
numSameIDs | [in] The number of elements that are the same on the source and target objects. These elements live on the same process in both the source and target objects. |
permuteToLIDs | [in] List of the elements that are permuted. They are listed by their local index (LID) in the destination object. |
permuteFromLIDs | [in] List of the elements that are permuted. They are listed by their local index (LID) in the source object. |
CM | [in] CombineMode to be used during copyAndPermute; may or may not be used by the particular object being called; behavior with respect to CombineMode may differ by object. |
|
protectedvirtualinherited |
Same as copyAndPermute, but do operations in space
.
|
protectedvirtualinherited |
Pack data and metadata for communication (sends).
Subclasses must reimplement this function. Its default implementation does nothing. Note that the <t>target object of the Export or Import, namely *this
, packs the source object's data.
exportLIDs.need_sync_host ()
and exportLIDs.need_sync_device()
are both false.source | [in] Source object for the redistribution. |
exportLIDs | [in] List of the entries (as local IDs in the source object) that Tpetra will send to other processes. |
exports | [out] On exit, the packed data to send. Implementations must reallocate this as needed (prefer reusing the existing allocation if possible), and may modify and/or sync this wherever they like. |
numPacketsPerLID | [out] On exit, the implementation of this method must do one of two things: either set numPacketsPerLID[i] to the number of packets to be packed for exportLIDs[i] and set constantNumPackets to zero, or set constantNumPackets to a nonzero value. If the latter, the implementation must not modify the entries of numPacketsPerLID . If the former, the implementation may sync numPacketsPerLID this wherever it likes, either to host or to device. The allocation belongs to DistObject, not to subclasses; don't be tempted to change this to pass by reference. |
constantNumPackets | [out] On exit, 0 if the number of packets per LID could differ, else (if nonzero) the number of packets per LID (which must be constant). |
|
protectedvirtualinherited |
Same as packAndPrepare, but in an execution space instance.
|
protectedvirtualinherited |
Perform any unpacking and combining after communication.
Subclasses must reimplement this function. Its default implementation does nothing. Note that the <t>target object of the Export or Import, namely *this
, unpacks the received data into itself, possibly modifying its entries.
importLIDs.need_sync_host ()
and importLIDs.need_sync_device()
are both false.importLIDs | [in] List of the entries (as LIDs in the destination object) we received from other processes. |
imports | [in/out] On input: Buffer of received data to unpack. DistObject promises nothing about where this is sync'd. Implementations may sync this wherever they like, either to host or to device. The allocation belongs to DistObject, not to subclasses; don't be tempted to change this to pass by reference. |
numPacketsPerLID | [in/out] On input: If constantNumPackets is zero, then numPacketsPerLID[i] contains the number of packets imported for importLIDs[i]. DistObject promises nothing about where this is sync'd. Implementations may sync this wherever they like, either to host or to device. The allocation belongs to DistObject, not to subclasses; don't be tempted to change this to pass by reference. |
constantNumPackets | [in] If nonzero, then the number of packets per LID is the same for all entries ("constant") and constantNumPackets is that number. If zero, then numPacketsPerLID[i] is the number of packets to unpack for LID importLIDs[i] . |
combineMode | [in] The CombineMode to use when combining the imported entries with existing entries. |
|
protectedvirtualinherited |
Reallocate imports_ if needed.
This unfortunately must be declared protected, for the same reason that imports_ is declared protected.
newSize | [in] New size of imports_. |
verbose | [in] Whether to print verbose debugging output to stderr on every (MPI) process in the communicator. |
prefix | [in] If verbose is true , then this is a nonnull prefix to print at the beginning of each line of verbose debugging output. Otherwise, not used. |
We don't need a "reallocExportsIfNeeded" method, because exports_
always gets passed into packAndPrepare() by nonconst reference. Thus, that method can resize the DualView without needing to call other DistObject methods.
Reimplemented in Tpetra::MultiVector< SC, LO, GO, NO >, Tpetra::MultiVector< Scalar, LO, GO, Node >, and Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node >.
|
protectedinherited |
The Map over which this object is distributed.
Definition at line 968 of file Tpetra_DistObject_decl.hpp.
|
protectedinherited |
Buffer into which packed data are imported (received from other processes).
Unfortunately, I had to declare these protected, because CrsMatrix uses them at one point. Please, nobody else use them.
Definition at line 981 of file Tpetra_DistObject_decl.hpp.
|
protectedinherited |
Number of packets to receive for each receive operation.
This array is used in Distributor::doPosts() (and doReversePosts()) when starting the ireceive operation.
This may be ignored in doTransfer() if constantNumPackets is nonzero, indicating a constant number of packets per LID. (For example, MultiVector sets the constantNumPackets output argument of packAndPrepare() to the number of columns in the multivector.)
Unfortunately, I had to declare this protected, because CrsMatrix uses it at one point. Please, nobody else use it.
Definition at line 1021 of file Tpetra_DistObject_decl.hpp.
|
protectedinherited |
Buffer from which packed data are exported (sent to other processes).
Unfortunately, I had to declare this protected, because CrsMatrix uses it at one point. Please, nobody else use it.
Definition at line 1028 of file Tpetra_DistObject_decl.hpp.
|
protectedinherited |
Number of packets to send for each send operation.
This array is used in Distributor::doPosts() (and doReversePosts()) for preparing for the send operation.
This may be ignored in doTransfer() if constantNumPackets is nonzero, indicating a constant number of packets per LID. (For example, MultiVector sets the constantNumPackets output argument of packAndPrepare() to the number of columns in the multivector.)
Unfortunately, I had to declare this protected, because CrsMatrix uses them at one point. Please, nobody else use it.
Definition at line 1043 of file Tpetra_DistObject_decl.hpp.