10.1. Description of the SUNMATRIX Modules
For problems that involve direct methods for solving linear systems,
the SUNDIALS packages not only operate on generic vectors, but also
on generic matrices (of type SUNMatrix
), through a set of
operations defined by the particular SUNMATRIX implementation.
Users can provide their own specific implementation of the
SUNMATRIX module, particularly in cases where they provide their
own N_Vector
and/or linear solver modules, and require matrices
that are compatible with those implementations. The generic
SUNMatrix
operations are described below, and descriptions of
the SUNMATRIX implementations provided with SUNDIALS follow.
The generic SUNMatrix
type has been modeled after the
objectoriented style of the generic N_Vector
type.
Specifically, a generic SUNMatrix
is a pointer to a structure
that has an implementationdependent content field containing
the description and actual data of the matrix, and an ops field
pointing to a structure with generic matrix operations.
The type SUNMatrix
is defined as:

typedef struct _generic_SUNMatrix *SUNMatrix
and the generic structure is defined as
struct _generic_SUNMatrix {
void *content;
struct _generic_SUNMatrix_Ops *ops;
};
Here, the _generic_SUNMatrix_Ops
structure is essentially a list of
function pointers to the various actual matrix operations, and is
defined as
struct _generic_SUNMatrix_Ops {
SUNMatrix_ID (*getid)(SUNMatrix);
SUNMatrix (*clone)(SUNMatrix);
void (*destroy)(SUNMatrix);
int (*zero)(SUNMatrix);
int (*copy)(SUNMatrix, SUNMatrix);
int (*scaleadd)(realtype, SUNMatrix, SUNMatrix);
int (*scaleaddi)(realtype, SUNMatrix);
int (*matvecsetup)(SUNMatrix);
int (*matvec)(SUNMatrix, N_Vector, N_Vector);
int (*space)(SUNMatrix, long int*, long int*);
};
The generic SUNMATRIX module defines and implements the matrix
operations acting on a SUNMatrix
. These routines are nothing but
wrappers for the matrix operations defined by a particular SUNMATRIX
implementation, which are accessed through the ops field of the
SUNMatrix
structure. To illustrate this point we show below the
implementation of a typical matrix operation from the generic
SUNMATRIX module, namely SUNMatZero
, which sets all values of a
matrix A
to zero, returning a flag denoting a successful/failed
operation:
int SUNMatZero(SUNMatrix A)
{
return((int) A>ops>zero(A));
}
§10.2 contains a complete list of all matrix operations defined by the generic SUNMATRIX module. A particular implementation of the SUNMATRIX module must:
Specify the content field of the
SUNMatrix
object.Define and implement a minimal subset of the matrix operations. See the documentation for each SUNDIALS package and/or linear solver to determine which SUNMATRIX operations they require.
Note that the names of these routines should be unique to that implementation in order to permit using more than one SUNMATRIX module (each with different
SUNMatrix
internal data representations) in the same code.Define and implement usercallable constructor and destructor routines to create and free a
SUNMatrix
with the new content field and with ops pointing to the new matrix operations.Optionally, define and implement additional usercallable routines acting on the newly defined
SUNMatrix
(e.g., a routine to print the content for debugging purposes).Optionally, provide accessor macros as needed for that particular implementation to be used to access different parts in the content field of the newly defined
SUNMatrix
.
To aid in the creation of custom SUNMATRIX modules the generic SUNMATRIX module
provides three utility functions SUNMatNewEmpty()
, SUNMatCopyOps()
,
and SUNMatFreeEmpty()
. When used in custom SUNMATRIX constructors and clone
routines these functions will ease the introduction of any new optional matrix
operations to the SUNMATRIX API by ensuring only required operations need to be
set and all operations are copied when cloning a matrix.

SUNMatrix SUNMatNewEmpty()
This function allocates a new generic
SUNMatrix
object and initializes its content pointer and the function pointers in the operations structure toNULL
. Return value:
If successful, this function returns a
SUNMatrix
object. If an error occurs when allocating the object, then this routine will returnNULL
.

int SUNMatCopyOps(SUNMatrix A, SUNMatrix B)
This function copies the function pointers in the
ops
structure ofA
into theops
structure ofB
. Arguments:
A – the matrix to copy operations from.
B – the matrix to copy operations to.
 Return value:
If successful, this function returns
0
. If either of the inputs areNULL
or theops
structure of either input isNULL
, then is function returns a nonzero value.

void SUNMatFreeEmpty(SUNMatrix A)
This routine frees the generic
SUNMatrix
object, under the assumption that any implementationspecific data that was allocated within the underlying content structure has already been freed. It will additionally test whether the ops pointer isNULL
, and, if it is not, it will free it as well. Arguments:
A – the SUNMatrix object to free
Each SUNMATRIX implementation included in SUNDIALS has a unique
identifier specified in enumeration and shown in
Table 10.1. It is recommended that a
usersupplied SUNMATRIX implementation use the SUNMATRIX_CUSTOM
identifier.
Matrix ID 
Matrix type 

