Skip to content

Commit b556049

Browse files
sunqmsunqm
authored
Improve docs and examples (pyscf#2359)
* Improve the example of UHF Hamiltonian customization. (pyscf#2350) * Add comments for tdscf modules (for issue pyscf#2246) * Improve hf.dispersion imports * Fix tdrhf
1 parent 1de6adc commit b556049

File tree

10 files changed

+97
-79
lines changed

10 files changed

+97
-79
lines changed

‎examples/scf/40-customizing_hamiltonian.py

Lines changed: 14 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -50,3 +50,17 @@
5050
# (particularly in the MO integral transformation) may ignore the customized
5151
# Hamiltonian if memory is not enough.
5252
mol.incore_anyway = True
53+
54+
#
55+
# In the case of spin-polarized system, the get_hcore method of UHF class
56+
# requires a tuple for the spin-up and spin-down one-electron Hamiltonian.
57+
# To define the number of spin-up and spin-down electrons, setting the mol.spin
58+
# attribute and assigning a tuple to mf.nelec both work.
59+
#
60+
mf = scf.UHF(mol)
61+
mf.get_hcore = lambda *args: (h1, h1)
62+
mf.get_ovlp = lambda *args: numpy.eye(n)
63+
mf._eri = ao2mo.restore(8, eri, n)
64+
# Setting mf.nelec below has the same effects to setting mol.spin = 2
65+
mf.nelec = (6, 4)
66+
mf.kernel()

‎pyscf/pbc/tdscf/krhf.py

Lines changed: 8 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -263,18 +263,24 @@ def vind(xys):
263263
dmov[i,k]+= reduce(numpy.dot, (orbv[k], dmy.T, orbo[k].T.conj()))
264264

265265
with lib.temporary_env(mf, exxdiv=None):
266-
v1ao = vresp(dmov, kshift)
266+
v1ao = vresp(dmov, kshift) # = <mb||nj> Xjb + <mj||nb> Yjb
267267
v1s = []
268268
for i in range(nz):
269269
dmx = z1xs[i]
270270
dmy = z1ys[i]
271271
v1xs = []
272272
v1ys = []
273273
for k in range(nkpts):
274+
# AX + BY
275+
# = <ib||aj> Xjb + <ij||ab> Yjb
276+
# = (<mb||nj> Xjb + <mj||nb> Yjb) Cmi* Cna
274277
v1x = reduce(numpy.dot, (orbo[k].T.conj(), v1ao[i,k], orbv[k]))
278+
# (B*)X + (A*)Y
279+
# = <ab||ij> Xjb + <aj||ib> Yjb
280+
# = (<mb||nj> Xjb + <mj||nb> Yjb) Cma* Cni
275281
v1y = reduce(numpy.dot, (orbv[k].T.conj(), v1ao[i,k], orbo[k])).T
276282
v1x+= e_ia[k] * dmx[k]
277-
v1y+= e_ia[k] * dmy[k]
283+
v1y+= e_ia[k].conj() * dmy[k]
278284
v1xs.append(v1x.ravel())
279285
v1ys.append(-v1y.ravel())
280286
# v1s += v1xs + v1ys

‎pyscf/pbc/tdscf/kuhf.py

Lines changed: 4 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -205,10 +205,10 @@ def vind(xys):
205205
for k in range(nkpts):
206206
dmx = reduce(numpy.dot, (orboa[k], xa[k] , orbva[k].conj().T))
207207
dmy = reduce(numpy.dot, (orbva[k], ya[k].T, orboa[k].conj().T))
208-
dmov[0,i,k] = dmx + dmy # AX + BY
208+
dmov[0,i,k] = dmx + dmy
209209
dmx = reduce(numpy.dot, (orbob[k], xb[k] , orbvb[k].conj().T))
210210
dmy = reduce(numpy.dot, (orbvb[k], yb[k].T, orbob[k].conj().T))
211-
dmov[1,i,k] = dmx + dmy # AX + BY
211+
dmov[1,i,k] = dmx + dmy
212212

213213
with lib.temporary_env(mf, exxdiv=None):
214214
dmov = dmov.reshape(2,nz,nkpts,nao,nao)
@@ -230,8 +230,8 @@ def vind(xys):
230230
v1yb = reduce(numpy.dot, (orbvb[k].conj().T, v1ao[1,i,k], orbob[k])).T
231231
v1xa+= e_ia_a[k] * xa[k]
232232
v1xb+= e_ia_b[k] * xb[k]
233-
v1ya+= e_ia_a[k] * ya[k]
234-
v1yb+= e_ia_b[k] * yb[k]
233+
v1ya+= e_ia_a[k].conj() * ya[k]
234+
v1yb+= e_ia_b[k].conj() * yb[k]
235235
v1xsa.append(v1xa.ravel())
236236
v1xsb.append(v1xb.ravel())
237237
v1ysa.append(-v1ya.ravel())

‎pyscf/pbc/tdscf/test/test_rhf.py

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -46,6 +46,7 @@ def setUpClass(cls):
4646

4747
cls.nstates = 5 # make sure first `nstates_test` states are converged
4848
cls.nstates_test = 2
49+
4950
@classmethod
5051
def tearDownClass(cls):
5152
cls.cell.stdout.close()
@@ -98,6 +99,7 @@ def setUpClass(cls):
9899

99100
cls.nstates = 5 # make sure first `nstates_test` states are converged
100101
cls.nstates_test = 2
102+
101103
@classmethod
102104
def tearDownClass(cls):
103105
cls.cell.stdout.close()

‎pyscf/scf/dispersion.py

Lines changed: 0 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -176,7 +176,3 @@ def get_dispersion(mf, disp=None, with_3body=None, verbose=None):
176176
return e_d4
177177
else:
178178
raise RuntimeError(f'dispersion correction: {disp_version} is not supported.')
179-
180-
# Inject to SCF class
181-
scf.hf.SCF.do_disp = check_disp
182-
scf.hf.SCF.get_dispersion = get_dispersion

‎pyscf/scf/hf.py

Lines changed: 4 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -33,6 +33,7 @@
3333
from pyscf.scf import diis
3434
from pyscf.scf import _vhf
3535
from pyscf.scf import chkfile
36+
from pyscf.scf import dispersion
3637
from pyscf.data import nist
3738
from pyscf import __config__
3839

@@ -1875,6 +1876,9 @@ def get_init_guess(self, mol=None, key='minao', **kwargs):
18751876
energy_elec = energy_elec
18761877
energy_tot = energy_tot
18771878

1879+
do_disp = dispersion.check_disp
1880+
get_dispersion = dispersion.get_dispersion
1881+
18781882
def energy_nuc(self):
18791883
return self.mol.enuc
18801884

‎pyscf/tdscf/dhf.py

Lines changed: 20 additions & 17 deletions
Original file line numberDiff line numberDiff line change
@@ -50,15 +50,13 @@ def gen_tda_operation(mf, fock_ao=None):
5050
orbo = mo_coeff[:,occidx]
5151

5252
if fock_ao is None:
53-
#dm0 = mf.make_rdm1(mo_coeff, mo_occ)
54-
#fock_ao = mf.get_hcore() + mf.get_veff(mol, dm0)
55-
foo = numpy.diag(mo_energy[occidx])
56-
fvv = numpy.diag(mo_energy[viridx])
53+
e_ia = hdiag = mo_energy[viridx] - mo_energy[occidx,None]
5754
else:
5855
fock = reduce(numpy.dot, (mo_coeff.conj().T, fock_ao, mo_coeff))
5956
foo = fock[occidx[:,None],occidx]
6057
fvv = fock[viridx[:,None],viridx]
61-
hdiag = (fvv.diagonal() - foo.diagonal()[:,None]).ravel()
58+
hdiag = fvv.diagonal() - foo.diagonal()[:,None]
59+
hdiag = hdiag.ravel()
6260

6361
mo_coeff = numpy.asarray(numpy.hstack((orbo,orbv)), order='F')
6462
vresp = mf.gen_response(hermi=0)
@@ -68,8 +66,11 @@ def vind(zs):
6866
dmov = lib.einsum('xov,qv,po->xpq', zs, orbv.conj(), orbo)
6967
v1ao = vresp(dmov)
7068
v1ov = lib.einsum('xpq,po,qv->xov', v1ao, orbo.conj(), orbv)
71-
v1ov += lib.einsum('xqs,sp->xqp', zs, fvv)
72-
v1ov -= lib.einsum('xpr,sp->xsr', zs, foo)
69+
if fock_ao is None:
70+
v1ov += numpy.einsum('xia,ia->xia', zs, e_ia)
71+
else:
72+
v1ov += lib.einsum('xqs,sp->xqp', zs, fvv)
73+
v1ov -= lib.einsum('xpr,sp->xsr', zs, foo)
7374
return v1ov.reshape(v1ov.shape[0], -1)
7475

7576
return vind, hdiag
@@ -490,9 +491,8 @@ def gen_tdhf_operation(mf, fock_ao=None):
490491
orbv = mo_coeff[:,viridx]
491492
orbo = mo_coeff[:,occidx]
492493

493-
foo = numpy.diag(mo_energy[occidx])
494-
fvv = numpy.diag(mo_energy[viridx])
495-
hdiag = fvv.diagonal() - foo.diagonal()[:,None]
494+
assert fock_ao is None
495+
e_ia = hdiag = mo_energy[viridx] - mo_energy[occidx,None]
496496
hdiag = numpy.hstack((hdiag.ravel(), -hdiag.ravel())).real
497497

498498
mo_coeff = numpy.asarray(numpy.hstack((orbo,orbv)), order='F')
@@ -501,17 +501,20 @@ def gen_tdhf_operation(mf, fock_ao=None):
501501
def vind(xys):
502502
xys = numpy.asarray(xys).reshape(-1,2,nocc,nvir)
503503
xs, ys = xys.transpose(1,0,2,3)
504-
# dms = AX + BY
505504
dms = lib.einsum('xov,qv,po->xpq', xs, orbv.conj(), orbo)
506505
dms += lib.einsum('xov,pv,qo->xpq', ys, orbv, orbo.conj())
507-
508-
v1ao = vresp(dms)
506+
v1ao = vresp(dms) # = <mb||nj> Xjb + <mj||nb> Yjb
507+
# A ~= <ib||aj>, B = <ij||ab>
508+
# AX + BY
509+
# = <ib||aj> Xjb + <ij||ab> Yjb
510+
# = (<mb||nj> Xjb + <mj||nb> Yjb) Cmi* Cna
509511
v1ov = lib.einsum('xpq,po,qv->xov', v1ao, orbo.conj(), orbv)
512+
# (B*)X + (A*)Y
513+
# = <ab||ij> Xjb + <aj||ib> Yjb
514+
# = (<mb||nj> Xjb + <mj||nb> Yjb) Cma* Cni
510515
v1vo = lib.einsum('xpq,qo,pv->xov', v1ao, orbo, orbv.conj())
511-
v1ov += lib.einsum('xqs,sp->xqp', xs, fvv) # AX
512-
v1ov -= lib.einsum('xpr,sp->xsr', xs, foo) # AX
513-
v1vo += lib.einsum('xqs,sp->xqp', ys, fvv.conj()) # (A*)Y
514-
v1vo -= lib.einsum('xpr,sp->xsr', ys, foo.conj()) # (A*)Y
516+
v1ov += numpy.einsum('xia,ia->xia', xs, e_ia) # AX
517+
v1vo += numpy.einsum('xia,ia->xia', ys, e_ia.conj()) # (A*)Y
515518

516519
# (AX, (-A*)Y)
517520
nz = xys.shape[0]

‎pyscf/tdscf/ghf.py

Lines changed: 19 additions & 22 deletions
Original file line numberDiff line numberDiff line change
@@ -67,16 +67,13 @@ def gen_tda_operation(mf, fock_ao=None, wfnsym=None):
6767
sym_forbid = (orbsym_in_d2h[occidx,None] ^ orbsym_in_d2h[viridx]) != wfnsym
6868

6969
if fock_ao is None:
70-
#dm0 = mf.make_rdm1(mo_coeff, mo_occ)
71-
#fock_ao = mf.get_hcore() + mf.get_veff(mol, dm0)
72-
foo = numpy.diag(mo_energy[occidx])
73-
fvv = numpy.diag(mo_energy[viridx])
70+
e_ia = hdiag = mo_energy[viridx] - mo_energy[occidx,None]
7471
else:
7572
fock = reduce(numpy.dot, (mo_coeff.conj().T, fock_ao, mo_coeff))
7673
foo = fock[occidx[:,None],occidx]
7774
fvv = fock[viridx[:,None],viridx]
75+
hdiag = fvv.diagonal() - foo.diagonal()[:,None]
7876

79-
hdiag = fvv.diagonal() - foo.diagonal()[:,None]
8077
if wfnsym is not None and mol.symmetry:
8178
hdiag[sym_forbid] = 0
8279
hdiag = hdiag.ravel().real
@@ -93,8 +90,11 @@ def vind(zs):
9390
dmov = lib.einsum('xov,qv,po->xpq', zs, orbv.conj(), orbo)
9491
v1ao = vresp(dmov)
9592
v1ov = lib.einsum('xpq,po,qv->xov', v1ao, orbo.conj(), orbv)
96-
v1ov += lib.einsum('xqs,sp->xqp', zs, fvv)
97-
v1ov -= lib.einsum('xpr,sp->xsr', zs, foo)
93+
if fock_ao is None:
94+
v1ov += numpy.einsum('xia,ia->xia', zs, e_ia)
95+
else:
96+
v1ov += lib.einsum('xqs,sp->xqp', zs, fvv)
97+
v1ov -= lib.einsum('xpr,sp->xsr', zs, foo)
9898
if wfnsym is not None and mol.symmetry:
9999
v1ov[:,sym_forbid] = 0
100100
return v1ov.reshape(v1ov.shape[0],-1)
@@ -483,15 +483,9 @@ def gen_tdhf_operation(mf, fock_ao=None, wfnsym=None):
483483
orbsym_in_d2h = numpy.asarray(orbsym) % 10 # convert to D2h irreps
484484
sym_forbid = (orbsym_in_d2h[occidx,None] ^ orbsym_in_d2h[viridx]) != wfnsym
485485

486-
#dm0 = mf.make_rdm1(mo_coeff, mo_occ)
487-
#fock_ao = mf.get_hcore() + mf.get_veff(mol, dm0)
488-
#fock = reduce(numpy.dot, (mo_coeff.T, fock_ao, mo_coeff))
489-
#foo = fock[occidx[:,None],occidx]
490-
#fvv = fock[viridx[:,None],viridx]
491-
foo = numpy.diag(mo_energy[occidx])
492-
fvv = numpy.diag(mo_energy[viridx])
486+
assert fock_ao is None
493487

494-
hdiag = fvv.diagonal() - foo.diagonal()[:,None]
488+
e_ia = hdiag = mo_energy[viridx] - mo_energy[occidx,None]
495489
if wfnsym is not None and mol.symmetry:
496490
hdiag[sym_forbid] = 0
497491
hdiag = numpy.hstack((hdiag.ravel(), -hdiag.ravel())).real
@@ -507,17 +501,20 @@ def vind(xys):
507501
xys[:,:,sym_forbid] = 0
508502

509503
xs, ys = xys.transpose(1,0,2,3)
510-
# dms = AX + BY
511504
dms = lib.einsum('xov,qv,po->xpq', xs, orbv.conj(), orbo)
512505
dms += lib.einsum('xov,pv,qo->xpq', ys, orbv, orbo.conj())
513-
514-
v1ao = vresp(dms)
506+
v1ao = vresp(dms) # = <mb||nj> Xjb + <mj||nb> Yjb
507+
# A ~= <ib||aj>, B = <ij||ab>
508+
# AX + BY
509+
# = <ib||aj> Xjb + <ij||ab> Yjb
510+
# = (<mb||nj> Xjb + <mj||nb> Yjb) Cmi* Cna
515511
v1ov = lib.einsum('xpq,po,qv->xov', v1ao, orbo.conj(), orbv)
512+
# (B*)X + (A*)Y
513+
# = <ab||ij> Xjb + <aj||ib> Yjb
514+
# = (<mb||nj> Xjb + <mj||nb> Yjb) Cma* Cni
516515
v1vo = lib.einsum('xpq,qo,pv->xov', v1ao, orbo, orbv.conj())
517-
v1ov += lib.einsum('xqs,sp->xqp', xs, fvv) # AX
518-
v1ov -= lib.einsum('xpr,sp->xsr', xs, foo) # AX
519-
v1vo += lib.einsum('xqs,sp->xqp', ys, fvv.conj()) # (A*)Y
520-
v1vo -= lib.einsum('xpr,sp->xsr', ys, foo.conj()) # (A*)Y
516+
v1ov += numpy.einsum('xia,ia->xia', xs, e_ia) # AX
517+
v1vo += numpy.einsum('xia,ia->xia', ys, e_ia.conj()) # (A*)Y
521518

522519
if wfnsym is not None and mol.symmetry:
523520
v1ov[:,sym_forbid] = 0

‎pyscf/tdscf/rhf.py

Lines changed: 23 additions & 24 deletions
Original file line numberDiff line numberDiff line change
@@ -70,16 +70,13 @@ def gen_tda_operation(mf, fock_ao=None, singlet=True, wfnsym=None):
7070
sym_forbid = (orbsym_in_d2h[occidx,None] ^ orbsym_in_d2h[viridx]) != wfnsym
7171

7272
if fock_ao is None:
73-
#dm0 = mf.make_rdm1(mo_coeff, mo_occ)
74-
#fock_ao = mf.get_hcore() + mf.get_veff(mol, dm0)
75-
foo = numpy.diag(mo_energy[occidx])
76-
fvv = numpy.diag(mo_energy[viridx])
73+
e_ia = hdiag = mo_energy[viridx] - mo_energy[occidx,None]
7774
else:
7875
fock = reduce(numpy.dot, (mo_coeff.conj().T, fock_ao, mo_coeff))
7976
foo = fock[occidx[:,None],occidx]
8077
fvv = fock[viridx[:,None],viridx]
78+
hdiag = fvv.diagonal() - foo.diagonal()[:,None]
8179

82-
hdiag = fvv.diagonal() - foo.diagonal()[:,None]
8380
if wfnsym is not None and mol.symmetry:
8481
hdiag[sym_forbid] = 0
8582
hdiag = hdiag.ravel()
@@ -97,8 +94,11 @@ def vind(zs):
9794
dmov = lib.einsum('xov,qv,po->xpq', zs*2, orbv.conj(), orbo)
9895
v1ao = vresp(dmov)
9996
v1ov = lib.einsum('xpq,po,qv->xov', v1ao, orbo.conj(), orbv)
100-
v1ov += lib.einsum('xqs,sp->xqp', zs, fvv)
101-
v1ov -= lib.einsum('xpr,sp->xsr', zs, foo)
97+
if fock_ao is None:
98+
v1ov += numpy.einsum('xia,ia->xia', zs, e_ia)
99+
else:
100+
v1ov += lib.einsum('xqs,sp->xqp', zs, fvv)
101+
v1ov -= lib.einsum('xpr,sp->xsr', zs, foo)
102102
if wfnsym is not None and mol.symmetry:
103103
v1ov[:,sym_forbid] = 0
104104
return v1ov.reshape(v1ov.shape[0],-1)
@@ -111,6 +111,8 @@ def get_ab(mf, mo_energy=None, mo_coeff=None, mo_occ=None):
111111
112112
A[i,a,j,b] = \delta_{ab}\delta_{ij}(E_a - E_i) + (ia||bj)
113113
B[i,a,j,b] = (ia||jb)
114+
115+
Ref: Chem Phys Lett, 256, 454
114116
'''
115117
if mo_energy is None: mo_energy = mf.mo_energy
116118
if mo_coeff is None: mo_coeff = mf.mo_coeff
@@ -884,8 +886,8 @@ def pickeig(w, v, nroots, envs):
884886
def gen_tdhf_operation(mf, fock_ao=None, singlet=True, wfnsym=None):
885887
'''Generate function to compute
886888
887-
[ A B][X]
888-
[-B -A][Y]
889+
[ A B ][X]
890+
[-B* -A*][Y]
889891
'''
890892
mol = mf.mol
891893
mo_coeff = mf.mo_coeff
@@ -908,15 +910,9 @@ def gen_tdhf_operation(mf, fock_ao=None, singlet=True, wfnsym=None):
908910
orbsym_in_d2h = numpy.asarray(orbsym) % 10 # convert to D2h irreps
909911
sym_forbid = (orbsym_in_d2h[occidx,None] ^ orbsym_in_d2h[viridx]) != wfnsym
910912

911-
#dm0 = mf.make_rdm1(mo_coeff, mo_occ)
912-
#fock_ao = mf.get_hcore() + mf.get_veff(mol, dm0)
913-
#fock = reduce(numpy.dot, (mo_coeff.T, fock_ao, mo_coeff))
914-
#foo = fock[occidx[:,None],occidx]
915-
#fvv = fock[viridx[:,None],viridx]
916-
foo = numpy.diag(mo_energy[occidx])
917-
fvv = numpy.diag(mo_energy[viridx])
913+
assert fock_ao is None
918914

919-
hdiag = fvv.diagonal() - foo.diagonal()[:,None]
915+
e_ia = hdiag = mo_energy[viridx] - mo_energy[occidx,None]
920916
if wfnsym is not None and mol.symmetry:
921917
hdiag[sym_forbid] = 0
922918
hdiag = numpy.hstack((hdiag.ravel(), -hdiag.ravel()))
@@ -932,18 +928,21 @@ def vind(xys):
932928
xys[:,:,sym_forbid] = 0
933929

934930
xs, ys = xys.transpose(1,0,2,3)
935-
# dms = AX + BY
936931
# *2 for double occupancy
937932
dms = lib.einsum('xov,qv,po->xpq', xs*2, orbv.conj(), orbo)
938933
dms += lib.einsum('xov,pv,qo->xpq', ys*2, orbv, orbo.conj())
939-
940-
v1ao = vresp(dms)
934+
v1ao = vresp(dms) # = <mb||nj> Xjb + <mj||nb> Yjb
935+
# A ~= <ib||aj>, B = <ij||ab>
936+
# AX + BY
937+
# = <ib||aj> Xjb + <ij||ab> Yjb
938+
# = (<mb||nj> Xjb + <mj||nb> Yjb) Cmi* Cna
941939
v1ov = lib.einsum('xpq,po,qv->xov', v1ao, orbo.conj(), orbv)
940+
# (B*)X + (A*)Y
941+
# = <ab||ij> Xjb + <aj||ib> Yjb
942+
# = (<mb||nj> Xjb + <mj||nb> Yjb) Cma* Cni
942943
v1vo = lib.einsum('xpq,qo,pv->xov', v1ao, orbo, orbv.conj())
943-
v1ov += lib.einsum('xqs,sp->xqp', xs, fvv) # AX
944-
v1ov -= lib.einsum('xpr,sp->xsr', xs, foo) # AX
945-
v1vo += lib.einsum('xqs,sp->xqp', ys, fvv) # AY
946-
v1vo -= lib.einsum('xpr,sp->xsr', ys, foo) # AY
944+
v1ov += numpy.einsum('xia,ia->xia', xs, e_ia) # AX
945+
v1vo += numpy.einsum('xia,ia->xia', ys, e_ia.conj()) # (A*)Y
947946

948947
if wfnsym is not None and mol.symmetry:
949948
v1ov[:,sym_forbid] = 0

‎pyscf/tdscf/uhf.py

Lines changed: 3 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -711,8 +711,8 @@ def pickeig(w, v, nroots, envs):
711711
def gen_tdhf_operation(mf, fock_ao=None, singlet=True, wfnsym=None):
712712
'''Generate function to compute
713713
714-
[ A B][X]
715-
[-B -A][Y]
714+
[ A B ][X]
715+
[-B* -A*][Y]
716716
'''
717717
mol = mf.mol
718718
mo_coeff = mf.mo_coeff
@@ -768,17 +768,14 @@ def vind(xys):
768768
xb = xs[:,nocca*nvira:].reshape(nz,noccb,nvirb)
769769
ya = ys[:,:nocca*nvira].reshape(nz,nocca,nvira)
770770
yb = ys[:,nocca*nvira:].reshape(nz,noccb,nvirb)
771-
# dms = AX + BY
772771
dmsa = lib.einsum('xov,qv,po->xpq', xa, orbva.conj(), orboa)
773772
dmsb = lib.einsum('xov,qv,po->xpq', xb, orbvb.conj(), orbob)
774773
dmsa += lib.einsum('xov,pv,qo->xpq', ya, orbva, orboa.conj())
775774
dmsb += lib.einsum('xov,pv,qo->xpq', yb, orbvb, orbob.conj())
776-
777775
v1ao = vresp(numpy.asarray((dmsa,dmsb)))
778-
779776
v1aov = lib.einsum('xpq,po,qv->xov', v1ao[0], orboa.conj(), orbva)
780-
v1avo = lib.einsum('xpq,qo,pv->xov', v1ao[0], orboa, orbva.conj())
781777
v1bov = lib.einsum('xpq,po,qv->xov', v1ao[1], orbob.conj(), orbvb)
778+
v1avo = lib.einsum('xpq,qo,pv->xov', v1ao[0], orboa, orbva.conj())
782779
v1bvo = lib.einsum('xpq,qo,pv->xov', v1ao[1], orbob, orbvb.conj())
783780

784781
v1ov = xs * e_ia # AX

0 commit comments

Comments
 (0)