CUDA reducción de muchas matrices pequeñas, de un tamaño desigual

Question

Me pregunto si alguien podría sugerir el mejor enfoque para el cálculo de la desviación media / estándar de un gran número de relativamente pequeño, pero matrices de diferentes tamaños en CUDA?

El ejemplo reducción paralela en el SDK funciona en una única matriz muy grande y parece que el tamaño es, convenientemente, un múltiplo del número de hilos por bloque, pero mi caso es bastante diferente:

Conceptualmente, I sin embargo tienen un gran número de objetos que contienen cada uno dos componentes, upper y lower y cada uno de estos componentes tiene un x y coordinar un y. es decir.

upper.x, lower.x, upper.y, lower.y

Cada una de estas matrices es aproximadamente 800 en longitud pero varía entre los objetos (no dentro de un objeto) por ejemplo.

Object1.lower.x = 1.1, 2.2, 3.3
Object1.lower.y = 4.4, 5.5, 6.6
Object1.upper.x = 7.7, 8.8, 9.9
Object1.upper.y = 1.1, 2.2, 3.3

Object2.lower.x = 1.0,  2.0,  3.0,  4.0, 5.0 
Object2.lower.y = 6.0,  7.0,  8.0,  9.0, 10.0
Object2.upper.x = 11.0, 12.0, 13.0, 14.0, 15.0 
Object2.upper.y = 16.0, 17.0, 18.0, 19.0, 20.0

Tenga en cuenta lo anterior es sólo mi manera de representar la matriz y mi datos no se almacenan en estructuras C ni nada de eso: los datos pueden organizarse en la forma que necesito. El punto es que, para cada matriz, la media, la desviación estándar y, finalmente, un histograma necesita ser computado y dentro de uno objeto particular, relaciones y diferencias entre las matrices necesitan ser calculados.

¿Cómo debo ir sobre el envío de estos datos al dispositivo de GPU y la organización de mi jerarquía de hilo bloque? Una idea que tenía era para rellenar cero todos mis arreglos para que sean de la misma longitud, y tener un grupo de bloques que trabajan en cada objeto, pero parece que hay todo tipo de problemas con ese método si funcionaría en absoluto.

Gracias de antemano

Answer 1

Si he entendido bien, desea reducir Object1.lower.x a un resultado, Object1.lower.y a otro resultado y así sucesivamente. Para cualquier objeto dado hay cuatro matrices para ser reducido, todas de igual longitud (para el objeto).

Hay muchos enfoques posibles para esto, uno de los factores que influyen sería el número total de objetos en el sistema. Vamos a suponer que el número es grande.

Para obtener un mejor rendimiento que desea un patrón de acceso a la memoria óptima y se quiere evitar la divergencia. Dado que el número de matrices congruentes es de cuatro, si se toma el enfoque ingenuo de hacer una matriz por hilo, a continuación, no sólo va a sufrir de tener acceso deficiente memoria, sino también la h / w tendrá que comprobar en cada iteración que se enrosca en la urdimbre necesidad de ejecutar el bucle -. aquellos que no será deshabilitado que puede ser ineficiente (especialmente si una matriz es mucho más largo que los otros, por ejemplo)

for (int i = 0 ; i < myarraylength ; i++)
    sum += myarray[i];

En cambio si se obtiene cada urdimbre para resumir una matriz entonces no sólo va a ser más eficiente, sino también su patrón de acceso a la memoria será mucho mejor, ya que los hilos adyacentes leerán los elementos adyacentes [1].

for (int i = tidwithinwarp ; i < warparraylength ; i += warpsize)
{
    mysum += warparray[i];
}
mysum = warpreduce(mysum);

También debe tomar la alineación de las matrices en consideración, preferentemente alinear en un límite de 64 bytes, aunque si está desarrollando para la capacidad de cómputo 1.2 o superior, entonces esto no es tan importante como en las GPU de mayor edad.

En este ejemplo se lanzará cuatro urdimbres por bloque, es decir, 128 hilos, y otros tantos bloques como usted tiene objetos.

[1] Usted dice que usted puede elegir cualquier disposición de memoria se quiere, a menudo puede ser útil para intercalar las matrices de manera que array [0] [0] está al lado de array [1] [0], ya que esto significa que los hilos adyacentes pueden operar en matrices adyacentes y obtener accesos coalescentes. Sin embargo dado que la longitud de las matrices no es constante esto es probablemente compleja, que requiere el relleno de las matrices más cortos.

Answer 2

Como seguimiento de la respuesta de Tom, me gustaría mencionar que reducción de urdimbre se puede implementar fácilmente por CUB .

Este es un ejemplo práctico:

#include <cub/cub.cuh>
#include <cuda.h>

#include "Utilities.cuh"

#include <iostream>

#define WARPSIZE    32
#define BLOCKSIZE   256

const int N = 1024;

/*************************/
/* WARP REDUCTION KERNEL */
/*************************/
__global__ void sum(const float * __restrict__ indata, float * __restrict__ outdata) {

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

    unsigned int warp_id = threadIdx.x / WARPSIZE;

    // --- Specialize WarpReduce for type float. 
    typedef cub::WarpReduce<float, WARPSIZE> WarpReduce;

    // --- Allocate WarpReduce shared memory for (N / WARPSIZE) warps
    __shared__ typename WarpReduce::TempStorage temp_storage[BLOCKSIZE / WARPSIZE];

    float result;
    if(tid < N) result = WarpReduce(temp_storage[warp_id]).Sum(indata[tid]);

    if(tid % WARPSIZE == 0) outdata[tid / WARPSIZE] = result;
}

/********/
/* MAIN */
/********/
int main() {

    // --- Allocate host side space for 
    float *h_data       = (float *)malloc(N * sizeof(float));
    float *h_result     = (float *)malloc((N / WARPSIZE) * sizeof(float));

    float *d_data;      gpuErrchk(cudaMalloc(&d_data, N * sizeof(float)));
    float *d_result;    gpuErrchk(cudaMalloc(&d_result, (N / WARPSIZE) * sizeof(float)));

    for (int i = 0; i < N; i++) h_data[i] = (float)i;

    gpuErrchk(cudaMemcpy(d_data, h_data, N * sizeof(float), cudaMemcpyHostToDevice));

    sum<<<iDivUp(N, BLOCKSIZE), BLOCKSIZE>>>(d_data, d_result);
    gpuErrchk(cudaPeekAtLastError());
    gpuErrchk(cudaDeviceSynchronize());

    gpuErrchk(cudaMemcpy(h_result, d_result, (N / WARPSIZE) * sizeof(float), cudaMemcpyDeviceToHost));

    std::cout << "output: ";
    for(int i = 0; i < (N / WARPSIZE); i++) std::cout << h_result[i] << " ";
    std::cout << std::endl;

    gpuErrchk(cudaFree(d_data));
    gpuErrchk(cudaFree(d_result));

    return 0;
}

En este ejemplo, se crea una matriz de longitud N y el resultado es la suma de elementos consecutivos 32. Así

result[0] = data[0] + ... + data[31];
result[1] = data[32] + ... + data[63];
....