From 96f5b7faaca2b1d2ad35372e19268cc737fc1e03 Mon Sep 17 00:00:00 2001 From: Matthew R Hermes Date: Thu, 19 Dec 2024 18:17:35 -0600 Subject: [PATCH] logoff commit lassi op_o1 prescreen linear dependencies not working; wrong energies --- my_pyscf/lassi/excitations.py | 3 ++- my_pyscf/lassi/lassi.py | 13 +++++++------ my_pyscf/lassi/op_o0.py | 4 +++- my_pyscf/lassi/op_o1/hams2ovlp.py | 2 +- my_pyscf/lassi/op_o1/stdm.py | 12 ++++++------ tests/lassi/test_4frag.py | 4 ++-- tests/lassi/test_opt57_slow.py | 4 ++-- 7 files changed, 23 insertions(+), 19 deletions(-) diff --git a/my_pyscf/lassi/excitations.py b/my_pyscf/lassi/excitations.py index d9105f7c..bc64003a 100644 --- a/my_pyscf/lassi/excitations.py +++ b/my_pyscf/lassi/excitations.py @@ -312,7 +312,8 @@ def get_ham_pq (self, h0, h1, h2, ci_p): nelec_frs = np.append (nelec_fs_p[:,None,:], nelec_frs_q, axis=1) with temporary_env (self, ncas_sub=norb_ref, mol=fcisolvers[0].mol): ham_pq, _, ovlp_pq = op[self.opt].ham (self, h1, h2, ci_fr, nelec_frs, soc=0, - orbsym=self.orbsym_ref, wfnsym=self.wfnsym_ref) + orbsym=self.orbsym_ref, + wfnsym=self.wfnsym_ref)[:3] t1 = self.log.timer ('get_ham_pq', *t0) return ham_pq + (h0*ovlp_pq) diff --git a/my_pyscf/lassi/lassi.py b/my_pyscf/lassi/lassi.py index 60632a5c..70b5b61f 100644 --- a/my_pyscf/lassi/lassi.py +++ b/my_pyscf/lassi/lassi.py @@ -356,7 +356,7 @@ def _eig_block (las, e0, h1, h2, ci_blk, nelec_blk, rootsym, soc, orbsym, wfnsym t0 = (lib.logger.process_clock (), lib.logger.perf_counter ()) if (las.verbose > lib.logger.INFO) and (o0_memcheck): ham_ref, s2_ref, ovlp_ref = op_o0.ham (las, h1, h2, ci_blk, nelec_blk, soc=soc, - orbsym=orbsym, wfnsym=wfnsym) + orbsym=orbsym, wfnsym=wfnsym)[:3] t0 = lib.logger.timer (las, 'LASSI diagonalizer rootsym {} CI algorithm'.format ( rootsym), *t0) @@ -364,8 +364,8 @@ def _eig_block (las, e0, h1, h2, ci_blk, nelec_blk, rootsym, soc, orbsym, wfnsym if soc: h1_sf = (h1[0:las.ncas,0:las.ncas] - h1[las.ncas:2*las.ncas,las.ncas:2*las.ncas]).real/2 - ham_blk, s2_blk, ovlp_blk = op[opt].ham (las, h1_sf, h2, ci_blk, nelec_blk, orbsym=orbsym, - wfnsym=wfnsym) + ham_blk, s2_blk, ovlp_blk, get_ovlp = op[opt].ham (las, h1_sf, h2, ci_blk, nelec_blk, + orbsym=orbsym, wfnsym=wfnsym) t0 = lib.logger.timer (las, 'LASSI diagonalizer rootsym {} TDM algorithm'.format ( rootsym), *t0) lib.logger.debug (las, @@ -388,8 +388,8 @@ def _eig_block (las, e0, h1, h2, ci_blk, nelec_blk, rootsym, soc, orbsym, wfnsym else: if (las.verbose > lib.logger.INFO): lib.logger.debug ( las, 'Insufficient memory to test against o0 LASSI algorithm') - ham_blk, s2_blk, ovlp_blk = op[opt].ham (las, h1, h2, ci_blk, nelec_blk, soc=soc, - orbsym=orbsym, wfnsym=wfnsym) + ham_blk, s2_blk, ovlp_blk, get_ovlp = op[opt].ham (las, h1, h2, ci_blk, nelec_blk, soc=soc, + orbsym=orbsym, wfnsym=wfnsym) t0 = lib.logger.timer (las, 'LASSI H build rootsym {}'.format (rootsym), *t0) log_debug = lib.logger.debug2 if las.nroots>10 else lib.logger.debug if np.iscomplexobj (ham_blk): @@ -430,7 +430,8 @@ def _eig_block (las, e0, h1, h2, ci_blk, nelec_blk, rootsym, soc, orbsym, wfnsym lc = 'checking if LASSI basis has lindeps: |ovlp| = {:.6e}'.format (ovlp_det) lib.logger.info (las, 'Caught error %s, %s', str (err), lc) if ovlp_det < LINDEP_THRESH: - raw2orth, orth2raw = citools.get_orth_basis (ci_blk, las.ncas_sub, nelec_blk) + raw2orth, orth2raw = citools.get_orth_basis (ci_blk, las.ncas_sub, nelec_blk, + _get_ovlp=get_ovlp) xhx = raw2orth (ham_blk).conj ().T #x = canonical_orth_(ovlp_blk, thr=LINDEP_THRESH) lib.logger.info (las, '%d/%d linearly independent model states', diff --git a/my_pyscf/lassi/op_o0.py b/my_pyscf/lassi/op_o0.py index c3ae5f51..dc8fe7c6 100644 --- a/my_pyscf/lassi/op_o0.py +++ b/my_pyscf/lassi/op_o0.py @@ -1,5 +1,6 @@ import sys import numpy as np +import functools from scipy import linalg from pyscf.fci import cistring from pyscf import fci, lib @@ -585,7 +586,8 @@ def contract_h (c, nel): ket = None ham_eff[i,:] = dotter (hket, nelec_ket, spinless2ss=spinless2ss, iket=i, oporder=2) - return ham_eff, s2_eff, ovlp_eff + _get_ovlp = functools.partial (get_ovlp, ci_fr, norb_f, nelec_frs) + return ham_eff, s2_eff, ovlp_eff, _get_ovlp def contract_ham_ci (las, h1, h2, ci_fr_ket, nelec_frs_ket, ci_fr_bra, nelec_frs_bra, soc=0, orbsym=None, wfnsym=None): diff --git a/my_pyscf/lassi/op_o1/hams2ovlp.py b/my_pyscf/lassi/op_o1/hams2ovlp.py index b3547219..044606f3 100644 --- a/my_pyscf/lassi/op_o1/hams2ovlp.py +++ b/my_pyscf/lassi/op_o1/hams2ovlp.py @@ -472,6 +472,6 @@ def ham (las, h1, h2, ci, nelec_frs, soc=0, nlas=None, _HamS2Ovlp_class=HamS2Ovl if las.verbose >= lib.logger.TIMER_LEVEL: lib.logger.info (las, 'LASSI Hamiltonian crunching profile:\n%s', outerprod.sprint_profile ()) - return ham, s2, ovlp + return ham, s2, ovlp, outerprod.get_ovlp diff --git a/my_pyscf/lassi/op_o1/stdm.py b/my_pyscf/lassi/op_o1/stdm.py index 633cbd4a..e271a278 100644 --- a/my_pyscf/lassi/op_o1/stdm.py +++ b/my_pyscf/lassi/op_o1/stdm.py @@ -425,19 +425,19 @@ def crunch_ovlp (self, bra, ket): def get_ovlp (self, rootidx=None): exc_null = self.exc_null - offs_lroots = self.offs_lroots + offs_lroots = self.offs_lroots.copy () nstates = self.nstates if rootidx is not None: rootidx = np.atleast_1d (rootidx) bra_null = np.isin (self.exc_null[:,0], rootidx) ket_null = np.isin (self.exc_null[:,1], rootidx) - idx_null = bra_null & ket_null - exc_null = exc_null[idx_null,:] - lroots = self.lroots[:,idx_null] + exc_null = exc_null[bra_null&ket_null,:] + lroots = self.lroots[:,rootidx] nprods = np.prod (lroots, axis=0) offs1 = np.cumsum (nprods) - offs0 = offs1 - nprods - offs_lroots = np.stack ([offs0, offs1], axis=1) + offs0 = offs1 - nprods + for i, iroot in enumerate (rootidx): + offs_lroots[iroot,:] = [offs0[i], offs1[i]] nstates = offs1[-1] ovlp = np.zeros ([nstates,]*2, dtype=self.dtype) for bra, ket in exc_null: diff --git a/tests/lassi/test_4frag.py b/tests/lassi/test_4frag.py index f9931f25..a76fa529 100644 --- a/tests/lassi/test_4frag.py +++ b/tests/lassi/test_4frag.py @@ -124,9 +124,9 @@ def test_stdm12s (self): def test_ham_s2_ovlp (self): h1, h2 = ham_2q (las, las.mo_coeff, veff_c=None, h2eff_sub=None)[1:] lbls = ('ham','s2','ovlp') - mats_o0 = op_o0.ham (las, h1, h2, las.ci, nelec_frs)#, orbsym=orbsym, wfnsym=wfnsym) + mats_o0 = op_o0.ham (las, h1, h2, las.ci, nelec_frs)[:3]#, orbsym=orbsym, wfnsym=wfnsym) fps_o0 = [lib.fp (mat) for mat in mats_o0] - mats_o1 = op_o1.ham (las, h1, h2, las.ci, nelec_frs)#, orbsym=orbsym, wfnsym=wfnsym) + mats_o1 = op_o1.ham (las, h1, h2, las.ci, nelec_frs)[:3]#, orbsym=orbsym, wfnsym=wfnsym) for lbl, mat, fp in zip (lbls, mats_o1, fps_o0): with self.subTest(opt=1, matrix=lbl): self.assertAlmostEqual (lib.fp (mat), fp, 9) diff --git a/tests/lassi/test_opt57_slow.py b/tests/lassi/test_opt57_slow.py index 46e9a6ba..0dd0da59 100644 --- a/tests/lassi/test_opt57_slow.py +++ b/tests/lassi/test_opt57_slow.py @@ -158,9 +158,9 @@ def test_ham_s2_ovlp (self): h1, h2 = ham_2q (las, las.mo_coeff, veff_c=None, h2eff_sub=None)[1:] lbls = ('ham','s2','ovlp') t0, w0 = lib.logger.process_clock (), lib.logger.perf_counter () - mats_o0 = op_o0.ham (las, h1, h2, las.ci, nelec_frs, orbsym=orbsym, wfnsym=wfnsym) + mats_o0 = op_o0.ham (las, h1, h2, las.ci, nelec_frs, orbsym=orbsym, wfnsym=wfnsym)[:3] t1, w1 = lib.logger.process_clock (), lib.logger.perf_counter () - mats_o1 = op_o1.ham (las, h1, h2, las.ci, nelec_frs, orbsym=orbsym, wfnsym=wfnsym) + mats_o1 = op_o1.ham (las, h1, h2, las.ci, nelec_frs, orbsym=orbsym, wfnsym=wfnsym)[:3] t2, w2 = lib.logger.process_clock (), lib.logger.perf_counter () #print (t1-t0, t2-t1) #print (w1-w0, w2-w1)