Algoritmo de Strassen

Algoritmo recursivo para la multiplicación de matrices

En álgebra lineal , el algoritmo de Strassen , llamado así por Volker Strassen , es un algoritmo para la multiplicación de matrices . Es más rápido que el algoritmo de multiplicación de matrices estándar para matrices grandes, con una mejor complejidad asintótica , aunque el algoritmo ingenuo suele ser mejor para matrices más pequeñas. El algoritmo de Strassen es más lento que los algoritmos más rápidos conocidos para matrices extremadamente grandes, pero estos algoritmos galácticos no son útiles en la práctica, ya que son mucho más lentos para matrices de tamaño práctico. Para matrices pequeñas existen algoritmos aún más rápidos.

El algoritmo de Strassen funciona para cualquier anillo , como plus/multiplicación, pero no para todos los semianillos , como min-plus o el álgebra booleana , donde el algoritmo ingenuo todavía funciona, y por eso se llama multiplicación de matrices combinatorias.

Historia

Volker Strassen publicó este algoritmo por primera vez en 1969 y demostró así que el algoritmo general de multiplicación de matrices no era óptimo. ^[1] La publicación del algoritmo de Strassen dio lugar a más investigaciones sobre la multiplicación de matrices que condujeron tanto a límites inferiores asintóticamente como a límites superiores computacionales mejorados. ${\displaystyle n^{3}}$

Algoritmo

La columna de la izquierda visualiza los cálculos necesarios para determinar el resultado de una multiplicación de matrices 2x2 . La multiplicación de matrices ingenua requiere una multiplicación por cada "1" de la columna de la izquierda. Cada una de las otras columnas (M1-M7) representa una sola de las 7 multiplicaciones del algoritmo de Strassen. La suma de las columnas M1-M7 da el mismo resultado que la multiplicación de matrices completa de la izquierda.

Sean , dos matrices cuadradas sobre un anillo , por ejemplo matrices cuyas entradas son números enteros o números reales. El objetivo de la multiplicación de matrices es calcular el producto matricial . La siguiente exposición del algoritmo supone que todas estas matrices tienen tamaños que son potencias de dos (es decir, ), pero esto es solo conceptualmente necesario - si las matrices , no son del tipo , las filas y columnas "faltantes" se pueden llenar con ceros para obtener matrices con tamaños de potencias de dos - aunque las implementaciones reales del algoritmo no hacen esto en la práctica. ${\displaystyle A}$ ${\displaystyle B}$ ${\displaystyle {\mathcal {R}}}$ ${\displaystyle C=AB}$ ${\displaystyle A,\,B,\,C\in \operatorname {Matr} _{2^{n}\times 2^{n}}({\mathcal {R}})}$ ${\displaystyle A}$ ${\displaystyle B}$ ${\displaystyle 2^{n}\times 2^{n}}$

El algoritmo de Strassen divide y en matrices de bloques de igual tamaño . ${\displaystyle A}$ ${\displaystyle B}$ ${\displaystyle C}$

${\displaystyle A={\begin{bmatrix}A_{11}&A_{12}\\A_{21}&A_{22}\end{bmatrix}},\quad B={\begin{bmatrix}B_{11}&B_{12}\\B_{21}&B_{22}\end{bmatrix}},\quad C={\begin{bmatrix}C_{11}&C_{12}\\C_{21}&C_{22}\end{bmatrix}},\quad }$

con . El algoritmo ingenuo sería: ${\displaystyle A_{ij},B_{ij},C_{ij}\in \operatorname {Mat} _{2^{n-1}\times 2^{n-1}}({\mathcal {R}})}$

${\displaystyle {\begin{bmatrix}C_{11}&C_{12}\\C_{21}&C_{22}\end{bmatrix}}={\begin{bmatrix}A_{11}{\color {red}\times }B_{11}+A_{12}{\color {red}\times }B_{21}\quad &A_{11}{\color {red}\times }B_{12}+A_{12}{\color {red}\times }B_{22}\\A_{21}{\color {red}\times }B_{11}+A_{22}{\color {red}\times }B_{21}\quad &A_{21}{\color {red}\times }B_{12}+A_{22}{\color {red}\times }B_{22}\end{bmatrix}}.}$

