(git:e7e05ae)
grid_api.F
Go to the documentation of this file.
1 !--------------------------------------------------------------------------------------------------!
2 ! CP2K: A general program to perform molecular dynamics simulations !
3 ! Copyright 2000-2024 CP2K developers group <https://cp2k.org> !
4 ! !
5 ! SPDX-License-Identifier: 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  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  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  border_mask = 0
179  IF (PRESENT(use_subpatch)) THEN
180  IF (use_subpatch) THEN
181  cpassert(PRESENT(subpatch_pattern))
182  border_mask = iand(63, not(subpatch_pattern)) ! invert last 6 bits
183  END IF
184  END IF
185 
186  orthorhombic = LOGICAL(rsgrid%desc%orthorhombic, c_bool)
187 
188  cpassert(lbound(pab, 1) == 1)
189  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  border_width=border_width)
195 
196  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  cpassert(is_contiguous(pab))
208  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  grid=c_loc(grid(1, 1, 1)))
238 
239  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  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  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  IF (radius == 0.0_dp) THEN
357  RETURN
358  END IF
359 
360  border_mask = 0
361  IF (PRESENT(use_subpatch)) THEN
362  IF (use_subpatch) THEN
363  cpassert(PRESENT(subpatch_pattern))
364  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  IF (PRESENT(compute_tau)) THEN
370  my_compute_tau = LOGICAL(compute_tau, c_bool)
371  ELSE
372  my_compute_tau = .false.
373  END IF
374 
375  IF (PRESENT(use_virial)) THEN
376  my_use_virial = use_virial
377  ELSE
378  my_use_virial = .false.
379  END IF
380 
381  IF (calculate_forces) THEN
382  cpassert(PRESENT(pab))
383  pab_cptr = c_loc(pab(1, 1))
384  forces(:, :) = 0.0_dp
385  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  IF (calculate_forces .AND. my_use_virial) THEN
392  virials(:, :, :) = 0.0_dp
393  virials_cptr = c_loc(virials(1, 1, 1))
394  ELSE
395  virials_cptr = c_null_ptr
396  END IF
397 
398  IF (calculate_forces .AND. PRESENT(hdab)) THEN
399  hdab_cptr = c_loc(hdab(1, 1, 1))
400  ELSE
401  hdab_cptr = c_null_ptr
402  END IF
403 
404  IF (calculate_forces .AND. PRESENT(hadb)) THEN
405  hadb_cptr = c_loc(hadb(1, 1, 1))
406  ELSE
407  hadb_cptr = c_null_ptr
408  END IF
409 
410  IF (calculate_forces .AND. my_use_virial .AND. PRESENT(a_hdab)) THEN
411  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  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  border_width=border_width)
422 
423  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  cpassert(is_contiguous(grid))
435  cpassert(is_contiguous(hab))
436  cpassert(is_contiguous(forces))
437  cpassert(is_contiguous(virials))
438  IF (PRESENT(pab)) THEN
439  cpassert(is_contiguous(pab))
440  END IF
441  IF (PRESENT(hdab)) THEN
442  cpassert(is_contiguous(hdab))
443  END IF
444  IF (PRESENT(a_hdab)) THEN
445  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  a_hdab=a_hdab_cptr)
479 
480  IF (PRESENT(force_a) .AND. c_associated(forces_cptr)) &
481  force_a = force_a + forces(:, 1)
482  IF (PRESENT(force_b) .AND. c_associated(forces_cptr)) &
483  force_b = force_b + forces(:, 2)
484  IF (PRESENT(my_virial_a) .AND. c_associated(virials_cptr)) &
485  my_virial_a = my_virial_a + virials(:, :, 1)
486  IF (PRESENT(my_virial_b) .AND. c_associated(virials_cptr)) &
487  my_virial_b = my_virial_b + virials(:, :, 2)
488 
489  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  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  cpassert(lbound(rsgrid%r, 1) == rsgrid%lb_local(1))
509  cpassert(ubound(rsgrid%r, 1) == rsgrid%ub_local(1))
510  cpassert(lbound(rsgrid%r, 2) == rsgrid%lb_local(2))
511  cpassert(ubound(rsgrid%r, 2) == rsgrid%ub_local(2))
512  cpassert(lbound(rsgrid%r, 3) == rsgrid%lb_local(3))
513  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  npts_global = rsgrid%desc%ub - rsgrid%desc%lb + 1
521 
522  ! Number of local grid points in each direction.
523  npts_local = rsgrid%ub_local - rsgrid%lb_local + 1
524 
525  ! Number of points the local grid is shifted wrt global grid.
526  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  DO i = 1, 3
530  IF (rsgrid%desc%perd(i) == 1) THEN
531  ! Periodic meaning the grid in this direction is entriely present on every processor.
532  cpassert(npts_local(i) == npts_global(i))
533  cpassert(shift_local(i) == 0)
534  ! No need for halo regions.
535  border_width(i) = 0
536  ELSE
537  ! Not periodic meaning the grid in this direction is distributed among processors.
538  cpassert(npts_local(i) <= npts_global(i))
539  ! Check bounds of grid section that is owned by this processor.
540  cpassert(rsgrid%lb_real(i) == rsgrid%lb_local(i) + rsgrid%desc%border)
541  cpassert(rsgrid%ub_real(i) == rsgrid%ub_local(i) - rsgrid%desc%border)
542  ! We have halo regions.
543  border_width(i) = rsgrid%desc%border
544  END IF
545  END DO
546  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  SUBROUTINE grid_create_basis_set(nset, nsgf, maxco, maxpgf, &
565  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  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  CALL timeset(routinen, handle)
601 
602  cpassert(SIZE(lmin) == nset)
603  cpassert(SIZE(lmin) == nset)
604  cpassert(SIZE(lmax) == nset)
605  cpassert(SIZE(npgf) == nset)
606  cpassert(SIZE(nsgf_set) == nset)
607  cpassert(SIZE(first_sgf, 2) == nset)
608  cpassert(SIZE(sphi, 1) == maxco .AND. SIZE(sphi, 2) == nsgf)
609  cpassert(SIZE(zet, 1) == maxpgf .AND. SIZE(zet, 2) == nset)
610  cpassert(.NOT. c_associated(basis_set%c_ptr))
611 
612 #if __GNUC__ >= 9
613  cpassert(is_contiguous(lmin))
614  cpassert(is_contiguous(lmax))
615  cpassert(is_contiguous(npgf))
616  cpassert(is_contiguous(nsgf_set))
617  cpassert(is_contiguous(my_first_sgf))
618  cpassert(is_contiguous(sphi))
619  cpassert(is_contiguous(zet))
620 #endif
621 
622  lmin_c = c_null_ptr
623  lmax_c = c_null_ptr
624  npgf_c = c_null_ptr
625  nsgf_set_c = c_null_ptr
626  first_sgf_c = c_null_ptr
627  sphi_c = c_null_ptr
628  zet_c = c_null_ptr
629 
630  ! Basis sets arrays can be empty, need to check before accessing the first element.
631  IF (nset > 0) THEN
632  lmin_c = c_loc(lmin(1))
633  lmax_c = c_loc(lmax(1))
634  npgf_c = c_loc(npgf(1))
635  nsgf_set_c = c_loc(nsgf_set(1))
636  END IF
637  IF (SIZE(first_sgf) > 0) THEN
638  my_first_sgf(:) = first_sgf(1, :) ! make a contiguous copy
639  first_sgf_c = c_loc(my_first_sgf(1))
640  END IF
641  IF (SIZE(sphi) > 0) THEN
642  sphi_c = c_loc(sphi(1, 1))
643  END IF
644  IF (SIZE(zet) > 0) THEN
645  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  basis_set=basis_set%c_ptr)
660  cpassert(c_associated(basis_set%c_ptr))
661 
662  CALL timestop(handle)
663  END SUBROUTINE grid_create_basis_set
664 
665 ! **************************************************************************************************
666 !> \brief Deallocates given basis set.
667 !> \param basis_set ...
668 !> \author Ole Schuett
669 ! **************************************************************************************************
670  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  CALL timeset(routinen, handle)
685 
686  cpassert(c_associated(basis_set%c_ptr))
687 
688  CALL grid_free_basis_set_c(basis_set%c_ptr)
689 
690  basis_set%c_ptr = c_null_ptr
691 
692  CALL timestop(handle)
693  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  SUBROUTINE grid_create_task_list(ntasks, natoms, nkinds, nblocks, &
721  block_offsets, atom_positions, atom_kinds, basis_sets, &
722  level_list, iatom_list, jatom_list, &
723  iset_list, jset_list, ipgf_list, jpgf_list, &
724  border_mask_list, block_num_list, &
725  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  INTEGER, ALLOCATABLE, DIMENSION(:, :), TARGET :: border_width, npts_global, npts_local, &
747  shift_local
748  LOGICAL(KIND=C_BOOL) :: orthorhombic
749  REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :, :), &
750  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  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  CALL timeset(routinen, handle)
799 
800  cpassert(SIZE(block_offsets) == nblocks)
801  cpassert(SIZE(atom_positions, 1) == 3 .AND. SIZE(atom_positions, 2) == natoms)
802  cpassert(SIZE(atom_kinds) == natoms)
803  cpassert(SIZE(basis_sets) == nkinds)
804  cpassert(SIZE(level_list) == ntasks)
805  cpassert(SIZE(iatom_list) == ntasks)
806  cpassert(SIZE(jatom_list) == ntasks)
807  cpassert(SIZE(iset_list) == ntasks)
808  cpassert(SIZE(jset_list) == ntasks)
809  cpassert(SIZE(ipgf_list) == ntasks)
810  cpassert(SIZE(jpgf_list) == ntasks)
811  cpassert(SIZE(border_mask_list) == ntasks)
812  cpassert(SIZE(block_num_list) == ntasks)
813  cpassert(SIZE(radius_list) == ntasks)
814  cpassert(SIZE(rab_list, 1) == 3 .AND. SIZE(rab_list, 2) == ntasks)
815 
816  ALLOCATE (basis_sets_c(nkinds))
817  DO ikind = 1, nkinds
818  basis_sets_c(ikind) = basis_sets(ikind)%c_ptr
819  END DO
820 
821  nlevels = SIZE(rs_grids)
822  cpassert(nlevels > 0)
823  orthorhombic = LOGICAL(rs_grids(1)%desc%orthorhombic, c_bool)
824 
825  ALLOCATE (npts_global(3, nlevels), npts_local(3, nlevels))
826  ALLOCATE (shift_local(3, nlevels), border_width(3, nlevels))
827  ALLOCATE (dh(3, 3, nlevels), dh_inv(3, 3, nlevels))
828  DO ilevel = 1, nlevels
829  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  border_width=border_width(:, ilevel))
835  cpassert(rsgrid%desc%orthorhombic .EQV. orthorhombic) ! should be the same for all levels
836  dh(:, :, ilevel) = rsgrid%desc%dh(:, :)
837  dh_inv(:, :, ilevel) = rsgrid%desc%dh_inv(:, :)
838  END associate
839  END DO
840 
841 #if __GNUC__ >= 9
842  cpassert(is_contiguous(block_offsets))
843  cpassert(is_contiguous(atom_positions))
844  cpassert(is_contiguous(atom_kinds))
845  cpassert(is_contiguous(basis_sets))
846  cpassert(is_contiguous(level_list))
847  cpassert(is_contiguous(iatom_list))
848  cpassert(is_contiguous(jatom_list))
849  cpassert(is_contiguous(iset_list))
850  cpassert(is_contiguous(jset_list))
851  cpassert(is_contiguous(ipgf_list))
852  cpassert(is_contiguous(jpgf_list))
853  cpassert(is_contiguous(border_mask_list))
854  cpassert(is_contiguous(block_num_list))
855  cpassert(is_contiguous(radius_list))
856  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  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  block_offsets_c = c_null_ptr
881  level_list_c = c_null_ptr
882  iatom_list_c = c_null_ptr
883  jatom_list_c = c_null_ptr
884  iset_list_c = c_null_ptr
885  jset_list_c = c_null_ptr
886  ipgf_list_c = c_null_ptr
887  jpgf_list_c = c_null_ptr
888  border_mask_list_c = c_null_ptr
889  block_num_list_c = c_null_ptr
890  radius_list_c = c_null_ptr
891  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  task_list=task_list%c_ptr)
923 
924  cpassert(c_associated(task_list%c_ptr))
925 
926  CALL timestop(handle)
927  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  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  CALL timeset(routinen, handle)
949 
950  IF (c_associated(task_list%c_ptr)) THEN
951  CALL grid_free_task_list_c(task_list%c_ptr)
952  END IF
953 
954  task_list%c_ptr = c_null_ptr
955 
956  CALL timestop(handle)
957  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  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  INTEGER, ALLOCATABLE, DIMENSION(:, :), TARGET :: npts_local
978  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  CALL timeset(routinen, handle)
994 
995  nlevels = SIZE(rs_grids)
996  cpassert(nlevels > 0)
997 
998  ALLOCATE (grids_c(nlevels))
999  ALLOCATE (npts_local(3, nlevels))
1000  DO ilevel = 1, nlevels
1001  associate(rsgrid => rs_grids(ilevel))
1002  npts_local(:, ilevel) = rsgrid%ub_local - rsgrid%lb_local + 1
1003  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  cpassert(c_associated(task_list%c_ptr))
1013  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  grids=c_loc(grids_c(1)))
1021 
1022  CALL timestop(handle)
1023  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  SUBROUTINE grid_integrate_task_list(task_list, compute_tau, calculate_forces, calculate_virial, &
1039  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  INTEGER, ALLOCATABLE, DIMENSION(:, :), TARGET :: npts_local
1056  TYPE(c_ptr) :: forces_c, virial_c
1057  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  CALL timeset(routinen, handle)
1078 
1079  nlevels = SIZE(rs_grids)
1080  cpassert(nlevels > 0)
1081 
1082  ALLOCATE (grids_c(nlevels))
1083  ALLOCATE (npts_local(3, nlevels))
1084  DO ilevel = 1, nlevels
1085  associate(rsgrid => rs_grids(ilevel))
1086  npts_local(:, ilevel) = rsgrid%ub_local - rsgrid%lb_local + 1
1087  grids_c(ilevel) = rsgrid%buffer%c_ptr
1088  END associate
1089  END DO
1090 
1091  IF (calculate_forces) THEN
1092  forces_c = c_loc(forces(1, 1))
1093  ELSE
1094  forces_c = c_null_ptr
1095  END IF
1096 
1097  IF (calculate_virial) THEN
1098  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  cpassert(is_contiguous(forces))
1107  cpassert(is_contiguous(virial))
1108 #endif
1109 
1110  cpassert(SIZE(forces, 1) == 3)
1111  cpassert(c_associated(task_list%c_ptr))
1112  cpassert(c_associated(hab_blocks%c_ptr))
1113  cpassert(c_associated(pab_blocks%c_ptr) .OR. .NOT. calculate_forces)
1114  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  virial=virial_c)
1126 
1127  CALL timestop(handle)
1128  END SUBROUTINE grid_integrate_task_list
1129 
1130 ! **************************************************************************************************
1131 !> \brief Initialize grid library
1132 !> \author Ole Schuett
1133 ! **************************************************************************************************
1134  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  CALL grid_library_init_c()
1141 
1142  END SUBROUTINE grid_library_init
1143 
1144 ! **************************************************************************************************
1145 !> \brief Finalize grid library
1146 !> \author Ole Schuett
1147 ! **************************************************************************************************
1148  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  CALL grid_library_finalize_c()
1155 
1156  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  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  apply_cutoff=logical(apply_cutoff, c_bool))
1182 
1183  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  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  output_unit=output_unit)
1211 
1212  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  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  CALL my_mpi_comm%set_handle(int(mpi_comm))
1228 
1229  CALL my_mpi_comm%sum(number)
1230  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  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  IF (output_unit <= 0) &
1246  RETURN
1247 
1248  ! Convert C char array into Fortran string.
1249  nchars = strlcpy_c2f(buffer, message)
1250 
1251  ! Print the message.
1252  WRITE (output_unit, fmt="(A)", advance="NO") buffer(1:nchars)
1253  END SUBROUTINE print_func
1254 
1255 END MODULE grid_api
print_func
static void print_func(char *message, int output_unit)
Wrapper for printf, passed to dbm_library_print_stats.
Definition: dbm_miniapp.c:28
grid_create_basis_set
void grid_create_basis_set(const int nset, const int nsgf, const int maxco, const int maxpgf, const int lmin[nset], const int lmax[nset], const int npgf[nset], const int nsgf_set[nset], const int first_sgf[nset], const double sphi[nsgf][maxco], const double zet[nset][maxpgf], grid_basis_set **basis_set_out)
Allocates a basis set which can be passed to grid_create_task_list. See grid_task_list....
Definition: grid_basis_set.c:18
grid_free_basis_set
void grid_free_basis_set(grid_basis_set *basis_set)
Deallocates given basis set.
Definition: grid_basis_set.c:58
apply_cutoff
void apply_cutoff(void *ptr)
Definition: grid_dgemm_context.c:484
grid_library_finalize
void grid_library_finalize(void)
Finalizes the grid library.
Definition: grid_library.c:83
grid_library_print_stats
void grid_library_print_stats(void(*mpi_sum_func)(long *, int), const int mpi_comm, void(*print_func)(char *, int), const int output_unit)
Prints statistics gathered by the grid library.
Definition: grid_library.c:155
grid_library_init
void grid_library_init(void)
Initializes the grid library.
Definition: grid_library.c:49
grid_library_set_config
void grid_library_set_config(const enum grid_backend backend, const bool validate, const bool apply_cutoff)
Configures the grid library.
Definition: grid_library.c:112
grid_api
Fortran API for the grid package, which is written in C.
Definition: grid_api.F:12
grid_api::grid_func_adbmdab_z
integer, parameter, public grid_func_adbmdab_z
Definition: grid_api.F:33
grid_api::grid_func_core_x
integer, parameter, public grid_func_core_x
Definition: grid_api.F:62
grid_api::grid_func_adbmdab_y
integer, parameter, public grid_func_adbmdab_y
Definition: grid_api.F:32
grid_api::grid_func_ardbmdarb_yx
integer, parameter, public grid_func_ardbmdarb_yx
Definition: grid_api.F:37
grid_api::grid_func_dab_z
integer, parameter, public grid_func_dab_z
Definition: grid_api.F:57
grid_api::grid_collocate_task_list
subroutine, public grid_collocate_task_list(task_list, ga_gb_function, pab_blocks, rs_grids)
Collocate all tasks of in given list onto given grids.
Definition: grid_api.F:968
grid_api::grid_func_dzdx
integer, parameter, public grid_func_dzdx
Definition: grid_api.F:51
grid_api::grid_func_ardbmdarb_zz
integer, parameter, public grid_func_ardbmdarb_zz
Definition: grid_api.F:42
grid_api::grid_backend_auto
integer, parameter, public grid_backend_auto
Definition: grid_api.F:66
grid_api::grid_backend_gpu
integer, parameter, public grid_backend_gpu
Definition: grid_api.F:70
grid_api::grid_free_task_list
subroutine, public grid_free_task_list(task_list)
Deallocates given task list, basis_sets have to be freed separately.
Definition: grid_api.F:935
grid_api::grid_func_dzdz
integer, parameter, public grid_func_dzdz
Definition: grid_api.F:54
grid_api::grid_func_dydz
integer, parameter, public grid_func_dydz
Definition: grid_api.F:50
grid_api::grid_func_adb_y
integer, parameter, public grid_func_adb_y
Definition: grid_api.F:59
grid_api::grid_func_dxdy
integer, parameter, public grid_func_dxdy
Definition: grid_api.F:49
grid_api::grid_func_dabpadb_y
integer, parameter, public grid_func_dabpadb_y
Definition: grid_api.F:44
grid_api::grid_func_ardbmdarb_xy
integer, parameter, public grid_func_ardbmdarb_xy
Definition: grid_api.F:35
grid_api::grid_func_dab_y
integer, parameter, public grid_func_dab_y
Definition: grid_api.F:56
grid_api::grid_backend_hip
integer, parameter, public grid_backend_hip
Definition: grid_api.F:71
grid_api::grid_create_task_list
subroutine, public grid_create_task_list(ntasks, natoms, nkinds, nblocks, block_offsets, atom_positions, atom_kinds, basis_sets, level_list, iatom_list, jatom_list, iset_list, jset_list, ipgf_list, jpgf_list, border_mask_list, block_num_list, radius_list, rab_list, rs_grids, task_list)
Allocates a task list which can be passed to grid_collocate_task_list.
Definition: grid_api.F:726
grid_api::grid_func_adb_z
integer, parameter, public grid_func_adb_z
Definition: grid_api.F:60
grid_api::grid_func_ardbmdarb_zx
integer, parameter, public grid_func_ardbmdarb_zx
Definition: grid_api.F:40
grid_api::grid_func_adb_x
integer, parameter, public grid_func_adb_x
Definition: grid_api.F:58
grid_api::grid_func_dxdx
integer, parameter, public grid_func_dxdx
Definition: grid_api.F:52
grid_api::grid_func_ardbmdarb_xx
integer, parameter, public grid_func_ardbmdarb_xx
Definition: grid_api.F:34
grid_api::grid_func_dadb
integer, parameter, public grid_func_dadb
Definition: grid_api.F:30
grid_api::grid_backend_dgemm
integer, parameter, public grid_backend_dgemm
Definition: grid_api.F:69
grid_api::grid_func_dydy
integer, parameter, public grid_func_dydy
Definition: grid_api.F:53
grid_api::grid_func_dabpadb_z
integer, parameter, public grid_func_dabpadb_z
Definition: grid_api.F:45
grid_api::grid_backend_cpu
integer, parameter, public grid_backend_cpu
Definition: grid_api.F:68
grid_api::grid_func_dabpadb_x
integer, parameter, public grid_func_dabpadb_x
Definition: grid_api.F:43
grid_api::grid_func_dx
integer, parameter, public grid_func_dx
Definition: grid_api.F:46
grid_api::grid_func_dz
integer, parameter, public grid_func_dz
Definition: grid_api.F:48
grid_api::grid_func_ardbmdarb_yz
integer, parameter, public grid_func_ardbmdarb_yz
Definition: grid_api.F:39
grid_api::grid_func_ab
integer, parameter, public grid_func_ab
Definition: grid_api.F:29
grid_api::integrate_pgf_product
subroutine, public integrate_pgf_product(la_max, zeta, la_min, lb_max, zetb, lb_min, ra, rab, rsgrid, hab, pab, o1, o2, radius, calculate_forces, force_a, force_b, compute_tau, use_virial, my_virial_a, my_virial_b, hdab, hadb, a_hdab, use_subpatch, subpatch_pattern)
low level function to compute matrix elements of primitive gaussian functions
Definition: grid_api.F:279
grid_api::grid_func_ardbmdarb_yy
integer, parameter, public grid_func_ardbmdarb_yy
Definition: grid_api.F:38
grid_api::grid_integrate_task_list
subroutine, public grid_integrate_task_list(task_list, compute_tau, calculate_forces, calculate_virial, pab_blocks, rs_grids, hab_blocks, forces, virial)
Integrate all tasks of in given list from given grids.
Definition: grid_api.F:1040
grid_api::grid_func_core_y
integer, parameter, public grid_func_core_y
Definition: grid_api.F:63
grid_api::grid_backend_ref
integer, parameter, public grid_backend_ref
Definition: grid_api.F:67
grid_api::grid_func_adbmdab_x
integer, parameter, public grid_func_adbmdab_x
Definition: grid_api.F:31
grid_api::grid_func_dab_x
integer, parameter, public grid_func_dab_x
Definition: grid_api.F:55
grid_api::collocate_pgf_product
subroutine, public collocate_pgf_product(la_max, zeta, la_min, lb_max, zetb, lb_min, ra, rab, scale, pab, o1, o2, rsgrid, ga_gb_function, radius, use_subpatch, subpatch_pattern)
low level collocation of primitive gaussian functions
Definition: grid_api.F:119
grid_api::grid_func_ardbmdarb_zy
integer, parameter, public grid_func_ardbmdarb_zy
Definition: grid_api.F:41
grid_api::grid_func_core_z
integer, parameter, public grid_func_core_z
Definition: grid_api.F:64
grid_api::grid_func_dy
integer, parameter, public grid_func_dy
Definition: grid_api.F:47
grid_api::grid_func_ardbmdarb_xz
integer, parameter, public grid_func_ardbmdarb_xz
Definition: grid_api.F:36
kinds
Defines the basic variable types.
Definition: kinds.F:23
kinds::dp
integer, parameter, public dp
Definition: kinds.F:34
message_passing
Interface to the message passing library MPI.
Definition: message_passing.F:23
offload_api
Fortran API for the offload package, which is written in C.
Definition: offload_api.F:12
realspace_grid_types
Definition: realspace_grid_types.F:18
string_utilities
Utilities for string manipulations.
Definition: string_utilities.F:16
string_utilities::strlcpy_c2f
integer function, public strlcpy_c2f(fstring, cstring)
Copy the content of a \0-terminated C-string to a finite-length Fortran string.
Definition: string_utilities.F:3463