tutorial - usar gpu en python

Matrices de orden C no cuadradas en cuBLAS(numba)

Question

Estoy tratando de usar las funciones cuBLAS en el paquete Numba de Anaconda y tengo un problema. Necesito que las matrices de entrada estén en orden C La salida puede estar en orden de Fortran.

Puedo ejecutar el script de ejemplo proporcionado con el paquete, here . El script tiene dos funciones, gemm_v1 y gemm_v2 . En gemm_v1 , el usuario tiene que crear las matrices de entrada en orden de Fortran. En gemm_v2 , se pueden pasar a la implementación cuda de GEMM y se pueden transponer en el dispositivo. Puedo conseguir que estos ejemplos trabajen con matrices cuadradas. Sin embargo, no puedo averiguar cómo hacer que gemm_v2 funcione con matrices de entrada no cuadradas. ¿Hay alguna manera de trabajar con matrices de entrada de orden C que no sean cuadradas?

Nota:
Idealmente, tanto las matrices de entrada como las de salida permanecerían en el dispositivo después de la llamada a GEMM para ser utilizadas en otros cálculos (esto es parte de un método iterativo).

Answer 1

El problema con este ejemplo es que solo funciona para matrices cuadradas. Si las matrices no son cuadradas, no puede calcular A^t*B^t debido a la falta de coincidencia de dimensión (suponiendo que las dimensiones fueran correctas para A*B ).

No tengo una instalación de cuBLAS en funcionamiento, así que es una especie de tiro en la oscuridad, pero me sorprendería mucho si cuBLAS funcionara de manera diferente a la BLAS habitual. BLAS espera que las matrices estén en orden mayor de columna (también conocido como orden de Fortran), pero también se pueden usar para matrices en orden mayor de fila (también conocido como orden C).

En mi opinión, lo que podría estar completamente equivocado, gemm_v2 no es la forma habitual / mejor de manejar la multiplicación de dos matrices de orden C, por ejemplo, porque si una multiplica dos matrices de orden C, también tendríamos una matriz de orden C como respuesta .

El truco para calcular el producto de dos matrices de orden C con la ayuda de gemm funcionaría de la siguiente manera:

Incluso si es probable que lo conozca, primero me gustaría desarrollar el orden de fila mayor (diseño de memoria-c) y el orden de columna mayor (diseño de memoria de fortran) para poder desarrollar mi responder.

Entonces, si tenemos una matriz A 2x3 (es decir, 2 filas y 3 columnas), y la almacenamos en alguna memoria continua obtenemos:

row-major-order(A) = A11, A12, A13, A21, A22, A23
col-major-order(A) = A11, A21, A12, A22, A13, A33

Eso significa que si obtenemos una memoria continua, que representa una matriz en el orden mayor de la fila, y si la interpretamos como una matriz en el orden mayor de la columna, obtendremos una matriz bastante diferente.

Sin embargo, si echamos un vistazo a la matriz transpuesta A^t podemos ver fácilmente:

row-major-order(A) = col-major-order(A^t)
col-major-order(A) = row-major-order(A^t)

Eso significa que, si quisiéramos obtener la matriz C en orden mayor de filas como resultado, la rutina blas debería escribir la matriz C transpuesta en orden mayor de columnas (después de todo esto no podemos cambiar) en esta misma memoria. Sin embargo, C^t=(AB)^t=B^t*A^t y B^t y A^t son las matrices originales reinterpretadas en orden de columna mayor.

Ahora, sea A una nxk nxk y B una kxm kxm, la rutina de llamada de gemm debería ser la siguiente:

gemm(''N'', ''N'', m, n, k, 1.0, B, m, A, k, 0.0, C, m)

Tenga en cuenta:

1. No tenemos que transponer las matrices A y B , porque se maneja reinterpretando C-order como Fortran-order.
2. Tenemos que cambiar los lugares de las matrices A y B para obtener C^t en orden de Fortran como resultado.
3. La matriz C resultante está en orden C (reinterpretándola de orden de Fortran a orden C de la que nos deshacemos de ^t ).