Tpetra parallel linear algebra Version of the Day
|
Declaration and definition of Tpetra::Details::readAndDealOutTriples, which reads a Matrix Market file or input stream on one process, and distributes the resulting sparse matrix entries to the other processes. More...
#include "TpetraCore_config.h"
#include "Tpetra_Details_PackTriples.hpp"
#include "Kokkos_ArithTraits.hpp"
#include "Teuchos_MatrixMarket_generic.hpp"
#include "Teuchos_CommHelpers.hpp"
#include <iostream>
#include <typeinfo>
Go to the source code of this file.
Classes | |
struct | Tpetra::Details::Impl::ReadLine< SC, GO, isComplex > |
Implementation of the readLine stand-alone function in this namespace (see below). More... | |
struct | Tpetra::Details::Impl::ReadLine< SC, GO, true > |
Complex-arithmetic partial specialization of ReadLine. More... | |
struct | Tpetra::Details::Impl::ReadLine< SC, GO, false > |
Real-arithmetic partial specialization of ReadLine. More... | |
Namespaces | |
namespace | Tpetra |
Namespace Tpetra contains the class and methods constituting the Tpetra library. | |
namespace | Tpetra::Details |
Nonmember function that computes a residual Computes R = B - A * X. | |
Functions | |
template<class OrdinalType , class RealType > | |
bool | Tpetra::Details::Impl::readComplexData (std::istream &istr, OrdinalType &rowIndex, OrdinalType &colIndex, RealType &realPart, RealType &imagPart, const std::size_t lineNumber, const bool tolerant) |
Read "<rowIndex> <colIndex> <realPart> <imagPart>" from a line. | |
template<class SC , class GO > | |
int | Tpetra::Details::Impl::readLine (std::function< int(const GO, const GO, const SC &)> processTriple, const std::string &line, const std::size_t lineNumber, const bool tolerant=false, std::ostream *errStrm=NULL, const bool debug=false) |
Take a line from the Matrix Market file or input stream, and process the sparse matrix entry in that line. | |
template<class SC , class GO > | |
int | Tpetra::Details::Impl::readTriples (std::istream &inputStream, std::size_t &curLineNum, std::size_t &numTriplesRead, std::function< int(const GO, const GO, const SC &)> processTriple, const std::size_t maxNumTriplesToRead, const bool tolerant=false, std::ostream *errStrm=NULL, const bool debug=false) |
Read at most numTriplesToRead triples from the given Matrix Market input stream, and pass along any resulting matrix entries to the given closure. | |
template<class SC , class GO > | |
int | Tpetra::Details::Impl::readAndSendOneBatchOfTriples (std::istream &inputStream, std::size_t &curLineNum, std::size_t &numEntRead, ::Teuchos::ArrayRCP< int > &sizeBuf, ::Teuchos::ArrayRCP< char > &msgBuf, std::vector< GO > &rowInds, std::vector< GO > &colInds, std::vector< SC > &vals, const std::size_t maxNumEntPerMsg, const int destRank, const ::Teuchos::Comm< int > &comm, const bool tolerant=false, std::ostream *errStrm=NULL, const bool debug=false) |
Read at most maxNumEntPerMsg sparse matrix entries from the input stream, and send them to the process with rank destRank. | |
template<class SC , class GO , class CommRequestPtr > | |
int | Tpetra::Details::Impl::recvOneBatchOfTriples (std::vector< GO > &rowInds, std::vector< GO > &colInds, std::vector< SC > &vals, int &numEnt, ::Teuchos::ArrayRCP< int > &sizeBuf, ::Teuchos::ArrayRCP< char > &msgBuf, CommRequestPtr &sizeReq, const int srcRank, const ::Teuchos::Comm< int > &comm, const bool tolerant=false, std::ostream *errStrm=NULL, const bool debug=false) |
Read at most maxNumEntPerMsg sparse matrix entries from the input stream, and send them to the process with rank destRank. | |
template<class SC , class GO > | |
int | Tpetra::Details::readAndDealOutTriples (std::istream &inputStream, std::size_t &curLineNum, std::size_t &totalNumEntRead, std::function< int(const GO, const GO, const SC &)> processTriple, const std::size_t maxNumEntPerMsg, const ::Teuchos::Comm< int > &comm, const bool tolerant=false, std::ostream *errStrm=NULL, const bool debug=false) |
On Process 0 in the given communicator, read sparse matrix entries (in chunks of at most maxNumEntPerMsg entries at a time) from the input stream, and "deal them out" to all other processes in the communicator. | |
Declaration and definition of Tpetra::Details::readAndDealOutTriples, which reads a Matrix Market file or input stream on one process, and distributes the resulting sparse matrix entries to the other processes.
Definition in file Tpetra_Details_ReadTriples.hpp.
bool Tpetra::Details::Impl::readComplexData | ( | std::istream & | istr, |
OrdinalType & | rowIndex, | ||
OrdinalType & | colIndex, | ||
RealType & | realPart, | ||
RealType & | imagPart, | ||
const std::size_t | lineNumber, | ||
const bool | tolerant | ||
) |
Read "<rowIndex> <colIndex> <realPart> <imagPart>" from a line.
Matrix Market files that store a sparse matrix with complex values do so with one sparse matrix entry per line. It is stored as space-delimited ASCII text: the row index, the column index, the real part, and the imaginary part, in that order. Both the row and column indices are 1-based. This function attempts to read one line from the given input stream istr, extract the row and column indices and the real and imaginary parts, and write them to the corresponding output variables.
istr | [in/out] Input stream from which to attempt to read one line. |
rowIndex | [out] On output: if successful, the row index read from the line. |
colIndex | [out] On output: if successful, the column index read from the line. |
realPart | [out] On output: if successful, the real part of the matrix entry's value read from the line. |
imagPart | [out] On output: if successful, the imaginary part of the matrix entry's value read from the line. |
lineNumber | [in] The current line number. Used only for diagnostic error messages. |
tolerant | [in] Whether to parse tolerantly. In tolerant mode, if this function fails in any way to read any of the data, it will return false without throwing an exception. Otherwise, this function will either throw an exception or return true. |
Definition at line 112 of file Tpetra_Details_ReadTriples.hpp.
int Tpetra::Details::Impl::readLine | ( | std::function< int(const GO, const GO, const SC &)> | processTriple, |
const std::string & | line, | ||
const std::size_t | lineNumber, | ||
const bool | tolerant = false , |
||
std::ostream * | errStrm = NULL , |
||
const bool | debug = false |
||
) |
Take a line from the Matrix Market file or input stream, and process the sparse matrix entry in that line.
The line must be a valid Matrix Market line, not a comment.
SC | The type of the value of each matrix entry. |
GO | The type of each (global) index of each matrix entry. |
processTriple | [in] Closure, generally with side effects, that takes in and stores off a sparse matrix entry. First argument is the (global) row index, second argument is the (global) column index, and third argument is the value of the entry. The closure must NOT do MPI communication. Return value is an error code, that is zero if and only if the closure succeeded. |
line | [in] The line from the Matrix Market file or input stream to read. |
lineNumber | [in] Current line number in the file or input stream. |
tolerant | [in] Whether to read tolerantly. |
errStrm | [in] If not NULL, print any error messages to this stream. |
debug | [in] If true, print debug messages to *errStrm . |
Definition at line 401 of file Tpetra_Details_ReadTriples.hpp.
int Tpetra::Details::Impl::readTriples | ( | std::istream & | inputStream, |
std::size_t & | curLineNum, | ||
std::size_t & | numTriplesRead, | ||
std::function< int(const GO, const GO, const SC &)> | processTriple, | ||
const std::size_t | maxNumTriplesToRead, | ||
const bool | tolerant = false , |
||
std::ostream * | errStrm = NULL , |
||
const bool | debug = false |
||
) |
Read at most numTriplesToRead triples from the given Matrix Market input stream, and pass along any resulting matrix entries to the given closure.
The line must be a valid Matrix Market line, not a comment.
SC | The type of the value of each matrix entry. |
GO | The type of each (global) index of each matrix entry. |
inputStream | [in/out] Input stream from which to read. |
curLineNum | [in/out] Current line number in the input stream. |
numTriplesRead | [out] On output: Number of matrix triples (row index, column index, value) successfully read from the input stream on this call to the function.x |
processTriple | [in] Closure, generally with side effects, that takes in and stores off a sparse matrix entry. First argument is the (global) row index, second argument is the (global) column index, and third argument is the value of the entry. The closure must NOT do MPI communication. Return value is an error code, that is zero if and only if the closure succeeded. |
maxNumTriplesToRead | [in] Maximum number of triples to read from the input stream on this call of the function. This is a strict upper bound for numTriplesRead (see above). |
tolerant | [in] Whether to read tolerantly. |
errStrm | [in] If not NULL, print any error messages to this stream. |
debug | [in] If true, print debug messages to *errStrm . |
Definition at line 444 of file Tpetra_Details_ReadTriples.hpp.
int Tpetra::Details::Impl::readAndSendOneBatchOfTriples | ( | std::istream & | inputStream, |
std::size_t & | curLineNum, | ||
std::size_t & | numEntRead, | ||
::Teuchos::ArrayRCP< int > & | sizeBuf, | ||
::Teuchos::ArrayRCP< char > & | msgBuf, | ||
std::vector< GO > & | rowInds, | ||
std::vector< GO > & | colInds, | ||
std::vector< SC > & | vals, | ||
const std::size_t | maxNumEntPerMsg, | ||
const int | destRank, | ||
const ::Teuchos::Comm< int > & | comm, | ||
const bool | tolerant = false , |
||
std::ostream * | errStrm = NULL , |
||
const bool | debug = false |
||
) |
Read at most maxNumEntPerMsg sparse matrix entries from the input stream, and send them to the process with rank destRank.
To be called only by the sending process.
SC | The type of the value of each matrix entry. |
GO | The type of each (global) index of each matrix entry. |
inputStream | [in/out] Input stream from which to read. |
curLineNum | [in/out] Current line number in the input stream. |
numEntRead | [out] Number of matrix entries successfully read from the input stream on this call of the function. |
sizeBuf | [in/out] Array of length 1, for sending message size. |
msgBuf | [in/out] Message buffer; to be resized as needed. |
rowInds | [out] Row indices read from the file. |
colInds | [out] Column indices read from the file. |
vals | [out] Matrix values read from the file. |
maxNumEntPerMsg | [in] Maximum number of sparse matrix entries to read from the input stream on this call to the function. |
destRank | [in] Rank of the process where to send the triples. |
comm | [in] Communicator to use for sending the triples. |
tolerant | [in] Whether to read tolerantly. |
errStrm | [in] If not NULL, print any error messages to this stream. |
Definition at line 562 of file Tpetra_Details_ReadTriples.hpp.
int Tpetra::Details::Impl::recvOneBatchOfTriples | ( | std::vector< GO > & | rowInds, |
std::vector< GO > & | colInds, | ||
std::vector< SC > & | vals, | ||
int & | numEnt, | ||
::Teuchos::ArrayRCP< int > & | sizeBuf, | ||
::Teuchos::ArrayRCP< char > & | msgBuf, | ||
CommRequestPtr & | sizeReq, | ||
const int | srcRank, | ||
const ::Teuchos::Comm< int > & | comm, | ||
const bool | tolerant = false , |
||
std::ostream * | errStrm = NULL , |
||
const bool | debug = false |
||
) |
Read at most maxNumEntPerMsg sparse matrix entries from the input stream, and send them to the process with rank destRank.
To be called only by the sending process.
SC | The type of the value of each matrix entry. |
GO | The type of each (global) index of each matrix entry. |
CommRequestPtr | The type of a (smart) pointer to a "communication request" returned by, e.g., ::Teuchos::ireceive. It must implement operator= and operator->, and the thing to which it points must implement void wait() . A model for this is ::Teuchos::RCP< ::Teuchos::CommRequest<int> >. |
rowInds | [out] Row indices to receive. Will be resized as needed. |
colInds | [out] Column indices to receive. Will be resized as needed. |
vals | [out] Matrix values to receive. Will be resized as needed. |
numEnt | [out] Number of matrix entries (triples) received. |
sizeBuf | [in/out] Array of length 1, for receiving message size (size in bytes of msgBuf ). |
msgBuf | [in/out] Message buffer; to be resized as needed. |
sizeReq | [in/out] Preposed receive request for message size. After waiting on this, you may read the contents of sizeBuf. This is a nonconst (smart) pointer reference, so that we can assign to it. A model for this is ::Teuchos::RCP< ::Teuchos::CommRequest<int> >&. |
srcRank | [in] Rank of the process from which to receive the matrix entries (triples). |
comm | [in] Communicator to use for receiving the triples. |
tolerant | [in] Whether to read tolerantly. |
errStrm | [in] If not NULL, print any error messages to this stream. |
Definition at line 812 of file Tpetra_Details_ReadTriples.hpp.