The Performance of GPU Stream-Assisted Asynchronous Multiplications for Deep Neural Networks

Introduction

Deep neural networks (DNNs) have become the computational workhorses behind everything from real-time language translation to autonomous driving. But beneath their remarkable capabilities lies an insatiable hunger—not for data alone, but for dense matrix multiplications. Every forward pass, every backward propagation, every gradient update ultimately boils down to tens of thousands of General Matrix Multiply (GEMM) operations, each one shoving billions of floating-point numbers through graphics processing units (GPUs) as fast as the silicon will allow.

For years, this arrangement worked spectacularly well. GPU vendors poured transistor budgets into specialized Tensor Cores and increased memory bandwidth, while software frameworks like cuBLAS and cuDNN squeezed near-theoretical peak performance from the hardware. But as models cross the trillion-parameter mark and inference latency requirements shrink to milliseconds, a subtle bottleneck has emerged: the GPU’s own internal scheduling model. The conventional approach—launching one large, synchronous GEMM, waiting for it to finish, then launching the next—leaves compute resources idle while data transfers or lower-priority kernels are waiting. The GPU is a massively parallel machine, yet its internal pipelines often stall because they are forced to operate in lockstep.

This is where CUDA streams enter the narrative. Streams allow a programmer to break a single GPU kernel into multiple independent sequences of operations that can overlap in execution, hiding latency and maximizing utilization. When applied to the core GEMM operations of a DNN, stream-assisted async multiplication offers a promising path to reclaim the lost cycles—but the path is not without its own pitfalls. Overlapping data transfers with computation, partitioning matrix dimensions across streams, and managing synchronization overhead all require careful engineering.

In this article, we will take a deep dive into the performance of GPU stream-assisted asynchronous multiplications for deep neural networks. We’ll start by examining exactly where and how matrix multiplications appear in modern DNNs, then dissect the GPU architecture to understand why synchronous execution leaves performance on the table. We’ll build a thorough understanding of CUDA streams, then present a detailed methodology for implementing asynchronous GEMMs. Through concrete code examples, performance analysis, and case studies of real DNN layers, we’ll quantify the benefits and address the pitfalls. By the end, you’ll have a clear roadmap for leveraging streams to accelerate deep learning workloads, whether you’re training the next GPT model or deploying an efficient real-time inference pipeline.

The Ubiquity of Matrix Multiplications in Deep Learning

Deep learning may seem like a zoo of different layer types, but at the hardware level almost everything reduces to matrix multiplication. Let’s walk through the major operations.

Fully Connected Layers

Consider a fully connected (FC) layer with input size (M) and output size (N) over a batch of (B) samples. The forward pass computes:

[ Y = X \cdot W ]

where (X) is a ([B \times M]) matrix and (W) is ([M \times N]). This is a direct GEMM of size (B \times M \times N). In practice, (M) and (N) can be thousands, and (B) often ranges from 1 (inference) to 128 or more (training). The backward pass requires three more GEMMs: one for the gradient w.r.t. weights ((dW = X^T \cdot dY)), one for the gradient w.r.t. input ((dX = dY \cdot W^T)), and sometimes an update step. For a network with hundreds of FC layers, the GEMM count multiplies.

Convolutional Layers

Convolutions are traditionally implemented via the im2col transformation, which flattens each convolutional patch into a column, producing a matrix (\text{im2col}(X)) of shape ([K \times C \cdot H_f \cdot W_f]) (where (K) is the number of patches) and then performing a GEMM with the filter matrix of shape ([C \cdot H_f \cdot W_f \times F]). The result is then reshaped. While newer algorithms like Winograd or direct convolution exist, the GEMM-based approach remains dominant in many frameworks because it leverages highly optimized cuBLAS.

Attention Mechanisms

Transformers have brought attention layers—specifically scaled dot-product attention—to the forefront. The core operation is:

[ \text{Attention}(Q, K, V) = \text{softmax}\left(\frac{Q K^T}{\sqrt{d_k}}\right) V ]

Here (Q), (K), (V) are matrices of shape ([B \times L \times d]) (sequence length (L), head dimension (d)). The multiplication (Q K^T) is a batched GEMM (batching over heads), and the subsequent multiplication with (V) is another. With multi-head attention, we have four GEMMs per layer (including the projection matrices for Q, K, V). For large language models with hundreds of layers and sequence lengths up to 8192, this dominates compute.