Esta construcción no reduce el número de multiplicaciones: todavía se necesitan 8 multiplicaciones de bloques de matrices para calcular las matrices, la misma cantidad de multiplicaciones necesarias cuando se utiliza la multiplicación de matrices estándar. ${\displaystyle C_{ij}}$

El algoritmo de Strassen define en cambio nuevos valores:

${\displaystyle {\begin{aligned}M_{1}&=(A_{11}+A_{22}){\color {red}\times }(B_{11}+B_{22});\\M_{2}&=(A_{21}+A_{22}){\color {red}\times }B_{11};\\M_{3}&=A_{11}{\color {red}\times }(B_{12}-B_{22});\\M_{4}&=A_{22}{\color {red}\times }(B_{21}-B_{11});\\M_{5}&=(A_{11}+A_{12}){\color {red}\times }B_{22};\\M_{6}&=(A_{21}-A_{11}){\color {red}\times }(B_{11}+B_{12});\\M_{7}&=(A_{12}-A_{22}){\color {red}\times }(B_{21}+B_{22}),\\\end{aligned}}}$

utilizando sólo 7 multiplicaciones (una para cada ) en lugar de 8. Ahora podemos expresar en términos de : ${\displaystyle M_{k}}$ ${\displaystyle C_{ij}}$ ${\displaystyle M_{k}}$

${\displaystyle {\begin{bmatrix}C_{11}&C_{12}\\C_{21}&C_{22}\end{bmatrix}}={\begin{bmatrix}M_{1}+M_{4}-M_{5}+M_{7}\quad &M_{3}+M_{5}\\M_{2}+M_{4}\quad &M_{1}-M_{2}+M_{3}+M_{6}\end{bmatrix}}.}$

Repetimos recursivamente este proceso de división hasta que las submatrices degeneren en números (elementos del anillo ). Si, como se mencionó anteriormente, la matriz original tenía un tamaño que no era una potencia de 2, entonces el producto resultante tendrá cero filas y columnas al igual que y , y estas se eliminarán en este punto para obtener la matriz (más pequeña) que realmente queríamos. ${\displaystyle {\mathcal {R}}}$ ${\displaystyle A}$ ${\displaystyle B}$ ${\displaystyle C}$

Las implementaciones prácticas del algoritmo de Strassen cambian a métodos estándar de multiplicación de matrices para submatrices lo suficientemente pequeñas, para las cuales esos algoritmos son más eficientes. El punto de cruce particular para el cual el algoritmo de Strassen es más eficiente depende de la implementación específica y el hardware. Autores anteriores habían estimado que el algoritmo de Strassen es más rápido para matrices con anchos de 32 a 128 para implementaciones optimizadas. ^[2] Sin embargo, se ha observado que este punto de cruce ha aumentado en los últimos años, y un estudio de 2010 encontró que incluso un solo paso del algoritmo de Strassen a menudo no es beneficioso en las arquitecturas actuales, en comparación con una multiplicación tradicional altamente optimizada, hasta que los tamaños de matriz exceden 1000 o más, e incluso para tamaños de matriz de varios miles, el beneficio es típicamente marginal en el mejor de los casos (alrededor del 10% o menos). ^[3] Un estudio más reciente (2016) observó beneficios para matrices tan pequeñas como 512 y un beneficio de alrededor del 20%. ^[4]

Formulario Winograd

Es posible reducir el número de adiciones de matrices utilizando la siguiente forma descubierta por Winograd:

