# A composable abstraction of hierarchical methods for matrix-vector product acceleration

Antoine Gicquel

► To cite this version:

Antoine Gicquel. A composable abstraction of hierarchical methods for matrix-vector product acceleration. COMPAS 2025 - Conférence francophone d'informatique en Parallélisme, Architecture et Système, Inria Bordeaux, Jun 2025, Bordeaux, France. hal-05185249

HAL Id: hal-05185249

https://inria.hal.science/hal-05185249v1

Submitted on 25 Jul 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# A composable abstraction of hierarchical methods for matrix-vector product acceleration

Antoine Gicquel

Inria Bordeaux (Concace team)

COMPAS, June 2025

<figure><img src="image_2.png" /></figure>
### Context

We focus on the **acceleration** of the **matrix-vector product**

$$y = Ax$$

where $A \in \mathbb{R}^{N \times N}$ and $x, y \in \mathbb{R}^N$

### Cost of this $N \times N$ matrix-vector product

* If $A$ is a **dense matrix**: $O(N^2)$

* If $A$ is a **(hierarchically) rank-structured matrix**: $O(N \log N)$ or $O(N)$ with accuracy $\epsilon$ via **hierarchical methods**

### Applications (rank-structured matrix)

$N$-body problem, Boundary Element Method (BEM), Volume Integral Methods (VIM),...
# Introduction

## Model problem

### Particle interactions (N-body problem)

<figure><img src="image_3.png" /></figure>

$$ \varphi(s_j) = - \sum_{i=1, i \neq j}^{N} \log |s_i - s_j| q_i \quad j \le N $$

$$ \mathbf{y} = \mathbf{A}\mathbf{x} \qquad \mathbf{A}_{ij} = -\delta_{ij}\log|s_i - s_j| $$

<figure><img src="image_4.png" /></figure>

Matrix A

### Low-rank approximation

* Algebraic approximation: $\mathbf{A}_{IJ} \approx \hat{\mathbf{U}}\hat{\mathbf{V}}^*$

* Analytical approximation (if $\mathbf{A}_{IJ} = k(s_i, s_j)$):

$$ k(x, y) \approx \sum_{\lambda} \varphi_{\lambda}(x) \gamma_{\lambda}(y) $$

<figure><img src="image_5.png" /></figure>
# Introduction

## Model problem

<figure><img src="image_6.png" /></figure>

<figure><img src="image_7.png" /></figure>

Flat partitioning
⇒ single-level tiling

<figure><img src="image_8.png" /></figure>

<figure><img src="image_9.png" /></figure>

Hierarchical partitioning
⇒ multi-level tiling
# Introduction

## State-of-the-art (Hierarchical methods)

<figure><img src="image_10.png" /></figure>

# Contributions

## Composable abstraction + generic matrix-vector product algorithm
Abstraction of a hierarchical method

<figure><img src="image_11.png" /></figure>

* **Reordering and partitioning:**

* **Subdivide matrix** $A$ into tiles $\{A_{ij}^{(\ell)}\}$ at different levels $\ell$

* **Admissibility criterion:**

* **Determines which tiles** $A_{ij}^{(\ell)}$ **are approximated (or admissible)**

* **Approximation strategy:**

* **Approximates the admissible tiles** $A_{ij}^{(\ell)} \approx \hat{V}_I^{(\ell)} \hat{A}_{IJ}^{(\ell)} \hat{W}_J^{(\ell)}$*

* **Basis mutualization:**

* **Factorizes the approximation bases** $\hat{V}_I^{(\ell)}$ and $\hat{W}_J^{(\ell)}$
Generic algorithm for the hierarchical matrix-vector product

<figure><img src="image_12.png" /></figure>
## Existing implementations

* **Analytical approximation (FMM):** PVFMM [Malhotra and Biros 2015], exaFMM [Barba and Yokota 2011], ...

* **Algebraic approximation:** hmat [Lizé 2014], HLIBpro [Kriemann 2012], Dense_HODLR, ...

## Our implementations

* **Software:** `composyx` (composability-oriented C++ library)

* **Composability at the algorithmic level (modular, maintainable, ...)**

* **Include various components of the abstraction:**
  * Different approximations, admissibility criteria, partitioning...

* **Flexibility according to the hardware/architecture:**
  * e.g., analytical approx. (high Flop/Byte), algebraic approx. (more storage), ...

* **Flexibility according to the application:**
  * e.g., N-body problem (analytical approx.), iterative solver (algebraic approx.), ...
# Numerical results

## Sequential implementation

* Particle interactions (2D): $A_{ij} = -\delta_{ij}\log \|s_i - s_j\|$

* Reference software: ScalFMM <fnref n="1" /> (FMM) and hmat-oss (H-matrix)

* Numerical validations

* Preliminary benchmarks <fnref n="2" />

<figure><img src="image_13.png" /></figure>

<figure><img src="image_14.png" /></figure>

<footnote n="1">¹For a fair comparison, optimizations based on kernel symmetries for the P2P and M2L operators were disabled.</footnote>

<footnote n="2">²plafrim/bora: 2x 18-core Cascade Lake Intel Xeon Skylake Gold 6240 @ 2.6 GHz</footnote>
# Parallelization

## Parallelization of hierarchical methods

* Extensive work (various approaches on different architectures [Lashuk et al. 2012], [Agullo et al. 2016], [Keyes et al. 2020])

* Little unification [Yokota et al. 2014]

## Question

We focused on composability at the **algorithmic level**. What about **parallelization** ?

* **Early results**

► Leverage the parallelization of
**distributed matrix multiplication**
directly in `composyx`

* **Perspectives**

► High-level abstraction: **Map Reduce**
model (cf: H. Dodelin's talk)

<figure><img src="image_15.png" /></figure>
# Conclusion

## Summary

* A **composable abstraction** of hierarchical methods

* A **generic algorithm** for the matrix-vector product

* Full working implementation in `composyx`

* **Composability** at the algorithmic level

* Single-level parallelization (match existing **distributed** GEMM/GEMV)

## Perspectives

* Full parallelization via **Map Reduce** model

* More applications to real cases
<figure><img src="image_16.png" /></figure>

Agullo, Emmanuel, Berenger Bramas, Olivier Coulaud, Eric Darve, Matthias Messner, and Toru Takahashi (2016). "Task-based FMM for heterogeneous architectures". In: *Concurrency and Computation: Practice and Experience* 28.9, pp. 2608-2629.

<figure><img src="image_17.png" /></figure>

Börm, Steffen, Lars Grasedyck, and Wolfgang Hackbusch (2003). "Introduction to hierarchical matrices with applications". In: *Engineering analysis with boundary elements* 27.5, pp. 405-422.

<figure><img src="image_18.png" /></figure>

Greengard, Leslie and Vladimir Rokhlin (1987). "A fast algorithm for particle simulations". In: *Journal of computational physics* 73.2, pp. 325-348.

<figure><img src="image_19.png" /></figure>

Keyes, David E, Hatem Ltaief, and George Turkiyyah (2020). "Hierarchical algorithms on hierarchical architectures". In: *Philosophical Transactions of the Royal Society A* 378.2166, p. 20190055.

<figure><img src="image_20.png" /></figure>

Lashuk, Ilya, Aparna Chandramowlishwaran, Harper Langston, Tuan-Anh Nguyen, Rahul Sampath, Aashay Shringarpure, Richard Vuduc, Lexing Ying, Denis Zorin, and George Biros (2012). "A massively parallel adaptive fast multipole method on heterogeneous architectures". In: *Commun. ACM* 55.5, pp. 101-109. DOI: 10.1145/2160718.2160740.