# Variational Quantum Eigensolver (VQE) by scipy minimize in cirq does not work

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]))
])
print(observable)
print(qc)

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.

np.random.seed(2023)
parameter = np.random.random(8)

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


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

np.random.seed(2023)

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'Circuit:\n{qc}')
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 = []
cost_history.append(cost(parameter))
min_result = minimize(cost, parameter, method="L-BFGS-B", callback=lambda x: cost_history.append(cost(x)))
print(min_result)
print(cost_history)

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

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


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?

• 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. Oct 13, 2023 at 2:17