${\displaystyle {\begin{bmatrix}a&b\\c&d\end{bmatrix}}{\begin{bmatrix}A&C\\B&D\end{bmatrix}}={\begin{bmatrix}t+b{\color {red}\times }B&w+v+(a+b-c-d){\color {red}\times }D\\w+u+d{\color {red}\times }(B+C-A-D)&w+u+v\end{bmatrix}}}$

dónde . ${\displaystyle t=a{\color {red}\times }A,\;u=(c-a){\color {red}\times }(C-D),\;v=(c+d){\color {red}\times }(C-A),\;w=t+(c+d-a){\color {red}\times }(A+D-C)}$

Esto reduce el número de adiciones y sustracciones de matrices de 18 a 15. El número de multiplicaciones de matrices sigue siendo 7 y la complejidad asintótica es la misma. ^[5]

Complejidad asintótica

El esquema del algoritmo anterior mostró que se puede lograr el éxito con sólo 7, en lugar de las 8 tradicionales, multiplicaciones matriz-matriz para los subbloques de la matriz. Por otro lado, hay que hacer sumas y restas de bloques, aunque esto no es un problema para la complejidad general: sumar matrices de tamaño grande requiere sólo operaciones, mientras que la multiplicación es sustancialmente más costosa (tradicionalmente, operaciones de suma o multiplicación). ${\displaystyle N/2}$ ${\displaystyle (N/2)^{2}}$ ${\displaystyle 2(N/2)^{3}}$

La pregunta entonces es cuántas operaciones exactamente se necesitan para los algoritmos de Strassen, y cómo se compara esto con la multiplicación de matrices estándar que requiere aproximadamente (donde ) operaciones aritméticas, es decir, una complejidad asintótica . ${\displaystyle 2N^{3}}$ ${\displaystyle N=2^{n}}$ ${\displaystyle \Theta (N^{3})}$

El número de adiciones y multiplicaciones requeridas en el algoritmo de Strassen se puede calcular de la siguiente manera: sea el número de operaciones para una matriz. Luego, mediante la aplicación recursiva del algoritmo de Strassen, vemos que , para alguna constante que depende del número de adiciones realizadas en cada aplicación del algoritmo. Por lo tanto , , es decir, la complejidad asintótica para multiplicar matrices de tamaño utilizando el algoritmo de Strassen es . Sin embargo, la reducción en el número de operaciones aritméticas tiene el precio de una estabilidad numérica algo reducida , ^[6] y el algoritmo también requiere significativamente más memoria en comparación con el algoritmo ingenuo. Ambas matrices iniciales deben tener sus dimensiones expandidas a la siguiente potencia de 2, lo que resulta en almacenar hasta cuatro veces más elementos, y las siete matrices auxiliares contienen cada una una cuarta parte de los elementos de las expandidas. ${\displaystyle f(n)}$ ${\displaystyle 2^{n}\times 2^{n}}$ ${\displaystyle f(n)=7f(n-1)+l4^{n}}$ ${\displaystyle l}$ ${\displaystyle f(n)=(7+o(1))^{n}}$ ${\displaystyle N=2^{n}}$ ${\displaystyle O([7+o(1)]^{n})=O(N^{\log _{2}7+o(1)})\approx O(N^{2.8074})}$

El algoritmo de Strassen debe compararse con la forma "ingenua" de realizar la multiplicación de matrices que requeriría 8 multiplicaciones de subbloques en lugar de 7. Esto daría lugar a la complejidad que se espera del enfoque estándar: . La comparación de estos dos algoritmos muestra que asintóticamente , el algoritmo de Strassen es más rápido: Existe un tamaño tal que las matrices que son más grandes se multiplican de manera más eficiente con el algoritmo de Strassen que con la forma "tradicional". Sin embargo, la afirmación asintótica no implica que el algoritmo de Strassen sea siempre más rápido incluso para matrices pequeñas, y en la práctica, de hecho, este no es el caso: Para matrices pequeñas, el costo de las adiciones adicionales de bloques de matriz supera el ahorro en el número de multiplicaciones. También hay otros factores que no se capturan en el análisis anterior, como la diferencia de costo en el hardware actual entre cargar datos desde la memoria a los procesadores frente al costo de realizar realmente operaciones en estos datos. Como consecuencia de este tipo de consideraciones, el algoritmo de Strassen se suele utilizar únicamente en matrices "grandes". Este tipo de efecto es aún más pronunciado con algoritmos alternativos como el de Coppersmith y Winograd : aunque asintóticamente es incluso más rápido, el punto de cruce es tan grande que el algoritmo no suele utilizarse en matrices que se encuentran en la práctica. ${\displaystyle O(8^{n})=O(N^{\log _{2}8})=O(N^{3})}$ ${\displaystyle N_{\text{threshold}}}$ ${\displaystyle N_{\text{threshold}}}$

