SE, LVN, LME types and integration

prop_eq_types_2.jl

using QuBase

import Base: start, next, done

abstract QuPropagatorMethod

abstract QuantumEquation

immutable SE{H<:QuBase.AbstractQuMatrix} <: QuantumEquation
    hamiltonian::H
    SE(hamiltonian) = new(hamiltonian)
end

SE{H<:QuBase.AbstractQuMatrix}(hamiltonian::H) = SE{H}(hamiltonian)

immutable LVN{H<:QuBase.AbstractQuMatrix} <: QuantumEquation
    hamiltonian::H
    LVN(hamiltonian) = new(hamiltonian)
end

LVN{H<:QuBase.AbstractQuMatrix}(hamiltonian::H) = LVN{H}(hamiltonian)

immutable QuPropagator{QPM<:QuPropagatorMethod, V<:Union(QuBase.AbstractQuVector,QuBase.AbstractQuMatrix)}
    eq
    init_state::V
    tlist
    method::QPM
    QuPropagator(eq, init_state, tlist, method) = new(eq, init_state, tlist, method)
end

QuPropagator{QPM<:QuPropagatorMethod, V<:QuBase.AbstractQuVector}(eq::SE, init_state::V, tlist, method::QPM) = QuPropagator{QPM,V}(eq, init_state, tlist, method)

QuPropagator{QPM<:QuPropagatorMethod, V<:QuBase.AbstractQuVector}(param::QuBase.AbstractQuMatrix, init_state::V, tlist, method::QPM) = QuPropagator{QPM,V}(SE(param), init_state, tlist, method)

immutable QuPropagatorState
    psi
    t
    t_state
end

function Base.start(prob::QuPropagator)
    init_state = prob.init_state
    t_state = start(prob.tlist)
    t,t_state = next(prob.tlist,t_state)
    return QuPropagatorState(init_state,t,t_state)
end

function op(qu_eq::LVN)
    return qu_eq.hamiltonian
end

function op(qu_eq::SE)
    return qu_eq.hamiltonian
end

function Base.next{QPM<:QuPropagatorMethod}(prob::QuPropagator{QPM}, qustate::QuPropagatorState)
    current_qustate = qustate.psi
    current_t = qustate.t
    t,t_state = next(prob.tlist, qustate.t_state)
    next_qustate = propagate(prob.method, t, current_t, current_qustate, prob.eq)
    return (t, next_qustate), QuPropagatorState(next_qustate, t, t_state)
end

Base.done(prob::QuPropagator, qustate::QuPropagatorState) = done(prob.tlist, qustate.t_state)

immutable QuEuler <: QuPropagatorMethod
end

function propagate(prob::QuEuler, t, current_t, current_qustate, qu_eq::QuantumEquation)
    dt = t - current_t
    return (eye(op(qu_eq))-im*op(qu_eq)*dt)*current_qustate
end

qv = normalize!(QuArray([0.5+0.1im, 0.2+0.2im]))
t = 0.:0.1:2*pi

se = SE(sigmax)
qv = normalize!(QuArray([0.5+0.1im, 0.2+0.2im]))
lvn = LVN(sigmax)
qe = QuPropagator(se, qv, t, QuEuler())
qe_m = QuPropagator(sigmax, qv, t, QuEuler())