Scan implementations often operate on predicates, truth values (0 or 1) computed by evaluating a condition. As mentioned at the beginning of the chapter, an exclusive scan of predicates implements stream compaction, a class of parallel problems where only the “interesting” elements of an input array are written to the output: where the predicate is 1 for an interesting element, its exclusive-scan value gives that element’s output index.
Because we already have a single-pass scan, compaction can be done in
a single pass too – which is how CUB’s DeviceSelect works. As an
example, let’s compact an array of int down to its odd values5. One thread block processes one tile
of b elements: it evaluates the predicate
(isOdd()) and inclusive-scans the resulting 0/1 flags in
shared memory. That scan gives each kept element its index among the
tile’s survivors, and its last element is the tile’s keep count. The
cooperative look-back of Section 13.3.4 then supplies the number of
elements kept by all earlier tiles – this tile’s base output index – and
each surviving element is scattered to
out[base + within-tile index]. The highest-numbered tile
writes the grand total.
The only differences from the decoupled-look-back scan are that we scan the 0/1 predicate flags rather than the input values, and scatter the survivors rather than writing a scanned value in place. Listing 13-17 gives the kernel.
template<class T>__host__ __device__ boolisOdd( T x ){ return x & 1;} template<class T>__global__ voidstreamCompact_odd_kernel( T *out, int *outCount, // total kept (written by the last tile) const T *in, volatile scanStatus *status, // one descriptor per tile (SCAN_X-initialized) uint32_t *tileCounter, // one global counter, 0-initialized uint32_t numTiles, size_t N ){ extern __shared__ int sPartials[]; // blockDim.x predicate flags __shared__ uint32_t s_tile; __shared__ int s_base; // number of elements kept by earlier tiles if ( threadIdx.x == 0 ) s_tile = atomicAdd( tileCounter, 1 ); __syncthreads(); const uint32_t tile = s_tile; const size_t gidx = (size_t) tile * blockDim.x + threadIdx.x; // // Evaluate the predicate, then inclusive-scan the 0/1 flags in shared // memory (Kogge-Stone). sPartials[blockDim.x-1] is the tile's keep count. // T value = (T) 0; int pred = 0; if ( gidx < N ) { value = in[gidx]; pred = isOdd( value ) ? 1 : 0; } sPartials[threadIdx.x] = pred; __syncthreads(); for ( int off = 1; off < blockDim.x; off <<= 1 ) { int add = ( threadIdx.x >= off ) ? sPartials[threadIdx.x - off] : 0; __syncthreads(); sPartials[threadIdx.x] += add; __syncthreads(); } const int aggregate = sPartials[blockDim.x - 1]; // // Cooperative look-back over the per-tile keep counts: s_base is the number // of elements kept by every earlier tile -- this tile's base output index. // scanCoopLookback<int>( status, tile, aggregate, s_base ); __syncthreads(); // // Scatter. sPartials[threadIdx.x] is the inclusive keep count, so // sPartials[threadIdx.x]-1 is this element's index among the tile's kept // elements; s_base offsets it into the global output. // if ( gidx < N && pred ) out[s_base + sPartials[threadIdx.x] - 1] = value; if ( tile == numTiles - 1 && threadIdx.x == 0 ) *outCount = s_base + aggregate;}streamCompact_odd_kernel() (source on GitHub)Listing 13-18 gives the host function; like the scan’s host code (Listing 13-15), it allocates the per-tile descriptor array and the tile counter, initializes them, and makes a single kernel launch.
template<class T>voidstreamCompact_odd( T *out, int *outCount, const T *in, size_t N, int b ){ cudaError_t status_cudart; scanStatus *gStatus = 0; uint32_t *tileCounter = 0; if ( N == 0 ) return; uint32_t numTiles = (uint32_t) ( ( N + b - 1 ) / b ); // // gStatus and tileCounter are transient per-call scratch. Allocate and // zero them stream-ordered on the default stream so repeated calls (e.g. // a timing loop) recycle the same pool memory instead of paying for a // synchronizing cudaMalloc/cudaFree on every invocation. // cuda(MallocAsync( &gStatus, numTiles * sizeof(scanStatus), 0 ) ); cuda(MemsetAsync( gStatus, 0, numTiles * sizeof(scanStatus), 0 ) ); // SCAN_X == 0 cuda(MallocAsync( &tileCounter, sizeof(uint32_t), 0 ) ); cuda(MemsetAsync( tileCounter, 0, sizeof(uint32_t), 0 ) ); streamCompact_odd_kernel<T><<<numTiles, b, b * sizeof(int)>>>( out, outCount, in, gStatus, tileCounter, numTiles, N ); Error_cudart: if ( gStatus ) cudaFreeAsync( gStatus, 0 ); if ( tileCounter ) cudaFreeAsync( tileCounter, 0 );}streamCompact_odd() host function. (source on GitHub)The code is easily modified to evaluate more-complicated predicates.↩︎