MatX: NumPy-style tensors with C++ control on the GPU

The expression model, written in C++

MatX is a C++20 library for general numerical computing on NVIDIA GPUs and CPUs, and the README states the idea in one line: NumPy-like tensor expressions, C++ control, NVIDIA GPU speed. In its own words, it is NumPy-style array programming with native C++ control and speed, so you create arbitrary-rank tensors, compose lazy expressions, and choose an optimized GPU or multithreaded CPU executor when the work runs.

The example in the README is ordinary numerical code:

#include <matx.h>

constexpr matx::index_t batch = 1024;
constexpr matx::index_t inputs = 256;
constexpr matx::index_t outputs = 128;

auto x = matx::make_tensor<float>({batch, inputs});
auto weights = matx::make_tensor<float>({inputs, outputs});
auto bias = matx::make_tensor<float>({outputs});
auto activations = matx::make_tensor<float>({batch, outputs});

// bias broadcasts naturally across the batch dimension.
auto layer = matx::tanh(matx::matmul(x, weights) + bias);

(activations = layer).run(matx::cudaExecutor{});

A batched matrix product, a broadcast add and an elementwise activation, and the bias is a vector that broadcasts across the batch dimension without any explicit reshaping. The last line is where the execution model appears: the expression is assigned to a tensor and then run with an executor, so nothing computes until the run call.

The mapping to NumPy is laid out in a table in the README, and it is a fair guide to the surface. Arbitrary-rank tensors, non-owning views, slicing, reshaping and permutation stand in for ndarray and views. zeros, ones, eye, linspace, random generators and custom generators stand in for array creation. Lazy elementwise operators, automatic broadcasting and batched transforms stand in for vectorized expressions. And reductions, BLAS, sparse operations, FFTs, decompositions and solvers cover the numerical routines.

The difference from NumPy is not the vocabulary, it is what happens at the run. The same expression can be executed by a CUDA executor, by the JIT executor, or by a multithreaded host executor, and the choice is the only thing that changes.

The executor is the switch, not the algorithm

The feature list puts it as fuse compatible numerical pipelines, eliminate launches and intermediate memory traffic by changing the executor, not the algorithm.

That is the design claim, and it is a strong one. In a typical tensor framework the pipeline is a sequence of kernel launches, with a temporary tensor allocated for each intermediate result, so a chain of small operations spends its time in launch overhead and memory traffic rather than arithmetic. MatX's claim is that the expression is a single description, and the backend decides how to execute it, so the intermediates can disappear.

The three executors mentioned in the README are a CUDA executor, a JIT executor, and a host executor, and the NumPy-like table makes the same point in one row: run the same expression with a CUDA, JIT, or host executor.

The host path is the part that makes this more than a GPU library. The README says you can target multithreaded host execution with NVPL, FFTW, OpenBLAS or BLIS, which are the CPU equivalents of the optimised GPU libraries. So the same numerical code has a path that runs on a machine without an NVIDIA GPU, and the executor choice is a deployment decision rather than a rewrite.

The FFT example is described in the README as a useful fusion example because it crosses an optimized-library boundary and makes the cost of intermediate tensors easy to measure:

#include <matx.h>

constexpr matx::index_t batches = 256;
constexpr matx::index_t samples = 4096;

auto signal      = matx::make_tensor<matx::fcomplex>({batches, samples});
auto calibration = matx::make_tensor<float>({batches, samples});
auto spectrum_db = matx::make_tensor<float>({batches, samples});

auto normalized_signal = signal / (matx::abs(signal) + 1.0e-6f);
auto pipeline = 20.0f * matx::log10(
  matx::abs(matx::fft(normalized_signal) * calibration) + 1.0e-6f);

(spectrum_db = pipeline).run(matx::CUDAJITExecutor{});

Normalise every complex sample, transform every batch, apply the frequency domain calibration, take the magnitude, convert to decibels. With -DMATX_EN_MATHDX=ON and a supported FFT shape, the README states that MatX lowers the compatible pipeline into one generated kernel, with no normalized-input or FFT temporary and no separate pre- or post-processing launch.

The benefit it claims is for memory-bound workloads, where eliminating launches and intermediate reads and writes lets a simple expression reach speed-of-light performance without hand-writing CUDA kernel plumbing.

JIT fusion across four vendor libraries

The JIT path is the part of the project with the most moving pieces, and the README documents exactly what it can join into a single generated kernel.

The mechanism is that CUDAJITExecutor can fuse supported transforms and their surrounding operators into a single runtime-generated kernel, and with -DMATX_EN_MATHDX=ON the compatible paths are listed.

cuBLASDx contributes compatible matrix multiplication and GEMM paths. cuSolverDx contributes Cholesky, inverse, and selected LU, QR and eigensolver projections. cuFFTDx contributes 1D FFTs and supported single-block 2D complex FFTs. And cuRANDDx contributes small floating-point and complex random expressions.

The qualifier on every row matters as much as the row itself. The cuSolverDx entry says selected LU, QR and eigensolver projections, and the cuFFTDx entry says supported single-block 2D complex FFTs. A library that fuses across a BLAS, a solver, a transform and a random number generator has a combinatorial space of combinations, and the documentation is careful to say which parts of that space are covered.

