VQE with ibmq machine

I’m trying to do VQE with the the qistkit.ibmq device.
I followed the video here:

and changed the quantum simulator to qiskit.ibmq
but i met some problems, but I’m not sure what results in it
is it the step_and_cost function?
now I excluded the max_iteration possibility by executing only once

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
# Put code here
import pennylane as qml
from pennylane import numpy as np
from pennylane import qchem

symbols = ["H", "H", "H"]
coordinates = np.array([[0.0102, 0.0442, 0.0], [0.9867, 1.6303, 0.0], [1.8720, -0.0085, 0.0]])

hamiltonian, qubits = qchem.molecular_hamiltonian(symbols, coordinates, charge = 1)

hf = qchem.hf_state(electrons = 2, orbitals = 6)

token_str = ""
with open('token.txt', 'r') as file:
    token_str = file.read().replace('\n', '')

num_wires = qubits
#dev = qml.device("default.qubit", wires = num_wires)
#dev = qml.device("qiskit.ibmq", wires = num_wires, backend="ibm_kyoto", ibmqx_token = token_str)
dev = qml.device("qiskit.ibmq", wires = num_wires, backend="ibm_kyoto")

@qml.qnode(dev)
def exp_energy(state):
    qml.BasisState(np.array(state), wires = range(num_wires))
    return qml.expval(hamiltonian)

def ansatz(params):
    qml.BasisState(hf, wires = range(num_wires))
    qml.DoubleExcitation(params[0], wires = [0, 1, 2, 3])
    qml.DoubleExcitation(params[1], wires = [0, 1, 4, 5])

@qml.qnode(dev)
def cost_function(params):
    ansatz(params)
    return qml.expval(hamiltonian)

opt = qml.GradientDescentOptimizer(stepsize = 0.4)
theta = np.array([0.0, 0.0], requires_grad = True)

energy = [cost_function(theta)]
angle = [theta]
max_iterations = 20

#for n in range(max_iterations):
#    print(len(angle))  
#    theta, prev_energy = opt.step_and_cost(cost_function, theta)
#

#    energy.append(cost_function(theta))
#    angle.append(theta)

#    if n % 2 == 0:
#        print(f"Step = {n}, Energy = {energy[-1]:.8f} Ha")

print(len(angle))  
theta, prev_energy = opt.step_and_cost(cost_function, theta)


energy.append(cost_function(theta))
angle.append(theta)
print(f"Step = {n}, Energy = {energy[-1]:.8f} Ha")

If you want help with diagnosing an error, please put the full error message below:

# Put full error message here

---------------------------------------------------------------------------
IBMBackendValueError                      Traceback (most recent call last)
Cell In[13], line 20
      8 #for n in range(max_iterations):
      9 #    print(len(angle))  
     10 #    theta, prev_energy = opt.step_and_cost(cost_function, theta)
   (...)
     16 #    if n % 2 == 0:
     17 #        print(f"Step = {n}, Energy = {energy[-1]:.8f} Ha")
     19 print(len(angle))  
---> 20 theta, prev_energy = opt.step_and_cost(cost_function, theta)
     23 energy.append(cost_function(theta))
     24 angle.append(theta)

File ~/.cache/pypoetry/virtualenvs/q1-LJbbFjuQ-py3.12/lib/python3.12/site-packages/pennylane/optimize/gradient_descent.py:64, in GradientDescentOptimizer.step_and_cost(self, objective_fn, grad_fn, *args, **kwargs)
     44 def step_and_cost(self, objective_fn, *args, grad_fn=None, **kwargs):
     45     """Update trainable arguments with one step of the optimizer and return the corresponding
     46     objective function value prior to the step.
     47 
   (...)
     61         If single arg is provided, list [array] is replaced by array.
     62     """
---> 64     g, forward = self.compute_grad(objective_fn, args, kwargs, grad_fn=grad_fn)
     65     new_args = self.apply_grad(g, args)
     67     if forward is None:

