diff --git a/src/MNH/ice4_compute_pdf.f90 b/src/MNH/ice4_compute_pdf.f90
index 2cc1f35ee01622b8abd1417feb8699214b4ac7f6..b95d226ed0f51a0908d1f745af69bf93631d9d57 100644
--- a/src/MNH/ice4_compute_pdf.f90
+++ b/src/MNH/ice4_compute_pdf.f90
@@ -10,7 +10,7 @@ SUBROUTINE ICE4_COMPUTE_PDF(KSIZE, HSUBG_AUCV, HSUBG_PR_PDF, &
IMPLICIT NONE
INTEGER, INTENT(IN) :: KSIZE
CHARACTER(LEN=4), INTENT(IN) :: HSUBG_AUCV ! Kind of Subgrid autoconversion method
-CHARACTER*80, INTENT(IN) :: HSUBG_PR_PDF ! pdf for subgrid precipitation
+CHARACTER(LEN=80), INTENT(IN) :: HSUBG_PR_PDF ! pdf for subgrid precipitation
REAL, DIMENSION(KSIZE), INTENT(IN) :: PRHODREF ! Reference density
REAL, DIMENSION(KSIZE), INTENT(IN) :: PRCT ! Cloud water m.r. at t
REAL, DIMENSION(KSIZE), INTENT(IN) :: PCF ! Cloud fraction
@@ -46,9 +46,9 @@ SUBROUTINE ICE4_COMPUTE_PDF(KSIZE, HSUBG_AUCV, HSUBG_PR_PDF, &
! ------------
!
!
-USE MODD_RAIN_ICE_DESCR
-USE MODD_RAIN_ICE_PARAM
-USE MODD_LUNIT_n, ONLY : TLUOUT
+USE MODD_RAIN_ICE_DESCR, ONLY: XRTMIN
+USE MODD_RAIN_ICE_PARAM, ONLY: XCRIAUTC
+!
USE MODE_MSG
!
#ifdef MNH_BITREP
@@ -63,7 +63,7 @@ IMPLICIT NONE
!
INTEGER, INTENT(IN) :: KSIZE
CHARACTER(LEN=4), INTENT(IN) :: HSUBG_AUCV ! Kind of Subgrid autoconversion method
-CHARACTER*80, INTENT(IN) :: HSUBG_PR_PDF ! pdf for subgrid precipitation
+CHARACTER(LEN=80), INTENT(IN) :: HSUBG_PR_PDF ! pdf for subgrid precipitation
REAL, DIMENSION(KSIZE), INTENT(IN) :: PRHODREF ! Reference density
REAL, DIMENSION(KSIZE), INTENT(IN) :: PRCT ! Cloud water m.r. at t
REAL, DIMENSION(KSIZE), INTENT(IN) :: PCF ! Cloud fraction
@@ -76,8 +76,13 @@ REAL, DIMENSION(KSIZE), INTENT(OUT) :: PHLC_LRC ! HLCLOUDS : LWC that is Low
! note that PRC = PHLC_HRC + PHLC_LRC
REAL, DIMENSION(KSIZE), INTENT(OUT) :: PRF ! Rain fraction
!
+!$acc declare present(KSIZE,HSUBG_AUCV,HSUBG_PR_PDF,PRHODREF,PRCT,PCF,PSIGMA_RC, &
+!$acc& PHLC_HCF,PHLC_LCF,PHLC_HRC,PHLC_LRC,PRF)
+!
!* 0.2 declaration of local variables
!
+INTEGER :: JI
+LOGICAL, DIMENSION(KSIZE) :: GWORK, GWORK2, GWORK3
REAL, DIMENSION(SIZE(PRHODREF)) :: ZRCRAUTC, & !RC value to begin rain formation =XCRIAUTC/RHODREF
ZHLC_RCMAX, & !HLCLOUDS : maximum value for RC in distribution
ZHLC_LRCLOCAL, & !HLCLOUDS : LWC that is Low LWC local in LCF
@@ -85,6 +90,11 @@ REAL, DIMENSION(SIZE(PRHODREF)) :: ZRCRAUTC, & !RC value to begin rain form
! note that ZRC/CF = ZHLC_HRCLOCAL+ ZHLC_LRCLOCAL
! = PHLC_HRC/HCF+ PHLC_LRC/LCF
REAL :: ZCOEFFRCM
+!
+!$acc declare create(GWORK,GWORK2,GWORK3,ZRCRAUTC,ZHLC_RCMAX,ZHLC_LRCLOCAL,ZHLC_HRCLOCAL,ZCOEFFRCM)
+!
+!$acc declare copyin(XCRIAUTC,XRTMIN)
+!
!-------------------------------------------------------------------------------
!
IF (MPPDB_INITIALIZED) THEN
@@ -96,16 +106,23 @@ IF (MPPDB_INITIALIZED) THEN
END IF
!
!Cloud water split between high and low content part is done according to autoconversion option
+!$acc kernels
ZRCRAUTC(:)=XCRIAUTC/PRHODREF(:) ! Autoconversion rc threshold
+!$acc end kernels
IF(HSUBG_AUCV=='NONE') THEN
+#if 0
+!PW: bug: crash at execution with PGI 18.10
+!$acc kernels
!Cloud water is entirely in low or high part
- WHERE(PRCT(:)>ZRCRAUTC(:))
+ GWORK(:) = PRCT(:)>ZRCRAUTC(:)
+ GWORK2(:) = PRCT(:)>XRTMIN(2)
+ WHERE(GWORK(:))
PHLC_HCF(:)=1.
PHLC_LCF(:)=0.
PHLC_HRC(:)=PRCT(:)
PHLC_LRC(:)=0.
PRF(:) =1.
- ELSEWHERE(PRCT(:)>XRTMIN(2))
+ ELSEWHERE(GWORK2(:))
PHLC_HCF(:)=0.
PHLC_LCF(:)=1.
PHLC_HRC(:)=0.
@@ -118,22 +135,51 @@ IF(HSUBG_AUCV=='NONE') THEN
PHLC_LRC(:)=0.
PRF(:) =0.
END WHERE
-
+!$acc end kernels
+#else
+!$acc kernels
+ DO JI = 1,KSIZE
+ IF (PRCT(JI)>ZRCRAUTC(JI)) THEN
+ PHLC_HCF(JI)=1.
+ PHLC_LCF(JI)=0.
+ PHLC_HRC(JI)=PRCT(JI)
+ PHLC_LRC(JI)=0.
+ PRF(JI) =1.
+ ELSE IF (PRCT(JI)>XRTMIN(2)) THEN
+ PHLC_HCF(JI)=0.
+ PHLC_LCF(JI)=1.
+ PHLC_HRC(JI)=0.
+ PHLC_LRC(JI)=PRCT(JI)
+ PRF(JI) =0.
+ ELSE
+ PHLC_HCF(JI)=0.
+ PHLC_LCF(JI)=0.
+ PHLC_HRC(JI)=0.
+ PHLC_LRC(JI)=0.
+ PRF(JI) =0.
+ END IF
+ END DO
+!$acc end kernels
+#endif
ELSEIF(HSUBG_AUCV=='CLFR') THEN
+!$acc kernels
!Cloud water is only in the cloudy part and entirely in low or high part
- WHERE(PCF(:) > 0. .AND. PRCT(:)>ZRCRAUTC(:)*PCF(:))
+ GWORK3(:) = PCF(:) > 0.
+ GWORK(:) = GWORK3(:) .AND. PRCT(:)>ZRCRAUTC(:)*PCF(:)
+ GWORK2(:) = GWORK3(:) .AND. PRCT(:)>XRTMIN(2)
+ WHERE(GWORK(:))
PHLC_HCF(:)=PCF(:)
PHLC_LCF(:)=0.
PHLC_HRC(:)=PRCT(:)
PHLC_LRC(:)=0.
PRF(:) =PCF(:)
- ELSEWHERE(PCF(:) > 0. .AND. PRCT(:)>XRTMIN(2))
+ ELSEWHERE(GWORK2(:))
PHLC_HCF(:)=0.
PHLC_LCF(:)=PCF(:)
PHLC_HRC(:)=0.0
PHLC_LRC(:)=PRCT(:)
PRF(:) =0.
- ELSEWHERE (PCF(:) > 0.)
+ ELSEWHERE (GWORK3(:))
PHLC_HCF(:)=0.
PHLC_LCF(:)=0.
PHLC_HRC(:)=0.
@@ -146,6 +192,7 @@ ELSEIF(HSUBG_AUCV=='CLFR') THEN
PHLC_LRC(:)=0.
PRF(:) =0.
END WHERE
+!$acc end kernels
ELSEIF(HSUBG_AUCV=='PDF ') THEN
!Cloud water is split between high and low part according to a PDF
@@ -155,15 +202,19 @@ ELSEIF(HSUBG_AUCV=='PDF ') THEN
! 'HLCISOTRIPDF' : isocele triangular PDF
! 'SIGM' : Redelsperger and Sommeria (1986)
IF(HSUBG_PR_PDF=='SIGM') THEN
+!$acc kernels
! Redelsperger and Sommeria (1986) but organised according to Turner (2011, 2012)
- WHERE (PRCT(:)>ZRCRAUTC(:)+PSIGMA_RC(:))
+ GWORK(:) = PRCT(:)> ZRCRAUTC(:)+PSIGMA_RC(:)
+ GWORK2(:) = PRCT(:)> (ZRCRAUTC(:)-PSIGMA_RC(:)) .AND. &
+ PRCT(:)<=(ZRCRAUTC(:)+PSIGMA_RC(:))
+ GWORK3(:) = PRCT(:)>XRTMIN(2) .AND. PCF(:)>0.
+ WHERE (GWORK(:))
PHLC_HCF(:)=1.
PHLC_LCF(:)=0.
PHLC_HRC(:)=PRCT(:)
PHLC_LRC(:)=0.
PRF(:) =1.
- ELSEWHERE(PRCT(:)> (ZRCRAUTC(:)-PSIGMA_RC(:)) .AND. &
- & PRCT(:)<=(ZRCRAUTC(:)+PSIGMA_RC(:)) )
+ ELSEWHERE(GWORK2(:))
PHLC_HCF(:)=(PRCT(:)+PSIGMA_RC(:)-ZRCRAUTC(:))/ &
&(2.*PSIGMA_RC(:))
PHLC_LCF(:)=MAX(0., PCF(:)-PHLC_HCF(:))
@@ -172,7 +223,7 @@ ELSEIF(HSUBG_AUCV=='PDF ') THEN
&(4.*PSIGMA_RC(:))
PHLC_LRC(:)=MAX(0., PRCT(:)-PHLC_HRC(:))
PRF(:) =PHLC_HCF(:)
- ELSEWHERE(PRCT(:)>XRTMIN(2) .AND. PCF(:)>0.)
+ ELSEWHERE(GWORK3(:))
PHLC_HCF(:)=0.
PHLC_LCF(:)=PCF(:)
PHLC_HRC(:)=0.
@@ -186,6 +237,7 @@ ELSEIF(HSUBG_AUCV=='PDF ') THEN
PRF(:) =0.
END WHERE
! Turner (2011, 2012)
+!$acc end kernels
ELSEIF(HSUBG_PR_PDF=='HLCRECTPDF' .OR. HSUBG_PR_PDF=='HLCISOTRIPDF' .OR. &
&HSUBG_PR_PDF=='HLCTRIANGPDF' .OR. HSUBG_PR_PDF=='HLCQUADRAPDF') THEN
! Calculate maximum value r_cM from PDF forms
@@ -196,24 +248,28 @@ ELSEIF(HSUBG_AUCV=='PDF ') THEN
ELSE IF(HSUBG_PR_PDF=='HLCQUADRAPDF') THEN
ZCOEFFRCM=4.
END IF
- WHERE(PRCT(:).GT.0. .AND. PCF(:).GT.0.)
+!$acc update device(ZCOEFFRCM)
+!$acc kernels
+ GWORK(:) = PRCT(:)>0. .AND. PCF(:)>0.
+ WHERE(GWORK(:))
ZHLC_RCMAX(:)=ZCOEFFRCM*PRCT(:)/PCF(:)
END WHERE
! Split available water and cloud fraction in two parts
! Calculate local mean values int he low and high parts for the 3 PDF forms:
+ GWORK(:) = PRCT(:)>0. .AND. PCF(:)>0. .AND. ZHLC_RCMAX(:)>ZRCRAUTC(:)
IF(HSUBG_PR_PDF=='HLCRECTPDF') THEN
- WHERE(PRCT(:).GT.0. .AND. PCF(:).GT.0. .AND. ZHLC_RCMAX(:).GT.ZRCRAUTC(:))
+ WHERE(GWORK(:))
ZHLC_LRCLOCAL(:)=0.5*ZRCRAUTC(:)
ZHLC_HRCLOCAL(:)=( ZHLC_RCMAX(:) + ZRCRAUTC(:))/2.0
END WHERE
ELSE IF(HSUBG_PR_PDF=='HLCTRIANGPDF') THEN
- WHERE(PRCT(:).GT.0. .AND. PCF(:).GT.0. .AND. ZHLC_RCMAX(:).GT.ZRCRAUTC(:))
+ WHERE(GWORK(:))
ZHLC_LRCLOCAL(:)=( ZRCRAUTC(:) *(3.0 * ZHLC_RCMAX(:) - 2.0 * ZRCRAUTC(:) ) ) &
/ (3.0 * (2.0 * ZHLC_RCMAX(:) - ZRCRAUTC(:) ) )
ZHLC_HRCLOCAL(:)=(ZHLC_RCMAX(:) + 2.0*ZRCRAUTC(:)) / 3.0
END WHERE
ELSE IF(HSUBG_PR_PDF=='HLCQUADRAPDF') THEN
- WHERE(PRCT(:).GT.0. .AND. PCF(:).GT.0. .AND. ZHLC_RCMAX(:).GT.ZRCRAUTC(:))
+ WHERE(GWORK(:))
#ifndef MNH_BITREP
ZHLC_LRCLOCAL(:)=(3.0 *ZRCRAUTC(:)**3 - 8.0 *ZRCRAUTC(:)**2 * ZHLC_RCMAX(:) &
+ 6.0*ZRCRAUTC(:) *ZHLC_RCMAX(:)**2 ) &
@@ -230,8 +286,9 @@ ELSEIF(HSUBG_AUCV=='PDF ') THEN
ZHLC_HRCLOCAL(:)=(ZHLC_RCMAX(:) + 3.0*ZRCRAUTC(:))/4.0
END WHERE
ELSE IF(HSUBG_PR_PDF=='HLCISOTRIPDF') THEN
- WHERE(PRCT(:).GT.0. .AND. PCF(:).GT.0. .AND. ZHLC_RCMAX(:).GT.ZRCRAUTC(:))
- WHERE((PRCT(:) / PCF(:)).LE.ZRCRAUTC(:))
+ GWORK2(:) = PRCT(:) <= ZRCRAUTC(:)*PCF(:)
+ WHERE(GWORK(:))
+ WHERE(GWORK2(:))
#ifndef MNH_BITREP
ZHLC_LRCLOCAL(:)=( (ZHLC_RCMAX(:))**3 &
-(12.0 * (ZHLC_RCMAX(:))*(ZRCRAUTC(:))**2) &
@@ -261,7 +318,8 @@ ELSEIF(HSUBG_AUCV=='PDF ') THEN
END WHERE
END IF
! Compare r_cM to r_cR to know if cloud water content is high enough to split in two parts or not
- WHERE (PRCT(:).GT.0. .AND. PCF(:).GT.0. .AND. ZHLC_RCMAX(:).GT.ZRCRAUTC(:))
+ GWORK2(:) = PRCT(:)>0. .AND. PCF(:)>0. .AND. ZHLC_RCMAX(:)<=ZRCRAUTC(:)
+ WHERE(GWORK(:))
! Calculate final values for LCF and HCF:
PHLC_LCF(:)=PCF(:) &
*(ZHLC_HRCLOCAL- &
@@ -272,7 +330,7 @@ ELSEIF(HSUBG_AUCV=='PDF ') THEN
! Calculate final values for LRC and HRC:
PHLC_LRC(:)=ZHLC_LRCLOCAL*PHLC_LCF(:)
PHLC_HRC(:)=MAX(0., PRCT(:)-PHLC_LRC(:))
- ELSEWHERE (PRCT(:).GT.0. .AND. PCF(:).GT.0. .AND. ZHLC_RCMAX(:).LE.ZRCRAUTC(:))
+ ELSEWHERE (GWORK2(:))
! Put all available cloud water and his fraction in the low part
PHLC_LCF(:)=PCF(:)
PHLC_HCF(:)=0.
@@ -285,6 +343,7 @@ ELSEIF(HSUBG_AUCV=='PDF ') THEN
PHLC_HRC(:)=0.
END WHERE
PRF(:)=PHLC_HCF(:) !Precipitation fraction
+!$acc end kernels
ELSE
CALL PRINT_MSG(NVERB_FATAL,'GEN','ICE4_COMPUTE_PDF','wrong HSUBG_PR_PDF case')
ENDIF
diff --git a/src/MNH/ice4_fast_rg.f90 b/src/MNH/ice4_fast_rg.f90
index c37a8cd88a163825dec691b975892fd636d6b894..25af952d6fd5f6d524afd113e92653b2651af20b 100644
--- a/src/MNH/ice4_fast_rg.f90
+++ b/src/MNH/ice4_fast_rg.f90
@@ -96,10 +96,13 @@ SUBROUTINE ICE4_FAST_RG(KSIZE, LDSOFT, LDCOMPUTE, KRR, &
!* 0. DECLARATIONS
! ------------
!
-USE MODD_CST
-USE MODD_RAIN_ICE_PARAM
-USE MODD_RAIN_ICE_DESCR
-USE MODD_PARAM_ICE, ONLY : LEVLIMIT, LNULLWETG, LWETGPOST, LCRFLIMIT
+USE MODD_CST, ONLY: XALPI,XALPW,XBETAI,XBETAW,XGAMW,XCI,XCL,XCPV,XESTT,XGAMI,XLMTT,XLVTT,XMD,XMV,XRV,XTT
+USE MODD_PARAM_ICE, ONLY: LCRFLIMIT,LEVLIMIT,LNULLWETG,LWETGPOST
+USE MODD_RAIN_ICE_DESCR, ONLY: XBS,XCEXVT,XCXG,XCXS,XDG,XRTMIN
+USE MODD_RAIN_ICE_PARAM, ONLY: NDRYLBDAG,NDRYLBDAR,NDRYLBDAS,X0DEPG,X1DEPG,XCOLEXIG,XCOLEXSG,XCOLIG,XCOLSG,XDRYINTP1G, &
+ XDRYINTP1R,XDRYINTP1S,XDRYINTP2G,XDRYINTP2R,XDRYINTP2S,XEX0DEPG,XEX1DEPG,XEXICFRR, &
+ XEXRCFRI,XFCDRYG,XFIDRYG,XFRDRYG,XFSDRYG,XICFRR,XKER_RDRYG,XKER_SDRYG,XLBRDRYG1, &
+ XLBRDRYG2,XLBRDRYG3,XLBSDRYG1,XLBSDRYG2,XLBSDRYG3,XRCFRI
!
#ifdef MNH_BITREP
USE MODI_BITREP
@@ -161,18 +164,35 @@ REAL, DIMENSION(KSIZE), INTENT(INOUT) :: PA_RH
REAL, DIMENSION(KSIZE), INTENT(INOUT) :: PB_RG
REAL, DIMENSION(KSIZE), INTENT(INOUT) :: PB_RH
!
+!$acc declare present(KSIZE,LDSOFT,LDCOMPUTE,KRR,PRHODREF,PLVFACT,PLSFACT,PPRES,PDV,PKA,PCJ,PCIT,PLBDAR,PLBDAS,PLBDAG, &
+!$acc& PT,PRVT,PRCT,PRRT,PRIT,PRST,PRGT,PRGSI,PRGSI_MR,LDWETG,PRICFRRG,PRRCFRIG,PRICFRR,PRCWETG,PRIWETG, &
+!$acc& PRRWETG,PRSWETG,PRCDRYG,PRIDRYG,PRRDRYG,PRSDRYG,PRWETGH,PRWETGH_MR,PRGMLTR,PRG_TEND, &
+!$acc& PA_TH,PA_RC,PA_RR,PA_RI,PA_RS,PA_RG,PA_RH,PB_RG,PB_RH)
+!
!* 0.2 declaration of local variables
!
INTEGER, PARAMETER :: IRCDRYG=1, IRIDRYG=2, IRIWETG=3, IRSDRYG=4, IRSWETG=5, IRRDRYG=6
!
-LOGICAL, DIMENSION(SIZE(PRHODREF)) :: GDRY, LLDRYG, GMASK
-INTEGER :: IGDRY
-REAL, DIMENSION(SIZE(PRHODREF)) :: ZVEC1, ZVEC2, ZVEC3
+INTEGER :: IGDRY
+INTEGER :: IDX, JJ
+INTEGER :: IRR !Workaround of PGI bug with OpenACC (at least up to 18.10 version)
+LOGICAL :: GCRFLIMIT,GDSOFT,GEVLIMIT,GNULLWETG,GWETGPOST !Workaround of PGI bug with OpenACC (at least up to 18.10 version)
+LOGICAL, DIMENSION(SIZE(PRHODREF)) :: GDRY, GLDRYG, GMASK
INTEGER, DIMENSION(SIZE(PRHODREF)) :: IVEC1, IVEC2
-REAL, DIMENSION(SIZE(PRHODREF)) :: ZZW, &
- ZRDRYG_INIT, & !Initial dry growth rate of the graupeln
- ZRWETG_INIT !Initial wet growth rate of the graupeln
-INTEGER :: JJ
+REAL, DIMENSION(SIZE(PRHODREF)) :: ZVEC1, ZVEC2, ZVEC3
+REAL, DIMENSION(SIZE(PRHODREF)) :: ZZW, &
+ ZRDRYG_INIT, & !Initial dry growth rate of the graupeln
+ ZRWETG_INIT !Initial wet growth rate of the graupeln
+!
+!$acc declare create(IGDRY,IDX,JJ,GDRY,GLDRYG,GMASK,IVEC1,IVEC2,ZVEC1,ZVEC2,ZVEC3,ZZW,ZRDRYG_INIT,ZRWETG_INIT)
+!
+!$acc declare copyin(XALPI,XALPW,XBETAI,XBETAW,XGAMW,XCI,XCL,XCPV,XESTT,XGAMI,XLMTT,XLVTT,XMD,XMV,XRV,XTT, &
+!$acc& LCRFLIMIT,LEVLIMIT,LNULLWETG,LWETGPOST, &
+!$acc& XBS,XCEXVT,XCXG,XCXS,XDG,XRTMIN, &
+!$acc& NDRYLBDAG,NDRYLBDAR,NDRYLBDAS,X0DEPG,X1DEPG,XCOLEXIG,XCOLEXSG,XCOLIG,XCOLSG,XDRYINTP1G, &
+!$acc& XDRYINTP1R,XDRYINTP1S,XDRYINTP2G,XDRYINTP2R,XDRYINTP2S,XEX0DEPG,XEX1DEPG,XEXICFRR, &
+!$acc& XEXRCFRI,XFCDRYG,XFIDRYG,XFRDRYG,XFSDRYG,XICFRR,XKER_RDRYG,XKER_SDRYG,XLBRDRYG1, &
+!$acc& XLBRDRYG2,XLBRDRYG3,XLBSDRYG1,XLBSDRYG2,XLBSDRYG3,XRCFRI)
!
!-------------------------------------------------------------------------------
!
@@ -216,10 +236,21 @@ IF (MPPDB_INITIALIZED) THEN
CALL MPPDB_CHECK1D(PB_RH,"ICE4_FAST_RG beg:PB_RH",PRECISION)
END IF
!
+!Workaround of PGI bug with OpenACC (at least up to 18.10 version)
+GCRFLIMIT = LCRFLIMIT
+GEVLIMIT = LEVLIMIT
+GNULLWETG = LNULLWETG
+GWETGPOST = LWETGPOST
+!$acc kernels
+GDSOFT = LDSOFT
+IRR = KRR
+!$acc end kernels
+!
!* 6.1 rain contact freezing
!
+!$acc kernels
GMASK(:)=PRIT(:)>XRTMIN(4) .AND. PRRT(:)>XRTMIN(3) .AND. LDCOMPUTE(:)
-IF(LDSOFT) THEN
+IF(GDSOFT) THEN
WHERE(.NOT. GMASK(:))
PRICFRRG(:)=0.
PRRCFRIG(:)=0.
@@ -247,7 +278,7 @@ ELSE
#endif
END WHERE
ZZW(:)=1.
- IF(LCRFLIMIT) THEN
+ IF(GCRFLIMIT) THEN
WHERE(GMASK(:))
!Comparison between heat to be released (to freeze rain) and heat sink (rain and ice temperature change)
!ZZW is the proportion of process that can take place
@@ -269,7 +300,7 @@ PA_TH(:) = PA_TH(:) + (PRRCFRIG(:) - PRICFRR(:))*(PLSFACT(:)-PLVFACT(:))
!
! Wet and dry collection of rc and ri on graupel
GMASK(:)=PRGT(:)>XRTMIN(6) .AND. PRCT(:)>XRTMIN(2) .AND. LDCOMPUTE(:)
-IF(LDSOFT) THEN
+IF(GDSOFT) THEN
WHERE(.NOT. GMASK(:))
PRG_TEND(:, IRCDRYG)=0.
END WHERE
@@ -285,7 +316,7 @@ ELSE
END WHERE
ENDIF
GMASK(:)=PRGT(:)>XRTMIN(6) .AND. PRIT(:)>XRTMIN(4) .AND. LDCOMPUTE(:)
-IF(LDSOFT) THEN
+IF(GDSOFT) THEN
WHERE(.NOT. GMASK(:))
PRG_TEND(:, IRIDRYG)=0.
PRG_TEND(:, IRIWETG)=0.
@@ -308,7 +339,7 @@ ENDIF
! Wet and dry collection of rs on graupel (6.2.1)
GDRY(:)=PRST(:)>XRTMIN(5) .AND. PRGT(:)>XRTMIN(6) .AND. LDCOMPUTE(:)
-IF(LDSOFT) THEN
+IF(GDSOFT) THEN
WHERE(.NOT. GDRY(:))
PRG_TEND(:, IRSDRYG)=0.
PRG_TEND(:, IRSWETG)=0.
@@ -321,8 +352,19 @@ ELSE
!
!* 6.2.3 select the (PLBDAG,PLBDAS) couplet
!
+#ifndef _OPENACC
ZVEC1(1:IGDRY)=PACK(PLBDAG(:), MASK=GDRY(:))
ZVEC2(1:IGDRY)=PACK(PLBDAS(:), MASK=GDRY(:))
+#else
+ IDX = 0
+ DO JJ = 1,KSIZE
+ IF (GDRY(JJ)) THEN
+ IDX = IDX+1
+ ZVEC1(IDX) = PLBDAG(JJ)
+ ZVEC2(IDX) = PLBDAS(JJ)
+ END IF
+ END DO
+#endif
!
!* 6.2.4 find the next lower indice for the PLBDAG and for the PLBDAS
! in the geometrical set of (Lbda_g,Lbda_s) couplet use to
@@ -357,7 +399,19 @@ ELSE
- XKER_SDRYG(IVEC1(JJ) ,IVEC2(JJ) )*(ZVEC2(JJ) - 1.0) ) &
*(ZVEC1(JJ) - 1.0)
END DO
+#ifndef _OPENACC
ZZW(:)=UNPACK(VECTOR=ZVEC3(1:IGDRY), MASK=GDRY(:), FIELD=0.0)
+#else
+ IDX = 0
+ DO JJ=1,KSIZE
+ IF (GDRY(JJ)) THEN
+ IDX = IDX + 1
+ ZZW(JJ) = ZVEC3(IDX)
+ ELSE
+ ZZW(JJ) = 0.
+ END IF
+ END DO
+#endif
!
WHERE(GDRY(:))
#ifndef MNH_BITREP
@@ -386,7 +440,7 @@ ENDIF
!* 6.2.6 accretion of raindrops on the graupeln
!
GDRY(:)=PRRT(:)>XRTMIN(3) .AND. PRGT(:)>XRTMIN(6) .AND. LDCOMPUTE(:)
-IF(LDSOFT) THEN
+IF(GDSOFT) THEN
WHERE(.NOT. GDRY(:))
PRG_TEND(:, IRRDRYG)=0.
END WHERE
@@ -398,8 +452,19 @@ ELSE
!
!* 6.2.8 select the (PLBDAG,PLBDAR) couplet
!
+#ifndef _OPENACC
ZVEC1(1:IGDRY)=PACK(PLBDAG(:), MASK=GDRY(:))
ZVEC2(1:IGDRY)=PACK(PLBDAR(:), MASK=GDRY(:))
+#else
+ IDX = 0
+ DO JJ = 1,KSIZE
+ IF (GDRY(JJ)) THEN
+ IDX = IDX+1
+ ZVEC1(IDX) = PLBDAG(JJ)
+ ZVEC2(IDX) = PLBDAR(JJ)
+ END IF
+ END DO
+#endif
!
!* 6.2.9 find the next lower indice for the PLBDAG and for the PLBDAR
! in the geometrical set of (Lbda_g,Lbda_r) couplet use to
@@ -434,7 +499,19 @@ ELSE
- XKER_RDRYG(IVEC1(JJ) ,IVEC2(JJ) )*(ZVEC2(JJ) - 1.0) ) &
*(ZVEC1(JJ) - 1.0)
END DO
+#ifndef _OPENACC
ZZW(:)=UNPACK(VECTOR=ZVEC3(1:IGDRY), MASK=GDRY, FIELD=0.)
+#else
+ IDX = 0
+ DO JJ=1,KSIZE
+ IF (GDRY(JJ)) THEN
+ IDX = IDX + 1
+ ZZW(JJ) = ZVEC3(IDX)
+ ELSE
+ ZZW(JJ) = 0.
+ END IF
+ END DO
+#endif
!
WHERE(GDRY(:))
#ifndef MNH_BITREP
@@ -464,7 +541,7 @@ GMASK(:)=PRGT(:)>XRTMIN(6) .AND. LDCOMPUTE(:)
WHERE(GMASK(:))
ZRWETG_INIT(:)=PRVT(:)*PPRES(:)/((XMV/XMD)+PRVT(:)) ! Vapor pressure
END WHERE
-IF(LEVLIMIT) THEN
+IF(GEVLIMIT) THEN
WHERE(GMASK(:))
#ifndef MNH_BITREP
ZRWETG_INIT(:)=MIN(ZRWETG_INIT(:), EXP(XALPI-XBETAI/PT(:)-XGAMI*ALOG(PT(:)))) ! min(ev, es_i(T))
@@ -496,14 +573,14 @@ END WHERE
LDWETG(:)=GMASK(:) .AND. &
&MAX(0., ZRDRYG_INIT(:)-PRG_TEND(:, IRIDRYG)-PRG_TEND(:, IRSDRYG))>= &
&MAX(0., ZRWETG_INIT(:)-PRG_TEND(:, IRIWETG)-PRG_TEND(:, IRSWETG))
-IF(LNULLWETG) THEN
+IF(GNULLWETG) THEN
LDWETG(:)=LDWETG(:) .AND. ZRDRYG_INIT(:)>0.
ELSE
LDWETG(:)=LDWETG(:) .AND. ZRWETG_INIT(:)>0.
ENDIF
-IF(.NOT. LWETGPOST) LDWETG(:)=LDWETG(:) .AND. PT(:)<XTT
+IF(.NOT. GWETGPOST) LDWETG(:)=LDWETG(:) .AND. PT(:)<XTT
-LLDRYG(:)=GMASK(:) .AND. PT(:)<XTT .AND. ZRDRYG_INIT(:)>0. .AND. &
+GLDRYG(:)=GMASK(:) .AND. PT(:)<XTT .AND. ZRDRYG_INIT(:)>0. .AND. &
&MAX(0., ZRDRYG_INIT(:)-PRG_TEND(:, IRIDRYG)-PRG_TEND(:, IRSDRYG))<&
&MAX(0., ZRWETG_INIT(:)-PRG_TEND(:, IRIWETG)-PRG_TEND(:, IRSWETG))
@@ -512,7 +589,7 @@ LLDRYG(:)=GMASK(:) .AND. PT(:)<XTT .AND. ZRDRYG_INIT(:)>0. .AND. &
! as a tendency (PRWETGH)
PRWETGH(:)=0.
PRWETGH_MR(:)=0.
-IF(KRR==7) THEN
+IF(IRR==7) THEN
WHERE(LDWETG(:))
!assume a linear percent of conversion of produced graupel into hail
PRWETGH(:)=(MAX(0., PRGSI(:)+PRICFRRG(:)+PRRCFRIG(:))+ZRWETG_INIT(:))*ZRDRYG_INIT(:)/(ZRWETG_INIT(:)+ZRDRYG_INIT(:))
@@ -535,7 +612,7 @@ PRCDRYG(:)=0.
PRIDRYG(:)=0.
PRRDRYG(:)=0.
PRSDRYG(:)=0.
-WHERE(LLDRYG(:))
+WHERE(GLDRYG(:))
PRCDRYG(:)=PRG_TEND(:, IRCDRYG)
PRRDRYG(:)=PRG_TEND(:, IRRDRYG)
PRIDRYG(:)=PRG_TEND(:, IRIDRYG)
@@ -562,7 +639,7 @@ PA_TH(:) = PA_TH(:) + (PRCDRYG(:)+PRRDRYG(:))*(PLSFACT(:)-PLVFACT(:))
!* 6.5 Melting of the graupeln
!
GMASK(:)=PRGT(:)>XRTMIN(6) .AND. PT(:)>XTT .AND. LDCOMPUTE(:)
-IF(LDSOFT) THEN
+IF(GDSOFT) THEN
WHERE(.NOT. GMASK(:))
PRGMLTR(:) = 0.
END WHERE
@@ -571,7 +648,7 @@ ELSE
WHERE(GMASK(:))
PRGMLTR(:) = PRVT(:)*PPRES(:)/((XMV/XMD)+PRVT(:)) ! Vapor pressure
END WHERE
- IF(LEVLIMIT) THEN
+ IF(GEVLIMIT) THEN
WHERE(GMASK(:))
#ifndef MNH_BITREP
PRGMLTR(:)=MIN(PRGMLTR(:), EXP(XALPW-XBETAW/PT(:)-XGAMW*ALOG(PT(:)))) ! min(ev, es_w(T))
@@ -605,6 +682,7 @@ ENDIF
PA_RR(:) = PA_RR(:) + PRGMLTR(:)
PA_RG(:) = PA_RG(:) - PRGMLTR(:)
PA_TH(:) = PA_TH(:) - PRGMLTR(:)*(PLSFACT(:)-PLVFACT(:))
+!$acc end kernels
!
IF (MPPDB_INITIALIZED) THEN
!Check all INOUT arrays
diff --git a/src/MNH/ice4_fast_rh.f90 b/src/MNH/ice4_fast_rh.f90
index 09d3860396d10a77b12c1a33001566c21a36a958..255d719d66635df53cc8b2d56034018487138a17 100644
--- a/src/MNH/ice4_fast_rh.f90
+++ b/src/MNH/ice4_fast_rh.f90
@@ -86,10 +86,13 @@ SUBROUTINE ICE4_FAST_RH(KSIZE, LDSOFT, LDCOMPUTE, LDWETG, &
!* 0. DECLARATIONS
! ------------
!
-USE MODD_CST
-USE MODD_RAIN_ICE_PARAM
-USE MODD_RAIN_ICE_DESCR
-USE MODD_PARAM_ICE, ONLY : LEVLIMIT, LNULLWETH, LWETHPOST, LCONVHG
+USE MODD_CST, ONLY: XALPI,XALPW,XBETAI,XBETAW,XGAMW,XCI,XCL,XCPV,XESTT,XGAMI,XLMTT,XLVTT,XMD,XMV,XRV,XTT
+USE MODD_PARAM_ICE, ONLY: LCONVHG,LEVLIMIT,LNULLWETH,LWETHPOST
+USE MODD_RAIN_ICE_DESCR, ONLY: XBG,XBS,XCEXVT,XCXG,XCXH,XCXS,XDH,XRTMIN
+USE MODD_RAIN_ICE_PARAM, ONLY: NWETLBDAG,NWETLBDAH,NWETLBDAR,NWETLBDAS,X0DEPH,X1DEPH,XCOLEXGH,XCOLEXIH,XCOLGH,XCOLIH,XCOLEXSH, &
+ XCOLSH,XEX0DEPH,XEX1DEPH,XFGWETH,XFRWETH,XFSWETH,XFWETH,XKER_GWETH,XKER_RWETH,XKER_SWETH, &
+ XLBGWETH1,XLBGWETH2,XLBGWETH3,XLBRWETH1,XLBRWETH2,XLBRWETH3,XLBSWETH1,XLBSWETH2,XLBSWETH3, &
+ XWETINTP1G,XWETINTP1H,XWETINTP1R,XWETINTP1S,XWETINTP2G,XWETINTP2H,XWETINTP2R,XWETINTP2S
!
#ifdef MNH_BITREP
USE MODI_BITREP
@@ -145,21 +148,41 @@ REAL, DIMENSION(KSIZE), INTENT(INOUT) :: PA_RS
REAL, DIMENSION(KSIZE), INTENT(INOUT) :: PA_RG
REAL, DIMENSION(KSIZE), INTENT(INOUT) :: PA_RH
!
+!$acc declare present(KSIZE,LDSOFT,LDCOMPUTE,LDWETG,PRHODREF,PLVFACT,PLSFACT,PPRES,PDV,PKA,PCJ, &
+!$acc& PLBDAS,PLBDAG,PLBDAR,PLBDAH,PT,PRVT,PRCT,PRRT,PRIT,PRST,PRGT,PRHT, &
+!$acc& PRCWETH,PRIWETH,PRSWETH,PRGWETH,PRRWETH, &
+!$acc& PRCDRYH,PRIDRYH,PRSDRYH,PRRDRYH,PRGDRYH,PRDRYHG,PRHMLTR, &
+!$acc& PRH_TEND,PA_TH,PA_RC,PA_RR,PA_RI,PA_RS,PA_RG,PA_RH)
+!
!* 0.2 declaration of local variables
!
INTEGER, PARAMETER :: IRCWETH=1, IRRWETH=2, IRIDRYH=3, IRIWETH=4, IRSDRYH=5, IRSWETH=6, IRGDRYH=7, IRGWETH=8
!
-LOGICAL, DIMENSION(SIZE(PRHODREF)) :: GHAIL, GWET, GMASK, LLWETH, LLDRYH
-INTEGER :: IHAIL, IGWET
-REAL, DIMENSION(SIZE(PRHODREF)) :: ZVEC1, ZVEC2, ZVEC3
+INTEGER :: IHAIL, IGWET
+INTEGER :: IDX, JJ
INTEGER, DIMENSION(SIZE(PRHODREF)) :: IVEC1, IVEC2
-REAL, DIMENSION(SIZE(PRHODREF)) :: ZZW, &
- ZRDRYH_INIT, ZRWETH_INIT, &
- ZRDRYHG
-INTEGER :: JJ
+LOGICAL :: GCONVHG,GDSOFT,GEVLIMIT,GNULLWETH,GWETHPOST !Workaround of PGI bug with OpenACC (at least up to 18.10 version)
+LOGICAL, DIMENSION(SIZE(PRHODREF)) :: GHAIL, GWET, GMASK, LLWETH, LLDRYH
+REAL, DIMENSION(SIZE(PRHODREF)) :: ZVEC1, ZVEC2, ZVEC3
+REAL, DIMENSION(SIZE(PRHODREF)) :: ZZW, &
+ ZRDRYH_INIT, ZRWETH_INIT, &
+ ZRDRYHG
+!
+!$acc declare create(IHAIL,IGWET,IDX,JJ,IVEC1,IVEC2,GHAIL,GWET,GMASK,LLWETH,LLDRYH, &
+!$acc& ZVEC1,ZVEC2,ZVEC3,ZZW,ZRDRYH_INIT,ZRWETH_INIT,ZRDRYHG)
+!
+!$acc declare copyin(XALPI,XALPW,XBETAI,XBETAW,XGAMW,XCI,XCL,XCPV,XESTT,XGAMI,XLMTT,XLVTT,XMD,XMV,XRV,XTT, &
+!$acc& LCONVHG,LEVLIMIT,LNULLWETH,LWETHPOST, &
+!$acc& XBG,XBS,XCEXVT,XCXG,XCXH,XCXS,XDH,XRTMIN, &
+!$acc& NWETLBDAG,NWETLBDAH,NWETLBDAR,NWETLBDAS,X0DEPH,X1DEPH,XCOLEXGH,XCOLEXIH,XCOLGH,XCOLIH,XCOLEXSH, &
+!$acc& XCOLSH,XEX0DEPH,XEX1DEPH,XFGWETH,XFRWETH,XFSWETH,XFWETH,XKER_GWETH,XKER_RWETH,XKER_SWETH, &
+!$acc& XLBGWETH1,XLBGWETH2,XLBGWETH3,XLBRWETH1,XLBRWETH2,XLBRWETH3,XLBSWETH1,XLBSWETH2,XLBSWETH3, &
+!$acc& XWETINTP1G,XWETINTP1H,XWETINTP1R,XWETINTP1S,XWETINTP2G,XWETINTP2H,XWETINTP2R,XWETINTP2S)
!
!-------------------------------------------------------------------------------
!
+CALL PRINT_MSG(NVERB_ERROR,'GEN','ICE4_FAST_RH','OPENACC: implemented but never tested')
+!
IF (MPPDB_INITIALIZED) THEN
!Check all IN arrays
CALL MPPDB_CHECK1D_LOG(LDCOMPUTE,"ICE4_FAST_RH beg:LDCOMPUTE")
@@ -195,10 +218,20 @@ IF (MPPDB_INITIALIZED) THEN
CALL MPPDB_CHECK1D(PA_RH,"ICE4_FAST_RH beg:PA_RH",PRECISION)
END IF
!
+!Workaround of PGI bug with OpenACC (at least up to 18.10 version)
+GCONVHG = LCONVHG
+GEVLIMIT = LEVLIMIT
+GNULLWETH = LNULLWETH
+GWETHPOST = LWETHPOST
+!$acc kernels
+GDSOFT = LDSOFT
+!$acc end kernels
+!
!* 7.2 compute the Wet and Dry growth of hail
!
+!$acc kernels
GMASK(:)=PRHT(:)>XRTMIN(7) .AND. PRCT(:)>XRTMIN(2) .AND. LDCOMPUTE(:)
-IF(LDSOFT) THEN
+IF(GDSOFT) THEN
WHERE(.NOT. GMASK(:))
PRH_TEND(:, IRCWETH)=0.
END WHERE
@@ -214,7 +247,7 @@ ELSE
END WHERE
ENDIF
GMASK(:)=PRHT(:)>XRTMIN(7) .AND. PRIT(:)>XRTMIN(4) .AND. LDCOMPUTE(:)
-IF(LDSOFT) THEN
+IF(GDSOFT) THEN
WHERE(.NOT. GMASK(:))
PRH_TEND(:, IRIWETH)=0.
PRH_TEND(:, IRIDRYH)=0.
@@ -234,12 +267,11 @@ ELSE
#endif
END WHERE
ENDIF
-
!
!* 7.2.1 accretion of aggregates on the hailstones
!
GWET(:) = PRHT(:)>XRTMIN(7) .AND. PRST(:)>XRTMIN(5) .AND. LDCOMPUTE(:)
-IF(LDSOFT) THEN
+IF(GDSOFT) THEN
WHERE(.NOT. GWET(:))
PRH_TEND(:, IRSWETH)=0.
PRH_TEND(:, IRSDRYH)=0.
@@ -252,8 +284,19 @@ ELSE
!
!* 7.2.3 select the (PLBDAH,PLBDAS) couplet
!
+#ifndef _OPENACC
ZVEC1(1:IGWET) = PACK( PLBDAH(:),MASK=GWET(:) )
ZVEC2(1:IGWET) = PACK( PLBDAS(:),MASK=GWET(:) )
+#else
+ IDX = 0
+ DO JJ = 1,KSIZE
+ IF (GWET(JJ)) THEN
+ IDX = IDX+1
+ ZVEC1(IDX) = PLBDAH(JJ)
+ ZVEC2(IDX) = PLBDAS(JJ)
+ END IF
+ END DO
+#endif
!
!* 7.2.4 find the next lower indice for the PLBDAG and for the PLBDAS
! in the geometrical set of (Lbda_h,Lbda_s) couplet use to
@@ -288,7 +331,19 @@ ELSE
- XKER_SWETH(IVEC1(JJ) ,IVEC2(JJ) )*(ZVEC2(JJ) - 1.0) ) &
* (ZVEC1(JJ) - 1.0)
END DO
+#ifndef _OPENACC
ZZW(:) = UNPACK( VECTOR=ZVEC3(1:IGWET),MASK=GWET,FIELD=0.0 )
+#else
+ IDX = 0
+ DO JJ=1,KSIZE
+ IF (GWET(JJ)) THEN
+ IDX = IDX + 1
+ ZZW(JJ) = ZVEC3(IDX)
+ ELSE
+ ZZW(JJ) = 0.
+ END IF
+ END DO
+#endif
!
WHERE(GWET(:))
#ifndef MNH_BITREP
@@ -315,7 +370,7 @@ ENDIF
!* 7.2.6 accretion of graupeln on the hailstones
!
GWET(:) = PRHT(:)>XRTMIN(7) .AND. PRGT(:)>XRTMIN(6) .AND. LDCOMPUTE(:)
-IF(LDSOFT) THEN
+IF(GDSOFT) THEN
WHERE(.NOT. GWET(:))
PRH_TEND(:, IRGWETH)=0.
PRH_TEND(:, IRGDRYH)=0.
@@ -328,8 +383,19 @@ ELSE
!
!* 7.2.8 select the (PLBDAH,PLBDAG) couplet
!
+#ifndef _OPENACC
ZVEC1(1:IGWET) = PACK( PLBDAH(:),MASK=GWET(:) )
ZVEC2(1:IGWET) = PACK( PLBDAG(:),MASK=GWET(:) )
+#else
+ IDX = 0
+ DO JJ = 1,KSIZE
+ IF (GWET(JJ)) THEN
+ IDX = IDX+1
+ ZVEC1(IDX) = PLBDAH(JJ)
+ ZVEC2(IDX) = PLBDAG(JJ)
+ END IF
+ END DO
+#endif
!
!* 7.2.9 find the next lower indice for the PLBDAH and for the PLBDAG
! in the geometrical set of (Lbda_h,Lbda_g) couplet use to
@@ -364,7 +430,19 @@ ELSE
- XKER_GWETH(IVEC1(JJ) ,IVEC2(JJ) )*(ZVEC2(JJ) - 1.0) ) &
* (ZVEC1(JJ) - 1.0)
END DO
+#ifndef _OPENACC
ZZW(:) = UNPACK( VECTOR=ZVEC3(1:IGWET),MASK=GWET,FIELD=0.0 )
+#else
+ IDX = 0
+ DO JJ=1,KSIZE
+ IF (GWET(JJ)) THEN
+ IDX = IDX + 1
+ ZZW(JJ) = ZVEC3(IDX)
+ ELSE
+ ZZW(JJ) = 0.
+ END IF
+ END DO
+#endif
!
WHERE(GWET(:))
#ifndef MNH_BITREP
@@ -398,7 +476,7 @@ ENDIF
!* 7.2.11 accretion of raindrops on the hailstones
!
GWET(:) = PRHT(:)>XRTMIN(7) .AND. PRRT(:)>XRTMIN(3) .AND. LDCOMPUTE(:)
-IF(LDSOFT) THEN
+IF(GDSOFT) THEN
WHERE(.NOT. GWET(:))
PRH_TEND(:, IRRWETH)=0.
END WHERE
@@ -409,8 +487,19 @@ ELSE
!
!* 7.2.12 select the (PLBDAH,PLBDAR) couplet
!
+#ifndef _OPENACC
ZVEC1(1:IGWET)=PACK(PLBDAH(:), MASK=GWET(:))
ZVEC2(1:IGWET)=PACK(PLBDAR(:), MASK=GWET(:))
+#else
+ IDX = 0
+ DO JJ = 1,KSIZE
+ IF (GWET(JJ)) THEN
+ IDX = IDX+1
+ ZVEC1(IDX) = PLBDAH(JJ)
+ ZVEC2(IDX) = PLBDAR(JJ)
+ END IF
+ END DO
+#endif
!
!* 7.2.13 find the next lower indice for the PLBDAH and for the PLBDAR
! in the geometrical set of (Lbda_h,Lbda_r) couplet use to
@@ -445,7 +534,19 @@ ELSE
- XKER_RWETH(IVEC1(JJ) ,IVEC2(JJ) )*(ZVEC2(JJ) - 1.0) ) &
*(ZVEC1(JJ) - 1.0)
END DO
+#ifndef _OPENACC
ZZW(:)=UNPACK(VECTOR=ZVEC3(1:IGWET), MASK=GWET, FIELD=0.)
+#else
+ IDX = 0
+ DO JJ=1,KSIZE
+ IF (GWET(JJ)) THEN
+ IDX = IDX + 1
+ ZZW(JJ) = ZVEC3(IDX)
+ ELSE
+ ZZW(JJ) = 0.
+ END IF
+ END DO
+#endif
!
WHERE(GWET(:))
#ifndef MNH_BITREP
@@ -476,7 +577,7 @@ ZRWETH_INIT(:)=0.
WHERE(GHAIL(:))
ZRWETH_INIT(:) = PRVT(:)*PPRES(:)/((XMV/XMD)+PRVT(:)) ! Vapor pressure
END WHERE
-IF(LEVLIMIT) THEN
+IF(GEVLIMIT) THEN
WHERE(GHAIL(:))
#ifndef MNH_BITREP
ZRWETH_INIT(:) = MIN(ZRWETH_INIT(:), EXP(XALPI-XBETAI/PT(:)-XGAMI*ALOG(PT(:)))) ! min(ev, es_i(T))
@@ -510,12 +611,12 @@ END WHERE
!Wet case
LLWETH(:)=GHAIL(:) .AND. MAX(0., ZRDRYH_INIT(:)-PRH_TEND(:, IRIDRYH)-PRH_TEND(:, IRSDRYH)-PRH_TEND(:, IRGDRYH))>= &
& MAX(0., ZRWETH_INIT(:)-PRH_TEND(:, IRIWETH)-PRH_TEND(:, IRSWETH)-PRH_TEND(:, IRGWETH))
-IF(LNULLWETH) THEN
+IF(GNULLWETH) THEN
LLWETH(:)=LLWETH(:) .AND. ZRDRYH_INIT(:)>0.
ELSE
LLWETH(:)=LLWETH(:) .AND. ZRWETH_INIT(:)>0.
ENDIF
-IF(.NOT. LWETHPOST) LLWETH(:)=LLWETH(:) .AND. PT(:)<XTT
+IF(.NOT. GWETHPOST) LLWETH(:)=LLWETH(:) .AND. PT(:)<XTT
LLDRYH(:)=GHAIL(:) .AND. PT(:)<XTT .AND. ZRDRYH_INIT(:)>0. .AND. &
& MAX(0., ZRDRYH_INIT(:)-PRH_TEND(:, IRIDRYH)-PRH_TEND(:, IRSDRYH))< &
& MAX(0., ZRWETH_INIT(:)-PRH_TEND(:, IRIWETH)-PRH_TEND(:, IRSWETH))
@@ -541,7 +642,7 @@ PRRDRYH(:) = 0.
PRGDRYH(:) = 0.
PRDRYHG(:) = 0.
ZRDRYHG(:)=0.
-IF(LCONVHG)THEN
+IF(GCONVHG)THEN
WHERE(LLDRYH(:))
ZRDRYHG(:)=ZRDRYH_INIT(:)*ZRWETH_INIT(:)/(ZRDRYH_INIT(:)+ZRWETH_INIT(:))
END WHERE
@@ -572,7 +673,7 @@ PA_TH(:) = PA_TH(:) + (PRCDRYH(:)+PRRDRYH(:))*(PLSFACT(:)-PLVFACT(:))
!* 7.5 Melting of the hailstones
!
GMASK(:)=PRHT(:)>XRTMIN(7) .AND. PT(:)>XTT .AND. LDCOMPUTE(:)
-IF(LDSOFT) THEN
+IF(GDSOFT) THEN
WHERE(.NOT. GMASK(:))
PRHMLTR(:) = 0.
END WHERE
@@ -581,7 +682,7 @@ ELSE
WHERE(GMASK(:))
PRHMLTR(:) = PRVT(:)*PPRES(:)/((XMV/XMD)+PRVT(:)) ! Vapor pressure
END WHERE
- IF(LEVLIMIT) THEN
+ IF(GEVLIMIT) THEN
WHERE(GMASK(:))
#ifndef MNH_BITREP
PRHMLTR(:)=MIN(PRHMLTR(:), EXP(XALPW-XBETAW/PT(:)-XGAMW*ALOG(PT(:)))) ! min(ev, es_w(T))
@@ -615,6 +716,7 @@ END IF
PA_RR(:) = PA_RR(:) + PRHMLTR(:)
PA_RH(:) = PA_RH(:) - PRHMLTR(:)
PA_TH(:) = PA_TH(:) - PRHMLTR(:)*(PLSFACT(:)-PLVFACT(:))
+!$acc end kernels
!
IF (MPPDB_INITIALIZED) THEN
!Check all INOUT arrays
diff --git a/src/MNH/ice4_fast_ri.f90 b/src/MNH/ice4_fast_ri.f90
index 7cf9af5412a1a2826eb245062f7eec01c0206c18..d5c969a0bf86564c51812c171a11583b6c7c3724 100644
--- a/src/MNH/ice4_fast_ri.f90
+++ b/src/MNH/ice4_fast_ri.f90
@@ -66,9 +66,8 @@ SUBROUTINE ICE4_FAST_RI(KSIZE, LDSOFT, LDCOMPUTE, &
!* 0. DECLARATIONS
! ------------
!
-USE MODD_CST
-USE MODD_RAIN_ICE_PARAM
-USE MODD_RAIN_ICE_DESCR
+USE MODD_RAIN_ICE_DESCR, ONLY: XDI,XLBEXI,XLBI,XRTMIN
+USE MODD_RAIN_ICE_PARAM, ONLY: X0DEPI,X2DEPI
!
#ifdef MNH_BITREP
USE MODI_BITREP
@@ -97,10 +96,15 @@ REAL, DIMENSION(KSIZE), INTENT(INOUT) :: PA_TH
REAL, DIMENSION(KSIZE), INTENT(INOUT) :: PA_RC
REAL, DIMENSION(KSIZE), INTENT(INOUT) :: PA_RI
!
+!$acc declare present(KSIZE,LDSOFT,LDCOMPUTE,PRHODREF,PLVFACT,PLSFACT,PAI,PCJ,PCIT, &
+!$acc& PSSI,PRCT,PRIT,PRCBERI,PA_TH,PA_RC,PA_RI)
+!
!* 0.2 declaration of local variables
!
-REAL, DIMENSION(SIZE(PRHODREF)) :: ZZW
+LOGICAL :: GDSOFT !Workaround of PGI bug with OpenACC (at least up to 18.10 version)
LOGICAL, DIMENSION(SIZE(PRHODREF)) :: GMASK
+REAL, DIMENSION(SIZE(PRHODREF)) :: ZZW
+!$acc declare create(GMASK,ZZW)
!
!-------------------------------------------------------------------------------
!
@@ -123,11 +127,16 @@ IF (MPPDB_INITIALIZED) THEN
CALL MPPDB_CHECK1D(PA_RI,"ICE4_FAST_RI beg:PA_RI",PRECISION)
END IF
!
+!$acc kernels
+GDSOFT = LDSOFT !Workaround of PGI bug with OpenACC (at least up to 18.10 version)
+!$acc end kernels
+!
!* 7.2 Bergeron-Findeisen effect: RCBERI
!
+!$acc kernels
GMASK(:)=PSSI(:)>0. .AND. PRCT(:)>XRTMIN(2) .AND. PRIT(:)>XRTMIN(4) .AND. &
&PCIT(:)>0. .AND. LDCOMPUTE(:)
-IF(LDSOFT) THEN
+IF(GDSOFT) THEN
WHERE(.NOT. GMASK(:))
PRCBERI(:) = 0.
END WHERE
@@ -148,6 +157,7 @@ ENDIF
PA_RC(:) = PA_RC(:) - PRCBERI(:)
PA_RI(:) = PA_RI(:) + PRCBERI(:)
PA_TH(:) = PA_TH(:) + PRCBERI(:)*(PLSFACT(:)-PLVFACT(:))
+!$acc end kernels
!
IF (MPPDB_INITIALIZED) THEN
!Check all INOUT arrays
diff --git a/src/MNH/ice4_fast_rs.f90 b/src/MNH/ice4_fast_rs.f90
index a9ad33f5000bcfde3a233ab2fc3704d9943087d9..6094f023666ab8ad98129614356d88590e424c6a 100644
--- a/src/MNH/ice4_fast_rs.f90
+++ b/src/MNH/ice4_fast_rs.f90
@@ -79,10 +79,14 @@ SUBROUTINE ICE4_FAST_RS(KSIZE, LDSOFT, LDCOMPUTE, &
!* 0. DECLARATIONS
! ------------
!
-USE MODD_CST
-USE MODD_RAIN_ICE_PARAM
-USE MODD_RAIN_ICE_DESCR
-USE MODD_PARAM_ICE, ONLY : LEVLIMIT, CSNOWRIMING
+USE MODD_CST, ONLY: XALPI,XALPW,XBETAI,XBETAW,XCI,XCL,XCPV,XESTT,XGAMI,XGAMW,XLMTT,XLVTT,XMD,XMV,XRV,XTT
+USE MODD_PARAM_ICE, ONLY: LEVLIMIT, CSNOWRIMING
+USE MODD_RAIN_ICE_DESCR, ONLY: XBS,XCEXVT,XCXS,XRTMIN
+USE MODD_RAIN_ICE_PARAM, ONLY: NACCLBDAR,NACCLBDAS,NGAMINC,X0DEPS,X1DEPS,XACCINTP1R,XACCINTP1S,XACCINTP2R,XACCINTP2S, &
+ XCRIMSG,XCRIMSS,XEX0DEPS,XEX1DEPS,XEXCRIMSG,XEXCRIMSS,XEXSRIMCG,XEXSRIMCG2,XFRACCSS, &
+ XFSACCRG,XFSCVMG,XGAMINC_RIM1,XGAMINC_RIM1,XGAMINC_RIM2,XGAMINC_RIM4,XKER_RACCS, &
+ XKER_RACCSS,XKER_SACCRG,XLBRACCS1,XLBRACCS2,XLBRACCS3,XLBSACCR1,XLBSACCR2,XLBSACCR3, &
+ XRIMINTP1,XRIMINTP2,XSRIMCG,XSRIMCG2,XSRIMCG3
!
#ifdef MNH_BITREP
USE MODI_BITREP
@@ -127,16 +131,34 @@ REAL, DIMENSION(KSIZE), INTENT(INOUT) :: PA_RR
REAL, DIMENSION(KSIZE), INTENT(INOUT) :: PA_RS
REAL, DIMENSION(KSIZE), INTENT(INOUT) :: PA_RG
!
+!$acc declare present(KSIZE,LDSOFT,LDCOMPUTE,PRHODREF,PLVFACT,PLSFACT,PPRES,PDV,PKA,PCJ, &
+!$acc& PLBDAR,PLBDAS,PT,PRVT,PRCT,PRRT,PRST,PRIAGGS,PRCRIMSS,PRCRIMSG,PRSRIMCG, &
+!$acc& PRRACCSS,PRRACCSG,PRSACCRG,PRSMLTG,PRCMLTSR,PRS_TEND,PA_TH, &
+!$acc& PA_RC,PA_RR,PA_RS,PA_RG)
+!
!* 0.2 declaration of local variables
!
INTEGER, PARAMETER :: IRCRIMS=1, IRCRIMSS=2, IRSRIMCG=3, IRRACCS=4, IRRACCSS=5, IRSACCRG=6
!
-LOGICAL, DIMENSION(SIZE(PRHODREF)) :: GRIM, GACC, GMASK
-INTEGER :: IGRIM, IGACC
-REAL, DIMENSION(SIZE(PRHODREF)) :: ZVEC1, ZVEC2, ZVEC3
+INTEGER :: IGRIM, IGACC
+INTEGER :: IDX, JJ
INTEGER, DIMENSION(SIZE(PRHODREF)) :: IVEC1, IVEC2
-REAL, DIMENSION(SIZE(PRHODREF)) :: ZZW, ZZW2, ZZW6, ZFREEZ_RATE
-INTEGER :: JJ
+LOGICAL :: GDSOFT,GEVLIMIT !Workaround of PGI bug with OpenACC (at least up to 18.10 version)
+LOGICAL, DIMENSION(SIZE(PRHODREF)) :: GRIM, GACC, GMASK
+REAL, DIMENSION(SIZE(PRHODREF)) :: ZVEC1, ZVEC2, ZVEC3
+REAL, DIMENSION(SIZE(PRHODREF)) :: ZZW, ZZW2, ZZW6, ZFREEZ_RATE
+!
+!$acc declare create(IGRIM,IGACC,IDX,JJ,IVEC1,IVEC2,GRIM,GACC,GMASK,ZVEC1,ZVEC2,ZVEC3,ZZW,ZZW2,ZZW6,ZFREEZ_RATE)
+!
+!$acc declare copyin(XALPI,XALPW,XBETAI,XBETAW,XCI,XCL,XCPV,XESTT,XGAMI,XGAMW,XLMTT,XLVTT,XMD,XMV,XRV,XTT, &
+!$acc& LEVLIMIT, CSNOWRIMING, &
+!$acc& XBS,XCEXVT,XCXS,XRTMIN, &
+!$acc& NACCLBDAR,NACCLBDAS,NGAMINC,X0DEPS,X1DEPS,XACCINTP1R,XACCINTP1S,XACCINTP2R,XACCINTP2S, &
+!$acc& XCRIMSG,XCRIMSS,XEX0DEPS,XEX1DEPS,XEXCRIMSG,XEXCRIMSS,XEXSRIMCG,XEXSRIMCG2,XFRACCSS, &
+!$acc& XFSACCRG,XFSCVMG,XGAMINC_RIM1,XGAMINC_RIM1,XGAMINC_RIM2,XGAMINC_RIM4,XKER_RACCS, &
+!$acc& XKER_RACCSS,XKER_SACCRG,XLBRACCS1,XLBRACCS2,XLBRACCS3,XLBSACCR1,XLBSACCR2,XLBSACCR3, &
+!$acc& XRIMINTP1,XRIMINTP2,XSRIMCG,XSRIMCG2,XSRIMCG3)
+!
!-------------------------------------------------------------------------------
!
IF (MPPDB_INITIALIZED) THEN
@@ -168,14 +190,19 @@ IF (MPPDB_INITIALIZED) THEN
CALL MPPDB_CHECK1D(PA_RG,"ICE4_FAST_RS beg:PA_RG",PRECISION)
END IF
!
+!$acc kernels
+GDSOFT = LDSOFT;GEVLIMIT = LEVLIMIT !Workaround of PGI bug with OpenACC (at least up to 18.10 version)
+!$acc end kernels
+!
!* 5.0 maximum freezing rate
!
+!$acc kernels
ZFREEZ_RATE(:)=0.
GMASK(:)=PRST(:)>XRTMIN(5) .AND. LDCOMPUTE(:)
WHERE(GMASK(:))
ZFREEZ_RATE(:)=PRVT(:)*PPRES(:)/((XMV/XMD)+PRVT(:)) ! Vapor pressure
END WHERE
-IF(LEVLIMIT) THEN
+IF(GEVLIMIT) THEN
WHERE(GMASK(:))
#ifndef MNH_BITREP
ZFREEZ_RATE(:)=MIN(ZFREEZ_RATE(:), EXP(XALPI-XBETAI/PT(:)-XGAMI*ALOG(PT(:)))) ! min(ev, es_i(T))
@@ -209,7 +236,7 @@ END WHERE
GRIM(:) = PRCT(:)>XRTMIN(2) .AND. PRST(:)>XRTMIN(5) .AND. LDCOMPUTE(:)
!
! Collection of cloud droplets by snow: this rate is used for riming (T<0) and for conversion/melting (T>0)
-IF(LDSOFT) THEN
+IF(GDSOFT) THEN
WHERE(.NOT. GRIM(:))
PRS_TEND(:, IRCRIMS)=0.
PRS_TEND(:, IRCRIMSS)=0.
@@ -225,7 +252,17 @@ ELSE
!
! 5.1.1 select the PLBDAS
!
+#ifndef _OPENACC
ZVEC1(1:IGRIM) = PACK( PLBDAS(:),MASK=GRIM(:) )
+#else
+ IDX = 0
+ DO JJ = 1,KSIZE
+ IF (GRIM(JJ)) THEN
+ IDX = IDX+1
+ ZVEC1(IDX) = PLBDAS(JJ)
+ END IF
+ END DO
+#endif
!
! 5.1.2 find the next lower indice for the PLBDAS in the geometrical
! set of Lbda_s used to tabulate some moments of the incomplete
@@ -245,7 +282,19 @@ ELSE
!
ZVEC1(1:IGRIM) = XGAMINC_RIM1( IVEC2(1:IGRIM)+1 )* ZVEC2(1:IGRIM) &
- XGAMINC_RIM1( IVEC2(1:IGRIM) )*(ZVEC2(1:IGRIM) - 1.0)
+#ifndef _OPENACC
ZZW(:) = UNPACK( VECTOR=ZVEC1(1:IGRIM),MASK=GRIM,FIELD=0.0 )
+#else
+ IDX = 0
+ DO JJ=1,KSIZE
+ IF (GRIM(JJ)) THEN
+ IDX = IDX + 1
+ ZZW(JJ) = ZVEC1(IDX)
+ ELSE
+ ZZW(JJ) = 0.
+ END IF
+ END DO
+#endif
!
! 5.1.4 riming of the small sized aggregates
!
@@ -266,11 +315,35 @@ ELSE
!
ZVEC1(1:IGRIM) = XGAMINC_RIM2( IVEC2(1:IGRIM)+1 )* ZVEC2(1:IGRIM) &
- XGAMINC_RIM2( IVEC2(1:IGRIM) )*(ZVEC2(1:IGRIM) - 1.0)
+#ifndef _OPENACC
ZZW(:) = UNPACK( VECTOR=ZVEC1(1:IGRIM),MASK=GRIM,FIELD=0.0 )
+#else
+ IDX = 0
+ DO JJ=1,KSIZE
+ IF (GRIM(JJ)) THEN
+ IDX = IDX + 1
+ ZZW(JJ) = ZVEC1(IDX)
+ ELSE
+ ZZW(JJ) = 0.
+ END IF
+ END DO
+#endif
ZVEC1(1:IGRIM) = XGAMINC_RIM4( IVEC2(1:IGRIM)+1 )* ZVEC2(1:IGRIM) &
- XGAMINC_RIM4( IVEC2(1:IGRIM) )*(ZVEC2(1:IGRIM) - 1.0)
+#ifndef _OPENACC
ZZW2(:) = UNPACK( VECTOR=ZVEC1(1:IGRIM),MASK=GRIM,FIELD=0.0)
+#else
+ IDX = 0
+ DO JJ=1,KSIZE
+ IF (GRIM(JJ)) THEN
+ IDX = IDX + 1
+ ZZW2(JJ) = ZVEC1(IDX)
+ ELSE
+ ZZW2(JJ) = 0.
+ END IF
+ END DO
+#endif
!
! 5.1.6 riming-conversion of the large sized aggregates into graupeln
!
@@ -322,7 +395,8 @@ WHERE(GRIM(:))
ZFREEZ_RATE(:) = MAX(0., ZFREEZ_RATE(:)-PRCRIMSG(:))
PRSRIMCG(:) = ZZW(:) * PRS_TEND(:, IRSRIMCG)
END WHERE
-WHERE(PRCRIMSG(:)<=0.)
+GMASK(:) = PRCRIMSG(:)<=0.
+WHERE(GMASK(:))
PRCRIMSG(:)=0.
PRSRIMCG(:)=0.
END WHERE
@@ -337,8 +411,9 @@ PA_TH(:) = PA_TH(:) + PRCRIMSG(:)*(PLSFACT(:)-PLVFACT(:))
!* 5.2 rain accretion onto the aggregates
!
GACC(:) = PRRT(:)>XRTMIN(3) .AND. PRST(:)>XRTMIN(5) .AND. LDCOMPUTE(:)
-IF(LDSOFT) THEN
- WHERE(.NOT. GACC(:))
+IF(GDSOFT) THEN
+ GMASK(:) = .NOT. GACC(:)
+ WHERE(GMASK(:))
PRS_TEND(:, IRRACCS)=0.
PRS_TEND(:, IRRACCSS)=0.
PRS_TEND(:, IRSACCRG)=0.
@@ -353,8 +428,19 @@ ELSE
!
! 5.2.1 select the (PLBDAS,PLBDAR) couplet
!
+#ifndef _OPENACC
ZVEC1(1:IGACC) = PACK( PLBDAS(:),MASK=GACC(:) )
ZVEC2(1:IGACC) = PACK( PLBDAR(:),MASK=GACC(:) )
+#else
+ IDX = 0
+ DO JJ = 1,KSIZE
+ IF (GACC(JJ)) THEN
+ IDX = IDX+1
+ ZVEC1(IDX) = PLBDAS(JJ)
+ ZVEC2(IDX) = PLBDAR(JJ)
+ END IF
+ END DO
+#endif
!
! 5.2.2 find the next lower indice for the PLBDAS and for the PLBDAR
! in the geometrical set of (Lbda_s,Lbda_r) couplet use to
@@ -389,7 +475,19 @@ ELSE
- XKER_RACCSS(IVEC1(JJ) ,IVEC2(JJ) )*(ZVEC2(JJ) - 1.0) ) &
* (ZVEC1(JJ) - 1.0)
END DO
+#ifndef _OPENACC
ZZW(:) = UNPACK( VECTOR=ZVEC3(1:IGACC),MASK=GACC,FIELD=0.0 )
+#else
+ IDX = 0
+ DO JJ=1,KSIZE
+ IF (GACC(JJ)) THEN
+ IDX = IDX + 1
+ ZZW(JJ) = ZVEC3(IDX)
+ ELSE
+ ZZW(JJ) = 0.
+ END IF
+ END DO
+#endif
!
! 5.2.4 raindrop accretion on the small sized aggregates
!
@@ -421,7 +519,19 @@ ELSE
- XKER_RACCS(IVEC1(JJ) ,IVEC2(JJ) )*(ZVEC2(JJ) - 1.0) ) &
* (ZVEC1(JJ) - 1.0)
END DO
+#ifndef _OPENACC
ZZW(:) = UNPACK( VECTOR=ZVEC3(1:IGACC),MASK=GACC(:),FIELD=0.0 )
+#else
+ IDX = 0
+ DO JJ=1,KSIZE
+ IF (GACC(JJ)) THEN
+ IDX = IDX + 1
+ ZZW(JJ) = ZVEC3(IDX)
+ ELSE
+ ZZW(JJ) = 0.
+ END IF
+ END DO
+#endif
WHERE(GACC(:))
PRS_TEND(:, IRRACCS) = ZZW(:)*ZZW6(:)
END WHERE
@@ -436,7 +546,19 @@ ELSE
- XKER_SACCRG(IVEC2(JJ) ,IVEC1(JJ) )*(ZVEC1(JJ) - 1.0) ) &
* (ZVEC2(JJ) - 1.0)
END DO
+#ifndef _OPENACC
ZZW(:) = UNPACK( VECTOR=ZVEC3(1:IGACC),MASK=GACC,FIELD=0.0 )
+#else
+ IDX = 0
+ DO JJ=1,KSIZE
+ IF (GACC(JJ)) THEN
+ IDX = IDX + 1
+ ZZW(JJ) = ZVEC3(IDX)
+ ELSE
+ ZZW(JJ) = 0.
+ END IF
+ END DO
+#endif
!
! 5.2.6 raindrop accretion-conversion of the large sized aggregates
! into graupeln
@@ -471,7 +593,8 @@ WHERE(GACC(:))
ZFREEZ_RATE(:) = MAX(0., ZFREEZ_RATE(:)-PRRACCSG(:))
PRSACCRG(:)=ZZW(:) * PRS_TEND(:, IRSACCRG)
END WHERE
-WHERE(PRRACCSG(:)<=0.)
+GMASK(:) = PRRACCSG(:)<=0.
+WHERE(GMASK(:))
PRRACCSG(:)=0.
PRSACCRG(:)=0.
END WHERE
@@ -487,7 +610,7 @@ PA_TH(:) = PA_TH(:) + PRRACCSG(:)*(PLSFACT(:)-PLVFACT(:))
!* 5.3 Conversion-Melting of the aggregates
!
GMASK(:)=PRST(:)>XRTMIN(5) .AND. PT(:)>XTT .AND. LDCOMPUTE(:)
-IF(LDSOFT) THEN
+IF(GDSOFT) THEN
WHERE(.NOT. GMASK(:))
PRSMLTG(:) = 0.
PRCMLTSR(:) = 0.
@@ -498,7 +621,7 @@ ELSE
WHERE(GMASK(:))
PRSMLTG(:) = PRVT(:)*PPRES(:)/((XMV/XMD)+PRVT(:)) ! Vapor pressure
END WHERE
- IF(LEVLIMIT) THEN
+ IF(GEVLIMIT) THEN
WHERE(GMASK(:))
#ifndef MNH_BITREP
PRSMLTG(:)=MIN(PRSMLTG(:), EXP(XALPW-XBETAW/PT(:)-XGAMW*ALOG(PT(:)))) ! min(ev, es_w(T))
@@ -540,6 +663,7 @@ PA_RS(:) = PA_RS(:) - PRSMLTG(:)
PA_RG(:) = PA_RG(:) + PRSMLTG(:)
PA_RC(:) = PA_RC(:) - PRCMLTSR(:)
PA_RR(:) = PA_RR(:) + PRCMLTSR(:)
+!$acc end kernels
!
IF (MPPDB_INITIALIZED) THEN
!Check all INOUT arrays
diff --git a/src/MNH/ice4_nucleation.f90 b/src/MNH/ice4_nucleation.f90
index ad0c38d14b092c48f403a3c83a43e4f1155e7674..d0098d3fdc10ed9df2d4a985f56e88d856596bc8 100644
--- a/src/MNH/ice4_nucleation.f90
+++ b/src/MNH/ice4_nucleation.f90
@@ -48,10 +48,10 @@ SUBROUTINE ICE4_NUCLEATION(KSIZE, LDSOFT, LDCOMPUTE, &
!* 0. DECLARATIONS
! ------------
!
-USE MODD_CST
-USE MODD_RAIN_ICE_PARAM
-USE MODD_RAIN_ICE_DESCR, ONLY : XRTMIN
-USE MODD_PARAM_ICE, ONLY : LFEEDBACKT
+USE MODD_CST, ONLY: XALPI,XALPW,XBETAI,XBETAW,XGAMI,XGAMW,XMD,XMV,XTT
+USE MODD_PARAM_ICE, ONLY: LFEEDBACKT
+USE MODD_RAIN_ICE_PARAM, ONLY: XALPHA1,XALPHA2,XBETA1,XBETA2,XMNU0,XNU10,XNU20
+USE MODD_RAIN_ICE_DESCR, ONLY: XRTMIN
!
USE MODE_MPPDB
!
@@ -79,13 +79,26 @@ REAL, DIMENSION(KSIZE), INTENT(INOUT) :: PB_TH
REAL, DIMENSION(KSIZE), INTENT(INOUT) :: PB_RV
REAL, DIMENSION(KSIZE), INTENT(INOUT) :: PB_RI
!
+!$acc declare present(KSIZE,LDSOFT,LDCOMPUTE,PTHT,PPABST,PRHODREF,PEXN,PLSFACT,PT,PRVT,PCIT,PRVHENI_MR,PB_TH,PB_RV,PB_RI)
+!
!* 0.2 declaration of local variables
!
-REAL, DIMENSION(KSIZE) :: ZW ! work array
+INTEGER :: JI
+LOGICAL :: GDSOFT, GFEEDBACKT !Workaround of PGI bug with OpenACC (at least up to 18.10 version)
LOGICAL, DIMENSION(KSIZE) :: GNEGT ! Test where to compute the HEN process
-REAL, DIMENSION(KSIZE) :: ZZW, & ! Work array
- ZUSW, & ! Undersaturation over water
- ZSSI ! Supersaturation over ice
+LOGICAL, DIMENSION(KSIZE) :: GWORK, GWORK2
+REAL, DIMENSION(KSIZE) :: ZW ! work array
+REAL, DIMENSION(KSIZE) :: ZZW, & ! Work array
+ ZUSW, & ! Undersaturation over water
+ ZSSI ! Supersaturation over ice
+!
+!$acc declare create(JI,GNEGT,GWORK,GWORK2,ZW,ZZW,ZUSW,ZSSI)
+!
+!$acc declare copyin(XALPI,XALPW,XBETAI,XBETAW,XGAMI,XGAMW,XMD,XMV,XTT, &
+!$acc& LFEEDBACKT, &
+!$acc& XALPHA1,XALPHA2,XBETA1,XBETA2,XMNU0,XNU10,XNU20, &
+!$acc& XRTMIN)
+!
!-------------------------------------------------------------------------------
IF (MPPDB_INITIALIZED) THEN
!Check all IN arrays
@@ -103,8 +116,15 @@ IF (MPPDB_INITIALIZED) THEN
CALL MPPDB_CHECK1D(PB_RI,"ICE4_NUCLEATION beg:PB_RI",PRECISION)
END IF
!
+!Workaround of PGI bug with OpenACC (at least up to 18.10 version)
+GFEEDBACKT = LFEEDBACKT
+!$acc kernels
+GDSOFT = LDSOFT
+!$acc end kernels
+!
+!$acc kernels
PRVHENI_MR(:)=0.
-IF(.NOT. LDSOFT) THEN
+IF(.NOT. GDSOFT) THEN
GNEGT(:)=PT(:)<XTT .AND. PRVT>XRTMIN(1) .AND. LDCOMPUTE(:)
PRVHENI_MR(:)=0.
ZSSI(:)=0.
@@ -136,20 +156,44 @@ IF(.NOT. LDSOFT) THEN
ZSSI(:)=MIN(ZSSI(:), ZUSW(:)) ! limitation of SSi according to SSw=0
END WHERE
ZZW(:)=0.
+ GWORK(:) = GNEGT(:) .AND. PT(:)<XTT-5.0 .AND. ZSSI(:)>0.0
+ GWORK2(:) = GNEGT(:) .AND. PT(:)<=XTT-2.0 .AND. PT(:)>=XTT-5.0 .AND. ZSSI(:)>0.0
#ifndef MNH_BITREP
- WHERE(GNEGT(:) .AND. PT(:)<XTT-5.0 .AND. ZSSI(:)>0.0 )
+#ifndef _OPENACC
+ WHERE(GWORK(:))
ZZW(:)=XNU20*EXP(XALPHA2*ZSSI(:)-XBETA2)
- ELSEWHERE(GNEGT(:) .AND. PT(:)<=XTT-2.0 .AND. PT(:)>=XTT-5.0 .AND. ZSSI(:)>0.0)
+ ELSEWHERE(GWORK2(:))
ZZW(:)=MAX(XNU20*EXP(-XBETA2 ), &
XNU10*EXP(-XBETA1*(PT(:)-XTT))*(ZSSI(:)/ZUSW(:))**XALPHA1)
END WHERE
#else
- WHERE(GNEGT(:) .AND. PT(:)<XTT-5.0 .AND. ZSSI(:)>0.0 )
+ DO JI=1,KSIZE
+ IF (GWORK(JI)) THEN
+ ZZW(JI)=XNU20*EXP(XALPHA2*ZSSI(JI)-XBETA2)
+ ELSE IF (GWORK2(JI)) THEN
+ ZZW(JI)=MAX(XNU20*(-XBETA2 ), &
+ XNU10*EXP(-XBETA1*(PT(JI)-XTT))*(ZSSI(JI)/ZUSW(JI)**XALPHA1)
+ END IF
+ END DO
+#endif
+#else
+#ifndef _OPENACC
+ WHERE(GWORK(:))
ZZW(:)=XNU20*BR_EXP(XALPHA2*ZSSI(:)-XBETA2)
- ELSEWHERE(GNEGT(:) .AND. PT(:)<=XTT-2.0 .AND. PT(:)>=XTT-5.0 .AND. ZSSI(:)>0.0)
+ ELSEWHERE(GWORK2(:))
ZZW(:)=MAX(XNU20*BR_EXP(-XBETA2 ), &
XNU10*BR_EXP(-XBETA1*(PT(:)-XTT))*BR_POW(ZSSI(:)/ZUSW(:),XALPHA1))
END WHERE
+#else
+ DO JI=1,KSIZE
+ IF (GWORK(JI)) THEN
+ ZZW(JI)=XNU20*BR_EXP(XALPHA2*ZSSI(JI)-XBETA2)
+ ELSE IF (GWORK2(JI)) THEN
+ ZZW(JI)=MAX(XNU20*BR_EXP(-XBETA2 ), &
+ XNU10*BR_EXP(-XBETA1*(PT(JI)-XTT))*BR_POW(ZSSI(JI)/ZUSW(JI),XALPHA1))
+ END IF
+ END DO
+#endif
#endif
WHERE(GNEGT(:))
ZZW(:)=ZZW(:)-PCIT(:)
@@ -163,7 +207,7 @@ IF(.NOT. LDSOFT) THEN
PRVHENI_MR(:)=MIN(PRVT(:), PRVHENI_MR(:))
END WHERE
!Limitation due to 0 crossing of temperature
- IF(LFEEDBACKT) THEN
+ IF(GFEEDBACKT) THEN
ZW(:)=0.
WHERE(GNEGT(:))
ZW(:)=MIN(PRVHENI_MR(:), &
@@ -180,6 +224,7 @@ IF(.NOT. LDSOFT) THEN
PB_RV(:)=PB_RV(:) - PRVHENI_MR(:)
PB_TH(:)=PB_TH(:) + PRVHENI_MR(:)*PLSFACT(:)
ENDIF
+!$acc end kernels
!
IF (MPPDB_INITIALIZED) THEN
!Check all OUT arrays
diff --git a/src/MNH/ice4_nucleation_wrapper.f90 b/src/MNH/ice4_nucleation_wrapper.f90
index a8a99c998e404fd1bd979f3bb120207672cf23dd..7961112de55ca03da1aa2a92f1a47d874ec63aa2 100644
--- a/src/MNH/ice4_nucleation_wrapper.f90
+++ b/src/MNH/ice4_nucleation_wrapper.f90
@@ -64,9 +64,13 @@ REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PRVT ! Water vapor m.r. at t
REAL, DIMENSION(KIT,KJT,KKT), INTENT(INOUT) :: PCIT ! Pristine ice n.c. at t
REAL, DIMENSION(KIT,KJT,KKT), INTENT(OUT) :: PRVHENI_MR ! Mixing ratio change due to the heterogeneous nucleation
!
+!$acc declare present(KIT,KJT,KKT,LDMASK,PTHT,PPABST,PRHODREF,PEXN,PLSFACT,PT,PRVT,PCIT,PRVHENI_MR)
+!
!* 0.2 declaration of local variables
!
INTEGER :: JL ! and PACK intrinsics
+LOGICAL :: GDSOFT
+LOGICAL, DIMENSION(:),ALLOCATABLE :: GLDCOMPUTE
LOGICAL, DIMENSION(SIZE(PRVT,1),SIZE(PRVT,2),SIZE(PRVT,3)) :: GNEGT ! Test where to compute the HEN process
INTEGER :: INEGT, INEGT_TMP
INTEGER, DIMENSION(:),ALLOCATABLE :: I1,I2,I3 ! Used to replace the COUNT
@@ -80,6 +84,11 @@ REAL, DIMENSION(:),ALLOCATABLE :: ZZT, & ! Temperature
ZLSFACT, &
ZRVHENI_MR, &
ZB_TH, ZB_RV, ZB_RI
+!
+!$acc declare create(JL,GDSOFT,GLDCOMPUTE,GNEGT,INEGT,INEGT_TMP,I1,I2,I3,ZZT,ZPRES,ZRVT,ZCIT,ZTHT, &
+!$acc& ZRHODREF,ZEXN,ZLSFACT,ZRVHENI_MR,ZB_TH,ZB_RV,ZB_RI)
+!
+!$acc declare copyin(XTT)
!-------------------------------------------------------------------------------
!
IF (MPPDB_INITIALIZED) THEN
@@ -102,6 +111,7 @@ END IF
GNEGT(:,:,:)=PT(:,:,:)<XTT .AND. LDMASK
INEGT = COUNT(GNEGT(:,:,:))
!
+ALLOCATE(GLDCOMPUTE(INEGT))
ALLOCATE(I1(INEGT),I2(INEGT),I3(INEGT))
ALLOCATE(ZZT(INEGT))
ALLOCATE(ZPRES(INEGT))
@@ -134,7 +144,9 @@ IF(INEGT>0) THEN
ZEXN(JL)=PEXN(I1(JL), I2(JL), I3(JL))
ZLSFACT(JL)=PLSFACT(I1(JL), I2(JL), I3(JL))
ENDDO
- CALL ICE4_NUCLEATION(INEGT, .FALSE., ZZT(:)<XTT, &
+ GDSOFT = .FALSE.
+ GLDCOMPUTE(:) = ZZT(:)<XTT
+ CALL ICE4_NUCLEATION(INEGT, GDSOFT, GLDCOMPUTE, &
ZTHT, ZPRES, ZRHODREF, ZEXN, ZLSFACT, ZZT, &
ZRVT, &
ZCIT, ZRVHENI_MR, ZB_TH, ZB_RV, ZB_RI)
@@ -142,6 +154,7 @@ IF(INEGT>0) THEN
PCIT(:,:,:)=UNPACK(ZCIT(:), MASK=GNEGT(:,:,:), FIELD=PCIT(:,:,:))
END IF
!
+DEALLOCATE(GLDCOMPUTE)
DEALLOCATE(I1,I2,I3)
DEALLOCATE(ZZT,ZPRES,ZRVT,ZCIT,ZTHT,ZRHODREF,ZEXN,ZLSFACT,ZRVHENI_MR,ZB_TH,ZB_RV,ZB_RI)
!
diff --git a/src/MNH/ice4_rainfr_vert.f90 b/src/MNH/ice4_rainfr_vert.f90
index f67a05aa32a5e6ea77ad65313f05de2524264044..0c349e68ca538f35b09acb48729d8503f14e7872 100644
--- a/src/MNH/ice4_rainfr_vert.f90
+++ b/src/MNH/ice4_rainfr_vert.f90
@@ -1,14 +1,14 @@
-!MNH_LIC Copyright 1994-2014 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC Copyright 1994-2018 CNRS, Meteo-France and Universite Paul Sabatier
!MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence
-!MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt
+!MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt
!MNH_LIC for details. version 1.
MODULE MODI_ICE4_RAINFR_VERT
INTERFACE
SUBROUTINE ICE4_RAINFR_VERT(KIB, KIE, KIT, KJB, KJE, KJT, KKB, KKE, KKT, KKL, PPRFR, PRR)
IMPLICIT NONE
-INTEGER, INTENT(IN) :: KIB, KIE, KIT, KJB, KJE, KJT, KKB, KKE, KKT, KKL
-REAL, DIMENSION(KIT,KJT,KKT), INTENT(OUT) :: PPRFR !Precipitation fraction
-REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PRR !Rain field
+INTEGER, INTENT(IN) :: KIB, KIE, KIT, KJB, KJE, KJT, KKB, KKE, KKT, KKL
+REAL, DIMENSION(KIT,KJT,KKT), INTENT(INOUT) :: PPRFR !Precipitation fraction
+REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PRR !Rain field
END SUBROUTINE ICE4_RAINFR_VERT
END INTERFACE
END MODULE MODI_ICE4_RAINFR_VERT
@@ -38,21 +38,30 @@ IMPLICIT NONE
!
!* 0.1 Declarations of dummy arguments :
!
-INTEGER, INTENT(IN) :: KIB, KIE, KIT, KJB, KJE, KJT, KKB, KKE, KKT, KKL
-REAL, DIMENSION(KIT,KJT,KKT), INTENT(OUT) :: PPRFR !Precipitation fraction
-REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PRR !Rain field
+INTEGER, INTENT(IN) :: KIB, KIE, KIT, KJB, KJE, KJT, KKB, KKE, KKT, KKL
+REAL, DIMENSION(KIT,KJT,KKT), INTENT(INOUT) :: PPRFR !Precipitation fraction
+REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PRR !Rain field
+!
+!$acc declare present(KIB, KIE, KIT, KJB, KJE, KJT, KKB, KKE, KKT, KKL, PPRFR, PRR)
!
!* 0.2 declaration of local variables
!
INTEGER :: JI, JJ, JK
!
+!$acc declare copyin(XRTMIN)
+!
!-------------------------------------------------------------------------------
!
-!Check all IN arrays
-CALL MPPDB_CHECK3D(PRR,"ICE4_RAINFR_VERT beg:PRR",PRECISION)
+IF (MPPDB_INITIALIZED) THEN
+ !Check all IN arrays
+ CALL MPPDB_CHECK3D(PRR,"ICE4_RAINFR_VERT beg:PRR",PRECISION)
+ !Check all INOUT arrays
+ CALL MPPDB_CHECK3D(PPRFR,"ICE4_RAINFR_VERT beg:PPRFR",PRECISION)
+END IF
!
-DO JI = KIB,KIE
- DO JJ = KJB, KJE
+!$acc kernels
+DO JJ = KJB, KJE
+ DO JI = KIB,KIE
PPRFR(JI,JJ,KKE)=0.
DO JK=KKE-KKL, KKB, -KKL
IF (PRR(JI,JJ,JK) .GT. XRTMIN(3)) THEN
@@ -66,8 +75,9 @@ DO JI = KIB,KIE
END DO
END DO
END DO
+!$acc end kernels
!
-!Check all OUT arrays
-CALL MPPDB_CHECK3D(PPRFR,"ICE4_RAINFR_VERT beg:PPRFR",PRECISION)
+!Check all INOUT arrays
+CALL MPPDB_CHECK3D(PPRFR,"ICE4_RAINFR_VERT end:PPRFR",PRECISION)
!
END SUBROUTINE ICE4_RAINFR_VERT
diff --git a/src/MNH/ice4_rimltc.f90 b/src/MNH/ice4_rimltc.f90
index 47bf2588521e12745e48a9557eded0ef73ad8358..0f6baa117764080c4563b3d0169f5f1c4e128d19 100644
--- a/src/MNH/ice4_rimltc.f90
+++ b/src/MNH/ice4_rimltc.f90
@@ -48,10 +48,8 @@ SUBROUTINE ICE4_RIMLTC(KSIZE, LDSOFT, LDCOMPUTE, &
!* 0. DECLARATIONS
! ------------
!
-USE MODD_CST
-USE MODD_RAIN_ICE_PARAM
-USE MODD_RAIN_ICE_DESCR
-USE MODD_PARAM_ICE, ONLY : LFEEDBACKT
+USE MODD_CST, ONLY: XTT
+USE MODD_PARAM_ICE, ONLY: LFEEDBACKT
!
USE MODE_MPPDB
!
@@ -73,10 +71,16 @@ REAL, DIMENSION(KSIZE), INTENT(INOUT) :: PB_TH
REAL, DIMENSION(KSIZE), INTENT(INOUT) :: PB_RC
REAL, DIMENSION(KSIZE), INTENT(INOUT) :: PB_RI
!
+!$acc declare present(KSIZE,LDSOFT,LDCOMPUTE,PEXN,PLVFACT,PLSFACT,PT, &
+!$acc& PTHT,PRIT,PRIMLTC_MR,PB_TH,PB_RC,PB_RI)
+!
!* 0.2 declaration of local variables
!
+LOGICAL :: GDSOFT, GFEEDBACKT !Workaround of PGI bug with OpenACC (at least up to 18.10 version)
LOGICAL, DIMENSION(KSIZE) :: GMASK
+!$acc declare create(GMASK)
!
+!$acc declare copyin(XTT,LFEEDBACKT)
!-------------------------------------------------------------------------------
IF (MPPDB_INITIALIZED) THEN
!Check all IN arrays
@@ -94,14 +98,18 @@ END IF
!
!* 7.1 cloud ice melting
!
+!$acc kernels
+GDSOFT=LDSOFT;GFEEDBACKT=LFEEDBACKT !Workaround of PGI bug with OpenACC (at least up to 18.10 version)
+!$acc end kernels
+!$acc kernels
PRIMLTC_MR(:)=0.
-IF(.NOT. LDSOFT) THEN
+IF(.NOT. GDSOFT) THEN
GMASK(:)=PRIT(:)>0. .AND. PT(:)>XTT .AND. LDCOMPUTE(:)
WHERE(GMASK(:))
PRIMLTC_MR(:)=PRIT(:)
END WHERE
- IF(LFEEDBACKT) THEN
+ IF(GFEEDBACKT) THEN
!Limitation due to 0 crossing of temperature
WHERE(GMASK(:))
PRIMLTC_MR(:)=MIN(PRIMLTC_MR(:), MAX(0., (PTHT(:)-XTT/PEXN(:)) / (PLSFACT(:)-PLVFACT(:))))
@@ -111,6 +119,7 @@ ENDIF
PB_RC(:) = PB_RC(:) + PRIMLTC_MR(:)
PB_RI(:) = PB_RI(:) - PRIMLTC_MR(:)
PB_TH(:) = PB_TH(:) - PRIMLTC_MR(:)*(PLSFACT(:)-PLVFACT(:))
+!$acc end kernels
!
IF (MPPDB_INITIALIZED) THEN
!Check all INOUT arrays
diff --git a/src/MNH/ice4_rrhong.f90 b/src/MNH/ice4_rrhong.f90
index 04be4fe593f03cbe7744074722dcb42cbe005fde..be585493a6c439a51a43f15b70af956219408679 100644
--- a/src/MNH/ice4_rrhong.f90
+++ b/src/MNH/ice4_rrhong.f90
@@ -48,10 +48,9 @@ SUBROUTINE ICE4_RRHONG(KSIZE, LDSOFT, LDCOMPUTE, &
!* 0. DECLARATIONS
! ------------
!
-USE MODD_CST
-USE MODD_RAIN_ICE_PARAM
-USE MODD_RAIN_ICE_DESCR
-USE MODD_PARAM_ICE, ONLY : LFEEDBACKT
+USE MODD_CST, ONLY: XTT
+USE MODD_RAIN_ICE_DESCR, ONLY: XRTMIN
+USE MODD_PARAM_ICE, ONLY: LFEEDBACKT
!
USE MODE_MPPDB
!
@@ -73,9 +72,16 @@ REAL, DIMENSION(KSIZE), INTENT(INOUT) :: PB_TH
REAL, DIMENSION(KSIZE), INTENT(INOUT) :: PB_RR
REAL, DIMENSION(KSIZE), INTENT(INOUT) :: PB_RG
!
+!$acc declare present(KSIZE,LDSOFT,LDCOMPUTE,PEXN,PLVFACT,PLSFACT,PT,PRRT,PTHT, &
+!$acc& PRRHONG_MR,PB_TH,PB_RR,PB_RG)
+!
!* 0.2 declaration of local variables
!
+LOGICAL :: GDSOFT, GFEEDBACKT !Workaround of PGI bug with OpenACC (at least up to 18.10 version)
LOGICAL, DIMENSION(SIZE(PRRT)) :: GMASK
+!$acc declare create(GMASK)
+!
+!$acc declare copyin(XTT,XRTMIN,LFEEDBACKT)
!
!-------------------------------------------------------------------------------
IF (MPPDB_INITIALIZED) THEN
@@ -94,13 +100,17 @@ END IF
!
!* 3.3 compute the spontaneous freezing source: RRHONG
!
+!$acc kernels
+GDSOFT=LDSOFT;GFEEDBACKT=LFEEDBACKT !Workaround of PGI bug with OpenACC (at least up to 18.10 version)
+!$acc end kernels
+!$acc kernels
PRRHONG_MR(:) = 0.
-IF(.NOT. LDSOFT) THEN
+IF(GDSOFT) THEN
GMASK(:)=PT(:)<XTT-35.0 .AND. PRRT(:)>XRTMIN(3) .AND. LDCOMPUTE(:)
WHERE(GMASK(:))
PRRHONG_MR(:) = PRRT(:)
ENDWHERE
- IF(LFEEDBACKT) THEN
+ IF(GFEEDBACKT) THEN
!Limitation due to -35 crossing of temperature
WHERE(GMASK(:))
PRRHONG_MR(:)=MIN(PRRHONG_MR(:), MAX(0., ((XTT-35.)/PEXN(:)-PTHT)/(PLSFACT(:)-PLVFACT(:))))
@@ -110,6 +120,7 @@ ENDIF
PB_RG(:) = PB_RG(:) + PRRHONG_MR(:)
PB_RR(:) = PB_RR(:) - PRRHONG_MR(:)
PB_TH(:) = PB_TH(:) + PRRHONG_MR(:)*(PLSFACT(:)-PLVFACT(:))
+!$acc end kernels
!
IF (MPPDB_INITIALIZED) THEN
!Check all INOUT arrays
diff --git a/src/MNH/ice4_rsrimcg_old.f90 b/src/MNH/ice4_rsrimcg_old.f90
index 5333ed3477ed42fd01723d55d90ef02604beed2a..c568907eb32c318c837e60b921fc931352db25b5 100644
--- a/src/MNH/ice4_rsrimcg_old.f90
+++ b/src/MNH/ice4_rsrimcg_old.f90
@@ -46,10 +46,10 @@ SUBROUTINE ICE4_RSRIMCG_OLD(KSIZE, LDSOFT, LDCOMPUTE, &
!* 0. DECLARATIONS
! ------------
!
-USE MODD_CST
-USE MODD_RAIN_ICE_PARAM
-USE MODD_RAIN_ICE_DESCR
-USE MODD_PARAM_ICE, ONLY : CSNOWRIMING
+USE MODD_CST, ONLY: XTT
+USE MODD_PARAM_ICE, ONLY: CSNOWRIMING
+USE MODD_RAIN_ICE_DESCR, ONLY: XRTMIN
+USE MODD_RAIN_ICE_PARAM, ONLY: NGAMINC,XEXSRIMCG,XGAMINC_RIM2,XRIMINTP1,XRIMINTP2,XSRIMCG
!
#ifdef MNH_BITREP
USE MODI_BITREP
@@ -73,14 +73,23 @@ REAL, DIMENSION(KSIZE), INTENT(OUT) :: PRSRIMCG_MR ! Mr change due to cl
REAL, DIMENSION(KSIZE), INTENT(INOUT) :: PB_RS
REAL, DIMENSION(KSIZE), INTENT(INOUT) :: PB_RG
!
+!$acc declare present(KSIZE,LDSOFT,LDCOMPUTE,PRHODREF,PLBDAS, &
+!$acc& PT,PRCT,PRST,PRSRIMCG_MR,PB_RS,PB_RG)
+!
!* 0.2 declaration of local variables
!
-LOGICAL, DIMENSION(SIZE(PRHODREF)) :: GRIM, GACC, GMASK
-INTEGER :: IGRIM, IGACC
-REAL, DIMENSION(SIZE(PRHODREF)) :: ZVEC1, ZVEC2, ZVEC3
+INTEGER :: IDX,JL
+INTEGER :: IGRIM
INTEGER, DIMENSION(SIZE(PRHODREF)) :: IVEC1, IVEC2
-REAL, DIMENSION(SIZE(PRHODREF)) :: ZZW, ZZW2, ZZW6
-INTEGER :: JJ
+LOGICAL :: GDSOFT !Workaround of PGI bug with OpenACC (at least up to 18.10 version)
+LOGICAL, DIMENSION(SIZE(PRHODREF)) :: GRIM
+REAL, DIMENSION(SIZE(PRHODREF)) :: ZVEC1, ZVEC2
+REAL, DIMENSION(SIZE(PRHODREF)) :: ZZW
+!
+!$acc declare create(GRIM, IGRIM, ZVEC1, ZVEC2, IVEC1, IVEC2, ZZW, JL)
+!
+!$acc declare copyin(XTT,CSNOWRIMING,XRTMIN,NGAMINC,XEXSRIMCG,XGAMINC_RIM2,XRIMINTP1,XRIMINTP2,XSRIMCG)
+!
!-------------------------------------------------------------------------------
!
IF (MPPDB_INITIALIZED) THEN
@@ -100,17 +109,34 @@ END IF
!
!* 5.1 cloud droplet riming of the aggregates
!
+!$acc kernels
+GDSOFT=LDSOFT !Workaround of PGI bug with OpenACC (at least up to 18.10 version)
PRSRIMCG_MR(:)=0.
+!$acc end kernels
!
-IF(.NOT. LDSOFT) THEN
+IF(.NOT. GDSOFT) THEN
+!$acc kernels
GRIM(:) = PRCT(:)>XRTMIN(2) .AND. PRST(:)>XRTMIN(5) .AND. LDCOMPUTE(:) .AND. PT(:)<XTT
IGRIM = COUNT(GRIM(:))
+!$acc end kernels
+!$acc update self(IGRIM)
!
IF(IGRIM>0 .AND. CSNOWRIMING=='OLD ') THEN
+!$acc kernels
!
! 5.1.1 select the PLBDAS
!
+#ifndef _OPENACC
ZVEC1(1:IGRIM) = PACK( PLBDAS(:),MASK=GRIM(:) )
+#else
+ IDX = 0
+ DO JL=1,KSIZE
+ IF (GRIM(JL)) THEN
+ IDX = IDX+1
+ ZVEC1(IDX) = PLBDAS(JL)
+ END IF
+ END DO
+#endif
!
! 5.1.2 find the next lower indice for the PLBDAS in the geometrical
! set of Lbda_s used to tabulate some moments of the incomplete
@@ -131,8 +157,19 @@ IF(.NOT. LDSOFT) THEN
!
ZVEC1(1:IGRIM) = XGAMINC_RIM2( IVEC2(1:IGRIM)+1 )* ZVEC2(1:IGRIM) &
- XGAMINC_RIM2( IVEC2(1:IGRIM) )*(ZVEC2(1:IGRIM) - 1.0)
+#ifndef _OPENACC
ZZW(:) = UNPACK( VECTOR=ZVEC1(1:IGRIM),MASK=GRIM,FIELD=0.0 )
-
+#else
+ IDX = 0
+ DO JL=1,KSIZE
+ IF (GRIM(JL)) THEN
+ IDX = IDX+1
+ ZZW(JL) = ZVEC1(IDX)
+ ELSE
+ ZZW(JL) = 0.
+ END IF
+ END DO
+#endif
!
! 5.1.6 riming-conversion of the large sized aggregates into graupeln
!
@@ -146,10 +183,13 @@ IF(.NOT. LDSOFT) THEN
* (1.0 - ZZW(:) )/PRHODREF(:)
END WHERE
PRSRIMCG_MR(:)=MIN(PRST(:), PRSRIMCG_MR(:))
+!$acc end kernels
END IF
ENDIF
+!$acc kernels
PB_RS(:) = PB_RS(:) - PRSRIMCG_MR(:)
PB_RG(:) = PB_RG(:) + PRSRIMCG_MR(:)
+!$acc end kernels
!
IF (MPPDB_INITIALIZED) THEN
!Check all INOUT arrays
diff --git a/src/MNH/ice4_sedimentation_split.f90 b/src/MNH/ice4_sedimentation_split.f90
index 74967f73782e2bc1008d0ef73519f59a9d1094e8..a60f53083f9dce9a8f1a716f799bd187717149d6 100644
--- a/src/MNH/ice4_sedimentation_split.f90
+++ b/src/MNH/ice4_sedimentation_split.f90
@@ -79,8 +79,8 @@ USE MODD_PARAM_ICE, ONLY: XSPLIT_MAXCFL
USE MODD_RAIN_ICE_DESCR, ONLY: XALPHAC,XALPHAC2,XCONC_LAND,XCONC_SEA,XCONC_URBAN,XLBC,XNUC,XNUC2
USE MODD_RAIN_ICE_PARAM, ONLY: XFSEDC
!
-USE MODE_MSG
USE MODE_MPPDB
+USE MODE_MSG
!
USE MODI_GAMMA
!
@@ -123,18 +123,20 @@ REAL, DIMENSION(KIT,KJT,KKT), OPTIONAL, INTENT(IN) :: PRHT ! Hail m.r.
REAL, DIMENSION(KIT,KJT,KKT), OPTIONAL, INTENT(INOUT) :: PRHS ! Hail m.r. source
REAL, DIMENSION(KIT,KJT,KKT,KRR), OPTIONAL, INTENT(OUT) :: PFPR ! upper-air precipitation fluxes
!
-!$acc declare copyin(KIB,KIE,KIT,KJB,KJE,KJT,KKB,KKE,KKTB,KKTE,KKT,KKL,PTSTEP,KRR,PVDEPOSC)
-!acc declare copyin(OSEDIC,ODEPOSC)
-!
-!$acc declare present(PDZZ,PRHODREF,PPABST,PTHT,PRHODJ,PRCS,PRCT,PRRS,PRRT,PRIS,PRIT,PRSS,PRST,PRGS,PRGT, &
-!$acc & PINPRC,PINDEP,PINPRR,PINPRI,PINPRS,PINPRG,PSEA,PTOWN,PINPRH,PRHT,PRHS,PFPR)
+!$acc declare present(KIB,KIE,KIT,KJB,KJE,KJT,KKB,KKE,KKTB,KKTE,KKT,KKL,PTSTEP,KRR,OSEDIC,ODEPOSC,PVDEPOSC,PDZZ, &
+!$acc& PRHODREF,PPABST,PTHT,PRHODJ,PRCS,PRCT,PRRS,PRRT,PRIS,PRIT,PRSS,PRST,PRGS,PRGT, &
+!$acc& PINPRC,PINDEP,PINPRR,PINPRI,PINPRS,PINPRG,PSEA,PTOWN,PINPRH,PRHT,PRHS,PFPR)
!
!* 0.2 declaration of local variables
!
!
-INTEGER :: JI,JJ,JK
-LOGICAL :: GPRESENT_PFPR, GPRESENT_PSEA
-REAL :: ZINVTSTEP
+INTEGER :: JI,JJ,JK
+INTEGER :: IRR !Workaround of PGI bug with OpenACC (at least up to 18.10 version)
+LOGICAL :: GDEPOSC, GSEDIC !Workaround of PGI bug with OpenACC (at least up to 18.10 version)
+LOGICAL :: GPRESENT_PFPR, GPRESENT_PSEA
+REAL :: ZINVTSTEP
+REAL, DIMENSION(SIZE(PRHODREF,1),SIZE(PRHODREF,2)) :: ZCONC_TMP ! Weighted concentration
+REAL, DIMENSION(SIZE(PRHODREF,1),SIZE(PRHODREF,2),KKTB:KKTE) :: ZW ! work array
REAL, DIMENSION(SIZE(PRHODREF,1),SIZE(PRHODREF,2),SIZE(PRHODREF,3)) :: ZCONC3D, & ! droplet condensation
& ZRAY, & ! Cloud Mean radius
& ZLBC, & ! XLBC weighted by sea fraction
@@ -146,41 +148,50 @@ REAL, DIMENSION(SIZE(PRHODREF,1),SIZE(PRHODREF,2),SIZE(PRHODREF,3)) :: ZCONC3D,
& ZRST, &
& ZRGT, &
& ZRHT
-REAL, DIMENSION(SIZE(PRHODREF,1),SIZE(PRHODREF,2),KKTB:KKTE) :: ZW ! work array
-REAL, DIMENSION(SIZE(PRHODREF,1),SIZE(PRHODREF,2)) :: ZCONC_TMP ! Weighted concentration
!
-!acc declare create(GPRESENT_PFPR,GPRESENT_PSEA)
-!$acc declare create(ZINVTSTEP,JI,JJ,JK, &
-!$acc ZCONC3D,ZRAY,ZLBC,ZFSEDC,ZPRCS,ZPRRS,ZPRIS,ZPRSS,ZPRGS,ZPRHS, &
-!$acc& ZW,ZRCT,ZRRT,ZRIT,ZRST,ZRGT,ZRHT,ZCONC_TMP) &
-!$acc& copyin(XALPHAC,XALPHAC2,XCONC_LAND,XCONC_SEA,XCONC_URBAN,XFSEDC,XLBC, &
-!$acc& XNUC,XNUC2,XRHOLW,XSPLIT_MAXCFL)
+!$acc declare create(JI,JJ,JK,ZINVTSTEP,ZCONC_TMP,ZW, &
+!$acc& ZCONC3D,ZRAY,ZLBC,ZFSEDC,ZPRCS,ZPRRS,ZPRIS,ZPRSS,ZPRGS,ZPRHS, &
+!$acc& ZRCT,ZRRT,ZRIT,ZRST,ZRGT,ZRHT)
+!
+!$acc declare copyin(XRHOLW, &
+!$acc XSPLIT_MAXCFL, &
+!$acc XALPHAC,XALPHAC2,XCONC_LAND,XCONC_SEA,XCONC_URBAN,XLBC,XNUC,XNUC2, &
+!$acc XFSEDC)
!
!-------------------------------------------------------------------------------
!
!-------------------------------------------------------------------------------
!
-!Check all IN arrays
-CALL MPPDB_CHECK3D(PDZZ,"ICE4_SEDIMENTATION_SPLIT beg:PDZZ",PRECISION)
-CALL MPPDB_CHECK3D(PRHODREF,"ICE4_SEDIMENTATION_SPLIT beg:PRHODREF",PRECISION)
-CALL MPPDB_CHECK3D(PPABST,"ICE4_SEDIMENTATION_SPLIT beg:PPABST",PRECISION)
-CALL MPPDB_CHECK3D(PTHT,"ICE4_SEDIMENTATION_SPLIT beg:PTHT",PRECISION)
-CALL MPPDB_CHECK3D(PRHODJ,"ICE4_SEDIMENTATION_SPLIT beg:PRHODJ",PRECISION)
-CALL MPPDB_CHECK3D(PRCT,"ICE4_SEDIMENTATION_SPLIT beg:PRCT",PRECISION)
-CALL MPPDB_CHECK3D(PRRT,"ICE4_SEDIMENTATION_SPLIT beg:PRRT",PRECISION)
-CALL MPPDB_CHECK3D(PRIT,"ICE4_SEDIMENTATION_SPLIT beg:PRIT",PRECISION)
-CALL MPPDB_CHECK3D(PRST,"ICE4_SEDIMENTATION_SPLIT beg:PRST",PRECISION)
-CALL MPPDB_CHECK3D(PRGT,"ICE4_SEDIMENTATION_SPLIT beg:PRGT",PRECISION)
-IF (PRESENT(PSEA)) CALL MPPDB_CHECK2D(PSEA,"ICE4_SEDIMENTATION_SPLIT beg:PSEA",PRECISION)
-IF (PRESENT(PTOWN)) CALL MPPDB_CHECK2D(PTOWN,"ICE4_SEDIMENTATION_SPLIT beg:PTOWN",PRECISION)
-IF (PRESENT(PRHT)) CALL MPPDB_CHECK3D(PRHT,"ICE4_SEDIMENTATION_SPLIT beg:PRHT",PRECISION)
-!Check all INOUT arrays
-CALL MPPDB_CHECK3D(PRCS,"ICE4_SEDIMENTATION_SPLIT beg:PRCS",PRECISION)
-CALL MPPDB_CHECK3D(PRRS,"ICE4_SEDIMENTATION_SPLIT beg:PRRS",PRECISION)
-CALL MPPDB_CHECK3D(PRIS,"ICE4_SEDIMENTATION_SPLIT beg:PRIS",PRECISION)
-CALL MPPDB_CHECK3D(PRSS,"ICE4_SEDIMENTATION_SPLIT beg:PRSS",PRECISION)
-CALL MPPDB_CHECK3D(PRGS,"ICE4_SEDIMENTATION_SPLIT beg:PRGS",PRECISION)
-IF (PRESENT(PRHS)) CALL MPPDB_CHECK3D(PRHS,"ICE4_SEDIMENTATION_SPLIT beg:PRHS",PRECISION)
+IF (MPPDB_INITIALIZED) THEN
+ !Check all IN arrays
+ CALL MPPDB_CHECK3D(PDZZ,"ICE4_SEDIMENTATION_SPLIT beg:PDZZ",PRECISION)
+ CALL MPPDB_CHECK3D(PRHODREF,"ICE4_SEDIMENTATION_SPLIT beg:PRHODREF",PRECISION)
+ CALL MPPDB_CHECK3D(PPABST,"ICE4_SEDIMENTATION_SPLIT beg:PPABST",PRECISION)
+ CALL MPPDB_CHECK3D(PTHT,"ICE4_SEDIMENTATION_SPLIT beg:PTHT",PRECISION)
+ CALL MPPDB_CHECK3D(PRHODJ,"ICE4_SEDIMENTATION_SPLIT beg:PRHODJ",PRECISION)
+ CALL MPPDB_CHECK3D(PRCT,"ICE4_SEDIMENTATION_SPLIT beg:PRCT",PRECISION)
+ CALL MPPDB_CHECK3D(PRRT,"ICE4_SEDIMENTATION_SPLIT beg:PRRT",PRECISION)
+ CALL MPPDB_CHECK3D(PRIT,"ICE4_SEDIMENTATION_SPLIT beg:PRIT",PRECISION)
+ CALL MPPDB_CHECK3D(PRST,"ICE4_SEDIMENTATION_SPLIT beg:PRST",PRECISION)
+ CALL MPPDB_CHECK3D(PRGT,"ICE4_SEDIMENTATION_SPLIT beg:PRGT",PRECISION)
+ IF (PRESENT(PSEA)) CALL MPPDB_CHECK2D(PSEA,"ICE4_SEDIMENTATION_SPLIT beg:PSEA",PRECISION)
+ IF (PRESENT(PTOWN)) CALL MPPDB_CHECK2D(PTOWN,"ICE4_SEDIMENTATION_SPLIT beg:PTOWN",PRECISION)
+ IF (PRESENT(PRHT)) CALL MPPDB_CHECK3D(PRHT,"ICE4_SEDIMENTATION_SPLIT beg:PRHT",PRECISION)
+ !Check all INOUT arrays
+ CALL MPPDB_CHECK3D(PRCS,"ICE4_SEDIMENTATION_SPLIT beg:PRCS",PRECISION)
+ CALL MPPDB_CHECK3D(PRRS,"ICE4_SEDIMENTATION_SPLIT beg:PRRS",PRECISION)
+ CALL MPPDB_CHECK3D(PRIS,"ICE4_SEDIMENTATION_SPLIT beg:PRIS",PRECISION)
+ CALL MPPDB_CHECK3D(PRSS,"ICE4_SEDIMENTATION_SPLIT beg:PRSS",PRECISION)
+ CALL MPPDB_CHECK3D(PRGS,"ICE4_SEDIMENTATION_SPLIT beg:PRGS",PRECISION)
+ IF (PRESENT(PRHS)) CALL MPPDB_CHECK3D(PRHS,"ICE4_SEDIMENTATION_SPLIT beg:PRHS",PRECISION)
+END IF
+!
+!$acc kernels
+!Workaround of PGI bug with OpenACC (at least up to 18.10 version)
+GDEPOSC = ODEPOSC
+GSEDIC = OSEDIC
+IRR = KRR
+!$acc end kernels
!
IF (PRESENT(PFPR)) THEN
GPRESENT_PFPR = .TRUE.
@@ -193,7 +204,6 @@ IF (PRESENT(PSEA)) THEN
ELSE
GPRESENT_PSEA = .FALSE.
END IF
-!acc update device(GPRESENT_PFPR,GPRESENT_PSEA)
!
! O. Initialization of for sedimentation
!
@@ -209,7 +219,7 @@ END IF
!* 1. Parameters for cloud sedimentation
!
!$acc kernels
-IF (OSEDIC) THEN
+IF (GSEDIC) THEN
ZRAY(:,:,:) = 0.
ZLBC(:,:,:) = XLBC(1)
ZFSEDC(:,:,:) = XFSEDC(1)
@@ -232,6 +242,7 @@ IF (OSEDIC) THEN
ZRAY(:,:,:) = MAX(1.,ZRAY(:,:,:))
ZLBC(:,:,:) = MAX(MIN(XLBC(1),XLBC(2)),ZLBC(:,:,:))
ENDIF
+!$acc end kernels
!
!* 2. compute the fluxes
!
@@ -240,17 +251,18 @@ ENDIF
! For optimization we consider each variable separately
!
! External tendecies
-IF (OSEDIC) ZPRCS(:,:,:) = PRCS(:,:,:)-PRCT(:,:,:)*ZINVTSTEP
+!$acc kernels
+IF (GSEDIC) THEN
+ ZPRCS(:,:,:) = PRCS(:,:,:)-PRCT(:,:,:)*ZINVTSTEP
+ENDIF
ZPRRS(:,:,:) = PRRS(:,:,:)-PRRT(:,:,:)*ZINVTSTEP
ZPRIS(:,:,:) = PRIS(:,:,:)-PRIT(:,:,:)*ZINVTSTEP
ZPRSS(:,:,:) = PRSS(:,:,:)-PRST(:,:,:)*ZINVTSTEP
ZPRGS(:,:,:) = PRGS(:,:,:)-PRGT(:,:,:)*ZINVTSTEP
-!$acc end kernels
-IF ( KRR == 7 ) THEN
-!$acc kernels
+IF ( IRR == 7 ) THEN
ZPRHS(:,:,:) = PRHS(:,:,:)-PRHT(:,:,:)*ZINVTSTEP
-!$acc end kernels
END IF
+!$acc end kernels
!
! mr values inside the time-splitting loop
!$acc kernels
@@ -259,12 +271,10 @@ ZRRT(:,:,:) = PRRT(:,:,:)
ZRIT(:,:,:) = PRIT(:,:,:)
ZRST(:,:,:) = PRST(:,:,:)
ZRGT(:,:,:) = PRGT(:,:,:)
-!$acc end kernels
-IF (KRR==7) THEN
-!$acc kernels
+IF (IRR==7) THEN
ZRHT(:,:,:) = PRHT(:,:,:)
-!$acc end kernels
END IF
+!$acc end kernels
!
!$acc kernels
ZW(:,:,KKTB:KKTE) =1./(PRHODREF(:,:,KKTB:KKTE)* PDZZ(:,:,KKTB:KKTE))
@@ -273,7 +283,7 @@ ZW(:,:,KKTB:KKTE) =1./(PRHODREF(:,:,KKTB:KKTE)* PDZZ(:,:,KKTB:KKTE))
!
!* 2.1 for cloud
!
-IF (OSEDIC) THEN
+IF (GSEDIC) THEN
CALL INTERNAL_SEDIM_SPLI(KIB,KIE,KIT,KJB,KJE,KJT, KKB, KKTB, KKTE, KKT, KKL, KRR, &
&XSPLIT_MAXCFL, &
&PRHODREF, ZW, PDZZ, PPABST, PTHT, PTSTEP, &
@@ -286,7 +296,7 @@ ENDIF
!* 2.1bis DROPLET DEPOSITION AT THE 1ST LEVEL ABOVE GROUND
!
!$acc kernels
-IF (ODEPOSC) THEN
+IF (GDEPOSC) THEN
PINDEP (:,:) = 0.
!$acc loop collapse(2) independent
DO JJ=KJB,KJE
@@ -339,7 +349,7 @@ END IF
!
!* 2.6 for hail
!
-IF (KRR==7) THEN
+IF (IRR==7) THEN
CALL INTERNAL_SEDIM_SPLI(KIB,KIE,KIT,KJB,KJE,KJT, KKB, KKTB, KKTE, KKT, KKL, KRR, &
&XSPLIT_MAXCFL, &
&PRHODREF, ZW, PDZZ, PPABST, PTHT, PTSTEP, &
@@ -351,25 +361,27 @@ ENDIF
!PW:TODO: to remove when everything on GPU in rain_ice_red
!$acc update self(PRCS,PRRS,PRIS,PRSS,PRGS,PRHS,PINPRC,PINDEP,PINPRR,PINPRI,PINPRS,PINPRG,PINPRH,PFPR)
!
-!Check all INOUT arrays
-CALL MPPDB_CHECK3D(PRCS,"ICE4_SEDIMENTATION_SPLIT end:PRCS",PRECISION)
-CALL MPPDB_CHECK3D(PRRS,"ICE4_SEDIMENTATION_SPLIT end:PRRS",PRECISION)
-CALL MPPDB_CHECK3D(PRIS,"ICE4_SEDIMENTATION_SPLIT end:PRIS",PRECISION)
-CALL MPPDB_CHECK3D(PRSS,"ICE4_SEDIMENTATION_SPLIT end:PRSS",PRECISION)
-CALL MPPDB_CHECK3D(PRGS,"ICE4_SEDIMENTATION_SPLIT end:PRGS",PRECISION)
-IF (PRESENT(PRHS)) CALL MPPDB_CHECK3D(PRHS,"ICE4_SEDIMENTATION_SPLIT end:PRHS",PRECISION)
-!Check all OUT arrays
-CALL MPPDB_CHECK2D(PINPRC,"ICE4_SEDIMENTATION_SPLIT end:PINPRC",PRECISION)
-CALL MPPDB_CHECK2D(PINDEP,"ICE4_SEDIMENTATION_SPLIT end:PINDEP",PRECISION)
-CALL MPPDB_CHECK2D(PINPRR,"ICE4_SEDIMENTATION_SPLIT end:PINPRR",PRECISION)
-CALL MPPDB_CHECK2D(PINPRI,"ICE4_SEDIMENTATION_SPLIT end:PINPRI",PRECISION)
-CALL MPPDB_CHECK2D(PINPRS,"ICE4_SEDIMENTATION_SPLIT end:PINPRS",PRECISION)
-CALL MPPDB_CHECK2D(PINPRG,"ICE4_SEDIMENTATION_SPLIT end:PINPRG",PRECISION)
-IF (PRESENT(PINPRH)) CALL MPPDB_CHECK2D(PINPRH,"ICE4_SEDIMENTATION_SPLIT end:PINPRH",PRECISION)
-IF (PRESENT(PFPR)) THEN
- DO JK=1,SIZE(PFPR,4)
- CALL MPPDB_CHECK3D(PFPR(:,:,:,JK),"ICE4_SEDIMENTATION_SPLIT end:PFPR(:,:,:,JK)",PRECISION)
- END DO
+IF (MPPDB_INITIALIZED) THEN
+ !Check all INOUT arrays
+ CALL MPPDB_CHECK3D(PRCS,"ICE4_SEDIMENTATION_SPLIT end:PRCS",PRECISION)
+ CALL MPPDB_CHECK3D(PRRS,"ICE4_SEDIMENTATION_SPLIT end:PRRS",PRECISION)
+ CALL MPPDB_CHECK3D(PRIS,"ICE4_SEDIMENTATION_SPLIT end:PRIS",PRECISION)
+ CALL MPPDB_CHECK3D(PRSS,"ICE4_SEDIMENTATION_SPLIT end:PRSS",PRECISION)
+ CALL MPPDB_CHECK3D(PRGS,"ICE4_SEDIMENTATION_SPLIT end:PRGS",PRECISION)
+ IF (PRESENT(PRHS)) CALL MPPDB_CHECK3D(PRHS,"ICE4_SEDIMENTATION_SPLIT end:PRHS",PRECISION)
+ !Check all OUT arrays
+ CALL MPPDB_CHECK2D(PINPRC,"ICE4_SEDIMENTATION_SPLIT end:PINPRC",PRECISION)
+ CALL MPPDB_CHECK2D(PINDEP,"ICE4_SEDIMENTATION_SPLIT end:PINDEP",PRECISION)
+ CALL MPPDB_CHECK2D(PINPRR,"ICE4_SEDIMENTATION_SPLIT end:PINPRR",PRECISION)
+ CALL MPPDB_CHECK2D(PINPRI,"ICE4_SEDIMENTATION_SPLIT end:PINPRI",PRECISION)
+ CALL MPPDB_CHECK2D(PINPRS,"ICE4_SEDIMENTATION_SPLIT end:PINPRS",PRECISION)
+ CALL MPPDB_CHECK2D(PINPRG,"ICE4_SEDIMENTATION_SPLIT end:PINPRG",PRECISION)
+ IF (PRESENT(PINPRH)) CALL MPPDB_CHECK2D(PINPRH,"ICE4_SEDIMENTATION_SPLIT end:PINPRH",PRECISION)
+ IF (PRESENT(PFPR)) THEN
+ DO JK=1,SIZE(PFPR,4)
+ CALL MPPDB_CHECK3D(PFPR(:,:,:,JK),"ICE4_SEDIMENTATION_SPLIT end:PFPR(:,:,:,JK)",PRECISION)
+ END DO
+ END IF
END IF
!
!
@@ -414,9 +426,10 @@ REAL, DIMENSION(KIT,KJT), INTENT(OUT) :: PINPRX ! instant precip
REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PPRXS ! external tendencie
REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN), OPTIONAL :: PRAY, PLBC, PFSEDC, PCONC3D
REAL, DIMENSION(KIT,KJT,KKT,KRR), INTENT(INOUT), OPTIONAL :: PFPR ! upper-air precipitation fluxes
-!$acc declare present(KIB,KIE,KIT,KJB,KJE,KJT,KKB,KKTB,KKTE,KKT,KKL,KRR,PMAXCFL,PTSTEP, &
-!$acc& PRHODREF,POORHODZ,PDZZ,PPABST,PTHT,PRXT,PRXS,PINPRX,PPRXS, &
-!$acc& PRAY,PLBC,PFSEDC,PCONC3D,PFPR) &
+!
+!$acc declare present(KIB,KIE,KIT,KJB,KJE,KJT,KKB,KKTB,KKTE,KKT,KKL,KRR,PMAXCFL, &
+!$acc& PRHODREF,POORHODZ,PDZZ,PPABST,PTHT,PTSTEP, PRXT,PRXS,PINPRX,PPRXS, &
+!$acc& PRAY,PLBC,PFSEDC,PCONC3D,PFPR) &
!$acc& copyin(KSPE)
!
!* 0.2 declaration of local variables
@@ -434,33 +447,37 @@ REAL, DIMENSION(SIZE(XRTMIN)) :: ZRSMIN
REAL, DIMENSION(SIZE(PRHODREF,1),SIZE(PRHODREF,2)) :: ZREMAINT ! Remaining time until the timestep end
REAL, DIMENSION(SIZE(PRHODREF,1),SIZE(PRHODREF,2),0:SIZE(PRHODREF,3)+1) :: ZWSED ! Sedimentation fluxes
!
-!$acc declare create(ISEDIM,I1,I2,I3,ZREMAINT,ZWSED,ZRSMIN, &
-!$acc& ZFSED,ZEXSED,ZMAX_TSTEP,ZMRCHANGE) &
-!$acc& copyin(XCPD,XP00,XRD,XEXSEDG,XEXSEDH,XEXSEDR,XEXSEDS,XFSEDG,XFSEDH,XFSEDI,XFSEDR,XFSEDS,XRTMIN, &
-!$acc& XCC,XCEXVT,XDC,XLBEXC,XEXCSEDI)
-!
-!Check all IN arrays
-CALL MPPDB_CHECK3D(PRHODREF,"INTERNAL_SEDIM_SPLI beg:PRHODREF",PRECISION)
-CALL MPPDB_CHECK3D(POORHODZ,"INTERNAL_SEDIM_SPLI beg:POORHODZ",PRECISION)
-CALL MPPDB_CHECK3D(PDZZ,"INTERNAL_SEDIM_SPLI beg:PDZZ",PRECISION)
-CALL MPPDB_CHECK3D(PPABST,"INTERNAL_SEDIM_SPLI beg:PPABST",PRECISION)
-CALL MPPDB_CHECK3D(PTHT,"INTERNAL_SEDIM_SPLI beg:PTHT",PRECISION)
-CALL MPPDB_CHECK3D(PPRXS,"INTERNAL_SEDIM_SPLI beg:PPRXS",PRECISION)
-IF (PRESENT(PRAY)) CALL MPPDB_CHECK3D(PRAY,"INTERNAL_SEDIM_SPLI beg:PRAY",PRECISION)
-IF (PRESENT(PLBC)) CALL MPPDB_CHECK3D(PLBC,"INTERNAL_SEDIM_SPLI beg:PLBC",PRECISION)
-IF (PRESENT(PFSEDC)) CALL MPPDB_CHECK3D(PFSEDC,"INTERNAL_SEDIM_SPLI beg:PFSEDC",PRECISION)
-IF (PRESENT(PCONC3D)) CALL MPPDB_CHECK3D(PCONC3D,"INTERNAL_SEDIM_SPLI beg:PCONC3D",PRECISION)
-!Check all INOUT arrays
-CALL MPPDB_CHECK3D(PRXT,"INTERNAL_SEDIM_SPLI beg:PRXT",PRECISION)
-CALL MPPDB_CHECK3D(PRXS,"INTERNAL_SEDIM_SPLI beg:PRXS",PRECISION)
-IF (PRESENT(PFPR)) THEN
- DO JI=1,KRR
- CALL MPPDB_CHECK3D(PFPR(:,:,:,JI),"INTERNAL_SEDIM_SPLI beg:",PRECISION)
- END DO
+!$acc declare create(IDX,ISEDIM,JI,JJ,JK,JL,I1,I2,I3,ZINVTSTEP,ZZWLBDC,ZRAY,ZZT,ZZWLBDA,ZZCC, &
+!$acc& ZFSED,ZEXSED,ZMRCHANGE,ZMAX_TSTEP,ZRSMIN,ZREMAINT,ZWSED)
+!
+!$acc declare copyin(XCPD,XP00,XRD, &
+!$acc& XCC,XCEXVT,XDC,XLBEXC,XRTMIN, &
+!$acc& XEXCSEDI,XEXSEDG,XEXSEDH,XEXSEDR,XEXSEDS,XFSEDG,XFSEDH,XFSEDI,XFSEDR,XFSEDS)
+!
+IF (MPPDB_INITIALIZED) THEN
+ !Check all IN arrays
+ CALL MPPDB_CHECK3D(PRHODREF,"INTERNAL_SEDIM_SPLI beg:PRHODREF",PRECISION)
+ CALL MPPDB_CHECK3D(POORHODZ,"INTERNAL_SEDIM_SPLI beg:POORHODZ",PRECISION)
+ CALL MPPDB_CHECK3D(PDZZ,"INTERNAL_SEDIM_SPLI beg:PDZZ",PRECISION)
+ CALL MPPDB_CHECK3D(PPABST,"INTERNAL_SEDIM_SPLI beg:PPABST",PRECISION)
+ CALL MPPDB_CHECK3D(PTHT,"INTERNAL_SEDIM_SPLI beg:PTHT",PRECISION)
+ CALL MPPDB_CHECK3D(PPRXS,"INTERNAL_SEDIM_SPLI beg:PPRXS",PRECISION)
+ IF (PRESENT(PRAY)) CALL MPPDB_CHECK3D(PRAY,"INTERNAL_SEDIM_SPLI beg:PRAY",PRECISION)
+ IF (PRESENT(PLBC)) CALL MPPDB_CHECK3D(PLBC,"INTERNAL_SEDIM_SPLI beg:PLBC",PRECISION)
+ IF (PRESENT(PFSEDC)) CALL MPPDB_CHECK3D(PFSEDC,"INTERNAL_SEDIM_SPLI beg:PFSEDC",PRECISION)
+ IF (PRESENT(PCONC3D)) CALL MPPDB_CHECK3D(PCONC3D,"INTERNAL_SEDIM_SPLI beg:PCONC3D",PRECISION)
+ !Check all INOUT arrays
+ CALL MPPDB_CHECK3D(PRXT,"INTERNAL_SEDIM_SPLI beg:PRXT",PRECISION)
+ CALL MPPDB_CHECK3D(PRXS,"INTERNAL_SEDIM_SPLI beg:PRXS",PRECISION)
+ IF (PRESENT(PFPR)) THEN
+ DO JI=1,KRR
+ CALL MPPDB_CHECK3D(PFPR(:,:,:,JI),"INTERNAL_SEDIM_SPLI beg:",PRECISION)
+ END DO
+ END IF
END IF
!
!-------------------------------------------------------------------------------
-IF (KSPE<2 .OR. KSPE>7) CALL PRINT_MSG(NVERB_FATAL,'GEN','ICE4_SEDIMENTATION_SPLIT','invalid species (KSPE variable)')
+IF (KSPE<2 .OR. KSPE>7) CALL PRINT_MSG(NVERB_FATAL,'GEN','INTERNAL_SEDIM_SPLIT','invalid species (KSPE variable)')
!
IF (PRESENT(PFPR)) THEN
GPRESENT_PFPR = .TRUE.
@@ -479,7 +496,7 @@ ZREMAINT(:,:) = PTSTEP
DO WHILE (ANY(ZREMAINT>0.))
!$acc kernels
ISEDIM = 0
-!$acc loop collapse(3) private(IDX,JI,JJ,JK) independent
+!acc loop collapse(3) private(IDX,JI,JJ,JK) independent
DO JK = KKTB,KKTE
DO JJ = KJB,KJE
DO JI = KIB,KIE
@@ -515,7 +532,7 @@ if (JK==-9999) print *,'PW: ISEDIM=',ISEDIM
!$acc kernels
!******* for cloud
ZWSED(:,:,:) = 0.
-!$acc loop independent private(JI,JJ,JK,JL,ZZWLBDC,ZRAY,ZZT,ZZWLBDA,ZZCC)
+!acc loop independent private(JI,JJ,JK,JL,ZZWLBDC,ZRAY,ZZT,ZZWLBDA,ZZCC)
DO JL=1, ISEDIM
JI=I1(JL)
JJ=I2(JL)
@@ -547,7 +564,7 @@ if (JK==-9999) print *,'PW: ISEDIM=',ISEDIM
!$acc kernels
! ******* for pristine ice
ZWSED(:,:,:) = 0.
-!$acc loop independent private(JI,JJ,JK,JL)
+!acc loop independent private(JI,JJ,JK,JL)
DO JL=1, ISEDIM
JI=I1(JL)
JJ=I2(JL)
@@ -586,7 +603,7 @@ if (JK==-9999) print *,'PW: ISEDIM=',ISEDIM
END SELECT
!
ZWSED(:,:,:) = 0.
-!$acc loop independent private(JI,JJ,JK,JL)
+!acc loop independent private(JI,JJ,JK,JL)
DO JL=1, ISEDIM
JI=I1(JL)
JJ=I2(JL)
@@ -633,16 +650,18 @@ if (JK==-9999) print *,'PW: ISEDIM=',ISEDIM
!
END DO
!
-!Check all INOUT arrays
-CALL MPPDB_CHECK3D(PRXT,"INTERNAL_SEDIM_SPLI end:PRXT",PRECISION)
-CALL MPPDB_CHECK3D(PRXS,"INTERNAL_SEDIM_SPLI end:PRXS",PRECISION)
-IF (PRESENT(PFPR)) THEN
- DO JI=1,KRR
- CALL MPPDB_CHECK3D(PFPR(:,:,:,JI),"INTERNAL_SEDIM_SPLI beg:",PRECISION)
- END DO
+IF (MPPDB_INITIALIZED) THEN
+ !Check all INOUT arrays
+ CALL MPPDB_CHECK3D(PRXT,"INTERNAL_SEDIM_SPLI end:PRXT",PRECISION)
+ CALL MPPDB_CHECK3D(PRXS,"INTERNAL_SEDIM_SPLI end:PRXS",PRECISION)
+ IF (PRESENT(PFPR)) THEN
+ DO JI=1,KRR
+ CALL MPPDB_CHECK3D(PFPR(:,:,:,JI),"INTERNAL_SEDIM_SPLI beg:",PRECISION)
+ END DO
+ END IF
+ !Check all OUT arrays
+ CALL MPPDB_CHECK2D(PINPRX,"INTERNAL_SEDIM_SPLI end:PINPRX",PRECISION)
END IF
-!Check all OUT arrays
-CALL MPPDB_CHECK2D(PINPRX,"INTERNAL_SEDIM_SPLI end:PINPRX",PRECISION)
!
END SUBROUTINE INTERNAL_SEDIM_SPLI
!
diff --git a/src/MNH/ice4_slow.f90 b/src/MNH/ice4_slow.f90
index 92713b432001e5194d1835323c15c1013a5346cd..c106a610795e1db9b4c35fa8c89c531df49a079d 100644
--- a/src/MNH/ice4_slow.f90
+++ b/src/MNH/ice4_slow.f90
@@ -67,9 +67,10 @@ SUBROUTINE ICE4_SLOW(KSIZE, LDSOFT, LDCOMPUTE, PRHODREF, PT, &
!* 0. DECLARATIONS
! ------------
!
-USE MODD_CST
-USE MODD_RAIN_ICE_PARAM
-USE MODD_RAIN_ICE_DESCR
+USE MODD_CST, ONLY: XTT
+USE MODD_RAIN_ICE_DESCR, ONLY: XCEXVT,XRTMIN
+USE MODD_RAIN_ICE_PARAM, ONLY: X0DEPG,X0DEPS,X1DEPG,X1DEPS,XACRIAUTI,XALPHA3,XBCRIAUTI,XBETA3,XCOLEXIS,XCRIAUTI, &
+ XEX0DEPG,XEX0DEPS,XEX1DEPG,XEX1DEPS,XEXIAGGS,XFIAGGS,XHON,XTEXAUTI,XTIMAUTI
!
#ifdef MNH_BITREP
USE MODI_BITREP
@@ -110,11 +111,25 @@ REAL, DIMENSION(KSIZE), INTENT(INOUT) :: PA_RI
REAL, DIMENSION(KSIZE), INTENT(INOUT) :: PA_RS
REAL, DIMENSION(KSIZE), INTENT(INOUT) :: PA_RG
!
+!$acc declare present(KSIZE,LDSOFT,LDCOMPUTE,PRHODREF,PT,PSSI,PLVFACT,PLSFACT, &
+!$acc& PRVT,PRCT,PRIT,PRST,PRGT,PLBDAS,PLBDAG,PAI,PCJ, &
+!$acc& PRCHONI,PRVDEPS,PRIAGGS,PRIAUTS,PRVDEPG, &
+!$acc& PA_TH,PA_RV,PA_RC,PA_RI,PA_RS,PA_RG )
+!
!* 0.2 declaration of local variables
!
-REAL, DIMENSION(SIZE(PRHODREF)) :: ZCRIAUTI
-REAL :: ZTIMAUTIC
+LOGICAL :: GDSOFT !Workaround of PGI bug with OpenACC (at least up to 18.10 version)
LOGICAL, DIMENSION(SIZE(PRHODREF)) :: GMASK
+REAL, DIMENSION(SIZE(PRHODREF)) :: ZCRIAUTI
+REAL :: ZTIMAUTIC
+!
+!$acc declare create(ZCRIAUTI,ZTIMAUTIC,GMASK)
+!
+!$acc declare copyin(XTT, &
+!$acc& XCEXVT,XRTMIN, &
+!$acc& X0DEPG,X0DEPS,X1DEPG,X1DEPS,XACRIAUTI,XALPHA3,XBCRIAUTI,XBETA3,XCOLEXIS,XCRIAUTI, &
+!$acc& XEX0DEPG,XEX0DEPS,XEX1DEPG,XEX1DEPS,XEXIAGGS,XFIAGGS,XHON,XTEXAUTI,XTIMAUTI )
+!
!-------------------------------------------------------------------------------
!
IF (MPPDB_INITIALIZED) THEN
@@ -148,10 +163,15 @@ IF (MPPDB_INITIALIZED) THEN
CALL MPPDB_CHECK1D(PA_RG,"ICE4_SLOW beg:PA_RG",PRECISION)
END IF
!
+!$acc kernels
+GDSOFT = LDSOFT !Workaround of PGI bug with OpenACC (at least up to 18.10 version)
+!$acc end kernels
+!
!* 3.2 compute the homogeneous nucleation source: RCHONI
!
+!$acc kernels
GMASK(:)=PT(:)<XTT-35.0 .AND. PRCT(:)>XRTMIN(2) .AND. LDCOMPUTE(:)
-IF(LDSOFT) THEN
+IF(GDSOFT) THEN
WHERE(.NOT. GMASK(:))
PRCHONI(:) = 0.
END WHERE
@@ -187,7 +207,7 @@ PA_TH(:) = PA_TH(:) + PRCHONI(:)*(PLSFACT(:)-PLVFACT(:))
!* 3.4.3 compute the deposition on r_s: RVDEPS
!
GMASK(:)=PRVT(:)>XRTMIN(1) .AND. PRST(:)>XRTMIN(5) .AND. LDCOMPUTE(:)
-IF(LDSOFT) THEN
+IF(GDSOFT) THEN
WHERE(.NOT. GMASK(:))
PRVDEPS(:) = 0.
END WHERE
@@ -209,7 +229,7 @@ PA_TH(:) = PA_TH(:) + PRVDEPS(:)*PLSFACT(:)
!* 3.4.4 compute the aggregation on r_s: RIAGGS
!
GMASK(:)=PRIT(:)>XRTMIN(4) .AND. PRST(:)>XRTMIN(5) .AND. LDCOMPUTE(:)
-IF(LDSOFT) THEN
+IF(GDSOFT) THEN
WHERE(.NOT. GMASK(:))
PRIAGGS(:) = 0.
END WHERE
@@ -235,7 +255,7 @@ PA_RI(:) = PA_RI(:) - PRIAGGS(:)
!* 3.4.5 compute the autoconversion of r_i for r_s production: RIAUTS
!
GMASK(:)=PRIT(:)>XRTMIN(4) .AND. LDCOMPUTE(:)
-IF(LDSOFT) THEN
+IF(GDSOFT) THEN
WHERE(.NOT. GMASK(:))
PRIAUTS(:) = 0.
END WHERE
@@ -263,7 +283,7 @@ PA_RI(:) = PA_RI(:) - PRIAUTS(:)
!
!
GMASK(:)=PRVT(:)>XRTMIN(1) .AND. PRGT(:)>XRTMIN(6) .AND. LDCOMPUTE(:)
-IF(LDSOFT) THEN
+IF(GDSOFT) THEN
WHERE(.NOT. GMASK(:))
PRVDEPG(:) = 0.
END WHERE
@@ -281,6 +301,7 @@ ENDIF
PA_RG(:) = PA_RG(:) + PRVDEPG(:)
PA_RV(:) = PA_RV(:) - PRVDEPG(:)
PA_TH(:) = PA_TH(:) + PRVDEPG(:)*PLSFACT(:)
+!$acc end kernels
!
IF (MPPDB_INITIALIZED) THEN
!Check all INOUT arrays
diff --git a/src/MNH/ice4_tendencies.f90 b/src/MNH/ice4_tendencies.f90
index 9a6a3f27ec58c79a065288c777113ce5bbe38090..915afbca2f7cbef6aebbebe5b62fcd4353477e6e 100644
--- a/src/MNH/ice4_tendencies.f90
+++ b/src/MNH/ice4_tendencies.f90
@@ -169,22 +169,22 @@ SUBROUTINE ICE4_TENDENCIES(KSIZE, KIB, KIE, KIT, KJB, KJE, KJT, KKB, KKE, KKT, K
!* 0. DECLARATIONS
! ------------
!
-USE MODD_CST
-USE MODD_RAIN_ICE_PARAM
-USE MODD_RAIN_ICE_DESCR
-
+USE MODD_CST, ONLY: XALPI,XBETAI,XCI,XCPV,XGAMI,XLSTT,XMD,XMV,XP00,XRV,XTT
+USE MODD_RAIN_ICE_DESCR, ONLY: XLBDAS_MAX,XLBEXG,XLBEXH,XLBEXR,XLBEXS,XLBG,XLBH,XLBR,XLBS,XRTMIN
+USE MODD_RAIN_ICE_PARAM, ONLY: XSCFAC
+!
+USE MODI_ICE4_COMPUTE_PDF
+USE MODI_ICE4_FAST_RG
+USE MODI_ICE4_FAST_RH
+USE MODI_ICE4_FAST_RI
+USE MODI_ICE4_FAST_RS
USE MODI_ICE4_NUCLEATION
-USE MODI_ICE4_RRHONG
+USE MODI_ICE4_RAINFR_VERT
USE MODI_ICE4_RIMLTC
+USE MODI_ICE4_RRHONG
USE MODI_ICE4_RSRIMCG_OLD
-USE MODI_ICE4_COMPUTE_PDF
-USE MODI_ICE4_RAINFR_VERT
USE MODI_ICE4_SLOW
USE MODI_ICE4_WARM
-USE MODI_ICE4_FAST_RS
-USE MODI_ICE4_FAST_RG
-USE MODI_ICE4_FAST_RH
-USE MODI_ICE4_FAST_RI
!
#ifdef MNH_BITREP
USE MODI_BITREP
@@ -299,18 +299,38 @@ REAL, DIMENSION(KSIZE), INTENT(OUT) :: PHLC_HRC
REAL, DIMENSION(KSIZE), INTENT(OUT) :: PHLC_LRC
REAL, DIMENSION(KIT,KJT,KKT), INTENT(OUT) :: PRAINFR ! Rain fraction
!
+!$acc declare present(KSIZE,KIB,KIE,KIT,KJB,KJE,KJT,KKB,KKE,KKT,KKL,KRR,LDSOFT, &
+!$acc& OWARM,HSUBG_RC_RR_ACCR,HSUBG_RR_EVAP,HSUBG_AUCV_RC,HSUBG_PR_PDF)
+!$acc declare present(LDCOMPUTE,PEXN,PRHODREF,PLVFACT,PLSFACT,LDMICRO,K1,K2,K3,PPRES,PCF,PCIT,PT,PTHT, &
+!$acc& PRVT,PRCT,PRRT,PRIT,PRST,PRGT,PRHT,PRRT3D,PSIGMA_RC,PRVHENI_MR,PRRHONG_MR,PRIMLTC_MR, &
+!$acc& PRSRIMCG_MR,PRCHONI,PRVDEPS,PRIAGGS,PRIAUTS,PRVDEPG,PRCAUTR,PRCACCR,PRREVAV, &
+!$acc& PRCRIMSS,PRCRIMSG,PRSRIMCG,PRRACCSS,PRRACCSG,PRSACCRG,PRSMLTG,PRCMLTSR,PRICFRRG,PRRCFRIG, &
+!$acc& PRICFRR,PRCWETG,PRIWETG,PRRWETG,PRSWETG,PRCDRYG,PRIDRYG,PRRDRYG,PRSDRYG, &
+!$acc& PRWETGH,PRWETGH_MR,PRGMLTR,PRCWETH,PRIWETH,PRSWETH,PRGWETH,PRRWETH, &
+!$acc& PRCDRYH,PRIDRYH,PRSDRYH,PRRDRYH,PRGDRYH,PRDRYHG,PRHMLTR,PRCBERI,PRS_TEND,PRG_TEND,PRH_TEND, &
+!$acc& PA_TH,PA_RV,PA_RC,PA_RR,PA_RI,PA_RS,PA_RG,PA_RH,PB_TH,PB_RV,PB_RC,PB_RR,PB_RI,PB_RS,PB_RG,PB_RH, &
+!$acc& PHLC_HCF,PHLC_LCF,PHLC_HRC,PHLC_LRC,PRAINFR)
+!
!* 0.2 declaration of local variables
!
-REAL, DIMENSION(KSIZE) :: ZRVT, ZRCT, ZRRT, ZRIT, ZRST, ZRGT, &
- & ZT, ZTHT, &
- & ZZW, &
- & ZSSI, ZKA, ZDV, ZAI, ZCJ, &
- & ZRF, &
- & ZLBDAR, ZLBDAS, ZLBDAG, ZLBDAH, ZLBDAR_RF, &
- & ZRGSI, ZRGSI_MR
+INTEGER :: IDX,JI,JJ,JK,JL
+INTEGER :: IRR, ISIZE
+LOGICAL, DIMENSION(KSIZE) :: LLWETG
+REAL, DIMENSION(KSIZE) :: ZRVT, ZRCT, ZRRT, ZRIT, ZRST, ZRGT, &
+ & ZT, ZTHT, &
+ & ZZW, &
+ & ZSSI, ZKA, ZDV, ZAI, ZCJ, &
+ & ZRF, &
+ & ZLBDAR, ZLBDAS, ZLBDAG, ZLBDAH, ZLBDAR_RF, &
+ & ZRGSI, ZRGSI_MR
REAL, DIMENSION(KIT,KJT,KKT) :: ZRRT3D
-INTEGER :: JL
-LOGICAL, DIMENSION(KSIZE) :: LLWETG
+!
+!$acc declare create(ZRVT,ZRCT,ZRRT,ZRIT,ZRST,ZRGT,ZT,ZTHT,ZZW,ZSSI,ZKA,ZDV,ZAI,ZCJ,ZRF, &
+!$acc& ZLBDAR,ZLBDAS,ZLBDAG,ZLBDAH,ZLBDAR_RF,ZRGSI,ZRGSI_MR,ZRRT3D,IDX,JI,JJ,JK,JL,LLWETG)
+!
+!$acc declare copyin(XALPI,XBETAI,XCI,XCPV,XGAMI,XLSTT,XMD,XMV,XP00,XRV,XTT, &
+!$acc& XSCFAC, &
+!$acc& XLBDAS_MAX,XLBEXG,XLBEXH,XLBEXR,XLBEXS,XLBG,XLBH,XLBR,XLBS,XRTMIN)
!
IF (MPPDB_INITIALIZED) THEN
!Check all IN arrays
@@ -358,6 +378,21 @@ IF (MPPDB_INITIALIZED) THEN
CALL MPPDB_CHECK2D(PRH_TEND,"ICE4_TENDENCIES beg:PRH_TEND",PRECISION)
END IF
!
+!$acc kernels
+IRR = KRR
+ISIZE = KSIZE
+!$acc end kernels
+#ifdef _OPENACC
+!Check if KRR and KSIZE have the same value on host and device
+IF (IRR/=KRR) THEN
+ CALL PRINT_MSG(NVERB_FATAL,'GEN','ICE4_TENDENCIES','OpenACC: KRR different between host and device')
+END IF
+IF (ISIZE/=KSIZE) THEN
+ CALL PRINT_MSG(NVERB_FATAL,'GEN','ICE4_TENDENCIES','OpenACC: KSIZE different between host and device')
+END IF
+#endif
+!
+!$acc kernels
PA_TH(:)=0.
PA_RV(:)=0.
PA_RC(:)=0.
@@ -383,6 +418,7 @@ ZRST(:)=PRST(:)
ZRGT(:)=PRGT(:)
ZTHT(:)=PTHT(:)
ZT(:)=PT(:)
+!$acc end kernels
!
!* 2. COMPUTES THE SLOW COLD PROCESS SOURCES
! --------------------------------------
@@ -390,10 +426,12 @@ CALL ICE4_NUCLEATION(KSIZE, LDSOFT, LDCOMPUTE, &
ZTHT, PPRES, PRHODREF, PEXN, PLSFACT, ZT, &
ZRVT, &
PCIT, PRVHENI_MR, PB_TH, PB_RV, PB_RI)
+!$acc kernels
ZRIT(:)=ZRIT(:) + PRVHENI_MR(:)
ZRVT(:)=ZRVT(:) - PRVHENI_MR(:)
ZTHT(:)=ZTHT(:) + PRVHENI_MR(:)*PLSFACT(:)
ZT(:) = ZTHT(:) * PEXN(:)
+!$acc end kernels
!
!* 3.3 compute the spontaneous freezing source: RRHONG
!
@@ -402,10 +440,12 @@ CALL ICE4_RRHONG(KSIZE, LDSOFT, LDCOMPUTE, &
&ZT, ZRRT, &
&ZTHT, &
&PRRHONG_MR, PB_TH, PB_RR, PB_RG)
+!$acc kernels
ZRGT(:) = ZRGT(:) + PRRHONG_MR(:)
ZRRT(:) = ZRRT(:) - PRRHONG_MR(:)
ZTHT(:) = ZTHT(:) + PRRHONG_MR(:)*(PLSFACT(:)-PLVFACT(:)) ! f(L_f*(RRHONG))
ZT(:) = ZTHT(:) * PEXN(:)
+!$acc end kernels
!
!* 7.1 cloud ice melting
!
@@ -414,6 +454,7 @@ CALL ICE4_RIMLTC(KSIZE, LDSOFT, LDCOMPUTE, &
&ZT, &
&ZTHT, ZRIT, &
&PRIMLTC_MR, PB_TH, PB_RC, PB_RI)
+!$acc kernels
ZRCT(:) = ZRCT(:) + PRIMLTC_MR(:)
ZRIT(:) = ZRIT(:) - PRIMLTC_MR(:)
ZTHT(:) = ZTHT(:) - PRIMLTC_MR(:)*(PLSFACT(:)-PLVFACT(:)) ! f(L_f*(-RIMLTC))
@@ -429,17 +470,21 @@ WHERE(ZRST(:)>0.)
ZLBDAS(:) = MIN(XLBDAS_MAX, XLBS*BR_POW(PRHODREF(:)*MAX(ZRST(:), XRTMIN(5)),XLBEXS))
#endif
END WHERE
+!$acc end kernels
CALL ICE4_RSRIMCG_OLD(KSIZE, LDSOFT, LDCOMPUTE, &
&PRHODREF, &
&ZLBDAS, &
&ZT, ZRCT, ZRST, &
&PRSRIMCG_MR, PB_RS, PB_RG)
+!$acc kernels
ZRST(:) = ZRST(:) - PRSRIMCG_MR(:)
ZRGT(:) = ZRGT(:) + PRSRIMCG_MR(:)
+!$acc end kernels
!
!* Derived fields
!
IF(KSIZE>0) THEN
+!$acc kernels
#ifndef MNH_BITREP
ZZW(:) = EXP(XALPI-XBETAI/ZT(:)-XGAMI*ALOG(ZT(:)))
#else
@@ -461,16 +506,37 @@ IF(KSIZE>0) THEN
ZCJ(JL) = XSCFAC*BR_POW(PRHODREF(JL),0.3) / SQRT(1.718E-5+0.0049E-5*(ZT(JL)-XTT))
#endif
ENDDO
+!$acc end kernels
!
!Cloud water split between high and low content part is done here
CALL ICE4_COMPUTE_PDF(KSIZE, HSUBG_AUCV_RC, HSUBG_PR_PDF,&
PRHODREF, ZRCT, PCF, PSIGMA_RC,&
PHLC_HCF, PHLC_LCF, PHLC_HRC, PHLC_LRC, ZRF)
!Diagnostic of precipitation fraction
- PRAINFR(:,:,:)=0.
- PRAINFR(:,:,:)=UNPACK(ZRF(:), MASK=LDMICRO(:,:,:), FIELD=PRAINFR(:,:,:))
+#ifndef _OPENACC
+ PRAINFR(:,:,:)=UNPACK(ZRF(:), MASK=LDMICRO(:,:,:), FIELD=0.)
ZRRT3D(:,:,:)=PRRT3D(:,:,:)-UNPACK(PRRHONG_MR(:), MASK=LDMICRO(:,:,:), FIELD=0.)
+#else
+!$acc kernels
+ IDX = 0
+ DO JK=1,KKT
+ DO JJ=1,KJT
+ DO JI=1,KIT
+ IF (LDMICRO(JI,JJ,JK)) THEN
+ IDX = IDX+1
+ PRAINFR(JI,JJ,JK) = ZRF(IDX)
+ ZRRT3D(JI,JJ,JK) = PRRT3D(JI,JJ,JK) - PRRHONG_MR(IDX)
+ ELSE
+ PRAINFR(JI,JJ,JK) = 0.
+ ZRRT3D(JI,JJ,JK) = PRRT3D(JI,JJ,JK)
+ END IF
+ END DO
+ END DO
+ END DO
+!$acc end kernels
+#endif
CALL ICE4_RAINFR_VERT(KIB, KIE, KIT, KJB, KJE, KJT, KKB, KKE, KKT, KKL, PRAINFR(:,:,:), ZRRT3D(:,:,:))
+!$acc kernels
DO JL=1,KSIZE
ZRF(JL)=PRAINFR(K1(JL), K2(JL), K3(JL))
END DO
@@ -511,7 +577,7 @@ IF(KSIZE>0) THEN
ZLBDAR_RF(:) = XLBR*BR_POW(PRHODREF(:) *MAX( ZRRT(:)/ZRF(:) , XRTMIN(3)),XLBEXR)
#endif
END WHERE
- IF(KRR==7) THEN
+ IF(IRR==7) THEN
ZLBDAH(:)=0.
WHERE(PRHT(:)>0.)
#ifndef MNH_BITREP
@@ -521,6 +587,7 @@ IF(KSIZE>0) THEN
#endif
END WHERE
ENDIF
+!$acc end kernels
ENDIF
!
!
@@ -550,9 +617,11 @@ IF(OWARM) THEN ! Check if the formation of the raindrops by the slow
&PRCAUTR, PRCACCR, PRREVAV, &
&PA_TH, PA_RV, PA_RC, PA_RR)
ELSE
+!$acc kernels
PRCAUTR(:)=0.
PRCACCR(:)=0.
PRREVAV(:)=0.
+!$acc end kernels
END IF
!
!-------------------------------------------------------------------------------
@@ -579,8 +648,10 @@ CALL ICE4_FAST_RS(KSIZE, LDSOFT, LDCOMPUTE, &
!* 5. COMPUTES THE FAST COLD PROCESS SOURCES FOR r_g
! ------------------------------------------------------
!
+!$acc kernels
ZRGSI(:) = PRVDEPG(:) + PRSMLTG(:) + PRRACCSG(:) + PRSACCRG(:) + PRCRIMSG(:) + PRSRIMCG(:)
ZRGSI_MR(:) = PRRHONG_MR(:) + PRSRIMCG_MR(:)
+!$acc end kernels
CALL ICE4_FAST_RG(KSIZE, LDSOFT, LDCOMPUTE, KRR, &
&PRHODREF, PLVFACT, PLSFACT, PPRES, &
&ZDV, ZKA, ZCJ, PCIT, &
@@ -610,6 +681,7 @@ IF (KRR==7) THEN
&PRH_TEND, &
&PA_TH, PA_RC, PA_RR, PA_RI, PA_RS, PA_RG, PA_RH)
ELSE
+!$acc kernels
PRCWETH(:)=0.
PRIWETH(:)=0.
PRSWETH(:)=0.
@@ -622,6 +694,7 @@ ELSE
PRGDRYH(:)=0.
PRDRYHG(:)=0.
PRHMLTR(:)=0.
+!$acc end kernels
END IF
!
!-------------------------------------------------------------------------------
diff --git a/src/MNH/ice4_warm.f90 b/src/MNH/ice4_warm.f90
index 659a44d880caba204b8825049d995d4f17ffe7e3..a0e5d366333cb1fb1733e03cc32b1773618b0e24 100644
--- a/src/MNH/ice4_warm.f90
+++ b/src/MNH/ice4_warm.f90
@@ -5,13 +5,13 @@
MODULE MODI_ICE4_WARM
INTERFACE
SUBROUTINE ICE4_WARM(KSIZE, LDSOFT, LDCOMPUTE, HSUBG_RC_RR_ACCR, HSUBG_RR_EVAP, &
- &PRHODREF, PLVFACT, PT, PPRES, PTHT, &
- &PLBDAR, PLBDAR_RF, PKA, PDV, PCJ, &
- &PHLC_LCF, PHLC_HCF, PHLC_LRC, PHLC_HRC, &
- &PCF, PRF, &
- &PRVT, PRCT, PRRT, &
- &PRCAUTR, PRCACCR, PRREVAV, &
- &PA_TH, PA_RV, PA_RC, PA_RR)
+ PRHODREF, PLVFACT, PT, PPRES, PTHT, &
+ PLBDAR, PLBDAR_RF, PKA, PDV, PCJ, &
+ PHLC_LCF, PHLC_HCF, PHLC_LRC, PHLC_HRC, &
+ PCF, PRF, &
+ PRVT, PRCT, PRRT, &
+ PRCAUTR, PRCACCR, PRREVAV, &
+ PA_TH, PA_RV, PA_RC, PA_RR)
IMPLICIT NONE
INTEGER, INTENT(IN) :: KSIZE
LOGICAL, INTENT(IN) :: LDSOFT
@@ -48,13 +48,13 @@ END SUBROUTINE ICE4_WARM
END INTERFACE
END MODULE MODI_ICE4_WARM
SUBROUTINE ICE4_WARM(KSIZE, LDSOFT, LDCOMPUTE, HSUBG_RC_RR_ACCR, HSUBG_RR_EVAP, &
- &PRHODREF, PLVFACT, PT, PPRES, PTHT, &
- &PLBDAR, PLBDAR_RF, PKA, PDV, PCJ, &
- &PHLC_LCF, PHLC_HCF, PHLC_LRC, PHLC_HRC, &
- &PCF, PRF, &
- &PRVT, PRCT, PRRT, &
- &PRCAUTR, PRCACCR, PRREVAV, &
- &PA_TH, PA_RV, PA_RC, PA_RR)
+ PRHODREF, PLVFACT, PT, PPRES, PTHT, &
+ PLBDAR, PLBDAR_RF, PKA, PDV, PCJ, &
+ PHLC_LCF, PHLC_HCF, PHLC_LRC, PHLC_HRC, &
+ PCF, PRF, &
+ PRVT, PRCT, PRRT, &
+ PRCAUTR, PRCACCR, PRREVAV, &
+ PA_TH, PA_RV, PA_RC, PA_RR)
!!
!!** PURPOSE
!! -------
@@ -72,9 +72,10 @@ SUBROUTINE ICE4_WARM(KSIZE, LDSOFT, LDCOMPUTE, HSUBG_RC_RR_ACCR, HSUBG_RR_EVAP,
!* 0. DECLARATIONS
! ------------
!
-USE MODD_CST
-USE MODD_RAIN_ICE_PARAM
-USE MODD_RAIN_ICE_DESCR
+USE MODD_CST, ONLY: XALPW,XBETAW,XCL,XCPD,XCPV,XGAMW,XLVTT,XMD,XMV,XRV,XTT
+USE MODD_RAIN_ICE_DESCR, ONLY: XCEXVT,XRTMIN
+USE MODD_RAIN_ICE_PARAM, ONLY: X0EVAR,X1EVAR,XCRIAUTC,XEX0EVAR,XEX1EVAR,XEXCACCR,XFCACCR,XTIMAUTC
+!
USE MODE_MSG
!
#ifdef MNH_BITREP
@@ -119,13 +120,20 @@ REAL, DIMENSION(KSIZE), INTENT(INOUT) :: PA_RV
REAL, DIMENSION(KSIZE), INTENT(INOUT) :: PA_RC
REAL, DIMENSION(KSIZE), INTENT(INOUT) :: PA_RR
!
+!$acc declare present(KSIZE,LDSOFT,LDCOMPUTE,HSUBG_RC_RR_ACCR,HSUBG_RR_EVAP,PRHODREF,PLVFACT,PT,PPRES,PTHT, &
+!$acc& PLBDAR,PLBDAR_RF,PKA,PDV,PCJ,PHLC_LCF,PHLC_HCF,PHLC_LRC,PHLC_HRC,PCF,PRF, &
+!$acc& PRVT,PRCT,PRRT,PRCAUTR,PRCACCR,PRREVAV,PA_TH,PA_RV,PA_RC,PA_RR )
+!
!* 0.2 declaration of local variables
!
-REAL, DIMENSION(SIZE(PRHODREF)) :: ZZW2, ZZW3, ZZW4
-REAL, DIMENSION(SIZE(PRHODREF)) :: ZUSW ! Undersaturation over water
-REAL, DIMENSION(SIZE(PRHODREF)) :: ZTHLT ! Liquid potential temperature
-REAL :: ZTIMAUTIC
+LOGICAL :: GDSOFT !Workaround of PGI bug with OpenACC (at least up to 18.10 version)
LOGICAL, DIMENSION(SIZE(PRHODREF)) :: GMASK, GMASK1, GMASK2
+REAL, DIMENSION(SIZE(PRHODREF)) :: ZZW2, ZZW3, ZZW4
+REAL, DIMENSION(SIZE(PRHODREF)) :: ZUSW ! Undersaturation over water
+REAL, DIMENSION(SIZE(PRHODREF)) :: ZTHLT ! Liquid potential temperature
+!
+!$acc declare create(ZZW2,ZZW3,ZZW4,ZUSW,ZTHLT,GMASK,GMASK1,GMASK2)
+!
!-------------------------------------------------------------------------------
!
IF (MPPDB_INITIALIZED) THEN
@@ -160,10 +168,15 @@ IF (MPPDB_INITIALIZED) THEN
CALL MPPDB_CHECK1D(PA_RR,"ICE4_WARM beg:PA_RR",PRECISION)
END IF
!
+!$acc kernels
+GDSOFT = LDSOFT !Workaround of PGI bug with OpenACC (at least up to 18.10 version)
+!$acc end kernels
+!
!* 4.2 compute the autoconversion of r_c for r_r production: RCAUTR
!
+!$acc kernels
GMASK(:)=PHLC_HRC(:)>XRTMIN(2) .AND. PHLC_HCF(:) .GT. 0. .AND. LDCOMPUTE(:)
-IF(LDSOFT) THEN
+IF(GDSOFT) THEN
WHERE(.NOT. GMASK(:))
PRCAUTR(:) = 0.
END WHERE
@@ -177,13 +190,12 @@ ENDIF
PA_RC(:) = PA_RC(:) - PRCAUTR(:)
PA_RR(:) = PA_RR(:) + PRCAUTR(:)
!
-!
!* 4.3 compute the accretion of r_c for r_r production: RCACCR
!
IF (HSUBG_RC_RR_ACCR=='NONE') THEN
!CLoud water and rain are diluted over the grid box
GMASK(:)=PRCT(:)>XRTMIN(2) .AND. PRRT(:)>XRTMIN(3) .AND. LDCOMPUTE(:)
- IF(LDSOFT) THEN
+ IF(GDSOFT) THEN
WHERE(.NOT. GMASK(:))
PRCACCR(:)=0.
END WHERE
@@ -213,7 +225,7 @@ ELSEIF (HSUBG_RC_RR_ACCR=='PRFR') THEN
GMASK(:)=PRCT(:)>XRTMIN(2) .AND. PRRT(:)>XRTMIN(3) .AND. LDCOMPUTE(:)
GMASK1(:)=GMASK(:) .AND. PHLC_HRC(:)>XRTMIN(2) .AND. PHLC_HCF(:)>0.
GMASK2(:)=GMASK(:) .AND. PHLC_LRC(:)>XRTMIN(2) .AND. PHLC_LCF(:)>0.
- IF(LDSOFT) THEN
+ IF(GDSOFT) THEN
WHERE(.NOT. (GMASK1(:) .OR. GMASK2(:)))
PRCACCR(:)=0.
END WHERE
@@ -251,6 +263,7 @@ ELSEIF (HSUBG_RC_RR_ACCR=='PRFR') THEN
ELSE
CALL PRINT_MSG(NVERB_FATAL,'GEN','ICE4_WARM','wrong HSUBG_RC_RR_ACCR case')
ENDIF
+!
PA_RC(:) = PA_RC(:) - PRCACCR(:)
PA_RR(:) = PA_RR(:) + PRCACCR(:)
!
@@ -258,7 +271,7 @@ PA_RR(:) = PA_RR(:) + PRCACCR(:)
!
IF (HSUBG_RR_EVAP=='NONE') THEN
GMASK(:)=PRRT(:)>XRTMIN(3) .AND. PRCT(:)<=XRTMIN(2) .AND. LDCOMPUTE(:)
- IF(LDSOFT) THEN
+ IF(GDSOFT) THEN
WHERE(.NOT. GMASK(:))
PRREVAV(:)=0.
END WHERE
@@ -305,7 +318,7 @@ ELSEIF (HSUBG_RR_EVAP=='CLFR' .OR. HSUBG_RR_EVAP=='PRFR') THEN
!et plusieurs versions (comme actuellement, en ciel clair, en ciel nuageux) de PKA, PDV, PCJ dans rain_ice
!On utiliserait la bonne version suivant l'option NONE, CLFR... dans l'évaporation et ailleurs
GMASK(:)=PRRT(:)>XRTMIN(3) .AND. ZZW4(:) > PCF(:) .AND. LDCOMPUTE(:)
- IF(LDSOFT) THEN
+ IF(GDSOFT) THEN
WHERE(.NOT. GMASK(:))
PRREVAV(:)=0.
END WHERE
@@ -355,6 +368,7 @@ END IF
PA_RR(:) = PA_RR(:) - PRREVAV(:)
PA_RV(:) = PA_RV(:) + PRREVAV(:)
PA_TH(:) = PA_TH(:) - PRREVAV(:)*PLVFACT(:)
+!$acc end kernels
!
IF (MPPDB_INITIALIZED) THEN
!Check all INOUT arrays
diff --git a/src/MNH/rain_ice_red.f90 b/src/MNH/rain_ice_red.f90
index 3dca1dbb515d038cb474497e39141bbac3ff52e7..22037c94acd9664aae4c1eb7f3e0496fbb5de93a 100644
--- a/src/MNH/rain_ice_red.f90
+++ b/src/MNH/rain_ice_red.f90
@@ -241,19 +241,20 @@ END MODULE MODI_RAIN_ICE_RED
!* 0. DECLARATIONS
! ------------
!
-USE MODD_PARAMETERS
-USE MODD_CST
-USE MODD_RAIN_ICE_DESCR
-USE MODD_RAIN_ICE_PARAM
-USE MODD_PARAM_ICE
-USE MODD_BUDGET
-USE MODD_LES
+USE MODD_PARAMETERS, ONLY: JPVEXT,XUNDEF
+USE MODD_CST, ONLY: XCI,XCL,XCPD,XCPV,XLSTT,XLVTT,XTT
+USE MODD_RAIN_ICE_DESCR, ONLY: XRTMIN
+USE MODD_PARAM_ICE, ONLY: CSUBG_PR_PDF,CSUBG_RC_RR_ACCR,CSUBG_RR_EVAP,LDEPOSC,LFEEDBACKT,LSEDIM_AFTER, &
+ NMAXITER,XMRSTEP,XTSTEP_TS,XVDEPOSC
+USE MODD_BUDGET, ONLY: LBU_ENABLE,LBUDGET_RC,LBUDGET_RR,LBUDGET_RI,LBUDGET_RS,LBUDGET_RG,LBUDGET_RH,LBUDGET_RV,LBUDGET_TH
+!
USE MODI_BUDGET
+USE MODI_ICE4_NUCLEATION_WRAPPER
USE MODI_ICE4_RAINFR_VERT
USE MODI_ICE4_SEDIMENTATION_STAT
USE MODI_ICE4_SEDIMENTATION_SPLIT
-USE MODI_ICE4_NUCLEATION_WRAPPER
USE MODI_ICE4_TENDENCIES
+!
USE MODE_FMWRIT
USE MODE_ll
USE MODE_MSG
@@ -322,6 +323,7 @@ REAL, DIMENSION(:,:,:), OPTIONAL, INTENT(INOUT) :: PRHS ! Hail m.r. source
REAL, DIMENSION(:,:), OPTIONAL, INTENT(OUT) :: PINPRH! Hail instant precip
REAL, DIMENSION(:,:,:,:), OPTIONAL, INTENT(OUT) :: PFPR ! upper-air precipitation fluxes
!
+!$acc declare present(OSEDIC,HSEDIM,HSUBG_AUCV_RC,OWARM,KKA,KKU,KKL,PTSTEP,KRR)
!$acc declare present(LDMICRO,PEXN,PDZZ,PRHODJ,PRHODREF,PEXNREF,PPABST,PCIT,PCLDFR,PTHT,PRVT, &
!$acc & PRCT,PRRT,PRIT,PRST,PRGT,PSIGS,PTHS,PRVS,PRCS,PRRS,PRIS,PRSS,PRGS, &
!$acc & PINPRC,PINDEP,PINPRR,PINPRR3D,PEVAP3D,PINPRS,PINPRG) &
@@ -343,13 +345,15 @@ INTEGER :: IMICRO ! Case r_x>0 locations
INTEGER, DIMENSION(COUNT(LDMICRO)) :: I1,I2,I3 ! Used to replace the COUNT
INTEGER :: JL ! and PACK intrinsics
!
+!$acc declare create(IIB,IIE,IIT,IJB,IJE,IJT,IKB,IKTB,IKT,IKE,IKTE,IMICRO,I1,I2,I3,JL)
+!
!Arrays for nucleation call outisde of LDMICRO points
REAL, DIMENSION(SIZE(PEXNREF,1),SIZE(PEXNREF,2),SIZE(PEXNREF,3)) :: ZW ! work array
REAL, DIMENSION(SIZE(PEXNREF,1),SIZE(PEXNREF,2),SIZE(PEXNREF,3)) :: ZT ! Temperature
REAL, DIMENSION(SIZE(PTHT,1),SIZE(PTHT,2),SIZE(PTHT,3)) :: &
& ZZ_RVHENI_MR, & ! heterogeneous nucleation mixing ratio change
& ZZ_RVHENI ! heterogeneous nucleation
-REAL, DIMENSION(SIZE(PTHT,1),SIZE(PTHT,2),SIZE(PTHT,3)) :: ZZ_LVFACT, ZZ_LSFACT
+REAL, DIMENSION(SIZE(PTHT,1),SIZE(PTHT,2),SIZE(PTHT,3)) :: ZZ_LVFACT, ZZ_LSFACT, ZLSFACT3D
!
!Diagnostics
REAL, DIMENSION(SIZE(PTHT,1),SIZE(PTHT,2),SIZE(PTHT,3)) :: ZRAINFR, &
@@ -358,7 +362,9 @@ REAL, DIMENSION(SIZE(PTHT,1),SIZE(PTHT,2),SIZE(PTHT,3)) :: ZRAINFR, &
& ZHLC_HRC3D,& ! HLCLOUDS cloud water content in high water content
& ZHLC_LRC3D ! HLCLOUDS cloud water content in low water content
REAL, DIMENSION(SIZE(PTHT,1),SIZE(PTHT,2)) :: ZINPRI ! Pristine ice instant precip
-!$acc declare create(ZINPRI)
+!
+!$acc declare create(ZW,ZT,ZZ_RVHENI_MR,ZZ_RVHENI,ZZ_LVFACT,ZZ_LSFACT,ZLSFACT3D, &
+!$acc & ZRAINFR,ZHLC_HCF3D,ZHLC_LCF3D,ZHLC_HRC3D,ZHLC_LRC3D,ZINPRI)
!
!Packed variables
REAL, DIMENSION(COUNT(LDMICRO)) :: ZRVT, & ! Water vapor m.r. at t
@@ -385,6 +391,9 @@ REAL, DIMENSION(COUNT(LDMICRO)) :: ZRVT, & ! Water vapor m.r. at t
& ZHLC_LRC ! HLCLOUDS : LWC that is Low LWC in grid
! note that ZRC = ZHLC_HRC + ZHLC_LRC
!
+!$acc declare create(ZRVT,ZRCT,ZRRT,ZRIT,ZRST,ZRGT,ZRHT,ZCIT,ZTHT,ZRHODREF,ZZT,ZPRES,ZEXN, &
+!$acc& ZLSFACT,ZLVFACT,ZSIGMA_RC,ZCF,ZHLC_HCF,ZHLC_LCF,ZHLC_HRC,ZHLC_LRC)
+!
!Output packed tendencies (for budgets only)
REAL, DIMENSION(COUNT(LDMICRO)) :: ZRVHENI_MR, & ! heterogeneous nucleation mixing ratio change
& ZRCHONI, & ! Homogeneous nucleation
@@ -413,6 +422,11 @@ REAL, DIMENSION(COUNT(LDMICRO)) :: ZRVHENI_MR, & ! heterogeneous nucleation mixi
& ZRCDRYH, ZRIDRYH, ZRSDRYH, ZRRDRYH, ZRGDRYH, & ! Wet growth of hailstone
& ZRDRYHG ! Conversion of hailstone into graupel
!
+!$acc declare create(ZRVHENI_MR,ZRCHONI,ZRRHONG_MR,ZRVDEPS,ZRIAGGS,ZRIAUTS,ZRVDEPG,ZRCAUTR,ZRCACCR,ZRREVAV,ZRIMLTC_MR, &
+!$acc& ZRCBERI,ZRHMLTR,ZRSMLTG,ZRCMLTSR,ZRRACCSS,ZRRACCSG,ZRSACCRG,ZRCRIMSS,ZRCRIMSG,ZRSRIMCG,ZRSRIMCG_MR, &
+!$acc& ZRICFRRG, ZRRCFRIG, ZRICFRR,ZRCWETG,ZRIWETG,ZRRWETG,ZRSWETG,ZRCDRYG,ZRIDRYG,ZRRDRYG,ZRSDRYG,ZRWETGH, &
+!$acc& ZRWETGH_MR,ZRGMLTR,ZRCWETH,ZRIWETH,ZRSWETH,ZRGWETH,ZRRWETH,ZRCDRYH,ZRIDRYH,ZRSDRYH,ZRRDRYH,ZRGDRYH,ZRDRYHG)
+!
!Output packed total mixing ratio change (for budgets only)
REAL, DIMENSION(COUNT(LDMICRO)) :: ZTOT_RVHENI, & ! heterogeneous nucleation mixing ratio change
& ZTOT_RCHONI, & ! Homogeneous nucleation
@@ -440,6 +454,13 @@ REAL, DIMENSION(COUNT(LDMICRO)) :: ZTOT_RVHENI, & ! heterogeneous nucleation mix
& ZTOT_RCDRYH, ZTOT_RIDRYH, ZTOT_RSDRYH, ZTOT_RRDRYH, ZTOT_RGDRYH, & ! Wet growth of hailstone
& ZTOT_RDRYHG ! Conversion of hailstone into graupel
!
+!$acc declare create(ZTOT_RVHENI,ZTOT_RCHONI,ZTOT_RRHONG,ZTOT_RVDEPS,ZTOT_RIAGGS,ZTOT_RIAUTS,ZTOT_RVDEPG,ZTOT_RCAUTR, &
+!$acc& ZTOT_RCACCR,ZTOT_RREVAV,ZTOT_RCRIMSS,ZTOT_RCRIMSG,ZTOT_RSRIMCG,ZTOT_RIMLTC,ZTOT_RCBERI,ZTOT_RHMLTR, &
+!$acc& ZTOT_RSMLTG,ZTOT_RCMLTSR,ZTOT_RRACCSS, ZTOT_RRACCSG, ZTOT_RSACCRG,ZTOT_RICFRRG, ZTOT_RRCFRIG, &
+!$acc& ZTOT_RICFRR,ZTOT_RCWETG,ZTOT_RIWETG,ZTOT_RRWETG,ZTOT_RSWETG,ZTOT_RCDRYG,ZTOT_RIDRYG,ZTOT_RRDRYG, &
+!$acc& ZTOT_RSDRYG,ZTOT_RWETGH,ZTOT_RGMLTR,ZTOT_RCWETH,ZTOT_RIWETH,ZTOT_RSWETH,ZTOT_RGWETH,ZTOT_RRWETH, &
+!$acc& ZTOT_RCDRYH,ZTOT_RIDRYH,ZTOT_RSDRYH,ZTOT_RRDRYH,ZTOT_RGDRYH,ZTOT_RDRYHG)
+!
!For time- or mixing-ratio- splitting
REAL, DIMENSION(COUNT(LDMICRO)) :: Z0RVT, & ! Water vapor m.r. at the beginig of the current loop
& Z0RCT, & ! Cloud water m.r. at the beginig of the current loop
@@ -451,6 +472,10 @@ REAL, DIMENSION(COUNT(LDMICRO)) :: Z0RVT, & ! Water vapor m.r. at the begi
& ZA_TH, ZA_RV, ZA_RC, ZA_RR, ZA_RI, ZA_RS, ZA_RG, ZA_RH, &
& ZB_TH, ZB_RV, ZB_RC, ZB_RR, ZB_RI, ZB_RS, ZB_RG, ZB_RH
!
+!$acc declare create(Z0RVT,Z0RCT,Z0RRT,Z0RIT,Z0RST,Z0RGT,Z0RHT, &
+!$acc& ZA_TH,ZA_RV,ZA_RC,ZA_RR,ZA_RI,ZA_RS,ZA_RG,ZA_RH, &
+!$acc& ZB_TH,ZB_RV,ZB_RC,ZB_RR,ZB_RI,ZB_RS,ZB_RG,ZB_RH )
+!
!To take into acount external tendencies inside the splitting
REAL, DIMENSION(COUNT(LDMICRO)) :: ZEXT_RV, & ! External tendencie for rv
ZEXT_RC, & ! External tendencie for rc
@@ -459,10 +484,11 @@ REAL, DIMENSION(COUNT(LDMICRO)) :: ZEXT_RV, & ! External tendencie for rv
ZEXT_RS, & ! External tendencie for rs
ZEXT_RG, & ! External tendencie for rg
ZEXT_RH, & ! External tendencie for rh
- ZEXT_TH, & ! External tendencie for th
- ZEXT_WW ! Working array
+ ZEXT_TH ! External tendencie for th
LOGICAL :: LEXT_TEND
!
+!$acc declare create(ZEXT_RV,ZEXT_RC,ZEXT_RR,ZEXT_RI,ZEXT_RS,ZEXT_RG,ZEXT_RH,ZEXT_TH,LEXT_TEND)
+!
INTEGER, DIMENSION(COUNT(LDMICRO)) :: IITER ! Number of iterations done (with real tendencies computation)
INTEGER :: INB_ITER_MAX ! Maximum number of iterations (with real tendencies computation)
REAL, DIMENSION(COUNT(LDMICRO)) :: ZTIME, & ! Current integration time (starts with 0 and ends with PTSTEP)
@@ -472,6 +498,7 @@ REAL, DIMENSION(COUNT(LDMICRO)) :: ZTIME, & ! Current integration time (start
LOGICAL, DIMENSION(COUNT(LDMICRO)) :: LLCOMPUTE ! Points where we must compute tendenceis
LOGICAL :: LSOFT ! Must we really compute tendencies or only adjust them to new T variables
LOGICAL, DIMENSION(SIZE(PRHODREF,1),SIZE(PRHODREF,2)):: GDEP
+LOGICAL, DIMENSION(SIZE(LDMICRO,1),SIZE(LDMICRO,2),SIZE(LDMICRO,3)) :: GDNOTMICRO ! = .NOT.LDMICRO
REAL :: ZTSTEP ! length of sub-timestep in case of time splitting
REAL :: ZINV_TSTEP ! Inverse ov PTSTEP
REAL, DIMENSION(COUNT(LDMICRO), 6) :: ZRS_TEND
@@ -482,8 +509,16 @@ REAL, DIMENSION(COUNT(LDMICRO), 8) :: ZRH_TEND
REAL, DIMENSION(SIZE(PTHT,1),SIZE(PTHT,2),SIZE(PTHT,3)) :: &
&ZW_RVS, ZW_RCS, ZW_RRS, ZW_RIS, ZW_RSS, ZW_RGS, ZW_RHS, ZW_THS
!
-!PW: to remove later when everything on GPU
-!$acc update self(PDZZ,PCLDFR,PTHS,PRVS,PRCS,PRIS)
+!$acc declare create(IITER,INB_ITER_MAX,ZTIME,ZMAXTIME,ZTIME_THRESHOLD,ZTIME_LASTCALL,LLCOMPUTE,LSOFT,GDEP,GDNOTMICRO,ZTSTEP, &
+!$acc & ZINV_TSTEP,ZRS_TEND,ZRG_TEND,ZRH_TEND,ZW_RVS,ZW_RCS,ZW_RRS,ZW_RIS,ZW_RSS,ZW_RGS,ZW_RHS,ZW_THS)
+!
+!acc declare copyin(JPVEXT,XUNDEF) !These are parameters
+!$acc declare copyin(XCI,XCL,XCPD,XCPV,XLSTT,XLVTT,XTT, &
+!$acc& XRTMIN, &
+!$acc& CSUBG_PR_PDF,CSUBG_RC_RR_ACCR,CSUBG_RR_EVAP,LDEPOSC,LFEEDBACKT,LSEDIM_AFTER, &
+!$acc& NMAXITER,XMRSTEP,XTSTEP_TS,XVDEPOSC, &
+!$acc& LBU_ENABLE,LBUDGET_RC,LBUDGET_RR,LBUDGET_RI,LBUDGET_RS,LBUDGET_RG,LBUDGET_RH,LBUDGET_RV,LBUDGET_TH)
+!
IF (MPPDB_INITIALIZED) THEN
!Check all IN arrays
CALL MPPDB_CHECK3D_LOG(LDMICRO,"RAIN_ICE_RED beg:LDMICRO")
@@ -526,6 +561,8 @@ END IF
! -----------------------
!
CALL GET_INDICE_ll (IIB,IJB,IIE,IJE)
+!$acc update device(IIB,IJB,IIE,IJE)
+!$acc kernels
IIT=SIZE(PDZZ,1)
IJT=SIZE(PDZZ,2)
IKB=KKA+JPVEXT*KKL
@@ -533,29 +570,60 @@ IKE=KKU-JPVEXT*KKL
IKT=SIZE(PDZZ,3)
IKTB=1+JPVEXT
IKTE=IKT-JPVEXT
+!$acc end kernels
+!$acc update self(IIT,IJT,IKB,IKE,IKT,IKTB,IKTE)
!
+!$acc kernels
ZINV_TSTEP=1./PTSTEP
LEXT_TEND=.TRUE.
!
+!Not necessary (done in ICE4_TENDENCIES when LSOFT=.FALSE.)
+!but useful for calls to MPPDB_CHECK
+ZRS_TEND(:,:) = 0.
+ZRG_TEND(:,:) = 0.
+ZRH_TEND(:,:) = 0.
+ZRCHONI(:) = 0.
+ZRVDEPS(:) = 0.
+ZRIAGGS(:) = 0.
+ZRIAUTS(:) = 0.
+ZRVDEPG(:) = 0.
+ZRCAUTR(:) = 0.
+ZRCACCR(:) = 0.
+ZRREVAV(:) = 0.
+ZRSMLTG(:) = 0.
+ZRCMLTSR(:)= 0.
+ZRICFRRG(:) = 0.
+ZRRCFRIG(:) = 0.
+ZRICFRR(:) = 0.
+ZRGMLTR(:) = 0.
+ZRHMLTR(:) = 0.
+ZRCBERI(:) = 0.
+!$acc end kernels
+!$acc update self(LEXT_TEND)
+!
!$acc kernels
ZT(:,:,:) = PTHT(:,:,:) * PEXN(:,:,:)
! LSFACT and LVFACT without exner
+!$acc end kernels
IF(KRR==7) THEN
+!$acc kernels
ZZ_LSFACT(:,:,:)=(XLSTT+(XCPV-XCI)*(ZT(:,:,:)-XTT)) &
/( XCPD + XCPV*PRVT(:,:,:) + XCL*(PRCT(:,:,:)+PRRT(:,:,:)) &
+ XCI*(PRIT(:,:,:)+PRST(:,:,:)+PRGT(:,:,:)+PRHT(:,:,:)))
ZZ_LVFACT(:,:,:)=(XLVTT+(XCPV-XCL)*(ZT(:,:,:)-XTT)) &
/( XCPD + XCPV*PRVT(:,:,:) + XCL*(PRCT(:,:,:)+PRRT(:,:,:)) &
+ XCI*(PRIT(:,:,:)+PRST(:,:,:)+PRGT(:,:,:)+PRHT(:,:,:)))
+!$acc end kernels
ELSE
+!$acc kernels
ZZ_LSFACT(:,:,:)=(XLSTT+(XCPV-XCI)*(ZT(:,:,:)-XTT)) &
/( XCPD + XCPV*PRVT(:,:,:) + XCL*(PRCT(:,:,:)+PRRT(:,:,:)) &
+ XCI*(PRIT(:,:,:)+PRST(:,:,:)+PRGT(:,:,:)))
ZZ_LVFACT(:,:,:)=(XLVTT+(XCPV-XCL)*(ZT(:,:,:)-XTT)) &
/( XCPD + XCPV*PRVT(:,:,:) + XCL*(PRCT(:,:,:)+PRRT(:,:,:)) &
+ XCI*(PRIT(:,:,:)+PRST(:,:,:)+PRGT(:,:,:)))
-ENDIF
!$acc end kernels
+ENDIF
!
!-------------------------------------------------------------------------------
!
@@ -568,8 +636,7 @@ IF(.NOT. LSEDIM_AFTER) THEN
!
IF(HSEDIM=='STAT') THEN
#ifdef _OPENACC
-PRINT *,'OPENACC: RAIN_ICE_RED::HSEDIM=STAT not yet implemented'
-CALL ABORT
+ CALL PRINT_MSG(NVERB_FATAL,'GEN','RAIN_ICE_RED','OpenACC: HSEDIM=STAT not yet implemented')
#endif
!SR: It *seems* that we must have two separate calls for ifort
IF(KRR==7) THEN
@@ -597,7 +664,6 @@ CALL ABORT
ELSEIF(HSEDIM=='SPLI') THEN
#ifdef _OPENACC
PRINT *,'OPENACC: RAIN_ICE_RED::HSEDIM=SPLI being implemented'
-! CALL ABORT
#endif
!SR: It *seems* that we must have two separate calls for ifort
IF(KRR==7) THEN
@@ -617,8 +683,9 @@ PRINT *,'OPENACC: RAIN_ICE_RED::HSEDIM=SPLI being implemented'
&PSEA, PTOWN, &
&PFPR=PFPR)
ENDIF
-!$acc update self(ZINPRI)
+!$acc kernels
PINPRS(:,:) = PINPRS(:,:) + ZINPRI(:,:)
+!$acc end kernels
!We correct negativities with conservation
!SPLI algorith uses a time-splitting. Inside the loop a temporary m.r. is used.
! It is initialized with the m.r. at T and is modified by two tendencies:
@@ -633,9 +700,7 @@ PRINT *,'OPENACC: RAIN_ICE_RED::HSEDIM=SPLI being implemented'
&PRIS, PRSS, PRGS, &
&PTHS, ZZ_LVFACT, ZZ_LSFACT, PRHS)
ELSE
- WRITE(*,*) ' STOP'
- WRITE(*,*) ' NO SEDIMENTATION SCHEME FOR HSEDIM=', HSEDIM
- CALL PRINT_MSG(NVERB_FATAL,'GEN','RAIN_ICE_RED','')
+ CALL PRINT_MSG(NVERB_FATAL,'GEN','RAIN_ICE_RED','no sedimentation scheme for HSEDIM='//TRIM(HSEDIM))
END IF
!
!* 2.2 budget storage
@@ -659,10 +724,16 @@ ENDIF
! optimization by looking for locations where
! the microphysical fields are larger than a minimal value only !!!
!
+#ifndef _OPENACC
IMICRO=0
IF(COUNT(LDMICRO)/=0) IMICRO=RAIN_ICE_COUNTJV(LDMICRO(:,:,:), IIT, IJT, IKT, SIZE(I1), I1(:), I2(:), I3(:))
+#else
+CALL RAIN_ICE_COUNTJV3D_DEVICE(LDMICRO(:,:,:),I1(:),I2(:),I3(:),IMICRO)
+!$acc update self(I1,I2,I3,IMICRO)
+#endif
!Packing
IF(IMICRO>0) THEN
+!$acc kernels
DO JL=1, IMICRO
ZRVT(JL) = PRVT(I1(JL),I2(JL),I3(JL))
ZRCT(JL) = PRCT(I1(JL),I2(JL),I3(JL))
@@ -677,6 +748,8 @@ IF(IMICRO>0) THEN
ZPRES(JL) = PPABST(I1(JL),I2(JL),I3(JL))
ZEXN(JL) = PEXN(I1(JL),I2(JL),I3(JL))
ENDDO
+!$acc end kernels
+!$acc kernels
IF(LEXT_TEND) THEN
DO JL=1, IMICRO
ZEXT_RV(JL) = PRVS(I1(JL),I2(JL),I3(JL)) - ZRVT(JL)*ZINV_TSTEP
@@ -689,13 +762,19 @@ IF(IMICRO>0) THEN
!The th tendency is not related to a mixing ratio change, there is no exn/exnref issue here
ENDDO
ENDIF
+!$acc end kernels
+!$acc kernels
IF(HSUBG_AUCV_RC=='PDF ' .AND. CSUBG_PR_PDF=='SIGM') THEN
DO JL=1, IMICRO
ZSIGMA_RC(JL) = PSIGS(I1(JL),I2(JL),I3(JL))*2.
! ZSIGMA_RC(JL) = MAX(PSIGS(I1(JL),I2(JL),I3(JL)) * 2., 1.E-12)
ENDDO
+ ELSE !useful when doing calls to MPPDB_CHECK
+ ZSIGMA_RC(:) = XUNDEF
ENDIF
+!$acc end kernels
IF(KRR==7) THEN
+!$acc kernels
DO JL=1, IMICRO
ZRHT(JL) = PRHT(I1(JL),I2(JL),I3(JL))
ENDDO
@@ -704,11 +783,15 @@ IF(IMICRO>0) THEN
ZEXT_RH(JL) = PRHS(I1(JL),I2(JL),I3(JL)) - ZRHT(JL)*ZINV_TSTEP
ENDDO
ENDIF
+!$acc end kernels
ELSE
+!$acc kernels
ZRHT(:)=0.
IF(LEXT_TEND) ZEXT_RH(:)=0.
+!$acc end kernels
ENDIF
IF(LBU_ENABLE) THEN
+!$acc kernels
ZTOT_RVHENI(:)=0.
ZTOT_RCHONI(:)=0.
ZTOT_RRHONG(:)=0.
@@ -754,6 +837,7 @@ IF(IMICRO>0) THEN
ZTOT_RRDRYH(:)=0.
ZTOT_RGDRYH(:)=0.
ZTOT_RDRYHG(:)=0.
+!$acc end kernels
ENDIF
ENDIF
!-------------------------------------------------------------------------------
@@ -763,16 +847,22 @@ ENDIF
!
!Maximum number of iterations
!We only count real iterations (those for which we *compute* tendencies)
+!$acc kernels
INB_ITER_MAX=NMAXITER
IF(XTSTEP_TS/=0.)THEN
INB_ITER_MAX=MAX(1, INT(PTSTEP/XTSTEP_TS)) !At least the number of iterations needed for the time-splitting
ZTSTEP=PTSTEP/INB_ITER_MAX
- INB_ITER_MAX=MAX(NMAXITER, INB_ITER_MAX) !Fot the case XMRSTEP/=0. at the same time
+ INB_ITER_MAX=MAX(NMAXITER, INB_ITER_MAX) !For the case XMRSTEP/=0. at the same time
ENDIF
IITER(:)=0
ZTIME(:)=0. ! Current integration time (all points may have a different integration time)
+!$acc end kernels
+
+!$acc update self(INB_ITER_MAX,ZTSTEP,IITER,ZTIME)
+
DO WHILE(ANY(ZTIME(:)<PTSTEP)) ! Loop to *really* compute tendencies
IF(XMRSTEP/=0.) THEN
+!$acc kernels
! In this case we need to remember the mixing ratios used to compute the tendencies
! because when mixing ratio has evolved more than a threshold, we must re-compute tendecies
Z0RVT(:)=ZRVT(:)
@@ -782,33 +872,50 @@ DO WHILE(ANY(ZTIME(:)<PTSTEP)) ! Loop to *really* compute tendencies
Z0RST(:)=ZRST(:)
Z0RGT(:)=ZRGT(:)
Z0RHT(:)=ZRHT(:)
+!$acc end kernels
ENDIF
IF(XTSTEP_TS/=0.) THEN
+!$acc kernels
! In this case we need to remember the time when tendencies were computed
! because when time has evolved more than a limit, we must re-compute tendecies
ZTIME_LASTCALL(:)=ZTIME(:)
+!$acc end kernels
ENDIF
+!$acc kernels
LLCOMPUTE(:)=ZTIME(:)<PTSTEP ! Compuation only for points for which integration time has not reached the timestep
LSOFT=.FALSE. ! We *really* compute the tendencies
+!$acc end kernels
+!$acc kernels
WHERE(LLCOMPUTE(:))
IITER(:)=IITER(:)+1
END WHERE
+!$acc end kernels
+
+
+!$acc update self(LLCOMPUTE)
+
DO WHILE(ANY(LLCOMPUTE(:))) ! Loop to adjust tendencies when we cross the 0°C or when a specie disappears
+!$acc kernels
ZZT(:) = ZTHT(:) * ZEXN(:)
+!$acc end kernels
IF(KRR==7) THEN
+!$acc kernels
ZLSFACT(:)=(XLSTT+(XCPV-XCI)*(ZZT(:)-XTT)) &
/( (XCPD + XCPV*ZRVT(:) + XCL*(ZRCT(:)+ZRRT(:)) &
+ XCI*(ZRIT(:)+ZRST(:)+ZRGT(:)+ZRHT(:)))*ZEXN(:) )
ZLVFACT(:)=(XLVTT+(XCPV-XCL)*(ZZT(:)-XTT)) &
/( (XCPD + XCPV*ZRVT(:) + XCL*(ZRCT(:)+ZRRT(:)) &
+ XCI*(ZRIT(:)+ZRST(:)+ZRGT(:)+ZRHT(:)))*ZEXN(:) )
+!$acc end kernels
ELSE
+!$acc kernels
ZLSFACT(:)=(XLSTT+(XCPV-XCI)*(ZZT(:)-XTT)) &
/( (XCPD + XCPV*ZRVT(:) + XCL*(ZRCT(:)+ZRRT(:)) &
+ XCI*(ZRIT(:)+ZRST(:)+ZRGT(:)))*ZEXN(:) )
ZLVFACT(:)=(XLVTT+(XCPV-XCL)*(ZZT(:)-XTT)) &
/( (XCPD + XCPV*ZRVT(:) + XCL*(ZRCT(:)+ZRRT(:)) &
+ XCI*(ZRIT(:)+ZRST(:)+ZRGT(:)))*ZEXN(:) )
+!$acc end kernels
ENDIF
!
!*** 4.1 Tendecies computation
@@ -1088,7 +1195,7 @@ IF(IMICRO>0) THEN
ZHLC_LRC3D(:,:,:) = UNPACK(ZHLC_LRC(:), MASK=LDMICRO(:,:,:), FIELD=ZW(:,:,:))
PCIT(:,:,:) = UNPACK(ZCIT(:), MASK=LDMICRO(:,:,:), FIELD=PCIT(:,:,:))
ELSE
- ZRAINFR(:,:,:)=0.
+! ZRAINFR(:,:,:)=0.
ZHLC_HCF3D(:,:,:)=0.
ZHLC_LCF3D(:,:,:)=0.
ZHLC_HRC3D(:,:,:)=0.
@@ -1104,8 +1211,10 @@ ENDIF
!* 6. COMPUTES THE SLOW COLD PROCESS SOURCES OUTSIDE OF LDMICRO POINTS
! ----------------------------------------------------------------
!
-CALL ICE4_NUCLEATION_WRAPPER(IIT, IJT, IKT, .NOT. LDMICRO, &
- PTHT, PPABST, PRHODREF, PEXN, ZZ_LSFACT/PEXN, ZT, &
+GDNOTMICRO = .NOT.LDMICRO
+ZLSFACT3D(:,:,:) = ZZ_LSFACT(:,:,:)/PEXN(:,:,:)
+CALL ICE4_NUCLEATION_WRAPPER(IIT, IJT, IKT, GDNOTMICRO, &
+ PTHT, PPABST, PRHODREF, PEXN, ZLSFACT3D, ZT, &
PRVT, &
PCIT, ZZ_RVHENI_MR)
!$acc update device(PRIS,PRVS,PTHS)
@@ -1155,7 +1264,6 @@ ZW_THS(:,:,:) = (ZW_RCS(:,:,:)+ZW_RRS(:,:,:) )*ZZ_LVF
& (ZW_RIS(:,:,:)+ZW_RSS(:,:,:)+ZW_RGS(:,:,:)+ZW_RHS(:,:,:))*ZZ_LSFACT(:,:,:)
!We apply these tendencies to the S variables
!$acc update device(PRVS,PRCS,PRRS,PRIS,PRSS,PRGS,PRHS,PTHS)
-!$acc kernels
ZW_RVS(:,:,:) = PRVS(:,:,:) + ZW_RVS(:,:,:)
ZW_RCS(:,:,:) = PRCS(:,:,:) + ZW_RCS(:,:,:)
ZW_RRS(:,:,:) = PRRS(:,:,:) + ZW_RRS(:,:,:)
@@ -1164,7 +1272,6 @@ ZW_RSS(:,:,:) = PRSS(:,:,:) + ZW_RSS(:,:,:)
ZW_RGS(:,:,:) = PRGS(:,:,:) + ZW_RGS(:,:,:)
IF(KRR==7) ZW_RHS(:,:,:) = PRHS(:,:,:) + ZW_RHS(:,:,:)
ZW_THS(:,:,:) = PTHS(:,:,:) + ZW_THS(:,:,:)
-!$acc end kernels
!We correct negativities with conservation
CALL CORRECT_NEGATIVITIES(KRR, ZW_RVS, ZW_RCS, ZW_RRS, &
&ZW_RIS, ZW_RSS, ZW_RGS, &
@@ -1587,9 +1694,7 @@ IF(LSEDIM_AFTER) THEN
&PRIS, PRSS, PRGS, &
&PTHS, ZZ_LVFACT, ZZ_LSFACT, PRHS)
ELSE
- WRITE(*,*) ' STOP'
- WRITE(*,*) ' NO SEDIMENTATION SCHEME FOR HSEDIM=', HSEDIM
- CALL PRINT_MSG(NVERB_FATAL,'GEN','RAIN_ICE_RED','')
+ CALL PRINT_MSG(NVERB_FATAL,'GEN','RAIN_ICE_RED','no sedimentation scheme for HSEDIM='//TRIM(HSEDIM))
END IF
!
!* 8.2 budget storage
@@ -1603,11 +1708,13 @@ IF(LSEDIM_AFTER) THEN
IF ( KRR == 7 .AND. LBUDGET_RH) &
CALL BUDGET (PRHS(:,:,:)*PRHODJ(:,:,:), 12, 'SEDI_BU_RRH')
!
- !sedimentation of rain fraction
- CALL ICE4_RAINFR_VERT(IIB, IIE, IIT, IJB, IJE, IJT, IKB, IKE, IKT, KKL, ZRAINFR, PRRS(:,:,:)*PTSTEP)
-
+!PW: ZRAINFR (only (IN)OUT variable) not used after
+! !sedimentation of rain fraction
+! CALL ICE4_RAINFR_VERT(IIB, IIE, IIT, IJB, IJE, IJT, IKB, IKE, IKT, KKL, ZRAINFR, PRRS(:,:,:)*PTSTEP)
ENDIF
!
+
+
!PW:to remove when all on GPU
!$acc update device(PCIT,PTHS,PRVS,PRCS,PRRS,PRIS,PRSS,PRGS,PINDEP,PINPRC,PINPRR,PINPRR3D,PEVAP3D,PINPRS,PINPRG)
IF (PRESENT(PINPRH)) THEN
@@ -1678,6 +1785,61 @@ CONTAINS
!
END FUNCTION RAIN_ICE_COUNTJV
!
+#ifdef _OPENACC
+ SUBROUTINE RAIN_ICE_COUNTJV3D_DEVICE(LTAB,I1,I2,I3,IC)
+!
+!* 0. DECLARATIONS
+! ------------
+!
+IMPLICIT NONE
+!
+LOGICAL, DIMENSION(:,:,:), INTENT(IN) :: LTAB ! Mask
+INTEGER, DIMENSION(:), INTENT(OUT) :: I1,I2,I3 ! Used to replace the COUNT and PACK
+INTEGER, INTENT(OUT) :: IC ! Count
+!$acc declare present(LTAB,I1,I2,I3,IC)
+!
+!* 0.2 declaration of local variables
+!
+INTEGER :: JI,JJ,JK,IDX
+!
+!-------------------------------------------------------------------------------
+!
+!$acc kernels present(LTAB,I1,I2,I3)
+
+!To allow comparisons... (I1/I2/I3 are not fully used)
+!Can be removed in production
+! I1(:) = -999
+! I2(:) = -999
+! I3(:) = -999
+
+
+IC = 0
+!Warning: if "independent" is set, content of I1, I2 and I3 can vary between 2
+! different runs of this subroutine BUT final result should be the same
+!Comment the following line + atomic directives to have consistent values for debugging
+!Warning: huge impact on performance
+!acc loop collapse(3) private(IDX) independent
+DO JK = 1,SIZE(LTAB,3)
+ DO JJ = 1,SIZE(LTAB,2)
+ DO JI = 1,SIZE(LTAB,1)
+ IF( LTAB(JI,JJ,JK) ) THEN
+!acc atomic capture
+ IC = IC +1
+ IDX = IC
+!acc end atomic
+ I1(IDX) = JI
+ I2(IDX) = JJ
+ I3(IDX) = JK
+ END IF
+ END DO
+ END DO
+END DO
+!$acc end kernels
+! print *,'PW: RAIN_ICE_COUNTJV3D_DEVICE: IC=',IC
+!
+END SUBROUTINE RAIN_ICE_COUNTJV3D_DEVICE
+#endif
+ !
!
SUBROUTINE CORRECT_NEGATIVITIES(KRR, PRV, PRC, PRR, &
&PRI, PRS, PRG, &
@@ -1690,44 +1852,76 @@ CONTAINS
REAL, DIMENSION(:,:,:), INTENT(IN) :: PLVFACT, PLSFACT
REAL, DIMENSION(:,:,:), OPTIONAL, INTENT(INOUT) :: PRH
!
+!$acc declare present(KRR,PRV,PRC,PRR,PRI,PRS,PRG,PTH,PLVFACT,PLSFACT,PRH)
+ !
+ LOGICAL, DIMENSION(SIZE(PRV,1), SIZE(PRV,2), SIZE(PRV,3)) :: GW
REAL, DIMENSION(SIZE(PRV,1), SIZE(PRV,2), SIZE(PRV,3)) :: ZW
+!$acc declare create(GW,ZW)
+ !
!
+ IF (MPPDB_INITIALIZED) THEN
+ !Check all IN arrays
+ CALL MPPDB_CHECK3D(PLVFACT,"CORRECT_NEGATIVITIES beg:PLVFACT",PRECISION)
+ CALL MPPDB_CHECK3D(PLSFACT,"CORRECT_NEGATIVITIES beg:PLSFACT",PRECISION)
+ !Check all INOUT arrays
+ CALL MPPDB_CHECK3D(PRV,"CORRECT_NEGATIVITIES beg:PRV",PRECISION)
+ CALL MPPDB_CHECK3D(PRC,"CORRECT_NEGATIVITIES beg:PRC",PRECISION)
+ CALL MPPDB_CHECK3D(PRR,"CORRECT_NEGATIVITIES beg:PRR",PRECISION)
+ CALL MPPDB_CHECK3D(PRI,"CORRECT_NEGATIVITIES beg:PRI",PRECISION)
+ CALL MPPDB_CHECK3D(PRS,"CORRECT_NEGATIVITIES beg:PRS",PRECISION)
+ CALL MPPDB_CHECK3D(PRG,"CORRECT_NEGATIVITIES beg:PRG",PRECISION)
+ IF(PRESENT(PRH)) CALL MPPDB_CHECK3D(PRH,"CORRECT_NEGATIVITIES beg:PRH",PRECISION)
+ CALL MPPDB_CHECK3D(PTH,"CORRECT_NEGATIVITIES beg:PTH",PRECISION)
+ END IF
+ !
+!$acc kernels
!We correct negativities with conservation
! 1) deal with negative values for mixing ratio, except for vapor
- WHERE(PRC(:,:,:)<0.)
+ GW(:,:,:) = PRC(:,:,:)<0.
+ WHERE(GW(:,:,:))
PRV(:,:,:)=PRV(:,:,:)+PRC(:,:,:)
PTH(:,:,:)=PTH(:,:,:)-PRC(:,:,:)*PLVFACT(:,:,:)
PRC(:,:,:)=0.
ENDWHERE
- WHERE(PRR(:,:,:)<0.)
+ GW(:,:,:) = PRR(:,:,:)<0.
+ WHERE(GW(:,:,:))
PRV(:,:,:)=PRV(:,:,:)+PRR(:,:,:)
PTH(:,:,:)=PTH(:,:,:)-PRR(:,:,:)*PLVFACT(:,:,:)
PRR(:,:,:)=0.
ENDWHERE
- WHERE(PRI(:,:,:)<0.)
+ GW(:,:,:) = PRI(:,:,:)<0.
+ WHERE(GW(:,:,:))
PRV(:,:,:)=PRV(:,:,:)+PRI(:,:,:)
PTH(:,:,:)=PTH(:,:,:)-PRI(:,:,:)*PLSFACT(:,:,:)
PRI(:,:,:)=0.
ENDWHERE
- WHERE(PRS(:,:,:)<0.)
+ GW(:,:,:) = PRS(:,:,:)<0.
+ WHERE(GW(:,:,:))
PRV(:,:,:)=PRV(:,:,:)+PRS(:,:,:)
PTH(:,:,:)=PTH(:,:,:)-PRS(:,:,:)*PLSFACT(:,:,:)
PRS(:,:,:)=0.
ENDWHERE
- WHERE(PRG(:,:,:)<0.)
+ GW(:,:,:) = PRG(:,:,:)<0.
+ WHERE(GW(:,:,:))
PRV(:,:,:)=PRV(:,:,:)+PRG(:,:,:)
PTH(:,:,:)=PTH(:,:,:)-PRG(:,:,:)*PLSFACT(:,:,:)
PRG(:,:,:)=0.
ENDWHERE
+!$acc end kernels
IF(KRR==7) THEN
- WHERE(PRH(:,:,:)<0.)
+!$acc kernels
+ GW(:,:,:) = PRH(:,:,:)<0.
+ WHERE(GW(:,:,:))
PRV(:,:,:)=PRV(:,:,:)+PRH(:,:,:)
PTH(:,:,:)=PTH(:,:,:)-PRH(:,:,:)*PLSFACT(:,:,:)
PRH(:,:,:)=0.
ENDWHERE
+!$acc end kernels
ENDIF
+!$acc kernels
! 2) deal with negative vapor mixing ratio
- WHERE(PRV(:,:,:)<0. .AND. PRC(:,:,:)+PRI(:,:,:)>0.)
+ GW(:,:,:) = PRV(:,:,:)<0. .AND. PRC(:,:,:)+PRI(:,:,:)>0.
+ WHERE(GW(:,:,:))
! for rc and ri, we keep ice fraction constant
ZW(:,:,:)=MIN(1., -PRV(:,:,:)/(PRC(:,:,:)+PRI(:,:,:))) ! Proportion of rc+ri to convert into rv
PTH(:,:,:)=PTH(:,:,:)-ZW(:,:,:)*(PRC(:,:,:)*PLVFACT(:,:,:)+PRI(:,:,:)*PLSFACT(:,:,:))
@@ -1735,33 +1929,52 @@ CONTAINS
PRC(:,:,:)=(1.-ZW(:,:,:))*PRC(:,:,:)
PRI(:,:,:)=(1.-ZW(:,:,:))*PRI(:,:,:)
ENDWHERE
- WHERE(PRV(:,:,:)<0. .AND. PRR(:,:,:)>0.)
+ GW(:,:,:) = PRV(:,:,:)<0. .AND. PRR(:,:,:)>0.
+ WHERE(GW(:,:,:))
ZW(:,:,:)=MIN(PRR(:,:,:), -PRV(:,:,:)) ! Quantity of rr to convert into rv
PRV(:,:,:)=PRV(:,:,:)+ZW(:,:,:)
PRR(:,:,:)=PRR(:,:,:)-ZW(:,:,:)
PTH(:,:,:)=PTH(:,:,:)-ZW(:,:,:)*PLVFACT(:,:,:)
ENDWHERE
- WHERE(PRV(:,:,:)<0. .AND. PRS(:,:,:)>0.)
+ GW(:,:,:) = PRV(:,:,:)<0. .AND. PRS(:,:,:)>0.
+ WHERE(GW(:,:,:))
ZW(:,:,:)=MIN(PRS(:,:,:), -PRV(:,:,:)) ! Quantity of rs to convert into rv
PRV(:,:,:)=PRV(:,:,:)+ZW(:,:,:)
PRS(:,:,:)=PRS(:,:,:)-ZW(:,:,:)
PTH(:,:,:)=PTH(:,:,:)-ZW(:,:,:)*PLSFACT(:,:,:)
ENDWHERE
- WHERE(PRV(:,:,:)<0. .AND. PRG(:,:,:)>0.)
+ GW(:,:,:) = PRV(:,:,:)<0. .AND. PRG(:,:,:)>0.
+ WHERE(GW(:,:,:))
ZW(:,:,:)=MIN(PRG(:,:,:), -PRV(:,:,:)) ! Quantity of rg to convert into rv
PRV(:,:,:)=PRV(:,:,:)+ZW(:,:,:)
PRG(:,:,:)=PRG(:,:,:)-ZW(:,:,:)
PTH(:,:,:)=PTH(:,:,:)-ZW(:,:,:)*PLSFACT(:,:,:)
ENDWHERE
+!$acc end kernels
IF(KRR==7) THEN
- WHERE(PRV(:,:,:)<0. .AND. PRH(:,:,:)>0.)
+!$acc kernels
+ GW(:,:,:) = PRV(:,:,:)<0. .AND. PRH(:,:,:)>0.
+ WHERE(GW(:,:,:))
ZW(:,:,:)=MIN(PRH(:,:,:), -PRV(:,:,:)) ! Quantity of rh to convert into rv
PRV(:,:,:)=PRV(:,:,:)+ZW(:,:,:)
PRH(:,:,:)=PRH(:,:,:)-ZW(:,:,:)
PTH(:,:,:)=PTH(:,:,:)-ZW(:,:,:)*PLSFACT(:,:,:)
ENDWHERE
+!$acc end kernels
ENDIF
!
+ IF (MPPDB_INITIALIZED) THEN
+ !Check all INOUT arrays
+ CALL MPPDB_CHECK3D(PRV,"CORRECT_NEGATIVITIES end:PRV",PRECISION)
+ CALL MPPDB_CHECK3D(PRC,"CORRECT_NEGATIVITIES end:PRC",PRECISION)
+ CALL MPPDB_CHECK3D(PRR,"CORRECT_NEGATIVITIES end:PRR",PRECISION)
+ CALL MPPDB_CHECK3D(PRI,"CORRECT_NEGATIVITIES end:PRI",PRECISION)
+ CALL MPPDB_CHECK3D(PRS,"CORRECT_NEGATIVITIES end:PRS",PRECISION)
+ CALL MPPDB_CHECK3D(PRG,"CORRECT_NEGATIVITIES end:PRG",PRECISION)
+ IF(PRESENT(PRH)) CALL MPPDB_CHECK3D(PRH,"CORRECT_NEGATIVITIES end:PRH",PRECISION)
+ CALL MPPDB_CHECK3D(PTH,"CORRECT_NEGATIVITIES end:PTH",PRECISION)
+ END IF
+ !
!
END SUBROUTINE CORRECT_NEGATIVITIES
!
diff --git a/src/MNH/resolved_cloud.f90 b/src/MNH/resolved_cloud.f90
index 47da8f9cc053c806f71ac60682e73670d62c8630..318f3fb4ba1aa101ca035f6e7678dd7c8d93f27e 100644
--- a/src/MNH/resolved_cloud.f90
+++ b/src/MNH/resolved_cloud.f90
@@ -443,6 +443,7 @@ REAL, DIMENSION(:,:), OPTIONAL, INTENT(IN) :: PTOWN ! Town fraction
!
! IN variables
!
+!$acc declare copyin(KRR,OSEDIC,HSUBG_AUCV,PTSTEP,OWARM)
!$acc declare present(PRHODJ) &
!$acc & copyin(PZZ,PRHODREF,PEXNREF,PPABST,PTHT,PSIGS,PSIGQSAT,PMFCONV,PTHM,PPABSM,PRCM,PW_ACT,PDTHRAD,&
!$acc & PCF_MF,PRC_MF,PRI_MF)
@@ -458,6 +459,8 @@ REAL, DIMENSION(:,:), OPTIONAL, INTENT(IN) :: PTOWN ! Town fraction
!
!$acc declare create(PSRCS)
!
+! Variables from modules
+!$acc declare copyin(CSEDIM)
!
!* 0.2 Declarations of local variables :
!
@@ -466,12 +469,14 @@ INTEGER :: IIB ! Define the physical domain
INTEGER :: IIE !
INTEGER :: IJB !
INTEGER :: IJE !
-INTEGER :: IKB !
+INTEGER :: IKA,IKB !
INTEGER :: IKE !
INTEGER :: IKU
+INTEGER :: IKL
INTEGER :: IINFO_ll ! return code of parallel routine
INTEGER :: JI,JJ,JK,JL
!
+!$acc declare create(IKA,IKU,IKL)
!
!
REAL, DIMENSION(SIZE(PZZ,1),SIZE(PZZ,2),SIZE(PZZ,3)):: ZDZZ
@@ -577,9 +582,18 @@ IF (MPPDB_INITIALIZED) THEN
END IF
!
CALL GET_INDICE_ll (IIB,IJB,IIE,IJE)
+!$acc kernels
+IKA=1
+!$acc end kernels
IKB=1+JPVEXT
IKE=SIZE(PZZ,3) - JPVEXT
IKU=SIZE(PZZ,3)
+!$acc update device(IKU)
+!
+!$acc kernels
+IKL = 1
+!$acc end kernels
+!$acc update host(IKL)
!
GWEST = LWEST_ll()
GEAST = LEAST_ll()
@@ -1186,7 +1200,7 @@ PRINT *,'OPENACC: RESOLVED_CLOUD::LRED not yet implemented'
!$acc update device(PCIT)
!$acc data copyin(PSEA,PTOWN)
- CALL RAIN_ICE_RED ( OSEDIC,CSEDIM, HSUBG_AUCV, OWARM,1,IKU,1, &
+ CALL RAIN_ICE_RED ( OSEDIC,CSEDIM, HSUBG_AUCV, OWARM,IKA,IKU,IKL, &
PTSTEP, KRR, LLMICRO, ZEXN, &
ZDZZ, PRHODJ, PRHODREF, PEXNREF, PPABST, PCIT,PCLDFR,&
PTHT, PRT(:,:,:,1), PRT(:,:,:,2), &
@@ -1292,7 +1306,7 @@ CALL ABORT
PRS(:,:,:,5)>ZRSMIN(5) .OR. &
PRS(:,:,:,6)>ZRSMIN(6) .OR. &
PRS(:,:,:,7)>ZRSMIN(7)
- CALL RAIN_ICE_RED ( OSEDIC,CSEDIM, HSUBG_AUCV, OWARM,1,IKU,1, &
+ CALL RAIN_ICE_RED ( OSEDIC,CSEDIM, HSUBG_AUCV, OWARM,IKA,IKU,IKL, &
PTSTEP, KRR, LLMICRO, ZEXN, &
ZDZZ, PRHODJ, PRHODREF, PEXNREF, PPABST, PCIT, PCLDFR,&
PTHT, PRT(:,:,:,1), PRT(:,:,:,2), &