UPP  V11.0.0
 All Data Structures Files Functions Pages
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 !
Definition: MASKS_mod.f:1
Definition: physcons.f:1
elemental real function, public fpvsnew(t)
Definition: UPP_PHYSICS.f:345
calcape() computes CAPE/CINS and other storm related variables.
Definition: UPP_PHYSICS.f:27