flake 8 happy?

adc-connect · Dec 14, 2021 · f3cd5f2 · f3cd5f2
1 parent d0f9950
commit f3cd5f2
Show file tree

Hide file tree

Showing 8 changed files with 109 additions and 167 deletions.
diff --git a/adcc/backends/pyscf.py b/adcc/backends/pyscf.py
@@ -48,11 +48,6 @@ def correlated_gradient(self, g1_ao, w_ao, g2_ao_1, g2_ao_2):
         Gradient["TEI"] = np.zeros((natoms, 3))
         Gradient["Total"] = np.zeros((natoms, 3))
 
-        # Iulia added variables; TODO: remove
-        #ovlp_iulia = np.zeros((natoms, 3))
-        #hcore_iulia = np.zeros((natoms, 3))
-        #eri_iulia = np.zeros((natoms, 3))
-
         # TODO: does RHF/UHF matter here?
         gradient = grad.RHF(self.scfres)
         hcore_deriv = gradient.hcore_generator()
@@ -71,13 +66,10 @@ def correlated_gradient(self, g1_ao, w_ao, g2_ao_1, g2_ao_2):
             Sx_a[:, k0:k1] = Sx[:, k0:k1]
             Sx_a += Sx_a.transpose(0, 2, 1)
             Gradient["S"][ia] += np.einsum("xpq,pq->x", Sx_a, w_ao)
-            #ovlp_iulia[ia] += np.einsum("xpq,pq->x", Sx_a, w_ao)
-
 
             # derivative of the core Hamiltonian
             Hx_a = hcore_deriv(ia)
             Gradient["T+V"][ia] += np.einsum("xpq,pq->x", Hx_a, g1_ao)
-            #hcore_iulia[ia] += np.einsum("xpq,pq->x", Hx_a, g1_ao)
 
             # derivatives of the ERIs
             ERIx_a = np.zeros_like(ERIx)
@@ -93,23 +85,6 @@ def correlated_gradient(self, g1_ao, w_ao, g2_ao_1, g2_ao_2):
             Gradient["TEI"][ia] -= np.einsum(
                 "pqrs,xpsqr->x", g2_ao_2, ERIx_a, optimize=True
             )
-            #eri_iulia[ia] += (
-            #      np.einsum("pqrs,xprqs->x", g2_ao_1, ERIx_a, optimize=True)
-            #    - np.einsum("pqrs,xpsqr->x", g2_ao_2, ERIx_a, optimize=True)
-            #)
-
-
-        # Iulia print TODO:remove
-        #print("2PDM_1 in AO:\n", g2_ao_1)
-        #print("Hcore contribution to the gradient:")
-        #print(hcore_iulia)
-        #print()
-        #print("Overlp contribution to the gradient:")
-        #print(ovlp_iulia)
-        #print()
-        #print("ERI contribution to the gradient:")
-        #print(eri_iulia)
-        #print()
 
         Gradient["N"] = gradient.grad_nuc()
         Gradient["OEI"] = Gradient["T+V"] + Gradient["S"]

diff --git a/adcc/gradients/TwoParticleDensityMatrix.py b/adcc/gradients/TwoParticleDensityMatrix.py
@@ -62,7 +62,7 @@ def __init__(self, spaces):
                 b.cccc, b.ococ, b.oooo,  b.cvcv,
                 b.ocov, b.cccv, b.cocv, b.ocoo,
                 b.ccco, b.occv, b.ccvv, b.ocvv,
-				b.vvvv,
+                b.vvvv,
             ]
         self._tensors = {}
 

diff --git a/adcc/gradients/__init__.py b/adcc/gradients/__init__.py
@@ -63,13 +63,13 @@ def nuclear_gradient(excitation_or_mp):
 
     with timer.record("orbital_response"):
         rhs = orbital_response_rhs(hf, g1a, g2a).evaluate()
-        l = orbital_response(hf, rhs)
+        lam = orbital_response(hf, rhs)
 
     # orbital-relaxed OPDM (without reference state)
     g1o = g1a.copy()
-    g1o.ov = 0.5 * l.ov
+    g1o.ov = 0.5 * lam.ov
     if hf.has_core_occupied_space:
-        g1o.cv = 0.5 * l.cv 
+        g1o.cv = 0.5 * lam.cv
     # orbital-relaxed OPDM (including reference state)
     g1 = g1o.copy()
     g1 += hf.density
