[doc] add gw_embedding to DOC

doc/documentation.rst

@@ -57,6 +57,7 @@ Module reference manual
     dft_managers
     dmft_cycle
     dmft_tools
+    gw_embedding
     io_tools
     postprocessing
     util

doc/input_output/DMFT_input/generate_doc_from_comments.sh

@@ -38,6 +38,7 @@ The aim of this section is to provide a comprehensive listing of all the input f
     general
     solver
     dft
+    gw
     advanced
 
 Below an exhaustive list containing all the parameters marked by section.
@@ -65,6 +66,8 @@ awk '/\[  solver  \]/{flag=1; c=0} flag; /\[ /&& ++c==2{flag=0}' matches_comment
 grep '::' solver.tmp | grep -o ':: \(.*\)' | sed 's/:: //g' | tr ' \n' '; ' > solver_list.tmp
 awk '/\[  dft  \]/{flag=1; c=0} flag; /\[ /&& ++c==2{flag=0}' matches_comments.txt | head -n -2 | tail -n +3 > dft.tmp
 grep '::' dft.tmp | grep -o ':: \(.*\)' | sed 's/:: //g' | tr ' \n' '; ' > dft_list.tmp
+awk '/\[  gw  \]/{flag=1; c=0} flag; /\[ /&& ++c==2{flag=0}' matches_comments.txt | head -n -2 | tail -n +3 > gw.tmp
+grep '::' gw.tmp | grep -o ':: \(.*\)' | sed 's/:: //g' | tr ' \n' '; ' > gw_list.tmp
 awk '/\[  advanced  \]/{flag=1; c=0} flag; /\[ /&& ++c==2{flag=0}' matches_comments.txt | head -n -2 | tail -n +3 > advanced.tmp
 grep '::' advanced.tmp | grep -o ':: \(.*\)' | sed 's/:: //g' | tr ' \n' '; ' > advanced_list.tmp
 
@@ -81,6 +84,11 @@ printf '**[dft]**' >> input.rst
 echo -e "\n"  >> input.rst
 cat dft_list.tmp >> input.rst
 
+echo -e "\n"  >> input.rst
+printf '**[gw]**' >> input.rst
+echo -e "\n"  >> input.rst
+cat gw_list.tmp >> input.rst
+
 echo -e "\n"  >> input.rst
 printf '**[advanced]**' >> input.rst
 echo -e "\n"  >> input.rst
@@ -139,6 +147,22 @@ echo -e "\n"  >> dft.rst
 cat dft.tmp >> dft.rst
 ##############
 
+###############
+cat > gw.rst << EOF
+
+[GW]: GW embedding inputs
+-------------------------
+
+List of parameters for the GW embedding calculation. The parameters are ignored unless gw_embedding=true, or interactions are read from cRPA.
+
+
+EOF
+
+#cat gw_list.tmp >> gw.rst
+echo -e "\n"  >> gw.rst
+cat gw.tmp >> gw.rst
+##############
+#
 ##############
 cat > advanced.rst << EOF
 [advanced]: Advanced inputs

python/solid_dmft/gw_embedding/bdft_converter.py

@@ -47,6 +47,26 @@
 HARTREE_EV = physical_constants['Hartree energy in eV'][0]
 
 def _get_dlr_from_IR(Gf_ir, ir_kernel, mesh_dlr_iw, dim=2):
+    r"""
+    Interpolate a given Gf from IR mesh to DLR mesh
+
+    Parameters
+    ----------
+    Gf_ir : np.ndarray
+        Green's function on IR mesh
+    ir_kernel : sparse_ir
+        IR kernel object
+    mesh_dlr_iw : MeshDLRImFreq
+        DLR mesh
+    dim : int, optional
+        dimension of the Green's function, defaults to 2
+
+    Returns
+    -------
+    Gf_dlr : BlockGf or Gf
+        Green's function on DLR mesh
+    """
+
     n_orb = Gf_ir.shape[-1]
     stats = 'f' if mesh_dlr_iw.statistic == 'Fermion' else 'b'
 
@@ -71,6 +91,30 @@ def _get_dlr_from_IR(Gf_ir, ir_kernel, mesh_dlr_iw, dim=2):
 
 
 def calc_Sigma_DC_gw(Wloc_dlr, Gloc_dlr, Vloc, verbose=False):
+    r"""
+    Calculate the double counting part of the self-energy from the screened Coulomb interaction
+
+    Parameters
+    ----------
+    Wloc_dlr : BlockGf or Gf with MeshDLR
+        screened Coulomb interaction
+    Gloc_dlr : BlockGf or Gf with MeshDLR
+        local Green's function
+    Vloc : np.ndarray
+        local Coulomb interaction
+    verbose : bool, optional
+        print additional information, defaults to False
+
+    Returns
+    -------
+    Sig_DC_dlr : BlockGf or Gf
+        double counting part of the self-energy
+    Sig_DC_hartree : np.ndarray
+        static Hartree part of the self-energy
+    Sig_DC_exchange : np.ndarray
+        static exchange part of the self-energy
+    """
+
     if isinstance(Gloc_dlr, BlockGf):
         Sig_DC_dlr_list = []
         Sig_DC_hartree_list = {}
@@ -129,7 +173,7 @@ def calc_Sigma_DC_gw(Wloc_dlr, Gloc_dlr, Vloc, verbose=False):
     return Sig_DC_dlr, Sig_DC_hartree, Sig_DC_exchange
 
 
-def calc_W_from_Gloc(Gloc_dlr: Gf | BlockGf, U: np.ndarray | dict) -> Gf | BlockGf:
+def calc_W_from_Gloc(Gloc_dlr, U):
     r"""
 
     Calculate Wijkl from given constant U tensor and Gf on DLRMesh

python/solid_dmft/gw_embedding/iaft.py

@@ -2,7 +2,7 @@
 import sparse_ir
 
 """
