diff --git a/src/MNH/ice4_nucleation_wrapper.f90 b/src/MNH/ice4_nucleation_wrapper.f90
index 8a2711fc7970a3b4b7cc2e0f8d7b3424bd537a42..662f7d67067af16f61606bdcad524f85bf7c5602 100644
--- a/src/MNH/ice4_nucleation_wrapper.f90
+++ b/src/MNH/ice4_nucleation_wrapper.f90
@@ -1,6 +1,6 @@
-!MNH_LIC Copyright 1994-2014 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC Copyright 1994-2019 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_NUCLEATION_WRAPPER
INTERFACE
@@ -64,20 +64,24 @@ REAL, DIMENSION(KIT,KJT,KKT), INTENT(OUT) :: PRVHENI_MR ! Mixing ratio chang
!
!* 0.2 declaration of local variables
!
-INTEGER :: JL ! and PACK intrinsics
-LOGICAL, DIMENSION(SIZE(PRVT,1),SIZE(PRVT,2),SIZE(PRVT,3)) :: GNEGT ! Test where to compute the HEN process
-INTEGER :: INEGT
-INTEGER, DIMENSION(COUNT(PT<XTT .AND. LDMASK)) :: I1,I2,I3 ! Used to replace the COUNT
-REAL, DIMENSION(COUNT(PT<XTT .AND. LDMASK)) :: ZZT, & ! Temperature
- ZPRES, & ! Pressure
- ZRVT, & ! Water vapor m.r. at t
- ZCIT, & ! Pristine ice conc. at t
- ZTHT, & ! Theta at t
- ZRHODREF, &
- ZEXN, &
- ZLSFACT, &
- ZRVHENI_MR, &
- ZB_TH, ZB_RV, ZB_RI
+INTEGER :: IDX, JI, JJ, JK
+INTEGER :: JL
+INTEGER :: INEGT, INEGT_TMP
+INTEGER, DIMENSION(:), ALLOCATABLE :: I1,I2,I3
+LOGICAL :: GDSOFT
+LOGICAL, DIMENSION(:), ALLOCATABLE :: GLDCOMPUTE
+LOGICAL, DIMENSION(KIT,KJT,KKT) :: GNEGT ! Test where to compute the HEN process
+REAL, DIMENSION(:), ALLOCATABLE :: ZZT, & ! Temperature
+ ZPRES, & ! Pressure
+ ZRVT, & ! Water vapor m.r. at t
+ ZCIT, & ! Pristine ice conc. at t
+ ZTHT, & ! Theta at t
+ ZRHODREF, &
+ ZEXN, &
+ ZLSFACT, &
+ ZRVHENI_MR, &
+ ZB_TH, ZB_RV, ZB_RI
+!
!-------------------------------------------------------------------------------
!
!
@@ -86,15 +90,31 @@ REAL, DIMENSION(COUNT(PT<XTT .AND. LDMASK)) :: ZZT, & ! Temperature
! the temperature is negative only !!!
!
GNEGT(:,:,:)=PT(:,:,:)<XTT .AND. LDMASK
-INEGT=0
-IF(COUNT(GNEGT)/=0) INEGT=ICE4_NUCLEATION_COUNTJV(GNEGT(:,:,:), KIT, KJT, KKT, SIZE(I1), I1(:), I2(:), I3(:))
+INEGT = COUNT(GNEGT(:,:,:))
+!
+ALLOCATE(GLDCOMPUTE(INEGT))
+ALLOCATE(I1(INEGT),I2(INEGT),I3(INEGT))
+ALLOCATE(ZZT(INEGT))
+ALLOCATE(ZPRES(INEGT))
+ALLOCATE(ZRVT(INEGT))
+ALLOCATE(ZCIT(INEGT))
+ALLOCATE(ZTHT(INEGT))
+ALLOCATE(ZRHODREF(INEGT))
+ALLOCATE(ZEXN(INEGT))
+ALLOCATE(ZLSFACT(INEGT))
+ALLOCATE(ZRVHENI_MR(INEGT))
+ALLOCATE(ZB_TH(INEGT))
+ALLOCATE(ZB_RV(INEGT))
+ALLOCATE(ZB_RI(INEGT))
!
ZB_TH(:) = 0.
ZB_RV(:) = 0.
ZB_RI(:) = 0.
!
+IF(INEGT>0) INEGT_TMP=ICE4_NUCLEATION_COUNTJV(GNEGT(:,:,:), KIT, KJT, KKT, SIZE(I1), I1(:), I2(:), I3(:))
+!
PRVHENI_MR(:,:,:)=0.
-IF(INEGT>=1) THEN
+IF(INEGT>0) THEN
DO JL=1, INEGT
ZRVT(JL)=PRVT(I1(JL), I2(JL), I3(JL))
ZCIT(JL)=PCIT(I1(JL), I2(JL), I3(JL))
@@ -105,15 +125,20 @@ IF(INEGT>=1) 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)
PRVHENI_MR(:,:,:)=UNPACK(ZRVHENI_MR(:), MASK=GNEGT(:,:,:), FIELD=0.0)
- PCIT(:,:,:)=UNPACK(ZCIT(:), MASK=GNEGT(:,:,:), FIELD=PCIT(:,:,:))
+ 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)
+!
CONTAINS
FUNCTION ICE4_NUCLEATION_COUNTJV(LTAB,KIT,KJT,KKT,KSIZE,I1,I2,I3) RESULT(IC)
IMPLICIT NONE
diff --git a/src/MNH/ice4_sedimentation_split.f90 b/src/MNH/ice4_sedimentation_split.f90
index 10dd37ae9af8494bc6be65bfcc7a973060c8ee0e..cb0a147d070b865e8e9391bacf7ef143c3727311 100644
--- a/src/MNH/ice4_sedimentation_split.f90
+++ b/src/MNH/ice4_sedimentation_split.f90
@@ -77,9 +77,9 @@ SUBROUTINE ICE4_SEDIMENTATION_SPLIT(KIB, KIE, KIT, KJB, KJE, KJT, KKB, KKE, KKTB
!* 0. DECLARATIONS
! ------------
!
-USE MODD_CST, ONLY: XCPD,XP00,XRD,XRHOLW
+USE MODD_CST, ONLY: XRHOLW
USE MODD_PARAM_ICE, ONLY: XSPLIT_MAXCFL
-USE MODD_RAIN_ICE_DESCR, ONLY: XALPHAC,XALPHAC2,XCONC_LAND,XCONC_SEA,XCONC_URBAN,XLBC,XNUC,XNUC2,XRTMIN
+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
@@ -128,57 +128,62 @@ REAL, DIMENSION(KIT,KJT,KKT,KRR), OPTIONAL, INTENT(OUT) :: PFPR ! upper-air
!* 0.2 declaration of local variables
!
!
-LOGICAL, DIMENSION(SIZE(PRHODREF,1),SIZE(PRHODREF,2),SIZE(PRHODREF,3)) &
- :: GSEDIM ! Test where to compute the SED processes
-LOGICAL, DIMENSION(SIZE(PRHODREF,1),SIZE(PRHODREF,2)):: GDEP
-INTEGER , DIMENSION(SIZE(GSEDIM)) :: I1,I2,I3 ! Used to replace the COUNT
-
+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
& ZFSEDC, &
& ZPRCS,ZPRRS,ZPRIS,ZPRSS,ZPRGS,ZPRHS, & ! Mixing ratios created during the time step
- & ZW, & ! work array
& ZRCT, &
& ZRRT, &
& ZRIT, &
& ZRST, &
& ZRGT, &
& ZRHT
-REAL, DIMENSION(SIZE(PRHODREF,1),SIZE(PRHODREF,2),0:SIZE(PRHODREF,3)+1) :: ZWSED ! sedimentation fluxes
-REAL, DIMENSION(SIZE(PRHODREF,1),SIZE(PRHODREF,2)) :: ZCONC_TMP ! Weighted concentration
-REAL, DIMENSION(SIZE(PRHODREF,1),SIZE(PRHODREF,2)) :: ZREMAINT ! Remaining time until the timestep end
-REAL :: ZINVTSTEP
-INTEGER :: ISEDIM ! ! Case number of sedimentation
-INTEGER :: JK
-REAL, DIMENSION(SIZE(XRTMIN)) :: ZRSMIN
+!
!-------------------------------------------------------------------------------
!
!-------------------------------------------------------------------------------
!
!
+GDEPOSC = ODEPOSC
+GSEDIC = OSEDIC
+IRR = KRR
+!
+IF (PRESENT(PFPR)) THEN
+ GPRESENT_PFPR = .TRUE.
+ELSE
+ GPRESENT_PFPR = .FALSE.
+END IF
+!
+IF (PRESENT(PSEA)) THEN
+ GPRESENT_PSEA = .TRUE.
+ELSE
+ GPRESENT_PSEA = .FALSE.
+END IF
+!
! O. Initialization of for sedimentation
!
ZINVTSTEP=1./PTSTEP
-ZRSMIN(:) = XRTMIN(:) * ZINVTSTEP
-IF (OSEDIC) PINPRC (:,:) = 0.
-IF (ODEPOSC) PINDEP (:,:) = 0.
-PINPRR (:,:) = 0.
-PINPRI (:,:) = 0.
-PINPRS (:,:) = 0.
-PINPRG (:,:) = 0.
-IF ( KRR == 7 ) PINPRH (:,:) = 0.
-IF (PRESENT(PFPR)) PFPR(:,:,:,:) = 0.
+IF (GPRESENT_PFPR) THEN
+ PFPR(:,:,:,:) = 0.
+END IF
!
!* 1. Parameters for cloud sedimentation
!
-IF (OSEDIC) THEN
+IF (GSEDIC) THEN
ZRAY(:,:,:) = 0.
ZLBC(:,:,:) = XLBC(1)
ZFSEDC(:,:,:) = XFSEDC(1)
ZCONC3D(:,:,:)= XCONC_LAND
ZCONC_TMP(:,:)= XCONC_LAND
- IF (PRESENT(PSEA)) THEN
+ IF (GPRESENT_PSEA) THEN
ZCONC_TMP(:,:)=PSEA(:,:)*XCONC_SEA+(1.-PSEA(:,:))*XCONC_LAND
DO JK=KKTB, KKTE
ZLBC(:,:,JK) = PSEA(:,:)*XLBC(2)+(1.-PSEA(:,:))*XLBC(1)
@@ -203,12 +208,16 @@ ENDIF
! For optimization we consider each variable separately
!
! External tendecies
-IF (OSEDIC) ZPRCS(:,:,:) = PRCS(:,:,:)-PRCT(:,:,:)*ZINVTSTEP
+IF (GSEDIC) THEN
+ ZPRCS(:,:,:) = PRCS(:,:,:)-PRCT(:,:,:)*ZINVTSTEP
+ENDIF
ZPRRS(:,:,:) = PRRS(:,:,:)-PRRT(:,:,:)*ZINVTSTEP
ZPRIS(:,:,:) = PRIS(:,:,:)-PRIT(:,:,:)*ZINVTSTEP
ZPRSS(:,:,:) = PRSS(:,:,:)-PRST(:,:,:)*ZINVTSTEP
ZPRGS(:,:,:) = PRGS(:,:,:)-PRGT(:,:,:)*ZINVTSTEP
-IF ( KRR == 7 ) ZPRHS(:,:,:) = PRHS(:,:,:)-PRHT(:,:,:)*ZINVTSTEP
+IF ( IRR == 7 ) THEN
+ ZPRHS(:,:,:) = PRHS(:,:,:)-PRHT(:,:,:)*ZINVTSTEP
+END IF
!
! mr values inside the time-splitting loop
ZRCT(:,:,:) = PRCT(:,:,:)
@@ -216,150 +225,85 @@ ZRRT(:,:,:) = PRRT(:,:,:)
ZRIT(:,:,:) = PRIT(:,:,:)
ZRST(:,:,:) = PRST(:,:,:)
ZRGT(:,:,:) = PRGT(:,:,:)
-IF (KRR==7) ZRHT(:,:,:) = PRHT(:,:,:)
+IF (IRR==7) THEN
+ ZRHT(:,:,:) = PRHT(:,:,:)
+END IF
!
-DO JK = KKTB , KKTE
- ZW(:,:,JK) =1./(PRHODREF(:,:,JK)* PDZZ(:,:,JK))
-END DO
+ZW(:,:,KKTB:KKTE) =1./(PRHODREF(:,:,KKTB:KKTE)* PDZZ(:,:,KKTB:KKTE))
!
!
!* 2.1 for cloud
!
-IF (OSEDIC) THEN
- ZREMAINT(:,:) = PTSTEP
- DO WHILE (ANY(ZREMAINT>0.))
- GSEDIM(:,:,:)=.FALSE.
- DO JK = KKTB , KKTE
- GSEDIM(KIB:KIE,KJB:KJE,JK) = &
- (ZRCT(KIB:KIE,KJB:KJE,JK)>XRTMIN(2) .OR. &
- ZPRCS(KIB:KIE,KJB:KJE,JK)>ZRSMIN(2)) .AND. &
- ZREMAINT(KIB:KIE,KJB:KJE)>0.
- ENDDO
- ISEDIM = ICE4_SEDIMENTATION_SPLIT_COUNTJV(GSEDIM(:,:,:),KIT,KJT,KKT,&
- &SIZE(I1),I1(:),I2(:),I3(:))
- CALL INTERNAL_SEDIM_SPLI(KIT, KJT, KKB, KKTB, KKTE, KKT, KKL, KRR, &
- &ISEDIM, GSEDIM, I1, I2, I3, XSPLIT_MAXCFL, ZREMAINT, &
+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, &
&2, &
- &ZRCT, PRCS, ZWSED, PINPRC, ZPRCS, &
- &ZRAY, ZLBC, ZFSEDC, ZCONC3D, PSEA, PTOWN, PFPR=PFPR)
- ENDDO
+ &ZRCT, PRCS, PINPRC, ZPRCS, &
+ &ZRAY, ZLBC, ZFSEDC, ZCONC3D, PFPR=PFPR)
ENDIF
!
!
!* 2.1bis DROPLET DEPOSITION AT THE 1ST LEVEL ABOVE GROUND
!
-IF (ODEPOSC) THEN
- GDEP(:,:) = .FALSE.
- GDEP(KIB:KIE,KJB:KJE) = PRCS(KIB:KIE,KJB:KJE,KKB) >0
- WHERE (GDEP)
- PRCS(:,:,KKB) = PRCS(:,:,KKB) - PVDEPOSC * PRCT(:,:,KKB) / PDZZ(:,:,KKB)
- PINPRC(:,:) = PINPRC(:,:) + PVDEPOSC * PRCT(:,:,KKB) * PRHODREF(:,:,KKB) /XRHOLW
- PINDEP(:,:) = PVDEPOSC * PRCT(:,:,KKB) * PRHODREF(:,:,KKB) /XRHOLW
- END WHERE
+IF (GDEPOSC) THEN
+ PINDEP (:,:) = 0.
+ DO JJ=KJB,KJE
+ DO JI=KIB,KIE
+ IF (PRCS(JI,JJ,KKB)>0.) THEN
+ PRCS(JI,JJ,KKB) = PRCS(JI,JJ,KKB) - PVDEPOSC * PRCT(JI,JJ,KKB) / PDZZ(JI,JJ,KKB)
+ PINPRC(JI,JJ) = PINPRC(JI,JJ) + PVDEPOSC * PRCT(JI,JJ,KKB) * PRHODREF(JI,JJ,KKB) /XRHOLW
+ PINDEP(JI,JJ) = PVDEPOSC * PRCT(JI,JJ,KKB) * PRHODREF(JI,JJ,KKB) /XRHOLW
+ END IF
+ END DO
+ END DO
END IF
!
!* 2.2 for rain
!
-ZREMAINT(:,:) = PTSTEP
-DO WHILE (ANY(ZREMAINT>0.))
- GSEDIM(:,:,:)=.FALSE.
- DO JK = KKTB , KKTE
- GSEDIM(KIB:KIE,KJB:KJE,JK) = &
- (ZRRT(KIB:KIE,KJB:KJE,JK)>XRTMIN(3) .OR. &
- ZPRRS(KIB:KIE,KJB:KJE,JK)>ZRSMIN(3)) .AND. &
- ZREMAINT(KIB:KIE,KJB:KJE)>0.
- ENDDO
- ISEDIM = ICE4_SEDIMENTATION_SPLIT_COUNTJV(GSEDIM(:,:,:),KIT,KJT,KKT,&
- &SIZE(I1),I1(:),I2(:),I3(:))
- CALL INTERNAL_SEDIM_SPLI(KIT, KJT, KKB, KKTB, KKTE, KKT, KKL, KRR, &
- &ISEDIM, GSEDIM, I1, I2, I3, XSPLIT_MAXCFL, ZREMAINT, &
+ CALL INTERNAL_SEDIM_SPLI(KIB,KIE,KIT,KJB,KJE,KJT, KKB, KKTB, KKTE, KKT, KKL, KRR, &
+ &XSPLIT_MAXCFL, &
&PRHODREF, ZW, PDZZ, PPABST, PTHT, PTSTEP, &
&3, &
- &ZRRT, PRRS, ZWSED, PINPRR, ZPRRS, &
- &PFPR=PFPR)
-ENDDO
+ &ZRRT, PRRS, PINPRR, ZPRRS, &
+ PFPR=PFPR)
!
!* 2.3 for pristine ice
!
-ZREMAINT(:,:) = PTSTEP
-DO WHILE (ANY(ZREMAINT>0.))
- GSEDIM(:,:,:)=.FALSE.
- DO JK = KKTB , KKTE
- GSEDIM(KIB:KIE,KJB:KJE,JK) = &
- (ZRIT(KIB:KIE,KJB:KJE,JK)>XRTMIN(4) .OR. &
- ZPRIS(KIB:KIE,KJB:KJE,JK)>ZRSMIN(4)) .AND. &
- ZREMAINT(KIB:KIE,KJB:KJE)>0.
- ENDDO
- ISEDIM = ICE4_SEDIMENTATION_SPLIT_COUNTJV(GSEDIM(:,:,:),KIT,KJT,KKT,&
- &SIZE(I1),I1(:),I2(:),I3(:))
- CALL INTERNAL_SEDIM_SPLI(KIT, KJT, KKB, KKTB, KKTE, KKT, KKL, KRR, &
- &ISEDIM, GSEDIM, I1, I2, I3, XSPLIT_MAXCFL, ZREMAINT, &
+ CALL INTERNAL_SEDIM_SPLI(KIB,KIE,KIT,KJB,KJE,KJT, KKB, KKTB, KKTE, KKT, KKL, KRR, &
+ &XSPLIT_MAXCFL, &
&PRHODREF, ZW, PDZZ, PPABST, PTHT, PTSTEP, &
&4, &
- &ZRIT, PRIS, ZWSED, PINPRI, ZPRIS, PFPR=PFPR)
-ENDDO
+ &ZRIT, PRIS, PINPRI, ZPRIS, &
+ PFPR=PFPR)
!
!* 2.4 for aggregates/snow
!
-ZREMAINT(:,:) = PTSTEP
-DO WHILE (ANY(ZREMAINT>0.))
- GSEDIM(:,:,:)=.FALSE.
- DO JK = KKTB , KKTE
- GSEDIM(KIB:KIE,KJB:KJE,JK) = &
- (ZRST(KIB:KIE,KJB:KJE,JK)>XRTMIN(5) .OR. &
- ZPRSS(KIB:KIE,KJB:KJE,JK)>ZRSMIN(5)) .AND. &
- ZREMAINT(KIB:KIE,KJB:KJE)>0.
- ENDDO
- ISEDIM = ICE4_SEDIMENTATION_SPLIT_COUNTJV(GSEDIM(:,:,:),KIT,KJT,KKT,&
- &SIZE(I1),I1(:),I2(:),I3(:))
- CALL INTERNAL_SEDIM_SPLI(KIT, KJT, KKB, KKTB, KKTE, KKT, KKL, KRR, &
- &ISEDIM, GSEDIM, I1, I2, I3, XSPLIT_MAXCFL, ZREMAINT, &
- &PRHODREF, ZW, PDZZ, PPABST, PTHT, PTSTEP, &
- &5, &
- &ZRST, PRSS, ZWSED, PINPRS, ZPRSS, PFPR=PFPR)
-ENDDO
+ CALL INTERNAL_SEDIM_SPLI(KIB,KIE,KIT,KJB,KJE,KJT, KKB, KKTB, KKTE, KKT, KKL, KRR, &
+ &XSPLIT_MAXCFL, &
+ &PRHODREF, ZW, PDZZ, PPABST, PTHT, PTSTEP, &
+ &5, &
+ &ZRST, PRSS, PINPRS, ZPRSS, &
+ PFPR=PFPR)
!
!* 2.5 for graupeln
!
-ZREMAINT(:,:) = PTSTEP
-DO WHILE (ANY(ZREMAINT>0.))
- GSEDIM(:,:,:)=.FALSE.
- DO JK = KKTB , KKTE
- GSEDIM(KIB:KIE,KJB:KJE,JK) = &
- (ZRGT(KIB:KIE,KJB:KJE,JK)>XRTMIN(6) .OR. &
- ZPRGS(KIB:KIE,KJB:KJE,JK)>ZRSMIN(6)) .AND. &
- ZREMAINT(KIB:KIE,KJB:KJE)>0.
- ENDDO
- ISEDIM = ICE4_SEDIMENTATION_SPLIT_COUNTJV(GSEDIM(:,:,:),KIT,KJT,KKT,&
- &SIZE(I1),I1(:),I2(:),I3(:))
- CALL INTERNAL_SEDIM_SPLI(KIT, KJT, KKB, KKTB, KKTE, KKT, KKL, KRR, &
- &ISEDIM, GSEDIM, I1, I2, I3, XSPLIT_MAXCFL, ZREMAINT, &
- &PRHODREF, ZW, PDZZ, PPABST, PTHT, PTSTEP, &
- &6, &
- &ZRGT, PRGS, ZWSED, PINPRG, ZPRGS, PFPR=PFPR)
-ENDDO
+ CALL INTERNAL_SEDIM_SPLI(KIB,KIE,KIT,KJB,KJE,KJT, KKB, KKTB, KKTE, KKT, KKL, KRR, &
+ &XSPLIT_MAXCFL, &
+ &PRHODREF, ZW, PDZZ, PPABST, PTHT, PTSTEP, &
+ &6, &
+ &ZRGT, PRGS, PINPRG, ZPRGS, &
+ PFPR=PFPR)
!
!* 2.6 for hail
!
-IF (KRR==7) THEN
- ZREMAINT(:,:) = PTSTEP
- DO WHILE (ANY(ZREMAINT>0.))
- GSEDIM(:,:,:)=.FALSE.
- DO JK = KKTB , KKTE
- GSEDIM(KIB:KIE,KJB:KJE,JK) = &
- (ZRHT(KIB:KIE,KJB:KJE,JK)>XRTMIN(7) .OR. &
- ZPRHS(KIB:KIE,KJB:KJE,JK)>ZRSMIN(7)) .AND. &
- ZREMAINT(KIB:KIE,KJB:KJE)>0.
- ENDDO
- ISEDIM = ICE4_SEDIMENTATION_SPLIT_COUNTJV(GSEDIM(:,:,:),KIT,KJT,KKT,&
- &SIZE(I1),I1(:),I2(:),I3(:))
- CALL INTERNAL_SEDIM_SPLI(KIT, KJT, KKB, KKTB, KKTE, KKT, KKL, KRR, &
- &ISEDIM, GSEDIM, I1, I2, I3, XSPLIT_MAXCFL, ZREMAINT, &
- &PRHODREF, ZW, PDZZ, PPABST, PTHT, PTSTEP, &
- &7, &
- &ZRHT, PRHS, ZWSED, PINPRH, ZPRHS, PFPR=PFPR)
- END DO
+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, &
+ &7, &
+ &ZRHT, PRHS, PINPRH, ZPRHS, &
+ PFPR=PFPR)
ENDIF
!
!
@@ -369,103 +313,127 @@ CONTAINS
!-------------------------------------------------------------------------------
!
!
- SUBROUTINE INTERNAL_SEDIM_SPLI(KIT, KJT, KKB, KKTB, KKTE, KKT, KKL, KRR, &
- &KSEDIM, LDSEDIM, I1, I2, I3, PMAXCFL, PREMAINT, &
- &PRHODREF, POORHODZ, PDZZ, PPABST, PTHT, PTSTEP, &
- &KSPE, &
- &PRXT, PRXS, PWSED, PINPRX, PPRXS, &
- &PRAY, PLBC, PFSEDC, PCONC3D, PSEA, PTOWN, PFPR)
- !
- !* 0. DECLARATIONS
- ! ------------
- !
- USE MODD_RAIN_ICE_DESCR
- USE MODD_RAIN_ICE_PARAM
- !
- IMPLICIT NONE
- !
- !* 0.1 Declarations of dummy arguments :
- !
- INTEGER, INTENT(IN) :: KIT, KJT, KKB, KKTB, KKTE, KKT, KKL, KRR
- INTEGER, INTENT(IN) :: KSEDIM
- LOGICAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: LDSEDIM
- INTEGER, DIMENSION(KSEDIM), INTENT(IN) :: I1, I2, I3
- REAL, INTENT(IN) :: PMAXCFL ! maximum CFL allowed
- REAL, DIMENSION(KIT,KJT), INTENT(INOUT) :: PREMAINT ! Time remaining until the end of the timestep
- REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PRHODREF ! Reference density
- REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: POORHODZ ! One Over (Rhodref times delta Z)
- REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PDZZ ! layer thikness (m)
- REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PPABST
- REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PTHT
- REAL, INTENT(IN) :: PTSTEP ! total timestep
- INTEGER, INTENT(IN) :: KSPE ! 1 for rc, 2 for rr...
- REAL, DIMENSION(KIT,KJT,KKT), INTENT(INOUT) :: PRXT ! mr of specy X
- REAL, DIMENSION(KIT,KJT,KKT), INTENT(INOUT) :: PRXS !Tendency of the specy KSPE
- REAL, DIMENSION(KIT,KJT,0:KKT+1), INTENT(OUT) :: PWSED ! sedimentation flux
- REAL, DIMENSION(KIT,KJT), INTENT(INOUT) :: PINPRX ! instant precip
- REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PPRXS ! external tendencie
- REAL, DIMENSION(KIT,KJT), INTENT(IN), OPTIONAL :: PSEA ! Sea Mask
- REAL, DIMENSION(KIT,KJT), INTENT(IN), OPTIONAL :: PTOWN ! Fraction that is town
- REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN), OPTIONAL :: PRAY, PLBC, PFSEDC, PCONC3D
- REAL, DIMENSION(KIT,KJT,KKT,KRR), OPTIONAL, INTENT(INOUT) :: PFPR ! upper-air precipitation fluxes
- !
- !* 0.2 declaration of local variables
- !
- !
- character(len=10) :: yspe ! String for error message
- INTEGER :: JK, JL, JI, JJ
- REAL :: ZINVTSTEP
- REAL :: ZZWLBDC, ZRAY, ZZT, ZZWLBDA, ZZCC
- REAL :: ZFSED, ZEXSED
- REAL, DIMENSION(KIT, KJT) :: ZMRCHANGE
- REAL, DIMENSION(KIT, KJT) :: ZMAX_TSTEP ! Maximum CFL in column
- REAL, DIMENSION(SIZE(XRTMIN)) :: ZRSMIN
- !
- !-------------------------------------------------------------------------------
- !
- !
- !* 1. Parameters for cloud sedimentation
- !
- !
- !* 2. compute the fluxes
- !
- !
- ZINVTSTEP = 1./PTSTEP
- ZRSMIN(:) = XRTMIN(:) * ZINVTSTEP
- IF(KSPE==2) THEN
- !******* for cloud
- PWSED(:,:,:) = 0.
- DO JL=1, KSEDIM
- JI=I1(JL)
- JJ=I2(JL)
- JK=I3(JL)
- IF(PRXT(JI,JJ,JK)>XRTMIN(KSPE)) THEN
- ZZWLBDC = PLBC(JI,JJ,JK) * PCONC3D(JI,JJ,JK) / &
- (PRHODREF(JI,JJ,JK) * PRXT(JI,JJ,JK))
- ZZWLBDC = ZZWLBDC**XLBEXC
- ZRAY = PRAY(JI,JJ,JK) / ZZWLBDC
- ZZT = PTHT(JI,JJ,JK) * (PPABST(JI,JJ,JK)/XP00)**(XRD/XCPD)
- ZZWLBDA = 6.6E-8*(101325./PPABST(JI,JJ,JK))*(ZZT/293.15)
- ZZCC = XCC*(1.+1.26*ZZWLBDA/ZRAY)
- PWSED(JI, JJ, JK) = PRHODREF(JI,JJ,JK)**(-XCEXVT +1 ) * &
- ZZWLBDC**(-XDC)*ZZCC*PFSEDC(JI,JJ,JK) * PRXT(JI,JJ,JK)
- ENDIF
- ENDDO
- ELSEIF(KSPE==4) THEN
- ! ******* for pristine ice
- PWSED(:,:,:) = 0.
- DO JL=1, KSEDIM
- JI=I1(JL)
- JJ=I2(JL)
- JK=I3(JL)
- IF(PRXT(JI, JJ, JK) .GT. MAX(XRTMIN(4), 1.0E-7)) THEN
- PWSED(JI, JJ, JK) = XFSEDI * PRXT(JI, JJ, JK) * &
- & PRHODREF(JI,JJ,JK)**(1.-XCEXVT) * & ! McF&H
- & MAX( 0.05E6,-0.15319E6-0.021454E6* &
- & ALOG(PRHODREF(JI,JJ,JK)*PRXT(JI,JJ,JK)) )**XEXCSEDI
- ENDIF
- ENDDO
- ELSE
+SUBROUTINE INTERNAL_SEDIM_SPLI(KIB,KIE,KIT,KJB,KJE,KJT,KKB,KKTB,KKTE,KKT,KKL,KRR, &
+ &PMAXCFL,PRHODREF,POORHODZ,PDZZ,PPABST,PTHT,PTSTEP, &
+ &KSPE,PRXT,PRXS,PINPRX,PPRXS, &
+ &PRAY,PLBC,PFSEDC,PCONC3D,PFPR)
+!
+!* 0. DECLARATIONS
+! ------------
+!
+USE MODD_CST, ONLY: XCPD,XP00,XRD
+USE MODD_RAIN_ICE_DESCR, ONLY: XCC,XCEXVT,XDC,XLBEXC,XRTMIN
+USE MODD_RAIN_ICE_PARAM, ONLY: XEXCSEDI,XEXSEDG,XEXSEDH,XEXSEDR,XEXSEDS,XFSEDG,XFSEDH,XFSEDI,XFSEDR,XFSEDS
+!
+IMPLICIT NONE
+!
+!* 0.1 Declarations of dummy arguments :
+!
+INTEGER, INTENT(IN) :: KIB,KIE,KIT, KJB,KJE,KJT, KKB, KKTB, KKTE, KKT, KKL, KRR
+REAL, INTENT(IN) :: PMAXCFL ! maximum CFL allowed
+REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PRHODREF ! Reference density
+REAL, DIMENSION(KIT,KJT,KKTB:KKTE), INTENT(IN) :: POORHODZ ! One Over (Rhodref times delta Z)
+REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PDZZ ! layer thikness (m)
+REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PPABST
+REAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PTHT
+REAL, INTENT(IN) :: PTSTEP ! total timestep
+INTEGER, INTENT(IN) :: KSPE ! 1 for rc, 2 for rr...
+REAL, DIMENSION(KIT,KJT,KKT), INTENT(INOUT) :: PRXT ! mr of specy X
+REAL, DIMENSION(KIT,KJT,KKT), INTENT(INOUT) :: PRXS !Tendency of the specy KSPE
+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
+!
+!* 0.2 declaration of local variables
+!
+character(len=10) :: yspe ! String for error message
+INTEGER :: IDX, ISEDIM
+INTEGER :: JI, JJ, JK, JL
+INTEGER, DIMENSION(KIT*KJT*KKT) :: I1,I2,I3 ! Used to replace the COUNT
+LOGICAL :: GPRESENT_PFPR
+REAL :: ZINVTSTEP
+REAL :: ZZWLBDC, ZRAY, ZZT, ZZWLBDA, ZZCC
+REAL :: ZFSED, ZEXSED
+REAL, DIMENSION(KIT, KJT) :: ZMRCHANGE
+REAL, DIMENSION(KIT, KJT) :: ZMAX_TSTEP ! Maximum CFL in column
+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
+!
+!-------------------------------------------------------------------------------
+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.
+ELSE
+ GPRESENT_PFPR = .FALSE.
+END IF
+!
+PINPRX(:,:) = 0.
+ZINVTSTEP=1./PTSTEP
+ZRSMIN(:) = XRTMIN(:) * ZINVTSTEP
+ZREMAINT(:,:) = PTSTEP
+!
+DO WHILE (ANY(ZREMAINT>0.))
+ ISEDIM = 0
+ DO JK = KKTB,KKTE
+ DO JJ = KJB,KJE
+ DO JI = KIB,KIE
+ IF( (PRXT (JI,JJ,JK)>XRTMIN(KSPE) .OR. &
+ PPRXS(JI,JJ,JK)>ZRSMIN(KSPE)) .AND. &
+ ZREMAINT(JI,JJ)>0. ) THEN
+ ISEDIM = ISEDIM + 1
+ IDX = ISEDIM
+ I1(IDX) = JI
+ I2(IDX) = JJ
+ I3(IDX) = JK
+ END IF
+ END DO
+ END DO
+ END DO
+ !
+ !
+ !* 1. Parameters for cloud sedimentation
+ !
+ !
+ !* 2. compute the fluxes
+ !
+ !
+ IF(KSPE==2) THEN
+ !******* for cloud
+ ZWSED(:,:,:) = 0.
+ DO JL=1, ISEDIM
+ JI=I1(JL)
+ JJ=I2(JL)
+ JK=I3(JL)
+ IF(PRXT(JI,JJ,JK)>XRTMIN(KSPE)) THEN
+ ZZWLBDC = PLBC(JI,JJ,JK) * PCONC3D(JI,JJ,JK) / &
+ (PRHODREF(JI,JJ,JK) * PRXT(JI,JJ,JK))
+ ZZWLBDC = ZZWLBDC**XLBEXC
+ ZRAY = PRAY(JI,JJ,JK) / ZZWLBDC
+ ZZT = PTHT(JI,JJ,JK) * (PPABST(JI,JJ,JK)/XP00)**(XRD/XCPD)
+ ZZWLBDA = 6.6E-8*(101325./PPABST(JI,JJ,JK))*(ZZT/293.15)
+ ZZCC = XCC*(1.+1.26*ZZWLBDA/ZRAY)
+ ZWSED(JI, JJ, JK) = PRHODREF(JI,JJ,JK)**(-XCEXVT +1 ) * &
+ ZZWLBDC**(-XDC)*ZZCC*PFSEDC(JI,JJ,JK) * PRXT(JI,JJ,JK)
+ ENDIF
+ ENDDO
+ ELSEIF(KSPE==4) THEN
+ ! ******* for pristine ice
+ ZWSED(:,:,:) = 0.
+ DO JL=1, ISEDIM
+ JI=I1(JL)
+ JJ=I2(JL)
+ JK=I3(JL)
+ IF(PRXT(JI, JJ, JK) .GT. MAX(XRTMIN(4), 1.0E-7)) THEN
+ ZWSED(JI, JJ, JK) = XFSEDI * PRXT(JI, JJ, JK) * &
+ & PRHODREF(JI,JJ,JK)**(1.-XCEXVT) * & ! McF&H
+ & MAX( 0.05E6,-0.15319E6-0.021454E6* &
+ & ALOG(PRHODREF(JI,JJ,JK)*PRXT(JI,JJ,JK)) )**XEXCSEDI
+ ENDIF
+ ENDDO
+ ELSE
! ******* for other species
SELECT CASE(KSPE)
CASE(3)
@@ -484,73 +452,43 @@ CONTAINS
write( yspe, '( I10 )' ) kspe
call Print_msg( NVERB_FATAL, 'GEN', 'ICE4_SEDIMENTATION_SPLIT', 'no sedimentation parameter for KSPE='//trim(yspe) )
END SELECT
- PWSED(:,:,:) = 0.
- DO JL=1, KSEDIM
- JI=I1(JL)
- JJ=I2(JL)
- JK=I3(JL)
- IF(PRXT(JI,JJ,JK)>XRTMIN(KSPE)) THEN
- PWSED(JI, JJ, JK) = ZFSED * PRXT(JI, JJ, JK)**ZEXSED * &
- PRHODREF(JI, JJ, JK)**(ZEXSED-XCEXVT)
- ENDIF
- ENDDO
- ENDIF
- ZMAX_TSTEP(:,:) = PREMAINT(:,:)
- DO JL=1, KSEDIM
+ !
+ ZWSED(:,:,:) = 0.
+ DO JL=1, ISEDIM
JI=I1(JL)
JJ=I2(JL)
JK=I3(JL)
- IF(PRXT(JI,JJ,JK)>XRTMIN(KSPE) .AND. PWSED(JI, JJ, JK)>1.E-20) THEN
- ZMAX_TSTEP(JI, JJ) = MIN(ZMAX_TSTEP(JI, JJ), PMAXCFL * PRHODREF(JI, JJ, JK) * &
- PRXT(JI, JJ, JK) * PDZZ(JI, JJ, JK) / PWSED(JI, JJ, JK))
+ IF(PRXT(JI,JJ,JK)>XRTMIN(KSPE)) THEN
+ ZWSED(JI, JJ, JK) = ZFSED * PRXT(JI, JJ, JK)**ZEXSED &
+ * PRHODREF(JI, JJ, JK)**(ZEXSED-XCEXVT)
ENDIF
ENDDO
- ZMRCHANGE(:,:) = 0.
- PREMAINT(:,:) = PREMAINT(:,:) - ZMAX_TSTEP(:,:)
+ ENDIF
+ ZMAX_TSTEP(:,:) = ZREMAINT(:,:)
+ DO JL=1, ISEDIM
+ JI=I1(JL)
+ JJ=I2(JL)
+ JK=I3(JL)
+ IF(PRXT(JI,JJ,JK)>XRTMIN(KSPE) .AND. ZWSED(JI, JJ, JK)>1.E-20) THEN
+ ZMAX_TSTEP(JI, JJ) = MIN(ZMAX_TSTEP(JI, JJ), PMAXCFL * PRHODREF(JI, JJ, JK) * &
+ PRXT(JI, JJ, JK) * PDZZ(JI, JJ, JK) / ZWSED(JI, JJ, JK))
+ ENDIF
+ ENDDO
+ ZREMAINT(:,:) = ZREMAINT(:,:) - ZMAX_TSTEP(:,:)
+ DO JK = KKTB , KKTE
+ ZMRCHANGE(:,:) = ZMAX_TSTEP(:,:) * POORHODZ(:,:,JK)*(ZWSED(:,:,JK+KKL)-ZWSED(:,:,JK))
+ PRXT(:,:,JK) = PRXT(:,:,JK) + ZMRCHANGE(:,:) + PPRXS(:,:,JK) * ZMAX_TSTEP(:,:)
+ PRXS(:,:,JK) = PRXS(:,:,JK) + ZMRCHANGE(:,:) * ZINVTSTEP
+ ENDDO
+ PINPRX(:,:) = PINPRX(:,:) + ZWSED(:,:,KKB) / XRHOLW * (ZMAX_TSTEP(:,:) * ZINVTSTEP)
+ IF (GPRESENT_PFPR) THEN
DO JK = KKTB , KKTE
- ZMRCHANGE(:,:) = ZMAX_TSTEP(:,:) * POORHODZ(:,:,JK)*(PWSED(:,:,JK+KKL)-PWSED(:,:,JK))
- PRXT(:,:,JK) = PRXT(:,:,JK) + ZMRCHANGE(:,:) + PPRXS(:,:,JK) * ZMAX_TSTEP(:,:)
- PRXS(:,:,JK) = PRXS(:,:,JK) + ZMRCHANGE(:,:) * ZINVTSTEP
+ PFPR(:,:,JK,KSPE) = PFPR(:,:,JK,KSPE) + ZWSED(:,:,JK) * (ZMAX_TSTEP(:,:) * ZINVTSTEP)
ENDDO
- PINPRX(:,:) = PINPRX(:,:) + ZWSED(:,:,KKB) / XRHOLW * (ZMAX_TSTEP(:,:) * ZINVTSTEP)
- IF (PRESENT(PFPR)) THEN
- DO JK = KKTB , KKTE
- PFPR(:,:,JK,KSPE) = PFPR(:,:,JK,KSPE) + ZWSED(:,:,JK) * (ZMAX_TSTEP(:,:) * ZINVTSTEP)
- ENDDO
- ENDIF
- !
- END SUBROUTINE INTERNAL_SEDIM_SPLI
- !
- FUNCTION ICE4_SEDIMENTATION_SPLIT_COUNTJV(LTAB,KIT,KJT,KKT,KSIZE,I1,I2,I3) RESULT(IC)
- !
- !* 0. DECLARATIONS
- ! ------------
- !
- IMPLICIT NONE
- !
- !* 0.2 declaration of local variables
- !
- INTEGER, INTENT(IN) :: KIT,KJT,KKT,KSIZE
- LOGICAL, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: LTAB ! Mask
- INTEGER, DIMENSION(KSIZE), INTENT(OUT) :: I1,I2,I3 ! Used to replace the COUNT and PACK
- INTEGER :: JI,JJ,JK,IC
- !
- !-------------------------------------------------------------------------------
- !
- IC = 0
- DO JK = 1,SIZE(LTAB,3)
- DO JJ = 1,SIZE(LTAB,2)
- DO JI = 1,SIZE(LTAB,1)
- IF( LTAB(JI,JJ,JK) ) THEN
- IC = IC +1
- I1(IC) = JI
- I2(IC) = JJ
- I3(IC) = JK
- END IF
- END DO
- END DO
- END DO
- !
- END FUNCTION ICE4_SEDIMENTATION_SPLIT_COUNTJV
- !
+ ENDIF
+!
+END DO
+!
+END SUBROUTINE INTERNAL_SEDIM_SPLI
+!
END SUBROUTINE ICE4_SEDIMENTATION_SPLIT