Ignore:
Timestamp:
Jul 29, 2024, 3:07:34 PM (8 weeks ago)
Author:
abarral
Message:

Put cvparam.h, fcg_gcssold.h, planete.h, tsoilnudge.h, YOECUMF.h into modules

File:
1 edited

Legend:

Unmodified
Added
Removed

  • LMDZ6/branches/Amaury_dev/libf/phylmd/conflx.F90

    r5112 r5142  
    1 
    21! $Header$
    32
    43SUBROUTINE conflx(dtime, pres_h, pres_f, t, q, con_t, con_q, pqhfl, w, d_t, &
    5     d_q, rain, snow, pmfu, pmfd, pen_u, pde_u, pen_d, pde_d, kcbot, kctop, &
    6     kdtop, pmflxr, pmflxs)
     4        d_q, rain, snow, pmfu, pmfd, pen_u, pde_u, pen_d, pde_d, kcbot, kctop, &
     5        kdtop, pmflxr, pmflxs)
    76
    87  USE dimphy
    … …  
    2019  ! Entree:
    2120  REAL dtime ! pas d'integration (s)
    22   REAL pres_h(klon, klev+1) ! pression half-level (Pa)
     21  REAL pres_h(klon, klev + 1) ! pression half-level (Pa)
    2322  REAL pres_f(klon, klev) ! pression full-level (Pa)
    2423  REAL t(klon, klev) ! temperature (K)
    … …  
    3938  REAL rain(klon) ! pluie (mm/s)
    4039  REAL snow(klon) ! neige (mm/s)
    41   REAL pmflxr(klon, klev+1)
    42   REAL pmflxs(klon, klev+1)
     40  REAL pmflxr(klon, klev + 1)
     41  REAL pmflxs(klon, klev + 1)
    4342  INTEGER kcbot(klon) ! niveau du bas de la convection
    4443  INTEGER kctop(klon) ! niveau du haut de la convection
    … …  
    5352  REAL d_t_bis(klon, klev)
    5453  REAL d_q_bis(klon, klev)
    55   REAL paprs(klon, klev+1)
     54  REAL paprs(klon, klev + 1)
    5655  REAL paprsf(klon, klev)
    5756  REAL zgeom(klon, klev)
    … …  
    6564  REAL zde_u(klon, klev)
    6665  REAL zde_d(klon, klev)
    67   REAL zmflxr(klon, klev+1)
    68   REAL zmflxs(klon, klev+1)
     66  REAL zmflxr(klon, klev + 1)
     67  REAL zmflxs(klon, klev + 1)
    6968  ! AA
    70 
    7169
    7270  INTEGER i, k
    … …  
    117115  DO k = 1, klev
    118116    DO i = 1, klon
    119       pt(i, k) = t(i, klev-k+1)
    120       pq(i, k) = q(i, klev-k+1)
    121       paprsf(i, k) = pres_f(i, klev-k+1)
    122       paprs(i, k) = pres_h(i, klev+1-k+1)
    123       pvervel(i, k) = w(i, klev+1-k)
    124       zcvgt(i, k) = con_t(i, klev-k+1)
    125       zcvgq(i, k) = con_q(i, klev-k+1)
    126 
    127       zdelta = max(0., sign(1.,rtt-pt(i,k)))
    128       zqsat = r2es*foeew(pt(i,k), zdelta)/paprsf(i, k)
     117      pt(i, k) = t(i, klev - k + 1)
     118      pq(i, k) = q(i, klev - k + 1)
     119      paprsf(i, k) = pres_f(i, klev - k + 1)
     120      paprs(i, k) = pres_h(i, klev + 1 - k + 1)
     121      pvervel(i, k) = w(i, klev + 1 - k)
     122      zcvgt(i, k) = con_t(i, klev - k + 1)
     123      zcvgq(i, k) = con_q(i, klev - k + 1)
     124
     125      zdelta = max(0., sign(1., rtt - pt(i, k)))
     126      zqsat = r2es * foeew(pt(i, k), zdelta) / paprsf(i, k)
    129127      zqsat = min(0.5, zqsat)
    130       zqsat = zqsat/(1.-retv*zqsat)
     128      zqsat = zqsat / (1. - retv * zqsat)
    131129      pqs(i, k) = zqsat
    132130    END DO
    133131  END DO
    134132  DO i = 1, klon
    135     paprs(i, klev+1) = pres_h(i, 1)
    136     zgeom(i, klev) = rd*pt(i, klev)/(0.5*(paprs(i,klev+1)+paprsf(i, &
    137       klev)))*(paprs(i,klev+1)-paprsf(i,klev))
     133    paprs(i, klev + 1) = pres_h(i, 1)
     134    zgeom(i, klev) = rd * pt(i, klev) / (0.5 * (paprs(i, klev + 1) + paprsf(i, &
     135            klev))) * (paprs(i, klev + 1) - paprsf(i, klev))
    138136  END DO
    139137  DO k = klev - 1, 1, -1
    140138    DO i = 1, klon
    141       zgeom(i, k) = zgeom(i, k+1) + rd*0.5*(pt(i,k+1)+pt(i,k))/paprs(i, k+1)* &
    142         (paprsf(i,k+1)-paprsf(i,k))
     139      zgeom(i, k) = zgeom(i, k + 1) + rd * 0.5 * (pt(i, k + 1) + pt(i, k)) / paprs(i, k + 1) * &
     140              (paprsf(i, k + 1) - paprsf(i, k))
    143141    END DO
    144142  END DO
    … …  
    147145
    148146  CALL flxmain(dtime, pt, pq, pqs, pqhfl, paprsf, paprs, zgeom, land, zcvgt, &
    149     zcvgq, pvervel, rain, snow, kcbot, kctop, kdtop, zmfu, zmfd, zen_u, &
    150     zde_u, zen_d, zde_d, d_t_bis, d_q_bis, zmflxr, zmflxs)
     147          zcvgq, pvervel, rain, snow, kcbot, kctop, kdtop, zmfu, zmfd, zen_u, &
     148          zde_u, zen_d, zde_d, d_t_bis, d_q_bis, zmflxr, zmflxs)
    151149
    152150  ! AA--------------------------------------------------------
    … …  
    158156  DO k = 1, klev
    159157    DO i = 1, klon
    160       d_q(i, klev+1-k) = dtime*d_q_bis(i, k)
    161       d_t(i, klev+1-k) = dtime*d_t_bis(i, k)
     158      d_q(i, klev + 1 - k) = dtime * d_q_bis(i, k)
     159      d_t(i, klev + 1 - k) = dtime * d_t_bis(i, k)
    162160    END DO
    163161  END DO
    … …  
    172170  DO k = 2, klev
    173171    DO i = 1, klon
    174       pmfu(i, klev+2-k) = zmfu(i, k)
    175       pmfd(i, klev+2-k) = zmfd(i, k)
     172      pmfu(i, klev + 2 - k) = zmfu(i, k)
     173      pmfd(i, klev + 2 - k) = zmfd(i, k)
    176174    END DO
    177175  END DO
    … …  
    179177  DO k = 1, klev
    180178    DO i = 1, klon
    181       pen_u(i, klev+1-k) = zen_u(i, k)
    182       pde_u(i, klev+1-k) = zde_u(i, k)
     179      pen_u(i, klev + 1 - k) = zen_u(i, k)
     180      pde_u(i, klev + 1 - k) = zde_u(i, k)
    183181    END DO
    184182  END DO
    … …  
    186184  DO k = 1, klev - 1
    187185    DO i = 1, klon
    188       pen_d(i, klev+1-k) = -zen_d(i, k+1)
    189       pde_d(i, klev+1-k) = -zde_d(i, k+1)
     186      pen_d(i, klev + 1 - k) = -zen_d(i, k + 1)
     187      pde_d(i, klev + 1 - k) = -zde_d(i, k + 1)
    190188    END DO
    191189  END DO
    … …  
    193191  DO k = 1, klev + 1
    194192    DO i = 1, klon
    195       pmflxr(i, klev+2-k) = zmflxr(i, k)
    196       pmflxs(i, klev+2-k) = zmflxs(i, k)
    197     END DO
    198   END DO
    199 
     193      pmflxr(i, klev + 2 - k) = zmflxr(i, k)
     194      pmflxs(i, klev + 2 - k) = zmflxs(i, k)
     195    END DO
     196  END DO
    200197
    201198END SUBROUTINE conflx
    202199! --------------------------------------------------------------------
    203200SUBROUTINE flxmain(pdtime, pten, pqen, pqsen, pqhfl, pap, paph, pgeo, ldland, &
    204     ptte, pqte, pvervel, prsfc, pssfc, kcbot, kctop, kdtop, & ! *
    205                                                               ! ldcum, ktype,
    206     pmfu, pmfd, pen_u, pde_u, pen_d, pde_d, dt_con, dq_con, pmflxr, pmflxs)
     201        ptte, pqte, pvervel, prsfc, pssfc, kcbot, kctop, kdtop, & ! *
     202        ! ldcum, ktype,
     203        pmfu, pmfd, pen_u, pde_u, pen_d, pde_d, dt_con, dq_con, pmflxr, pmflxs)
    207204  USE dimphy
     205  USE lmdz_YOECUMF
     206
    208207  IMPLICIT NONE
    209208  ! ------------------------------------------------------------------
    210209  include "YOMCST.h"
    211210  include "YOETHF.h"
    212   include "YOECUMF.h"
    213211  ! ----------------------------------------------------------------
    214212  REAL pten(klon, klev), pqen(klon, klev), pqsen(klon, klev)
    … …  
    216214  REAL pqte(klon, klev)
    217215  REAL pvervel(klon, klev)
    218   REAL pgeo(klon, klev), pap(klon, klev), paph(klon, klev+1)
     216  REAL pgeo(klon, klev), pap(klon, klev), paph(klon, klev + 1)
    219217  REAL pqhfl(klon)
    220218
    … …  
    234232  REAL zdqpbl(klon), zdqcv(klon), zdhpbl(klon)
    235233  REAL zrfl(klon)
    236   REAL pmflxr(klon, klev+1)
    237   REAL pmflxs(klon, klev+1)
     234  REAL pmflxr(klon, klev + 1)
     235  REAL pmflxs(klon, klev + 1)
    238236  INTEGER ilab(klon, klev), ictop0(klon)
    239237  LOGICAL llo1
    … …  
    275273  ! ----------------------------------------------------------------------
    276274  CALL flxini(pten, pqen, pqsen, pgeo, paph, zgeoh, ztenh, zqenh, zqsenh, &
    277     ptu, pqu, ptd, pqd, pmfd, zmfds, zmfdq, zdmfdp, pmfu, zmfus, zmfuq, &
    278     zdmfup, zdpmel, plu, plude, ilab, pen_u, pde_u, pen_d, pde_d)
     275          ptu, pqu, ptd, pqd, pmfd, zmfds, zmfdq, zdmfdp, pmfu, zmfus, zmfuq, &
     276          zdmfup, zdpmel, plu, plude, ilab, pen_u, pde_u, pen_d, pde_d)
    279277  ! ---------------------------------------------------------------------
    280278  ! determiner les valeurs au niveau de base de la tour convective
    … …  
    290288  k = 1
    291289  DO i = 1, klon
    292     zdqcv(i) = pqte(i, k)*(paph(i,k+1)-paph(i,k))
     290    zdqcv(i) = pqte(i, k) * (paph(i, k + 1) - paph(i, k))
    293291    zdhpbl(i) = 0.0
    294292    zdqpbl(i) = 0.0
    … …  
    297295  DO k = 2, klev
    298296    DO i = 1, klon
    299       zdqcv(i) = zdqcv(i) + pqte(i, k)*(paph(i,k+1)-paph(i,k))
     297      zdqcv(i) = zdqcv(i) + pqte(i, k) * (paph(i, k + 1) - paph(i, k))
    300298      IF (k>=kcbot(i)) THEN
    301         zdqpbl(i) = zdqpbl(i) + pqte(i, k)*(paph(i,k+1)-paph(i,k))
    302         zdhpbl(i) = zdhpbl(i) + (rcpd*ptte(i,k)+rlvtt*pqte(i,k))*(paph(i,k+1) &
    303           -paph(i,k))
     299        zdqpbl(i) = zdqpbl(i) + pqte(i, k) * (paph(i, k + 1) - paph(i, k))
     300        zdhpbl(i) = zdhpbl(i) + (rcpd * ptte(i, k) + rlvtt * pqte(i, k)) * (paph(i, k + 1) &
     301                - paph(i, k))
    304302      END IF
    305303    END DO
    … …  
    308306  DO i = 1, klon
    309307    ktype(i) = 2
    310     IF (zdqcv(i)>max(0.,-1.5*pqhfl(i)*rg)) ktype(i) = 1
     308    IF (zdqcv(i)>max(0., -1.5 * pqhfl(i) * rg)) ktype(i) = 1
    311309    ! cc         if (zdqcv(i).GT.MAX(0.,-1.1*pqhfl(i)*RG)) ktype(i) = 1
    312310  END DO
    … …  
    319317    ikb = kcbot(i)
    320318    zqumqe = pqu(i, ikb) + plu(i, ikb) - zqenh(i, ikb)
    321     zdqmin = max(0.01*zqenh(i,ikb), 1.E-10)
     319    zdqmin = max(0.01 * zqenh(i, ikb), 1.E-10)
    322320    IF (zdqpbl(i)>0. .AND. zqumqe>zdqmin .AND. ldcum(i)) THEN
    323       zmfub(i) = zdqpbl(i)/(rg*max(zqumqe,zdqmin))
     321      zmfub(i) = zdqpbl(i) / (rg * max(zqumqe, zdqmin))
    324322    ELSE
    325323      zmfub(i) = 0.01
    … …  
    327325    END IF
    328326    IF (ktype(i)==2) THEN
    329       zdh = rcpd*(ptu(i,ikb)-ztenh(i,ikb)) + rlvtt*zqumqe
    330       zdh = rg*max(zdh, 1.0E5*zdqmin)
    331       IF (zdhpbl(i)>0. .AND. ldcum(i)) zmfub(i) = zdhpbl(i)/zdh
     327      zdh = rcpd * (ptu(i, ikb) - ztenh(i, ikb)) + rlvtt * zqumqe
     328      zdh = rg * max(zdh, 1.0E5 * zdqmin)
     329      IF (zdhpbl(i)>0. .AND. ldcum(i)) zmfub(i) = zdhpbl(i) / zdh
    332330    END IF
    333     zmfmax = (paph(i,ikb)-paph(i,ikb-1))/(rg*pdtime)
     331    zmfmax = (paph(i, ikb) - paph(i, ikb - 1)) / (rg * pdtime)
    334332    zmfub(i) = min(zmfub(i), zmfmax)
    335333    zentr(i) = entrscv
    … …  
    345343  DO i = 1, klon
    346344    ikb = kcbot(i)
    347     zhcbase(i) = rcpd*ptu(i, ikb) + zgeoh(i, ikb) + rlvtt*pqu(i, ikb)
     345    zhcbase(i) = rcpd * ptu(i, ikb) + zgeoh(i, ikb) + rlvtt * pqu(i, ikb)
    348346    ictop0(i) = kcbot(i) - 1
    349347  END DO
    350348
    351   zalvdcp = rlvtt/rcpd
     349  zalvdcp = rlvtt / rcpd
    352350  DO k = klev - 1, 3, -1
    353351    DO i = 1, klon
    354       zhsat = rcpd*ztenh(i, k) + zgeoh(i, k) + rlvtt*zqsenh(i, k)
    355       zgam = r5les*zalvdcp*zqsenh(i, k)/((1.-retv*zqsenh(i,k))*(ztenh(i, &
    356         k)-r4les)**2)
    357       zzz = rcpd*ztenh(i, k)*0.608
    358       zhhat = zhsat - (zzz+zgam*zzz)/(1.+zgam*zzz/rlvtt)*max(zqsenh(i,k)- &
    359         zqenh(i,k), 0.)
     352      zhsat = rcpd * ztenh(i, k) + zgeoh(i, k) + rlvtt * zqsenh(i, k)
     353      zgam = r5les * zalvdcp * zqsenh(i, k) / ((1. - retv * zqsenh(i, k)) * (ztenh(i, &
     354              k) - r4les)**2)
     355      zzz = rcpd * ztenh(i, k) * 0.608
     356      zhhat = zhsat - (zzz + zgam * zzz) / (1. + zgam * zzz / rlvtt) * max(zqsenh(i, k) - &
     357              zqenh(i, k), 0.)
    360358      IF (k<ictop0(i) .AND. zhcbase(i)>zhhat) ictop0(i) = k
    361359    END DO
    … …  
    365363
    366364  CALL flxasc(pdtime, ztenh, zqenh, pten, pqen, pqsen, pgeo, zgeoh, pap, &
    367     paph, pqte, pvervel, ldland, ldcum, ktype, ilab, ptu, pqu, plu, pmfu, &
    368     zmfub, zentr, zmfus, zmfuq, zmful, plude, zdmfup, kcbot, kctop, ictop0, &
    369     kcum, pen_u, pde_u)
     365          paph, pqte, pvervel, ldland, ldcum, ktype, ilab, ptu, pqu, plu, pmfu, &
     366          zmfub, zentr, zmfus, zmfuq, zmful, plude, zdmfup, kcbot, kctop, ictop0, &
     367          kcum, pen_u, pde_u)
    370368  IF (kcum==0) GO TO 1000
    371369
    … …  
    395393    ! determiner le LFS (level of free sinking: niveau de plonge libre)
    396394    CALL flxdlfs(ztenh, zqenh, zgeoh, paph, ptu, pqu, ldcum, kcbot, kctop, &
    397       zmfub, zrfl, ptd, pqd, pmfd, zmfds, zmfdq, zdmfdp, kdtop, lddraf)
     395            zmfub, zrfl, ptd, pqd, pmfd, zmfds, zmfdq, zdmfdp, kdtop, lddraf)
    398396
    399397    ! calculer le panache descendant
    400398    CALL flxddraf(ztenh, zqenh, zgeoh, paph, zrfl, ptd, pqd, pmfd, zmfds, &
    401       zmfdq, zdmfdp, lddraf, pen_d, pde_d)
     399            zmfdq, zdmfdp, lddraf, pen_d, pde_d)
    402400
    403401    ! calculer de nouveau le flux de masse entrant a travers la base
    … …  
    410408        zeps = 0.
    411409        IF (llo1) zeps = cmfdeps
    412         zqumqe = pqu(i, ikb) + plu(i, ikb) - zeps*pqd(i, ikb) - &
    413           (1.-zeps)*zqenh(i, ikb)
    414         zdqmin = max(0.01*zqenh(i,ikb), 1.E-10)
    415         zmfmax = (paph(i,ikb)-paph(i,ikb-1))/(rg*pdtime)
     410        zqumqe = pqu(i, ikb) + plu(i, ikb) - zeps * pqd(i, ikb) - &
     411                (1. - zeps) * zqenh(i, ikb)
     412        zdqmin = max(0.01 * zqenh(i, ikb), 1.E-10)
     413        zmfmax = (paph(i, ikb) - paph(i, ikb - 1)) / (rg * pdtime)
    416414        IF (zdqpbl(i)>0. .AND. zqumqe>zdqmin .AND. ldcum(i) .AND. &
    417             zmfub(i)<zmfmax) THEN
    418           zmfub1(i) = zdqpbl(i)/(rg*max(zqumqe,zdqmin))
     415                zmfub(i)<zmfmax) THEN
     416          zmfub1(i) = zdqpbl(i) / (rg * max(zqumqe, zdqmin))
    419417        ELSE
    420418          zmfub1(i) = zmfub(i)
    421419        END IF
    422420        IF (ktype(i)==2) THEN
    423           zdh = rcpd*(ptu(i,ikb)-zeps*ptd(i,ikb)-(1.-zeps)*ztenh(i,ikb)) + &
    424             rlvtt*zqumqe
    425           zdh = rg*max(zdh, 1.0E5*zdqmin)
    426           IF (zdhpbl(i)>0. .AND. ldcum(i)) zmfub1(i) = zdhpbl(i)/zdh
     421          zdh = rcpd * (ptu(i, ikb) - zeps * ptd(i, ikb) - (1. - zeps) * ztenh(i, ikb)) + &
     422                  rlvtt * zqumqe
     423          zdh = rg * max(zdh, 1.0E5 * zdqmin)
     424          IF (zdhpbl(i)>0. .AND. ldcum(i)) zmfub1(i) = zdhpbl(i) / zdh
    427425        END IF
    428         IF (.NOT. ((ktype(i)==1 .OR. ktype(i)==2) .AND. abs(zmfub1(i)-zmfub(i &
    429           ))<0.2*zmfub(i))) zmfub1(i) = zmfub(i)
     426        IF (.NOT. ((ktype(i)==1 .OR. ktype(i)==2) .AND. abs(zmfub1(i) - zmfub(i &
     427                ))<0.2 * zmfub(i))) zmfub1(i) = zmfub(i)
    430428      END IF
    431429    END DO
    … …  
    433431      DO i = 1, klon
    434432        IF (lddraf(i)) THEN
    435           zfac = zmfub1(i)/max(zmfub(i), 1.E-10)
    436           pmfd(i, k) = pmfd(i, k)*zfac
    437           zmfds(i, k) = zmfds(i, k)*zfac
    438           zmfdq(i, k) = zmfdq(i, k)*zfac
    439           zdmfdp(i, k) = zdmfdp(i, k)*zfac
    440           pen_d(i, k) = pen_d(i, k)*zfac
    441           pde_d(i, k) = pde_d(i, k)*zfac
     433          zfac = zmfub1(i) / max(zmfub(i), 1.E-10)
     434          pmfd(i, k) = pmfd(i, k) * zfac
     435          zmfds(i, k) = zmfds(i, k) * zfac
     436          zmfdq(i, k) = zmfdq(i, k) * zfac
     437          zdmfdp(i, k) = zdmfdp(i, k) * zfac
     438          pen_d(i, k) = pen_d(i, k) * zfac
     439          pde_d(i, k) = pde_d(i, k) * zfac
    442440        END IF
    443441      END DO
    … …  
    453451  ! -----------------------------------------------------------------------
    454452  CALL flxasc(pdtime, ztenh, zqenh, pten, pqen, pqsen, pgeo, zgeoh, pap, &
    455     paph, pqte, pvervel, ldland, ldcum, ktype, ilab, ptu, pqu, plu, pmfu, &
    456     zmfub, zentr, zmfus, zmfuq, zmful, plude, zdmfup, kcbot, kctop, ictop0, &
    457     kcum, pen_u, pde_u)
     453          paph, pqte, pvervel, ldland, ldcum, ktype, ilab, ptu, pqu, plu, pmfu, &
     454          zmfub, zentr, zmfus, zmfuq, zmful, plude, zdmfup, kcbot, kctop, ictop0, &
     455          kcum, pen_u, pde_u)
    458456
    459457  ! -----------------------------------------------------------------------
    … …  
    462460  ! -----------------------------------------------------------------------
    463461  CALL flxflux(pdtime, pqen, pqsen, ztenh, zqenh, pap, paph, ldland, zgeoh, &
    464     kcbot, kctop, lddraf, kdtop, ktype, ldcum, pmfu, pmfd, zmfus, zmfds, &
    465     zmfuq, zmfdq, zmful, plude, zdmfup, zdmfdp, pten, prsfc, pssfc, zdpmel, &
    466     itopm2, pmflxr, pmflxs)
     462          kcbot, kctop, lddraf, kdtop, ktype, ldcum, pmfu, pmfd, zmfus, zmfds, &
     463          zmfuq, zmfdq, zmful, plude, zdmfup, zdmfdp, pten, prsfc, pssfc, zdpmel, &
     464          itopm2, pmflxr, pmflxs)
    467465
    468466  ! ----------------------------------------------------------------------
    … …  
    470468  ! ----------------------------------------------------------------------
    471469  CALL flxdtdq(pdtime, itopm2, paph, ldcum, pten, zmfus, zmfds, zmfuq, zmfdq, &
    472     zmful, zdmfup, zdmfdp, zdpmel, dt_con, dq_con)
    473 
    474 1000 CONTINUE
     470          zmful, zdmfup, zdmfdp, zdpmel, dt_con, dq_con)
     471
     472  1000 CONTINUE
    475473
    476474END SUBROUTINE flxmain
    477475SUBROUTINE flxini(pten, pqen, pqsen, pgeo, paph, pgeoh, ptenh, pqenh, pqsenh, &
    478     ptu, pqu, ptd, pqd, pmfd, pmfds, pmfdq, pdmfdp, pmfu, pmfus, pmfuq, &
    479     pdmfup, pdpmel, plu, plude, klab, pen_u, pde_u, pen_d, pde_d)
     476        ptu, pqu, ptd, pqd, pmfd, pmfds, pmfdq, pdmfdp, pmfu, pmfus, pmfuq, &
     477        pdmfup, pdpmel, plu, plude, klab, pen_u, pde_u, pen_d, pde_d)
    480478  USE dimphy
    481479  IMPLICIT NONE
    … …  
    493491  REAL pgeo(klon, klev) ! geopotentiel (g * metre)
    494492  REAL pgeoh(klon, klev) ! geopotentiel aux demi-niveaux
    495   REAL paph(klon, klev+1) ! pression aux demi-niveaux
     493  REAL paph(klon, klev + 1) ! pression aux demi-niveaux
    496494  REAL ptenh(klon, klev) ! temperature aux demi-niveaux
    497495  REAL pqenh(klon, klev) ! humidite aux demi-niveaux
    … …  
    532530
    533531    DO i = 1, klon
    534       pgeoh(i, k) = pgeo(i, k) + (pgeo(i,k-1)-pgeo(i,k))*0.5
    535       ptenh(i, k) = (max(rcpd*pten(i,k-1)+pgeo(i,k-1),rcpd*pten(i,k)+pgeo(i, &
    536         k))-pgeoh(i,k))/rcpd
    537       pqsenh(i, k) = pqsen(i, k-1)
     532      pgeoh(i, k) = pgeo(i, k) + (pgeo(i, k - 1) - pgeo(i, k)) * 0.5
     533      ptenh(i, k) = (max(rcpd * pten(i, k - 1) + pgeo(i, k - 1), rcpd * pten(i, k) + pgeo(i, &
     534              k)) - pgeoh(i, k)) / rcpd
     535      pqsenh(i, k) = pqsen(i, k - 1)
    538536      llflag(i) = .TRUE.
    539537    END DO
    540538
    541539    iCALL = 0
    542     CALL flxadjtq(paph(1,k), ptenh(1,k), pqsenh(1,k), llflag, icall)
    543 
    544     DO i = 1, klon
    545       pqenh(i, k) = min(pqen(i,k-1), pqsen(i,k-1)) + &
    546         (pqsenh(i,k)-pqsen(i,k-1))
    547       pqenh(i, k) = max(pqenh(i,k), 0.)
    548     END DO
    549 
    550   END DO
    551 
    552   DO i = 1, klon
    553     ptenh(i, klev) = (rcpd*pten(i,klev)+pgeo(i,klev)-pgeoh(i,klev))/rcpd
     540    CALL flxadjtq(paph(1, k), ptenh(1, k), pqsenh(1, k), llflag, icall)
     541
     542    DO i = 1, klon
     543      pqenh(i, k) = min(pqen(i, k - 1), pqsen(i, k - 1)) + &
     544              (pqsenh(i, k) - pqsen(i, k - 1))
     545      pqenh(i, k) = max(pqenh(i, k), 0.)
     546    END DO
     547
     548  END DO
     549
     550  DO i = 1, klon
     551    ptenh(i, klev) = (rcpd * pten(i, klev) + pgeo(i, klev) - pgeoh(i, klev)) / rcpd
    554552    pqenh(i, klev) = pqen(i, klev)
    555553    ptenh(i, 1) = pten(i, 1)
    … …  
    560558  DO k = klev - 1, 2, -1
    561559    DO i = 1, klon
    562       zzs = max(rcpd*ptenh(i,k)+pgeoh(i,k), rcpd*ptenh(i,k+1)+pgeoh(i,k+1))
    563       ptenh(i, k) = (zzs-pgeoh(i,k))/rcpd
     560      zzs = max(rcpd * ptenh(i, k) + pgeoh(i, k), rcpd * ptenh(i, k + 1) + pgeoh(i, k + 1))
     561      ptenh(i, k) = (zzs - pgeoh(i, k)) / rcpd
    564562    END DO
    565563  END DO
    … …  
    596594  END DO
    597595
    598 
    599596END SUBROUTINE flxini
    600597SUBROUTINE flxbase(ptenh, pqenh, pgeoh, paph, ptu, pqu, plu, ldcum, kcbot, &
    601     klab)
     598        klab)
    602599  USE dimphy
    603600  IMPLICIT NONE
    … …  
    617614  ! ----------------------------------------------------------------
    618615  REAL ptenh(klon, klev), pqenh(klon, klev)
    619   REAL pgeoh(klon, klev), paph(klon, klev+1)
     616  REAL pgeoh(klon, klev), paph(klon, klev + 1)
    620617
    621618  REAL ptu(klon, klev), pqu(klon, klev), plu(klon, klev)
    … …  
    643640    is = 0
    644641    DO i = 1, klon
    645       IF (klab(i,k+1)==1) is = is + 1
     642      IF (klab(i, k + 1)==1) is = is + 1
    646643      llflag(i) = .FALSE.
    647       IF (klab(i,k+1)==1) llflag(i) = .TRUE.
     644      IF (klab(i, k + 1)==1) llflag(i) = .TRUE.
    648645    END DO
    649646    IF (is==0) GO TO 290
    … …  
    651648    DO i = 1, klon
    652649      IF (llflag(i)) THEN
    653         pqu(i, k) = pqu(i, k+1)
    654         ptu(i, k) = ptu(i, k+1) + (pgeoh(i,k+1)-pgeoh(i,k))/rcpd
    655         zbuo = ptu(i, k)*(1.+retv*pqu(i,k)) - ptenh(i, k)*(1.+retv*pqenh(i,k) &
    656           ) + 0.5
     650        pqu(i, k) = pqu(i, k + 1)
     651        ptu(i, k) = ptu(i, k + 1) + (pgeoh(i, k + 1) - pgeoh(i, k)) / rcpd
     652        zbuo = ptu(i, k) * (1. + retv * pqu(i, k)) - ptenh(i, k) * (1. + retv * pqenh(i, k) &
     653                ) + 0.5
    657654        IF (zbuo>0.) klab(i, k) = 1
    658655        zqold(i) = pqu(i, k)
    … …  
    661658
    662659    iCALL = 1
    663     CALL flxadjtq(paph(1,k), ptu(1,k), pqu(1,k), llflag, icall)
    664 
    665     DO i = 1, klon
    666       IF (llflag(i) .AND. pqu(i,k)/=zqold(i)) THEN
     660    CALL flxadjtq(paph(1, k), ptu(1, k), pqu(1, k), llflag, icall)
     661
     662    DO i = 1, klon
     663      IF (llflag(i) .AND. pqu(i, k)/=zqold(i)) THEN
    667664        klab(i, k) = 2
    668665        plu(i, k) = plu(i, k) + zqold(i) - pqu(i, k)
    669         zbuo = ptu(i, k)*(1.+retv*pqu(i,k)) - ptenh(i, k)*(1.+retv*pqenh(i,k) &
    670           ) + 0.5
     666        zbuo = ptu(i, k) * (1. + retv * pqu(i, k)) - ptenh(i, k) * (1. + retv * pqenh(i, k) &
     667                ) + 0.5
    671668        IF (zbuo>0.) kcbot(i) = k
    672669        IF (zbuo>0.) ldcum(i) = .TRUE.
    … …  
    674671    END DO
    675672
    676 290 END DO
    677 
     673  290 END DO
    678674
    679675END SUBROUTINE flxbase
    680676SUBROUTINE flxasc(pdtime, ptenh, pqenh, pten, pqen, pqsen, pgeo, pgeoh, pap, &
    681     paph, pqte, pvervel, ldland, ldcum, ktype, klab, ptu, pqu, plu, pmfu, &
    682     pmfub, pentr, pmfus, pmfuq, pmful, plude, pdmfup, kcbot, kctop, kctop0, &
    683     kcum, pen_u, pde_u)
     677        paph, pqte, pvervel, ldland, ldcum, ktype, klab, ptu, pqu, plu, pmfu, &
     678        pmfub, pentr, pmfus, pmfuq, pmful, plude, pdmfup, kcbot, kctop, kctop0, &
     679        kcum, pen_u, pde_u)
    684680  USE dimphy
     681  USE lmdz_YOECUMF
     682
    685683  IMPLICIT NONE
    686684  ! ----------------------------------------------------------------------
    … …  
    690688  include "YOMCST.h"
    691689  include "YOETHF.h"
    692   include "YOECUMF.h"
    693690
    694691  REAL pdtime
    … …  
    696693  REAL pqen(klon, klev), pqenh(klon, klev), pqsen(klon, klev)
    697694  REAL pgeo(klon, klev), pgeoh(klon, klev)
    698   REAL pap(klon, klev), paph(klon, klev+1)
     695  REAL pap(klon, klev), paph(klon, klev + 1)
    699696  REAL pqte(klon, klev)
    700697  REAL pvervel(klon, klev) ! vitesse verticale en Pa/s
    … …  
    735732  DO k = klev, 3, -1
    736733    DO i = 1, klon
    737       IF (pvervel(i,k)<zwmax(i)) THEN
     734      IF (pvervel(i, k)<zwmax(i)) THEN
    738735        zwmax(i) = pvervel(i, k)
    739736        klwmin(i) = k
    … …  
    758755      pdmfup(i, k) = 0.
    759756      IF (.NOT. ldcum(i) .OR. ktype(i)==3) klab(i, k) = 0
    760       IF (.NOT. ldcum(i) .AND. paph(i,k)<4.E4) kctop0(i) = k
     757      IF (.NOT. ldcum(i) .AND. paph(i, k)<4.E4) kctop0(i) = k
    761758    END DO
    762759  END DO
    … …  
    766763      zdland(i) = 3.0E4
    767764      zdphi = pgeoh(i, kctop0(i)) - pgeoh(i, kcbot(i))
    768       IF (ptu(i,kctop0(i))>=ztglace) zdland(i) = zdphi
     765      IF (ptu(i, kctop0(i))>=ztglace) zdland(i) = zdphi
    769766      zdland(i) = max(3.0E4, zdland(i))
    770767      zdland(i) = min(5.0E4, zdland(i))
    … …  
    782779    END IF
    783780    pmfu(i, klev) = pmfub(i)
    784     pmfus(i, klev) = pmfub(i)*(rcpd*ptu(i,klev)+pgeoh(i,klev))
    785     pmfuq(i, klev) = pmfub(i)*pqu(i, klev)
     781    pmfus(i, klev) = pmfub(i) * (rcpd * ptu(i, klev) + pgeoh(i, klev))
     782    pmfuq(i, klev) = pmfub(i) * pqu(i, klev)
    786783  END DO
    787784
    … …  
    797794  DO k = klev - 1, 3, -1
    798795
    799     IF (lmfmid .AND. k<klev-1) THEN
     796    IF (lmfmid .AND. k<klev - 1) THEN
    800797      DO i = 1, klon
    801         IF (.NOT. ldcum(i) .AND. klab(i,k+1)==0 .AND. &
    802             pqen(i,k)>0.9*pqsen(i,k) .AND. pap(i,k)/paph(i,klev+1)>0.4) THEN
    803           ptu(i, k+1) = pten(i, k) + (pgeo(i,k)-pgeoh(i,k+1))/rcpd
    804           pqu(i, k+1) = pqen(i, k)
    805           plu(i, k+1) = 0.0
    806           zzzmb = max(cmfcmin, -pvervel(i,k)/rg)
    807           zmfmax = (paph(i,k)-paph(i,k-1))/(rg*pdtime)
     798        IF (.NOT. ldcum(i) .AND. klab(i, k + 1)==0 .AND. &
     799                pqen(i, k)>0.9 * pqsen(i, k) .AND. pap(i, k) / paph(i, klev + 1)>0.4) THEN
     800          ptu(i, k + 1) = pten(i, k) + (pgeo(i, k) - pgeoh(i, k + 1)) / rcpd
     801          pqu(i, k + 1) = pqen(i, k)
     802          plu(i, k + 1) = 0.0
     803          zzzmb = max(cmfcmin, -pvervel(i, k) / rg)
     804          zmfmax = (paph(i, k) - paph(i, k - 1)) / (rg * pdtime)
    808805          pmfub(i) = min(zzzmb, zmfmax)
    809           pmfu(i, k+1) = pmfub(i)
    810           pmfus(i, k+1) = pmfub(i)*(rcpd*ptu(i,k+1)+pgeoh(i,k+1))
    811           pmfuq(i, k+1) = pmfub(i)*pqu(i, k+1)
    812           pmful(i, k+1) = 0.0
    813           pdmfup(i, k+1) = 0.0
     806          pmfu(i, k + 1) = pmfub(i)
     807          pmfus(i, k + 1) = pmfub(i) * (rcpd * ptu(i, k + 1) + pgeoh(i, k + 1))
     808          pmfuq(i, k + 1) = pmfub(i) * pqu(i, k + 1)
     809          pmful(i, k + 1) = 0.0
     810          pdmfup(i, k + 1) = 0.0
    814811          kcbot(i) = k
    815           klab(i, k+1) = 1
     812          klab(i, k + 1) = 1
    816813          ktype(i) = 3
    817814          pentr(i) = entrmid
    … …  
    822819    is = 0
    823820    DO i = 1, klon
    824       is = is + klab(i, k+1)
    825       IF (klab(i,k+1)==0) klab(i, k) = 0
     821      is = is + klab(i, k + 1)
     822      IF (klab(i, k + 1)==0) klab(i, k) = 0
    826823      llflag(i) = .FALSE.
    827       IF (klab(i,k+1)>0) llflag(i) = .TRUE.
     824      IF (klab(i, k + 1)>0) llflag(i) = .TRUE.
    828825    END DO
    829826    IF (is==0) GO TO 480
    … …  
    834831      pen_u(i, k) = 0.0
    835832      pde_u(i, k) = 0.0
    836       zrho(i) = paph(i, k+1)/(rd*ptenh(i,k+1))
     833      zrho(i) = paph(i, k + 1) / (rd * ptenh(i, k + 1))
    837834      zpbot(i) = paph(i, kcbot(i))
    838835      zptop(i) = paph(i, kctop0(i))
    … …  
    841838    DO i = 1, klon
    842839      IF (ldcum(i)) THEN
    843         zdprho = (paph(i,k+1)-paph(i,k))/(rg*zrho(i))
    844         zentr = pentr(i)*pmfu(i, k+1)*zdprho
     840        zdprho = (paph(i, k + 1) - paph(i, k)) / (rg * zrho(i))
     841        zentr = pentr(i) * pmfu(i, k + 1) * zdprho
    845842        llo1 = k < kcbot(i)
    846843        IF (llo1) pde_u(i, k) = zentr
    847         zpmid = 0.5*(zpbot(i)+zptop(i))
    848         llo2 = llo1 .AND. ktype(i) == 2 .AND. (zpbot(i)-paph(i,k)<0.2E5 .OR. &
    849           paph(i,k)>zpmid)
     844        zpmid = 0.5 * (zpbot(i) + zptop(i))
     845        llo2 = llo1 .AND. ktype(i) == 2 .AND. (zpbot(i) - paph(i, k)<0.2E5 .OR. &
     846                paph(i, k)>zpmid)
    850847        IF (llo2) pen_u(i, k) = zentr
    851848        llo2 = llo1 .AND. (ktype(i)==1 .OR. ktype(i)==3) .AND. &
    852           (k>=max(klwmin(i),kctop0(i)+2) .OR. pap(i,k)>zpmid)
     849                (k>=max(klwmin(i), kctop0(i) + 2) .OR. pap(i, k)>zpmid)
    853850        IF (llo2) pen_u(i, k) = zentr
    854851        llo1 = pen_u(i, k) > 0. .AND. (ktype(i)==1 .OR. ktype(i)==2)
    855852        IF (llo1) THEN
    856           fact = 1. + 3.*(1.-min(1.,(zpbot(i)-pap(i,k))/1.5E4))
    857           zentr = zentr*fact
    858           pen_u(i, k) = pen_u(i, k)*fact
    859           pde_u(i, k) = pde_u(i, k)*fact
     853          fact = 1. + 3. * (1. - min(1., (zpbot(i) - pap(i, k)) / 1.5E4))
     854          zentr = zentr * fact
     855          pen_u(i, k) = pen_u(i, k) * fact
     856          pde_u(i, k) = pde_u(i, k) * fact
    860857        END IF
    861         IF (llo2 .AND. pqenh(i,k+1)>1.E-5) pen_u(i, k) = zentr + &
    862           max(pqte(i,k), 0.)/pqenh(i, k+1)*zrho(i)*zdprho
     858        IF (llo2 .AND. pqenh(i, k + 1)>1.E-5) pen_u(i, k) = zentr + &
     859                max(pqte(i, k), 0.) / pqenh(i, k + 1) * zrho(i) * zdprho
    863860      END IF
    864861    END DO
    … …  
    871868      IF (llflag(i)) THEN
    872869        IF (k<kcbot(i)) THEN
    873           zmftest = pmfu(i, k+1) + pen_u(i, k) - pde_u(i, k)
    874           zmfmax = min(zmftest, (paph(i,k)-paph(i,k-1))/(rg*pdtime))
    875           pen_u(i, k) = max(pen_u(i,k)-max(0.0,zmftest-zmfmax), 0.0)
     870          zmftest = pmfu(i, k + 1) + pen_u(i, k) - pde_u(i, k)
     871          zmfmax = min(zmftest, (paph(i, k) - paph(i, k - 1)) / (rg * pdtime))
     872          pen_u(i, k) = max(pen_u(i, k) - max(0.0, zmftest - zmfmax), 0.0)
    876873        END IF
    877         pde_u(i, k) = min(pde_u(i,k), 0.75*pmfu(i,k+1))
     874        pde_u(i, k) = min(pde_u(i, k), 0.75 * pmfu(i, k + 1))
    878875        ! calculer le flux de masse du niveau k a partir de celui du k+1
    879         pmfu(i, k) = pmfu(i, k+1) + pen_u(i, k) - pde_u(i, k)
     876        pmfu(i, k) = pmfu(i, k + 1) + pen_u(i, k) - pde_u(i, k)
    880877        ! calculer les valeurs Su, Qu et l du niveau k dans le panache
    881878        ! montant
    882         zqeen = pqenh(i, k+1)*pen_u(i, k)
    883         zseen = (rcpd*ptenh(i,k+1)+pgeoh(i,k+1))*pen_u(i, k)
    884         zscde = (rcpd*ptu(i,k+1)+pgeoh(i,k+1))*pde_u(i, k)
    885         zqude = pqu(i, k+1)*pde_u(i, k)
    886         plude(i, k) = plu(i, k+1)*pde_u(i, k)
    887         zmfusk = pmfus(i, k+1) + zseen - zscde
    888         zmfuqk = pmfuq(i, k+1) + zqeen - zqude
    889         zmfulk = pmful(i, k+1) - plude(i, k)
    890         plu(i, k) = zmfulk*(1./max(cmfcmin,pmfu(i,k)))
    891         pqu(i, k) = zmfuqk*(1./max(cmfcmin,pmfu(i,k)))
    892         ptu(i, k) = (zmfusk*(1./max(cmfcmin,pmfu(i,k)))-pgeoh(i,k))/rcpd
    893         ptu(i, k) = max(100., ptu(i,k))
    894         ptu(i, k) = min(400., ptu(i,k))
     879        zqeen = pqenh(i, k + 1) * pen_u(i, k)
     880        zseen = (rcpd * ptenh(i, k + 1) + pgeoh(i, k + 1)) * pen_u(i, k)
     881        zscde = (rcpd * ptu(i, k + 1) + pgeoh(i, k + 1)) * pde_u(i, k)
     882        zqude = pqu(i, k + 1) * pde_u(i, k)
     883        plude(i, k) = plu(i, k + 1) * pde_u(i, k)
     884        zmfusk = pmfus(i, k + 1) + zseen - zscde
     885        zmfuqk = pmfuq(i, k + 1) + zqeen - zqude
     886        zmfulk = pmful(i, k + 1) - plude(i, k)
     887        plu(i, k) = zmfulk * (1. / max(cmfcmin, pmfu(i, k)))
     888        pqu(i, k) = zmfuqk * (1. / max(cmfcmin, pmfu(i, k)))
     889        ptu(i, k) = (zmfusk * (1. / max(cmfcmin, pmfu(i, k))) - pgeoh(i, k)) / rcpd
     890        ptu(i, k) = max(100., ptu(i, k))
     891        ptu(i, k) = min(400., ptu(i, k))
    895892        zqold(i) = pqu(i, k)
    896893      ELSE
    … …  
    904901
    905902    iCALL = 1
    906     CALL flxadjtq(paph(1,k), ptu(1,k), pqu(1,k), llflag, icall)
    907 
    908     DO i = 1, klon
    909       IF (llflag(i) .AND. pqu(i,k)/=zqold(i)) THEN
     903    CALL flxadjtq(paph(1, k), ptu(1, k), pqu(1, k), llflag, icall)
     904
     905    DO i = 1, klon
     906      IF (llflag(i) .AND. pqu(i, k)/=zqold(i)) THEN
    910907        klab(i, k) = 2
    911908        plu(i, k) = plu(i, k) + zqold(i) - pqu(i, k)
    912         zbuo = ptu(i, k)*(1.+retv*pqu(i,k)) - ptenh(i, k)*(1.+retv*pqenh(i,k) &
    913           )
    914         IF (klab(i,k+1)==1) zbuo = zbuo + 0.5
    915         IF (zbuo>0. .AND. pmfu(i,k)>=0.1*pmfub(i)) THEN
     909        zbuo = ptu(i, k) * (1. + retv * pqu(i, k)) - ptenh(i, k) * (1. + retv * pqenh(i, k) &
     910                )
     911        IF (klab(i, k + 1)==1) zbuo = zbuo + 0.5
     912        IF (zbuo>0. .AND. pmfu(i, k)>=0.1 * pmfub(i)) THEN
    916913          kctop(i) = k
    917914          ldcum(i) = .TRUE.
    … …  
    919916          IF (ldland(i)) zdnoprc = zdland(i)
    920917          zprcon = cprcon
    921           IF ((zpbot(i)-paph(i,k))<zdnoprc) zprcon = 0.0
    922           zlnew = plu(i, k)/(1.+zprcon*(pgeoh(i,k)-pgeoh(i,k+1)))
    923           pdmfup(i, k) = max(0., (plu(i,k)-zlnew)*pmfu(i,k))
     918          IF ((zpbot(i) - paph(i, k))<zdnoprc) zprcon = 0.0
     919          zlnew = plu(i, k) / (1. + zprcon * (pgeoh(i, k) - pgeoh(i, k + 1)))
     920          pdmfup(i, k) = max(0., (plu(i, k) - zlnew) * pmfu(i, k))
    924921          plu(i, k) = zlnew
    925922        ELSE
    … …  
    931928    DO i = 1, klon
    932929      IF (llflag(i)) THEN
    933         pmful(i, k) = plu(i, k)*pmfu(i, k)
    934         pmfus(i, k) = (rcpd*ptu(i,k)+pgeoh(i,k))*pmfu(i, k)
    935         pmfuq(i, k) = pqu(i, k)*pmfu(i, k)
    936       END IF
    937     END DO
    938 
    939 480 END DO
     930        pmful(i, k) = plu(i, k) * pmfu(i, k)
     931        pmfus(i, k) = (rcpd * ptu(i, k) + pgeoh(i, k)) * pmfu(i, k)
     932        pmfuq(i, k) = pqu(i, k) * pmfu(i, k)
     933      END IF
     934    END DO
     935
     936  480 END DO
    940937  ! ----------------------------------------------------------------------
    941938  ! DETERMINE CONVECTIVE FLUXES ABOVE NON-BUOYANCY LEVEL
    … …  
    945942  ! ----------------------------------------------------------------------
    946943  DO i = 1, klon
    947     IF (kctop(i)==klev-1) ldcum(i) = .FALSE.
     944    IF (kctop(i)==klev - 1) ldcum(i) = .FALSE.
    948945    kcbot(i) = max(kcbot(i), kctop(i))
    949946  END DO
    … …  
    961958    IF (ldcum(i)) THEN
    962959      k = kctop(i) - 1
    963       pde_u(i, k) = (1.-cmfctop)*pmfu(i, k+1)
    964       plude(i, k) = pde_u(i, k)*plu(i, k+1)
    965       pmfu(i, k) = pmfu(i, k+1) - pde_u(i, k)
     960      pde_u(i, k) = (1. - cmfctop) * pmfu(i, k + 1)
     961      plude(i, k) = pde_u(i, k) * plu(i, k + 1)
     962      pmfu(i, k) = pmfu(i, k + 1) - pde_u(i, k)
    966963      zlnew = plu(i, k)
    967       pdmfup(i, k) = max(0., (plu(i,k)-zlnew)*pmfu(i,k))
     964      pdmfup(i, k) = max(0., (plu(i, k) - zlnew) * pmfu(i, k))
    968965      plu(i, k) = zlnew
    969       pmfus(i, k) = (rcpd*ptu(i,k)+pgeoh(i,k))*pmfu(i, k)
    970       pmfuq(i, k) = pqu(i, k)*pmfu(i, k)
    971       pmful(i, k) = plu(i, k)*pmfu(i, k)
    972       plude(i, k-1) = pmful(i, k)
     966      pmfus(i, k) = (rcpd * ptu(i, k) + pgeoh(i, k)) * pmfu(i, k)
     967      pmfuq(i, k) = pqu(i, k) * pmfu(i, k)
     968      pmful(i, k) = plu(i, k) * pmfu(i, k)
     969      plude(i, k - 1) = pmful(i, k)
    973970    END IF
    974971  END DO
    975972
    976 800 CONTINUE
     973  800 CONTINUE
    977974
    978975END SUBROUTINE flxasc
    979976SUBROUTINE flxflux(pdtime, pqen, pqsen, ptenh, pqenh, pap, paph, ldland, &
    980     pgeoh, kcbot, kctop, lddraf, kdtop, ktype, ldcum, pmfu, pmfd, pmfus, &
    981     pmfds, pmfuq, pmfdq, pmful, plude, pdmfup, pdmfdp, pten, prfl, psfl, &
    982     pdpmel, ktopm2, pmflxr, pmflxs)
     977        pgeoh, kcbot, kctop, lddraf, kdtop, ktype, ldcum, pmfu, pmfd, pmfus, &
     978        pmfds, pmfuq, pmfdq, pmful, plude, pdmfup, pdmfdp, pten, prfl, psfl, &
     979        pdpmel, ktopm2, pmflxr, pmflxs)
    983980  USE dimphy
    984981  USE lmdz_print_control, ONLY: prt_level
     982  USE lmdz_YOECUMF
     983
    985984  IMPLICIT NONE
    986985  ! ----------------------------------------------------------------------
    … …  
    990989  include "YOMCST.h"
    991990  include "YOETHF.h"
    992   include "YOECUMF.h"
    993991
    994992  REAL cevapcu(klon, klev)
    … …  
    996994  REAL pqen(klon, klev), pqenh(klon, klev), pqsen(klon, klev)
    997995  REAL pten(klon, klev), ptenh(klon, klev)
    998   REAL paph(klon, klev+1), pgeoh(klon, klev)
     996  REAL paph(klon, klev + 1), pgeoh(klon, klev)
    999997
    1000998  REAL pap(klon, klev)
    … …  
    10111009  REAL pdmfdp(klon, klev), maxpdmfdp(klon, klev)
    10121010  REAL prfl(klon), psfl(klon)
    1013   REAL pmflxr(klon, klev+1), pmflxs(klon, klev+1)
     1011  REAL pmflxr(klon, klev + 1), pmflxs(klon, klev + 1)
    10141012  INTEGER kcbot(klon), kctop(klon), ktype(klon)
    10151013  LOGICAL ldland(klon), ldcum(klon)
    … …  
    10281026  DO k = 1, klev
    10291027    DO i = 1, klon
    1030       cevapcu(i, k) = 1.93E-6*261.*sqrt(1.E3/(38.3*0.293)*sqrt(0.5*(paph(i,k) &
    1031         +paph(i,k+1))/paph(i,klev+1)))*0.5/rg
     1028      cevapcu(i, k) = 1.93E-6 * 261. * sqrt(1.E3 / (38.3 * 0.293) * sqrt(0.5 * (paph(i, k) &
     1029              + paph(i, k + 1)) / paph(i, klev + 1))) * 0.5 / rg
    10321030    END DO
    10331031  END DO
    … …  
    10351033  ! SPECIFY CONSTANTS
    10361034
    1037   zcons1 = rcpd/(rlmlt*rg*pdtime)
    1038   zcons2 = 1./(rg*pdtime)
     1035  zcons1 = rcpd / (rlmlt * rg * pdtime)
     1036  zcons2 = 1. / (rg * pdtime)
    10391037  zcucov = 0.05
    10401038  ztmelp2 = rtt + 2.
    … …  
    10521050  DO k = ktopm2, klev
    10531051    DO i = 1, klon
    1054       IF (ldcum(i) .AND. k>=kctop(i)-1) THEN
    1055         pmfus(i, k) = pmfus(i, k) - pmfu(i, k)*(rcpd*ptenh(i,k)+pgeoh(i,k))
    1056         pmfuq(i, k) = pmfuq(i, k) - pmfu(i, k)*pqenh(i, k)
     1052      IF (ldcum(i) .AND. k>=kctop(i) - 1) THEN
     1053        pmfus(i, k) = pmfus(i, k) - pmfu(i, k) * (rcpd * ptenh(i, k) + pgeoh(i, k))
     1054        pmfuq(i, k) = pmfuq(i, k) - pmfu(i, k) * pqenh(i, k)
    10571055        zdp = 1.5E4
    10581056        IF (ldland(i)) zdp = 3.E4
    … …  
    10621060        ! evaporee dans l'environnement)
    10631061
    1064         IF (paph(i,kcbot(i))-paph(i,kctop(i))>=zdp .AND. pqen(i,k-1)>0.8* &
    1065           pqsen(i,k-1)) pdmfup(i, k-1) = pdmfup(i, k-1) + plude(i, k-1)
     1062        IF (paph(i, kcbot(i)) - paph(i, kctop(i))>=zdp .AND. pqen(i, k - 1)>0.8 * &
     1063                pqsen(i, k - 1)) pdmfup(i, k - 1) = pdmfup(i, k - 1) + plude(i, k - 1)
    10661064
    10671065        IF (lddraf(i) .AND. k>=kdtop(i)) THEN
    1068           pmfds(i, k) = pmfds(i, k) - pmfd(i, k)*(rcpd*ptenh(i,k)+pgeoh(i,k))
    1069           pmfdq(i, k) = pmfdq(i, k) - pmfd(i, k)*pqenh(i, k)
     1066          pmfds(i, k) = pmfds(i, k) - pmfd(i, k) * (rcpd * ptenh(i, k) + pgeoh(i, k))
     1067          pmfdq(i, k) = pmfdq(i, k) - pmfd(i, k) * pqenh(i, k)
    10701068        ELSE
    10711069          pmfd(i, k) = 0.
    10721070          pmfds(i, k) = 0.
    10731071          pmfdq(i, k) = 0.
    1074           pdmfdp(i, k-1) = 0.
     1072          pdmfdp(i, k - 1) = 0.
    10751073        END IF
    10761074      ELSE
    … …  
    10791077        pmfuq(i, k) = 0.
    10801078        pmful(i, k) = 0.
    1081         pdmfup(i, k-1) = 0.
    1082         plude(i, k-1) = 0.
     1079        pdmfup(i, k - 1) = 0.
     1080        plude(i, k - 1) = 0.
    10831081        pmfd(i, k) = 0.
    10841082        pmfds(i, k) = 0.
    10851083        pmfdq(i, k) = 0.
    1086         pdmfdp(i, k-1) = 0.
     1084        pdmfdp(i, k - 1) = 0.
    10871085      END IF
    10881086    END DO
    … …  
    10931091      IF (ldcum(i) .AND. k>kcbot(i)) THEN
    10941092        ikb = kcbot(i)
    1095         zzp = ((paph(i,klev+1)-paph(i,k))/(paph(i,klev+1)-paph(i,ikb)))
     1093        zzp = ((paph(i, klev + 1) - paph(i, k)) / (paph(i, klev + 1) - paph(i, ikb)))
    10961094        IF (ktype(i)==3) zzp = zzp**2
    1097         pmfu(i, k) = pmfu(i, ikb)*zzp
    1098         pmfus(i, k) = pmfus(i, ikb)*zzp
    1099         pmfuq(i, k) = pmfuq(i, ikb)*zzp
    1100         pmful(i, k) = pmful(i, ikb)*zzp
     1095        pmfu(i, k) = pmfu(i, ikb) * zzp
     1096        pmfus(i, k) = pmfus(i, ikb) * zzp
     1097        pmfuq(i, k) = pmfuq(i, ikb) * zzp
     1098        pmful(i, k) = pmful(i, ikb) * zzp
    11011099      END IF
    11021100    END DO
    … …  
    11161114    DO i = 1, klon
    11171115      IF (ldcum(i)) THEN
    1118         IF (pmflxs(i,k)>0.0 .AND. pten(i,k)>ztmelp2) THEN
    1119           zfac = zcons1*(paph(i,k+1)-paph(i,k))
    1120           zsnmlt = min(pmflxs(i,k), zfac*(pten(i,k)-ztmelp2))
     1116        IF (pmflxs(i, k)>0.0 .AND. pten(i, k)>ztmelp2) THEN
     1117          zfac = zcons1 * (paph(i, k + 1) - paph(i, k))
     1118          zsnmlt = min(pmflxs(i, k), zfac * (pten(i, k) - ztmelp2))
    11211119          pdpmel(i, k) = zsnmlt
    1122           ztmsmlt = pten(i, k) - zsnmlt/zfac
    1123           zdelta = max(0., sign(1.,rtt-ztmsmlt))
    1124           zqsat = r2es*foeew(ztmsmlt, zdelta)/pap(i, k)
     1120          ztmsmlt = pten(i, k) - zsnmlt / zfac
     1121          zdelta = max(0., sign(1., rtt - ztmsmlt))
     1122          zqsat = r2es * foeew(ztmsmlt, zdelta) / pap(i, k)
    11251123          zqsat = min(0.5, zqsat)
    1126           zqsat = zqsat/(1.-retv*zqsat)
     1124          zqsat = zqsat / (1. - retv * zqsat)
    11271125          pqsen(i, k) = zqsat
    11281126        END IF
    1129         IF (pten(i,k)>rtt) THEN
    1130           pmflxr(i, k+1) = pmflxr(i, k) + pdmfup(i, k) + pdmfdp(i, k) + &
    1131             pdpmel(i, k)
    1132           pmflxs(i, k+1) = pmflxs(i, k) - pdpmel(i, k)
     1127        IF (pten(i, k)>rtt) THEN
     1128          pmflxr(i, k + 1) = pmflxr(i, k) + pdmfup(i, k) + pdmfdp(i, k) + &
     1129                  pdpmel(i, k)
     1130          pmflxs(i, k + 1) = pmflxs(i, k) - pdpmel(i, k)
    11331131        ELSE
    1134           pmflxs(i, k+1) = pmflxs(i, k) + pdmfup(i, k) + pdmfdp(i, k)
    1135           pmflxr(i, k+1) = pmflxr(i, k)
     1132          pmflxs(i, k + 1) = pmflxs(i, k) + pdmfup(i, k) + pdmfdp(i, k)
     1133          pmflxr(i, k + 1) = pmflxr(i, k)
    11361134        END IF
    11371135        ! si la precipitation est negative, on ajuste le plux du
    11381136        ! panache descendant pour eliminer la negativite
    1139         IF ((pmflxr(i,k+1)+pmflxs(i,k+1))<0.0) THEN
     1137        IF ((pmflxr(i, k + 1) + pmflxs(i, k + 1))<0.0) THEN
    11401138          pdmfdp(i, k) = -pmflxr(i, k) - pmflxs(i, k) - pdmfup(i, k)
    1141           pmflxr(i, k+1) = 0.0
    1142           pmflxs(i, k+1) = 0.0
     1139          pmflxr(i, k + 1) = 0.0
     1140          pmflxs(i, k + 1) = 0.0
    11431141          pdpmel(i, k) = 0.0
    11441142        END IF
    … …  
    11741172        zrfl = pmflxr(i, k) + pmflxs(i, k)
    11751173        IF (zrfl>1.0E-20) THEN
    1176           zrnew = (max(0.,sqrt(zrfl/zcucov)-cevapcu(i, &
    1177             k)*(paph(i,k+1)-paph(i,k))*max(0.,pqsen(i,k)-pqen(i,k))))**2* &
    1178             zcucov
    1179           zrmin = zrfl - zcucov*max(0., 0.8*pqsen(i,k)-pqen(i,k))*zcons2*( &
    1180             paph(i,k+1)-paph(i,k))
     1174          zrnew = (max(0., sqrt(zrfl / zcucov) - cevapcu(i, &
     1175                  k) * (paph(i, k + 1) - paph(i, k)) * max(0., pqsen(i, k) - pqen(i, k))))**2 * &
     1176                  zcucov
     1177          zrmin = zrfl - zcucov * max(0., 0.8 * pqsen(i, k) - pqen(i, k)) * zcons2 * (&
     1178                  paph(i, k + 1) - paph(i, k))
    11811179          zrnew = max(zrnew, zrmin)
    11821180          zrfln = max(zrnew, 0.)
    1183           zdrfl = min(0., zrfln-zrfl)
     1181          zdrfl = min(0., zrfln - zrfl)
    11841182          ! jq At least the amount of precipiation needed to feed the
    11851183          ! downdraft
    … …  
    11871185          ! can't
    11881186          ! jq be evaporated (surely the evaporation can't be positive):
    1189           zdrfl = max(zdrfl, min(-pmflxr(i,k)-pmflxs(i,k)-maxpdmfdp(i, &
    1190             k),0.0))
     1187          zdrfl = max(zdrfl, min(-pmflxr(i, k) - pmflxs(i, k) - maxpdmfdp(i, &
     1188                  k), 0.0))
    11911189          ! jq End of insertion
    11921190
    1193           zdenom = 1.0/max(1.0E-20, pmflxr(i,k)+pmflxs(i,k))
    1194           IF (pten(i,k)>rtt) THEN
     1191          zdenom = 1.0 / max(1.0E-20, pmflxr(i, k) + pmflxs(i, k))
     1192          IF (pten(i, k)>rtt) THEN
    11951193            zpdr = pdmfdp(i, k)
    11961194            zpds = 0.0
    … …  
    11991197            zpds = pdmfdp(i, k)
    12001198          END IF
    1201           pmflxr(i, k+1) = pmflxr(i, k) + zpdr + pdpmel(i, k) + &
    1202             zdrfl*pmflxr(i, k)*zdenom
    1203           pmflxs(i, k+1) = pmflxs(i, k) + zpds - pdpmel(i, k) + &
    1204             zdrfl*pmflxs(i, k)*zdenom
     1199          pmflxr(i, k + 1) = pmflxr(i, k) + zpdr + pdpmel(i, k) + &
     1200                  zdrfl * pmflxr(i, k) * zdenom
     1201          pmflxs(i, k + 1) = pmflxs(i, k) + zpds - pdpmel(i, k) + &
     1202                  zdrfl * pmflxs(i, k) * zdenom
    12051203          pdmfup(i, k) = pdmfup(i, k) + zdrfl
    12061204        ELSE
    1207           pmflxr(i, k+1) = 0.0
    1208           pmflxs(i, k+1) = 0.0
     1205          pmflxr(i, k + 1) = 0.0
     1206          pmflxs(i, k + 1) = 0.0
    12091207          pdmfdp(i, k) = 0.0
    12101208          pdpmel(i, k) = 0.0
    12111209        END IF
    1212         IF (pmflxr(i,k)+pmflxs(i,k)<-1.E-26 .AND. prt_level>=1) WRITE (*, *) &
    1213           'precip. < 1e-16 ', pmflxr(i, k) + pmflxs(i, k)
    1214       END IF
    1215     END DO
    1216   END DO
    1217 
    1218   DO i = 1, klon
    1219     prfl(i) = pmflxr(i, klev+1)
    1220     psfl(i) = pmflxs(i, klev+1)
    1221   END DO
    1222 
     1210        IF (pmflxr(i, k) + pmflxs(i, k)<-1.E-26 .AND. prt_level>=1) WRITE (*, *) &
     1211                'precip. < 1e-16 ', pmflxr(i, k) + pmflxs(i, k)
     1212      END IF
     1213    END DO
     1214  END DO
     1215
     1216  DO i = 1, klon
     1217    prfl(i) = pmflxr(i, klev + 1)
     1218    psfl(i) = pmflxs(i, klev + 1)
     1219  END DO
    12231220
    12241221END SUBROUTINE flxflux
    12251222SUBROUTINE flxdtdq(pdtime, ktopm2, paph, ldcum, pten, pmfus, pmfds, pmfuq, &
    1226     pmfdq, pmful, pdmfup, pdmfdp, pdpmel, dt_con, dq_con)
     1223        pmfdq, pmful, pdmfup, pdmfdp, pdpmel, dt_con, dq_con)
    12271224  USE dimphy
     1225  USE lmdz_YOECUMF
     1226
    12281227  IMPLICIT NONE
    12291228  ! ----------------------------------------------------------------------
    … …  
    12321231  include "YOMCST.h"
    12331232  include "YOETHF.h"
    1234   include "YOECUMF.h"
    12351233  ! -----------------------------------------------------------------
    12361234  LOGICAL llo1
    12371235
    1238   REAL pten(klon, klev), paph(klon, klev+1)
     1236  REAL pten(klon, klev), paph(klon, klev + 1)
    12391237  REAL pmfus(klon, klev), pmfuq(klon, klev), pmful(klon, klev)
    12401238  REAL pmfds(klon, klev), pmfdq(klon, klev)
    … …  
    12541252    DO i = 1, klon
    12551253      IF (ldcum(i)) THEN
    1256         llo1 = (pten(i,k)-rtt) > 0.
     1254        llo1 = (pten(i, k) - rtt) > 0.
    12571255        zalv = rlstt
    12581256        IF (llo1) zalv = rlvtt
    1259         zdtdt = rg/(paph(i,k+1)-paph(i,k))/rcpd*(pmfus(i,k+1)-pmfus(i,k)+ &
    1260           pmfds(i,k+1)-pmfds(i,k)-rlmlt*pdpmel(i,k)-zalv*(pmful(i, &
    1261           k+1)-pmful(i,k)-pdmfup(i,k)-pdmfdp(i,k)))
     1257        zdtdt = rg / (paph(i, k + 1) - paph(i, k)) / rcpd * (pmfus(i, k + 1) - pmfus(i, k) + &
     1258                pmfds(i, k + 1) - pmfds(i, k) - rlmlt * pdpmel(i, k) - zalv * (pmful(i, &
     1259                k + 1) - pmful(i, k) - pdmfup(i, k) - pdmfdp(i, k)))
    12621260        dt_con(i, k) = zdtdt
    1263         zdqdt = rg/(paph(i,k+1)-paph(i,k))*(pmfuq(i,k+1)-pmfuq(i,k)+pmfdq(i,k &
    1264           +1)-pmfdq(i,k)+pmful(i,k+1)-pmful(i,k)-pdmfup(i,k)-pdmfdp(i,k))
     1261        zdqdt = rg / (paph(i, k + 1) - paph(i, k)) * (pmfuq(i, k + 1) - pmfuq(i, k) + pmfdq(i, k &
     1262                + 1) - pmfdq(i, k) + pmful(i, k + 1) - pmful(i, k) - pdmfup(i, k) - pdmfdp(i, k))
    12651263        dq_con(i, k) = zdqdt
    12661264      END IF
    … …  
    12711269  DO i = 1, klon
    12721270    IF (ldcum(i)) THEN
    1273       llo1 = (pten(i,k)-rtt) > 0.
     1271      llo1 = (pten(i, k) - rtt) > 0.
    12741272      zalv = rlstt
    12751273      IF (llo1) zalv = rlvtt
    1276       zdtdt = -rg/(paph(i,k+1)-paph(i,k))/rcpd*(pmfus(i,k)+pmfds(i,k)+rlmlt* &
    1277         pdpmel(i,k)-zalv*(pmful(i,k)+pdmfup(i,k)+pdmfdp(i,k)))
     1274      zdtdt = -rg / (paph(i, k + 1) - paph(i, k)) / rcpd * (pmfus(i, k) + pmfds(i, k) + rlmlt * &
     1275              pdpmel(i, k) - zalv * (pmful(i, k) + pdmfup(i, k) + pdmfdp(i, k)))
    12781276      dt_con(i, k) = zdtdt
    1279       zdqdt = -rg/(paph(i,k+1)-paph(i,k))*(pmfuq(i,k)+pmfdq(i,k)+pmful(i,k)+ &
    1280         pdmfup(i,k)+pdmfdp(i,k))
     1277      zdqdt = -rg / (paph(i, k + 1) - paph(i, k)) * (pmfuq(i, k) + pmfdq(i, k) + pmful(i, k) + &
     1278              pdmfup(i, k) + pdmfdp(i, k))
    12811279      dq_con(i, k) = zdqdt
    12821280    END IF
    12831281  END DO
    12841282
    1285 
    12861283END SUBROUTINE flxdtdq
    12871284SUBROUTINE flxdlfs(ptenh, pqenh, pgeoh, paph, ptu, pqu, ldcum, kcbot, kctop, &
    1288     pmfub, prfl, ptd, pqd, pmfd, pmfds, pmfdq, pdmfdp, kdtop, lddraf)
     1285        pmfub, prfl, ptd, pqd, pmfd, pmfds, pmfdq, pdmfdp, kdtop, lddraf)
    12891286  USE dimphy
     1287  USE lmdz_YOECUMF
     1288
    12901289  IMPLICIT NONE
    12911290
    … …  
    13071306  include "YOMCST.h"
    13081307  include "YOETHF.h"
    1309   include "YOECUMF.h"
    13101308
    13111309  REAL ptenh(klon, klev)
    13121310  REAL pqenh(klon, klev)
    1313   REAL pgeoh(klon, klev), paph(klon, klev+1)
     1311  REAL pgeoh(klon, klev), paph(klon, klev + 1)
    13141312  REAL ptu(klon, klev), pqu(klon, klev)
    13151313  REAL pmfub(klon)
    … …  
    13541352      zqenwb(i, k) = pqenh(i, k)
    13551353      llo2(i) = ldcum(i) .AND. prfl(i) > 0. .AND. .NOT. lddraf(i) .AND. &
    1356         (k<kcbot(i) .AND. k>kctop(i))
     1354              (k<kcbot(i) .AND. k>kctop(i))
    13571355      IF (llo2(i)) is = is + 1
    13581356    END DO
    … …  
    13601358
    13611359    iCALL = 2
    1362     CALL flxadjtq(paph(1,k), ztenwb(1,k), zqenwb(1,k), llo2, icall)
     1360    CALL flxadjtq(paph(1, k), ztenwb(1, k), zqenwb(1, k), llo2, icall)
    13631361
    13641362    ! ----------------------------------------------------------------------
    … …  
    13691367    DO i = 1, klon
    13701368      IF (llo2(i)) THEN
    1371         zttest = 0.5*(ptu(i,k)+ztenwb(i,k))
    1372         zqtest = 0.5*(pqu(i,k)+zqenwb(i,k))
    1373         zbuo = zttest*(1.+retv*zqtest) - ptenh(i, k)*(1.+retv*pqenh(i,k))
     1369        zttest = 0.5 * (ptu(i, k) + ztenwb(i, k))
     1370        zqtest = 0.5 * (pqu(i, k) + zqenwb(i, k))
     1371        zbuo = zttest * (1. + retv * zqtest) - ptenh(i, k) * (1. + retv * pqenh(i, k))
    13741372        zcond(i) = pqenh(i, k) - zqenwb(i, k)
    1375         zmftop = -cmfdeps*pmfub(i)
    1376         IF (zbuo<0. .AND. prfl(i)>10.*zmftop*zcond(i)) THEN
     1373        zmftop = -cmfdeps * pmfub(i)
     1374        IF (zbuo<0. .AND. prfl(i)>10. * zmftop * zcond(i)) THEN
    13771375          kdtop(i) = k
    13781376          lddraf(i) = .TRUE.
    … …  
    13801378          pqd(i, k) = zqtest
    13811379          pmfd(i, k) = zmftop
    1382           pmfds(i, k) = pmfd(i, k)*(rcpd*ptd(i,k)+pgeoh(i,k))
    1383           pmfdq(i, k) = pmfd(i, k)*pqd(i, k)
    1384           pdmfdp(i, k-1) = -0.5*pmfd(i, k)*zcond(i)
    1385           prfl(i) = prfl(i) + pdmfdp(i, k-1)
     1380          pmfds(i, k) = pmfd(i, k) * (rcpd * ptd(i, k) + pgeoh(i, k))
     1381          pmfdq(i, k) = pmfd(i, k) * pqd(i, k)
     1382          pdmfdp(i, k - 1) = -0.5 * pmfd(i, k) * zcond(i)
     1383          prfl(i) = prfl(i) + pdmfdp(i, k - 1)
    13861384        END IF
    13871385      END IF
    13881386    END DO
    13891387
    1390 290 END DO
    1391 
     1388  290 END DO
    13921389
    13931390END SUBROUTINE flxdlfs
    13941391SUBROUTINE flxddraf(ptenh, pqenh, pgeoh, paph, prfl, ptd, pqd, pmfd, pmfds, &
    1395     pmfdq, pdmfdp, lddraf, pen_d, pde_d)
     1392        pmfdq, pdmfdp, lddraf, pen_d, pde_d)
    13961393  USE dimphy
     1394  USE lmdz_YOECUMF
     1395
    13971396  IMPLICIT NONE
    13981397
    … …  
    14141413  include "YOMCST.h"
    14151414  include "YOETHF.h"
    1416   include "YOECUMF.h"
    14171415
    14181416  REAL ptenh(klon, klev), pqenh(klon, klev)
    1419   REAL pgeoh(klon, klev), paph(klon, klev+1)
     1417  REAL pgeoh(klon, klev), paph(klon, klev + 1)
    14201418
    14211419  REAL ptd(klon, klev), pqd(klon, klev)
    … …  
    14431441    is = 0
    14441442    DO i = 1, klon
    1445       llo2(i) = lddraf(i) .AND. pmfd(i, k-1) < 0.
     1443      llo2(i) = lddraf(i) .AND. pmfd(i, k - 1) < 0.
    14461444      IF (llo2(i)) is = is + 1
    14471445    END DO
    … …  
    14501448    DO i = 1, klon
    14511449      IF (llo2(i)) THEN
    1452         zentr = entrdd*pmfd(i, k-1)*rd*ptenh(i, k-1)/(rg*paph(i,k-1))* &
    1453           (paph(i,k)-paph(i,k-1))
     1450        zentr = entrdd * pmfd(i, k - 1) * rd * ptenh(i, k - 1) / (rg * paph(i, k - 1)) * &
     1451                (paph(i, k) - paph(i, k - 1))
    14541452        pen_d(i, k) = zentr
    14551453        pde_d(i, k) = zentr
    … …  
    14621460        IF (llo2(i)) THEN
    14631461          pen_d(i, k) = 0.
    1464           pde_d(i, k) = pmfd(i, itopde)*(paph(i,k)-paph(i,k-1))/ &
    1465             (paph(i,klev+1)-paph(i,itopde))
     1462          pde_d(i, k) = pmfd(i, itopde) * (paph(i, k) - paph(i, k - 1)) / &
     1463                  (paph(i, klev + 1) - paph(i, itopde))
    14661464        END IF
    14671465      END DO
    … …  
    14701468    DO i = 1, klon
    14711469      IF (llo2(i)) THEN
    1472         pmfd(i, k) = pmfd(i, k-1) + pen_d(i, k) - pde_d(i, k)
    1473         zseen = (rcpd*ptenh(i,k-1)+pgeoh(i,k-1))*pen_d(i, k)
    1474         zqeen = pqenh(i, k-1)*pen_d(i, k)
    1475         zsdde = (rcpd*ptd(i,k-1)+pgeoh(i,k-1))*pde_d(i, k)
    1476         zqdde = pqd(i, k-1)*pde_d(i, k)
    1477         zmfdsk = pmfds(i, k-1) + zseen - zsdde
    1478         zmfdqk = pmfdq(i, k-1) + zqeen - zqdde
    1479         pqd(i, k) = zmfdqk*(1./min(-cmfcmin,pmfd(i,k)))
    1480         ptd(i, k) = (zmfdsk*(1./min(-cmfcmin,pmfd(i,k)))-pgeoh(i,k))/rcpd
    1481         ptd(i, k) = min(400., ptd(i,k))
    1482         ptd(i, k) = max(100., ptd(i,k))
     1470        pmfd(i, k) = pmfd(i, k - 1) + pen_d(i, k) - pde_d(i, k)
     1471        zseen = (rcpd * ptenh(i, k - 1) + pgeoh(i, k - 1)) * pen_d(i, k)
     1472        zqeen = pqenh(i, k - 1) * pen_d(i, k)
     1473        zsdde = (rcpd * ptd(i, k - 1) + pgeoh(i, k - 1)) * pde_d(i, k)
     1474        zqdde = pqd(i, k - 1) * pde_d(i, k)
     1475        zmfdsk = pmfds(i, k - 1) + zseen - zsdde
     1476        zmfdqk = pmfdq(i, k - 1) + zqeen - zqdde
     1477        pqd(i, k) = zmfdqk * (1. / min(-cmfcmin, pmfd(i, k)))
     1478        ptd(i, k) = (zmfdsk * (1. / min(-cmfcmin, pmfd(i, k))) - pgeoh(i, k)) / rcpd
     1479        ptd(i, k) = min(400., ptd(i, k))
     1480        ptd(i, k) = max(100., ptd(i, k))
    14831481        zcond(i) = pqd(i, k)
    14841482      END IF
    … …  
    14861484
    14871485    iCALL = 2
    1488     CALL flxadjtq(paph(1,k), ptd(1,k), pqd(1,k), llo2, icall)
     1486    CALL flxadjtq(paph(1, k), ptd(1, k), pqd(1, k), llo2, icall)
    14891487
    14901488    DO i = 1, klon
    14911489      IF (llo2(i)) THEN
    14921490        zcond(i) = zcond(i) - pqd(i, k)
    1493         zbuo = ptd(i, k)*(1.+retv*pqd(i,k)) - ptenh(i, k)*(1.+retv*pqenh(i,k) &
    1494           )
    1495         llo1 = zbuo < 0. .AND. (prfl(i)-pmfd(i,k)*zcond(i)>0.)
     1491        zbuo = ptd(i, k) * (1. + retv * pqd(i, k)) - ptenh(i, k) * (1. + retv * pqenh(i, k) &
     1492                )
     1493        llo1 = zbuo < 0. .AND. (prfl(i) - pmfd(i, k) * zcond(i)>0.)
    14961494        IF (.NOT. llo1) pmfd(i, k) = 0.0
    1497         pmfds(i, k) = (rcpd*ptd(i,k)+pgeoh(i,k))*pmfd(i, k)
    1498         pmfdq(i, k) = pqd(i, k)*pmfd(i, k)
    1499         zdmfdp = -pmfd(i, k)*zcond(i)
    1500         pdmfdp(i, k-1) = zdmfdp
     1495        pmfds(i, k) = (rcpd * ptd(i, k) + pgeoh(i, k)) * pmfd(i, k)
     1496        pmfdq(i, k) = pqd(i, k) * pmfd(i, k)
     1497        zdmfdp = -pmfd(i, k) * zcond(i)
     1498        pdmfdp(i, k - 1) = zdmfdp
    15011499        prfl(i) = prfl(i) + zdmfdp
    15021500      END IF
    15031501    END DO
    15041502
    1505 180 END DO
     1503  180 END DO
    15061504
    15071505END SUBROUTINE flxddraf
    … …  
    15301528  include "FCTTRE.h"
    15311529
    1532   z5alvcp = r5les*rlvtt/rcpd
    1533   z5alscp = r5ies*rlstt/rcpd
    1534   zalvdcp = rlvtt/rcpd
    1535   zalsdcp = rlstt/rcpd
    1536 
     1530  z5alvcp = r5les * rlvtt / rcpd
     1531  z5alscp = r5ies * rlstt / rcpd
     1532  zalvdcp = rlvtt / rcpd
     1533  zalsdcp = rlstt / rcpd
    15371534
    15381535  DO i = 1, klon
    … …  
    15421539  DO i = 1, klon
    15431540    IF (ldflag(i)) THEN
    1544       zdelta = max(0., sign(1.,rtt-pt(i)))
    1545       zcvm5 = z5alvcp*(1.-zdelta) + zdelta*z5alscp
    1546       zldcp = zalvdcp*(1.-zdelta) + zdelta*zalsdcp
    1547       zqsat = r2es*foeew(pt(i), zdelta)/pp(i)
     1541      zdelta = max(0., sign(1., rtt - pt(i)))
     1542      zcvm5 = z5alvcp * (1. - zdelta) + zdelta * z5alscp
     1543      zldcp = zalvdcp * (1. - zdelta) + zdelta * zalsdcp
     1544      zqsat = r2es * foeew(pt(i), zdelta) / pp(i)
    15481545      zqsat = min(0.5, zqsat)
    1549       zcor = 1./(1.-retv*zqsat)
    1550       zqsat = zqsat*zcor
    1551       zcond(i) = (pq(i)-zqsat)/(1.+foede(pt(i),zdelta,zcvm5,zqsat,zcor))
     1546      zcor = 1. / (1. - retv * zqsat)
     1547      zqsat = zqsat * zcor
     1548      zcond(i) = (pq(i) - zqsat) / (1. + foede(pt(i), zdelta, zcvm5, zqsat, zcor))
    15521549      IF (kcall==1) zcond(i) = max(zcond(i), 0.)
    15531550      IF (kcall==2) zcond(i) = min(zcond(i), 0.)
    1554       pt(i) = pt(i) + zldcp*zcond(i)
     1551      pt(i) = pt(i) + zldcp * zcond(i)
    15551552      pq(i) = pq(i) - zcond(i)
    15561553    END IF
    … …  
    15651562  DO i = 1, klon
    15661563    IF (ldflag(i) .AND. zcond(i)/=0.) THEN
    1567       zdelta = max(0., sign(1.,rtt-pt(i)))
    1568       zcvm5 = z5alvcp*(1.-zdelta) + zdelta*z5alscp
    1569       zldcp = zalvdcp*(1.-zdelta) + zdelta*zalsdcp
    1570       zqsat = r2es*foeew(pt(i), zdelta)/pp(i)
     1564      zdelta = max(0., sign(1., rtt - pt(i)))
     1565      zcvm5 = z5alvcp * (1. - zdelta) + zdelta * z5alscp
     1566      zldcp = zalvdcp * (1. - zdelta) + zdelta * zalsdcp
     1567      zqsat = r2es * foeew(pt(i), zdelta) / pp(i)
    15711568      zqsat = min(0.5, zqsat)
    1572       zcor = 1./(1.-retv*zqsat)
    1573       zqsat = zqsat*zcor
    1574       zcond1 = (pq(i)-zqsat)/(1.+foede(pt(i),zdelta,zcvm5,zqsat,zcor))
    1575       pt(i) = pt(i) + zldcp*zcond1
     1569      zcor = 1. / (1. - retv * zqsat)
     1570      zqsat = zqsat * zcor
     1571      zcond1 = (pq(i) - zqsat) / (1. + foede(pt(i), zdelta, zcvm5, zqsat, zcor))
     1572      pt(i) = pt(i) + zldcp * zcond1
    15761573      pq(i) = pq(i) - zcond1
    15771574    END IF
    15781575  END DO
    15791576
    1580 230 CONTINUE
     1577  230 CONTINUE
    15811578
    15821579END SUBROUTINE flxadjtq
    15831580SUBROUTINE flxsetup
     1581  USE lmdz_YOECUMF
     1582
    15841583  IMPLICIT NONE
    15851584
    15861585  ! THIS ROUTINE DEFINES DISPOSABLE PARAMETERS FOR MASSFLUX SCHEME
    1587 
    1588   include "YOECUMF.h"
    15891586
    15901587  entrpen = 1.0E-4 ! ENTRAINMENT RATE FOR PENETRATIVE CONVECTION
    … …  
    16051602  lmfdudv = .TRUE.
    16061603
    1607 
    16081604END SUBROUTINE flxsetup
Note: See TracChangeset for help on using the changeset viewer.