From 6bab98c86c89db733e64b34b70a4c8c05a713df8 Mon Sep 17 00:00:00 2001
From: Philippe WAUTELET <philippe.wautelet@aero.obs-mip.fr>
Date: Tue, 19 Dec 2023 14:52:32 +0100
Subject: [PATCH] Philippe 19/12/2023: OpenACC: developments for KTEST/003_KW78
 (work in progress)

---
 src/MNH/fast_terms.f90     | 57 ++++++++++++++++++++++++++++++++++++--
 src/MNH/ini_modeln.f90     |  2 +-
 src/MNH/pressurez.f90      | 52 ++++++++++++++++++++++++++++------
 src/MNH/resolved_cloud.f90 |  2 +-
 4 files changed, 101 insertions(+), 12 deletions(-)

diff --git a/src/MNH/fast_terms.f90 b/src/MNH/fast_terms.f90
index 58308a500..f618e6aea 100644
--- a/src/MNH/fast_terms.f90
+++ b/src/MNH/fast_terms.f90
@@ -163,11 +163,17 @@ USE MODD_LUNIT_n, ONLY: TLUOUT
 USE MODD_PARAMETERS
 
 use mode_budget,     only: Budget_store_end, Budget_store_init
+#ifdef MNH_OPENACC
+USE MODE_MNH_ZWORK,   ONLY: MNH_MEM_GET, MNH_MEM_POSITION_PIN, MNH_MEM_RELEASE
+#endif
 use mode_mppdb
 #ifdef MNH_OPENACC
 use mode_msg
 #endif
 
+#if defined(MNH_BITREP) || defined(MNH_BITREP_OMP)
+USE MODI_BITREP
+#endif
 USE MODI_CONDENS
 USE MODI_GET_HALO
 !
@@ -214,13 +220,25 @@ REAL, DIMENSION(:,:,:),   INTENT(OUT)   :: PCLDFR  ! Cloud fraction
 !
 REAL  :: ZEPS  ! Mv/Md
 REAL  :: ZTL,ZTHETA,ZRVS
-REAL, DIMENSION(SIZE(PRHODJ,1),SIZE(PRHODJ,2),SIZE(PRHODJ,3)) &
+#ifndef MNH_OPENACC
+LOGICAL, DIMENSION(:,:,:), ALLOCATABLE :: GWORK
+REAL,    DIMENSION(:,:,:), ALLOCATABLE &
+                         :: ZEXNS,&      ! guess of the Exner function at t+1
+                            ZT,   &      ! guess of the temperature at t+1
+                            ZCPH, &      ! guess of the CPh for the mixing
+                            ZLV,  &      ! guess of the Lv at t+1
+                            ZW1,ZW2,ZW3  ! Work arrays for intermediate
+                                         ! fields
+#else
+LOGICAL, DIMENSION(:,:,:), POINTER, CONTIGUOUS :: GWORK
+REAL,    DIMENSION(:,:,:), POINTER, CONTIGUOUS &
                          :: ZEXNS,&      ! guess of the Exner function at t+1
                             ZT,   &      ! guess of the temperature at t+1
                             ZCPH, &      ! guess of the CPh for the mixing
                             ZLV,  &      ! guess of the Lv at t+1
                             ZW1,ZW2,ZW3  ! Work arrays for intermediate
                                          ! fields
+#endif
 !
 INTEGER             :: IRESP      ! Return code of FM routines
 INTEGER             :: ILENG      ! Length of comment string in LFIFM file
@@ -229,7 +247,6 @@ INTEGER             :: ILENCH     ! Length of comment string in LFIFM file
 INTEGER             :: JK         ! Var for vertical DO loops
 INTEGER             :: JITER,ITERMAX  ! iterative loop for first order adjustment
 INTEGER             :: ILUOUT     ! Logical unit of output listing 
-LOGICAL,DIMENSION(SIZE(PRCS,1),SIZE(PRCS,2),SIZE(PRCS,3))::GWORK
 !-------------------------------------------------------------------------------
 #ifdef MNH_OPENACC
 call Print_msg( NVERB_WARNING, 'GEN', 'FAST_TERMS', 'OpenACC: being implemented' )
@@ -264,10 +281,38 @@ ELSE
   ITERMAX=1
 END IF
 
+#ifndef MNH_OPENACC
+ALLOCATE( GWORK(SIZE( PRCS,   1 ), SIZE( PRCS,   2 ), SIZE( PRCS,   3 )) )
+ALLOCATE( ZEXNS(SIZE( PRHODJ ,1 ), SIZE( PRHODJ, 2 ), SIZE( PRHODJ, 3 )) )
+ALLOCATE( ZT   (SIZE( PRHODJ ,1 ), SIZE( PRHODJ, 2 ), SIZE( PRHODJ, 3 )) )
+ALLOCATE( ZCPH (SIZE( PRHODJ ,1 ), SIZE( PRHODJ, 2 ), SIZE( PRHODJ, 3 )) )
+ALLOCATE( ZLV  (SIZE( PRHODJ ,1 ), SIZE( PRHODJ, 2 ), SIZE( PRHODJ, 3 )) )
+ALLOCATE( ZW1  (SIZE( PRHODJ ,1 ), SIZE( PRHODJ, 2 ), SIZE( PRHODJ, 3 )) )
+ALLOCATE( ZW2  (SIZE( PRHODJ ,1 ), SIZE( PRHODJ, 2 ), SIZE( PRHODJ, 3 )) )
+ALLOCATE( ZW3  (SIZE( PRHODJ ,1 ), SIZE( PRHODJ, 2 ), SIZE( PRHODJ, 3 )) )
+#else
+!Pin positions in the pools of MNH memory
+CALL MNH_MEM_POSITION_PIN( 'FAST_TERMS' )
+
+CALL MNH_MEM_GET( GWORK, SIZE( PRCS,   1 ), SIZE( PRCS,   2 ), SIZE( PRCS,   3 ) )
+CALL MNH_MEM_GET( ZEXNS, SIZE( PRHODJ, 1 ), SIZE( PRHODJ, 2 ), SIZE( PRHODJ, 3 ) )
+CALL MNH_MEM_GET( ZT,    SIZE( PRHODJ, 1 ), SIZE( PRHODJ, 2 ), SIZE( PRHODJ, 3 ) )
+CALL MNH_MEM_GET( ZCPH,  SIZE( PRHODJ, 1 ), SIZE( PRHODJ, 2 ), SIZE( PRHODJ, 3 ) )
+CALL MNH_MEM_GET( ZLV,   SIZE( PRHODJ, 1 ), SIZE( PRHODJ, 2 ), SIZE( PRHODJ, 3 ) )
+CALL MNH_MEM_GET( ZW1,   SIZE( PRHODJ, 1 ), SIZE( PRHODJ, 2 ), SIZE( PRHODJ, 3 ) )
+CALL MNH_MEM_GET( ZW2,   SIZE( PRHODJ, 1 ), SIZE( PRHODJ, 2 ), SIZE( PRHODJ, 3 ) )
+CALL MNH_MEM_GET( ZW3,   SIZE( PRHODJ, 1 ), SIZE( PRHODJ, 2 ), SIZE( PRHODJ, 3 ) )
+#endif
+
 if ( lbudget_rv ) call Budget_store_init( tbudgets(NBUDGET_RV), 'COND', prvs(:, :, :) * prhodj(:, :, :) )
 if ( lbudget_rc ) call Budget_store_init( tbudgets(NBUDGET_RC), 'COND', prcs(:, :, :) * prhodj(:, :, :) )
 if ( lbudget_th ) call Budget_store_init( tbudgets(NBUDGET_TH), 'COND', pths(:, :, :) * prhodj(:, :, :) )
 
+!$acc data present( PRHODJ, PSIGS, PPABST, PRVT, PRCT, PRVS, PRCS, PCF_MF, PRC_MF, PTHS, PSRCS, PCLDFR ) &
+!$acc&     present( GWORK, ZEXNS, ZT, ZCPH, ZLV, ZW1, ZW2, ZW3 )
+
+!$acc data present( PRRS ) if( KRR==3 )
+
 !-------------------------------------------------------------------------------
 !
 !*       2.     COMPUTE QUANTITIES WITH THE GUESS OF THE FUTURE INSTANT
@@ -465,6 +510,14 @@ ELSE
 !$acc end kernels
   END IF
 ENDIF
+
+!$acc end data
+!$acc end data
+
+#ifdef MNH_OPENACC
+!Release all memory allocated with MNH_MEM_GET calls since last call to MNH_MEM_POSITION_PIN
+CALL MNH_MEM_RELEASE( 'FAST_TERMS' )
+#endif
 !
 !
 !
diff --git a/src/MNH/ini_modeln.f90 b/src/MNH/ini_modeln.f90
index 77d897d77..9fb03b3af 100644
--- a/src/MNH/ini_modeln.f90
+++ b/src/MNH/ini_modeln.f90
@@ -1046,7 +1046,7 @@ ALLOCATE(XPHIT(IIU,IJU,IKU))
 ALLOCATE(XRHODREF(IIU,IJU,IKU))
 ALLOCATE(XTHVREF(IIU,IJU,IKU))
 ALLOCATE(XEXNREF(IIU,IJU,IKU))
-!$acc enter data create( XRHODREF, XTHVREF, XEXNREF )
+!$acc enter data create( XPHIT, XRHODREF, XTHVREF, XEXNREF )
 ALLOCATE(XRHODJ(IIU,IJU,IKU))
 !$acc enter data create(XRHODJ)
 IF (CEQNSYS=='DUR' .AND. LUSERV) THEN
diff --git a/src/MNH/pressurez.f90 b/src/MNH/pressurez.f90
index 6bcc60244..69ed18b44 100644
--- a/src/MNH/pressurez.f90
+++ b/src/MNH/pressurez.f90
@@ -461,7 +461,8 @@ IF (MPPDB_INITIALIZED) THEN
 END IF
 !$acc data present( PRHODJ, PDXX, PDYY, PDZZ, PDZX, PDZY, PTHT, PRT, PRHODREF, PTHVREF, PRVREF, PEXNREF ) &
 !$acc &    present( PRUS, PRVS, PRWS, PPABST ) &
-!$acc &    present( PRHOT, PAF, PBF, PCF, PTRIGSX, PTRIGSY, KIFAXX, KIFAXY, PBFB, PBF_SXP2_YP1_Z )
+!$acc &    present( PRHOT, PAF, PBF, PCF, PTRIGSX, PTRIGSY, KIFAXX, KIFAXY, PBFB, PBF_SXP2_YP1_Z ) &
+!$acc &    present( XPHIT )
 
 !-------------------------------------------------------------------------------
 !
@@ -517,7 +518,10 @@ CALL MNH_MEM_GET( ZMXM_PRHODJ, IIU, IJU, IKU )
 CALL MNH_MEM_GET( ZMYM_PRHODJ, IIU, IJU, IKU )
 CALL MNH_MEM_GET( ZMZM_PRHODJ, IIU, IJU, IKU )
 #endif
-!
+
+!$acc data present( ZDV_SOURCE, ZTHETAV, ZPHIT, ZPABS_S, ZPABS_N, ZPABS_E, ZPABS_W ) &
+!$acc&     present( ZPRHODJ, ZGZ_M_W, ZMXM_PRHODJ, ZMYM_PRHODJ, ZMZM_PRHODJ )
+
 !$acc kernels
 ZPABS_S(:,:) = 0.
 ZPABS_N(:,:) = 0.
@@ -687,12 +691,23 @@ ELSEIF(CEQNSYS=='LHE') THEN
   CALL PRINT_MSG( NVERB_FATAL, 'GEN', 'PRESSUREZ', 'OpenACC: CEQNSYS=LHE not yet ported' )
 #endif
   IF ( .NOT. LOCEAN) THEN
-    ZPHIT(:,:,:)= ((PPABST(:,:,:)/XP00)**(XRD/XCPD)-PEXNREF(:,:,:))   &
-                 * XCPD * PTHVREF(:,:,:)
+  !$acc kernels
+  !$mnh_do_concurrent ( JI=1:IIU,JJ=1:IJU,JK=1:IKU )
+#if !defined(MNH_BITREP) && !defined(MNH_BITREP_OMP)
+    ZPHIT(JI,JJ,JK)= ((PPABST(JI,JJ,JK)/XP00)**(XRD/XCPD)-PEXNREF(JI,JJ,JK))   &
+                 * XCPD * PTHVREF(JI,JJ,JK)
+#else
+    ZPHIT(JI,JJ,JK)= ( BR_POW( PPABST(JI,JJ,JK) / XP00, XRD / XCPD ) - PEXNREF(JI,JJ,JK) ) &
+                 * XCPD * PTHVREF(JI,JJ,JK)
+#endif
+  !$mnh_end_do()
+  !$acc end kernels
   ELSE
     ! Field at T- DT for LHE anelastic approx
     ! not used in ocean case (flat LHE)
+    !$acc kernels
     ZPHIT(:,:,:)=0.
+    !$acc end kernels
   END IF
 !
 END IF
