wrftladj/module_cu_du_tl.F

   1 !        Generated by TAPENADE     (INRIA, Tropics team)
   2 !  Tapenade 3.7 (r4786) - 21 Feb 2013 15:53
   3 !
   4 !  Differentiation of ducu in forward (tangent) mode (with options r8):
   5 !   variations   of useful results: raincv rthcuten pratec rqvcuten
   6 !   with respect to varying inputs: th raincv t rthcuten pcps qv
   7 !                z pratec rqvcuten rho dz8w
   8 MODULE g_module_cu_du
   9    USE module_wrf_error
  10   REAL    , PARAMETER :: cincap = -10.
  11   REAL    , PARAMETER :: capemin = 10.
  12   REAL    , PARAMETER :: dpthmin = 1000.
  13   REAL    , PARAMETER :: alpha = 0.00002
  14   REAL    , PARAMETER :: eps = 0.5
  15   REAL    , PARAMETER :: Vfall = 5.
  16 !--------------------------------------------------------------------
  17 CONTAINS
  18
  19 SUBROUTINE DUCU_D(ids, ide, jds, jde, kds, kde, ims, ime, jms, jme, kms&
  20 &  , kme, its, ite, jts, jte, kts, kte, dt, ktau, dx, rho, rhod, raincv, &
  21 &  raincvd, nca, pratec, pratecd, u, v, th, thd, t, td, w, dz8w, dz8wd, z&
  22 &  , zd, pcps, pcpsd, pi, w0avg, cp, rd, rv, g, xlv, ep2, svp1, svp2, &
  23 &  svp3, svpt0, stepcu, cu_act_flag, warm_rain, cutop, cubot, qv, qvd, &
  24 &  rthcuten, rthcutend, rqvcuten, rqvcutend)
  25   IMPLICIT NONE
  26 !
  27 !-------------------------------------------------------------
  28   INTEGER, INTENT(IN) :: ids, ide, jds, jde, kds, kde, ims, ime, jms, &
  29 &  jme, kms, kme, its, ite, jts, jte, kts, kte
  30   INTEGER, INTENT(IN) :: stepcu
  31   LOGICAL, INTENT(IN) :: warm_rain
  32   REAL, INTENT(IN) :: xlv
  33   REAL, INTENT(IN) :: cp, rd, rv, g, ep2
  34   REAL, INTENT(IN) :: svp1, svp2, svp3, svpt0
  35   INTEGER, INTENT(IN) :: ktau
  36   REAL, DIMENSION(ims:ime, kms:kme, jms:jme), INTENT(IN) :: u, v, w, th&
  37 &  , t, qv, dz8w, z, pcps, rho, pi
  38   REAL, DIMENSION(ims:ime, kms:kme, jms:jme), INTENT(IN) :: thd, td, qvd&
  39 &  , dz8wd, zd, pcpsd, rhod
  40 !
  41   REAL, DIMENSION(ims:ime, kms:kme, jms:jme), INTENT(INOUT) :: w0avg
  42   REAL, INTENT(IN) :: dt, dx
  43 !
  44   REAL, DIMENSION(ims:ime, jms:jme), INTENT(INOUT) :: raincv, pratec
  45   REAL, DIMENSION(ims:ime, jms:jme), INTENT(INOUT) :: raincvd, pratecd
  46   REAL, DIMENSION(ims:ime, jms:jme), INTENT(INOUT) :: nca
  47   REAL, DIMENSION(ims:ime, jms:jme), INTENT(OUT) :: cubot, cutop
  48   LOGICAL, DIMENSION(ims:ime, jms:jme), INTENT(INOUT) :: cu_act_flag
  49 !
  50 ! Optional arguments
  51 !
  52   REAL, DIMENSION(ims:ime, kms:kme, jms:jme), OPTIONAL, INTENT(INOUT) ::&
  53 &  rthcuten, rqvcuten
  54   REAL, DIMENSION(ims:ime, kms:kme, jms:jme), OPTIONAL, INTENT(INOUT) ::&
  55 &  rthcutend, rqvcutend
  56 !
  57 ! LOCAL VARS
  58   LOGICAL :: flag_qr, flag_qi, flag_qs
  59   REAL, DIMENSION(kts:kte) :: u1d, v1d, t1d, th1d, dz1d, z1d, qv1d, p1d&
  60 &  , rho1d, w0avg1d
  61   REAL, DIMENSION(kts:kte) :: t1dd, th1dd, dz1dd, z1dd, qv1dd, p1dd, &
  62 &  rho1dd
  63   REAL, DIMENSION(kts:kte) :: dqvdt, dthdt
  64   REAL, DIMENSION(kts:kte) :: dqvdtd, dthdtd
  65   REAL :: pprate, tst, tv, prs, rhoe, w0, scr1, dxsq, rthcumax
  66   REAL :: pprated
  67   INTEGER :: i, j, k, i_start, i_end, j_start, j_end, sz, ntst, icldck
  68 !
  69   dxsq = dx*dx
  70   ntst = stepcu
  71   icldck = MOD(ktau, ntst)
  72   IF (icldck .EQ. 0 .OR. ktau .EQ. 1) THEN
  73 !
  74 !  Keep away from specified and relaxation zone (should be for just specified and nested bc)
  75     sz = 1
  76     IF (ids + sz .LT. its) THEN
  77       i_start = its
  78     ELSE
  79       i_start = ids + sz
  80     END IF
  81     IF (ide - 1 - sz .GT. ite) THEN
  82       i_end = ite
  83     ELSE
  84       i_end = ide - 1 - sz
  85     END IF
  86     IF (jds + sz .LT. jts) THEN
  87       j_start = jts
  88     ELSE
  89       j_start = jds + sz
  90     END IF
  91     IF (jde - 1 - sz .GT. jte) THEN
  92       j_end = jte
  93     ELSE
  94       j_end = jde - 1 - sz
  95     END IF
  96     qv1dd = 0.0_8
  97     dthdtd = 0.0_8
  98     th1dd = 0.0_8
  99     rho1dd = 0.0_8
 100     p1dd = 0.0_8
 101     z1dd = 0.0_8
 102     dz1dd = 0.0_8
 103     t1dd = 0.0_8
 104     dqvdtd = 0.0_8
 105 !
 106     DO j=j_start,j_end
 107       DO i=i_start,i_end
 108         DO k=kts,kte
 109           dqvdtd(k) = 0.0_8
 110           dqvdt(k) = 0.
 111           dthdtd(k) = 0.0_8
 112           dthdt(k) = 0.
 113         END DO
 114         raincvd(i, j) = 0.0_8
 115         raincv(i, j) = 0.
 116         pratecd(i, j) = 0.0_8
 117         pratec(i, j) = 0.
 118         cutop(i, j) = kts
 119         cubot(i, j) = kte + 1
 120 !
 121 ! assign vars from 3D to 1D
 122         DO k=kts,kte
 123           u1d(k) = u(i, k, j)
 124           v1d(k) = v(i, k, j)
 125           t1dd(k) = td(i, k, j)
 126           t1d(k) = t(i, k, j)
 127           th1dd(k) = thd(i, k, j)
 128           th1d(k) = th(i, k, j)
 129           rho1dd(k) = rhod(i, k, j)
 130           rho1d(k) = rho(i, k, j)
 131           qv1dd(k) = qvd(i, k, j)
 132           qv1d(k) = qv(i, k, j)
 133           IF (qv1d(k) .LT. 1.e-08) THEN
 134             qv1dd(k) = 0.0_8
 135             qv1d(k) = 1.e-08
 136           END IF
 137           p1dd(k) = pcpsd(i, k, j)
 138           p1d(k) = pcps(i, k, j)
 139           w0avg1d(k) = w0avg(i, k, j)
 140           dz1dd(k) = dz8wd(i, k, j)
 141           dz1d(k) = dz8w(i, k, j)
 142           z1dd(k) = zd(i, k, j)
 143           z1d(k) = z(i, k, j)
 144         END DO
 145         CALL DUCU1D_D(i, j, u1d, v1d, t1d, t1dd, qv1d, qv1dd, p1d, p1dd&
 146 &                , dz1d, dz1dd, z1d, z1dd, w0avg1d, dt, dx, dxsq, rho1d, &
 147 &                rho1dd, th1d, th1dd, xlv, cp, rd, rv, g, ep2, svp1, svp2&
 148 &                , svp3, svpt0, dqvdt, dqvdtd, dthdt, dthdtd, pprate, &
 149 &                pprated, nca, ntst, cutop, cubot, ids, ide, jds, jde, &
 150 &                kds, kde, ims, ime, jms, jme, kms, kme, its, ite, jts, &
 151 &                jte, kts, kte)
 152         IF (PRESENT(rthcuten) .AND. PRESENT(rqvcuten)) THEN
 153           DO k=kts,kte
 154             rthcutend(i, k, j) = dthdtd(k)
 155             rthcuten(i, k, j) = dthdt(k)
 156             rqvcutend(i, k, j) = dqvdtd(k)
 157             rqvcuten(i, k, j) = dqvdt(k)
 158           END DO
 159           pratecd(i, j) = pprated
 160           pratec(i, j) = pprate
 161           raincvd(i, j) = dt*pprated
 162           raincv(i, j) = pprate*dt
 163         END IF
 164       END DO
 165     END DO
 166   END IF
 167 END SUBROUTINE DUCU_D
 168
 169 !        Generated by TAPENADE     (INRIA, Tropics team)
 170 !  Tapenade 3.7 (r4786) - 21 Feb 2013 15:53
 171 !
 172 !  Differentiation of ducu1d in forward (tangent) mode (with options r8):
 173 !   variations   of useful results: dthdt pprate dqvdt
 174 !   with respect to varying inputs: dthdt p0 z t0 th0 rhoe dzq
 175 !                dqvdt qv0
 176 ! ****************************************************************************
 177 !-----------------------------------------------------------
 178 SUBROUTINE DUCU1D_D(i, j, u0, v0, t0, t0d, qv0, qv0d, p0, p0d, dzq, dzqd&
 179 &  , z, zd, w0avg1d, delt, dx, dxsq, rhoe, rhoed, th0, th0d, xlv, cp, rd&
 180 &  , rv, g, ep2, svp1, svp2, svp3, svpt0, dqvdt, dqvdtd, dthdt, dthdtd, &
 181 &  pprate, pprated, nca, ntst, cutop, cubot, ids, ide, jds, jde, kds, kde&
 182 &  , ims, ime, jms, jme, kms, kme, its, ite, jts, jte, kts, kte)
 183   IMPLICIT NONE
 184 !-----------------------------------------------------------------------
 185 !-----------------------------------------------------------
 186   INTEGER, INTENT(IN) :: ids, ide, jds, jde, kds, kde, ims, ime, jms, &
 187 &  jme, kms, kme, its, ite, jts, jte, kts, kte, i, j, ntst
 188 !
 189   REAL, DIMENSION(kts:kte), INTENT(IN) :: u0, v0, t0, th0, qv0, p0, rhoe&
 190 &  , dzq, z, w0avg1d
 191   REAL, DIMENSION(kts:kte), INTENT(IN) :: t0d, th0d, qv0d, p0d, rhoed, &
 192 &  dzqd, zd
 193 !
 194   REAL, INTENT(IN) :: delt, dx, dxsq
 195 !
 196   REAL, INTENT(IN) :: xlv, cp, rd, rv, g
 197   REAL, INTENT(IN) :: ep2, svp1, svp2, svp3, svpt0
 198 !
 199   REAL, DIMENSION(kts:kte), INTENT(INOUT) :: dqvdt, dthdt
 200   REAL, DIMENSION(kts:kte), INTENT(INOUT) :: dqvdtd, dthdtd
 201   REAL, DIMENSION(ims:ime, jms:jme), INTENT(INOUT) :: nca
 202   REAL, DIMENSION(ims:ime, jms:jme), INTENT(OUT) :: cubot, cutop
 203   REAL, INTENT(OUT) :: pprate
 204   REAL, INTENT(OUT) :: pprated
 205 !
 206 !...DEFINE LOCAL VARIABLES...
 207 !
 208   REAL, DIMENSION(kts:kte) :: cond, h, hs, qs, x
 209   REAL, DIMENSION(kts:kte) :: condd, hd, hsd, qsd, xd
 210   REAL :: buoy, cape, cin, dh, dq, dt, dtm, ep, es, evap, hp, mp, qp, &
 211 &  qsp, rrk, rrkp, tadp, tdp, zb, zg, zi, zt
 212   REAL :: dhd, dqd, dtd, esd, evapd, hpd, mpd, qpd, qspd, rrkd, rrkpd
 213   INTEGER :: ipos, isat, k, kb, ki, kt
 214   REAL :: arg1
 215   REAL :: arg1d
 216   REAL :: arg10
 217   REAL :: arg2
 218   REAL :: abs1
 219   REAL :: abs0
 220   hd = 0.0_8
 221   qsd = 0.0_8
 222   xd = 0.0_8
 223   hsd = 0.0_8
 224 !
 225 !...DEFINE PROFILES
 226   DO k=kts,kte
 227     hd(k) = cp*t0d(k) + g*zd(k) + xlv*qv0d(k)
 228     h(k) = cp*t0(k) + g*z(k) + xlv*qv0(k)
 229     arg1d = (svp2*t0d(k)*(t0(k)-svp3)-svp2*(t0(k)-svpt0)*t0d(k))/(t0(k)-&
 230 &      svp3)**2
 231     arg1 = svp2*(t0(k)-svpt0)/(t0(k)-svp3)
 232     esd = 1000.*svp1*arg1d*EXP(arg1)
 233     es = 1000.*svp1*EXP(arg1)
 234     qsd(k) = (ep2*esd*(p0(k)-es)-ep2*es*(p0d(k)-esd))/(p0(k)-es)**2
 235     qs(k) = ep2*es/(p0(k)-es)
 236     hsd(k) = cp*t0d(k) + g*zd(k) + xlv*qsd(k)
 237     hs(k) = cp*t0(k) + g*z(k) + xlv*qs(k)
 238     xd(k) = (xlv**2*qsd(k)*cp*rv*t0(k)**2-xlv**2*qs(k)*cp*rv*(t0d(k)*t0(&
 239 &      k)+t0(k)*t0d(k)))/(cp*rv*t0(k)*t0(k))**2
 240     x(k) = xlv*xlv*qs(k)/(cp*rv*t0(k)*t0(k))
 241     dthdtd(k) = 0.0_8
 242     dthdt(k) = 0.
 243     dqvdtd(k) = 0.0_8
 244     dqvdt(k) = 0.
 245   END DO
 246   pprate = 0.
 247   zg = z(1) - 0.5*dzq(1)
 248   pprated = 0.0_8
 249 !
 250 !...LOOP OVER PARCELS
 251 loop_origin:DO ki=kts,kte
 252     hpd = hd(ki)
 253     hp = h(ki)
 254     qpd = qv0d(ki)
 255     qp = qv0(ki)
 256     mpd = alpha*(rhoed(ki)*dzq(ki)+rhoe(ki)*dzqd(ki))
 257     mp = alpha*rhoe(ki)*dzq(ki)
 258     zi = z(ki)
 259     buoy = 0.
 260     cape = 0.
 261     cin = 0.
 262     dtm = 0.
 263     isat = 0
 264     ipos = 0
 265     kt = 0
 266     kb = 0
 267     zt = 0.
 268     zb = 0.
 269     cond = 0.
 270     condd = 0.0_8
 271 !
 272 !...LIFT PARCEL
 273 loop_lift:DO k=ki+1,kte
 274       tadp = t0(ki) + g/cp*(z(ki)-z(k))
 275       ep = p0(k)*qv0(ki)/(ep2+qv0(ki))
 276       arg10 = 0.001*ep/svp1
 277       arg2 = 0.001*ep/svp1
 278       tdp = (svpt0-svp3/svp2*ALOG(arg10))/(1.-1./svp2*ALOG(arg2))
 279       IF (tadp .GE. tdp) THEN
 280 !         unsaturated
 281         IF (isat .EQ. 1) PRINT*, i, j, &
 282 &                         'sounding warning: unsat above sat'
 283         dt = tadp - t0(k)
 284         condd(k) = 0.0_8
 285         cond(k) = 0.
 286       ELSE
 287 !         saturated
 288         IF (isat .EQ. 0) THEN
 289           kb = k
 290           zb = z(k) - 0.5*dzq(k)
 291         END IF
 292         isat = 1
 293         dhd = hpd - hsd(k)
 294         dh = hp - hs(k)
 295         dt = dh/cp/(1.+x(k))
 296         qspd = qsd(k) + ((dhd*x(k)/xlv+dh*xd(k)/xlv)*(1.+x(k))-dh*x(k)*&
 297 &          xd(k)/xlv)/(1.+x(k))**2
 298         qsp = qs(k) + dh/xlv*x(k)/(1.+x(k))
 299 !...CONDENSATE PRODUCED
 300         condd(k) = mpd*(qp-qsp) + mp*(qpd-qspd)
 301         cond(k) = mp*(qp-qsp)
 302         qpd = qspd
 303         qp = qsp
 304       END IF
 305       buoy = buoy + g*dt*dzq(k)/t0(k)
 306       IF (cape .LT. buoy) THEN
 307         cape = buoy
 308       ELSE
 309         cape = cape
 310       END IF
 311       IF (buoy .GE. cincap) THEN
 312         IF (cin .GT. buoy) THEN
 313           cin = buoy
 314         ELSE
 315           cin = cin
 316         END IF
 317       END IF
 318       IF (dt .GE. 0.) THEN
 319         kt = k
 320         zt = z(k) + 0.5*dzq(k)
 321       ELSE IF (dt .LT. 0. .AND. dtm .GE. 0.) THEN
 322         IF (dt .GE. 0.) THEN
 323           abs0 = dt
 324         ELSE
 325           abs0 = -dt
 326         END IF
 327         IF (dtm .GE. 0.) THEN
 328           abs1 = dtm
 329         ELSE
 330           abs1 = -dtm
 331         END IF
 332 ! cloud top is level closest to parcel temperature
 333         IF (abs0 .LT. abs1) THEN
 334           kt = k
 335           zt = z(k) + 0.5*dzq(k)
 336         END IF
 337       END IF
 338       dtm = dt
 339 ! continue lifting until buoyancy is gone
 340       IF (buoy .LT. cincap) THEN
 341         GOTO 100
 342       ELSE
 343         IF (buoy .GT. 0.) ipos = 1
 344 !         positive area detected
 345         IF (k .EQ. 1) THEN
 346           kt = k
 347           zt = z(k) + 0.5*dzq(k)
 348           zi = z(ki)
 349           PRINT*, 'sounding warning: cloud top at model top'
 350         END IF
 351       END IF
 352     END DO loop_lift
 353 !
 354 !...CHECK FOR CLOUD
 355  100 IF (isat .NE. 0) THEN
 356 !       no cloud from lifting - no convection
 357       IF (zt - zb .GT. dpthmin) THEN
 358 !       not more than one cloud level - no convection
 359         IF (ipos .NE. 0) THEN
 360 !       no buoyancy in cloud - no convection
 361           IF (cape .GT. capemin) THEN
 362 !       not enough cape
 363 !
 364 !...IF CHECK FOR CLOUD SUCCESSFUL
 365 !
 366 !...DETRAINMENT
 367             dhd = hpd - hsd(kt)
 368             dh = hp - hs(kt)
 369             dtd = (dhd*(1.+x(kt))/cp-dh*xd(kt)/cp)/(1.+x(kt))**2
 370             dt = dh/cp/(1.+x(kt))
 371             dqd = qsd(kt) + ((dhd*x(kt)/xlv+dh*xd(kt)/xlv)*(1.+x(kt))-dh&
 372 &              *x(kt)*xd(kt)/xlv)/(1.+x(kt))**2 - qv0d(kt)
 373             dq = qs(kt) + dh/xlv*x(kt)/(1.+x(kt)) - qv0(kt)
 374             dthdtd(kt) = dthdtd(kt) + (((mpd*dt+mp*dtd)*th0(kt)/t0(kt)+&
 375 &              mp*dt*(th0d(kt)*t0(kt)-th0(kt)*t0d(kt))/t0(kt)**2)*rhoe(kt&
 376 &              )*dzq(kt)-mp*th0(kt)*dt*(rhoed(kt)*dzq(kt)+rhoe(kt)*dzqd(&
 377 &              kt))/t0(kt))/(rhoe(kt)*dzq(kt))**2
 378             dthdt(kt) = dthdt(kt) + mp*(th0(kt)/t0(kt))*dt/(rhoe(kt)*dzq&
 379 &              (kt))
 380             dqvdtd(kt) = dqvdtd(kt) + ((mpd*dq+mp*dqd)*rhoe(kt)*dzq(kt)-&
 381 &              mp*dq*(rhoed(kt)*dzq(kt)+rhoe(kt)*dzqd(kt)))/(rhoe(kt)*dzq&
 382 &              (kt))**2
 383             dqvdt(kt) = dqvdt(kt) + mp*dq/(rhoe(kt)*dzq(kt))
 384 !
 385 !...SUBSIDENCE
 386 loop_subsidence:DO k=kt-1,ki,-1
 387               dthdtd(k) = dthdtd(k) + ((mpd*(th0(k+1)-th0(k))+mp*(th0d(k&
 388 &                +1)-th0d(k)))*rhoe(k)*dzq(k)-mp*(th0(k+1)-th0(k))*(rhoed&
 389 &                (k)*dzq(k)+rhoe(k)*dzqd(k)))/(rhoe(k)*dzq(k))**2
 390               dthdt(k) = dthdt(k) + mp*(th0(k+1)-th0(k))/(rhoe(k)*dzq(k)&
 391 &                )
 392               dqvdtd(k) = dqvdtd(k) + ((mpd*(qv0(k+1)-qv0(k))+mp*(qv0d(k&
 393 &                +1)-qv0d(k)))*rhoe(k)*dzq(k)-mp*(qv0(k+1)-qv0(k))*(rhoed&
 394 &                (k)*dzq(k)+rhoe(k)*dzqd(k)))/(rhoe(k)*dzq(k))**2
 395               dqvdt(k) = dqvdt(k) + mp*(qv0(k+1)-qv0(k))/(rhoe(k)*dzq(k)&
 396 &                )
 397             END DO loop_subsidence
 398 !
 399 !...RAINFALL AND EVAPORATION
 400             rrkp = 0.
 401             rrkpd = 0.0_8
 402 loop_rainfall:DO k=kt,1,-1
 403               rrkd = rrkpd + condd(k)
 404               rrk = rrkp + cond(k)
 405               evapd = eps*((dzqd(k)*rrkp+dzq(k)*rrkpd)*(qs(k)-qv0(k))/&
 406 &                vfall+dzq(k)*rrkp*(qsd(k)-qv0d(k))/vfall)
 407               evap = dzq(k)*rrkp/vfall*eps*(qs(k)-qv0(k))
 408 ! restrict evap to below cloud base
 409               IF (k .GE. kb) THEN
 410                 evap = 0.
 411                 evapd = 0.0_8
 412               END IF
 413               IF (rrk .GT. evap) THEN
 414                 evap = evap
 415               ELSE
 416                 evapd = rrkd
 417                 evap = rrk
 418               END IF
 419               rrkd = rrkd - evapd
 420               rrk = rrk - evap
 421               dqvdtd(k) = dqvdtd(k) + (evapd*rhoe(k)*dzq(k)-evap*(rhoed(&
 422 &                k)*dzq(k)+rhoe(k)*dzqd(k)))/(rhoe(k)*dzq(k))**2
 423               dqvdt(k) = dqvdt(k) + evap/(rhoe(k)*dzq(k))
 424               dthdtd(k) = dthdtd(k) - (xlv*((th0d(kt)*t0(kt)-th0(kt)*t0d&
 425 &                (kt))*evap/t0(kt)**2+th0(kt)*evapd/t0(kt))*rhoe(k)*dzq(k&
 426 &                )/cp-xlv*th0(kt)*evap*(rhoed(k)*dzq(k)+rhoe(k)*dzqd(k))/&
 427 &                (cp*t0(kt)))/(rhoe(k)*dzq(k))**2
 428               dthdt(k) = dthdt(k) - xlv/cp*(th0(kt)/t0(kt))*evap/(rhoe(k&
 429 &                )*dzq(k))
 430               rrkpd = rrkd
 431               rrkp = rrk
 432             END DO loop_rainfall
 433             pprated = pprated + rrkpd
 434             pprate = pprate + rrkp
 435           END IF
 436         END IF
 437       END IF
 438     END IF
 439   END DO loop_origin
 440 END SUBROUTINE DUCU1D_D
 441
 442 !        Generated by TAPENADE     (INRIA, Tropics team)
 443 !  Tapenade 3.6 (r4343) - 10 Feb 2012 10:52
 444 !
 445 !  Differentiation of ducuinit in forward (tangent) mode:
 446 !   variations   of useful results: rqccuten rthcuten w0avg rqvcuten
 447 !   with respect to varying inputs: rqccuten rthcuten w0avg rqvcuten
 448 SUBROUTINE DUCUINIT_D(rthcuten, rthcutend, rqvcuten, rqvcutend, rqccuten&
 449 &  , rqccutend, rqrcuten, rqicuten, rqscuten, nca, w0avg, w0avgd, p_qc, &
 450 &  p_qr, svp1, svp2, svp3, svpt0, p_first_scalar, restart, &
 451 &  allowed_to_read, ids, ide, jds, jde, kds, kde, ims, ime, jms, jme, kms&
 452 &  , kme, its, ite, jts, jte, kts, kte)
 453   IMPLICIT NONE
 454 !--------------------------------------------------------------------
 455   LOGICAL, INTENT(IN) :: restart, allowed_to_read
 456   INTEGER, INTENT(IN) :: ids, ide, jds, jde, kds, kde, ims, ime, jms, &
 457 &  jme, kms, kme, its, ite, jts, jte, kts, kte
 458   INTEGER, INTENT(IN) :: p_qc, p_qr, p_first_scalar
 459   REAL, DIMENSION(ims:ime, kms:kme, jms:jme), INTENT(OUT) :: rthcuten, &
 460 &  rqvcuten, rqccuten, rqrcuten, rqicuten, rqscuten
 461   REAL, DIMENSION(ims:ime, kms:kme, jms:jme), INTENT(OUT) :: rthcutend, &
 462 &  rqvcutend, rqccutend
 463   REAL, DIMENSION(ims:ime, kms:kme, jms:jme), INTENT(OUT) :: w0avg
 464   REAL, DIMENSION(ims:ime, kms:kme, jms:jme), INTENT(OUT) :: w0avgd
 465   REAL, DIMENSION(ims:ime, jms:jme), INTENT(INOUT) :: nca
 466   INTEGER :: i, j, k, itf, jtf, ktf
 467   REAL, INTENT(IN) :: svp1, svp2, svp3, svpt0
 468   IF (jte .GT. jde - 1) THEN
 469     jtf = jde - 1
 470   ELSE
 471     jtf = jte
 472   END IF
 473   IF (kte .GT. kde - 1) THEN
 474     ktf = kde - 1
 475   ELSE
 476     ktf = kte
 477   END IF
 478   IF (ite .GT. ide - 1) THEN
 479     itf = ide - 1
 480   ELSE
 481     itf = ite
 482   END IF
 483   IF (.NOT.restart) THEN
 484     DO j=jts,jtf
 485       DO k=kts,ktf
 486         DO i=its,itf
 487           rthcutend(i, k, j) = 0.0
 488           rthcuten(i, k, j) = 0.
 489           rqvcutend(i, k, j) = 0.0
 490           rqvcuten(i, k, j) = 0.
 491         END DO
 492       END DO
 493     END DO
 494     IF (p_qc .GE. p_first_scalar) THEN
 495       DO j=jts,jtf
 496         DO k=kts,ktf
 497           DO i=its,itf
 498             rqccutend(i, k, j) = 0.0
 499             rqccuten(i, k, j) = 0.
 500           END DO
 501         END DO
 502       END DO
 503     END IF
 504     DO j=jts,jtf
 505       DO i=its,itf
 506         nca(i, j) = -100.
 507       END DO
 508     END DO
 509     DO j=jts,jtf
 510       DO k=kts,ktf
 511         DO i=its,itf
 512           w0avgd(i, k, j) = 0.0
 513           w0avg(i, k, j) = 0.
 514         END DO
 515       END DO
 516     END DO
 517   END IF
 518 END SUBROUTINE DUCUINIT_D
 519
 520 END MODULE g_module_cu_du