DifferentialEquations solvers
In this page, we will benchmark several solvers provided by DifferentialEquations.jl for solving time evolution in hierarchical equations of motion approach.
using OrdinaryDiffEq ## or "using DifferentialEquations"
using BenchmarkTools
using HierarchicalEOM
HierarchicalEOM.versioninfo()
__
/ \
__ __ __ \__/ __
| | | | / \ / \
| | | | ______ ______ \__/_ _\__/
| |___| |/ __ \ / __ \ / ' \/ \
| ___ | |__) | / \ | _ _ |
| | | | ____/| ( ) | / \ / \ |
| | | | |____ | \__/ | | | | | |
|__| |__|\______) \______/|__| |_| |_|
Julia framework for Hierarchical Equations of Motion
≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡
Copyright © NCKU-QFORT 2022 and later.
Lead developer : Yi-Te Huang
Other developers:
Simon Cross, Neill Lambert, Po-Chen Kuo and Shen-Liang Yang
Package information:
====================================
HierarchicalEOM Ver. 1.6.0
LinearSolve Ver. 2.22.1
OrdinaryDiffEq Ver. 6.80.1
FastExpm Ver. 1.1.0
JLD2 Ver. 0.4.48
System information:
====================================
Julia Version: 1.10.4
OS : Linux (x86_64-linux-gnu)
CPU : 4 × AMD EPYC 7763 64-Core Processor
Memory : 15.606 GB
WORD_SIZE: 64
LIBM : libopenlibm
LLVM : libLLVM-15.0.7 (ORCJIT, znver3)
BLAS : libopenblas64_.so (ilp64)Here, we use the example of driven systems and dynamical decoupling:
Γ = 0.0005
W = 0.005
kT = 0.05
N = 3
tier = 6
amp = 0.50
delay = 20
tlist = 0:0.4:400
σz = [1 0; 0 -1]
σx = [0 1; 1 0]
H0 = 0.0 * σz
ρ0 = 0.5 * [1 1; 1 1];
param = (amp, delay, σx)
function pulse(V, Δ, t)
τ = 0.5 * π / V
period = τ + Δ
if (t % period) < τ
return V
else
return 0
end
end
function H_D(p::Tuple, t)
V, Δ, σx = p
return pulse(V, Δ, t) * σx
end
bath = Boson_DrudeLorentz_Pade(σz, Γ, W, kT, N)
M = M_Boson(H0, tier, bath);Preparing block matrices for HEOM Liouvillian superoperator (using 1 threads)...
Processing: 1%|▎ | ETA: 0:00:21 (99.26 ms/it)
Processing: 100%|████████████████████| Time: 0:00:00 ( 0.96 ms/it)
Constructing matrix...[DONE]ODE Solver List
(click here to see the full solver list provided by DifferentialEquations.jl)
For any extra solver options, we can add it in the function evolution with keyword arguments. These keyword arguments will be directly pass to the solvers in DifferentialEquations (click here to see the documentation for the common solver options)
Furthermore, since we are using BenchmarkTools (@benchmark) in the following, we set verbose = false to disable the output message.
DP5 (Default solver)
Dormand-Prince's 5/4 Runge-Kutta method. (free 4th order interpolant)
@benchmark evolution(M, ρ0, tlist, H_D, param; abstol = 1e-12, reltol = 1e-12, verbose = false)BenchmarkTools.Trial: 6 samples with 1 evaluation.
Range (min … max): 770.823 ms … 1.097 s ┊ GC (min … max): 2.05% … 25.52%
Time (median): 798.421 ms ┊ GC (median): 2.09%
Time (mean ± σ): 849.998 ms ± 124.155 ms ┊ GC (mean ± σ): 7.11% ± 9.57%
█ █ ██ █ █
█▁█▁██▁▁▁▁▁▁▁▁▁█▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁█ ▁
771 ms Histogram: frequency by time 1.1 s <
Memory estimate: 763.00 MiB, allocs estimate: 1276588.RK4
The canonical Runge-Kutta Order 4 method. Uses a defect control for adaptive stepping using maximum error over the whole interval.
@benchmark evolution(M, ρ0, tlist, H_D, param; solver = RK4(), abstol = 1e-12, reltol = 1e-12, verbose = false)BenchmarkTools.Trial: 3 samples with 1 evaluation.
Range (min … max): 1.775 s … 1.793 s ┊ GC (min … max): 1.87% … 1.92%
Time (median): 1.781 s ┊ GC (median): 1.91%
Time (mean ± σ): 1.783 s ± 9.557 ms ┊ GC (mean ± σ): 1.90% ± 0.03%
█ █ █
█▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁█▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁█ ▁
1.77 s Histogram: frequency by time 1.79 s <
Memory estimate: 1.60 GiB, allocs estimate: 2691713.Tsit5
Tsitouras 5/4 Runge-Kutta method. (free 4th order interpolant).
@benchmark evolution(M, ρ0, tlist, H_D, param; solver = Tsit5(), abstol = 1e-12, reltol = 1e-12, verbose = false)BenchmarkTools.Trial: 5 samples with 1 evaluation.
Range (min … max): 996.222 ms … 1.245 s ┊ GC (min … max): 2.26% … 24.10%
Time (median): 1.022 s ┊ GC (median): 2.20%
Time (mean ± σ): 1.061 s ± 103.061 ms ┊ GC (mean ± σ): 7.36% ± 9.78%
█
▇▁▁▁▁▁█▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▇ ▁
996 ms Histogram: frequency by time 1.24 s <
Memory estimate: 939.73 MiB, allocs estimate: 1604252.Vern7
Verner's “Most Efficient” 7/6 Runge-Kutta method. (lazy 7th order interpolant).
@benchmark evolution(M, ρ0, tlist, H_D, param; solver = Vern7(), abstol = 1e-12, reltol = 1e-12, verbose = false)BenchmarkTools.Trial: 4 samples with 1 evaluation.
Range (min … max): 1.331 s … 1.376 s ┊ GC (min … max): 2.59% … 2.56%
Time (median): 1.345 s ┊ GC (median): 2.60%
Time (mean ± σ): 1.349 s ± 19.585 ms ┊ GC (mean ± σ): 2.62% ± 0.18%
█ █ █ █
█▁▁▁▁▁▁▁▁▁▁▁█▁▁▁▁▁▁▁▁▁█▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁█ ▁
1.33 s Histogram: frequency by time 1.38 s <
Memory estimate: 1.25 GiB, allocs estimate: 2211323.Vern9
Verner's “Most Efficient” 9/8 Runge-Kutta method. (lazy 9th order interpolant)
@benchmark evolution(M, ρ0, tlist, H_D, param; solver = Vern9(), abstol = 1e-12, reltol = 1e-12, verbose = false)BenchmarkTools.Trial: 3 samples with 1 evaluation.
Range (min … max): 2.057 s … 2.087 s ┊ GC (min … max): 2.54% … 2.62%
Time (median): 2.069 s ┊ GC (median): 2.64%
Time (mean ± σ): 2.071 s ± 15.168 ms ┊ GC (mean ± σ): 2.64% ± 0.11%
█ █ █
█▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁█▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁█ ▁
2.06 s Histogram: frequency by time 2.09 s <
Memory estimate: 1.89 GiB, allocs estimate: 3363360.This page was generated using Literate.jl.