Finding Problem Shape#
Note
This chapter is part of the agent-targeted Agent BLAS Performance Guide and is also published as Markdown for
agent consumption, in the markdown/ folder of the HTML documentation tree.
Describe the problem in the terms the library reasons about: what is multiplied, how much, where it lives, how many independent instances at once. Output of this chapter: an operator answer per row below, plus the fixed/free record from Axes That Can Be Relaxed.
Describing A GEMM#
Only three operators have no default; the rest resolve silently to something.
Operator |
Required |
What you have to decide |
|---|---|---|
|
Yes |
The per-tile problem, not the whole matrix. All three must be greater than 0. |
|
Yes |
That this is a GEMM. Note the spelling — there is no |
|
Yes |
The target architecture, e.g. CMake and descriptor spellings differ on three-digit architectures: CMake/nvcc |
|
Yes, to execute |
GEMM is block-level only — |
|
Defaults |
Input and accumulation precision. One-, two- and three-argument forms are all legal. Defaults to |
|
Defaults |
Real or complex. Defaults to |
|
Defaults |
Row- or column-majorness of each matrix, written with Defaults to |
|
Defaults |
How many of the block’s threads do the GEMM. Flat size must be within Defaults to |
|
Optional |
Cluster extents in CTAs; takes effect with
On supported paths, clustering can enable TMA multicast and 2-CTA MMA. The tile grid and launch must match the selected cluster; see the Hopper and Blackwell cards. |
|
Optional |
Only if your data is strided. Warning Setting this disables the suggested-layout database and the suggested accumulator — both fall
back silently. Prefer passing the leading dimension to |
Tuning-only operators include Alignment, StaticBlockDim, EnableInputStreaming, WithPipeline, and
RequiredMantissaBits. TRSM descriptors have additional operators and defaults that select the solve.
Composing The Descriptor#
One operator per answer, summed with +, captured with decltype. Order does not matter; what matters is
that every line is a decision you made rather than a default you inherited. The result is a type — everything
downstream hangs off it.
#include <cuda_fp16.h>
#include <cublasdx.hpp>
using namespace cublasdx; // operators live in namespace cublasdx
using GEMM = decltype(Size<128, 128, 32>() // per-tile M, N, K
+ Precision<__half, __half, float>() // A, B in half; accumulate in float
+ Type<type::real>() // real, not complex
+ Function<function::MM>() // this is a GEMM
+ Arrangement<row_major, col_major, col_major>()
+ SM<900, arch_specific>() // sm_90a
+ Block() // block-level operation
+ BlockDim<256>()); // threads spent on the GEMM
Partial descriptors are legal and reusable — this is how one axis is swept:
using problem = decltype(Size<128, 128, 32>() + Precision<__half, __half, float>() +
Type<type::real>() + Function<function::MM>() +
Arrangement<row_major, col_major, col_major>() + Block());
using on_ampere = decltype(problem() + SM<800>() + BlockDim<256>());
using on_hopper = decltype(problem() + SM<900, arch_specific>() + BlockDim<256>());
Axes That Can Be Relaxed#
Algorithmic-level relaxations are usually worth more than anything in the later chapters, and only you can make them — they trade accuracy, memory or interface constraints for speed. Record each axis as fixed by the application or free to experiment with; that record is the search space every later chapter operates in.
Axis |
Question to settle |
|---|---|
Computation precision |
Can it be lowered? Inputs may be The single highest-leverage question here. Experiment, and report the precision loss — a relaxation that silently degrades accuracy is not a win. |
Thread count |
Is it fixed by the surrounding kernel, or can it be swept? |
Tile shape |
Is the tile fixed by the surrounding problem, or can it be swept? |
Emulation |
Can an emulated-precision scheme be applied with |
Shared memory budget |
How much shared memory is free for the GEMM? Lowering precision lets more pipeline stages or a bigger tile fit in the same budget. |
Register budget |
Is the surrounding kernel register-heavy? Accumulator registers grow with |
Emulated Precision (RequiredMantissaBits)#
RequiredMantissaBits <X> turns a float or double pipelined GEMM into an integer-slice emulation: A and
B are split into ceil((X + 1) / 8) int8 slices under a shared exponent per row of A and per column of B, the
slice products run on the integer tensor cores into an int32 accumulator, and the result is recombined into the
descriptor’s floating-point type. It is an accuracy/throughput trade: X is the number of mantissa bits you ask the
result to preserve, so measure the error against a native reference before adopting it.
The descriptor constraints enforced by this library (each is a static_assert in blas_description.hpp):
exactly one
RequiredMantissaBitsoperator;WithPipeline()andBlock()are present, andFunction<function::MM>();Type<type::real>();Precision<float, float, float>()orPrecision<double, double, double>()(all three equal, nothing mixed).
At pipeline creation (suggest_emulation_pipeline, reached through suggest_pipeline / make_pipeline):
AandBmust be rank-2 tensors whose element type equals the descriptor’s precision (both accumulator modes work; the slice loop drives its internal integer pipeline in reusable mode either way);global
Kis bounded:pipeline_error_code::emulation_k_dimension_too_largeis returned whenKexceeds 33025, the largestKfor which theint32accumulation cannot overflow; splitKacross several pipelines and add the partial results for larger problems;the divisibility rules of Shapes That Do Not Divide The Tile apply to the sliced tensors as well.
The operator is a no-op request on some architectures: where the native float / double pipeline already
carries X bits, suggest_pipeline keeps it (see the operator reference). Runnable starting points are
example/cublasdx/16_dgemm_emulation (the block-level building blocks), 20_dgemm_emulation_pipeline and
21_sgemm_emulation_pipeline. Lifetime and stream rules of the emulation pipeline (it allocates and launches
slicing kernels when it is created) are in Pipeline Lifetime, Streams And Launch.
Shapes That Do Not Divide The Tile#
A pipeline has no partial-tile path. suggest_pipeline / make_pipeline return
pipeline_error_code::tile_divisibility (global dimensions are not divisible by the BLAS tile) unless M,
N and K are each a whole multiple of the Size tile, and pipeline_error_code::cluster_grid_divisibility
unless the tile grid is a whole number of clusters. The block-level API has no global shape at all: execute
works on shared-memory (or register) tensors whose shapes must be static and equal to Size (All layout shapes
must be static, only strides can be dynamic), so a ragged edge is the caller’s tile loop to handle. Three strategies:
Strategy |
Use when |
Cost |
|---|---|---|
Pad the allocation |
You own the allocation (or can copy once), the extra rows and columns of |
Simplest: no extra code and no extra compile time. Wasted work is the padded fraction of the problem. |
Pad the inputs, mask the store |
Only the write of |
One extra shared-memory hop in the epilogue: write the accumulator to a shared |
Bulk plus remainder descriptors |
Inputs cannot be padded, or the waste is unacceptable. |
General, but one more descriptor and kernel instantiation per distinct remainder, so more compile time.
The remainder |
Decision rule: pad if you can; mask the store when the output extent is contractual; split only when padding the inputs is impossible and the remainder is a few fixed, instruction-legal sizes.
Sketch of the split for M (N is symmetric; for K, either zero-pad or run the tail pipeline with
beta = 1 into the same C):
// Base: your descriptor without Size (precision, Arrangement<row_major, col_major, col_major>, Block, BlockDim,
// WithPipeline, SM, ...)
constexpr int tile_m = 128, tail_m = 64; // tail_m == M % tile_m, known at compile time
using Bulk = decltype(Base() + cublasdx::Size<tile_m, 256, 64>());
using Tail = decltype(Base() + cublasdx::Size<tail_m, 256, 64>());
int m_bulk = M - tail_m; // a whole number of tile_m rows
auto b = cublasdx::make_gmem_tensor<cublasdx::col_major>(b_ptr, K, N, ldb);
auto a0 = cublasdx::make_gmem_tensor<cublasdx::row_major>(a_ptr, m_bulk, K, lda);
auto p0 = cublasdx::suggest_pipeline<Bulk>(a0, b); // check p0, as always
kernel<Bulk><<<grid(m_bulk, N, Bulk{}), p0->get_block_dim(), smem0, stream>>>(p0->get_device_handle(), c_ptr, ldc);
auto a1 = cublasdx::make_gmem_tensor<cublasdx::row_major>(a_ptr + size_t(m_bulk) * lda, tail_m, K, lda);
auto p1 = cublasdx::suggest_pipeline<Tail>(a1, b);
kernel<Tail><<<grid(tail_m, N, Tail{}), p1->get_block_dim(), smem1, stream>>>(p1->get_device_handle(),
c_ptr + m_bulk, ldc); // col-major C
Each launch is an ordinary divisible problem, so everything else in the guide applies to it unchanged; the sub-tensor
base pointers must still satisfy the declared Alignment.