Desarrollo Cuda

Una ordenación matricial en GPU

2012-06-03T17:41:00.000-07:00

Hola a todos,

hoy quiero ilustraros un método de ordenación que implementé para la ordenación de vectores en GPU, la verdad es que primero lo ideé y como suele pasar en estas cosas luego busqué en internet y efectivamente el método ya existía y para mi sorpresa desde hace muchos años, en este caso aplicado a lás máquinas matriciales.

La solución de la que hablo está basada en el método de ordenación paralelo Parallel Neighbor-Sort propuesto por Nico Habermann(1), tambien conocido como odd-even sort, se diseño inicialmente para computación en arrays de procesadores a principios de los 70. Se trata realizar comparaciones de los elementos del vector dos a dos, siendo estos vecinos, realizando paralelamente Elementos/2 comparaciones en cada ejecución del kernel, las llamadas se intercalarán variando la paridad del elemento comparado, repitiendo este proceso hasta que no se produzca ningún cambio, esto determinará el final de la ordenación.

Es un método sencillo pero muy eficaz en sistemas multiprocesador con capacidades de paralelismo.

A continuación os ilustro el kernel para la ordenación en CUDA, está ordenando 3 vectoresunas reglas especificas demás complejo del que ya os hablaré, sirva ahora de ejemplo:

__global__ void ordenar_especies(int par,int *cambios,int *contar, int M, int D, float *Puntos, int *Nivel, int *Orden, float *Valores) {

__shared__ int cuenta;

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

int elem = ((3*M*3)+M)/2;

float puntoA, puntoB, aux = 0.0f;

int ordenA, ordenB = 0;

int nivelA, nivelB = 0;

int cambiar_puntos = 0;

if (tid == 0){

cambios = 0;

if (contar[0] == 0){

cuenta = 0;

}

tid = (tid*2)+par;

if (tid < elem) {

puntoA = tex1Dfetch(valoresTex,tid);

puntoB = tex1Dfetch(valoresTex,tid+1);

ordenA = tex1Dfetch(ordenTex,tid);

ordenB = tex1Dfetch(ordenTex,tid+1);

nivelA = tex1Dfetch(nivelTex,tid);

nivelB = tex1Dfetch(nivelTex,tid+1);

if (ordenB > 0) {

if (ordenA == 0) {

cambiar_puntos = 1;

} else {

if (nivelB > nivelA){

cambiar_puntos = 1;

} else if ((nivelB == nivelA) && (puntoB > puntoA)){

cambiar_puntos = 1;

}

if (cambiar_puntos == 1) {

for (int x=0;x<D;++x){

aux = Puntos[tid+(x*elem)];

Puntos[tid+(x*elem)] = Puntos[tid+1+(x*elem)];

Puntos[tid+1+(x*elem)] = aux;

}

Valores[tid] = puntoB;

Valores[tid+1] = puntoA;

Orden[tid] = ordenB;

Orden[tid+1] = ordenA;

Nivel[tid] = nivelB;

Nivel[tid+1] = nivelA;

cambios[0] = 1;

}

if (contar[0] == 0){

if ((ordenA == 0) && (ordenB == 0)) {

__syncthreads();

cuenta = cuenta + 2;

}

else if ((ordenA == 0) || (ordenB == 0)) {

__syncthreads();

cuenta = cuenta + 1;

}

if (tid == (elem-1)){ //el ultimo que ponga la cuenta de su bloque

if (contar[0] == 0){

__syncthreads();

contar[blockIdx.x] = elem-cuenta;

}

Hasta la próxima!!
____________________
(1)Parallel Neighbor-Sort (or the Glory of the Induction Principle). 5285 Port Royal Rd Sprigfield VA : National Technical Information Service, US Department of Commerce, 1972. AD0759248.

Bienvenido Mr. Kepler

2012-05-12T16:35:00.000-07:00

En esta entrega compartiré con vosotros la nueva noticia: por fin tengo disponible en la plataforma del departamento un dispositivo con arquitectuta Kepler.

Según informa nVidia:

Kepler es la arquitectura de GPU más avanzada de nuestra historia. Optimizada para DirectX 11, Kepler no es solo rápida, también es supereficiente.

cSegún parece se Consigue 3 veces más velocidad de procesamiento con NVIDIA® Kepler, aseguran que es la arquitectura de alta computación (HPC) más rápida y eficiente del mundo. Está dotada de tecnología y funciones de computación avanzadas, puede utilizarse en una extensa variedad de aplicaciones de cálculo científico y pone los sistemas de computación híbridos al alcance de mayor número de programadores e investigadores.

Sus principales novedades son:

SMXProporciona mayor velocidad de procesamiento y eficiencia gracias a su innovador multiprocesador de streaming, que permite dedicar más espacio a los núcleos de procesamiento que a la lógica de control.

GPU dinámica con Kepler
Simplifica la programación en la GPU ya que facilita la aceleración de bucles anidados paralelos, lo que significa que una GPU puede iniciar nuevos subprocesos de forma dinámica por sí misma, sin necesidad de volver a la CPU.

Hyper-Q de KeplerReduce el tiempo de inactividad de la CPU al permitir que múltiples núcleos de ésta utilicen una misma GPU Kepler, lo que mejora drásticamente la programabilidad y la eficiencia.

Estas tres mejoras en la arqutectura prometen ser un salto cualitativo y cuantotativo en la aplicación de GPGPU.

Yo por el momento cuento con una modesta GeForce GTX 660, las pruebas que he realizado en principio no me han dejado muy impresionado. Supongo que CUDA 5 permitirá exprimir a la perfección esta arquitectura.

Ya os contaré como se comporta esta nueva arquitectura con CUAD 5

Hasta entoces, feclices desarrollos!!!

Generación de números aleatorios en CUDA

2012-05-04T13:45:00.000-07:00

Uno de los caballos de batalla más empleado en sistemas de optimización, minería de datos, sistemas expertos probabilísticos y simulación, es la utilización de números aleatorios, los números aleatorios nos permiten generar sobre la marcha poblaciones dentro de un espacio muestral controlado, que nos van a permitir, aplicando ciertas heurísticas, realizar inferencias y/o llegar a resultados. Normalmente estos métodos están clasificados dentro de lo que se conoce como "métodos Montecarlo" o basados en Montecarlo.

Por la intensa aplicación de estos métodos, se hace necesario disponer de mecanismos para la generación de números aleatorios que sea eficiente y con poco costo de procesamiento.

Existen varias soluciones para abordar el problema de la generación de números pseudoaleatorios con CUDA. Así, la librería CURAND provee al desarrollador de un importante surtido de funciones de alto nivel para la obtención de números aleatorios.

La librería CURAND contiene dos APIs: una para uso en código host, y la otra para uso en código de dispositivo.

Host API: para utilizar este API es necesario incluir curand.h en la lista de librerías en nuestro programa. Esta librería oculta al programador toda la problemática relacionada con el estado del generador de números aleatorios y su progresión de estado dentro del kernel de dispositivo. El funcionamiento básico es el siguiente:

1º) Se realiza la inicialización del estado del generador, para ello se utiliza la función curandCreateGenerator().

2º) Se especifican las opciones del generador, indicando además la semilla a utilizar. CURAND ofrece 2 tipos de generadores: CURAND_RNG_XORWOW, es un generador de números pseudoaleatorios que ha sido implementado usando el algoritmo XORWOW (32) introducido por George Marsaglia es una implementación de un generador del tipo xorshif también propuesto por él. Es un algoritmo muy rápido, calcula el siguiente elemento de la secuencia aplicando repetidamente un OR exclusivo sobre el resultado de aplicar un desplazamiento de bit al propio número. El código en lenguaje C publicado por Marsaglia es el siguiente:

unsigned long xor128(){

static unsigned long x=123456789,y=362436069,z=521288629,w=88675123;

unsigned long t;

t=(xˆ(x<<11));x=y;y=z;z=w;

return( w=(wˆ(w>>19))ˆ(tˆ(t>>8)) ) }

El Segundo tipo de generador que ofrece CURAND es un tipo de generador números quasialeatorios: CURAND_RNG_SOBOL32, que permite generar secuencias de 32-bit y hasta 20.000 dimensiones. Es un algoritmo algo más complejo ideado por I.M. Sobol en 1967, que permite crear secuencias uniformes de baja discrepancia utilizando muchas dimensiones.

Una vez seleccionado el tipo de RNG[1] realizaremos la llamada a la función de la librería: curandSet-PseudoRandomGeneratorSeed() indicando la semilla y tipo de generador deseado.

3º) Reservar el espacio de memoria en el dispositivo.

4º) Invocar una función de generación de las disponibles pasando como parámetro el puntero a la memoria del dispositivo donde se generará la secuencia. Las distintas funciones de generación son:

- curandSetPseudoRandomGeneratorSeed(), permite generar números pseudoaleatorios y quasialeatorios cada elemento es de 32 bits, sin signo, donde todos sus bits son aleatorios.

- curandGenerateUniform(), genera una secuencia uniforme-mente distribuida de números en coma flotante comprendidos entre (0.0, 1.0].

- curandGenerateNormal(), permite obtener una secuencia de distribución normal, dados los valores de la media y la desviación estándar.

Para obtener secuencias con valores utilizando doble precisión se cuenta con las funciones homologas a las dos últimas: curandGenerateUniformDouble(), para la distribución uniforme y curandGenerateNormaDouble() para la normal.

5º) Una vez terminados los procesos deseados, liberamos los recursos solicitados mediante el uso de la instrucción curandDestroyGenerator().

Device API: para utilizar esta modalidad, deberemos incluir la librería curand_kernel.h en nuestro programa. Como es evidente, las funciones aquí definidas se invocarán en los kernels de dispositivo o en funciones de dispositivo, diseñadas para ser ejecutadas por los kernel. En general, el uso de esta API permite una mejor optimización del código y ofrece mejores tiempos de respuesta. La librería permite, al igual que su versión para host, generar secuencias de números pseudoaleatorios y quasialeatorios.

Las funciones principales son las siguientes: curand_init(), para iniciar la secuencia, la función curand(), que permite obtener un número pseudoaleatorio entero sin signo. curand_uniform(), obtiene un número en coma flotante de una distribución uniforme en el intervalo (0.0, 1.0]. curand_normal(), obtiene un número en coma flotante de una distribución normal en el intervalo (0.0, 1.0] utilizando un valor de 0 para la media y de 1 para la desviación estándar.

También se dispone de versiones de estas dos últimas funciones para obtener secuencias con números de doble precisión: curand_uniform_double() y curand_normal_double(). Adicionalmente contamos con dos funciones que permiten obtener pares de elementos utilizando la estructura float2 y double2, que componen de un campo x y otro y, para obtener secuencias de distribución normal: curand_normal2() y curand_normal2_double() respectiva-mente.

El valor que se pasa a estas funciones es el estado del generador, que puede ser de dos tipos, dependiendo del generador que deseemos utilizar: curandState si deseamos utilizar un RNG pseudoaleatorio, o curandStateSobol32 si deseamos obtener valores de la secuencia quasialeatorios. Para el primer caso, en el que se trabaja con números pseudoaleatorios hay que pasar una semilla a la función, mientras que en el caso de quasialeatorios hay que pasar como argumento los vectores de dirección y el desplazamiento. Finalmente cabe destacar una función para pasar por alto elementos de la secuencia, que es lo mismo que llamar varias veces a la función curand() pero mucho más rápido: skipahead(), en esta función indicamos el número de elementos a descartar y el estado.

Para la inicialización del generador de números aleatorios en CUDA, utilizamos la función curand_init(). Como particularidad hay que utilizar una estructura de datos que mantenga el estado del generador y lo actualice durante el proceso de generación de números. Necesitamos mantener y pasar esa estructura en cada kernel que quiera generar números aleatorios:

curandState* devStates;

cudaMalloc ( &devStates, Elementos*sizeof( curandState ) );

/* inicializacion de numeros aletorios */

randomize<<< ( Elementos + ( hilos – 1 )) / hilos, hilos >>> ( devStates, time(NULL) );

En el código anterior, en primer lugar se declara la estructura de estados de tipo curandState ya que se van a utilizar números pseudoaleatorios.

En la siguiente sentencia se reserva la memoria que necesita el dispositivo para actualizar los estados a cada hilo del kernel. Finalmente la última instrucción consiste en una llamada al kernel reandomize (la hemos llamado igual que la función existente en muchos lenguajes de alto nivel para la inicialización de la generación de números aleatorios). Como semilla se emplea el valor del tiempo actual del sistema, preparando así los estados para la obtención de números aleatorios.

__global__ void randomize( curandState *state, unsigned long semilla )

{

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

curand_init ( semilla, id, 0, &state[id] );

}

El kernel randomize se encarga de iniciar el estado para cada hilo al llamar a curand_init. Separamos esta tarea en un kernel independiente con el objetivo de aumentar el rendimiento global del algoritmo. Por un lado para reducir el uso de la pila de memoria por hilo, ya que la generación de los estados de números aleatorios, en el peor de los casos puede llegar a consumir hasta 16KB. Por otro lado al separar esta tarea de la obtención de números aleatorios curand() permitimos al kernel, que está manejando la obtención, disponer de más velocidad y aumentar sus posibilidades de concurrencia al requerir menos recursos del procesador. Esto permite manejar más bloques por procesador.

Para obtener números aleatorios dentro del kernel utilizaremos la función curand_uniform( &localState ), que como hemos indicando anteriormente obtienen un elemento de la secuencia (en nuestro caso pseudoaleatoria) de elementos en coma flotante dentro del intervalo (0.0, 1.0].

[1] Las siglas en inglés de: Random Number Generator (Generador de Números Aleatorios)

CUDA 5.0 beta

2012-04-29T17:06:00.000-07:00

Buenas Noticias, parece que el lanzamiento de las primeras GPU con arquitectura Kepler, estan acelerando la salida del nuevo SDK CUDA y sus drivers.

Ya está disponible la versión 5.0 Preview. Al parecer, vienen con nuevas herramientas y mejoras que harán más facil la implementación de aceleración de programas y procesos mediante GPU, y como no podía ser de otra manera, se incluyen los aceleradores basados en la arquitectura de cálculo NVIDIA Kepler™, la arquitectura de procesamiento más rápida, más eficiente y de mayor rendimiento que se haya inventado jamás.

Entre las principales características de CUDA 5 se incluyen (citando fuentes de nVidia):

Paralelismo dinámico: lleva la aceleración en la GPU a nuevos algoritmos
Los procesos de la GPU pueden generar dinámicamente nuevos procesos, lo que permite a la GPU adaptarse a los datos. Al reducir las interacciones con la CPU, el paralelismo dinámico simplifica enormemente la programación paralela. Además, posibilita la aceleración en la GPU de un número más amplio de algoritmos conocidos, tales como los utilizados en muchas aplicaciones de dinámica de fluidos computacional y métodos de refinamiento adaptable de mallas.

Llamadas a otras librerías desde la GPU: el origen de un ecosistema de desarrolladores
Una nueva librería BLAS de CUDA permite a los desarrolladores utilizar el paralelismo dinámico para hacer llamadas a sus propias librerías desde la GPU. El desarrollador puede diseñar API que permitan a otros desarrolladores ampliar la funcionalidad de sus kernels e implementar retrollamadas (callbacks) en la GPU para adaptar la funcionalidad de otras librerías externas que admitan llamadas desde la GPU. La posibilidad de “vincular objetos” proporciona una forma eficiente y conocida de desarrollar aplicaciones de GPU grandes ya que otorga a los desarrolladores la capacidad de compilar múltiples archivos CUDA de origen en diferentes archivos de objetos y vincularlos a librerías y aplicaciones de mayor tamaño.

Soporte de GPUDirect para RDMA: minimiza los cuellos de botella en la memoria del sistema
La tecnología GPUDirect proporciona comunicación directa entre las GPU y otros dispositivos PCI-E, y acceso directo a la memoria entre las tarjetas de red y la GPU. También reduce considerablemente la latencia de MPISendRecv entre los nodos de GPU de un cluster y mejora el rendimiento general de las aplicaciones.

NVIDIA Nsight Eclipse Edition: una forma rápida y fácil de generar código CUDA
Con NVIDIA Nsight Eclipse Edition, los programadores pueden programar, depurar y analizar el rendimiento de las aplicaciones de GPU dentro del conocido entorno de desarrollo integrado basado en Eclipse para las plataformas Linux y Mac OS X.Un editor de CUDA integrado y los ejemplos de código CUDA ayudan a agilizar la programación, y la refactorización automática del código facilita el traslado de bucles basados en la CPU a kernels CUDA.Un sistema de análisis integrado proporciona análisis automatizado del rendimiento y una guía de instrucciones detalladas para eliminar los cuellos de botella del código. Por último, la diferenciación visual de la sintaxis ayuda a distinguir mejor entre el código de la GPU y el de la CPU.

Estoy deseando probarlo!! Un saludo a todos.

Paralelizar en CPU vs GPU

2012-04-21T13:24:00.000-07:00

Hoy vamos a tratar dos cuestiones clave a la hora de abordar una implementación en GPU:

los que estén acostumbrados a implementar algoritmos paralelos en CPU tendrán un acceso algo mas llano a la programación GPGPU, pero hay ciertas diferencias a tener en cuenta: la paralelización en CPU difiere considerablemente de la paralelización sobre GPU. La diferencia fundamental se encuentra en el número de procesadores con el que contamos en una y otra arquitectura, si bien en CPU contamos con pocos procesadores pero muy rápidos, en GPU disponemos de muchos procesadores algo más lentos.

La consecuencia directa de esta diferencia implica adoptar un enfoque diferente, orientado a la unidad mínima de procesamiento independiente, evitando en la medida de lo posible, la dependencia entre estas unidades, esto lo tenemos que reflejar tanto en el diseño de las estructuras de datos, como en las tareas que realicemos sobre estas.

Además en GPU podemos adoptar 2 enfoques diferentes:

A. Enfoque orientado al bloque

Cuando implementamos un kernel en CUDA, éste puede diseñarse de tal manera que deba existir cierta comunicación entre los hilos que ejecutan las tareas. De alguna manera estos hilos colaboran entre sí y pueden pasarse información entre ellos, estableciendo mecanismos de sincronización para acceder a la memoria local donde residen las variables que son comunes al proceso. En este enfoque existe una limitación, la comunicación o colaboración sólo es posible entre hilos que pertenecen al mismo bloque de ejecución. La ventaja principal es que es posible realizar una comunicación muy rápida mediante la memoria local, ya que su tiempo de latencia es mínimo (en comparación con el acceso a memoria global).

La implementación de este enfoque requiere entonces de dos fases, una primera en la que todos los hilos de cada bloque procesan en colaboración obteniendo un resultado por bloque, y una segunda fase en la que son procesados los resultados de cada bloque. Esta tarea se puede ejecutar en otro kernel de GPU o bien directamente en CPU.

B. Enfoque orientado al hilo

Este enfoque implica una estructura más plana, en la que para una determinada ejecución existe independencia total entre los elementos de cómputo y la tarea que se aplica. Por ejemplo, una tarea que asigne un número aleatorio a cada elemento de un vector puede ser un caso típico.

Aunque inicialmente puede parecer que sí tienen dependencia dada la naturaleza pseudo-aleatoria de estos números, es posible hacer una descomposición de tareas de forma que se pueda dar la independencia de procesos.

Nuestro algoritmo de reducción de vectores propuesto en post anteriores, está implementado utilizando este enfoque. Así, para desligar la dependencia se han dividido las ejecuciones del kernel en niveles, que son ejecutados secuencialmente pero para cada ejecución todas las operaciones son independientes.

El enfoque orientado al hilo permite de una manera más sencilla alinear las lecturas de memoria con cada uno de los hilos que hacen uso de esa memoria, reduciendo el “lag” que produce la latencia de la memoria sobre todo cuando esa relación es de uno a uno.

En resumen: dependiendo del tipo de problema una estrategia será más beneficiosa que la otra, aunque en muchos casos lo idóneo será una mezcla de todas.

Un ejemplo básico, parte III: Análisis del rendimiento

2012-04-06T13:16:00.000-07:00

Para realizar las pruebas se ha utilizado la plataforma "Davinci" que es una computadora gestionada por el Grupo de Investigación TIC-146 Supercomputación: Algoritmos de la Universidad de Almería. Este servidor multicore dedicado a computación GPGPU cuenta con:

 2 Intel Xeon Quad Core E5640 2.67 GHz 64 bits

 GeForce GTX 295 (2x240 cores)

 2 Tesla C2050 (Fermi) (448 cores)

 GeForce GTX 480 (480 cores)

Todas las GPU son Compatibles con CUDA, si bien sólo las Tesla C2050 y la GeForce GTX 480 ofrecen las prestaciones de cómputo 2.0 necesarias para aprovechar al máximo las prestaciones de CUDA 4.

El sistema operativo de este servidor es un Linux Version de kernel 2.6.26-2-amd64, con una distribución Debian GNU/Linux 5.0.8 (lenny). La versión de Cuda instalada es la 4.0.

Los programas se han desarrollado en lenguaje C y CUDA, utilizando los compiladores gcc (The GNU Compiler Collection) y nvgc (el compilador CUDA de nVidia).

Las pruebas se han realizado sobre diferentes tamaños de vectores partiendo de 10.000 elementos, y con incrementos de ese mismo número hasta llegar a los 1.000.000 de elementos del vector.

A. Reducción de matrices

En las pruebas no se ha tenido en cuenta el tiempo de introducir valores aleatorios a los elementos del vector, ya que es común para los tres escenarios. Los tres escenarios comparados son los siguientes: por un lado la versión secuencial de la suma ejecutada en la CPU, nuestro algoritmo de reducción de vectores propuesto, y la reducción de vectores que la nueva librería Thrust® ofrece y que simplifica la programación paralela para CUDA.

Cada uno de los test se ha realizado 10 veces consecutivas y se ha calculado la media aritmética de los resultados obtenidos.

Para medir los tiempos de respuesta de cada uno de los algoritmos hemos utilizado la gestión de eventos de CUDA. Dentro del tiempo de cómputo para las versiones paralelas, se ha incluido el tiempo de copia de memoria principal a memoria de dispositivo y viceversa. La siguiente tabla muestra los resultados obtenidos en el test:

N	CPU	GPU1	Thrust
10000	0.0000382	0.00007563	0.00031225
20000	0.0000758	0.00010429	0.00031417
30000	0.0001138	0.00013749	0.00033243
40000	0.0001512	0.00016749	0.00035026
50000	0.0001884	0.00019928	0.00036693
60000	0.0002266	0.00022946	0.0003838
70000	0.0002640	0.00025606	0.00040174
80000	0.0003020	0.00028915	0.00041772
90000	0.0003392	0.00031836	0.00043793
100000	0.0003770	0.00035073	0.00046362
110000	0.0004146	0.00038296	0.00049638
120000	0.0004536	0.00041684	0.00051473
130000	0.0004896	0.00044498	0.00056173
140000	0.0005274	0.00046391	0.00058575
150000	0.0005678	0.00048342	0.0005885
160000	0.0006090	0.00050197	0.00060389
170000	0.0006442	0.00051199	0.00061768
180000	0.0006782	0.00052612	0.00063185
190000	0.0007172	0.00054233	0.00064476
200000	0.0007548	0.00055549	0.00066341

Tabla 1. Tiempo (en segundos) de respuesta para los algoritmos

La tabla anterior muestra en la columna N el número de elementos del vector tratado, la columna CPU muestra los tiempos de respuesta obtenidos para el algoritmo en serie sobre la CPU, la columna GPU1 se refiere al algoritmo de reducción propuesto y finalmente Thrust muestra los tiempos de respuesta obtenidos utilizando las librerías Thrust®.

Como se puede observar el rendimiento obtenido con nuestra propuesta de algoritmo es bastante superior al rendimiento obtenido por el algoritmo de reducción de la librería Thrust de nVidia.

Figura 1. Comparativa rendimientos obtenidos

Se realizó una segunda comparativa entre CPU y GPU1 con un crecimiento exponencial de elementos de la forma 1024ⁱ para i = 1 hasta 12 con los siguientes resultados:

N	CPU	GPU1	SpeedUp
1024	0.0000042	0.00004873	0.08618921
2048	0.0000082	0.00004486	0.18279091
4096	0.000016	0.00005397	0.296461
8192	0.0000318	0.00006936	0.45847751
16384	0.0000624	0.00009197	0.67848211
32768	0.000124	0.00014606	0.84896618
65536	0.0002474	0.00024655	1.00344758
131072	0.000494	0.00045282	1.09094121
262144	0.000992	0.00070173	1.41364912
524288	0.0019882	0.00116573	1.70554073
1048576	0.0039892	0.00213745	1.86633605
2097152	0.0109062	0.00441923	2.46789599

Tabla 2. SpeedUp respecto a GPU1

En la siguiente figura se puede observar cómo la pendiente de la tendencia creciente exponencial de la GPU1 es minorada por el factor de SpeedUp obtenido en cada tramo:

Figura 2. Evolución del SpeedUp al aumentar N

B. Determinación del Punto Muerto Computacional[1]

Trataremos de determinar el punto donde es indiferente realizar los cálculos en CPU respecto de la GPU ó

GPU respecto de la CPU porque el rendimiento sea el mismo. En las pruebas anteriores, hemos podido observar cómo en los tramos donde el número de elementos que componen el vector a reducir es bajo, se daba la situación que el algoritmo en serie ejecutado en la CPU obtenía mejores resultados que sus homólogos paralelizados sobre GPU. Evidentemente debemos de tener esto en cuenta a la hora de decidir en donde debemos ejecutar cada trabajo.

Dicho punto coincidirá con el punto de corte de las funciones que representan el tiempo de ejecución que emplea cada algoritmo en procesar exactamente el mismo trabajo (incluidos en este caso para los algoritmos de GPU, los tiempos de proceso consumidos en la transferencia de memoria principal a memoria del dispositivo y una vez realizado el cálculo, vuelta a recibir los datos a memoria principal) y obtener exactamente el mismo resultado.

Si realizamos un zoom sobre los datos tratados apreciamos

Figura 3. Punto de corte función CPU y GPU1

un punto de corte entre las funciones de tiempo de respuesta sobre el tamaño N del vector. Así, en la Figura 18 podemos observar cómo ese punto de corte está cerca de los 70000 elementos. Podemos decir que en este punto es indiferente la utilización de un algoritmo u otro en términos de tiempos de ejecución. Por debajo de este valor sería mejor computar sobre la CPU, mientras que por encima de ese valor sería mejor hacerlo en GPU.

Evidentemente, esta visión es algo ingenua, y se deberían tener otros factores en cuenta a la hora de tomar esta decisión.

- Coste de utilización: es el coste que implica la utilización de cada escenario. En muchos casos (como en grandes centros de computación) se calcula como la suma del coste de adquisición más el coste en potencia (watios) por euro de cada recurso.

- Coste de oportunidad[2]: es un concepto más sutil pero muy interesante a la hora de tomar decisiones económicas. Aplicado a nuestro caso, se puede definir como el coste que supone la decisión de computar un algoritmo en un recurso (que es escaso) y no computar o descartar la computación de otro algoritmo en ese mismo recurso.

Teniendo en cuenta los costes anteriores se podrán tomar decisiones por ejemplo: de utilizar una GPU para realizar un determinado proceso aunque el tiempo de respuesta sea similar o peor al obtenido en CPU, debido a que el coste por watio de la GPU es muy inferior, y a que el coste del trabajo que puede realizarse en la CPU sea mucho mayor.

C. Determinación del Factor de Desenrollado

Ya se habló en post anteriores de la importancia de este factor, y de cómo ha sido un caballo de batalla en los compiladores desde los inicios de la informática.

El objetivo fundamental que perseguimos con esta técnica es compensar los efectos negativos en el rendimiento que producen los tiempos de latencia que se producen en los accesos a memoria.

De forma intuitiva, podemos expresar este efecto como el que se produce en las líneas de autobuses urbanos: un autobús puede transportar 50 pasajeros a la vez, tarda prácticamente lo mismo si trasporta 1 o si trasporta a 50. El problema está en que la compañía de autobuses espera optimizar su beneficio, y establece unos horarios de salida y llegada, por eso, aunque un pasajero suba al autobús, deberá esperar cierto tiempo hasta ser trasportado a su destino. En nuestro caso los pasajeros son los datos que necesitamos tratar, la parada de origen es la memoria global del dispositivo, y la parada final, los registros locales del hilo que ejecuta un kernel. La arquitectura CUDA, está diseñada para esperar en la parada de origen entre 400 y 600 ciclos de reloj, además hay que tener en cuenta el número de pasajeros que se permite trasportar en cada viaje.

Como podemos deducir, para optimizar el proceso, cada vez que realicemos un acceso a memoria global del dispositivo, lo haremos de tal manera que llenemos el autobús, y todos esos datos trasportados, los procesaremos en el mismo hilo.

Para nuestro algoritmo, hemos realizado una serie de pruebas partiendo de un factor de desenrollado de 2 hasta 10, 25, 50 y 100, la siguiente tabla muestra los tiempos de respuesta obtenidos en cada caso.

UF	GPU1
2	0.00010693
3	0.00008918
4	0.00008308
5	0.00007948
6	0.00008024
7	0.00007638
8	0.00007826
9	0.00008013
10	0.00007762
25	0.0000807
50	0.00008548
100	0.00011343

Tabla 3. Tiempos de respuesta para distintos valores de UF

Hemos determinado de forma empírica que el mejor factor es el de 7. Para otros algoritmos puede variar ya que todo esto depende del tipo de dato que se esté transportando, del tamaño de los registros necesarios en el kernel, etc.

Figura 4. Selección del mejor factor de Unroll

[1] Este término ha sido acuñado para el presente trabajo como analogía al término económico “Punto Muerto” o “Umbral de rentabilidad” donde los ingresos son iguales a los gastos

[2] Es el término económico utilizado para referirnos al coste de la mejor opción no realizada. Este término fue inventado por Friedrich von Wieser en su Theorie der gesellschaftlichen Wirtschaft (Teoría de la economía social). El coste de la oportunidad es aquello a lo que renunciamos cuando tomamos una decisión económica.

Un ejemplo básico, parte II: Implementación del Kernel

2012-03-23T12:59:00.000-07:00

Como os prometí en el anterior post, hoy veremos la implementación del kernel para la reducción de vectores:

Nuestro kernel básicamente consiste en realizar la suma de una serie de elementos y guardarla en el primero de ellos. Cada ejecución se encargará de una posición de memoria diferente, con lo que no es necesario preocuparnos de implementar ningún mecanismo de sincronización. El número de elementos a sumar viene dado por el factor de desenrollado (UF).

Los parámetros de entrada del Algoritmo son: Datos: puntero a la estructura que contiene los datos a sumar, Elementos: número de elementos que se suman, Impar: indica si el número de elementos es impar (para evitar la realización del cálculo del módulo, dentro del kernel, ya que es computacionalmente costoso), y finalmente UF, que se utilizará para indicar el número de elementos a sumar dentro de cada hilo.

*Algoritmo 3:* kernel de reducción de vectores
1	tid ← threadIdx.x + blockIdx.x blockDim.x*
2	t ← (ene/UF)
3	si tid < t hacer
4	para i←1 hasta UF hacer
5	a[tid] ← a[tid] + a[(tid+(ti))]*
6	Si tid = 0 y impar > 0 hacer
7	para i←1 hasta impar+1 hacer
8		a[tid] ← a[tid] + a[(ene-i)]

Algoritmo 1 Kernel de reducción de vectores

1.-La primera línea del algoritmo, “tid ← threadIdx.x + blockIdx.x * blockDim.x”, es la típica asignación para el cálculo de la posición global de ejecución del kernel. Así se identifica exactamente al hilo que se está referenciando en cada uno de los kernels y es la forma que tenemos para organizar el proceso en paralelo. Como en nuestro caso estamos trabajando únicamente con una dimensión, sólo se referencia al valor del campo “x” de cada variable (threadIdx, blockIdx y blockDim), pero podrían ser hasta tridimensionales, para ello se utilizan los campos “y” y “z” para referirnos a la segunda y tercera dimensión respectivamente.

2.- Dado que el número real de cálculos que tenemos que ejecutar no siempre coincide con el total de hilos ejecutados en los Kernels, debemos controlar esta situación para que el dispositivo no haga trabajo inútil (o redundante). Por este motivo la segunda línea del algoritmo calcula el número total de threads a emplear y en la tercera sólo se permite la ejecución del proceso de reducción a los hilos que están dentro del rango necesario.

4.-El siguiente paso define al bucle que empleamos para la suma de elementos del vector, que viene definido por el factor de desenrollado que empleamos en la llamada.

5.-La siguiente línea realiza la suma de elementos y la acumula en el primer elemento.

6.-Finalmente, tenemos una condición especial: sólo si se trata del primer hilo y el número de elementos a procesar es impar, adicionamos esos elementos impares también al primero.

Las operaciones realizadas en el siguiente ejemplo son: (3+6=9), (8+5=13), (4+2=6), y en la siguiente iteración (9+13+6=28).

En la siguiente Figura se puede observar gráficamente el proceso, para un factor UF de 2:

Figura 1. Reducción con UF = 2

Las propiedades fundamentales de este método de reducción que acabamos de proponer son:

· No hay relación entre bloques ni hilos del mismo stream, con lo que no es necesario establecer ningún tipo de mecanismo para sincronizar el algoritmo en ningún punto. Es un algoritmo orientado al hilo y no al bloque.

· Adecúa la petición de recursos del sistema en función al estado de su ejecución.

· Sólo se realiza un desplazamiento de datos de memoria central a memoria de dispositivo por el total de elementos a sumar, y además solo se necesita desplazar el primer elemento a memoria principal para obtener el resultado de la reducción.

· Con el control de UF se reduce el efecto de la latencia de memoria global del dispositivo, no siendo por tanto necesario realizar las lecturas alineadas con el orden de los hilos.

Ya hemos presentado un algoritmo para la reducción de vectores, que como sabréis, es una de las operaciones más típicas en multitud de escenarios, tanto en análisis de datos, como en optimización. El próximo post presentará los datos obtenidos en las pruebas.

Un ejemplo básico, parte I: la llamada al kernel

2012-03-09T12:47:00.000-08:00

Vamos a diseñar un algoritmo que permita realizar llamadas a nuestro kernel de tal forma que conociendo N número de elementos a procesar o tamaño del bucle en versión serie, utilice los mínimos recursos (de procesador gráfico) minimizando también las transferencias de memoria. Nos centraremos en un típico caso que en el que se emplea la programación paralela: Reducción de vectores. En la literatura hay multitud de ejemplos y diversos enfoques, nVidia también realiza varias implementaciones hasta 7 diferentes, van ilustrando los diferentes factores para mejorar el rendimiento. Yo implementé mi propio método, que ha resultado ser más rápido que las implementaciones de nVidia y de la librería Trust en varios escenarios. En el presente post y el siguiente lo presento:

Algoritmo para la llamada al kernel

Primero necesitamos obtener las características del dispositivo donde se van a ejecutar los programas, eso lo podemos obtener en tiempo de ejecución interrogando las capacidades de cada dispositivo. Las tres variables cuyo valor necesitamos saber son: N número de elementos a procesar, UF factor de desenrollado (Unroll factor), que determinará la velocidad de segmentación de los datos y el número de instrucciones serie que ejecuta el kernel en cada hilo, y datos que es el puntero al bloque de datos a procesar.

Una vez determinados los valores de estas variables podremos calcular el factor de ocupación. El consumo máximo de tiempo se emplea en la migración de los datos desde memoria principal a la memoria del dispositivo, por lo que realizaremos sólo un proceso de carga y otro de descarga. El tamaño de bloques necesarios para la ejecución se irá reduciendo de forma logarítmica en cada iteración, ya que este algoritmo está orientado a las operaciones de reducción de vectores.

Para acumular los valores obtenidos se sobrescribe el vector de entrada con los resultados que se van obteniendo en cada cálculo, evitando así crear dos estructuras de datos, una de entrada y otra de salida. Además para cada llamada se aprovecha la información residente en la memoria del dispositivo, reduciendo la cantidad de datos transferida entre la memoria central y la del dispositivo.

*Algoritmo 2:* Algoritmo de llamada a un kernel
1	Obtener_caracteristicasGPU (max_bloques, max_hilos)
2	Determinar_factor_ocupacion(N,max_bloques, max_hilos, registros_kernel, nhilos)
3	UF ← 7
4	auxN ← N
5	huerfanos ← (auxN modulo UF)
6	Datos_a_dispositivo(datos_dispositivo, datos)
7	Para i ← 1 hasta (log(N) / Log(UF))
8		Ejecutar en paralelo Reduce<<<(auxN+nhilos-1)/nhilos,nhilos> >>(datos_dispositivo, auxN, huerfanos, UF)
9		auxN ← auxN / UF
10		huerfanos ← (auxN modulo UF)
11	Datos_a_memoria_principal(datos, datos_dispositivo)

Algoritmo 1.Llamada al Kernel

Pasaremos a comentar cada una de las líneas de l algoritmo:

1.- Obtener_caracteristicas GPU (max_bloques, max_hilos): Indica los límites de la arquitectura seleccionada, que se utilizarán para aplicar las reglas que obtendrán el mejor factor de ocupación. Adicionalmente se obtienen otras cualidades de los dispositivos como la capacidad computacional (Compute capability) que indicará cuales son los dispositivos que permiten aprovechar las nuevas capacidades de CUDA 4, caracterizados por tener una capacidad igual o superior a 2.0.

2.- Si se puede estimar el tamaño N de nuestro proceso, y esto se puede traducir en el número de hilos necesarios para invocar al kernel, se podrá determinar el mejor factor de ocupación. Este cálculo óptimo de la relación hilos/bloques se puede realizar fácilmente porque al compilar el kernel se puede conocer el tamaño de de los registros necesarios para cada hilo, y consecuentemente, se puede determinar el número máximo de hilos por bloque y de bloques por procesador.

3.- En la línea 3 hay que proporcionar un valor de Unroll Factor (UF). Para ello se deberán de realizar varias pruebas para obtener el mejor valor para este parámetro, ya que dependerá mucho de la estructura de datos en memoria que estemos utilizando, el tipo de memoria que utilicemos para su gestión y la coalescencia de los accesos que realicemos. En nuestras pruebas hemos utilizado la memoria global del dispositivo, y los accesos no han sido alineados con la ejecución de los kernels, Como veremos más adelante para nuestro algoritmo el mejor valor ha sido 7.

4.- En la línea 4 utilizamos una variable para guardar el tamaño de procesamiento remanente. El valor de esta variable se va reduciendo logarítmicamente conforme avanza el proceso.

5.- La variable "huérfanos" se utiliza para determinar cuándo se produce un elemento huérfano, el deberá de tratarse de forma especial.

6.- La línea 6 copia los datos desde memoria principal a la memoria global del dispositivo para ser procesados por el kernel.

7.- El bucle realiza una serie de llamadas al kernel en las que se ajustan los bloques que demandamos al dispositivo en función del remanente de ejecuciones. El número de iteraciones se corresponde con el cociente entre el logaritmo de N y el logaritmo de UF-.

11.- Una vez concluidas las ejecuciones se descargará el dato deseado a memoria principal para obtener el resultado final. La utilización del las llamadas descendentes, la vamos a entender a continuación en cuanto se exponga el Kernel.

Os dejo analizando el algoritmo, en el próximo post, veremos la implementación del Kernel. Hasta pronto

En busca del Kernel perdido

2012-03-02T10:48:00.000-08:00

Con todo lo expuesto en los post anteriores debemos confeccionar nuestra lista de prioridades a la hora de diseñar los kernel que ejecuten las GPU. Resumiendo, debemos tener en cuenta los siguientes factores:

· Factor de “unroll”

· ILP: Nivel de Paralelización por instrucción

· Factor de ocupación de la GPU, que se compone de:

o Nº de bloques que se necesitan

o Nº de hilos por cada bloque

o Registros utilizados en el kernel

· Transferencias de memoria CPU <-> GPU

· Tipos de memoria a utilizar, Global, compartida,…

Teniendo todo esto en cuenta, nuestro kernel debe tener un número de instrucciones mínimo que ocupe la latencia de la arquitectura, y se intentará invocar al kernel con el tamaño justo de bloques e hilos por bloques que maximicen la ocupación del tamaño requerido. Se minimizará el número de transferencias de memoria entre CPU y GPU. Se intentará evitar las latencias de acceso a memoria global mediante accesos ordenados en memoria compartida.

En el próximo post, ilustraré con un ejemplo todo esto que hemos comentado

Los diversos usos de la memoria

2012-02-24T10:40:00.000-08:00

Es muy importante saber que ventajas e inconvenientes tiene el uso de cada tipo de memoria de dispositivo, antes de realizar cualquier diseño, deberías tener claro cómo CUDA gestiona cada tipo de memoria, y cual es el coste o penalización que supone su utilización.

La siguiente Figura muestra los tipos de memoria en las GPUs y el ámbito de acceso de cada elemento

Figura 1. Los diferentes tipos de memoria en la GPU

Cada una de ellas tiene unos propósitos concretos que pasamos a detallar: Los registros y la memoria local son espacios de memoria accesibles únicamente por el hilo al que pertenecen. Por otro lado, a la memoria compartida pueden acceder todos los hilos del mismo bloque, utilizando para ello mecanismos de sincronización que mantienen el orden en las lecturas y escrituras. Al resto de memorias se puede acceder desde toda la malla que ejecuta un kernel. Se distinguen en que la memoria global no dispone de memoria caché, la memoria de constantes es de sólo lectura, y la memoria de texturas dispone de una caché espacial.

El cuello de botella principal lo encontramos en la transferencia de memoria entre la CPU y la GPU, de modo que la penalización nos la encontramos en la latencia del acceso a la memoria global (entre 400 y 600 ciclos). Como solución posible se plantea el alineamiento de la memoria global para poder hacer lecturas ordenadas según la posición de cada hilo dentro de la tanda, se puede utilizar la memoria compartida para hacer la ordenación dentro de cada bloque. No siempre es posible hacer esta ordenación y el problema se complica cuando aumentan las dimensiones (matrices por ejemplo).

Entorno respetuoso de ejecución

2012-02-10T10:36:00.000-08:00

A simple vista podríamos pensar que si un Grid solo puede ejecutar el mismo kernel y este es controlado solo por un dispositivo el cual solo es accesible desde el mismo thread del host, no es necesario tener en consideración otras interacciones entre distintos procesos. Hasta la versión 3.2 de CUDA esto era así,

Figura 1. CUDA 4.0 revoluciona la computación paralela en GPGPU

pero con la aparición de CUDA 4.0, se han eliminado ésta y muchas otras barreras: en cuanto a la plataforma GPUDirect™ se pueden realizar transferencias de memoria directamente de una GPU a otra sin pasar por memoria host; arquitectura 64-bit nativa; gestión de clúster; soporte multi-GPU en la que se pueden compartir múltiples GPUs a través de múltiples hilos (de host), un solo hilo (de host) puede acceder a todas las GPUs. En cuanto al modelo de programación: direccionamiento de memoria unificado, ahora se podrá contemplar la memoria de sistema y la memoria de cada una de las GPUs como una sola unidad de memoria simplificando así la migración de programas existentes.

Además el nuevo toolkit viene acompañado de librerías OpenSource y plantillas de primitivas que ayudan y facilitan el manejo de vectores, matrices, etc. Un nuevo sistema de análisis de rendimiento permitirá obtener todas las métricas necesarias para afinar nuestros desarrollos.

Por eso es importante diseñar no sólo kernels, sino también las rutinas de host que los llaman, para que hagan un uso razonable de los recursos, por ejemplo modificando el tamaño de los bloques e hilos en las llamadas dinámicamente, en lugar de dejar un valor fijo para todas las llamadas.

No es “Speed Up” todo lo que reluce

2012-01-28T10:18:00.000-08:00

Cuando realizamos optimizaciones en nuestro código, debemos ser conscientes de las limitaciones que impone la arquitectura, en muchos casos lo que puede parecer una mejora acaba siendo un cuello de botella.

Aunque la técnica de “unroll” o desenrollado de bucles descrita en el post anterior parezca la “Panacea[1]”, hay que ser cautelosos, porque no siempre será así. Es necesario hacer un estudio detallado de cada caso para aplicar factores de desenrollado apropiados a cada uno de los algoritmos que se pretenden computar, ya que un desenrollado de bucles agresivo, puede aumentar demasiado las tareas de planificación (scheduling) y aumentar demasiado el uso de registro en el cuerpo del bucle, que en las CPU puede implicar una migración de los valores de los registros a memoria y esto puede desembocar en una degradación del rendimiento del programa. Por otro lado, si se utilizan factores de desenrollado muy altos, el tamaño del cuerpo del bucle resultante podría desbordar la caché de instrucciones, seguido de pérdidas de caché, y en consecuencia reduciéndose el rendimiento.

En programación GPGPU nos encontramos además con una serie de restricciones que limitan el rendimiento de los programas. Antes de proceder a la implementación deberemos estudiar el impacto que supone el desenrollado de bucles para los programas GPGPU en el contexto de las restricciones de esos recursos y hallar el factor de “unroll” óptimo para ese programa y dispositivo.

ILP y Factor de Ocupación (OF): como ya se ha visto, el nivel de paralelismo por instrucción juega un papel importante a la hora de aumentar el rendimiento, ya que el tamaño de la caché de instrucciones “I-cache” es limitado. Es importante no sobrecargar el cuerpo de la función con demasiadas instrucciones que provoquen un fallo de caché y la consiguiente penalización en el rendimiento. Por otro lado en las GPU el factor de ocupación puede ser igual o más importante que el anterior a la hora de ocultar las latencias de acceso a memoria y de flujo de instrucciones. Debido a la interrelación que hay entre ILP y OF, se hace muy complicado estimar los valores óptimos. A continuación se muestra un caso concreto que sirve de ejemplo.

En el apartado B de esta sección se expuso brevemente el modelo de programación de CUDA compuesto por Mallas (Grids), Bloques (Blocks) e Hilos (Threads), ahora vamos a compararlo con la arquitectura de las GPUs para comprender mejor el problema del factor de ocupación. Una GPU está compuesta por x SM “stream multiprocessor”, cada uno dispone de c SP “stream processor”, en cada uno hay w Warps[2], que a su vez pueden ejecutar t Threads. La clave para comprender el problema del factor de ocupación reside en la diferencia real que existe entre la capacidad computacional de la arquitectura que utilicemos y la capacidad que solicitamos en la llamada al kernel. La Figura 1 muestra las capacidades para distintos dispositivos con arquitectura CUDA.

Device Name	GeForce GTX 480	Tesla C2050	GeForce GTX 295
Compute capability	2.0	2.0	1.3
Memory Informat.
Total global mem	1.536 MB	2.687 MB	896 MB
Total constant Mem	64 Kb	64 Kb	64 Kb
Max mem pitch	2048 MB	2048 MB	2048 MB
MP Information
Multiprocessor count	15	14	30
Shared mem per mp	48	48	16
Registers per mp	32768	32768	16384
Threads in warp	32	32	32
Warps/MP	48	48	32
Blocks/MP	8	8	8
Max threads/block	1024	1024	512
Max thread dim	(1024,1024, 64)	(1024,1024, 64)	(512, 512, 64)
Max grid dim	(65535, 65535, 65535)	(65535, 65535, 65535)	(65535, 65535, 1)
Warps/block in SP	6	6	4
threads in block	192	192	128
total pararell threads	2880	2688	3840

Figura 1. Capacidades de algunos dispositivos

De todos los valores de la tabla ahora mismo nos interesan: “Multiprocessor count”, que es el número de multiprocesadores que tiene el dispositivo, “Blocks/MP” bloques por multiprocesador, “”Warps/MP” que son agrupaciones de hilos dentro de cada bloque y “Threads in warp” que son el número de hilos que ejecutan en cada warp. Sin tener nada más en cuenta, podemos calcular el número total de hilos que cada dispositivo es capaz de ejecutar de forma concurrente: por ejemplo para el dispositivo GeForce GTX 480 tenemos 15 procesadores x 8 bloques x 6 warps x 48 hilos, en total 34.560 hilos concurrentes ejecutando el mismo kernel. Evidentemente la arquitectura nos está imponiendo además ciertas restricciones: “Max threads per block” está indicando un número (1024) mayor al teórico. Si en un bloque hay 6 warps y 32 hilos (192 hilos en total) ¿cómo es posible alcanzar esa cifra? CUDA puede utilizar hasta 8 bloques por procesador, es decir, que si utilizamos 192 hilos por bloque, dispondremos de 8 bloques por procesador, pero si necesitamos 256 hilos, dispondremos de 6 bloques, para 512 hilos solo 3 bloques, para 1024 solo un bloque por procesador. Cuando se solicitan los recursos en las llamadas a los kernels, CUDA realiza una segmentación y planificación de ejecución de todos los bloques solicitados, ahí es donde se pueden producir holguras si no hemos diseñado bien nuestra llamada, reduciéndose el factor de ocupación.

Otra cuestión a considerar, también muy importante es el uso de los registros que necesita nuestro Kernel, a mayor uso de registros, menor disponibilidad de bloques por procesador. Estas dos variables son las que van a modificar el factor de ocupación que, como hemos podido inferir, es el aprovechamiento físico de la arquitectura para cada tanda de instrucciones. El primero se maneja afinando los parámetros de llamada al núcleo y el segundo, mediante la creación de kernels lo más sencillos que se pueda (sin perder de vista el tiempo de latencia). Podemos saber el uso de registros que hace nuestro kernel añadiendo el parámetro (31) – Xptxas –v en la línea de compilación.

nvcc -Xptxas –v ejemplo.cu

ptxas info : Compiling entry function 'acos_main'

ptxas info : Used 4 registers, 60+56 bytes lmem, 44+40 bytes smem,

20 bytes cmem[1],
12 bytes cmem[14]

Figura 2. Obteniendo información de nuestro kernel

[1] De panacea universal: Remedio que buscaban los antiguos alquimistas para curar todas las enfermedades

[2] El término “Warp” puede asimilarse al de “tanda” ya que hace referencia al conjunto de hilos de un procesador que se ejecutan en el mismo instante por pertenecer a la misma programación del scheduler.

Unroll y Vencerás

2012-01-14T09:59:00.000-08:00

Este, es el primero de una serie de 5 post en los cuales quiero presentar las claves básicas para desarrollar kernels CUDA de calidad.

Desde hace décadas, uno de los caballos de batalla en el mundo de los compiladores para CPU es el análisis de bucles y la reducción de iteraciones con el objetivo de optimizar el uso de recursos y mejorar los tiempos de ejecución. El resultado suele ser un bucle de menor número de repeticiones con instrucciones paralelas que son independientes unas de otras, en muchos casos esta simple acción puede reducir drásticamente el tiempo de procesamiento. El número de instrucciones serializadas, dentro del nuevo bucle coincide con el factor de desenrollado “unroll factor”. Este valor es clave a la hora de afinar y sintonizar los mejores valores para la paralelización como veremos más adelante. La Figura 1 muestra un ejemplo sencillo de esta técnica para UF = 4.

/* Antes del desenrollado */

for (i = 0; i < N; ++i) {

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

}

/* Despues del desenrollado */

for (i = 0; i < N - (4 - 1); i += 4) {

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

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

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

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

}

/* Bucle para el resto */

for (; i < N; ++i) {

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

}

Figura 1. Ejemplo de "unroll"

Los beneficios inmediatos de este proceso son:

1. Se reduce el contador de instrucciones dinámicas, se reducen el número de comparaciones y operaciones ramificadas, para la misma carga de trabajo.

2. Aumentan las posibles combinaciones para la planificación, debido a que surgen más instrucciones independientes, aumentando así el nivel de paralelización por instrucción (ILP).

3. Aumenta la oportunidad para explotar la herencia de los registros y memoria local cuando los bucles externos se desenrollan y los internos se fusionan en el mismo juego de instrucciones iguales para todos.

En las arquitecturas GPGPU esta técnica también es habitualmente usada, aunque con otro enfoque, la diferencia de arquitecturas y modelos de programación necesita de una nueva forma de contemplar el impacto del desenrollado de bucles en el rendimiento en los programas para GPGPU.

A diferencia de los compiladores para CPU, no es tan fácil poder realizar optimizaciones automáticas de código usando esta técnica, por lo que el programador es el que deberá diseñar un buen código paralelizado en GPU. El modelo de programación CUDA y las restricciones de los dispositivos nos dan la guía para realizar la paralelización de bucles, como no es el objetivo de este post exponer los pormenores de la arquitectura y su forma de uso, sólo haremos una breve reseña con objeto de aportar claridad a esta exposición, para ampliar conceptos puede referirse el lector a las guías que acompañan el toolkit de CUDA.

Básicamente podemos decir que un Kernel[1] equivale a una iteración del cuerpo de un bucle (en programación secuencial). La GPU tiene la habilidad de ejecutar de forma paralela el mismo kernel tantas veces como se defina en su llamada. Cada kernel corre en un hilo, y existen unas variables locales a cada uno de los kernels, pudiendo así cada kernel disponer de la información del índice o contador de iteración equivalente en caso de los bucles secuenciales. Los hilos (Thread), están organizados en bloques (Block), y los bloques en mallas (Grid). Los hilos pueden llegar a tener hasta un máximo de 3 dimensiones, y los bloques y las mallas hasta 2 o 3 dependiendo de la arquitectura. La siguiente línea de código muestra como un kernel puede obtener su identificador único (utilizando una dimensión en cada nivel):

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

Cómo pueden compartir y acceder a la memoria cada uno de los hilos y su ámbito se describe en el post "Los diversos usos de la memoria".

El tamaño de la malla a utilizar y del bloque de hilos, es una elección que se deja al programador en función del problema a resolver. Esta decisión tendrá impacto en el rendimiento final del programa, que como veremos a continuación, está ligada al problema y sujeta a las restricciones que impone la arquitectura GPU y los recursos propios del dispositivo que se utilice.

[1] El Kernel es una porción de código que es invocada desde un proceso en CPU (host) y se ejecuta en un dispositivo GPU.

CUDA vs OpenCL. ¿Cual es mejor?

2011-12-30T09:36:00.000-08:00

A la hora de seleccionar una plataforma u otra, deberemos tener en mente qué es exactamente lo que queremos implementar y qué mejoras deseamos obtener.

Por un lado la disponibilidad de plataformas puede ser el único indicador para realizar la elección: CUDA solo está disponible para la arquitectura de nVidia. Así que, esta cuestión es la más importante.

Si puedes hacer una elección entre ambos, la siguiente cuestión es el objetivo de la paralelización, casi siempre suele ser el incremento del speedUp con respecto a la versión serie, o de CPU, si este es el caso, considera utilizar CUDA ya que es más eficiente tanto a la hora de compilar, obteniendo un código mejorado y optimizado, como a la hora de gestionar en tiempo de ejecución transferencias de memoria, estas dos cuestiones conducen a un comportamiento mejor con respecto a su código homólogo en OpenCL.

Puedes ampliar información a este respecto en este paper, en este otro estudio también se expone la mejora en la utilización de CUDA frente a OpenCL.

Un aspecto a destacar sobre OpenCL es que se trataría de un lenguaje portable, válido para cualquier arquitectura, si has leído los enlaces anteriores te habrás dado cuenta que no es oro todo lo que reluce. Desafortunadamente, es necesario realizar ciertos ajustes en el código cuando deseamos cambiar de plataforma, pero no son muy grandes, esto se debe a la capa de plataforma, que es la parte específica que gestiona el interface con cada arquitectura.

Si aceptas mi consejo, utiliza CUDA para mejorar el rendimiento y aplicaciones para la vida real, y mantén OpenCL cerca de ti, trata de experimentar y portar a este lenguaje tus desarrollos, algún día OpenCL será el estandar, aunque para ese entonces llegarán otros gigantes imponiendo su criterio. El grupo Khronos ya sabe de estos casos: ¿OpenGL?

Si el rendimiento no lo es todo, podrías directamente utilizar OpenCL, sobre todo si no cuentas con plataformas de nVidia, aunque también están disponibles las plataformas específicas de amd: su famoso ATI Stream o APP SDK

Espero que esto te ayude a tomar una decisión, cualquier camino que tomes seguro que te conducirá al éxito. Hasta el próximo post.

Estado del Arte III (CUDA)

2011-12-15T13:14:00.001-08:00

CUDA[1]es una arquitectura de cálculo paralelo de nVidia que aprovecha la gran potencia de la GPU (unidad de procesamiento gráfico) para proporcionar un incremento extraordinario del rendimiento del sistema. CUDA ha sido recibida con entusiasmo por la comunidad científica. Por ejemplo, se está utilizando para acelerar AMBER, un simulador de dinámica molecular empleado por más de 60.000 investigadores del ámbito académico y farmacéutico de todo el mundo para acelerar el descubrimiento de nuevos medicamentos. En el mercado financiero, Numerix y CompatibL introdujeron soporte de CUDA para una nueva aplicación de cálculo de riesgo de contraparte y, como resultado, se ha multiplicado por 18 la velocidad de la aplicación.

CUDA Toolkit es el conjunto de herramientas y librerías que comprenden el entorno software de la arquitectura. Los componentes básicos son: El compilador C/C++, herramientas para depurar, librerías matemáticas optimizadas para GPU, Primitivas de rendimiento optimizado para GPU, todo ello acompañado de una extensa documentación y abundantes ejemplos. El paquete está disponible para los sistemas operativos más habituales: Windows, Linux y MacOS. Además el desarrollador cuenta con un enorme repertorio de servicios en línea, Webminars, whitepapers y numerosas publicaciones electrónicas y como no, los foros de la comunidad de desarrolladores.

Versiones del toolkit: desde 2006 hasta la fecha han aparecido varias versiones de CUDA que han supuesto un avance en cada una de ellas, la Figura siguiente muestra dicha evolución a lo largo del tiempo.

Evolución de CUDA

En mayo de 2011 se ha publicado la versión 4.0 del toolkit, el cual supone grandes avances en la facilidad de uso y programación, entre otras destacan las siguientes novedades: Se simplifica la portabilidad de aplicaciones: se pueden compartir múltiples GPUs a través de múltiples hilos, se pueden utilizar todas las GPUs de forma concurrente desde el mismo hilo del host. Reserva de memoria de host más rápida, implementación de las funciones C++ new/delete y soporte para funciones virtuales, soporte para código ensamblador inline PTX, librerías para la gestión de vectores y matrices optimizadas, y más.

El siguiente post, tratará de argumentar la elección de CUDA como plataforma de desarrollo frente a otras alternativas existentes.

[1] son las siglas de Compute Unified Device Architecture

Estado del Arte II (Software)

2011-12-10T12:56:00.001-08:00

A. Librerías para GPGPU

Han ido de la mano a la aparición del hardware y la presentación de las nuevas librerías gráficas y SDK para el manejo del nuevo hardware. Por un lado los fabricantes han sacado al mercado soluciones específicas para su arquitectura, por otro lado, nos encontramos con la iniciativa del grupo Khronos para crear un lenguaje común para desarrollo de aplicaciones que exploten las tecnologías GPGPU.

Así AMD en 2006 comenzaría con su anuncio de CTM™ ("Close To Metal"). Se trata de un SDK para uso de su arquitectura FireStream recién estrenada por esa misma época. Esta interfaz proporciona a los desarrolladores un completo acceso al conjunto de instrucciones nativas y la memoria de los elementos de cálculo masivamente paralelos de las GPUs de AMD. Con CTM, las GPUs se convierten en poderosas arquitecturas abiertas programables, como las CPUs de hoy. Al liberar la arquitectura, CTM ofrece a los desarrolladores de hardware, un acceso de bajo nivel, determinista y repetible necesario para desarrollar las herramientas esenciales como compiladores, depuradores, bibliotecas de matemáticas y plataformas de aplicaciones.

En diciembre de 2007 y como continuación con su compromiso con un entorno completo de desarrollo de software que libere el sorprendente poder informático de las GPUs de AMD, AMD presentó ATI Stream SDK. Con ATI Stream SDK, AMD añadió un nuevo lenguaje de alto nivel, denominado ATI Brook+. CTM evolucionó a ATI CAL (nivel de abstracción de cálculo), el nivel de API compatible con Brook+. La aparición de Brook+ significó que AMD fue capaz de ofrecer una completa pila de desarrollo completamente gratis para el desarrollo de ATI Stream.

En junio de 2008, AMD y varios líderes del sector de GPGPUs, así como otras tecnologías de aceleradores, formaron el grupo de trabajo OpenCL bajo la denominación The Khronos Group. Khronos es artífice de otras especificaciones bien conocidas como OpenGL, y realizó la apuesta adecuada al enfocarse en la especificación OpenCL.

El 9 de diciembre de 2008 The Khronos Group anunció la puesta en el mercado de la especificación OpenCL 1.0. Mientras tanto, AMD también integró bibliotecas de rutinas para CAL en la suite del controlador ATI Catalyst, liberando las capacidades de aceleración de ATI Stream ya incorporadas en millones de tarjetas gráficas ATI Radeon. A principios de 2009, inmediatamente después de la presentación de la especificación OpenCL 1.0, AMD anunció su intención de adoptar rápidamente el estándar de programación OpenCL 1.0 e integrar un compilador y rutinas compatibles en su ATI Stream SDK.

nVidia desde la presentación en 2006 de la GeForce 8800 GTX comenzó la senda de CUDA, que se propone de un entorno software de librerías especializadas para extraer el máximo rendimiento al hardware. Evidentemente este entorno implica una dependencia del fabricante para obtener nuevas características y prestaciones, se discutirá todo esto y más en el siguiente post.

Una iniciativa que recientemente ha terminado su andadura es Sh o Libsh , se trata de un lenguaje de alto nivel para programación en GPU, con una sintaxis similar a C++, es el resultado de la investigación del laboratorio de gráficos por ordenador de la universidad de Waterloo. Concretamente el lenguaje fue concebido por Michael McCool e implementado por varios miembros del laboratorio.

En el próximo post hablaremos más detalladamente de la evolución de CUDA, hasta entonces!

Estado del Arte I (Hardware)

2011-12-08T12:40:00.006-08:00

Hoy trataré de hacer una panorámica sobre el estado del arte en paralelización GPGPU tanto en software como en hardware.

Actualmente los fabricantes de tarjetas graficas producen para o están agrupados en 4 rúbricas: AMD[1]/ATI, nVidia, Intel e IBM. Estos cuatro gigantes dominan en la actualidad el mercado de las GPGPUs[2] cada uno de ellos propone su propia arquitectura que va normalmente acompañada de soluciones software especificas para explotar al máximo las capacidades concretas de cada una. Los dos primeros con tendencia a desarrollar dispositivos GPU puros, mientras que los dos segundos se decantan por soluciones hibridas CPU+GPU, solo el futuro revelará cual ha sido la mejor estrategia.

1. AMD: En noviembre de 2006, AMD revolucionó las GPGPU con la presentación de CTM™ ("Close To Metal"), AMD anunció también la primera generación de AMD FireStream como la primera solución de procesamiento a través de hardware comercialmente disponible. Se trataba de un procesador basado en PCI Express con 1 GB de memory GDDR3 para computación con características de procesamiento de alto rendimiento optimizado. En 2007, AMD presentó también su solución de procesamiento de segunda generación, la AMD FireStream 9170. FireStream 9170 revolucionó la industria de GPGPU al ser la primera solución de procesamiento que ofreció compatibilidad de puntos flotantes en hardware.

Esta GPU alcanza más de 5 GFLOPS/W en precisión normal, y 1 GFLOPS/W en operaciones con doble precisión, 2 GB de memoria GDDR3, 320 núcleos. En agosto de 2008, AMD presentó el primer producto de su línea de soluciones de procesamiento de tercera generación, la AMD FireStream 9250. FireStream 9250 fue la primera solución de procesamiento que rompió la barrera de 1 TFLOPS y además, lo hizo en el espacio de una única ranura PCI-Express. Las soluciones de otros proveedores logran menos rendimiento a un mayor consumo y en mayores tamaños. En el primer trimestre de 2009, AMD anunció AMD FireStream 9270, que satisface las necesidades de los usuarios de alto rendimiento. FireStream 9270 proporciona 1,2 TFLOPS de increíble rendimiento a un consumo típico de 160 W. Destacan sus 800 núcleos de procesamiento, 2 GB de GDDR5, un ancho de banda de memoria de 108,8 GB/s, con un consumo medio de 160 watios.

2. nVidia.

En noviembre de 2006, nVidia presentó al mercado la GPU GeForce 8800 GTX con soporte para DirectX 10, supone la adopción de una nueva arquitectura basada en el procesador. Dotada de 138 procesadores y 768MB de memoria de video con un ancho de banda de memoria de 86.4 GB/s.

Esta nueva arquitectura vendría acompañada de la presentación de CUDA®. A partir de este momento, nVidia trabajará duro para obtener GPU’s para computación de propósito general, la Figura siguiente, muestra la evolución de sus GPGPU y la proyección futura de la compañía.

Evolución de la arquitectura de nVidida

La gama Tesla aparece en 2007, la familia de productos Tesla se ha diseñado para el cálculo paralelo y ofrece funciones de computación exclusivas:
Rendimiento en operaciones de coma flotante de doble precisión: 515 GigaFlops en productos Tesla C2050, M2050 y S2070, Comunicación PCIe más rápida: con dos motores DMA para comunicación PCIe bidireccional. Mayor rendimiento en aplicaciones técnicas con grandes volúmenes de datos: mayor memoria en placa (3 GB y 6 GB).

En 2009 aparece la tecnología Fermi que convertirá el coprocesamiento distribuido entre la GPU y la CPU en algo generalizado porque es capaz de manejar todo tipo de aplicaciones de cálculo en cualquier ámbito de la ciencia o la industria. Fermi incorpora algunas innovaciones importantes: 512 núcleos CUDA, Tecnología Parallel DataCache™ de NVIDIA, Motor GigaThread™ 3.0 de NVIDIA, Soporte completo de ECC.

Se prevé que para finales de 2011 o primeros de 2012, estén disponibles las primeras GPU con tecnología Kepler. Kepler es el sobrenombre utilizado por nVidia para referirse a su núcleo gráfico de nueva generación y, al igual que las GPUs Southern Islands de AMD, estas soluciones serán fabricadas utilizando el nodo de producción 28nm. nVidia prometió que la GPU Kepler y su sucesora, Maxwell, incluirían espacio de memoria virtual y soporte preferente, además de una serie de tecnologías que mejorarán la capacidad de la GPU de procesar datos sin la ayuda del procesador del sistema. Para un horizonte algo más lejano 2013-2014, la compañía tiene previsto lanzar la primera arquitectura basada en 22/20nm, su nombre en clave es "Maxwell".

En respuesta a la tendencia de la competencia en el diseño de soluciones para computación de altas prestaciones basadas en CPU+GPU, se manejan proyectos que pretender sacar al mercado en 2013 – 2015 sus nombres en clave son “Denver” y “Echelon”.

3. Intel

Larrabee es el nombre en código para un chip GPGPU que Intel está desarrollando separadamente de su actual línea de aceleradores de gráficos integrados. El chip iba a ser lanzado en 2010 como el núcleo de una tarjeta gráfica 3D para los consumidores, pero estos planes fueron cancelados debido a los retrasos y las primeras decepcionantes cifras de rendimiento. Larrabee ahora será lanzado como una plataforma para la investigación y el desarrollo de los gráficos por ordenador y HPC. Una versión futura de Larrabee eventualmente podría llegar a integrarse en una tarjeta de gráficos de consumo, pero Intel no ha discutido los planes específicos.

4. IBM

Con el nombre de Proyecto Cell , junto a Sony y Toshiba, el gigante azul lanzó una tecnología híbrida para ofrecer Computación de Alto Rendimiento a la industria y la investigación. Con mayor o menor fortuna esta arquitectura no está liderando actualmente las tendencias, y la propia firma está montado soluciones de nVidia para sus servidores de altas prestaciones.

En el próximo post, analizaremos la situación actual del software para GPGPU.

[1] Advanced Micro Devices, Inc.

[2] General Porpouse Graphics Processor Unit

Bienvenidos a desarrollo CUDA

2011-12-01T12:11:00.000-08:00

Bienvenidos!

Con este blog quiero compartir mi experiencia paralelizando algoritmos mediante el uso de GPU's. Trataré de presentar aportes cortos que permitan obtener soluciones concretas.

La computación sobre Unidades de Procesamiento Gráfico (GPU), es sin lugar a dudas, uno de los avances tecnológicos más importantes de los últimos años, que aporta un nuevo paradigma en el campo de la computación en paralelo y la computación de altas prestaciones. Este tipo de hardware (GPU), inicialmente diseñado para representar imágenes en el ordenador, tiene unas características peculiares que ofrecen una capacidad de cómputo inusitada. A diferencia de los actuales procesadores centrales (CPU) en los cuales podemos encontrar configuraciones de unos cuantos procesadores con 2, 4 u 8 núcleos de proceso, las GPU ofrecen configuraciones de hasta 30 procesadores con hasta 8 núcleos de procesamientos cada uno, permitiendo ejecutar hasta 1536 hilos cada procesador.

Este potencial inicialmente era explotado por librerías gráficas de renderizado y trazadores de rayos para obtener imágenes realistas aplicando geometría computacional. El objetivo fundamental era obtener ciertos los valores para cada punto (en el plano o en el espacio para sistemas tridimensionales) como son: visible o no visible, color y profundidad de color. Así, con objeto de producir imágenes en movimiento se diseñaron sistemas que permitiesen realizar el cálculo de esos valores para cada punto de forma paralela y obtener el resultado de la imagen en tiempos tales que permitan la producción de una ilusión óptica de movimiento al presentar imágenes con una frecuencia mínima de 24 imágenes por segundo. Se trataba de poder calcular todos los valores de la matriz de puntos que forman la imagen en paralelo y de forma independiente del procesador central el cual estaría realizando otro tipo de tareas. Este conjunto de técnicas fueron denominados como Aceleración Gráfica por Hardware. Para poder explotar estas capacidades del nuevo hardware, deberá haber asociado a éste, un controlador que permita un repertorio específico de instrucciones, que serán llamadas utilizando un determinado API[1]. Términos tales como OpenGL® o DitectX® empezaron a sonar como librerías gráficas que eran capaces de explotar las nuevas capacidades de las GPUs.

Con el paso del tiempo cada vez fueron mayores las prestaciones de estas GPUs, aumentando tanto la capacidad de memoria como de los procesadores. Esto permitió aumentar el tipo de operaciones a realizar en dichos dispositivos de modo que además de poder realizar operaciones de dibujo, también se empezaron a realizar operaciones de codificación y decodificación, lo que permitió liberar a los procesadores centrales de tales tareas.

El empleo de procesadores gráficos para la resolución de problemas, comienza en el mundo académico en un principio dando forma gráfica al problema a solucionar, ya que las primitivas y repertorio de instrucciones con el que se contaba eran muy básicos y demasiado orientados al mundo de los gráficos. Actualmente se ha conseguido evolucionar estas librerías para proporcionar soporte a sistemas de desarrollo de propósito general que permiten aplicar con gran éxito la computación sobre GPU a cualquier tipo de problema sea o no de tipo gráfico.

[1] Application Programming Interface