Suma de matrices (método 2)

Código para sumar una matriz en CUDA

#include "stdio.h"
#define columnas 300
#define filas 300
__global__ void add(int *a, int *b, int *c) {

int x = blockIdx.x * blockDim.x + threadIdx.x;
int y = blockIdx.y * blockDim.y + threadIdx.y;
int i = (columnas * y) + x;

c[i] = a[i] + b[i];
}

int main() {
int cont = 0;
int i, j;
// matrices en host
int a[filas][columnas], b[filas][columnas], c[filas][columnas];

// matrices en GPGPU
int *dev_a, *dev_b, *dev_c;

cudaMalloc((void **) &dev_a, filas * columnas * sizeof(int));
cudaMalloc((void **) &dev_b, filas * columnas * sizeof(int));
cudaMalloc((void **) &dev_c, filas * columnas * sizeof(int));

/* inicializando variables con datos foo*/
for (i = 0; i < filas; i++) {
cont = 0;
for (j = 0; j < columnas; j++) {
a[i][j] = cont;
b[i][j] = cont;
cont++;
}
}
cudaMemcpy(dev_a, a, filas * columnas * sizeof(int),cudaMemcpyHostToDevice);
cudaMemcpy(dev_b, b, filas * columnas * sizeof(int),cudaMemcpyHostToDevice);

// definiendo grid
dim3 grid(columnas, filas);

// grid del tamaño de la matriz, con un thread por bloque
add<<<grid, 1>>>(dev_a, dev_b, dev_c);

cudaMemcpy(c, dev_c, filas * columnas * sizeof(int), cudaMemcpyDeviceToHost);

// imprimiendo
for (int y = 0; y < filas; y++)
{
for (int x = 0; x < columnas; x++) {
printf("[%d][%d]=%d ", y, x, c[y][x]);
}
printf("\n");
}
return 0;
}

Notas importantes del código:

• El código anterior realiza la suma de una matriz cuadrada de 300 x 300.
• La diferencia principal con el código del método 1, es que el grid de bloques mapea la matriz definida, y el número de threads por bloque es sólo de 1.
• Los demás cálculos son equivalentes.

Suma de Matrices (método 1)

Código para sumar una matriz en CUDA

#include <stdio.h>
#include <stdlib.h>
#include <cuda_runtime.h>
#include <cuda.h>
#include <math.h>

#define T 10 // max threads x bloque
#define N 300


__global__ void sumaMatrices(int *m1, int *m2, int *m3) {

int col = blockIdx.x * blockDim.x + threadIdx.x;
int fil = blockIdx.y * blockDim.y + threadIdx.y;

int indice = fil * N + col;


if (col < N && fil < N) {
// debido a que en los últimos bloques no se realizan todos los threads
m3[indice] = m1[indice] + m2[indice];
}
}

int main(int argc, char** argv) {

int m1[N][N];
int m2[N][N];
int m3[N][N];
int i, j;
int c = 0;

/* inicializando variables con datos foo*/
for (i = 0; i < N; i++) {
c = 0;
for (j = 0; j < N; j++) {
m1[i][j] = c;
m2[i][j] = c;
c++;
}
}

int *dm1, *dm2, *dm3;

cudaMalloc((void**) &dm1, N * N * sizeof(int));
cudaMalloc((void**) &dm2, N * N * sizeof(int));
cudaMalloc((void**) &dm3, N * N * sizeof(int));

// copiando memoria a la GPGPU
cudaMemcpy(dm1, m1, N * N * sizeof(int), cudaMemcpyHostToDevice);
cudaMemcpy(dm2, m2, N * N * sizeof(int), cudaMemcpyHostToDevice);

// cada bloque en dimensión x y y tendrá un tamaño de T Threads
dim3 dimThreadsBloque(T, T);

// Calculando el número de bloques en 1D
float BFloat = (float) N / (float) T;
int B = (int) ceil(BFloat);

// El grid tendrá B número de bloques en x y y
dim3 dimBloques(B, B);

// Llamando a ejecutar el kernel
sumaMatrices<<<dimBloques, dimThreadsBloque>>>(dm1, dm2, dm3);

// copiando el resultado a la memoria Host
cudaMemcpy(m3, dm3, N * N * sizeof(int), cudaMemcpyDeviceToHost);
//cudaMemcpy(m2, dm2, N * N * sizeof(int), cudaMemcpyDeviceToHost);

cudaFree(dm1);
cudaFree(dm2);
cudaFree(dm3);

printf("\n");

for (i = 0; i < N; i++) {
for (j = 0; j < N; j++) {
printf(" [%d,%d]=%d", i, j, m3[i][j]);

}
printf("\n\n");

}
printf("\nB = %d", B);
printf("\n%d, %d",dimBloques.x, dimBloques.y);
printf("\n%d, %d",dimThreadsBloque.x, dimThreadsBloque.y);
return (EXIT_SUCCESS);
}

