100 DAYS OF GPU Course preview
Day 72Module 8
in-technical-review

Tensor cores with WMMA

If you time a WMMA kernel against an old FP32 kernel, the ratio combines three changes. The arithmetic moved from CUDA cores to tensor cores, the inputs shrank from FP32 to FP16, and the memory path was probably rewritten on the way. A number that mixes three changes cannot show the effect of one change.

This course already has a measured hand-written FP32 kernel: day 44's float4 kernel ran 2048-cubed in 4.219 ms, 67.7 percent of cublasGemmEx at FP32. By the end of this page you can write a WMMA matmul, say which of the three changes your speedup came from, and name the two things the API's fragment layout forbids you to do.

A warp owns a tile, and no thread owns an element

WMMA is one header, mma.h, one class and four calls in nvcuda::wmma. The unit of work is a fragment: "an overloaded class containing a section of a matrix distributed across all threads in the warp".

You declare three of them, fill or load two, multiply, and store the third. The whole warp makes every call together, so this is warp-level programming of the kind day 21 introduced.

It is not thread-level programming with a wider instruction.

The shape this page uses is 16x16x16 with __half inputs and a float accumulator, one row of the supported-combinations table in the programming guide. That table also allows a __half accumulator at the same shape, which is the exercise, and 32x8x16 and 8x32x16 for the same types. On a compute capability 7.5 card, __half is the only floating-point input the table offers; bf16, tf32 and double all need 8.0 or higher, which is day 71's format table restated as an API restriction.

The rule that shapes everything else is one sentence from the same page: "The mapping of matrix elements into fragment internal storage is unspecified and subject to change in future architectures." A fragment is not an array you can index by row and column.

Two members are documented, num_elements and x[num_elements], and they are for uniform element-wise work, scaling the whole tile by alpha. To read one element you store the fragment to memory and read it back.

You can measure the layout instead of guessing it. wmma_layout.cu multiplies two constant fragments, then has each lane overwrite its own accumulator slots with lane * 10 + slot and stores the tile, so the printed 16x16 grid says which lane held which element.

Executed on the Compiler Explorer Tesla T4, exit code 0, transcript in code/day72-wmma/evidence/ce-execute-2026-09-01.txt: lane 0's eight slots land at rows 0 and 8, columns 0, 1, 8 and 9. Lane 1 takes columns 2, 3, 10 and 11 of the same two rows.

Every lane holds eight of the 256 elements, scattered in two row groups. Nothing in the API promises that, and an sm_80 card is free to answer differently.

The compiler agrees about the count. Built at -arch=sm_75 with -ptx, that kernel's multiply comes out as one PTX instruction, wmma.mma.sync.aligned.row.row.m16n16k16.f32.f32, with eight destination registers per lane, and each store_matrix_sync as one wmma.store.d.sync.aligned line. Both are in code/day72-wmma/evidence/compile-2026-09-01.txt, on the compile side of that file, where nothing ran.

Diagram: which lane holds which element of a 16x16 accumulator tile. A 16 by 16 grid of cells on the left, one cell per element of the accumulator tile, with 32 lane colours; a strip of 32 lane boxes on the right, each labelled with the eight cells it owns. Band 1: the whole tile, coloured by owning lane. Caption "256 elements, 32 lanes, 8 elements each" Band 2: lane 0's eight cells highlighted. Caption "rows 0 and 8, columns 0, 1, 8 and 9, measured on a Tesla T4" Band 3: the same tile after an sm_80 relabelling, greyed with a question mark. Caption "unspecified by the API, so a different card may answer differently" Alt text: "One warp's 16 by 16 accumulator tile split across 32 lanes, eight elements each. On a T4, lane zero holds rows zero and eight at columns zero, one, eight and nine. The API promises none of it."

The habit this API takes away

Every kernel up to here let a thread identify its output. Day 16's tiled matmul computes row and col from threadIdx, guards them against the matrix edge, and writes one element. WMMA does not expose that mapping.

