From a033b14b818a25723b1af2b222b3df8519b65556 Mon Sep 17 00:00:00 2001 From: Qiming Sun Date: Mon, 26 Aug 2024 20:42:41 -0700 Subject: [PATCH] Update tdghf, tdgks, tddhf, tddks --- pyscf/tdscf/dhf.py | 32 +++++---- pyscf/tdscf/dks.py | 4 -- pyscf/tdscf/ghf.py | 117 ++++++++++++++++++++------------- pyscf/tdscf/gks.py | 32 ++++----- pyscf/tdscf/rhf.py | 8 +-- pyscf/tdscf/rks.py | 4 +- pyscf/tdscf/test/test_tdrhf.py | 23 +++++++ pyscf/tdscf/test/test_tdrks.py | 12 ---- pyscf/tdscf/test/test_tduhf.py | 16 ++++- pyscf/tdscf/test/test_tduks.py | 9 +-- pyscf/tdscf/uhf.py | 2 - pyscf/tdscf/uks.py | 5 +- 12 files changed, 145 insertions(+), 119 deletions(-) diff --git a/pyscf/tdscf/dhf.py b/pyscf/tdscf/dhf.py index 2aab4e6a8c..6e93d11135 100644 --- a/pyscf/tdscf/dhf.py +++ b/pyscf/tdscf/dhf.py @@ -23,12 +23,11 @@ from functools import reduce import numpy from pyscf import lib -from pyscf import dft +from pyscf import scf from pyscf import ao2mo from pyscf.lib import logger from pyscf.tdscf import rhf -from pyscf.tdscf._lr_eig import eig as lr_eig -from pyscf.data import nist +from pyscf.tdscf._lr_eig import eigh as lr_eigh, eig as lr_eig from pyscf import __config__ OUTPUT_THRESHOLD = getattr(__config__, 'tdscf_uhf_get_nto_threshold', 0.3) @@ -121,7 +120,7 @@ def add_hf_(a, b, hyb=1): b = b - numpy.einsum('jaib->iajb', eri_mo[:nocc,nocc:,:nocc,nocc:]) * hyb return a, b - if isinstance(mf, dft.KohnShamDFT): + if isinstance(mf, scf.hf.KohnShamDFT): from pyscf.dft import xc_deriv ni = mf._numint ni.libxc.test_deriv_order(mf.xc, 2, raise_error=True) @@ -399,6 +398,9 @@ def nuc_grad_method(self): @lib.with_doc(rhf.TDA.__doc__) class TDA(TDBase): + + singlet = None + def gen_vind(self, mf=None): '''Generate function to compute Ax''' if mf is None: @@ -450,14 +452,10 @@ def pickeig(w, v, nroots, envs): if x0 is None: x0 = self.init_guess(self._scf, self.nstates) - # FIXME: Is it correct to call davidson1 for complex integrals? - self.converged, self.e, x1 = \ - lib.davidson1(vind, x0, precond, - tol=self.conv_tol, - nroots=nstates, lindep=self.lindep, - max_cycle=self.max_cycle, - max_space=self.max_space, pick=pickeig, - verbose=log) + self.converged, self.e, x1 = lr_eigh( + vind, x0, precond, tol_residual=self.conv_tol, lindep=self.lindep, + nroots=nstates, pick=pickeig, max_cycle=self.max_cycle, + max_memory=self.max_memory, verbose=log) nocc = (self._scf.mo_occ>0).sum() nmo = self._scf.mo_occ.size @@ -527,7 +525,8 @@ def vind(xys): class TDHF(TDBase): - conv_tol = 1e-5 + conv_tol = 1e-4 + singlet = None @lib.with_doc(gen_tdhf_operation.__doc__) def gen_vind(self, mf=None): @@ -569,9 +568,9 @@ def pickeig(w, v, nroots, envs): x0 = self.init_guess(self._scf, self.nstates) self.converged, w, x1 = lr_eig( - vind, x0, precond, tol_residual=self.conv_tol, nroots=nstates, - lindep=self.lindep, max_cycle=self.max_cycle, - max_space=self.max_space, pick=pickeig, verbose=log) + vind, x0, precond, tol_residual=self.conv_tol, lindep=self.lindep, + nroots=nstates, pick=pickeig, max_cycle=self.max_cycle, + max_memory=self.max_memory, verbose=log) nocc = (self._scf.mo_occ>0).sum() nmo = self._scf.mo_occ.size @@ -593,7 +592,6 @@ def norm_xy(z): self._finalize() return self.e, self.xy -from pyscf import scf scf.dhf.DHF.TDA = scf.dhf.RDHF.TDA = lib.class_as_method(TDA) scf.dhf.DHF.TDHF = scf.dhf.RDHF.TDHF = lib.class_as_method(TDHF) diff --git a/pyscf/tdscf/dks.py b/pyscf/tdscf/dks.py index 1f7bdf0563..fec0ad1e85 100644 --- a/pyscf/tdscf/dks.py +++ b/pyscf/tdscf/dks.py @@ -17,12 +17,8 @@ # -import numpy from pyscf import lib from pyscf.tdscf import dhf -from pyscf.data import nist -from pyscf.dft.rks import KohnShamDFT -from pyscf import __config__ class TDA(dhf.TDA): diff --git a/pyscf/tdscf/ghf.py b/pyscf/tdscf/ghf.py index 69f70ea13b..52071f7468 100644 --- a/pyscf/tdscf/ghf.py +++ b/pyscf/tdscf/ghf.py @@ -25,14 +25,13 @@ from functools import reduce import numpy from pyscf import lib -from pyscf import dft -from pyscf.dft import numint +from pyscf import scf from pyscf import ao2mo from pyscf import symm from pyscf.lib import logger from pyscf.tdscf import rhf +from pyscf.tdscf._lr_eig import eigh as lr_eigh, eig as lr_eig from pyscf.scf import ghf_symm -from pyscf.data import nist from pyscf import __config__ OUTPUT_THRESHOLD = getattr(__config__, 'tdscf_rhf_get_nto_threshold', 0.3) @@ -62,9 +61,7 @@ def gen_tda_operation(mf, fock_ao=None, wfnsym=None): if isinstance(wfnsym, str): wfnsym = symm.irrep_name2id(mol.groupname, wfnsym) wfnsym = wfnsym % 10 # convert to D2h subgroup - orbsym = ghf_symm.get_orbsym(mol, mo_coeff) - orbsym_in_d2h = numpy.asarray(orbsym) % 10 # convert to D2h irreps - sym_forbid = (orbsym_in_d2h[occidx,None] ^ orbsym_in_d2h[viridx]) != wfnsym + sym_forbid = _get_x_sym_table(mf) != wfnsym if fock_ao is None: e_ia = hdiag = mo_energy[viridx] - mo_energy[occidx,None] @@ -102,6 +99,14 @@ def vind(zs): return vind, hdiag gen_tda_hop = gen_tda_operation +def _get_x_sym_table(mf): + '''Irrep (up to D2h symmetry) of each coefficient in X[nocc,nvir]''' + mol = mf.mol + mo_occ = mf.mo_occ + orbsym = ghf_symm.get_orbsym(mol, mf.mo_coeff) + orbsym = numpy.asarray(orbsym) % 10 # convert to D2h irreps + return orbsym[mo_occ==1,None] ^ orbsym[mo_occ==0] + def get_ab(mf, mo_energy=None, mo_coeff=None, mo_occ=None): r'''A and B matrices for TDDFT response function. @@ -151,7 +156,7 @@ def add_hf_(a, b, hyb=1): b -= numpy.einsum('jaib->iajb', eri_mo[:nocc,nocc:,:nocc,nocc:]) * hyb return a, b - if isinstance(mf, dft.KohnShamDFT): + if isinstance(mf, scf.hf.KohnShamDFT): from pyscf.dft import xc_deriv ni = mf._numint ni.libxc.test_deriv_order(mf.xc, 2, raise_error=True) @@ -379,7 +384,7 @@ def gen_vind(self, mf=None): mf = self._scf return gen_tda_hop(mf, wfnsym=self.wfnsym) - def init_guess(self, mf, nstates=None, wfnsym=None): + def init_guess(self, mf, nstates=None, wfnsym=None, return_symmetry=False): if nstates is None: nstates = self.nstates if wfnsym is None: wfnsym = self.wfnsym @@ -387,27 +392,36 @@ def init_guess(self, mf, nstates=None, wfnsym=None): mo_occ = mf.mo_occ occidx = numpy.where(mo_occ==1)[0] viridx = numpy.where(mo_occ==0)[0] - e_ia = mo_energy[viridx] - mo_energy[occidx,None] - - if wfnsym is not None and mf.mol.symmetry: - if isinstance(wfnsym, str): - wfnsym = symm.irrep_name2id(mf.mol.groupname, wfnsym) - wfnsym = wfnsym % 10 # convert to D2h subgroup - orbsym = ghf_symm.get_orbsym(mf.mol, mf.mo_coeff) - orbsym_in_d2h = numpy.asarray(orbsym) % 10 # convert to D2h irreps - e_ia[(orbsym_in_d2h[occidx,None] ^ orbsym_in_d2h[viridx]) != wfnsym] = 1e99 - + e_ia = (mo_energy[viridx] - mo_energy[occidx,None]).ravel() nov = e_ia.size nstates = min(nstates, nov) - e_ia = e_ia.ravel() - e_threshold = numpy.sort(e_ia)[nstates-1] + + if (wfnsym is not None or return_symmetry) and mf.mol.symmetry: + x_sym = _get_x_sym_table(mf).ravel() + if wfnsym is not None: + if isinstance(wfnsym, str): + wfnsym = symm.irrep_name2id(mf.mol.groupname, wfnsym) + wfnsym = wfnsym % 10 # convert to D2h subgroup + e_ia[x_sym != wfnsym] = 1e99 + nov_allowed = numpy.count_nonzero(x_sym == wfnsym) + nstates = min(nstates, nov_allowed) + + e_threshold = numpy.partition(e_ia, nstates-1)[nstates-1] e_threshold += self.deg_eia_thresh idx = numpy.where(e_ia <= e_threshold)[0] x0 = numpy.zeros((idx.size, nov)) for i, j in enumerate(idx): x0[i, j] = 1 # Koopmans' excitations - return x0 + + if return_symmetry: + if mf.mol.symmetry: + x0sym = x_sym[idx] + else: + x0sym = None + return x0, x0sym + else: + return x0 def kernel(self, x0=None, nstates=None): '''TDA diagonalization solver @@ -419,6 +433,7 @@ def kernel(self, x0=None, nstates=None): nstates = self.nstates else: self.nstates = nstates + mol = self.mol log = logger.Logger(self.stdout, self.verbose) @@ -429,17 +444,18 @@ def pickeig(w, v, nroots, envs): idx = numpy.where(w > self.positive_eig_threshold)[0] return w[idx], v[:,idx], idx + x0sym = None if x0 is None: - x0 = self.init_guess(self._scf, self.nstates) + x0, x0sym = self.init_guess( + self._scf, self.nstates, return_symmetry=True) + elif mol.symmetry: + x_sym = _get_x_sym_table(self._scf).ravel() + x0sym = [rhf._guess_wfnsym_id(self, x_sym, x) for x in x0] - # FIXME: Is it correct to call davidson1 for complex integrals - self.converged, self.e, x1 = \ - lib.davidson1(vind, x0, precond, - tol=self.conv_tol, - nroots=nstates, lindep=self.lindep, - max_cycle=self.max_cycle, - max_space=self.max_space, pick=pickeig, - verbose=log) + self.converged, self.e, x1 = lr_eigh( + vind, x0, precond, tol_residual=self.conv_tol, lindep=self.lindep, + nroots=nstates, x0sym=x0sym, pick=pickeig, max_cycle=self.max_cycle, + max_memory=self.max_memory, verbose=log) nocc = (self._scf.mo_occ>0).sum() nmo = self._scf.mo_occ.size @@ -479,9 +495,7 @@ def gen_tdhf_operation(mf, fock_ao=None, wfnsym=None): if isinstance(wfnsym, str): wfnsym = symm.irrep_name2id(mol.groupname, wfnsym) wfnsym = wfnsym % 10 # convert to D2h subgroup - orbsym = ghf_symm.get_orbsym(mol, mo_coeff) - orbsym_in_d2h = numpy.asarray(orbsym) % 10 # convert to D2h irreps - sym_forbid = (orbsym_in_d2h[occidx,None] ^ orbsym_in_d2h[viridx]) != wfnsym + sym_forbid = _get_x_sym_table(mf) != wfnsym assert fock_ao is None @@ -530,6 +544,7 @@ def vind(xys): class TDHF(TDBase): + conv_tol = 1e-4 singlet = None @lib.with_doc(gen_tdhf_operation.__doc__) @@ -538,10 +553,15 @@ def gen_vind(self, mf=None): mf = self._scf return gen_tdhf_operation(mf, wfnsym=self.wfnsym) - def init_guess(self, mf, nstates=None, wfnsym=None): - x0 = TDA.init_guess(self, mf, nstates, wfnsym) - y0 = numpy.zeros_like(x0) - return numpy.asarray(numpy.block([[x0, y0], [y0, x0.conj()]])) + def init_guess(self, mf, nstates=None, wfnsym=None, return_symmetry=False): + if return_symmetry: + x0, x0sym = TDA.init_guess(self, mf, nstates, wfnsym, return_symmetry) + y0 = numpy.zeros_like(x0) + return numpy.hstack([x0, y0]), x0sym + else: + x0 = TDA.init_guess(self, mf, nstates, wfnsym, return_symmetry) + y0 = numpy.zeros_like(x0) + return numpy.hstack([x0, y0]) def kernel(self, x0=None, nstates=None): '''TDHF diagonalization with non-Hermitian eigenvalue solver @@ -553,6 +573,7 @@ def kernel(self, x0=None, nstates=None): nstates = self.nstates else: self.nstates = nstates + mol = self.mol log = logger.Logger(self.stdout, self.verbose) @@ -566,16 +587,19 @@ def pickeig(w, v, nroots, envs): # FIXME: Should the amplitudes be real? It also affects x2c-tdscf return lib.linalg_helper._eigs_cmplx2real(w, v, realidx, ensure_real) + x0sym = None if x0 is None: - x0 = self.init_guess(self._scf, self.nstates) - - self.converged, w, x1 = \ - lib.davidson_nosym1(vind, x0, precond, - tol=self.conv_tol, - nroots=nstates, lindep=self.lindep, - max_cycle=self.max_cycle, - max_space=self.max_space, pick=pickeig, - verbose=log) + x0, x0sym = self.init_guess( + self._scf, self.nstates, return_symmetry=True) + elif mol.symmetry: + x_sym = y_sym = _get_x_sym_table(self._scf).ravel() + x_sym = numpy.append(x_sym, y_sym) + x0sym = [rhf._guess_wfnsym_id(self, x_sym, x) for x in x0] + + self.converged, w, x1 = lr_eig( + vind, x0, precond, tol_residual=self.conv_tol, lindep=self.lindep, + nroots=nstates, x0sym=x0sym, pick=pickeig, max_cycle=self.max_cycle, + max_memory=self.max_memory, verbose=log) nocc = (self._scf.mo_occ>0).sum() nmo = self._scf.mo_occ.size @@ -598,6 +622,5 @@ def norm_xy(z): RPA = TDGHF = TDHF -from pyscf import scf scf.ghf.GHF.TDA = lib.class_as_method(TDA) scf.ghf.GHF.TDHF = lib.class_as_method(TDHF) diff --git a/pyscf/tdscf/gks.py b/pyscf/tdscf/gks.py index 461fc60635..3d372a6a08 100644 --- a/pyscf/tdscf/gks.py +++ b/pyscf/tdscf/gks.py @@ -20,10 +20,8 @@ import numpy from pyscf import lib from pyscf import symm -from pyscf.tdscf import ghf -from pyscf.scf import ghf_symm -from pyscf.scf import _response_functions # noqa -from pyscf.data import nist +from pyscf.tdscf import ghf, rhf +from pyscf.tdscf._lr_eig import eigh as lr_eigh from pyscf.dft.rks import KohnShamDFT from pyscf import __config__ @@ -39,6 +37,7 @@ class TDDFT(ghf.TDHF): class CasidaTDDFT(TDDFT, TDA): '''Solve the Casida TDDFT formula (A-B)(A+B)(X+Y) = (X+Y)w^2 ''' + init_guess = TDA.init_guess def gen_vind(self, mf=None): @@ -62,9 +61,7 @@ def gen_vind(self, mf=None): if isinstance(wfnsym, str): wfnsym = symm.irrep_name2id(mol.groupname, wfnsym) wfnsym = wfnsym % 10 # convert to D2h subgroup - orbsym = ghf_symm.get_orbsym(mol, mo_coeff) - orbsym_in_d2h = numpy.asarray(orbsym) % 10 # convert to D2h irreps - sym_forbid = (orbsym_in_d2h[occidx,None] ^ orbsym_in_d2h[viridx]) != wfnsym + sym_forbid = ghf._get_x_sym_table(mf) != wfnsym e_ia = (mo_energy[viridx].reshape(-1,1) - mo_energy[occidx]).T if wfnsym is not None and mol.symmetry: @@ -101,6 +98,7 @@ def kernel(self, x0=None, nstates=None): '''TDDFT diagonalization solver ''' cpu0 = (lib.logger.process_clock(), lib.logger.perf_counter()) + mol = self.mol mf = self._scf if mf._numint.libxc.is_hybrid_xc(mf.xc): raise RuntimeError('%s cannot be used with hybrid functional' @@ -121,16 +119,18 @@ def pickeig(w, v, nroots, envs): idx = numpy.where(w > self.positive_eig_threshold)[0] return w[idx], v[:,idx], idx + x0sym = None if x0 is None: - x0 = self.init_guess(self._scf, self.nstates) - - self.converged, w2, x1 = \ - lib.davidson1(vind, x0, precond, - tol=self.conv_tol, - nroots=nstates, lindep=self.lindep, - max_cycle=self.max_cycle, - max_space=self.max_space, pick=pickeig, - verbose=log) + x0, x0sym = self.init_guess( + self._scf, self.nstates, return_symmetry=True) + elif mol.symmetry: + x_sym = ghf._get_x_sym_table(self._scf).ravel() + x0sym = [rhf._guess_wfnsym_id(self, x_sym, x) for x in x0] + + self.converged, w2, x1 = lr_eigh( + vind, x0, precond, tol_residual=self.conv_tol, lindep=self.lindep, + nroots=nstates, x0sym=x0sym, pick=pickeig, max_cycle=self.max_cycle, + max_memory=self.max_memory, verbose=log) mo_energy = self._scf.mo_energy mo_occ = self._scf.mo_occ diff --git a/pyscf/tdscf/rhf.py b/pyscf/tdscf/rhf.py index e77da0a5bd..8dd0bf9fee 100644 --- a/pyscf/tdscf/rhf.py +++ b/pyscf/tdscf/rhf.py @@ -32,9 +32,9 @@ from pyscf import symm from pyscf.lib import logger from pyscf.scf import hf_symm -from pyscf.scf import _response_functions +from pyscf.scf import _response_functions # noqa from pyscf.data import nist -from pyscf.tdscf._lr_eig import eig as lr_eig, eigh as lr_eigh +from pyscf.tdscf._lr_eig import eigh as lr_eigh, eig as lr_eig from pyscf import __config__ OUTPUT_THRESHOLD = getattr(__config__, 'tdscf_rhf_get_nto_threshold', 0.3) @@ -655,7 +655,7 @@ def __call__(self, mol_or_geom, **kwargs): class TDBase(lib.StreamObject): - conv_tol = getattr(__config__, 'tdscf_rhf_TDA_conv_tol', 1e-9) + conv_tol = getattr(__config__, 'tdscf_rhf_TDA_conv_tol', 1e-5) nstates = getattr(__config__, 'tdscf_rhf_TDA_nstates', 3) singlet = getattr(__config__, 'tdscf_rhf_TDA_singlet', True) lindep = getattr(__config__, 'tdscf_rhf_TDA_lindep', 1e-12) @@ -804,8 +804,6 @@ class TDA(TDBase): normalized to 1 for UHF/UKS methods. In the TDA calculation, Y = 0. ''' - conv_tol = 1e-5 - def gen_vind(self, mf=None): '''Generate function to compute Ax''' if mf is None: diff --git a/pyscf/tdscf/rks.py b/pyscf/tdscf/rks.py index 2e53a21e1a..5b22e06449 100644 --- a/pyscf/tdscf/rks.py +++ b/pyscf/tdscf/rks.py @@ -26,8 +26,6 @@ from pyscf import symm from pyscf.tdscf import rhf from pyscf.tdscf._lr_eig import eigh as lr_eigh -from pyscf.scf import hf_symm -from pyscf.data import nist from pyscf.dft.rks import KohnShamDFT from pyscf import __config__ @@ -47,7 +45,7 @@ def nuc_grad_method(self): class CasidaTDDFT(TDDFT, TDA): '''Solve the Casida TDDFT formula (A-B)(A+B)(X+Y) = (X+Y)w^2 ''' - conv_tol = 1e-5 + init_guess = TDA.init_guess def gen_vind(self, mf=None): diff --git a/pyscf/tdscf/test/test_tdrhf.py b/pyscf/tdscf/test/test_tdrhf.py index 6f0296af9b..776006e44e 100644 --- a/pyscf/tdscf/test/test_tdrhf.py +++ b/pyscf/tdscf/test/test_tdrhf.py @@ -64,6 +64,29 @@ def test_tdhf_triplet(self): ref = [10.8919234, 10.8919234, 12.63440705] self.assertAlmostEqual(abs(e[:len(ref)] * 27.2114 - ref).max(), 0, 5) + def test_symmetry_init_guess(self): + mol = gto.M(atom='N 0 0 0; N 0 0 1.2', basis='631g', symmetry='D2h') + mf = mol.RHF.run() + td = mf.TDA().run(nstates=1) + self.assertAlmostEqual(td.e[0], 0.22349707455528, 7) + mo_coeff = mf.mo_coeff + mo_occ = mf.mo_occ + orbsym = tdscf.rhf.hf_symm.get_orbsym(mol, mo_coeff) + x_sym = tdscf.rhf.symm.direct_prod(orbsym[mo_occ==2], orbsym[mo_occ==0], mol.groupname) + wfnsym = tdscf.rhf._analyze_wfnsym(td, x_sym, td.xy[0][0]) + self.assertEqual(wfnsym, 'Au') + + mol = gto.M(atom='N 0 0 0; N 0 0 1.2', basis='631g', symmetry=True) + mf = mol.RHF.run() + td = mf.TDA().run(nstates=1, singlet=False) + self.assertAlmostEqual(td.e[0], 0.14147328219131602, 7) + mo_coeff = mf.mo_coeff + mo_occ = mf.mo_occ + orbsym = tdscf.rhf.hf_symm.get_orbsym(mol, mo_coeff) + x_sym = tdscf.rhf.symm.direct_prod(orbsym[mo_occ==2], orbsym[mo_occ==0], mol.groupname) + wfnsym = tdscf.rhf._analyze_wfnsym(td, x_sym, td.xy[0][0]) + self.assertEqual(wfnsym, 'A1u') + if __name__ == "__main__": print("Full Tests for rhf-TDA and rhf-TDHF") diff --git a/pyscf/tdscf/test/test_tdrks.py b/pyscf/tdscf/test/test_tdrks.py index 691f49bf5d..7bc7a228d4 100644 --- a/pyscf/tdscf/test/test_tdrks.py +++ b/pyscf/tdscf/test/test_tdrks.py @@ -473,18 +473,6 @@ def test_custom_rsh(self): ref = [16.14837289, 28.01968627, 49.00854076] self.assertAlmostEqual(abs(e_td*nist.HARTREE2EV - ref).max(), 0, 4) - def test_symmetry_init_guess(self): - mol = gto.M(atom='N 0 0 0; N 0 0 1.2', basis='631g', output='/dev/null', symmetry='D2h') - mf = mol.RHF.run() - td = mf.TDA().run(nstates=1) - self.assertAlmostEqual(td.e[0], 0.22349707455528, 7) - mo_coeff = mf.mo_coeff - mo_occ = mf.mo_occ - orbsym = rhf.hf_symm.get_orbsym(mol, mo_coeff) - x_sym = rhf.symm.direct_prod(orbsym[mo_occ==2], orbsym[mo_occ==0], mol.groupname) - wfnsym = rhf._analyze_wfnsym(td, x_sym, td.xy[0][0]) - self.assertEqual(wfnsym, 'Au') - if __name__ == "__main__": print("Full Tests for TD-RKS") unittest.main() diff --git a/pyscf/tdscf/test/test_tduhf.py b/pyscf/tdscf/test/test_tduhf.py index 8769765f87..08f8734bab 100644 --- a/pyscf/tdscf/test/test_tduhf.py +++ b/pyscf/tdscf/test/test_tduhf.py @@ -15,7 +15,7 @@ # import unittest -from pyscf import lib, gto, scf, tdscf +from pyscf import lib, gto, scf, tdscf, symm def setUpModule(): global mol, mol1, mf, mf1 @@ -76,6 +76,20 @@ def test_tdhf_triplet(self): ref = [3.31267103, 18.4954748, 20.84935404, 21.54808392] self.assertAlmostEqual(abs(e * 27.2114 - ref).max(), 0, 4) + def test_symmetry_init_guess(self): + mol = gto.M(atom='N 0 0 0; N 0 0 1.2', basis='631g', symmetry=True, verbose=0) + mf = mol.UHF.run() + td = mf.TDA().run(nstates=1) + self.assertAlmostEqual(td.e[0], 0.14147328219131602, 7) + mo_coeff = mf.mo_coeff + mo_occa, mo_occb = mf.mo_occ + orbsyma, orbsymb = scf.uhf_symm.get_orbsym(mol, mo_coeff) + x_syma = symm.direct_prod(orbsyma[mo_occa==1], orbsyma[mo_occa==0], mol.groupname) + x_symb = symm.direct_prod(orbsymb[mo_occb==1], orbsymb[mo_occb==0], mol.groupname) + wfnsyma = tdscf.rhf._analyze_wfnsym(td, x_syma, td.xy[0][0][0]) + wfnsymb = tdscf.rhf._analyze_wfnsym(td, x_symb, td.xy[0][0][1]) + self.assertAlmostEqual(wfnsyma, 'A1u') + self.assertAlmostEqual(wfnsymb, 'A1u') if __name__ == "__main__": print("Full Tests for uhf-TDA and uhf-TDHF") diff --git a/pyscf/tdscf/test/test_tduks.py b/pyscf/tdscf/test/test_tduks.py index 9f0c921d72..5c7aa31851 100644 --- a/pyscf/tdscf/test/test_tduks.py +++ b/pyscf/tdscf/test/test_tduks.py @@ -18,7 +18,7 @@ import unittest import numpy -from pyscf import lib, gto, scf, dft +from pyscf import lib, gto, scf, dft, symm from pyscf import tdscf def diagonalize(a, b, nroots=4): @@ -470,13 +470,6 @@ def test_tddft_with_wfnsym(self): self.assertAlmostEqual(lib.fp(es), 0.15403661700414412, 6) td.analyze() - def test_symmetry_init_guess(self): - mol = gto.M(atom='N 0 0 0; N 0 0 1.2', basis='631g', symmetry=True, verbose=0) - td = mol.UHF.run().TDA().run(nstates=1) - self.assertAlmostEqual(td.e[0], 0.14147328219131602, 7) - # TODO: verify symmetry of td.x == A1u, close to triplet state - - if __name__ == "__main__": print("Full Tests for TD-UKS") unittest.main() diff --git a/pyscf/tdscf/uhf.py b/pyscf/tdscf/uhf.py index f04b716fb9..abdad81c31 100644 --- a/pyscf/tdscf/uhf.py +++ b/pyscf/tdscf/uhf.py @@ -24,7 +24,6 @@ from pyscf import ao2mo from pyscf.lib import logger from pyscf.scf import uhf_symm -from pyscf.scf import _response_functions from pyscf.tdscf import rhf from pyscf.tdscf._lr_eig import eigh as lr_eigh, eig as lr_eig from pyscf.data import nist @@ -619,7 +618,6 @@ def nuc_grad_method(self): @lib.with_doc(rhf.TDA.__doc__) class TDA(TDBase): - conv_tol = 1e-5 singlet = None def gen_vind(self, mf=None): diff --git a/pyscf/tdscf/uks.py b/pyscf/tdscf/uks.py index c2af5e34a8..25c49e2ffc 100644 --- a/pyscf/tdscf/uks.py +++ b/pyscf/tdscf/uks.py @@ -21,10 +21,8 @@ from pyscf import symm from pyscf import lib from pyscf.lib import logger -from pyscf.tdscf import rks, uhf, rhf +from pyscf.tdscf import uhf, rhf from pyscf.tdscf._lr_eig import eigh as lr_eigh -from pyscf.scf import uhf_symm -from pyscf.data import nist from pyscf.dft.rks import KohnShamDFT from pyscf import __config__ @@ -46,7 +44,6 @@ class CasidaTDDFT(TDDFT, TDA): '''Solve the Casida TDDFT formula (A-B)(A+B)(X+Y) = (X+Y)w^2 ''' - conv_tol = 1e-5 init_guess = TDA.init_guess def gen_vind(self, mf=None):