File ~/.cache/pypoetry/virtualenvs/q1-LJbbFjuQ-py3.12/lib/python3.12/site-packages/pennylane/optimize/gradient_descent.py:122, in GradientDescentOptimizer.compute_grad(objective_fn, args, kwargs, grad_fn)
    104 r"""Compute gradient of the objective function at the given point and return it along with
    105 the objective function forward pass (if available).
    106 
   (...)
    119     will not be evaluted and instead ``None`` will be returned.
    120 """
    121 g = get_gradient(objective_fn) if grad_fn is None else grad_fn
--> 122 grad = g(*args, **kwargs)
    123 forward = getattr(g, "forward", None)
    125 num_trainable_args = sum(getattr(arg, "requires_grad", False) for arg in args)

File ~/.cache/pypoetry/virtualenvs/q1-LJbbFjuQ-py3.12/lib/python3.12/site-packages/pennylane/_grad.py:165, in grad.__call__(self, *args, **kwargs)
    162     self._forward = self._fun(*args, **kwargs)
    163     return ()
--> 165 grad_value, ans = grad_fn(*args, **kwargs)  # pylint: disable=not-callable
    166 self._forward = ans
    168 return grad_value

File ~/.cache/pypoetry/virtualenvs/q1-LJbbFjuQ-py3.12/lib/python3.12/site-packages/autograd/wrap_util.py:20, in unary_to_nary.<locals>.nary_operator.<locals>.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)

File ~/.cache/pypoetry/virtualenvs/q1-LJbbFjuQ-py3.12/lib/python3.12/site-packages/pennylane/_grad.py:191, in grad._grad_with_forward(fun, x)
    185 if vspace(ans).size != 1:
    186     raise TypeError(
    187         "Grad only applies to real scalar-output functions. "
    188         "Try jacobian, elementwise_grad or holomorphic_grad."
    189     )
--> 191 grad_value = vjp(vspace(ans).ones())
    192 return grad_value, ans

File ~/.cache/pypoetry/virtualenvs/q1-LJbbFjuQ-py3.12/lib/python3.12/site-packages/autograd/core.py:14, in make_vjp.<locals>.vjp(g)
---> 14 def vjp(g): return backward_pass(g, end_node)

File ~/.cache/pypoetry/virtualenvs/q1-LJbbFjuQ-py3.12/lib/python3.12/site-packages/autograd/core.py:21, in backward_pass(g, end_node)
     19 for node in toposort(end_node):
     20     outgrad = outgrads.pop(node)
---> 21     ingrads = node.vjp(outgrad[0])
     22     for parent, ingrad in zip(node.parents, ingrads):
     23         outgrads[parent] = add_outgrads(outgrads.get(parent), ingrad)

File ~/.cache/pypoetry/virtualenvs/q1-LJbbFjuQ-py3.12/lib/python3.12/site-packages/autograd/core.py:67, in defvjp.<locals>.vjp_argnums.<locals>.<lambda>(g)
     64         raise NotImplementedError(
     65             "VJP of {} wrt argnum 0 not defined".format(fun.__name__))
     66     vjp = vjpfun(ans, *args, **kwargs)
---> 67     return lambda g: (vjp(g),)
     68 elif L == 2:
     69     argnum_0, argnum_1 = argnums

File ~/.cache/pypoetry/virtualenvs/q1-LJbbFjuQ-py3.12/lib/python3.12/site-packages/pennylane/interfaces/autograd.py:199, in vjp.<locals>.grad_fn(dy)
    196 def grad_fn(dy):
    197     """Returns the vector-Jacobian product with given
    198     parameter values and output gradient dy"""
--> 199     vjps = jpc.compute_vjp(tapes, dy)
    200     return tuple(
    201         qml.math.to_numpy(v, max_depth=1) if isinstance(v, ArrayBox) else v for v in vjps
    202     )

File ~/.cache/pypoetry/virtualenvs/q1-LJbbFjuQ-py3.12/lib/python3.12/site-packages/pennylane/interfaces/jacobian_products.py:293, in TransformJacobianProducts.compute_vjp(self, tapes, dy)
    290     logger.debug("compute_vjp called with (%s, %s)", tapes, dy)
    292 if self._cache_full_jacobian:
--> 293     jacs = self.compute_jacobian(tapes)
    294     return _compute_vjps(jacs, dy, tapes)
    296 vjp_tapes, processing_fn = qml.gradients.batch_vjp(
    297     tapes, dy, self._gradient_transform, gradient_kwargs=self._gradient_kwargs
    298 )

File ~/.cache/pypoetry/virtualenvs/q1-LJbbFjuQ-py3.12/lib/python3.12/site-packages/pennylane/interfaces/jacobian_products.py:330, in TransformJacobianProducts.compute_jacobian(self, tapes)
    326 partial_gradient_fn = partial(self._gradient_transform, **self._gradient_kwargs)
    327 jac_tapes, batch_post_processing = qml.transforms.map_batch_transform(
    328     partial_gradient_fn, tapes
    329 )
--> 330 results = self._inner_execute(jac_tapes)
    331 jacs = tuple(batch_post_processing(results))
    332 self._cache[tapes] = jacs

File ~/.cache/pypoetry/virtualenvs/q1-LJbbFjuQ-py3.12/lib/python3.12/site-packages/pennylane/interfaces/execution.py:261, in _make_inner_execute.<locals>.inner_execute(tapes, **_)
    259 if numpy_only:
    260     tapes = tuple(qml.transforms.convert_to_numpy_parameters(t) for t in tapes)
--> 261 return cached_device_execution(tapes)

File ~/.cache/pypoetry/virtualenvs/q1-LJbbFjuQ-py3.12/lib/python3.12/site-packages/pennylane/interfaces/execution.py:383, in cache_execute.<locals>.wrapper(tapes, **kwargs)
    378         return (res, []) if return_tuple else res
    380 else:
    381     # execute all unique tapes that do not exist in the cache
    382     # convert to list as new device interface returns a tuple
--> 383     res = list(fn(tuple(execution_tapes.values()), **kwargs))
    385 final_res = []
    387 for i, tape in enumerate(tapes):

File /usr/lib/python3.12/contextlib.py:81, in ContextDecorator.__call__.<locals>.inner(*args, **kwds)
     78 @wraps(func)
     79 def inner(*args, **kwds):
     80     with self._recreate_cm():
---> 81         return func(*args, **kwds)

File ~/.cache/pypoetry/virtualenvs/q1-LJbbFjuQ-py3.12/lib/python3.12/site-packages/pennylane_qiskit/ibmq.py:91, in IBMQDevice.batch_execute(self, circuits)
     90 def batch_execute(self, circuits):  # pragma: no cover, pylint:disable=arguments-differ
---> 91     res = super().batch_execute(circuits, timeout=self.timeout_secs)
     92     if self.tracker.active:
     93         self._track_run()

File ~/.cache/pypoetry/virtualenvs/q1-LJbbFjuQ-py3.12/lib/python3.12/site-packages/pennylane_qiskit/qiskit_device.py:498, in QiskitDevice.batch_execute(self, circuits, timeout)
    495 compiled_circuits = self.compile_circuits(circuits)
    497 # Send the batch of circuit objects using backend.run
--> 498 self._current_job = self.backend.run(compiled_circuits, shots=self.shots, **self.run_args)
    500 try:
    501     result = self._current_job.result(timeout=timeout)

File ~/.cache/pypoetry/virtualenvs/q1-LJbbFjuQ-py3.12/lib/python3.12/site-packages/qiskit_ibm_provider/ibm_backend.py:423, in IBMBackend.run(self, circuits, dynamic, job_tags, init_circuit, init_num_resets, header, shots, memory, meas_level, meas_return, rep_delay, init_qubits, use_measure_esp, noise_model, seed_simulator, **run_config)
    421 if not isinstance(circuits, List):
    422     circuits = [circuits]
--> 423 self._check_circuits_attributes(circuits)
    425 if (
    426     use_measure_esp
    427     and getattr(self.configuration(), "measure_esp_enabled", False) is False
    428 ):
    429     raise IBMBackendValueError(
    430         "ESP readout not supported on this device. Please make sure the flag "
    431         "'use_measure_esp' is unset or set to 'False'."
    432     )

File ~/.cache/pypoetry/virtualenvs/q1-LJbbFjuQ-py3.12/lib/python3.12/site-packages/qiskit_ibm_provider/ibm_backend.py:825, in IBMBackend._check_circuits_attributes(self, circuits)
    819 """Check that circuits can be executed on backend.
    820 Raises:
    821     IBMBackendValueError:
    822         - If one of the circuits contains more qubits than on the backend."""
    824 if len(circuits) > self._max_circuits:
--> 825     raise IBMBackendValueError(
    826         f"Number of circuits, {len(circuits)} exceeds the "
    827         f"maximum for this backend, {self._max_circuits})"
    828     )
    829 for circ in circuits:
    830     if isinstance(circ, QuantumCircuit):

IBMBackendValueError: 'Number of circuits, 804 exceeds the maximum for this backend, 300)'

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

Name: PennyLane
Version: 0.34.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: /home/usagi/.cache/pypoetry/virtualenvs/q1-LJbbFjuQ-py3.12/lib/python3.12/site-packages
Requires: appdirs, autograd, autoray, cachetools, networkx, numpy, pennylane-lightning, requests, rustworkx, scipy, semantic-version, toml, typing-extensions
Required-by: PennyLane-Lightning, PennyLane-qiskit, q1

Platform info: Linux-5.15.133.1-microsoft-standard-WSL2-x86_64-with-glibc2.35
Python version: 3.12.1
Numpy version: 1.26.3
Scipy version: 1.11.4
Installed devices:

  • default.gaussian (PennyLane-0.34.0)
  • default.mixed (PennyLane-0.34.0)
  • default.qubit (PennyLane-0.34.0)
  • default.qubit.autograd (PennyLane-0.34.0)
  • default.qubit.jax (PennyLane-0.34.0)
  • default.qubit.legacy (PennyLane-0.34.0)
  • default.qubit.tf (PennyLane-0.34.0)
  • default.qubit.torch (PennyLane-0.34.0)
  • default.qutrit (PennyLane-0.34.0)
  • null.qubit (PennyLane-0.34.0)
  • lightning.qubit (PennyLane-Lightning-0.34.0)
  • qiskit.aer (PennyLane-qiskit-0.34.0)
  • qiskit.basicaer (PennyLane-qiskit-0.34.0)
  • qiskit.ibmq (PennyLane-qiskit-0.34.0)
  • qiskit.ibmq.circuit_runner (PennyLane-qiskit-0.34.0)
  • qiskit.ibmq.sampler (PennyLane-qiskit-0.34.0)
  • qiskit.remote (PennyLane-qiskit-0.34.0)

Hello,
the main problem is, you need to run 814 circuits, and the limit is 300.
I’ve tried to split the circuits into multiple jobs, but it doesn’t work.
I’ve tried to reduce the circuit number by running more of the required gates in a circuit, but that doesn’t help either.
It is just that this code needs many circuits to work, for the given molecule.(more than maximum 300).
I recommend you run the code on a local simulator, since it doesn’t have a circuit limit.

I also came up with this sample, give it a try. It uses code from the tutorial “A brief overview of VQE”, from the quantum chemistry tutorials in Pennylane.
It runs for h2, but also exceeds the limit for h3+. (it takes 318 circuits)
But you can do well with a smaller molecule.
Lastly, if you want to execute on a real quantum computer (ibm_kyoto), it will take about 20 seconds per job, so beware not to expend all your time.

from jax import numpy as np
import jax
jax.config.update("jax_platform_name", "cpu")
jax.config.update('jax_enable_x64', True)
!pip install --upgrade pennylane
!pip install pennylane-qiskit
import pennylane as qml
#choose a molecule from the dataset
dataset = qml.data.load('qchem', molname="H2")[0]
H, qubits = dataset.hamiltonian, len(dataset.hamiltonian.wires)

print("Number of qubits = ", qubits)
print("The Hamiltonian is ", H)
#choose a backend
#dev = qml.device("default.qubit", wires=qubits)
#or
#dev = qml.device('qiskit.ibmq', wires=qubits, backend='ibmq_qasm_simulator', ibmqx_token="")
electrons = 2
hf = qml.qchem.hf_state(electrons, qubits)
print(hf)
@qml.qnode(dev)
def circuit(param, wires):
    qml.BasisState(hf, wires=wires)
    qml.DoubleExcitation(param, wires=[0, 1, 2, 3])
    return qml.expval(H)
def cost_fn(param):
    return circuit(param, wires=range(qubits))
import optax

max_iterations = 20
conv_tol = 1e-06

opt = optax.sgd(learning_rate=0.4)
theta = np.array(0.)

# store the values of the cost function
energy = [cost_fn(theta)]

# store the values of the circuit parameter
angle = [theta]

opt_state = opt.init(theta)

for n in range(max_iterations):

    gradient = jax.grad(cost_fn)(theta)
    updates, opt_state = opt.update(gradient, opt_state)
    theta = optax.apply_updates(theta, updates)

    angle.append(theta)
    energy.append(cost_fn(theta))

    conv = np.abs(energy[-1] - energy[-2])

    if n % 2 == 0:
        print(f"Step = {n},  Energy = {energy[-1]:.8f} Ha")

    if conv <= conv_tol:
        break

print("\n" f"Final value of the ground-state energy = {energy[-1]:.8f} Ha")
print("\n" f"Optimal value of the circuit parameter = {angle[-1]:.4f}")
import matplotlib.pyplot as plt

fig = plt.figure()
fig.set_figheight(5)
fig.set_figwidth(12)

# Full configuration interaction (FCI) energy computed classically
E_fci = -1.136189454088

# Add energy plot on column 1
ax1 = fig.add_subplot(121)
ax1.plot(range(n + 2), energy, "go", ls="dashed")
ax1.plot(range(n + 2), np.full(n + 2, E_fci), color="red")
ax1.set_xlabel("Optimization step", fontsize=13)
ax1.set_ylabel("Energy (Hartree)", fontsize=13)
ax1.text(0.5, -1.1176, r"$E_\mathrm{HF}$", fontsize=15)
ax1.text(0, -1.1357, r"$E_\mathrm{FCI}$", fontsize=15)
plt.xticks(fontsize=12)
plt.yticks(fontsize=12)

# Add angle plot on column 2
ax2 = fig.add_subplot(122)
ax2.plot(range(n + 2), angle, "go", ls="dashed")
ax2.set_xlabel("Optimization step", fontsize=13)
ax2.set_ylabel("Gate parameter $\\theta$ (rad)", fontsize=13)
plt.xticks(fontsize=12)
plt.yticks(fontsize=12)

plt.subplots_adjust(wspace=0.3, bottom=0.2)
plt.show()

Hey @chocolatesouffle0609! Welcome to the forum :slight_smile:

Was @Alex’s response helpful?

Hi @isaacdevlugt!
I think I can improve my code. Could you help me, please?
So far I read the demo “Measurement optimization” from Josh Izaac.
I wish we could simulate H2O on a real backend, using measurement optimization. The issue is, for the H2O molecule, the optimization is still using more than 300 circuits/job . In the demo there are references 2,3,4 pointing to techniques that further reduce this number.
Are there any code examples using full commutativity [2], unitary partitioning [3], and Fermionic basis rotation grouping[4]? These are needed so we can have below 300 circuits per job.

Hey @Alex,

This is probably a little too much detail for what I’m able to offer on my side :sweat_smile:. Sounds like there’s some research to do but that’s life! Let us know if you have any troubles with PennyLane along the way :slight_smile:

1 Like