# Why does sparse linear algebra have a low arithmetic intensity?

I often see the terms "low arithmetic intensity" and "memory-bound" associated with sparse matrix operations. However, my intuition is that a sparse matrix operation should be less memory-bound, if anything. After all, for a highly sparse matrix, you only have to transfer a fraction of the data required for a dense matrix operation (i.e. only the non-zero, actually meaningful entries of the matrix have to be stored and transferred).

So my question is: why are sparse matrix operations usually considered memory-bound?

I suppose there must be some level of "denseness" at which a sparse matrix computation becomes compute-bound.

• If one of the operands is dense and requires much memory, they often can't use CPU caching of memory very well. This can take latency from getting data from like a few clock cycles to 10s or even 100s of them. – mathreadler May 21 '18 at 18:54

BLAS1-operations, BLAS2-operations, and sparse-operations share the same curse of low arithmetic intensity, that they perform $O(1)$ flops for each memory read (contrast this to a BLAS3-operation like gemm, which performs $O(N^3)$ flops over $O(N^2)$ reads and only becomes more and more arithmetic-intensive/compute-bound for large $N$).

However, sparse matvecs have even further complications, because of indirect addressing. This adds (at the very least) additional memory pressure to load the integer indices. Even worse, it results in much poorer cache utilization, since a cache typically fetches a large run of data nearby any single read, but a sparse matvec will generally only use just the one value. (BLAS1/BLAS2-operations are low arithmetic intensity, but at least cache-friendly).

There are sparse matrix formats that attempt to mitigate some of this cache-unfriendliness, typically they are like blocked variants of the traditional formats (CSR/CSC/etc), where each "nonzero" of the matrix is instead a small dense matrix. This lets you cut out some of the indirect addressing and improve the locality of access. They are a good fit for high-order FEM's, where many dof's are often collocated on the same geometrical entity and share the same pattern. But they are not a good fit for low-order methods (there's just not enough local similarity in such a matrix to effectively "build-up" a large block).

All this said, the upshot to "low arithmetic intensity" is "low arithmetic" .. sparsity-exploiting methods are generally fast even though they're not particularly efficient at the hardware level. They're just innately doing "less work"... not necessarily a bad thing from a user perspective.

• Regarding indirect addressing, I made some benchmarks once in which sparse matvec with the CSR format was appreciably faster when it had been reordered using the Cuthill-McKee algorithm. Same number of FLOPs, but reducing the bandwidth of the matrix made for greater locality of reference and thus better cache utilization. – Daniel Shapero May 21 '18 at 17:43
• Caches in this context probably work more against performance than for it. Caching helps performance when there is re-use of data, but in a matrix-vector product there is virtually no data re-use to speak of. When we request data in our code, the CPU gives us an entire cache line whether we want it or not. It's up to us to make full use of that cache line. For codes with a lot of random access that line gets wasted and eventually evicted from cache. If we can regularize access though we use that cache line - but our overall performance is still bound by the system's memory bandwidth. – Reid.Atcheson May 21 '18 at 21:51

This really depends on the operations you are including in your question. If you took the sparse equivalent of any level 1 BLAS or level 2 BLAS algorithm, then yes they are memory bound (not compute bound), but that is also true of level 1 BLAS and level 2 BLAS for dense matrices. To answer this we need to understand what makes an algorithm compute bound.

Compute bound algorithms are those for which their performance relies largely on the capability of the underlying system to perform arithmetic, often measured in FLOPS.

For an algorithm to have this property it must do significantly more arithmetic operations than memory operations, because memory operations have very high latency. In the context of linear algebra algorithms this means there must be significant data re-use in the course of the algorithm. Data re-use means we can avoid resorting to main memory. Level 3 BLAS operations often have this property.

For example to do a sparse matrix-vector product you scan through its list of nonzeros only once. There is no way to cache any values in any way that could help performance. The same is true of dense matrix-vector product for that matter. These are not examples of compute bound algorithms.

Unfortunately most compute architectures today are very good at arithmetic and only decent at memory performance. Thus A lot of sparse linear algebra algorithms are specifically designed to reduce to dense linear algebra algorithms from level 3 BLAS at its most demanding points. Multifrontal factorization is an example of an algorithm like this.

So the overall answer to you is that: Sparse linear algebra is memory bound for many operations, but so is dense linear algebra. We can choose/design our algorithms to match the capabilities of the system at hand however, so this is not necessarily a bad thing.