Changeset 3634 for palm/trunk/SOURCE/fft_xy_mod.f90

Timestamp:

Dec 18, 2018 12:31:28 PM (5 years ago)

Author:

knoop

Message:

OpenACC port for SPEC

File:

• palm/trunk/SOURCE/fft_xy_mod.f90 (modified) (9 diffs)

palm/trunk/SOURCE/fft_xy_mod.f90

@@ -25 +25 @@
 ! -----------------
 ! $Id$
+! OpenACC port for SPEC
+! 
+! 3241 2018-09-12 15:02:00Z raasch
 ! preprocessor switches for variables that are required on NEC only
 ! 
@@ -156 +159 @@
         ONLY:  fft_method, message_string
         
+    USE cuda_fft_interfaces
+        
     USE indices,                                                               &
         ONLY:  nx, ny, nz
         
-#if defined( __fftw )
+#if defined( __cuda_fft )
+    USE ISO_C_BINDING
+#elif defined( __fftw )
     USE, INTRINSIC ::  ISO_C_BINDING
 #endif
@@ -210 +217 @@
     REAL(wp), DIMENSION(:), ALLOCATABLE, SAVE ::  trig_yf  !<
     
+#elif defined( __cuda_fft )
+    INTEGER(C_INT), SAVE ::  plan_xf  !<
+    INTEGER(C_INT), SAVE ::  plan_xi  !<
+    INTEGER(C_INT), SAVE ::  plan_yf  !<
+    INTEGER(C_INT), SAVE ::  plan_yi  !<
+
 #endif
 
@@ -298 +311 @@
        ENDIF
 
+#if defined( _OPENACC ) && defined( __cuda_fft )
+       fft_method = 'system-specific'
+#endif
+
        IF ( fft_method == 'system-specific' )  THEN
 
@@ -346 +363 @@
           CALL ZDFFTM( 0, ny+1, nz1, sqr_dny, work_y, ny+4, work_y, ny+4,      &
                        trig_yb, worky, 0 )
+#elif defined( __cuda_fft )
+          CALL CUFFTPLAN1D( plan_xf, nx+1, CUFFT_D2Z, (nyn_x-nys_x+1) * (nzt_x-nzb_x+1) )
+          CALL CUFFTPLAN1D( plan_xi, nx+1, CUFFT_Z2D, (nyn_x-nys_x+1) * (nzt_x-nzb_x+1) )
+          CALL CUFFTPLAN1D( plan_yf, ny+1, CUFFT_D2Z, (nxr_y-nxl_y+1) * (nzt_y-nzb_y+1) )
+          CALL CUFFTPLAN1D( plan_yi, ny+1, CUFFT_Z2D, (nxr_y-nxl_y+1) * (nzt_y-nzb_y+1) )
 #else
           message_string = 'no system-specific fft-call available'
@@ -425 +447 @@
 #elif defined( __nec )
        REAL(wp), DIMENSION(6*(nx+1)) ::  work2  !<
+#elif defined( __cuda_fft )
+       COMPLEX(dp), DIMENSION(0:(nx+1)/2,nys_x:nyn_x,nzb_x:nzt_x) ::           &
+          ar_tmp  !<
+       !$ACC DECLARE CREATE(ar_tmp)
 #endif
 
@@ -726 +752 @@
              ENDDO
              !$OMP END PARALLEL