Recurrent Layers and Others

Even layers like LSTMs or GRUs are essentially composed of several GEMMs inside a loop. The point is clear: DNNs are GEMM factories. Any improvement in GEMM throughput directly translates to faster training and inference.

GPU Architecture and the Synchronous GEMM Bottleneck

To understand why asynchronous execution matters, we must first understand how a GPU executes a single GEMM kernel.

The SIMT Model

A GPU consists of multiple Streaming Multiprocessors (SMs), each containing a set of CUDA cores, shared memory, register files, and warp schedulers. Work is organized into warps (groups of 32 threads) that execute instructions in lockstep—SIMT (Single Instruction, Multiple Threads). The GPU scheduler issues warps from a pool of thread blocks (CTAs) that are assigned to SMs.

Memory Hierarchy

• Global memory: Large (tens of GB), high latency (~300-800 cycles), accessible by all threads.
• Shared memory: Small (tens of KB per SM), low latency (~5-10 cycles), on-chip, shared within a thread block.
• Registers: Fastest, but limited per thread (typically 255 32-bit registers per thread on modern architectures).

A typical GEMM kernel (e.g., cuBLAS’s own) works like this:

1. Load tiles of input matrices from global memory to shared memory (coalesced access).
2. Compute partial products in shared memory / registers.
3. Accumulate results and write back to global memory.

This is a memory-bound operation for small matrices, but for large matrices it becomes compute-bound, limited by the number of floating-point operations per second (FLOPS) and the availability of Tensor Cores.

The Synchronous Bottleneck

When you launch a GEMM kernel synchronously (e.g., cublasSgemm followed by cudaDeviceSynchronize), the CPU issues the kernel launch to the GPU command queue. The GPU starts executing it on as many SMs as possible. However, during execution:

• If the kernel is memory-bound, many warps stall waiting for global memory loads, leaving SMs underutilized. The GPU may hide some latency by switching warps, but there is a limit.
• If the kernel is compute-bound, the Tensor Cores are fully utilized, but the launch itself introduces overhead (typically a few microseconds for the driver to handle the command). For very large GEMMs this overhead is negligible, but for the many small GEMMs in attention or backward passes, overhead becomes significant.
• Worse, no other kernel can run on the same GPU while this kernel is executing unless you explicitly overlap them with streams. In synchronous mode, consecutive GEMMs are serialized. The GPU sits idle between the end of one kernel and the launch of the next (driver overhead, kernel setup). For a chain of 100 GEMMs, these tiny idle periods add up.

Tensor Cores and Synchronous Limits

NVIDIA’s Tensor Cores perform fused multiply-add on 4x4 matrices in one cycle, achieving massive throughput. However, Tensor Cores require data in a specific layout (e.g., row-major floats or half-precision). The cuBLAS library handles this. But Tensor Core usage does not break the synchronous model; the whole kernel still runs to completion before the next kernel begins.

CUDA Streams: A Primer

A CUDA stream is a sequence of operations that are executed in order on the GPU. Operations from different streams can be interleaved, and in many cases can run concurrently. Streams are the fundamental mechanism for expressing asynchronous, overlapping work.

Stream Creation and Usage

cudaStream_t stream1, stream2;
cudaStreamCreate(&stream1);
cudaStreamCreate(&stream2);

You then pass the stream to any asynchronous CUDA API:

cudaMemcpyAsync(dst, src, size, cudaMemcpyHostToDevice, stream1);
kernel<<<grid, block, 0, stream1>>>(args);
cublasSgemm(handle, ... , stream1);

By default, operations use the default stream (stream 0), which is synchronous with the host and with other streams unless you enable stream-ordered semantics. All operations in the same stream are ordered; operations in different streams may run concurrently if the hardware supports it.

Overlapping Data Transfers and Computation

The classic use case is to overlap memory copies with kernel execution. For example, while one stream copies data from host to device, another stream can execute a kernel on already-resident data.

// Stream 0: copy input A
cudaMemcpyAsync(d_A, h_A, size, cudaMemcpyHostToDevice, stream0);
// Stream 1: compute with B (already on device)
kernel_B<<<grid, block, 0, stream1>>>(d_B, d_C);

The GPU can simultaneously move data over the PCIe/ NVLink bus and compute on one or more SMs. This hides the latency of data transfers.

Concurrent Kernel Execution

Modern GPUs (compute capability 3.5+) can run multiple kernels concurrently on different SMs, as long as resource constraints (registers, shared memory, thread blocks) allow. For example, if one kernel uses half the SMs, another kernel can use the other half. This is where streams shine for DNNs: instead of one large synchronous GEMM, you can split a large GEMM into smaller sub-GEMMs, each launched on its own stream, and they can run concurrently, increasing GPU utilization.

Events and Synchronization

To coordinate streams, CUDA provides events:

cudaEvent_t event;
cudaEventCreate(&event);
cudaEventRecord(event, stream0); // record after all ops in stream0
cudaStreamWaitEvent(stream1, event, 0); // stream1 waits for event

This allows fine-grained synchronization without host intervention.

Stream-Assisted Asynchronous Multiplication: Methodology

The core idea is to decompose a single large GEMM into multiple independent sub-GEMMs that can be executed concurrently on different streams, and to overlap the data transfers (if any) with computation.

Partitioning the GEMM

A GEMM C = A * B (with dimensions M x K and K x N) can be partitioned in several ways:

1. Split along M dimension: Divide rows of A and C into chunks. Each stream computes a sub-matrix of C = A_chunk * B.
2. Split along N dimension: Divide columns of B and C. Each stream computes C = A * B_chunk.
3. Split along K dimension: This is more complex because it requires reduction (accumulating partial sums). It is less common for stream-based concurrency because it introduces dependencies.

The simplest is splitting along M or N, as each sub-GEMM is independent.

For a batch dimension (e.g., batched GEMM in attention), we can split the batch across streams.

Assignment to Streams

Suppose we have a GEMM of size (M=8192, K=4096, N=4096) and we create 4 streams. We split M into 4 chunks of 2048 rows each. Each stream launches a cuBLAS cublasGemmEx (or cublasSgemm) on its sub-matrix. The cuBLAS handle must be associated with a stream (via cublasSetStream). The four kernels can run concurrently if the GPU has enough resources.

Overlapping Memory Transfers

In many inference scenarios, the weights (B matrix) are already on the device, but the input (A) may be streamed from host. In that case, we can use double buffering: while one stream computes on chunk i of A, another stream transfers chunk i+1. This requires careful synchronization.

Synchronization Strategy

After all sub-GEMMs complete, we must combine the results (if necessary) or proceed to the next layer. We can record events in each stream after its kernel, then wait on those events from the next stream that needs the full result. For example, if the next layer requires the entire C matrix, we can have a single synchronization barrier that waits for all streams.

Implementation Sketch

// Assume 4 streams
cudaStream_t streams[4];
for (int i=0; i<4; i++) cudaStreamCreate(&streams[i]);

// Partition M into 4 chunks
int chunkM = M / 4;
cublasSetStream(handle, streams[i]);

for (int i=0; i<4; i++) {
    float* A_chunk = d_A + i * chunkM * K;  // row-major
    float* C_chunk = d_C + i * chunkM * N;
    cublasSgemm(handle, CUBLAS_OP_N, CUBLAS_OP_N,
                chunkM, N, K, &alpha,
                A_chunk, K, d_B, N, &beta,
                C_chunk, N);
}
// Sync all streams
cudaDeviceSynchronize();

However, this naive approach might not achieve concurrency due to resource contention. The GPU schedules thread blocks from all kernels. If the kernel uses many registers or shared memory, the number of resident blocks per SM is limited. With 4 kernels trying to occupy the same SMs, they may serialize. We need to tune the chunk size and kernel parameters.

Implementation Details and Code Examples

Let’s dive deeper into practical implementation.

Using cuBLAS with Streams

cuBLAS provides two ways to associate a stream with a handle:

• cublasSetStream(cublasHandle_t handle, cudaStream_t streamId) sets the stream for all subsequent cuBLAS calls on that handle.
• Alternatively, each cuBLAS routine accepts a stream argument in newer versions.

We must create separate cuBLAS handles per stream or reuse one handle and switch streams before each call. Creating multiple handles is safer to avoid race conditions, but uses more memory.

Example: Splitting a large GEMM across 2 streams

cublasHandle_t handle1, handle2;
cublasCreate(&handle1);
cublasCreate(&handle2);
cudaStream_t s1, s2;
cudaStreamCreate(&s1);
cudaStreamCreate(&s2);
cublasSetStream(handle1, s1);
cublasSetStream(handle2, s2);

int M=8192, K=4096, N=4096;
int halfM = M/2;
float alpha=1.0f, beta=0.0f;

// Stream 1 computes top half
cublasSgemm(handle1, CUBLAS_OP_N, CUBLAS_OP_N,
            halfM, N, K, &alpha,
            d_A, K, d_B, N, &beta,
            d_C, N);

// Stream 2 computes bottom half
cublasSgemm(handle2, CUBLAS_OP_N, CUBLAS_OP_N,
            halfM, N, K, &alpha,
            d_A + halfM*K, K, d_B, N, &beta,
            d_C + halfM*N, N);

// Wait for both
cudaStreamSynchronize(s1);
cudaStreamSynchronize(s2);

To maximize concurrency, we may want to launch both kernels quickly from the host. The host can launch both without waiting, and the GPU scheduler will try to execute them together.

Custom GEMM Kernel with Streams

Instead of relying on cuBLAS, we can write a custom GEMM that is deliberately underutilized per thread block to allow more kernels to run concurrently. For example, a naive GEMM using fewer registers per thread leaves more room for other kernels’ thread blocks. But performance may suffer.

A balanced approach is to use cuBLAS for the inner GEMMs but adjust the block size. However, we cannot control cuBLAS internals. An alternative is to use CUTLASS, NVIDIA’s open-source template library for GEMM, which allows tuning tile sizes and launch parameters. With CUTLASS, we can create a kernel that uses fewer SMs (by specifying a smaller grid) and then launch multiple instances on different streams.

CUTLASS example (simplified):

// Assuming CUTLASS::Gemm is defined
using Gemm = cutlass::gemm::device::Gemm<float, cutlass::layout::RowMajor>;
Gemm gemm_op;

// For each stream, set up parameters with sub-matrix offsets
for (int i=0; i<num_streams; i++) {
    Gemm::Arguments args(
        {chunkM, N, K},  // problem size
        {A_chunk, K},    // tensor ref
        {d_B, N},
        {C_chunk, N},
        {alpha, beta},
        1  // splitK = 1 (no reduction)
    );
    // Launch on stream
    gemm_op(args, stream[i]);
}

Overlapping Data Transfers

A critical use case is when input data arrives from host (e.g., real-time inference). The data for the next layer can be transferred while the current layer computes. This requires a pipeline:

Pipeline design:

• Stream 0: Transfer A (layer1) from host to device.
• Stream 1: Compute layer1 with A, then transfer A for layer2 (from host).
• Stream 2: Compute layer2, etc.

Each stream has a double buffer. Events ensure that transfers don’t overwrite buffers still in use.

// Pseudo-code for inference with two streams
cudaStream_t transferStream, computeStream;
// ... create events

while (input_available) {
    // Transfer next input to buffer[ping]
    cudaMemcpyAsync(d_input[ping], h_input, size,
                    cudaMemcpyHostToDevice, transferStream);
    // Compute with previous buffer[pong] (if ready)
    cudaStreamWaitEvent(computeStream, transferDoneEvent, 0);
    kernel<<<grid, block, 0, computeStream>>>(d_input[pong], ...);
    // Swap ping-pong
    ping = 1 - ping;
}

This hides the PCIe transfer latency behind computation.

Performance Analysis: Metrics and Benchmarks

To evaluate stream-assisted async multiplication, we need to define the right metrics.

Latency vs Throughput

• Latency: time to complete a single forward pass for a given batch. With streams, we may increase latency if we split work inefficiently, but overall throughput (samples per second) can increase because the GPU is utilized more continuously.
• Throughput: number of GEMMs completed per second. This is the key metric for batch training and inference.

Occupancy and Warp Utilization