There is no if (row < m && col < n) you can write inside a fragment multiply, because the lane does not know which row and column it is holding and the answer is not the same on the next architecture.

Two consequences follow, and both add work to a GEMM. Ragged edges have to be handled outside the tile, by padding the matrix, by zero-filling a staging tile in shared memory, or by a separate clean-up kernel; this program sidesteps the question by using sizes that divide 16, and tilesDivide is a static_assert, not a comment.

A fused epilogue that depends on the element's position, a bias per column say, cannot read that position from the fragment. Element-wise work that is the same for every element is fine, which is what x[] is for.

// Day 16's habit. Inside a WMMA kernel there is no `row` and no `col`
// to test: the fragment does not tell a lane which element it holds.
if (row < m && col < n) {
    c[row * n + col] = acc;
}

Why this program calls cuBLAS twice

Full program in wmma_matmul.cu, under code/day72-wmma/. It computes one square C = A * B four ways at 256, 1024 and 2048: a WMMA kernel fed straight from global memory, a WMMA kernel fed from staged shared-memory tiles, cublasGemmEx with FP32 in and out, and cublasGemmEx with f16 inputs and an f32 output.

Both cuBLAS calls are there because a percentage needs the right denominator. The FP32 call is day 44's, byte for byte, so the wall clock on this page and the wall clock in day 44's evidence describe the same work. The f16 call runs the same dtypes as the WMMA kernels, and it is what the %f16row column divides by. Score a f16 kernel against an FP32 library call and you have measured a precision change and a hardware change together, which is the mistake in the hook.

The inputs are small integers, so correctness is exact and any mismatch is a bug. Day 44's fill is unchanged, values in [-3, 3] and [-2, 2]. They survive conversion to __half bit for bit, every product is an integer at most 6, and the largest dot product a case reaches is 6 * 2048 = 12,288, inside the range where a float represents consecutive integers. The accuracy story that needs non-representable inputs is day 71's; this page is about the units, so it removes rounding from the comparison entirely.

The inner loop is the whole API in five lines:

__global__ void matmulWmmaGlobal(const __half* __restrict__ a,
                                 const __half* __restrict__ b,
                                 AccT* __restrict__ c, size_t n, size_t k) {
    const unsigned int warp = threadIdx.x / 32;
    const size_t row0 =
        (static_cast<size_t>(blockIdx.y) * kWarpsPerBlock + warp) * kWmmaM;
    const size_t col0 = blockIdx.x * static_cast<size_t>(kWmmaN);

    wmma::fragment<wmma::matrix_a, kWmmaM, kWmmaN, kWmmaK, __half,
                   wmma::row_major>
        aFrag;
    wmma::fragment<wmma::matrix_b, kWmmaM, kWmmaN, kWmmaK, __half,
                   wmma::row_major>
        bFrag;
    wmma::fragment<wmma::accumulator, kWmmaM, kWmmaN, kWmmaK, AccT> accFrag;
    wmma::fill_fragment(accFrag, static_cast<AccT>(0.0f));

    for (size_t p = 0; p < k; p += kWmmaK) {
        wmma::load_matrix_sync(aFrag, a + row0 * k + p,
                               static_cast<unsigned int>(k));
        wmma::load_matrix_sync(bFrag, b + p * n + col0,
                               static_cast<unsigned int>(n));
        wmma::mma_sync(accFrag, aFrag, bFrag, accFrag);
    }

    wmma::store_matrix_sync(c + row0 * n + col0, accFrag,
                            static_cast<unsigned int>(n), wmma::mem_row_major);
}

That kernel reloads both operands from global memory on every K step, which is day 43's memory path carrying this day's arithmetic. matmulWmmaStaged keeps the same fragment calls and feeds them from a 32x16 A tile and a 16x32 B tile in shared memory, so a 2 x 2 arrangement of warps reads each staged fragment twice. Timing is CUDA events, three warm-ups and the mean of ten runs, copies excluded, and every path is checked against the CPU reference before anything is timed.

