! $Id: nuage.F90 5153 2024-07-31 16:20:03Z fhourdin $ SUBROUTINE nuage(paprs, pplay, t, pqlwp, picefra, pclc, pcltau, pclemi, pch, pcl, pcm, & pct, pctlwp, ok_aie, mass_solu_aero, mass_solu_aero_pi, bl95_b0, bl95_b1, distcltop, & temp_cltop, cldtaupi, re, fl) USE dimphy USE lmdz_lscp_tools, ONLY: icefrac_lscp USE icefrac_lsc_mod ! computes ice fraction (JBM 3/14) USE lmdz_lscp_ini, ONLY: iflag_t_glace USE phys_local_var_mod, ONLY: ptconv USE lmdz_clesphys USE lmdz_nuage_params ! JBM 3/14 USE lmdz_yomcst IMPLICIT NONE ! ====================================================================== ! Auteur(s): Z.X. Li (LMD/CNRS) date: 19930910 ! Objet: Calculer epaisseur optique et emmissivite des nuages ! ====================================================================== ! Arguments: ! t-------input-R-temperature ! pqlwp---input-R-eau liquide nuageuse dans l'atmosphere (kg/kg) ! picefra--inout-R-fraction de glace dans les nuages (-) ! pclc----input-R-couverture nuageuse pour le rayonnement (0 a 1) ! ok_aie--input-L-apply aerosol indirect effect or not ! mass_solu_aero-----input-R-total mass concentration for all soluble ! aerosols[ug/m^3] ! mass_solu_aero_pi--input-R-dito, pre-industrial value ! bl95_b0-input-R-a parameter, may be varied for tests (s-sea, l-land) ! bl95_b1-input-R-a parameter, may be varied for tests ( -"- ) ! cldtaupi-output-R-pre-industrial value of cloud optical thickness, ! needed for the diagnostics of the aerosol indirect ! radiative forcing (see radlwsw) ! re------output-R-Cloud droplet effective radius multiplied by fl [um] ! fl------output-R-Denominator to re, introduced to avoid problems in ! the averaging of the output. fl is the fraction of liquid ! water clouds within a grid cell ! pcltau--output-R-epaisseur optique des nuages ! pclemi--output-R-emissivite des nuages (0 a 1) ! ====================================================================== REAL paprs(klon, klev + 1), pplay(klon, klev) REAL t(klon, klev) REAL pclc(klon, klev) REAL pqlwp(klon, klev), picefra(klon, klev) REAL pcltau(klon, klev), pclemi(klon, klev) REAL pct(klon), pctlwp(klon), pch(klon), pcl(klon), pcm(klon) REAL distcltop(klon, klev) REAL temp_cltop(klon, klev) LOGICAL lo REAL cetahb, cetamb PARAMETER (cetahb = 0.45, cetamb = 0.80) INTEGER i, k REAL zflwp, zradef, zfice(klon), zmsac REAL radius, rad_chaud ! JBM (3/14) parameters already defined in nuage.h: ! REAL rad_froid, rad_chau1, rad_chau2 ! PARAMETER (rad_chau1=13.0, rad_chau2=9.0, rad_froid=35.0) ! cc PARAMETER (rad_chaud=15.0, rad_froid=35.0) ! sintex initial PARAMETER (rad_chaud=10.0, rad_froid=30.0) REAL coef, coef_froi, coef_chau PARAMETER (coef_chau = 0.13, coef_froi = 0.09) REAL seuil_neb PARAMETER (seuil_neb = 0.001) ! JBM (3/14) nexpo is replaced by exposant_glace ! REAL nexpo ! exponentiel pour glace/eau ! PARAMETER (nexpo=6.) REAL, PARAMETER :: t_glace_min_old = 258. INTEGER, PARAMETER :: exposant_glace_old = 6 ! jq for the aerosol indirect effect ! jq introduced by Johannes Quaas (quaas@lmd.jussieu.fr), 27/11/2003 ! jq LOGICAL ok_aie ! Apply AIE or not? REAL mass_solu_aero(klon, klev) ! total mass concentration for all soluble aerosols[ug m-3] REAL mass_solu_aero_pi(klon, klev) ! - " - pre-industrial value REAL cdnc(klon, klev) ! cloud droplet number concentration [m-3] REAL re(klon, klev) ! cloud droplet effective radius [um] REAL cdnc_pi(klon, klev) ! cloud droplet number concentration [m-3] (pi value) REAL re_pi(klon, klev) ! cloud droplet effective radius [um] (pi value) REAL fl(klon, klev) ! xliq * rneb (denominator to re; fraction of liquid water clouds ! within the grid cell) REAL bl95_b0, bl95_b1 ! Parameter in B&L 95-Formula REAL cldtaupi(klon, klev) ! pre-industrial cloud opt thickness for diag REAl dzfice(klon) ! jq-end ! cc PARAMETER (nexpo=1) ! Calculer l'epaisseur optique et l'emmissivite des nuages DO k = 1, klev IF (iflag_t_glace==0) THEN DO i = 1, klon zfice(i) = 1.0 - (t(i, k) - t_glace_min_old) / (273.13 - t_glace_min_old) zfice(i) = min(max(zfice(i), 0.0), 1.0) zfice(i) = zfice(i)**exposant_glace_old ENDDO ELSE ! of IF (iflag_t_glace.EQ.0) ! JBM: icefrac_lsc is now a function contained in icefrac_lsc_mod ! zfice(i) = icefrac_lsc(t(i,k), t_glace_min, & ! t_glace_max, exposant_glace) IF (ok_new_lscp) THEN CALL icefrac_lscp(klon, t(:, k), iflag_ice_thermo, distcltop(:, k), temp_cltop(:, k), zfice(:), dzfice(:)) ELSE CALL icefrac_lsc(klon, t(:, k), pplay(:, k) / paprs(:, 1), zfice(:)) ENDIF IF (ok_new_lscp .AND. ok_icefra_lscp) THEN ! EV: take the ice fraction directly from the lscp code ! consistent only for non convective grid points ! critical for mixed phase clouds DO i = 1, klon IF (.NOT. ptconv(i, k)) THEN zfice(i) = picefra(i, k) ENDIF ENDDO ENDIF ENDIF DO i = 1, klon rad_chaud = rad_chau1 IF (k<=3) rad_chaud = rad_chau2 pclc(i, k) = max(pclc(i, k), seuil_neb) zflwp = 1000. * pqlwp(i, k) / rg / pclc(i, k) * (paprs(i, k) - paprs(i, k + 1)) IF (ok_aie) THEN ! Formula "D" of Boucher and Lohmann, Tellus, 1995 cdnc(i, k) = 10.**(bl95_b0 + bl95_b1 * log(max(mass_solu_aero(i, k), & 1.E-4)) / log(10.)) * 1.E6 !-m-3 ! Cloud droplet number concentration (CDNC) is restricted ! to be within [20, 1000 cm^3] cdnc(i, k) = min(1000.E6, max(20.E6, cdnc(i, k))) cdnc_pi(i, k) = 10.**(bl95_b0 + bl95_b1 * log(max(mass_solu_aero_pi(i, k), & 1.E-4)) / log(10.)) * 1.E6 !-m-3 cdnc_pi(i, k) = min(1000.E6, max(20.E6, cdnc_pi(i, k))) ! air density: pplay(i,k) / (RD * zT(i,k)) ! factor 1.1: derive effective radius from volume-mean radius ! factor 1000 is the water density ! _chaud means that this is the CDR for liquid water clouds rad_chaud = 1.1 * ((pqlwp(i, k) * pplay(i, k) / (rd * t(i, k))) / (4. / 3. * rpi * 1000. & * cdnc(i, k)))**(1. / 3.) ! Convert to um. CDR shall be at least 3 um. rad_chaud = max(rad_chaud * 1.E6, 3.) ! For output diagnostics ! Cloud droplet effective radius [um] ! we multiply here with f * xl (fraction of liquid water ! clouds in the grid cell) to avoid problems in the ! averaging of the output. ! In the output of IOIPSL, derive the real cloud droplet ! effective radius as re/fl fl(i, k) = pclc(i, k) * (1. - zfice(i)) re(i, k) = rad_chaud * fl(i, k) ! Pre-industrial cloud opt thickness ! "radius" is calculated as rad_chaud above (plus the ! ice cloud contribution) but using cdnc_pi instead of ! cdnc. radius = max(1.1E6 * ((pqlwp(i, k) * pplay(i, k) / (rd * t(i, k))) / (4. / 3. * rpi * & 1000. * cdnc_pi(i, k)))**(1. / 3.), 3.) * (1. - zfice(i)) + rad_froid * zfice(i) cldtaupi(i, k) = 3.0 / 2.0 * zflwp / radius END IF ! ok_aie radius = rad_chaud * (1. - zfice(i)) + rad_froid * zfice(i) coef = coef_chau * (1. - zfice(i)) + coef_froi * zfice(i) pcltau(i, k) = 3.0 / 2.0 * zflwp / radius pclemi(i, k) = 1.0 - exp(-coef * zflwp) lo = (pclc(i, k)<=seuil_neb) IF (lo) pclc(i, k) = 0.0 IF (lo) pcltau(i, k) = 0.0 IF (lo) pclemi(i, k) = 0.0 IF (.NOT. ok_aie) cldtaupi(i, k) = pcltau(i, k) END DO END DO ! cc DO k = 1, klev ! cc DO i = 1, klon ! cc t(i,k) = t(i,k) ! cc pclc(i,k) = MAX( 1.e-5 , pclc(i,k) ) ! cc lo = pclc(i,k) .GT. (2.*1.e-5) ! cc zflwp = pqlwp(i,k)*1000.*(paprs(i,k)-paprs(i,k+1)) ! cc . /(rg*pclc(i,k)) ! cc zradef = 10.0 + (1.-sigs(k))*45.0 ! cc pcltau(i,k) = 1.5 * zflwp / zradef ! cc zfice=1.0-MIN(MAX((t(i,k)-263.)/(273.-263.),0.0),1.0) ! cc zmsac = 0.13*(1.0-zfice) + 0.08*zfice ! cc pclemi(i,k) = 1.-EXP(-zmsac*zflwp) ! cc if (.NOT.lo) pclc(i,k) = 0.0 ! cc if (.NOT.lo) pcltau(i,k) = 0.0 ! cc if (.NOT.lo) pclemi(i,k) = 0.0 ! cc ENDDO ! cc ENDDO ! ccccc PRINT*, 'pas de nuage dans le rayonnement' ! ccccc DO k = 1, klev ! ccccc DO i = 1, klon ! ccccc pclc(i,k) = 0.0 ! ccccc pcltau(i,k) = 0.0 ! ccccc pclemi(i,k) = 0.0 ! ccccc ENDDO ! ccccc ENDDO ! COMPUTE CLOUD LIQUID PATH AND TOTAL CLOUDINESS DO i = 1, klon pct(i) = 1.0 pch(i) = 1.0 pcm(i) = 1.0 pcl(i) = 1.0 pctlwp(i) = 0.0 END DO DO k = klev, 1, -1 DO i = 1, klon pctlwp(i) = pctlwp(i) + pqlwp(i, k) * (paprs(i, k) - paprs(i, k + 1)) / rg pct(i) = pct(i) * (1.0 - pclc(i, k)) IF (pplay(i, k)<=cetahb * paprs(i, 1)) pch(i) = pch(i) * (1.0 - pclc(i, k)) IF (pplay(i, k)>cetahb * paprs(i, 1) .AND. pplay(i, k)<=cetamb * paprs(i, 1)) & pcm(i) = pcm(i) * (1.0 - pclc(i, k)) IF (pplay(i, k)>cetamb * paprs(i, 1)) pcl(i) = pcl(i) * (1.0 - pclc(i, k)) END DO END DO DO i = 1, klon pct(i) = 1. - pct(i) pch(i) = 1. - pch(i) pcm(i) = 1. - pcm(i) pcl(i) = 1. - pcl(i) END DO END SUBROUTINE nuage SUBROUTINE diagcld1(paprs, pplay, rain, snow, kbot, ktop, diafra, dialiq) USE dimphy USE lmdz_yomcst IMPLICIT NONE ! Laurent Li (LMD/CNRS), le 12 octobre 1998 ! (adaptation du code ECMWF) ! Dans certains cas, le schema pronostique des nuages n'est ! pas suffisament performant. On a donc besoin de diagnostiquer ! ces nuages. Je dois avouer que c'est une frustration. ! Arguments d'entree: REAL paprs(klon, klev + 1) ! pression (Pa) a inter-couche REAL pplay(klon, klev) ! pression (Pa) au milieu de couche REAL t(klon, klev) ! temperature (K) REAL q(klon, klev) ! humidite specifique (Kg/Kg) REAL rain(klon) ! pluie convective (kg/m2/s) REAL snow(klon) ! neige convective (kg/m2/s) INTEGER ktop(klon) ! sommet de la convection INTEGER kbot(klon) ! bas de la convection ! Arguments de sortie: REAL diafra(klon, klev) ! fraction nuageuse diagnostiquee REAL dialiq(klon, klev) ! eau liquide nuageuse ! Constantes ajustables: REAL canva, canvb, canvh PARAMETER (canva = 2.0, canvb = 0.3, canvh = 0.4) REAL cca, ccb, ccc PARAMETER (cca = 0.125, ccb = 1.5, ccc = 0.8) REAL ccfct, ccscal PARAMETER (ccfct = 0.400) PARAMETER (ccscal = 1.0E+11) REAL cetahb, cetamb PARAMETER (cetahb = 0.45, cetamb = 0.80) REAL cclwmr PARAMETER (cclwmr = 1.E-04) REAL zepscr PARAMETER (zepscr = 1.0E-10) ! Variables locales: INTEGER i, k REAL zcc(klon) ! Initialisation: DO k = 1, klev DO i = 1, klon diafra(i, k) = 0.0 dialiq(i, k) = 0.0 END DO END DO DO i = 1, klon ! Calculer la fraction nuageuse zcc(i) = 0.0 IF ((rain(i) + snow(i))>0.) THEN zcc(i) = cca * log(max(zepscr, (rain(i) + snow(i)) * ccscal)) - ccb zcc(i) = min(ccc, max(0.0, zcc(i))) END IF END DO DO i = 1, klon ! pour traiter les enclumes diafra(i, ktop(i)) = max(diafra(i, ktop(i)), zcc(i) * ccfct) IF ((zcc(i)>=canvh) .AND. (pplay(i, ktop(i))<=cetahb * paprs(i, & 1))) diafra(i, ktop(i)) = max(diafra(i, ktop(i)), max(zcc(& i) * ccfct, canva * (zcc(i) - canvb))) dialiq(i, ktop(i)) = cclwmr * diafra(i, ktop(i)) END DO DO k = 1, klev ! nuages convectifs (sauf enclumes) DO i = 1, klon IF (k=kbot(i)) THEN diafra(i, k) = max(diafra(i, k), zcc(i) * ccfct) dialiq(i, k) = cclwmr * diafra(i, k) END IF END DO END DO END SUBROUTINE diagcld1 SUBROUTINE diagcld2(paprs, pplay, t, q, diafra, dialiq) USE dimphy USE lmdz_yoethf USE lmdz_yomcst IMPLICIT NONE INCLUDE "FCTTRE.h" ! Arguments d'entree: REAL paprs(klon, klev + 1) ! pression (Pa) a inter-couche REAL pplay(klon, klev) ! pression (Pa) au milieu de couche REAL t(klon, klev) ! temperature (K) REAL q(klon, klev) ! humidite specifique (Kg/Kg) ! Arguments de sortie: REAL diafra(klon, klev) ! fraction nuageuse diagnostiquee REAL dialiq(klon, klev) ! eau liquide nuageuse REAL cetamb PARAMETER (cetamb = 0.80) REAL cloia, cloib, cloic, cloid PARAMETER (cloia = 1.0E+02, cloib = -10.00, cloic = -0.6, cloid = 5.0) ! cc PARAMETER (CLOIA=1.0E+02, CLOIB=-10.00, CLOIC=-0.9, CLOID=5.0) REAL rgammas PARAMETER (rgammas = 0.05) REAL crhl PARAMETER (crhl = 0.15) ! cc PARAMETER (CRHL=0.70) REAL t_coup PARAMETER (t_coup = 234.0) ! Variables locales: INTEGER i, k, kb, invb(klon) REAL zqs, zrhb, zcll, zdthmin(klon), zdthdp REAL zdelta, zcor ! Initialisation: DO k = 1, klev DO i = 1, klon diafra(i, k) = 0.0 dialiq(i, k) = 0.0 END DO END DO DO i = 1, klon invb(i) = klev zdthmin(i) = 0.0 END DO DO k = 2, klev / 2 - 1 DO i = 1, klon zdthdp = (t(i, k) - t(i, k + 1)) / (pplay(i, k) - pplay(i, k + 1)) - & rd * 0.5 * (t(i, k) + t(i, k + 1)) / rcpd / paprs(i, k + 1) zdthdp = zdthdp * cloia IF (pplay(i, k)>cetamb * paprs(i, 1) .AND. zdthdp0.0 .AND. zrhb