Iterative quantum phase estimation - real device

I'm trying to implement the iterative quantum phase estimation on a real (IBM) quantum computer. I'm using the code below. When I run this code on a simulator the results are the expected ones, but when running on a real device the results don't follow any pattern.

import matplotlib.pyplot as plt
# QML
from pennylane import numpy as np
from qiskit import *
from qiskit.visualization import plot_histogram
from qiskit.tools.monitor import job_monitor

shots=32000

# Key with the maximum probability - maior in Portuguese
def maior(dic):
m=list(dic)[0]
for n in list(dic):
if dic[n]>dic[m]:
m=n
return m

from key import tok
from qiskit import IBMQ #2
IBMQ.save_account(tok, overwrite=True)

backend = provider.get_backend('ibmq_toronto') #4

# ## Iterative Quantum Phase Estimation Algorithm
def get_circuit_phase(t,
QC,
clbits,
qubits,
ancilla,
backend=None,
):
mycircuit=QuantumCircuit(2)
mycircuit.cx(0,1)
mycircuit.rx(2*t,0)
mycircuit.rz(2*t,1)
mycircuit.cx(0,1)
mycircuit.cy(0,1)
mycircuit.ry(2*t,0)
mycircuit.cy(0,1)
#print(mycircuit.draw())

# Circuit -> controlled gate
CU=mycircuit.to_gate().control(1)

res = []
phase = -2 * np.pi
factor = 0
iterations = 3

# generate the qubit list on which the Unitary is applied
qargs = [ancilla]
for q in qubits:
qargs.append(q)

exponent = 2 ** (iterations - 1)
for it in range(iterations):
# start
QC.reset(ancilla)
QC.h(ancilla)
inv_phase = phase * factor

QC.p(inv_phase, ancilla)

# add the controlled Unitary of iteration it

# need to add exponential amount of matrices
for _ in range(int(exponent)):
QC = QC.compose(CU, qubits=qargs)
exponent /= 2

QC.h(ancilla)
QC.measure(ancilla, clbits[it])

if backend == None:  # simulating
backend=Aer.get_backend("qasm_simulator")

t_qpe = transpile(QC, backend,optimization_level=3)
job = backend.run(t_qpe, shots=shots)
job_monitor(job)

counts = job.result().get_counts(QC)

#  mai is the key with the most probability.
mai=maior(counts)

# Save the bit
res.append(int(mai[3-it-1]))

# if bit measured is 1
if mai[3-it-1] == "1":
factor += 1 / 2  # add the phase factor

factor = factor / 2  # shift each towards one weight right

# phase has now been stored in the clbits
# returning its binary representation

# need to reverse as LSB is stored at the zeroth index and
# not the last
res = res[::-1]

# find decimal phase
dec = 0
weight = 1 / 2
for k in res:
dec += (weight) * k
weight /= 2

return dec

tau=[]
phase=[]
for t in range(0,30,1):
nq = 3    # number of qubits
m = 3    # number of classical bits
q = QuantumRegister(nq,'q')
c = ClassicalRegister(m,'c')
qc = QuantumCircuit(q,c)
qc.h(0)
t=t/5
tau.append(t)
x=get_circuit_phase(t,
QC=qc,
clbits=[0,1,2],
qubits=[1,2],
ancilla=[0], backend=backend
)
phase.append(x)
plt.plot(tau,phase,label="- state")
plt.plot(tau,phase,label="+ state")
plt.xlabel("$$tau$$")
plt.ylabel("$$theta$$")
plt.legend()
plt.show()


I hope anyone can help with this issue, my sincere thanks,

Gabriela Oliveira.

• Please note that not all features available on the simulator can be used on real QPU. For example conditioning on classical register which is employed in middle-of-circuit measurement is one of such. May 15 at 7:22