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

Size<M, N, K>

Yes

The per-tile problem, not the whole matrix. All three must be greater than 0.

Function<function::MM>

Yes

That this is a GEMM. Note the spelling — there is no function::gemm.

SM<Arch>

Yes

The target architecture, e.g. SM<800> for sm_80. An optional second parameter takes arch_specific or family_specific (the a / f suffixes), which are only legal for Arch >= 900 and are what unlock warp specialization and the newest MMA instructions.

CMake and descriptor spellings differ on three-digit architectures: CMake/nvcc 100a maps to SM<1000, arch_specific>(), not SM<100>. Keep the build architecture and descriptor synchronized; the performance-guide CMake does this through GUIDE_ARCH.

Block()

Yes, to execute

GEMM is block-level only — Thread() is TRSM-only. Without Block() you have a description you can query, but not run.

Precision<PA, PB, PC>

Defaults

Input and accumulation precision. One-, two- and three-argument forms are all legal.

Defaults to Precision<float, float, float> — so a descriptor that never mentions precision is an fp32 one, whatever your data is.

Type<type::real>

Defaults

Real or complex. Defaults to type::real.

Arrangement<A, B, C>

Defaults

Row- or column-majorness of each matrix, written with cublasdx::row_major / col_major. This is the basis for layout construction and instruction selection, so for good performance it should reflect how the data actually sits in global memory.

Defaults to <row_major, col_major, col_major>.

BlockDim<X, Y, Z>

Defaults

How many of the block’s threads do the GEMM. Flat size must be within [32, 1024], and the launch rule is unchecked.

Defaults to BLAS::suggested_block_dim, the tuning database’s choice.

BlocksPerCluster<M, N, K>

Optional

Cluster extents in CTAs; takes effect with WithPipeline(). Size remains the per-CTA tile.

N and K default to 1; K must be 1 and M * N must not exceed 16. Without the operator, the extents are 1 x 1 x 1 and clustered execution is disabled.

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.

LeadingDimension<LDA, LDB, LDC>

Optional

Only if your data is strided. LDA and LDB are both required when the operator is used; only LDC is optional.

Warning

Setting this disables the suggested-layout database and the suggested accumulator — both fall back silently. Prefer passing the leading dimension to get_layout_smem_*() at the call site.

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 float while the multiply runs in TF32 or FP16.

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 RequiredMantissaBits? See Emulated Precision (RequiredMantissaBits).

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 M * N, so a smaller tile may be necessary.

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 RequiredMantissaBits operator;

  • WithPipeline() and Block() are present, and Function<function::MM>();

  • Type<type::real>();

  • Precision<float, float, float>() or Precision<double, double, double>() (all three equal, nothing mixed).

At pipeline creation (suggest_emulation_pipeline, reached through suggest_pipeline / make_pipeline):

  • A and B must 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 K is bounded: pipeline_error_code::emulation_k_dimension_too_large is returned when K exceeds 33025, the largest K for which the int32 accumulation cannot overflow; split K across 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 C may be written, and padded K can be zero-filled.

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 C must stay inside the logical extent, but A and B can be padded (zeros in the K padding).

One extra shared-memory hop in the epilogue: write the accumulator to a shared C tile (partition_and_store), synchronize, then copy out with a per-element bounds check. There is no dynamic-shape C on the accumulator API of this library, so the check is yours.

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 Size is a compile-time value, and it must itself be legal for the instruction path (for example M a multiple of 64 for wgmma); a remainder that is not, or that is only known at run time, falls back to padding.

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.