The behaviour when a combination is not supported is the part worth trusting. The README says JIT kernel fusion is experimental and intentionally rejects unsupported combinations rather than silently changing the execution model.

That sentence describes a specific and sensible engineering choice. A fusion system that quietly falls back to a different execution path for a case it does not handle would produce results that look identical and behave differently, and the failure would surface as a performance regression weeks later, in a benchmark, with nothing to connect it back to the change. Refusing the combination means the failure is immediate and obvious. The cost is that the supported set has to be read before the feature is relied on, and the fusion guide is where the current support matrix, the constraints and the compile-cache behaviour are documented.

The FFT benchmark and the numbers behind it

The README points at a reproducible benchmark example, examples/fft_benchmark.cu, which compares handwritten CUDA with cuFFT against MatX with the cudaExecutor, and with MathDx enabled it adds the JIT case.

The comparison is set up as three implementations of the same pipeline. First, handwritten CUDA normalization and post-processing kernels around a batched cuFFT plan. Second, one MatX expression with cudaExecutor, which uses cuFFT plus generated normalization and post-processing kernels. Third, the exact same MatX expression with CUDAJITExecutor, which fuses both sides of the FFT through cuFFTDx.

The claims that follow are specific. The implementation code shrinks by 92 percent, with 5 MatX lines replacing 59 lines of CUDA and cuFFT. The concise cudaExecutor path lands within about 3 percent of the handwritten implementation. And the same five line expression runs 1.61 times faster by changing only the executor.

The headline in the README puts the same result more bluntly: measured on an NVIDIA DGX Spark, 5 lines of MatX replace 59 lines of CUDA plus cuFFT and ran 1.61 times faster with JIT enabled. A table follows with median time and throughput columns for the implementations, and the benchmark repository also publishes a Python driver, black_scholes_benchmark.py.

The caveat is stated in the README itself, immediately after the number: actual performance depends on the operation, shape, data type, and target GPU. That is the sentence to keep in mind when reading any of these figures, including the ones in this article. A 1.61 times improvement for one FFT shape on one DGX Spark is a demonstration that the mechanism works, not a general statement about the library. The examples directory supports the point with a much wider spread of workloads than the FFT alone, including black_scholes, convolution, eigen decomposition, QR, interpolation, a recursive filter, a resampler, a MVDR beamformer, a Welch power spectrum, and a simple radar pipeline, so the honest evaluation is to find the example closest to your own workload and read its numbers.

The CPU path and the examples that show the scope

Two things in the repository broaden the picture beyond the GPU.

The first is the CPU backend list. The README names NVPL, FFTW, OpenBLAS and BLIS for multithreaded host execution, which are the CPU counterparts of the optimised GPU libraries the JIT path fuses. The same expression therefore has a path for a machine with no NVIDIA GPU, and the executor is what selects it.

The second is the examples directory, and its contents are a map of the library's actual capabilities. There are the linear algebra and solver examples, eigenExample.cu, qr.cu and cgsolve.cu. There are the signal processing ones, fft_benchmark.cu, fft_conv.cu, pwelch.cu, resample.cu, recursive_filter.cu and interpolate.cu. There are the machine learning ones, black_scholes.cu, conv2d.cu, convolution.cu and kernel_fusion.cu, which is the example that most directly demonstrates the fusion feature. And there are the domain ones, mvdr_beamformer.cu for beamforming and simple_radar_pipeline.cu for a radar processing chain, with a sarbp/ directory for synthetic aperture radar.

The spread is the point. The fusion mechanism is not tied to one workload, but the cases where it helps are the ones with many small operations between the large ones, which is what a radar pipeline or a filter bank looks like.

The rest of the tree shows a project with real infrastructure around the library. There is a CMakeLists.txt with a cmake/ directory, a docker/ directory, a bench/ directory, a test/ directory, a devcontainer configuration, a CITATION.cff file, a make_matx_release.sh script and an add_license.sh helper, plus an AGENTS.md and a CODE_OF_CONDUCT.md.

The licence is BSD 3-Clause, the default branch is main, the last push was on 2026-10-17, and the releases are v0.9.4 from 2025-10-27, v1.0.0 from 2026-03-04 and v1.1.0 from 2026-08-07.

Editorial conclusion

Adopt MatX if you are writing C++ that has to run on an NVIDIA GPU and you want array style expressions without giving up C++ control flow, and if the pipeline you care about can be expressed as compatible transforms. Do not adopt it for a workload that is not numerical, and do not plan around the JIT path until you have checked the fusion support matrix, because the feature is experimental and deliberately rejects combinations it cannot handle. Verify first that your shapes and data types fall inside the supported matrix, that the MathDx build option is what you want, and that the performance claim in the README, measured on one DGX Spark configuration, is revisited on your own GPU.

Frequently asked questions

Official sources

1. License: BSD-3-Clause
2. NVIDIA/MatX on GitHub
3. Project website
4. README
5. Releases