Philippe 18/07/2019: OpenACC: add support for HTURBLEN=DELT

src/MNH/turb.f90

@@ -712,9 +712,6 @@ SELECT CASE (HTURBLEN)
 !           -------------------
 !
   CASE ('DELT')
-#ifdef _OPENACC
-    call Print_msg( NVERB_WARNING, 'GEN', 'TURB', 'OpenACC: HTURBLEN=DELT not yet tested' )
-#endif
     CALL DELT(KKA,KKU,KKL,IKB, IKE,IKTB, IKTE,ORMC01,HTURBDIM,PDXX, PDYY,PZZ,PDIRCOSZW,PLEM)
 !
 !*      3.4 Deardorff mixing length
@@ -1617,13 +1614,11 @@ ELSE
 !*         3.2 Delta mixing length
 !           -------------------
   CASE ('DELT')
-!     CALL DELT(ZLM_CLOUD)
     CALL DELT(KKA,KKU,KKL,IKB, IKE,IKTB, IKTE,ORMC01,HTURBDIM,PDXX, PDYY,PZZ,PDIRCOSZW,ZLM_CLOUD)
 !
 !*         3.3 Deardorff mixing length
 !           -----------------------
   CASE ('DEAR')
-!     CALL DEAR(ZLM_CLOUD)
     CALL DEAR(KKA,KKU,KKL,KRR, KRRI, IKB, IKE,IKTB, IKTE, &
                 ORMC01,HTURBDIM,PDXX, PDYY, PDZZ,PZZ,PDIRCOSZW,PTHLT,PTHVREF,PTKET,PSRCT,PRT,&
                 ZLOCPEXNM,ZATHETA, ZAMOIST, ZLM_CLOUD)
