Juan 14/09/2020 : MNH/rain_ice_fast_rg.f90 , OpenACC optimisation -> DO...

Juan 14/09/2020 : MNH/rain_ice_fast_rg.f90 , OpenACC optimisation -> DO CONCURRENT for BR_POW parallelization

src/MNH/rain_ice_fast_rg.f90

@@ -101,6 +101,7 @@ REAL,    DIMENSION(:),   ALLOCATABLE :: ZVEC1,ZVEC2,ZVEC3 ! Work vectors for int
 REAL,    DIMENSION(:),   ALLOCATABLE :: ZVECLBDAG, ZVECLBDAR, ZVECLBDAS
 REAL,    DIMENSION(:,:), ALLOCATABLE :: ZZW1              ! Work arrays
 !
+INTEGER                              :: JLU
 !-------------------------------------------------------------------------------
 !
 ! IN variables
@@ -159,6 +160,8 @@ IF (MPPDB_INITIALIZED) THEN
   CALL MPPDB_CHECK(XKER_SDRYG,"RAIN_ICE_FAST_RG beg:XKER_RDRYG")
 END IF
 !
+JLU = size(PRHODREF) 
+!
 ALLOCATE( I1   (size(PRHODREF)) )
 ALLOCATE( GWORK(size(PRHODREF)) )
 ALLOCATE( ZZW  (size(PRHODREF)) )
@@ -170,28 +173,36 @@ ALLOCATE( ZZW1 (size(PRHODREF),7) )
 !*       6.1    rain contact freezing
 !
 !$acc kernels
-  GWORK(:) = PRIT(:)>XRTMIN(4) .AND. PRRT(:)>XRTMIN(3) .AND. PRIS(:)>0.0 .AND. PRRS(:)>0.0
-  WHERE( GWORK(:) )
+GWORK(:) = PRIT(:)>XRTMIN(4) .AND. PRRT(:)>XRTMIN(3) .AND. PRIS(:)>0.0 .AND. PRRS(:)>0.0
 #ifndef MNH_BITREP
+  WHERE( GWORK(:) )
     ZZW1(:,3) = MIN( PRIS(:),XICFRR * PRIT(:)                & ! RICFRRG
                                     * PLBDAR(:)**XEXICFRR    &
                                     * PRHODREF(:)**(-XCEXVT) )
     ZZW1(:,4) = MIN( PRRS(:),XRCFRI * PCIT(:)                & ! RRCFRIG
                                     * PLBDAR(:)**XEXRCFRI    &
                                     * PRHODREF(:)**(-XCEXVT-1.) )
-#else
-    ZZW1(:,3) = MIN( PRIS(:),XICFRR * PRIT(:)                     & ! RICFRRG
-                                    * BR_POW(PLBDAR(:),XEXICFRR)  &
-                                    * BR_POW(PRHODREF(:),-XCEXVT) )
-    ZZW1(:,4) = MIN( PRRS(:),XRCFRI * PCIT(:)                     & ! RRCFRIG
-                                    * BR_POW(PLBDAR(:),XEXRCFRI)  &
-                                    * BR_POW(PRHODREF(:),-XCEXVT-1.) )
-#endif
     PRIS(:) = PRIS(:) - ZZW1(:,3)
     PRRS(:) = PRRS(:) - ZZW1(:,4)
     PRGS(:) = PRGS(:) + ZZW1(:,3)+ZZW1(:,4)
     PTHS(:) = PTHS(:) + ZZW1(:,4)*(PLSFACT(:)-PLVFACT(:)) ! f(L_f*RRCFRIG)
-  END WHERE
+ END WHERE
+#else
+ DO CONCURRENT ( JL=1:JLU )
+    IF ( GWORK(JL) ) THEN
+       ZZW1(JL,3) = MIN( PRIS(JL),XICFRR * PRIT(JL)                     & ! RICFRRG
+            * BR_POW(PLBDAR(JL),XEXICFRR)  &
+            * BR_POW(PRHODREF(JL),-XCEXVT) )
+       ZZW1(JL,4) = MIN( PRRS(JL),XRCFRI * PCIT(JL)                     & ! RRCFRIG
+            * BR_POW(PLBDAR(JL),XEXRCFRI)  &
+            * BR_POW(PRHODREF(JL),-XCEXVT-1.) )
+       PRIS(JL) = PRIS(JL) - ZZW1(JL,3)
+       PRRS(JL) = PRRS(JL) - ZZW1(JL,4)
+       PRGS(JL) = PRGS(JL) + ZZW1(JL,3)+ZZW1(JL,4)
+       PTHS(JL) = PTHS(JL) + ZZW1(JL,4)*(PLSFACT(JL)-PLVFACT(JL)) ! f(L_f*RRCFRIG)
+    END IF
+ END DO ! CONCURRENT
+#endif 
 !$acc end kernels
 IF (MPPDB_INITIALIZED) THEN
    CALL MPPDB_CHECK(PRRS,"RAIN_ICE_FAST_RG 6.1:PRRS")
@@ -451,18 +462,12 @@ END IF
 !
   PRWETG(:) = 0.0
   GWORK(:) = PRGT(:)>XRTMIN(6)
+#ifndef MNH_BITREP  
   WHERE( GWORK(:) )
-#ifndef MNH_BITREP
     ZZW1(:,5) = MIN( PRIS(:),                                    &
                 ZZW1(:,2) / (XCOLIG*EXP(XCOLEXIG*(PZT(:)-XTT)) ) ) ! RIWETG
     ZZW1(:,6) = MIN( PRSS(:),                                    &
                 ZZW1(:,3) / (XCOLSG*EXP(XCOLEXSG*(PZT(:)-XTT)) ) ) ! RSWETG
-#else
-    ZZW1(:,5) = MIN( PRIS(:),                                    &
-                ZZW1(:,2) / (XCOLIG*BR_EXP(XCOLEXIG*(PZT(:)-XTT)) ) ) ! RIWETG
-    ZZW1(:,6) = MIN( PRSS(:),                                    &
-                ZZW1(:,3) / (XCOLSG*BR_EXP(XCOLEXSG*(PZT(:)-XTT)) ) ) ! RSWETG
-#endif
 !
     ZZW(:) = PRVT(:)*PPRES(:)/((XMV/XMD)+PRVT(:)) ! Vapor pressure
     ZZW(:) =   PKA(:)*(XTT-PZT(:)) +                              &
