Merge pull request #12 from ChrisRackauckas/update_bench

Update DifferentialEquations.jl benchmarks
Nicholaswogan · Dec 13, 2021 · 6b37966 · 6b37966
2 parents a389ae8 + b05b7b2
commit 6b37966
Show file tree

Hide file tree

Showing 2 changed files with 39 additions and 23 deletions.
diff --git a/benchmark/README.md b/benchmark/README.md
@@ -14,8 +14,8 @@ Also, Scipy might be faster than NumbaLSODA for stiff problems with > ~500 ODEs
 | Scipy LSODA  | 595.029 ms  | 174x                   |
 | Scipy RK45   | 1744.073 ms | 510x                   |
 | Scipy DOP853 | 1052.442 ms | 308x                   |
-| Julia Tsit5  | 2.612 ms    | 0.76x                  |
-| Julia Vern8  | 1.803 ms    | 0.53x                  |
+| Julia Tsit5  | 1.340 ms    | 0.39x                  |
+| Julia Vern8  | 1.803 ms    | 0.21x                  |
 
 ## Stiff (Rober)
 
@@ -25,9 +25,9 @@ Also, Scipy might be faster than NumbaLSODA for stiff problems with > ~500 ODEs
 | Scipy LSODA        | 26.483 ms  | 120x                   |
 | Scipy BDF          | 158.911 ms | 719x                   |
 | Scipy Radau        | 167.977 ms | 760x                   |
-| Julia TRBDF2       | 5.765 ms   | 26.08x                 |
-| Julia CVODE_BDF    | 1.115 ms   | 5.045x                 |
-| Julia LSODA        | 0.264 ms   | 1.194x                 |
-| Julia Rodas5       | 0.762 ms   | 3.45x                  |
-| Julia Rodas5 + StaticArrays | 0.173 ms   | 0.78x                 |
-| Julia Rodas5 + StaticArrays + Analytical Jacobian       | 0.113 ms   | 0.511x                 |
+| Julia TRBDF2       | 2.202 ms   | 9.96x                 |
+| Julia CVODE_BDF    | 0.628 ms   | 2.84x                 |
+| Julia LSODA        | 0.254 ms   | 1.15x                 |
+| Julia Rodas5       | 0.227 ms   | 1.03x                  |
+| Julia Rodas5 + StaticArrays | 0.0423 ms   | 0.19x                 |
+| Julia Rodas5 + StaticArrays + Analytical Jacobian       | 0.0420 ms ms   | 0.19x                 |
diff --git a/benchmark/benchmark.jl b/benchmark/benchmark.jl
@@ -1,5 +1,4 @@
-
-using OrdinaryDiffEq, StaticArrays, BenchmarkTools, Sundials, LSODA, ModelingToolkit
+using OrdinaryDiffEq, StaticArrays, BenchmarkTools, Sundials, LSODA
 
 println("Lorenz")
 # Lorenz
@@ -15,27 +14,27 @@ u0 = SA[1.0,0.0,0.0]
 p  = SA[10.0,28.0,8/3]
 tspan = (0.0,100.0)
 prob = ODEProblem(lorenz_static,u0,tspan,p)
-@btime solve(prob,Tsit5(),saveat=0.1,reltol=1.0e-8,abstol=1.0e-8) # 2.612 ms (30 allocations: 59.22 KiB)
-@btime solve(prob,Vern8(),saveat=0.1,reltol=1.0e-8,abstol=1.0e-8); # 1.803 ms (31 allocations: 65.70 KiB)
+@btime solve(prob,Tsit5(),saveat=0.1,reltol=1.0e-8,abstol=1.0e-8) # 1.340 ms (27 allocations: 64.55 KiB)
+@btime solve(prob,Vern8(),saveat=0.1,reltol=1.0e-8,abstol=1.0e-8); # 718.000 μs (28 allocations: 71.52 KiB)
 
 println("\nRober")
 # rober
 function rober(du,u,p,t)
     y₁,y₂,y₃ = u
     k₁,k₂,k₃ = p
     @inbounds begin