@@ -84,48 +84,48 @@ def nuclear_gradient(excitation_or_mp):
             delta_IJ = hf.density.cc
             g2_hf = TwoParticleDensityMatrix(hf)
 
-            g2_hf.oooo = -0.25 * ( einsum("ik,jl->ijkl", delta_ij, delta_ij))
-            g2_hf.cccc = -0.5 * ( einsum("IK,JL->IJKL", delta_IJ, delta_IJ))
-            g2_hf.ococ = -1.0 * ( einsum("ik,JL->iJkL", delta_ij, delta_IJ))
+            g2_hf.oooo = -0.25 * einsum("ik,jl->ijkl", delta_ij, delta_ij)
+            g2_hf.cccc = -0.5 * einsum("IK,JL->IJKL", delta_IJ, delta_IJ)
+            g2_hf.ococ = -1.0 * einsum("ik,JL->iJkL", delta_ij, delta_IJ)
 
             g2_oresp = TwoParticleDensityMatrix(hf)
-            g2_oresp.cccc = einsum("IK,JL->IJKL", delta_IJ, g1o.cc+delta_IJ)
-            g2_oresp.ococ = ( 
-                 einsum("ik,JL->iJkL", delta_ij, g1o.cc+2*delta_IJ)
+            g2_oresp.cccc = einsum("IK,JL->IJKL", delta_IJ, g1o.cc + delta_IJ)
+            g2_oresp.ococ = (
+                + einsum("ik,JL->iJkL", delta_ij, g1o.cc + 2.0 * delta_IJ)
                 + einsum("ik,JL->iJkL", g1o.oo, delta_IJ)
             )
             g2_oresp.oooo = 0.25 * (
-                einsum("ik,jl->ijkl", delta_ij, g1o.oo+delta_ij)
+                einsum("ik,jl->ijkl", delta_ij, g1o.oo + delta_ij)
             )
             g2_oresp.ovov = einsum("ij,ab->iajb", delta_ij, g1o.vv)
             g2_oresp.cvcv = einsum("IJ,ab->IaJb", delta_IJ, g1o.vv)
-            g2_oresp.ocov = 2*einsum("ik,Ja->iJka", delta_ij, g1o.cv)
-            g2_oresp.cccv = 2*einsum("IK,Ja->IJKa", delta_IJ, g1o.cv)
-            g2_oresp.ooov = 2*einsum("ik,ja->ijka", delta_ij, g1o.ov)
-            g2_oresp.cocv = 2*einsum("IK,ja->IjKa", delta_IJ, g1o.ov)
-            g2_oresp.ocoo = 2*einsum("ik,Jl->iJkl", delta_ij, g1o.co)
-            g2_oresp.ccco = 2*einsum("IK,Jl->IJKl", delta_IJ, g1o.co)
+            g2_oresp.ocov = 2 * einsum("ik,Ja->iJka", delta_ij, g1o.cv)
+            g2_oresp.cccv = 2 * einsum("IK,Ja->IJKa", delta_IJ, g1o.cv)
+            g2_oresp.ooov = 2 * einsum("ik,ja->ijka", delta_ij, g1o.ov)
+            g2_oresp.cocv = 2 * einsum("IK,ja->IjKa", delta_IJ, g1o.ov)
+            g2_oresp.ocoo = 2 * einsum("ik,Jl->iJkl", delta_ij, g1o.co)
+            g2_oresp.ccco = 2 * einsum("IK,Jl->IJKl", delta_IJ, g1o.co)
 
             # scale for contraction with integrals
             g2a.oovv *= 0.5
             g2a.ccvv *= 0.5
             g2a.occv *= 2.0
-            g2a.vvvv *= 0.25 # Scaling twice! 
-            g2a.vvvv *= 0.25 # it's the only way it works...
+            g2a.vvvv *= 0.25  # Scaling twice!
+            g2a.vvvv *= 0.25  # it's the only way it works...
 
             g2_total = evaluate(g2_hf + g2a + g2_oresp)
         else:
             delta_ij = hf.density.oo
             g2_hf = TwoParticleDensityMatrix(hf)
             g2_hf.oooo = 0.25 * (- einsum("li,jk->ijkl", delta_ij, delta_ij)
                                  + einsum("ki,jl->ijkl", delta_ij, delta_ij))
-    
+
             g2_oresp = TwoParticleDensityMatrix(hf)
             g2_oresp.oooo = einsum("ij,kl->kilj", delta_ij, g1o.oo)
             g2_oresp.ovov = einsum("ij,ab->iajb", delta_ij, g1o.vv)
             g2_oresp.ooov = (- einsum("kj,ia->ijka", delta_ij, g1o.ov)
                              + einsum("ki,ja->ijka", delta_ij, g1o.ov))
-    
+
             # scale for contraction with integrals
             g2a.oovv *= 0.5
             g2_total = evaluate(g2_hf + g2a + g2_oresp)

diff --git a/adcc/gradients/amplitude_response.py b/adcc/gradients/amplitude_response.py
@@ -54,17 +54,14 @@ def t2bar_oovv_adc2(exci, g1a_adc0):
     )
     return t2bar
 
