Sparse matrices and storage formats

The problem: storing coefficients that are mostly zero

In a linear system $Ax=b$, each row of $A$ records how the unknowns contribute to one equation. Large systems can have many coefficients equal to zero. Storing every zero uses memory without contributing to the equations. The storage discussion in Templates for the Solution of Linear Systems, “Survey of Sparse Matrix Storage Formats” starts from this difficulty: retain the useful coefficients and enough information to recover their positions.

A second problem is exchanging those coefficients between programs. NIST’s Matrix Market provides test matrices for comparing numerical algorithms. Its exchange-format documentation supplies the published example used below. These sources document computational and exchange problems; they do not establish who first invented COO, CSR, or CSC. The numerical baseline is a published format illustration, not measured application data or a reconstruction of an original physical system. No physical units are specified.

For matrix notation and multiplication, see matrix factorization and decomposition. Here $A\in\mathbb R^{m\times n}$ means a real matrix with $m$ rows and $n$ columns. We number rows and columns from zero: $a_{ij}$ is the entry at row $i$, column $j$.

Published worked data: a matrix to exchange

The real, general $5\times5$ example in NIST’s “Matrix Market Exchange Formats” section is

All entries are preserved. The source’s coordinate file uses one-based indices; our arrays subtract one from each row and column index. Subsequent array conversions and products are modern worked derivations.

There are 25 positions, of which eight are nonzero. Its density is the nonzero fraction, $8/25=0.32$; its zero fraction is $17/25=0.68$. The small matrix is convenient for tracing a format, rather than a demonstration of large-scale performance.

What does sparse mean?

A sparse matrix has enough zeros that representing or processing its nonzero structure can be useful. There is no universal density threshold: index sizes, dimensions, and operations matter. A dense representation stores every position. In COO, CSR, and CSC, positions without a stored entry represent zero. See SciPy’s sparse-array introduction.

The sparsity pattern describes which positions are nonzero. Changing a number at an existing position and inserting a new position are different storage tasks. Also distinguish the mathematical number of nonzeros from the number of stored entries, $s$: software can retain explicit zeros or repeated coordinates. In this baseline both counts equal eight. SciPy’s nnz attribute counts stored entries, including explicit zeros.

COO: attach coordinates to every value

Coordinate format (COO) stores three parallel arrays: row, col, and data. Entry number $k$ contributes data[k] at (row[k], col[k]). Coordinates need not arrive in sorted order. This is the triplet representation described by SciPy’s COO documentation.

Keeping NIST’s file entry order, but using zero-based coordinates, gives:

shape = (5, 5)
row   = [0,   1,     2, 0,     3,    3,     3,  4]
col   = [0,   1,     2, 3,     1,    3,     4,  4]
data  = [1, 10.5, 0.015, 6, 250.5, -280, 33.32, 12]

Read vertically: slot $k=4$ means $a_{3,1}=250.5$. Slot $k=3$ returns to row 0; COO does not require all entries of a row to be adjacent. To reconstruct the matrix, start with zeros and add each triplet’s value at its coordinate.

Always retain shape. The largest stored indices cannot reveal trailing empty rows or columns; an entirely empty coordinate list cannot reveal either dimension.

Assembly and repeated coordinates

Assembly combines contributions into one matrix. Separate pieces of a model can contribute to the same coefficient. Under an additive convention, supplementary invented triplets (1, 2, 3) and (1, 2, 4) combine to $a_{1,2}=7$, rather than one overwriting the other. SciPy’s COO conversion to CSR/CSC sums duplicate coordinates; its documentation explicitly connects this to finite-element assembly. This is a library convention to check when importing other representations, not permission to assume every interface handles duplicates identically. SciPy COO notes

COO is convenient while collecting triplets. Finding one particular row in an unindexed COO list requires examining its row coordinates; recording both coordinates for every slot also costs space. We can derive a different representation by grouping the baseline’s entries by row.

CSR: replace repeated row numbers with boundaries

Compressed Sparse Row (CSR) stores data, column indices indices, and row boundaries indptr. It is also called compressed row storage (CRS). Entries of a row occupy one contiguous segment. See Netlib’s CRS definition.

First order the baseline’s triplets by row, then by column:

The counts are [2, 1, 1, 3, 1]. Start at zero and successively add them: $0,2,3,4,7,8$. These cumulative counts give the row boundaries:

shape   = (5, 5)
indptr  = [0, 2, 3, 4, 7, 8]
indices = [0, 3,    1,     2,     1,    3,     4,  4]
data    = [1, 6, 10.5, 0.015, 250.5, -280, 33.32, 12]

For row $i$, slots run from indptr[i] inclusive to indptr[i+1] exclusive, often written [indptr[i], indptr[i+1]). Each slot’s column is indices[k]. Row 3 occupies slots $4,5,6$: these recover $(3,1,250.5)$, $(3,3,-280)$, and $(3,4,33.32)$. The final pointer 8 marks the end; it is not a ninth stored coefficient. SciPy’s CSR array constructor uses this slice convention.

Empty rows and validity