@@ -1708,13 +1703,18 @@ SUBROUTINE DELT(KKA,KKU,KKL,KKB, KKE,KKTB, KKTE,ORMC01,HTURBDIM,PDXX, PDYY,PZZ,P
 !!    -------------
 !!      Original   01/05
 !!
+!  P. Wautelet 18/07/2019: add OpenACC directives
 !-------------------------------------------------------------------------------
 !
 !*       0.    DECLARATIONS
 !              ------------
 use modd_conf, only: l2d
 
+#ifndef _OPENACC
 USE MODI_SHUMAN
+#else
+USE MODI_SHUMAN_DEVICE
+#endif
 
 implicit none
 !
@@ -1737,6 +1737,8 @@ REAL, DIMENSION(:,:),   INTENT(IN)      ::  PDIRCOSZW
 ! Director Cosinus along x, y and z directions at surface w-point
 REAL, DIMENSION(:,:,:), INTENT(OUT)   :: PLM
 !
+!$acc declare present(PDXX, PDYY, PZZ, PDIRCOSZW, PLM)
+!
 !*       0.2   Declarations of local variables
 !
 
@@ -1744,38 +1746,73 @@ integer :: ji, jj, jk
 REAL                :: ZALPHA       ! proportionnality constant between Dz/2 and
 !                                   ! BL89 mixing length near the surface
 REAL                :: ZD           ! distance to the surface
+#ifdef _OPENACC
+real, dimension(:,:,:), allocatable :: ztmp1_device, ztmp2_device
+!$acc declare create(ztmp1_device, ztmp2_device)
+#endif
 !
 !-------------------------------------------------------------------------------
 !
 #ifdef _OPENACC
-    call Print_msg( NVERB_FATAL, 'GEN', 'TURB', 'OpenACC: DELT not yet implemented' )
+allocate( ztmp1_device( size( pdxx, 1 ), size( pdxx, 2 ), size( pdxx, 3 ) ) )
+allocate( ztmp2_device( size( pdxx, 1 ), size( pdxx, 2 ), size( pdxx, 3 ) ) )
 #endif
+
+!$acc kernels
 DO JK = KKTB,KKTE ! 1D turbulence scheme
   PLM(:,:,JK) = PZZ(:,:,JK+KKL) - PZZ(:,:,JK)
 END DO
 PLM(:,:,KKU) = PLM(:,:,KKE)
 PLM(:,:,KKA) = PZZ(:,:,KKB) - PZZ(:,:,KKA)
+!$acc end kernels
 IF ( HTURBDIM /= '1DIM' ) THEN  ! 3D turbulence scheme
-  IF ( L2D) THEN
+  IF ( L2D ) THEN
+#ifndef _OPENACC
     PLM(:,:,:) = SQRT( PLM(:,:,:)*MXF(PDXX(:,:,:)) )
+#else
+    CALL MXF_DEVICE( PDXX, ZTMP1_DEVICE )
+!$acc kernels
+    PLM(:,:,:) = SQRT( PLM(:,:,:) * ZTMP1_DEVICE )
+!$acc end kernels
+#endif
   ELSE
+#ifndef _OPENACC
+#ifndef MNH_BITREP
     PLM(:,:,:) = (PLM(:,:,:)*MXF(PDXX(:,:,:))*MYF(PDYY(:,:,:)) ) ** (1./3.)
+#else
+    PLM(:,:,:) = BR_POW( PLM(:, :, : ) * MXF( PDXX(:, :, : ) ) * MYF( PDYY(:, :, : ) ), 1. / 3. )
+#endif
+#else
+    CALL MXF_DEVICE( PDXX, ZTMP1_DEVICE )
+    CALL MYF_DEVICE( PDYY, ZTMP2_DEVICE )
+!$acc kernels
+#ifndef MNH_BITREP
+    PLM(:,:,:) = ( PLM(:,:,:) * ZTMP1_DEVICE * ZTMP2_DEVICE ) ** (1./3.)
+#else
+    PLM(:,:,:) = BR_POW( PLM(:,:,:) * ZTMP1_DEVICE * ZTMP2_DEVICE, 1./3. )
+#endif
+!$acc end kernels
+#endif
   END IF
 END IF
 !
 !  mixing length limited by the distance normal to the surface
 !  (with the same factor as for BL89)
 !
+!$acc kernels
 IF (.NOT. ORMC01) THEN
-  ZALPHA=0.5**(-1.5)
+#ifndef MNH_BITREP
+  ZALPHA = 0.5**(-1.5)
+#else
+  ZALPHA = BR_POW( 0.5, -1.5 )
+#endif
   !
   DO JJ=1,SIZE(PLM,2)
     DO JI=1,SIZE(PLM,1)
       DO JK=KKTB,KKTE
-        ZD=ZALPHA*(0.5*(PZZ(JI,JJ,JK)+PZZ(JI,JJ,JK+KKL))&
-        -PZZ(JI,JJ,KKB)) *PDIRCOSZW(JI,JJ)
-        IF ( PLM(JI,JJ,JK)>ZD) THEN
-          PLM(JI,JJ,JK)=ZD
+        ZD = ZALPHA * ( 0.5 * ( PZZ(JI, JJ, JK ) + PZZ(JI, JJ, JK+KKL ) ) - PZZ(JI, JJ, KKB ) ) * PDIRCOSZW(JI, JJ )
+        IF ( PLM(JI,JJ,JK) > ZD ) THEN
+          PLM(JI,JJ,JK) = ZD
         ELSE
           EXIT
         ENDIF
@@ -1786,6 +1823,7 @@ END IF
 !
 PLM(:,:,KKA) = PLM(:,:,KKB  )
 PLM(:,:,KKU  ) = PLM(:,:,KKE)
+!$acc end kernels
 !
 END SUBROUTINE DELT
 
@@ -1797,7 +1835,7 @@ SUBROUTINE DEAR(KKA,KKU,KKL,KRR, KRRI, KKB, KKE,KKTB, KKTE, &
                 PLOCPEXNM,PATHETA, PAMOIST, PLM)
 !###################
 !!
-!!****  *DELT* routine to compute mixing length for DEARdorff case
+!!****  *DEAR* routine to compute mixing length for DEARdorff case
 !
 !!    AUTHOR
 !!    ------