Next: Sparse Matrices References and Further Reading, Previous: Sparse Matrices Conversion Between Sparse and Dense, Up: Sparse Matrices [Index]
The following example program builds a 5-by-4 sparse matrix and prints it in triplet, compressed column, and compressed row format. The matrix which is constructed is The output of the program is
printing all matrix elements: A(0,0) = 0 A(0,1) = 0 A(0,2) = 3.1 A(0,3) = 4.6 A(1,0) = 1 . . . A(4,0) = 4.1 A(4,1) = 0 A(4,2) = 0 A(4,3) = 0 matrix in triplet format (i,j,Aij): (0, 2, 3.1) (0, 3, 4.6) (1, 0, 1.0) (1, 2, 7.2) (3, 0, 2.1) (3, 1, 2.9) (3, 3, 8.5) (4, 0, 4.1) matrix in compressed column format: i = [ 1, 3, 4, 3, 0, 1, 0, 3, ] p = [ 0, 3, 4, 6, 8, ] d = [ 1, 2.1, 4.1, 2.9, 3.1, 7.2, 4.6, 8.5, ] matrix in compressed row format: i = [ 2, 3, 0, 2, 0, 1, 3, 0, ] p = [ 0, 2, 4, 4, 7, 8, ] d = [ 3.1, 4.6, 1, 7.2, 2.1, 2.9, 8.5, 4.1, ]
We see in the compressed column output, the data array stores each column contiguously, the array i stores the row index of the corresponding data element, and the array p stores the index of the start of each column in the data array. Similarly, for the compressed row output, the data array stores each row contiguously, the array i stores the column index of the corresponding data element, and the p array stores the index of the start of each row in the data array.
#include <stdio.h>
#include <stdlib.h>
#include <gsl/gsl_spmatrix.h>
int
main()
{
gsl_spmatrix *A = gsl_spmatrix_alloc(5, 4); /* triplet format */
gsl_spmatrix *B, *C;
size_t i, j;
/* build the sparse matrix */
gsl_spmatrix_set(A, 0, 2, 3.1);
gsl_spmatrix_set(A, 0, 3, 4.6);
gsl_spmatrix_set(A, 1, 0, 1.0);
gsl_spmatrix_set(A, 1, 2, 7.2);
gsl_spmatrix_set(A, 3, 0, 2.1);
gsl_spmatrix_set(A, 3, 1, 2.9);
gsl_spmatrix_set(A, 3, 3, 8.5);
gsl_spmatrix_set(A, 4, 0, 4.1);
printf("printing all matrix elements:\n");
for (i = 0; i < 5; ++i)
for (j = 0; j < 4; ++j)
printf("A(%zu,%zu) = %g\n", i, j,
gsl_spmatrix_get(A, i, j));
/* print out elements in triplet format */
printf("matrix in triplet format (i,j,Aij):\n");
gsl_spmatrix_fprintf(stdout, A, "%.1f");
/* convert to compressed column format */
B = gsl_spmatrix_ccs(A);
printf("matrix in compressed column format:\n");
printf("i = [ ");
for (i = 0; i < B->nz; ++i)
printf("%zu, ", B->i[i]);
printf("]\n");
printf("p = [ ");
for (i = 0; i < B->size2 + 1; ++i)
printf("%zu, ", B->p[i]);
printf("]\n");
printf("d = [ ");
for (i = 0; i < B->nz; ++i)
printf("%g, ", B->data[i]);
printf("]\n");
/* convert to compressed row format */
C = gsl_spmatrix_crs(A);
printf("matrix in compressed row format:\n");
printf("i = [ ");
for (i = 0; i < C->nz; ++i)
printf("%zu, ", C->i[i]);
printf("]\n");
printf("p = [ ");
for (i = 0; i < C->size1 + 1; ++i)
printf("%zu, ", C->p[i]);
printf("]\n");
printf("d = [ ");
for (i = 0; i < C->nz; ++i)
printf("%g, ", C->data[i]);
printf("]\n");
gsl_spmatrix_free(A);
gsl_spmatrix_free(B);
gsl_spmatrix_free(C);
return 0;
}