Issues using step_and_cost as part of an optimization method

PennyLane Help
1

Hello! If applicable, put your complete code example down below. Make sure that your code:

  • is 100% self-contained — someone can copy-paste exactly what is here and run it to
    reproduce the behaviour you are observing
  • includes comments
import pennylane as qml 
from pennylane import qchem
from pennylane import numpy as np
import matplotlib.pyplot as plt

#Defining an algorithm class for Grover search
#This class will have added functionality to optimize the number of calls to the oracle to experimentally verify the optimal value of sqrt(n)

class Grover:
#------------
#    INITIALIZATION, Input number of qubits, shots, and the bit string to be marked for search
#------------    
    def __init__(self, qubits, shots, oracle_state):
        #num of qubits
        self.qubits = qubits
        #num of shots on backend
        self.shots = shots
        
        #marked state defining the oracle - input in the form of an integer
        #integer is converted to bit string which is converted to binary vector for later use in FlipState method
        self.num = oracle_state
        if oracle_state >= 2**self.qubits:
            raise ValueError(f"Insufficient qubits. Input value less than {2**self.qubits}")
        v = []
        B = bin(oracle_state)[2:].zfill(self.qubits)
        for i in range(len(B)):
            v.append(int(B[i]))
        self.oracle = v
        
        #optimal number of oracle+diffuser cycles - this value is used as the default to reproduce the standard Grover search
        self.cycles = int(np.floor(np.sqrt(self.qubits)))

#------------
#    RUN METHOD, This is the core method which is used as input to all other methods
#    Default run uses the optimal number of oracle cycles
#    User can specify alternate number of cycles if desired
#------------    
    def run(self, cycles = None):
        #cycles is an optional input if the user wants to implement Grover with a nonstandard number of oracle calls
        if cycles is None:
            cycles = self.cycles
        #the run method takes in an integer to define the number of cycles - during optimization a tensor is passed to run - this converts the tensor input to one the method can use
        if type(cycles) == np.tensor:
            cycles = int(np.floor(cycles.item()))
        
        #setup for Grover algorithm with initial state [1,1,...,1] (this is more optimal as the diffuser does not need to be conjugated by x gates)
        initial_state = []
        for i in range(self.qubits):
            initial_state.append(1)
            
        #oracle+diffusion
        def OD():
        #oracle
            qml.FlipSign(self.oracle, wires = range(self.qubits))
        #diffusion
            for i in range(self.qubits):
                qml.Hadamard(i)
            qml.FlipSign(initial_state, wires = range(self.qubits))
            for i in range(self.qubits):
                qml.Hadamard(i)
                
        #repeater oracle+diffuser 
        def Oracle_Diffusion(n):
            for i in range(n):
                OD()
                
    #run grover
        dev = qml.device("default.qubit", wires = self.qubits, shots = self.shots)
        @qml.qnode(dev, interface="autograd")
        def G(q):
        #initial state prep
            for i in range(self.qubits):
                qml.PauliX(i)
                qml.Hadamard(i)
        #repeated component
            Oracle_Diffusion(cycles)
        #measurement
            return qml.probs(wires = range(self.qubits))
        #the output of the method is a 1D tensor of probabilities for each bit string given the number of cycles specified
        return (G(cycles))
#------------
#    PLOT METHOD, visualizes probabilities outputted by run method
#    Default run uses the optimal number of oracle cycles
#    User can specify alternate number of cycles if desired
#------------    
    
    def plot(self, cycles = None):
        #same as for run
        if cycles is None:
            cycles = self.cycles
        if type(cycles) == np.tensor:
            cycles = int(np.floor(cycles.item()))
            
        #store y values as output of run method
        y = self.run(cycles)
        bit_strings = [f"{x:0{self.qubits}b}" for x in range(len(y))]
        plt.bar(bit_strings, y, color = "#212121")
        plt.xticks(rotation="vertical")
        plt.xlabel("State label")
        plt.ylabel("Probability Amplitude")
        plt.title("States probabilities amplitudes")
        plt.show()
        
