Library Overview#

GLASS is a header-only CUDA C++ library of composable __device__ primitives for small, block-local linear algebra and robotics math.

What GLASS is#

GLASS functions are __device__ helpers that operate on data in shared or device memory. It began as a small set of hand-rolled SIMT subroutines tuned for very small matrices — sizes where the launch and dispatch overhead of a vendor library would dominate the actual work — and has since grown into a unified single-block linear-algebra surface that also wraps NVIDIA’s state-of-the-art device-side libraries (CUB, cuBLASDx, cuSOLVERDx) under the same calling convention.

The intent is one predictable API vocabulary across execution scopes: pure SIMT for small problems, optional vendor kernels where they help, and warp/thread forms where packing more problems into a block is the better mapping.

The single-block execution model#

Every GLASS function assumes it runs within one CUDA thread block. The caller is responsible for launching one block per independent data item:

my_kernel<<<num_items, 256>>>(A, B, C, m, n, k);

Inside the kernel, all threads of the block cooperate on a single problem. This design enables composable GPU kernels for applications such as model-predictive control and rigid-body dynamics, where many small independent linear-algebra problems run in parallel — one per block.

Interfaces#

GLASS exposes three execution scopes — block, warp, and thread — in a dependency-free implementation family and, where NVIDIA provides a suitable device routine, an optional vendor family. This is a small matrix of explicit spellings rather than one flat list of interchangeable interfaces:

Interface

Scope

What it is

Header

glass::block:: (Block)

block

Hand-rolled SIMT, threadIdx.{x,y,z} / blockDim.* (no dependencies). The contract tier — bit-exact, thread-count invariant for deterministic-order ops, never re-dispatched

glass.cuh

glass::warp:: (Warp)

warp

Single-warp SIMT via __shfl_*_sync — selected L1/L2/L3 ops, no __syncthreads / shared

inline in the base L1/L2/L3 headers

glass::thread:: (Thread)

thread

Sequential, thread-per-problem — compile-time sizes, usually register-resident around N≤7 but correct and measured beyond it; branch-free ops only

inline in the base L1/L2/L3 headers

glass::nvidia::block:: (Nvidia)

block

CUB (L1) + cuBLASDx (L2/L3, batched) + cuSOLVERDx (LAPACK) — compile-time sizes only; plus glass::nvidia::warp:: CUB WarpReduce L1 reductions (one full 32-lane warp per problem)

glass-nvidia.cuh

glass::nvidia::warp:: (Nvidia warp)

warp

CUB WarpReduce L1 reductions, one full warp per problem and explicit per-warp scratch

glass-nvidia.cuh

glass::nvidia::thread:: (Nvidia thread)

thread

cuSOLVERDx 0.4+ LAPACK, one packed compile-time problem per CUDA thread; no dynamic shared scratch or block barrier

glass-nvidia.cuh

Bare glass::op is the measured-default face: the same block-scope calling contract, body chosen per (op, size, dtype) by glass::dispatch_body() (glass-dispatch.cuh). Measured cells may use a warp-0 or thread-0 implementation behind a wrapper; operations with no moved cell remain the same entity as glass::block::. The calling contract stays block-scoped in every case. Pin glass::block:: where determinism is load-bearing; see Namespaces, suffixes, and flags.

The interfaces intentionally overlap but are not interchangeable inventories. The API reference lists the exact operations and overloads. In general, glass::warp:: mirrors selected vector, dense, factor/solve, and fused families and requires a full warp. glass::thread:: is a smaller branch-free, compile-time subset; pivoted operations and reduction-strategy variants are deliberately absent.

Note

glass::cgrps:: (header glass-cgrps.cuh) is a convenience cooperative-groups alias of the Block interface — the same SIMT loop indexed via a g.thread_rank() / g.size() handle, with identical numerics. Use it from cooperative-groups code or to tile an arbitrary sub-block group; it is not a separately-tuned backend.

