H2 energy plot computed using SWAPRZ variational form

This notebook demonstrates using Qiskit Aqua Chemistry to plot graphs of the ground state energy of the Hydrogen (H2) molecule over a range of inter-atomic distances using VQE and SWAPRZ. It is compared to the same energies as computed by the ExactEigensolver

This notebook populates a dictionary, that is a progammatic representation of an input file, in order to drive the Qiskit Aqua Chemistry stack. Such a dictionary can be manipulated programmatically and this is indeed the case here where we alter the molecule supplied to the driver in each loop.

This notebook has been written to use the PYQUANTE chemistry driver. See the PYQUANTE chemistry driver readme if you need to install the external PyQuante2 library that this driver requires.


In [1]:
import numpy as np
import pylab
from qiskit_aqua_chemistry import AquaChemistry

# Input dictionary to configure qischem for the chemistry problem.
aqua_chemistry_dict = {
    'problem': {'random_seed': 50},
    'driver': {'name': 'PYQUANTE'},
    'PYQUANTE': {'atoms': '', 'basis': 'sto3g'},
    'operator': {'name': 'hamiltonian', 'qubit_mapping': 'jordan_wigner',
                 'two_qubit_reduction': False},
    'algorithm': {'name': ''},
    'optimizer': {'name': 'COBYLA', 'maxiter': 10000 },
    'variational_form': {'name': 'SWAPRZ'},
    'initial_state': {'name': 'HartreeFock'}
}
molecule = 'H .0 .0 -{0}; H .0 .0 {0}'
algorithms = ['VQE', 'ExactEigensolver']

start = 0.5  # Start distance
by    = 0.5  # How much to increase distance by
steps = 20   # Number of steps to increase by
energies = np.empty([len(algorithms), steps+1])
hf_energies = np.empty(steps+1)
distances = np.empty(steps+1)
eval_counts = np.empty(steps+1)

print('Processing step __', end='')
for i in range(steps+1):
    print('\b\b{:2d}'.format(i), end='', flush=True)
    d = start + i*by/steps
    aqua_chemistry_dict['PYQUANTE']['atoms'] = molecule.format(d/2) 
    for j in range(len(algorithms)):
        aqua_chemistry_dict['algorithm']['name'] = algorithms[j] 
        solver = AquaChemistry()
        result = solver.run(aqua_chemistry_dict)
        energies[j][i] = result['energy']
        hf_energies[i] = result['hf_energy']
        if algorithms[j] == 'VQE':
            eval_counts[i] = result['algorithm_retvals']['eval_count']
    distances[i] = d
print(' --- complete')

print('Distances: ', distances)
print('Energies:', energies)
print('Hartree-Fock energies:', hf_energies)
print('VQE num evaluations:', eval_counts)


Processing step 20 --- complete
Distances:  [0.5   0.525 0.55  0.575 0.6   0.625 0.65  0.675 0.7   0.725 0.75  0.775
 0.8   0.825 0.85  0.875 0.9   0.925 0.95  0.975 1.   ]
Energies: [[-1.05515973 -1.07591361 -1.09262987 -1.10591801 -1.11628597 -1.12416089
  -1.12990475 -1.1338262  -1.13618943 -1.13722136 -1.13711706 -1.13604435
  -1.13414767 -1.1315512  -1.12836188 -1.12467173 -1.12056028 -1.11609624
  -1.11133942 -1.10634211 -1.10115033]
 [-1.05515974 -1.07591361 -1.09262987 -1.10591802 -1.11628599 -1.12416089
  -1.12990476 -1.1338262  -1.13618944 -1.13722136 -1.13711707 -1.13604436
  -1.13414767 -1.13155121 -1.12836188 -1.12467175 -1.12056028 -1.11609624
  -1.11133943 -1.10634212 -1.10115034]]
Hartree-Fock energies: [-1.04299622 -1.0630621  -1.0790507  -1.09157046 -1.10112822 -1.10814997
 -1.11299652 -1.11597525 -1.11734902 -1.11734325 -1.11615145 -1.11393966
 -1.1108504  -1.10700581 -1.10251056 -1.09745432 -1.09191405 -1.08595588
 -1.07963694 -1.07300677 -1.06610866]
VQE num evaluations: [ 685.  687.  707.  717.  666.  755.  828.  668.  750.  786.  645.  875.
  649.  788.  832. 2379.  938.  875.  816.  917.  757.]

In [2]:
pylab.plot(distances, hf_energies, label='Hartree-Fock')
for j in range(len(algorithms)):
    pylab.plot(distances, energies[j], label=algorithms[j])
pylab.xlabel('Interatomic distance')
pylab.ylabel('Energy')
pylab.title('H2 Ground State Energy')
pylab.legend(loc='upper right')


Out[2]:
<matplotlib.legend.Legend at 0x10e44e4a8>

In [3]:
pylab.plot(distances, np.subtract(hf_energies, energies[1]), label='Hartree-Fock')
pylab.plot(distances, np.subtract(energies[0], energies[1]), label='VQE')
pylab.xlabel('Interatomic distance')
pylab.ylabel('Energy')
pylab.yscale('log')
pylab.title('Energy difference from ExactEigensolver')
pylab.legend(loc='center right')


Out[3]:
<matplotlib.legend.Legend at 0x10df285c0>

In [4]:
pylab.plot(distances, eval_counts, '-o', color=[0.8500, 0.3250, 0.0980], label='VQE')
pylab.xlabel('Interatomic distance')
pylab.ylabel('Evaluations')
pylab.title('VQE number of evaluations')
pylab.legend(loc='upper left')


Out[4]:
<matplotlib.legend.Legend at 0x10e4dfb70>

In [ ]: