41.1 Overview

These routines provide support for constructing and manipulating sparse matrices in GSL, using an API similar to gsl_matrix. The basic structure is called gsl_spmatrix. There are three supported storage formats for sparse matrices: the triplet, compressed column storage (CCS) and compressed row storage (CRS) formats. The triplet format stores triplets (i,j,x) for each non-zero element of the matrix. This notation means that the (i,j) element of the matrix A is A_{ij} = x. Compressed column storage stores each column of non-zero values in the sparse matrix in a continuous memory block, keeping pointers to the beginning of each column in that memory block, and storing the row indices of each non-zero element. Compressed row storage stores each row of non-zero values in a continuous memory block, keeping pointers to the beginning of each row in the block and storing the column indices of each non-zero element. The triplet format is ideal for adding elements to the sparse matrix structure while it is being constructed, while the compressed storage formats are better suited for matrix-matrix multiplication or linear solvers.

The gsl_spmatrix structure is defined as

typedef struct
{
  size_t size1;
  size_t size2;
  size_t *i;
  double *data;
  size_t *p;
  size_t nzmax;
  size_t nz;
  gsl_spmatrix_tree *tree_data;
  void *work;
  size_t sptype;
} gsl_spmatrix;

This defines a size1-by-size2 sparse matrix. The number of non-zero elements currently in the matrix is given by nz. For the triplet representation, i, p, and data are arrays of size nz which contain the row indices, column indices, and element value, respectively. So if data[k] = A(i,j), then i = i[k] and j = p[k].

For compressed column storage, i and data are arrays of size nz containing the row indices and element values, identical to the triplet case. p is an array of size size2 + 1 where p[j] points to the index in data of the start of column j. Thus, if data[k] = A(i,j), then i = i[k] and p[j] <= k < p[j+1].

For compressed row storage, i and data are arrays of size nz containing the column indices and element values, identical to the triplet case. p is an array of size size1 + 1 where p[i] points to the index in data of the start of row i. Thus, if data[k] = A(i,j), then j = i[k] and p[i] <= k < p[i+1].

The parameter tree_data is a binary tree structure used in the triplet representation, specifically a balanced AVL tree. This speeds up element searches and duplicate detection during the matrix assembly process. The parameter work is additional workspace needed for various operations like converting from triplet to compressed storage. sptype indicates the type of storage format being used (triplet, CCS or CRS).

The compressed storage format defined above makes it very simple to interface with sophisticated external linear solver libraries which accept compressed storage input. The user can simply pass the arrays i, p, and data as the inputs to external libraries.