VQE Application in Quantum Chemistry Computing
Overview
Quantum chemistry refers to solving the numerical values of the time-dependent or time-independent Schrödinger equations by using the basic theory and method of quantum mechanics. Quantum chemical simulation on high-performance computers has become an important method to study the physical and chemical properties of materials. However, the exact solution of the Schrödinger equation has exponential complexity, which severely constrains the scale of the chemical system that can be simulated. The development of quantum computing in recent years provides a feasible way to solve this problem. It is expected that the Schrödinger equation can be solved with high accuracy on quantum computers under the complexity of polynomials.
Peruzzo et al. first applied the VQE and unitary coupled-cluster theory to quantum chemistry simulation in 2014 to solve the ground state energy of He-H+. The VQE is a hybrid quantum-classical algorithm and is widely used in chemical simulation based on quantum algorithms. This tutorial describes how to use the VQE to solve the ground-state energy of a molecular system.
This tutorial consists of the following parts:
Introduction to the quantum chemistry
VQE application
Using MindQuantum to perform VQE simulation with efficient and automatic derivation
This document applies to the CPU environment. You can obtain the complete executable sample code at https://gitee.com/mindspore/mindquantum/blob/r0.6/tutorials/source/7.vqe_for_quantum_chemistry.py.
Environment Preparation
In this tutorial, the following environments need to be installed:
NumPy
SciPy
PySCF
OpenFermion
OpenFermion-PySCF
The preceding dependencies can be installed by running the
pip
command.
Importing Dependencies
Import the modules on which this tutorial depends.
import numpy as np
from openfermion.chem import MolecularData
from openfermionpyscf import run_pyscf
import mindquantum as mq
from mindquantum import Hamiltonian
from mindquantum.core import X, RX
from mindquantum.core import Circuit
from mindquantum.algorithm.nisq.chem import generate_uccsd
import mindspore as ms
import mindspore.context as context
from mindspore import Parameter
from mindspore.common.initializer import initializer
context.set_context(mode=context.PYNATIVE_MODE, device_target="CPU")
Quantum Chemistry Computing Method
The core of quantum chemistry is to solve the Schrödinger equation. In general, the solution of time-dependent Schrödinger equation is more complex, so Born-Oppenheimer approximation (BO approximation) is introduced. In BO approximation, the mass of the nucleus is far greater than that of electrons, and the velocity of the nucleus is far lower than that of electrons. Therefore, the nucleus and electrons can be separated from each other, and the time-independent electron motion equation (also called the time-independent Schrödinger equation) can be obtained as follows:
They are electron kinetic energy, electron-electron potential energy and electron-nuclear potential energy.
There are many numerical algorithms that can be used to solve the time-independent Schrödinger equation. This tutorial introduces one of these methods: the wave function. Wave function directly solves the eigenfunction and eigenenergy of a given molecular Hamiltonian. At present, there are a large number of open-source software packages, such as PySCF, which can be implemented. Here is a simple example: lithium hydride molecules, using the OpenFermion and OpenFermion-PySCF plug-ins. First, define the molecular structure:
dist = 1.5
geometry = [
["Li", [0.0, 0.0, 0.0 * dist]],
["H", [0.0, 0.0, 1.0 * dist]],
]
basis = "sto3g"
spin = 0
print("Geometry: \n", geometry)
Geometry:
[['Li', [0.0, 0.0, 0.0]], ['H', [0.0, 0.0, 1.5]]]
The code above defines a Li-H key with a length of 1.5Å molecules. The STO-3G basis set is used for computing. Then, OpenFermion-PySCF is used to call PySCF to perform Hartree-Fock (HF), coupled-cluster with singles and doubles (CCSD), and full configuration interaction (FCI) computing. These three methods belong to the wave function. Before starting the computing, first make a brief introduction to these methods.
Wave Function
One of the methods to solve the time-independent Schrödinger equation is the Hartree-Fock (HF) method, which was proposed by Hartree et al. in the 1930s and is the basic method in quantum chemistry computing. The HF method introduces a single determinant approximation, that is, a wave function of the
Where
The improvement of the HF method can be derived from the wave function expansion theorem. The wave function expansion theorem can be expressed as follows: if
You can obtain the configuration interaction (CI) method:
Second Quantization
Under the second quantization expression, the Hamiltonian of the system has the following form:
The excited-state wave function can be expressed conveniently by using a second quantization expression method:
An improvement to the CI method is the coupled-cluster theory (CC). Exponential operators are introduced to CC:
The coupled-cluster operator
The effect of electron correlation is to reduce the total energy, so the ground state energy of HF is slightly higher than that of CCSD and FCI. In addition, it is easy to find that the computing volume of FCI is much greater than that of CCSD and HF. The MolecularData
function encapsulated by OpenFermion and the run_pyscf
function encapsulated by OpenFermion-PySCF are used for demonstration.
molecule_of = MolecularData(
geometry,
basis,
multiplicity=2 * spin + 1
)
molecule_of = run_pyscf(
molecule_of,
run_scf=1,
run_ccsd=1,
run_fci=1
)
print("Hartree-Fock energy: %20.16f Ha" % (molecule_of.hf_energy))
print("CCSD energy: %20.16f Ha" % (molecule_of.ccsd_energy))
print("FCI energy: %20.16f Ha" % (molecule_of.fci_energy))
Hartree-Fock energy: -7.8633576215351200 Ha
CCSD energy: -7.8823529091527051 Ha
FCI energy: -7.8823622867987249 Ha
In the preceding example, HF, CCSD, and FCI are used to compute the total energy. If you collect statistics on the runtime, you will find that molecule_file
file (molecule_of.filename
).
molecule_of.save()
molecule_file = molecule_of.filename
print(molecule_file)
One of the major obstacles to quantum chemistry is the volume of computation. As the system size (electron number and atomic number) increases, the time required for solving the FCI wave function and ground state energy increases by about
Variational Quantum Eigensolver (VQE)
The VQE is a hybrid quantum-classical algorithm. It uses the variational principle to solve the ground state wave function. The optimization of variational parameters is carried out on the classical computer.
Variational Principle
The variational principle may be expressed in the following form:
In the preceding formula,
Initial State Preparation
The
The above formula builds a bridge from quantum chemical wave function to quantum computing:
The following code builds an HF initial state wave function corresponding to the LiH molecule. In Jordan-Wigner transformation,
hartreefock_wfn_circuit = Circuit([X.on(i) for i in range(molecule_of.n_electrons)])
print(hartreefock_wfn_circuit)
q0: ──X──
q1: ──X──
q2: ──X──
q3: ──X──
We can build a probe wave function in the following form:
Wave Function Ansatz
The coupled-cluster theory mentioned above is a very efficient wave function ansatz. To use it on a quantum computer, you need to make the following modifications:
UCC is short for unitary coupled-cluster theory.
The generate_uccsd
function in the circuit module of MindQuantum can be used to read the computing result saved in molecule_file
, build the UCCSD wave function by one click, and obtain the corresponding quantum circuit.
ansatz_circuit, \
init_amplitudes, \
ansatz_parameter_names, \
hamiltonian_QubitOp, \
n_qubits, n_electrons = generate_uccsd(molecule_file, th=-1)
ccsd:-7.882352909152705.
fci:-7.882362286798725.
generate_uccsd
packs functions related to the unitary coupled-cluster, including multiple steps such as deriving a molecular Hamiltonian, building a unitary coupled-cluster ansatz operator, and extracting a coupled-cluster coefficient computed by CCSD. This function reads the molecule by entering its file path. The parameter th
indicates the to-be-updated gradient threshold of a parameter in the quantum circuit. In the section Building a Unitary Coupled-Cluster Ansatz Step by Step, we will demonstrate how to use the related interfaces of MindQuantum to complete the steps. A complete quantum circuit includes an HF initial state and a UCCSD ansatz, as shown in the following code:
total_circuit = hartreefock_wfn_circuit + ansatz_circuit
total_circuit.summary()
print("Number of parameters: %d" % (len(ansatz_parameter_names)))
==============================Circuit Summary==============================
|Total number of gates : 12612. |
|Parameter gates : 640. |
|with 44 parameters are : p40, p9, p8, p3, p32, p28, p15, p4, p18, p22... |
|Number qubit of circuit: 12 |
===========================================================================
Number of parameters: 44
For the LiH molecule, the UCCSD wave function ansatz includes 44 variational parameters. The total number of quantum bit gates of the circuit is 12612, and a total of 12 quantum bits are needed for simulation.
VQE Procedure
The procedure for solving the molecular ground state by using the VQE is as follows:
Prepare the HF initial state:
.Define the wave function ansatz, such as UCCSD.
Convert the wave function into a variational quantum circuit.
Initialize the variational parameters, for example, set all parameters to 0.
Obtain the energy
of the molecular Hamiltonian under the set of variational parameters and the derivative of the energy about the parameters by means of multiple measurements on the quantum computer.Use optimization algorithms, such as gradient descent and BFGS, to update variational parameters on classical computers.
Transfer the new variational parameters to the quantum circuit for updating.
Repeat steps 5 to 7 until the convergence criteria are met.
End.
In step 5, the derivative
from mindquantum.simulator import Simulator
sim = Simulator('projectq', total_circuit.n_qubits)
molecule_pqc = sim.get_expectation_with_grad(Hamiltonian(hamiltonian_QubitOp),
total_circuit)
You can obtain the energy molecule_pqc
.
Next, steps 5 to 7 in VQE optimization need to be performed, that is, parameterized quantum circuits need to be optimized. Based on the MindSpore framework, you can use the parameterized quantum circuit operator molecule_pqc
to build a neural network model, and then optimize the variational parameters by using a method similar to training the neural network.
from mindquantum.framework import MQAnsatzOnlyLayer
class PQCNet(ms.nn.Cell):
def __init__(self, pqc):
super(PQCNet, self).__init__()
self.pqc = pqc
self.net = MQAnsatzOnlyLayer(self.pqc, weight='Zeros')
def construct(self):
return self.net()
molecule_pqcnet = PQCNet(molecule_pqc)
Here, we manually build a basic PQCNet
as a model example. This model can be used similar to a conventional machine learning model, for example, optimizing weights and calculating derivatives. A better choice is to use MQAnsatzOnlyLayer
encapsulated in MindQuantum, which will be demonstrated later.
The built PQCNet
uses the "Zeros"
keyword to initialize all variational parameters to 0. The computing result of CCSD or second order Møller-Plesset perturbation theory (MP2) can also be used as the initial value of the variational parameters of unitary coupled-clusters. In this case,
initial_energy = molecule_pqcnet()
print("Initial energy: %20.16f" % (initial_energy.asnumpy()))
Initial energy: -7.8633575439453125
Finally, the Adam optimizer of MindSpore is used for optimization. The learning rate is set to
optimizer = ms.nn.Adagrad(molecule_pqcnet.trainable_params(), learning_rate=4e-2)
train_pqcnet = ms.nn.TrainOneStepCell(molecule_pqcnet, optimizer)
eps = 1.e-8
energy_diff = eps * 1000
energy_last = initial_energy.asnumpy() + energy_diff
iter_idx = 0
while abs(energy_diff) > eps:
energy_i = train_pqcnet().asnumpy()
if iter_idx % 5 == 0:
print("Step %3d energy %20.16f" % (iter_idx, float(energy_i)))
energy_diff = energy_last - energy_i
energy_last = energy_i
iter_idx += 1
print("Optimization completed at step %3d" % (iter_idx - 1))
print("Optimized energy: %20.16f" % (energy_i))
print("Optimized amplitudes: \n", molecule_pqcnet.net.weight.asnumpy())
Step 0 energy -7.8633575439453125
Step 5 energy -7.8726239204406738
Step 10 energy -7.8821778297424316
Step 15 energy -7.8822836875915527
Step 20 energy -7.8823199272155762
Step 25 energy -7.8823370933532715
Step 30 energy -7.8823437690734863
Step 35 energy -7.8618836402893066
Step 40 energy -7.8671770095825195
Step 45 energy -7.8751692771911621
Step 50 energy -7.8822755813598633
Step 55 energy -7.8812966346740723
Step 60 energy -7.8823189735412598
Step 65 energy -7.8823523521423340
Optimization completed at step 67
Optimized energy: -7.8823528289794922
Optimized amplitudes:
[ 2.3980068e-04 1.8912849e-03 3.5044324e-02 1.6005965e-02
-1.9985158e-07 9.0940151e-04 1.6222824e-05 1.4160988e-02
-1.1072063e-07 9.0867787e-04 1.3825165e-05 1.4166672e-02
-5.4699212e-04 4.2679289e-04 2.8641545e-03 5.3817011e-02
2.3320253e-04 1.7034533e-07 6.6684343e-08 -2.7686235e-07
7.2332718e-08 1.2834757e-05 -1.0439425e-04 7.1826143e-08
3.6483241e-06 6.1677817e-08 3.1003920e-06 7.9770159e-04
-5.4951470e-02 3.0904056e-03 -4.4321241e-05 8.5840838e-07
-1.9589644e-08 -4.9430941e-08 8.6163556e-07 -2.5008637e-07
2.1493735e-08 -4.6331229e-06 3.0904033e-03 9.5311613e-08
-4.8755901e-08 2.0483398e-08 -3.9453280e-06 3.7235476e-04]
It can be seen that the computing result of unitary coupled-cluster is very close to that of FCI, and has good accuracy.
Building a Unitary Coupled-Cluster Ansatz Step by Step
In the preceding part, the generate_uccsd
is used to build all the content required for designing the unitary coupled-cluster. In this section, the steps are split, we get the coupled-cluster operator, the corresponding quantum circuit and the initial guess of the variational parameters from the classical CCSD results.
First, import some extra dependencies, including the related functions of the HiQfermion module in MindQuantum.
from mindquantum.algorithm.nisq.chem import Transform
from mindquantum.algorithm.nisq.chem import get_qubit_hamiltonian
from mindquantum.algorithm.nisq.chem import uccsd_singlet_generator, uccsd_singlet_get_packed_amplitudes
from mindquantum.core import TimeEvolution
from mindquantum.framework import MQAnsatzOnlyLayer
The molecule Hamiltonian uses get_qubit_hamiltonian
to read the previous computing result. The result is as follows:
hamiltonian_QubitOp = get_qubit_hamiltonian(molecule_of)
The unitary coupled-cluster operator uccsd_singlet_generator
. Provide the total number of quantum bits (total number of spin orbits) and the total number of electrons, and set anti_hermitian=True
.
ucc_fermion_ops = uccsd_singlet_generator(
molecule_of.n_qubits, molecule_of.n_electrons, anti_hermitian=True)
The ucc_fermion_ops
built in the previous step is parameterized. Use the Jordan-Wigner transformation to map the Fermi excitation operator to the Pauli operator:
ucc_qubit_ops = Transform(ucc_fermion_ops).jordan_wigner()
Next, we need to obtain the quantum circuit corresponding to the unitary operator TimeEvolution
can generate the circuit corresponding to TimeEvolution
is used, ucc_qubit_ops
already contains the complex number factor ucc_qubit_ops
by ucc_qubit_ops
.
ansatz_circuit = TimeEvolution(ucc_qubit_ops.imag, 1.0).circuit
ansatz_parameter_names = ansatz_circuit.params_name
ansatz_parameter_names
is used to record the parameter names in the circuit. So far, we have obtained the contents required by the VQE quantum circuit, including the Hamiltonian hamiltonian_QubitOp
and the parameterized wave function ansatz ansatz_circuit
. By referring to the preceding steps, we can obtain a complete state preparation circuit. hartreefock_wfn_circuit
mentioned above is used as the Hartree-Fock reference state:
total_circuit = hartreefock_wfn_circuit + ansatz_circuit
total_circuit.summary()
======================================Circuit Summary======================================
|Total number of gates : 12612. |
|Parameter gates : 640. |
|with 44 parameters are : d1_3, d2_26, d2_6, d2_1, d2_2, d2_14, d1_1, s_1, d2_16, d2_11...|
|Number qubit of circuit: 12 |
===========================================================================================
Next, you need to provide a reasonable initial value for the variational parameter. The PQCNet
built in the preceding text uses 0 as the initial guess by default, which is feasible in most cases. However, using CCSD’s computational data as a starting point for UCC may be better. Use the uccsd_singlet_get_packed_amplitudes
function to extract CCSD parameters from molecule_of
.
init_amplitudes_ccsd = uccsd_singlet_get_packed_amplitudes(
molecule_of.ccsd_single_amps, molecule_of.ccsd_double_amps, molecule_of.n_qubits, molecule_of.n_electrons)
init_amplitudes_ccsd = [init_amplitudes_ccsd[param_i] for param_i in ansatz_parameter_names]
MQAnsatzOnlyLayer
can be used to easily obtain a machine learning model based on a variational quantum circuit by using a parameter and a quantum circuit:
sim = Simulator('projectq', total_circuit.n_qubits)
grad_ops = sim.get_expectation_with_grad(Hamiltonian(hamiltonian_QubitOp.real),
total_circuit)
molecule_pqcnet = MQAnsatzOnlyLayer(grad_ops)
init_amplitudes_ccsd
(coupled-cluster coefficient computed by CCSD) is used as an initial variational parameter:
molecule_pqcnet.weight = Parameter(ms.Tensor(init_amplitudes_ccsd, molecule_pqcnet.weight.dtype))
initial_energy = molecule_pqcnet()
print("Initial energy: %20.16f" % (initial_energy.asnumpy()))
Initial energy: -7.8173098564147949
In this example, CCSD’s initial guess does not provide a better starting point. You can test and explore more molecules and more types of initial values (such as initial guesses of random numbers). Finally, the VQE is optimized. The optimizer still uses Adam, and the convergence standard remains unchanged. The code used for optimization is basically the same as that described in the preceding sections. You only need to update the corresponding variables.
optimizer = ms.nn.Adagrad(molecule_pqcnet.trainable_params(), learning_rate=4e-2)
train_pqcnet = ms.nn.TrainOneStepCell(molecule_pqcnet, optimizer)
print("eps: ", eps)
energy_diff = eps * 1000
energy_last = initial_energy.asnumpy() + energy_diff
iter_idx = 0
while abs(energy_diff) > eps:
energy_i = train_pqcnet().asnumpy()
if iter_idx % 5 == 0:
print("Step %3d energy %20.16f" % (iter_idx, float(energy_i)))
energy_diff = energy_last - energy_i
energy_last = energy_i
iter_idx += 1
print("Optimization completed at step %3d" % (iter_idx - 1))
print("Optimized energy: %20.16f" % (energy_i))
print("Optimized amplitudes: \n", molecule_pqcnet.weight.asnumpy())
eps: 1e-08
Step 0 energy -7.8173098564147949
Step 5 energy -7.8740763664245605
Step 10 energy -7.8818783760070801
Step 15 energy -7.8821649551391602
Step 20 energy -7.8822622299194336
Step 25 energy -7.8823084831237793
Step 30 energy -7.8823180198669434
Step 35 energy -7.8737111091613770
Step 40 energy -7.8724455833435059
Step 45 energy -7.8801403045654297
Step 50 energy -7.8821926116943359
Step 55 energy -7.8818311691284180
Step 60 energy -7.8823456764221191
Optimization completed at step 64
Optimized energy: -7.8823523521423340
Optimized amplitudes:
[-2.4216002e-04 1.8924323e-03 -3.4653045e-02 1.5943546e-02
3.6362690e-07 9.0936717e-04 -1.7181528e-05 1.4154296e-02
-4.4650793e-08 9.0864423e-04 -2.6399141e-06 1.4159971e-02
5.4558384e-04 4.2672374e-04 -2.8494308e-03 5.3833455e-02
2.3033506e-04 1.2578158e-06 3.3855862e-08 7.3955505e-08
-5.2005623e-07 2.9746575e-08 1.2325607e-08 1.1919828e-05
-1.0492613e-04 7.9503102e-04 3.8478893e-06 5.9738107e-07
-5.4855812e-02 3.0889052e-03 7.9252044e-05 -1.5384763e-06
-1.5373821e-06 -3.0784176e-07 -3.5303248e-08 1.7360321e-08
4.4359115e-07 -4.9067144e-06 3.0889027e-03 1.3888703e-07
-1.6715177e-08 6.3234533e-09 -7.5149819e-07 3.7140178e-04]
Summary
In this case, the ground state energy of the LiH molecule is obtained by using the quantum neural network in two methods. In the first method, we use the generate_uccsd
function packaged by MindQuantum to generate a quantum neural network that can solve this problem. In the second method, we build a similar quantum neural network step by step. The final results are the same.