An empty row still has two boundaries; they are equal. In a hypothetical modification that removes row 3’s entries while keeping the shape, the counts become [2, 1, 1, 0, 1], so indptr = [0, 2, 3, 4, 4, 5]. The repeated 4 encodes an empty segment, not an error.

For the tightly packed, zero-based representation used here, check:

• data and indices both have length $s$; indptr has length $m+1$.
• indptr[0] = 0, indptr[m] = s, and pointers never decrease.
• Every column index is an integer from $0$ through $n-1$.

These conditions follow from partitioning the slots into row segments. Canonical CSR additionally has sorted column indices within each row and no duplicate coordinates, matching SciPy’s canonical-format definition. Grouping by row alone does not guarantee either property.

CSC: compress columns with the same idea

Compressed Sparse Column (CSC) groups by column and stores row indices. It is also called compressed column storage (CCS). The roles of rows and columns are exchanged relative to CSR. Netlib’s CCS definition

For our baseline, the columns contain respectively $1,2,1,2,2$ entries. Cumulative counts give $0,1,3,4,6,8$:

shape   = (5, 5)
indptr  = [0, 1, 3, 4, 6, 8]
indices = [0,    1,     3,     2, 0,    3,     3,  4]
data    = [1, 10.5, 250.5, 0.015, 6, -280, 33.32, 12]

Column $j$ occupies [indptr[j], indptr[j+1]). Thus column 1 uses slots 1 and 2: row indices 1 and 3 carry values 10.5 and 250.5. Notice that indices means rows here and columns in CSR. SciPy’s CSC constructor documents these roles.

The validity checks exchange $m$ and $n$: there are $n+1$ pointers and row indices must lie in $[0,m)$. Canonical CSC sorts row indices within each column and removes duplicate coordinates. An empty column has equal adjacent pointers, just as an empty row does in CSR. SciPy CSC notes

Transposition versus format conversion

The transpose $A^T$ exchanges row and column coordinates. Consequently, CSR arrays for an $m\times n$ matrix can be interpreted as CSC arrays for its $n\times m$ transpose: row groups become column groups, and column indices become row indices. This follows directly from the definitions; Netlib describes the same relationship in its CCS discussion.

That interpretation represents $A^T$, not $A$. To convert CSR of $A$ to CSC of the same $A$, regroup the entries. Our baseline makes the distinction visible: CSR starts data = [1, 6, …], whereas CSC starts [1, 10.5, …].

Explore the published baseline

Choose a format and a row or column to trace. The small k labels identify positions in data, not coefficient magnitudes. COO keeps the published entry order; CSR and CSC use sorted inner indices. The removal scenarios are hypothetical experiments, retain the $5\times5$ shape, and omit zeros rather than explicitly storing them. Restore the baseline to return to the published entries. The arrays above remain available without JavaScript.

Computation: the same product, different traversal

The following products and operation counts are modern derivations from the array definitions, not timings from the sources. For a supplementary input $x=(1,1,1,1,1)^T$, each output is a row sum:

CSR gathers one output at a time

In $y_i=\sum_j a_{ij}x_j$, row $i$ determines one output. Its segment tells us exactly which terms to visit:

for i = 0, ..., m-1:
    y[i] = 0
    for k = indptr[i], ..., indptr[i+1]-1:
        y[i] += data[k] * x[indices[k]]

For row 3 this evaluates $250.5x_1-280x_3+33.32x_4$. An empty segment leaves its output zero. Each stored slot is visited once, with $m$ row initializations: $O(m+s)$ work. Here $O(\cdot)$ describes how the work grows, ignoring fixed constants. Unsorted indices and additive duplicates still give the algebraic product with this loop; they need not be canonical first.

CSC scatters each column’s contribution

The equivalent identity $Ax=\sum_j x_j a_j$, where $a_j$ denotes column $j$, says to scale each column by $x_j$ and add it to the output:

initialize y[0:m] to zero
for j = 0, ..., n-1:
    for k = indptr[j], ..., indptr[j+1]-1:
        y[indices[k]] += data[k] * x[j]

Column 1 adds $10.5x_1$ to $y_1$ and $250.5x_1$ to $y_3$. This visits $s$ slots, $n$ columns, and initializes $m$ outputs: $O(m+n+s)$ work for this particular loop. For $A^Tx$, the CSC column segments instead compute individual output dot products. An additive COO loop, y[row[k]] += data[k] * x[col[k]], also computes $Ax$ in $O(m+s)$ work after initialization. Formats affect access and organization even when asymptotic counts look similar.

Parallel CSR row computations write distinct output entries. Parallel CSC columns may update the same output and need coordinated accumulation. These observations follow from the loops; they do not establish a universal speed ranking. Machine behavior, row lengths, and implementation still matter.

Conversion: count, reserve, and place

To group $s$ valid COO triplets into CSR without sorting within rows:

1. Count triplets per row in an array of length $m$.
2. Form cumulative counts to create indptr.
3. Copy each row’s starting boundary to a temporary next array.
4. For each (i, j, v), put j and v at slot next[i], then advance next[i].

