Multiplicación de matrices optimizada en C

10 minutos de lectura

avatar de usuario
Pedro

Estoy tratando de comparar diferentes métodos para la multiplicación de matrices. El primero es el método normal:

do
{
    for (j = 0; j < i; j++)
    {
        for (k = 0; k < i; k++)
        {
            suma = 0;
            for (l = 0; l < i; l++)
                suma += MatrixA[j][l]*MatrixB[l][k];
                MatrixR[j][k] = suma;
            }
        }
    }
    c++;
} while (c<iteraciones);

La segunda consiste en transponer primero la matriz B y luego hacer la multiplicación por filas:

int f, co;
for (f = 0; f < i; f++) {
    for ( co = 0; co < i; co++) {
        MatrixB[f][co] = MatrixB[co][f];
    }
}

c = 0;
do
{
    for (j = 0; j < i; j++)
    {
        for (k = 0; k < i; k++)
        {
            suma = 0;
            for (l = 0; l < i; l++)
                suma += MatrixA[j][l]*MatrixB[k][l];
                MatrixR[j][k] = suma;
            }
        }
     }
     c++;
} while (c<iteraciones);

Se supone que el segundo método es mucho más rápido, porque estamos accediendo a ranuras de memoria contiguas, pero no obtengo una mejora significativa en el rendimiento. ¿Estoy haciendo algo mal?

Puedo publicar el código completo, pero creo que no es necesario.

  • A menos que esté implementando sus propias rutinas de multiplicación de matrices como un ejercicio de aprendizaje, debería considerar seriamente el uso de una biblioteca existente, revisada y optimizada como BLAS o LAPACK.

    – maéricos

    21 de febrero de 2013 a las 16:33


  • El primer fragmento tiene 3 { y 4 }. Mi impresión es que lo más íntimo } no se quiere en absoluto, y la asignación MatrixR[j][k] = suma; no es parte de lo más íntimo for bucle, a pesar de la sangría (por lo que está mal sangrado, debe estar en el mismo nivel que suma = 0;).

    –Jonathan Leffler

    5 de agosto de 2018 a las 5:47


  • Esta respuesta podría ser útil: stackoverflow.com/a/54546544/3234205

    –Jianyu Huang

    6 de febrero de 2019 a las 4:34


Lo que todo programador debe saber sobre la memoria (enlace pdf) por Ulrich Drepper tiene muchas buenas ideas sobre la eficiencia de la memoria, pero en particular, usa la multiplicación de matrices como un ejemplo de cómo conocer la memoria y usar ese conocimiento puede acelerar este proceso. Consulte el apéndice A.1 de su artículo y lea la sección 6.2.1. La Tabla 6.2 del documento muestra que podría obtener un tiempo de ejecución del 10 % a partir del tiempo de una implementación ingenua para una matriz de 1000×1000.

De acuerdo, su código final es bastante peludo y usa muchas cosas específicas del sistema y ajustes en tiempo de compilación, pero aun así, si De Verdad necesita velocidad, leer ese documento y leer su implementación definitivamente vale la pena.

  • Si realmente necesita velocidad, simplemente use una implementación BLAS correcta y no desarrolle sus rutinas de álgebra lineal usted mismo.

    – Alejandro C.

    3 de junio de 2012 a las 22:19


  • Ulrich Drepper es mi héroe.

    – Ciro Santilli Путлер Капут 六四事

    13 de marzo de 2017 a las 9:46

avatar de usuario
peter joot

Hacer esto bien no es trivial. Se recomienda encarecidamente utilizar una biblioteca BLAS existente.

Si realmente te inclinas a rodar tu propia multiplicación de matrices, mosaico de bucle es una optimización que es de particular importancia para matrices grandes. El mosaico debe ajustarse al tamaño de la memoria caché para garantizar que la memoria caché no esté siendo golpeada continuamente, lo que ocurrirá con una implementación ingenua. Una vez medí una diferencia de rendimiento de 12x al combinar una matriz multiplicada con tamaños de matriz seleccionados para consumir múltiplos de mi caché (alrededor del ’97, por lo que el caché probablemente era pequeño).

Los algoritmos de mosaico en bucle suponen que se utiliza una matriz lineal contigua de elementos, en lugar de filas o columnas de punteros. Con una opción de almacenamiento de este tipo, su esquema de indexación determina qué dimensión cambia más rápido y usted es libre de decidir si el acceso a filas o columnas tendrá el mejor rendimiento de caché.

Hay una lote de literatura sobre el tema. Las siguientes referencias, especialmente los libros de Banerjee, pueden ser útiles:

