(git:b279b6b)
soft_basis_set.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 !> none
11 !> \author MI (08.01.2004)
12 ! **************************************************************************************************
13 MODULE soft_basis_set
14 
15  USE ao_util, ONLY: exp_radius
16  USE basis_set_types, ONLY: copy_gto_basis_set,&
17  deallocate_gto_basis_set,&
18  get_gto_basis_set,&
19  gto_basis_set_type,&
20  init_cphi_and_sphi
21  USE kinds, ONLY: default_string_length,&
22  dp
23  USE memory_utilities, ONLY: reallocate
24  USE orbital_pointers, ONLY: indco,&
25  nco,&
26  ncoset,&
27  nso
28  USE orbital_symbols, ONLY: cgf_symbol,&
29  sgf_symbol
30 #include "../base/base_uses.f90"
31 
32  IMPLICIT NONE
33 
34  PRIVATE
35 
36 ! *** Global parameters (only in this module)
37 
38  CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'soft_basis_set'
39 
40  INTEGER, PARAMETER :: max_name_length = 60
41 
42 ! *** Public subroutines ***
43 
44  PUBLIC :: create_soft_basis
45 
46 CONTAINS
47 
48 ! **************************************************************************************************
49 !> \brief create the soft basis from a GTO basis
50 !> \param orb_basis ...
51 !> \param soft_basis ...
52 !> \param eps_fit ...
53 !> \param rc ...
54 !> \param paw_atom ...
55 !> \param paw_type_forced ...
56 !> \param gpw_r3d_rs_type_forced ...
57 !> \version 1.0
58 ! **************************************************************************************************
59  SUBROUTINE create_soft_basis(orb_basis, soft_basis, eps_fit, rc, paw_atom, &
60  paw_type_forced, gpw_r3d_rs_type_forced)
61 
62  TYPE(gto_basis_set_type), POINTER :: orb_basis, soft_basis
63  REAL(dp), INTENT(IN) :: eps_fit, rc
64  LOGICAL, INTENT(OUT) :: paw_atom
65  LOGICAL, INTENT(IN) :: paw_type_forced, gpw_r3d_rs_type_forced
66 
67  CHARACTER(LEN=default_string_length) :: bsname
68  INTEGER :: ico, ipgf, ipgf_s, iset, iset_s, ishell, lshell, lshell_old, m, maxco, maxpgf, &
69  maxpgf_s, maxshell, maxshell_s, ncgf, nset, nset_s, nsgf
70  INTEGER, ALLOCATABLE, DIMENSION(:) :: iset_s2h
71  INTEGER, DIMENSION(:), POINTER :: lmax, lmin, npgf, nshell
72  INTEGER, DIMENSION(:, :), POINTER :: l, n
73  LOGICAL :: my_gpw_r3d_rs_type_forced
74  REAL(kind=dp) :: minzet, radius
75  REAL(kind=dp), DIMENSION(:, :), POINTER :: zet
76  REAL(kind=dp), DIMENSION(:, :, :), POINTER :: gcc
77 
78  NULLIFY (gcc, l, lmax, lmin, n, npgf, nshell, zet)
79  paw_atom = .false.
80  my_gpw_r3d_rs_type_forced = gpw_r3d_rs_type_forced
81  IF (paw_type_forced) THEN
82  paw_atom = .true.
83  my_gpw_r3d_rs_type_forced = .false.
84  END IF
85  IF (.NOT. my_gpw_r3d_rs_type_forced) THEN
86  CALL get_gto_basis_set(gto_basis_set=orb_basis, name=bsname, &
87  maxpgf=maxpgf, maxshell=maxshell, nset=nset)
88 
89  soft_basis%name = trim(bsname)//"_soft"
90 
91  CALL reallocate(npgf, 1, nset)
92  CALL reallocate(nshell, 1, nset)
93  CALL reallocate(lmax, 1, nset)
94  CALL reallocate(lmin, 1, nset)
95 
96  CALL reallocate(n, 1, maxshell, 1, nset)
97  CALL reallocate(l, 1, maxshell, 1, nset)
98 
99  CALL reallocate(zet, 1, maxpgf, 1, nset)
100  CALL reallocate(gcc, 1, maxpgf, 1, maxshell, 1, nset)
101 
102  ALLOCATE (iset_s2h(nset))
103 
104  iset_s2h = 0
105  iset_s = 0
106  maxpgf_s = 0
107  maxshell_s = 0
108 
109  DO iset = 1, nset ! iset
110  minzet = orb_basis%zet(orb_basis%npgf(iset), iset)
111  DO ipgf = orb_basis%npgf(iset) - 1, 1, -1
112  IF (orb_basis%zet(ipgf, iset) < minzet) THEN
113  minzet = orb_basis%zet(ipgf, iset)
114  END IF
115  END DO
116  radius = exp_radius(orb_basis%lmax(iset), minzet, eps_fit, 1.0_dp)
117 
118  ! The soft basis contains this set
119  iset_s = iset_s + 1
120  nshell(iset_s) = orb_basis%nshell(iset)
121  lmax(iset_s) = orb_basis%lmax(iset)
122  lmin(iset_s) = orb_basis%lmin(iset)
123 
124  iset_s2h(iset_s) = iset
125 
126  DO ishell = 1, nshell(iset_s)
127  n(ishell, iset_s) = orb_basis%n(ishell, iset)
128  l(ishell, iset_s) = orb_basis%l(ishell, iset)
129  END DO
130 
131  IF (nshell(iset_s) > maxshell_s) THEN
132  maxshell_s = nshell(iset_s)
133  END IF
134 
135  IF (radius < rc) THEN
136  ! The soft basis does not contain this set
137  ! For the moment I keep the set as a dummy set
138  ! with no exponents, in order to have the right number of contractions
139  ! In a second time it can be taken away, by creating a pointer
140  ! which connects the remaining sets to the right contraction index
141  paw_atom = .true.
142  npgf(iset_s) = 0
143  cycle
144  END IF
145 
146  ipgf_s = 0
147  DO ipgf = 1, orb_basis%npgf(iset) ! ipgf
148  IF (orb_basis%zet(ipgf, iset) > 100.0_dp) THEN
149  ! The soft basis does not contain this exponent
150  paw_atom = .true.
151  cycle
152  END IF
153 
154  radius = exp_radius(orb_basis%lmax(iset), orb_basis%zet(ipgf, iset), &
155  eps_fit, 1.0_dp)
156  IF (radius < rc) THEN
157  ! The soft basis does not contain this exponent
158  paw_atom = .true.
159  cycle
160  END IF
161 
162  ! The soft basis contains this exponent
163  ipgf_s = ipgf_s + 1
164  zet(ipgf_s, iset_s) = orb_basis%zet(ipgf, iset)
165 
166  lshell_old = orb_basis%l(1, iset)
167  radius = exp_radius(lshell_old, zet(ipgf_s, iset_s), eps_fit, 1.0_dp)
168 
169  DO ishell = 1, nshell(iset_s)
170  lshell = orb_basis%l(ishell, iset)
171  IF (lshell == lshell_old) THEN
172  ELSE
173  lshell_old = lshell
174  radius = exp_radius(lshell_old, zet(ipgf_s, iset_s), eps_fit, 1.0_dp)
175  END IF
176  IF (radius < rc) THEN
177  gcc(ipgf_s, ishell, iset_s) = 0.0_dp
178  paw_atom = .true.
179  ELSE
180  gcc(ipgf_s, ishell, iset_s) = orb_basis%gcc(ipgf, ishell, iset)
181  END IF
182  END DO
183  END DO
184  npgf(iset_s) = ipgf_s
185  IF (ipgf_s > maxpgf_s) THEN
186  maxpgf_s = ipgf_s
187  END IF
188  END DO
189  nset_s = iset_s
190  IF (paw_atom) THEN
191  soft_basis%nset = nset_s
192  CALL reallocate(soft_basis%lmax, 1, nset_s)
193  CALL reallocate(soft_basis%lmin, 1, nset_s)
194  CALL reallocate(soft_basis%npgf, 1, nset_s)
195  CALL reallocate(soft_basis%nshell, 1, nset_s)
196  CALL reallocate(soft_basis%n, 1, maxshell_s, 1, nset_s)
197  CALL reallocate(soft_basis%l, 1, maxshell_s, 1, nset_s)
198  CALL reallocate(soft_basis%zet, 1, maxpgf_s, 1, nset_s)
199  CALL reallocate(soft_basis%gcc, 1, maxpgf_s, 1, maxshell_s, 1, nset_s)
200 
201  ! *** Copy the basis set information into the data structure ***
202 
203  DO iset = 1, nset_s
204  soft_basis%lmax(iset) = lmax(iset)
205  soft_basis%lmin(iset) = lmin(iset)
206  soft_basis%npgf(iset) = npgf(iset)
207  soft_basis%nshell(iset) = nshell(iset)
208  DO ishell = 1, nshell(iset)
209  soft_basis%n(ishell, iset) = n(ishell, iset)
210  soft_basis%l(ishell, iset) = l(ishell, iset)
211  DO ipgf = 1, npgf(iset)
212  soft_basis%gcc(ipgf, ishell, iset) = gcc(ipgf, ishell, iset)
213  END DO
214  END DO
215  DO ipgf = 1, npgf(iset)
216  soft_basis%zet(ipgf, iset) = zet(ipgf, iset)
217  END DO
218  END DO
219 
220  ! *** Initialise the depending soft_basis pointers ***
221  CALL reallocate(soft_basis%set_radius, 1, nset_s)
222  CALL reallocate(soft_basis%pgf_radius, 1, maxpgf_s, 1, nset_s)
223  CALL reallocate(soft_basis%first_cgf, 1, maxshell_s, 1, nset_s)
224  CALL reallocate(soft_basis%first_sgf, 1, maxshell_s, 1, nset_s)
225  CALL reallocate(soft_basis%last_cgf, 1, maxshell_s, 1, nset_s)
226  CALL reallocate(soft_basis%last_sgf, 1, maxshell_s, 1, nset_s)
227  CALL reallocate(soft_basis%ncgf_set, 1, nset_s)
228  CALL reallocate(soft_basis%nsgf_set, 1, nset_s)
229 
230  maxco = 0
231  ncgf = 0
232  nsgf = 0
233 
234  DO iset = 1, nset_s
235  soft_basis%ncgf_set(iset) = 0
236  soft_basis%nsgf_set(iset) = 0
237  DO ishell = 1, nshell(iset)
238  lshell = soft_basis%l(ishell, iset)
239  soft_basis%first_cgf(ishell, iset) = ncgf + 1
240  ncgf = ncgf + nco(lshell)
241  soft_basis%last_cgf(ishell, iset) = ncgf
242  soft_basis%ncgf_set(iset) = &
243  soft_basis%ncgf_set(iset) + nco(lshell)
244  soft_basis%first_sgf(ishell, iset) = nsgf + 1
245  nsgf = nsgf + nso(lshell)
246  soft_basis%last_sgf(ishell, iset) = nsgf
247  soft_basis%nsgf_set(iset) = &
248  soft_basis%nsgf_set(iset) + nso(lshell)
249  END DO
250  maxco = max(maxco, npgf(iset)*ncoset(lmax(iset)))
251  END DO
252  soft_basis%ncgf = ncgf
253  soft_basis%nsgf = nsgf
254 
255  CALL reallocate(soft_basis%cphi, 1, maxco, 1, ncgf)
256  CALL reallocate(soft_basis%sphi, 1, maxco, 1, nsgf)
257  CALL reallocate(soft_basis%scon, 1, maxco, 1, nsgf)
258  CALL reallocate(soft_basis%lx, 1, ncgf)
259  CALL reallocate(soft_basis%ly, 1, ncgf)
260  CALL reallocate(soft_basis%lz, 1, ncgf)
261  CALL reallocate(soft_basis%m, 1, nsgf)
262  CALL reallocate(soft_basis%norm_cgf, 1, ncgf)
263  ALLOCATE (soft_basis%cgf_symbol(ncgf))
264  ALLOCATE (soft_basis%sgf_symbol(nsgf))
265 
266  ncgf = 0
267  nsgf = 0
268 
269  DO iset = 1, nset_s
270  DO ishell = 1, nshell(iset)
271  lshell = soft_basis%l(ishell, iset)
272  DO ico = ncoset(lshell - 1) + 1, ncoset(lshell)
273  ncgf = ncgf + 1
274  soft_basis%lx(ncgf) = indco(1, ico)
275  soft_basis%ly(ncgf) = indco(2, ico)
276  soft_basis%lz(ncgf) = indco(3, ico)
277  soft_basis%cgf_symbol(ncgf) = &
278  cgf_symbol(n(ishell, iset), (/soft_basis%lx(ncgf), &
279  soft_basis%ly(ncgf), &
280  soft_basis%lz(ncgf)/))
281  END DO
282  DO m = -lshell, lshell
283  nsgf = nsgf + 1
284  soft_basis%m(nsgf) = m
285  soft_basis%sgf_symbol(nsgf) = &
286  sgf_symbol(n(ishell, iset), lshell, m)
287  END DO
288  END DO
289  END DO
290 
291  ! *** Normalization factor of the contracted Gaussians ***
292  soft_basis%norm_type = orb_basis%norm_type
293  soft_basis%norm_cgf = orb_basis%norm_cgf
294  ! *** Initialize the transformation matrices ***
295  CALL init_cphi_and_sphi(soft_basis)
296  END IF
297 
298  DEALLOCATE (gcc, l, lmax, lmin, n, npgf, nshell, zet, iset_s2h)
299  END IF
300 
301  IF (.NOT. paw_atom) THEN
302  CALL deallocate_gto_basis_set(soft_basis)
303  CALL copy_gto_basis_set(orb_basis, soft_basis)
304  CALL get_gto_basis_set(gto_basis_set=orb_basis, name=bsname)
305  soft_basis%name = trim(bsname)//"_soft"
306  END IF
307 
308  END SUBROUTINE create_soft_basis
309 
310 END MODULE soft_basis_set
ao_util
All kind of helpful little routines.
Definition: ao_util.F:14
ao_util::exp_radius
real(kind=dp) function, public exp_radius(l, alpha, threshold, prefactor, epsabs, epsrel, rlow)
The radius of a primitive Gaussian function for a given threshold is calculated. g(r) = prefactor*r**...
Definition: ao_util.F:96
basis_set_types
Definition: basis_set_types.F:15
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::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::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
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
memory_utilities
Utility routines for the memory handling.
Definition: memory_utilities.F:14
orbital_pointers
Provides Cartesian and spherical orbital pointers and indices.
Definition: orbital_pointers.F:15
orbital_pointers::nco
integer, dimension(:), allocatable, public nco
Definition: orbital_pointers.F:42
orbital_pointers::ncoset
integer, dimension(:), allocatable, public ncoset
Definition: orbital_pointers.F:42
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
soft_basis_set
Definition: soft_basis_set.F:13
soft_basis_set::create_soft_basis
subroutine, public create_soft_basis(orb_basis, soft_basis, eps_fit, rc, paw_atom, paw_type_forced, gpw_r3d_rs_type_forced)
create the soft basis from a GTO basis
Definition: soft_basis_set.F:61