From a07d27d292757e5d5bb06971ce556af5fa8c9f20 Mon Sep 17 00:00:00 2001 From: Alexander Hampel Date: Thu, 23 May 2024 10:01:03 -0400 Subject: [PATCH] [doc] add gw_embedding to DOC --- doc/documentation.rst | 1 + .../DMFT_input/generate_doc_from_comments.sh | 24 +++++ .../solid_dmft/gw_embedding/bdft_converter.py | 46 +++++++++- python/solid_dmft/gw_embedding/iaft.py | 10 +-- .../solid_dmft/gw_embedding/qp_evs_to_eig.py | 88 ++++++++++--------- requirements.txt | 1 - 6 files changed, 121 insertions(+), 49 deletions(-) diff --git a/doc/documentation.rst b/doc/documentation.rst index 4502a333..02d02b57 100644 --- a/doc/documentation.rst +++ b/doc/documentation.rst @@ -57,6 +57,7 @@ Module reference manual dft_managers dmft_cycle dmft_tools + gw_embedding io_tools postprocessing util diff --git a/doc/input_output/DMFT_input/generate_doc_from_comments.sh b/doc/input_output/DMFT_input/generate_doc_from_comments.sh index 3ee3f120..a6d69b17 100755 --- a/doc/input_output/DMFT_input/generate_doc_from_comments.sh +++ b/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 diff --git a/python/solid_dmft/gw_embedding/bdft_converter.py b/python/solid_dmft/gw_embedding/bdft_converter.py index adc9a2b6..c2cb8f2b 100644 --- a/python/solid_dmft/gw_embedding/bdft_converter.py +++ b/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 diff --git a/python/solid_dmft/gw_embedding/iaft.py b/python/solid_dmft/gw_embedding/iaft.py index 568d2141..062d37dc 100644 --- a/python/solid_dmft/gw_embedding/iaft.py +++ b/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. diff --git a/python/solid_dmft/gw_embedding/qp_evs_to_eig.py b/python/solid_dmft/gw_embedding/qp_evs_to_eig.py index 547a0c1b..18350c57 100644 --- a/python/solid_dmft/gw_embedding/qp_evs_to_eig.py +++ b/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
") - 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
') + 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() diff --git a/requirements.txt b/requirements.txt index c5c035a1..00dec4af 100644 --- a/requirements.txt +++ b/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