The accumulator type is one line, because the exercise changes it:

using AccT = float;

Results

Re-verified on the same Tesla T4 with driver 580.173.02 and CUDA 13.0 (V13.0.88). Both shipped programs again passed every correctness gate, and the compiler resource table stayed byte-for-byte stable.

At 2048, staged WMMA improved from 7.908 to 5.018 ms and now exceeds the 1.5x staging-speedup prediction. Both hardware transcripts are listed in front matter; the table below is the CUDA 13.0 shipped-program run. The FP16-accumulator exercise remains a separate CUDA 12.6 result.

size path ms GFLOP/s %f16row
256 matmulWmmaGlobal 0.035 952.5 55.8
256 matmulWmmaStaged 0.036 942.4 55.2
256 cublasGemmEx fp32 0.049 679.8 39.8
256 cublasGemmEx f16 0.020 1706.4 100.0
1024 matmulWmmaGlobal 1.571 1366.7 9.6
1024 matmulWmmaStaged 1.122 1914.3 13.5
1024 cublasGemmEx fp32 0.827 2596.4 18.3
1024 cublasGemmEx f16 0.151 14224.6 100.0
2048 matmulWmmaGlobal 10.998 1562.1 5.2
2048 matmulWmmaStaged 5.018 3423.4 11.3
2048 cublasGemmEx fp32 3.628 4735.0 15.6
2048 cublasGemmEx f16 0.567 30316.6 100.0

The transcript settles five predictions.

  1. Refuted. matmulWmmaStaged took 5.018 ms at 2048, slower than day 44's 4.219 ms target. The optimization ladder has a feed problem to solve.
  2. Held under CUDA 13.0. Staging improved 10.998 ms to 5.018 ms, comfortably past the predicted 1.5x. CUDA 12.6 had missed at 1.431x.
  3. Held. The 2048 FP16 cuBLAS call took 0.567 ms against 3.628 ms for its FP32 call, comfortably above the predicted 2x floor.
  4. Refuted. The staged kernel reached 11.3 percent of the FP16 cuBLAS row at 2048, below the predicted 15 to 40 percent.
  5. Held. Every path matched the CPU reference at every element and size, and both programs exited 0.

The CUDA 12.6 appended exercise run also settled the accumulator question: with AccT = __half, every path still matched at 256, 1024 and 2048. There was no first failing size. Its 2048 WMMA rows were 12.309 ms and 4.5 percent for global feed, then 7.168 ms and 7.7 percent for staged feed.

A percentage compares two programs, so print the denominator beside it. A WMMA kernel scored against an FP32 baseline reports two changes as one.

Run it yourself

The 16x16x16 shape and __half inputs require compute capability 7.0 or newer. wmma_layout.cu also fits within Compiler Explorer's 20-second cap. The embed uses compiler nvcc126u2 (NVCC 12.6.2) at -arch=sm_75: https://godbolt.org/z/MY1f4s4jW (minted and run 2026-09-01).

The measurement program takes longer. wmma_matmul.cu spends minutes in a double-precision CPU reference at 2048, so it belongs on a local or remote GPU session rather than a browser tab:

nvcc -std=c++17 -O3 -arch=sm_75 -o wmma_matmul wmma_matmul.cu -lcublas

Exercise

Rebuild wmma_matmul.cu with using AccT = __half;, run it, and report the first case size at which it stops passing and why. Then take the %f16row column for both WMMA kernels at 2048 from the shipped build's own run, say how far each sits from day 44's 67.7 percent, and say why those two percentages are not comparable.

Time: 25 to 40 minutes. Submit: the size where the run first fails, the failing path's printed element with the value it got and the value it wanted, the two %f16row numbers, and two sentences on the denominators.

