I tried to make a cirq program calculating an eigenvalue of the observable by VQE. Inspired by the qulacs VQE tutorial, I defined a cost function from the expectation value of the observable with hardware efficient ansatz.

import cirq
import sympy
import numpy as np

def hea_multilayer(qubits, depth):
    n_qubits = len(qubits)
    n_params = 2*n_qubits*(depth + 1)
    theta = sympy.symbols(f'a:{n_params}')
    for j in range(depth):
        for i in range(n_qubits):
            yield cirq.ry(theta[2*j*n_qubits + i])(qubits[i])
            yield cirq.rz(theta[2*j*n_qubits + i + n_qubits])(qubits[i])

        for i in range(n_qubits - 1):
            yield cirq.CNOT(qubits[i], qubits[i+1])

    for i in range(n_qubits):
        yield cirq.ry(theta[2*depth*n_qubits + i])(qubits[i])
        yield cirq.rz(theta[2*depth*n_qubits + n_qubits + i])(qubits[i])

n_qubits = 2
depth = 1
qubits = cirq.LineQubit.range(n_qubits)
qc = cirq.Circuit()
qc.append(hea_multilayer(qubits, depth))

observable = cirq.PauliSum.from_pauli_strings([
    cirq.PauliString(1.0, cirq.I(qubits[0]), cirq.I(qubits[1])),
    cirq.PauliString(1.0, cirq.X(qubits[0]), cirq.X(qubits[1])),
    cirq.PauliString(-1.0, cirq.Y(qubits[0]), cirq.Y(qubits[1])),
    cirq.PauliString(1.0, cirq.Z(qubits[0]), cirq.Z(qubits[1]))

sim = cirq.Simulator()
def cost(parameter):
    subs = {}
    for i in range(2*n_qubits*(depth+1)):
        subs[f'a{i}'] = parameter[i]

    result = sim.simulate_expectation_values(qc, observable, subs)
    return result[0].real

After that, I minimized this cost function using scipy minimize function.

parameter = np.random.random(8)

from scipy.optimize import minimize
cost_history = []
min_result = minimize(cost, parameter, method="BFGS", callback=lambda x: cost_history.append(cost(x)))

However, this minimization did not work with the following message:

  message: Desired error not necessarily achieved due to precision loss.
  success: False
   status: 2
      fun: 0.2531989514827728
        x: [ 3.220e-01  8.904e-01  5.881e-01 -1.969e-01 -6.673e-01
             6.296e-01  2.209e-02  2.421e-01]
      nit: 7
      jac: [ 0.000e+00  0.000e+00  0.000e+00  0.000e+00  0.000e+00
             2.000e+00  0.000e+00  0.000e+00]
 hess_inv: [[ 1.000e+00  0.000e+00 ...  0.000e+00  0.000e+00]
            [ 0.000e+00  1.000e+00 ...  0.000e+00  0.000e+00]
            [ 0.000e+00  0.000e+00 ...  2.805e-06 -5.530e-06]
            [ 0.000e+00  0.000e+00 ... -5.530e-06  1.104e-05]]
     nfev: 993
     njev: 109
[0.36777588725090027, 0.25321047008037567, 0.2532104179263115, 0.25321032106876373, 0.2532094195485115, 0.25320935994386673, 0.2532089278101921, 0.2531989514827728]

I implemented this program using qiskit, and compared the result to the qiskit VQE solver:

from qiskit.primitives import Estimator
from qiskit.quantum_info import SparsePauliOp
from qiskit.circuit.library import EfficientSU2
import numpy as np


qc = EfficientSU2(2, reps = 1, insert_barriers=True, flatten=True)

observable = SparsePauliOp(["II", "XX", "YY", "ZZ"], coeffs = [1, 1, -1, 1])

estimator = Estimator()

parameter = np.random.random(8)
job = estimator.run(qc, observable, parameter_values = parameter)
print(f'Parameter: {parameter}')
print(f'Job result: {job.result()}')

def cost(parameter):
    job = estimator.run(qc, observable, parameter_values=parameter)
    return job.result().values[0]

from scipy.optimize import minimize
cost_history = []
min_result = minimize(cost, parameter, method="L-BFGS-B", callback=lambda x: cost_history.append(cost(x)))

from qiskit.algorithms.minimum_eigensolvers import VQE
from qiskit.algorithms.optimizers import L_BFGS_B

cost_history = []
vqe = VQE(estimator=estimator, ansatz=qc, optimizer=L_BFGS_B(), initial_point=parameter)
result = vqe.compute_minimum_eigenvalue(observable)

And my result is similar to the qiskit VQE solver result.

Job result: EstimatorResult(values=array([1.43214656]), metadata=[{}])
  success: True
   status: 0
      fun: 5.907418998418734e-11
        x: [ 2.681e-01  1.710e+00  8.852e-01  1.201e-02 -1.934e-01
             1.231e+00  3.499e-01  1.055e+00]
      nit: 6
      jac: [-2.187e-06 -9.281e-06 -2.776e-06 -4.619e-06  6.217e-06
            -9.692e-06 -1.099e-06 -1.099e-06]
     nfev: 63
     njev: 7
 hess_inv: <8x8 LbfgsInvHessProduct with dtype=float64>
[1.432146558079639, 0.2513704685300491, 0.15156282475910865, 0.0019124923978873776, 3.6796306554398583e-05, 1.1992205672939349e-08, 5.907418998418734e-11]
/home/wleelinux/sources/qiskit_expect/hea_estimator.py:31: DeprecationWarning: ``qiskit.algorithms`` has been migrated to an independent package: https://github.com/qiskit-community/qiskit-algorithms. The ``qiskit.algorithms`` import path is deprecated as of qiskit-terra 0.25.0 and will be removed no earlier than 3 months after the release date. Please run ``pip install qiskit_algorithms`` and use ``import qiskit_algorithms`` instead.
  from qiskit.algorithms.minimum_eigensolvers import VQE
{   'aux_operators_evaluated': None,
    'cost_function_evals': 63,
    'eigenvalue': 2.642741581126984e-11,
    'optimal_circuit': <qiskit.circuit.library.n_local.efficient_su2.EfficientSU2 object at 0x7ff1ac903590>,
    'optimal_parameters': {   ParameterVectorElement(θ[0]): -0.08608117349727416,
                              ParameterVectorElement(θ[1]): 1.9438040339495675,
                              ParameterVectorElement(θ[2]): 1.1497492983826083,
                              ParameterVectorElement(θ[3]): -0.1815167532084752,
                              ParameterVectorElement(θ[4]): -0.23455036012933897,
                              ParameterVectorElement(θ[5]): 1.413331393301791,
                              ParameterVectorElement(θ[6]): 0.7537967436651564,
                              ParameterVectorElement(θ[7]): 1.4589817775219438},
    'optimal_point': array([-0.08608117,  1.94380403,  1.1497493 , -0.18151675, -0.23455036,
        1.41333139,  0.75379674,  1.45898178]),
    'optimal_value': 2.642741581126984e-11,
    'optimizer_evals': None,
    'optimizer_result': <qiskit.algorithms.optimizers.optimizer.OptimizerResult object at 0x7ff1a92233d0>,
    'optimizer_time': 0.2594926357269287}

Why my cirq code failed? Is it wrong to use simulate_expectation_values method and should I implement the expectation value from the measure of the quantum circuit in cirq?

  • $\begingroup$ Inspired another cirq vqe implementation (github.com/mafaldaramoa/VQE/blob/main/CIRQ_VQE.ipynb), I changed the minimization method to Nelder-Mead, and the optimization terminated successfully. The minimum was -1.2246891856193542e-07 at x: [ 9.571e-03 2.241e+00 8.869e-01 1.190e-01 1.241e-01 9.869e-01 -4.678e-03 7.675e-01]. I think that the calculation of derivative in cirq would have an error. $\endgroup$ Oct 13, 2023 at 2:17


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Browse other questions tagged or ask your own question.