+
 def t2bar_oovv_cvs_adc2(exci, g1a_adc0):
     mp = exci.ground_state
     hf = mp.reference_state
-    u = exci.excitation_vector
     df_ia = mp.df(b.ov)
-    t2bar = 0.5*(
+    t2bar = 0.5 * (
         - einsum("ijcb,ac->ijab", hf.oovv, g1a_adc0.vv).antisymmetrise((2, 3))
-        ) / (
-            direct_sum("ia+jb->ijab", df_ia, df_ia).symmetrise((0, 1))
-        )
-
+    ) / direct_sum("ia+jb->ijab", df_ia, df_ia).symmetrise((0, 1))
     return t2bar
 
 
@@ -78,6 +75,7 @@ def ampl_relaxed_dms_adc1(exci):
     g2a.ovov = - 1.0 * einsum("ja,ib->iajb", u.ph, u.ph)
     return g1a, g2a
 
+
 def ampl_relaxed_dms_adc0(exci):
     hf = exci.reference_state
     u = exci.excitation_vector
@@ -89,6 +87,7 @@ def ampl_relaxed_dms_adc0(exci):
     g1a.vv = + 1.0 * einsum("ia,ib->ab", u.ph, u.ph)
     return g1a, g2a
 
+
 def ampl_relaxed_dms_cvs_adc0(exci):
     hf = exci.reference_state
     u = exci.excitation_vector
@@ -100,6 +99,7 @@ def ampl_relaxed_dms_cvs_adc0(exci):
     g1a.vv = + 1.0 * einsum("Ia,Ib->ab", u.ph, u.ph)
     return g1a, g2a
 
+
 def ampl_relaxed_dms_cvs_adc1(exci):
     hf = exci.reference_state
     u = exci.excitation_vector
@@ -118,6 +118,7 @@ def ampl_relaxed_dms_cvs_adc1(exci):
     g1a.co = - 1.0 * einsum('JbKc,ibKc->Ji', g2a.cvcv, hf.ovcv) / fco
     return g1a, g2a
 
+
 def ampl_relaxed_dms_cvs_adc2(exci):
     hf = exci.reference_state
     mp = exci.ground_state
@@ -139,10 +140,9 @@ def ampl_relaxed_dms_cvs_adc2(exci):
         - 0.5 * einsum('kIab,kJab->IJ', t2ocvv, t2ocvv)
     )
 
