Working with sparse data in RandBLAS

Sparse matrix data structures

The common interface for our sparse matrix types

enum class RandBLAS::sparse_data::IndexBase : int

Values:

enumerator Zero

enumerator One

template<typename SpMat> concept SparseMatrix

#include <base.hh>

An object $\texttt{M}$ of type $\texttt{SpMat}$ has the following attributes.

	type	description
$\texttt{M.n_rows}$	$\texttt{const int64_t}$	number of rows
$\texttt{M.n_cols}$	$\texttt{const int64_t}$	number of columns
$\texttt{M.nnz}$	$\texttt{int64_t}$	number of structural nonzeros
$\texttt{M.vals}$	$\texttt{SpMat::scalar_t *}$	pointer to values of structural nonzeros
$\texttt{M.own_memory}$	$\texttt{bool}$	A flag indicating if memory attached to $\texttt{M}$ should be deallocated when $\texttt{M}$ is deleted. This flag is set automatically based on the type of constructor used for $\texttt{M}.$

Memory-owning constructors

$\texttt{SpMat}$ must have a constructor for an empty matrix of given dimensions. Conformant implementations of this constructor look like the following.
SpMat(int64_t n_rows, int64_t n_cols) 
 : n_rows(n_rows), n_cols(n_cols), nnz(0), vals(nullptr), own_memory(true) {
     // class-specific code ...
 };
If we construct $\texttt{SpMat M(m, n)},$ then we can’t store data in $\texttt{M}$ until a function call of the form $\texttt{M.reserve(nnz)}.$ Here’s an outline of a conformant implementation of this function.
void reserve(int64_t nnz) {
    assert this->own_memory;
    this->nnz = nnz;
    this->vals = new SpMat::scalar_t[nnz];
    // ... class-specific code ...
}
The destructor of $\texttt{M}$ is responsible for deallocating $\texttt{M.vals}$ and other attached data. A conformant implemnentation of the destructor will look like the following.
~SpMat() {
    if (own_memory) {
        if (vals != nullptr) delete [] vals;
        // ... class-specific code ...
    }
}        

View constructors

This concept doesn’t place requirements on constructors for sparse matrix views of existing data. However, all of RandBLAS’ sparse matrix classes offer such constructors. See individual classes’ documentation for details.

COOMatrix

template<typename T, SignedInteger sint_t = int64_t> struct COOMatrix

Let ${\mathbf{A}}$ denote a sparse matrix with ${\texttt{nnz}}$ structural nonzeros. Its COO representation consists of declared dimensions, ${\texttt{n_rows}}$ and ${\texttt{n_cols}}$ , as well as a triplet of arrays ${(\texttt{vals},\texttt{rows},\texttt{cols})}$ where

\[\mathbf{A}_{\texttt{rows}[\ell],\texttt{cols}[\ell]} = \texttt{vals}[\ell] \quad\text{for all}\quad \ell \in \{0,\ldots,\texttt{nnz}-1\}.\]

This type conforms to the SparseMatrix concept.

Public Types

using scalar_t = T : Real scalar type used for structural nonzeros in this matrix.

using index_t = sint_t : Signed integer type used in the rows and cols array members.

Public Functions

inline COOMatrix(int64_t n_rows, int64_t n_cols)

Standard constructor. Initializes n_rows and n_cols at the provided values. The vals, rows, and cols members are set to null pointers; nnz is set to Zero, index_base is set to zero, and COOMatrix::own_memory is set to true.

This constructor is intended for use with reserve_coo(int64_t nnz, COOMatrix &A).

inline COOMatrix(int64_t n_rows, int64_t n_cols, int64_t nnz, T *vals, sint_t *rows, sint_t *cols, bool compute_sort_type = true, IndexBase index_base = IndexBase::Zero): Expert constructor. Arguments passed to this function are used to initialize members of the same names; COOMatrix::own_memory is set to false. If compute_sort_type is true, then the sort member will be computed by inspecting the contents of (rows, cols). If compute_sort_type is false, then the sort member is set to None.

Public Members

const int64_t n_rows: The number of rows in this sparse matrix.

const int64_t n_cols: The number of columns in this sparse matrix.

bool own_memory

If true, then RandBLAS has permission to allocate and attach memory to the reference members of this matrix (vals, rows, and cols). If true at destruction time, then delete [] will be called on each non-null reference member of this matrix.

RandBLAS only writes to this member at construction time.

int64_t nnz: The number of structral nonzeros in this matrix.

IndexBase index_base: A flag to indicate whether (rows, cols) are interpreted with zero-based or one-based indexing.

T *vals

Reference to an array that holds values for the structural nonzeros of this matrix.

If non-null, this must have length at least nnz.

sint_t *rows

