Unexpected behaviour with qml.QubitUnitary

PennyLane Help
1

I want to apply a custom unitary operation, for that I am building the matrix and then using qml.QubitUnitary to use it as a circuit operation.

Building the matrix

import pennylane as qml
from pennylane import numpy as np

one = np.diag([1,0,0,0])
two = np.diag([0,1,0,0])
three = np.diag([0,0,1,0])
four = np.diag([0,0,0,1])
numbers = [one,two,three,four]

def R_R(theta_R):
    res = np.zeros((8,8),dtype=complex)
    control_number = 0
    for i in range(4):
        if i == control_number:
            continue
        res += np.kron(np.eye(2),numbers[3-i])
    res += np.kron(qml.RY(2*theta_R,wires=1).matrix(),numbers[3-control_number])
    return res

def R_G(theta_G):
    res = np.zeros((8,8),dtype=complex)
    control_number = 1
    for i in range(4):
        if i == control_number:
            continue
        res += np.kron(np.eye(2),numbers[3-i])
    res += np.kron(qml.RY(2*theta_G,wires=1).matrix(),numbers[3-control_number])
    return res

def R_B(theta_B):
    res = np.zeros((8,8),dtype=complex)
    control_number = 2
    for i in range(4):
        if i == control_number:
            continue
        res += np.kron(np.eye(2),numbers[3-i])
    res += np.kron(qml.RY(2*theta_B,wires=1).matrix(),numbers[3-control_number])
    return res

def R_prime(theta_R,theta_G,theta_B):
    # return np.linalg.multi_dot([R_R(theta_R),R_G(theta_G),R_B(theta_B)])
    return R_R(theta_R) @ R_G(theta_G) @ R_B(theta_B)

def R_i_mat(theta_R,theta_G,theta_B,i,wires):
    U = R_prime(theta_R,theta_G,theta_B)
    i_dens_mat = np.zeros((2**(len(wires)-3),2**(len(wires)-3)),dtype=complex)
    i_dens_mat[i,i] = 1
    R_i = np.kron(U,i_dens_mat)
    for k in range(2**(len(wires)-3)):
        if k == i:
            continue
        k_dens_mat = np.zeros((2**(len(wires)-3),2**(len(wires)-3)),dtype=complex)
        k_dens_mat[k,k] = 1
        R_i += np.kron(np.eye(2**(3)),k_dens_mat)
    return R_i

H = qml.Hadamard(wires=0).matrix()
def initial_hadamard(n):
    res = np.eye(2)
    for i in range(n-1):
        res = np.kron(res,H)
    return res

def embedding_matrix(image,n):
    mat = initial_hadamard(n) 
    i = 0
    for j in range(2):
        for k in range(2):
            mat = R_i_mat(image[0,j,k],image[1,j,k],image[2,j,k],i,range(n)) @ mat
            i += 1
    return mat

# set seed for reproducibility
np.random.seed(42)
image_test = np.random.rand(3,2,2)
print(image_test)

CASE I: matrix multiplication

Now if use matrix multiplication from numpy and apply the embedding matrix to the |0> state I get the expected output (which I do not copy because it is a 32 dimensional array, you may run the code on your own computer)

# set seed for reproducibility
np.random.seed(42)
image_test = np.random.rand(3,2,2)
print(image_test)

EM = embedding_matrix(image_test,3+2)

# create initial state
initial_state = np.zeros(2**5,dtype=complex)
initial_state[0] = 1

# correct result
print(EM@initial_state)

CASE II: generate matrix within qnode

However, if I use qml.QubitUnitary to apply the custum unitary EM as a circuit operation, I get a different result:

dev = qml.device('default.qubit',wires=5)

@qml.qnode(dev)
def EM_circuit():
    U = embedding_matrix(image_test,3+2)
    qml.QubitUnitary(U,wires=range(5))
    return qml.state()

print(EM_circuit())

CASE III: pass custom unitary matrix as argument

