Note

Run interactively in the IBM Quantum lab.

# Electronic structure¶

## Introduction¶

The molecular Hamiltonian is

$\mathcal{H} = - \sum_I \frac{\nabla_{R_I}^2}{M_I} - \sum_i \frac{\nabla_{r_i}^2}{m_e} - \sum_I\sum_i \frac{Z_I e^2}{|R_I-r_i|} + \sum_i \sum_{j>i} \frac{e^2}{|r_i-r_j|} + \sum_I\sum_{J>I} \frac{Z_I Z_J e^2}{|R_I-R_J|}$

Because the nuclei are much heavier than the electrons they do not move on the same time scale and therefore, the behavior of nuclei and electronc can be decoupled. This is the Born-Oppenheimer approximation.

Therefore, one can first tackle the electronic problem with nuclear coordinate entering only as parameters. The energy levels of the electrons in the molecule can be found by solving the non-relativistic time independant Schroedinger equation,

$\mathcal{H}_{\text{el}} |\Psi_{n}\rangle = E_{n} |\Psi_{n}\rangle$

where

$\mathcal{H}_{\text{el}} = - \sum_i \frac{\nabla_{r_i}^2}{m_e} - \sum_I\sum_i \frac{Z_I e^2}{|R_I-r_i|} + \sum_i \sum_{j>i} \frac{e^2}{|r_i-r_j|}.$

In particular the ground state energy is given by:

\begin{align}\begin{aligned} E_0 = \frac{\langle \Psi_0 | H_{\text{el}} | \Psi_0 \rangle}{\langle \Psi_0 | \Psi_0 \rangle}\\where :math:\Psi_0 is the ground state of the system.\end{aligned}\end{align}

However, the dimensionality of this problem grows exponentially with the number of degrees of freedom. To tackle this issue we would like to prepare $$\Psi_0$$ on a quantum computer and measure the Hamiltonian expectation value (or $$E_0$$) directly.

So how do we do that concretely?

## The Hartree-Fock initial state¶

A good starting point for solving this problem is the Hartree-Fock (HF) method. This method approximates a N-body problem into N one-body problems where each electron evolves in the mean-field of the others. Classically solving the HF equations is efficient and leads to the exact exchange energy but does not include any electron correlation. Therefore, it is usually a good starting point to start adding correlation.

The Hamtiltonian can then be re-expressed in the basis of the solutions of the HF method, also called Molecular Orbitals (MOs):

\begin{align}\begin{aligned} \hat{H}_{elec}=\sum_{pq} h_{pq} \hat{a}^{\dagger}_p \hat{a}_q + \frac{1}{2} \sum_{pqrs} h_{pqrs} \hat{a}^{\dagger}_p \hat{a}^{\dagger}_q \hat{a}_r \hat{a}_s\\with the 1-body integrals\end{aligned}\end{align}
\begin{align}\begin{aligned} h_{pq} = \int \phi^*_p(r) \left( -\frac{1}{2} \nabla^2 - \sum_{I} \frac{Z_I}{R_I- r} \right) \phi_q(r)\\and 2-body integrals\end{aligned}\end{align}
$h_{pqrs} = \int \frac{\phi^*_p(r_1) \phi^*_q(r_2) \phi_r(r_2) \phi_s(r_1)}{|r_1-r_2|}.$

The MOs ($$\phi_u$$) can be occupied or virtual (unoccupied). One MO can contain 2 electrons. However, in what follows we actually work with Spin Orbitals which are associated with a spin up ($$\alpha$$) of spin down ($$\beta$$) electron. Thus Spin Orbitals can contain one electron or be unoccupied.

We now show how to concreatly realise these steps with Qiskit.

Qiskit is interfaced with different classical codes which are able to find the HF solutions. Interfacing between Qiskit and the following codes is already available: * Gaussian * Psi4 * PyQuante * PySCF

In the following we set up a PySCF driver, for the hydrogen molecule at equilibrium bond length (0.735 angstrom) in the singlet state and with no charge.

[1]:

from qiskit.chemistry.drivers import PySCFDriver, UnitsType, Molecule
molecule = Molecule(geometry=[['H', [0., 0., 0.]],
['H', [0., 0., 0.735]]],
charge=0, multiplicity=1)
driver = PySCFDriver(molecule = molecule, unit=UnitsType.ANGSTROM, basis='sto3g')


For further information about the drivers see https://qiskit.org/documentation/apidoc/qiskit.chemistry.drivers.html

## The mapping from fermions to qubits¶

The Hamiltonian given in the previous section is expressed in terms of fermionic operators. To encode the problem into the state of a quantum computer, these operators must be mapped to spin operators (indeed the qubits follow spin statistics).

There exist different mapping types with different properties. Qiskit already supports the following mappings: * The Jordan-Wigner ‘jordan_wigner’ mapping (über das paulische äquivalenzverbot. In The Collected Works of Eugene Paul Wigner (pp. 109-129). Springer, Berlin, Heidelberg (1993)). * The Parity ‘parity’ (The Journal of chemical physics, 137(22), 224109 (2012)) * The Bravyi-Kitaev ‘bravyi_kitaev’ (Annals of Physics, 298(1), 210-226 (2002)) * The Bravyi-Kitaev Super Fast ‘bksf’ (Annals of Physics, 298(1), 210-226 (2002))

The Jordan-Wigner mapping is particularly interesting as it maps each Spin Orbital to a qubit (as shown on the Figure above).

Here we set up an object which contains all the information about any transformation of the fermionic Hamiltonian to the qubits Hamiltonian. In this example we simply ask for the Jordan-Wigner mapping.

[2]:

from qiskit.chemistry.transformations import (FermionicTransformation,
FermionicTransformationType,
FermionicQubitMappingType)

fermionic_transformation = FermionicTransformation(
transformation=FermionicTransformationType.FULL,
qubit_mapping=FermionicQubitMappingType.JORDAN_WIGNER,
two_qubit_reduction=False,
freeze_core=False)


If we now transform this Hamiltonian for the given driver defined above we get our qubit operator:

[3]:

qubit_op, _ = fermionic_transformation.transform(driver)
print(qubit_op)
print(fermionic_transformation.molecule_info)

SummedOp([
-0.8105479805373264 * IIII,
0.17218393261915527 * IIIZ,
-0.22575349222402502 * IIZI,
0.17218393261915532 * IZII,
-0.22575349222402497 * ZIII,
0.12091263261776641 * IIZZ,
0.16892753870087926 * IZIZ,
0.04523279994605789 * XXYY,
0.04523279994605789 * YYYY,
0.04523279994605789 * XXXX,
0.04523279994605789 * YYXX,
0.1661454325638243 * ZIIZ,
0.1661454325638243 * IZZI,
0.17464343068300467 * ZIZI,
0.12091263261776641 * ZZII
])
{'num_particles': [1, 1], 'num_orbitals': 4, 'two_qubit_reduction': False, 'z2_symmetries': <qiskit.aqua.operators.legacy.weighted_pauli_operator.Z2Symmetries object at 0x7f6894047eb0>}


In the minimal (STO-3G) basis set 4 qubits are required. We could even lower the qubit count by using the Parity mapping which allows to get rid of to qubits by symmetry considerations.

[4]:

fermionic_transformation_2 = FermionicTransformation(
transformation=FermionicTransformationType.FULL,
qubit_mapping=FermionicQubitMappingType.PARITY,
two_qubit_reduction=True,
freeze_core=False)
qubit_op_2, _ = fermionic_transformation_2.transform(driver)
print(qubit_op_2)

SummedOp([
-1.0523732457728587 * II,
0.3979374248431802 * IZ,
-0.3979374248431802 * ZI,
-0.011280104256235324 * ZZ,
0.18093119978423147 * XX
])


This time only 2 qubits are needed.

Another possibility is to use the Particle-Hole tranformation (Physical Review A, 98(2), 022322 (2018)). This shifts the vacuum state to a state lying in the N-particle Fock space. In this representation the HF (reference) state has a null energy and the optimization procedure is more faster.

[5]:

fermionic_transformation_3 = FermionicTransformation(
transformation=FermionicTransformationType.PARTICLE_HOLE,
qubit_mapping=FermionicQubitMappingType.JORDAN_WIGNER,
two_qubit_reduction=False,
freeze_core=False)
qubit_op_3, _ = fermionic_transformation_3.transform(driver)
print(qubit_op_3)

SummedOp([
1.026420010665659 * IIII,
0.17218393261915532 * IIIZ,
-0.22575349222402522 * IIZI,
0.17218393261915532 * IZII,
-0.22575349222402522 * ZIII,
0.16892753870087926 * IZIZ,
0.04523279994605789 * YYYY,
0.04523279994605789 * XXYY,
0.04523279994605789 * YYXX,
0.04523279994605789 * XXXX,
0.12091263261776641 * IIZZ,
0.1661454325638243 * IZZI,
0.1661454325638243 * ZIIZ,
0.17464343068300467 * ZIZI,
0.12091263261776641 * ZZII
])


The list of available mappings and transformations are

[6]:

print('*Transformations')
for fer_transform in FermionicTransformationType:
print(fer_transform)

print('\n*Mappings')
for fer_mapping in FermionicQubitMappingType:
print(fer_mapping)

*Transformations
FermionicTransformationType.FULL
FermionicTransformationType.PARTICLE_HOLE

*Mappings
FermionicQubitMappingType.JORDAN_WIGNER
FermionicQubitMappingType.PARITY
FermionicQubitMappingType.BRAVYI_KITAEV


Now that the Hamiltonian is ready, it can be used in a quantum algorithm to find information about the electronic structure of the corresponding molecule. Check out our tutorials on Ground State Calculation and Excited States Calculation to learn more about how to do that in Qiskit!

[7]:

import qiskit.tools.jupyter
%qiskit_version_table


### Version Information

Qiskit SoftwareVersion
Qiskit0.23.1
Terra0.16.1
Aer0.7.1
Ignis0.5.1
Aqua0.8.1
IBM Q Provider0.11.1
System information
Python3.8.6 (default, Oct 28 2020, 13:08:18) [GCC 7.5.0]
OSLinux
CPUs2
Memory (Gb)6.791393280029297
Mon Nov 30 18:47:36 2020 UTC

### This code is a part of Qiskit

[ ]: