fcdmft.gw.mol.ugw_exact module#
- class fcdmft.gw.mol.ugw_exact.UGWExact(mf, auxbasis=None)[source]#
Bases:
UGWAC- Attributes:
- nmo
- nocc
Methods
__call__(*args, **kwargs)Update the attributes of the current object.
ao2mo([mo_coeff])Transform density-fitting integral from AO to MO.
apply(fn, *args, **kwargs)Apply the fn to rest arguments: return fn(*args, **kwargs). The return value of method set is the object itself. This allows a series of functions/methods to be executed in pipe.
check_sanity()Check input of class/object attributes, check whether a class method is overwritten.
copy()Returns a shallow copy
energy_tot([nw])Calculate GW total energy using Galitskii-Migdal formula.
get_ef([mo_energy])Get Fermi level.
get_frozen_mask()Get boolean mask for the unrestricted reference orbitals.
get_sigma_exchange(mo_coeff)Get exchange self-energy (EXX).
initialize_df([auxbasis])Initialize density fitting.
kernel()Do one-shot spin-unrestricted GW calculation using analytical continuation.
loop_ao2mo([mo_coeff, spin, ijslicea, ijsliceb])Transform density-fitting integral from AO to MO by block.
make_diag_dos(omega, eta)Get density of states using diagonal self-energy.
make_gf(omega, eta[, fullsigma, mode])Get exact dynamical Green's function and self-energy.
make_rdm1([nw, mode, ao_repr])Get GW one-particle density matrix.
post_kernel(envs)A hook to be run after the main body of the kernel function.
pre_kernel(envs)A hook to be run before the main body of kernel function is executed.
run(*args, **kwargs)Call the kernel function of current object.
set(*args, **kwargs)Update the attributes of the current object.
setup_evaluation_grid([fallback_freqs, ...])Set up self-energy grid, aka freqs2, aka gw.freqs.
view(cls)New view of object with the same attributes.
dump_flags
get_nmo
get_nocc
- energy_tot(nw=60)[source]#
Calculate GW total energy using Galitskii-Migdal formula. V. M. Galitskii and A. B. Migdal, Zh. E ksp. Teor. Fiz. 34, 139~1958! @Sov. Phys. JETP 139, 96 ~1958!# Working equation: equation A5 in Phys. Rev. B 88, 075105
- Parameters:
- nw
int, optional number for imaginary frequency grids for integration, by default 60
- nw
- Returns:
- e_totdouble
GW total energy
- e_hfdouble
HF total energy
- e_cdouble
GW correlation energy
- make_diag_dos(omega, eta)#
Get density of states using diagonal self-energy. Equation 7 in doi.org/10.1021/acs.jctc.2c00617
- make_gf(omega, eta, fullsigma=True, mode='linear')#
Get exact dynamical Green’s function and self-energy.
Dynamical self-energy is evaluated as equation 78 in doi.org/10.1021/ct300648t
Two modes for solving Dyson equation “dyson” for using inverse Dyson equation. “linear” for G = G0 + G0 Sigma G0, as equation 16 in doi.org/10.1021/acs.jctc.0c01264
- Parameters:
- gwUGWExact
gw object, provides attributes: nocc, nmo, _scf.mo_energy, exci, rho, vk, vxc
- omegadouble or complex
array frequency grids
- etadouble
broadening parameter
- fullsigmabool, optional
calculate off-diagonal elements, by default True
- mode
str, optional mode for Dyson equation, ‘linear’ or ‘dyson’, by default ‘linear’
- Returns:
- fcdmft.gw.mol.ugw_exact.diagonalize_phrpa(nocc, mo_energy, Lpq, RPAE=False)[source]#
Diagonalize particle-hole RPA matrix. The phRPA equation is solved for the response function for the screening interaction, which is defined in equation 66-68 in doi.org/10.1021/ct300648t The phRPA equation is reformulated as (A-B)(A+B)|X+Y>= omega^2|X+Y>, which is equation 15 in 10.1063/1.477483 RPAE equations are defined as equation 2 and 3 in 10.1103/PhysRevA.85.042507
- fcdmft.gw.mol.ugw_exact.get_sigma(nocc, mo_energy, mo_energy_prev, exci, rho, eta=1e-05, fullsigma=False, mode='b')[source]#
Get the real part of the GW correlation self-energy. Equation 83 in doi.org/10.1021/ct300648t mode ‘a’ and ‘b’ correspond to equation 10 and 11 in doi.org/10.1103/PhysRevB.76.165106
- Parameters:
- nocc
int1darray the number of occupied orbitals
- mo_energydouble 2d
array orbital energy
- mo_energy_prevdouble 2d
array orbital energy in previous iteration
- excidouble 1d
array phRPA excitation energy
- rhodouble 4d
array transition density
- etadouble, optional
broadening parameter, by default 1.0e-5
- fullsigmabool, optional
calculate off-diagonal elements, by default False
- mode
str, optional mode for off-diagonal elements, by default “b”
- nocc
- Returns:
- sigmadouble 3d
array real part of the GW static correlation self-energy
- sigmadouble 3d
- fcdmft.gw.mol.ugw_exact.get_sigma_derivative(nocc, mo_energy, mo_energy_prev, exci, rho, eta=1e-05)[source]#
Get the first-order derivative of the self-energy to the frequency. Equation.84 in doi.org/10.1021/ct300648t
- Parameters:
- Returns:
- derivativedouble 2d
ndarray first-order derivative of the correlation self-energy to the frequency
- derivativedouble 2d
- fcdmft.gw.mol.ugw_exact.get_transition_density(nocc, xpy, Lpq)[source]#
Calculate the transition density. Equation 85 in doi/abs/10.1021/ct300648t
- fcdmft.gw.mol.ugw_exact.make_diag_dos(gw, omega, eta)[source]#
Get density of states using diagonal self-energy. Equation 7 in doi.org/10.1021/acs.jctc.2c00617
- fcdmft.gw.mol.ugw_exact.make_gf(gw, omega, eta, fullsigma=True, mode='linear')[source]#
Get exact dynamical Green’s function and self-energy.
Dynamical self-energy is evaluated as equation 78 in doi.org/10.1021/ct300648t
Two modes for solving Dyson equation “dyson” for using inverse Dyson equation. “linear” for G = G0 + G0 Sigma G0, as equation 16 in doi.org/10.1021/acs.jctc.0c01264
- Parameters:
- gwUGWExact
gw object, provides attributes: nocc, nmo, _scf.mo_energy, exci, rho, vk, vxc
- omegadouble or complex
array frequency grids
- etadouble
broadening parameter
- fullsigmabool, optional
calculate off-diagonal elements, by default True
- mode
str, optional mode for Dyson equation, ‘linear’ or ‘dyson’, by default ‘linear’
- Returns: