Kokkos Node API and Local Linear Algebra Kernels Version of the Day
Public Member Functions
TestSparseOps< SparseOpsType > Class Template Reference

Helper functions for tests of Kokkos sparse kernels. More...

#include <DefaultSparseOps_TestSparseOps.hpp>

List of all members.

Public Member Functions

void makeDenseTestProblem (Teuchos::RCP< dense_matrix_type > &A_out, Teuchos::RCP< dense_matrix_type > &L_out, Teuchos::RCP< dense_matrix_type > &U_out, Teuchos::Array< ordinal_type > &pivots, const ordinal_type numRows) const
 Compute dense test matrices for sparse triangular solve.
void readFile (size_t &numRows, size_t &numCols, Teuchos::ArrayRCP< const size_t > &rowptr, Teuchos::ArrayRCP< const ordinal_type > &colind, Teuchos::ArrayRCP< const scalar_type > &values, const std::string &filename) const
 Read a CSR-format sparse matrix from a Matrix Market file.
Teuchos::RCP< SparseOpsType > makeSparseOpsFromFile (ordinal_type &numRows, ordinal_type &numCols, const Teuchos::RCP< node_type > &node, const std::string &filename, const Teuchos::EUplo uplo=Teuchos::UNDEF_TRI, const Teuchos::EDiag diag=Teuchos::NON_UNIT_DIAG) const
 Initialize and return a sparse kernels object from a Matrix Market file.
Teuchos::RCP< SparseOpsType > makeSparseOps (const Teuchos::RCP< node_type > &node, const Teuchos::ArrayRCP< const size_t > &ptr, const Teuchos::ArrayRCP< const ordinal_type > &ind, const Teuchos::ArrayRCP< const scalar_type > &val, const Teuchos::EUplo uplo=Teuchos::UNDEF_TRI, const Teuchos::EDiag diag=Teuchos::NON_UNIT_DIAG) const
 Initialize and return a sparse kernels object.
Teuchos::RCP< SparseOpsType > denseTriToSparseOps (const dense_matrix_type &A, const Teuchos::RCP< node_type > &node, const Teuchos::EUplo uplo, const Teuchos::EDiag diag) const
 Convert the given dense triangular matrix to a sparse kernels object.
Teuchos::RCP< SparseOpsType > denseToSparseOps (const dense_matrix_type &A, const Teuchos::RCP< node_type > &node) const
 Convert the given dense matrix to a sparse kernels object.
Teuchos::RCP
< Kokkos::MultiVector
< scalar_type, node_type > > 
makeMultiVector (const Teuchos::RCP< node_type > &node, const ordinal_type numRows, const ordinal_type numCols) const
magnitude_type maxRelativeError (const Teuchos::RCP< const Kokkos::MultiVector< scalar_type, node_type > > &X, const Teuchos::RCP< const Kokkos::MultiVector< scalar_type, node_type > > &Y, const Teuchos::RCP< Kokkos::MultiVector< scalar_type, node_type > > &Z) const
 Return the maximum relative error between the two multivectors.
void testSparseOps (const Teuchos::RCP< node_type > &node, const ordinal_type N, const ordinal_type numVecs, const magnitude_type tol) const
 Test sparse matrix-(multi)vector multiply and sparse triangular solve.
void benchmarkSparseMatVecFromFile (std::vector< std::pair< std::string, double > > &results, const std::string &filename, const std::string &label, const Teuchos::RCP< node_type > &node, const ordinal_type numVecs, const int numTrials) const
 Benchmark sparse mat-vec with a matrix read from a Matrix Market file.

Dense matrix, BLAS, and LAPACK typedefs

Teuchos' BLAS and LAPACK wrappers do accept an "OrdinalType" template parameter, but they only work for OrdinalType=int (unless you've built your BLAS and LAPACK with 64-bit integer indices). Thus, in turn, we use int for SerialDenseMatrix's ordinal_type for compatibility.

The BLAS and LAPACK typedefs are private because they are implementation details.

typedef
Teuchos::SerialDenseMatrix
< int, scalar_type > 
dense_matrix_type
 The type of local dense matrices in column-major order.

Detailed Description

template<class SparseOpsType>
class TestSparseOps< SparseOpsType >

Helper functions for tests of Kokkos sparse kernels.

Template Parameters:
SparseOpsTypeAny class which implements the Kokkos sparse kernels interface, for which DefaultHostSparseOps is a working example and EmptySparseKernel (in kokkos/LinAlg/examples/KokkosExamples_EmptySparseKernelClass.hpp) is a stub.

Definition at line 69 of file DefaultSparseOps_TestSparseOps.hpp.


Member Typedef Documentation

template<class SparseOpsType>
typedef Teuchos::SerialDenseMatrix<int, scalar_type> TestSparseOps< SparseOpsType >::dense_matrix_type

The type of local dense matrices in column-major order.

Definition at line 93 of file DefaultSparseOps_TestSparseOps.hpp.


Member Function Documentation

template<class SparseOpsType>
void TestSparseOps< SparseOpsType >::makeDenseTestProblem ( Teuchos::RCP< dense_matrix_type > &  A_out,
Teuchos::RCP< dense_matrix_type > &  L_out,
Teuchos::RCP< dense_matrix_type > &  U_out,
Teuchos::Array< ordinal_type > &  pivots,
const ordinal_type  numRows 
) const [inline]

Compute dense test matrices for sparse triangular solve.

To make nonsingular lower and upper triangular matrices for testing sparse triangular solve, we start with a dense random matrix A and compute its LU factorization using LAPACK. Random matrices tend to be well conditioned, so the L and U factors will also be well conditioned. We output the dense matrices here so that you can test sparse routines by comparing their results with those of using the BLAS and LAPACK.

Parameters:
A_out[out] A numRows by numRows random matrix, of which L_out, U_out is the LU factorization with permutation pivots.
L_out[out] The L factor in the LU factorization of A_out.
U_out[out] The U factor in the LU factorization of A_out.
pivots[out] The permutation array in the LU factorization of A_out. Row i in the matrix A was interchanged with row pivots[i].
numRows[in] The number of rows and columns in A_out.

If you're solving AX=B using the LU factorization, you first have to apply the row permutation to B. You can do this in the same way that the reference _GETRS implementation does, using the _LASWP BLAS routine. Here's how the Fortran calls it (where _ is replaced by D for the special case of scalar_type=double):

CALL DLASWP( NRHS, B, LDB, 1, N, IPIV, 1 )

Then, after the triangular solves, apply the reverse-direction permutation (last parameter is -1 instead of 1) to the solution vector X:

CALL DLASWP( NRHS, X, LDX, 1, N, IPIV, -1 )

Definition at line 139 of file DefaultSparseOps_TestSparseOps.hpp.

template<class SparseOpsType>
void TestSparseOps< SparseOpsType >::readFile ( size_t &  numRows,
size_t &  numCols,
Teuchos::ArrayRCP< const size_t > &  rowptr,
Teuchos::ArrayRCP< const ordinal_type > &  colind,
Teuchos::ArrayRCP< const scalar_type > &  values,
const std::string &  filename 
) const [inline]

Read a CSR-format sparse matrix from a Matrix Market file.

CSR stands for "compressed sparse row." It's the desired input format for Kokkos' sparse kernels interface.

Parameters:
numRows[out] The number of rows in the sparse matrix.
numCols[out] The number of columns in the sparse matrix.
rowptr[out] The first of the three CSR arrays. For row r (zero-based row indices), rowptr[r] .. rowptr[r+1]-1 give the index range of colind and values for the entries of that row.
colind[out] The second of the three CSR arrays; the column indices.
values[out] The third of the three CSR arrays; the values of the matrix.
filename[in] The name of the Matrix Market - format file to read.

Definition at line 418 of file DefaultSparseOps_TestSparseOps.hpp.

template<class SparseOpsType>
Teuchos::RCP<SparseOpsType> TestSparseOps< SparseOpsType >::makeSparseOpsFromFile ( ordinal_type &  numRows,
ordinal_type &  numCols,
const Teuchos::RCP< node_type > &  node,
const std::string &  filename,
const Teuchos::EUplo  uplo = Teuchos::UNDEF_TRI,
const Teuchos::EDiag  diag = Teuchos::NON_UNIT_DIAG 
) const [inline]

Initialize and return a sparse kernels object from a Matrix Market file.

A Kokkos sparse kernels object turns a sparse matrix (represented in compressed sparse row format, more or less) into an opaque implementation of sparse matrix-(multi)vector multiply and sparse triangular solve.

