From c8b549123b0b8c8ff3eed0cf4844263c63dd1457 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?S=C3=A9bastien=20Riette?= <sebastien.riette@meteo.fr> Date: Tue, 18 Jan 2022 15:49:38 +0100 Subject: [PATCH] =?UTF-8?q?S=C3=A9bastien=20Riette=2018/01/2022=20ICE3=20a?= =?UTF-8?q?djustment?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit OCND2 included in common code (from AROME version) OCND2 specific code moved in common/micro or arome/dead_code directories High/Low cloud included for condensation (from Méso-NH version) Harmonie SPP code moved outside of physics condensation and ice_adjust optimised by Ryad compute_frac_ice as an elemental procedure --- docs/TODO | 15 +- src/arome/{micro => dead_code}/aro_icecld.F90 | 0 src/arome/{micro => dead_code}/aro_tiwmx.F90 | 0 src/arome/dead_code/mode_qsatmx_tab.F90 | 27 + .../{micro => dead_code}/modi_aro_icecld.F90 | 0 .../{micro => dead_code}/modi_aro_tiwmx.F90 | 0 src/arome/ext/apl_arome.F90 | 10 +- .../micro/externals => ext}/aro_adjust.F90 | 117 +++- src/arome/ext/aro_adjust.h | 59 ++ src/arome/{micro => ext}/modd_spp_type.F90 | 0 src/arome/gmkpack_ignored_files | 11 +- src/arome/micro/condensation.F90 | 497 ----------------- src/arome/micro/ice_adjust.F90 | 390 ------------- src/arome/micro/modi_condensation.F90 | 57 -- src/arome/micro/modi_icecloud.F90 | 26 - src/arome/micro/modi_tiwmx.F90 | 11 - src/arome/micro/rain_ice_old.F90 | 4 +- src/arome/turb/shallow_mf.F90 | 4 +- src/arome/turb/th_r_from_thl_rt_1d.F90 | 3 +- .../micro}/condensation.F90 | 388 +++++++------ .../internals => common/micro}/ice_adjust.F90 | 236 +++++--- src/common/micro/mode_compute_frac_ice.F90 | 91 +-- .../micro/mode_icecloud.F90} | 6 +- .../micro/mode_tiwmx.F90} | 6 +- .../micro/mode_tiwmx_tab.F90} | 31 +- .../micro}/modi_condensation.F90 | 30 +- .../micro/modi_ice_adjust.F90 | 80 +-- src/common/micro/rain_ice.F90 | 18 +- src/mesonh/micro/condensation.f90 | 513 ----------------- src/mesonh/micro/ice_adjust.f90 | 524 ------------------ src/mesonh/turb/shallow_mf.f90 | 6 +- 31 files changed, 694 insertions(+), 2466 deletions(-) rename src/arome/{micro => dead_code}/aro_icecld.F90 (100%) rename src/arome/{micro => dead_code}/aro_tiwmx.F90 (100%) create mode 100644 src/arome/dead_code/mode_qsatmx_tab.F90 rename src/arome/{micro => dead_code}/modi_aro_icecld.F90 (100%) rename src/arome/{micro => dead_code}/modi_aro_tiwmx.F90 (100%) rename src/arome/{modset_Ryad/mpa/micro/externals => ext}/aro_adjust.F90 (79%) create mode 100644 src/arome/ext/aro_adjust.h rename src/arome/{micro => ext}/modd_spp_type.F90 (100%) delete mode 100644 src/arome/micro/condensation.F90 delete mode 100644 src/arome/micro/ice_adjust.F90 delete mode 100644 src/arome/micro/modi_condensation.F90 delete mode 100644 src/arome/micro/modi_icecloud.F90 delete mode 100644 src/arome/micro/modi_tiwmx.F90 rename src/{arome/modset_Ryad/mpa/micro/internals => common/micro}/condensation.F90 (59%) rename src/{arome/modset_Ryad/mpa/micro/internals => common/micro}/ice_adjust.F90 (64%) rename src/{arome/micro/icecloud.F90 => common/micro/mode_icecloud.F90} (98%) rename src/{arome/micro/modd_tiwmx.F90 => common/micro/mode_tiwmx.F90} (97%) rename src/{arome/micro/tiwmx_tab.F90 => common/micro/mode_tiwmx_tab.F90} (85%) rename src/{arome/modset_Ryad/mpa/micro/module => common/micro}/modi_condensation.F90 (73%) rename src/{arome => common}/micro/modi_ice_adjust.F90 (56%) delete mode 100644 src/mesonh/micro/condensation.f90 delete mode 100644 src/mesonh/micro/ice_adjust.f90 diff --git a/docs/TODO b/docs/TODO index a6a30bb9b..55b6e2dc3 100644 --- a/docs/TODO +++ b/docs/TODO @@ -21,12 +21,16 @@ Merge pb: Pas introduits (pour l'instant?) dans version common. Ryad a fait des tests pour regarder impact des allocatable sur CPU => temps * 2 Code à nettoyer quelque soit l'option retenue + Dernier code de Ryad: /home/gmap/mrpm/khatib/public/modset/mods_ice4_nucleation_wrapper.tgz et/ou /home/gmap/mrpm/khatib/public/modset/ice4_nucleation_wrapper.f90 Etape 2: array syntax -> loop - en profiter pour supprimer args PA/PB des routines appelées depuis ice4_tendencies, comme pour nucleation - regarder si pcompute et llcompute sont toujours tous deux nécessaires dans les mode_ice4* avec le passage en do - si possible, modifier ice4_sedimentation_split* dans le même esprit que stat - transformer sedimentation_split_momentum comme sedimentation_split +- une distinction KSIZE/KPROMA est à faire dans la microphysique. Pour une meilleure compréhension du code, il faut partir + du principe de déclarer tous les tableaux en KPROMA et de les utiliser que jusqu'à KSIZE +- essayer l'outil développé par Juan Pb identifiés à corriger plus tard: - deposition devrait être déplacée dans ice4_tendencies @@ -36,12 +40,19 @@ Pb identifiés à corriger plus tard: sont relatifs à la boucle IMICRO et que les calculs faits en debut de routine ne concernent qu'une partie des points => à corriger -- avec découpage en KPROMA, le code ne produit plus les memes resultats - seules les options oper ont été testées, il manque des test pour sedim_after, nmaxiter, xmrstep, xtstep, autoconv, rainfr Répertoire arome/ext contient les codes non PHYEX qu'il faut modifier dans le pack pour qu'il puisse être compilé. -Ce répertoire devra être vidé à la fin du phasage, les modifications nécessaires ayant été fournies par ailleurs +Ce répertoire devra être vidé à la fin du phasage, les modifications nécessaires ayadevront avoir été fournies par ailleurs Budgets/DDH - Le code dans budget_DDH devra être transféré dans mode_budget - les routines arome specifiques aux budgets sont dans mpa/micro, il faudrait les mettre ailleurs + +SPP +- modd_spp_type est pour l'instant dans mpa/micro/externals mais n'est pas de la microphysique + +Gradients/shuman: +- essayer de mettre des abort dans les routines arome (shuman doit suffire) + +compute_frac_ice: prendre la solution retenue à la fin de la discussion entre LAERO et GMAP diff --git a/src/arome/micro/aro_icecld.F90 b/src/arome/dead_code/aro_icecld.F90 similarity index 100% rename from src/arome/micro/aro_icecld.F90 rename to src/arome/dead_code/aro_icecld.F90 diff --git a/src/arome/micro/aro_tiwmx.F90 b/src/arome/dead_code/aro_tiwmx.F90 similarity index 100% rename from src/arome/micro/aro_tiwmx.F90 rename to src/arome/dead_code/aro_tiwmx.F90 diff --git a/src/arome/dead_code/mode_qsatmx_tab.F90 b/src/arome/dead_code/mode_qsatmx_tab.F90 new file mode 100644 index 000000000..01d697b19 --- /dev/null +++ b/src/arome/dead_code/mode_qsatmx_tab.F90 @@ -0,0 +1,27 @@ +MODULE MODE_QSATMX_TAB +IMPLICIT NONE +CONTAINS +FUNCTION QSATMX_TAB(P,T,FICE) + + USE PARKIND1, ONLY : JPRB + USE MODD_CST ,ONLY : XEPSILO + USE MODE_TIWMX, ONLY : ESATI,ESATW + + IMPLICIT NONE + + REAL :: QSATMX_TAB + REAL, INTENT(IN) :: P,T,FICE + + REAL :: ZES + + ZES = ESATI(T)*FICE + ESATW(T)*(1.-FICE) + IF(ZES >= P)THEN ! temp > boiling point, condensation not possible. + ! Then this function lacks physical meaning, + ! here set to one + QSATMX_TAB = 1. + ELSE + QSATMX_TAB = XEPSILO*ZES/(P-ZES) !r + ENDIF + +END FUNCTION QSATMX_TAB +END MODULE MODE_QSATMX_TAB diff --git a/src/arome/micro/modi_aro_icecld.F90 b/src/arome/dead_code/modi_aro_icecld.F90 similarity index 100% rename from src/arome/micro/modi_aro_icecld.F90 rename to src/arome/dead_code/modi_aro_icecld.F90 diff --git a/src/arome/micro/modi_aro_tiwmx.F90 b/src/arome/dead_code/modi_aro_tiwmx.F90 similarity index 100% rename from src/arome/micro/modi_aro_tiwmx.F90 rename to src/arome/dead_code/modi_aro_tiwmx.F90 diff --git a/src/arome/ext/apl_arome.F90 b/src/arome/ext/apl_arome.F90 index 0601ffe03..298f6c90d 100644 --- a/src/arome/ext/apl_arome.F90 +++ b/src/arome/ext/apl_arome.F90 @@ -1645,10 +1645,12 @@ IF (LMICRO) THEN CALL ARO_ADJUST (KLEV,IKU,IKL,KFDIA,KLEV,NRR,& & NGFL_EZDIAG, & - & CFRAC_ICE_ADJUST, LOSUBG_COND, LOSIGMAS, CMICRO, LOCND2, & + & CFRAC_ICE_ADJUST, CCONDENS, CLAMBDA3, LOSUBG_COND, & + & LOSIGMAS, CMICRO, LOCND2, CSUBG_MF_PDF, & & ZDT,VSIGQSAT,ZZZ_F_,& & ZRHODJM__(:,1:KLEV),& & ZEXNREFM_,& + & ZRHODREFM__(:,1:KLEV),& & ZPABSM__(:,1:KLEV),& & ZTHM__(:,1:KLEV),& & ZRM_,ZSIGM_,& @@ -1656,6 +1658,8 @@ IF (LMICRO) THEN & ZRI_MF_,ZCF_MF_,& & ZTHS__(:,1:KLEV),ZRS_,& & ZSRCS__(:,1:KLEV),ZNEBMNH_,& + & ZHLC_HRC__(:,1:KLEV), ZHLC_HCF__(:,1:KLEV),& + & ZHLI_HRI__(:,1:KLEV), ZHLI_HCF__(:,1:KLEV),& & PGP2DSPP,PEZDIAG,& & YDDDH, YDMODEL%YRML_DIAG%YRLDDH, YDMODEL%YRML_DIAG%YRMDDH) @@ -3292,8 +3296,8 @@ IF (LMICRO) THEN CALL ARO_RAIN_ICE (NPROMICRO,KLEV,IKU,IKL,KFDIA,KLEV,NRR,KSTEP+1,NSPLITR,NGFL_EZDIAG,& & LOSUBG_COND, CSUBG_AUCV_RC, CSUBG_AUCV_RI, LOSEDIC,CSEDIM, CMICRO, ZDT,ZDZZ_ ,& & ZRHODJM__(:,1:KLEV),ZRHODREFM__(:,1:KLEV), ZEXNREFM_, ZPABSM__(:,1:KLEV),& - & ZHLC_HRC__(:,1:KLEV), ZHLC_HCF__(KIDIA:KFDIA,1:KLEV),& - & ZHLI_HRI__(:,1:KLEV), ZHLI_HCF__(KIDIA:KFDIA,1:KLEV),& + & ZHLC_HRC__(:,1:KLEV), ZHLC_HCF__(:,1:KLEV),& + & ZHLI_HRI__(:,1:KLEV), ZHLI_HCF__(:,1:KLEV),& & ZTHM__(:,1:KLEV),ZRM_, ZSIGS__(:,1:KLEV), ZNEBMNH_, ZTHS__(:,1:KLEV),ZRS_,& & ZEVAP_, ZCIT_,LOWARM,ZSEA_,ZTOWN_, LOCND2,LGRSN,& & ZINPRR_NOTINCR_, ZINPRS_NOTINCR_, ZINPRG_NOTINCR_, ZINPRH_NOTINCR_,ZPFPR_,& diff --git a/src/arome/modset_Ryad/mpa/micro/externals/aro_adjust.F90 b/src/arome/ext/aro_adjust.F90 similarity index 79% rename from src/arome/modset_Ryad/mpa/micro/externals/aro_adjust.F90 rename to src/arome/ext/aro_adjust.F90 index 637b31c5d..213a974dd 100644 --- a/src/arome/modset_Ryad/mpa/micro/externals/aro_adjust.F90 +++ b/src/arome/ext/aro_adjust.F90 @@ -1,13 +1,14 @@ ! ######spl SUBROUTINE ARO_ADJUST(KKA,KKU,KKL,KLON,KLEV, KRR, & NGFL_EZDIAG, & - HFRAC_ICE, OSUBG_COND, OSIGMAS, & - CMICRO, OCND2, & + HFRAC_ICE, HCONDENS, HLAMBDA3, OSUBG_COND, & + OSIGMAS, CMICRO, OCND2, HSUBG_MF_PDF, & PTSTEP, PSIGQSAT, & - PZZF, PRHODJ, PEXNREF,& + PZZF, PRHODJ, PEXNREF, PRHODREF,& PPABSM, PTHT, PRT, PSIGS, & PMFCONV, PRC_MF, PRI_MF, PCF_MF, & PTHS, PRS, PSRCS, PCLDFR,& + PHLC_HRC, PHLC_HCF, PHLI_HRI, PHLI_HCF,& PGP2DSPP,PEZDIAG, & YDDDH,YDLDDH,YDMDDH) USE PARKIND1, ONLY : JPRB @@ -87,11 +88,10 @@ USE MODD_CONF USE MODD_CST USE MODD_PARAMETERS USE MODD_RAIN_ICE_DESCR -USE MODD_BUDGET +USE MODD_BUDGET, ONLY: TBUDGETDATA, NBUDGET_RI USE MODD_SPP_TYPE ! USE MODI_ICE_ADJUST -USE MODI_BUDGET USE SPP_MOD, ONLY : YSPP_CONFIG,YSPP ! USE DDH_MIX , ONLY : TYP_DDH @@ -113,12 +113,15 @@ INTEGER, INTENT(IN) :: KLEV !Number of vertical levels INTEGER, INTENT(IN) :: KRR ! Number of moist variables INTEGER, INTENT(IN) :: NGFL_EZDIAG ! Diagnostic array dimension CHARACTER*1, INTENT(IN) :: HFRAC_ICE +CHARACTER*80, INTENT(IN) :: HCONDENS +CHARACTER*4, INTENT(IN) :: HLAMBDA3 ! formulation for lambda3 coeff LOGICAL, INTENT(IN) :: OSUBG_COND ! Switch for Subgrid Cond. LOGICAL, INTENT(IN) :: OSIGMAS ! Switch for Sigma_s: ! use values computed in CONDENSATION ! or that from turbulence scheme CHARACTER (LEN=4), INTENT(IN) :: CMICRO ! Microphysics scheme LOGICAL, INTENT(IN) :: OCND2 +CHARACTER*80, INTENT(IN) :: HSUBG_MF_PDF REAL, INTENT(IN) :: PTSTEP ! Time step REAL, INTENT(IN) :: PSIGQSAT ! coeff applied to qsat variance contribution ! @@ -126,6 +129,7 @@ REAL, INTENT(IN) :: PSIGQSAT ! coeff applied to qsat varia REAL, DIMENSION(KLON,1,KLEV), INTENT(IN) :: PZZF ! Height (z) REAL, DIMENSION(KLON,1,KLEV), INTENT(IN) :: PRHODJ !Dry density * Jacobian REAL, DIMENSION(KLON,1,KLEV), INTENT(IN) :: PEXNREF ! Reference Exner function +REAL, DIMENSION(KLON,1,KLEV), INTENT(IN) :: PRHODREF ! ! REAL, DIMENSION(KLON,1,KLEV), INTENT(IN) :: PPABSM ! abs. pressure at time t-dt @@ -144,12 +148,18 @@ REAL, DIMENSION(KLON,1,KLEV), INTENT(OUT) :: PSRCS ! Second-order flux ! s'rc'/2Sigma_s2 at time t+1 ! multiplied by Lambda_3 REAL, DIMENSION(KLON,1,KLEV), INTENT(INOUT) :: PCLDFR! Cloud fraction +! +REAL, DIMENSION(KLON,1,KLEV), INTENT(OUT) :: PHLC_HRC +REAL, DIMENSION(KLON,1,KLEV), INTENT(OUT) :: PHLC_HCF +REAL, DIMENSION(KLON,1,KLEV), INTENT(OUT) :: PHLI_HRI +REAL, DIMENSION(KLON,1,KLEV), INTENT(OUT) :: PHLI_HCF +! REAL, DIMENSION(KLON,YSPP%N2D), TARGET, INTENT(INOUT) :: PGP2DSPP REAL, DIMENSION(KLON,KLEV,NGFL_EZDIAG), INTENT(INOUT) :: PEZDIAG ! -TYPE(TYP_DDH), INTENT(INOUT) :: YDDDH -TYPE(TLDDH), INTENT(IN) :: YDLDDH -TYPE(TMDDH), INTENT(IN) :: YDMDDH +TYPE(TYP_DDH), INTENT(INOUT), TARGET :: YDDDH +TYPE(TLDDH), INTENT(IN), TARGET :: YDLDDH +TYPE(TMDDH), INTENT(IN), TARGET :: YDMDDH ! !* 0.2 Declarations of local variables : @@ -171,6 +181,10 @@ REAL :: ZRATIO ! ZMASSTOT / ZMASSCOR REAL :: ZCOR(KLON) ! for the correction of negative rv ! TYPE(TSPP_CONFIG_MPA) :: YSPP_PSIGQSAT,YSPP_ICE_CLD_WGT +REAL :: ZMU, ZVAL +REAL, DIMENSION(KLON,1) :: ZSIGQSAT, ZICE_CLD_WGT +INTEGER :: JI +TYPE(TBUDGETDATA), DIMENSION(NBUDGET_RI) :: YLBUDGET !NBUDGET_RI is the one with the highest number ! ! REAL(KIND=JPRB) :: ZHOOK_HANDLE @@ -208,6 +222,44 @@ IF ( YSPP_CONFIG%LSPP ) THEN ENDIF +! Compute perturbations +! Perturb PSIGQSAT +IF (YSPP_PSIGQSAT%LPERT) THEN + IF (YSPP_PSIGQSAT%LLNN_MEAN1.OR.YSPP_PSIGQSAT%LLNN_MEAN1_SELF) THEN + ZMU = -0.5 * (YSPP_PSIGQSAT%CMPERT * YSPP_PSIGQSAT%SDEV)**2 + ELSE + ZMU = 0. + ENDIF + DO JI=1,KLON + ZVAL = & + PSIGQSAT*EXP(ZMU+YSPP_PSIGQSAT%CMPERT*YSPP_PSIGQSAT%PGP2DSPP(JI)) + ZSIGQSAT(JI,1) = MAX(YSPP_PSIGQSAT%CLIP(1),MIN(ZVAL,YSPP_PSIGQSAT%CLIP(2))) + ENDDO +ELSE + ZSIGQSAT(:,1) = PSIGQSAT +ENDIF + +! Perturb ICE_CLD_WGT +IF (YSPP_ICE_CLD_WGT%LPERT) THEN + IF (YSPP_ICE_CLD_WGT%LLNN_MEAN1.OR.YSPP_ICE_CLD_WGT%LLNN_MEAN1_SELF) THEN + ZMU = -0.5 * (YSPP_ICE_CLD_WGT%CMPERT * YSPP_ICE_CLD_WGT%SDEV)**2 + ELSE + ZMU = 0. + ENDIF + DO JI=1,KLON + ! Awaiting HARMONIE-AROME physics changes + ! ZVAL = & + ! XFRMIN(21)* EXP(ZMU+YSPP_ICE_CLD_WGT%CMPERT*YSPP_ICE_CLD_WGT%PGP2DSPP(JI)) + ZVAL = & + 1.5* EXP(ZMU+YSPP_ICE_CLD_WGT%CMPERT*YSPP_ICE_CLD_WGT%PGP2DSPP(JI)) + ZICE_CLD_WGT(JI,1) = & + MAX(YSPP_ICE_CLD_WGT%CLIP(1),MIN(ZVAL,YSPP_ICE_CLD_WGT%CLIP(2))) + ENDDO +ELSE +! ZICE_CLD_WGT(:) = XFRMIN(21) + ZICE_CLD_WGT(:,1) = 1.5 +ENDIF + HBUNAME='DEPI' ! !* 2. TRANSFORMATION INTO PHYSICAL TENDENCIES @@ -300,12 +352,19 @@ SELECT CASE ( CMICRO ) ! cloud DO JLON=1,KLON LL(JLON) = (PRS(JLON,1,JLEV,1) <0.) .AND. (PRS(JLON,1,JLEV,2)> 0.) + IF (LL(JLON)) THEN +#ifdef REPRO48 + ZCOR(JLON)=PRS(JLON,1,JLEV,2) +#else + ZCOR(JLON)=MIN(-PRS(JLON,1,JLEV,1),PRS(JLON,1,JLEV,2)) +#endif + ENDIF ENDDO DO JLON=1,KLON IF (LL(JLON)) THEN - PRS(JLON,1,JLEV,1) = PRS(JLON,1,JLEV,1) + PRS(JLON,1,JLEV,2) - PTHS(JLON,1,JLEV) = PTHS(JLON,1,JLEV) - PRS(JLON,1,JLEV,2) * ZLV(JLON) / ZCPH(JLON) / PEXNREF(JLON,1,JLEV) - PRS(JLON,1,JLEV,2) = 0. + PRS(JLON,1,JLEV,1) = PRS(JLON,1,JLEV,1) + ZCOR(JLON) + PTHS(JLON,1,JLEV) = PTHS(JLON,1,JLEV) - ZCOR(JLON) * ZLV(JLON) / ZCPH(JLON) / PEXNREF(JLON,1,JLEV) + PRS(JLON,1,JLEV,2) = PRS(JLON,1,JLEV,2) - ZCOR(JLON) ENDIF ENDDO @@ -343,6 +402,12 @@ END SELECT !IF (LBUDGET_RH) CALL BUDGET (PRS(:,:,:,7) * PRHODJ(:,:,:),12,'NEGA_BU_RRH',YDDDH) !IF (LBUDGET_TH) CALL BUDGET (PTHS(:,:,:) * PRHODJ(:,:,:),4,'NEGA_BU_RTH',YDDDH) +DO JRR=1, NBUDGET_RI + YLBUDGET(JRR)%NBUDGET=JRR + YLBUDGET(JRR)%YDDDH=>YDDDH + YLBUDGET(JRR)%YDLDDH=>YDLDDH + YLBUDGET(JRR)%YDMDDH=>YDMDDH +ENDDO ! !------------------------------------------------------------------------------- @@ -364,11 +429,11 @@ ZZZ(:,:,:) = PZZF(:,:,:) ! IF (KRR==6) THEN CALL ICE_ADJUST ( KKA=KKA,KKU=KKU,KKL=KKL,KRR=KRR,& - & HFRAC_ICE=HFRAC_ICE, HBUNAME=HBUNAME, & + & HFRAC_ICE=HFRAC_ICE, HCONDENS=HCONDENS, HLAMBDA3=HLAMBDA3, HBUNAME=HBUNAME, & & OSUBG_COND=OSUBG_COND, OSIGMAS=OSIGMAS, & - & OCND2=OCND2, & - & PTSTEP=ZTWOTSTEP,PSIGQSAT=PSIGQSAT, & - & PRHODJ=PRHODJ ,PEXNREF=PEXNREF, & + & OCND2=OCND2, HSUBG_MF_PDF=HSUBG_MF_PDF, & + & PTSTEP=ZTWOTSTEP,PSIGQSAT=ZSIGQSAT, & + & PRHODJ=PRHODJ ,PEXNREF=PEXNREF, PRHODREF=PRHODREF, & & PSIGS=PSIGS, PMFCONV=PMFCONV, PPABST=PPABSM, PZZ=ZZZ, & & PEXN=PEXNREF, PCF_MF=PCF_MF,PRC_MF=PRC_MF,PRI_MF=PRI_MF, & & PRV=ZRS(:,:,:,1), PRC=ZRS(:,:,:,2), & @@ -378,16 +443,17 @@ IF (KRR==6) THEN & PRI=ZRS(:,:,:,4), PRIS=PRS(:,:,:,4), & & PRS=ZRS(:,:,:,5), & & PRG=ZRS(:,:,:,6), & - & YSPP_PSIGQSAT=YSPP_PSIGQSAT, & - & YSPP_ICE_CLD_WGT=YSPP_ICE_CLD_WGT, & - & YDDDH=YDDDH, YDLDDH=YDLDDH, YDMDDH=YDMDDH) + & PHLC_HRC=PHLC_HRC(:,:,:), PHLC_HCF=PHLC_HCF(:,:,:), & + & PHLI_HRI=PHLI_HRI(:,:,:), PHLI_HCF=PHLI_HCF(:,:,:), & + & PICE_CLD_WGT=ZICE_CLD_WGT(:,:), & + & TBUDGETS=YLBUDGET, KBUDGETS=SIZE(YLBUDGET)) ELSE CALL ICE_ADJUST ( KKA=KKA,KKU=KKU,KKL=KKL,KRR=KRR,& - & HFRAC_ICE=HFRAC_ICE, HBUNAME=HBUNAME, & + & HFRAC_ICE=HFRAC_ICE, HCONDENS=HCONDENS, HLAMBDA3=HLAMBDA3, HBUNAME=HBUNAME, & & OSUBG_COND=OSUBG_COND, OSIGMAS=OSIGMAS, & - & OCND2=OCND2, & - & PTSTEP=ZTWOTSTEP,PSIGQSAT=PSIGQSAT, & - & PRHODJ=PRHODJ ,PEXNREF=PEXNREF, & + & OCND2=OCND2, HSUBG_MF_PDF=HSUBG_MF_PDF, & + & PTSTEP=ZTWOTSTEP,PSIGQSAT=ZSIGQSAT, & + & PRHODJ=PRHODJ ,PEXNREF=PEXNREF, PRHODREF=PRHODREF, & & PSIGS=PSIGS, PMFCONV=PMFCONV, PPABST=PPABSM, PZZ=ZZZ, & & PEXN=PEXNREF, PCF_MF=PCF_MF,PRC_MF=PRC_MF,PRI_MF=PRI_MF, & & PRV=ZRS(:,:,:,1), PRC=ZRS(:,:,:,2), & @@ -398,9 +464,10 @@ ELSE & PRS=ZRS(:,:,:,5), & & PRG=ZRS(:,:,:,6), & & PRH=ZRS(:,:,:,7), & - & YSPP_PSIGQSAT=YSPP_PSIGQSAT, & - & YSPP_ICE_CLD_WGT=YSPP_ICE_CLD_WGT, & - & YDDDH=YDDDH, YDLDDH=YDLDDH, YDMDDH=YDMDDH) + & PHLC_HRC=PHLC_HRC(:,:,:), PHLC_HCF=PHLC_HCF(:,:,:), & + & PHLI_HRI=PHLI_HRI(:,:,:), PHLI_HCF=PHLI_HCF(:,:,:), & + & PICE_CLD_WGT=ZICE_CLD_WGT(:,:), & + & TBUDGETS=YLBUDGET, KBUDGETS=SIZE(YLBUDGET)) ENDIF CALL CLEAR_SPP_TYPE(YSPP_PSIGQSAT) diff --git a/src/arome/ext/aro_adjust.h b/src/arome/ext/aro_adjust.h new file mode 100644 index 000000000..6f2856602 --- /dev/null +++ b/src/arome/ext/aro_adjust.h @@ -0,0 +1,59 @@ +INTERFACE + SUBROUTINE ARO_ADJUST(KKA,KKU,KKL,KLON,KLEV, KRR,& + & NGFL_EZDIAG, & + & HFRAC_ICE, HCONDENS, HLAMBDA3, OSUBG_COND, & + & OSIGMAS, CMICRO, OCND2, HSUBG_MF_PDF,& + & PTSTEP, PSIGQSAT, PZZF, PRHODJ, PEXNREF, PRHODREF,& + & PPABSM, PTHT, PRT, PSIGS,& + & PMFCONV, PRC_MF, PRI_MF, PCF_MF,& + & PTHS, PRS, PSRCS, PCLDFR, & + & PHLC_HRC, PHLC_HCF, PHLI_HRI, PHLI_HCF, & + & PGP2DSPP,PEZDIAG, & + & YDDDH,YDLDDH,YDMDDH) +USE PARKIND1 ,ONLY : JPIM ,JPRB +USE SPP_MOD, ONLY : YSPP +USE DDH_MIX, ONLY : TYP_DDH +USE YOMLDDH, ONLY : TLDDH +USE YOMMDDH, ONLY : TMDDH +INTEGER(KIND=JPIM), INTENT(IN) :: KKA +INTEGER(KIND=JPIM), INTENT(IN) :: KKU +INTEGER(KIND=JPIM), INTENT(IN) :: KKL +INTEGER(KIND=JPIM), INTENT(IN) :: KLON +INTEGER(KIND=JPIM), INTENT(IN) :: KLEV +INTEGER(KIND=JPIM), INTENT(IN) :: KRR +INTEGER(KIND=JPIM), INTENT(IN) :: NGFL_EZDIAG +CHARACTER*1, INTENT(IN) :: HFRAC_ICE +CHARACTER(LEN=80), INTENT(IN) :: HCONDENS +CHARACTER*4, INTENT(IN) :: HLAMBDA3 +LOGICAL, INTENT(IN) :: OSUBG_COND +LOGICAL, INTENT(IN) :: OSIGMAS +CHARACTER(LEN=4), INTENT(IN) :: CMICRO +LOGICAL, INTENT(IN) :: OCND2 +CHARACTER(LEN=80), INTENT(IN) :: HSUBG_MF_PDF +REAL(KIND=JPRB), INTENT(IN) :: PTSTEP +REAL(KIND=JPRB), INTENT(IN) :: PSIGQSAT +REAL(KIND=JPRB), DIMENSION(KLON,1,KLEV), INTENT(IN) :: PZZF +REAL(KIND=JPRB), DIMENSION(KLON,1,KLEV), INTENT(IN) :: PRHODJ +REAL(KIND=JPRB), DIMENSION(KLON,1,KLEV), INTENT(IN) :: PEXNREF +REAL(KIND=JPRB), DIMENSION(KLON,1,KLEV), INTENT(IN) :: PRHODREF +REAL(KIND=JPRB), DIMENSION(KLON,1,KLEV), INTENT(IN) :: PPABSM +REAL(KIND=JPRB), DIMENSION(KLON,1,KLEV), INTENT(IN) :: PTHT +REAL(KIND=JPRB), DIMENSION(KLON,1,KLEV,KRR), INTENT(INOUT) :: PRT +REAL(KIND=JPRB), DIMENSION(KLON,1,KLEV), INTENT(IN) :: PSIGS +REAL(KIND=JPRB), DIMENSION(KLON,1,KLEV), INTENT(IN) :: PMFCONV +REAL(KIND=JPRB), DIMENSION(KLON,1,KLEV), INTENT(IN) :: PRC_MF,PRI_MF,PCF_MF +REAL(KIND=JPRB), DIMENSION(KLON,1,KLEV), INTENT(INOUT) :: PTHS +REAL(KIND=JPRB), DIMENSION(KLON,1,KLEV,KRR), INTENT(INOUT) :: PRS +REAL(KIND=JPRB), DIMENSION(KLON,1,KLEV), INTENT(OUT) :: PSRCS +REAL(KIND=JPRB), DIMENSION(KLON,1,KLEV), INTENT(INOUT) :: PCLDFR +REAL(KIND=JPRB), DIMENSION(KLON,1,KLEV), INTENT(OUT) :: PHLC_HRC +REAL(KIND=JPRB), DIMENSION(KLON,1,KLEV), INTENT(OUT) :: PHLC_HCF +REAL(KIND=JPRB), DIMENSION(KLON,1,KLEV), INTENT(OUT) :: PHLI_HRI +REAL(KIND=JPRB), DIMENSION(KLON,1,KLEV), INTENT(OUT) :: PHLI_HCF +REAL(KIND=JPRB), DIMENSION(KLON,YSPP%N2D), INTENT(INOUT) :: PGP2DSPP +REAL(KIND=JPRB), DIMENSION(KLON,KLEV,NGFL_EZDIAG), INTENT(INOUT) :: PEZDIAG +TYPE(TYP_DDH) , INTENT(INOUT) :: YDDDH +TYPE(TLDDH) , INTENT(IN) :: YDLDDH +TYPE(TMDDH) , INTENT(IN) :: YDMDDH +END SUBROUTINE ARO_ADJUST +END INTERFACE diff --git a/src/arome/micro/modd_spp_type.F90 b/src/arome/ext/modd_spp_type.F90 similarity index 100% rename from src/arome/micro/modd_spp_type.F90 rename to src/arome/ext/modd_spp_type.F90 diff --git a/src/arome/gmkpack_ignored_files b/src/arome/gmkpack_ignored_files index 7621cb09a..f4e9b065d 100644 --- a/src/arome/gmkpack_ignored_files +++ b/src/arome/gmkpack_ignored_files @@ -91,4 +91,13 @@ phyex/turb/compute_frac_ice3d.F90 phyex/turb/modi_compute_frac_ice.F90 phyex/turb/modi_compute_frac_ice1d.F90 phyex/turb/modi_compute_frac_ice3d.F90 - +phyex/micro/aro_tiwmx.F90 +phyex/micro/modd_tiwmx.F90 +phyex/micro/modi_aro_tiwmx.F90 +phyex/micro/aro_icecld.F90 +phyex/micro/modi_aro_icecld.F90 +phyex/micro/icecloud.F90 +phyex/micro/modi_icecloud.F90 +phyex/micro/tiwmx_tab.F90 +phyex/micro/modi_tiwmx.F90 +phyex/micro/modd_spp_type.F90 diff --git a/src/arome/micro/condensation.F90 b/src/arome/micro/condensation.F90 deleted file mode 100644 index a3033817e..000000000 --- a/src/arome/micro/condensation.F90 +++ /dev/null @@ -1,497 +0,0 @@ -! ######spl - SUBROUTINE CONDENSATION( KIU, KJU, KKU, KIB, KIE, KJB, KJE, KKB, KKE, KKL, & - HFRAC_ICE, & - PPABS, PZZ, PT, PRV, PRC, PRI, PRS, PRG, PSIGS, PMFCONV, PCLDFR, PSIGRC, OUSERI,& - OSIGMAS, OCND2, PSIGQSAT, YSPP_PSIGQSAT,YSPP_ICE_CLD_WGT, & - PLV, PLS, PCPH ) - USE PARKIND1, ONLY : JPRB - USE MODD_SPP_TYPE - USE YOMHOOK , ONLY : LHOOK, DR_HOOK -! ################################################################################ -! -!! -!! PURPOSE -!! ------- -!!** Routine to diagnose cloud fraction, liquid and ice condensate mixing ratios -!! and s'rl'/sigs^2 -!! -!! -!!** METHOD -!! ------ -!! Based on the large-scale fields of temperature, water vapor, and possibly -!! liquid and solid condensate, the conserved quantities r_t and h_l are constructed -!! and then fractional cloudiness, liquid and solid condensate is diagnosed. -!! -!! The total variance is parameterized as the sum of stratiform/turbulent variance -!! and a convective variance. -!! The turbulent variance is parameterized as a function of first-order moments, and -!! the convective variance is modelled as a function of the convective mass flux -!! (units kg/s m^2) as provided by the mass flux convection scheme. -!! -!! Nota: if the host model does not use prognostic values for liquid and solid condensate -!! or does not provide a convective mass flux, put all these values to zero. -!! -!! -!! EXTERNAL -!! -------- -!! INI_CST -!! -!! IMPLICIT ARGUMENTS -!! ------------------ -!! Module MODD_CST : contains physical constants -!! -!! REFERENCE -!! --------- -!! Chaboureau J.P. and P. Bechtold (J. Atmos. Sci. 2002) -!! -!! AUTHOR -!! ------ -!! P. BECHTOLD * Laboratoire d'Aerologie * -!! -!! MODIFICATIONS -!! ------------- -!! Original: 31.1.2002 -!! modified : 21.3.2002 -!! S.Malardel : 05.2006 : Correction sur le calcul de la fonction de -!! Bougeault F2 -!! W. de Rooy: 06-06-2010: Modification in the statistical cloud scheme -!! more specifically adding a variance term -!! following ideas of Lenderink & Siebesma 2002 -!! and adding a height dependence -!! S. Riette, 18 May 2010 : PSIGQSAT is added -!! S. Riette, 11 Oct 2011 : MIN function in PDF for continuity -!! modification of minimum value for Rc+Ri to create cloud and minimum value for sigma -!! Use of guess point as a starting point instead of liquid point -!! Better computation of ZCPH and dRsat/dT -!! Set ZCOND to zero if PCLDFR==0 -!! Safety limitation to .99*Pressure for saturation vapour pressure -!! 2012-02 Y. Seity, add possibility to run with reversed vertical levels -!! 2014-11 K.I Ivarsson add possibility to run with OCND2 option -!! 2016-11 S. Riette: use HFRAC_ICE, output adjusted state -!! 2020-12 U. Andrae : Introduce SPP for HARMONIE-AROME -!------------------------------------------------------------------------------- -! -!* 0. DECLARATIONS -! ------------ -! -USE MODD_CST -USE MODD_PARAMETERS -USE MODE_COMPUTE_FRAC_ICE, ONLY : COMPUTE_FRAC_ICE -USE MODD_TIWMX -USE MODI_ICECLOUD -USE MODI_TIWMX -! -IMPLICIT NONE -! -!* 0.1 Declarations of dummy arguments : -! -! -INTEGER, INTENT(IN) :: KIU ! horizontal dimension in x -INTEGER, INTENT(IN) :: KJU ! horizontal dimension in y -INTEGER, INTENT(IN) :: KKU ! vertical dimension -INTEGER, INTENT(IN) :: KIB ! value of the first point in x -INTEGER, INTENT(IN) :: KIE ! value of the last point in x -INTEGER, INTENT(IN) :: KJB ! value of the first point in y -INTEGER, INTENT(IN) :: KJE ! value of the last point in y -INTEGER, INTENT(IN) :: KKB ! value of the first point in z -INTEGER, INTENT(IN) :: KKE ! value of the last point in z -INTEGER, INTENT(IN) :: KKL ! +1 if grid goes from ground to atmosphere top, -1 otherwise -CHARACTER*1, INTENT(IN) :: HFRAC_ICE -REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PPABS ! pressure (Pa) -REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PZZ ! height of model levels (m) -REAL, DIMENSION(KIU,KJU,KKU), INTENT(INOUT) :: PT ! grid scale T (K) -REAL, DIMENSION(KIU,KJU,KKU), INTENT(INOUT) :: PRV ! grid scale water vapor mixing ratio (kg/kg) -LOGICAL, INTENT(IN) :: OUSERI ! logical switch to compute both - ! liquid and solid condensate (OUSERI=.TRUE.) - ! or only solid condensate (OUSERI=.FALSE.) -LOGICAL, INTENT(IN) :: OCND2 ! logical switch to sparate liquid and ice - ! more rigid (DEFALT value : .FALSE.) -LOGICAL, INTENT(IN) :: OSIGMAS! use present global Sigma_s values - ! or that from turbulence scheme -REAL, INTENT(IN) :: PSIGQSAT ! use an extra "qsat" variance contribution (OSIGMAS case) - ! multiplied by PSIGQSAT -TYPE(TSPP_CONFIG_MPA), INTENT(IN) :: YSPP_PSIGQSAT,YSPP_ICE_CLD_WGT -REAL, DIMENSION(KIU,KJU,KKU), INTENT(INOUT) :: PRC ! grid scale r_c mixing ratio (kg/kg) -REAL, DIMENSION(KIU,KJU,KKU), INTENT(INOUT) :: PRI ! grid scale r_i (kg/kg) -REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PRS ! grid scale mixing ration of snow (kg/kg) -REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PRG ! grid scale mixing ration of graupel (kg/kg) -REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PSIGS ! Sigma_s from turbulence scheme -REAL, DIMENSION(:,:,:), INTENT(IN) :: PMFCONV! convective mass flux (kg /s m^2) -REAL, DIMENSION(KIU,KJU,KKU), INTENT(OUT) :: PCLDFR ! cloud fraction -REAL, DIMENSION(KIU,KJU,KKU), INTENT(OUT) :: PSIGRC ! s r_c / sig_s^2 -REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(IN) :: PLV -REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(IN) :: PLS -REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(IN) :: PCPH -! -! -!* 0.2 Declarations of local variables : -! -INTEGER :: JI, JJ, JK, JKP, JKM, IKTB, IKTE ! loop index -REAL, DIMENSION(KIU,KJU,KKU) :: ZTLK, ZRT ! work arrays for T_l and total water mixing ratio -REAL, DIMENSION(KIU,KJU,KKU) :: ZL ! length scale -REAL, DIMENSION(KIU,KJU,KKU) :: ZFRAC ! Ice fraction -INTEGER, DIMENSION(KIU,KJU) :: ITPL ! top levels of troposphere -REAL, DIMENSION(KIU,KJU) :: ZTMIN ! minimum Temp. related to ITPL -! -REAL, DIMENSION(KIU,KJU,KKU) :: ZLV, ZLS, ZCPD -REAL :: ZTEMP, ZPV, ZQSL, ZPIV, ZQSI, ZCOND, ZLVS ! thermodynamics -REAL :: ZLL, DZZ, ZZZ ! used for length scales -REAL :: ZAH, ZA, ZB, ZSBAR, ZQ1, ZSIGMA, ZDRW, ZDTL, ZSIG_CONV ! related to computation of Sig_s -REAL :: ZRIOLD - -! related to OCND2 noise check : -REAL :: ZRCOLD,ZTSP,ZRSP, ZRSW, ZRFRAC,ZRSDIF ,ZRSB -! related to OCND2 ice cloud calulation : -REAL :: ZDUM1,ZDUM2,ZDUM3,ZDUM4,ZPRIFACT -REAL, DIMENSION(KIU,KJU,KKU) :: TCLD -REAL :: ZDZ(KIB:KIE), ZARDUM(KIE-KIB+1),ZCLDUM(KIE-KIB+1) -REAL, DIMENSION(KIB:KIE) :: ZSIGQSAT,ZICE_CLD_WGT -! end OCND2 - -! -!* 0.3 Definition of constants : -! -!------------------------------------------------------------------------------- -! -REAL :: ZL0 = 600. ! tropospheric length scale -REAL :: ZCSIGMA = 0.2 ! constant in sigma_s parameterization -REAL :: ZCSIG_CONV = 0.30E-2 ! scaling factor for ZSIG_CONV as function of mass flux -! - -REAL, DIMENSION(-22:11) :: ZSRC_1D =(/ & - 0. , 0. , 2.0094444E-04, 0.316670E-03, & - 4.9965648E-04, 0.785956E-03 , 1.2341294E-03, 0.193327E-02, & - 3.0190963E-03, 0.470144E-02 , 7.2950651E-03, 0.112759E-01, & - 1.7350994E-02, 0.265640E-01 , 4.0427860E-02, 0.610997E-01, & - 9.1578111E-02, 0.135888E+00 , 0.1991484 , 0.230756E+00, & - 0.2850565 , 0.375050E+00 , 0.5000000 , 0.691489E+00, & - 0.8413813 , 0.933222E+00 , 0.9772662 , 0.993797E+00, & - 0.9986521 , 0.999768E+00 , 0.9999684 , 0.999997E+00, & - 1.0000000 , 1.000000 /) -INTEGER :: INQ1 -REAL :: ZINC,ZMU,ZVAL -! -!------------------------------------------------------------------------------- -! -! -REAL(KIND=JPRB) :: ZHOOK_HANDLE -IF (LHOOK) CALL DR_HOOK('CONDENSATION',0,ZHOOK_HANDLE) - -IKTB=1+JPVEXT -IKTE=KKU-JPVEXT - -PCLDFR(:,:,:) = 0. ! Initialize values -PSIGRC(:,:,:) = 0. ! Initialize values -ZPRIFACT = 1. ! Initialize value -ZCLDUM=-1. ! Initialize value -IF(OCND2)ZPRIFACT = 0. - -! -! Perturb PSIGQSAT -! - -IF (YSPP_PSIGQSAT%LPERT) THEN - IF (YSPP_PSIGQSAT%LLNN_MEAN1.OR.YSPP_PSIGQSAT%LLNN_MEAN1_SELF) THEN - ZMU = -0.5_JPRB * (YSPP_PSIGQSAT%CMPERT * YSPP_PSIGQSAT%SDEV)**2 - ELSE - ZMU = 0._JPRB - ENDIF - DO JI=KIB,KIE - ZVAL = & - PSIGQSAT*EXP(ZMU+YSPP_PSIGQSAT%CMPERT*YSPP_PSIGQSAT%PGP2DSPP(JI)) - ZSIGQSAT(JI) = MAX(YSPP_PSIGQSAT%CLIP(1),MIN(ZVAL,YSPP_PSIGQSAT%CLIP(2))) - ENDDO -ELSE - ZSIGQSAT(KIB:KIE) = PSIGQSAT -ENDIF - -! -! Perturb ICE_CLD_WGT -! - -IF (YSPP_ICE_CLD_WGT%LPERT) THEN - IF (YSPP_ICE_CLD_WGT%LLNN_MEAN1.OR.YSPP_ICE_CLD_WGT%LLNN_MEAN1_SELF) THEN - ZMU = -0.5_JPRB * (YSPP_ICE_CLD_WGT%CMPERT * YSPP_ICE_CLD_WGT%SDEV)**2 - ELSE - ZMU = 0._JPRB - ENDIF - DO JI=KIB,KIE -! Awaiting HARMONIE-AROME physics changes -! ZVAL = & -! XFRMIN(21)* EXP(ZMU+YSPP_ICE_CLD_WGT%CMPERT*YSPP_ICE_CLD_WGT%PGP2DSPP(JI)) - ZVAL = & - 1.5* EXP(ZMU+YSPP_ICE_CLD_WGT%CMPERT*YSPP_ICE_CLD_WGT%PGP2DSPP(JI)) - ZICE_CLD_WGT(JI) = & - MAX(YSPP_ICE_CLD_WGT%CLIP(1),MIN(ZVAL,YSPP_ICE_CLD_WGT%CLIP(2))) - ENDDO -ELSE -! ZICE_CLD_WGT(KIB:KIE) = XFRMIN(21) - ZICE_CLD_WGT(KIB:KIE) = 1.5 -ENDIF - -! -! -!------------------------------------------------------------------------------- -! store total water mixing ratio -DO JK=IKTB,IKTE - DO JJ=KJB,KJE - DO JI=KIB,KIE - ZRT(JI,JJ,JK) = PRV(JI,JJ,JK) + PRC(JI,JJ,JK) + PRI(JI,JJ,JK)*ZPRIFACT - END DO - END DO -END DO -!------------------------------------------------------------------------------- -! Preliminary calculations -! latent heat of vaporisation/sublimation -IF(PRESENT(PLV) .AND. PRESENT(PLS)) THEN - ZLV(:,:,:)=PLV(:,:,:) - ZLS(:,:,:)=PLS(:,:,:) -ELSE - DO JK=IKTB,IKTE - DO JJ=KJB,KJE - DO JI=KIB,KIE - ZTEMP = PT(JI,JJ,JK) - ! latent heat of vaporisation/sublimation - ZLV(JI,JJ,JK) = XLVTT + ( XCPV - XCL ) * ( ZTEMP - XTT ) - ZLS(JI,JJ,JK) = XLSTT + ( XCPV - XCI ) * ( ZTEMP - XTT ) - ENDDO - ENDDO - ENDDO -ENDIF -IF(PRESENT(PCPH)) THEN - ZCPD(:,:,:)=PCPH(:,:,:) -ELSE - DO JK=IKTB,IKTE - DO JJ=KJB,KJE - DO JI=KIB,KIE - ZCPD(JI,JJ,JK) = XCPD + XCPV*PRV(JI,JJ,JK) + XCL*PRC(JI,JJ,JK) + XCI*PRI(JI,JJ,JK) + & - XCI*(PRS(JI,JJ,JK) + PRG(JI,JJ,JK) ) - ENDDO - ENDDO - ENDDO -ENDIF -! Preliminary calculations needed for computing the "turbulent part" of Sigma_s -IF ( .NOT. OSIGMAS ) THEN - DO JK=IKTB,IKTE - DO JJ=KJB,KJE - DO JI=KIB,KIE - ZTEMP = PT(JI,JJ,JK) - ! store temperature at saturation - ZTLK(JI,JJ,JK) = ZTEMP - ZLV(JI,JJ,JK)*PRC(JI,JJ,JK)/ZCPD(JI,JJ,JK) & - - ZLS(JI,JJ,JK)*PRI(JI,JJ,JK)/ZCPD(JI,JJ,JK)*ZPRIFACT - END DO - END DO - END DO - ! Determine tropopause/inversion height from minimum temperature - ITPL(:,:) = KIB+1 - ZTMIN(:,:) = 400. - DO JK = IKTB+1,IKTE-1 - DO JJ=KJB,KJE - DO JI=KIB,KIE - IF ( PT(JI,JJ,JK) < ZTMIN(JI,JJ) ) THEN - ZTMIN(JI,JJ) = PT(JI,JJ,JK) - ITPL(JI,JJ) = JK - ENDIF - END DO - END DO - END DO - ! Set the mixing length scale - ZL(:,:,KKB) = 20. - DO JK = KKB+KKL,KKE,KKL - DO JJ=KJB,KJE - DO JI=KIB,KIE - ! free troposphere - ZL(JI,JJ,JK) = ZL0 - ZZZ = PZZ(JI,JJ,JK) - PZZ(JI,JJ,KKB) - JKP = ITPL(JI,JJ) - ! approximate length for boundary-layer - IF ( ZL0 > ZZZ ) ZL(JI,JJ,JK) = ZZZ - ! gradual decrease of length-scale near and above tropopause - IF ( ZZZ > 0.9*(PZZ(JI,JJ,JKP)-PZZ(JI,JJ,KKB)) ) & - ZL(JI,JJ,JK) = .6 * ZL(JI,JJ,JK-KKL) - END DO - END DO - END DO -END IF -!------------------------------------------------------------------------------- -! -! -!Ice fraction -ZFRAC(:,:,:) = 0. -WHERE(PRC(:,:,:)+PRI(:,:,:) > 1.E-20) - ZFRAC(:,:,:) = PRI(:,:,:) / (PRC(:,:,:)+PRI(:,:,:)) -ENDWHERE -CALL COMPUTE_FRAC_ICE(SIZE(ZFRAC, 1), SIZE(ZFRAC, 2), SIZE(ZFRAC, 3), HFRAC_ICE, ZFRAC, PT) -IF(.NOT. OUSERI) ZFRAC(:,:,:)=0. -IF(OCND2) ZFRAC=0. -! -DO JK=IKTB,IKTE - JKP=MAX(MIN(JK+KKL,IKTE),IKTB) - JKM=MAX(MIN(JK-KKL,IKTE),IKTB) - DO JJ=KJB,KJE - IF (OCND2) THEN - ZDZ(KIB:KIE) = PZZ(KIB:KIE,JJ,JKP) - PZZ(KIB:KIE,JJ,JKP+1) - CALL ICECLOUD(KIE-KIB+1,PPABS(KIB,JJ,JK),PZZ(KIB,JJ,JK),ZDZ(KIB), & - & PT(KIB,JJ,JK),PRV(KIB,JJ,JK),1.,-1., & - & ZCLDUM,1.,TCLD(KIB,JJ,JK), & - & ZARDUM,ZARDUM,ZARDUM,ZARDUM) - ENDIF - DO JI=KIB,KIE - ! latent heats - ZTEMP = PT(JI,JJ,JK) - ! saturated water vapor mixing ratio over liquid water - IF (OCND2) THEN - ZPV = MIN(ESATW(ZTEMP), .99*PPABS(JI,JJ,JK)) - ELSE - ZPV = MIN(EXP( XALPW - XBETAW / ZTEMP - XGAMW * LOG( ZTEMP ) ), .99*PPABS(JI,JJ,JK)) - ENDIF - ZQSL = XRD / XRV * ZPV / ( PPABS(JI,JJ,JK) - ZPV ) - - ! saturated water vapor mixing ratio over ice - IF (OCND2) THEN - ZPIV = MIN(ESATI(ZTEMP), .99*PPABS(JI,JJ,JK)) - ELSE - ZPIV = MIN(EXP( XALPI - XBETAI / ZTEMP - XGAMI * LOG( ZTEMP ) ), .99*PPABS(JI,JJ,JK)) - ENDIF - ZQSI = XRD / XRV * ZPIV / ( PPABS(JI,JJ,JK) - ZPIV ) - - ! interpolate between liquid and solid as function of temperature - ZQSL = (1. - ZFRAC(JI,JJ,JK)) * ZQSL + ZFRAC(JI,JJ,JK) * ZQSI - ZLVS = (1. - ZFRAC(JI,JJ,JK)) * ZLV(JI,JJ,JK) + & - & ZFRAC(JI,JJ,JK) * ZLS(JI,JJ,JK) - - ! coefficients a and b - ZAH = ZLVS * ZQSL / ( XRV * ZTEMP**2 ) * (XRV * ZQSL / XRD + 1.) - ZA = 1. / ( 1. + ZLVS/ZCPD(JI,JJ,JK) * ZAH ) - ZB = ZAH * ZA - - ZSBAR = ZA * ( ZRT(JI,JJ,JK) - ZQSL + & - & ZAH * ZLVS * (PRC(JI,JJ,JK)+PRI(JI,JJ,JK)*ZPRIFACT) / ZCPD(JI,JJ,JK)) - - ! switch to take either present computed value of SIGMAS - ! or that of Meso-NH turbulence scheme - IF ( OSIGMAS ) THEN - IF (ZSIGQSAT(JI)/=0.) THEN - ZSIGMA = SQRT((2*PSIGS(JI,JJ,JK))**2 + (ZSIGQSAT(JI)*ZQSL*ZA)**2) - ELSE - ZSIGMA = 2*PSIGS(JI,JJ,JK) - END IF - ELSE - ! parameterize Sigma_s with first_order closure - DZZ = PZZ(JI,JJ,JKP) - PZZ(JI,JJ,JKM) - ZDRW = ZRT(JI,JJ,JKP) - ZRT(JI,JJ,JKM) - ZDTL = ZTLK(JI,JJ,JKP) - ZTLK(JI,JJ,JKM) + XG/ZCPD(JI,JJ,JK) * DZZ - ZLL = ZL(JI,JJ,JK) - ! standard deviation due to convection - ZSIG_CONV =0. - IF( SIZE(PMFCONV) /= 0) & - ZSIG_CONV = ZCSIG_CONV * PMFCONV(JI,JJ,JK) / ZA - ! zsigma should be of order 4.e-4 in lowest 5 km of atmosphere - ZSIGMA = SQRT( MAX( 1.E-25, ZCSIGMA * ZCSIGMA * ZLL*ZLL/(DZZ*DZZ)*(& - ZA*ZA*ZDRW*ZDRW - 2.*ZA*ZB*ZDRW*ZDTL + ZB*ZB*ZDTL*ZDTL) + & - ZSIG_CONV * ZSIG_CONV ) ) - END IF - ZSIGMA= MAX( 1.E-10, ZSIGMA ) - - ! normalized saturation deficit - ZQ1 = ZSBAR/ZSIGMA - - ! cloud fraction - PCLDFR(JI,JJ,JK) = MAX( 0., MIN(1.,0.5+0.36*ATAN(1.55*ZQ1)) ) - - ! total condensate - IF (ZQ1 > 0. .AND. ZQ1 <= 2 ) THEN - ZCOND = MIN(EXP(-1.)+.66*ZQ1+.086*ZQ1**2, 2.) ! We use the MIN function for continuity - ELSE IF (ZQ1 > 2.) THEN - ZCOND = ZQ1 - ELSE - ZCOND = EXP( 1.2*ZQ1-1. ) - ENDIF - ZCOND = ZCOND * ZSIGMA - - IF ( ZCOND < 1.E-12 ) THEN - ZCOND = 0. - PCLDFR(JI,JJ,JK) = 0. - ENDIF - IF (PCLDFR(JI,JJ,JK)==0.) THEN - ZCOND=0. - ENDIF - - ZRCOLD=PRC(JI,JJ,JK) - ZRIOLD=PRI(JI,JJ,JK) - - IF(.NOT. OCND2) THEN - PRC(JI,JJ,JK) = (1.-ZFRAC(JI,JJ,JK)) * ZCOND ! liquid condensate - PRI(JI,JJ,JK) = ZFRAC(JI,JJ,JK) * ZCOND ! solid condensate - ELSE - PRC(JI,JJ,JK) = (1.-ZFRAC(JI,JJ,JK)) * ZCOND ! liquid condensate - ! -! This check is mainly for noise reduction : -! ------------------------- - IF(ABS(ZRCOLD-PRC(JI,JJ,JK))>1.0E-12 .AND. ESATW(ZTEMP) < PPABS(JI,JJ,JK)*0.5)THEN - ZRCOLD=PRC(JI,JJ,JK) - ZTSP= TIWMX_TAB(PPABS(JI,JJ,JK),PT(JI,JJ,JK),PRV(JI,JJ,JK),0.,ZRSP,ZRSW,0.1) - - ZRFRAC = PRV(JI,JJ,JK) - ZCOND + ZRCOLD - IF( PRV(JI,JJ,JK) < ZRSW )THEN ! sub - saturation over water: - ! Avoid drying of cloudwater leading to supersaturation with - ! respect to water -! - ZRSDIF= MIN(0.,ZRSP-ZRFRAC) - ELSE ! super - saturation over water: - ! Avoid depostition of water leading to sub-saturation with - ! respect to water - ! ZRSDIF= MAX(0.,ZRSP-ZRFRAC) - ZRSDIF= MAX(0.,ZRSP*PCLDFR(JI,JJ,JK) - ZRFRAC) - ENDIF - ZRSB= ZCOND - ZRSDIF - - PRC(JI,JJ,JK) = ZRSB - ENDIF - ! end check - - ! compute separate ice cloud: - ZZZ = ZDZ(JI) - - ZDUM1 = MIN(1.0,20.* PRC(JI,JJ,JK)*SQRT(ZZZ)/ZQSL) ! clould liquid water - ! factor - - ZDUM2 = (0.8*PCLDFR(JI,JJ,JK)+0.2)*MIN(1.,ZDUM1 + ZDUM4*PCLDFR(JI,JJ,JK)) - ! water cloud, use 'statistical' cloud, but reduce it in case of low liquid content - - ZDUM3 = MAX(0.,TCLD(JI,JJ,JK)-PCLDFR(JI,JJ,JK)) ! pure ice cloud part - - ZDUM4 = PRI(JI,JJ,JK) + PRS(JI,JJ,JK)*0.5 + PRG(JI,JJ,JK)*0.25 - IF (JK==IKTB) ZDUM4 = PRI(JI,JJ,JK) - - ZDUM4 = MAX(0.,MIN(1.,ZICE_CLD_WGT(JI)*ZDUM4*SQRT(ZZZ)/ZQSI)) ! clould ice+solid - ! precip. water factor - PCLDFR(JI,JJ,JK) = MIN(1., ZDUM2 + (0.9*ZDUM3+0.1)*ZDUM4) ! Rad cloud - ! Reduce ice cloud part in case of low ice water content - ENDIF ! End OCND2 - - PT(JI,JJ,JK) = PT(JI,JJ,JK) + ((PRC(JI,JJ,JK)-ZRCOLD)*ZLV(JI,JJ,JK) + & - &(PRI(JI,JJ,JK)-ZRIOLD)*ZLS(JI,JJ,JK) ) & - & /ZCPD(JI,JJ,JK) - PRV(JI,JJ,JK) = ZRT(JI,JJ,JK) - PRC(JI,JJ,JK) - PRI(JI,JJ,JK)*ZPRIFACT - -! s r_c/ sig_s^2 -! PSIGRC(JI,JJ,JK) = PCLDFR(JI,JJ,JK) ! use simple Gaussian relation -! -! multiply PSRCS by the lambda3 coefficient -! -! PSIGRC(JI,JJ,JK) = 2.*PCLDFR(JI,JJ,JK) * MIN( 3. , MAX(1.,1.-ZQ1) ) -! in the 3D case lambda_3 = 1. - INQ1 = MIN( MAX(-22,FLOOR(MIN(100., MAX(-100., 2*ZQ1))) ), 10) !inner min/max prevents sigfpe when 2*zq1 does not fit into an int - ZINC = 2.*ZQ1 - INQ1 - - PSIGRC(JI,JJ,JK) = MIN(1.,(1.-ZINC)*ZSRC_1D(INQ1)+ZINC*ZSRC_1D(INQ1+1)) - - PSIGRC(JI,JJ,JK) = PSIGRC(JI,JJ,JK)* MIN( 3. , MAX(1.,1.-ZQ1) ) - - END DO - END DO -END DO -! -IF (LHOOK) CALL DR_HOOK('CONDENSATION',1,ZHOOK_HANDLE) -END SUBROUTINE CONDENSATION diff --git a/src/arome/micro/ice_adjust.F90 b/src/arome/micro/ice_adjust.F90 deleted file mode 100644 index 8e23954f6..000000000 --- a/src/arome/micro/ice_adjust.F90 +++ /dev/null @@ -1,390 +0,0 @@ -! ########################################################################## - SUBROUTINE ICE_ADJUST (KKA, KKU, KKL, KRR, HFRAC_ICE, & - HBUNAME, OSUBG_COND, OSIGMAS, OCND2, & - PTSTEP, PSIGQSAT, & - PRHODJ, PEXNREF, PSIGS, PMFCONV, PPABST, PZZ, & - PEXN, PCF_MF, PRC_MF, PRI_MF, & - PRV, PRC, PRVS, PRCS, PTH, PTHS, PSRCS, PCLDFR, & - PRR, PRI, PRIS, PRS, PRG, & - PRH, POUT_RV, POUT_RC, POUT_RI, POUT_TH, & - YSPP_PSIGQSAT, YSPP_ICE_CLD_WGT, & - YDDDH, YDLDDH, YDMDDH) - USE PARKIND1, ONLY : JPRB - USE MODD_SPP_TYPE - USE YOMHOOK , ONLY : LHOOK, DR_HOOK - USE DDH_MIX , ONLY : TYP_DDH - USE YOMLDDH , ONLY : TLDDH - USE YOMMDDH , ONLY : TMDDH -! ######################################################################### -! -!!**** *ICE_ADJUST* - compute the ajustment of water vapor in mixed-phase -!! clouds -!! -!! PURPOSE -!! ------- -!! The purpose of this routine is to compute the fast microphysical sources -!! through a saturation ajustement procedure in case of mixed-phase clouds. -!! -!! -!!** METHOD -!! ------ -!! Langlois, Tellus, 1973 for the cloudless version. -!! When cloud water is taken into account, refer to book 1 of the -!! documentation. -!! -!! -!! -!! EXTERNAL -!! -------- -!! None -!! -!! -!! IMPLICIT ARGUMENTS -!! ------------------ -!! Module MODD_CST -!! XP00 ! Reference pressure -!! XMD,XMV ! Molar mass of dry air and molar mass of vapor -!! XRD,XRV ! Gaz constant for dry air, gaz constant for vapor -!! XCPD,XCPV ! Cpd (dry air), Cpv (vapor) -!! XCL ! Cl (liquid) -!! XCI ! Ci (ice) -!! XTT ! Triple point temperature -!! XLVTT ! Vaporization heat constant -!! XLSTT ! Sublimation heat constant -!! XALPW,XBETAW,XGAMW ! Constants for saturation vapor over liquid -!! ! pressure function -!! XALPI,XBETAI,XGAMI ! Constants for saturation vapor over ice -!! ! pressure function -!! Module MODD_CONF -!! CCONF -!! Module MODD_BUDGET: -!! NBUMOD -!! CBUTYPE -!! NBUPROCCTR -!! LBU_RTH -!! LBU_RRV -!! LBU_RRC -!! LBU_RRI -!! -!! -!! REFERENCE -!! --------- -!! Book 1 and Book2 of documentation ( routine ICE_ADJUST ) -!! Langlois, Tellus, 1973 -!! -!! AUTHOR -!! ------ -!! J.-P. Pinty * Laboratoire d'Aerologie* -!! -!! -!! MODIFICATIONS -!! ------------- -!! Original 06/12/96 -!! M. Tomasini 27/11/00 Change CND and DEP fct of the T instead of rc and ri -!! Avoid the sub- and super-saturation before the ajustment -!! Avoid rc>0 if T<T00 before the ajustment -!! P Bechtold 12/02/02 change subgrid condensation -!! JP Pinty 29/11/02 add ICE2 and IC4 cases -!! (P. Jabouille) 27/05/04 safety test for case where esw/i(T)> pabs (~Z>40km) -!! J.Pergaud and S.Malardel Add EDKF case -!! S. Riette ice for EDKF -!! 2012-02 Y. Seity, add possibility to run with reversed vertical levels -!! 2016-07 S. Riette: adjustement is now realized on state variables (PRV, PRC, PRI, PTH) -!! whereas tendencies are still applied on S variables. -!! This modification allows to call ice_adjust on T variable -!! or to call it on S variables -!! 2016-11 S. Riette: all-or-nothing adjustment now uses condensation -!! 2020-12 U. Andrae : Introduce SPP for HARMONIE-AROME -!! -!------------------------------------------------------------------------------- -! -!* 0. DECLARATIONS -! ------------ -! -USE MODD_PARAMETERS -USE MODD_CST -USE MODD_CONF -USE MODD_BUDGET -! -USE MODI_CONDENSATION -USE MODI_BUDGET_DDH -USE MODE_FMWRIT -! -IMPLICIT NONE -! -! -!* 0.1 Declarations of dummy arguments : -! -! -INTEGER, INTENT(IN) :: KKA !near ground array index -INTEGER, INTENT(IN) :: KKU !uppest atmosphere array index -INTEGER, INTENT(IN) :: KKL !vert. levels type 1=MNH -1=ARO -INTEGER, INTENT(IN) :: KRR ! Number of moist variables -CHARACTER*1, INTENT(IN) :: HFRAC_ICE -CHARACTER*4, INTENT(IN) :: HBUNAME ! Name of the budget -LOGICAL, INTENT(IN) :: OSUBG_COND ! Switch for Subgrid - ! Condensation -LOGICAL :: OSIGMAS ! Switch for Sigma_s: - ! use values computed in CONDENSATION - ! or that from turbulence scheme -LOGICAL :: OCND2 ! logical switch to sparate liquid - ! and ice - ! more rigid (DEFALT value : .FALSE.) -REAL, INTENT(IN) :: PTSTEP ! Double Time step - ! (single if cold start) -REAL, INTENT(IN) :: PSIGQSAT ! coeff applied to qsat variance contribution -! -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRHODJ ! Dry density * Jacobian -REAL, DIMENSION(:,:,:), INTENT(IN) :: PEXNREF ! Reference Exner function -! -REAL, DIMENSION(:,:,:), INTENT(IN) :: PSIGS ! Sigma_s at time t -REAL, DIMENSION(:,:,:), INTENT(IN) :: PMFCONV ! convective mass flux -REAL, DIMENSION(:,:,:), INTENT(IN) :: PPABST ! Absolute Pressure at t -REAL, DIMENSION(:,:,:), INTENT(IN) :: PZZ ! height of model layer -REAL, DIMENSION(:,:,:), INTENT(IN) :: PEXN ! Exner function -! -REAL, DIMENSION(:,:,:), INTENT(IN) :: PCF_MF! Convective Mass Flux Cloud fraction -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRC_MF! Convective Mass Flux liquid mixing ratio -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRI_MF! Convective Mass Flux ice mixing ratio -! -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRV ! Water vapor m.r. to adjust -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRC ! Cloud water m.r. to adjust -REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PRVS ! Water vapor m.r. source -REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PRCS ! Cloud water m.r. source -REAL, DIMENSION(:,:,:), INTENT(IN) :: PTH ! Theta to adjust -REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PTHS ! Theta source -REAL, DIMENSION(:,:,:), INTENT(OUT) :: PSRCS ! Second-order flux - ! s'rc'/2Sigma_s2 at time t+1 - ! multiplied by Lambda_3 -REAL, DIMENSION(:,:,:), INTENT(OUT) :: PCLDFR ! Cloud fraction -! -REAL, DIMENSION(:,:,:), INTENT(INOUT):: PRIS ! Cloud ice m.r. at t+1 -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRR ! Rain water m.r. to adjust -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRI ! Cloud ice m.r. to adjust -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRS ! Aggregate m.r. to adjust -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRG ! Graupel m.r. to adjust -REAL, DIMENSION(:,:,:), OPTIONAL, INTENT(IN) :: PRH ! Hail m.r. to adjust -REAL, DIMENSION(:,:,:), OPTIONAL, INTENT(OUT) :: POUT_RV ! Adjusted value -REAL, DIMENSION(:,:,:), OPTIONAL, INTENT(OUT) :: POUT_RC ! Adjusted value -REAL, DIMENSION(:,:,:), OPTIONAL, INTENT(OUT) :: POUT_RI ! Adjusted value -REAL, DIMENSION(:,:,:), OPTIONAL, INTENT(OUT) :: POUT_TH ! Adjusted value -TYPE(TSPP_CONFIG_MPA), INTENT(IN) :: YSPP_PSIGQSAT,YSPP_ICE_CLD_WGT -TYPE(TYP_DDH) , INTENT(INOUT) :: YDDDH -TYPE(TLDDH) , INTENT(IN) :: YDLDDH -TYPE(TMDDH) , INTENT(IN) :: YDMDDH -! -!* 0.2 Declarations of local variables : -! -! -REAL :: ZEPS ! Mv/Md -REAL :: ZT00,ZT0 ! Min and max temperature for the mixed phase liquid and solid water - ! for the coeff CND of the barycentric mixing ratio -REAL, DIMENSION(SIZE(PEXNREF,1),SIZE(PEXNREF,2),SIZE(PEXNREF,3)) & - :: ZT, & ! adjusted temperature - ZRV, ZRC, ZRI, & ! adjusted state - ZCPH, & ! guess of the CPh for the mixing - ZLV, & ! guess of the Lv at t+1 - ZLS, & ! guess of the Ls at t+1 - ZW1,ZW2, & ! Work arrays for intermediate fields - ZH2O ! Fraction of all types of water (kg/kg) (for OCND2) -! -INTEGER :: IIU,IJU,IKU! dimensions of dummy arrays -INTEGER :: IIB,IJB ! Horz index values of the first inner mass points -INTEGER :: IIE,IJE ! Horz index values of the last inner mass points -INTEGER :: IKB ! K index value of the first inner mass point -INTEGER :: IKE ! K index value of the last inner mass point -INTEGER :: JITER,ITERMAX ! iterative loop for first order adjustment -! -REAL, DIMENSION(SIZE(PEXNREF,1),SIZE(PEXNREF,2),SIZE(PEXNREF,3)) :: ZSIGS -! -!------------------------------------------------------------------------------- -! -!* 1. PRELIMINARIES -! ------------- -! -REAL(KIND=JPRB) :: ZHOOK_HANDLE -IF (LHOOK) CALL DR_HOOK('ICE_ADJUST',0,ZHOOK_HANDLE) -IIU = SIZE(PEXNREF,1) -IJU = SIZE(PEXNREF,2) -IKU = SIZE(PEXNREF,3) -IIB = 1 + JPHEXT -IIE = IIU - JPHEXT -IJB = 1 + JPHEXT -IJE = IJU - JPHEXT -IKB=KKA+JPVEXT*KKL -IKE=KKU-JPVEXT*KKL -! -ITERMAX=1 -! -!------------------------------------------------------------------------------- -! -!* 2. COMPUTE QUANTITIES WITH THE GUESS OF THE FUTURE INSTANT -! ------------------------------------------------------- -! -! -! beginning of the iterative loop (to compute the adjusted state) -ZRV(:,:,:)=PRV(:,:,:) -ZRC(:,:,:)=PRC(:,:,:) -ZRI(:,:,:)=PRI(:,:,:) -ZT(:,:,:)=PTH(:,:,:) * PEXN(:,:,:) -! -DO JITER =1,ITERMAX - ! - !* 2.3 compute the latent heat of vaporization Lv(T*) at t+1 - ! and the latent heat of sublimation Ls(T*) at t+1 - ! - ZLV(:,:,:) = XLVTT + ( XCPV - XCL ) * ( ZT(:,:,:) -XTT ) - ZLS(:,:,:) = XLSTT + ( XCPV - XCI ) * ( ZT(:,:,:) -XTT ) - ! - !* 2.4 compute the specific heat for moist air (Cph) at t+1 - ! - IF ( KRR == 7 ) THEN - ZCPH(:,:,:) = XCPD + XCPV * ZRV(:,:,:) & - + XCL * (ZRC(:,:,:) + PRR(:,:,:)) & - + XCI * (ZRI(:,:,:) + PRS(:,:,:) + PRG(:,:,:) + PRH(:,:,:)) - ELSE IF( KRR == 6 ) THEN - ZCPH(:,:,:) = XCPD + XCPV * ZRV(:,:,:) & - + XCL * (ZRC(:,:,:) + PRR(:,:,:)) & - + XCI * (ZRI(:,:,:) + PRS(:,:,:) + PRG(:,:,:)) - ELSE IF( KRR == 5 ) THEN - ZCPH(:,:,:) = XCPD + XCPV * ZRV(:,:,:) & - + XCL * (ZRC(:,:,:) + PRR(:,:,:)) & - + XCI * (ZRI(:,:,:) + PRS(:,:,:)) - ELSE IF( KRR == 3 ) THEN - ZCPH(:,:,:) = XCPD + XCPV * ZRV(:,:,:) & - + XCL * (ZRC(:,:,:) + PRR(:,:,:)) - ELSE IF( KRR == 2 ) THEN - ZCPH(:,:,:) = XCPD + XCPV * ZRV(:,:,:) & - + XCL * ZRC(:,:,:) - END IF - ! - IF ( OSUBG_COND ) THEN - ! - !* 3. SUBGRID CONDENSATION SCHEME - ! --------------------------- - ! - ! PSRC= s'rci'/Sigma_s^2 - ! ZT, ZRV, ZRC and ZRI are INOUT -! ZW8 = 4.*ZW2 ! + 0.5_JPRB*PRSS*PTSTEP + 0.25_JPRB*PRGS*PTSTEP ! include snow and graupel to cloudcover -! ZW8(:,:,1) = ZW2(:,:,1) ! but don't do this for the lowest model -! ZW8(:,:,2) = ZW2(:,:,2) ! level (Hail not concidered yet) - - CALL CONDENSATION(IIU, IJU, IKU, IIB, IIE, IJB, IJE, IKB, IKE, KKL, & - HFRAC_ICE, & - PPABST, PZZ, ZT, ZRV, ZRC, ZRI, PRS, PRG, PSIGS, PMFCONV, PCLDFR, & - PSRCS, .TRUE., OSIGMAS, & - OCND2, PSIGQSAT, YSPP_PSIGQSAT,YSPP_ICE_CLD_WGT, & - PLV=ZLV, PLS=ZLS, PCPH=ZCPH) - ELSE - ! - !* 4. ALL OR NOTHING CONDENSATION SCHEME - ! FOR MIXED-PHASE CLOUD - ! ----------------------------------------------- - ZSIGS=0. - CALL CONDENSATION(IIU, IJU, IKU, IIB, IIE, IJB, IJE, IKB, IKE, KKL, & - HFRAC_ICE, & - PPABST, PZZ, ZT, ZRV, ZRC, ZRI, PRS, PRG, ZSIGS, PMFCONV, PCLDFR, & - PSRCS, .TRUE., OSIGMAS=.TRUE., & - OCND2=OCND2, PSIGQSAT=0., & - YSPP_PSIGQSAT=YSPP_PSIGQSAT, & - YSPP_ICE_CLD_WGT=YSPP_ICE_CLD_WGT, & - PLV=ZLV, PLS=ZLS, PCPH=ZCPH) - ENDIF -ENDDO ! end of the iterative loop -! -!* 5. COMPUTE THE SOURCES AND STORES THE CLOUD FRACTION -! ------------------------------------------------- -! -! -!* 5.0 compute the variation of mixing ratio -! - ! Rc - Rc* -ZW1(:,:,:) = (ZRC(:,:,:) - PRC(:,:,:)) / PTSTEP ! Pcon = ---------- - ! 2 Delta t - -ZW2(:,:,:) = (ZRI(:,:,:) - PRI(:,:,:)) / PTSTEP ! idem ZW1 but for Ri -! -!* 5.1 compute the sources -! -WHERE( ZW1(:,:,:) < 0.0 ) - ZW1(:,:,:) = MAX ( ZW1(:,:,:), -PRCS(:,:,:) ) -ELSEWHERE - ZW1(:,:,:) = MIN ( ZW1(:,:,:), PRVS(:,:,:) ) -END WHERE -PRVS(:,:,:) = PRVS(:,:,:) - ZW1(:,:,:) -PRCS(:,:,:) = PRCS(:,:,:) + ZW1(:,:,:) -PTHS(:,:,:) = PTHS(:,:,:) + & - ZW1(:,:,:) * ZLV(:,:,:) / (ZCPH(:,:,:) * PEXNREF(:,:,:)) -! -WHERE( ZW2(:,:,:) < 0.0 ) - ZW2(:,:,:) = MAX ( ZW2(:,:,:), -PRIS(:,:,:) ) -ELSEWHERE - ZW2(:,:,:) = MIN ( ZW2(:,:,:), PRVS(:,:,:) ) -END WHERE -PRVS(:,:,:) = PRVS(:,:,:) - ZW2(:,:,:) -PRIS(:,:,:) = PRIS(:,:,:) + ZW2(:,:,:) -PTHS(:,:,:) = PTHS(:,:,:) + & - ZW2(:,:,:) * ZLS(:,:,:) / (ZCPH(:,:,:) * PEXNREF(:,:,:)) -! -! -!* 5.2 compute the cloud fraction PCLDFR -! -IF ( .NOT. OSUBG_COND ) THEN - WHERE (PRCS(:,:,:) + PRIS(:,:,:) > 1.E-12 / PTSTEP) - PCLDFR(:,:,:) = 1. - ELSEWHERE - PCLDFR(:,:,:) = 0. - ENDWHERE - IF ( SIZE(PSRCS,3) /= 0 ) THEN - PSRCS(:,:,:) = PCLDFR(:,:,:) - END IF -ELSE - !We limit PRC_MF+PRI_MF to PRVS*PTSTEP to avoid negative humidity - ZW1(:,:,:)=PRC_MF(:,:,:)/PTSTEP - ZW2(:,:,:)=PRI_MF(:,:,:)/PTSTEP - WHERE(ZW1(:,:,:)+ZW2(:,:,:)>PRVS(:,:,:)) - ZW1(:,:,:)=ZW1(:,:,:)*PRVS(:,:,:)/(ZW1(:,:,:)+ZW2(:,:,:)) - ZW2(:,:,:)=PRVS(:,:,:)-ZW1(:,:,:) - ENDWHERE - PCLDFR(:,:,:)=MIN(1.,PCLDFR(:,:,:)+PCF_MF(:,:,:)) - PRCS(:,:,:)=PRCS(:,:,:)+ZW1(:,:,:) - PRIS(:,:,:)=PRIS(:,:,:)+ZW2(:,:,:) - PRVS(:,:,:)=PRVS(:,:,:)-(ZW1(:,:,:)+ZW2(:,:,:)) - PTHS(:,:,:) = PTHS(:,:,:) + & - (ZW1 * ZLV(:,:,:) + ZW2 * ZLS(:,:,:)) / ZCPH(:,:,:) & - / PEXNREF(:,:,:) - IF(PRESENT(POUT_RV) .OR. PRESENT(POUT_RC) .OR. & - &PRESENT(POUT_RI) .OR. PRESENT(POUT_TH)) THEN - ZW1(:,:,:)=PRC_MF(:,:,:) - ZW2(:,:,:)=PRI_MF(:,:,:) - WHERE(ZW1(:,:,:)+ZW2(:,:,:)>ZRV(:,:,:)) - ZW1(:,:,:)=ZW1(:,:,:)*ZRV(:,:,:)/(ZW1(:,:,:)+ZW2(:,:,:)) - ZW2(:,:,:)=ZRV(:,:,:)-ZW1(:,:,:) - ENDWHERE - ZRC(:,:,:)=ZRC(:,:,:)+ZW1(:,:,:) - ZRI(:,:,:)=ZRI(:,:,:)+ZW2(:,:,:) - ZRV(:,:,:)=ZRV(:,:,:)-(ZW1(:,:,:)+ZW2(:,:,:)) - ZT(:,:,:) = ZT(:,:,:) + & - (ZW1 * ZLV(:,:,:) + ZW2 * ZLS(:,:,:)) / ZCPH(:,:,:) - IF(PRESENT(POUT_RV)) POUT_RV=ZRV - IF(PRESENT(POUT_RC)) POUT_RC=ZRC - IF(PRESENT(POUT_RI)) POUT_RI=ZRI - IF(PRESENT(POUT_TH)) POUT_TH=ZT / PEXN(:,:,:) - ENDIF -ENDIF -! -! -! -!* 6. STORE THE BUDGET TERMS -! ---------------------- -! -IF (LBUDGET_RV) CALL BUDGET_DDH (PRVS(:,:,:) * PRHODJ(:,:,:),6,HBUNAME//'_BU_RRV',YDDDH, YDLDDH, YDMDDH) -IF (LBUDGET_RC) CALL BUDGET_DDH (PRCS(:,:,:) * PRHODJ(:,:,:),7,HBUNAME//'_BU_RRC',YDDDH, YDLDDH, YDMDDH) -IF (LBUDGET_RI) CALL BUDGET_DDH (PRIS(:,:,:) * PRHODJ(:,:,:),9,HBUNAME//'_BU_RRI',YDDDH, YDLDDH, YDMDDH) -IF (LBUDGET_TH) CALL BUDGET_DDH (PTHS(:,:,:) * PRHODJ(:,:,:),4,HBUNAME//'_BU_RTH',YDDDH, YDLDDH, YDMDDH) -! -!------------------------------------------------------------------------------ -! -! -IF (LHOOK) CALL DR_HOOK('ICE_ADJUST',1,ZHOOK_HANDLE) -END SUBROUTINE ICE_ADJUST diff --git a/src/arome/micro/modi_condensation.F90 b/src/arome/micro/modi_condensation.F90 deleted file mode 100644 index fa0a449ce..000000000 --- a/src/arome/micro/modi_condensation.F90 +++ /dev/null @@ -1,57 +0,0 @@ -! ######spl - MODULE MODI_CONDENSATION -! ######################## -! -INTERFACE -! - SUBROUTINE CONDENSATION( KIU, KJU, KKU, KIB, KIE, KJB, KJE, KKB, KKE, KKL,& - HFRAC_ICE, & - PPABS, PZZ, PT, PRV, PRC, PRI, PRS, PRG, PSIGS, PMFCONV, PCLDFR, PSIGRC, OUSERI,& - OSIGMAS, OCND2, PSIGQSAT, & - YSPP_PSIGQSAT, YSPP_ICE_CLD_WGT,& - PLV, PLS, PCPH) -! -USE MODD_SPP_TYPE -! -INTEGER, INTENT(IN) :: KIU ! horizontal dimension in x -INTEGER, INTENT(IN) :: KJU ! horizontal dimension in y -INTEGER, INTENT(IN) :: KKU ! vertical dimension -INTEGER, INTENT(IN) :: KIB ! value of the first point in x -INTEGER, INTENT(IN) :: KIE ! value of the last point in x -INTEGER, INTENT(IN) :: KJB ! value of the first point in y -INTEGER, INTENT(IN) :: KJE ! value of the last point in y -INTEGER, INTENT(IN) :: KKB ! value of the first point in z -INTEGER, INTENT(IN) :: KKE ! value of the last point in z -INTEGER, INTENT(IN) :: KKL ! +1 if grid goes from ground to atmosphere top, -1 otherwise -CHARACTER*1, INTENT(IN) :: HFRAC_ICE -REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PPABS ! pressure (Pa) -REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PZZ ! height of model levels (m) -REAL, DIMENSION(KIU,KJU,KKU), INTENT(INOUT) :: PT ! grid scale T (K) -REAL, DIMENSION(KIU,KJU,KKU), INTENT(INOUT) :: PRV ! grid scale water vapor mixing ratio (kg/kg) -LOGICAL, INTENT(IN) :: OUSERI ! logical switch to compute both - ! liquid and solid condensate (OUSERI=.TRUE.) - ! or only solid condensate (OUSERI=.FALSE.) -LOGICAL, INTENT(IN) :: OSIGMAS! use present global Sigma_s values - ! or that from turbulence scheme -LOGICAL, INTENT(IN) :: OCND2 ! logical switch to sparate liquid and ice - ! more rigid (DEFALT value : .FALSE.) - -REAL, INTENT(IN) :: PSIGQSAT ! use an extra "qsat" variance contribution (OSIGMAS case) - ! multiplied by PSIGQSAT -TYPE(TSPP_CONFIG_MPA), INTENT(IN) :: YSPP_PSIGQSAT,YSPP_ICE_CLD_WGT -REAL, DIMENSION(KIU,KJU,KKU), INTENT(INOUT) :: PRC ! grid scale r_c mixing ratio (kg/kg) -REAL, DIMENSION(KIU,KJU,KKU), INTENT(INOUT) :: PRI ! grid scale r_i (kg/kg) -REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PRS ! grid scale mixing ration of snow (kg/kg) -REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PRG ! grid scale mixing ration of graupel (kg/kg) -REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PSIGS ! Sigma_s from turbulence scheme -REAL, DIMENSION(:,:,:), INTENT(IN) :: PMFCONV! convective mass flux (kg /s m^2) -REAL, DIMENSION(KIU,KJU,KKU), INTENT(OUT) :: PCLDFR ! cloud fraction -REAL, DIMENSION(KIU,KJU,KKU), INTENT(OUT) :: PSIGRC ! s r_c / sig_s^2 -REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(IN) :: PLV -REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(IN) :: PLS -REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(IN) :: PCPH -END SUBROUTINE CONDENSATION -! -END INTERFACE -! -END MODULE MODI_CONDENSATION diff --git a/src/arome/micro/modi_icecloud.F90 b/src/arome/micro/modi_icecloud.F90 deleted file mode 100644 index eede63c72..000000000 --- a/src/arome/micro/modi_icecloud.F90 +++ /dev/null @@ -1,26 +0,0 @@ -! ######spl - MODULE MODI_ICECLOUD -! #################### -! -INTERFACE -SUBROUTINE ICECLOUD & - & ( NP,PP,PZ,PDZ,PT,PR,PTSTEP,PPBLH,PWCLD,XW2D, & - & SIFRC,SSIO,SSIU,W2D,RSI) -INTEGER, INTENT(IN) :: NP -REAL, INTENT(IN) :: PP(NP) -REAL, INTENT(IN) :: PZ(NP) -REAL, INTENT(IN) :: PDZ(NP) -REAL, INTENT(IN) :: PT(NP) -REAL, INTENT(IN) :: PR(NP) -REAL, INTENT(IN) :: PTSTEP -REAL, INTENT(IN) :: PPBLH -REAL, INTENT(IN) :: PWCLD(NP) -REAL, INTENT(IN) :: XW2D -REAL, INTENT(OUT) :: SIFRC(NP) -REAL, INTENT(OUT) :: SSIO(NP) -REAL, INTENT(OUT) :: SSIU(NP) -REAL, INTENT(OUT) :: W2D(NP) -REAL, INTENT(OUT) :: RSI(NP) -END SUBROUTINE ICECLOUD -END INTERFACE -END MODULE MODI_ICECLOUD diff --git a/src/arome/micro/modi_tiwmx.F90 b/src/arome/micro/modi_tiwmx.F90 deleted file mode 100644 index 137fdc80d..000000000 --- a/src/arome/micro/modi_tiwmx.F90 +++ /dev/null @@ -1,11 +0,0 @@ -! ######spl -MODULE MODI_TIWMX -! #################### -! -INTERFACE -FUNCTION TIWMX_TAB(P,T,QR,FICE,QRSN,RS,EPS) - REAL, INTENT(IN) :: P,T,QR,FICE,EPS - REAL, INTENT(OUT) :: QRSN,RS -END FUNCTION TIWMX_TAB -END INTERFACE -END MODULE MODI_TIWMX diff --git a/src/arome/micro/rain_ice_old.F90 b/src/arome/micro/rain_ice_old.F90 index 9a83426d0..97bb21a53 100644 --- a/src/arome/micro/rain_ice_old.F90 +++ b/src/arome/micro/rain_ice_old.F90 @@ -172,8 +172,8 @@ USE MODD_LES USE MODI_BUDGET_DDH USE MODI_GAMMA USE MODD_TIWMX -USE MODI_ICECLOUD -USE MODI_TIWMX +USE MODE_ICECLOUD, ONLY : ICECLOUD +USE MODE_TIWMX_TAB, ONLY : TIWMX_TAB USE DDH_MIX, ONLY : TYP_DDH USE YOMLDDH, ONLY : TLDDH USE YOMMDDH, ONLY : TMDDH diff --git a/src/arome/turb/shallow_mf.F90 b/src/arome/turb/shallow_mf.F90 index 03c55548d..d73839a6b 100644 --- a/src/arome/turb/shallow_mf.F90 +++ b/src/arome/turb/shallow_mf.F90 @@ -162,6 +162,7 @@ REAL, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) :: ZRSAT_UP ! Rsat in updraft LOGICAL :: GENTR_DETR ! flag to recompute entrainment, detrainment and mass flux INTEGER :: IKB ! near ground physical index INTEGER :: IKE ! uppest atmosphere physical index +INTEGER, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) :: IERR !------------------------------------------------------------------------ !!! 1. Initialisation @@ -186,8 +187,7 @@ ZFRAC_ICE(:,:) = 0. WHERE(PRM(:,:,2)+PRM(:,:,4) > 1.E-20) ZFRAC_ICE(:,:) = PRM(:,:,4) / (PRM(:,:,2)+PRM(:,:,4)) ENDWHERE -CALL COMPUTE_FRAC_ICE(SIZE(ZFRAC_ICE, 1), SIZE(ZFRAC_ICE, 2), & - HFRAC_ICE,ZFRAC_ICE(:,:),PTHM(:,:)*PEXNM(:,:)) +CALL COMPUTE_FRAC_ICE(HFRAC_ICE,ZFRAC_ICE(:,:),PTHM(:,:)*PEXNM(:,:), IERR(:,:)) ! Conservative variables at t-dt CALL THL_RT_FROM_TH_R_MF(KRR,KRRL,KRRI, & diff --git a/src/arome/turb/th_r_from_thl_rt_1d.F90 b/src/arome/turb/th_r_from_thl_rt_1d.F90 index 6add0bb64..f553f4585 100644 --- a/src/arome/turb/th_r_from_thl_rt_1d.F90 +++ b/src/arome/turb/th_r_from_thl_rt_1d.F90 @@ -78,6 +78,7 @@ REAL, DIMENSION(SIZE(PTHL,1)) :: ZDRSATODT,ZDRSATODTW,ZDRSATODTI REAL, DIMENSION(SIZE(PTHL,1)) :: ZFOESW, ZFOESI REAL, DIMENSION(SIZE(PTHL,1)) :: ZLOGT, Z99PP, Z1PRT REAL(KIND=JPRB) :: ZVAR1, ZVAR2, ZTPOW2, ZDELT +INTEGER, DIMENSION(SIZE(PTHL,1)) :: IERR REAL(KIND=JPRB) :: ZHOOK_HANDLE !---------------------------------------------------------------------------- @@ -127,7 +128,7 @@ DO II=1,JITER PFRAC_ICE(J) = PRI(J) / (PRL(J)+PRI(J)) ENDIF ENDDO - CALL COMPUTE_FRAC_ICE(SIZE(PFRAC_ICE, 1), HFRAC_ICE,PFRAC_ICE(:),ZT(:)) + CALL COMPUTE_FRAC_ICE(HFRAC_ICE,PFRAC_ICE(:),ZT(:), IERR(:)) !Computation of Rvsat and dRsat/dT !In this version QSAT, QSATI, DQSAT and DQASATI functions are not used diff --git a/src/arome/modset_Ryad/mpa/micro/internals/condensation.F90 b/src/common/micro/condensation.F90 similarity index 59% rename from src/arome/modset_Ryad/mpa/micro/internals/condensation.F90 rename to src/common/micro/condensation.F90 index 70ff1aa90..2e8112fc2 100644 --- a/src/arome/modset_Ryad/mpa/micro/internals/condensation.F90 +++ b/src/common/micro/condensation.F90 @@ -1,12 +1,17 @@ +!MNH_LIC Copyright 2002-2021 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 for details. version 1. +!----------------------------------------------------------------- ! ######spl SUBROUTINE CONDENSATION( KIU, KJU, KKU, KIB, KIE, KJB, KJE, KKB, KKE, KKL, & - HFRAC_ICE, & - PPABS, PZZ, PT, PRV_IN, PRV_OUT, PRC_IN, PRC_OUT, PRI_IN, PRI_OUT, PRS, PRG, PSIGS, PMFCONV, PCLDFR, PSIGRC, OUSERI,& - OSIGMAS, OCND2, PSIGQSAT, YSPP_PSIGQSAT,YSPP_ICE_CLD_WGT, & - PLV, PLS, PCPH ) - USE PARKIND1, ONLY : JPRB - USE MODD_SPP_TYPE - USE YOMHOOK , ONLY : LHOOK, DR_HOOK + HFRAC_ICE, HCONDENS, HLAMBDA3, & + PPABS, PZZ, PRHODREF, PT, PRV_IN, PRV_OUT, PRC_IN, PRC_OUT, PRI_IN, PRI_OUT, & + PRS, PRG, PSIGS, PMFCONV, PCLDFR, PSIGRC, OUSERI, & + OSIGMAS, OCND2, PSIGQSAT, & + PLV, PLS, PCPH, & + PHLC_HRC, PHLC_HCF, PHLI_HRI, PHLI_HCF, & + PICE_CLD_WGT) ! ################################################################################ ! !! @@ -14,8 +19,8 @@ !! ------- !!** Routine to diagnose cloud fraction, liquid and ice condensate mixing ratios !! and s'rl'/sigs^2 -!! -!! +!! +!! !!** METHOD !! ------ !! Based on the large-scale fields of temperature, water vapor, and possibly @@ -67,20 +72,24 @@ !! Safety limitation to .99*Pressure for saturation vapour pressure !! 2012-02 Y. Seity, add possibility to run with reversed vertical levels !! 2014-11 K.I Ivarsson add possibility to run with OCND2 option +!! 2016 S.Riette Change INQ1 !! 2016-11 S. Riette: use HFRAC_ICE, output adjusted state !! 2020-12 U. Andrae : Introduce SPP for HARMONIE-AROME !! R. El Khatib 24-Aug-2021 Optimizations +!! 2021-01: SPP computations moved in aro_adjust (AROME/HARMONIE) !------------------------------------------------------------------------------- ! !* 0. DECLARATIONS ! ------------ ! +USE PARKIND1, ONLY : JPRB +USE YOMHOOK , ONLY : LHOOK, DR_HOOK USE MODD_CST USE MODD_PARAMETERS -USE MODI_COMPUTE_FRAC_ICE -USE MODD_TIWMX -USE MODI_ICECLOUD -USE MODI_TIWMX +USE MODD_RAIN_ICE_PARAM, ONLY : XCRIAUTC, XCRIAUTI, XACRIAUTI, XBCRIAUTI +USE MODE_COMPUTE_FRAC_ICE, ONLY : COMPUTE_FRAC_ICE +USE MODE_TIWMX, ONLY : ESATW, ESATI +USE MODE_ICECLOUD, ONLY : ICECLOUD ! IMPLICIT NONE ! @@ -97,35 +106,43 @@ INTEGER, INTENT(IN) :: KJE ! value of the last point INTEGER, INTENT(IN) :: KKB ! value of the first point in z INTEGER, INTENT(IN) :: KKE ! value of the last point in z INTEGER, INTENT(IN) :: KKL ! +1 if grid goes from ground to atmosphere top, -1 otherwise -CHARACTER*1, INTENT(IN) :: HFRAC_ICE +CHARACTER(LEN=1), INTENT(IN) :: HFRAC_ICE +CHARACTER(LEN=4), INTENT(IN) :: HCONDENS +CHARACTER(LEN=*), INTENT(IN) :: HLAMBDA3 ! formulation for lambda3 coeff REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PPABS ! pressure (Pa) REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PZZ ! height of model levels (m) +REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PRHODREF REAL, DIMENSION(KIU,KJU,KKU), INTENT(INOUT) :: PT ! grid scale T (K) REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PRV_IN ! grid scale water vapor mixing ratio (kg/kg) in input REAL, DIMENSION(KIU,KJU,KKU), INTENT(OUT) :: PRV_OUT! grid scale water vapor mixing ratio (kg/kg) in output +REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PRC_IN ! grid scale r_c mixing ratio (kg/kg) in input +REAL, DIMENSION(KIU,KJU,KKU), INTENT(OUT) :: PRC_OUT! grid scale r_c mixing ratio (kg/kg) in output +REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PRI_IN ! grid scale r_i (kg/kg) in input +REAL, DIMENSION(KIU,KJU,KKU), INTENT(OUT) :: PRI_OUT! grid scale r_i (kg/kg) in output LOGICAL, INTENT(IN) :: OUSERI ! logical switch to compute both ! liquid and solid condensate (OUSERI=.TRUE.) ! or only solid condensate (OUSERI=.FALSE.) -LOGICAL, INTENT(IN) :: OCND2 ! logical switch to sparate liquid and ice - ! more rigid (DEFALT value : .FALSE.) LOGICAL, INTENT(IN) :: OSIGMAS! use present global Sigma_s values ! or that from turbulence scheme -REAL, INTENT(IN) :: PSIGQSAT ! use an extra "qsat" variance contribution (OSIGMAS case) +LOGICAL, INTENT(IN) :: OCND2 ! logical switch to sparate liquid and ice + ! more rigid (DEFALT value : .FALSE.) +REAL, DIMENSION(KIU,KJU), INTENT(IN) :: PSIGQSAT ! use an extra "qsat" variance contribution (OSIGMAS case) ! multiplied by PSIGQSAT -TYPE(TSPP_CONFIG_MPA), INTENT(IN) :: YSPP_PSIGQSAT,YSPP_ICE_CLD_WGT -REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PRC_IN ! grid scale r_c mixing ratio (kg/kg) in input -REAL, DIMENSION(KIU,KJU,KKU), INTENT(OUT) :: PRC_OUT! grid scale r_c mixing ratio (kg/kg) in output -REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PRI_IN ! grid scale r_i (kg/kg) in input -REAL, DIMENSION(KIU,KJU,KKU), INTENT(OUT) :: PRI_OUT! grid scale r_i (kg/kg) in output REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PRS ! grid scale mixing ration of snow (kg/kg) REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PRG ! grid scale mixing ration of graupel (kg/kg) REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PSIGS ! Sigma_s from turbulence scheme REAL, DIMENSION(:,:,:), INTENT(IN) :: PMFCONV! convective mass flux (kg /s m^2) REAL, DIMENSION(KIU,KJU,KKU), INTENT(OUT) :: PCLDFR ! cloud fraction REAL, DIMENSION(KIU,KJU,KKU), INTENT(OUT) :: PSIGRC ! s r_c / sig_s^2 -REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(IN) :: PLV -REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(IN) :: PLS -REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(IN) :: PCPH + +REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(IN) :: PLV ! Latent heat L_v +REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(IN) :: PLS ! Latent heat L_s +REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(IN) :: PCPH ! Specific heat C_ph +REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(OUT) :: PHLC_HRC +REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(OUT) :: PHLC_HCF ! cloud fraction +REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(OUT) :: PHLI_HRI +REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(OUT) :: PHLI_HCF +REAL, DIMENSION(KIU,KJU), OPTIONAL, INTENT(IN) :: PICE_CLD_WGT ! ! !* 0.2 Declarations of local variables : @@ -137,32 +154,37 @@ INTEGER, DIMENSION(KIU,KJU) :: ITPL ! top levels of troposphere REAL, DIMENSION(KIU,KJU) :: ZTMIN ! minimum Temp. related to ITPL ! REAL, DIMENSION(KIU,KJU,KKU) :: ZLV, ZLS, ZCPD -REAL :: ZLVS ! thermodynamics +REAL :: ZGCOND, ZAUTC, ZAUTI, ZGAUV, ZGAUC, ZGAUI, ZGAUTC, ZGAUTI, ZCRIAUTI ! Used for Gaussian PDF integration +REAL :: ZLVS ! thermodynamics +REAL, DIMENSION(KIB:KIE) :: ZPV, ZPIV, ZQSL, ZQSI ! thermodynamics REAL :: ZLL, DZZ, ZZZ ! used for length scales -REAL :: ZAH, ZDRW, ZDTL, ZSIG_CONV ! related to computation of Sig_s - -REAL, DIMENSION(KIB:KIE) :: ZTEMP, ZQSL, ZQSI, ZA, ZB, ZSIGMA, ZQ1, ZCOND, ZSBAR, ZPV, ZPIV, ESATW_T, ZDUM4 +REAL :: ZAH, ZDRW, ZDTL, ZSIG_CONV ! related to computation of Sig_s +REAL, DIMENSION(KIB:KIE) :: ZA, ZB, ZSBAR, ZSIGMA, ZQ1 ! related to computation of Sig_s +REAL, DIMENSION(KIB:KIE) :: ZCOND REAL, DIMENSION(KIB:KIE) :: ZFRAC ! Ice fraction +INTEGER :: INQ1 +REAL :: ZINC ! related to OCND2 noise check : -REAL :: ZRSP, ZRSW, ZRFRAC,ZRSDIF, ZRCOLD +REAL :: ZRSP, ZRSW, ZRFRAC, ZRSDIF, ZRCOLD ! related to OCND2 ice cloud calulation : -REAL :: ZDUM1,ZDUM2,ZDUM3,ZPRIFACT +REAL, DIMENSION(KIB:KIE) :: ESATW_T +REAL :: ZDUM1,ZDUM2,ZDUM3,ZDUM4,ZPRIFACT REAL, DIMENSION(KIU,KJU,KKU) :: TCLD REAL :: ZDZ(KIB:KIE), ZARDUM(KIE-KIB+1),ZCLDUM(KIE-KIB+1) -REAL, DIMENSION(KIB:KIE) :: ZSIGQSAT,ZICE_CLD_WGT ! end OCND2 - +REAL(KIND=JPRB) :: ZHOOK_HANDLE +INTEGER, DIMENSION(KIU) :: IERR ! !* 0.3 Definition of constants : ! !------------------------------------------------------------------------------- ! -REAL :: ZL0 = 600. ! tropospheric length scale -REAL :: ZCSIGMA = 0.2 ! constant in sigma_s parameterization -REAL :: ZCSIG_CONV = 0.30E-2 ! scaling factor for ZSIG_CONV as function of mass flux +REAL,PARAMETER :: ZL0 = 600. ! tropospheric length scale +REAL,PARAMETER :: ZCSIGMA = 0.2 ! constant in sigma_s parameterization +REAL,PARAMETER :: ZCSIG_CONV = 0.30E-2 ! scaling factor for ZSIG_CONV as function of mass flux ! -REAL, DIMENSION(-22:11) :: ZSRC_1D =(/ & +REAL, DIMENSION(-22:11),PARAMETER :: ZSRC_1D =(/ & 0. , 0. , 2.0094444E-04, 0.316670E-03, & 4.9965648E-04, 0.785956E-03 , 1.2341294E-03, 0.193327E-02, & 3.0190963E-03, 0.470144E-02 , 7.2950651E-03, 0.112759E-01, & @@ -172,13 +194,10 @@ REAL, DIMENSION(-22:11) :: ZSRC_1D =(/ & 0.8413813 , 0.933222E+00 , 0.9772662 , 0.993797E+00, & 0.9986521 , 0.999768E+00 , 0.9999684 , 0.999997E+00, & 1.0000000 , 1.000000 /) -INTEGER :: INQ1 -REAL :: ZINC,ZMU,ZVAL ! !------------------------------------------------------------------------------- ! ! -REAL(KIND=JPRB) :: ZHOOK_HANDLE IF (LHOOK) CALL DR_HOOK('CONDENSATION',0,ZHOOK_HANDLE) IKTB=1+JPVEXT @@ -189,50 +208,6 @@ PSIGRC(:,:,:) = 0. ! Initialize values ZPRIFACT = 1. ! Initialize value ZCLDUM=-1. ! Initialize value IF(OCND2)ZPRIFACT = 0. - -! -! Perturb PSIGQSAT -! - -IF (YSPP_PSIGQSAT%LPERT) THEN - IF (YSPP_PSIGQSAT%LLNN_MEAN1.OR.YSPP_PSIGQSAT%LLNN_MEAN1_SELF) THEN - ZMU = -0.5_JPRB * (YSPP_PSIGQSAT%CMPERT * YSPP_PSIGQSAT%SDEV)**2 - ELSE - ZMU = 0._JPRB - ENDIF - DO JI=KIB,KIE - ZVAL = & - PSIGQSAT*EXP(ZMU+YSPP_PSIGQSAT%CMPERT*YSPP_PSIGQSAT%PGP2DSPP(JI)) - ZSIGQSAT(JI) = MAX(YSPP_PSIGQSAT%CLIP(1),MIN(ZVAL,YSPP_PSIGQSAT%CLIP(2))) - ENDDO -ELSE - ZSIGQSAT(KIB:KIE) = PSIGQSAT -ENDIF - -! -! Perturb ICE_CLD_WGT -! - -IF (YSPP_ICE_CLD_WGT%LPERT) THEN - IF (YSPP_ICE_CLD_WGT%LLNN_MEAN1.OR.YSPP_ICE_CLD_WGT%LLNN_MEAN1_SELF) THEN - ZMU = -0.5_JPRB * (YSPP_ICE_CLD_WGT%CMPERT * YSPP_ICE_CLD_WGT%SDEV)**2 - ELSE - ZMU = 0._JPRB - ENDIF - DO JI=KIB,KIE -! Awaiting HARMONIE-AROME physics changes -! ZVAL = & -! XFRMIN(21)* EXP(ZMU+YSPP_ICE_CLD_WGT%CMPERT*YSPP_ICE_CLD_WGT%PGP2DSPP(JI)) - ZVAL = & - 1.5* EXP(ZMU+YSPP_ICE_CLD_WGT%CMPERT*YSPP_ICE_CLD_WGT%PGP2DSPP(JI)) - ZICE_CLD_WGT(JI) = & - MAX(YSPP_ICE_CLD_WGT%CLIP(1),MIN(ZVAL,YSPP_ICE_CLD_WGT%CLIP(2))) - ENDDO -ELSE -! ZICE_CLD_WGT(KIB:KIE) = XFRMIN(21) - ZICE_CLD_WGT(KIB:KIE) = 1.5 -ENDIF - ! ! !------------------------------------------------------------------------------- @@ -280,7 +255,7 @@ IF ( .NOT. OSIGMAS ) THEN DO JI=KIB,KIE ! store temperature at saturation ZTLK(JI,JJ,JK) = PT(JI,JJ,JK) - ZLV(JI,JJ,JK)*PRC_IN(JI,JJ,JK)/ZCPD(JI,JJ,JK) & - - ZLS(JI,JJ,JK)*PRI_IN(JI,JJ,JK)/ZCPD(JI,JJ,JK)*ZPRIFACT + - ZLS(JI,JJ,JK)*PRI_IN(JI,JJ,JK)/ZCPD(JI,JJ,JK)*ZPRIFACT END DO END DO END DO @@ -319,7 +294,7 @@ END IF ! DO JK=IKTB,IKTE JKP=MAX(MIN(JK+KKL,IKTE),IKTB) - JKM=MAX(MIN(JK-KKL,IKTE),IKTB) + JKM=MAX(MIN(JK-KKL,IKTE),IKTB) DO JJ=KJB,KJE IF (OCND2) THEN ZDZ(KIB:KIE) = PZZ(KIB:KIE,JJ,JKP) - PZZ(KIB:KIE,JJ,JKP+1) @@ -327,16 +302,16 @@ DO JK=IKTB,IKTE & PT(KIB,JJ,JK),PRV_IN(KIB,JJ,JK),1.,-1., & & ZCLDUM,1.,TCLD(KIB,JJ,JK), & & ZARDUM,ZARDUM,ZARDUM,ZARDUM) - ! latent heats - ! saturated water vapor mixing ratio over liquid water and ice - DO JI=KIB,KIE - ESATW_T(JI)=ESATW(PT(JI,JJ,JK)) - ZPV(JI) = MIN(ESATW_T(JI), .99*PPABS(JI,JJ,JK)) - ZPIV(JI) = MIN(ESATI(PT(JI,JJ,JK)), .99*PPABS(JI,JJ,JK)) - END DO + ! latent heats + ! saturated water vapor mixing ratio over liquid water and ice + DO JI=KIB,KIE + ESATW_T(JI)=ESATW(PT(JI,JJ,JK)) + ZPV(JI) = MIN(ESATW_T(JI), .99*PPABS(JI,JJ,JK)) + ZPIV(JI) = MIN(ESATI(PT(JI,JJ,JK)), .99*PPABS(JI,JJ,JK)) + END DO ELSE - ! latent heats - ! saturated water vapor mixing ratio over liquid water and ice + ! latent heats + ! saturated water vapor mixing ratio over liquid water and ice DO JI=KIB,KIE ZPV(JI) = MIN(EXP( XALPW - XBETAW / PT(JI,JJ,JK) - XGAMW * LOG( PT(JI,JJ,JK) ) ), .99*PPABS(JI,JJ,JK)) ZPIV(JI) = MIN(EXP( XALPI - XBETAI / PT(JI,JJ,JK) - XGAMI * LOG( PT(JI,JJ,JK) ) ), .99*PPABS(JI,JJ,JK)) @@ -350,7 +325,7 @@ DO JK=IKTB,IKTE ZFRAC(JI) = PRI_IN(JI,JJ,JK) / (PRC_IN(JI,JJ,JK)+PRI_IN(JI,JJ,JK)) ENDIF END DO - CALL COMPUTE_FRAC_ICE(HFRAC_ICE, ZFRAC(:), PT(:,JJ,JK)) + CALL COMPUTE_FRAC_ICE(HFRAC_ICE, ZFRAC(:), PT(:,JJ,JK), IERR) !error code IERR cannot be checked here to not break vectorization ENDIF DO JI=KIB,KIE ZQSL(JI) = XRD / XRV * ZPV(JI) / ( PPABS(JI,JJ,JK) - ZPV(JI) ) @@ -372,8 +347,8 @@ DO JK=IKTB,IKTE ! or that of Meso-NH turbulence scheme IF ( OSIGMAS ) THEN DO JI=KIB,KIE - IF (ZSIGQSAT(JI)/=0.) THEN - ZSIGMA(JI) = SQRT((2*PSIGS(JI,JJ,JK))**2 + (ZSIGQSAT(JI)*ZQSL(JI)*ZA(JI))**2) + IF (PSIGQSAT(JI,JJ)/=0.) THEN + ZSIGMA(JI) = SQRT((2*PSIGS(JI,JJ,JK))**2 + (PSIGQSAT(JI,JJ)*ZQSL(JI)*ZA(JI))**2) ELSE ZSIGMA(JI) = 2*PSIGS(JI,JJ,JK) END IF @@ -385,9 +360,9 @@ DO JK=IKTB,IKTE ZDRW = ZRT(JI,JJ,JKP) - ZRT(JI,JJ,JKM) ZDTL = ZTLK(JI,JJ,JKP) - ZTLK(JI,JJ,JKM) + XG/ZCPD(JI,JJ,JK) * DZZ ZLL = ZL(JI,JJ,JK) - ! standard deviation due to convection + ! standard deviation due to convection ZSIG_CONV =0. - IF( SIZE(PMFCONV) /= 0) ZSIG_CONV = ZCSIG_CONV * PMFCONV(JI,JJ,JK) / ZA(JI) + IF( SIZE(PMFCONV) /= 0) ZSIG_CONV = ZCSIG_CONV * PMFCONV(JI,JJ,JK) / ZA(JI) ! zsigma should be of order 4.e-4 in lowest 5 km of atmosphere ZSIGMA(JI) = SQRT( MAX( 1.E-25, ZCSIGMA * ZCSIGMA * ZLL*ZLL/(DZZ*DZZ)*(& ZA(JI)*ZA(JI)*ZDRW*ZDRW - 2.*ZA(JI)*ZB(JI)*ZDRW*ZDTL + ZB(JI)*ZB(JI)*ZDTL*ZDTL) + & @@ -396,27 +371,103 @@ DO JK=IKTB,IKTE END IF DO JI=KIB,KIE ZSIGMA(JI)= MAX( 1.E-10, ZSIGMA(JI) ) + ! normalized saturation deficit ZQ1(JI) = ZSBAR(JI)/ZSIGMA(JI) - ! total condensate - IF (ZQ1(JI) > 0. .AND. ZQ1(JI) <= 2 ) THEN - ZCOND(JI) = MIN(EXP(-1.)+.66*ZQ1(JI)+.086*ZQ1(JI)**2, 2.) ! We use the MIN function for continuity - ELSE IF (ZQ1(JI) > 2.) THEN - ZCOND(JI) = ZQ1(JI) - ELSE - ZCOND(JI) = EXP( 1.2*ZQ1(JI)-1. ) + END DO + IF(HCONDENS == 'GAUS') THEN + DO JI=KIB,KIE + ! Gaussian Probability Density Function around ZQ1 + ! Computation of ZG and ZGAM(=erf(ZG)) + ZGCOND = -ZQ1(JI)/SQRT(2.) + + !Approximation of erf function for Gaussian distribution + ZGAUV = 1 - SIGN(1., ZGCOND) * SQRT(1-EXP(-4*ZGCOND**2/XPI)) + + !Computation Cloud Fraction + PCLDFR(JI,JJ,JK) = MAX( 0., MIN(1.,0.5*ZGAUV)) + + !Computation of condensate + ZCOND(JI) = (EXP(-ZGCOND**2)-ZGCOND*SQRT(XPI)*ZGAUV)*ZSIGMA(JI)/SQRT(2.*XPI) + ZCOND(JI) = MAX(ZCOND(JI), 0.) + + PSIGRC(JI,JJ,JK) = PCLDFR(JI,JJ,JK) + END DO + !Computation warm/cold Cloud Fraction and content in high water content part + IF(PRESENT(PHLC_HCF) .AND. PRESENT(PHLC_HRC))THEN + DO JI=KIB,KIE + IF(1-ZFRAC(JI) > 1.E-20)THEN + ZAUTC = (ZSBAR(JI) - XCRIAUTC/(PRHODREF(JI,JJ,JK)*(1-ZFRAC(JI))))/ZSIGMA(JI) + ZGAUTC = -ZAUTC/SQRT(2.) + !Approximation of erf function for Gaussian distribution + ZGAUC = 1 - SIGN(1., ZGAUTC) * SQRT(1-EXP(-4*ZGAUTC**2/XPI)) + PHLC_HCF(JI,JJ,JK) = MAX( 0., MIN(1.,0.5*ZGAUC)) + PHLC_HRC(JI,JJ,JK) = (1-ZFRAC(JI))*(EXP(-ZGAUTC**2)-ZGAUTC*SQRT(XPI)*ZGAUC)*ZSIGMA(JI)/SQRT(2.*XPI) + PHLC_HRC(JI,JJ,JK) = PHLC_HRC(JI,JJ,JK) + XCRIAUTC/PRHODREF(JI,JJ,JK) * PHLC_HCF(JI,JJ,JK) + PHLC_HRC(JI,JJ,JK) = MAX(PHLC_HRC(JI,JJ,JK), 0.) + ELSE + PHLC_HCF(JI,JJ,JK)=0. + PHLC_HRC(JI,JJ,JK)=0. + ENDIF + END DO ENDIF - ZCOND(JI) = ZCOND(JI) * ZSIGMA(JI) - ! cloud fraction - IF ( ZCOND(JI) < 1.E-12 ) THEN - PCLDFR(JI,JJ,JK) = 0. - ELSE - PCLDFR(JI,JJ,JK) = MAX( 0., MIN(1.,0.5+0.36*ATAN(1.55*ZQ1(JI))) ) + + IF(PRESENT(PHLI_HCF) .AND. PRESENT(PHLI_HRI))THEN + DO JI=KIB,KIE + IF(ZFRAC(JI) > 1.E-20)THEN + ZCRIAUTI=MIN(XCRIAUTI,10**(XACRIAUTI*(PT(JI,JJ,JK)-XTT)+XBCRIAUTI)) + ZAUTI = (ZSBAR(JI) - ZCRIAUTI/ZFRAC(JI))/ZSIGMA(JI) + ZGAUTI = -ZAUTI/SQRT(2.) + !Approximation of erf function for Gaussian distribution + ZGAUI = 1 - SIGN(1., ZGAUTI) * SQRT(1-EXP(-4*ZGAUTI**2/XPI)) + PHLI_HCF(JI,JJ,JK) = MAX( 0., MIN(1.,0.5*ZGAUI)) + PHLI_HRI(JI,JJ,JK) = ZFRAC(JI)*(EXP(-ZGAUTI**2)-ZGAUTI*SQRT(XPI)*ZGAUI)*ZSIGMA(JI)/SQRT(2.*XPI) + PHLI_HRI(JI,JJ,JK) = PHLI_HRI(JI,JJ,JK) + ZCRIAUTI*PHLI_HCF(JI,JJ,JK) + PHLI_HRI(JI,JJ,JK) = MAX(PHLI_HRI(JI,JJ,JK), 0.) + ELSE + PHLI_HCF(JI,JJ,JK)=0. + PHLI_HRI(JI,JJ,JK)=0. + ENDIF + END DO ENDIF - IF (PCLDFR(JI,JJ,JK)==0.) THEN - ZCOND(JI)=0. + + ELSEIF(HCONDENS == 'CB02')THEN + DO JI=KIB,KIE + !Total condensate + IF (ZQ1(JI) > 0. .AND. ZQ1(JI) <= 2) THEN + ZCOND(JI) = MIN(EXP(-1.)+.66*ZQ1(JI)+.086*ZQ1(JI)**2, 2.) ! We use the MIN function for continuity + ELSE IF (ZQ1(JI) > 2.) THEN + ZCOND(JI) = ZQ1(JI) + ELSE + ZCOND(JI) = EXP( 1.2*ZQ1(JI)-1. ) + ENDIF + ZCOND(JI) = ZCOND(JI) * ZSIGMA(JI) + + !Cloud fraction + IF (ZCOND(JI) < 1.E-12) THEN + PCLDFR(JI,JJ,JK) = 0. + ELSE + PCLDFR(JI,JJ,JK) = MAX( 0., MIN(1.,0.5+0.36*ATAN(1.55*ZQ1(JI))) ) + ENDIF + IF (PCLDFR(JI,JJ,JK)==0.) THEN + ZCOND(JI)=0. + ENDIF + + INQ1 = MIN( MAX(-22,FLOOR(MIN(100., MAX(-100., 2*ZQ1(JI)))) ), 10) !inner min/max prevents sigfpe when 2*zq1 does not fit into an int + ZINC = 2.*ZQ1(JI) - INQ1 + + PSIGRC(JI,JJ,JK) = MIN(1.,(1.-ZINC)*ZSRC_1D(INQ1)+ZINC*ZSRC_1D(INQ1+1)) + END DO + IF(PRESENT(PHLC_HCF) .AND. PRESENT(PHLC_HRC))THEN + PHLC_HCF(:,JJ,JK)=0. + PHLC_HRC(:,JJ,JK)=0. ENDIF - END DO + IF(PRESENT(PHLI_HCF) .AND. PRESENT(PHLI_HRI))THEN + PHLI_HCF(:,JJ,JK)=0. + PHLI_HRI(:,JJ,JK)=0. + ENDIF + END IF !HCONDENS + IF(.NOT. OCND2) THEN DO JI=KIB,KIE PRC_OUT(JI,JJ,JK) = (1.-ZFRAC(JI)) * ZCOND(JI) ! liquid condensate @@ -428,72 +479,69 @@ DO JK=IKTB,IKTE END DO ELSE DO JI=KIB,KIE - PRC_OUT(JI,JJ,JK) = (1.-ZFRAC(JI)) * ZCOND(JI) ! liquid condensate - ! -! This check is mainly for noise reduction : -! ------------------------- - IF(ABS(PRC_IN(JI,JJ,JK)-PRC_OUT(JI,JJ,JK))>1.0E-12 .AND. ESATW_T(JI) < PPABS(JI,JJ,JK)*0.5)THEN - ZRCOLD = PRC_OUT(JI,JJ,JK) - ZRFRAC = PRV_IN(JI,JJ,JK) - ZCOND(JI) + PRC_OUT(JI,JJ,JK) - IF( PRV_IN(JI,JJ,JK) < ZRSW )THEN ! sub - saturation over water: - ! Avoid drying of cloudwater leading to supersaturation with - ! respect to water - ZRSDIF= MIN(0.,ZRSP-ZRFRAC) - ELSE ! super - saturation over water: - ! Avoid depostition of water leading to sub-saturation with - ! respect to water - ! ZRSDIF= MAX(0.,ZRSP-ZRFRAC) - ZRSDIF= MAX(0.,ZRSP*PCLDFR(JI,JJ,JK) - ZRFRAC) - ENDIF - PRC_OUT(JI,JJ,JK) = ZCOND(JI) - ZRSDIF - ELSE - ZRCOLD = PRC_IN(JI,JJ,JK) - ENDIF - ! end check + PRC_OUT(JI,JJ,JK) = (1.-ZFRAC(JI)) * ZCOND(JI) ! liquid condensate + ! +! This check is mainly for noise reduction : +! ------------------------- + IF(ABS(PRC_IN(JI,JJ,JK)-PRC_OUT(JI,JJ,JK))>1.0E-12 .AND. ESATW_T(JI) < PPABS(JI,JJ,JK)*0.5)THEN + ZRCOLD = PRC_OUT(JI,JJ,JK) + ZRFRAC = PRV_IN(JI,JJ,JK) - ZCOND(JI) + PRC_OUT(JI,JJ,JK) + IF( PRV_IN(JI,JJ,JK) < ZRSW )THEN ! sub - saturation over water: + ! Avoid drying of cloudwater leading to supersaturation with + ! respect to water + ZRSDIF= MIN(0.,ZRSP-ZRFRAC) + ELSE ! super - saturation over water: + ! Avoid depostition of water leading to sub-saturation with + ! respect to water + ! ZRSDIF= MAX(0.,ZRSP-ZRFRAC) + ZRSDIF= MAX(0.,ZRSP*PCLDFR(JI,JJ,JK) - ZRFRAC) + ENDIF + PRC_OUT(JI,JJ,JK) = ZCOND(JI) - ZRSDIF + ELSE + ZRCOLD = PRC_IN(JI,JJ,JK) + ENDIF + ! end check - ! compute separate ice cloud: + ! compute separate ice cloud: + ZDUM1 = MIN(1.0,20.* PRC_OUT(JI,JJ,JK)*SQRT(ZDZ(JI))/ZQSL(JI)) ! clould liquid water + ! factor - ZDUM1 = MIN(1.0,20.* PRC_OUT(JI,JJ,JK)*SQRT(ZDZ(JI))/ZQSL(JI)) ! clould liquid water - ! factor - ZDUM3 = MAX(0.,TCLD(JI,JJ,JK)-PCLDFR(JI,JJ,JK)) ! pure ice cloud part + ZDUM3 = MAX(0.,TCLD(JI,JJ,JK)-PCLDFR(JI,JJ,JK)) ! pure ice cloud part - IF (JK==IKTB) THEN - ZDUM4(JI) = PRI_IN(JI,JJ,JK) - ELSE - ZDUM4(JI) = PRI_IN(JI,JJ,JK) + PRS(JI,JJ,JK)*0.5 + PRG(JI,JJ,JK)*0.25 - ENDIF + IF (JK==IKTB) THEN + ZDUM4 = PRI_IN(JI,JJ,JK) + ELSE + ZDUM4 = PRI_IN(JI,JJ,JK) + PRS(JI,JJ,JK)*0.5 + PRG(JI,JJ,JK)*0.25 + ENDIF - ZDUM4(JI) = MAX(0.,MIN(1.,ZICE_CLD_WGT(JI)*ZDUM4(JI)*SQRT(ZDZ(JI))/ZQSI(JI))) ! clould ice+solid - ! precip. water factor + ZDUM4 = MAX(0.,MIN(1.,PICE_CLD_WGT(JI,JJ)*ZDUM4*SQRT(ZDZ(JI))/ZQSI(JI))) ! clould ice+solid + ! precip. water factor - ZDUM2 = (0.8*PCLDFR(JI,JJ,JK)+0.2)*MIN(1.,ZDUM1 + ZDUM4(JI)*PCLDFR(JI,JJ,JK)) - ! water cloud, use 'statistical' cloud, but reduce it in case of low liquid content + ZDUM2 = (0.8*PCLDFR(JI,JJ,JK)+0.2)*MIN(1.,ZDUM1 + ZDUM4*PCLDFR(JI,JJ,JK)) + ! water cloud, use 'statistical' cloud, but reduce it in case of low liquid content - PCLDFR(JI,JJ,JK) = MIN(1., ZDUM2 + (0.9*ZDUM3+0.1)*ZDUM4(JI)) ! Rad cloud - ! Reduce ice cloud part in case of low ice water content - PRI_OUT(JI,JJ,JK) = PRI_IN(JI,JJ,JK) + PCLDFR(JI,JJ,JK) = MIN(1., ZDUM2 + (0.9*ZDUM3+0.1)*ZDUM4) ! Rad cloud + ! Reduce ice cloud part in case of low ice water content + PRI_OUT(JI,JJ,JK) = PRI_IN(JI,JJ,JK) PT(JI,JJ,JK) = PT(JI,JJ,JK) + ((PRC_OUT(JI,JJ,JK)-ZRCOLD)*ZLV(JI,JJ,JK) + & &(PRI_OUT(JI,JJ,JK)-PRI_IN(JI,JJ,JK))*ZLS(JI,JJ,JK) ) & & /ZCPD(JI,JJ,JK) PRV_OUT(JI,JJ,JK) = ZRT(JI,JJ,JK) - PRC_OUT(JI,JJ,JK) - PRI_OUT(JI,JJ,JK)*ZPRIFACT END DO - ENDIF ! End OCND2 - DO JI=KIB,KIE - -! s r_c/ sig_s^2 -! PSIGRC(JI,JJ,JK) = PCLDFR(JI,JJ,JK) ! use simple Gaussian relation -! -! multiply PSRCS by the lambda3 coefficient -! -! PSIGRC(JI,JJ,JK) = 2.*PCLDFR(JI,JJ,JK) * MIN( 3. , MAX(1.,1.-ZQ1(JI)) ) -! in the 3D case lambda_3 = 1. - INQ1 = MIN( MAX(-22,FLOOR(MIN(100., MAX(-100., 2*ZQ1(JI)))) ), 10) !inner min/max prevents sigfpe when 2*zq1 does not fit into an int - ZINC = 2.*ZQ1(JI) - INQ1 - - PSIGRC(JI,JJ,JK) = MIN(1.,(1.-ZINC)*ZSRC_1D(INQ1)+ZINC*ZSRC_1D(INQ1+1)) + END IF ! End OCND2 + IF(HLAMBDA3=='CB')THEN + DO JI=KIB,KIE + ! s r_c/ sig_s^2 + ! PSIGRC(JI,JJ,JK) = PCLDFR(JI,JJ,JK) ! use simple Gaussian relation + ! + ! multiply PSRCS by the lambda3 coefficient + ! + ! PSIGRC(JI,JJ,JK) = 2.*PCLDFR(JI,JJ,JK) * MIN( 3. , MAX(1.,1.-ZQ1(JI)) ) + ! in the 3D case lambda_3 = 1. - PSIGRC(JI,JJ,JK) = PSIGRC(JI,JJ,JK)* MIN( 3. , MAX(1.,1.-ZQ1(JI)) ) - END DO + PSIGRC(JI,JJ,JK) = PSIGRC(JI,JJ,JK)* MIN( 3. , MAX(1.,1.-ZQ1(JI)) ) + END DO + END IF END DO END DO ! diff --git a/src/arome/modset_Ryad/mpa/micro/internals/ice_adjust.F90 b/src/common/micro/ice_adjust.F90 similarity index 64% rename from src/arome/modset_Ryad/mpa/micro/internals/ice_adjust.F90 rename to src/common/micro/ice_adjust.F90 index 41d3cb26e..355cb63cc 100644 --- a/src/arome/modset_Ryad/mpa/micro/internals/ice_adjust.F90 +++ b/src/common/micro/ice_adjust.F90 @@ -1,20 +1,21 @@ +!MNH_LIC Copyright 1996-2021 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 for details. version 1. +!----------------------------------------------------------------- ! ########################################################################## - SUBROUTINE ICE_ADJUST (KKA, KKU, KKL, KRR, HFRAC_ICE, & - HBUNAME, OSUBG_COND, OSIGMAS, OCND2, & + SUBROUTINE ICE_ADJUST (KKA, KKU, KKL, KRR, HFRAC_ICE, HCONDENS, HLAMBDA3,& + HBUNAME, OSUBG_COND, OSIGMAS, OCND2, HSUBG_MF_PDF,& PTSTEP, PSIGQSAT, & - PRHODJ, PEXNREF, PSIGS, PMFCONV, PPABST, PZZ, & + PRHODJ, PEXNREF, PRHODREF, PSIGS, PMFCONV, & + PPABST, PZZ, & PEXN, PCF_MF, PRC_MF, PRI_MF, & PRV, PRC, PRVS, PRCS, PTH, PTHS, PSRCS, PCLDFR, & - PRR, PRI, PRIS, PRS, PRG, & - PRH, POUT_RV, POUT_RC, POUT_RI, POUT_TH, & - YSPP_PSIGQSAT, YSPP_ICE_CLD_WGT, & - YDDDH, YDLDDH, YDMDDH) - USE PARKIND1, ONLY : JPRB - USE MODD_SPP_TYPE - USE YOMHOOK , ONLY : LHOOK, DR_HOOK - USE DDH_MIX , ONLY : TYP_DDH - USE YOMLDDH , ONLY : TLDDH - USE YOMMDDH , ONLY : TMDDH + PRR, PRI, PRIS, PRS, PRG, PRH, & + POUT_RV, POUT_RC, POUT_RI, POUT_TH, & + PHLC_HRC, PHLC_HCF, PHLI_HRI, PHLI_HCF, & + TBUDGETS, KBUDGETS, & + PICE_CLD_WGT) ! ######################################################################### ! !!**** *ICE_ADJUST* - compute the ajustment of water vapor in mixed-phase @@ -60,7 +61,6 @@ !! Module MODD_BUDGET: !! NBUMOD !! CBUTYPE -!! NBUPROCCTR !! LBU_RTH !! LBU_RRV !! LBU_RRC @@ -89,11 +89,14 @@ !! J.Pergaud and S.Malardel Add EDKF case !! S. Riette ice for EDKF !! 2012-02 Y. Seity, add possibility to run with reversed vertical levels +!! J.Escobar : 15/09/2015 : WENO5 & JPHEXT <> 1 !! 2016-07 S. Riette: adjustement is now realized on state variables (PRV, PRC, PRI, PTH) !! whereas tendencies are still applied on S variables. !! This modification allows to call ice_adjust on T variable !! or to call it on S variables !! 2016-11 S. Riette: all-or-nothing adjustment now uses condensation +! P. Wautelet 05/2016-04/2018: new data structures and calls for I/O +! P. Wautelet 02/2020: use the new data structures and subroutines for budgets !! 2020-12 U. Andrae : Introduce SPP for HARMONIE-AROME !! R. El Khatib 24-Aug-2021 Optimizations !! @@ -102,14 +105,19 @@ !* 0. DECLARATIONS ! ------------ ! -USE MODD_PARAMETERS -USE MODD_CST -USE MODD_CONF -USE MODD_BUDGET +USE PARKIND1, ONLY : JPRB +USE YOMHOOK , ONLY : LHOOK, DR_HOOK +USE MODD_BUDGET, ONLY: TBUDGETDATA, LBU_ENABLE, & + & LBUDGET_TH, LBUDGET_RV, LBUDGET_RC, LBUDGET_RI, & + & NBUDGET_TH, NBUDGET_RV, NBUDGET_RC, NBUDGET_RI +USE MODD_CST, ONLY: XLVTT, XLSTT, XCPV, XCL, XCI, XTT, XCPD, XCPV +USE MODD_PARAMETERS, ONLY: JPVEXT +USE MODD_RAIN_ICE_PARAM, ONLY : XCRIAUTC, XCRIAUTI, XACRIAUTI, XBCRIAUTI +! +USE MODE_BUDGET, ONLY: BUDGET_STORE_INIT, BUDGET_STORE_END +USE MODE_ll, ONLY: GET_INDICE_ll ! USE MODI_CONDENSATION -USE MODI_BUDGET -USE MODE_FMWRIT ! IMPLICIT NONE ! @@ -121,8 +129,10 @@ INTEGER, INTENT(IN) :: KKA !near ground array index INTEGER, INTENT(IN) :: KKU !uppest atmosphere array index INTEGER, INTENT(IN) :: KKL !vert. levels type 1=MNH -1=ARO INTEGER, INTENT(IN) :: KRR ! Number of moist variables -CHARACTER*1, INTENT(IN) :: HFRAC_ICE -CHARACTER*4, INTENT(IN) :: HBUNAME ! Name of the budget +CHARACTER(LEN=1), INTENT(IN) :: HFRAC_ICE +CHARACTER(LEN=80), INTENT(IN) :: HCONDENS +CHARACTER(LEN=4), INTENT(IN) :: HLAMBDA3 ! formulation for lambda3 coeff +CHARACTER(LEN=4), INTENT(IN) :: HBUNAME ! Name of the budget LOGICAL, INTENT(IN) :: OSUBG_COND ! Switch for Subgrid ! Condensation LOGICAL :: OSIGMAS ! Switch for Sigma_s: @@ -131,12 +141,14 @@ LOGICAL :: OSIGMAS ! Switch for Sigma_s: LOGICAL :: OCND2 ! logical switch to sparate liquid ! and ice ! more rigid (DEFALT value : .FALSE.) +CHARACTER(LEN=80), INTENT(IN) :: HSUBG_MF_PDF REAL, INTENT(IN) :: PTSTEP ! Double Time step ! (single if cold start) -REAL, INTENT(IN) :: PSIGQSAT ! coeff applied to qsat variance contribution +REAL, DIMENSION(:,:), INTENT(IN) :: PSIGQSAT ! coeff applied to qsat variance contribution ! REAL, DIMENSION(:,:,:), CONTIGUOUS, INTENT(IN) :: PRHODJ ! Dry density * Jacobian REAL, DIMENSION(:,:,:), CONTIGUOUS, INTENT(IN) :: PEXNREF ! Reference Exner function +REAL, DIMENSION(:,:,:), CONTIGUOUS, INTENT(IN) :: PRHODREF ! REAL, DIMENSION(:,:,:), CONTIGUOUS, INTENT(IN) :: PSIGS ! Sigma_s at time t REAL, DIMENSION(:,:,:), CONTIGUOUS, INTENT(IN) :: PMFCONV ! convective mass flux @@ -169,25 +181,26 @@ REAL, DIMENSION(:,:,:), CONTIGUOUS, OPTIONAL, INTENT(OUT) :: POUT_RV ! Adjuste REAL, DIMENSION(:,:,:), CONTIGUOUS, OPTIONAL, INTENT(OUT) :: POUT_RC ! Adjusted value REAL, DIMENSION(:,:,:), CONTIGUOUS, OPTIONAL, INTENT(OUT) :: POUT_RI ! Adjusted value REAL, DIMENSION(:,:,:), CONTIGUOUS, OPTIONAL, INTENT(OUT) :: POUT_TH ! Adjusted value -TYPE(TSPP_CONFIG_MPA), INTENT(IN) :: YSPP_PSIGQSAT,YSPP_ICE_CLD_WGT -TYPE(TYP_DDH) , INTENT(INOUT) :: YDDDH -TYPE(TLDDH) , INTENT(IN) :: YDLDDH -TYPE(TMDDH) , INTENT(IN) :: YDMDDH +REAL, DIMENSION(:,:,:), CONTIGUOUS, OPTIONAL, INTENT(OUT) :: PHLC_HRC +REAL, DIMENSION(:,:,:), CONTIGUOUS, OPTIONAL, INTENT(OUT) :: PHLC_HCF +REAL, DIMENSION(:,:,:), CONTIGUOUS, OPTIONAL, INTENT(OUT) :: PHLI_HRI +REAL, DIMENSION(:,:,:), CONTIGUOUS, OPTIONAL, INTENT(OUT) :: PHLI_HCF +TYPE(TBUDGETDATA), DIMENSION(KBUDGETS), INTENT(INOUT) :: TBUDGETS +INTEGER, INTENT(IN) :: KBUDGETS +REAL, DIMENSION(:,:), CONTIGUOUS, OPTIONAL, INTENT(IN) :: PICE_CLD_WGT ! !* 0.2 Declarations of local variables : ! ! -REAL :: ZEPS ! Mv/Md -REAL :: ZT00,ZT0 ! Min and max temperature for the mixed phase liquid and solid water - ! for the coeff CND of the barycentric mixing ratio REAL :: ZW1,ZW2 ! intermediate fields REAL, DIMENSION(SIZE(PEXNREF,1),SIZE(PEXNREF,2),SIZE(PEXNREF,3)) & :: ZT, & ! adjusted temperature ZRV, ZRC, ZRI, & ! adjusted state ZCPH, & ! guess of the CPh for the mixing ZLV, & ! guess of the Lv at t+1 - ZLS, & ! guess of the Ls at t+1 - ZH2O ! Fraction of all types of water (kg/kg) (for OCND2) + ZLS ! guess of the Ls at t+1 +REAL :: ZCRIAUT, & ! Autoconversion thresholds + ZHCF, ZHR ! INTEGER :: IIU,IJU,IKU! dimensions of dummy arrays INTEGER :: IIB,IJB ! Horz index values of the first inner mass points @@ -197,27 +210,30 @@ INTEGER :: IKE ! K index value of the last inner mass point INTEGER :: JITER,ITERMAX ! iterative loop for first order adjustment INTEGER :: JI, JJ, JK ! -REAL, DIMENSION(SIZE(PEXNREF,1),SIZE(PEXNREF,2),SIZE(PEXNREF,3)) :: ZSIGS -! +REAL, DIMENSION(SIZE(PEXNREF,1),SIZE(PEXNREF,2),SIZE(PEXNREF,3)) :: ZSIGS, ZSRCS +REAL, DIMENSION(SIZE(PEXNREF,1),SIZE(PEXNREF,2)) :: ZSIGQSAT REAL(KIND=JPRB) :: ZHOOK_HANDLE +! !------------------------------------------------------------------------------- ! !* 1. PRELIMINARIES ! ------------- ! IF (LHOOK) CALL DR_HOOK('ICE_ADJUST',0,ZHOOK_HANDLE) +! IIU = SIZE(PEXNREF,1) IJU = SIZE(PEXNREF,2) IKU = SIZE(PEXNREF,3) -IIB = 1 + JPHEXT -IIE = IIU - JPHEXT -IJB = 1 + JPHEXT -IJE = IJU - JPHEXT +CALL GET_INDICE_ll(IIB,IJB,IIE,IJE,IIU,IJU) IKB=KKA+JPVEXT*KKL IKE=KKU-JPVEXT*KKL ! ITERMAX=1 ! +IF(LBUDGET_TH) CALL BUDGET_STORE_INIT(TBUDGETS(NBUDGET_TH), TRIM(HBUNAME), PTHS(:, :, :)*PRHODJ(:, :, :)) +IF(LBUDGET_RV) CALL BUDGET_STORE_INIT(TBUDGETS(NBUDGET_RV), TRIM(HBUNAME), PRVS(:, :, :)*PRHODJ(:, :, :)) +IF(LBUDGET_RC) CALL BUDGET_STORE_INIT(TBUDGETS(NBUDGET_RC), TRIM(HBUNAME), PRCS(:, :, :)*PRHODJ(:, :, :)) +IF(LBUDGET_RI) CALL BUDGET_STORE_INIT(TBUDGETS(NBUDGET_RI), TRIM(HBUNAME), PRIS(:, :, :)*PRHODJ(:, :, :)) !------------------------------------------------------------------------------- ! !* 2. COMPUTE QUANTITIES WITH THE GUESS OF THE FUTURE INSTANT @@ -250,28 +266,25 @@ DO JITER =1,ITERMAX ! compute with updated values CALL ITERATION(ZRV,ZRC,ZRI,ZRV,ZRC,ZRI) ENDIF - ENDDO ! end of the iterative loop ! !* 5. COMPUTE THE SOURCES AND STORES THE CLOUD FRACTION ! ------------------------------------------------- ! -! ZW8 = 4.*ZW2 ! + 0.5_JPRB*PRSS*PTSTEP + 0.25_JPRB*PRGS*PTSTEP ! include snow and graupel to cloudcover -! ZW8(:,:,1) = ZW2(:,:,1) ! but don't do this for the lowest model -! ZW8(:,:,2) = ZW2(:,:,2) ! level (Hail not concidered yet) -! -!* 5.0 compute the variation of mixing ratio ! DO JK=1,IKU DO JJ=1,IJU DO JI=1,IIU + ! + !* 5.0 compute the variation of mixing ratio + ! ! Rc - Rc* ZW1 = (ZRC(JI,JJ,JK) - PRC(JI,JJ,JK)) / PTSTEP ! Pcon = ---------- ! 2 Delta t ZW2 = (ZRI(JI,JJ,JK) - PRI(JI,JJ,JK)) / PTSTEP ! idem ZW1 but for Ri -! -!* 5.1 compute the sources -! + ! + !* 5.1 compute the sources + ! IF( ZW1 < 0.0 ) THEN ZW1 = MAX ( ZW1, -PRCS(JI,JJ,JK) ) ELSE @@ -281,7 +294,7 @@ DO JK=1,IKU PRCS(JI,JJ,JK) = PRCS(JI,JJ,JK) + ZW1 PTHS(JI,JJ,JK) = PTHS(JI,JJ,JK) + & ZW1 * ZLV(JI,JJ,JK) / (ZCPH(JI,JJ,JK) * PEXNREF(JI,JJ,JK)) -! + ! IF( ZW2 < 0.0 ) THEN ZW2 = MAX ( ZW2, -PRIS(JI,JJ,JK) ) ELSE @@ -292,10 +305,9 @@ DO JK=1,IKU PTHS(JI,JJ,JK) = PTHS(JI,JJ,JK) + & ZW2 * ZLS(JI,JJ,JK) / (ZCPH(JI,JJ,JK) * PEXNREF(JI,JJ,JK)) ENDDO -! -! -!* 5.2 compute the cloud fraction PCLDFR -! + ! + !* 5.2 compute the cloud fraction PCLDFR + ! IF ( .NOT. OSUBG_COND ) THEN DO JI=1,IIU IF (PRCS(JI,JJ,JK) + PRIS(JI,JJ,JK) > 1.E-12 / PTSTEP) THEN @@ -307,7 +319,7 @@ DO JK=1,IKU PSRCS(JI,JJ,JK) = PCLDFR(JI,JJ,JK) END IF ENDDO - ELSE + ELSE !OSUBG_COND case DO JI=1,IIU !We limit PRC_MF+PRI_MF to PRVS*PTSTEP to avoid negative humidity ZW1=PRC_MF(JI,JJ,JK)/PTSTEP @@ -322,9 +334,64 @@ DO JK=1,IKU PRVS(JI,JJ,JK)=PRVS(JI,JJ,JK)-(ZW1+ZW2) PTHS(JI,JJ,JK) = PTHS(JI,JJ,JK) + & (ZW1 * ZLV(JI,JJ,JK) + ZW2 * ZLS(JI,JJ,JK)) / ZCPH(JI,JJ,JK) / PEXNREF(JI,JJ,JK) + ! + IF(PRESENT(PHLC_HRC) .AND. PRESENT(PHLC_HCF)) THEN + ZCRIAUT=XCRIAUTC/PRHODREF(JI,JJ,JK) + IF(HSUBG_MF_PDF=='NONE')THEN + IF(ZW1*PTSTEP>PCF_MF(JI,JJ,JK) * ZCRIAUT) THEN + PHLC_HRC(JI,JJ,JK)=PHLC_HRC(JI,JJ,JK)+ZW1*PTSTEP + PHLC_HCF(JI,JJ,JK)=MIN(1.,PHLC_HCF(JI,JJ,JK)+PCF_MF(JI,JJ,JK)) + ENDIF + ELSEIF(HSUBG_MF_PDF=='TRIANGLE')THEN + !ZHCF is the precipitating part of the *cloud* and not of the grid cell + IF(ZW1*PTSTEP>PCF_MF(JI,JJ,JK)*ZCRIAUT) THEN + ZHCF=1.-.5*(ZCRIAUT*PCF_MF(JI,JJ,JK) / MAX(1.E-20, ZW1*PTSTEP))**2 + ZHR=ZW1*PTSTEP-(ZCRIAUT*PCF_MF(JI,JJ,JK))**3 / & + &(3*MAX(1.E-20, ZW1*PTSTEP)**2) + ELSEIF(2.*ZW1*PTSTEP<=PCF_MF(JI,JJ,JK) * ZCRIAUT) THEN + ZHCF=0. + ZHR=0. + ELSE + ZHCF=(2.*ZW1*PTSTEP-ZCRIAUT*PCF_MF(JI,JJ,JK))**2 / & + &(2.*MAX(1.E-20, ZW1*PTSTEP)**2) + ZHR=(4.*(ZW1*PTSTEP)**3-3.*ZW1*PTSTEP*(ZCRIAUT*PCF_MF(JI,JJ,JK))**2+& + (ZCRIAUT*PCF_MF(JI,JJ,JK))**3) / & + &(3*MAX(1.E-20, ZW1*PTSTEP)**2) + ENDIF + ZHCF=ZHCF*PCF_MF(JI,JJ,JK) !to retrieve the part of the grid cell + PHLC_HCF(JI,JJ,JK)=MIN(1.,PHLC_HCF(JI,JJ,JK)+ZHCF) !total part of the grid cell that is precipitating + PHLC_HRC(JI,JJ,JK)=PHLC_HRC(JI,JJ,JK)+ZHR + ENDIF + ENDIF + IF(PRESENT(PHLI_HRI) .AND. PRESENT(PHLI_HCF)) THEN + ZCRIAUT=MIN(XCRIAUTI,10**(XACRIAUTI*(ZT(JI,JJ,JK)-XTT)+XBCRIAUTI)) + IF(HSUBG_MF_PDF=='NONE')THEN + IF(ZW2*PTSTEP>PCF_MF(JI,JJ,JK) * ZCRIAUT) THEN + PHLI_HRI(JI,JJ,JK)=PHLI_HRI(JI,JJ,JK)+ZW2*PTSTEP + PHLI_HCF(JI,JJ,JK)=MIN(1.,PHLI_HCF(JI,JJ,JK)+PCF_MF(JI,JJ,JK)) + ENDIF + ELSEIF(HSUBG_MF_PDF=='TRIANGLE')THEN + !ZHCF is the precipitating part of the *cloud* and not of the grid cell + IF(ZW2*PTSTEP>PCF_MF(JI,JJ,JK)*ZCRIAUT) THEN + ZHCF=1.-.5*(ZCRIAUT*PCF_MF(JI,JJ,JK) / (ZW2*PTSTEP))**2 + ZHR=ZW2*PTSTEP-(ZCRIAUT*PCF_MF(JI,JJ,JK))**3/(3*(ZW2*PTSTEP)**2) + ELSEIF(2.*ZW2*PTSTEP<=PCF_MF(JI,JJ,JK) * ZCRIAUT) THEN + ZHCF=0. + ZHR=0. + ELSE + ZHCF=(2.*ZW2*PTSTEP-ZCRIAUT*PCF_MF(JI,JJ,JK))**2 / (2.*(ZW2*PTSTEP)**2) + ZHR=(4.*(ZW2*PTSTEP)**3-3.*ZW2*PTSTEP*(ZCRIAUT*PCF_MF(JI,JJ,JK))**2+& + (ZCRIAUT*PCF_MF(JI,JJ,JK))**3)/(3*(ZW2*PTSTEP)**2) + ENDIF + ZHCF=ZHCF*PCF_MF(JI,JJ,JK) !to retrieve the part of the grid cell + PHLI_HCF(JI,JJ,JK)=MIN(1.,PHLI_HCF(JI,JJ,JK)+ZHCF) !total part of the grid cell that is precipitating + PHLI_HRI(JI,JJ,JK)=PHLI_HRI(JI,JJ,JK)+ZHR + ENDIF + ENDIF ENDDO + ! IF(PRESENT(POUT_RV) .OR. PRESENT(POUT_RC) .OR. & - &PRESENT(POUT_RI) .OR. PRESENT(POUT_TH)) THEN + &PRESENT(POUT_RI) .OR. PRESENT(POUT_TH)) THEN DO JI=1,IIU ZW1=PRC_MF(JI,JJ,JK) ZW2=PRI_MF(JI,JJ,JK) @@ -337,30 +404,30 @@ DO JK=1,IKU ZRV(JI,JJ,JK)=ZRV(JI,JJ,JK)-(ZW1+ZW2) ZT(JI,JJ,JK) = ZT(JI,JJ,JK) + & (ZW1 * ZLV(JI,JJ,JK) + ZW2 * ZLS(JI,JJ,JK)) / ZCPH(JI,JJ,JK) - IF(PRESENT(POUT_RV)) POUT_RV(JI,JJ,JK)=ZRV(JI,JJ,JK) - IF(PRESENT(POUT_RC)) POUT_RC(JI,JJ,JK)=ZRC(JI,JJ,JK) - IF(PRESENT(POUT_RI)) POUT_RI(JI,JJ,JK)=ZRI(JI,JJ,JK) - IF(PRESENT(POUT_TH)) POUT_TH(JI,JJ,JK)=ZT(JI,JJ,JK) / PEXN(JI,JJ,JK) ENDDO ENDIF - ENDIF + ENDIF !OSUBG_COND ENDDO ENDDO ! +IF(PRESENT(POUT_RV)) POUT_RV=ZRV +IF(PRESENT(POUT_RC)) POUT_RC=ZRC +IF(PRESENT(POUT_RI)) POUT_RI=ZRI +IF(PRESENT(POUT_TH)) POUT_TH=ZT / PEXN(:,:,:) +! ! !* 6. STORE THE BUDGET TERMS ! ---------------------- ! -IF (LBUDGET_RV) CALL BUDGET (PRVS(:,:,:) * PRHODJ(:,:,:),6,HBUNAME//'_BU_RRV',YDDDH, YDLDDH, YDMDDH) -IF (LBUDGET_RC) CALL BUDGET (PRCS(:,:,:) * PRHODJ(:,:,:),7,HBUNAME//'_BU_RRC',YDDDH, YDLDDH, YDMDDH) -IF (LBUDGET_RI) CALL BUDGET (PRIS(:,:,:) * PRHODJ(:,:,:),9,HBUNAME//'_BU_RRI',YDDDH, YDLDDH, YDMDDH) -IF (LBUDGET_TH) CALL BUDGET (PTHS(:,:,:) * PRHODJ(:,:,:),4,HBUNAME//'_BU_RTH',YDDDH, YDLDDH, YDMDDH) -! +IF(LBUDGET_TH) CALL BUDGET_STORE_END(TBUDGETS(NBUDGET_TH), TRIM(HBUNAME), PTHS(:, :, :)*PRHODJ(:, :, :)) +IF(LBUDGET_RV) CALL BUDGET_STORE_END(TBUDGETS(NBUDGET_RV), TRIM(HBUNAME), PRVS(:, :, :)*PRHODJ(:, :, :)) +IF(LBUDGET_RC) CALL BUDGET_STORE_END(TBUDGETS(NBUDGET_RC), TRIM(HBUNAME), PRCS(:, :, :)*PRHODJ(:, :, :)) +IF(LBUDGET_RI) CALL BUDGET_STORE_END(TBUDGETS(NBUDGET_RI), TRIM(HBUNAME), PRIS(:, :, :)*PRHODJ(:, :, :)) !------------------------------------------------------------------------------ ! ! IF (LHOOK) CALL DR_HOOK('ICE_ADJUST',1,ZHOOK_HANDLE) - +! CONTAINS SUBROUTINE ITERATION(PRV_IN,PRC_IN,PRI_IN,PRV_OUT,PRC_OUT,PRI_OUT) @@ -370,7 +437,7 @@ REAL, DIMENSION(:,:,:), CONTIGUOUS, INTENT(IN) :: PRI_IN ! Cloud ice m.r. to a REAL, DIMENSION(:,:,:), CONTIGUOUS, INTENT(OUT) :: PRV_OUT ! Water vapor m.r. to adjust in output REAL, DIMENSION(:,:,:), CONTIGUOUS, INTENT(OUT) :: PRC_OUT ! Cloud water m.r. to adjust in output REAL, DIMENSION(:,:,:), CONTIGUOUS, INTENT(OUT) :: PRI_OUT ! Cloud ice m.r. to adjust in output - +! !* 2.4 compute the specific heat for moist air (Cph) at t+1 SELECT CASE(KRR) @@ -401,27 +468,34 @@ IF ( OSUBG_COND ) THEN ! ! PSRC= s'rci'/Sigma_s^2 ! ZT is INOUT - - CALL CONDENSATION(IIU, IJU, IKU, IIB, IIE, IJB, IJE, IKB, IKE, KKL, & - HFRAC_ICE, & - PPABST, PZZ, ZT, PRV_IN, PRV_OUT, PRC_IN, PRC_OUT, PRI_IN, PRI_OUT, PRS, PRG, PSIGS, PMFCONV, PCLDFR, & + CALL CONDENSATION(IIU, IJU, IKU, IIB, IIE, IJB, IJE, IKB, IKE, KKL, & + HFRAC_ICE, HCONDENS, HLAMBDA3, & + PPABST, PZZ, PRHODREF, ZT, PRV_IN, PRV_OUT, PRC_IN, PRC_OUT, PRI_IN, PRI_OUT, & + PRS, PRG, PSIGS, PMFCONV, PCLDFR, & PSRCS, .TRUE., OSIGMAS, & - OCND2, PSIGQSAT, YSPP_PSIGQSAT,YSPP_ICE_CLD_WGT, & - PLV=ZLV, PLS=ZLS, PCPH=ZCPH) + OCND2, PSIGQSAT, & + PLV=ZLV, PLS=ZLS, PCPH=ZCPH, & + PHLC_HRC=PHLC_HRC, PHLC_HCF=PHLC_HCF, PHLI_HRI=PHLI_HRI, PHLI_HCF=PHLI_HCF,& + PICE_CLD_WGT=PICE_CLD_WGT) ELSE ! !* 4. ALL OR NOTHING CONDENSATION SCHEME ! FOR MIXED-PHASE CLOUD ! ----------------------------------------------- - ZSIGS(:,:,:)=0. - CALL CONDENSATION(IIU, IJU, IKU, IIB, IIE, IJB, IJE, IKB, IKE, KKL, & - HFRAC_ICE, & - PPABST, PZZ, ZT, PRV_IN, PRV_OUT, PRC_IN, PRC_OUT, PRI_IN, PRI_OUT, PRS, PRG, ZSIGS, PMFCONV, PCLDFR, & - PSRCS, .TRUE., OSIGMAS=.TRUE., & - OCND2=OCND2, PSIGQSAT=0., & - YSPP_PSIGQSAT=YSPP_PSIGQSAT, & - YSPP_ICE_CLD_WGT=YSPP_ICE_CLD_WGT, & - PLV=ZLV, PLS=ZLS, PCPH=ZCPH) + ! + ZSIGS(:,:,:)=0. + ZSIGQSAT(:,:)=0. + !We use ZSRCS because in Méso-NH, PSRCS can be a zero-length array in this case + !ZT is INOUT + CALL CONDENSATION(IIU, IJU, IKU, IIB, IIE, IJB, IJE, IKB, IKE, KKL, & + HFRAC_ICE, HCONDENS, HLAMBDA3, & + PPABST, PZZ, PRHODREF, ZT, PRV_IN, PRV_OUT, PRC_IN, PRC_OUT, PRI_IN, PRI_OUT, & + PRS, PRG, ZSIGS, PMFCONV, PCLDFR, & + ZSRCS, .TRUE., OSIGMAS=.TRUE., & + OCND2=OCND2, PSIGQSAT=ZSIGQSAT, & + PLV=ZLV, PLS=ZLS, PCPH=ZCPH, & + PHLC_HRC=PHLC_HRC, PHLC_HCF=PHLC_HCF, PHLI_HRI=PHLI_HRI, PHLI_HCF=PHLI_HCF,& + PICE_CLD_WGT=PICE_CLD_WGT) ENDIF END SUBROUTINE ITERATION diff --git a/src/common/micro/mode_compute_frac_ice.F90 b/src/common/micro/mode_compute_frac_ice.F90 index 6ebe4a28d..ba2dffd7f 100644 --- a/src/common/micro/mode_compute_frac_ice.F90 +++ b/src/common/micro/mode_compute_frac_ice.F90 @@ -37,110 +37,57 @@ !! P. Wautelet 10/04/2019: replace ABORT and STOP calls by Print_msg !! R. El Khatib 24-Aug-2021 Optimization by cache re-use + assume data is contiguous !! S. Riette Jan-2022 Merge the 3 procedures in one module + array shape declaration +! R. El Khatib / S. Riette Jan-2022 written as an elemental subroutine ! ! -USE PARKIND1, ONLY : JPRB -USE YOMHOOK , ONLY : LHOOK, DR_HOOK ! -INTERFACE COMPUTE_FRAC_ICE - MODULE PROCEDURE COMPUTE_FRAC_ICE1D, COMPUTE_FRAC_ICE2D, COMPUTE_FRAC_ICE3D -END INTERFACE COMPUTE_FRAC_ICE -CONTAINS -!! ========== -!! 3D version -!! ========== -! -SUBROUTINE COMPUTE_FRAC_ICE3D(KIT, KJT, KKT, HFRAC_ICE, PFRAC_ICE, PT) -IMPLICIT NONE -INTEGER, INTENT(IN) :: KIT, KJT, KKT -CHARACTER(LEN=1), INTENT(IN) :: HFRAC_ICE ! scheme to use -REAL, TARGET, DIMENSION(KIT,KJT,KKT), INTENT(IN) :: PT ! Temperature -REAL, TARGET, DIMENSION(KIT,KJT,KKT), INTENT(INOUT) :: PFRAC_ICE ! Ice fraction (1 for ice only, 0 for liquid only) -REAL, POINTER, DIMENSION(:) :: ZT, ZFRAC_ICE -REAL(KIND=JPRB) :: ZHOOK_HANDLE -! -IF (LHOOK) CALL DR_HOOK('COMPUTE_FRAC_ICE3D',0,ZHOOK_HANDLE) -! -ZT(1:SIZE(PT))=>PT -ZFRAC_ICE(1:SIZE(PFRAC_ICE))=>PFRAC_ICE -CALL COMPUTE_FRAC_ICE1D(SIZE(PT), HFRAC_ICE,ZFRAC_ICE,ZT) -! -IF (LHOOK) CALL DR_HOOK('COMPUTE_FRAC_ICE3D',1,ZHOOK_HANDLE) -END SUBROUTINE COMPUTE_FRAC_ICE3D -! -! -!! ========== -!! 2D version -!! ========== -! -SUBROUTINE COMPUTE_FRAC_ICE2D(KIT, KJT, HFRAC_ICE, PFRAC_ICE, PT) +!****************** Don't use drHook !!! + + + IMPLICIT NONE -INTEGER, INTENT(IN) :: KIT, KJT -CHARACTER(LEN=1), INTENT(IN) :: HFRAC_ICE ! scheme to use -REAL, TARGET, DIMENSION(KIT,KJT), INTENT(IN) :: PT ! Temperature -REAL, TARGET, DIMENSION(KIT,KJT), INTENT(INOUT) :: PFRAC_ICE ! Ice fraction (1 for ice only, 0 for liquid only) -REAL, POINTER, DIMENSION(:) :: ZT, ZFRAC_ICE -REAL(KIND=JPRB) :: ZHOOK_HANDLE -! -IF (LHOOK) CALL DR_HOOK('COMPUTE_FRAC_ICE2D',0,ZHOOK_HANDLE) -! -ZT(1:SIZE(PT))=>PT -ZFRAC_ICE(1:SIZE(PFRAC_ICE))=>PFRAC_ICE -CALL COMPUTE_FRAC_ICE1D(SIZE(PT), HFRAC_ICE,ZFRAC_ICE,ZT) -! -IF (LHOOK) CALL DR_HOOK('COMPUTE_FRAC_ICE2D',1,ZHOOK_HANDLE) -END SUBROUTINE COMPUTE_FRAC_ICE2D -! -! -!! ========== -!! 1D version -!! ========== -! -SUBROUTINE COMPUTE_FRAC_ICE1D(KIT, HFRAC_ICE, PFRAC_ICE, PT) +CONTAINS + +ELEMENTAL SUBROUTINE COMPUTE_FRAC_ICE(HFRAC_ICE, PFRAC_ICE, PT, KERR) !! -------------------------------------------------------------------------- ! 0. DECLARATIONS ! ------------ ! USE MODD_NEB, ONLY : XTMINMIX, XTMAXMIX USE MODD_CST, ONLY : XTT -USE MODE_MSG, ONLY : PRINT_MSG, NVERB_FATAL ! IMPLICIT NONE ! ! !* 0.1 declarations of arguments ! -INTEGER, INTENT(IN) :: KIT -CHARACTER(LEN=1), INTENT(IN) :: HFRAC_ICE ! scheme to use -REAL, TARGET, DIMENSION(KIT), INTENT(IN) :: PT ! Temperature -REAL, TARGET, DIMENSION(KIT), INTENT(INOUT) :: PFRAC_ICE ! Ice fraction (1 for ice only, 0 for liquid only) +CHARACTER(LEN=1), INTENT(IN) :: HFRAC_ICE ! scheme to use +REAL, INTENT(IN) :: PT ! Temperature +REAL, INTENT(INOUT) :: PFRAC_ICE ! Ice fraction (1 for ice only, 0 for liquid only) +INTEGER, INTENT(OUT) :: KERR ! Error code in return ! ! 0.2 declaration of local variables ! -REAL(KIND=JPRB) :: ZHOOK_HANDLE -! -IF (LHOOK) CALL DR_HOOK('COMPUTE_FRAC_ICE1D',0,ZHOOK_HANDLE) !------------------------------------------------------------------------ ! 1. Compute FRAC_ICE ! +KERR=0 SELECT CASE(HFRAC_ICE) CASE ('T') !using Temperature - PFRAC_ICE(:) = MAX( 0., MIN(1., (( XTMAXMIX - PT(:) ) / ( XTMAXMIX - XTMINMIX )) ) ) ! freezing interval + PFRAC_ICE = MAX( 0., MIN(1., (( XTMAXMIX - PT ) / ( XTMAXMIX - XTMINMIX )) ) ) ! freezing interval CASE ('O') !using Temperature with old formulae - PFRAC_ICE(:) = MAX( 0., MIN(1., (( XTT - PT(:) ) / 40.) ) ) ! freezing interval + PFRAC_ICE = MAX( 0., MIN(1., (( XTT - PT ) / 40.) ) ) ! freezing interval CASE ('N') !No ice - PFRAC_ICE(:) = 0. + PFRAC_ICE = 0. CASE ('S') !Same as previous ! (almost) nothing to do - PFRAC_ICE(:) = MAX( 0., MIN(1., PFRAC_ICE(:) ) ) + PFRAC_ICE = MAX( 0., MIN(1., PFRAC_ICE ) ) CASE DEFAULT - CALL PRINT_MSG(NVERB_FATAL, 'GEN', 'COMPUTE_FRAC_ICE', 'invalid option for HFRAC_ICE='//HFRAC_ICE) + KERR=1 END SELECT ! -IF (LHOOK) CALL DR_HOOK('COMPUTE_FRAC_ICE1D',1,ZHOOK_HANDLE) -! -END SUBROUTINE COMPUTE_FRAC_ICE1D +END SUBROUTINE COMPUTE_FRAC_ICE ! END MODULE MODE_COMPUTE_FRAC_ICE diff --git a/src/arome/micro/icecloud.F90 b/src/common/micro/mode_icecloud.F90 similarity index 98% rename from src/arome/micro/icecloud.F90 rename to src/common/micro/mode_icecloud.F90 index abd63cfde..43b394d21 100644 --- a/src/arome/micro/icecloud.F90 +++ b/src/common/micro/mode_icecloud.F90 @@ -1,3 +1,6 @@ +MODULE MODE_ICECLOUD +IMPLICIT NONE +CONTAINS SUBROUTINE ICECLOUD & ! Input : & ( NP,PP,PZ,PDZ,PT,PR,PTSTEP,PPBLH,PWCLD,XW2D, & @@ -7,7 +10,7 @@ SUBROUTINE ICECLOUD & USE PARKIND1, ONLY : JPRB USE YOMHOOK , ONLY : LHOOK, DR_HOOK USE MODD_CST,ONLY : XCPD,XCPV,XLVTT,XLSTT,XG,XRD,XEPSILO - USE MODD_TIWMX + USE MODE_TIWMX, ONLY: ESATW, ESATI USE MODI_TIWMX IMPLICIT NONE !----------------------------------------------------------------------- @@ -202,3 +205,4 @@ ENDDO IF (LHOOK) CALL DR_HOOK('ICECLOUD',1,ZHOOK_HANDLE) END SUBROUTINE ICECLOUD +END MODULE MODE_ICECLOUD diff --git a/src/arome/micro/modd_tiwmx.F90 b/src/common/micro/mode_tiwmx.F90 similarity index 97% rename from src/arome/micro/modd_tiwmx.F90 rename to src/common/micro/mode_tiwmx.F90 index 4ff4a9836..dbfbdc712 100644 --- a/src/arome/micro/modd_tiwmx.F90 +++ b/src/common/micro/mode_tiwmx.F90 @@ -1,9 +1,9 @@ !@no_insert_drhook ! ######spl - MODULE MODD_TIWMX + MODULE MODE_TIWMX ! ############### ! -!!**** *MODD_TIWMX* - +!!**** *MODE_TIWMX* - !! !! PURPOSE !! ------- @@ -110,4 +110,4 @@ CONTAINS REAL,INTENT(IN) :: TT REDIN = REDINTAB(NINT(XNDEGR*TT)) END FUNCTION -END MODULE MODD_TIWMX +END MODULE MODE_TIWMX diff --git a/src/arome/micro/tiwmx_tab.F90 b/src/common/micro/mode_tiwmx_tab.F90 similarity index 85% rename from src/arome/micro/tiwmx_tab.F90 rename to src/common/micro/mode_tiwmx_tab.F90 index cc04ed7c7..e8e42bcd6 100644 --- a/src/arome/micro/tiwmx_tab.F90 +++ b/src/common/micro/mode_tiwmx_tab.F90 @@ -1,3 +1,6 @@ +MODULE MODE_TIWMX_TAB +IMPLICIT NONE +CONTAINS FUNCTION TIWMX_TAB(P,T,QR,FICE,QRSN,RS,EPS) ! Purpose: (*) @@ -46,7 +49,7 @@ FUNCTION TIWMX_TAB(P,T,QR,FICE,QRSN,RS,EPS) ! 1.1 MODULES USED USE MODD_CST, ONLY : XEPSILO, XCPD, XLSTT, XLVTT - USE MODD_TIWMX, ONLY : ESATI, ESATW, DESDTI, DESDTW + USE MODE_TIWMX, ONLY : ESATI, ESATW, DESDTI, DESDTW USE YOMHOOK , ONLY : LHOOK, DR_HOOK USE PARKIND1, ONLY : JPRB @@ -103,28 +106,4 @@ FUNCTION TIWMX_TAB(P,T,QR,FICE,QRSN,RS,EPS) IF (LHOOK) CALL DR_HOOK('TIWMX_TAB',1,ZHOOK_HANDLE) END FUNCTION TIWMX_TAB - - -FUNCTION QSATMX_TAB(P,T,FICE) - - USE PARKIND1, ONLY : JPRB - USE MODD_CST ,ONLY : XEPSILO - USE MODD_TIWMX, ONLY : ESATI,ESATW - - IMPLICIT NONE - - REAL :: QSATMX_TAB - REAL, INTENT(IN) :: P,T,FICE - - REAL :: ZES - - ZES = ESATI(T)*FICE + ESATW(T)*(1.-FICE) - IF(ZES >= P)THEN ! temp > boiling point, condensation not possible. - ! Then this function lacks physical meaning, - ! here set to one - QSATMX_TAB = 1. - ELSE - QSATMX_TAB = XEPSILO*ZES/(P-ZES) !r - ENDIF - -END FUNCTION QSATMX_TAB +END MODULE MODE_TIWMX_TAB diff --git a/src/arome/modset_Ryad/mpa/micro/module/modi_condensation.F90 b/src/common/micro/modi_condensation.F90 similarity index 73% rename from src/arome/modset_Ryad/mpa/micro/module/modi_condensation.F90 rename to src/common/micro/modi_condensation.F90 index b46a1960a..7163e8208 100644 --- a/src/arome/modset_Ryad/mpa/micro/module/modi_condensation.F90 +++ b/src/common/micro/modi_condensation.F90 @@ -5,11 +5,12 @@ INTERFACE ! SUBROUTINE CONDENSATION( KIU, KJU, KKU, KIB, KIE, KJB, KJE, KKB, KKE, KKL,& - HFRAC_ICE, & - PPABS, PZZ, PT, PRV_IN, PRV_OUT, PRC_IN, PRC_OUT, PRI_IN, PRI_OUT, PRS, PRG, PSIGS, PMFCONV, PCLDFR, PSIGRC, OUSERI,& + HFRAC_ICE, HCONDENS, HLAMBDA3, & + PPABS, PZZ, PRHODREF, PT, PRV_IN, PRV_OUT, PRC_IN, PRC_OUT, PRI_IN, PRI_OUT, & + PRS, PRG, PSIGS, PMFCONV, PCLDFR, PSIGRC, OUSERI,& OSIGMAS, OCND2, PSIGQSAT, & - YSPP_PSIGQSAT, YSPP_ICE_CLD_WGT,& - PLV, PLS, PCPH) + PLV, PLS, PCPH, & + PHLC_HRC, PHLC_HCF, PHLI_HRI, PHLI_HCF, PICE_CLD_WGT) ! USE MODD_SPP_TYPE ! @@ -23,9 +24,12 @@ INTEGER, INTENT(IN) :: KJE ! value of the last point INTEGER, INTENT(IN) :: KKB ! value of the first point in z INTEGER, INTENT(IN) :: KKE ! value of the last point in z INTEGER, INTENT(IN) :: KKL ! +1 if grid goes from ground to atmosphere top, -1 otherwise -CHARACTER*1, INTENT(IN) :: HFRAC_ICE +CHARACTER(LEN=1), INTENT(IN) :: HFRAC_ICE +CHARACTER(LEN=4), INTENT(IN) :: HCONDENS +CHARACTER(LEN=*), INTENT(IN) :: HLAMBDA3 ! formulation for lambda3 coeff REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PPABS ! pressure (Pa) REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PZZ ! height of model levels (m) +REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PRHODREF REAL, DIMENSION(KIU,KJU,KKU), INTENT(INOUT) :: PT ! grid scale T (K) REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PRV_IN ! grid scale water vapor mixing ratio (kg/kg) in input REAL, DIMENSION(KIU,KJU,KKU), INTENT(OUT) :: PRV_OUT! grid scale water vapor mixing ratio (kg/kg) in output @@ -36,10 +40,8 @@ LOGICAL, INTENT(IN) :: OSIGMAS! use present global Sigma ! or that from turbulence scheme LOGICAL, INTENT(IN) :: OCND2 ! logical switch to sparate liquid and ice ! more rigid (DEFALT value : .FALSE.) - -REAL, INTENT(IN) :: PSIGQSAT ! use an extra "qsat" variance contribution (OSIGMAS case) +REAL, DIMENSION(KIU,KJU), INTENT(IN) :: PSIGQSAT ! use an extra "qsat" variance contribution (OSIGMAS case) ! multiplied by PSIGQSAT -TYPE(TSPP_CONFIG_MPA), INTENT(IN) :: YSPP_PSIGQSAT,YSPP_ICE_CLD_WGT REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PRC_IN ! grid scale r_c mixing ratio (kg/kg) in input REAL, DIMENSION(KIU,KJU,KKU), INTENT(OUT) :: PRC_OUT! grid scale r_c mixing ratio (kg/kg) in output REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PRI_IN ! grid scale r_i (kg/kg) in input @@ -50,9 +52,15 @@ REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PSIGS ! Sigma_s from turbulence REAL, DIMENSION(:,:,:), INTENT(IN) :: PMFCONV! convective mass flux (kg /s m^2) REAL, DIMENSION(KIU,KJU,KKU), INTENT(OUT) :: PCLDFR ! cloud fraction REAL, DIMENSION(KIU,KJU,KKU), INTENT(OUT) :: PSIGRC ! s r_c / sig_s^2 -REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(IN) :: PLV -REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(IN) :: PLS -REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(IN) :: PCPH + +REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(IN) :: PLV ! Latent heat L_v +REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(IN) :: PLS ! Latent heat L_s +REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(IN) :: PCPH ! Specific heat C_ph +REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(OUT) :: PHLC_HRC +REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(OUT) :: PHLC_HCF ! cloud fraction +REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(OUT) :: PHLI_HRI +REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(OUT) :: PHLI_HCF +REAL, DIMENSION(KIU,KJU), OPTIONAL, INTENT(IN) :: PICE_CLD_WGT END SUBROUTINE CONDENSATION ! END INTERFACE diff --git a/src/arome/micro/modi_ice_adjust.F90 b/src/common/micro/modi_ice_adjust.F90 similarity index 56% rename from src/arome/micro/modi_ice_adjust.F90 rename to src/common/micro/modi_ice_adjust.F90 index 0bcf1be1b..e949147f5 100644 --- a/src/arome/micro/modi_ice_adjust.F90 +++ b/src/common/micro/modi_ice_adjust.F90 @@ -4,52 +4,59 @@ ! INTERFACE ! - SUBROUTINE ICE_ADJUST (KKA, KKU, KKL, KRR, HFRAC_ICE, & - HBUNAME, OSUBG_COND, OSIGMAS, OCND2, & + SUBROUTINE ICE_ADJUST (KKA, KKU, KKL, KRR, HFRAC_ICE, HCONDENS, HLAMBDA3,& + HBUNAME, OSUBG_COND, OSIGMAS, OCND2, HSUBG_MF_PDF,& PTSTEP, PSIGQSAT, & - PRHODJ, PEXNREF, PSIGS, PMFCONV, PPABST, PZZ, & + PRHODJ, PEXNREF, PRHODREF, PSIGS, PMFCONV, & + PPABST, PZZ, & PEXN, PCF_MF, PRC_MF, PRI_MF, & - PRV, PRC, PRVS, PRCS, PTH, PTHS, PSRCS, PCLDFR , & - PRR, PRI, PRIS, PRS, PRG, & - PRH, POUT_RV, POUT_RC, POUT_RI, POUT_TH, & - YSPP_PSIGQSAT, YSPP_ICE_CLD_WGT, & - YDDDH, YDLDDH, YDMDDH ) + PRV, PRC, PRVS, PRCS, PTH, PTHS, PSRCS, PCLDFR, & + PRR, PRI, PRIS, PRS, PRG, PRH, & + POUT_RV, POUT_RC, POUT_RI, POUT_TH, & + PHLC_HRC, PHLC_HCF, PHLI_HRI, PHLI_HCF, & + TBUDGETS, KBUDGETS, & + PICE_CLD_WGT) +USE MODD_BUDGET, ONLY: TBUDGETDATA +IMPLICIT NONE ! -USE MODD_SPP_TYPE ! -USE DDH_MIX, ONLY : TYP_DDH -USE YOMLDDH, ONLY : TLDDH -USE YOMMDDH, ONLY : TMDDH +!* 0.1 Declarations of dummy arguments : ! -INTEGER, INTENT(IN) :: KKA !near ground array index -INTEGER, INTENT(IN) :: KKU !uppest atmosphere array index -INTEGER, INTENT(IN) :: KKL !vert. levels type 1=MNH -1=ARO +! +INTEGER, INTENT(IN) :: KKA !near ground array index +INTEGER, INTENT(IN) :: KKU !uppest atmosphere array index +INTEGER, INTENT(IN) :: KKL !vert. levels type 1=MNH -1=ARO INTEGER, INTENT(IN) :: KRR ! Number of moist variables -CHARACTER*1, INTENT(IN) :: HFRAC_ICE -CHARACTER*4, INTENT(IN) :: HBUNAME ! Name of the budget -LOGICAL, INTENT(IN) :: OSUBG_COND ! Switch for Subgrid +CHARACTER(LEN=1), INTENT(IN) :: HFRAC_ICE +CHARACTER(LEN=80), INTENT(IN) :: HCONDENS +CHARACTER(LEN=4), INTENT(IN) :: HLAMBDA3 ! formulation for lambda3 coeff +CHARACTER(LEN=4), INTENT(IN) :: HBUNAME ! Name of the budget +LOGICAL, INTENT(IN) :: OSUBG_COND ! Switch for Subgrid ! Condensation -LOGICAL :: OSIGMAS ! Switch for Sigma_s: +LOGICAL :: OSIGMAS ! Switch for Sigma_s: ! use values computed in CONDENSATION ! or that from turbulence scheme -LOGICAL, INTENT(IN) :: OCND2 ! logical switch to separate liquid +LOGICAL :: OCND2 ! logical switch to sparate liquid ! and ice + ! more rigid (DEFALT value : .FALSE.) +CHARACTER(LEN=80), INTENT(IN) :: HSUBG_MF_PDF REAL, INTENT(IN) :: PTSTEP ! Double Time step ! (single if cold start) -REAL, INTENT(IN) :: PSIGQSAT ! coeff applied to qsat variance contribution +REAL, DIMENSION(:,:), INTENT(IN) :: PSIGQSAT ! coeff applied to qsat variance contribution ! REAL, DIMENSION(:,:,:), INTENT(IN) :: PRHODJ ! Dry density * Jacobian REAL, DIMENSION(:,:,:), INTENT(IN) :: PEXNREF ! Reference Exner function +REAL, DIMENSION(:,:,:), INTENT(IN) :: PRHODREF ! REAL, DIMENSION(:,:,:), INTENT(IN) :: PSIGS ! Sigma_s at time t REAL, DIMENSION(:,:,:), INTENT(IN) :: PMFCONV ! convective mass flux -REAL, DIMENSION(:,:,:), INTENT(IN) :: PPABST ! Absolute Pressure at t +REAL, DIMENSION(:,:,:), INTENT(IN) :: PPABST ! Absolute Pressure at t REAL, DIMENSION(:,:,:), INTENT(IN) :: PZZ ! height of model layer REAL, DIMENSION(:,:,:), INTENT(IN) :: PEXN ! Exner function ! -REAL, DIMENSION(:,:,:), INTENT(IN) :: PCF_MF! Convective Mass Flux Cloud fraction -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRI_MF! Convective Mass Flux ice mixing ratio +REAL, DIMENSION(:,:,:), INTENT(IN) :: PCF_MF! Convective Mass Flux Cloud fraction REAL, DIMENSION(:,:,:), INTENT(IN) :: PRC_MF! Convective Mass Flux liquid mixing ratio +REAL, DIMENSION(:,:,:), INTENT(IN) :: PRI_MF! Convective Mass Flux ice mixing ratio ! REAL, DIMENSION(:,:,:), INTENT(IN) :: PRV ! Water vapor m.r. to adjust REAL, DIMENSION(:,:,:), INTENT(IN) :: PRC ! Cloud water m.r. to adjust @@ -60,22 +67,25 @@ REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PTHS ! Theta source REAL, DIMENSION(:,:,:), INTENT(OUT) :: PSRCS ! Second-order flux ! s'rc'/2Sigma_s2 at time t+1 ! multiplied by Lambda_3 -REAL, DIMENSION(:,:,:), INTENT(OUT) :: PCLDFR ! Cloud fraction -REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PRIS ! Cloud ice m.r. at t+1 -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRR ! Rain water m.r. to adjust -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRI ! Cloud ice m.r. to adjust -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRS ! Aggregate m.r. to adjust -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRG ! Graupel m.r. to adjust +REAL, DIMENSION(:,:,:), INTENT(OUT) :: PCLDFR ! Cloud fraction +! +REAL, DIMENSION(:,:,:), INTENT(INOUT):: PRIS ! Cloud ice m.r. at t+1 +REAL, DIMENSION(:,:,:), INTENT(IN) :: PRR ! Rain water m.r. to adjust +REAL, DIMENSION(:,:,:), INTENT(IN) :: PRI ! Cloud ice m.r. to adjust +REAL, DIMENSION(:,:,:), INTENT(IN) :: PRS ! Aggregate m.r. to adjust +REAL, DIMENSION(:,:,:), INTENT(IN) :: PRG ! Graupel m.r. to adjust REAL, DIMENSION(:,:,:), OPTIONAL, INTENT(IN) :: PRH ! Hail m.r. to adjust REAL, DIMENSION(:,:,:), OPTIONAL, INTENT(OUT) :: POUT_RV ! Adjusted value REAL, DIMENSION(:,:,:), OPTIONAL, INTENT(OUT) :: POUT_RC ! Adjusted value REAL, DIMENSION(:,:,:), OPTIONAL, INTENT(OUT) :: POUT_RI ! Adjusted value REAL, DIMENSION(:,:,:), OPTIONAL, INTENT(OUT) :: POUT_TH ! Adjusted value -TYPE(TSPP_CONFIG_MPA), INTENT(IN) :: YSPP_PSIGQSAT,YSPP_ICE_CLD_WGT -TYPE(TYP_DDH) , INTENT(INOUT) :: YDDDH -TYPE(TLDDH) , INTENT(IN) :: YDLDDH -TYPE(TMDDH) , INTENT(IN) :: YDMDDH -! +REAL, DIMENSION(:,:,:), OPTIONAL, INTENT(OUT) :: PHLC_HRC +REAL, DIMENSION(:,:,:), OPTIONAL, INTENT(OUT) :: PHLC_HCF +REAL, DIMENSION(:,:,:), OPTIONAL, INTENT(OUT) :: PHLI_HRI +REAL, DIMENSION(:,:,:), OPTIONAL, INTENT(OUT) :: PHLI_HCF +TYPE(TBUDGETDATA), DIMENSION(KBUDGETS), INTENT(INOUT) :: TBUDGETS +INTEGER, INTENT(IN) :: KBUDGETS +REAL, DIMENSION(:,:), OPTIONAL, INTENT(IN) :: PICE_CLD_WGT ! END SUBROUTINE ICE_ADJUST ! diff --git a/src/common/micro/rain_ice.F90 b/src/common/micro/rain_ice.F90 index b75c0de45..0fa4b7bd7 100644 --- a/src/common/micro/rain_ice.F90 +++ b/src/common/micro/rain_ice.F90 @@ -179,7 +179,7 @@ USE PARKIND1, ONLY : JPRB USE YOMHOOK , ONLY : LHOOK, DR_HOOK -USE MODD_BUDGET, ONLY: TBUDGETDATA, LBU_ENABLE, & +USE MODD_BUDGET, ONLY: TBUDGETDATA, LBU_ENABLE, & & LBUDGET_TH, LBUDGET_RV, LBUDGET_RC, LBUDGET_RR, LBUDGET_RI, LBUDGET_RS, LBUDGET_RG, LBUDGET_RH, & & NBUDGET_TH, NBUDGET_RV, NBUDGET_RC, NBUDGET_RR, NBUDGET_RI, NBUDGET_RS, NBUDGET_RG, NBUDGET_RH USE MODD_CST, ONLY: XCI, XCL, XCPD, XCPV, XLSTT, XLVTT, XTT, XRHOLW @@ -198,7 +198,7 @@ USE MODD_FIELDS_ADDRESS, ONLY : & ! common fields adress & IRH ! Hail USE MODE_BUDGET, ONLY: BUDGET_STORE_ADD, BUDGET_STORE_INIT, BUDGET_STORE_END -USE MODE_ll +USE MODE_ll, ONLY: GET_INDICE_ll USE MODE_MSG, ONLY: PRINT_MSG, NVERB_FATAL USE MODE_ICE4_RAINFR_VERT, ONLY: ICE4_RAINFR_VERT @@ -444,14 +444,12 @@ IF(OCND2) THEN END IF IF(KPROMA /= KSIZE) THEN CALL PRINT_MSG(NVERB_FATAL, 'GEN', 'RAIN_ICE', 'For now, KPROMA must be equal to KSIZE, see code for explanation') - ! 2 issues - ! * Microphyscs was optimized by introducing chunks of KPROMA size - ! Thus, in ice4_tendencies, the 1D array represent only a fraction of the points where microphisical species are present - ! We cannot rebuild the entire 3D arrays in the subroutine, so we cannot call ice4_rainfr_vert in it - ! A solution would be to suppress optimisation in this case by setting KPROMA=KSIZE in rain_ice - ! Another solution would be to compute column by column? - ! Another one would be to cut tendencies in 3 parts: before rainfr_vert, rainfr_vert, after rainfr_vert - ! * When chuncks are used, result is different + ! Microphyscs was optimized by introducing chunks of KPROMA size + ! Thus, in ice4_tendencies, the 1D array represent only a fraction of the points where microphisical species are present + ! We cannot rebuild the entire 3D arrays in the subroutine, so we cannot call ice4_rainfr_vert in it + ! A solution would be to suppress optimisation in this case by setting KPROMA=KSIZE in rain_ice + ! Another solution would be to compute column by column? + ! Another one would be to cut tendencies in 3 parts: before rainfr_vert, rainfr_vert, after rainfr_vert ENDIF ! !* 1. COMPUTE THE LOOP BOUNDS diff --git a/src/mesonh/micro/condensation.f90 b/src/mesonh/micro/condensation.f90 deleted file mode 100644 index d0c88e2a3..000000000 --- a/src/mesonh/micro/condensation.f90 +++ /dev/null @@ -1,513 +0,0 @@ -!MNH_LIC Copyright 2002-2021 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 for details. version 1. -!----------------------------------------------------------------- -! ######spl - MODULE MODI_CONDENSATION -! ######################## -! -INTERFACE -! - SUBROUTINE CONDENSATION( KIU, KJU, KKU, KIB, KIE, KJB, KJE, KKB, KKE, KKL,& - HFRAC_ICE, HCONDENS, HLAMBDA3, & - PPABS, PZZ, PRHODREF, PT, PRV, PRC, PRI, PRS, PRG, PSIGS, PMFCONV, PCLDFR, PSIGRC, OUSERI,& - OSIGMAS, PSIGQSAT, PLV, PLS, PCPH, PHLC_HRC, PHLC_HCF, PHLI_HRI, PHLI_HCF) -! -INTEGER, INTENT(IN) :: KIU ! horizontal dimension in x -INTEGER, INTENT(IN) :: KJU ! horizontal dimension in y -INTEGER, INTENT(IN) :: KKU ! vertical dimension -INTEGER, INTENT(IN) :: KIB ! value of the first point in x -INTEGER, INTENT(IN) :: KIE ! value of the last point in x -INTEGER, INTENT(IN) :: KJB ! value of the first point in y -INTEGER, INTENT(IN) :: KJE ! value of the last point in y -INTEGER, INTENT(IN) :: KKB ! value of the first point in z -INTEGER, INTENT(IN) :: KKE ! value of the last point in z -INTEGER, INTENT(IN) :: KKL ! +1 if grid goes from ground to atmosphere top, -1 otherwise -CHARACTER(len=1), INTENT(IN) :: HFRAC_ICE -CHARACTER(len=4), INTENT(IN) :: HCONDENS -CHARACTER(len=*), INTENT(IN) :: HLAMBDA3 ! formulation for lambda3 coeff -REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PPABS ! pressure (Pa) -REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PZZ ! height of model levels (m) -REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PRHODREF -REAL, DIMENSION(KIU,KJU,KKU), INTENT(INOUT) :: PT ! grid scale T (K) -REAL, DIMENSION(KIU,KJU,KKU), INTENT(INOUT) :: PRV ! grid scale water vapor mixing ratio (kg/kg) -REAL, DIMENSION(KIU,KJU,KKU), INTENT(INOUT) :: PRC ! grid scale r_c mixing ratio (kg/kg) -REAL, DIMENSION(KIU,KJU,KKU), INTENT(INOUT) :: PRI ! grid scale r_i (kg/kg) -REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PRS ! grid scale mixing ration of snow (kg/kg) -REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PRG ! grid scale mixing ration of graupel (kg/kg) -REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PSIGS ! Sigma_s from turbulence scheme -REAL, DIMENSION(:,:,:), INTENT(IN) :: PMFCONV! convective mass flux (kg /s m^2) -REAL, DIMENSION(KIU,KJU,KKU), INTENT(OUT) :: PCLDFR ! cloud fraction -REAL, DIMENSION(KIU,KJU,KKU), INTENT(OUT) :: PSIGRC ! s r_c / sig_s^2 -LOGICAL, INTENT(IN) :: OUSERI ! logical switch to compute both - ! liquid and solid condensate (OUSERI=.TRUE.) - ! or only solid condensate (OUSERI=.FALSE.) -LOGICAL, INTENT(IN) :: OSIGMAS! use present global Sigma_s values - ! or that from turbulence scheme -REAL, INTENT(IN) :: PSIGQSAT ! use an extra "qsat" variance contribution (OSIGMAS case) -REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(IN) :: PLV -REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(IN) :: PLS -REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(IN) :: PCPH -REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(OUT) :: PHLC_HRC !cloud water content in precipitating part -REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(OUT) :: PHLC_HCF !precipitating part -REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(OUT) :: PHLI_HRI ! -REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(OUT) :: PHLI_HCF ! - -END SUBROUTINE CONDENSATION -! -END INTERFACE -! -END MODULE MODI_CONDENSATION -! ######spl - SUBROUTINE CONDENSATION( KIU, KJU, KKU, KIB, KIE, KJB, KJE, KKB, KKE, KKL, & - HFRAC_ICE, HCONDENS, HLAMBDA3, & - PPABS, PZZ, PRHODREF, PT, PRV, PRC, PRI, PRS, PRG, PSIGS, PMFCONV, PCLDFR, PSIGRC, OUSERI,& - OSIGMAS, PSIGQSAT, PLV, PLS, PCPH, PHLC_HRC, PHLC_HCF, PHLI_HRI, PHLI_HCF ) -! ################################################################################ -! -!! -!! PURPOSE -!! ------- -!!** Routine to diagnose cloud fraction, liquid and ice condensate mixing ratios -!! and s'rl'/sigs^2 -!! -!! -!!** METHOD -!! ------ -!! Based on the large-scale fields of temperature, water vapor, and possibly -!! liquid and solid condensate, the conserved quantities r_t and h_l are constructed -!! and then fractional cloudiness, liquid and solid condensate is diagnosed. -!! -!! The total variance is parameterized as the sum of stratiform/turbulent variance -!! and a convective variance. -!! The turbulent variance is parameterized as a function of first-order moments, and -!! the convective variance is modelled as a function of the convective mass flux -!! (units kg/s m^2) as provided by the mass flux convection scheme. -!! -!! Nota: if the host model does not use prognostic values for liquid and solid condensate -!! or does not provide a convective mass flux, put all these values to zero. -!! -!! -!! EXTERNAL -!! -------- -!! INI_CST -!! -!! IMPLICIT ARGUMENTS -!! ------------------ -!! Module MODD_CST : contains physical constants -!! -!! REFERENCE -!! --------- -!! Chaboureau J.P. and P. Bechtold (J. Atmos. Sci. 2002) -!! -!! AUTHOR -!! ------ -!! P. BECHTOLD * Laboratoire d'Aerologie * -!! -!! MODIFICATIONS -!! ------------- -!! Original: 31.1.2002 -!! modified : 21.3.2002 -!! S.Malardel : 05.2006 : Correction sur le calcul de la fonction de -!! Bougeault F2 -!! W. de Rooy: 06-06-2010: Modification in the statistical cloud scheme -!! more specifically adding a variance term -!! following ideas of Lenderink & Siebesma 2002 -!! and adding a height dependence -!! S. Riette, 18 May 2010 : PSIGQSAT is added -!! S. Riette, 11 Oct 2011 : MIN function in PDF for continuity -!! modification of minimum value for Rc+Ri to create cloud and minimum value for sigma -!! Use of guess point as a starting point instead of liquid point -!! Better computation of ZCPH and dRsat/dT -!! Set ZCOND to zero if PCLDFR==0 -!! Safety limitation to .99*Pressure for saturation vapour pressure -!! 2012-02 Y. Seity, add possibility to run with reversed vertical levels -!! 2015 C.Lac Change min value of ZSIGMA to be in agreement with AROME -!! 2016 G.Delautier Restore min value of ZSIGMA (instability) -!! 2016 S.Riette Change INQ1 -!! 2016-11 S. Riette: use HFRAC_ICE, output adjusted state -!------------------------------------------------------------------------------- -! -!* 0. DECLARATIONS -! ------------ -! -USE MODD_CST -USE MODD_PARAMETERS -USE MODD_RAIN_ICE_PARAM, ONLY : XCRIAUTC, XCRIAUTI, XACRIAUTI, XBCRIAUTI -! -use mode_msg -! -USE MODE_COMPUTE_FRAC_ICE, ONLY : COMPUTE_FRAC_ICE -! -! -IMPLICIT NONE -! -!* 0.1 Declarations of dummy arguments : -! -! -INTEGER, INTENT(IN) :: KIU ! horizontal dimension in x -INTEGER, INTENT(IN) :: KJU ! horizontal dimension in y -INTEGER, INTENT(IN) :: KKU ! vertical dimension -INTEGER, INTENT(IN) :: KIB ! value of the first point in x -INTEGER, INTENT(IN) :: KIE ! value of the last point in x -INTEGER, INTENT(IN) :: KJB ! value of the first point in y -INTEGER, INTENT(IN) :: KJE ! value of the last point in y -INTEGER, INTENT(IN) :: KKB ! value of the first point in z -INTEGER, INTENT(IN) :: KKE ! value of the last point in z -INTEGER, INTENT(IN) :: KKL ! +1 if grid goes from ground to atmosphere top, -1 otherwise -CHARACTER(len=1), INTENT(IN) :: HFRAC_ICE -CHARACTER(len=4), INTENT(IN) :: HCONDENS -CHARACTER(len=*), INTENT(IN) :: HLAMBDA3 ! formulation for lambda3 coeff -REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PPABS ! pressure (Pa) -REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PZZ ! height of model levels (m) -REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PRHODREF -REAL, DIMENSION(KIU,KJU,KKU), INTENT(INOUT) :: PT ! grid scale T (K) -REAL, DIMENSION(KIU,KJU,KKU), INTENT(INOUT) :: PRV ! grid scale water vapor mixing ratio (kg/kg) -REAL, DIMENSION(KIU,KJU,KKU), INTENT(INOUT) :: PRC ! grid scale r_c mixing ratio (kg/kg) -REAL, DIMENSION(KIU,KJU,KKU), INTENT(INOUT) :: PRI ! grid scale r_i (kg/kg) -REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PRS ! grid scale mixing ration of snow (kg/kg) -REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PRG ! grid scale mixing ration of graupel (kg/kg) -REAL, DIMENSION(KIU,KJU,KKU), INTENT(IN) :: PSIGS ! Sigma_s from turbulence scheme -REAL, DIMENSION(:,:,:), INTENT(IN) :: PMFCONV! convective mass flux (kg /s m^2) -REAL, DIMENSION(KIU,KJU,KKU), INTENT(OUT) :: PCLDFR ! cloud fraction -REAL, DIMENSION(KIU,KJU,KKU), INTENT(OUT) :: PSIGRC ! s r_c / sig_s^2 - -REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(IN) :: PLV ! Latent heat L_v -REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(IN) :: PLS ! Latent heat L_s -REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(IN) :: PCPH ! Specific heat C_ph -REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(OUT) :: PHLC_HRC -REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(OUT) :: PHLC_HCF ! cloud fraction -REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(OUT) :: PHLI_HRI -REAL, DIMENSION(KIU,KJU,KKU), OPTIONAL, INTENT(OUT) :: PHLI_HCF -LOGICAL, INTENT(IN) :: OUSERI ! logical switch to compute both - ! liquid and solid condensate (OUSERI=.TRUE.) - ! or only solid condensate (OUSERI=.FALSE.) -LOGICAL, INTENT(IN) :: OSIGMAS! use present global Sigma_s values - ! or that from turbulence scheme -REAL, INTENT(IN) :: PSIGQSAT ! use an extra "qsat" variance contribution (OSIGMAS case) - -! -! -!* 0.2 Declarations of local variables : -! -INTEGER :: JI, JJ, JK, JKP, JKM, IKTB, IKTE ! loop index -REAL, DIMENSION(KIU,KJU,KKU) :: ZTLK, ZRT ! work arrays for T_l and total water mixing ratio -REAL, DIMENSION(KIU,KJU,KKU) :: ZL ! length scale -REAL, DIMENSION(KIU,KJU,KKU) :: ZFRAC ! Ice fraction -REAL, DIMENSION(KIU,KJU,KKU) :: ZCRIAUTI ! -INTEGER, DIMENSION(KIU,KJU) :: ITPL ! top levels of troposphere -REAL, DIMENSION(KIU,KJU) :: ZTMIN ! minimum Temp. related to ITPL -! -REAL, DIMENSION(KIU,KJU,KKU) :: ZLV, ZLS, ZCPD -REAL, DIMENSION(KIU,KJU,KKU) :: ZCOND -REAL :: ZGCOND, ZSBAR, ZSBARC, ZQ1, ZAUTC, ZAUTI, ZGAUV, ZGAUC, ZGAUI, ZGAUTC, ZGAUTI ! Used for integration in Gaussian Probability Density Function -REAL :: ZTEMP, ZPV, ZQSL, ZPIV, ZQSI, ZLVS ! thermodynamics -REAL :: ZLL, DZZ, ZZZ ! used for length scales -REAL :: ZAH, ZA, ZB, ZSIGMA, ZDRW, ZDTL, ZSIG_CONV ! related to computation of Sig_s -REAL :: ZRCOLD, ZRIOLD -INTEGER :: INQ1 -REAL :: ZINC -! -!* 0.3 Definition of constants : -! -!------------------------------------------------------------------------------- -! -REAL,PARAMETER :: ZL0 = 600. ! tropospheric length scale -REAL,PARAMETER :: ZCSIGMA = 0.2 ! constant in sigma_s parameterization -REAL,PARAMETER :: ZCSIG_CONV = 0.30E-2 ! scaling factor for ZSIG_CONV as function of mass flux -! - -REAL, DIMENSION(-22:11),PARAMETER :: ZSRC_1D =(/ & - 0. , 0. , 2.0094444E-04, 0.316670E-03, & - 4.9965648E-04, 0.785956E-03 , 1.2341294E-03, 0.193327E-02, & - 3.0190963E-03, 0.470144E-02 , 7.2950651E-03, 0.112759E-01, & - 1.7350994E-02, 0.265640E-01 , 4.0427860E-02, 0.610997E-01, & - 9.1578111E-02, 0.135888E+00 , 0.1991484 , 0.230756E+00, & - 0.2850565 , 0.375050E+00 , 0.5000000 , 0.691489E+00, & - 0.8413813 , 0.933222E+00 , 0.9772662 , 0.993797E+00, & - 0.9986521 , 0.999768E+00 , 0.9999684 , 0.999997E+00, & - 1.0000000 , 1.000000 /) -! -!------------------------------------------------------------------------------- -! -! - -IKTB=1+JPVEXT -IKTE=KKU-JPVEXT - -PCLDFR(:,:,:) = 0. ! Initialize values -PSIGRC(:,:,:) = 0. ! Initialize values -! -! -!------------------------------------------------------------------------------- -! store total water mixing ratio -DO JK=IKTB,IKTE - DO JJ=KJB,KJE - DO JI=KIB,KIE - ZRT(JI,JJ,JK) = PRV(JI,JJ,JK) + PRC(JI,JJ,JK) + PRI(JI,JJ,JK) - END DO - END DO -END DO -!------------------------------------------------------------------------------- -! Preliminary calculations -! latent heat of vaporisation/sublimation -IF(PRESENT(PLV) .AND. PRESENT(PLS)) THEN - ZLV(:,:,:)=PLV(:,:,:) - ZLS(:,:,:)=PLS(:,:,:) -ELSE - DO JK=IKTB,IKTE - DO JJ=KJB,KJE - DO JI=KIB,KIE - ZTEMP = PT(JI,JJ,JK) - ! latent heat of vaporisation/sublimation - ZLV(JI,JJ,JK) = XLVTT + ( XCPV - XCL ) * ( ZTEMP - XTT ) - ZLS(JI,JJ,JK) = XLSTT + ( XCPV - XCI ) * ( ZTEMP - XTT ) - ENDDO - ENDDO - ENDDO -ENDIF -IF(PRESENT(PCPH)) THEN - ZCPD(:,:,:)=PCPH(:,:,:) -ELSE - DO JK=IKTB,IKTE - DO JJ=KJB,KJE - DO JI=KIB,KIE - ZCPD(JI,JJ,JK) = XCPD + XCPV*PRV(JI,JJ,JK) + XCL*PRC(JI,JJ,JK) + XCI*PRI(JI,JJ,JK) + & - XCI*(PRS(JI,JJ,JK) + PRG(JI,JJ,JK) ) - ENDDO - ENDDO - ENDDO -ENDIF -!------------------------------------------------------------------------------- -! Preliminary calculations needed for computing the "turbulent part" of Sigma_s -IF ( .NOT. OSIGMAS ) THEN - DO JK=IKTB,IKTE - DO JJ=KJB,KJE - DO JI=KIB,KIE - ZTEMP = PT(JI,JJ,JK) - ! store temperature at saturation - ZTLK(JI,JJ,JK) = ZTEMP - ZLV(JI,JJ,JK)*PRC(JI,JJ,JK)/ZCPD(JI,JJ,JK) & - - ZLS(JI,JJ,JK)*PRI(JI,JJ,JK)/ZCPD(JI,JJ,JK) - END DO - END DO - END DO - ! Determine tropopause/inversion height from minimum temperature - ITPL(:,:) = KIB+1 - ZTMIN(:,:) = 400. - DO JK = IKTB+1,IKTE-1 - DO JJ=KJB,KJE - DO JI=KIB,KIE - IF ( PT(JI,JJ,JK) < ZTMIN(JI,JJ) ) THEN - ZTMIN(JI,JJ) = PT(JI,JJ,JK) - ITPL(JI,JJ) = JK - ENDIF - END DO - END DO - END DO - ! Set the mixing length scale - ZL(:,:,KKB) = 20. - DO JK = KKB+KKL,KKE,KKL - DO JJ=KJB,KJE - DO JI=KIB,KIE - ! free troposphere - ZL(JI,JJ,JK) = ZL0 - ZZZ = PZZ(JI,JJ,JK) - PZZ(JI,JJ,KKB) - JKP = ITPL(JI,JJ) - ! approximate length for boundary-layer - IF ( ZL0 > ZZZ ) ZL(JI,JJ,JK) = ZZZ - ! gradual decrease of length-scale near and above tropopause - IF ( ZZZ > 0.9*(PZZ(JI,JJ,JKP)-PZZ(JI,JJ,KKB)) ) & - ZL(JI,JJ,JK) = .6 * ZL(JI,JJ,JK-KKL) - END DO - END DO - END DO -END IF -!------------------------------------------------------------------------------- -! -! -!Ice fraction -ZFRAC(:,:,:) = 0. -IF (OUSERI) THEN - WHERE(PRC(:,:,:)+PRI(:,:,:) > 1.E-20) - ZFRAC(:,:,:) = PRI(:,:,:) / (PRC(:,:,:)+PRI(:,:,:)) - ENDWHERE - CALL COMPUTE_FRAC_ICE(SIZE(ZFRAC, 1), SIZE(ZFRAC, 2), SIZE(ZFRAC, 3), HFRAC_ICE, ZFRAC, PT) -ENDIF -! -DO JK=IKTB,IKTE - JKP=MAX(MIN(JK+KKL,IKTE),IKTB) - JKM=MAX(MIN(JK-KKL,IKTE),IKTB) - DO JJ=KJB,KJE - DO JI=KIB,KIE - ! latent heats - ZTEMP = PT(JI,JJ,JK) - ! saturated water vapor mixing ratio over liquid water - ZPV = MIN(EXP( XALPW - XBETAW / ZTEMP - XGAMW * LOG( ZTEMP ) ), .99*PPABS(JI,JJ,JK)) - ZQSL = XRD / XRV * ZPV / ( PPABS(JI,JJ,JK) - ZPV ) - - ! saturated water vapor mixing ratio over ice - ZPIV = MIN(EXP( XALPI - XBETAI / ZTEMP - XGAMI * LOG( ZTEMP ) ), .99*PPABS(JI,JJ,JK)) - ZQSI = XRD / XRV * ZPIV / ( PPABS(JI,JJ,JK) - ZPIV ) - - ! interpolate between liquid and solid as function of temperature - ZQSL = (1. - ZFRAC(JI,JJ,JK)) * ZQSL + ZFRAC(JI,JJ,JK) * ZQSI - ZLVS = (1. - ZFRAC(JI,JJ,JK)) * ZLV(JI,JJ,JK) + & - & ZFRAC(JI,JJ,JK) * ZLS(JI,JJ,JK) - - ! coefficients a and b - ZAH = ZLVS * ZQSL / ( XRV * ZTEMP**2 ) * (XRV * ZQSL / XRD + 1.) - ZA = 1. / ( 1. + ZLVS/ZCPD(JI,JJ,JK) * ZAH ) - ZB = ZAH * ZA - - ZSBAR = ZA * ( ZRT(JI,JJ,JK) - ZQSL + & - & ZAH * ZLVS * (PRC(JI,JJ,JK)+PRI(JI,JJ,JK)) / ZCPD(JI,JJ,JK)) - - ! switch to take either present computed value of SIGMAS - ! or that of Meso-NH turbulence scheme - IF ( OSIGMAS ) THEN - IF (PSIGQSAT/=0.) THEN - ZSIGMA = SQRT((2*PSIGS(JI,JJ,JK))**2 + (PSIGQSAT*ZQSL*ZA)**2) - ELSE - ZSIGMA = 2*PSIGS(JI,JJ,JK) - END IF - ELSE - ! parameterize Sigma_s with first_order closure - DZZ = PZZ(JI,JJ,JKP) - PZZ(JI,JJ,JKM) - ZDRW = ZRT(JI,JJ,JKP) - ZRT(JI,JJ,JKM) - ZDTL = ZTLK(JI,JJ,JKP) - ZTLK(JI,JJ,JKM) + XG/ZCPD(JI,JJ,JK) * DZZ - ZLL = ZL(JI,JJ,JK) - ! standard deviation due to convection - ZSIG_CONV =0. - IF( SIZE(PMFCONV) /= 0) & - ZSIG_CONV = ZCSIG_CONV * PMFCONV(JI,JJ,JK) / ZA - ! zsigma should be of order 4.e-4 in lowest 5 km of atmosphere - ZSIGMA = SQRT( MAX( 1.E-25, ZCSIGMA * ZCSIGMA * ZLL*ZLL/(DZZ*DZZ)*(& - ZA*ZA*ZDRW*ZDRW - 2.*ZA*ZB*ZDRW*ZDTL + ZB*ZB*ZDTL*ZDTL) + & - ZSIG_CONV * ZSIG_CONV ) ) - END IF - ZSIGMA= MAX( 1.E-10, ZSIGMA ) -! ZSIGMA= MAX( 1.E-12, ZSIGMA ) - - ! normalized saturation deficit - ZQ1 = ZSBAR/ZSIGMA - - IF(HCONDENS == 'GAUS')THEN - ! Gaussian Probability Density Function around ZQ1 - ! Computation of ZG and ZGAM(=erf(ZG)) - ZGCOND = -ZQ1/SQRT(2.) - - !Approximation of erf function for Gaussian distribution - ZGAUV = 1 - SIGN(1., ZGCOND) * SQRT(1-EXP(-4*ZGCOND**2/XPI)) - - !Computation Cloud Fraction - PCLDFR(JI,JJ,JK) = MAX( 0., MIN(1.,0.5*ZGAUV)) - - !Computation of condensate - ZCOND(JI,JJ,JK) = (EXP(-ZGCOND**2)-ZGCOND*SQRT(XPI)*ZGAUV)*ZSIGMA/SQRT(2.*XPI) - ZCOND(JI,JJ,JK) = MAX(ZCOND(JI,JJ,JK), 0.) - - PSIGRC(JI,JJ,JK) = PCLDFR(JI,JJ,JK) - - !Computation warm/cold Cloud Fraction and content in high water content part - IF(PRESENT(PHLC_HCF) .AND. PRESENT(PHLC_HRC))THEN - IF(1-ZFRAC(JI,JJ,JK) > 1.E-20)THEN - ZAUTC = (ZSBAR - XCRIAUTC/(PRHODREF(JI,JJ,JK)*(1-ZFRAC(JI,JJ,JK))))/ZSIGMA - ZGAUTC = -ZAUTC/SQRT(2.) - !Approximation of erf function for Gaussian distribution - ZGAUC = 1 - SIGN(1., ZGAUTC) * SQRT(1-EXP(-4*ZGAUTC**2/XPI)) - PHLC_HCF(JI,JJ,JK) = MAX( 0., MIN(1.,0.5*ZGAUC)) - PHLC_HRC(JI,JJ,JK) = (1-ZFRAC(JI,JJ,JK))*(EXP(-ZGAUTC**2)-ZGAUTC*SQRT(XPI)*ZGAUC)*ZSIGMA/SQRT(2.*XPI) - PHLC_HRC(JI,JJ,JK) = PHLC_HRC(JI,JJ,JK) + XCRIAUTC/PRHODREF(JI,JJ,JK) * PHLC_HCF(JI,JJ,JK) - PHLC_HRC(JI,JJ,JK) = MAX(PHLC_HRC(JI,JJ,JK), 0.) - ELSE - PHLC_HCF(JI,JJ,JK)=0. - PHLC_HRC(JI,JJ,JK)=0. - ENDIF - ENDIF - - IF(PRESENT(PHLI_HCF) .AND. PRESENT(PHLI_HRI))THEN - IF(ZFRAC(JI,JJ,JK) > 1.E-20)THEN - ZCRIAUTI(JI,JJ,JK)=MIN(XCRIAUTI,10**(XACRIAUTI*(PT(JI,JJ,JK)-XTT)+XBCRIAUTI)) - ZAUTI = (ZSBAR - ZCRIAUTI(JI,JJ,JK)/ZFRAC(JI,JJ,JK))/ZSIGMA - ZGAUTI = -ZAUTI/SQRT(2.) - !Approximation of erf function for Gaussian distribution - ZGAUI = 1 - SIGN(1., ZGAUTI) * SQRT(1-EXP(-4*ZGAUTI**2/XPI)) - PHLI_HCF(JI,JJ,JK) = MAX( 0., MIN(1.,0.5*ZGAUI)) - PHLI_HRI(JI,JJ,JK) = ZFRAC(JI,JJ,JK)*(EXP(-ZGAUTI**2)-ZGAUTI*SQRT(XPI)*ZGAUI)*ZSIGMA/SQRT(2.*XPI) - PHLI_HRI(JI,JJ,JK) = PHLI_HRI(JI,JJ,JK) + ZCRIAUTI(JI,JJ,JK)*PHLI_HCF(JI,JJ,JK) - PHLI_HRI(JI,JJ,JK) = MAX(PHLI_HRI(JI,JJ,JK), 0.) - ELSE - PHLI_HCF(JI,JJ,JK)=0. - PHLI_HRI(JI,JJ,JK)=0. - ENDIF - ENDIF - - ELSEIF(HCONDENS == 'CB02')THEN - !Cloud fraction - PCLDFR(JI,JJ,JK) = MAX( 0., MIN(1.,0.5+0.36*ATAN(1.55*ZQ1)) ) - - !Total condensate - IF (ZQ1 > 0. .AND. ZQ1 <= 2) THEN - ZCOND(JI,JJ,JK) = MIN(EXP(-1.)+.66*ZQ1+.086*ZQ1**2, 2.) ! We use the MIN function for continuity - ELSE IF (ZQ1 > 2.) THEN - ZCOND(JI,JJ,JK) = ZQ1 - ELSE - ZCOND(JI,JJ,JK) = EXP( 1.2*ZQ1-1. ) - ENDIF - ZCOND(JI,JJ,JK) = ZCOND(JI,JJ,JK) * ZSIGMA - - INQ1 = MIN( MAX(-22,FLOOR(MIN(100., MAX(-100., 2*ZQ1))) ), 10) !inner min/max prevents sigfpe when 2*zq1 does not fit into an int - ZINC = 2.*ZQ1 - INQ1 - - PSIGRC(JI,JJ,JK) = MIN(1.,(1.-ZINC)*ZSRC_1D(INQ1)+ZINC*ZSRC_1D(INQ1+1)) - - IF(PRESENT(PHLC_HCF) .AND. PRESENT(PHLC_HRC))THEN - PHLC_HCF(JI,JJ,JK)=0. - PHLC_HRC(JI,JJ,JK)=0. - ENDIF - IF(PRESENT(PHLI_HCF) .AND. PRESENT(PHLI_HRI))THEN - PHLI_HCF(JI,JJ,JK)=0. - PHLI_HRI(JI,JJ,JK)=0. - ENDIF - ENDIF - - IF ( ZCOND(JI,JJ,JK) < 1.E-12 ) THEN - ZCOND(JI,JJ,JK) = 0. - PCLDFR(JI,JJ,JK) = 0. - ENDIF - IF (PCLDFR(JI,JJ,JK)==0.) THEN - ZCOND(JI,JJ,JK)=0. - ENDIF - - ZRCOLD=PRC(JI,JJ,JK) - ZRIOLD=PRI(JI,JJ,JK) - - PRC(JI,JJ,JK) = (1.-ZFRAC(JI,JJ,JK)) * ZCOND(JI,JJ,JK) ! liquid condensate - PRI(JI,JJ,JK) = ZFRAC(JI,JJ,JK) * ZCOND(JI,JJ,JK) ! solid condensate - - PT(JI,JJ,JK) = PT(JI,JJ,JK) + ((PRC(JI,JJ,JK)-ZRCOLD)*ZLV(JI,JJ,JK) + & - &(PRI(JI,JJ,JK)-ZRIOLD)*ZLS(JI,JJ,JK) ) & - & /ZCPD(JI,JJ,JK) - PRV(JI,JJ,JK) = ZRT(JI,JJ,JK) - PRC(JI,JJ,JK) - PRI(JI,JJ,JK) - -! s r_c/ sig_s^2 -! PSIGRC(JI,JJ,JK) = PCLDFR(JI,JJ,JK) ! use simple Gaussian relation -! -! multiply PSRCS by the lambda3 coefficient -! -! PSIGRC(JI,JJ,JK) = 2.*PCLDFR(JI,JJ,JK) * MIN( 3. , MAX(1.,1.-ZQ1) ) -! in the 3D case lambda_3 = 1. - - IF(HLAMBDA3=='CB')THEN - PSIGRC(JI,JJ,JK) = PSIGRC(JI,JJ,JK)* MIN( 3. , MAX(1.,1.-ZQ1) ) - ELSEIF(HLAMBDA3=='NONE') THEN - ELSE - call Print_msg( NVERB_FATAL, 'GEN', 'CONDENSATION', 'invalid value for HLAMBDA3: ' // TRIM( HLAMBDA3 ) ) - ENDIF - - END DO - END DO -END DO -! -END SUBROUTINE CONDENSATION diff --git a/src/mesonh/micro/ice_adjust.f90 b/src/mesonh/micro/ice_adjust.f90 deleted file mode 100644 index 8f5a8b35e..000000000 --- a/src/mesonh/micro/ice_adjust.f90 +++ /dev/null @@ -1,524 +0,0 @@ -!MNH_LIC Copyright 1996-2021 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 for details. version 1. -!----------------------------------------------------------------- -! ###################### - MODULE MODI_ICE_ADJUST -! ###################### -! -INTERFACE -! - SUBROUTINE ICE_ADJUST (KKA, KKU, KKL, KRR, HFRAC_ICE, HCONDENS, HLAMBDA3,& - HBUNAME, OSUBG_COND, OSIGMAS, HSUBG_MF_PDF, & - PTSTEP, PSIGQSAT, & - PRHODJ, PEXNREF, PRHODREF, PSIGS, PMFCONV, & - PPABST, PZZ, & - PEXN, PCF_MF, PRC_MF, PRI_MF, & - PRV, PRC, PRVS, PRCS, PTH, PTHS, PSRCS, PCLDFR, & - PRR, PRI, PRIS, PRS, PRG, PRH, & - POUT_RV, POUT_RC, POUT_RI, POUT_TH, & - PHLC_HRC, PHLC_HCF, PHLI_HRI, PHLI_HCF ) -! -INTEGER, INTENT(IN) :: KKA !near ground array index -INTEGER, INTENT(IN) :: KKU !uppest atmosphere array index -INTEGER, INTENT(IN) :: KKL !vert. levels type 1=MNH -1=ARO -INTEGER, INTENT(IN) :: KRR ! Number of moist variables -CHARACTER(len=1), INTENT(IN) :: HFRAC_ICE -CHARACTER(len=80), INTENT(IN) :: HCONDENS -CHARACTER(len=4), INTENT(IN) :: HLAMBDA3 ! formulation for lambda3 coeff -CHARACTER(len=*), INTENT(IN) :: HBUNAME ! Name of the budget -LOGICAL, INTENT(IN) :: OSUBG_COND ! Switch for Subgrid - ! Condensation -LOGICAL :: OSIGMAS ! Switch for Sigma_s: - ! use values computed in CONDENSATION - ! or that from turbulence scheme -CHARACTER(len=*), INTENT(IN) :: HSUBG_MF_PDF -REAL, INTENT(IN) :: PTSTEP ! Double Time step - ! (single if cold start) -REAL, INTENT(IN) :: PSIGQSAT ! coeff applied to qsat variance contribution -! -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRHODJ ! Dry density * Jacobian -REAL, DIMENSION(:,:,:), INTENT(IN) :: PEXNREF ! Reference Exner function -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRHODREF -! -REAL, DIMENSION(:,:,:), INTENT(IN) :: PSIGS ! Sigma_s at time t -REAL, DIMENSION(:,:,:), INTENT(IN) :: PMFCONV ! convective mass flux -REAL, DIMENSION(:,:,:), INTENT(IN) :: PPABST ! Absolute Pressure at t -REAL, DIMENSION(:,:,:), INTENT(IN) :: PZZ ! height of model layer -REAL, DIMENSION(:,:,:), INTENT(IN) :: PEXN ! Exner function -! -REAL, DIMENSION(:,:,:), INTENT(IN) :: PCF_MF! Convective Mass Flux Cloud fraction -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRI_MF! Convective Mass Flux ice mixing ratio -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRC_MF! Convective Mass Flux liquid mixing ratio -! -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRV ! Water vapor m.r. to adjust -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRC ! Cloud water m.r. to adjust -REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PRVS ! Water vapor m.r. source -REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PRCS ! Cloud water m.r. source -REAL, DIMENSION(:,:,:), INTENT(IN) :: PTH ! Theta to adjust -REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PTHS ! Theta source -REAL, DIMENSION(:,:,:), INTENT(OUT) :: PSRCS ! Second-order flux - ! s'rc'/2Sigma_s2 at time t+1 - ! multiplied by Lambda_3 -REAL, DIMENSION(:,:,:), INTENT(OUT) :: PCLDFR ! Cloud fraction -REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PRIS ! Cloud ice m.r. at t+1 -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRR ! Rain water m.r. to adjust -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRI ! Cloud ice m.r. to adjust -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRS ! Aggregate m.r. to adjust -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRG ! Graupel m.r. to adjust -REAL, DIMENSION(:,:,:), OPTIONAL, INTENT(IN) :: PRH ! Hail m.r. to adjust -REAL, DIMENSION(:,:,:), OPTIONAL, INTENT(OUT) :: POUT_RV ! Adjusted value -REAL, DIMENSION(:,:,:), OPTIONAL, INTENT(OUT) :: POUT_RC ! Adjusted value -REAL, DIMENSION(:,:,:), OPTIONAL, INTENT(OUT) :: POUT_RI ! Adjusted value -REAL, DIMENSION(:,:,:), OPTIONAL, INTENT(OUT) :: POUT_TH ! Adjusted value -REAL, DIMENSION(:,:,:), OPTIONAL, INTENT(OUT) :: PHLC_HRC -REAL, DIMENSION(:,:,:), OPTIONAL, INTENT(OUT) :: PHLC_HCF -REAL, DIMENSION(:,:,:), OPTIONAL, INTENT(OUT) :: PHLI_HRI -REAL, DIMENSION(:,:,:), OPTIONAL, INTENT(OUT) :: PHLI_HCF -! -! -END SUBROUTINE ICE_ADJUST -! -END INTERFACE -! -END MODULE MODI_ICE_ADJUST - -! ########################################################################## - SUBROUTINE ICE_ADJUST (KKA, KKU, KKL, KRR, HFRAC_ICE, HCONDENS, HLAMBDA3,& - HBUNAME, OSUBG_COND, OSIGMAS, HSUBG_MF_PDF, & - PTSTEP, PSIGQSAT, & - PRHODJ, PEXNREF, PRHODREF, PSIGS, PMFCONV, & - PPABST, PZZ, & - PEXN, PCF_MF, PRC_MF, PRI_MF, & - PRV, PRC, PRVS, PRCS, PTH, PTHS, PSRCS, PCLDFR, & - PRR, PRI, PRIS, PRS, PRG, PRH, & - POUT_RV, POUT_RC, POUT_RI, POUT_TH, & - PHLC_HRC, PHLC_HCF, PHLI_HRI, PHLI_HCF ) -! ######################################################################### -! -!!**** *ICE_ADJUST* - compute the ajustment of water vapor in mixed-phase -!! clouds -!! -!! PURPOSE -!! ------- -!! The purpose of this routine is to compute the fast microphysical sources -!! through a saturation ajustement procedure in case of mixed-phase clouds. -!! -!! -!!** METHOD -!! ------ -!! Langlois, Tellus, 1973 for the cloudless version. -!! When cloud water is taken into account, refer to book 1 of the -!! documentation. -!! -!! -!! -!! EXTERNAL -!! -------- -!! None -!! -!! -!! IMPLICIT ARGUMENTS -!! ------------------ -!! Module MODD_CST -!! XP00 ! Reference pressure -!! XMD,XMV ! Molar mass of dry air and molar mass of vapor -!! XRD,XRV ! Gaz constant for dry air, gaz constant for vapor -!! XCPD,XCPV ! Cpd (dry air), Cpv (vapor) -!! XCL ! Cl (liquid) -!! XCI ! Ci (ice) -!! XTT ! Triple point temperature -!! XLVTT ! Vaporization heat constant -!! XLSTT ! Sublimation heat constant -!! XALPW,XBETAW,XGAMW ! Constants for saturation vapor over liquid -!! ! pressure function -!! XALPI,XBETAI,XGAMI ! Constants for saturation vapor over ice -!! ! pressure function -!! Module MODD_CONF -!! CCONF -!! Module MODD_BUDGET: -!! NBUMOD -!! CBUTYPE -!! LBU_RTH -!! LBU_RRV -!! LBU_RRC -!! LBU_RRI -!! -!! -!! REFERENCE -!! --------- -!! Book 1 and Book2 of documentation ( routine ICE_ADJUST ) -!! Langlois, Tellus, 1973 -!! -!! AUTHOR -!! ------ -!! J.-P. Pinty * Laboratoire d'Aerologie* -!! -!! -!! MODIFICATIONS -!! ------------- -!! Original 06/12/96 -!! M. Tomasini 27/11/00 Change CND and DEP fct of the T instead of rc and ri -!! Avoid the sub- and super-saturation before the ajustment -!! Avoid rc>0 if T<T00 before the ajustment -!! P Bechtold 12/02/02 change subgrid condensation -!! JP Pinty 29/11/02 add ICE2 and IC4 cases -!! (P. Jabouille) 27/05/04 safety test for case where esw/i(T)> pabs (~Z>40km) -!! J.Pergaud and S.Malardel Add EDKF case -!! S. Riette ice for EDKF -!! 2012-02 Y. Seity, add possibility to run with reversed vertical levels -!! J.Escobar : 15/09/2015 : WENO5 & JPHEXT <> 1 -!! 2016-07 S. Riette: adjustement is now realized on state variables (PRV, PRC, PRI, PTH) -!! whereas tendencies are still applied on S variables. -!! This modification allows to call ice_adjust on T variable -!! or to call it on S variables -!! 2016-11 S. Riette: all-or-nothing adjustment now uses condensation -! P. Wautelet 05/2016-04/2018: new data structures and calls for I/O -! P. Wautelet 02/2020: use the new data structures and subroutines for budgets -!------------------------------------------------------------------------------- -! -!* 0. DECLARATIONS -! ------------ -! -use modd_budget, only: lbudget_th, lbudget_rv, lbudget_rc, lbudget_ri, & - NBUDGET_TH, NBUDGET_RV, NBUDGET_RC, NBUDGET_RI, & - tbudgets -USE MODD_CONF -USE MODD_CST -USE MODD_PARAMETERS - -use mode_budget, only: Budget_store_init, Budget_store_end -USE MODD_RAIN_ICE_PARAM, ONLY : XCRIAUTC, XCRIAUTI, XACRIAUTI, XBCRIAUTI -use mode_tools_ll, only: GET_INDICE_ll - -USE MODI_CONDENSATION -USE MODI_GET_HALO - -IMPLICIT NONE -! -! -!* 0.1 Declarations of dummy arguments : -! -! -INTEGER, INTENT(IN) :: KKA !near ground array index -INTEGER, INTENT(IN) :: KKU !uppest atmosphere array index -INTEGER, INTENT(IN) :: KKL !vert. levels type 1=MNH -1=ARO -INTEGER, INTENT(IN) :: KRR ! Number of moist variables -CHARACTER(len=1), INTENT(IN) :: HFRAC_ICE -CHARACTER(len=80), INTENT(IN) :: HCONDENS -CHARACTER(len=4), INTENT(IN) :: HLAMBDA3 ! formulation for lambda3 coeff -CHARACTER(len=*), INTENT(IN) :: HBUNAME ! Name of the budget -LOGICAL, INTENT(IN) :: OSUBG_COND ! Switch for Subgrid - ! Condensation -LOGICAL :: OSIGMAS ! Switch for Sigma_s: - ! use values computed in CONDENSATION - ! or that from turbulence scheme -CHARACTER(len=*), INTENT(IN) :: HSUBG_MF_PDF -REAL, INTENT(IN) :: PTSTEP ! Double Time step - ! (single if cold start) -REAL, INTENT(IN) :: PSIGQSAT ! coeff applied to qsat variance contribution -! -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRHODJ ! Dry density * Jacobian -REAL, DIMENSION(:,:,:), INTENT(IN) :: PEXNREF ! Reference Exner function -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRHODREF -! -REAL, DIMENSION(:,:,:), INTENT(IN) :: PSIGS ! Sigma_s at time t -REAL, DIMENSION(:,:,:), INTENT(IN) :: PMFCONV ! convective mass flux -REAL, DIMENSION(:,:,:), INTENT(IN) :: PPABST ! Absolute Pressure at t -REAL, DIMENSION(:,:,:), INTENT(IN) :: PZZ ! height of model layer -REAL, DIMENSION(:,:,:), INTENT(IN) :: PEXN ! Exner function -! -REAL, DIMENSION(:,:,:), INTENT(IN) :: PCF_MF! Convective Mass Flux Cloud fraction -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRC_MF! Convective Mass Flux liquid mixing ratio -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRI_MF! Convective Mass Flux ice mixing ratio -! -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRV ! Water vapor m.r. to adjust -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRC ! Cloud water m.r. to adjust -REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PRVS ! Water vapor m.r. source -REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PRCS ! Cloud water m.r. source -REAL, DIMENSION(:,:,:), INTENT(IN) :: PTH ! Theta to adjust -REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PTHS ! Theta source -REAL, DIMENSION(:,:,:), INTENT(OUT) :: PSRCS ! Second-order flux - ! s'rc'/2Sigma_s2 at time t+1 - ! multiplied by Lambda_3 -REAL, DIMENSION(:,:,:), INTENT(OUT) :: PCLDFR ! Cloud fraction -! -REAL, DIMENSION(:,:,:), INTENT(INOUT):: PRIS ! Cloud ice m.r. at t+1 -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRR ! Rain water m.r. to adjust -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRI ! Cloud ice m.r. to adjust -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRS ! Aggregate m.r. to adjust -REAL, DIMENSION(:,:,:), INTENT(IN) :: PRG ! Graupel m.r. to adjust -REAL, DIMENSION(:,:,:), OPTIONAL, INTENT(IN) :: PRH ! Hail m.r. to adjust -REAL, DIMENSION(:,:,:), OPTIONAL, INTENT(OUT) :: POUT_RV ! Adjusted value -REAL, DIMENSION(:,:,:), OPTIONAL, INTENT(OUT) :: POUT_RC ! Adjusted value -REAL, DIMENSION(:,:,:), OPTIONAL, INTENT(OUT) :: POUT_RI ! Adjusted value -REAL, DIMENSION(:,:,:), OPTIONAL, INTENT(OUT) :: POUT_TH ! Adjusted value -REAL, DIMENSION(:,:,:), OPTIONAL, INTENT(OUT) :: PHLC_HRC -REAL, DIMENSION(:,:,:), OPTIONAL, INTENT(OUT) :: PHLC_HCF -REAL, DIMENSION(:,:,:), OPTIONAL, INTENT(OUT) :: PHLI_HRI -REAL, DIMENSION(:,:,:), OPTIONAL, INTENT(OUT) :: PHLI_HCF -! -!* 0.2 Declarations of local variables : -! -! -REAL, DIMENSION(SIZE(PEXNREF,1),SIZE(PEXNREF,2),SIZE(PEXNREF,3)) & - :: ZT, & ! adjusted temperature - ZRV, ZRC, ZRI, & ! adjusted state - ZCPH, & ! guess of the CPh for the mixing - ZLV, & ! guess of the Lv at t+1 - ZLS, & ! guess of the Ls at t+1 - ZW1,ZW2, & ! Work arrays for intermediate fields - ZCRIAUT, & ! Autoconversion thresholds - ZHCF, ZHR -! -INTEGER :: IIU,IJU,IKU! dimensions of dummy arrays -INTEGER :: IIB,IJB ! Horz index values of the first inner mass points -INTEGER :: IIE,IJE ! Horz index values of the last inner mass points -INTEGER :: IKB ! K index value of the first inner mass point -INTEGER :: IKE ! K index value of the last inner mass point -INTEGER :: JITER,ITERMAX ! iterative loop for first order adjustment -! -REAL, DIMENSION(SIZE(PEXNREF,1),SIZE(PEXNREF,2),SIZE(PEXNREF,3)) :: ZSIGS,ZSRCS -! -!------------------------------------------------------------------------------- -! -!* 1. PRELIMINARIES -! ------------- -! -if ( lbudget_th ) call Budget_store_init( tbudgets(NBUDGET_TH), trim( hbuname ), pths(:, :, :) * prhodj(:, :, :) ) -if ( lbudget_rv ) call Budget_store_init( tbudgets(NBUDGET_RV), trim( hbuname ), prvs(:, :, :) * prhodj(:, :, :) ) -if ( lbudget_rc ) call Budget_store_init( tbudgets(NBUDGET_RC), trim( hbuname ), prcs(:, :, :) * prhodj(:, :, :) ) -if ( lbudget_ri ) call Budget_store_init( tbudgets(NBUDGET_RI), trim( hbuname ), pris(:, :, :) * prhodj(:, :, :) ) - -IIU = SIZE(PEXNREF,1) -IJU = SIZE(PEXNREF,2) -IKU = SIZE(PEXNREF,3) -CALL GET_INDICE_ll (IIB,IJB,IIE,IJE) -IKB=KKA+JPVEXT*KKL -IKE=KKU-JPVEXT*KKL -! -ITERMAX=1 -! -!------------------------------------------------------------------------------- -! -!* 2. COMPUTE QUANTITIES WITH THE GUESS OF THE FUTURE INSTANT -! ------------------------------------------------------- -! -! -! beginning of the iterative loop (to compute the adjusted state) -ZRV(:,:,:)=PRV(:,:,:) -ZRC(:,:,:)=PRC(:,:,:) -ZRI(:,:,:)=PRI(:,:,:) -ZT(:,:,:)=PTH(:,:,:) * PEXN(:,:,:) -! -DO JITER =1,ITERMAX - ! - !* 2.3 compute the latent heat of vaporization Lv(T*) at t+1 - ! and the latent heat of sublimation Ls(T*) at t+1 - ! - ZLV(:,:,:) = XLVTT + ( XCPV - XCL ) * ( ZT(:,:,:) -XTT ) - ZLS(:,:,:) = XLSTT + ( XCPV - XCI ) * ( ZT(:,:,:) -XTT ) - ! - !* 2.4 compute the specific heat for moist air (Cph) at t+1 - ! - IF ( KRR == 7 ) THEN - ZCPH(:,:,:) = XCPD + XCPV * ZRV(:,:,:) & - + XCL * (ZRC(:,:,:) + PRR(:,:,:)) & - + XCI * (ZRI(:,:,:) + PRS(:,:,:) + PRG(:,:,:) + PRH(:,:,:)) - ELSE IF( KRR == 6 ) THEN - ZCPH(:,:,:) = XCPD + XCPV * ZRV(:,:,:) & - + XCL * (ZRC(:,:,:) + PRR(:,:,:)) & - + XCI * (ZRI(:,:,:) + PRS(:,:,:) + PRG(:,:,:)) - ELSE IF( KRR == 5 ) THEN - ZCPH(:,:,:) = XCPD + XCPV * ZRV(:,:,:) & - + XCL * (ZRC(:,:,:) + PRR(:,:,:)) & - + XCI * (ZRI(:,:,:) + PRS(:,:,:)) - ELSE IF( KRR == 3 ) THEN - ZCPH(:,:,:) = XCPD + XCPV * ZRV(:,:,:) & - + XCL * (ZRC(:,:,:) + PRR(:,:,:)) - ELSE IF( KRR == 2 ) THEN - ZCPH(:,:,:) = XCPD + XCPV * ZRV(:,:,:) & - + XCL * ZRC(:,:,:) - END IF - ! - IF ( OSUBG_COND ) THEN - ! - !* 3. SUBGRID CONDENSATION SCHEME - ! --------------------------- - ! - ! PSRC= s'rci'/Sigma_s^2 - ! ZT, ZRV, ZRC and ZRI are INOUT - CALL CONDENSATION(IIU, IJU, IKU, IIB, IIE, IJB, IJE, IKB, IKE, KKL, & - HFRAC_ICE, HCONDENS, HLAMBDA3, & - PPABST, PZZ, PRHODREF, ZT, ZRV, ZRC, ZRI, PRS, PRG, PSIGS, PMFCONV, PCLDFR, & - PSRCS, .TRUE., OSIGMAS, & - PSIGQSAT, PLV=ZLV, PLS=ZLS, PCPH=ZCPH, PHLC_HRC=PHLC_HRC, PHLC_HCF=PHLC_HCF, PHLI_HRI=PHLI_HRI, PHLI_HCF=PHLI_HCF) - ELSE - ! - !* 4. ALL OR NOTHING CONDENSATION SCHEME - ! FOR MIXED-PHASE CLOUD - ! ----------------------------------------------- - ! - ! - ! ZT, ZRV, ZRC and ZRI are INOUT - ! - !CALL ADJUST_LANGLOIS(IIU, IJU, IKU, IIB, IIE, IJB, IJE, IKB, IKE, KKL, & - ! PPABST, ZT, ZRV, ZRC, ZRI, ZLV, ZLS, ZCPH) HFRAC_ICE must be implemented in Langlois before using it again - ZSIGS=0. - CALL CONDENSATION(IIU, IJU, IKU, IIB, IIE, IJB, IJE, IKB, IKE, KKL, & - HFRAC_ICE, HCONDENS, HLAMBDA3, & - PPABST, PZZ, PRHODREF, ZT, ZRV, ZRC, ZRI, PRS, PRG, ZSIGS, PMFCONV, PCLDFR, & - ZSRCS, .TRUE., OSIGMAS=.TRUE., & - PSIGQSAT=0., PLV=ZLV, PLS=ZLS, PHLC_HRC=PHLC_HRC, PHLC_HCF=PHLC_HCF, PHLI_HRI=PHLI_HRI, PHLI_HCF=PHLI_HCF) - ENDIF -ENDDO ! end of the iterative loop -! -!* 5. COMPUTE THE SOURCES AND STORES THE CLOUD FRACTION -! ------------------------------------------------- -! -! -!* 5.0 compute the variation of mixing ratio -! - ! Rc - Rc* -ZW1(:,:,:) = (ZRC(:,:,:) - PRC(:,:,:)) / PTSTEP ! Pcon = ---------- - ! 2 Delta t - -ZW2(:,:,:) = (ZRI(:,:,:) - PRI(:,:,:)) / PTSTEP ! idem ZW1 but for Ri -! -!* 5.1 compute the sources -! -WHERE( ZW1(:,:,:) < 0.0 ) - ZW1(:,:,:) = MAX ( ZW1(:,:,:), -PRCS(:,:,:) ) -ELSEWHERE - ZW1(:,:,:) = MIN ( ZW1(:,:,:), PRVS(:,:,:) ) -END WHERE -PRVS(:,:,:) = PRVS(:,:,:) - ZW1(:,:,:) -PRCS(:,:,:) = PRCS(:,:,:) + ZW1(:,:,:) -PTHS(:,:,:) = PTHS(:,:,:) + & - ZW1(:,:,:) * ZLV(:,:,:) / (ZCPH(:,:,:) * PEXNREF(:,:,:)) -! -WHERE( ZW2(:,:,:) < 0.0 ) - ZW2(:,:,:) = MAX ( ZW2(:,:,:), -PRIS(:,:,:) ) -ELSEWHERE - ZW2(:,:,:) = MIN ( ZW2(:,:,:), PRVS(:,:,:) ) -END WHERE -PRVS(:,:,:) = PRVS(:,:,:) - ZW2(:,:,:) -PRIS(:,:,:) = PRIS(:,:,:) + ZW2(:,:,:) -PTHS(:,:,:) = PTHS(:,:,:) + & - ZW2(:,:,:) * ZLS(:,:,:) / (ZCPH(:,:,:) * PEXNREF(:,:,:)) -! -! -!* 5.2 compute the cloud fraction PCLDFR -! -IF ( .NOT. OSUBG_COND ) THEN - WHERE (PRCS(:,:,:) + PRIS(:,:,:) > 1.E-12 / PTSTEP) - PCLDFR(:,:,:) = 1. - ELSEWHERE - PCLDFR(:,:,:) = 0. - ENDWHERE - IF ( SIZE(PSRCS,3) /= 0 ) THEN - PSRCS(:,:,:) = PCLDFR(:,:,:) - END IF -ELSE - !We limit PRC_MF+PRI_MF to PRVS*PTSTEP to avoid negative humidity - ZW1(:,:,:)=PRC_MF(:,:,:)/PTSTEP - ZW2(:,:,:)=PRI_MF(:,:,:)/PTSTEP - WHERE(ZW1(:,:,:)+ZW2(:,:,:)>PRVS(:,:,:)) - ZW1(:,:,:)=ZW1(:,:,:)*PRVS(:,:,:)/(ZW1(:,:,:)+ZW2(:,:,:)) - ZW2(:,:,:)=PRVS(:,:,:)-ZW1(:,:,:) - ENDWHERE - IF(PRESENT(PHLC_HRC) .AND. PRESENT(PHLC_HCF)) THEN - ZCRIAUT(:,:,:)=XCRIAUTC/PRHODREF - IF(HSUBG_MF_PDF=='NONE')THEN - WHERE(ZW1(:,:,:)*PTSTEP>PCF_MF * ZCRIAUT) - PHLC_HRC(:,:,:)=PHLC_HRC(:,:,:)+ZW1(:,:,:)*PTSTEP - PHLC_HCF(:,:,:)=MIN(1.,PHLC_HCF(:,:,:)+PCF_MF(:,:,:)) - ENDWHERE - ELSEIF(HSUBG_MF_PDF=='TRIANGLE')THEN - !ZHCF is the precipitating part of the *cloud* and not of the grid cell - WHERE(ZW1(:,:,:)*PTSTEP>PCF_MF*ZCRIAUT(:,:,:)) - ZHCF(:,:,:)=1.-.5*(ZCRIAUT(:,:,:)*PCF_MF(:,:,:) / MAX(1.E-20, ZW1(:,:,:)*PTSTEP))**2 - ZHR(:,:,:)=ZW1(:,:,:)*PTSTEP-(ZCRIAUT(:,:,:)*PCF_MF(:,:,:))**3 / & - &(3*MAX(1.E-20, ZW1(:,:,:)*PTSTEP)**2) - ELSEWHERE(2.*ZW1(:,:,:)*PTSTEP<=PCF_MF * ZCRIAUT(:,:,:)) - ZHCF(:,:,:)=0. - ZHR(:,:,:)=0. - ELSEWHERE - ZHCF(:,:,:)=(2.*ZW1(:,:,:)*PTSTEP-ZCRIAUT(:,:,:)*PCF_MF(:,:,:))**2 / & - &(2.*MAX(1.E-20, ZW1(:,:,:)*PTSTEP)**2) - ZHR(:,:,:)=(4.*(ZW1(:,:,:)*PTSTEP)**3-3.*ZW1(:,:,:)*PTSTEP*(ZCRIAUT(:,:,:)*PCF_MF(:,:,:))**2+& - (ZCRIAUT(:,:,:)*PCF_MF(:,:,:))**3) / & - &(3*MAX(1.E-20, ZW1(:,:,:)*PTSTEP)**2) - ENDWHERE - ZHCF(:,:,:)=ZHCF(:,:,:)*PCF_MF(:,:,:) !to retrieve the part of the grid cell - PHLC_HCF(:,:,:)=MIN(1.,PHLC_HCF(:,:,:)+ZHCF(:,:,:)) !total part of the grid cell that is precipitating - PHLC_HRC(:,:,:)=PHLC_HRC(:,:,:)+ZHR(:,:,:) - ENDIF - ENDIF - IF(PRESENT(PHLI_HRI) .AND. PRESENT(PHLI_HCF)) THEN - ZCRIAUT(:,:,:)=MIN(XCRIAUTI,10**(XACRIAUTI*(ZT(:,:,:)-XTT)+XBCRIAUTI)) - IF(HSUBG_MF_PDF=='NONE')THEN - WHERE(ZW2(:,:,:)*PTSTEP>PCF_MF * ZCRIAUT(:,:,:)) - PHLI_HRI(:,:,:)=PHLI_HRI(:,:,:)+ZW2(:,:,:)*PTSTEP - PHLI_HCF(:,:,:)=MIN(1.,PHLI_HCF(:,:,:)+PCF_MF(:,:,:)) - ENDWHERE - ELSEIF(HSUBG_MF_PDF=='TRIANGLE')THEN - !ZHCF is the precipitating part of the *cloud* and not of the grid cell - WHERE(ZW2(:,:,:)*PTSTEP>PCF_MF*ZCRIAUT) - ZHCF(:,:,:)=1.-.5*(ZCRIAUT*PCF_MF(:,:,:) / (ZW2(:,:,:)*PTSTEP))**2 - ZHR(:,:,:)=ZW2(:,:,:)*PTSTEP-(ZCRIAUT*PCF_MF(:,:,:))**3/(3*(ZW2(:,:,:)*PTSTEP)**2) - ELSEWHERE(2.*ZW2(:,:,:)*PTSTEP<=PCF_MF * ZCRIAUT) - ZHCF(:,:,:)=0. - ZHR(:,:,:)=0. - ELSEWHERE - ZHCF(:,:,:)=(2.*ZW2(:,:,:)*PTSTEP-ZCRIAUT*PCF_MF(:,:,:))**2 / (2.*(ZW2(:,:,:)*PTSTEP)**2) - ZHR(:,:,:)=(4.*(ZW2(:,:,:)*PTSTEP)**3-3.*ZW2(:,:,:)*PTSTEP*(ZCRIAUT*PCF_MF(:,:,:))**2+& - (ZCRIAUT*PCF_MF(:,:,:))**3)/(3*(ZW2(:,:,:)*PTSTEP)**2) - ENDWHERE - ZHCF(:,:,:)=ZHCF(:,:,:)*PCF_MF(:,:,:) !to retrieve the part of the grid cell - PHLI_HCF(:,:,:)=MIN(1.,PHLI_HCF(:,:,:)+ZHCF(:,:,:)) !total part of the grid cell that is precipitating - PHLI_HRI(:,:,:)=PHLI_HRI(:,:,:)+ZHR(:,:,:) - ENDIF - ENDIF - PCLDFR(:,:,:)=MIN(1.,PCLDFR(:,:,:)+PCF_MF(:,:,:)) - PRCS(:,:,:)=PRCS(:,:,:)+ZW1(:,:,:) - PRIS(:,:,:)=PRIS(:,:,:)+ZW2(:,:,:) - PRVS(:,:,:)=PRVS(:,:,:)-(ZW1(:,:,:)+ZW2(:,:,:)) - PTHS(:,:,:) = PTHS(:,:,:) + & - (ZW1 * ZLV(:,:,:) + ZW2 * ZLS(:,:,:)) / ZCPH(:,:,:) & - / PEXNREF(:,:,:) - IF(PRESENT(POUT_RV) .OR. PRESENT(POUT_RC) .OR. & - &PRESENT(POUT_RI) .OR. PRESENT(POUT_TH)) THEN - ZW1(:,:,:)=PRC_MF(:,:,:) - ZW2(:,:,:)=PRI_MF(:,:,:) - WHERE(ZW1(:,:,:)+ZW2(:,:,:)>ZRV(:,:,:)) - ZW1(:,:,:)=ZW1(:,:,:)*ZRV(:,:,:)/(ZW1(:,:,:)+ZW2(:,:,:)) - ZW2(:,:,:)=ZRV(:,:,:)-ZW1(:,:,:) - ENDWHERE - ZRC(:,:,:)=ZRC(:,:,:)+ZW1(:,:,:) - ZRI(:,:,:)=ZRI(:,:,:)+ZW2(:,:,:) - ZRV(:,:,:)=ZRV(:,:,:)-(ZW1(:,:,:)+ZW2(:,:,:)) - ZT(:,:,:) = ZT(:,:,:) + & - (ZW1 * ZLV(:,:,:) + ZW2 * ZLS(:,:,:)) / ZCPH(:,:,:) - ENDIF -ENDIF -! -IF(PRESENT(POUT_RV)) POUT_RV=ZRV -IF(PRESENT(POUT_RC)) POUT_RC=ZRC -IF(PRESENT(POUT_RI)) POUT_RI=ZRI -IF(PRESENT(POUT_TH)) POUT_TH=ZT / PEXN(:,:,:) -! -! -!* 6. STORE THE BUDGET TERMS -! ---------------------- -! -if ( lbudget_th ) call Budget_store_end( tbudgets(NBUDGET_TH), trim( hbuname ), pths(:, :, :) * prhodj(:, :, :) ) -if ( lbudget_rv ) call Budget_store_end( tbudgets(NBUDGET_RV), trim( hbuname ), prvs(:, :, :) * prhodj(:, :, :) ) -if ( lbudget_rc ) call Budget_store_end( tbudgets(NBUDGET_RC), trim( hbuname ), prcs(:, :, :) * prhodj(:, :, :) ) -if ( lbudget_ri ) call Budget_store_end( tbudgets(NBUDGET_RI), trim( hbuname ), pris(:, :, :) * prhodj(:, :, :) ) - -!------------------------------------------------------------------------------ -! -! -END SUBROUTINE ICE_ADJUST diff --git a/src/mesonh/turb/shallow_mf.f90 b/src/mesonh/turb/shallow_mf.f90 index 6bf7075e2..23dbbec39 100644 --- a/src/mesonh/turb/shallow_mf.f90 +++ b/src/mesonh/turb/shallow_mf.f90 @@ -298,7 +298,8 @@ REAL, DIMENSION(SIZE(PTHM,1)) :: ZRESOL_NORM, ZRESOL_GRID,& ! normalized grid ZLUP, ZPLAW ! Test if the ascent continue, if LCL or ETL is reached LOGICAL :: GLMIX - INTEGER :: JI,JJ,JK ! loop counter +INTEGER :: JI,JJ,JK ! loop counter +INTEGER, DIMENSION(SIZE(PTHM,1),SIZE(PTHM,2)) :: IERR !------------------------------------------------------------------------ !!! 1. Initialisation @@ -322,8 +323,7 @@ IF (SIZE(PRM,3).GE.4) THEN ZFRAC_ICE(:,:) = PRM(:,:,4) / (PRM(:,:,2)+PRM(:,:,4)) ENDWHERE ENDIF -CALL COMPUTE_FRAC_ICE(SIZE(ZFRAC_ICE, 1), SIZE(ZFRAC_ICE, 2), & - &HFRAC_ICE,ZFRAC_ICE(:,:),PTHM(:,:)*PEXNM(:,:)) +CALL COMPUTE_FRAC_ICE(HFRAC_ICE,ZFRAC_ICE(:,:),PTHM(:,:)*PEXNM(:,:), IERR(:,:)) ! Conservative variables at t-dt CALL THL_RT_FROM_TH_R_MF(KRR,KRRL,KRRI, & -- GitLab