Benchmark on dense random symmetric positive definite matrices
using RandomLinearAlgebraSolvers, Random, LinearAlgebra, SparseArrays
Random.seed!(1234)MersenneTwister(1234)Benchmark of methods to solve Ax = b with A an n x n randomly generated symmetric positive definite matrix, and b = A * xref.
using DataFrames, Printf, SolverBenchmark, StoppingN = 5 # number of problems
n = 100 # size of A: n x n100Names of solvers:
names = [:RandomizedKaczmarz, :RandomizedBlockKaczmarz, :RandomizedCD, :RandomizedCD2, :RandomizedNewton]5-element Vector{Symbol}:
:RandomizedKaczmarz
:RandomizedBlockKaczmarz
:RandomizedCD
:RandomizedCD2
:RandomizedNewton#Initialization of the DataFrame for n problems.
stats = Dict(name => DataFrame(
:id => 1:N,
:nvar => zeros(Int64, N),
:status => [:Unknown for i = 1:N],
:time => NaN*ones(N),
:iter => zeros(Int64, N),
:score => NaN*ones(N)) for name in names)Dict{Symbol, DataFrames.DataFrame} with 5 entries:
:RandomizedBlockKaczmarz => 5×6 DataFrame…
:RandomizedCD => 5×6 DataFrame…
:RandomizedCD2 => 5×6 DataFrame…
:RandomizedKaczmarz => 5×6 DataFrame…
:RandomizedNewton => 5×6 DataFrame…for i=1:N
Ar = rand(n,n); A = 1/sqrt(n) * Ar'*Ar + Matrix{Float64}(I, n, n)
xref = 100 * rand(n)
b = A * xref
x0 = zeros(size(A,2))
la_stop = RLAStopping(A, b, max_iter = 100000, rtol = sqrt(eps()), atol = sqrt(eps()))
for name in names
#solve the problem
reinit!(la_stop, rstate = true, x = x0, res = similar(b))
la_stop.meta.start_time = time()
@time eval(name)(la_stop, r = 80, is_zero_start=true)
sol = la_stop.current_state.x
#update the stats from the Stopping
stats[name].nvar[i] = n
stats[name].status[i] = status(la_stop)
stop_has_time = (la_stop.current_state.current_time != nothing)
stats[name].time[i] = stop_has_time ? la_stop.current_state.current_time - la_stop.meta.start_time : time() - la_stop.meta.start_time
stats[name].iter[i] = la_stop.meta.nb_of_stop
stats[name].score[i] = norm(la_stop.current_state.current_score, Inf)
end
end 0.709734 seconds (1.10 M allocations: 268.649 MiB, 18.27% gc time, 0.35% compilation time)
0.172047 seconds (1.58 k allocations: 18.731 MiB)
0.056013 seconds (800.01 k allocations: 12.207 MiB)
0.271580 seconds (1.11 M allocations: 29.017 MiB, 4.81% gc time, 70.01% compilation time)
1.940081 seconds (2.74 M allocations: 463.573 MiB, 3.91% gc time, 64.18% compilation time)
0.660372 seconds (1.10 M allocations: 268.557 MiB, 12.16% gc time)
0.179197 seconds (1.39 k allocations: 16.460 MiB, 15.98% gc time)
0.058878 seconds (800.01 k allocations: 12.207 MiB)
0.079600 seconds (800.01 k allocations: 12.207 MiB)
0.680473 seconds (383.89 k allocations: 334.943 MiB, 4.01% gc time)
0.677140 seconds (1.10 M allocations: 268.557 MiB, 14.52% gc time)
0.150991 seconds (1.49 k allocations: 17.595 MiB)
0.059845 seconds (800.01 k allocations: 12.207 MiB)
0.086932 seconds (800.01 k allocations: 12.207 MiB, 19.75% gc time)
0.660140 seconds (384.09 k allocations: 335.198 MiB, 3.27% gc time)
0.677236 seconds (1.10 M allocations: 268.557 MiB, 14.75% gc time)
0.174273 seconds (1.58 k allocations: 18.731 MiB)
0.058495 seconds (800.01 k allocations: 12.207 MiB)
0.077586 seconds (800.01 k allocations: 12.207 MiB)
0.681826 seconds (385.89 k allocations: 335.926 MiB, 4.65% gc time)
0.666673 seconds (1.10 M allocations: 268.557 MiB, 12.06% gc time)
0.174085 seconds (1.39 k allocations: 16.460 MiB, 16.02% gc time)
0.057900 seconds (800.01 k allocations: 12.207 MiB)
0.077279 seconds (800.01 k allocations: 12.207 MiB)
0.671613 seconds (384.69 k allocations: 335.453 MiB, 4.07% gc time)for name in names
@show name
@show stats[name]
endname = :RandomizedKaczmarz
stats[name] = 5×6 DataFrame
Row │ id nvar status time iter score
│ Int64 Int64 Symbol Float64 Int64 Float64
─────┼──────────────────────────────────────────────────────────
1 │ 1 100 IterationLimit 0.709741 100001 11.4904
2 │ 2 100 IterationLimit 0.660373 100001 10.3227
3 │ 3 100 IterationLimit 0.677143 100001 7.35573
4 │ 4 100 IterationLimit 0.677238 100001 9.91928
5 │ 5 100 IterationLimit 0.666675 100001 9.25658
name = :RandomizedBlockKaczmarz
stats[name] = 5×6 DataFrame
Row │ id nvar status time iter score
│ Int64 Int64 Symbol Float64 Int64 Float64
─────┼────────────────────────────────────────────────────
1 │ 1 100 Optimal 0.166995 33 1.21199e-8
2 │ 2 100 Optimal 0.174489 29 1.30985e-8
3 │ 3 100 Optimal 0.14652 31 8.51651e-9
4 │ 4 100 Optimal 0.169391 33 6.2355e-9
5 │ 5 100 Optimal 0.169176 29 4.05817e-9
name = :RandomizedCD
stats[name] = 5×6 DataFrame
Row │ id nvar status time iter score
│ Int64 Int64 Symbol Float64 Int64 Float64
─────┼──────────────────────────────────────────────────────────
1 │ 1 100 IterationLimit 0.056015 100001 334.969
2 │ 2 100 IterationLimit 0.0588779 100001 266.787
3 │ 3 100 IterationLimit 0.0598452 100001 321.47
4 │ 4 100 IterationLimit 0.058495 100001 304.792
5 │ 5 100 IterationLimit 0.0579019 100001 330.454
name = :RandomizedCD2
stats[name] = 5×6 DataFrame
Row │ id nvar status time iter score
│ Int64 Int64 Symbol Float64 Int64 Float64
─────┼──────────────────────────────────────────────────────────
1 │ 1 100 IterationLimit 0.271581 100001 28121.7
2 │ 2 100 IterationLimit 0.0796008 100001 29849.0
3 │ 3 100 IterationLimit 0.0869341 100001 27224.9
4 │ 4 100 IterationLimit 0.0775869 100001 29506.1
5 │ 5 100 IterationLimit 0.0772791 100001 29165.6
name = :RandomizedNewton
stats[name] = 5×6 DataFrame
Row │ id nvar status time iter score
│ Int64 Int64 Symbol Float64 Int64 Float64
─────┼─────────────────────────────────────────────────
1 │ 1 100 Optimal 1.93947 1010 NaN
2 │ 2 100 Optimal 0.679822 1009 NaN
3 │ 3 100 Optimal 0.659448 1010 NaN
4 │ 4 100 Optimal 0.681213 1010 NaN
5 │ 5 100 Optimal 0.670995 1010 NaNor run a performance profile:
using Plots
gr()
cost(df) = (df.status .!= :Optimal) * Inf + df.time
p = performance_profile(stats, cost)