(git:97501a3)
Loading...
Searching...
No Matches
basis_set_types.F
Go to the documentation of this file.
1!--------------------------------------------------------------------------------------------------!
2! CP2K: A general program to perform molecular dynamics simulations !
3! Copyright 2000-2025 CP2K developers group <https://cp2k.org> !
4! !
5! SPDX-License-Identifier: GPL-2.0-or-later !
6!--------------------------------------------------------------------------------------------------!
7
8! **************************************************************************************************
9!> \par History
10!> - 02.2004 flexible normalization of basis sets [jgh]
11!> - 07.2014 Add a set of contraction coefficient that only work on active
12!> functions
13!> \author Matthias Krack (04.07.2000)
14! **************************************************************************************************
16
17 USE ai_coulomb, ONLY: coulomb2
19 cite_reference
30 USE input_val_types, ONLY: val_get,&
32 USE kinds, ONLY: default_path_length,&
34 dp
35 USE mathconstants, ONLY: dfac,&
36 pi
39 USE orbital_pointers, ONLY: coset,&
40 indco,&
42 nco,&
43 ncoset,&
44 nso,&
45 nsoset
46 USE orbital_symbols, ONLY: cgf_symbol,&
50 USE sto_ng, ONLY: get_sto_ng
54 USE util, ONLY: sort
55#include "../base/base_uses.f90"
56
57 IMPLICIT NONE
58
59 PRIVATE
60
61 ! Global parameters (only in this module)
62
63 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'basis_set_types'
64
65 ! basis set sort criteria
66 INTEGER, PARAMETER, PUBLIC :: basis_sort_default = 0, &
68
69! **************************************************************************************************
70 ! Define the Gaussian-type orbital basis set type
71
73 !MK PRIVATE
74 CHARACTER(LEN=default_string_length) :: name = ""
75 CHARACTER(LEN=default_string_length) :: aliases = ""
76 REAL(kind=dp) :: kind_radius = 0.0_dp
77 REAL(kind=dp) :: short_kind_radius = 0.0_dp
78 INTEGER :: norm_type = -1
79 INTEGER :: ncgf = -1, nset = -1, nsgf = -1
80 CHARACTER(LEN=12), DIMENSION(:), POINTER :: cgf_symbol => null()
81 CHARACTER(LEN=6), DIMENSION(:), POINTER :: sgf_symbol => null()
82 REAL(kind=dp), DIMENSION(:), POINTER :: norm_cgf => null(), set_radius => null()
83 INTEGER, DIMENSION(:), POINTER :: lmax => null(), lmin => null(), &
84 lx => null(), ly => null(), lz => null(), &
85 m => null(), ncgf_set => null(), &
86 npgf => null(), nsgf_set => null(), nshell => null()
87 REAL(kind=dp), DIMENSION(:, :), POINTER :: cphi => null(), pgf_radius => null(), sphi => null(), &
88 scon => null(), zet => null()
89 INTEGER, DIMENSION(:, :), POINTER :: first_cgf => null(), first_sgf => null(), l => null(), &
90 last_cgf => null(), last_sgf => null(), n => null()
91 REAL(kind=dp), DIMENSION(:, :, :), POINTER :: gcc => null()
92 END TYPE gto_basis_set_type
93
95 TYPE(gto_basis_set_type), POINTER :: gto_basis_set => null()
97
98! **************************************************************************************************
99 ! Define the Slater-type orbital basis set type
100
102 PRIVATE
103 CHARACTER(LEN=default_string_length) :: name = ""
104 INTEGER :: nshell = -1
105 CHARACTER(LEN=6), DIMENSION(:), POINTER :: symbol => null()
106 INTEGER, DIMENSION(:), POINTER :: nq => null(), lq => null()
107 REAL(kind=dp), DIMENSION(:), POINTER :: zet => null()
108 END TYPE sto_basis_set_type
109
110! **************************************************************************************************
112 MODULE PROCEDURE read_gto_basis_set1, read_gto_basis_set2
113 END INTERFACE
114! **************************************************************************************************
115
116 ! Public subroutines
117 PUBLIC :: allocate_gto_basis_set, &
131
132 PUBLIC :: allocate_sto_basis_set, &
138
139 ! Public data types
140 PUBLIC :: gto_basis_set_p_type, &
143
144CONTAINS
145
146! **************************************************************************************************
147!> \brief ...
148!> \param gto_basis_set ...
149! **************************************************************************************************
150 SUBROUTINE allocate_gto_basis_set(gto_basis_set)
151
152 ! Allocate a Gaussian-type orbital (GTO) basis set data set.
153
154 ! - Creation (26.10.2000,MK)
155
156 TYPE(gto_basis_set_type), POINTER :: gto_basis_set
157
158 CALL deallocate_gto_basis_set(gto_basis_set)
159
160 ALLOCATE (gto_basis_set)
161
162 END SUBROUTINE allocate_gto_basis_set
163
164! **************************************************************************************************
165!> \brief ...
166!> \param gto_basis_set ...
167! **************************************************************************************************
168 SUBROUTINE deallocate_gto_basis_set(gto_basis_set)
169
170 ! Deallocate a Gaussian-type orbital (GTO) basis set data set.
171
172 ! - Creation (03.11.2000,MK)
173
174 TYPE(gto_basis_set_type), POINTER :: gto_basis_set
175
176 IF (ASSOCIATED(gto_basis_set)) THEN
177 IF (ASSOCIATED(gto_basis_set%cgf_symbol)) DEALLOCATE (gto_basis_set%cgf_symbol)
178 IF (ASSOCIATED(gto_basis_set%sgf_symbol)) DEALLOCATE (gto_basis_set%sgf_symbol)
179 IF (ASSOCIATED(gto_basis_set%norm_cgf)) DEALLOCATE (gto_basis_set%norm_cgf)
180 IF (ASSOCIATED(gto_basis_set%set_radius)) DEALLOCATE (gto_basis_set%set_radius)
181 IF (ASSOCIATED(gto_basis_set%lmax)) DEALLOCATE (gto_basis_set%lmax)
182 IF (ASSOCIATED(gto_basis_set%lmin)) DEALLOCATE (gto_basis_set%lmin)
183 IF (ASSOCIATED(gto_basis_set%lx)) DEALLOCATE (gto_basis_set%lx)
184 IF (ASSOCIATED(gto_basis_set%ly)) DEALLOCATE (gto_basis_set%ly)
185 IF (ASSOCIATED(gto_basis_set%lz)) DEALLOCATE (gto_basis_set%lz)
186 IF (ASSOCIATED(gto_basis_set%m)) DEALLOCATE (gto_basis_set%m)
187 IF (ASSOCIATED(gto_basis_set%ncgf_set)) DEALLOCATE (gto_basis_set%ncgf_set)
188 IF (ASSOCIATED(gto_basis_set%npgf)) DEALLOCATE (gto_basis_set%npgf)
189 IF (ASSOCIATED(gto_basis_set%nsgf_set)) DEALLOCATE (gto_basis_set%nsgf_set)
190 IF (ASSOCIATED(gto_basis_set%nshell)) DEALLOCATE (gto_basis_set%nshell)
191 IF (ASSOCIATED(gto_basis_set%cphi)) DEALLOCATE (gto_basis_set%cphi)
192 IF (ASSOCIATED(gto_basis_set%pgf_radius)) DEALLOCATE (gto_basis_set%pgf_radius)
193 IF (ASSOCIATED(gto_basis_set%sphi)) DEALLOCATE (gto_basis_set%sphi)
194 IF (ASSOCIATED(gto_basis_set%scon)) DEALLOCATE (gto_basis_set%scon)
195 IF (ASSOCIATED(gto_basis_set%zet)) DEALLOCATE (gto_basis_set%zet)
196 IF (ASSOCIATED(gto_basis_set%first_cgf)) DEALLOCATE (gto_basis_set%first_cgf)
197 IF (ASSOCIATED(gto_basis_set%first_sgf)) DEALLOCATE (gto_basis_set%first_sgf)
198 IF (ASSOCIATED(gto_basis_set%l)) DEALLOCATE (gto_basis_set%l)
199 IF (ASSOCIATED(gto_basis_set%last_cgf)) DEALLOCATE (gto_basis_set%last_cgf)
200 IF (ASSOCIATED(gto_basis_set%last_sgf)) DEALLOCATE (gto_basis_set%last_sgf)
201 IF (ASSOCIATED(gto_basis_set%n)) DEALLOCATE (gto_basis_set%n)
202 IF (ASSOCIATED(gto_basis_set%gcc)) DEALLOCATE (gto_basis_set%gcc)
203 DEALLOCATE (gto_basis_set)
204 END IF
205 END SUBROUTINE deallocate_gto_basis_set
206
207! **************************************************************************************************
208!> \brief ...
209!> \param basis_set_in ...
210!> \param basis_set_out ...
211! **************************************************************************************************
212 SUBROUTINE copy_gto_basis_set(basis_set_in, basis_set_out)
213
214 ! Copy a Gaussian-type orbital (GTO) basis set data set.
215
216 TYPE(gto_basis_set_type), INTENT(IN) :: basis_set_in
217 TYPE(gto_basis_set_type), POINTER :: basis_set_out
218
219 INTEGER :: maxco, maxpgf, maxshell, ncgf, nset, nsgf
220
221 CALL allocate_gto_basis_set(basis_set_out)
222
223 basis_set_out%name = basis_set_in%name
224 basis_set_out%aliases = basis_set_in%aliases
225 basis_set_out%kind_radius = basis_set_in%kind_radius
226 basis_set_out%norm_type = basis_set_in%norm_type
227 basis_set_out%nset = basis_set_in%nset
228 basis_set_out%ncgf = basis_set_in%ncgf
229 basis_set_out%nsgf = basis_set_in%nsgf
230 nset = basis_set_in%nset
231 ncgf = basis_set_in%ncgf
232 nsgf = basis_set_in%nsgf
233 ALLOCATE (basis_set_out%cgf_symbol(ncgf))
234 ALLOCATE (basis_set_out%sgf_symbol(nsgf))
235 basis_set_out%cgf_symbol = basis_set_in%cgf_symbol
236 basis_set_out%sgf_symbol = basis_set_in%sgf_symbol
237 ALLOCATE (basis_set_out%norm_cgf(ncgf))
238 basis_set_out%norm_cgf = basis_set_in%norm_cgf
239 ALLOCATE (basis_set_out%set_radius(nset))
240 basis_set_out%set_radius = basis_set_in%set_radius
241 ALLOCATE (basis_set_out%lmax(nset), basis_set_out%lmin(nset), basis_set_out%npgf(nset), basis_set_out%nshell(nset))
242 basis_set_out%lmax = basis_set_in%lmax
243 basis_set_out%lmin = basis_set_in%lmin
244 basis_set_out%npgf = basis_set_in%npgf
245 basis_set_out%nshell = basis_set_in%nshell
246 ALLOCATE (basis_set_out%lx(ncgf), basis_set_out%ly(ncgf), basis_set_out%lz(ncgf), basis_set_out%m(nsgf))
247 basis_set_out%lx = basis_set_in%lx
248 basis_set_out%ly = basis_set_in%ly
249 basis_set_out%lz = basis_set_in%lz
250 basis_set_out%m = basis_set_in%m
251 ALLOCATE (basis_set_out%ncgf_set(nset), basis_set_out%nsgf_set(nset))
252 basis_set_out%ncgf_set = basis_set_in%ncgf_set
253 basis_set_out%nsgf_set = basis_set_in%nsgf_set
254 maxco = SIZE(basis_set_in%cphi, 1)
255 ALLOCATE (basis_set_out%cphi(maxco, ncgf), basis_set_out%sphi(maxco, nsgf), basis_set_out%scon(maxco, nsgf))
256 basis_set_out%cphi = basis_set_in%cphi
257 basis_set_out%sphi = basis_set_in%sphi
258 basis_set_out%scon = basis_set_in%scon
259 maxpgf = maxval(basis_set_in%npgf)
260 ALLOCATE (basis_set_out%pgf_radius(maxpgf, nset), basis_set_out%zet(maxpgf, nset))
261 basis_set_out%pgf_radius = basis_set_in%pgf_radius
262 basis_set_out%zet = basis_set_in%zet
263 maxshell = maxval(basis_set_in%nshell)
264 ALLOCATE (basis_set_out%first_cgf(maxshell, nset), basis_set_out%first_sgf(maxshell, nset))
265 ALLOCATE (basis_set_out%last_cgf(maxshell, nset), basis_set_out%last_sgf(maxshell, nset))
266 basis_set_out%first_cgf = basis_set_in%first_cgf
267 basis_set_out%first_sgf = basis_set_in%first_sgf
268 basis_set_out%last_cgf = basis_set_in%last_cgf
269 basis_set_out%last_sgf = basis_set_in%last_sgf
270 ALLOCATE (basis_set_out%n(maxshell, nset), basis_set_out%l(maxshell, nset))
271 basis_set_out%n = basis_set_in%n
272 basis_set_out%l = basis_set_in%l
273 ALLOCATE (basis_set_out%gcc(maxpgf, maxshell, nset))
274 basis_set_out%gcc = basis_set_in%gcc
275
276 END SUBROUTINE copy_gto_basis_set
277
278! **************************************************************************************************
279!> \brief ...
280!> \param basis_set ...
281!> \param pbasis ...
282!> \param lmax ...
283! **************************************************************************************************
284 SUBROUTINE create_primitive_basis_set(basis_set, pbasis, lmax)
285
286 ! Create a primitives only basis set
287
288 TYPE(gto_basis_set_type), INTENT(IN) :: basis_set
289 TYPE(gto_basis_set_type), POINTER :: pbasis
290 INTEGER, INTENT(IN), OPTIONAL :: lmax
291
292 INTEGER :: i, ico, ip, ipgf, iset, ishell, l, lm, &
293 lshell, m, maxco, mpgf, nc, ncgf, ns, &
294 nset, nsgf
295 INTEGER, ALLOCATABLE, DIMENSION(:) :: nindex, nprim
296 REAL(kind=dp) :: zet0
297 REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :) :: zet, zeta
298
299 mpgf = sum(basis_set%npgf)
300 lm = maxval(basis_set%lmax)
301 ALLOCATE (zet(mpgf, 0:lm), zeta(mpgf, lm + 1), nindex(mpgf), nprim(0:lm))
302 zet = 0.0_dp
303 zeta = 0.0_dp
304 DO l = 0, lm
305 ip = 0
306 DO iset = 1, basis_set%nset
307 IF (basis_set%lmin(iset) <= l .AND. basis_set%lmax(iset) >= l) THEN
308 DO ipgf = 1, basis_set%npgf(iset)
309 ip = ip + 1
310 zet(ip, l) = basis_set%zet(ipgf, iset)
311 END DO
312 END IF
313 END DO
314 nprim(l) = ip
315 END DO
316
317 ! sort exponents
318 DO l = 0, lm
319 zet(1:nprim(l), l) = -zet(1:nprim(l), l)
320 CALL sort(zet(1:nprim(l), l), nprim(l), nindex)
321 ! remove duplicates
322 ip = 0
323 zet0 = 0.0_dp
324 DO i = 1, nprim(l)
325 IF (abs(zet0 - zet(i, l)) > 1.e-6_dp) THEN
326 ip = ip + 1
327 zeta(ip, l + 1) = zet(i, l)
328 END IF
329 END DO
330 nprim(l) = ip
331 !
332 zeta(1:ip, l + 1) = -zeta(1:ip, l + 1)
333 END DO
334
335 CALL allocate_gto_basis_set(pbasis)
336
337 IF (PRESENT(lmax)) THEN
338 ! if requested, reduce max l val
339 lm = min(lm, lmax)
340 END IF
341
342 IF (len_trim(basis_set%name) + 10 > default_string_length) THEN
343 cpwarn("The name of the primitive basis set will be truncated.")
344 END IF
345 pbasis%name = trim(basis_set%name)//"_primitive"
346 pbasis%kind_radius = basis_set%kind_radius
347 pbasis%short_kind_radius = basis_set%short_kind_radius
348 pbasis%norm_type = basis_set%norm_type
349 nset = lm + 1
350 pbasis%nset = nset
351 ALLOCATE (pbasis%lmax(nset), pbasis%lmin(nset), pbasis%npgf(nset), pbasis%nshell(nset))
352 DO iset = 1, nset
353 pbasis%lmax(iset) = iset - 1
354 pbasis%lmin(iset) = iset - 1
355 pbasis%npgf(iset) = nprim(iset - 1)
356 pbasis%nshell(iset) = nprim(iset - 1)
357 END DO
358 pbasis%ncgf = 0
359 pbasis%nsgf = 0
360 DO l = 0, lm
361 pbasis%ncgf = pbasis%ncgf + nprim(l)*((l + 1)*(l + 2))/2
362 pbasis%nsgf = pbasis%nsgf + nprim(l)*(2*l + 1)
363 END DO
364 mpgf = maxval(nprim)
365 ALLOCATE (pbasis%zet(mpgf, nset))
366 pbasis%zet(1:mpgf, 1:nset) = zeta(1:mpgf, 1:nset)
367
368 ALLOCATE (pbasis%l(mpgf, nset), pbasis%n(mpgf, nset))
369 DO iset = 1, nset
370 DO ip = 1, nprim(iset - 1)
371 pbasis%l(ip, iset) = iset - 1
372 pbasis%n(ip, iset) = iset + ip - 1
373 END DO
374 END DO
375
376 ALLOCATE (pbasis%cgf_symbol(pbasis%ncgf))
377 ALLOCATE (pbasis%lx(pbasis%ncgf))
378 ALLOCATE (pbasis%ly(pbasis%ncgf))
379 ALLOCATE (pbasis%lz(pbasis%ncgf))
380 ALLOCATE (pbasis%m(pbasis%nsgf))
381 ALLOCATE (pbasis%sgf_symbol(pbasis%nsgf))
382 ALLOCATE (pbasis%ncgf_set(nset), pbasis%nsgf_set(nset))
383
384 ncgf = 0
385 nsgf = 0
386 DO iset = 1, nset
387 l = iset - 1
388 pbasis%ncgf_set(iset) = nprim(l)*((l + 1)*(l + 2))/2
389 pbasis%nsgf_set(iset) = nprim(l)*(2*l + 1)
390 DO ishell = 1, pbasis%nshell(iset)
391 lshell = pbasis%l(ishell, iset)
392 DO ico = ncoset(lshell - 1) + 1, ncoset(lshell)
393 ncgf = ncgf + 1
394 pbasis%lx(ncgf) = indco(1, ico)
395 pbasis%ly(ncgf) = indco(2, ico)
396 pbasis%lz(ncgf) = indco(3, ico)
397 pbasis%cgf_symbol(ncgf) = &
398 cgf_symbol(pbasis%n(ishell, iset), (/pbasis%lx(ncgf), pbasis%ly(ncgf), pbasis%lz(ncgf)/))
399 END DO
400 DO m = -lshell, lshell
401 nsgf = nsgf + 1
402 pbasis%m(nsgf) = m
403 pbasis%sgf_symbol(nsgf) = sgf_symbol(pbasis%n(ishell, iset), lshell, m)
404 END DO
405 END DO
406 END DO
407 cpassert(ncgf == pbasis%ncgf)
408 cpassert(nsgf == pbasis%nsgf)
409
410 ALLOCATE (pbasis%gcc(mpgf, mpgf, nset))
411 pbasis%gcc = 0.0_dp
412 DO iset = 1, nset
413 DO i = 1, mpgf
414 pbasis%gcc(i, i, iset) = 1.0_dp
415 END DO
416 END DO
417
418 ALLOCATE (pbasis%first_cgf(mpgf, nset))
419 ALLOCATE (pbasis%first_sgf(mpgf, nset))
420 ALLOCATE (pbasis%last_cgf(mpgf, nset))
421 ALLOCATE (pbasis%last_sgf(mpgf, nset))
422 nc = 0
423 ns = 0
424 maxco = 0
425 DO iset = 1, nset
426 DO ishell = 1, pbasis%nshell(iset)
427 lshell = pbasis%l(ishell, iset)
428 pbasis%first_cgf(ishell, iset) = nc + 1
429 nc = nc + nco(lshell)
430 pbasis%last_cgf(ishell, iset) = nc
431 pbasis%first_sgf(ishell, iset) = ns + 1
432 ns = ns + nso(lshell)
433 pbasis%last_sgf(ishell, iset) = ns
434 END DO
435 maxco = max(maxco, pbasis%npgf(iset)*ncoset(pbasis%lmax(iset)))
436 END DO
437
438 ALLOCATE (pbasis%norm_cgf(ncgf))
439 ALLOCATE (pbasis%cphi(maxco, ncgf))
440 pbasis%cphi = 0.0_dp
441 ALLOCATE (pbasis%sphi(maxco, nsgf))
442 pbasis%sphi = 0.0_dp
443 ALLOCATE (pbasis%scon(maxco, ncgf))
444 pbasis%scon = 0.0_dp
445 ALLOCATE (pbasis%set_radius(nset))
446 ALLOCATE (pbasis%pgf_radius(mpgf, nset))
447 pbasis%pgf_radius = 0.0_dp
448
449 CALL init_orb_basis_set(pbasis)
450
451 DEALLOCATE (zet, zeta, nindex, nprim)
452
453 END SUBROUTINE create_primitive_basis_set
454
455! **************************************************************************************************
456!> \brief ...
457!> \param basis_set ...
458!> \param basis_set_add ...
459! **************************************************************************************************
460 SUBROUTINE combine_basis_sets(basis_set, basis_set_add)
461
462 ! Combine two Gaussian-type orbital (GTO) basis sets.
463
464 TYPE(gto_basis_set_type), INTENT(INOUT) :: basis_set
465 TYPE(gto_basis_set_type), INTENT(IN) :: basis_set_add
466
467 CHARACTER(LEN=12), DIMENSION(:), POINTER :: cgf_symbol
468 CHARACTER(LEN=6), DIMENSION(:), POINTER :: sgf_symbol
469 INTEGER :: iset, ishell, lshell, maxco, maxpgf, &
470 maxshell, nc, ncgf, ncgfn, ncgfo, ns, &
471 nset, nsetn, nseto, nsgf, nsgfn, nsgfo
472
473 IF (len_trim(basis_set%name) + len_trim(basis_set_add%name) > default_string_length) THEN
474 cpwarn("The name of the combined GTO basis set will be truncated.")
475 END IF
476 basis_set%name = trim(basis_set%name)//trim(basis_set_add%name)
477 basis_set%nset = basis_set%nset + basis_set_add%nset
478 basis_set%ncgf = basis_set%ncgf + basis_set_add%ncgf
479 basis_set%nsgf = basis_set%nsgf + basis_set_add%nsgf
480 nset = basis_set%nset
481 ncgf = basis_set%ncgf
482 nsgf = basis_set%nsgf
483
484 nsetn = basis_set_add%nset
485 nseto = nset - nsetn
486 CALL reallocate(basis_set%set_radius, 1, nset) ! to be defined later
487 CALL reallocate(basis_set%lmax, 1, nset)
488 CALL reallocate(basis_set%lmin, 1, nset)
489 CALL reallocate(basis_set%npgf, 1, nset)
490 CALL reallocate(basis_set%nshell, 1, nset)
491 basis_set%lmax(nseto + 1:nset) = basis_set_add%lmax(1:nsetn)
492 basis_set%lmin(nseto + 1:nset) = basis_set_add%lmin(1:nsetn)
493 basis_set%npgf(nseto + 1:nset) = basis_set_add%npgf(1:nsetn)
494 basis_set%nshell(nseto + 1:nset) = basis_set_add%nshell(1:nsetn)
495 CALL reallocate(basis_set%ncgf_set, 1, nset)
496 CALL reallocate(basis_set%nsgf_set, 1, nset)
497 basis_set%ncgf_set(nseto + 1:nset) = basis_set_add%ncgf_set(1:nsetn)
498 basis_set%nsgf_set(nseto + 1:nset) = basis_set_add%nsgf_set(1:nsetn)
499
500 nsgfn = basis_set_add%nsgf
501 nsgfo = nsgf - nsgfn
502 ncgfn = basis_set_add%ncgf
503 ncgfo = ncgf - ncgfn
504
505 ALLOCATE (cgf_symbol(ncgf), sgf_symbol(nsgf))
506 cgf_symbol(1:ncgfo) = basis_set%cgf_symbol(1:ncgfo)
507 cgf_symbol(ncgfo + 1:ncgf) = basis_set_add%cgf_symbol(1:ncgfn)
508 sgf_symbol(1:nsgfo) = basis_set%sgf_symbol(1:nsgfo)
509 sgf_symbol(nsgfo + 1:nsgf) = basis_set_add%sgf_symbol(1:nsgfn)
510 DEALLOCATE (basis_set%cgf_symbol, basis_set%sgf_symbol)
511 ALLOCATE (basis_set%cgf_symbol(ncgf), basis_set%sgf_symbol(nsgf))
512 basis_set%cgf_symbol = cgf_symbol
513 basis_set%sgf_symbol = sgf_symbol
514 DEALLOCATE (cgf_symbol, sgf_symbol)
515
516 CALL reallocate(basis_set%lx, 1, ncgf)
517 CALL reallocate(basis_set%ly, 1, ncgf)
518 CALL reallocate(basis_set%lz, 1, ncgf)
519 CALL reallocate(basis_set%m, 1, nsgf)
520 basis_set%lx(ncgfo + 1:ncgf) = basis_set_add%lx(1:ncgfn)
521 basis_set%ly(ncgfo + 1:ncgf) = basis_set_add%ly(1:ncgfn)
522 basis_set%lz(ncgfo + 1:ncgf) = basis_set_add%lz(1:ncgfn)
523 basis_set%m(nsgfo + 1:nsgf) = basis_set_add%m(1:nsgfn)
524
525 maxpgf = maxval(basis_set%npgf)
526 CALL reallocate(basis_set%zet, 1, maxpgf, 1, nset)
527 nc = SIZE(basis_set_add%zet, 1)
528 DO iset = 1, nsetn
529 basis_set%zet(1:nc, nseto + iset) = basis_set_add%zet(1:nc, iset)
530 END DO
531
532 maxshell = maxval(basis_set%nshell)
533 CALL reallocate(basis_set%l, 1, maxshell, 1, nset)
534 CALL reallocate(basis_set%n, 1, maxshell, 1, nset)
535 nc = SIZE(basis_set_add%l, 1)
536 DO iset = 1, nsetn
537 basis_set%l(1:nc, nseto + iset) = basis_set_add%l(1:nc, iset)
538 basis_set%n(1:nc, nseto + iset) = basis_set_add%n(1:nc, iset)
539 END DO
540
541 CALL reallocate(basis_set%first_cgf, 1, maxshell, 1, nset)
542 CALL reallocate(basis_set%first_sgf, 1, maxshell, 1, nset)
543 CALL reallocate(basis_set%last_cgf, 1, maxshell, 1, nset)
544 CALL reallocate(basis_set%last_sgf, 1, maxshell, 1, nset)
545 nc = 0
546 ns = 0
547 DO iset = 1, nset
548 DO ishell = 1, basis_set%nshell(iset)
549 lshell = basis_set%l(ishell, iset)
550 basis_set%first_cgf(ishell, iset) = nc + 1
551 nc = nc + nco(lshell)
552 basis_set%last_cgf(ishell, iset) = nc
553 basis_set%first_sgf(ishell, iset) = ns + 1
554 ns = ns + nso(lshell)
555 basis_set%last_sgf(ishell, iset) = ns
556 END DO
557 END DO
558
559 CALL reallocate(basis_set%gcc, 1, maxpgf, 1, maxshell, 1, nset)
560 nc = SIZE(basis_set_add%gcc, 1)
561 ns = SIZE(basis_set_add%gcc, 2)
562 DO iset = 1, nsetn
563 basis_set%gcc(1:nc, 1:ns, nseto + iset) = basis_set_add%gcc(1:nc, 1:ns, iset)
564 END DO
565
566 ! these arrays are determined later using initialization calls
567 CALL reallocate(basis_set%norm_cgf, 1, ncgf)
568 maxco = max(SIZE(basis_set%cphi, 1), SIZE(basis_set_add%cphi, 1))
569 CALL reallocate(basis_set%cphi, 1, maxco, 1, ncgf)
570 CALL reallocate(basis_set%sphi, 1, maxco, 1, nsgf)
571 CALL reallocate(basis_set%scon, 1, maxco, 1, nsgf)
572 CALL reallocate(basis_set%pgf_radius, 1, maxpgf, 1, nset)
573
574 END SUBROUTINE combine_basis_sets
575
576! **************************************************************************************************
577!> \brief ...
578!> \param gto_basis_set ...
579!> \param name ...
580!> \param aliases ...
581!> \param norm_type ...
582!> \param kind_radius ...
583!> \param ncgf ...
584!> \param nset ...
585!> \param nsgf ...
586!> \param cgf_symbol ...
587!> \param sgf_symbol ...
588!> \param norm_cgf ...
589!> \param set_radius ...
590!> \param lmax ...
591!> \param lmin ...
592!> \param lx ...
593!> \param ly ...
594!> \param lz ...
595!> \param m ...
596!> \param ncgf_set ...
597!> \param npgf ...
598!> \param nsgf_set ...
599!> \param nshell ...
600!> \param cphi ...
601!> \param pgf_radius ...
602!> \param sphi ...
603!> \param scon ...
604!> \param zet ...
605!> \param first_cgf ...
606!> \param first_sgf ...
607!> \param l ...
608!> \param last_cgf ...
609!> \param last_sgf ...
610!> \param n ...
611!> \param gcc ...
612!> \param maxco ...
613!> \param maxl ...
614!> \param maxpgf ...
615!> \param maxsgf_set ...
616!> \param maxshell ...
617!> \param maxso ...
618!> \param nco_sum ...
619!> \param npgf_sum ...
620!> \param nshell_sum ...
621!> \param maxder ...
622!> \param short_kind_radius ...
623!> \param npgf_seg_sum number of primitives in "segmented contraction format"
624! **************************************************************************************************
625 SUBROUTINE get_gto_basis_set(gto_basis_set, name, aliases, norm_type, kind_radius, ncgf, &
626 nset, nsgf, cgf_symbol, sgf_symbol, norm_cgf, set_radius, lmax, lmin, lx, ly, lz, &
627 m, ncgf_set, npgf, nsgf_set, nshell, cphi, pgf_radius, sphi, scon, zet, first_cgf, first_sgf, l, &
628 last_cgf, last_sgf, n, gcc, maxco, maxl, maxpgf, maxsgf_set, maxshell, maxso, nco_sum, &
629 npgf_sum, nshell_sum, maxder, short_kind_radius, npgf_seg_sum)
630
631 ! Get informations about a Gaussian-type orbital (GTO) basis set.
632
633 ! - Creation (10.01.2002,MK)
634
635 TYPE(gto_basis_set_type), INTENT(IN) :: gto_basis_set
636 CHARACTER(LEN=default_string_length), &
637 INTENT(OUT), OPTIONAL :: name, aliases
638 INTEGER, INTENT(OUT), OPTIONAL :: norm_type
639 REAL(kind=dp), INTENT(OUT), OPTIONAL :: kind_radius
640 INTEGER, INTENT(OUT), OPTIONAL :: ncgf, nset, nsgf
641 CHARACTER(LEN=12), DIMENSION(:), OPTIONAL, POINTER :: cgf_symbol
642 CHARACTER(LEN=6), DIMENSION(:), OPTIONAL, POINTER :: sgf_symbol
643 REAL(kind=dp), DIMENSION(:), OPTIONAL, POINTER :: norm_cgf, set_radius
644 INTEGER, DIMENSION(:), OPTIONAL, POINTER :: lmax, lmin, lx, ly, lz, m, ncgf_set, &
645 npgf, nsgf_set, nshell
646 REAL(kind=dp), DIMENSION(:, :), OPTIONAL, POINTER :: cphi, pgf_radius, sphi, scon, zet
647 INTEGER, DIMENSION(:, :), OPTIONAL, POINTER :: first_cgf, first_sgf, l, last_cgf, &
648 last_sgf, n
649 REAL(kind=dp), DIMENSION(:, :, :), OPTIONAL, &
650 POINTER :: gcc
651 INTEGER, INTENT(OUT), OPTIONAL :: maxco, maxl, maxpgf, maxsgf_set, &
652 maxshell, maxso, nco_sum, npgf_sum, &
653 nshell_sum
654 INTEGER, INTENT(IN), OPTIONAL :: maxder
655 REAL(kind=dp), INTENT(OUT), OPTIONAL :: short_kind_radius
656 INTEGER, INTENT(OUT), OPTIONAL :: npgf_seg_sum
657
658 INTEGER :: iset, nder
659
660 IF (PRESENT(name)) name = gto_basis_set%name
661 IF (PRESENT(aliases)) aliases = gto_basis_set%aliases
662 IF (PRESENT(norm_type)) norm_type = gto_basis_set%norm_type
663 IF (PRESENT(kind_radius)) kind_radius = gto_basis_set%kind_radius
664 IF (PRESENT(short_kind_radius)) short_kind_radius = gto_basis_set%short_kind_radius
665 IF (PRESENT(ncgf)) ncgf = gto_basis_set%ncgf
666 IF (PRESENT(nset)) nset = gto_basis_set%nset
667 IF (PRESENT(nsgf)) nsgf = gto_basis_set%nsgf
668 IF (PRESENT(cgf_symbol)) cgf_symbol => gto_basis_set%cgf_symbol
669 IF (PRESENT(sgf_symbol)) sgf_symbol => gto_basis_set%sgf_symbol
670 IF (PRESENT(norm_cgf)) norm_cgf => gto_basis_set%norm_cgf
671 IF (PRESENT(set_radius)) set_radius => gto_basis_set%set_radius
672 IF (PRESENT(lmax)) lmax => gto_basis_set%lmax
673 IF (PRESENT(lmin)) lmin => gto_basis_set%lmin
674 IF (PRESENT(lx)) lx => gto_basis_set%lx
675 IF (PRESENT(ly)) ly => gto_basis_set%ly
676 IF (PRESENT(lz)) lz => gto_basis_set%lz
677 IF (PRESENT(m)) m => gto_basis_set%m
678 IF (PRESENT(ncgf_set)) ncgf_set => gto_basis_set%ncgf_set
679 IF (PRESENT(npgf)) npgf => gto_basis_set%npgf
680 IF (PRESENT(nsgf_set)) nsgf_set => gto_basis_set%nsgf_set
681 IF (PRESENT(nshell)) nshell => gto_basis_set%nshell
682 IF (PRESENT(cphi)) cphi => gto_basis_set%cphi
683 IF (PRESENT(pgf_radius)) pgf_radius => gto_basis_set%pgf_radius
684 IF (PRESENT(sphi)) sphi => gto_basis_set%sphi
685 IF (PRESENT(scon)) scon => gto_basis_set%scon
686 IF (PRESENT(zet)) zet => gto_basis_set%zet
687 IF (PRESENT(first_cgf)) first_cgf => gto_basis_set%first_cgf
688 IF (PRESENT(first_sgf)) first_sgf => gto_basis_set%first_sgf
689 IF (PRESENT(l)) l => gto_basis_set%l
690 IF (PRESENT(last_cgf)) last_cgf => gto_basis_set%last_cgf
691 IF (PRESENT(last_sgf)) last_sgf => gto_basis_set%last_sgf
692 IF (PRESENT(n)) n => gto_basis_set%n
693 IF (PRESENT(gcc)) gcc => gto_basis_set%gcc
694 IF (PRESENT(maxco)) THEN
695 maxco = 0
696 IF (PRESENT(maxder)) THEN
697 nder = maxder
698 ELSE
699 nder = 0
700 END IF
701 DO iset = 1, gto_basis_set%nset
702 maxco = max(maxco, gto_basis_set%npgf(iset)* &
703 ncoset(gto_basis_set%lmax(iset) + nder))
704 END DO
705 END IF
706 IF (PRESENT(maxl)) THEN
707 maxl = -1
708 DO iset = 1, gto_basis_set%nset
709 maxl = max(maxl, gto_basis_set%lmax(iset))
710 END DO
711 END IF
712 IF (PRESENT(maxpgf)) THEN
713 maxpgf = 0
714 DO iset = 1, gto_basis_set%nset
715 maxpgf = max(maxpgf, gto_basis_set%npgf(iset))
716 END DO
717 END IF
718 IF (PRESENT(maxsgf_set)) THEN
719 maxsgf_set = 0
720 DO iset = 1, gto_basis_set%nset
721 maxsgf_set = max(maxsgf_set, gto_basis_set%nsgf_set(iset))
722 END DO
723 END IF
724 IF (PRESENT(maxshell)) THEN ! MAXVAL on structure component avoided
725 maxshell = 0
726 DO iset = 1, gto_basis_set%nset
727 maxshell = max(maxshell, gto_basis_set%nshell(iset))
728 END DO
729 END IF
730 IF (PRESENT(maxso)) THEN
731 maxso = 0
732 DO iset = 1, gto_basis_set%nset
733 maxso = max(maxso, gto_basis_set%npgf(iset)* &
734 nsoset(gto_basis_set%lmax(iset)))
735 END DO
736 END IF
737
738 IF (PRESENT(nco_sum)) THEN
739 nco_sum = 0
740 DO iset = 1, gto_basis_set%nset
741 nco_sum = nco_sum + gto_basis_set%npgf(iset)* &
742 ncoset(gto_basis_set%lmax(iset))
743 END DO
744 END IF
745 IF (PRESENT(npgf_sum)) npgf_sum = sum(gto_basis_set%npgf)
746 IF (PRESENT(nshell_sum)) nshell_sum = sum(gto_basis_set%nshell)
747 IF (PRESENT(npgf_seg_sum)) THEN
748 npgf_seg_sum = 0
749 DO iset = 1, gto_basis_set%nset
750 npgf_seg_sum = npgf_seg_sum + gto_basis_set%npgf(iset)*gto_basis_set%nshell(iset)
751 END DO
752 END IF
753
754 END SUBROUTINE get_gto_basis_set
755
756! **************************************************************************************************
757!> \brief ...
758!> \param gto_basis_set ...
759! **************************************************************************************************
760 SUBROUTINE init_aux_basis_set(gto_basis_set)
761
762 ! Initialise a Gaussian-type orbital (GTO) basis set data set.
763
764 ! - Creation (06.12.2000,MK)
765
766 TYPE(gto_basis_set_type), POINTER :: gto_basis_set
767
768 CHARACTER(len=*), PARAMETER :: routinen = 'init_aux_basis_set'
769
770 INTEGER :: handle
771
772! -------------------------------------------------------------------------
773
774 IF (.NOT. ASSOCIATED(gto_basis_set)) RETURN
775
776 CALL timeset(routinen, handle)
777
778 SELECT CASE (gto_basis_set%norm_type)
779 CASE (0)
780 ! No normalisation requested
781 CASE (1)
782 CALL init_norm_cgf_aux_2(gto_basis_set)
783 CASE (2)
784 ! WARNING this was never tested
785 CALL init_norm_cgf_aux(gto_basis_set)
786 CASE DEFAULT
787 cpabort("Normalization method not specified")
788 END SELECT
789
790 ! Initialise the transformation matrices "pgf" -> "cgf"
791 CALL init_cphi_and_sphi(gto_basis_set)
792
793 CALL timestop(handle)
794
795 END SUBROUTINE init_aux_basis_set
796
797! **************************************************************************************************
798!> \brief ...
799!> \param gto_basis_set ...
800! **************************************************************************************************
801 SUBROUTINE init_cphi_and_sphi(gto_basis_set)
802
803 ! Initialise the matrices for the transformation of primitive Cartesian
804 ! Gaussian-type functions to contracted Cartesian (cphi) and spherical
805 ! (sphi) Gaussian-type functions.
806
807 ! - Creation (20.09.2000,MK)
808
809 TYPE(gto_basis_set_type), INTENT(INOUT) :: gto_basis_set
810
811 CHARACTER(len=*), PARAMETER :: routinen = 'init_cphi_and_sphi'
812
813 INTEGER :: first_cgf, first_sgf, handle, icgf, ico, &
814 ipgf, iset, ishell, l, last_sgf, lmax, &
815 lmin, n, n1, n2, ncgf, nn, nn1, nn2, &
816 npgf, nsgf
817
818! -------------------------------------------------------------------------
819! Build the Cartesian transformation matrix "cphi"
820
821 CALL timeset(routinen, handle)
822
823 gto_basis_set%cphi = 0.0_dp
824 DO iset = 1, gto_basis_set%nset
825 n = ncoset(gto_basis_set%lmax(iset))
826 DO ishell = 1, gto_basis_set%nshell(iset)
827 DO icgf = gto_basis_set%first_cgf(ishell, iset), &
828 gto_basis_set%last_cgf(ishell, iset)
829 ico = coset(gto_basis_set%lx(icgf), &
830 gto_basis_set%ly(icgf), &
831 gto_basis_set%lz(icgf))
832 DO ipgf = 1, gto_basis_set%npgf(iset)
833 gto_basis_set%cphi(ico, icgf) = gto_basis_set%norm_cgf(icgf)* &
834 gto_basis_set%gcc(ipgf, ishell, iset)
835 ico = ico + n
836 END DO
837 END DO
838 END DO
839 END DO
840
841 ! Build the spherical transformation matrix "sphi"
842
843 n = SIZE(gto_basis_set%cphi, 1)
844
845 gto_basis_set%sphi = 0.0_dp
846 IF (n .GT. 0) THEN
847 lmax = -1
848 ! Ensure proper setup of orbtramat
849 DO iset = 1, gto_basis_set%nset
850 DO ishell = 1, gto_basis_set%nshell(iset)
851 lmax = max(lmax, gto_basis_set%l(ishell, iset))
852 END DO
853 END DO
854 CALL init_spherical_harmonics(lmax, -1)
855
856 DO iset = 1, gto_basis_set%nset
857 DO ishell = 1, gto_basis_set%nshell(iset)
858 l = gto_basis_set%l(ishell, iset)
859 first_cgf = gto_basis_set%first_cgf(ishell, iset)
860 first_sgf = gto_basis_set%first_sgf(ishell, iset)
861 ncgf = nco(l)
862 nsgf = nso(l)
863 CALL dgemm("N", "T", n, nsgf, ncgf, &
864 1.0_dp, gto_basis_set%cphi(1, first_cgf), n, &
865 orbtramat(l)%c2s(1, 1), nsgf, &
866 0.0_dp, gto_basis_set%sphi(1, first_sgf), n)
867 END DO
868 END DO
869 END IF
870
871 ! Build the reduced transformation matrix "scon"
872 ! This matrix transforms from cartesian primitifs to spherical contracted functions
873 ! "scon" only includes primitifs (lmin -> lmax), whereas sphi is (0 -> lmax)
874
875 n = SIZE(gto_basis_set%scon, 1)
876
877 gto_basis_set%scon = 0.0_dp
878 IF (n .GT. 0) THEN
879 DO iset = 1, gto_basis_set%nset
880 lmin = gto_basis_set%lmin(iset)
881 lmax = gto_basis_set%lmax(iset)
882 npgf = gto_basis_set%npgf(iset)
883 nn = ncoset(lmax) - ncoset(lmin - 1)
884 DO ishell = 1, gto_basis_set%nshell(iset)
885 first_sgf = gto_basis_set%first_sgf(ishell, iset)
886 last_sgf = gto_basis_set%last_sgf(ishell, iset)
887 DO ipgf = 1, npgf
888 nn1 = (ipgf - 1)*ncoset(lmax) + ncoset(lmin - 1) + 1
889 nn2 = ipgf*ncoset(lmax)
890 n1 = (ipgf - 1)*nn + 1
891 n2 = ipgf*nn
892 gto_basis_set%scon(n1:n2, first_sgf:last_sgf) = gto_basis_set%sphi(nn1:nn2, first_sgf:last_sgf)
893 END DO
894 END DO
895 END DO
896 END IF
897
898 CALL timestop(handle)
899
900 END SUBROUTINE init_cphi_and_sphi
901
902! **************************************************************************************************
903!> \brief ...
904!> \param gto_basis_set ...
905! **************************************************************************************************
906 SUBROUTINE init_norm_cgf_aux(gto_basis_set)
907
908 ! Initialise the normalization factors of the contracted Cartesian Gaussian
909 ! functions, if the Gaussian functions represent charge distributions.
910
911 ! - Creation (07.12.2000,MK)
912
913 TYPE(gto_basis_set_type), INTENT(INOUT) :: gto_basis_set
914
915 INTEGER :: icgf, ico, ipgf, iset, ishell, jco, &
916 jpgf, ll, lmax, lmin, lx, ly, lz, n, &
917 npgfa
918 REAL(kind=dp) :: fnorm, gcca, gccb
919 REAL(kind=dp), ALLOCATABLE, DIMENSION(:) :: ff
920 REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :) :: gaa
921 REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :, :) :: vv
922 REAL(kind=dp), DIMENSION(:), POINTER :: rpgfa, zeta
923
924! -------------------------------------------------------------------------
925
926 n = 0
927 ll = 0
928 DO iset = 1, gto_basis_set%nset
929 n = max(n, gto_basis_set%npgf(iset)*ncoset(gto_basis_set%lmax(iset)))
930 ll = max(ll, gto_basis_set%lmax(iset))
931 END DO
932
933 ALLOCATE (gaa(n, n))
934 ALLOCATE (vv(ncoset(ll), ncoset(ll), ll + ll + 1))
935 ALLOCATE (ff(0:ll + ll))
936
937 DO iset = 1, gto_basis_set%nset
938 lmax = gto_basis_set%lmax(iset)
939 lmin = gto_basis_set%lmin(iset)
940 n = ncoset(lmax)
941 npgfa = gto_basis_set%npgf(iset)
942 rpgfa => gto_basis_set%pgf_radius(1:npgfa, iset)
943 zeta => gto_basis_set%zet(1:npgfa, iset)
944 CALL coulomb2(lmax, npgfa, zeta, rpgfa, lmin, &
945 lmax, npgfa, zeta, rpgfa, lmin, &
946 (/0.0_dp, 0.0_dp, 0.0_dp/), 0.0_dp, gaa, vv, ff(0:))
947 DO ishell = 1, gto_basis_set%nshell(iset)
948 DO icgf = gto_basis_set%first_cgf(ishell, iset), &
949 gto_basis_set%last_cgf(ishell, iset)
950 lx = gto_basis_set%lx(icgf)
951 ly = gto_basis_set%ly(icgf)
952 lz = gto_basis_set%lz(icgf)
953 ico = coset(lx, ly, lz)
954 fnorm = 0.0_dp
955 DO ipgf = 1, npgfa
956 gcca = gto_basis_set%gcc(ipgf, ishell, iset)
957 jco = coset(lx, ly, lz)
958 DO jpgf = 1, npgfa
959 gccb = gto_basis_set%gcc(jpgf, ishell, iset)
960 fnorm = fnorm + gcca*gccb*gaa(ico, jco)
961 jco = jco + n
962 END DO
963 ico = ico + n
964 END DO
965 gto_basis_set%norm_cgf(icgf) = 1.0_dp/sqrt(fnorm)
966 END DO
967 END DO
968 END DO
969
970 DEALLOCATE (vv, ff)
971 DEALLOCATE (gaa)
972
973 END SUBROUTINE init_norm_cgf_aux
974
975! **************************************************************************************************
976!> \brief ...
977!> \param gto_basis_set ...
978! **************************************************************************************************
979 ELEMENTAL SUBROUTINE init_norm_cgf_aux_2(gto_basis_set)
980
981 ! Initialise the normalization factors of the auxiliary Cartesian Gaussian
982 ! functions (Kim-Gordon polarization basis) Norm = 1.
983
984 ! - Creation (07.12.2000,GT)
985
986 TYPE(gto_basis_set_type), INTENT(INOUT) :: gto_basis_set
987
988 INTEGER :: icgf, iset, ishell
989
990 DO iset = 1, gto_basis_set%nset
991 DO ishell = 1, gto_basis_set%nshell(iset)
992 DO icgf = gto_basis_set%first_cgf(ishell, iset), &
993 gto_basis_set%last_cgf(ishell, iset)
994 gto_basis_set%norm_cgf(icgf) = 1.0_dp
995 END DO
996 END DO
997 END DO
998
999 END SUBROUTINE init_norm_cgf_aux_2
1000
1001! **************************************************************************************************
1002!> \brief Initialise the normalization factors of the contracted Cartesian Gaussian functions.
1003!> \param gto_basis_set ...
1004!> \author MK
1005! **************************************************************************************************
1006 ELEMENTAL SUBROUTINE init_norm_cgf_orb(gto_basis_set)
1007
1008 TYPE(gto_basis_set_type), INTENT(INOUT) :: gto_basis_set
1009
1010 INTEGER :: icgf, ipgf, iset, ishell, jpgf, l, lx, &
1011 ly, lz
1012 REAL(kind=dp) :: expzet, fnorm, gcca, gccb, prefac, zeta, &
1013 zetb
1014
1015 DO iset = 1, gto_basis_set%nset
1016 DO ishell = 1, gto_basis_set%nshell(iset)
1017
1018 l = gto_basis_set%l(ishell, iset)
1019
1020 expzet = 0.5_dp*real(2*l + 3, dp)
1021
1022 fnorm = 0.0_dp
1023
1024 DO ipgf = 1, gto_basis_set%npgf(iset)
1025 gcca = gto_basis_set%gcc(ipgf, ishell, iset)
1026 zeta = gto_basis_set%zet(ipgf, iset)
1027 DO jpgf = 1, gto_basis_set%npgf(iset)
1028 gccb = gto_basis_set%gcc(jpgf, ishell, iset)
1029 zetb = gto_basis_set%zet(jpgf, iset)
1030 fnorm = fnorm + gcca*gccb/(zeta + zetb)**expzet
1031 END DO
1032 END DO
1033
1034 fnorm = 0.5_dp**l*pi**1.5_dp*fnorm
1035
1036 DO icgf = gto_basis_set%first_cgf(ishell, iset), &
1037 gto_basis_set%last_cgf(ishell, iset)
1038 lx = gto_basis_set%lx(icgf)
1039 ly = gto_basis_set%ly(icgf)
1040 lz = gto_basis_set%lz(icgf)
1041 prefac = dfac(2*lx - 1)*dfac(2*ly - 1)*dfac(2*lz - 1)
1042 gto_basis_set%norm_cgf(icgf) = 1.0_dp/sqrt(prefac*fnorm)
1043 END DO
1044
1045 END DO
1046 END DO
1047
1048 END SUBROUTINE init_norm_cgf_orb
1049
1050! **************************************************************************************************
1051!> \brief Initialise the normalization factors of the contracted Cartesian Gaussian
1052!> functions used for frozen density representation.
1053!> \param gto_basis_set ...
1054!> \author GT
1055! **************************************************************************************************
1056 ELEMENTAL SUBROUTINE init_norm_cgf_orb_den(gto_basis_set)
1057
1058 TYPE(gto_basis_set_type), INTENT(INOUT) :: gto_basis_set
1059
1060 INTEGER :: icgf, ipgf, iset, ishell, l
1061 REAL(kind=dp) :: expzet, gcca, prefac, zeta
1062
1063 DO iset = 1, gto_basis_set%nset
1064 DO ishell = 1, gto_basis_set%nshell(iset)
1065 l = gto_basis_set%l(ishell, iset)
1066 expzet = 0.5_dp*real(2*l + 3, dp)
1067 prefac = (1.0_dp/pi)**1.5_dp
1068 DO ipgf = 1, gto_basis_set%npgf(iset)
1069 gcca = gto_basis_set%gcc(ipgf, ishell, iset)
1070 zeta = gto_basis_set%zet(ipgf, iset)
1071 gto_basis_set%gcc(ipgf, ishell, iset) = prefac*zeta**expzet*gcca
1072 END DO
1073 DO icgf = gto_basis_set%first_cgf(ishell, iset), &
1074 gto_basis_set%last_cgf(ishell, iset)
1075 gto_basis_set%norm_cgf(icgf) = 1.0_dp
1076 END DO
1077 END DO
1078 END DO
1079
1080 END SUBROUTINE init_norm_cgf_orb_den
1081
1082! **************************************************************************************************
1083!> \brief Initialise a Gaussian-type orbital (GTO) basis set data set.
1084!> \param gto_basis_set ...
1085!> \author MK
1086! **************************************************************************************************
1087 SUBROUTINE init_orb_basis_set(gto_basis_set)
1088
1089 TYPE(gto_basis_set_type), POINTER :: gto_basis_set
1090
1091 CHARACTER(len=*), PARAMETER :: routinen = 'init_orb_basis_set'
1092
1093 INTEGER :: handle
1094
1095! -------------------------------------------------------------------------
1096
1097 IF (.NOT. ASSOCIATED(gto_basis_set)) RETURN
1098
1099 CALL timeset(routinen, handle)
1100
1101 SELECT CASE (gto_basis_set%norm_type)
1102 CASE (0)
1103 ! No normalisation requested
1104 CASE (1)
1105 CALL init_norm_cgf_orb_den(gto_basis_set)
1106 CASE (2)
1107 ! Normalise the primitive Gaussian functions
1108 CALL normalise_gcc_orb(gto_basis_set)
1109 ! Compute the normalization factors of the contracted Gaussian-type
1110 ! functions
1111 CALL init_norm_cgf_orb(gto_basis_set)
1112 CASE (3)
1113 CALL init_norm_cgf_orb(gto_basis_set)
1114 CASE DEFAULT
1115 cpabort("Normalization method not specified")
1116 END SELECT
1117
1118 ! Initialise the transformation matrices "pgf" -> "cgf"
1119
1120 CALL init_cphi_and_sphi(gto_basis_set)
1121
1122 CALL timestop(handle)
1123
1124 END SUBROUTINE init_orb_basis_set
1125
1126! **************************************************************************************************
1127!> \brief Normalise the primitive Cartesian Gaussian functions. The normalization
1128!> factor is included in the Gaussian contraction coefficients.
1129!> \param gto_basis_set ...
1130!> \author MK
1131! **************************************************************************************************
1132 SUBROUTINE normalise_gcc_orb(gto_basis_set)
1133
1134 TYPE(gto_basis_set_type), POINTER :: gto_basis_set
1135
1136 INTEGER :: ipgf, iset, ishell, l
1137 REAL(kind=dp) :: expzet, gcca, prefac, zeta
1138
1139 DO iset = 1, gto_basis_set%nset
1140 DO ishell = 1, gto_basis_set%nshell(iset)
1141 l = gto_basis_set%l(ishell, iset)
1142 expzet = 0.25_dp*real(2*l + 3, dp)
1143 prefac = 2.0_dp**l*(2.0_dp/pi)**0.75_dp
1144 DO ipgf = 1, gto_basis_set%npgf(iset)
1145 gcca = gto_basis_set%gcc(ipgf, ishell, iset)
1146 zeta = gto_basis_set%zet(ipgf, iset)
1147 gto_basis_set%gcc(ipgf, ishell, iset) = prefac*zeta**expzet*gcca
1148 END DO
1149 END DO
1150 END DO
1151
1152 END SUBROUTINE normalise_gcc_orb
1153
1154! **************************************************************************************************
1155!> \brief Read a Gaussian-type orbital (GTO) basis set from the database file.
1156!> \param element_symbol ...
1157!> \param basis_set_name ...
1158!> \param gto_basis_set ...
1159!> \param para_env ...
1160!> \param dft_section ...
1161!> \author MK
1162! **************************************************************************************************
1163 SUBROUTINE read_gto_basis_set1(element_symbol, basis_set_name, gto_basis_set, &
1164 para_env, dft_section)
1165
1166 CHARACTER(LEN=*), INTENT(IN) :: element_symbol, basis_set_name
1167 TYPE(gto_basis_set_type), INTENT(INOUT) :: gto_basis_set
1168 TYPE(mp_para_env_type), POINTER :: para_env
1169 TYPE(section_vals_type), POINTER :: dft_section
1170
1171 CHARACTER(LEN=240) :: line
1172 CHARACTER(LEN=242) :: line2
1173 CHARACTER(len=default_path_length) :: basis_set_file_name, tmp
1174 CHARACTER(LEN=default_path_length), DIMENSION(:), &
1175 POINTER :: cbasis
1176 CHARACTER(LEN=LEN(basis_set_name)) :: bsname
1177 CHARACTER(LEN=LEN(basis_set_name)+2) :: bsname2
1178 CHARACTER(LEN=LEN(element_symbol)) :: symbol
1179 CHARACTER(LEN=LEN(element_symbol)+2) :: symbol2
1180 INTEGER :: i, ibasis, ico, ipgf, irep, iset, ishell, lshell, m, maxco, maxl, maxpgf, &
1181 maxshell, nbasis, ncgf, nmin, nset, nsgf, sort_method, strlen1, strlen2
1182 INTEGER, DIMENSION(:), POINTER :: lmax, lmin, npgf, nshell
1183 INTEGER, DIMENSION(:, :), POINTER :: l, n
1184 LOGICAL :: basis_found, found, match
1185 REAL(KIND=dp), DIMENSION(:, :), POINTER :: zet
1186 REAL(KIND=dp), DIMENSION(:, :, :), POINTER :: gcc
1187 TYPE(cp_parser_type) :: parser
1188
1189 line = ""
1190 line2 = ""
1191 symbol = ""
1192 symbol2 = ""
1193 bsname = ""
1194 bsname2 = ""
1195
1196 nbasis = 1
1197
1198 gto_basis_set%name = basis_set_name
1199 gto_basis_set%aliases = basis_set_name
1200 CALL section_vals_val_get(dft_section, "BASIS_SET_FILE_NAME", &
1201 n_rep_val=nbasis)
1202 ALLOCATE (cbasis(nbasis))
1203 DO ibasis = 1, nbasis
1204 CALL section_vals_val_get(dft_section, "BASIS_SET_FILE_NAME", &
1205 i_rep_val=ibasis, c_val=cbasis(ibasis))
1206 basis_set_file_name = cbasis(ibasis)
1207 tmp = basis_set_file_name
1208 CALL uppercase(tmp)
1209 IF (index(tmp, "MOLOPT") .NE. 0) CALL cite_reference(vandevondele2007)
1210 END DO
1211
1212 ! Search for the requested basis set in the basis set file
1213 ! until the basis set is found or the end of file is reached
1214
1215 basis_found = .false.
1216 basis_loop: DO ibasis = 1, nbasis
1217 IF (basis_found) EXIT basis_loop
1218 basis_set_file_name = cbasis(ibasis)
1219 CALL parser_create(parser, basis_set_file_name, para_env=para_env)
1220
1221 bsname = basis_set_name
1222 symbol = element_symbol
1223 irep = 0
1224
1225 tmp = basis_set_name
1226 CALL uppercase(tmp)
1227 IF (index(tmp, "MOLOPT") .NE. 0) CALL cite_reference(vandevondele2007)
1228
1229 nset = 0
1230 maxshell = 0
1231 maxpgf = 0
1232 maxco = 0
1233 ncgf = 0
1234 nsgf = 0
1235 gto_basis_set%nset = nset
1236 gto_basis_set%ncgf = ncgf
1237 gto_basis_set%nsgf = nsgf
1238 CALL reallocate(gto_basis_set%lmax, 1, nset)
1239 CALL reallocate(gto_basis_set%lmin, 1, nset)
1240 CALL reallocate(gto_basis_set%npgf, 1, nset)
1241 CALL reallocate(gto_basis_set%nshell, 1, nset)
1242 CALL reallocate(gto_basis_set%n, 1, maxshell, 1, nset)
1243 CALL reallocate(gto_basis_set%l, 1, maxshell, 1, nset)
1244 CALL reallocate(gto_basis_set%zet, 1, maxpgf, 1, nset)
1245 CALL reallocate(gto_basis_set%gcc, 1, maxpgf, 1, maxshell, 1, nset)
1246 CALL reallocate(gto_basis_set%set_radius, 1, nset)
1247 CALL reallocate(gto_basis_set%pgf_radius, 1, maxpgf, 1, nset)
1248 CALL reallocate(gto_basis_set%first_cgf, 1, maxshell, 1, nset)
1249 CALL reallocate(gto_basis_set%first_sgf, 1, maxshell, 1, nset)
1250 CALL reallocate(gto_basis_set%last_cgf, 1, maxshell, 1, nset)
1251 CALL reallocate(gto_basis_set%last_sgf, 1, maxshell, 1, nset)
1252 CALL reallocate(gto_basis_set%ncgf_set, 1, nset)
1253 CALL reallocate(gto_basis_set%nsgf_set, 1, nset)
1254 CALL reallocate(gto_basis_set%cphi, 1, maxco, 1, ncgf)
1255 CALL reallocate(gto_basis_set%sphi, 1, maxco, 1, nsgf)
1256 CALL reallocate(gto_basis_set%scon, 1, maxco, 1, nsgf)
1257 CALL reallocate(gto_basis_set%lx, 1, ncgf)
1258 CALL reallocate(gto_basis_set%ly, 1, ncgf)
1259 CALL reallocate(gto_basis_set%lz, 1, ncgf)
1260 CALL reallocate(gto_basis_set%m, 1, nsgf)
1261 CALL reallocate(gto_basis_set%norm_cgf, 1, ncgf)
1262
1263 IF (tmp .NE. "NONE") THEN
1264 search_loop: DO
1265
1266 CALL parser_search_string(parser, trim(bsname), .true., found, line)
1267 IF (found) THEN
1268 CALL uppercase(symbol)
1269 CALL uppercase(bsname)
1270
1271 match = .false.
1272 CALL uppercase(line)
1273 ! Check both the element symbol and the basis set name
1274 line2 = " "//line//" "
1275 symbol2 = " "//trim(symbol)//" "
1276 bsname2 = " "//trim(bsname)//" "
1277 strlen1 = len_trim(symbol2) + 1
1278 strlen2 = len_trim(bsname2) + 1
1279
1280 IF ((index(line2, symbol2(:strlen1)) > 0) .AND. &
1281 (index(line2, bsname2(:strlen2)) > 0)) match = .true.
1282 IF (match) THEN
1283 ! copy all names into aliases field
1284 i = index(line2, symbol2(:strlen1))
1285 i = i + 1 + index(line2(i + 1:), " ")
1286 gto_basis_set%aliases = line2(i:)
1287
1288 NULLIFY (gcc, l, lmax, lmin, n, npgf, nshell, zet)
1289 ! Read the basis set information
1290 CALL parser_get_object(parser, nset, newline=.true.)
1291
1292 CALL reallocate(npgf, 1, nset)
1293 CALL reallocate(nshell, 1, nset)
1294 CALL reallocate(lmax, 1, nset)
1295 CALL reallocate(lmin, 1, nset)
1296 CALL reallocate(n, 1, 1, 1, nset)
1297
1298 maxl = 0
1299 maxpgf = 0
1300 maxshell = 0
1301
1302 DO iset = 1, nset
1303 CALL parser_get_object(parser, n(1, iset), newline=.true.)
1304 CALL parser_get_object(parser, lmin(iset))
1305 CALL parser_get_object(parser, lmax(iset))
1306 CALL parser_get_object(parser, npgf(iset))
1307 maxl = max(maxl, lmax(iset))
1308 IF (npgf(iset) > maxpgf) THEN
1309 maxpgf = npgf(iset)
1310 CALL reallocate(zet, 1, maxpgf, 1, nset)
1311 CALL reallocate(gcc, 1, maxpgf, 1, maxshell, 1, nset)
1312 END IF
1313 nshell(iset) = 0
1314 DO lshell = lmin(iset), lmax(iset)
1315 nmin = n(1, iset) + lshell - lmin(iset)
1316 CALL parser_get_object(parser, ishell)
1317 nshell(iset) = nshell(iset) + ishell
1318 IF (nshell(iset) > maxshell) THEN
1319 maxshell = nshell(iset)
1320 CALL reallocate(n, 1, maxshell, 1, nset)
1321 CALL reallocate(l, 1, maxshell, 1, nset)
1322 CALL reallocate(gcc, 1, maxpgf, 1, maxshell, 1, nset)
1323 END IF
1324 DO i = 1, ishell
1325 n(nshell(iset) - ishell + i, iset) = nmin + i - 1
1326 l(nshell(iset) - ishell + i, iset) = lshell
1327 END DO
1328 END DO
1329 DO ipgf = 1, npgf(iset)
1330 CALL parser_get_object(parser, zet(ipgf, iset), newline=.true.)
1331 DO ishell = 1, nshell(iset)
1332 CALL parser_get_object(parser, gcc(ipgf, ishell, iset))
1333 END DO
1334 END DO
1335 END DO
1336
1337 ! Maximum angular momentum quantum number of the atomic kind
1338
1339 CALL init_orbital_pointers(maxl)
1340
1341 ! Allocate the global variables
1342
1343 gto_basis_set%nset = nset
1344 CALL reallocate(gto_basis_set%lmax, 1, nset)
1345 CALL reallocate(gto_basis_set%lmin, 1, nset)
1346 CALL reallocate(gto_basis_set%npgf, 1, nset)
1347 CALL reallocate(gto_basis_set%nshell, 1, nset)
1348 CALL reallocate(gto_basis_set%n, 1, maxshell, 1, nset)
1349 CALL reallocate(gto_basis_set%l, 1, maxshell, 1, nset)
1350 CALL reallocate(gto_basis_set%zet, 1, maxpgf, 1, nset)
1351 CALL reallocate(gto_basis_set%gcc, 1, maxpgf, 1, maxshell, 1, nset)
1352
1353 ! Copy the basis set information into the data structure
1354
1355 DO iset = 1, nset
1356 gto_basis_set%lmax(iset) = lmax(iset)
1357 gto_basis_set%lmin(iset) = lmin(iset)
1358 gto_basis_set%npgf(iset) = npgf(iset)
1359 gto_basis_set%nshell(iset) = nshell(iset)
1360 DO ishell = 1, nshell(iset)
1361 gto_basis_set%n(ishell, iset) = n(ishell, iset)
1362 gto_basis_set%l(ishell, iset) = l(ishell, iset)
1363 DO ipgf = 1, npgf(iset)
1364 gto_basis_set%gcc(ipgf, ishell, iset) = gcc(ipgf, ishell, iset)
1365 END DO
1366 END DO
1367 DO ipgf = 1, npgf(iset)
1368 gto_basis_set%zet(ipgf, iset) = zet(ipgf, iset)
1369 END DO
1370 END DO
1371
1372 ! Initialise the depending atomic kind information
1373
1374 CALL reallocate(gto_basis_set%set_radius, 1, nset)
1375 CALL reallocate(gto_basis_set%pgf_radius, 1, maxpgf, 1, nset)
1376 CALL reallocate(gto_basis_set%first_cgf, 1, maxshell, 1, nset)
1377 CALL reallocate(gto_basis_set%first_sgf, 1, maxshell, 1, nset)
1378 CALL reallocate(gto_basis_set%last_cgf, 1, maxshell, 1, nset)
1379 CALL reallocate(gto_basis_set%last_sgf, 1, maxshell, 1, nset)
1380 CALL reallocate(gto_basis_set%ncgf_set, 1, nset)
1381 CALL reallocate(gto_basis_set%nsgf_set, 1, nset)
1382
1383 maxco = 0
1384 ncgf = 0
1385 nsgf = 0
1386
1387 DO iset = 1, nset
1388 gto_basis_set%ncgf_set(iset) = 0
1389 gto_basis_set%nsgf_set(iset) = 0
1390 DO ishell = 1, nshell(iset)
1391 lshell = gto_basis_set%l(ishell, iset)
1392 gto_basis_set%first_cgf(ishell, iset) = ncgf + 1
1393 ncgf = ncgf + nco(lshell)
1394 gto_basis_set%last_cgf(ishell, iset) = ncgf
1395 gto_basis_set%ncgf_set(iset) = &
1396 gto_basis_set%ncgf_set(iset) + nco(lshell)
1397 gto_basis_set%first_sgf(ishell, iset) = nsgf + 1
1398 nsgf = nsgf + nso(lshell)
1399 gto_basis_set%last_sgf(ishell, iset) = nsgf
1400 gto_basis_set%nsgf_set(iset) = &
1401 gto_basis_set%nsgf_set(iset) + nso(lshell)
1402 END DO
1403 maxco = max(maxco, npgf(iset)*ncoset(lmax(iset)))
1404 END DO
1405
1406 gto_basis_set%ncgf = ncgf
1407 gto_basis_set%nsgf = nsgf
1408
1409 CALL reallocate(gto_basis_set%cphi, 1, maxco, 1, ncgf)
1410 CALL reallocate(gto_basis_set%sphi, 1, maxco, 1, nsgf)
1411 CALL reallocate(gto_basis_set%scon, 1, maxco, 1, nsgf)
1412 CALL reallocate(gto_basis_set%lx, 1, ncgf)
1413 CALL reallocate(gto_basis_set%ly, 1, ncgf)
1414 CALL reallocate(gto_basis_set%lz, 1, ncgf)
1415 CALL reallocate(gto_basis_set%m, 1, nsgf)
1416 CALL reallocate(gto_basis_set%norm_cgf, 1, ncgf)
1417 ALLOCATE (gto_basis_set%cgf_symbol(ncgf))
1418
1419 ALLOCATE (gto_basis_set%sgf_symbol(nsgf))
1420
1421 ncgf = 0
1422 nsgf = 0
1423
1424 DO iset = 1, nset
1425 DO ishell = 1, nshell(iset)
1426 lshell = gto_basis_set%l(ishell, iset)
1427 DO ico = ncoset(lshell - 1) + 1, ncoset(lshell)
1428 ncgf = ncgf + 1
1429 gto_basis_set%lx(ncgf) = indco(1, ico)
1430 gto_basis_set%ly(ncgf) = indco(2, ico)
1431 gto_basis_set%lz(ncgf) = indco(3, ico)
1432 gto_basis_set%cgf_symbol(ncgf) = &
1433 cgf_symbol(n(ishell, iset), (/gto_basis_set%lx(ncgf), &
1434 gto_basis_set%ly(ncgf), &
1435 gto_basis_set%lz(ncgf)/))
1436 END DO
1437 DO m = -lshell, lshell
1438 nsgf = nsgf + 1
1439 gto_basis_set%m(nsgf) = m
1440 gto_basis_set%sgf_symbol(nsgf) = &
1441 sgf_symbol(n(ishell, iset), lshell, m)
1442 END DO
1443 END DO
1444 END DO
1445
1446 DEALLOCATE (gcc, l, lmax, lmin, n, npgf, nshell, zet)
1447
1448 basis_found = .true.
1449 EXIT search_loop
1450 END IF
1451 ELSE
1452 EXIT search_loop
1453 END IF
1454 END DO search_loop
1455 ELSE
1456 match = .false.
1457 ALLOCATE (gto_basis_set%cgf_symbol(ncgf))
1458 ALLOCATE (gto_basis_set%sgf_symbol(nsgf))
1459 END IF
1460
1461 CALL parser_release(parser)
1462
1463 END DO basis_loop
1464
1465 IF (tmp .NE. "NONE") THEN
1466 IF (.NOT. basis_found) THEN
1467 basis_set_file_name = ""
1468 DO ibasis = 1, nbasis
1469 basis_set_file_name = trim(basis_set_file_name)//"<"//trim(cbasis(ibasis))//"> "
1470 END DO
1471 CALL cp_abort(__location__, &
1472 "The requested basis set <"//trim(bsname)// &
1473 "> for element <"//trim(symbol)//"> was not "// &
1474 "found in the basis set files "// &
1475 trim(basis_set_file_name))
1476 END IF
1477 END IF
1478 DEALLOCATE (cbasis)
1479
1480 CALL section_vals_val_get(dft_section, "SORT_BASIS", i_val=sort_method)
1481 CALL sort_gto_basis_set(gto_basis_set, sort_method)
1482
1483 END SUBROUTINE read_gto_basis_set1
1484
1485! **************************************************************************************************
1486!> \brief Read a Gaussian-type orbital (GTO) basis set from the database file.
1487!> \param element_symbol ...
1488!> \param basis_type ...
1489!> \param gto_basis_set ...
1490!> \param basis_section ...
1491!> \param irep ...
1492!> \param dft_section ...
1493!> \author MK
1494! **************************************************************************************************
1495 SUBROUTINE read_gto_basis_set2(element_symbol, basis_type, gto_basis_set, &
1496 basis_section, irep, dft_section)
1497
1498 CHARACTER(LEN=*), INTENT(IN) :: element_symbol
1499 CHARACTER(LEN=*), INTENT(INOUT) :: basis_type
1500 TYPE(gto_basis_set_type), INTENT(INOUT) :: gto_basis_set
1501 TYPE(section_vals_type), OPTIONAL, POINTER :: basis_section
1502 INTEGER :: irep
1503 TYPE(section_vals_type), OPTIONAL, POINTER :: dft_section
1504
1505 CHARACTER(len=20*default_string_length) :: line_att
1506 CHARACTER(LEN=240) :: line
1507 CHARACTER(LEN=242) :: line2
1508 CHARACTER(LEN=default_path_length) :: bsname, bsname2
1509 CHARACTER(LEN=LEN(element_symbol)) :: symbol
1510 CHARACTER(LEN=LEN(element_symbol)+2) :: symbol2
1511 INTEGER :: i, ico, ipgf, iset, ishell, lshell, m, &
1512 maxco, maxl, maxpgf, maxshell, ncgf, &
1513 nmin, nset, nsgf, sort_method
1514 INTEGER, DIMENSION(:), POINTER :: lmax, lmin, npgf, nshell
1515 INTEGER, DIMENSION(:, :), POINTER :: l, n
1516 LOGICAL :: is_ok
1517 REAL(KIND=dp), DIMENSION(:, :), POINTER :: zet
1518 REAL(KIND=dp), DIMENSION(:, :, :), POINTER :: gcc
1519 TYPE(cp_sll_val_type), POINTER :: list
1520 TYPE(val_type), POINTER :: val
1521
1522 line = ""
1523 line2 = ""
1524 symbol = ""
1525 symbol2 = ""
1526 bsname = ""
1527 bsname2 = ""
1528
1529 gto_basis_set%name = " "
1530 gto_basis_set%aliases = " "
1531
1532 bsname = " "
1533 symbol = element_symbol
1534
1535 nset = 0
1536 maxshell = 0
1537 maxpgf = 0
1538 maxco = 0
1539 ncgf = 0
1540 nsgf = 0
1541 gto_basis_set%nset = nset
1542 gto_basis_set%ncgf = ncgf
1543 gto_basis_set%nsgf = nsgf
1544 CALL reallocate(gto_basis_set%lmax, 1, nset)
1545 CALL reallocate(gto_basis_set%lmin, 1, nset)
1546 CALL reallocate(gto_basis_set%npgf, 1, nset)
1547 CALL reallocate(gto_basis_set%nshell, 1, nset)
1548 CALL reallocate(gto_basis_set%n, 1, maxshell, 1, nset)
1549 CALL reallocate(gto_basis_set%l, 1, maxshell, 1, nset)
1550 CALL reallocate(gto_basis_set%zet, 1, maxpgf, 1, nset)
1551 CALL reallocate(gto_basis_set%gcc, 1, maxpgf, 1, maxshell, 1, nset)
1552 CALL reallocate(gto_basis_set%set_radius, 1, nset)
1553 CALL reallocate(gto_basis_set%pgf_radius, 1, maxpgf, 1, nset)
1554 CALL reallocate(gto_basis_set%first_cgf, 1, maxshell, 1, nset)
1555 CALL reallocate(gto_basis_set%first_sgf, 1, maxshell, 1, nset)
1556 CALL reallocate(gto_basis_set%last_cgf, 1, maxshell, 1, nset)
1557 CALL reallocate(gto_basis_set%last_sgf, 1, maxshell, 1, nset)
1558 CALL reallocate(gto_basis_set%ncgf_set, 1, nset)
1559 CALL reallocate(gto_basis_set%nsgf_set, 1, nset)
1560 CALL reallocate(gto_basis_set%cphi, 1, maxco, 1, ncgf)
1561 CALL reallocate(gto_basis_set%sphi, 1, maxco, 1, nsgf)
1562 CALL reallocate(gto_basis_set%scon, 1, maxco, 1, nsgf)
1563 CALL reallocate(gto_basis_set%lx, 1, ncgf)
1564 CALL reallocate(gto_basis_set%ly, 1, ncgf)
1565 CALL reallocate(gto_basis_set%lz, 1, ncgf)
1566 CALL reallocate(gto_basis_set%m, 1, nsgf)
1567 CALL reallocate(gto_basis_set%norm_cgf, 1, ncgf)
1568
1569 basis_type = ""
1570 CALL section_vals_val_get(basis_section, "_SECTION_PARAMETERS_", i_rep_section=irep, c_val=basis_type)
1571 IF (basis_type == "Orbital") basis_type = "ORB"
1572
1573 NULLIFY (list, val)
1574 CALL section_vals_list_get(basis_section, "_DEFAULT_KEYWORD_", i_rep_section=irep, list=list)
1575 CALL uppercase(symbol)
1576 CALL uppercase(bsname)
1577
1578 NULLIFY (gcc, l, lmax, lmin, n, npgf, nshell, zet)
1579 ! Read the basis set information
1580 is_ok = cp_sll_val_next(list, val)
1581 IF (.NOT. is_ok) cpabort("Error reading the Basis set from input file!")
1582 CALL val_get(val, c_val=line_att)
1583 READ (line_att, *) nset
1584
1585 CALL reallocate(npgf, 1, nset)
1586 CALL reallocate(nshell, 1, nset)
1587 CALL reallocate(lmax, 1, nset)
1588 CALL reallocate(lmin, 1, nset)
1589 CALL reallocate(n, 1, 1, 1, nset)
1590
1591 maxl = 0
1592 maxpgf = 0
1593 maxshell = 0
1594
1595 DO iset = 1, nset
1596 is_ok = cp_sll_val_next(list, val)
1597 IF (.NOT. is_ok) cpabort("Error reading the Basis set from input file!")
1598 CALL val_get(val, c_val=line_att)
1599 READ (line_att, *) n(1, iset)
1600 CALL remove_word(line_att)
1601 READ (line_att, *) lmin(iset)
1602 CALL remove_word(line_att)
1603 READ (line_att, *) lmax(iset)
1604 CALL remove_word(line_att)
1605 READ (line_att, *) npgf(iset)
1606 CALL remove_word(line_att)
1607 maxl = max(maxl, lmax(iset))
1608 IF (npgf(iset) > maxpgf) THEN
1609 maxpgf = npgf(iset)
1610 CALL reallocate(zet, 1, maxpgf, 1, nset)
1611 CALL reallocate(gcc, 1, maxpgf, 1, maxshell, 1, nset)
1612 END IF
1613 nshell(iset) = 0
1614 DO lshell = lmin(iset), lmax(iset)
1615 nmin = n(1, iset) + lshell - lmin(iset)
1616 READ (line_att, *) ishell
1617 CALL remove_word(line_att)
1618 nshell(iset) = nshell(iset) + ishell
1619 IF (nshell(iset) > maxshell) THEN
1620 maxshell = nshell(iset)
1621 CALL reallocate(n, 1, maxshell, 1, nset)
1622 CALL reallocate(l, 1, maxshell, 1, nset)
1623 CALL reallocate(gcc, 1, maxpgf, 1, maxshell, 1, nset)
1624 END IF
1625 DO i = 1, ishell
1626 n(nshell(iset) - ishell + i, iset) = nmin + i - 1
1627 l(nshell(iset) - ishell + i, iset) = lshell
1628 END DO
1629 END DO
1630 IF (len_trim(line_att) /= 0) &
1631 cpabort("Error reading the Basis from input file!")
1632 DO ipgf = 1, npgf(iset)
1633 is_ok = cp_sll_val_next(list, val)
1634 IF (.NOT. is_ok) cpabort("Error reading the Basis set from input file!")
1635 CALL val_get(val, c_val=line_att)
1636 READ (line_att, *) zet(ipgf, iset), (gcc(ipgf, ishell, iset), ishell=1, nshell(iset))
1637 END DO
1638 END DO
1639
1640 ! Maximum angular momentum quantum number of the atomic kind
1641
1642 CALL init_orbital_pointers(maxl)
1643
1644 ! Allocate the global variables
1645
1646 gto_basis_set%nset = nset
1647 CALL reallocate(gto_basis_set%lmax, 1, nset)
1648 CALL reallocate(gto_basis_set%lmin, 1, nset)
1649 CALL reallocate(gto_basis_set%npgf, 1, nset)
1650 CALL reallocate(gto_basis_set%nshell, 1, nset)
1651 CALL reallocate(gto_basis_set%n, 1, maxshell, 1, nset)
1652 CALL reallocate(gto_basis_set%l, 1, maxshell, 1, nset)
1653 CALL reallocate(gto_basis_set%zet, 1, maxpgf, 1, nset)
1654 CALL reallocate(gto_basis_set%gcc, 1, maxpgf, 1, maxshell, 1, nset)
1655
1656 ! Copy the basis set information into the data structure
1657
1658 DO iset = 1, nset
1659 gto_basis_set%lmax(iset) = lmax(iset)
1660 gto_basis_set%lmin(iset) = lmin(iset)
1661 gto_basis_set%npgf(iset) = npgf(iset)
1662 gto_basis_set%nshell(iset) = nshell(iset)
1663 DO ishell = 1, nshell(iset)
1664 gto_basis_set%n(ishell, iset) = n(ishell, iset)
1665 gto_basis_set%l(ishell, iset) = l(ishell, iset)
1666 DO ipgf = 1, npgf(iset)
1667 gto_basis_set%gcc(ipgf, ishell, iset) = gcc(ipgf, ishell, iset)
1668 END DO
1669 END DO
1670 DO ipgf = 1, npgf(iset)
1671 gto_basis_set%zet(ipgf, iset) = zet(ipgf, iset)
1672 END DO
1673 END DO
1674
1675 ! Initialise the depending atomic kind information
1676
1677 CALL reallocate(gto_basis_set%set_radius, 1, nset)
1678 CALL reallocate(gto_basis_set%pgf_radius, 1, maxpgf, 1, nset)
1679 CALL reallocate(gto_basis_set%first_cgf, 1, maxshell, 1, nset)
1680 CALL reallocate(gto_basis_set%first_sgf, 1, maxshell, 1, nset)
1681 CALL reallocate(gto_basis_set%last_cgf, 1, maxshell, 1, nset)
1682 CALL reallocate(gto_basis_set%last_sgf, 1, maxshell, 1, nset)
1683 CALL reallocate(gto_basis_set%ncgf_set, 1, nset)
1684 CALL reallocate(gto_basis_set%nsgf_set, 1, nset)
1685
1686 maxco = 0
1687 ncgf = 0
1688 nsgf = 0
1689
1690 DO iset = 1, nset
1691 gto_basis_set%ncgf_set(iset) = 0
1692 gto_basis_set%nsgf_set(iset) = 0
1693 DO ishell = 1, nshell(iset)
1694 lshell = gto_basis_set%l(ishell, iset)
1695 gto_basis_set%first_cgf(ishell, iset) = ncgf + 1
1696 ncgf = ncgf + nco(lshell)
1697 gto_basis_set%last_cgf(ishell, iset) = ncgf
1698 gto_basis_set%ncgf_set(iset) = &
1699 gto_basis_set%ncgf_set(iset) + nco(lshell)
1700 gto_basis_set%first_sgf(ishell, iset) = nsgf + 1
1701 nsgf = nsgf + nso(lshell)
1702 gto_basis_set%last_sgf(ishell, iset) = nsgf
1703 gto_basis_set%nsgf_set(iset) = &
1704 gto_basis_set%nsgf_set(iset) + nso(lshell)
1705 END DO
1706 maxco = max(maxco, npgf(iset)*ncoset(lmax(iset)))
1707 END DO
1708
1709 gto_basis_set%ncgf = ncgf
1710 gto_basis_set%nsgf = nsgf
1711
1712 CALL reallocate(gto_basis_set%cphi, 1, maxco, 1, ncgf)
1713 CALL reallocate(gto_basis_set%sphi, 1, maxco, 1, nsgf)
1714 CALL reallocate(gto_basis_set%scon, 1, maxco, 1, nsgf)
1715 CALL reallocate(gto_basis_set%lx, 1, ncgf)
1716 CALL reallocate(gto_basis_set%ly, 1, ncgf)
1717 CALL reallocate(gto_basis_set%lz, 1, ncgf)
1718 CALL reallocate(gto_basis_set%m, 1, nsgf)
1719 CALL reallocate(gto_basis_set%norm_cgf, 1, ncgf)
1720 ALLOCATE (gto_basis_set%cgf_symbol(ncgf))
1721
1722 ALLOCATE (gto_basis_set%sgf_symbol(nsgf))
1723
1724 ncgf = 0
1725 nsgf = 0
1726
1727 DO iset = 1, nset
1728 DO ishell = 1, nshell(iset)
1729 lshell = gto_basis_set%l(ishell, iset)
1730 DO ico = ncoset(lshell - 1) + 1, ncoset(lshell)
1731 ncgf = ncgf + 1
1732 gto_basis_set%lx(ncgf) = indco(1, ico)
1733 gto_basis_set%ly(ncgf) = indco(2, ico)
1734 gto_basis_set%lz(ncgf) = indco(3, ico)
1735 gto_basis_set%cgf_symbol(ncgf) = &
1736 cgf_symbol(n(ishell, iset), (/gto_basis_set%lx(ncgf), &
1737 gto_basis_set%ly(ncgf), &
1738 gto_basis_set%lz(ncgf)/))
1739 END DO
1740 DO m = -lshell, lshell
1741 nsgf = nsgf + 1
1742 gto_basis_set%m(nsgf) = m
1743 gto_basis_set%sgf_symbol(nsgf) = &
1744 sgf_symbol(n(ishell, iset), lshell, m)
1745 END DO
1746 END DO
1747 END DO
1748
1749 DEALLOCATE (gcc, l, lmax, lmin, n, npgf, nshell, zet)
1750
1751 IF (PRESENT(dft_section)) THEN
1752 CALL section_vals_val_get(dft_section, "SORT_BASIS", i_val=sort_method)
1753 CALL sort_gto_basis_set(gto_basis_set, sort_method)
1754 END IF
1755
1756 END SUBROUTINE read_gto_basis_set2
1757
1758! **************************************************************************************************
1759!> \brief Set the components of Gaussian-type orbital (GTO) basis set data set.
1760!> \param gto_basis_set ...
1761!> \param name ...
1762!> \param aliases ...
1763!> \param norm_type ...
1764!> \param kind_radius ...
1765!> \param ncgf ...
1766!> \param nset ...
1767!> \param nsgf ...
1768!> \param cgf_symbol ...
1769!> \param sgf_symbol ...
1770!> \param norm_cgf ...
1771!> \param set_radius ...
1772!> \param lmax ...
1773!> \param lmin ...
1774!> \param lx ...
1775!> \param ly ...
1776!> \param lz ...
1777!> \param m ...
1778!> \param ncgf_set ...
1779!> \param npgf ...
1780!> \param nsgf_set ...
1781!> \param nshell ...
1782!> \param cphi ...
1783!> \param pgf_radius ...
1784!> \param sphi ...
1785!> \param scon ...
1786!> \param zet ...
1787!> \param first_cgf ...
1788!> \param first_sgf ...
1789!> \param l ...
1790!> \param last_cgf ...
1791!> \param last_sgf ...
1792!> \param n ...
1793!> \param gcc ...
1794!> \param short_kind_radius ...
1795!> \author MK
1796! **************************************************************************************************
1797 SUBROUTINE set_gto_basis_set(gto_basis_set, name, aliases, norm_type, kind_radius, ncgf, &
1798 nset, nsgf, cgf_symbol, sgf_symbol, norm_cgf, set_radius, lmax, &
1799 lmin, lx, ly, lz, m, ncgf_set, npgf, nsgf_set, nshell, &
1800 cphi, pgf_radius, sphi, scon, zet, first_cgf, first_sgf, l, &
1801 last_cgf, last_sgf, n, gcc, short_kind_radius)
1802
1803 TYPE(gto_basis_set_type), INTENT(INOUT) :: gto_basis_set
1804 CHARACTER(LEN=default_string_length), INTENT(IN), &
1805 OPTIONAL :: name, aliases
1806 INTEGER, INTENT(IN), OPTIONAL :: norm_type
1807 REAL(kind=dp), INTENT(IN), OPTIONAL :: kind_radius
1808 INTEGER, INTENT(IN), OPTIONAL :: ncgf, nset, nsgf
1809 CHARACTER(LEN=12), DIMENSION(:), OPTIONAL, POINTER :: cgf_symbol
1810 CHARACTER(LEN=6), DIMENSION(:), OPTIONAL, POINTER :: sgf_symbol
1811 REAL(kind=dp), DIMENSION(:), OPTIONAL, POINTER :: norm_cgf, set_radius
1812 INTEGER, DIMENSION(:), OPTIONAL, POINTER :: lmax, lmin, lx, ly, lz, m, ncgf_set, &
1813 npgf, nsgf_set, nshell
1814 REAL(kind=dp), DIMENSION(:, :), OPTIONAL, POINTER :: cphi, pgf_radius, sphi, scon, zet
1815 INTEGER, DIMENSION(:, :), OPTIONAL, POINTER :: first_cgf, first_sgf, l, last_cgf, &
1816 last_sgf, n
1817 REAL(kind=dp), DIMENSION(:, :, :), OPTIONAL, &
1818 POINTER :: gcc
1819 REAL(kind=dp), INTENT(IN), OPTIONAL :: short_kind_radius
1820
1821 IF (PRESENT(name)) gto_basis_set%name = name
1822 IF (PRESENT(aliases)) gto_basis_set%aliases = aliases
1823 IF (PRESENT(norm_type)) gto_basis_set%norm_type = norm_type
1824 IF (PRESENT(kind_radius)) gto_basis_set%kind_radius = kind_radius
1825 IF (PRESENT(short_kind_radius)) gto_basis_set%short_kind_radius = short_kind_radius
1826 IF (PRESENT(ncgf)) gto_basis_set%ncgf = ncgf
1827 IF (PRESENT(nset)) gto_basis_set%nset = nset
1828 IF (PRESENT(nsgf)) gto_basis_set%nsgf = nsgf
1829 IF (PRESENT(cgf_symbol)) gto_basis_set%cgf_symbol(:) = cgf_symbol(:)
1830 IF (PRESENT(sgf_symbol)) gto_basis_set%sgf_symbol(:) = sgf_symbol(:)
1831 IF (PRESENT(norm_cgf)) gto_basis_set%norm_cgf(:) = norm_cgf(:)
1832 IF (PRESENT(set_radius)) gto_basis_set%set_radius(:) = set_radius(:)
1833 IF (PRESENT(lmax)) gto_basis_set%lmax(:) = lmax(:)
1834 IF (PRESENT(lmin)) gto_basis_set%lmin(:) = lmin(:)
1835 IF (PRESENT(lx)) gto_basis_set%lx(:) = lx(:)
1836 IF (PRESENT(ly)) gto_basis_set%ly(:) = ly(:)
1837 IF (PRESENT(lz)) gto_basis_set%lz(:) = lz(:)
1838 IF (PRESENT(m)) gto_basis_set%m(:) = m(:)
1839 IF (PRESENT(ncgf_set)) gto_basis_set%ncgf_set(:) = ncgf_set(:)
1840 IF (PRESENT(npgf)) gto_basis_set%npgf(:) = npgf(:)
1841 IF (PRESENT(nsgf_set)) gto_basis_set%nsgf_set(:) = nsgf_set(:)
1842 IF (PRESENT(nshell)) gto_basis_set%nshell(:) = nshell(:)
1843 IF (PRESENT(cphi)) gto_basis_set%cphi(:, :) = cphi(:, :)
1844 IF (PRESENT(pgf_radius)) gto_basis_set%pgf_radius(:, :) = pgf_radius(:, :)
1845 IF (PRESENT(sphi)) gto_basis_set%sphi(:, :) = sphi(:, :)
1846 IF (PRESENT(scon)) gto_basis_set%scon(:, :) = scon(:, :)
1847 IF (PRESENT(zet)) gto_basis_set%zet(:, :) = zet(:, :)
1848 IF (PRESENT(first_cgf)) gto_basis_set%first_cgf(:, :) = first_cgf(:, :)
1849 IF (PRESENT(first_sgf)) gto_basis_set%first_sgf(:, :) = first_sgf(:, :)
1850 IF (PRESENT(l)) l(:, :) = gto_basis_set%l(:, :)
1851 IF (PRESENT(last_cgf)) gto_basis_set%last_cgf(:, :) = last_cgf(:, :)
1852 IF (PRESENT(last_sgf)) gto_basis_set%last_sgf(:, :) = last_sgf(:, :)
1853 IF (PRESENT(n)) gto_basis_set%n(:, :) = n(:, :)
1854 IF (PRESENT(gcc)) gto_basis_set%gcc(:, :, :) = gcc(:, :, :)
1855
1856 END SUBROUTINE set_gto_basis_set
1857
1858! **************************************************************************************************
1859!> \brief Write a Gaussian-type orbital (GTO) basis set data set to the output unit.
1860!> \param gto_basis_set ...
1861!> \param output_unit ...
1862!> \param header ...
1863!> \author MK
1864! **************************************************************************************************
1865 SUBROUTINE write_gto_basis_set(gto_basis_set, output_unit, header)
1866
1867 TYPE(gto_basis_set_type), INTENT(IN) :: gto_basis_set
1868 INTEGER, INTENT(in) :: output_unit
1869 CHARACTER(len=*), OPTIONAL :: header
1870
1871 INTEGER :: ipgf, iset, ishell
1872
1873 IF (output_unit > 0) THEN
1874
1875 IF (PRESENT(header)) THEN
1876 WRITE (unit=output_unit, fmt="(/,T6,A,T41,A40)") &
1877 trim(header), trim(gto_basis_set%name)
1878 END IF
1879
1880 WRITE (unit=output_unit, fmt="(/,(T8,A,T71,I10))") &
1881 "Number of orbital shell sets: ", &
1882 gto_basis_set%nset, &
1883 "Number of orbital shells: ", &
1884 sum(gto_basis_set%nshell(:)), &
1885 "Number of primitive Cartesian functions: ", &
1886 sum(gto_basis_set%npgf(:)), &
1887 "Number of Cartesian basis functions: ", &
1888 gto_basis_set%ncgf, &
1889 "Number of spherical basis functions: ", &
1890 gto_basis_set%nsgf, &
1891 "Norm type: ", &
1892 gto_basis_set%norm_type
1893
1894 WRITE (unit=output_unit, fmt="(/,T6,A,T41,A40,/,/,T25,A)") &
1895 "GTO basis set information for", trim(gto_basis_set%name), &
1896 "Set Shell n l Exponent Coefficient"
1897
1898 DO iset = 1, gto_basis_set%nset
1899 WRITE (unit=output_unit, fmt="(A)") ""
1900 DO ishell = 1, gto_basis_set%nshell(iset)
1901 WRITE (unit=output_unit, &
1902 fmt="(T25,I3,4X,I4,4X,I2,2X,I2,(T51,2F15.6))") &
1903 iset, ishell, &
1904 gto_basis_set%n(ishell, iset), &
1905 gto_basis_set%l(ishell, iset), &
1906 (gto_basis_set%zet(ipgf, iset), &
1907 gto_basis_set%gcc(ipgf, ishell, iset), &
1908 ipgf=1, gto_basis_set%npgf(iset))
1909 END DO
1910 END DO
1911
1912 END IF
1913
1914 END SUBROUTINE write_gto_basis_set
1915
1916! **************************************************************************************************
1917
1918! **************************************************************************************************
1919!> \brief Write a Gaussian-type orbital (GTO) basis set data set to the output unit.
1920!> \param orb_basis_set ...
1921!> \param output_unit ...
1922!> \param header ...
1923!> \author MK
1924! **************************************************************************************************
1925 SUBROUTINE write_orb_basis_set(orb_basis_set, output_unit, header)
1926
1927 TYPE(gto_basis_set_type), INTENT(IN) :: orb_basis_set
1928 INTEGER, INTENT(in) :: output_unit
1929 CHARACTER(len=*), OPTIONAL :: header
1930
1931 INTEGER :: icgf, ico, ipgf, iset, ishell
1932
1933 IF (output_unit > 0) THEN
1934 IF (PRESENT(header)) THEN
1935 WRITE (unit=output_unit, fmt="(/,T6,A,T41,A40)") &
1936 trim(header), trim(orb_basis_set%name)
1937 END IF
1938
1939 WRITE (unit=output_unit, fmt="(/,(T8,A,T71,I10))") &
1940 "Number of orbital shell sets: ", &
1941 orb_basis_set%nset, &
1942 "Number of orbital shells: ", &
1943 sum(orb_basis_set%nshell(:)), &
1944 "Number of primitive Cartesian functions: ", &
1945 sum(orb_basis_set%npgf(:)), &
1946 "Number of Cartesian basis functions: ", &
1947 orb_basis_set%ncgf, &
1948 "Number of spherical basis functions: ", &
1949 orb_basis_set%nsgf, &
1950 "Norm type: ", &
1951 orb_basis_set%norm_type
1952
1953 WRITE (unit=output_unit, fmt="(/,T8,A,/,/,T25,A)") &
1954 "Normalised Cartesian orbitals:", &
1955 "Set Shell Orbital Exponent Coefficient"
1956
1957 icgf = 0
1958
1959 DO iset = 1, orb_basis_set%nset
1960 DO ishell = 1, orb_basis_set%nshell(iset)
1961 WRITE (unit=output_unit, fmt="(A)") ""
1962 DO ico = 1, nco(orb_basis_set%l(ishell, iset))
1963 icgf = icgf + 1
1964 WRITE (unit=output_unit, &
1965 fmt="(T25,I3,4X,I4,3X,A12,(T51,2F15.6))") &
1966 iset, ishell, orb_basis_set%cgf_symbol(icgf), &
1967 (orb_basis_set%zet(ipgf, iset), &
1968 orb_basis_set%norm_cgf(icgf)* &
1969 orb_basis_set%gcc(ipgf, ishell, iset), &
1970 ipgf=1, orb_basis_set%npgf(iset))
1971 END DO
1972 END DO
1973 END DO
1974 END IF
1975
1976 END SUBROUTINE write_orb_basis_set
1977
1978! **************************************************************************************************
1979!> \brief ...
1980!> \param sto_basis_set ...
1981! **************************************************************************************************
1982 SUBROUTINE allocate_sto_basis_set(sto_basis_set)
1983
1984 TYPE(sto_basis_set_type), POINTER :: sto_basis_set
1985
1986 CALL deallocate_sto_basis_set(sto_basis_set)
1987
1988 ALLOCATE (sto_basis_set)
1989
1990 END SUBROUTINE allocate_sto_basis_set
1991
1992! **************************************************************************************************
1993!> \brief ...
1994!> \param sto_basis_set ...
1995! **************************************************************************************************
1996 SUBROUTINE deallocate_sto_basis_set(sto_basis_set)
1997
1998 TYPE(sto_basis_set_type), POINTER :: sto_basis_set
1999
2000! -------------------------------------------------------------------------
2001
2002 IF (ASSOCIATED(sto_basis_set)) THEN
2003 IF (ASSOCIATED(sto_basis_set%symbol)) THEN
2004 DEALLOCATE (sto_basis_set%symbol)
2005 END IF
2006 IF (ASSOCIATED(sto_basis_set%nq)) THEN
2007 DEALLOCATE (sto_basis_set%nq)
2008 END IF
2009 IF (ASSOCIATED(sto_basis_set%lq)) THEN
2010 DEALLOCATE (sto_basis_set%lq)
2011 END IF
2012 IF (ASSOCIATED(sto_basis_set%zet)) THEN
2013 DEALLOCATE (sto_basis_set%zet)
2014 END IF
2015
2016 DEALLOCATE (sto_basis_set)
2017 END IF
2018 END SUBROUTINE deallocate_sto_basis_set
2019
2020! **************************************************************************************************
2021!> \brief ...
2022!> \param sto_basis_set ...
2023!> \param name ...
2024!> \param nshell ...
2025!> \param symbol ...
2026!> \param nq ...
2027!> \param lq ...
2028!> \param zet ...
2029!> \param maxlq ...
2030!> \param numsto ...
2031! **************************************************************************************************
2032 SUBROUTINE get_sto_basis_set(sto_basis_set, name, nshell, symbol, nq, lq, zet, maxlq, numsto)
2033
2034 TYPE(sto_basis_set_type), INTENT(IN) :: sto_basis_set
2035 CHARACTER(LEN=default_string_length), &
2036 INTENT(OUT), OPTIONAL :: name
2037 INTEGER, INTENT(OUT), OPTIONAL :: nshell
2038 CHARACTER(LEN=6), DIMENSION(:), OPTIONAL, POINTER :: symbol
2039 INTEGER, DIMENSION(:), OPTIONAL, POINTER :: nq, lq
2040 REAL(kind=dp), DIMENSION(:), OPTIONAL, POINTER :: zet
2041 INTEGER, INTENT(OUT), OPTIONAL :: maxlq, numsto
2042
2043 INTEGER :: iset
2044
2045 IF (PRESENT(name)) name = sto_basis_set%name
2046 IF (PRESENT(nshell)) nshell = sto_basis_set%nshell
2047 IF (PRESENT(symbol)) symbol => sto_basis_set%symbol
2048 IF (PRESENT(nq)) nq => sto_basis_set%nq
2049 IF (PRESENT(lq)) lq => sto_basis_set%lq
2050 IF (PRESENT(zet)) zet => sto_basis_set%zet
2051 IF (PRESENT(maxlq)) THEN
2052 maxlq = maxval(sto_basis_set%lq(1:sto_basis_set%nshell))
2053 END IF
2054 IF (PRESENT(numsto)) THEN
2055 numsto = 0
2056 DO iset = 1, sto_basis_set%nshell
2057 numsto = numsto + 2*sto_basis_set%lq(iset) + 1
2058 END DO
2059 END IF
2060
2061 END SUBROUTINE get_sto_basis_set
2062
2063! **************************************************************************************************
2064!> \brief ...
2065!> \param sto_basis_set ...
2066!> \param name ...
2067!> \param nshell ...
2068!> \param symbol ...
2069!> \param nq ...
2070!> \param lq ...
2071!> \param zet ...
2072! **************************************************************************************************
2073 SUBROUTINE set_sto_basis_set(sto_basis_set, name, nshell, symbol, nq, lq, zet)
2074
2075 TYPE(sto_basis_set_type), INTENT(INOUT) :: sto_basis_set
2076 CHARACTER(LEN=default_string_length), INTENT(IN), &
2077 OPTIONAL :: name
2078 INTEGER, INTENT(IN), OPTIONAL :: nshell
2079 CHARACTER(LEN=6), DIMENSION(:), OPTIONAL, POINTER :: symbol
2080 INTEGER, DIMENSION(:), OPTIONAL, POINTER :: nq, lq
2081 REAL(kind=dp), DIMENSION(:), OPTIONAL, POINTER :: zet
2082
2083 INTEGER :: ns
2084
2085 IF (PRESENT(name)) sto_basis_set%name = name
2086 IF (PRESENT(nshell)) sto_basis_set%nshell = nshell
2087 IF (PRESENT(symbol)) THEN
2088 ns = SIZE(symbol)
2089 IF (ASSOCIATED(sto_basis_set%symbol)) DEALLOCATE (sto_basis_set%symbol)
2090 ALLOCATE (sto_basis_set%symbol(1:ns))
2091 sto_basis_set%symbol(:) = symbol(:)
2092 END IF
2093 IF (PRESENT(nq)) THEN
2094 ns = SIZE(nq)
2095 CALL reallocate(sto_basis_set%nq, 1, ns)
2096 sto_basis_set%nq = nq(:)
2097 END IF
2098 IF (PRESENT(lq)) THEN
2099 ns = SIZE(lq)
2100 CALL reallocate(sto_basis_set%lq, 1, ns)
2101 sto_basis_set%lq = lq(:)
2102 END IF
2103 IF (PRESENT(zet)) THEN
2104 ns = SIZE(zet)
2105 CALL reallocate(sto_basis_set%zet, 1, ns)
2106 sto_basis_set%zet = zet(:)
2107 END IF
2108
2109 END SUBROUTINE set_sto_basis_set
2110
2111! **************************************************************************************************
2112!> \brief ...
2113!> \param element_symbol ...
2114!> \param basis_set_name ...
2115!> \param sto_basis_set ...
2116!> \param para_env ...
2117!> \param dft_section ...
2118! **************************************************************************************************
2119 SUBROUTINE read_sto_basis_set(element_symbol, basis_set_name, sto_basis_set, para_env, dft_section)
2120
2121 ! Read a Slater-type orbital (STO) basis set from the database file.
2122
2123 CHARACTER(LEN=*), INTENT(IN) :: element_symbol, basis_set_name
2124 TYPE(sto_basis_set_type), INTENT(INOUT) :: sto_basis_set
2125 TYPE(mp_para_env_type), POINTER :: para_env
2126 TYPE(section_vals_type), POINTER :: dft_section
2127
2128 CHARACTER(LEN=10) :: nlsym
2129 CHARACTER(LEN=2) :: lsym
2130 CHARACTER(LEN=240) :: line
2131 CHARACTER(LEN=242) :: line2
2132 CHARACTER(len=default_path_length) :: basis_set_file_name, tmp
2133 CHARACTER(LEN=default_path_length), DIMENSION(:), &
2134 POINTER :: cbasis
2135 CHARACTER(LEN=LEN(basis_set_name)) :: bsname
2136 CHARACTER(LEN=LEN(basis_set_name)+2) :: bsname2
2137 CHARACTER(LEN=LEN(element_symbol)) :: symbol
2138 CHARACTER(LEN=LEN(element_symbol)+2) :: symbol2
2139 INTEGER :: ibasis, irep, iset, nbasis, nq, nset, &
2140 strlen1, strlen2
2141 LOGICAL :: basis_found, found, match
2142 REAL(kind=dp) :: zet
2143 TYPE(cp_parser_type) :: parser
2144
2145 line = ""
2146 line2 = ""
2147 symbol = ""
2148 symbol2 = ""
2149 bsname = ""
2150 bsname2 = ""
2151
2152 nbasis = 1
2153
2154 sto_basis_set%name = basis_set_name
2155 CALL section_vals_val_get(dft_section, "BASIS_SET_FILE_NAME", &
2156 n_rep_val=nbasis)
2157 ALLOCATE (cbasis(nbasis))
2158 DO ibasis = 1, nbasis
2159 CALL section_vals_val_get(dft_section, "BASIS_SET_FILE_NAME", &
2160 i_rep_val=ibasis, c_val=cbasis(ibasis))
2161 basis_set_file_name = cbasis(ibasis)
2162 tmp = basis_set_file_name
2163 CALL uppercase(tmp)
2164 END DO
2165
2166 ! Search for the requested basis set in the basis set file
2167 ! until the basis set is found or the end of file is reached
2168
2169 basis_found = .false.
2170 basis_loop: DO ibasis = 1, nbasis
2171 IF (basis_found) EXIT basis_loop
2172 basis_set_file_name = cbasis(ibasis)
2173 CALL parser_create(parser, basis_set_file_name, para_env=para_env)
2174
2175 bsname = basis_set_name
2176 symbol = element_symbol
2177 irep = 0
2178
2179 tmp = basis_set_name
2180 CALL uppercase(tmp)
2181
2182 IF (tmp .NE. "NONE") THEN
2183 search_loop: DO
2184
2185 CALL parser_search_string(parser, trim(bsname), .true., found, line)
2186 IF (found) THEN
2187 CALL uppercase(symbol)
2188 CALL uppercase(bsname)
2189
2190 match = .false.
2191 CALL uppercase(line)
2192 ! Check both the element symbol and the basis set name
2193 line2 = " "//line//" "
2194 symbol2 = " "//trim(symbol)//" "
2195 bsname2 = " "//trim(bsname)//" "
2196 strlen1 = len_trim(symbol2) + 1
2197 strlen2 = len_trim(bsname2) + 1
2198
2199 IF ((index(line2, symbol2(:strlen1)) > 0) .AND. &
2200 (index(line2, bsname2(:strlen2)) > 0)) match = .true.
2201 IF (match) THEN
2202 ! Read the basis set information
2203 CALL parser_get_object(parser, nset, newline=.true.)
2204 sto_basis_set%nshell = nset
2205
2206 CALL reallocate(sto_basis_set%nq, 1, nset)
2207 CALL reallocate(sto_basis_set%lq, 1, nset)
2208 CALL reallocate(sto_basis_set%zet, 1, nset)
2209 ALLOCATE (sto_basis_set%symbol(nset))
2210
2211 DO iset = 1, nset
2212 CALL parser_get_object(parser, nq, newline=.true.)
2213 CALL parser_get_object(parser, lsym)
2214 CALL parser_get_object(parser, zet)
2215 sto_basis_set%nq(iset) = nq
2216 sto_basis_set%zet(iset) = zet
2217 WRITE (nlsym, "(I2,A)") nq, trim(lsym)
2218 sto_basis_set%symbol(iset) = trim(nlsym)
2219 SELECT CASE (trim(lsym))
2220 CASE ("S", "s")
2221 sto_basis_set%lq(iset) = 0
2222 CASE ("P", "p")
2223 sto_basis_set%lq(iset) = 1
2224 CASE ("D", "d")
2225 sto_basis_set%lq(iset) = 2
2226 CASE ("F", "f")
2227 sto_basis_set%lq(iset) = 3
2228 CASE ("G", "g")
2229 sto_basis_set%lq(iset) = 4
2230 CASE ("H", "h")
2231 sto_basis_set%lq(iset) = 5
2232 CASE ("I", "i", "J", "j")
2233 sto_basis_set%lq(iset) = 6
2234 CASE ("K", "k")
2235 sto_basis_set%lq(iset) = 7
2236 CASE ("L", "l")
2237 sto_basis_set%lq(iset) = 8
2238 CASE ("M", "m")
2239 sto_basis_set%lq(iset) = 9
2240 CASE DEFAULT
2241 CALL cp_abort(__location__, &
2242 "The requested basis set <"//trim(bsname)// &
2243 "> for element <"//trim(symbol)//"> has an invalid component: ")
2244 END SELECT
2245 END DO
2246
2247 basis_found = .true.
2248 EXIT search_loop
2249 END IF
2250 ELSE
2251 EXIT search_loop
2252 END IF
2253 END DO search_loop
2254 ELSE
2255 match = .false.
2256 END IF
2257
2258 CALL parser_release(parser)
2259
2260 END DO basis_loop
2261
2262 IF (tmp .NE. "NONE") THEN
2263 IF (.NOT. basis_found) THEN
2264 basis_set_file_name = ""
2265 DO ibasis = 1, nbasis
2266 basis_set_file_name = trim(basis_set_file_name)//"<"//trim(cbasis(ibasis))//"> "
2267 END DO
2268 CALL cp_abort(__location__, &
2269 "The requested basis set <"//trim(bsname)// &
2270 "> for element <"//trim(symbol)//"> was not "// &
2271 "found in the basis set files "// &
2272 trim(basis_set_file_name))
2273 END IF
2274 END IF
2275 DEALLOCATE (cbasis)
2276
2277 END SUBROUTINE read_sto_basis_set
2278
2279! **************************************************************************************************
2280!> \brief ...
2281!> \param sto_basis_set ...
2282!> \param gto_basis_set ...
2283!> \param ngauss ...
2284!> \param ortho ...
2285! **************************************************************************************************
2286 SUBROUTINE create_gto_from_sto_basis(sto_basis_set, gto_basis_set, ngauss, ortho)
2287
2288 TYPE(sto_basis_set_type), INTENT(IN) :: sto_basis_set
2289 TYPE(gto_basis_set_type), POINTER :: gto_basis_set
2290 INTEGER, INTENT(IN), OPTIONAL :: ngauss
2291 LOGICAL, INTENT(IN), OPTIONAL :: ortho
2292
2293 INTEGER, PARAMETER :: maxng = 6
2294
2295 CHARACTER(LEN=default_string_length) :: name, sng
2296 INTEGER :: ipgf, iset, maxl, ng, nset, nshell
2297 INTEGER, DIMENSION(:), POINTER :: lq, nq
2298 LOGICAL :: do_ortho
2299 REAL(kind=dp), DIMENSION(:), POINTER :: zet
2300 REAL(kind=dp), DIMENSION(maxng) :: gcc, zetg
2301
2302 ng = 6
2303 IF (PRESENT(ngauss)) ng = ngauss
2304 IF (ng > maxng) cpabort("Too many Gaussian primitives requested")
2305 do_ortho = .false.
2306 IF (PRESENT(ortho)) do_ortho = ortho
2307
2308 CALL allocate_gto_basis_set(gto_basis_set)
2309
2310 CALL get_sto_basis_set(sto_basis_set, name=name, nshell=nshell, nq=nq, &
2311 lq=lq, zet=zet)
2312
2313 maxl = maxval(lq)
2314 CALL init_orbital_pointers(maxl)
2315
2316 CALL integer_to_string(ng, sng)
2317 gto_basis_set%name = trim(name)//"_STO-"//trim(sng)//"G"
2318
2319 nset = nshell
2320 gto_basis_set%nset = nset
2321 CALL reallocate(gto_basis_set%lmax, 1, nset)
2322 CALL reallocate(gto_basis_set%lmin, 1, nset)
2323 CALL reallocate(gto_basis_set%npgf, 1, nset)
2324 CALL reallocate(gto_basis_set%nshell, 1, nset)
2325 CALL reallocate(gto_basis_set%n, 1, 1, 1, nset)
2326 CALL reallocate(gto_basis_set%l, 1, 1, 1, nset)
2327 CALL reallocate(gto_basis_set%zet, 1, ng, 1, nset)
2328 CALL reallocate(gto_basis_set%gcc, 1, ng, 1, 1, 1, nset)
2329
2330 DO iset = 1, nset
2331 CALL get_sto_ng(zet(iset), ng, nq(iset), lq(iset), zetg, gcc)
2332 gto_basis_set%lmax(iset) = lq(iset)
2333 gto_basis_set%lmin(iset) = lq(iset)
2334 gto_basis_set%npgf(iset) = ng
2335 gto_basis_set%nshell(iset) = 1
2336 gto_basis_set%n(1, iset) = lq(iset) + 1
2337 gto_basis_set%l(1, iset) = lq(iset)
2338 DO ipgf = 1, ng
2339 gto_basis_set%gcc(ipgf, 1, iset) = gcc(ipgf)
2340 gto_basis_set%zet(ipgf, iset) = zetg(ipgf)
2341 END DO
2342 END DO
2343
2344 CALL process_gto_basis(gto_basis_set, do_ortho, nset, maxl)
2345
2346 END SUBROUTINE create_gto_from_sto_basis
2347
2348! **************************************************************************************************
2349!> \brief ...
2350!> \param gto_basis_set ...
2351!> \param do_ortho ...
2352!> \param nset ...
2353!> \param maxl ...
2354! **************************************************************************************************
2355 SUBROUTINE process_gto_basis(gto_basis_set, do_ortho, nset, maxl)
2356
2357 TYPE(gto_basis_set_type), POINTER :: gto_basis_set
2358 LOGICAL, INTENT(IN), OPTIONAL :: do_ortho
2359 INTEGER, INTENT(IN) :: nset, maxl
2360
2361 INTEGER :: i1, i2, ico, iset, jset, l, lshell, m, &
2362 maxco, ncgf, ng, ngs, np, nsgf
2363 INTEGER, DIMENSION(0:10) :: mxf
2364 REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :) :: gal, zal, zll
2365
2366 ng = gto_basis_set%npgf(1)
2367 DO iset = 1, nset
2368 IF ((ng /= gto_basis_set%npgf(iset)) .AND. do_ortho) &
2369 cpabort("different number of primitves")
2370 END DO
2371
2372 IF (do_ortho) THEN
2373 mxf = 0
2374 DO iset = 1, nset
2375 l = gto_basis_set%l(1, iset)
2376 mxf(l) = mxf(l) + 1
2377 END DO
2378 m = maxval(mxf)
2379 IF (m > 1) THEN
2380 ALLOCATE (gal(ng, nset), zal(ng, nset), zll(m*ng, 0:maxl))
2381 DO iset = 1, nset
2382 zal(1:ng, iset) = gto_basis_set%zet(1:ng, iset)
2383 gal(1:ng, iset) = gto_basis_set%gcc(1:ng, 1, iset)
2384 END DO
2385 CALL reallocate(gto_basis_set%zet, 1, m*ng, 1, nset)
2386 CALL reallocate(gto_basis_set%gcc, 1, m*ng, 1, 1, 1, nset)
2387 DO iset = 1, nset
2388 l = gto_basis_set%l(1, iset)
2389 gto_basis_set%npgf(iset) = ng*mxf(l)
2390 END DO
2391 gto_basis_set%zet = 0.0_dp
2392 gto_basis_set%gcc = 0.0_dp
2393 zll = 0.0_dp
2394 mxf = 0
2395 DO iset = 1, nset
2396 l = gto_basis_set%l(1, iset)
2397 mxf(l) = mxf(l) + 1
2398 i1 = mxf(l)*ng - ng + 1
2399 i2 = mxf(l)*ng
2400 zll(i1:i2, l) = zal(1:ng, iset)
2401 gto_basis_set%gcc(i1:i2, 1, iset) = gal(1:ng, iset)
2402 END DO
2403 DO iset = 1, nset
2404 l = gto_basis_set%l(1, iset)
2405 gto_basis_set%zet(:, iset) = zll(:, l)
2406 END DO
2407 DO iset = 1, nset
2408 l = gto_basis_set%l(1, iset)
2409 DO jset = 1, iset - 1
2410 IF (gto_basis_set%l(1, iset) == l) THEN
2411 m = mxf(l)*ng
2412 CALL orthofun(gto_basis_set%zet(1:m, iset), gto_basis_set%gcc(1:m, 1, iset), &
2413 gto_basis_set%gcc(1:m, 1, jset), l)
2414 END IF
2415 END DO
2416 END DO
2417 DEALLOCATE (gal, zal, zll)
2418 END IF
2419 END IF
2420
2421 ngs = maxval(gto_basis_set%npgf(1:nset))
2422 CALL reallocate(gto_basis_set%set_radius, 1, nset)
2423 CALL reallocate(gto_basis_set%pgf_radius, 1, ngs, 1, nset)
2424 CALL reallocate(gto_basis_set%first_cgf, 1, 1, 1, nset)
2425 CALL reallocate(gto_basis_set%first_sgf, 1, 1, 1, nset)
2426 CALL reallocate(gto_basis_set%last_cgf, 1, 1, 1, nset)
2427 CALL reallocate(gto_basis_set%last_sgf, 1, 1, 1, nset)
2428 CALL reallocate(gto_basis_set%ncgf_set, 1, nset)
2429 CALL reallocate(gto_basis_set%nsgf_set, 1, nset)
2430
2431 maxco = 0
2432 ncgf = 0
2433 nsgf = 0
2434
2435 DO iset = 1, nset
2436 gto_basis_set%ncgf_set(iset) = 0
2437 gto_basis_set%nsgf_set(iset) = 0
2438 lshell = gto_basis_set%l(1, iset)
2439 gto_basis_set%first_cgf(1, iset) = ncgf + 1
2440 ncgf = ncgf + nco(lshell)
2441 gto_basis_set%last_cgf(1, iset) = ncgf
2442 gto_basis_set%ncgf_set(iset) = &
2443 gto_basis_set%ncgf_set(iset) + nco(lshell)
2444 gto_basis_set%first_sgf(1, iset) = nsgf + 1
2445 nsgf = nsgf + nso(lshell)
2446 gto_basis_set%last_sgf(1, iset) = nsgf
2447 gto_basis_set%nsgf_set(iset) = &
2448 gto_basis_set%nsgf_set(iset) + nso(lshell)
2449 ngs = gto_basis_set%npgf(iset)
2450 maxco = max(maxco, ngs*ncoset(lshell))
2451 END DO
2452
2453 gto_basis_set%ncgf = ncgf
2454 gto_basis_set%nsgf = nsgf
2455
2456 CALL reallocate(gto_basis_set%cphi, 1, maxco, 1, ncgf)
2457 CALL reallocate(gto_basis_set%sphi, 1, maxco, 1, nsgf)
2458 CALL reallocate(gto_basis_set%scon, 1, maxco, 1, nsgf)
2459 CALL reallocate(gto_basis_set%lx, 1, ncgf)
2460 CALL reallocate(gto_basis_set%ly, 1, ncgf)
2461 CALL reallocate(gto_basis_set%lz, 1, ncgf)
2462 CALL reallocate(gto_basis_set%m, 1, nsgf)
2463 CALL reallocate(gto_basis_set%norm_cgf, 1, ncgf)
2464 ALLOCATE (gto_basis_set%cgf_symbol(ncgf))
2465 ALLOCATE (gto_basis_set%sgf_symbol(nsgf))
2466
2467 ncgf = 0
2468 nsgf = 0
2469
2470 DO iset = 1, nset
2471 lshell = gto_basis_set%l(1, iset)
2472 np = lshell + 1
2473 DO ico = ncoset(lshell - 1) + 1, ncoset(lshell)
2474 ncgf = ncgf + 1
2475 gto_basis_set%lx(ncgf) = indco(1, ico)
2476 gto_basis_set%ly(ncgf) = indco(2, ico)
2477 gto_basis_set%lz(ncgf) = indco(3, ico)
2478 gto_basis_set%cgf_symbol(ncgf) = &
2479 cgf_symbol(np, (/gto_basis_set%lx(ncgf), &
2480 gto_basis_set%ly(ncgf), &
2481 gto_basis_set%lz(ncgf)/))
2482 END DO
2483 DO m = -lshell, lshell
2484 nsgf = nsgf + 1
2485 gto_basis_set%m(nsgf) = m
2486 gto_basis_set%sgf_symbol(nsgf) = sgf_symbol(np, lshell, m)
2487 END DO
2488 END DO
2489
2490 gto_basis_set%norm_type = -1
2491
2492 END SUBROUTINE process_gto_basis
2493
2494! **************************************************************************************************
2495!> \brief ...
2496!> \param zet ...
2497!> \param co ...
2498!> \param cr ...
2499!> \param l ...
2500! **************************************************************************************************
2501 SUBROUTINE orthofun(zet, co, cr, l)
2502 REAL(kind=dp), DIMENSION(:), INTENT(IN) :: zet
2503 REAL(kind=dp), DIMENSION(:), INTENT(INOUT) :: co, cr
2504 INTEGER, INTENT(IN) :: l
2505
2506 REAL(kind=dp) :: ss
2507
2508 CALL aovlp(l, zet, cr, cr, ss)
2509 cr(:) = cr(:)/sqrt(ss)
2510 CALL aovlp(l, zet, co, cr, ss)
2511 co(:) = co(:) - ss*cr(:)
2512 CALL aovlp(l, zet, co, co, ss)
2513 co(:) = co(:)/sqrt(ss)
2514
2515 END SUBROUTINE orthofun
2516
2517! **************************************************************************************************
2518!> \brief ...
2519!> \param l ...
2520!> \param zet ...
2521!> \param ca ...
2522!> \param cb ...
2523!> \param ss ...
2524! **************************************************************************************************
2525 SUBROUTINE aovlp(l, zet, ca, cb, ss)
2526 INTEGER, INTENT(IN) :: l
2527 REAL(kind=dp), DIMENSION(:), INTENT(IN) :: zet, ca, cb
2528 REAL(kind=dp), INTENT(OUT) :: ss
2529
2530 INTEGER :: i, j, m
2531 REAL(kind=dp) :: ab, ai, aj, s00, sss
2532
2533!
2534! use init_norm_cgf_orb
2535!
2536 m = SIZE(zet)
2537 ss = 0.0_dp
2538 DO i = 1, m
2539 ai = (2.0_dp*zet(i)/pi)**0.75_dp
2540 DO j = 1, m
2541 aj = (2.0_dp*zet(j)/pi)**0.75_dp
2542 ab = 1._dp/(zet(i) + zet(j))
2543 s00 = ai*aj*(pi*ab)**1.50_dp
2544 IF (l == 0) THEN
2545 sss = s00
2546 ELSEIF (l == 1) THEN
2547 sss = s00*ab*0.5_dp
2548 ELSE
2549 cpabort("aovlp lvalue")
2550 END IF
2551 ss = ss + sss*ca(i)*cb(j)
2552 END DO
2553 END DO
2554
2555 END SUBROUTINE aovlp
2556
2557! **************************************************************************************************
2558!> \brief ...
2559!> \param z ...
2560!> \param ne ...
2561!> \param n ...
2562!> \param l ...
2563!> \return ...
2564! **************************************************************************************************
2565 PURE FUNCTION srules(z, ne, n, l)
2566 ! Slater rules
2567 INTEGER, INTENT(IN) :: z
2568 INTEGER, DIMENSION(:, :), INTENT(IN) :: ne
2569 INTEGER, INTENT(IN) :: n, l
2570 REAL(dp) :: srules
2571
2572 REAL(dp), DIMENSION(7), PARAMETER :: &
2573 xns = (/1.0_dp, 2.0_dp, 3.0_dp, 3.7_dp, 4.0_dp, 4.2_dp, 4.4_dp/)
2574
2575 INTEGER :: i, l1, l2, m, m1, m2, nn
2576 REAL(dp) :: s
2577
2578 s = 0.0_dp
2579 ! The complete shell
2580 l1 = min(l + 1, 4)
2581 nn = min(n, 7)
2582 IF (l1 == 1) l2 = 2
2583 IF (l1 == 2) l2 = 1
2584 IF (l1 == 3) l2 = 4
2585 IF (l1 == 4) l2 = 3
2586 ! Rule a) no contribution from shells further out
2587 ! Rule b) 0.35 (1s 0.3) from each other electron in the same shell
2588 IF (n == 1) THEN
2589 m = ne(1, 1)
2590 s = s + 0.3_dp*real(m - 1, dp)
2591 ELSE
2592 m = ne(l1, nn) + ne(l2, nn)
2593 s = s + 0.35_dp*real(m - 1, dp)
2594 END IF
2595 ! Rule c) if (s,p) shell 0.85 from each electron with n-1, and 1.0
2596 ! from all electrons further in
2597 IF (l1 + l2 == 3) THEN
2598 IF (nn > 1) THEN
2599 m1 = ne(1, nn - 1) + ne(2, nn - 1) + ne(3, nn - 1) + ne(4, nn - 1)
2600 m2 = 0
2601 DO i = 1, nn - 2
2602 m2 = m2 + ne(1, i) + ne(2, i) + ne(3, i) + ne(4, i)
2603 END DO
2604 s = s + 0.85_dp*real(m1, dp) + 1._dp*real(m2, dp)
2605 END IF
2606 ELSE
2607 ! Rule d) if (d,f) shell 1.0 from each electron inside
2608 m = 0
2609 DO i = 1, nn - 1
2610 m = m + ne(1, i) + ne(2, i) + ne(3, i) + ne(4, i)
2611 END DO
2612 s = s + 1._dp*real(m, dp)
2613 END IF
2614 ! Slater exponent is (Z-S)/NS
2615 srules = (real(z, dp) - s)/xns(nn)
2616 END FUNCTION srules
2617
2618! **************************************************************************************************
2619!> \brief sort basis sets w.r.t. radius
2620!> \param basis_set ...
2621!> \param sort_method ...
2622! **************************************************************************************************
2623 SUBROUTINE sort_gto_basis_set(basis_set, sort_method)
2624 TYPE(gto_basis_set_type), INTENT(INOUT) :: basis_set
2625 INTEGER, INTENT(IN) :: sort_method
2626
2627 CHARACTER(LEN=12), DIMENSION(:), POINTER :: cgf_symbol
2628 CHARACTER(LEN=6), DIMENSION(:), POINTER :: sgf_symbol
2629 INTEGER :: ic, ic_max, icgf, icgf_new, icgf_old, ico, is, is_max, iset, isgf, isgf_new, &
2630 isgf_old, ishell, lshell, maxco, maxpgf, maxshell, mm, nc, ncgf, ns, nset
2631 INTEGER, ALLOCATABLE, DIMENSION(:) :: sort_index
2632 INTEGER, ALLOCATABLE, DIMENSION(:, :) :: icgf_set, isgf_set
2633 INTEGER, DIMENSION(:), POINTER :: lx, ly, lz, m, npgf
2634 REAL(dp), ALLOCATABLE, DIMENSION(:) :: tmp
2635 REAL(dp), DIMENSION(:), POINTER :: set_radius
2636 REAL(dp), DIMENSION(:, :), POINTER :: zet
2637 REAL(kind=dp), DIMENSION(:), POINTER :: norm_cgf
2638 REAL(kind=dp), DIMENSION(:, :), POINTER :: cphi, scon, sphi
2639
2640 NULLIFY (set_radius, zet)
2641
2642 IF (sort_method == basis_sort_default) RETURN
2643
2644 CALL get_gto_basis_set(gto_basis_set=basis_set, &
2645 nset=nset, &
2646 maxshell=maxshell, &
2647 maxpgf=maxpgf, &
2648 maxco=maxco, &
2649 ncgf=ncgf, &
2650 npgf=npgf, &
2651 set_radius=set_radius, &
2652 zet=zet)
2653
2654 ALLOCATE (sort_index(nset))
2655 ALLOCATE (tmp(nset))
2656 SELECT CASE (sort_method)
2657 CASE (basis_sort_zet)
2658 DO iset = 1, nset
2659 tmp(iset) = minval(basis_set%zet(:npgf(iset), iset))
2660 END DO
2661 CASE DEFAULT
2662 cpabort("Request basis sort criterion not implemented.")
2663 END SELECT
2664
2665 CALL sort(tmp(1:nset), nset, sort_index)
2666
2667 ic_max = 0
2668 is_max = 0
2669 DO iset = 1, nset
2670 ic = 0
2671 is = 0
2672 DO ishell = 1, basis_set%nshell(iset)
2673 DO ico = 1, nco(basis_set%l(ishell, iset))
2674 ic = ic + 1
2675 IF (ic > ic_max) ic_max = ic
2676 END DO
2677 lshell = basis_set%l(ishell, iset)
2678 DO mm = -lshell, lshell
2679 is = is + 1
2680 IF (is > is_max) is_max = is
2681 END DO
2682 END DO
2683 END DO
2684
2685 icgf = 0
2686 isgf = 0
2687 ALLOCATE (icgf_set(nset, ic_max))
2688 icgf_set(:, :) = 0
2689 ALLOCATE (isgf_set(nset, is_max))
2690 isgf_set(:, :) = 0
2691
2692 DO iset = 1, nset
2693 ic = 0
2694 is = 0
2695 DO ishell = 1, basis_set%nshell(iset)
2696 DO ico = 1, nco(basis_set%l(ishell, iset))
2697 icgf = icgf + 1
2698 ic = ic + 1
2699 icgf_set(iset, ic) = icgf
2700 END DO
2701 lshell = basis_set%l(ishell, iset)
2702 DO mm = -lshell, lshell
2703 isgf = isgf + 1
2704 is = is + 1
2705 isgf_set(iset, is) = isgf
2706 END DO
2707 END DO
2708 END DO
2709
2710 ALLOCATE (cgf_symbol(SIZE(basis_set%cgf_symbol)))
2711 ALLOCATE (norm_cgf(SIZE(basis_set%norm_cgf)))
2712 ALLOCATE (lx(SIZE(basis_set%lx)))
2713 ALLOCATE (ly(SIZE(basis_set%ly)))
2714 ALLOCATE (lz(SIZE(basis_set%lz)))
2715 ALLOCATE (cphi(SIZE(basis_set%cphi, 1), SIZE(basis_set%cphi, 2)))
2716 cphi = 0.0_dp
2717 ALLOCATE (sphi(SIZE(basis_set%sphi, 1), SIZE(basis_set%sphi, 2)))
2718 sphi = 0.0_dp
2719 ALLOCATE (scon(SIZE(basis_set%scon, 1), SIZE(basis_set%scon, 2)))
2720 scon = 0.0_dp
2721
2722 ALLOCATE (sgf_symbol(SIZE(basis_set%sgf_symbol)))
2723 ALLOCATE (m(SIZE(basis_set%m)))
2724
2725 icgf_new = 0
2726 isgf_new = 0
2727 DO iset = 1, nset
2728 DO ic = 1, ic_max
2729 icgf_old = icgf_set(sort_index(iset), ic)
2730 IF (icgf_old == 0) cycle
2731 icgf_new = icgf_new + 1
2732 norm_cgf(icgf_new) = basis_set%norm_cgf(icgf_old)
2733 lx(icgf_new) = basis_set%lx(icgf_old)
2734 ly(icgf_new) = basis_set%ly(icgf_old)
2735 lz(icgf_new) = basis_set%lz(icgf_old)
2736 cphi(:, icgf_new) = basis_set%cphi(:, icgf_old)
2737 cgf_symbol(icgf_new) = basis_set%cgf_symbol(icgf_old)
2738 END DO
2739 DO is = 1, is_max
2740 isgf_old = isgf_set(sort_index(iset), is)
2741 IF (isgf_old == 0) cycle
2742 isgf_new = isgf_new + 1
2743 m(isgf_new) = basis_set%m(isgf_old)
2744 sphi(:, isgf_new) = basis_set%sphi(:, isgf_old)
2745 scon(:, isgf_new) = basis_set%scon(:, isgf_old)
2746 sgf_symbol(isgf_new) = basis_set%sgf_symbol(isgf_old)
2747 END DO
2748 END DO
2749
2750 DEALLOCATE (basis_set%cgf_symbol)
2751 basis_set%cgf_symbol => cgf_symbol
2752 DEALLOCATE (basis_set%norm_cgf)
2753 basis_set%norm_cgf => norm_cgf
2754 DEALLOCATE (basis_set%lx)
2755 basis_set%lx => lx
2756 DEALLOCATE (basis_set%ly)
2757 basis_set%ly => ly
2758 DEALLOCATE (basis_set%lz)
2759 basis_set%lz => lz
2760 DEALLOCATE (basis_set%cphi)
2761 basis_set%cphi => cphi
2762 DEALLOCATE (basis_set%sphi)
2763 basis_set%sphi => sphi
2764 DEALLOCATE (basis_set%scon)
2765 basis_set%scon => scon
2766
2767 DEALLOCATE (basis_set%m)
2768 basis_set%m => m
2769 DEALLOCATE (basis_set%sgf_symbol)
2770 basis_set%sgf_symbol => sgf_symbol
2771
2772 basis_set%lmax = basis_set%lmax(sort_index)
2773 basis_set%lmin = basis_set%lmin(sort_index)
2774 basis_set%npgf = basis_set%npgf(sort_index)
2775 basis_set%nshell = basis_set%nshell(sort_index)
2776 basis_set%ncgf_set = basis_set%ncgf_set(sort_index)
2777 basis_set%nsgf_set = basis_set%nsgf_set(sort_index)
2778
2779 basis_set%n(:, :) = basis_set%n(:, sort_index)
2780 basis_set%l(:, :) = basis_set%l(:, sort_index)
2781 basis_set%zet(:, :) = basis_set%zet(:, sort_index)
2782
2783 basis_set%gcc(:, :, :) = basis_set%gcc(:, :, sort_index)
2784 basis_set%set_radius(:) = basis_set%set_radius(sort_index)
2785 basis_set%pgf_radius(:, :) = basis_set%pgf_radius(:, sort_index)
2786
2787 nc = 0
2788 ns = 0
2789 DO iset = 1, nset
2790 DO ishell = 1, basis_set%nshell(iset)
2791 lshell = basis_set%l(ishell, iset)
2792 basis_set%first_cgf(ishell, iset) = nc + 1
2793 nc = nc + nco(lshell)
2794 basis_set%last_cgf(ishell, iset) = nc
2795 basis_set%first_sgf(ishell, iset) = ns + 1
2796 ns = ns + nso(lshell)
2797 basis_set%last_sgf(ishell, iset) = ns
2798 END DO
2799 END DO
2800
2801 END SUBROUTINE
2802
2803END MODULE basis_set_types
static void dgemm(const char transa, const char transb, const int m, const int n, const int k, const double alpha, const double *a, const int lda, const double *b, const int ldb, const double beta, double *c, const int ldc)
Convenient wrapper to hide Fortran nature of dgemm_, swapping a and b.
Calculation of Coulomb integrals over Cartesian Gaussian-type functions (electron repulsion integrals...
Definition ai_coulomb.F:41
subroutine, public coulomb2(la_max, npgfa, zeta, rpgfa, la_min, lc_max, npgfc, zetc, rpgfc, lc_min, rac, rac2, vac, v, f, maxder, vac_plus)
Calculation of the primitive two-center Coulomb integrals over Cartesian Gaussian-type functions.
Definition ai_coulomb.F:86
subroutine, public get_gto_basis_set(gto_basis_set, name, aliases, norm_type, kind_radius, ncgf, nset, nsgf, cgf_symbol, sgf_symbol, norm_cgf, set_radius, lmax, lmin, lx, ly, lz, m, ncgf_set, npgf, nsgf_set, nshell, cphi, pgf_radius, sphi, scon, zet, first_cgf, first_sgf, l, last_cgf, last_sgf, n, gcc, maxco, maxl, maxpgf, maxsgf_set, maxshell, maxso, nco_sum, npgf_sum, nshell_sum, maxder, short_kind_radius, npgf_seg_sum)
...
integer, parameter, public basis_sort_zet
subroutine, public process_gto_basis(gto_basis_set, do_ortho, nset, maxl)
...
subroutine, public set_gto_basis_set(gto_basis_set, name, aliases, norm_type, kind_radius, ncgf, nset, nsgf, cgf_symbol, sgf_symbol, norm_cgf, set_radius, lmax, lmin, lx, ly, lz, m, ncgf_set, npgf, nsgf_set, nshell, cphi, pgf_radius, sphi, scon, zet, first_cgf, first_sgf, l, last_cgf, last_sgf, n, gcc, short_kind_radius)
Set the components of Gaussian-type orbital (GTO) basis set data set.
subroutine, public deallocate_gto_basis_set(gto_basis_set)
...
pure real(dp) function, public srules(z, ne, n, l)
...
subroutine, public write_orb_basis_set(orb_basis_set, output_unit, header)
Write a Gaussian-type orbital (GTO) basis set data set to the output unit.
subroutine, public sort_gto_basis_set(basis_set, sort_method)
sort basis sets w.r.t. radius
subroutine, public deallocate_sto_basis_set(sto_basis_set)
...
subroutine, public init_aux_basis_set(gto_basis_set)
...
subroutine, public allocate_gto_basis_set(gto_basis_set)
...
subroutine, public copy_gto_basis_set(basis_set_in, basis_set_out)
...
subroutine, public init_cphi_and_sphi(gto_basis_set)
...
subroutine, public combine_basis_sets(basis_set, basis_set_add)
...
subroutine, public write_gto_basis_set(gto_basis_set, output_unit, header)
Write a Gaussian-type orbital (GTO) basis set data set to the output unit.
subroutine, public allocate_sto_basis_set(sto_basis_set)
...
subroutine, public create_primitive_basis_set(basis_set, pbasis, lmax)
...
subroutine, public create_gto_from_sto_basis(sto_basis_set, gto_basis_set, ngauss, ortho)
...
subroutine, public set_sto_basis_set(sto_basis_set, name, nshell, symbol, nq, lq, zet)
...
subroutine, public read_sto_basis_set(element_symbol, basis_set_name, sto_basis_set, para_env, dft_section)
...
subroutine, public init_orb_basis_set(gto_basis_set)
Initialise a Gaussian-type orbital (GTO) basis set data set.
integer, parameter, public basis_sort_default
collects all references to literature in CP2K as new algorithms / method are included from literature...
integer, save, public vandevondele2007
logical function, public cp_sll_val_next(iterator, el_att)
returns true if the actual element is valid (i.e. iterator ont at end) moves the iterator to the next...
Utility routines to read data from files. Kept as close as possible to the old parser because.
subroutine, public parser_search_string(parser, string, ignore_case, found, line, begin_line, search_from_begin_of_file)
Search a string pattern in a file defined by its logical unit number "unit". A case sensitive search ...
Utility routines to read data from files. Kept as close as possible to the old parser because.
subroutine, public parser_release(parser)
releases the parser
subroutine, public parser_create(parser, file_name, unit_nr, para_env, end_section_label, separator_chars, comment_char, continuation_char, quote_char, section_char, parse_white_lines, initial_variables, apply_preprocessing)
Start a parser run. Initial variables allow to @SET stuff before opening the file.
objects that represent the structure of input sections and the data contained in an input section
subroutine, public section_vals_list_get(section_vals, keyword_name, i_rep_section, list)
returns the requested list
subroutine, public section_vals_val_get(section_vals, keyword_name, i_rep_section, i_rep_val, n_rep_val, val, l_val, i_val, r_val, c_val, l_vals, i_vals, r_vals, c_vals, explicit)
returns the requested value
a wrapper for basic fortran types.
subroutine, public val_get(val, has_l, has_i, has_r, has_lc, has_c, l_val, l_vals, i_val, i_vals, r_val, r_vals, c_val, c_vals, len_c, type_of_var, enum)
returns the stored values
Defines the basic variable types.
Definition kinds.F:23
integer, parameter, public dp
Definition kinds.F:34
integer, parameter, public default_string_length
Definition kinds.F:57
integer, parameter, public default_path_length
Definition kinds.F:58
An array-based list which grows on demand. When the internal array is full, a new array of twice the ...
Definition list.F:24
Definition of mathematical constants and functions.
real(kind=dp), parameter, public pi
real(kind=dp), dimension(-1:2 *maxfac+1), parameter, public dfac
Utility routines for the memory handling.
Interface to the message passing library MPI.
Provides Cartesian and spherical orbital pointers and indices.
subroutine, public init_orbital_pointers(maxl)
Initialize or update the orbital pointers.
integer, dimension(:, :, :), allocatable, public co
integer, dimension(:), allocatable, public nco
integer, dimension(:), allocatable, public nsoset
integer, dimension(:), allocatable, public ncoset
integer, dimension(:, :, :), allocatable, public coset
integer, dimension(:, :), allocatable, public indco
integer, dimension(:), allocatable, public nso
orbital_symbols
character(len=12) function, public cgf_symbol(n, lxyz)
Build a Cartesian orbital symbol (orbital labels for printing).
character(len=6) function, public sgf_symbol(n, l, m)
Build a spherical orbital symbol (orbital labels for printing).
Calculation of the spherical harmonics and the corresponding orbital transformation matrices.
subroutine, public init_spherical_harmonics(maxl, output_unit)
Initialize or update the orbital transformation matrices.
type(orbtramat_type), dimension(:), pointer, public orbtramat
subroutine, public get_sto_ng(zeta, n, nq, lq, alpha, coef)
return STO-NG parameters; INPUT: zeta (Slater exponent) n (Expansion length) nq (principle quantum nu...
Definition sto_ng.F:54
Utilities for string manipulations.
subroutine, public integer_to_string(inumber, string)
Converts an integer number to a string. The WRITE statement will return an error message,...
character(len=1), parameter, public newline
subroutine, public remove_word(string)
remove a word from a string (words are separated by white spaces)
elemental subroutine, public uppercase(string)
Convert all lower case characters in a string to upper case.
All kind of helpful little routines.
Definition util.F:14
represent a single linked list that stores pointers to the elements
a type to have a wrapper that stores any basic fortran type
stores all the informations relevant to an mpi environment