Many operations offer runtime (size as a function argument) and compile-time (size as a template argument) overloads. Reduction operations additionally offer _lowmem (no scratch, thread 0 accumulates) and _fast (warp-shuffle plus shared-memory inter-warp reduction) suffixed forms — e.g. glass::reduce_lowmem / glass::reduce_fast — keeping namespace = scope.

Higher-level solvers build on these primitives (and are likewise single-block): glass::bdmv (block-tridiagonal matvec) and glass::pcg (preconditioned conjugate gradient) for the block-tridiagonal SPD systems of trajectory optimization / MPC — see Block-tridiagonal Solves.

Warning

Factorizations do not check their input by default. CHECK defaults to false, so potrf / posv on a non-SPD matrix (or a non-pivoted ldlt hitting a zero pivot) silently produces NaN/Inf — there is no error return. To detect failure, instantiate with CHECK=true and pass an s_fail flag (and, for ldlt, the optional s_inertia pivot-sign counts); the reporting path compiles out entirely when CHECK is off. See examples/10_ldlt_solve.cu for the pattern on both a good and a zero-pivot matrix.

Choosing the right backend#

First pick the execution scope: block, warp, or thread. For a block-scoped implementation, three questions narrow the choice:

  1. Are sizes known at compile time?

  2. Is the matrix large enough that vendor-tuned tensor-core kernels matter?

  3. Can you launch with the thread count the backend wants?

Scenario

Use

Reason

Sizes only known at runtime

glass::gemm(m, n, k, ...)

Pure-SIMT, accepts dynamic args

Compile-time sizes, small matrices (≤ ~8×8), simple kernel

glass::gemm<float, M, N, K>(...)

Compiler unrolls inner loops; ~1 µs/op overhead is hard to beat for tiny sizes

Compile-time sizes, larger matrices, tensor-core hardware

glass::nvidia::block::gemm<float, M, N, K>(...)

cuBLASDx generates SM-specific tensor-core code

Compile-time sizes inside a kernel using a different thread count

glass::nvidia::block::gemm<float, M, N, K, TC>(...) with DEFINE_NVIDIA_GEMM_BLOCKDIM(M,N,K,TC)

Pins cuBLASDx’s BlockDim<TC,1,1>; lets you launch with any thread count ≥ TC

Need a transposed B / row-major storage in the NVIDIA path

glass::nvidia::block::gemm<...,LA,LB,LC> with DEFINE_NVIDIA_GEMM_BLOCKDIM_LAYOUT(...)

cuBLASDx Arrangement; no SIMT fallback needed

Linear solve Mx = b for SPD M

glass::nvidia::block::posv<float, N, NRHS>(...)

cuSOLVERDx fused factor + solve; faster than chol+trsm at N ≥ 8

General linear solve (non-SPD)

glass::nvidia::block::gesv_no_pivot<float, N, NRHS>(...)

cuSOLVERDx LU + solve

Least-squares / over- or under-determined

glass::nvidia::block::gels<float, M, N, NRHS>(...)

cuSOLVERDx QR (or LQ) + solve

BATCH independent GEMMs of the same shape, amortize launch

glass::nvidia::block::gemm_batched<...,BATCH,TC>

Single block, all batches active via threadIdx.y

When not to use glass::nvidia:::

  • Sizes only known at runtime (the templates require compile-time M, N, K).

  • You can’t add a DEFINE_NVIDIA_GEMM* macro for the size you need (the macro instantiation cost grows fast if you want every conceivable triple).

  • You’re on an SM cuBLASDx doesn’t tune for — it falls back to a generic config, and the pure-SIMT compile-time path is often competitive there.

The glass::nvidia::block::gemm<> / gemv<> / row_strided_* / gemm_batched_1d<> primary templates auto-dispatch at compile time: small shapes route to SIMT automatically without any DEFINE macro (a constexpr selection — nothing is decided at runtime). See Backend Dispatch for the full decision logic.

Next steps#

  • Installation — set up the headers and the optional MathDx backend.

  • Quickstart — a minimal end-to-end kernel.

  • Concepts — backend dispatch, TRAILING_SYNC, tuning, and batched-1D APIs.