Line data Source code
1 : !--------------------------------------------------------------------------------------------------!
2 : ! CP2K: A general program to perform molecular dynamics simulations !
3 : ! Copyright 2000-2025 CP2K developers group <https://cp2k.org> !
4 : ! !
5 : ! SPDX-License-Identifier: GPL-2.0-or-later !
6 : !--------------------------------------------------------------------------------------------------!
7 : ! **************************************************************************************************
8 : MODULE qs_local_rho_types
9 :
10 : USE kinds, ONLY: dp
11 : USE mathconstants, ONLY: fourpi,&
12 : pi
13 : USE memory_utilities, ONLY: reallocate
14 : USE qs_cneo_types, ONLY: deallocate_rhoz_cneo_set,&
15 : rhoz_cneo_type
16 : USE qs_grid_atom, ONLY: grid_atom_type
17 : USE qs_harmonics_atom, ONLY: harmonics_atom_type
18 : USE qs_rho0_types, ONLY: deallocate_rho0_atom,&
19 : deallocate_rho0_mpole,&
20 : rho0_atom_type,&
21 : rho0_mpole_type
22 : USE qs_rho_atom_types, ONLY: deallocate_rho_atom_set,&
23 : rho_atom_type
24 : #include "./base/base_uses.f90"
25 :
26 : IMPLICIT NONE
27 :
28 : PRIVATE
29 :
30 : ! *** Global parameters (only in this module)
31 :
32 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_local_rho_types'
33 :
34 : ! *** Define rhoz and local_rho types ***
35 :
36 : ! **************************************************************************************************
37 : TYPE rhoz_type
38 : REAL(dp) :: one_atom = -1.0_dp
39 : REAL(dp), DIMENSION(:), POINTER :: r_coef => NULL()
40 : REAL(dp), DIMENSION(:), POINTER :: dr_coef => NULL()
41 : REAL(dp), DIMENSION(:), POINTER :: vr_coef => NULL()
42 : END TYPE rhoz_type
43 :
44 : ! **************************************************************************************************
45 : TYPE local_rho_type
46 : TYPE(rho_atom_type), DIMENSION(:), POINTER :: rho_atom_set => NULL()
47 : TYPE(rho0_mpole_type), POINTER :: rho0_mpole => NULL()
48 : TYPE(rho0_atom_type), DIMENSION(:), POINTER :: rho0_atom_set => NULL()
49 : TYPE(rhoz_type), DIMENSION(:), POINTER :: rhoz_set => NULL()
50 : TYPE(rhoz_cneo_type), DIMENSION(:), POINTER :: rhoz_cneo_set => NULL()
51 : REAL(dp) :: rhoz_tot = -1.0_dp, &
52 : rhoz_cneo_tot = -1.0_dp
53 : END TYPE local_rho_type
54 :
55 : ! Public Types
56 : PUBLIC :: local_rho_type, rhoz_type
57 :
58 : ! Public Subroutine
59 : PUBLIC :: allocate_rhoz, calculate_rhoz, &
60 : get_local_rho, local_rho_set_create, &
61 : local_rho_set_release, set_local_rho
62 :
63 : CONTAINS
64 :
65 : ! **************************************************************************************************
66 : !> \brief ...
67 : !> \param rhoz_set ...
68 : !> \param nkind ...
69 : ! **************************************************************************************************
70 2176 : SUBROUTINE allocate_rhoz(rhoz_set, nkind)
71 :
72 : TYPE(rhoz_type), DIMENSION(:), POINTER :: rhoz_set
73 : INTEGER :: nkind
74 :
75 : INTEGER :: ikind
76 :
77 2176 : IF (ASSOCIATED(rhoz_set)) THEN
78 0 : CALL deallocate_rhoz(rhoz_set)
79 : END IF
80 :
81 10816 : ALLOCATE (rhoz_set(nkind))
82 :
83 6464 : DO ikind = 1, nkind
84 4288 : NULLIFY (rhoz_set(ikind)%r_coef)
85 4288 : NULLIFY (rhoz_set(ikind)%dr_coef)
86 6464 : NULLIFY (rhoz_set(ikind)%vr_coef)
87 : END DO
88 :
89 2176 : END SUBROUTINE allocate_rhoz
90 :
91 : ! **************************************************************************************************
92 : !> \brief ...
93 : !> \param rhoz ...
94 : !> \param grid_atom ...
95 : !> \param alpha ...
96 : !> \param zeff ...
97 : !> \param natom ...
98 : !> \param rhoz_tot ...
99 : !> \param harmonics ...
100 : ! **************************************************************************************************
101 4280 : SUBROUTINE calculate_rhoz(rhoz, grid_atom, alpha, zeff, natom, rhoz_tot, harmonics)
102 :
103 : TYPE(rhoz_type) :: rhoz
104 : TYPE(grid_atom_type) :: grid_atom
105 : REAL(dp), INTENT(IN) :: alpha
106 : REAL(dp) :: zeff
107 : INTEGER :: natom
108 : REAL(dp), INTENT(INOUT) :: rhoz_tot
109 : TYPE(harmonics_atom_type) :: harmonics
110 :
111 : INTEGER :: ir, na, nr
112 : REAL(dp) :: c1, c2, c3, prefactor1, prefactor2, &
113 : prefactor3, sum
114 :
115 4280 : nr = grid_atom%nr
116 4280 : na = grid_atom%ng_sphere
117 4280 : CALL reallocate(rhoz%r_coef, 1, nr)
118 4280 : CALL reallocate(rhoz%dr_coef, 1, nr)
119 4280 : CALL reallocate(rhoz%vr_coef, 1, nr)
120 :
121 4280 : c1 = alpha/pi
122 4280 : c2 = c1*c1*c1*fourpi
123 4280 : c3 = SQRT(alpha)
124 4280 : prefactor1 = zeff*SQRT(c2)
125 4280 : prefactor2 = -2.0_dp*alpha
126 4280 : prefactor3 = -zeff*SQRT(fourpi)
127 :
128 4280 : sum = 0.0_dp
129 223200 : DO ir = 1, nr
130 218920 : c1 = -alpha*grid_atom%rad2(ir)
131 218920 : rhoz%r_coef(ir) = -EXP(c1)*prefactor1
132 218920 : IF (ABS(rhoz%r_coef(ir)) < 1.0E-30_dp) THEN
133 135758 : rhoz%r_coef(ir) = 0.0_dp
134 135758 : rhoz%dr_coef(ir) = 0.0_dp
135 : ELSE
136 83162 : rhoz%dr_coef(ir) = prefactor2*rhoz%r_coef(ir)
137 : END IF
138 218920 : rhoz%vr_coef(ir) = prefactor3*erf(grid_atom%rad(ir)*c3)/grid_atom%rad(ir)
139 223200 : sum = sum + rhoz%r_coef(ir)*grid_atom%wr(ir)
140 : END DO
141 4280 : rhoz%one_atom = sum*harmonics%slm_int(1)
142 4280 : rhoz_tot = rhoz_tot + natom*rhoz%one_atom
143 :
144 4280 : END SUBROUTINE calculate_rhoz
145 :
146 : ! **************************************************************************************************
147 : !> \brief ...
148 : !> \param rhoz_set ...
149 : ! **************************************************************************************************
150 2176 : SUBROUTINE deallocate_rhoz(rhoz_set)
151 :
152 : TYPE(rhoz_type), DIMENSION(:), POINTER :: rhoz_set
153 :
154 : INTEGER :: ikind, nkind
155 :
156 2176 : nkind = SIZE(rhoz_set)
157 :
158 6464 : DO ikind = 1, nkind
159 4288 : IF (ASSOCIATED(rhoz_set(ikind)%r_coef)) &
160 4280 : DEALLOCATE (rhoz_set(ikind)%r_coef)
161 4288 : IF (ASSOCIATED(rhoz_set(ikind)%dr_coef)) &
162 4280 : DEALLOCATE (rhoz_set(ikind)%dr_coef)
163 4288 : IF (ASSOCIATED(rhoz_set(ikind)%vr_coef)) &
164 6456 : DEALLOCATE (rhoz_set(ikind)%vr_coef)
165 : END DO
166 :
167 2176 : DEALLOCATE (rhoz_set)
168 :
169 2176 : END SUBROUTINE deallocate_rhoz
170 :
171 : ! **************************************************************************************************
172 : !> \brief ...
173 : !> \param local_rho_set ...
174 : !> \param rho_atom_set ...
175 : !> \param rho0_atom_set ...
176 : !> \param rho0_mpole ...
177 : !> \param rhoz_set ...
178 : !> \param rhoz_cneo_set ...
179 : ! **************************************************************************************************
180 314021 : SUBROUTINE get_local_rho(local_rho_set, rho_atom_set, rho0_atom_set, rho0_mpole, rhoz_set, &
181 : rhoz_cneo_set)
182 :
183 : TYPE(local_rho_type), POINTER :: local_rho_set
184 : TYPE(rho_atom_type), DIMENSION(:), OPTIONAL, &
185 : POINTER :: rho_atom_set
186 : TYPE(rho0_atom_type), DIMENSION(:), OPTIONAL, &
187 : POINTER :: rho0_atom_set
188 : TYPE(rho0_mpole_type), OPTIONAL, POINTER :: rho0_mpole
189 : TYPE(rhoz_type), DIMENSION(:), OPTIONAL, POINTER :: rhoz_set
190 : TYPE(rhoz_cneo_type), DIMENSION(:), OPTIONAL, &
191 : POINTER :: rhoz_cneo_set
192 :
193 314021 : IF (PRESENT(rho_atom_set)) rho_atom_set => local_rho_set%rho_atom_set
194 314021 : IF (PRESENT(rho0_atom_set)) rho0_atom_set => local_rho_set%rho0_atom_set
195 314021 : IF (PRESENT(rho0_mpole)) rho0_mpole => local_rho_set%rho0_mpole
196 314021 : IF (PRESENT(rhoz_set)) rhoz_set => local_rho_set%rhoz_set
197 314021 : IF (PRESENT(rhoz_cneo_set)) rhoz_cneo_set => local_rho_set%rhoz_cneo_set
198 :
199 314021 : END SUBROUTINE get_local_rho
200 :
201 : ! **************************************************************************************************
202 : !> \brief ...
203 : !> \param local_rho_set ...
204 : ! **************************************************************************************************
205 9680 : SUBROUTINE local_rho_set_create(local_rho_set)
206 :
207 : TYPE(local_rho_type), POINTER :: local_rho_set
208 :
209 9680 : ALLOCATE (local_rho_set)
210 :
211 : NULLIFY (local_rho_set%rho_atom_set)
212 : NULLIFY (local_rho_set%rho0_atom_set)
213 : NULLIFY (local_rho_set%rho0_mpole)
214 : NULLIFY (local_rho_set%rhoz_set)
215 : NULLIFY (local_rho_set%rhoz_cneo_set)
216 :
217 9680 : local_rho_set%rhoz_tot = 0.0_dp
218 9680 : local_rho_set%rhoz_cneo_tot = 0.0_dp
219 :
220 9680 : END SUBROUTINE local_rho_set_create
221 :
222 : ! **************************************************************************************************
223 : !> \brief ...
224 : !> \param local_rho_set ...
225 : ! **************************************************************************************************
226 9680 : SUBROUTINE local_rho_set_release(local_rho_set)
227 :
228 : TYPE(local_rho_type), POINTER :: local_rho_set
229 :
230 9680 : IF (ASSOCIATED(local_rho_set)) THEN
231 9680 : IF (ASSOCIATED(local_rho_set%rho_atom_set)) THEN
232 3232 : CALL deallocate_rho_atom_set(local_rho_set%rho_atom_set)
233 : END IF
234 :
235 9680 : IF (ASSOCIATED(local_rho_set%rho0_atom_set)) THEN
236 2176 : CALL deallocate_rho0_atom(local_rho_set%rho0_atom_set)
237 : END IF
238 :
239 9680 : IF (ASSOCIATED(local_rho_set%rho0_mpole)) THEN
240 2176 : CALL deallocate_rho0_mpole(local_rho_set%rho0_mpole)
241 : END IF
242 :
243 9680 : IF (ASSOCIATED(local_rho_set%rhoz_set)) THEN
244 2176 : CALL deallocate_rhoz(local_rho_set%rhoz_set)
245 : END IF
246 :
247 9680 : IF (ASSOCIATED(local_rho_set%rhoz_cneo_set)) THEN
248 8 : CALL deallocate_rhoz_cneo_set(local_rho_set%rhoz_cneo_set)
249 : END IF
250 :
251 9680 : DEALLOCATE (local_rho_set)
252 : END IF
253 :
254 9680 : END SUBROUTINE local_rho_set_release
255 :
256 : ! **************************************************************************************************
257 : !> \brief ...
258 : !> \param local_rho_set ...
259 : !> \param rho_atom_set ...
260 : !> \param rho0_atom_set ...
261 : !> \param rho0_mpole ...
262 : !> \param rhoz_set ...
263 : !> \param rhoz_cneo_set ...
264 : ! **************************************************************************************************
265 3174 : SUBROUTINE set_local_rho(local_rho_set, rho_atom_set, rho0_atom_set, rho0_mpole, &
266 : rhoz_set, rhoz_cneo_set)
267 :
268 : TYPE(local_rho_type), POINTER :: local_rho_set
269 : TYPE(rho_atom_type), DIMENSION(:), OPTIONAL, &
270 : POINTER :: rho_atom_set
271 : TYPE(rho0_atom_type), DIMENSION(:), OPTIONAL, &
272 : POINTER :: rho0_atom_set
273 : TYPE(rho0_mpole_type), OPTIONAL, POINTER :: rho0_mpole
274 : TYPE(rhoz_type), DIMENSION(:), OPTIONAL, POINTER :: rhoz_set
275 : TYPE(rhoz_cneo_type), DIMENSION(:), OPTIONAL, &
276 : POINTER :: rhoz_cneo_set
277 :
278 3174 : IF (PRESENT(rho_atom_set)) THEN
279 998 : IF (ASSOCIATED(local_rho_set%rho_atom_set)) THEN
280 0 : CALL deallocate_rho_atom_set(local_rho_set%rho_atom_set)
281 : END IF
282 998 : local_rho_set%rho_atom_set => rho_atom_set
283 : END IF
284 :
285 3174 : IF (PRESENT(rho0_atom_set)) THEN
286 2176 : IF (ASSOCIATED(local_rho_set%rho0_atom_set)) THEN
287 0 : CALL deallocate_rho0_atom(local_rho_set%rho0_atom_set)
288 : END IF
289 2176 : local_rho_set%rho0_atom_set => rho0_atom_set
290 : END IF
291 :
292 3174 : IF (PRESENT(rho0_mpole)) THEN
293 2176 : IF (ASSOCIATED(local_rho_set%rho0_mpole)) THEN
294 0 : CALL deallocate_rho0_mpole(local_rho_set%rho0_mpole)
295 : END IF
296 2176 : local_rho_set%rho0_mpole => rho0_mpole
297 : END IF
298 :
299 3174 : IF (PRESENT(rhoz_set)) THEN
300 2176 : IF (ASSOCIATED(local_rho_set%rhoz_set)) THEN
301 0 : CALL deallocate_rhoz(local_rho_set%rhoz_set)
302 : END IF
303 2176 : local_rho_set%rhoz_set => rhoz_set
304 : END IF
305 :
306 3174 : IF (PRESENT(rhoz_cneo_set)) THEN
307 2176 : IF (ASSOCIATED(local_rho_set%rhoz_cneo_set)) THEN
308 0 : CALL deallocate_rhoz_cneo_set(local_rho_set%rhoz_cneo_set)
309 : END IF
310 2176 : local_rho_set%rhoz_cneo_set => rhoz_cneo_set
311 : END IF
312 :
313 3174 : END SUBROUTINE set_local_rho
314 :
315 0 : END MODULE qs_local_rho_types
|
