Partial trace for sequential measurement simulation

I’m trying to simulate a sequential circuit process that uses PNR detection, a BS gate, and a displacement gate at each step. After each step I trace out the detected mode (mode 1) and use the state in the non-detected mode (mode 2) as the input state in mode 1 for the next step and repeat. For instance at the start of the sequence:

eng = sf.Engine("fock", backend_options={"cutoff_dim": 15})
prog = sf.Program(2) # 2 mode circuit
with prog.context as q:
     Dgate(0.5) | q[1]  # displace vacuum
     BSgate(np.pi/6, 0) | (q[0], q[1]) 
     MeasureFock() | q[0]
result = eng.run(prog)

I then take the output of mode 2 by doing:

dm = result.state.reduced_dm([1])

and then perform the next step in the sequence as so:

prog = sf.Program(2) # 2 mode circuit
with prog.context as q:
     Dgate(1) | q[0]  # displace vacuum
     DensityMatrix(dm) | q[1]
     BSgate(np.pi/3, 0) | (q[0], q[1]) 
     MeasureFock() | q[0]
result = eng.run(prog)

My question is: is using DensityMatrix(dm) | q[0] the correct way to go about doing this or is there a more proper way to accomplish this in strawberry fields?

Hey @Amanuel!

It seems fine to me :slight_smile:. Let’s double check with @Sebastian_Duque_Mesa!

1 Like

Hi @Amanuel — As @isaacdevlugt said, this works!

For completeness: you can always concatenate programs, i.e., continue one program with another. For example, your code can be rewritten as

eng = sf.Engine("fock", backend_options={"cutoff_dim": 15})

prog1 = sf.Program(2) # 2 mode circuit
with prog1.context as q:
     Dgate(0.5) | q[1] 
     BSgate(np.pi/6, 0) | (q[0], q[1]) 
     MeasureFock() | q[0] # measurements reset quantum registers to vacuum!

result = eng.run(prog)
dm1 = result.state.reduced_dm([1])

# continue previous program by initializing 
# a new program with the previous one
prog2 = sf.Program(prog1)
with prog2.context as q:
     Dgate(1) | q[0]
     BSgate(np.pi/3, 0) | (q[0], q[1]) 
     MeasureFock() | q[0]

result = eng.run(prog)
dm2 = result.state.reduced_dm([1])

Remember that measurements reset quantum registers. If there is some structure in your circuit architecture, you can put all this in a for loop.

1 Like

Thanks for the reply @Sebastian_Duque_Mesa,

I forgot to ask this in the original post as well but I wanted to know if strawberry fields also allows one to calculate the probability of detection events when performing a measurement?

For example if one does

with prog.context as q:
     MeasureHomodyne(phi=0) | q[0]

I know we can get the resulting x value with

result = eng.run(prog)
x = result.samples[0][0] # homodyne sample

However is there a way to also receive the probability of measuring this particular sample as well? I can see a way of doing this by using something like x_quad_values but I was wondering if there was a simpler way.

Hey @Amanuel,

There is currently no method other than using x_quad_values to compute the probabilities of a sample — or at least not that I am aware of.

1 Like