-    du[1] = -k₁*y₁+k₃*y₂*y₃
-    du[2] =  k₁*y₁-k₂*y₂^2-k₃*y₂*y₃
-    du[3] =  k₂*y₂^2
+        du[1] = -k₁*y₁+k₃*y₂*y₃
+        du[2] =  k₁*y₁-k₂*y₂^2-k₃*y₂*y₃
+        du[3] =  k₂*y₂^2
     end
     nothing
 end
-prob = ODEProblem(rober,[1.0,0.0,0.0],(0.0,1e5),[0.04,3e7,1e4]) 
+prob = ODEProblem(rober,[1.0,0.0,0.0],(0.0,1e5),[0.04,3e7,1e4])
 
-@btime solve(prob,Rodas5(),reltol=1.0e-8,abstol=1.0e-8, saveat = 1000) # 761.578 μs (633 allocations: 53.05 KiB)
-@btime solve(prob,TRBDF2(),reltol=1.0e-8,abstol=1.0e-8, saveat = 1000) # 5.765 ms (11415 allocations: 524.14 KiB)
-@btime solve(prob,CVODE_BDF(),reltol=1.0e-8,abstol=1.0e-8, saveat = 1000) # 1.115 ms (6954 allocations: 295.83 KiB)
-@btime solve(prob,lsoda(),reltol=1.0e-8,abstol=1.0e-8, saveat = 1000) # 263.684 μs (2169 allocations: 185.09 KiB)
+@btime solve(prob,Rodas5(),reltol=1.0e-8,abstol=1.0e-8, saveat = 1000) # 227.200 μs (1288 allocations: 52.56 KiB)
+@btime solve(prob,TRBDF2(),reltol=1.0e-8,abstol=1.0e-8, saveat = 1000) # 2.202 ms (3424 allocations: 194.41 KiB)
+@btime solve(prob,CVODE_BDF(),reltol=1.0e-8,abstol=1.0e-8, saveat = 1000) # 627.500 μs (7183 allocations: 214.33 KiB)
+@btime solve(prob,lsoda(),reltol=1.0e-8,abstol=1.0e-8, saveat = 1000) # 254.400 μs (2370 allocations: 127.09 KiB)
 
 # rober
 function rober_static(u,p,t)
@@ -47,8 +46,25 @@ function rober_static(u,p,t)
     SA[du1,du2,du3]
 end
 prob = ODEProblem{false}(rober_static,SA[1.0,0.0,0.0],(0.0,1e5),SA[0.04,3e7,1e4])
-@btime solve(prob,Rodas5(),reltol=1.0e-8,abstol=1.0e-8, saveat = 1000) # 172.953 μs (580 allocations: 53.36 KiB)
+@btime solve(prob,Rodas5(),reltol=1.0e-8,abstol=1.0e-8, saveat = 1000) # 42.300 μs (29 allocations: 9.75 KiB)
 
-prob2 = ODEProblem{false}(modelingtoolkitize(prob),SA[1.0,0.0,0.0],(0.0,1e5),SA[0.04,3e7,1e4],jac=true)
-@btime solve(prob2,Rodas5(), reltol=1.0e-8, abstol=1.0e-8, saveat = 1000); # 113.475 μs (400 allocations: 51.19 KiB)
+function rober_jac(u,p,t)
+  y₁,y₂,y₃ = u
+  k₁,k₂,k₃ = p
+  J11 = k₁ * -1
+  J21 = k₁
+  J31 = 0
+  J12 = y₃ * k₃
+  J22 = y₂ * k₂ * -2 + y₃ * k₃ * -1
+  J32 = y₂ * 2 * k₂
+  J13 = k₃ * y₂
+  J23 = k₃ * y₂ * -1
+  J33 = 0
+  SA[J11 J12 J13
+     J21 J22 J23
+     J31 J32 J33]
+end
 
+ff = ODEFunction(rober_static, jac=rober_jac)
+prob2 = ODEProblem{false}(ff,SA[1.0,0.0,0.0],(0.0,1e5),SA[0.04,3e7,1e4])
+@btime solve(prob2,Rodas5(), reltol=1.0e-8, abstol=1.0e-8, saveat = 1000); # 42.000 μs (30 allocations: 9.70 KiB)