#------------BROKEN METHOD
#    OPT METHOD, utilizes cycle input capability to experimentally recover optimal number of cycles
#------------BROKEN METHOD           
    def opt(self):
        #initialize optimizer using integer step size as cycles must always be an integer
        opt = qml.GradientDescentOptimizer(stepsize=1)
        #initialize training parameter
        theta = np.array(0.0, requires_grad=True)
        #cost function to be trained
        #sends theta -> run method and selects the amplitude of marked state from output tensor
        #value is negative to permit minimization
        Cost_vec = [-self.run(theta)[self.num]]
        angle = [theta]
    
        max_iterations = 100
        convergence_tolerance = 1e-06
    
        for n in range(max_iterations):
            #cost function
            cost = -self.run(theta)[self.num]
            #BROKEN LINE - step and cost does not appear to be working with this setup
            theta, prev_prob = opt.step_and_cost(cost, theta)
            
            #store subsequent optimization steps in the corresponding vectors
            Cost_vec.append(self.run(theta)[self.num])
            angle.append(theta)
            
            #convergence condition left over from VQE demo where opt code was first sourced - doesn't apply to this problem
            #conv = np.abs(Prob[-1] - prev_prob)
              
            if n % 2 == 0:
                print(f"Step = {n}, probability = {Prob[-1]:.8f}")
              
            #if conv <= convergence_tolerance:
            #    break

        print("\n" f"Final amplitude of the target state = {Prob[-1]:.8f}")
        print("\n" f"Optimal value of the circuit parameter = {angle[-1]:.4f}")

TypeError Traceback (most recent call last)
/tmp/ipykernel_70/2615489177.py in <cell line: 1>()
----> 1 G.opt()

/tmp/ipykernel_70/3335889083.py in opt(self)
89
90 cost = -self.run(theta)[self.num]
—> 91 theta, prev_prob = opt.step_and_cost(cost, theta)
92
93 Prob.append(self.run(theta)[self.num])

/opt/conda/envs/pennylane/lib/python3.9/site-packages/pennylane/optimize/gradient_descent.py in step_and_cost(self, objective_fn, grad_fn, *args, **kwargs)
57 “”"
58
—> 59 g, forward = self.compute_grad(objective_fn, args, kwargs, grad_fn=grad_fn)
60 new_args = self.apply_grad(g, args)
61

/opt/conda/envs/pennylane/lib/python3.9/site-packages/pennylane/optimize/gradient_descent.py in compute_grad(objective_fn, args, kwargs, grad_fn)
115 “”"
116 g = get_gradient(objective_fn) if grad_fn is None else grad_fn
→ 117 grad = g(*args, **kwargs)
118 forward = getattr(g, “forward”, None)
119

/opt/conda/envs/pennylane/lib/python3.9/site-packages/pennylane/_grad.py in call(self, *args, **kwargs)
113 return ()
114
→ 115 grad_value, ans = grad_fn(*args, **kwargs)
116 self._forward = ans
117

/opt/conda/envs/pennylane/lib/python3.9/site-packages/autograd/wrap_util.py in nary_f(*args, **kwargs)
18 else:
19 x = tuple(args[i] for i in argnum)
—> 20 return unary_operator(unary_f, x, *nary_op_args, **nary_op_kwargs)
21 return nary_f
22 return nary_operator

/opt/conda/envs/pennylane/lib/python3.9/site-packages/pennylane/_grad.py in _grad_with_forward(fun, x)
131 difference being that it returns both the gradient and the forward pass
132 value.“”"
→ 133 vjp, ans = _make_vjp(fun, x)
134
135 if not vspace(ans).size == 1:

/opt/conda/envs/pennylane/lib/python3.9/site-packages/autograd/core.py in make_vjp(fun, x)
8 def make_vjp(fun, x):
9 start_node = VJPNode.new_root()
—> 10 end_value, end_node = trace(start_node, fun, x)
11 if end_node is None:
12 def vjp(g): return vspace(x).zeros()

/opt/conda/envs/pennylane/lib/python3.9/site-packages/autograd/tracer.py in trace(start_node, fun, x)
8 with trace_stack.new_trace() as t:
9 start_box = new_box(x, t, start_node)
—> 10 end_box = fun(start_box)
11 if isbox(end_box) and end_box._trace == start_box._trace:
12 return end_box._value, end_box._node

/opt/conda/envs/pennylane/lib/python3.9/site-packages/autograd/wrap_util.py in unary_f(x)
13 else:
14 subargs = subvals(args, zip(argnum, x))
—> 15 return fun(*subargs, **kwargs)
16 if isinstance(argnum, int):
17 x = args[argnum]

TypeError: ‘tensor’ object is not callable

And, finally, make sure to include the versions of your packages. Specifically, show us the output of qml.about().

Name: PennyLane
Version: 0.28.0
Summary: PennyLane is a Python quantum machine learning library by Xanadu Inc.
Home-page: GitHub - PennyLaneAI/pennylane: PennyLane is a cross-platform Python library for differentiable programming of quantum computers. Train a quantum computer the same way as a neural network.
Author:
Author-email:
License: Apache License 2.0
Location: /opt/conda/envs/pennylane/lib/python3.9/site-packages
Requires: appdirs, autograd, autoray, cachetools, networkx, numpy, pennylane-lightning, requests, retworkx, scipy, semantic-version, toml
Required-by: PennyLane-Cirq, PennyLane-Lightning, PennyLane-qiskit, pennylane-qulacs, PennyLane-SF

