🔢 mathematics

Algorithm librla: A library of randomized linear algebra routines

The paper introduces \texttt{librla}, a stable and efficient randomized linear algebra library implemented in MATLAB, Python, and Julia for intermediate-sized matrices, offering flexible factorization options (QR, SVD, interpolative decomposition) based on fixed rank or tolerance and supporting both matrix and matrix-vector multiplication access.

Original authors: Adrianna Gillman, Zydrunas Gimbutas

Published 2026-08-06

📖 4 min read🧠 Deep dive

CC BY 4.0

Original authors: Adrianna Gillman, Zydrunas Gimbutas

Original paper licensed under CC BY 4.0 (http://creativecommons.org/licenses/by/4.0/). ✨ This is an AI-generated explanation of the paper below. It is not written or endorsed by the authors. For technical accuracy, refer to the original paper. Read full disclaimer

Imagine you are trying to organize a massive, chaotic library where every book is a giant spreadsheet of numbers. Some of these spreadsheets are so huge they would fill a stadium, but most are just big enough to fit in a very large backpack. In the world of computer science and math, these are called matrices. Often, these giant spreadsheets are full of "noise" or redundant information, meaning they actually hold a much simpler, smaller story underneath. The challenge is finding that hidden story without getting lost in the millions of numbers. This is where "randomized linear algebra" comes in. Think of it like trying to guess the flavor of a giant pot of soup. Instead of tasting every single spoonful (which takes forever), you take a few random sips. If you pick your sips wisely, you can figure out the main ingredients and the overall taste with just a handful of samples. This technique is crucial for everything from training artificial intelligence to compressing images and solving complex physics problems, because it lets computers work with huge data sets much faster.

Now, meet librla, a new software toolkit designed by Adrianna Gillman and Zydrunas Gimbutas to make this "soup-tasting" easier, faster, and more reliable. Before this library, researchers often had to juggle different tools that were either too slow, prone to crashing, or only worked in one specific computer language. The authors built librla to be a universal translator that works smoothly in three popular languages: Python, MATLAB, and Julia. It is specifically designed for "medium-sized" matrices—those with dimensions up to roughly 10,000—which are the sweet spot for many real-world applications like reducing complex models or compressing data.

The library offers three main ways to simplify these data sets, which the authors call factorizations. You can think of these as different ways to summarize a long novel. The first is QR, which is like creating a clean, organized outline of the story. The second is SVD (Singular Value Decomposition), which breaks the story down into its most important themes and characters. The third is Interpolative Decomposition, which is like picking a few key sentences from the book that can be used to reconstruct the rest of the text. What makes librla special is that it lets you choose how you want to summarize: you can tell it, "Give me a summary with exactly 15 key points," or "Give me a summary that is accurate enough for my needs," letting the software decide how many points it needs to get the job done.

The paper shows that librla is not just a new way of doing things, but a faster and more stable one. In tests, the library's Python version was compared to existing tools in PyTorch and SciPy. The results showed that librla was often significantly faster—sometimes up to 34 times faster than SciPy for certain tasks—while producing results with similar accuracy. It also handles "matrix-free" problems, meaning it can work even if you don't have the full spreadsheet in memory, but only have a way to calculate the result of multiplying the data by a number. This is like being able to taste the soup without ever seeing the pot.

The authors also explored how to make these summaries even better. They found that while taking a few extra random samples helps a little, using a technique called "power iteration" (which is like re-reading the most important pages of the book to make sure you understood them) makes a huge difference. Combining extra samples with power iteration gave the best results, allowing the software to capture the true "spectrum" or hidden structure of the data much more accurately.

Ultimately, this paper presents librla as a robust, efficient, and user-friendly tool that fills a gap in the software world. It doesn't claim to solve every math problem in existence, but it provides a stable and fast foundation for handling the medium-sized data challenges that appear frequently in science and engineering. By making these powerful randomized techniques accessible in multiple languages and offering flexible options for accuracy and speed, librla helps researchers and developers turn their messy, giant data sets into clear, manageable insights without the headache of unstable code.

Technical Summary of librla: A Library of Randomized Linear Algebra Routines

Problem Statement

While randomized linear algebra has become a vital tool for developing algorithms in neural networks, fast direct solvers, reduced order modeling, and data compression, a significant gap exists in available software. Current factorization libraries often suffer from stability issues, manifesting as invalid factorizations or failures to return results. Furthermore, dependencies on low-level legacy libraries (e.g., Fortran) hinder the portability and usability of new high-level implementations. There is a lack of a unified, stable, and efficient library that supports multiple high-level languages (MATLAB, Python, Julia) and offers flexible factorization options (fixed rank vs. tolerance) for intermediate-sized matrices (dimensions up to roughly 10,000).

Methodology

The librla library is built upon the orth_sketch algorithm, which creates an orthogonal basis for the range of a linear operator via randomized sampling. The core methodology involves:

Randomized Range Finding: The algorithm samples the range of a matrix $A$ (or operator) by applying it to a random matrix $\Omega$ with entries sampled uniformly from $[-1, 1]$ .
Orthogonalization: An orthogonal basis $Q$ is formed by computing a pivoted QR factorization of $A\Omega$ .
Enhancements:
- Power Iteration: Optional power iterations ( $q$ ) are applied to sharpen the spectral decay of the sketched operator, improving accuracy for matrices with slow singular value decay.
- Oversampling: Extra samples ( $p$ ) are added to the requested rank to improve expected error and reduce the probability of large deviations from the true spectrum.
Factorization Modes:
- Fixed Rank Mode: The user specifies a target rank $k$ . The sample matrix $\Omega$ has $b = k + p$ columns.
- Tolerance Mode: The user specifies a desired accuracy (tolerance). The algorithm iteratively increases the block size of $\Omega$ (starting at 42 columns, growing geometrically by a factor of 4) until the ratio of the smallest to largest diagonal entries of the resulting $R$ matrix falls below the tolerance. If the tolerance cannot be met before the block size reaches the matrix dimension, the routine falls back to a deterministic full pivoted QR factorization.
Matrix-Free Capability: The library supports "matrix-free" operations where the user provides a subroutine for matrix-vector multiplication (and its conjugate transpose) rather than the explicit matrix.
Implementation: The library is written natively in Python, Julia, and MATLAB. It is designed to exploit Level 3 BLAS operators to ensure efficient execution via low-level optimized code.

Key Contributions

The paper presents librla as the first software package that is simultaneously stable, efficient, and available in MATLAB, Python, and Julia. Its specific contributions include:

Unified Factorization Options: It provides three distinct factorization routines:
- qr_sketch: Randomized QR factorization.
- svd_sketch: Randomized Singular Value Decomposition (SVD).
- id_sketch: Interpolative Decomposition (skeletonization).
Flexible Stopping Criteria: Users can choose between a fixed rank or a desired tolerance, a feature not uniformly available in existing high-level packages.
Robustness and Stability: The library addresses stability issues common in other libraries by ensuring valid factorizations and providing fallback mechanisms (e.g., to deterministic QR) when tolerance criteria cannot be met.
Interpolative Decomposition Variants: The id_sketch routine offers three methods for constructing the interpolation matrix $T$ : a fast triangular solve, a minimum-norm solution via SVD (to handle ill-conditioning), and a least-squares approach using original columns for higher accuracy.

Results

Numerical experiments were conducted on a MacBook Pro (M4 Max, 64 GB RAM) comparing librla against existing Python libraries (PyTorch's svd_lowrank and SciPy's interp_decomp).

Performance Comparison:
- vs. PyTorch: For fixed-rank SVDs on Hilbert matrices, decaying spectrum matrices, and Gaussian mixture models, librla demonstrated comparable timing to PyTorch (e.g., 0.0018s vs 0.002s for a 2000x1000 Hilbert matrix).
- vs. SciPy: librla showed significant speedups over SciPy's interpolative decomposition. For a 2000x1000 Hilbert matrix, librla was 14.6 times faster; for a 4000x2000 matrix, it was 33.9 times faster. Similar speedups (12.8x to 20x) were observed for tolerance-based approximations.
Accuracy: The relative errors for librla were comparable to or better than the reference libraries.
Application Examples:
- Climate Data: Applied to NOAA ERSST v5 data, the results confirmed that combining oversampling and power iterations yields the best approximation of singular values.
- Image Compression: Applied to image data, librla successfully reconstructed images using both SVD and interpolative decomposition. The SVD approach produced images with less "smearing" compared to the interpolative decomposition at the same rank.

Significance and Claims

The authors claim that librla fills a critical gap in the software domain by providing a robust, multi-language solution for randomized linear algebra. The library is significant because:

Portability: It eliminates the dependency on legacy low-level libraries (like Fortran) that often hinder the adoption of randomized methods in high-level languages.
Functionality: It offers a broader range of options (QR, SVD, ID) and control parameters (rank vs. tolerance, power iterations, oversampling) than current available packages.
Efficiency: By leveraging Level 3 BLAS, the library achieves performance comparable to optimized C/Fortran implementations while maintaining the readability and accessibility of high-level code.
Stability: It provides a stable alternative where existing libraries may fail or produce invalid factorizations, particularly in tolerance-based modes.

The paper concludes that librla is suitable for intermediate-sized matrices arising in applications such as reduced order modeling, fast direct solvers, and data compression, but is not intended for matrices that are too large to fit in memory.