# How to check if a statevector is made of product of EPR pairs in Qiskit?

I have a statevector of 10 qubits. I want to check if the qubits 0,1 and 2,3 are EPR pairs/Bell pairs. The other qubits (4,5,6,7,8,9) have been measured.

There are four possible 2-qubit states that are Bell pairs. I created circuits for all possible tensor products of these 4 states with the following code (total 16 states):

num_qubits = 10

#measurement_string contains the measurement values of qubits 4,5,6,7,8,9
#example:measurement_string  = 001010 means q4 was measured = 0,q5 = 0, q6 = 1, q7 = 0,q8 = 1, q9 = 0

def epr_product_1(measurement_string):
#first entangled pair
qc = QuantumCircuit(num_qubits)
qc.h(0)
qc.cx(0,1)

#second entangled pair
qc.h(1)
qc.cx(2,3)

#qubits 4,5,6,7,8,9 were measured as either 0 or 1, create the circuit accordingly
for index,s in enumerate(measurement_string):
if s == '1':
qc.x(index + 4)

return Statevector.from_instruction(qc)

def epr_product_2(measurement_string):
#first entangled pair
qc = QuantumCircuit(num_qubits)

qc.x(0)
qc.h(0)

qc.cx(0,1)

#second entangled pair
qc.x(2)
qc.h(2)

qc.cx(2,3)

for index,s in enumerate(measurement_string):
if s == '1':
qc.x(index + 4)

return Statevector.from_instruction(qc)

def epr_product_1_2(measurement_string):
#firdt entangled pair
qc = QuantumCircuit(num_qubits)

qc.h(0)

qc.cx(0,1)

#second entangled pair
qc.x(2)
qc.h(2)

qc.cx(2,3)

for index,s in enumerate(measurement_string):
if s == '1':
qc.x(index + 4)

return Statevector.from_instruction(qc)

def epr_product_2_1(measurement_string):
#first entangled pair
qc = QuantumCircuit(num_qubits)

#second entangled pair
qc.x(0)
qc.h(0)

qc.cx(0,1)

qc.h(2)

qc.cx(2,3)

for index,s in enumerate(measurement_string):
if s == '1':
qc.x(index + 4)

return Statevector.from_instruction(qc)

def epr_product_3(measurement_string):
#first entangled pair
qc = QuantumCircuit(num_qubits)

qc.x(1)
qc.h(0)

qc.cx(0,1)

#second entangled pair
qc.x(3)
qc.h(2)

qc.cx(2,3)

for index,s in enumerate(measurement_string):
if s == '1':
qc.x(index + 4)

return Statevector.from_instruction(qc)

def epr_product_4(measurement_string):
#first entangled pair
qc = QuantumCircuit(num_qubits)

qc.x(1)
qc.h(0)

qc.z(0)
qc.z(1)

#second entangled pair
qc.x(3)
qc.h(2)

qc.z(2)
qc.z(3)

for index,s in enumerate(measurement_string):
if s == '1':
qc.x(index + 4)

return Statevector.from_instruction(qc)

def epr_product_3_4(measurement_string):
#first entangled pair
qc = QuantumCircuit(num_qubits)

qc.x(1)
qc.h(0)

qc.cx(0,1)

#second entangled pair
qc.x(3)
qc.h(2)

qc.z(2)
qc.z(3)

for index,s in enumerate(measurement_string):
if s == '1':
qc.x(index + 4)

return Statevector.from_instruction(qc)

def epr_product_4_3(measurement_string):
#first entangled pair
qc = QuantumCircuit(num_qubits)

qc.x(1)
qc.h(0)
qc.z(0)
qc.z(1)

qc.x(3)
qc.h(2)
qc.cx(2,3)

for index,s in enumerate(measurement_string):
if s == '1':
qc.x(index + 4)

return Statevector.from_instruction(qc)

def epr_product_1_3(measurement_string):
#first entangled pair
qc = QuantumCircuit(num_qubits)

qc.h(0)

qc.cx(0,1)

#second entangled pair
qc.x(3)
qc.h(2)

qc.cx(2,3)

for index,s in enumerate(measurement_string):
if s == '1':
qc.x(index + 4)

return Statevector.from_instruction(qc)

def epr_product_3_1(measurement_string):
#first entangled pair
qc = QuantumCircuit(num_qubits)

#second entangled pair
qc.x(1)
qc.h(0)

# Apply a CNOT
qc.cx(0,1)

qc.h(2)

qc.cx(2, 3)

for index,s in enumerate(measurement_string):
if s == '1':
qc.x(index + 4)

return Statevector.from_instruction(qc)

def epr_product_2_3(measurement_string):
#first entangled pair
qc = QuantumCircuit(num_qubits)

qc.x(0)
qc.h(0)

qc.cx(0,1)

#second entangled pair
qc.x(3)
qc.h(2)