-Fourier transform on the imaginary axis based on IR basis and the sparse sampling technique.  
+Fourier transform on the imaginary axis based on IR basis and the sparse sampling technique.
 """
 
 
@@ -97,7 +97,7 @@ def __str__(self):
                "*******************************".format(self.prec, self.beta, self.lmbda, self.nt_f, self.nw_f,
                                                         self.nt_b, self.nw_b)
 
-    def wn_mesh(self, stats: str, ir_notation: bool = True):
+    def wn_mesh(self, stats, ir_notation= True):
         """
         Return Matsubara frequency indices.
         :param stats: str
@@ -117,7 +117,7 @@ def wn_mesh(self, stats: str, ir_notation: bool = True):
             wn_mesh = (wn_mesh-1)//2 if stats == 'f' else wn_mesh//2
         return wn_mesh
 
-    def tau_to_w(self, Ot, stats: str):
+    def tau_to_w(self, Ot, stats):
         """
         Fourier transform from imaginary-time axis to Matsubara-frequency axis
         :param Ot: numpy.ndarray
@@ -172,7 +172,7 @@ def w_to_tau(self, Ow, stats):
         Ot = Ot.reshape((Ttw.shape[0],) + Ow_shape[1:])
         return Ot
 
-    def w_interpolate(self, Ow, wn_mesh_interp, stats: str, ir_notation: bool = True):
+    def w_interpolate(self, Ow, wn_mesh_interp, stats, ir_notation=True):
         """
         Interpolate a dynamic object to arbitrary points on the Matsubara axis.
 
@@ -213,7 +213,7 @@ def w_interpolate(self, Ow, wn_mesh_interp, stats: str, ir_notation: bool = True
         Ow_interp = Ow_interp.reshape((wn_indices.shape[0],) + Ow_shape[1:])
         return Ow_interp
 
-    def tau_interpolate(self, Ot, tau_mesh_interp, stats: str):
+    def tau_interpolate(self, Ot, tau_mesh_interp, stats):
         """
          Interpolate a dynamic object to arbitrary points on the imaginary-time axis.
 

python/solid_dmft/gw_embedding/qp_evs_to_eig.py

@@ -1,47 +1,51 @@
 import sys
 
-import numpy as np
 from h5 import HDFArchive
 from scipy.constants import physical_constants
 
-###
-# This script read bdft output and dump g0w0 eigenvalues to si.eig for wannier90 interpolation
-###
-
-# first arg is the h5 to use
-if len(sys.argv) < 2:
-    print ("Usage: python qp_evs_to_eig.py <h5>")
-    quit()
-print('h5 archive:', str(sys.argv[1]))
-bdft_output = str(sys.argv[1])
-
-Hartree_eV = physical_constants['Hartree energy in eV'][0]
-
-###### params ######
-# bdft output is defined by "ouptut" in "evgw0" block.
-# number of bands used in wannier90.
-# It should be consistent with  "num_bands" in si.win
-nbnd_pp = 3
-# number of bands to include in the beginning
-excl = 20
-# iteration of evGW0
-it = None
-# seed for w90
-seed = 'svo'
-
-############
-############
-
-# Read evGW0 eigenvalues
-with HDFArchive(bdft_output, 'r') as ar:
-    if not it:
-        it = ar['scf']['final_iter']
-    eig = ar[f'scf/iter{it}/E_ska'].real * Hartree_eV
-
-
-# Write eigenvalues in the format of Quantum Espresso .eig file
-ns, nkpts, nbnd_gw = eig.shape
-with open(f'{seed}.eig', 'w') as f:
-    for k in range(1, nkpts+1):
-        for i in range(excl+1, excl+nbnd_pp+1):
-            f.write("   {}   {}   {}\n".format(i-excl, k, eig[0, k-1, i-1]))
+
+def extract_qp_eig():
+    r"""
+    This script read bdft output and dump g0w0 eigenvalues to si.eig for wannier90 interpolation
+    """
+
+    # first arg is the h5 to use
+    if len(sys.argv) < 2:
+        print('Usage: python qp_evs_to_eig.py <h5>')
+        quit()
+    print('h5 archive:', str(sys.argv[1]))
+    bdft_output = str(sys.argv[1])
+
+    Hartree_eV = physical_constants['Hartree energy in eV'][0]
+
+    ###### params ######
+    # bdft output is defined by "ouptut" in "evgw0" block.
+    # number of bands used in wannier90.
+    # It should be consistent with  "num_bands" in si.win
+    nbnd_pp = 3
+    # number of bands to include in the beginning
+    excl = 20
+    # iteration of evGW0
+    it = None
+    # seed for w90
+    seed = 'svo'
+
+    ############
+    ############
+
+    # Read evGW0 eigenvalues
+    with HDFArchive(bdft_output, 'r') as ar:
+        if not it:
+            it = ar['scf']['final_iter']
+        eig = ar[f'scf/iter{it}/E_ska'].real * Hartree_eV
+
+    # Write eigenvalues in the format of Quantum Espresso .eig file
+    ns, nkpts, nbnd_gw = eig.shape
+    with open(f'{seed}.eig', 'w') as f:
+        for k in range(1, nkpts + 1):
+            for i in range(excl + 1, excl + nbnd_pp + 1):
+                f.write('   {}   {}   {}\n'.format(i - excl, k, eig[0, k - 1, i - 1]))
+
+
+if __name__ == '__main__':
+    extract_qp_eig()

requirements.txt

@@ -1,7 +1,6 @@
 # Required python packages for this application
 numpy
 scipy
-argparse
 scikit-image
 numpydoc
 # From python 3.11 and newer, toml is part of standard implementation