I was creating a simple implementation of the Hadamard test when I came across the following, what seems like, strange behavior.
Consider the following snippet, which is part of the computating the expected value of the Pauli-$X$ operator with respect to the state $|\psi\rangle = \frac{1}{2}|0\rangle + \frac{\sqrt{3}}{2}|1\rangle$. After measuring the ancillary, we compute $Prob(|0\rangle) - Prob(|1\rangle) = Re(\langle \psi |X|\psi\rangle)$ (see https://en.wikipedia.org/wiki/Hadamard_test_(quantum_computation))
Circuit 1:
reg = QuantumRegister(2)
classical = ClassicalRegister(1)
qc = QuantumCircuit(reg, classical)
def op():
r = QuantumRegister(1)
gate = QuantumCircuit(r)
gate.append(Operator(Pauli(label='X')), [*r])
return gate.to_gate()
qc.initialize([1/2, np.sqrt(3)/2], [reg[1]])
qc.barrier()
qc.h(reg[0])
qc.append(op().control(),[*reg])
qc.h(reg[0])
qc.barrier()
qc.measure(reg[0], classical)
job = execute(qc, Aer.get_backend('qasm_simulator'), shots=1024)
print(job.result().get_counts())
The result from QASM is $\{'0': 499, '1': 525\}$, which is not correct. I verified with Qiskit's snapshot expectation that expected value is $\frac{\sqrt{3}}{2}$. This was also verified through density matrix computation with numpy.
However, if I replace this with the $CX$ as opposed to the Pauli-$X$ from the quantum info package of qiskit.
Circuit 2
reg = QuantumRegister(2)
classical = ClassicalRegister(1)
qc = QuantumCircuit(reg, classical)
qc.initialize([1/2, np.sqrt(3)/2], [reg[1]])
qc.barrier()
qc.h(reg[0])
qc.cx(reg[0], reg[1])
qc.h(reg[0])
qc.barrier()
qc.measure(reg[0], classical)
job = execute(qc, Aer.get_backend('qasm_simulator'), shots=1024)
print(job.result().get_counts())
I got the following: $\{'0': 962, '1': 62\}$, which taking $(962-62)/1024$ produces an approximation to the correct result.
Next I printed out the transpilation of Circui1
t = transpile(qc, Aer.get_backend('qasm_simulator'), basis_gates=['u1','u2','u3','cx'])
t.draw(output='mpl')
I verified that this circuit does in fact produce the same result as Circuit 1. However, I then did a density matrix computation using numpy, which also produced the same result. However, this revealed that there is an incorrect phase difference, not a global phase difference. If we view the bipartite system as $reg[0] \otimes reg[1]$ where $reg[0]$ is the ancillary control bit, we get a unitary matrix for the gates in between the second set of barriers above (not including the Hadamards) like the following: $\begin{bmatrix} -1 & 0 &0 & 0 \\ 0 & -1 &0 & 0 \\ 0 & 0 &0 & -i \\ 0 & 0 &-i & 0 \\ \end{bmatrix}$
I verified that my computation in numpy was correct through multiple test; as a sanity check I compared my the matrices in the numpy simulation against those printed out by qiskit in the transpilation. So this circuit is not the same as $CX$.
Any ideas as to why I am seeing this?
I thought this might be a Qiskit issue, but I posted it here in case I am missing something, which most likely the case.