julia> versioninfo()
Julia Version 0.7.0
Commit a4cb80f3ed (2018-08-08 06:46 UTC)
Platform Info:
OS: macOS (x86_64-apple-darwin14.5.0)
CPU: Intel(R) Core(TM) i7-4770HQ CPU @ 2.20GHz
WORD_SIZE: 64
LIBM: libopenlibm
LLVM: libLLVM-6.0.0 (ORCJIT, haswell)
julia> using ExponentialUtilities, Expokit, LazySets, SparseArrays
julia> # define an nxn tridiagonal matrix
A(n) = sparse(diagm(0 => fill(0.05, n), -1 => fill(-1, n-1), 1 => fill(-1, n-1)))
A (generic function with 1 method)
julia> # step size and initial set
δ = 0.1
0.1
julia> X0(n) = Ball2(ones(n), 0.1)
X0 (generic function with 1 method)
julia> # input coefficients matrix (nx2 matrix with coefficients from -1 to 1)
b(n) = vcat(range(-1, stop=1, length=n))
b (generic function with 1 method)
julia> B(n) = [b(n) b(n)]
B (generic function with 1 method)
julia> U = BallInf(zeros(2), 1.2)
BallInf{Float64}([0.0, 0.0], 1.2)
julia> # lazy matrix exponential
eAδ(n) = SparseMatrixExp(A(n) * δ)
eAδ (generic function with 1 method)
julia> # set that we want to overapproximate with an interval hull
Y(n) = ConvexHull(eAδ(n) * X0(n) ⊕ (δ * B(n) * U), X0(n))
Y (generic function with 1 method)
julia> A(10)
10×10 SparseMatrixCSC{Float64,Int64} with 28 stored entries:
[1 , 1] = 0.05
[2 , 1] = -1.0
[1 , 2] = -1.0
[2 , 2] = 0.05
[3 , 2] = -1.0
[2 , 3] = -1.0
[3 , 3] = 0.05
[4 , 3] = -1.0
[3 , 4] = -1.0
[4 , 4] = 0.05
[5 , 4] = -1.0
[4 , 5] = -1.0
[5 , 5] = 0.05
[6 , 5] = -1.0
[5 , 6] = -1.0
[6 , 6] = 0.05
[7 , 6] = -1.0
[6 , 7] = -1.0
[7 , 7] = 0.05
[8 , 7] = -1.0
[7 , 8] = -1.0
[8 , 8] = 0.05
[9 , 8] = -1.0
[8 , 9] = -1.0
[9 , 9] = 0.05
[10, 9] = -1.0
[9 , 10] = -1.0
[10, 10] = 0.05
julia> expmv(1.0, A(10), b(10))
10-element Array{Float64,1}:
-0.8536046924069195
0.3365829070231674
-0.2437063579211246
-0.003670541256143919
-0.02637707499315963
0.026377074993159617
0.003670541256143919
0.2437063579211246
-0.3365829070231674
0.8536046924069195
julia> @which expmv(1.0, A(10), b(10))
expmv(t::Number, A, vec::Array{T,1}) where T in Expokit at /Users/forets/.julia/packages/Expokit/fQbKs/src/expmv.jl:33
julia> @time expmv(1.0, A(10), b(10))
0.000237 seconds (259 allocations: 21.219 KiB)
10-element Array{Float64,1}:
-0.8536046924069195
0.3365829070231674
-0.2437063579211246
-0.003670541256143919
-0.02637707499315963
0.026377074993159617
0.003670541256143919
0.2437063579211246
-0.3365829070231674
0.8536046924069195
julia> @time expmv(1.0, A(10), b(10))
0.005259 seconds (259 allocations: 21.219 KiB)
10-element Array{Float64,1}:
-0.8536046924069195
0.3365829070231674
-0.2437063579211246
-0.003670541256143919
-0.02637707499315963
0.026377074993159617
0.003670541256143919
0.2437063579211246
-0.3365829070231674
0.8536046924069195
julia> @time expmv(1.0, A(10), b(10))
0.000577 seconds (259 allocations: 21.219 KiB)
10-element Array{Float64,1}:
-0.8536046924069195
0.3365829070231674
-0.2437063579211246
-0.003670541256143919
-0.02637707499315963
0.026377074993159617
0.003670541256143919
0.2437063579211246
-0.3365829070231674
0.8536046924069195
julia> @time expmv(1.0, A(10), b(10))
0.001631 seconds (259 allocations: 21.219 KiB)
10-element Array{Float64,1}:
-0.8536046924069195
0.3365829070231674
-0.2437063579211246
-0.003670541256143919
-0.02637707499315963
0.026377074993159617
0.003670541256143919
0.2437063579211246
-0.3365829070231674
0.8536046924069195
julia> @btime expmv(1.0, $A(10), $b(10))
ERROR: LoadError: UndefVarError: @btime not defined
in expression starting at REPL[17]:1
julia> using BenchmarkTools
julia> expv(1.0, A(10), b(10))
10-element Array{Float64,1}:
-0.8536046924069198
0.3365829070231669
-0.24370635792112402
-0.0036705412561437558
-0.026377074993159316
0.02637707499315932
0.0036705412561437558
0.24370635792112405
-0.336582907023167
0.8536046924069197
julia> expv(1.0, A(10), b(10))
10-element Array{Float64,1}:
-0.8536046924069198
0.3365829070231669
-0.24370635792112402
-0.0036705412561437558
-0.026377074993159316
0.02637707499315932
0.0036705412561437558
0.24370635792112405
-0.336582907023167
0.8536046924069197
julia> a10 = A(10); b10 = b(10);
julia> @which expv(1.0, a10, b10)
expv(t, A, b) in ExponentialUtilities at /Users/forets/.julia/packages/ExponentialUtilities/tWv1t/src/krylov_phiv.jl:37
julia> eutil = expv(1.0, a10, b10)
10-element Array{Float64,1}:
-0.8536046924069198
0.3365829070231669
-0.24370635792112402
-0.0036705412561437558
-0.026377074993159316
0.02637707499315932
0.0036705412561437558
0.24370635792112405
-0.336582907023167
0.8536046924069197
julia> ekit = expmv(1.0, a10, b10)
10-element Array{Float64,1}:
-0.8536046924069195
0.3365829070231674
-0.2437063579211246
-0.003670541256143919
-0.02637707499315963
0.026377074993159617
0.003670541256143919
0.2437063579211246
-0.3365829070231674
0.8536046924069195
julia> norm(eutil - ekit)
1.2044079031002385e-15
julia> @btime expv(1.0, $a10, $b10)
13.987 μs (77 allocations: 8.22 KiB)
10-element Array{Float64,1}:
-0.8536046924069198
0.3365829070231669
-0.24370635792112402
-0.0036705412561437558
-0.026377074993159316
0.02637707499315932
0.0036705412561437558
0.24370635792112405
-0.336582907023167
0.8536046924069197
julia> @btime expmv(1.0, $a10, $b10)
27.846 μs (85 allocations: 14.30 KiB)
10-element Array{Float64,1}:
-0.8536046924069195
0.3365829070231674
-0.2437063579211246
-0.003670541256143919
-0.02637707499315963
0.026377074993159617
0.003670541256143919
0.2437063579211246
-0.3365829070231674
0.8536046924069195
julia> A100 = A(100); b100 = b(100);
julia> @btime expmv(1.0, $A100, $b100);
430.845 μs (104 allocations: 267.77 KiB)
julia> @btime expv(1.0, $A100, $b100);
135.513 μs (153 allocations: 65.23 KiB)
julia> 430e-6 / 135e-6
3.185185185185185
julia> @btime expv(1.0, $A100, $b100);^C
julia> A1000 = A(1000); b1000 = b(1000);
julia> @btime expmv(1.0, $A100, $b100);
405.285 μs (104 allocations: 267.77 KiB)
julia> @btime expv(1.0, $A100, $b100);
135.430 μs (153 allocations: 65.23 KiB)
julia> @btime expmv(1.0, $A1000, $b1000);
704.306 μs (71 allocations: 515.02 KiB)
julia> @btime expv(1.0, $A1000, $b1000);
448.999 μs (153 allocations: 304.42 KiB)
julia> A10000 = A(10000); b10000 = b(10000);
julia> @btime expmv(1.0, $A1000, $b1000);
666.230 μs (71 allocations: 515.02 KiB)
julia> @btime expmv(1.0, $A10000, $b10000);
5.308 ms (105 allocations: 2.84 MiB)
julia> @btime expv(1.0, $A10000, $b10000);
3.664 ms (156 allocations: 2.63 MiB)
julia> x = expmv(1.0, A10000, b10000); y = expv(1.0, A10000, b10000);
julia> norm(x-y)
1.0537327108923815e-13
julia> 5.3 / 3.66
1.448087431693989
Last active
November 12, 2018 22:26
-
-
Save mforets/d7d2f6afc29c1ea649661ac4ebd67b26 to your computer and use it in GitHub Desktop.
expmv and expv benchs
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment