3

我不是一个有任何能力的程序员。只是有人对 CUDA 感到好奇，所以我正在做一些阅读。我遇到了一个使用推力做移动平均线的例子：

简单移动平均线推力示例

示例，例如，运行并且大部分正常工作。然而，从某种意义上说，它只进行一次移动平均操作，这是微不足道的。

我会怎么说352这些并行的移动平均操作，都在同一个数据流上操作？在我看来，程序流程可能是：

  1. 生成数据并将其发送到一个 CUDA 核心。（与现有代码相同，但认为长度为1000or10000而不是30）
  2. 将它从它所在的 CUDA 核心复制到351我的 GTX 465 中的所有其他 CUDA 核心
  3. 告诉每个 CUDA 核心要平均多少数据项。( 4, 5, 6,..., 352, 353, 354)
  4. 告诉设备在每个内核中并行运行平均值
  5. 读回每个核心的结果

我知道这段代码

// compute SMA using standard summation
simple_moving_average(data, w, averages);

让这一切发生，但我如何让 Thrust 并行执行其中的许多操作？

我在这里感兴趣的是股票数据之类的东西。如果我正在查看 GOOG 价格，我会使用所有内核将其放入 GPU 中并留在那里。然后，我可以自由地进行大量处理，而无需再加载数据，只需从每个核心读取结果。注意：我可能不想在所有内核中都使用 GOOG。有些核心可能是 GOOG，有些核心可能带有其他符号，但我稍后会介绍。我只是在想，如果每个核心都有足够的空间，我不想要全局内存中的库存数据。

我认为这对于 CUDA 和 Thrust 来说非常简单？

4

2 回答 2

3

以下是使用 arrayfire 执行此操作的可能方法：请注意，我不隶属于该库。
我很确定这也可以通过推力来完成，但我发现使用 arrayfire 更简单。如果图书馆是免费的，为什么我不能用它来代替推力？

在 arrayfire 中，您可以使用矩阵并行运行多个 SMA 操作：

unsigned n_SMAs = 1000;   // # of SMA indicators to evaluate 
unsigned len = 2000;      // # of stock prices per indicator
unsigned w = 6; // window size

// generate stock prices: [0..10] 
af::array data = af::randu(n_SMAs, len) * 10;

// compute inclusive prefix sums along colums of the matrix
af::array s = af::accum(data, 1);

// compute the average
af::array avg = (s.cols(w, af::end) - s.cols(0, af::end - w)) / w;
af::eval(avg);

std::cout << avg.dims() << "\n" << avg << "\n";

让我知道这是否是您正在寻找的。这就是我理解您的问题的方式：并行计算多个 SMA 指标

于 2012-09-14T11:49:11.887 回答
1

我的理解是您对以下两种情况感兴趣：

  1. 您有很长的项目序列，并且您希望通过对不同数量的项目进行平均来计算一定数量的平均值，即，对移动平均窗口使用不同的长度。这是我从你原来的问题中理解的。
  2. 您有一系列序列，连续存储在内存中，并且您希望使用大小为 的固定平均窗口并行平均它们2 * RADIUS + 1。这就是@asm 提出的 ArrayFire 代码的作用——你已经接受了。

而不是使用 CUDA Thrust，我认为编写自己的 CUDA 内核来执行上述操作会更容易。下面是一个完整的示例，其操作方式与@asm 提出的 ArrayFire 代码相同，因此涵盖了案例 #2。修改它以涵盖案例#1 将很简单。

#include <thrust/device_vector.h>

#define RADIUS        3
#define BLOCK_SIZE_X  8
#define BLOCK_SIZE_Y  8

/*******************/
/* iDivUp FUNCTION */
/*******************/
int iDivUp(int a, int b){ return ((a % b) != 0) ? (a / b + 1) : (a / b); }