Check: the rebuilt program's own correctness gate does the grading. Every path is compared against the double-precision CPU reference before anything is timed, and a failure prints the path, the size, the flat index, the row, the column, got and want, then returns EXIT_FAILURE. A build that passes all three sizes is a valid result; confirm the banner says accumulator type: __half, then report that no case failed.

Hint 1

The inputs are integers and the products are integers, so the only question is which integers the accumulator can hold exactly. Work out the largest dot product each size can produce, then ask what a 16-bit float does with it.

Hint 2

A __half has 10 stored mantissa bits, so it holds every integer up to 2048 exactly, then every second one to 4096, then every fourth to 8192. The largest dot product at size N is 6N. Which of 256, 1024 and 2048 crosses the line, and does crossing it force a wrong answer or only allow one?

Solution

The measured answer is that no size fails. The banner confirms AccT is __half, and all four paths match at 256, 1024 and 2048. Size 256 tops out at 1536, inside __half's consecutive-integer range.

The larger cases cross that range, but crossing it only permits rounding; this input pattern's partial sums still land on representable values. Magnitude alone therefore could not predict a failure.

The two percentages are not comparable because they have different denominators and different dtypes: day 44's 67.7 is a hand-written FP32 kernel against an FP32 library call on CUDA cores, and %f16row is a f16 tensor-core kernel against a f16 library call. The generalisation is the one the module keeps repeating: an accumulator is an accuracy decision and a storage format is a bandwidth decision, and a percentage is only a measurement when both sides of the ratio hold them fixed.

Pitfalls

36 compiler errors starting at your namespace alias. Building for a card below compute capability 7.0 gives <source>(49): error: name must be a namespace name pointing at namespace wmma = nvcuda::wmma;, because mma.h defines nothing there. The verbatim capture of that build, at -arch=sm_60, is in code/day72-wmma/evidence/compile-2026-09-01.txt. Fix: target sm_70 or higher.

You built for sm_70, it runs on your Turing card, and it is half speed. The Turing tuning guide says it plainly: "Volta binaries using Tensor Cores will only be able to reach half of Turing's Tensor Core peak performance. Recompiling the binary specifically for Turing would allow it to reach the peak performance."

The sm_70 build in the evidence file exits 0 with no diagnostics. Name the arch you are actually on.

The kernel hangs and no tool names a line. Every WMMA call must be reached by all 32 lanes: "these operations are allowed in conditional code only if the condition evaluates identically across the entire warp, otherwise the code execution is likely to hang". This is day 14's barrier rule applied to a whole warp, since a fragment call inside a if (row < n) guard divides a warp exactly the way a matrix edge does.

Padding your shared tile broke the fragment load. Day 15's fix for a bank conflict is to pad a tile row by one element. load_matrix_sync will not have it: ldm "must be a multiple of 8 for __half element type or multiple of 4 for float element type", and the pointer must be 256-bit aligned. Pad by 8 halves, or swizzle, which is the route CUTLASS takes.

You read frag.x[3] to find element (0, 3). The mapping is unspecified, and the captured transcript is not a specification. Uniform element-wise work over x[] is documented and safe; addressing a particular matrix element means storing the fragment and reading memory. mma.sync is the level where the register layout becomes yours to control, and it is the next day's subject.

Your "tensor-core speedup" is a dtype change. A f16 WMMA kernel against an FP32 baseline moves the arithmetic units and halves the input width at once. Report both denominators, the way this page's table does, or say which change you did not isolate. Mixed precision and the accumulator type are day 71's subject for exactly this reason, and cuBLAS picks its own path per dtype.

Go deeper

Next

Day 73 uses mma.sync and ldmatrix, where you place every value in the right lane's registers yourself and the layout this page could only print becomes something you specify. The f16 m16n8k8 shape and ldmatrix both require sm_75, not sm_80. Day 80 then takes whichever path measures better here and builds the capstone GEMM on it.