@@ -472,17 +477,36 @@ END IF
 ! compute RWETG
 !
     PRWETG(:)=MAX( 0.0,                                               &
-#ifndef MNH_BITREP
                  ( ZZW(:) * ( X0DEPG*       PLBDAG(:)**XEX0DEPG +     &
                               X1DEPG*PCJ(:)*PLBDAG(:)**XEX1DEPG ) +   &
-#else
-                 ( ZZW(:) * ( X0DEPG*       BR_POW(PLBDAG(:),XEX0DEPG) +     &
-                              X1DEPG*PCJ(:)*BR_POW(PLBDAG(:),XEX1DEPG) ) +   &
-#endif
                  ( ZZW1(:,5)+ZZW1(:,6) ) *                            &
                  ( PRHODREF(:)*(XLMTT+(XCI-XCL)*(XTT-PZT(:)))   ) ) / &
                             ( PRHODREF(:)*(XLMTT-XCL*(XTT-PZT(:))) )   )
-  END WHERE
+ END WHERE
+#else
+ DO CONCURRENT ( JL=1:JLU )
+    IF ( GWORK(JL) ) THEN
+       ZZW1(JL,5) = MIN( PRIS(JL),                                    &
+            ZZW1(JL,2) / (XCOLIG*BR_EXP(XCOLEXIG*(PZT(JL)-XTT)) ) ) ! RIWETG
+       ZZW1(JL,6) = MIN( PRSS(JL),                                    &
+            ZZW1(JL,3) / (XCOLSG*BR_EXP(XCOLEXSG*(PZT(JL)-XTT)) ) ) ! RSWETG
+       !
+       ZZW(JL) = PRVT(JL)*PPRES(JL)/((XMV/XMD)+PRVT(JL)) ! Vapor pressure
+       ZZW(JL) =   PKA(JL)*(XTT-PZT(JL)) +                              &
+            ( PDV(JL)*(XLVTT + ( XCPV - XCL ) * ( PZT(JL) - XTT )) &
+            *(XESTT-ZZW(JL))/(XRV*PZT(JL))           )
+       !
+       ! compute RWETG
+       !
+       PRWETG(JL)=MAX( 0.0,                                               &
+            ( ZZW(JL) * ( X0DEPG*       BR_POW(PLBDAG(JL),XEX0DEPG) +     &
+            X1DEPG*PCJ(JL)*BR_POW(PLBDAG(JL),XEX1DEPG) ) +   &
+            ( ZZW1(JL,5)+ZZW1(JL,6) ) *                            &
+            ( PRHODREF(JL)*(XLMTT+(XCI-XCL)*(XTT-PZT(JL)))   ) ) / &
+            ( PRHODREF(JL)*(XLMTT-XCL*(XTT-PZT(JL))) )   )
+    END IF
+ END DO ! CONCURRENT
+#endif 
 !
 !*       6.4    Select Wet or Dry case
 !
@@ -610,7 +634,8 @@ END IF
 !*       6.5    Melting of the graupeln
 !
 !$acc kernels
-   GWORK(:) = PRGT(:)>XRTMIN(6) .AND. PRGS(:)>0.0 .AND. PZT(:)>XTT
+  GWORK(:) = PRGT(:)>XRTMIN(6) .AND. PRGS(:)>0.0 .AND. PZT(:)>XTT
+#ifndef MNH_BITREP  
    WHERE( GWORK(:) )
     ZZW(:) = PRVT(:)*PPRES(:)/((XMV/XMD)+PRVT(:)) ! Vapor pressure
     ZZW(:) =  PKA(:)*(XTT-PZT(:)) +                                 &
@@ -620,20 +645,37 @@ END IF
 ! compute RGMLTR
 !
     ZZW(:)  = MIN( PRGS(:), MAX( 0.0,( -ZZW(:) *                     &
-#ifndef MNH_BITREP
                            ( X0DEPG*       PLBDAG(:)**XEX0DEPG +     &
                              X1DEPG*PCJ(:)*PLBDAG(:)**XEX1DEPG ) -   &
-#else
-                           ( X0DEPG*       BR_POW(PLBDAG(:),XEX0DEPG) +     &
-                             X1DEPG*PCJ(:)*BR_POW(PLBDAG(:),XEX1DEPG) ) -   &
-#endif
                                      ( ZZW1(:,1)+ZZW1(:,4) ) *       &
                               ( PRHODREF(:)*XCL*(XTT-PZT(:))) ) /    &
                                              ( PRHODREF(:)*XLMTT ) ) )
     PRRS(:) = PRRS(:) + ZZW(:)
     PRGS(:) = PRGS(:) - ZZW(:)
     PTHS(:) = PTHS(:) - ZZW(:)*(PLSFACT(:)-PLVFACT(:)) ! f(L_f*(-RGMLTR))
-  END WHERE
+ END WHERE
+#else
+ DO CONCURRENT ( JL=1:JLU )
+    IF ( GWORK(JL) ) THEN
+       ZZW(JL) = PRVT(JL)*PPRES(JL)/((XMV/XMD)+PRVT(JL)) ! Vapor pressure
+       ZZW(JL) =  PKA(JL)*(XTT-PZT(JL)) +                                 &
+            ( PDV(JL)*(XLVTT + ( XCPV - XCL ) * ( PZT(JL) - XTT )) &
+            *(XESTT-ZZW(JL))/(XRV*PZT(JL))             )
+       !
+       ! compute RGMLTR
+       !
+       ZZW(JL)  = MIN( PRGS(JL), MAX( 0.0,( -ZZW(JL) *                     &
+            ( X0DEPG*       BR_POW(PLBDAG(JL),XEX0DEPG) +     &
+            X1DEPG*PCJ(JL)*BR_POW(PLBDAG(JL),XEX1DEPG) ) -   &
+            ( ZZW1(JL,1)+ZZW1(JL,4) ) *       &
+            ( PRHODREF(JL)*XCL*(XTT-PZT(JL))) ) /    &
+            ( PRHODREF(JL)*XLMTT ) ) )
+       PRRS(JL) = PRRS(JL) + ZZW(JL)
+       PRGS(JL) = PRGS(JL) - ZZW(JL)
+       PTHS(JL) = PTHS(JL) - ZZW(JL)*(PLSFACT(JL)-PLVFACT(JL)) ! f(L_f*(-RGMLTR))
+    END IF
+ END DO ! CONCURRENT
+#endif 
 !$acc end kernels
 IF (MPPDB_INITIALIZED) THEN
    CALL MPPDB_CHECK(PRRS,"RAIN_ICE_FAST_RG 6.5:PRRS")