Line data Source code
1 : !--------------------------------------------------------------------------------------------------!
2 : ! CP2K: A general program to perform molecular dynamics simulations !
3 : ! Copyright 2000-2026 CP2K developers group <https://cp2k.org> !
4 : ! !
5 : ! SPDX-License-Identifier: GPL-2.0-or-later !
6 : !--------------------------------------------------------------------------------------------------!
7 :
8 : #ifdef __GAUXC
9 : #include "gauxc/gauxc_config.f"
10 : #endif
11 :
12 : MODULE xc_gauxc_cache
13 :
14 : USE iso_c_binding, ONLY: c_double
15 : USE kinds, ONLY: default_path_length,&
16 : default_string_length
17 : USE message_passing, ONLY: mp_comm_self,&
18 : mp_para_env_type
19 : USE particle_types, ONLY: particle_type
20 : USE qs_kind_types, ONLY: qs_kind_type
21 : USE xc_gauxc_interface, ONLY: &
22 : cp_gauxc_basisset_type, cp_gauxc_grid_type, cp_gauxc_integrator_type, &
23 : cp_gauxc_molecule_type, cp_gauxc_status_type, gauxc_check_status, gauxc_create_basisset, &
24 : gauxc_create_grid, gauxc_create_integrator, gauxc_create_molecule, gauxc_destroy_basisset, &
25 : gauxc_destroy_grid, gauxc_destroy_integrator, gauxc_destroy_molecule
26 : #include "../base/base_uses.f90"
27 :
28 : IMPLICIT NONE
29 : PRIVATE
30 :
31 : TYPE cp_gauxc_cache_params
32 : INTEGER :: natom = -1
33 : INTEGER :: nspins = -1
34 : INTEGER :: batch_size = -1
35 : REAL(c_double) :: device_runtime_fill_fraction = 0.0_c_double
36 : CHARACTER(LEN=default_string_length) :: xc_fun_name = ""
37 : CHARACTER(LEN=default_string_length) :: grid_type = ""
38 : CHARACTER(LEN=default_string_length) :: radial_quadrature = ""
39 : CHARACTER(LEN=default_string_length) :: pruning_scheme = ""
40 : CHARACTER(LEN=default_string_length) :: lb_exec_space = ""
41 : CHARACTER(LEN=default_string_length) :: int_exec_space = ""
42 : CHARACTER(LEN=default_path_length) :: model_eval_name = ""
43 : CHARACTER(len=default_string_length) :: lwd_kernel = ""
44 : LOGICAL :: use_mpi_runtime = .FALSE.
45 : LOGICAL :: use_gradient_mpi_runtime = .FALSE.
46 : LOGICAL :: use_self_runtime = .FALSE.
47 : LOGICAL :: use_gradient_self_runtime = .FALSE.
48 : LOGICAL :: use_fd_gradient = .FALSE.
49 : LOGICAL :: use_gauxc_model = .FALSE.
50 : END TYPE cp_gauxc_cache_params
51 :
52 : TYPE, extends(cp_gauxc_cache_params) :: cp_gauxc_cache_type
53 : TYPE(cp_gauxc_molecule_type) :: molecule = cp_gauxc_molecule_type()
54 : TYPE(cp_gauxc_basisset_type) :: basisset = cp_gauxc_basisset_type()
55 : TYPE(cp_gauxc_grid_type) :: grid = cp_gauxc_grid_type()
56 : TYPE(cp_gauxc_grid_type) :: gradient_grid = cp_gauxc_grid_type()
57 : TYPE(cp_gauxc_integrator_type) :: integrator = cp_gauxc_integrator_type()
58 : TYPE(cp_gauxc_integrator_type) :: gradient_integrator = cp_gauxc_integrator_type()
59 : INTEGER :: mpi_comm = -1
60 : LOGICAL :: is_init = .FALSE.
61 : REAL(c_double), DIMENSION(:, :), ALLOCATABLE :: last_positions
62 : CONTAINS
63 :
64 : PROCEDURE, PUBLIC :: needs_gradient_grid => gauxc_cache_type_needs_gradient_grid
65 : PROCEDURE, PUBLIC :: max_l => gauxc_cache_type_max_l
66 : END TYPE cp_gauxc_cache_type
67 :
68 : PUBLIC :: &
69 : cp_gauxc_cache_type, &
70 : gauxc_cache_init, &
71 : gauxc_cache_is_valid, &
72 : gauxc_cache_release, cp_gauxc_cache_params
73 :
74 : CONTAINS
75 :
76 : ! **************************************************************************************************
77 : !> \brief ...
78 : !> \param this ...
79 : !> \return ...
80 : ! **************************************************************************************************
81 44 : FUNCTION gauxc_cache_type_needs_gradient_grid(this) RESULT(res)
82 : class(cp_gauxc_cache_type) :: this
83 : LOGICAL :: res
84 :
85 44 : res = this%use_gradient_self_runtime
86 44 : END FUNCTION gauxc_cache_type_needs_gradient_grid
87 :
88 : ! **************************************************************************************************
89 : !> \brief ...
90 : !> \param this ...
91 : !> \return ...
92 : ! **************************************************************************************************
93 0 : FUNCTION gauxc_cache_type_max_l(this) RESULT(res)
94 : class(cp_gauxc_cache_type) :: this
95 : INTEGER :: res
96 :
97 0 : res = this%basisset%max_l
98 0 : END FUNCTION gauxc_cache_type_max_l
99 :
100 : ! **************************************************************************************************
101 : !> \brief ...
102 : !> \param cache ...
103 : !> \param params ...
104 : !> \param para_env ...
105 : !> \param particle_set ...
106 : !> \param qs_kind_set ...
107 : !> \param status ...
108 : ! **************************************************************************************************
109 378 : SUBROUTINE gauxc_cache_init( &
110 : cache, &
111 : params, &
112 : para_env, &
113 : particle_set, &
114 : qs_kind_set, &
115 : status)
116 :
117 : TYPE(cp_gauxc_cache_type), INTENT(inout) :: cache
118 : TYPE(cp_gauxc_cache_params), INTENT(inout) :: params
119 : TYPE(mp_para_env_type), INTENT(in), POINTER :: para_env
120 : TYPE(particle_type), DIMENSION(:), INTENT(in), &
121 : POINTER :: particle_set
122 : TYPE(qs_kind_type), DIMENSION(:), INTENT(in), &
123 : POINTER :: qs_kind_set
124 : TYPE(cp_gauxc_status_type), INTENT(inout) :: status
125 :
126 : #ifdef __GAUXC
127 378 : IF (gauxc_cache_is_valid( &
128 : cache, &
129 : params, &
130 : para_env, &
131 : particle_set)) THEN
132 334 : status%status%code = 0
133 334 : RETURN
134 : END IF
135 :
136 44 : CALL gauxc_cache_release(cache)
137 :
138 : cache%molecule = gauxc_create_molecule( &
139 : particle_set, &
140 44 : status)
141 44 : IF (status%status%code /= 0) THEN
142 0 : CALL cleanup_after_molecule()
143 0 : RETURN
144 : END IF
145 :
146 : cache%basisset = gauxc_create_basisset( &
147 : qs_kind_set, &
148 : particle_set, &
149 44 : status)
150 44 : IF (status%status%code /= 0) THEN
151 0 : CALL cleanup_after_basisset()
152 0 : RETURN
153 : END IF
154 :
155 : params%use_fd_gradient = &
156 44 : params%use_fd_gradient .AND. (cache%max_l() > 3)
157 : params%use_gradient_mpi_runtime = &
158 44 : params%use_gradient_mpi_runtime .AND. .NOT. params%use_fd_gradient
159 : params%use_gradient_self_runtime = &
160 44 : params%use_gradient_self_runtime .AND. .NOT. params%use_gradient_mpi_runtime
161 :
162 44 : IF (params%use_fd_gradient .AND. para_env%mepos == 0) THEN
163 : CALL cp_warn( &
164 : __LOCATION__, &
165 : "Using finite-difference GauXC XC gradients for METHOD GAPW/GAPW_XC with "// &
166 : "basis functions beyond f shells. The upstream analytical GauXC gradient path is "// &
167 0 : "not yet reliable for this case.")
168 : END IF
169 44 : IF (params%use_self_runtime .OR. .NOT. params%use_gauxc_model) THEN
170 : ! SKALA currently needs a replicated molecular runtime for reproducible
171 : ! open-shell densities across CP2K MPI ranks.
172 : ! Conventional GauXC also uses a replicated runtime here: CP2K has
173 : ! already allreduced the dense AO density matrix on every rank.
174 : cache%grid = gauxc_create_grid( &
175 : cache%molecule, &
176 : cache%basisset, &
177 : params%grid_type, &
178 : params%radial_quadrature, &
179 : params%pruning_scheme, &
180 : params%lb_exec_space, &
181 : params%batch_size, &
182 : params%device_runtime_fill_fraction, &
183 : status, &
184 : mpi_comm=mp_comm_self%get_handle(), &
185 30 : force_new_runtime=.TRUE.)
186 : ELSE
187 : ! Use the QS force-evaluation communicator rather than GauXC's global
188 : ! runtime. In mixed CDFT the diabatic states can be built on disjoint
189 : ! MPI subgroups; using the global communicator would make GauXC wait
190 : ! for ranks that are working on another state.
191 : cache%grid = gauxc_create_grid( &
192 : cache%molecule, &
193 : cache%basisset, &
194 : params%grid_type, &
195 : params%radial_quadrature, &
196 : params%pruning_scheme, &
197 : params%lb_exec_space, &
198 : params%batch_size, &
199 : params%device_runtime_fill_fraction, &
200 : status, &
201 14 : mpi_comm=para_env%get_handle())
202 : END IF
203 44 : IF (status%status%code /= 0) THEN
204 0 : CALL cleanup_after_grid()
205 0 : RETURN
206 : END IF
207 :
208 : cache%integrator = gauxc_create_integrator( &
209 : TRIM(params%xc_fun_name), &
210 : cache%grid, &
211 : params%int_exec_space, &
212 : params%lwd_kernel, &
213 : params%nspins, &
214 44 : status)
215 44 : IF (status%status%code /= 0) THEN
216 0 : CALL cleanup_after_integrator()
217 0 : RETURN
218 : END IF
219 :
220 44 : IF (params%use_gradient_self_runtime) THEN
221 : ! Upstream GauXC does not yet support OneDFT/SKALA nuclear gradients
222 : ! on an MPI runtime. Keep the energy/VXC path on the normal MPI
223 : ! runtime and use an isolated runtime only for the replicated gradient.
224 : cache%gradient_grid = gauxc_create_grid( &
225 : cache%molecule, &
226 : cache%basisset, &
227 : params%grid_type, &
228 : params%radial_quadrature, &
229 : params%pruning_scheme, &
230 : params%lb_exec_space, &
231 : params%batch_size, &
232 : params%device_runtime_fill_fraction, &
233 : status, &
234 : mpi_comm=mp_comm_self%get_handle(), &
235 0 : force_new_runtime=.TRUE.)
236 0 : IF (status%status%code /= 0) THEN
237 0 : CALL cleanup_after_gradient_grid()
238 0 : RETURN
239 : END IF
240 : cache%gradient_integrator = gauxc_create_integrator( &
241 : TRIM(params%xc_fun_name), &
242 : cache%gradient_grid, &
243 : params%int_exec_space, &
244 : params%lwd_kernel, &
245 : params%nspins, &
246 0 : status)
247 0 : IF (status%status%code /= 0) THEN
248 0 : CALL cleanup_after_gradient_integrator()
249 0 : RETURN
250 : END IF
251 : END IF
252 :
253 44 : cache%natom = params%natom
254 44 : cache%nspins = params%nspins
255 44 : cache%batch_size = params%batch_size
256 44 : cache%device_runtime_fill_fraction = REAL(params%device_runtime_fill_fraction, c_double)
257 44 : cache%xc_fun_name = TRIM(params%xc_fun_name)
258 44 : cache%grid_type = TRIM(params%grid_type)
259 44 : cache%radial_quadrature = TRIM(params%radial_quadrature)
260 44 : cache%pruning_scheme = TRIM(params%pruning_scheme)
261 44 : cache%lb_exec_space = TRIM(params%lb_exec_space)
262 44 : cache%int_exec_space = TRIM(params%int_exec_space)
263 44 : cache%model_eval_name = TRIM(params%model_eval_name)
264 44 : cache%lwd_kernel = TRIM(params%lwd_kernel)
265 44 : cache%use_gradient_mpi_runtime = params%use_gradient_mpi_runtime
266 44 : cache%use_mpi_runtime = params%use_mpi_runtime
267 44 : cache%use_gradient_self_runtime = params%use_gradient_self_runtime
268 44 : cache%use_self_runtime = params%use_self_runtime
269 44 : cache%use_fd_gradient = params%use_fd_gradient
270 44 : IF (ALLOCATED(cache%last_positions)) THEN
271 0 : DEALLOCATE (cache%last_positions)
272 : END IF
273 132 : ALLOCATE (cache%last_positions(3, params%natom))
274 :
275 : BLOCK
276 : INTEGER :: iatom
277 164 : DO iatom = 1, params%natom
278 524 : cache%last_positions(:, iatom) = REAL(particle_set(iatom)%r(:), c_double)
279 : END DO
280 : END BLOCK
281 :
282 44 : cache%mpi_comm = para_env%get_handle()
283 44 : cache%is_init = .TRUE.
284 :
285 : #else
286 : MARK_USED(params)
287 : MARK_USED(para_env)
288 : MARK_USED(particle_set)
289 : MARK_USED(qs_kind_set)
290 : MARK_USED(status)
291 : cache%is_init = .FALSE.
292 : #endif
293 :
294 : CONTAINS
295 : ! **************************************************************************************************
296 : !> \brief ...
297 : ! **************************************************************************************************
298 0 : SUBROUTINE cleanup_after_gradient_integrator()
299 0 : CALL gauxc_destroy_integrator(cache%gradient_integrator, status)
300 0 : CALL gauxc_check_status(status)
301 0 : CALL cleanup_after_gradient_grid()
302 0 : END SUBROUTINE cleanup_after_gradient_integrator
303 :
304 : ! **************************************************************************************************
305 : !> \brief ...
306 : ! **************************************************************************************************
307 0 : SUBROUTINE cleanup_after_gradient_grid()
308 0 : CALL gauxc_destroy_grid(cache%gradient_grid, status)
309 0 : CALL gauxc_check_status(status)
310 0 : CALL cleanup_after_integrator()
311 0 : END SUBROUTINE cleanup_after_gradient_grid
312 :
313 : ! **************************************************************************************************
314 : !> \brief ...
315 : ! **************************************************************************************************
316 0 : SUBROUTINE cleanup_after_integrator()
317 0 : CALL gauxc_destroy_integrator(cache%integrator, status)
318 0 : CALL gauxc_check_status(status)
319 0 : CALL cleanup_after_grid()
320 0 : END SUBROUTINE cleanup_after_integrator
321 :
322 : ! **************************************************************************************************
323 : !> \brief ...
324 : ! **************************************************************************************************
325 0 : SUBROUTINE cleanup_after_grid()
326 0 : CALL gauxc_destroy_grid(cache%grid, status)
327 0 : CALL gauxc_check_status(status)
328 0 : CALL cleanup_after_basisset()
329 0 : END SUBROUTINE cleanup_after_grid
330 :
331 : ! **************************************************************************************************
332 : !> \brief ...
333 : ! **************************************************************************************************
334 0 : SUBROUTINE cleanup_after_basisset()
335 0 : CALL gauxc_destroy_basisset(cache%basisset, status)
336 0 : CALL gauxc_check_status(status)
337 0 : CALL cleanup_after_molecule()
338 0 : END SUBROUTINE cleanup_after_basisset
339 :
340 : ! **************************************************************************************************
341 : !> \brief ...
342 : ! **************************************************************************************************
343 0 : SUBROUTINE cleanup_after_molecule()
344 0 : CALL gauxc_destroy_molecule(cache%molecule, status)
345 0 : CALL gauxc_check_status(status)
346 0 : CALL cleanup_final()
347 0 : END SUBROUTINE cleanup_after_molecule
348 :
349 : ! **************************************************************************************************
350 : !> \brief ...
351 : ! **************************************************************************************************
352 0 : SUBROUTINE cleanup_final()
353 0 : IF (ALLOCATED(cache%last_positions)) DEALLOCATE (cache%last_positions)
354 0 : cache%is_init = .FALSE.
355 0 : END SUBROUTINE cleanup_final
356 :
357 : END SUBROUTINE gauxc_cache_init
358 :
359 : ! **************************************************************************************************
360 : !> \brief Release all GauXC objects in a cache and deallocate arrays
361 : !> \param cache ...
362 : ! **************************************************************************************************
363 80 : SUBROUTINE gauxc_cache_release(cache)
364 : TYPE(cp_gauxc_cache_type), INTENT(INOUT) :: cache
365 :
366 : #ifdef __GAUXC
367 : TYPE(cp_gauxc_status_type) :: status
368 :
369 80 : IF (cache%is_init) THEN
370 44 : IF (cache%needs_gradient_grid()) THEN
371 0 : CALL gauxc_destroy_integrator(cache%gradient_integrator, status)
372 0 : CALL gauxc_check_status(status)
373 0 : CALL gauxc_destroy_grid(cache%gradient_grid, status)
374 0 : CALL gauxc_check_status(status)
375 : END IF
376 44 : CALL gauxc_destroy_integrator(cache%integrator, status)
377 44 : CALL gauxc_check_status(status)
378 44 : CALL gauxc_destroy_grid(cache%grid, status)
379 44 : CALL gauxc_check_status(status)
380 44 : CALL gauxc_destroy_basisset(cache%basisset, status)
381 44 : CALL gauxc_check_status(status)
382 44 : CALL gauxc_destroy_molecule(cache%molecule, status)
383 44 : CALL gauxc_check_status(status)
384 : END IF
385 80 : IF (ALLOCATED(cache%last_positions)) DEALLOCATE (cache%last_positions)
386 80 : cache%is_init = .FALSE.
387 : #else
388 : MARK_USED(cache)
389 : IF (ALLOCATED(cache%last_positions)) DEALLOCATE (cache%last_positions)
390 : cache%is_init = .FALSE.
391 : #endif
392 80 : END SUBROUTINE gauxc_cache_release
393 :
394 : ! **************************************************************************************************
395 : !> \brief Check if the GauXC cache is valid for the given parameters
396 : !> \param cache ...
397 : !> \param params ...
398 : !> \param para_env ...
399 : !> \param particle_set ...
400 : !> \return ...
401 : !> \retval cache_valid ...
402 : ! **************************************************************************************************
403 378 : FUNCTION gauxc_cache_is_valid( &
404 : cache, &
405 : params, &
406 : para_env, &
407 : particle_set) &
408 : RESULT(cache_valid)
409 :
410 : TYPE(cp_gauxc_cache_type), INTENT(IN) :: cache
411 : TYPE(cp_gauxc_cache_params), INTENT(in) :: params
412 : TYPE(mp_para_env_type), INTENT(in), POINTER :: para_env
413 : TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
414 : LOGICAL :: use_gradient_self_runtime, use_gradient_mpi_runtime, use_fd_gradient
415 : LOGICAL :: cache_valid
416 :
417 : #ifdef __GAUXC
418 : INTEGER :: iatom
419 : REAL(c_double), PARAMETER :: pos_tol = 1.0E-8_c_double
420 :
421 378 : cache_valid = .FALSE.
422 378 : IF (.NOT. cache%is_init) RETURN
423 342 : IF (cache%natom /= params%natom) RETURN
424 342 : IF (cache%nspins /= params%nspins) RETURN
425 342 : IF (cache%batch_size /= params%batch_size) RETURN
426 342 : IF (ABS(cache%device_runtime_fill_fraction - REAL(params%device_runtime_fill_fraction, c_double)) > &
427 : 1.0E-10_c_double) RETURN
428 342 : IF (TRIM(cache%xc_fun_name) /= TRIM(params%xc_fun_name)) RETURN
429 342 : IF (TRIM(cache%grid_type) /= TRIM(params%grid_type)) RETURN
430 342 : IF (TRIM(cache%radial_quadrature) /= TRIM(params%radial_quadrature)) RETURN
431 342 : IF (TRIM(cache%pruning_scheme) /= TRIM(params%pruning_scheme)) RETURN
432 342 : IF (TRIM(cache%lb_exec_space) /= TRIM(params%lb_exec_space)) RETURN
433 342 : IF (TRIM(cache%int_exec_space) /= TRIM(params%int_exec_space)) RETURN
434 342 : IF (TRIM(cache%model_eval_name) /= TRIM(params%model_eval_name)) RETURN
435 342 : IF (TRIM(cache%lwd_kernel) /= TRIM(params%lwd_kernel)) RETURN
436 342 : IF (cache%use_self_runtime .NEQV. params%use_self_runtime) RETURN
437 342 : IF (cache%use_mpi_runtime .NEQV. params%use_mpi_runtime) RETURN
438 :
439 : ! These three are modified during cache init, so we must emulate these changes here
440 342 : use_gradient_self_runtime = params%use_gradient_self_runtime
441 342 : use_gradient_mpi_runtime = params%use_gradient_mpi_runtime
442 342 : use_fd_gradient = params%use_fd_gradient
443 :
444 : use_fd_gradient = &
445 342 : use_fd_gradient .AND. (cache%max_l() > 3)
446 : use_gradient_mpi_runtime = &
447 342 : use_gradient_mpi_runtime .AND. .NOT. use_fd_gradient
448 : use_gradient_self_runtime = &
449 342 : use_gradient_self_runtime .AND. .NOT. use_gradient_mpi_runtime
450 :
451 342 : IF (cache%use_gradient_mpi_runtime .NEQV. use_gradient_mpi_runtime) RETURN
452 342 : IF (cache%use_fd_gradient .NEQV. use_fd_gradient) RETURN
453 342 : IF (cache%use_gradient_self_runtime .NEQV. use_gradient_self_runtime) RETURN
454 342 : IF (cache%mpi_comm /= para_env%get_handle()) RETURN
455 :
456 342 : IF (.NOT. ALLOCATED(cache%last_positions)) RETURN
457 342 : IF (SIZE(cache%last_positions, 2) /= params%natom) RETURN
458 1202 : DO iatom = 1, params%natom
459 3464 : IF (ANY(ABS(cache%last_positions(:, iatom) - &
460 334 : REAL(particle_set(iatom)%r(:), c_double)) > pos_tol)) RETURN
461 : END DO
462 378 : cache_valid = .TRUE.
463 : #else
464 : MARK_USED(cache)
465 : MARK_USED(particle_set)
466 : MARK_USED(params)
467 : MARK_USED(use_gradient_self_runtime)
468 : MARK_USED(use_gradient_mpi_runtime)
469 : MARK_USED(use_fd_gradient)
470 : MARK_USED(para_env)
471 : cache_valid = .FALSE.
472 : #endif
473 : END FUNCTION gauxc_cache_is_valid
474 :
475 0 : END MODULE xc_gauxc_cache
|
