(git:08825b9)
Loading...
Searching...
No Matches
fft_tools.F
Go to the documentation of this file.
1!--------------------------------------------------------------------------------------------------!
2! CP2K: A general program to perform molecular dynamics simulations !
3! Copyright 2000-2026 CP2K developers group <https://cp2k.org> !
4! !
5! SPDX-License-Identifier: GPL-2.0-or-later !
6!--------------------------------------------------------------------------------------------------!
7
8! **************************************************************************************************
9!> \par History
10!> JGH (30-Nov-2000): ESSL FFT Library added
11!> JGH (05-Jan-2001): Added SGI library FFT
12!> JGH (14-Jan-2001): Added parallel 3d FFT
13!> JGH (10-Feb-2006): New interface type
14!> JGH (31-Mar-2008): Remove local allocates and reshapes (performance)
15!> Possible problems can be related with setting arrays
16!> not to zero
17!> Some interfaces could be further simplified by avoiding
18!> an initial copy. However, this assumes contiguous arrays
19!> IAB (15-Oct-2008): Moved mp_cart_sub calls out of cube_tranpose_* and into
20!> fft_scratch type, reducing number of calls dramatically
21!> IAB (05-Dec-2008): Moved all other non-essential MPI calls into scratch type
22!> IAB (09-Jan-2009): Added fft_plan_type to store FFT data, including cached FFTW plans
23!> IAB (13-Feb-2009): Extended plan caching to serial 3D FFT (fft3d_s)
24!> IAB (09-Oct-2009): Added OpenMP directives to parallel 3D FFT
25!> (c) The Numerical Algorithms Group (NAG) Ltd, 2008-2009 on behalf of the HECToR project
26!> HFP (17-Oct-2024): Thread-safety insurance (CPASSERT), and OMP ATOMIC (tick_fft_pool)
27!> \author JGH
28! **************************************************************************************************
30 USE iso_c_binding, ONLY: c_f_pointer,&
31 c_loc,&
32 c_ptr,&
33 c_size_t
35 USE fft_lib, ONLY: &
36 fft_1d, fft_3d, fft_alloc, fft_create_plan_1d, fft_create_plan_3d, fft_dealloc, &
38 USE fft_plan, ONLY: fft_plan_type
39 USE kinds, ONLY: dp,&
40 dp_size,&
41 sp
42 USE mathconstants, ONLY: z_zero
43 USE message_passing, ONLY: mp_cart_type,&
50
51!$ USE OMP_LIB, ONLY: omp_get_max_threads, omp_get_thread_num, omp_in_parallel
52
53#include "../base/base_uses.f90"
54
55 IMPLICIT NONE
56
57 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'fft_tools'
58
59 ! Types for the pool of scratch data needed in FFT routines
60 ! keep the subroutine "is_equal" up-to-date
61 ! needs a default initialization
63 INTEGER :: nx = 0, ny = 0, nz = 0
64 INTEGER :: lmax = 0, mmax = 0, nmax = 0
65 INTEGER :: mx1 = 0, mx2 = 0, mx3 = 0
66 INTEGER :: my1 = 0, my2 = 0, my3 = 0
67 INTEGER :: mz1 = 0, mz2 = 0, mz3 = 0
68 INTEGER :: mcz1 = 0, mcz2 = 0, mcy3 = 0, mcx2 = 0
69 INTEGER :: lg = 0, mg = 0
70 INTEGER :: nbx = 0, nbz = 0
71 INTEGER :: nmray = 0, nyzray = 0
72 TYPE(mp_cart_type) :: rs_group = mp_cart_type()
73 INTEGER, DIMENSION(2) :: g_pos = 0, r_pos = 0, r_dim = 0
74 INTEGER :: numtask = 0
75 END TYPE fft_scratch_sizes
76
78 INTEGER :: fft_scratch_id = 0
79 INTEGER :: tf_type = -1
80 LOGICAL :: in_use = .true.
81 TYPE(mp_comm_type) :: group = mp_comm_type()
82 INTEGER, DIMENSION(3) :: nfft = -1
83 ! to be used in cube_transpose_* routines
84 TYPE(mp_cart_type), DIMENSION(2) :: cart_sub_comm = mp_cart_type()
85 INTEGER, DIMENSION(2) :: dim = -1, pos = -1
86 ! to be used in fft3d_s
87 COMPLEX(KIND=dp), DIMENSION(:, :, :), POINTER, CONTIGUOUS &
88 :: ziptr => null(), zoptr => null()
89 ! to be used in fft3d_ps : block distribution
90 COMPLEX(KIND=dp), DIMENSION(:, :), CONTIGUOUS, POINTER &
91 :: p1buf => null(), p2buf => null(), p3buf => null(), p4buf => null(), &
92 p5buf => null(), p6buf => null(), p7buf => null()
93 ! to be used in fft3d_ps : plane distribution
94 COMPLEX(KIND=dp), DIMENSION(:, :), POINTER, CONTIGUOUS &
95 :: r1buf => null(), r2buf => null()
96 COMPLEX(KIND=dp), DIMENSION(:, :, :), POINTER, CONTIGUOUS &
97 :: tbuf => null()
98 ! to be used in fft3d_pb
99 COMPLEX(KIND=dp), DIMENSION(:, :), POINTER, CONTIGUOUS &
100 :: a1buf => null(), a2buf => null(), a3buf => null(), &
101 a4buf => null(), a5buf => null(), a6buf => null()
102 ! to be used in communication routines
103 INTEGER, DIMENSION(:), CONTIGUOUS, POINTER :: scount => null(), rcount => null(), sdispl => null(), rdispl => null()
104 INTEGER, DIMENSION(:, :), CONTIGUOUS, POINTER :: pgcube => null()
105 INTEGER, DIMENSION(:), CONTIGUOUS, POINTER :: xzcount => null(), yzcount => null(), xzdispl => null(), yzdispl => null()
106 INTEGER :: in = 0, mip = -1
107 REAL(kind=dp) :: rsratio = 1.0_dp
108 COMPLEX(KIND=dp), DIMENSION(:), POINTER, CONTIGUOUS &
109 :: xzbuf => null(), yzbuf => null()
110 COMPLEX(KIND=sp), DIMENSION(:), POINTER, CONTIGUOUS &
111 :: xzbuf_sgl => null(), yzbuf_sgl => null()
112 COMPLEX(KIND=dp), DIMENSION(:, :), POINTER, CONTIGUOUS &
113 :: rbuf1 => null(), rbuf2 => null(), rbuf3 => null(), rbuf4 => null(), &
114 rbuf5 => null(), rbuf6 => null(), rr => null()
115 COMPLEX(KIND=sp), DIMENSION(:, :), POINTER, CONTIGUOUS &
116 :: ss => null(), tt => null()
117 INTEGER, DIMENSION(:, :), POINTER, CONTIGUOUS :: pgrid => null()
118 INTEGER, DIMENSION(:), POINTER, CONTIGUOUS :: xcor => null(), zcor => null(), pzcoord => null()
120 TYPE(fft_plan_type), DIMENSION(6) :: fft_plan = fft_plan_type()
121 INTEGER :: last_tick = -1
122 END TYPE fft_scratch_type
123
125 TYPE(fft_scratch_type), POINTER :: fft_scratch => null()
126 TYPE(fft_scratch_pool_type), POINTER :: fft_scratch_next => null()
127 END TYPE fft_scratch_pool_type
128
129 INTEGER, SAVE :: init_fft_pool = 0
130 ! the clock for fft pool. Allows to identify the least recently used scratch
131 INTEGER, SAVE :: tick_fft_pool = 0
132 ! limit the number of scratch pools to fft_pool_scratch_limit.
133 INTEGER, SAVE :: fft_pool_scratch_limit = 15
135 ! END of types for the pool of scratch data needed in FFT routines
136
137 PRIVATE
138 PUBLIC :: init_fft, fft3d, finalize_fft
139 PUBLIC :: fft_alloc, fft_dealloc
140 PUBLIC :: init_fft_scratch_pool
142 PUBLIC :: fwfft, bwfft
144 PUBLIC :: fft_radix_next_odd
145
146 INTEGER, PARAMETER :: fwfft = +1, bwfft = -1
147 INTEGER, PARAMETER :: fft_radix_closest = 493, fft_radix_next = 494
148 INTEGER, PARAMETER :: fft_radix_allowed = 495, fft_radix_disallowed = 496
149 INTEGER, PARAMETER :: fft_radix_next_odd = 497
150
151 REAL(kind=dp), PARAMETER :: ratio_sparse_alltoall = 0.5_dp
152
153 ! these saved variables are FFT globals
154 LOGICAL, SAVE :: alltoall_sgl = .false.
155
156 ! these are only needed for pw_gpu (-D__OFFLOAD)
159 PUBLIC :: yz_to_x, x_to_yz, xz_to_yz, yz_to_xz
161
162 INTERFACE fft3d
163 MODULE PROCEDURE fft3d_s, fft3d_ps, fft3d_pb
164 END INTERFACE
165
166! **************************************************************************************************
167
168CONTAINS
169
170! **************************************************************************************************
171!> \brief ...
172!> \param fftlib ...
173!> \param alltoall ...
174!> \param pool_limit ...
175!> \param wisdom_file ...
176!> \param plan_style ...
177!> \author JGH
178! **************************************************************************************************
179 SUBROUTINE init_fft(fftlib, alltoall, pool_limit, wisdom_file, &
180 plan_style)
181
182 CHARACTER(LEN=*), INTENT(IN) :: fftlib
183 LOGICAL, INTENT(IN) :: alltoall
184 INTEGER, INTENT(IN) :: pool_limit
185 CHARACTER(LEN=*), INTENT(IN) :: wisdom_file
186 INTEGER, INTENT(IN) :: plan_style
187
188 alltoall_sgl = alltoall
189 fft_pool_scratch_limit = pool_limit
190
191 CALL release_fft_scratch_pool()
192
193 CALL fft_do_init(fftlib, plan_style, wisdom_file)
194
195 ! setup the FFT scratch pool, if one is associated, clear first
197
198 END SUBROUTINE init_fft
199
200! **************************************************************************************************
201!> \brief does whatever is needed to finalize the current fft setup
202!> \param para_env ...
203!> \param wisdom_file ...
204!> \param export_wisdom Export the FFTW wisdom to wisdom_file (default: no export)
205!> \par History
206!> 10.2007 created [Joost VandeVondele]
207! **************************************************************************************************
208 SUBROUTINE finalize_fft(para_env, wisdom_file, export_wisdom)
209 CLASS(mp_comm_type) :: para_env
210 CHARACTER(LEN=*), INTENT(IN) :: wisdom_file
211 LOGICAL, INTENT(IN), OPTIONAL :: export_wisdom
212
213! release the FFT scratch pool
214
215 CALL release_fft_scratch_pool()
216
217 ! finalize fft libs
218
219 CALL fft_do_cleanup(wisdom_file, para_env%is_source(), export_wisdom)
220
221 END SUBROUTINE finalize_fft
222
223! **************************************************************************************************
224!> \brief Determine the allowed lengths of FFT's '''
225!> \param radix_in ...
226!> \param radix_out ...
227!> \param operation ...
228!> \par History
229!> new library structure (JGH)
230!> \author Ari Seitsonen
231! **************************************************************************************************
232 SUBROUTINE fft_radix_operations(radix_in, radix_out, operation)
233
234 INTEGER, INTENT(IN) :: radix_in
235 INTEGER, INTENT(OUT) :: radix_out
236 INTEGER, INTENT(IN) :: operation
237
238 INTEGER, PARAMETER :: number_of_small_radices = 81
239 INTEGER, DIMENSION(number_of_small_radices), PARAMETER :: small_radices = [2, 4, 6, 8, 9, 12,&
240 15, 16, 18, 20, 24, 25, 27, 30, 32, 36, 40, 45, 48, 54, 60, 64, 72, 75, 80, 81, 90, 96, &
241 100, 108, 120, 125, 128, 135, 144, 150, 160, 162, 180, 192, 200, 216, 225, 240, 243, 256, &
242 270, 288, 300, 320, 324, 360, 375, 384, 400, 405, 432, 450, 480, 486, 500, 512, 540, 576, &
243 600, 625, 640, 648, 675, 720, 729, 750, 768, 800, 810, 864, 900, 960, 972, 1000, 1024]
244
245 INTEGER :: i, iloc, ldata, trial_number
246 INTEGER, ALLOCATABLE, DIMENSION(:) :: data
247
248 ldata = 1024
249 ALLOCATE (DATA(ldata))
250 DATA = -1
251
252 DATA(1:number_of_small_radices) = small_radices(1:number_of_small_radices)
253
254 ! First run: try smaller radices
255 ! Second run: Fetch extended list of radices and try these
256 DO trial_number = 1, 2
257 IF (trial_number == 2) CALL fft_get_lengths(DATA, ldata)
258 iloc = 0
259 DO i = 1, ldata
260 IF (DATA(i) == radix_in) THEN
261 iloc = i
262 EXIT
263 ELSE
264 IF (operation == fft_radix_allowed) THEN
265 cycle
266 ELSE IF (DATA(i) > radix_in) THEN
267 iloc = i
268 EXIT
269 END IF
270 END IF
271 END DO
272
273 IF (iloc == 0) THEN
274 IF (trial_number == 1) THEN
275 ! We haven't found any suitable size in the list of smaller radices
276 ! because the requested radix is too large (large cell or large cutoff)
277 ! So, try again with the larger radices
278 cycle
279 ELSE
280 IF (operation == fft_radix_closest) THEN
281 radix_out = maxval(data)
282 ELSE
283 CALL cp_abort(__location__, &
284 "No suitable FFT lengths found! Either enable FFTW or reduce the grid size!")
285 END IF
286 RETURN
287 END IF
288
289 END IF
290
291 IF (operation == fft_radix_allowed) THEN
292 IF (DATA(iloc) == radix_in) THEN
293 radix_out = fft_radix_allowed
294 ELSE
295 radix_out = fft_radix_disallowed
296 END IF
297
298 ELSE IF (operation == fft_radix_closest) THEN
299 IF (DATA(iloc) == radix_in) THEN
300 radix_out = DATA(iloc)
301 ELSE
302 IF (abs(DATA(iloc - 1) - radix_in) <= &
303 abs(DATA(iloc) - radix_in)) THEN
304 radix_out = DATA(iloc - 1)
305 ELSE
306 radix_out = DATA(iloc)
307 END IF
308 END IF
309
310 ELSE IF (operation == fft_radix_next) THEN
311 radix_out = DATA(iloc)
312
313 ELSE IF (operation == fft_radix_next_odd) THEN
314 DO i = iloc, ldata
315 IF (mod(DATA(i), 2) == 1) THEN
316 radix_out = DATA(i)
317 EXIT
318 END IF
319 END DO
320 IF (mod(radix_out, 2) == 0) THEN
321 ! No larger odd radix found
322 IF (trial_number == 1) THEN
323 cycle
324 ELSE
325 radix_out = radix_in + 1 - mod(radix_in, 2)
326 END IF
327 END IF
328
329 END IF
330 EXIT
331 END DO
332
333 DEALLOCATE (data)
334
335 END SUBROUTINE fft_radix_operations
336
337! **************************************************************************************************
338!> \brief Performs m 1-D forward FFT-s of size n.
339!> \param sign ...
340!> \param n size of the FFT
341!> \param m number of FFT-s
342!> \param trans_in transposition state of the input array
343!> \param trans_out transposition state of the output array
344!> \param ldx_in ...
345!> \param ldx_out ...
346!> \param zin input array
347!> \param zout output array
348!> \param scale scaling factor
349!> \param stat status of the operation, non-zero code indicates an error
350! **************************************************************************************************
351 SUBROUTINE fft_1d_many(sign, n, m, trans_in, trans_out, ldx_in, ldx_out, zin, zout, scale, stat)
352 INTEGER, INTENT(in) :: sign, n, m
353 LOGICAL, INTENT(in) :: trans_in, trans_out
354 INTEGER, INTENT(in) :: ldx_in, ldx_out
355 COMPLEX(kind=dp), DIMENSION(*), INTENT(inout) :: zin, zout
356 REAL(kind=dp), INTENT(in) :: scale
357 INTEGER, INTENT(out) :: stat
358
359 CHARACTER(len=*), PARAMETER :: routinen = 'fft_1d_many'
360
361 COMPLEX(KIND=dp), CONTIGUOUS, DIMENSION(:), &
362 POINTER :: ztmp
363 INTEGER :: handle
364 TYPE(fft_plan_type) :: fft_plan
365
366 CALL timeset(routinen, handle)
367
368 IF (sign /= fwfft .AND. sign /= bwfft) THEN
369 cpabort("Unknown FFT-direction!")
370 END IF
371
372 NULLIFY (ztmp)
373 IF (trans_in) THEN
374 CALL fft_alloc(ztmp, [ldx_in*n])
375 ELSE
376 CALL fft_alloc(ztmp, [ldx_in*m])
377 END IF
378 CALL fft_create_plan_1d(fft_plan, sign, trans_in, trans_out, ldx_in, ldx_out, n, m, ztmp, zout)
379 CALL fft_dealloc(ztmp)
380 CALL fft_1d(fft_plan, zin, zout, scale, stat)
382
383 CALL timestop(handle)
384 END SUBROUTINE fft_1d_many
385
386! **************************************************************************************************
387!> \brief Calls the 3D-FFT function from the initialized library
388!> \param fsign ...
389!> \param n ...
390!> \param zin ...
391!> \param zout ...
392!> \param status ...
393!> \param debug ...
394!> \par History
395!> none
396!> \author JGH
397! **************************************************************************************************
398 SUBROUTINE fft3d_s(fsign, n, zin, zout, status, debug)
399
400 INTEGER, INTENT(IN) :: fsign
401 INTEGER, DIMENSION(:), INTENT(INOUT) :: n
402 COMPLEX(KIND=dp), CONTIGUOUS, DIMENSION(:, :, :), &
403 INTENT(INOUT) :: zin
404 COMPLEX(KIND=dp), CONTIGUOUS, DIMENSION(:, :, :), &
405 INTENT(INOUT), OPTIONAL, TARGET :: zout
406 INTEGER, INTENT(OUT), OPTIONAL :: status
407 LOGICAL, INTENT(IN), OPTIONAL :: debug
408
409 CHARACTER(len=*), PARAMETER :: routineN = 'fft3d_s'
410
411 COMPLEX(KIND=dp), CONTIGUOUS, DIMENSION(:, :, :), &
412 POINTER :: zoptr
413 COMPLEX(KIND=dp), DIMENSION(1, 1, 1), TARGET :: zdum
414 INTEGER :: handle, ld(3), lo(3), output_unit, sign, &
415 stat
416 LOGICAL :: fft_in_place, test
417 REAL(KIND=dp) :: in_sum, norm, out_sum
418 TYPE(fft_scratch_type), POINTER :: fft_scratch
419
420 CALL timeset(routinen, handle)
421 output_unit = cp_logger_get_default_io_unit()
422
423 IF (fsign == fwfft) THEN
424 norm = 1.0_dp/real(product(n), kind=dp)
425 ELSE IF (fsign == bwfft) THEN
426 norm = 1.0_dp
427 ELSE
428 cpabort("Unknown FFT direction!")
429 END IF
430
431 IF (PRESENT(debug)) THEN
432 test = debug
433 ELSE
434 test = .false.
435 END IF
436
437 IF (PRESENT(zout)) THEN
438 fft_in_place = .false.
439 ELSE
440 fft_in_place = .true.
441 END IF
442
443 IF (test) THEN
444 in_sum = sum(abs(zin))
445 END IF
446
447 ld(1) = SIZE(zin, 1)
448 ld(2) = SIZE(zin, 2)
449 ld(3) = SIZE(zin, 3)
450
451 IF (n(1) /= ld(1) .OR. n(2) /= ld(2) .OR. n(3) /= ld(3)) THEN
452 cpabort("Size and dimension (zin) have to be the same.")
453 END IF
454
455 sign = fsign
456 CALL get_fft_scratch(fft_scratch, tf_type=400, n=n)
457
458 IF (fft_in_place) THEN
459 zoptr => zdum
460 IF (fsign == fwfft) THEN
461 CALL fft_3d(fft_scratch%fft_plan(1), norm, zin, zoptr, stat)
462 ELSE
463 CALL fft_3d(fft_scratch%fft_plan(2), norm, zin, zoptr, stat)
464 END IF
465 ELSE
466 IF (fsign == fwfft) THEN
467 CALL fft_3d(fft_scratch%fft_plan(3), norm, zin, zout, stat)
468 ELSE
469 CALL fft_3d(fft_scratch%fft_plan(4), norm, zin, zout, stat)
470 END IF
471 END IF
472
473 CALL release_fft_scratch(fft_scratch)
474
475 IF (PRESENT(zout)) THEN
476 lo(1) = SIZE(zout, 1)
477 lo(2) = SIZE(zout, 2)
478 lo(3) = SIZE(zout, 3)
479 IF (n(1) /= lo(1) .OR. n(2) /= lo(2) .OR. n(3) /= lo(3)) THEN
480 cpabort("Size and dimension (zout) have to be the same.")
481 END IF
482 END IF
483
484 IF (PRESENT(status)) THEN
485 status = stat
486 END IF
487
488 IF (test .AND. output_unit > 0) THEN
489 IF (PRESENT(zout)) THEN
490 out_sum = sum(abs(zout))
491 WRITE (output_unit, '(A)') " Out of place 3D FFT (local) : fft3d_s"
492 WRITE (output_unit, '(A,T60,3I7)') " Transform lengths ", n
493 WRITE (output_unit, '(A,T60,3I7)') " Input array dimensions ", ld
494 WRITE (output_unit, '(A,T60,3I7)') " Output array dimensions ", lo
495 WRITE (output_unit, '(A,T61,E20.14)') " Sum of input data ", in_sum
496 WRITE (output_unit, '(A,T61,E20.14)') " Sum of output data ", out_sum
497 ELSE
498 out_sum = sum(abs(zin))
499 WRITE (output_unit, '(A)') " In place 3D FFT (local) : fft3d_s"
500 WRITE (output_unit, '(A,T60,3I7)') " Transform lengths ", n
501 WRITE (output_unit, '(A,T60,3I7)') " Input/output array dimensions ", ld
502 WRITE (output_unit, '(A,T61,E20.14)') " Sum of input data ", in_sum
503 WRITE (output_unit, '(A,T61,E20.14)') " Sum of output data ", out_sum
504 END IF
505 END IF
506
507 CALL timestop(handle)
508
509 END SUBROUTINE fft3d_s
510
511! **************************************************************************************************
512!> \brief ...
513!> \param fsign ...
514!> \param n ...
515!> \param cin ...
516!> \param gin ...
517!> \param rs_group ...
518!> \param yzp ...
519!> \param nyzray ...
520!> \param bo ...
521!> \param status ...
522!> \param debug ...
523! **************************************************************************************************
524 SUBROUTINE fft3d_ps(fsign, n, cin, gin, rs_group, yzp, nyzray, &
525 bo, status, debug)
526
527 INTEGER, INTENT(IN) :: fsign
528 INTEGER, DIMENSION(:), INTENT(IN) :: n
529 COMPLEX(KIND=dp), CONTIGUOUS, DIMENSION(:, :, :), &
530 INTENT(INOUT) :: cin
531 COMPLEX(KIND=dp), CONTIGUOUS, DIMENSION(:, :), &
532 INTENT(INOUT) :: gin
533 TYPE(mp_cart_type), INTENT(IN) :: rs_group
534 INTEGER, CONTIGUOUS, DIMENSION(:, :, 0:), &
535 INTENT(IN) :: yzp
536 INTEGER, CONTIGUOUS, DIMENSION(0:), INTENT(IN) :: nyzray
537 INTEGER, CONTIGUOUS, DIMENSION(:, :, 0:, :), &
538 INTENT(IN) :: bo
539 INTEGER, INTENT(OUT), OPTIONAL :: status
540 LOGICAL, INTENT(IN), OPTIONAL :: debug
541
542 CHARACTER(len=*), PARAMETER :: routineN = 'fft3d_ps'
543
544 COMPLEX(KIND=dp), CONTIGUOUS, DIMENSION(:, :), &
545 POINTER :: pbuf, qbuf, rbuf, sbuf
546 COMPLEX(KIND=dp), CONTIGUOUS, DIMENSION(:, :, :), &
547 POINTER :: tbuf
548 INTEGER :: g_pos, handle, lg, lmax, mcx2, mcz1, mcz2, mg, mmax, mx1, mx2, my1, mz2, n1, n2, &
549 nmax, numtask, nx, ny, nz, output_unit, r_dim(2), r_pos(2), rp, sign, stat
550 INTEGER, ALLOCATABLE, DIMENSION(:) :: p2p
551 LOGICAL :: test
552 REAL(KIND=dp) :: norm, sum_data
553 TYPE(fft_scratch_sizes) :: fft_scratch_size
554 TYPE(fft_scratch_type), POINTER :: fft_scratch
555
556 CALL timeset(routinen, handle)
557 output_unit = cp_logger_get_default_io_unit()
558
559 IF (PRESENT(debug)) THEN
560 test = debug
561 ELSE
562 test = .false.
563 END IF
564
565 g_pos = rs_group%mepos
566 numtask = rs_group%num_pe
567 r_dim = rs_group%num_pe_cart
568 r_pos = rs_group%mepos_cart
569
570 IF (fsign == fwfft) THEN
571 norm = 1.0_dp/real(product(n), kind=dp)
572 ELSE IF (fsign == bwfft) THEN
573 norm = 1.0_dp
574 ELSE
575 cpabort("Unknown FFT direction!")
576 END IF
577
578 sign = fsign
579
580 lg = SIZE(gin, 1)
581 mg = SIZE(gin, 2)
582
583 nx = SIZE(cin, 1)
584 ny = SIZE(cin, 2)
585 nz = SIZE(cin, 3)
586
587 IF (mg == 0) THEN
588 mmax = 1
589 ELSE
590 mmax = mg
591 END IF
592 lmax = max(lg, (nx*ny*nz)/mmax + 1)
593
594 ALLOCATE (p2p(0:numtask - 1))
595
596 CALL rs_group%rank_compare(rs_group, p2p)
597
598 rp = p2p(g_pos)
599 mx1 = bo(2, 1, rp, 1) - bo(1, 1, rp, 1) + 1
600 my1 = bo(2, 2, rp, 1) - bo(1, 2, rp, 1) + 1
601 mx2 = bo(2, 1, rp, 2) - bo(1, 1, rp, 2) + 1
602 mz2 = bo(2, 3, rp, 2) - bo(1, 3, rp, 2) + 1
603
604 n1 = maxval(bo(2, 1, :, 1) - bo(1, 1, :, 1) + 1)
605 n2 = maxval(bo(2, 2, :, 1) - bo(1, 2, :, 1) + 1)
606 nmax = max((2*n2)/numtask, 2)*mx2*mz2
607 nmax = max(nmax, n1*maxval(nyzray))
608 n1 = maxval(bo(2, 1, :, 2))
609 n2 = maxval(bo(2, 3, :, 2))
610
611 fft_scratch_size%nx = nx
612 fft_scratch_size%ny = ny
613 fft_scratch_size%nz = nz
614 fft_scratch_size%lmax = lmax
615 fft_scratch_size%mmax = mmax
616 fft_scratch_size%mx1 = mx1
617 fft_scratch_size%mx2 = mx2
618 fft_scratch_size%my1 = my1
619 fft_scratch_size%mz2 = mz2
620 fft_scratch_size%lg = lg
621 fft_scratch_size%mg = mg
622 fft_scratch_size%nbx = n1
623 fft_scratch_size%nbz = n2
624 mcz1 = maxval(bo(2, 3, :, 1) - bo(1, 3, :, 1) + 1)
625 mcx2 = maxval(bo(2, 1, :, 2) - bo(1, 1, :, 2) + 1)
626 mcz2 = maxval(bo(2, 3, :, 2) - bo(1, 3, :, 2) + 1)
627 fft_scratch_size%mcz1 = mcz1
628 fft_scratch_size%mcx2 = mcx2
629 fft_scratch_size%mcz2 = mcz2
630 fft_scratch_size%nmax = nmax
631 fft_scratch_size%nmray = maxval(nyzray)
632 fft_scratch_size%nyzray = nyzray(g_pos)
633 fft_scratch_size%rs_group = rs_group
634 fft_scratch_size%g_pos = g_pos
635 fft_scratch_size%r_pos = r_pos
636 fft_scratch_size%r_dim = r_dim
637 fft_scratch_size%numtask = numtask
638
639 IF (test) THEN
640 IF (g_pos == 0 .AND. output_unit > 0) THEN
641 WRITE (output_unit, '(A)') " Parallel 3D FFT : fft3d_ps"
642 WRITE (output_unit, '(A,T60,3I7)') " Transform lengths ", n
643 WRITE (output_unit, '(A,T67,2I7)') " Array dimensions (gin) ", lg, mg
644 WRITE (output_unit, '(A,T60,3I7)') " Array dimensions (cin) ", nx, ny, nz
645 END IF
646 END IF
647
648 IF (r_dim(2) > 1) THEN
649
650 !
651 ! real space is distributed over x and y coordinate
652 ! we have two stages of communication
653 !
654
655 IF (r_dim(1) == 1) THEN
656 cpabort("This processor distribution is not supported.")
657 END IF
658 CALL get_fft_scratch(fft_scratch, tf_type=300, n=n, fft_sizes=fft_scratch_size)
659
660 IF (sign == fwfft) THEN
661 ! cin -> gin
662
663 IF (test) THEN
664 sum_data = abs(sum(cin))
665 CALL rs_group%sum(sum_data)
666 IF (g_pos == 0 .AND. output_unit > 0) THEN
667 WRITE (output_unit, '(A)') " Two step communication algorithm "
668 WRITE (output_unit, '(A,T60,3I7)') " Transform Z ", n(3), mx1*my1
669 WRITE (output_unit, '(A,T60,3I7)') " Transform Y ", n(2), mx2*mz2
670 WRITE (output_unit, '(A,T67,2I7)') " Transform X ", n(1), nyzray(g_pos)
671 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(1) ", sum_data
672 END IF
673 END IF
674
675 pbuf => fft_scratch%p1buf
676 qbuf => fft_scratch%p2buf
677
678 ! FFT along z
679 CALL fft_1d(fft_scratch%fft_plan(1), cin, qbuf, norm, stat)
680
681 rbuf => fft_scratch%p3buf
682
683 IF (test) THEN
684 sum_data = abs(sum(qbuf))
685 CALL rs_group%sum(sum_data)
686 IF (g_pos == 0 .AND. output_unit > 0) THEN
687 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(2) T", sum_data
688 END IF
689 END IF
690
691 ! Exchange data ( transpose of matrix )
692 CALL cube_transpose_2(qbuf, bo(:, :, :, 1), bo(:, :, :, 2), rbuf, fft_scratch)
693
694 IF (test) THEN
695 sum_data = abs(sum(rbuf))
696 CALL rs_group%sum(sum_data)
697 IF (g_pos == 0 .AND. output_unit > 0) THEN
698 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(3) T", sum_data
699 END IF
700 END IF
701
702 pbuf => fft_scratch%p4buf
703
704 ! FFT along y
705 CALL fft_1d(fft_scratch%fft_plan(2), rbuf, pbuf, 1.0_dp, stat)
706
707 qbuf => fft_scratch%p5buf
708
709 IF (test) THEN
710 sum_data = abs(sum(pbuf))
711 CALL rs_group%sum(sum_data)
712 IF (g_pos == 0 .AND. output_unit > 0) THEN
713 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(4) TS", sum_data
714 END IF
715 END IF
716
717 ! Exchange data ( transpose of matrix ) and sort
718 CALL xz_to_yz(pbuf, rs_group, r_dim, g_pos, p2p, yzp, nyzray, &
719 bo(:, :, :, 2), qbuf, fft_scratch)
720
721 IF (test) THEN
722 sum_data = abs(sum(qbuf))
723 CALL rs_group%sum(sum_data)
724 IF (g_pos == 0 .AND. output_unit > 0) THEN
725 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(5) TS", sum_data
726 END IF
727 END IF
728
729 ! FFT along x
730 CALL fft_1d(fft_scratch%fft_plan(3), qbuf, gin, 1.0_dp, stat)
731
732 IF (test) THEN
733 sum_data = abs(sum(gin))
734 CALL rs_group%sum(sum_data)
735 IF (g_pos == 0 .AND. output_unit > 0) THEN
736 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(6) ", sum_data
737 END IF
738 END IF
739
740 ELSE IF (sign == bwfft) THEN
741 ! gin -> cin
742
743 IF (test) THEN
744 sum_data = abs(sum(gin))
745 CALL rs_group%sum(sum_data)
746 IF (g_pos == 0 .AND. output_unit > 0) THEN
747 WRITE (output_unit, '(A)') " Two step communication algorithm "
748 WRITE (output_unit, '(A,T67,2I7)') " Transform X ", n(1), nyzray(g_pos)
749 WRITE (output_unit, '(A,T60,3I7)') " Transform Y ", n(2), mx2*mz2
750 WRITE (output_unit, '(A,T60,3I7)') " Transform Z ", n(3), mx1*my1
751 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(1) ", sum_data
752 END IF
753 END IF
754
755 pbuf => fft_scratch%p7buf
756
757 ! FFT along x
758 CALL fft_1d(fft_scratch%fft_plan(4), gin, pbuf, norm, stat)
759
760 qbuf => fft_scratch%p4buf
761
762 IF (test) THEN
763 sum_data = abs(sum(pbuf))
764 CALL rs_group%sum(sum_data)
765 IF (g_pos == 0 .AND. output_unit > 0) THEN
766 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(2) TS", sum_data
767 END IF
768 END IF
769
770 ! Exchange data ( transpose of matrix ) and sort
771 CALL yz_to_xz(pbuf, rs_group, r_dim, g_pos, p2p, yzp, nyzray, &
772 bo(:, :, :, 2), qbuf, fft_scratch)
773
774 IF (test) THEN
775 sum_data = abs(sum(qbuf))
776 CALL rs_group%sum(sum_data)
777 IF (g_pos == 0 .AND. output_unit > 0) THEN
778 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(3) TS", sum_data
779 END IF
780 END IF
781
782 rbuf => fft_scratch%p3buf
783
784 ! FFT along y
785 CALL fft_1d(fft_scratch%fft_plan(5), qbuf, rbuf, 1.0_dp, stat)
786
787 pbuf => fft_scratch%p2buf
788
789 IF (test) THEN
790 sum_data = abs(sum(rbuf))
791 CALL rs_group%sum(sum_data)
792 IF (g_pos == 0 .AND. output_unit > 0) THEN
793 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(4) T", sum_data
794 END IF
795 END IF
796
797 ! Exchange data ( transpose of matrix )
798 CALL cube_transpose_1(rbuf, bo(:, :, :, 2), bo(:, :, :, 1), pbuf, fft_scratch)
799
800 IF (test) THEN
801 sum_data = abs(sum(pbuf))
802 CALL rs_group%sum(sum_data)
803 IF (g_pos == 0 .AND. output_unit > 0) THEN
804 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(5) T", sum_data
805 END IF
806 END IF
807
808 qbuf => fft_scratch%p1buf
809
810 ! FFT along z
811 CALL fft_1d(fft_scratch%fft_plan(6), pbuf, cin, 1.0_dp, stat)
812
813 IF (test) THEN
814 sum_data = abs(sum(cin))
815 CALL rs_group%sum(sum_data)
816 IF (g_pos == 0 .AND. output_unit > 0) THEN
817 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(6) ", sum_data
818 END IF
819 END IF
820
821 ELSE
822
823 cpabort("Illegal fsign parameter.")
824
825 END IF
826
827 CALL release_fft_scratch(fft_scratch)
828
829 ELSE
830
831 !
832 ! real space is only distributed over x coordinate
833 ! we have one stage of communication, after the transform of
834 ! direction x
835 !
836
837 CALL get_fft_scratch(fft_scratch, tf_type=200, n=n, fft_sizes=fft_scratch_size)
838
839 sbuf => fft_scratch%r1buf
840 tbuf => fft_scratch%tbuf
841
842 sbuf = z_zero
843 tbuf = z_zero
844
845 IF (sign == fwfft) THEN
846 ! cin -> gin
847
848 IF (test) THEN
849 sum_data = abs(sum(cin))
850 CALL rs_group%sum(sum_data)
851 IF (g_pos == 0 .AND. output_unit > 0) THEN
852 WRITE (output_unit, '(A)') " One step communication algorithm "
853 WRITE (output_unit, '(A,T60,3I7)') " Transform YZ ", n(2), n(3), nx
854 WRITE (output_unit, '(A,T67,2I7)') " Transform X ", n(1), nyzray(g_pos)
855 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(1) ", sum_data
856 END IF
857 END IF
858
859 ! FFT along y and z
860 CALL fft_1d(fft_scratch%fft_plan(1), cin, sbuf, 1._dp, stat)
861 CALL fft_1d(fft_scratch%fft_plan(2), sbuf, tbuf, 1._dp, stat)
862
863 IF (test) THEN
864 sum_data = abs(sum(tbuf))
865 CALL rs_group%sum(sum_data)
866 IF (g_pos == 0 .AND. output_unit > 0) THEN
867 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(2) TS", sum_data
868 END IF
869 END IF
870
871 ! Exchange data ( transpose of matrix ) and sort
872 CALL yz_to_x(tbuf, rs_group, g_pos, p2p, yzp, nyzray, &
873 bo(:, :, :, 2), sbuf, fft_scratch)
874
875 IF (test) THEN
876 sum_data = abs(sum(sbuf))
877 CALL rs_group%sum(sum_data)
878 IF (g_pos == 0 .AND. output_unit > 0) THEN
879 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(3) TS", sum_data
880 END IF
881 END IF
882 ! FFT along x
883 CALL fft_1d(fft_scratch%fft_plan(3), sbuf, gin, norm, stat)
884
885 IF (test) THEN
886 sum_data = abs(sum(gin))
887 CALL rs_group%sum(sum_data)
888 IF (g_pos == 0 .AND. output_unit > 0) THEN
889 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(4) ", sum_data
890 END IF
891 END IF
892
893 ELSE IF (sign == bwfft) THEN
894 ! gin -> cin
895
896 IF (test) THEN
897 sum_data = abs(sum(gin))
898 CALL rs_group%sum(sum_data)
899 IF (g_pos == 0 .AND. output_unit > 0) THEN
900 WRITE (output_unit, '(A)') " One step communication algorithm "
901 WRITE (output_unit, '(A,T67,2I7)') " Transform X ", n(1), nyzray(g_pos)
902 WRITE (output_unit, '(A,T60,3I7)') " Transform YZ ", n(2), n(3), nx
903 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(1) ", sum_data
904 END IF
905 END IF
906
907 ! FFT along x
908 CALL fft_1d(fft_scratch%fft_plan(4), gin, sbuf, norm, stat)
909
910 IF (test) THEN
911 sum_data = abs(sum(sbuf))
912 CALL rs_group%sum(sum_data)
913 IF (g_pos == 0 .AND. output_unit > 0) THEN
914 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(2) TS", sum_data
915 END IF
916 END IF
917
918 ! Exchange data ( transpose of matrix ) and sort
919 CALL x_to_yz(sbuf, rs_group, g_pos, p2p, yzp, nyzray, &
920 bo(:, :, :, 2), tbuf, fft_scratch)
921
922 IF (test) THEN
923 sum_data = abs(sum(tbuf))
924 CALL rs_group%sum(sum_data)
925 IF (g_pos == 0 .AND. output_unit > 0) THEN
926 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(3) TS", sum_data
927 END IF
928 END IF
929
930 ! FFT along y and z
931 CALL fft_1d(fft_scratch%fft_plan(5), tbuf, sbuf, 1._dp, stat)
932 CALL fft_1d(fft_scratch%fft_plan(6), sbuf, cin, 1._dp, stat)
933
934 IF (test) THEN
935 sum_data = abs(sum(cin))
936 CALL rs_group%sum(sum_data)
937 IF (g_pos == 0 .AND. output_unit > 0) THEN
938 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(4) ", sum_data
939 END IF
940 END IF
941 ELSE
942 cpabort("Illegal fsign parameter.")
943 END IF
944
945 CALL release_fft_scratch(fft_scratch)
946
947 END IF
948
949 DEALLOCATE (p2p)
950
951 IF (PRESENT(status)) THEN
952 status = stat
953 END IF
954 CALL timestop(handle)
955
956 END SUBROUTINE fft3d_ps
957
958! **************************************************************************************************
959!> \brief ...
960!> \param fsign ...
961!> \param n ...
962!> \param zin ...
963!> \param gin ...
964!> \param group ...
965!> \param bo ...
966!> \param status ...
967!> \param debug ...
968! **************************************************************************************************
969 SUBROUTINE fft3d_pb(fsign, n, zin, gin, group, bo, status, debug)
970
971 INTEGER, INTENT(IN) :: fsign
972 INTEGER, DIMENSION(3), INTENT(IN) :: n
973 COMPLEX(KIND=dp), CONTIGUOUS, DIMENSION(:, :, :), &
974 INTENT(INOUT) :: zin
975 COMPLEX(KIND=dp), CONTIGUOUS, DIMENSION(:, :), &
976 INTENT(INOUT) :: gin
977 TYPE(mp_cart_type), INTENT(IN) :: group
978 INTEGER, CONTIGUOUS, DIMENSION(:, :, 0:, :), &
979 INTENT(IN) :: bo
980 INTEGER, INTENT(OUT), OPTIONAL :: status
981 LOGICAL, INTENT(IN), OPTIONAL :: debug
982
983 CHARACTER(len=*), PARAMETER :: routineN = 'fft3d_pb'
984
985 COMPLEX(KIND=dp), CONTIGUOUS, DIMENSION(:, :), &
986 POINTER :: abuf, bbuf
987 INTEGER :: handle, lg(2), lz(3), mcx2, mcy3, mcz1, &
988 mcz2, mx1, mx2, mx3, my1, my2, my3, &
989 my_pos, mz1, mz2, mz3, output_unit, &
990 sign, stat
991 INTEGER, DIMENSION(2) :: dim
992 LOGICAL :: test
993 REAL(KIND=dp) :: norm, sum_data
994 TYPE(fft_scratch_sizes) :: fft_scratch_size
995 TYPE(fft_scratch_type), POINTER :: fft_scratch
996
997!------------------------------------------------------------------------------
998! "Real Space" 1) xyZ or 1) xYZ
999! 2) xYz or not used
1000! "G Space" 3) Xyz or 3) XYz
1001!
1002! There is one communicator (2-dimensional) for all distributions
1003! np = n1 * n2, where np is the total number of processors
1004! If n2 = 1, we have the second case and only one transpose step is needed
1005!
1006! Assignment of dimensions to axis for different steps
1007! First case: 1) n1=x; n2=y
1008! 2) n1=x; n2=z
1009! 3) n1=y; n2=z
1010! Second case 1) n1=x
1011! 3) n1=z
1012!
1013! The more general case with two communicators for the initial and final
1014! distribution is not covered.
1015!------------------------------------------------------------------------------
1016
1017 CALL timeset(routinen, handle)
1018 output_unit = cp_logger_get_default_io_unit()
1019
1020 dim = group%num_pe_cart
1021 my_pos = group%mepos
1022
1023 IF (PRESENT(debug)) THEN
1024 test = debug
1025 ELSE
1026 test = .false.
1027 END IF
1028
1029 IF (fsign == fwfft) THEN
1030 norm = 1.0_dp/real(product(n), kind=dp)
1031 ELSE IF (fsign == bwfft) THEN
1032 norm = 1.0_dp
1033 ELSE
1034 cpabort("Unknown FFT direction!")
1035 END IF
1036
1037 sign = fsign
1038
1039 IF (test) THEN
1040 lg(1) = SIZE(gin, 1)
1041 lg(2) = SIZE(gin, 2)
1042 lz(1) = SIZE(zin, 1)
1043 lz(2) = SIZE(zin, 2)
1044 lz(3) = SIZE(zin, 3)
1045 IF (my_pos == 0 .AND. output_unit > 0) THEN
1046 WRITE (output_unit, '(A)') " Parallel 3D FFT : fft3d_pb"
1047 WRITE (output_unit, '(A,T60,3I7)') " Transform lengths ", n
1048 WRITE (output_unit, '(A,T67,2I7)') " Array dimensions (gin) ", lg
1049 WRITE (output_unit, '(A,T60,3I7)') " Array dimensions (cin) ", lz
1050 END IF
1051 END IF
1052
1053 mx1 = bo(2, 1, my_pos, 1) - bo(1, 1, my_pos, 1) + 1
1054 my1 = bo(2, 2, my_pos, 1) - bo(1, 2, my_pos, 1) + 1
1055 mz1 = bo(2, 3, my_pos, 1) - bo(1, 3, my_pos, 1) + 1
1056 mx2 = bo(2, 1, my_pos, 2) - bo(1, 1, my_pos, 2) + 1
1057 my2 = bo(2, 2, my_pos, 2) - bo(1, 2, my_pos, 2) + 1
1058 mz2 = bo(2, 3, my_pos, 2) - bo(1, 3, my_pos, 2) + 1
1059 mx3 = bo(2, 1, my_pos, 3) - bo(1, 1, my_pos, 3) + 1
1060 my3 = bo(2, 2, my_pos, 3) - bo(1, 2, my_pos, 3) + 1
1061 mz3 = bo(2, 3, my_pos, 3) - bo(1, 3, my_pos, 3) + 1
1062 fft_scratch_size%mx1 = mx1
1063 fft_scratch_size%mx2 = mx2
1064 fft_scratch_size%mx3 = mx3
1065 fft_scratch_size%my1 = my1
1066 fft_scratch_size%my2 = my2
1067 fft_scratch_size%my3 = my3
1068 fft_scratch_size%mz1 = mz1
1069 fft_scratch_size%mz2 = mz2
1070 fft_scratch_size%mz3 = mz3
1071 mcz1 = maxval(bo(2, 3, :, 1) - bo(1, 3, :, 1) + 1)
1072 mcx2 = maxval(bo(2, 1, :, 2) - bo(1, 1, :, 2) + 1)
1073 mcz2 = maxval(bo(2, 3, :, 2) - bo(1, 3, :, 2) + 1)
1074 mcy3 = maxval(bo(2, 2, :, 3) - bo(1, 2, :, 3) + 1)
1075 fft_scratch_size%mcz1 = mcz1
1076 fft_scratch_size%mcx2 = mcx2
1077 fft_scratch_size%mcz2 = mcz2
1078 fft_scratch_size%mcy3 = mcy3
1079 fft_scratch_size%rs_group = group
1080 fft_scratch_size%g_pos = my_pos
1081 fft_scratch_size%numtask = dim(1)*dim(2)
1082
1083 IF (dim(1) > 1 .AND. dim(2) > 1) THEN
1084
1085 !
1086 ! First case; two stages of communication
1087 !
1088
1089 CALL get_fft_scratch(fft_scratch, tf_type=100, n=n, fft_sizes=fft_scratch_size)
1090
1091 IF (sign == fwfft) THEN
1092 ! Stage 1 -> 3
1093
1094 bbuf => fft_scratch%a2buf
1095
1096 IF (test) THEN
1097 sum_data = abs(sum(zin))
1098 CALL group%sum(sum_data)
1099 IF (my_pos == 0 .AND. output_unit > 0) THEN
1100 WRITE (output_unit, '(A)') " Two step communication algorithm "
1101 WRITE (output_unit, '(A,T67,2I7)') " Transform Z ", n(3), mx1*my1
1102 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(1) ", sum_data
1103 END IF
1104 END IF
1105
1106 ! FFT along z
1107 CALL fft_1d(fft_scratch%fft_plan(1), zin, bbuf, norm, stat)
1108
1109 abuf => fft_scratch%a3buf
1110
1111 IF (test) THEN
1112 sum_data = abs(sum(bbuf))
1113 CALL group%sum(sum_data)
1114 IF (my_pos == 0 .AND. output_unit > 0) THEN
1115 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(2) T", sum_data
1116 END IF
1117 END IF
1118
1119 CALL cube_transpose_2(bbuf, bo(:, :, :, 1), bo(:, :, :, 2), abuf, fft_scratch)
1120
1121 bbuf => fft_scratch%a4buf
1122
1123 IF (test) THEN
1124 sum_data = abs(sum(abuf))
1125 CALL group%sum(sum_data)
1126 IF (my_pos == 0 .AND. output_unit > 0) THEN
1127 WRITE (output_unit, '(A,T67,2I7)') " Transform Y ", n(2), mx2*mz2
1128 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(3) ", sum_data
1129 END IF
1130 END IF
1131
1132 ! FFT along y
1133 CALL fft_1d(fft_scratch%fft_plan(2), abuf, bbuf, 1.0_dp, stat)
1134
1135 abuf => fft_scratch%a5buf
1136
1137 IF (test) THEN
1138 sum_data = abs(sum(bbuf))
1139 CALL group%sum(sum_data)
1140 IF (my_pos == 0 .AND. output_unit > 0) THEN
1141 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(4) T", sum_data
1142 END IF
1143 END IF
1144
1145 CALL cube_transpose_4(bbuf, bo(:, :, :, 2), bo(:, :, :, 3), abuf, fft_scratch)
1146
1147 IF (test) THEN
1148 sum_data = abs(sum(abuf))
1149 CALL group%sum(sum_data)
1150 IF (my_pos == 0 .AND. output_unit > 0) THEN
1151 WRITE (output_unit, '(A,T67,2I7)') " Transform X ", n(1), my3*mz3
1152 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(5) ", sum_data
1153 END IF
1154 END IF
1155
1156 ! FFT along x
1157 CALL fft_1d(fft_scratch%fft_plan(3), abuf, gin, 1.0_dp, stat)
1158
1159 IF (test) THEN
1160 sum_data = abs(sum(gin))
1161 CALL group%sum(sum_data)
1162 IF (my_pos == 0 .AND. output_unit > 0) THEN
1163 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(6) ", sum_data
1164 END IF
1165 END IF
1166
1167 ELSE IF (sign == bwfft) THEN
1168 ! Stage 3 -> 1
1169
1170 bbuf => fft_scratch%a5buf
1171
1172 IF (test) THEN
1173 sum_data = abs(sum(gin))
1174 CALL group%sum(sum_data)
1175 IF (my_pos == 0 .AND. output_unit > 0) THEN
1176 WRITE (output_unit, '(A)') " Two step communication algorithm "
1177 WRITE (output_unit, '(A,T67,2I7)') " Transform X ", n(1), my3*mz3
1178 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(1) ", sum_data
1179 END IF
1180 END IF
1181
1182 ! FFT along x
1183 CALL fft_1d(fft_scratch%fft_plan(4), gin, bbuf, 1.0_dp, stat)
1184
1185 abuf => fft_scratch%a4buf
1186
1187 IF (test) THEN
1188 sum_data = abs(sum(bbuf))
1189 CALL group%sum(sum_data)
1190 IF (my_pos == 0 .AND. output_unit > 0) THEN
1191 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(2) T", sum_data
1192 END IF
1193 END IF
1194
1195 CALL cube_transpose_3(bbuf, bo(:, :, :, 3), bo(:, :, :, 2), abuf, fft_scratch)
1196
1197 bbuf => fft_scratch%a3buf
1198
1199 IF (test) THEN
1200 sum_data = abs(sum(abuf))
1201 CALL group%sum(sum_data)
1202 IF (my_pos == 0 .AND. output_unit > 0) THEN
1203 WRITE (output_unit, '(A,T67,2I7)') " Transform Y ", n(2), mx2*mz2
1204 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(3) ", sum_data
1205 END IF
1206 END IF
1207
1208 ! FFT along y
1209 CALL fft_1d(fft_scratch%fft_plan(5), abuf, bbuf, 1.0_dp, stat)
1210
1211 abuf => fft_scratch%a2buf
1212
1213 IF (test) THEN
1214 sum_data = abs(sum(bbuf))
1215 CALL group%sum(sum_data)
1216 IF (my_pos == 0 .AND. output_unit > 0) THEN
1217 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(4) T", sum_data
1218 END IF
1219 END IF
1220
1221 CALL cube_transpose_1(bbuf, bo(:, :, :, 2), bo(:, :, :, 1), abuf, fft_scratch)
1222
1223 IF (test) THEN
1224 sum_data = abs(sum(abuf))
1225 CALL group%sum(sum_data)
1226 IF (my_pos == 0 .AND. output_unit > 0) THEN
1227 WRITE (output_unit, '(A,T67,2I7)') " Transform Z ", n(3), mx1*my1
1228 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(5) ", sum_data
1229 END IF
1230 END IF
1231
1232 ! FFT along z
1233 CALL fft_1d(fft_scratch%fft_plan(6), abuf, zin, norm, stat)
1234
1235 IF (test) THEN
1236 sum_data = abs(sum(zin))
1237 CALL group%sum(sum_data)
1238 IF (my_pos == 0 .AND. output_unit > 0) THEN
1239 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(6) ", sum_data
1240 END IF
1241 END IF
1242
1243 ELSE
1244 cpabort("Illegal fsign parameter.")
1245 END IF
1246
1247 CALL release_fft_scratch(fft_scratch)
1248
1249 ELSE IF (dim(2) == 1) THEN
1250
1251 !
1252 ! Second case; one stage of communication
1253 !
1254
1255 CALL get_fft_scratch(fft_scratch, tf_type=101, n=n, fft_sizes=fft_scratch_size)
1256
1257 IF (sign == fwfft) THEN
1258 ! Stage 1 -> 3
1259
1260 IF (test) THEN
1261 sum_data = abs(sum(zin))
1262 CALL group%sum(sum_data)
1263 IF (my_pos == 0 .AND. output_unit > 0) THEN
1264 WRITE (output_unit, '(A)') " one step communication algorithm "
1265 WRITE (output_unit, '(A,T67,2I7)') " Transform Z ", n(3), mx1*my1
1266 WRITE (output_unit, '(A,T67,2I7)') " Transform Y ", n(2), mx1*mz1
1267 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(1) ", sum_data
1268 END IF
1269 END IF
1270
1271 abuf => fft_scratch%a3buf
1272 bbuf => fft_scratch%a4buf
1273 ! FFT along z and y
1274 CALL fft_1d(fft_scratch%fft_plan(1), zin, abuf, norm, stat)
1275 CALL fft_1d(fft_scratch%fft_plan(2), abuf, bbuf, 1.0_dp, stat)
1276
1277 abuf => fft_scratch%a5buf
1278
1279 IF (test) THEN
1280 sum_data = abs(sum(bbuf))
1281 CALL group%sum(sum_data)
1282 IF (my_pos == 0 .AND. output_unit > 0) THEN
1283 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(2) T", sum_data
1284 END IF
1285 END IF
1286
1287 CALL cube_transpose_6(bbuf, group, bo(:, :, :, 1), bo(:, :, :, 3), abuf, fft_scratch)
1288
1289 IF (test) THEN
1290 sum_data = abs(sum(abuf))
1291 CALL group%sum(sum_data)
1292 IF (my_pos == 0 .AND. output_unit > 0) THEN
1293 WRITE (output_unit, '(A,T67,2I7)') " Transform X ", n(1), my3*mz3
1294 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(3) ", sum_data
1295 END IF
1296 END IF
1297
1298 ! FFT along x
1299 CALL fft_1d(fft_scratch%fft_plan(3), abuf, gin, 1.0_dp, stat)
1300
1301 IF (test) THEN
1302 sum_data = abs(sum(gin))
1303 CALL group%sum(sum_data)
1304 IF (my_pos == 0 .AND. output_unit > 0) THEN
1305 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(4) ", sum_data
1306 END IF
1307 END IF
1308
1309 ELSE IF (sign == bwfft) THEN
1310 ! Stage 3 -> 1
1311
1312 IF (test) THEN
1313 sum_data = abs(sum(gin))
1314 CALL group%sum(sum_data)
1315 IF (my_pos == 0 .AND. output_unit > 0) THEN
1316 WRITE (output_unit, '(A)') " one step communication algorithm "
1317 WRITE (output_unit, '(A,T67,2I7)') " Transform X ", n(1), my3*mz3
1318 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(1) ", sum_data
1319 END IF
1320 END IF
1321
1322 bbuf => fft_scratch%a5buf
1323
1324 ! FFT along x
1325 CALL fft_1d(fft_scratch%fft_plan(4), gin, bbuf, 1.0_dp, stat)
1326
1327 abuf => fft_scratch%a4buf
1328
1329 IF (test) THEN
1330 sum_data = abs(sum(bbuf))
1331 CALL group%sum(sum_data)
1332 IF (my_pos == 0 .AND. output_unit > 0) THEN
1333 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(2) T", sum_data
1334 END IF
1335 END IF
1336
1337 CALL cube_transpose_5(bbuf, group, bo(:, :, :, 3), bo(:, :, :, 1), abuf, fft_scratch)
1338
1339 bbuf => fft_scratch%a3buf
1340
1341 IF (test) THEN
1342 sum_data = abs(sum(abuf))
1343 CALL group%sum(sum_data)
1344 IF (my_pos == 0 .AND. output_unit > 0) THEN
1345 WRITE (output_unit, '(A,T67,2I7)') " Transform Y ", n(2), mx1*mz1
1346 WRITE (output_unit, '(A,T67,2I7)') " Transform Z ", n(3), mx1*my1
1347 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(3) ", sum_data
1348 END IF
1349 END IF
1350
1351 ! FFT along y
1352 CALL fft_1d(fft_scratch%fft_plan(5), abuf, bbuf, 1.0_dp, stat)
1353
1354 ! FFT along z
1355 CALL fft_1d(fft_scratch%fft_plan(6), bbuf, zin, norm, stat)
1356
1357 IF (test) THEN
1358 sum_data = abs(sum(zin))
1359 CALL group%sum(sum_data)
1360 IF (my_pos == 0 .AND. output_unit > 0) THEN
1361 WRITE (output_unit, '(A,T61,E20.14)') " Sum of data(4) ", sum_data
1362 END IF
1363 END IF
1364
1365 ELSE
1366 cpabort("Illegal fsign parameter.")
1367 END IF
1368
1369 CALL release_fft_scratch(fft_scratch)
1370
1371 ELSE
1372
1373 cpabort("Partition not implemented.")
1374
1375 END IF
1376
1377 IF (PRESENT(status)) THEN
1378 status = stat
1379 END IF
1380
1381 CALL timestop(handle)
1382
1383 END SUBROUTINE fft3d_pb
1384
1385! **************************************************************************************************
1386!> \brief ...
1387!> \param sb ...
1388!> \param group ...
1389!> \param my_pos ...
1390!> \param p2p ...
1391!> \param yzp ...
1392!> \param nray ...
1393!> \param bo ...
1394!> \param tb ...
1395!> \param fft_scratch ...
1396!> \par History
1397!> 15. Feb. 2006 : single precision all_to_all
1398!> \author JGH (14-Jan-2001)
1399! **************************************************************************************************
1400 SUBROUTINE x_to_yz(sb, group, my_pos, p2p, yzp, nray, bo, tb, fft_scratch)
1401
1402 COMPLEX(KIND=dp), CONTIGUOUS, DIMENSION(:, :), &
1403 INTENT(IN) :: sb
1404
1405 CLASS(mp_comm_type), INTENT(IN) :: group
1406 INTEGER, INTENT(IN) :: my_pos
1407 INTEGER, CONTIGUOUS, DIMENSION(0:), INTENT(IN) :: p2p
1408 INTEGER, CONTIGUOUS, DIMENSION(:, :, 0:), &
1409 INTENT(IN) :: yzp
1410 INTEGER, CONTIGUOUS, DIMENSION(0:), INTENT(IN) :: nray
1411 INTEGER, CONTIGUOUS, DIMENSION(:, :, 0:), &
1412 INTENT(IN) :: bo
1413 COMPLEX(KIND=dp), CONTIGUOUS, DIMENSION(:, :, :), &
1414 INTENT(INOUT) :: tb
1415 TYPE(fft_scratch_type), INTENT(IN) :: fft_scratch
1416
1417 CHARACTER(len=*), PARAMETER :: routinen = 'x_to_yz'
1418
1419 COMPLEX(KIND=dp), CONTIGUOUS, DIMENSION(:, :), &
1420 POINTER :: rr
1421 COMPLEX(KIND=sp), CONTIGUOUS, DIMENSION(:, :), &
1422 POINTER :: ss, tt
1423 INTEGER :: handle, ip, ir, ix, ixx, iy, iz, mpr, &
1424 nm, np, nr, nx
1425 INTEGER, CONTIGUOUS, DIMENSION(:), POINTER :: rcount, rdispl, scount, sdispl
1426
1427 CALL timeset(routinen, handle)
1428
1429 np = SIZE(p2p)
1430 scount => fft_scratch%scount
1431 rcount => fft_scratch%rcount
1432 sdispl => fft_scratch%sdispl
1433 rdispl => fft_scratch%rdispl
1434
1435 IF (alltoall_sgl) THEN
1436 ss => fft_scratch%ss
1437 tt => fft_scratch%tt
1438 ss(:, :) = cmplx(sb(:, :), kind=sp)
1439 tt(:, :) = 0._sp
1440 ELSE
1441 rr => fft_scratch%rr
1442 END IF
1443
1444 mpr = p2p(my_pos)
1445 nm = maxval(nray(0:np - 1))
1446 nr = nray(my_pos)
1447!$OMP PARALLEL DO DEFAULT(NONE), &
1448!$OMP PRIVATE(ix,nx), &
1449!$OMP SHARED(np,p2p,bo,nr,scount,sdispl)
1450 DO ip = 0, np - 1
1451 ix = p2p(ip)
1452 nx = bo(2, 1, ix) - bo(1, 1, ix) + 1
1453 scount(ip) = nr*nx
1454 sdispl(ip) = nr*(bo(1, 1, ix) - 1)
1455 END DO
1456!$OMP END PARALLEL DO
1457 nx = bo(2, 1, mpr) - bo(1, 1, mpr) + 1
1458!$OMP PARALLEL DO DEFAULT(NONE), &
1459!$OMP PRIVATE(nr), &
1460!$OMP SHARED(np,nray,nx,rcount,rdispl,nm)
1461 DO ip = 0, np - 1
1462 nr = nray(ip)
1463 rcount(ip) = nr*nx
1464 rdispl(ip) = nm*nx*ip
1465 END DO
1466!$OMP END PARALLEL DO
1467 IF (alltoall_sgl) THEN
1468 CALL group%alltoall(ss, scount, sdispl, tt, rcount, rdispl)
1469 ELSE
1470 CALL group%alltoall(sb, scount, sdispl, rr, rcount, rdispl)
1471 END IF
1472
1473 nx = bo(2, 1, mpr) - bo(1, 1, mpr) + 1
1474!$OMP PARALLEL DO DEFAULT(NONE) COLLAPSE(2) &
1475!$OMP PRIVATE(ixx,ir,iy,iz,ix) &
1476!$OMP SHARED(np,nray,nx,alltoall_sgl,yzp,tt,rr,tb)
1477 DO ip = 0, np - 1
1478 DO ix = 1, nx
1479 ixx = nray(ip)*(ix - 1)
1480 IF (alltoall_sgl) THEN
1481 DO ir = 1, nray(ip)
1482 iy = yzp(1, ir, ip)
1483 iz = yzp(2, ir, ip)
1484 tb(iy, iz, ix) = tt(ir + ixx, ip)
1485 END DO
1486 ELSE
1487 DO ir = 1, nray(ip)
1488 iy = yzp(1, ir, ip)
1489 iz = yzp(2, ir, ip)
1490 tb(iy, iz, ix) = rr(ir + ixx, ip)
1491 END DO
1492 END IF
1493 END DO
1494 END DO
1495!$OMP END PARALLEL DO
1496
1497 CALL timestop(handle)
1498
1499 END SUBROUTINE x_to_yz
1500
1501! **************************************************************************************************
1502!> \brief ...
1503!> \param tb ...
1504!> \param group ...
1505!> \param my_pos ...
1506!> \param p2p ...
1507!> \param yzp ...
1508!> \param nray ...
1509!> \param bo ...
1510!> \param sb ...
1511!> \param fft_scratch ...
1512!> \par History
1513!> 15. Feb. 2006 : single precision all_to_all
1514!> \author JGH (14-Jan-2001)
1515! **************************************************************************************************
1516 SUBROUTINE yz_to_x(tb, group, my_pos, p2p, yzp, nray, bo, sb, fft_scratch)
1517
1518 COMPLEX(KIND=dp), CONTIGUOUS, DIMENSION(:, :, :), &
1519 INTENT(IN) :: tb
1520
1521 CLASS(mp_comm_type), INTENT(IN) :: group
1522 INTEGER, INTENT(IN) :: my_pos
1523 INTEGER, CONTIGUOUS, DIMENSION(0:), INTENT(IN) :: p2p
1524 INTEGER, CONTIGUOUS, DIMENSION(:, :, 0:), &
1525 INTENT(IN) :: yzp
1526 INTEGER, CONTIGUOUS, DIMENSION(0:), INTENT(IN) :: nray
1527 INTEGER, DIMENSION(:, :, 0:), INTENT(IN) :: bo
1528 COMPLEX(KIND=dp), CONTIGUOUS, DIMENSION(:, :), &
1529 INTENT(INOUT) :: sb
1530 TYPE(fft_scratch_type), INTENT(IN) :: fft_scratch
1531
1532 CHARACTER(len=*), PARAMETER :: routinen = 'yz_to_x'
1533
1534 COMPLEX(KIND=dp), CONTIGUOUS, DIMENSION(:, :), &
1535 POINTER :: rr
1536 COMPLEX(KIND=sp), CONTIGUOUS, DIMENSION(:, :), &
1537 POINTER :: ss, tt
1538 INTEGER :: handle, ip, ir, ix, ixx, iy, iz, mpr, &
1539 nm, np, nr, nx
1540 INTEGER, CONTIGUOUS, DIMENSION(:), POINTER :: rcount, rdispl, scount, sdispl
1541
1542 CALL timeset(routinen, handle)
1543
1544 np = SIZE(p2p)
1545 mpr = p2p(my_pos)
1546 scount => fft_scratch%scount
1547 rcount => fft_scratch%rcount
1548 sdispl => fft_scratch%sdispl
1549 rdispl => fft_scratch%rdispl
1550
1551 IF (alltoall_sgl) THEN
1552 ss => fft_scratch%ss
1553 tt => fft_scratch%tt
1554 ss = 0._sp
1555 tt = 0._sp
1556 ELSE
1557 rr => fft_scratch%rr
1558 END IF
1559
1560 nx = bo(2, 1, mpr) - bo(1, 1, mpr) + 1
1561!$OMP PARALLEL DO DEFAULT(NONE) COLLAPSE(2) &
1562!$OMP PRIVATE(ip, ixx, ir, iy, iz, ix) &
1563!$OMP SHARED(np,nray,nx,alltoall_sgl,yzp,tb,tt,rr)
1564 DO ip = 0, np - 1
1565 DO ix = 1, nx
1566 ixx = nray(ip)*(ix - 1)
1567 IF (alltoall_sgl) THEN
1568 DO ir = 1, nray(ip)
1569 iy = yzp(1, ir, ip)
1570 iz = yzp(2, ir, ip)
1571 tt(ir + ixx, ip) = cmplx(tb(iy, iz, ix), kind=sp)
1572 END DO
1573 ELSE
1574 DO ir = 1, nray(ip)
1575 iy = yzp(1, ir, ip)
1576 iz = yzp(2, ir, ip)
1577 rr(ir + ixx, ip) = tb(iy, iz, ix)
1578 END DO
1579 END IF
1580 END DO
1581 END DO
1582!$OMP END PARALLEL DO
1583 nm = maxval(nray(0:np - 1))
1584 nr = nray(my_pos)
1585!$OMP PARALLEL DO DEFAULT(NONE), &
1586!$OMP PRIVATE(ix,nx), &
1587!$OMP SHARED(np,p2p,bo,rcount,rdispl,nr)
1588 DO ip = 0, np - 1
1589 ix = p2p(ip)
1590 nx = bo(2, 1, ix) - bo(1, 1, ix) + 1
1591 rcount(ip) = nr*nx
1592 rdispl(ip) = nr*(bo(1, 1, ix) - 1)
1593 END DO
1594!$OMP END PARALLEL DO
1595 nx = bo(2, 1, mpr) - bo(1, 1, mpr) + 1
1596!$OMP PARALLEL DO DEFAULT(NONE), &
1597!$OMP PRIVATE(nr), &
1598!$OMP SHARED(np,nray,scount,sdispl,nx,nm)
1599 DO ip = 0, np - 1
1600 nr = nray(ip)
1601 scount(ip) = nr*nx
1602 sdispl(ip) = nm*nx*ip
1603 END DO
1604!$OMP END PARALLEL DO
1605
1606 IF (alltoall_sgl) THEN
1607 CALL group%alltoall(tt, scount, sdispl, ss, rcount, rdispl)
1608 sb = ss
1609 ELSE
1610 CALL group%alltoall(rr, scount, sdispl, sb, rcount, rdispl)
1611 END IF
1612
1613 CALL timestop(handle)
1614
1615 END SUBROUTINE yz_to_x
1616
1617! **************************************************************************************************
1618!> \brief ...
1619!> \param sb ...
1620!> \param group ...
1621!> \param dims ...
1622!> \param my_pos ...
1623!> \param p2p ...
1624!> \param yzp ...
1625!> \param nray ...
1626!> \param bo ...
1627!> \param tb ...
1628!> \param fft_scratch ...
1629!> \par History
1630!> 15. Feb. 2006 : single precision all_to_all
1631!> \author JGH (18-Jan-2001)
1632! **************************************************************************************************
1633 SUBROUTINE yz_to_xz(sb, group, dims, my_pos, p2p, yzp, nray, bo, tb, fft_scratch)
1634
1635 COMPLEX(KIND=dp), CONTIGUOUS, DIMENSION(:, :), &
1636 INTENT(IN) :: sb
1637
1638 CLASS(mp_comm_type), INTENT(IN) :: group
1639 INTEGER, DIMENSION(2), INTENT(IN) :: dims
1640 INTEGER, INTENT(IN) :: my_pos
1641 INTEGER, CONTIGUOUS, DIMENSION(0:), INTENT(IN) :: p2p
1642 INTEGER, CONTIGUOUS, DIMENSION(:, :, 0:), INTENT(IN) :: yzp
1643 INTEGER, CONTIGUOUS, DIMENSION(0:), INTENT(IN) :: nray
1644 INTEGER, CONTIGUOUS, DIMENSION(:, :, 0:), INTENT(IN) :: bo
1645 COMPLEX(KIND=dp), DIMENSION(:, :), INTENT(INOUT), CONTIGUOUS :: tb
1646 TYPE(fft_scratch_type), INTENT(INOUT) :: fft_scratch
1647
1648 CHARACTER(len=*), PARAMETER :: routinen = 'yz_to_xz'
1649
1650 COMPLEX(KIND=dp), DIMENSION(:), POINTER, CONTIGUOUS :: xzbuf, yzbuf
1651 COMPLEX(KIND=sp), DIMENSION(:), POINTER, CONTIGUOUS :: xzbuf_sgl, yzbuf_sgl
1652 INTEGER :: handle, icrs, ip, ipl, ipr, ir, ix, iz, &
1653 jj, jx, jy, jz, myx, myz, np, npx, &
1654 npz, nx, nz, rs_pos
1655 INTEGER, DIMENSION(:), POINTER, CONTIGUOUS :: pzcoord, rcount, rdispl, scount, sdispl, &
1656 xcor, zcor
1657 INTEGER, DIMENSION(:, :), CONTIGUOUS, POINTER :: pgrid
1658
1659 CALL timeset(routinen, handle)
1660
1661 np = SIZE(p2p)
1662
1663 rs_pos = p2p(my_pos)
1664
1665 IF (alltoall_sgl) THEN
1666 yzbuf_sgl => fft_scratch%yzbuf_sgl
1667 xzbuf_sgl => fft_scratch%xzbuf_sgl
1668 ELSE
1669 yzbuf => fft_scratch%yzbuf
1670 xzbuf => fft_scratch%xzbuf
1671 END IF
1672 npx = dims(1)
1673 npz = dims(2)
1674 pgrid => fft_scratch%pgrid
1675 xcor => fft_scratch%xcor
1676 zcor => fft_scratch%zcor
1677 pzcoord => fft_scratch%pzcoord
1678 scount => fft_scratch%scount
1679 rcount => fft_scratch%rcount
1680 sdispl => fft_scratch%sdispl
1681 rdispl => fft_scratch%rdispl
1682
1683 nx = SIZE(sb, 2)
1684
1685! If the send and recv counts are not already cached, then
1686! calculate and store them
1687 IF (fft_scratch%in == 0) THEN
1688
1689 scount = 0
1690
1691 DO ix = 0, npx - 1
1692 ip = pgrid(ix, 0)
1693 xcor(bo(1, 1, ip):bo(2, 1, ip)) = ix
1694 END DO
1695 DO iz = 0, npz - 1
1696 ip = pgrid(0, iz)
1697 zcor(bo(1, 3, ip):bo(2, 3, ip)) = iz
1698 END DO
1699 DO jx = 1, nx
1700 IF (alltoall_sgl) THEN
1701 DO ir = 1, nray(my_pos)
1702 jy = yzp(1, ir, my_pos)
1703 jz = yzp(2, ir, my_pos)
1704 ip = pgrid(xcor(jx), zcor(jz))
1705 scount(ip) = scount(ip) + 1
1706 END DO
1707 ELSE
1708 DO ir = 1, nray(my_pos)
1709 jy = yzp(1, ir, my_pos)
1710 jz = yzp(2, ir, my_pos)
1711 ip = pgrid(xcor(jx), zcor(jz))
1712 scount(ip) = scount(ip) + 1
1713 END DO
1714 END IF
1715 END DO
1716
1717 CALL group%alltoall(scount, rcount, 1)
1718 fft_scratch%yzcount = scount
1719 fft_scratch%xzcount = rcount
1720
1721 ! Work out the correct displacements in the buffers
1722 sdispl(0) = 0
1723 rdispl(0) = 0
1724 DO ip = 1, np - 1
1725 sdispl(ip) = sdispl(ip - 1) + scount(ip - 1)
1726 rdispl(ip) = rdispl(ip - 1) + rcount(ip - 1)
1727 END DO
1728
1729 fft_scratch%yzdispl = sdispl
1730 fft_scratch%xzdispl = rdispl
1731
1732 icrs = 0
1733 DO ip = 0, np - 1
1734 IF (scount(ip) /= 0) icrs = icrs + 1
1735 IF (rcount(ip) /= 0) icrs = icrs + 1
1736 END DO
1737 CALL group%sum(icrs)
1738 fft_scratch%rsratio = real(icrs, kind=dp)/(real(2*np, kind=dp)*real(np, kind=dp))
1739
1740 fft_scratch%in = 1
1741 ELSE
1742 scount = fft_scratch%yzcount
1743 rcount = fft_scratch%xzcount
1744 sdispl = fft_scratch%yzdispl
1745 rdispl = fft_scratch%xzdispl
1746 END IF
1747
1748! Do the actual packing
1749!$OMP PARALLEL DO DEFAULT(NONE), &
1750!$OMP PRIVATE(ipl,jj,nx,ir,jx,jy,jz),&
1751!$OMP SHARED(np,p2p,pzcoord,bo,nray,yzp,zcor),&
1752!$OMP SHARED(yzbuf,sb,scount,sdispl,my_pos),&
1753!$OMP SHARED(yzbuf_sgl,alltoall_sgl)
1754 DO ip = 0, np - 1
1755 IF (scount(ip) == 0) cycle
1756 ipl = p2p(ip)
1757 jj = 0
1758 nx = bo(2, 1, ipl) - bo(1, 1, ipl) + 1
1759 DO ir = 1, nray(my_pos)
1760 jz = yzp(2, ir, my_pos)
1761 IF (zcor(jz) == pzcoord(ipl)) THEN
1762 jj = jj + 1
1763 jy = yzp(1, ir, my_pos)
1764 IF (alltoall_sgl) THEN
1765 DO jx = 0, nx - 1
1766 yzbuf_sgl(sdispl(ip) + jj + jx*scount(ip)/nx) = cmplx(sb(ir, jx + bo(1, 1, ipl)), kind=sp)
1767 END DO
1768 ELSE
1769 DO jx = 0, nx - 1
1770 yzbuf(sdispl(ip) + jj + jx*scount(ip)/nx) = sb(ir, jx + bo(1, 1, ipl))
1771 END DO
1772 END IF
1773 END IF
1774 END DO
1775 END DO
1776!$OMP END PARALLEL DO
1777
1778 IF (alltoall_sgl) THEN
1779 CALL group%alltoall(yzbuf_sgl, scount, sdispl, xzbuf_sgl, rcount, rdispl)
1780 ELSE
1781 IF (fft_scratch%rsratio < ratio_sparse_alltoall) THEN
1782 CALL sparse_alltoall(yzbuf, scount, sdispl, xzbuf, rcount, rdispl, group)
1783 ELSE
1784 CALL group%alltoall(yzbuf, scount, sdispl, xzbuf, rcount, rdispl)
1785 END IF
1786 END IF
1787
1788 myx = fft_scratch%sizes%r_pos(1)
1789 myz = fft_scratch%sizes%r_pos(2)
1790 nz = bo(2, 3, rs_pos) - bo(1, 3, rs_pos) + 1
1791
1792!$OMP PARALLEL DO DEFAULT(NONE), &
1793!$OMP PRIVATE(ipr,jj,ir,jx,jy,jz),&
1794!$OMP SHARED(tb,np,p2p,bo,rs_pos,nray),&
1795!$OMP SHARED(yzp,alltoall_sgl,zcor,myz),&
1796!$OMP SHARED(xzbuf,xzbuf_sgl,nz,rdispl)
1797 DO ip = 0, np - 1
1798 ipr = p2p(ip)
1799 jj = 0
1800 DO jx = 0, bo(2, 1, rs_pos) - bo(1, 1, rs_pos)
1801 DO ir = 1, nray(ip)
1802 jz = yzp(2, ir, ip)
1803 IF (alltoall_sgl) THEN
1804 IF (zcor(jz) == myz) THEN
1805 jj = jj + 1
1806 jy = yzp(1, ir, ip)
1807 jz = jz - bo(1, 3, rs_pos) + 1
1808 tb(jy, jz + jx*nz) = xzbuf_sgl(jj + rdispl(ipr))
1809 END IF
1810 ELSE
1811 IF (zcor(jz) == myz) THEN
1812 jj = jj + 1
1813 jy = yzp(1, ir, ip)
1814 jz = jz - bo(1, 3, rs_pos) + 1
1815 tb(jy, jz + jx*nz) = xzbuf(jj + rdispl(ipr))
1816 END IF
1817 END IF
1818 END DO
1819 END DO
1820 END DO
1821!$OMP END PARALLEL DO
1822
1823 CALL timestop(handle)
1824
1825 END SUBROUTINE yz_to_xz
1826
1827! **************************************************************************************************
1828!> \brief ...
1829!> \param sb ...
1830!> \param group ...
1831!> \param dims ...
1832!> \param my_pos ...
1833!> \param p2p ...
1834!> \param yzp ...
1835!> \param nray ...
1836!> \param bo ...
1837!> \param tb ...
1838!> \param fft_scratch ...
1839!> \par History
1840!> 15. Feb. 2006 : single precision all_to_all
1841!> \author JGH (19-Jan-2001)
1842! **************************************************************************************************
1843 SUBROUTINE xz_to_yz(sb, group, dims, my_pos, p2p, yzp, nray, bo, tb, fft_scratch)
1844
1845 COMPLEX(KIND=dp), CONTIGUOUS, DIMENSION(:, :), &
1846 INTENT(IN) :: sb
1847
1848 CLASS(mp_comm_type), INTENT(IN) :: group
1849 INTEGER, DIMENSION(2), INTENT(IN) :: dims
1850 INTEGER, INTENT(IN) :: my_pos
1851 INTEGER, CONTIGUOUS, DIMENSION(0:), INTENT(IN) :: p2p
1852 INTEGER, CONTIGUOUS, DIMENSION(:, :, 0:), INTENT(IN) :: yzp
1853 INTEGER, CONTIGUOUS, DIMENSION(0:), INTENT(IN) :: nray
1854 INTEGER, CONTIGUOUS, DIMENSION(:, :, 0:), INTENT(IN) :: bo
1855 COMPLEX(KIND=dp), DIMENSION(:, :), INTENT(INOUT), CONTIGUOUS :: tb
1856 TYPE(fft_scratch_type), INTENT(INOUT) :: fft_scratch
1857
1858 CHARACTER(len=*), PARAMETER :: routinen = 'xz_to_yz'
1859
1860 COMPLEX(KIND=dp), DIMENSION(:), POINTER, CONTIGUOUS :: xzbuf, yzbuf
1861 COMPLEX(KIND=sp), DIMENSION(:), POINTER, CONTIGUOUS :: xzbuf_sgl, yzbuf_sgl
1862 INTEGER :: handle, icrs, ip, ipl, ir, ix, ixx, iz, &
1863 jj, jx, jy, jz, mp, myx, myz, np, npx, &
1864 npz, nx, nz
1865 INTEGER, DIMENSION(:), POINTER, CONTIGUOUS :: pzcoord, rcount, rdispl, scount, sdispl, &
1866 xcor, zcor
1867 INTEGER, DIMENSION(:, :), CONTIGUOUS, POINTER :: pgrid
1868
1869 CALL timeset(routinen, handle)
1870
1871 np = SIZE(p2p)
1872
1873 IF (alltoall_sgl) THEN
1874 yzbuf_sgl => fft_scratch%yzbuf_sgl
1875 xzbuf_sgl => fft_scratch%xzbuf_sgl
1876 ELSE
1877 yzbuf => fft_scratch%yzbuf
1878 xzbuf => fft_scratch%xzbuf
1879 END IF
1880 npx = dims(1)
1881 npz = dims(2)
1882 pgrid => fft_scratch%pgrid
1883 xcor => fft_scratch%xcor
1884 zcor => fft_scratch%zcor
1885 pzcoord => fft_scratch%pzcoord
1886 scount => fft_scratch%scount
1887 rcount => fft_scratch%rcount
1888 sdispl => fft_scratch%sdispl
1889 rdispl => fft_scratch%rdispl
1890
1891! If the send and recv counts are not already cached, then
1892! calculate and store them
1893 IF (fft_scratch%in == 0) THEN
1894
1895 rcount = 0
1896 nx = maxval(bo(2, 1, :))
1897
1898 DO ix = 0, npx - 1
1899 ip = pgrid(ix, 0)
1900 xcor(bo(1, 1, ip):bo(2, 1, ip)) = ix
1901 END DO
1902 DO iz = 0, npz - 1
1903 ip = pgrid(0, iz)
1904 zcor(bo(1, 3, ip):bo(2, 3, ip)) = iz
1905 END DO
1906 DO jx = 1, nx
1907 DO ir = 1, nray(my_pos)
1908 jy = yzp(1, ir, my_pos)
1909 jz = yzp(2, ir, my_pos)
1910 ip = pgrid(xcor(jx), zcor(jz))
1911 rcount(ip) = rcount(ip) + 1
1912 END DO
1913 END DO
1914
1915 CALL group%alltoall(rcount, scount, 1)
1916 fft_scratch%xzcount = scount
1917 fft_scratch%yzcount = rcount
1918
1919 ! Work out the correct displacements in the buffers
1920 sdispl(0) = 0
1921 rdispl(0) = 0
1922 DO ip = 1, np - 1
1923 sdispl(ip) = sdispl(ip - 1) + scount(ip - 1)
1924 rdispl(ip) = rdispl(ip - 1) + rcount(ip - 1)
1925 END DO
1926
1927 fft_scratch%xzdispl = sdispl
1928 fft_scratch%yzdispl = rdispl
1929
1930 icrs = 0
1931 DO ip = 0, np - 1
1932 IF (scount(ip) /= 0) icrs = icrs + 1
1933 IF (rcount(ip) /= 0) icrs = icrs + 1
1934 END DO
1935 CALL group%sum(icrs)
1936 fft_scratch%rsratio = real(icrs, kind=dp)/(real(2*np, kind=dp)*real(np, kind=dp))
1937
1938 fft_scratch%in = 1
1939 ELSE
1940 scount = fft_scratch%xzcount
1941 rcount = fft_scratch%yzcount
1942 sdispl = fft_scratch%xzdispl
1943 rdispl = fft_scratch%yzdispl
1944 END IF
1945
1946! Now do the actual packing
1947 myx = fft_scratch%sizes%r_pos(1)
1948 myz = fft_scratch%sizes%r_pos(2)
1949 mp = p2p(my_pos)
1950 nz = bo(2, 3, mp) - bo(1, 3, mp) + 1
1951 nx = bo(2, 1, mp) - bo(1, 1, mp) + 1
1952
1953!$OMP PARALLEL DO DEFAULT(NONE), &
1954!$OMP PRIVATE(jj,ipl,ir,jx,jy,jz,ixx),&
1955!$OMP SHARED(np,p2p,nray,yzp,zcor,myz,bo,mp),&
1956!$OMP SHARED(alltoall_sgl,nx,scount,sdispl),&
1957!$OMP SHARED(xzbuf,xzbuf_sgl,sb,nz)
1958 DO ip = 0, np - 1
1959 jj = 0
1960 ipl = p2p(ip)
1961 DO ir = 1, nray(ip)
1962 jz = yzp(2, ir, ip)
1963 IF (zcor(jz) == myz) THEN
1964 jj = jj + 1
1965 jy = yzp(1, ir, ip)
1966 jz = yzp(2, ir, ip) - bo(1, 3, mp) + 1
1967 IF (alltoall_sgl) THEN
1968 DO jx = 0, nx - 1
1969 ixx = jj + jx*scount(ipl)/nx
1970 xzbuf_sgl(ixx + sdispl(ipl)) = cmplx(sb(jy, jz + jx*nz), kind=sp)
1971 END DO
1972 ELSE
1973 DO jx = 0, nx - 1
1974 ixx = jj + jx*scount(ipl)/nx
1975 xzbuf(ixx + sdispl(ipl)) = sb(jy, jz + jx*nz)
1976 END DO
1977 END IF
1978 END IF
1979 END DO
1980 END DO
1981!$OMP END PARALLEL DO
1982
1983 IF (alltoall_sgl) THEN
1984 CALL group%alltoall(xzbuf_sgl, scount, sdispl, yzbuf_sgl, rcount, rdispl)
1985 ELSE
1986 IF (fft_scratch%rsratio < ratio_sparse_alltoall) THEN
1987 CALL sparse_alltoall(xzbuf, scount, sdispl, yzbuf, rcount, rdispl, group)
1988 ELSE
1989 CALL group%alltoall(xzbuf, scount, sdispl, yzbuf, rcount, rdispl)
1990 END IF
1991 END IF
1992
1993!$OMP PARALLEL DO DEFAULT(NONE), &
1994!$OMP PRIVATE(ipl,jj,nx,ir,jx,jy,jz),&
1995!$OMP SHARED(p2p,pzcoord,bo,nray,my_pos,yzp),&
1996!$OMP SHARED(rcount,rdispl,tb,yzbuf,zcor),&
1997!$OMP SHARED(yzbuf_sgl,alltoall_sgl,np)
1998 DO ip = 0, np - 1
1999 IF (rcount(ip) == 0) cycle
2000 ipl = p2p(ip)
2001 jj = 0
2002 nx = bo(2, 1, ipl) - bo(1, 1, ipl) + 1
2003 DO ir = 1, nray(my_pos)
2004 jz = yzp(2, ir, my_pos)
2005 IF (zcor(jz) == pzcoord(ipl)) THEN
2006 jj = jj + 1
2007 jy = yzp(1, ir, my_pos)
2008 IF (alltoall_sgl) THEN
2009 DO jx = 0, nx - 1
2010 tb(ir, jx + bo(1, 1, ipl)) = yzbuf_sgl(rdispl(ip) + jj + jx*rcount(ip)/nx)
2011 END DO
2012 ELSE
2013 DO jx = 0, nx - 1
2014 tb(ir, jx + bo(1, 1, ipl)) = yzbuf(rdispl(ip) + jj + jx*rcount(ip)/nx)
2015 END DO
2016 END IF
2017 END IF
2018 END DO
2019 END DO
2020!$OMP END PARALLEL DO
2021
2022 CALL timestop(handle)
2023
2024 END SUBROUTINE xz_to_yz
2025
2026! **************************************************************************************************
2027!> \brief ...
2028!> \param cin ...
2029!> \param boin ...
2030!> \param boout ...
2031!> \param sout ...
2032!> \param fft_scratch ...
2033!> \par History
2034!> none
2035!> \author JGH (20-Jan-2001)
2036! **************************************************************************************************
2037 SUBROUTINE cube_transpose_1(cin, boin, boout, sout, fft_scratch)
2038
2039 COMPLEX(KIND=dp), CONTIGUOUS, DIMENSION(:, :), &
2040 INTENT(IN) :: cin
2041 INTEGER, CONTIGUOUS, DIMENSION(:, :, 0:), &
2042 INTENT(IN) :: boin, boout
2043 COMPLEX(KIND=dp), CONTIGUOUS, DIMENSION(:, :), &
2044 INTENT(OUT) :: sout
2045 TYPE(fft_scratch_type), INTENT(IN) :: fft_scratch
2046
2047 CHARACTER(len=*), PARAMETER :: routinen = 'cube_transpose_1'
2048
2049 COMPLEX(KIND=dp), CONTIGUOUS, DIMENSION(:, :), &
2050 POINTER :: rbuf
2051 INTEGER :: handle, ip, ipl, ir, is, ixy, iz, mip, &
2052 mz, np, nx, ny, nz
2053 INTEGER, CONTIGUOUS, DIMENSION(:), POINTER :: rcount, rdispl, scount, sdispl
2054 INTEGER, CONTIGUOUS, DIMENSION(:, :), POINTER :: pgrid
2055 INTEGER, DIMENSION(2) :: dim, pos
2056
2057 CALL timeset(routinen, handle)
2058
2059 mip = fft_scratch%mip
2060 dim = fft_scratch%dim
2061 pos = fft_scratch%pos
2062 scount => fft_scratch%scount
2063 rcount => fft_scratch%rcount
2064 sdispl => fft_scratch%sdispl
2065 rdispl => fft_scratch%rdispl
2066 pgrid => fft_scratch%pgcube
2067 np = dim(2)
2068
2069 nx = boin(2, 1, mip) - boin(1, 1, mip) + 1
2070 nz = boin(2, 3, mip) - boin(1, 3, mip) + 1
2071
2072!$OMP PARALLEL DO DEFAULT(NONE), &
2073!$OMP PRIVATE(ipl,ny), &
2074!$OMP SHARED(np,pgrid,boout,scount,sdispl,nx,nz)
2075 DO ip = 0, np - 1
2076 ipl = pgrid(ip, 2)
2077 ny = boout(2, 2, ipl) - boout(1, 2, ipl) + 1
2078 scount(ip) = nx*nz*ny
2079 sdispl(ip) = nx*nz*(boout(1, 2, ipl) - 1)
2080 END DO
2081!$OMP END PARALLEL DO
2082 ny = boout(2, 2, mip) - boout(1, 2, mip) + 1
2083 mz = maxval(boin(2, 3, :) - boin(1, 3, :) + 1)
2084!$OMP PARALLEL DO DEFAULT(NONE), &
2085!$OMP PRIVATE(ipl,nz), &
2086!$OMP SHARED(np,pgrid,boin,nx,ny,rcount,rdispl,mz)
2087 DO ip = 0, np - 1
2088 ipl = pgrid(ip, 2)
2089 nz = boin(2, 3, ipl) - boin(1, 3, ipl) + 1
2090 rcount(ip) = nx*nz*ny
2091 rdispl(ip) = nx*ny*mz*ip
2092 END DO
2093!$OMP END PARALLEL DO
2094
2095 rbuf => fft_scratch%rbuf1
2096
2097 CALL fft_scratch%cart_sub_comm(2)%alltoall(cin, scount, sdispl, rbuf, rcount, rdispl)
2098
2099!$OMP PARALLEL DO DEFAULT(NONE) COLLAPSE(2) &
2100!$OMP PRIVATE(ip,ipl,nz,iz,is,ir) &
2101!$OMP SHARED(nx,ny,np,pgrid,boin,sout,rbuf)
2102 DO ixy = 1, nx*ny
2103 DO ip = 0, np - 1
2104 ipl = pgrid(ip, 2)
2105 nz = boin(2, 3, ipl) - boin(1, 3, ipl) + 1
2106 DO iz = 1, nz
2107 is = boin(1, 3, ipl) + iz - 1
2108 ir = iz + nz*(ixy - 1)
2109 sout(is, ixy) = rbuf(ir, ip)
2110 END DO
2111 END DO
2112 END DO
2113!$OMP END PARALLEL DO
2114
2115 CALL timestop(handle)
2116
2117 END SUBROUTINE cube_transpose_1
2118
2119! **************************************************************************************************
2120!> \brief ...
2121!> \param cin ...
2122!> \param boin ...
2123!> \param boout ...
2124!> \param sout ...
2125!> \param fft_scratch ...
2126! **************************************************************************************************
2127 SUBROUTINE cube_transpose_2(cin, boin, boout, sout, fft_scratch)
2128
2129 COMPLEX(KIND=dp), CONTIGUOUS, DIMENSION(:, :), &
2130 INTENT(IN) :: cin
2131 INTEGER, CONTIGUOUS, DIMENSION(:, :, 0:), &
2132 INTENT(IN) :: boin, boout
2133 COMPLEX(KIND=dp), CONTIGUOUS, DIMENSION(:, :), &
2134 INTENT(OUT) :: sout
2135 TYPE(fft_scratch_type), INTENT(IN) :: fft_scratch
2136
2137 CHARACTER(len=*), PARAMETER :: routinen = 'cube_transpose_2'
2138
2139 COMPLEX(KIND=dp), CONTIGUOUS, DIMENSION(:, :), &
2140 POINTER :: rbuf
2141 INTEGER :: handle, ip, ipl, ir, ixy, iz, mip, mz, &
2142 np, nx, ny, nz
2143 INTEGER, CONTIGUOUS, DIMENSION(:), POINTER :: rcount, rdispl, scount, sdispl
2144 INTEGER, CONTIGUOUS, DIMENSION(:, :), POINTER :: pgrid
2145 INTEGER, DIMENSION(2) :: dim, pos
2146 TYPE(mp_comm_type) :: sub_group
2147
2148 CALL timeset(routinen, handle)
2149
2150 sub_group = fft_scratch%cart_sub_comm(2)
2151 mip = fft_scratch%mip
2152 dim = fft_scratch%dim
2153 pos = fft_scratch%pos
2154 scount => fft_scratch%scount
2155 rcount => fft_scratch%rcount
2156 sdispl => fft_scratch%sdispl
2157 rdispl => fft_scratch%rdispl
2158 pgrid => fft_scratch%pgcube
2159 np = dim(2)
2160
2161 nx = boin(2, 1, mip) - boin(1, 1, mip) + 1
2162 ny = boin(2, 2, mip) - boin(1, 2, mip) + 1
2163 mz = maxval(boout(2, 3, :) - boout(1, 3, :) + 1)
2164
2165 rbuf => fft_scratch%rbuf2
2166
2167!$OMP PARALLEL DEFAULT(NONE), &
2168!$OMP PRIVATE(ip,ipl,nz,iz,ir), &
2169!$OMP SHARED(nx,ny,np,pgrid,boout,rbuf,cin,scount,sdispl,mz)
2170!$OMP DO COLLAPSE(2)
2171 DO ixy = 1, nx*ny
2172 DO ip = 0, np - 1
2173 ipl = pgrid(ip, 2)
2174 nz = boout(2, 3, ipl) - boout(1, 3, ipl) + 1
2175 DO iz = boout(1, 3, ipl), boout(2, 3, ipl)
2176 ir = iz - boout(1, 3, ipl) + 1 + (ixy - 1)*nz
2177 rbuf(ir, ip) = cin(iz, ixy)
2178 END DO
2179 END DO
2180 END DO
2181!$OMP END DO
2182!$OMP DO
2183 DO ip = 0, np - 1
2184 ipl = pgrid(ip, 2)
2185 nz = boout(2, 3, ipl) - boout(1, 3, ipl) + 1
2186 scount(ip) = nx*ny*nz
2187 sdispl(ip) = nx*ny*mz*ip
2188 END DO
2189!$OMP END DO
2190!$OMP END PARALLEL
2191 nz = boout(2, 3, mip) - boout(1, 3, mip) + 1
2192!$OMP PARALLEL DO DEFAULT(NONE), &
2193!$OMP PRIVATE(ipl,ny), &
2194!$OMP SHARED(np,pgrid,boin,nx,nz,rcount,rdispl)
2195 DO ip = 0, np - 1
2196 ipl = pgrid(ip, 2)
2197 ny = boin(2, 2, ipl) - boin(1, 2, ipl) + 1
2198 rcount(ip) = nx*ny*nz
2199 rdispl(ip) = nx*nz*(boin(1, 2, ipl) - 1)
2200 END DO
2201!$OMP END PARALLEL DO
2202
2203 CALL sub_group%alltoall(rbuf, scount, sdispl, sout, rcount, rdispl)
2204
2205 CALL timestop(handle)
2206
2207 END SUBROUTINE cube_transpose_2
2208
2209! **************************************************************************************************
2210!> \brief ...
2211!> \param cin ...
2212!> \param boin ...
2213!> \param boout ...
2214!> \param sout ...
2215!> \param fft_scratch ...
2216! **************************************************************************************************
2217 SUBROUTINE cube_transpose_3(cin, boin, boout, sout, fft_scratch)
2218
2219 COMPLEX(KIND=dp), CONTIGUOUS, DIMENSION(:, :), &
2220 INTENT(IN) :: cin
2221 INTEGER, CONTIGUOUS, DIMENSION(:, :, 0:), &
2222 INTENT(IN) :: boin, boout
2223 COMPLEX(KIND=dp), CONTIGUOUS, DIMENSION(:, :), &
2224 INTENT(OUT) :: sout
2225 TYPE(fft_scratch_type), INTENT(IN) :: fft_scratch
2226
2227 CHARACTER(len=*), PARAMETER :: routinen = 'cube_transpose_3'
2228
2229 COMPLEX(KIND=dp), CONTIGUOUS, DIMENSION(:, :), &
2230 POINTER :: rbuf
2231 INTEGER :: handle, ip, ipl, ir, is, ixz, iy, lb, &
2232 mip, my, my_id, np, num_threads, nx, &
2233 ny, nz, ub
2234 INTEGER, CONTIGUOUS, DIMENSION(:), POINTER :: rcount, rdispl, scount, sdispl
2235 INTEGER, CONTIGUOUS, DIMENSION(:, :), POINTER :: pgrid
2236 INTEGER, DIMENSION(2) :: dim, pos
2237 TYPE(mp_comm_type) :: sub_group
2238
2239 CALL timeset(routinen, handle)
2240
2241 sub_group = fft_scratch%cart_sub_comm(1)
2242 mip = fft_scratch%mip
2243 dim = fft_scratch%dim
2244 pos = fft_scratch%pos
2245 np = dim(1)
2246 scount => fft_scratch%scount
2247 rcount => fft_scratch%rcount
2248 sdispl => fft_scratch%sdispl
2249 rdispl => fft_scratch%rdispl
2250 pgrid => fft_scratch%pgcube
2251
2252 ny = boin(2, 2, mip) - boin(1, 2, mip) + 1
2253 nz = boin(2, 3, mip) - boin(1, 3, mip) + 1
2254!$OMP PARALLEL DO DEFAULT(NONE), &
2255!$OMP PRIVATE(ipl, nx), &
2256!$OMP SHARED(np,pgrid,boout,ny,nz,scount,sdispl)
2257 DO ip = 0, np - 1
2258 ipl = pgrid(ip, 1)
2259 nx = boout(2, 1, ipl) - boout(1, 1, ipl) + 1
2260 scount(ip) = nx*nz*ny
2261 sdispl(ip) = ny*nz*(boout(1, 1, ipl) - 1)
2262 END DO
2263!$OMP END PARALLEL DO
2264 nx = boout(2, 1, mip) - boout(1, 1, mip) + 1
2265 my = maxval(boin(2, 2, :) - boin(1, 2, :) + 1)
2266!$OMP PARALLEL DO DEFAULT(NONE), &
2267!$OMP PRIVATE(ipl, ny), &
2268!$OMP SHARED(np,pgrid,boin,nx,nz,my,rcount,rdispl)
2269 DO ip = 0, np - 1
2270 ipl = pgrid(ip, 1)
2271 ny = boin(2, 2, ipl) - boin(1, 2, ipl) + 1
2272 rcount(ip) = nx*nz*ny
2273 rdispl(ip) = nx*my*nz*ip
2274 END DO
2275!$OMP END PARALLEL DO
2276
2277 rbuf => fft_scratch%rbuf3
2278 num_threads = 1
2279 my_id = 0
2280!$OMP PARALLEL DEFAULT(NONE), &
2281!$OMP PRIVATE(NUM_THREADS, my_id, lb, ub) &
2282!$OMP SHARED(rbuf)
2283!$ num_threads = MIN(omp_get_max_threads(), SIZE(rbuf, 2))
2284!$ my_id = omp_get_thread_num()
2285 IF (my_id < num_threads) THEN
2286 lb = (SIZE(rbuf, 2)*my_id)/num_threads
2287 ub = (SIZE(rbuf, 2)*(my_id + 1))/num_threads - 1
2288 rbuf(:, lb:ub) = 0.0_dp
2289 END IF
2290!$OMP END PARALLEL
2291
2292 CALL sub_group%alltoall(cin, scount, sdispl, rbuf, rcount, rdispl)
2293
2294!$OMP PARALLEL DO DEFAULT(NONE) COLLAPSE(2) &
2295!$OMP PRIVATE(ip,ipl,ny,iy,is,ir) &
2296!$OMP SHARED(nx,nz,np,pgrid,boin,rbuf,sout)
2297 DO ixz = 1, nx*nz
2298 DO ip = 0, np - 1
2299 ipl = pgrid(ip, 1)
2300 ny = boin(2, 2, ipl) - boin(1, 2, ipl) + 1
2301 DO iy = 1, ny
2302 is = boin(1, 2, ipl) + iy - 1
2303 ir = iy + ny*(ixz - 1)
2304 sout(is, ixz) = rbuf(ir, ip)
2305 END DO
2306 END DO
2307 END DO
2308!$OMP END PARALLEL DO
2309
2310 CALL timestop(handle)
2311
2312 END SUBROUTINE cube_transpose_3
2313
2314! **************************************************************************************************
2315!> \brief ...
2316!> \param cin ...
2317!> \param boin ...
2318!> \param boout ...
2319!> \param sout ...
2320!> \param fft_scratch ...
2321! **************************************************************************************************
2322 SUBROUTINE cube_transpose_4(cin, boin, boout, sout, fft_scratch)
2323
2324 COMPLEX(KIND=dp), CONTIGUOUS, DIMENSION(:, :), &
2325 INTENT(IN) :: cin
2326 INTEGER, CONTIGUOUS, DIMENSION(:, :, 0:), &
2327 INTENT(IN) :: boin, boout
2328 COMPLEX(KIND=dp), CONTIGUOUS, DIMENSION(:, :), &
2329 INTENT(OUT) :: sout
2330 TYPE(fft_scratch_type), INTENT(IN) :: fft_scratch
2331
2332 CHARACTER(len=*), PARAMETER :: routinen = 'cube_transpose_4'
2333
2334 COMPLEX(KIND=dp), CONTIGUOUS, DIMENSION(:, :), &
2335 POINTER :: rbuf
2336 INTEGER :: handle, ip, ipl, ir, iy, izx, lb, mip, &
2337 my, my_id, np, num_threads, nx, ny, &
2338 nz, ub
2339 INTEGER, CONTIGUOUS, DIMENSION(:), POINTER :: rcount, rdispl, scount, sdispl
2340 INTEGER, CONTIGUOUS, DIMENSION(:, :), POINTER :: pgrid
2341 INTEGER, DIMENSION(2) :: dim, pos
2342 TYPE(mp_comm_type) :: sub_group
2343
2344 CALL timeset(routinen, handle)
2345
2346 sub_group = fft_scratch%cart_sub_comm(1)
2347 mip = fft_scratch%mip
2348 dim = fft_scratch%dim
2349 pos = fft_scratch%pos
2350 np = dim(1)
2351 scount => fft_scratch%scount
2352 rcount => fft_scratch%rcount
2353 sdispl => fft_scratch%sdispl
2354 rdispl => fft_scratch%rdispl
2355 pgrid => fft_scratch%pgcube
2356
2357 nx = boin(2, 1, mip) - boin(1, 1, mip) + 1
2358 nz = boin(2, 3, mip) - boin(1, 3, mip) + 1
2359 my = maxval(boout(2, 2, :) - boout(1, 2, :) + 1)
2360
2361 rbuf => fft_scratch%rbuf4
2362 num_threads = 1
2363 my_id = 0
2364!$OMP PARALLEL DEFAULT(NONE), &
2365!$OMP PRIVATE(NUM_THREADS,my_id,lb,ub,ip,ipl,ny,iy,ir), &
2366!$OMP SHARED(rbuf,nz,nx,np,pgrid,boout,cin,my,scount,sdispl)
2367!$ num_threads = MIN(omp_get_max_threads(), SIZE(rbuf, 2))
2368!$ my_id = omp_get_thread_num()
2369 IF (my_id < num_threads) THEN
2370 lb = (SIZE(rbuf, 2)*my_id)/num_threads
2371 ub = (SIZE(rbuf, 2)*(my_id + 1))/num_threads - 1
2372 rbuf(:, lb:ub) = 0.0_dp
2373 END IF
2374!$OMP BARRIER
2375
2376!$OMP DO COLLAPSE(2)
2377 DO izx = 1, nz*nx
2378 DO ip = 0, np - 1
2379 ipl = pgrid(ip, 1)
2380 ny = boout(2, 2, ipl) - boout(1, 2, ipl) + 1
2381 DO iy = boout(1, 2, ipl), boout(2, 2, ipl)
2382 ir = iy - boout(1, 2, ipl) + 1 + (izx - 1)*ny
2383 rbuf(ir, ip) = cin(iy, izx)
2384 END DO
2385 END DO
2386 END DO
2387!$OMP END DO
2388!$OMP DO
2389 DO ip = 0, np - 1
2390 ipl = pgrid(ip, 1)
2391 ny = boout(2, 2, ipl) - boout(1, 2, ipl) + 1
2392 scount(ip) = nx*ny*nz
2393 sdispl(ip) = nx*nz*my*ip
2394 END DO
2395!$OMP END DO
2396!$OMP END PARALLEL
2397 ny = boout(2, 2, mip) - boout(1, 2, mip) + 1
2398!$OMP PARALLEL DO DEFAULT(NONE), &
2399!$OMP PRIVATE(ipl,nx), &
2400!$OMP SHARED(np,pgrid,boin,rcount,rdispl,ny,nz)
2401 DO ip = 0, np - 1
2402 ipl = pgrid(ip, 1)
2403 nx = boin(2, 1, ipl) - boin(1, 1, ipl) + 1
2404 rcount(ip) = nx*ny*nz
2405 rdispl(ip) = ny*nz*(boin(1, 1, ipl) - 1)
2406 END DO
2407!$OMP END PARALLEL DO
2408
2409 CALL sub_group%alltoall(rbuf, scount, sdispl, sout, rcount, rdispl)
2410
2411 CALL timestop(handle)
2412
2413 END SUBROUTINE cube_transpose_4
2414
2415! **************************************************************************************************
2416!> \brief ...
2417!> \param cin ...
2418!> \param group ...
2419!> \param boin ...
2420!> \param boout ...
2421!> \param sout ...
2422!> \param fft_scratch ...
2423! **************************************************************************************************
2424 SUBROUTINE cube_transpose_5(cin, group, boin, boout, sout, fft_scratch)
2425
2426 COMPLEX(KIND=dp), CONTIGUOUS, DIMENSION(:, :), &
2427 INTENT(IN) :: cin
2428
2429 CLASS(mp_comm_type), INTENT(IN) :: group
2430 INTEGER, CONTIGUOUS, DIMENSION(:, :, 0:), INTENT(IN) :: boin, boout
2431 COMPLEX(KIND=dp), DIMENSION(:, :), INTENT(OUT), CONTIGUOUS :: sout
2432 TYPE(fft_scratch_type), INTENT(IN) :: fft_scratch
2433
2434 CHARACTER(len=*), PARAMETER :: routinen = 'cube_transpose_5'
2435
2436 COMPLEX(KIND=dp), DIMENSION(:, :), POINTER, CONTIGUOUS :: rbuf
2437 INTEGER :: handle, ip, ir, is, ixz, iy, lb, mip, &
2438 my, my_id, np, num_threads, nx, ny, &
2439 nz, ub
2440 INTEGER, DIMENSION(:), POINTER, CONTIGUOUS :: rcount, rdispl, scount, sdispl
2441
2442 CALL timeset(routinen, handle)
2443
2444 np = fft_scratch%sizes%numtask
2445 mip = fft_scratch%mip
2446 scount => fft_scratch%scount
2447 rcount => fft_scratch%rcount
2448 sdispl => fft_scratch%sdispl
2449 rdispl => fft_scratch%rdispl
2450
2451 ny = boin(2, 2, mip) - boin(1, 2, mip) + 1
2452 nz = boin(2, 3, mip) - boin(1, 3, mip) + 1
2453!$OMP PARALLEL DO DEFAULT(NONE), &
2454!$OMP PRIVATE(nx), &
2455!$OMP SHARED(np,boout,ny,nz,scount,sdispl)
2456 DO ip = 0, np - 1
2457 nx = boout(2, 1, ip) - boout(1, 1, ip) + 1
2458 scount(ip) = nx*nz*ny
2459 sdispl(ip) = ny*nz*(boout(1, 1, ip) - 1)
2460 END DO
2461!$OMP END PARALLEL DO
2462 nx = boout(2, 1, mip) - boout(1, 1, mip) + 1
2463 my = maxval(boin(2, 2, :) - boin(1, 2, :) + 1)
2464!$OMP PARALLEL DO DEFAULT(NONE), &
2465!$OMP PRIVATE(ny), &
2466!$OMP SHARED(np,boin,nx,nz,rcount,rdispl,my)
2467 DO ip = 0, np - 1
2468 ny = boin(2, 2, ip) - boin(1, 2, ip) + 1
2469 rcount(ip) = nx*nz*ny
2470 rdispl(ip) = nx*my*nz*ip
2471 END DO
2472!$OMP END PARALLEL DO
2473
2474 rbuf => fft_scratch%rbuf5
2475 num_threads = 1
2476 my_id = 0
2477!$OMP PARALLEL DEFAULT(NONE), &
2478!$OMP PRIVATE(NUM_THREADS, my_id, lb, ub), &
2479!$OMP SHARED(rbuf)
2480!$ num_threads = MIN(omp_get_max_threads(), SIZE(rbuf, 2))
2481!$ my_id = omp_get_thread_num()
2482 IF (my_id < num_threads) THEN
2483 lb = (SIZE(rbuf, 2)*my_id)/num_threads
2484 ub = (SIZE(rbuf, 2)*(my_id + 1))/num_threads - 1
2485 rbuf(:, lb:ub) = 0.0_dp
2486 END IF
2487!$OMP END PARALLEL
2488
2489 CALL group%alltoall(cin, scount, sdispl, rbuf, rcount, rdispl)
2490
2491!$OMP PARALLEL DO DEFAULT(NONE) COLLAPSE(2) &
2492!$OMP PRIVATE(ip,ny,iy,is,ir) &
2493!$OMP SHARED(nx,nz,np,boin,sout,rbuf)
2494 DO ixz = 1, nx*nz
2495 DO ip = 0, np - 1
2496 ny = boin(2, 2, ip) - boin(1, 2, ip) + 1
2497 DO iy = 1, ny
2498 is = boin(1, 2, ip) + iy - 1
2499 ir = iy + ny*(ixz - 1)
2500 sout(is, ixz) = rbuf(ir, ip)
2501 END DO
2502 END DO
2503 END DO
2504!$OMP END PARALLEL DO
2505
2506 CALL timestop(handle)
2507
2508 END SUBROUTINE cube_transpose_5
2509
2510! **************************************************************************************************
2511!> \brief ...
2512!> \param cin ...
2513!> \param group ...
2514!> \param boin ...
2515!> \param boout ...
2516!> \param sout ...
2517!> \param fft_scratch ...
2518! **************************************************************************************************
2519 SUBROUTINE cube_transpose_6(cin, group, boin, boout, sout, fft_scratch)
2520
2521 COMPLEX(KIND=dp), CONTIGUOUS, DIMENSION(:, :), &
2522 INTENT(IN) :: cin
2523
2524 CLASS(mp_comm_type), INTENT(IN) :: group
2525 INTEGER, CONTIGUOUS, DIMENSION(:, :, 0:), INTENT(IN) :: boin, boout
2526 COMPLEX(KIND=dp), DIMENSION(:, :), INTENT(OUT), CONTIGUOUS :: sout
2527 TYPE(fft_scratch_type), INTENT(IN) :: fft_scratch
2528
2529 CHARACTER(len=*), PARAMETER :: routinen = 'cube_transpose_6'
2530
2531 COMPLEX(KIND=dp), DIMENSION(:, :), POINTER, CONTIGUOUS :: rbuf
2532 INTEGER :: handle, ip, ir, iy, izx, lb, mip, my, &
2533 my_id, np, num_threads, nx, ny, nz, ub
2534 INTEGER, DIMENSION(:), POINTER, CONTIGUOUS :: rcount, rdispl, scount, sdispl
2535
2536 CALL timeset(routinen, handle)
2537
2538 np = fft_scratch%sizes%numtask
2539 mip = fft_scratch%mip
2540 scount => fft_scratch%scount
2541 rcount => fft_scratch%rcount
2542 sdispl => fft_scratch%sdispl
2543 rdispl => fft_scratch%rdispl
2544
2545 nx = boin(2, 1, mip) - boin(1, 1, mip) + 1
2546 nz = boin(2, 3, mip) - boin(1, 3, mip) + 1
2547 my = maxval(boout(2, 2, :) - boout(1, 2, :) + 1)
2548
2549 rbuf => fft_scratch%rbuf5
2550 num_threads = 1
2551 my_id = 0
2552!$OMP PARALLEL DEFAULT(NONE), &
2553!$OMP PRIVATE(NUM_THREADS,my_id,lb,ub,ip,ny,iy,ir), &
2554!$OMP SHARED(rbuf,nx,nz,np,boout,cin,my,scount,sdispl)
2555!$ num_threads = MIN(omp_get_max_threads(), SIZE(rbuf, 2))
2556!$ my_id = omp_get_thread_num()
2557 IF (my_id < num_threads) THEN
2558 lb = (SIZE(rbuf, 2)*my_id)/num_threads
2559 ub = (SIZE(rbuf, 2)*(my_id + 1))/num_threads - 1
2560 rbuf(:, lb:ub) = 0.0_dp
2561 END IF
2562!$OMP BARRIER
2563
2564!$OMP DO COLLAPSE(2)
2565 DO izx = 1, nz*nx
2566 DO ip = 0, np - 1
2567 ny = boout(2, 2, ip) - boout(1, 2, ip) + 1
2568 DO iy = boout(1, 2, ip), boout(2, 2, ip)
2569 ir = iy - boout(1, 2, ip) + 1 + (izx - 1)*ny
2570 rbuf(ir, ip) = cin(iy, izx)
2571 END DO
2572 END DO
2573 END DO
2574!$OMP END DO
2575!$OMP DO
2576 DO ip = 0, np - 1
2577 ny = boout(2, 2, ip) - boout(1, 2, ip) + 1
2578 scount(ip) = nx*ny*nz
2579 sdispl(ip) = nx*nz*my*ip
2580 END DO
2581!$OMP END DO
2582!$OMP END PARALLEL
2583 ny = boout(2, 2, mip) - boout(1, 2, mip) + 1
2584!$OMP PARALLEL DO DEFAULT(NONE), &
2585!$OMP PRIVATE(nx), &
2586!$OMP SHARED(np,boin,rcount,rdispl,nz,ny)
2587 DO ip = 0, np - 1
2588 nx = boin(2, 1, ip) - boin(1, 1, ip) + 1
2589 rcount(ip) = nx*ny*nz
2590 rdispl(ip) = ny*nz*(boin(1, 1, ip) - 1)
2591 END DO
2592!$OMP END PARALLEL DO
2593
2594 CALL group%alltoall(rbuf, scount, sdispl, sout, rcount, rdispl)
2595
2596 CALL timestop(handle)
2597
2598 END SUBROUTINE cube_transpose_6
2599
2600! **************************************************************************************************
2601!> \brief ...
2602! **************************************************************************************************
2604
2605 CALL release_fft_scratch_pool()
2606
2607 ! Allocate first scratch and mark it as used
2608 ALLOCATE (fft_scratch_first)
2609 ALLOCATE (fft_scratch_first%fft_scratch)
2610 ! this is a very special scratch, it seems, we always keep it 'most - recent' so we will never delete it
2611 fft_scratch_first%fft_scratch%last_tick = huge(fft_scratch_first%fft_scratch%last_tick)
2612
2614
2615 END SUBROUTINE init_fft_scratch_pool
2616
2617! **************************************************************************************************
2618!> \brief ...
2619!> \param fft_scratch ...
2620! **************************************************************************************************
2621 SUBROUTINE deallocate_fft_scratch_type(fft_scratch)
2622 TYPE(fft_scratch_type), INTENT(INOUT) :: fft_scratch
2623
2624#if defined(__OFFLOAD) && !defined(__NO_OFFLOAD_PW)
2625 INTEGER :: ierr
2626 COMPLEX(KIND=dp), POINTER :: dummy_ptr_z
2627#endif
2628
2629 ! deallocate structures
2630 IF (ASSOCIATED(fft_scratch%ziptr)) THEN
2631 CALL fft_dealloc(fft_scratch%ziptr)
2632 END IF
2633 IF (ASSOCIATED(fft_scratch%zoptr)) THEN
2634 CALL fft_dealloc(fft_scratch%zoptr)
2635 END IF
2636 IF (ASSOCIATED(fft_scratch%p1buf)) THEN
2637 CALL fft_dealloc(fft_scratch%p1buf)
2638 END IF
2639 IF (ASSOCIATED(fft_scratch%p2buf)) THEN
2640#if defined(__OFFLOAD) && !defined(__NO_OFFLOAD_PW)
2641 dummy_ptr_z => fft_scratch%p2buf(1, 1)
2642 ierr = offload_free_pinned_mem(c_loc(dummy_ptr_z))
2643#else
2644 CALL fft_dealloc(fft_scratch%p2buf)
2645#endif
2646 END IF
2647 IF (ASSOCIATED(fft_scratch%p3buf)) THEN
2648#if defined(__OFFLOAD) && !defined(__NO_OFFLOAD_PW)
2649 dummy_ptr_z => fft_scratch%p3buf(1, 1)
2650 ierr = offload_free_pinned_mem(c_loc(dummy_ptr_z))
2651#else
2652 CALL fft_dealloc(fft_scratch%p3buf)
2653#endif
2654 END IF
2655 IF (ASSOCIATED(fft_scratch%p4buf)) THEN
2656#if defined(__OFFLOAD) && !defined(__NO_OFFLOAD_PW)
2657 dummy_ptr_z => fft_scratch%p4buf(1, 1)
2658 ierr = offload_free_pinned_mem(c_loc(dummy_ptr_z))
2659#else
2660 CALL fft_dealloc(fft_scratch%p4buf)
2661#endif
2662 END IF
2663 IF (ASSOCIATED(fft_scratch%p5buf)) THEN
2664#if defined(__OFFLOAD) && !defined(__NO_OFFLOAD_PW)
2665 dummy_ptr_z => fft_scratch%p5buf(1, 1)
2666 ierr = offload_free_pinned_mem(c_loc(dummy_ptr_z))
2667#else
2668 CALL fft_dealloc(fft_scratch%p5buf)
2669#endif
2670 END IF
2671 IF (ASSOCIATED(fft_scratch%p6buf)) THEN
2672 CALL fft_dealloc(fft_scratch%p6buf)
2673 END IF
2674 IF (ASSOCIATED(fft_scratch%p7buf)) THEN
2675#if defined(__OFFLOAD) && !defined(__NO_OFFLOAD_PW)
2676 dummy_ptr_z => fft_scratch%p7buf(1, 1)
2677 ierr = offload_free_pinned_mem(c_loc(dummy_ptr_z))
2678#else
2679 CALL fft_dealloc(fft_scratch%p7buf)
2680#endif
2681 END IF
2682 IF (ASSOCIATED(fft_scratch%r1buf)) THEN
2683#if defined(__OFFLOAD) && !defined(__NO_OFFLOAD_PW)
2684 dummy_ptr_z => fft_scratch%r1buf(1, 1)
2685 ierr = offload_free_pinned_mem(c_loc(dummy_ptr_z))
2686#else
2687 CALL fft_dealloc(fft_scratch%r1buf)
2688#endif
2689 END IF
2690 IF (ASSOCIATED(fft_scratch%r2buf)) THEN
2691 CALL fft_dealloc(fft_scratch%r2buf)
2692 END IF
2693 IF (ASSOCIATED(fft_scratch%tbuf)) THEN
2694#if defined(__OFFLOAD) && !defined(__NO_OFFLOAD_PW)
2695 dummy_ptr_z => fft_scratch%tbuf(1, 1, 1)
2696 ierr = offload_free_pinned_mem(c_loc(dummy_ptr_z))
2697#else
2698 CALL fft_dealloc(fft_scratch%tbuf)
2699#endif
2700 END IF
2701 IF (ASSOCIATED(fft_scratch%a1buf)) THEN
2702 CALL fft_dealloc(fft_scratch%a1buf)
2703 END IF
2704 IF (ASSOCIATED(fft_scratch%a2buf)) THEN
2705 CALL fft_dealloc(fft_scratch%a2buf)
2706 END IF
2707 IF (ASSOCIATED(fft_scratch%a3buf)) THEN
2708 CALL fft_dealloc(fft_scratch%a3buf)
2709 END IF
2710 IF (ASSOCIATED(fft_scratch%a4buf)) THEN
2711 CALL fft_dealloc(fft_scratch%a4buf)
2712 END IF
2713 IF (ASSOCIATED(fft_scratch%a5buf)) THEN
2714 CALL fft_dealloc(fft_scratch%a5buf)
2715 END IF
2716 IF (ASSOCIATED(fft_scratch%a6buf)) THEN
2717 CALL fft_dealloc(fft_scratch%a6buf)
2718 END IF
2719 IF (ASSOCIATED(fft_scratch%scount)) THEN
2720 DEALLOCATE (fft_scratch%scount, fft_scratch%rcount, &
2721 fft_scratch%sdispl, fft_scratch%rdispl)
2722 END IF
2723 IF (ASSOCIATED(fft_scratch%rr)) THEN
2724 DEALLOCATE (fft_scratch%rr)
2725 END IF
2726 IF (ASSOCIATED(fft_scratch%xzbuf)) THEN
2727 DEALLOCATE (fft_scratch%xzbuf)
2728 END IF
2729 IF (ASSOCIATED(fft_scratch%yzbuf)) THEN
2730 DEALLOCATE (fft_scratch%yzbuf)
2731 END IF
2732 IF (ASSOCIATED(fft_scratch%xzbuf_sgl)) THEN
2733 DEALLOCATE (fft_scratch%xzbuf_sgl)
2734 END IF
2735 IF (ASSOCIATED(fft_scratch%yzbuf_sgl)) THEN
2736 DEALLOCATE (fft_scratch%yzbuf_sgl)
2737 END IF
2738 IF (ASSOCIATED(fft_scratch%ss)) THEN
2739 DEALLOCATE (fft_scratch%ss)
2740 END IF
2741 IF (ASSOCIATED(fft_scratch%tt)) THEN
2742 DEALLOCATE (fft_scratch%tt)
2743 END IF
2744 IF (ASSOCIATED(fft_scratch%pgrid)) THEN
2745 DEALLOCATE (fft_scratch%pgrid)
2746 END IF
2747 IF (ASSOCIATED(fft_scratch%pgcube)) THEN
2748 DEALLOCATE (fft_scratch%pgcube)
2749 END IF
2750 IF (ASSOCIATED(fft_scratch%xcor)) THEN
2751 DEALLOCATE (fft_scratch%xcor, fft_scratch%zcor)
2752 END IF
2753 IF (ASSOCIATED(fft_scratch%pzcoord)) THEN
2754 DEALLOCATE (fft_scratch%pzcoord)
2755 END IF
2756 IF (ASSOCIATED(fft_scratch%xzcount)) THEN
2757 DEALLOCATE (fft_scratch%xzcount, fft_scratch%yzcount)
2758 DEALLOCATE (fft_scratch%xzdispl, fft_scratch%yzdispl)
2759 fft_scratch%in = 0
2760 fft_scratch%rsratio = 1._dp
2761 END IF
2762 IF (ASSOCIATED(fft_scratch%rbuf1)) THEN
2763 DEALLOCATE (fft_scratch%rbuf1)
2764 END IF
2765 IF (ASSOCIATED(fft_scratch%rbuf2)) THEN
2766 DEALLOCATE (fft_scratch%rbuf2)
2767 END IF
2768 IF (ASSOCIATED(fft_scratch%rbuf3)) THEN
2769 DEALLOCATE (fft_scratch%rbuf3)
2770 END IF
2771 IF (ASSOCIATED(fft_scratch%rbuf4)) THEN
2772 DEALLOCATE (fft_scratch%rbuf4)
2773 END IF
2774 IF (ASSOCIATED(fft_scratch%rbuf5)) THEN
2775 DEALLOCATE (fft_scratch%rbuf5)
2776 END IF
2777 IF (ASSOCIATED(fft_scratch%rbuf6)) THEN
2778 DEALLOCATE (fft_scratch%rbuf6)
2779 END IF
2780
2781 IF (fft_scratch%cart_sub_comm(1) /= mp_comm_null) THEN
2782 CALL fft_scratch%cart_sub_comm(1)%free()
2783 END IF
2784 IF (fft_scratch%cart_sub_comm(2) /= mp_comm_null) THEN
2785 CALL fft_scratch%cart_sub_comm(2)%free()
2786 END IF
2787
2788 CALL fft_destroy_plan(fft_scratch%fft_plan(1))
2789 CALL fft_destroy_plan(fft_scratch%fft_plan(2))
2790 CALL fft_destroy_plan(fft_scratch%fft_plan(3))
2791 CALL fft_destroy_plan(fft_scratch%fft_plan(4))
2792 CALL fft_destroy_plan(fft_scratch%fft_plan(5))
2793 CALL fft_destroy_plan(fft_scratch%fft_plan(6))
2794
2795 END SUBROUTINE deallocate_fft_scratch_type
2796
2797! **************************************************************************************************
2798!> \brief ...
2799! **************************************************************************************************
2800 SUBROUTINE release_fft_scratch_pool()
2801 TYPE(fft_scratch_pool_type), POINTER :: fft_scratch, fft_scratch_current
2802
2803!$ CPASSERT(.NOT. omp_in_parallel() .OR. 0 == omp_get_thread_num())
2804
2805 IF (init_fft_pool == 0) NULLIFY (fft_scratch_first)
2806
2807 fft_scratch => fft_scratch_first
2808 DO
2809 IF (ASSOCIATED(fft_scratch)) THEN
2810 fft_scratch_current => fft_scratch
2811 fft_scratch => fft_scratch_current%fft_scratch_next
2812 NULLIFY (fft_scratch_current%fft_scratch_next)
2813
2814 CALL deallocate_fft_scratch_type(fft_scratch_current%fft_scratch)
2815
2816 DEALLOCATE (fft_scratch_current%fft_scratch)
2817 DEALLOCATE (fft_scratch_current)
2818 ELSE
2819 EXIT
2820 END IF
2821 END DO
2822
2823 init_fft_pool = 0
2824
2825 END SUBROUTINE release_fft_scratch_pool
2826
2827! **************************************************************************************************
2828!> \brief ...
2829! **************************************************************************************************
2830 SUBROUTINE resize_fft_scratch_pool()
2831
2832 INTEGER :: last_tick, nscratch
2833 TYPE(fft_scratch_pool_type), POINTER :: fft_scratch_current, fft_scratch_old
2834
2835 nscratch = 0
2836
2837 last_tick = huge(last_tick)
2838 NULLIFY (fft_scratch_old)
2839
2840 ! start at the global pool, count, and find a deletion candidate
2841 fft_scratch_current => fft_scratch_first
2842 DO
2843 IF (ASSOCIATED(fft_scratch_current)) THEN
2844 nscratch = nscratch + 1
2845 ! is this a candidate for deletion (i.e. least recently used, and not in use)
2846 IF (.NOT. fft_scratch_current%fft_scratch%in_use) THEN
2847 IF (fft_scratch_current%fft_scratch%last_tick < last_tick) THEN
2848 last_tick = fft_scratch_current%fft_scratch%last_tick
2849 fft_scratch_old => fft_scratch_current
2850 END IF
2851 END IF
2852 fft_scratch_current => fft_scratch_current%fft_scratch_next
2853 ELSE
2854 EXIT
2855 END IF
2856 END DO
2857
2858 ! we should delete a scratch
2859 IF (nscratch > fft_pool_scratch_limit) THEN
2860 ! note that we never deallocate the first (special) element of the list
2861 IF (ASSOCIATED(fft_scratch_old)) THEN
2862 fft_scratch_current => fft_scratch_first
2863 DO
2864 IF (ASSOCIATED(fft_scratch_current)) THEN
2865 ! should we delete the next in the list?
2866 IF (ASSOCIATED(fft_scratch_current%fft_scratch_next, fft_scratch_old)) THEN
2867 ! fix the linked list
2868 fft_scratch_current%fft_scratch_next => fft_scratch_old%fft_scratch_next
2869
2870 ! deallocate the element
2871 CALL deallocate_fft_scratch_type(fft_scratch_old%fft_scratch)
2872 DEALLOCATE (fft_scratch_old%fft_scratch)
2873 DEALLOCATE (fft_scratch_old)
2874
2875 ELSE
2876 fft_scratch_current => fft_scratch_current%fft_scratch_next
2877 END IF
2878 ELSE
2879 EXIT
2880 END IF
2881 END DO
2882
2883 ELSE
2884 cpwarn("The number of the scratches exceeded the limit, but none could be deallocated")
2885 END IF
2886 END IF
2887
2888 END SUBROUTINE resize_fft_scratch_pool
2889
2890! **************************************************************************************************
2891!> \brief ...
2892!> \param fft_scratch ...
2893!> \param tf_type ...
2894!> \param n ...
2895!> \param fft_sizes ...
2896! **************************************************************************************************
2897 SUBROUTINE get_fft_scratch(fft_scratch, tf_type, n, fft_sizes)
2898 TYPE(fft_scratch_type), POINTER :: fft_scratch
2899 INTEGER, INTENT(IN) :: tf_type
2900 INTEGER, DIMENSION(:), INTENT(IN) :: n
2901 TYPE(fft_scratch_sizes), INTENT(IN), &
2902 OPTIONAL :: fft_sizes
2903
2904 CHARACTER(len=*), PARAMETER :: routinen = 'get_fft_scratch'
2905
2906 INTEGER :: coord(2), dim(2), handle, i, ix, iz, lg, lmax, m1, m2, &
2907 mcx2, mcy3, mcz1, mcz2, mg, mmax, mx1, mx2, my1, my3, mz1, mz2, mz3, &
2908 nbx, nbz, nm, nmax, nmray, np, nx, ny, nyzray, nz, pos(2)
2909 INTEGER, DIMENSION(3) :: pcoord
2910 LOGICAL :: equal
2911 LOGICAL, DIMENSION(2) :: dims
2912 TYPE(fft_scratch_pool_type), POINTER :: fft_scratch_current, &
2913 fft_scratch_last, &
2914 fft_scratch_new
2915#if defined(__OFFLOAD) && !defined(__NO_OFFLOAD_PW)
2916 INTEGER :: ierr
2917 INTEGER(KIND=C_SIZE_T) :: length
2918 TYPE(c_ptr) :: cptr_r1buf, cptr_tbuf, &
2919 cptr_p2buf, cptr_p3buf, cptr_p4buf, cptr_p5buf, cptr_p7buf
2920#endif
2921 CALL timeset(routinen, handle)
2922
2923 ! this is the place to check that the scratch_pool does not grow without limits
2924 ! before we add a new scratch check the size of the pool and release some of the list if needed
2925 CALL resize_fft_scratch_pool()
2926
2927 ! get the required scratch
2928!$OMP ATOMIC
2930 fft_scratch_current => fft_scratch_first
2931 DO
2932 IF (ASSOCIATED(fft_scratch_current)) THEN
2933 IF (fft_scratch_current%fft_scratch%in_use) THEN
2934 fft_scratch_last => fft_scratch_current
2935 fft_scratch_current => fft_scratch_current%fft_scratch_next
2936 cycle
2937 END IF
2938 IF (tf_type /= fft_scratch_current%fft_scratch%tf_type) THEN
2939 fft_scratch_last => fft_scratch_current
2940 fft_scratch_current => fft_scratch_current%fft_scratch_next
2941 cycle
2942 END IF
2943 IF (.NOT. all(n == fft_scratch_current%fft_scratch%nfft)) THEN
2944 fft_scratch_last => fft_scratch_current
2945 fft_scratch_current => fft_scratch_current%fft_scratch_next
2946 cycle
2947 END IF
2948 IF (PRESENT(fft_sizes)) THEN
2949 IF (fft_sizes%rs_group /= fft_scratch_current%fft_scratch%group) THEN
2950 fft_scratch_last => fft_scratch_current
2951 fft_scratch_current => fft_scratch_current%fft_scratch_next
2952 cycle
2953 END IF
2954 CALL is_equal(fft_sizes, fft_scratch_current%fft_scratch%sizes, equal)
2955 IF (.NOT. equal) THEN
2956 fft_scratch_last => fft_scratch_current
2957 fft_scratch_current => fft_scratch_current%fft_scratch_next
2958 cycle
2959 END IF
2960 END IF
2961 ! Success
2962 fft_scratch => fft_scratch_current%fft_scratch
2963 fft_scratch_current%fft_scratch%in_use = .true.
2964 EXIT
2965 ELSE
2966 ! We cannot find the scratch type in this pool
2967 ! Generate a new scratch set
2968 ALLOCATE (fft_scratch_new)
2969 ALLOCATE (fft_scratch_new%fft_scratch)
2970
2971 IF (tf_type /= 400) THEN
2972 fft_scratch_new%fft_scratch%sizes = fft_sizes
2973 np = fft_sizes%numtask
2974 ALLOCATE (fft_scratch_new%fft_scratch%scount(0:np - 1), fft_scratch_new%fft_scratch%rcount(0:np - 1), &
2975 fft_scratch_new%fft_scratch%sdispl(0:np - 1), fft_scratch_new%fft_scratch%rdispl(0:np - 1), &
2976 fft_scratch_new%fft_scratch%pgcube(0:np - 1, 2))
2977 END IF
2978
2979 SELECT CASE (tf_type)
2980 CASE DEFAULT
2981 cpabort("Invalid scratch type.")
2982 CASE (100) ! fft3d_pb: full cube distribution
2983 cpassert(PRESENT(fft_sizes))
2984 mx1 = fft_sizes%mx1
2985 my1 = fft_sizes%my1
2986 mx2 = fft_sizes%mx2
2987 mz2 = fft_sizes%mz2
2988 my3 = fft_sizes%my3
2989 mz3 = fft_sizes%mz3
2990 CALL fft_alloc(fft_scratch_new%fft_scratch%a1buf, [mx1*my1, n(3)])
2991 CALL fft_alloc(fft_scratch_new%fft_scratch%a2buf, [n(3), mx1*my1])
2992 CALL fft_alloc(fft_scratch_new%fft_scratch%a3buf, [mx2*mz2, n(2)])
2993 CALL fft_alloc(fft_scratch_new%fft_scratch%a4buf, [n(2), mx2*mz2])
2994 CALL fft_alloc(fft_scratch_new%fft_scratch%a5buf, [my3*mz3, n(1)])
2995 CALL fft_alloc(fft_scratch_new%fft_scratch%a6buf, [n(1), my3*mz3])
2996 fft_scratch_new%fft_scratch%group = fft_sizes%rs_group
2997
2998 dim = fft_sizes%rs_group%num_pe_cart
2999 pos = fft_sizes%rs_group%mepos_cart
3000 fft_scratch_new%fft_scratch%mip = fft_sizes%rs_group%mepos
3001 fft_scratch_new%fft_scratch%dim = dim
3002 fft_scratch_new%fft_scratch%pos = pos
3003 mcz1 = fft_sizes%mcz1
3004 mcx2 = fft_sizes%mcx2
3005 mcz2 = fft_sizes%mcz2
3006 mcy3 = fft_sizes%mcy3
3007 ALLOCATE (fft_scratch_new%fft_scratch%rbuf1(mx2*my1*mcz2, 0:dim(2) - 1))
3008 ALLOCATE (fft_scratch_new%fft_scratch%rbuf2(mx1*my1*mcz2, 0:dim(2) - 1))
3009 ALLOCATE (fft_scratch_new%fft_scratch%rbuf3(mx2*mz3*mcy3, 0:dim(1) - 1))
3010 ALLOCATE (fft_scratch_new%fft_scratch%rbuf4(mx2*mz2*mcy3, 0:dim(1) - 1))
3011
3012 dims = [.true., .false.]
3013 CALL fft_scratch_new%fft_scratch%cart_sub_comm(1)%from_sub(fft_sizes%rs_group, dims)
3014 dims = [.false., .true.]
3015 CALL fft_scratch_new%fft_scratch%cart_sub_comm(2)%from_sub(fft_sizes%rs_group, dims)
3016
3017 !initialise pgcube
3018 DO i = 0, dim(1) - 1
3019 coord = [i, pos(2)]
3020 CALL fft_sizes%rs_group%rank_cart(coord, fft_scratch_new%fft_scratch%pgcube(i, 1))
3021 END DO
3022 DO i = 0, dim(2) - 1
3023 coord = [pos(1), i]
3024 CALL fft_sizes%rs_group%rank_cart(coord, fft_scratch_new%fft_scratch%pgcube(i, 2))
3025 END DO
3026
3027 !set up fft plans
3028 CALL fft_create_plan_1d(fft_scratch_new%fft_scratch%fft_plan(1), fwfft, .true., .false., &
3029 mx1*my1, n(3), n(3), mx1*my1, &
3030 fft_scratch_new%fft_scratch%a1buf, fft_scratch_new%fft_scratch%a2buf)
3031 CALL fft_create_plan_1d(fft_scratch_new%fft_scratch%fft_plan(2), fwfft, .true., .false., &
3032 mx2*mz2, n(2), n(2), mx2*mz2, &
3033 fft_scratch_new%fft_scratch%a3buf, fft_scratch_new%fft_scratch%a4buf)
3034 CALL fft_create_plan_1d(fft_scratch_new%fft_scratch%fft_plan(3), fwfft, .true., .false., &
3035 my3*mz3, n(1), n(1), my3*mz3, &
3036 fft_scratch_new%fft_scratch%a5buf, fft_scratch_new%fft_scratch%a6buf)
3037 CALL fft_create_plan_1d(fft_scratch_new%fft_scratch%fft_plan(4), bwfft, .false., .true., &
3038 n(1), my3*mz3, n(1), my3*mz3, &
3039 fft_scratch_new%fft_scratch%a6buf, fft_scratch_new%fft_scratch%a5buf)
3040 CALL fft_create_plan_1d(fft_scratch_new%fft_scratch%fft_plan(5), bwfft, .false., .true., &
3041 n(2), mx2*mz2, n(2), mx2*mz2, &
3042 fft_scratch_new%fft_scratch%a4buf, fft_scratch_new%fft_scratch%a3buf)
3043 CALL fft_create_plan_1d(fft_scratch_new%fft_scratch%fft_plan(6), bwfft, .false., .true., &
3044 n(3), mx1*my1, n(3), mx1*my1, &
3045 fft_scratch_new%fft_scratch%a2buf, fft_scratch_new%fft_scratch%a1buf)
3046
3047 CASE (101) ! fft3d_pb: full cube distribution (dim 1)
3048 cpassert(PRESENT(fft_sizes))
3049 mx1 = fft_sizes%mx1
3050 my1 = fft_sizes%my1
3051 mz1 = fft_sizes%mz1
3052 my3 = fft_sizes%my3
3053 mz3 = fft_sizes%mz3
3054 CALL fft_alloc(fft_scratch_new%fft_scratch%a1buf, [mx1*my1, n(3)])
3055 CALL fft_alloc(fft_scratch_new%fft_scratch%a2buf, [n(3), mx1*my1])
3056 fft_scratch_new%fft_scratch%group = fft_sizes%rs_group
3057 CALL fft_alloc(fft_scratch_new%fft_scratch%a3buf, [mx1*mz1, n(2)])
3058 CALL fft_alloc(fft_scratch_new%fft_scratch%a4buf, [n(2), mx1*mz1])
3059 CALL fft_alloc(fft_scratch_new%fft_scratch%a5buf, [my3*mz3, n(1)])
3060 CALL fft_alloc(fft_scratch_new%fft_scratch%a6buf, [n(1), my3*mz3])
3061
3062 dim = fft_sizes%rs_group%num_pe_cart
3063 pos = fft_sizes%rs_group%mepos_cart
3064 fft_scratch_new%fft_scratch%mip = fft_sizes%rs_group%mepos
3065 fft_scratch_new%fft_scratch%dim = dim
3066 fft_scratch_new%fft_scratch%pos = pos
3067 mcy3 = fft_sizes%mcy3
3068 ALLOCATE (fft_scratch_new%fft_scratch%rbuf5(mx1*mz3*mcy3, 0:dim(1) - 1))
3069 ALLOCATE (fft_scratch_new%fft_scratch%rbuf6(mx1*mz1*mcy3, 0:dim(1) - 1))
3070
3071 !set up fft plans
3072 CALL fft_create_plan_1d(fft_scratch_new%fft_scratch%fft_plan(1), fwfft, .true., .false., &
3073 mx1*my1, n(3), n(3), mx1*my1, &
3074 fft_scratch_new%fft_scratch%a1buf, fft_scratch_new%fft_scratch%a3buf)
3075 CALL fft_create_plan_1d(fft_scratch_new%fft_scratch%fft_plan(2), fwfft, .true., .false., &
3076 mx1*mz1, n(2), n(2), mx1*mz1, &
3077 fft_scratch_new%fft_scratch%a3buf, fft_scratch_new%fft_scratch%a4buf)
3078 CALL fft_create_plan_1d(fft_scratch_new%fft_scratch%fft_plan(3), fwfft, .true., .false., &
3079 my3*mz3, n(1), n(1), my3*mz3, &
3080 fft_scratch_new%fft_scratch%a5buf, fft_scratch_new%fft_scratch%a6buf)
3081 CALL fft_create_plan_1d(fft_scratch_new%fft_scratch%fft_plan(4), bwfft, .false., .true., &
3082 n(1), my3*mz3, n(1), my3*mz3, &
3083 fft_scratch_new%fft_scratch%a6buf, fft_scratch_new%fft_scratch%a5buf)
3084 CALL fft_create_plan_1d(fft_scratch_new%fft_scratch%fft_plan(5), bwfft, .false., .true., &
3085 n(2), mx1*mz1, n(2), mx1*mz1, &
3086 fft_scratch_new%fft_scratch%a4buf, fft_scratch_new%fft_scratch%a3buf)
3087 CALL fft_create_plan_1d(fft_scratch_new%fft_scratch%fft_plan(6), bwfft, .false., .true., &
3088 n(3), mx1*my1, n(3), mx1*my1, &
3089 fft_scratch_new%fft_scratch%a3buf, fft_scratch_new%fft_scratch%a1buf)
3090
3091 CASE (200) ! fft3d_ps: plane distribution
3092 cpassert(PRESENT(fft_sizes))
3093 nx = fft_sizes%nx
3094 ny = fft_sizes%ny
3095 nz = fft_sizes%nz
3096 mx2 = fft_sizes%mx2
3097 lmax = fft_sizes%lmax
3098 mmax = fft_sizes%mmax
3099 lg = fft_sizes%lg
3100 mg = fft_sizes%mg
3101 np = fft_sizes%numtask
3102 nmray = fft_sizes%nmray
3103 nyzray = fft_sizes%nyzray
3104#if defined(__OFFLOAD) && !defined(__NO_OFFLOAD_PW)
3105 length = int(2*dp_size*max(mmax, 1)*max(lmax, 1), kind=c_size_t)
3106 ierr = offload_malloc_pinned_mem(cptr_r1buf, length)
3107 cpassert(ierr == 0)
3108 CALL c_f_pointer(cptr_r1buf, fft_scratch_new%fft_scratch%r1buf, [max(mmax, 1), max(lmax, 1)])
3109 length = int(2*dp_size*max(ny, 1)*max(nz, 1)*max(nx, 1), kind=c_size_t)
3110 ierr = offload_malloc_pinned_mem(cptr_tbuf, length)
3111 cpassert(ierr == 0)
3112 CALL c_f_pointer(cptr_tbuf, fft_scratch_new%fft_scratch%tbuf, [max(ny, 1), max(nz, 1), max(nx, 1)])
3113#else
3114 CALL fft_alloc(fft_scratch_new%fft_scratch%r1buf, [mmax, lmax])
3115 CALL fft_alloc(fft_scratch_new%fft_scratch%tbuf, [ny, nz, nx])
3116#endif
3117 fft_scratch_new%fft_scratch%group = fft_sizes%rs_group
3118 CALL fft_alloc(fft_scratch_new%fft_scratch%r2buf, [lg, mg])
3119 nm = nmray*mx2
3120 IF (alltoall_sgl) THEN
3121 ALLOCATE (fft_scratch_new%fft_scratch%ss(mmax, lmax))
3122 ALLOCATE (fft_scratch_new%fft_scratch%tt(nm, 0:np - 1))
3123 ELSE
3124 ALLOCATE (fft_scratch_new%fft_scratch%rr(nm, 0:np - 1))
3125 END IF
3126
3127 !set up fft plans
3128 CALL fft_create_plan_1d(fft_scratch_new%fft_scratch%fft_plan(1), fwfft, .true., .false., &
3129 nx*ny, nz, nz, nx*ny, &
3130 fft_scratch_new%fft_scratch%tbuf, fft_scratch_new%fft_scratch%r1buf)
3131 CALL fft_create_plan_1d(fft_scratch_new%fft_scratch%fft_plan(2), fwfft, .true., .false., &
3132 nx*nz, ny, ny, nx*nz, &
3133 fft_scratch_new%fft_scratch%r1buf, fft_scratch_new%fft_scratch%tbuf)
3134 CALL fft_create_plan_1d(fft_scratch_new%fft_scratch%fft_plan(3), fwfft, .true., .false., &
3135 nyzray, n(1), n(1), nyzray, &
3136 fft_scratch_new%fft_scratch%r1buf, fft_scratch_new%fft_scratch%r2buf)
3137 CALL fft_create_plan_1d(fft_scratch_new%fft_scratch%fft_plan(4), bwfft, .false., .true., &
3138 n(1), nyzray, n(1), nyzray, &
3139 fft_scratch_new%fft_scratch%r2buf, fft_scratch_new%fft_scratch%r1buf)
3140 CALL fft_create_plan_1d(fft_scratch_new%fft_scratch%fft_plan(5), bwfft, .false., .true., &
3141 ny, nx*nz, ny, nx*nz, &
3142 fft_scratch_new%fft_scratch%tbuf, fft_scratch_new%fft_scratch%r1buf)
3143 CALL fft_create_plan_1d(fft_scratch_new%fft_scratch%fft_plan(6), bwfft, .false., .true., &
3144 nz, nx*ny, nz, nx*ny, &
3145 fft_scratch_new%fft_scratch%r1buf, fft_scratch_new%fft_scratch%tbuf)
3146
3147 CASE (300) ! fft3d_ps: block distribution
3148 cpassert(PRESENT(fft_sizes))
3149 mx1 = fft_sizes%mx1
3150 mx2 = fft_sizes%mx2
3151 my1 = fft_sizes%my1
3152 mz2 = fft_sizes%mz2
3153 mcx2 = fft_sizes%mcx2
3154 lg = fft_sizes%lg
3155 mg = fft_sizes%mg
3156 nmax = fft_sizes%nmax
3157 nmray = fft_sizes%nmray
3158 nyzray = fft_sizes%nyzray
3159 m1 = fft_sizes%r_dim(1)
3160 m2 = fft_sizes%r_dim(2)
3161 nbx = fft_sizes%nbx
3162 nbz = fft_sizes%nbz
3163 CALL fft_alloc(fft_scratch_new%fft_scratch%p1buf, [mx1*my1, n(3)])
3164 CALL fft_alloc(fft_scratch_new%fft_scratch%p6buf, [lg, mg])
3165#if defined(__OFFLOAD) && !defined(__NO_OFFLOAD_PW)
3166 length = int(2*dp_size*max(n(3), 1)*max(mx1*my1, 1), kind=c_size_t)
3167 ierr = offload_malloc_pinned_mem(cptr_p2buf, length)
3168 cpassert(ierr == 0)
3169 CALL c_f_pointer(cptr_p2buf, fft_scratch_new%fft_scratch%p2buf, [max(n(3), 1), max(mx1*my1, 1)])
3170 length = int(2*dp_size*max(mx2*mz2, 1)*max(n(2), 1), kind=c_size_t)
3171 ierr = offload_malloc_pinned_mem(cptr_p3buf, length)
3172 cpassert(ierr == 0)
3173 CALL c_f_pointer(cptr_p3buf, fft_scratch_new%fft_scratch%p3buf, [max(mx2*mz2, 1), max(n(2), 1)])
3174 length = int(2*dp_size*max(n(2), 1)*max(mx2*mz2, 1), kind=c_size_t)
3175 ierr = offload_malloc_pinned_mem(cptr_p4buf, length)
3176 cpassert(ierr == 0)
3177 CALL c_f_pointer(cptr_p4buf, fft_scratch_new%fft_scratch%p4buf, [max(n(2), 1), max(mx2*mz2, 1)])
3178 length = int(2*dp_size*max(nyzray, 1)*max(n(1), 1), kind=c_size_t)
3179 ierr = offload_malloc_pinned_mem(cptr_p5buf, length)
3180 cpassert(ierr == 0)
3181 CALL c_f_pointer(cptr_p5buf, fft_scratch_new%fft_scratch%p5buf, [max(nyzray, 1), max(n(1), 1)])
3182 length = int(2*dp_size*max(mg, 1)*max(lg, 1), kind=c_size_t)
3183 ierr = offload_malloc_pinned_mem(cptr_p7buf, length)
3184 cpassert(ierr == 0)
3185 CALL c_f_pointer(cptr_p7buf, fft_scratch_new%fft_scratch%p7buf, [max(mg, 1), max(lg, 1)])
3186#else
3187 CALL fft_alloc(fft_scratch_new%fft_scratch%p2buf, [n(3), mx1*my1])
3188 CALL fft_alloc(fft_scratch_new%fft_scratch%p3buf, [mx2*mz2, n(2)])
3189 CALL fft_alloc(fft_scratch_new%fft_scratch%p4buf, [n(2), mx2*mz2])
3190 CALL fft_alloc(fft_scratch_new%fft_scratch%p5buf, [nyzray, n(1)])
3191 CALL fft_alloc(fft_scratch_new%fft_scratch%p7buf, [mg, lg])
3192#endif
3193 IF (alltoall_sgl) THEN
3194 ALLOCATE (fft_scratch_new%fft_scratch%yzbuf_sgl(mg*lg))
3195 ALLOCATE (fft_scratch_new%fft_scratch%xzbuf_sgl(n(2)*mx2*mz2))
3196 ELSE
3197 ALLOCATE (fft_scratch_new%fft_scratch%yzbuf(mg*lg))
3198 ALLOCATE (fft_scratch_new%fft_scratch%xzbuf(n(2)*mx2*mz2))
3199 END IF
3200 ALLOCATE (fft_scratch_new%fft_scratch%pgrid(0:m1 - 1, 0:m2 - 1))
3201 ALLOCATE (fft_scratch_new%fft_scratch%xcor(nbx))
3202 ALLOCATE (fft_scratch_new%fft_scratch%zcor(nbz))
3203 ALLOCATE (fft_scratch_new%fft_scratch%pzcoord(0:np - 1))
3204 ALLOCATE (fft_scratch_new%fft_scratch%xzcount(0:np - 1), &
3205 fft_scratch_new%fft_scratch%yzcount(0:np - 1))
3206 ALLOCATE (fft_scratch_new%fft_scratch%xzdispl(0:np - 1), &
3207 fft_scratch_new%fft_scratch%yzdispl(0:np - 1))
3208 fft_scratch_new%fft_scratch%group = fft_sizes%rs_group
3209
3210 dim = fft_sizes%rs_group%num_pe_cart
3211 pos = fft_sizes%rs_group%mepos_cart
3212 fft_scratch_new%fft_scratch%mip = fft_sizes%rs_group%mepos
3213 fft_scratch_new%fft_scratch%dim = dim
3214 fft_scratch_new%fft_scratch%pos = pos
3215 mcz1 = fft_sizes%mcz1
3216 mcz2 = fft_sizes%mcz2
3217 ALLOCATE (fft_scratch_new%fft_scratch%rbuf1(mx2*my1*mcz2, 0:dim(2) - 1))
3218 ALLOCATE (fft_scratch_new%fft_scratch%rbuf2(mx1*my1*mcz2, 0:dim(2) - 1))
3219
3220 dims = [.false., .true.]
3221 CALL fft_scratch_new%fft_scratch%cart_sub_comm(2)%from_sub(fft_sizes%rs_group, dims)
3222
3223 !initialise pgcube
3224 DO i = 0, dim(2) - 1
3225 coord = [pos(1), i]
3226 CALL fft_sizes%rs_group%rank_cart(coord, fft_scratch_new%fft_scratch%pgcube(i, 2))
3227 END DO
3228
3229 !initialise pgrid
3230 DO ix = 0, m1 - 1
3231 DO iz = 0, m2 - 1
3232 coord = [ix, iz]
3233 CALL fft_sizes%rs_group%rank_cart(coord, fft_scratch_new%fft_scratch%pgrid(ix, iz))
3234 END DO
3235 END DO
3236
3237 !initialise pzcoord
3238 DO i = 0, np - 1
3239 CALL fft_sizes%rs_group%coords(i, pcoord)
3240 fft_scratch_new%fft_scratch%pzcoord(i) = pcoord(2)
3241 END DO
3242
3243 !set up fft plans
3244 ! The scratch buffers must cover the geometry declared in each call
3245 ! below, because the FFTW planner may write the buffers during
3246 ! planning. The plans pair the buffers such that each of p1-p6 is
3247 ! checked once; p7 shares the size of p6.
3248 cpassert(SIZE(fft_scratch_new%fft_scratch%p1buf) >= mx1*my1*n(3))
3249 CALL fft_create_plan_1d(fft_scratch_new%fft_scratch%fft_plan(1), fwfft, .true., .false., &
3250 mx1*my1, n(3), n(3), mx1*my1, &
3251 fft_scratch_new%fft_scratch%p1buf, fft_scratch_new%fft_scratch%p2buf)
3252 cpassert(SIZE(fft_scratch_new%fft_scratch%p3buf) >= mx2*mz2*n(2))
3253 CALL fft_create_plan_1d(fft_scratch_new%fft_scratch%fft_plan(2), fwfft, .true., .false., &
3254 mx2*mz2, n(2), n(2), mx2*mz2, &
3255 fft_scratch_new%fft_scratch%p3buf, fft_scratch_new%fft_scratch%p4buf)
3256 cpassert(SIZE(fft_scratch_new%fft_scratch%p5buf) >= nyzray*n(1))
3257 CALL fft_create_plan_1d(fft_scratch_new%fft_scratch%fft_plan(3), fwfft, .true., .false., &
3258 nyzray, n(1), n(1), nyzray, &
3259 fft_scratch_new%fft_scratch%p5buf, fft_scratch_new%fft_scratch%p6buf)
3260 cpassert(SIZE(fft_scratch_new%fft_scratch%p6buf) >= n(1)*nyzray)
3261 CALL fft_create_plan_1d(fft_scratch_new%fft_scratch%fft_plan(4), bwfft, .false., .true., &
3262 n(1), nyzray, n(1), nyzray, &
3263 fft_scratch_new%fft_scratch%p6buf, fft_scratch_new%fft_scratch%p7buf)
3264 cpassert(SIZE(fft_scratch_new%fft_scratch%p4buf) >= n(2)*mx2*mz2)
3265 CALL fft_create_plan_1d(fft_scratch_new%fft_scratch%fft_plan(5), bwfft, .false., .true., &
3266 n(2), mx2*mz2, n(2), mx2*mz2, &
3267 fft_scratch_new%fft_scratch%p4buf, fft_scratch_new%fft_scratch%p3buf)
3268 cpassert(SIZE(fft_scratch_new%fft_scratch%p2buf) >= n(3)*mx1*my1)
3269 CALL fft_create_plan_1d(fft_scratch_new%fft_scratch%fft_plan(6), bwfft, .false., .true., &
3270 n(3), mx1*my1, n(3), mx1*my1, &
3271 fft_scratch_new%fft_scratch%p2buf, fft_scratch_new%fft_scratch%p1buf)
3272
3273 CASE (400) ! serial FFT
3274 np = 0
3275 CALL fft_alloc(fft_scratch_new%fft_scratch%ziptr, n)
3276 CALL fft_alloc(fft_scratch_new%fft_scratch%zoptr, n)
3277
3278 !in place plans
3279 CALL fft_create_plan_3d(fft_scratch_new%fft_scratch%fft_plan(1), .true., fwfft, n, &
3280 fft_scratch_new%fft_scratch%ziptr, fft_scratch_new%fft_scratch%ziptr)
3281 CALL fft_create_plan_3d(fft_scratch_new%fft_scratch%fft_plan(2), .true., bwfft, n, &
3282 fft_scratch_new%fft_scratch%ziptr, fft_scratch_new%fft_scratch%ziptr)
3283 ! out of place plans
3284 CALL fft_create_plan_3d(fft_scratch_new%fft_scratch%fft_plan(3), .false., fwfft, n, &
3285 fft_scratch_new%fft_scratch%ziptr, fft_scratch_new%fft_scratch%zoptr)
3286 CALL fft_create_plan_3d(fft_scratch_new%fft_scratch%fft_plan(4), .false., bwfft, n, &
3287 fft_scratch_new%fft_scratch%ziptr, fft_scratch_new%fft_scratch%zoptr)
3288
3289 END SELECT
3290
3291 NULLIFY (fft_scratch_new%fft_scratch_next)
3292 fft_scratch_new%fft_scratch%fft_scratch_id = &
3293 fft_scratch_last%fft_scratch%fft_scratch_id + 1
3294 fft_scratch_new%fft_scratch%in_use = .true.
3295 fft_scratch_new%fft_scratch%nfft = n
3296 fft_scratch_last%fft_scratch_next => fft_scratch_new
3297 fft_scratch_new%fft_scratch%tf_type = tf_type
3298 fft_scratch => fft_scratch_new%fft_scratch
3299 EXIT
3300
3301 END IF
3302 END DO
3303
3304!$OMP ATOMIC READ
3305 fft_scratch%last_tick = tick_fft_pool
3306
3307 CALL timestop(handle)
3308
3309 END SUBROUTINE get_fft_scratch
3310
3311! **************************************************************************************************
3312!> \brief ...
3313!> \param fft_scratch ...
3314! **************************************************************************************************
3315 SUBROUTINE release_fft_scratch(fft_scratch)
3316
3317 TYPE(fft_scratch_type), POINTER :: fft_scratch
3318
3319 INTEGER :: scratch_id
3320 TYPE(fft_scratch_pool_type), POINTER :: fft_scratch_current
3321
3322 scratch_id = fft_scratch%fft_scratch_id
3323
3324 fft_scratch_current => fft_scratch_first
3325 DO
3326 IF (ASSOCIATED(fft_scratch_current)) THEN
3327 IF (scratch_id == fft_scratch_current%fft_scratch%fft_scratch_id) THEN
3328 fft_scratch%in_use = .false.
3329 NULLIFY (fft_scratch)
3330 EXIT
3331 END IF
3332 fft_scratch_current => fft_scratch_current%fft_scratch_next
3333 ELSE
3334 ! We cannot find the scratch type in this pool
3335 cpabort("Invalid scratch type.")
3336 EXIT
3337 END IF
3338 END DO
3339
3340 END SUBROUTINE release_fft_scratch
3341
3342! **************************************************************************************************
3343!> \brief ...
3344!> \param rs ...
3345!> \param scount ...
3346!> \param sdispl ...
3347!> \param rq ...
3348!> \param rcount ...
3349!> \param rdispl ...
3350!> \param group ...
3351! **************************************************************************************************
3352 SUBROUTINE sparse_alltoall(rs, scount, sdispl, rq, rcount, rdispl, group)
3353 COMPLEX(KIND=dp), DIMENSION(:), POINTER :: rs
3354 INTEGER, DIMENSION(:), POINTER :: scount, sdispl
3355 COMPLEX(KIND=dp), DIMENSION(:), POINTER :: rq
3356 INTEGER, DIMENSION(:), POINTER :: rcount, rdispl
3357
3358 CLASS(mp_comm_type), INTENT(IN) :: group
3359
3360 COMPLEX(KIND=dp), DIMENSION(:), POINTER :: msgin, msgout
3361 INTEGER :: ip, n, nr, ns, pos
3362 TYPE(mp_request_type), ALLOCATABLE, DIMENSION(:) :: rreq, sreq
3363
3364 CALL group%sync()
3365 n = group%num_pe
3366 pos = group%mepos
3367 ALLOCATE (sreq(0:n - 1))
3368 ALLOCATE (rreq(0:n - 1))
3369 nr = 0
3370 DO ip = 0, n - 1
3371 IF (rcount(ip) == 0) cycle
3372 IF (ip == pos) cycle
3373 msgout => rq(rdispl(ip) + 1:rdispl(ip) + rcount(ip))
3374 CALL group%irecv(msgout, ip, rreq(nr))
3375 nr = nr + 1
3376 END DO
3377 ns = 0
3378 DO ip = 0, n - 1
3379 IF (scount(ip) == 0) cycle
3380 IF (ip == pos) cycle
3381 msgin => rs(sdispl(ip) + 1:sdispl(ip) + scount(ip))
3382 CALL group%isend(msgin, ip, sreq(ns))
3383 ns = ns + 1
3384 END DO
3385 IF (rcount(pos) /= 0) THEN
3386 IF (rcount(pos) /= scount(pos)) cpabort("Invalid count.")
3387 rq(rdispl(pos) + 1:rdispl(pos) + rcount(pos)) = rs(sdispl(pos) + 1:sdispl(pos) + scount(pos))
3388 END IF
3389 CALL mp_waitall(sreq(0:ns - 1))
3390 CALL mp_waitall(rreq(0:nr - 1))
3391 DEALLOCATE (sreq)
3392 DEALLOCATE (rreq)
3393 CALL group%sync()
3394
3395 END SUBROUTINE sparse_alltoall
3396
3397! **************************************************************************************************
3398!> \brief test data structures for equality. It is assumed that if they are
3399!> different for one mpi task they are different for all (??)
3400!> \param fft_size_1 ...
3401!> \param fft_size_2 ...
3402!> \param equal ...
3403! **************************************************************************************************
3404 SUBROUTINE is_equal(fft_size_1, fft_size_2, equal)
3405 TYPE(fft_scratch_sizes) :: fft_size_1, fft_size_2
3406 LOGICAL :: equal
3407
3408 equal = .true.
3409
3410 equal = equal .AND. fft_size_1%nx == fft_size_2%nx
3411 equal = equal .AND. fft_size_1%ny == fft_size_2%ny
3412 equal = equal .AND. fft_size_1%nz == fft_size_2%nz
3413
3414 equal = equal .AND. fft_size_1%lmax == fft_size_2%lmax
3415 equal = equal .AND. fft_size_1%mmax == fft_size_2%mmax
3416 equal = equal .AND. fft_size_1%nmax == fft_size_2%nmax
3417
3418 equal = equal .AND. fft_size_1%mx1 == fft_size_2%mx1
3419 equal = equal .AND. fft_size_1%mx2 == fft_size_2%mx2
3420 equal = equal .AND. fft_size_1%mx3 == fft_size_2%mx3
3421
3422 equal = equal .AND. fft_size_1%my1 == fft_size_2%my1
3423 equal = equal .AND. fft_size_1%my2 == fft_size_2%my2
3424 equal = equal .AND. fft_size_1%my3 == fft_size_2%my3
3425
3426 equal = equal .AND. fft_size_1%mcz1 == fft_size_2%mcz1
3427 equal = equal .AND. fft_size_1%mcx2 == fft_size_2%mcx2
3428 equal = equal .AND. fft_size_1%mcz2 == fft_size_2%mcz2
3429 equal = equal .AND. fft_size_1%mcy3 == fft_size_2%mcy3
3430
3431 equal = equal .AND. fft_size_1%lg == fft_size_2%lg
3432 equal = equal .AND. fft_size_1%mg == fft_size_2%mg
3433
3434 equal = equal .AND. fft_size_1%nbx == fft_size_2%nbx
3435 equal = equal .AND. fft_size_1%nbz == fft_size_2%nbz
3436
3437 equal = equal .AND. fft_size_1%nmray == fft_size_2%nmray
3438 equal = equal .AND. fft_size_1%nyzray == fft_size_2%nyzray
3439
3440 equal = equal .AND. fft_size_1%rs_group == fft_size_2%rs_group
3441
3442 equal = equal .AND. all(fft_size_1%g_pos == fft_size_2%g_pos)
3443 equal = equal .AND. all(fft_size_1%r_pos == fft_size_2%r_pos)
3444 equal = equal .AND. all(fft_size_1%r_dim == fft_size_2%r_dim)
3445
3446 equal = equal .AND. fft_size_1%numtask == fft_size_2%numtask
3447
3448 END SUBROUTINE is_equal
3449
3450END MODULE fft_tools
various routines to log and control the output. The idea is that decisions about where to log should ...
integer function, public cp_logger_get_default_io_unit(logger)
returns the unit nr for the ionode (-1 on all other processors) skips as well checks if the procs cal...
subroutine, public fft_3d(plan, scale, zin, zout, stat)
...
Definition fft_lib.F:172
subroutine, public fft_destroy_plan(plan)
...
Definition fft_lib.F:256
subroutine, public fft_1d(plan, zin, zout, scale, stat)
...
Definition fft_lib.F:276
subroutine, public fft_do_cleanup(wisdom_file, ionode, export_wisdom)
...
Definition fft_lib.F:99
subroutine, public fft_create_plan_1d(plan, fsign, trans_in, trans_out, ldx_in, ldx_out, n, m, zin, zout)
...
Definition fft_lib.F:213
subroutine, public fft_do_init(fftlib, plan_style, wisdom_file)
...
Definition fft_lib.F:76
subroutine, public fft_get_lengths(data, max_length)
...
Definition fft_lib.F:120
subroutine, public fft_create_plan_3d(plan, fft_in_place, fsign, n, zin, zout)
...
Definition fft_lib.F:139
Type to store data about a (1D or 3D) FFT, including FFTW plan.
Definition fft_plan.F:18
subroutine, public fft_radix_operations(radix_in, radix_out, operation)
Determine the allowed lengths of FFT's '''.
Definition fft_tools.F:233
integer, save init_fft_pool
Definition fft_tools.F:129
subroutine, public yz_to_x(tb, group, my_pos, p2p, yzp, nray, bo, sb, fft_scratch)
...
Definition fft_tools.F:1517
integer, parameter, public bwfft
Definition fft_tools.F:146
subroutine, public yz_to_xz(sb, group, dims, my_pos, p2p, yzp, nray, bo, tb, fft_scratch)
...
Definition fft_tools.F:1634
subroutine, public init_fft_scratch_pool()
...
Definition fft_tools.F:2604
subroutine, public finalize_fft(para_env, wisdom_file, export_wisdom)
does whatever is needed to finalize the current fft setup
Definition fft_tools.F:209
subroutine, public cube_transpose_1(cin, boin, boout, sout, fft_scratch)
...
Definition fft_tools.F:2038
subroutine, public fft_1d_many(sign, n, m, trans_in, trans_out, ldx_in, ldx_out, zin, zout, scale, stat)
Performs m 1-D forward FFT-s of size n.
Definition fft_tools.F:352
integer, parameter, public fft_radix_closest
Definition fft_tools.F:147
subroutine, public release_fft_scratch(fft_scratch)
...
Definition fft_tools.F:3316
subroutine, public x_to_yz(sb, group, my_pos, p2p, yzp, nray, bo, tb, fft_scratch)
...
Definition fft_tools.F:1401
integer, save tick_fft_pool
Definition fft_tools.F:131
subroutine, public init_fft(fftlib, alltoall, pool_limit, wisdom_file, plan_style)
...
Definition fft_tools.F:181
subroutine, public get_fft_scratch(fft_scratch, tf_type, n, fft_sizes)
...
Definition fft_tools.F:2898
integer, parameter, public fwfft
Definition fft_tools.F:146
integer, parameter, public fft_radix_next_odd
Definition fft_tools.F:149
subroutine, public xz_to_yz(sb, group, dims, my_pos, p2p, yzp, nray, bo, tb, fft_scratch)
...
Definition fft_tools.F:1844
subroutine, public cube_transpose_2(cin, boin, boout, sout, fft_scratch)
...
Definition fft_tools.F:2128
type(fft_scratch_pool_type), pointer, save fft_scratch_first
Definition fft_tools.F:134
integer, save fft_pool_scratch_limit
Definition fft_tools.F:133
integer, parameter, public fft_radix_next
Definition fft_tools.F:147
Defines the basic variable types.
Definition kinds.F:23
integer, parameter, public dp_size
Definition kinds.F:36
integer, parameter, public dp
Definition kinds.F:34
integer, parameter, public sp
Definition kinds.F:33
Definition of mathematical constants and functions.
complex(kind=dp), parameter, public z_zero
Interface to the message passing library MPI.
type(mp_comm_type), parameter, public mp_comm_null
Fortran API for the offload package, which is written in C.
Definition offload_api.F:12
integer function, public offload_free_pinned_mem(buffer)
free pinned memory
Definition offload_api.F:75
integer function, public offload_malloc_pinned_mem(buffer, length)
allocate pinned memory.
Definition offload_api.F:52