da_cumulus.inc

1 subroutine da_cumulus (ZCB, TCB, QCB, PCB, PK, TE,                &
2                      Z, T, Q, LCB, LCT, PCT, ZCT, kts, kte)
3 
4    !-----------------------------------------------------------------------
5    ! Purpose: TBD
6    !-----------------------------------------------------------------------
7 
8    implicit none
9 
10    integer, intent(in)                      :: kts, kte
11    real, intent(inout)                      :: ZCB, TCB, QCB, PCB
12    real, intent(in),   dimension(kts:kte)   :: PK, TE
13    real, intent(out),  dimension(kts:kte)   :: Z, T, Q
14    integer, intent(out)                     :: LCB, LCT
15    real, intent(out)                        :: PCT, ZCT
16 
17    integer   :: k, IA, L, NCB
18    real      :: CP, R, HL, EM, ET, P
19    real      :: TLL, QLL, PLL, ZLL, TBAR, PBAR, QBAR
20    real      :: DP, DZ, DDT, DT
21 
22    CP=1004.
23    R=2000./7.
24    HL=2.49E06
25    DT=0.1
26    IA=1000
27 
28    do k = kts, kte
29       Z(k) = 0.0
30       T(k) = 0.0
31       Q(k) = 0.0
32    end do
33 
34    EM=gravity*ZCB+CP*TCB+HL*QCB
35 
36    NCB=kts
37 
38    if (PK(kte) > PCB) then
39       NCB=kte
40    end if
41 
42    do L=kte-1,kts,-1
43       if (PK(L) > PCB) then
44          NCB=L+1
45          exit
46       end if
47    end do
48 
49    do L=NCB,kte
50       P=PK(L)
51       do K=1,IA
52          if (L == NCB) then
53             TLL=TCB
54             QLL=QCB
55             PLL=PCB
56             ZLL=ZCB
57          else
58             TLL=T(L-1)
59             QLL=Q(L-1)
60             PLL=PK(L-1)
61             ZLL=Z(L-1)
62          end if
63 
64          T(L)=TLL-(K*DT)
65 
66          call da_qfrmrh(P, T(L), 100., Q(L))
67 
68          TBAR=0.5*(T(L)+TLL)
69          QBAR=0.5*(Q(L)+QLL)
70          PBAR=0.5*(P+PLL)
71          DP=PLL-P
72          DZ=(R*TBAR*(1.+0.61*QBAR)*DP)/(gravity*PBAR)
73          Z(L)=ZLL+DZ
74          ET=gravity*Z(L)+CP*T(L)+HL*Q(L)
75          if ((ET-EM) <= 0.) exit
76       end do
77    end do
78 
79    LCT=NCB
80 
81    do K=kte,NCB+1,-1
82       DDT=T(K)-TE(K)
83 
84       if (DDT >= 0.) then
85          LCT=K
86          exit
87       end if
88    end do
89 
90    LCB=LCT
91 
92    do K=NCB,kte
93       DDT=T(K)-TE(K)
94       if (DDT >= 0.) then
95          LCB=K
96          exit
97       end if
98    end do
99 
100    PCT=PK(LCT)
101    ZCT=Z(LCT)
102    PCB=PK(LCB)
103    ZCB=Z(LCB)
104 
105 end subroutine da_cumulus
106 
107