Complejidad de rango o bilineal

La complejidad bilineal o rango de una función bilineal es un concepto importante en la complejidad asintótica de la multiplicación de matrices. El rango de una función bilineal sobre un cuerpo F se define como (lo que es un abuso de notación ) ${\displaystyle \phi :\mathbf {A} \times \mathbf {B} \rightarrow \mathbf {C} }$

${\displaystyle R(\phi /\mathbf {F} )=\min \left\{r\left|\exists f_{i}\in \mathbf {A} ^{*},g_{i}\in \mathbf {B} ^{*},w_{i}\in \mathbf {C} ,\forall \mathbf {a} \in \mathbf {A} ,\mathbf {b} \in \mathbf {B} ,\phi (\mathbf {a} ,\mathbf {b} )=\sum _{i=1}^{r}f_{i}(\mathbf {a} )g_{i}(\mathbf {b} )w_{i}\right.\right\}}$

En otras palabras, el rango de una función bilineal es la longitud de su cálculo bilineal más corto. ^[7] La ​​existencia del algoritmo de Strassen muestra que el rango de la multiplicación de matrices no es mayor que siete. Para ver esto, expresemos este algoritmo (junto con el algoritmo estándar) como un cálculo bilineal de este tipo. En el caso de las matrices, los espacios duales A * y B * consisten en funciones en el cuerpo F inducidas por un producto escalar de doble punto , (es decir, en este caso la suma de todas las entradas de un producto de Hadamard .) ${\displaystyle 2\times 2}$

Se puede demostrar que el número total de multiplicaciones elementales necesarias para la multiplicación de matrices está estrechamente ligado asintóticamente al rango , es decir , o más específicamente, dado que las constantes son conocidas, . Una propiedad útil del rango es que es submultiplicativo para productos tensoriales , y esto permite demostrar que la multiplicación de matrices se puede lograr con no más que multiplicaciones elementales para cualquier . (Este producto tensorial -fold del mapa de multiplicación de matrices consigo mismo —una -ésima potencia tensorial— se realiza mediante el paso recursivo en el algoritmo mostrado). ${\displaystyle L}$ ${\displaystyle R}$ ${\displaystyle L=\Theta (R)}$ ${\displaystyle R/2\leq L\leq R}$ ${\displaystyle 2^{n}\times 2^{n}\times 2^{n}}$ ${\displaystyle 7n}$ ${\displaystyle n}$ ${\displaystyle n}$ ${\displaystyle 2\times 2\times 2}$ ${\displaystyle n}$

Comportamiento de la caché

El algoritmo de Strassen no tiene en cuenta la memoria caché . El análisis de su algoritmo de comportamiento de la memoria caché ha demostrado que incurre en

${\displaystyle \Theta \left(1+{\frac {n^{2}}{b}}+{\frac {n^{\log _{2}7}}{b{\sqrt {M}}}}\right)}$

La caché falla durante su ejecución, asumiendo una caché idealizada de tamaño (es decir, con líneas de longitud ). ^{[8] : 13} ${\displaystyle M}$ ${\displaystyle M/b}$ ${\displaystyle b}$

Consideraciones de implementación

