Wrong evolution when using Newton Propagator?

[alastair-marshall]

Maybe I'm using the propagate function wrongly but here's a not very minimal working example:

using QuantumPropagators
using LinearAlgebra
ρ0 = [1 0;0 0] .+ 0.0im
ρ0vec = reshape(ρ0, (4))
const sx = [0 1;1 0] .+ 0.0im
super_sx = kron(eye2, sx) - kron(sx', eye2)

T = 1.0
dt = 0.1
N = floor(Int, T/dt)
K = 1 # number of controls
drive = rand(K, N) .* 0.0 .+ 1 * pi

function get_genfunc(drive, sx)
    gen_func = (tlist, i; kwargs...) -> drive[i] * sx 
end

gen_func = get_genfunc(drive, super_sx)
# should be at the half points but for now just at one end
tlist = collect(0:dt:T-dt)
wrk = initpropwrk(ρ0vec, tlist, nothing; method=:newton)

test = propagate(ρ0vec, gen_func, tlist; wrk = wrk)
test ≈ ρ0vec # false

# working version using ExponentialUtilities
using ExponentialUtilities

function test_prop(state, gen_func, tlist)
    ev = copy(state)

    for i = 1:length(tlist)
        H = gen_func(tlist, i)
        ev .= expv(-1.0im * dt, H, ev)
    end
    ev
end

test2 = test_prop(ρ0vec, gen_func, tlist)
test2 ≈ ρ0vec # true

[goerz]

Maybe I misunderstand the code, but isn’t test the propagated state from the original ρ0vec for Newton? The propagate routine is not in-place

[goerz]

Apart from that I think propagate might be doing one less time step than test_prop. It make one step per interval, not one step per time grid point

[alastair-marshall]

Yes it is but the evolution should just be a Rabi period so they should end up the same, it's definitely not a good example

[alastair-marshall]

Ah okay I can quickly check that

[alastair-marshall]

Apart from that I think propagate might be doing one less time step than test_prop. It make one step per interval, not one step per time grid point

Yup thats it