Tpetra Matrix/Vector Services Version of the Day
Public Types | Static Public Member Functions
Tpetra::MatrixMarket::Reader< SparseMatrixType > Class Template Reference

Matrix Market file reader for CrsMatrix and MultiVector. More...

#include <MatrixMarket_Tpetra.hpp>

List of all members.

Public Types

typedef SparseMatrixType sparse_matrix_type
 This class' template parameter; a specialization of CrsMatrix.
typedef
SparseMatrixType::scalar_type 
scalar_type
 Type of the entries of the sparse matrix.
typedef
SparseMatrixType::local_ordinal_type 
local_ordinal_type
 Only used to define map_type.
typedef
SparseMatrixType::global_ordinal_type 
global_ordinal_type
 Type of indices as read from the Matrix Market file.
typedef SparseMatrixType::node_type node_type
 The Kokkos Node type.
typedef MultiVector
< scalar_type,
local_ordinal_type,
global_ordinal_type, node_type
multivector_type
 The MultiVector type associated with SparseMatrixType.

Static Public Member Functions

static sparse_matrix_ptr readSparseFile (const std::string &filename, const RCP< const Comm< int > > &pComm, const RCP< node_type > &pNode, const bool callFillComplete=true, const bool tolerant=false, const bool debug=false)
 Read sparse matrix from the given Matrix Market file.
static sparse_matrix_ptr readSparse (std::istream &in, const RCP< const Comm< int > > &pComm, const RCP< node_type > &pNode, const bool callFillComplete=true, const bool tolerant=false, const bool debug=false)
 Read sparse matrix from the given Matrix Market input stream.
static RCP< multivector_typereadDenseFile (const std::string &filename, const RCP< const comm_type > &pComm, const RCP< node_type > &pNode, RCP< const map_type > &pMap, const bool tolerant=false, const bool debug=false)
 Read dense matrix (as a MultiVector) from the given Matrix Market file.
static RCP< multivector_typereadDense (std::istream &in, const RCP< const comm_type > &pComm, const RCP< node_type > &pNode, RCP< const map_type > &pMap, const bool tolerant=false, const bool debug=false)
 Read dense matrix (as a MultiVector) from the given Matrix Market input stream.

Detailed Description

template<class SparseMatrixType>
class Tpetra::MatrixMarket::Reader< SparseMatrixType >

Matrix Market file reader for CrsMatrix and MultiVector.

Author:
Mark Hoemmen

The Matrix Market (see their web site for details) defines a human-readable ASCII text file format for interchange of sparse and dense matrices. This class defines methods for reading sparse and dense matrices from a Matrix Market file or input stream.

All methods of this class assume that the file is only openable resp. the input stream is only readable, on the MPI process with Rank 0 (with respect to the MPI communicator over which the given CrsMatrix or MultiVector is to be distributed).

We define the MultiVector type returned by readDense() and readDenseFile() using the scalar_type, local_ordinal_type, global_ordinal_type, and node_type typedefs in SparseMatrixType. This ensures that the multivectors returned by those methods have a type compatible with the CrsMatrix sparse matrices returned by readSparse() and readSparseFile(). We do this because the typical use case of Matrix Market files in Trilinos is to test sparse matrix methods, which usually involves reading a sparse matrix A and perhaps also a dense right-hand side b. Also, this lets you use CrsMatrix objects with non-default LocalMatOps template parameters. (If we templated on Scalar, LocalOrdinal, GlobalOrdinal, and Node, we would also have to template on LocalMatOps in order to deal with CrsMatrix types with nondefault LocalMatOps. That would tie Reader to CrsMatrix anyway, since MultiVector is not templated on LocalMatOps. As a result, we might as well just template on the CrsMatrix type, in order to use arbitrary LocalMatOps types without additional code.)

Template Parameters:
SparseMatrixTypeA specialization of Tpetra::CrsMatrix.

Templating on the specialization of CrsMatrix means that the Reader expects matrix data of a type compatible with the CrsMatrix's scalar_type. In general, Matrix Market files may contain data of integer, real, or complex type. However, the reader methods have to return a CrsMatrix of a specific type, so we require that you declare a Reader with the CrsMatrix type that you want and that you expect the file(s) to contain.