[Ban93] Banerjee, Utpal, Loop Transformations for Restructuring Compilers: the Foundations, Kluwer Academic Publishers, Norwell, MA, 1993.

[Ban94] Banerjee, Utpal, Paralelización de bucles, Kluwer Academic Publishers, Norwell, MA, 1994.

[BGS93] Bacon, David F., Susan L. Graham y Oliver Sharp, Compiler Transformations for High-Performance Computing, Computer Science Division, University of California, Berkeley, Calif., Technical Report No UCB/CSD-93-781.

[LRW91] Lam, Monica S., Edward E. Rothberg y Michael E Wolf. The Cache Performance and Optimizations of Blocked Algorithms, en la 4ª Conferencia Internacional sobre Soporte Arquitectónico para Lenguajes de Programación, celebrada en Santa Clara, California, abril de 1991, 63-74.

[LW91] Lam, Monica S. y Michael E Wolf. Una teoría de transformación de bucles y un algoritmo para maximizar el paralelismo, en IEEE Transactions on Parallel and Distributed Systems, 1991, 2(4):452-471.

[PW86] Padua, David A. y Michael J. Wolfe, Optimizaciones avanzadas del compilador para supercomputadoras, en Comunicaciones de la ACM, 29(12):1184-1201, 1986.

[Wolfe89] Wolfe, Michael J. Optimización de supercompiladores para supercomputadoras, The MIT Press, Cambridge, MA, 1989.

[Wolfe96] Wolfe, Michael J., Compiladores de alto rendimiento para computación en paralelo, Addison-Wesley, CA, 1996.

ATENCIÓN: Tienes un BUG en tu segunda implementación

for (f = 0; f < i; f++) {
    for (co = 0; co < i; co++) {
        MatrixB[f][co] = MatrixB[co][f];
    }
}

Cuando haces f=0, c=1

        MatrixB[0][1] = MatrixB[1][0];

tu sobrescribes MatrixB[0][1] y perder ese valor! Cuando el ciclo llega a f=1, c=0

        MatrixB[1][0] = MatrixB[0][1];

el valor copiado es el mismo que ya estaba ahí.

  • ¡Oh! Ya veo, ¿debería usar una matriz tmp? Eso aumentará el tiempo aún más que el primer método.

    – Pedro

    15 de diciembre de 2009 a las 14:41

  • Puede transponer la matriz con una variable temporal: for (f=0; f<i-1; f++) for (co=f+1; co<i; co++) swap(MatrixB[f]+co, MatrixB[co]+f);

    – pmg

    16 de diciembre de 2009 a las 8:58

avatar de usuario
pablo rodriguez

Si la matriz no es lo suficientemente grande o no repites las operaciones un número elevado de veces no verás diferencias apreciables.

Si la matriz es, digamos, de 1.000×1.000 empezarás a ver mejoras, pero yo diría que si está por debajo de 100×100 no debes preocuparte.

Además, cualquier ‘mejora’ puede ser del orden de milisegundos, a menos que esté trabajando con matrices extremadamente grandes o repitiendo la operación miles de veces.

Finalmente, si cambias la computadora que estás usando por una más rápida, ¡las diferencias serán aún más pequeñas!

No debes escribir la multiplicación de matrices. Debe depender de bibliotecas externas. En particular, debe utilizar el GEMM rutina de la BLAS Biblioteca. GEMM a menudo proporciona las siguientes optimizaciones

Bloqueo

La multiplicación eficiente de matrices se basa en bloquear su matriz y realizar varias multiplicaciones bloqueadas más pequeñas. Idealmente, el tamaño de cada bloque se elige para que encaje bien en el caché. mejorando mucho el rendimiento.

Afinación

El tamaño de bloque ideal depende de la jerarquía de memoria subyacente (¿qué tan grande es el caché?). Como resultado, las bibliotecas deben ajustarse y compilarse para cada máquina específica. Esto lo hace, entre otros, el ATLAS implementación de BLAS.

Optimización del nivel de ensamblaje

La multiplicación de matrices es tan común que los desarrolladores la optimizarán manualmente. En particular, esto se hace en GotoBLAS.

Computación heterogénea (GPU)

Matrix Multiply es muy intensivo en FLOP/computación, lo que lo convierte en un candidato ideal para ejecutarse en GPU. cuBLAS y MAGMA son buenos candidatos para esto.

En resumen, el álgebra lineal densa es un tema bien estudiado. La gente dedica su vida a la mejora de estos algoritmos. Deberías usar su trabajo; Los hará felices.

  • Es decir solamente verdadero si las matrices son lo suficientemente grandes como para justificar el uso de algoritmos más complejos implementados en BLAS (AFAIK, no hay una ruta rápida para matrices pequeñas). Pero +1 por citar las rutinas GEMM (“multiplicación general de matrices”) y explicar por qué ayuda.

    – mctylr

    23 mayo 2014 a las 14:21

