Last spring, I was a TA for CS179, Prof. Al Barr's GPU programming course. ¹ I was lucky enough to give a lectures on shared memory and matrix transpose, and while preparing it I couldn't find a single resource that answered all the questions students actually ask about it. So I wrote this.

The matrix transpose operation shows up inside attention layers, inside convolution layers, and more generally in matrix-heavy applications. So it's worth doing fast. But the reason it's the best way to learn about the GPU memory model is that a transpose does zero arithmetic. You're not adding or multiplying anything; you're only moving data from one place to another. Therefore, what determines how fast the kernel runs is how efficiently it uses memory bandwidth. Additionally, a transpose reads every element once and writes every element once, which is exactly what a memory copy does. An optimal transpose should thus run as fast as a plain memcpy of the same matrix.

#CUDA Memory Hierarchy

Each level of the CUDA thread hierarchy runs on a level of the hardware and can share a level of memory.

Global memory is the GPU's off-chip DRAM. It is gigabytes in size and accessible to every thread in the grid (and to the host), but it's slow. Shared memory sits on-chip inside each SM, so it's roughly an order of magnitude faster, but it's usually only tens of kilobytes per SM and private to a single thread block. Both are accessed per warp, and each has its optimal access pattern.

#Tiling

A big matrix won't fit in shared memory, so we break it into 64×64 tiles. Each tile is assigned to one thread block, and every block transposes its tile independently. Each block has 1024 threads, which is 32 warps of 32 threads. And 1024 threads covering $64 \times 64 = 4096$ elements means each thread handles 4 values.

The threads are arranged 64 wide and 16 tall, blockDim = (64, 16). Thread (tx, ty) owns column tx of the tile and handles rows ty, ty + 16, ty + 32 and ty + 48.

#Coalesced Memory Access

Global memory is read and written in aligned chunks of 128-byte transactions. ² When a warp issues a load or store, the hardware looks at the 32 addresses and works out the minimum number of transactions needed to serve them all. So the best case scenario would be that thread 0 accesses 0x0, thread 1 accesses 0x4, thread 2 accesses 0x8, … 32 consecutive floats, 128 bytes in total. That's one transaction. This is a coalesced access. The worst case scenario is when every thread accesses an address far from the others, say with a stride of a whole matrix row. Now the warp needs up to 32 separate transactions, and of each 128 bytes fetched you use only 4 which is 32× the bandwidth for the same data.

So the optimal access pattern is that every global read and every global write should be coalesced.

#The Naive Transpose Kernel

__global__ void naiveTransposeKernel(const float *input, float *output, int n) {
    int col = blockIdx.x * 64 + threadIdx.x;
    int row = blockIdx.y * 64 + threadIdx.y;

    for (int k = 0; k < 64; k += 16)
        output[col * n + (row + k)] = input[(row + k) * n + col];
}

Observe that the read is good, but the write is completely uncoalesced. If we flip it, reading columns and writing rows, then the read is uncoalesced instead. With only global memory, one side of a transpose is always strided. It is obvious that we cannot coalesce both the read and the write without an intermediate step.

1. Load from global memory row by row. Each warp reads 32 adjacent floats.
2. Store those values into a shared-memory tile in the same row order.
3. Read back from shared memory transposed, so what was a column of the tile becomes a row.
4. Write that row to global memory. Each warp writes 32 adjacent floats.

The load side is 4 coalesced loads per warp, and $32 \times 4 = 128$ coalesced transactions bring in the whole tile. The store side is the same. Between steps 2 and 3 we need a __syncthreads(). The column a warp reads in step 3 contains values that other warps stored in step 2, so every warp must finish storing before any warp starts reading.

#Bank conflicts

Bank conflicts play the same role that coalescing plays for global memory. Shared memory is split into 32 banks, interleaved every 4 bytes: bytes 0–3 are in bank 0, bytes 4–7 in bank 1, …, and the 33rd float wraps back to bank 0. When a warp accesses shared memory:

