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