(git:d18deda)
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, &
137 srules
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 DEFAULT
1107 cpabort("Normalization method not specified")
1108 END SELECT
1109
1110 ! Initialise the transformation matrices "pgf" -> "cgf"
1111
1112 CALL init_cphi_and_sphi(gto_basis_set)
1113
1114 CALL timestop(handle)
1115
1116 END SUBROUTINE init_orb_basis_set
1117
1118! **************************************************************************************************
1119!> \brief Normalise the primitive Cartesian Gaussian functions. The normalization
1120!> factor is included in the Gaussian contraction coefficients.
1121!> \param gto_basis_set ...
1122!> \author MK
1123! **************************************************************************************************
1124 SUBROUTINE normalise_gcc_orb(gto_basis_set)
1125
1126 TYPE(gto_basis_set_type), POINTER :: gto_basis_set
1127
1128 INTEGER :: ipgf, iset, ishell, l
1129 REAL(kind=dp) :: expzet, gcca, prefac, zeta
1130
1131 DO iset = 1, gto_basis_set%nset
1132 DO ishell = 1, gto_basis_set%nshell(iset)
1133 l = gto_basis_set%l(ishell, iset)
1134 expzet = 0.25_dp*real(2*l + 3, dp)
1135 prefac = 2.0_dp**l*(2.0_dp/pi)**0.75_dp
1136 DO ipgf = 1, gto_basis_set%npgf(iset)
1137 gcca = gto_basis_set%gcc(ipgf, ishell, iset)
1138 zeta = gto_basis_set%zet(ipgf, iset)
1139 gto_basis_set%gcc(ipgf, ishell, iset) = prefac*zeta**expzet*gcca
1140 END DO
1141 END DO
1142 END DO
1143
1144 END SUBROUTINE normalise_gcc_orb
1145
1146! **************************************************************************************************
1147!> \brief Read a Gaussian-type orbital (GTO) basis set from the database file.
1148!> \param element_symbol ...
1149!> \param basis_set_name ...
1150!> \param gto_basis_set ...
1151!> \param para_env ...
1152!> \param dft_section ...
1153!> \author MK
1154! **************************************************************************************************
1155 SUBROUTINE read_gto_basis_set1(element_symbol, basis_set_name, gto_basis_set, &
1156 para_env, dft_section)
1157
1158 CHARACTER(LEN=*), INTENT(IN) :: element_symbol, basis_set_name
1159 TYPE(gto_basis_set_type), INTENT(INOUT) :: gto_basis_set
1160 TYPE(mp_para_env_type), POINTER :: para_env
1161 TYPE(section_vals_type), POINTER :: dft_section
1162
1163 CHARACTER(LEN=240) :: line
1164 CHARACTER(LEN=242) :: line2
1165 CHARACTER(len=default_path_length) :: basis_set_file_name, tmp
1166 CHARACTER(LEN=default_path_length), DIMENSION(:), &
1167 POINTER :: cbasis
1168 CHARACTER(LEN=LEN(basis_set_name)) :: bsname
1169 CHARACTER(LEN=LEN(basis_set_name)+2) :: bsname2
1170 CHARACTER(LEN=LEN(element_symbol)) :: symbol
1171 CHARACTER(LEN=LEN(element_symbol)+2) :: symbol2
1172 INTEGER :: i, ibasis, ico, ipgf, irep, iset, ishell, lshell, m, maxco, maxl, maxpgf, &
1173 maxshell, nbasis, ncgf, nmin, nset, nsgf, sort_method, strlen1, strlen2
1174 INTEGER, DIMENSION(:), POINTER :: lmax, lmin, npgf, nshell
1175 INTEGER, DIMENSION(:, :), POINTER :: l, n
1176 LOGICAL :: basis_found, found, match
1177 REAL(KIND=dp), DIMENSION(:, :), POINTER :: zet
1178 REAL(KIND=dp), DIMENSION(:, :, :), POINTER :: gcc
1179 TYPE(cp_parser_type) :: parser
1180
1181 line = ""
1182 line2 = ""
1183 symbol = ""
1184 symbol2 = ""
1185 bsname = ""
1186 bsname2 = ""
1187
1188 nbasis = 1
1189
1190 gto_basis_set%name = basis_set_name
1191 gto_basis_set%aliases = basis_set_name
1192 CALL section_vals_val_get(dft_section, "BASIS_SET_FILE_NAME", &
1193 n_rep_val=nbasis)
1194 ALLOCATE (cbasis(nbasis))
1195 DO ibasis = 1, nbasis
1196 CALL section_vals_val_get(dft_section, "BASIS_SET_FILE_NAME", &
1197 i_rep_val=ibasis, c_val=cbasis(ibasis))
1198 basis_set_file_name = cbasis(ibasis)
1199 tmp = basis_set_file_name
1200 CALL uppercase(tmp)
1201 IF (index(tmp, "MOLOPT") .NE. 0) CALL cite_reference(vandevondele2007)
1202 END DO
1203
1204 ! Search for the requested basis set in the basis set file
1205 ! until the basis set is found or the end of file is reached
1206
1207 basis_found = .false.
1208 basis_loop: DO ibasis = 1, nbasis
1209 IF (basis_found) EXIT basis_loop
1210 basis_set_file_name = cbasis(ibasis)
1211 CALL parser_create(parser, basis_set_file_name, para_env=para_env)
1212
1213 bsname = basis_set_name
1214 symbol = element_symbol
1215 irep = 0
1216
1217 tmp = basis_set_name
1218 CALL uppercase(tmp)
1219 IF (index(tmp, "MOLOPT") .NE. 0) CALL cite_reference(vandevondele2007)
1220
1221 nset = 0
1222 maxshell = 0
1223 maxpgf = 0
1224 maxco = 0
1225 ncgf = 0
1226 nsgf = 0
1227 gto_basis_set%nset = nset
1228 gto_basis_set%ncgf = ncgf
1229 gto_basis_set%nsgf = nsgf
1230 CALL reallocate(gto_basis_set%lmax, 1, nset)
1231 CALL reallocate(gto_basis_set%lmin, 1, nset)
1232 CALL reallocate(gto_basis_set%npgf, 1, nset)
1233 CALL reallocate(gto_basis_set%nshell, 1, nset)
1234 CALL reallocate(gto_basis_set%n, 1, maxshell, 1, nset)
1235 CALL reallocate(gto_basis_set%l, 1, maxshell, 1, nset)
1236 CALL reallocate(gto_basis_set%zet, 1, maxpgf, 1, nset)
1237 CALL reallocate(gto_basis_set%gcc, 1, maxpgf, 1, maxshell, 1, nset)
1238 CALL reallocate(gto_basis_set%set_radius, 1, nset)
1239 CALL reallocate(gto_basis_set%pgf_radius, 1, maxpgf, 1, nset)
1240 CALL reallocate(gto_basis_set%first_cgf, 1, maxshell, 1, nset)
1241 CALL reallocate(gto_basis_set%first_sgf, 1, maxshell, 1, nset)
1242 CALL reallocate(gto_basis_set%last_cgf, 1, maxshell, 1, nset)
1243 CALL reallocate(gto_basis_set%last_sgf, 1, maxshell, 1, nset)
1244 CALL reallocate(gto_basis_set%ncgf_set, 1, nset)
1245 CALL reallocate(gto_basis_set%nsgf_set, 1, nset)
1246 CALL reallocate(gto_basis_set%cphi, 1, maxco, 1, ncgf)
1247 CALL reallocate(gto_basis_set%sphi, 1, maxco, 1, nsgf)
1248 CALL reallocate(gto_basis_set%scon, 1, maxco, 1, nsgf)
1249 CALL reallocate(gto_basis_set%lx, 1, ncgf)
1250 CALL reallocate(gto_basis_set%ly, 1, ncgf)
1251 CALL reallocate(gto_basis_set%lz, 1, ncgf)
1252 CALL reallocate(gto_basis_set%m, 1, nsgf)
1253 CALL reallocate(gto_basis_set%norm_cgf, 1, ncgf)
1254
1255 IF (tmp .NE. "NONE") THEN
1256 search_loop: DO
1257
1258 CALL parser_search_string(parser, trim(bsname), .true., found, line)
1259 IF (found) THEN
1260 CALL uppercase(symbol)
1261 CALL uppercase(bsname)
1262
1263 match = .false.
1264 CALL uppercase(line)
1265 ! Check both the element symbol and the basis set name
1266 line2 = " "//line//" "
1267 symbol2 = " "//trim(symbol)//" "
1268 bsname2 = " "//trim(bsname)//" "
1269 strlen1 = len_trim(symbol2) + 1
1270 strlen2 = len_trim(bsname2) + 1
1271
1272 IF ((index(line2, symbol2(:strlen1)) > 0) .AND. &
1273 (index(line2, bsname2(:strlen2)) > 0)) match = .true.
1274 IF (match) THEN
1275 ! copy all names into aliases field
1276 i = index(line2, symbol2(:strlen1))
1277 i = i + 1 + index(line2(i + 1:), " ")
1278 gto_basis_set%aliases = line2(i:)
1279
1280 NULLIFY (gcc, l, lmax, lmin, n, npgf, nshell, zet)
1281 ! Read the basis set information
1282 CALL parser_get_object(parser, nset, newline=.true.)
1283
1284 CALL reallocate(npgf, 1, nset)
1285 CALL reallocate(nshell, 1, nset)
1286 CALL reallocate(lmax, 1, nset)
1287 CALL reallocate(lmin, 1, nset)
1288 CALL reallocate(n, 1, 1, 1, nset)
1289
1290 maxl = 0
1291 maxpgf = 0
1292 maxshell = 0
1293
1294 DO iset = 1, nset
1295 CALL parser_get_object(parser, n(1, iset), newline=.true.)
1296 CALL parser_get_object(parser, lmin(iset))
1297 CALL parser_get_object(parser, lmax(iset))
1298 CALL parser_get_object(parser, npgf(iset))
1299 maxl = max(maxl, lmax(iset))
1300 IF (npgf(iset) > maxpgf) THEN
1301 maxpgf = npgf(iset)
1302 CALL reallocate(zet, 1, maxpgf, 1, nset)
1303 CALL reallocate(gcc, 1, maxpgf, 1, maxshell, 1, nset)
1304 END IF
1305 nshell(iset) = 0
1306 DO lshell = lmin(iset), lmax(iset)
1307 nmin = n(1, iset) + lshell - lmin(iset)
1308 CALL parser_get_object(parser, ishell)
1309 nshell(iset) = nshell(iset) + ishell
1310 IF (nshell(iset) > maxshell) THEN
1311 maxshell = nshell(iset)
1312 CALL reallocate(n, 1, maxshell, 1, nset)
1313 CALL reallocate(l, 1, maxshell, 1, nset)
1314 CALL reallocate(gcc, 1, maxpgf, 1, maxshell, 1, nset)
1315 END IF
1316 DO i = 1, ishell
1317 n(nshell(iset) - ishell + i, iset) = nmin + i - 1
1318 l(nshell(iset) - ishell + i, iset) = lshell
1319 END DO
1320 END DO
1321 DO ipgf = 1, npgf(iset)
1322 CALL parser_get_object(parser, zet(ipgf, iset), newline=.true.)
1323 DO ishell = 1, nshell(iset)
1324 CALL parser_get_object(parser, gcc(ipgf, ishell, iset))
1325 END DO
1326 END DO
1327 END DO
1328
1329 ! Maximum angular momentum quantum number of the atomic kind
1330
1331 CALL init_orbital_pointers(maxl)
1332
1333 ! Allocate the global variables
1334
1335 gto_basis_set%nset = nset
1336 CALL reallocate(gto_basis_set%lmax, 1, nset)
1337 CALL reallocate(gto_basis_set%lmin, 1, nset)
1338 CALL reallocate(gto_basis_set%npgf, 1, nset)
1339 CALL reallocate(gto_basis_set%nshell, 1, nset)
1340 CALL reallocate(gto_basis_set%n, 1, maxshell, 1, nset)
1341 CALL reallocate(gto_basis_set%l, 1, maxshell, 1, nset)
1342 CALL reallocate(gto_basis_set%zet, 1, maxpgf, 1, nset)
1343 CALL reallocate(gto_basis_set%gcc, 1, maxpgf, 1, maxshell, 1, nset)
1344
1345 ! Copy the basis set information into the data structure
1346
1347 DO iset = 1, nset
1348 gto_basis_set%lmax(iset) = lmax(iset)
1349 gto_basis_set%lmin(iset) = lmin(iset)
1350 gto_basis_set%npgf(iset) = npgf(iset)
1351 gto_basis_set%nshell(iset) = nshell(iset)
1352 DO ishell = 1, nshell(iset)
1353 gto_basis_set%n(ishell, iset) = n(ishell, iset)
1354 gto_basis_set%l(ishell, iset) = l(ishell, iset)
1355 DO ipgf = 1, npgf(iset)
1356 gto_basis_set%gcc(ipgf, ishell, iset) = gcc(ipgf, ishell, iset)
1357 END DO
1358 END DO
1359 DO ipgf = 1, npgf(iset)
1360 gto_basis_set%zet(ipgf, iset) = zet(ipgf, iset)
1361 END DO
1362 END DO
1363
1364 ! Initialise the depending atomic kind information
1365
1366 CALL reallocate(gto_basis_set%set_radius, 1, nset)
1367 CALL reallocate(gto_basis_set%pgf_radius, 1, maxpgf, 1, nset)
1368 CALL reallocate(gto_basis_set%first_cgf, 1, maxshell, 1, nset)
1369 CALL reallocate(gto_basis_set%first_sgf, 1, maxshell, 1, nset)
1370 CALL reallocate(gto_basis_set%last_cgf, 1, maxshell, 1, nset)
1371 CALL reallocate(gto_basis_set%last_sgf, 1, maxshell, 1, nset)
1372 CALL reallocate(gto_basis_set%ncgf_set, 1, nset)
1373 CALL reallocate(gto_basis_set%nsgf_set, 1, nset)
1374
1375 maxco = 0
1376 ncgf = 0
1377 nsgf = 0
1378
1379 DO iset = 1, nset
1380 gto_basis_set%ncgf_set(iset) = 0
1381 gto_basis_set%nsgf_set(iset) = 0
1382 DO ishell = 1, nshell(iset)
1383 lshell = gto_basis_set%l(ishell, iset)
1384 gto_basis_set%first_cgf(ishell, iset) = ncgf + 1
1385 ncgf = ncgf + nco(lshell)
1386 gto_basis_set%last_cgf(ishell, iset) = ncgf
1387 gto_basis_set%ncgf_set(iset) = &
1388 gto_basis_set%ncgf_set(iset) + nco(lshell)
1389 gto_basis_set%first_sgf(ishell, iset) = nsgf + 1
1390 nsgf = nsgf + nso(lshell)
1391 gto_basis_set%last_sgf(ishell, iset) = nsgf
1392 gto_basis_set%nsgf_set(iset) = &
1393 gto_basis_set%nsgf_set(iset) + nso(lshell)
1394 END DO
1395 maxco = max(maxco, npgf(iset)*ncoset(lmax(iset)))
1396 END DO
1397
1398 gto_basis_set%ncgf = ncgf
1399 gto_basis_set%nsgf = nsgf
1400
1401 CALL reallocate(gto_basis_set%cphi, 1, maxco, 1, ncgf)
1402 CALL reallocate(gto_basis_set%sphi, 1, maxco, 1, nsgf)
1403 CALL reallocate(gto_basis_set%scon, 1, maxco, 1, nsgf)
1404 CALL reallocate(gto_basis_set%lx, 1, ncgf)
1405 CALL reallocate(gto_basis_set%ly, 1, ncgf)
1406 CALL reallocate(gto_basis_set%lz, 1, ncgf)
1407 CALL reallocate(gto_basis_set%m, 1, nsgf)
1408 CALL reallocate(gto_basis_set%norm_cgf, 1, ncgf)
1409 ALLOCATE (gto_basis_set%cgf_symbol(ncgf))
1410
1411 ALLOCATE (gto_basis_set%sgf_symbol(nsgf))
1412
1413 ncgf = 0
1414 nsgf = 0
1415
1416 DO iset = 1, nset
1417 DO ishell = 1, nshell(iset)
1418 lshell = gto_basis_set%l(ishell, iset)
1419 DO ico = ncoset(lshell - 1) + 1, ncoset(lshell)
1420 ncgf = ncgf + 1
1421 gto_basis_set%lx(ncgf) = indco(1, ico)
1422 gto_basis_set%ly(ncgf) = indco(2, ico)
1423 gto_basis_set%lz(ncgf) = indco(3, ico)
1424 gto_basis_set%cgf_symbol(ncgf) = &
1425 cgf_symbol(n(ishell, iset), (/gto_basis_set%lx(ncgf), &
1426 gto_basis_set%ly(ncgf), &
1427 gto_basis_set%lz(ncgf)/))
1428 END DO
1429 DO m = -lshell, lshell
1430 nsgf = nsgf + 1
1431 gto_basis_set%m(nsgf) = m
1432 gto_basis_set%sgf_symbol(nsgf) = &
1433 sgf_symbol(n(ishell, iset), lshell, m)
1434 END DO
1435 END DO
1436 END DO
1437
1438 DEALLOCATE (gcc, l, lmax, lmin, n, npgf, nshell, zet)
1439
1440 basis_found = .true.
1441 EXIT search_loop
1442 END IF
1443 ELSE
1444 EXIT search_loop
1445 END IF
1446 END DO search_loop
1447 ELSE
1448 match = .false.
1449 ALLOCATE (gto_basis_set%cgf_symbol(ncgf))
1450 ALLOCATE (gto_basis_set%sgf_symbol(nsgf))
1451 END IF
1452
1453 CALL parser_release(parser)
1454
1455 END DO basis_loop
1456
1457 IF (tmp .NE. "NONE") THEN
1458 IF (.NOT. basis_found) THEN
1459 basis_set_file_name = ""
1460 DO ibasis = 1, nbasis
1461 basis_set_file_name = trim(basis_set_file_name)//"<"//trim(cbasis(ibasis))//"> "
1462 END DO
1463 CALL cp_abort(__location__, &
1464 "The requested basis set <"//trim(bsname)// &
1465 "> for element <"//trim(symbol)//"> was not "// &
1466 "found in the basis set files "// &
1467 trim(basis_set_file_name))
1468 END IF
1469 END IF
1470 DEALLOCATE (cbasis)
1471
1472 CALL section_vals_val_get(dft_section, "SORT_BASIS", i_val=sort_method)
1473 CALL sort_gto_basis_set(gto_basis_set, sort_method)
1474
1475 END SUBROUTINE read_gto_basis_set1
1476
1477! **************************************************************************************************
1478!> \brief Read a Gaussian-type orbital (GTO) basis set from the database file.
1479!> \param element_symbol ...
1480!> \param basis_type ...
1481!> \param gto_basis_set ...
1482!> \param basis_section ...
1483!> \param irep ...
1484!> \param dft_section ...
1485!> \author MK
1486! **************************************************************************************************
1487 SUBROUTINE read_gto_basis_set2(element_symbol, basis_type, gto_basis_set, &
1488 basis_section, irep, dft_section)
1489
1490 CHARACTER(LEN=*), INTENT(IN) :: element_symbol
1491 CHARACTER(LEN=*), INTENT(INOUT) :: basis_type
1492 TYPE(gto_basis_set_type), INTENT(INOUT) :: gto_basis_set
1493 TYPE(section_vals_type), OPTIONAL, POINTER :: basis_section
1494 INTEGER :: irep
1495 TYPE(section_vals_type), OPTIONAL, POINTER :: dft_section
1496
1497 CHARACTER(len=20*default_string_length) :: line_att
1498 CHARACTER(LEN=240) :: line
1499 CHARACTER(LEN=242) :: line2
1500 CHARACTER(LEN=default_path_length) :: bsname, bsname2
1501 CHARACTER(LEN=LEN(element_symbol)) :: symbol
1502 CHARACTER(LEN=LEN(element_symbol)+2) :: symbol2
1503 INTEGER :: i, ico, ipgf, iset, ishell, lshell, m, &
1504 maxco, maxl, maxpgf, maxshell, ncgf, &
1505 nmin, nset, nsgf, sort_method
1506 INTEGER, DIMENSION(:), POINTER :: lmax, lmin, npgf, nshell
1507 INTEGER, DIMENSION(:, :), POINTER :: l, n
1508 LOGICAL :: is_ok
1509 REAL(KIND=dp), DIMENSION(:, :), POINTER :: zet
1510 REAL(KIND=dp), DIMENSION(:, :, :), POINTER :: gcc
1511 TYPE(cp_sll_val_type), POINTER :: list
1512 TYPE(val_type), POINTER :: val
1513
1514 line = ""
1515 line2 = ""
1516 symbol = ""
1517 symbol2 = ""
1518 bsname = ""
1519 bsname2 = ""
1520
1521 gto_basis_set%name = " "
1522 gto_basis_set%aliases = " "
1523
1524 bsname = " "
1525 symbol = element_symbol
1526
1527 nset = 0
1528 maxshell = 0
1529 maxpgf = 0
1530 maxco = 0
1531 ncgf = 0
1532 nsgf = 0
1533 gto_basis_set%nset = nset
1534 gto_basis_set%ncgf = ncgf
1535 gto_basis_set%nsgf = nsgf
1536 CALL reallocate(gto_basis_set%lmax, 1, nset)
1537 CALL reallocate(gto_basis_set%lmin, 1, nset)
1538 CALL reallocate(gto_basis_set%npgf, 1, nset)
1539 CALL reallocate(gto_basis_set%nshell, 1, nset)
1540 CALL reallocate(gto_basis_set%n, 1, maxshell, 1, nset)
1541 CALL reallocate(gto_basis_set%l, 1, maxshell, 1, nset)
1542 CALL reallocate(gto_basis_set%zet, 1, maxpgf, 1, nset)
1543 CALL reallocate(gto_basis_set%gcc, 1, maxpgf, 1, maxshell, 1, nset)
1544 CALL reallocate(gto_basis_set%set_radius, 1, nset)
1545 CALL reallocate(gto_basis_set%pgf_radius, 1, maxpgf, 1, nset)
1546 CALL reallocate(gto_basis_set%first_cgf, 1, maxshell, 1, nset)
1547 CALL reallocate(gto_basis_set%first_sgf, 1, maxshell, 1, nset)
1548 CALL reallocate(gto_basis_set%last_cgf, 1, maxshell, 1, nset)
1549 CALL reallocate(gto_basis_set%last_sgf, 1, maxshell, 1, nset)
1550 CALL reallocate(gto_basis_set%ncgf_set, 1, nset)
1551 CALL reallocate(gto_basis_set%nsgf_set, 1, nset)
1552 CALL reallocate(gto_basis_set%cphi, 1, maxco, 1, ncgf)
1553 CALL reallocate(gto_basis_set%sphi, 1, maxco, 1, nsgf)
1554 CALL reallocate(gto_basis_set%scon, 1, maxco, 1, nsgf)
1555 CALL reallocate(gto_basis_set%lx, 1, ncgf)
1556 CALL reallocate(gto_basis_set%ly, 1, ncgf)
1557 CALL reallocate(gto_basis_set%lz, 1, ncgf)
1558 CALL reallocate(gto_basis_set%m, 1, nsgf)
1559 CALL reallocate(gto_basis_set%norm_cgf, 1, ncgf)
1560
1561 basis_type = ""
1562 CALL section_vals_val_get(basis_section, "_SECTION_PARAMETERS_", i_rep_section=irep, c_val=basis_type)
1563 IF (basis_type == "Orbital") basis_type = "ORB"
1564
1565 NULLIFY (list, val)
1566 CALL section_vals_list_get(basis_section, "_DEFAULT_KEYWORD_", i_rep_section=irep, list=list)
1567 CALL uppercase(symbol)
1568 CALL uppercase(bsname)
1569
1570 NULLIFY (gcc, l, lmax, lmin, n, npgf, nshell, zet)
1571 ! Read the basis set information
1572 is_ok = cp_sll_val_next(list, val)
1573 IF (.NOT. is_ok) cpabort("Error reading the Basis set from input file!")
1574 CALL val_get(val, c_val=line_att)
1575 READ (line_att, *) nset
1576
1577 CALL reallocate(npgf, 1, nset)
1578 CALL reallocate(nshell, 1, nset)
1579 CALL reallocate(lmax, 1, nset)
1580 CALL reallocate(lmin, 1, nset)
1581 CALL reallocate(n, 1, 1, 1, nset)
1582
1583 maxl = 0
1584 maxpgf = 0
1585 maxshell = 0
1586
1587 DO iset = 1, nset
1588 is_ok = cp_sll_val_next(list, val)
1589 IF (.NOT. is_ok) cpabort("Error reading the Basis set from input file!")
1590 CALL val_get(val, c_val=line_att)
1591 READ (line_att, *) n(1, iset)
1592 CALL remove_word(line_att)
1593 READ (line_att, *) lmin(iset)
1594 CALL remove_word(line_att)
1595 READ (line_att, *) lmax(iset)
1596 CALL remove_word(line_att)
1597 READ (line_att, *) npgf(iset)
1598 CALL remove_word(line_att)
1599 maxl = max(maxl, lmax(iset))
1600 IF (npgf(iset) > maxpgf) THEN
1601 maxpgf = npgf(iset)
1602 CALL reallocate(zet, 1, maxpgf, 1, nset)
1603 CALL reallocate(gcc, 1, maxpgf, 1, maxshell, 1, nset)
1604 END IF
1605 nshell(iset) = 0
1606 DO lshell = lmin(iset), lmax(iset)
1607 nmin = n(1, iset) + lshell - lmin(iset)
1608 READ (line_att, *) ishell
1609 CALL remove_word(line_att)
1610 nshell(iset) = nshell(iset) + ishell
1611 IF (nshell(iset) > maxshell) THEN
1612 maxshell = nshell(iset)
1613 CALL reallocate(n, 1, maxshell, 1, nset)
1614 CALL reallocate(l, 1, maxshell, 1, nset)
1615 CALL reallocate(gcc, 1, maxpgf, 1, maxshell, 1, nset)
1616 END IF
1617 DO i = 1, ishell
1618 n(nshell(iset) - ishell + i, iset) = nmin + i - 1
1619 l(nshell(iset) - ishell + i, iset) = lshell
1620 END DO
1621 END DO
1622 IF (len_trim(line_att) /= 0) &
1623 cpabort("Error reading the Basis from input file!")
1624 DO ipgf = 1, npgf(iset)
1625 is_ok = cp_sll_val_next(list, val)
1626 IF (.NOT. is_ok) cpabort("Error reading the Basis set from input file!")
1627 CALL val_get(val, c_val=line_att)
1628 READ (line_att, *) zet(ipgf, iset), (gcc(ipgf, ishell, iset), ishell=1, nshell(iset))
1629 END DO
1630 END DO
1631
1632 ! Maximum angular momentum quantum number of the atomic kind
1633
1634 CALL init_orbital_pointers(maxl)
1635
1636 ! Allocate the global variables
1637
1638 gto_basis_set%nset = nset
1639 CALL reallocate(gto_basis_set%lmax, 1, nset)
1640 CALL reallocate(gto_basis_set%lmin, 1, nset)
1641 CALL reallocate(gto_basis_set%npgf, 1, nset)
1642 CALL reallocate(gto_basis_set%nshell, 1, nset)
1643 CALL reallocate(gto_basis_set%n, 1, maxshell, 1, nset)
1644 CALL reallocate(gto_basis_set%l, 1, maxshell, 1, nset)
1645 CALL reallocate(gto_basis_set%zet, 1, maxpgf, 1, nset)
1646 CALL reallocate(gto_basis_set%gcc, 1, maxpgf, 1, maxshell, 1, nset)
1647
1648 ! Copy the basis set information into the data structure
1649
1650 DO iset = 1, nset
1651 gto_basis_set%lmax(iset) = lmax(iset)
1652 gto_basis_set%lmin(iset) = lmin(iset)
1653 gto_basis_set%npgf(iset) = npgf(iset)
1654 gto_basis_set%nshell(iset) = nshell(iset)
1655 DO ishell = 1, nshell(iset)
1656 gto_basis_set%n(ishell, iset) = n(ishell, iset)
1657 gto_basis_set%l(ishell, iset) = l(ishell, iset)
1658 DO ipgf = 1, npgf(iset)
1659 gto_basis_set%gcc(ipgf, ishell, iset) = gcc(ipgf, ishell, iset)
1660 END DO
1661 END DO
1662 DO ipgf = 1, npgf(iset)
1663 gto_basis_set%zet(ipgf, iset) = zet(ipgf, iset)
1664 END DO
1665 END DO
1666
1667 ! Initialise the depending atomic kind information
1668
1669 CALL reallocate(gto_basis_set%set_radius, 1, nset)
1670 CALL reallocate(gto_basis_set%pgf_radius, 1, maxpgf, 1, nset)
1671 CALL reallocate(gto_basis_set%first_cgf, 1, maxshell, 1, nset)
1672 CALL reallocate(gto_basis_set%first_sgf, 1, maxshell, 1, nset)
1673 CALL reallocate(gto_basis_set%last_cgf, 1, maxshell, 1, nset)
1674 CALL reallocate(gto_basis_set%last_sgf, 1, maxshell, 1, nset)
1675 CALL reallocate(gto_basis_set%ncgf_set, 1, nset)
1676 CALL reallocate(gto_basis_set%nsgf_set, 1, nset)
1677
1678 maxco = 0
1679 ncgf = 0
1680 nsgf = 0
1681
1682 DO iset = 1, nset
1683 gto_basis_set%ncgf_set(iset) = 0
1684 gto_basis_set%nsgf_set(iset) = 0
1685 DO ishell = 1, nshell(iset)
1686 lshell = gto_basis_set%l(ishell, iset)
1687 gto_basis_set%first_cgf(ishell, iset) = ncgf + 1
1688 ncgf = ncgf + nco(lshell)
1689 gto_basis_set%last_cgf(ishell, iset) = ncgf
1690 gto_basis_set%ncgf_set(iset) = &
1691 gto_basis_set%ncgf_set(iset) + nco(lshell)
1692 gto_basis_set%first_sgf(ishell, iset) = nsgf + 1
1693 nsgf = nsgf + nso(lshell)
1694 gto_basis_set%last_sgf(ishell, iset) = nsgf
1695 gto_basis_set%nsgf_set(iset) = &
1696 gto_basis_set%nsgf_set(iset) + nso(lshell)
1697 END DO
1698 maxco = max(maxco, npgf(iset)*ncoset(lmax(iset)))
1699 END DO
1700
1701 gto_basis_set%ncgf = ncgf
1702 gto_basis_set%nsgf = nsgf
1703
1704 CALL reallocate(gto_basis_set%cphi, 1, maxco, 1, ncgf)
1705 CALL reallocate(gto_basis_set%sphi, 1, maxco, 1, nsgf)
1706 CALL reallocate(gto_basis_set%scon, 1, maxco, 1, nsgf)
1707 CALL reallocate(gto_basis_set%lx, 1, ncgf)
1708 CALL reallocate(gto_basis_set%ly, 1, ncgf)
1709 CALL reallocate(gto_basis_set%lz, 1, ncgf)
1710 CALL reallocate(gto_basis_set%m, 1, nsgf)
1711 CALL reallocate(gto_basis_set%norm_cgf, 1, ncgf)
1712 ALLOCATE (gto_basis_set%cgf_symbol(ncgf))
1713
1714 ALLOCATE (gto_basis_set%sgf_symbol(nsgf))
1715
1716 ncgf = 0
1717 nsgf = 0
1718
1719 DO iset = 1, nset
1720 DO ishell = 1, nshell(iset)
1721 lshell = gto_basis_set%l(ishell, iset)
1722 DO ico = ncoset(lshell - 1) + 1, ncoset(lshell)
1723 ncgf = ncgf + 1
1724 gto_basis_set%lx(ncgf) = indco(1, ico)
1725 gto_basis_set%ly(ncgf) = indco(2, ico)
1726 gto_basis_set%lz(ncgf) = indco(3, ico)
1727 gto_basis_set%cgf_symbol(ncgf) = &
1728 cgf_symbol(n(ishell, iset), (/gto_basis_set%lx(ncgf), &
1729 gto_basis_set%ly(ncgf), &
1730 gto_basis_set%lz(ncgf)/))
1731 END DO
1732 DO m = -lshell, lshell
1733 nsgf = nsgf + 1
1734 gto_basis_set%m(nsgf) = m
1735 gto_basis_set%sgf_symbol(nsgf) = &
1736 sgf_symbol(n(ishell, iset), lshell, m)
1737 END DO
1738 END DO
1739 END DO
1740
1741 DEALLOCATE (gcc, l, lmax, lmin, n, npgf, nshell, zet)
1742
1743 IF (PRESENT(dft_section)) THEN
1744 CALL section_vals_val_get(dft_section, "SORT_BASIS", i_val=sort_method)
1745 CALL sort_gto_basis_set(gto_basis_set, sort_method)
1746 END IF
1747
1748 END SUBROUTINE read_gto_basis_set2
1749
1750! **************************************************************************************************
1751!> \brief Set the components of Gaussian-type orbital (GTO) basis set data set.
1752!> \param gto_basis_set ...
1753!> \param name ...
1754!> \param aliases ...
1755!> \param norm_type ...
1756!> \param kind_radius ...
1757!> \param ncgf ...
1758!> \param nset ...
1759!> \param nsgf ...
1760!> \param cgf_symbol ...
1761!> \param sgf_symbol ...
1762!> \param norm_cgf ...
1763!> \param set_radius ...
1764!> \param lmax ...
1765!> \param lmin ...
1766!> \param lx ...
1767!> \param ly ...
1768!> \param lz ...
1769!> \param m ...
1770!> \param ncgf_set ...
1771!> \param npgf ...
1772!> \param nsgf_set ...
1773!> \param nshell ...
1774!> \param cphi ...
1775!> \param pgf_radius ...
1776!> \param sphi ...
1777!> \param scon ...
1778!> \param zet ...
1779!> \param first_cgf ...
1780!> \param first_sgf ...
1781!> \param l ...
1782!> \param last_cgf ...
1783!> \param last_sgf ...
1784!> \param n ...
1785!> \param gcc ...
1786!> \param short_kind_radius ...
1787!> \author MK
1788! **************************************************************************************************
1789 SUBROUTINE set_gto_basis_set(gto_basis_set, name, aliases, norm_type, kind_radius, ncgf, &
1790 nset, nsgf, cgf_symbol, sgf_symbol, norm_cgf, set_radius, lmax, &
1791 lmin, lx, ly, lz, m, ncgf_set, npgf, nsgf_set, nshell, &
1792 cphi, pgf_radius, sphi, scon, zet, first_cgf, first_sgf, l, &
1793 last_cgf, last_sgf, n, gcc, short_kind_radius)
1794
1795 TYPE(gto_basis_set_type), INTENT(INOUT) :: gto_basis_set
1796 CHARACTER(LEN=default_string_length), INTENT(IN), &
1797 OPTIONAL :: name, aliases
1798 INTEGER, INTENT(IN), OPTIONAL :: norm_type
1799 REAL(kind=dp), INTENT(IN), OPTIONAL :: kind_radius
1800 INTEGER, INTENT(IN), OPTIONAL :: ncgf, nset, nsgf
1801 CHARACTER(LEN=12), DIMENSION(:), OPTIONAL, POINTER :: cgf_symbol
1802 CHARACTER(LEN=6), DIMENSION(:), OPTIONAL, POINTER :: sgf_symbol
1803 REAL(kind=dp), DIMENSION(:), OPTIONAL, POINTER :: norm_cgf, set_radius
1804 INTEGER, DIMENSION(:), OPTIONAL, POINTER :: lmax, lmin, lx, ly, lz, m, ncgf_set, &
1805 npgf, nsgf_set, nshell
1806 REAL(kind=dp), DIMENSION(:, :), OPTIONAL, POINTER :: cphi, pgf_radius, sphi, scon, zet
1807 INTEGER, DIMENSION(:, :), OPTIONAL, POINTER :: first_cgf, first_sgf, l, last_cgf, &
1808 last_sgf, n
1809 REAL(kind=dp), DIMENSION(:, :, :), OPTIONAL, &
1810 POINTER :: gcc
1811 REAL(kind=dp), INTENT(IN), OPTIONAL :: short_kind_radius
1812
1813 IF (PRESENT(name)) gto_basis_set%name = name
1814 IF (PRESENT(aliases)) gto_basis_set%aliases = aliases
1815 IF (PRESENT(norm_type)) gto_basis_set%norm_type = norm_type
1816 IF (PRESENT(kind_radius)) gto_basis_set%kind_radius = kind_radius
1817 IF (PRESENT(short_kind_radius)) gto_basis_set%short_kind_radius = short_kind_radius
1818 IF (PRESENT(ncgf)) gto_basis_set%ncgf = ncgf
1819 IF (PRESENT(nset)) gto_basis_set%nset = nset
1820 IF (PRESENT(nsgf)) gto_basis_set%nsgf = nsgf
1821 IF (PRESENT(cgf_symbol)) gto_basis_set%cgf_symbol(:) = cgf_symbol(:)
1822 IF (PRESENT(sgf_symbol)) gto_basis_set%sgf_symbol(:) = sgf_symbol(:)
1823 IF (PRESENT(norm_cgf)) gto_basis_set%norm_cgf(:) = norm_cgf(:)
1824 IF (PRESENT(set_radius)) gto_basis_set%set_radius(:) = set_radius(:)
1825 IF (PRESENT(lmax)) gto_basis_set%lmax(:) = lmax(:)
1826 IF (PRESENT(lmin)) gto_basis_set%lmin(:) = lmin(:)
1827 IF (PRESENT(lx)) gto_basis_set%lx(:) = lx(:)
1828 IF (PRESENT(ly)) gto_basis_set%ly(:) = ly(:)
1829 IF (PRESENT(lz)) gto_basis_set%lz(:) = lz(:)
1830 IF (PRESENT(m)) gto_basis_set%m(:) = m(:)
1831 IF (PRESENT(ncgf_set)) gto_basis_set%ncgf_set(:) = ncgf_set(:)
1832 IF (PRESENT(npgf)) gto_basis_set%npgf(:) = npgf(:)
1833 IF (PRESENT(nsgf_set)) gto_basis_set%nsgf_set(:) = nsgf_set(:)
1834 IF (PRESENT(nshell)) gto_basis_set%nshell(:) = nshell(:)
1835 IF (PRESENT(cphi)) gto_basis_set%cphi(:, :) = cphi(:, :)
1836 IF (PRESENT(pgf_radius)) gto_basis_set%pgf_radius(:, :) = pgf_radius(:, :)
1837 IF (PRESENT(sphi)) gto_basis_set%sphi(:, :) = sphi(:, :)
1838 IF (PRESENT(scon)) gto_basis_set%scon(:, :) = scon(:, :)
1839 IF (PRESENT(zet)) gto_basis_set%zet(:, :) = zet(:, :)
1840 IF (PRESENT(first_cgf)) gto_basis_set%first_cgf(:, :) = first_cgf(:, :)
1841 IF (PRESENT(first_sgf)) gto_basis_set%first_sgf(:, :) = first_sgf(:, :)
1842 IF (PRESENT(l)) l(:, :) = gto_basis_set%l(:, :)
1843 IF (PRESENT(last_cgf)) gto_basis_set%last_cgf(:, :) = last_cgf(:, :)
1844 IF (PRESENT(last_sgf)) gto_basis_set%last_sgf(:, :) = last_sgf(:, :)
1845 IF (PRESENT(n)) gto_basis_set%n(:, :) = n(:, :)
1846 IF (PRESENT(gcc)) gto_basis_set%gcc(:, :, :) = gcc(:, :, :)
1847
1848 END SUBROUTINE set_gto_basis_set
1849
1850! **************************************************************************************************
1851!> \brief Write a Gaussian-type orbital (GTO) basis set data set to the output unit.
1852!> \param gto_basis_set ...
1853!> \param output_unit ...
1854!> \param header ...
1855!> \author MK
1856! **************************************************************************************************
1857 SUBROUTINE write_gto_basis_set(gto_basis_set, output_unit, header)
1858
1859 TYPE(gto_basis_set_type), INTENT(IN) :: gto_basis_set
1860 INTEGER, INTENT(in) :: output_unit
1861 CHARACTER(len=*), OPTIONAL :: header
1862
1863 INTEGER :: ipgf, iset, ishell
1864
1865 IF (output_unit > 0) THEN
1866
1867 IF (PRESENT(header)) THEN
1868 WRITE (unit=output_unit, fmt="(/,T6,A,T41,A40)") &
1869 trim(header), trim(gto_basis_set%name)
1870 END IF
1871
1872 WRITE (unit=output_unit, fmt="(/,(T8,A,T71,I10))") &
1873 "Number of orbital shell sets: ", &
1874 gto_basis_set%nset, &
1875 "Number of orbital shells: ", &
1876 sum(gto_basis_set%nshell(:)), &
1877 "Number of primitive Cartesian functions: ", &
1878 sum(gto_basis_set%npgf(:)), &
1879 "Number of Cartesian basis functions: ", &
1880 gto_basis_set%ncgf, &
1881 "Number of spherical basis functions: ", &
1882 gto_basis_set%nsgf, &
1883 "Norm type: ", &
1884 gto_basis_set%norm_type
1885
1886 WRITE (unit=output_unit, fmt="(/,T6,A,T41,A40,/,/,T25,A)") &
1887 "GTO basis set information for", trim(gto_basis_set%name), &
1888 "Set Shell n l Exponent Coefficient"
1889
1890 DO iset = 1, gto_basis_set%nset
1891 WRITE (unit=output_unit, fmt="(A)") ""
1892 DO ishell = 1, gto_basis_set%nshell(iset)
1893 WRITE (unit=output_unit, &
1894 fmt="(T25,I3,4X,I4,4X,I2,2X,I2,(T51,2F15.6))") &
1895 iset, ishell, &
1896 gto_basis_set%n(ishell, iset), &
1897 gto_basis_set%l(ishell, iset), &
1898 (gto_basis_set%zet(ipgf, iset), &
1899 gto_basis_set%gcc(ipgf, ishell, iset), &
1900 ipgf=1, gto_basis_set%npgf(iset))
1901 END DO
1902 END DO
1903
1904 END IF
1905
1906 END SUBROUTINE write_gto_basis_set
1907
1908! **************************************************************************************************
1909
1910! **************************************************************************************************
1911!> \brief Write a Gaussian-type orbital (GTO) basis set data set to the output unit.
1912!> \param orb_basis_set ...
1913!> \param output_unit ...
1914!> \param header ...
1915!> \author MK
1916! **************************************************************************************************
1917 SUBROUTINE write_orb_basis_set(orb_basis_set, output_unit, header)
1918
1919 TYPE(gto_basis_set_type), INTENT(IN) :: orb_basis_set
1920 INTEGER, INTENT(in) :: output_unit
1921 CHARACTER(len=*), OPTIONAL :: header
1922
1923 INTEGER :: icgf, ico, ipgf, iset, ishell
1924
1925 IF (output_unit > 0) THEN
1926 IF (PRESENT(header)) THEN
1927 WRITE (unit=output_unit, fmt="(/,T6,A,T41,A40)") &
1928 trim(header), trim(orb_basis_set%name)
1929 END IF
1930
1931 WRITE (unit=output_unit, fmt="(/,(T8,A,T71,I10))") &
1932 "Number of orbital shell sets: ", &
1933 orb_basis_set%nset, &
1934 "Number of orbital shells: ", &
1935 sum(orb_basis_set%nshell(:)), &
1936 "Number of primitive Cartesian functions: ", &
1937 sum(orb_basis_set%npgf(:)), &
1938 "Number of Cartesian basis functions: ", &
1939 orb_basis_set%ncgf, &
1940 "Number of spherical basis functions: ", &
1941 orb_basis_set%nsgf, &
1942 "Norm type: ", &
1943 orb_basis_set%norm_type
1944
1945 WRITE (unit=output_unit, fmt="(/,T8,A,/,/,T25,A)") &
1946 "Normalised Cartesian orbitals:", &
1947 "Set Shell Orbital Exponent Coefficient"
1948
1949 icgf = 0
1950
1951 DO iset = 1, orb_basis_set%nset
1952 DO ishell = 1, orb_basis_set%nshell(iset)
1953 WRITE (unit=output_unit, fmt="(A)") ""
1954 DO ico = 1, nco(orb_basis_set%l(ishell, iset))
1955 icgf = icgf + 1
1956 WRITE (unit=output_unit, &
1957 fmt="(T25,I3,4X,I4,3X,A12,(T51,2F15.6))") &
1958 iset, ishell, orb_basis_set%cgf_symbol(icgf), &
1959 (orb_basis_set%zet(ipgf, iset), &
1960 orb_basis_set%norm_cgf(icgf)* &
1961 orb_basis_set%gcc(ipgf, ishell, iset), &
1962 ipgf=1, orb_basis_set%npgf(iset))
1963 END DO
1964 END DO
1965 END DO
1966 END IF
1967
1968 END SUBROUTINE write_orb_basis_set
1969
1970! **************************************************************************************************
1971!> \brief ...
1972!> \param sto_basis_set ...
1973! **************************************************************************************************
1974 SUBROUTINE allocate_sto_basis_set(sto_basis_set)
1975
1976 TYPE(sto_basis_set_type), POINTER :: sto_basis_set
1977
1978 CALL deallocate_sto_basis_set(sto_basis_set)
1979
1980 ALLOCATE (sto_basis_set)
1981
1982 END SUBROUTINE allocate_sto_basis_set
1983
1984! **************************************************************************************************
1985!> \brief ...
1986!> \param sto_basis_set ...
1987! **************************************************************************************************
1988 SUBROUTINE deallocate_sto_basis_set(sto_basis_set)
1989
1990 TYPE(sto_basis_set_type), POINTER :: sto_basis_set
1991
1992! -------------------------------------------------------------------------
1993
1994 IF (ASSOCIATED(sto_basis_set)) THEN
1995 IF (ASSOCIATED(sto_basis_set%symbol)) THEN
1996 DEALLOCATE (sto_basis_set%symbol)
1997 END IF
1998 IF (ASSOCIATED(sto_basis_set%nq)) THEN
1999 DEALLOCATE (sto_basis_set%nq)
2000 END IF
2001 IF (ASSOCIATED(sto_basis_set%lq)) THEN
2002 DEALLOCATE (sto_basis_set%lq)
2003 END IF
2004 IF (ASSOCIATED(sto_basis_set%zet)) THEN
2005 DEALLOCATE (sto_basis_set%zet)
2006 END IF
2007
2008 DEALLOCATE (sto_basis_set)
2009 END IF
2010 END SUBROUTINE deallocate_sto_basis_set
2011
2012! **************************************************************************************************
2013!> \brief ...
2014!> \param sto_basis_set ...
2015!> \param name ...
2016!> \param nshell ...
2017!> \param symbol ...
2018!> \param nq ...
2019!> \param lq ...
2020!> \param zet ...
2021!> \param maxlq ...
2022!> \param numsto ...
2023! **************************************************************************************************
2024 SUBROUTINE get_sto_basis_set(sto_basis_set, name, nshell, symbol, nq, lq, zet, maxlq, numsto)
2025
2026 TYPE(sto_basis_set_type), INTENT(IN) :: sto_basis_set
2027 CHARACTER(LEN=default_string_length), &
2028 INTENT(OUT), OPTIONAL :: name
2029 INTEGER, INTENT(OUT), OPTIONAL :: nshell
2030 CHARACTER(LEN=6), DIMENSION(:), OPTIONAL, POINTER :: symbol
2031 INTEGER, DIMENSION(:), OPTIONAL, POINTER :: nq, lq
2032 REAL(kind=dp), DIMENSION(:), OPTIONAL, POINTER :: zet
2033 INTEGER, INTENT(OUT), OPTIONAL :: maxlq, numsto
2034
2035 INTEGER :: iset
2036
2037 IF (PRESENT(name)) name = sto_basis_set%name
2038 IF (PRESENT(nshell)) nshell = sto_basis_set%nshell
2039 IF (PRESENT(symbol)) symbol => sto_basis_set%symbol
2040 IF (PRESENT(nq)) nq => sto_basis_set%nq
2041 IF (PRESENT(lq)) lq => sto_basis_set%lq
2042 IF (PRESENT(zet)) zet => sto_basis_set%zet
2043 IF (PRESENT(maxlq)) THEN
2044 maxlq = maxval(sto_basis_set%lq(1:sto_basis_set%nshell))
2045 END IF
2046 IF (PRESENT(numsto)) THEN
2047 numsto = 0
2048 DO iset = 1, sto_basis_set%nshell
2049 numsto = numsto + 2*sto_basis_set%lq(iset) + 1
2050 END DO
2051 END IF
2052
2053 END SUBROUTINE get_sto_basis_set
2054
2055! **************************************************************************************************
2056!> \brief ...
2057!> \param sto_basis_set ...
2058!> \param name ...
2059!> \param nshell ...
2060!> \param symbol ...
2061!> \param nq ...
2062!> \param lq ...
2063!> \param zet ...
2064! **************************************************************************************************
2065 SUBROUTINE set_sto_basis_set(sto_basis_set, name, nshell, symbol, nq, lq, zet)
2066
2067 TYPE(sto_basis_set_type), INTENT(INOUT) :: sto_basis_set
2068 CHARACTER(LEN=default_string_length), INTENT(IN), &
2069 OPTIONAL :: name
2070 INTEGER, INTENT(IN), OPTIONAL :: nshell
2071 CHARACTER(LEN=6), DIMENSION(:), OPTIONAL, POINTER :: symbol
2072 INTEGER, DIMENSION(:), OPTIONAL, POINTER :: nq, lq
2073 REAL(kind=dp), DIMENSION(:), OPTIONAL, POINTER :: zet
2074
2075 INTEGER :: ns
2076
2077 IF (PRESENT(name)) sto_basis_set%name = name
2078 IF (PRESENT(nshell)) sto_basis_set%nshell = nshell
2079 IF (PRESENT(symbol)) THEN
2080 ns = SIZE(symbol)
2081 IF (ASSOCIATED(sto_basis_set%symbol)) DEALLOCATE (sto_basis_set%symbol)
2082 ALLOCATE (sto_basis_set%symbol(1:ns))
2083 sto_basis_set%symbol(:) = symbol(:)
2084 END IF
2085 IF (PRESENT(nq)) THEN
2086 ns = SIZE(nq)
2087 CALL reallocate(sto_basis_set%nq, 1, ns)
2088 sto_basis_set%nq = nq(:)
2089 END IF
2090 IF (PRESENT(lq)) THEN
2091 ns = SIZE(lq)
2092 CALL reallocate(sto_basis_set%lq, 1, ns)
2093 sto_basis_set%lq = lq(:)
2094 END IF
2095 IF (PRESENT(zet)) THEN
2096 ns = SIZE(zet)
2097 CALL reallocate(sto_basis_set%zet, 1, ns)
2098 sto_basis_set%zet = zet(:)
2099 END IF
2100
2101 END SUBROUTINE set_sto_basis_set
2102
2103! **************************************************************************************************
2104!> \brief ...
2105!> \param element_symbol ...
2106!> \param basis_set_name ...
2107!> \param sto_basis_set ...
2108!> \param para_env ...
2109!> \param dft_section ...
2110! **************************************************************************************************
2111 SUBROUTINE read_sto_basis_set(element_symbol, basis_set_name, sto_basis_set, para_env, dft_section)
2112
2113 ! Read a Slater-type orbital (STO) basis set from the database file.
2114
2115 CHARACTER(LEN=*), INTENT(IN) :: element_symbol, basis_set_name
2116 TYPE(sto_basis_set_type), INTENT(INOUT) :: sto_basis_set
2117 TYPE(mp_para_env_type), POINTER :: para_env
2118 TYPE(section_vals_type), POINTER :: dft_section
2119
2120 CHARACTER(LEN=10) :: nlsym
2121 CHARACTER(LEN=2) :: lsym
2122 CHARACTER(LEN=240) :: line
2123 CHARACTER(LEN=242) :: line2
2124 CHARACTER(len=default_path_length) :: basis_set_file_name, tmp
2125 CHARACTER(LEN=default_path_length), DIMENSION(:), &
2126 POINTER :: cbasis
2127 CHARACTER(LEN=LEN(basis_set_name)) :: bsname
2128 CHARACTER(LEN=LEN(basis_set_name)+2) :: bsname2
2129 CHARACTER(LEN=LEN(element_symbol)) :: symbol
2130 CHARACTER(LEN=LEN(element_symbol)+2) :: symbol2
2131 INTEGER :: ibasis, irep, iset, nbasis, nq, nset, &
2132 strlen1, strlen2
2133 LOGICAL :: basis_found, found, match
2134 REAL(kind=dp) :: zet
2135 TYPE(cp_parser_type) :: parser
2136
2137 line = ""
2138 line2 = ""
2139 symbol = ""
2140 symbol2 = ""
2141 bsname = ""
2142 bsname2 = ""
2143
2144 nbasis = 1
2145
2146 sto_basis_set%name = basis_set_name
2147 CALL section_vals_val_get(dft_section, "BASIS_SET_FILE_NAME", &
2148 n_rep_val=nbasis)
2149 ALLOCATE (cbasis(nbasis))
2150 DO ibasis = 1, nbasis
2151 CALL section_vals_val_get(dft_section, "BASIS_SET_FILE_NAME", &
2152 i_rep_val=ibasis, c_val=cbasis(ibasis))
2153 basis_set_file_name = cbasis(ibasis)
2154 tmp = basis_set_file_name
2155 CALL uppercase(tmp)
2156 END DO
2157
2158 ! Search for the requested basis set in the basis set file
2159 ! until the basis set is found or the end of file is reached
2160
2161 basis_found = .false.
2162 basis_loop: DO ibasis = 1, nbasis
2163 IF (basis_found) EXIT basis_loop
2164 basis_set_file_name = cbasis(ibasis)
2165 CALL parser_create(parser, basis_set_file_name, para_env=para_env)
2166
2167 bsname = basis_set_name
2168 symbol = element_symbol
2169 irep = 0
2170
2171 tmp = basis_set_name
2172 CALL uppercase(tmp)
2173
2174 IF (tmp .NE. "NONE") THEN
2175 search_loop: DO
2176
2177 CALL parser_search_string(parser, trim(bsname), .true., found, line)
2178 IF (found) THEN
2179 CALL uppercase(symbol)
2180 CALL uppercase(bsname)
2181
2182 match = .false.
2183 CALL uppercase(line)
2184 ! Check both the element symbol and the basis set name
2185 line2 = " "//line//" "
2186 symbol2 = " "//trim(symbol)//" "
2187 bsname2 = " "//trim(bsname)//" "
2188 strlen1 = len_trim(symbol2) + 1
2189 strlen2 = len_trim(bsname2) + 1
2190
2191 IF ((index(line2, symbol2(:strlen1)) > 0) .AND. &
2192 (index(line2, bsname2(:strlen2)) > 0)) match = .true.
2193 IF (match) THEN
2194 ! Read the basis set information
2195 CALL parser_get_object(parser, nset, newline=.true.)
2196 sto_basis_set%nshell = nset
2197
2198 CALL reallocate(sto_basis_set%nq, 1, nset)
2199 CALL reallocate(sto_basis_set%lq, 1, nset)
2200 CALL reallocate(sto_basis_set%zet, 1, nset)
2201 ALLOCATE (sto_basis_set%symbol(nset))
2202
2203 DO iset = 1, nset
2204 CALL parser_get_object(parser, nq, newline=.true.)
2205 CALL parser_get_object(parser, lsym)
2206 CALL parser_get_object(parser, zet)
2207 sto_basis_set%nq(iset) = nq
2208 sto_basis_set%zet(iset) = zet
2209 WRITE (nlsym, "(I2,A)") nq, trim(lsym)
2210 sto_basis_set%symbol(iset) = trim(nlsym)
2211 SELECT CASE (trim(lsym))
2212 CASE ("S", "s")
2213 sto_basis_set%lq(iset) = 0
2214 CASE ("P", "p")
2215 sto_basis_set%lq(iset) = 1
2216 CASE ("D", "d")
2217 sto_basis_set%lq(iset) = 2
2218 CASE ("F", "f")
2219 sto_basis_set%lq(iset) = 3
2220 CASE ("G", "g")
2221 sto_basis_set%lq(iset) = 4
2222 CASE ("H", "h")
2223 sto_basis_set%lq(iset) = 5
2224 CASE ("I", "i", "J", "j")
2225 sto_basis_set%lq(iset) = 6
2226 CASE ("K", "k")
2227 sto_basis_set%lq(iset) = 7
2228 CASE ("L", "l")
2229 sto_basis_set%lq(iset) = 8
2230 CASE ("M", "m")
2231 sto_basis_set%lq(iset) = 9
2232 CASE DEFAULT
2233 CALL cp_abort(__location__, &
2234 "The requested basis set <"//trim(bsname)// &
2235 "> for element <"//trim(symbol)//"> has an invalid component: ")
2236 END SELECT
2237 END DO
2238
2239 basis_found = .true.
2240 EXIT search_loop
2241 END IF
2242 ELSE
2243 EXIT search_loop
2244 END IF
2245 END DO search_loop
2246 ELSE
2247 match = .false.
2248 END IF
2249
2250 CALL parser_release(parser)
2251
2252 END DO basis_loop
2253
2254 IF (tmp .NE. "NONE") THEN
2255 IF (.NOT. basis_found) THEN
2256 basis_set_file_name = ""
2257 DO ibasis = 1, nbasis
2258 basis_set_file_name = trim(basis_set_file_name)//"<"//trim(cbasis(ibasis))//"> "
2259 END DO
2260 CALL cp_abort(__location__, &
2261 "The requested basis set <"//trim(bsname)// &
2262 "> for element <"//trim(symbol)//"> was not "// &
2263 "found in the basis set files "// &
2264 trim(basis_set_file_name))
2265 END IF
2266 END IF
2267 DEALLOCATE (cbasis)
2268
2269 END SUBROUTINE read_sto_basis_set
2270
2271! **************************************************************************************************
2272!> \brief ...
2273!> \param sto_basis_set ...
2274!> \param gto_basis_set ...
2275!> \param ngauss ...
2276!> \param ortho ...
2277! **************************************************************************************************
2278 SUBROUTINE create_gto_from_sto_basis(sto_basis_set, gto_basis_set, ngauss, ortho)
2279
2280 TYPE(sto_basis_set_type), INTENT(IN) :: sto_basis_set
2281 TYPE(gto_basis_set_type), POINTER :: gto_basis_set
2282 INTEGER, INTENT(IN), OPTIONAL :: ngauss
2283 LOGICAL, INTENT(IN), OPTIONAL :: ortho
2284
2285 INTEGER, PARAMETER :: maxng = 6
2286
2287 CHARACTER(LEN=default_string_length) :: name, sng
2288 INTEGER :: i1, i2, ico, ipgf, iset, jset, l, &
2289 lshell, m, maxco, maxl, ncgf, ng, ngs, &
2290 np, nset, nsgf, nshell
2291 INTEGER, DIMENSION(0:10) :: mxf
2292 INTEGER, DIMENSION(:), POINTER :: lq, nq
2293 LOGICAL :: do_ortho
2294 REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :) :: gal, zal, zll
2295 REAL(kind=dp), DIMENSION(:), POINTER :: zet
2296 REAL(kind=dp), DIMENSION(maxng) :: gcc, zetg
2297
2298 ng = 6
2299 IF (PRESENT(ngauss)) ng = ngauss
2300 IF (ng > maxng) cpabort("Too many Gaussian primitives requested")
2301 do_ortho = .false.
2302 IF (PRESENT(ortho)) do_ortho = ortho
2303
2304 CALL allocate_gto_basis_set(gto_basis_set)
2305
2306 CALL get_sto_basis_set(sto_basis_set, name=name, nshell=nshell, nq=nq, &
2307 lq=lq, zet=zet)
2308
2309 maxl = maxval(lq)
2310 CALL init_orbital_pointers(maxl)
2311
2312 CALL integer_to_string(ng, sng)
2313 gto_basis_set%name = trim(name)//"_STO-"//trim(sng)//"G"
2314
2315 nset = nshell
2316 gto_basis_set%nset = nset
2317 CALL reallocate(gto_basis_set%lmax, 1, nset)
2318 CALL reallocate(gto_basis_set%lmin, 1, nset)
2319 CALL reallocate(gto_basis_set%npgf, 1, nset)
2320 CALL reallocate(gto_basis_set%nshell, 1, nset)
2321 CALL reallocate(gto_basis_set%n, 1, 1, 1, nset)
2322 CALL reallocate(gto_basis_set%l, 1, 1, 1, nset)
2323 CALL reallocate(gto_basis_set%zet, 1, ng, 1, nset)
2324 CALL reallocate(gto_basis_set%gcc, 1, ng, 1, 1, 1, nset)
2325
2326 DO iset = 1, nset
2327 CALL get_sto_ng(zet(iset), ng, nq(iset), lq(iset), zetg, gcc)
2328 gto_basis_set%lmax(iset) = lq(iset)
2329 gto_basis_set%lmin(iset) = lq(iset)
2330 gto_basis_set%npgf(iset) = ng
2331 gto_basis_set%nshell(iset) = 1
2332 gto_basis_set%n(1, iset) = lq(iset) + 1
2333 gto_basis_set%l(1, iset) = lq(iset)
2334 DO ipgf = 1, ng
2335 gto_basis_set%gcc(ipgf, 1, iset) = gcc(ipgf)
2336 gto_basis_set%zet(ipgf, iset) = zetg(ipgf)
2337 END DO
2338 END DO
2339
2340 IF (do_ortho) THEN
2341 mxf = 0
2342 DO iset = 1, nset
2343 l = gto_basis_set%l(1, iset)
2344 mxf(l) = mxf(l) + 1
2345 END DO
2346 m = maxval(mxf)
2347 IF (m > 1) THEN
2348 ALLOCATE (gal(ng, nset), zal(ng, nset), zll(m*ng, 0:maxl))
2349 DO iset = 1, nset
2350 zal(1:ng, iset) = gto_basis_set%zet(1:ng, iset)
2351 gal(1:ng, iset) = gto_basis_set%gcc(1:ng, 1, iset)
2352 END DO
2353 CALL reallocate(gto_basis_set%zet, 1, m*ng, 1, nset)
2354 CALL reallocate(gto_basis_set%gcc, 1, m*ng, 1, 1, 1, nset)
2355 DO iset = 1, nset
2356 l = gto_basis_set%l(1, iset)
2357 gto_basis_set%npgf(iset) = ng*mxf(l)
2358 END DO
2359 gto_basis_set%zet = 0.0_dp
2360 gto_basis_set%gcc = 0.0_dp
2361 zll = 0.0_dp
2362 mxf = 0
2363 DO iset = 1, nset
2364 l = gto_basis_set%l(1, iset)
2365 mxf(l) = mxf(l) + 1
2366 i1 = mxf(l)*ng - ng + 1
2367 i2 = mxf(l)*ng
2368 zll(i1:i2, l) = zal(1:ng, iset)
2369 gto_basis_set%gcc(i1:i2, 1, iset) = gal(1:ng, iset)
2370 END DO
2371 DO iset = 1, nset
2372 l = gto_basis_set%l(1, iset)
2373 gto_basis_set%zet(:, iset) = zll(:, l)
2374 END DO
2375 DO iset = 1, nset
2376 l = gto_basis_set%l(1, iset)
2377 DO jset = 1, iset - 1
2378 IF (gto_basis_set%l(1, iset) == l) THEN
2379 m = mxf(l)*ng
2380 CALL orthofun(gto_basis_set%zet(1:m, iset), gto_basis_set%gcc(1:m, 1, iset), &
2381 gto_basis_set%gcc(1:m, 1, jset), l)
2382 END IF
2383 END DO
2384 END DO
2385 DEALLOCATE (gal, zal, zll)
2386 END IF
2387 END IF
2388
2389 ngs = maxval(gto_basis_set%npgf(1:nset))
2390 CALL reallocate(gto_basis_set%set_radius, 1, nset)
2391 CALL reallocate(gto_basis_set%pgf_radius, 1, ngs, 1, nset)
2392 CALL reallocate(gto_basis_set%first_cgf, 1, 1, 1, nset)
2393 CALL reallocate(gto_basis_set%first_sgf, 1, 1, 1, nset)
2394 CALL reallocate(gto_basis_set%last_cgf, 1, 1, 1, nset)
2395 CALL reallocate(gto_basis_set%last_sgf, 1, 1, 1, nset)
2396 CALL reallocate(gto_basis_set%ncgf_set, 1, nset)
2397 CALL reallocate(gto_basis_set%nsgf_set, 1, nset)
2398
2399 maxco = 0
2400 ncgf = 0
2401 nsgf = 0
2402
2403 DO iset = 1, nset
2404 gto_basis_set%ncgf_set(iset) = 0
2405 gto_basis_set%nsgf_set(iset) = 0
2406 lshell = gto_basis_set%l(1, iset)
2407 gto_basis_set%first_cgf(1, iset) = ncgf + 1
2408 ncgf = ncgf + nco(lshell)
2409 gto_basis_set%last_cgf(1, iset) = ncgf
2410 gto_basis_set%ncgf_set(iset) = &
2411 gto_basis_set%ncgf_set(iset) + nco(lshell)
2412 gto_basis_set%first_sgf(1, iset) = nsgf + 1
2413 nsgf = nsgf + nso(lshell)
2414 gto_basis_set%last_sgf(1, iset) = nsgf
2415 gto_basis_set%nsgf_set(iset) = &
2416 gto_basis_set%nsgf_set(iset) + nso(lshell)
2417 ngs = gto_basis_set%npgf(iset)
2418 maxco = max(maxco, ngs*ncoset(lshell))
2419 END DO
2420
2421 gto_basis_set%ncgf = ncgf
2422 gto_basis_set%nsgf = nsgf
2423
2424 CALL reallocate(gto_basis_set%cphi, 1, maxco, 1, ncgf)
2425 CALL reallocate(gto_basis_set%sphi, 1, maxco, 1, nsgf)
2426 CALL reallocate(gto_basis_set%scon, 1, maxco, 1, nsgf)
2427 CALL reallocate(gto_basis_set%lx, 1, ncgf)
2428 CALL reallocate(gto_basis_set%ly, 1, ncgf)
2429 CALL reallocate(gto_basis_set%lz, 1, ncgf)
2430 CALL reallocate(gto_basis_set%m, 1, nsgf)
2431 CALL reallocate(gto_basis_set%norm_cgf, 1, ncgf)
2432 ALLOCATE (gto_basis_set%cgf_symbol(ncgf))
2433 ALLOCATE (gto_basis_set%sgf_symbol(nsgf))
2434
2435 ncgf = 0
2436 nsgf = 0
2437
2438 DO iset = 1, nset
2439 lshell = gto_basis_set%l(1, iset)
2440 np = lshell + 1
2441 DO ico = ncoset(lshell - 1) + 1, ncoset(lshell)
2442 ncgf = ncgf + 1
2443 gto_basis_set%lx(ncgf) = indco(1, ico)
2444 gto_basis_set%ly(ncgf) = indco(2, ico)
2445 gto_basis_set%lz(ncgf) = indco(3, ico)
2446 gto_basis_set%cgf_symbol(ncgf) = &
2447 cgf_symbol(np, (/gto_basis_set%lx(ncgf), &
2448 gto_basis_set%ly(ncgf), &
2449 gto_basis_set%lz(ncgf)/))
2450 END DO
2451 DO m = -lshell, lshell
2452 nsgf = nsgf + 1
2453 gto_basis_set%m(nsgf) = m
2454 gto_basis_set%sgf_symbol(nsgf) = sgf_symbol(np, lshell, m)
2455 END DO
2456 END DO
2457
2458 gto_basis_set%norm_type = -1
2459
2460 END SUBROUTINE create_gto_from_sto_basis
2461
2462! **************************************************************************************************
2463!> \brief ...
2464!> \param zet ...
2465!> \param co ...
2466!> \param cr ...
2467!> \param l ...
2468! **************************************************************************************************
2469 SUBROUTINE orthofun(zet, co, cr, l)
2470 REAL(kind=dp), DIMENSION(:), INTENT(IN) :: zet
2471 REAL(kind=dp), DIMENSION(:), INTENT(INOUT) :: co, cr
2472 INTEGER, INTENT(IN) :: l
2473
2474 REAL(kind=dp) :: ss
2475
2476 CALL aovlp(l, zet, cr, cr, ss)
2477 cr(:) = cr(:)/sqrt(ss)
2478 CALL aovlp(l, zet, co, cr, ss)
2479 co(:) = co(:) - ss*cr(:)
2480 CALL aovlp(l, zet, co, co, ss)
2481 co(:) = co(:)/sqrt(ss)
2482
2483 END SUBROUTINE orthofun
2484
2485! **************************************************************************************************
2486!> \brief ...
2487!> \param l ...
2488!> \param zet ...
2489!> \param ca ...
2490!> \param cb ...
2491!> \param ss ...
2492! **************************************************************************************************
2493 SUBROUTINE aovlp(l, zet, ca, cb, ss)
2494 INTEGER, INTENT(IN) :: l
2495 REAL(kind=dp), DIMENSION(:), INTENT(IN) :: zet, ca, cb
2496 REAL(kind=dp), INTENT(OUT) :: ss
2497
2498 INTEGER :: i, j, m
2499 REAL(kind=dp) :: ab, ai, aj, s00, sss
2500
2501!
2502! use init_norm_cgf_orb
2503!
2504 m = SIZE(zet)
2505 ss = 0.0_dp
2506 DO i = 1, m
2507 ai = (2.0_dp*zet(i)/pi)**0.75_dp
2508 DO j = 1, m
2509 aj = (2.0_dp*zet(j)/pi)**0.75_dp
2510 ab = 1._dp/(zet(i) + zet(j))
2511 s00 = ai*aj*(pi*ab)**1.50_dp
2512 IF (l == 0) THEN
2513 sss = s00
2514 ELSEIF (l == 1) THEN
2515 sss = s00*ab*0.5_dp
2516 ELSE
2517 cpabort("aovlp lvalue")
2518 END IF
2519 ss = ss + sss*ca(i)*cb(j)
2520 END DO
2521 END DO
2522
2523 END SUBROUTINE aovlp
2524
2525! **************************************************************************************************
2526!> \brief ...
2527!> \param z ...
2528!> \param ne ...
2529!> \param n ...
2530!> \param l ...
2531!> \return ...
2532! **************************************************************************************************
2533 PURE FUNCTION srules(z, ne, n, l)
2534 ! Slater rules
2535 INTEGER, INTENT(IN) :: z
2536 INTEGER, DIMENSION(:, :), INTENT(IN) :: ne
2537 INTEGER, INTENT(IN) :: n, l
2538 REAL(dp) :: srules
2539
2540 REAL(dp), DIMENSION(7), PARAMETER :: &
2541 xns = (/1.0_dp, 2.0_dp, 3.0_dp, 3.7_dp, 4.0_dp, 4.2_dp, 4.4_dp/)
2542
2543 INTEGER :: i, l1, l2, m, m1, m2, nn
2544 REAL(dp) :: s
2545
2546 s = 0.0_dp
2547 ! The complete shell
2548 l1 = min(l + 1, 4)
2549 nn = min(n, 7)
2550 IF (l1 == 1) l2 = 2
2551 IF (l1 == 2) l2 = 1
2552 IF (l1 == 3) l2 = 4
2553 IF (l1 == 4) l2 = 3
2554 ! Rule a) no contribution from shells further out
2555 ! Rule b) 0.35 (1s 0.3) from each other electron in the same shell
2556 IF (n == 1) THEN
2557 m = ne(1, 1)
2558 s = s + 0.3_dp*real(m - 1, dp)
2559 ELSE
2560 m = ne(l1, nn) + ne(l2, nn)
2561 s = s + 0.35_dp*real(m - 1, dp)
2562 END IF
2563 ! Rule c) if (s,p) shell 0.85 from each electron with n-1, and 1.0
2564 ! from all electrons further in
2565 IF (l1 + l2 == 3) THEN
2566 IF (nn > 1) THEN
2567 m1 = ne(1, nn - 1) + ne(2, nn - 1) + ne(3, nn - 1) + ne(4, nn - 1)
2568 m2 = 0
2569 DO i = 1, nn - 2
2570 m2 = m2 + ne(1, i) + ne(2, i) + ne(3, i) + ne(4, i)
2571 END DO
2572 s = s + 0.85_dp*real(m1, dp) + 1._dp*real(m2, dp)
2573 END IF
2574 ELSE
2575 ! Rule d) if (d,f) shell 1.0 from each electron inside
2576 m = 0
2577 DO i = 1, nn - 1
2578 m = m + ne(1, i) + ne(2, i) + ne(3, i) + ne(4, i)
2579 END DO
2580 s = s + 1._dp*real(m, dp)
2581 END IF
2582 ! Slater exponent is (Z-S)/NS
2583 srules = (real(z, dp) - s)/xns(nn)
2584 END FUNCTION srules
2585
2586! **************************************************************************************************
2587!> \brief sort basis sets w.r.t. radius
2588!> \param basis_set ...
2589!> \param sort_method ...
2590! **************************************************************************************************
2591 SUBROUTINE sort_gto_basis_set(basis_set, sort_method)
2592 TYPE(gto_basis_set_type), INTENT(INOUT) :: basis_set
2593 INTEGER, INTENT(IN) :: sort_method
2594
2595 CHARACTER(LEN=12), DIMENSION(:), POINTER :: cgf_symbol
2596 CHARACTER(LEN=6), DIMENSION(:), POINTER :: sgf_symbol
2597 INTEGER :: ic, ic_max, icgf, icgf_new, icgf_old, ico, is, is_max, iset, isgf, isgf_new, &
2598 isgf_old, ishell, lshell, maxco, maxpgf, maxshell, mm, nc, ncgf, ns, nset
2599 INTEGER, ALLOCATABLE, DIMENSION(:) :: sort_index
2600 INTEGER, ALLOCATABLE, DIMENSION(:, :) :: icgf_set, isgf_set
2601 INTEGER, DIMENSION(:), POINTER :: lx, ly, lz, m, npgf
2602 REAL(dp), ALLOCATABLE, DIMENSION(:) :: tmp
2603 REAL(dp), DIMENSION(:), POINTER :: set_radius
2604 REAL(dp), DIMENSION(:, :), POINTER :: zet
2605 REAL(kind=dp), DIMENSION(:), POINTER :: norm_cgf
2606 REAL(kind=dp), DIMENSION(:, :), POINTER :: cphi, scon, sphi
2607
2608 NULLIFY (set_radius, zet)
2609
2610 IF (sort_method == basis_sort_default) RETURN
2611
2612 CALL get_gto_basis_set(gto_basis_set=basis_set, &
2613 nset=nset, &
2614 maxshell=maxshell, &
2615 maxpgf=maxpgf, &
2616 maxco=maxco, &
2617 ncgf=ncgf, &
2618 npgf=npgf, &
2619 set_radius=set_radius, &
2620 zet=zet)
2621
2622 ALLOCATE (sort_index(nset))
2623 ALLOCATE (tmp(nset))
2624 SELECT CASE (sort_method)
2625 CASE (basis_sort_zet)
2626 DO iset = 1, nset
2627 tmp(iset) = minval(basis_set%zet(:npgf(iset), iset))
2628 END DO
2629 CASE DEFAULT
2630 cpabort("Request basis sort criterion not implemented.")
2631 END SELECT
2632
2633 CALL sort(tmp(1:nset), nset, sort_index)
2634
2635 ic_max = 0
2636 is_max = 0
2637 DO iset = 1, nset
2638 ic = 0
2639 is = 0
2640 DO ishell = 1, basis_set%nshell(iset)
2641 DO ico = 1, nco(basis_set%l(ishell, iset))
2642 ic = ic + 1
2643 IF (ic > ic_max) ic_max = ic
2644 END DO
2645 lshell = basis_set%l(ishell, iset)
2646 DO mm = -lshell, lshell
2647 is = is + 1
2648 IF (is > is_max) is_max = is
2649 END DO
2650 END DO
2651 END DO
2652
2653 icgf = 0
2654 isgf = 0
2655 ALLOCATE (icgf_set(nset, ic_max))
2656 icgf_set(:, :) = 0
2657 ALLOCATE (isgf_set(nset, is_max))
2658 isgf_set(:, :) = 0
2659
2660 DO iset = 1, nset
2661 ic = 0
2662 is = 0
2663 DO ishell = 1, basis_set%nshell(iset)
2664 DO ico = 1, nco(basis_set%l(ishell, iset))
2665 icgf = icgf + 1
2666 ic = ic + 1
2667 icgf_set(iset, ic) = icgf
2668 END DO
2669 lshell = basis_set%l(ishell, iset)
2670 DO mm = -lshell, lshell
2671 isgf = isgf + 1
2672 is = is + 1
2673 isgf_set(iset, is) = isgf
2674 END DO
2675 END DO
2676 END DO
2677
2678 ALLOCATE (cgf_symbol(SIZE(basis_set%cgf_symbol)))
2679 ALLOCATE (norm_cgf(SIZE(basis_set%norm_cgf)))
2680 ALLOCATE (lx(SIZE(basis_set%lx)))
2681 ALLOCATE (ly(SIZE(basis_set%ly)))
2682 ALLOCATE (lz(SIZE(basis_set%lz)))
2683 ALLOCATE (cphi(SIZE(basis_set%cphi, 1), SIZE(basis_set%cphi, 2)))
2684 cphi = 0.0_dp
2685 ALLOCATE (sphi(SIZE(basis_set%sphi, 1), SIZE(basis_set%sphi, 2)))
2686 sphi = 0.0_dp
2687 ALLOCATE (scon(SIZE(basis_set%scon, 1), SIZE(basis_set%scon, 2)))
2688 scon = 0.0_dp
2689
2690 ALLOCATE (sgf_symbol(SIZE(basis_set%sgf_symbol)))
2691 ALLOCATE (m(SIZE(basis_set%m)))
2692
2693 icgf_new = 0
2694 isgf_new = 0
2695 DO iset = 1, nset
2696 DO ic = 1, ic_max
2697 icgf_old = icgf_set(sort_index(iset), ic)
2698 IF (icgf_old == 0) cycle
2699 icgf_new = icgf_new + 1
2700 norm_cgf(icgf_new) = basis_set%norm_cgf(icgf_old)
2701 lx(icgf_new) = basis_set%lx(icgf_old)
2702 ly(icgf_new) = basis_set%ly(icgf_old)
2703 lz(icgf_new) = basis_set%lz(icgf_old)
2704 cphi(:, icgf_new) = basis_set%cphi(:, icgf_old)
2705 cgf_symbol(icgf_new) = basis_set%cgf_symbol(icgf_old)
2706 END DO
2707 DO is = 1, is_max
2708 isgf_old = isgf_set(sort_index(iset), is)
2709 IF (isgf_old == 0) cycle
2710 isgf_new = isgf_new + 1
2711 m(isgf_new) = basis_set%m(isgf_old)
2712 sphi(:, isgf_new) = basis_set%sphi(:, isgf_old)
2713 scon(:, isgf_new) = basis_set%scon(:, isgf_old)
2714 sgf_symbol(isgf_new) = basis_set%sgf_symbol(isgf_old)
2715 END DO
2716 END DO
2717
2718 DEALLOCATE (basis_set%cgf_symbol)
2719 basis_set%cgf_symbol => cgf_symbol
2720 DEALLOCATE (basis_set%norm_cgf)
2721 basis_set%norm_cgf => norm_cgf
2722 DEALLOCATE (basis_set%lx)
2723 basis_set%lx => lx
2724 DEALLOCATE (basis_set%ly)
2725 basis_set%ly => ly
2726 DEALLOCATE (basis_set%lz)
2727 basis_set%lz => lz
2728 DEALLOCATE (basis_set%cphi)
2729 basis_set%cphi => cphi
2730 DEALLOCATE (basis_set%sphi)
2731 basis_set%sphi => sphi
2732 DEALLOCATE (basis_set%scon)
2733 basis_set%scon => scon
2734
2735 DEALLOCATE (basis_set%m)
2736 basis_set%m => m
2737 DEALLOCATE (basis_set%sgf_symbol)
2738 basis_set%sgf_symbol => sgf_symbol
2739
2740 basis_set%lmax = basis_set%lmax(sort_index)
2741 basis_set%lmin = basis_set%lmin(sort_index)
2742 basis_set%npgf = basis_set%npgf(sort_index)
2743 basis_set%nshell = basis_set%nshell(sort_index)
2744 basis_set%ncgf_set = basis_set%ncgf_set(sort_index)
2745 basis_set%nsgf_set = basis_set%nsgf_set(sort_index)
2746
2747 basis_set%n(:, :) = basis_set%n(:, sort_index)
2748 basis_set%l(:, :) = basis_set%l(:, sort_index)
2749 basis_set%zet(:, :) = basis_set%zet(:, sort_index)
2750
2751 basis_set%gcc(:, :, :) = basis_set%gcc(:, :, sort_index)
2752 basis_set%set_radius(:) = basis_set%set_radius(sort_index)
2753 basis_set%pgf_radius(:, :) = basis_set%pgf_radius(:, sort_index)
2754
2755 nc = 0
2756 ns = 0
2757 DO iset = 1, nset
2758 DO ishell = 1, basis_set%nshell(iset)
2759 lshell = basis_set%l(ishell, iset)
2760 basis_set%first_cgf(ishell, iset) = nc + 1
2761 nc = nc + nco(lshell)
2762 basis_set%last_cgf(ishell, iset) = nc
2763 basis_set%first_sgf(ishell, iset) = ns + 1
2764 ns = ns + nso(lshell)
2765 basis_set%last_sgf(ishell, iset) = ns
2766 END DO
2767 END DO
2768
2769 END SUBROUTINE
2770
2771END 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 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