UPP (develop)
Loading...
Searching...
No Matches
SLP_new.f
1 SUBROUTINE memslp(TPRES,QPRES,FIPRES)
2!
3!$$$ SUBPROGRAM DOCUMENTATION BLOCK
4! . . .
5! SUBROUTINE: MEMSLP MEMBRANE SLP REDUCTION
6!
7! ABSTRACT: THIS ROUTINE COMPUTES THE SEA LEVEL PRESSURE
8! REDUCTION USING THE MESINGER RELAXATION
9! METHOD FOR SIGMA COORDINATES.
10! A BY-PRODUCT IS THE
11! SET OF VALUES FOR THE UNDERGROUND TEMPERATURES
12! ON THE SPECIFIED PRESSURE LEVELS
13!
14! PROGRAM HISTORY LOG:
15! 99-09-23 T BLACK - REWRITTEN FROM ROUTINE SLP (ETA
16! COORDINATES)
17! 02-07-26 H CHUANG - PARALLIZE AND MODIFIED FOR WRF A/C GRIDS
18! ALSO REDUCE S.O.R. COEFF FROM 1.75 to 1.25
19! BECAUSE THERE WAS NUMERICAL INSTABILITY
20! 02-08-21 H CHUANG - MODIFIED TO ALWAYS USE OLD TTV FOR RELAXATION
21! SO THAT THERE WAS BIT REPRODUCIBILITY BETWEEN
22! USING ONE AND MULTIPLE TASKS
23! 11-04-29 H CHUANG - FIX GFS GIBSING BY USING LM-1 STATE VARIABLES
24! TO DERIVE SLP HYDROSTATICALLY
25! 13-12-06 H CHUANG - REMOVE EXTRA SMOOTHING OF SLP ITSELF
26! CHANGES TO AVOID RELAXATION FOR ABOVE G GIBSING
27! ARE COMMENTED OUT FOR NOW
28! 19-10-30 Bo CUI - REMOVE "GOTO" STATEMENT
29! 21-07-26 W Meng - Restrict computation from undefined grids
30! 21-07-07 J Meng - 2D DECOMPOSITION
31! 21-09-25 W Meng - Further modification for restricting computation
32! from undefined grids.
33!
34! USAGE: CALL SLPSIG FROM SUBROUITNE ETA2P
35!
36! INPUT ARGUMENT LIST:
37! PD - SFC PRESSURE MINUS PTOP
38! FIS - SURFACE GEOPOTENTIAL
39! T - TEMPERATURE
40! Q - SPECIFIC HUMIDITY
41! FI - GEOPOTENTIAL
42! PT - TOP PRESSURE OF DOMAIN
43!
44! OUTPUT ARGUMENT LIST:
45! PSLP - THE FINAL REDUCED SEA LEVEL PRESSURE ARRAY
46!
47! SUBPROGRAMS CALLED:
48! UNIQUE:
49! NONE
50!
51!-----------------------------------------------------------------------
52 use vrbls3d, only: pint, zint, t, q
53 use vrbls2d, only: pslp, fis
54 use masks, only: lmh
55 use params_mod, only: overrc, ad05, cft0, g, rd, d608, h1, kslpd
56 use ctlblk_mod, only: jend, jsta, spval, spl, num_procs, mpi_comm_comp, lsmp1, &
57 jsta_m, jend_m, lm, im, jsta_2l, jend_2u, lsm, jm,&
58 im_jm, iend, ista, ista_m, iend_m, ista_2l, iend_2u
59!- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
60 implicit none
61!
62 include "mpif.h"
63!-----------------------------------------------------------------------
64 integer,PARAMETER :: NFILL=0,nrlx1=500,nrlx2=100
65 real,parameter:: def_of_mountain=2.0
66!-----------------------------------------------------------------------
67 real,dimension(ISTA_2L:IEND_2U,JSTA_2L:JEND_2U,LSM),intent(in) :: QPRES
68 real,dimension(ISTA_2L:IEND_2U,JSTA_2L:JEND_2U,LSM),intent(inout) :: TPRES,FIPRES
69 REAL :: TTV(ISTA_2L:IEND_2U,JSTA_2L:JEND_2U),TNEW(ISTA_2L:IEND_2U,JSTA_2L:JEND_2U) &
70 , p1(ista_2l:iend_2u,jsta_2l:jend_2u),htm2d(ista_2l:iend_2u,jsta_2l:jend_2u)
71 REAL :: HTMO(ISTA_2L:IEND_2U,JSTA_2L:JEND_2U,LSM)
72 real :: P2,TLYR,GZ1,GZ2,SPLL,PSFC,PCHK,SLOPE,TVRTC,DIS,TVRT,tem
73!-----------------------------------------------------------------------
74!-----------------------------------------------------------------------
75 INTEGER :: KMNTM(LSM),IMNT(IM_JM,LSM),JMNT(IM_JM,LSM) &
76 , lmho(ista_2l:iend_2u,jsta_2l:jend_2u)
77 INTEGER :: IHE(JM),IHW(JM),IVE(JM),IVW(JM),IHS(JM),IHN(JM)
78 integer ii,jj,I,J,L,N,LLMH,KM,KS,IHH2,KOUNT,KMN,NRLX,LHMNT, &
79 lmhij,lmap1,kmm,lp,lxxx,ierr
80! dong
81 real a1,a2,a3,a4,a5,a6,a7,a8
82!-----------------------------------------------------------------------
83 LOGICAL :: DONE(ISTA_2L:IEND_2U,JSTA_2L:JEND_2U)
84!-----------------------------------------------------------------------
85!***
86!*** CALCULATE THE I-INDEX EAST-WEST INCREMENTS
87!***
88!
89 ii = (iend-ista)/2
90 jj = (jend-jsta)/2
91 DO j=1,jm
92 ihe(j) = 1
93 ihw(j) = -1
94 ihs(j) = -1
95 ihn(j) = 1
96 ive(j) = mod(j,2)
97 ivw(j) = ive(j)-1
98 ENDDO
99! print*,'relaxation coeff= ',OVERRC
100!-----------------------------------------------------------------------
101!***
102!*** INITIALIZE ARRAYS. LOAD SLP ARRAY WITH SURFACE PRESSURE.
103!***
104!$omp parallel do private(i,j,llmh)
105 DO j=jsta,jend
106 DO i=ista,iend
107 llmh = nint(lmh(i,j))
108 pslp(i,j) = pint(i,j,llmh+1)
109! dong
110! TTV(I,J) = 0.
111 ttv(i,j) = spval
112 tnew(i,j) = spval
113
114
115 lmho(i,j) = lsm
116 done(i,j) = .false.
117 ENDDO
118 ENDDO
119!
120!--------------------------------------------------------------------
121!***
122!*** CREATE A 3-D "HEIGHT MASK" FOR THE SPECIFIED PRESSURE LEVELS
123!*** (1 => ABOVE GROUND) AND A 2-D INDICATOR ARRAY THAT SAYS
124!*** WHICH PRESSURE LEVEL IS THE LOWEST ONE ABOVE THE GROUND
125!***
126 DO l=1,lsm
127 spll = spl(l)
128!
129!$omp parallel do private(j,i,psfc,pchk)
130 DO j=jsta,jend
131 DO i=ista,iend
132
133 htmo(i,j,l)=1.
134 if(pslp(i,j)<spval) then
135
136 psfc = pslp(i,j)
137 pchk = psfc
138 IF(nfill > 0) THEN
139 pchk = pint(i,j,nint(lmh(i,j))+1-nfill)
140 ENDIF
141 IF(fis(i,j) < 1.) pchk = psfc
142!
143 IF(spll < pchk) THEN
144 htmo(i,j,l) = 1.
145 ELSE
146 htmo(i,j,l) = 0.
147 IF(l > 1 .AND. htmo(i,j,l-1) > 0.5) lmho(i,j) = l-1
148 ENDIF
149 IF(l == lsm .AND. htmo(i,j,l) > 0.5) lmho(i,j) = lsm
150!
151! test new idea of filtering above-ground pressure levels for Gibsing
152! IF(L==LSM.AND.HTMO(I,J,L)>0.5)THEN
153! IF(FIS(I,J)>0.)THEN
154! LMHO(I,J)=LSM
155! ELSE
156! LMHO(I,J)=LSM-2
157! HTMO(I,J,LSM)=0.
158! HTMO(I,J,LSM-1)=0.
159! END IF
160! END IF
161! if(i==ii.and.j==jj)print*,'Debug: HTMO= ',HTMO(I,J,L)
162
163 endif !if pslp
164
165 ENDDO
166 ENDDO
167!
168 ENDDO
169! if(jj>=jsta.and.jj<=jend) print*,'Debug: LMHO=',LMHO(ii,jj)
170!--------------------------------------------------------------------
171!***
172!*** WE REACH THIS LINE IF WE WANT THE MESINGER ETA SLP REDUCTION
173!*** BASED ON RELAXATION TEMPERATURES. THE FIRST STEP IS TO
174!*** FIND THE HIGHEST LAYER CONTAINING MOUNTAINS.
175!***
176 lhmnt = lsm
177 loop210: DO l=lsm,1,-1
178
179 DO j=jsta,jend
180 DO i=ista,iend
181 if(pslp(i,j)<spval) then
182 IF(htmo(i,j,l) < 0.5) cycle loop210
183 endif
184 ENDDO
185 ENDDO
186 lhmnt = l+1
187 EXIT loop210
188 ENDDO loop210
189 !210 continue
190 !220 continue
191
192! print*,'Debug in SLP: LHMNT=',LHMNT
193
194 if ( num_procs > 1 ) then
195 CALL mpi_allreduce &
196 (lhmnt,lxxx,1,mpi_integer,mpi_min,mpi_comm_comp,ierr)
197 lhmnt = lxxx
198 end if
199
200 IF(lhmnt == lsmp1) GO TO 325
201
202! print*,'Debug in SLP: LHMNT A ALLREDUCE=',LHMNT
203!***
204!*** NOW GATHER THE ADDRESSES OF ALL THE UNDERGROUND POINTS.
205!***
206!!$omp parallel do private(kmn,kount)
207 DO l=lhmnt,lsm
208 kmn = 0
209 kmntm(l) = 0
210 kount = 0
211! DO 240 J=JSTA_M2,JEND_M2
212 DO j=jsta_m,jend_m
213 DO i=ista_m,iend_m
214 if(pslp(i,j)<spval) then
215 kount = kount + 1
216 imnt(kount,l) = 0
217 jmnt(kount,l) = 0
218 IF(htmo(i,j,l) > 0.5) cycle
219 kmn = kmn + 1
220 imnt(kmn,l) = i
221 jmnt(kmn,l) = j
222 endif
223 enddo
224 enddo
225 kmntm(l) = kmn
226 enddo
227!
228!
229!*** CREATE A TEMPORARY TV ARRAY, AND FOLLOW BY SEQUENTIAL
230!*** OVERRELAXATION, DOING NRLX PASSES.
231!
232! IF(NTSD==1)THEN
233 nrlx = nrlx2
234! ELSE
235! NRLX=NRLX2
236! ENDIF
237!
238!!$omp parallel do private(i,j,ttv,tem,kmma (Can this loop be threaded?))
239 DO l=lhmnt,lsm
240!
241!$omp parallel do private(i,j)
242 DO j=jsta,jend
243 DO i=ista,iend
244! dong
245! if (QPRES(I,J,LSM) < spval) then
246 !if(PSLP(I,J)<spval) then
247 ttv(i,j) = tpres(i,j,l)
248 htm2d(i,j) = htmo(i,j,l)
249 !end if ! spval if
250! end if ! spval if
251! IF(TTV(I,J)<150. .and. TTV(I,J)>325.0)print* &
252! ,'abnormal IC for T relaxation',i,j,TTV(I,J)
253 enddo
254 enddo
255!
256!*** FOR GRID BOXES NEXT TO MOUNTAINS, COMPUTE TV TO USE AS
257!*** BOUNDARY CONDITIONS FOR THE RELAXATION UNDERGROUND
258!
259 CALL exch(htm2d(ista_2l,jsta_2l)) !ONLY NEED TO EXCHANGE ONE ROW FOR A/C GRID
260! DO J=JSTA_M2,JEND_M2
261!$omp parallel do private(i,j,tem)
262 DO j=jsta_m,jend_m
263 DO i=ista_m,iend_m
264
265 if(pslp(i,j)<spval) then
266
267!HC IF(HTM2D(I,J,L)>0.5.AND.
268!HC 1 HTM2D(I+IHW(J),J-1,L)*HTM2D(I+IHE(J),J-1,L)
269!HC 2 *HTM2D(I+IHW(J),J+1,L)*HTM2D(I+IHE(J),J+1,L)
270!HC 3 *HTM2D(I-1 ,J ,L)*HTM2D(I+1 ,J ,L)
271!HC 4 *HTM2D(I ,J-2,L)*HTM2D(I ,J+2,L)<0.5)THEN
272!HC MODIFICATION FOR C AND A GRIDS
273
274 tem = htm2d(i-1,j)*htm2d(i+1,j)*htm2d(i,j-1)*htm2d(i,j+1) &
275 * htm2d(i-1,j-1)*htm2d(i+1,j-1)*htm2d(i-1,j+1)*htm2d(i+1,j+1)
276 IF(htm2d(i,j) > 0.5 .AND. tem < 0.5) then
277 ttv(i,j) = tpres(i,j,l)*(1.+0.608*qpres(i,j,l))
278 ENDIF
279! if(i==ii.and.j==jj)print*,'Debug:L,TTV B SMOO= ',l,TTV(I,J)
280 end if ! spval
281 ENDDO
282 ENDDO
283!
284 kmm = kmntm(l)
285
286! print*,'Debug:L,KMM=',L,KMM
287!
288 DO n=1,nrlx
289 CALL exch(ttv(ista_2l,jsta_2l))
290!!$omp parallel do private(i,j,km,a1,a2,a3,a4,a5,a6,a7,a8)
291 DO km=1,kmm
292 i = imnt(km,l)
293 j = jmnt(km,l)
294
295 if(pslp(i,j)<spval) then
296
297!HC TTV(I,J)=AD05*(4.*(TTV(I+IHW(J),J-1)+TTV(I+IHE(J),J-1)
298!HC 1 +TTV(I+IHW(J),J+1)+TTV(I+IHE(J),J+1))
299!HC 2 +TTV(I-1,J) +TTV(I+1,J)
300!HC 3 +TTV(I,J-2) +TTV(I,J+2))
301!HC 4 -CFT0*TTV(I,J)
302!HC MODIFICATION FOR C AND A GRIDS
303! eight point relaxation using updated TTV to the lower and left
304! TTV(I,J)=AD05*(4.*(TTV(I-1,J)+TTV(I+1,J)
305! 1 +TTV(I,J-1)+TTV(I,J+1))
306! 2 +TTV(I-1,J-1)+TTV(I+1,J-1)
307! 3 +TTV(I-1,J+1)+TTV(I+1,J+1))
308! 4 -CFT0*TTV(I,J)
309! eight point relaxation using old TTV
310 a1=ttv(i-1,j)
311 a2=ttv(i+1,j)
312 a3=ttv(i,j-1)
313 a4=ttv(i,j+1)
314 a5=ttv(i-1,j-1)
315 a6=ttv(i+1,j-1)
316 a7=ttv(i-1,j+1)
317 a8=ttv(i+1,j+1)
318! if ((a1-spval) <= 1e-10) a1=TTV(I,J)
319! if ((a2-spval) <= 1e-10) a2=TTV(I,J)
320! if ((a3-spval) <= 1e-10) a3=TTV(I,J)
321! if ((a4-spval) <= 1e-10) a4=TTV(I,J)
322! if ((a5-spval) <= 1e-10) a5=TTV(I,J)
323! if ((a6-spval) <= 1e-10) a6=TTV(I,J)
324! if ((a7-spval) <= 1e-10) a7=TTV(I,J)
325! if ((a8-spval) <= 1e-10) a8=TTV(I,J)
326
327 if ((a1 < spval) .and. &
328 (a2 < spval) .and. &
329 (a3 < spval) .and. &
330 (a4 < spval) .and. &
331 (a5 < spval) .and. &
332 (a6 < spval) .and. &
333 (a7 < spval) .and. &
334 (a8 < spval) .and. (ttv(i,j) < spval)) then
335
336! TNEW(I,J) = AD05*(4.*(a1 +a2 +a3 &
337! +a4) +a5 +a6 &
338! +a7+a8)-TTV(I,J)*CFT0
339
340 tnew(i,j) = ad05*(4.*(ttv(i-1,j) +ttv(i+1,j) +ttv(i,j-1) &
341 +ttv(i,j+1)) +ttv(i-1,j-1) +ttv(i+1,j-1) &
342 +ttv(i-1,j+1)+ttv(i+1,j+1))-ttv(i,j)*cft0
343 else
344 tnew(i,j) = ttv(i,j)
345 end if ! spval
346
347! four point relaxation using old TTV
348! TNEW(I,J)=TTV(I,J)+1.0*((TTV(I-1,J)+TTV(I+1,J)
349! 1 +TTV(I,J-1)+TTV(I,J+1)-4.0*TTV(I,J))/4.0)
350! four point relaxation using updated TTV to the lower and left
351! TTV(I,J)=TTV(I,J)+1.0*((TTV(I-1,J)+TTV(I+1,J)
352! 1 +TTV(I,J-1)+TTV(I,J+1)-4.0*TTV(I,J))/4.0)
353!
354! if(i==ii.and.j==jj)print*,'Debug: L,TTV A S'
355! 1,l,TTV(I,J),N
356! 1,l,TNEW(I,J),N
357
358 end if ! spval
359
360 enddo
361!
362!!$omp parallel do private(i,j,km)
363 DO km=1,kmm
364 i = imnt(km,l)
365 j = jmnt(km,l)
366 if(pslp(i,j)<spval .and. tnew(i,j)< spval/100.) then
367 ttv(i,j) = tnew(i,j)
368 end if ! spval
369 END DO
370 END DO ! NRLX loop
371!
372!!$omp parallel do private(i,j,km)
373 DO km=1,kmm
374 i = imnt(km,l)
375 j = jmnt(km,l)
376
377 if(pslp(i,j)<spval) then
378
379! dong try to fix missing value for hgtprs at 1000 mb
380 tpres(i,j,l) = ttv(i,j)
381 end if ! spval
382
383! if (QPRES(I,J,L) < 1000) TPRES(I,J,L) = TTV(I,J)
384! if (QPRES(I,J,L) < 1000) TPRES(I,J,L) = 1
385
386 END DO
387 enddo ! end of l loop
388!----------------------------------------------------------------
389!***
390!*** CALCULATE THE SEA LEVEL PRESSURE AS PER THE NEW SCHEME.
391!*** INTEGRATE THE HYDROSTATIC EQUATION DOWNWARD FROM THE
392!*** GROUND THROUGH EACH OUTPUT PRESSURE LEVEL (WHERE TV
393!*** IS NOW KNOWN) TO FIND GZ AT THE NEXT MIDPOINT BETWEEN
394!*** PRESSURE LEVELS. WHEN GZ=0 IS REACHED, SOLVE FOR THE
395!*** PRESSURE.
396!***
397!
398!*** BEFORE APPLYING RELAXATION FOR UNDERGROUND POINTS,
399!*** FIRST FIND GRID POINTS AT/NEAR/BELOW SEA LEVEL AND DERIVE
400!*** SEA LEVEL PRESSURE TO AVOID MEMBRANE RELAXATION
401!*** AT THESE GRID POINTS. E.G. HURRICANE CENTER NEAR COAST
402!
403 kount = 0
404 DO j=jsta,jend
405 DO i=ista,iend
406
407 if(pslp(i,j)<spval) then
408! P1(I,J)=SPL(NINT(LMH(I,J)))
409! DONE(I,J)=.FALSE.
410
411 IF(abs(fis(i,j)) < 1.) THEN
412 pslp(i,j) = pint(i,j,nint(lmh(i,j))+1)
413 done(i,j) = .true.
414 kount = kount + 1
415! if(i==ii.and.j==jj)print*,'Debug:DONE,PSLP A S1=' &
416! ,done(i,j),PSLP(I,J)
417 ELSE IF(fis(i,j) < -1.0) THEN
418 DO l=lm,1,-1
419 IF(zint(i,j,l) > 0.)THEN
420! PSLP(I,J)=PINT(I,J,L)/EXP(-ZINT(I,J,L)*G &
421! /(RD*T(I,J,L)*(Q(I,J,L)*D608+1.0)))
422
423 tem = 0.5*(t(i,j,l)+t(i,j,l-1))*(1.0+0.5*d608*(q(i,j,l)+q(i,j,l-1)))
424 pslp(i,j) = pint(i,j,l-1)/exp(-zint(i,j,l-1)*g/(rd*tem))
425 done(i,j) = .true.
426! if(i==ii.and.j==jj)print* &
427! ,'Debug:DONE,PINT,PSLP A S1=' &
428! ,done(i,j),PINT(I,J,L),PSLP(I,J)
429 exit
430 END IF
431 END DO
432 ENDIF
433
434 end if ! spval
435
436 ENDDO
437 ENDDO
438!
439 kmm = kmntm(lsm)
440!!$omp parallel do private(gz1,gz2,i,j,lmap1,p1,p2),shared(pslp)
441loop320:DO km=1,kmm
442 i = imnt(km,lsm)
443 j = jmnt(km,lsm)
444
445 if(pslp(i,j)<spval) then
446
447 IF(done(i,j)) cycle
448 lmhij = lmho(i,j)
449 gz1 = fipres(i,j,lmhij)
450 p1(i,j) = spl(lmhij)
451!
452 lmap1 = lmhij+1
453 DO l=lmap1,lsm
454 IF(gz1<spval .AND. p1(i,j)<spval .AND. tpres(i,j,l)<spval .AND. tpres(i,j,l-1)<spval) THEN
455 p2 = spl(l)
456 tlyr = 0.5*(tpres(i,j,l)+tpres(i,j,l-1))
457 gz2 = gz1 + rd*tlyr*log(p1(i,j)/p2)
458 fipres(i,j,l) = gz2
459! if(i==ii.and.j==jj)print*,'Debug:L,FI A S2=',L,GZ2
460 IF(gz2 <= 0.)THEN
461 pslp(i,j) = p1(i,j)/exp(-gz1/(rd*tpres(i,j,l-1)))
462! if(i==ii.and.j==jj)print*,'Debug:PSLP A S2=',PSLP(I,J)
463 done(i,j) = .true.
464 kount = kount + 1
465 cycle loop320
466 ENDIF
467 p1(i,j) = p2
468 gz1 = gz2
469 ENDIF
470 ENDDO
471!HC EXPERIMENT
472 lp = lsm
473 slope = -6.6e-4
474 IF(tpres(i,j,lp)<spval .AND. fipres(i,j,lp)<spval .AND. spl(lp)<spval ) THEN
475 tlyr = tpres(i,j,lp)-0.5*fipres(i,j,lp)*slope
476 pslp(i,j) = spl(lp)/exp(-fipres(i,j,lp)/(rd*tlyr))
477 ENDIF
478 done(i,j) = .true.
479! if(i==ii.and.j==jj)print*,'Debug:spl,FI,TLYR,PSLPA3=' &
480! ,spl(lp),FIPRES(I,J,LP),TLYR,PSLP(I,J)
481!HC EXPERIMENT
482 end if ! spval
483ENDDO loop320
484 320 CONTINUE
485!
486!*** WHEN SEA LEVEL IS BELOW THE LOWEST OUTPUT PRESSURE LEVEL,
487!*** SOLVE THE HYDROSTATIC EQUATION BY CHOOSING A TEMPERATURE
488!*** AT THE MIDPOINT OF THE LAYER BETWEEN THAT LOWEST PRESSURE
489!*** LEVEL AND THE GROUND BY EXTRAPOLATING DOWNWARD FROM T ON
490!*** THE LOWEST PRESSURE LEVEL USING THE DT/DFI BETWEEN THE
491!*** LOWEST PRESSURE LEVEL AND THE ONE ABOVE IT.
492!
493! TOTAL=(IM-2)*(JM-4)
494!
495!HC DO 340 LP=LSM,1,-1
496! IF(KOUNT==TOTAL)GO TO 350
497!HC MODIFICATION FOR SMALL HILL HIGH PRESSURE SITUATION
498!HC IF SURFACE PRESSURE IS CLOSER TO SEA LEVEL THAN LWOEST
499!HC OUTPUT PRESSURE LEVEL, USE SURFACE PRESSURE TO DO EXTRAPOLATION
500
501 325 CONTINUE
502 lp = lsm
503 DO j=jsta,jend
504 DO i=ista,iend
505
506 if(pslp(i,j)<spval) then
507
508! if(i==ii.and.j==jj)print*,'Debug: with 330 loop'
509 IF(done(i,j)) cycle
510
511! if(i==ii.and.j==jj)print*,'Debug: still within 330 loop'
512!HC Comment out the following line for situation with terrain
513!HC at boundary (ie FIPRES<0)
514!HC because they were not counted as undergound point for 8 pt
515!HC relaxation
516!HC IF(FIPRES(I,J,LP)<0.)GO TO 330
517! IF(FIPRES(I,J,LP)<0.)THEN
518! DO LP=LSM,1,-1
519! IF (FIPRES(I,J) <= 0)
520
521! IF(FIPRES(I,J,LP)<0..OR.DONE(I,J))GO TO 330
522! SLOPE=(TPRES(I,J,LP)-TPRES(I,J,LP-1))
523! & /(FIPRES(I,J,LP)-FIPRES(I,J,LP-1))
524
525 slope = -6.6e-4
526
527 IF(pint(i,j,nint(lmh(i,j))+1) > spl(lp))THEN
528 llmh = nint(lmh(i,j))
529 IF(t(i,j,llmh)<spval .AND. q(i,j,llmh)<spval .AND. &
530 zint(i,j,llmh)<spval .AND. zint(i,j,llmh+1)<spval .AND. &
531 pint(i,j,llmh+1)<spval) THEN
532 tvrt = t(i,j,llmh)*(h1+d608*q(i,j,llmh))
533 dis = zint(i,j,llmh+1)-zint(i,j,llmh)+0.5*zint(i,j,llmh+1)
534 tlyr = tvrt-dis*g*slope
535 pslp(i,j) = pint(i,j,llmh+1)*exp(zint(i,j,llmh+1)*g &
536 /(rd*tlyr))
537 ENDIF
538! if(i==ii.and.j==jj)print*,'Debug:PSFC,zsfc,TLYR,PSLPA3='
539! 1,PINT(I,J,LLMH+1),ZINT(I,J,LLMH+1),TLYR,PSLP(I,J)
540 ELSE
541 IF(tpres(i,j,lp)<spval .AND. fipres(i,j,lp)<spval .AND. spl(lp)<spval ) THEN
542 tlyr=tpres(i,j,lp)-0.5*fipres(i,j,lp)*slope
543 pslp(i,j)=spl(lp)/exp(-fipres(i,j,lp)/(rd*tlyr))
544 ENDIF
545! if(i==ii.and.j==jj)print*,'Debug:spl,FI,TLYR,PSLPA3=' &
546! ,spl(lp),FIPRES(I,J,LP),TLYR,PSLP(I,J)
547 END IF
548 done(i,j) = .true.
549 kount = kount + 1
550 end if ! spval
551
552 enddo
553 enddo
554!HC 340 CONTINUE
555!
556! 350 CONTINUE
557!--------------------------------------------------------------------
558 RETURN
559 END
exch
subroutine exch(a)
exch() Subroutine that exchanges one halo row.
Definition EXCH.f:27