qojulia / QuantumOptics.jl

Library for the numerical simulation of closed as well as open quantum systems.
http://qojulia.org
Other
543 stars 110 forks source link

Support for measurements #296

Open seadra opened 3 years ago

seadra commented 3 years ago

It would be great to have POVM measurements, as well as measurements with post-selection. Both are available in QuantumInformation.jl, Yao.jl and recently added to QuTip (git only) as well.

david-pl commented 3 years ago

That sounds like an interesting feature. The issue is that I'm not too knowledgeable when it comes to POVM or PVM measurements. Could you provide some more details?

seadra commented 3 years ago

Thanks! POVM measurements are projective measurements with respect to a given set of operators M_i which satisfy the completeness relation \sum_i M_i^\dagger M_i = I (Section 2.2.6 in Nielsen & Chuang). In QuTip, this is implemented in measurement.py as _measurement_statistics_povm_dm and returns the array of probabilities tr(M_i^\daggerM_i'ρ), Eq. 2.143 in Nielsen & Chuang (for kets, it is def _measurement_statistics_povm_ket).

Measurement with post-selection returns the state of the system after a specified measurement M_i (Eqs. 2.131 and 2.147 in Nielsen and Chuang), essentially performing a wave-function collapse using Born rule. In QuTip, it's implemented in measure_observable (although it does a little more than just calculating the collapsed state, by diagonalization and finding eigenvalues to additionally report the result of the measurement as well). For simulating quantum circuits which include measurement boxes anywhere along the circuit, this would be the necessary function to have.

david-pl commented 3 years ago

You say that the measurement statistics should return the probabilities, but the functions you reference in QuTiP return a set of projected states as well. Also, given a vector M that contains the POVM operators and a state, the probabilities are literally just

probs = [expect(m'*m, state) for m=M]

for both Ket and Operator states. So that doesn't really need to be implemented. Returning the states as well is not hard, but probably worth to just have implemented. I suppose that is what you want on top of the probabilities. I'd do something like this:

function povm_statistics(state, ops::Vector)
    ps = zeros(length(ops))
    projected_states = Vector{typeof(state)}(undef, length(ops))
    for i=1:length(ops)
        p, state_projected = _projected_state(state, ops[i])
        ps[i] = p
        projected_states[i] = state_projected
    end
    return ps, projected_states
end
function _projected_state(state::Ket, op)
    state_proj = copy(state)
    QuantumOpticsBase.mul!(state_proj, op, state)
    p = dot(state_proj.data, state_proj.data)
    rmul!(state_proj.data, 1/sqrt(p))
    return p, state_proj
end
function _projected_state(state::Operator, op)
    state_proj = copy(state)
    tmp = copy(state)
    QuantumOpticsBase.mul!(tmp, op', state)
    QuantumOpticsBase.mul!(state_proj, tmp, op)
    p = tr(state_proj)
    rmul!(state_proj.data, 1/p)
    return p, state_proj
end

Measuring an observable can be implemented using this as well. However, it would require building the POVM by actually writing down a bunch of matrices that are just projectors. You can actually do it without building up those matrices, which is beneficial when working with many qubits:

function measure_observable(state, op; kwargs...)
    d, v = eigen(op.data; kwargs...)
    probs = zeros(length(d))
    for i=1:length(d)
        probs[i] = _prob(state, v[:,i])
    end
    cumprobs = cumsum(probs./sum(probs))
    i = findfirst(cumprobs.>rand())
    state_projected = _make_state(state, v[:,i])
    return d[i], state_projected
end
_prob(state::Ket, v) = abs2(dot(v, state.data))
_prob(state::Operator, v) = real(dot(v, state.data, v))
_make_state(state::Ket, v) = Ket(state.basis, v)
_make_state(state::Operator, v) = Operator(state.basis_l, state.basis_r, v*v')

Is that what you need? Could you test this in your use-case?

seadra commented 3 years ago

POVMMeasurement in QuantumInformation.jl actually only returns a diagonal matrix of probabilities. Yes, it's a one liner (assuming you don't verify the given set of operators is complete and skip other sanity checks), just like Pauli matrices, but IMHO it's one of the fundamental concepts in QIP that probably deserves to have its function.

For my purposes, I'm currently using QuTip. but thanks anyway for the code. I'll be looking forward to use QuantumOptics.jl in the future. I just wanted to point out that it would be nice to have measurements out-of-the-box in QuantumOptics.jl since that's an essential and very practical concept in quantum mechanics.

seadra commented 3 years ago

For reference, measurement API is now included in the recently released version of QuTiP 4.6 (changelog).