Reference to an array that holds row indices for the structural nonzeros of this matrix.

If non-null, this must have length at least nnz.

sint_t *cols

Reference to an array that holds column indicies for the structural nonzeros of this matrix.

If non-null, this must have length at least nnz.

NonzeroSort sort: A flag to indicate if the data in (vals, rows, cols) is sorted in a CSC-like order, a CSR-like order, or neither order.

template<typename T, SignedInteger sint_t> void RandBLAS::sparse_data::reserve_coo(int64_t nnz, COOMatrix<T, sint_t> &M)

This function requires that M.own_memory is true and that M.vals, M.rows, and M.cols are all null. If any of these conditions aren’t satisfied then RandBLAS will raise an error.

If no error is raised, then M.nnz is overwritten by nnz, M.vals is redirected to a new length-nnz array of type T, and (M.rows, M.cols) are redirected to new length-nnz arrays of type sint_t.

enum class RandBLAS::sparse_data::NonzeroSort : char

Indicates whether the (vals, rows, cols) data of a COO-format sparse matrix are known to be sorted in CSC order, CSR order, or neither of those orders.

Values:

enumerator CSC

enumerator CSR

enumerator None

CSRMatrix

template<typename T, SignedInteger sint_t = int64_t> struct CSRMatrix

Let ${\mathbf{A}}$ denote a sparse matrix with ${\texttt{nnz}}$ structural nonzeros. Its CSR representation consists of declared dimensions, ${\texttt{n_rows}}$ and ${\texttt{n_cols}}$ , and a triplet of arrays ${(\texttt{vals},\texttt{rowptr},\texttt{colidxs}).}$

The ${\texttt{i}^{\text{th}}}$ row of ${\mathbf{A}}$ has ${\texttt{rowptr[i+1] - rowptr[i]}}$ structural nonzeros. The ${\texttt{k}^{\text{th}}}$ structural nonzero in row ${\texttt{i}}$ appears in column ${\texttt{colidxs[rowptr[i] + k]}}$ and is equal to ${\texttt{vals[rowptr[i] + k]}.}$

This type conforms to the SparseMatrix concept.

Public Types

using scalar_t = T : Real scalar type used for structural nonzeros in this matrix.

using index_t = sint_t : Signed integer type used in the rowptr and colidxs array members.

Public Functions

inline CSRMatrix(int64_t n_rows, int64_t n_cols)

Standard constructor. Initializes n_rows and n_cols at the provided values. The vals, rowptr, and colidxs members are set to null pointers; nnz is set to zero, index_base is set to Zero, and CSRMatrix::own_memory is set to true.

This constructor is intended for use with reserve_csr(int64_t nnz, CSRMatrix &A).

inline CSRMatrix(int64_t n_rows, int64_t n_cols, int64_t nnz, T *vals, sint_t *rowptr, sint_t *colidxs, IndexBase index_base = IndexBase::Zero): Expert constructor. Arguments passed to this function are used to initialize members of the same names; CSRMatrix::own_memory is set to false.

Public Members

const int64_t n_rows: The number of rows in this sparse matrix.

const int64_t n_cols: The number of columns in this sparse matrix.

bool own_memory

If true, then RandBLAS has permission to allocate and attach memory to the reference members of this matrix (vals, rowptr, colidxs). If true at destruction time, then delete [] will be called on each non-null reference member of this matrix.

RandBLAS only writes to this member at construction time.

int64_t nnz: The number of structral nonzeros in this matrix.

IndexBase index_base: A flag to indicate whether colidxs is interpreted with zero-based or one-based indexing.

T *vals

Reference to an array that holds values for the structural nonzeros of this matrix.

If non-null, this must have length at least nnz.

sint_t *rowptr

Reference to a pointer offset array for the CSR format.

If non-null, then must have length at least ${\texttt{n_rows + 1}}$ .

sint_t *colidxs

Reference to a column index array in the CSR format, interpreted in ${\texttt{index_base}}$ .

If non-null, then must have length at least nnz.

template<typename T, SignedInteger sint_t> void RandBLAS::sparse_data::reserve_csr(int64_t nnz, CSRMatrix<T, sint_t> &M)

This function requires that M.own_memory is true, that M.colidxs is null, and that M.vals is null. If any of these conditions are not met then this function will raise an error.

Special logic applies to M.rowptr because its documented length requirement is determined by the const variable M.n_rows.

If M.rowptr is non-null then it is left unchanged, and it is presumed to point to an array of length at least M.n_rows + 1.
If M.rowptr is null, then it will be redirected to a new array of type sint_t and length (M.n_rows + 1).

From there, M.nnz is overwritten by nnz, and the reference members M.colidxs and M.vals are redirected to new arrays of length nnz (of types sint_t and T, respectively).

