Line data Source code
1 : !--------------------------------------------------------------------------------------------------!
2 : ! CP2K: A general program to perform molecular dynamics simulations !
3 : ! Copyright 2000-2024 CP2K developers group <https://cp2k.org> !
4 : ! !
5 : ! SPDX-License-Identifier: GPL-2.0-or-later !
6 : !--------------------------------------------------------------------------------------------------!
7 :
8 : ! **************************************************************************************************
9 : !> \brief tensor index and mapping to DBM index
10 : !> \author Patrick Seewald
11 : ! **************************************************************************************************
12 : MODULE dbt_index
13 : USE dbt_allocate_wrap, ONLY: allocate_any
14 : USE kinds, ONLY: int_8
15 : #include "../base/base_uses.f90"
16 : #:include "dbt_macros.fypp"
17 :
18 : IMPLICIT NONE
19 : PRIVATE
20 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'dbt_index'
21 :
22 : PUBLIC :: &
23 : combine_tensor_index, &
24 : combine_pgrid_index, &
25 : create_nd_to_2d_mapping, &
26 : destroy_nd_to_2d_mapping, &
27 : get_2d_indices_tensor, &
28 : get_2d_indices_pgrid, &
29 : dbt_get_mapping_info, &
30 : get_nd_indices_tensor, &
31 : get_nd_indices_pgrid, &
32 : nd_to_2d_mapping, &
33 : ndims_mapping, &
34 : split_tensor_index, &
35 : split_pgrid_index, &
36 : ndims_mapping_row, &
37 : ndims_mapping_column, &
38 : dbt_inverse_order, &
39 : permute_index
40 :
41 : TYPE nd_to_2d_mapping
42 : INTEGER :: ndim_nd = -1
43 : INTEGER :: ndim1_2d = -1
44 : INTEGER :: ndim2_2d = -1
45 :
46 : INTEGER, DIMENSION(:), ALLOCATABLE :: dims_nd
47 : INTEGER(KIND=int_8), DIMENSION(2) :: dims_2d = -1
48 : INTEGER, DIMENSION(:), ALLOCATABLE :: dims1_2d
49 : INTEGER, DIMENSION(:), ALLOCATABLE :: dims2_2d
50 :
51 : INTEGER, DIMENSION(:), ALLOCATABLE :: map1_2d
52 : INTEGER, DIMENSION(:), ALLOCATABLE :: map2_2d
53 : INTEGER, DIMENSION(:), ALLOCATABLE :: map_nd
54 :
55 : INTEGER :: base = -1
56 : LOGICAL :: col_major = .FALSE.
57 : END TYPE nd_to_2d_mapping
58 :
59 : CONTAINS
60 :
61 : ! **************************************************************************************************
62 : !> \brief Create all data needed to quickly map between nd index and 2d index.
63 : !> \param map index mapping data
64 : !> \param dims nd sizes
65 : !> \param map1_2d which nd-indices map to first matrix index and in which order
66 : !> \param map2_2d which nd-indices map to second matrix index and in which order
67 : !> \param base base index (1 for Fortran-style, 0 for C-style, default is 1)
68 : !> \param col_major whether index should be column major order
69 : !> (.TRUE. for Fortran-style, .FALSE. for C-style, default is .TRUE.).
70 : !> \author Patrick Seewald
71 : ! **************************************************************************************************
72 11810116 : SUBROUTINE create_nd_to_2d_mapping(map, dims, map1_2d, map2_2d, base, col_major)
73 : TYPE(nd_to_2d_mapping), INTENT(OUT) :: map
74 : INTEGER, DIMENSION(:), INTENT(IN) :: dims, map1_2d, map2_2d
75 : INTEGER, INTENT(IN), OPTIONAL :: base
76 : LOGICAL, INTENT(IN), OPTIONAL :: col_major
77 :
78 : INTEGER :: i
79 :
80 2952529 : IF (PRESENT(col_major)) THEN
81 771146 : map%col_major = col_major
82 : ELSE
83 2181383 : map%col_major = .TRUE.
84 : END IF
85 :
86 2952529 : IF (PRESENT(base)) THEN
87 771146 : map%base = base
88 : ELSE
89 2181383 : map%base = 1
90 : END IF
91 :
92 2952529 : map%ndim1_2d = SIZE(map1_2d)
93 2952529 : map%ndim2_2d = SIZE(map2_2d)
94 2952529 : map%ndim_nd = SIZE(dims)
95 :
96 12562273 : ALLOCATE (map%map1_2d, source=map1_2d)
97 12720666 : ALLOCATE (map%map2_2d, source=map2_2d)
98 16425352 : ALLOCATE (map%dims_nd, source=dims)
99 13314430 : ALLOCATE (map%dims1_2d, source=dims(map1_2d))
100 13631216 : ALLOCATE (map%dims2_2d, source=dims(map2_2d))
101 :
102 8857587 : ALLOCATE (map%map_nd(map%ndim_nd))
103 13314430 : map%map_nd(map1_2d) = (/(i, i=1, SIZE(map1_2d))/)
104 13631216 : map%map_nd(map2_2d) = (/(i + SIZE(map1_2d), i=1, SIZE(map2_2d))/)
105 :
106 16425352 : map%dims_2d = [PRODUCT(INT(map%dims1_2d, KIND=int_8)), PRODUCT(INT(map%dims2_2d, KIND=int_8))]
107 :
108 2952529 : END SUBROUTINE create_nd_to_2d_mapping
109 :
110 : ! **************************************************************************************************
111 : !> \brief
112 : !> \author Patrick Seewald
113 : ! **************************************************************************************************
114 3368271 : SUBROUTINE destroy_nd_to_2d_mapping(map)
115 : TYPE(nd_to_2d_mapping), INTENT(INOUT) :: map
116 :
117 3368271 : DEALLOCATE (map%dims1_2d)
118 3368271 : DEALLOCATE (map%dims2_2d)
119 3368271 : DEALLOCATE (map%map1_2d)
120 3368271 : DEALLOCATE (map%map2_2d)
121 3368271 : DEALLOCATE (map%map_nd)
122 3368271 : DEALLOCATE (map%dims_nd)
123 3368271 : END SUBROUTINE destroy_nd_to_2d_mapping
124 :
125 : ! **************************************************************************************************
126 : !> \brief
127 : !> \author Patrick Seewald
128 : ! **************************************************************************************************
129 6904744 : PURE FUNCTION ndims_mapping(map)
130 : TYPE(nd_to_2d_mapping), INTENT(IN) :: map
131 : INTEGER :: ndims_mapping
132 :
133 6904744 : ndims_mapping = map%ndim_nd
134 6904744 : END FUNCTION
135 :
136 : ! **************************************************************************************************
137 : !> \brief how many tensor dimensions are mapped to matrix row
138 : !> \author Patrick Seewald
139 : ! **************************************************************************************************
140 12887566 : PURE FUNCTION ndims_mapping_row(map)
141 : TYPE(nd_to_2d_mapping), INTENT(IN) :: map
142 : INTEGER :: ndims_mapping_row
143 12887566 : ndims_mapping_row = map%ndim1_2d
144 12887566 : END FUNCTION
145 :
146 : ! **************************************************************************************************
147 : !> \brief how many tensor dimensions are mapped to matrix column
148 : !> \author Patrick Seewald
149 : ! **************************************************************************************************
150 12886978 : PURE FUNCTION ndims_mapping_column(map)
151 : TYPE(nd_to_2d_mapping), INTENT(IN) :: map
152 : INTEGER :: ndims_mapping_column
153 12886978 : ndims_mapping_column = map%ndim2_2d
154 12886978 : END FUNCTION
155 :
156 : ! **************************************************************************************************
157 : !> \brief get mapping info
158 : !> \param map index mapping data
159 : !> \param ndim_nd number of dimensions
160 : !> \param ndim1_2d number of dimensions that map to first 2d index
161 : !> \param ndim2_2d number of dimensions that map to first 2d index
162 : !> \param dims_2d 2d dimensions
163 : !> \param dims_nd nd dimensions
164 : !> \param dims1_2d dimensions that map to first 2d index
165 : !> \param dims2_2d dimensions that map to second 2d index
166 : !> \param map1_2d indices that map to first 2d index
167 : !> \param map2_2d indices that map to second 2d index
168 : !> \param map_nd inverse of [map1_2d, map2_2d]
169 : !> \param base base index
170 : !> \param col_major is index in column major order
171 : !> \author Patrick Seewald
172 : ! **************************************************************************************************
173 8478612 : PURE SUBROUTINE dbt_get_mapping_info(map, ndim_nd, ndim1_2d, ndim2_2d, dims_2d_i8, &
174 5260330 : dims_2d, dims_nd, dims1_2d, dims2_2d, &
175 9201340 : map1_2d, map2_2d, map_nd, base, col_major)
176 : TYPE(nd_to_2d_mapping), INTENT(IN) :: map
177 : INTEGER, INTENT(OUT), OPTIONAL :: ndim_nd, ndim1_2d, ndim2_2d
178 : INTEGER(KIND=int_8), DIMENSION(2), INTENT(OUT), OPTIONAL :: dims_2d_i8
179 : INTEGER, DIMENSION(2), INTENT(OUT), OPTIONAL :: dims_2d
180 : INTEGER, DIMENSION(ndims_mapping(map)), &
181 : INTENT(OUT), OPTIONAL :: dims_nd
182 : INTEGER, DIMENSION(ndims_mapping_row(map)), INTENT(OUT), &
183 : OPTIONAL :: dims1_2d
184 : INTEGER, DIMENSION(ndims_mapping_column(map)), INTENT(OUT), &
185 : OPTIONAL :: dims2_2d
186 : INTEGER, DIMENSION(ndims_mapping_row(map)), INTENT(OUT), &
187 : OPTIONAL :: map1_2d
188 : INTEGER, DIMENSION(ndims_mapping_column(map)), INTENT(OUT), &
189 : OPTIONAL :: map2_2d
190 : INTEGER, DIMENSION(ndims_mapping(map)), &
191 : INTENT(OUT), OPTIONAL :: map_nd
192 : INTEGER, INTENT(OUT), OPTIONAL :: base
193 : LOGICAL, INTENT(OUT), OPTIONAL :: col_major
194 :
195 8478612 : IF (PRESENT(ndim_nd)) ndim_nd = map%ndim_nd
196 8478612 : IF (PRESENT(ndim1_2d)) ndim1_2d = map%ndim1_2d
197 8478612 : IF (PRESENT(ndim2_2d)) ndim2_2d = map%ndim2_2d
198 11159182 : IF (PRESENT(dims_2d_i8)) dims_2d_i8(:) = map%dims_2d(:)
199 11553160 : IF (PRESENT(dims_2d)) dims_2d(:) = INT(map%dims_2d(:))
200 8478612 : IF (PRESENT(dims_nd)) THEN
201 12089248 : dims_nd(:) = map%dims_nd(:)
202 : END IF
203 8478612 : IF (PRESENT(dims1_2d)) THEN
204 1996449 : dims1_2d(:) = map%dims1_2d
205 : END IF
206 8478612 : IF (PRESENT(dims2_2d)) THEN
207 2275654 : dims2_2d(:) = map%dims2_2d
208 : END IF
209 8478612 : IF (PRESENT(map1_2d)) THEN
210 10713931 : map1_2d(:) = map%map1_2d
211 : END IF
212 8478612 : IF (PRESENT(map2_2d)) THEN
213 10479437 : map2_2d(:) = map%map2_2d
214 : END IF
215 8478612 : IF (PRESENT(map_nd)) THEN
216 0 : map_nd(:) = map%map_nd(:)
217 : END IF
218 8478612 : IF (PRESENT(base)) THEN
219 0 : base = map%base
220 : END IF
221 8478612 : IF (PRESENT(col_major)) THEN
222 0 : col_major = map%col_major
223 : END IF
224 :
225 22940282 : END SUBROUTINE dbt_get_mapping_info
226 :
227 : ! **************************************************************************************************
228 : !> \brief transform nd index to flat index
229 : !> \param ind_in nd index
230 : !> \param dims nd dimensions
231 : !> \param ind_out flat index
232 : !> \author Patrick Seewald
233 : ! **************************************************************************************************
234 145505408 : PURE FUNCTION combine_tensor_index(ind_in, dims) RESULT(ind_out)
235 : INTEGER, DIMENSION(:), INTENT(IN) :: ind_in, dims
236 : INTEGER(KIND=int_8) :: ind_out
237 : INTEGER :: i_dim
238 :
239 145505408 : ind_out = ind_in(SIZE(dims))
240 221849254 : DO i_dim = SIZE(dims) - 1, 1, -1
241 221849254 : ind_out = (ind_out - 1)*dims(i_dim) + ind_in(i_dim)
242 : END DO
243 :
244 145505408 : END FUNCTION
245 :
246 : ! **************************************************************************************************
247 : !> \brief transform nd index to flat index
248 : !> \param ind_in nd index
249 : !> \param dims nd dimensions
250 : !> \param ind_out flat index
251 : !> \author Patrick Seewald
252 : ! **************************************************************************************************
253 29229884 : PURE FUNCTION combine_pgrid_index(ind_in, dims) RESULT(ind_out)
254 : INTEGER, DIMENSION(:), INTENT(IN) :: ind_in, dims
255 : INTEGER :: ind_out
256 :
257 : INTEGER :: i_dim
258 :
259 29229884 : ind_out = ind_in(1)
260 42679197 : DO i_dim = 2, SIZE(dims)
261 42679197 : ind_out = ind_out*dims(i_dim) + ind_in(i_dim)
262 : END DO
263 29229884 : END FUNCTION
264 :
265 : ! **************************************************************************************************
266 : !> \brief transform flat index to nd index
267 : !> \param ind_in flat index
268 : !> \param dims nd dimensions
269 : !> \param ind_out nd index
270 : !> \author Patrick Seewald
271 : ! **************************************************************************************************
272 112711755 : PURE FUNCTION split_tensor_index(ind_in, dims) RESULT(ind_out)
273 : INTEGER(KIND=int_8), INTENT(IN) :: ind_in
274 : INTEGER, DIMENSION(:), INTENT(IN) :: dims
275 : INTEGER, DIMENSION(SIZE(dims)) :: ind_out
276 :
277 : INTEGER(KIND=int_8) :: tmp
278 : INTEGER :: i_dim
279 :
280 112711755 : tmp = ind_in
281 272779298 : DO i_dim = 1, SIZE(dims)
282 160067543 : ind_out(i_dim) = INT(MOD(tmp - 1, INT(dims(i_dim), int_8)) + 1)
283 272779298 : tmp = (tmp - 1)/dims(i_dim) + 1
284 : END DO
285 :
286 112711755 : END FUNCTION
287 :
288 : ! **************************************************************************************************
289 : !> \brief transform flat index to nd index
290 : !> \param ind_in flat index
291 : !> \param dims nd dimensions
292 : !> \param ind_out nd index
293 : !> \author Patrick Seewald
294 : ! **************************************************************************************************
295 2255436 : PURE FUNCTION split_pgrid_index(ind_in, dims) RESULT(ind_out)
296 : INTEGER, INTENT(IN) :: ind_in
297 : INTEGER, DIMENSION(:), INTENT(IN) :: dims
298 : INTEGER, DIMENSION(SIZE(dims)) :: ind_out
299 :
300 : INTEGER :: tmp
301 : INTEGER :: i_dim
302 :
303 2255436 : tmp = ind_in
304 5130350 : DO i_dim = SIZE(dims), 1, -1
305 2874914 : ind_out(i_dim) = MOD(tmp, dims(i_dim))
306 5130350 : tmp = tmp/dims(i_dim)
307 : END DO
308 2255436 : END FUNCTION
309 :
310 : ! **************************************************************************************************
311 : !> \brief transform nd index to 2d index, using info from index mapping.
312 : !> \param map index mapping
313 : !> \param ind_in nd index
314 : !> \param ind_out 2d index
315 : !> \author Patrick Seewald
316 : ! **************************************************************************************************
317 60076986 : PURE FUNCTION get_2d_indices_tensor(map, ind_in) RESULT(ind_out)
318 : TYPE(nd_to_2d_mapping), INTENT(IN) :: map
319 : INTEGER, DIMENSION(map%ndim_nd), INTENT(IN) :: ind_in
320 : INTEGER(KIND=int_8), DIMENSION(2) :: ind_out
321 : INTEGER :: i
322 : INTEGER, DIMENSION(${maxrank}$) :: ind_tmp
323 :
324 130336540 : DO i = 1, map%ndim1_2d
325 130336540 : ind_tmp(i) = ind_in(map%map1_2d(i))
326 : END DO
327 60076986 : ind_out(1) = combine_tensor_index(ind_tmp(:map%ndim1_2d), map%dims1_2d)
328 :
329 159779796 : DO i = 1, map%ndim2_2d
330 159779796 : ind_tmp(i) = ind_in(map%map2_2d(i))
331 : END DO
332 60076986 : ind_out(2) = combine_tensor_index(ind_tmp(:map%ndim2_2d), map%dims2_2d)
333 60076986 : END FUNCTION
334 :
335 : ! **************************************************************************************************
336 : !> \brief transform nd index to 2d index, using info from index mapping.
337 : !> \param map index mapping
338 : !> \param ind_in nd index
339 : !> \param ind_out 2d index
340 : !> \author Patrick Seewald
341 : ! **************************************************************************************************
342 0 : PURE FUNCTION get_2d_indices_pgrid(map, ind_in) RESULT(ind_out)
343 : TYPE(nd_to_2d_mapping), INTENT(IN) :: map
344 : INTEGER, DIMENSION(map%ndim_nd), INTENT(IN) :: ind_in
345 : INTEGER, DIMENSION(2) :: ind_out
346 : INTEGER :: i
347 : INTEGER, DIMENSION(${maxrank}$) :: ind_tmp
348 :
349 0 : DO i = 1, map%ndim1_2d
350 0 : ind_tmp(i) = ind_in(map%map1_2d(i))
351 : END DO
352 0 : ind_out(1) = combine_pgrid_index(ind_tmp(:map%ndim1_2d), map%dims1_2d)
353 :
354 0 : DO i = 1, map%ndim2_2d
355 0 : ind_tmp(i) = ind_in(map%map2_2d(i))
356 : END DO
357 0 : ind_out(2) = combine_pgrid_index(ind_tmp(:map%ndim2_2d), map%dims2_2d)
358 0 : END FUNCTION
359 :
360 : ! **************************************************************************************************
361 : !> \brief transform 2d index to nd index, using info from index mapping.
362 : !> \param map index mapping
363 : !> \param ind_in 2d index
364 : !> \param ind_out nd index
365 : !> \author Patrick Seewald
366 : ! **************************************************************************************************
367 34894236 : PURE FUNCTION get_nd_indices_tensor(map, ind_in) RESULT(ind_out)
368 : TYPE(nd_to_2d_mapping), INTENT(IN) :: map
369 : INTEGER(KIND=int_8), DIMENSION(2), INTENT(IN) :: ind_in
370 : INTEGER, DIMENSION(map%ndim_nd) :: ind_out
371 : INTEGER, DIMENSION(${maxrank}$) :: ind_tmp
372 : INTEGER :: i
373 :
374 34894236 : ind_tmp(:map%ndim1_2d) = split_tensor_index(ind_in(1), map%dims1_2d)
375 :
376 79589164 : DO i = 1, map%ndim1_2d
377 79589164 : ind_out(map%map1_2d(i)) = ind_tmp(i)
378 : END DO
379 :
380 34894236 : ind_tmp(:map%ndim2_2d) = split_tensor_index(ind_in(2), map%dims2_2d)
381 :
382 89323784 : DO i = 1, map%ndim2_2d
383 89323784 : ind_out(map%map2_2d(i)) = ind_tmp(i)
384 : END DO
385 :
386 34894236 : END FUNCTION
387 :
388 : ! **************************************************************************************************
389 : !> \brief transform 2d index to nd index, using info from index mapping.
390 : !> \param map index mapping
391 : !> \param ind_in 2d index
392 : !> \param ind_out nd index
393 : !> \author Patrick Seewald
394 : ! **************************************************************************************************
395 728012 : PURE FUNCTION get_nd_indices_pgrid(map, ind_in) RESULT(ind_out)
396 : TYPE(nd_to_2d_mapping), INTENT(IN) :: map
397 : INTEGER, DIMENSION(2), INTENT(IN) :: ind_in
398 : INTEGER, DIMENSION(map%ndim_nd) :: ind_out
399 :
400 3162424 : ind_out(map%map1_2d) = split_pgrid_index(ind_in(1), map%dims1_2d)
401 3399924 : ind_out(map%map2_2d) = split_pgrid_index(ind_in(2), map%dims2_2d)
402 :
403 728012 : END FUNCTION
404 :
405 : ! **************************************************************************************************
406 : !> \brief Invert order
407 : !> \author Patrick Seewald
408 : ! **************************************************************************************************
409 2662852 : PURE FUNCTION dbt_inverse_order(order)
410 : INTEGER, DIMENSION(:), INTENT(IN) :: order
411 : INTEGER, DIMENSION(SIZE(order)) :: dbt_inverse_order
412 :
413 : INTEGER :: i
414 :
415 19579138 : dbt_inverse_order(order) = (/(i, i=1, SIZE(order))/)
416 : END FUNCTION
417 :
418 : ! **************************************************************************************************
419 : !> \brief reorder tensor index (no data)
420 : !> \author Patrick Seewald
421 : ! **************************************************************************************************
422 5661756 : SUBROUTINE permute_index(map_in, map_out, order)
423 : TYPE(nd_to_2d_mapping), INTENT(IN) :: map_in
424 : TYPE(nd_to_2d_mapping), INTENT(OUT) :: map_out
425 : INTEGER, DIMENSION(ndims_mapping(map_in)), &
426 : INTENT(IN) :: order
427 :
428 : INTEGER :: ndim_nd
429 2830878 : INTEGER, DIMENSION(ndims_mapping_row(map_in)) :: map1_2d, map1_2d_reorder
430 2830878 : INTEGER, DIMENSION(ndims_mapping_column(map_in)) :: map2_2d, map2_2d_reorder
431 1415439 : INTEGER, DIMENSION(ndims_mapping(map_in)) :: dims_nd, dims_reorder
432 :
433 1415439 : CALL dbt_get_mapping_info(map_in, ndim_nd, dims_nd=dims_nd, map1_2d=map1_2d, map2_2d=map2_2d)
434 :
435 5196834 : dims_reorder(order) = dims_nd
436 :
437 3452193 : map1_2d_reorder(:) = order(map1_2d)
438 3160080 : map2_2d_reorder(:) = order(map2_2d)
439 :
440 1415439 : CALL create_nd_to_2d_mapping(map_out, dims_reorder, map1_2d_reorder, map2_2d_reorder)
441 1415439 : END SUBROUTINE
442 :
443 0 : END MODULE dbt_index
|