Time Evolution and Quantum System Dynamics

Time Evolution Solutions

Solution

QuantumToolbox utilizes the powerful DifferentialEquation.jl to simulate different kinds of quantum system dynamics. Thus, we will first look at the data structure used for returning the solution (sol) from DifferentialEquation.jl. The solution stores all the crucial data needed for analyzing and plotting the results of a simulation. A generic structure TimeEvolutionSol contains the following properties for storing simulation data:

Fields (Attributes)Description
sol.timesThe time list of the evolution.
sol.statesThe list of result states.
sol.expectThe expectation values corresponding to each time point in sol.times.
sol.algThe algorithm which is used during the solving process.
sol.abstolThe absolute tolerance which is used during the solving process.
sol.reltolThe relative tolerance which is used during the solving process.
sol.retcode (or sol.converged)The returned status from the solver.

Accessing data in solutions

To understand how to access the data in solution, we will use an example as a guide, although we do not worry about the simulation details at this stage. The Schrödinger equation solver (sesolve) used in this example returns TimeEvolutionSol:

H = 0.5 * sigmax()
ψ0 = basis(2, 0)
e_ops = [
    proj(basis(2, 0)),
    proj(basis(2, 1)),
    basis(2, 0) * basis(2, 1)'
]
tlist = LinRange(0, 10, 100)

To see what is contained inside the solution, we can use the print function:

print(sol)
Solution of time evolution
(return code: Success)
--------------------------
num_states = 1
num_expect = 3
ODE alg.: OrdinaryDiffEqTsit5.Tsit5{typeof(OrdinaryDiffEqCore.trivial_limiter!), typeof(OrdinaryDiffEqCore.trivial_limiter!), Static.False}(OrdinaryDiffEqCore.trivial_limiter!, OrdinaryDiffEqCore.trivial_limiter!, static(false))
abstol = 1.0e-8
reltol = 1.0e-6

It tells us the number of expectation values are computed and the number of states are stored. Now we have all the information needed to analyze the simulation results. To access the data for the three expectation values, one can do:

expt1 = real(sol.expect[1,:])
expt2 = real(sol.expect[2,:])

Recall that Julia uses Fortran-style indexing that begins with one (i.e., [1,:] represents the 1-st observable, where : represents all values corresponding to tlist).

Together with the array of times at which these expectation values are calculated:

we can plot the resulting expectation values:

using CairoMakie
CairoMakie.enable_only_mime!(MIME"image/svg+xml"())

fig = Figure()
ax = Axis(fig[1, 1])
lines!(ax, times, expt1, label = L"P_00")
lines!(ax, times, expt2, label = L"P_11")
lines!(ax, times, expt3, label = L"P_01")

fig
Example block output

State vectors, or density matrices, are accessed in a similar manner:

sol.states
1-element Vector{QuantumObject{Vector{ComplexF64}, KetQuantumObject, 1}}:
 Quantum Object:   type=Ket   dims=[2]   size=(2,)
2-element Vector{ComplexF64}:
 0.28366218546240907 + 0.0im
                 0.0 + 0.9589242745990754im

Here, the solution contains only one (final) state. Because the states will be saved depend on the keyword argument saveat in kwargs. If e_ops is specified, the default value of saveat=[tlist[end]] (only save the final state), otherwise, saveat=tlist (saving the states corresponding to tlist). One can also specify e_ops and saveat separately.

Some other solvers can have other output.

Multiple trajectories solution

This part is still under construction, please visit API first.