UPP (develop)
Loading...
Searching...
No Matches
FDLVL.f
Go to the documentation of this file.
1
47!--------------------------------------------------------------------------
57!--------------------------------------------------------------------------
58
59 SUBROUTINE fdlvl(ITYPE,TFD,QFD,UFD,VFD,PFD,ICINGFD)
60
61!
62!
63 use vrbls3d, only: zmid, t, q, pmid, icing_gfip, uh, vh
64 use vrbls2d, only: fis
65 use masks, only: lmh
66 use params_mod, only: gi, g
67 use ctlblk_mod, only: jsta, jend, spval, jsta_2l, jend_2u, lm, jsta_m, &
68 jend_m, htfd, nfd, im, jm, nbin_du, &
69 modelname, ista, iend, ista_2l, iend_2u, ista_m, iend_m
70 use gridspec_mod, only: gridtype
71!- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
72 implicit none
73!
74! SET NUMBER OF FD LEVELS.
75!jw integer,intent(in) :: NFD ! coming from calling subroutine
76!
77! DECLARE VARIABLES
78!
79 integer,intent(in) :: ITYPE(NFD)
80!jw real,intent(in) :: HTFD(NFD)
81 real,dimension(ISTA:IEND,JSTA:JEND,NFD),intent(out) :: TFD,QFD,UFD,VFD,PFD,ICINGFD
82!
83 INTEGER LVL(NFD),LHL(NFD)
84 INTEGER IVE(JM),IVW(JM)
85 REAL DZABV(NFD), DZABH(NFD)
86 LOGICAL DONEH, DONEV
87!jw
88 integer I,J,JVS,JVN,IE,IW,JN,JS,JNT,L,LLMH,IFD,N
89 integer ISTART,ISTOP,JSTART,JSTOP
90 real htt,htsfc,httuv,dz,rdz,delt,delq,delu,delv,z1,z2,htabv,htabh,htsfcv
91!
92! SET FD LEVEL HEIGHTS IN METERS.
93! DATA HTFD / 30.E0,50.E0,80.E0,100.E0,305.E0,457.E0,610.E0,914.E0,1524.E0, &
94! 1829.E0,2134.E0,2743.E0,3658.E0,4572.E0,6000.E0/
95!
96!****************************************************************
97! START FDLVL HERE
98!
99! INITIALIZE ARRAYS.
100!
101!$omp parallel do
102 DO ifd = 1,nfd
103 DO j=jsta,jend
104 DO i=ista,iend
105 tfd(i,j,ifd) = spval
106 qfd(i,j,ifd) = spval
107 ufd(i,j,ifd) = spval
108 vfd(i,j,ifd) = spval
109 pfd(i,j,ifd) = spval
110 icingfd(i,j,ifd) = spval
111 ENDDO
112 ENDDO
113 ENDDO
114
115 IF(gridtype == 'E') THEN
116 jvn = 1
117 jvs = -1
118 do j=jsta,jend
119 ive(j) = mod(j,2)
120 ivw(j) = ive(j)-1
121 enddo
122 END IF
123
124 IF(gridtype /= 'A')THEN
125 CALL exch(fis(ista_2l:iend_2u,jsta_2l:jend_2u))
126 DO l=1,lm
127 CALL exch(zmid(ista_2l:iend_2u,jsta_2l:jend_2u,l))
128 END DO
129 istart = ista_m
130 istop = iend_m
131 jstart = jsta_m
132 jstop = jend_m
133 ELSE
134 istart = ista
135 istop = iend
136 jstart = jsta
137 jstop = jend
138 END IF
139 DO ifd = 1, nfd
140!
141! MSL FD LEVELS
142!
143 IF (itype(ifd)==1) THEN
144! write(6,*) 'computing above MSL'
145!
146! LOOP OVER HORIZONTAL GRID.
147!
148 DO j=jstart,jstop
149 DO i=istart,istop
150 htsfc = fis(i,j)*gi
151 llmh = nint(lmh(i,j))
152! IFD = 1
153!
154! LOCATE VERTICAL INDICES OF T,Q,U,V, LEVEL JUST
155! ABOVE EACH FD LEVEL.
156!
157! DO 22 IFD = 1, NFD
158 doneh=.false.
159 donev=.false.
160 DO l = lm,1,-1
161 htt = zmid(i,j,l)
162 IF(gridtype == 'E') THEN
163 ie = i+ive(j)
164 iw = i+ivw(j)
165 jn = j+jvn
166 js = j+jvs
167 httuv = 0.25*(zmid(iw,j,l) &
168 + zmid(ie,j,l)+zmid(i,jn,l)+zmid(i,js,l))
169 ELSE IF(gridtype=='B')THEN
170 ie = i+1
171 iw = i
172 jn = j+1
173 js = j
174 httuv = 0.25*(zmid(iw,j,l) &
175 + zmid(ie,j,l)+zmid(i,jn,l)+zmid(ie,jn,l))
176 ELSE
177 httuv = htt
178 END IF
179
180 IF (.NOT. doneh .AND. htt>htfd(ifd)) THEN
181 lhl(ifd) = l
182 dzabh(ifd) = htt-htfd(ifd)
183 doneh = .true.
184! THIS SHOULD SET BELOW GROUND VALUES TO SPVAL
185 IF(htsfc > htfd(ifd)) THEN
186!mp
187 lhl(ifd) = lm+1 ! CHUANG: changed to lm+1
188!mp
189 ENDIF
190! THIS SHOULD SET BELOW GROUND VALUES TO SPVAL
191! IFD = IFD + 1
192! IF (IFD>NFD) GOTO 30
193 END IF
194
195 IF (.NOT. donev .AND. httuv>htfd(ifd)) THEN
196 lvl(ifd) = l
197 dzabv(ifd) = httuv-htfd(ifd)
198 donev=.true.
199! THIS SHOULD SET BELOW GROUND VALUES TO SPVAL
200 IF(htsfc>htfd(ifd)) THEN
201!mp
202 lvl(ifd)=lm+1 ! CHUANG: changed to lm+1
203!mp
204 ENDIF
205! THIS SHOULD SET BELOW GROUND VALUES TO SPVAL
206! IFD = IFD + 1
207! IF (IFD>NFD) GOTO 30
208 ENDIF
209
210 IF(doneh .AND. donev) exit
211 enddo ! end of l loop
212! 22 CONTINUE
213!
214! COMPUTE T, Q, U, AND V AT FD LEVELS.
215!
216! DO 40 IFD = 1,NFD
217
218 l = lhl(ifd)
219 IF (l < lm) THEN
220 dz = zmid(i,j,l)-zmid(i,j,l+1)
221 rdz = 1./dz
222 delt = t(i,j,l)-t(i,j,l+1)
223 delq = q(i,j,l)-q(i,j,l+1)
224 tfd(i,j,ifd) = t(i,j,l) - delt*rdz*dzabh(ifd)
225 qfd(i,j,ifd) = q(i,j,l) - delq*rdz*dzabh(ifd)
226 pfd(i,j,ifd) = pmid(i,j,l) - (pmid(i,j,l)-pmid(i,j,l+1))*rdz*dzabh(ifd)
227 icingfd(i,j,ifd) = icing_gfip(i,j,l) - &
228 (icing_gfip(i,j,l)-icing_gfip(i,j,l+1))*rdz*dzabh(ifd)
229 ELSEIF (l == lm) THEN
230 tfd(i,j,ifd) = t(i,j,l)
231 qfd(i,j,ifd) = q(i,j,l)
232 pfd(i,j,ifd) = pmid(i,j,l)
233 icingfd(i,j,ifd) = icing_gfip(i,j,l)
234 ENDIF
235
236 l = lvl(ifd)
237 IF (l < lm) THEN
238 IF(gridtype == 'E')THEN
239 ie = i+ive(j)
240 iw = i+ivw(j)
241 jn = j+jvn
242 js = j+jvs
243 z1 = 0.25*(zmid(iw,j,l) &
244 + zmid(ie,j,l)+zmid(i,jn,l)+zmid(i,js,l))
245 z2 = 0.25*(zmid(iw,j,l+1) &
246 + zmid(ie,j,l+1)+zmid(i,jn,l+1)+zmid(i,js,l+1))
247 dz = z1-z2
248
249 ELSE IF(gridtype=='B')THEN
250 ie =i+1
251 iw = i
252 jn = j+1
253 js = j
254 z1 = 0.25*(zmid(iw,j,l) &
255 + zmid(ie,j,l)+zmid(i,jn,l)+zmid(ie,jn,l))
256 z2 = 0.25*(zmid(iw,j,l+1) &
257 + zmid(ie,j,l+1)+zmid(i,jn,l+1)+zmid(ie,jn,l+1))
258 dz = z1-z2
259 ELSE
260 dz = zmid(i,j,l)-zmid(i,j,l+1)
261 END IF
262 rdz = 1./dz
263 delu = uh(i,j,l) - uh(i,j,l+1)
264 delv = vh(i,j,l) - vh(i,j,l+1)
265 ufd(i,j,ifd) = uh(i,j,l) - delu*rdz*dzabv(ifd)
266 vfd(i,j,ifd) = vh(i,j,l) - delv*rdz*dzabv(ifd)
267 ELSEIF (l==lm) THEN
268 ufd(i,j,ifd)=uh(i,j,l)
269 vfd(i,j,ifd)=vh(i,j,l)
270 ENDIF
271! 40 CONTINUE
272!
273! COMPUTE FD LEVEL T, Q, U, AND V AT NEXT K.
274!
275 enddo ! end of i loop
276 enddo ! end of j loop
277! END OF MSL FD LEVELS
278 ELSE
279! write(6,*) 'computing above AGL'
280!
281! AGL FD LEVELS
282!
283!
284! LOOP OVER HORIZONTAL GRID.
285!
286 DO j=jstart,jstop
287 DO i=istart,istop
288 htsfc = fis(i,j)*gi
289 IF(gridtype == 'E') THEN
290 ie = i+ive(j)
291 iw = i+ivw(j)
292 jn = j+jvn
293 js = j+jvs
294 htsfcv = (fis(iw,j)+fis(ie,j)+fis(i,jn)+fis(i,js))*(0.25/g)
295 ELSE IF(gridtype == 'B')THEN
296 ie = i+1
297 iw = i
298 jn = j+1
299 js = j
300 htsfcv = (fis(iw,j)+fis(ie,j)+fis(i,jn)+fis(ie,jn))*(0.25/g)
301 END IF
302 llmh = nint(lmh(i,j))
303! IFD = 1
304!
305! LOCATE VERTICAL INDICES OF T,U,V, LEVEL JUST
306! ABOVE EACH FD LEVEL.
307!
308! DO 222 IFD = 1, NFD
309 doneh=.false.
310 donev=.false.
311 DO l = llmh,1,-1
312 htabh = zmid(i,j,l)-htsfc
313! if(i==245.and.j==813)print*,'Debug FDL HTABH= ',htabh,zmid(i,j,l),htsfc
314 IF(gridtype=='E')THEN
315 htabv = 0.25*(zmid(iw,j,l) &
316 + zmid(ie,j,l)+zmid(i,jn,l)+zmid(i,js,l))-htsfcv
317 ELSE IF(gridtype=='B')THEN
318 htabv = 0.25*(zmid(iw,j,l) &
319 + zmid(ie,j,l)+zmid(i,jn,l)+zmid(ie,jn,l))-htsfcv
320 ELSE
321 htabv = htabh
322 END IF
323
324 IF (.NOT. doneh .AND. htabh>htfd(ifd)) THEN
325 lhl(ifd) = l
326 dzabh(ifd) = htabh-htfd(ifd)
327 doneh=.true.
328! IFD = IFD + 1
329! IF (IFD>NFD) GOTO 230
330 ENDIF
331
332 IF (.NOT. donev .AND. htabv>htfd(ifd)) THEN
333 lvl(ifd) = l
334 dzabv(ifd) = htabv-htfd(ifd)
335 donev = .true.
336! IFD = IFD + 1
337! IF (IFD>NFD) GOTO 230
338 ENDIF
339 IF(doneh .AND. donev) exit
340 enddo ! end of l loop
341!
342! COMPUTE T, Q, U, AND V AT FD LEVELS.
343!
344! 222 CONTINUE
345!
346! DO 240 IFD = 1,NFD
347 l = lhl(ifd)
348 IF (l<lm) THEN
349 dz = zmid(i,j,l)-zmid(i,j,l+1)
350 rdz = 1./dz
351 delt = t(i,j,l)-t(i,j,l+1)
352 delq = q(i,j,l)-q(i,j,l+1)
353 tfd(i,j,ifd) = t(i,j,l) - delt*rdz*dzabh(ifd)
354 qfd(i,j,ifd) = q(i,j,l) - delq*rdz*dzabh(ifd)
355 pfd(i,j,ifd) = pmid(i,j,l) - (pmid(i,j,l)-pmid(i,j,l+1))*rdz*dzabh(ifd)
356 icingfd(i,j,ifd) = icing_gfip(i,j,l) - &
357 (icing_gfip(i,j,l)-icing_gfip(i,j,l+1))*rdz*dzabh(ifd)
358 ELSE
359 tfd(i,j,ifd) = t(i,j,l)
360 qfd(i,j,ifd) = q(i,j,l)
361 pfd(i,j,ifd) = pmid(i,j,l)
362 icingfd(i,j,ifd) = icing_gfip(i,j,l)
363 ENDIF
364
365 l = lvl(ifd)
366 IF (l < lm) THEN
367 IF(gridtype == 'E')THEN
368 ie = i+ive(j)
369 iw = i+ivw(j)
370 jn = j+jvn
371 js = j+jvs
372 z1 = 0.25*(zmid(iw,j,l) &
373 + zmid(ie,j,l)+zmid(i,jn,l)+zmid(i,js,l))
374 z2 = 0.25*(zmid(iw,j,l+1) &
375 + zmid(ie,j,l+1)+zmid(i,jn,l+1)+zmid(i,js,l+1))
376 dz = z1-z2
377 ELSE IF(gridtype=='B')THEN
378 ie = i+1
379 iw = i
380 jn = j+1
381 js = j
382 z1 = 0.25*(zmid(iw,j,l) &
383 + zmid(ie,j,l)+zmid(i,jn,l)+zmid(ie,jn,l))
384 z2 = 0.25*(zmid(iw,j,l+1) &
385 + zmid(ie,j,l+1)+zmid(i,jn,l+1)+zmid(ie,jn,l+1))
386 dz = z1-z2
387 ELSE
388 dz = zmid(i,j,l)-zmid(i,j,l+1)
389 END IF
390 rdz = 1./dz
391 delu = uh(i,j,l)-uh(i,j,l+1)
392 delv = vh(i,j,l)-vh(i,j,l+1)
393 ufd(i,j,ifd) = uh(i,j,l) - delu*rdz*dzabv(ifd)
394 vfd(i,j,ifd) = vh(i,j,l) - delv*rdz*dzabv(ifd)
395 ELSE
396 ufd(i,j,ifd) = uh(i,j,l)
397 vfd(i,j,ifd) = vh(i,j,l)
398 ENDIF
399! 240 CONTINUE
400!
401! COMPUTE FD LEVEL T, U, AND V AT NEXT K.
402!
403 enddo ! end of i loop
404 enddo ! end of j loop
405! END OF AGL FD LEVELS
406 ENDIF
407 enddo ! end of IFD loop
408
409! safety check to avoid tiny QFD values
410 !krf: need ncar and nmm wrf cores in this check as well?
411 IF(modelname=='RAPR' .OR. modelname=='NCAR' .OR. modelname=='NMM') THEN !
412 DO 420 ifd = 1,nfd
413 DO j=jsta,jend
414 DO i=ista,iend
415 if(qfd(i,j,ifd) < 1.0e-8) qfd(i,j,ifd)=0.0
416 ENDDO
417 ENDDO
418420 CONTINUE
419 endif
420!
421! END OF ROUTINE.
422!
423 RETURN
424 END
425
467 SUBROUTINE fdlvl_uv(ITYPE,NFD,HTFD,UFD,VFD)
468!
469!
470 use vrbls3d, only: zmid, pmid, uh, vh
471 use vrbls2d, only: fis
472 use masks, only: lmh
473 use params_mod, only: gi, g
474 use ctlblk_mod, only: jsta, jend, spval, jsta_2l, jend_2u, lm, jsta_m, &
475 jend_m, im, jm, modelname, &
476 ista, iend, ista_2l, iend_2u, ista_m, iend_m
477 use gridspec_mod, only: gridtype
478!- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
479 implicit none
480!
481! DECLARE VARIABLES
482!
483 integer,intent(in) :: ITYPE(NFD)
484 integer,intent(in) :: NFD ! coming from calling subroutine
485 real,intent(in) :: HTFD(NFD)
486 real,dimension(ISTA_2L:IEND_2U,JSTA_2L:JEND_2U,NFD),intent(out) :: UFD,VFD
487!
488 INTEGER LVL(NFD)
489 INTEGER IVE(JM),IVW(JM)
490 REAL DZABV(NFD)
491!jw
492 integer I,J,JVS,JVN,IE,IW,JN,JS,L,LLMH,IFD,N
493 integer ISTART,ISTOP,JSTART,JSTOP
494 real htt,htsfc,httuv,dz,rdz,delu,delv,z1,z2,htabv,htabh,htsfcv
495!
496!****************************************************************
497! START FDLVL_UV HERE
498!
499! INITIALIZE ARRAYS.
500!
501!$omp parallel do
502 DO ifd = 1,nfd
503 DO j=jsta,jend
504 DO i=ista,iend
505 ufd(i,j,ifd) = spval
506 vfd(i,j,ifd) = spval
507 ENDDO
508 ENDDO
509 ENDDO
510
511 IF(gridtype == 'E') THEN
512 jvn = 1
513 jvs = -1
514 do j=jsta,jend
515 ive(j) = mod(j,2)
516 ivw(j) = ive(j)-1
517 enddo
518 END IF
519
520 IF(gridtype /= 'A')THEN
521 CALL exch(fis(ista_2l:iend_2u,jsta_2l:jend_2u))
522 DO l=1,lm
523 CALL exch(zmid(ista_2l:iend_2u,jsta_2l:jend_2u,l))
524 END DO
525 istart = ista_m
526 istop = iend_m
527 jstart = jsta_m
528 jstop = jend_m
529 ELSE
530 istart = ista
531 istop = iend
532 jstart = jsta
533 jstop = jend
534 END IF
535 DO ifd = 1, nfd
536!
537! MSL FD LEVELS
538!
539 IF (itype(ifd) == 1) THEN
540! write(6,*) 'computing above MSL'
541!
542! LOOP OVER HORIZONTAL GRID.
543!
544 DO j=jstart,jstop
545 DO i=istart,istop
546 htsfc = fis(i,j)*gi
547 llmh = nint(lmh(i,j))
548!
549! LOCATE VERTICAL INDICES OF U,V, LEVEL JUST
550! ABOVE EACH FD LEVEL.
551!
552 DO l = lm,1,-1
553 htt = zmid(i,j,l)
554 IF(gridtype == 'E') THEN
555 ie = i+ive(j)
556 iw = i+ivw(j)
557 jn = j+jvn
558 js = j+jvs
559 httuv = 0.25*(zmid(iw,j,l) &
560 + zmid(ie,j,l)+zmid(i,jn,l)+zmid(i,js,l))
561 ELSE IF(gridtype=='B')THEN
562 ie = i+1
563 iw = i
564 jn = j+1
565 js = j
566 httuv = 0.25*(zmid(iw,j,l) &
567 + zmid(ie,j,l)+zmid(i,jn,l)+zmid(ie,jn,l))
568 ELSE
569 httuv = htt
570 END IF
571
572 IF (httuv > htfd(ifd)) THEN
573 lvl(ifd) = l
574 dzabv(ifd) = httuv-htfd(ifd)
575! THIS SHOULD SET BELOW GROUND VALUES TO SPVAL
576 IF(htsfc > htfd(ifd)) THEN
577!mp
578 lvl(ifd)=lm+1 ! CHUANG: changed to lm+1
579!mp
580 ENDIF
581! THIS SHOULD SET BELOW GROUND VALUES TO SPVAL
582
583 exit
584 ENDIF
585 enddo ! end of l loop
586!
587! COMPUTE U V AT FD LEVELS.
588!
589 l = lvl(ifd)
590 IF (l < lm) THEN
591 IF(gridtype == 'E')THEN
592 ie = i+ive(j)
593 iw = i+ivw(j)
594 jn = j+jvn
595 js = j+jvs
596 z1 = 0.25*(zmid(iw,j,l) &
597 + zmid(ie,j,l)+zmid(i,jn,l)+zmid(i,js,l))
598 z2 = 0.25*(zmid(iw,j,l+1) &
599 + zmid(ie,j,l+1)+zmid(i,jn,l+1)+zmid(i,js,l+1))
600 dz = z1-z2
601
602 ELSE IF(gridtype=='B')THEN
603 ie =i+1
604 iw = i
605 jn = j+1
606 js = j
607 z1 = 0.25*(zmid(iw,j,l) &
608 + zmid(ie,j,l)+zmid(i,jn,l)+zmid(ie,jn,l))
609 z2 = 0.25*(zmid(iw,j,l+1) &
610 + zmid(ie,j,l+1)+zmid(i,jn,l+1)+zmid(ie,jn,l+1))
611 dz = z1-z2
612 ELSE
613 dz = zmid(i,j,l)-zmid(i,j,l+1)
614 END IF
615 rdz = 1./dz
616 delu = uh(i,j,l) - uh(i,j,l+1)
617 delv = vh(i,j,l) - vh(i,j,l+1)
618 ufd(i,j,ifd) = uh(i,j,l) - delu*rdz*dzabv(ifd)
619 vfd(i,j,ifd) = vh(i,j,l) - delv*rdz*dzabv(ifd)
620 ELSEIF (l == lm) THEN
621 ufd(i,j,ifd)=uh(i,j,l)
622 vfd(i,j,ifd)=vh(i,j,l)
623 ELSE ! Underground
624 ufd(i,j,ifd)=uh(i,j,lm)
625 vfd(i,j,ifd)=vh(i,j,lm)
626 ENDIF
627!
628 enddo ! end of i loop
629 enddo ! end of j loop
630! END OF MSL FD LEVELS
631 ELSE
632! write(6,*) 'computing above AGL'
633!
634! AGL FD LEVELS
635!
636!
637! LOOP OVER HORIZONTAL GRID.
638!
639 DO j=jstart,jstop
640 DO i=istart,istop
641 htsfc = fis(i,j)*gi
642 IF(gridtype == 'E') THEN
643 ie = i+ive(j)
644 iw = i+ivw(j)
645 jn = j+jvn
646 js = j+jvs
647 htsfcv = (fis(iw,j)+fis(ie,j)+fis(i,jn)+fis(i,js))*(0.25/g)
648 ELSE IF(gridtype == 'B')THEN
649 ie = i+1
650 iw = i
651 jn = j+1
652 js = j
653 htsfcv = (fis(iw,j)+fis(ie,j)+fis(i,jn)+fis(ie,jn))*(0.25/g)
654 END IF
655 llmh = nint(lmh(i,j))
656!
657! LOCATE VERTICAL INDICES OF U,V, LEVEL JUST
658! ABOVE EACH FD LEVEL.
659!
660 DO l = llmh,1,-1
661 htabh = zmid(i,j,l)-htsfc
662 IF(gridtype=='E')THEN
663 htabv = 0.25*(zmid(iw,j,l) &
664 + zmid(ie,j,l)+zmid(i,jn,l)+zmid(i,js,l))-htsfcv
665 ELSE IF(gridtype=='B')THEN
666 htabv = 0.25*(zmid(iw,j,l) &
667 + zmid(ie,j,l)+zmid(i,jn,l)+zmid(ie,jn,l))-htsfcv
668 ELSE
669 htabv = htabh
670 END IF
671
672 IF (htabv > htfd(ifd)) THEN
673 lvl(ifd) = l
674 dzabv(ifd) = htabv-htfd(ifd)
675! IFD = IFD + 1
676 exit
677 ENDIF
678 enddo ! end of l loop
679!
680! COMPUTE U V AT FD LEVELS.
681!
682 l = lvl(ifd)
683 IF (l < lm) THEN
684 IF(gridtype == 'E')THEN
685 ie = i+ive(j)
686 iw = i+ivw(j)
687 jn = j+jvn
688 js = j+jvs
689 z1 = 0.25*(zmid(iw,j,l) &
690 + zmid(ie,j,l)+zmid(i,jn,l)+zmid(i,js,l))
691 z2 = 0.25*(zmid(iw,j,l+1) &
692 + zmid(ie,j,l+1)+zmid(i,jn,l+1)+zmid(i,js,l+1))
693 dz = z1-z2
694 ELSE IF(gridtype=='B')THEN
695 ie = i+1
696 iw = i
697 jn = j+1
698 js = j
699 z1 = 0.25*(zmid(iw,j,l) &
700 + zmid(ie,j,l)+zmid(i,jn,l)+zmid(ie,jn,l))
701 z2 = 0.25*(zmid(iw,j,l+1) &
702 + zmid(ie,j,l+1)+zmid(i,jn,l+1)+zmid(ie,jn,l+1))
703 dz = z1-z2
704 ELSE
705 dz = zmid(i,j,l)-zmid(i,j,l+1)
706 END IF
707 rdz = 1./dz
708 delu = uh(i,j,l)-uh(i,j,l+1)
709 delv = vh(i,j,l)-vh(i,j,l+1)
710 ufd(i,j,ifd) = uh(i,j,l) - delu*rdz*dzabv(ifd)
711 vfd(i,j,ifd) = vh(i,j,l) - delv*rdz*dzabv(ifd)
712 ELSE
713 ufd(i,j,ifd) = uh(i,j,l)
714 vfd(i,j,ifd) = vh(i,j,l)
715 ENDIF
716!
717! COMPUTE FD LEVEL T, U, AND V AT NEXT K.
718!
719 enddo ! end of i loop
720 enddo ! end of j loop
721! END OF AGL FD LEVELS
722 ENDIF
723 enddo ! end of IFD loop
724
725 RETURN
726 END
727
798 SUBROUTINE fdlvl_mass(ITYPE,NFD,PTFD,HTFD,NIN,QIN,QTYPE,QFD)
799 use vrbls3d, only: t,q,zmid,pmid,pint,zint
800 use vrbls2d, only: fis
801 use masks, only: lmh
802 use params_mod, only: gi, g, gamma,pq0, a2, a3, a4, rhmin,rgamog
803 use ctlblk_mod, only: jsta, jend, spval, jsta_2l, jend_2u, lm, jsta_m, &
804 jend_m, im, jm,global,modelname, &
805 ista, iend, ista_2l, iend_2u, ista_m, iend_m
806 use gridspec_mod, only: gridtype
807 use physcons_post,only: con_fvirt, con_rog, con_eps, con_epsm1
808 use upp_physics, only: fpvsnew
809!- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
810 implicit none
811!
812! SET NUMBER OF FD LEVELS.
813!
814! DECLARE VARIABLES
815!
816 real,parameter:: zshul=75.,tvshul=290.66
817
818 integer,intent(in) :: ITYPE(NFD)
819 integer,intent(in) :: NFD ! coming from calling subroutine
820 real, intent(in) :: PTFD(NFD)
821 real,intent(in) :: HTFD(NFD)
822 integer,intent(in) :: NIN
823 real,intent(in) :: QIN(ISTA:IEND,JSTA:JEND,LM,NIN)
824 character, intent(in) :: QTYPE(NIN)
825 real,intent(out) :: QFD(ISTA:IEND,JSTA:JEND,NFD,NIN)
826
827!
828 INTEGER LHL(NFD)
829 REAL DZABH(NFD)
830!jw
831 integer I,J,L,LLMH,IFD,N
832 integer ISTART,ISTOP,JSTART,JSTOP
833 real htt,htsfc,dz,rdz,delq,htabh
834
835 real :: tvu,tvd,gammas,part,ES,QSAT,RHL,PL,ZL,TL,QL
836 real :: TVRL,TVRBLO,TBLO,QBLO
837!
838!****************************************************************
839! START FDLVL_MASS HERE
840!
841! INITIALIZE ARRAYS.
842!
843!$omp parallel do
844 DO n=1,nin
845 DO ifd = 1,nfd
846 DO j=jsta,jend
847 DO i=ista,iend
848 qfd(i,j,ifd,n) = spval
849 ENDDO
850 ENDDO
851 ENDDO
852 ENDDO
853
854 IF(gridtype /= 'A')THEN
855 istart = ista_m
856 istop = iend_m
857 jstart = jsta_m
858 jstop = jend_m
859 ELSE
860 istart = ista
861 istop = iend
862 jstart = jsta
863 jstop = jend
864 END IF
865
866 DO ifd = 1, nfd
867
868!
869! MSL FD LEVELS
870!
871 IF (itype(ifd) == 1) THEN
872! write(6,*) 'computing above MSL'
873!
874! LOOP OVER HORIZONTAL GRID.
875!
876 DO j=jstart,jstop
877 DO i=istart,istop
878 htsfc = fis(i,j)*gi
879 llmh = nint(lmh(i,j))
880!
881! LOCATE VERTICAL INDICES OF Q, LEVEL JUST
882! ABOVE EACH FD LEVEL.
883!
884 DO l = lm,1,-1
885 htt = zmid(i,j,l)
886
887 IF (htt > htfd(ifd)) THEN
888 lhl(ifd) = l
889 dzabh(ifd) = htt-htfd(ifd)
890! THIS SHOULD SET BELOW GROUND VALUES TO SPVAL
891 IF(htsfc > htfd(ifd)) THEN
892!mp
893 lhl(ifd) = lm+1 ! CHUANG: changed to lm+1
894!mp
895 ENDIF
896! THIS SHOULD SET BELOW GROUND VALUES TO SPVAL
897
898 exit
899 END IF
900
901 ENDDO ! end of L loop
902!
903! COMPUTE Q AT FD LEVELS.
904!
905 l = lhl(ifd)
906 IF (l < lm) THEN
907 dz = zmid(i,j,l)-zmid(i,j,l+1)
908 rdz = 1./dz
909 DO n = 1, nin
910 if(qin(i,j,l,n)<spval) then
911 qfd(i,j,ifd,n)=qin(i,j,l+1,n)
912 elseif(qin(i,j,l+1,n)<spval) then
913 qfd(i,j,ifd,n)=qin(i,j,l,n)
914 else
915 qfd(i,j,ifd,n) = qin(i,j,l,n) - &
916 (qin(i,j,l,n)-qin(i,j,l+1,n))*rdz*dzabh(ifd)
917 endif
918 ENDDO
919 ELSEIF (l == lm) THEN
920 DO n = 1, nin
921 qfd(i,j,ifd,n) = qin(i,j,l,n)
922 ENDDO
923 ELSE ! Underground
924 DO n = 1, nin
925 ! Deduce T and Q differently by different models
926 IF(modelname == 'GFS')THEN ! GFS deduce T using Shuell
927 if(qtype(n) == "T" .or. qtype(n) == "Q") then
928 tvu = t(i,j,lm) * (1.+con_fvirt*q(i,j,lm))
929 if(zmid(i,j,lm) > zshul) then
930 tvd = tvu + gamma*zmid(i,j,lm)
931 if(tvd > tvshul) then
932 if(tvu > tvshul) then
933 tvd = tvshul - 5.e-3*(tvu-tvshul)*(tvu-tvshul)
934 else
935 tvd = tvshul
936 endif
937 endif
938 gammas = (tvu-tvd)/zmid(i,j,lm)
939 else
940 gammas = 0.
941 endif
942 part = con_rog*(log(ptfd(ifd))-log(pmid(i,j,lm)))
943 part = zmid(i,j,lm) - tvu*part/(1.+0.5*gammas*part)
944 part = t(i,j,lm) - gamma*(part-zmid(i,j,lm))
945
946 if(qtype(n) == "T") qfd(i,j,ifd,n) = part
947
948 if(qtype(n) == "Q") then
949
950! Compute RH at lowest model layer because Iredell and Chuang decided to compute
951! underground GFS Q to maintain RH
952 es = min(fpvsnew(t(i,j,lm)), pmid(i,j,lm))
953 qsat = con_eps*es/(pmid(i,j,lm)+con_epsm1*es)
954 rhl = q(i,j,lm)/qsat
955! compute saturation water vapor at isobaric level
956 es = min(fpvsnew(part), ptfd(ifd))
957 qsat = con_eps*es/(ptfd(ifd)+con_epsm1*es)
958! Q at isobaric level is computed by maintaining constant RH
959 qfd(i,j,ifd,n) = rhl*qsat
960 endif
961 endif
962
963 ELSE
964 if(qtype(n) == "T" .or. qtype(n) == "Q") then
965 pl = pint(i,j,lm-1)
966 zl = zint(i,j,lm-1)
967 tl = 0.5*(t(i,j,lm-2)+t(i,j,lm-1))
968 ql = 0.5*(q(i,j,lm-2)+q(i,j,lm-1))
969
970 qsat = pq0/pl*exp(a2*(tl-a3)/(tl-a4))
971 rhl = ql/qsat
972!
973 IF(rhl > 1.)THEN
974 rhl = 1.
975 ql = rhl*qsat
976 ENDIF
977!
978 IF(rhl < rhmin)THEN
979 rhl = rhmin
980 ql = rhl*qsat
981 ENDIF
982!
983 tvrl = tl*(1.+0.608*ql)
984 tvrblo = tvrl*(ptfd(ifd)/pl)**rgamog
985 tblo = tvrblo/(1.+0.608*ql)
986
987 qsat = pq0/ptfd(ifd)*exp(a2*(tblo-a3)/(tblo-a4))
988 if(qtype(n) == "T") qfd(i,j,ifd,n) = tblo
989 qblo = rhl*qsat
990 if(qtype(n) == "Q") qfd(i,j,ifd,n) = max(1.e-12,qblo)
991 endif
992 END IF ! endif loop for deducing T and Q differently for GFS
993
994 if(qtype(n) == "K") qfd(i,j,ifd,n)= max(0.0,0.5*(qin(i,j,lm,n)+qin(i,j,lm-1,n))) ! TKE
995 END DO
996
997 ENDIF ! Underground
998
999!
1000! COMPUTE FD LEVEL Q AT NEXT K.
1001!
1002 enddo ! end of i loop
1003 enddo ! end of j loop
1004! END OF MSL FD LEVELS
1005 ELSE
1006! write(6,*) 'computing above AGL'
1007!
1008! AGL FD LEVELS
1009!
1010!
1011! LOOP OVER HORIZONTAL GRID.
1012!
1013 DO j=jstart,jstop
1014 DO i=istart,istop
1015 htsfc = fis(i,j)*gi
1016 llmh = nint(lmh(i,j))
1017!
1018! LOCATE VERTICAL INDICES OF Q, LEVEL JUST
1019! ABOVE EACH FD LEVEL.
1020!
1021 DO l = llmh,1,-1
1022 htabh = zmid(i,j,l)-htsfc
1023
1024 IF ( htabh > htfd(ifd)) THEN
1025 lhl(ifd) = l
1026 dzabh(ifd) = htabh-htfd(ifd)
1027
1028 exit
1029 ENDIF
1030 enddo ! end of l loop
1031!
1032! COMPUTE Q AT FD LEVELS.
1033!
1034 l = lhl(ifd)
1035 IF (l < lm) THEN
1036 dz = zmid(i,j,l)-zmid(i,j,l+1)
1037 rdz = 1./dz
1038 DO n = 1, nin
1039 if(qin(i,j,l,n)<spval) then
1040 qfd(i,j,ifd,n)=qin(i,j,l+1,n)
1041 elseif(qin(i,j,l+1,n)<spval) then
1042 qfd(i,j,ifd,n)=qin(i,j,l,n)
1043 else
1044 qfd(i,j,ifd,n) = qin(i,j,l,n) - &
1045 (qin(i,j,l,n)-qin(i,j,l+1,n))*rdz*dzabh(ifd)
1046 endif
1047 ENDDO
1048 ELSE
1049 DO n = 1, nin
1050 qfd(i,j,ifd,n) = qin(i,j,l,n)
1051 ENDDO
1052 ENDIF
1053
1054!
1055! COMPUTE FD LEVEL Q AT NEXT K.
1056!
1057 enddo ! end of i loop
1058 enddo ! end of j loop
1059! END OF AGL FD LEVELS
1060 ENDIF
1061 enddo ! end of IFD loop
1062
1063!
1064! END OF ROUTINE.
1065!
1066 RETURN
1067 END
subroutine exch(a)
exch() Subroutine that exchanges one halo row.
Definition EXCH.f:27
subroutine fdlvl(itype, tfd, qfd, ufd, vfd, pfd, icingfd)
fdlvl() Subroutine that computes T, Q, U, V, P, ICING on the flight levels (FD).
Definition FDLVL.f:60
subroutine fdlvl_uv(itype, nfd, htfd, ufd, vfd)
Computes FD level for u,v.
Definition FDLVL.f:468
subroutine fdlvl_mass(itype, nfd, ptfd, htfd, nin, qin, qtype, qfd)
Computes FD level for mass variables.
Definition FDLVL.f:799