Every input slot is placed once. The counts and cumulative sums take $O(m+s)$ work in total, with $O(m)$ temporary bookkeeping beyond the output. This is a derivation for the stated algorithm, not a promise about every library conversion. Sorting and combining duplicates require additional steps. Applying comparison sorting to $s$ triplets by (row, column) first gives an $O(s\log s+m)$ route to canonical row ordering; sum equal coordinates, and optionally remove resulting zeros. CSC uses column counts and the analogous $n$-dependent bounds. Preserve the original shape throughout.

Storage cost and choosing a format

Let each value occupy $b_v$ bytes and each index or pointer occupy $b_i$ bytes. Counting array entries gives these payload estimates, excluding object headers, alignment, spare capacity, shape metadata, and conversion workspaces:

For the published baseline, assume eight-byte values and four-byte integers. Dense uses $25\cdot8=200$ bytes, COO $8\cdot16=128$, and either compressed format $8\cdot12+6\cdot4=120$. These are illustrative arithmetic, not measured allocations.

CSR saves payload space relative to this dense model exactly when

For a supplementary hypothetical dense $2\times2$ matrix, $s=4$ with the same widths: dense takes 32 bytes and CSR takes $4\cdot12+3\cdot4=60$. Sparse storage can cost more. Likewise, a very tall matrix with few entries can pay heavily for CSR’s row pointers; COO avoids that dimension-sized array, and CSC may have fewer pointers.

The row/column access tradeoffs and cost of changing compressed structure are documented in SciPy CSR notes and CSC notes. DOK stores a dictionary keyed by coordinates; LIL stores lists per row. See SciPy’s format overview for these and block or diagonal formats. Choose a solver’s expected format when handing data to it, and account for conversion if operations repeat.

Misconceptions and edge cases

• “Compressed means values lose precision.” CSR and CSC compress position information by grouping entries. The format itself does not quantize values; the chosen numeric type controls their precision.
• “Every stored value is nonzero.” An explicit zero occupies a slot. Invented duplicate contributions $+3$ and $-3$ at the same coordinate can also cancel. Summing duplicates and removing stored zeros are separate operations. SciPy’s discussion of implicit zeros and duplicates
• “Every operation stays sparse.” Adding a nonzero constant to every position makes formerly zero positions nonzero. Factorization can also introduce fill-in, new nonzero positions during elimination; see the coupling example in Schur complement.
• “An empty matrix needs no metadata.” For an all-zero $m\times n$ matrix, COO has empty triplet arrays; CSR has $m+1$ zero pointers and CSC $n+1$. All still retain the shape. These follow from the segment definitions.

Practice problems

The following exercises are original teaching exercises. Problems 1 and 4 use the published baseline; other numerical data are supplementary invented examples.

1. Hand calculation: recover a column

Use the baseline CSC arrays to recover column 3, and compute that column’s contribution when $x_3=2$.

2. Hand calculation: preserve empty dimensions

Encode $B$ in CSR and CSC:

3. Strategy: choose a construction workflow

A program collects contributions in arbitrary order, sometimes repeating coordinates, then repeatedly evaluates $Ax$. What storage workflow would you try, and what would you measure?

4. Logic: prove the transpose interpretation

Explain why baseline CSR arrays represent $A^T$ when reinterpreted as CSC. Why do those arrays fail to represent the same $A$ in CSC?

5. Logic: distinguish storage from coefficients

An additive COO representation of shape $(2,3)$ has triplets (0,1,3), (0,1,-3), and (1,2,0). How many stored entries and mathematical nonzeros are there before and after summing duplicates and removing explicit zeros?

6. Application: estimate a memory budget

A $1000\times1000$ matrix has 5000 stored entries, no duplicates, and no explicit zeros. Estimate dense, COO, and CSR payload sizes with eight-byte values and four-byte indices. Does this establish which representation computes fastest?

Further learning

• NIST: Matrix Market file formats: connect coordinate triplets to an exchange file, and compare its one-based indices with the arrays here.
• Netlib: sparse matrix storage survey: explore how block, diagonal, and other structure can motivate additional formats.
• SciPy: sparse arrays tutorial: practice construction, conversion, and the distinction between stored and implicit zeros.

Related pages

• Matrix factorization and decomposition
• LU decomposition
• Cholesky decomposition
• Schur complement
• Linear programming

References

1. NIST. Matrix Market. The repository’s purpose: test data for comparing numerical linear algebra algorithms.
2. NIST. Matrix Market: Text File Formats, “Matrix Market Exchange Formats,” real general $5\times5$ example. Source of the eight baseline coefficients and the file’s one-based convention.
3. Barrett, R., et al. Templates for the Solution of Linear Systems: Building Blocks for Iterative Methods. Netlib online edition: storage survey, compressed row storage, and compressed column storage.
4. SciPy documentation. coo_array. Coordinate storage, duplicate handling, and stored-entry counts.
5. SciPy documentation. csr_array. Row segments, canonical form, and access tradeoffs.
6. SciPy documentation. csc_array. Column segments, canonical form, and access tradeoffs.
7. SciPy documentation. Sparse Arrays tutorial. Format choices and explicit versus implicit zeros.