(git:b195825)
dct.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: GPL-2.0-or-later !
6 !--------------------------------------------------------------------------------------------------!
7 
8 ! **************************************************************************************************
9 !> \brief the type I Discrete Cosine Transform (DCT-I)
10 !> \par History
11 !> 07.2014 created [Hossein Bani-Hashemian]
12 !> 11.2015 dealt with periodic grids [Hossein Bani-Hashemian]
13 !> 03.2016 dct in one or two directions [Hossein Bani-Hashemian]
14 !> \author Mohammad Hossein Bani-Hashemian
15 ! **************************************************************************************************
16 MODULE dct
17 
18  USE kinds, ONLY: dp
19  USE message_passing, ONLY: mp_comm_type,&
20  mp_request_null,&
21  mp_request_type,&
22  mp_waitall
23  USE pw_grid_types, ONLY: pw_grid_type
24  USE pw_grids, ONLY: pw_grid_create,&
25  pw_grid_setup
26  USE pw_types, ONLY: pw_r3d_rs_type
27 #include "../base/base_uses.f90"
28 
29  IMPLICIT NONE
30  PRIVATE
31  CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'dct'
32 
33  TYPE :: dct_type
34  INTEGER, DIMENSION(:), POINTER :: dests_expand => null()
35  INTEGER, DIMENSION(:), POINTER :: srcs_expand => null()
36  INTEGER, DIMENSION(:), POINTER :: flipg_stat => null()
37  INTEGER, DIMENSION(:), POINTER :: dests_shrink => null()
38  INTEGER :: srcs_shrink = 0
39  INTEGER, DIMENSION(:, :, :), POINTER :: recv_msgs_bnds => null()
40  INTEGER, DIMENSION(2, 3) :: dct_bounds = 0
41  INTEGER, DIMENSION(2, 3) :: dct_bounds_local = 0
42  INTEGER, DIMENSION(2, 3) :: bounds_shftd = 0
43  INTEGER, DIMENSION(2, 3) :: bounds_local_shftd = 0
44  END TYPE dct_type
45 
46  TYPE dct_msg_type
47  PRIVATE
48  REAL(dp), DIMENSION(:, :, :), ALLOCATABLE :: msg
49  END TYPE dct_msg_type
50 
51  PUBLIC dct_type, &
52  dct_type_init, &
53  dct_type_release, &
54  setup_dct_pw_grids, &
55  pw_shrink, &
56  pw_expand
57 
58  INTEGER, PARAMETER, PRIVATE :: NOT_FLIPPED = 0, &
59  ud_flipped = 1, &
60  lr_flipped = 2, &
61  bf_flipped = 3, &
62  rotated = 4
63 
64  INTEGER, PARAMETER, PUBLIC :: neumannxyz = 111, &
65  neumannxy = 110, &
66  neumannxz = 101, &
67  neumannyz = 011, &
68  neumannx = 100, &
69  neumanny = 010, &
70  neumannz = 001
71 
72 CONTAINS
73 
74 ! **************************************************************************************************
75 !> \brief Initializes a dct_type
76 !> \param pw_grid the original plane wave grid
77 !> \param neumann_directions directions in which dct should be performed
78 !> \param dct_env dct_type to be initialized
79 !> \par History
80 !> 08.2014 created [Hossein Bani-Hashemian]
81 !> \author Mohammad Hossein Bani-Hashemian
82 ! **************************************************************************************************
83  SUBROUTINE dct_type_init(pw_grid, neumann_directions, dct_env)
84 
85  TYPE(pw_grid_type), INTENT(IN), POINTER :: pw_grid
86  INTEGER, INTENT(IN) :: neumann_directions
87  TYPE(dct_type), INTENT(INOUT) :: dct_env
88 
89  CHARACTER(len=*), PARAMETER :: routinen = 'dct_type_init'
90 
91  INTEGER :: handle, maxn_sendrecv
92 
93  CALL timeset(routinen, handle)
94 
95  SELECT CASE (neumann_directions)
96  CASE (neumannxyz, neumannxy)
97  maxn_sendrecv = 4
98  CASE (neumannx, neumanny, neumannxz, neumannyz)
99  maxn_sendrecv = 2
100  CASE (neumannz)
101  maxn_sendrecv = 1
102  CASE DEFAULT
103  cpabort("Invalid combination of Neumann and periodic conditions.")
104  END SELECT
105 
106  ALLOCATE (dct_env%flipg_stat(maxn_sendrecv))
107  ALLOCATE (dct_env%dests_expand(maxn_sendrecv), dct_env%srcs_expand(maxn_sendrecv))
108  ALLOCATE (dct_env%dests_shrink(maxn_sendrecv))
109  ALLOCATE (dct_env%recv_msgs_bnds(2, 3, maxn_sendrecv))
110 
111  CALL set_dests_srcs_pid(pw_grid, neumann_directions, &
112  dct_env%dests_expand, dct_env%srcs_expand, dct_env%flipg_stat, &
113  dct_env%dests_shrink, dct_env%srcs_shrink)
114  CALL expansion_bounds(pw_grid, neumann_directions, &
115  dct_env%srcs_expand, dct_env%flipg_stat, &
116  dct_env%bounds_shftd, dct_env%bounds_local_shftd, &
117  dct_env%recv_msgs_bnds, dct_env%dct_bounds, &
118  dct_env%dct_bounds_local)
119 
120  CALL timestop(handle)
121 
122  END SUBROUTINE dct_type_init
123 
124 ! **************************************************************************************************
125 !> \brief Releases a dct_type
126 !> \param dct_env dct_type to be released
127 !> \par History
128 !> 03.2016 created [Hossein Bani-Hashemian]
129 !> \author Mohammad Hossein Bani-Hashemian
130 ! **************************************************************************************************
131  SUBROUTINE dct_type_release(dct_env)
132 
133  TYPE(dct_type), INTENT(INOUT) :: dct_env
134 
135  CHARACTER(len=*), PARAMETER :: routinen = 'dct_type_release'
136 
137  INTEGER :: handle
138 
139  CALL timeset(routinen, handle)
140 
141  IF (ASSOCIATED(dct_env%dests_shrink)) DEALLOCATE (dct_env%dests_shrink)
142  IF (ASSOCIATED(dct_env%dests_expand)) DEALLOCATE (dct_env%dests_expand)
143  IF (ASSOCIATED(dct_env%srcs_expand)) DEALLOCATE (dct_env%srcs_expand)
144  IF (ASSOCIATED(dct_env%flipg_stat)) DEALLOCATE (dct_env%flipg_stat)
145  IF (ASSOCIATED(dct_env%recv_msgs_bnds)) DEALLOCATE (dct_env%recv_msgs_bnds)
146 
147  CALL timestop(handle)
148 
149  END SUBROUTINE dct_type_release
150 
151 ! **************************************************************************************************
152 !> \brief sets up an extended pw_grid for Discrete Cosine Transform (DCT)
153 !> calculations
154 !> \param pw_grid the original plane wave grid
155 !> \param cell_hmat cell hmat
156 !> \param neumann_directions directions in which dct should be performed
157 !> \param dct_pw_grid DCT plane-wave grid
158 !> \par History
159 !> 07.2014 created [Hossein Bani-Hashemian]
160 !> \author Mohammad Hossein Bani-Hashemian
161 ! **************************************************************************************************
162  SUBROUTINE setup_dct_pw_grids(pw_grid, cell_hmat, neumann_directions, dct_pw_grid)
163 
164  TYPE(pw_grid_type), INTENT(IN), POINTER :: pw_grid
165  REAL(dp), DIMENSION(3, 3), INTENT(IN) :: cell_hmat
166  INTEGER, INTENT(IN) :: neumann_directions
167  TYPE(pw_grid_type), INTENT(INOUT), POINTER :: dct_pw_grid
168 
169  CHARACTER(LEN=*), PARAMETER :: routinen = 'setup_dct_pw_grids'
170 
171  INTEGER :: blocked, handle, maxn_sendrecv, &
172  srcs_shrink
173  INTEGER, DIMENSION(2, 3) :: bounds_local_new, bounds_local_shftd, &
174  bounds_new, bounds_shftd
175  INTEGER, DIMENSION(:), POINTER :: dests_expand, dests_shrink, flipg_stat, &
176  srcs_expand
177  INTEGER, DIMENSION(:, :, :), POINTER :: recv_msgs_bnds
178  REAL(kind=dp), DIMENSION(3) :: scfac
179  REAL(kind=dp), DIMENSION(3, 3) :: hmat2
180 
181  CALL timeset(routinen, handle)
182 
183  SELECT CASE (neumann_directions)
184  CASE (neumannxyz)
185  maxn_sendrecv = 4
186  scfac = (/2.0_dp, 2.0_dp, 2.0_dp/)
187  CASE (neumannxy)
188  maxn_sendrecv = 4
189  scfac = (/2.0_dp, 2.0_dp, 1.0_dp/)
190  CASE (neumannxz)
191  maxn_sendrecv = 2
192  scfac = (/2.0_dp, 1.0_dp, 2.0_dp/)
193  CASE (neumannyz)
194  maxn_sendrecv = 2
195  scfac = (/1.0_dp, 2.0_dp, 2.0_dp/)
196  CASE (neumannx)
197  maxn_sendrecv = 2
198  scfac = (/2.0_dp, 1.0_dp, 1.0_dp/)
199  CASE (neumanny)
200  maxn_sendrecv = 2
201  scfac = (/1.0_dp, 2.0_dp, 1.0_dp/)
202  CASE (neumannz)
203  maxn_sendrecv = 1
204  scfac = (/1.0_dp, 1.0_dp, 2.0_dp/)
205  CASE DEFAULT
206  cpabort("Invalid combination of Neumann and periodic conditions.")
207  END SELECT
208 
209  ALLOCATE (flipg_stat(maxn_sendrecv))
210  ALLOCATE (dests_expand(maxn_sendrecv), srcs_expand(maxn_sendrecv), dests_shrink(maxn_sendrecv))
211  ALLOCATE (recv_msgs_bnds(2, 3, maxn_sendrecv))
212 
213  CALL set_dests_srcs_pid(pw_grid, neumann_directions, dests_expand, srcs_expand, flipg_stat, &
214  dests_shrink, srcs_shrink)
215  CALL expansion_bounds(pw_grid, neumann_directions, srcs_expand, flipg_stat, &
216  bounds_shftd, bounds_local_shftd, recv_msgs_bnds, bounds_new, bounds_local_new)
217  CALL pw_grid_create(dct_pw_grid, pw_grid%para%rs_group, local=.false.)
218 
219  hmat2 = 0.0_dp
220  hmat2(1, 1) = scfac(1)*cell_hmat(1, 1)
221  hmat2(2, 2) = scfac(2)*cell_hmat(2, 2)
222  hmat2(3, 3) = scfac(3)*cell_hmat(3, 3)
223 
224  ! uses bounds_local_new that is 2*n-2 in size....this is only rarely fft-able by fftsg, and needs fftw3,
225  ! where it might use sizes that are not declared available in fft_get_radix.
226 
227  IF (pw_grid%para%blocked) THEN
228  blocked = 1
229  ELSE IF (pw_grid%para%ray_distribution) THEN
230  blocked = 0
231  END IF
232 
233  CALL pw_grid_setup(hmat2, dct_pw_grid, &
234  bounds=bounds_new, &
235  rs_dims=pw_grid%para%rs_dims, &
236  blocked=blocked, &
237  bounds_local=bounds_local_new)
238 
239  DEALLOCATE (flipg_stat, dests_expand, srcs_expand, dests_shrink, recv_msgs_bnds)
240 
241  CALL timestop(handle)
242 
243  END SUBROUTINE setup_dct_pw_grids
244 
245 ! **************************************************************************************************
246 !> \brief Finds the process ids for mpi_isend destiations and mpi_irecv sources
247 !> for expanding and shrinking a pw_r3d_rs_type data
248 !> \param pw_grid the original plane wave grid
249 !> \param neumann_directions directions in which dct should be performed
250 !> \param dests_expand list of the destination processes (pw_expand)
251 !> \param srcs_expand list of the source processes (pw_expand)
252 !> \param flipg_stat flipping status for the received data chunks (pw_expand)
253 !> \param dests_shrink list of the destination processes (pw_shrink)
254 !> \param srcs_shrink list of the source processes (pw_shrink)
255 !> \par History
256 !> 07.2014 created [Hossein Bani-Hashemian]
257 !> \author Mohammad Hossein Bani-Hashemian
258 ! **************************************************************************************************
259  SUBROUTINE set_dests_srcs_pid(pw_grid, neumann_directions, dests_expand, srcs_expand, &
260  flipg_stat, dests_shrink, srcs_shrink)
261 
262  TYPE(pw_grid_type), INTENT(IN), POINTER :: pw_grid
263  INTEGER, INTENT(IN) :: neumann_directions
264  INTEGER, DIMENSION(:), INTENT(INOUT), POINTER :: dests_expand, srcs_expand, flipg_stat, &
265  dests_shrink
266  INTEGER, INTENT(OUT) :: srcs_shrink
267 
268  CHARACTER(LEN=*), PARAMETER :: routinen = 'set_dests_srcs_pid'
269 
270  INTEGER :: group_size, handle, i, j, k, &
271  maxn_sendrecv, rs_dim1, rs_dim2, &
272  rs_mpo, tmp_size1, tmp_size2
273  INTEGER, ALLOCATABLE, DIMENSION(:) :: src_pos1_onesdd, src_pos2_onesdd, &
274  tmp1_arr, tmp2_arr
275  INTEGER, ALLOCATABLE, DIMENSION(:, :) :: dests_shrink_all, src_pos1, src_pos2, &
276  srcs_coord, srcs_expand_all
277  INTEGER, DIMENSION(2) :: rs_dims, rs_pos
278  TYPE(mp_comm_type) :: rs_group
279 
280  CALL timeset(routinen, handle)
281 
282 ! example: 3x4 process grid
283 ! XYZ or XY
284 ! rs_dim1 = 3 --> src_pos1 = [1 3 -2]
285 ! [2 -3 -1]
286 ! rs_dim2 = 4 --> src_pos2 = [1 3 -4 -2]
287 ! [2 4 -3 -1]
288 ! => (1,1) receives from (1,1) ; (1,2) ; (2,1) ; (2,2) | flipping status 0 0 0 0
289 ! (2,4) receives from (3,2) ; (3,1) ; (3,2) ; (3,1) | flipping status 2 2 4 4
290 ! and so on ...
291 ! or equivalently
292 ! => 0 receives from 0 ; 1 ; 4 ; 5 | flipping status 0 0 0 0
293 ! 7 receives from 9 ; 8 ; 9 ; 8 | flipping status 2 2 4 4
294 ! from srcs_coord :
295 ! => rs_mpo = 0 -> rs_pos = 0,0 -> srcs_coord = [ 1 1 2 2] -> 0(0) 1(0) 4(0) 5(0)
296 ! [ 1 2 1 2]
297 ! rs_mpo = 7 -> rs_pos = 1,3 -> srcs_coord = [ 3 3 -3 -3] -> 9(2) 8(2) 9(4) 8(4)
298 ! [-2 -1 -2 -1]
299 ! schematically :
300 ! ij : coordinates in a 2D process grid (starting from 1 just for demonstration)
301 ! () : to be flipped from left to right
302 ! [] : to be flipped from up to down
303 ! {} : to be rotated 180 degrees
304 ! 11 | 12 | 13 | 14 11 12 | 13 14 | (14) (13) | (12) (11)
305 ! ----------------- ==> 21 22 | 23 24 | (24) (23) | (22) (21)
306 ! 21 | 22 | 23 | 24 ---------------------------------------------
307 ! ----------------- 31 32 | 33 34 | (34) (33) | (32) (31)
308 ! 31 | 32 | 33 | 34 [31] [32] | [33] [34] | {34} {33} | {32} {31}
309 ! ---------------------------------------------
310 ! [21] [22] | [23] [24] | {24} {23} | {22} {21}
311 ! [11] [12] | [13] [14] | {14} {13} | {12} {11}
312 ! one(two)-sided :
313 ! YZ or Y
314 ! rs_dim1 = 3 --> src_pos1 = [1 2 3]
315 ! rs_dim2 = 4 --> src_pos2 = [1 3 -4 -2]
316 ! [2 4 -3 -1]
317 ! XZ or X
318 ! rs_dim1 = 3 --> src_pos1 = [1 3 -2]
319 ! [2 -3 -1]
320 ! rs_dim2 = 4 --> src_pos2 = [1 2 3 4]
321 ! Z
322 ! rs_dim1 = 3 --> src_pos1 = [1 2 3]
323 ! rs_dim2 = 4 --> src_pos2 = [1 2 3 4]
324 
325  rs_group = pw_grid%para%rs_group
326  rs_mpo = pw_grid%para%rs_mpo
327  group_size = pw_grid%para%group_size
328  rs_dims = pw_grid%para%rs_dims
329 
330  rs_pos = pw_grid%para%rs_pos
331  rs_dim1 = rs_dims(1); rs_dim2 = rs_dims(2)
332 
333 ! prepare srcs_coord
334  SELECT CASE (neumann_directions)
335  CASE (neumannxyz, neumannxy)
336  maxn_sendrecv = 4
337  ALLOCATE (srcs_coord(2, maxn_sendrecv))
338 
339  IF (mod(rs_dim1, 2) .EQ. 0) THEN
340  tmp_size1 = rs_dim1
341  ELSE
342  tmp_size1 = rs_dim1 + 1
343  END IF
344  ALLOCATE (tmp1_arr(tmp_size1), src_pos1(2, 0:rs_dim1 - 1))
345 
346  IF (mod(rs_dim1, 2) .EQ. 0) THEN
347  tmp1_arr(:) = (/(i, i=1, rs_dim1)/)
348  src_pos1(:, :) = reshape((/tmp1_arr, -tmp1_arr(tmp_size1:1:-1)/), (/2, rs_dim1/))
349  ELSE
350  tmp1_arr(:) = (/(i, i=1, rs_dim1), -rs_dim1/)
351  src_pos1(:, :) = reshape((/tmp1_arr, -tmp1_arr(tmp_size1 - 2:1:-1)/), (/2, rs_dim1/))
352  END IF
353 !---
354  IF (mod(rs_dim2, 2) .EQ. 0) THEN
355  tmp_size2 = rs_dim2
356  ELSE
357  tmp_size2 = rs_dim2 + 1
358  END IF
359  ALLOCATE (tmp2_arr(tmp_size2), src_pos2(2, 0:rs_dim2 - 1))
360 
361  IF (mod(rs_dim2, 2) .EQ. 0) THEN
362  tmp2_arr(:) = (/(i, i=1, rs_dim2)/)
363  src_pos2(:, :) = reshape((/tmp2_arr, -tmp2_arr(tmp_size2:1:-1)/), (/2, rs_dim2/))
364  ELSE
365  tmp2_arr(:) = (/(i, i=1, rs_dim2), -rs_dim2/)
366  src_pos2(:, :) = reshape((/tmp2_arr, -tmp2_arr(tmp_size2 - 2:1:-1)/), (/2, rs_dim2/))
367  END IF
368 !---
369  srcs_coord(:, 1) = (/src_pos1(1, rs_pos(1)), src_pos2(1, rs_pos(2))/)
370  srcs_coord(:, 2) = (/src_pos1(1, rs_pos(1)), src_pos2(2, rs_pos(2))/)
371  srcs_coord(:, 3) = (/src_pos1(2, rs_pos(1)), src_pos2(1, rs_pos(2))/)
372  srcs_coord(:, 4) = (/src_pos1(2, rs_pos(1)), src_pos2(2, rs_pos(2))/)
373  CASE (neumannxz, neumannx)
374  maxn_sendrecv = 2
375  ALLOCATE (srcs_coord(2, maxn_sendrecv))
376 
377  IF (mod(rs_dim1, 2) .EQ. 0) THEN
378  tmp_size1 = rs_dim1
379  ELSE
380  tmp_size1 = rs_dim1 + 1
381  END IF
382  ALLOCATE (tmp1_arr(tmp_size1), src_pos1(2, 0:rs_dim1 - 1))
383 
384  IF (mod(rs_dim1, 2) .EQ. 0) THEN
385  tmp1_arr(:) = (/(i, i=1, rs_dim1)/)
386  src_pos1(:, :) = reshape((/tmp1_arr, -tmp1_arr(tmp_size1:1:-1)/), (/2, rs_dim1/))
387  ELSE
388  tmp1_arr(:) = (/(i, i=1, rs_dim1), -rs_dim1/)
389  src_pos1(:, :) = reshape((/tmp1_arr, -tmp1_arr(tmp_size1 - 2:1:-1)/), (/2, rs_dim1/))
390  END IF
391 !---
392  ALLOCATE (src_pos2_onesdd(0:rs_dim2 - 1))
393  src_pos2_onesdd(:) = (/(i, i=1, rs_dim2)/)
394 !---
395  srcs_coord(:, 1) = (/src_pos1(1, rs_pos(1)), src_pos2_onesdd(rs_pos(2))/)
396  srcs_coord(:, 2) = (/src_pos1(2, rs_pos(1)), src_pos2_onesdd(rs_pos(2))/)
397  CASE (neumannyz, neumanny)
398  maxn_sendrecv = 2
399  ALLOCATE (srcs_coord(2, maxn_sendrecv))
400 
401  ALLOCATE (src_pos1_onesdd(0:rs_dim1 - 1))
402  src_pos1_onesdd(:) = (/(i, i=1, rs_dim1)/)
403 !---
404  IF (mod(rs_dim2, 2) .EQ. 0) THEN
405  tmp_size2 = rs_dim2
406  ELSE
407  tmp_size2 = rs_dim2 + 1
408  END IF
409  ALLOCATE (tmp2_arr(tmp_size2), src_pos2(2, 0:rs_dim2 - 1))
410 
411  IF (mod(rs_dim2, 2) .EQ. 0) THEN
412  tmp2_arr(:) = (/(i, i=1, rs_dim2)/)
413  src_pos2(:, :) = reshape((/tmp2_arr, -tmp2_arr(tmp_size2:1:-1)/), (/2, rs_dim2/))
414  ELSE
415  tmp2_arr(:) = (/(i, i=1, rs_dim2), -rs_dim2/)
416  src_pos2(:, :) = reshape((/tmp2_arr, -tmp2_arr(tmp_size2 - 2:1:-1)/), (/2, rs_dim2/))
417  END IF
418 !---
419  srcs_coord(:, 1) = (/src_pos1_onesdd(rs_pos(1)), src_pos2(1, rs_pos(2))/)
420  srcs_coord(:, 2) = (/src_pos1_onesdd(rs_pos(1)), src_pos2(2, rs_pos(2))/)
421  CASE (neumannz)
422  maxn_sendrecv = 1
423  ALLOCATE (srcs_coord(2, maxn_sendrecv))
424  ALLOCATE (src_pos1_onesdd(0:rs_dim1 - 1))
425  ALLOCATE (src_pos2_onesdd(0:rs_dim2 - 1))
426 
427  src_pos1_onesdd(:) = (/(i, i=1, rs_dim1)/)
428 !---
429  src_pos2_onesdd(:) = (/(i, i=1, rs_dim2)/)
430 !---
431  srcs_coord(:, 1) = (/src_pos1_onesdd(rs_pos(1)), src_pos2_onesdd(rs_pos(2))/)
432  END SELECT
433 
434 ! default flipping status
435  flipg_stat = not_flipped
436 
437  DO k = 1, maxn_sendrecv
438 ! convert srcs_coord to pid
439  CALL pw_grid%para%rs_group%rank_cart(abs(srcs_coord(:, k)) - 1, srcs_expand(k))
440 ! find out the flipping status
441  IF ((srcs_coord(1, k) .GT. 0) .AND. (srcs_coord(2, k) .GT. 0)) THEN
442  flipg_stat(k) = not_flipped
443  ELSE IF ((srcs_coord(1, k) .LT. 0) .AND. (srcs_coord(2, k) .GT. 0)) THEN
444  flipg_stat(k) = ud_flipped
445  ELSE IF ((srcs_coord(1, k) .GT. 0) .AND. (srcs_coord(2, k) .LT. 0)) THEN
446  flipg_stat(k) = lr_flipped
447  ELSE
448  flipg_stat(k) = rotated
449  END IF
450  END DO
451 
452 ! let all the nodes know about each others srcs_expand list
453  ALLOCATE (srcs_expand_all(maxn_sendrecv, group_size))
454  CALL rs_group%allgather(srcs_expand, srcs_expand_all)
455 ! now scan the srcs_expand_all list and check if I am on the srcs_expand list of the other nodes
456 ! if that is the case then I am obliged to send data to those nodes (the nodes are on my dests_expand list)
457  k = 1
458  DO i = 1, group_size
459  DO j = 1, maxn_sendrecv
460  IF (srcs_expand_all(j, i) .EQ. rs_mpo) THEN
461  dests_expand(k) = i - 1
462  k = k + 1
463  END IF
464  END DO
465  END DO
466 
467 ! find srcs and dests for the reverse procedure :
468 ! initialize dests_shrink and srcs_shrink with invalid process id
469  dests_shrink = -1
470  srcs_shrink = -1
471 ! scan the flipping status of the data that I am supposed to receive
472 ! if the flipping status for a process is NOT_FLIPPED that means I will have to resend
473 ! data to that process in the reverse procedure (the process is on my dests_shrink list)
474  DO i = 1, maxn_sendrecv
475  IF (flipg_stat(i) .EQ. not_flipped) dests_shrink(i) = srcs_expand(i)
476  END DO
477 
478 ! let all the nodes know about each others dests_shrink list
479  ALLOCATE (dests_shrink_all(maxn_sendrecv, group_size))
480  CALL rs_group%allgather(dests_shrink, dests_shrink_all)
481 ! now scan the dests_shrink_all list and check if I am on the dests_shrink list of any other node
482 ! if that is the case then I'll receive data from that node (the node is on my srcs_shrink list)
483 ! note that in the shrinking procedure I will receive from only one node
484  DO i = 1, group_size
485  DO j = 1, maxn_sendrecv
486  IF (dests_shrink_all(j, i) .EQ. rs_mpo) THEN
487  srcs_shrink = i - 1
488  EXIT
489  END IF
490  END DO
491  END DO
492 
493  CALL timestop(handle)
494 
495  END SUBROUTINE set_dests_srcs_pid
496 
497 ! **************************************************************************************************
498 !> \brief expands a pw_r3d_rs_type data to an evenly symmetric pw_r3d_rs_type data that is 8 times
499 !> larger than the original one:
500 !> the even symmetry for a 1D sequence of length n is defined as:
501 !> 1 2 3 ... n-2 n-1 n --> 1 2 3 ... n-2 n-1 n n-1 n-2 ... 3 2
502 !> and is generalized to 3D by applying the same rule in all three directions
503 !>
504 !> \param neumann_directions directions in which dct should be performed
505 !> \param recv_msgs_bnds bounds of the messages to be received
506 !> \param dests_expand list of the destination processes
507 !> \param srcs_expand list of the source processes
508 !> \param flipg_stat flipping status for the received data chunks
509 !> \param bounds_shftd bounds of the original grid shifted to have g0 in the middle of the cell
510 !> \param pw_in the original plane wave data
511 !> \param pw_expanded the pw data after expansion
512 !> \par History
513 !> 07.2014 created [Hossein Bani-Hashemian]
514 !> \author Mohammad Hossein Bani-Hashemian
515 ! **************************************************************************************************
516  SUBROUTINE pw_expand(neumann_directions, recv_msgs_bnds, dests_expand, srcs_expand, &
517  flipg_stat, bounds_shftd, pw_in, pw_expanded)
518 
519  INTEGER, INTENT(IN) :: neumann_directions
520  INTEGER, DIMENSION(:, :, :), INTENT(IN), POINTER :: recv_msgs_bnds
521  INTEGER, DIMENSION(:), INTENT(IN), POINTER :: dests_expand, srcs_expand, flipg_stat
522  INTEGER, DIMENSION(2, 3), INTENT(IN) :: bounds_shftd
523  TYPE(pw_r3d_rs_type), INTENT(IN) :: pw_in
524  TYPE(pw_r3d_rs_type), INTENT(INOUT) :: pw_expanded
525 
526  CHARACTER(LEN=*), PARAMETER :: routinen = 'pw_expand'
527 
528  INTEGER :: group_size, handle, i, ind, lb1, lb1_loc, lb1_new, lb2, lb2_loc, lb2_new, lb3, &
529  lb3_loc, lb3_new, loc, maxn_sendrecv, rs_mpo, ub1, ub1_loc, ub1_new, ub2, ub2_loc, &
530  ub2_new, ub3, ub3_loc, ub3_new
531  INTEGER, ALLOCATABLE, DIMENSION(:) :: dest_hist, src_hist
532  INTEGER, ALLOCATABLE, DIMENSION(:, :, :) :: pcs_bnds
533  INTEGER, DIMENSION(2, 3) :: bounds_local_new
534  REAL(dp), ALLOCATABLE, DIMENSION(:, :, :) :: catd, catd_flipdbf, cr3d_xpndd
535  TYPE(dct_msg_type), ALLOCATABLE, DIMENSION(:) :: pcs, recv_msgs
536  TYPE(mp_comm_type) :: rs_group
537  TYPE(mp_request_type), ALLOCATABLE, DIMENSION(:) :: recv_reqs, send_reqs
538  TYPE(pw_grid_type), POINTER :: pw_grid
539 
540  CALL timeset(routinen, handle)
541 
542  pw_grid => pw_in%pw_grid
543  rs_group = pw_grid%para%rs_group
544  rs_mpo = pw_grid%para%my_pos
545  group_size = pw_grid%para%group_size
546 
547  bounds_local_new = pw_expanded%pw_grid%bounds_local
548 
549  SELECT CASE (neumann_directions)
550  CASE (neumannxyz, neumannxy)
551  maxn_sendrecv = 4
552  CASE (neumannx, neumanny, neumannxz, neumannyz)
553  maxn_sendrecv = 2
554  CASE (neumannz)
555  maxn_sendrecv = 1
556  END SELECT
557 
558  ALLOCATE (recv_reqs(maxn_sendrecv), send_reqs(maxn_sendrecv))
559  ALLOCATE (dest_hist(maxn_sendrecv), src_hist(maxn_sendrecv))
560  ALLOCATE (pcs_bnds(2, 3, maxn_sendrecv))
561 
562  ALLOCATE (pcs(maxn_sendrecv), recv_msgs(maxn_sendrecv))
563 
564  send_reqs = mp_request_null
565  recv_reqs = mp_request_null
566 
567  src_hist = -1 ! keeps the history of sources
568  dest_hist = -1 ! keeps the history of destinations
569 
570  DO i = 1, maxn_sendrecv
571 ! no need to send to myself or to the destination that I have already sent to
572  IF ((dests_expand(i) .NE. rs_mpo) .AND. .NOT. any(dest_hist .EQ. dests_expand(i))) THEN
573  CALL rs_group%isend(pw_in%array, dests_expand(i), send_reqs(i))
574  END IF
575  dest_hist(i) = dests_expand(i)
576  END DO
577 
578  DO i = 1, maxn_sendrecv
579  lb1 = recv_msgs_bnds(1, 1, i)
580  ub1 = recv_msgs_bnds(2, 1, i)
581  lb2 = recv_msgs_bnds(1, 2, i)
582  ub2 = recv_msgs_bnds(2, 2, i)
583  lb3 = recv_msgs_bnds(1, 3, i)
584  ub3 = recv_msgs_bnds(2, 3, i)
585 ! no need to receive from myself
586  IF (srcs_expand(i) .EQ. rs_mpo) THEN
587  ALLOCATE (recv_msgs(i)%msg(lb1:ub1, lb2:ub2, lb3:ub3))
588  recv_msgs(i)%msg(:, :, :) = pw_in%array
589 ! if I have already received data from the source, just use the one from the last time
590  ELSE IF (any(src_hist .EQ. srcs_expand(i))) THEN
591  loc = minloc(abs(src_hist - srcs_expand(i)), 1)
592  lb1_loc = recv_msgs_bnds(1, 1, loc)
593  ub1_loc = recv_msgs_bnds(2, 1, loc)
594  lb2_loc = recv_msgs_bnds(1, 2, loc)
595  ub2_loc = recv_msgs_bnds(2, 2, loc)
596  lb3_loc = recv_msgs_bnds(1, 3, loc)
597  ub3_loc = recv_msgs_bnds(2, 3, loc)
598  ALLOCATE (recv_msgs(i)%msg(lb1_loc:ub1_loc, lb2_loc:ub2_loc, lb3_loc:ub3_loc))
599  recv_msgs(i)%msg(:, :, :) = recv_msgs(loc)%msg
600  ELSE
601  ALLOCATE (recv_msgs(i)%msg(lb1:ub1, lb2:ub2, lb3:ub3))
602  CALL rs_group%irecv(recv_msgs(i)%msg, srcs_expand(i), recv_reqs(i))
603  CALL recv_reqs(i)%wait()
604  END IF
605  src_hist(i) = srcs_expand(i)
606  END DO
607 ! cleanup mpi_request to prevent memory leak
608  CALL mp_waitall(send_reqs(:))
609 
610 ! flip the received data according on the flipping status
611  DO i = 1, maxn_sendrecv
612  SELECT CASE (flipg_stat(i))
613  CASE (not_flipped)
614  lb1 = recv_msgs_bnds(1, 1, i)
615  ub1 = recv_msgs_bnds(2, 1, i)
616  lb2 = recv_msgs_bnds(1, 2, i)
617  ub2 = recv_msgs_bnds(2, 2, i)
618  lb3 = recv_msgs_bnds(1, 3, i)
619  ub3 = recv_msgs_bnds(2, 3, i)
620  ALLOCATE (pcs(i)%msg(lb1:ub1, lb2:ub2, lb3:ub3))
621  pcs(i)%msg(:, :, :) = recv_msgs(i)%msg
622  CASE (ud_flipped)
623  CALL flipud(recv_msgs(i)%msg, pcs(i)%msg, bounds_shftd)
624  CASE (lr_flipped)
625  CALL fliplr(recv_msgs(i)%msg, pcs(i)%msg, bounds_shftd)
626  CASE (bf_flipped)
627  CALL flipbf(recv_msgs(i)%msg, pcs(i)%msg, bounds_shftd)
628  CASE (rotated)
629  CALL rot180(recv_msgs(i)%msg, pcs(i)%msg, bounds_shftd)
630  END SELECT
631  END DO
632 ! concatenate the received (flipped) data store the result as catd
633 ! need the bounds of the four pieces for concatenation
634  DO i = 1, maxn_sendrecv
635  pcs_bnds(1, 1, i) = lbound(pcs(i)%msg, 1)
636  pcs_bnds(2, 1, i) = ubound(pcs(i)%msg, 1)
637  pcs_bnds(1, 2, i) = lbound(pcs(i)%msg, 2)
638  pcs_bnds(2, 2, i) = ubound(pcs(i)%msg, 2)
639  pcs_bnds(1, 3, i) = lbound(pcs(i)%msg, 3)
640  pcs_bnds(2, 3, i) = ubound(pcs(i)%msg, 3)
641  END DO
642 
643  lb1_new = bounds_local_new(1, 1); ub1_new = bounds_local_new(2, 1)
644  lb2_new = bounds_local_new(1, 2); ub2_new = bounds_local_new(2, 2)
645  lb3_new = bounds_local_new(1, 3); ub3_new = bounds_local_new(2, 3)
646 
647  SELECT CASE (neumann_directions)
648  CASE (neumannxyz, neumannxz, neumannyz, neumannz)
649  ind = int(0.5*(ub3_new + lb3_new + 1))
650  ALLOCATE (catd(lb1_new:ub1_new, lb2_new:ub2_new, lb3_new:ind - 1))
651  CASE (neumannxy, neumannx, neumanny)
652  ALLOCATE (catd(lb1_new:ub1_new, lb2_new:ub2_new, lb3_new:ub3_new))
653  END SELECT
654 
655  DO i = 1, maxn_sendrecv
656  catd(pcs_bnds(1, 1, i):pcs_bnds(2, 1, i), &
657  pcs_bnds(1, 2, i):pcs_bnds(2, 2, i), &
658  pcs_bnds(1, 3, i):pcs_bnds(2, 3, i)) = pcs(i)%msg
659  END DO
660 
661 ! flip catd from back to front
662  CALL flipbf(catd, catd_flipdbf, bounds_shftd)
663 ! concatenate catd and catd_flipdbf to get cr3d_xpndd
664  ALLOCATE (cr3d_xpndd(lb1_new:ub1_new, lb2_new:ub2_new, lb3_new:ub3_new))
665  SELECT CASE (neumann_directions)
666  CASE (neumannxyz, neumannxz, neumannyz, neumannz)
667  cr3d_xpndd(:, :, lb3_new:ind - 1) = catd
668  cr3d_xpndd(:, :, ind:ub3_new) = catd_flipdbf
669  CASE (neumannxy, neumannx, neumanny)
670  cr3d_xpndd(:, :, :) = catd
671  END SELECT
672 
673  pw_expanded%array = cr3d_xpndd
674 
675  DO i = 1, maxn_sendrecv
676  DEALLOCATE (pcs(i)%msg)
677  DEALLOCATE (recv_msgs(i)%msg)
678  END DO
679  DEALLOCATE (pcs, recv_msgs)
680  DEALLOCATE (catd, catd_flipdbf, cr3d_xpndd)
681 
682  CALL timestop(handle)
683 
684  END SUBROUTINE pw_expand
685 
686 ! **************************************************************************************************
687 !> \brief shrinks an evenly symmetric pw_r3d_rs_type data to a pw_r3d_rs_type data that is 8
688 !> times smaller (the reverse procedure of pw_expand).
689 !>
690 !> \param neumann_directions directions in which dct should be performed
691 !> \param dests_shrink list of the destination processes
692 !> \param srcs_shrink list of the source processes
693 !> \param bounds_local_shftd local bounds of the original grid after shifting
694 !> \param pw_in the original plane wave data
695 !> \param pw_shrinked the shrunk plane wave data
696 !> \par History
697 !> 07.2014 created [Hossein Bani-Hashemian]
698 !> \author Mohammad Hossein Bani-Hashemian
699 ! **************************************************************************************************
700  SUBROUTINE pw_shrink(neumann_directions, dests_shrink, srcs_shrink, bounds_local_shftd, &
701  pw_in, pw_shrinked)
702 
703  INTEGER, INTENT(IN) :: neumann_directions
704  INTEGER, DIMENSION(:), INTENT(IN), POINTER :: dests_shrink
705  INTEGER, INTENT(INOUT) :: srcs_shrink
706  INTEGER, DIMENSION(2, 3), INTENT(IN) :: bounds_local_shftd
707  TYPE(pw_r3d_rs_type), INTENT(IN) :: pw_in
708  TYPE(pw_r3d_rs_type), INTENT(INOUT) :: pw_shrinked
709 
710  CHARACTER(LEN=*), PARAMETER :: routinen = 'pw_shrink'
711 
712  INTEGER :: group_size, handle, i, lb1_orig, lb1_xpnd, lb2_orig, lb2_xpnd, lb3_orig, &
713  lb3_xpnd, maxn_sendrecv, rs_mpo, send_lb1, send_lb2, send_lb3, send_ub1, send_ub2, &
714  send_ub3, tag, ub1_orig, ub1_xpnd, ub2_orig, ub2_xpnd, ub3_orig, ub3_xpnd
715  INTEGER, ALLOCATABLE, DIMENSION(:, :, :) :: bounds_local_all
716  INTEGER, DIMENSION(2, 3) :: bounds_local_xpnd
717  REAL(dp), ALLOCATABLE, DIMENSION(:, :, :) :: cr3d, send_crmsg
718  TYPE(mp_comm_type) :: rs_group
719  TYPE(pw_grid_type), POINTER :: pw_grid_orig
720 
721  CALL timeset(routinen, handle)
722 
723  pw_grid_orig => pw_shrinked%pw_grid
724  rs_group = pw_grid_orig%para%rs_group
725  rs_mpo = pw_grid_orig%para%my_pos
726  group_size = pw_grid_orig%para%group_size
727  bounds_local_xpnd = pw_in%pw_grid%bounds_local
728  tag = 1
729 
730  SELECT CASE (neumann_directions)
731  CASE (neumannxyz, neumannxy)
732  maxn_sendrecv = 4
733  CASE (neumannx, neumanny, neumannxz, neumannyz)
734  maxn_sendrecv = 2
735  CASE (neumannz)
736  maxn_sendrecv = 1
737  END SELECT
738 
739 ! cosine transform is a real transform. The cosine transform of a 3D data must be real and 3D.
740  lb1_xpnd = bounds_local_xpnd(1, 1)
741  ub1_xpnd = bounds_local_xpnd(2, 1)
742  lb2_xpnd = bounds_local_xpnd(1, 2)
743  ub2_xpnd = bounds_local_xpnd(2, 2)
744  lb3_xpnd = bounds_local_xpnd(1, 3)
745  ub3_xpnd = bounds_local_xpnd(2, 3)
746  ALLOCATE (cr3d(lb1_xpnd:ub1_xpnd, lb2_xpnd:ub2_xpnd, lb3_xpnd:ub3_xpnd))
747  cr3d(:, :, :) = pw_in%array
748 
749 ! let all the nodes know about each others shifted local bounds
750  ALLOCATE (bounds_local_all(2, 3, group_size))
751  CALL rs_group%allgather(bounds_local_shftd, bounds_local_all)
752 
753  DO i = 1, maxn_sendrecv
754 ! no need to send to myself or to an invalid destination (pid = -1)
755  IF ((dests_shrink(i) .NE. rs_mpo) .AND. (dests_shrink(i) .NE. -1)) THEN
756  send_lb1 = bounds_local_all(1, 1, dests_shrink(i) + 1)
757  send_ub1 = bounds_local_all(2, 1, dests_shrink(i) + 1)
758  send_lb2 = bounds_local_all(1, 2, dests_shrink(i) + 1)
759  send_ub2 = bounds_local_all(2, 2, dests_shrink(i) + 1)
760  send_lb3 = bounds_local_all(1, 3, dests_shrink(i) + 1)
761  send_ub3 = bounds_local_all(2, 3, dests_shrink(i) + 1)
762 
763  ALLOCATE (send_crmsg(send_lb1:send_ub1, send_lb2:send_ub2, send_lb3:send_ub3))
764  send_crmsg(:, :, :) = cr3d(send_lb1:send_ub1, send_lb2:send_ub2, send_lb3:send_ub3)
765  CALL rs_group%send(send_crmsg, dests_shrink(i), tag)
766  DEALLOCATE (send_crmsg)
767  END IF
768  END DO
769 
770  lb1_orig = bounds_local_shftd(1, 1)
771  ub1_orig = bounds_local_shftd(2, 1)
772  lb2_orig = bounds_local_shftd(1, 2)
773  ub2_orig = bounds_local_shftd(2, 2)
774  lb3_orig = bounds_local_shftd(1, 3)
775  ub3_orig = bounds_local_shftd(2, 3)
776 
777 ! no need to receive from myself
778  IF (srcs_shrink .EQ. rs_mpo) THEN
779  pw_shrinked%array = cr3d(lb1_orig:ub1_orig, lb2_orig:ub2_orig, lb3_orig:ub3_orig)
780  ELSE IF (srcs_shrink .EQ. -1) THEN
781 ! the source is invalid ... do nothing
782  ELSE
783  CALL rs_group%recv(pw_shrinked%array, srcs_shrink, tag)
784  END IF
785 
786  DEALLOCATE (bounds_local_all)
787  DEALLOCATE (cr3d)
788  CALL timestop(handle)
789 
790  END SUBROUTINE pw_shrink
791 
792 ! **************************************************************************************************
793 !> \brief flips a 3d (real dp) array up to down (the way needed to expand data
794 !> as explained in the description of the afore-defined subroutines)
795 !> \param cr3d_in input array
796 !> \param cr3d_out output array
797 !> \param bounds global lower and upper bounds
798 !> \par History
799 !> 07.2014 created [Hossein Bani-Hashemian]
800 !> \author Mohammad Hossein Bani-Hashemian
801 ! **************************************************************************************************
802  SUBROUTINE flipud(cr3d_in, cr3d_out, bounds)
803 
804  REAL(dp), ALLOCATABLE, DIMENSION(:, :, :), &
805  INTENT(IN) :: cr3d_in
806  REAL(dp), ALLOCATABLE, DIMENSION(:, :, :), &
807  INTENT(INOUT) :: cr3d_out
808  INTEGER, DIMENSION(2, 3), INTENT(IN) :: bounds
809 
810  CHARACTER(LEN=*), PARAMETER :: routinen = 'flipud'
811 
812  INTEGER :: handle, i, lb1, lb1_glbl, lb1_new, lb2, lb2_glbl, lb2_new, lb3, lb3_glbl, &
813  lb3_new, ub1, ub1_glbl, ub1_new, ub2, ub2_glbl, ub2_new, ub3, ub3_glbl, ub3_new
814  INTEGER, ALLOCATABLE, DIMENSION(:) :: indx
815  INTEGER, DIMENSION(2, 3) :: bndsl, bndsl_new
816 
817  CALL timeset(routinen, handle)
818 
819  lb1 = lbound(cr3d_in, 1); ub1 = ubound(cr3d_in, 1)
820  lb2 = lbound(cr3d_in, 2); ub2 = ubound(cr3d_in, 2)
821  lb3 = lbound(cr3d_in, 3); ub3 = ubound(cr3d_in, 3)
822 
823  lb1_glbl = bounds(1, 1); ub1_glbl = bounds(2, 1)
824  lb2_glbl = bounds(1, 2); ub2_glbl = bounds(2, 2)
825  lb3_glbl = bounds(1, 3); ub3_glbl = bounds(2, 3)
826 
827  bndsl = reshape((/lb1, ub1, lb2, ub2, lb3, ub3/), (/2, 3/))
828  bndsl_new = flipud_bounds_local(bndsl, bounds)
829 
830  lb1_new = bndsl_new(1, 1); ub1_new = bndsl_new(2, 1)
831  lb2_new = bndsl_new(1, 2); ub2_new = bndsl_new(2, 2)
832  lb3_new = bndsl_new(1, 3); ub3_new = bndsl_new(2, 3)
833 
834  ALLOCATE (cr3d_out(lb1_new:ub1_new, lb2_new:ub2_new, lb3_new:ub3_new))
835  cr3d_out = 0.0_dp
836 
837 ! set the data at the missing grid points (in a periodic grid) equal to the data at
838 ! the last existing grid points
839  ALLOCATE (indx(ub1_new - lb1_new + 1))
840  indx(:) = (/(i, i=2*(ub1_glbl + 1) - lb1_new, 2*(ub1_glbl + 1) - ub1_new, -1)/)
841  IF (lb1_new .EQ. ub1_glbl + 1) indx(1) = indx(2)
842  cr3d_out(lb1_new:ub1_new, lb2_new:ub2_new, lb3_new:ub3_new) = cr3d_in(indx, :, :)
843 
844  CALL timestop(handle)
845 
846  END SUBROUTINE flipud
847 
848 ! **************************************************************************************************
849 !> \brief flips a 3d (real dp) array left to right (the way needed to expand data
850 !> as explained in the description of the afore-defined subroutines)
851 !> \param cr3d_in input array
852 !> \param cr3d_out output array
853 !> \param bounds global lower and upper bounds
854 !> \par History
855 !> 07.2014 created [Hossein Bani-Hashemian]
856 !> \author Mohammad Hossein Bani-Hashemian
857 ! **************************************************************************************************
858  SUBROUTINE fliplr(cr3d_in, cr3d_out, bounds)
859 
860  REAL(dp), ALLOCATABLE, DIMENSION(:, :, :), &
861  INTENT(IN) :: cr3d_in
862  REAL(dp), ALLOCATABLE, DIMENSION(:, :, :), &
863  INTENT(INOUT) :: cr3d_out
864  INTEGER, DIMENSION(2, 3), INTENT(IN) :: bounds
865 
866  CHARACTER(LEN=*), PARAMETER :: routinen = 'fliplr'
867 
868  INTEGER :: handle, i, lb1, lb1_glbl, lb1_new, lb2, lb2_glbl, lb2_new, lb3, lb3_glbl, &
869  lb3_new, ub1, ub1_glbl, ub1_new, ub2, ub2_glbl, ub2_new, ub3, ub3_glbl, ub3_new
870  INTEGER, ALLOCATABLE, DIMENSION(:) :: indy
871  INTEGER, DIMENSION(2, 3) :: bndsl, bndsl_new
872 
873  CALL timeset(routinen, handle)
874 
875  lb1 = lbound(cr3d_in, 1); ub1 = ubound(cr3d_in, 1)
876  lb2 = lbound(cr3d_in, 2); ub2 = ubound(cr3d_in, 2)
877  lb3 = lbound(cr3d_in, 3); ub3 = ubound(cr3d_in, 3)
878 
879  lb1_glbl = bounds(1, 1); ub1_glbl = bounds(2, 1)
880  lb2_glbl = bounds(1, 2); ub2_glbl = bounds(2, 2)
881  lb3_glbl = bounds(1, 3); ub3_glbl = bounds(2, 3)
882 
883  bndsl = reshape((/lb1, ub1, lb2, ub2, lb3, ub3/), (/2, 3/))
884  bndsl_new = fliplr_bounds_local(bndsl, bounds)
885 
886  lb1_new = bndsl_new(1, 1); ub1_new = bndsl_new(2, 1)
887  lb2_new = bndsl_new(1, 2); ub2_new = bndsl_new(2, 2)
888  lb3_new = bndsl_new(1, 3); ub3_new = bndsl_new(2, 3)
889 
890  ALLOCATE (cr3d_out(lb1_new:ub1_new, lb2_new:ub2_new, lb3_new:ub3_new))
891  cr3d_out = 0.0_dp
892 
893 ! set the data at the missing grid points (in a periodic grid) equal to the data at
894 ! the last existing grid points
895  ALLOCATE (indy(ub2_new - lb2_new + 1))
896  indy(:) = (/(i, i=2*(ub2_glbl + 1) - lb2_new, 2*(ub2_glbl + 1) - ub2_new, -1)/)
897  IF (lb2_new .EQ. ub2_glbl + 1) indy(1) = indy(2)
898  cr3d_out(lb1_new:ub1_new, lb2_new:ub2_new, lb3_new:ub3_new) = cr3d_in(:, indy, :)
899 
900  CALL timestop(handle)
901 
902  END SUBROUTINE fliplr
903 
904 ! **************************************************************************************************
905 !> \brief flips a 3d (real dp) array back to front (the way needed to expand data
906 !> as explained in the description of the afore-defined subroutines)
907 !> \param cr3d_in input array
908 !> \param cr3d_out output array
909 !> \param bounds global lower and upper bounds
910 !> \par History
911 !> 07.2014 created [Hossein Bani-Hashemian]
912 !> \author Mohammad Hossein Bani-Hashemian
913 ! **************************************************************************************************
914  SUBROUTINE flipbf(cr3d_in, cr3d_out, bounds)
915 
916  REAL(dp), ALLOCATABLE, DIMENSION(:, :, :), &
917  INTENT(IN) :: cr3d_in
918  REAL(dp), ALLOCATABLE, DIMENSION(:, :, :), &
919  INTENT(INOUT) :: cr3d_out
920  INTEGER, DIMENSION(2, 3), INTENT(IN) :: bounds
921 
922  CHARACTER(LEN=*), PARAMETER :: routinen = 'flipbf'
923 
924  INTEGER :: handle, i, lb1, lb1_glbl, lb1_new, lb2, lb2_glbl, lb2_new, lb3, lb3_glbl, &
925  lb3_new, ub1, ub1_glbl, ub1_new, ub2, ub2_glbl, ub2_new, ub3, ub3_glbl, ub3_new
926  INTEGER, ALLOCATABLE, DIMENSION(:) :: indz
927  INTEGER, DIMENSION(2, 3) :: bndsl, bndsl_new
928 
929  CALL timeset(routinen, handle)
930 
931  lb1 = lbound(cr3d_in, 1); ub1 = ubound(cr3d_in, 1)
932  lb2 = lbound(cr3d_in, 2); ub2 = ubound(cr3d_in, 2)
933  lb3 = lbound(cr3d_in, 3); ub3 = ubound(cr3d_in, 3)
934 
935  lb1_glbl = bounds(1, 1); ub1_glbl = bounds(2, 1)
936  lb2_glbl = bounds(1, 2); ub2_glbl = bounds(2, 2)
937  lb3_glbl = bounds(1, 3); ub3_glbl = bounds(2, 3)
938 
939  bndsl = reshape((/lb1, ub1, lb2, ub2, lb3, ub3/), (/2, 3/))
940  bndsl_new = flipbf_bounds_local(bndsl, bounds)
941 
942  lb1_new = bndsl_new(1, 1); ub1_new = bndsl_new(2, 1)
943  lb2_new = bndsl_new(1, 2); ub2_new = bndsl_new(2, 2)
944  lb3_new = bndsl_new(1, 3); ub3_new = bndsl_new(2, 3)
945 
946  ALLOCATE (cr3d_out(lb1_new:ub1_new, lb2_new:ub2_new, lb3_new:ub3_new))
947  cr3d_out = 0.0_dp
948 
949 ! set the data at the missing grid points (in a periodic grid) equal to the data at
950 ! the last existing grid points
951  ALLOCATE (indz(ub3_new - lb3_new + 1))
952  indz(:) = (/(i, i=2*(ub3_glbl + 1) - lb3_new, 2*(ub3_glbl + 1) - ub3_new, -1)/)
953  IF (lb3_new .EQ. ub3_glbl + 1) indz(1) = indz(2)
954  cr3d_out(lb1_new:ub1_new, lb2_new:ub2_new, lb3_new:ub3_new) = cr3d_in(:, :, indz)
955 
956  CALL timestop(handle)
957 
958  END SUBROUTINE flipbf
959 
960 ! **************************************************************************************************
961 !> \brief rotates a 3d (real dp) array by 180 degrees (the way needed to expand data
962 !> as explained in the description of the afore-defined subroutines)
963 !> \param cr3d_in input array
964 !> \param cr3d_out output array
965 !> \param bounds global lower and upper bounds
966 !> \par History
967 !> 07.2014 created [Hossein Bani-Hashemian]
968 !> \author Mohammad Hossein Bani-Hashemian
969 ! **************************************************************************************************
970  SUBROUTINE rot180(cr3d_in, cr3d_out, bounds)
971 
972  REAL(dp), ALLOCATABLE, DIMENSION(:, :, :), &
973  INTENT(IN) :: cr3d_in
974  REAL(dp), ALLOCATABLE, DIMENSION(:, :, :), &
975  INTENT(INOUT) :: cr3d_out
976  INTEGER, DIMENSION(2, 3), INTENT(IN) :: bounds
977 
978  CHARACTER(LEN=*), PARAMETER :: routinen = 'rot180'
979 
980  INTEGER :: handle, i, lb1, lb1_glbl, lb1_new, lb2, lb2_glbl, lb2_new, lb3, lb3_glbl, &
981  lb3_new, ub1, ub1_glbl, ub1_new, ub2, ub2_glbl, ub2_new, ub3, ub3_glbl, ub3_new
982  INTEGER, ALLOCATABLE, DIMENSION(:) :: indx, indy
983  INTEGER, DIMENSION(2, 3) :: bndsl, bndsl_new
984 
985  CALL timeset(routinen, handle)
986 
987  lb1 = lbound(cr3d_in, 1); ub1 = ubound(cr3d_in, 1)
988  lb2 = lbound(cr3d_in, 2); ub2 = ubound(cr3d_in, 2)
989  lb3 = lbound(cr3d_in, 3); ub3 = ubound(cr3d_in, 3)
990 
991  lb1_glbl = bounds(1, 1); ub1_glbl = bounds(2, 1)
992  lb2_glbl = bounds(1, 2); ub2_glbl = bounds(2, 2)
993  lb3_glbl = bounds(1, 3); ub3_glbl = bounds(2, 3)
994 
995  bndsl = reshape((/lb1, ub1, lb2, ub2, lb3, ub3/), (/2, 3/))
996  bndsl_new = rot180_bounds_local(bndsl, bounds)
997 
998  lb1_new = bndsl_new(1, 1); ub1_new = bndsl_new(2, 1)
999  lb2_new = bndsl_new(1, 2); ub2_new = bndsl_new(2, 2)
1000  lb3_new = bndsl_new(1, 3); ub3_new = bndsl_new(2, 3)
1001 
1002  ALLOCATE (cr3d_out(lb1_new:ub1_new, lb2_new:ub2_new, lb3_new:ub3_new))
1003  cr3d_out = 0.0_dp
1004 
1005 ! set the data at the missing grid points (in a periodic grid) equal to the data at
1006 ! the last existing grid points
1007  ALLOCATE (indx(ub1_new - lb1_new + 1), indy(ub2_new - lb2_new + 1))
1008  indx(:) = (/(i, i=2*(ub1_glbl + 1) - lb1_new, 2*(ub1_glbl + 1) - ub1_new, -1)/)
1009  indy(:) = (/(i, i=2*(ub2_glbl + 1) - lb2_new, 2*(ub2_glbl + 1) - ub2_new, -1)/)
1010  IF (lb1_new .EQ. ub1_glbl + 1) indx(1) = indx(2)
1011  IF (lb2_new .EQ. ub2_glbl + 1) indy(1) = indy(2)
1012  cr3d_out(lb1_new:ub1_new, lb2_new:ub2_new, lb3_new:ub3_new) = cr3d_in(indx, indy, :)
1013 
1014  CALL timestop(handle)
1015 
1016  END SUBROUTINE rot180
1017 
1018 ! **************************************************************************************************
1019 !> \brief calculates the global and local bounds of the expanded data
1020 !> \param pw_grid original plane-wave grid
1021 !> \param neumann_directions directions in which dct should be performed
1022 !> \param srcs_expand list of the source processes (pw_expand)
1023 !> \param flipg_stat flipping status for the received data chunks (pw_expand)
1024 !> \param bounds_shftd bounds of the original grid shifted to have g0 in the middle of the cell
1025 !> \param bounds_local_shftd local bounds of the original grid after shifting
1026 !> \param recv_msgs_bnds bounds of the messages to be received (pw_expand)
1027 !> \param bounds_new new global lower and upper bounds
1028 !> \param bounds_local_new new local lower and upper bounds
1029 !> \par History
1030 !> 07.2014 created [Hossein Bani-Hashemian]
1031 !> \author Mohammad Hossein Bani-Hashemian
1032 ! **************************************************************************************************
1033  SUBROUTINE expansion_bounds(pw_grid, neumann_directions, srcs_expand, flipg_stat, &
1034  bounds_shftd, bounds_local_shftd, &
1035  recv_msgs_bnds, bounds_new, bounds_local_new)
1036 
1037  TYPE(pw_grid_type), INTENT(IN), POINTER :: pw_grid
1038  INTEGER, INTENT(IN) :: neumann_directions
1039  INTEGER, DIMENSION(:), INTENT(IN), POINTER :: srcs_expand, flipg_stat
1040  INTEGER, DIMENSION(2, 3), INTENT(OUT) :: bounds_shftd, bounds_local_shftd
1041  INTEGER, DIMENSION(:, :, :), INTENT(INOUT), &
1042  POINTER :: recv_msgs_bnds
1043  INTEGER, DIMENSION(2, 3), INTENT(OUT) :: bounds_new, bounds_local_new
1044 
1045  CHARACTER(LEN=*), PARAMETER :: routinen = 'expansion_bounds'
1046 
1047  INTEGER :: group_size, handle, i, lb1_new, lb2_new, &
1048  lb3_new, loc, maxn_sendrecv, rs_mpo, &
1049  ub1_new, ub2_new, ub3_new
1050  INTEGER, ALLOCATABLE, DIMENSION(:) :: src_hist
1051  INTEGER, ALLOCATABLE, DIMENSION(:, :, :) :: bounds_local_all, bounds_local_new_all, &
1052  pcs_bnds
1053  INTEGER, DIMENSION(2, 3) :: bounds, bounds_local
1054  INTEGER, DIMENSION(3) :: npts_new, shf_yesno, shift
1055  TYPE(mp_comm_type) :: rs_group
1056 
1057  CALL timeset(routinen, handle)
1058 
1059  rs_group = pw_grid%para%rs_group
1060  rs_mpo = pw_grid%para%my_pos
1061  group_size = pw_grid%para%group_size
1062  bounds = pw_grid%bounds
1063  bounds_local = pw_grid%bounds_local
1064 
1065  SELECT CASE (neumann_directions)
1066  CASE (neumannxyz)
1067  maxn_sendrecv = 4
1068  shf_yesno = (/1, 1, 1/)
1069  CASE (neumannxy)
1070  maxn_sendrecv = 4
1071  shf_yesno = (/1, 1, 0/)
1072  CASE (neumannxz)
1073  maxn_sendrecv = 2
1074  shf_yesno = (/1, 0, 1/)
1075  CASE (neumannyz)
1076  maxn_sendrecv = 2
1077  shf_yesno = (/0, 1, 1/)
1078  CASE (neumannx)
1079  maxn_sendrecv = 2
1080  shf_yesno = (/1, 0, 0/)
1081  CASE (neumanny)
1082  maxn_sendrecv = 2
1083  shf_yesno = (/0, 1, 0/)
1084  CASE (neumannz)
1085  maxn_sendrecv = 1
1086  shf_yesno = (/0, 0, 1/)
1087  CASE DEFAULT
1088  cpabort("Unknown neumann direction")
1089  shf_yesno = (/0, 0, 0/)
1090  END SELECT
1091 
1092  ALLOCATE (pcs_bnds(2, 3, maxn_sendrecv))
1093  ALLOCATE (src_hist(maxn_sendrecv))
1094 
1095  ! Note that this is not easily FFT-able ... needed anyway, so link in FFTW.
1096  npts_new = 2*pw_grid%npts
1097  shift = -npts_new/2
1098  shift = shift - bounds(1, :)
1099  bounds_shftd(:, 1) = bounds(:, 1) + shf_yesno(1)*shift(1)
1100  bounds_shftd(:, 2) = bounds(:, 2) + shf_yesno(2)*shift(2)
1101  bounds_shftd(:, 3) = bounds(:, 3) + shf_yesno(3)*shift(3)
1102  bounds_local_shftd(:, 1) = bounds_local(:, 1) + shf_yesno(1)*shift(1)
1103  bounds_local_shftd(:, 2) = bounds_local(:, 2) + shf_yesno(2)*shift(2)
1104  bounds_local_shftd(:, 3) = bounds_local(:, 3) + shf_yesno(3)*shift(3)
1105 
1106 ! let all the nodes know about each others local shifted bounds
1107  ALLOCATE (bounds_local_all(2, 3, group_size))
1108  CALL rs_group%allgather(bounds_local_shftd, bounds_local_all)
1109 
1110  src_hist = -1 ! keeps the history of sources
1111 
1112  DO i = 1, maxn_sendrecv
1113 ! no need to receive from myself
1114  IF (srcs_expand(i) .EQ. rs_mpo) THEN
1115  recv_msgs_bnds(1, 1, i) = bounds_local_shftd(1, 1)
1116  recv_msgs_bnds(2, 1, i) = bounds_local_shftd(2, 1)
1117  recv_msgs_bnds(1, 2, i) = bounds_local_shftd(1, 2)
1118  recv_msgs_bnds(2, 2, i) = bounds_local_shftd(2, 2)
1119  recv_msgs_bnds(1, 3, i) = bounds_local_shftd(1, 3)
1120  recv_msgs_bnds(2, 3, i) = bounds_local_shftd(2, 3)
1121 ! if I have already received data from the source, just use the one from the last time
1122  ELSE IF (any(src_hist .EQ. srcs_expand(i))) THEN
1123  loc = minloc(abs(src_hist - srcs_expand(i)), 1)
1124  recv_msgs_bnds(1, 1, i) = bounds_local_all(1, 1, srcs_expand(loc) + 1)
1125  recv_msgs_bnds(2, 1, i) = bounds_local_all(2, 1, srcs_expand(loc) + 1)
1126  recv_msgs_bnds(1, 2, i) = bounds_local_all(1, 2, srcs_expand(loc) + 1)
1127  recv_msgs_bnds(2, 2, i) = bounds_local_all(2, 2, srcs_expand(loc) + 1)
1128  recv_msgs_bnds(1, 3, i) = bounds_local_all(1, 3, srcs_expand(loc) + 1)
1129  recv_msgs_bnds(2, 3, i) = bounds_local_all(2, 3, srcs_expand(loc) + 1)
1130  ELSE
1131  recv_msgs_bnds(1, 1, i) = bounds_local_all(1, 1, srcs_expand(i) + 1)
1132  recv_msgs_bnds(2, 1, i) = bounds_local_all(2, 1, srcs_expand(i) + 1)
1133  recv_msgs_bnds(1, 2, i) = bounds_local_all(1, 2, srcs_expand(i) + 1)
1134  recv_msgs_bnds(2, 2, i) = bounds_local_all(2, 2, srcs_expand(i) + 1)
1135  recv_msgs_bnds(1, 3, i) = bounds_local_all(1, 3, srcs_expand(i) + 1)
1136  recv_msgs_bnds(2, 3, i) = bounds_local_all(2, 3, srcs_expand(i) + 1)
1137  END IF
1138  src_hist(i) = srcs_expand(i)
1139  END DO
1140 
1141 ! flip the received data based on the flipping status
1142  DO i = 1, maxn_sendrecv
1143  SELECT CASE (flipg_stat(i))
1144  CASE (not_flipped)
1145  pcs_bnds(:, :, i) = recv_msgs_bnds(:, :, i)
1146  CASE (ud_flipped)
1147  pcs_bnds(:, :, i) = flipud_bounds_local(recv_msgs_bnds(:, :, i), bounds_shftd)
1148  CASE (lr_flipped)
1149  pcs_bnds(:, :, i) = fliplr_bounds_local(recv_msgs_bnds(:, :, i), bounds_shftd)
1150  CASE (bf_flipped)
1151  pcs_bnds(:, :, i) = flipbf_bounds_local(recv_msgs_bnds(:, :, i), bounds_shftd)
1152  CASE (rotated)
1153  pcs_bnds(:, :, i) = rot180_bounds_local(recv_msgs_bnds(:, :, i), bounds_shftd)
1154  END SELECT
1155  END DO
1156 
1157  lb1_new = minval(pcs_bnds(1, 1, :)); ub1_new = maxval(pcs_bnds(2, 1, :))
1158  lb2_new = minval(pcs_bnds(1, 2, :)); ub2_new = maxval(pcs_bnds(2, 2, :))
1159  lb3_new = minval(pcs_bnds(1, 3, :)); ub3_new = maxval(pcs_bnds(2, 3, :))
1160 
1161 ! calculate the new local and global bounds
1162  bounds_local_new(1, 1) = minval(pcs_bnds(1, 1, :))
1163  bounds_local_new(2, 1) = maxval(pcs_bnds(2, 1, :))
1164  bounds_local_new(1, 2) = minval(pcs_bnds(1, 2, :))
1165  bounds_local_new(2, 2) = maxval(pcs_bnds(2, 2, :))
1166  bounds_local_new(1, 3) = minval(pcs_bnds(1, 3, :))
1167  SELECT CASE (neumann_directions)
1168  CASE (neumannxyz, neumannxz, neumannyz, neumannz)
1169  bounds_local_new(2, 3) = 2*(maxval(pcs_bnds(2, 3, :)) + 1) - bounds_local_new(1, 3) - 1
1170  CASE (neumannxy, neumannx, neumanny)
1171  bounds_local_new(2, 3) = maxval(pcs_bnds(2, 3, :))
1172  END SELECT
1173 
1174  ALLOCATE (bounds_local_new_all(2, 3, group_size))
1175  CALL rs_group%allgather(bounds_local_new, bounds_local_new_all)
1176  bounds_new(1, 1) = minval(bounds_local_new_all(1, 1, :))
1177  bounds_new(2, 1) = maxval(bounds_local_new_all(2, 1, :))
1178  bounds_new(1, 2) = minval(bounds_local_new_all(1, 2, :))
1179  bounds_new(2, 2) = maxval(bounds_local_new_all(2, 2, :))
1180  bounds_new(1, 3) = minval(bounds_local_new_all(1, 3, :))
1181  bounds_new(2, 3) = maxval(bounds_local_new_all(2, 3, :))
1182 
1183  DEALLOCATE (bounds_local_all, bounds_local_new_all)
1184 
1185  CALL timestop(handle)
1186 
1187  END SUBROUTINE expansion_bounds
1188 
1189 ! **************************************************************************************************
1190 !> \brief precalculates the local bounds of a 3d array after applying flipud
1191 !> \param bndsl_in current local lower and upper bounds
1192 !> \param bounds global lower and upper bounds
1193 !> \return new local lower and upper bounds
1194 !> \par History
1195 !> 07.2014 created [Hossein Bani-Hashemian]
1196 !> \author Mohammad Hossein Bani-Hashemian
1197 ! **************************************************************************************************
1198  PURE FUNCTION flipud_bounds_local(bndsl_in, bounds) RESULT(bndsl_out)
1199 
1200  INTEGER, DIMENSION(2, 3), INTENT(IN) :: bndsl_in, bounds
1201  INTEGER, DIMENSION(2, 3) :: bndsl_out
1202 
1203  bndsl_out(1, 1) = 2*(bounds(2, 1) + 1) - bndsl_in(2, 1)
1204  bndsl_out(2, 1) = 2*(bounds(2, 1) + 1) - bndsl_in(1, 1)
1205  IF (bndsl_out(1, 1) .EQ. bounds(2, 1) + 2) bndsl_out(1, 1) = bndsl_out(1, 1) - 1
1206  IF (bndsl_out(2, 1) .EQ. 2*(bounds(2, 1) + 1) - bounds(1, 1)) bndsl_out(2, 1) = bndsl_out(2, 1) - 1
1207 
1208  bndsl_out(1, 2) = bndsl_in(1, 2)
1209  bndsl_out(2, 2) = bndsl_in(2, 2)
1210 
1211  bndsl_out(1, 3) = bndsl_in(1, 3)
1212  bndsl_out(2, 3) = bndsl_in(2, 3)
1213 
1214  END FUNCTION flipud_bounds_local
1215 
1216 ! **************************************************************************************************
1217 !> \brief precalculates the local bounds of a 3d array after applying fliplr
1218 !> \param bndsl_in current local lower and upper bounds
1219 !> \param bounds global lower and upper bounds
1220 !> \return new local lower and upper bounds
1221 !> \par History
1222 !> 07.2014 created [Hossein Bani-Hashemian]
1223 !> \author Mohammad Hossein Bani-Hashemian
1224 ! **************************************************************************************************
1225  PURE FUNCTION fliplr_bounds_local(bndsl_in, bounds) RESULT(bndsl_out)
1226 
1227  INTEGER, DIMENSION(2, 3), INTENT(IN) :: bndsl_in, bounds
1228  INTEGER, DIMENSION(2, 3) :: bndsl_out
1229 
1230  bndsl_out(1, 1) = bndsl_in(1, 1)
1231  bndsl_out(2, 1) = bndsl_in(2, 1)
1232 
1233  bndsl_out(1, 2) = 2*(bounds(2, 2) + 1) - bndsl_in(2, 2)
1234  bndsl_out(2, 2) = 2*(bounds(2, 2) + 1) - bndsl_in(1, 2)
1235  IF (bndsl_out(1, 2) .EQ. bounds(2, 2) + 2) bndsl_out(1, 2) = bndsl_out(1, 2) - 1
1236  IF (bndsl_out(2, 2) .EQ. 2*(bounds(2, 2) + 1) - bounds(1, 2)) bndsl_out(2, 2) = bndsl_out(2, 2) - 1
1237 
1238  bndsl_out(1, 3) = bndsl_in(1, 3)
1239  bndsl_out(2, 3) = bndsl_in(2, 3)
1240 
1241  END FUNCTION fliplr_bounds_local
1242 
1243 ! **************************************************************************************************
1244 !> \brief precalculates the local bounds of a 3d array after applying flipbf
1245 !> \param bndsl_in current local lower and upper bounds
1246 !> \param bounds global lower and upper bounds
1247 !> \return new local lower and upper bounds
1248 !> \par History
1249 !> 07.2014 created [Hossein Bani-Hashemian]
1250 !> \author Mohammad Hossein Bani-Hashemian
1251 ! **************************************************************************************************
1252  PURE FUNCTION flipbf_bounds_local(bndsl_in, bounds) RESULT(bndsl_out)
1253 
1254  INTEGER, DIMENSION(2, 3), INTENT(IN) :: bndsl_in, bounds
1255  INTEGER, DIMENSION(2, 3) :: bndsl_out
1256 
1257  bndsl_out(1, 1) = bndsl_in(1, 1)
1258  bndsl_out(2, 1) = bndsl_in(2, 1)
1259 
1260  bndsl_out(1, 2) = bndsl_in(1, 2)
1261  bndsl_out(2, 2) = bndsl_in(2, 2)
1262 
1263  bndsl_out(1, 3) = 2*(bounds(2, 3) + 1) - bndsl_in(2, 3)
1264  bndsl_out(2, 3) = 2*(bounds(2, 3) + 1) - bndsl_in(1, 3)
1265  IF (bndsl_out(1, 3) .EQ. bounds(2, 3) + 2) bndsl_out(1, 3) = bndsl_out(1, 3) - 1
1266  IF (bndsl_out(2, 3) .EQ. 2*(bounds(2, 3) + 1) - bounds(1, 3)) bndsl_out(2, 3) = bndsl_out(2, 3) - 1
1267 
1268  END FUNCTION flipbf_bounds_local
1269 
1270 ! **************************************************************************************************
1271 !> \brief precalculates the local bounds of a 3d array after applying rot180
1272 !> \param bndsl_in current local lower and upper bounds
1273 !> \param bounds global lower and upper bounds
1274 !> \return new local lower and upper bounds
1275 !> \par History
1276 !> 07.2014 created [Hossein Bani-Hashemian]
1277 !> \author Mohammad Hossein Bani-Hashemian
1278 ! **************************************************************************************************
1279  PURE FUNCTION rot180_bounds_local(bndsl_in, bounds) RESULT(bndsl_out)
1280 
1281  INTEGER, DIMENSION(2, 3), INTENT(IN) :: bndsl_in, bounds
1282  INTEGER, DIMENSION(2, 3) :: bndsl_out
1283 
1284  bndsl_out(1, 1) = 2*(bounds(2, 1) + 1) - bndsl_in(2, 1)
1285  bndsl_out(2, 1) = 2*(bounds(2, 1) + 1) - bndsl_in(1, 1)
1286  IF (bndsl_out(1, 1) .EQ. bounds(2, 1) + 2) bndsl_out(1, 1) = bndsl_out(1, 1) - 1
1287  IF (bndsl_out(2, 1) .EQ. 2*(bounds(2, 1) + 1) - bounds(1, 1)) bndsl_out(2, 1) = bndsl_out(2, 1) - 1
1288 
1289  bndsl_out(1, 2) = 2*(bounds(2, 2) + 1) - bndsl_in(2, 2)
1290  bndsl_out(2, 2) = 2*(bounds(2, 2) + 1) - bndsl_in(1, 2)
1291  IF (bndsl_out(1, 2) .EQ. bounds(2, 2) + 2) bndsl_out(1, 2) = bndsl_out(1, 2) - 1
1292  IF (bndsl_out(2, 2) .EQ. 2*(bounds(2, 2) + 1) - bounds(1, 2)) bndsl_out(2, 2) = bndsl_out(2, 2) - 1
1293 
1294  bndsl_out(1, 3) = bndsl_in(1, 3)
1295  bndsl_out(2, 3) = bndsl_in(2, 3)
1296 
1297  END FUNCTION rot180_bounds_local
1298 
1299 END MODULE dct
1300 
dct
the type I Discrete Cosine Transform (DCT-I)
Definition: dct.F:16
dct::dct_type_release
subroutine, public dct_type_release(dct_env)
Releases a dct_type.
Definition: dct.F:132
dct::neumannx
integer, parameter, public neumannx
Definition: dct.F:64
dct::neumannxy
integer, parameter, public neumannxy
Definition: dct.F:64
dct::dct_type_init
subroutine, public dct_type_init(pw_grid, neumann_directions, dct_env)
Initializes a dct_type.
Definition: dct.F:84
dct::neumannxz
integer, parameter, public neumannxz
Definition: dct.F:64
dct::pw_shrink
subroutine, public pw_shrink(neumann_directions, dests_shrink, srcs_shrink, bounds_local_shftd, pw_in, pw_shrinked)
shrinks an evenly symmetric pw_r3d_rs_type data to a pw_r3d_rs_type data that is 8 times smaller (the...
Definition: dct.F:702
dct::neumannxyz
integer, parameter, public neumannxyz
Definition: dct.F:64
dct::neumannz
integer, parameter, public neumannz
Definition: dct.F:64
dct::neumannyz
integer, parameter, public neumannyz
Definition: dct.F:64
dct::setup_dct_pw_grids
subroutine, public setup_dct_pw_grids(pw_grid, cell_hmat, neumann_directions, dct_pw_grid)
sets up an extended pw_grid for Discrete Cosine Transform (DCT) calculations
Definition: dct.F:163
dct::neumanny
integer, parameter, public neumanny
Definition: dct.F:64
dct::pw_expand
subroutine, public pw_expand(neumann_directions, recv_msgs_bnds, dests_expand, srcs_expand, flipg_stat, bounds_shftd, pw_in, pw_expanded)
expands a pw_r3d_rs_type data to an evenly symmetric pw_r3d_rs_type data that is 8 times larger than ...
Definition: dct.F:518
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
message_passing::mp_request_null
type(mp_request_type), parameter, public mp_request_null
Definition: message_passing.F:814
pw_grid_types
Definition: pw_grid_types.F:13
pw_grids
This module defines the grid data type and some basic operations on it.
Definition: pw_grids.F:36
pw_grids::pw_grid_setup
subroutine, public pw_grid_setup(cell_hmat, pw_grid, grid_span, cutoff, bounds, bounds_local, npts, spherical, odd, fft_usage, ncommensurate, icommensurate, blocked, ref_grid, rs_dims, iounit)
sets up a pw_grid
Definition: pw_grids.F:286
pw_grids::pw_grid_create
subroutine, public pw_grid_create(pw_grid, pe_group, local)
Initialize a PW grid with all defaults.
Definition: pw_grids.F:93
pw_types
Definition: pw_types.F:24