MODULE lmdz_lscp_poprecip !---------------------------------------------------------------- ! Module for the process-oriented treament of precipitation ! that are called in LSCP ! Authors: Atelier Nuage (G. Riviere, L. Raillard, M. Wimmer, ! N. Dutrievoz, E. Vignon, A. Borella, et al.) ! Jan. 2024 IMPLICIT NONE CONTAINS !---------------------------------------------------------------- ! Computes the processes-oriented precipitation formulations for ! evaporation and sublimation ! SUBROUTINE poprecip_evapsub( & klon, dtime, iftop, paprsdn, paprsup, pplay, temp, tempupnew, qvap, & qprecip, precipfracclr, precipfraccld, & rain, rainclr, raincld, snow, snowclr, snowcld, dqreva, dqssub & ) USE lmdz_lscp_ini, ONLY : prt_level, lunout USE lmdz_lscp_ini, ONLY : coef_eva, coef_eva_i, expo_eva, expo_eva_i, thresh_precip_frac USE lmdz_lscp_ini, ONLY : RCPD, RLSTT, RLVTT, RLMLT, RVTMP2, RTT, RD, RG USE lmdz_lscp_tools, ONLY : calc_qsat_ecmwf IMPLICIT NONE INTEGER, INTENT(IN) :: klon !--number of horizontal grid points [-] REAL, INTENT(IN) :: dtime !--time step [s] LOGICAL, INTENT(IN) :: iftop !--if top of the column REAL, INTENT(IN), DIMENSION(klon) :: paprsdn !--pressure at the bottom interface of the layer [Pa] REAL, INTENT(IN), DIMENSION(klon) :: paprsup !--pressure at the top interface of the layer [Pa] REAL, INTENT(IN), DIMENSION(klon) :: pplay !--pressure in the middle of the layer [Pa] REAL, INTENT(INOUT), DIMENSION(klon) :: temp !--current temperature [K] REAL, INTENT(INOUT), DIMENSION(klon) :: tempupnew !--updated temperature of the overlying layer [K] REAL, INTENT(INOUT), DIMENSION(klon) :: qvap !--current water vapor specific humidity (includes evaporated qi and ql) [kg/kg] REAL, INTENT(INOUT), DIMENSION(klon) :: qprecip !--specific humidity in the precipitation falling from the upper layer [kg/kg] REAL, INTENT(INOUT), DIMENSION(klon) :: precipfracclr !--fraction of precipitation in the clear sky IN THE LAYER ABOVE [-] REAL, INTENT(INOUT), DIMENSION(klon) :: precipfraccld !--fraction of precipitation in the cloudy air IN THE LAYER ABOVE [-] REAL, INTENT(INOUT), DIMENSION(klon) :: rain !--flux of rain gridbox-mean coming from the layer above [kg/s/m2] REAL, INTENT(INOUT), DIMENSION(klon) :: rainclr !--flux of rain gridbox-mean in clear sky coming from the layer above [kg/s/m2] REAL, INTENT(IN), DIMENSION(klon) :: raincld !--flux of rain gridbox-mean in cloudy air coming from the layer above [kg/s/m2] REAL, INTENT(INOUT), DIMENSION(klon) :: snow !--flux of snow gridbox-mean coming from the layer above [kg/s/m2] REAL, INTENT(INOUT), DIMENSION(klon) :: snowclr !--flux of snow gridbox-mean in clear sky coming from the layer above [kg/s/m2] REAL, INTENT(IN), DIMENSION(klon) :: snowcld !--flux of snow gridbox-mean in cloudy air coming from the layer above [kg/s/m2] REAL, INTENT(OUT), DIMENSION(klon) :: dqreva !--rain tendency due to evaporation [kg/kg/s] REAL, INTENT(OUT), DIMENSION(klon) :: dqssub !--snow tendency due to sublimation [kg/kg/s] ! integer for interating over klon INTEGER :: i ! saturation values REAL, DIMENSION(klon) :: qzero, qsat, dqsat, qsatl, dqsatl, qsati, dqsati ! fluxes tendencies because of evaporation REAL :: flevapmax, flevapl, flevapi, flevaptot ! specific humidity tendencies because of evaporation REAL :: dqevapl, dqevapi ! specific heat constant REAL :: cpair, cpw qzero(:) = 0.0 dqreva(:) = 0.0 dqssub(:) = 0.0 dqevapl=0.0 dqevapi=0.0 ! Calculation of saturation specific humidity ! depending on temperature: CALL calc_qsat_ecmwf(klon,temp(:),qzero(:),pplay(:),RTT,0,.false.,qsat(:),dqsat(:)) ! wrt liquid water CALL calc_qsat_ecmwf(klon,temp(:),qzero(:),pplay(:),RTT,1,.false.,qsatl(:),dqsatl(:)) ! wrt ice CALL calc_qsat_ecmwf(klon,temp(:),qzero(:),pplay(:),RTT,2,.false.,qsati(:),dqsati(:)) ! First step consists in "thermalizing" the layer: ! as the flux of precip from layer above "advects" some heat (as the precip is at the temperature ! of the overlying layer) we recalculate a mean temperature that both the air and the precip in the ! layer have. IF (iftop) THEN DO i = 1, klon qprecip(i) = 0. ENDDO ELSE DO i = 1, klon ! no condensed water so cp=cp(vapor+dry air) ! RVTMP2=rcpv/rcpd-1 cpair=RCPD*(1.0+RVTMP2*qvap(i)) cpw=RCPD*RVTMP2 ! qprecip has to be thermalized with ! layer's air so that precipitation at the ground has the ! same temperature as the lowermost layer ! we convert the flux into a specific quantity qprecip qprecip(i) = (rain(i)+snow(i))*dtime/((paprsdn(i)-paprsup(i))/RG) ! t(i,k+1)+d_t(i,k+1): new temperature of the overlying layer temp(i) = ( (tempupnew(i))*qprecip(i)*cpw + cpair*temp(i) ) & / (cpair + qprecip(i)*cpw) ENDDO ENDIF DO i = 1, klon ! if precipitation from the layer above IF ( ( rain(i) + snow(i) ) .GT. 0. ) THEN ! Evaporation of liquid precipitation coming from above ! dP/dz=beta*(1-q/qsat)*(P**expo_eva) (lines 1-2) ! multiplying by dz = - dP / g / rho (line 3-4) ! formula from Sundqvist 1988, Klemp & Wilhemson 1978 ! LTP: evaporation only in the clear sky part flevapl = precipfracclr(i) * coef_eva * (1.0 - qvap(i) / qsatl(i)) & * ( rainclr(i) / MAX(thresh_precip_frac, precipfracclr(i)) ) ** expo_eva & * temp(i) * RD / pplay(i) & * ( paprsdn(i) - paprsup(i) ) / RG ! evaporation is limited by 0 and by the total water amount in ! the precipitation flevapl = MAX(0.0, MIN(flevapl, rainclr(i))) ! sublimation of the solid precipitation coming from above ! (same formula as for liquid precip) flevapi = precipfracclr(i) * coef_eva_i * (1.0 - qvap(i) / qsati(i)) & * ( snowclr(i) / MAX(thresh_precip_frac, precipfracclr(i)) ) ** expo_eva_i & * temp(i) * RD / pplay(i) & * ( paprsdn(i) - paprsup(i) ) / RG ! sublimation is limited by 0 and by the total water amount in ! the precipitation ! TODO: change max when we will allow for vapor deposition in supersaturated regions flevapi = MAX(0.0, MIN(flevapi, snowclr(i))) ! Evaporation limit: we ensure that the layer's fraction below ! the clear sky does not reach saturation. In this case, we ! redistribute the maximum flux flevapmax conserving the ratio liquid/ice ! Max evaporation is computed not to saturate the clear sky precip fraction ! (i.e., the fraction where evaporation occurs) ! It is expressed as a max flux flevapmax ! flevapmax = MAX(0.0, ( qsat(i) - qvap(i) ) * precipfracclr(i)) & * ( paprsdn(i) - paprsup(i) ) / RG / dtime flevaptot = flevapl + flevapi IF ( flevaptot .GT. flevapmax ) THEN flevapl = flevapmax * flevapl / flevaptot flevapi = flevapmax * flevapi / flevaptot ENDIF ! New solid and liquid precipitation fluxes after evap and sublimation dqevapl = flevapl / ( paprsdn(i) - paprsup(i) ) * RG * dtime dqevapi = flevapi / ( paprsdn(i) - paprsup(i) ) * RG * dtime ! vapor is updated after evaporation/sublimation (it is increased) qvap(i) = qvap(i) + dqevapl + dqevapi ! qprecip is the total condensed water in the precip flux (it is decreased) qprecip(i) = qprecip(i) - dqevapl - dqevapi ! air and precip temperature (i.e., gridbox temperature) ! is updated due to latent heat cooling temp(i) = temp(i) & - dqevapl * RLVTT / RCPD & / ( 1.0 + RVTMP2 * ( qvap(i) + qprecip(i) ) ) & - dqevapi * RLSTT / RCPD & / ( 1.0 + RVTMP2 * ( qvap(i) + qprecip(i) ) ) ! New values of liquid and solid precipitation rainclr(i) = rainclr(i) - flevapl snowclr(i) = snowclr(i) - flevapi ! if there is no more precip fluxes, the precipitation fraction in clear ! sky is set to 0 IF ( ( rainclr(i) + snowclr(i) ) .LE. 0. ) precipfracclr(i) = 0. ! calculation of the total fluxes rain(i) = rainclr(i) + raincld(i) snow(i) = snowclr(i) + snowcld(i) ELSE ! if no precip, we reinitialize the cloud fraction used for the precip to 0 precipfraccld(i) = 0. precipfracclr(i) = 0. ENDIF ! ( ( rain(i) + snow(i) ) .GT. 0. ) ! write output tendencies for rain and snow dqssub(i) = -dqevapi/dtime dqreva(i) = -dqevapl/dtime ENDDO ! loop on klon END SUBROUTINE poprecip_evapsub !---------------------------------------------------------------- ! Computes the processes-oriented precipitation formulations for ! - autoconversion (auto) via a deposition process ! - aggregation (agg) ! - riming (rim) ! - collection (coll) ! - melting (melt) ! - freezing (free) ! SUBROUTINE poprecip_postcld( & klon, dtime, paprsdn, paprsup, pplay, ctot_vol, ptconv, & temp, qvap, qliq, qice, icefrac, cldfra, & precipfracclr, precipfraccld, & rain, rainclr, raincld, snow, snowclr, snowcld, & dqrauto,dqrcol,dqrmelt,dqrfreez,dqsauto,dqsagg,dqsrim,dqsmelt,dqsfreez) USE lmdz_lscp_ini, ONLY : prt_level, lunout USE lmdz_lscp_ini, ONLY : RCPD, RLSTT, RLVTT, RLMLT, RVTMP2, RTT, RD, RG USE lmdz_lscp_tools, ONLY : calc_qsat_ecmwf USE lmdz_lscp_ini, ONLY : cld_lc_con, cld_tau_con, cld_expo_con, seuil_neb, & cld_lc_lsc, cld_tau_lsc, cld_expo_lsc, rain_int_min, & thresh_precip_frac, gamma_col, gamma_agg, gamma_rim, & rho_rain, rho_snow, r_rain, r_snow, Eff_rain_liq, & Eff_snow_ice, Eff_snow_liq, tau_auto_snow_min, & tau_auto_snow_max, thresh_precip_frac, eps, & iflag_cloudth_vert, iflag_rain_incloud_vol IMPLICIT NONE INTEGER, INTENT(IN) :: klon !--number of horizontal grid points [-] REAL, INTENT(IN) :: dtime !--time step [s] REAL, INTENT(IN), DIMENSION(klon) :: paprsdn !--pressure at the bottom interface of the layer [Pa] REAL, INTENT(IN), DIMENSION(klon) :: paprsup !--pressure at the top interface of the layer [Pa] REAL, INTENT(IN), DIMENSION(klon) :: pplay !--pressure in the middle of the layer [Pa] REAL, INTENT(IN), DIMENSION(klon) :: ctot_vol !-- LOGICAL, INTENT(IN), DIMENSION(klon) :: ptconv !-- REAL, INTENT(INOUT), DIMENSION(klon) :: temp !--current temperature [K] REAL, INTENT(INOUT), DIMENSION(klon) :: qvap !--current water vapor specific humidity [kg/kg] REAL, INTENT(INOUT), DIMENSION(klon) :: qliq !--current liquid water specific humidity [kg/kg] REAL, INTENT(INOUT), DIMENSION(klon) :: qice !--current ice water specific humidity [kg/kg] REAL, INTENT(IN), DIMENSION(klon) :: icefrac !-- REAL, INTENT(IN), DIMENSION(klon) :: cldfra !-- REAL, INTENT(INOUT), DIMENSION(klon) :: precipfracclr !--fraction of precipitation in the clear sky IN THE LAYER ABOVE [-] REAL, INTENT(INOUT), DIMENSION(klon) :: precipfraccld !--fraction of precipitation in the cloudy air IN THE LAYER ABOVE [-] !--NB. at the end of the routine, becomes the fraction of precip !--in the current layer REAL, INTENT(INOUT), DIMENSION(klon) :: rain !--flux of rain gridbox-mean coming from the layer above [kg/s/m2] REAL, INTENT(INOUT), DIMENSION(klon) :: rainclr !--flux of rain gridbox-mean in clear sky coming from the layer above [kg/s/m2] REAL, INTENT(INOUT), DIMENSION(klon) :: raincld !--flux of rain gridbox-mean in cloudy air coming from the layer above [kg/s/m2] REAL, INTENT(INOUT), DIMENSION(klon) :: snow !--flux of snow gridbox-mean coming from the layer above [kg/s/m2] REAL, INTENT(INOUT), DIMENSION(klon) :: snowclr !--flux of snow gridbox-mean in clear sky coming from the layer above [kg/s/m2] REAL, INTENT(INOUT), DIMENSION(klon) :: snowcld !--flux of snow gridbox-mean in cloudy air coming from the layer above [kg/s/m2] REAL, INTENT(OUT), DIMENSION(klon) :: dqrcol !-- rain tendendy due to collection by rain of liquid cloud droplets [kg/kg/s] REAL, INTENT(OUT), DIMENSION(klon) :: dqsagg !-- snow tendency due to collection of lcoud ice by aggregation [kg/kg/s] REAL, INTENT(OUT), DIMENSION(klon) :: dqrauto !-- rain tendency due to autoconversion of cloud liquid [kg/kg/s] REAL, INTENT(OUT), DIMENSION(klon) :: dqsauto !-- snow tendency due to autoconversion of cloud ice [kg/kg/s] REAL, INTENT(OUT), DIMENSION(klon) :: dqsrim !-- snow tendency due to riming [kg/kg/s] REAL, INTENT(OUT), DIMENSION(klon) :: dqsmelt !-- snow tendency due to melting [kg/kg/s] REAL, INTENT(OUT), DIMENSION(klon) :: dqrmelt !-- rain tendency due to melting [kg/kg/s] REAL, INTENT(OUT), DIMENSION(klon) :: dqsfreez !-- snow tendency due to freezing [kg/kg/s] REAL, INTENT(OUT), DIMENSION(klon) :: dqrfreez !-- rain tendency due to freezing [kg/kg/s] !--Local variables INTEGER :: i REAL :: hum_to_flux REAL :: dcldfra REAL :: precipfractot REAL :: dprecipfracclr, dprecipfraccld REAL :: drainclr, dsnowclr REAL :: draincld, dsnowcld REAL :: eff_cldfra REAL :: coef_col, coef_agg, coef_tmp, qrain REAL :: qthresh_auto_rain, tau_auto_rain, expo_auto_rain REAL :: qthresh_auto_snow, tau_auto_snow, expo_auto_snow REAL :: dqlcol ! loss of liquid cloud content due to collection by rain [kg/kg/s] REAL :: dqiagg ! loss of ice cloud content due to collection by aggregation [kg/kg/s] REAL :: dqlauto ! loss of liquid cloud content due to autoconversion to rain [kg/kg/s] REAL :: dqiauto ! loss of ice cloud content due to autoconversion to snow [kg/kg/s] REAL :: dqlrim ! loss of liquid cloud content due to riming on snow[kg/kg/s] !--Initialisation of variables dqrcol(:) = 0. dqsagg(:) = 0. dqsauto(:) = 0. dqrauto(:) = 0. dqsrim(:) = 0. dqrmelt(:) = 0. dqsmelt(:) = 0. dqrfreez(:) = 0. dqsfreez(:) = 0. DO i = 1, klon ! variables initialisation dqlrim = 0.0 dqlcol = 0.0 dqiagg = 0.0 dqiauto = 0.0 dqlauto = 0.0 !------------------------------------------------------------ !-- PRECIPITATION FRACTIONS UPDATE !------------------------------------------------------------ !--The goal of this routine is to reattribute precipitation fractions !--and fluxes to clear or cloudy air, depending on the variation of !--the cloud fraction on the vertical dimension. We assume a !--maximum-random overlap of the cloud cover (see Jakob and Klein, 2000, !--and LTP thesis, 2021) !--NB. in fact, we assume a maximum-random overlap of the total precip. frac !--Initialisation !--hum_to_flux: coef to convert a specific quantity to a flux !-- hum_to_flux = rho * dz/dt = 1 / g * dP/dt hum_to_flux = ( paprsdn(i) - paprsup(i) ) / RG / dtime precipfractot = precipfracclr(i) + precipfraccld(i) !--Instead of using the cloud cover which was use in LTP thesis, we use the !--total precip. fraction to compute the maximum-random overlap. This is !--because all the information of the cloud cover is embedded into !--precipfractot, and this allows for taking into account the potential !--reduction of the precipitation fraction because either the flux is too !--small (see barrier at the end of poprecip_postcld) or the flux is completely !--evaporated (see barrier at the end of poprecip_precld) !--NB. precipfraccld(i) is here the cloud fraction of the layer above precipfractot = 1. - ( 1. - precipfractot ) * & ( 1. - MAX( cldfra(i), precipfraccld(i) ) ) & / ( 1. - MIN( precipfraccld(i), 1. - eps ) ) !--precipfraccld(i) is here the cloud fraction of the layer above dcldfra = cldfra(i) - precipfraccld(i) !--Tendency of the clear-sky precipitation fraction. We add a MAX on the !--calculation of the current CS precip. frac. dprecipfracclr = MAX( 0., ( precipfractot - cldfra(i) ) ) - precipfracclr(i) !--Tendency of the cloudy precipitation fraction. We add a MAX on the !--calculation of the current CS precip. frac. !dprecipfraccld = MAX( dcldfra , - precipfraccld(i) ) !--We remove it, because cldfra is guaranteed to be > O (the MAX is activated !--if cldfra < 0) dprecipfraccld = dcldfra !--If the cloud extends IF ( dprecipfraccld .GT. 0. ) THEN !--If there is no CS precip, nothing happens. !--If there is, we reattribute some of the CS precip flux !--to the cloud precip flux, proportionnally to the !--decrease of the CS precip fraction IF ( precipfracclr(i) .LE. 0. ) THEN drainclr = 0. dsnowclr = 0. ELSE drainclr = dprecipfracclr / precipfracclr(i) * rainclr(i) dsnowclr = dprecipfracclr / precipfracclr(i) * snowclr(i) ENDIF !--If the cloud narrows ELSEIF ( dprecipfraccld .LT. 0. ) THEN !--We reattribute some of the cloudy precip flux !--to the CS precip flux, proportionnally to the !--decrease of the cloud precip fraction draincld = dprecipfraccld / precipfraccld(i) * raincld(i) dsnowcld = dprecipfraccld / precipfraccld(i) * snowcld(i) drainclr = - draincld dsnowclr = - dsnowcld !--If the cloud stays the same or if there is no cloud above and !--in the current layer, nothing happens ELSE drainclr = 0. dsnowclr = 0. ENDIF !--We add the tendencies precipfraccld(i) = precipfraccld(i) + dprecipfraccld precipfracclr(i) = precipfracclr(i) + dprecipfracclr rainclr(i) = rainclr(i) + drainclr snowclr(i) = snowclr(i) + dsnowclr raincld(i) = raincld(i) - drainclr snowcld(i) = snowcld(i) - dsnowclr ! if vertical heterogeneity is taken into account, we use ! the "true" volume fraction instead of a modified ! surface fraction (which is larger and artificially ! reduces the in-cloud water). IF ( ( iflag_cloudth_vert .GE. 3 ) .AND. ( iflag_rain_incloud_vol .EQ. 1 ) ) THEN eff_cldfra = ctot_vol(i) ELSE eff_cldfra = cldfra(i) ENDIF ! Start precipitation growth processes !--If the cloud is big enough, the precipitation processes activate IF ( cldfra(i) .GE. seuil_neb ) THEN !--------------------------------------------------------- !-- COLLECTION AND AGGREGATION !--------------------------------------------------------- !--Collection: processus through which rain collects small liquid droplets !--in suspension, and add it to the rain flux !--Aggregation: same for snow (precip flux) and ice crystals (in suspension) !--Those processes are treated before autoconversion because we do not !--want to collect/aggregate the newly formed fluxes, which already !--"saw" the cloud as they come from it !--gamma_col: tuning coefficient [-] !--rho_rain: volumic mass of rain [kg/m3] !--r_rain: size of the rain droplets [m] !--Eff_rain_liq: efficiency of the collection process [-] (between 0 and 1) !--dqlcol is a gridbox-mean quantity, as is qliq and raincld. They are !--divided by respectively eff_cldfra, eff_cldfra and precipfraccld to !--get in-cloud mean quantities. The two divisions by eff_cldfra are !--then simplified. coef_col = gamma_col * 3. / 4. / rho_rain / r_rain * Eff_rain_liq IF ((raincld(i) .GT. 0.) .AND. (coef_col .GT. 0.)) THEN !--Explicit version !dqlcol = - coef_col * qliq(i) * raincld(i) / precipfraccld(i) *dtime !--Semi-implicit version !dqlcol = qliq(i) * ( 1. / ( 1. + coef_col * raincld(i) / precipfraccld(i)*dtime ) - 1. ) !--Implicit version !qrain = raincld(i) / hum_to_flux !coef_tmp = coef_col * dtime *( qrain / precipfraccld(i) + qliq(i) / eff_cldfra ) !dqlcol = qliq(i) * ( 1. / ( 1. + 0.5 * ( coef_tmp - 1. + SQRT( & ! ( 1. - coef_tmp )**2. + 4. * coef_col * dtime *qrain / precipfraccld(i) ) & ! ) ) - 1. ) !dqlcol=max(dqlcol,-qliq(i)) ! Exact version dqlcol=qliq(i)*(exp(-dtime * coef_col * raincld(i) / precipfraccld(i))-1.) ENDIF !--Same as for aggregation coef_agg=gamma_agg * 3. / 4. / rho_snow / r_snow * Eff_snow_ice IF ((snowcld(i) .GT. 0.) .AND. (coef_agg .GT. 0.)) THEN !--Explicit version !dqiagg = - coef_agg & ! * qice(i) * snowcld(i) / precipfraccld(i) * dtime ! Exact version dqiagg=qice(i)*(exp(-dtime * coef_agg * snowcld(i) / precipfraccld(i))-1.) ENDIF !--Barriers so that the processes do not consume more liquid/ice than !--available. dqlcol = MAX( - qliq(i), dqlcol ) dqiagg = MAX( - qice(i), dqiagg ) !--Add tendencies qliq(i) = qliq(i) + dqlcol qice(i) = qice(i) + dqiagg raincld(i) = raincld(i) - dqlcol * hum_to_flux snowcld(i) = snowcld(i) - dqiagg * hum_to_flux !--------------------------------------------------------- !-- AUTOCONVERSION !--------------------------------------------------------- ! TODO IF ( ptconv(i) ) THEN ! if convective point qthresh_auto_rain = cld_lc_con qthresh_auto_snow = cld_lc_con tau_auto_rain = cld_tau_con tau_auto_snow = tau_auto_snow_max & + ( tau_auto_snow_min - tau_auto_snow_max ) * ( 1. - icefrac(i) ) expo_auto_rain = cld_expo_con expo_auto_snow = cld_expo_con ELSE qthresh_auto_rain = cld_lc_lsc qthresh_auto_snow = cld_lc_lsc tau_auto_rain = cld_tau_lsc tau_auto_snow = tau_auto_snow_max & + ( tau_auto_snow_min - tau_auto_snow_max ) * ( 1. - icefrac(i) ) expo_auto_rain = cld_expo_lsc expo_auto_snow = cld_expo_lsc ENDIF ! Liquid water quantity to remove according to (Sundqvist, 1978) ! dqliq/dt=-qliq/tau*(1-exp(-qcin/clw)**2) !......................................................... ! we first treat the second term (with exponential) in an explicit way ! and then treat the first term (-q/tau) in an exact way dqlauto = - qliq(i) * ( 1. - exp( - dtime / tau_auto_rain * ( 1. - exp( & - ( qliq(i) / eff_cldfra / qthresh_auto_rain ) ** expo_auto_rain ) ) ) ) dqiauto = - qice(i) * ( 1. - exp( - dtime / tau_auto_snow * ( 1. - exp( & - ( qice(i) / eff_cldfra / qthresh_auto_snow ) ** expo_auto_snow ) ) ) ) dqlauto = MAX( - qliq(i), dqlauto ) dqiauto = MAX( - qice(i), dqiauto ) qliq(i) = qliq(i) + dqlauto qice(i) = qice(i) + dqiauto raincld(i) = raincld(i) - dqlauto * hum_to_flux snowcld(i) = snowcld(i) - dqiauto * hum_to_flux ! FOLLOWING PROCESSES IMPLY A PHASE CHANGE SO A TEMPERATURE ! ADJUSTMENT !--------------------------------------------------------- !-- RIMING !--------------------------------------------------------- dqlrim=0.0 ! remplacer la premiere ligne par "coef_rim" ? IF (snowcld(i) .GT. 0.) THEN dqlrim = - gamma_rim * 3. / 4. / rho_snow / r_snow * Eff_snow_liq & * qliq(i) * snowcld(i) / precipfraccld(i) * dtime ENDIF dqlrim = MAX( - qliq(i), dqlrim ) qliq(i) = qliq(i) + dqlrim snowcld(i) = snowcld(i) - dqlrim * hum_to_flux ENDIF ! rneb .GE. seuil_neb !--------------------------------------------------------- !-- FREEZING !--------------------------------------------------------- !dqrfree_max = MINOUMAX(0., ( RTT - temp(i) ) / RLMLT * RCPD / ( 1. + RVTMP2 * ( qtot(i) + qprecip(i) ) )) !--------------------------------------------------------- !-- MELTING !--------------------------------------------------------- !flux = velocity * N_0 * 4. / 3. * PI * r_snow**3. * rho_snow !IF ( ( snowclr(i) + snowcld(i) ) .GT. 0. ) THEN ! dqsmelt_max = MIN(0., ( RTT - temp(i) ) / RLMLT * RCPD / ( 1. + RVTMP2 * ( qtot(i) + qprecip(i) ) )) ! dsnowtotmelt_max = dqsmelt_max * hum_to_flux(i) ! ! dsnowtotmelt = - nb_snowflake * 4. * PI * mol_diff_vap * snowflake_capa / RLMLT * coef_ventil & ! * MAX(0., ticebulb - RTT) & ! * ( paprsdn(i) - paprsup(i) ) / RG ! ! max bec. negative values ! dsnowtotmelt = MAX(dsnowtotmelt, dsnowmelt_max) ! dsnowclrmelt = dsnowtotmelt * snowclr(i) / ( snowclr(i) + snowcld(i) ) ! dsnowcldmelt = dsnowtotmelt - dsnowclrmelt ! ! update of rainfall and snowfall due to melting ! rainclr(i) = rainclr(i) - dsnowclrmelt(i) ! raincld(i) = raincld(i) - dsnowcldmelt(i) ! snowclr(i) = snowclr(i) + dsnowclrmelt(i) ! snowcld(i) = snowcld(i) + dsnowcldmelt(i) !ENDIF ! !! Latent heat of melting with precipitation thermalization !zt(i)=zt(i)-zifl(i)*zmelt*(RG*dtime)/(paprs(i,k)-paprs(i,k+1)) & !*RLMLT/RCPD/(1.0+RVTMP2*(zq(i)+zmqc(i))) !! MISE A JOUR DES FRACTIONS PRECIP CLD et CS ! LTP: limit of surface cloud fraction covered by precipitation when the local intensity of the flux is below rain_int_min precipfracclr(i) = MIN( precipfracclr(i), ( rainclr(i) + snowclr(i) ) / rain_int_min ) precipfraccld(i) = MIN( precipfraccld(i), ( raincld(i) + snowcld(i) ) / rain_int_min ) rain(i) = rainclr(i) + raincld(i) snow(i) = snowclr(i) + snowcld(i) ! write output tendencies for rain and snow dqsrim(i) = -dqlrim/dtime dqrcol(i) = -dqlcol/dtime dqsagg(i) = -dqiagg/dtime dqsauto(i) = -dqiauto/dtime dqrauto(i) = -dqlauto/dtime ENDDO END SUBROUTINE poprecip_postcld END MODULE lmdz_lscp_poprecip