Объединить два ядра CUDA в одно

Question

Я использую CUDA для вычисления out = C(b(A(in))), где функции A и C являются сверточными, а b является поэлементной функцией.Вот пример:

#define N 1000

__device__ float b(float d_in){return min(d_in + 10.0f, 100.0f);}
__global__ void bA(float *d_in, float *d_out){
    const int x = threadIdx.x + blockIdx.x * blockDim.x;
    if (x >= N)  return;

    // replicate boundary
    int x_left  = max(x-1, 0); int x_right = min(x+1, N-1);

    d_out[x] = b( d_in[x_left] + d_in[x] + d_in[x_right] );
}
__global__ void C(float *d_in, float *d_out){
    const int x = threadIdx.x + blockIdx.x * blockDim.x;
    if (x >= N)  return;

    // replicate boundary        
    int x_left  = max(x-1, 0); int x_right = min(x+1, N-1);

    d_out[x] = d_in[x_left] + d_in[x] + d_in[x_right];
}
void myfunc(float *d_data, float *d_temp){
    dim3 threads(256);
    dim3 blocks( (N + threads.x - 1) / threads.x ); // divide up

    // kernels that I would like to merge into one:
    bA<<<blocks, threads>>>(d_data, d_temp);
    C <<<blocks, threads>>>(d_temp, d_data);
}

Для таких вычислений нужна дополнительная переменная d_temp, которая мне не нужна.Поэтому я хотел бы объединить эти ядра в одно, то есть одно ядро ​​для вычисления C(b(A(in))).

Одна трудность состоит в том, как я могу сохранить временные результаты из b(A(in)), а затем выполнить функцию свертки C()?Я пытался использовать разделяемую память, но потеря при загрузке временного результата b(A(in)) в разделяемую память.Например:

#define BLOCK_SIZE 32

__global__ void CbA(float *d_in, float *d_out){
    const int x = threadIdx.x + blockIdx.x * blockDim.x;
    if (x >= N)  return;

    // replicate boundary
    int x_left  = max(x-1, 0); int x_right = min(x+1, N-1);

    // temp result for b(A(in))
    float temp = b( d_in[x_left] + d_in[x] + d_in[x_right] );

    // shared memory for convolution (stencil size of 3)
    __shared__ float shmem[BLOCK_SIZE+2];

    // load center part to shared memory
    shmem[threadIdx.x+1] = temp;

    // but how to load boundary parts from temp to shmem?
    // ...

    __syncthreads();

    // perform function C()
    // ...
}

Любые советы или советы приветствуются.

Answer 1

Сначала комментарий о

// load center part to shared memory
shmem[threadIdx.x+1] = temp;

Я бы назвал это сохранением в общей памяти ...

Кроме того, несколько идей:

Используйте первый и последний поток в блоке, чтобы вычислять только b(A(in))

Конечно, вам нужно учитывать это при расчете x (const int x = threadIdx.x + blockIdx.x * (blockDim.x-2);) и вызывать ваше ядро ​​с большим количеством потоков./blocks.
При выполнении команды C() у вас будет по два потока на каждый простаивающий блок.Но это не должно иметь большого влияния.
Вот ядро.Это легче понять, если попытаться визуализировать последовательность вычислений.

__global__ void CbA(float *d_in, float *d_out)
{
  const int x = threadIdx.x + blockIdx.x * (blockDim.x - 2);
  if (x >= N) return;
  int x_left  = max(x-1, 0); int x_right = min(x+1, N-1);
  float temp = b( d_in[x_left] + d_in[x] + d_in[x_right] );
  __shared__ float shmem[BLOCK_SIZE]; // = 256
  shmem[threadIdx.x] = temp;
  __syncthreads();
  if (threadIdx.x > 0 && threadIdx.x < blockDim.x-1)
    d_out[x-1] = shmem[threadIdx.x-1] + d_in[threadIdx.x] + d_in[threadIdx.x+1];
}

Пусть один поток в блоке также выполняет b(A()) для "граничных частей" блока

Но тогда вы будете использовать только 1 из 32 потоков для этого расчета для каждого блока.Наихудшим случаем будет отношение 1/32 для всего SM на время дополнительных вычислений.

...
// but how to load boundary parts from temp to shmem?
if (threadIdx.x == 0)
{
  {
    const int x = 0 + blockIdx.x * blockDim.x;
    int x_left = max(x-1, 0); int x_right = min(x+1, N-1);
    float temp = b( d_in[x_left] + d_in[x] + d_in[x_right] );
    shmem[0] = temp;
  }
  {
    const int x = blockDim.x-1 + blockIdx.x * blockDim.x;
    int x_left = max(x-1, 0); int x_right = min(x+1, N-1);
    float temp = b( d_in[x_left] + d_in[x] + d_in[x_right] );
    shmem[blockDim.x-1] = temp;
  }
}
// perform function C()
...

Избегайте использования общей памяти

(по крайней мере, в вашем упрощенном примере) Значение temp является результатом очень простого вычисления.Возможно, лучше всего рассчитать все значения, которые вам нужно выполнить C() в потоке локально в этом потоке.

__global__ void CbA(float *d_in, float *d_out)
{
  const int x = threadIdx.x + blockIdx.x * blockDim.x;
  if (x >= N)  return;

  float temp[3];
  for (int i(0); i < 3; ++i)
  {
    int x_left  = max(x-1-1+i, 0); int x_right = min(x+1-1+i, N-1);
    temp[i] = b( d_in[x_left] + d_in[x-1+i] + d_in[x_right] );
  }

  // perform function C()
  ...
}