Funciones __device__

En ejercicios anteriores hemos aprendido como usar la palabra reservada __global__ para marcar una función como código que el host puede llamar y que origina una invocación de un kernel paralelo en la GPGPU.

Con una función __global__, cada thread CUDA sigue su propio patrón de ejecución en forma serial. En CUDA, los kernels consisten en mayormente en código C/C++ que pueden ser muy rápidos.

Como programadores en paralelo, tendremos la necesidad de abstraer y encapsular el kernel en funciones. La palabra reservada __device__ nos permite marcar una función que es llamada desde threads ejecutándose en la GPGPU. Por ejemplo:

__device__ float mi_funcion_device(float x)
{
   return x + 1;
}

Aunque las funciones __device__ son similares a las funciones __global__ en que pueden ser ejecutadas por threads CUDA, de hecho se comportan más como funciones normales en C. A diferencia de funciones __global__, las funciones __device__ no pueden ser configuradas con <<B, T>>, y no están sujetas a ninguna restricción especial en los tipos de sus parámetros o sus resultados. El código del host no puede llamar a las funciones __device__ directamente; si se desea acceder a una función __device__, necesitamos escribir una función __global__ que la llame.

Como se podría esperar, las funciones __device__ pueden llamar a otras funciones __device__:

__device__ float mi_segunda_funcion_device(float y)
{
   return mi_funcion_device(y)/2;
}

Siempre y cuando no se llamen a si mismas:

__device__ int float mi_funcion_device_ilegal_recursiva(int x)
{
   if (x == 0) return 1;
   return x * mi_funcion_device_ilegal_recursiva(x-1);
}

El siguiente código muestra como se podría usar la funciones __device__ para empaquetar varios bits de código cuando se desarrolla un kernel CUDA.

#include <stdlib.h>
#include <stdio.h>

__device__ int get_global_index(void)
{
  return blockIdx.x * blockDim.x + threadIdx.x;
}

__device__ int get_constant(void)
{
  return 7;
}

__global__ void kernel1(int *array)
{
  int index = get_global_index();
  array[index] = get_constant();
}

__global__ void kernel2(int *array)
{
  int index = get_global_index();
  array[index] = get_global_index();
}

int main(void)
{
  int num_elements = 256;
  int num_bytes = num_elements * sizeof(int);

  int *device_array = 0;
  int *host_array = 0;

  // reserva memoria
  host_array = (int*)malloc(num_bytes);
  cudaMalloc((void**)&device_array, num_bytes);

  int block_size = 128;
  int grid_size = num_elements / block_size;

  // lanza kernel1 e inspecciona sus resultados
  kernel1<<<grid_size,block_size>>>(device_array);
  cudaMemcpy(host_array, device_array, num_bytes, cudaMemcpyDeviceToHost);

  printf("resultados de kernel1:\n");
  for(int i = 0; i < num_elements; ++i)
  {
    printf("%d ", host_array[i]);
  }
  printf("\n\n");

  // lanza kernel2 e inspecciona sus resultados
  kernel2<<<grid_size,block_size>>>(device_array);
  cudaMemcpy(host_array, device_array, num_bytes, cudaMemcpyDeviceToHost);

  printf("resultados de kernel2:\n");
  for(int i = 0; i < num_elements; ++i)
  {
    printf("%d ", host_array[i]);
  }
  printf("\n\n");

  // liberar memoria
  free(host_array);
  cudaFree(device_array);
  return 0;
}

Jugando con dimensiones en grids y bloques

En CUDA es posible tener estructura de 1D, 2D y 3D para bloques y grids de threads, lo cual se explicó en Identificación de un thread. Para aclarar la forma de utilizar estas posibles configuraciones se resume en lo siguiente, usando configuraciones 1D y 2D para bloques y 1D, 2D y 3D para threads:

1. Array 1D de bloques donde cada bloque tiene un array 1D de threads.

UniqueBlockIndex = blockIdx.x;
UniqueThreadIndex = blockIdx.x * blockDim.x + threadIdx.x;
   
2. Array 1D de bloques donde cada bloque tiene un array 2D de threads.


UniqueBlockIndex = blockIdx.x;
UniqueThreadIndex = blockIdx.x * blockDim.x * blockDim.y + threadIdx.y * blockDim.x + threadIdx.x;
                      
3. Array 1D de bloques donde cada bloque tiene un array 3D de threads.


UniqueBlockIndex = blockIdx.x;
UniqueThreadIndex = blockIdx.x * blockDim.x * blockDim.y + threadIdx.y * blockDim.x + threadIdx.x;
                      
4. Array 2D de bloques donde cada bloque tiene un array 1D de threads.


UniqueBlockIndex = blockIdx.y * gridDim.x + blockIdx.x;
UniqueThreadIndex = UniqueBlockIndex * blockDim.x + threadIdx.x;
                                                         
5. Array 2D de bloques donde cada bloque tiene un array 2D de threads.


UniqueBlockIndex = blockIdx.y * gridDim.x + blockIdx.x;
UniqueThreadIndex = UniqueBlockIndex * blockDim.y * blockDim.x + threadIdx.y * blockDim.x + threadIdx.x;
                      
6. Array 2D de bloques donde cada bloque tiene un array 3D de threads.


UniqueBlockIndex = blockIdx.y * gridDim.z * + blockIdx.x; UniqueThreadIndex = UniqueBlockIndex * blockDim.z * blockDim.y * blockDim.x + threadIdx.z * blockDim.y * blockDim.z + threadIdz.y * blockDim.x + threadIdx.x;