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: BSD-3-Clause !
6 : !--------------------------------------------------------------------------------------------------!
7 :
8 : ! **************************************************************************************************
9 : !> \brief Fortran API for the grid package, which is written in C.
10 : !> \author Ole Schuett
11 : ! **************************************************************************************************
12 : MODULE grid_api
13 : USE ISO_C_BINDING, ONLY: &
14 : C_ASSOCIATED, C_BOOL, C_CHAR, C_DOUBLE, C_FUNLOC, C_FUNPTR, C_INT, C_LOC, C_LONG, &
15 : C_NULL_PTR, C_PTR
16 : USE kinds, ONLY: dp
17 : USE message_passing, ONLY: mp_comm_type
18 : USE offload_api, ONLY: offload_buffer_type
19 : USE realspace_grid_types, ONLY: realspace_grid_type
20 : USE string_utilities, ONLY: strlcpy_c2f
21 : #include "../base/base_uses.f90"
22 :
23 : IMPLICIT NONE
24 :
25 : PRIVATE
26 :
27 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'grid_api'
28 :
29 : INTEGER, PARAMETER, PUBLIC :: GRID_FUNC_AB = 100
30 : INTEGER, PARAMETER, PUBLIC :: GRID_FUNC_DADB = 200
31 : INTEGER, PARAMETER, PUBLIC :: GRID_FUNC_ADBmDAB_X = 301
32 : INTEGER, PARAMETER, PUBLIC :: GRID_FUNC_ADBmDAB_Y = 302
33 : INTEGER, PARAMETER, PUBLIC :: GRID_FUNC_ADBmDAB_Z = 303
34 : INTEGER, PARAMETER, PUBLIC :: GRID_FUNC_ARDBmDARB_XX = 411
35 : INTEGER, PARAMETER, PUBLIC :: GRID_FUNC_ARDBmDARB_XY = 412
36 : INTEGER, PARAMETER, PUBLIC :: GRID_FUNC_ARDBmDARB_XZ = 413
37 : INTEGER, PARAMETER, PUBLIC :: GRID_FUNC_ARDBmDARB_YX = 421
38 : INTEGER, PARAMETER, PUBLIC :: GRID_FUNC_ARDBmDARB_YY = 422
39 : INTEGER, PARAMETER, PUBLIC :: GRID_FUNC_ARDBmDARB_YZ = 423
40 : INTEGER, PARAMETER, PUBLIC :: GRID_FUNC_ARDBmDARB_ZX = 431
41 : INTEGER, PARAMETER, PUBLIC :: GRID_FUNC_ARDBmDARB_ZY = 432
42 : INTEGER, PARAMETER, PUBLIC :: GRID_FUNC_ARDBmDARB_ZZ = 433
43 : INTEGER, PARAMETER, PUBLIC :: GRID_FUNC_DABpADB_X = 501
44 : INTEGER, PARAMETER, PUBLIC :: GRID_FUNC_DABpADB_Y = 502
45 : INTEGER, PARAMETER, PUBLIC :: GRID_FUNC_DABpADB_Z = 503
46 : INTEGER, PARAMETER, PUBLIC :: GRID_FUNC_DX = 601
47 : INTEGER, PARAMETER, PUBLIC :: GRID_FUNC_DY = 602
48 : INTEGER, PARAMETER, PUBLIC :: GRID_FUNC_DZ = 603
49 : INTEGER, PARAMETER, PUBLIC :: GRID_FUNC_DXDY = 701
50 : INTEGER, PARAMETER, PUBLIC :: GRID_FUNC_DYDZ = 702
51 : INTEGER, PARAMETER, PUBLIC :: GRID_FUNC_DZDX = 703
52 : INTEGER, PARAMETER, PUBLIC :: GRID_FUNC_DXDX = 801
53 : INTEGER, PARAMETER, PUBLIC :: GRID_FUNC_DYDY = 802
54 : INTEGER, PARAMETER, PUBLIC :: GRID_FUNC_DZDZ = 803
55 : INTEGER, PARAMETER, PUBLIC :: GRID_FUNC_DAB_X = 901
56 : INTEGER, PARAMETER, PUBLIC :: GRID_FUNC_DAB_Y = 902
57 : INTEGER, PARAMETER, PUBLIC :: GRID_FUNC_DAB_Z = 903
58 : INTEGER, PARAMETER, PUBLIC :: GRID_FUNC_ADB_X = 904
59 : INTEGER, PARAMETER, PUBLIC :: GRID_FUNC_ADB_Y = 905
60 : INTEGER, PARAMETER, PUBLIC :: GRID_FUNC_ADB_Z = 906
61 :
62 : INTEGER, PARAMETER, PUBLIC :: GRID_FUNC_CORE_X = 1001
63 : INTEGER, PARAMETER, PUBLIC :: GRID_FUNC_CORE_Y = 1002
64 : INTEGER, PARAMETER, PUBLIC :: GRID_FUNC_CORE_Z = 1003
65 :
66 : INTEGER, PARAMETER, PUBLIC :: GRID_BACKEND_AUTO = 10
67 : INTEGER, PARAMETER, PUBLIC :: GRID_BACKEND_REF = 11
68 : INTEGER, PARAMETER, PUBLIC :: GRID_BACKEND_CPU = 12
69 : INTEGER, PARAMETER, PUBLIC :: GRID_BACKEND_DGEMM = 13
70 : INTEGER, PARAMETER, PUBLIC :: GRID_BACKEND_GPU = 14
71 : INTEGER, PARAMETER, PUBLIC :: GRID_BACKEND_HIP = 15
72 :
73 : PUBLIC :: grid_library_init, grid_library_finalize
74 : PUBLIC :: grid_library_set_config, grid_library_print_stats
75 : PUBLIC :: collocate_pgf_product, integrate_pgf_product
76 : PUBLIC :: grid_basis_set_type, grid_create_basis_set, grid_free_basis_set
77 : PUBLIC :: grid_task_list_type, grid_create_task_list, grid_free_task_list
78 : PUBLIC :: grid_collocate_task_list, grid_integrate_task_list
79 :
80 : TYPE grid_basis_set_type
81 : PRIVATE
82 : TYPE(C_PTR) :: c_ptr = C_NULL_PTR
83 : END TYPE grid_basis_set_type
84 :
85 : TYPE grid_task_list_type
86 : PRIVATE
87 : TYPE(C_PTR) :: c_ptr = C_NULL_PTR
88 : END TYPE grid_task_list_type
89 :
90 : CONTAINS
91 :
92 : ! **************************************************************************************************
93 : !> \brief low level collocation of primitive gaussian functions
94 : !> \param la_max ...
95 : !> \param zeta ...
96 : !> \param la_min ...
97 : !> \param lb_max ...
98 : !> \param zetb ...
99 : !> \param lb_min ...
100 : !> \param ra ...
101 : !> \param rab ...
102 : !> \param scale ...
103 : !> \param pab ...
104 : !> \param o1 ...
105 : !> \param o2 ...
106 : !> \param rsgrid ...
107 : !> \param ga_gb_function ...
108 : !> \param radius ...
109 : !> \param use_subpatch ...
110 : !> \param subpatch_pattern ...
111 : !> \author Ole Schuett
112 : ! **************************************************************************************************
113 1144613 : SUBROUTINE collocate_pgf_product(la_max, zeta, la_min, &
114 : lb_max, zetb, lb_min, &
115 : ra, rab, scale, pab, o1, o2, &
116 : rsgrid, &
117 : ga_gb_function, radius, &
118 : use_subpatch, subpatch_pattern)
119 :
120 : INTEGER, INTENT(IN) :: la_max
121 : REAL(KIND=dp), INTENT(IN) :: zeta
122 : INTEGER, INTENT(IN) :: la_min, lb_max
123 : REAL(KIND=dp), INTENT(IN) :: zetb
124 : INTEGER, INTENT(IN) :: lb_min
125 : REAL(KIND=dp), DIMENSION(3), INTENT(IN), TARGET :: ra, rab
126 : REAL(KIND=dp), INTENT(IN) :: scale
127 : REAL(KIND=dp), DIMENSION(:, :), POINTER :: pab
128 : INTEGER, INTENT(IN) :: o1, o2
129 : TYPE(realspace_grid_type) :: rsgrid
130 : INTEGER, INTENT(IN) :: ga_gb_function
131 : REAL(KIND=dp), INTENT(IN) :: radius
132 : LOGICAL, OPTIONAL :: use_subpatch
133 : INTEGER, INTENT(IN), OPTIONAL :: subpatch_pattern
134 :
135 : INTEGER :: border_mask
136 : INTEGER, DIMENSION(3), TARGET :: border_width, npts_global, npts_local, &
137 : shift_local
138 : LOGICAL(KIND=C_BOOL) :: orthorhombic
139 1144613 : REAL(KIND=dp), DIMENSION(:, :, :), POINTER :: grid
140 : INTERFACE
141 : SUBROUTINE grid_cpu_collocate_pgf_product_c(orthorhombic, &
142 : border_mask, func, &
143 : la_max, la_min, lb_max, lb_min, &
144 : zeta, zetb, rscale, dh, dh_inv, ra, rab, &
145 : npts_global, npts_local, shift_local, border_width, &
146 : radius, o1, o2, n1, n2, pab, &
147 : grid) &
148 : BIND(C, name="grid_cpu_collocate_pgf_product")
149 : IMPORT :: C_PTR, C_INT, C_DOUBLE, C_BOOL
150 : LOGICAL(KIND=C_BOOL), VALUE :: orthorhombic
151 : INTEGER(KIND=C_INT), VALUE :: border_mask
152 : INTEGER(KIND=C_INT), VALUE :: func
153 : INTEGER(KIND=C_INT), VALUE :: la_max
154 : INTEGER(KIND=C_INT), VALUE :: la_min
155 : INTEGER(KIND=C_INT), VALUE :: lb_max
156 : INTEGER(KIND=C_INT), VALUE :: lb_min
157 : REAL(KIND=C_DOUBLE), VALUE :: zeta
158 : REAL(KIND=C_DOUBLE), VALUE :: zetb
159 : REAL(KIND=C_DOUBLE), VALUE :: rscale
160 : TYPE(C_PTR), VALUE :: dh
161 : TYPE(C_PTR), VALUE :: dh_inv
162 : TYPE(C_PTR), VALUE :: ra
163 : TYPE(C_PTR), VALUE :: rab
164 : TYPE(C_PTR), VALUE :: npts_global
165 : TYPE(C_PTR), VALUE :: npts_local
166 : TYPE(C_PTR), VALUE :: shift_local
167 : TYPE(C_PTR), VALUE :: border_width
168 : REAL(KIND=C_DOUBLE), VALUE :: radius
169 : INTEGER(KIND=C_INT), VALUE :: o1
170 : INTEGER(KIND=C_INT), VALUE :: o2
171 : INTEGER(KIND=C_INT), VALUE :: n1
172 : INTEGER(KIND=C_INT), VALUE :: n2
173 : TYPE(C_PTR), VALUE :: pab
174 : TYPE(C_PTR), VALUE :: grid
175 : END SUBROUTINE grid_cpu_collocate_pgf_product_c
176 : END INTERFACE
177 :
178 1144613 : border_mask = 0
179 1144613 : IF (PRESENT(use_subpatch)) THEN
180 65606 : IF (use_subpatch) THEN
181 54770 : CPASSERT(PRESENT(subpatch_pattern))
182 54770 : border_mask = IAND(63, NOT(subpatch_pattern)) ! invert last 6 bits
183 : END IF
184 : END IF
185 :
186 1144613 : orthorhombic = LOGICAL(rsgrid%desc%orthorhombic, C_BOOL)
187 :
188 1144613 : CPASSERT(LBOUND(pab, 1) == 1)
189 1144613 : CPASSERT(LBOUND(pab, 2) == 1)
190 :
191 : CALL get_rsgrid_properties(rsgrid, npts_global=npts_global, &
192 : npts_local=npts_local, &
193 : shift_local=shift_local, &
194 1144613 : border_width=border_width)
195 :
196 1144613 : grid(1:, 1:, 1:) => rsgrid%r(:, :, :) ! pointer assignment
197 :
198 : #if __GNUC__ >= 9
199 : CPASSERT(IS_CONTIGUOUS(rsgrid%desc%dh))
200 : CPASSERT(IS_CONTIGUOUS(rsgrid%desc%dh_inv))
201 : CPASSERT(IS_CONTIGUOUS(ra))
202 : CPASSERT(IS_CONTIGUOUS(rab))
203 : CPASSERT(IS_CONTIGUOUS(npts_global))
204 : CPASSERT(IS_CONTIGUOUS(npts_local))
205 : CPASSERT(IS_CONTIGUOUS(shift_local))
206 : CPASSERT(IS_CONTIGUOUS(border_width))
207 1144613 : CPASSERT(IS_CONTIGUOUS(pab))
208 1144613 : CPASSERT(IS_CONTIGUOUS(grid))
209 : #endif
210 :
211 : ! For collocating a single pgf product we use the optimized cpu backend.
212 :
213 : CALL grid_cpu_collocate_pgf_product_c(orthorhombic=orthorhombic, &
214 : border_mask=border_mask, &
215 : func=ga_gb_function, &
216 : la_max=la_max, &
217 : la_min=la_min, &
218 : lb_max=lb_max, &
219 : lb_min=lb_min, &
220 : zeta=zeta, &
221 : zetb=zetb, &
222 : rscale=scale, &
223 : dh=C_LOC(rsgrid%desc%dh(1, 1)), &
224 : dh_inv=C_LOC(rsgrid%desc%dh_inv(1, 1)), &
225 : ra=C_LOC(ra(1)), &
226 : rab=C_LOC(rab(1)), &
227 : npts_global=C_LOC(npts_global(1)), &
228 : npts_local=C_LOC(npts_local(1)), &
229 : shift_local=C_LOC(shift_local(1)), &
230 : border_width=C_LOC(border_width(1)), &
231 : radius=radius, &
232 : o1=o1, &
233 : o2=o2, &
234 : n1=SIZE(pab, 1), &
235 : n2=SIZE(pab, 2), &
236 : pab=C_LOC(pab(1, 1)), &
237 1144613 : grid=C_LOC(grid(1, 1, 1)))
238 :
239 1144613 : END SUBROUTINE collocate_pgf_product
240 :
241 : ! **************************************************************************************************
242 : !> \brief low level function to compute matrix elements of primitive gaussian functions
243 : !> \param la_max ...
244 : !> \param zeta ...
245 : !> \param la_min ...
246 : !> \param lb_max ...
247 : !> \param zetb ...
248 : !> \param lb_min ...
249 : !> \param ra ...
250 : !> \param rab ...
251 : !> \param rsgrid ...
252 : !> \param hab ...
253 : !> \param pab ...
254 : !> \param o1 ...
255 : !> \param o2 ...
256 : !> \param radius ...
257 : !> \param calculate_forces ...
258 : !> \param force_a ...
259 : !> \param force_b ...
260 : !> \param compute_tau ...
261 : !> \param use_virial ...
262 : !> \param my_virial_a ...
263 : !> \param my_virial_b ...
264 : !> \param hdab Derivative with respect to the primitive on the left.
265 : !> \param hadb Derivative with respect to the primitive on the right.
266 : !> \param a_hdab ...
267 : !> \param use_subpatch ...
268 : !> \param subpatch_pattern ...
269 : ! **************************************************************************************************
270 882452 : SUBROUTINE integrate_pgf_product(la_max, zeta, la_min, &
271 : lb_max, zetb, lb_min, &
272 : ra, rab, rsgrid, &
273 : hab, pab, o1, o2, &
274 : radius, &
275 : calculate_forces, force_a, force_b, &
276 : compute_tau, &
277 : use_virial, my_virial_a, &
278 : my_virial_b, hdab, hadb, a_hdab, use_subpatch, subpatch_pattern)
279 :
280 : INTEGER, INTENT(IN) :: la_max
281 : REAL(KIND=dp), INTENT(IN) :: zeta
282 : INTEGER, INTENT(IN) :: la_min, lb_max
283 : REAL(KIND=dp), INTENT(IN) :: zetb
284 : INTEGER, INTENT(IN) :: lb_min
285 : REAL(KIND=dp), DIMENSION(3), INTENT(IN), TARGET :: ra, rab
286 : TYPE(realspace_grid_type), INTENT(IN) :: rsgrid
287 : REAL(KIND=dp), DIMENSION(:, :), POINTER :: hab
288 : REAL(KIND=dp), DIMENSION(:, :), OPTIONAL, POINTER :: pab
289 : INTEGER, INTENT(IN) :: o1, o2
290 : REAL(KIND=dp), INTENT(IN) :: radius
291 : LOGICAL, INTENT(IN) :: calculate_forces
292 : REAL(KIND=dp), DIMENSION(3), INTENT(INOUT), &
293 : OPTIONAL :: force_a, force_b
294 : LOGICAL, INTENT(IN), OPTIONAL :: compute_tau, use_virial
295 : REAL(KIND=dp), DIMENSION(3, 3), OPTIONAL :: my_virial_a, my_virial_b
296 : REAL(KIND=dp), DIMENSION(:, :, :), OPTIONAL, &
297 : POINTER :: hdab, hadb
298 : REAL(KIND=dp), DIMENSION(:, :, :, :), OPTIONAL, &
299 : POINTER :: a_hdab
300 : LOGICAL, OPTIONAL :: use_subpatch
301 : INTEGER, INTENT(IN), OPTIONAL :: subpatch_pattern
302 :
303 : INTEGER :: border_mask
304 : INTEGER, DIMENSION(3), TARGET :: border_width, npts_global, npts_local, &
305 : shift_local
306 : LOGICAL :: my_use_virial
307 : LOGICAL(KIND=C_BOOL) :: my_compute_tau, orthorhombic
308 : REAL(KIND=dp), DIMENSION(3, 2), TARGET :: forces
309 : REAL(KIND=dp), DIMENSION(3, 3, 2), TARGET :: virials
310 882452 : REAL(KIND=dp), DIMENSION(:, :, :), POINTER :: grid
311 : TYPE(C_PTR) :: a_hdab_cptr, forces_cptr, hadb_cptr, &
312 : hdab_cptr, pab_cptr, virials_cptr
313 : INTERFACE
314 : SUBROUTINE grid_cpu_integrate_pgf_product_c(orthorhombic, compute_tau, &
315 : border_mask, &
316 : la_max, la_min, lb_max, lb_min, &
317 : zeta, zetb, dh, dh_inv, ra, rab, &
318 : npts_global, npts_local, shift_local, border_width, &
319 : radius, o1, o2, n1, n2, grid, hab, pab, &
320 : forces, virials, hdab, hadb, a_hdab) &
321 : BIND(C, name="grid_cpu_integrate_pgf_product")
322 : IMPORT :: C_PTR, C_INT, C_DOUBLE, C_BOOL
323 : LOGICAL(KIND=C_BOOL), VALUE :: orthorhombic
324 : LOGICAL(KIND=C_BOOL), VALUE :: compute_tau
325 : INTEGER(KIND=C_INT), VALUE :: border_mask
326 : INTEGER(KIND=C_INT), VALUE :: la_max
327 : INTEGER(KIND=C_INT), VALUE :: la_min
328 : INTEGER(KIND=C_INT), VALUE :: lb_max
329 : INTEGER(KIND=C_INT), VALUE :: lb_min
330 : REAL(KIND=C_DOUBLE), VALUE :: zeta
331 : REAL(KIND=C_DOUBLE), VALUE :: zetb
332 : TYPE(C_PTR), VALUE :: dh
333 : TYPE(C_PTR), VALUE :: dh_inv
334 : TYPE(C_PTR), VALUE :: ra
335 : TYPE(C_PTR), VALUE :: rab
336 : TYPE(C_PTR), VALUE :: npts_global
337 : TYPE(C_PTR), VALUE :: npts_local
338 : TYPE(C_PTR), VALUE :: shift_local
339 : TYPE(C_PTR), VALUE :: border_width
340 : REAL(KIND=C_DOUBLE), VALUE :: radius
341 : INTEGER(KIND=C_INT), VALUE :: o1
342 : INTEGER(KIND=C_INT), VALUE :: o2
343 : INTEGER(KIND=C_INT), VALUE :: n1
344 : INTEGER(KIND=C_INT), VALUE :: n2
345 : TYPE(C_PTR), VALUE :: grid
346 : TYPE(C_PTR), VALUE :: hab
347 : TYPE(C_PTR), VALUE :: pab
348 : TYPE(C_PTR), VALUE :: forces
349 : TYPE(C_PTR), VALUE :: virials
350 : TYPE(C_PTR), VALUE :: hdab
351 : TYPE(C_PTR), VALUE :: hadb
352 : TYPE(C_PTR), VALUE :: a_hdab
353 : END SUBROUTINE grid_cpu_integrate_pgf_product_c
354 : END INTERFACE
355 :
356 882452 : IF (radius == 0.0_dp) THEN
357 0 : RETURN
358 : END IF
359 :
360 882452 : border_mask = 0
361 882452 : IF (PRESENT(use_subpatch)) THEN
362 836419 : IF (use_subpatch) THEN
363 48560 : CPASSERT(PRESENT(subpatch_pattern))
364 48560 : border_mask = IAND(63, NOT(subpatch_pattern)) ! invert last 6 bits
365 : END IF
366 : END IF
367 :
368 : ! When true then 0.5 * (nabla x_a).(v(r) nabla x_b) is computed.
369 882452 : IF (PRESENT(compute_tau)) THEN
370 3612 : my_compute_tau = LOGICAL(compute_tau, C_BOOL)
371 : ELSE
372 : my_compute_tau = .FALSE.
373 : END IF
374 :
375 882452 : IF (PRESENT(use_virial)) THEN
376 386575 : my_use_virial = use_virial
377 : ELSE
378 : my_use_virial = .FALSE.
379 : END IF
380 :
381 882452 : IF (calculate_forces) THEN
382 355207 : CPASSERT(PRESENT(pab))
383 355207 : pab_cptr = C_LOC(pab(1, 1))
384 355207 : forces(:, :) = 0.0_dp
385 355207 : forces_cptr = C_LOC(forces(1, 1))
386 : ELSE
387 : pab_cptr = C_NULL_PTR
388 : forces_cptr = C_NULL_PTR
389 : END IF
390 :
391 882452 : IF (calculate_forces .AND. my_use_virial) THEN
392 80078 : virials(:, :, :) = 0.0_dp
393 80078 : virials_cptr = C_LOC(virials(1, 1, 1))
394 : ELSE
395 : virials_cptr = C_NULL_PTR
396 : END IF
397 :
398 882452 : IF (calculate_forces .AND. PRESENT(hdab)) THEN
399 2714 : hdab_cptr = C_LOC(hdab(1, 1, 1))
400 : ELSE
401 : hdab_cptr = C_NULL_PTR
402 : END IF
403 :
404 882452 : IF (calculate_forces .AND. PRESENT(hadb)) THEN
405 1806 : hadb_cptr = C_LOC(hadb(1, 1, 1))
406 : ELSE
407 : hadb_cptr = C_NULL_PTR
408 : END IF
409 :
410 882452 : IF (calculate_forces .AND. my_use_virial .AND. PRESENT(a_hdab)) THEN
411 36 : a_hdab_cptr = C_LOC(a_hdab(1, 1, 1, 1))
412 : ELSE
413 : a_hdab_cptr = C_NULL_PTR
414 : END IF
415 :
416 882452 : orthorhombic = LOGICAL(rsgrid%desc%orthorhombic, C_BOOL)
417 :
418 : CALL get_rsgrid_properties(rsgrid, npts_global=npts_global, &
419 : npts_local=npts_local, &
420 : shift_local=shift_local, &
421 882452 : border_width=border_width)
422 :
423 882452 : grid(1:, 1:, 1:) => rsgrid%r(:, :, :) ! pointer assignment
424 :
425 : #if __GNUC__ >= 9
426 : CPASSERT(IS_CONTIGUOUS(rsgrid%desc%dh))
427 : CPASSERT(IS_CONTIGUOUS(rsgrid%desc%dh_inv))
428 : CPASSERT(IS_CONTIGUOUS(ra))
429 : CPASSERT(IS_CONTIGUOUS(rab))
430 : CPASSERT(IS_CONTIGUOUS(npts_global))
431 : CPASSERT(IS_CONTIGUOUS(npts_local))
432 : CPASSERT(IS_CONTIGUOUS(shift_local))
433 : CPASSERT(IS_CONTIGUOUS(border_width))
434 882452 : CPASSERT(IS_CONTIGUOUS(grid))
435 882452 : CPASSERT(IS_CONTIGUOUS(hab))
436 : CPASSERT(IS_CONTIGUOUS(forces))
437 : CPASSERT(IS_CONTIGUOUS(virials))
438 882452 : IF (PRESENT(pab)) THEN
439 390210 : CPASSERT(IS_CONTIGUOUS(pab))
440 : END IF
441 882452 : IF (PRESENT(hdab)) THEN
442 24835 : CPASSERT(IS_CONTIGUOUS(hdab))
443 : END IF
444 882452 : IF (PRESENT(a_hdab)) THEN
445 23029 : CPASSERT(IS_CONTIGUOUS(a_hdab))
446 : END IF
447 : #endif
448 :
449 : CALL grid_cpu_integrate_pgf_product_c(orthorhombic=orthorhombic, &
450 : compute_tau=my_compute_tau, &
451 : border_mask=border_mask, &
452 : la_max=la_max, &
453 : la_min=la_min, &
454 : lb_max=lb_max, &
455 : lb_min=lb_min, &
456 : zeta=zeta, &
457 : zetb=zetb, &
458 : dh=C_LOC(rsgrid%desc%dh(1, 1)), &
459 : dh_inv=C_LOC(rsgrid%desc%dh_inv(1, 1)), &
460 : ra=C_LOC(ra(1)), &
461 : rab=C_LOC(rab(1)), &
462 : npts_global=C_LOC(npts_global(1)), &
463 : npts_local=C_LOC(npts_local(1)), &
464 : shift_local=C_LOC(shift_local(1)), &
465 : border_width=C_LOC(border_width(1)), &
466 : radius=radius, &
467 : o1=o1, &
468 : o2=o2, &
469 : n1=SIZE(hab, 1), &
470 : n2=SIZE(hab, 2), &
471 : grid=C_LOC(grid(1, 1, 1)), &
472 : hab=C_LOC(hab(1, 1)), &
473 : pab=pab_cptr, &
474 : forces=forces_cptr, &
475 : virials=virials_cptr, &
476 : hdab=hdab_cptr, &
477 : hadb=hadb_cptr, &
478 882452 : a_hdab=a_hdab_cptr)
479 :
480 882452 : IF (PRESENT(force_a) .AND. C_ASSOCIATED(forces_cptr)) &
481 1407548 : force_a = force_a + forces(:, 1)
482 882452 : IF (PRESENT(force_b) .AND. C_ASSOCIATED(forces_cptr)) &
483 1407548 : force_b = force_b + forces(:, 2)
484 882452 : IF (PRESENT(my_virial_a) .AND. C_ASSOCIATED(virials_cptr)) &
485 1041014 : my_virial_a = my_virial_a + virials(:, :, 1)
486 882452 : IF (PRESENT(my_virial_b) .AND. C_ASSOCIATED(virials_cptr)) &
487 1041014 : my_virial_b = my_virial_b + virials(:, :, 2)
488 :
489 882452 : END SUBROUTINE integrate_pgf_product
490 :
491 : ! **************************************************************************************************
492 : !> \brief Helper routines for getting rsgrid properties and asserting underlying assumptions.
493 : !> \param rsgrid ...
494 : !> \param npts_global ...
495 : !> \param npts_local ...
496 : !> \param shift_local ...
497 : !> \param border_width ...
498 : !> \author Ole Schuett
499 : ! **************************************************************************************************
500 2080253 : SUBROUTINE get_rsgrid_properties(rsgrid, npts_global, npts_local, shift_local, border_width)
501 : TYPE(realspace_grid_type), INTENT(IN) :: rsgrid
502 : INTEGER, DIMENSION(:) :: npts_global, npts_local, shift_local, &
503 : border_width
504 :
505 : INTEGER :: i
506 :
507 : ! See rs_grid_create() in ./src/pw/realspace_grid_types.F.
508 4160506 : CPASSERT(LBOUND(rsgrid%r, 1) == rsgrid%lb_local(1))
509 4160506 : CPASSERT(UBOUND(rsgrid%r, 1) == rsgrid%ub_local(1))
510 4160506 : CPASSERT(LBOUND(rsgrid%r, 2) == rsgrid%lb_local(2))
511 4160506 : CPASSERT(UBOUND(rsgrid%r, 2) == rsgrid%ub_local(2))
512 4160506 : CPASSERT(LBOUND(rsgrid%r, 3) == rsgrid%lb_local(3))
513 4160506 : CPASSERT(UBOUND(rsgrid%r, 3) == rsgrid%ub_local(3))
514 :
515 : ! While the rsgrid code assumes that the grid starts at rsgrid%lb,
516 : ! the collocate code assumes that the grid starts at (1,1,1) in Fortran, or (0,0,0) in C.
517 : ! So, a point rp(:) gets the following grid coordinates MODULO(rp(:)/dr(:),npts_global(:))
518 :
519 : ! Number of global grid points in each direction.
520 8321012 : npts_global = rsgrid%desc%ub - rsgrid%desc%lb + 1
521 :
522 : ! Number of local grid points in each direction.
523 8321012 : npts_local = rsgrid%ub_local - rsgrid%lb_local + 1
524 :
525 : ! Number of points the local grid is shifted wrt global grid.
526 8321012 : shift_local = rsgrid%lb_local - rsgrid%desc%lb
527 :
528 : ! Convert rsgrid%desc%border and rsgrid%desc%perd into the more convenient border_width array.
529 8321012 : DO i = 1, 3
530 8321012 : IF (rsgrid%desc%perd(i) == 1) THEN
531 : ! Periodic meaning the grid in this direction is entriely present on every processor.
532 6240323 : CPASSERT(npts_local(i) == npts_global(i))
533 6240323 : CPASSERT(shift_local(i) == 0)
534 : ! No need for halo regions.
535 6240323 : border_width(i) = 0
536 : ELSE
537 : ! Not periodic meaning the grid in this direction is distributed among processors.
538 436 : CPASSERT(npts_local(i) <= npts_global(i))
539 : ! Check bounds of grid section that is owned by this processor.
540 436 : CPASSERT(rsgrid%lb_real(i) == rsgrid%lb_local(i) + rsgrid%desc%border)
541 436 : CPASSERT(rsgrid%ub_real(i) == rsgrid%ub_local(i) - rsgrid%desc%border)
542 : ! We have halo regions.
543 436 : border_width(i) = rsgrid%desc%border
544 : END IF
545 : END DO
546 2080253 : END SUBROUTINE get_rsgrid_properties
547 :
548 : ! **************************************************************************************************
549 : !> \brief Allocates a basis set which can be passed to grid_create_task_list.
550 : !> \param nset ...
551 : !> \param nsgf ...
552 : !> \param maxco ...
553 : !> \param maxpgf ...
554 : !> \param lmin ...
555 : !> \param lmax ...
556 : !> \param npgf ...
557 : !> \param nsgf_set ...
558 : !> \param first_sgf ...
559 : !> \param sphi ...
560 : !> \param zet ...
561 : !> \param basis_set ...
562 : !> \author Ole Schuett
563 : ! **************************************************************************************************
564 14335 : SUBROUTINE grid_create_basis_set(nset, nsgf, maxco, maxpgf, &
565 14335 : lmin, lmax, npgf, nsgf_set, first_sgf, sphi, zet, &
566 : basis_set)
567 : INTEGER, INTENT(IN) :: nset, nsgf, maxco, maxpgf
568 : INTEGER, DIMENSION(:), INTENT(IN), TARGET :: lmin, lmax, npgf, nsgf_set
569 : INTEGER, DIMENSION(:, :), INTENT(IN) :: first_sgf
570 : REAL(KIND=dp), DIMENSION(:, :), INTENT(IN), TARGET :: sphi, zet
571 : TYPE(grid_basis_set_type), INTENT(INOUT) :: basis_set
572 :
573 : CHARACTER(LEN=*), PARAMETER :: routineN = 'grid_create_basis_set'
574 :
575 : INTEGER :: handle
576 28670 : INTEGER, DIMENSION(nset), TARGET :: my_first_sgf
577 : TYPE(C_PTR) :: first_sgf_c, lmax_c, lmin_c, npgf_c, &
578 : nsgf_set_c, sphi_c, zet_c
579 : INTERFACE
580 : SUBROUTINE grid_create_basis_set_c(nset, nsgf, maxco, maxpgf, &
581 : lmin, lmax, npgf, nsgf_set, first_sgf, sphi, zet, &
582 : basis_set) &
583 : BIND(C, name="grid_create_basis_set")
584 : IMPORT :: C_PTR, C_INT
585 : INTEGER(KIND=C_INT), VALUE :: nset
586 : INTEGER(KIND=C_INT), VALUE :: nsgf
587 : INTEGER(KIND=C_INT), VALUE :: maxco
588 : INTEGER(KIND=C_INT), VALUE :: maxpgf
589 : TYPE(C_PTR), VALUE :: lmin
590 : TYPE(C_PTR), VALUE :: lmax
591 : TYPE(C_PTR), VALUE :: npgf
592 : TYPE(C_PTR), VALUE :: nsgf_set
593 : TYPE(C_PTR), VALUE :: first_sgf
594 : TYPE(C_PTR), VALUE :: sphi
595 : TYPE(C_PTR), VALUE :: zet
596 : TYPE(C_PTR) :: basis_set
597 : END SUBROUTINE grid_create_basis_set_c
598 : END INTERFACE
599 :
600 14335 : CALL timeset(routineN, handle)
601 :
602 14335 : CPASSERT(SIZE(lmin) == nset)
603 14335 : CPASSERT(SIZE(lmin) == nset)
604 14335 : CPASSERT(SIZE(lmax) == nset)
605 14335 : CPASSERT(SIZE(npgf) == nset)
606 14335 : CPASSERT(SIZE(nsgf_set) == nset)
607 14335 : CPASSERT(SIZE(first_sgf, 2) == nset)
608 14335 : CPASSERT(SIZE(sphi, 1) == maxco .AND. SIZE(sphi, 2) == nsgf)
609 14335 : CPASSERT(SIZE(zet, 1) == maxpgf .AND. SIZE(zet, 2) == nset)
610 14335 : CPASSERT(.NOT. C_ASSOCIATED(basis_set%c_ptr))
611 :
612 : #if __GNUC__ >= 9
613 14335 : CPASSERT(IS_CONTIGUOUS(lmin))
614 14335 : CPASSERT(IS_CONTIGUOUS(lmax))
615 14335 : CPASSERT(IS_CONTIGUOUS(npgf))
616 14335 : CPASSERT(IS_CONTIGUOUS(nsgf_set))
617 : CPASSERT(IS_CONTIGUOUS(my_first_sgf))
618 14335 : CPASSERT(IS_CONTIGUOUS(sphi))
619 14335 : CPASSERT(IS_CONTIGUOUS(zet))
620 : #endif
621 :
622 14335 : lmin_c = C_NULL_PTR
623 14335 : lmax_c = C_NULL_PTR
624 14335 : npgf_c = C_NULL_PTR
625 14335 : nsgf_set_c = C_NULL_PTR
626 14335 : first_sgf_c = C_NULL_PTR
627 14335 : sphi_c = C_NULL_PTR
628 14335 : zet_c = C_NULL_PTR
629 :
630 : ! Basis sets arrays can be empty, need to check before accessing the first element.
631 14335 : IF (nset > 0) THEN
632 14327 : lmin_c = C_LOC(lmin(1))
633 14327 : lmax_c = C_LOC(lmax(1))
634 14327 : npgf_c = C_LOC(npgf(1))
635 14327 : nsgf_set_c = C_LOC(nsgf_set(1))
636 : END IF
637 43005 : IF (SIZE(first_sgf) > 0) THEN
638 47142 : my_first_sgf(:) = first_sgf(1, :) ! make a contiguous copy
639 14327 : first_sgf_c = C_LOC(my_first_sgf(1))
640 : END IF
641 43005 : IF (SIZE(sphi) > 0) THEN
642 14325 : sphi_c = C_LOC(sphi(1, 1))
643 : END IF
644 43005 : IF (SIZE(zet) > 0) THEN
645 14325 : zet_c = C_LOC(zet(1, 1))
646 : END IF
647 :
648 : CALL grid_create_basis_set_c(nset=nset, &
649 : nsgf=nsgf, &
650 : maxco=maxco, &
651 : maxpgf=maxpgf, &
652 : lmin=lmin_c, &
653 : lmax=lmax_c, &
654 : npgf=npgf_c, &
655 : nsgf_set=nsgf_set_c, &
656 : first_sgf=first_sgf_c, &
657 : sphi=sphi_c, &
658 : zet=zet_c, &
659 14335 : basis_set=basis_set%c_ptr)
660 14335 : CPASSERT(C_ASSOCIATED(basis_set%c_ptr))
661 :
662 14335 : CALL timestop(handle)
663 14335 : END SUBROUTINE grid_create_basis_set
664 :
665 : ! **************************************************************************************************
666 : !> \brief Deallocates given basis set.
667 : !> \param basis_set ...
668 : !> \author Ole Schuett
669 : ! **************************************************************************************************
670 14335 : SUBROUTINE grid_free_basis_set(basis_set)
671 : TYPE(grid_basis_set_type), INTENT(INOUT) :: basis_set
672 :
673 : CHARACTER(LEN=*), PARAMETER :: routineN = 'grid_free_basis_set'
674 :
675 : INTEGER :: handle
676 : INTERFACE
677 : SUBROUTINE grid_free_basis_set_c(basis_set) &
678 : BIND(C, name="grid_free_basis_set")
679 : IMPORT :: C_PTR
680 : TYPE(C_PTR), VALUE :: basis_set
681 : END SUBROUTINE grid_free_basis_set_c
682 : END INTERFACE
683 :
684 14335 : CALL timeset(routineN, handle)
685 :
686 14335 : CPASSERT(C_ASSOCIATED(basis_set%c_ptr))
687 :
688 14335 : CALL grid_free_basis_set_c(basis_set%c_ptr)
689 :
690 14335 : basis_set%c_ptr = C_NULL_PTR
691 :
692 14335 : CALL timestop(handle)
693 14335 : END SUBROUTINE grid_free_basis_set
694 :
695 : ! **************************************************************************************************
696 : !> \brief Allocates a task list which can be passed to grid_collocate_task_list.
697 : !> \param ntasks ...
698 : !> \param natoms ...
699 : !> \param nkinds ...
700 : !> \param nblocks ...
701 : !> \param block_offsets ...
702 : !> \param atom_positions ...
703 : !> \param atom_kinds ...
704 : !> \param basis_sets ...
705 : !> \param level_list ...
706 : !> \param iatom_list ...
707 : !> \param jatom_list ...
708 : !> \param iset_list ...
709 : !> \param jset_list ...
710 : !> \param ipgf_list ...
711 : !> \param jpgf_list ...
712 : !> \param border_mask_list ...
713 : !> \param block_num_list ...
714 : !> \param radius_list ...
715 : !> \param rab_list ...
716 : !> \param rs_grids ...
717 : !> \param task_list ...
718 : !> \author Ole Schuett
719 : ! **************************************************************************************************
720 13410 : SUBROUTINE grid_create_task_list(ntasks, natoms, nkinds, nblocks, &
721 13410 : block_offsets, atom_positions, atom_kinds, basis_sets, &
722 13410 : level_list, iatom_list, jatom_list, &
723 13410 : iset_list, jset_list, ipgf_list, jpgf_list, &
724 13410 : border_mask_list, block_num_list, &
725 13410 : radius_list, rab_list, rs_grids, task_list)
726 :
727 : INTEGER, INTENT(IN) :: ntasks, natoms, nkinds, nblocks
728 : INTEGER, DIMENSION(:), INTENT(IN), TARGET :: block_offsets
729 : REAL(KIND=dp), DIMENSION(:, :), INTENT(IN), TARGET :: atom_positions
730 : INTEGER, DIMENSION(:), INTENT(IN), TARGET :: atom_kinds
731 : TYPE(grid_basis_set_type), DIMENSION(:), &
732 : INTENT(IN), TARGET :: basis_sets
733 : INTEGER, DIMENSION(:), INTENT(IN), TARGET :: level_list, iatom_list, jatom_list, &
734 : iset_list, jset_list, ipgf_list, &
735 : jpgf_list, border_mask_list, &
736 : block_num_list
737 : REAL(KIND=dp), DIMENSION(:), INTENT(IN), TARGET :: radius_list
738 : REAL(KIND=dp), DIMENSION(:, :), INTENT(IN), TARGET :: rab_list
739 : TYPE(realspace_grid_type), DIMENSION(:), &
740 : INTENT(IN) :: rs_grids
741 : TYPE(grid_task_list_type), INTENT(INOUT) :: task_list
742 :
743 : CHARACTER(LEN=*), PARAMETER :: routineN = 'grid_create_task_list'
744 :
745 : INTEGER :: handle, ikind, ilevel, nlevels
746 13410 : INTEGER, ALLOCATABLE, DIMENSION(:, :), TARGET :: border_width, npts_global, npts_local, &
747 13410 : shift_local
748 : LOGICAL(KIND=C_BOOL) :: orthorhombic
749 : REAL(KIND=dp), ALLOCATABLE, DIMENSION(:, :, :), &
750 13410 : TARGET :: dh, dh_inv
751 : TYPE(C_PTR) :: block_num_list_c, block_offsets_c, border_mask_list_c, iatom_list_c, &
752 : ipgf_list_c, iset_list_c, jatom_list_c, jpgf_list_c, jset_list_c, level_list_c, &
753 : rab_list_c, radius_list_c
754 13410 : TYPE(C_PTR), ALLOCATABLE, DIMENSION(:), TARGET :: basis_sets_c
755 : INTERFACE
756 : SUBROUTINE grid_create_task_list_c(orthorhombic, &
757 : ntasks, nlevels, natoms, nkinds, nblocks, &
758 : block_offsets, atom_positions, atom_kinds, basis_sets, &
759 : level_list, iatom_list, jatom_list, &
760 : iset_list, jset_list, ipgf_list, jpgf_list, &
761 : border_mask_list, block_num_list, &
762 : radius_list, rab_list, &
763 : npts_global, npts_local, shift_local, &
764 : border_width, dh, dh_inv, task_list) &
765 : BIND(C, name="grid_create_task_list")
766 : IMPORT :: C_PTR, C_INT, C_BOOL
767 : LOGICAL(KIND=C_BOOL), VALUE :: orthorhombic
768 : INTEGER(KIND=C_INT), VALUE :: ntasks
769 : INTEGER(KIND=C_INT), VALUE :: nlevels
770 : INTEGER(KIND=C_INT), VALUE :: natoms
771 : INTEGER(KIND=C_INT), VALUE :: nkinds
772 : INTEGER(KIND=C_INT), VALUE :: nblocks
773 : TYPE(C_PTR), VALUE :: block_offsets
774 : TYPE(C_PTR), VALUE :: atom_positions
775 : TYPE(C_PTR), VALUE :: atom_kinds
776 : TYPE(C_PTR), VALUE :: basis_sets
777 : TYPE(C_PTR), VALUE :: level_list
778 : TYPE(C_PTR), VALUE :: iatom_list
779 : TYPE(C_PTR), VALUE :: jatom_list
780 : TYPE(C_PTR), VALUE :: iset_list
781 : TYPE(C_PTR), VALUE :: jset_list
782 : TYPE(C_PTR), VALUE :: ipgf_list
783 : TYPE(C_PTR), VALUE :: jpgf_list
784 : TYPE(C_PTR), VALUE :: border_mask_list
785 : TYPE(C_PTR), VALUE :: block_num_list
786 : TYPE(C_PTR), VALUE :: radius_list
787 : TYPE(C_PTR), VALUE :: rab_list
788 : TYPE(C_PTR), VALUE :: npts_global
789 : TYPE(C_PTR), VALUE :: npts_local
790 : TYPE(C_PTR), VALUE :: shift_local
791 : TYPE(C_PTR), VALUE :: border_width
792 : TYPE(C_PTR), VALUE :: dh
793 : TYPE(C_PTR), VALUE :: dh_inv
794 : TYPE(C_PTR) :: task_list
795 : END SUBROUTINE grid_create_task_list_c
796 : END INTERFACE
797 :
798 13410 : CALL timeset(routineN, handle)
799 :
800 13410 : CPASSERT(SIZE(block_offsets) == nblocks)
801 13410 : CPASSERT(SIZE(atom_positions, 1) == 3 .AND. SIZE(atom_positions, 2) == natoms)
802 13410 : CPASSERT(SIZE(atom_kinds) == natoms)
803 13410 : CPASSERT(SIZE(basis_sets) == nkinds)
804 13410 : CPASSERT(SIZE(level_list) == ntasks)
805 13410 : CPASSERT(SIZE(iatom_list) == ntasks)
806 13410 : CPASSERT(SIZE(jatom_list) == ntasks)
807 13410 : CPASSERT(SIZE(iset_list) == ntasks)
808 13410 : CPASSERT(SIZE(jset_list) == ntasks)
809 13410 : CPASSERT(SIZE(ipgf_list) == ntasks)
810 13410 : CPASSERT(SIZE(jpgf_list) == ntasks)
811 13410 : CPASSERT(SIZE(border_mask_list) == ntasks)
812 13410 : CPASSERT(SIZE(block_num_list) == ntasks)
813 13410 : CPASSERT(SIZE(radius_list) == ntasks)
814 13410 : CPASSERT(SIZE(rab_list, 1) == 3 .AND. SIZE(rab_list, 2) == ntasks)
815 :
816 40230 : ALLOCATE (basis_sets_c(nkinds))
817 37467 : DO ikind = 1, nkinds
818 37467 : basis_sets_c(ikind) = basis_sets(ikind)%c_ptr
819 : END DO
820 :
821 13410 : nlevels = SIZE(rs_grids)
822 13410 : CPASSERT(nlevels > 0)
823 13410 : orthorhombic = LOGICAL(rs_grids(1)%desc%orthorhombic, C_BOOL)
824 :
825 53640 : ALLOCATE (npts_global(3, nlevels), npts_local(3, nlevels))
826 40230 : ALLOCATE (shift_local(3, nlevels), border_width(3, nlevels))
827 53640 : ALLOCATE (dh(3, 3, nlevels), dh_inv(3, 3, nlevels))
828 66598 : DO ilevel = 1, nlevels
829 13410 : ASSOCIATE (rsgrid => rs_grids(ilevel))
830 : CALL get_rsgrid_properties(rsgrid=rsgrid, &
831 : npts_global=npts_global(:, ilevel), &
832 : npts_local=npts_local(:, ilevel), &
833 : shift_local=shift_local(:, ilevel), &
834 53188 : border_width=border_width(:, ilevel))
835 53188 : CPASSERT(rsgrid%desc%orthorhombic .EQV. orthorhombic) ! should be the same for all levels
836 691444 : dh(:, :, ilevel) = rsgrid%desc%dh(:, :)
837 744632 : dh_inv(:, :, ilevel) = rsgrid%desc%dh_inv(:, :)
838 : END ASSOCIATE
839 : END DO
840 :
841 : #if __GNUC__ >= 9
842 13410 : CPASSERT(IS_CONTIGUOUS(block_offsets))
843 13410 : CPASSERT(IS_CONTIGUOUS(atom_positions))
844 13410 : CPASSERT(IS_CONTIGUOUS(atom_kinds))
845 13410 : CPASSERT(IS_CONTIGUOUS(basis_sets))
846 13410 : CPASSERT(IS_CONTIGUOUS(level_list))
847 13410 : CPASSERT(IS_CONTIGUOUS(iatom_list))
848 13410 : CPASSERT(IS_CONTIGUOUS(jatom_list))
849 13410 : CPASSERT(IS_CONTIGUOUS(iset_list))
850 13410 : CPASSERT(IS_CONTIGUOUS(jset_list))
851 13410 : CPASSERT(IS_CONTIGUOUS(ipgf_list))
852 13410 : CPASSERT(IS_CONTIGUOUS(jpgf_list))
853 13410 : CPASSERT(IS_CONTIGUOUS(border_mask_list))
854 13410 : CPASSERT(IS_CONTIGUOUS(block_num_list))
855 13410 : CPASSERT(IS_CONTIGUOUS(radius_list))
856 13410 : CPASSERT(IS_CONTIGUOUS(rab_list))
857 : CPASSERT(IS_CONTIGUOUS(npts_global))
858 : CPASSERT(IS_CONTIGUOUS(npts_local))
859 : CPASSERT(IS_CONTIGUOUS(shift_local))
860 : CPASSERT(IS_CONTIGUOUS(border_width))
861 : CPASSERT(IS_CONTIGUOUS(dh))
862 : CPASSERT(IS_CONTIGUOUS(dh_inv))
863 : #endif
864 :
865 13410 : IF (ntasks > 0) THEN
866 : block_offsets_c = C_LOC(block_offsets(1))
867 : level_list_c = C_LOC(level_list(1))
868 : iatom_list_c = C_LOC(iatom_list(1))
869 : jatom_list_c = C_LOC(jatom_list(1))
870 : iset_list_c = C_LOC(iset_list(1))
871 : jset_list_c = C_LOC(jset_list(1))
872 : ipgf_list_c = C_LOC(ipgf_list(1))
873 : jpgf_list_c = C_LOC(jpgf_list(1))
874 : border_mask_list_c = C_LOC(border_mask_list(1))
875 : block_num_list_c = C_LOC(block_num_list(1))
876 : radius_list_c = C_LOC(radius_list(1))
877 : rab_list_c = C_LOC(rab_list(1, 1))
878 : ELSE
879 : ! Without tasks the lists are empty and there is no first element to call C_LOC on.
880 227 : block_offsets_c = C_NULL_PTR
881 227 : level_list_c = C_NULL_PTR
882 227 : iatom_list_c = C_NULL_PTR
883 227 : jatom_list_c = C_NULL_PTR
884 227 : iset_list_c = C_NULL_PTR
885 227 : jset_list_c = C_NULL_PTR
886 227 : ipgf_list_c = C_NULL_PTR
887 227 : jpgf_list_c = C_NULL_PTR
888 227 : border_mask_list_c = C_NULL_PTR
889 227 : block_num_list_c = C_NULL_PTR
890 227 : radius_list_c = C_NULL_PTR
891 227 : rab_list_c = C_NULL_PTR
892 : END IF
893 :
894 : !If task_list%c_ptr is already allocated, then its memory will be reused or freed.
895 : CALL grid_create_task_list_c(orthorhombic=orthorhombic, &
896 : ntasks=ntasks, &
897 : nlevels=nlevels, &
898 : natoms=natoms, &
899 : nkinds=nkinds, &
900 : nblocks=nblocks, &
901 : block_offsets=block_offsets_c, &
902 : atom_positions=C_LOC(atom_positions(1, 1)), &
903 : atom_kinds=C_LOC(atom_kinds(1)), &
904 : basis_sets=C_LOC(basis_sets_c(1)), &
905 : level_list=level_list_c, &
906 : iatom_list=iatom_list_c, &
907 : jatom_list=jatom_list_c, &
908 : iset_list=iset_list_c, &
909 : jset_list=jset_list_c, &
910 : ipgf_list=ipgf_list_c, &
911 : jpgf_list=jpgf_list_c, &
912 : border_mask_list=border_mask_list_c, &
913 : block_num_list=block_num_list_c, &
914 : radius_list=radius_list_c, &
915 : rab_list=rab_list_c, &
916 : npts_global=C_LOC(npts_global(1, 1)), &
917 : npts_local=C_LOC(npts_local(1, 1)), &
918 : shift_local=C_LOC(shift_local(1, 1)), &
919 : border_width=C_LOC(border_width(1, 1)), &
920 : dh=C_LOC(dh(1, 1, 1)), &
921 : dh_inv=C_LOC(dh_inv(1, 1, 1)), &
922 13410 : task_list=task_list%c_ptr)
923 :
924 13410 : CPASSERT(C_ASSOCIATED(task_list%c_ptr))
925 :
926 13410 : CALL timestop(handle)
927 26820 : END SUBROUTINE grid_create_task_list
928 :
929 : ! **************************************************************************************************
930 : !> \brief Deallocates given task list, basis_sets have to be freed separately.
931 : !> \param task_list ...
932 : !> \author Ole Schuett
933 : ! **************************************************************************************************
934 8016 : SUBROUTINE grid_free_task_list(task_list)
935 : TYPE(grid_task_list_type), INTENT(INOUT) :: task_list
936 :
937 : CHARACTER(LEN=*), PARAMETER :: routineN = 'grid_free_task_list'
938 :
939 : INTEGER :: handle
940 : INTERFACE
941 : SUBROUTINE grid_free_task_list_c(task_list) &
942 : BIND(C, name="grid_free_task_list")
943 : IMPORT :: C_PTR
944 : TYPE(C_PTR), VALUE :: task_list
945 : END SUBROUTINE grid_free_task_list_c
946 : END INTERFACE
947 :
948 8016 : CALL timeset(routineN, handle)
949 :
950 8016 : IF (C_ASSOCIATED(task_list%c_ptr)) THEN
951 8016 : CALL grid_free_task_list_c(task_list%c_ptr)
952 : END IF
953 :
954 8016 : task_list%c_ptr = C_NULL_PTR
955 :
956 8016 : CALL timestop(handle)
957 8016 : END SUBROUTINE grid_free_task_list
958 :
959 : ! **************************************************************************************************
960 : !> \brief Collocate all tasks of in given list onto given grids.
961 : !> \param task_list ...
962 : !> \param ga_gb_function ...
963 : !> \param pab_blocks ...
964 : !> \param rs_grids ...
965 : !> \author Ole Schuett
966 : ! **************************************************************************************************
967 193857 : SUBROUTINE grid_collocate_task_list(task_list, ga_gb_function, pab_blocks, rs_grids)
968 : TYPE(grid_task_list_type), INTENT(IN) :: task_list
969 : INTEGER, INTENT(IN) :: ga_gb_function
970 : TYPE(offload_buffer_type), INTENT(IN) :: pab_blocks
971 : TYPE(realspace_grid_type), DIMENSION(:), &
972 : INTENT(IN) :: rs_grids
973 :
974 : CHARACTER(LEN=*), PARAMETER :: routineN = 'grid_collocate_task_list'
975 :
976 : INTEGER :: handle, ilevel, nlevels
977 193857 : INTEGER, ALLOCATABLE, DIMENSION(:, :), TARGET :: npts_local
978 193857 : TYPE(C_PTR), ALLOCATABLE, DIMENSION(:), TARGET :: grids_c
979 : INTERFACE
980 : SUBROUTINE grid_collocate_task_list_c(task_list, func, nlevels, &
981 : npts_local, pab_blocks, grids) &
982 : BIND(C, name="grid_collocate_task_list")
983 : IMPORT :: C_PTR, C_INT, C_BOOL
984 : TYPE(C_PTR), VALUE :: task_list
985 : INTEGER(KIND=C_INT), VALUE :: func
986 : INTEGER(KIND=C_INT), VALUE :: nlevels
987 : TYPE(C_PTR), VALUE :: npts_local
988 : TYPE(C_PTR), VALUE :: pab_blocks
989 : TYPE(C_PTR), VALUE :: grids
990 : END SUBROUTINE grid_collocate_task_list_c
991 : END INTERFACE
992 :
993 193857 : CALL timeset(routineN, handle)
994 :
995 193857 : nlevels = SIZE(rs_grids)
996 193857 : CPASSERT(nlevels > 0)
997 :
998 581571 : ALLOCATE (grids_c(nlevels))
999 581571 : ALLOCATE (npts_local(3, nlevels))
1000 961279 : DO ilevel = 1, nlevels
1001 193857 : ASSOCIATE (rsgrid => rs_grids(ilevel))
1002 3069688 : npts_local(:, ilevel) = rsgrid%ub_local - rsgrid%lb_local + 1
1003 1534844 : grids_c(ilevel) = rsgrid%buffer%c_ptr
1004 : END ASSOCIATE
1005 : END DO
1006 :
1007 : #if __GNUC__ >= 9
1008 : CPASSERT(IS_CONTIGUOUS(npts_local))
1009 : CPASSERT(IS_CONTIGUOUS(grids_c))
1010 : #endif
1011 :
1012 193857 : CPASSERT(C_ASSOCIATED(task_list%c_ptr))
1013 193857 : CPASSERT(C_ASSOCIATED(pab_blocks%c_ptr))
1014 :
1015 : CALL grid_collocate_task_list_c(task_list=task_list%c_ptr, &
1016 : func=ga_gb_function, &
1017 : nlevels=nlevels, &
1018 : npts_local=C_LOC(npts_local(1, 1)), &
1019 : pab_blocks=pab_blocks%c_ptr, &
1020 193857 : grids=C_LOC(grids_c(1)))
1021 :
1022 193857 : CALL timestop(handle)
1023 387714 : END SUBROUTINE grid_collocate_task_list
1024 :
1025 : ! **************************************************************************************************
1026 : !> \brief Integrate all tasks of in given list from given grids.
1027 : !> \param task_list ...
1028 : !> \param compute_tau ...
1029 : !> \param calculate_forces ...
1030 : !> \param calculate_virial ...
1031 : !> \param pab_blocks ...
1032 : !> \param rs_grids ...
1033 : !> \param hab_blocks ...
1034 : !> \param forces ...
1035 : !> \param virial ...
1036 : !> \author Ole Schuett
1037 : ! **************************************************************************************************
1038 181955 : SUBROUTINE grid_integrate_task_list(task_list, compute_tau, calculate_forces, calculate_virial, &
1039 181955 : pab_blocks, rs_grids, hab_blocks, forces, virial)
1040 : TYPE(grid_task_list_type), INTENT(IN) :: task_list
1041 : LOGICAL, INTENT(IN) :: compute_tau, calculate_forces, &
1042 : calculate_virial
1043 : TYPE(offload_buffer_type), INTENT(IN) :: pab_blocks
1044 : TYPE(realspace_grid_type), DIMENSION(:), &
1045 : INTENT(IN) :: rs_grids
1046 : TYPE(offload_buffer_type), INTENT(INOUT) :: hab_blocks
1047 : REAL(KIND=dp), DIMENSION(:, :), INTENT(INOUT), &
1048 : TARGET :: forces
1049 : REAL(KIND=dp), DIMENSION(3, 3), INTENT(INOUT), &
1050 : TARGET :: virial
1051 :
1052 : CHARACTER(LEN=*), PARAMETER :: routineN = 'grid_integrate_task_list'
1053 :
1054 : INTEGER :: handle, ilevel, nlevels
1055 181955 : INTEGER, ALLOCATABLE, DIMENSION(:, :), TARGET :: npts_local
1056 : TYPE(C_PTR) :: forces_c, virial_c
1057 181955 : TYPE(C_PTR), ALLOCATABLE, DIMENSION(:), TARGET :: grids_c
1058 : INTERFACE
1059 : SUBROUTINE grid_integrate_task_list_c(task_list, compute_tau, natoms, &
1060 : nlevels, npts_local, &
1061 : pab_blocks, grids, hab_blocks, forces, virial) &
1062 : BIND(C, name="grid_integrate_task_list")
1063 : IMPORT :: C_PTR, C_INT, C_BOOL
1064 : TYPE(C_PTR), VALUE :: task_list
1065 : LOGICAL(KIND=C_BOOL), VALUE :: compute_tau
1066 : INTEGER(KIND=C_INT), VALUE :: natoms
1067 : INTEGER(KIND=C_INT), VALUE :: nlevels
1068 : TYPE(C_PTR), VALUE :: npts_local
1069 : TYPE(C_PTR), VALUE :: pab_blocks
1070 : TYPE(C_PTR), VALUE :: grids
1071 : TYPE(C_PTR), VALUE :: hab_blocks
1072 : TYPE(C_PTR), VALUE :: forces
1073 : TYPE(C_PTR), VALUE :: virial
1074 : END SUBROUTINE grid_integrate_task_list_c
1075 : END INTERFACE
1076 :
1077 181955 : CALL timeset(routineN, handle)
1078 :
1079 181955 : nlevels = SIZE(rs_grids)
1080 181955 : CPASSERT(nlevels > 0)
1081 :
1082 545865 : ALLOCATE (grids_c(nlevels))
1083 545865 : ALLOCATE (npts_local(3, nlevels))
1084 901859 : DO ilevel = 1, nlevels
1085 181955 : ASSOCIATE (rsgrid => rs_grids(ilevel))
1086 2879616 : npts_local(:, ilevel) = rsgrid%ub_local - rsgrid%lb_local + 1
1087 1439808 : grids_c(ilevel) = rsgrid%buffer%c_ptr
1088 : END ASSOCIATE
1089 : END DO
1090 :
1091 181955 : IF (calculate_forces) THEN
1092 : forces_c = C_LOC(forces(1, 1))
1093 : ELSE
1094 159192 : forces_c = C_NULL_PTR
1095 : END IF
1096 :
1097 181955 : IF (calculate_virial) THEN
1098 3585 : virial_c = C_LOC(virial(1, 1))
1099 : ELSE
1100 : virial_c = C_NULL_PTR
1101 : END IF
1102 :
1103 : #if __GNUC__ >= 9
1104 : CPASSERT(IS_CONTIGUOUS(npts_local))
1105 : CPASSERT(IS_CONTIGUOUS(grids_c))
1106 181955 : CPASSERT(IS_CONTIGUOUS(forces))
1107 : CPASSERT(IS_CONTIGUOUS(virial))
1108 : #endif
1109 :
1110 181955 : CPASSERT(SIZE(forces, 1) == 3)
1111 181955 : CPASSERT(C_ASSOCIATED(task_list%c_ptr))
1112 181955 : CPASSERT(C_ASSOCIATED(hab_blocks%c_ptr))
1113 181955 : CPASSERT(C_ASSOCIATED(pab_blocks%c_ptr) .OR. .NOT. calculate_forces)
1114 181955 : CPASSERT(C_ASSOCIATED(pab_blocks%c_ptr) .OR. .NOT. calculate_virial)
1115 :
1116 : CALL grid_integrate_task_list_c(task_list=task_list%c_ptr, &
1117 : compute_tau=LOGICAL(compute_tau, C_BOOL), &
1118 : natoms=SIZE(forces, 2), &
1119 : nlevels=nlevels, &
1120 : npts_local=C_LOC(npts_local(1, 1)), &
1121 : pab_blocks=pab_blocks%c_ptr, &
1122 : grids=C_LOC(grids_c(1)), &
1123 : hab_blocks=hab_blocks%c_ptr, &
1124 : forces=forces_c, &
1125 181955 : virial=virial_c)
1126 :
1127 181955 : CALL timestop(handle)
1128 363910 : END SUBROUTINE grid_integrate_task_list
1129 :
1130 : ! **************************************************************************************************
1131 : !> \brief Initialize grid library
1132 : !> \author Ole Schuett
1133 : ! **************************************************************************************************
1134 8530 : SUBROUTINE grid_library_init()
1135 : INTERFACE
1136 : SUBROUTINE grid_library_init_c() BIND(C, name="grid_library_init")
1137 : END SUBROUTINE grid_library_init_c
1138 : END INTERFACE
1139 :
1140 8530 : CALL grid_library_init_c()
1141 :
1142 8530 : END SUBROUTINE grid_library_init
1143 :
1144 : ! **************************************************************************************************
1145 : !> \brief Finalize grid library
1146 : !> \author Ole Schuett
1147 : ! **************************************************************************************************
1148 8530 : SUBROUTINE grid_library_finalize()
1149 : INTERFACE
1150 : SUBROUTINE grid_library_finalize_c() BIND(C, name="grid_library_finalize")
1151 : END SUBROUTINE grid_library_finalize_c
1152 : END INTERFACE
1153 :
1154 8530 : CALL grid_library_finalize_c()
1155 :
1156 8530 : END SUBROUTINE grid_library_finalize
1157 :
1158 : ! **************************************************************************************************
1159 : !> \brief Configures the grid library
1160 : !> \param backend : backend to be used for collocate/integrate, possible values are REF, CPU, GPU
1161 : !> \param validate : if set to true, compare the results of all backend to the reference backend
1162 : !> \param apply_cutoff : apply a spherical cutoff before collocating or integrating. Only relevant for CPU backend
1163 : !> \author Ole Schuett
1164 : ! **************************************************************************************************
1165 8648 : SUBROUTINE grid_library_set_config(backend, validate, apply_cutoff)
1166 : INTEGER, INTENT(IN) :: backend
1167 : LOGICAL, INTENT(IN) :: validate, apply_cutoff
1168 :
1169 : INTERFACE
1170 : SUBROUTINE grid_library_set_config_c(backend, validate, apply_cutoff) &
1171 : BIND(C, name="grid_library_set_config")
1172 : IMPORT :: C_INT, C_BOOL
1173 : INTEGER(KIND=C_INT), VALUE :: backend
1174 : LOGICAL(KIND=C_BOOL), VALUE :: validate
1175 : LOGICAL(KIND=C_BOOL), VALUE :: apply_cutoff
1176 : END SUBROUTINE grid_library_set_config_c
1177 : END INTERFACE
1178 :
1179 : CALL grid_library_set_config_c(backend=backend, &
1180 : validate=LOGICAL(validate, C_BOOL), &
1181 8648 : apply_cutoff=LOGICAL(apply_cutoff, C_BOOL))
1182 :
1183 8648 : END SUBROUTINE grid_library_set_config
1184 :
1185 : ! **************************************************************************************************
1186 : !> \brief Print grid library statistics
1187 : !> \param mpi_comm ...
1188 : !> \param output_unit ...
1189 : !> \author Ole Schuett
1190 : ! **************************************************************************************************
1191 8648 : SUBROUTINE grid_library_print_stats(mpi_comm, output_unit)
1192 : TYPE(mp_comm_type) :: mpi_comm
1193 : INTEGER, INTENT(IN) :: output_unit
1194 :
1195 : INTERFACE
1196 : SUBROUTINE grid_library_print_stats_c(mpi_sum_func, mpi_comm, print_func, output_unit) &
1197 : BIND(C, name="grid_library_print_stats")
1198 : IMPORT :: C_FUNPTR, C_INT
1199 : TYPE(C_FUNPTR), VALUE :: mpi_sum_func
1200 : INTEGER(KIND=C_INT), VALUE :: mpi_comm
1201 : TYPE(C_FUNPTR), VALUE :: print_func
1202 : INTEGER(KIND=C_INT), VALUE :: output_unit
1203 : END SUBROUTINE grid_library_print_stats_c
1204 : END INTERFACE
1205 :
1206 : ! Since Fortran units and mpi groups can't be used from C, we pass function pointers instead.
1207 : CALL grid_library_print_stats_c(mpi_sum_func=C_FUNLOC(mpi_sum_func), &
1208 : mpi_comm=mpi_comm%get_handle(), &
1209 : print_func=C_FUNLOC(print_func), &
1210 8648 : output_unit=output_unit)
1211 :
1212 8648 : END SUBROUTINE grid_library_print_stats
1213 :
1214 : ! **************************************************************************************************
1215 : !> \brief Callback to run mpi_sum on a Fortran MPI communicator.
1216 : !> \param number ...
1217 : !> \param mpi_comm ...
1218 : !> \author Ole Schuett
1219 : ! **************************************************************************************************
1220 3459200 : SUBROUTINE mpi_sum_func(number, mpi_comm) BIND(C, name="grid_api_mpi_sum_func")
1221 : INTEGER(KIND=C_LONG), INTENT(INOUT) :: number
1222 : INTEGER(KIND=C_INT), INTENT(IN), VALUE :: mpi_comm
1223 :
1224 : TYPE(mp_comm_type) :: my_mpi_comm
1225 :
1226 : ! Convert the handle to the default integer kind and convert it to the communicator type
1227 3459200 : CALL my_mpi_comm%set_handle(INT(mpi_comm))
1228 :
1229 3459200 : CALL my_mpi_comm%sum(number)
1230 3459200 : END SUBROUTINE mpi_sum_func
1231 :
1232 : ! **************************************************************************************************
1233 : !> \brief Callback to write to a Fortran output unit.
1234 : !> \param message ...
1235 : !> \param output_unit ...
1236 : !> \author Ole Schuett
1237 : ! **************************************************************************************************
1238 119166 : SUBROUTINE print_func(message, output_unit) BIND(C, name="grid_api_print_func")
1239 : CHARACTER(LEN=1, KIND=C_CHAR), INTENT(IN) :: message(*)
1240 : INTEGER(KIND=C_INT), INTENT(IN), VALUE :: output_unit
1241 :
1242 : CHARACTER(LEN=1000) :: buffer
1243 : INTEGER :: nchars
1244 :
1245 119166 : IF (output_unit <= 0) &
1246 59359 : RETURN
1247 :
1248 : ! Convert C char array into Fortran string.
1249 59807 : nchars = strlcpy_c2f(buffer, message)
1250 :
1251 : ! Print the message.
1252 59807 : WRITE (output_unit, FMT="(A)", ADVANCE="NO") buffer(1:nchars)
1253 : END SUBROUTINE print_func
1254 :
1255 0 : END MODULE grid_api
|