qc.cx(2,3)

for index,s in enumerate(measurement_string):
if s == '1':
qc.x(index + 4)

return Statevector.from_instruction(qc)

def epr_product_3_2(measurement_string):
#first entangled pair
qc = QuantumCircuit(num_qubits)

qc.x(1)
qc.h(0)

qc.cx(0,1)

qc.x(2)
qc.h(2)

qc.cx(2,3)

for index,s in enumerate(measurement_string):
if s == '1':
qc.x(index + 4)

return Statevector.from_instruction(qc)

def epr_product_1_4(measurement_string):
#first entangled pair
qc = QuantumCircuit(num_qubits)

qc.h(0)

qc.cx(0,1)

#second entangled pair
qc.x(3)
qc.h(2)

qc.z(2)
qc.z(3)

for index,s in enumerate(measurement_string):
if s == '1':
qc.x(index + 4)

return Statevector.from_instruction(qc)

def epr_product_4_1(measurement_string):
#first entangled pair
qc = QuantumCircuit(num_qubits)

qc.x(1)
qc.h(0)
qc.z(0)
qc.z(1)

qc.h(2)

qc.cx(2,3)

for index,s in enumerate(measurement_string):
if s == '1':
qc.x(index + 4)

return Statevector.from_instruction(qc)

def epr_product_2_4(measurement_string):
#first entangled pair
qc = QuantumCircuit(num_qubits)

qc.x(0)
qc.h(0)

qc.cx(0,1)

#second entangled pair
qc.x(3)
qc.h(2)

qc.z(2)
qc.z(3)

for index,s in enumerate(measurement_string):
if s == '1':
qc.x(index + 4)

return Statevector.from_instruction(qc)

def epr_product_4_2(measurement_string):
#first entangled pair
qc = QuantumCircuit(num_qubits)

qc.x(1)
qc.h(0)

qc.z(0)
qc.z(1)

qc.x(2)
qc.h(2)

qc.cx(2,3)

for index,s in enumerate(measurement_string):
if s == '1':
qc.x(index + 4)

return Statevector.from_instruction(qc)


Then I compared the statevectors obtained from the above circuit with the statevector that I have(called 'my_state_vec' in code) using the equiv function:

def check_epr_pairs():

measurement = '010101'

epr_products = []
epr_products.append(epr_product_1(measurement))
epr_products.append(epr_product_2(measurement))
epr_products.append(epr_product_3(measurement))

epr_products.append(epr_product_4(measurement))
epr_products.append(epr_product_1_2(measurement))
epr_products.append(epr_product_1_3(measurement))

epr_products.append(epr_product_1_4(measurement))
epr_products.append(epr_product_2_1(measurement))
epr_products.append(epr_product_2_3(measurement))

epr_products.append(epr_product_2_4(measurement))
epr_products.append(epr_product_3_1(measurement))
epr_products.append(epr_product_3_2(measurement))

epr_products.append(epr_product_3_4(measurement))
epr_products.append(epr_product_4_1(measurement))
epr_products.append(epr_product_4_2(measurement))
epr_products.append(epr_product_4_3(measurement))

for ep in epr_products:
#Check if statevectors are equivalent
if my_state_vec.equiv(ep):
print("Statevec equivalent with", ep)


I have theoretical reason to believe that my_state_vec should be such that qubits 0,1 and 2,3 are EPR pairs. But I find that it is not equivalent to any of the EPR product statevectors.

The statevector produced by the qubits 0,1,2,3 is: $$\begin{bmatrix}-0.25\\ -0.25\\ 0.25\\ 0.25\\ 0.25\\ 0.25\\ 0.25\\ 0.25\\ 0.25\\ -0.25\\ -0.25\\ 0.25\\ -0.25\\ 0.25\\ -0.25\\ 0.25\end{bmatrix}$$

Is there a better way of checking for EPR pair products?

• A high-level suggestion. To see if qubits A and B form an EPR pair (1) compute density matrix of the two-qubit system AB. It should have rank 1, i.e. be a pure state. (2) compute the density matrix of A alone, it should be an identity matrix if A and B form a Bell pair. Jul 26, 2022 at 19:49
• @NikitaNemkov I tried this suggestion. Step (1) checked out and AB is a pure state. Step (2), I get the density matrix of A as a diagonal matrix (all diagonal entries = 0.5). So this suggests, AB don't form a bell pair? Aug 18, 2022 at 0:44
• I also updated the question (added the resultant statevector of the 4 qubits) Aug 18, 2022 at 0:53
• I didn't look at your statevector, but diagonal matrix with entries 0.5 sounds exactly right. Actually, a density matrix can not be an identity because it must have unit trace. So, saying that $\rho$ is an identity is often a jargon for $\rho=\mathbb{I}/2^n$. Aug 18, 2022 at 17:35