La descripción anterior indica que las matrices son cuadradas y que el tamaño es una potencia de dos, y que se debe utilizar relleno si es necesario. Esta restricción permite que las matrices se dividan por la mitad, de forma recursiva, hasta que se alcance el límite de la multiplicación escalar. La restricción simplifica la explicación y el análisis de la complejidad, pero en realidad no es necesaria; ^[9] y, de hecho, rellenar la matriz como se describe aumentará el tiempo de cálculo y puede eliminar fácilmente los ahorros de tiempo relativamente limitados obtenidos al utilizar el método en primer lugar.

Una buena implementación observará lo siguiente:

• No es necesario ni deseable utilizar el algoritmo de Strassen hasta el límite de los escalares. En comparación con la multiplicación de matrices convencional, el algoritmo agrega una carga de trabajo considerable en las sumas y restas; por lo tanto, por debajo de un cierto tamaño, será mejor utilizar la multiplicación convencional. Así, por ejemplo, no es necesario rellenar a , ya que podría subdividirse en matrices y, luego, se puede utilizar la multiplicación convencional en ese nivel. ${\displaystyle O(n^{2})}$ ${\displaystyle 1600\times 1600}$ ${\displaystyle 2048\times 2048}$ ${\displaystyle 25\times 25}$
• El método puede aplicarse a matrices cuadradas de cualquier dimensión. ^[3] Si la dimensión es par, se dividen por la mitad como se describe. Si la dimensión es impar, primero se aplica un relleno de ceros en una fila y una columna. Dicho relleno se puede aplicar sobre la marcha y de forma diferida, y las filas y columnas adicionales se descartan a medida que se forma el resultado. Por ejemplo, supongamos que las matrices son . Se pueden dividir de modo que la parte superior izquierda sea y la inferior derecha sea . Siempre que las operaciones lo requieran, las dimensiones de se rellenan con ceros primero. Nótese, por ejemplo, que el producto solo se usa en la fila inferior de la salida, por lo que solo se requiere que tenga filas de alto; y, por lo tanto, el factor izquierdo usado para generarlo solo necesita tener filas de alto; en consecuencia, no hay necesidad de rellenar esa suma a filas; solo es necesario rellenar a columnas para que coincida con . ${\displaystyle 199\times 199}$ ${\displaystyle 100\times 100}$ ${\displaystyle 99\times 99}$ ${\displaystyle 99}$ ${\displaystyle 100}$ ${\displaystyle M_{2}}$ ${\displaystyle 99}$ ${\displaystyle A_{21}+A_{22}}$ ${\displaystyle 99}$ ${\displaystyle 100}$ ${\displaystyle A_{22}}$ ${\displaystyle 100}$ ${\displaystyle A_{21}}$

Además, no es necesario que las matrices sean cuadradas. Las matrices no cuadradas se pueden dividir por la mitad utilizando los mismos métodos, lo que da lugar a matrices no cuadradas más pequeñas. Si las matrices son suficientemente no cuadradas, valdrá la pena reducir la operación inicial a más productos cuadrados, utilizando métodos simples que son esencialmente , por ejemplo: ${\displaystyle O(n^{2})}$

• Un producto de tamaño se puede realizar como 20 operaciones separadas, ordenadas para formar el resultado; ${\displaystyle [2N\times N]\ast [N\times 10N]}$ ${\displaystyle [N\times N]\ast [N\times N]}$
• Un producto de tamaño se puede realizar como 10 operaciones separadas, sumadas para formar el resultado. ${\displaystyle [N\times 10N]\ast [10N\times N]}$ ${\displaystyle [N\times N]\ast [N\times N]}$

Estas técnicas harán que la implementación sea más complicada, en comparación con el simple relleno hasta un cuadrado de potencia de dos; sin embargo, es una suposición razonable que cualquiera que emprenda una implementación de Strassen, en lugar de la multiplicación convencional, colocará una mayor prioridad en la eficiencia computacional que en la simplicidad de la implementación.

En la práctica, el algoritmo de Strassen se puede implementar para lograr un mejor rendimiento que la multiplicación convencional incluso para matrices tan pequeñas como , para matrices que no son en absoluto cuadradas y sin requerir espacio de trabajo más allá de los buffers que ya se necesitan para una multiplicación convencional de alto rendimiento. ^[4] ${\displaystyle 500\times 500}$

Véase también

• Complejidad computacional de operaciones matemáticas
• Eliminación de Gauss-Jordan
• Complejidad computacional de la multiplicación de matrices
• Curva de orden Z
• Algoritmo de Karatsuba , para multiplicar números enteros de n dígitos en lugar de en el tiempo ${\displaystyle O(n^{\log _{2}3})}$ ${\displaystyle O(n^{2})}$
  • Un algoritmo de multiplicación complejo similar multiplica dos números complejos utilizando 3 multiplicaciones reales en lugar de 4.
• Algoritmo Toom-Cook , una generalización más rápida del algoritmo Karatsuba que permite la descomposición recursiva de dividir y vencer en más de 2 bloques a la vez

Referencias

1. Strassen, Volker (1969). "La eliminación gaussiana no es óptima". Numer. Math . 13 (4): 354–356. doi :10.1007/BF02165411. S2CID  121656251.
2. Skiena, Steven S. (1998), "§8.2.3 Multiplicación de matrices", The Algorithm Design Manual , Berlín, Nueva York: Springer-Verlag , ISBN 978-0-387-94860-7.
3. D'Alberto, Paolo; Nicolau, Alexandru (2005). Uso de la recursión para mejorar el rendimiento de ATLAS (PDF) . Sexto Simposio Internacional sobre Computación de Alto Rendimiento.
4. Huang, Jianyu; Smith, Tyler M.; Henry, Greg M.; van de Geijn, Robert A. (13 de noviembre de 2016). Strassen's Algorithm Reloaded. SC16: La conferencia internacional sobre computación de alto rendimiento, redes, almacenamiento y análisis. IEEE Press. págs. 690–701. doi :10.1109/SC.2016.58. ISBN. 9781467388153. Recuperado el 1 de noviembre de 2022 .
5. Knuth (1997), pág. 500.
6. Webb, Miller (1975). "Complejidad computacional y estabilidad numérica". SIAM J. Comput . 4 (2): 97–107. doi :10.1137/0204009.
7. Burgués; cláusulas; Shokrollahi (1997). Teoría de la complejidad algebraica . Springer-Verlag. ISBN 3-540-60582-7.
8. Frigo, M.; Leiserson, CE ; Prokop, H .; Ramachandran, S. (1999). Algoritmos ajenos a la caché (PDF) . Proc. Simposio IEEE sobre Fundamentos de la Ciencia de la Computación (FOCS). págs. 285–297.
9. Higham, Nicholas J. (1990). "Explotación de la multiplicación rápida de matrices dentro del nivel 3 de BLAS" (PDF) . ACM Transactions on Mathematical Software . 16 (4): 352–368. doi :10.1145/98267.98290. hdl : 1813/6900 . S2CID  5715053.

• Thomas H. Cormen , Charles E. Leiserson , Ronald L. Rivest y Clifford Stein . Introducción a los algoritmos , segunda edición. MIT Press y McGraw-Hill, 2001. ISBN 0-262-03293-7 . Capítulo 28: Sección 28.2: Algoritmo de Strassen para la multiplicación de matrices, págs. 735–741. 
• Knuth, Donald (1997). El arte de la programación informática, algoritmos seminuméricos . Vol. II (3.ª ed.). Addison-Wesley. ISBN 0-201-89684-2.

Enlaces externos

• Weisstein, Eric W. "Fórmulas de Strassen". MundoMatemático .(también incluye fórmulas para la inversión rápida de matrices )
• Tyler J. Earnest, Algoritmo de Strassen sobre el motor de banda ancha celular