Bug when implementing Controlled Arbitrary Unitary operations when num_qubits > 3

[stfnmangini]

Informations

• Qiskit Aer version: 0.10.4
• Python version: 3.8.12
• Operating system: MacOS Monterey

What is the current behavior?

The implementation of controlled operation for arbitrary unitaries, for example generated via qiskit.quantum_info.random_unitary, is incorrect when the number of qubits in the unitary to be controlled is greater than three.

To my understanding, there may be a problem with the circuit implementation of the random unitary when you ask it to be controlled, as some undesired approximation is happening somewhere. In fact, the problem is not present when trying to implement a controlled version of a circuit with an explicit form.

Steps to reproduce the problem

Pick a random unitary via qiskit.quantum_info.random_unitary and ask Qiskit to implement it on a quantum circuit.

import numpy as np
import qiskit as qk
from qiskit.quantum_info import random_unitary
from qiskit.circuit.random import random_circuit
uny_bknd = qk.Aer.get_backend('unitary_simulator')

num_qubits = 3
qc = qk.QuantumCircuit(num_qubits)

rand_u = random_unitary(2**num_qubits)
qc.unitary(rand_u, range(num_qubits))

# Inspect the unitary implemented by the circuit
u_circ = qk.execute(qc, uny_bknd).result().get_unitary()

# Sanity check
print("Random Unitary = Unitary from circuit? ->", np.all(np.isclose(np.asarray(rand_u), np.asarray(u_circ))))
# ---> TRUE

qc.decompose().draw()

Then create a quantum circuit with a controlled version of the unitary defined above

qc_ctrl = qc.control(1)

# Unitary implemented by the quantum circuit
cu_circ = qk.execute(qc_ctrl, uny_bknd).result().get_unitary()

# Theory definition of controlled operation (|0x0| x Id + |1x1| x U)
proj0 = np.array([[1, 0], [0, 0]])
proj1 = np.array([[0, 0], [0, 1]])
cu_theory = np.kron(np.eye(2**num_qubits), proj0) + np.kron(np.asarray(u_circ), proj1) # order to cope with Qiskit little-endian ordering

print("Circuit = Theory? ->", np.all(np.isclose(cu_circ, cu_theory)))
# ---> FALSE

qc_ctrl.decompose().draw()

The unitary implemented by the controlled circuit is not equal to the theoretical one when num_qubits >= 3, while is it correct for 1 and 2 qubits. My guess is that when Qiskit is asked to perform a controlled version, it only finds an approximate, rather than exact, circuital implementation of the unitary, without explicitly mentioning it.

This does not happen if instead of an arbitrary unitary you use an explicit circuit, for example:

num_qubits = 3
qc = qk.QuantumCircuit(num_qubits)

# Generate an *explicit* random circuit, already in gate-based form
qc.append(random_circuit(num_qubits, depth=2), range(num_qubits))
qc = qk.transpile(qc, uny_bknd) #otherwise compilation error on "non-unitary" gates being not controllable 

u_circ = qk.execute(qc, uny_bknd).result().get_unitary()

qc.decompose().draw()

# CONTROLLED CIRCUIT
qc_ctrl = qc.control(1)
cu_circ = qk.execute(qc_ctrl, uny_bknd).result().get_unitary()

# Definition of controlled operation
proj0 = np.array([[1, 0], [0, 0]])
proj1 = np.array([[0, 0], [0, 1]])
cu_theory = np.kron(np.eye(2**num_qubits), proj0) + np.kron(np.asarray(u_circ), proj1) # order to cope with Qiskit little-endian ordering

print("Circuit = Theory? ->", np.all(np.isclose(cu_circ, cu_theory)))
# ----> TRUE

qc_ctrl.decompose().draw()

Version Information

Qiskit Software	Version
qiskit-terra	0.20.1
qiskit-aer	0.10.4
qiskit-ignis	0.7.0
qiskit-ibmq-provider	0.19.1
qiskit	0.36.1
Python version	3.8.12

What is the expected behavior?

The controlled unitary should be correct.

Suggested solutions

Inspect what happens when calling controlled on a unitary gate rather than an explicit circuit.

[jakelishman]

This may well be related to the approximation_degree setting of the transpiler in Terra (I think this bug might be related to Terra, not to Aer?), which is under quite a bit of discussion about design philosophy - see Qiskit/qiskit-terra#8043, for example. If you transpile with approximation_degree = 0, I think that should have the behaviour you expect.

[stfnmangini]

I agree on the problem being somewhere in the transpiler, probably using a low approximation_degree. However I cannot obtain the exact unitary even if transpiling the circuit, playing with the optimization_level and approximation_degree. I feel there is something going bad when going through the control function, apparently overriding the user commands.

As for the appropriate page, you're probably right, I didn't know where to open this issue as I'm not sure of its origin 😅

[jakelishman]

When you're setting approximation_degree, are you doing it in a separate transpile pass, or trying to pass it through execute? (We're in discussion about completely deprecating execute, because it was only ever meant to be a quick "just run it" method, and doesn't allow full access to the transpiler.) The full path is usually:

from qiskit import Aer, transpile
circuit = ...
backend = Aer.get_backend("aer_simulator")
transpiled = transpile(circuit, backend, **transpile_options)
backend.run(transpiled).result().get_unitary()

As a faster check, you can also do qiskit.quantum_info.Operator(circuit) to get the unitary of small-ish circuits without going to Aer.

At any rate, I'll transfer this over to Terra (no big deal).