CSCMatrix

template<typename T, SignedInteger sint_t = int64_t> struct CSCMatrix

Let ${\mathbf{A}}$ denote a sparse matrix with ${\texttt{nnz}}$ structural nonzeros. Its CSC representation consists of declared dimensions, ${\texttt{n_rows}}$ and ${\texttt{n_cols}}$ , and a triplet of arrays ${(\texttt{vals},\texttt{rowidxs},\texttt{colptr}).}$

The ${\texttt{j}^{\text{th}}}$ column of ${\mathbf{A}}$ has ${\texttt{colptr[j+1] - colptr[j]}}$ structural nonzeros. The ${\texttt{k}^{\text{th}}}$ structural nonzero in column ${\texttt{j}}$ appears in row ${\texttt{rowidxs[colptr[j] + k]}}$ and is equal to ${\texttt{vals[colptr[j] + k]}.}$

This type conforms to the SparseMatrix concept.

Public Types

using scalar_t = T : Real scalar type used for structural nonzeros in this matrix.

using index_t = sint_t : Signed integer type used in the rowptr and colidxs array members.

Public Functions

inline CSCMatrix(int64_t n_rows, int64_t n_cols)

Standard constructor. Initializes n_rows and n_cols at the provided values. The vals, rowidxs, and colptr members are null-initialized; ${\texttt{nnz}}$ is set to zero, ${\texttt{index_base}}$ is set to Zero, and CSCMatrix::own_memory is set to true.

This constructor is intended for use with reserve_csc(int64_t nnz, CSCMatrix &A).

inline CSCMatrix(int64_t n_rows, int64_t n_cols, int64_t nnz, T *vals, sint_t *rowidxs, sint_t *colptr, IndexBase index_base = IndexBase::Zero): Expert constructor. Arguments passed to this function are used to initialize members of the same names; CSCMatrix::own_memory is set to false.

Public Members

const int64_t n_rows: The number of rows in this sparse matrix.

const int64_t n_cols: The number of columns in this sparse matrix.

bool own_memory

If true, then RandBLAS has permission to allocate and attach memory to the reference members of this matrix (vals, rowidxs, colptr). If true at destruction time, then delete [] will be called on each non-null reference member of this matrix.

RandBLAS only writes to this member at construction time.

int64_t nnz: The number of structral nonzeros in this matrix.

IndexBase index_base: A flag to indicate whether colidxs is interpreted with zero-based or one-based indexing.

T *vals

Reference to an array that holds values for the structural nonzeros of this matrix.

If non-null, this must have length at least nnz.

sint_t *rowidxs

Reference to a row index array in the CSC format, interpreted in ${\texttt{index_base}}$ .

If non-null, then must have length at least nnz.

sint_t *colptr

Reference to a pointer offset array for the CSC format.

If non-null, then must have length at least ${\texttt{n_cols + 1}}$ .

template<typename T, SignedInteger sint_t> void RandBLAS::sparse_data::reserve_csc(int64_t nnz, CSCMatrix<T, sint_t> &M)

This function requires that M.own_memory is true, that M.rowidxs is null, and that M.vals is null. If any of these conditions are not met then this function will raise an error.

Special logic applies to M.colptr because its documented length requirement is determined by the const variable M.n_cols.

If M.colptr is non-null then it is left unchanged, and it is presumed to point to an array of length at least M.n_cols + 1.
If M.colptr is null, then it will be redirected to a new array of type sint_t and length (M.n_cols + 1).

From there, M.nnz is overwritten by nnz, and the reference members M.rowidxs and M.vals are redirected to new arrays of length nnz (of types sint_t and T, respectively).

Operations with sparse matrices

Sketching

$\mathbf{B} = \alpha \cdot \operatorname{op}(\operatorname{submat}(\mathbf{S}))\cdot \operatorname{op}(\mathbf{A}) + \beta \cdot \mathbf{B}$

template<SparseMatrix SpMat, typename DenseSkOp, typename T = DenseSkOp::scalar_t> inline void RandBLAS::sketch_sparse(blas::Layout layout, blas::Op opS, blas::Op opA, int64_t d, int64_t n, int64_t m, T alpha, DenseSkOp &S, int64_t ro_s, int64_t co_s, SpMat &A, T beta, T *B, int64_t ldb)

Sketch from the left in an SpMM-like operation

\[\operatorname{mat}(B) = \alpha \cdot \underbrace{\operatorname{op}(\operatorname{submat}(\mathbf{S}))}_{d \times m} \cdot \underbrace{\operatorname{op}(\mathbf{A})}_{m \times n} + \beta \cdot \underbrace{\operatorname{mat}(B)}_{d \times n}, \tag{$\star$}\]