+
+          ENDIF
+
+#elif defined( __cuda_fft )
+
+          IF ( forward_fft )  THEN
+
+             !$ACC HOST_DATA USE_DEVICE(ar, ar_tmp)
+             CALL CUFFTEXECD2Z( plan_xf, ar, ar_tmp )
+             !$ACC END HOST_DATA
+
+             !$ACC PARALLEL LOOP COLLAPSE(2) PRIVATE(i,j,k) &
+             !$ACC PRESENT(ar, ar_tmp)
+             DO  k = nzb_x, nzt_x
+                DO  j = nys_x, nyn_x
+
+                   DO  i = 0, (nx+1)/2
+                      ar(i,j,k)      = REAL( ar_tmp(i,j,k), KIND=wp )  * dnx
+                   ENDDO
+
+                   DO  i = 1, (nx+1)/2 - 1
+                      ar(nx+1-i,j,k) = AIMAG( ar_tmp(i,j,k) ) * dnx
+                   ENDDO
+
+                ENDDO
+             ENDDO
+
+          ELSE
+
+             !$ACC PARALLEL LOOP COLLAPSE(2) PRIVATE(i,j,k) &
+             !$ACC PRESENT(ar, ar_tmp)
+             DO  k = nzb_x, nzt_x
+                DO  j = nys_x, nyn_x
+
+                   ar_tmp(0,j,k) = CMPLX( ar(0,j,k), 0.0_wp, KIND=wp )
+
+                   DO  i = 1, (nx+1)/2 - 1
+                      ar_tmp(i,j,k) = CMPLX( ar(i,j,k), ar(nx+1-i,j,k),        &
+                                             KIND=wp )
+                   ENDDO
+                   ar_tmp((nx+1)/2,j,k) = CMPLX( ar((nx+1)/2,j,k), 0.0_wp,     &
+                                                 KIND=wp )
+
+                ENDDO
+             ENDDO
+
+             !$ACC HOST_DATA USE_DEVICE(ar, ar_tmp)
+             CALL CUFFTEXECZ2D( plan_xi, ar_tmp, ar )
+             !$ACC END HOST_DATA
 
           ENDIF
@@ -1001 +1076 @@
 #elif defined( __nec )
        REAL(wp), DIMENSION(6*(ny+1)) ::  work2  !<
+#elif defined( __cuda_fft )
+       COMPLEX(dp), DIMENSION(0:(ny+1)/2,nxl_y:nxr_y,nzb_y:nzt_y) ::           &
+          ar_tmp  !<
+       !$ACC DECLARE CREATE(ar_tmp)
 #endif
 
@@ -1278 +1357 @@
 
           ENDIF
+#elif defined( __cuda_fft )
+
+          IF ( forward_fft )  THEN
+
+             !$ACC HOST_DATA USE_DEVICE(ar, ar_tmp)
+             CALL CUFFTEXECD2Z( plan_yf, ar, ar_tmp )
+             !$ACC END HOST_DATA
+
+             !$ACC PARALLEL LOOP COLLAPSE(2) PRIVATE(i,j,k) &
+             !$ACC PRESENT(ar, ar_tmp)
+             DO  k = nzb_y, nzt_y
+                DO  i = nxl_y, nxr_y
+
+                   DO  j = 0, (ny+1)/2
+                      ar(j,i,k)      = REAL( ar_tmp(j,i,k), KIND=wp )  * dny
+                   ENDDO
+
+                   DO  j = 1, (ny+1)/2 - 1
+                      ar(ny+1-j,i,k) = AIMAG( ar_tmp(j,i,k) ) * dny
+                   ENDDO
+
+                ENDDO
+             ENDDO
+
+          ELSE
+
+             !$ACC PARALLEL LOOP COLLAPSE(2) PRIVATE(i,j,k) &
+             !$ACC PRESENT(ar, ar_tmp)
+             DO  k = nzb_y, nzt_y
+                DO  i = nxl_y, nxr_y
+
+                   ar_tmp(0,i,k) = CMPLX( ar(0,i,k), 0.0_wp, KIND=wp )
+
+                   DO  j = 1, (ny+1)/2 - 1
+                      ar_tmp(j,i,k) = CMPLX( ar(j,i,k), ar(ny+1-j,i,k),        &
+                                             KIND=wp )
+                   ENDDO
+                   ar_tmp((ny+1)/2,i,k) = CMPLX( ar((ny+1)/2,i,k), 0.0_wp,     &
+                                                 KIND=wp )
+
+                ENDDO
+             ENDDO
+
+             !$ACC HOST_DATA USE_DEVICE(ar, ar_tmp)
+             CALL CUFFTEXECZ2D( plan_yi, ar_tmp, ar )
+             !$ACC END HOST_DATA
+
+          ENDIF
+
 #endif