performance - suma - ¿Por qué MATLAB es tan rápido en la multiplicación de matrices?
multiplicacion de matrices matlab codigo (14)
Estoy haciendo algunos puntos de referencia con CUDA, C ++, C # y Java, y el uso de MATLAB para la verificación y generación de matriz. Pero cuando me multiplico con MATLAB, 2048x2048 e incluso matrices más grandes se multiplican casi al instante.
1024x1024 2048x2048 4096x4096
--------- --------- ---------
CUDA C (ms) 43.11 391.05 3407.99
C++ (ms) 6137.10 64369.29 551390.93
C# (ms) 10509.00 300684.00 2527250.00
Java (ms) 9149.90 92562.28 838357.94
MATLAB (ms) 75.01 423.10 3133.90
Solo CUDA es competitivo, pero pensé que al menos C ++ será algo cercano y no 60 veces más lento.
Entonces mi pregunta es: ¿Cómo lo está haciendo MATLAB tan rápido?
Código C ++:
float temp = 0;
timer.start();
for(int j = 0; j < rozmer; j++)
{
for (int k = 0; k < rozmer; k++)
{
temp = 0;
for (int m = 0; m < rozmer; m++)
{
temp = temp + matice1[j][m] * matice2[m][k];
}
matice3[j][k] = temp;
}
}
timer.stop();
Editar: Tampoco sé qué pensar sobre los resultados C #. El algoritmo es lo mismo que C ++ y Java, pero hay un salto gigante 2048 desde 1024.
Edit2: resultados actualizados de MATLAB y 4096x4096
¿Verificaron que todas las implementaciones usaran optimizaciones multihilo para el algoritmo? ¿Y usaron el mismo algoritmo de multiplicación?
Realmente lo dudo.
Matlab no es intrínsecamente rápido, probablemente haya utilizado implementaciones lentas.
Al hacer la multiplicación de matrices, se usa un método de multiplicación ingenuo que toma tiempo de O(n^3) .
Existe un algoritmo de multiplicación de matriz que toma O(n^2.4) . Lo que significa que en n=2000 su algoritmo requiere ~ 100 veces más computación que el mejor algoritmo.
Debería consultar la página de wikipedia para ver la multiplicación de matrices para obtener más información sobre las formas eficientes de implementarla.
Aquí están mis resultados usando MATLAB R2011a + Parallel Computing Toolbox en una máquina con Tesla C2070:
>> A = rand(1024); gA = gpuArray(A);
% warm up by executing the operations a couple of times, and then:
>> tic, C = A * A; toc
Elapsed time is 0.075396 seconds.
>> tic, gC = gA * gA; toc
Elapsed time is 0.008621 seconds.
MATLAB utiliza bibliotecas altamente optimizadas para la multiplicación de matrices, razón por la cual la simple multiplicación de matrices de MATLAB es muy rápida. La versión de gpuArray usa MAGMA .
Actualice usando R2014a en una máquina con Tesla K20c, y las nuevas funciones timeit y gputimeit :
>> A = rand(1024); gA = gpuArray(A);
>> timeit(@()A*A)
ans =
0.0324
>> gputimeit(@()gA*gA)
ans =
0.0022
Debe tener cuidado al hacer comparaciones justas con C ++. ¿Puede publicar el código de C ++ que muestra los bucles interiores centrales que está utilizando para la multiplicación de matrices? Sobre todo, me preocupa tu diseño de memoria y si estás haciendo las cosas de forma inútil.
He escrito la multiplicación de matrices C ++ que es tan rápida como la de Matlab, pero me tomó un poco de cuidado. (EDITAR: antes de que Matlab utilizara GPU para esto).
Puede estar virtualmente garantizado que Matlab está desperdiciando muy pocos ciclos en estas funciones "integradas". Mi pregunta es, ¿dónde estás perdiendo ciclos? (Sin ofender)
Dependiendo de su versión de Matlab, creo que ya podría estar usando su GPU.
Otra cosa; Matlab realiza un seguimiento de muchas propiedades de su matriz; ya sea diagonal, hermetista, etc., y especializa sus algoritmos basados en eso. Tal vez su especialización se basa en la matriz cero que está pasando, o algo así? Tal vez está almacenando en caché las llamadas a funciones repetidas, lo que estropea los tiempos. ¿Tal vez optimiza productos de matriz no utilizados repetidos?
Para evitar que sucedan estas cosas, use una matriz de números aleatorios y asegúrese de forzar la ejecución imprimiendo el resultado en una pantalla o disco o algo así.
El fuerte contraste no solo se debe a la sorprendente optimización de Matlab (como ya se discutió en muchas otras respuestas), sino también a la forma en que formuló la matriz como un objeto.
Parece que hiciste matriz una lista de listas? Una lista de listas contiene punteros a listas que luego contienen sus elementos de matriz. Las ubicaciones de las listas contenidas se asignan arbitrariamente. A medida que recorre su primer índice (¿número de fila?), El tiempo de acceso a la memoria es muy significativo. En comparación, ¿por qué no intentas implementar matrix como una lista / vector único usando el siguiente método?
#include <vector>
struct matrix {
matrix(int x, int y) : n_row(x), n_col(y), M(x * y) {}
int n_row;
int n_col;
std::vector<double> M;
double &operator()(int i, int j);
};
Y
double &matrix::operator()(int i, int j) {
return M[n_col * i + j];
}
Se debe usar el mismo algoritmo de multiplicación para que el número de flop sea el mismo. (n ^ 3 para matrices cuadradas de tamaño n)
Te pido que tengas tiempo para que el resultado sea comparable a lo que tenías antes (en la misma máquina). Con la comparación, usted mostrará exactamente cuán significativo puede ser el tiempo de acceso a la memoria.
El uso de dobles y una matriz sólida en lugar de tres conduce mi código C # a casi los mismos resultados que C ++ / Java (con tu código: 1024 - fue un poco más rápido, 2048 - fue de aproximadamente 140 y 4096 - fue de aproximadamente 22 minutos)
1024x1024 2048x2048 4096x4096 --------- --------- --------- your C++ (ms) 6137.10 64369.29 551390.93 my C# (ms) 9730.00 90875.00 1062156.00
aquí está mi código:
const int rozmer = 1024;
double[][] matice1 = new double[rozmer * 3][];
Random rnd = new Random();
public Form1()
{
InitializeComponent();
System.Threading.Thread thr = new System.Threading.Thread(new System.Threading.ThreadStart(() =>
{
string res = "";
Stopwatch timer = new Stopwatch();
timer.Start();
double temp = 0;
int r2 = rozmer * 2;
for (int i = 0; i < rozmer*3; i++)
{
if (matice1[i] == null)
{
matice1[i] = new double[rozmer];
{
for (int e = 0; e < rozmer; e++)
{
matice1[i][e] = rnd.NextDouble();
}
}
}
}
timer.Stop();
res += timer.ElapsedMilliseconds.ToString();
int j = 0; int k = 0; int m = 0;
timer.Reset();
timer.Start();
for (j = 0; j < rozmer; j++)
{
for (k = 0; k < rozmer; k++)
{
temp = 0;
for (m = 0; m < rozmer; m++)
{
temp = temp + matice1[j][m] * matice1[m + rozmer][k];
}
matice1[j + r2][k] = temp;
}
}
timer.Stop();
this.Invoke((Action)delegate
{
this.Text = res + " : " + timer.ElapsedMilliseconds.ToString();
});
}));
thr.Start();
}
Es lento en C ++ porque no está utilizando multihilo. Básicamente, si A = BC, donde todas son matrices, la primera fila de A se puede calcular independientemente de la 2da fila, etc. Si A, B y C son todas n por matrices, puedes acelerar la multiplicación por un factor de n ^ 2, como
a_ {i, j} = suma_ {k} b_ {i, k} c_ {k, j}
Si usa, digamos, Eigen [ http://eigen.tuxfamily.org/dox/GettingStarted.html ], el multiproceso está incorporado y el número de subprocesos es ajustable.
Este tipo de pregunta es recurrente y debe responderse con más claridad que "Matlab usa bibliotecas altamente optimizadas" o "Matlab usa el MKL" por primera vez en .
Historia:
La multiplicación de matrices (junto con la matriz-vector, la multiplicación vector-vector y muchas de las descomposiciones de la matriz) es (son) los problemas más importantes en el algrebra lineal. Los ingenieros han estado resolviendo estos problemas con las computadoras desde los primeros días.
No soy un experto en la historia, pero aparentemente en aquel entonces, todo el mundo simplemente reescribió su versión de Fortran con bucles simples. Luego apareció cierta estandarización, con la identificación de "kernels" (rutinas básicas) que la mayoría de los problemas de álgebra lineal necesitaban para ser resueltos. Estas operaciones básicas se estandarizaron en una especificación llamada: Subprogramas básicos de álgebra lineal (BLAS). Los ingenieros podrían llamar a estas rutinas BLAS estándar y bien probadas en su código, haciendo su trabajo mucho más fácil.
BLAS:
BLAS evolucionó del nivel 1 (la primera versión que definía operaciones vectoriales escalares y vectoriales) al nivel 2 (operaciones vector-matriz) al nivel 3 (operaciones matriz-matriz) y proporcionó más y más "núcleos" más estandarizados. y más de las operaciones fundamentales de álgebra lineal. Las implementaciones originales de Fortran 77 todavía están disponibles en el sitio web de Netlib .
Hacia un mejor rendimiento:
Por lo tanto, a lo largo de los años (especialmente entre los lanzamientos de BLAS nivel 1 y nivel 2: principios de los 80), el hardware cambió, con el advenimiento de operaciones vectoriales y jerarquías de caché. Estas evoluciones permitieron aumentar sustancialmente el rendimiento de las subrutinas BLAS. Luego vinieron diferentes proveedores con su implementación de rutinas BLAS que eran cada vez más eficientes.
No conozco todas las implementaciones históricas (yo no nací ni era un niño en ese entonces), pero dos de las más notables aparecieron a principios de la década de 2000: Intel MKL y GotoBLAS. Su Matlab usa Intel MKL, que es un BLAS muy bueno y optimizado, y eso explica el gran rendimiento que ve.
Detalles técnicos sobre la multiplicación de Matrix:
Entonces, ¿por qué Matlab (el MKL) es tan rápido en dgemm (multiplicación de matriz matricial general de doble precisión)? En términos simples: porque usa vectorización y buen almacenamiento en caché de datos. En términos más complejos: vea el article proporcionado por Jonathan Moore.
Básicamente, cuando realiza su multiplicación en el código de C ++ que proporcionó, no es para nada amigable con el caché. Como sospecho que creó una matriz de punteros para filas de arrays, sus accesos en su bucle interno a la columna k-th de "matice2": matice2[m][k] son muy lentos. De hecho, cuando accedes a matice2[0][k] , debes obtener el k-ésimo elemento de la matriz 0 de tu matriz. Luego, en la siguiente iteración, debe acceder a matice2[1][k] , que es el k-ésimo elemento de otra matriz (la matriz 1). Luego, en la siguiente iteración accede a otra matriz, y así sucesivamente ... Dado que toda la matriz matice2 no puede caber en las memorias caché más altas (tiene 8*1024*1024 bytes de tamaño), el programa debe recuperar el elemento deseado de la fuente principal memoria, perdiendo mucho tiempo.
Si acabara de transponer la matriz, para que los accesos estuvieran en direcciones de memoria contiguas, su código ya se ejecutaría mucho más rápido porque ahora el compilador puede cargar filas completas en la memoria caché al mismo tiempo. Solo prueba esta versión modificada:
timer.start();
float temp = 0;
//transpose matice2
for (int p = 0; p < rozmer; p++)
{
for (int q = 0; q < rozmer; q++)
{
tempmat[p][q] = matice2[q][p];
}
}
for(int j = 0; j < rozmer; j++)
{
for (int k = 0; k < rozmer; k++)
{
temp = 0;
for (int m = 0; m < rozmer; m++)
{
temp = temp + matice1[j][m] * tempmat[k][m];
}
matice3[j][k] = temp;
}
}
timer.stop();
Para que pueda ver cómo la localidad de caché simplemente aumentó el rendimiento de su código bastante sustancialmente. Ahora las implementaciones de dgemm real explotan eso a un nivel muy extenso: realizan la multiplicación en bloques de la matriz definida por el tamaño del TLB (buffer de traducción lookaside, larga historia corta: lo que se puede almacenar en la memoria caché), de modo que se transmiten al procesa exactamente la cantidad de datos que puede procesar. El otro aspecto es la vectorización, usan las instrucciones vectorizadas del procesador para un rendimiento óptimo de la instrucción, lo que realmente no se puede hacer desde el código C ++ multiplataforma.
Finalmente, las personas que afirman que es debido al algoritmo de Strassen o Coppersmith-Winograd están equivocados, estos dos algoritmos no son implementables en la práctica, debido a las consideraciones de hardware mencionadas anteriormente.
La respuesta es que las bibliotecas LAPACK y BLAS hacen que MATLAB sea deslumbrantemente rápido en las operaciones de la matriz, no en ningún código de propiedad de la gente de MATLAB.
Use las bibliotecas LAPACK y / o BLAS en su código C ++ para las operaciones matriciales y debería obtener un rendimiento similar al de MATLAB. Estas bibliotecas deben estar disponibles gratuitamente en cualquier sistema moderno y las partes se desarrollaron durante décadas en la academia. Tenga en cuenta que hay varias implementaciones, incluidas algunas de código cerrado como Intel MKL .
Una discusión sobre cómo BLAS obtiene un alto rendimiento está disponible aquí.
Por cierto, es un dolor grave en mi experiencia llamar bibliotecas LAPACK directamente desde c (pero vale la pena). Necesita leer la documentación MUY precisa.
La respuesta general a "¿Por qué Matlab es más rápido haciendo xxx que otros programas?" Es que matlab tiene muchas funciones integradas y optimizadas.
Los otros programas que se utilizan a menudo no tienen estas funciones para que las personas apliquen sus propias soluciones creativas, que son sorprendentemente más lentas que el código optimizado profesionalmente.
Esto se puede interpretar de dos maneras:
1) La forma común / teórica: Matlab no es significativamente más rápido, solo está haciendo el benchmark incorrecto
2) La forma más realista: para esto, Matlab es más rápido en la práctica porque los lenguajes como c ++ se usan con demasiada facilidad de formas ineficaces.
MATLAB utiliza una implementación altamente optimizada de LAPACK de Intel, conocida como Intel Math Kernel Library (Intel MKL), específicamente la función dgemm . La velocidad Esta biblioteca aprovecha las características del procesador, incluidas las instrucciones SIMD y los procesadores multi-core. No documentan qué algoritmo específico usan. Si tuviera que llamar a Intel MKL desde C ++, debería ver un rendimiento similar.
No estoy seguro de qué biblioteca utiliza MATLAB para la multiplicación de GPU, pero probablemente algo así como nVidia CUBLAS .
Matlab incorporó LAPACK hace algún tiempo, así que supongo que su multiplicación de matriz usa algo al menos tan rápido. El código fuente y la documentación de LAPACK están disponibles.
También puede consultar el artículo de Goto y Van De Geijn "Anatomía de la multiplicación de matrices de alto rendimiento" en article
Esta es la razón por la cual MATLAB no realiza una multiplicación de matrices ingenua haciendo un bucle sobre cada elemento individual de la forma en que lo hizo en su código C ++.
Por supuesto, supongo que acaba de usar C=A*B lugar de escribir una función de multiplicación usted mismo.