Platform info: Linux-5.4.209-116.367.amzn2.x86_64-x86_64-with-glibc2.31
Python version: 3.9.15
Numpy version: 1.23.5
Scipy version: 1.10.0
Installed devices:

  • default.gaussian (PennyLane-0.28.0)
  • default.mixed (PennyLane-0.28.0)
  • default.qubit (PennyLane-0.28.0)
  • default.qubit.autograd (PennyLane-0.28.0)
  • default.qubit.jax (PennyLane-0.28.0)
  • default.qubit.tf (PennyLane-0.28.0)
  • default.qubit.torch (PennyLane-0.28.0)
  • default.qutrit (PennyLane-0.28.0)
  • null.qubit (PennyLane-0.28.0)
  • cirq.mixedsimulator (PennyLane-Cirq-0.28.0)
  • cirq.pasqal (PennyLane-Cirq-0.28.0)
  • cirq.qsim (PennyLane-Cirq-0.28.0)
  • cirq.qsimh (PennyLane-Cirq-0.28.0)
  • cirq.simulator (PennyLane-Cirq-0.28.0)
  • lightning.qubit (PennyLane-Lightning-0.28.2)
  • strawberryfields.fock (PennyLane-SF-0.20.1)
  • strawberryfields.gaussian (PennyLane-SF-0.20.1)
  • strawberryfields.gbs (PennyLane-SF-0.20.1)
  • strawberryfields.remote (PennyLane-SF-0.20.1)
  • strawberryfields.tf (PennyLane-SF-0.20.1)
  • qiskit.aer (PennyLane-qiskit-0.28.0)
  • qiskit.basicaer (PennyLane-qiskit-0.28.0)
  • qiskit.ibmq (PennyLane-qiskit-0.28.0)
  • qiskit.ibmq.circuit_runner (PennyLane-qiskit-0.28.0)
  • qiskit.ibmq.sampler (PennyLane-qiskit-0.28.0)
  • qulacs.simulator (pennylane-qulacs-0.28.0)
    ​
2

Hey @Joan!

There’s quite a bit going on, and I think it would best serve us to work with a much simpler example that highlights the key parts to any circuit optimization in PennyLane. The code I’m about to show is in the qubit rotation demo.

In your code, there’s a parameterized circuit whose parameters need to be tuned to optimize some cost function. Let’s first define a parameterized circuit:

import pennylane as qml
from pennylane import numpy as np

dev = qml.device("default.qubit", wires=1)

@qml.qnode(dev)
def circuit(params):
    qml.RX(params[0], wires=0)
    qml.RY(params[1], wires=0)
    return qml.expval(qml.PauliZ(0))

This circuit takes a vector of parameters — length 2 — and returns an expectation value. Now, we need to define a separate function that we want to optimize — the cost function. It could literally be anything so long as it returns a scalar output (e.g., mean square loss, KL divergence, etc. — any loss function). Here, let’s just use the circuit’s output as the cost function:

def cost(x):
    return circuit(x)

init_params = np.array([0.011, 0.012], requires_grad=True)
print(cost(init_params))

0.9998675058299389

Now we can go ahead and optimize the cost function. To do that, we define an optimizer then update the parameters with the optimizer’s step method, which takes the loss function as the first argument and the arguments of the loss function next. step_and_cost returns both the parameters and the loss function value.

# initialise the optimizer
opt = qml.GradientDescentOptimizer(stepsize=0.4)

# set the number of steps
steps = 5
# set the initial parameter values
params = init_params

for i in range(steps):
    # update the circuit parameters
    params = opt.step(cost, params)

    print("Cost after step {:5d}: {: .7f}".format(i + 1, cost(params)))

print("Optimized rotation angles: {}".format(params))

It sounds like you’re aware of most of this already, but it helps to simplify what you’re doing if you get stuck. If what you’re trying to do is more complex than what’s above, then I can definitely help! Just fill me in on the details :slight_smile:

3

Hey Isaac,

Thanks for your response. I am familiar with the qubit rotation tutorial - it doesn’t seem to make a difference whether I use opt.step or opt.step_and_cost.

It may be more helpful to compare how I am trying to optimize my problem to the opt module in the VQE demo, as this was my starting point.