• If all 32 threads hit different banks, the access completes at full speed, all 32 at once.
• If several threads hit the same bank at different addresses, those accesses are serialized. A 2-way conflict costs two passes.
• If several threads hit the same address, it's a broadcast, and that's free.

Now apply this to our tile. We have __shared__ float tile[64][64], element tile[r][c] lives in bank $(64r + c) \bmod 32 = c \bmod 32$. the 32 threads of a warp have consecutive c, so they hit 32 different banks. We have no conflict. But, the 32 threads have consecutive r and the same c, so they all hit the same bank. That's a 32-way conflict on every read.

The fix is so easy. We can just pad each row by one float: __shared__ float tile[64][65]. Now tile[r][c] lives in bank $(65r + c) \bmod 32 = (r + c) \bmod 32$, so walking down a column steps through the banks one by one. Row accesses are still conflict-free, and column accesses are now conflict-free too.

#The Optimal kernel

__global__ void optimalTransposeKernel(const float *input, float *output, int n) {
    // One padding column so that reading a column is bank-conflict free.
    __shared__ float tile[64][65];

    int tx = threadIdx.x, ty = threadIdx.y;

    // Steps 1-2: coalesced read of input tile (blockIdx.y, blockIdx.x).
    int inCol = blockIdx.x * 64 + tx;
    int inRow = blockIdx.y * 64 + ty;
    for (int k = 0; k < 64; k += 16)
        tile[ty + k][tx] = input[(inRow + k) * n + inCol];

    __syncthreads();

    // Steps 3-4: that tile lands at (blockIdx.x, blockIdx.y) of the output.
    int outCol = blockIdx.y * 64 + tx;
    int outRow = blockIdx.x * 64 + ty;
    for (int k = 0; k < 64; k += 16)
        output[(outRow + k) * n + outCol] = tile[tx][ty + k];
}

// Launch: one 64x16 block per 64x64 tile.
dim3 block(64, 16);
dim3 grid(n / 64, n / 64);
optimalTransposeKernel<<<grid, block>>>(d_input, d_output, n);

In both global accesses, input[(inRow + k) * n + inCol] and output[(outRow + k) * n + outCol], the index that varies across a warp is a consecutive column (inCol or outCol), so both are coalesced. In shared memory, the row-wise store tile[ty + k][tx] puts the warp's 32 threads in banks $(r + tx) \bmod 32$ (no conflict). The column-wise read tile[tx][ty + k] puts them in banks $(tx + c) \bmod 32$, which is also conflict-free, thanks to the padding column.

#Profiling with Nsight Compute

I profiled three kernels in Nsight Compute on a $4096 \times 4096$ float matrix (67 MB) on an RTX A5000 (768 GB/s peak DRAM bandwidth): the naive kernel, the shared-memory kernel without padding, and the optimal kernel.

The naive kernel takes 961 µs. Look at the L1 → L2 link where 558 GB/s leaves L1 while only about 70 GB/s reaches DRAM. That is 8× more traffic than needed, because of the strided writes.

The shared-memory kernels take about 193 µs, roughly 5× faster. Shared memory is now busy, the L1 hit rate drops to 0% because every global address is touched exactly once, and L1 → L2 traffic matches DRAM traffic. Here, padding makes no difference here because the kernel is bound by DRAM, which hides the bank conflicts.

The optimal kernel moves about 687 GB/s through DRAM, close to 90% of the A5000's peak, which is about what a device-to-device copy achieves. It moves the same bytes a copy would, at the rate a copy would, which is the bound we wanted to reach at the start.

#References

1. [1]
   See Week 2 lecture notes: https://courses.cms.caltech.edu/cs179/2026_lectures/cs179_2026_lec06.pdf
2. [2]
   On recent NVIDIA GPUs the real unit is a 32-byte sector, and a 128-byte cache line is four sectors. The 128-byte picture is the right mental model for a warp of 32 floats: a perfectly coalesced float access touches exactly one line (four sectors), and a fully strided one touches 32 different lines.