En la realización de este reporte se aborda el desempeño para algoritmos que realizan operaciones con matrices y que son muy usados en el área de Algebra Lineal, dichos algoritmos son multiplicación de matrices de la forma Naive, multiplicación por Strassen, eliminación Gaussiana y eliminación de Gauss-Jordan. Se verá cómo es el desempeño para cada uno tanto con mediciones de tiempo como el número de operaciones que ejecuta cada uno a partir de un análisis teórico. Se hacen algunas aclaraciones de por qué se decide usar un factor de \(n\) ligeramente diferente para Strassen y se muestran las tablas de desempeño y total de operaciones para cada caso.
Todas las implementaciones fueron hechas usando Julia y usando el paquete de BenchmarkTools para medidas de desempeño.
El reporte esta estructurado de la siguiente forma
En esta sección se muestran los algoritmos implementados, de forma que se realizaron 2 códigos para cada categoría, la implementación Naive o simple y el algoritmo de strassen para multiplicaciones de matrices. Para el otro caso se implemetó la eliminación gaussiana y Gauss-Jordan. Ya que la eliminación gaussiana deja las matrices de forma escalonada reducida y hace sustituciones hacia atrás, mientras que Gauss-Jordan solo deja 1 en la diagonal a traves de restas en las filas de acuerdo al pivote.
function multiplicarMatrizCuadrada(A, B, C, n)
for i in 1:n
for j in 1:n
C[i,j] = 0
for k in 1:n
C[i,j] +=A[i,k]*B[k,j]
end
end
end
return C
end multiplicarMatrizCuadrada (generic function with 1 method)El algoritmo de Strassen solo está definido para matrices que tengan tamaño 2^n como lo indica en el paper original¹ por lo que las pruebas de desempeño en la siguiente sección se usarán matrices que tengan esos tamaños.
function strassen(A, B)
n = size(A, 1)
if n <= 64
return A * B
end
# Seperar matrices en cuadrantes
half = n ÷ 2
A11 = A[1:half, 1:half]
A12 = A[1:half, half+1:end]
A21 = A[half+1:end, 1:half]
A22 = A[half+1:end, half+1:end]
B11 = B[1:half, 1:half]
B12 = B[1:half, half+1:end]
B21 = B[half+1:end, 1:half]
B22 = B[half+1:end, half+1:end]
# Realizar productos de acuerdo al algoritmo de strassen
P1 = strassen(A11 + A22, B11 + B22)
P2 = strassen(A21 + A22, B11)
P3 = strassen(A11, B12 - B22)
P4 = strassen(A22, B21 - B11)
P5 = strassen(A11 + A12, B22)
P6 = strassen(A21 - A11, B11 + B12)
P7 = strassen(A12 - A22, B21 + B22)
# Realizar las operaciones en los cuadrantes
# de la matriz resultante
C11 = P1 + P4 - P5 + P7
C12 = P3 + P5
C21 = P2 + P4
C22 = P1 - P2 + P3 + P6
# Combinar los resultados en la matriz resultante
C = zeros(n, n)
C[1:half, 1:half] = C11
C[1:half, half+1:end] = C12
C[half+1:end, 1:half] = C21
C[half+1:end, half+1:end] = C22
return C
end
# A = rand(128, 128)
# B = rand(128, 128)
# C = strassen(A, B)strassen (generic function with 1 method)La implementación de eliminación gaussiana se uso de base el pseudocodigo encontrado en el libro de Métdos Numéricos²
function gaussiana(A,b,n)
# n = size(A, 1)
for k in 1:n-1
for i in k+1:n
factor = A[i,k] / A[k,k]
for j in k:n
A[i,j] = A[i,j] - factor * A[k,j]
end
b[i] = b[i] - factor * b[k]
end
end
#sustitucion hacia atraas
x = zeros(n)
x[n] = b[n] / A[n, n]
for i in n-1:-1:1
sum = b[i]
for j in i+1:n
sum = sum - A[i, j] * x[j]
end
x[i] = sum / A[i, i]
end
return x
end
# A = [3 -0.1 -0.2; 0.1 7 -0.3; 0.3 -0.2 10]
# b = [7.85, -19.3, 71.4]
# gaussiana(A,b)gaussiana (generic function with 2 methods)function gauss_jordan(A, b, n)
for i in 1:n
pivot = A[i, i]
A[i, :] = A[i, :] / pivot
b[i] = b[i] / pivot
for j in 1:n
if j != i
factor = A[j, i]
A[j, :] = A[j, :] - factor * A[i, :]
b[j] = b[j] - factor * b[i]
end
end
end
return b
end
# A = [3 -0.1 -0.2; 0.1 7 -0.3; 0.3 -0.2 10]
# b = [7.85, -19.3, 71.4]
# gauss_jordan(A, b, size(A, 1))gauss_jordan (generic function with 2 methods)En esta sección se estará midiendo el desempeño para los cuatro algoritmos a través de la librería de benchmarkTools para obtener mejores mediciones de tiempo. Solo en el caso del algoritmo de strassen se estarán usando tamaños de \(n=128, 512, 1024\) en lugar de \(n=100,300,1000\) debido a la naturaleza del algoritmo que no permite tamaños que no estén en factores de base 2.
En el caso de eliminaciones se corroboró que las matrices aleatorias generadas tuvieran solución y cumplieran con la característica que su determinante fuera diferente de cero para asegurar que es invertible la matriz A y por lo tanto asegurar que tiene solución el sistema.
También se comprobó que la solución proporcionada para el vector \(x\) tuviera un error mínimo comparado con el generado por la operación \ (que es el operador usado en Julia para matrices), encontrando valores muy pequeños del orden de 10^(-12) asegurando que los algoritmos son consistentes.
Todas las medidas de desempeño usaron la función tune!() para poder dar mejores estimaciones de tiempo.
using BenchmarkToolsn = 100
b = @benchmarkable multiplicarMatrizCuadrada(A, B, C, $n) setup=(A=rand(1:100, $n, $n); B=rand(1:100, $n, $n); C=zeros($n,$n))
tune!(b)
run(b)BenchmarkTools.Trial: 1609 samples with 1 evaluation per sample. Range (min … max): 2.560 ms … 4.243 ms ┊ GC (min … max): 0.00% … 0.00% Time (median): 2.875 ms ┊ GC (median): 0.00% Time (mean ± σ): 2.957 ms ± 250.917 μs ┊ GC (mean ± σ): 0.00% ± 0.00% ▂▄▇▆▃▇█▃▁▁ ▄▁▄▃▄▇▇██████████▇▇▆▇▆▇▇▇█▇▇▆▅▆▄▇▄▄▅▅▅▄▄▃▄▅▄▃▃▃▃▃▃▃▃▂▃▃▂▃▂▃ ▄ 2.56 ms Histogram: frequency by time 3.68 ms < Memory estimate: 0 bytes, allocs estimate: 0.
n=128
b = @benchmarkable strassen(A, B) setup=(A=rand(1:100, $n, $n); B=rand(1:100, $n, $n))
tune!(b)
run(b)BenchmarkTools.Trial: 5110 samples with 1 evaluation per sample. Range (min … max): 653.699 μs … 3.313 ms ┊ GC (min … max): 0.00% … 72.24% Time (median): 702.138 μs ┊ GC (median): 0.00% Time (mean ± σ): 815.508 μs ± 271.428 μs ┊ GC (mean ± σ): 6.75% ± 11.98% █▄▄▅▄▄▃▁▁ ▁▁ ▁ █████████████▇▇▇▇▇▅▆▇▇▇▆▅▅▅▄▄▆███▇███▇▇▇█▇▇▇▇▇▆▆▆▆▆▅▆▅▆▅▄▅▅▄▅ █ 654 μs Histogram: log(frequency) by time 1.77 ms < Memory estimate: 1.13 MiB, allocs estimate: 101.
using LinearAlgebra
n=100
A = rand(Float64, n, n)
while det(A) == 0
A = rand(Float64, n, n)
end
b = rand(Float64, n)
# x = gaussiana(A, b)
# x2 = A \ b
# println("Factor de Error: ", norm(A * x - b))
# println("Factor de Error: ", norm(A * x2 - b))
# print(norm(A * x - b)-norm(A * x2 - b))
b = @benchmarkable gaussiana($A, b) setup=(b = rand(Float64, $n))
tune!(b)
run(b)BenchmarkTools.Trial: 1193 samples with 1 evaluation per sample. Range (min … max): 3.900 ms … 5.633 ms ┊ GC (min … max): 0.00% … 0.00% Time (median): 4.088 ms ┊ GC (median): 0.00% Time (mean ± σ): 4.188 ms ± 291.297 μs ┊ GC (mean ± σ): 0.00% ± 0.00% ▅ █▆▅█▂▁ ▁ █████████▇▅▆▆███▇▆▅▅▄▆▄▄▄▃▃▃▃▃▂▃▂▃▂▂▂▂▃▂▁▁▂▃▃▂▂▂▃▂▂▃▃▂▂▁▂▃▃ ▄ 3.9 ms Histogram: frequency by time 5.28 ms < Memory estimate: 928 bytes, allocs estimate: 2.
A = rand(Float64, n, n)
while det(A) == 0
A = rand(Float64, n, n)
end
# x = gauss_jordan(A, b,n)
# x2 = A \ b
# println("Factor de error: ", norm(A * x - b))
# println("Factor de error: ", norm(A * x2 - b))
b = @benchmarkable gauss_jordan($A, B) setup=(B = rand(Float64, $n))
tune!(b)
run(b)BenchmarkTools.Trial: 928 samples with 1 evaluation per sample. Range (min … max): 4.220 ms … 54.857 ms ┊ GC (min … max): 0.00% … 90.63% Time (median): 4.636 ms ┊ GC (median): 12.97% Time (mean ± σ): 5.381 ms ± 3.766 ms ┊ GC (mean ± σ): 15.59% ± 5.97% ▇█▆▃▁ ███████▅▅▅▁▄▁▁▄▁▁▄▁▁▁▁▁▁▁▁▁▁▄▁▁▄▁▁▁▁▁▁▁▁▁▄▁▁▁▁▄▁▁▁▁▄▁▅▆▅▄▅ ▇ 4.22 ms Histogram: log(frequency) by time 23.6 ms < Memory estimate: 35.30 MiB, allocs estimate: 79605.
n=300300b = @benchmarkable multiplicarMatrizCuadrada(A, B, C, $n) setup=(A=rand(1:100, $n, $n); B=rand(1:100, $n, $n); C=zeros($n,$n))
tune!(b)
run(b)BenchmarkTools.Trial: 63 samples with 1 evaluation per sample. Range (min … max): 76.488 ms … 85.424 ms ┊ GC (min … max): 0.00% … 0.00% Time (median): 78.010 ms ┊ GC (median): 0.00% Time (mean ± σ): 78.862 ms ± 2.203 ms ┊ GC (mean ± σ): 0.00% ± 0.00% █▂ ▂ █▆▆▆████▄███▁▄▄▆▆▁▁▄▆▄▁▁▁▄▁▄▄▄▄▁▁█▁▁▁▄▁▄▄▁▁▆▁▁▁▄▁▁▁▁▁▁▁▁▁▁▄ ▁ 76.5 ms Histogram: frequency by time 85.2 ms < Memory estimate: 0 bytes, allocs estimate: 0.
n=512
b = @benchmarkable strassen(A, B) setup=(A=rand(1:100, $n, $n); B=rand(1:100, $n, $n))
tune!(b)
run(b)BenchmarkTools.Trial: 72 samples with 1 evaluation per sample. Range (min … max): 40.152 ms … 185.792 ms ┊ GC (min … max): 2.36% … 67.21% Time (median): 47.352 ms ┊ GC (median): 2.63% Time (mean ± σ): 66.107 ms ± 44.573 ms ┊ GC (mean ± σ): 31.27% ± 25.59% █▂ ████▅▄▄▄▄▁▁▃▄▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▃▁▁▁▁▁▆▃▁▃▃▁▁▃ ▁ 40.2 ms Histogram: frequency by time 182 ms < Memory estimate: 95.14 MiB, allocs estimate: 5589.
n=300
A = rand(Float64, n, n)
while det(A) == 0
A = rand(Float64, n, n)
end
b = @benchmarkable gaussiana($A, B) setup=(B = rand(Float64, $n))
tune!(b)
run(b)BenchmarkTools.Trial: 20 samples with 1 evaluation per sample. Range (min … max): 247.601 ms … 279.131 ms ┊ GC (min … max): 0.00% … 0.00% Time (median): 253.197 ms ┊ GC (median): 0.00% Time (mean ± σ): 255.660 ms ± 8.987 ms ┊ GC (mean ± σ): 0.00% ± 0.00% █▃▃ ▃ ███▁▁▁▇▁▇▁█▇▁▇▁▁▁▇▁▁▁▇▁▁▇▁▁▁▁▇▁▇▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▇▁▁▁▁▁▇ ▁ 248 ms Histogram: frequency by time 279 ms < Memory estimate: 2.44 KiB, allocs estimate: 3.
A = rand(Float64, n, n)
while det(A) == 0
A = rand(Float64, n, n)
end
b = @benchmarkable gauss_jordan($A, B) setup=(B = rand(Float64, $n))
tune!(b)
run(b)BenchmarkTools.Trial: 13 samples with 1 evaluation per sample. Range (min … max): 350.985 ms … 530.231 ms ┊ GC (min … max): 27.80% … 17.52% Time (median): 360.742 ms ┊ GC (median): 27.89% Time (mean ± σ): 385.925 ms ± 50.920 ms ┊ GC (mean ± σ): 27.21% ± 3.10% █ ▃ ▃ ▇█▇█▇▁▇▁▁▁▁▁▁▁▁▁▁▇▁▁▁▁▁▁▁█▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▇ ▁ 351 ms Histogram: frequency by time 530 ms < Memory estimate: 856.20 MiB, allocs estimate: 1078206.
# size(A,1)n = 1000
b = @benchmarkable multiplicarMatrizCuadrada(A, B, C, $n) setup=(A=rand(1:100, $n, $n); B=rand(1:100, $n, $n); C=zeros($n,$n))
tune!(b)
run(b)BenchmarkTools.Trial: 2 samples with 1 evaluation per sample. Range (min … max): 3.452 s … 3.455 s ┊ GC (min … max): 0.00% … 0.00% Time (median): 3.453 s ┊ GC (median): 0.00% Time (mean ± σ): 3.453 s ± 1.903 ms ┊ GC (mean ± σ): 0.00% ± 0.00% █ █ █▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁█ ▁ 3.45 s Histogram: frequency by time 3.45 s < Memory estimate: 0 bytes, allocs estimate: 0.
n=1024
b = @benchmarkable strassen(A, B) setup=(A=rand(1:100, $n, $n); B=rand(1:100, $n, $n))
tune!(b)
run(b)BenchmarkTools.Trial: 13 samples with 1 evaluation per sample. Range (min … max): 309.369 ms … 563.416 ms ┊ GC (min … max): 2.55% … 1.03% Time (median): 339.075 ms ┊ GC (median): 3.45% Time (mean ± σ): 400.824 ms ± 88.966 ms ┊ GC (mean ± σ): 15.11% ± 13.98% █ ▁▁▁▁▁ ▁▁ ▁ ▁ ▁ ▁ █▁▁█████▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁██▁▁▁█▁▁▁▁█▁█▁▁▁▁▁▁▁▁▁▁▁▁▁▁█ ▁ 309 ms Histogram: frequency by time 563 ms < Memory estimate: 724.00 MiB, allocs estimate: 39203.
n=1000
A = rand(Float64, n, n)
while det(A) == 0
A = rand(Float64, n, n)
end
b = @benchmarkable gaussiana($A, B) setup=(B = rand(Float64, $n))
tune!(b)
run(b)BenchmarkTools.Trial: 5 samples with 1 evaluation per sample. Range (min … max): 1.131 s … 1.222 s ┊ GC (min … max): 0.00% … 0.00% Time (median): 1.175 s ┊ GC (median): 0.00% Time (mean ± σ): 1.180 s ± 36.526 ms ┊ GC (mean ± σ): 0.00% ± 0.00% █ █ █ █ █ █▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁█▁▁▁▁▁▁█▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁█▁▁▁▁▁▁█ ▁ 1.13 s Histogram: frequency by time 1.22 s < Memory estimate: 7.88 KiB, allocs estimate: 3.
A = rand(Float64, n, n)
while det(A) == 0
A = rand(Float64, n, n)
end
b = @benchmarkable gauss_jordan($A, B) setup=(B = rand(Float64, $n))
tune!(b)
run(b)BenchmarkTools.Trial: 1 sample with 1 evaluation per sample. Single result which took 13.019 s (10.35% GC) to evaluate, with a memory estimate of 30.27 GiB, over 11994006 allocations.
Se muestran a continuación los tiempos medidos para cada algoritmo de acuerdo al valor median que arroja benchamark, ya que en su documentación menciona que es menos proponeso a dar valores atípicos³.
| Algoritmo | n=100 | n=300 | n=1000 |
|---|---|---|---|
| Multiplicación Simple | 2.875 ms | 78.010 ms | 3.453 s |
| Elimnación Gaussiana | 4.088 ms | 253.197 ms | 1.175 s |
| Eliminación Gauss-Jordan | 4.636 ms | 360.742 ms | 13.019 s |
| Algoritmo | n=128 | n=512 | n=1024 |
|---|---|---|---|
| Strassen | 702.138 μs | 47.352 ms | 339.075 ms |
Para cada algoritmo hacemos un análisis del total de operaciones aritméticas y las resumiremos en una tabla. Usaremos la notación \(T(n)\) para expresar el número de operaciones requeridas.
Se tiene que el número de multiplicaciones para cada elemento de la matriz C se obtiene como \(n^2\times n\) ya que en los loops mas internos son para recorrer las matrices A y B y el for externo es para recorrer todos los elementos de C.
Ahora para el caso de las sumas se tiene que se son \(n-1\) ya que se suman todos los elementos con todos menos excepto para la primer posición que sería sumar una posición consigo mismo. Y considerando que se hace para cada posición de C se tiene que \(n^2(n-1)\)
Tenemos entonces que el total de operaciones aritméticas está dado por
\[T(n) = n^3 + n^2(n-1)\] \[T(n) = 2n^3 + n^2 \]
Para este método se tiene que el algoritmo divide las matrices quedando en submatrices de tamaños de \(\frac{n}{2}\), luego hace 7 operaciones recursivas quedando como \(7T(\frac{n}{2})\) y realiza un total de 18 operaciones entre sumas y restas para submatriz.
Por lo que se puede expresar la siguiente relación de recurrencia
\[ T(n) = 7T(\frac{n}{2}) + 18(\frac{n}{2})^2\]
En la que el término de \(18(\frac{n}{2})^2\) hace referencia al tamaño total de la submatriz.
Es importante aclarar que el algoritmo al encontrar tamaños de \(n \leq 64\) hace una multiplicación estándar, por lo que para el caso de \(n=100\) el número total de operaciones aritméticas está dado por el siguiente desarrollo
\[T(100) = 7T(50) + 18(50)^2\]
\[T(50) = 7T(25) + 18(25)^2\]
\[T(25) = 2(25)^3 - 25^2\]
Analizamos para cada linea de codigo el total de operaciones aritmeticas
factor = A[i,k] / A[k,k]
Se contabiliza de forma que tiende a \((n-1)+(n-2)+\dots+1 = \frac{n(n-1)}{2}\)
A[i,j] = A[i,j] - factor * A[k,j]
Se puede expresar como \[\sum_{k=1}^{n-1}(n-k)(n-k+1)\] \[\sum_{k=1}^{n-1}(n^2 -2nk +n -k^2 -k)\]
Si desarrollamos cada término de acuerdo a su serie de convergencia tenemos lo siguientes
\[(n^3 -n^2) -(n^3 - n^2) +(n^2 - n)-(\frac{2n^3 -3n^2+n}{6})-(\frac{n^2-n}{2}) = \]
Simplificando tenemos que
\[\frac{n^2-n}{2} - \frac{n^3}{3} +\frac{n^2}{2} - \frac{n}{6} = \] \[-\frac{n^3}{3} - \frac{2n}{3} \]
Por lo que para una \(n\) muy grande se reduce únicamente al término de mayor grado \[ \frac{n^3}{3}\]
b[i] = b[i] - factor * b[k]
Se ejecuta \((n-1)+(n-2)+\dots+1 = \frac{n(n-1)}{2}\)
x[n] = b[n] / A[n, n]
Se ejecuta \(n\)
sum = sum - A[i, j] * x[j]
Se ejecuta \[\sum_{i=1}^{n-1}(n-i)=\frac{n(n-1)}{2}\]
Sumando todos los términos anteriores tenemos
\[ \frac{n(n-1)}{2} + \frac{n^3}{3}+ \frac{n(n-1)}{2} + n \frac{n(n-1)}{2} = \] \[ \frac{n^3}{3}+ \frac{3(n^2-n)}{2} - \frac{3n}{2}+ n = \] \[ \frac{n^3}{3}+ \frac{3n^2}{2} + \frac{n}{2}\]
Siendo que para una \(n\) muy grande el total de operaciones se tomen en cuenta de acuerdo al término de mayor orden quedando como \[ T(n) = \frac{n^3}{3} + \frac{3n^2}{2} + \frac{n}{2}\]
Calculamos para cada linea de código - A[i, :] = A[i, :] / pivot y b[i] = b[i] / pivot Se tienen \(n+1\) divisiones por cada fila por lo que se tiene un total de \[n(n+1)\]
A[j, :] = A[j, :] - factor * A[i, :]
Se puede expresar como
\[ n(n-1)n = n^3 - n^2\]
b[j] = b[j] - factor * b[i]
Se puede expresar como
\[ n(n-1) = n^2 - n\]
Sumando todos los terminos tenemos que el total de operaciones es\[T(n) = (n^2 +n) +(n^3 - n^2) + (n^2 -n)\] \[T(n) = n^3 + n^2\]
Se tienen entonces las siguientes expresiones para calcular el total de operaciones para cada algoritmo
En la siguiente tabla se muestran el número de operaciones para cada algoritmo dependiendo el tamaño de \(n\), de acuerdo a las expresiones obtenidas anteriormente
| Algoritmo | n=100 | n=300 | n=1000 |
|---|---|---|---|
| Multiplicación Simple | 2010000 | 54090000 | 2001000000 |
| Elimnación Gaussiana | 348383.33 | 9135150 | 1001500500 |
| Eliminación Gauss-Jordan | 1010000 | 27090000 | 1001000000 |
| Algoritmo | n=128 | n=512 | n=1024 |
|---|---|---|---|
| Strassen | 1880064 | 95367168 | 672288768 |
A continuación observamos que después de calcular tiempos y operaciones hay evidencias de que los algoritmos que tienen mayor costo computacional son los de multiplicación simple y la eliminación de gauss-Jordan, ya que éstos tienen en general un mayor número de operaciones a ejecutar y sus tiempos de ejecución fueron superiores contra Strassen y eliminación gaussiana respectivamente.
Por otro lado a pesar de que para el algoritmo de Strassen se hizo experimentos con una \(n\) ligeramente mayor para cada caso, se tuvieron tiempos de ejecución menores, y muchas menos operaciones calculadas para el peor caso.
Otra observación importante fue que al momento de realizar los cálculos de operaciones encontramos que el número de cifras significativas eran las que se obtenían a partir del término de mayor grado por lo que nuestras expresiones anteriormente encontradas se pueden reducir al término de mayor grado que era una \(n^3\) o usando la notación \(O(n^3)\) exceptuando Strassen ya que ese es una relación de recurrencia.
| Algoritmo | t | n=100 | t | n=300 | t | n=1000 |
|---|---|---|---|---|---|---|
| Multiplicación Simple | 2.875 ms | 2010000 | 78.010 ms | 54090000 | 3.453 s | 2001000000 |
| Elimnación Gaussiana | 4.088 ms | 348383.33 | 253.197 ms | 9135150 | 1.175 s | 1001500500 |
| Eliminación Gauss-Jordan | 4.636 ms | 1010000 | 360.742 ms | 27090000 | 13.019 s | 1001000000 |
| Algoritmo | t | n=128 | t | n=512 | t | n=1024 |
|---|---|---|---|---|---|---|
| Strassen | 702.138 us | 1880064 | 47.352 ms | 95367168 | 339.075 ms | 672288768 |
Con los algoritmos vistos notamos que el costo computacional es muy alto en cuanto a accesos de memoria, ya que como incluso en las pruebas de desempeño para una n=1000 nos reportaba un consumo alto memoria y de Garbage Collector, suponiendo que el tamaño de \(n\) crezca indefinidamente llevaría al desborde de memoria haciendo que el sistema se sature y quede sin espacio para ejecutar otros procesos. El hecho de tener que acceder a grandes segmentos de memoria hace lento al algoritmo y por lo tanto no es óptimo para tamaños de problemas grandes, en estos casos es recomendable ya sea generando mejores algoritmos o encontrar alternativas en la paralelización del mismo algoritmo, aunque eso implica otros componentes de hardware como GPUs.
El hecho de plantear el uso de matrices dispersas como se menciona en el libro de Algorithms for sparse linear systems⁴ tiene que tener ciertas propiedade las matrices como el hecho de que tenga algunos elementos cero para reducir la cantidad de operaciones, ya que al multiplicar columas o filas por elementos ceros no es necesario llevar a cabo la evaluación, haciendo que no sea necesaria la asignación de memoria para esas posiciones en cero, y reacomodando las matrices para que se puedan hacer las operaciones de forma correcta. Pero, para determinar los costos es necesario tener en cuenta cómo están acomodados los elementos cero que ayudarán a reducir elo número de operaciones.
[1] Strassen, V. (1969). Gaussian elimination is not optimal. Numerische mathematik, 13(4), 354-356.
[2] Chapra, S. C., Canale, R. P., Ruiz, R. S. G., Mercado, V. H. I., Díaz, E. M., & Benites, G. E. (2011). Métodos numéricos para ingenieros (Vol. 5, pp. 154-196). New York, NY, USA: McGraw-Hill.
[3] https://juliaci.github.io/BenchmarkTools.jl/stable/manual/
[4] Scott, J., & Tůma, M. (2023). Algorithms for sparse linear systems (p. 242). Springer Nature.