Quantum Computing Stack Exchange is a question and answer site for engineers, scientists, programmers, and computing professionals interested in quantum computing. It only takes a minute to sign up.
Sign up to join this community
How can I get a subspace of a quantum circuit? More precisely, I'm dealing with quantum circuit with data qubits ('q') and ancilla qubits ('anc'), such as $(q_0,q_1,...,q_n,anc_0,..anc_m)$.
After some quantum operations, I'd like to get a unitary matrix representation of the data qubits. Is this possible?
What you are looking for is the partial trace of a density matrix, $|\psi\rangle\langle\psi|$, where $|\psi\rangle=|q_0,q_1,...,ac_0,ac_1,...\rangle$. This matrix is surely a unitary matrix, and it can be viewed as a tensor product of all the density matrices of each qubit.
To get a subspace of a quantum circuit, you have to partial trace over all the irrelevant qubits, in your case, $Tr_{ac}(|\psi\rangle\langle\psi|)$.
If you are using qiskit, see the lower code as an example:
from qiskit import QuantumCircuit,QuantumRegister
from qiskit.quantum_info import DensityMatrix,partial_trace
import numpy as np
qr=QuantumRegister(2)
circ=QuantumCircuit(qr)
circ.h(qr[0])
circ.cx(qr[0],qr[1])
DM=DensityMatrix.from_instruction(circ)
print(DM.data)
PT=partial_trace(DM,[0])
print(PT.data)
If you are not using the python package, and suppose you already got the density matrix of the entire quantum circuit, see the lower code(qiskit is again used to generate the density matrix, but I provide another function for partial trace in this place):
def partial_trace(rho,qubit2keep):
num_qubit=int(np.log2(rho.shape[0]))
for i in range(len(qubit2keep)):
qubit2keep[i]=num_qubit-1-qubit2keep[i]
qubit_axis=[(i,num_qubit+i) for i in range(num_qubit)
if i not in qubit2keep]
minus_factor=[(i,2*i) for i in range(len(qubit_axis))]
minus_qubit_axis=[(q[0]-m[0],q[1]-m[1])
for q, m in zip(qubit_axis,minus_factor)]
rho_res=np.reshape(rho,[2,2]*num_qubit)
qubit_left=num_qubit-len(qubit_axis)
for i,j in minus_qubit_axis:
rho_res=np.trace(rho_res,axis1=i,axis2=j)
if qubit_left>1:
rho_res=np.reshape(rho_res,[2**qubit_left]*2)
return rho_res
from qiskit import QuantumCircuit,QuantumRegister
from qiskit.quantum_info import DensityMatrix
import numpy as np
qr=QuantumRegister(2)
circ=QuantumCircuit(qr)
circ.h(qr[0])
circ.cx(qr[0],qr[1])
DM=DensityMatrix.from_instruction(circ)
print(DM.data)
PT=partial_trace(DM.data,[0])
print(PT)
(Note that the partial_trace function is a rewritten of this one, the default sequence in this code and mine are reversed.)
The two codes give identical results. In PT=partial_trace(DM.data,[0]), [0] denotes
the qubit to be traced out, maybe some test should be helpful for you to trace out the qubits that you don't need correctly.