Notas importantes:

• Realiza la suma de una matriz cuadrada de enteros de 300 x 300 (90,000 elementos en total).
• Define cada bloque con Threads de 2D con valores de 10, es decir 100 threads por bloque.
• Define un grid de bloques en 2D. Para que alcancen los elementos de la matriz el grid será de 30 x 30 (900 bloques).
• Los elementos de la matriz se toman como un único vector. C se encarga de hacerlo parecer como una matriz. Por lo tanto es necesario calcular el desplazamiento que se genera al aumentar la fila del elemento que se esté calculando:

int col = blockIdx.x * blockDim.x + threadIdx.x;
int fil = blockIdx.y * blockDim.y + threadIdx.y;

int indice = fil * N + col;

SumaVectores.cu

El siguiente código en CUDA-C realiza una suma básica de dos vectores, y almacena el resultado en un tercer vector.

#include <stdio.h>
#include <stdlib.h>
#include <cuda_runtime.h>
#include <cuda.h>
#include <math.h>

#define T 1024 // max threads x bloque
#define N 100000

__global__ void sumaVector(int *v1, int *v2, int *v3) {

int tid = blockIdx.x * blockDim.x + threadIdx.x;

if (tid < N) {
// debido a que en el último bloque no se realizan todos los threads
v3[tid] = v1[tid] + v2[tid];
}
}

int main(int argc, char** argv) {

int v1[N];
int v2[N];
int v3[N];
int i;

/* inicializando variables con datos foo*/
for (i = 0; i < N; i++) {
v1[i] = i;
v2[i] = i;
v3[i] = 0;
}

for (i = 0; i < N; i++) {
printf(" %d %d \t", v1[i], v2[i]);
if (i % 20 == 0)
printf("\n");

}

int *dv1, *dv2, *dv3;

cudaMalloc((void**) &dv1, N * sizeof(int));
cudaMalloc((void**) &dv2, N * sizeof(int));
cudaMalloc((void**) &dv3, N * sizeof(int));

/* copiando memoria a la GPGPU*/
cudaMemcpy(dv1, v1, N * sizeof(int), cudaMemcpyHostToDevice);
cudaMemcpy(dv2, v2, N * sizeof(int), cudaMemcpyHostToDevice);

/* Calculando el número de bloques*/
float BFloat = (float) N / (float) T;
int B = (int) ceil(BFloat);

/* Llamando a ejecutar el kernel */
sumaVector<<<B, T>>>(dv1, dv2, dv3);

/* copiando el resultado a la memoria Host */
cudaMemcpy(v3, dv3, N * sizeof(int), cudaMemcpyDeviceToHost);

cudaFree(dv1);
cudaFree(dv2);
cudaFree(dv3);

printf("\n");

for (i = 0; i < N; i++) {
printf(" %d=%d", i, v3[i]);
if (i % 40 == 0)
printf("\n");

}

return (EXIT_SUCCESS);
}

Características importantes a tomar en cuenta

• Este cálculo se realiza en dimensión 1D.
• Debido a que el número de threads por dimensión en un bloque es 1024, se define ese valor como T. 
• N representa el tamaño del vector, en este caso se definió en 100 mil.
• Debido al tamaño de N, es necesario dividir el trabajo en diferentes bloques, por lo que se realiza el cálculo para la variable B (en este caso seria 100000/1024 =97.65, convertido a 98 bloques).
• Cada bloque es reservado completamente, aunque parte del último bloque no se utilice.

Arquitectura de una GPGPU NVIDIA Tesla C2075

Número máximo de threads por bloque: 1024

Tamaño máximo de las dimensiones x, y, z de un bloque de threads: 1024 x 1024 x 64

Tamaño máximo de cada dimensión del grid de bloques de threads: 65535 x 65535 x 65535

NVIDIA Nsight Eclipse Edition

Nsight Eclipse es un IDE basado en Eclipse que permite editar, construir y debuggear aplicaciones en CUDA-C. Viene incluido por default en la versión 5.5 de Cuda Toolkit.

Sólo es necesario escribir en la consola

$nsight

Y la aplicación se iniciará:

Características de NVIDIA Nsight Eclipse Edition

Instalar Fedora 18 & CUDA 5.5

A continuación se describen los pasos necesarios para instalar CUDA 5.5 en Fedora 18 kernel version 3.6.10-4.fc18.x86_64

Instalar paquetes necesarios

sudo yum install kernel-devel
sudo yum install gcc-c++
yum update audit

Repositorio CUDA

Descargar el repositorio para Fedora 18 desde  CUDA Donwload e instalarlo en una terminal:

sudo rpm -Uhv cuda-repo-fedora18-5.5-0.x86_64.rpm

Instalar el controlador propietario de Video

El controlador instalado por default "nouveau" es incompatible con el Toolkit de CUDA, por lo que es necesario reemplazarlo:

sudo yum remove xorg-x11-drv-nouveau
sudo yum install nvidia-settings nvidia-kmod xorg-x11-drv-nvidia

Para prevenir que el controlador nouveau se active accidentalmente después, se debe editar el archivo /etc/default/grub. En la línea GRUB_CMDLINE_LINUX_DEFAULT añadir:

"rdblacklist=nouvear nouveau.modeset=0"

Reiniciar el sistema

Ahora ya es posible reiniciar el sistema operativo, debiendo trabajar correctamente.

Instalar CUDA Toolkit

En una terminal instalar CUDA Toolkit mediante:

sudo yum install cuda

* Descargará aproximadamente 700 Mb de datos en la instalación

Añadir variables de ambiente

Editar el archivo .bashrc del home, añadiendo:

export CUDA_HOME=/usr/local/cuda
export LD_LIBRARY_PATH=${CUDA_HOME}/lib64

PATH=${CUDA_HOME]/bin:${PATH}
export PATH

Comprobar la instalación

La instalación debe estar completa, para comprobar es posible realizar un ejemplo básico de Hello world y compilarlo.

O también dentro de /usr/local/cuda/samples compilar usando Make, y correr /usr/local/cuda/samples/bin/x86_64/linux/release/.deviceQuery.

Definiendo las estructuras Grid/Block para un kernel

Para llamar a ejecutar un kernel es necesario proveer una configuración de ejecución, esto es, las dimensiones del grid y del bloque con que se ejecutará dicho kernel. Esta información se proporciona mediante dos estructuras claves:

• Número de bloques en cada dimensión
• Threads por bloque en cada dimensión

Sintaxis de llamada a un kernel:

      nombreKernel <<< B, T >>> (arg1, arg2, ... argN);

donde:

• B es una estructura que define el número de bloques en el grid para cada dimensión (1D o 2D)
• T es una estructura que define el número de threads en un bloque para cada dimensión (1D, 2D o 3D).

Si se desea definir una estructura de 1D, se puede usar un entero para B y T en:

       nombreKernel <<< B, T >>> (arg1, arg2, ... argN);

donde:

• B es un entero que define un grid de dimensión 1D
• T es un entero que define un bloque 1D de ese tamaño.

por ejemplo:

       nombreKernel <<< 1, 100 >>> (arg1, arg2, ... argN);