Header-only · C++20 · MPI + ScaLAPACK, without the boilerplate.
Documentation · Quick start · Solvers · Standalone solvers
ezp is a lightweight, header-only C++ wrapper around a curated set of distributed linear-system solvers.
It hides the BLACS, ScaLAPACK, MPI and low-level C/Fortran minutiae behind an interface that feels just like solving a system sequentially.
solver.solve(A, B); // that's the whole ideaYou may be interested if you...
- 🧩 want to solve linear systems on distributed clusters without memorising
ScaLAPACKcalling conventions. - 🛠️ want to avoid the hassles of learning how to call various low-level
C/Fortransolvers and managing their buffers by hand. - 🌐 want an object-oriented,
MPI-aware API that resembles theLAPACKexperience you already know. - 🚀 want to separate the solver from the application through standalone solver executables.
| Feature | |
|---|---|
| ✨ |
Ergonomic interface — prepare |
| 📦 |
Drop-in, header-only library — no build step, no ABI to manage; just #include
|
| 🧮 |
Full precision coverage — all four DSZC precisions, with 32-bit and 64-bit indexing |
| 🖥️ |
Standalone solver binaries — invoke the solver from any MPI-capable caller, including Python via mpi4py
|
| 🧠 | Build once, solve many — factorization is retained and reused across right-hand sides |
| ✅ | Validated with randomly generated inputs — an extensive test suite backs every solver |
The library implements all precisions with both 32-bit and 64-bit integer indexing, and supports both row- and column-major process grids.
| availability | type of matrix | operation | solver | package |
|---|---|---|---|---|
| ✅ | general (partial pivoting) | simple | PxGESV | ScaLAPACK |
| ✅ | general (partial pivoting) | expert | PxGESVX | ScaLAPACK |
| ✅ | symmetric/Hermitian positive definite | simple | PxPOSV | ScaLAPACK |
| ✅ | symmetric/Hermitian positive definite | expert | PxPOSVX | ScaLAPACK |
| ✅ | general band (partial pivoting) | simple | PxGBSV | ScaLAPACK |
| ✅ | general band (no pivoting) | simple | PxDBSV | ScaLAPACK |
| ✅ | symmetric/Hermitian positive definite band | simple | PxPBSV | ScaLAPACK |
| availability | type of matrix | operation | solver | package |
|---|---|---|---|---|
| ✅ | sparse (one- or zero-indexing CSR format) | direct | PARDISO | MKL |
| ✅ | sparse (one-indexing COO format) | direct | MUMPS | MUMPS |
| ✅ | sparse (zero-indexing CSR format) | iterative | Lis | Lis |
Expert drivers additionally provide equilibration, iterative refinement, condition estimation and determinant computation where applicable.
ezp requires a C++20-compatible compiler and the following drivers:
- an implementation of
LAPACK/BLAS, such asOpenBLAS,MKL,AOCL, ... - an implementation of
ScaLAPACK; - an implementation of
MPI, such asOpenMPI,MPICH,Intel® MPI, ...
The easiest setup is the Intel® oneAPI toolkit, which supplies MPI, oneMKL (with a complete ScaLAPACK) in one shot:
git clone --recurse-submodules --depth 1 https://github.com/TLCFEM/ezp.git
mkdir ezp/build && cd ezp/build
. /opt/intel/oneapi/setvars.sh && cmake -DEZP_TEST=ON .. && cmake --build . --config ReleaseUsing system libraries is also possible — see the compilation guide for openSUSE, Fedora and other distros.
It is assumed that the root node (rank 0) prepares the left-hand side
All BLACS and ScaLAPACK details are hidden: prepare the matrices (on the root node), wrap them, and call solve.
The result highly resembles how one would solve a linear system sequentially.
#include <ezp/pgesv.hpp>
#include <iomanip>
#include <iostream>
using namespace ezp;
int main() {
// get the current blacs environment
const auto env = get_env<int>();
constexpr auto N = 6, NRHS = 2;
// storage for the matrices A and B
std::vector<double> A, B;
// helper functor to convert 2D indices to 1D indices
// it's likely the matrices are provided by some other subsystem
const auto IDX = par_dgesv<int>::indexer{N};
if(0 == env.rank()) {
// the matrices are only initialized on the root process
A.resize(N * N, 0.);
B.resize(N * NRHS, 1.);
for(auto I = 0; I < N; ++I) A[IDX(I, I)] = static_cast<double>(I + 1);
}
// create a parallel solver
// it takes the number of rows and columns of the process grid as arguments,
// or let the library automatically determine them
auto solver = par_dgesv<int>();
// wrap the data in full_mat objects:
// number of rows, number of columns, and a pointer to the data
// on non-root processes the pointer is nullptr, as the vectors are empty
const auto info = solver.solve({N, N, A.data()}, {N, NRHS, B.data()});
// reuse the stored factorization to compute the determinant
const auto det = solver.det({N, N, A.data()});
if(0 == env.rank() && 0 == info) {
std::cout << std::setprecision(10) << "Determinant: " << det << '\n';
std::cout << "Solution:\n";
for(const auto i : B) std::cout << i << '\n';
}
return info;
}More runnable examples live in examples/, covering every solver and input format.
MPI can dynamically spawn worker processes via MPI_Comm_spawn.
ezp exploits this to ship its solvers as standalone executables, cleanly separating the main application from the solver.
This is particularly handy when the main application runs a different parallelism model, has latency constraints, or simply wants a clear boundary around the solver.
Because the boundary is a plain MPI interface, callers can be written in any MPI-compatible language — the bundled runner.py demonstrates doing it from Python with mpi4py.
Enable them at configure time with -DEZP_STANDALONE=ON; see the standalone guide.
Full guides, per-solver references and the API docs are published at https://tlcfem.github.io/ezp/:
- Motivation and Compilation
- Usage, Custom matrices, MPI management, Standalone solvers
- Per-solver reference for dense and sparse solvers
Important
ezp is distributed under the terms of the GNU General Public License v3.0.
See LICENSE for details.
Copyright © 2025–2026 Theodore Chang