Line data Source code
1 : !--------------------------------------------------------------------------------------------------!
2 : ! CP2K: A general program to perform molecular dynamics simulations !
3 : ! Copyright 2000-2024 CP2K developers group <https://cp2k.org> !
4 : ! !
5 : ! SPDX-License-Identifier: GPL-2.0-or-later !
6 : !--------------------------------------------------------------------------------------------------!
7 :
8 : ! **************************************************************************************************
9 : !> \brief Definition of the xTB parameter types.
10 : !> \author JGH (10.2018)
11 : ! **************************************************************************************************
12 : ! To be done:
13 : ! 1) Ewald defaults options for GMAX, ALPHA, RCUT
14 : ! 2) QM/MM debugging of forces -- done
15 : ! 3) Periodic displacement field (debugging)
16 : ! 4) Check for RTP and EMD
17 : ! 5) Wannier localization
18 : ! 6) Charge Mixing methods: Broyden/Pulay (more debugging needed, also add to DFTB)
19 : ! **************************************************************************************************
20 : MODULE xtb_types
21 :
22 : USE cp_log_handling, ONLY: cp_get_default_logger,&
23 : cp_logger_type
24 : USE cp_output_handling, ONLY: cp_p_file,&
25 : cp_print_key_finished_output,&
26 : cp_print_key_should_output,&
27 : cp_print_key_unit_nr
28 : USE input_section_types, ONLY: section_vals_type
29 : USE kinds, ONLY: default_string_length,&
30 : dp
31 : #include "./base/base_uses.f90"
32 :
33 : IMPLICIT NONE
34 :
35 : PRIVATE
36 :
37 : ! *** Global parameters ***
38 :
39 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'xtb_types'
40 :
41 : ! **************************************************************************************************
42 : TYPE xtb_atom_type
43 : ! PRIVATE
44 : CHARACTER(LEN=default_string_length) :: typ
45 : CHARACTER(LEN=default_string_length) :: aname
46 : CHARACTER(LEN=2) :: symbol
47 : LOGICAL :: defined
48 : INTEGER :: z !atomic number
49 : REAL(KIND=dp) :: zeff !effective core charge
50 : INTEGER :: natorb !number of orbitals
51 : INTEGER :: lmax !max angular momentum
52 : !
53 : REAL(KIND=dp) :: rcut !cutoff radius for sr-Coulomb
54 : REAL(KIND=dp) :: rcov !covalent radius
55 : REAL(KIND=dp) :: electronegativity !electronegativity
56 : !
57 : REAL(KIND=dp) :: kx !scaling for halogen term
58 : !
59 : REAL(KIND=dp) :: eta !Atomic Hubbard parameter
60 : REAL(KIND=dp) :: xgamma !charge derivative of eta
61 : REAL(KIND=dp) :: alpha !exponential scaling parameter for repulsion potential
62 : REAL(KIND=dp) :: zneff !effective core charge for repulsion potential
63 : ! shell specific parameters
64 : INTEGER :: nshell !number of orbital shells
65 : INTEGER, DIMENSION(5) :: nval ! n-quantum number of shell i
66 : INTEGER, DIMENSION(5) :: lval ! l-quantum number of shell i
67 : INTEGER, DIMENSION(5) :: occupation ! occupation of shell i
68 : REAL(KIND=dp), DIMENSION(5) :: kpoly
69 : REAL(KIND=dp), DIMENSION(5) :: kappa
70 : REAL(KIND=dp), DIMENSION(5) :: hen
71 : REAL(KIND=dp), DIMENSION(5) :: zeta
72 : ! AO to shell pointer
73 : INTEGER, DIMENSION(25) :: nao, lao
74 : ! Upper limit of Mulliken charge
75 : REAL(KIND=dp) :: chmax
76 : END TYPE xtb_atom_type
77 :
78 : ! *** Public data types ***
79 :
80 : PUBLIC :: xtb_atom_type, get_xtb_atom_param, set_xtb_atom_param, write_xtb_atom_param
81 : PUBLIC :: allocate_xtb_atom_param, deallocate_xtb_atom_param
82 :
83 : CONTAINS
84 :
85 : ! **************************************************************************************************
86 : !> \brief ...
87 : !> \param xtb_parameter ...
88 : ! **************************************************************************************************
89 592 : SUBROUTINE allocate_xtb_atom_param(xtb_parameter)
90 :
91 : TYPE(xtb_atom_type), POINTER :: xtb_parameter
92 :
93 592 : IF (ASSOCIATED(xtb_parameter)) &
94 0 : CALL deallocate_xtb_atom_param(xtb_parameter)
95 :
96 592 : ALLOCATE (xtb_parameter)
97 :
98 592 : xtb_parameter%defined = .FALSE.
99 592 : xtb_parameter%aname = ""
100 592 : xtb_parameter%symbol = ""
101 592 : xtb_parameter%typ = "NONE"
102 592 : xtb_parameter%z = -1
103 592 : xtb_parameter%zeff = -1.0_dp
104 592 : xtb_parameter%natorb = 0
105 592 : xtb_parameter%lmax = -1
106 592 : xtb_parameter%rcut = 0.0_dp
107 592 : xtb_parameter%rcov = 0.0_dp
108 592 : xtb_parameter%electronegativity = 0.0_dp
109 592 : xtb_parameter%kx = -100.0_dp
110 592 : xtb_parameter%eta = 0.0_dp
111 592 : xtb_parameter%xgamma = 0.0_dp
112 592 : xtb_parameter%alpha = 0.0_dp
113 592 : xtb_parameter%zneff = 0.0_dp
114 592 : xtb_parameter%nshell = 0
115 3552 : xtb_parameter%nval = 0
116 3552 : xtb_parameter%lval = 0
117 3552 : xtb_parameter%occupation = 0
118 3552 : xtb_parameter%kpoly = 0.0_dp
119 3552 : xtb_parameter%kappa = 0.0_dp
120 3552 : xtb_parameter%hen = 0.0_dp
121 3552 : xtb_parameter%zeta = 0.0_dp
122 15392 : xtb_parameter%nao = 0
123 15392 : xtb_parameter%lao = 0
124 592 : xtb_parameter%chmax = 0.0_dp
125 :
126 592 : END SUBROUTINE allocate_xtb_atom_param
127 :
128 : ! **************************************************************************************************
129 : !> \brief ...
130 : !> \param xtb_parameter ...
131 : ! **************************************************************************************************
132 592 : SUBROUTINE deallocate_xtb_atom_param(xtb_parameter)
133 :
134 : TYPE(xtb_atom_type), POINTER :: xtb_parameter
135 :
136 592 : CPASSERT(ASSOCIATED(xtb_parameter))
137 592 : DEALLOCATE (xtb_parameter)
138 :
139 592 : END SUBROUTINE deallocate_xtb_atom_param
140 :
141 : ! **************************************************************************************************
142 : !> \brief ...
143 : !> \param xtb_parameter ...
144 : !> \param symbol ...
145 : !> \param aname ...
146 : !> \param typ ...
147 : !> \param defined ...
148 : !> \param z ...
149 : !> \param zeff ...
150 : !> \param natorb ...
151 : !> \param lmax ...
152 : !> \param nao ...
153 : !> \param lao ...
154 : !> \param rcut ...
155 : !> \param rcov ...
156 : !> \param kx ...
157 : !> \param eta ...
158 : !> \param xgamma ...
159 : !> \param alpha ...
160 : !> \param zneff ...
161 : !> \param nshell ...
162 : !> \param nval ...
163 : !> \param lval ...
164 : !> \param kpoly ...
165 : !> \param kappa ...
166 : !> \param hen ...
167 : !> \param zeta ...
168 : !> \param occupation ...
169 : !> \param electronegativity ...
170 : !> \param chmax ...
171 : ! **************************************************************************************************
172 41652001 : SUBROUTINE get_xtb_atom_param(xtb_parameter, symbol, aname, typ, defined, z, zeff, natorb, lmax, nao, lao, &
173 : rcut, rcov, kx, eta, xgamma, alpha, zneff, nshell, nval, lval, kpoly, kappa, &
174 : hen, zeta, occupation, electronegativity, chmax)
175 :
176 : TYPE(xtb_atom_type), POINTER :: xtb_parameter
177 : CHARACTER(LEN=2), INTENT(OUT), OPTIONAL :: symbol
178 : CHARACTER(LEN=default_string_length), &
179 : INTENT(OUT), OPTIONAL :: aname, typ
180 : LOGICAL, INTENT(OUT), OPTIONAL :: defined
181 : INTEGER, INTENT(OUT), OPTIONAL :: z
182 : REAL(KIND=dp), INTENT(OUT), OPTIONAL :: zeff
183 : INTEGER, INTENT(OUT), OPTIONAL :: natorb, lmax
184 : INTEGER, DIMENSION(25), INTENT(OUT), OPTIONAL :: nao, lao
185 : REAL(KIND=dp), INTENT(OUT), OPTIONAL :: rcut, rcov, kx, eta, xgamma, alpha, zneff
186 : INTEGER, INTENT(OUT), OPTIONAL :: nshell
187 : INTEGER, DIMENSION(5), INTENT(OUT), OPTIONAL :: nval, lval
188 : REAL(KIND=dp), DIMENSION(5), INTENT(OUT), OPTIONAL :: kpoly, kappa, hen, zeta
189 : INTEGER, DIMENSION(5), INTENT(OUT), OPTIONAL :: occupation
190 : REAL(KIND=dp), INTENT(OUT), OPTIONAL :: electronegativity, chmax
191 :
192 41652001 : CPASSERT(ASSOCIATED(xtb_parameter))
193 :
194 41652001 : IF (PRESENT(symbol)) symbol = xtb_parameter%symbol
195 41652001 : IF (PRESENT(aname)) aname = xtb_parameter%aname
196 41652001 : IF (PRESENT(typ)) typ = xtb_parameter%typ
197 41652001 : IF (PRESENT(defined)) defined = xtb_parameter%defined
198 41652001 : IF (PRESENT(z)) z = xtb_parameter%z
199 41652001 : IF (PRESENT(zeff)) zeff = xtb_parameter%zeff
200 41652001 : IF (PRESENT(natorb)) natorb = xtb_parameter%natorb
201 41652001 : IF (PRESENT(lmax)) lmax = xtb_parameter%lmax
202 90082351 : IF (PRESENT(nao)) nao = xtb_parameter%nao
203 260283251 : IF (PRESENT(lao)) lao = xtb_parameter%lao
204 : !
205 41652001 : IF (PRESENT(rcut)) rcut = xtb_parameter%rcut
206 41652001 : IF (PRESENT(rcov)) rcov = xtb_parameter%rcov
207 41652001 : IF (PRESENT(kx)) kx = xtb_parameter%kx
208 41652001 : IF (PRESENT(electronegativity)) electronegativity = xtb_parameter%electronegativity
209 41652001 : IF (PRESENT(eta)) eta = xtb_parameter%eta
210 41652001 : IF (PRESENT(xgamma)) xgamma = xtb_parameter%xgamma
211 41652001 : IF (PRESENT(alpha)) alpha = xtb_parameter%alpha
212 41652001 : IF (PRESENT(zneff)) zneff = xtb_parameter%zneff
213 41652001 : IF (PRESENT(nshell)) nshell = xtb_parameter%nshell
214 41652001 : IF (PRESENT(nval)) nval = xtb_parameter%nval
215 41794571 : IF (PRESENT(lval)) lval = xtb_parameter%lval
216 41937901 : IF (PRESENT(occupation)) occupation = xtb_parameter%occupation
217 51187421 : IF (PRESENT(kpoly)) kpoly = xtb_parameter%kpoly
218 127257341 : IF (PRESENT(kappa)) kappa = xtb_parameter%kappa
219 51329991 : IF (PRESENT(hen)) hen = xtb_parameter%hen
220 41652001 : IF (PRESENT(zeta)) zeta = xtb_parameter%zeta
221 41652001 : IF (PRESENT(chmax)) chmax = xtb_parameter%chmax
222 :
223 41652001 : END SUBROUTINE get_xtb_atom_param
224 :
225 : ! **************************************************************************************************
226 : !> \brief ...
227 : !> \param xtb_parameter ...
228 : !> \param aname ...
229 : !> \param typ ...
230 : !> \param defined ...
231 : !> \param z ...
232 : !> \param zeff ...
233 : !> \param natorb ...
234 : !> \param lmax ...
235 : !> \param nao ...
236 : !> \param lao ...
237 : !> \param rcut ...
238 : !> \param rcov ...
239 : !> \param kx ...
240 : !> \param eta ...
241 : !> \param xgamma ...
242 : !> \param alpha ...
243 : !> \param zneff ...
244 : !> \param nshell ...
245 : !> \param nval ...
246 : !> \param lval ...
247 : !> \param kpoly ...
248 : !> \param kappa ...
249 : !> \param hen ...
250 : !> \param zeta ...
251 : !> \param electronegativity ...
252 : !> \param occupation ...
253 : !> \param chmax ...
254 : ! **************************************************************************************************
255 0 : SUBROUTINE set_xtb_atom_param(xtb_parameter, aname, typ, defined, z, zeff, natorb, lmax, nao, lao, &
256 : rcut, rcov, kx, eta, xgamma, alpha, zneff, nshell, nval, lval, kpoly, kappa, &
257 : hen, zeta, electronegativity, occupation, chmax)
258 :
259 : TYPE(xtb_atom_type), POINTER :: xtb_parameter
260 : CHARACTER(LEN=default_string_length), INTENT(IN), &
261 : OPTIONAL :: aname, typ
262 : LOGICAL, INTENT(IN), OPTIONAL :: defined
263 : INTEGER, INTENT(IN), OPTIONAL :: z
264 : REAL(KIND=dp), INTENT(IN), OPTIONAL :: zeff
265 : INTEGER, INTENT(IN), OPTIONAL :: natorb, lmax
266 : INTEGER, DIMENSION(25), INTENT(IN), OPTIONAL :: nao, lao
267 : REAL(KIND=dp), INTENT(IN), OPTIONAL :: rcut, rcov, kx, eta, xgamma, alpha, zneff
268 : INTEGER, INTENT(IN), OPTIONAL :: nshell
269 : INTEGER, DIMENSION(5), INTENT(IN), OPTIONAL :: nval, lval
270 : REAL(KIND=dp), DIMENSION(5), INTENT(IN), OPTIONAL :: kpoly, kappa, hen, zeta
271 : REAL(KIND=dp), INTENT(IN), OPTIONAL :: electronegativity
272 : INTEGER, DIMENSION(5), INTENT(IN), OPTIONAL :: occupation
273 : REAL(KIND=dp), INTENT(IN), OPTIONAL :: chmax
274 :
275 0 : CPASSERT(ASSOCIATED(xtb_parameter))
276 :
277 0 : IF (PRESENT(aname)) xtb_parameter%aname = aname
278 0 : IF (PRESENT(typ)) xtb_parameter%typ = typ
279 0 : IF (PRESENT(defined)) xtb_parameter%defined = defined
280 0 : IF (PRESENT(z)) xtb_parameter%z = z
281 0 : IF (PRESENT(zeff)) xtb_parameter%zeff = zeff
282 0 : IF (PRESENT(natorb)) xtb_parameter%natorb = natorb
283 0 : IF (PRESENT(lmax)) xtb_parameter%lmax = lmax
284 0 : IF (PRESENT(nao)) xtb_parameter%nao = nao
285 0 : IF (PRESENT(lao)) xtb_parameter%lao = lao
286 : !
287 0 : IF (PRESENT(rcut)) xtb_parameter%rcut = rcut
288 0 : IF (PRESENT(rcov)) xtb_parameter%rcov = rcov
289 0 : IF (PRESENT(kx)) xtb_parameter%kx = kx
290 0 : IF (PRESENT(electronegativity)) xtb_parameter%electronegativity = electronegativity
291 0 : IF (PRESENT(eta)) xtb_parameter%eta = eta
292 0 : IF (PRESENT(xgamma)) xtb_parameter%xgamma = xgamma
293 0 : IF (PRESENT(alpha)) xtb_parameter%alpha = alpha
294 0 : IF (PRESENT(zneff)) xtb_parameter%zneff = zneff
295 0 : IF (PRESENT(nshell)) xtb_parameter%nshell = nshell
296 0 : IF (PRESENT(nval)) xtb_parameter%nval = nval
297 0 : IF (PRESENT(lval)) xtb_parameter%lval = lval
298 0 : IF (PRESENT(occupation)) xtb_parameter%occupation = occupation
299 0 : IF (PRESENT(kpoly)) xtb_parameter%kpoly = kpoly
300 0 : IF (PRESENT(kappa)) xtb_parameter%kappa = kappa
301 0 : IF (PRESENT(hen)) xtb_parameter%hen = hen
302 0 : IF (PRESENT(zeta)) xtb_parameter%zeta = zeta
303 0 : IF (PRESENT(chmax)) xtb_parameter%chmax = chmax
304 :
305 0 : END SUBROUTINE set_xtb_atom_param
306 :
307 : ! **************************************************************************************************
308 : !> \brief ...
309 : !> \param xtb_parameter ...
310 : !> \param subsys_section ...
311 : ! **************************************************************************************************
312 592 : SUBROUTINE write_xtb_atom_param(xtb_parameter, subsys_section)
313 :
314 : TYPE(xtb_atom_type), POINTER :: xtb_parameter
315 : TYPE(section_vals_type), POINTER :: subsys_section
316 :
317 : CHARACTER(LEN=default_string_length) :: aname, bb
318 : INTEGER :: i, io_unit, m, natorb, nshell
319 : INTEGER, DIMENSION(5) :: lval, nval, occupation
320 : LOGICAL :: defined
321 : REAL(dp) :: zeff
322 : REAL(KIND=dp) :: alpha, en, eta, xgamma, zneff
323 : REAL(KIND=dp), DIMENSION(5) :: hen, kappa, kpoly, zeta
324 : TYPE(cp_logger_type), POINTER :: logger
325 :
326 592 : NULLIFY (logger)
327 592 : logger => cp_get_default_logger()
328 592 : IF (ASSOCIATED(xtb_parameter) .AND. &
329 : BTEST(cp_print_key_should_output(logger%iter_info, subsys_section, &
330 : "PRINT%KINDS/POTENTIAL"), cp_p_file)) THEN
331 :
332 : io_unit = cp_print_key_unit_nr(logger, subsys_section, "PRINT%KINDS", &
333 0 : extension=".Log")
334 :
335 0 : IF (io_unit > 0) THEN
336 0 : CALL get_xtb_atom_param(xtb_parameter, aname=aname, defined=defined, zeff=zeff, natorb=natorb)
337 0 : CALL get_xtb_atom_param(xtb_parameter, nshell=nshell, lval=lval, nval=nval, occupation=occupation)
338 0 : CALL get_xtb_atom_param(xtb_parameter, kpoly=kpoly, kappa=kappa, hen=hen, zeta=zeta)
339 0 : CALL get_xtb_atom_param(xtb_parameter, electronegativity=en, xgamma=xgamma, eta=eta, alpha=alpha, zneff=zneff)
340 :
341 0 : bb = " "
342 0 : WRITE (UNIT=io_unit, FMT="(/,A,T67,A14)") " xTB parameters: ", TRIM(aname)
343 0 : IF (defined) THEN
344 0 : m = 5 - nshell
345 0 : WRITE (UNIT=io_unit, FMT="(T16,A,T71,F10.2)") "Effective core charge:", zeff
346 0 : WRITE (UNIT=io_unit, FMT="(T16,A,T71,I10)") "Number of orbitals:", natorb
347 0 : WRITE (UNIT=io_unit, FMT="(T16,A,T41,A,5(A4,I1,I2,A1))") "Basis set [nl]", bb(1:8*m), &
348 0 : (" [", nval(i), lval(i), "]", i=1, nshell)
349 0 : WRITE (UNIT=io_unit, FMT="(T16,A,T41,A,5F8.3)") "Slater Exponent", bb(1:8*m), (zeta(i), i=1, nshell)
350 0 : WRITE (UNIT=io_unit, FMT="(T16,A,T41,A,5I8)") "Ref. occupation", bb(1:8*m), (occupation(i), i=1, nshell)
351 0 : WRITE (UNIT=io_unit, FMT="(T16,A,T41,A,5F8.3)") "Energy levels [au]", bb(1:8*m), (hen(i), i=1, nshell)
352 0 : WRITE (UNIT=io_unit, FMT="(T16,A,T41,A,5F8.3)") "Kpoly", bb(1:8*m), (kpoly(i), i=1, nshell)
353 0 : WRITE (UNIT=io_unit, FMT="(T16,A,T71,F10.3)") "Electronegativity", en
354 0 : WRITE (UNIT=io_unit, FMT="(T16,A,T71,F10.3)") "Mataga-Nishimoto constant (eta)", eta
355 0 : WRITE (UNIT=io_unit, FMT="(T16,A,T41,A,5F8.3)") "Mataga-Nishimoto scaling kappa", bb(1:8*m), (kappa(i), i=1, nshell)
356 0 : WRITE (UNIT=io_unit, FMT="(T16,A,T71,F10.3)") "3rd Order constant", xgamma
357 0 : WRITE (UNIT=io_unit, FMT="(T16,A,T61,2F10.3)") "Repulsion potential [Z,alpha]", zneff, alpha
358 : ELSE
359 0 : WRITE (UNIT=io_unit, FMT="(T55,A)") "Parameters are not defined"
360 : END IF
361 : END IF
362 0 : CALL cp_print_key_finished_output(io_unit, logger, subsys_section, "PRINT%KINDS")
363 : END IF
364 :
365 592 : END SUBROUTINE write_xtb_atom_param
366 :
367 0 : END MODULE xtb_types
368 :
|
