Line data Source code
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 : !> \brief contains the structure
10 : !> \par History
11 : !> 11.2003 created [fawzi]
12 : !> \author fawzi
13 : ! **************************************************************************************************
14 : MODULE xc_rho_set_types
15 : USE cp_array_utils, ONLY: cp_3d_r_cp_type
16 : USE kinds, ONLY: dp
17 : USE pw_grid_types, ONLY: pw_grid_type
18 : USE pw_methods, ONLY: pw_copy, &
19 : pw_transfer
20 : USE pw_pool_types, ONLY: &
21 : pw_pool_type
22 : USE pw_spline_utils, ONLY: pw_spline_scale_deriv
23 : USE pw_types, ONLY: &
24 : pw_c1d_gs_type, &
25 : pw_r3d_rs_type
26 : USE xc_input_constants, ONLY: xc_deriv_pw, &
27 : xc_deriv_spline2, &
28 : xc_deriv_spline3, &
29 : xc_rho_no_smooth
30 : USE xc_rho_cflags_types, ONLY: xc_rho_cflags_equal, &
31 : xc_rho_cflags_setall, &
32 : xc_rho_cflags_type
33 : USE xc_util, ONLY: xc_pw_gradient, &
34 : xc_pw_laplace, &
35 : xc_pw_smooth
36 : #include "../base/base_uses.f90"
37 :
38 : IMPLICIT NONE
39 : PRIVATE
40 : LOGICAL, PRIVATE, PARAMETER :: debug_this_module = .FALSE.
41 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'xc_rho_set_types'
42 :
43 : PUBLIC :: xc_rho_set_type
44 : PUBLIC :: xc_rho_set_create, xc_rho_set_release, &
45 : xc_rho_set_update, xc_rho_set_get, xc_rho_set_recover_pw
46 :
47 : ! **************************************************************************************************
48 : !> \brief represent a density, with all the representation and data needed
49 : !> to perform a functional evaluation
50 : !> \param local_bounds the part of the 3d array on which the functional is
51 : !> computed
52 : !> \param owns which components are owned by this structure (and should
53 : !> be deallocated
54 : !> \param has which components are present and up to date
55 : !> \param rho the density
56 : !> \param drho the gradient of the density (x,y and z direction)
57 : !> \param norm_drho the norm of the gradient of the density
58 : !> \param rhoa , rhob: spin alpha and beta parts of the density in the LSD case
59 : !> \param drhoa , drhob: gradient of the spin alpha and beta parts of the density
60 : !> in the LSD case (x,y and z direction)
61 : !> \param norm_drhoa , norm_drhob: norm of the gradient of rhoa and rhob
62 : !> \param rho_ 1_3: rho^(1.0_dp/3.0_dp)
63 : !> \param rhoa_ 1_3, rhob_1_3: rhoa^(1.0_dp/3.0_dp), rhob^(1.0_dp/3.0_dp)
64 : !> \param tau the kinetic (KohnSham) part of rho
65 : !> \param tau_a the kinetic (KohnSham) part of rhoa
66 : !> \param tau_b the kinetic (KohnSham) part of rhob
67 : !> \note
68 : !> the use of 3d arrays is the result of trying to use only basic
69 : !> types (to be generic and independent from the method), and
70 : !> avoiding copies using the actual structure.
71 : !> only the part defined by local bounds is guaranteed to be present,
72 : !> and it is the only meaningful part.
73 : !> \par History
74 : !> 11.2003 created [fawzi & thomas]
75 : !> 12.2008 added laplace parts [mguidon]
76 : !> \author fawzi & thomas
77 : ! **************************************************************************************************
78 : TYPE xc_rho_set_type
79 : INTEGER, DIMENSION(2, 3) :: local_bounds = -1
80 : REAL(kind=dp) :: rho_cutoff = EPSILON(0.0_dp), drho_cutoff = EPSILON(0.0_dp), tau_cutoff = EPSILON(0.0_dp)
81 : TYPE(xc_rho_cflags_type) :: owns = xc_rho_cflags_type(), has = xc_rho_cflags_type()
82 : ! for spin restricted systems
83 : REAL(KIND=dp), DIMENSION(:, :, :), CONTIGUOUS, POINTER :: rho => NULL()
84 : TYPE(cp_3d_r_cp_type), DIMENSION(3) :: drho = cp_3d_r_cp_type()
85 : REAL(KIND=dp), DIMENSION(:, :, :), CONTIGUOUS, POINTER :: norm_drho => NULL()
86 : REAL(KIND=dp), DIMENSION(:, :, :), CONTIGUOUS, POINTER :: rho_1_3 => NULL()
87 : REAL(kind=dp), DIMENSION(:, :, :), CONTIGUOUS, POINTER :: tau => NULL()
88 : ! for UNrestricted systems
89 : REAL(KIND=dp), DIMENSION(:, :, :), CONTIGUOUS, POINTER :: rhoa => NULL(), rhob => NULL()
90 : TYPE(cp_3d_r_cp_type), DIMENSION(3) :: drhoa = cp_3d_r_cp_type(), &
91 : drhob = cp_3d_r_cp_type()
92 : REAL(KIND=dp), DIMENSION(:, :, :), CONTIGUOUS, POINTER :: norm_drhoa => NULL(), norm_drhob => NULL()
93 : REAL(kind=dp), DIMENSION(:, :, :), CONTIGUOUS, POINTER :: rhoa_1_3 => NULL(), rhob_1_3 => NULL()
94 : REAL(kind=dp), DIMENSION(:, :, :), CONTIGUOUS, POINTER :: tau_a => NULL(), tau_b => NULL()
95 : REAL(kind=dp), DIMENSION(:, :, :), CONTIGUOUS, POINTER :: laplace_rho => NULL(), laplace_rhoa => NULL(), &
96 : laplace_rhob => NULL()
97 : END TYPE xc_rho_set_type
98 :
99 : CONTAINS
100 :
101 : ! **************************************************************************************************
102 : !> \brief allocates and does (minimal) initialization of a rho_set
103 : !> \param rho_set the structure to allocate
104 : !> \param local_bounds ...
105 : !> \param rho_cutoff ...
106 : !> \param drho_cutoff ...
107 : !> \param tau_cutoff ...
108 : ! **************************************************************************************************
109 251487 : SUBROUTINE xc_rho_set_create(rho_set, local_bounds, rho_cutoff, drho_cutoff, &
110 : tau_cutoff)
111 : TYPE(xc_rho_set_type), INTENT(INOUT) :: rho_set
112 : INTEGER, DIMENSION(2, 3), INTENT(in) :: local_bounds
113 : REAL(kind=dp), INTENT(in), OPTIONAL :: rho_cutoff, drho_cutoff, tau_cutoff
114 :
115 251487 : IF (PRESENT(rho_cutoff)) rho_set%rho_cutoff = rho_cutoff
116 251487 : IF (PRESENT(drho_cutoff)) rho_set%drho_cutoff = drho_cutoff
117 251487 : IF (PRESENT(tau_cutoff)) rho_set%tau_cutoff = tau_cutoff
118 2514870 : rho_set%local_bounds = local_bounds
119 251487 : CALL xc_rho_cflags_setall(rho_set%owns, .TRUE.)
120 251487 : CALL xc_rho_cflags_setall(rho_set%has, .FALSE.)
121 251487 : END SUBROUTINE xc_rho_set_create
122 :
123 : ! **************************************************************************************************
124 : !> \brief releases the given rho_set
125 : !> \param rho_set the structure to release
126 : !> \param pw_pool the plae where to give back the arrays
127 : ! **************************************************************************************************
128 251487 : SUBROUTINE xc_rho_set_release(rho_set, pw_pool)
129 : TYPE(xc_rho_set_type), INTENT(INOUT) :: rho_set
130 : TYPE(pw_pool_type), OPTIONAL, POINTER :: pw_pool
131 :
132 : INTEGER :: i
133 :
134 251487 : IF (PRESENT(pw_pool)) THEN
135 120637 : IF (ASSOCIATED(pw_pool)) THEN
136 120637 : CALL xc_rho_set_clean(rho_set, pw_pool)
137 : ELSE
138 0 : CPABORT("pw_pool must be associated")
139 : END IF
140 : END IF
141 :
142 1005948 : rho_set%local_bounds(1, :) = -HUGE(0) ! we want to crash...
143 1005948 : rho_set%local_bounds(1, :) = HUGE(0)
144 251487 : IF (rho_set%owns%rho .AND. ASSOCIATED(rho_set%rho)) THEN
145 113715 : DEALLOCATE (rho_set%rho)
146 : ELSE
147 137772 : NULLIFY (rho_set%rho)
148 : END IF
149 251487 : IF (rho_set%owns%rho_spin) THEN
150 109286 : IF (ASSOCIATED(rho_set%rhoa)) THEN
151 17057 : DEALLOCATE (rho_set%rhoa)
152 : END IF
153 109286 : IF (ASSOCIATED(rho_set%rhob)) THEN
154 17057 : DEALLOCATE (rho_set%rhob)
155 : END IF
156 : ELSE
157 142201 : NULLIFY (rho_set%rhoa, rho_set%rhob)
158 : END IF
159 251487 : IF (rho_set%owns%rho_1_3 .AND. ASSOCIATED(rho_set%rho_1_3)) THEN
160 5233 : DEALLOCATE (rho_set%rho_1_3)
161 : ELSE
162 246254 : NULLIFY (rho_set%rho_1_3)
163 : END IF
164 251487 : IF (rho_set%owns%rho_spin) THEN
165 109286 : IF (ASSOCIATED(rho_set%rhoa_1_3)) THEN
166 2468 : DEALLOCATE (rho_set%rhoa_1_3)
167 : END IF
168 109286 : IF (ASSOCIATED(rho_set%rhob_1_3)) THEN
169 2468 : DEALLOCATE (rho_set%rhob_1_3)
170 : END IF
171 : ELSE
172 142201 : NULLIFY (rho_set%rhoa_1_3, rho_set%rhob_1_3)
173 : END IF
174 251487 : IF (rho_set%owns%drho) THEN
175 472880 : DO i = 1, 3
176 472880 : IF (ASSOCIATED(rho_set%drho(i)%array)) THEN
177 170292 : DEALLOCATE (rho_set%drho(i)%array)
178 : END IF
179 : END DO
180 : ELSE
181 533068 : DO i = 1, 3
182 533068 : NULLIFY (rho_set%drho(i)%array)
183 : END DO
184 : END IF
185 251487 : IF (rho_set%owns%drho_spin) THEN
186 430232 : DO i = 1, 3
187 322674 : IF (ASSOCIATED(rho_set%drhoa(i)%array)) THEN
188 27504 : DEALLOCATE (rho_set%drhoa(i)%array)
189 : END IF
190 430232 : IF (ASSOCIATED(rho_set%drhob(i)%array)) THEN
191 27504 : DEALLOCATE (rho_set%drhob(i)%array)
192 : END IF
193 : END DO
194 : ELSE
195 575716 : DO i = 1, 3
196 575716 : NULLIFY (rho_set%drhoa(i)%array, rho_set%drhob(i)%array)
197 : END DO
198 : END IF
199 251487 : IF (rho_set%owns%laplace_rho .AND. ASSOCIATED(rho_set%laplace_rho)) THEN
200 202 : DEALLOCATE (rho_set%laplace_rho)
201 : ELSE
202 251285 : NULLIFY (rho_set%laplace_rho)
203 : END IF
204 :
205 251487 : IF (rho_set%owns%norm_drho .AND. ASSOCIATED(rho_set%norm_drho)) THEN
206 88012 : DEALLOCATE (rho_set%norm_drho)
207 : ELSE
208 163475 : NULLIFY (rho_set%norm_drho)
209 : END IF
210 251487 : IF (rho_set%owns%laplace_rho_spin) THEN
211 106078 : IF (ASSOCIATED(rho_set%laplace_rhoa)) THEN
212 50 : DEALLOCATE (rho_set%laplace_rhoa)
213 : END IF
214 106078 : IF (ASSOCIATED(rho_set%laplace_rhob)) THEN
215 50 : DEALLOCATE (rho_set%laplace_rhob)
216 : END IF
217 : ELSE
218 145409 : NULLIFY (rho_set%laplace_rhoa, rho_set%laplace_rhob)
219 : END IF
220 :
221 251487 : IF (rho_set%owns%norm_drho_spin) THEN
222 107558 : IF (ASSOCIATED(rho_set%norm_drhoa)) THEN
223 10047 : DEALLOCATE (rho_set%norm_drhoa)
224 : END IF
225 107558 : IF (ASSOCIATED(rho_set%norm_drhob)) THEN
226 10047 : DEALLOCATE (rho_set%norm_drhob)
227 : END IF
228 : ELSE
229 143929 : NULLIFY (rho_set%norm_drhoa, rho_set%norm_drhob)
230 : END IF
231 251487 : IF (rho_set%owns%tau .AND. ASSOCIATED(rho_set%tau)) THEN
232 996 : DEALLOCATE (rho_set%tau)
233 : ELSE
234 250491 : NULLIFY (rho_set%tau)
235 : END IF
236 251487 : IF (rho_set%owns%tau_spin) THEN
237 106028 : IF (ASSOCIATED(rho_set%tau_a)) THEN
238 2 : DEALLOCATE (rho_set%tau_a)
239 : END IF
240 106028 : IF (ASSOCIATED(rho_set%tau_b)) THEN
241 2 : DEALLOCATE (rho_set%tau_b)
242 : END IF
243 : ELSE
244 145459 : NULLIFY (rho_set%tau_a, rho_set%tau_b)
245 : END IF
246 251487 : END SUBROUTINE xc_rho_set_release
247 :
248 : ! **************************************************************************************************
249 : !> \brief returns the various attributes of rho_set
250 : !> \param rho_set the object you want info about
251 : !> \param can_return_null if true the object returned can be null,
252 : !> if false (the default) it stops with an error if a requested
253 : !> component is not associated
254 : !> \param rho ...
255 : !> \param drho ...
256 : !> \param norm_drho ...
257 : !> \param rhoa ...
258 : !> \param rhob ...
259 : !> \param norm_drhoa ...
260 : !> \param norm_drhob ...
261 : !> \param rho_1_3 ...
262 : !> \param rhoa_1_3 ...
263 : !> \param rhob_1_3 ...
264 : !> \param laplace_rho ...
265 : !> \param laplace_rhoa ...
266 : !> \param laplace_rhob ...
267 : !> \param drhoa ...
268 : !> \param drhob ...
269 : !> \param rho_cutoff ...
270 : !> \param drho_cutoff ...
271 : !> \param tau_cutoff ...
272 : !> \param tau ...
273 : !> \param tau_a ...
274 : !> \param tau_b ...
275 : !> \param local_bounds ...
276 : ! **************************************************************************************************
277 1010048 : SUBROUTINE xc_rho_set_get(rho_set, can_return_null, rho, drho, norm_drho, &
278 : rhoa, rhob, norm_drhoa, norm_drhob, rho_1_3, rhoa_1_3, &
279 : rhob_1_3, laplace_rho, laplace_rhoa, laplace_rhob, drhoa, drhob, rho_cutoff, &
280 : drho_cutoff, tau_cutoff, tau, tau_a, tau_b, local_bounds)
281 : TYPE(xc_rho_set_type), INTENT(IN) :: rho_set
282 : LOGICAL, INTENT(in), OPTIONAL :: can_return_null
283 : REAL(KIND=dp), DIMENSION(:, :, :), OPTIONAL, &
284 : POINTER :: rho
285 : TYPE(cp_3d_r_cp_type), DIMENSION(3), OPTIONAL :: drho
286 : REAL(KIND=dp), DIMENSION(:, :, :), OPTIONAL, &
287 : POINTER :: norm_drho, rhoa, rhob, norm_drhoa, &
288 : norm_drhob, rho_1_3, rhoa_1_3, &
289 : rhob_1_3, laplace_rho, laplace_rhoa, &
290 : laplace_rhob
291 : TYPE(cp_3d_r_cp_type), DIMENSION(3), OPTIONAL :: drhoa, drhob
292 : REAL(kind=dp), INTENT(out), OPTIONAL :: rho_cutoff, drho_cutoff, tau_cutoff
293 : REAL(KIND=dp), DIMENSION(:, :, :), OPTIONAL, &
294 : POINTER :: tau, tau_a, tau_b
295 : INTEGER, DIMENSION(2, 3), INTENT(OUT), OPTIONAL :: local_bounds
296 :
297 : INTEGER :: i
298 : LOGICAL :: my_can_return_null
299 :
300 1010048 : my_can_return_null = .FALSE.
301 1010048 : IF (PRESENT(can_return_null)) my_can_return_null = can_return_null
302 :
303 1010048 : IF (PRESENT(rho)) THEN
304 264139 : rho => rho_set%rho
305 264139 : CPASSERT(my_can_return_null .OR. ASSOCIATED(rho))
306 : END IF
307 1010048 : IF (PRESENT(drho)) THEN
308 691044 : DO i = 1, 3
309 518283 : drho(i)%array => rho_set%drho(i)%array
310 691044 : CPASSERT(my_can_return_null .OR. ASSOCIATED(rho_set%drho(i)%array))
311 : END DO
312 : END IF
313 1010048 : IF (PRESENT(norm_drho)) THEN
314 365400 : norm_drho => rho_set%norm_drho
315 365400 : CPASSERT(my_can_return_null .OR. ASSOCIATED(norm_drho))
316 : END IF
317 1010048 : IF (PRESENT(laplace_rho)) THEN
318 12232 : laplace_rho => rho_set%laplace_rho
319 12232 : CPASSERT(my_can_return_null .OR. ASSOCIATED(laplace_rho))
320 : END IF
321 1010048 : IF (PRESENT(rhoa)) THEN
322 144711 : rhoa => rho_set%rhoa
323 144711 : CPASSERT(my_can_return_null .OR. ASSOCIATED(rhoa))
324 : END IF
325 1010048 : IF (PRESENT(rhob)) THEN
326 144711 : rhob => rho_set%rhob
327 144711 : CPASSERT(my_can_return_null .OR. ASSOCIATED(rhob))
328 : END IF
329 1010048 : IF (PRESENT(drhoa)) THEN
330 575044 : DO i = 1, 3
331 431283 : drhoa(i)%array => rho_set%drhoa(i)%array
332 575044 : CPASSERT(my_can_return_null .OR. ASSOCIATED(rho_set%drhoa(i)%array))
333 : END DO
334 : END IF
335 1010048 : IF (PRESENT(drhob)) THEN
336 575044 : DO i = 1, 3
337 431283 : drhob(i)%array => rho_set%drhob(i)%array
338 575044 : CPASSERT(my_can_return_null .OR. ASSOCIATED(rho_set%drhob(i)%array))
339 : END DO
340 : END IF
341 1010048 : IF (PRESENT(laplace_rhoa)) THEN
342 3434 : laplace_rhoa => rho_set%laplace_rhoa
343 3434 : CPASSERT(my_can_return_null .OR. ASSOCIATED(laplace_rhoa))
344 : END IF
345 1010048 : IF (PRESENT(laplace_rhob)) THEN
346 3434 : laplace_rhob => rho_set%laplace_rhob
347 3434 : CPASSERT(my_can_return_null .OR. ASSOCIATED(laplace_rhob))
348 : END IF
349 1010048 : IF (PRESENT(norm_drhoa)) THEN
350 193691 : norm_drhoa => rho_set%norm_drhoa
351 193691 : CPASSERT(my_can_return_null .OR. ASSOCIATED(norm_drhoa))
352 : END IF
353 1010048 : IF (PRESENT(norm_drhob)) THEN
354 193691 : norm_drhob => rho_set%norm_drhob
355 193691 : CPASSERT(my_can_return_null .OR. ASSOCIATED(norm_drhob))
356 : END IF
357 1010048 : IF (PRESENT(rho_1_3)) THEN
358 21023 : rho_1_3 => rho_set%rho_1_3
359 21023 : CPASSERT(my_can_return_null .OR. ASSOCIATED(rho_1_3))
360 : END IF
361 1010048 : IF (PRESENT(rhoa_1_3)) THEN
362 3274 : rhoa_1_3 => rho_set%rhoa_1_3
363 3274 : CPASSERT(my_can_return_null .OR. ASSOCIATED(rhoa_1_3))
364 : END IF
365 1010048 : IF (PRESENT(rhob_1_3)) THEN
366 3274 : rhob_1_3 => rho_set%rhob_1_3
367 3274 : CPASSERT(my_can_return_null .OR. ASSOCIATED(rhob_1_3))
368 : END IF
369 1010048 : IF (PRESENT(tau)) THEN
370 14170 : tau => rho_set%tau
371 14170 : CPASSERT(my_can_return_null .OR. ASSOCIATED(tau))
372 : END IF
373 1010048 : IF (PRESENT(tau_a)) THEN
374 3786 : tau_a => rho_set%tau_a
375 3786 : CPASSERT(my_can_return_null .OR. ASSOCIATED(tau_a))
376 : END IF
377 1010048 : IF (PRESENT(tau_b)) THEN
378 3786 : tau_b => rho_set%tau_b
379 3786 : CPASSERT(my_can_return_null .OR. ASSOCIATED(tau_b))
380 : END IF
381 1010048 : IF (PRESENT(rho_cutoff)) rho_cutoff = rho_set%rho_cutoff
382 1010048 : IF (PRESENT(drho_cutoff)) drho_cutoff = rho_set%drho_cutoff
383 1010048 : IF (PRESENT(tau_cutoff)) tau_cutoff = rho_set%tau_cutoff
384 2585728 : IF (PRESENT(local_bounds)) local_bounds = rho_set%local_bounds
385 1010048 : END SUBROUTINE xc_rho_set_get
386 :
387 : #:mute
388 : #:def recover(variable)
389 : #! Determine component of the actual data
390 : #:set var_long_pw = (variable+"(i)" if variable.startswith("drho") else variable)
391 : #:set var_long_rho = (variable+"(i)%array" if variable.startswith("drho") else variable)
392 : #! Determine the flag name
393 : #! Remove spin states and potential underscore
394 : #:set is_1_3 = variable.endswith("_1_3")
395 : #:set var_base = variable.strip("_13")
396 : #:set is_spin = var_base.endswith("a") or var_base.endswith("b")
397 : #:set var_base = var_base.strip("ab_")
398 : #:set var_cflags = (var_base if not is_spin else var_base+"_spin")
399 : #:set var_cflags = (var_cflags if not is_1_3 else var_cflags+"_1_3")
400 : IF (PRESENT(${variable}$)) THEN
401 : #:if variable.startswith("drho")
402 : DO i = 1, 3
403 : #:else
404 : NULLIFY (${var_long_pw}$)
405 : ALLOCATE (${var_long_pw}$)
406 : #:endif
407 : CALL xc_rho_set_recover_pw_low(${var_long_pw}$, rho_set%${var_long_rho}$, pw_grid, pw_pool#{if variable =="drho"}#, rho_set%drhoa(i)%array, rho_set%drhob(i)%array#{endif}#)
408 : #:if not variable.startswith("drho")
409 : NULLIFY (rho_set%${var_long_rho}$)
410 : #:else
411 : END DO
412 : #:endif
413 : owns_data = #{if variable =="drho"}#.TRUE.#{else}#rho_set%owns%${var_cflags}$#{endif}#
414 : END IF
415 : #:enddef
416 : #:endmute
417 :
418 : ! **************************************************************************************************
419 : !> \brief Shifts association of the requested array to a pw grid
420 : !> Requires that the corresponding component of rho_set is associated
421 : !> If owns_data returns TRUE, the caller has to allocate the data later
422 : !> It is allowed to task for only one component per call.
423 : !> In case of drho, the array is allocated if not internally available and calculated from drhoa and drhob.
424 : !> \param rho_set the object you want info about
425 : !> \param pw_grid ...
426 : !> \param pw_pool ...
427 : !> \param owns_data ...
428 : !> \param rho ...
429 : !> \param drho ...
430 : !> \param norm_drho ...
431 : !> \param rhoa ...
432 : !> \param rhob ...
433 : !> \param norm_drhoa ...
434 : !> \param norm_drhob ...
435 : !> \param rho_1_3 ...
436 : !> \param rhoa_1_3 ...
437 : !> \param rhob_1_3 ...
438 : !> \param laplace_rho ...
439 : !> \param laplace_rhoa ...
440 : !> \param laplace_rhob ...
441 : !> \param drhoa ...
442 : !> \param drhob ...
443 : !> \param tau ...
444 : !> \param tau_a ...
445 : !> \param tau_b ...
446 : ! **************************************************************************************************
447 504980 : SUBROUTINE xc_rho_set_recover_pw(rho_set, pw_grid, pw_pool, owns_data, rho, drho, norm_drho, &
448 : rhoa, rhob, norm_drhoa, norm_drhob, rho_1_3, rhoa_1_3, &
449 : rhob_1_3, laplace_rho, laplace_rhoa, laplace_rhob, drhoa, drhob, &
450 : tau, tau_a, tau_b)
451 : TYPE(xc_rho_set_type) :: rho_set
452 : TYPE(pw_r3d_rs_type), DIMENSION(3), OPTIONAL, INTENT(OUT) :: drho, drhoa, drhob
453 : TYPE(pw_r3d_rs_type), OPTIONAL, POINTER :: rho, norm_drho, rhoa, rhob, norm_drhoa, &
454 : norm_drhob, rho_1_3, rhoa_1_3, &
455 : rhob_1_3, laplace_rho, laplace_rhoa, &
456 : laplace_rhob, tau, tau_a, tau_b
457 : TYPE(pw_grid_type), POINTER, INTENT(IN) :: pw_grid
458 : TYPE(pw_pool_type), POINTER, INTENT(IN) :: pw_pool
459 : LOGICAL, INTENT(OUT) :: owns_data
460 :
461 : INTEGER :: i
462 :
463 : #:for variable in ["rho", "drho", "norm_drho", "rhoa", "rhob", "norm_drhoa", "norm_drhob", "rho_1_3", "rhoa_1_3", "rhob_1_3", "laplace_rho", "laplace_rhoa", "laplace_rhob", "drhoa", "drhob", "tau", "tau_a", "tau_b"]
464 432375 : $:recover(variable)
465 : #:endfor
466 :
467 86475 : END SUBROUTINE
468 :
469 259425 : SUBROUTINE xc_rho_set_recover_pw_low(rho_pw, rho, pw_grid, pw_pool, rhoa, rhob)
470 : TYPE(pw_r3d_rs_type), INTENT(OUT) :: rho_pw
471 : REAL(KIND=dp), DIMENSION(:, :, :), POINTER, CONTIGUOUS :: rho
472 : TYPE(pw_grid_type), POINTER, INTENT(IN) :: pw_grid
473 : TYPE(pw_pool_type), POINTER, INTENT(IN) :: pw_pool
474 : REAL(KIND=dp), DIMENSION(:, :, :), POINTER, OPTIONAL :: rhoa, rhob
475 :
476 259425 : IF (ASSOCIATED(rho)) THEN
477 219897 : CALL rho_pw%create(pw_grid=pw_grid, array_ptr=rho)
478 219897 : NULLIFY (rho)
479 39528 : ELSE IF (PRESENT(rhoa) .AND. PRESENT(rhob)) THEN
480 39528 : CALL pw_pool%create_pw(rho_pw)
481 39528 : !$OMP PARALLEL WORKSHARE DEFAULT(NONE) SHARED(rho_pw,rhoa,rhob)
482 : rho_pw%array(:, :, :) = rhoa(:, :, :) + rhob(:, :, :)
483 : !$OMP END PARALLEL WORKSHARE
484 : ELSE
485 : CALL cp_abort(__LOCATION__, "Either component or its spin parts (if applicable) "// &
486 0 : "have to be associated in rho_set!")
487 : END IF
488 :
489 259425 : END SUBROUTINE
490 :
491 : ! **************************************************************************************************
492 : !> \brief cleans (releases) most of the data stored in the rho_set giving back
493 : !> what it can to the pw_pool
494 : !> \param rho_set the rho_set to be cleaned
495 : !> \param pw_pool place to give back 3d arrays,...
496 : !> \author Fawzi Mohamed
497 : ! **************************************************************************************************
498 272546 : SUBROUTINE xc_rho_set_clean(rho_set, pw_pool)
499 : TYPE(xc_rho_set_type), INTENT(INOUT) :: rho_set
500 : TYPE(pw_pool_type), POINTER :: pw_pool
501 :
502 : INTEGER :: idir
503 :
504 272546 : IF (rho_set%owns%rho) THEN
505 244916 : CALL pw_pool%give_back_cr3d(rho_set%rho)
506 : ELSE
507 27630 : NULLIFY (rho_set%rho)
508 : END IF
509 272546 : IF (rho_set%owns%rho_1_3) THEN
510 154103 : CALL pw_pool%give_back_cr3d(rho_set%rho_1_3)
511 : ELSE
512 118443 : NULLIFY (rho_set%rho_1_3)
513 : END IF
514 272546 : IF (rho_set%owns%drho) THEN
515 783560 : DO idir = 1, 3
516 783560 : CALL pw_pool%give_back_cr3d(rho_set%drho(idir)%array)
517 : END DO
518 : ELSE
519 306624 : DO idir = 1, 3
520 306624 : NULLIFY (rho_set%drho(idir)%array)
521 : END DO
522 : END IF
523 272546 : IF (rho_set%owns%norm_drho) THEN
524 214392 : CALL pw_pool%give_back_cr3d(rho_set%norm_drho)
525 : ELSE
526 58154 : NULLIFY (rho_set%norm_drho)
527 : END IF
528 272546 : IF (rho_set%owns%laplace_rho) THEN
529 146137 : CALL pw_pool%give_back_cr3d(rho_set%laplace_rho)
530 : ELSE
531 126409 : NULLIFY (rho_set%laplace_rho)
532 : END IF
533 272546 : IF (rho_set%owns%tau) THEN
534 145459 : CALL pw_pool%give_back_cr3d(rho_set%tau)
535 : ELSE
536 127087 : NULLIFY (rho_set%tau)
537 : END IF
538 272546 : IF (rho_set%owns%rho_spin) THEN
539 173089 : CALL pw_pool%give_back_cr3d(rho_set%rhoa)
540 173089 : CALL pw_pool%give_back_cr3d(rho_set%rhob)
541 : ELSE
542 99457 : NULLIFY (rho_set%rhoa, rho_set%rhob)
543 : END IF
544 272546 : IF (rho_set%owns%rho_spin_1_3) THEN
545 147113 : CALL pw_pool%give_back_cr3d(rho_set%rhoa_1_3)
546 147113 : CALL pw_pool%give_back_cr3d(rho_set%rhob_1_3)
547 : ELSE
548 125433 : NULLIFY (rho_set%rhoa_1_3, rho_set%rhob_1_3)
549 : END IF
550 272546 : IF (rho_set%owns%drho_spin) THEN
551 640188 : DO idir = 1, 3
552 480141 : CALL pw_pool%give_back_cr3d(rho_set%drhoa(idir)%array)
553 640188 : CALL pw_pool%give_back_cr3d(rho_set%drhob(idir)%array)
554 : END DO
555 : ELSE
556 449996 : DO idir = 1, 3
557 449996 : NULLIFY (rho_set%drhoa(idir)%array, rho_set%drhob(idir)%array)
558 : END DO
559 : END IF
560 272546 : IF (rho_set%owns%laplace_rho_spin) THEN
561 145703 : CALL pw_pool%give_back_cr3d(rho_set%laplace_rhoa)
562 145703 : CALL pw_pool%give_back_cr3d(rho_set%laplace_rhob)
563 : ELSE
564 126843 : NULLIFY (rho_set%laplace_rhoa, rho_set%laplace_rhob)
565 : END IF
566 272546 : IF (rho_set%owns%norm_drho_spin) THEN
567 162387 : CALL pw_pool%give_back_cr3d(rho_set%norm_drhoa)
568 162387 : CALL pw_pool%give_back_cr3d(rho_set%norm_drhob)
569 : ELSE
570 110159 : NULLIFY (rho_set%norm_drhoa, rho_set%norm_drhob)
571 : END IF
572 272546 : IF (rho_set%owns%tau_spin) THEN
573 145459 : CALL pw_pool%give_back_cr3d(rho_set%tau_a)
574 145459 : CALL pw_pool%give_back_cr3d(rho_set%tau_b)
575 : ELSE
576 127087 : NULLIFY (rho_set%tau_a, rho_set%tau_b)
577 : END IF
578 :
579 272546 : CALL xc_rho_cflags_setall(rho_set%has, .FALSE.)
580 272546 : CALL xc_rho_cflags_setall(rho_set%owns, .FALSE.)
581 :
582 272546 : END SUBROUTINE xc_rho_set_clean
583 :
584 : ! **************************************************************************************************
585 : !> \brief updates the given rho set with the density given by
586 : !> rho_r (and rho_g). The rho set will contain the components specified
587 : !> in needs
588 : !> \param rho_set the rho_set to update
589 : !> \param rho_r the new density (in r space)
590 : !> \param rho_g the new density (in g space, needed for some
591 : !> derivatives)
592 : !> \param tau ...
593 : !> \param needs the components of rho that are needed
594 : !> \param xc_deriv_method_id ...
595 : !> \param xc_rho_smooth_id ...
596 : !> \param pw_pool pool for the allocation of pw and cr3d
597 : ! **************************************************************************************************
598 151909 : SUBROUTINE xc_rho_set_update(rho_set, rho_r, rho_g, tau, needs, &
599 : xc_deriv_method_id, xc_rho_smooth_id, pw_pool)
600 : TYPE(xc_rho_set_type), INTENT(INOUT) :: rho_set
601 : TYPE(pw_r3d_rs_type), DIMENSION(:), INTENT(IN) :: rho_r
602 : TYPE(pw_c1d_gs_type), DIMENSION(:), POINTER :: rho_g
603 : TYPE(pw_r3d_rs_type), DIMENSION(:), POINTER, INTENT(IN) :: tau
604 : TYPE(xc_rho_cflags_type), INTENT(in) :: needs
605 : INTEGER, INTENT(IN) :: xc_deriv_method_id, xc_rho_smooth_id
606 : TYPE(pw_pool_type), POINTER :: pw_pool
607 :
608 : REAL(KIND=dp), PARAMETER :: f13 = (1.0_dp/3.0_dp)
609 :
610 : INTEGER :: i, idir, ispin, j, k, nspins
611 : LOGICAL :: gradient_f, my_rho_g_local, &
612 : needs_laplace, needs_rho_g
613 : REAL(kind=dp) :: rho_cutoff
614 455727 : TYPE(pw_r3d_rs_type), DIMENSION(2) :: laplace_rho_r
615 1367181 : TYPE(pw_r3d_rs_type), DIMENSION(3, 2) :: drho_r
616 : TYPE(pw_c1d_gs_type) :: my_rho_g, tmp_g
617 455727 : TYPE(pw_r3d_rs_type), DIMENSION(2) :: my_rho_r
618 :
619 1519090 : IF (ANY(rho_set%local_bounds /= pw_pool%pw_grid%bounds_local)) &
620 0 : CPABORT("pw_pool cr3d have different size than expected")
621 151909 : nspins = SIZE(rho_r)
622 1519090 : rho_set%local_bounds = rho_r(1)%pw_grid%bounds_local
623 151909 : rho_cutoff = 0.5*rho_set%rho_cutoff
624 :
625 151909 : my_rho_g_local = .FALSE.
626 : ! some checks
627 121021 : SELECT CASE (nspins)
628 : CASE (1)
629 121021 : CPASSERT(.NOT. needs%rho_spin)
630 121021 : CPASSERT(.NOT. needs%drho_spin)
631 121021 : CPASSERT(.NOT. needs%norm_drho_spin)
632 121021 : CPASSERT(.NOT. needs%rho_spin_1_3)
633 121021 : CPASSERT(.NOT. needs%tau_spin)
634 121021 : CPASSERT(.NOT. needs%laplace_rho_spin)
635 : CASE (2)
636 30888 : CPASSERT(.NOT. needs%rho)
637 30888 : CPASSERT(.NOT. needs%drho)
638 30888 : CPASSERT(.NOT. needs%rho_1_3)
639 30888 : CPASSERT(.NOT. needs%tau)
640 30888 : CPASSERT(.NOT. needs%laplace_rho)
641 : CASE default
642 151909 : CPABORT("Unknown number of spin states")
643 : END SELECT
644 :
645 151909 : CALL xc_rho_set_clean(rho_set, pw_pool=pw_pool)
646 :
647 151909 : needs_laplace = (needs%laplace_rho .OR. needs%laplace_rho_spin)
648 : gradient_f = (needs%drho_spin .OR. needs%norm_drho_spin .OR. &
649 : needs%drho .OR. needs%norm_drho .OR. &
650 151909 : needs_laplace)
651 : needs_rho_g = ((xc_deriv_method_id == xc_deriv_spline3 .OR. &
652 : xc_deriv_method_id == xc_deriv_spline2 .OR. &
653 151909 : xc_deriv_method_id == xc_deriv_pw)) .AND. (gradient_f .OR. needs_laplace)
654 151909 : IF ((gradient_f .AND. needs_laplace) .AND. &
655 : (xc_deriv_method_id /= xc_deriv_pw)) THEN
656 : CALL cp_abort(__LOCATION__, &
657 : "MGGA functionals that require the Laplacian are "// &
658 0 : "only compatible with 'XC_DERIV PW' and 'XC_SMOOTH_RHO NONE'")
659 : END IF
660 :
661 151909 : IF (needs_rho_g) THEN
662 81659 : CALL pw_pool%create_pw(tmp_g)
663 : END IF
664 334706 : DO ispin = 1, nspins
665 182797 : CALL pw_pool%create_pw(my_rho_r(ispin))
666 : ! introduce a smoothing kernel on the density
667 182797 : IF (xc_rho_smooth_id == xc_rho_no_smooth) THEN
668 182299 : IF (needs_rho_g) THEN
669 99399 : IF (ASSOCIATED(rho_g)) THEN
670 84215 : my_rho_g_local = .FALSE.
671 84215 : my_rho_g = rho_g(ispin)
672 : END IF
673 : END IF
674 :
675 182299 : CALL pw_copy(rho_r(ispin), my_rho_r(ispin))
676 : ELSE
677 498 : CALL xc_pw_smooth(rho_r(ispin), my_rho_r(ispin), xc_rho_smooth_id)
678 : END IF
679 :
680 334706 : IF (gradient_f) THEN ! calculate the grad of rho
681 : ! normally when you need the gradient you need the whole gradient
682 : ! (for the partial integration)
683 : ! deriv rho
684 408620 : DO idir = 1, 3
685 408620 : CALL pw_pool%create_pw(drho_r(idir, ispin))
686 : END DO
687 102155 : IF (needs_rho_g) THEN
688 99465 : IF (.NOT. ASSOCIATED(my_rho_g%pw_grid)) THEN
689 15250 : my_rho_g_local = .TRUE.
690 15250 : CALL pw_pool%create_pw(my_rho_g)
691 15250 : CALL pw_transfer(my_rho_r(ispin), my_rho_g)
692 : END IF
693 84215 : IF (.NOT. my_rho_g_local .AND. (xc_deriv_method_id == xc_deriv_spline2 .OR. &
694 : xc_deriv_method_id == xc_deriv_spline3)) THEN
695 7122 : CALL pw_pool%create_pw(my_rho_g)
696 7122 : my_rho_g_local = .TRUE.
697 7122 : CALL pw_copy(rho_g(ispin), my_rho_g)
698 : END IF
699 : END IF
700 102155 : IF (needs%laplace_rho .OR. needs%laplace_rho_spin) THEN
701 1468 : CALL pw_pool%create_pw(laplace_rho_r(ispin))
702 1468 : CALL xc_pw_laplace(my_rho_g, pw_pool, xc_deriv_method_id, laplace_rho_r(ispin), tmp_g=tmp_g)
703 : END IF
704 102155 : CALL xc_pw_gradient(my_rho_r(ispin), my_rho_g, tmp_g, drho_r(:, ispin), xc_deriv_method_id)
705 :
706 102155 : IF (needs_rho_g) THEN
707 99465 : IF (my_rho_g_local) THEN
708 22372 : my_rho_g_local = .FALSE.
709 22372 : CALL pw_pool%give_back_pw(my_rho_g)
710 : END IF
711 : END IF
712 :
713 102155 : IF (xc_deriv_method_id /= xc_deriv_pw) THEN
714 9866 : CALL pw_spline_scale_deriv(drho_r(:, ispin))
715 : END IF
716 :
717 : END IF
718 :
719 : END DO
720 :
721 151909 : IF (ASSOCIATED(tmp_g%pw_grid)) THEN
722 81659 : CALL pw_pool%give_back_pw(tmp_g)
723 : END IF
724 :
725 121021 : SELECT CASE (nspins)
726 : CASE (1)
727 121021 : IF (needs%rho_1_3) THEN
728 9870 : CALL pw_pool%create_cr3d(rho_set%rho_1_3)
729 9870 : !$OMP PARALLEL DO DEFAULT(NONE) PRIVATE(i,j,k) SHARED(rho_set,my_rho_r)
730 : DO k = rho_set%local_bounds(1, 3), rho_set%local_bounds(2, 3)
731 : DO j = rho_set%local_bounds(1, 2), rho_set%local_bounds(2, 2)
732 : DO i = rho_set%local_bounds(1, 1), rho_set%local_bounds(2, 1)
733 : rho_set%rho_1_3(i, j, k) = MAX(my_rho_r(1)%array(i, j, k), 0.0_dp)**f13
734 : END DO
735 : END DO
736 : END DO
737 9870 : rho_set%owns%rho_1_3 = .TRUE.
738 9870 : rho_set%has%rho_1_3 = .TRUE.
739 : END IF
740 121021 : IF (needs%rho) THEN
741 121021 : rho_set%rho => my_rho_r(1)%array
742 121021 : NULLIFY (my_rho_r(1)%array)
743 121021 : rho_set%owns%rho = .TRUE.
744 121021 : rho_set%has%rho = .TRUE.
745 : END IF
746 121021 : IF (needs%norm_drho) THEN
747 65239 : CALL pw_pool%create_cr3d(rho_set%norm_drho)
748 65239 : !$OMP PARALLEL DO DEFAULT(NONE) PRIVATE(i,j,k) SHARED(rho_set,drho_r)
749 : DO k = rho_set%local_bounds(1, 3), rho_set%local_bounds(2, 3)
750 : DO j = rho_set%local_bounds(1, 2), rho_set%local_bounds(2, 2)
751 : DO i = rho_set%local_bounds(1, 1), rho_set%local_bounds(2, 1)
752 : rho_set%norm_drho(i, j, k) = SQRT( &
753 : drho_r(1, 1)%array(i, j, k)**2 + &
754 : drho_r(2, 1)%array(i, j, k)**2 + &
755 : drho_r(3, 1)%array(i, j, k)**2)
756 : END DO
757 : END DO
758 : END DO
759 65239 : rho_set%owns%norm_drho = .TRUE.
760 65239 : rho_set%has%norm_drho = .TRUE.
761 : END IF
762 121021 : IF (needs%laplace_rho) THEN
763 880 : rho_set%laplace_rho => laplace_rho_r(1)%array
764 880 : NULLIFY (laplace_rho_r(1)%array)
765 880 : rho_set%owns%laplace_rho = .TRUE.
766 880 : rho_set%has%laplace_rho = .TRUE.
767 : END IF
768 :
769 121021 : IF (needs%drho) THEN
770 250492 : DO idir = 1, 3
771 187869 : rho_set%drho(idir)%array => drho_r(idir, 1)%array
772 250492 : NULLIFY (drho_r(idir, 1)%array)
773 : END DO
774 62623 : rho_set%owns%drho = .TRUE.
775 62623 : rho_set%has%drho = .TRUE.
776 : END IF
777 : CASE (2)
778 30888 : IF (needs%rho_spin_1_3) THEN
779 1666 : CALL pw_pool%create_cr3d(rho_set%rhoa_1_3)
780 : !assume that the bounds are the same?
781 1666 : !$OMP PARALLEL DO DEFAULT(NONE) PRIVATE(i,j,k) SHARED(rho_set,my_rho_r)
782 : DO k = rho_set%local_bounds(1, 3), rho_set%local_bounds(2, 3)
783 : DO j = rho_set%local_bounds(1, 2), rho_set%local_bounds(2, 2)
784 : DO i = rho_set%local_bounds(1, 1), rho_set%local_bounds(2, 1)
785 : rho_set%rhoa_1_3(i, j, k) = MAX(my_rho_r(1)%array(i, j, k), 0.0_dp)**f13
786 : END DO
787 : END DO
788 : END DO
789 1666 : CALL pw_pool%create_cr3d(rho_set%rhob_1_3)
790 : !assume that the bounds are the same?
791 1666 : !$OMP PARALLEL DO DEFAULT(NONE) PRIVATE(i,j,k) SHARED(rho_set,my_rho_r)
792 : DO k = rho_set%local_bounds(1, 3), rho_set%local_bounds(2, 3)
793 : DO j = rho_set%local_bounds(1, 2), rho_set%local_bounds(2, 2)
794 : DO i = rho_set%local_bounds(1, 1), rho_set%local_bounds(2, 1)
795 : rho_set%rhob_1_3(i, j, k) = MAX(my_rho_r(2)%array(i, j, k), 0.0_dp)**f13
796 : END DO
797 : END DO
798 : END DO
799 1666 : rho_set%owns%rho_spin_1_3 = .TRUE.
800 1666 : rho_set%has%rho_spin_1_3 = .TRUE.
801 : END IF
802 30888 : IF (needs%norm_drho) THEN
803 :
804 17396 : CALL pw_pool%create_cr3d(rho_set%norm_drho)
805 17396 : !$OMP PARALLEL DO DEFAULT(NONE) PRIVATE(i,j,k) SHARED(rho_set,drho_r)
806 : DO k = rho_set%local_bounds(1, 3), rho_set%local_bounds(2, 3)
807 : DO j = rho_set%local_bounds(1, 2), rho_set%local_bounds(2, 2)
808 : DO i = rho_set%local_bounds(1, 1), rho_set%local_bounds(2, 1)
809 : rho_set%norm_drho(i, j, k) = SQRT( &
810 : (drho_r(1, 1)%array(i, j, k) + drho_r(1, 2)%array(i, j, k))**2 + &
811 : (drho_r(2, 1)%array(i, j, k) + drho_r(2, 2)%array(i, j, k))**2 + &
812 : (drho_r(3, 1)%array(i, j, k) + drho_r(3, 2)%array(i, j, k))**2)
813 : END DO
814 : END DO
815 : END DO
816 :
817 17396 : rho_set%owns%norm_drho = .TRUE.
818 17396 : rho_set%has%norm_drho = .TRUE.
819 : END IF
820 30888 : IF (needs%rho_spin) THEN
821 :
822 30888 : rho_set%rhoa => my_rho_r(1)%array
823 30888 : NULLIFY (my_rho_r(1)%array)
824 :
825 30888 : rho_set%rhob => my_rho_r(2)%array
826 30888 : NULLIFY (my_rho_r(2)%array)
827 :
828 30888 : rho_set%owns%rho_spin = .TRUE.
829 30888 : rho_set%has%rho_spin = .TRUE.
830 : END IF
831 30888 : IF (needs%norm_drho_spin) THEN
832 :
833 18458 : CALL pw_pool%create_cr3d(rho_set%norm_drhoa)
834 18458 : !$OMP PARALLEL DO DEFAULT(NONE) PRIVATE(i,j,k) SHARED(rho_set,drho_r)
835 : DO k = rho_set%local_bounds(1, 3), rho_set%local_bounds(2, 3)
836 : DO j = rho_set%local_bounds(1, 2), rho_set%local_bounds(2, 2)
837 : DO i = rho_set%local_bounds(1, 1), rho_set%local_bounds(2, 1)
838 : rho_set%norm_drhoa(i, j, k) = SQRT( &
839 : drho_r(1, 1)%array(i, j, k)**2 + &
840 : drho_r(2, 1)%array(i, j, k)**2 + &
841 : drho_r(3, 1)%array(i, j, k)**2)
842 : END DO
843 : END DO
844 : END DO
845 :
846 18458 : CALL pw_pool%create_cr3d(rho_set%norm_drhob)
847 18458 : rho_set%owns%norm_drho_spin = .TRUE.
848 18458 : !$OMP PARALLEL DO DEFAULT(NONE) PRIVATE(i,j,k) SHARED(rho_set,drho_r)
849 : DO k = rho_set%local_bounds(1, 3), rho_set%local_bounds(2, 3)
850 : DO j = rho_set%local_bounds(1, 2), rho_set%local_bounds(2, 2)
851 : DO i = rho_set%local_bounds(1, 1), rho_set%local_bounds(2, 1)
852 : rho_set%norm_drhob(i, j, k) = SQRT( &
853 : drho_r(1, 2)%array(i, j, k)**2 + &
854 : drho_r(2, 2)%array(i, j, k)**2 + &
855 : drho_r(3, 2)%array(i, j, k)**2)
856 : END DO
857 : END DO
858 : END DO
859 :
860 18458 : rho_set%owns%norm_drho_spin = .TRUE.
861 18458 : rho_set%has%norm_drho_spin = .TRUE.
862 : END IF
863 30888 : IF (needs%laplace_rho_spin) THEN
864 294 : rho_set%laplace_rhoa => laplace_rho_r(1)%array
865 294 : NULLIFY (laplace_rho_r(1)%array)
866 :
867 294 : rho_set%laplace_rhob => laplace_rho_r(2)%array
868 294 : NULLIFY (laplace_rho_r(2)%array)
869 :
870 294 : rho_set%owns%laplace_rho_spin = .TRUE.
871 294 : rho_set%has%laplace_rho_spin = .TRUE.
872 : END IF
873 182797 : IF (needs%drho_spin) THEN
874 64472 : DO idir = 1, 3
875 48354 : rho_set%drhoa(idir)%array => drho_r(idir, 1)%array
876 48354 : NULLIFY (drho_r(idir, 1)%array)
877 48354 : rho_set%drhob(idir)%array => drho_r(idir, 2)%array
878 64472 : NULLIFY (drho_r(idir, 2)%array)
879 : END DO
880 16118 : rho_set%owns%drho_spin = .TRUE.
881 16118 : rho_set%has%drho_spin = .TRUE.
882 : END IF
883 : END SELECT
884 : ! post cleanup
885 334706 : DO ispin = 1, nspins
886 182797 : IF (needs%laplace_rho .OR. needs%laplace_rho_spin) THEN
887 1468 : CALL pw_pool%give_back_pw(laplace_rho_r(ispin))
888 : END IF
889 883097 : DO idir = 1, 3
890 731188 : CALL pw_pool%give_back_pw(drho_r(idir, ispin))
891 : END DO
892 : END DO
893 334706 : DO ispin = 1, nspins
894 334706 : CALL pw_pool%give_back_pw(my_rho_r(ispin))
895 : END DO
896 :
897 : ! tau part
898 151909 : IF (needs%tau .OR. needs%tau_spin) THEN
899 3588 : CPASSERT(ASSOCIATED(tau))
900 7932 : DO ispin = 1, nspins
901 156253 : CPASSERT(ASSOCIATED(tau(ispin)%array))
902 : END DO
903 : END IF
904 151909 : IF (needs%tau) THEN
905 2832 : rho_set%tau => tau(1)%array
906 2832 : rho_set%owns%tau = .FALSE.
907 2832 : rho_set%has%tau = .TRUE.
908 : END IF
909 151909 : IF (needs%tau_spin) THEN
910 756 : rho_set%tau_a => tau(1)%array
911 756 : rho_set%tau_b => tau(2)%array
912 756 : rho_set%owns%tau_spin = .FALSE.
913 756 : rho_set%has%tau_spin = .TRUE.
914 : END IF
915 :
916 151909 : CPASSERT(xc_rho_cflags_equal(rho_set%has, needs))
917 :
918 303818 : END SUBROUTINE xc_rho_set_update
919 :
920 0 : END MODULE xc_rho_set_types
|