We didn't find any of the alternatives to this approach acceptable. One possibility would have been to have the reader methods return a "container that can hold anything," like a boost::any. However, then you would have to know all five template arguments of the CrsMatrix in order to get the actual CrsMatrix object out. C++ doesn't have algebraic data types (see the Wikipedia entry for a good definition) that are disjoint unions of different types. Thus, we couldn't have had the readers return a CrsMatrix with scalar_type = "int or double or complex<double>." While you can implement such a type in C++ (see e.g., boost::variant), it would not be interchangeable for its component types. This is because it may not have the same memory layout (e.g., copying an array of boost::variant<int, double, complex<double> > bitwise into an array of int may not work).

Definition at line 165 of file MatrixMarket_Tpetra.hpp.


Member Typedef Documentation

template<class SparseMatrixType >
Tpetra::MatrixMarket::Reader< SparseMatrixType >::sparse_matrix_type

This class' template parameter; a specialization of CrsMatrix.

Definition at line 169 of file MatrixMarket_Tpetra.hpp.

template<class SparseMatrixType >
Tpetra::MatrixMarket::Reader< SparseMatrixType >::scalar_type

Type of the entries of the sparse matrix.

Definition at line 174 of file MatrixMarket_Tpetra.hpp.

template<class SparseMatrixType >
Tpetra::MatrixMarket::Reader< SparseMatrixType >::local_ordinal_type

Only used to define map_type.

Definition at line 177 of file MatrixMarket_Tpetra.hpp.

template<class SparseMatrixType >
Tpetra::MatrixMarket::Reader< SparseMatrixType >::global_ordinal_type

Type of indices as read from the Matrix Market file.

Indices of the sparse matrix are read in as global ordinals, since Matrix Market files represent the whole matrix and don't have a notion of distribution.

Definition at line 185 of file MatrixMarket_Tpetra.hpp.

template<class SparseMatrixType >
Tpetra::MatrixMarket::Reader< SparseMatrixType >::node_type

The Kokkos Node type.

Definition at line 189 of file MatrixMarket_Tpetra.hpp.

template<class SparseMatrixType >
Tpetra::MatrixMarket::Reader< SparseMatrixType >::multivector_type

The MultiVector type associated with SparseMatrixType.

Definition at line 196 of file MatrixMarket_Tpetra.hpp.


Member Function Documentation

template<class SparseMatrixType >
static sparse_matrix_ptr Tpetra::MatrixMarket::Reader< SparseMatrixType >::readSparseFile ( const std::string &  filename,
const RCP< const Comm< int > > &  pComm,
const RCP< node_type > &  pNode,
const bool  callFillComplete = true,
const bool  tolerant = false,
const bool  debug = false 
) [inline, static]

Read sparse matrix from the given Matrix Market file.

Open the given file on MPI Rank 0 (with respect to the given communicator). The file should contain Matrix Market "coordinate" format sparse matrix data. Read that data on Rank 0, and distribute it to all processors. Return the resulting distributed CrsMatrix.

Note:
This is a collective operation. Only Rank 0 opens the file and reads data from it, but all ranks participate and wait for the final result.
Parameters:
filename[in] Name of the Matrix Market file.
pComm[in] Communicator containing all processor(s) over which the sparse matrix will be distributed.
pNode[in] Kokkos Node object.
callFillComplete[in] Whether to call fillComplete() on the Tpetra::CrsMatrix, after adding all the entries read in from the input stream.
tolerant[in] Whether to read the data tolerantly from the file.
debug[in] Whether to produce copious status output useful for Tpetra developers, but probably not useful for anyone else.

Definition at line 1185 of file MatrixMarket_Tpetra.hpp.

template<class SparseMatrixType >
static sparse_matrix_ptr Tpetra::MatrixMarket::Reader< SparseMatrixType >::readSparse ( std::istream &  in,
const RCP< const Comm< int > > &  pComm,
const RCP< node_type > &  pNode,
const bool  callFillComplete = true,
const bool  tolerant = false,
const bool  debug = false 
) [inline, static]

Read sparse matrix from the given Matrix Market input stream.

The given input stream need only be readable by MPI Rank 0 (with respect to the given communicator). The input stream should contain Matrix Market "coordinate" format sparse matrix data. Read that data on Rank 0, and distribute it to all processors. Return the resulting distributed CrsMatrix.

Note:
This is a collective operation. Only Rank 0 reads data from the input stream, but all ranks participate and wait for the final result.
Parameters:
in[in] The input stream from which to read.
pComm[in] Communicator containing all processor(s) over which the sparse matrix will be distributed.
pNode[in] Kokkos Node object.
callFillComplete[in] Whether to call fillComplete() on the Tpetra::CrsMatrix, after adding all the entries read in from the input stream. (Not calling fillComplete() may be useful if you want to change the matrix after reading it from a file.)
tolerant[in] Whether to read the data tolerantly from the file.
debug[in] Whether to produce copious status output useful for Tpetra developers, but probably not useful for anyone else.

