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