/********************/
/* CUDA ERROR CHECK */
/********************/
#define gpuErrchk(ans) { gpuAssert((ans), __FILE__, __LINE__); }
inline void gpuAssert(cudaError_t code, const char *file, int line, bool abort=true)
{
   if (code != cudaSuccess) 
   {
      fprintf(stderr,"GPUassert: %s %s %d\n", cudaGetErrorString(code), file, line);
      if (abort) exit(code);
   }
}

/**********/
/* KERNEL */
/**********/
__global__ void moving_average(unsigned int *in, unsigned int *out, unsigned int M, unsigned int N) {

    __shared__ unsigned int temp[BLOCK_SIZE_Y][BLOCK_SIZE_X + 2 * RADIUS];

    unsigned int gindexx = threadIdx.x + blockIdx.x * blockDim.x;
    unsigned int gindexy = threadIdx.y + blockIdx.y * blockDim.y;
    unsigned int gindex  = gindexy * N + gindexx;

    unsigned int lindexx = threadIdx.x + RADIUS;
    unsigned int lindexy = threadIdx.y;

    // --- Read input elements into shared memory
    temp[lindexy][lindexx] = ((gindexx < N)&&(gindexy < M))? in[gindex] : 0;
    if (threadIdx.x < RADIUS) {
        temp[lindexy][threadIdx.x] = ((gindexx >= RADIUS)&&(gindexx < (N + RADIUS))&&(gindexy < M)) ? in[gindex - RADIUS] : 0;
        temp[lindexy][threadIdx.x + (RADIUS + min(BLOCK_SIZE_X, N - blockIdx.x * BLOCK_SIZE_X))] = (((gindexx + min(BLOCK_SIZE_X, N - blockIdx.x * BLOCK_SIZE_X)) < N)&&(gindexy < M))? in[gindexy * N + gindexx + min(BLOCK_SIZE_X, N - blockIdx.x * BLOCK_SIZE_X)] : 0;
        if ((threadIdx.y == 0)&&(gindexy < M)&&((gindexx + BLOCK_SIZE_X) < N)&&(gindexy < M)) printf("Inside 2 - tidx = %i; bidx = %i; tidy = %i; bidy = %i; lindexx = %i; temp = %i\n", threadIdx.x, blockIdx.x, threadIdx.y, blockIdx.y, threadIdx.x + (RADIUS + BLOCK_SIZE_X), temp[lindexy][threadIdx.x + (RADIUS + BLOCK_SIZE_X)]);
    }

    __syncthreads();

    // --- Apply the stencil
    unsigned int result = 0;
    for (int offset = -RADIUS ; offset <= RADIUS ; offset++) {
        result += temp[lindexy][lindexx + offset];
    }

    // --- Store the result
    out[gindexy * N + gindexx] = result;
}

/********/
/* MAIN */
/********/
int main() {

    const unsigned int M        = 2;
    const unsigned int N        = 4 + 2 * RADIUS;

    const unsigned int constant = 3;

    thrust::device_vector<unsigned int> d_in(M * N, constant);
    thrust::device_vector<unsigned int> d_out(M * N);

    dim3 GridSize(iDivUp(N, BLOCK_SIZE_X), iDivUp(M, BLOCK_SIZE_Y));
    dim3 BlockSize(BLOCK_SIZE_X, BLOCK_SIZE_Y);
    moving_average<<<GridSize, BlockSize>>>(thrust::raw_pointer_cast(d_in.data()), thrust::raw_pointer_cast(d_out.data()), M, N);
    gpuErrchk(cudaPeekAtLastError());
    gpuErrchk(cudaDeviceSynchronize());

    thrust::host_vector<unsigned int> h_out = d_out;

    for (int j=0; j<M; j++) {
        for (int i=0; i<N; i++)
            printf("Element j = %i; i = %i; h_out = %i\n", j, i, h_out[N*j+i]);
    }

    return 0;

}
于 2014-11-15T15:42:20.850 回答