diff --git a/pyscf/pbc/tdscf/krhf.py b/pyscf/pbc/tdscf/krhf.py index 30e9557f8d..afad46984f 100644 --- a/pyscf/pbc/tdscf/krhf.py +++ b/pyscf/pbc/tdscf/krhf.py @@ -45,6 +45,11 @@ def get_ab(mf, kshift=0): B[i,a,j,b] = (ia||jb) Ref: Chem Phys Lett, 256, 454 + + Kwargs: + kshift : integer + The index of the k-point that represents the transition between + k-points in the excitation coefficients. ''' cell = mf.cell mo_energy = scf.addons.mo_energy_with_exxdiv_none(mf) @@ -155,9 +160,17 @@ def add_hf_(a, b, hyb=1): return a, b class KTDBase(TDBase): + ''' + Attributes: + kshift_lst : list of integers + Each element in the list is the index of the k-point that + represents the transition between k-points in the excitation + coefficients. + ''' + conv_tol = getattr(__config__, 'pbc_tdscf_rhf_TDA_conv_tol', 1e-4) - _keys = {'kconserv', 'kshift_lst'} + _keys = {'kshift_lst'} def __init__(self, mf, kshift_lst=None): assert isinstance(mf, scf.khf.KSCF) @@ -165,8 +178,6 @@ def __init__(self, mf, kshift_lst=None): warn_pbc2d_eri(mf) if kshift_lst is None: kshift_lst = [0] - - self.kconserv = get_kconserv_ria(mf.cell, mf.kpts) self.kshift_lst = kshift_lst def dump_flags(self, verbose=None): @@ -222,12 +233,16 @@ def get_ab(self, mf=None, kshift=0): if mf is None: mf = self._scf return get_ab(mf, kshift) - def gen_vind(self, mf, kshift): - # exxdiv corrections are kept in hdiag while excluding them when calling - # the contractions between two-electron integrals and X/Y amplitudes. - # See also the relevant comments in function pbc.tdscf.rhf.TDA.gen_vind + def gen_vind(self, mf, kshift=0): + '''Compute Ax + + Kwargs: + kshift : integer + The index of the k-point that represents the transition between + k-points in the excitation coefficients. + ''' singlet = self.singlet - kconserv = self.kconserv[kshift] + kconserv = get_kconserv_ria(mf.cell, mf.kpts)[kshift] mo_coeff = mf.mo_coeff mo_occ = mf.mo_occ @@ -273,7 +288,7 @@ def init_guess(self, mf, kshift, nstates=None): mo_energy = mf.mo_energy mo_occ = mf.mo_occ - kconserv = self.kconserv[kshift] + kconserv = get_kconserv_ria(mf.cell, mf.kpts)[kshift] e_ia = numpy.concatenate( [x.reshape(-1) for x in _get_e_ia(mo_energy, mo_occ, kconserv)] ) @@ -315,7 +330,7 @@ def pickeig(w, v, nroots, envs): self.e = [] self.xy = [] for i,kshift in enumerate(self.kshift_lst): - kconserv = self.kconserv[kshift] + kconserv = get_kconserv_ria(mf.cell, mf.kpts)[kshift] vind, hdiag = self.gen_vind(self._scf, kshift) precond = self.get_precond(hdiag) @@ -467,7 +482,7 @@ def norm_xy(z, kconserv): self.e = [] self.xy = [] for i,kshift in enumerate(self.kshift_lst): - kconserv = self.kconserv[kshift] + kconserv = get_kconserv_ria(mf.cell, mf.kpts)[kshift] vind, hdiag = self.gen_vind(self._scf, kshift) precond = self.get_precond(hdiag) diff --git a/pyscf/pbc/tdscf/krks.py b/pyscf/pbc/tdscf/krks.py index b38c8d309a..7ce5b6f9a0 100644 --- a/pyscf/pbc/tdscf/krks.py +++ b/pyscf/pbc/tdscf/krks.py @@ -57,47 +57,3 @@ def _rebuild_df(td): dft.kroks.KROKS.TDA = None dft.kroks.KROKS.TDHF = None dft.kroks.KROKS.TDDFT = None - - -if __name__ == '__main__': - from pyscf.pbc import gto - from pyscf.pbc import dft - cell = gto.Cell() - cell.unit = 'B' - cell.atom = ''' - C 0. 0. 0. - C 1.68506879 1.68506879 1.68506879 - ''' - cell.a = ''' - 0. 3.37013758 3.37013758 - 3.37013758 0. 3.37013758 - 3.37013758 3.37013758 0. - ''' - cell.basis = 'gth-szv' - cell.pseudo = 'gth-pade' - cell.mesh = [25]*3 - cell.build() - - mf = dft.KRKS(cell, cell.make_kpts([2,1,1])) - #mf.with_df = df.MDF(cell, cell.make_kpts([2,1,1])) - #mf.with_df.auxbasis = 'weigend' - #mf.with_df._cderi = 'eri3d-mdf.h5' - #mf.with_df.build(with_j3c=False) - mf.xc = 'lda,' - mf.kernel() -#mesh=12 -10.3077341607895 -#mesh=5 -10.3086623157515 - - td = TDDFT(mf) - td.nstates = 5 - td.verbose = 5 - print(td.kernel()[0][0] * 27.2114) -#mesh=12 [ 6.08108297 6.10231481 6.10231478 6.38355803 6.38355804] -#MDF mesh=5 [ 6.07919157 6.10251718 6.10253961 6.37202499 6.37565246] - - td = TDA(mf) - td.singlet = False - td.verbose = 5 - print(td.kernel()[0][0] * 27.2114) -#mesh=12 [ 4.01539192 5.1750807 5.17508071] -#MDF mesh=5 [ 4.01148649 5.18043397 5.18043459] diff --git a/pyscf/pbc/tdscf/kuhf.py b/pyscf/pbc/tdscf/kuhf.py index 86eb31c2d6..30e476635b 100644 --- a/pyscf/pbc/tdscf/kuhf.py +++ b/pyscf/pbc/tdscf/kuhf.py @@ -24,21 +24,26 @@ from pyscf.tdscf._lr_eig import eigh as lr_eigh, eig as lr_eig from pyscf.pbc import scf from pyscf.pbc.tdscf.krhf import KTDBase, _get_e_ia -from pyscf.pbc.lib.kpts_helper import is_gamma_point +from pyscf.pbc.lib.kpts_helper import is_gamma_point, get_kconserv_ria from pyscf.pbc.scf import _response_functions # noqa from pyscf import __config__ REAL_EIG_THRESHOLD = getattr(__config__, 'pbc_tdscf_uhf_TDDFT_pick_eig_threshold', 1e-3) class TDA(KTDBase): - conv_tol = getattr(__config__, 'pbc_tdscf_rhf_TDA_conv_tol', 1e-7) def get_ab(self, mf=None, kshift=0): raise NotImplementedError - def gen_vind(self, mf, kshift): - '''Compute Ax''' - kconserv = self.kconserv[kshift] + def gen_vind(self, mf, kshift=0): + '''Compute Ax + + Kwargs: + kshift : integer + The index of the k-point that represents the transition between + k-points in the excitation coefficients. + ''' + kconserv = get_kconserv_ria(mf.cell, mf.kpts)[kshift] mo_coeff = mf.mo_coeff mo_occ = mf.mo_occ nkpts = len(mo_occ[0]) @@ -98,7 +103,7 @@ def init_guess(self, mf, kshift, nstates=None): mo_energy = mf.mo_energy mo_occ = mf.mo_occ - kconserv = self.kconserv[kshift] + kconserv = get_kconserv_ria(mf.cell, mf.kpts)[kshift] e_ia_a = _get_e_ia(mo_energy[0], mo_occ[0], kconserv) e_ia_b = _get_e_ia(mo_energy[1], mo_occ[1], kconserv) e_ia = numpy.hstack([x.ravel() for x in (e_ia_a + e_ia_b)]) @@ -136,7 +141,7 @@ def pickeig(w, v, nroots, envs): self.e = [] self.xy = [] for i,kshift in enumerate(self.kshift_lst): - kconserv = self.kconserv[kshift] + kconserv = get_kconserv_ria(mf.cell, mf.kpts)[kshift] vind, hdiag = self.gen_vind(self._scf, kshift) precond = self.get_precond(hdiag) @@ -164,8 +169,6 @@ def pickeig(w, v, nroots, envs): class TDHF(KTDBase): - conv_tol = 1e-5 - def get_ab(self, mf=None, kshift=0): raise NotImplementedError @@ -281,7 +284,7 @@ def pickeig(w, v, nroots, envs): self.e = [] self.xy = [] for i,kshift in enumerate(self.kshift_lst): - kconserv = self.kconserv[kshift] + kconserv = get_kconserv_ria(mf.cell, mf.kpts)[kshift] vind, hdiag = self.gen_vind(self._scf, kshift) precond = self.get_precond(hdiag) diff --git a/pyscf/pbc/tdscf/kuks.py b/pyscf/pbc/tdscf/kuks.py index 1f76d12c8f..aee05f0f61 100644 --- a/pyscf/pbc/tdscf/kuks.py +++ b/pyscf/pbc/tdscf/kuks.py @@ -39,59 +39,3 @@ def kernel(self, x0=None): dft.kuks.KUKS.TDA = lib.class_as_method(KTDA) dft.kuks.KUKS.TDHF = None dft.kuks.KUKS.TDDFT = lib.class_as_method(TDDFT) - - -if __name__ == '__main__': - from pyscf.pbc import gto - from pyscf.pbc import dft - from pyscf.pbc import df - cell = gto.Cell() - cell.unit = 'B' - cell.atom = ''' - C 0. 0. 0. - C 1.68506879 1.68506879 1.68506879 - ''' - cell.a = ''' - 0. 3.37013758 3.37013758 - 3.37013758 0. 3.37013758 - 3.37013758 3.37013758 0. - ''' - - cell.basis = 'gth-szv' - cell.pseudo = 'gth-pade' - cell.mesh = [37]*3 - cell.build() - mf = dft.KUKS(cell, cell.make_kpts([2,1,1])).set(exxdiv=None, xc='b88,p86') - #mf.with_df = df.MDF(cell, cell.make_kpts([2,1,1])) - #mf.with_df.auxbasis = 'weigend' - #mf.with_df._cderi = 'eri3d-mdf.h5' - #mf.with_df.build(with_j3c=False) - mf.run() - - td = TDDFT(mf) - td.verbose = 5 - td.nstates = 5 - print(td.kernel()[0][0] * 27.2114) - - mf.xc = 'lda,vwn' - mf.run() - td = TDA(mf) - td.verbose = 5 - td.nstates = 5 - print(td.kernel()[0][0] * 27.2114) - - cell.spin = 2 - mf = dft.KUKS(cell, cell.make_kpts([2,1,1])).set(exxdiv=None, xc='b88,p86') - mf.run() - - td = TDDFT(mf) - td.verbose = 5 - td.nstates = 5 - print(td.kernel()[0][0] * 27.2114) - - mf.xc = 'lda,vwn' - mf.run() - td = TDA(mf) - td.verbose = 5 - td.nstates = 5 - print(td.kernel()[0][0] * 27.2114) diff --git a/pyscf/tdscf/test/test_tduks.py b/pyscf/tdscf/test/test_tduks.py index 5c7aa31851..3bd42317af 100644 --- a/pyscf/tdscf/test/test_tduks.py +++ b/pyscf/tdscf/test/test_tduks.py @@ -133,7 +133,6 @@ def test_tddft_b3lyp(self): def test_tddft_camb3lyp(self): mf = mol1.UKS(xc='camb3lyp').run() td = mf.TDDFT() - td.conv_tol = 1e-9 es = td.kernel(nstates=4)[0] a,b = td.get_ab() e_ref = diagonalize(a, b, 5)