Occupancy is the ratio of active warps to maximum warps per SM. High occupancy helps hide memory latency. When launching multiple kernels concurrently, occupancy per kernel may drop because each kernel uses fewer SMs. However, total occupancy across all kernels may be higher if the kernels use different resources.

For example, a single large GEMM may occupy 80% of SMs with 100% occupancy per SM. Split into 4 sub-GEMMs, each may use 20% of SMs, but total occupancy might be 80% again (if no resource collisions). But if each sub-GEMM uses a lot of shared memory, the number of concurrent thread blocks per SM is limited, so total occupancy may increase less.

We can use NVIDIA Nsight Compute to profile occupancy and memory throughput.

Overlap Efficiency

Overlap efficiency measures how well two operations (e.g., copy and compute) run simultaneously. Ideally, the sum of their execution times should be less than the time of sequential execution. For overlapping data transfers with compute, a typical efficiency is 70-90% on modern GPUs with NVLink.

For concurrent kernel execution, the efficiency depends on resource contention. We can measure: speedup = time(sequential) / time(concurrent). With two streams, speedup can be close to 2 for memory-bound kernels if they use different memory channels. For compute-bound kernels, speedup is limited by SM count and power.

Baseline Comparison

We compare:

• Synchronous baseline: one cuBLAS call per GEMM, serial.
• Streamed: split GEMM into N streams, no overlapping of transfers.
• Overlapped: streams plus data transfer overlapping.

Experimental setup:

• GPU: NVIDIA A100 (108 SMs, 40GB HBM2e)
• GEMM: float32, dimensions (M=N=K=8192)
• Stream count: 1, 2, 4, 8, 16

Results (hypothetical):

For large matrices, splitting into streams actually hurts performance due to kernel launch overhead and reduced Tensor Core utilization (each sub-GEMM may be too small to fully occupy Tensor Cores). However, for small batch GEMMs (e.g., M=64, N=64, K=1024), the story reverses.

Small GEMM example:

Here, the ability to run 4 sub-GEMMs concurrently on different SMs nearly doubles throughput, because each small GEMM alone cannot saturate the GPU.

Case Studies: Applying to DNN Layers

Fully Connected Layer (Inference)

In inference with a batch size of 1, each FC layer is a matrix-vector multiply (MxN with M=1). This is a tiny GEMM. Synchronous cuBLAS will have high launch overhead relative to compute. Using streams, we can batch multiple layers’ weight matrices into a single stream pipeline? Actually, each layer depends on the previous layer’s output, so we cannot parallelize layers across streams. But we can parallelize within a layer if the weight matrix is large. For batch size 1, the GEMM is small; streams may not help. However, if we have multiple input samples (batch size >1), we can split the batch across streams for each layer.

Convolutional Layer via im2col

Consider a convolution with input (batch=64, channels=3, height=224, width=224), filters (64, 3, 3, 3). After im2col, we get (64, 27, 27, 27) patches? Actually, careful calculation yields a GEMM of shape (batch _ output_spatial) x (filter_size _ input_channels) - say (64*484=30976) x (27) by (64) filters. That’s a large GEMM. Splitting across batch dimension makes sense: each stream handles a subset of the batch. The compute and memory requirements per stream shrink, allowing concurrent execution and overlapping of im2col (which is also a kernel).

Attention Mechanism

In multi-head attention with batch=64, heads=12, seq_len=512, d=64. The QK^T multiplication is a batched GEMM of many small matrices: each head has (64, 512, 64) x (64, 512, 64)^T → (64, 512, 512). That’s 12 such GEMMs. They are independent across heads. We can assign each head to a different stream, achieving concurrent kernel execution. This is a perfect scenario for streams.

Implementation: For each head, launch cublasGemmStridedBatchedEx with a batch of size 64 (the batch dimension), but we can split the batch across streams as well. However, each individual GEMM (within a head) is 512x64 x 64x512 = 512x512. That’s medium-sized (512^2 = 262K elements). On an A100, such a GEMM takes about 20 µs. With 12 heads, sequential would be 240 µs. With 12 streams, if the GPU can accommodate 12 concurrent kernels, the time could drop to ~20 µs, a 12x speedup. In practice, limited by SM count (108 SMs), we can run maybe 4-8 concurrent kernels efficiently. Still, significant gain.

Backward Pass

The backward pass through a fully connected layer requires three GEMMs: dW, dX, and sometimes weight update. dW and dX can be independent (dW= X^T _ dY; dX= dY _ W^T). They can be launched on separate streams concurrently. This parallelizes the backward computation within a single layer, reducing the per-layer time.

Pitfalls and Challenges

Stream Synchronization Overhead

Every cudaStreamSynchronize or cudaEventSynchronize introduces a host-device round trip that can cost tens of microseconds. Overuse can negate gains. Use events and cudaStreamWaitEvent to avoid host synchronization.

Memory Contention

When multiple kernels access global memory simultaneously, they compete for memory bandwidth. If all sub-GEMMs access the same memory controller, performance degrades. Using different memory regions (e.g., different row stripes) can help, but total bandwidth is shared.

Diminishing Returns with Small Matrices

For very small GEMMs (e.g., <256 elements), kernel launch overhead dominates. Streams exacerbate this because you now have multiple small launches. In such cases, it’s better to use batch GEMM (cublasGemmBatched) which groups many small GEMMs into one kernel.

Debugging and Profiling

Concurrent asynchronous operations are notoriously hard to debug. Race conditions may not manifest every run. Use cuda-gdb with -g flags and careful event logging. Profiling with Nsight Systems can show timeline visualization of streams.

Resource Limits on Concurrent Kernels

The GPU can run at most a limited number of concurrent kernels. On older architectures, only 32 concurrent kernels; on Ampere, up to 128. However, each kernel consumes hardware resources (warp schedulers, shared memory). The number of concurrent kernels is also limited by SM occupancy.

Advanced Techniques

Dynamic Stream Management

Instead of statically assigning streams, we can have a thread pool on the CPU that launches kernels as resources become available. For example, a dispatch queue that monitors active kernel count and launches new ones only when occupancy is below threshold.

Multi-Stream with Priority

CUDA streams can have priorities (low, normal, high). High-priority streams get preference in scheduling. For time-critical layers (e.g., early in network), we can assign higher priority.

Combining with CUDA Graphs

CUDA Graphs capture a sequence of operations (including kernel launches, memory copies, events) into a graph that can be launched repeatedly with minimal overhead. For DNN layers with fixed dimensions (e.g., during inference), we can build a graph that includes streamed GEMMs and synchronization. This eliminates per-launch overhead, making stream-assisted multiplication almost free.

Example: Create a graph with two streams for a layer’s forward and backward GEMMs, then replay the graph each iteration. The graph launch overhead is a few microseconds, regardless of number of kernels.

Future Directions: MIG and NVLink

Multi-Instance GPU (MIG) partitions the GPU into isolated instances. Streams across instances work independently. For large models that fit on one MIG slice, stream concurrency is limited to that slice. However, NVLink allows multiple GPUs to work together. Streams can be used to overlap communication (via NVLink peer access) with computation on local GPU.

Conclusion

GPU stream-assisted asynchronous multiplication is not a silver bullet, but it is a powerful tool in the deep learning engineer’s optimization arsenal. As we have seen, the benefits are most pronounced for small-to-medium GEMMs, batched operations like attention, and when overlapping data transfers with computation. The pitfalls—synchronization overhead, resource contention, and debugging complexity—require careful measurement and tuning.

The key takeaway is that the GPU’s massive parallelism is not automatically realized by simply launching larger kernels. By breaking work into finer-grained, independent streams, we can fill the gaps left by memory latency, kernel launch overhead, and underutilized SMs. The future of DNN acceleration lies in software that thinks in terms of concurrent, overlapping operations rather than monolithic, synchronous blocks. Stream-assisted asynchronous multiplication is a step in that direction.

For practitioners, we recommend starting with a profiler to identify opportunities: measure whether your GEMMs are memory-bound or compute-bound, whether kernel launch overhead is significant, and whether data transfers are a bottleneck. Then experiment with splitting along batch or head dimensions, and use CUDA Graphs to amortize launch overhead. The results can be striking, bringing us closer to the dream of seamless, low-latency deep learning.

This article was written with the intent to provide a comprehensive, hands-on guide to GPU stream-assisted multiplications. All code snippets are illustrative and may require adaptation to your specific environment. Remember, the best optimizations come from understanding your workload and measuring religiously.