Sparse matrices are typically represented by storing only their non-zero entries. Along with the actual matrix values, their row and column information needs to be explicitly stored. For example, the sparse matrix - - | 1 1 0 0 0 0 0 | | 0 0 0 0 0 0 0 | | 0 0 3 3 3 0 0 | | 0 6 0 0 0 4 0 | | 2 0 0 0 0 5 5 | - - may be stored column-by-column by the three arrays: Column pointers (colptrs) : 0 2 4 5 6 7 9 10 Row indices (rowindices) : 0 4 0 3 2 2 2 3 4 4 Non-zero values (val) : 1 2 1 6 3 3 3 4 5 5 (In the above, note that all indices start at 0 (C-style) ) Assume the sparse matrix has m rows and n columns. "colptrs" is an integer array of size n+1. This array indexes into the "rowindices" and "val" arrays. In particular, the i-th column starts at colptrs(i) and ends at colptrs(i+1)-1. For programming convenience, colptrs(n+1) is set to nz, where nz is the number of nonzeros in the matrix. "rowindices" is an integer array of size nz (the number of nonzeros in the matrix), and it contains the row indices of the matrix elements. For example, rowindices[colptrs[i]] gives the row index of the 1st non-zero element in the i-th column. "val" is an array of size nz (the number of nonzeros in the matrix), and it contains the matrix values. Whether val is an integer, float or double (or whatever) array depends on the type of your matrix elements. For example, val[colptrs[i]] gives the value of the 1st non-zero element in the i-th column. Since many of you will be exchanging sparse matrices, I would like you to input and output the m x n sparse matrices (to files) in the following form: 5 7 10 % m, n, nz % Note that m is number of rows % and n is the number of columns 0 2 4 5 6 7 9 10 % colptrs 0 4 0 3 2 2 2 3 4 4 % rowindices 1.0 2.0 1.0 6.0 3.0 3.0 3.0 4.0 5.0 5.0 % val The following is an example of using this sparse matrix format: void at_times_x(int n, int nz, int *colptrs, int *rowindices, double *val, double *x, double *y) /* Subroutine at_times_x computes y = A^T * x, A is mxn and stored as a sparse matrix, and y and x are dense vectors */ { double yi; int i,j; for(i = 0;i < n;i++){ yi = 0.0; for(j = colptrs[i];j < colptrs[i+1];j++){ yi += val[j] * x[rowindices[j]]; } y[i] = yi; } }