The CUDA Handbook ▸ Chapter 11. Streaming Workloads

11.1 Device Memory

If the input and output data is in device memory, optimizing a low-density computation such as SAXPY is a matter of optimizing the global memory access. Besides alignment and coalescing constraints that inform performance, CUDA kernels are sensitive to the number of blocks and threads per block. The globalRead, globalWrite, globalCopy() and globalCopy2() applications (in the memory/ subdirectory of the source code) generate reports for the bandwidths achieved for a variety of operand sizes, threadblock sizes, and loop unroll factors. A sample report generated by globalCopy2() (which follows a memory access pattern similar to SAXPY: two reads and one write per loop iteration) is given in Listing 11-3.

If we reference the globalCopy2.cu application from Chapter 5 (Listing 5-8), running it on a GeForce RTX 3060 gets us the output of Listing 11-3 for 4-byte operands. The top row (unroll factor of 1) corresponds to the naïve implementation (similar to Listing 11-2). On this Ampere-class GPU, the naïve kernel already delivers 308 GiB/s—more than 90% of the board’s theoretical peak bandwidth—so unrolling the loop makes little difference: the fastest configuration, an unroll factor of 4, reaches just 310 GiB/s.

Because the naïve copy already runs so close to the memory ceiling, the micro-optimizations that mattered on older hardware—hand-unrolling the loop, or coaxing the compiler with #pragma unroll—have little room left to work. Recasting the kernel to perform SAXPY rather than a plain copy does not change the picture either: Listing 11-4 gives that kernel, which is implemented in the stream1Device.cu application (cudahandbook/streaming/), and measured on its own it sustains about 306 GiB/s.

Operand size: 4 bytesInput size: 16M operandsUnroll	32	64	128	256	512	maxBW	maxThreads1	240.62	294.75	306.58	307.14	308.14	308.14	5122	296.61	308.86	307.10	308.39	308.25	308.86	643	304.43	308.75	306.34	307.86	308.32	308.75	644	309.65	310.34	305.88	307.53	307.76	310.34	645	309.59	305.59	307.89	308.40	307.86	309.59	326	308.62	309.66	306.89	307.06	306.79	309.66	647	309.57	307.77	306.83	306.42	305.06	309.57	328	309.35	307.23	307.43	306.40	307.52	309.35	329	308.65	304.30	307.55	306.21	305.86	308.65	3210	308.19	306.39	307.74	305.73	303.02	308.19	3211	307.64	306.77	306.42	304.63	299.44	307.64	3212	307.86	304.92	305.84	303.70	294.41	307.86	3213	307.20	304.89	307.01	306.51	0.00	307.20	3214	306.63	303.20	305.51	304.76	0.00	306.63	3215	305.42	302.45	304.81	301.91	0.00	305.42	3216	306.14	301.85	303.51	298.18	0.00	306.14	32

template<const int n> __device__ voidsaxpy_unrolled(     float *out,     const float *px,     const float *py,     size_t N,     float alpha ){    float x[n], y[n];    size_t i;    for ( i = n*blockIdx.x*blockDim.x+threadIdx.x;           i < N-n*blockDim.x*gridDim.x;           i += n*blockDim.x*gridDim.x ) {        for ( int j = 0; j < n; j++ ) {            size_t index = i+j*blockDim.x;            x[j] = px[index];            y[j] = py[index];        }        for ( int j = 0; j < n; j++ ) {            size_t index = i+j*blockDim.x;            out[index] = alpha*x[j]+y[j];        }    }    // to avoid the (index<N) conditional in the inner loop,     // we left off some work at the end    for ( int j = 0; j < n; j++ ) {        for ( int j = 0; j < n; j++ ) {            size_t index = i+j*blockDim.x;            if ( index<N ) {                x[j] = px[index];                y[j] = py[index];            }        }        for ( int j = 0; j < n; j++ ) {            size_t index = i+j*blockDim.x;            if ( index<N ) out[index] = alpha*x[j]+y[j];        }    }} __global__ voidsaxpyGPU( float *out, const float *px, const float *py, size_t N, float alpha ){    saxpy_unrolled<4>( out, px, py, N, alpha );}

The stream1Device.cu application reports the total wall clock time needed to transfer data from pageable system memory to device memory, operate on the data with the kernel of Listing 11-4, and transfer the data back. On a test system with an AMD Ryzen 7 7700X running Ubuntu on a GeForce RTX 3060, the output of this application is as follows:

Measuring times with 128M floats (use --N to specify number of Mfloats)
Memcpy( host->device ): 53.62 ms (20024.40 MB/s)
Kernel processing     : 4.90 ms (328632.09 MB/s)
Memcpy (device->host ): 29.26 ms (18348.57 MB/s)

Total time (wall clock): 87.80 ms (18343.83 MB/s)

The kernel takes a tiny amount of the overall execution time – about 6% of the wall clock time. The other 94% of time is spent transferring data to and from the GPU! For transfer-bound workloads like this one, if some or all of the data being operated on is in host memory, the best way to optimize the application is to improve CPU/GPU overlap and transfer performance.