1 |
jmc |
1.31 |
C $Header: /u/gcmpack/MITgcm/pkg/fizhi/fizhi_moist.F,v 1.30 2005/03/09 23:25:18 molod Exp $ |
2 |
molod |
1.1 |
C $Name: $ |
3 |
molod |
1.2 |
|
4 |
molod |
1.13 |
#include "FIZHI_OPTIONS.h" |
5 |
molod |
1.5 |
subroutine moistio (ndmoist,istrip,npcs, |
6 |
molod |
1.4 |
. lowlevel,midlevel,nltop,nsubmin,nsubmax,Lup, |
7 |
molod |
1.21 |
. pz,plz,plze,dpres,pkht,pkl,uz,vz,tz,qz,bi,bj,ntracerin,ptracer, |
8 |
|
|
. qqz,dumoist,dvmoist,dtmoist,dqmoist,cumfric, |
9 |
molod |
1.4 |
. im,jm,lm,ptop, |
10 |
molod |
1.1 |
. iras,rainlsp,rainconv,snowfall, |
11 |
|
|
. nswcld,cldtot_sw,cldras_sw,cldlsp_sw,nswlz,swlz, |
12 |
|
|
. nlwcld,cldtot_lw,cldras_lw,cldlsp_lw,nlwlz,lwlz, |
13 |
molod |
1.2 |
. lpnt,myid) |
14 |
molod |
1.1 |
|
15 |
molod |
1.8 |
implicit none |
16 |
|
|
|
17 |
molod |
1.2 |
#ifdef ALLOW_DIAGNOSTICS |
18 |
molod |
1.8 |
#include "SIZE.h" |
19 |
molod |
1.25 |
#include "DIAGNOSTICS_SIZE.h" |
20 |
|
|
#include "DIAGNOSTICS.h" |
21 |
molod |
1.2 |
#endif |
22 |
molod |
1.1 |
|
23 |
|
|
c Input Variables |
24 |
|
|
c --------------- |
25 |
molod |
1.8 |
integer im,jm,lm |
26 |
molod |
1.5 |
integer ndmoist,istrip,npcs |
27 |
molod |
1.21 |
integer bi,bj,ntracerin,ptracer |
28 |
molod |
1.4 |
integer lowlevel,midlevel,nltop,nsubmin,nsubmax,Lup |
29 |
molod |
1.13 |
_RL pz(im,jm),plz(im,jm,lm),plze(im,jm,lm+1),dpres(im,jm,lm) |
30 |
|
|
_RL pkht(im,jm,lm+1),pkl(im,jm,lm) |
31 |
molod |
1.21 |
_RL tz(im,jm,lm),qz(im,jm,lm,ntracerin) |
32 |
|
|
_RL uz(im,jm,lm),vz(im,jm,lm) |
33 |
molod |
1.13 |
_RL qqz(im,jm,lm) |
34 |
|
|
_RL dumoist(im,jm,lm),dvmoist(im,jm,lm) |
35 |
molod |
1.21 |
_RL dtmoist(im,jm,lm),dqmoist(im,jm,lm,ntracerin) |
36 |
|
|
logical cumfric |
37 |
molod |
1.13 |
_RL ptop |
38 |
molod |
1.5 |
integer iras |
39 |
molod |
1.13 |
_RL rainlsp(im,jm),rainconv(im,jm),snowfall(im,jm) |
40 |
molod |
1.5 |
integer nswcld,nswlz |
41 |
molod |
1.13 |
_RL cldlsp_sw(im,jm,lm),cldras_sw(im,jm,lm) |
42 |
|
|
_RL cldtot_sw(im,jm,lm),swlz(im,jm,lm) |
43 |
molod |
1.5 |
integer nlwcld,nlwlz |
44 |
molod |
1.13 |
_RL cldlsp_lw(im,jm,lm),cldras_lw(im,jm,lm) |
45 |
|
|
_RL cldtot_lw(im,jm,lm),lwlz(im,jm,lm) |
46 |
molod |
1.5 |
logical lpnt |
47 |
|
|
integer myid |
48 |
molod |
1.1 |
|
49 |
|
|
c Local Variables |
50 |
|
|
c --------------- |
51 |
molod |
1.4 |
integer ncrnd,nsecf |
52 |
molod |
1.1 |
|
53 |
molod |
1.13 |
_RL fracqq, dum |
54 |
molod |
1.4 |
integer snowcrit |
55 |
molod |
1.1 |
parameter (fracqq = 0.1) |
56 |
molod |
1.14 |
_RL one |
57 |
|
|
parameter (one=1.) |
58 |
molod |
1.1 |
|
59 |
molod |
1.13 |
_RL cldsr(im,jm,lm) |
60 |
|
|
_RL srcld(istrip,lm) |
61 |
molod |
1.1 |
|
62 |
molod |
1.13 |
_RL plev |
63 |
|
|
_RL cldnow,cldlsp_mem,cldlsp,cldras_mem,cldras |
64 |
|
|
_RL watnow,watmin,cldmin |
65 |
|
|
_RL cldprs(im,jm),cldtmp(im,jm) |
66 |
|
|
_RL cldhi (im,jm),cldlow(im,jm) |
67 |
|
|
_RL cldmid(im,jm),totcld(im,jm) |
68 |
|
|
|
69 |
|
|
_RL CLDLS(im,jm,lm) , CPEN(im,jm,lm) |
70 |
|
|
_RL tmpimjm(im,jm) |
71 |
|
|
_RL lsp_new(im,jm) |
72 |
|
|
_RL conv_new(im,jm) |
73 |
|
|
_RL snow_new(im,jm) |
74 |
molod |
1.1 |
|
75 |
molod |
1.13 |
_RL qqcolmin(im,jm) |
76 |
|
|
_RL qqcolmax(im,jm) |
77 |
molod |
1.1 |
integer levpbl(im,jm) |
78 |
|
|
|
79 |
|
|
c Gathered Arrays for Variable Cloud Base |
80 |
|
|
c --------------------------------------- |
81 |
molod |
1.13 |
_RL raincgath(im*jm) |
82 |
|
|
_RL pigather(im*jm) |
83 |
|
|
_RL thgather(im*jm,lm) |
84 |
|
|
_RL shgather(im*jm,lm) |
85 |
|
|
_RL pkzgather(im*jm,lm) |
86 |
|
|
_RL pkegather(im*jm,lm+1) |
87 |
|
|
_RL plzgather(im*jm,lm) |
88 |
|
|
_RL plegather(im*jm,lm+1) |
89 |
|
|
_RL dpgather(im*jm,lm) |
90 |
|
|
_RL tmpgather(im*jm,lm) |
91 |
|
|
_RL deltgather(im*jm,lm) |
92 |
|
|
_RL delqgather(im*jm,lm) |
93 |
molod |
1.21 |
_RL ugather(im*jm,lm,ntracerin+2-ptracer) |
94 |
|
|
_RL delugather(im*jm,lm,ntracerin+2-ptracer) |
95 |
molod |
1.13 |
_RL deltrnev(im*jm,lm) |
96 |
|
|
_RL delqrnev(im*jm,lm) |
97 |
molod |
1.1 |
|
98 |
|
|
integer nindeces(lm) |
99 |
|
|
integer pblindex(im*jm) |
100 |
|
|
integer levgather(im*jm) |
101 |
|
|
|
102 |
|
|
c Stripped Arrays |
103 |
|
|
c --------------- |
104 |
molod |
1.13 |
_RL saveth (istrip,lm) |
105 |
|
|
_RL saveq (istrip,lm) |
106 |
molod |
1.21 |
_RL saveu (istrip,lm,ntracerin+2-ptracer) |
107 |
|
|
_RL usubcl (istrip, ntracerin+2-ptracer) |
108 |
molod |
1.13 |
|
109 |
|
|
_RL ple(istrip,lm+1) |
110 |
|
|
_RL dp(istrip,lm) |
111 |
|
|
_RL TL(ISTRIP,lm) , SHL(ISTRIP,lm) |
112 |
|
|
_RL PL(ISTRIP,lm) , PLK(ISTRIP,lm) |
113 |
|
|
_RL PLKE(ISTRIP,lm+1) |
114 |
|
|
_RL TH(ISTRIP,lm) ,CVTH(ISTRIP,lm) |
115 |
|
|
_RL CVQ(ISTRIP,lm) |
116 |
molod |
1.21 |
_RL UL(ISTRIP,lm,ntracerin+2-ptracer) |
117 |
|
|
_RL cvu(istrip,lm,ntracerin+2-ptracer) |
118 |
molod |
1.13 |
_RL CLMAXO(ISTRIP,lm),CLBOTH(ISTRIP,lm) |
119 |
|
|
_RL CLSBTH(ISTRIP,lm) |
120 |
|
|
_RL TMP1(ISTRIP,lm), TMP2(ISTRIP,lm) |
121 |
|
|
_RL TMP3(ISTRIP,lm), TMP4(ISTRIP,lm+1) |
122 |
|
|
_RL TMP5(ISTRIP,lm+1) |
123 |
molod |
1.1 |
integer ITMP1(ISTRIP,lm), ITMP2(ISTRIP,lm) |
124 |
|
|
|
125 |
molod |
1.13 |
_RL PRECIP(ISTRIP), PCNET(ISTRIP) |
126 |
|
|
_RL SP(ISTRIP), PREP(ISTRIP) |
127 |
|
|
_RL PCPEN (ISTRIP,lm) |
128 |
molod |
1.1 |
integer pbl(istrip),depths(lm) |
129 |
|
|
|
130 |
molod |
1.13 |
_RL cldlz(istrip,lm), cldwater(im,jm,lm) |
131 |
|
|
_RL rhfrac(istrip), rhmin, pup, ppbl, rhcrit(istrip,lm) |
132 |
|
|
_RL offset, alpha, rasmax |
133 |
molod |
1.1 |
|
134 |
|
|
logical first |
135 |
|
|
logical lras |
136 |
molod |
1.13 |
_RL clfrac (istrip,lm) |
137 |
|
|
_RL cldmas (istrip,lm) |
138 |
|
|
_RL detrain(istrip,lm) |
139 |
|
|
_RL psubcld (istrip), psubcldg (im,jm) |
140 |
|
|
_RL psubcld_cnt(istrip), psubcldgc(im,jm) |
141 |
|
|
_RL rnd(lm/2) |
142 |
molod |
1.1 |
DATA FIRST /.TRUE./ |
143 |
|
|
|
144 |
molod |
1.4 |
integer imstp,nsubcl,nlras |
145 |
molod |
1.16 |
integer i,j,iloop,indx,indgath,l,nn,num,numdeps,nt |
146 |
molod |
1.13 |
_RL tmstp,tminv,sday,grav,alhl,cp,elocp,gamfac |
147 |
|
|
_RL rkappa,p0kappa,p0kinv,ptopkap,pcheck |
148 |
|
|
_RL tice,getcon,pi |
149 |
molod |
1.21 |
integer ntracer,ntracedim, ntracex |
150 |
molod |
1.1 |
|
151 |
|
|
C ********************************************************************** |
152 |
|
|
C **** INITIALIZATION **** |
153 |
|
|
C ********************************************************************** |
154 |
|
|
|
155 |
molod |
1.21 |
C Add U and V components to tracer array for cumulus friction |
156 |
|
|
|
157 |
|
|
if(cumfric) then |
158 |
|
|
ntracer = ntracerin + 2 |
159 |
|
|
else |
160 |
|
|
ntracer = ntracerin |
161 |
|
|
endif |
162 |
molod |
1.18 |
ntracedim= max(ntracer-ptracer,1) |
163 |
molod |
1.21 |
ntracex= ntracer-ptracer |
164 |
molod |
1.1 |
IMSTP = nsecf(NDMOIST) |
165 |
|
|
TMSTP = FLOAT(IMSTP) |
166 |
|
|
TMINV = 1. / TMSTP |
167 |
|
|
|
168 |
|
|
C Minimum Large-Scale Cloud Fraction at rhcrit |
169 |
|
|
alpha = 0.80 |
170 |
molod |
1.2 |
C Difference in fraction between SR and LS Threshold |
171 |
molod |
1.1 |
offset = 0.10 |
172 |
molod |
1.2 |
C Large-Scale Relative Humidity Threshold in PBL |
173 |
molod |
1.1 |
rhmin = 0.90 |
174 |
|
|
C Maximum Cloud Fraction associated with RAS |
175 |
|
|
rasmax = 1.00 |
176 |
|
|
|
177 |
|
|
nn = 3*3600.0/tmstp + 1 |
178 |
|
|
C Threshold for Cloud Fraction Memory |
179 |
|
|
cldmin = rasmax*(1.0-tmstp/3600.)**nn |
180 |
|
|
C Threshold for Cloud Liquid Water Memory |
181 |
|
|
watmin = 1.0e-8 |
182 |
|
|
|
183 |
|
|
SDAY = GETCON('SDAY') |
184 |
|
|
GRAV = GETCON('GRAVITY') |
185 |
|
|
ALHL = GETCON('LATENT HEAT COND') |
186 |
|
|
CP = GETCON('CP') |
187 |
|
|
ELOCP = GETCON('LATENT HEAT COND') / GETCON('CP') |
188 |
|
|
GAMFAC = GETCON('LATENT HEAT COND') * GETCON('EPS') * ELOCP |
189 |
|
|
. / GETCON('RGAS') |
190 |
|
|
RKAPPA = GETCON('KAPPA') |
191 |
|
|
P0KAPPA = 1000.0**RKAPPA |
192 |
|
|
P0KINV = 1. / P0KAPPA |
193 |
|
|
PTOPKAP = PTOP**RKAPPA |
194 |
|
|
tice = getcon('FREEZING-POINT') |
195 |
|
|
PI = 4.*atan(1.) |
196 |
|
|
|
197 |
molod |
1.4 |
c Determine Total number of Random Clouds to Check |
198 |
molod |
1.1 |
c --------------------------------------------- |
199 |
molod |
1.24 |
ncrnd = (lm-nltop+1)/2 |
200 |
molod |
1.1 |
|
201 |
molod |
1.19 |
if(first .and. myid.eq.1 .and. bi.eq.1 ) then |
202 |
molod |
1.1 |
print * |
203 |
molod |
1.4 |
print *,'Top Level Allowed for Convection : ',nltop |
204 |
|
|
print *,' Highest Sub-Cloud Level: ',nsubmax |
205 |
|
|
print *,' Lowest Sub-Cloud Level: ',nsubmin |
206 |
molod |
1.1 |
print *,' Total Number of Random Clouds: ',ncrnd |
207 |
|
|
print * |
208 |
|
|
first = .false. |
209 |
|
|
endif |
210 |
|
|
|
211 |
|
|
c And now find PBL depth - the level where qq = fracqq * qq at surface |
212 |
|
|
c -------------------------------------------------------------------- |
213 |
|
|
do j = 1,jm |
214 |
|
|
do i = 1,im |
215 |
|
|
qqcolmin(i,j) = qqz(i,j,lm)*fracqq |
216 |
|
|
qqcolmax(i,j) = qqz(i,j,lm) |
217 |
|
|
levpbl(i,j) = lm+1 |
218 |
|
|
enddo |
219 |
|
|
enddo |
220 |
|
|
|
221 |
|
|
DO L = lm-1,1,-1 |
222 |
|
|
DO j = 1,jm |
223 |
|
|
DO i = 1,im |
224 |
|
|
IF((qqz(i,j,l).gt.qqcolmax(i,j)) |
225 |
|
|
1 .and.(levpbl(i,j).eq.lm+1))then |
226 |
|
|
qqcolmax(i,j) = qqz(i,j,l) |
227 |
|
|
qqcolmin(i,j) = fracqq*qqcolmax(i,j) |
228 |
|
|
endif |
229 |
|
|
if((qqz(i,j,l).lt.qqcolmin(i,j)) |
230 |
|
|
1 .and.(levpbl(i,j).eq.lm+1))levpbl(i,j)=L+1 |
231 |
|
|
enddo |
232 |
|
|
enddo |
233 |
|
|
enddo |
234 |
|
|
|
235 |
|
|
do j = 1,jm |
236 |
|
|
do i = 1,im |
237 |
|
|
if(levpbl(i,j).gt.nsubmin) levpbl(i,j) = nsubmin |
238 |
|
|
if(levpbl(i,j).lt.nsubmax) levpbl(i,j) = nsubmax |
239 |
|
|
enddo |
240 |
|
|
enddo |
241 |
|
|
|
242 |
|
|
|
243 |
|
|
c Set up the array of indeces of subcloud levels for the gathering |
244 |
|
|
c ---------------------------------------------------------------- |
245 |
molod |
1.16 |
indx = 0 |
246 |
molod |
1.1 |
do L = nsubmin,nltop,-1 |
247 |
|
|
do j = 1,jm |
248 |
|
|
do i = 1,im |
249 |
|
|
if(levpbl(i,j).eq.L) then |
250 |
molod |
1.16 |
indx = indx + 1 |
251 |
|
|
pblindex(indx) = (j-1)*im + i |
252 |
molod |
1.1 |
endif |
253 |
|
|
enddo |
254 |
|
|
enddo |
255 |
|
|
enddo |
256 |
|
|
|
257 |
molod |
1.16 |
do indx = 1,im*jm |
258 |
|
|
levgather(indx) = levpbl(pblindex(indx),1) |
259 |
|
|
pigather(indx) = pz(pblindex(indx),1) |
260 |
|
|
pkegather(indx,lm+1) = pkht(pblindex(indx),1,lm+1) |
261 |
|
|
plegather(indx,lm+1) = plze(pblindex(indx),1,lm+1) |
262 |
molod |
1.1 |
enddo |
263 |
|
|
|
264 |
|
|
do L = 1,lm |
265 |
molod |
1.16 |
do indx = 1,im*jm |
266 |
|
|
thgather(indx,L) = tz(pblindex(indx),1,L) |
267 |
|
|
shgather(indx,L) = qz(pblindex(indx),1,L,1) |
268 |
|
|
pkegather(indx,L) = pkht(pblindex(indx),1,L) |
269 |
|
|
pkzgather(indx,L) = pkl(pblindex(indx),1,L) |
270 |
|
|
plegather(indx,L) = plze(pblindex(indx),1,L) |
271 |
|
|
plzgather(indx,L) = plz(pblindex(indx),1,L) |
272 |
|
|
dpgather(indx,L) = dpres(pblindex(indx),1,L) |
273 |
molod |
1.1 |
enddo |
274 |
|
|
enddo |
275 |
molod |
1.21 |
C General Tracers |
276 |
|
|
C---------------- |
277 |
|
|
do nt = 1,ntracerin-ptracer |
278 |
|
|
do L = 1,lm |
279 |
|
|
do indx = 1,im*jm |
280 |
|
|
ugather(indx,L,nt) = qz(pblindex(indx),1,L,nt+ptracer) |
281 |
|
|
enddo |
282 |
|
|
enddo |
283 |
|
|
enddo |
284 |
|
|
|
285 |
|
|
if(cumfric) then |
286 |
|
|
C Cumulus Friction - load u and v wind components into tracer array |
287 |
|
|
C------------------------------------------------------------------ |
288 |
|
|
do L = 1,lm |
289 |
|
|
do indx = 1,im*jm |
290 |
|
|
ugather(indx,L,ntracerin-ptracer+1) = uz(pblindex(indx),1,L) |
291 |
|
|
ugather(indx,L,ntracerin-ptracer+2) = vz(pblindex(indx),1,L) |
292 |
|
|
enddo |
293 |
|
|
enddo |
294 |
|
|
|
295 |
|
|
endif |
296 |
molod |
1.1 |
|
297 |
|
|
c bump the counter for number of calls to convection |
298 |
|
|
c -------------------------------------------------- |
299 |
|
|
iras = iras + 1 |
300 |
|
|
if( iras.ge.1e9 ) iras = 1 |
301 |
|
|
|
302 |
|
|
c select the 'random' cloud detrainment levels for RAS |
303 |
|
|
c ---------------------------------------------------- |
304 |
|
|
call rndcloud(iras,ncrnd,rnd,myid) |
305 |
|
|
|
306 |
|
|
do l=1,lm |
307 |
|
|
do j=1,jm |
308 |
|
|
do i=1,im |
309 |
molod |
1.21 |
dumoist(i,j,l) = 0. |
310 |
|
|
dvmoist(i,j,l) = 0. |
311 |
molod |
1.1 |
dtmoist(i,j,l) = 0. |
312 |
molod |
1.21 |
do nt = 1,ntracerin |
313 |
molod |
1.1 |
dqmoist(i,j,l,nt) = 0. |
314 |
|
|
enddo |
315 |
|
|
enddo |
316 |
|
|
enddo |
317 |
|
|
enddo |
318 |
|
|
|
319 |
|
|
C*********************************************************************** |
320 |
|
|
C **** LOOP OVER NPCS PEICES **** |
321 |
|
|
C ********************************************************************** |
322 |
|
|
|
323 |
|
|
DO 1000 NN = 1,NPCS |
324 |
|
|
|
325 |
|
|
C ********************************************************************** |
326 |
|
|
C **** VARIABLE INITIALIZATION **** |
327 |
|
|
C ********************************************************************** |
328 |
|
|
|
329 |
|
|
CALL STRIP ( pigather, SP ,im*jm,ISTRIP,1 ,NN ) |
330 |
|
|
CALL STRIP ( pkzgather, PLK ,im*jm,ISTRIP,lm,NN ) |
331 |
molod |
1.6 |
CALL STRIP ( pkegather, PLKE ,im*jm,ISTRIP,lm+1,NN ) |
332 |
|
|
CALL STRIP ( plzgather, PL ,im*jm,ISTRIP,lm,NN ) |
333 |
|
|
CALL STRIP ( plegather, PLE ,im*jm,ISTRIP,lm+1,NN ) |
334 |
|
|
CALL STRIP ( dpgather, dp ,im*jm,ISTRIP,lm,NN ) |
335 |
molod |
1.1 |
CALL STRIP ( thgather, TH ,im*jm,ISTRIP,lm,NN ) |
336 |
|
|
CALL STRIP ( shgather, SHL ,im*jm,ISTRIP,lm,NN ) |
337 |
|
|
CALL STRINT( levgather, pbl ,im*jm,ISTRIP,1 ,NN ) |
338 |
|
|
|
339 |
molod |
1.21 |
do nt = 1,ntracer-ptracer |
340 |
|
|
call strip ( ugather(1,1,nt), ul(1,1,nt),im*jm,istrip,lm,nn ) |
341 |
|
|
enddo |
342 |
molod |
1.1 |
|
343 |
|
|
C ********************************************************************** |
344 |
|
|
C **** SETUP FOR RAS CUMULUS PARAMETERIZATION **** |
345 |
|
|
C ********************************************************************** |
346 |
|
|
|
347 |
|
|
DO L = 1,lm |
348 |
|
|
DO I = 1,ISTRIP |
349 |
|
|
TH(I,L) = TH(I,L) * P0KAPPA |
350 |
|
|
CLMAXO(I,L) = 0. |
351 |
|
|
CLBOTH(I,L) = 0. |
352 |
|
|
cldmas(I,L) = 0. |
353 |
|
|
detrain(I,L) = 0. |
354 |
|
|
ENDDO |
355 |
|
|
ENDDO |
356 |
|
|
|
357 |
|
|
do L = 1,lm |
358 |
|
|
depths(L) = 0 |
359 |
|
|
enddo |
360 |
|
|
|
361 |
|
|
numdeps = 0 |
362 |
|
|
do L = nsubmin,nltop,-1 |
363 |
|
|
nindeces(L) = 0 |
364 |
|
|
do i = 1,istrip |
365 |
|
|
if(pbl(i).eq.L) nindeces(L) = nindeces(L) + 1 |
366 |
|
|
enddo |
367 |
|
|
if(nindeces(L).gt.0) then |
368 |
|
|
numdeps = numdeps + 1 |
369 |
|
|
depths(numdeps) = L |
370 |
|
|
endif |
371 |
|
|
enddo |
372 |
|
|
|
373 |
|
|
C Initiate a do-loop around RAS for the number of different |
374 |
|
|
C sub-cloud layer depths in this strip |
375 |
|
|
C --If all subcloud depths are the same, execute loop once |
376 |
|
|
C Otherwise loop over different subcloud layer depths |
377 |
|
|
|
378 |
|
|
num = 1 |
379 |
|
|
DO iloop = 1,numdeps |
380 |
|
|
|
381 |
|
|
nsubcl = depths(iloop) |
382 |
|
|
|
383 |
|
|
c Compute sub-cloud values for Temperature and Spec.Hum. |
384 |
|
|
c ------------------------------------------------------ |
385 |
|
|
DO 600 I=num,num+nindeces(nsubcl)-1 |
386 |
|
|
TMP1(I,2) = 0. |
387 |
|
|
TMP1(I,3) = 0. |
388 |
|
|
600 CONTINUE |
389 |
|
|
|
390 |
|
|
NLRAS = NSUBCL - NLTOP + 1 |
391 |
|
|
DO 601 L=NSUBCL,lm |
392 |
|
|
DO 602 I=num,num+nindeces(nsubcl)-1 |
393 |
|
|
TMP1(I,2) = TMP1(I,2) + (PLE(I,L+1)-PLE(I,L))*TH (I,L)/sp(i) |
394 |
|
|
TMP1(I,3) = TMP1(I,3) + (PLE(I,L+1)-PLE(I,L))*SHL(I,L)/sp(i) |
395 |
|
|
602 CONTINUE |
396 |
|
|
601 CONTINUE |
397 |
|
|
DO 603 I=num,num+nindeces(nsubcl)-1 |
398 |
|
|
TMP1(I,4) = 1. / ( (PLE(I,lm+1)-PLE(I,NSUBCL))/sp(I) ) |
399 |
|
|
TH(I,NSUBCL) = TMP1(I,2)*TMP1(I,4) |
400 |
|
|
SHL(I,NSUBCL) = TMP1(I,3)*TMP1(I,4) |
401 |
|
|
603 CONTINUE |
402 |
|
|
|
403 |
|
|
c Save initial value of tracers and compute sub-cloud value |
404 |
|
|
c --------------------------------------------------------- |
405 |
molod |
1.21 |
DO NT = 1,ntracer-ptracer |
406 |
|
|
do L = 1,lm |
407 |
|
|
do i = num,num+nindeces(nsubcl)-1 |
408 |
|
|
saveu(i,L,nt) = ul(i,L,nt) |
409 |
|
|
enddo |
410 |
|
|
enddo |
411 |
|
|
DO I=num,num+nindeces(nsubcl)-1 |
412 |
|
|
TMP1(I,2) = 0. |
413 |
|
|
ENDDO |
414 |
|
|
DO L=NSUBCL,lm |
415 |
|
|
DO I=num,num+nindeces(nsubcl)-1 |
416 |
|
|
TMP1(I,2) = TMP1(I,2)+(PLE(I,L+1)-PLE(I,L))*UL(I,L,NT)/sp(i) |
417 |
|
|
ENDDO |
418 |
|
|
ENDDO |
419 |
|
|
DO I=num,num+nindeces(nsubcl)-1 |
420 |
|
|
UL(I,NSUBCL,NT) = TMP1(I,2)*TMP1(I,4) |
421 |
|
|
usubcl(i,nt) = ul(i,nsubcl,nt) |
422 |
|
|
ENDDO |
423 |
|
|
ENDDO |
424 |
molod |
1.1 |
|
425 |
|
|
c Compute Pressure Arrays for RAS |
426 |
|
|
c ------------------------------- |
427 |
|
|
DO 111 L=1,lm |
428 |
|
|
DO 112 I=num,num+nindeces(nsubcl)-1 |
429 |
|
|
TMP4(I,L) = PLE(I,L) |
430 |
|
|
112 CONTINUE |
431 |
|
|
111 CONTINUE |
432 |
|
|
DO I=num,num+nindeces(nsubcl)-1 |
433 |
|
|
TMP5(I,1) = PTOPKAP / P0KAPPA |
434 |
|
|
ENDDO |
435 |
|
|
DO L=2,lm |
436 |
|
|
DO I=num,num+nindeces(nsubcl)-1 |
437 |
molod |
1.6 |
TMP5(I,L) = PLKE(I,L)*P0KINV |
438 |
molod |
1.1 |
ENDDO |
439 |
|
|
ENDDO |
440 |
|
|
DO I=num,num+nindeces(nsubcl)-1 |
441 |
|
|
TMP4(I,lm+1) = PLE (I,lm+1) |
442 |
molod |
1.6 |
TMP5(I,lm+1) = PLKE(I,lm+1)*P0KINV |
443 |
molod |
1.1 |
ENDDO |
444 |
|
|
DO 113 I=num,num+nindeces(nsubcl)-1 |
445 |
|
|
TMP4(I,NSUBCL+1) = PLE (I,lm+1) |
446 |
molod |
1.6 |
TMP5(I,NSUBCL+1) = PLKE(I,lm+1)*P0KINV |
447 |
molod |
1.1 |
113 CONTINUE |
448 |
|
|
|
449 |
|
|
do i=num,num+nindeces(nsubcl)-1 |
450 |
|
|
C Temperature at top of sub-cloud layer |
451 |
|
|
tmp2(i,1) = TH(i,NSUBCL) * PLKE(i,NSUBCL)/P0KAPPA |
452 |
|
|
C Pressure at top of sub-cloud layer |
453 |
|
|
tmp2(i,2) = tmp4(i,nsubcl) |
454 |
|
|
enddo |
455 |
|
|
|
456 |
|
|
C CHANGED THIS: no RH requirement for RAS |
457 |
|
|
c call vqsat ( tmp2(num,1),tmp2(num,2),tmp2(num,3), |
458 |
|
|
c . dum,.false.,nindeces(nsubcl) ) |
459 |
|
|
c do i=num,num+nindeces(nsubcl)-1 |
460 |
|
|
c rh = SHL(I,NSUBCL) / tmp2(i,3) |
461 |
|
|
c if (rh .le. 0.85) then |
462 |
|
|
c rhfrac(i) = 0. |
463 |
|
|
c else if (rh .ge. 0.95) then |
464 |
|
|
c rhfrac(i) = 1. |
465 |
|
|
c else |
466 |
|
|
c rhfrac(i) = (rh-0.85)*10. |
467 |
|
|
c endif |
468 |
|
|
c enddo |
469 |
|
|
do i=num,num+nindeces(nsubcl)-1 |
470 |
|
|
rhfrac(i) = 1. |
471 |
|
|
enddo |
472 |
|
|
|
473 |
|
|
C Compute RH threshold for Large-scale condensation |
474 |
|
|
C Used in Slingo-Ritter clouds as well - define offset between SR and LS |
475 |
|
|
|
476 |
|
|
C Top level of atan func above this rh_threshold = rhmin |
477 |
|
|
pup = 600. |
478 |
|
|
do i=num,num+nindeces(nsubcl)-1 |
479 |
molod |
1.2 |
do L = nsubcl, lm |
480 |
|
|
rhcrit(i,L) = 1. |
481 |
|
|
enddo |
482 |
|
|
do L = 1, nsubcl-1 |
483 |
molod |
1.6 |
pcheck = pl(i,L) |
484 |
molod |
1.2 |
if (pcheck .le. pup) then |
485 |
|
|
rhcrit(i,L) = rhmin |
486 |
|
|
else |
487 |
molod |
1.6 |
ppbl = pl(i,nsubcl) |
488 |
molod |
1.2 |
rhcrit(i,L) = rhmin + (1.-rhmin)/(19.) * |
489 |
|
|
. ((atan( (2.*(pcheck-pup)/(ppbl-pup)-1.) * |
490 |
molod |
1.1 |
. tan(20.*pi/21.-0.5*pi) ) |
491 |
|
|
. + 0.5*pi) * 21./pi - 1.) |
492 |
molod |
1.2 |
endif |
493 |
|
|
enddo |
494 |
molod |
1.1 |
enddo |
495 |
|
|
|
496 |
|
|
c Save Initial Values of Temperature and Specific Humidity |
497 |
|
|
c -------------------------------------------------------- |
498 |
|
|
do L = 1,lm |
499 |
|
|
do i = num,num+nindeces(nsubcl)-1 |
500 |
|
|
saveth(i,L) = th (i,L) |
501 |
|
|
saveq (i,L) = shl(i,L) |
502 |
|
|
PCPEN (i,L) = 0. |
503 |
|
|
CLFRAC(i,L) = 0. |
504 |
|
|
enddo |
505 |
|
|
enddo |
506 |
|
|
|
507 |
|
|
CALL RAS ( NN,istrip,nindeces(nsubcl),NLRAS,NLTOP,lm,TMSTP |
508 |
molod |
1.21 |
1, UL(num,1,1),ntracedim,ntracex,TH(num,NLTOP),SHL(num,NLTOP) |
509 |
molod |
1.1 |
2, TMP4(num,NLTOP), TMP5(num,NLTOP),rnd, ncrnd, PCPEN(num,NLTOP) |
510 |
|
|
3, CLBOTH(num,NLTOP), CLFRAC(num,NLTOP) |
511 |
|
|
4, cldmas(num,nltop), detrain(num,nltop) |
512 |
|
|
8, cp,grav,rkappa,alhl,rhfrac(num),rasmax ) |
513 |
|
|
|
514 |
|
|
c Compute Diagnostic CLDMAS in RAS Subcloud Layers |
515 |
|
|
c ------------------------------------------------ |
516 |
|
|
do L=nsubcl,lm |
517 |
|
|
do I=num,num+nindeces(nsubcl)-1 |
518 |
molod |
1.6 |
dum = dp(i,L)/(ple(i,lm+1)-ple(i,nsubcl)) |
519 |
molod |
1.1 |
cldmas(i,L) = cldmas(i,L-1) - dum*cldmas(i,nsubcl-1) |
520 |
|
|
enddo |
521 |
|
|
enddo |
522 |
|
|
|
523 |
|
|
c Update Theta and Moisture due to RAS |
524 |
|
|
c ------------------------------------ |
525 |
|
|
DO L=1,nsubcl |
526 |
|
|
DO I=num,num+nindeces(nsubcl)-1 |
527 |
|
|
CVTH(I,L) = (TH (I,L) - saveth(i,l)) |
528 |
|
|
CVQ (I,L) = (SHL(I,L) - saveq (i,l)) |
529 |
|
|
ENDDO |
530 |
|
|
ENDDO |
531 |
|
|
DO L=nsubcl+1,lm |
532 |
|
|
DO I=num,num+nindeces(nsubcl)-1 |
533 |
|
|
CVTH(I,L) = cvth(i,nsubcl) |
534 |
|
|
CVQ (I,L) = cvq (i,nsubcl) |
535 |
|
|
ENDDO |
536 |
|
|
ENDDO |
537 |
|
|
|
538 |
|
|
DO L=nsubcl+1,lm |
539 |
|
|
DO I=num,num+nindeces(nsubcl)-1 |
540 |
|
|
TH (I,L) = saveth(i,l) + cvth(i,l) |
541 |
|
|
SHL(I,L) = saveq (i,l) + cvq (i,l) |
542 |
|
|
ENDDO |
543 |
|
|
ENDDO |
544 |
|
|
DO L=1,lm |
545 |
|
|
DO I=num,num+nindeces(nsubcl)-1 |
546 |
|
|
CVTH(I,L) = CVTH(I,L) *P0KINV*SP(I)*tminv |
547 |
|
|
CVQ (I,L) = CVQ (I,L) *SP(I)*tminv |
548 |
|
|
ENDDO |
549 |
|
|
ENDDO |
550 |
|
|
|
551 |
|
|
c Compute Tracer Tendency due to RAS |
552 |
|
|
c ---------------------------------- |
553 |
molod |
1.21 |
do nt = 1,ntracer-ptracer |
554 |
|
|
DO L=1,nsubcl-1 |
555 |
|
|
DO I=num,num+nindeces(nsubcl)-1 |
556 |
|
|
CVU(I,L,nt) = ( UL(I,L,nt)-saveu(i,l,nt) )*sp(i)*tminv |
557 |
|
|
ENDDO |
558 |
|
|
ENDDO |
559 |
|
|
DO L=nsubcl,lm |
560 |
|
|
DO I=num,num+nindeces(nsubcl)-1 |
561 |
|
|
if( usubcl(i,nt).ne.0.0 ) then |
562 |
|
|
cvu(i,L,nt) = ( ul(i,nsubcl,nt)-usubcl(i,nt) ) * |
563 |
|
|
. ( saveu(i,L,nt)/usubcl(i,nt) )*sp(i)*tminv |
564 |
|
|
else |
565 |
|
|
cvu(i,L,nt) = 0.0 |
566 |
|
|
endif |
567 |
|
|
ENDDO |
568 |
|
|
ENDDO |
569 |
|
|
enddo |
570 |
molod |
1.1 |
|
571 |
|
|
c Compute Diagnostic PSUBCLD (Subcloud Layer Pressure) |
572 |
|
|
c ---------------------------------------------------- |
573 |
|
|
do i=num,num+nindeces(nsubcl)-1 |
574 |
|
|
lras = .false. |
575 |
|
|
do L=nltop,nsubcl |
576 |
|
|
if( cvq(i,L).ne.0.0 ) lras = .true. |
577 |
|
|
enddo |
578 |
|
|
psubcld (i) = 0.0 |
579 |
|
|
psubcld_cnt(i) = 0.0 |
580 |
|
|
if( lras ) then |
581 |
|
|
psubcld (i) = sp(i)+ptop-ple(i,nsubcl) |
582 |
|
|
psubcld_cnt(i) = 1.0 |
583 |
|
|
endif |
584 |
|
|
enddo |
585 |
|
|
|
586 |
|
|
|
587 |
|
|
C End of subcloud layer depth loop (iloop) |
588 |
|
|
|
589 |
|
|
num = num+nindeces(nsubcl) |
590 |
|
|
|
591 |
|
|
ENDDO |
592 |
|
|
|
593 |
|
|
C ********************************************************************** |
594 |
|
|
C **** TENDENCY UPDATES **** |
595 |
|
|
C **** (Keep 'Gathered' tendencies in 'gather' arrays now) **** |
596 |
|
|
C ********************************************************************** |
597 |
|
|
|
598 |
|
|
call paste( CVTH,deltgather,istrip,im*jm,lm,NN ) |
599 |
|
|
call paste( CVQ,delqgather,istrip,im*jm,lm,NN ) |
600 |
molod |
1.21 |
do nt = 1,ntracer-ptracer |
601 |
|
|
call paste( CVU(1,1,nt),delugather(1,1,nt),istrip,im*jm,lm,NN ) |
602 |
|
|
enddo |
603 |
molod |
1.1 |
|
604 |
|
|
C ********************************************************************** |
605 |
|
|
C And now paste some arrays for filling diagnostics |
606 |
|
|
C (use pkegather to hold detrainment and tmpgather for cloud mass flux) |
607 |
|
|
C ********************************************************************** |
608 |
|
|
|
609 |
|
|
if(icldmas .gt.0) call paste( cldmas,tmpgather,istrip,im*jm,lm,NN) |
610 |
|
|
if(idtrain .gt.0) call paste(detrain,pkegather,istrip,im*jm,lm,NN) |
611 |
|
|
if(ipsubcld.gt.0) then |
612 |
|
|
call paste(psubcld ,psubcldg ,istrip,im*jm,1,NN) |
613 |
|
|
call paste(psubcld_cnt,psubcldgc,istrip,im*jm,1,NN) |
614 |
|
|
endif |
615 |
|
|
|
616 |
|
|
C ********************************************************************* |
617 |
|
|
C **** RE-EVAPORATION OF PENETRATING CONVECTIVE RAIN **** |
618 |
|
|
C ********************************************************************* |
619 |
|
|
|
620 |
|
|
CALL STRIP ( thgather,TH ,im*jm,ISTRIP,lm,NN) |
621 |
|
|
CALL STRIP ( shgather,SHL,im*jm,ISTRIP,lm,NN) |
622 |
|
|
DO L=1,lm |
623 |
|
|
DO I=1,ISTRIP |
624 |
|
|
TH(I,L) = TH(I,L) + CVTH(I,L)*tmstp/SP(I) |
625 |
|
|
SHL(I,L) = SHL(I,L) + CVQ(I,L)*tmstp/SP(I) |
626 |
|
|
TL(I,L) = TH(I,L)*PLK(I,L) |
627 |
|
|
saveth(I,L) = th(I,L) |
628 |
|
|
saveq (I,L) = SHL(I,L) |
629 |
|
|
ENDDO |
630 |
|
|
ENDDO |
631 |
|
|
|
632 |
molod |
1.6 |
CALL RNEVP (NN,ISTRIP,lm,TL,SHL,PCPEN,PL,CLFRAC,SP,DP,PLKE, |
633 |
molod |
1.1 |
. PLK,TH,TMP1,TMP2,TMP3,ITMP1,ITMP2,PCNET,PRECIP, |
634 |
molod |
1.14 |
. CLSBTH,TMSTP,one,cp,grav,alhl,gamfac,cldlz,rhcrit,offset,alpha) |
635 |
molod |
1.1 |
|
636 |
|
|
C ********************************************************************** |
637 |
|
|
C **** TENDENCY UPDATES **** |
638 |
|
|
C ********************************************************************** |
639 |
|
|
|
640 |
|
|
DO L=1,lm |
641 |
|
|
|
642 |
|
|
DO I =1,ISTRIP |
643 |
|
|
TMP1(I,L) = sp(i) * (SHL(I,L)-saveq(I,L)) * tminv |
644 |
|
|
ENDDO |
645 |
|
|
CALL PSTBMP(TMP1(1,L),delqgather(1,L),ISTRIP,im*jm,1,NN) |
646 |
|
|
|
647 |
|
|
DO I =1,ISTRIP |
648 |
|
|
TMP1(I,L) = sp(i) * ((TL(I,L)/PLK(I,L))-saveth(i,l)) * tminv |
649 |
|
|
ENDDO |
650 |
|
|
CALL PSTBMP(TMP1(1,L),deltgather(1,L),ISTRIP,im*jm,1,NN) |
651 |
|
|
|
652 |
|
|
C Paste rain evap tendencies into arrays for diagnostic output |
653 |
|
|
c ------------------------------------------------------------ |
654 |
|
|
if(idtls.gt.0)then |
655 |
|
|
DO I =1,ISTRIP |
656 |
|
|
TMP1(I,L) = ((TL(I,L)/PLK(I,L))-saveth(i,l))*plk(i,l)*sday*tminv |
657 |
|
|
ENDDO |
658 |
|
|
call paste(tmp1(1,L),deltrnev(1,L),istrip,im*jm,1,NN) |
659 |
|
|
endif |
660 |
|
|
|
661 |
|
|
if(idqls.gt.0)then |
662 |
|
|
DO I =1,ISTRIP |
663 |
|
|
TMP1(I,L) = (SHL(I,L)-saveq(I,L)) * 1000. * sday * tminv |
664 |
|
|
ENDDO |
665 |
|
|
call paste(tmp1(1,L),delqrnev(1,L),istrip,im*jm,1,NN) |
666 |
|
|
endif |
667 |
|
|
|
668 |
|
|
ENDDO |
669 |
|
|
|
670 |
|
|
C ********************************************************************* |
671 |
|
|
C Add Non-Precipitating Clouds where the relative |
672 |
|
|
C humidity is less than 100% |
673 |
|
|
C Apply Cloud Top Entrainment Instability |
674 |
|
|
C ********************************************************************* |
675 |
|
|
|
676 |
|
|
do L=1,lm |
677 |
|
|
do i=1,istrip |
678 |
|
|
srcld(i,L) = -clsbth(i,L) |
679 |
|
|
enddo |
680 |
|
|
enddo |
681 |
|
|
|
682 |
|
|
call srclouds (saveth,saveq,plk,pl,plke,clsbth,cldlz,istrip,lm, |
683 |
|
|
. rhcrit,offset,alpha) |
684 |
|
|
|
685 |
|
|
do L=1,lm |
686 |
|
|
do i=1,istrip |
687 |
|
|
srcld(i,L) = srcld(i,L)+clsbth(i,L) |
688 |
|
|
enddo |
689 |
|
|
enddo |
690 |
|
|
|
691 |
|
|
C ********************************************************************* |
692 |
|
|
C **** PASTE CLOUD AMOUNTS **** |
693 |
|
|
C ********************************************************************* |
694 |
|
|
|
695 |
|
|
call paste ( srcld, cldsr,istrip,im*jm,lm,nn ) |
696 |
|
|
call paste ( cldlz,cldwater,istrip,im*jm,lm,nn ) |
697 |
|
|
call paste ( clsbth, cldls,istrip,im*jm,lm,nn ) |
698 |
|
|
call paste ( clboth, cpen ,istrip,im*jm,lm,nn ) |
699 |
|
|
|
700 |
|
|
c compute Total Accumulated Precip for Landsurface Model |
701 |
|
|
c ------------------------------------------------------ |
702 |
|
|
do i = 1,istrip |
703 |
|
|
C Initialize Rainlsp, Rainconv and Snowfall |
704 |
|
|
tmp1(i,1) = 0.0 |
705 |
|
|
tmp1(i,2) = 0.0 |
706 |
|
|
tmp1(i,3) = 0.0 |
707 |
|
|
enddo |
708 |
|
|
|
709 |
|
|
do i = 1,istrip |
710 |
|
|
prep(i) = PRECIP(I) + PCNET(I) |
711 |
|
|
tmp1(i,1) = PRECIP(I) |
712 |
|
|
tmp1(i,2) = pcnet(i) |
713 |
|
|
enddo |
714 |
|
|
c |
715 |
|
|
c check whether there is snow |
716 |
|
|
c------------------------------------------------------- |
717 |
|
|
c snow algorthm: |
718 |
|
|
c if temperature profile from the surface level to 700 mb |
719 |
|
|
c uniformaly c below zero, then precipitation (total) is |
720 |
|
|
c snowfall. Else there is no snow. |
721 |
|
|
c------------------------------------------------------- |
722 |
|
|
|
723 |
|
|
do i = 1,istrip |
724 |
|
|
snowcrit=0 |
725 |
|
|
do l=lup,lm |
726 |
|
|
if (saveth(i,l)*plk(i,l).le. tice ) then |
727 |
|
|
snowcrit=snowcrit+1 |
728 |
|
|
endif |
729 |
|
|
enddo |
730 |
|
|
if (snowcrit .eq. (lm-lup+1)) then |
731 |
|
|
tmp1(i,3) = prep(i) |
732 |
|
|
tmp1(i,1)=0.0 |
733 |
|
|
tmp1(i,2)=0.0 |
734 |
|
|
endif |
735 |
|
|
enddo |
736 |
|
|
|
737 |
|
|
CALL paste (tmp1(1,1), lsp_new,ISTRIP,im*jm,1,NN) |
738 |
|
|
CALL paste (tmp1(1,2),conv_new,ISTRIP,im*jm,1,NN) |
739 |
|
|
CALL paste (tmp1(1,3),snow_new,ISTRIP,im*jm,1,NN) |
740 |
|
|
|
741 |
|
|
if(iprecon.gt.0) then |
742 |
|
|
CALL paste (pcnet,raincgath,ISTRIP,im*jm,1,NN) |
743 |
|
|
endif |
744 |
|
|
|
745 |
|
|
C ********************************************************************* |
746 |
|
|
C **** End Major Stripped Region **** |
747 |
|
|
C ********************************************************************* |
748 |
|
|
|
749 |
|
|
1000 CONTINUE |
750 |
|
|
|
751 |
|
|
C Large Scale Rainfall, Conv rain, and snowfall |
752 |
|
|
c --------------------------------------------- |
753 |
|
|
call back2grd ( lsp_new,pblindex, lsp_new,im*jm) |
754 |
|
|
call back2grd (conv_new,pblindex,conv_new,im*jm) |
755 |
|
|
call back2grd (snow_new,pblindex,snow_new,im*jm) |
756 |
|
|
|
757 |
|
|
if(iprecon.gt.0) then |
758 |
|
|
call back2grd (raincgath,pblindex,raincgath,im*jm) |
759 |
|
|
endif |
760 |
|
|
|
761 |
|
|
c Subcloud Layer Pressure |
762 |
|
|
c ----------------------- |
763 |
|
|
if(ipsubcld.gt.0) then |
764 |
|
|
call back2grd (psubcldg ,pblindex,psubcldg ,im*jm) |
765 |
|
|
call back2grd (psubcldgc,pblindex,psubcldgc,im*jm) |
766 |
|
|
endif |
767 |
|
|
|
768 |
|
|
do L = 1,lm |
769 |
|
|
C Delta theta,q, convective, max and ls clouds |
770 |
|
|
c -------------------------------------------- |
771 |
|
|
call back2grd (deltgather(1,L),pblindex, dtmoist(1,1,L) ,im*jm) |
772 |
|
|
call back2grd (delqgather(1,L),pblindex, dqmoist(1,1,L,1),im*jm) |
773 |
|
|
call back2grd ( cpen(1,1,L),pblindex, cpen(1,1,L) ,im*jm) |
774 |
|
|
call back2grd ( cldls(1,1,L),pblindex, cldls(1,1,L) ,im*jm) |
775 |
|
|
call back2grd (cldwater(1,1,L),pblindex,cldwater(1,1,L) ,im*jm) |
776 |
|
|
call back2grd ( pkzgather(1,L),pblindex, pkzgather(1,L) ,im*jm) |
777 |
|
|
|
778 |
|
|
C Diagnostics: |
779 |
|
|
c ------------ |
780 |
|
|
if(icldmas.gt.0)call back2grd(tmpgather(1,L),pblindex, |
781 |
|
|
. tmpgather(1,L),im*jm) |
782 |
|
|
if(idtrain.gt.0)call back2grd(pkegather(1,L),pblindex, |
783 |
|
|
. pkegather(1,L),im*jm) |
784 |
|
|
if(idtls.gt.0)call back2grd(deltrnev(1,L),pblindex, |
785 |
|
|
. deltrnev(1,L),im*jm) |
786 |
|
|
if(idqls.gt.0)call back2grd(delqrnev(1,L),pblindex, |
787 |
|
|
. delqrnev(1,L),im*jm) |
788 |
|
|
if(icldnp.gt.0)call back2grd(cldsr(1,1,L),pblindex, |
789 |
|
|
. cldsr(1,1,L),im*jm) |
790 |
|
|
enddo |
791 |
|
|
|
792 |
molod |
1.21 |
c General Tracers |
793 |
|
|
c --------------- |
794 |
|
|
do nt = 1,ntracerin-ptracer |
795 |
|
|
do L = 1,lm |
796 |
|
|
call back2grd (delugather(1,L,nt),pblindex, |
797 |
|
|
. dqmoist(1,1,L,ptracer+nt),im*jm) |
798 |
|
|
enddo |
799 |
|
|
enddo |
800 |
|
|
|
801 |
|
|
if(cumfric) then |
802 |
|
|
|
803 |
|
|
c U and V for cumulus friction |
804 |
|
|
c ---------------------------- |
805 |
|
|
do L = 1,lm |
806 |
|
|
call back2grd (delugather(1,L,ntracerin-ptracer+1),pblindex, |
807 |
|
|
. dumoist(1,1,L),im*jm) |
808 |
|
|
call back2grd (delugather(1,L,ntracerin-ptracer+2),pblindex, |
809 |
|
|
. dvmoist(1,1,L),im*jm) |
810 |
|
|
enddo |
811 |
|
|
|
812 |
|
|
C Remove pi-weighting for u and v tendencies |
813 |
|
|
do j = 1,jm |
814 |
|
|
do i = 1,im |
815 |
|
|
tmpimjm(i,j) = 1./pz(i,j) |
816 |
|
|
enddo |
817 |
|
|
enddo |
818 |
|
|
do L = 1,lm |
819 |
|
|
do j = 1,jm |
820 |
|
|
do i = 1,im |
821 |
|
|
dumoist(i,j,L) = dumoist(i,j,L) * tmpimjm(i,j) |
822 |
molod |
1.22 |
dvmoist(i,j,L) = dvmoist(i,j,L) * tmpimjm(i,j) |
823 |
molod |
1.21 |
enddo |
824 |
|
|
enddo |
825 |
|
|
enddo |
826 |
molod |
1.1 |
|
827 |
molod |
1.22 |
|
828 |
molod |
1.21 |
endif |
829 |
molod |
1.1 |
|
830 |
|
|
C ********************************************************************** |
831 |
|
|
C BUMP DIAGNOSTICS |
832 |
|
|
C ********************************************************************** |
833 |
|
|
|
834 |
molod |
1.16 |
|
835 |
molod |
1.1 |
c Sub-Cloud Layer |
836 |
|
|
c ------------------------- |
837 |
|
|
if( ipsubcld.ne.0 ) then |
838 |
|
|
do j = 1,jm |
839 |
|
|
do i = 1,im |
840 |
molod |
1.3 |
qdiag(i,j,ipsubcld,bi,bj) = qdiag(i,j,ipsubcld,bi,bj) + |
841 |
|
|
. psubcldg (i,j) |
842 |
|
|
qdiag(i,j,ipsubcldc,bi,bj) = qdiag(i,j,ipsubcldc,bi,bj) + |
843 |
|
|
. psubcldgc(i,j) |
844 |
molod |
1.1 |
enddo |
845 |
|
|
enddo |
846 |
|
|
endif |
847 |
|
|
|
848 |
|
|
c Non-Precipitating Cloud Fraction |
849 |
|
|
c -------------------------------- |
850 |
|
|
if( icldnp.ne.0 ) then |
851 |
|
|
do L = 1,lm |
852 |
|
|
do j = 1,jm |
853 |
|
|
do i = 1,im |
854 |
molod |
1.3 |
qdiag(i,j,icldnp+L-1,bi,bj) = qdiag(i,j,icldnp+L-1,bi,bj) + |
855 |
|
|
. cldsr(i,j,L) |
856 |
molod |
1.1 |
enddo |
857 |
|
|
enddo |
858 |
|
|
enddo |
859 |
|
|
ncldnp = ncldnp + 1 |
860 |
|
|
endif |
861 |
|
|
|
862 |
|
|
c Moist Processes Heating Rate |
863 |
|
|
c ---------------------------- |
864 |
|
|
if(imoistt.gt.0) then |
865 |
|
|
do L = 1,lm |
866 |
molod |
1.16 |
do j = 1,jm |
867 |
|
|
do i = 1,im |
868 |
|
|
indgath = (j-1)*im + i |
869 |
|
|
qdiag(i,j,imoistt+L-1,bi,bj) = qdiag(i,j,imoistt+L-1,bi,bj) + |
870 |
|
|
. (dtmoist(i,j,L)*sday*pkzgather(indgath,L)/pz(i,j)) |
871 |
|
|
enddo |
872 |
molod |
1.1 |
enddo |
873 |
|
|
enddo |
874 |
|
|
endif |
875 |
|
|
|
876 |
|
|
c Moist Processes Moistening Rate |
877 |
|
|
c ------------------------------- |
878 |
|
|
if(imoistq.gt.0) then |
879 |
|
|
do L = 1,lm |
880 |
|
|
do j = 1,jm |
881 |
|
|
do i = 1,im |
882 |
molod |
1.3 |
qdiag(i,j,imoistq+L-1,bi,bj) = qdiag(i,j,imoistq+L-1,bi,bj) + |
883 |
molod |
1.1 |
. (dqmoist(i,j,L,1)*sday*1000.0/pz(i,j)) |
884 |
|
|
enddo |
885 |
|
|
enddo |
886 |
|
|
enddo |
887 |
|
|
endif |
888 |
|
|
|
889 |
|
|
c Cloud Mass Flux |
890 |
|
|
c --------------- |
891 |
|
|
if(icldmas.gt.0) then |
892 |
|
|
do L = 1,lm |
893 |
molod |
1.16 |
do j = 1,jm |
894 |
|
|
do i = 1,im |
895 |
|
|
indgath = (j-1)*im + i |
896 |
|
|
qdiag(i,j,icldmas+L-1,bi,bj) = qdiag(i,j,icldmas+L-1,bi,bj) + |
897 |
|
|
. tmpgather(indgath,L) |
898 |
|
|
enddo |
899 |
molod |
1.1 |
enddo |
900 |
|
|
enddo |
901 |
|
|
endif |
902 |
|
|
|
903 |
|
|
c Detrained Cloud Mass Flux |
904 |
|
|
c ------------------------- |
905 |
|
|
if(idtrain.gt.0) then |
906 |
|
|
do L = 1,lm |
907 |
molod |
1.16 |
do j = 1,jm |
908 |
|
|
do i = 1,im |
909 |
|
|
indgath = (j-1)*im + i |
910 |
|
|
qdiag(i,j,idtrain+L-1,bi,bj) = qdiag(i,j,idtrain+L-1,bi,bj) + |
911 |
|
|
. pkegather(indgath,L) |
912 |
|
|
enddo |
913 |
molod |
1.1 |
enddo |
914 |
|
|
enddo |
915 |
|
|
endif |
916 |
|
|
|
917 |
|
|
c Grid-Scale Condensational Heating Rate |
918 |
|
|
c -------------------------------------- |
919 |
|
|
if(idtls.gt.0) then |
920 |
|
|
do L = 1,lm |
921 |
molod |
1.16 |
do j = 1,jm |
922 |
|
|
do i = 1,im |
923 |
|
|
indgath = (j-1)*im + i |
924 |
|
|
qdiag(i,j,idtls+L-1,bi,bj) = qdiag(i,j,idtls+L-1,bi,bj) + |
925 |
|
|
. deltrnev(indgath,L) |
926 |
|
|
enddo |
927 |
molod |
1.1 |
enddo |
928 |
|
|
enddo |
929 |
|
|
endif |
930 |
|
|
|
931 |
|
|
c Grid-Scale Condensational Moistening Rate |
932 |
|
|
c ----------------------------------------- |
933 |
|
|
if(idqls.gt.0) then |
934 |
|
|
do L = 1,lm |
935 |
molod |
1.16 |
do j = 1,jm |
936 |
|
|
do i = 1,im |
937 |
|
|
indgath = (j-1)*im + i |
938 |
|
|
qdiag(i,j,idqls+L-1,bi,bj) = qdiag(i,j,idqls+L-1,bi,bj) + |
939 |
|
|
. delqrnev(indgath,L) |
940 |
|
|
enddo |
941 |
molod |
1.1 |
enddo |
942 |
|
|
enddo |
943 |
|
|
endif |
944 |
|
|
|
945 |
|
|
c Total Precipitation |
946 |
|
|
c ------------------- |
947 |
|
|
if(ipreacc.gt.0) then |
948 |
|
|
do j = 1,jm |
949 |
|
|
do i = 1,im |
950 |
molod |
1.3 |
qdiag(i,j,ipreacc,bi,bj) = qdiag(i,j,ipreacc,bi,bj) |
951 |
molod |
1.1 |
. + ( lsp_new(I,j) |
952 |
|
|
. + snow_new(I,j) |
953 |
|
|
. + conv_new(i,j) ) *sday*tminv |
954 |
|
|
enddo |
955 |
|
|
enddo |
956 |
|
|
endif |
957 |
|
|
|
958 |
|
|
c Convective Precipitation |
959 |
|
|
c ------------------------ |
960 |
|
|
if(iprecon.gt.0) then |
961 |
molod |
1.16 |
do j = 1,jm |
962 |
|
|
do i = 1,im |
963 |
|
|
indgath = (j-1)*im + i |
964 |
|
|
qdiag(i,j,iprecon,bi,bj) = qdiag(i,j,iprecon,bi,bj) + |
965 |
|
|
. raincgath(indgath)*sday*tminv |
966 |
|
|
enddo |
967 |
molod |
1.1 |
enddo |
968 |
|
|
endif |
969 |
|
|
|
970 |
|
|
C ********************************************************************** |
971 |
|
|
C **** Fill Rainfall and Snowfall Arrays for Land Surface Model **** |
972 |
|
|
C **** Note: Precip Rates work when DT(turb)<DT(moist) **** |
973 |
|
|
C ********************************************************************** |
974 |
|
|
|
975 |
|
|
do j = 1,jm |
976 |
|
|
do i = 1,im |
977 |
|
|
rainlsp (i,j) = rainlsp (i,j) + lsp_new(i,j)*tminv |
978 |
|
|
rainconv(i,j) = rainconv(i,j) + conv_new(i,j)*tminv |
979 |
|
|
snowfall(i,j) = snowfall(i,j) + snow_new(i,j)*tminv |
980 |
|
|
enddo |
981 |
|
|
enddo |
982 |
|
|
|
983 |
|
|
C ********************************************************************** |
984 |
|
|
C *** Compute Time-averaged Quantities for Radiation *** |
985 |
|
|
C *** CPEN => Cloud Fraction from RAS *** |
986 |
|
|
C *** CLDLS => Cloud Fraction from RNEVP *** |
987 |
|
|
C ********************************************************************** |
988 |
|
|
|
989 |
|
|
do j = 1,jm |
990 |
|
|
do i = 1,im |
991 |
|
|
cldhi (i,j) = 0. |
992 |
|
|
cldmid(i,j) = 0. |
993 |
|
|
cldlow(i,j) = 0. |
994 |
|
|
cldtmp(i,j) = 0. |
995 |
|
|
cldprs(i,j) = 0. |
996 |
|
|
tmpimjm(i,j) = 0. |
997 |
|
|
enddo |
998 |
|
|
enddo |
999 |
|
|
|
1000 |
|
|
c Set Moist-Process Memory Coefficient |
1001 |
|
|
c ------------------------------------ |
1002 |
|
|
cldras_mem = 1.0-tmstp/ 3600.0 |
1003 |
|
|
cldlsp_mem = 1.0-tmstp/(3600.0*3) |
1004 |
|
|
|
1005 |
|
|
do L = 1,lm |
1006 |
|
|
do i = 1,im*jm |
1007 |
molod |
1.6 |
plev = pl(i,L) |
1008 |
molod |
1.1 |
|
1009 |
|
|
c Compute Time-averaged Cloud and Water Amounts for Longwave Radiation |
1010 |
|
|
c -------------------------------------------------------------------- |
1011 |
|
|
watnow = cldwater(i,1,L) |
1012 |
|
|
if( plev.le.500.0 ) then |
1013 |
|
|
cldras = min( max( cldras_lw(i,1,L)*cldras_mem,cpen(i,1,L)),1.0) |
1014 |
|
|
else |
1015 |
|
|
cldras = 0.0 |
1016 |
|
|
endif |
1017 |
|
|
cldlsp = min( max( cldlsp_lw(i,1,L)*cldlsp_mem,cldls(i,1,L)),1.0) |
1018 |
|
|
|
1019 |
|
|
if( cldras.lt.cldmin ) cldras = 0.0 |
1020 |
|
|
if( cldlsp.lt.cldmin ) cldlsp = 0.0 |
1021 |
|
|
|
1022 |
|
|
cldnow = max( cldlsp,cldras ) |
1023 |
|
|
|
1024 |
|
|
lwlz(i,1,L) = ( nlwlz*lwlz(i,1,L) + watnow)/(nlwlz +1) |
1025 |
|
|
cldtot_lw(i,1,L) = (nlwcld*cldtot_lw(i,1,L) + cldnow)/(nlwcld+1) |
1026 |
|
|
cldlsp_lw(i,1,L) = (nlwcld*cldlsp_lw(i,1,L) + cldlsp)/(nlwcld+1) |
1027 |
|
|
cldras_lw(i,1,L) = (nlwcld*cldras_lw(i,1,L) + cldras)/(nlwcld+1) |
1028 |
|
|
|
1029 |
|
|
|
1030 |
|
|
c Compute Time-averaged Cloud and Water Amounts for Shortwave Radiation |
1031 |
|
|
c --------------------------------------------------------------------- |
1032 |
|
|
watnow = cldwater(i,1,L) |
1033 |
|
|
if( plev.le.500.0 ) then |
1034 |
|
|
cldras = min( max(cldras_sw(i,1,L)*cldras_mem, cpen(i,1,L)),1.0) |
1035 |
|
|
else |
1036 |
|
|
cldras = 0.0 |
1037 |
|
|
endif |
1038 |
|
|
cldlsp = min( max(cldlsp_sw(i,1,L)*cldlsp_mem,cldls(i,1,L)),1.0) |
1039 |
|
|
|
1040 |
|
|
if( cldras.lt.cldmin ) cldras = 0.0 |
1041 |
|
|
if( cldlsp.lt.cldmin ) cldlsp = 0.0 |
1042 |
|
|
|
1043 |
|
|
cldnow = max( cldlsp,cldras ) |
1044 |
|
|
|
1045 |
|
|
swlz(i,1,L) = ( nswlz*swlz(i,1,L) + watnow)/(nswlz +1) |
1046 |
|
|
cldtot_sw(i,1,L) = (nswcld*cldtot_sw(i,1,L) + cldnow)/(nswcld+1) |
1047 |
|
|
cldlsp_sw(i,1,L) = (nswcld*cldlsp_sw(i,1,L) + cldlsp)/(nswcld+1) |
1048 |
|
|
cldras_sw(i,1,L) = (nswcld*cldras_sw(i,1,L) + cldras)/(nswcld+1) |
1049 |
|
|
|
1050 |
|
|
|
1051 |
|
|
c Compute Instantaneous Low-Mid-High Maximum Overlap Cloud Fractions |
1052 |
|
|
c ---------------------------------------------------------------------- |
1053 |
|
|
|
1054 |
|
|
if( L.lt.midlevel ) cldhi (i,1) = max( cldnow,cldhi (i,1) ) |
1055 |
|
|
if( L.ge.midlevel .and. |
1056 |
|
|
. L.lt.lowlevel ) cldmid(i,1) = max( cldnow,cldmid(i,1) ) |
1057 |
|
|
if( L.ge.lowlevel ) cldlow(i,1) = max( cldnow,cldlow(i,1) ) |
1058 |
|
|
|
1059 |
|
|
c Compute Cloud-Top Temperature and Pressure |
1060 |
|
|
c ------------------------------------------ |
1061 |
|
|
cldtmp(i,1) = cldtmp(i,1) + cldnow*pkzgather(i,L) |
1062 |
|
|
. * ( tz(i,1,L) + dtmoist(i,1,L)*tmstp/pz(i,1) ) |
1063 |
|
|
cldprs(i,1) = cldprs(i,1) + cldnow*plev |
1064 |
|
|
tmpimjm(i,1) = tmpimjm(i,1) + cldnow |
1065 |
|
|
|
1066 |
|
|
enddo |
1067 |
|
|
enddo |
1068 |
|
|
|
1069 |
molod |
1.21 |
c Compute Instantaneous Total 2-D Cloud Fraction |
1070 |
molod |
1.1 |
c ---------------------------------------------- |
1071 |
|
|
do j = 1,jm |
1072 |
|
|
do i = 1,im |
1073 |
|
|
totcld(i,j) = 1.0 - (1.-cldhi (i,j)) |
1074 |
|
|
. * (1.-cldmid(i,j)) |
1075 |
|
|
. * (1.-cldlow(i,j)) |
1076 |
|
|
enddo |
1077 |
|
|
enddo |
1078 |
|
|
|
1079 |
|
|
|
1080 |
|
|
C ********************************************************************** |
1081 |
|
|
C *** Fill Cloud Top Pressure and Temperature Diagnostic *** |
1082 |
|
|
C ********************************************************************** |
1083 |
|
|
|
1084 |
|
|
if(icldtmp.gt.0) then |
1085 |
|
|
do j = 1,jm |
1086 |
|
|
do i = 1,im |
1087 |
|
|
if( cldtmp(i,j).gt.0.0 ) then |
1088 |
molod |
1.3 |
qdiag(i,j,icldtmp,bi,bj) = qdiag(i,j,icldtmp,bi,bj) + |
1089 |
molod |
1.1 |
. cldtmp(i,j)*totcld(i,j)/tmpimjm(i,j) |
1090 |
molod |
1.3 |
qdiag(i,j,icttcnt,bi,bj) = qdiag(i,j,icttcnt,bi,bj) + |
1091 |
|
|
. totcld(i,j) |
1092 |
molod |
1.1 |
endif |
1093 |
|
|
enddo |
1094 |
|
|
enddo |
1095 |
|
|
endif |
1096 |
|
|
|
1097 |
|
|
if(icldprs.gt.0) then |
1098 |
|
|
do j = 1,jm |
1099 |
|
|
do i = 1,im |
1100 |
|
|
if( cldprs(i,j).gt.0.0 ) then |
1101 |
molod |
1.3 |
qdiag(i,j,icldprs,bi,bj) = qdiag(i,j,icldprs,bi,bj) + |
1102 |
molod |
1.1 |
. cldprs(i,j)*totcld(i,j)/tmpimjm(i,j) |
1103 |
molod |
1.3 |
qdiag(i,j,ictpcnt,bi,bj) = qdiag(i,j,ictpcnt,bi,bj) + |
1104 |
|
|
. totcld(i,j) |
1105 |
molod |
1.1 |
endif |
1106 |
|
|
enddo |
1107 |
|
|
enddo |
1108 |
|
|
endif |
1109 |
|
|
|
1110 |
|
|
C ********************************************************************** |
1111 |
|
|
C **** INCREMENT COUNTERS **** |
1112 |
|
|
C ********************************************************************** |
1113 |
|
|
|
1114 |
|
|
nlwlz = nlwlz + 1 |
1115 |
|
|
nswlz = nswlz + 1 |
1116 |
|
|
|
1117 |
|
|
nlwcld = nlwcld + 1 |
1118 |
|
|
nswcld = nswcld + 1 |
1119 |
|
|
|
1120 |
molod |
1.16 |
#ifdef ALLOW_DIAGNOSTICS |
1121 |
|
|
if( (bi.eq.1) .and. (bj.eq.1) ) then |
1122 |
molod |
1.1 |
nmoistt = nmoistt + 1 |
1123 |
|
|
nmoistq = nmoistq + 1 |
1124 |
|
|
npreacc = npreacc + 1 |
1125 |
|
|
nprecon = nprecon + 1 |
1126 |
|
|
|
1127 |
|
|
ncldmas = ncldmas + 1 |
1128 |
|
|
ndtrain = ndtrain + 1 |
1129 |
|
|
|
1130 |
|
|
ndtls = ndtls + 1 |
1131 |
|
|
ndqls = ndqls + 1 |
1132 |
molod |
1.22 |
|
1133 |
molod |
1.16 |
endif |
1134 |
|
|
#endif |
1135 |
molod |
1.1 |
|
1136 |
|
|
RETURN |
1137 |
|
|
END |
1138 |
molod |
1.16 |
SUBROUTINE RAS( NN, LNG, LENC, K, NLTOP, nlayr, DT |
1139 |
molod |
1.21 |
*, UOI, ntracedim, ntracer, POI, QOI, PRS, PRJ, rnd, ncrnd |
1140 |
|
|
*, RAINS, CLN, CLF, cldmas, detrain |
1141 |
|
|
*, cp,grav,rkappa,alhl,rhfrac,rasmax ) |
1142 |
molod |
1.1 |
C |
1143 |
|
|
C********************************************************************* |
1144 |
|
|
C********************* SUBROUTINE RAS ***************************** |
1145 |
|
|
C********************** 16 MARCH 1988 ****************************** |
1146 |
|
|
C********************************************************************* |
1147 |
|
|
C |
1148 |
molod |
1.8 |
implicit none |
1149 |
|
|
|
1150 |
molod |
1.9 |
C Argument List |
1151 |
molod |
1.16 |
integer nn,lng,lenc,k,nltop,nlayr |
1152 |
molod |
1.21 |
integer ntracedim, ntracer |
1153 |
molod |
1.9 |
integer ncrnd |
1154 |
molod |
1.13 |
_RL dt |
1155 |
molod |
1.21 |
_RL UOI(lng,nlayr,ntracedim), POI(lng,K) |
1156 |
molod |
1.16 |
_RL QOI(lng,K), PRS(lng,K+1), PRJ(lng,K+1) |
1157 |
molod |
1.13 |
_RL rnd(ncrnd) |
1158 |
molod |
1.16 |
_RL RAINS(lng,K), CLN(lng,K), CLF(lng,K) |
1159 |
|
|
_RL cldmas(lng,K), detrain(lng,K) |
1160 |
|
|
_RL cp,grav,rkappa,alhl,rhfrac(lng),rasmax |
1161 |
molod |
1.9 |
|
1162 |
|
|
C Local Variables |
1163 |
molod |
1.16 |
_RL TCU(lng,K), QCU(lng,K) |
1164 |
molod |
1.21 |
_RL ucu(lng,K,ntracedim) |
1165 |
molod |
1.16 |
_RL ALF(lng,K), BET(lng,K), GAM(lng,K) |
1166 |
|
|
*, ETA(lng,K), HOI(lng,K) |
1167 |
|
|
*, PRH(lng,K), PRI(lng,K) |
1168 |
|
|
_RL HST(lng,K), QOL(lng,K), GMH(lng,K) |
1169 |
|
|
|
1170 |
|
|
_RL TX1(lng), TX2(lng), TX3(lng), TX4(lng), TX5(lng) |
1171 |
|
|
*, TX6(lng), TX7(lng), TX8(lng), TX9(lng) |
1172 |
molod |
1.21 |
*, TX11(lng), TX12(lng), TX13(lng), TX14(lng,ntracedim) |
1173 |
molod |
1.16 |
*, TX15(lng) |
1174 |
|
|
*, WFN(lng) |
1175 |
|
|
integer IA1(lng), IA2(lng), IA3(lng) |
1176 |
|
|
_RL cloudn(lng), pcu(lng) |
1177 |
molod |
1.1 |
|
1178 |
molod |
1.8 |
integer krmin,icm |
1179 |
molod |
1.13 |
_RL rknob, cmb2pa |
1180 |
molod |
1.8 |
PARAMETER (KRMIN=01) |
1181 |
|
|
PARAMETER (ICM=1000) |
1182 |
|
|
PARAMETER (CMB2PA=100.0) |
1183 |
|
|
PARAMETER (rknob = 10.) |
1184 |
molod |
1.9 |
|
1185 |
|
|
integer IC(ICM), IRND(icm) |
1186 |
molod |
1.16 |
_RL cmass(lng,K) |
1187 |
molod |
1.9 |
LOGICAL SETRAS |
1188 |
molod |
1.24 |
integer ifound |
1189 |
|
|
_RL temp |
1190 |
|
|
_RL thbef(lng,K) |
1191 |
molod |
1.9 |
|
1192 |
|
|
integer i,L,nc,ib,nt |
1193 |
molod |
1.8 |
integer km1,kp1,kprv,kcr,kfx,ncmx |
1194 |
molod |
1.13 |
_RL p00, crtmsf, frac, rasblf |
1195 |
molod |
1.8 |
|
1196 |
|
|
do L = 1,k |
1197 |
|
|
do I = 1,LENC |
1198 |
|
|
rains(i,l) = 0. |
1199 |
|
|
enddo |
1200 |
|
|
enddo |
1201 |
molod |
1.1 |
|
1202 |
|
|
p00 = 1000. |
1203 |
|
|
crtmsf = 0. |
1204 |
|
|
|
1205 |
|
|
C The numerator here is the fraction of the subcloud layer mass flux |
1206 |
|
|
C allowed to entrain into the cloud |
1207 |
|
|
|
1208 |
|
|
CCC FRAC = 1./dt |
1209 |
molod |
1.24 |
CCC FRAC = 0.5/dt |
1210 |
molod |
1.1 |
FRAC = 0.5/dt |
1211 |
|
|
|
1212 |
|
|
KM1 = K - 1 |
1213 |
|
|
KP1 = K + 1 |
1214 |
|
|
C we want the ras adjustment time scale to be one hour (indep of dt) |
1215 |
|
|
RASBLF = 1./3600. |
1216 |
|
|
C |
1217 |
|
|
KPRV = KM1 |
1218 |
|
|
C Removed KRMAX parameter |
1219 |
|
|
KCR = MIN(KM1,nlayr-2) |
1220 |
molod |
1.24 |
KFX = KM1 - KCR |
1221 |
molod |
1.1 |
NCMX = KFX + NCRND |
1222 |
|
|
C |
1223 |
|
|
IF (KFX .GT. 0) THEN |
1224 |
|
|
DO NC=1,KFX |
1225 |
|
|
IC(NC) = K - NC |
1226 |
|
|
ENDDO |
1227 |
|
|
ENDIF |
1228 |
|
|
C |
1229 |
|
|
IF (NCRND .GT. 0) THEN |
1230 |
|
|
DO I=1,ncrnd |
1231 |
|
|
IRND(I) = (RND(I)-0.0005)*(KCR-KRMIN+1) |
1232 |
|
|
IRND(I) = IRND(I) + KRMIN |
1233 |
|
|
ENDDO |
1234 |
|
|
C |
1235 |
|
|
DO NC=1,NCRND |
1236 |
|
|
IC(KFX+NC) = IRND(NC) |
1237 |
|
|
ENDDO |
1238 |
|
|
ENDIF |
1239 |
|
|
C |
1240 |
|
|
DO 100 NC=1,NCMX |
1241 |
|
|
C |
1242 |
|
|
IF (NC .EQ. 1 ) THEN |
1243 |
|
|
SETRAS = .TRUE. |
1244 |
|
|
ELSE |
1245 |
|
|
SETRAS = .FALSE. |
1246 |
|
|
ENDIF |
1247 |
|
|
IB = IC(NC) |
1248 |
molod |
1.16 |
|
1249 |
molod |
1.1 |
c Initialize Cloud Fraction Array |
1250 |
|
|
c ------------------------------- |
1251 |
|
|
do i = 1,lenc |
1252 |
|
|
cloudn(i) = 0.0 |
1253 |
|
|
enddo |
1254 |
|
|
|
1255 |
molod |
1.16 |
CALL CLOUD(nn,lng, LENC, K, NLTOP, nlayr, IB, RASBLF,SETRAS,FRAC |
1256 |
molod |
1.1 |
*, CP, ALHL, RKAPPA, GRAV, P00, CRTMSF |
1257 |
molod |
1.21 |
*, POI, QOI, UOI, ntracedim, Ntracer, PRS, PRJ |
1258 |
molod |
1.1 |
*, PCU, CLOUDN, TCU, QCU, UCU, CMASS |
1259 |
|
|
*, ALF, BET, GAM, PRH, PRI, HOI, ETA |
1260 |
|
|
*, HST, QOL, GMH |
1261 |
|
|
*, TX1, TX2, TX3, TX4, TX5, TX6, TX7, TX8, TX9 |
1262 |
|
|
*, WFN, TX11, TX12, TX13, TX14, TX15 |
1263 |
|
|
*, IA1,IA2,IA3,rhfrac) |
1264 |
|
|
|
1265 |
|
|
C Compute fraction of grid box into which rain re-evap occurs (clf) |
1266 |
|
|
c ----------------------------------------------------------------- |
1267 |
|
|
do i = 1,lenc |
1268 |
|
|
|
1269 |
|
|
c mass in detrainment layer |
1270 |
|
|
c ------------------------- |
1271 |
|
|
tx1(i) = cmb2pa * (prs(i,ib+1) - prs(i,ib))/(grav*dt) |
1272 |
|
|
|
1273 |
|
|
c ratio of detraining cloud mass to mass in detrainment layer |
1274 |
|
|
c ----------------------------------------------------------- |
1275 |
|
|
tx1(i) = rhfrac(i)*rknob * cmass(i,ib) / tx1(i) |
1276 |
|
|
if(cmass(i,K).gt.0.) clf(i,ib) = clf(i,ib) + tx1(i) |
1277 |
|
|
if( clf(i,ib).gt.1.) clf(i,ib) = 1. |
1278 |
|
|
enddo |
1279 |
|
|
|
1280 |
|
|
c Compute Total Cloud Mass Flux |
1281 |
|
|
c ***************************** |
1282 |
|
|
do L=ib,k |
1283 |
|
|
do i=1,lenc |
1284 |
|
|
cmass(i,L) = rhfrac(i)*cmass(i,L) * dt |
1285 |
|
|
enddo |
1286 |
|
|
enddo |
1287 |
|
|
|
1288 |
|
|
do L=ib,k |
1289 |
|
|
do i=1,lenc |
1290 |
|
|
cldmas(i,L) = cldmas(i,L) + cmass(i,L) |
1291 |
|
|
enddo |
1292 |
|
|
enddo |
1293 |
|
|
|
1294 |
|
|
do i=1,lenc |
1295 |
|
|
detrain(i,ib) = detrain(i,ib) + cmass(i,ib) |
1296 |
|
|
enddo |
1297 |
|
|
|
1298 |
|
|
DO L=IB,K |
1299 |
|
|
DO I=1,LENC |
1300 |
molod |
1.24 |
thbef(I,L) = POI(I,L) |
1301 |
molod |
1.1 |
POI(I,L) = POI(I,L) + TCU(I,L) * DT * rhfrac(i) |
1302 |
|
|
QOI(I,L) = QOI(I,L) + QCU(I,L) * DT * rhfrac(i) |
1303 |
|
|
ENDDO |
1304 |
|
|
ENDDO |
1305 |
molod |
1.21 |
DO NT=1,Ntracer |
1306 |
|
|
DO L=IB,K |
1307 |
|
|
DO I=1,LENC |
1308 |
|
|
UOI(I,L+nltop-1,NT)=UOI(I,L+nltop-1,NT)+UCU(I,L,NT)*DT*rhfrac(i) |
1309 |
|
|
ENDDO |
1310 |
|
|
ENDDO |
1311 |
|
|
ENDDO |
1312 |
molod |
1.1 |
DO I=1,LENC |
1313 |
|
|
rains(I,ib) = rains(I,ib) + PCU(I)*dt * rhfrac(i) |
1314 |
|
|
ENDDO |
1315 |
|
|
|
1316 |
molod |
1.24 |
do i = 1,lenc |
1317 |
|
|
ifound = 0 |
1318 |
|
|
do L = 1,k |
1319 |
|
|
if(tcu(i,L).ne.0.)ifound = ifound + 1 |
1320 |
|
|
enddo |
1321 |
|
|
if(ifound.ne.0) then |
1322 |
|
|
c print *,i,' made a cloud detraining at ',ib |
1323 |
|
|
do L = 1,k |
1324 |
|
|
temp = TCU(I,L) * DT * rhfrac(i) |
1325 |
|
|
c write(6,122)L,thbef(i,L),poi(i,L),temp |
1326 |
|
|
enddo |
1327 |
|
|
endif |
1328 |
|
|
enddo |
1329 |
|
|
|
1330 |
molod |
1.1 |
100 CONTINUE |
1331 |
molod |
1.24 |
|
1332 |
|
|
122 format(' ',i3,' TH B ',e10.3,' TH A ',e10.3,' DTH ',e10.3) |
1333 |
molod |
1.1 |
|
1334 |
|
|
c Fill Convective Cloud Fractions based on 3-D Rain Amounts |
1335 |
|
|
c --------------------------------------------------------- |
1336 |
|
|
do L=k-1,1,-1 |
1337 |
|
|
do i=1,lenc |
1338 |
|
|
tx1(i) = 100*(prs(i,L+1)-prs(i,L))/grav |
1339 |
|
|
cln(i,L) = min(1600*rains(i,L)/tx1(i),rasmax ) |
1340 |
|
|
enddo |
1341 |
|
|
enddo |
1342 |
|
|
|
1343 |
|
|
RETURN |
1344 |
|
|
END |
1345 |
|
|
subroutine rndcloud (iras,nrnd,rnd,myid) |
1346 |
|
|
implicit none |
1347 |
|
|
integer n,iras,nrnd,myid |
1348 |
molod |
1.13 |
_RL random_numbx |
1349 |
jmc |
1.29 |
c _RL rnd(nrnd) |
1350 |
|
|
_RL rnd(*) |
1351 |
molod |
1.1 |
integer irm |
1352 |
|
|
parameter (irm = 1000) |
1353 |
molod |
1.13 |
_RL random(irm) |
1354 |
jmc |
1.29 |
integer i,mcheck,iseed,indx |
1355 |
molod |
1.1 |
logical first |
1356 |
|
|
data first /.true./ |
1357 |
|
|
integer iras0 |
1358 |
|
|
data iras0 /0/ |
1359 |
|
|
save random, iras0 |
1360 |
|
|
|
1361 |
jmc |
1.29 |
if(nrnd.eq.0)then |
1362 |
molod |
1.1 |
do i = 1,nrnd |
1363 |
|
|
rnd(i) = 0 |
1364 |
|
|
enddo |
1365 |
molod |
1.15 |
if(first .and. myid.eq.1) print *,' NO RANDOM CLOUDS IN RAS ' |
1366 |
molod |
1.1 |
go to 100 |
1367 |
|
|
endif |
1368 |
|
|
|
1369 |
|
|
mcheck = mod(iras-1,irm/nrnd) |
1370 |
|
|
|
1371 |
jmc |
1.29 |
c print *,' RNDCLOUD: first ',first,' iras ',iras,' iras0 ',iras0 |
1372 |
|
|
c print *,' RNDCLOUD: irm,nrnd,mcheck=',irm,nrnd,mcheck |
1373 |
|
|
|
1374 |
|
|
if ( iras.eq.iras0 ) then |
1375 |
|
|
C- Not the 1rst tile: we are all set (already done for the 1rst tile): |
1376 |
|
|
c ----------------------------------------------------------------------- |
1377 |
|
|
indx = (iras-1)*nrnd |
1378 |
|
|
|
1379 |
molod |
1.1 |
c First Time In From a Continuing RESTART (IRAS.GT.1) or Reading a New RESTART |
1380 |
jmc |
1.29 |
c -- or -- |
1381 |
|
|
c Multiple Time In But have Used Up all 1000 numbers (MCHECK.EQ.0) |
1382 |
molod |
1.1 |
c ---------------------------------------------------------------------------- |
1383 |
jmc |
1.29 |
elseif ( first.and.(iras.gt.1) .or. mcheck.eq.0 ) then |
1384 |
|
|
iseed = (iras-1-mcheck)*nrnd |
1385 |
molod |
1.1 |
call random_seedx(iseed) |
1386 |
|
|
do i = 1,irm |
1387 |
molod |
1.12 |
random(i) = random_numbx(iseed) |
1388 |
molod |
1.1 |
enddo |
1389 |
molod |
1.16 |
indx = (iras-1)*nrnd |
1390 |
molod |
1.1 |
|
1391 |
jmc |
1.29 |
if( myid.eq.1 ) print *, 'Creating Rand Numb Array in RNDCLOUD' |
1392 |
|
|
& ,', iseed=', iseed |
1393 |
|
|
if( myid.eq.1 ) print *, 'IRAS: ',iras,' IRAS0: ',iras0, |
1394 |
|
|
& ' indx: ', mod(indx,irm) |
1395 |
molod |
1.1 |
|
1396 |
|
|
c Multiple Time In But have NOT Used Up all 1000 numbers (MCHECK.NE.0) |
1397 |
|
|
c -------------------------------------------------------------------- |
1398 |
|
|
else |
1399 |
molod |
1.16 |
indx = (iras-1)*nrnd |
1400 |
molod |
1.1 |
endif |
1401 |
|
|
|
1402 |
molod |
1.16 |
indx = mod(indx,irm) |
1403 |
jmc |
1.29 |
if( indx+nrnd.gt.irm ) then |
1404 |
|
|
c if( myid.eq.1 .AND. iras.ne.iras0 ) print *, |
1405 |
|
|
c & 'reach end of Rand Numb Array in RNDCLOUD',indx,irm-nrnd |
1406 |
|
|
indx=irm-nrnd |
1407 |
|
|
endif |
1408 |
|
|
|
1409 |
molod |
1.1 |
do n = 1,nrnd |
1410 |
molod |
1.16 |
rnd(n) = random(indx+n) |
1411 |
molod |
1.1 |
enddo |
1412 |
|
|
|
1413 |
|
|
100 continue |
1414 |
|
|
first = .false. |
1415 |
|
|
iras0 = iras |
1416 |
jmc |
1.29 |
|
1417 |
molod |
1.1 |
return |
1418 |
|
|
end |
1419 |
molod |
1.12 |
function random_numbx(iseed) |
1420 |
molod |
1.1 |
implicit none |
1421 |
molod |
1.12 |
integer iseed |
1422 |
|
|
real *8 seed,port_rand |
1423 |
molod |
1.13 |
_RL random_numbx |
1424 |
molod |
1.10 |
#ifdef CRAY |
1425 |
molod |
1.13 |
_RL ranf |
1426 |
molod |
1.1 |
random_numbx = ranf() |
1427 |
molod |
1.12 |
#else |
1428 |
molod |
1.10 |
#ifdef SGI |
1429 |
molod |
1.13 |
_RL rand |
1430 |
molod |
1.1 |
random_numbx = rand() |
1431 |
jmc |
1.26 |
#else |
1432 |
jmc |
1.27 |
seed = -1.d0 |
1433 |
jmc |
1.26 |
random_numbx = port_rand(seed) |
1434 |
molod |
1.1 |
#endif |
1435 |
molod |
1.12 |
#endif |
1436 |
molod |
1.1 |
return |
1437 |
|
|
end |
1438 |
|
|
subroutine random_seedx (iseed) |
1439 |
|
|
implicit none |
1440 |
|
|
integer iseed |
1441 |
molod |
1.28 |
real *8 port_rand |
1442 |
molod |
1.10 |
#ifdef CRAY |
1443 |
molod |
1.1 |
call ranset (iseed) |
1444 |
jmc |
1.27 |
#else |
1445 |
molod |
1.10 |
#ifdef SGI |
1446 |
molod |
1.1 |
integer*4 seed |
1447 |
|
|
seed = iseed |
1448 |
|
|
call srand (seed) |
1449 |
jmc |
1.27 |
#else |
1450 |
|
|
real*8 tmpRdN |
1451 |
|
|
real*8 seed |
1452 |
|
|
seed = iseed |
1453 |
|
|
tmpRdN = port_rand(seed) |
1454 |
|
|
#endif |
1455 |
molod |
1.1 |
#endif |
1456 |
|
|
return |
1457 |
|
|
end |
1458 |
molod |
1.16 |
SUBROUTINE CLOUD(nn,lng, LENC, K, NLTOP, nlayr, IC, RASALF |
1459 |
molod |
1.1 |
*, SETRAS, FRAC |
1460 |
|
|
*, CP, ALHL, RKAP, GRAV, P00, CRTMSF |
1461 |
molod |
1.21 |
*, POI, QOI, UOI, ntracedim, Ntracer, PRS, PRJ |
1462 |
molod |
1.1 |
*, PCU, CLN, TCU, QCU, UCU, CMASS |
1463 |
|
|
*, ALF, BET, GAM, PRH, PRI, HOL, ETA |
1464 |
|
|
*, HST, QOL, GMH |
1465 |
|
|
*, TX1, TX2, TX3, TX4, TX5, TX6, TX7, TX8, ALM |
1466 |
|
|
*, WFN, AKM, QS1, CLF, UHT, WLQ |
1467 |
|
|
*, IA, I1, I2,rhfrac) |
1468 |
|
|
C |
1469 |
|
|
C********************************************************************* |
1470 |
|
|
C******************** Relaxed Arakawa-Schubert *********************** |
1471 |
|
|
C********************* Plug Compatible Version ********************** |
1472 |
|
|
C************************ SUBROUTINE CLOUD *************************** |
1473 |
|
|
C************************* 23 JULY 1992 *************************** |
1474 |
|
|
C********************************************************************* |
1475 |
|
|
C********************************************************************* |
1476 |
|
|
C********************************************************************* |
1477 |
|
|
C************************** Developed By ***************************** |
1478 |
|
|
C************************** ***************************** |
1479 |
|
|
C************************ Shrinivas Moorthi ************************** |
1480 |
|
|
C************************ and ************************** |
1481 |
|
|
C************************ Max J. Suarez ***************************** |
1482 |
|
|
C************************ ***************************** |
1483 |
|
|
C******************** Laboratory for Atmospheres ********************* |
1484 |
|
|
C****************** NASA/GSFC, Greenbelt, MD 20771 ******************* |
1485 |
|
|
C********************************************************************* |
1486 |
|
|
C********************************************************************* |
1487 |
|
|
C |
1488 |
|
|
C The calculations of Moorthi and Suarez (1992, MWR) are |
1489 |
|
|
C contained in the CLOUD routine. |
1490 |
|
|
C It is probably advisable, at least initially, to treat CLOUD |
1491 |
|
|
C as a black box that computes the single cloud adjustments. RAS, |
1492 |
|
|
C on the other hand, can be tailored to each GCMs configuration |
1493 |
|
|
C (ie, number and placement of levels, nature of boundary layer, |
1494 |
|
|
C time step and frequency with which RAS is called). |
1495 |
|
|
C |
1496 |
|
|
C |
1497 |
|
|
C Input: |
1498 |
|
|
C ------ |
1499 |
|
|
C |
1500 |
molod |
1.16 |
C lng : The inner dimension of update and input arrays. |
1501 |
molod |
1.1 |
C |
1502 |
|
|
C LENC : The run: the number of soundings processes in a single call. |
1503 |
molod |
1.16 |
C RAS works on the first LENC of the lng soundings |
1504 |
molod |
1.1 |
C passed. This allows working on pieces of the world |
1505 |
|
|
C say for multitasking, without declaring temporary arrays |
1506 |
|
|
C and copying the data to and from them. This is an f77 |
1507 |
|
|
C version. An F90 version would have to allow more |
1508 |
|
|
C flexibility in the argument declarations. Obviously |
1509 |
molod |
1.16 |
C (LENC<=lng). |
1510 |
molod |
1.1 |
C |
1511 |
|
|
C K : Number of vertical layers (increasing downwards). |
1512 |
|
|
C Need not be the same as the number of layers in the |
1513 |
|
|
C GCM, since it is the outer dimension. The bottom layer |
1514 |
|
|
C (K) is the subcloud layer. |
1515 |
|
|
C |
1516 |
|
|
C IC : Detrainment level to check for presence of convection |
1517 |
|
|
C |
1518 |
|
|
C RASALF : Relaxation parameter (< 1.) for present cloud-type |
1519 |
|
|
C |
1520 |
|
|
C SETRAS : Logical parameter to control re-calculation of |
1521 |
|
|
C saturation specific humidity and mid level P**kappa |
1522 |
|
|
C |
1523 |
|
|
C FRAC : Fraction of the PBL (layer K) mass allowed to be used |
1524 |
|
|
C by a cloud-type in time DT |
1525 |
|
|
C |
1526 |
|
|
C CP : Specific heat at constant pressure |
1527 |
|
|
C |
1528 |
|
|
C ALHL : Latent Heat of condensation |
1529 |
|
|
C |
1530 |
|
|
C RKAP : R/Cp, where R is the gas constant |
1531 |
|
|
C |
1532 |
|
|
C GRAV : Acceleration due to gravity |
1533 |
|
|
C |
1534 |
|
|
C P00 : A reference pressure in hPa, useually 1000 hPa |
1535 |
|
|
C |
1536 |
|
|
C CRTMSF : Critical value of mass flux above which cloudiness at |
1537 |
|
|
C the detrainment layer of that cloud-type is assumed. |
1538 |
|
|
C Affects only cloudiness calculation. |
1539 |
|
|
C |
1540 |
molod |
1.16 |
C POI : 2D array of dimension (lng,K) containing potential |
1541 |
molod |
1.1 |
C temperature. Updated but not initialized by RAS. |
1542 |
|
|
C |
1543 |
molod |
1.16 |
C QOI : 2D array of dimension (lng,K) containing specific |
1544 |
molod |
1.1 |
C humidity. Updated but not initialized by RAS. |
1545 |
|
|
C |
1546 |
molod |
1.16 |
C UOI : 3D array of dimension (lng,K,NTRACER) containing tracers |
1547 |
molod |
1.1 |
C Updated but not initialized by RAS. |
1548 |
|
|
C |
1549 |
molod |
1.16 |
C PRS : 2D array of dimension (lng,K+1) containing pressure |
1550 |
molod |
1.1 |
C in hPa at the interfaces of K-layers from top of the |
1551 |
|
|
C atmosphere to the bottom. Not modified. |
1552 |
|
|
C |
1553 |
molod |
1.16 |
C PRJ : 2D array of dimension (lng,K+1) containing (PRS/P00) ** |
1554 |
molod |
1.1 |
C RKAP. i.e. Exner function at layer edges. Not modified. |
1555 |
|
|
C |
1556 |
molod |
1.16 |
C rhfrac : 1D array of dimension (lng) containing a rel.hum. scaling |
1557 |
molod |
1.1 |
C fraction. Not modified. |
1558 |
|
|
C |
1559 |
|
|
C Output: |
1560 |
|
|
C ------- |
1561 |
|
|
C |
1562 |
molod |
1.16 |
C PCU : 1D array of length lng containing accumulated |
1563 |
molod |
1.1 |
C precipitation in mm/sec. |
1564 |
|
|
C |
1565 |
molod |
1.16 |
C CLN : 2D array of dimension (lng,K) containing cloudiness |
1566 |
molod |
1.1 |
C Note: CLN is bumped but NOT initialized |
1567 |
|
|
C |
1568 |
molod |
1.16 |
C TCU : 2D array of dimension (lng,K) containing accumulated |
1569 |
molod |
1.1 |
C convective heating (K/sec). |
1570 |
|
|
C |
1571 |
molod |
1.16 |
C QCU : 2D array of dimension (lng,K) containing accumulated |
1572 |
molod |
1.1 |
C convective drying (kg/kg/sec). |
1573 |
|
|
C |
1574 |
molod |
1.16 |
C CMASS : 2D array of dimension (lng,K) containing the |
1575 |
molod |
1.1 |
C cloud mass flux (kg/sec). Filled from cloud top |
1576 |
|
|
C to base. |
1577 |
|
|
C |
1578 |
|
|
C Temporaries: |
1579 |
|
|
C |
1580 |
|
|
C ALF, BET, GAM, ETA, PRH, PRI, HOI, HST, QOL, GMH are temporary |
1581 |
|
|
C 2D real arrays of dimension of at least (LENC,K) where LENC is |
1582 |
|
|
C the horizontal dimension over which convection is invoked. |
1583 |
|
|
C |
1584 |
|
|
C |
1585 |
|
|
C TX1, TX2, TX3, TX4, TX5, TX6, TX7, TX8, TX9, AKM, QS1, CLF, UHT |
1586 |
|
|
C VHT, WLQ WFN are temporary real arrays of length at least LENC |
1587 |
|
|
C |
1588 |
|
|
C IA, I1, and I2 are temporary integer arrays of length LENC |
1589 |
|
|
C |
1590 |
|
|
C |
1591 |
|
|
C************************************************************************ |
1592 |
molod |
1.9 |
implicit none |
1593 |
|
|
C Argument List declarations |
1594 |
molod |
1.21 |
integer nn,lng,LENC,K,NLTOP,nlayr,ic,ntracedim, ntracer |
1595 |
molod |
1.13 |
_RL rasalf |
1596 |
molod |
1.9 |
LOGICAL SETRAS |
1597 |
molod |
1.13 |
_RL frac, cp, alhl, rkap, grav, p00, crtmsf |
1598 |
molod |
1.16 |
_RL POI(lng,K),QOI(lng,K),PRS(lng,K+1),PRJ(lng,K+1) |
1599 |
molod |
1.21 |
_RL uoi(lng,nlayr,ntracedim) |
1600 |
molod |
1.16 |
_RL PCU(LENC), CLN(lng) |
1601 |
molod |
1.23 |
_RL TCU(lng,K),QCU(lng,K),ucu(lng,k,ntracedim),CMASS(lng,K) |
1602 |
molod |
1.16 |
_RL ALF(lng,K), BET(lng,K), GAM(lng,K), PRH(lng,K), PRI(lng,K) |
1603 |
molod |
1.13 |
_RL HOL(LENC,K), ETA(LENC,K), HST(LENC,K), QOL(LENC,K) |
1604 |
|
|
_RL GMH(LENC,K) |
1605 |
|
|
_RL TX1(LENC), TX2(LENC), TX3(LENC), TX4(LENC) |
1606 |
|
|
_RL TX5(LENC), TX6(LENC), TX7(LENC), TX8(LENC) |
1607 |
|
|
_RL ALM(LENC), WFN(LENC), AKM(LENC), QS1(LENC) |
1608 |
|
|
_RL WLQ(LENC), CLF(LENC) |
1609 |
molod |
1.21 |
_RL uht(lng,ntracedim) |
1610 |
molod |
1.9 |
integer IA(LENC), I1(LENC),I2(LENC) |
1611 |
molod |
1.16 |
_RL rhfrac(lng) |
1612 |
molod |
1.1 |
|
1613 |
molod |
1.9 |
C Local Variables |
1614 |
molod |
1.13 |
_RL daylen,half,one,zero,cmb2pa,rhmax |
1615 |
molod |
1.1 |
PARAMETER (DAYLEN=86400.0, HALF=0.5, ONE=1.0, ZERO=0.0) |
1616 |
|
|
PARAMETER (CMB2PA=100.0) |
1617 |
|
|
PARAMETER (RHMAX=0.9999) |
1618 |
molod |
1.13 |
_RL rkapp1,onebcp,albcp,onebg,cpbg,twobal |
1619 |
molod |
1.1 |
C |
1620 |
molod |
1.9 |
integer nt,km1,ic1,i,L,len1,len2,isav,len11,ii |
1621 |
molod |
1.16 |
integer lena,lena1,lenb |
1622 |
|
|
_RL tem,tem1 |
1623 |
molod |
1.1 |
|
1624 |
|
|
c Explicit Inline Directives |
1625 |
|
|
c -------------------------- |
1626 |
molod |
1.10 |
#ifdef CRAY |
1627 |
|
|
#ifdef f77 |
1628 |
molod |
1.1 |
cfpp$ expand (qsat) |
1629 |
|
|
#endif |
1630 |
|
|
#endif |
1631 |
|
|
|
1632 |
|
|
RKAPP1 = 1.0 + RKAP |
1633 |
|
|
ONEBCP = 1.0 / CP |
1634 |
|
|
ALBCP = ALHL * ONEBCP |
1635 |
|
|
ONEBG = 1.0 / GRAV |
1636 |
|
|
CPBG = CP * ONEBG |
1637 |
|
|
TWOBAL = 2.0 / ALHL |
1638 |
molod |
1.21 |
|
1639 |
molod |
1.1 |
C |
1640 |
|
|
KM1 = K - 1 |
1641 |
|
|
IC1 = IC + 1 |
1642 |
|
|
C |
1643 |
molod |
1.9 |
C SETTING ALF, BET, GAM, PRH, AND PRI : DONE ONLY WHEN SETRAS=.T. |
1644 |
molod |
1.1 |
C |
1645 |
|
|
|
1646 |
|
|
IF (SETRAS) THEN |
1647 |
|
|
|
1648 |
|
|
DO 2050 L=1,K |
1649 |
|
|
DO 2030 I=1,LENC |
1650 |
|
|
PRH(I,L) = (PRJ(I,L+1)*PRS(I,L+1) - PRJ(I,L)*PRS(I,L)) |
1651 |
|
|
* / ((PRS(I,L+1)-PRS(I,L)) * RKAPP1) |
1652 |
|
|
2030 CONTINUE |
1653 |
|
|
2050 CONTINUE |
1654 |
|
|
|
1655 |
|
|
DO 2070 L=1,K |
1656 |
|
|
DO 2060 I=1,LENC |
1657 |
|
|
TX5(I) = POI(I,L) * PRH(I,L) |
1658 |
|
|
TX1(I) = (PRS(I,L) + PRS(I,L+1)) * 0.5 |
1659 |
|
|
TX3(I) = TX5(I) |
1660 |
|
|
CALL QSAT(TX3(I), TX1(I), TX2(I), TX4(I), .TRUE.) |
1661 |
|
|
ALF(I,L) = TX2(I) - TX4(I) * TX5(I) |
1662 |
|
|
BET(I,L) = TX4(I) * PRH(I,L) |
1663 |
|
|
GAM(I,L) = 1.0 / ((1.0 + TX4(I)*ALBCP) * PRH(I,L)) |
1664 |
|
|
PRI(I,L) = (CP/CMB2PA) / (PRS(I,L+1) - PRS(I,L)) |
1665 |
|
|
2060 CONTINUE |
1666 |
|
|
2070 CONTINUE |
1667 |
|
|
|
1668 |
|
|
ENDIF |
1669 |
|
|
C |
1670 |
|
|
C |
1671 |
|
|
DO 10 L=1,K |
1672 |
molod |
1.16 |
DO 10 I=1,lng |
1673 |
molod |
1.1 |
TCU(I,L) = 0.0 |
1674 |
|
|
QCU(I,L) = 0.0 |
1675 |
|
|
CMASS(I,L) = 0.0 |
1676 |
|
|
10 CONTINUE |
1677 |
|
|
|
1678 |
molod |
1.21 |
do nt = 1,ntracer |
1679 |
|
|
do L=1,K |
1680 |
|
|
do I=1,LENC |
1681 |
|
|
ucu(I,L,nt) = 0.0 |
1682 |
|
|
enddo |
1683 |
|
|
enddo |
1684 |
|
|
enddo |
1685 |
molod |
1.1 |
C |
1686 |
|
|
DO 30 I=1,LENC |
1687 |
|
|
TX1(I) = PRJ(I,K+1) * POI(I,K) |
1688 |
|
|
QS1(I) = ALF(I,K) + BET(I,K)*POI(I,K) |
1689 |
|
|
QOL(I,K) = MIN(QS1(I)*RHMAX,QOI(I,K)) |
1690 |
|
|
|
1691 |
|
|
HOL(I,K) = TX1(I)*CP + QOL(I,K)*ALHL |
1692 |
|
|
ETA(I,K) = ZERO |
1693 |
|
|
TX2(I) = (PRJ(I,K+1) - PRJ(I,K)) * POI(I,K) * CP |
1694 |
|
|
30 CONTINUE |
1695 |
|
|
C |
1696 |
|
|
IF (IC .LT. KM1) THEN |
1697 |
|
|
DO 3703 L=KM1,IC1,-1 |
1698 |
|
|
DO 50 I=1,LENC |
1699 |
|
|
QS1(I) = ALF(I,L) + BET(I,L)*POI(I,L) |
1700 |
|
|
QOL(I,L) = MIN(QS1(I)*RHMAX,QOI(I,L)) |
1701 |
|
|
C |
1702 |
|
|
TEM1 = TX2(I) + PRJ(I,L+1) * POI(I,L) * CP |
1703 |
|
|
HOL(I,L) = TEM1 + QOL(I,L )* ALHL |
1704 |
|
|
HST(I,L) = TEM1 + QS1(I) * ALHL |
1705 |
|
|
|
1706 |
|
|
TX1(I) = (PRJ(I,L+1) - PRJ(I,L)) * POI(I,L) |
1707 |
|
|
ETA(I,L) = ETA(I,L+1) + TX1(I)*CPBG |
1708 |
|
|
TX2(I) = TX2(I) + TX1(I)*CP |
1709 |
|
|
50 CONTINUE |
1710 |
|
|
C |
1711 |
|
|
3703 CONTINUE |
1712 |
|
|
ENDIF |
1713 |
|
|
|
1714 |
|
|
|
1715 |
|
|
DO 70 I=1,LENC |
1716 |
|
|
HOL(I,IC) = TX2(I) |
1717 |
|
|
QS1(I) = ALF(I,IC) + BET(I,IC)*POI(I,IC) |
1718 |
|
|
QOL(I,IC) = MIN(QS1(I)*RHMAX,QOI(I,IC)) |
1719 |
|
|
c |
1720 |
|
|
TEM1 = TX2(I) + PRJ(I,IC1) * POI(I,IC) * CP |
1721 |
|
|
HOL(I,IC) = TEM1 + QOL(I,IC) * ALHL |
1722 |
|
|
HST(I,IC) = TEM1 + QS1(I) * ALHL |
1723 |
|
|
C |
1724 |
|
|
TX3(I ) = (PRJ(I,IC1) - PRH(I,IC)) * POI(I,IC) |
1725 |
|
|
ETA(I,IC) = ETA(I,IC1) + CPBG * TX3(I) |
1726 |
|
|
70 CONTINUE |
1727 |
|
|
C |
1728 |
|
|
DO 130 I=1,LENC |
1729 |
|
|
TX2(I) = HOL(I,K) - HST(I,IC) |
1730 |
|
|
TX1(I) = ZERO |
1731 |
|
|
|
1732 |
|
|
130 CONTINUE |
1733 |
|
|
C |
1734 |
|
|
C ENTRAINMENT PARAMETER |
1735 |
|
|
C |
1736 |
|
|
DO 160 L=IC,KM1 |
1737 |
|
|
DO 160 I=1,LENC |
1738 |
|
|
TX1(I) = TX1(I) + (HST(I,IC) - HOL(I,L)) * (ETA(I,L) - ETA(I,L+1)) |
1739 |
|
|
160 CONTINUE |
1740 |
|
|
C |
1741 |
|
|
LEN1 = 0 |
1742 |
|
|
LEN2 = 0 |
1743 |
|
|
ISAV = 0 |
1744 |
|
|
DO 195 I=1,LENC |
1745 |
|
|
IF (TX1(I) .GT. ZERO .AND. TX2(I) .GT. ZERO |
1746 |
|
|
. .AND. rhfrac(i).ne.0.0 ) THEN |
1747 |
|
|
LEN1 = LEN1 + 1 |
1748 |
|
|
IA(LEN1) = I |
1749 |
|
|
ALM(LEN1) = TX2(I) / TX1(I) |
1750 |
|
|
ENDIF |
1751 |
|
|
195 CONTINUE |
1752 |
|
|
C |
1753 |
|
|
LEN2 = LEN1 |
1754 |
|
|
if (IC1 .lt. K) then |
1755 |
|
|
DO 196 I=1,LENC |
1756 |
|
|
IF (TX2(I) .LE. 0.0 .AND. (HOL(I,K) .GT. HST(I,IC1)) |
1757 |
|
|
. .AND. rhfrac(i).ne.0.0 ) THEN |
1758 |
|
|
LEN2 = LEN2 + 1 |
1759 |
|
|
IA(LEN2) = I |
1760 |
|
|
ALM(LEN2) = 0.0 |
1761 |
|
|
ENDIF |
1762 |
|
|
196 CONTINUE |
1763 |
|
|
endif |
1764 |
|
|
C |
1765 |
|
|
IF (LEN2 .EQ. 0) THEN |
1766 |
|
|
DO 5010 I=1,LENC*K |
1767 |
|
|
HST(I,1) = 0.0 |
1768 |
|
|
QOL(I,1) = 0.0 |
1769 |
|
|
5010 CONTINUE |
1770 |
|
|
DO 5020 I=1,LENC |
1771 |
|
|
PCU(I) = 0.0 |
1772 |
|
|
5020 CONTINUE |
1773 |
|
|
RETURN |
1774 |
|
|
ENDIF |
1775 |
|
|
LEN11 = LEN1 + 1 |
1776 |
|
|
C |
1777 |
|
|
C NORMALIZED MASSFLUX |
1778 |
|
|
C |
1779 |
|
|
DO 250 I=1,LEN2 |
1780 |
|
|
ETA(I,K) = 1.0 |
1781 |
|
|
II = IA(I) |
1782 |
|
|
TX2(I) = 0.5 * (PRS(II,IC) + PRS(II,IC1)) |
1783 |
|
|
TX4(I) = PRS(II,K) |
1784 |
|
|
250 CONTINUE |
1785 |
|
|
C |
1786 |
|
|
DO 252 I=LEN11,LEN2 |
1787 |
|
|
WFN(I) = 0.0 |
1788 |
|
|
II = IA(I) |
1789 |
|
|
IF (HST(II,IC1) .LT. HST(II,IC)) THEN |
1790 |
|
|
TX6(I) = (HST(II,IC1)-HOL(II,K))/(HST(II,IC1)-HST(II,IC)) |
1791 |
|
|
ELSE |
1792 |
|
|
TX6(I) = 0.0 |
1793 |
|
|
ENDIF |
1794 |
|
|
TX2(I) = 0.5 * (PRS(II,IC1)+PRS(II,IC1+1)) * (1.0-TX6(I)) |
1795 |
|
|
* + TX2(I) * TX6(I) |
1796 |
|
|
252 CONTINUE |
1797 |
|
|
C |
1798 |
|
|
CALL ACRITN(LEN2, TX2, TX4, TX3) |
1799 |
|
|
C |
1800 |
|
|
DO 260 L=KM1,IC,-1 |
1801 |
|
|
DO 255 I=1,LEN2 |
1802 |
|
|
TX1(I) = ETA(IA(I),L) |
1803 |
|
|
255 CONTINUE |
1804 |
|
|
DO 260 I=1,LEN2 |
1805 |
|
|
ETA(I,L) = 1.0 + ALM(I) * TX1(I) |
1806 |
|
|
260 CONTINUE |
1807 |
|
|
C |
1808 |
|
|
C CLOUD WORKFUNCTION |
1809 |
|
|
C |
1810 |
|
|
IF (LEN1 .GT. 0) THEN |
1811 |
|
|
DO 270 I=1,LEN1 |
1812 |
|
|
II = IA(I) |
1813 |
|
|
WFN(I) = - GAM(II,IC) * (PRJ(II,IC1) - PRH(II,IC)) |
1814 |
|
|
* * HST(II,IC) * ETA(I,IC1) |
1815 |
|
|
270 CONTINUE |
1816 |
|
|
ENDIF |
1817 |
|
|
C |
1818 |
|
|
DO 290 I=1,LEN2 |
1819 |
|
|
II = IA(I) |
1820 |
|
|
TX1(I) = HOL(II,K) |
1821 |
|
|
290 CONTINUE |
1822 |
|
|
C |
1823 |
|
|
IF (IC1 .LE. KM1) THEN |
1824 |
|
|
|
1825 |
|
|
DO 380 L=KM1,IC1,-1 |
1826 |
|
|
DO 380 I=1,LEN2 |
1827 |
|
|
II = IA(I) |
1828 |
|
|
TEM = TX1(I) + (ETA(I,L) - ETA(I,L+1)) * HOL(II,L) |
1829 |
|
|
C |
1830 |
|
|
PCU(I) = PRJ(II,L+1) - PRH(II,L) |
1831 |
|
|
TEM1 = ETA(I,L+1) * PCU(I) |
1832 |
|
|
TX1(I) = TX1(I)*PCU(I) |
1833 |
|
|
C |
1834 |
|
|
PCU(I) = PRH(II,L) - PRJ(II,L) |
1835 |
|
|
TEM1 = (TEM1 + ETA(I,L) * PCU(I)) * HST(II,L) |
1836 |
|
|
TX1(I) = TX1(I) + TEM*PCU(I) |
1837 |
|
|
C |
1838 |
|
|
WFN(I) = WFN(I) + (TX1(I) - TEM1) * GAM(II,L) |
1839 |
|
|
TX1(I) = TEM |
1840 |
|
|
380 CONTINUE |
1841 |
|
|
ENDIF |
1842 |
|
|
C |
1843 |
|
|
LENA = 0 |
1844 |
|
|
IF (LEN1 .GT. 0) THEN |
1845 |
|
|
DO 512 I=1,LEN1 |
1846 |
|
|
II = IA(I) |
1847 |
|
|
WFN(I) = WFN(I) + TX1(I) * GAM(II,IC)*(PRJ(II,IC1)-PRH(II,IC)) |
1848 |
|
|
* - TX3(I) |
1849 |
|
|
IF (WFN(I) .GT. 0.0) THEN |
1850 |
|
|
LENA = LENA + 1 |
1851 |
|
|
I1(LENA) = IA(I) |
1852 |
|
|
I2(LENA) = I |
1853 |
|
|
TX1(LENA) = WFN(I) |
1854 |
|
|
TX2(LENA) = QS1(IA(I)) |
1855 |
|
|
TX6(LENA) = 1.0 |
1856 |
|
|
ENDIF |
1857 |
|
|
512 CONTINUE |
1858 |
|
|
ENDIF |
1859 |
|
|
LENB = LENA |
1860 |
|
|
DO 515 I=LEN11,LEN2 |
1861 |
|
|
WFN(I) = WFN(I) - TX3(I) |
1862 |
|
|
IF (WFN(I) .GT. 0.0 .AND. TX6(I) .GT. 0.0) THEN |
1863 |
|
|
LENB = LENB + 1 |
1864 |
|
|
I1(LENB) = IA(I) |
1865 |
|
|
I2(LENB) = I |
1866 |
|
|
TX1(LENB) = WFN(I) |
1867 |
|
|
TX2(LENB) = QS1(IA(I)) |
1868 |
|
|
TX4(LENB) = TX6(I) |
1869 |
|
|
ENDIF |
1870 |
|
|
515 CONTINUE |
1871 |
|
|
C |
1872 |
|
|
IF (LENB .LE. 0) THEN |
1873 |
|
|
DO 5030 I=1,LENC*K |
1874 |
|
|
HST(I,1) = 0.0 |
1875 |
|
|
QOL(I,1) = 0.0 |
1876 |
|
|
5030 CONTINUE |
1877 |
|
|
DO 5040 I=1,LENC |
1878 |
|
|
PCU(I) = 0.0 |
1879 |
|
|
5040 CONTINUE |
1880 |
|
|
RETURN |
1881 |
|
|
ENDIF |
1882 |
|
|
|
1883 |
|
|
C |
1884 |
|
|
DO 516 I=1,LENB |
1885 |
|
|
WFN(I) = TX1(I) |
1886 |
|
|
QS1(I) = TX2(I) |
1887 |
|
|
516 CONTINUE |
1888 |
|
|
C |
1889 |
|
|
DO 520 L=IC,K |
1890 |
|
|
DO 517 I=1,LENB |
1891 |
|
|
TX1(I) = ETA(I2(I),L) |
1892 |
|
|
517 CONTINUE |
1893 |
|
|
DO 520 I=1,LENB |
1894 |
|
|
ETA(I,L) = TX1(I) |
1895 |
|
|
520 CONTINUE |
1896 |
|
|
C |
1897 |
|
|
LENA1 = LENA + 1 |
1898 |
|
|
C |
1899 |
|
|
DO 510 I=1,LENA |
1900 |
|
|
II = I1(I) |
1901 |
|
|
TX8(I) = HST(II,IC) - HOL(II,IC) |
1902 |
|
|
510 CONTINUE |
1903 |
|
|
DO 530 I=LENA1,LENB |
1904 |
|
|
II = I1(I) |
1905 |
|
|
TX6(I) = TX4(I) |
1906 |
|
|
TEM = TX6(I) * (HOL(II,IC)-HOL(II,IC1)) + HOL(II,IC1) |
1907 |
|
|
TX8(I) = HOL(II,K) - TEM |
1908 |
|
|
|
1909 |
|
|
TEM1 = TX6(I) * (QOL(II,IC)-QOL(II,IC1)) + QOL(II,IC1) |
1910 |
|
|
TX5(I) = TEM - TEM1 * ALHL |
1911 |
|
|
QS1(I) = TEM1 + TX8(I)*(ONE/ALHL) |
1912 |
|
|
TX3(I) = HOL(II,IC) |
1913 |
|
|
530 CONTINUE |
1914 |
|
|
C |
1915 |
|
|
C |
1916 |
|
|
DO 620 I=1,LENB |
1917 |
|
|
II = I1(I) |
1918 |
|
|
WLQ(I) = QOL(II,K) - QS1(I) * ETA(I,IC) |
1919 |
|
|
TX7(I) = HOL(II,K) |
1920 |
|
|
620 CONTINUE |
1921 |
molod |
1.21 |
DO NT=1,Ntracer |
1922 |
|
|
DO 621 I=1,LENB |
1923 |
|
|
II = I1(I) |
1924 |
|
|
UHT(I,NT) = UOI(II,K+nltop-1,NT)-UOI(II,IC+nltop-1,NT) * ETA(I,IC) |
1925 |
|
|
621 CONTINUE |
1926 |
|
|
ENDDO |
1927 |
molod |
1.1 |
C |
1928 |
|
|
DO 635 L=KM1,IC,-1 |
1929 |
|
|
DO 630 I=1,LENB |
1930 |
|
|
II = I1(I) |
1931 |
|
|
TEM = ETA(I,L) - ETA(I,L+1) |
1932 |
|
|
WLQ(I) = WLQ(I) + TEM * QOL(II,L) |
1933 |
|
|
630 CONTINUE |
1934 |
|
|
635 CONTINUE |
1935 |
molod |
1.21 |
DO NT=1,Ntracer |
1936 |
|
|
DO L=KM1,IC,-1 |
1937 |
|
|
DO I=1,LENB |
1938 |
|
|
II = I1(I) |
1939 |
|
|
TEM = ETA(I,L) - ETA(I,L+1) |
1940 |
|
|
UHT(I,NT) = UHT(I,NT) + TEM * UOI(II,L+nltop-1,NT) |
1941 |
|
|
ENDDO |
1942 |
|
|
ENDDO |
1943 |
|
|
ENDDO |
1944 |
molod |
1.1 |
C |
1945 |
|
|
C CALCULATE GS AND PART OF AKM (THAT REQUIRES ETA) |
1946 |
|
|
C |
1947 |
|
|
DO 690 I=1,LENB |
1948 |
|
|
II = I1(I) |
1949 |
|
|
c TX7(I) = HOL(II,K) |
1950 |
|
|
TEM = (POI(II,KM1) - POI(II,K)) / (PRH(II,K) - PRH(II,KM1)) |
1951 |
|
|
HOL(I,K) = TEM * (PRJ(II,K)-PRH(II,KM1))*PRH(II,K)*PRI(II,K) |
1952 |
|
|
HOL(I,KM1) = TEM * (PRH(II,K)-PRJ(II,K))*PRH(II,KM1)*PRI(II,KM1) |
1953 |
|
|
AKM(I) = ZERO |
1954 |
|
|
TX2(I) = 0.5 * (PRS(II,IC) + PRS(II,IC1)) |
1955 |
|
|
690 CONTINUE |
1956 |
|
|
|
1957 |
|
|
IF (IC1 .LE. KM1) THEN |
1958 |
|
|
DO 750 L=KM1,IC1,-1 |
1959 |
|
|
DO 750 I=1,LENB |
1960 |
|
|
II = I1(I) |
1961 |
|
|
TEM = (POI(II,L-1) - POI(II,L)) * ETA(I,L) |
1962 |
|
|
* / (PRH(II,L) - PRH(II,L-1)) |
1963 |
|
|
C |
1964 |
|
|
HOL(I,L) = TEM * (PRJ(II,L)-PRH(II,L-1)) * PRH(II,L) |
1965 |
|
|
* * PRI(II,L) + HOL(I,L) |
1966 |
|
|
HOL(I,L-1) = TEM * (PRH(II,L)-PRJ(II,L)) * PRH(II,L-1) |
1967 |
|
|
* * PRI(II,L-1) |
1968 |
|
|
C |
1969 |
|
|
AKM(I) = AKM(I) - HOL(I,L) |
1970 |
|
|
* * (ETA(I,L) * (PRH(II,L)-PRJ(II,L)) + |
1971 |
|
|
* ETA(I,L+1) * (PRJ(II,L+1)-PRH(II,L))) / PRH(II,L) |
1972 |
|
|
750 CONTINUE |
1973 |
|
|
ENDIF |
1974 |
|
|
C |
1975 |
|
|
C |
1976 |
|
|
CALL RNCL(LENB, TX2, TX1, CLF) |
1977 |
|
|
|
1978 |
|
|
DO 770 I=1,LENB |
1979 |
|
|
TX2(I) = (ONE - TX1(I)) * WLQ(I) |
1980 |
|
|
WLQ(I) = TX1(I) * WLQ(I) |
1981 |
|
|
C |
1982 |
|
|
TX1(I) = HOL(I,IC) |
1983 |
|
|
770 CONTINUE |
1984 |
|
|
DO 790 I=LENA1, LENB |
1985 |
|
|
II = I1(I) |
1986 |
|
|
TX1(I) = TX1(I) + (TX5(I)-TX3(I)+QOL(II,IC)*ALHL)*(PRI(II,IC)/CP) |
1987 |
|
|
790 CONTINUE |
1988 |
|
|
|
1989 |
|
|
DO 800 I=1,LENB |
1990 |
|
|
HOL(I,IC) = TX1(I) - TX2(I) * ALBCP * PRI(I1(I),IC) |
1991 |
|
|
800 CONTINUE |
1992 |
|
|
|
1993 |
|
|
IF (LENA .GT. 0) THEN |
1994 |
|
|
DO 810 I=1,LENA |
1995 |
|
|
II = I1(I) |
1996 |
|
|
AKM(I) = AKM(I) - ETA(I,IC1) * (PRJ(II,IC1) - PRH(II,IC)) |
1997 |
|
|
* * TX1(I) / PRH(II,IC) |
1998 |
|
|
810 CONTINUE |
1999 |
|
|
ENDIF |
2000 |
|
|
c |
2001 |
|
|
C CALCULATE GH |
2002 |
|
|
C |
2003 |
|
|
DO 830 I=1,LENB |
2004 |
|
|
II = I1(I) |
2005 |
|
|
TX3(I) = QOL(II,KM1) - QOL(II,K) |
2006 |
|
|
GMH(I,K) = HOL(I,K) + TX3(I) * PRI(II,K) * (ALBCP) |
2007 |
|
|
|
2008 |
|
|
AKM(I) = AKM(I) + GAM(II,KM1)*(PRJ(II,K)-PRH(II,KM1)) |
2009 |
|
|
* * GMH(I,K) |
2010 |
|
|
TX3(I) = zero |
2011 |
|
|
830 CONTINUE |
2012 |
|
|
C |
2013 |
|
|
IF (IC1 .LE. KM1) THEN |
2014 |
|
|
DO 840 L=KM1,IC1,-1 |
2015 |
|
|
DO 840 I=1,LENB |
2016 |
|
|
II = I1(I) |
2017 |
|
|
TX2(I) = TX3(I) |
2018 |
|
|
TX3(I) = (QOL(II,L-1) - QOL(II,L)) * ETA(I,L) |
2019 |
|
|
TX2(I) = TX2(I) + TX3(I) |
2020 |
|
|
C |
2021 |
|
|
GMH(I,L) = HOL(I,L) + TX2(I) * PRI(II,L) * (ALBCP*HALF) |
2022 |
|
|
840 CONTINUE |
2023 |
|
|
C |
2024 |
|
|
C |
2025 |
|
|
ENDIF |
2026 |
|
|
DO 850 I=LENA1,LENB |
2027 |
|
|
TX3(I) = TX3(I) + TWOBAL |
2028 |
|
|
* * (TX7(I) - TX8(I) - TX5(I) - QOL(I1(I),IC)*ALHL) |
2029 |
|
|
850 CONTINUE |
2030 |
|
|
DO 860 I=1,LENB |
2031 |
|
|
GMH(I,IC) = TX1(I) + PRI(I1(I),IC) * ONEBCP |
2032 |
|
|
* * (TX3(I)*(ALHL*HALF) + ETA(I,IC) * TX8(I)) |
2033 |
|
|
860 CONTINUE |
2034 |
|
|
C |
2035 |
|
|
C CALCULATE HC PART OF AKM |
2036 |
|
|
C |
2037 |
|
|
IF (IC1 .LE. KM1) THEN |
2038 |
|
|
DO 870 I=1,LENB |
2039 |
|
|
TX1(I) = GMH(I,K) |
2040 |
|
|
870 CONTINUE |
2041 |
|
|
DO 3725 L=KM1,IC1,-1 |
2042 |
|
|
DO 880 I=1,LENB |
2043 |
|
|
II = I1(I) |
2044 |
|
|
TX1(I) = TX1(I) + (ETA(I,L) - ETA(I,L+1)) * GMH(I,L) |
2045 |
|
|
TX2(I) = GAM(II,L-1) * (PRJ(II,L) - PRH(II,L-1)) |
2046 |
|
|
880 CONTINUE |
2047 |
|
|
C |
2048 |
|
|
IF (L .EQ. IC1) THEN |
2049 |
|
|
DO 890 I=LENA1,LENB |
2050 |
|
|
TX2(I) = ZERO |
2051 |
|
|
890 CONTINUE |
2052 |
|
|
ENDIF |
2053 |
|
|
DO 900 I=1,LENB |
2054 |
|
|
II = I1(I) |
2055 |
|
|
AKM(I) = AKM(I) + TX1(I) * |
2056 |
|
|
* (TX2(I) + GAM(II,L)*(PRH(II,L)-PRJ(II,L))) |
2057 |
|
|
900 CONTINUE |
2058 |
|
|
3725 CONTINUE |
2059 |
|
|
ENDIF |
2060 |
|
|
C |
2061 |
|
|
DO 920 I=LENA1,LENB |
2062 |
|
|
II = I1(I) |
2063 |
|
|
TX2(I) = 0.5 * (PRS(II,IC) + PRS(II,IC1)) |
2064 |
|
|
* + 0.5*(PRS(II,IC+2) - PRS(II,IC)) * (ONE-TX6(I)) |
2065 |
|
|
c |
2066 |
|
|
TX1(I) = PRS(II,IC1) |
2067 |
|
|
TX5(I) = 0.5 * (PRS(II,IC1) + PRS(II,IC+2)) |
2068 |
|
|
C |
2069 |
|
|
IF ((TX2(I) .GE. TX1(I)) .AND. (TX2(I) .LT. TX5(I))) THEN |
2070 |
|
|
TX6(I) = ONE - (TX2(I) - TX1(I)) / (TX5(I) - TX1(I)) |
2071 |
|
|
C |
2072 |
|
|
TEM = PRI(II,IC1) / PRI(II,IC) |
2073 |
|
|
HOL(I,IC1) = HOL(I,IC1) + HOL(I,IC) * TEM |
2074 |
|
|
HOL(I,IC) = ZERO |
2075 |
|
|
C |
2076 |
|
|
GMH(I,IC1) = GMH(I,IC1) + GMH(I,IC) * TEM |
2077 |
|
|
GMH(I,IC) = ZERO |
2078 |
|
|
ELSEIF (TX2(I) .LT. TX1(I)) THEN |
2079 |
|
|
TX6(I) = 1.0 |
2080 |
|
|
ELSE |
2081 |
|
|
TX6(I) = 0.0 |
2082 |
|
|
ENDIF |
2083 |
|
|
920 CONTINUE |
2084 |
|
|
C |
2085 |
|
|
C |
2086 |
|
|
DO I=1,LENC |
2087 |
|
|
PCU(I) = 0.0 |
2088 |
|
|
ENDDO |
2089 |
|
|
|
2090 |
|
|
DO 970 I=1,LENB |
2091 |
|
|
II = I1(I) |
2092 |
|
|
IF (AKM(I) .LT. ZERO .AND. WLQ(I) .GE. 0.0) THEN |
2093 |
|
|
WFN(I) = - TX6(I) * WFN(I) * RASALF / AKM(I) |
2094 |
|
|
ELSE |
2095 |
|
|
WFN(I) = ZERO |
2096 |
|
|
ENDIF |
2097 |
|
|
TEM = (PRS(II,K+1)-PRS(II,K))*(CMB2PA*FRAC) |
2098 |
|
|
WFN(I) = MIN(WFN(I), TEM) |
2099 |
|
|
C |
2100 |
|
|
C compute cloud amount |
2101 |
|
|
C |
2102 |
|
|
CC TX1(I) = CLN(II) |
2103 |
|
|
CC IF (WFN(I) .GT. CRTMSF) TX1(I) = TX1(I) + CLF(I) |
2104 |
|
|
CC IF (TX1(I) .GT. ONE) TX1(I) = ONE |
2105 |
|
|
C |
2106 |
|
|
C PRECIPITATION |
2107 |
|
|
C |
2108 |
|
|
PCU(II) = WLQ(I) * WFN(I) * ONEBG |
2109 |
|
|
C |
2110 |
|
|
C CUMULUS FRICTION AT THE BOTTOM LAYER |
2111 |
|
|
C |
2112 |
|
|
TX4(I) = WFN(I) * (1.0/ALHL) |
2113 |
|
|
TX5(I) = WFN(I) * ONEBCP |
2114 |
|
|
970 CONTINUE |
2115 |
|
|
C |
2116 |
|
|
C compute cloud mass flux for diagnostic output |
2117 |
|
|
C |
2118 |
|
|
DO L = IC,K |
2119 |
|
|
DO I=1,LENB |
2120 |
|
|
II = I1(I) |
2121 |
|
|
if(L.lt.K)then |
2122 |
|
|
CMASS(II,L) = ETA(I,L+1) * WFN(I) * ONEBG |
2123 |
|
|
else |
2124 |
|
|
CMASS(II,L) = WFN(I) * ONEBG |
2125 |
|
|
endif |
2126 |
|
|
ENDDO |
2127 |
|
|
ENDDO |
2128 |
|
|
C |
2129 |
|
|
CC DO 975 I=1,LENB |
2130 |
|
|
CC II = I1(I) |
2131 |
|
|
CC CLN(II) = TX1(I) |
2132 |
|
|
CC975 CONTINUE |
2133 |
|
|
C |
2134 |
|
|
C THETA AND Q CHANGE DUE TO CLOUD TYPE IC |
2135 |
|
|
C |
2136 |
|
|
|
2137 |
|
|
c TEMA = 0.0 |
2138 |
|
|
c TEMB = 0.0 |
2139 |
|
|
DO 990 L=IC,K |
2140 |
|
|
DO 980 I=1,LENB |
2141 |
|
|
II = I1(I) |
2142 |
|
|
TEM = (GMH(I,L) - HOL(I,L)) * TX4(I) |
2143 |
|
|
TEM1 = HOL(I,L) * TX5(I) |
2144 |
|
|
C |
2145 |
|
|
TCU(II,L) = TEM1 / PRH(II,L) |
2146 |
|
|
QCU(II,L) = TEM |
2147 |
|
|
980 CONTINUE |
2148 |
|
|
|
2149 |
|
|
c I = I1(IP1) |
2150 |
|
|
c |
2151 |
|
|
c TEM = (PRS(I,L+1)-PRS(I,L)) * (ONEBG*100.0) |
2152 |
|
|
c TEMA = TEMA + TCU(I,L) * PRH(I,L) * TEM * (CP/ALHL) |
2153 |
|
|
c TEMB = TEMB + QCU(I,L) * TEM |
2154 |
|
|
C |
2155 |
|
|
990 CONTINUE |
2156 |
|
|
C |
2157 |
|
|
c Compute Tracer Tendencies |
2158 |
|
|
c ------------------------- |
2159 |
molod |
1.21 |
do nt = 1,ntracer |
2160 |
molod |
1.18 |
c |
2161 |
molod |
1.1 |
c Tracer Tendency at the Bottom Layer |
2162 |
|
|
c ----------------------------------- |
2163 |
molod |
1.21 |
DO 995 I=1,LENB |
2164 |
|
|
II = I1(I) |
2165 |
|
|
TEM = half*TX5(I) * PRI(II,K) |
2166 |
|
|
TX1(I) = ( UOI(II,KM1+nltop-1,nt) - UOI(II,K+nltop-1,nt)) |
2167 |
|
|
ucu(II,K,nt) = TEM * TX1(I) |
2168 |
|
|
995 CONTINUE |
2169 |
molod |
1.18 |
c |
2170 |
molod |
1.1 |
c Tracer Tendency at all other Levels |
2171 |
|
|
c ----------------------------------- |
2172 |
molod |
1.21 |
DO 1020 L=KM1,IC1,-1 |
2173 |
|
|
DO 1010 I=1,LENB |
2174 |
|
|
II = I1(I) |
2175 |
|
|
TEM = half*TX5(I) * PRI(II,L) |
2176 |
|
|
TEM1 = TX1(I) |
2177 |
|
|
TX1(I) = (UOI(II,L-1+nltop-1,nt)-UOI(II,L+nltop-1,nt)) * ETA(I,L) |
2178 |
|
|
TX3(I) = (TX1(I) + TEM1) * TEM |
2179 |
|
|
1010 CONTINUE |
2180 |
|
|
DO 1020 I=1,LENB |
2181 |
|
|
II = I1(I) |
2182 |
|
|
ucu(II,L,nt) = TX3(I) |
2183 |
|
|
1020 CONTINUE |
2184 |
molod |
1.22 |
|
2185 |
molod |
1.21 |
DO 1030 I=1,LENB |
2186 |
|
|
II = I1(I) |
2187 |
|
|
IF (TX6(I) .GE. 1.0) THEN |
2188 |
|
|
TEM = half*TX5(I) * PRI(II,IC) |
2189 |
|
|
ELSE |
2190 |
|
|
TEM = 0.0 |
2191 |
|
|
ENDIF |
2192 |
|
|
TX1(I) = (TX1(I) + UHT(I,nt) + UHT(I,nt)) * TEM |
2193 |
|
|
1030 CONTINUE |
2194 |
|
|
DO 1040 I=1,LENB |
2195 |
|
|
II = I1(I) |
2196 |
|
|
ucu(II,IC,nt) = TX1(I) |
2197 |
|
|
1040 CONTINUE |
2198 |
|
|
|
2199 |
|
|
enddo |
2200 |
molod |
1.1 |
C |
2201 |
|
|
C PENETRATIVE CONVECTION CALCULATION OVER |
2202 |
|
|
C |
2203 |
|
|
|
2204 |
|
|
RETURN |
2205 |
|
|
END |
2206 |
molod |
1.16 |
SUBROUTINE RNCL(lng, PL, RNO, CLF) |
2207 |
molod |
1.1 |
C |
2208 |
|
|
C********************************************************************* |
2209 |
|
|
C********************** Relaxed Arakawa-Schubert ********************* |
2210 |
|
|
C************************ SUBROUTINE RNCL ************************ |
2211 |
|
|
C**************************** 23 July 1992 *************************** |
2212 |
|
|
C********************************************************************* |
2213 |
molod |
1.9 |
implicit none |
2214 |
|
|
C Argument List declarations |
2215 |
molod |
1.16 |
integer lng |
2216 |
|
|
_RL PL(lng), RNO(lng), CLF(lng) |
2217 |
molod |
1.1 |
|
2218 |
molod |
1.9 |
C Local Variables |
2219 |
molod |
1.13 |
_RL p5,p8,pt8,pt2,pfac,p4,p6,p7,p9,cucld,cfac |
2220 |
molod |
1.1 |
PARAMETER (P5=500.0, P8=800.0, PT8=0.8, PT2=0.2) |
2221 |
|
|
PARAMETER (PFAC=PT2/(P8-P5)) |
2222 |
|
|
PARAMETER (P4=400.0, P6=401.0) |
2223 |
|
|
PARAMETER (P7=700.0, P9=900.0) |
2224 |
|
|
PARAMETER (CUCLD=0.5,CFAC=CUCLD/(P6-P4)) |
2225 |
molod |
1.9 |
|
2226 |
|
|
integer i |
2227 |
molod |
1.1 |
C |
2228 |
molod |
1.16 |
DO 10 I=1,lng |
2229 |
molod |
1.1 |
rno(i) = 1.0 |
2230 |
|
|
ccc if( pl(i).le.400.0 ) rno(i) = max( 0.75, 1.0-0.0025*(400.0-pl(i)) ) |
2231 |
|
|
|
2232 |
|
|
ccc IF ( PL(I).GE.P7 .AND. PL(I).LE.P9 ) THEN |
2233 |
|
|
ccc RNO(I) = ((P9-PL(I))/(P9-P7)) **2 |
2234 |
|
|
ccc ELSE IF (PL(I).GT.P9) THEN |
2235 |
|
|
ccc RNO(I) = 0. |
2236 |
|
|
ccc ENDIF |
2237 |
|
|
|
2238 |
|
|
CLF(I) = CUCLD |
2239 |
|
|
C |
2240 |
|
|
CARIESIF (PL(I) .GE. P5 .AND. PL(I) .LE. P8) THEN |
2241 |
|
|
CARIES RNO(I) = (P8-PL(I))*PFAC + PT8 |
2242 |
|
|
CARIESELSEIF (PL(I) .GT. P8 ) THEN |
2243 |
|
|
CARIES RNO(I) = PT8 |
2244 |
|
|
CARIESENDIF |
2245 |
|
|
CARIES |
2246 |
|
|
IF (PL(I) .GE. P4 .AND. PL(I) .LE. P6) THEN |
2247 |
|
|
CLF(I) = (P6-PL(I))*CFAC |
2248 |
|
|
ELSEIF (PL(I) .GT. P6 ) THEN |
2249 |
|
|
CLF(I) = 0.0 |
2250 |
|
|
ENDIF |
2251 |
|
|
10 CONTINUE |
2252 |
|
|
C |
2253 |
|
|
RETURN |
2254 |
|
|
END |
2255 |
molod |
1.16 |
SUBROUTINE ACRITN ( lng,PL,PLB,ACR ) |
2256 |
molod |
1.1 |
|
2257 |
|
|
C********************************************************************* |
2258 |
|
|
C********************** Relaxed Arakawa-Schubert ********************* |
2259 |
|
|
C************************** SUBROUTINE ACRIT ********************* |
2260 |
|
|
C****************** modified August 28, 1996 L.Takacs ************ |
2261 |
|
|
C**** ***** |
2262 |
|
|
C**** Note: Data obtained from January Mean After-Analysis ***** |
2263 |
|
|
C**** from 4x5 46-layer GEOS Assimilation ***** |
2264 |
|
|
C**** ***** |
2265 |
|
|
C********************************************************************* |
2266 |
molod |
1.9 |
implicit none |
2267 |
|
|
C Argument List declarations |
2268 |
molod |
1.16 |
integer lng |
2269 |
|
|
_RL PL(lng), PLB(lng), ACR(lng) |
2270 |
molod |
1.1 |
|
2271 |
molod |
1.9 |
C Local variables |
2272 |
|
|
integer lma |
2273 |
molod |
1.1 |
parameter (lma=18) |
2274 |
molod |
1.13 |
_RL p(lma) |
2275 |
|
|
_RL a(lma) |
2276 |
molod |
1.9 |
integer i,L |
2277 |
molod |
1.13 |
_RL temp |
2278 |
molod |
1.1 |
|
2279 |
|
|
data p / 93.81, 111.65, 133.46, 157.80, 186.51, |
2280 |
|
|
. 219.88, 257.40, 301.21, 352.49, 409.76, |
2281 |
|
|
. 471.59, 535.04, 603.33, 672.79, 741.12, |
2282 |
|
|
. 812.52, 875.31, 930.20/ |
2283 |
|
|
|
2284 |
|
|
data a / 3.35848, 3.13645, 2.48072, 2.08277, 1.53364, |
2285 |
|
|
. 1.01971, .65846, .45867, .38687, .31002, |
2286 |
|
|
. .25574, .20347, .17254, .15260, .16756, |
2287 |
|
|
. .09916, .10360, .05880/ |
2288 |
|
|
|
2289 |
|
|
|
2290 |
|
|
do L=1,lma-1 |
2291 |
molod |
1.16 |
do i=1,lng |
2292 |
molod |
1.1 |
if( pl(i).ge.p(L) .and. |
2293 |
|
|
. pl(i).le.p(L+1)) then |
2294 |
|
|
temp = ( pl(i)-p(L) )/( p(L+1)-p(L) ) |
2295 |
|
|
acr(i) = a(L+1)*temp + a(L)*(1-temp) |
2296 |
|
|
endif |
2297 |
|
|
enddo |
2298 |
|
|
enddo |
2299 |
|
|
|
2300 |
molod |
1.16 |
do i=1,lng |
2301 |
molod |
1.1 |
if( pl(i).lt.p(1) ) acr(i) = a(1) |
2302 |
|
|
if( pl(i).gt.p(lma) ) acr(i) = a(lma) |
2303 |
|
|
enddo |
2304 |
|
|
|
2305 |
molod |
1.16 |
do i=1,lng |
2306 |
molod |
1.1 |
acr(i) = acr(i) * (plb(i)-pl(i)) |
2307 |
|
|
enddo |
2308 |
|
|
|
2309 |
|
|
RETURN |
2310 |
|
|
END |
2311 |
molod |
1.6 |
SUBROUTINE RNEVP(NN,IRUN,NLAY,TL,QL,RAIN,PL,CLFRAC,SP,DP,PLKE, |
2312 |
molod |
1.1 |
1 PLK,TH,TEMP1,TEMP2,TEMP3,ITMP1,ITMP2,RCON,RLAR,CLSBTH,tmscl, |
2313 |
|
|
2 tmfrc,cp,gravity,alhl,gamfac,cldlz,RHCRIT,offset,alpha) |
2314 |
|
|
|
2315 |
molod |
1.9 |
implicit none |
2316 |
|
|
C Argument List declarations |
2317 |
|
|
integer nn,irun,nlay |
2318 |
molod |
1.13 |
_RL TL(IRUN,NLAY),QL(IRUN,NLAY),RAIN(IRUN,NLAY), |
2319 |
molod |
1.9 |
. PL(IRUN,NLAY),CLFRAC(IRUN,NLAY),SP(IRUN),TEMP1(IRUN,NLAY), |
2320 |
|
|
. TEMP2(IRUN,NLAY),PLKE(IRUN,NLAY+1), |
2321 |
|
|
. RCON(IRUN),RLAR(IRUN),DP(IRUN,NLAY),PLK(IRUN,NLAY),TH(IRUN,NLAY), |
2322 |
|
|
. TEMP3(IRUN,NLAY) |
2323 |
|
|
integer ITMP1(IRUN,NLAY),ITMP2(IRUN,NLAY) |
2324 |
molod |
1.13 |
_RL CLSBTH(IRUN,NLAY) |
2325 |
|
|
_RL tmscl,tmfrc,cp,gravity,alhl,gamfac,offset,alpha |
2326 |
|
|
_RL cldlz(irun,nlay) |
2327 |
|
|
_RL rhcrit(irun,nlay) |
2328 |
molod |
1.9 |
C |
2329 |
|
|
C Local Variables |
2330 |
molod |
1.13 |
_RL zm1p04,zero,two89,zp44,zp01,half,zp578,one,thousand,z3600 |
2331 |
|
|
_RL zp1,zp001 |
2332 |
molod |
1.1 |
PARAMETER (ZM1P04 = -1.04E-4 ) |
2333 |
|
|
PARAMETER (ZERO = 0.) |
2334 |
|
|
PARAMETER (TWO89= 2.89E-5) |
2335 |
|
|
PARAMETER ( ZP44= 0.44) |
2336 |
|
|
PARAMETER ( ZP01= 0.01) |
2337 |
|
|
PARAMETER ( ZP1 = 0.1 ) |
2338 |
|
|
PARAMETER ( ZP001= 0.001) |
2339 |
|
|
PARAMETER ( HALF= 0.5) |
2340 |
|
|
PARAMETER ( ZP578 = 0.578 ) |
2341 |
|
|
PARAMETER ( ONE = 1.) |
2342 |
|
|
PARAMETER ( THOUSAND = 1000.) |
2343 |
|
|
PARAMETER ( Z3600 = 3600.) |
2344 |
|
|
C |
2345 |
molod |
1.13 |
_RL EVP9(IRUN,NLAY) |
2346 |
|
|
_RL water(irun),crystal(irun) |
2347 |
|
|
_RL watevap(irun),iceevap(irun) |
2348 |
|
|
_RL fracwat,fracice, tice,rh,fact,dum |
2349 |
|
|
_RL rainmax(irun) |
2350 |
|
|
_RL getcon,rphf,elocp,cpog,relax |
2351 |
|
|
_RL exparg,arearat,rpow |
2352 |
molod |
1.9 |
|
2353 |
|
|
integer i,L,n,nlaym1,irnlay,irnlm1 |
2354 |
molod |
1.1 |
|
2355 |
|
|
c Explicit Inline Directives |
2356 |
|
|
c -------------------------- |
2357 |
molod |
1.10 |
#ifdef CRAY |
2358 |
|
|
#ifdef f77 |
2359 |
molod |
1.1 |
cfpp$ expand (qsat) |
2360 |
|
|
#endif |
2361 |
|
|
#endif |
2362 |
|
|
|
2363 |
|
|
tice = getcon('FREEZING-POINT') |
2364 |
|
|
|
2365 |
|
|
fracwat = 0.70 |
2366 |
|
|
fracice = 0.01 |
2367 |
|
|
|
2368 |
|
|
NLAYM1 = NLAY - 1 |
2369 |
|
|
IRNLAY = IRUN*NLAY |
2370 |
|
|
IRNLM1 = IRUN*(NLAY-1) |
2371 |
|
|
|
2372 |
|
|
RPHF = Z3600/tmscl |
2373 |
|
|
|
2374 |
|
|
ELOCP = alhl/cp |
2375 |
|
|
CPOG = cp/gravity |
2376 |
|
|
|
2377 |
|
|
DO I = 1,IRUN |
2378 |
|
|
RLAR(I) = 0. |
2379 |
|
|
water(i) = 0. |
2380 |
|
|
crystal(i) = 0. |
2381 |
|
|
ENDDO |
2382 |
|
|
|
2383 |
|
|
do L = 1,nlay |
2384 |
|
|
do i = 1,irun |
2385 |
|
|
EVP9(i,L) = 0. |
2386 |
|
|
TEMP1(i,L) = 0. |
2387 |
|
|
TEMP2(i,L) = 0. |
2388 |
|
|
TEMP3(i,L) = 0. |
2389 |
|
|
CLSBTH(i,L) = 0. |
2390 |
|
|
cldlz(i,L) = 0. |
2391 |
|
|
enddo |
2392 |
|
|
enddo |
2393 |
|
|
|
2394 |
|
|
C RHO(ZERO) / RHO FOR TERMINAL VELOCITY APPROX. |
2395 |
|
|
c --------------------------------------------- |
2396 |
|
|
DO L = 1,NLAY |
2397 |
|
|
DO I = 1,IRUN |
2398 |
|
|
TEMP2(I,L) = PL(I,L)*ZP001 |
2399 |
|
|
TEMP2(I,L) = SQRT(TEMP2(I,L)) |
2400 |
|
|
ENDDO |
2401 |
|
|
ENDDO |
2402 |
|
|
|
2403 |
|
|
C INVERSE OF MASS IN EACH LAYER |
2404 |
|
|
c ----------------------------- |
2405 |
|
|
DO L = 1,NLAY |
2406 |
|
|
DO I = 1,IRUN |
2407 |
molod |
1.6 |
TEMP3(I,L) = GRAVITY*ZP01 / DP(I,L) |
2408 |
molod |
1.1 |
ENDDO |
2409 |
|
|
ENDDO |
2410 |
|
|
|
2411 |
|
|
C DO LOOP FOR MOISTURE EVAPORATION ABILITY AND CONVEC EVAPORATION. |
2412 |
|
|
c ---------------------------------------------------------------- |
2413 |
|
|
DO 100 L=1,NLAY |
2414 |
|
|
|
2415 |
|
|
DO I = 1,IRUN |
2416 |
|
|
TEMP1(I,3) = TL(I,L) |
2417 |
|
|
TEMP1(I,4) = QL(I,L) |
2418 |
|
|
ENDDO |
2419 |
|
|
|
2420 |
|
|
DO 50 N=1,2 |
2421 |
|
|
IF(N.EQ.1)RELAX=HALF |
2422 |
|
|
IF(N.GT.1)RELAX=ONE |
2423 |
|
|
|
2424 |
|
|
DO I = 1,IRUN |
2425 |
|
|
call qsat ( temp1(i,3),pl(i,L),temp1(i,2),temp1(i,6),.true. ) |
2426 |
|
|
TEMP1(I,5)=TEMP1(I,2)-TEMP1(I,4) |
2427 |
|
|
TEMP1(I,6)=TEMP1(I,6)*ELOCP |
2428 |
|
|
TEMP1(I,5)=TEMP1(I,5)/(ONE+TEMP1(I,6)) |
2429 |
|
|
TEMP1(I,4)=TEMP1(I,4)+TEMP1(I,5)*RELAX |
2430 |
|
|
TEMP1(I,3)=TEMP1(I,3)-TEMP1(I,5)*ELOCP*RELAX |
2431 |
|
|
ENDDO |
2432 |
|
|
50 CONTINUE |
2433 |
|
|
|
2434 |
|
|
DO I = 1,IRUN |
2435 |
|
|
EVP9(I,L) = (TEMP1(I,4) - QL(I,L))/TEMP3(I,L) |
2436 |
|
|
C convective detrained water |
2437 |
|
|
cldlz(i,L) = rain(i,L)*temp3(i,L) |
2438 |
|
|
if( tl(i,L).gt.tice-20.) then |
2439 |
|
|
water(i) = water(i) + rain(i,L) |
2440 |
|
|
else |
2441 |
|
|
crystal(i) = crystal(i) + rain(i,L) |
2442 |
|
|
endif |
2443 |
|
|
ENDDO |
2444 |
|
|
|
2445 |
|
|
C********************************************************************** |
2446 |
|
|
C FOR CONVECTIVE PRECIP, FIND THE "EVAPORATION EFFICIENCY" USING * |
2447 |
|
|
C KESSLERS PARAMETERIZATION * |
2448 |
|
|
C********************************************************************** |
2449 |
|
|
|
2450 |
|
|
DO 20 I=1,IRUN |
2451 |
|
|
|
2452 |
|
|
iceevap(i) = 0. |
2453 |
|
|
watevap(i) = 0. |
2454 |
|
|
|
2455 |
|
|
if( (evp9(i,L).gt.0.) .and. (crystal(i).gt.0.) ) then |
2456 |
|
|
iceevap(I) = EVP9(I,L)*fracice |
2457 |
|
|
IF(iceevap(i).GE.crystal(i)) iceevap(i) = crystal(i) |
2458 |
|
|
EVP9(I,L)=EVP9(I,L)-iceevap(I) |
2459 |
|
|
crystal(i) = crystal(i) - iceevap(i) |
2460 |
|
|
endif |
2461 |
|
|
|
2462 |
|
|
C and now warm precipitate |
2463 |
|
|
if( (evp9(i,L).gt.0.) .and. (water(i).gt.0.) ) then |
2464 |
|
|
exparg = ZM1P04*tmscl*((water(i)*RPHF*TEMP2(I,L))**ZP578) |
2465 |
|
|
AREARAT = ONE-(EXP(EXPARG)) |
2466 |
|
|
watevap(I) = EVP9(I,L)*AREARAT*fracwat |
2467 |
|
|
IF(watevap(I).GE.water(i)) watevap(I) = water(i) |
2468 |
|
|
EVP9(I,L)=EVP9(I,L)-watevap(I) |
2469 |
|
|
water(i) = water(i) - watevap(i) |
2470 |
|
|
endif |
2471 |
|
|
|
2472 |
|
|
QL(I,L) = QL(I,L)+(iceevap(i)+watevap(i))*TEMP3(I,L) |
2473 |
|
|
TL(I,L) = TL(I,L)-(iceevap(i)+watevap(i))*TEMP3(I,L)*ELOCP |
2474 |
|
|
|
2475 |
|
|
20 CONTINUE |
2476 |
|
|
|
2477 |
|
|
100 CONTINUE |
2478 |
|
|
|
2479 |
|
|
do i = 1,irun |
2480 |
|
|
rcon(i) = water(i) + crystal(i) |
2481 |
|
|
enddo |
2482 |
|
|
|
2483 |
|
|
C********************************************************************** |
2484 |
|
|
C Large Scale Precip |
2485 |
|
|
C********************************************************************** |
2486 |
|
|
|
2487 |
|
|
DO 200 L=1,NLAY |
2488 |
|
|
DO I = 1,IRUN |
2489 |
|
|
rainmax(i) = rhcrit(i,L)*evp9(i,L) + |
2490 |
|
|
. ql(i,L)*(rhcrit(i,L)-1.)/temp3(i,L) |
2491 |
|
|
|
2492 |
|
|
if (rainmax(i).LE.0.0) then |
2493 |
|
|
call qsat( tl(i,L),pl(i,L),rh,dum,.false.) |
2494 |
|
|
rh = ql(i,L)/rh |
2495 |
|
|
|
2496 |
|
|
if( rhcrit(i,L).eq.1.0 ) then |
2497 |
|
|
fact = 1.0 |
2498 |
|
|
else |
2499 |
|
|
fact = min( 1.0, alpha + (1.0-alpha)*( rh-rhcrit(i,L)) / |
2500 |
|
|
1 (1.0-rhcrit(i,L)) ) |
2501 |
|
|
endif |
2502 |
|
|
|
2503 |
|
|
C Do not allow clouds above 10 mb |
2504 |
|
|
if( pl(i,L).ge.10.0 ) CLSBTH(I,L) = fact |
2505 |
|
|
RLAR(I) = RLAR(I)-rainmax(I) |
2506 |
|
|
QL(I,L) = QL(I,L)+rainmax(I)*TEMP3(I,L) |
2507 |
|
|
TL(I,L) = TL(I,L)-rainmax(I)*TEMP3(I,L)*ELOCP |
2508 |
|
|
C Large-scale water |
2509 |
|
|
cldlz(i,L) = cldlz(i,L) - rainmax(i)*temp3(i,L) |
2510 |
|
|
ENDIF |
2511 |
|
|
ENDDO |
2512 |
|
|
|
2513 |
|
|
DO I=1,IRUN |
2514 |
|
|
IF((RLAR(I).GT.0.0).AND.(rainmax(I).GT.0.0))THEN |
2515 |
|
|
RPOW=(RLAR(I)*RPHF*TEMP2(I,L))**ZP578 |
2516 |
|
|
EXPARG = ZM1P04*tmscl*RPOW |
2517 |
|
|
AREARAT = ONE-(EXP(EXPARG)) |
2518 |
|
|
TEMP1(I,7) = rainmax(I)*AREARAT |
2519 |
|
|
IF(TEMP1(I,7).GE.RLAR(I)) TEMP1(I,7) = RLAR(I) |
2520 |
|
|
RLAR(I) = RLAR(I)-TEMP1(I,7) |
2521 |
|
|
QL(I,L) = QL(I,L)+TEMP1(I,7)*TEMP3(I,L) |
2522 |
|
|
TL(I,L) = TL(I,L)-TEMP1(I,7)*TEMP3(I,L)*ELOCP |
2523 |
|
|
ENDIF |
2524 |
|
|
ENDDO |
2525 |
|
|
|
2526 |
|
|
200 CONTINUE |
2527 |
|
|
|
2528 |
|
|
RETURN |
2529 |
|
|
END |
2530 |
|
|
subroutine srclouds (th,q,plk,pl,plke,cloud,cldwat,irun,irise, |
2531 |
|
|
1 rhc,offset,alpha) |
2532 |
|
|
C*********************************************************************** |
2533 |
|
|
C |
2534 |
|
|
C PURPOSE: |
2535 |
|
|
C ======== |
2536 |
|
|
C Compute non-precipitating cloud fractions |
2537 |
|
|
C based on Slingo and Ritter (1985). |
2538 |
|
|
C Remove cloudiness where conditionally unstable. |
2539 |
|
|
C |
2540 |
|
|
C INPUT: |
2541 |
|
|
C ====== |
2542 |
|
|
C th ......... Potential Temperature (irun,irise) |
2543 |
|
|
C q .......... Specific Humidity (irun,irise) |
2544 |
|
|
C plk ........ P**Kappa at mid-layer (irun,irise) |
2545 |
|
|
C pl ......... Pressure at mid-layer (irun,irise) |
2546 |
|
|
C plke ....... P**Kappa at edge (irun,irise+1) |
2547 |
|
|
C irun ....... Horizontal dimension |
2548 |
|
|
C irise ...... Vertical dimension |
2549 |
|
|
C |
2550 |
|
|
C OUTPUT: |
2551 |
|
|
C ======= |
2552 |
|
|
C cloud ...... Cloud Fraction (irun,irise) |
2553 |
|
|
C |
2554 |
|
|
C*********************************************************************** |
2555 |
|
|
|
2556 |
|
|
implicit none |
2557 |
|
|
integer irun,irise |
2558 |
|
|
|
2559 |
molod |
1.13 |
_RL th(irun,irise) |
2560 |
|
|
_RL q(irun,irise) |
2561 |
|
|
_RL plk(irun,irise) |
2562 |
|
|
_RL pl(irun,irise) |
2563 |
|
|
_RL plke(irun,irise+1) |
2564 |
|
|
|
2565 |
|
|
_RL cloud(irun,irise) |
2566 |
|
|
_RL cldwat(irun,irise) |
2567 |
|
|
_RL qs(irun,irise) |
2568 |
|
|
|
2569 |
|
|
_RL cp, alhl, getcon, akap |
2570 |
|
|
_RL ratio, temp, elocp |
2571 |
|
|
_RL rhcrit,rh,dum |
2572 |
molod |
1.9 |
integer i,L |
2573 |
molod |
1.1 |
|
2574 |
molod |
1.13 |
_RL rhc(irun,irise) |
2575 |
|
|
_RL offset,alpha |
2576 |
molod |
1.1 |
|
2577 |
|
|
c Explicit Inline Directives |
2578 |
|
|
c -------------------------- |
2579 |
molod |
1.10 |
#ifdef CRAY |
2580 |
|
|
#ifdef f77 |
2581 |
molod |
1.1 |
cfpp$ expand (qsat) |
2582 |
|
|
#endif |
2583 |
|
|
#endif |
2584 |
|
|
|
2585 |
|
|
cp = getcon('CP') |
2586 |
|
|
alhl = getcon('LATENT HEAT COND') |
2587 |
|
|
elocp = alhl/cp |
2588 |
|
|
akap = getcon('KAPPA') |
2589 |
|
|
|
2590 |
|
|
do L = 1,irise |
2591 |
|
|
do i = 1,irun |
2592 |
|
|
temp = th(i,L)*plk(i,L) |
2593 |
|
|
call qsat ( temp,pl(i,L),qs(i,L),dum,.false. ) |
2594 |
|
|
enddo |
2595 |
|
|
enddo |
2596 |
|
|
|
2597 |
|
|
do L = 2,irise |
2598 |
|
|
do i = 1,irun |
2599 |
|
|
rh = q(i,L)/qs(i,L) |
2600 |
|
|
|
2601 |
|
|
rhcrit = rhc(i,L) - offset |
2602 |
|
|
ratio = alpha*(rh-rhcrit)/offset |
2603 |
|
|
|
2604 |
|
|
if(cloud(i,L).eq. 0.0 .and. ratio.gt.0.0 ) then |
2605 |
|
|
cloud(i,L) = min( ratio,1.0 ) |
2606 |
|
|
endif |
2607 |
|
|
|
2608 |
|
|
enddo |
2609 |
|
|
enddo |
2610 |
|
|
|
2611 |
|
|
c Reduce clouds from conditionally unstable layer |
2612 |
|
|
c ----------------------------------------------- |
2613 |
|
|
call ctei ( th,q,cloud,cldwat,pl,plk,plke,irun,irise ) |
2614 |
|
|
|
2615 |
|
|
return |
2616 |
|
|
end |
2617 |
|
|
|
2618 |
|
|
subroutine ctei ( th,q,cldfrc,cldwat,pl,plk,plke,im,lm ) |
2619 |
|
|
implicit none |
2620 |
|
|
integer im,lm |
2621 |
molod |
1.13 |
_RL th(im,lm),q(im,lm),plke(im,lm+1),cldwat(im,lm) |
2622 |
|
|
_RL plk(im,lm),pl(im,lm),cldfrc(im,lm) |
2623 |
molod |
1.1 |
integer i,L |
2624 |
molod |
1.13 |
_RL getcon,cp,alhl,elocp,cpoel,t,p,s,qs,dqsdt,dq |
2625 |
|
|
_RL k,krd,kmm,f |
2626 |
molod |
1.1 |
|
2627 |
|
|
cp = getcon('CP') |
2628 |
|
|
alhl = getcon('LATENT HEAT COND') |
2629 |
|
|
cpoel = cp/alhl |
2630 |
|
|
elocp = alhl/cp |
2631 |
|
|
|
2632 |
|
|
do L=lm,2,-1 |
2633 |
|
|
do i=1,im |
2634 |
|
|
dq = q(i,L)+cldwat(i,L)-q(i,L-1)-cldwat(i,L-1) |
2635 |
|
|
if( dq.eq.0.0 ) dq = 1.0e-20 |
2636 |
|
|
k = 1.0 + cpoel*plke(i,L)*( th(i,L)-th(i,L-1) ) / dq |
2637 |
|
|
|
2638 |
|
|
t = th(i,L)*plk(i,L) |
2639 |
|
|
p = pl(i,L) |
2640 |
|
|
call qsat ( t,p,qs,dqsdt,.true. ) |
2641 |
|
|
|
2642 |
|
|
krd = ( cpoel*t*(1+elocp*dqsdt) )/( 1 + 1.608*dqsdt*t ) |
2643 |
|
|
|
2644 |
|
|
kmm = ( 1+elocp*dqsdt )*( 1 + 0.392*cpoel*t ) |
2645 |
|
|
. / ( 2+(1+1.608*cpoel*t)*elocp*dqsdt ) |
2646 |
|
|
|
2647 |
|
|
s = ( (k-krd)/(kmm-krd) ) |
2648 |
|
|
f = 1.0 - min( 1.0, max(0.0,1.0-exp(-s)) ) |
2649 |
|
|
|
2650 |
|
|
cldfrc(i,L) = cldfrc(i,L)*f |
2651 |
|
|
cldwat(i,L) = cldwat(i,L)*f |
2652 |
|
|
|
2653 |
|
|
enddo |
2654 |
|
|
enddo |
2655 |
|
|
|
2656 |
|
|
return |
2657 |
|
|
end |
2658 |
|
|
|
2659 |
|
|
subroutine back2grd(gathered,indeces,scattered,irun) |
2660 |
|
|
implicit none |
2661 |
|
|
integer i,irun,indeces(irun) |
2662 |
molod |
1.13 |
_RL gathered(irun),scattered(irun) |
2663 |
|
|
_RL temp(irun) |
2664 |
molod |
1.1 |
do i = 1,irun |
2665 |
|
|
temp(indeces(i)) = gathered(i) |
2666 |
|
|
enddo |
2667 |
|
|
do i = 1,irun |
2668 |
|
|
scattered(i) = temp(i) |
2669 |
|
|
enddo |
2670 |
|
|
return |
2671 |
|
|
end |