Ignore:
Timestamp:
Jul 26, 2024, 12:23:19 PM (4 months ago)
Author:
abarral
Message:

Correct arguments order in abort_gcm
Merge r5085 r5097 r5109 r5124 r5125 r5126 r5127
Replace calls to get_ioipsl* by IOIPSL in phylmdiso/
Symlink inlandsis and lmdz_simu_airs into phylmdiso as it's needed for the compilation
Remove now unused key from makelmdz_fcm

Location:
LMDZ6/branches/Amaury_dev
Files:
17 edited
4 copied

Legend:

Unmodified
Added
Removed

  • LMDZ6/branches/Amaury_dev

  • LMDZ6/branches/Amaury_dev/libf/phylmdiso/change_srf_frac_mod.F90

    r5112 r5132  
    33
    44MODULE change_srf_frac_mod
    5 
     5USE lmdz_abort_physic, ONLY: abort_physic
    66  IMPLICIT NONE
    77

  • LMDZ6/branches/Amaury_dev/libf/phylmdiso/cv30_routines.F90

    r5117 r5132  
    946946          )
    947947  USE lmdz_print_control, ONLY: lunout
     948  USE lmdz_abort_physic, ONLY: abort_physic
    948949#ifdef ISO
    949950    USE infotrac_phy, ONLY: ntraciso=>ntiso
    … …  
    11321133         )
    11331134    ! epmax_cape: ajout arguments
     1135USE lmdz_abort_physic, ONLY: abort_physic
    11341136#ifdef ISO
    11351137USE infotrac_phy, ONLY: ntraciso=>ntiso
    … …  
    63336335                ,cape,ep,hp,icb,inb,clw,nk,t,h,lv &
    63346336                ,epmax_diag)
     6337        USE lmdz_abort_physic, ONLY: abort_physic
    63356338        IMPLICIT NONE
    63366339

  • LMDZ6/branches/Amaury_dev/libf/phylmdiso/cv3_routines.F90

    r5117 r5132  
    12661266          )
    12671267  USE lmdz_print_control, ONLY: lunout
     1268
     1269USE lmdz_abort_physic, ONLY: abort_physic
    12681270#ifdef ISO
    12691271    USE infotrac_phy, ONLY: ntraciso=>ntiso
    … …  
    14641466         )
    14651467  USE lmdz_print_control, ONLY: prt_level
     1468  USE lmdz_abort_physic, ONLY: abort_physic
    14661469#ifdef ISO
    14671470USE infotrac_phy, ONLY: ntraciso=>ntiso
    … …  
    35843587                   )
    35853588  USE lmdz_print_control, ONLY: prt_level, lunout
     3589  USE lmdz_abort_physic, ONLY: abort_physic
    35863590#ifdef ISO
    35873591    USE infotrac_phy, ONLY: ntraciso=>ntiso

  • LMDZ6/branches/Amaury_dev/libf/phylmdiso/cv3a_compress.F90

    r5117 r5132  
    4040#endif
    4141#endif
    42 
     42  USE lmdz_abort_physic, ONLY: abort_physic
    4343  IMPLICIT NONE
    4444

  • LMDZ6/branches/Amaury_dev/libf/phylmdiso/cv_routines.F90

    r5117 r5132  
    1 
    21! $Id$
    32
    … …  
    3837  include "cvparam.h"
    3938  INTEGER nd
    40   CHARACTER (LEN=20) :: modname = 'cv_routines'
    41   CHARACTER (LEN=80) :: abort_message
     39  CHARACTER (LEN = 20) :: modname = 'cv_routines'
     40  CHARACTER (LEN = 80) :: abort_message
    4241
    4342  ! noff: integer limit for convection (nd-noff)
    … …  
    7170  delta = 0.01 ! cld
    7271
    73 
    7472END SUBROUTINE cv_param
    7573
    … …  
    9694  include "cvparam.h"
    9795
    98 
    9996  DO k = 1, nlp
    10097    DO i = 1, len
    101       lv(i, k) = lv0 - clmcpv*(t(i,k)-t0)
    102       cpn(i, k) = cpd*(1.0-q(i,k)) + cpv*q(i, k)
    103       cpx(i, k) = cpd*(1.0-q(i,k)) + cl*q(i, k)
    104       tv(i, k) = t(i, k)*(1.0+q(i,k)*epsim1)
     98      lv(i, k) = lv0 - clmcpv * (t(i, k) - t0)
     99      cpn(i, k) = cpd * (1.0 - q(i, k)) + cpv * q(i, k)
     100      cpx(i, k) = cpd * (1.0 - q(i, k)) + cl * q(i, k)
     101      tv(i, k) = t(i, k) * (1.0 + q(i, k) * epsim1)
    105102    END DO
    106103  END DO
    … …  
    113110  DO k = 2, nlp
    114111    DO i = 1, len
    115       gz(i, k) = gz(i, k-1) + hrd*(tv(i,k-1)+tv(i,k))*(p(i,k-1)-p(i,k))/ph(i, &
    116         k)
     112      gz(i, k) = gz(i, k - 1) + hrd * (tv(i, k - 1) + tv(i, k)) * (p(i, k - 1) - p(i, k)) / ph(i, &
     113              k)
    117114    END DO
    118115  END DO
    … …  
    123120  DO k = 1, nlp
    124121    DO i = 1, len
    125       h(i, k) = gz(i, k) + cpn(i, k)*t(i, k)
    126       hm(i, k) = gz(i, k) + cpx(i, k)*(t(i,k)-t(i,1)) + lv(i, k)*q(i, k)
    127     END DO
    128   END DO
    129 
     122      h(i, k) = gz(i, k) + cpn(i, k) * t(i, k)
     123      hm(i, k) = gz(i, k) + cpx(i, k) * (t(i, k) - t(i, 1)) + lv(i, k) * q(i, k)
     124    END DO
     125  END DO
    130126
    131127END SUBROUTINE cv_prelim
    132128
    133129SUBROUTINE cv_feed(len, nd, t, q, qs, p, hm, gz, nk, icb, icbmax, iflag, tnk, &
    134     qnk, gznk, plcl)
     130        qnk, gznk, plcl)
    135131  IMPLICIT NONE
    136132
    … …  
    169165  DO k = 2, nlp
    170166    DO i = 1, len
    171       IF ((hm(i,k)<work(i)) .AND. (hm(i,k)<hm(i,k-1))) THEN
     167      IF ((hm(i, k)<work(i)) .AND. (hm(i, k)<hm(i, k - 1))) THEN
    172168        work(i) = hm(i, k)
    173169        ihmin(i) = k
    … …  
    193189  DO k = minorig + 1, nl
    194190    DO i = 1, len
    195       IF ((hm(i,k)>work(i)) .AND. (k<=ihmin(i))) THEN
     191      IF ((hm(i, k)>work(i)) .AND. (k<=ihmin(i))) THEN
    196192        work(i) = hm(i, k)
    197193        nk(i) = k
    … …  
    204200  ! -------------------------------------------------------------------
    205201  DO i = 1, len
    206     IF (((t(i,nk(i))<250.0) .OR. (q(i,nk(i))<=0.0) .OR. (p(i,ihmin(i))< &
    207       400.0)) .AND. (iflag(i)==0)) iflag(i) = 7
     202    IF (((t(i, nk(i))<250.0) .OR. (q(i, nk(i))<=0.0) .OR. (p(i, ihmin(i))< &
     203            400.0)) .AND. (iflag(i)==0)) iflag(i) = 7
    208204  END DO
    209205  ! -------------------------------------------------------------------
    … …  
    218214    qsnk(i) = qs(i, nk(i))
    219215
    220     rh(i) = qnk(i)/qsnk(i)
     216    rh(i) = qnk(i) / qsnk(i)
    221217    rh(i) = min(1.0, rh(i))
    222     chi(i) = tnk(i)/(1669.0-122.0*rh(i)-tnk(i))
    223     plcl(i) = pnk(i)*(rh(i)**chi(i))
     218    chi(i) = tnk(i) / (1669.0 - 122.0 * rh(i) - tnk(i))
     219    plcl(i) = pnk(i) * (rh(i)**chi(i))
    224220    IF (((plcl(i)<200.0) .OR. (plcl(i)>=2000.0)) .AND. (iflag(i)==0)) iflag(i &
    225       ) = 8
     221            ) = 8
    226222  END DO
    227223  ! -------------------------------------------------------------------
    … …  
    234230  DO k = minorig, nl
    235231    DO i = 1, len
    236       IF ((k>=(nk(i)+1)) .AND. (p(i,k)<plcl(i))) icb(i) = min(icb(i), k)
     232      IF ((k>=(nk(i) + 1)) .AND. (p(i, k)<plcl(i))) icb(i) = min(icb(i), k)
    237233    END DO
    238234  END DO
    … …  
    249245  END DO
    250246
    251 
    252247END SUBROUTINE cv_feed
    253248
    254249SUBROUTINE cv_undilute1(len, nd, t, q, qs, gz, p, nk, icb, icbmax, tp, tvp, &
    255     clw)
     250        clw)
    256251  IMPLICIT NONE
    257252
    … …  
    292287
    293288  DO i = 1, len
    294     ah0(i) = (cpd*(1.-qnk(i))+cl*qnk(i))*tnk(i) + qnk(i)*(lv0-clmcpv*(tnk(i)- &
    295       273.15)) + gznk(i)
    296     cpp(i) = cpd*(1.-qnk(i)) + qnk(i)*cpv
     289    ah0(i) = (cpd * (1. - qnk(i)) + cl * qnk(i)) * tnk(i) + qnk(i) * (lv0 - clmcpv * (tnk(i) - &
     290            273.15)) + gznk(i)
     291    cpp(i) = cpd * (1. - qnk(i)) + qnk(i) * cpv
    297292  END DO
    298293
    … …  
    301296  DO k = minorig, icbmax - 1
    302297    DO i = 1, len
    303       tp(i, k) = tnk(i) - (gz(i,k)-gznk(i))/cpp(i)
    304       tvp(i, k) = tp(i, k)*(1.+qnk(i)*epsi)
     298      tp(i, k) = tnk(i) - (gz(i, k) - gznk(i)) / cpp(i)
     299      tvp(i, k) = tp(i, k) * (1. + qnk(i) * epsi)
    305300    END DO
    306301  END DO
    … …  
    311306    tg = ticb(i)
    312307    qg = qs(i, icb(i))
    313     alv = lv0 - clmcpv*(ticb(i)-t0)
     308    alv = lv0 - clmcpv * (ticb(i) - t0)
    314309
    315310    ! First iteration.
    316311
    317     s = cpd + alv*alv*qg/(rrv*ticb(i)*ticb(i))
    318     s = 1./s
    319     ahg = cpd*tg + (cl-cpd)*qnk(i)*ticb(i) + alv*qg + gzicb(i)
    320     tg = tg + s*(ah0(i)-ahg)
     312    s = cpd + alv * alv * qg / (rrv * ticb(i) * ticb(i))
     313    s = 1. / s
     314    ahg = cpd * tg + (cl - cpd) * qnk(i) * ticb(i) + alv * qg + gzicb(i)
     315    tg = tg + s * (ah0(i) - ahg)
    321316    tg = max(tg, 35.0)
    322317    tc = tg - t0
    323318    denom = 243.5 + tc
    324319    IF (tc>=0.0) THEN
    325       es = 6.112*exp(17.67*tc/denom)
     320      es = 6.112 * exp(17.67 * tc / denom)
    326321    ELSE
    327       es = exp(23.33086-6111.72784/tg+0.15215*log(tg))
     322      es = exp(23.33086 - 6111.72784 / tg + 0.15215 * log(tg))
    328323    END IF
    329     qg = eps*es/(p(i,icb(i))-es*(1.-eps))
     324    qg = eps * es / (p(i, icb(i)) - es * (1. - eps))
    330325
    331326    ! Second iteration.
    332327
    333     s = cpd + alv*alv*qg/(rrv*ticb(i)*ticb(i))
    334     s = 1./s
    335     ahg = cpd*tg + (cl-cpd)*qnk(i)*ticb(i) + alv*qg + gzicb(i)
    336     tg = tg + s*(ah0(i)-ahg)
     328    s = cpd + alv * alv * qg / (rrv * ticb(i) * ticb(i))
     329    s = 1. / s
     330    ahg = cpd * tg + (cl - cpd) * qnk(i) * ticb(i) + alv * qg + gzicb(i)
     331    tg = tg + s * (ah0(i) - ahg)
    337332    tg = max(tg, 35.0)
    338333    tc = tg - t0
    339334    denom = 243.5 + tc
    340335    IF (tc>=0.0) THEN
    341       es = 6.112*exp(17.67*tc/denom)
     336      es = 6.112 * exp(17.67 * tc / denom)
    342337    ELSE
    343       es = exp(23.33086-6111.72784/tg+0.15215*log(tg))
     338      es = exp(23.33086 - 6111.72784 / tg + 0.15215 * log(tg))
    344339    END IF
    345     qg = eps*es/(p(i,icb(i))-es*(1.-eps))
    346 
    347     alv = lv0 - clmcpv*(ticb(i)-273.15)
    348     tp(i, icb(i)) = (ah0(i)-(cl-cpd)*qnk(i)*ticb(i)-gz(i,icb(i))-alv*qg)/cpd
     340    qg = eps * es / (p(i, icb(i)) - es * (1. - eps))
     341
     342    alv = lv0 - clmcpv * (ticb(i) - 273.15)
     343    tp(i, icb(i)) = (ah0(i) - (cl - cpd) * qnk(i) * ticb(i) - gz(i, icb(i)) - alv * qg) / cpd
    349344    clw(i, icb(i)) = qnk(i) - qg
    350     clw(i, icb(i)) = max(0.0, clw(i,icb(i)))
    351     rg = qg/(1.-qnk(i))
    352     tvp(i, icb(i)) = tp(i, icb(i))*(1.+rg*epsi)
     345    clw(i, icb(i)) = max(0.0, clw(i, icb(i)))
     346    rg = qg / (1. - qnk(i))
     347    tvp(i, icb(i)) = tp(i, icb(i)) * (1. + rg * epsi)
    353348  END DO
    354349
    355350  DO k = minorig, icbmax
    356351    DO i = 1, len
    357       tvp(i, k) = tvp(i, k) - tp(i, k)*qnk(i)
    358     END DO
    359   END DO
    360 
     352      tvp(i, k) = tvp(i, k) - tp(i, k) * qnk(i)
     353    END DO
     354  END DO
    361355
    362356END SUBROUTINE cv_undilute1
    … …  
    383377  INTEGER i
    384378
    385 
    386379  DO i = 1, len
    387380    IF ((cbmf(i)==0.0) .AND. (iflag(i)==0) .AND. (tvp(i, &
    388       icb(i))<=(tv(i,icb(i))-dtmax))) iflag(i) = 4
    389   END DO
    390 
     381            icb(i))<=(tv(i, icb(i)) - dtmax))) iflag(i) = 4
     382  END DO
    391383
    392384END SUBROUTINE cv_trigger
    393385
    394386SUBROUTINE cv_compress(len, nloc, ncum, nd, iflag1, nk1, icb1, cbmf1, plcl1, &
    395     tnk1, qnk1, gznk1, t1, q1, qs1, u1, v1, gz1, h1, lv1, cpn1, p1, ph1, tv1, &
    396     tp1, tvp1, clw1, iflag, nk, icb, cbmf, plcl, tnk, qnk, gznk, t, q, qs, u, &
    397     v, gz, h, lv, cpn, p, ph, tv, tp, tvp, clw, dph)
     387        tnk1, qnk1, gznk1, t1, q1, qs1, u1, v1, gz1, h1, lv1, cpn1, p1, ph1, tv1, &
     388        tp1, tvp1, clw1, iflag, nk, icb, cbmf, plcl, tnk, qnk, gznk, t, q, qs, u, &
     389        v, gz, h, lv, cpn, p, ph, tv, tp, tvp, clw, dph)
    398390  USE lmdz_print_control, ONLY: lunout
     391  USE lmdz_abort_physic, ONLY: abort_physic
    399392  IMPLICIT NONE
    400393
    … …  
    407400  REAL t1(len, nd), q1(len, nd), qs1(len, nd), u1(len, nd), v1(len, nd)
    408401  REAL gz1(len, nd), h1(len, nd), lv1(len, nd), cpn1(len, nd)
    409   REAL p1(len, nd), ph1(len, nd+1), tv1(len, nd), tp1(len, nd)
     402  REAL p1(len, nd), ph1(len, nd + 1), tv1(len, nd), tp1(len, nd)
    410403  REAL tvp1(len, nd), clw1(len, nd)
    411404
    … …  
    415408  REAL t(nloc, nd), q(nloc, nd), qs(nloc, nd), u(nloc, nd), v(nloc, nd)
    416409  REAL gz(nloc, nd), h(nloc, nd), lv(nloc, nd), cpn(nloc, nd)
    417   REAL p(nloc, nd), ph(nloc, nd+1), tv(nloc, nd), tp(nloc, nd)
     410  REAL p(nloc, nd), ph(nloc, nd + 1), tv(nloc, nd), tp(nloc, nd)
    418411  REAL tvp(nloc, nd), clw(nloc, nd)
    419412  REAL dph(nloc, nd)
    … …  
    421414  ! local variables:
    422415  INTEGER i, k, nn
    423   CHARACTER (LEN=20) :: modname = 'cv_compress'
    424   CHARACTER (LEN=80) :: abort_message
    425 
     416  CHARACTER (LEN = 20) :: modname = 'cv_compress'
     417  CHARACTER (LEN = 80) :: abort_message
    426418
    427419  DO k = 1, nl + 1
    … …  
    472464  DO k = 1, nl
    473465    DO i = 1, ncum
    474       dph(i, k) = ph(i, k) - ph(i, k+1)
    475     END DO
    476   END DO
    477 
     466      dph(i, k) = ph(i, k) - ph(i, k + 1)
     467    END DO
     468  END DO
    478469
    479470END SUBROUTINE cv_compress
    480471
    481472SUBROUTINE cv_undilute2(nloc, ncum, nd, icb, nk, tnk, qnk, gznk, t, q, qs, &
    482     gz, p, dph, h, tv, lv, inb, inb1, tp, tvp, clw, hp, ep, sigp, frac)
     473        gz, p, dph, h, tv, lv, inb, inb1, tp, tvp, clw, hp, ep, sigp, frac)
    483474  IMPLICIT NONE
    484475
    … …  
    537528  ! ***  Calculate certain parcel quantities, including static energy   ***
    538529
    539 
    540   DO i = 1, ncum
    541     ah0(i) = (cpd*(1.-qnk(i))+cl*qnk(i))*tnk(i) + qnk(i)*(lv0-clmcpv*(tnk(i)- &
    542       t0)) + gznk(i)
     530  DO i = 1, ncum
     531    ah0(i) = (cpd * (1. - qnk(i)) + cl * qnk(i)) * tnk(i) + qnk(i) * (lv0 - clmcpv * (tnk(i) - &
     532            t0)) + gznk(i)
    543533  END DO
    544534
    … …  
    546536  ! ***  Find lifted parcel quantities above cloud base    ***
    547537
    548 
    549538  DO k = minorig + 1, nl
    550539    DO i = 1, ncum
    551       IF (k>=(icb(i)+1)) THEN
     540      IF (k>=(icb(i) + 1)) THEN
    552541        tg = t(i, k)
    553542        qg = qs(i, k)
    554         alv = lv0 - clmcpv*(t(i,k)-t0)
     543        alv = lv0 - clmcpv * (t(i, k) - t0)
    555544
    556545        ! First iteration.
    557546
    558         s = cpd + alv*alv*qg/(rrv*t(i,k)*t(i,k))
    559         s = 1./s
    560         ahg = cpd*tg + (cl-cpd)*qnk(i)*t(i, k) + alv*qg + gz(i, k)
    561         tg = tg + s*(ah0(i)-ahg)
     547        s = cpd + alv * alv * qg / (rrv * t(i, k) * t(i, k))
     548        s = 1. / s
     549        ahg = cpd * tg + (cl - cpd) * qnk(i) * t(i, k) + alv * qg + gz(i, k)
     550        tg = tg + s * (ah0(i) - ahg)
    562551        tg = max(tg, 35.0)
    563552        tc = tg - t0
    564553        denom = 243.5 + tc
    565554        IF (tc>=0.0) THEN
    566           es = 6.112*exp(17.67*tc/denom)
     555          es = 6.112 * exp(17.67 * tc / denom)
    567556        ELSE
    568           es = exp(23.33086-6111.72784/tg+0.15215*log(tg))
     557          es = exp(23.33086 - 6111.72784 / tg + 0.15215 * log(tg))
    569558        END IF
    570         qg = eps*es/(p(i,k)-es*(1.-eps))
     559        qg = eps * es / (p(i, k) - es * (1. - eps))
    571560
    572561        ! Second iteration.
    573562
    574         s = cpd + alv*alv*qg/(rrv*t(i,k)*t(i,k))
    575         s = 1./s
    576         ahg = cpd*tg + (cl-cpd)*qnk(i)*t(i, k) + alv*qg + gz(i, k)
    577         tg = tg + s*(ah0(i)-ahg)
     563        s = cpd + alv * alv * qg / (rrv * t(i, k) * t(i, k))
     564        s = 1. / s
     565        ahg = cpd * tg + (cl - cpd) * qnk(i) * t(i, k) + alv * qg + gz(i, k)
     566        tg = tg + s * (ah0(i) - ahg)
    578567        tg = max(tg, 35.0)
    579568        tc = tg - t0
    580569        denom = 243.5 + tc
    581570        IF (tc>=0.0) THEN
    582           es = 6.112*exp(17.67*tc/denom)
     571          es = 6.112 * exp(17.67 * tc / denom)
    583572        ELSE
    584           es = exp(23.33086-6111.72784/tg+0.15215*log(tg))
     573          es = exp(23.33086 - 6111.72784 / tg + 0.15215 * log(tg))
    585574        END IF
    586         qg = eps*es/(p(i,k)-es*(1.-eps))
    587 
    588         alv = lv0 - clmcpv*(t(i,k)-t0)
     575        qg = eps * es / (p(i, k) - es * (1. - eps))
     576
     577        alv = lv0 - clmcpv * (t(i, k) - t0)
    589578        ! PRINT*,'cpd dans convect2 ',cpd
    590579        ! PRINT*,'tp(i,k),ah0(i),cl,cpd,qnk(i),t(i,k),gz(i,k),alv,qg,cpd'
    591580        ! PRINT*,tp(i,k),ah0(i),cl,cpd,qnk(i),t(i,k),gz(i,k),alv,qg,cpd
    592         tp(i, k) = (ah0(i)-(cl-cpd)*qnk(i)*t(i,k)-gz(i,k)-alv*qg)/cpd
     581        tp(i, k) = (ah0(i) - (cl - cpd) * qnk(i) * t(i, k) - gz(i, k) - alv * qg) / cpd
    593582        ! if (.NOT.cpd.gt.1000.) THEN
    594583        ! PRINT*,'CPD=',cpd
    … …  
    596585        ! END IF
    597586        clw(i, k) = qnk(i) - qg
    598         clw(i, k) = max(0.0, clw(i,k))
    599         rg = qg/(1.-qnk(i))
    600         tvp(i, k) = tp(i, k)*(1.+rg*epsi)
     587        clw(i, k) = max(0.0, clw(i, k))
     588        rg = qg / (1. - qnk(i))
     589        tvp(i, k) = tp(i, k) * (1. + rg * epsi)
    601590      END IF
    602591    END DO
    … …  
    611600  DO k = minorig + 1, nl
    612601    DO i = 1, ncum
    613       IF (k>=(nk(i)+1)) THEN
     602      IF (k>=(nk(i) + 1)) THEN
    614603        tca = tp(i, k) - t0
    615604        IF (tca>=0.0) THEN
    616605          elacrit = elcrit
    617606        ELSE
    618           elacrit = elcrit*(1.0-tca/tlcrit)
     607          elacrit = elcrit * (1.0 - tca / tlcrit)
    619608        END IF
    620609        elacrit = max(elacrit, 0.0)
    621         ep(i, k) = 1.0 - elacrit/max(clw(i,k), 1.0E-8)
    622         ep(i, k) = max(ep(i,k), 0.0)
    623         ep(i, k) = min(ep(i,k), 1.0)
     610        ep(i, k) = 1.0 - elacrit / max(clw(i, k), 1.0E-8)
     611        ep(i, k) = max(ep(i, k), 0.0)
     612        ep(i, k) = min(ep(i, k), 1.0)
    624613        sigp(i, k) = sigs
    625614      END IF
    … …  
    634623  DO k = minorig + 1, nl
    635624    DO i = 1, ncum
    636       IF (k>=(icb(i)+1)) THEN
    637         tvp(i, k) = tvp(i, k)*(1.0-qnk(i)+ep(i,k)*clw(i,k))
     625      IF (k>=(icb(i) + 1)) THEN
     626        tvp(i, k) = tvp(i, k) * (1.0 - qnk(i) + ep(i, k) * clw(i, k))
    638627        ! PRINT*,'i,k,tvp(i,k),qnk(i),ep(i,k),clw(i,k)'
    639628        ! PRINT*, i,k,tvp(i,k),qnk(i),ep(i,k),clw(i,k)
    … …  
    642631  END DO
    643632  DO i = 1, ncum
    644     tvp(i, nlp) = tvp(i, nl) - (gz(i,nlp)-gz(i,nl))/cpd
     633    tvp(i, nlp) = tvp(i, nl) - (gz(i, nlp) - gz(i, nl)) / cpd
    645634  END DO
    646635
    … …  
    720709    DO i = 1, ncum
    721710      IF (cape(i)<0.0) lcape(i) = .FALSE.
    722       IF ((k>=(icb(i)+1)) .AND. lcape(i)) THEN
    723         by = (tvp(i,k)-tv(i,k))*dph(i, k)/p(i, k)
    724         byp(i) = (tvp(i,k+1)-tv(i,k+1))*dph(i, k+1)/p(i, k+1)
     711      IF ((k>=(icb(i) + 1)) .AND. lcape(i)) THEN
     712        by = (tvp(i, k) - tv(i, k)) * dph(i, k) / p(i, k)
     713        byp(i) = (tvp(i, k + 1) - tv(i, k + 1)) * dph(i, k + 1) / p(i, k + 1)
    725714        cape(i) = cape(i) + by
    726715        IF (by>=0.0) inb1(i) = k + 1
    … …  
    736725    defrac = capem(i) - cape(i)
    737726    defrac = max(defrac, 0.001)
    738     frac(i) = -cape(i)/defrac
     727    frac(i) = -cape(i) / defrac
    739728    frac(i) = min(frac(i), 1.0)
    740729    frac(i) = max(frac(i), 0.0)
    … …  
    746735
    747736  ! initialization:
    748   DO i = 1, ncum*nlp
     737  DO i = 1, ncum * nlp
    749738    hp(i, 1) = h(i, 1)
    750739  END DO
    … …  
    753742    DO i = 1, ncum
    754743      IF ((k>=icb(i)) .AND. (k<=inb(i))) THEN
    755         hp(i, k) = h(i, nk(i)) + (lv(i,k)+(cpd-cpv)*t(i,k))*ep(i, k)*clw(i, k &
    756           )
    757       END IF
    758     END DO
    759   END DO
    760 
     744        hp(i, k) = h(i, nk(i)) + (lv(i, k) + (cpd - cpv) * t(i, k)) * ep(i, k) * clw(i, k &
     745                )
     746      END IF
     747    END DO
     748  END DO
    761749
    762750END SUBROUTINE cv_undilute2
    763751
    764752SUBROUTINE cv_closure(nloc, ncum, nd, nk, icb, tv, tvp, p, ph, dph, plcl, &
    765     cpn, iflag, cbmf)
     753        cpn, iflag, cbmf)
    766754  IMPLICIT NONE
    767755
    … …  
    770758  INTEGER nk(nloc), icb(nloc)
    771759  REAL tv(nloc, nd), tvp(nloc, nd), p(nloc, nd), dph(nloc, nd)
    772   REAL ph(nloc, nd+1) ! caution nd instead ndp1 to be consistent...
     760  REAL ph(nloc, nd + 1) ! caution nd instead ndp1 to be consistent...
    773761  REAL plcl(nloc), cpn(nloc, nd)
    774762
    … …  
    804792  DO i = 1, ncum
    805793    dtpbl(i) = 0.0
    806     tvpplcl(i) = tvp(i, icb(i)-1) - rrd*tvp(i, icb(i)-1)*(p(i,icb(i)-1)-plcl( &
    807       i))/(cpn(i,icb(i)-1)*p(i,icb(i)-1))
    808     tvaplcl(i) = tv(i, icb(i)) + (tvp(i,icb(i))-tvp(i,icb(i)+1))*(plcl(i)-p(i &
    809       ,icb(i)))/(p(i,icb(i))-p(i,icb(i)+1))
     794    tvpplcl(i) = tvp(i, icb(i) - 1) - rrd * tvp(i, icb(i) - 1) * (p(i, icb(i) - 1) - plcl(&
     795            i)) / (cpn(i, icb(i) - 1) * p(i, icb(i) - 1))
     796    tvaplcl(i) = tv(i, icb(i)) + (tvp(i, icb(i)) - tvp(i, icb(i) + 1)) * (plcl(i) - p(i &
     797            , icb(i))) / (p(i, icb(i)) - p(i, icb(i) + 1))
    810798  END DO
    811799
    … …  
    819807  DO k = minorig, icbmax
    820808    DO i = 1, ncum
    821       IF ((k>=nk(i)) .AND. (k<=(icb(i)-1))) THEN
    822         dtpbl(i) = dtpbl(i) + (tvp(i,k)-tv(i,k))*dph(i, k)
    823       END IF
    824     END DO
    825   END DO
    826   DO i = 1, ncum
    827     dtpbl(i) = dtpbl(i)/(ph(i,nk(i))-ph(i,icb(i)))
     809      IF ((k>=nk(i)) .AND. (k<=(icb(i) - 1))) THEN
     810        dtpbl(i) = dtpbl(i) + (tvp(i, k) - tv(i, k)) * dph(i, k)
     811      END IF
     812    END DO
     813  END DO
     814  DO i = 1, ncum
     815    dtpbl(i) = dtpbl(i) / (ph(i, nk(i)) - ph(i, icb(i)))
    828816    dtmin(i) = tvpplcl(i) - tvaplcl(i) + dtmax + dtpbl(i)
    829817  END DO
    … …  
    835823  DO i = 1, ncum
    836824    work(i) = cbmf(i)
    837     cbmf(i) = max(0.0, (1.0-damp)*cbmf(i)+0.1*alpha*dtmin(i))
     825    cbmf(i) = max(0.0, (1.0 - damp) * cbmf(i) + 0.1 * alpha * dtmin(i))
    838826    IF ((work(i)==0.0) .AND. (cbmf(i)==0.0)) THEN
    839827      iflag(i) = 3
    … …  
    841829  END DO
    842830
    843 
    844831END SUBROUTINE cv_closure
    845832
    846833SUBROUTINE cv_mixing(nloc, ncum, nd, icb, nk, inb, inb1, ph, t, q, qs, u, v, &
    847     h, lv, qnk, hp, tv, tvp, ep, clw, cbmf, m, ment, qent, uent, vent, nent, &
    848     sij, elij)
     834        h, lv, qnk, hp, tv, tvp, ep, clw, cbmf, m, ment, qent, uent, vent, nent, &
     835        sij, elij)
    849836  IMPLICIT NONE
    850837
    … …  
    856843  INTEGER icb(nloc), inb(nloc), inb1(nloc), nk(nloc)
    857844  REAL cbmf(nloc), qnk(nloc)
    858   REAL ph(nloc, nd+1)
     845  REAL ph(nloc, nd + 1)
    859846  REAL t(nloc, nd), q(nloc, nd), qs(nloc, nd), lv(nloc, nd)
    860847  REAL u(nloc, nd), v(nloc, nd), h(nloc, nd), hp(nloc, nd)
    … …  
    880867  ! =====================================================================
    881868
    882   DO i = 1, ncum*nlp
     869  DO i = 1, ncum * nlp
    883870    nent(i, 1) = 0
    884871    m(i, 1) = 0.0
    … …  
    906893  DO j = minorig + 1, nl
    907894    DO i = 1, ncum
    908       IF ((j>=(icb(i)+1)) .AND. (j<=inb(i))) THEN
     895      IF ((j>=(icb(i) + 1)) .AND. (j<=inb(i))) THEN
    909896        k = min(j, inb1(i))
    910         dbo = abs(tv(i,k+1)-tvp(i,k+1)-tv(i,k-1)+tvp(i,k-1)) + &
    911           entp*0.04*(ph(i,k)-ph(i,k+1))
     897        dbo = abs(tv(i, k + 1) - tvp(i, k + 1) - tv(i, k - 1) + tvp(i, k - 1)) + &
     898                entp * 0.04 * (ph(i, k) - ph(i, k + 1))
    912899        work(i) = work(i) + dbo
    913         m(i, j) = cbmf(i)*dbo
     900        m(i, j) = cbmf(i) * dbo
    914901      END IF
    915902    END DO
    … …  
    917904  DO k = minorig + 1, nl
    918905    DO i = 1, ncum
    919       IF ((k>=(icb(i)+1)) .AND. (k<=inb(i))) THEN
    920         m(i, k) = m(i, k)/work(i)
     906      IF ((k>=(icb(i) + 1)) .AND. (k<=inb(i))) THEN
     907        m(i, k) = m(i, k) / work(i)
    921908      END IF
    922909    END DO
    … …  
    930917  ! =====================================================================
    931918
    932 
    933919  DO i = minorig + 1, nl
    934920    DO j = minorig + 1, nl
    935921      DO ij = 1, ncum
    936         IF ((i>=(icb(ij)+1)) .AND. (j>=icb(ij)) .AND. (i<=inb(ij)) .AND. (j<= &
    937             inb(ij))) THEN
    938           qti = qnk(ij) - ep(ij, i)*clw(ij, i)
    939           bf2 = 1. + lv(ij, j)*lv(ij, j)*qs(ij, j)/(rrv*t(ij,j)*t(ij,j)*cpd)
    940           anum = h(ij, j) - hp(ij, i) + (cpv-cpd)*t(ij, j)*(qti-q(ij,j))
    941           denom = h(ij, i) - hp(ij, i) + (cpd-cpv)*(q(ij,i)-qti)*t(ij, j)
     922        IF ((i>=(icb(ij) + 1)) .AND. (j>=icb(ij)) .AND. (i<=inb(ij)) .AND. (j<= &
     923                inb(ij))) THEN
     924          qti = qnk(ij) - ep(ij, i) * clw(ij, i)
     925          bf2 = 1. + lv(ij, j) * lv(ij, j) * qs(ij, j) / (rrv * t(ij, j) * t(ij, j) * cpd)
     926          anum = h(ij, j) - hp(ij, i) + (cpv - cpd) * t(ij, j) * (qti - q(ij, j))
     927          denom = h(ij, i) - hp(ij, i) + (cpd - cpv) * (q(ij, i) - qti) * t(ij, j)
    942928          dei = denom
    943929          IF (abs(dei)<0.01) dei = 0.01
    944           sij(ij, i, j) = anum/dei
     930          sij(ij, i, j) = anum / dei
    945931          sij(ij, i, i) = 1.0
    946           altem = sij(ij, i, j)*q(ij, i) + (1.-sij(ij,i,j))*qti - qs(ij, j)
    947           altem = altem/bf2
    948           cwat = clw(ij, j)*(1.-ep(ij,j))
     932          altem = sij(ij, i, j) * q(ij, i) + (1. - sij(ij, i, j)) * qti - qs(ij, j)
     933          altem = altem / bf2
     934          cwat = clw(ij, j) * (1. - ep(ij, j))
    949935          stemp = sij(ij, i, j)
    950936          IF ((stemp<0.0 .OR. stemp>1.0 .OR. altem>cwat) .AND. j>i) THEN
    951             anum = anum - lv(ij, j)*(qti-qs(ij,j)-cwat*bf2)
    952             denom = denom + lv(ij, j)*(q(ij,i)-qti)
     937            anum = anum - lv(ij, j) * (qti - qs(ij, j) - cwat * bf2)
     938            denom = denom + lv(ij, j) * (q(ij, i) - qti)
    953939            IF (abs(denom)<0.01) denom = 0.01
    954             sij(ij, i, j) = anum/denom
    955             altem = sij(ij, i, j)*q(ij, i) + (1.-sij(ij,i,j))*qti - qs(ij, j)
    956             altem = altem - (bf2-1.)*cwat
     940            sij(ij, i, j) = anum / denom
     941            altem = sij(ij, i, j) * q(ij, i) + (1. - sij(ij, i, j)) * qti - qs(ij, j)
     942            altem = altem - (bf2 - 1.) * cwat
    957943          END IF
    958           IF (sij(ij,i,j)>0.0 .AND. sij(ij,i,j)<0.9) THEN
    959             qent(ij, i, j) = sij(ij, i, j)*q(ij, i) + (1.-sij(ij,i,j))*qti
    960             uent(ij, i, j) = sij(ij, i, j)*u(ij, i) + &
    961               (1.-sij(ij,i,j))*u(ij, nk(ij))
    962             vent(ij, i, j) = sij(ij, i, j)*v(ij, i) + &
    963               (1.-sij(ij,i,j))*v(ij, nk(ij))
     944          IF (sij(ij, i, j)>0.0 .AND. sij(ij, i, j)<0.9) THEN
     945            qent(ij, i, j) = sij(ij, i, j) * q(ij, i) + (1. - sij(ij, i, j)) * qti
     946            uent(ij, i, j) = sij(ij, i, j) * u(ij, i) + &
     947                    (1. - sij(ij, i, j)) * u(ij, nk(ij))
     948            vent(ij, i, j) = sij(ij, i, j) * v(ij, i) + &
     949                    (1. - sij(ij, i, j)) * v(ij, nk(ij))
    964950            elij(ij, i, j) = altem
    965             elij(ij, i, j) = max(0.0, elij(ij,i,j))
    966             ment(ij, i, j) = m(ij, i)/(1.-sij(ij,i,j))
     951            elij(ij, i, j) = max(0.0, elij(ij, i, j))
     952            ment(ij, i, j) = m(ij, i) / (1. - sij(ij, i, j))
    967953            nent(ij, i) = nent(ij, i) + 1
    968954          END IF
    969           sij(ij, i, j) = max(0.0, sij(ij,i,j))
    970           sij(ij, i, j) = min(1.0, sij(ij,i,j))
     955          sij(ij, i, j) = max(0.0, sij(ij, i, j))
     956          sij(ij, i, j) = min(1.0, sij(ij, i, j))
    971957        END IF
    972958      END DO
    … …  
    979965
    980966    DO ij = 1, ncum
    981       IF ((i>=(icb(ij)+1)) .AND. (i<=inb(ij)) .AND. (nent(ij,i)==0)) THEN
     967      IF ((i>=(icb(ij) + 1)) .AND. (i<=inb(ij)) .AND. (nent(ij, i)==0)) THEN
    982968        ment(ij, i, i) = m(ij, i)
    983         qent(ij, i, i) = q(ij, nk(ij)) - ep(ij, i)*clw(ij, i)
     969        qent(ij, i, i) = q(ij, nk(ij)) - ep(ij, i) * clw(ij, i)
    984970        uent(ij, i, i) = u(ij, nk(ij))
    985971        vent(ij, i, i) = v(ij, nk(ij))
    … …  
    999985  ! =====================================================================
    1000986
    1001   CALL zilch(bsum, ncum*nlp)
     987  CALL zilch(bsum, ncum * nlp)
    1002988  DO ij = 1, ncum
    1003989    lwork(ij) = .FALSE.
    … …  
    1007993    num1 = 0
    1008994    DO ij = 1, ncum
    1009       IF ((i>=icb(ij)+1) .AND. (i<=inb(ij))) num1 = num1 + 1
     995      IF ((i>=icb(ij) + 1) .AND. (i<=inb(ij))) num1 = num1 + 1
    1010996    END DO
    1011997    IF (num1<=0) GO TO 789
    1012998
    1013999    DO ij = 1, ncum
    1014       IF ((i>=icb(ij)+1) .AND. (i<=inb(ij))) THEN
    1015         lwork(ij) = (nent(ij,i)/=0)
    1016         qp1 = q(ij, nk(ij)) - ep(ij, i)*clw(ij, i)
    1017         anum = h(ij, i) - hp(ij, i) - lv(ij, i)*(qp1-qs(ij,i))
    1018         denom = h(ij, i) - hp(ij, i) + lv(ij, i)*(q(ij,i)-qp1)
     1000      IF ((i>=icb(ij) + 1) .AND. (i<=inb(ij))) THEN
     1001        lwork(ij) = (nent(ij, i)/=0)
     1002        qp1 = q(ij, nk(ij)) - ep(ij, i) * clw(ij, i)
     1003        anum = h(ij, i) - hp(ij, i) - lv(ij, i) * (qp1 - qs(ij, i))
     1004        denom = h(ij, i) - hp(ij, i) + lv(ij, i) * (q(ij, i) - qp1)
    10191005        IF (abs(denom)<0.01) denom = 0.01
    1020         scrit(ij) = anum/denom
    1021         alt = qp1 - qs(ij, i) + scrit(ij)*(q(ij,i)-qp1)
     1006        scrit(ij) = anum / denom
     1007        alt = qp1 - qs(ij, i) + scrit(ij) * (q(ij, i) - qp1)
    10221008        IF (scrit(ij)<0.0 .OR. alt<0.0) scrit(ij) = 1.0
    10231009        asij(ij) = 0.0
    … …  
    10291015      num2 = 0
    10301016      DO ij = 1, ncum
    1031         IF ((i>=icb(ij)+1) .AND. (i<=inb(ij)) .AND. (j>=icb( &
    1032           ij)) .AND. (j<=inb(ij)) .AND. lwork(ij)) num2 = num2 + 1
     1017        IF ((i>=icb(ij) + 1) .AND. (i<=inb(ij)) .AND. (j>=icb(&
     1018                ij)) .AND. (j<=inb(ij)) .AND. lwork(ij)) num2 = num2 + 1
    10331019      END DO
    10341020      IF (num2<=0) GO TO 783
    10351021
    10361022      DO ij = 1, ncum
    1037         IF ((i>=icb(ij)+1) .AND. (i<=inb(ij)) .AND. (j>=icb( &
    1038             ij)) .AND. (j<=inb(ij)) .AND. lwork(ij)) THEN
    1039           IF (sij(ij,i,j)>0.0 .AND. sij(ij,i,j)<0.9) THEN
     1023        IF ((i>=icb(ij) + 1) .AND. (i<=inb(ij)) .AND. (j>=icb(&
     1024                ij)) .AND. (j<=inb(ij)) .AND. lwork(ij)) THEN
     1025          IF (sij(ij, i, j)>0.0 .AND. sij(ij, i, j)<0.9) THEN
    10401026            IF (j>i) THEN
    1041               smid = min(sij(ij,i,j), scrit(ij))
     1027              smid = min(sij(ij, i, j), scrit(ij))
    10421028              sjmax = smid
    10431029              sjmin = smid
    1044               IF (smid<smin(ij) .AND. sij(ij,i,j+1)<smid) THEN
     1030              IF (smid<smin(ij) .AND. sij(ij, i, j + 1)<smid) THEN
    10451031                smin(ij) = smid
    1046                 sjmax = min(sij(ij,i,j+1), sij(ij,i,j), scrit(ij))
    1047                 sjmin = max(sij(ij,i,j-1), sij(ij,i,j))
     1032                sjmax = min(sij(ij, i, j + 1), sij(ij, i, j), scrit(ij))
     1033                sjmin = max(sij(ij, i, j - 1), sij(ij, i, j))
    10481034                sjmin = min(sjmin, scrit(ij))
    10491035              END IF
    10501036            ELSE
    1051               sjmax = max(sij(ij,i,j+1), scrit(ij))
    1052               smid = max(sij(ij,i,j), scrit(ij))
     1037              sjmax = max(sij(ij, i, j + 1), scrit(ij))
     1038              smid = max(sij(ij, i, j), scrit(ij))
    10531039              sjmin = 0.0
    1054               IF (j>1) sjmin = sij(ij, i, j-1)
     1040              IF (j>1) sjmin = sij(ij, i, j - 1)
    10551041              sjmin = max(sjmin, scrit(ij))
    10561042            END IF
    1057             delp = abs(sjmax-smid)
    1058             delm = abs(sjmin-smid)
    1059             asij(ij) = asij(ij) + (delp+delm)*(ph(ij,j)-ph(ij,j+1))
    1060             ment(ij, i, j) = ment(ij, i, j)*(delp+delm)*(ph(ij,j)-ph(ij,j+1))
     1043            delp = abs(sjmax - smid)
     1044            delm = abs(sjmin - smid)
     1045            asij(ij) = asij(ij) + (delp + delm) * (ph(ij, j) - ph(ij, j + 1))
     1046            ment(ij, i, j) = ment(ij, i, j) * (delp + delm) * (ph(ij, j) - ph(ij, j + 1))
    10611047          END IF
    10621048        END IF
    10631049      END DO
    1064 783 END DO
     1050    783 END DO
    10651051    DO ij = 1, ncum
    1066       IF ((i>=icb(ij)+1) .AND. (i<=inb(ij)) .AND. lwork(ij)) THEN
     1052      IF ((i>=icb(ij) + 1) .AND. (i<=inb(ij)) .AND. lwork(ij)) THEN
    10671053        asij(ij) = max(1.0E-21, asij(ij))
    1068         asij(ij) = 1.0/asij(ij)
     1054        asij(ij) = 1.0 / asij(ij)
    10691055        bsum(ij, i) = 0.0
    10701056      END IF
    … …  
    10721058    DO j = minorig, nl + 1
    10731059      DO ij = 1, ncum
    1074         IF ((i>=icb(ij)+1) .AND. (i<=inb(ij)) .AND. (j>=icb( &
    1075             ij)) .AND. (j<=inb(ij)) .AND. lwork(ij)) THEN
    1076           ment(ij, i, j) = ment(ij, i, j)*asij(ij)
     1060        IF ((i>=icb(ij) + 1) .AND. (i<=inb(ij)) .AND. (j>=icb(&
     1061                ij)) .AND. (j<=inb(ij)) .AND. lwork(ij)) THEN
     1062          ment(ij, i, j) = ment(ij, i, j) * asij(ij)
    10771063          bsum(ij, i) = bsum(ij, i) + ment(ij, i, j)
    10781064        END IF
    … …  
    10801066    END DO
    10811067    DO ij = 1, ncum
    1082       IF ((i>=icb(ij)+1) .AND. (i<=inb(ij)) .AND. (bsum(ij, &
    1083           i)<1.0E-18) .AND. lwork(ij)) THEN
     1068      IF ((i>=icb(ij) + 1) .AND. (i<=inb(ij)) .AND. (bsum(ij, &
     1069              i)<1.0E-18) .AND. lwork(ij)) THEN
    10841070        nent(ij, i) = 0
    10851071        ment(ij, i, i) = m(ij, i)
    1086         qent(ij, i, i) = q(ij, nk(ij)) - ep(ij, i)*clw(ij, i)
     1072        qent(ij, i, i) = q(ij, nk(ij)) - ep(ij, i) * clw(ij, i)
    10871073        uent(ij, i, i) = u(ij, nk(ij))
    10881074        vent(ij, i, i) = v(ij, nk(ij))
    … …  
    10911077      END IF
    10921078    END DO
    1093 789 END DO
    1094 
     1079  789 END DO
    10951080
    10961081END SUBROUTINE cv_mixing
    10971082
    10981083SUBROUTINE cv_unsat(nloc, ncum, nd, inb, t, q, qs, gz, u, v, p, ph, h, lv, &
    1099     ep, sigp, clw, m, ment, elij, iflag, mp, qp, up, vp, wt, water, evap)
     1084        ep, sigp, clw, m, ment, elij, iflag, mp, qp, up, vp, wt, water, evap)
    11001085  IMPLICIT NONE
    1101 
    11021086
    11031087  include "cvthermo.h"
    … …  
    11091093  REAL t(nloc, nd), q(nloc, nd), qs(nloc, nd)
    11101094  REAL gz(nloc, nd), u(nloc, nd), v(nloc, nd)
    1111   REAL p(nloc, nd), ph(nloc, nd+1), h(nloc, nd)
     1095  REAL p(nloc, nd), ph(nloc, nd + 1), h(nloc, nd)
    11121096  REAL lv(nloc, nd), ep(nloc, nd), sigp(nloc, nd), clw(nloc, nd)
    11131097  REAL m(nloc, nd), ment(nloc, nd, nd), elij(nloc, nd, nd)
    … …  
    11491133  DO k = 2, nl + 1
    11501134    DO i = 1, ncum
    1151       qp(i, k) = q(i, k-1)
    1152       up(i, k) = u(i, k-1)
    1153       vp(i, k) = v(i, k-1)
     1135      qp(i, k) = q(i, k - 1)
     1136      up(i, k) = u(i, k - 1)
     1137      vp(i, k) = v(i, k - 1)
    11541138    END DO
    11551139  END DO
    … …  
    11631147  ! ***                and condensed water flux                    ***
    11641148
    1165 
    11661149  DO i = 1, ncum
    11671150    jtt(i) = 2
    1168     IF (ep(i,inb(i))<=0.0001) iflag(i) = 2
     1151    IF (ep(i, inb(i))<=0.0001) iflag(i) = 2
    11691152    IF (iflag(i)==0) THEN
    11701153      lwork(i) = .TRUE.
    … …  
    11761159  ! ***                    Begin downdraft loop                    ***
    11771160
    1178 
    11791161  CALL zilch(wdtrain, ncum)
    11801162  DO i = nl + 1, 1, -1
    … …  
    11911173    DO ij = 1, ncum
    11921174      IF ((i<=inb(ij)) .AND. (lwork(ij))) THEN
    1193         wdtrain(ij) = g*ep(ij, i)*m(ij, i)*clw(ij, i)
     1175        wdtrain(ij) = g * ep(ij, i) * m(ij, i) * clw(ij, i)
    11941176      END IF
    11951177    END DO
    … …  
    11991181        DO ij = 1, ncum
    12001182          IF ((i<=inb(ij)) .AND. (lwork(ij))) THEN
    1201             awat = elij(ij, j, i) - (1.-ep(ij,i))*clw(ij, i)
     1183            awat = elij(ij, j, i) - (1. - ep(ij, i)) * clw(ij, i)
    12021184            awat = max(0.0, awat)
    1203             wdtrain(ij) = wdtrain(ij) + g*awat*ment(ij, j, i)
     1185            wdtrain(ij) = wdtrain(ij) + g * awat * ment(ij, j, i)
    12041186          END IF
    12051187        END DO
    … …  
    12221204        ! rain   ***
    12231205
    1224         IF (t(ij,i)>273.0) THEN
     1206        IF (t(ij, i)>273.0) THEN
    12251207          coeff = coeffr
    12261208          wt(ij, i) = omtrain
    12271209        END IF
    1228         qsm = 0.5*(q(ij,i)+qp(ij,i+1))
    1229         afac = coeff*ph(ij, i)*(qs(ij,i)-qsm)/(1.0E4+2.0E3*ph(ij,i)*qs(ij,i))
     1210        qsm = 0.5 * (q(ij, i) + qp(ij, i + 1))
     1211        afac = coeff * ph(ij, i) * (qs(ij, i) - qsm) / (1.0E4 + 2.0E3 * ph(ij, i) * qs(ij, i))
    12301212        afac = max(afac, 0.0)
    12311213        sigt = sigp(ij, i)
    12321214        sigt = max(0.0, sigt)
    12331215        sigt = min(1.0, sigt)
    1234         b6 = 100.*(ph(ij,i)-ph(ij,i+1))*sigt*afac/wt(ij, i)
    1235         c6 = (water(ij,i+1)*wt(ij,i+1)+wdtrain(ij)/sigd)/wt(ij, i)
    1236         revap = 0.5*(-b6+sqrt(b6*b6+4.*c6))
    1237         evap(ij, i) = sigt*afac*revap
    1238         water(ij, i) = revap*revap
     1216        b6 = 100. * (ph(ij, i) - ph(ij, i + 1)) * sigt * afac / wt(ij, i)
     1217        c6 = (water(ij, i + 1) * wt(ij, i + 1) + wdtrain(ij) / sigd) / wt(ij, i)
     1218        revap = 0.5 * (-b6 + sqrt(b6 * b6 + 4. * c6))
     1219        evap(ij, i) = sigt * afac * revap
     1220        water(ij, i) = revap * revap
    12391221
    12401222        ! ***  Calculate precipitating downdraft mass flux under     ***
    … …  
    12421224
    12431225        IF (i>1) THEN
    1244           dhdp = (h(ij,i)-h(ij,i-1))/(p(ij,i-1)-p(ij,i))
     1226          dhdp = (h(ij, i) - h(ij, i - 1)) / (p(ij, i - 1) - p(ij, i))
    12451227          dhdp = max(dhdp, 10.0)
    1246           mp(ij, i) = 100.*ginv*lv(ij, i)*sigd*evap(ij, i)/dhdp
    1247           mp(ij, i) = max(mp(ij,i), 0.0)
     1228          mp(ij, i) = 100. * ginv * lv(ij, i) * sigd * evap(ij, i) / dhdp
     1229          mp(ij, i) = max(mp(ij, i), 0.0)
    12481230
    12491231          ! ***   Add small amount of inertia to downdraft              ***
    12501232
    1251           fac = 20.0/(ph(ij,i-1)-ph(ij,i))
    1252           mp(ij, i) = (fac*mp(ij,i+1)+mp(ij,i))/(1.+fac)
     1233          fac = 20.0 / (ph(ij, i - 1) - ph(ij, i))
     1234          mp(ij, i) = (fac * mp(ij, i + 1) + mp(ij, i)) / (1. + fac)
    12531235
    12541236          ! ***      Force mp to decrease linearly to zero
    … …  
    12571239          ! ***
    12581240
    1259           IF (p(ij,i)>(0.949*p(ij,1))) THEN
     1241          IF (p(ij, i)>(0.949 * p(ij, 1))) THEN
    12601242            jtt(ij) = max(jtt(ij), i)
    1261             mp(ij, i) = mp(ij, jtt(ij))*(p(ij,1)-p(ij,i))/ &
    1262               (p(ij,1)-p(ij,jtt(ij)))
     1243            mp(ij, i) = mp(ij, jtt(ij)) * (p(ij, 1) - p(ij, i)) / &
     1244                    (p(ij, 1) - p(ij, jtt(ij)))
    12631245          END IF
    12641246        END IF
    … …  
    12701252            qstm = qs(ij, 1)
    12711253          ELSE
    1272             qstm = qs(ij, i-1)
     1254            qstm = qs(ij, i - 1)
    12731255          END IF
    1274           IF (mp(ij,i)>mp(ij,i+1)) THEN
    1275             rat = mp(ij, i+1)/mp(ij, i)
    1276             qp(ij, i) = qp(ij, i+1)*rat + q(ij, i)*(1.0-rat) + &
    1277               100.*ginv*sigd*(ph(ij,i)-ph(ij,i+1))*(evap(ij,i)/mp(ij,i))
    1278             up(ij, i) = up(ij, i+1)*rat + u(ij, i)*(1.-rat)
    1279             vp(ij, i) = vp(ij, i+1)*rat + v(ij, i)*(1.-rat)
     1256          IF (mp(ij, i)>mp(ij, i + 1)) THEN
     1257            rat = mp(ij, i + 1) / mp(ij, i)
     1258            qp(ij, i) = qp(ij, i + 1) * rat + q(ij, i) * (1.0 - rat) + &
     1259                    100. * ginv * sigd * (ph(ij, i) - ph(ij, i + 1)) * (evap(ij, i) / mp(ij, i))
     1260            up(ij, i) = up(ij, i + 1) * rat + u(ij, i) * (1. - rat)
     1261            vp(ij, i) = vp(ij, i + 1) * rat + v(ij, i) * (1. - rat)
    12801262          ELSE
    1281             IF (mp(ij,i+1)>0.0) THEN
    1282               qp(ij, i) = (gz(ij,i+1)-gz(ij,i)+qp(ij,i+1)*(lv(ij,i+1)+t(ij, &
    1283                 i+1)*(cl-cpd))+cpd*(t(ij,i+1)-t(ij, &
    1284                 i)))/(lv(ij,i)+t(ij,i)*(cl-cpd))
    1285               up(ij, i) = up(ij, i+1)
    1286               vp(ij, i) = vp(ij, i+1)
     1263            IF (mp(ij, i + 1)>0.0) THEN
     1264              qp(ij, i) = (gz(ij, i + 1) - gz(ij, i) + qp(ij, i + 1) * (lv(ij, i + 1) + t(ij, &
     1265                      i + 1) * (cl - cpd)) + cpd * (t(ij, i + 1) - t(ij, &
     1266                      i))) / (lv(ij, i) + t(ij, i) * (cl - cpd))
     1267              up(ij, i) = up(ij, i + 1)
     1268              vp(ij, i) = vp(ij, i + 1)
    12871269            END IF
    12881270          END IF
    1289           qp(ij, i) = min(qp(ij,i), qstm)
    1290           qp(ij, i) = max(qp(ij,i), 0.0)
     1271          qp(ij, i) = min(qp(ij, i), qstm)
     1272          qp(ij, i) = max(qp(ij, i), 0.0)
    12911273        END IF
    12921274      END IF
    12931275    END DO
    1294 899 END DO
    1295 
     1276  899 END DO
    12961277
    12971278END SUBROUTINE cv_unsat
    12981279
    12991280SUBROUTINE cv_yield(nloc, ncum, nd, nk, icb, inb, delt, t, q, u, v, gz, p, &
    1300     ph, h, hp, lv, cpn, ep, clw, frac, m, mp, qp, up, vp, wt, water, evap, &
    1301     ment, qent, uent, vent, nent, elij, tv, tvp, iflag, wd, qprime, tprime, &
    1302     precip, cbmf, ft, fq, fu, fv, ma, qcondc)
     1281        ph, h, hp, lv, cpn, ep, clw, frac, m, mp, qp, up, vp, wt, water, evap, &
     1282        ment, qent, uent, vent, nent, elij, tv, tvp, iflag, wd, qprime, tprime, &
     1283        precip, cbmf, ft, fq, fu, fv, ma, qcondc)
    13031284  IMPLICIT NONE
    13041285
    … …  
    13131294  REAL t(nloc, nd), q(nloc, nd), u(nloc, nd), v(nloc, nd)
    13141295  REAL gz(nloc, nd)
    1315   REAL p(nloc, nd), ph(nloc, nd+1), h(nloc, nd)
     1296  REAL p(nloc, nd), ph(nloc, nd + 1), h(nloc, nd)
    13161297  REAL hp(nloc, nd), lv(nloc, nd)
    13171298  REAL cpn(nloc, nd), ep(nloc, nd), clw(nloc, nd), frac(nloc)
    … …  
    13431324  ! -- initializations:
    13441325
    1345   delti = 1.0/delt
     1326  delti = 1.0 / delt
    13461327
    13471328  DO i = 1, ncum
    … …  
    13551336      fv(i, k) = 0.0
    13561337      fq(i, k) = 0.0
    1357       lvcp(i, k) = lv(i, k)/cpn(i, k)
     1338      lvcp(i, k) = lv(i, k) / cpn(i, k)
    13581339      qcondc(i, k) = 0.0 ! cld
    13591340      qcond(i, k) = 0.0 ! cld
    … …  
    13701351      ! c     &                /(rowl*g)
    13711352      ! c            precip(i)=precip(i)*delt/86400.
    1372       precip(i) = wt(i, 1)*sigd*water(i, 1)*86400/g
     1353      precip(i) = wt(i, 1) * sigd * water(i, 1) * 86400 / g
    13731354    END IF
    13741355  END DO
    … …  
    13791360
    13801361  DO i = 1, ncum
    1381     wd(i) = betad*abs(mp(i,icb(i)))*0.01*rrd*t(i, icb(i))/(sigd*p(i,icb(i)))
    1382     qprime(i) = 0.5*(qp(i,1)-q(i,1))
    1383     tprime(i) = lv0*qprime(i)/cpd
     1362    wd(i) = betad * abs(mp(i, icb(i))) * 0.01 * rrd * t(i, icb(i)) / (sigd * p(i, icb(i)))
     1363    qprime(i) = 0.5 * (qp(i, 1) - q(i, 1))
     1364    tprime(i) = lv0 * qprime(i) / cpd
    13841365  END DO
    13851366
    … …  
    13881369
    13891370  DO i = 1, ncum
    1390     work(i) = 0.01/(ph(i,1)-ph(i,2))
     1371    work(i) = 0.01 / (ph(i, 1) - ph(i, 2))
    13911372    am(i) = 0.0
    13921373  END DO
    … …  
    13991380  END DO
    14001381  DO i = 1, ncum
    1401     IF ((g*work(i)*am(i))>=delti) iflag(i) = 1
    1402     ft(i, 1) = ft(i, 1) + g*work(i)*am(i)*(t(i,2)-t(i,1)+(gz(i,2)-gz(i, &
    1403       1))/cpn(i,1))
    1404     ft(i, 1) = ft(i, 1) - lvcp(i, 1)*sigd*evap(i, 1)
    1405     ft(i, 1) = ft(i, 1) + sigd*wt(i, 2)*(cl-cpd)*water(i, 2)*(t(i,2)-t(i,1))* &
    1406       work(i)/cpn(i, 1)
    1407     fq(i, 1) = fq(i, 1) + g*mp(i, 2)*(qp(i,2)-q(i,1))*work(i) + &
    1408       sigd*evap(i, 1)
    1409     fq(i, 1) = fq(i, 1) + g*am(i)*(q(i,2)-q(i,1))*work(i)
    1410     fu(i, 1) = fu(i, 1) + g*work(i)*(mp(i,2)*(up(i,2)-u(i,1))+am(i)*(u(i, &
    1411       2)-u(i,1)))
    1412     fv(i, 1) = fv(i, 1) + g*work(i)*(mp(i,2)*(vp(i,2)-v(i,1))+am(i)*(v(i, &
    1413       2)-v(i,1)))
     1382    IF ((g * work(i) * am(i))>=delti) iflag(i) = 1
     1383    ft(i, 1) = ft(i, 1) + g * work(i) * am(i) * (t(i, 2) - t(i, 1) + (gz(i, 2) - gz(i, &
     1384            1)) / cpn(i, 1))
     1385    ft(i, 1) = ft(i, 1) - lvcp(i, 1) * sigd * evap(i, 1)
     1386    ft(i, 1) = ft(i, 1) + sigd * wt(i, 2) * (cl - cpd) * water(i, 2) * (t(i, 2) - t(i, 1)) * &
     1387            work(i) / cpn(i, 1)
     1388    fq(i, 1) = fq(i, 1) + g * mp(i, 2) * (qp(i, 2) - q(i, 1)) * work(i) + &
     1389            sigd * evap(i, 1)
     1390    fq(i, 1) = fq(i, 1) + g * am(i) * (q(i, 2) - q(i, 1)) * work(i)
     1391    fu(i, 1) = fu(i, 1) + g * work(i) * (mp(i, 2) * (up(i, 2) - u(i, 1)) + am(i) * (u(i, &
     1392            2) - u(i, 1)))
     1393    fv(i, 1) = fv(i, 1) + g * work(i) * (mp(i, 2) * (vp(i, 2) - v(i, 1)) + am(i) * (v(i, &
     1394            2) - v(i, 1)))
    14141395  END DO
    14151396  DO j = 2, nl
    14161397    DO i = 1, ncum
    14171398      IF (j<=inb(i)) THEN
    1418         fq(i, 1) = fq(i, 1) + g*work(i)*ment(i, j, 1)*(qent(i,j,1)-q(i,1))
    1419         fu(i, 1) = fu(i, 1) + g*work(i)*ment(i, j, 1)*(uent(i,j,1)-u(i,1))
    1420         fv(i, 1) = fv(i, 1) + g*work(i)*ment(i, j, 1)*(vent(i,j,1)-v(i,1))
     1399        fq(i, 1) = fq(i, 1) + g * work(i) * ment(i, j, 1) * (qent(i, j, 1) - q(i, 1))
     1400        fu(i, 1) = fu(i, 1) + g * work(i) * ment(i, j, 1) * (uent(i, j, 1) - u(i, 1))
     1401        fv(i, 1) = fv(i, 1) + g * work(i) * ment(i, j, 1) * (vent(i, j, 1) - v(i, 1))
    14211402      END IF
    14221403    END DO
    … …  
    14421423    DO k = i + 1, nl + 1
    14431424      DO ij = 1, ncum
    1444         IF ((i>=nk(ij)) .AND. (i<=inb(ij)) .AND. (k<=(inb(ij)+1))) THEN
     1425        IF ((i>=nk(ij)) .AND. (i<=inb(ij)) .AND. (k<=(inb(ij) + 1))) THEN
    14451426          amp1(ij) = amp1(ij) + m(ij, k)
    14461427        END IF
    … …  
    14511432      DO j = i + 1, nl + 1
    14521433        DO ij = 1, ncum
    1453           IF ((j<=(inb(ij)+1)) .AND. (i<=inb(ij))) THEN
     1434          IF ((j<=(inb(ij) + 1)) .AND. (i<=inb(ij))) THEN
    14541435            amp1(ij) = amp1(ij) + ment(ij, k, j)
    14551436          END IF
    … …  
    14691450    DO ij = 1, ncum
    14701451      IF (i<=inb(ij)) THEN
    1471         dpinv = 0.01/(ph(ij,i)-ph(ij,i+1))
    1472         cpinv = 1.0/cpn(ij, i)
    1473 
    1474         ft(ij, i) = ft(ij, i) + g*dpinv*(amp1(ij)*(t(ij,i+1)-t(ij, &
    1475           i)+(gz(ij,i+1)-gz(ij,i))*cpinv)-ad(ij)*(t(ij,i)-t(ij, &
    1476           i-1)+(gz(ij,i)-gz(ij,i-1))*cpinv)) - sigd*lvcp(ij, i)*evap(ij, i)
    1477         ft(ij, i) = ft(ij, i) + g*dpinv*ment(ij, i, i)*(hp(ij,i)-h(ij,i)+t(ij &
    1478           ,i)*(cpv-cpd)*(q(ij,i)-qent(ij,i,i)))*cpinv
    1479         ft(ij, i) = ft(ij, i) + sigd*wt(ij, i+1)*(cl-cpd)*water(ij, i+1)*(t( &
    1480           ij,i+1)-t(ij,i))*dpinv*cpinv
    1481         fq(ij, i) = fq(ij, i) + g*dpinv*(amp1(ij)*(q(ij,i+1)-q(ij, &
    1482           i))-ad(ij)*(q(ij,i)-q(ij,i-1)))
    1483         fu(ij, i) = fu(ij, i) + g*dpinv*(amp1(ij)*(u(ij,i+1)-u(ij, &
    1484           i))-ad(ij)*(u(ij,i)-u(ij,i-1)))
    1485         fv(ij, i) = fv(ij, i) + g*dpinv*(amp1(ij)*(v(ij,i+1)-v(ij, &
    1486           i))-ad(ij)*(v(ij,i)-v(ij,i-1)))
     1452        dpinv = 0.01 / (ph(ij, i) - ph(ij, i + 1))
     1453        cpinv = 1.0 / cpn(ij, i)
     1454
     1455        ft(ij, i) = ft(ij, i) + g * dpinv * (amp1(ij) * (t(ij, i + 1) - t(ij, &
     1456                i) + (gz(ij, i + 1) - gz(ij, i)) * cpinv) - ad(ij) * (t(ij, i) - t(ij, &
     1457                i - 1) + (gz(ij, i) - gz(ij, i - 1)) * cpinv)) - sigd * lvcp(ij, i) * evap(ij, i)
     1458        ft(ij, i) = ft(ij, i) + g * dpinv * ment(ij, i, i) * (hp(ij, i) - h(ij, i) + t(ij &
     1459                , i) * (cpv - cpd) * (q(ij, i) - qent(ij, i, i))) * cpinv
     1460        ft(ij, i) = ft(ij, i) + sigd * wt(ij, i + 1) * (cl - cpd) * water(ij, i + 1) * (t(&
     1461                ij, i + 1) - t(ij, i)) * dpinv * cpinv
     1462        fq(ij, i) = fq(ij, i) + g * dpinv * (amp1(ij) * (q(ij, i + 1) - q(ij, &
     1463                i)) - ad(ij) * (q(ij, i) - q(ij, i - 1)))
     1464        fu(ij, i) = fu(ij, i) + g * dpinv * (amp1(ij) * (u(ij, i + 1) - u(ij, &
     1465                i)) - ad(ij) * (u(ij, i) - u(ij, i - 1)))
     1466        fv(ij, i) = fv(ij, i) + g * dpinv * (amp1(ij) * (v(ij, i + 1) - v(ij, &
     1467                i)) - ad(ij) * (v(ij, i) - v(ij, i - 1)))
    14871468      END IF
    14881469    END DO
    … …  
    14901471      DO ij = 1, ncum
    14911472        IF (i<=inb(ij)) THEN
    1492           awat = elij(ij, k, i) - (1.-ep(ij,i))*clw(ij, i)
     1473          awat = elij(ij, k, i) - (1. - ep(ij, i)) * clw(ij, i)
    14931474          awat = max(awat, 0.0)
    1494           fq(ij, i) = fq(ij, i) + g*dpinv*ment(ij, k, i)*(qent(ij,k,i)-awat-q &
    1495             (ij,i))
    1496           fu(ij, i) = fu(ij, i) + g*dpinv*ment(ij, k, i)*(uent(ij,k,i)-u(ij,i &
    1497             ))
    1498           fv(ij, i) = fv(ij, i) + g*dpinv*ment(ij, k, i)*(vent(ij,k,i)-v(ij,i &
    1499             ))
     1475          fq(ij, i) = fq(ij, i) + g * dpinv * ment(ij, k, i) * (qent(ij, k, i) - awat - q &
     1476                  (ij, i))
     1477          fu(ij, i) = fu(ij, i) + g * dpinv * ment(ij, k, i) * (uent(ij, k, i) - u(ij, i &
     1478                  ))
     1479          fv(ij, i) = fv(ij, i) + g * dpinv * ment(ij, k, i) * (vent(ij, k, i) - v(ij, i &
     1480                  ))
    15001481          ! (saturated updrafts resulting from mixing)               ! cld
    1501           qcond(ij, i) = qcond(ij, i) + (elij(ij,k,i)-awat) ! cld
     1482          qcond(ij, i) = qcond(ij, i) + (elij(ij, k, i) - awat) ! cld
    15021483          nqcond(ij, i) = nqcond(ij, i) + 1. ! cld
    15031484        END IF
    … …  
    15071488      DO ij = 1, ncum
    15081489        IF ((i<=inb(ij)) .AND. (k<=inb(ij))) THEN
    1509           fq(ij, i) = fq(ij, i) + g*dpinv*ment(ij, k, i)*(qent(ij,k,i)-q(ij,i &
    1510             ))
    1511           fu(ij, i) = fu(ij, i) + g*dpinv*ment(ij, k, i)*(uent(ij,k,i)-u(ij,i &
    1512             ))
    1513           fv(ij, i) = fv(ij, i) + g*dpinv*ment(ij, k, i)*(vent(ij,k,i)-v(ij,i &
    1514             ))
     1490          fq(ij, i) = fq(ij, i) + g * dpinv * ment(ij, k, i) * (qent(ij, k, i) - q(ij, i &
     1491                  ))
     1492          fu(ij, i) = fu(ij, i) + g * dpinv * ment(ij, k, i) * (uent(ij, k, i) - u(ij, i &
     1493                  ))
     1494          fv(ij, i) = fv(ij, i) + g * dpinv * ment(ij, k, i) * (vent(ij, k, i) - v(ij, i &
     1495                  ))
    15151496        END IF
    15161497      END DO
    … …  
    15181499    DO ij = 1, ncum
    15191500      IF (i<=inb(ij)) THEN
    1520         fq(ij, i) = fq(ij, i) + sigd*evap(ij, i) + g*(mp(ij,i+1)*(qp(ij, &
    1521           i+1)-q(ij,i))-mp(ij,i)*(qp(ij,i)-q(ij,i-1)))*dpinv
    1522         fu(ij, i) = fu(ij, i) + g*(mp(ij,i+1)*(up(ij,i+1)-u(ij, &
    1523           i))-mp(ij,i)*(up(ij,i)-u(ij,i-1)))*dpinv
    1524         fv(ij, i) = fv(ij, i) + g*(mp(ij,i+1)*(vp(ij,i+1)-v(ij, &
    1525           i))-mp(ij,i)*(vp(ij,i)-v(ij,i-1)))*dpinv
     1501        fq(ij, i) = fq(ij, i) + sigd * evap(ij, i) + g * (mp(ij, i + 1) * (qp(ij, &
     1502                i + 1) - q(ij, i)) - mp(ij, i) * (qp(ij, i) - q(ij, i - 1))) * dpinv
     1503        fu(ij, i) = fu(ij, i) + g * (mp(ij, i + 1) * (up(ij, i + 1) - u(ij, &
     1504                i)) - mp(ij, i) * (up(ij, i) - u(ij, i - 1))) * dpinv
     1505        fv(ij, i) = fv(ij, i) + g * (mp(ij, i + 1) * (vp(ij, i + 1) - v(ij, &
     1506                i)) - mp(ij, i) * (vp(ij, i) - v(ij, i - 1))) * dpinv
    15261507        ! (saturated downdrafts resulting from mixing)               ! cld
    15271508        DO k = i + 1, inb(ij) ! cld
    … …  
    15301511        END DO ! cld
    15311512        ! (particular case: no detraining level is found)            ! cld
    1532         IF (nent(ij,i)==0) THEN ! cld
    1533           qcond(ij, i) = qcond(ij, i) + (1.-ep(ij,i))*clw(ij, i) ! cld
     1513        IF (nent(ij, i)==0) THEN ! cld
     1514          qcond(ij, i) = qcond(ij, i) + (1. - ep(ij, i)) * clw(ij, i) ! cld
    15341515          nqcond(ij, i) = nqcond(ij, i) + 1. ! cld
    15351516        END IF ! cld
    1536         IF (nqcond(ij,i)/=0.) THEN ! cld
    1537           qcond(ij, i) = qcond(ij, i)/nqcond(ij, i) ! cld
     1517        IF (nqcond(ij, i)/=0.) THEN ! cld
     1518          qcond(ij, i) = qcond(ij, i) / nqcond(ij, i) ! cld
    15381519        END IF ! cld
    15391520      END IF
    15401521    END DO
    1541 1500 END DO
     1522  1500 END DO
    15421523
    15431524  ! *** Adjust tendencies at top of convection layer to reflect  ***
    … …  
    15461527  DO ij = 1, ncum
    15471528    fqold = fq(ij, inb(ij))
    1548     fq(ij, inb(ij)) = fq(ij, inb(ij))*(1.-frac(ij))
    1549     fq(ij, inb(ij)-1) = fq(ij, inb(ij)-1) + frac(ij)*fqold*((ph(ij, &
    1550       inb(ij))-ph(ij,inb(ij)+1))/(ph(ij,inb(ij)-1)-ph(ij, &
    1551       inb(ij))))*lv(ij, inb(ij))/lv(ij, inb(ij)-1)
     1529    fq(ij, inb(ij)) = fq(ij, inb(ij)) * (1. - frac(ij))
     1530    fq(ij, inb(ij) - 1) = fq(ij, inb(ij) - 1) + frac(ij) * fqold * ((ph(ij, &
     1531            inb(ij)) - ph(ij, inb(ij) + 1)) / (ph(ij, inb(ij) - 1) - ph(ij, &
     1532            inb(ij)))) * lv(ij, inb(ij)) / lv(ij, inb(ij) - 1)
    15521533    ftold = ft(ij, inb(ij))
    1553     ft(ij, inb(ij)) = ft(ij, inb(ij))*(1.-frac(ij))
    1554     ft(ij, inb(ij)-1) = ft(ij, inb(ij)-1) + frac(ij)*ftold*((ph(ij, &
    1555       inb(ij))-ph(ij,inb(ij)+1))/(ph(ij,inb(ij)-1)-ph(ij, &
    1556       inb(ij))))*cpn(ij, inb(ij))/cpn(ij, inb(ij)-1)
     1534    ft(ij, inb(ij)) = ft(ij, inb(ij)) * (1. - frac(ij))
     1535    ft(ij, inb(ij) - 1) = ft(ij, inb(ij) - 1) + frac(ij) * ftold * ((ph(ij, &
     1536            inb(ij)) - ph(ij, inb(ij) + 1)) / (ph(ij, inb(ij) - 1) - ph(ij, &
     1537            inb(ij)))) * cpn(ij, inb(ij)) / cpn(ij, inb(ij) - 1)
    15571538    fuold = fu(ij, inb(ij))
    1558     fu(ij, inb(ij)) = fu(ij, inb(ij))*(1.-frac(ij))
    1559     fu(ij, inb(ij)-1) = fu(ij, inb(ij)-1) + frac(ij)*fuold*((ph(ij, &
    1560       inb(ij))-ph(ij,inb(ij)+1))/(ph(ij,inb(ij)-1)-ph(ij,inb(ij))))
     1539    fu(ij, inb(ij)) = fu(ij, inb(ij)) * (1. - frac(ij))
     1540    fu(ij, inb(ij) - 1) = fu(ij, inb(ij) - 1) + frac(ij) * fuold * ((ph(ij, &
     1541            inb(ij)) - ph(ij, inb(ij) + 1)) / (ph(ij, inb(ij) - 1) - ph(ij, inb(ij))))
    15611542    fvold = fv(ij, inb(ij))
    1562     fv(ij, inb(ij)) = fv(ij, inb(ij))*(1.-frac(ij))
    1563     fv(ij, inb(ij)-1) = fv(ij, inb(ij)-1) + frac(ij)*fvold*((ph(ij, &
    1564       inb(ij))-ph(ij,inb(ij)+1))/(ph(ij,inb(ij)-1)-ph(ij,inb(ij))))
     1543    fv(ij, inb(ij)) = fv(ij, inb(ij)) * (1. - frac(ij))
     1544    fv(ij, inb(ij) - 1) = fv(ij, inb(ij) - 1) + frac(ij) * fvold * ((ph(ij, &
     1545            inb(ij)) - ph(ij, inb(ij) + 1)) / (ph(ij, inb(ij) - 1) - ph(ij, inb(ij))))
    15651546  END DO
    15661547
    … …  
    15731554    vav(ij) = 0.0
    15741555    DO i = 1, inb(ij)
    1575       ents(ij) = ents(ij) + (cpn(ij,i)*ft(ij,i)+lv(ij,i)*fq(ij,i))*(ph(ij,i)- &
    1576         ph(ij,i+1))
    1577       uav(ij) = uav(ij) + fu(ij, i)*(ph(ij,i)-ph(ij,i+1))
    1578       vav(ij) = vav(ij) + fv(ij, i)*(ph(ij,i)-ph(ij,i+1))
     1556      ents(ij) = ents(ij) + (cpn(ij, i) * ft(ij, i) + lv(ij, i) * fq(ij, i)) * (ph(ij, i) - &
     1557              ph(ij, i + 1))
     1558      uav(ij) = uav(ij) + fu(ij, i) * (ph(ij, i) - ph(ij, i + 1))
     1559      vav(ij) = vav(ij) + fv(ij, i) * (ph(ij, i) - ph(ij, i + 1))
    15791560    END DO
    15801561  END DO
    15811562  DO ij = 1, ncum
    1582     ents(ij) = ents(ij)/(ph(ij,1)-ph(ij,inb(ij)+1))
    1583     uav(ij) = uav(ij)/(ph(ij,1)-ph(ij,inb(ij)+1))
    1584     vav(ij) = vav(ij)/(ph(ij,1)-ph(ij,inb(ij)+1))
     1563    ents(ij) = ents(ij) / (ph(ij, 1) - ph(ij, inb(ij) + 1))
     1564    uav(ij) = uav(ij) / (ph(ij, 1) - ph(ij, inb(ij) + 1))
     1565    vav(ij) = vav(ij) / (ph(ij, 1) - ph(ij, inb(ij) + 1))
    15851566  END DO
    15861567  DO ij = 1, ncum
    15871568    DO i = 1, inb(ij)
    1588       ft(ij, i) = ft(ij, i) - ents(ij)/cpn(ij, i)
    1589       fu(ij, i) = (1.-cu)*(fu(ij,i)-uav(ij))
    1590       fv(ij, i) = (1.-cu)*(fv(ij,i)-vav(ij))
     1569      ft(ij, i) = ft(ij, i) - ents(ij) / cpn(ij, i)
     1570      fu(ij, i) = (1. - cu) * (fu(ij, i) - uav(ij))
     1571      fv(ij, i) = (1. - cu) * (fv(ij, i) - vav(ij))
    15911572    END DO
    15921573  END DO
    … …  
    15941575  DO k = 1, nl + 1
    15951576    DO i = 1, ncum
    1596       IF ((q(i,k)+delt*fq(i,k))<0.0) iflag(i) = 10
    1597     END DO
    1598   END DO
    1599 
     1577      IF ((q(i, k) + delt * fq(i, k))<0.0) iflag(i) = 10
     1578    END DO
     1579  END DO
    16001580
    16011581  DO i = 1, ncum
    … …  
    16241604  DO k = nl, 1, -1
    16251605    DO i = 1, ncum
    1626       ma(i, k) = ma(i, k+1) + m(i, k)
     1606      ma(i, k) = ma(i, k + 1) + m(i, k)
    16271607    END DO
    16281608  END DO
    … …  
    16461626      ax(ij, i) = 0. ! cld
    16471627      DO j = icb(ij), i ! cld
    1648         ax(ij, i) = ax(ij, i) + rrd*(tvp(ij,j)-tv(ij,j)) & ! cld
    1649           *(ph(ij,j)-ph(ij,j+1))/p(ij, j) ! cld
     1628        ax(ij, i) = ax(ij, i) + rrd * (tvp(ij, j) - tv(ij, j)) & ! cld
     1629                * (ph(ij, j) - ph(ij, j + 1)) / p(ij, j) ! cld
    16501630      END DO ! cld
    1651       IF (ax(ij,i)>0.0) THEN ! cld
    1652         wa(ij, i) = sqrt(2.*ax(ij,i)) ! cld
     1631      IF (ax(ij, i)>0.0) THEN ! cld
     1632        wa(ij, i) = sqrt(2. * ax(ij, i)) ! cld
    16531633      END IF ! cld
    16541634    END DO ! cld
    16551635    DO i = 1, nl ! cld
    1656       IF (wa(ij,i)>0.0) &          ! cld
    1657         siga(ij, i) = mac(ij, i)/wa(ij, i) & ! cld
    1658         *rrd*tvp(ij, i)/p(ij, i)/100./delta ! cld
    1659       siga(ij, i) = min(siga(ij,i), 1.0) ! cld
    1660       qcondc(ij, i) = siga(ij, i)*clw(ij, i)*(1.-ep(ij,i)) & ! cld
    1661         +(1.-siga(ij,i))*qcond(ij, i) ! cld
     1636      IF (wa(ij, i)>0.0) &          ! cld
     1637              siga(ij, i) = mac(ij, i) / wa(ij, i) & ! cld
     1638                      * rrd * tvp(ij, i) / p(ij, i) / 100. / delta ! cld
     1639      siga(ij, i) = min(siga(ij, i), 1.0) ! cld
     1640      qcondc(ij, i) = siga(ij, i) * clw(ij, i) * (1. - ep(ij, i)) & ! cld
     1641              + (1. - siga(ij, i)) * qcond(ij, i) ! cld
    16621642    END DO ! cld
    16631643  END DO ! cld
    16641644
    1665 
    16661645END SUBROUTINE cv_yield
    16671646
    16681647SUBROUTINE cv_uncompress(nloc, len, ncum, nd, idcum, iflag, precip, cbmf, ft, &
    1669     fq, fu, fv, ma, qcondc, iflag1, precip1, cbmf1, ft1, fq1, fu1, fv1, ma1, &
    1670     qcondc1)
     1648        fq, fu, fv, ma, qcondc, iflag1, precip1, cbmf1, ft1, fq1, fu1, fv1, ma1, &
     1649        qcondc1)
    16711650  IMPLICIT NONE
    16721651
    … …  
    17091688  END DO
    17101689
    1711 
    17121690END SUBROUTINE cv_uncompress
    17131691

  • LMDZ6/branches/Amaury_dev/libf/phylmdiso/cva_driver.F90

    r5117 r5132  
    5353  USE lmdz_print_control, ONLY: prt_level, lunout
    5454  USE add_phys_tend_mod, ONLY: fl_cor_ebil
     55  USE lmdz_abort_physic, ONLY: abort_physic
    5556#ifdef ISO
    5657  USE infotrac_phy, ONLY: ntraciso=>ntiso,niso,niso,itZonIso,nzone

  • LMDZ6/branches/Amaury_dev/libf/phylmdiso/isotopes_mod.F90

    r5117 r5132  
    432432
    433433SUBROUTINE getinp_s(nam, val, def, lDisp)
    434    USE ioipsl_getincom, ONLY: getin
     434   USE IOIPSL, ONLY: getin
    435435   USE lmdz_phys_mpi_data, ONLY:  is_mpi_root
    436436   USE lmdz_phys_omp_data, ONLY:  is_omp_root
    … …  
    451451
    452452SUBROUTINE getinp_i(nam, val, def, lDisp)
    453    USE ioipsl_getincom, ONLY: getin
     453   USE IOIPSL, ONLY: getin
    454454   USE lmdz_phys_mpi_data, ONLY:  is_mpi_root
    455455   USE lmdz_phys_omp_data, ONLY:  is_omp_root
    … …  
    470470
    471471SUBROUTINE getinp_r(nam, val, def, lDisp)
    472    USE ioipsl_getincom, ONLY: getin
     472   USE IOIPSL, ONLY: getin
    473473   USE lmdz_phys_mpi_data, ONLY:  is_mpi_root
    474474   USE lmdz_phys_omp_data, ONLY:  is_omp_root
    … …  
    489489
    490490SUBROUTINE getinp_l(nam, val, def, lDisp)
    491    USE ioipsl_getincom, ONLY: getin
     491   USE IOIPSL, ONLY: getin
    492492   USE lmdz_phys_mpi_data, ONLY:  is_mpi_root
    493493   USE lmdz_phys_omp_data, ONLY:  is_omp_root

  • LMDZ6/branches/Amaury_dev/libf/phylmdiso/isotopes_routines_mod.F90

    r5117 r5132  
    44MODULE isotopes_routines_mod
    55  USE infotrac_phy, ONLY: niso, ntraciso=>ntiso, index_trac=>itZonIso, ntraceurs_zone=>nzone
     6  USE lmdz_abort_physic, ONLY: abort_physic
    67IMPLICIT NONE
    78

  • LMDZ6/branches/Amaury_dev/libf/phylmdiso/isotopes_verif_mod.F90

    r5117 r5132  
    99!#endif
    1010USE infotrac_phy, ONLY: ntraciso=>ntiso, niso, itZonIso, nzone
     11USE lmdz_abort_physic, ONLY: abort_physic
    1112IMPLICIT NONE
    1213save

  • LMDZ6/branches/Amaury_dev/libf/phylmdiso/isotrac_mod.F90

    r5117 r5132  
    66  USE lmdz_readTracFiles, ONLY: delPhase
    77  USE isotopes_mod,      ONLY: ridicule, get_in
     8  USE lmdz_abort_physic, ONLY: abort_physic
    89
    910  IMPLICIT NONE

  • LMDZ6/branches/Amaury_dev/libf/phylmdiso/limit_read_mod.F90

    r5117 r5132  
    1111!  limit_read_sst     : return sea ice temperature   
    1212!  limit_read_tot     : read limit.nc and store the fields in local modules variables
     13
     14  USE lmdz_abort_physic, ONLY: abort_physic
    1315
    1416  IMPLICIT NONE

  • LMDZ6/branches/Amaury_dev/libf/phylmdiso/lmdz_lscp_old.F90

    r5117 r5132  
    44
    55MODULE lmdz_lscp_old
     6  USE lmdz_abort_physic, ONLY: abort_physic
    67CONTAINS
    78SUBROUTINE fisrtilp(dtime,paprs,pplay,t,q,ptconv,ratqs, &

  • LMDZ6/branches/Amaury_dev/libf/phylmdiso/phyaqua_mod.F90

    r5117 r5132  
    44MODULE phyaqua_mod
    55  ! Routines complementaires pour la physique planetaire.
     6  USE lmdz_abort_physic, ONLY: abort_physic
    67  IMPLICIT NONE
    78

  • LMDZ6/branches/Amaury_dev/libf/phylmdiso/phyetat0_mod.F90

    r5117 r5132  
    22
    33MODULE phyetat0_mod
     4  USE lmdz_abort_physic, ONLY: abort_physic
    45
    56  PRIVATE

  • LMDZ6/branches/Amaury_dev/libf/phylmdiso/phyredem.F90

    r5117 r5132  
    5151  USE time_phylmdz_mod, ONLY: annee_ref, day_end, itau_phy, pdtphys
    5252  USE config_ocean_skin_m, ONLY: activate_ocean_skin
     53  USE lmdz_abort_physic, ONLY: abort_physic
    5354
    5455  IMPLICIT none

  • LMDZ6/branches/Amaury_dev/libf/phylmdiso/physiq_mod.F90

    r5128 r5132  
    111111    USE phys_output_write_mod
    112112
    113     USE lmdz_cppkeys_wrapper, ONLY: CPPKEY_DUST
     113    USE lmdz_cppkeys_wrapper, ONLY: CPPKEY_DUST, CPPKEY_STRATAER
    114114
    115115!!!!!!!!!!!!!!!!!! "USE" section for CPP keys !!!!!!!!!!!!!!!!!!!!!!!!
    … …  
    14941494
    14951495IF (CPPKEY_STRATAER) THEN
     1496#ifdef ISO
     1497CALL abort_gcm("physiq_mod", "StratAer isn't ISO-compatible for now, 07/24",1)
     1498#else
    14961499       CALL strataer_init
    14971500       CALL strataer_nuc_init
    14981501       CALL strataer_emiss_init
     1502#endif
    14991503END IF
    15001504
    … …  
    16201624       CALL getin_p('iflag_phytrac',iflag_phytrac)
    16211625IF (CPPKEY_DUST) THEN
    1622        IF (iflag_phytrac.EQ.0) THEN
     1626       IF (iflag_phytrac==0) THEN
    16231627         WRITE(lunout,*) 'In order to run with SPLA, iflag_phytrac will be forced to 1'
    16241628         iflag_phytrac = 1
    … …  
    57765780
    57775781#ifdef CPP_RRTM
    5778 IF (CPPKEY_STRATER) THEN
     5782IF (CPPKEY_STRATAER) THEN
     5783#ifdef ISO
     5784CALL abort_gcm("physiq_mod", "StratAer isn't ISO-compatible for now, 07/24",1)
     5785#else
    57795786       !--compute stratospheric mask
    57805787       CALL stratosphere_mask(missing_val, pphis, t_seri, pplay, latitude_deg)
    57815788       !--interactive strat aerosols
    57825789       CALL calcaerosolstrato_rrtm(pplay,t_seri,paprs,debut)
     5790#endif
    57835791END IF
    57845792#endif
    … …  
    64116419    IF (ok_qch4) THEN
    64126420!      d_q_ch4: H2O source from CH4 in MMR/s (mass mixing ratio/s or kg H2O/kg air/s)
    6413 IF (CPPKEY_STRATER) THEN
     6421IF (CPPKEY_STRATAER) THEN
     6422#ifdef ISO
     6423CALL abort_gcm("physiq_mod", "StratAer isn't ISO-compatible for now, 07/24",1)
     6424#else
    64146425       CALL stratH2O_methox(debut,paprs,d_q_ch4)
     6426#endif
    64156427ELSE
    64166428!      ECMWF routine METHOX
    … …  
    64356447
    64366448
    6437 IF (CPPKEY_STRATER) THEN
     6449IF (CPPKEY_STRATAER) THEN
     6450#ifdef ISO
     6451CALL abort_gcm("physiq_mod", "StratAer isn't ISO-compatible for now, 07/24",1)
     6452#else
    64386453    IF (ok_qemiss) THEN
    64396454       flh2o=1
    … …  
    64786493       END SELECT ! emission scenario (flag_emit)
    64796494    ENDIF
     6495#endif
    64806496END IF
    64816497
    … …  
    69326948    ENDDO
    69336949
    6934 IF (CPPKEY_STRATER) THEN
     6950IF (CPPKEY_STRATAER) THEN
     6951#ifdef ISO
     6952CALL abort_gcm("physiq_mod", "StratAer isn't ISO-compatible for now, 07/24",1)
     6953#else
    69356954    IF (ok_qemiss) THEN
    69366955       DO k = 1, klev
    … …  
    69386957       ENDDO
    69396958    ENDIF
     6959#endif
    69406960END IF
    69416961    IF (ok_qch4) THEN
Note: See TracChangeset for help on using the changeset viewer.