avatar de usuario
Andreas Brick

Las grandes mejoras que obtenga dependerán de:

  1. El tamaño del caché
  2. El tamaño de una línea de caché
  3. El grado de asociatividad de la memoria caché.

Para tamaños de matriz pequeños y procesadores modernos, es muy probable que los datos de ambos MatrixA y MatrixB se mantendrá casi por completo en el caché después de tocarlo por primera vez.

  • Es decir solamente verdadero si las matrices son lo suficientemente grandes como para justificar el uso de algoritmos más complejos implementados en BLAS (AFAIK, no hay una ruta rápida para matrices pequeñas). Pero +1 por citar las rutinas GEMM (“multiplicación general de matrices”) y explicar por qué ayuda.

    – mctylr

    23 mayo 2014 a las 14:21

avatar de usuario
San Jacinto

Solo algo para que pruebe (pero esto solo marcaría la diferencia para matrices grandes): separe su lógica de suma de la lógica de multiplicación en el ciclo interno de la siguiente manera:

for (k = 0; k < i; k++)
{
    int sums[i];//I know this size declaration is illegal in C. consider 
            //this pseudo-code.
    for (l = 0; l < i; l++)
        sums[l] = MatrixA[j][l]*MatrixB[k][l];

    int suma = 0;
    for(int s = 0; s < i; s++)
       suma += sums[s];
}

Esto se debe a que terminas estancando tu tubería cuando escribes a suma. De acuerdo, gran parte de esto se soluciona en el cambio de nombre de registro y cosas por el estilo, pero con mi comprensión limitada del hardware, si quisiera exprimir cada gramo de rendimiento del código, lo haría porque ahora no tiene que hacerlo. detenga la canalización para esperar una escritura en suma. Dado que la multiplicación es más costosa que la suma, desea dejar que la máquina la paralelice tanto como sea posible, por lo que guardar sus puestos para la suma significa que pasará menos tiempo esperando en el ciclo de suma que en el ciclo de multiplicación.

Esta es solo mi lógica. Otros con más conocimiento en el área pueden estar en desacuerdo.

  • Comentario 1: si le preocupa la canalización de instrucciones, considere desenrollar los bucles. Reducir el número de bifurcaciones acelerará mucho más la ejecución que separar la multiplicación de la suma (algunos procesadores tienen instrucciones simples que multiplican y suman).

    – Thomas Matthews

    15 de diciembre de 2009 a las 19:33

  • Comentario 2: Mover el sums matriz y el suma variables fuera del alcance de la for lazo. Considere también mover las variables de índice l y s también. Cuando están dentro del exterior for bucle que crearán y destruirán para cada iteración del exterior for lazo. La eliminación de instrucciones innecesarias también acelerará la ejecución.

    – Thomas Matthews

    15 de diciembre de 2009 a las 19:36

  • De hecho, estaba asumiendo que el compilador desenrollaría los bucles porque eso es bastante común, pero creo que es poco probable que almacene en caché las multiplicaciones por separado, ya que es difícil determinar la semántica de una reducción. Tal vez estoy desgastado y cada compilador lo hace. para su segundo comentario, usted es incorrecto. la semántica en el código muestra que la memoria se está saliendo de la pila, lo que significa que se reserva al comienzo de la función, una vez.

    – San Jacinto

    15 de diciembre de 2009 a las 19:51

  • Además, con respecto a tu primer comentario, creo que te estás perdiendo el punto. No es la atomicidad del mult lo que nos preocupa, es el tiempo que lleva ejecutar un mult. Los procesadores modernos aún pueden tardar decenas de veces más en hacer un mult que en un add. Por lo tanto, si detiene la canalización hasta que tenga un valor para poner en suma, es posible que esté esperando 40 ciclos o más más de lo necesario en cada iteración del ciclo, en lugar de detenerse solo una o dos veces (si desenrolla el ciclo , o el compilador lo hace).

    – San Jacinto

    16 de diciembre de 2009 a las 0:31

¿Ha sido útil esta solución?

Esta web utiliza cookies propias y de terceros para su correcto funcionamiento y para fines analíticos y para mostrarte publicidad relacionada con sus preferencias en base a un perfil elaborado a partir de tus hábitos de navegación. Al hacer clic en el botón Aceptar, acepta el uso de estas tecnologías y el procesamiento de tus datos para estos propósitos. Configurar y más información
Privacidad