I believe that the problem I am running into is at the interaction of GradientDescentOptimizer’s methods and my run() module.

In your example code, you call cost(init_params), and this returns a zero-dimensional tensor:
tensor(0.99986751, requires_grad=True)

In my code, when I initialize an instance of Grover G, then G.run(theta) returns a 1D tensor where each entry corresponds to the measurement probability of a particular bit string after applying some number theta of oracle+diffusion modules.

For my optimization problem, I define the 3 qubit grover with 100 shots that searches for the bit string 000 via

G = Grover(3,100,0)
def cost(x):
return -G.run(x)[G.num]

If I define

theta = np.array(0.0, requires_grad=True)
then cost(theta) gives
tensor(-0.111, requires_grad=True)

This outputs a zero-dimensional tensor corresponding to the measurement probability of the bit string which has been marked for search after zero calls to the oracle. This is the quantity I want to maximize (or minimize in this case as I have made it negative to be compatible with gradient descent).

What does step or step_and_cost expect for inputs? I have tried to ensure that cost(theta) and theta are both tensors with requires_grad=True - but is it possible that I am feeding the opt methods something they don’t know how to work with?

I added
def cost(x)
return -G.run(x)[G.num]
to the opt method as I thought it was cleaner, so now I have

def opt(self):
#initialize optimizer using integer step size as cycles must always be an integer
opt = qml.GradientDescentOptimizer(stepsize=1)
#initialize training parameter
#dev = qml.device(“default.qubit”, wires = self.qubits, shots = self.shots)
#@qml.qnode(dev, interface=“autograd”)
def cost(x):
return -self.run(x)[self.num]
theta = np.array(0.0, requires_grad=True)
#cost function to be trained
#sends theta → run method and selects the amplitude of marked state from output tensor
#value is negative to permit minimization
#Cost_vec = [-self.run(theta)[self.num]]
Cost_vec = [cost(theta)]
angle = [theta]
max_iterations = 5
#convergence_tolerance = 1e-06
for n in range(max_iterations):
#BROKEN LINE - step and cost does not appear to be working with this settup
theta = opt.step(cost, theta)
#store subsequent optimization steps in the corresponding vectors
Cost_vec.append(cost(theta))
angle.append(theta)
#convergence condition left over from VQE demo where opt code was first sourced - doesn’t apply to this problem
#conv = np.abs(Prob[-1] - prev_prob)
if n % 2 == 0:
print(f"Step = {n}, probability = {Cost_vec[-1]:.8f}“)
#if conv <= convergence_tolerance:
# break
print(”\n" f"Final amplitude of the target state = {Prob[-1]:.8f}“)
print(”\n" f"Optimal value of the circuit parameter = {angle[-1]:.4f}")

Which produces a different error: TypeError: ‘ArrayBox’ object cannot be interpreted as an integer

TypeError Traceback (most recent call last)
/tmp/ipykernel_51/2615489177.py in <cell line: 1>()
----> 1 G.opt()

/tmp/ipykernel_51/3346500594.py in opt(self)
121 for n in range(max_iterations):
122 #BROKEN LINE - step and cost does not appear to be working with this settup
→ 123 theta = opt.step(cost, theta)
124
125 #store subsequent optimization steps in the corresponding vectors

/opt/conda/envs/pennylane/lib/python3.9/site-packages/pennylane/optimize/gradient_descent.py in step(self, objective_fn, grad_fn, *args, **kwargs)
86 “”"
87
—> 88 g, _ = self.compute_grad(objective_fn, args, kwargs, grad_fn=grad_fn)
89 new_args = self.apply_grad(g, args)
90

/opt/conda/envs/pennylane/lib/python3.9/site-packages/pennylane/optimize/gradient_descent.py in compute_grad(objective_fn, args, kwargs, grad_fn)
115 “”"
116 g = get_gradient(objective_fn) if grad_fn is None else grad_fn
→ 117 grad = g(*args, **kwargs)
118 forward = getattr(g, “forward”, None)
119

/opt/conda/envs/pennylane/lib/python3.9/site-packages/pennylane/_grad.py in call(self, *args, **kwargs)
113 return ()
114
→ 115 grad_value, ans = grad_fn(*args, **kwargs)
116 self._forward = ans
117

/opt/conda/envs/pennylane/lib/python3.9/site-packages/autograd/wrap_util.py in nary_f(*args, **kwargs)
18 else:
19 x = tuple(args[i] for i in argnum)
—> 20 return unary_operator(unary_f, x, *nary_op_args, **nary_op_kwargs)
21 return nary_f
22 return nary_operator

/opt/conda/envs/pennylane/lib/python3.9/site-packages/pennylane/_grad.py in _grad_with_forward(fun, x)
131 difference being that it returns both the gradient and the forward pass
132 value.“”"
→ 133 vjp, ans = _make_vjp(fun, x)
134
135 if not vspace(ans).size == 1:

/opt/conda/envs/pennylane/lib/python3.9/site-packages/autograd/core.py in make_vjp(fun, x)
8 def make_vjp(fun, x):
9 start_node = VJPNode.new_root()
—> 10 end_value, end_node = trace(start_node, fun, x)
11 if end_node is None:
12 def vjp(g): return vspace(x).zeros()

/opt/conda/envs/pennylane/lib/python3.9/site-packages/autograd/tracer.py in trace(start_node, fun, x)
8 with trace_stack.new_trace() as t:
9 start_box = new_box(x, t, start_node)
—> 10 end_box = fun(start_box)
11 if isbox(end_box) and end_box._trace == start_box._trace:
12 return end_box._value, end_box._node

/opt/conda/envs/pennylane/lib/python3.9/site-packages/autograd/wrap_util.py in unary_f(x)
13 else:
14 subargs = subvals(args, zip(argnum, x))
—> 15 return fun(*subargs, **kwargs)
16 if isinstance(argnum, int):
17 x = args[argnum]

/tmp/ipykernel_51/3360673386.py in cost(x)
1 def cost(x):
----> 2 return -G.run(x)[G.num]

/tmp/ipykernel_51/3346500594.py in run(self, cycles)
73 return qml.probs(wires = range(self.qubits))
74 #the output of the method is a 1D tensor of probabilities for each bit string given the number of cycles specified
—> 75 return (G(cycles))
76 #------------
77 # PLOT METHOD, visualizes probabilities outputted by run method

/opt/conda/envs/pennylane/lib/python3.9/site-packages/pennylane/qnode.py in call(self, *args, **kwargs)
798
799 # construct the tape
→ 800 self.construct(args, kwargs)
801
802 cache = self.execute_kwargs.get(“cache”, False)

/opt/conda/envs/pennylane/lib/python3.9/site-packages/pennylane/qnode.py in construct(self, args, kwargs)
709 “”“Call the quantum function with a tape context, ensuring the operations get queued.”“”
710
→ 711 self._tape = make_qscript(self.func)(*args, **kwargs)
712 self._tape._queue_category = “_ops”
713 self._qfunc_output = self.tape._qfunc_output

/opt/conda/envs/pennylane/lib/python3.9/site-packages/pennylane/tape/qscript.py in wrapper(*args, **kwargs)
1344 def wrapper(*args, **kwargs):
1345 with AnnotatedQueue() as q:
→ 1346 result = fn(*args, **kwargs)
1347
1348 qscript = QuantumScript.from_queue(q)

/tmp/ipykernel_51/3346500594.py in G(q)
69 qml.Hadamard(i)
70 #repeated component
—> 71 Oracle_Diffusion(cycles)
72 #measurement
73 return qml.probs(wires = range(self.qubits))

/tmp/ipykernel_51/3346500594.py in Oracle_Diffusion(n)
57 #repeater oracle+diffuser
58 def Oracle_Diffusion(n):
—> 59 for i in range(n):
60 OD()
61

TypeError: ‘ArrayBox’ object cannot be interpreted as an integer

4

I did a bit more digging, and it seems as if when I pass a tensor like theta → autograd it converts the input to the ArrayBox data type.

This data type is then incompatible with my range(theta) as range needs an integer data type
#repeater oracle+diffuser
def Oracle_Diffusion(n):
for i in range(n):
OD()

I ran into this issue prior to bringing in gradient descent - and I thought I had fixed it by adding this line to the run method

def run(self, cycles = None):
#cycles is an optional input if the user wants to implement Grover with a nonstandard number of oracle calls
if cycles is None:
cycles = self.cycles
#the run method takes in an integer to define the number of cycles - during optimization a tensor is passed to run - this converts the tensor input to one the method can use
if type(cycles) == np.tensor:
cycles = int(np.floor(cycles.item()))

But apparently the issue is still there???

5

Hi @Joan,

ArrayBox issues usually come up when your cost function is updating a global variable or when you do things that Autograd doesn’t like (details here in the docs).

I will take a look at the code you shared to see if I can reproduce the issue and share some more specific insights. Stay tuned! :grinning:

6

Hey @CatalinaAlbornoz - great to hear from you!

Looking forward to learning more, the good news is that if I am right about the cause of the error then it is an artifact of my unusual framing of Grover as a variational algorithm that can only take integer theta

I was building this code as a conceptual bridge for my students to go from non-variational to variational algorithms - I don’t expect the VQA code I am working on next will have this issue since their variational parameters will be more conventional

7

Hi @Joan,

I wasn’t able to find an answer to your issue here so I will forward your question to my colleagues here who may have some insights.

8

Hey @Joan,

I think you can just get away with casting a larger net to make sure cycles is an int. Right now, checking if it’s a numpy array won’t catch subsequent instances when it’s converted to an ArrayBox. Let me know if this works:

class Grover:
    def __init__(self, qubits, shots, oracle_state):
        # num of qubits
        self.qubits = qubits
        # num of shots on backend
        self.shots = shots

        # marked state defining the oracle - input in the form of an integer
        # integer is converted to bit string which is converted to binary vector for later use in FlipState method
        self.num = oracle_state
        if oracle_state >= 2**self.qubits:
            raise ValueError(
                f"Insufficient qubits. Input value less than {2**self.qubits}"
            )
        v = []
        B = bin(oracle_state)[2:].zfill(self.qubits)
        for i in range(len(B)):
            v.append(int(B[i]))
        self.oracle = v

        # optimal number of oracle+diffuser cycles - this value is used as the default to reproduce the standard Grover search
        self.cycles = int(np.floor(np.sqrt(self.qubits)))

    def run(self, cycles=None):
        # cycles is an optional input if the user wants to implement Grover with a nonstandard number of oracle calls
        if cycles is None:
            cycles = self.cycles
        # the run method takes in an integer to define the number of cycles - during optimization a tensor is passed to run - this converts the tensor input to one the method can use
        elif type(cycles) != int:
            cycles = int(np.floor(cycles))

        # setup for Grover algorithm with initial state [1,1,...,1] (this is more optimal as the diffuser does not need to be conjugated by x gates)
        initial_state = []
        for i in range(self.qubits):
            initial_state.append(1)

        # oracle+diffusion
        def OD():
            # oracle
            qml.FlipSign(self.oracle, wires=range(self.qubits))
            # diffusion
            for i in range(self.qubits):
                qml.Hadamard(i)
            qml.FlipSign(initial_state, wires=range(self.qubits))
            for i in range(self.qubits):
                qml.Hadamard(i)

        # repeater oracle+diffuser
        def Oracle_Diffusion(n):
            for i in range(n):
                OD()

        # run grover
        dev = qml.device("default.qubit", wires=self.qubits, shots=self.shots)

        @qml.qnode(dev, interface="autograd")
        def G(q):
            # initial state prep
            for i in range(self.qubits):
                qml.PauliX(i)
                qml.Hadamard(i)
            # repeated component
            Oracle_Diffusion(cycles)
            # measurement
            return qml.probs(wires=range(self.qubits))

        # the output of the method is a 1D tensor of probabilities for each bit string given the number of cycles specified
        return G(cycles)

    def cost(self, x):
        return -self.run(x)[self.num]

    def opt(self):
        # initialize optimizer using integer step size as cycles must always be an integer
        opt = qml.GradientDescentOptimizer(stepsize=1)
        # initialize training parameter
        theta = np.array(0.0, requires_grad=True)
        # cost function to be trained
        # sends theta -> run method and selects the amplitude of marked state from output tensor
        # value is negative to permit minimization
        Cost_vec = [-self.run(theta)[self.num]]
        angle = [theta]

        max_iterations = 100
        convergence_tolerance = 1e-06

        for n in range(max_iterations):
            # cost function
            # BROKEN LINE - step and cost does not appear to be working with this setup
            theta, prev_prob = opt.step_and_cost(cost, theta)
            # store subsequent optimization steps in the corresponding vectors
            Cost_vec.append(self.run(theta)[self.num])
            angle.append(theta)

NB:

        if cycles is None:
            cycles = self.cycles
        # the run method takes in an integer to define the number of cycles - during optimization a tensor is passed to run - this converts the tensor input to one the method can use
        elif type(cycles) != int:

            cycles = int(np.floor(cycles))
9

That did something!

Now opt() runs, but it seems not to be passing new theta values to opt.step
I also tried running step_and_cost and got the same error

Here is my updated code:

#Defining an algorithm class for Grover search
#This class will have addded functionality to optimize the number of calls to the oracle to experimentally verify the optimal value of sqrt(n)

class Grover:
#------------

INITIALIZATION, Input number of qubits, shots, and the bit string to be marked for search

#------------
def init(self, qubits, shots, oracle_state):
#num of qubits
self.qubits = qubits
#num of shots on backend
self.shots = shots

    #marked state defining the oracle - input in the form of an integer
    #integer is converted to bit string which is converted to binary vector for later use in FlipState method
    self.num = oracle_state
    if oracle_state >= 2**self.qubits:
        raise ValueError(f"Insufficient qubits. Input value less than {2**self.qubits}")
    v = []
    B = bin(oracle_state)[2:].zfill(self.qubits)
    for i in range(len(B)):
        v.append(int(B[i]))
    self.oracle = v
    
    #optimal number of oracle+diffuser cycles - this value is used as the default to reproduce the standard Grover search
    self.cycles = int(np.floor(np.sqrt(self.qubits)))

#------------

RUN METHOD, This is the core method which is used as input to all other methods

Default run uses the optimal number of oracle cycles

User can specify alternate number of cycles if desired

#------------
def run(self, cycles = None):
#cycles is an optional input if the user wants to implement Grover with a nonstandard number of oracle calls
if cycles is None:
cycles = self.cycles
#the run method takes in an integer to define the number of cycles - during optimization a tensor is passed to run - this converts the tensor input to one the method can use
#if type(cycles) == np.tensor:
# cycles = int(np.floor(cycles.item()))
elif type(cycles) != int:
cycles = int(np.floor(cycles))

    #setup for Grover algorithm with initial state [1,1,...,1] (this is more optimal as the diffuser does not need to be conjugated by x gates)
    initial_state = []
    for i in range(self.qubits):
        initial_state.append(1)
        
    #oracle+diffusion
    def OD():
    #oracle
        qml.FlipSign(self.oracle, wires = range(self.qubits))
    #diffusion
        for i in range(self.qubits):
            qml.Hadamard(i)
        qml.FlipSign(initial_state, wires = range(self.qubits))
        for i in range(self.qubits):
            qml.Hadamard(i)
            
    #repeater oracle+diffuser 
    def Oracle_Diffusion(n):
        if type(n) == np.tensor:
            n = int(np.floor(n.item()))
        for i in range(int(n)):
            OD()
            
#run grover
    dev = qml.device("default.qubit", wires = self.qubits, shots = self.shots)
    @qml.qnode(dev, interface="autograd")
    def G(x):
        if type(x) == np.tensor:
            x = int(np.floor(x.item()))
    #initial state prep
        for i in range(self.qubits):
            qml.PauliX(i)
            qml.Hadamard(i)
    #repeated component
        Oracle_Diffusion(x)
    #measurement
        return qml.probs(wires = range(self.qubits))
    #the output of the method is a 1D tensor of probabilities for each bit string given the number of cycles specified
    return (G(cycles))

#------------

PLOT METHOD, visualizes probabilities outputted by run method

Default run uses the optimal number of oracle cycles

User can specify alternate number of cycles if desired

#------------

def plot(self, cycles = None):
    #same as for run
    if cycles is None:
        cycles = self.cycles
    if type(cycles) == np.tensor:
        cycles = int(np.floor(cycles.item()))
        
    #store y values as output of run method
    y = self.run(cycles)
    bit_strings = [f"{x:0{self.qubits}b}" for x in range(len(y))]
    plt.bar(bit_strings, y, color = "#212121")
    plt.xticks(rotation="vertical")
    plt.xlabel("State label")
    plt.ylabel("Probability Amplitude")
    plt.title("States probabilities amplitudes")
    plt.show()

#------------

OPT METHOD, utilizes cycle input capability to experimentally recover optimal number of cycles

#------------
def opt(self):
#initialize optimizer using integer step size as cycles must always be an integer
opt = qml.GradientDescentOptimizer(stepsize=1)
#initialize training parameter
#dev = qml.device(“default.qubit”, wires = self.qubits, shots = self.shots)
#@qml.qnode(dev, interface=“autograd”)
def cost(x):
return -self.run(x)[self.num]
theta = np.array(0.0, requires_grad=True)
#cost function to be trained
#sends theta → run method and selects the amplitude of marked state from output tensor
#value is negative to permit minimization
#Cost_vec = [-self.run(theta)[self.num]]
Cost_vec = [cost(theta)]
angle = [theta]

    max_iterations = 10
    #convergence_tolerance = 1e-06

    for n in range(max_iterations):
        #BROKEN LINE - step and cost does not appear to be working with this settup
        #theta = opt.step(cost, theta)
        theta, prev_cost = opt.step_and_cost(cost, theta)
        
        #store subsequent optimization steps in the corresponding vectors
        Cost_vec.append(cost(theta))
        angle.append(theta)
        
        #convergence condition left over from VQE demo where opt code was first sourced - doesn't apply to this problem
        #conv = np.abs(Prob[-1] - prev_prob)
          
        print(f"Step = {n}, probability = {Cost_vec[-1]:.8f}")
          
        #if conv <= convergence_tolerance:
        #    break

    print("\n" f"Final amplitude of the target state = {Cost_vec[-1]:.8f}")
    print("\n" f"Optimal value of the circuit parameter = {angle[-1]:.4f}")

Here is what I expect from ten iterations compared to the method output

image
image452×532 5.9 KB

10

I think the issue is (probably two-fold) in that that you’re creating a QNode every time you call run, which is (essentially) what is being differentiated. Can you try creating the QNode once by either (1) defining it outside of the Grover class or (2) in the Grover class making sure that the QNode is created when an instance of Grover is created? Let me know if that works!

11

So this is something I played with in the early stages of Grover’s development and I couldn’t seem to get it to work without having @qml.qnode(dev) immediately preceding the circuit I wanted to run - that is why it put it inside the run method

If I put qnode anywhere else I get errors - is there another way to do this?

Is qnode created when I define dev = qml.device(“default.qubit”, wires=self.qubits, shots = self.shots)?

or when I write @qml.node(dev)?

If it is the former, then I can store dev as self.dev in the following way (this is a simpler bit of code code designed to implement a superposition circuit within this OOP framework I am learning)

class Sup:

def __init__(self, qubits, shots):
    self.qubits = qubits
    self.shots = shots
    self.dev = qml.device("default.qubit", wires=self.qubits, shots = self.shots)
    
def run(self):
    @qml.qnode(self.dev, interface="autograd")
    def H(q):
        for i in range(self.qubits):
            qml.Hadamard(i)
        return [qml.expval(qml.PauliZ(i)) for i in range(self.qubits)]
    return (H(self.qubits))

But beyond this I’m not sure how to create a node without placing it immediately above the circuit I am passing to it

12

Is qnode created when I define dev = qml.device(“default.qubit”, wires=self.qubits, shots = self.shots)? or when I write @qml.node(dev)?

A QNode is created when you decorate a quantum function with @qml.qnode. You should be able to put it outside the class and have things work (option 1) or as a class attribute (option 2) :smile:

13

I think my issue is I dont really understand Qnode - I have tried placing it as a class attribute in the following example code

class Sup:

def __init__(self, qubits, shots):
    self.qubits = qubits
    self.shots = shots
    self.dev = qml.device("default.qubit", wires=self.qubits, shots = self.shots)
    self.node = qml.qnode(self.dev)
    
def run(self):
    #wrapper
    @self.node
    #circuit
    def H(q):
        for i in range(self.qubits):
            qml.Hadamard(i)
        return [qml.expval(qml.PauliZ(i)) for i in range(self.qubits)]
        #return qml.sample()
    return (H(self.qubits))

But I don’t think this is functionally any different from what I had before.

I tried removing it from the class and placing it outside
A = sup(3,100) #no device or qnode inside the class
dev = qml.device(“default.qubit”, wires=A.qubits, shots = A.shots)
@qml.node(dev)
A.run()

but this gives ‘invalid syntax’ error

I have tried several other placements of qnode all with the same effect

Can you give an example of how the qnode can be made a class attribute such that the run method still works?

I was playing with node = qml.QNode(circuit, dev) also but havent had any luck

14

Sure thing! Here’s a minimal example:

import pennylane as qml
from pennylane import numpy as np

class MyClass:
    def __init__(self, n_qubits):
        self.n_qubits = n_qubits
        device = qml.device("default.qubit", wires=n_qubits)

        self.quantum_node = qml.QNode(self.quantum_function, device)

    def quantum_function(self, params):
        for i in range(len(params)):
            qml.RX(params[i], wires=i)
        return [qml.expval(qml.PauliZ(i)) for i in range(len(params))]

    def cost(self, params):
        return -np.sum(self.quantum_node(params))

    def do_optimization(self):
        params = np.random.uniform(0, np.pi, size=(self.n_qubits,), requires_grad=True)
        opt = qml.GradientDescentOptimizer(0.1)
        for n in range(10):
            params = opt.step(self.cost, params)
            print(self.cost(params))

obj = MyClass(4)
obj.do_optimization()

-3.6838202359463894
-3.737987269714827
-3.783698057758141
-3.8220120753125184
-3.8539404434260787
-3.880416877392254
-3.9022816835753353
-3.920275812324963
-3.9350420219145708
-3.9471306216104813

Hope this helps!

1 Like
15

Fantastic thank you!

1 Like
16

My pleasure! Glad I could help :slight_smile: