/[lmdze]/trunk/phylmd/CV30_routines/cv30_undilute1.f
ViewVC logotype

Contents of /trunk/phylmd/CV30_routines/cv30_undilute1.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 195 - (show annotations)
Wed May 18 17:56:44 2016 UTC (8 years, 1 month ago) by guez
Original Path: trunk/Sources/phylmd/CV30_routines/cv30_undilute1.f
File size: 6883 byte(s)
In cv30_feed, iflag1 is 0 on entry so we can simplify the test for
iflag1 = 7.

In cv30_feed, for the computation of icb, replaced sequential search
(with a useless end of loop on k) by a call to locate.

In CV30 routines, replaced len, nloc, nd, na by klon or
klev. Philosophy: no more generality than actually necessary.

Converted as many variables as possible to named constants in
cv30_param_m and downgraded pbcrit, ptcrit, dtovsh, dpbase, dttrig,
tau, delta to local objects in procedures. spfac, betad and omtrain
are useless and removed.

Instead of filling the array sigp with the constant spfac in
cv30_undilute2, just made sigp a constant in cv30_unsat.

In cv_driver, define as allocatable variables that are only
used on the range (ncum, nl).

1 module cv30_undilute1_m
2
3 implicit none
4
5 contains
6
7 SUBROUTINE cv30_undilute1(t, q, qs, gz, plcl, p, nk, icb, tp, tvp, clw, icbs)
8
9 ! UNDILUTE (ADIABATIC) UPDRAFT / 1st part
10 ! (up through ICB for convect4, up through ICB + 1 for convect3)
11 ! Calculates the lifted parcel virtual temperature at nk, the
12 ! actual temperature, and the adiabatic liquid water content.
13
14 ! Equivalent de TLIFT entre NK et ICB+1 inclus
15
16 ! Differences with convect4:
17 ! - specify plcl in input
18 ! - icbs is the first level above LCL (may differ from icb)
19 ! - in the iterations, used x(icbs) instead x(icb)
20 ! - many minor differences in the iterations
21 ! - tvp is computed in only one time
22 ! - icbs: first level above Plcl (IMIN de TLIFT) in output
23 ! - if icbs=icb, compute also tp(icb+1), tvp(icb+1) & clw(icb+1)
24
25 use cv30_param_m, only: minorig, nl
26 use cv_thermo_m, only: cl, clmcpv, cpd, cpv, eps, lv0, rrv
27 USE dimphy, ONLY: klev, klon
28
29 ! inputs:
30 integer, intent(in):: nk(klon), icb(klon)
31 real, intent(in):: t(klon, klev)
32 real, intent(in):: q(klon, klev), qs(klon, klev), gz(klon, klev)
33 real, intent(in):: p(klon, klev)
34 real, intent(in):: plcl(klon) ! convect3
35
36 ! outputs:
37 real tp(klon, klev), tvp(klon, klev), clw(klon, klev)
38
39 ! local variables:
40 integer i, k
41 integer icb1(klon), icbs(klon), icbsmax2 ! convect3
42 real tg, qg, alv, s, ahg, tc, denom, es
43 real ah0(klon), cpp(klon)
44 real tnk(klon), qnk(klon), gznk(klon), ticb(klon), gzicb(klon)
45 real qsicb(klon) ! convect3
46 real cpinv(klon) ! convect3
47
48 !-------------------------------------------------------------------
49
50 ! Calculates the lifted parcel virtual temperature at nk,
51 ! the actual temperature, and the adiabatic
52 ! liquid water content. The procedure is to solve the equation.
53 ! cp*tp+L*qp+phi=cp*tnk+L*qnk+gznk.
54
55 do i=1, klon
56 tnk(i)=t(i, nk(i))
57 qnk(i)=q(i, nk(i))
58 gznk(i)=gz(i, nk(i))
59 end do
60
61 ! *** Calculate certain parcel quantities, including static energy ***
62
63 do i=1, klon
64 ah0(i)=(cpd*(1.-qnk(i))+cl*qnk(i))*tnk(i) &
65 +qnk(i)*(lv0-clmcpv*(tnk(i)-273.15))+gznk(i)
66 cpp(i)=cpd*(1.-qnk(i))+qnk(i)*cpv
67 cpinv(i)=1./cpp(i)
68 end do
69
70 ! *** Calculate lifted parcel quantities below cloud base ***
71
72 do i=1, klon !convect3
73 icb1(i)=MAX(icb(i), 2) !convect3
74 icb1(i)=MIN(icb(i), nl) !convect3
75 ! if icb is below LCL, start loop at ICB+1:
76 ! (icbs est le premier niveau au-dessus du LCL)
77 icbs(i)=icb1(i) !convect3
78 if (plcl(i) < p(i, icb1(i))) then
79 icbs(i)=MIN(icbs(i)+1, nl) !convect3
80 endif
81 enddo !convect3
82
83 do i=1, klon !convect3
84 ticb(i)=t(i, icbs(i)) !convect3
85 gzicb(i)=gz(i, icbs(i)) !convect3
86 qsicb(i)=qs(i, icbs(i)) !convect3
87 enddo !convect3
88
89 ! Re-compute icbsmax (icbsmax2): !convect3
90 ! !convect3
91 icbsmax2=2 !convect3
92 do i=1, klon !convect3
93 icbsmax2=max(icbsmax2, icbs(i)) !convect3
94 end do
95
96 ! initialization outputs:
97
98 do k=1, icbsmax2 ! convect3
99 do i=1, klon ! convect3
100 tp(i, k) = 0.0 ! convect3
101 tvp(i, k) = 0.0 ! convect3
102 clw(i, k) = 0.0 ! convect3
103 enddo ! convect3
104 enddo ! convect3
105
106 ! tp and tvp below cloud base:
107
108 do k=minorig, icbsmax2-1
109 do i=1, klon
110 tp(i, k)=tnk(i)-(gz(i, k)-gznk(i))*cpinv(i)
111 tvp(i, k)=tp(i, k)*(1.+qnk(i)/eps-qnk(i)) !whole thing (convect3)
112 end do
113 end do
114
115 ! *** Find lifted parcel quantities above cloud base ***
116
117 do i=1, klon
118 tg=ticb(i)
119 qg=qsicb(i) ! convect3
120 !debug alv=lv0-clmcpv*(ticb(i)-t0)
121 alv=lv0-clmcpv*(ticb(i)-273.15)
122
123 ! First iteration.
124
125 s=cpd*(1.-qnk(i))+cl*qnk(i) &
126 +alv*alv*qg/(rrv*ticb(i)*ticb(i)) ! convect3
127 s=1./s
128
129 ahg=cpd*tg+(cl-cpd)*qnk(i)*tg+alv*qg+gzicb(i) ! convect3
130 tg=tg+s*(ah0(i)-ahg)
131
132 !debug tc=tg-t0
133 tc=tg-273.15
134 denom=243.5+tc
135 denom=MAX(denom, 1.0) ! convect3
136
137 es=6.112*exp(17.67*tc/denom)
138 qg=eps*es/(p(i, icbs(i))-es*(1.-eps))
139
140 ! Second iteration.
141
142 ahg=cpd*tg+(cl-cpd)*qnk(i)*tg+alv*qg+gzicb(i) ! convect3
143 tg=tg+s*(ah0(i)-ahg)
144
145 !debug tc=tg-t0
146 tc=tg-273.15
147 denom=243.5+tc
148 denom=MAX(denom, 1.0) ! convect3
149
150 es=6.112*exp(17.67*tc/denom)
151
152 qg=eps*es/(p(i, icbs(i))-es*(1.-eps))
153
154 alv=lv0-clmcpv*(ticb(i)-273.15)
155
156 ! convect3: no approximation:
157 tp(i, icbs(i))=(ah0(i)-gz(i, icbs(i))-alv*qg) &
158 /(cpd+(cl-cpd)*qnk(i))
159
160 clw(i, icbs(i))=qnk(i)-qg
161 clw(i, icbs(i))=max(0.0, clw(i, icbs(i)))
162
163 ! convect3: (qg utilise au lieu du vrai mixing ratio rg)
164 tvp(i, icbs(i))=tp(i, icbs(i))*(1.+qg/eps-qnk(i)) !whole thing
165
166 end do
167
168 ! The following is only for convect3:
169
170 ! * icbs is the first level above the LCL:
171 ! if plcl<p(icb), then icbs=icb+1
172 ! if plcl>p(icb), then icbs=icb
173
174 ! * the routine above computes tvp from minorig to icbs (included).
175
176 ! * to compute buoybase (in cv30_trigger.F), both tvp(icb) and tvp(icb+1)
177 ! must be known. This is the case if icbs=icb+1, but not if icbs=icb.
178
179 ! * therefore, in the case icbs=icb, we compute tvp at level icb+1
180 ! (tvp at other levels will be computed in cv30_undilute2.F)
181
182 do i=1, klon
183 ticb(i)=t(i, icb(i)+1)
184 gzicb(i)=gz(i, icb(i)+1)
185 qsicb(i)=qs(i, icb(i)+1)
186 enddo
187
188 do i=1, klon
189 tg=ticb(i)
190 qg=qsicb(i) ! convect3
191 !debug alv=lv0-clmcpv*(ticb(i)-t0)
192 alv=lv0-clmcpv*(ticb(i)-273.15)
193
194 ! First iteration.
195
196 s=cpd*(1.-qnk(i))+cl*qnk(i) &
197 +alv*alv*qg/(rrv*ticb(i)*ticb(i)) ! convect3
198 s=1./s
199
200 ahg=cpd*tg+(cl-cpd)*qnk(i)*tg+alv*qg+gzicb(i) ! convect3
201 tg=tg+s*(ah0(i)-ahg)
202
203 !debug tc=tg-t0
204 tc=tg-273.15
205 denom=243.5+tc
206 denom=MAX(denom, 1.0) ! convect3
207
208 es=6.112*exp(17.67*tc/denom)
209
210 qg=eps*es/(p(i, icb(i)+1)-es*(1.-eps))
211
212 ! Second iteration.
213
214 ahg=cpd*tg+(cl-cpd)*qnk(i)*tg+alv*qg+gzicb(i) ! convect3
215 tg=tg+s*(ah0(i)-ahg)
216
217 !debug tc=tg-t0
218 tc=tg-273.15
219 denom=243.5+tc
220 denom=MAX(denom, 1.0) ! convect3
221
222 es=6.112*exp(17.67*tc/denom)
223
224 qg=eps*es/(p(i, icb(i)+1)-es*(1.-eps))
225
226 alv=lv0-clmcpv*(ticb(i)-273.15)
227
228 ! convect3: no approximation:
229 tp(i, icb(i)+1)=(ah0(i)-gz(i, icb(i)+1)-alv*qg) &
230 /(cpd+(cl-cpd)*qnk(i))
231
232 clw(i, icb(i)+1)=qnk(i)-qg
233 clw(i, icb(i)+1)=max(0.0, clw(i, icb(i)+1))
234
235 ! convect3: (qg utilise au lieu du vrai mixing ratio rg)
236 tvp(i, icb(i)+1)=tp(i, icb(i)+1)*(1.+qg/eps-qnk(i)) !whole thing
237 end do
238
239 end SUBROUTINE cv30_undilute1
240
241 end module cv30_undilute1_m

  ViewVC Help
Powered by ViewVC 1.1.21