Add algebra benchmarks

Author: janbruedigam

benchmark/adjacency_matrix.jl

@@ -0,0 +1,75 @@
+A = [
+    0 1 0 1 1 0 1 0 1 0
+    1 0 1 0 0 0 0 0 0 0
+    0 1 0 0 0 0 0 0 0 0
+    1 0 0 0 1 0 0 0 0 0
+    1 0 0 1 0 1 0 0 0 0
+    0 0 0 0 1 0 0 0 0 0
+    1 0 0 0 0 0 0 1 0 0
+    0 0 0 0 0 0 1 0 1 1
+    1 0 0 0 0 0 0 1 0 0
+    0 0 0 0 0 0 0 1 0 0
+]
+
+B = [
+    0 1 0 0 0 1 0 0 0
+    1 0 0 0 0 1 0 0 0
+    0 0 0 0 0 1 1 0 1
+    0 0 0 0 0 0 1 1 0
+    0 0 0 0 0 0 0 1 1
+    1 1 1 0 0 0 1 0 0
+    0 0 1 1 0 1 0 1 1
+    0 0 0 1 1 0 1 0 1
+    0 0 1 0 1 0 1 1 0
+]
+
+C = [
+    0 1 1 0 0 0
+    1 0 0 1 0 0
+    1 0 0 1 1 0
+    0 1 1 0 0 1
+    0 0 1 0 0 1
+    0 0 0 1 1 0
+]
+
+D = [
+    0 1 1 0 1 0
+    1 0 0 1 0 1
+    1 0 0 1 0 0
+    0 1 1 0 0 0
+    1 0 0 0 0 1
+    0 1 0 0 1 0]
+
+ZAA = zeros(Int64, 10, 10)
+ZAB = zeros(Int64, 10, 9)
+ZAC = zeros(Int64, 10, 6)
+ZAD = zeros(Int64, 10, 6)
+ZBB = zeros(Int64, 9, 9)
+ZBC = zeros(Int64, 9, 6)
+ZBD = zeros(Int64, 9, 6)
+ZCC = zeros(Int64, 6, 6)
+ZCD = zeros(Int64, 6, 6)
+
+ZBA = ZAB'
+ZCA = ZAC'
+ZDA = ZAD'
+ZCB = ZBC'
+ZDB = ZBD'
+ZDC = ZCD'
+
+
+# Graph 1 is disconnected, 2-3-4-5 are connected, 6 is disconnected, 7 is disconnected
+
+A = [
+    A   ZAA ZAB ZAC ZAA ZAA ZAD
+    ZAA A   ZAB ZAC ZAA ZAA ZAD
+    ZBA ZBA B   ZBC ZBA ZBA ZBD
+    ZCA ZCA ZCB C   ZCA ZCA ZCD
+    ZAA ZAA ZAB ZAC A   ZAA ZAD
+    ZAA ZAA ZAB ZAC ZAA A   ZAD
+    ZDA ZDA ZDB ZDC ZDA ZDA D
+]
+
+A[13,21] = A[21,13] = 1
+A[16,30] = A[30,16] = 1
+A[20,36] = A[36,20] = 1

benchmark/example_benchmark.jl

@@ -1,108 +1,33 @@
 using GraphBasedSystems
 using LinearAlgebra
 
-A = [
-    0 1 0 1 1 0 1 0 1 0
-    1 0 1 0 0 0 0 0 0 0
-    0 1 0 0 0 0 0 0 0 0
-    1 0 0 0 1 0 0 0 0 0
-    1 0 0 1 0 1 0 0 0 0
-    0 0 0 0 1 0 0 0 0 0
-    1 0 0 0 0 0 0 1 0 0
-    0 0 0 0 0 0 1 0 1 1
-    1 0 0 0 0 0 0 1 0 0
-    0 0 0 0 0 0 0 1 0 0
-]
+include("adjacency_matrix.jl")
 
-B = [
-    0 1 0 0 0 1 0 0 0
-    1 0 0 0 0 1 0 0 0
-    0 0 0 0 0 1 1 0 1
-    0 0 0 0 0 0 1 1 0
-    0 0 0 0 0 0 0 1 1
-    1 1 1 0 0 0 1 0 0
-    0 0 1 1 0 1 0 1 1
-    0 0 0 1 1 0 1 0 1
-    0 0 1 0 1 0 1 1 0
-]
-
-C = [
-    0 1 1 0 0 0
-    1 0 0 1 0 0
-    1 0 0 1 1 0
-    0 1 1 0 0 1
-    0 0 1 0 0 1
-    0 0 0 1 1 0
-]
-
-D = [
-    0 1 1 0 1 0
-    1 0 0 1 0 1
-    1 0 0 1 0 0
-    0 1 1 0 0 0
-    1 0 0 0 0 1
-    0 1 0 0 1 0]
-
-ZAA = zeros(Int64,10,10)
-ZAB = zeros(Int64,10,9)
-ZAC = zeros(Int64,10,6)
-ZAD = zeros(Int64,10,6)
-ZBB = zeros(Int64,9,9)
-ZBC = zeros(Int64,9,6)
-ZBD = zeros(Int64,9,6)
-ZCC = zeros(Int64,6,6)
-ZCD = zeros(Int64,6,6)
-
-ZBA = ZAB'
-ZCA = ZAC'
-ZDA = ZAD'
-ZCB = ZBC'
-ZDB = ZBD'
-ZDC = ZCD'
-
-
-# Graph 1 is disconnected, 2-3-4-5 are connected, 6 is disconnected, 7 is disconnected
-
-A = [
-    A   ZAA ZAB ZAC ZAA ZAA ZAD
-    ZAA A   ZAB ZAC ZAA ZAA ZAD
-    ZBA ZBA B   ZBC ZBA ZBA ZBD
-    ZCA ZCA ZCB C   ZCA ZCA ZCD
-    ZAA ZAA ZAB ZAC A   ZAA ZAD
-    ZAA ZAA ZAB ZAC ZAA A   ZAD
-    ZDA ZDA ZDB ZDC ZDA ZDA D
-]
-
-A[13,21] = A[21,13] = 1
-A[16,30] = A[30,16] = 1
-A[20,36] = A[36,20] = 1
 
 function initialize!_posdef!(system::System{N}) where N 
     initialize!(system,rand)
     for i=1:N
-        system.matrix_entries[i,i].value += 1000*I 
+        system.matrix_entries[i,i].value += 10*I 
     end
+    # display(isposdef(full_matrix(system)))
 end
 
 
 system = System{Float64}(A, ones(Int,size(A)[1])*3)
 systemldlt = System{Float64}(A, ones(Int,size(A)[1])*3, symmetric=true)
 systemllt = System{Float64}(A, ones(Int,size(A)[1])*3, symmetric=true)
 
-SUITE["ldu"] = @benchmarkable ldu_solve!($system) samples=2 setup=(initialize!($system))
-SUITE["lu"] = @benchmarkable lu_solve!($system) samples=2 setup=(initialize!($system))
-SUITE["ldlt"] = @benchmarkable ldlt_solve!($systemldlt) samples=2 setup=(initialize!($systemldlt))
-SUITE["llt"] = @benchmarkable llt_solve!($systemllt) samples=2 setup=(initialize!_posdef!($systemllt))
-
-# A = [
-#     0 1 1 1 1
-#     1 0 1 1 1
-#     1 1 0 1 1
-#     1 1 1 0 1
-#     1 1 1 1 0
-# ]
-
-
-
-
-
+SUITE["sparse_ldu"] = @benchmarkable ldu_solve!($system) setup=(initialize!($system))
+SUITE["sparse_lu"] = @benchmarkable lu_solve!($system) setup=(initialize!($system))
+SUITE["dense_lu"] = @benchmarkable lu(F)\f setup=(initialize!($system);F=full_matrix($system);f=full_vector($system))
+SUITE["sparse_ldlt"] = @benchmarkable ldlt_solve!($systemldlt) setup=(initialize!($systemldlt))
+SUITE["dense_ldlt"] = @benchmarkable bunchkaufman(F)\f setup=(initialize!($systemldlt);F=full_matrix($systemldlt);f=full_vector($systemldlt))
+SUITE["sparse_llt"] = @benchmarkable llt_solve!($systemllt) setup=(initialize!_posdef!($systemllt))
+SUITE["dense_llt"] = @benchmarkable cholesky(F)\f setup=(initialize!_posdef!($systemllt);F=full_matrix($systemllt);f=full_vector($systemllt))
+
+SUITE["sparse_add"] = @benchmarkable +($system,$system) setup=(initialize!($system))
+SUITE["dense_add"] = @benchmarkable (+(F,F);+(f+f)) setup=(initialize!($system);F=full_matrix($system);f=full_vector($system))
+SUITE["sparse_mul"] = @benchmarkable *($system,$system) setup=(initialize!($system))
+SUITE["dense_mul"] = @benchmarkable *(F,F) setup=(initialize!($system);F=full_matrix($system))
+SUITE["sparse_solve"] = @benchmarkable \($system,Bmat) setup=(initialize!($system);Bmat=deepcopy(system.matrix_entries);initialize!($system))
+SUITE["dense_solve"] = @benchmarkable \(F1,F2) setup=(initialize!($system);F2=full_matrix($system);initialize!($system);F1=full_matrix($system))

src/solvers/llt.jl

@@ -59,10 +59,7 @@ function llt_backsubstitution_d!(vector, diagonal, diagonal_inverse)
     if diagonal_inverse.isinverted
         vector.value = diagonal_inverse.value * vector.value
     else
-        invdiagonal = inv(diagonal.value)
-        diagonal_inverse.value = invdiagonal
-        diagonal_inverse.isinverted = true
-        vector.value = diagonal_inverse.value * vector.value
+        vector.value = diagonal.value \ vector.value
     end
     return
 end