Ignore:
Timestamp:
May 10, 2019, 12:17:35 PM (6 years ago)
Author:
jyg
Message:

Implementation of the ejection of liquid precipitation from the adiabatic ascents.
New flags:
+cvflag_prec_eject: logical

n -> old code, y -> new code

+ejectliq: real; possible values 0. & 1.

  1. -> no liquid precipitation is ejected
  2. -> all liquid precipitation is ejected

+ejectice: real; any value between 0. and 1.

fraction of solid precipitation ejected at each level

Note that the adiabatic ascent mass flux decrease due to precipitation ejection is not taken into account.

Attempts to do it led to water conservation violation.

Location:
LMDZ6/trunk/libf/phylmd
Files:
14 edited

Legend:

Unmodified
Added
Removed

  • LMDZ6/trunk/libf/phylmd/concvl.F90

    r3197 r3496  
    1313!RomP >>>
    1414!!     .             da,phi,mp,dd_t,dd_q,lalim_conv,wght_th)
    15                   da, phi, mp, phii, d1a, dam, sij, clw, elij, &     ! RomP
     15                  da, phi, mp, phii, d1a, dam, sij, qta, clw, elij, &! RomP
    1616                  dd_t, dd_q, lalim_conv, wght_th, &                 ! RomP
    1717                  evap, ep, epmlmMm, eplaMm, &                       ! RomP
    18                   wdtrainA, wdtrainM, wght, qtc, sigt, &
     18                  wdtrainA, wdtrainS, wdtrainM, wght, qtc, sigt, &
    1919                  tau_cld_cv, coefw_cld_cv, &                           ! RomP+RL, AJ
    2020!RomP <<<
    … …  
    8080! eplaMm-----output-R
    8181! wdtrainA---output-R
     82! wdtrainS---output-R
    8283! wdtrainM---output-R
    8384! wght-------output-R
    … …  
    134135  REAL, DIMENSION(klon,klev),   INTENT(OUT)     :: d1a, dam
    135136  REAL, DIMENSION(klon,klev,klev),INTENT(OUT)   :: sij, elij
     137  REAL, DIMENSION(klon,klev),   INTENT(OUT)     :: qta
    136138  REAL, DIMENSION(klon,klev),   INTENT(OUT)     :: clw
    137139  REAL, DIMENSION(klon,klev),   INTENT(OUT)     :: dd_t, dd_q
    … …  
    139141  REAL, DIMENSION(klon,klev),   INTENT(OUT)     :: eplaMm
    140142  REAL, DIMENSION(klon,klev,klev), INTENT(OUT)  :: epmlmMm
    141   REAL, DIMENSION(klon,klev),   INTENT(OUT)     :: wdtrainA, wdtrainM
     143  REAL, DIMENSION(klon,klev),   INTENT(OUT)     :: wdtrainA, wdtrainS, wdtrainM
    142144! RomP <<<
    143145  REAL, DIMENSION(klon,klev),   INTENT(OUT)     :: wght                       !RL
    … …  
    437439!!                   evap,ep,epmlmMm,eplaMm,                              ! RomP
    438440                    da, phi, mp, phii, d1a, dam, sij, wght, &           ! RomP+RL
    439                     clw, elij, evap, ep, epmlmMm, eplaMm, &             ! RomP+RL
    440                     wdtrainA, wdtrainM, qtc, sigt, &
     441                    qta, clw, elij, evap, ep, epmlmMm, eplaMm, &        ! RomP+RL
     442                    wdtrainA, wdtrainS, wdtrainM, qtc, sigt, &
    441443                    tau_cld_cv, coefw_cld_cv, &                         ! RomP,AJ
    442444!AC!+!RomP+jyg

  • LMDZ6/trunk/libf/phylmd/cv3_mixscale.F90

    r1992 r3496  
    1313  include "cv3param.h"
    1414
    15   INTEGER nloc, ncum, na
     15!inputs:
     16  INTEGER, INTENT (IN)                               :: ncum, na, nloc
     17  REAL, DIMENSION (nloc, na), INTENT (IN)            :: m
     18!input/outputs:
     19  REAL, DIMENSION (nloc, na, na), INTENT (INOUT)     :: ment
     20
     21!local variables:
    1622  INTEGER i, j, il
    17   REAL ment(nloc, na, na), m(nloc, na)
    1823
    19   DO j = 1, nl
    20     DO i = 1, nl
    21       DO il = 1, ncum
    22         ment(il, i, j) = m(il, i)*ment(il, i, j)
     24    DO j = 1, nl
     25      DO i = 1, nl
     26        DO il = 1, ncum
     27          ment(il, i, j) = m(il, i)*ment(il, i, j)
     28        END DO
    2329      END DO
    2430    END DO
    25   END DO
    2631
    2732

  • LMDZ6/trunk/libf/phylmd/cv3_routines.F90

    r3492 r3496  
    3535
    3636  include "cv3param.h"
     37  include "cvflag.h"
    3738  include "conema3.h"
    3839
    … …  
    125126     tlcrit=-55.0
    126127     CALL getin_p('tlcrit',tlcrit)
     128     ejectliq=0.
     129     CALL getin_p('ejectliq',ejectliq)
     130     ejectice=0.
     131     CALL getin_p('ejectice',ejectice)
     132     cvflag_prec_eject = .FALSE.
     133     CALL getin_p('cvflag_prec_eject',cvflag_prec_eject)
    127134     qsat_depends_on_qt = .FALSE.
    128135     CALL getin_p('qsat_depends_on_qt',qsat_depends_on_qt)
     136     adiab_ascent_mass_flux_depends_on_ejectliq = .FALSE.
     137     CALL getin_p('adiab_ascent_mass_flux_depends_on_ejectliq',adiab_ascent_mass_flux_depends_on_ejectliq)
    129138
    130139    WRITE (*, *) 't_top_max=', t_top_max
    … …  
    172181
    173182  include "cv3param.h"
     183  include "cvflag.h"
    174184
    175185!inputs:
    … …  
    292302  USE mod_phys_lmdz_transfert_para, ONLY : bcast
    293303  USE add_phys_tend_mod, ONLY: fl_cor_ebil
     304  USE print_control_mod, ONLY: prt_level
    294305  IMPLICIT NONE
    295306
    … …  
    519530    END DO
    520531  ENDIF
     532  IF (prt_level .GE. 10) THEN
     533    print *,'cv3_feed : iflag(1), pfeed(1), plcl(1), wghti(1,k) ', &
     534                        iflag(1), pfeed(1), plcl(1), (wghti(1,k),k=1,10)
     535  ENDIF
    521536
    522537! -------------------------------------------------------------------
    … …  
    11081123                         tnk, qnk, gznk, hnk, t, q, qs, gz, &
    11091124                         p, ph, h, tv, lv, lf, pbase, buoybase, plcl, &
    1110                          inb, tp, tvp, clw, hp, ep, sigp, buoy, frac)
     1125                         inb, tp, tvp, clw, hp, ep, sigp, buoy, &
     1126                         frac_a, frac_s, qpreca, qta)
    11111127  USE print_control_mod, ONLY: prt_level
    11121128  IMPLICIT NONE
    … …  
    11561172  REAL, DIMENSION (nloc, nd), INTENT (OUT)           :: ep, sigp, hp
    11571173  REAL, DIMENSION (nloc, nd), INTENT (OUT)           :: buoy
    1158   REAL, DIMENSION (nloc, nd), INTENT (OUT)           :: frac   ! ice fraction
     1174  REAL, DIMENSION (nloc, nd), INTENT (OUT)           :: frac_a, frac_s
     1175  REAL, DIMENSION (nloc, nd), INTENT (OUT)           :: qpreca
     1176  REAL, DIMENSION (nloc, nd), INTENT (OUT)           :: qta
    11591177
    11601178!local variables:
    … …  
    11621180  REAL smallestreal
    11631181  REAL tg, qg, dqgdT, ahg, alv, alf, s, tc, es, esi, denom, rg, tca, elacrit
     1182  REAL                                               :: phinu2p
    11641183  REAL als
    11651184  REAL                                               :: qsat_new, snew
    11661185  REAL, DIMENSION (nloc,nd)                          :: qi
     1186  REAL, DIMENSION (nloc,nd)                          :: ha    ! moist static energy of adiabatic ascents
     1187                                                              ! taking into account precip ejection
     1188  REAL, DIMENSION (nloc,nd)                          :: hla   ! liquid water static energy of adiabatic ascents
     1189                                                              ! taking into account precip ejection
     1190  REAL, DIMENSION (nloc,nd)                          :: qcld  ! specific cloud water
    11671191  REAL, DIMENSION (nloc,nd)                          :: qhsat    ! specific humidity at saturation
    11681192  REAL, DIMENSION (nloc,nd)                          :: dqhsatdT ! dqhsat/dT
     1193  REAL, DIMENSION (nloc,nd)                          :: frac  ! ice fraction function of envt temperature
     1194  REAL, DIMENSION (nloc,nd)                          :: qps   ! specific solid precipitation
     1195  REAL, DIMENSION (nloc,nd)                          :: qpl   ! specific liquid precipitation
    11691196  REAL, DIMENSION (nloc)                             :: ah0, cape, capem, byp
    11701197  LOGICAL, DIMENSION (nloc)                          :: lcape
    11711198  INTEGER, DIMENSION (nloc)                          :: iposit
     1199  REAL                                               :: denomm1
    11721200  REAL                                               :: by, defrac, pden, tbis
    11731201  REAL                                               :: fracg
    … …  
    11961224  END DO
    11971225
     1226
    11981227! =====================================================================
    11991228! --- FIND THE REST OF THE LIFTED PARCEL TEMPERATURES
    … …  
    12111240             qnk(i)*(lv0-clmcpv*(tnk(i)-273.15)) + gznk(i)
    12121241  END DO
     1242!
     1243!  Ice fraction
     1244!
     1245  IF (cvflag_ice) THEN
     1246    DO k = minorig, nl
     1247      DO i = 1, ncum
     1248          frac(i, k) = (Tx - t(i,k))/(Tx - Tm)
     1249          frac(i, k) = min(max(frac(i,k),0.0), 1.0)
     1250      END DO
     1251    END DO
     1252! Below cloud base, set ice fraction to cloud base value
     1253    DO k = 1, nl
     1254      DO i = 1, ncum
     1255        IF (k<icb(i)) THEN
     1256          frac(i,k) = frac(i,icb(i))
     1257        END IF
     1258      END DO
     1259    END DO
     1260  ELSE
     1261    DO k = 1, nl
     1262      DO i = 1, ncum
     1263          frac(i,k) = 0.
     1264      END DO
     1265    END DO
     1266  ENDIF ! (cvflag_ice)
     1267
    12131268
    12141269  DO k = minorig, nl
    12151270    DO i = 1,ncum
     1271      ha(i,k) = ah0(i)
     1272      hla(i,k) = hnk(i)
     1273      qta(i,k) = qnk(i)
     1274      qpreca(i,k) = 0.
     1275      frac_a(i,k) = 0.
     1276      frac_s(i,k) = frac(i,k)
     1277      qpl(i,k) = 0.
     1278      qps(i,k) = 0.
    12161279      qhsat(i,k) = qs(i,k)
     1280      qcld(i,k) = max(qta(i,k)-qhsat(i,k),0.)
     1281      IF (k <= icb(i)+1) THEN
     1282        qhsat(i,k) = qnk(i)-clw(i,k)
     1283        qcld(i,k) = clw(i,k)
     1284      ENDIF
    12171285    ENDDO
    12181286  ENDDO
    … …  
    13471415!----------------------------------------------------------------------------
    13481416!
    1349 !  Ice fraction
    1350 !
    1351     DO k = minorig + 1, nl
    1352       DO i = 1, ncum
    1353           frac(i, k) = (Tx - t(i,k))/(Tx - Tm)
    1354           frac(i, k) = min(max(frac(i,k),0.0), 1.0)
    1355       END DO
    1356     END DO
    1357 ! Below cloud base, set ice fraction to cloud base value
    1358     DO k = 1, nl
    1359       DO i = 1, ncum
    1360         IF (k<icb(i)) THEN
    1361           frac(i,k) = frac(i,icb(i))
    1362         END IF
    1363       END DO
    1364     END DO
    1365 
    13661417    DO k = minorig + 1, nl
    13671418      DO i = 1,ncum
    … …  
    13761427          IF (k>=(icbs(i)+1)) THEN                                ! convect3
    13771428            tg = tp(i, k)
    1378             IF (tg .gt. Tx) THEN
     1429            IF (tg .gt. Tx .OR. .NOT.cvflag_ice) THEN
    13791430              es = 6.112*exp(17.67*(tg - 273.15)/(tg + 243.5 - 273.15))
    13801431              qg = eps*es/(p(i,k)-es*(1.-eps))
    … …  
    13861437            ENDIF
    13871438            IF (qsat_depends_on_qt) THEN
    1388               dqgdT = dqgdT*(1.-qnk(i))/(1.-qg)**2
    1389               qg = qg*(1.-qnk(i))/(1.-qg)           
     1439              dqgdT = dqgdT*(1.-qta(i,k-1))/(1.-qg)**2
     1440              qg = qg*(1.-qta(i,k-1))/(1.-qg)           
    13901441            ENDIF
    1391             ahg = (cpd + (cl-cpd)*qnk(i))*tg + lv(i,k)*qg - &
    1392                   lf(i,k)*frac(i,k)*(qnk(i) - qg) + gz(i,k)
     1442            ahg = (cpd + (cl-cpd)*qta(i,k-1))*tg + lv(i,k)*qg - &
     1443                  lf(i,k)*frac(i,k)*(qta(i,k-1) - qg) + gz(i,k)
    13931444            Tp(i,k) = tg + (ah0(i) - ahg)/ &
    1394                     (cpd + (cl-cpd)*qnk(i) + (lv(i,k)+frac(i,k)*lf(i,k))*dqgdT)
     1445                    (cpd + (cl-cpd)*qta(i,k-1) + (lv(i,k)+frac(i,k)*lf(i,k))*dqgdT)
    13951446!!   print *,'undilute2 iterations k, Tp(i,k), ah0(i), ahg ', &
    13961447!!                                 k, Tp(i,k), ah0(i), ahg
    … …  
    14011452        IF (k>=(icbs(i)+1)) THEN                                ! convect3
    14021453          tg = tp(i, k)
    1403           IF (tg .gt. Tx) THEN
     1454          IF (tg .gt. Tx .OR. .NOT.cvflag_ice) THEN
    14041455            es = 6.112*exp(17.67*(tg - 273.15)/(tg + 243.5 - 273.15))
    14051456            qg = eps*es/(p(i,k)-es*(1.-eps))
    … …  
    14081459            qg = eps*esi/(p(i,k)-esi*(1.-eps))
    14091460          ENDIF
    1410             IF (qsat_depends_on_qt) THEN
    1411               qg = qg*(1.-qnk(i))/(1.-qg)           
    1412             ENDIF
    1413           clw(i, k) = qnk(i) - qg
     1461          IF (qsat_depends_on_qt) THEN
     1462            qg = qg*(1.-qta(i,k-1))/(1.-qg)           
     1463          ENDIF
     1464          qhsat(i,k) = qg
     1465        END IF ! (k>=(icbs(i)+1))
     1466      END DO ! i = 1, ncum
     1467      DO i = 1, ncum
     1468        IF (k>=(icbs(i)+1)) THEN                                ! convect3
     1469          clw(i, k) = qta(i,k-1) - qhsat(i,k)
    14141470          clw(i, k) = max(0.0, clw(i,k))
    1415           tvp(i, k) = max(0., tp(i,k)*(1.+qg/eps-qnk(i)))
     1471          tvp(i, k) = max(0., tp(i,k)*(1.+qhsat(i,k)/eps-qta(i,k-1)))
    14161472! print*,tvp(i,k),'tvp'
    14171473          IF (clw(i,k)<1.E-11) THEN
    … …  
    14211477        END IF ! (k>=(icbs(i)+1))
    14221478      END DO ! i = 1, ncum
     1479!
     1480      IF (cvflag_prec_eject) THEN
     1481        DO i = 1, ncum
     1482          IF (k>=(icbs(i)+1)) THEN                                ! convect3
     1483!  Specific precipitation (liquid and solid) and ice content
     1484!  before ejection of precipitation                                                     !!jygprl
     1485            elacrit = elcrit*min(max(1.-(tp(i,k)-T0)/Tlcrit, 0.), 1.)                   !!jygprl
     1486!!!!            qcld(i,k) = min(clw(i,k), elacrit)                                          !!jygprl
     1487            qcld(i,k) = min(clw(i,k), elacrit*(1.-qta(i,k-1))/(1.-elacrit))             !!jygprl
     1488            phinu2p = qhsat(i,k-1) + qcld(i,k-1) - (qhsat(i,k) + qcld(i,k))                     !!jygprl
     1489            qpl(i,k) = qpl(i,k-1) + (1.-frac(i,k))*phinu2p                            !!jygprl
     1490            qps(i,k) = qps(i,k-1) + frac(i,k)     *phinu2p                            !!jygprl
     1491            qi(i,k) = (1.-ejectliq)*clw(i,k)*frac(i,k) + &                            !!jygprl
     1492                     ejectliq*(qps(i,k-1) + frac(i,k)*(phinu2p+qcld(i,k)))            !!jygprl
     1493!!
     1494!  =====================================================================================
     1495!  Ejection of precipitation from adiabatic ascents if requested (cvflag_prec_eject=True):
     1496!  Compute the steps of total water (qta), of moist static energy (ha), of specific
     1497!  precipitation (qpl and qps) and of specific cloud water (qcld) associated with precipitation
     1498!   ejection.
     1499!  =====================================================================================
     1500! 
     1501!   Verif
     1502            qpreca(i,k) = ejectliq*qpl(i,k) + ejectice*qps(i,k)                                   !!jygprl
     1503            frac_a(i,k) = ejectice*qps(i,k)/max(qpreca(i,k),smallestreal)                         !!jygprl
     1504            frac_s(i,k) = (1.-ejectliq)*frac(i,k) + &                                             !!jygprl
     1505               ejectliq*(1. - (qpl(i,k)+(1.-frac(i,k))*qcld(i,k))/max(clw(i,k),smallestreal))     !!jygprl
     1506!         
     1507            denomm1 = 1./(1. - qpreca(i,k))
     1508!         
     1509            qta(i,k) = qta(i,k-1) - &
     1510                      qpreca(i,k)*(1.-qta(i,k-1))*denomm1
     1511            ha(i,k)  = ha(i,k-1) + &
     1512                      ( qpreca(i,k)*(-(1.-qta(i,k-1))*(cl-cpd)*tp(i,k) + &
     1513                                  lv(i,k)*qhsat(i,k) - lf(i,k)*(frac_s(i,k)*qcld(i,k)+qps(i,k))) + &
     1514                        lf(i,k)*ejectice*qps(i,k))*denomm1
     1515            hla(i,k) = hla(i,k-1) + &
     1516                      ( qpreca(i,k)*(-(1.-qta(i,k-1))*(cpv-cpd)*tp(i,k) - &
     1517                                  lv(i,k)*((1.-frac_s(i,k))*qcld(i,k)+qpl(i,k)) - &
     1518                                  (lv(i,k)+lf(i,k))*(frac_s(i,k)*qcld(i,k)+qps(i,k))) + &
     1519                        lv(i,k)*ejectliq*qpl(i,k) + (lv(i,k)+lf(i,k))*ejectice*qps(i,k))*denomm1
     1520            qpl(i,k) = qpl(i,k)*(1.-ejectliq)*denomm1
     1521            qps(i,k) = qps(i,k)*(1.-ejectice)*denomm1
     1522            qcld(i,k) = qcld(i,k)*denomm1
     1523            qhsat(i,k) = qhsat(i,k)*(1.-qta(i,k))/(1.-qta(i,k-1))
     1524         END IF ! (k>=(icbs(i)+1))
     1525        END DO ! i = 1, ncum
     1526      ENDIF  ! (cvflag_prec_eject)
     1527!
    14231528    END DO ! k = minorig + 1, nl
    14241529!
    14251530!----------------------------------------------------------------------------
    14261531!
    1427   ELSE ! (icvflag_Tpa == 2) ELSE IF(icvflag_Tpa == 1)
     1532  ELSE IF (icvflag_Tpa == 0) THEN! (icvflag_Tpa == 2) ELSE IF(icvflag_Tpa == 1)
    14281533!
    14291534!----------------------------------------------------------------------------
    … …  
    15871692!----------------------------------------------------------------------------
    15881693!
    1589   ENDIF ! (icvflag_Tpa == 2) ELSEIF (icvflag_Tpa == 1) ELSE
     1694  ENDIF ! (icvflag_Tpa == 2) ELSEIF (icvflag_Tpa == 1) ELSE (icvflag_Tpa == 0)
    15901695!
    15911696!----------------------------------------------------------------------------
    … …  
    18721977  IF (cvflag_ice) THEN
    18731978!
     1979  IF (cvflag_prec_eject) THEN
     1980!!    DO k = minorig + 1, nl
     1981!!      DO i = 1, ncum
     1982!!        IF ((k>=icb(i)) .AND. (k<=inb(i))) THEN
     1983!!          frac_s(i,k) = qi(i,k)/max(clw(i,k),smallestreal)   
     1984!!          frac_s(i,k) = 1. - (qpl(i,k)+(1.-frac_s(i,k))*qcld(i,k))/max(clw(i,k),smallestreal)   
     1985!!        END IF
     1986!!      END DO
     1987!!    END DO
     1988  ELSE    ! (cvflag_prec_eject)
    18741989    DO k = minorig + 1, nl
    18751990      DO i = 1, ncum
    … …  
    18771992!jyg< frac computation moved to beginning of cv3_undilute2.
    18781993!     kept here for compatibility test with CMip6 version
    1879           frac(i, k) = 1. - (t(i,k)-243.15)/(263.15-243.15)
    1880           frac(i, k) = min(max(frac(i,k),0.0), 1.0)
    1881           hp(i, k) = hnk(i) + (lv(i,k)+(cpd-cpv)*t(i,k)+frac(i,k)*lf(i,k))* &
    1882                               ep(i, k)*clw(i, k)
     1994          frac_s(i, k) = 1. - (t(i,k)-243.15)/(263.15-243.15)
     1995          frac_s(i, k) = min(max(frac_s(i,k),0.0), 1.0)
    18831996        END IF
    18841997      END DO
    18851998    END DO
    1886 ! Below cloud base, set ice fraction to cloud base value
    1887 !!    DO k = 1, nl
    1888 !!      DO i = 1, ncum
    1889 !!        IF (k<icb(i)) THEN
    1890 !!          frac(i,k) = frac(i,icb(i))
    1891 !!        END IF
    1892 !!      END DO
    1893 !!    END DO
     1999  ENDIF  ! (cvflag_prec_eject) ELSE
     2000    DO k = minorig + 1, nl
     2001      DO i = 1, ncum
     2002        IF ((k>=icb(i)) .AND. (k<=inb(i))) THEN
     2003!!          hp(i, k) = hnk(i) + (lv(i,k)+(cpd-cpv)*t(i,k)+frac_s(i,k)*lf(i,k))* &     !!jygprl
     2004!!                              ep(i, k)*clw(i, k)                                    !!jygprl
     2005          hp(i, k) = hla(i,k-1) + (lv(i,k)+(cpd-cpv)*t(i,k)+frac_s(i,k)*lf(i,k))* &   !!jygprl
     2006                              ep(i, k)*clw(i, k)                                      !!jygprl
     2007        END IF
     2008      END DO
     2009    END DO
    18942010!
    18952011  ELSE   ! (cvflag_ice)
    … …  
    25762692SUBROUTINE cv3_unsat(nloc, ncum, nd, na, ntra, icb, inb, iflag, &
    25772693                     t, rr, rs, gz, u, v, tra, p, ph, &
    2578                      th, tv, lv, lf, cpn, ep, sigp, clw, &
     2694                     th, tv, lv, lf, cpn, ep, sigp, clw, frac_s, qpreca, frac_a, qta , &                       !!jygprl
    25792695                     m, ment, elij, delt, plcl, coef_clos, &
    25802696                     mp, rp, up, vp, trap, wt, water, evap, fondue, ice, &
    25812697                     faci, b, sigd, &
    2582                      wdtrainA, wdtrainM)                                      ! RomP
     2698                     wdtrainA, wdtrainS, wdtrainM)                                      ! RomP
    25832699  USE print_control_mod, ONLY: prt_level, lunout
    25842700  IMPLICIT NONE
    … …  
    25982714  REAL, DIMENSION (nloc, na), INTENT (IN)            :: gz
    25992715  REAL, DIMENSION (nloc, nd), INTENT (IN)            :: u, v
    2600   REAL tra(nloc, nd, ntra)
    2601   REAL p(nloc, nd), ph(nloc, nd+1)
    2602   REAL, DIMENSION (nloc, na), INTENT (IN)            :: ep, sigp, clw
     2716  REAL, DIMENSION (nloc, nd, ntra), INTENT(IN)       :: tra
     2717  REAL, DIMENSION (nloc, nd), INTENT (IN)            :: p
     2718  REAL, DIMENSION (nloc, nd+1), INTENT (IN)          :: ph
     2719  REAL, DIMENSION (nloc, na), INTENT (IN)            :: ep, sigp, clw   !adiab ascent shedding
     2720  REAL, DIMENSION (nloc, na), INTENT (IN)            :: frac_s          !ice fraction in adiab ascent shedding !!jygprl
     2721  REAL, DIMENSION (nloc, na), INTENT (IN)            :: qpreca          !adiab ascent precip                   !!jygprl
     2722  REAL, DIMENSION (nloc, na), INTENT (IN)            :: frac_a          !ice fraction in adiab ascent precip   !!jygprl
     2723  REAL, DIMENSION (nloc, na), INTENT (IN)            :: qta             !adiab ascent specific total water     !!jygprl
    26032724  REAL, DIMENSION (nloc, na), INTENT (IN)            :: th, tv, lv, cpn
    26042725  REAL, DIMENSION (nloc, na), INTENT (IN)            :: lf
    … …  
    26132734  REAL, DIMENSION (nloc, na), INTENT (OUT)           :: mp, rp, up, vp
    26142735  REAL, DIMENSION (nloc, na), INTENT (OUT)           :: water, evap, wt
    2615   REAL, DIMENSION (nloc, na), INTENT (OUT)           :: ice, fondue, faci
     2736  REAL, DIMENSION (nloc, na), INTENT (OUT)           :: ice, fondue
     2737  REAL, DIMENSION (nloc, na), INTENT (OUT)           :: faci            ! ice fraction in precipitation
    26162738  REAL, DIMENSION (nloc, na, ntra), INTENT (OUT)     :: trap
    26172739  REAL, DIMENSION (nloc, na), INTENT (OUT)           :: b
    … …  
    26212743! Distinction des wdtrain
    26222744! Pa = wdtrainA     Pm = wdtrainM
    2623   REAL, DIMENSION (nloc, na), INTENT (OUT)           :: wdtrainA, wdtrainM
     2745  REAL, DIMENSION (nloc, na), INTENT (OUT)           :: wdtrainA, wdtrainS, wdtrainM
    26242746
    26252747!local variables
    26262748  INTEGER i, j, k, il, num1, ndp1
     2749  REAL smallestreal
    26272750  REAL tinv, delti, coef
    26282751  REAL awat, afac, afac1, afac2, bfac
    … …  
    26312754  REAL ampmax, thaw
    26322755  REAL tevap(nloc)
    2633   REAL lvcp(nloc, na), lfcp(nloc, na)
    2634   REAL h(nloc, na), hm(nloc, na)
    2635   REAL frac(nloc, na)
    2636   REAL fraci(nloc, na), prec(nloc, na)
     2756  REAL, DIMENSION (nloc, na)      :: lvcp, lfcp
     2757  REAL, DIMENSION (nloc, na)      :: h, hm
     2758  REAL, DIMENSION (nloc, na)      :: ma
     2759  REAL, DIMENSION (nloc, na)      :: frac          ! ice fraction in precipitation source
     2760  REAL, DIMENSION (nloc, na)      :: fraci         ! provisionnal ice fraction in precipitation
     2761  REAL, DIMENSION (nloc, na)      :: prec
    26372762  REAL wdtrain(nloc)
    26382763  LOGICAL lwork(nloc), mplus(nloc)
    … …  
    26412766! ------------------------------------------------------
    26422767IF (prt_level .GE. 10) print *,' ->cv3_unsat, iflag(1) ', iflag(1)
     2768
     2769smallestreal=tiny(smallestreal)
    26432770
    26442771! =============================
    … …  
    26602787!! RomP >>>
    26612788wdtrainA(:,:) = 0.
     2789wdtrainS(:,:) = 0.
    26622790wdtrainM(:,:) = 0.
    26632791!! RomP <<<
    … …  
    27152843  END DO
    27162844
     2845!
     2846! Get adiabatic ascent mass flux
     2847!
     2848!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     2849  IF (adiab_ascent_mass_flux_depends_on_ejectliq) THEN
     2850!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     2851!!! Warning : this option leads to water conservation violation
     2852!!!           Expert only
     2853!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     2854    DO il = 1, ncum
     2855      ma(il, nlp) = 0.
     2856      ma(il, 1)   = 0.
     2857    END DO
     2858
     2859  DO i = nl, 2, -1
     2860      DO il = 1, ncum
     2861        ma(il, i) = ma(il, i+1)*(1.-qta(il,i))/(1.-qta(il,i-1)) + m(il, i)
     2862      END DO
     2863  END DO
     2864!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     2865  ELSE ! (adiab_ascent_mass_flux_depends_on_ejectliq)
     2866!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     2867    DO il = 1, ncum
     2868      ma(il, nlp) = 0.
     2869      ma(il, 1)   = 0.
     2870    END DO
     2871
     2872  DO i = nl, 2, -1
     2873      DO il = 1, ncum
     2874        ma(il, i) = ma(il, i+1) + m(il, i)
     2875      END DO
     2876  END DO
     2877
     2878  ENDIF ! (adiab_ascent_mass_flux_depends_on_ejectliq) ELSE
     2879!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    27172880
    27182881! ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
    … …  
    27392902! ***              calculate detrained precipitation             ***
    27402903
    2741     DO il = 1, ncum
    2742       IF (i<=inb(il) .AND. lwork(il)) THEN
    2743         IF (cvflag_grav) THEN
    2744           wdtrain(il) = grav*ep(il, i)*m(il, i)*clw(il, i)
    2745           wdtrainA(il, i) = wdtrain(il)/grav                        !   Pa   RomP
    2746         ELSE
    2747           wdtrain(il) = 10.0*ep(il, i)*m(il, i)*clw(il, i)
    2748           wdtrainA(il, i) = wdtrain(il)/10.                         !   Pa   RomP
    2749         END IF
    2750       END IF
    2751     END DO
     2904
     2905    DO il = 1, ncum                                                   
     2906      IF (i<=inb(il) .AND. lwork(il)) THEN                           
     2907        wdtrain(il) = grav*ep(il, i)*m(il, i)*clw(il, i)           
     2908        wdtrainS(il, i) = wdtrain(il)/grav                                            !   Ps   jyg
     2909!!        wdtrainA(il, i) = wdtrain(il)/grav                                          !   Ps   RomP
     2910      END IF                                                         
     2911    END DO                                                           
    27522912
    27532913    IF (i>1) THEN
    … …  
    27572917            awat = elij(il, j, i) - (1.-ep(il,i))*clw(il, i)
    27582918            awat = max(awat, 0.0)
    2759             IF (cvflag_grav) THEN
    2760               wdtrain(il) = wdtrain(il) + grav*awat*ment(il, j, i)
    2761               wdtrainM(il, i) = wdtrain(il)/grav - wdtrainA(il, i)  !   Pm  RomP
    2762             ELSE
    2763               wdtrain(il) = wdtrain(il) + 10.0*awat*ment(il, j, i)
    2764               wdtrainM(il, i) = wdtrain(il)/10. - wdtrainA(il, i)   !   Pm  RomP
    2765             END IF
     2919            wdtrain(il) = wdtrain(il) + grav*awat*ment(il, j, i)
     2920            wdtrainM(il, i) = wdtrain(il)/grav - wdtrainS(il, i)    !   Pm  jyg
     2921!!            wdtrainM(il, i) = wdtrain(il)/grav - wdtrainA(il, i)  !   Pm  RomP
    27662922          END IF
    27672923        END DO
    … …  
    27692925    END IF
    27702926
     2927    IF (cvflag_prec_eject) THEN
     2928!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     2929      IF (adiab_ascent_mass_flux_depends_on_ejectliq) THEN
     2930!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     2931!!! Warning : this option leads to water conservation violation
     2932!!!           Expert only
     2933!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     2934          IF ( i > 1) THEN
     2935            DO il = 1, ncum
     2936              IF (i<=inb(il) .AND. lwork(il)) THEN
     2937                wdtrainA(il,i) = ma(il, i+1)*(qta(il, i-1)-qta(il,i))/(1. - qta(il, i-1))    !   Pa   jygprl
     2938                wdtrain(il) = wdtrain(il) + grav*wdtrainA(il,i)
     2939              END IF
     2940            END DO
     2941          ENDIF  ! ( i > 1)
     2942!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     2943      ELSE ! (adiab_ascent_mass_flux_depends_on_ejectliq)
     2944!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     2945          IF ( i > 1) THEN
     2946            DO il = 1, ncum
     2947              IF (i<=inb(il) .AND. lwork(il)) THEN
     2948                wdtrainA(il,i) = ma(il, i+1)*(qta(il, i-1)-qta(il,i))                        !   Pa   jygprl
     2949                wdtrain(il) = wdtrain(il) + grav*wdtrainA(il,i)
     2950              END IF
     2951            END DO
     2952          ENDIF  ! ( i > 1)
     2953
     2954      ENDIF ! (adiab_ascent_mass_flux_depends_on_ejectliq) ELSE
     2955!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     2956    ENDIF  ! (cvflag_prec_eject)
     2957
    27712958
    27722959! ***    find rain water and evaporation using provisional   ***
    … …  
    27742961
    27752962
     2963    IF (cvflag_ice) THEN                                                                                !!jygprl
     2964      DO il = 1, ncum                                                                                   !!jygprl
     2965        IF (i<=inb(il) .AND. lwork(il)) THEN                                                            !!jygprl
     2966          frac(il, i) = (frac_a(il,i)*wdtrainA(il,i)+frac_s(il,i)*(wdtrainS(il,i)+wdtrainM(il,i))) / &  !!jygprl
     2967                        max(wdtrainA(il,i)+wdtrainS(il,i)+wdtrainM(il,i),smallestreal)                  !!jygprl
     2968          fraci(il, i) = frac(il, i)                                                                    !!jygprl
     2969        END IF                                                                                          !!jygprl
     2970      END DO                                                                                            !!jygprl
     2971    END IF                                                                                              !!jygprl
     2972
    27762973    DO il = 1, ncum
    27772974      IF (i<=inb(il) .AND. lwork(il)) THEN
    … …  
    27792976        wt(il, i) = 45.0
    27802977
    2781         IF (cvflag_ice) THEN
    2782           frac(il, inb(il)) = 1. - (t(il,inb(il))-243.15)/(263.15-243.15)
    2783           frac(il, inb(il)) = min(max(frac(il,inb(il)),0.), 1.)
    2784           fraci(il, inb(il)) = frac(il, inb(il))
    2785         ELSE
    2786           CONTINUE
    2787         END IF
    27882978
    27892979        IF (i<inb(il)) THEN
    … …  
    28022992          rp(il, i) = 0.5*(rp(il,i)+rr(il,i))
    28032993        END IF
    2804         fraci(il, i) = 1. - (t(il,i)-243.15)/(263.15-243.15)
    2805         fraci(il, i) = min(max(fraci(il,i),0.0), 1.0)
     2994!!        fraci(il, i) = 1. - (t(il,i)-243.15)/(263.15-243.15)
     2995!!        fraci(il, i) = min(max(fraci(il,i),0.0), 1.0)
    28062996        rp(il, i) = max(rp(il,i), 0.0)
    28072997        rp(il, i) = amin1(rp(il,i), rs(il,i))
    … …  
    32303420                     t, rr, t_wake, rr_wake, s_wake, u, v, tra, &
    32313421                     gz, p, ph, h, hp, lv, lf, cpn, th, th_wake, &
    3232                      ep, clw, m, tp, mp, rp, up, vp, trap, &
     3422                     ep, clw, qpreca, m, tp, mp, rp, up, vp, trap, &
    32333423                     wt, water, ice, evap, fondue, faci, b, sigd, &
    32343424                     ment, qent, hent, iflag_mix, uent, vent, &
    … …  
    32403430!!                     tls, tps,                             ! useless . jyg
    32413431                     qcondc, wd, &
    3242                      ftd, fqd, qnk, qtc, sigt, tau_cld_cv, coefw_cld_cv)
     3432                     ftd, fqd, qta, qtc, sigt, tau_cld_cv, coefw_cld_cv)
    32433433
    32443434    USE print_control_mod, ONLY: lunout, prt_level
    … …  
    32803470      REAL, DIMENSION (nloc, na, na, ntra), INTENT (IN)  :: traent
    32813471      REAL, DIMENSION (nloc, nd), INTENT (IN)            :: tv, tvp, wghti
    3282       REAL,INTENT(IN)                                    :: tau_cld_cv, coefw_cld_cv
     3472      REAL, DIMENSION (nloc, nd), INTENT (IN)            :: qta
     3473      REAL, DIMENSION (nloc, na),INTENT(IN)              :: qpreca
     3474      REAL, INTENT(IN)                                   :: tau_cld_cv, coefw_cld_cv
    32833475!
    32843476!input/output:
    … …  
    33093501      REAL                                               :: ax, bx, cx, dx, ex
    33103502      REAL                                               :: cpinv, rdcp, dpinv
     3503      REAL                                               :: sigaq
    33113504      REAL, DIMENSION (nloc)                             ::  awat
    33123505      REAL, DIMENSION (nloc, nd)                         :: lvcp, lfcp              ! , mke ! unused . jyg
    … …  
    33263519      REAL, DIMENSION (nloc)                             :: sument
    33273520      REAL, DIMENSION (nloc, nd)                         :: sigment, qtment             ! cld
    3328       REAL, DIMENSION (nloc)                             :: qnk
    33293521      REAL sumdq !jyg
    33303522!
    … …  
    34373629  END DO
    34383630
     3631! - Adiabatic ascent mass flux "ma" and cloud base mass flux "cbmf"
     3632!-----------------------------------------------------------------
     3633!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     3634  IF (adiab_ascent_mass_flux_depends_on_ejectliq) THEN
     3635!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     3636!!! Warning : this option leads to water conservation violation
     3637!!!           Expert only
     3638!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     3639  DO il = 1, ncum
     3640    ma(il, nlp) = 0.
     3641    ma(il, 1)   = 0.
     3642  END DO
     3643  DO k = nl, 2, -1
     3644    DO il = 1, ncum
     3645      ma(il, k) = ma(il, k+1)*(1.-qta(il, k))/(1.-qta(il, k-1)) + m(il, k)
     3646      cbmf(il) = max(cbmf(il), ma(il,k))
     3647    END DO
     3648  END DO
     3649  DO k = 2,nl
     3650    DO il = 1, ncum
     3651      IF (k <icb(il)) THEN
     3652        ma(il, k) = ma(il, k-1) + wghti(il,k-1)*cbmf(il)
     3653      ENDIF
     3654    END DO
     3655  END DO
     3656!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     3657  ELSE ! (adiab_ascent_mass_flux_depends_on_ejectliq)
     3658!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     3659!! Line kept for compatibility with earlier versions
    34393660  DO k = 2, nl
    34403661    DO il = 1, ncum
    … …  
    34453666  END DO
    34463667
     3668  DO il = 1, ncum
     3669    ma(il, nlp) = 0.
     3670    ma(il, 1)   = 0.
     3671  END DO
     3672  DO k = nl, 2, -1
     3673    DO il = 1, ncum
     3674      ma(il, k) = ma(il, k+1) + m(il, k)
     3675    END DO
     3676  END DO
     3677  DO k = 2,nl
     3678    DO il = 1, ncum
     3679      IF (k <icb(il)) THEN
     3680        ma(il, k) = ma(il, k-1) + wghti(il,k-1)*cbmf(il)
     3681      ENDIF
     3682    END DO
     3683  END DO
     3684
     3685  ENDIF ! (adiab_ascent_mass_flux_depends_on_ejectliq) ELSE
     3686!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     3687!
    34473688!    print*,'cv3_yield avant ft'
    34483689! am is the part of cbmf taken from the first level
    … …  
    35813822!***    Compute convective mass fluxes upwd and dnwd      ***
    35823823
     3824!
     3825! =================================================
     3826!              upward fluxes                      |
     3827! ------------------------------------------------
     3828!
    35833829upwd(:,:) = 0.
    35843830up_to(:,:) = 0.
    35853831up_from(:,:) = 0.
    3586 dnwd(:,:) = 0.
    3587 dn_to(:,:) = 0.
    3588 dn_from(:,:) = 0.
    3589 !
    3590 ! =================================================
    3591 !              upward fluxes                      |
    3592 ! ------------------------------------------------
     3832!
     3833!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     3834  IF (adiab_ascent_mass_flux_depends_on_ejectliq) THEN
     3835!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     3836!! The decrease of the adiabatic ascent mass flux due to ejection of precipitation
     3837!! is taken into account.
     3838!! WARNING : in the present version, taking into account the mass-flux decrease due to
     3839!! precipitation ejection leads to water conservation violation.
     3840!
     3841! - Upward mass flux of mixed draughts
     3842!---------------------------------------
     3843DO i = 2, nl
     3844  DO j = 1, i-1
     3845    DO il = 1, ncum
     3846      IF (i<=inb(il)) THEN
     3847        up_to(il,i) = up_to(il,i) + ment(il,j,i)
     3848      ENDIF
     3849    ENDDO
     3850  ENDDO
     3851ENDDO
     3852!
     3853DO j = 3, nl
     3854  DO i = 2, j-1
     3855    DO il = 1, ncum
     3856      IF (j<=inb(il)) THEN
     3857        up_from(il,i) = up_from(il,i) + ment(il,i,j)
     3858      ENDIF
     3859    ENDDO
     3860  ENDDO
     3861ENDDO
     3862!
     3863! The difference between upwd(il,i) and upwd(il,i-1) is due to updrafts ending in layer
     3864!(i-1) (theses drafts cross interface (i-1) but not interface(i)) and to updrafts starting
     3865!from layer (i-1) (theses drafts cross interface (i) but not interface(i-1)):
     3866!
     3867DO i = 2, nlp
     3868  DO il = 1, ncum
     3869    IF (i<=inb(il)+1) THEN
     3870      upwd(il,i) = max(0., upwd(il,i-1) - up_to(il,i-1) + up_from(il,i-1))
     3871    ENDIF
     3872  ENDDO
     3873ENDDO
     3874!
     3875! - Total upward mass flux
     3876!---------------------------
     3877DO i = 2, nlp
     3878  DO il = 1, ncum
     3879    IF (i<=inb(il)+1) THEN
     3880      upwd(il,i) = upwd(il,i) + ma(il,i)
     3881    ENDIF
     3882  ENDDO
     3883ENDDO
     3884!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     3885  ELSE ! (adiab_ascent_mass_flux_depends_on_ejectliq)
     3886!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     3887!! The decrease of the adiabatic ascent mass flux due to ejection of precipitation
     3888!! is not taken into account.
     3889!
     3890! - Upward mass flux
     3891!-------------------
    35933892DO i = 2, nl
    35943893  DO il = 1, ncum
    … …  
    36133912  ENDDO
    36143913ENDDO
    3615 !!DO i = 2, nl
    3616 !!  DO j = i+1, nl          !! Permuter les boucles i et j
     3914!
    36173915DO j = 3, nl
    36183916  DO i = 2, j-1
    … …  
    36363934  ENDDO
    36373935ENDDO
     3936
     3937
     3938  ENDIF ! (adiab_ascent_mass_flux_depends_on_ejectliq) ELSE
     3939!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     3940
    36383941!
    36393942! =================================================
    36403943!              downward fluxes                    |
    36413944! ------------------------------------------------
     3945dnwd(:,:) = 0.
     3946dn_to(:,:) = 0.
     3947dn_from(:,:) = 0.
    36423948DO i = 1, nl
    36433949  DO j = i+1, nl
    … …  
    36503956ENDDO
    36513957!
    3652 !!DO i = 2, nl
    3653 !!  DO j = 1, i-1          !! Permuter les boucles i et j
    36543958DO j = 1, nl
    36553959  DO i = j+1, nl
    … …  
    44234727!!!!
    44244728!!!!      ENDDO
     4729
     4730!!  DO i = 1, nlp
     4731!!    DO il = 1, ncum
     4732!!      ma(il, i) = 0
     4733!!    END DO
     4734!!  END DO
     4735!!
     4736!!  DO i = 1, nl
     4737!!    DO j = i, nl
     4738!!      DO il = 1, ncum
     4739!!        ma(il, i) = ma(il, i) + m(il, j)
     4740!!      END DO
     4741!!    END DO
     4742!!  END DO
     4743
     4744!jyg<  (loops stop at nl)
     4745!!  DO i = nl + 1, nd
     4746!!    DO il = 1, ncum
     4747!!      ma(il, i) = 0.
     4748!!    END DO
     4749!!  END DO
     4750!>jyg
     4751
     4752!!  DO i = 1, nl
     4753!!    DO il = 1, ncum
     4754!!      IF (i<=(icb(il)-1)) THEN
     4755!!        ma(il, i) = 0
     4756!!      END IF
     4757!!    END DO
     4758!!  END DO
     4759
    44254760!-----------------------------------------------------------
    44264761        ENDIF !(.NOT.ok_optim_yield)                      !|
    … …  
    44474782!>jyg
    44484783
    4449   DO i = 1, nlp
    4450     DO il = 1, ncum
    4451       ma(il, i) = 0
    4452     END DO
    4453   END DO
    4454 
    4455   DO i = 1, nl
    4456     DO j = i, nl
    4457       DO il = 1, ncum
    4458         ma(il, i) = ma(il, i) + m(il, j)
    4459       END DO
    4460     END DO
    4461   END DO
    4462 
    4463 !jyg<  (loops stop at nl)
    4464 !!  DO i = nl + 1, nd
    4465 !!    DO il = 1, ncum
    4466 !!      ma(il, i) = 0.
    4467 !!    END DO
    4468 !!  END DO
    4469 !>jyg
    4470 
    4471   DO i = 1, nl
    4472     DO il = 1, ncum
    4473       IF (i<=(icb(il)-1)) THEN
    4474         ma(il, i) = 0
    4475       END IF
    4476     END DO
    4477   END DO
    44784784
    44794785! cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
    … …  
    45624868! 14/01/15 AJ delta n'a rien à faire là...                                                 
    45634869    DO il = 1, ncum                                                  ! cld
    4564       IF (wa(il,i)>0.0 .AND. iflag(il)<=1) &                         ! cld
     4870!!      IF (wa(il,i)>0.0 .AND. iflag(il)<=1) &                         ! cld
     4871!!        siga(il, i) = mac(il, i)/(coefw_cld_cv*wa(il, i)) &          ! cld
     4872!!        *rrd*tvp(il, i)/p(il, i)/100.                                ! cld
     4873!!
     4874!!      siga(il, i) = min(siga(il,i), 1.0)                             ! cld
     4875      sigaq = 0.
     4876      IF (wa(il,i)>0.0 .AND. iflag(il)<=1)  THEN                     ! cld
    45654877        siga(il, i) = mac(il, i)/(coefw_cld_cv*wa(il, i)) &          ! cld
    4566         *rrd*tvp(il, i)/p(il, i)/100.                                ! cld
    4567 
    4568       siga(il, i) = min(siga(il,i), 1.0)                             ! cld
     4878                     *rrd*tvp(il, i)/p(il, i)/100.                   ! cld
     4879        siga(il, i) = min(siga(il,i), 1.0)                           ! cld
     4880        sigaq = siga(il,i)*qta(il,i-1)                               ! cld
     4881      ENDIF
    45694882
    45704883! IM cf. FH
    … …  
    45784891        sigment(il,i)=sument(il)*tau_cld_cv/(ph(il,i)-ph(il,i+1))    ! cld
    45794892        sigment(il, i) = min(1.e-4+sigment(il,i), 1.0 - siga(il,i))  ! cld
    4580         qtc(il, i) = (siga(il,i)*qnk(il)+sigment(il,i)*qtment(il,i)) & ! cld
     4893!!        qtc(il, i) = (siga(il,i)*qta(il,i-1)+sigment(il,i)*qtment(il,i)) & ! cld
     4894        qtc(il, i) = (sigaq+sigment(il,i)*qtment(il,i)) & ! cld
    45814895                     /(siga(il,i)+sigment(il,i))                     ! cld
    45824896        sigt(il,i) = sigment(il, i) + siga(il, i)
    45834897
    4584 !        qtc(il, i) = siga(il,i)*qnk(il)+(1.-siga(il,i))*qtment(il,i) ! cld
     4898!        qtc(il, i) = siga(il,i)*qta(il,i-1)+(1.-siga(il,i))*qtment(il,i) ! cld
    45854899!     print*,'BIGAUSSIAN CONV',siga(il,i),sigment(il,i),qtc(il,i) 
    45864900               
    … …  
    48715185      do k=1,nl
    48725186        do i=1,ncum
    4873           hp(i,k)=h(i,k)
    4874         enddo
     5187          hp(i,k)=h(i,k)
     5188        enddo
    48755189      enddo
    48765190

  • LMDZ6/trunk/libf/phylmd/cv3a_uncompress.F90

    r2481 r3496  
    1010                           asupmaxmin, &
    1111                           da, phi, mp, phi2, d1a, dam, sigij, &                ! RomP+AC+jyg
    12                            clw, elij, evap, ep, epmlmMm, eplaMm, &              ! RomP
    13                            wdtrainA, wdtrainM, &                                ! RomP
     12                           qta, clw, elij, evap, ep, epmlmMm, eplaMm, &         ! RomP+jyg
     13                           wdtrainA, wdtrainS, wdtrainM, &                      ! RomP
    1414                           qtc, sigt,          &
    1515                           epmax_diag, & ! epmax_cape
    … …  
    2424                           asupmaxmin1, &
    2525                           da1, phi1, mp1, phi21, d1a1, dam1, sigij1, &         ! RomP+AC+jyg
    26                            clw1, elij1, evap1, ep1, epmlmMm1, eplaMm1, &        ! RomP
    27                            wdtrainA1, wdtrainM1, &                              ! RomP
     26                           qta1, clw1, elij1, evap1, ep1, epmlmMm1, eplaMm1, &  ! RomP+jyg
     27                           wdtrainA1, wdtrainS1, wdtrainM1, &                   ! RomP
    2828                           qtc1, sigt1, &
    2929                           epmax_diag1) ! epmax_cape
    … …  
    7575  REAL, DIMENSION (nloc, nd), INTENT (IN)            :: d1a, dam               !RomP
    7676  REAL, DIMENSION (nloc, nd, nd), INTENT (IN)        :: sigij                  !RomP
     77  REAL, DIMENSION (nloc, nd), INTENT (IN)            :: qta                    !jyg
    7778  REAL, DIMENSION (nloc, nd), INTENT (IN)            :: clw                    !RomP
    7879  REAL, DIMENSION (nloc, nd, nd), INTENT (IN)        :: elij                   !RomP
    … …  
    8182  REAL, DIMENSION (nloc, nd), INTENT (IN)            :: eplamM                 !RomP+jyg
    8283  REAL, DIMENSION (nloc, nd), INTENT (IN)            :: qtc, sigt              !RomP
    83   REAL, DIMENSION (nloc, nd), INTENT (IN)            :: wdtrainA, wdtrainM     !RomP
     84  REAL, DIMENSION (nloc, nd), INTENT (IN)            :: wdtrainA, wdtrainS, wdtrainM     !RomP
    8485
    8586  ! outputs:
    … …  
    111112  REAL, DIMENSION (len, nd), INTENT (OUT)            :: d1a1, dam1 !RomP       !RomP
    112113  REAL, DIMENSION (len, nd, nd), INTENT (OUT)        :: sigij1                 !RomP
     114  REAL, DIMENSION (len, nd), INTENT (OUT)            :: qta1                   !jyg
    113115  REAL, DIMENSION (len, nd), INTENT (OUT)            :: clw1                   !RomP
    114116  REAL, DIMENSION (len, nd, nd), INTENT (OUT)        :: elij1                  !RomP
    … …  
    117119  REAL, DIMENSION (len, nd), INTENT (OUT)            :: eplamM1                !RomP+jyg
    118120  REAL, DIMENSION (len, nd), INTENT (OUT)            :: qtc1, sigt1            !RomP
    119   REAL, DIMENSION (len, nd), INTENT (OUT)            :: wdtrainA1, wdtrainM1   !RomP
     121  REAL, DIMENSION (len, nd), INTENT (OUT)            :: wdtrainA1, wdtrainS1, wdtrainM1   !RomP
    120122
    121123
    … …  
    175177        d1a1(idcum(i), k) = d1a(i, k) !RomP
    176178        dam1(idcum(i), k) = dam(i, k) !RomP
     179        qta1(idcum(i), k) = qta(i, k) !jyg
    177180        clw1(idcum(i), k) = clw(i, k) !RomP
    178181        evap1(idcum(i), k) = evap(i, k) !RomP
    … …  
    180183        eplamM1(idcum(i), k) = eplamM(i, k) !RomP+jyg
    181184        wdtrainA1(idcum(i), k) = wdtrainA(i, k) !RomP
     185        wdtrainS1(idcum(i), k) = wdtrainS(i, k) !RomP
    182186        wdtrainM1(idcum(i), k) = wdtrainM(i, k) !RomP
    183187        qtc1(idcum(i), k) = qtc(i, k)
    … …  
    300304      d1a1(:, 1:nl) = d1a(:, 1:nl)            !RomP
    301305      dam1(:, 1:nl) = dam(:, 1:nl)            !RomP
     306      qta1(:, 1:nl) = qta(:, 1:nl)            !jyg
    302307      clw1(:, 1:nl) = clw(:, 1:nl)            !RomP
    303308      evap1(:, 1:nl) = evap(:, 1:nl)          !RomP
    … …  
    305310      eplamM1(:, 1:nl) = eplamM(:, 1:nl)       !RomP+jyg
    306311      wdtrainA1(:, 1:nl) = wdtrainA(:, 1:nl)  !RomP
     312      wdtrainS1(:, 1:nl) = wdtrainS(:, 1:nl)  !RomP
    307313      wdtrainM1(:, 1:nl) = wdtrainM(:, 1:nl)  !RomP
    308314      qtc1(:, 1:nl) = qtc(:, 1:nl)

  • LMDZ6/trunk/libf/phylmd/cv3p2_closure.F90

    r2502 r3496  
    2323  include "cvthermo.h"
    2424  include "cv3param.h"
     25  include "cvflag.h"
    2526  include "YOMCST2.h"
    2627  include "YOMCST.h"

  • LMDZ6/trunk/libf/phylmd/cv3p_mixing.F90

    r2905 r3496  
    11SUBROUTINE cv3p_mixing(nloc, ncum, nd, na, ntra, icb, nk, inb, &
    2                        ph, t, rr, rs, u, v, tra, h, lv, lf, frac, qnk, &
     2                       ph, t, rr, rs, u, v, tra, h, lv, lf, frac, qta, &
    33                       unk, vnk, hp, tv, tvp, ep, clw, sig, &
    44                       Ment, Qent, hent, uent, vent, nent, &
    … …  
    2929  INTEGER, DIMENSION (nloc), INTENT (IN)             :: icb, inb, nk
    3030  REAL, DIMENSION (nloc, nd), INTENT (IN)            :: sig
    31   REAL, DIMENSION (nloc), INTENT (IN)                :: qnk, unk, vnk
     31  REAL, DIMENSION (nloc), INTENT (IN)                :: unk, vnk
     32  REAL, DIMENSION (nloc, nd), INTENT (IN)            :: qta
    3233  REAL, DIMENSION (nloc, nd+1), INTENT (IN)          :: ph
    3334  REAL, DIMENSION (nloc, nd), INTENT (IN)            :: t, rr, rs
    … …  
    173174                           .AND. (j<=inb(il))) THEN
    174175
    175             rti = qnk(il) - ep(il, i)*clw(il, i)
     176!!            rti = qnk(il) - ep(il, i)*clw(il, i)
     177            rti = qta(il,i-1) - ep(il, i)*clw(il, i)
    176178            bf2 = 1. + lv(il, j)*lv(il, j)*rs(il, j)/(rrv*t(il,j)*t(il,j)*cpd)
    177179!jyg(from aj)<
    … …  
    219221            Sij(il, i, j) = amax1(0.0, Sij(il,i,j))
    220222            Sij(il, i, j) = amin1(1.0, Sij(il,i,j))
     223          ELSE IF (j > i) THEN
     224            IF (prt_level >= 10) THEN
     225              print *,'cv3p_mixing i, j, Sij given by the no-precip eq. ', i, j, Sij(il,i,j)
     226            ENDIF
    221227          END IF ! new
    222228        END DO
    … …  
    248254!!!       Ment(il,i,i)=m(il,i)
    249255        Ment(il, i, i) = 1.
    250         Qent(il, i, i) = qnk(il) - ep(il, i)*clw(il, i)
     256!!        Qent(il, i, i) = qnk(il) - ep(il, i)*clw(il, i)
     257        Qent(il, i, i) = qta(il,i-1) - ep(il, i)*clw(il, i)
    251258        uent(il, i, i) = unk(il)
    252259        vent(il, i, i) = vnk(il)
    … …  
    332339      IF (i>=icb(il) .AND. i<=inb(il)) THEN
    333340        lwork(il) = (nent(il,i)/=0)
    334         rti = qnk(il) - ep(il, i)*clw(il, i)
     341!!        rti = qnk(il) - ep(il, i)*clw(il, i)
     342        rti = qta(il,i-1) - ep(il, i)*clw(il, i)
    335343!jyg<
    336344        IF (cvflag_ice) THEN
    … …  
    462470            lwork(il)) THEN
    463471          IF (Sij(il,i,j)>0.0) THEN
    464             rti = qnk(il) - ep(il, i)*clw(il, i)
     472!!            rti = qnk(il) - ep(il, i)*clw(il, i)
     473            rti = qta(il,i-1) - ep(il, i)*clw(il, i)
    465474            Qmixmax(il) = Qmix(Sjmax(il))
    466475            Qmixmin(il) = Qmix(Sjmin(il))
    … …  
    590599              lwork(il)) THEN
    591600            IF (Sij(il,i,j)>0.0) THEN
    592               rti = qnk(il) - ep(il, i)*clw(il, i)
     601!!              rti = qnk(il) - ep(il, i)*clw(il, i)
     602              rti = qta(il,i-1) - ep(il, i)*clw(il, i)
    593603!!!             Ment(il,i,i) = m(il,i)*abs(Qmixmax(il)*(1.-Sjmax(il))
    594604              Ment(il, i, i) = abs(Qmixmax(il)*(1.-Sjmax(il))+Rmixmax(il) - &
    … …  
    659669! cc        Ment(il,i,i)=m(il,i)
    660670        Ment(il, i, i) = 1.
    661         Qent(il, i, i) = qnk(il) - ep(il, i)*clw(il, i)
     671!!        Qent(il, i, i) = qnk(il) - ep(il, i)*clw(il, i)
     672        Qent(il, i, i) = qta(il,i-1) - ep(il, i)*clw(il, i)
    662673        uent(il, i, i) = unk(il)
    663674        vent(il, i, i) = vnk(il)

  • LMDZ6/trunk/libf/phylmd/cv3param.h

    r3492 r3496  
    77!------------------------------------------------------------
    88
    9       logical ok_homo_tend
    10       logical ok_optim_yield
    11       logical ok_entrain
    12       logical ok_convstop
    13       logical ok_intermittent
    14       logical qsat_depends_on_qt
     9      integer flag_epKEorig, flag_wb
     10      integer cv_flag_feed
    1511      integer noff, minorig, nl, nlp, nlm
    16       integer cv_flag_feed
    17       integer flag_epKEorig,flag_wb
    1812      real sigdz, spfac
    1913      real pbcrit, ptcrit
    … …  
    2822      real delta
    2923      real betad
     24      real ejectliq
     25      real ejectice
    3026
    3127      COMMON /cv3param/ sigdz, spfac &
    … …  
    4036                      ,wbmax &
    4137                      ,delta, betad  &
     38                      ,ejectliq, ejectice &
    4239                      ,flag_epKEorig &
    4340                      ,flag_wb, cv_flag_feed &
    44                       ,noff, minorig, nl, nlp, nlm  &
    45                       ,ok_convstop, ok_intermittent &
    46                       ,ok_optim_yield &
    47                       ,ok_entrain &
    48                       ,ok_homo_tend &
    49                       ,qsat_depends_on_qt
     41                      ,noff, minorig, nl, nlp, nlm
    5042!$OMP THREADPRIVATE(/cv3param/)
    5143

  • LMDZ6/trunk/libf/phylmd/cva_driver.F90

    r3435 r3496  
    2525!!                      elij1,evap1,ep1,epmlmMm1,eplaMm1, &                ! RomP
    2626                      da1, phi1, mp1, phi21, d1a1, dam1, sigij1, wghti1, & ! RomP, RL
    27                       clw1, elij1, evap1, ep1, epmlmMm1, eplaMm1, &        ! RomP, RL
    28                       wdtrainA1, wdtrainM1, qtc1, sigt1, tau_cld_cv, &
     27                      qta1, clw1, elij1, evap1, ep1, epmlmMm1, eplaMm1, &  ! RomP, RL
     28                      wdtrainA1, wdtrainS1, wdtrainM1, qtc1, sigt1, tau_cld_cv, &     !!jygprl
    2929                      coefw_cld_cv, &                                      ! RomP, AJ
    3030                      epmax_diag1)  ! epmax_cape
    … …  
    124124!                                      of dimension ND, defined at same grid levels as T, Q, QS and P.
    125125
    126 ! wdtrainA1     Real           Output   precipitation detrained from adiabatic draught;
     126! wdtrainA1     Real           Output   precipitation ejected from adiabatic draught;
     127!                                         should be used in tracer transport (cvltr)
     128! wdtrainS1     Real           Output   precipitation detrained from shedding of adiabatic draught;
    127129!                                         used in tracer transport (cvltr)
    128130! wdtrainM1     Real           Output   precipitation detrained from mixed draughts;
    … …  
    248250
    249251! RomP >>>
    250   REAL, DIMENSION (len, nd), INTENT (OUT)            :: wdtrainA1, wdtrainM1 ! precipitation sources (extensive)
     252  REAL, DIMENSION (len, nd), INTENT (OUT)            :: wdtrainA1, wdtrainS1, wdtrainM1 ! precipitation sources (extensive)
    251253  REAL, DIMENSION (len, nd), INTENT (OUT)            :: mp1  ! unsat. mass flux (staggered grid)
    252254  REAL, DIMENSION (len, nd), INTENT (OUT)            :: da1  ! detrained mass flux of adiab. asc. air (extensive)
    … …  
    258260  REAL, DIMENSION (len, nd, nd), INTENT (OUT)        :: sigij1 ! mass fraction of env. air in mixed draughts (intensive)
    259261  REAL, DIMENSION (len, nd, nd), INTENT (OUT)        :: elij1! cond. water per unit mass of mixed draughts (intensive)
     262  REAL, DIMENSION (len, nd), INTENT (OUT)            :: qta1 ! total water per unit mass of the adiab. asc. (intensive)
    260263  REAL, DIMENSION (len, nd), INTENT (OUT)            :: clw1 ! cond. water per unit mass of the adiab. asc. (intensive)
    261264!JYG,RL
    … …  
    467470  REAL tv_wake(nloc, nd)
    468471  REAL clw(nloc, nd)
     472  REAL, DIMENSION(nloc, nd)    :: qta, qpreca                       !!jygprl
    469473  REAL dph(nloc, nd)
    470474  REAL pbase(nloc), buoybase(nloc), th(nloc, nd)
    … …  
    477481  REAL cin(nloc)
    478482  REAL m(nloc, nd)
     483  REAL mm(nloc, nd)
    479484  REAL ment(nloc, nd, nd), sigij(nloc, nd, nd)
    480485  REAL qent(nloc, nd, nd)
    … …  
    494499  REAL, DIMENSION(len,nd)     :: wt, water, evap
    495500  REAL, DIMENSION(len,nd)     :: ice, fondue, b
    496   REAL, DIMENSION(len,nd)     :: frac, faci
     501  REAL, DIMENSION(len,nd)     :: frac_a, frac_s, faci               !!jygprl
    497502  REAL ft(nloc, nd), fq(nloc, nd)
    498503  REAL ftd(nloc, nd), fqd(nloc, nd)
    … …  
    523528 
    524529! RomP >>>
    525   REAL wdtrainA(nloc, nd), wdtrainM(nloc, nd)
     530  REAL wdtrainA(nloc, nd), wdtrainS(nloc, nd), wdtrainM(nloc, nd)   !!jygprl
    526531  REAL da(len, nd), phi(len, nd, nd)
    527532  REAL epmlmMm(nloc, nd, nd), eplaMm(nloc, nd)
    … …  
    615620  tvp(:, :) = 0. !ym missing init, need to have a look by developpers
    616621  tv(:, :) = 0. !ym missing init, need to have a look by developpers
    617    
     622
    618623  DO il = 1, len
    619624    cin1(il) = -100000.
    … …  
    636641  qtc1(:, :) = 0.
    637642  wdtrainA1(:, :) = 0.
     643  wdtrainS1(:, :) = 0.
    638644  wdtrainM1(:, :) = 0.
    639645  da1(:, :) = 0.
    … …  
    646652  sigij1(:, :, :) = 0.
    647653  elij1(:, :, :) = 0.
     654  qta1(:,:) = 0.
    648655  clw1(:,:) = 0.
    649656  wghti1(:,:) = 0.
    … …  
    906913                         p, ph, h, tv, lv, lf, pbase, buoybase, plcl, &
    907914                         inb, tp, tvp, clw, hp, ep, sigp, buoy, &
    908                          frac)
     915                         frac_a, frac_s, qpreca, qta)                        !!jygprl
    909916    END IF
    910917
    … …  
    915922                        tnk, qnk, gznk, t, q, qs, gz, &
    916923                        p, dph, h, tv, lv, &
    917                         inb, inbis, tp, tvp, clw, hp, ep, sigp, frac)
     924                        inb, inbis, tp, tvp, clw, hp, ep, sigp, frac_s)
    918925    END IF
    919926
    … …  
    923930             PRINT *, 'cva_driver -> cv3_epmax_cape'
    924931    call cv3_epmax_fn_cape(nloc,ncum,nd &
    925                 , ep,hp,icb,inb,clw,nk,t,h,hnk,lv,lf,frac &
     932                , ep,hp,icb,inb,clw,nk,t,h,hnk,lv,lf,frac_s &
    926933                , pbase, p, ph, tv, buoy, sig, w0,iflag &
    927934                , epmax_diag)
    … …  
    941948             PRINT *, 'cva_driver -> cv3p_mixing'
    942949        CALL cv3p_mixing(nloc, ncum, nd, nd, ntra, icb, nk, inb, &           ! na->nd
    943                          ph, t, q, qs, u, v, tra, h, lv, lf, frac, qnk, &
     950!!                         ph, t, q, qs, u, v, tra, h, lv, lf, frac, qnk, &
     951                         ph, t, q, qs, u, v, tra, h, lv, lf, frac_s, qta, &      !!jygprl
    944952                         unk, vnk, hp, tv, tvp, ep, clw, sig, &
    945953                         ment, qent, hent, uent, vent, nent, &
    … …  
    10211029             PRINT *, 'cva_driver -> cv3_mixing'
    10221030        CALL cv3_mixing(nloc, ncum, nd, nd, ntra, icb, nk, inb, &             ! na->nd
    1023                         ph, t, q, qs, u, v, tra, h, lv, lf, frac, qnk, &
     1031                        ph, t, q, qs, u, v, tra, h, lv, lf, frac_s, qnk, &
    10241032                        unk, vnk, hp, tv, tvp, ep, clw, m, sig, &
    10251033                        ment, qent, uent, vent, nent, sigij, elij, ments, qents, traent)
    10261034        CALL zilch(hent, nloc*nd*nd)
    10271035      ELSE
    1028         CALL cv3_mixscale(nloc, ncum, nd, ment, m)
     1036!!jyg:  Essais absurde pour voir
     1037!!        mm(:,1) = 0.
     1038!!        DO  i = 2,nd
     1039!!          mm(:,i) = m(:,i)*(1.-qta(:,i-1))
     1040!!        ENDDO
     1041        mm(:,:) = m(:,:)
     1042        CALL cv3_mixscale(nloc, ncum, nd, ment, mm)
    10291043        IF (debut) THEN
    10301044          PRINT *, ' cv3_mixscale-> '
    … …  
    10621076                     t_wake, q_wake, qs_wake, gz, u, v, tra, p, ph, &
    10631077                     th_wake, tv_wake, lv_wake, lf_wake, cpn_wake, &
    1064                      ep, sigp, clw, &
     1078                     ep, sigp, clw, frac_s, qpreca, frac_a, qta, &                    !!jygprl
    10651079                     m, ment, elij, delt, plcl, coef_clos, &
    10661080                     mp, qp, up, vp, trap, wt, water, evap, fondue, ice, &
    10671081                     faci, b, sigd, &
    1068                      wdtrainA, wdtrainM)                                       ! RomP
     1082!!                     wdtrainA, wdtrainM)                                       ! RomP
     1083                     wdtrainA, wdtrainS, wdtrainM)                               !!jygprl
    10691084!
    10701085      IF (prt_level >= 10) THEN
    … …  
    10751090           evap(igout,k), fondue(igout,k)
    10761091        ENDDO
    1077         Print *, 'cva_driver after cv3_unsat: wdtrainA, wdtrainM '
     1092        Print *, 'cva_driver after cv3_unsat: wdtrainA, wdtrainS, wdtrainM '     !!jygprl
    10781093        DO k = 1,nd
    1079         write (6, '(i4,2(1x,e13.6))'), &
    1080            k, wdtrainA(igout,k), wdtrainM(igout,k)
     1094        write (6, '(i4,3(1x,e13.6))'), &
     1095           k, wdtrainA(igout,k), wdtrainS(igout,k), wdtrainM(igout,k)            !!jygprl
    10811096        ENDDO
    10821097      ENDIF
    … …  
    11121127                     t, q, t_wake, q_wake, s_wake, u, v, tra, &
    11131128                     gz, p, ph, h, hp, lv, lf, cpn, th, th_wake, &
    1114                      ep, clw, m, tp, mp, qp, up, vp, trap, &
     1129                     ep, clw, qpreca, m, tp, mp, qp, up, vp, trap, &
    11151130                     wt, water, ice, evap, fondue, faci, b, sigd, &
    11161131                     ment, qent, hent, iflag_mix, uent, vent, &
    … …  
    11211136!!                     tls, tps, &                            ! useless . jyg
    11221137                     qcondc, wd, &
    1123                      ftd, fqd, qnk, qtc, sigt, tau_cld_cv, coefw_cld_cv)
     1138!!                     ftd, fqd, qnk, qtc, sigt, tau_cld_cv, coefw_cld_cv)
     1139                     ftd, fqd, qta, qtc, sigt, tau_cld_cv, coefw_cld_cv)         !!jygprl
     1140!
     1141!         Test conseravtion de l'eau
    11241142!
    11251143      IF (debut) THEN
    … …  
    11421160                     t, q, u, v, &
    11431161                     gz, p, ph, h, hp, lv, cpn, &
    1144                      ep, clw, frac, m, mp, qp, up, vp, &
     1162                     ep, clw, frac_s, m, mp, qp, up, vp, &
    11451163                     wt, water, evap, &
    11461164                     ment, qent, uent, vent, nent, elij, &
    … …  
    11871205                           asupmaxmin, &
    11881206                           da, phi, mp, phi2, d1a, dam, sigij, &         ! RomP
    1189                            clw, elij, evap, ep, epmlmMm, eplaMm, &       ! RomP
    1190                            wdtrainA, wdtrainM, &                         ! RomP
     1207                           qta, clw, elij, evap, ep, epmlmMm, eplaMm, &  ! RomP
     1208                           wdtrainA, wdtrainS, wdtrainM, &                         ! RomP
    11911209                           qtc, sigt, epmax_diag, & ! epmax_cape
    11921210                           iflag1, kbas1, ktop1, &
    … …  
    11991217                           Plim11, plim21, asupmax1, supmax01, &
    12001218                           asupmaxmin1, &
    1201                            da1, phi1, mp1, phi21, d1a1, dam1, sigij1,  & ! RomP
    1202                            clw1, elij1, evap1, ep1, epmlmMm1, eplaMm1, & ! RomP
    1203                            wdtrainA1, wdtrainM1,                       & ! RomP
     1219                           da1, phi1, mp1, phi21, d1a1, dam1, sigij1,  &       ! RomP
     1220                           qta1, clw1, elij1, evap1, ep1, epmlmMm1, eplaMm1, & ! RomP
     1221                           wdtrainA1, wdtrainS1, wdtrainM1,                  & ! RomP
    12041222                           qtc1, sigt1, epmax_diag1) ! epmax_cape
    12051223!   

  • LMDZ6/trunk/libf/phylmd/cvflag.h

    r3492 r3496  
    44      logical cvflag_grav
    55      logical cvflag_ice
     6      logical ok_optim_yield
     7      logical ok_entrain
     8      logical ok_homo_tend
     9      logical ok_convstop
     10      logical ok_intermittent
     11      logical cvflag_prec_eject
     12      logical qsat_depends_on_qt
     13      logical adiab_ascent_mass_flux_depends_on_ejectliq
    614      integer icvflag_Tpa
    715
    8       COMMON /cvflag/ icvflag_Tpa, cvflag_grav, cvflag_ice
     16      COMMON /cvflag/ icvflag_Tpa, &
     17                      cvflag_grav, cvflag_ice, &
     18                      ok_optim_yield, &
     19                      ok_entrain, &
     20                      ok_homo_tend, &
     21                      ok_convstop, ok_intermittent, &
     22                      cvflag_prec_eject, &
     23                      qsat_depends_on_qt, &
     24                      adiab_ascent_mass_flux_depends_on_ejectliq
    925!$OMP THREADPRIVATE(/cvflag/)

  • LMDZ6/trunk/libf/phylmd/phys_local_var_mod.F90

    r3489 r3496  
    416416      REAL,ALLOCATABLE,SAVE,DIMENSION(:,:)  :: pmflxr, pmflxs
    417417!$OMP THREADPRIVATE(pmflxr, pmflxs)
    418       REAL,ALLOCATABLE,SAVE,DIMENSION(:,:)  :: wdtrainA, wdtrainM
    419 !$OMP THREADPRIVATE(wdtrainA, wdtrainM)
     418      REAL,ALLOCATABLE,SAVE,DIMENSION(:,:)  :: wdtrainA, wdtrainS, wdtrainM
     419!$OMP THREADPRIVATE(wdtrainA, wdtrainS, wdtrainM)
    420420      REAL,ALLOCATABLE,SAVE,DIMENSION(:,:)  :: da, mp
    421421!$OMP THREADPRIVATE(da, mp)
    … …  
    428428      REAL,ALLOCATABLE,SAVE,DIMENSION(:,:)  :: ev
    429429!$OMP THREADPRIVATE(ev)
     430      REAL,ALLOCATABLE,SAVE,DIMENSION(:,:)  :: qtaa
     431!$OMP THREADPRIVATE(qtaa)
    430432      REAL,ALLOCATABLE,SAVE,DIMENSION(:,:)  :: clw
    431433!$OMP THREADPRIVATE(clw)
    … …  
    791793!  Deep convective variables used in phytrac
    792794      ALLOCATE(pmflxr(klon, klev+1), pmflxs(klon, klev+1))
    793       ALLOCATE(wdtrainA(klon,klev),wdtrainM(klon,klev))
     795      ALLOCATE(wdtrainA(klon,klev),wdtrainS(klon,klev),wdtrainM(klon,klev))
    794796      ALLOCATE(dnwd(klon, klev), upwd(klon, klev) )
    795797      ALLOCATE(ep(klon,klev))                          ! epmax_cape
    … …  
    801803      ALLOCATE(ev(klon,klev) )
    802804      ALLOCATE(elij(klon,klev,klev) )
     805      ALLOCATE(qtaa(klon,klev) )
    803806      ALLOCATE(clw(klon,klev) )
    804807      ALLOCATE(epmlmMm(klon,klev,klev), eplaMm(klon,klev) )
    … …  
    10821085
    10831086      DEALLOCATE(pmflxr, pmflxs)
    1084       DEALLOCATE(wdtrainA, wdtrainM)
     1087      DEALLOCATE(wdtrainA, wdtrainS, wdtrainM)
    10851088      DEALLOCATE(upwd, dnwd)
    10861089      DEALLOCATE(ep)
    … …  
    10921095      DEALLOCATE(ev )
    10931096      DEALLOCATE(elij )
     1097      DEALLOCATE(qtaa )
    10941098      DEALLOCATE(clw )
    10951099      DEALLOCATE(epmlmMm, eplaMm )

  • LMDZ6/trunk/libf/phylmd/phys_output_ctrlout_mod.F90

    r3486 r3496  
    15261526  TYPE(ctrl_out), SAVE :: o_clwcon = ctrl_out((/ 4, 10, 10, 10, 10, 10, 11, 11, 11, 11/), &
    15271527    'clwcon', 'Convective Cloud Liquid water content', 'kg/kg', (/ ('', i=1, 10) /))
     1528  TYPE(ctrl_out), SAVE :: o_Mipsh = ctrl_out((/ 4, 10, 10, 10, 10, 10, 11, 11, 11, 11/), &
     1529    'Mipsh', 'mass flux shed from adiab. ascents', 'kg/m2/s', (/ ('', i=1, 10) /))
    15281530  TYPE(ctrl_out), SAVE :: o_Ma = ctrl_out((/ 4, 10, 10, 10, 10, 10, 11, 11, 11, 11/), &
    1529     'Ma', 'undilute adiab updraft', 'kg/m2/s', (/ ('', i=1, 10) /))
     1531    'Ma', 'undilute adiab updraft mass flux', 'kg/m2/s', (/ ('', i=1, 10) /))
    15301532  TYPE(ctrl_out), SAVE :: o_dnwd = ctrl_out((/ 4, 10, 10, 10, 10, 10, 11, 11, 11, 11/), &
    15311533    'dnwd', 'saturated downdraft', 'kg/m2/s', (/ ('', i=1, 10) /))
    … …  
    15911593  TYPE(ctrl_out), SAVE :: o_wdtrainA = ctrl_out((/ 4, 5, 10,  4, 10, 10, 11, 11, 11, 11 /), &
    15921594    'wdtrainA', 'precipitation from AA', '-', (/ ('', i=1, 10) /))
     1595  TYPE(ctrl_out), SAVE :: o_wdtrainS = ctrl_out((/ 4, 5, 10,  4, 10, 10, 11, 11, 11, 11 /), &
     1596    'wdtrainS', 'precipitation from shedding of AA', '-', (/ ('', i=1, 10) /))
    15931597  TYPE(ctrl_out), SAVE :: o_wdtrainM = ctrl_out((/ 4, 5, 10,  4, 10, 10, 11, 11, 11, 11 /), &
    15941598    'wdtrainM', 'precipitation from mixture', '-', (/ ('', i=1, 10) /))
    15951599  TYPE(ctrl_out), SAVE :: o_Vprecip = ctrl_out((/ 10, 10, 10, 10, 10, 10, 11, 11, 11, 11/), &
    15961600    'Vprecip', 'precipitation vertical profile', '-', (/ ('', i=1, 10) /))
     1601  TYPE(ctrl_out), SAVE :: o_qtaa = ctrl_out((/ 10, 10, 10, 10, 10, 10, 11, 11, 11, 11/), &
     1602    'qtaa', 'specific total water in adiabatic ascents', 'kg/kg', (/ ('', i=1, 10) /))
     1603  TYPE(ctrl_out), SAVE :: o_clwaa = ctrl_out((/ 10, 10, 10, 10, 10, 10, 11, 11, 11, 11/), &
     1604    'Clwaa', 'specific condensed water in adiabatic ascents', 'kg/kg', (/ ('', i=1, 10) /))
    15971605  TYPE(ctrl_out), SAVE :: o_ftd = ctrl_out((/ 4, 5, 10, 10, 10, 10, 11, 11, 11, 11/), &
    15981606    'ftd', 'tend temp due aux descentes precip', '-', (/ ('', i=1, 10) /))

  • LMDZ6/trunk/libf/phylmd/phys_output_write_mod.F90

    r3489 r3496  
    7373         o_uwat, o_vwat, &
    7474         o_ptop, o_fbase, o_plcl, o_plfc, &
    75          o_wbeff, o_convoccur, o_cape_max, o_upwd, o_ep,o_epmax_diag, o_Ma, &
     75         o_wbeff, o_convoccur, o_cape_max, o_upwd, o_ep,o_epmax_diag, &
     76         o_Mipsh, o_Ma, &
    7677         o_dnwd, o_dnwd0, o_ftime_deepcv, o_ftime_con, o_mc, &
    7778         o_prw, o_prlw, o_prsw, o_s_pblh, o_s_pblt, o_s_lcl, &
    … …  
    8889         o_wake_s, o_wake_deltat, o_wake_deltaq, &
    8990         o_wake_omg, o_dtwak, o_dqwak, o_dqwak2d, o_Vprecip, &
    90          o_ftd, o_fqd, o_wdtrainA, o_wdtrainM, &
     91         o_qtaa, o_Clwaa, &
     92         o_ftd, o_fqd, o_wdtrainA, o_wdtrainS, o_wdtrainM, &
    9193         o_n2, o_s2, o_proba_notrig, &
    9294         o_random_notrig, o_ale_bl_stat, &
    … …  
    232234         delta_tsurf, &
    233235         wstar, cape, ema_pcb, ema_pct, &
    234          ema_cbmf, Ma, fm_therm, ale_bl, alp_bl, ale, &
     236         ema_cbmf, Mipsh, Ma, fm_therm, ale_bl, alp_bl, ale, &
    235237         alp, cin, wake_pe, wake_dens, wake_s, wake_deltat, &
    236238         wake_deltaq, ftd, fqd, ale_bl_trig, albsol1, &
    … …  
    273275         kh         ,kh_x       ,kh_w       , &
    274276         cv_gen, wake_h, &
    275          wake_omg, d_t_wake, d_q_wake, Vprecip, &
    276          wdtrainA, wdtrainM, n2, s2, proba_notrig, &
     277         wake_omg, d_t_wake, d_q_wake, Vprecip, qtaa, Clw, &
     278         wdtrainA, wdtrainS, wdtrainM, n2, s2, proba_notrig, &
    277279         random_notrig, &
    278280         alp_bl_det, alp_bl_fluct_m, alp_bl_conv, &
    … …  
    10581060
    10591061          CALL histwrite_phy(o_upwd, upwd)
     1062          CALL histwrite_phy(o_Mipsh, Mipsh)
    10601063          CALL histwrite_phy(o_Ma, Ma)
    10611064          CALL histwrite_phy(o_dnwd, dnwd)
    … …  
    12421245          ! etendue a iflag_con=3 (jyg)
    12431246          CALL histwrite_phy(o_Vprecip, Vprecip)
     1247          CALL histwrite_phy(o_qtaa, qtaa)
     1248          CALL histwrite_phy(o_clwaa, clw)
    12441249          CALL histwrite_phy(o_wdtrainA, wdtrainA)
     1250          CALL histwrite_phy(o_wdtrainS, wdtrainS)
    12451251          CALL histwrite_phy(o_wdtrainM, wdtrainM)
    12461252       ENDIF !(iflag_con.EQ.3.or.iflag_con.EQ.30)

  • LMDZ6/trunk/libf/phylmd/phys_state_var_mod.F90

    r3479 r3496  
    202202      REAL,ALLOCATABLE,SAVE :: ema_pcb(:), ema_pct(:)
    203203!$OMP THREADPRIVATE(ema_pcb,ema_pct)
    204       REAL,ALLOCATABLE,SAVE :: Ma(:,:)        ! undilute upward mass flux
     204      REAL,ALLOCATABLE,SAVE :: Mipsh(:,:)     ! mass flux shed from  adiab. ascents
     205!$OMP THREADPRIVATE(Mipsh)
     206      REAL,ALLOCATABLE,SAVE :: Ma(:,:)       ! undilute upward mass flux
    205207!$OMP THREADPRIVATE(Ma)
    206208      REAL,ALLOCATABLE,SAVE :: qcondc(:,:)    ! in-cld water content from convect
    … …  
    541543      ALLOCATE(ema_pcb(klon), ema_pct(klon))
    542544!
     545      ALLOCATE(Mipsh(klon,klev))
    543546      ALLOCATE(Ma(klon,klev))
    544547      ALLOCATE(qcondc(klon,klev))
    … …  
    699702      deallocate(ema_cbmf)
    700703      deallocate(ema_pcb, ema_pct)
    701       deallocate(Ma, qcondc)
     704      deallocate(Mipsh, Ma, qcondc)
    702705      deallocate(wd, sigd)
    703706      deallocate(cin, ALE, ALP)

  • LMDZ6/trunk/libf/phylmd/physiq_mod.F90

    r3491 r3496  
    172172       !    Deep convective variables used in phytrac
    173173       pmflxr, pmflxs,  &
    174        wdtrainA, wdtrainM,  &
     174       wdtrainA, wdtrainS, wdtrainM,  &
    175175       upwd, dnwd, &
    176176       ep,  &
    … …  
    182182       ev, &
    183183       elij, &
     184       qtaa, &
    184185       clw, &
    185186       epmlmMm, eplaMm, &
    … …  
    559560    ! Variables li\'ees \`a la poche froide (jyg)
    560561
    561     REAL mip(klon,klev)  ! mass flux shed by the adiab ascent at each level
     562!!    REAL mipsh(klon,klev)  ! mass flux shed by the adiab ascent at each level
     563!!      Moved to phys_state_var_mod
    562564    !
    563565    REAL wape_prescr, fip_prescr
    … …  
    26292631    pmflxs(:,:) = 0.
    26302632    wdtrainA(:,:) = 0.
     2633    wdtrainS(:,:) = 0.
    26312634    wdtrainM(:,:) = 0.
    26322635    upwd(:,:) = 0.
    … …  
    26442647    elij(:,:,:)=0.
    26452648    ev(:,:)=0.
     2649    qtaa(:,:)=0.
    26462650    clw(:,:)=0.
    26472651    sij(:,:,:)=0.
    … …  
    27802784               rain_con, snow_con, ibas_con, itop_con, sigd, &
    27812785               ema_cbmf,plcl,plfc,wbeff,convoccur,upwd,dnwd,dnwd0, &
    2782                Ma,mip,Vprecip,cape,cin,tvp,Tconv,iflagctrl, &
     2786               Ma,mipsh,Vprecip,cape,cin,tvp,Tconv,iflagctrl, &
    27832787               pbase,bbase,dtvpdt1,dtvpdq1,dplcldt,dplcldr,qcondc,wd, &
    27842788                                ! RomP >>>
    27852789                                !!     .        pmflxr,pmflxs,da,phi,mp,
    27862790                                !!     .        ftd,fqd,lalim_conv,wght_th)
    2787                pmflxr,pmflxs,da,phi,mp,phi2,d1a,dam,sij,clw,elij, &
     2791               pmflxr,pmflxs,da,phi,mp,phi2,d1a,dam,sij,qtaa,clw,elij, &
    27882792               ftd,fqd,lalim_conv,wght_th, &
    27892793               ev, ep,epmlmMm,eplaMm, &
    2790                wdtrainA,wdtrainM,wght_cvfd,qtc_cv,sigt_cv, &
     2794               wdtrainA, wdtrainS, wdtrainM,wght_cvfd,qtc_cv,sigt_cv, &
    27912795               tau_cld_cv,coefw_cld_cv,epmax_diag)
    27922796
Note: See TracChangeset for help on using the changeset viewer.