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.