Interestingly enough, if I don’t generate the matrix EM inside the qnode, but rather pass it as an argument, I get the expected result:

@qml.qnode(dev)
def U_circuit(U):
    qml.QubitUnitary(U,wires=range(5))
    return qml.state()

print(U_circuit(EM)) # correct result as in case 1 (matrix multiplication)

I don’t understand the difference between the last two cases. Besides, the way of doing it that does not lead to the correct result fits better the pipeline in which I want to integrate this operation. Am I doing something wrong or is this a bug?

Thank you very much in advance for your help!

PS: qml.about()

Name: PennyLane
Version: 0.36.0
Summary: PennyLane is a cross-platform Python library for quantum computing, quantum machine learning, and quantum chemistry. Train a quantum computer the same way as a neural network.
Home-page: https://github.com/PennyLaneAI/pennylane
Author: 
Author-email: 
License: Apache License 2.0
Location: 
Requires: appdirs, autograd, autoray, cachetools, networkx, numpy, pennylane-lightning, requests, rustworkx, scipy, semantic-version, toml, typing-extensions
Required-by: PennyLane_Lightning

Platform info:           Windows-11-10.0.22621-SP0
Python version:          3.12.3
Numpy version:           1.26.4
Scipy version:           1.13.0
Installed devices:
- default.clifford (PennyLane-0.36.0)
- default.gaussian (PennyLane-0.36.0)
- default.mixed (PennyLane-0.36.0)
- default.qubit (PennyLane-0.36.0)
- default.qubit.autograd (PennyLane-0.36.0)
- default.qubit.jax (PennyLane-0.36.0)
- default.qubit.legacy (PennyLane-0.36.0)
- default.qubit.tf (PennyLane-0.36.0)
- default.qubit.torch (PennyLane-0.36.0)
- default.qutrit (PennyLane-0.36.0)
- default.qutrit.mixed (PennyLane-0.36.0)
- null.qubit (PennyLane-0.36.0)
- lightning.qubit (PennyLane_Lightning-0.37.0)
2

Hi @ARO-GZ ,

Welcome to the Forum and thank you for your excellent post! The way you wrote it helped me a lot to find the issue.

I’ll start with a short answer and the fix for your problem, followed by a detailed explanation.

Short answer:
You’re accidentally adding gates to your circuit in Case II.
You can fix it with one line of code!
with qml.QueuingManager.stop_recording():

See how to use it below:

dev = qml.device('default.qubit',wires=5)

@qml.qnode(dev)
def EM_circuit():
    with qml.QueuingManager.stop_recording():
      U = embedding_matrix(image_test,3+2)
    qml.QubitUnitary(U,wires=range(5))
    return qml.state()

print(EM_circuit())

Detailed explanation:

PennyLane works by queuing up all quantum operations in a tape that then runs them on the device of your choice. This is what allows you to create sub-circuits outside of your QNode (such as layers) and then easily add them to your circuit.

Example: In the example below, I’m creating a function called layer which allows me to define the layer once and then call it within my QNode as many times as I want.

In this example we get the layers of Hadamard gates properly added to the circuit.

n_wires = 5
dev = qml.device('default.qubit',wires=n_wires)

def layer(n_wires):
  for w in range(n_wires):
    qml.Hadamard(wires=w)

@qml.qnode(dev)
def circuit(n_layers):
    for i in range(n_layers):
      layer(n_wires)
    return qml.probs()

n_layers = 3
qml.draw_mpl(circuit, style='pennylane')(n_layers);

image
image619×619 18.1 KB

The issue in your case is that you have qml.RY within your functions R_R, R_G, and R_B so they’re getting added to the circuit without you realizing it.

I hope this helps you move along with your project and prevent this from causing issues in your code again!

PS: we’ve now released version v0.37.0 of PennyLane! I recommend that you always keep your code updated to run with the latest version :smiley:

3

Thank you very much @CatalinaAlbornoz for the quick reply and detailed explanation!! Your answer did solved my issue.

1 Like