(git:374b731)
Loading...
Searching...
No Matches
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! **************************************************************************************************
16MODULE dct
17
18 USE kinds, ONLY: dp
19 USE message_passing, ONLY: mp_comm_type,&
24 USE pw_grids, ONLY: pw_grid_create,&
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, &
55 pw_shrink, &
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
72CONTAINS
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
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
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)
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))
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)
667 cr3d_xpndd(:, :, lb3_new:ind - 1) = catd
668 cr3d_xpndd(:, :, ind:ub3_new) = catd_flipdbf
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
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)
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
1299END MODULE dct
1300
the type I Discrete Cosine Transform (DCT-I)
Definition dct.F:16
subroutine, public dct_type_release(dct_env)
Releases a dct_type.
Definition dct.F:132
integer, parameter, public neumannx
Definition dct.F:64
integer, parameter, public neumannxy
Definition dct.F:64
subroutine, public dct_type_init(pw_grid, neumann_directions, dct_env)
Initializes a dct_type.
Definition dct.F:84
integer, parameter, public neumannxz
Definition dct.F:64
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
integer, parameter, public neumannxyz
Definition dct.F:64
integer, parameter, public neumannz
Definition dct.F:64
integer, parameter, public neumannyz
Definition dct.F:64
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
integer, parameter, public neumanny
Definition dct.F:64
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
Defines the basic variable types.
Definition kinds.F:23
integer, parameter, public dp
Definition kinds.F:34
Interface to the message passing library MPI.
type(mp_request_type), parameter, public mp_request_null
This module defines the grid data type and some basic operations on it.
Definition pw_grids.F:36
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
subroutine, public pw_grid_create(pw_grid, pe_group, local)
Initialize a PW grid with all defaults.
Definition pw_grids.F:93