@@ -904,8 +919,20 @@ IF(CEQNSYS=='MAE' .OR. CEQNSYS=='DUR') THEN
   !$acc end kernels
 #endif
 ELSEIF(CEQNSYS=='LHE') THEN
+#ifndef MNH_OPENACC
   PRUS = PRUS - MXM(PRHODJ) * ZDV_SOURCE
   PRWS = PRWS - MZM(PRHODJ) * GZ_M_W(1,IKU,1,ZPHIT,PDZZ)
+#else
+  CALL MXM_DEVICE(PRHODJ, ZMXM_PRHODJ)
+  CALL MZM_DEVICE(PRHODJ, ZMZM_PRHODJ)
+  CALL GZ_M_W_DEVICE(1,IKU,1,ZPHIT,PDZZ,ZGZ_M_W)
+  !$acc kernels
+  !dir$ concurrent
+  PRUS(:,:,:) = PRUS(:,:,:) - ZMXM_PRHODJ(:,:,:) * ZDV_SOURCE(:,:,:)
+  !dir$ concurrent
+  PRWS(:,:,:) = PRWS(:,:,:) - ZMZM_PRHODJ(:,:,:) * ZGZ_M_W(:,:,:)
+  !$acc end kernels
+#endif
 END IF
 !
 IF(.NOT. L2D) THEN
@@ -1096,11 +1123,19 @@ IF ((ZMAX_ll > 1.E-12) .AND. KTCOUNT >0 ) THEN
     IF (.NOT. LOCEAN) THEN
       !$acc kernels
       ! Deep atmosphere case : computing of PI fluctuation ; ZPHI0 (computed in P_ABS routine) is added
-      XPHIT(:,:,:) = (ZPHIT+ZPHI0)/(XCPD*PTHVREF)
-      ! Absolute Pressure 
+      XPHIT(:,:,:) = (ZPHIT(:,:,:)+ZPHI0)/(XCPD*PTHVREF(:,:,:))
+      ! Absolute Pressure
+#if !defined(MNH_BITREP) && !defined(MNH_BITREP_OMP)
       PPABST(:,:,:)=XP00*(XPHIT(:,:,:)+PEXNREF)**(XCPD/XRD)
+#else
+      PPABST(:,:,:)=XP00 * BR_POW( XPHIT(:,:,:) + PEXNREF(:,:,:), XCPD / XRD )
+#endif
       ! Computing press fluctuation
+#if !defined(MNH_BITREP) && !defined(MNH_BITREP_OMP)
       XPHIT(:,:,:) = PPABST(:,:,:) - XP00*PEXNREF**(XCPD/XRD)
+#else
+      XPHIT(:,:,:) = PPABST(:,:,:) - XP00 * BR_POW( PEXNREF, XCPD / XRD )
+#endif
      !$acc end kernels
     ELSE
      !$acc kernels
@@ -1173,6 +1208,9 @@ IF ((ZMAX_ll > 1.E-12) .AND. KTCOUNT >0 ) THEN
 !
 END IF
 
+!$acc end data
+!$acc end data
+
 IF (MPPDB_INITIALIZED) THEN
   !Check all INOUT arrays
   CALL MPPDB_CHECK(PRUS,"PRESSUREZ end:PRUS")
@@ -1185,8 +1223,6 @@ END IF
 CALL MNH_MEM_RELEASE( 'PRESSUREZ' )
 #endif
 
-!$acc end data
-
 !-------------------------------------------------------------------------------
 !
 END SUBROUTINE PRESSUREZ
diff --git a/src/MNH/resolved_cloud.f90 b/src/MNH/resolved_cloud.f90
index 071c83e1f..a59a4eaff 100644
--- a/src/MNH/resolved_cloud.f90
+++ b/src/MNH/resolved_cloud.f90
@@ -1322,7 +1322,7 @@ IF (MPPDB_INITIALIZED) THEN
   CALL MPPDB_CHECK(PHLI_HCF,"RESOLVED_CLOUD end:PHLI_HCF")
   !Check all OUT arrays
   CALL MPPDB_CHECK(PSRCS,  "RESOLVED_CLOUD end:PSRCS")
-  CALL MPPDB_CHECK(PRAINFR,"RESOLVED_CLOUD end:PRAINFR")
+  IF (HCLOUD(1:3)=='ICE') CALL MPPDB_CHECK(PRAINFR,"RESOLVED_CLOUD end:PRAINFR")
 END IF
 
 !$acc end data
-- 
GitLab