propagate with efficiency on w

efficiency_w.jl

function propagate(prob::QuKrylov, op, dt, current_qustate)
    basis_size = prob.options[:NC]
    N = min(basis_size, length(current_qustate))
    v = Array(typeof(current_qustate),0)
    push!(v,zeros(current_qustate))
    push!(v,current_qustate)
    # w = Array(typeof(current_qustate),0)
    w = Array(typeof(current_qustate),1)
    alpha = Array(Complex{Float64},0)
    beta = Array(Complex{Float64},0)
    push!(beta,0.)
    for i=2:N
        # push!(w, op*v[i])
        w[1] = op*v[i]
        # push!(alpha, w[i-1]'*v[i])
        push!(alpha, w[1]'*v[i])
        # w[i-1] = w[i-1]-alpha[i-1]*v[i]-beta[i-1]*v[i-1]
        w[1] = w[1]-alpha[i-1]*v[i]-beta[i-1]*v[i-1]
        # push!(beta, norm(w[i-1]))
        push!(beta, norm(w[1]))
        # push!(v, w[i-1]/beta[i])
        push!(v, w[1]/beta[i])
    end
    # push!(w, op*v[end])
    w[1] = op*v[end]
    # push!(alpha, w[end]'*v[end])
    push!(alpha, w[1]'*v[end])
    deleteat!(v,1)
    H_k = full(Tridiagonal(beta[2:end], alpha, beta[2:end]))
    U_k = expm(-im*dt*H_k)
    next_state = zeros(current_qustate)
    for j=1:N
        next_state = next_state + v[j]*U_k[j,1]
    end
    return next_state
end