iphone - resueltos - Multiplicación de matriz 4x4 rápida en C
multiplicacion de matrices online (5)
Estoy tratando de encontrar una implementación optimizada de C o Assembler de una función que multiplique dos matrices 4x4 entre sí. La plataforma es un iPhone o iPod basado en ARM6 o ARM7.
Actualmente, estoy usando un enfoque bastante estándar, solo un pequeño ciclo desenrollado.
#define O(y,x) (y + (x<<2)) static inline void Matrix4x4MultiplyBy4x4 (float *src1, float *src2, float *dest) { *(dest+O(0,0)) = (*(src1+O(0,0)) * *(src2+O(0,0))) + (*(src1+O(0,1)) * *(src2+O(1,0))) + (*(src1+O(0,2)) * *(src2+O(2,0))) + (*(src1+O(0,3)) * *(src2+O(3,0))); *(dest+O(0,1)) = (*(src1+O(0,0)) * *(src2+O(0,1))) + (*(src1+O(0,1)) * *(src2+O(1,1))) + (*(src1+O(0,2)) * *(src2+O(2,1))) + (*(src1+O(0,3)) * *(src2+O(3,1))); *(dest+O(0,2)) = (*(src1+O(0,0)) * *(src2+O(0,2))) + (*(src1+O(0,1)) * *(src2+O(1,2))) + (*(src1+O(0,2)) * *(src2+O(2,2))) + (*(src1+O(0,3)) * *(src2+O(3,2))); *(dest+O(0,3)) = (*(src1+O(0,0)) * *(src2+O(0,3))) + (*(src1+O(0,1)) * *(src2+O(1,3))) + (*(src1+O(0,2)) * *(src2+O(2,3))) + (*(src1+O(0,3)) * *(src2+O(3,3))); *(dest+O(1,0)) = (*(src1+O(1,0)) * *(src2+O(0,0))) + (*(src1+O(1,1)) * *(src2+O(1,0))) + (*(src1+O(1,2)) * *(src2+O(2,0))) + (*(src1+O(1,3)) * *(src2+O(3,0))); *(dest+O(1,1)) = (*(src1+O(1,0)) * *(src2+O(0,1))) + (*(src1+O(1,1)) * *(src2+O(1,1))) + (*(src1+O(1,2)) * *(src2+O(2,1))) + (*(src1+O(1,3)) * *(src2+O(3,1))); *(dest+O(1,2)) = (*(src1+O(1,0)) * *(src2+O(0,2))) + (*(src1+O(1,1)) * *(src2+O(1,2))) + (*(src1+O(1,2)) * *(src2+O(2,2))) + (*(src1+O(1,3)) * *(src2+O(3,2))); *(dest+O(1,3)) = (*(src1+O(1,0)) * *(src2+O(0,3))) + (*(src1+O(1,1)) * *(src2+O(1,3))) + (*(src1+O(1,2)) * *(src2+O(2,3))) + (*(src1+O(1,3)) * *(src2+O(3,3))); *(dest+O(2,0)) = (*(src1+O(2,0)) * *(src2+O(0,0))) + (*(src1+O(2,1)) * *(src2+O(1,0))) + (*(src1+O(2,2)) * *(src2+O(2,0))) + (*(src1+O(2,3)) * *(src2+O(3,0))); *(dest+O(2,1)) = (*(src1+O(2,0)) * *(src2+O(0,1))) + (*(src1+O(2,1)) * *(src2+O(1,1))) + (*(src1+O(2,2)) * *(src2+O(2,1))) + (*(src1+O(2,3)) * *(src2+O(3,1))); *(dest+O(2,2)) = (*(src1+O(2,0)) * *(src2+O(0,2))) + (*(src1+O(2,1)) * *(src2+O(1,2))) + (*(src1+O(2,2)) * *(src2+O(2,2))) + (*(src1+O(2,3)) * *(src2+O(3,2))); *(dest+O(2,3)) = (*(src1+O(2,0)) * *(src2+O(0,3))) + (*(src1+O(2,1)) * *(src2+O(1,3))) + (*(src1+O(2,2)) * *(src2+O(2,3))) + (*(src1+O(2,3)) * *(src2+O(3,3))); *(dest+O(3,0)) = (*(src1+O(3,0)) * *(src2+O(0,0))) + (*(src1+O(3,1)) * *(src2+O(1,0))) + (*(src1+O(3,2)) * *(src2+O(2,0))) + (*(src1+O(3,3)) * *(src2+O(3,0))); *(dest+O(3,1)) = (*(src1+O(3,0)) * *(src2+O(0,1))) + (*(src1+O(3,1)) * *(src2+O(1,1))) + (*(src1+O(3,2)) * *(src2+O(2,1))) + (*(src1+O(3,3)) * *(src2+O(3,1))); *(dest+O(3,2)) = (*(src1+O(3,0)) * *(src2+O(0,2))) + (*(src1+O(3,1)) * *(src2+O(1,2))) + (*(src1+O(3,2)) * *(src2+O(2,2))) + (*(src1+O(3,3)) * *(src2+O(3,2))); *(dest+O(3,3)) = (*(src1+O(3,0)) * *(src2+O(0,3))) + (*(src1+O(3,1)) * *(src2+O(1,3))) + (*(src1+O(3,2)) * *(src2+O(2,3))) + (*(src1+O(3,3)) * *(src2+O(3,3))); };
¿Me beneficiaría usar el algoritmo Strassen o Coppersmith-Winograd?
¿Estás seguro de que tu código desenrollado es más rápido que el enfoque basado en bucle explícito? ¡Tenga en cuenta que los compiladores suelen ser mejores que los humanos que realizan optimizaciones!
De hecho, apostaría a que hay más posibilidades de que un compilador emita automáticamente instrucciones SIMD desde un ciclo bien escrito que desde una serie de declaraciones "no relacionadas" ...
También podría especificar los tamaños de matrices en la declaración del argumento. Luego, podría usar la sintaxis del paréntesis normal para acceder a los elementos, y también podría ser una buena pista para que el compilador también realice sus optimizaciones.
¿Son estas matrices arbitrarias o tienen simetrías? Si es así, esas simetrías pueden explotarse a menudo para mejorar el rendimiento (por ejemplo, en matrices de rotación).
Además, estoy de acuerdo con Fortran anteriormente, y ejecutaría algunas pruebas de tiempo para verificar que su código desenrollado a mano sea más rápido de lo que puede crear un compilador de optimización. Por lo menos, puede simplificar su código.
Pablo
No, el algoritmo de Strassen o Coppersmith-Winograd no haría mucha diferencia aquí. Comienzan a pagar por matrices más grandes solamente.
Si su multiplicación de matriz es realmente un cuello de botella, podría reescribir el algoritmo usando las instrucciones de NEON SIMD. Eso solo ayudaría con ARMv7 ya que ARMv6 no tiene esta extensión.
Esperaría una aceleración factor 3 sobre el código C compilado para su caso.
EDITAR: Puede encontrar una buena implementación en ARM-NEON aquí: http://code.google.com/p/math-neon/
Para su código C, hay dos cosas que puede hacer para acelerar el código:
No alinee la función. La multiplicación de la matriz genera bastante código a medida que se desenrolla, y el ARM solo tiene un caché de instrucciones muy pequeño. Una inversión excesiva puede hacer que su código sea más lento porque la CPU estará ocupada cargando código en la memoria caché en lugar de ejecutarlo.
Use la palabra clave restrict para indicar al compilador que los punteros de origen y de destino no se superponen en la memoria. Actualmente, el compilador se ve obligado a volver a cargar cada valor fuente de la memoria cada vez que se escribe un resultado, ya que debe asumir que el origen y el destino pueden superponerse o incluso apuntar a la misma memoria.
Sólo curioseando. Me pregunto por qué las personas todavía ofuscan su código voluntariamente. C ya es difícil de leer, no es necesario agregarlo.
static inline void Matrix4x4MultiplyBy4x4 (float src1[4][4], float src2[4][4], float dest[4][4])
{
dest[0][0] = src1[0][0] * src2[0][0] + src1[0][1] * src2[1][0] + src1[0][2] * src2[2][0] + src1[0][3] * src2[3][0];
dest[0][1] = src1[0][0] * src2[0][1] + src1[0][1] * src2[1][1] + src1[0][2] * src2[2][1] + src1[0][3] * src2[3][1];
dest[0][2] = src1[0][0] * src2[0][2] + src1[0][1] * src2[1][2] + src1[0][2] * src2[2][2] + src1[0][3] * src2[3][2];
dest[0][3] = src1[0][0] * src2[0][3] + src1[0][1] * src2[1][3] + src1[0][2] * src2[2][3] + src1[0][3] * src2[3][3];
dest[1][0] = src1[1][0] * src2[0][0] + src1[1][1] * src2[1][0] + src1[1][2] * src2[2][0] + src1[1][3] * src2[3][0];
dest[1][1] = src1[1][0] * src2[0][1] + src1[1][1] * src2[1][1] + src1[1][2] * src2[2][1] + src1[1][3] * src2[3][1];
dest[1][2] = src1[1][0] * src2[0][2] + src1[1][1] * src2[1][2] + src1[1][2] * src2[2][2] + src1[1][3] * src2[3][2];
dest[1][3] = src1[1][0] * src2[0][3] + src1[1][1] * src2[1][3] + src1[1][2] * src2[2][3] + src1[1][3] * src2[3][3];
dest[2][0] = src1[2][0] * src2[0][0] + src1[2][1] * src2[1][0] + src1[2][2] * src2[2][0] + src1[2][3] * src2[3][0];
dest[2][1] = src1[2][0] * src2[0][1] + src1[2][1] * src2[1][1] + src1[2][2] * src2[2][1] + src1[2][3] * src2[3][1];
dest[2][2] = src1[2][0] * src2[0][2] + src1[2][1] * src2[1][2] + src1[2][2] * src2[2][2] + src1[2][3] * src2[3][2];
dest[2][3] = src1[2][0] * src2[0][3] + src1[2][1] * src2[1][3] + src1[2][2] * src2[2][3] + src1[2][3] * src2[3][3];
dest[3][0] = src1[3][0] * src2[0][0] + src1[3][1] * src2[1][0] + src1[3][2] * src2[2][0] + src1[3][3] * src2[3][0];
dest[3][1] = src1[3][0] * src2[0][1] + src1[3][1] * src2[1][1] + src1[3][2] * src2[2][1] + src1[3][3] * src2[3][1];
dest[3][2] = src1[3][0] * src2[0][2] + src1[3][1] * src2[1][2] + src1[3][2] * src2[2][2] + src1[3][3] * src2[3][2];
dest[3][3] = src1[3][0] * src2[0][3] + src1[3][1] * src2[1][3] + src1[3][2] * src2[2][3] + src1[3][3] * src2[3][3];
};
Su producto tradicional completamente desenrollado probablemente sea bastante rápido.
Su matriz es demasiado pequeña para superar la administración de una multiplicación de Strassen en su forma tradicional con índices explícitos y códigos de particionamiento; es probable que pierda cualquier efecto en la optimización de esa sobrecarga.
Pero si quieres rápido, usaría instrucciones SIMD si están disponibles. Me sorprendería si los chips ARM en estos días no los tienen. Si lo hacen, puedes administrar todos los productos en fila / columna en una sola instrucción; si el SIMD tiene 8 anchos, puede administrar 2 repeticiones de fila / columna en una sola instrucción. Configurar los operandos para hacer esa instrucción puede requerir un poco de baile; Las instrucciones SIMD recogerán fácilmente sus filas (valores adyacentes), pero no recogerán las columnas (no contiguas). Y puede tomar algún esfuerzo calcular la suma de los productos en una fila / columna.