Parameters:
numRows[out] Number of rows in the sparse matrix. SparseOpsType doesn't currently require a method that tells you the number of rows in the sparse matrix, so we output it here for later use.
numCols[out] Number of columns in the sparse matrix. SparseOpsType doesn't currently require a method that tells you the number of columns in the sparse matrix, so we output it here for later use.
node[in/out] Kokkos Node instance, which the constructors of graph_type and SparseOpsType require.
filename[in] Name of a Matrix Market sparse matrix file.
uplo[in] Whether the matrix is lower triangular (LOWER_TRI), upper triangular (UPPER_TRI), or neither (UNDEF_TRI). The latter is the default.
diag[in] Whether the matrix has an implicitly stored unit diagonal (UNIT_DIAG) or not (NON_UNIT_DIAG). The latter is the default. Kokkos' convention is to ignore this for sparse matrix-vector multiply, and respect it only for triangular solves.

Definition at line 491 of file DefaultSparseOps_TestSparseOps.hpp.

template<class SparseOpsType>
Teuchos::RCP<SparseOpsType> TestSparseOps< SparseOpsType >::makeSparseOps ( const Teuchos::RCP< node_type > &  node,
const Teuchos::ArrayRCP< const size_t > &  ptr,
const Teuchos::ArrayRCP< const ordinal_type > &  ind,
const Teuchos::ArrayRCP< const scalar_type > &  val,
const Teuchos::EUplo  uplo = Teuchos::UNDEF_TRI,
const Teuchos::EDiag  diag = Teuchos::NON_UNIT_DIAG 
) const [inline]

Initialize and return a sparse kernels object.

A Kokkos sparse kernels object turns a sparse matrix (represented in compressed sparse row format, more or less) into an opaque implementation of sparse matrix-(multi)vector multiply and sparse triangular solve.

Parameters:
node[in/out] Kokkos Node instance, which the constructors of graph_type and SparseOpsType require.
ptr[in] The first of the three CSR arrays. For row r (zero-based row indices), ptr[r] .. ptr[r+1]-1 give the index range of ind and val for the entries of that row.
ind[in] The second of the three CSR arrays; the column indices.
val[in] The third of the three CSR arrays; the values of the matrix.
uplo[in] Whether the matrix is lower triangular (LOWER_TRI), upper triangular (UPPER_TRI), or neither (UNDEF_TRI). The latter is the default.
diag[in] Whether the matrix has an implicitly stored unit diagonal (UNIT_DIAG) or not (NON_UNIT_DIAG). The latter is the default. Kokkos' convention is to ignore this for sparse matrix-vector multiply, and respect it only for triangular solves.

After calling this method, you can set ptr, ind, and val to null. This may free memory if the SparseOpsType copies into its own internal format instead of just using the original arrays.

Definition at line 542 of file DefaultSparseOps_TestSparseOps.hpp.

template<class SparseOpsType>
Teuchos::RCP<SparseOpsType> TestSparseOps< SparseOpsType >::denseTriToSparseOps ( const dense_matrix_type A,
const Teuchos::RCP< node_type > &  node,
const Teuchos::EUplo  uplo,
const Teuchos::EDiag  diag 
) const [inline]

Convert the given dense triangular matrix to a sparse kernels object.

Parameters:
A[in] The dense (lower or upper) triangular matrix.
node[in/out] The Kokkos Node instance.
uplo[in] Whether A is lower (LOWER_TRI) or upper (UPPER_TRI) triangular.
diag[in] Whether A has an implicit unit diagonal (UNIT_DIAG) or not (NON_UNIT_DIAG).

Definition at line 582 of file DefaultSparseOps_TestSparseOps.hpp.

template<class SparseOpsType>
Teuchos::RCP<SparseOpsType> TestSparseOps< SparseOpsType >::denseToSparseOps ( const dense_matrix_type A,
const Teuchos::RCP< node_type > &  node 
) const [inline]

Convert the given dense matrix to a sparse kernels object.

Parameters:
A[in] The dense matrix.
node[in/out] The Kokkos Node instance.

Definition at line 654 of file DefaultSparseOps_TestSparseOps.hpp.

template<class SparseOpsType>
Teuchos::RCP<Kokkos::MultiVector<scalar_type, node_type> > TestSparseOps< SparseOpsType >::makeMultiVector ( const Teuchos::RCP< node_type > &  node,
const ordinal_type  numRows,
const ordinal_type  numCols 
) const [inline]

Return an initialized Kokkos::MultiVector, filled with zeros.

Parameters:
node[in] The Kokkos Node instance.
numRows[in] The number of rows in the MultiVector.
numCols[in] The number of columns in the MultiVector.

Definition at line 712 of file DefaultSparseOps_TestSparseOps.hpp.

template<class SparseOpsType>
magnitude_type TestSparseOps< SparseOpsType >::maxRelativeError ( const Teuchos::RCP< const Kokkos::MultiVector< scalar_type, node_type > > &  X,
const Teuchos::RCP< const Kokkos::MultiVector< scalar_type, node_type > > &  Y,
const Teuchos::RCP< Kokkos::MultiVector< scalar_type, node_type > > &  Z 
) const [inline]

Return the maximum relative error between the two multivectors.

We define "maximum relative error" between X and Y as \[ \| X_i - Y_i \|-2 / \| X_i \|_2, \] where $X_i$ indicates the i-th column of X.

Parameters:
X[in] The "correct" multivector, against which to compare Y.
Y[in] The "test" multivector, which we hope is equal or close to X.
Z[in] A "scratch" multivector to use for storing X - Y.

Definition at line 746 of file DefaultSparseOps_TestSparseOps.hpp.

template<class SparseOpsType>
void TestSparseOps< SparseOpsType >::testSparseOps ( const Teuchos::RCP< node_type > &  node,
const ordinal_type  N,
const ordinal_type  numVecs,
const magnitude_type  tol 
) const [inline]

Test sparse matrix-(multi)vector multiply and sparse triangular solve.

Parameters:
node[in/out] The Kokkos Node instance.
N[in] The number of rows (and columns) in the sparse matrices to test.
numVecs[in] The number of columns in the multivectors to test.
tol[in] Tolerance for relative errors.

Test methodology ================

Summary: Compute a dense LU factorization of a random matrix A. Store its L and U factors as sparse matrices. Use the factors to test sparse matrix-vector multiply and sparse triangular solve.

1. Make dense ${L}$ and ${U}$ matrices in the way described above. ${L}$ is upper triangular with unit diagonal, and ${U}$ is lower triangular. 2. Compute a random multivector X. Give X > 1 column, to make sure that the routine can handle multivectors correctly. 3. Compute ${Y} = {L} ({U} X)$ using the BLAS. 4. Convert ${L}$ to L and ${U}$ to U, where L and U are "sparse" matrices. We put "sparse" in quotes because they are actually dense, but stored sparsely. 5. Test sparse matrix-(multi)vector multiply by computing $Y = L (U X)$ and making sure that $\|Y - {Y}\| / \|{Y}\| $ for some reasonable tolerance $$. 6. Test sparse triangular solve: a. Compute ${Z} = {L}^{-1} {Y}$ and $Z = L^{-1} {Y}$, and compare the results. b. Compute ${W} = {U}^{-1} {Z}$ and $W = L^{-1} {Z}$, and compare the results.

Discussion ==========

Random square matrices tend to be well conditioned on average (this statement sounds vague, but can be made rigorous), so the L and U factors are likely to exist and be well conditioned on average. Of course, outliers are possible, so this test isn't 100% guaranteed to succeed, but success is very likely and passing falsely is even less likely.

Definition at line 824 of file DefaultSparseOps_TestSparseOps.hpp.

template<class SparseOpsType>
void TestSparseOps< SparseOpsType >::benchmarkSparseMatVecFromFile ( std::vector< std::pair< std::string, double > > &  results,
const std::string &  filename,
const std::string &  label,
const Teuchos::RCP< node_type > &  node,
const ordinal_type  numVecs,
const int  numTrials 
) const [inline]

Benchmark sparse mat-vec with a matrix read from a Matrix Market file.

Parameters:
results[out] Timing results. You can also use TimeMonitor::report() or TimeMonitor::summarize() to display results.
filename[in] Name of the Matrix Market sparse matrix file.
label[in] Label to distinguish the SparseOpsType, in case you are running multiple benchmarks with different SparseOpsType types.
node[in/out] Kokkos Node instance.
numVecsNumber of columns in the multivectors to benchmark.
numTrialsNumber of trials. We time a loop around all the trials to increase accuracy and smooth out variance.

Definition at line 1105 of file DefaultSparseOps_TestSparseOps.hpp.


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