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