where $\alpha$ and $\beta$ are real scalars, $\operatorname{op}(\mathbf{X})$ either returns a matrix $\mathbf{X}$ or its transpose, $\mathbf{A}$ is a sparse matrix, and $\mathbf{S}$ is a dense sketching operator.

$\mathbf{B} = \alpha \cdot \operatorname{op}(\mathbf{A})\cdot \operatorname{op}(\operatorname{submat}(\mathbf{S})) + \beta \cdot \mathbf{B}$

template<SparseMatrix SpMat, typename DenseSkOp, typename T = DenseSkOp::scalar_t> inline void RandBLAS::sketch_sparse(blas::Layout layout, blas::Op opA, blas::Op opS, int64_t m, int64_t d, int64_t n, T alpha, SpMat &A, DenseSkOp &S, int64_t ro_s, int64_t co_s, T beta, T *B, int64_t ldb)

Sketch from the right in an SpMM-like operation

\[\operatorname{mat}(B) = \alpha \cdot \underbrace{\operatorname{op}(\mathbf{A})}_{m \times n} \cdot \underbrace{\operatorname{op}(\operatorname{submat}(\mathbf{S}))}_{n \times d} + \beta \cdot \underbrace{\operatorname{mat}(B)}_{m \times d}, \tag{$\star$}\]

where $\alpha$ and $\beta$ are real scalars, $\operatorname{op}(\mathbf{X})$ either returns a matrix $\mathbf{X}$ or its transpose, $\mathbf{A}$ is a sparse matrix, and $\mathbf{S}$ is a dense sketching operator.

Deterministic operations

$\mathbf{C} = \alpha \cdot \operatorname{op}(\mathbf{A})\cdot \operatorname{op}(\mathbf{B}) + \beta \cdot \mathbf{C},$ with sparse $\mathbf{A}$

template<SparseMatrix SpMat, typename T = SpMat::scalar_t> inline void RandBLAS::spmm(blas::Layout layout, blas::Op opA, blas::Op opB, int64_t m, int64_t n, int64_t k, T alpha, SpMat &A, const T *B, int64_t ldb, T beta, T *C, int64_t ldc)

Perform an SPMM-like operation, multiplying a dense matrix on the left with a sparse matrix:

\[\operatorname{mat}(C) = \alpha \cdot \underbrace{\operatorname{op}(\mathbf{A})}_{m \times k} \cdot \underbrace{\operatorname{op}(\operatorname{mat}(B))}_{k \times n} + \beta \cdot \underbrace{\operatorname{mat}(C)}_{m \times n}, \tag{$\star$}\]

where $\alpha$ and $\beta$ are real scalars, $\operatorname{op}(\mathbf{X})$ either returns a matrix $\mathbf{X}$ or its transpose, and $\mathbf{A}$ is a sparse matrix.

$\mathbf{C} = \alpha \cdot \operatorname{op}(\mathbf{A})\cdot \operatorname{op}(\mathbf{B}) + \beta \cdot \mathbf{C},$ with sparse $\mathbf{B}$

template<SparseMatrix SpMat, typename T = SpMat::scalar_t> inline void RandBLAS::spmm(blas::Layout layout, blas::Op opA, blas::Op opB, int64_t m, int64_t n, int64_t k, T alpha, const T *A, int64_t lda, SpMat &B, T beta, T *C, int64_t ldc)

Perform an SPMM-like operation, multiplying a dense matrix on the right with a (submatrix of a) sparse matrix:

\[\operatorname{mat}(C) = \alpha \cdot \underbrace{\operatorname{op}(\operatorname{mat}(A))}_{m \times k} \cdot \underbrace{\operatorname{op}(\mathbf{B})}_{k \times n} + \beta \cdot \underbrace{\operatorname{mat}(C)}_{m \times n}, \tag{$\star$}\]

where $\alpha$ and $\beta$ are real scalars, $\operatorname{op}(\mathbf{X})$ either returns a matrix $\mathbf{X}$ or its transpose, and $\mathbf{B}$ is a sparse matrix.

	type	description
\(\texttt{M.n_rows}\)	\(\texttt{const int64_t}\)	number of rows
\(\texttt{M.n_cols}\)	\(\texttt{const int64_t}\)	number of columns
\(\texttt{M.nnz}\)	\(\texttt{int64_t}\)	number of structural nonzeros
\(\texttt{M.vals}\)	\(\texttt{SpMat::scalar_t *}\)	pointer to values of structural nonzeros
\(\texttt{M.own_memory}\)	\(\texttt{bool}\)	A flag indicating if memory attached to \(\texttt{M}\) should be deallocated when \(\texttt{M}\) is deleted. This flag is set automatically based on the type of constructor used for \(\texttt{M}.\)