Definition at line 1231 of file MatrixMarket_Tpetra.hpp.

template<class SparseMatrixType >
static RCP<multivector_type> Tpetra::MatrixMarket::Reader< SparseMatrixType >::readDenseFile ( const std::string &  filename,
const RCP< const comm_type > &  pComm,
const RCP< node_type > &  pNode,
RCP< const map_type > &  pMap,
const bool  tolerant = false,
const bool  debug = false 
) [inline, static]

Read dense matrix (as a MultiVector) from the given Matrix Market file.

Open the given file on MPI Rank 0 (with respect to the given communicator). The file should contain Matrix Market "array" format dense matrix data. Read that data on Rank 0, and distribute it to all processors. Return the resulting distributed MultiVector.

See documentation of readDense() for details.

Parameters:
filename[in] Name of the Matrix Market file from which to read. Both the filename and the file itself are only accessed on Rank 0 of the given communicator.
pComm[in] Communicator containing all process(es) over which the dense matrix will be distributed.
pNode[in] Kokkos Node object.
pMap[in/out] On input: if nonnull, the map describing how to distribute the vector (not modified). In this case, the map's communicator and node must equal pComm resp. pNode. If null on input, then on output, a sensible (contiguous and uniformly distributed over the given communicator) map describing the distribution of the output multivector.
tolerant[in] Whether to read the data tolerantly from the file.
debug[in] Whether to produce copious status output useful for Tpetra developers, but probably not useful for anyone else.

Definition at line 1728 of file MatrixMarket_Tpetra.hpp.

template<class SparseMatrixType >
static RCP<multivector_type> Tpetra::MatrixMarket::Reader< SparseMatrixType >::readDense ( std::istream &  in,
const RCP< const comm_type > &  pComm,
const RCP< node_type > &  pNode,
RCP< const map_type > &  pMap,
const bool  tolerant = false,
const bool  debug = false 
) [inline, static]

Read dense matrix (as a MultiVector) from the given Matrix Market input stream.

The given input stream need only be readable by MPI Rank 0 (with respect to the given communicator). The input stream should contain Matrix Market "array" format dense matrix data. Read that data on Rank 0, and distribute it to all processors. Return the resulting distributed MultiVector.

Unlike readSparse(), this method allows callers to supply a Map over which to distribute the resulting MultiVector. The Map argument is optional; if null, we construct our own reasonable Map. We let users supply their own Map, because a common case in Tpetra is to read in or construct a sparse matrix first, and then create dense (multi)vectors distributed with the sparse matrix's domain or range Map.

Note:
This is a collective operation. Only Rank 0 opens the file and reads data from it, but all ranks participate and wait for the final result.
"Tolerant" parsing mode means something different for dense matrices than it does for sparse matrices. Since Matrix Market dense matrix files don't store indices with each value, unlike sparse matrices, we can't determine the matrix dimensions from the data alone. Thus, we require the metadata to include a valid number of rows and columns. "Tolerant" in the dense case refers to the data; in tolerant mode, the number of stored matrix entries may be more or less than the number reported by the metadata (number of rows times number of columns). If more, the extra data are ignored; if less, the remainder is filled in with zeros.
On Map compatibility: Suppose that you write a multivector X to a file. Then, you read it back in as a different multivector Y distributed over the same communicator, but with a Map constructed by the input routine (i.e., a null Map on input to readDenseFile() or readDense()). In that case, the only properties shared by the maps of X and Y are that they have the same communicator and the same number of GIDs. The Maps need not necessarily be compatible (in the sense of isCompatible()), and they certainly need not necessarily be the same Map (in the sense of isSameAs()).
Parameters:
in[in] The input stream from which to read. The stream is only accessed on Rank 0 of the given communicator.
pComm[in] Communicator containing all process(es) over which the dense matrix will be distributed.
pNode[in] Kokkos Node object.
pMap[in/out] On input: if nonnull, the map describing how to distribute the vector (not modified). In this case, the map's communicator and node must equal pComm resp. pNode. If null on input, then on output, a sensible (contiguous and uniformly distributed over the given communicator) map describing the distribution of the output multivector.
tolerant[in] Whether to read the data tolerantly from the stream.
debug[in] Whether to produce copious status output useful for Tpetra developers, but probably not useful for anyone else.

Definition at line 1812 of file MatrixMarket_Tpetra.hpp.


The documentation for this class was generated from the following file:
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines