<jupyter_start><jupyter_text>Contour deformation In the context of GW method, contour deformation (CD) technique is used in conjunction with resolution of identity (RI) to reduce the formal scaling of the self-energy calculation. Compared to widely used analytic continuation approach it provides a means to evaluate self-energy  directly on the real axis without employing Pade approximants or non-linear least squares fit and potentially offering superior accuracy. Here, we provide a brief outline of the theory behind CD and give an example of the self-energy calculation within CD without invoking RI in order to facilitate comparison with the results prsented above. Detailed discussion of the CD can be found in the following papers:1. <NAME>., <NAME>., <NAME>., & <NAME>. (2018). Core-Level Binding Energies from GW : An Efficient Full-Frequency Approach within a Localized Basis. Journal of Chemical Theory and Computation, 14(9), 4856–4869. https://doi.org/10.1021/acs.jctc.8b004582. <NAME>., <NAME>., <NAME>., <NAME>., <NAME>., <NAME>., & <NAME>. (2011). Electronic properties of interfaces and defects from many-body perturbation theory: Recent developments and applications. Physica Status Solidi (B), 248(2), 275–289. https://doi.org/10.1002/pssb.201046094 CD is used to recast the convolution in the GW expression of self-energy as a difference between two integrals, one which can be performed analytically whereas the other can be evaluated numerically on a relatively small grid. This is achieved by closing the inegration contour as shown below [2]:![Integration contour used to evaluate $\Sigma(\omega)$](CD_scheme.jpg)$$\Sigma(r_1,r_2, \omega) = \frac{i}{2\pi} \int_{-\infty}^{+\infty} e^{i\omega^{\prime} \eta} G(r_1, r_2, \omega + \omega^{\prime}) W(r_1, r_2, \omega^{\prime}) d\omega^{\prime}\\ = \frac{i}{2\pi} \oint_{\Gamma} G(r_1, r_2, \omega + z) W(r_1, r_2, z) dz -  \frac{1}{2\pi} \int_{-\infty}^{+\infty} G(r_1, r_2, \omega + i\omega^{\prime}) W(r_1, r_2, i\omega^{\prime}) d\omega^{\prime}$$Depending on the $\omega$ value the lower-left and the upper-right loops of the contour can enclose one or several poles of the zero-order Green's function whereas the poles of the screened Coulomb interaction never fall within the contour. This allowes to evaluate the countour integral as a sum of corresponding residues with apropriate signs (note that the upper-right loop is traversed counter-clockwise, while the lower-left loop is traversed clockwise). The imaginary axis contribution is calculated using Gauss-Legendre grid. Importantly, the intgrals over the arches vanish iff the screened Coulomb interaction does not contain the exchange contribution.<jupyter_code>import psi4
import numpy as np
import scipy as sp
from matplotlib import pyplot as plt
%matplotlib inline
from IPython.core.display import display, HTML

display(HTML("<style>.container {width:95% !important;}</style>"))


psi4.set_options({'basis' : 'cc-pvdz', 'd_convergence' : 1e-7,'scf_type' : 'out_of_core', 'dft_spherical_points' : 974, 'dft_radial_points' : 150 })

sih4 = psi4.geometry("""
                       Si  0.0000  0.0000  0.0000
                       H   0.8544 -0.8544  0.8544
                       H  -0.8544  0.8544  0.8544
                       H  -0.8544 -0.8544 -0.8544
                       H   0.8544  0.8544 -0.8544
                       symmetry c1
                       units angstrom
""")



psi4.set_output_file('sih4_ccpvdz.out')

scf_e, scf_wfn = psi4.energy('PBE', return_wfn=True)


print("DFT energy is %16.10f" % scf_e)
epsilon = np.asarray(scf_wfn.epsilon_a())
print(epsilon*psi4.constants.hartree2ev)<jupyter_output>DFT energy is  -291.7041409214
[-1.78015161e+03 -1.37583520e+02 -9.43390440e+01 -9.43390436e+01
 -9.43390431e+01 -1.34588824e+01 -8.45064857e+00 -8.45064828e+00
 -8.45064799e+00  6.10841118e-01  6.10841192e-01  6.10841524e-01
  1.48415008e+00  6.57736831e+00  6.57736917e+00  6.57736959e+00
  7.81254643e+00  7.81254793e+00  9.19553485e+00  1.12297170e+01
  1.12297176e+01  1.12297180e+01  1.49820206e+01  2.23924790e+01
  2.23924791e+01  2.23924799e+01  3.18168512e+01  3.18168512e+01
  3.18168512e+01  3.94598408e+01  3.94598409e+01  3.94598409e+01
  4.47366037e+01  4.47366043e+01  4.78594674e+01  5.21424175e+01
  5.21424181e+01  5.21424191e+01]<jupyter_text>``` SCF Total Energy (Ha):    -291.7041423085   (MOLGW) ```<jupyter_code>import GW
gw_par = {'no_qp' : 9, 'nv_qp' : 0, 'nomega_sigma' : 501, 'step_sigma' : 0.01, 'gl_npoint' : 200, 'debug' : False, 'low_mem' : True }
gw_sih4_dz_cd1 = GW.GW_DFT(scf_wfn, sih4, gw_par)
gw_sih4_dz_cd1.print_summary()<jupyter_output>Number of basis functions:  38
occ/virt: 9/29
Attempting to create RI basis set for CC-PVDZ (RIFIT)... 
Auxiliary basis set has been generated!
Number of auxiliary basis functions:  132
Fraction of HF exchange is  0.000
Running in production mode!
Shape of the omega_grid_all is  (9, 501)
Caculating GW self-energy via contour deformation
Calculation of the integral term requires    0.057 Gb
Calculation of the residue term requires     0.144 Gb
Using low-memory algorithm
Finished calculating self-energy
Performing one-shot G0W0
SigX - Vxc
[-3.35157679 -1.05362364 -0.74918378 -0.74918379 -0.74918381 -0.23312109
 -0.16336648 -0.16336648 -0.16336649]
Perfoming graphic solution of the inverse Dyson equation
Done!
E^lin, eV  E^graph, eV  Z 
 -1837.458186   -1844.206682       0.745068
  -150.388473    -176.635681       0.541100
  -105.332418    -103.499110       0.811066
  -105.357625    -103.499110       0.811066
  -105.385443    -103.499110       0.811066
   -16.727115     -16.367367       0[...]<jupyter_text>``` GW eigenvalues (eV) RI           E0         SigX-Vxc       SigC          Z        E_qp^lin     E_qp^graph   1  -1780.564829   -91.201582    14.293043     0.745186 -1837.875974 -1844.625575   2   -137.607724   -28.671569     5.027952     0.553859  -150.702954  -151.679909   3    -94.335614   -20.386438     6.804644     0.811169  -105.352751  -103.499211   4    -94.335614   -20.386438     6.804644     0.811169  -105.352751  -103.499211   5    -94.335614   -20.386438     6.804644     0.811169  -105.352751  -103.499211   6    -13.460219    -6.343842     2.093938     0.768966   -16.728252   -16.369167   7     -8.451237    -4.446523     0.319568     0.881765   -12.090241   -12.020075   8     -8.451237    -4.446523     0.319568     0.881765   -12.090241   -12.020075   9     -8.451237    -4.446523     0.319568     0.881765   -12.090241   -12.020075```<jupyter_code>gw_par = {'no_qp' : 9, 'nv_qp' : 0, 'nomega_sigma' : 501, 'step_sigma' : 0.01, 'analytic_W': True, 'gl_npoint' : 200, 'debug' : False, 'low_mem' : False }
gw_sih4_dz_cd2 = GW.GW_DFT(scf_wfn, sih4, gw_par)
gw_sih4_dz_cd2.print_summary()<jupyter_output>Number of basis functions:  38
occ/virt: 9/29
Attempting to create RI basis set for CC-PVDZ (RIFIT)... 
Auxiliary basis set has been generated!
Number of auxiliary basis functions:  132
Fraction of HF exchange is  0.000
Running in production mode!
Shape of the omega_grid_all is  (9, 501)
Caculating GW self-energy via contour deformation
Analytic W has been requested; performing RPA calculation
Shape of omega tensor is  (38, 38, 261)
Calculation of the integral term requires    0.057 Gb
Calculation of the residue term requires     5.467 Gb
Finished calculating self-energy
Performing one-shot G0W0
SigX - Vxc
[-3.35157679 -1.05362364 -0.74918378 -0.74918379 -0.74918381 -0.23312109
 -0.16336648 -0.16336648 -0.16336649]
Perfoming graphic solution of the inverse Dyson equation
Done!
E^lin, eV  E^graph, eV  Z 
 -1837.461131   -1844.206671       0.745068
  -150.389764    -176.635691       0.541093
  -105.334685    -103.499112       0.811066
  -105.359892    -103.499112       0.811066
  -105.38[...]
{'fraction_non_alphanumeric': 0.10969813744380219, 'fraction_numerical': 0.2, 'mean_word_length': 3.195043103448276, 'pattern_counts': {'":': 0, '<': 23, '<?xml version=': 0, '>': 23, 'https://': 2, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '3', 'dir': 'jupyter-structured-clean-dedup', 'id': '603551', 'license': 'MIT', 'n_tokens_mistral': 3807, 'n_tokens_neox': 2881, 'n_words': 827, 'path': 'tests/CDGW_DFT4cp2k-SiH4.ipynb', 'repo_name': 'aakunitsa/GW-approximation', 'size': '7842'}