Index: LMDZ5/trunk/libf/phylmd/fisrtilp.F90 =================================================================== --- LMDZ5/trunk/libf/phylmd/fisrtilp.F90 (revision 2499) +++ LMDZ5/trunk/libf/phylmd/fisrtilp.F90 (revision 2500) @@ -21,4 +21,10 @@ ! Objet: condensation et precipitation stratiforme. ! schema de nuage + ! Fusion de fisrt (physique sursaturation, P. LeVan K. Laval) + ! et ilp (il pleut, L. Li) + ! Principales parties: + ! P1> Evaporation de la precipitation (qui vient du niveau k+1) + ! P2> Formation du nuage (en k) + ! P3> Formation de la precipitation (en k) !====================================================================== !====================================================================== @@ -28,5 +34,5 @@ ! - ! Arguments: + ! Principaux inputs: ! REAL dtime ! intervalle du temps (s) @@ -35,4 +41,7 @@ REAL t(klon,klev) ! temperature (K) REAL q(klon,klev) ! humidite specifique (kg/kg) + ! + ! Principaux outputs: + ! REAL d_t(klon,klev) ! incrementation de la temperature (K) REAL d_q(klon,klev) ! incrementation de la vapeur d'eau @@ -46,4 +55,7 @@ REAL prfl(klon,klev+1) ! flux d'eau precipitante aux interfaces (kg/m2/s) REAL psfl(klon,klev+1) ! flux d'eau precipitante aux interfaces (kg/m2/s) + ! + ! Autres arguments + ! REAL ztv(klon,klev) REAL zqta(klon,klev),fraca(klon,klev) @@ -155,11 +167,10 @@ ! Pour la conversion eau-neige REAL zlh_solid(klon), zm_solid - !IM - !ym INTEGER klevm1 !--------------------------------------------------------------- ! ! Fonctions en ligne: ! - REAL fallvs,fallvc ! vitesse de chute pour crystaux de glace + REAL fallvs,fallvc ! Vitesse de chute pour cristaux de glace + ! (Heymsfield & Donner, 1990) REAL zzz include "YOETHF.h" @@ -278,15 +289,14 @@ zfrac_lessi = 0. - !AA---------------------------------------------------------- + !AA================================================================== ! ncoreczq=0 - ! Boucle verticale (du haut vers le bas) - ! - !IM : klevm1 - !ym klevm1=klev-1 + ! BOUCLE VERTICALE (DU HAUT VERS LE BAS) + ! DO k = klev, 1, -1 ! - !AA---------------------------------------------------------- - ! + !AA=============================================================== + ! + ! Initialisation temperature et vapeur DO i = 1, klon zt(i)=t(i,k) @@ -312,14 +322,7 @@ ENDDO ENDIF - ! - ! - ! Calculer l'evaporation de la precipitation - ! - - - ! Calculer l'evaporation de la precipitation - ! - - + ! ---------------------------------------------------------------- + ! P1> Debut evaporation de la precipitation + ! ---------------------------------------------------------------- IF (evap_prec) THEN DO i = 1, klon @@ -328,4 +331,5 @@ IF (zrfl(i)+zifl(i).GT.0.) THEN !>AJ + ! Calcul du qsat IF (thermcep) THEN zdelta=MAX(0.,SIGN(1.,RTT-zt(i))) @@ -350,22 +354,23 @@ IF (zrfl(i)+zifl(i).GT.0.) THEN !>AJ + ! Evap max pour ne pas saturer la fraction sous le nuage zqev = MAX (0.0, (zqs(i)-zq(i))*zneb(i) ) + ! Calcul de l'evaporation du flux de precip herite + ! d'au-dessus zqevt = coef_eva * (1.0-zq(i)/zqs(i)) * SQRT(zrfl(i)) & * (paprs(i,k)-paprs(i,k+1))/pplay(i,k)*zt(i)*RD/RG zqevt = MAX(0.0,MIN(zqevt,zrfl(i))) & * RG*dtime/(paprs(i,k)-paprs(i,k+1)) + ! Seuil pour ne pas saturer la fraction sous le nuage zqev = MIN (zqev, zqevt) + ! Nouveau flux de precip zrfln(i) = zrfl(i) - zqev*(paprs(i,k)-paprs(i,k+1)) & /RG/dtime - - ! pour la glace, on ré-évapore toute la précip dans la - ! couche du dessous - ! la glace venant de la couche du dessus est simplement - ! dans la couche du dessous. - + ! Aucun flux liquide pour T < t_coup IF (zt(i) .LT. t_coup.and.reevap_ice) zrfln(i)=0. - + ! Nouvelle vapeur zq(i) = zq(i) - (zrfln(i)-zrfl(i)) & * (RG/(paprs(i,k)-paprs(i,k+1)))*dtime + ! Nouvelle temperature (chaleur latente) zt(i) = zt(i) + (zrfln(i)-zrfl(i)) & * (RG/(paprs(i,k)-paprs(i,k+1)))*dtime & @@ -384,18 +389,19 @@ !>AJ !JAM !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - ! Modification de l'évaporation avec la glace - ! Différentiation entre précipitation liquide et solide - ! On suppose que coef_evai=2*coef_eva + ! Modification de l'evaporation avec la glace + ! Differentiation entre precipitation liquide et solide !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + ! Evap max pour ne pas saturer la fraction sous le nuage zqev0 = MAX (0.0, (zqs(i)-zq(i))*zneb(i) ) ! zqev0 = MAX (0.0, zqs(i)-zq(i) ) !JAM !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - ! On différencie qsat pour l'eau et la glace + ! On differencie qsat pour l'eau et la glace ! Si zdelta=1. --> glace ! Si zdelta=0. --> eau liquide !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - + + ! Calcul du qsat par rapport a l'eau liquide qsl= R2ES*FOEEW(zt(i),0.)/pplay(i,k) qsl= MIN(0.5,qsl) @@ -403,4 +409,8 @@ qsl= qsl*zcor + ! Calcul de l'evaporation du flux de precip herite + ! d'au-dessus + ! Formulation en racine du flux de precip + ! (Klemp & Wilhelmson, 1978; Sundqvist, 1988) zqevt = 1.*coef_eva*(1.0-zq(i)/qsl)*SQRT(zrfl(i)) & *(paprs(i,k)-paprs(i,k+1))/pplay(i,k)*zt(i)*RD/RG @@ -408,4 +418,6 @@ *RG*dtime/(paprs(i,k)-paprs(i,k+1)) + + ! Calcul du qsat par rapport a la glace qsi= R2ES*FOEEW(zt(i),1.)/pplay(i,k) qsi= MIN(0.5,qsi) @@ -413,4 +425,6 @@ qsi= qsi*zcor + ! Calcul de la sublimation du flux de precip solide herite + ! d'au-dessus zqevti = 1.*coef_eva*(1.0-zq(i)/qsi)*SQRT(zifl(i)) & *(paprs(i,k)-paprs(i,k+1))/pplay(i,k)*zt(i)*RD/RG @@ -420,5 +434,9 @@ !JAM!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - ! Vérification sur l'évaporation + ! Verification sur l'evaporation + ! On s'assure qu'on ne sature pas + ! la fraction sous le nuage sinon on + ! repartit zqev0 en gardant la proportion + ! liquide / glace !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! @@ -436,16 +454,11 @@ ENDIF ENDIF - + ! Nouveaux flux de precip liquide et solide zrfln(i) = Max(0.,zrfl(i) - zqev*(paprs(i,k)-paprs(i,k+1)) & /RG/dtime) zifln(i) = Max(0.,zifl(i) - zqevi*(paprs(i,k)-paprs(i,k+1)) & /RG/dtime) - - ! Pour la glace, on révapore toute la précip dans la couche du dessous - ! la glace venant de la couche du dessus est simplement dans la couche - ! du dessous. - - ! IF (zt(i) .LT. t_coup.and.reevap_ice) zrfln(i)=0. -! print*,zrfl(i),zrfln(i),zqevt,zqevti,RLMLT,'fluxdeprecip' + + ! Mise a jour de la vapeur, temperature et flux de precip zq(i) = zq(i) - (zrfln(i)+zifln(i)-zrfl(i)-zifl(i)) & * (RG/(paprs(i,k)-paprs(i,k+1)))*dtime @@ -466,7 +479,9 @@ zmelt = ((zt(i)-273.15)/(ztfondue-273.15)) ! jyg zmelt = MIN(MAX(zmelt,0.),1.) + ! Fusion de la glace zrfl(i)=zrfl(i)+zmelt*zifl(i) zifl(i)=zifl(i)*(1.-zmelt) ! print*,zt(i),'octavio1' + ! Chaleur latente de fusion zt(i)=zt(i)-zifl(i)*zmelt*(RG*dtime)/(paprs(i,k)-paprs(i,k+1)) & *RLMLT/RCPD/(1.0+RVTMP2*zq(i)) @@ -481,4 +496,7 @@ ENDIF ! (.NOT. ice_thermo) + ! ---------------------------------------------------------------- + ! Fin evaporation de la precipitation + ! ---------------------------------------------------------------- ENDIF ! (evap_prec) ! @@ -518,4 +536,7 @@ ! endif + ! ---------------------------------------------------------------- + ! P2> Formation du nuage + ! ---------------------------------------------------------------- IF (cpartiel) THEN @@ -534,4 +555,5 @@ if (iflag_pdf.eq.0) then + ! version creneau de (Li, 1998) do i=1,klon zdelq = min(ratqs(i,k),0.99) * zq(i) @@ -575,5 +597,5 @@ end if -!CR: variation de qsat avec T en présence de glace ou non +!CR: variation de qsat avec T en presence de glace ou non !initialisations do i=1,klon @@ -587,4 +609,6 @@ !Boucle iterative: ATTENTION, l'option -1 n'est plus activable ici if (iflag_fisrtilp_qsat.ge.0) then + ! Iteration pour condensation avec variation de qsat(T) + ! ----------------------------------------------------- do iter=1,iflag_fisrtilp_qsat+1 @@ -603,4 +627,5 @@ endif endif + ! Calcul des PDF lognormales zcvm5 = R5LES*RLVTT*(1.-zdelta) + R5IES*RLSTT*zdelta zcvm5 = zcvm5 /RCPD/(1.0+RVTMP2*zq(i)) @@ -647,12 +672,5 @@ else - -!calcul de la fraction de glace -!CR: on utilise la nouvelle fonction de JBM pour l ancien calcul -! zfice(i) = icefrac_lsc(Tbef(i), t_glace_min, & -! t_glace_max, exposant_glace) -! zfice(i) = 1.0 - (Tbef(i)-ztglace) / (RTT-ztglace) -! zfice(i) = MIN(MAX(zfice(i),0.0),1.0) -! zfice(i) = zfice(i)**nexpo + ! Iteration pour convergence avec qsat(T) if (iflag_t_glace.eq.1) then CALL icefrac_lsc(klon,zt(:),pplay(:,k)/paprs(:,1),zfice(:)) @@ -699,5 +717,7 @@ - enddo + enddo ! iter=1,iflag_fisrtilp_qsat+1 + ! Fin d'iteration pour condensation avec variation de qsat(T) + ! ----------------------------------------------------------- endif @@ -707,5 +727,8 @@ ! if (iflag_fisrtilp_qsat.eq.-1) then -!CR: ATTENTION option fausse mais a existe: pour la re-activer, prendre iflag_fisrtilp_qsat=0 et activer les lignes suivantes: +!------------------------------------------ +!CR: ATTENTION option fausse mais a existe: +! pour la re-activer, prendre iflag_fisrtilp_qsat=0 et +! activer les lignes suivantes: IF (1.eq.0) THEN DO i=1,klon @@ -724,9 +747,14 @@ rhcl(i,k)=(zqs(i)+zq(i)-zdelq)/2./zqs(i) ENDIF - ENDDO - ENDIF + ENDDO + ENDIF +!------------------------------------------ ! ELSE + ! Calcul de l'eau in-cloud (zqn), + ! moyenne dans la maille (zcond), + ! fraction nuageuse (rneb) et + ! humidite relative ciel-clair (rhcl) DO i=1,klon IF (rneb(i,k) .LE. 0.0) THEN @@ -749,21 +777,6 @@ ! ENDIF - ! do i=1,klon - ! IF (rneb(i,k) .LE. 0.0) zqn(i) = 0.0 - ! IF (rneb(i,k) .GE. 1.0) zqn(i) = zq(i) - ! rneb(i,k) = MAX(0.0,MIN(1.0,rneb(i,k))) - !c zcond(i) = MAX(0.0,zqn(i)-zqs(i))*rneb(i,k)/(1.+zdqs(i)) - !c On ne divise pas par 1+zdqs pour forcer a avoir l'eau predite par - !c la convection. - !c ATTENTION !!! Il va falloir verifier tout ca. - ! zcond(i) = MAX(0.0,zqn(i)-zqs(i))*rneb(i,k) - !c print*,'ZDQS ',zdqs(i) - !c--Olivier - ! rhcl(i,k)=(zqs(i)+zq(i)-zdelq)/2./zqs(i) - ! IF (rneb(i,k) .LE. 0.0) rhcl(i,k)=zq(i)/zqs(i) - ! IF (rneb(i,k) .GE. 1.0) rhcl(i,k)=1.0 - !c--fin - ! ENDDO - ELSE + ELSE ! de IF (cpartiel) + ! Cas "tout ou rien" DO i = 1, klon IF (zq(i).GT.zqs(i)) THEN @@ -775,5 +788,9 @@ ENDDO ENDIF - ! + ! ---------------------------------------------------------------- + ! Fin de formation du nuage + ! ---------------------------------------------------------------- + ! + ! Mise a jour vapeur d'eau DO i = 1, klon zq(i) = zq(i) - zcond(i) @@ -781,4 +798,6 @@ ENDDO !AJ< + ! Chaleur latente apres formation nuage + ! ------------------------------------- IF (.NOT. ice_thermo) THEN if (iflag_fisrtilp_qsat.lt.1) then @@ -825,7 +844,12 @@ ENDIF !>AJ + ! ---------------------------------------------------------------- + ! P3> Formation des precipitations + ! ---------------------------------------------------------------- ! ! Partager l'eau condensee en precipitation et eau liquide nuageuse ! + + ! Initialisation de zoliq (eau condensee moyenne dans la maille) DO i = 1, klon IF (rneb(i,k).GT.0.0) THEN @@ -858,4 +882,11 @@ ENDIF ENDIF + + ! Calcul de radliq (eau condensee pour le rayonnement) + ! Iteration pour realiser une moyenne de l'eau nuageuse lors de la precip + ! Remarque: ce n'est donc pas l'eau restante en fin de precip mais une + ! eau moyenne restante dans le nuage sur la duree du pas de temps qui est + ! transmise au rayonnement; + ! ---------------------------------------------------------------- DO i = 1, klon IF (rneb(i,k).GT.0.0) THEN @@ -878,6 +909,6 @@ ztot = 0.0 ELSE - ! quantite d'eau a eliminer: zchau - ! meme chose pour la glace: zfroi + ! quantite d'eau a eliminer: zchau (Sundqvist, 1978) + ! meme chose pour la glace: zfroi (Zender & Kiehl, 1997) if (ptconv(i,k)) then zcl =cld_lc_con @@ -916,4 +947,5 @@ ENDDO ! i = 1,klon ENDDO ! n = 1,ninter + ! ---------------------------------------------------------------- ! IF (.NOT. ice_thermo) THEN @@ -1028,5 +1060,5 @@ ELSE ! JAM************************************************* - ! Revoir partie ci-dessous: à quoi servent psfl et prfl? + ! Revoir partie ci-dessous: a quoi servent psfl et prfl? ! ***************************************************** @@ -1040,4 +1072,7 @@ ENDDO ENDIF + ! ---------------------------------------------------------------- + ! Fin de formation des precipitations + ! ---------------------------------------------------------------- ! ! @@ -1110,9 +1145,9 @@ ENDDO ! - !AA---------------------------------------------------------- - ! FIN DE BOUCLE SUR K + !AA=============================================================== + ! FIN DE LA BOUCLE VERTICALE end DO ! - !AA----------------------------------------------------------- + !AA================================================================== ! ! Pluie ou neige au sol selon la temperature de la 1ere couche