At full occupancy, corrShared_kernel() is limited not by
arithmetic but by the latency of its shared-memory reads: each thread
runs a single dependent chain of accumulations, so it spends most of its
time waiting on loads rather than computing. Two changes, applied in
sequence, move and then address that bottleneck. Both are implemented in
the corrShared4Dp4a.cuh and
corrShared4Dp4aIlp.cuh header files, and together they
replace the scalar 4×-unrolled kernel of the first edition, which
measured no faster than the plain shared-memory version on modern
hardware.
The first change exposes instruction-level parallelism. Instead of computing one output column, each thread computes two, keeping an independent set of accumulators for each:
unsigned SumI[2], SumISq[2], SumIT[2];
for ( int c = 0; c < 2; c++ )
SumI[c] = SumISq[c] = SumIT[c] = 0;
...
for ( int j = 0; j < hTemplate; j++ ) {
for ( int i = 0; i < wTemplate/4; i++ ) {
for ( int c = 0; c < 2; c++ ) {
unsigned char I = /* column c's pixel */;
unsigned char T = /* template pixel */;
SumI[c] += I;
SumISq[c] += I*I;
SumIT[c] += I*T;
}
}
}
Because the two columns’ accumulation chains are independent, the compiler can interleave them: the second column’s loads are in flight while the first column’s arithmetic runs. This hides the shared-load latency and flips the kernel from latency-bound to arithmetic-bound—worth about 2.4× on a GeForce RTX 3060. Two columns is the measured sweet spot on this hardware; going wider costs registers and, past two, occupancy, causing throughput to fall off again.
Only once the kernel is arithmetic-bound does optimizing the
arithmetic improve performance. The inner loop reduces over four
adjacent bytes at a time, making the computation ideally suited to the
__dp4a integer dot-product instruction introduced with the
Pascal architecture (Section 8.2.6): in a single instruction, the four
8-bit lanes in each of the two 32-bit input operands are multiplied
separately, then added together and the sum of those four products is
returned in a 32-bit accumulator. Packing four image pixels into a
32-bit word I4 and four template pixels into
T4, the four scalar iterations collapse to three
instructions:
SumI[c] = __dp4a( I4, 0x01010101u, SumI[c] ); // sum of the four pixels
SumISq[c] = __dp4a( I4, I4, SumISq[c] ); // sum of squares
SumIT[c] = __dp4a( I4, T4, SumIT[c] ); // sum of image*template
Because the pixels are unsigned char, the unsigned
overload of __dp4a is required: under the signed
interpretation, any pixel above 127 would be read as negative and the
products would come out wrong. Packing the image pixels takes some
care—consecutive columns slide the template window by one byte, so each
thread’s four pixels begin at an unaligned offset and must be assembled
from the two aligned shared-memory words that straddle them with a
funnel shift (__funnelshift_r); the template, whose width
is a multiple of four, packs directly. Only SumIT is a
paired dot product, so image and template bytes must land in
matching lanes; the sums of pixels and of squares are order-independent,
which means a byte-order slip hides in SumIT alone and
survives any test that checks only the other two sums.
On its own, this DP4A rewrite (corrShared4Dp4a.cuh) is
no faster than the scalar shared kernel: with only a short dependent
chain per thread, the cheaper arithmetic simply hides in the same
latency shadow that limited the original—reducing instruction count does
not help when instructions are not the bottleneck. It is kept in the
sample as a teaching intermediate. Combined with the two-column ILP
(corrShared4Dp4aIlp.cuh), though, DP4A delivers about 6×
the throughput of corrShared, at 100% occupancy with no
register spilling. Neither optimization is worth much alone—DP4A needs
the parallelism to escape the latency shadow, and the parallelism needs
DP4A to win big once the kernel is throughput-bound—but together they
deliver about 7x faster performance.