How to run noisy circuits in `Yao.jl`

Specific Domains Quantum
1

I want to run simulations of noisy circuits (e.g. with a random Pauli error [I, X, Y, Z] with weights [1-3p, p, p, p] after each gate) in Yao.jl.

My first thought was to use full density matrix simulations, but it seems that there is no support for time evolution of DensityMatrixes.

My second thought was running a BatchedArrayReg through a UnitaryChannel, i.e. like this

# initialize the computational zero state on all states in the register
state = zeros(ComplexF32, 2^2, 3)
state[1,:] .= 1.
reg = BatchedArrayReg(state, 3)

# create a circuit corresponding to the depolarizing channel
noisechannel = UnitaryChannel([igate(1), X, Y, Z], [0.4, 0.2, 0.2, 0.2])
circuit = put(2, 1=>noisechannel)

# and apply the circuit to the register
apply!(reg, circuit)

But it seems that apply!() samples one of the gates in the noisechannel and then applies the same gate to all states in the batch. Is there a way, that a different gate gets applied to all states in the register?

And tangentially related: What is the use of the BatchedArrayRegister when all states in the batch evolve under exactly the same circuit?

2

We have not yet released the density matrix backend yet. Please check this pr: Density matrix implementation by GiggleLiu · Pull Request #394 · QuantumBFS/Yao.jl · GitHub . We can tag a new version at this weekend.

The batched register is for a batch of pure states, rather than the noisy circuit simulation.

3

Thanks for pointing me to that PR!

Out of curiosity: So the point of a batched register is to run multiple, different states through the same circuit, but not to run one state through a probabilistic circuit?

4

Yes, batch is for improving the performance of GPU simulation in certain applications, like some variational algorithms.

5

Hey, please check the newly released Yao@0.8. We have limited support to density matrix now.

julia> using Yao

julia> noisechannel = UnitaryChannel([put(2, 1=>igate(1)), put(2, 1=>X), put(2, 1=>Y), put(2, 1=>Z)], [0.4, 0.2, 0.2, 0.2])
nqubits: 2
unitary_channel
├─ [0.4] put on (1)
│  └─ igate(1)
├─ [0.2] put on (1)
│  └─ X
├─ [0.2] put on (1)
│  └─ Y
└─ [0.2] put on (1)
   └─ Z


julia> apply(density_matrix(zero_state(2)), noisechannel)
DensityMatrix{2, ComplexF64, Matrix{ComplexF64}}(ComplexF64[0.6000000000000001 + 0.0im 0.0 + 0.0im 0.0 + 0.0im 0.0 + 0.0im; 0.0 + 0.0im 0.4 + 0.0im 0.0 + 0.0im 0.0 + 0.0im; 0.0 + 0.0im 0.0 + 0.0im 0.0 + 0.0im 0.0 + 0.0im; 0.0 + 0.0im 0.0 + 0.0im 0.0 + 0.0im 0.0 + 0.0im])

julia> measure(r; nshots=3)
3-element Vector{DitStr{2, 2, Int64}}:
 01 ₍₂₎
 00 ₍₂₎
 01 ₍₂₎

julia> expect(EasyBuild.heisenberg(2), r)
0.40000000000000013 + 0.0im

julia> fidelity(r, r)
1.0

NOTE:
Since a UnitaryChannel block does not have a matrix representation, it is not composible with other blocks. It must be the outer most block. For example

julia> g = put(2, 1=>UnitaryChannel([igate(1), X, Y, Z], [0.4, 0.2, 0.2, 0.2]))
nqubits: 2
put on (1)
└─ unitary_channel
   ├─ [0.4] igate(1)
   ├─ [0.2] X
   ├─ [0.2] Y
   └─ [0.2] Z


julia> mat(g)
ERROR: `UnitaryChannel` does not have a matrix representation!
Stacktrace:
 [1] error(s::String)
   @ Base ./error.jl:35
 [2] mat(#unused#::Type{ComplexF64}, x::UnitaryChannel{2, Vector{Float64}})
   @ YaoBlocks ~/.julia/dev/Yao/lib/YaoBlocks/src/composite/unitary_channel.jl:64
 [3] mat(#unused#::Type{ComplexF64}, pb::PutBlock{2, 1, UnitaryChannel{2, Vector{Float64}}})
   @ YaoBlocks ~/.julia/dev/Yao/lib/YaoBlocks/src/composite/put_block.jl:86
 [4] mat(x::PutBlock{2, 1, UnitaryChannel{2, Vector{Float64}}})
   @ YaoBlocks ~/.julia/dev/Yao/lib/YaoBlocks/src/abstract_block.jl:120
 [5] top-level scope
   @ REPL[22]:1

Let me know if you want more features!

6

Nice! For now it seems this does pretty much everything I need.

I was also going to ask for gradient support, but then I realized that quantum backpropagation as implemented in Yao.jl only works because evolution with unitary gates is invertible. But evolution with a unitary channel is not trivially invertible, so I guess gradients are not that easily possible

7

That is true. We can add the AD rule to ChainRules if you really need it. Or if you already know how to differentiate, this is where you can contribute:

Just specify an rrule for the the channel apply function. Then you can use it in AD engines like Zygote.

(NOTE: Zygote’s memory management and stability can be an issue though.)