-    g1a.oo = ( 
+    g1a.oo = (
         - 1.0 * einsum("jKba,iKba->ij", u.pphh, u.pphh)
         - 2.0 * einsum("ikab,jkab->ij", t2bar, t2oovv).symmetrise((0, 1))
-        #- 1.0 * einsum('jkab,ikab->ij', t2oovv, t2bar)
         - 0.5 * einsum('iKab,jKab->ij', t2ocvv, t2ocvv)
         - 0.5 * einsum('ikab,jkab->ij', t2oovv, t2oovv)
     )
@@ -157,7 +157,6 @@ def ampl_relaxed_dms_cvs_adc2(exci):
         + 1.0 * einsum("Ia,Ib->ab", u.ph, u.ph)
         + 2.0 * einsum('jIcb,jIca->ab', u.pphh, u.pphh)
         + 2.0 * einsum('ijac,ijbc->ab', t2bar, t2oovv).symmetrise((0, 1))
-        #+ 1.0 * einsum('ijbc,ijac->ab', t2bar, t2oovv)
         + 0.5 * einsum('IJac,IJbc->ab', t2ccvv, t2ccvv)
         + 0.5 * einsum('ijac,ijbc->ab', t2oovv, t2oovv)
         + 1.0 * einsum('iJac,iJbc->ab', t2ocvv, t2ocvv)
@@ -169,32 +168,29 @@ def ampl_relaxed_dms_cvs_adc2(exci):
 
     # The factor 1/sqrt(2) is needed because of the scaling used in adcc
     # for the ph-pphh blocks.
-    g2a.occv  = (1 / sqrt(2)) * (
-         2.0 * einsum('Ib,kJba->kJIa', u.ph, u.pphh)
-        #- einsum('Ib,kJab->kJIa', u.ph, u.pphh)
+    g2a.occv = (1 / sqrt(2)) * (
+        2.0 * einsum('Ib,kJba->kJIa', u.ph, u.pphh)
     )
 
     g2a.oovv = (
         + 1.0 * einsum('ijcb,ca->ijab', t2oovv, g1a_cvs0.vv).antisymmetrise((2, 3))
-        #- 0.5 * einsum('ijca,cb->ijab', t2oovv, g1a_cvs0.vv)
         - 1.0 * t2oovv
         - 2.0 * t2bar
     )
 
     # The factor 2/sqrt(2) is necessary because of the way
     # that the ph-pphh is scaled.
-    g2a.ovvv = (2 / sqrt(2) ) * (
+    g2a.ovvv = (2 / sqrt(2)) * (
         einsum('Ja,iJcb->iabc', u.ph, u.pphh)
     )
-
-    g2a.ccvv = -t2ccvv 
-
-    g2a.ocvv = -t2ocvv
+
+    g2a.ccvv = - 1.0 * t2ccvv
+    g2a.ocvv = - 1.0 * t2ocvv
 
     # This is the OC block of the orbital response
     # Lagrange multipliers (lambda):
     g1a.co = (
-        - 1.0 * einsum('JbKc,ibKc->Ji', g2a.cvcv, hf.ovcv) 
+        - 1.0 * einsum('JbKc,ibKc->Ji', g2a.cvcv, hf.ovcv)
         - 0.5 * einsum('JKab,iKab->Ji', g2a.ccvv, hf.ocvv)
         + 1.0 * einsum('kJLa,ikLa->Ji', g2a.occv, hf.oocv)
         + 0.5 * einsum('kJab,ikab->Ji', g2a.ocvv, hf.oovv)
@@ -205,9 +201,9 @@ def ampl_relaxed_dms_cvs_adc2(exci):
         + 0.5 * einsum('iabc,Jabc->Ji', g2a.ovvv, hf.cvvv)
     ) / fco
 
-
     return g1a, g2a
 
+
 def ampl_relaxed_dms_cvs_adc2x(exci):
     hf = exci.reference_state
     mp = exci.ground_state
@@ -229,10 +225,9 @@ def ampl_relaxed_dms_cvs_adc2x(exci):
         - 0.5 * einsum('kIab,kJab->IJ', t2ocvv, t2ocvv)
     )
 
-    g1a.oo = ( 
+    g1a.oo = (
         - 1.0 * einsum("jKba,iKba->ij", u.pphh, u.pphh)
         - 2.0 * einsum("ikab,jkab->ij", t2bar, t2oovv).symmetrise((0, 1))
-        #- 1.0 * einsum('jkab,ikab->ij', t2oovv, t2bar) # use symmetrize?
         - 0.5 * einsum('iKab,jKab->ij', t2ocvv, t2ocvv)
         - 0.5 * einsum('ikab,jkab->ij', t2oovv, t2oovv)
     )
@@ -247,7 +242,6 @@ def ampl_relaxed_dms_cvs_adc2x(exci):
         + 1.0 * einsum("Ia,Ib->ab", u.ph, u.ph)
         + 2.0 * einsum('jIcb,jIca->ab', u.pphh, u.pphh)
         + 2.0 * einsum('ijac,ijbc->ab', t2bar, t2oovv).symmetrise((0, 1))
-        #+ 1.0 * einsum('ijbc,ijac->ab', t2bar, t2oovv)
         + 0.5 * einsum('IJac,IJbc->ab', t2ccvv, t2ccvv)
         + 0.5 * einsum('ijac,ijbc->ab', t2oovv, t2oovv)
         + 1.0 * einsum('iJac,iJbc->ab', t2ocvv, t2ocvv)
@@ -256,49 +250,39 @@ def ampl_relaxed_dms_cvs_adc2x(exci):
     g2a.cvcv = (
         - 1.0 * einsum("Ja,Ib->IaJb", u.ph, u.ph)
         - 1.0 * einsum('kIbc,kJac->IaJb', u.pphh, u.pphh)
-        #- 0.5 * einsum('kIcb,kJca->IaJb', u.pphh, u.pphh)
         + 1.0 * einsum('kIcb,kJac->IaJb', u.pphh, u.pphh)
-        #+ 0.5 * einsum('kIbc,kJca->IaJb', u.pphh, u.pphh)
     )
 
     # The factor 1/sqrt(2) is needed because of the scaling used in adcc
     # for the ph-pphh blocks.
-    g2a.occv  = (1 / sqrt(2)) * (
-          2.0 * einsum('Ib,kJba->kJIa', u.ph, u.pphh)
-        #- einsum('Ib,kJab->kJIa', u.ph, u.pphh)
+    g2a.occv = (1 / sqrt(2)) * (
+        2.0 * einsum('Ib,kJba->kJIa', u.ph, u.pphh)
     )
 
     g2a.oovv = (
         + 1.0 * einsum('ijcb,ca->ijab', t2oovv, g1a_cvs0.vv).antisymmetrise((2, 3))
-        #- 0.5 * einsum('ijca,cb->ijab', t2oovv, g1a_cvs0.vv)
         - 1.0 * t2oovv
         - 2.0 * t2bar
     )
 
-    # The factor 2/sqrt(2) is necessary because of 
+    # The factor 2/sqrt(2) is necessary because of
     # the way that the ph-pphh is scaled
     g2a.ovvv = (2 / sqrt(2)) * (
         einsum('Ja,iJcb->iabc', u.ph, u.pphh)
     )
 
     g2a.ovov = 1.0 * (
         - einsum("iKbc,jKac->iajb", u.pphh, u.pphh)
-        #- einsum("iKcb,jKca->iajb", u.pphh, u.pphh)
         + einsum("iKcb,jKac->iajb", u.pphh, u.pphh)
-        #+ einsum("iKbc,jKca->iajb", u.pphh, u.pphh)
     )
-
-    g2a.ccvv = -t2ccvv 
-
-    g2a.ocvv = -t2ocvv
 
+    g2a.ccvv = - 1.0 * t2ccvv
+    g2a.ocvv = - 1.0 * t2ocvv
     g2a.ococ = 1.0 * einsum("iJab,kLab->iJkL", u.pphh, u.pphh)
-
     g2a.vvvv = 2.0 * einsum("iJcd,iJab->abcd", u.pphh, u.pphh)
 
-    # These are the OC multipliers:
     g1a.co = (
-        - 1.0 * einsum('JbKc,ibKc->Ji', g2a.cvcv, hf.ovcv) 
+        - 1.0 * einsum('JbKc,ibKc->Ji', g2a.cvcv, hf.ovcv)
         - 0.5 * einsum('JKab,iKab->Ji', g2a.ccvv, hf.ocvv)
         + 1.0 * einsum('kJLa,ikLa->Ji', g2a.occv, hf.oocv)
         + 0.5 * einsum('kJab,ikab->Ji', g2a.ocvv, hf.oovv)
@@ -307,7 +291,7 @@ def ampl_relaxed_dms_cvs_adc2x(exci):
         + 0.5 * einsum('iKab,JKab->Ji', g2a.ocvv, hf.ccvv)
         - 0.5 * einsum('ikab,kJab->Ji', g2a.oovv, hf.ocvv)
         + 0.5 * einsum('iabc,Jabc->Ji', g2a.ovvv, hf.cvvv)
-        + 1.0 * einsum('kJmL,ikmL->Ji', g2a.ococ, hf.oooc) 
+        + 1.0 * einsum('kJmL,ikmL->Ji', g2a.ococ, hf.oooc)
         - 1.0 * einsum('iKlM,JKMl->Ji', g2a.ococ, hf.ccco)
         + 1.0 * einsum('iakb,kbJa->Ji', g2a.ovov, hf.ovcv)
     ) / fco
@@ -369,7 +353,7 @@ def ampl_relaxed_dms_mp2(mp):
     "cvs-adc0": ampl_relaxed_dms_cvs_adc0,
     "cvs-adc1": ampl_relaxed_dms_cvs_adc1,
     "cvs-adc2": ampl_relaxed_dms_cvs_adc2,
-    "cvs-adc2x": ampl_relaxed_dms_cvs_adc2x, 
+    "cvs-adc2x": ampl_relaxed_dms_cvs_adc2x,
 }