cuda - GPU 上的广义滑动窗口计算

Question

这是一些在两个 3D 矩阵 X 和 Y 上实现滑动窗口计算的 Python 代码。

import numpy

def sliding_dot( X,Y ) :

    assert X.ndim == Y.ndim == 3
    iw,ih,id = X.shape
    fw,fh,fd = Y.shape

    assert id == fd
    assert fw < iw and fh < ih

    ow,oh = iw-fw+1,ih-fh+1
    out = numpy.zeros( [ow,oh] )

    for x in xrange(ow) :
        for y in xrange(oh) :
            window = X[x:x+fw,y:y+fh,:]
            out[x,y] = numpy.dot( window.flatten(),Y.flatten() )

    return out

#################    

A_dims = (640,480,32)
B_dims = (6,6,32)

A = numpy.random.rand(*A_dims)
B = numpy.random.rand(*B_dims)

sliding_dot(A,B)

一般来说，Y 在第一和第二维度上总是比 X 小得多，但在第三维度上它们是相等的。

请注意，我们可以将 numpy.dot() 替换为 Y 和窗口的任何函数。这与卷积有点不同，因为 Y 仅沿 X 的第一和第二维度滑动。我正在寻找一种有效的策略来使用 CUDA 有效地实现这种滑动窗口计算。有人想给我一些方向吗？干杯!

更新 :您可以在下面的回答中观看我在其他用户的帮助下完成优化过程的工作。

Answer 1

在像 CUDA 这样的架构中，尝试设计一个可以适应几乎任何您可能想要的操作的“通用”实现将是一个巨大的权衡。对于您的具体点积示例，这是一个典型的归约操作，这是一个非常有用的实现:

__constant__ int ldaX[3];
__constant__ int ldaY[3];
__constant__ int dimX[3];
__constant__ int dimY[3];

template<typename real,int blocksize>
__global__ void sliding_k(const real *X, const real *Y, real *out)
{
    __shared__ volatile real buffer[blocksize];

    int tid = threadIdx.x;
    int gid = blockIdx.x * gridDim.y + blockIdx.y;

    real value = (real)0;
    int xpos = (blockIdx.y * ldaX[2]) + (blockIdx.x * ldaX[1]);
    int ypos = 0;
    for(int i=0; i<dimY[0]; i++) {
        for(int jk=tid; jk<ldaY[1]; jk+=blocksize) {
            value += X[xpos+jk] * Y[ypos+jk];
        }
        xpos += ldaX[1];
        ypos += ldaY[1];
    }

    buffer[tid] = value;
    __syncthreads();

# pragma unroll
    for(int i=(tid+32); ((tid<32)&&(i<blocksize)); i+=32)
        buffer[tid] += buffer[i];

    if (tid < 16) buffer[tid] += buffer[tid + 16];
    if (tid < 8)  buffer[tid] += buffer[tid + 8];
    if (tid < 4)  buffer[tid] += buffer[tid + 4];
    if (tid < 2)  buffer[tid] += buffer[tid + 2];
    if (tid == 0) out[gid] = buffer[0] + buffer[1];
}

您可以将您喜欢的任何类型的归约运算符替换为点积使用的浮点乘加/求和运算，并且代码应该可以正常工作。每个窗口计算由单个 block 执行。有足够的并行工作来证明在这个窗口大小下每个窗口一个 block 是合理的。这允许合并的全局内存访问，并且在 Fermi 卡上，大量的 L1 缓存命中。

这里我只在代码中建立了一个假设，即源数组和窗口数组的第三维是相等的。这允许内部两个循环“融合”成一个操作，因为它们共享公共(public)内存布局。使用引用代码的改进版本在 Python 中运行测试工具，主机代码用 PyCUDA 编写，我得到以下信息:

In [15]: %timeit -n3 -r3 out2=sliding_cuda(A,B)
3 loops, best of 3: 49.8 ms per loop

In [16]: %timeit -n3 -r3 out=sliding_dot(A,B)
3 loops, best of 3: 2.18 s per loop

In [17]: (numpy.abs(out2-out)/numpy.abs(out)).max()
Out[17]: 4.2921323635558404e-15

当在 3GHz Phenom II 和 GTX470 上运行时，在 635x475 2D 网格上使用 64 个线程 block ——即。使用可分页的主机内存分配，包括模块加载、设置和内存传输，速度提高了大约 50 倍。在不包括内存传输和设置开销的情况下，内核本身比 Python 快大约 100 倍。请注意，这是一个 double 版本 - Python 默认使用 double 浮点运算。