Skip to content

Commit fb6b6a3

Browse files
authored
Add nonadiabatic coupling vector using TDDFT-ris (#451)
* add TDDFT-ris-nacv * add unit tests * add density-fitting support * fix some typos * fix some typos * debug finite difference ee * * add finite difference. * add unit tests. * add files * merge with master ris * Review the codes, according to comments. * fix some typos * change difficult-converged unit test to get_ab * add a unit test for df tddft-ris * change the basis * change the basis set * fix a bug
1 parent 25f91d3 commit fb6b6a3

34 files changed

+3046
-106
lines changed

examples/34-tdhf-nacv.py renamed to examples/34-tddft-nacv.py

Lines changed: 55 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -14,7 +14,7 @@
1414
# limitations under the License.
1515

1616
'''
17-
Nonadiabatic coupling vectors between ground and excited states for RHF
17+
Nonadiabatic coupling vectors for TDRHF and TDRKS
1818
'''
1919

2020
# This example will gives the derivative coupling (DC),
@@ -24,6 +24,7 @@
2424
import pyscf
2525
import gpu4pyscf
2626
from gpu4pyscf.scf import hf
27+
from gpu4pyscf.dft import rks
2728

2829
atom = '''
2930
O 0.0000000000 -0.0000000000 0.1174000000
@@ -39,8 +40,8 @@
3940
td = mf.TDA().set(nstates=5) # TDHF is OK
4041
td.kernel() # [ 9.21540892 10.99036172 11.83380819 13.62301694 15.06349085]
4142

42-
nac = td.NAC()
43-
nac.state=(0,1) # same as (1,0) 0 means ground state, 1 means the first excited state
43+
nac = td.nac_method()
44+
nac.states=(0,1) # same as (1,0) 0 means ground state, 1 means the first excited state
4445
nac.kernel()
4546
'''
4647
--------- TDA nonadiabatic derivative coupling for state 0 and 1----------
@@ -66,6 +67,57 @@
6667
----------------------------------------------
6768
'''
6869

70+
print('-----------------------------------------------------')
71+
print("Non-adiabatic coupling matrix element (NACME) between ground and first excited state")
72+
print(nac.de)
73+
print('-----------------------------------------------------')
74+
print("NACME between ground and first excited state scaled by E (/E_ex)")
75+
print(nac.de_scaled)
76+
print('-----------------------------------------------------')
77+
print("NACME between ground and first excited state with ETF (electron translation factor)")
78+
# Without including the contribution of the electron translation factor (ETF), for some molecules,
79+
# the non-adiabatic coupling matrix element (NACME) may lack translational invariance,
80+
# which can further lead to errors in subsequent calculations such as MD simulations.
81+
# In this case, it is necessary to use the NACME that takes the ETF into account.
82+
print(nac.de_etf)
83+
print('-----------------------------------------------------')
84+
print("NACME between ground and first excited state with ETF (electron translation factor) scaled by E (/E_ex)")
85+
print(nac.de_etf_scaled)
86+
87+
88+
mf = rks.RKS(mol, xc='b3lyp') # -76.4203783335521
89+
mf.kernel()
90+
91+
td = mf.TDA().set(nstates=5) # TDHF is OK
92+
td.kernel() # [ 7.63727447 9.47865422 10.00032863 11.95971483 14.06564139]
93+
94+
nac = td.nac_method()
95+
nac.states=(1,2) # same as (1,2) 1 means the first excited state, 2 means the second excited state
96+
nac.kernel()
97+
"""
98+
--------- TDA nonadiabatic derivative coupling for states 1 and 2----------
99+
x y z
100+
0 O -0.1134916788 -0.0000000000 0.0000000000
101+
1 H 0.0639158824 0.0000000000 0.0424664269
102+
2 H 0.0639158824 -0.0000000000 -0.0424664269
103+
--------- TDA nonadiabatic derivative coupling for states 1 and 2 after E scaled (divided by E)----------
104+
x y z
105+
0 O -1.6771477392 -0.0000000000 0.0000000000
106+
1 H 0.9445307246 0.0000000000 0.6275567747
107+
2 H 0.9445307246 -0.0000000000 -0.6275567747
108+
--------- TDA nonadiabatic derivative coupling for states 1 and 2 with ETF----------
109+
x y z
110+
0 O -0.1180169882 0.0000000000 0.0000000000
111+
1 H 0.0590085046 0.0000000000 0.0438508910
112+
2 H 0.0590085046 -0.0000000000 -0.0438508910
113+
--------- TDA nonadiabatic derivative coupling for states 1 and 2 with ETF after E scaled (divided by E)----------
114+
x y z
115+
0 O -1.7440214738 0.0000000000 0.0000000000
116+
1 H 0.8720108929 0.0000000000 0.6480159906
117+
2 H 0.8720108929 -0.0000000000 -0.6480159906
118+
----------------------------------------------
119+
"""
120+
69121
print('-----------------------------------------------------')
70122
print("Non-adiabatic coupling matrix element (NACME) between ground and first excited state")
71123
print(nac.de)

examples/37-tddft_ris_gradient.py

Lines changed: 86 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,86 @@
1+
#!/usr/bin/env python
2+
# Copyright 2021-2025 The PySCF Developers. All Rights Reserved.
3+
#
4+
# Licensed under the Apache License, Version 2.0 (the "License");
5+
# you may not use this file except in compliance with the License.
6+
# You may obtain a copy of the License at
7+
#
8+
# http://www.apache.org/licenses/LICENSE-2.0
9+
#
10+
# Unless required by applicable law or agreed to in writing, software
11+
# distributed under the License is distributed on an "AS IS" BASIS,
12+
# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
13+
# See the License for the specific language governing permissions and
14+
# limitations under the License.
15+
16+
'''
17+
Gradient for TDDFT-RIS
18+
'''
19+
20+
import pyscf
21+
import numpy as np
22+
import gpu4pyscf
23+
from gpu4pyscf.dft import rks
24+
from gpu4pyscf.tdscf import ris
25+
26+
atom = '''
27+
O 0.0000000000 -0.0000000000 0.1174000000
28+
H -0.7570000000 -0.0000000000 -0.4696000000
29+
H 0.7570000000 0.0000000000 -0.4696000000
30+
'''
31+
32+
mol = pyscf.M(atom=atom, basis='def2tzvp')
33+
34+
mf = rks.RKS(mol, xc='pbe0') # -76.3773133945678
35+
mf.kernel()
36+
37+
td = mf.TDA()
38+
td.nstates=5
39+
td.kernel() # [ 7.81949919 9.71029362 10.13398432 12.10163229 13.93675959] (eV)
40+
41+
g = td.nuc_grad_method()
42+
g.state=1
43+
g.kernel()
44+
"""
45+
--------- TDA gradients for state 1 ----------
46+
x y z
47+
0 O 0.0000000000 0.0000000000 -0.0949023769
48+
1 H 0.0627472634 -0.0000000000 0.0474538726
49+
2 H -0.0627472634 -0.0000000000 0.0474538726
50+
----------------------------------------------
51+
"""
52+
53+
td_ris = ris.TDA(mf=mf.to_gpu(), nstates=5, spectra=False, single=False, gram_schmidt=True, Ktrunc=0.0)
54+
td_ris.conv_tol = 1.0E-4
55+
td_ris.kernel() # [ 7.56352157 9.65590899 10.1236409 12.09873137 13.921372 ] (eV)
56+
57+
g_ris = td_ris.nuc_grad_method()
58+
g_ris.state=1
59+
g_ris.kernel()
60+
"""
61+
--------- TDA gradients for state 1 ----------
62+
x y z
63+
0 O 0.0000000106 -0.0000000000 -0.0969465222
64+
1 H 0.0674194961 0.0000000000 0.0484759423
65+
2 H -0.0674195066 0.0000000000 0.0484759501
66+
----------------------------------------------
67+
"""
68+
69+
print("defference for excitation energy between TDA and TDA-ris (in eV)")
70+
print(td.e*27.21138602 - td_ris.energies.get())
71+
print()
72+
"""
73+
[0.25597762 0.05438464 0.01034341 0.00290092 0.01538759]
74+
"""
75+
print("defference for gradient between TDA and TDA-ris (in Hartree/Bohr)")
76+
print(g.de - g_ris.de)
77+
"""
78+
[[-1.05589088e-08 -1.97977506e-15 2.04414538e-03]
79+
[-4.67223270e-03 9.83637092e-16 -1.02206969e-03]
80+
[ 4.67224325e-03 1.00292678e-15 -1.02207751e-03]]
81+
"""
82+
print("norm of the diff")
83+
print(np.linalg.norm(g.de - g_ris.de))
84+
"""
85+
0.007065933384199997
86+
"""

examples/38-tddft_ris_nacv.py

Lines changed: 130 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,130 @@
1+
#!/usr/bin/env python
2+
# Copyright 2021-2025 The PySCF Developers. All Rights Reserved.
3+
#
4+
# Licensed under the Apache License, Version 2.0 (the "License");
5+
# you may not use this file except in compliance with the License.
6+
# You may obtain a copy of the License at
7+
#
8+
# http://www.apache.org/licenses/LICENSE-2.0
9+
#
10+
# Unless required by applicable law or agreed to in writing, software
11+
# distributed under the License is distributed on an "AS IS" BASIS,
12+
# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
13+
# See the License for the specific language governing permissions and
14+
# limitations under the License.
15+
16+
'''
17+
NACV for TDDFT-RIS
18+
'''
19+
20+
# This example will gives the derivative coupling (DC),
21+
# also known as NACME (non-adiabatic coupling matrix element)
22+
# between ground and excited states.
23+
24+
import numpy as np
25+
import pyscf
26+
import gpu4pyscf
27+
from gpu4pyscf.dft import rks
28+
from gpu4pyscf.tdscf import ris
29+
30+
atom = '''
31+
O 0.0000000000 -0.0000000000 0.1174000000
32+
H -0.7570000000 -0.0000000000 -0.4696000000
33+
H 0.7570000000 0.0000000000 -0.4696000000
34+
'''
35+
36+
mol = pyscf.M(atom=atom, basis='def2tzvp')
37+
38+
mf = rks.RKS(mol, xc='pbe0') # -76.3773133945678
39+
mf.conv_tol = 1e-10
40+
mf.kernel()
41+
42+
td = mf.TDA().set(nstates=5)
43+
td.kernel() # [ 7.81949919 9.71029362 10.13398432 12.10163229 13.93675959]
44+
45+
nac = td.nac_method()
46+
nac.states=(1,2)
47+
nac.kernel()
48+
"""
49+
--------- TDA nonadiabatic derivative coupling for states 1 and 2----------
50+
x y z
51+
0 O -0.0975602441 -0.0000000000 -0.0000000000
52+
1 H 0.0548213338 -0.0000000000 0.0360881697
53+
2 H 0.0548213338 0.0000000000 -0.0360881697
54+
--------- TDA nonadiabatic derivative coupling for states 1 and 2 after E scaled (divided by E)----------
55+
x y z
56+
0 O -1.4040391809 -0.0000000000 -0.0000000000
57+
1 H 0.7889617464 -0.0000000000 0.5193632370
58+
2 H 0.7889617464 0.0000000000 -0.5193632370
59+
--------- TDA nonadiabatic derivative coupling for states 1 and 2 with ETF----------
60+
x y z
61+
0 O -0.0965494688 -0.0000000000 -0.0000000000
62+
1 H 0.0482746550 -0.0000000000 0.0378920920
63+
2 H 0.0482746550 0.0000000000 -0.0378920920
64+
--------- TDA nonadiabatic derivative coupling for states 1 and 2 with ETF after E scaled (divided by E)----------
65+
x y z
66+
0 O -1.3894925990 -0.0000000000 -0.0000000000
67+
1 H 0.6947451570 -0.0000000000 0.5453244023
68+
2 H 0.6947451570 0.0000000000 -0.5453244023
69+
----------------------------------------------
70+
"""
71+
72+
td_ris = ris.TDA(mf=mf.to_gpu(), nstates=5, spectra=False, single=False, gram_schmidt=True, Ktrunc=0.0)
73+
td_ris.conv_tol = 1.0E-6
74+
td_ris.kernel() # [ 7.56352157 9.65590898 10.12364072 12.09873136 13.921372 ]
75+
print(td_ris.energies.get())
76+
77+
nac_ris = td_ris.nac_method()
78+
nac_ris.states=(1,2)
79+
nac_ris.kernel()
80+
"""
81+
--------- TDA nonadiabatic derivative coupling for states 1 and 2----------
82+
x y z
83+
0 O 0.1009731844 0.0000000000 -0.0000000014
84+
1 H -0.0575662857 -0.0000000000 -0.0380920213
85+
2 H -0.0575662879 0.0000000000 0.0380920227
86+
--------- TDA nonadiabatic derivative coupling for states 1 and 2 after E scaled (divided by E)----------
87+
x y z
88+
0 O 1.3131508452 0.0000000000 -0.0000000181
89+
1 H -0.7486464564 -0.0000000000 -0.4953846935
90+
2 H -0.7486464849 0.0000000000 0.4953847117
91+
--------- TDA nonadiabatic derivative coupling for states 1 and 2 with ETF----------
92+
x y z
93+
0 O 0.1006324741 0.0000000000 -0.0000000015
94+
1 H -0.0503161533 -0.0000000000 -0.0395091578
95+
2 H -0.0503161556 0.0000000000 0.0395091592
96+
--------- TDA nonadiabatic derivative coupling for states 1 and 2 with ETF after E scaled (divided by E)----------
97+
x y z
98+
0 O 1.3087199253 0.0000000000 -0.0000000189
99+
1 H -0.6543588731 -0.0000000000 -0.5138144769
100+
2 H -0.6543589035 0.0000000000 0.5138144958
101+
----------------------------------------------
102+
"""
103+
104+
print("defference for excitation energy between TDA and TDA-ris (in eV)")
105+
print(td.e*27.21138602 - td_ris.energies.get())
106+
print()
107+
"""
108+
[0.25597762 0.05438464 0.01034359 0.00290093 0.01538759]
109+
"""
110+
print("CIS derivative coupling without ETF")
111+
print(np.abs(nac.de_scaled) - np.abs(nac_ris.de_scaled))
112+
print("norm of difference", np.linalg.norm(np.abs(nac.de_scaled) - np.abs(nac_ris.de_scaled)))
113+
print()
114+
"""
115+
[[ 9.08883357e-02 7.20409145e-15 -1.80700507e-08]
116+
[ 4.03152900e-02 1.15468620e-14 2.39785435e-02]
117+
[ 4.03152615e-02 1.77621884e-16 2.39785253e-02]]
118+
0.11252232092869598
119+
"""
120+
print("difference for CIS derivative coupling with ETF")
121+
print(np.abs(nac.de_etf_scaled) - np.abs(nac_ris.de_etf_scaled))
122+
print("norm of difference", np.linalg.norm(np.abs(nac.de_etf_scaled) - np.abs(nac_ris.de_etf_scaled)))
123+
print()
124+
"""
125+
[[ 8.07726737e-02 1.38040250e-14 -1.88986440e-08]
126+
[ 4.03862838e-02 1.32300610e-14 3.15099254e-02]
127+
[ 4.03862535e-02 -6.54111691e-16 3.15099065e-02]]
128+
0.10849919731015174
129+
"""
130+

gpu4pyscf/df/grad/tdrks_ris.py

Lines changed: 39 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,39 @@
1+
# Copyright 2021-2025 The PySCF Developers. All Rights Reserved.
2+
#
3+
# Licensed under the Apache License, Version 2.0 (the "License");
4+
# you may not use this file except in compliance with the License.
5+
# You may obtain a copy of the License at
6+
#
7+
# http://www.apache.org/licenses/LICENSE-2.0
8+
#
9+
# Unless required by applicable law or agreed to in writing, software
10+
# distributed under the License is distributed on an "AS IS" BASIS,
11+
# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
12+
# See the License for the specific language governing permissions and
13+
# limitations under the License.
14+
15+
from gpu4pyscf.df import int3c2e, df
16+
from gpu4pyscf.df.grad import tdrhf as tdrhf_grad_df
17+
from gpu4pyscf.tdscf import ris
18+
from gpu4pyscf.grad import tdrks_ris as tdrks_ris_grad
19+
from gpu4pyscf import __config__
20+
21+
class Gradients(tdrks_ris_grad.Gradients):
22+
from gpu4pyscf.lib.utils import to_gpu, device
23+
24+
_keys = {'with_df', 'auxbasis_response'}
25+
def __init__(self, td):
26+
# Whether to include the response of DF auxiliary basis when computing
27+
# nuclear gradients of J/K matrices
28+
tdrks_ris_grad.Gradients.__init__(self, td)
29+
30+
auxbasis_response = True
31+
get_jk = tdrhf_grad_df.get_jk
32+
33+
def check_sanity(self):
34+
assert isinstance(self.base._scf, df.df_jk._DFHF)
35+
assert isinstance(self.base, ris.TDDFT) or isinstance(self.base, ris.TDA)
36+
37+
get_veff = tdrhf_grad_df.get_veff
38+
39+
Grad = Gradients

gpu4pyscf/df/nac/__init__.py

Lines changed: 15 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,15 @@
1+
# Copyright 2021-2024 The PySCF Developers. All Rights Reserved.
2+
#
3+
# Licensed under the Apache License, Version 2.0 (the "License");
4+
# you may not use this file except in compliance with the License.
5+
# You may obtain a copy of the License at
6+
#
7+
# http://www.apache.org/licenses/LICENSE-2.0
8+
#
9+
# Unless required by applicable law or agreed to in writing, software
10+
# distributed under the License is distributed on an "AS IS" BASIS,
11+
# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
12+
# See the License for the specific language governing permissions and
13+
# limitations under the License.
14+
15+
from . import tdrhf, tdrks, tdrks_ris

0 commit comments

Comments
 (0)