lassi hci op_o1 3frag syntax debug

my_pyscf/lassi/op_o1/hci.py

@@ -41,9 +41,7 @@ def make_exc_tables (self, hopping_index):
         exc['1s1c'] = np.append (exc_1s1c, exc_1s1cT, axis=0)
         return exc
 
-    def _fn_row_has_spin (self, _crunch_fn):
-        return _crunch_fn.__name__ in ('_crunch_1c_', '_crunch_1c1d_', '_crunch_1s1c_',
-                                       '_crunch_2c_')
+    interaction_has_spin = ('_1c_', '_1c1d_', '_1s1c_', '_2c_')
 
     def _init_vecs (self):
         hci_fr_pabq = []
@@ -159,10 +157,10 @@ def _crunch_1c1d_(self, bra, ket, i, j, k, s1):
         p_i = self.ints[i].get_1_p (bra, ket, s1)
         h_j = self.ints[j].get_1_h (bra, ket, s1)
         if iad or jad:
-            d1s_k = self.ints[k].get_1_dm1 (bra, ket)
+            d1s_kk = self.ints[k].get_1_dm1 (bra, ket)
             d1_kk = d1s_kk.sum (0)
             h1_ij  = np.tensordot (h2j, d1_kk, axes=2)
-            h1_ij -= np.tensordot (d1s_kk[s1], h2k, axes=((1,2),(2,1)))
+            h1_ij -= np.tensordot (d1s_kk[s1], h2k, axes=((0,1),(2,1)))
             if iad:
                 h_ = np.dot (h1_ij, h_j)
                 hci_f_ab[i] += fac * self.ints[i].contract_h10 (s1, h_, None, ket)
@@ -173,11 +171,11 @@ def _crunch_1c1d_(self, bra, ket, i, j, k, s1):
             h2j = np.tensordot (p_i, h2j, axes=1)
             h2j = np.tensordot (h_j, h2j, axes=1)
             h_ = h2j # opposite-spin; Coulomb effect only
-            hci_f_ab[k] += fac * self.ints[k].contract_h11 (2-(2*s1), h_, ket)
+            hci_f_ab[k] += fac * self.ints[k].contract_h11 (3*(1-s1), h_, ket)
             h2k = np.tensordot (p_i, h2k, axes=1)
             h2k = np.tensordot (h2k, h_j, axes=1)
             h_ -= h2k.T # same-spin; Coulomb and exchange
-            hci_f_ab[k] += fac * self.ints[k].contract_h11 (2*s1, h_, ket)
+            hci_f_ab[k] += fac * self.ints[k].contract_h11 (3*s1, h_, ket)
         dt, dw = logger.process_clock () - t0, logger.perf_counter () - w0
         self.dt_1c1d, self.dw_1c1d = self.dt_1c1d + dt, self.dw_1c1d + dw
         self._put_vecs_(bra, ket, hci_f_ab, i,j,k)
@@ -230,6 +228,8 @@ def _crunch_1s1c_(self, bra, ket, i, j, k, s1):
         s11 = s1
         s12 = 1-s1
         s2 = 2-s1
+        nelec_f_bra = self.nelec_rf[self.rootaddr[bra]]
+        nelec_f_ket = self.nelec_rf[self.rootaddr[ket]]
         fac = -1 # a'bb'a -> a'ab'b signi
         fac *= fermion_des_shuffle (nelec_f_bra, (i, j, k), i)
         fac *= fermion_des_shuffle (nelec_f_ket, (i, j, k), j)
@@ -247,7 +247,7 @@ def _crunch_1s1c_(self, bra, ket, i, j, k, s1):
                 hci_f_ab[j] += fac * self.ints[j].contract_h01 (s12, h_, None, ket)
         if kad:
             h_ = np.tensordot (h2, h_j, axes=1)
-            h_ = np.tensordot (h2, p_i, axes=1).T
+            h_ = np.tensordot (p_i, h_, axes=1).T
             hci_f_ab[k] += fac * self.ints[k].contract_h11 (s2, h_, ket)
         dt, dw = logger.process_clock () - t0, logger.perf_counter () - w0
         self.dt_1s1c, self.dw_1s1c = self.dt_1s1c + dt, self.dw_1s1c + dw

my_pyscf/lassi/op_o1/stdm.py

@@ -263,7 +263,7 @@ def make_exc_tables (self, hopping_index):
         # This will only be used when we are unable to restrict ourselves to the lower triangle
         idx = idx_2e & (nspin_index==3) & (ncharge_index==2) & (np.amax(spin_index,axis=0)==2)
         if nfrags > 2: exc['1s1c_T'] = np.vstack (
-            list (np.where (idx)) + [findf[-1][idx], findf[1][idx], findf[0][idx]]
+            list (np.where (idx)) + [findf[-2][idx], findf[0][idx], findf[-1][idx]]
         ).T
 
         # Symmetric two-electron interactions: lower triangle only
@@ -296,19 +296,22 @@ def make_exc_tables (self, hopping_index):
 
         return exc
 
+    interaction_has_spin = ('_1c_', '_1c1d_', '_2c_')
+
     def mask_exc_table_(self, exc, lbl, mask_bra_space=None, mask_ket_space=None):
         # Part 1: restrict to the caller-specified rectangle OR its transpose
         idx1  = mask_exc_table (exc, col=0, mask_space=mask_bra_space)
         idx1 &= mask_exc_table (exc, col=1, mask_space=mask_ket_space)
-        idx2  = mask_exc_table (exc, col=1, mask_space=mask_bra_space)
-        idx2 &= mask_exc_table (exc, col=0, mask_space=mask_ket_space)
-        exc = exc[idx1|idx2]
+        #idx2  = mask_exc_table (exc, col=1, mask_space=mask_bra_space)
+        #idx2 &= mask_exc_table (exc, col=0, mask_space=mask_ket_space)
+        exc = exc[idx1]#|idx2]
         # Part 2: identify interactions which are equivalent except for the overlap
         # factor of spectator fragments. Reduce the exc table only to the unique
         # interactions and populate self.nonuniq_exc with the corresponding
         # nonunique images.
         if lbl=='null': return exc
-        excp = exc[:,:-1] if lbl in ('1c', '1c1d', '2c') else exc
+        ulblu = '_' + lbl + '_'
+        excp = exc[:,:-1] if ulblu in self.interaction_has_spin else exc
         fprint = []
         for row in excp:
             frow = []
@@ -823,7 +826,7 @@ def _crunch_2c_(self, bra, ket, i, j, k, l, s2lt):
         self._put_D2_(bra, ket, d2, i, j, k, l)
 
     def _fn_row_has_spin (self, _crunch_fn):
-        return _crunch_fn.__name__ in ('_crunch_1c_', '_crunch_1c1d_', '_crunch_2c_')
+        return any ((i in _crunch_fn.__name__ for i in self.interaction_has_spin))
 
     def _crunch_env_(self, _crunch_fn, *row):
         if self._fn_row_has_spin (_crunch_fn):