Nota generada a partir de liga
Notas para contenedor de docker:
Comando de docker para ejecución de la nota de forma local:
nota: cambiar <ruta a mi directorio>
por la ruta de directorio que se desea mapear a /datos
dentro del contenedor de docker.
docker run --gpus all --rm -v <ruta a mi directorio>:/datos --name jupyterlab_nvidia_cuda_c_container -p 8888:8888 -d palmoreck/jupyterlab_nvidia_cuda_c:1.1.0_10.2
password para jupyterlab: qwerty
Detener el contenedor de docker:
docker stop jupyterlab_nvidia_cuda_c_container
Documentación de la imagen de docker palmoreck/jupyterlab_nvidia_cuda_c:1.1.0_10.2
en liga.
Esta nota utiliza métodos vistos en 1.5.Integracion_numerica
Nota: si se desean ejecutar los ejemplos que se presentan a continuación de forma local, es necesario tener una tarjeta gráfica NVIDIA.
Para ejecutar el notebook si están usando google colab:
Ver googlecolab github para información sobre google colab.
Consiste en extensiones al lenguaje C y en una runtime library. Ver 2.3.CUDA para más información.
__global__ void mifun(int param){
...
}
y siempre es tipo void
(no hay return
).
<<< >>>
y en la primera entrada se coloca el número de CUDA blocks y en la segunda entrada el número de CUDA threads:__global__ void mifun(int param){
...
}
int main(){
int par=0;
mifun<<<N,5>>> (par); //N bloques de 5 threads
}
* Los nombres de host y device hacen referencia a la CPU y a la GPU.
In [1]:
%%file hello_world.cu
#include<stdio.h>
__global__ void func(void){
}
int main(void){
func<<<1,1>>>(); //1 bloque de 1 thread
printf("Hello world!\n");
return 0;
}
Compilación:
In [2]:
%%bash
nvcc --compiler-options -Wall hello_world.cu -o hello_world.out
Ejecución:
In [3]:
%%bash
./hello_world.out
Comentarios:
La función main
se ejecuta en la CPU.
func
es un kernel y es ejecutada por los CUDA threads en el device (GPU). Obsérvese que tal función inicia con la sintaxis __global__
. En este caso el CUDA thread que fue lanzado no realiza ninguna acción pues el cuerpo del kernel está vacío.
El kernel sólo puede tener un return
tipo void: __global__ void func
por lo que el kernel debe regresar sus resultados a través de sus argumentos.
nvcc
es un wrapper para el compilador de programas escritos en C
. El compilador instalado en el contenedor de docker descrito al inicio de ésta nota es gcc
.
La extensión del archivo debe ser .cu
aunque esto puede modificarse al compilar con nvcc
:
$nvcc -x cu hello_world.c -o hello_world.out
-arch=sm_11
en la línea de nvcc
. En este caso se le indica al compilador que compile el programa para un compute capability de $1.1$. Ver liga y liga2 para más sobre esto.Los CUDA threads son divididos los CUDA blocks y éstos se encuentran en un grid. En el lanzamiento del kernel se debe especificar al hardware cuántos CUDA blocks tendrá nuestro grid y cuántos CUDA threads estarán en cada bloque.
In [4]:
%%file hello_world_2.cu
#include<stdio.h>
__global__ void func(void){
printf("Hello world del bloque %d del thread %d!\n", blockIdx.x, threadIdx.x);
}
int main(void){
func<<<2,3>>>(); //2 bloques de 3 threads cada uno
cudaDeviceSynchronize();
printf("Hola del cpu thread\n");
return 0;
}
In [5]:
%%bash
nvcc --compiler-options -Wall hello_world_2.cu -o hello_world_2.out
In [6]:
%%bash
./hello_world_2.out
Comentarios:
El llamado a la ejecución del kernel se realizó en el host (CPU) y se lanzaron $2$ bloques (primera posición en la sintaxis <<<>>>), cada uno con $3$ threads.
Se utiliza la función cudaDeviceSynchronize para que el cpu-thread espere la finalización de la ejecución del kernel.
En el ejemplo anterior, las variables blockIdx
y threadIdx
hacen referencia a los id's que tienen los bloques y los threads: el id del bloque dentro del grid y el id del thread dentro del bloque. La parte .x
de las variables: blockIdx.x
y threadIdx.x
refieren a la primera coordenada del bloque en el grid y a la primera coordenada del thread en en el bloque.
La elección del número de bloques en un grid o el número de threads en un bloque no corresponde a alguna disposición del hardware, esto es, si se lanza un kernel con <<< 1, 3 >>>
no implica que la GPU tenga en su hardware un bloque o 3 threads. Asimismo, las coordenadas que se obtienen vía blockIdx
o threadIdx
son meras abstracciones, no corresponden a algún ordenamiento en el hardware de la GPU.
Todos los threads de un bloque ejecutan el kernel por lo que se tienen tantas copias del kernel como número de bloques sean lanzados. Aquí encontramos que en una GPU se tiene el modelo Single Instruction Multiple Threads SIMT.
En una GPU podemos definir el grid de bloques y el bloque de threads utilizando el tipo de dato dim3
el cual también es parte de CUDA C:
In [7]:
%%file hello_world_3.cu
#include<stdio.h>
__global__ void func(void){
printf("Hello world del bloque %d del thread %d!\n", blockIdx.y, threadIdx.z);
}
int main(void){
dim3 dimGrid(1,2,1); //2 bloques en el grid
dim3 dimBlock(1,1,3); //3 threads por bloque
func<<<dimGrid,dimBlock>>>();
cudaDeviceSynchronize();
printf("Hola del cpu thread\n");
return 0;
}
In [8]:
%%bash
nvcc --compiler-options -Wall hello_world_3.cu -o hello_world_3.out
In [9]:
%%bash
./hello_world_3.out
In [10]:
%%file thread_idxs.cu
#include<stdio.h>
__global__ void func(void){
if(threadIdx.x==0 && threadIdx.y==0 && threadIdx.z==0){
printf("blockIdx.x:%d\n",blockIdx.x);
}
printf("thread idx.x:%d\n",threadIdx.x);
printf("thread idx.y:%d\n",threadIdx.y);
printf("thread idx.z:%d\n",threadIdx.z);
}
int main(void){
dim3 dimGrid(1,1,1); //1 bloque en el grid
dim3 dimBlock(1,3,1); //3 threads por bloque
func<<<dimGrid,dimBlock>>>();
cudaDeviceSynchronize();
return 0;
}
In [11]:
%%bash
nvcc --compiler-options -Wall thread_idxs.cu -o thread_idxs.out
In [12]:
%%bash
./thread_idxs.out
In [13]:
%%file block_idxs.cu
#include<stdio.h>
__global__ void func(void){
printf("blockIdx.x:%d\n",blockIdx.x);
printf("blockIdx.y:%d\n",blockIdx.y);
printf("blockIdx.z:%d\n",blockIdx.z);
}
int main(void){
dim3 dimGrid(1,2,2); //4 bloques en el grid
dim3 dimBlock(1,1,1); //1 thread por bloque
func<<<dimGrid,dimBlock>>>();
cudaDeviceSynchronize();
return 0;
}
In [14]:
%%bash
nvcc --compiler-options -Wall block_idxs.cu -o block_idxs.out
In [15]:
%%bash
./block_idxs.out
Podemos usar la variable blockDim
para cada coordenada x, y
o z
y obtener la dimensión de los bloques:
In [16]:
%%file block_dims.cu
#include<stdio.h>
__global__ void func(void){
if(threadIdx.x==0 && threadIdx.y==0 && threadIdx.z==0 && blockIdx.z==1){
printf("blockDim.x:%d\n",blockDim.x);
printf("blockDim.y:%d\n",blockDim.y);
printf("blockDim.z:%d\n",blockDim.z);
}
}
int main(void){
dim3 dimGrid(2,2,2); //8 bloques en el grid
dim3 dimBlock(3,1,2); //6 threads por bloque
func<<<dimGrid,dimBlock>>>();
cudaDeviceSynchronize();
return 0;
}
In [17]:
%%bash
nvcc --compiler-options -Wall block_dims.cu -o block_dims.out
In [18]:
%%bash
./block_dims.out
Para alojar memoria en el device se utiliza el llamado a cudaMalloc y para transferir datos del host al device o viceversa se llama a lafunción cudaMemcpy con respectivos parámetros como cudaMemcpyHostToDevice
o cudaMemcpyDeviceToHost
.
Para desalojar memoria del device se utiliza el llamado a cudaFree
N bloques de 1 thread
In [19]:
%%file suma_vectorial.cu
#include<stdio.h>
#define N 10
__global__ void suma_vect(int *a, int *b, int *c){
int block_id_x = blockIdx.x;
if(block_id_x<N) //aquí se asume que el valor de N
//es menor al número máximo de bloques que se pueden lanzar
//si fuese mayor, hay que hacer un ajuste
c[block_id_x] = a[block_id_x]+b[block_id_x];
}
int main(void){
int a[N], b[N],c[N];
int *device_a, *device_b, *device_c;
int i;
dim3 dimGrid(N,1,1); //N bloques en el grid
dim3 dimBlock(1,1,1); //1 threads por bloque
//alojando en device
cudaMalloc((void **)&device_a, sizeof(int)*N);
cudaMalloc((void **)&device_b, sizeof(int)*N);
cudaMalloc((void **)&device_c, sizeof(int)*N);
//llenando los arreglos con datos dummy:
for(i=0;i<N;i++){
a[i]=i;
b[i]=i*i;
}
//copiamos arreglos a, b a la GPU
cudaMemcpy(device_a,a,N*sizeof(int), cudaMemcpyHostToDevice);
cudaMemcpy(device_b,b,N*sizeof(int), cudaMemcpyHostToDevice);
//mandamos a llamar a suma_vect:
suma_vect<<<dimGrid,dimBlock>>>(device_a,device_b,device_c); //N bloques de 1 thread
cudaDeviceSynchronize();
//copia del resultado al arreglo c:
cudaMemcpy(c,device_c,N*sizeof(int),cudaMemcpyDeviceToHost);
for(i=0;i<N;i++)
printf("%d+%d = %d\n",a[i],b[i],c[i]);
cudaFree(device_a);
cudaFree(device_b);
cudaFree(device_c);
return 0;
}
In [20]:
%%bash
nvcc --compiler-options -Wall suma_vectorial.cu -o suma_vectorial.out
In [21]:
%%bash
./suma_vectorial.out
Comentarios:
int *device_a, *device_b, *device_c;
pero estos apuntadores no apuntan a una dirección de memoria en el device pues aunque NVIDIA añadió el feature de Unified Memory (un espacio de memoria accesible para el host y el device) aquí no se está usando tal feature. Más bien se están utilizando los apuntadores anteriores para apuntar a un struct) de C en el que uno de sus tipos de datos es una dirección de memoria en el device.
Obsérvese el uso de (void **)
por la definición de la función cudaMalloc
.
Obsérvese que en el programa anterior se coloca en comentario que se asume que $N$ el número de datos en el arreglo es menor al número de bloques que es posible lanzar. Esto como veremos más adelante es importante considerar pues aunque en un device se pueden lanzar muchos bloques y muchos threads, se tienen límites en el número de éstos que es posible lanzar.
Al instalar el CUDA toolkit en sus máquinas o bien si utilizan el contenedor de docker (descrito al inicio de la nota) se instala la línea de comando nvprof para perfilamiento. Se puede ejecutar con:
In [22]:
%%bash
nvprof --normalized-time-unit s ./suma_vectorial.out
Comentarios:
Las unidades en las que se reporta son s: second, ms: millisecond, us: microsecond, ns: nanosecond.
En la documentación de NVIDIA se menciona que nvprof
será reemplazada próximamente por NVIDIA Nsight Compute y NVIDIA Nsight Systems.
En el ejemplo anterior se lanzaron $N$ bloques con $1$ thread cada uno y a continuación se lanza $1$ bloque con $N$ threads.
In [23]:
%%file suma_vectorial_2.cu
#include<stdio.h>
#define N 10
__global__ void suma_vect(int *a, int *b, int *c){
int thread_id_x = threadIdx.x;
if(thread_id_x<N) //aquí se asume que el valor de N
//es menor al número máximo de threads que se pueden lanzar
//si fuese mayor, hay que hacer un ajuste
c[thread_id_x] = a[thread_id_x]+b[thread_id_x];
}
int main(void){
int *device_a, *device_b, *device_c;
int i;
dim3 dimGrid(1,1,1); //1 bloques en el grid
dim3 dimBlock(N,1,1); //N threads por bloque
//alojando en device con Unified Memory
cudaMallocManaged(&device_a, sizeof(int)*N);
cudaMallocManaged(&device_b, sizeof(int)*N);
cudaMallocManaged(&device_c, sizeof(int)*N);
//llenando los arreglos:
for(i=0;i<N;i++){
device_a[i]=i;
device_b[i]=i*i;
}
suma_vect<<<dimGrid,dimBlock>>>(device_a,device_b,device_c); //1 bloque con N threads
cudaDeviceSynchronize();
for(i=0;i<N;i++)
printf("%d+%d = %d\n",device_a[i],device_b[i],device_c[i]);
cudaFree(device_a);
cudaFree(device_b);
cudaFree(device_c);
return 0;
}
In [24]:
%%bash
nvcc --compiler-options -Wall suma_vectorial_2.cu -o suma_vectorial_2.out
In [25]:
%%bash
nvprof --normalized-time-unit s ./suma_vectorial_2.out
Comentarios:
El programa anterior utiliza la Unified Memory con la función cudaMallocManaged. La Unified Memory es un feature que se añadió a CUDA desde las arquitecturas de Kepler y Maxwell pero que ha ido mejorando (por ejemplo añadiendo page faulting and migration) en las arquitecturas siguientes a la de Kepler: la arquitectura Pascal y Volta. Por esto en el output anterior de nvprof aparece una sección de page fault.
Obsérvese que en el programa anterior se comenta que se asume que $N$ el número de datos en el arreglo es menor al número de threads que es posible lanzar. Esto como veremos más adelante es importante considerar pues aunque en el device se pueden lanzar muchos bloques y muchos threads, se tienen límites en el número de éstos que es posible lanzar.
En realidad no tenemos que realizarlo para el ejemplo de suma_vectorial.cu
anterior. Por ejemplo:
In [26]:
%%file suma_vectorial_3.cu
#include<stdio.h>
#define N 10
__global__ void init(int *a, int *b){
int thread_id_x = threadIdx.x;
a[thread_id_x]=thread_id_x;
b[thread_id_x]=thread_id_x*thread_id_x;
}
__global__ void suma_vect(int *a, int *b, int *c){
int thread_id_x = threadIdx.x;
if(thread_id_x<N) //aquí se asume que el valor de N
//es menor al número máximo de threads que se pueden lanzar
//si fuese mayor, hay que hacer un ajuste
c[thread_id_x] = a[thread_id_x]+b[thread_id_x];
}
int main(void){
int *device_a, *device_b, *device_c;
int i;
dim3 dimGrid(1,1,1); //1 bloques en el grid
dim3 dimBlock(N,1,1); //N threads por bloque
//alojando en device con Unified Memory
cudaMallocManaged(&device_a, sizeof(int)*N);
cudaMallocManaged(&device_b, sizeof(int)*N);
cudaMallocManaged(&device_c, sizeof(int)*N);
//llenando los arreglos:
init<<<dimGrid,dimBlock>>>(device_a,device_b);
cudaDeviceSynchronize();
//llamando al kernel suma_vect
suma_vect<<<dimGrid,dimBlock>>>(device_a,device_b,device_c); //1 bloque con N threads
cudaDeviceSynchronize();
for(i=0;i<N;i++)
printf("%d+%d = %d\n",device_a[i],device_b[i],device_c[i]);
cudaFree(device_a);
cudaFree(device_b);
cudaFree(device_c);
return 0;
}
In [27]:
%%bash
nvcc --compiler-options -Wall suma_vectorial_3.cu -o suma_vectorial_3.out
In [28]:
%%bash
nvprof --normalized-time-unit s ./suma_vectorial_3.out
Un device está compuesto por arreglos de streaming multiprocessors SM's (también denotados como MP's) y en cada SM encontramos un número (determinado por la arquitectura del device) de streaming processors SP's que comparten el caché y unidades de control (que están dentro de cada SM):
En el dibujo anterior se muestran las SM's en color rojo y los SP's en morado. Hay dos SM's por cada bloque anaranjado y ocho SP's por cada SM. Así, una GPU es una máquina multicore. Aunque cada SM ejecuta las instrucciones de forma independiente a otra SM, comparten la memoria global.
Los bloques de threads son asignados a cada SM por el CUDA runtime system, puede asignar más de un bloque a una SM pero hay un límite de bloques que pueden ser asignados a cada SM. Ver maximum number of blocks per multiprocessor.
Comentarios:
Para responder lo anterior se puede utilizar el siguiente programa que está basado en liga y cudaDeviceProp Struct Reference:
In [29]:
%%file device_properties.cu
#include<stdio.h>
int main(void){
cudaDeviceProp properties;
int count;
int i;
cudaGetDeviceCount(&count);
for(i=0;i<count;i++){
printf("----------------------\n");
cudaGetDeviceProperties(&properties, i);
printf("----device %d ----\n",i);
printf("Device Name: %s\n", properties.name);
printf("Compute capability: %d.%d\n", properties.major, properties.minor);
printf("Clock rate: %d\n", properties.clockRate);
printf("Unified memory: %d\n", properties.unifiedAddressing);
printf(" ---Memory Information for device %d (results on bytes)---\n", i);
printf("Total global mem: %ld\n", properties.totalGlobalMem);
printf("Total constant Mem: %ld\n", properties.totalConstMem);
printf("Shared memory per thread block: %ld\n", properties.sharedMemPerBlock);
printf("Shared memory per SM: %ld\n",properties.sharedMemPerMultiprocessor );
printf(" ---MP Information for device %d ---\n", i);
printf("SM count: %d\n", properties.multiProcessorCount);
printf("Threads in warp: %d\n", properties.warpSize);
printf("Max threads per SM: %d\n", properties.maxThreadsPerMultiProcessor);
printf("Max warps per SM: %d\n",properties.maxThreadsPerMultiProcessor/properties.warpSize);
printf("Max threads per block: %d\n", properties.maxThreadsPerBlock);
printf("Max thread dimensions: (%d, %d, %d)\n", properties.maxThreadsDim[0], properties.maxThreadsDim[1], properties.maxThreadsDim[2]);
printf("Max grid dimensions: (%d, %d, %d)\n", properties.maxGridSize[0], properties.maxGridSize[1], properties.maxGridSize[2]);
}
return 0;
}
In [30]:
%%bash
nvcc --compiler-options -Wall device_properties.cu -o device_properties.out
In [20]:
%%bash
./device_properties.out
Comentarios:
También en la documentación oficial de NVIDIA dentro de compute-capabilities se pueden revisar los valores anteriores y muchos más.
Obsérvese del output anterior que el sistema tiene dos devices con las mismas capacidades.
En un device encontramos diferentes tipos de memoria: global, constante, shared y texture. En esta nota únicamente trabajamos con la memoria global. Tenemos funciones en CUDA para poder comunicar/coordinar a los threads en un bloque por medio de la shared memory. Ver por ejemplo Using Shared Memory in CUDA C/C++ para un pequeño post del $2013$ sobre shared memory.
Los bloques de threads que son asignados a una SM son divididos en warps que es la unidad de thread scheduling* que tiene el CUDA run time system. El output anterior indica que son divisiones de $32$ threads: un warp consta de $32$ threads.
* El thread scheduling se puede pensar a la funcionalidad del hardware para seleccionar una instrucción del programa y asginar su ejecución por los threads en un warp (SIMT). Otro ejemplo es tener una instrucción que indica que se debe realizar lectura o escritura, entonces el hardware del device utiliza un warp de threads para tal operación mientras selecciona un warp de threads distinto para seleccionar otra instrucción diferente a la de I/O.
El número máximo de threads que pueden iniciarse de forma simultánea o en un instante por SM es de $2048$ o bien $2048/32 = 64$ warps.
El output anterior muestra los límites para número de bloques en las tres dimensiones de un grid y el número de threads en las tres dimensiones en un bloque.
Un bloque puede tener como máximo $1024$ threads en cualquier configuración: por ejemplo $(1024,1,1), (32,1,32), (4,4,64)$.
Por los puntos anteriores si lanzamos bloques de $1024$ threads entonces sólo $2$ bloques pueden residir en una SM en un instante. Con esta configuración alcanzaríamos $1024/32=32$ warps por cada bloque y como lanzamos $2$ bloques alcanzaríamos $64$ warps (que es el máximo de warps por SM que podemos tener en un instante). Otra configuración para alcanzar el máximo número de warps en un instante, es considerar $4$ bloques de $512$ threads pues tendríamos $512/32=16$ warps por bloque y en total serían $16*4$ (warps $\times$ bloques) $=64$ warps. Entre los datos que hay que elegir en los programas de CUDA C se tienen las configuraciones en el número de threads y el número de bloques a lanzar. La idea es alcanzar o rebasar el máximo número de warps en cada SM que soporta nuestro device en un instante.
Por ejemplo para el dibujo en el que se asumió que el CUDA runtime system había asignado $3$ bloques a cada SM, se tendría una división de cada bloque en un warp de $32$ threads como sigue:
Los programas de CUDA C tienen la opción de elegir el número de threads y de bloques a ser lanzados. En la referencia "Parallel Computing for Data Science. With Examples in R, C++ and CUDA" de N. Matloff se enlistan algunas consideraciones para elegir tales parámetros:
Given that scheduling is done on a warp basis, block size should be a multiple of the warp size (32).
One wants to utilize all the SMs. If one sets the block size too large, not all will be used, as a block cannot be split across SM's.
..., barrier synchronization can be done effectively only at the block level. The larger the block, the more the barrier delay, so one might want smaller blocks.
On the other hand, if one is using shared memory, this can only be done at the block level, and efficient use may indicate using a larger block.
Two threads doing unrelated work, or the same work but with many if/elses, would cause a lot of thread divergence if they were in the same block. In some cases, it may be known in advance which threads will do the "ifs" and which will do the "elses", in which case they should be placed in different blocks if possible.
A commonly-cited rule of thumb is to have between $128$ and $256$ threads per block.
En el uso de CUDA se recomienda que:
Users escriban código de CUDA C simple.
Utilicen las librerías ya hechas por NVIDIA o terceros para mantener simplicidad y eficiencia en el código.
Lo anterior para disminuir el tiempo y la cantidad de código que users tengan que hacer (o rehacer) y puesto que dominar la programación de CUDA C requiere una buena inversión de tiempo.
Así, tenemos a Thrust una template library basada en la Standard Template Library (STL) de C++ construída por NVIDIA que de acuerdo a su documentación:
"Thrust provides a rich collection of data parallel primitives such as scan, sort, and reduce, which can be composed together to implement complex algorithms with concise, readable source code. By describing your computation in terms of these high-level abstractions you provide Thrust with the freedom to select the most efficient implementation automatically. As a result, Thrust can be utilized in rapid prototyping of CUDA applications, where programmer productivity matters most, as well as in production, where robustness and absolute performance are crucial."
Thrust
tiene la opción de utilizarse con OpenMP, Thread Building Blocks (TBB) y con CUDA-C++. Ver por ejemplo Device Backends para conocer cómo cambiar entre OpenMP y CUDA-C++, lo cual se realiza en la compilación y sin hacer cambios en el código!. A los sistemas de software que tienen este feature se les llama Heterogeneous computing.
Si se instala el CUDA toolkit o se utiliza la imagen de docker descrita al inicio de la nota, los headers en la librería template de Thrust
estarán disponibles para su uso.
En el siguiente ejemplo de la regla del rectángulo compuesta se utiliza:
Los headers:
y se hace explícito el uso de la política de ejecucion thrust::device.
Referencias para el programa siguiente se encuentran en thrust inside user written kernels y cuda how to sum all elements of an array into one number within the gpu de stackoverflow.
Los siguientes tiempos se calcularon con una GPU de las siguientes características:
In [31]:
%%bash
./device_properties.out
$n=10^3$ subintervalos
In [1]:
%%file Rcf.cu
#include<stdio.h>
#include <thrust/reduce.h>
#include <thrust/execution_policy.h>
__global__ void Rcf(double *data, double a, double h_hat, int n, double *sum ) {
/*
Compute numerical approximation using rectangle or mid-point method in
an interval.
Nodes are generated via formula: x_i = a+(i+1/2)h_hat for i=0,1,...,n-1 and h_hat=(b-a)/n
Args:
a (int): left point of interval
n (int): number of subintervals
data (double): array that will hold values evaluated in function
Returns:
sum (double): pointer to result
*/
double x=0.0;
if(threadIdx.x<=n-1){
x=a+(threadIdx.x+1/2.0)*h_hat;
data[threadIdx.x]=std::exp(-std::pow(x,2));
}
if(threadIdx.x==0){
*sum = thrust::reduce(thrust::device, data , data + n, (double)0, thrust::plus<double>());
}
}
int main(int argc, char *argv[]){
double suma=0.0;
double *d_data;
double *d_suma;
double a=0.0, b=1.0;
double h_hat;
int n=1e3; //número de subintervalos
double objetivo=0.7468241328124271;
double time_spent;
clock_t begin,end;
cudaMalloc((void **)&d_data,sizeof(double)*n);
cudaMalloc((void**)&d_suma,sizeof(double));
h_hat=(b-a)/n;
begin=clock();
Rcf<<<1,n>>>(d_data, a,h_hat,n,d_suma); //1 bloque de n threads
cudaDeviceSynchronize();
end=clock();
time_spent = (double)(end - begin) / CLOCKS_PER_SEC;
cudaMemcpy(&suma, d_suma, sizeof(double), cudaMemcpyDeviceToHost);
suma=h_hat*suma;
cudaFree(d_data) ;
cudaFree(d_suma) ;
printf("Integral de %f a %f = %1.15e\n", a,b,suma);
printf("Error relativo de la solución: %1.15e\n", fabs(suma-objetivo)/fabs(objetivo));
printf("Tiempo de cálculo en la gpu %.5f\n", time_spent);
return 0;
}
In [33]:
%%bash
nvcc --compiler-options -Wall Rcf.cu -o Rcf.out
In [34]:
%%bash
./Rcf.out
In [35]:
%%bash
nvprof --normalized-time-unit s ./Rcf.out
$n=1025$ subintervalos
In [36]:
%%file Rcf2.cu
#include<stdio.h>
#include <thrust/reduce.h>
#include <thrust/execution_policy.h>
__global__ void Rcf(double *data, double a, double h_hat, int n, double *sum ) {
/*
Compute numerical approximation using rectangle or mid-point method in
an interval.
Nodes are generated via formula: x_i = a+(i+1/2)h_hat for i=0,1,...,n-1 and h_hat=(b-a)/n
Args:
a (int): left point of interval
n (int): number of subintervals
data (double): array that will hold values evaluated in function
Returns:
sum (double): pointer to result
*/
double x=0.0;
if(threadIdx.x<=n-1){
x=a+(threadIdx.x+1/2.0)*h_hat;
data[threadIdx.x]=std::exp(-std::pow(x,2));
}
if(threadIdx.x==0){
*sum = thrust::reduce(thrust::device, data , data + n, (double)0, thrust::plus<double>());
}
}
int main(int argc, char *argv[]){
double suma=0.0;
double *d_data;
double *d_suma;
double a=0.0, b=1.0;
double h_hat;
int n=1025; //número de subintervalos
double objetivo=0.7468241328124271;
double time_spent;
clock_t begin,end;
cudaMalloc((void **)&d_data,sizeof(double)*n);
cudaMalloc((void**)&d_suma,sizeof(double));
h_hat=(b-a)/n;
begin=clock();
Rcf<<<1,n>>>(d_data, a,h_hat,n,d_suma); //1 bloque de n threads
cudaDeviceSynchronize();
end=clock();
time_spent = (double)(end - begin) / CLOCKS_PER_SEC;
cudaMemcpy(&suma, d_suma, sizeof(double), cudaMemcpyDeviceToHost);
suma=h_hat*suma;
cudaFree(d_data) ;
cudaFree(d_suma) ;
printf("Integral de %f a %f = %1.15e\n", a,b,suma);
printf("Error relativo de la solución: %1.15e\n", fabs(suma-objetivo)/fabs(objetivo));
printf("Tiempo de cálculo en la gpu %.5f\n", time_spent);
return 0;
}
In [37]:
%%bash
nvcc --compiler-options -Wall Rcf2.cu -o Rcf2.out
In [38]:
%%bash
./Rcf2.out
Obsérvese error relativo de $100\%$
¿Cómo lo arreglamos?
In [39]:
%%file Rcf3.cu
#include<stdio.h>
#include <thrust/reduce.h>
#include <thrust/execution_policy.h>
__global__ void Rcf(double *data, double a, double h_hat, int n, double *sum) {
/*
Compute numerical approximation using rectangle or mid-point method in
an interval.
Nodes are generated via formula: x_i = a+(i+1/2)h_hat for i=0,1,...,n-1 and h_hat=(b-a)/n
Args:
a (int): left point of interval
n (int): number of subintervals
data (double): array that will hold values evaluated in function
Returns:
sum (double): pointer to result
*/
double x=0.0;
int stride=0;
if(threadIdx.x<=n-1){
x=a+(threadIdx.x+1/2.0)*h_hat;
data[threadIdx.x]=std::exp(-std::pow(x,2));
}
if(threadIdx.x==0){
stride=blockDim.x;
x=a+(threadIdx.x+stride+1/2.0)*h_hat;
data[threadIdx.x+stride]=std::exp(-std::pow(x,2));
*sum = thrust::reduce(thrust::device, data , data + n, (double)0, thrust::plus<double>());
}
}
int main(int argc, char *argv[]){
double suma=0.0;
double *d_data;
double *d_suma;
double a=0.0, b=1.0;
double h_hat;
int n_threads_per_block=1024;
int n_bloques=2;
int n=1025;//número de subintervalos
double objetivo=0.7468241328124271;
double time_spent;
clock_t begin,end;
cudaMalloc((void **)&d_data,sizeof(double)*n);
cudaMalloc((void**)&d_suma,sizeof(double));
h_hat=(b-a)/n;
begin=clock();
Rcf<<<n_bloques,n_threads_per_block>>>(d_data, a,h_hat,n,d_suma);
cudaDeviceSynchronize();
end=clock();
time_spent = (double)(end - begin) / CLOCKS_PER_SEC;
cudaMemcpy(&suma, d_suma, sizeof(double), cudaMemcpyDeviceToHost);
suma=h_hat*suma;
cudaFree(d_data) ;
cudaFree(d_suma) ;
printf("Integral de %f a %f = %1.15e\n", a,b,suma);
printf("Error relativo de la solución: %1.15e\n", fabs(suma-objetivo)/fabs(objetivo));
printf("Tiempo de cálculo en la gpu %.5f\n", time_spent);
return 0;
}
In [40]:
%%bash
nvcc --compiler-options -Wall Rcf3.cu -o Rcf3.out
In [41]:
%%bash
./Rcf3.out
Pero en la propuesta anterior lanzamos $2*1024$ (bloques $\times$ número de threads) $=2048$ threads y sólo ocupamos $1025$ threads. Entonces podemos cambiar el código anterior para aprovechar los $2048$ threads como sigue:
In [42]:
%%file Rcf4.cu
#include<stdio.h>
#include <thrust/reduce.h>
#include <thrust/execution_policy.h>
__global__ void Rcf(double *data, double a, double h_hat, int n, double *sum) {
/*
Compute numerical approximation using rectangle or mid-point method in
an interval.
Nodes are generated via formula: x_i = a+(i+1/2)h_hat for i=0,1,...,n-1 and h_hat=(b-a)/n
Args:
a (int): left point of interval
n (int): number of subintervals
data (double): array that will hold values evaluated in function
Returns:
sum (double): pointer to result
*/
double x=0.0;
int stride=0;
int i;
stride=blockDim.x;
for(i=threadIdx.x;i<=n-1;i+=stride){
if(i<=n-1){
x=a+(i+1/2.0)*h_hat;
data[i]=std::exp(-std::pow(x,2));
}
}
if(threadIdx.x==0){
*sum = thrust::reduce(thrust::device, data , data + n, (double)0, thrust::plus<double>());
}
}
int main(int argc, char *argv[]){
double suma=0.0;
double *d_data;
double *d_suma;
double a=0.0, b=1.0;
double h_hat;
int n_threads_per_block=1024;
int n_bloques=2;
int n=n_threads_per_block*n_bloques;//número de subintervalos
double objetivo=0.7468241328124271;
double time_spent;
clock_t begin,end;
cudaMalloc((void **)&d_data,sizeof(double)*n);
cudaMalloc((void**)&d_suma,sizeof(double));
h_hat=(b-a)/n;
begin=clock();
Rcf<<<n_bloques,n_threads_per_block>>>(d_data, a,h_hat,n,d_suma);
cudaDeviceSynchronize();
end=clock();
time_spent = (double)(end - begin) / CLOCKS_PER_SEC;
cudaMemcpy(&suma, d_suma, sizeof(double), cudaMemcpyDeviceToHost);
suma=h_hat*suma;
cudaFree(d_data) ;
cudaFree(d_suma) ;
printf("Integral de %f a %f = %1.15e\n", a,b,suma);
printf("Error relativo de la solución: %1.15e\n", fabs(suma-objetivo)/fabs(objetivo));
printf("Tiempo de cálculo en la gpu %.5f\n", time_spent);
return 0;
}
In [43]:
%%bash
nvcc --compiler-options -Wall Rcf4.cu -o Rcf4.out
In [44]:
%%bash
./Rcf4.out
Pero todavía podemos hacer aún mejor sin el for:
In [45]:
%%file Rcf5.cu
#include<stdio.h>
#include <thrust/reduce.h>
#include <thrust/execution_policy.h>
__global__ void Rcf(double *data, double a, double h_hat, int n, double *sum ) {
/*
Compute numerical approximation using rectangle or mid-point method in
an interval.
Nodes are generated via formula: x_i = a+(i+1/2)h_hat for i=0,1,...,n-1 and h_hat=(b-a)/n
Args:
a (int): left point of interval
n (int): number of subintervals
data (double): array that will hold values evaluated in function
Returns:
sum (double): pointer to result
*/
double x=0.0;
int idx;
idx = blockIdx.x * blockDim.x + threadIdx.x;
if(idx<=n-1){
x=a+(idx+1/2.0)*h_hat;
data[idx]=std::exp(-std::pow(x,2));
}
if(idx==0){
*sum = thrust::reduce(thrust::device, data , data + n, (double)0, thrust::plus<double>());
}
}
int main(int argc, char *argv[]){
double suma=0.0;
double *d_data;
double *d_suma;
double a=0.0, b=1.0;
double h_hat;
int n_threads_per_block=1024;
int n_bloques=2;
double objetivo=0.7468241328124271;
int n=n_bloques*n_threads_per_block;//número de subintervalos
double time_spent;
clock_t begin,end;
cudaMalloc((void **)&d_data,sizeof(double)*n);
cudaMalloc((void**)&d_suma,sizeof(double));
h_hat=(b-a)/n;
begin = clock();
Rcf<<<n_bloques,n_threads_per_block>>>(d_data, a,h_hat,n,d_suma);
cudaDeviceSynchronize();
end = clock();
time_spent = (double)(end - begin) / CLOCKS_PER_SEC;
cudaMemcpy(&suma, d_suma, sizeof(double), cudaMemcpyDeviceToHost);
suma=h_hat*suma;
cudaFree(d_data) ;
cudaFree(d_suma) ;
printf("Integral de %f a %f = %1.15e\n", a,b,suma);
printf("Error relativo de la solución: %1.15e\n", fabs(suma-objetivo)/fabs(objetivo));
printf("Tiempo de cálculo en la gpu %.5f\n", time_spent);
return 0;
}
In [46]:
%%bash
nvcc --compiler-options -Wall Rcf5.cu -o Rcf5.out
In [47]:
%%bash
./Rcf5.out
In [48]:
%%bash
nvprof --normalized-time-unit s ./Rcf5.out
Para una visualización sobre la construcción del índice en el kernel utilizando blockDim.x*blockIdx.x + threadIdx.x
ver An Even Easier Introduction to CUDA.
¿Más nodos?
Para este caso, incrementamos el número de bloques:
In [67]:
%%file Rcf6.cu
#include<stdio.h>
#include <thrust/reduce.h>
#include <thrust/execution_policy.h>
__global__ void Rcf(double *data, double a, double h_hat, int n, double *sum ) {
/*
Compute numerical approximation using rectangle or mid-point method in
an interval.
Nodes are generated via formula: x_i = a+(i+1/2)h_hat for i=0,1,...,n-1 and h_hat=(b-a)/n
Args:
a (int): left point of interval
n (int): number of subintervals
data (double): array that will hold values evaluated in function
Returns:
sum (double): pointer to result
*/
double x=0.0;
int idx;
idx = blockIdx.x * blockDim.x + threadIdx.x;
if(idx<=n-1){
x=a+(idx+1/2.0)*h_hat;
data[idx]=std::exp(-std::pow(x,2));
}
if(idx==0){
*sum = thrust::reduce(thrust::device, data , data + n, (double)0, thrust::plus<double>());
}
}
int main(int argc, char *argv[]){
double suma=0.0;
double *d_data;
double *d_suma;
double a=0.0, b=1.0;
double h_hat;
int n_threads_per_block=1024;
int n_bloques=0; //inicializamos en 0
double objetivo=0.7468241328124271;
int n=0; //inicializamos en 0
double time_spent;
clock_t begin,end;
cudaDeviceProp properties;
cudaGetDeviceProperties(&properties, 0);
n_bloques = 256 * properties.multiProcessorCount; //elegimos un múltiplo del número de SM's
//y además que sea menor a 10^4 (para tener 10^6 aprox subintervalos)
n = n_bloques*n_threads_per_block;//número de subintervalos
cudaMalloc((void **)&d_data,sizeof(double)*n);
cudaMalloc((void**)&d_suma,sizeof(double));
h_hat=(b-a)/n;
begin = clock();
Rcf<<<n_bloques,n_threads_per_block>>>(d_data, a,h_hat,n,d_suma);
cudaDeviceSynchronize();
end = clock();
time_spent = (double)(end - begin) / CLOCKS_PER_SEC;
cudaMemcpy(&suma, d_suma, sizeof(double), cudaMemcpyDeviceToHost);
suma=h_hat*suma;
cudaFree(d_data) ;
cudaFree(d_suma) ;
printf("Integral de %f a %f = %1.15e\n", a,b,suma);
printf("Error relativo de la solución: %1.15e\n", fabs(suma-objetivo)/fabs(objetivo));
printf("Tiempo de cálculo en la gpu %.5f\n", time_spent);
return 0;
}
In [68]:
%%bash
nvcc --compiler-options -Wall Rcf6.cu -o Rcf6.out
In [69]:
%%bash
./Rcf6.out
In [70]:
%%bash
nvprof --normalized-time-unit s ./Rcf6.out
$n=16777216 \approx 10^7$ subintervalos?
In [90]:
%%file Rcf7.cu
#include<stdio.h>
#include <thrust/reduce.h>
#include <thrust/execution_policy.h>
__global__ void Rcf(double *data, double a, double h_hat, int n, double *sum ) {
/*
Compute numerical approximation using rectangle or mid-point method in
an interval.
Nodes are generated via formula: x_i = a+(i+1/2)h_hat for i=0,1,...,n-1 and h_hat=(b-a)/n
Args:
a (int): left point of interval
n (int): number of subintervals
data (double): array that will hold values evaluated in function
Returns:
sum (double): pointer to result
*/
double x=0.0;
int idx;
idx = blockIdx.x * blockDim.x + threadIdx.x;
if(idx<=n-1){
x=a+(idx+1/2.0)*h_hat;
data[idx]=std::exp(-std::pow(x,2));
}
if(idx==0){
*sum = thrust::reduce(thrust::device, data , data + n, (double)0, thrust::plus<double>());
}
}
int main(int argc, char *argv[]){
double suma=0.0;
double *d_data;
double *d_suma;
double a=0.0, b=1.0;
double h_hat;
int n_threads_per_block=1024;
int n_bloques=0; //inicializamos en 0
double objetivo=0.7468241328124271;
int n=0; //inicializamos en 0
double time_spent;
clock_t begin,end;
cudaDeviceProp properties;
cudaGetDeviceProperties(&properties, 0);
n_bloques = 4096 * properties.multiProcessorCount; //elegimos un múltiplo del número de SM's
//y además que sea menor a 10^8 (para tener 10^7 aprox subintervalos)
n = n_bloques*n_threads_per_block;//número de subintervalos
cudaMalloc((void **)&d_data,sizeof(double)*n);
cudaMalloc((void**)&d_suma,sizeof(double));
h_hat=(b-a)/n;
begin = clock();
Rcf<<<n_bloques,n_threads_per_block>>>(d_data, a,h_hat,n,d_suma);
cudaDeviceSynchronize();
end = clock();
time_spent = (double)(end - begin) / CLOCKS_PER_SEC;
cudaMemcpy(&suma, d_suma, sizeof(double), cudaMemcpyDeviceToHost);
suma=h_hat*suma;
cudaFree(d_data) ;
cudaFree(d_suma) ;
printf("Integral de %f a %f = %1.15e\n", a,b,suma);
printf("Error relativo de la solución: %1.15e\n", fabs(suma-objetivo)/fabs(objetivo));
printf("Tiempo de cálculo en la gpu %.5f\n", time_spent);
return 0;
}
In [80]:
%%bash
nvcc --compiler-options -Wall Rcf7.cu -o Rcf7.out
In [81]:
%%bash
./Rcf7.out
In [75]:
%%bash
nvprof --normalized-time-unit s ./Rcf7.out
$n=117440512 \approx 10^8$ subintervalos?
Este caso se ejecutó en una GPU con las siguientes características:
----------------------
----device 0 ----
Device Name: Tesla K20Xm
Compute capability: 3.5
Clock rate: 732000
Unified memory: 1
---Memory Information for device 0 (results on bytes)---
Total global mem: 5977800704
Total constant Mem: 65536
Shared memory per thread block: 49152
Shared memory per SM: 49152
---MP Information for device 0 ---
SM count: 14
Threads in warp: 32
Max threads per SM: 2048
Max warps per SM: 64
Max threads per block: 1024
Max thread dimensions: (1024, 1024, 64)
Max grid dimensions: (2147483647, 65535, 65535)
----------------------
----device 1 ----
Device Name: Tesla K20Xm
Compute capability: 3.5
Clock rate: 732000
Unified memory: 1
---Memory Information for device 1 (results on bytes)---
Total global mem: 5977800704
Total constant Mem: 65536
Shared memory per thread block: 49152
Shared memory per SM: 49152
---MP Information for device 1 ---
SM count: 14
Threads in warp: 32
Max threads per SM: 2048
Max warps per SM: 64
Max threads per block: 1024
Max thread dimensions: (1024, 1024, 64)
Max grid dimensions: (2147483647, 65535, 65535)
In [1]:
%%file Rcf8.cu
#include<stdio.h>
#include <thrust/reduce.h>
#include <thrust/execution_policy.h>
__global__ void Rcf(double *data, double a, double h_hat, long int n, double *sum ) {
/*
Compute numerical approximation using rectangle or mid-point method in
an interval.
Nodes are generated via formula: x_i = a+(i+1/2)h_hat for i=0,1,...,n-1 and h_hat=(b-a)/n
Args:
a (int): left point of interval
n (int): number of subintervals
data (double): array that will hold values evaluated in function
Returns:
sum (double): pointer to result
*/
double x=0.0;
int idx;
idx = blockIdx.x * blockDim.x + threadIdx.x;
if(idx<=n-1){
x=a+(idx+1/2.0)*h_hat;
data[idx]=std::exp(-std::pow(x,2));
}
//esta alternativa no funciona... si el número de bloques rebasa el límite:
//int num_threads=gridDim.x * blockDim.x;
//int stride = num_threads;
//for(i=idx; i<=n-1; i+=stride){
// if(idx<=n-1){
// x=a+(idx+1/2.0)*h_hat;
// data[idx]=std::exp(-std::pow(x,2));
// }
//}
if(idx==0){
*sum = thrust::reduce(thrust::device, data , data + n, (double)0, thrust::plus<double>());
}
}
int main(int argc, char *argv[]){
double suma=0.0;
double *d_data;
double *d_suma;
double a=0.0, b=1.0;
double h_hat;
int n_threads_per_block=1024;
long int n_bloques=0; //inicializamos en 0
double objetivo=0.7468241328124271;
long int n=0; //inicializamos en 0
double time_spent;
clock_t begin,end;
cudaDeviceProp properties;
cudaGetDeviceProperties(&properties, 0);
n_bloques = 8192 * properties.multiProcessorCount; //elegimos un múltiplo del número de SM's
//y además que sea menor a 10^9 (para tener 10^8 aprox subintervalos)
n = n_bloques*n_threads_per_block;//número de subintervalos
dim3 dimGrid(n_bloques,1,1);
dim3 dimBlock(n_threads_per_block,1,1);
cudaMalloc((void **)&d_data,sizeof(double)*n);
cudaMalloc((void**)&d_suma,sizeof(double));
h_hat=(b-a)/n;
begin = clock();
Rcf<<<dimGrid,dimBlock>>>(d_data, a,h_hat,n,d_suma);
cudaDeviceSynchronize();
end = clock();
time_spent = (double)(end - begin) / CLOCKS_PER_SEC;
cudaMemcpy(&suma, d_suma, sizeof(double), cudaMemcpyDeviceToHost);
suma=h_hat*suma;
cudaFree(d_data) ;
cudaFree(d_suma) ;
printf("Integral de %f a %f = %1.15e\n", a,b,suma);
printf("Error relativo de la solución: %1.15e\n", fabs(suma-objetivo)/fabs(objetivo));
printf("Tiempo de cálculo en la gpu %.5f\n", time_spent);
return 0;
}
In [5]:
%%bash
nvcc --compiler-options -Wall Rcf8.cu -o Rcf8.out
In [6]:
%%bash
./Rcf8.out
In [7]:
%%bash
nvprof --normalized-time-unit s ./Rcf8.out
Obs: en la programación con CUDA-C es importante checar posibles errores de alojamiento de memoria. Una forma es con los tipos cudaError_t y cudaSuccess
. Ver liga:
In [1]:
%%file Rcf9.cu
#include<stdio.h>
#include <thrust/reduce.h>
#include <thrust/execution_policy.h>
__global__ void Rcf(double *data, double a, double h_hat, long int n, double *sum ) {
/*
Compute numerical approximation using rectangle or mid-point method in
an interval.
Nodes are generated via formula: x_i = a+(i+1/2)h_hat for i=0,1,...,n-1 and h_hat=(b-a)/n
Args:
a (int): left point of interval
n (int): number of subintervals
data (double): array that will hold values evaluated in function
Returns:
sum (double): pointer to result
*/
double x=0.0;
int idx;
int num_threads=gridDim.x * blockDim.x;
int stride = num_threads;
int i;
idx = blockIdx.x * blockDim.x + threadIdx.x;
for(i=idx; i<=n-1; i+=stride){
if(idx<=n-1){
x=a+(idx+1/2.0)*h_hat;
data[idx]=std::exp(-std::pow(x,2));
}
}
if(idx==0){
*sum = thrust::reduce(thrust::device, data , data + n, (double)0, thrust::plus<double>());
}
}
cudaError_t fun_checa_error(cudaError_t result) {
if (result != cudaSuccess) {
fprintf(stderr, "Error: %s\n", cudaGetErrorString(result));
}
return result;
}
int main(int argc, char *argv[]){
double suma=0.0;
double *d_data;
double *d_suma;
double a=0.0, b=1.0;
double h_hat;
int n_threads_per_block=1024;
long int n_bloques=0; //inicializamos en 0
double objetivo=0.7468241328124271;
long int n=0; //inicializamos en 0
double time_spent;
clock_t begin,end;
cudaDeviceProp properties;
cudaGetDeviceProperties(&properties, 0);
n_bloques = 32768 * properties.multiProcessorCount; //elegimos un múltiplo del número de SM's
//y además que sea menor a 10^9 (para tener 10^8 aprox subintervalos)
n = n_bloques*n_threads_per_block;//número de subintervalos
dim3 dimGrid(n_bloques,1,1);
dim3 dimBlock(n_threads_per_block,1,1);
fun_checa_error(cudaMalloc((void **)&d_data,sizeof(double)*n));
fun_checa_error(cudaMalloc((void**)&d_suma,sizeof(double)));
h_hat=(b-a)/n;
begin = clock();
Rcf<<<dimGrid,dimBlock>>>(d_data, a,h_hat,n,d_suma);
cudaDeviceSynchronize();
end = clock();
time_spent = (double)(end - begin) / CLOCKS_PER_SEC;
fun_checa_error(cudaMemcpy(&suma, d_suma, sizeof(double), cudaMemcpyDeviceToHost));
suma=h_hat*suma;
cudaFree(d_data);
cudaFree(d_suma);
printf("Integral de %f a %f = %1.15e\n", a,b,suma);
printf("Error relativo de la solución: %1.15e\n", fabs(suma-objetivo)/fabs(objetivo));
printf("Tiempo de cálculo en la gpu %.5f\n", time_spent);
return 0;
}
In [2]:
%%bash
nvcc --compiler-options -Wall Rcf9.cu -o Rcf9.out
In [3]:
%%bash
./Rcf9.out
Preguntas de comprehensión:
1)¿Qué es y en qué consiste CUDA C?
2)¿Qué es un kernel?
3)¿Qué pieza de CUDA se encarga de asignar los bloques de cuda-threads a las SM’s?
4)¿Qué características (recursos compartidos, dimensiones, forma de agendar la ejecución en threads) tienen los bloques que se asignan a una SM al lanzarse y ejecutarse un kernel?
5)¿Qué es un warp?
6)Menciona los tipos de memorias que existen en las GPU’s.
7) Supón que tienes una tarjeta GT200 cuyas características son:
* Máximo número de threads que soporta una SM en un mismo instante en el tiempo: 1024
* Máximo número de threads en un bloque: 512
* Máximo número de bloques por SM: 8
* Número de SM’s que tiene esta GPU: 30
Responde:
a) ¿Cuál es la máxima cantidad de threads que puede soportar esta GPU en un mismo instante en el tiempo?
b) ¿Cuál es la máxima cantidad de warps por SM que puede soportar esta GPU en un mismo instante en el tiempo?
c) ¿Cuáles configuraciones de bloques y threads siguientes aprovechan la máxima cantidad de warps en una SM de esta GPU para un mismo instante en el tiempo?
1.Una configuración del tipo: bloques de 64 threads y 16 bloques.
2.Una configuración del tipo: bloques de 1024 threads y 1 bloque.
3.Una configuración del tipo: bloques de 256 threads y 4 bloques.
4.Una configuración del tipo: bloques de 512 threads y 8 bloques.
*Debes considerar las restricciones/características de la GPU dadas para responder pues algunas configuraciones infringen las mismas. No estamos considerando registers o shared memory.
Referencias
N. Matloff, Parallel Computing for Data Science. With Examples in R, C++ and CUDA, 2014.
Para más sobre Unified Memory revisar:
Es importante el manejo de errores por ejemplo en el alojamiento de memoria en la GPU. En este caso es útil revisar:
En stackoverflow encontramos a personas desarrolladoras de CUDA que resuelven preguntas muy útiles para continuar con el aprendizaje de CUDA C. Por ejemplo:
Otros sistemas de software para el Heterogeneous computing son:
OpenCl. Ver NVIDIA OpenCL SDK Code Samples para ejemplos con NVIDIA GPU's.
Rth. Ver también rdrr.io matloff/Rth.