Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

Implement the OrdinaryDiffEq interface for Kets #16

Closed
wants to merge 2 commits into from
Closed

Implement the OrdinaryDiffEq interface for Kets #16

wants to merge 2 commits into from

Commits on Apr 7, 2021

  1. Implement the OrdinaryDiffEq interface for Kets 

    The following now works

```julia
using QuantumOptics
using DifferentialEquations

ℋ = SpinBasis(1//2)

σx = sigmax(ℋ)

↓ = s =  spindown(ℋ)

schrod(ψ,p,t) = im * σx * ψ

t₀, t₁ = (0.0, pi)
Δt = 0.1

prob = ODEProblem(schrod, ↓, (t₀, t₁))
sol = solve(prob,Tsit5())
```

It works for Bras as well.
It works for in-place operations and in some situations it is
faster than the standard `timeevolution.schroedinger`.

```julia
ℋ = SpinBasis(20//1)
↓ = spindown(ℋ)
t₀, t₁ = (0.0, pi)
const σx = sigmax(ℋ)
const iσx = im * σx
schrod!(dψ,ψ,p,t) = mul!(dψ, iσx, ψ)
prob! = ODEProblem(schrod!, ↓, (t₀, t₁))

julia> @benchmark sol = solve($prob!,DP5(),save_everystep=false)
BenchmarkTools.Trial:
  memory estimate:  22.67 KiB
  allocs estimate:  178
  --------------
  minimum time:     374.463 μs (0.00% GC)
  median time:      397.327 μs (0.00% GC)
  mean time:        406.738 μs (0.37% GC)
  maximum time:     4.386 ms (89.76% GC)
  --------------
  samples:          10000
  evals/sample:     1

julia> @benchmark timeevolution.schroedinger([$t₀,$t₁], $↓, $σx)
BenchmarkTools.Trial:
  memory estimate:  23.34 KiB
  allocs estimate:  161
  --------------
  minimum time:     748.106 μs (0.00% GC)
  median time:      774.601 μs (0.00% GC)
  mean time:        786.933 μs (0.14% GC)
  maximum time:     4.459 ms (80.46% GC)
  --------------
  samples:          6350
  evals/sample:     1
```
    @Krastanov
    Krastanov committed Apr 7, 2021
    Configuration menu
    Copy the full SHA
    e08608b View commit details
    Browse the repository at this point in the history

  2. Implement OrdinaryDiffEq interface for dense operators. 

    The following works now:

```julia
ℋ = SpinBasis(20//1)

const σx = sigmax(ℋ)
const iσx = im * σx
const σ₋ = sigmam(ℋ)
const σ₊ = σ₋'
const mhalfσ₊σ₋ = -σ₊*σ₋/2

↓ = spindown(ℋ)
ρ = dm(↓)

lind(ρ,p,t) = - iσx * ρ + ρ * iσx + σ₋*ρ*σ₊ + mhalfσ₊σ₋ * ρ + ρ * mhalfσ₊σ₋

t₀, t₁ = (0.0, pi)
Δt = 0.1

prob = ODEProblem(lind, ρ, (t₀, t₁))
sol = solve(prob,Tsit5())
```

Works in-place as well.

It is slightly slower than `timeevolution.master`:

```julia
function makelind!()
    tmp = zero(ρ) # this is the global rho
    function lind!(dρ,ρ,p,t) # TODO this can be much better with a good Tullio kernel
        mul!(tmp, ρ, σ₊)
        mul!(dρ, σ₋, ρ)
        mul!(dρ,    ρ, mhalfσ₊σ₋, true, true)
        mul!(dρ, mhalfσ₊σ₋,    ρ, true, true)
        mul!(dρ,  iσx,    ρ, -ComplexF64(1),   ComplexF64(1))
        mul!(dρ,    ρ,  iσx,  true,   true)
        return dρ
    end
end
lind! = makelind!()
prob! = ODEProblem(lind!, ρ, (t₀, t₁))

julia> @benchmark sol = solve($prob!,DP5(),save_everystep=false)
BenchmarkTools.Trial:
  memory estimate:  408.94 KiB
  allocs estimate:  213
  --------------
  minimum time:     126.334 ms (0.00% GC)
  median time:      127.359 ms (0.00% GC)
  mean time:        127.876 ms (0.00% GC)
  maximum time:     138.660 ms (0.00% GC)
  --------------
  samples:          40
  evals/sample:     1

julia> @benchmark timeevolution.master([$t₀,$t₁], $ρ, $σx, [$σ₋])
BenchmarkTools.Trial:
  memory estimate:  497.91 KiB
  allocs estimate:  210
  --------------
  minimum time:     97.902 ms (0.00% GC)
  median time:      98.469 ms (0.00% GC)
  mean time:        98.655 ms (0.00% GC)
  maximum time:     104.850 ms (0.00% GC)
  --------------
  samples:          51
  evals/sample:     1
```
    @Krastanov
    Krastanov committed Apr 7, 2021
    Configuration menu
    Copy the full SHA
    4ca6500 View commit details
    Browse the repository at this point in the history