UPP (upp-srw-2.2.0)
Loading...
Searching...
No Matches
BNDLYR.f
Go to the documentation of this file.
1
3!
48 SUBROUTINE bndlyr(PBND,TBND,QBND,RHBND,UBND,VBND, &
49 WBND,OMGBND,PWTBND,QCNVBND,LVLBND)
50
51!
52!
53 use vrbls3d, only: pint, q, uh, vh, pmid, t, omga, wh, cwm
54 use masks, only: lmh
55 use params_mod, only: d00, gi, pq0, a2, a3, a4
56 use ctlblk_mod, only: jsta_2l, jend_2u, lm, jsta, jend, modelname, &
57 jsta_m, jend_m, im, nbnd, spval, ista_2l, iend_2u, ista_m, iend_m, ista, iend
58 use physcons_post, only: con_rd, con_rv, con_eps, con_epsm1
59 use gridspec_mod, only: gridtype
60 use upp_physics, only: fpvsnew
61!
62 implicit none
63!
64! DECLARE VARIABLES.
65!
66 real,PARAMETER :: DPBND=30.e2
67 integer, dimension(ista:iend,jsta:jend,NBND),intent(inout) :: LVLBND
68 real, dimension(ista:iend,jsta:jend,NBND),intent(inout) :: PBND,TBND, &
69 qbnd,rhbnd,ubnd,vbnd,wbnd,omgbnd,pwtbnd,qcnvbnd
70
71 REAL Q1D(ISTA_2L:IEND_2U,JSTA_2L:JEND_2U),V1D(ISTA_2L:IEND_2U,JSTA_2L:JEND_2U), &
72 u1d(ista_2l:iend_2u,jsta_2l:jend_2u),qcnv1d(ista_2l:iend_2u,jsta_2l:jend_2u)
73!
74 REAL, ALLOCATABLE :: PBINT(:,:,:),QSBND(:,:,:)
75 REAL, ALLOCATABLE :: PSUM(:,:,:), QCNVG(:,:,:)
76 REAL, ALLOCATABLE :: PVSUM(:,:,:),NSUM(:,:,:)
77!
78 integer I,J,L,IE,IW,LL,LV,LBND
79 real DP,QSAT,PV1,PV2,PMV,RPSUM,RPVSUM,PMIN,PM,DELP,PMINV,DELPV
80 real es
81!
82!*****************************************************************************
83! START BNDLYR HERE
84!
85 ALLOCATE (pbint(ista_2l:iend_2u,jsta_2l:jend_2u,nbnd+1))
86 ALLOCATE (qsbnd(ista_2l:iend_2u,jsta_2l:jend_2u,nbnd))
87 ALLOCATE (psum(ista_2l:iend_2u,jsta_2l:jend_2u,nbnd))
88 ALLOCATE (qcnvg(ista_2l:iend_2u,jsta_2l:jend_2u,lm))
89 ALLOCATE (pvsum(ista_2l:iend_2u,jsta_2l:jend_2u,nbnd))
90 ALLOCATE (nsum(ista_2l:iend_2u,jsta_2l:jend_2u,nbnd))
91!
92! LOOP OVER HORIZONTAL GRID. AT EACH MASS POINT COMPUTE
93! PRESSURE AT THE INTERFACE OF EACH BOUNDARY LAYER.
94!
95!$omp parallel do private(i,j)
96 DO j=jsta,jend
97 DO i=ista,iend
98 pbint(i,j,1) = pint(i,j,nint(lmh(i,j))+1)
99 ENDDO
100 ENDDO
101!
102 DO lbnd=2,nbnd+1
103!$omp parallel do private(i,j)
104 DO j=jsta,jend
105 DO i=ista,iend
106 pbint(i,j,lbnd) = pbint(i,j,lbnd-1) - dpbnd
107 ENDDO
108 ENDDO
109 ENDDO
110
111! COMPUTE MOISTURE CONVERGENCE FOR EVERY LEVEL
112 DO l=1,lm
113!$omp parallel do private(i,j)
114 DO j=jsta_2l,jend_2u
115 DO i=ista_2l,iend_2u
116 q1d(i,j) = q(i,j,l)
117 u1d(i,j) = uh(i,j,l)
118 v1d(i,j) = vh(i,j,l)
119 ENDDO
120 ENDDO
121 CALL calmcvg(q1d,u1d,v1d,qcnv1d)
122!$omp parallel do private(i,j)
123 DO j=jsta,jend
124 DO i=ista,iend
125 qcnvg(i,j,l)=qcnv1d(i,j)
126 ENDDO
127 ENDDO
128 ENDDO
129
130!
131! LOOP OVER HORIZONTAL. AT EACH MASS POINT COMPUTE
132! MASS WEIGHTED LAYER MEAN P, T, Q, U, V, OMEGA,
133! WAND PRECIPITABLE WATER IN EACH BOUNDARY LAYER FROM THE SURFACE UP.
134!
135!!$omp+ private(dp,pm,qsat)
136!!$omp parallel do private(i,j,lbnd,l,ie,iw,dp,pm,qsat,pv1,pv2,pmv)
137 DO lbnd=1,nbnd
138!$omp parallel do private(i,j)
139 DO j=jsta,jend
140 DO i=ista,iend
141 pbnd(i,j,lbnd) = d00
142 tbnd(i,j,lbnd) = d00
143 qbnd(i,j,lbnd) = d00
144 qsbnd(i,j,lbnd) = d00
145 rhbnd(i,j,lbnd) = d00
146 ubnd(i,j,lbnd) = d00
147 vbnd(i,j,lbnd) = d00
148 wbnd(i,j,lbnd) = d00
149 omgbnd(i,j,lbnd) = d00
150 lvlbnd(i,j,lbnd) = 0
151 nsum(i,j,lbnd) = 0
152 psum(i,j,lbnd) = d00
153 pvsum(i,j,lbnd) = d00
154 pwtbnd(i,j,lbnd) = d00
155 qcnvbnd(i,j,lbnd)= d00
156 ENDDO
157 ENDDO
158!
159!!$omp parallel do private(i,j,l,dp,pm,es,qsat)
160 DO l=1,lm
161!$omp parallel do private(i,j,dp,pm,es,qsat)
162 DO j=jsta,jend
163 DO i=ista,iend
164!
165 pm = pmid(i,j,l)
166 IF(pm<spval)THEN
167 IF((pbint(i,j,lbnd) >= pm).AND. &
168 (pbint(i,j,lbnd+1) <= pm)) THEN
169 dp = pint(i,j,l+1) - pint(i,j,l)
170 psum(i,j,lbnd) = psum(i,j,lbnd) + dp
171 nsum(i,j,lbnd) = nsum(i,j,lbnd) + 1
172 lvlbnd(i,j,lbnd) = lvlbnd(i,j,lbnd) + l
173 tbnd(i,j,lbnd) = tbnd(i,j,lbnd) + t(i,j,l)*dp
174 qbnd(i,j,lbnd) = qbnd(i,j,lbnd) + q(i,j,l)*dp
175 omgbnd(i,j,lbnd) = omgbnd(i,j,lbnd) + omga(i,j,l)*dp
176 IF(gridtype == 'A')THEN
177 ubnd(i,j,lbnd) = ubnd(i,j,lbnd) + uh(i,j,l)*dp
178 vbnd(i,j,lbnd) = vbnd(i,j,lbnd) + vh(i,j,l)*dp
179 END IF
180 wbnd(i,j,lbnd) = wbnd(i,j,lbnd) + wh(i,j,l)*dp
181 qcnvbnd(i,j,lbnd) = qcnvbnd(i,j,lbnd) + qcnvg(i,j,l)*dp
182 pwtbnd(i,j,lbnd) = pwtbnd(i,j,lbnd) &
183 + ( q(i,j,l)+cwm(i,j,l))*dp*gi
184 IF(modelname == 'GFS')THEN
185 es = min(fpvsnew(t(i,j,l)),pm)
186 qsat = con_eps*es/(pm+con_epsm1*es)
187 ELSE
188 qsat = pq0/pm*exp(a2*(t(i,j,l)-a3)/(t(i,j,l)-a4))
189 END IF
190 qsbnd(i,j,lbnd) = qsbnd(i,j,lbnd) + qsat*dp
191 ENDIF
192 ELSE !undeined grids
193 pbnd(i,j,lbnd)=spval
194 tbnd(i,j,lbnd)=spval
195 ubnd(i,j,lbnd)=spval
196 vbnd(i,j,lbnd)=spval
197 wbnd(i,j,lbnd)=spval
198 omgbnd(i,j,lbnd)=spval
199 qcnvbnd(i,j,lbnd)=spval
200 pwtbnd(i,j,lbnd)=spval
201 qbnd(i,j,lbnd)=spval
202 qsbnd(i,j,lbnd)=spval
203 rhbnd(i,j,lbnd)=spval
204 ENDIF
205 ENDDO
206 ENDDO
207 ENDDO
208
209
210 IF(gridtype=='E')THEN
211 CALL exch(pint(ista_2l:iend_2u,jsta_2l:jend_2u,1))
212 DO l=1,lm
213 CALL exch(pint(ista_2l:iend_2u,jsta_2l:jend_2u,l+1))
214!$omp parallel do private(i,j,ie,iw,dp,pv1,pv2,pmv)
215 DO j=jsta_m,jend_m
216 DO i=ista_m,iend_m
217 ie = i+mod(j,2)
218 iw = i+mod(j,2)-1
219 pv1 = 0.25*(pint(iw,j,l) + pint(ie,j,l) &
220 +pint(i,j+1,l) + pint(i,j-1,l))
221 pv2 = 0.25*(pint(iw,j,l+1) + pint(ie,j,l+1) &
222 +pint(i,j+1,l+1) + pint(i,j-1,l+1))
223 dp = pv2-pv1
224 pmv = 0.5*(pv1+pv2)
225 IF((pbint(iw,j,lbnd)>=pmv).AND. &
226 (pbint(iw,j,lbnd+1)<=pmv)) THEN
227 pvsum(i,j,lbnd) = pvsum(i,j,lbnd) + dp
228 ubnd(i,j,lbnd) = ubnd(i,j,lbnd) + dp* uh(i,j,l)
229 vbnd(i,j,lbnd) = vbnd(i,j,lbnd) + dp*vh(i,j,l)
230 ENDIF
231!
232 ENDDO
233 ENDDO
234 ENDDO
235 ELSE IF (gridtype=='B')THEN
236 CALL exch(pint(ista_2l:iend_2u,jsta_2l:jend_2u,1))
237 DO l=1,lm
238 CALL exch(pint(ista_2l:iend_2u,jsta_2l:jend_2u,l+1))
239!$omp parallel do private(i,j,ie,iw,dp,pv1,pv2,pmv)
240 DO j=jsta_m,jend_m
241 DO i=ista_m,iend_m
242 ie = i+1
243 iw = i
244 pv1 = 0.25*(pint(iw,j,l) + pint(ie,j,l) &
245 +pint(iw,j+1,l) + pint(ie,j+1,l))
246 pv2 = 0.25*(pint(iw,j,l+1) + pint(ie,j,l+1) &
247 +pint(iw,j+1,l+1) + pint(ie,j+1,l+1))
248 dp = pv2-pv1
249 pmv = 0.5*(pv1+pv2)
250 IF((pbint(iw,j,lbnd)>=pmv).AND. &
251 (pbint(iw,j,lbnd+1)<=pmv)) THEN
252 pvsum(i,j,lbnd) = pvsum(i,j,lbnd)+dp
253 ubnd(i,j,lbnd) = ubnd(i,j,lbnd)+uh(i,j,l)*dp
254 vbnd(i,j,lbnd) = vbnd(i,j,lbnd)+vh(i,j,l)*dp
255 ENDIF
256!
257 ENDDO
258 ENDDO
259 ENDDO
260 END IF
261
262 ENDDO ! end of lbnd loop
263!
264!!$omp+ private(rpsum)
265!$omp parallel do private(i,j,lbnd,rpsum,rpvsum)
266 DO lbnd=1,nbnd
267 DO j=jsta,jend
268 DO i=ista,iend
269 IF(psum(i,j,lbnd)/=0..AND.tbnd(i,j,lbnd)<spval)THEN
270 rpsum = 1./psum(i,j,lbnd)
271 lvlbnd(i,j,lbnd)= lvlbnd(i,j,lbnd)/nsum(i,j,lbnd)
272 pbnd(i,j,lbnd) = (pbint(i,j,lbnd)+pbint(i,j,lbnd+1))*0.5
273 tbnd(i,j,lbnd) = tbnd(i,j,lbnd)*rpsum
274 qbnd(i,j,lbnd) = qbnd(i,j,lbnd)*rpsum
275 qsbnd(i,j,lbnd) = qsbnd(i,j,lbnd)*rpsum
276 omgbnd(i,j,lbnd)= omgbnd(i,j,lbnd)*rpsum
277 IF(gridtype=='A')THEN
278 ubnd(i,j,lbnd) = ubnd(i,j,lbnd)*rpsum
279 vbnd(i,j,lbnd) = vbnd(i,j,lbnd)*rpsum
280 END IF
281 wbnd(i,j,lbnd) = wbnd(i,j,lbnd)*rpsum
282 IF(qcnvbnd(i,j,lbnd)<spval) &
283 qcnvbnd(i,j,lbnd) = qcnvbnd(i,j,lbnd)*rpsum
284 ENDIF
285 ENDDO
286 ENDDO
287
288 IF(gridtype=='E' .or. gridtype=='B')THEN
289 DO j=jsta_m,jend_m
290 DO i=ista_m,iend_m
291 IF(pvsum(i,j,lbnd)/=0.)THEN
292 rpvsum = 1./pvsum(i,j,lbnd)
293 ubnd(i,j,lbnd) = ubnd(i,j,lbnd)*rpvsum
294 vbnd(i,j,lbnd) = vbnd(i,j,lbnd)*rpvsum
295 ENDIF
296 ENDDO
297 ENDDO
298 END IF
299 ENDDO
300!
301! IF NO ETA MID LAYER PRESSURES FELL WITHIN A BND LYR,
302! FIND THE CLOSEST LAYER TO THE BND LYR AND ASSIGN THE VALUES THERE
303!
304!!$omp+ private(delp,dp,l,pm,pmin,qsat)
305!$omp parallel do private(i,j,lbnd,l,ll,ie,iw,pminv,delp,dp,pm,pmin,es, &
306!$omp& qsat,pmv,delpv,lv)
307 DO lbnd=1,nbnd
308 DO j=jsta,jend
309 DO i=ista,iend
310 IF(psum(i,j,lbnd)==0..AND.pbnd(i,j,lbnd)<spval)THEN
311 l = lm
312 pmin = 9999999.
313 pbnd(i,j,lbnd) = (pbint(i,j,lbnd)+pbint(i,j,lbnd+1))*0.5
314!
315 DO ll=1,lm
316 pm = pmid(i,j,ll)
317 delp = abs(pm-pbnd(i,j,lbnd))
318 IF(delp<pmin)THEN
319 pmin = delp
320 l = ll
321 ENDIF
322 ENDDO
323!
324 dp = pint(i,j,l+1)-pint(i,j,l)
325 pm = pmid(i,j,l)
326 lvlbnd(i,j,lbnd) = l
327 tbnd(i,j,lbnd) = t(i,j,l)
328 qbnd(i,j,lbnd) = q(i,j,l)
329 IF(gridtype == 'A')THEN
330 ubnd(i,j,lbnd) = uh(i,j,l)
331 vbnd(i,j,lbnd) = vh(i,j,l)
332 END IF
333 wbnd(i,j,lbnd) = wh(i,j,l)
334 qcnvbnd(i,j,lbnd) = qcnvg(i,j,l)
335 IF(modelname == 'GFS' .OR. modelname == 'FV3R')THEN
336 es = fpvsnew(t(i,j,l))
337 es = min(es,pm)
338 qsat = con_eps*es/(pm+con_epsm1*es)
339 ELSE
340 qsat=pq0/pm*exp(a2*(t(i,j,l)-a3)/(t(i,j,l)-a4))
341 END IF
342 qsbnd(i,j,lbnd) = qsat
343 omgbnd(i,j,lbnd) = omga(i,j,l)
344 pwtbnd(i,j,lbnd) = (q(i,j,l)+cwm(i,j,l))*dp*gi
345 ENDIF
346!
347! RH, BOUNDS CHECK
348!
349 IF(qsbnd(i,j,lbnd)/=0..AND.qbnd(i,j,lbnd)<spval)THEN
350 rhbnd(i,j,lbnd) = qbnd(i,j,lbnd)/qsbnd(i,j,lbnd)
351 IF (rhbnd(i,j,lbnd)>1.0) THEN
352 rhbnd(i,j,lbnd) = 1.0
353 qbnd(i,j,lbnd) = rhbnd(i,j,lbnd)*qsbnd(i,j,lbnd)
354 ENDIF
355 IF (rhbnd(i,j,lbnd)<0.01) THEN
356 rhbnd(i,j,lbnd) = 0.01
357 qbnd(i,j,lbnd) = rhbnd(i,j,lbnd)*qsbnd(i,j,lbnd)
358 ENDIF
359 ENDIF
360 ENDDO
361 ENDDO
362!
363 IF(gridtype == 'E')THEN
364 DO j=jsta_m,jend_m
365 DO i=ista_m,iend_m
366 IF(pvsum(i,j,lbnd)==0.)THEN
367 lv = lm
368 pminv = 9999999.
369 ie = i+mod(j,2)
370 iw = i+mod(j,2)-1
371!
372! PINT HALOS UPDATED ALREADY
373!
374 DO ll=1,lm
375 pmv = 0.125*(pint(iw,j,ll) + pint(ie,j,ll) + &
376 pint(i,j+1,ll) + pint(i,j-1,ll) + &
377 pint(iw,j,ll+1) + pint(ie,j,ll+1) + &
378 pint(i,j+1,ll+1) + pint(i,j-1,ll+1))
379 delpv = abs(pmv-pbnd(i,j,lbnd))
380 IF(delpv<pminv)THEN
381 pminv = delpv
382 lv = ll
383 ENDIF
384 ENDDO
385!
386 ubnd(i,j,lbnd) = uh(i,j,lv)
387 vbnd(i,j,lbnd) = vh(i,j,lv)
388 ENDIF
389 ENDDO
390 ENDDO
391! END IF
392
393 ELSE IF(gridtype=='B')THEN
394 DO j=jsta_m,jend_m
395 DO i=ista_m,iend_m
396 IF(pvsum(i,j,lbnd)==0.)THEN
397 lv=lm
398 pminv=9999999.
399 ie=i+1
400 iw=i
401!
402! PINT HALOS UPDATED ALREADY
403!
404 DO ll=1,lm
405 pmv=0.125*(pint(iw,j,ll)+pint(ie,j,ll)+ &
406 pint(iw,j+1,ll)+pint(ie,j+1,ll)+ &
407 pint(iw,j,ll+1)+pint(ie,j,ll+1)+ &
408 pint(iw,j+1,ll+1)+pint(ie,j+1,ll+1))
409 delpv=abs(pmv-pbnd(i,j,lbnd))
410 IF(delpv<pminv)THEN
411 pminv=delpv
412 lv=ll
413 ENDIF
414 ENDDO
415
416 ubnd(i,j,lbnd) = uh(i,j,lv)
417 vbnd(i,j,lbnd) = vh(i,j,lv)
418 ENDIF
419 ENDDO
420 ENDDO
421 END IF
422 ENDDO ! end of lbnd loop
423!
424 DEALLOCATE (pbint, qsbnd, psum, pvsum, qcnvg, nsum)
425!
426! END OF ROUTINE
427!
428 RETURN
429 END
430!
subroutine bndlyr(pbnd, tbnd, qbnd, rhbnd, ubnd, vbnd, wbnd, omgbnd, pwtbnd, qcnvbnd, lvlbnd)
Computes constant mass mean fields.
Definition BNDLYR.f:50
subroutine calmcvg(q1d, u1d, v1d, qcnvg)
Subroutine that computes moisture convergence.
Definition CALMCVG.f:44