SUNMATRIX_BAND 
Band \(M \times M\) matrix 
SUNMATRIX_CUSPARSE 
CUDA sparse CSR matrix 
SUNMATRIX_CUSTOM 
Userprovided custom matrix 
SUNMATRIX_DENSE 
Dense \(M \times N\) matrix 
SUNMATRIX_GINKGO 
SUNMatrix wraper for Ginkgo matrices 
SUNMATRIX_MAGMADENSE 
Dense \(M \times N\) matrix 
SUNMATRIX_ONEMKLDENSE 
oneMKL dense \(M \times N\) matrix 
SUNMATRIX_SLUNRLOC 
SUNMatrix wrapper for SuperLU_DIST SuperMatrix 
SUNMATRIX_SPARSE 
Sparse (CSR or CSC) \(M\times N\) matrix 
10.2. Description of the SUNMATRIX operations
For each of the SUNMatrix
operations, we give the name, usage
of the function, and a description of its mathematical operations
below.

SUNMatrix_ID SUNMatGetID(SUNMatrix A)
Returns the type identifier for the matrix A. It is used to determine the matrix implementation type (e.g. dense, banded, sparse,…) from the abstract
SUNMatrix
interface. This is used to assess compatibility with SUNDIALSprovided linear solver implementations. Returned values are given in Table 10.1Usage:
id = SUNMatGetID(A);

SUNMatrix SUNMatClone(SUNMatrix A)
Creates a new
SUNMatrix
of the same type as an existing matrix A and sets the ops field. It does not copy the matrix values, but rather allocates storage for the new matrix.Usage:
B = SUNMatClone(A);

void SUNMatDestroy(SUNMatrix A)
Destroys the
SUNMatrix
A and frees memory allocated for its internal data.Usage:
SUNMatDestroy(A);

int SUNMatSpace(SUNMatrix A, long int *lrw, long int *liw)
Returns the storage requirements for the matrix A. lrw contains the number of realtype words and liw contains the number of integer words. The return value denotes success/failure of the operation.
This function is advisory only, for use in determining a user’s total space requirements; it could be a dummy function in a usersupplied
SUNMatrix
module if that information is not of interest.Usage:
retval = SUNMatSpace(A, &lrw, &liw);

int SUNMatZero(SUNMatrix A)
Zeros all entries of the
SUNMatrix
A. The return value is an integer flag denoting success/failure of the operation:\[A_{i,j} = 0, \quad i=1,\ldots,m, \; j=1,\ldots,n.\]Usage:
retval = SUNMatZero(A);

int SUNMatCopy(SUNMatrix A, SUNMatrix B)
Performs the operation B gets A for all entries of the matrices A and B. The return value is an integer flag denoting success/failure of the operation:
\[B_{i,j} = A_{i,j}, \quad i=1,\ldots,m, \; j=1,\ldots,n.\]Usage:
retval = SUNMatCopy(A,B);

int SUNMatScaleAdd(realtype c, SUNMatrix A, SUNMatrix B)
Performs the operation A gets cA + B. The return value is an integer flag denoting success/failure of the operation:
\[A_{i,j} = cA_{i,j} + B_{i,j}, \quad i=1,\ldots,m, \; j=1,\ldots,n.\]Usage:
retval = SUNMatScaleAdd(c, A, B);

int SUNMatScaleAddI(realtype c, SUNMatrix A)
Performs the operation A gets cA + I. The return value is an integer flag denoting success/failure of the operation:
\[A_{i,j} = cA_{i,j} + \delta_{i,j}, \quad i,j=1,\ldots,n.\]Usage:
retval = SUNMatScaleAddI(c, A);

int SUNMatMatvecSetup(SUNMatrix A)
Performs any setup necessary to perform a matrixvector product. The return value is an integer flag denoting success/failure of the operation. It is useful for SUNMatrix implementations which need to prepare the matrix itself, or communication structures before performing the matrixvector product.
Usage:
retval = SUNMatMatvecSetup(A);

int SUNMatMatvec(SUNMatrix A, N_Vector x, N_Vector y)
Performs the matrixvector product y gets Ax. It should only be called with vectors x and y that are compatible with the matrix A – both in storage type and dimensions. The return value is an integer flag denoting success/failure of the operation:
\[y_i = \sum_{j=1}^n A_{i,j} x_j, \quad i=1,\ldots,m.\]Usage:
retval = SUNMatMatvec(A, x, y);
10.2.1. SUNMatrix return codes
The functions provided to SUNMatrix modules within the SUNDIALSprovided SUNMatrix implementations utilize a common set of return codes, listed below. These adhere to a common pattern: 0 indicates success, a negative value indicates a failure. Aside from this pattern, the actual values of each error code are primarily to provide additional information to the user in case of a SUNMatrix failure.
SUNMAT_SUCCESS
(0) – successful callSUNMAT_ILL_INPUT
(1) – an illegal input has been provided to the functionSUNMAT_MEM_FAIL
(2) – failed memory access or allocationSUNMAT_OPERATION_FAIL
(3) – a SUNMatrix operation returned nonzeroSUNMAT_MATVEC_SETUP_REQUIRED
(4) – theSUNMatMatvecSetup()
routine needs to be called prior to callingSUNMatMatvec()