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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 201 - (show annotations)
Mon Jun 6 17:42:15 2016 UTC (7 years, 10 months ago) by guez
File size: 6390 byte(s)
Removed intermediary objects of cv_thermo_m, access suphec_m
directly. Procedure cv_thermo disappeared, all objects are named
constants.

In cv_driver and below, limited extents of arrays to what is needed.

lv, cpn and th in cv30_compress were set at level nl + 1 but lv1, cpn1
and th1 are not defined at this level. This did not lead to an error
because values at nl + 1 were not used.

Removed test on ok_sync in phystokenc because it is not read at run
time. Printing min and max of output NetCDF variables is heavy and
archaic.

Used histwrite_phy in phytrac.

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

  ViewVC Help
Powered by ViewVC 1.1.21