module_sf_myjsfc.F
References to this file elsewhere.
1 !----------------------------------------------------------------------
2 !
3 MODULE MODULE_SF_MYJSFC
4 !
5 !----------------------------------------------------------------------
6 !
7 USE MODULE_MODEL_CONSTANTS
8 USE MODULE_DM, ONLY : WRF_DM_MAXVAL
9 !
10 !----------------------------------------------------------------------
11 !
12 ! REFERENCES: Janjic (2002), NCEP Office Note 437
13 !
14 ! ABSTRACT:
15 ! MYJSFC GENERATES THE SURFACE EXCHANGE COEFFICIENTS FOR VERTICAL
16 ! TURBULENT EXCHANGE BASED UPON MONIN_OBUKHOV THEORY WITH
17 ! VARIOUS REFINEMENTS.
18 !
19 !----------------------------------------------------------------------
20 !
21 INTEGER :: ITRMX=5 ! Iteration count for sfc layer computations
22 !
23 REAL,PARAMETER :: VKARMAN=0.4
24 REAL,PARAMETER :: CAPA=R_D/CP,ELOCP=2.72E6/CP,RCAP=1./CAPA
25 REAL,PARAMETER :: GOCP02=G/CP*2.,GOCP10=G/CP*10.
26 ! REAL,PARAMETER :: EPSU2=1.E-4,EPSUST=0.07,EPSZT=1.E-5
27 ! ECMWF sets lower limit on ln(Z0h)=-20 --> EPSZT=2.e-9
28 REAL,PARAMETER :: EPSU2=1.E-6,EPSUST=1.e-9,EPSZT=1.E-28
29 REAL,PARAMETER :: A2S=17.2693882,A3S=273.16,A4S=35.86
30 REAL,PARAMETER :: SEAFC=0.98,PQ0SEA=PQ0*SEAFC
31 REAL,PARAMETER :: BETA=1./273.,CZIL=0.1,EXCML=0.0001,EXCMS=0.0001 &
32 !old REAL,PARAMETER :: BETA=1./273.,CZIL=0.1,EXCML=0.001,EXCMS=0.001 &
33 & ,GLKBR=10.,GLKBS=30.,PI=3.1415926 &
34 & ,QVISC=2.1E-5,RIC=0.505,SMALL=0.35 &
35 & ,SQPR=0.84,SQSC=0.84,SQVISC=258.2,TVISC=2.1E-5 &
36 & ,USTC=0.7,USTR=0.225,VISC=1.5E-5,FH=1.01 &
37 & ,WWST=1.2,ZTFC=1.,TOPOFAC=9.0e-6
38 !
39 REAL,PARAMETER :: BTG=BETA*G,CZIV=SMALL*GLKBS &
40 & ,GRRS=GLKBR/GLKBS &
41 & ,RTVISC=1./TVISC,RVISC=1./VISC &
42 & ,ZQRZT=SQSC/SQPR &
43 & ,FZQ1=RTVISC*QVISC*ZQRZT &
44 & ,FZQ2=RTVISC*QVISC*ZQRZT &
45 & ,FZT1=RVISC *TVISC*SQPR &
46 & ,FZT2=CZIV*GRRS*TVISC*SQPR &
47 & ,FZU1=CZIV*VISC &
48 & ,PIHF=0.5*PI &
49 & ,RQVISC=1./QVISC &
50 & ,USTFC=0.018/G &
51 & ,WWST2=WWST*WWST &
52 & ,ZILFC=-CZIL*VKARMAN*SQVISC
53 !
54 !----------------------------------------------------------------------
55 INTEGER, PARAMETER :: KZTM=10001,KZTM2=KZTM-2
56 !
57 REAL :: DZETA1,DZETA2,FH01,FH02,ZTMAX1,ZTMAX2,ZTMIN1,ZTMIN2
58 !
59 REAL,DIMENSION(KZTM) :: PSIH1,PSIH2,PSIM1,PSIM2
60 !
61 !----------------------------------------------------------------------
62 !
63 CONTAINS
64 !
65 !----------------------------------------------------------------------
66 SUBROUTINE MYJSFC(ITIMESTEP,HT,DZ &
67 & ,PMID,PINT,TH,T,QV,QC,U,V,Q2 &
68 & ,TSK,QSFC,THZ0,QZ0,UZ0,VZ0 &
69 & ,LOWLYR,XLAND &
70 & ,USTAR,ZNT,Z0BASE,PBLH,MAVAIL,RMOL &
71 & ,AKHS,AKMS &
72 & ,CHS,CHS2,CQS2,HFX,QFX,FLX_LH,FLHC,FLQC &
73 & ,QGH,CPM,CT &
74 & ,U10,V10,T02,TH02,TSHLTR,TH10,Q02,QSHLTR,Q10,PSHLTR &
75 & ,IDS,IDE,JDS,JDE,KDS,KDE &
76 & ,IMS,IME,JMS,JME,KMS,KME &
77 & ,ITS,ITE,JTS,JTE,KTS,KTE)
78 !----------------------------------------------------------------------
79 !
80 IMPLICIT NONE
81 !
82 !----------------------------------------------------------------------
83 INTEGER,INTENT(IN) :: IDS,IDE,JDS,JDE,KDS,KDE &
84 & ,IMS,IME,JMS,JME,KMS,KME &
85 & ,ITS,ITE,JTS,JTE,KTS,KTE
86 !
87 INTEGER,INTENT(IN) :: ITIMESTEP
88 !
89 INTEGER,DIMENSION(IMS:IME,JMS:JME),INTENT(IN) :: LOWLYR
90 !
91 REAL,DIMENSION(IMS:IME,JMS:JME),INTENT(IN) :: HT,MAVAIL,TSK &
92 & ,XLAND,Z0BASE
93 !
94 REAL,DIMENSION(IMS:IME,KMS:KME,JMS:JME),INTENT(IN) :: DZ &
95 & ,PMID,PINT &
96 & ,Q2,QC,QV &
97 & ,T,TH &
98 & ,U,V
99 !
100 REAL,DIMENSION(IMS:IME,JMS:JME),INTENT(OUT) :: FLX_LH,HFX,PSHLTR &
101 & ,QFX,Q10,QSHLTR &
102 & ,TH10,TSHLTR,T02 &
103 & ,U10,V10,TH02,Q02
104 !
105 REAL,DIMENSION(IMS:IME,JMS:JME),INTENT(INOUT) :: AKHS,AKMS &
106 & ,PBLH,QSFC
107 !
108 REAL,DIMENSION(IMS:IME,JMS:JME),INTENT(INOUT) :: QZ0,RMOL,THZ0 &
109 & ,USTAR,UZ0,VZ0 &
110 & ,ZNT
111 !
112 REAL,DIMENSION(IMS:IME,JMS:JME),INTENT(OUT) :: CHS,CHS2,CQS2 &
113 & ,CPM,CT,FLHC,FLQC &
114 & ,QGH
115 !----------------------------------------------------------------------
116 !***
117 !*** LOCAL VARIABLES
118 !***
119 INTEGER :: I,J,K,KFLIP,LMH,LPBL,NTSD
120 !
121 REAL :: A,APESFC,B,BTGX,CWMLOW &
122 & ,DQDT,DTDIF,DTDT,DUDT,DVDT &
123 & ,FIS &
124 & ,P02P,P10P,PLOW,PSFC,PTOP,QLOW,QS02,QS10 &
125 & ,RAPA,RAPA02,RAPA10,RATIOMX,RDZ,SEAMASK,SM &
126 & ,T02P,T10P,TEM,TH02P,TH10P,THLOW,THELOW,THM &
127 & ,TLOW,TZ0,ULOW,VLOW,ZSL
128 !
129 REAL,DIMENSION(KTS:KTE) :: CWMK,PK,Q2K,QK,THEK,THK,TK,UK,VK
130 !
131 REAL,DIMENSION(KTS:KTE-1) :: EL,ELM
132 !
133 REAL,DIMENSION(KTS:KTE+1) :: ZHK
134 !
135 REAL,DIMENSION(ITS:ITE,JTS:JTE) :: THSK
136 !
137 REAL,DIMENSION(ITS:ITE,KTS:KTE+1,JTS:JTE) :: ZINT
138 !
139 !----------------------------------------------------------------------
140 !**********************************************************************
141 !----------------------------------------------------------------------
142 !
143 !*** MAKE PREPARATIONS
144 !
145 !----------------------------------------------------------------------
146 DO J=JTS,JTE
147 DO K=KTS,KTE+1
148 DO I=ITS,ITE
149 ZINT(I,K,J)=0.
150 ENDDO
151 ENDDO
152 ENDDO
153 !
154 DO J=JTS,JTE
155 DO I=ITS,ITE
156 ZINT(I,KTE+1,J)=HT(I,J) ! Z at bottom of lowest sigma layer
157 PBLH(I,J)=-1.
158 !
159 !!!!!!!!!
160 !!!!!! UNCOMMENT THESE LINES IF USING ETA COORDINATES
161 !!!!!!!!!
162 !!!!!! ZINT(I,KTE+1,J)=1.E-4 ! Z of bottom of lowest eta layer
163 !!!!!! ZHK(KTE+1)=1.E-4 ! Z of bottom of lowest eta layer
164 !
165 ENDDO
166 ENDDO
167 !
168 DO J=JTS,JTE
169 DO K=KTE,KTS,-1
170 KFLIP=KTE+1-K
171 DO I=ITS,ITE
172 ZINT(I,K,J)=ZINT(I,K+1,J)+DZ(I,KFLIP,J)
173 ENDDO
174 ENDDO
175 ENDDO
176 !
177 NTSD=ITIMESTEP
178 !
179 #if ( NMM_CORE == 1 )
180 if(NTSD+1.eq.1) then
181 #else
182 IF(NTSD==1)THEN
183 #endif
184 !tgs IF(NTSD==1)THEN
185 DO J=JTS,JTE
186 DO I=ITS,ITE
187 USTAR(I,J)=0.1
188 FIS=HT(I,J)*G
189 SM=XLAND(I,J)-1.
190 !!! Z0(I,J)=SM*Z0SEA+(1.-SM)*(Z0(I,J)*Z0MAX+FIS*FCM+Z0LAND)
191 !!! ZNT(I,J)=SM*Z0SEA+(1.-SM)*(ZNT(I,J)*Z0MAX+FIS*FCM+Z0LAND)
192 ENDDO
193 ENDDO
194 ENDIF
195 !
196 !!!! IF(NTSD==1)THEN
197 DO J=JTS,JTE
198 DO I=ITS,ITE
199 CT(I,J)=0.
200 ENDDO
201 ENDDO
202 !!!! ENDIF
203 !
204 !----------------------------------------------------------------------
205 setup_integration: DO J=JTS,JTE
206 !----------------------------------------------------------------------
207 !
208 DO I=ITS,ITE
209 !
210 !*** LOWEST LAYER ABOVE GROUND MUST BE FLIPPED
211 !
212 LMH=KTE-LOWLYR(I,J)+1
213 !
214 PTOP=PINT(I,KTE+1,J) ! KTE+1=KME
215 PSFC=PINT(I,LOWLYR(I,J),J)
216 ! Define THSK here (for first timestep mostly)
217 THSK(I,J)=TSK(I,J)/(PSFC*1.E-5)**CAPA
218 !
219 !*** CONVERT LAND MASK (1 FOR SEA; 0 FOR LAND)
220 !
221 SEAMASK=XLAND(I,J)-1.
222 !
223 !*** FILL 1-D VERTICAL ARRAYS
224 !*** AND FLIP DIRECTION SINCE MYJ SCHEME
225 !*** COUNTS DOWNWARD FROM THE DOMAIN'S TOP
226 !
227 DO K=KTE,KTS,-1
228 KFLIP=KTE+1-K
229 THK(K)=TH(I,KFLIP,J)
230 TK(K)=T(I,KFLIP,J)
231 RATIOMX=QV(I,KFLIP,J)
232 QK(K)=RATIOMX/(1.+RATIOMX)
233 PK(K)=PMID(I,KFLIP,J)
234 CWMK(K)=QC(I,KFLIP,J)
235 THEK(K)=(CWMK(K)*(-ELOCP/TK(K))+1.)*THK(K)
236 Q2K(K)=2.*Q2(I,KFLIP,J)
237 !
238 !
239 !*** COMPUTE THE HEIGHTS OF THE LAYER INTERFACES
240 !
241 ZHK(K)=ZINT(I,K,J)
242 !
243 ENDDO
244 ZHK(KTE+1)=HT(I,J) ! Z at bottom of lowest sigma layer
245 !
246 DO K=KTE,KTS,-1
247 KFLIP=KTE+1-K
248 UK(K)=U(I,KFLIP,J)
249 VK(K)=V(I,KFLIP,J)
250 ENDDO
251 !
252 !*** FIND THE HEIGHT OF THE PBL
253 !
254 LPBL=LMH
255 DO K=LMH-1,1,-1
256 IF(Q2K(K)<=EPSQ2*FH) THEN
257 LPBL=K
258 GO TO 110
259 ENDIF
260 ENDDO
261 !
262 LPBL=1
263 !
264 !-----------------------------------------------------------------------
265 !--------------THE HEIGHT OF THE PBL------------------------------------
266 !-----------------------------------------------------------------------
267 !
268 110 PBLH(I,J)=ZHK(LPBL)-ZHK(LMH+1)
269 !
270 !----------------------------------------------------------------------
271 !***
272 !*** FIND THE SURFACE EXCHANGE COEFFICIENTS
273 !***
274 !----------------------------------------------------------------------
275 PLOW=PK(LMH)
276 TLOW=TK(LMH)
277 THLOW=THK(LMH)
278 THELOW=THEK(LMH)
279 QLOW=QK(LMH)
280 CWMLOW=CWMK(LMH)
281 ULOW=UK(LMH)
282 VLOW=VK(LMH)
283 ZSL=(ZHK(LMH)-ZHK(LMH+1))*0.5
284 APESFC=(PSFC*1.E-5)**CAPA
285 !tgs - in ARW THZ0 is not initialized when MYJSFC is called first time
286 #if ( NMM_CORE == 1 )
287 if(NTSD+1.eq.1) then
288 #else
289 IF(NTSD==1)THEN
290 #endif
291 ! if(itimestep.le.1) then
292 TZ0=TSK(I,J)
293 else
294 TZ0=THZ0(I,J)*APESFC
295 endif
296 !
297 CALL SFCDIF(NTSD,SEAMASK,THSK(I,J),QSFC(I,J),PSFC &
298 & ,UZ0(I,J),VZ0(I,J),TZ0,THZ0(I,J),QZ0(I,J) &
299 & ,USTAR(I,J),ZNT(I,J),Z0BASE(I,J),CT(I,J),RMOL(I,J) &
300 & ,AKMS(I,J),AKHS(I,J),PBLH(I,J),MAVAIL(I,J) &
301 & ,CHS(I,J),CHS2(I,J),CQS2(I,J) &
302 & ,HFX(I,J),QFX(I,J),FLX_LH(I,J) &
303 & ,FLHC(I,J),FLQC(I,J),QGH(I,J),CPM(I,J) &
304 & ,ULOW,VLOW,TLOW,THLOW,THELOW,QLOW,CWMLOW &
305 & ,ZSL,PLOW &
306 & ,U10(I,J),V10(I,J),TSHLTR(I,J),TH10(I,J) &
307 & ,QSHLTR(I,J),Q10(I,J),PSHLTR(I,J) &
308 & ,IDS,IDE,JDS,JDE,KDS,KDE &
309 & ,IMS,IME,JMS,JME,KMS,KME &
310 & ,ITS,ITE,JTS,JTE,KTS,KTE,i,j,ZHK(LMH+1))
311 !
312 !*** REMOVE SUPERATURATION AT 2M AND 10M
313 !
314 RAPA=APESFC
315 TH02P=TSHLTR(I,J)
316 TH10P=TH10(I,J)
317 TH02(I,J)=TSHLTR(I,J)
318 !
319 RAPA02=RAPA-GOCP02/TH02P
320 RAPA10=RAPA-GOCP10/TH10P
321 !
322 T02P=TH02P*RAPA02
323 T10P=TH10P*RAPA10
324 ! 1 may 06 tgs T02(I,J) = T02P
325 T02(I,J) = TH02(I,J)*APESFC
326 !
327 P02P=(RAPA02**RCAP)*1.E5
328 P10P=(RAPA10**RCAP)*1.E5
329 !
330 QS02=PQ0/P02P*EXP(A2*(T02P-A3)/(T02P-A4))
331 QS10=PQ0/P10P*EXP(A2*(T10P-A3)/(T10P-A4))
332 !
333 IF(QSHLTR(I,J)>QS02)QSHLTR(I,J)=QS02
334 IF(Q10 (I,J)>QS10)Q10 (I,J)=QS10
335 Q02(I,J)=QSHLTR(I,J)/(1.-QSHLTR(I,J))
336 !----------------------------------------------------------------------
337 !
338 ENDDO
339 !
340 !----------------------------------------------------------------------
341 ENDDO setup_integration
342 !----------------------------------------------------------------------
343
344 END SUBROUTINE MYJSFC
345 !XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX
346 !----------------------------------------------------------------------
347 SUBROUTINE SFCDIF(NTSD,SEAMASK,THS,QS,PSFC &
348 & ,UZ0,VZ0,TZ0,THZ0,QZ0 &
349 & ,USTAR,Z0,Z0BASE,CT,RLMO,AKMS,AKHS,PBLH,WETM &
350 & ,CHS,CHS2,CQS2,HFX,QFX,FLX_LH,FLHC,FLQC,QGH,CPM &
351 & ,ULOW,VLOW,TLOW,THLOW,THELOW,QLOW,CWMLOW &
352 & ,ZSL,PLOW &
353 & ,U10,V10,TH02,TH10,Q02,Q10,PSHLTR &
354 & ,IDS,IDE,JDS,JDE,KDS,KDE &
355 & ,IMS,IME,JMS,JME,KMS,KME &
356 & ,ITS,ITE,JTS,JTE,KTS,KTE,i,j,ZSFC)
357 ! ****************************************************************
358 ! * *
359 ! * SURFACE LAYER *
360 ! * *
361 ! ****************************************************************
362 !----------------------------------------------------------------------
363 !
364 IMPLICIT NONE
365 !
366 !----------------------------------------------------------------------
367 INTEGER,INTENT(IN) :: IDS,IDE,JDS,JDE,KDS,KDE &
368 & ,IMS,IME,JMS,JME,KMS,KME &
369 & ,ITS,ITE,JTS,JTE,KTS,KTE,i,j
370 !
371 INTEGER,INTENT(IN) :: NTSD
372 !
373 REAL,INTENT(IN) :: CWMLOW,PBLH,PLOW,QLOW,PSFC,SEAMASK,ZSFC &
374 & ,THELOW,THLOW,THS,TLOW,TZ0,ULOW,VLOW,WETM,ZSL &
375 & ,Z0BASE
376 !
377 REAL,INTENT(OUT) :: CHS,CHS2,CPM,CQS2,CT,FLHC,FLQC,FLX_LH,HFX &
378 & ,PSHLTR,Q02,Q10,QFX,QGH,RLMO,TH02,TH10,U10,V10
379 !
380 REAL,INTENT(INOUT) :: AKHS,AKMS,QZ0,THZ0,USTAR,UZ0,VZ0,Z0,QS
381 !----------------------------------------------------------------------
382 !***
383 !*** LOCAL VARIABLES
384 !***
385 INTEGER :: ITR,K
386 !
387 REAL :: A,B,BTGH,BTGX,CXCHL,CXCHS,DTHV,DU2,ELFC,FCT &
388 & ,HLFLX,HSFLX,HV,PSH02,PSH10,PSHZ,PSHZL,PSM10,PSMZ,PSMZL &
389 & ,RDZ,RDZT,RIB,RLMA,RLMN,RLMP &
390 & ,RLOGT,RLOGU,RWGH,RZ,RZST,RZSU,SIMH,SIMM,TEM,THM &
391 & ,UMFLX,USTARK,VMFLX,WGHT,WGHTT,WGHTQ,WSTAR2 &
392 & ,X,XLT,XLT4,XLU,XLU4,XT,XT4,XU,XU4,ZETALT,ZETALU &
393 & ,ZETAT,ZETAU,ZQ,ZSLT,ZSLU,ZT,ZU,TOPOTERM,ZZIL
394 !
395 !*** DIAGNOSTICS
396 !
397 REAL :: AKHS02,AKHS10,AKMS02,AKMS10,EKMS10,QSAT10,QSAT2 &
398 & ,RLNT02,RLNT10,RLNU10,SIMH02,SIMH10,SIMM10,T02,T10 &
399 & ,TERM1,RLOW,U10E,V10E,WSTAR,XLT02,XLT024,XLT10 &
400 & ,XLT104,XLU10,XLU104,XU10,XU104,ZT02,ZT10,ZTAT02,ZTAT10 &
401 & ,ZTAU,ZTAU10,ZU10,ZUUZ
402 !----------------------------------------------------------------------
403 !**********************************************************************
404 !----------------------------------------------------------------------
405 RDZ=1./ZSL
406 CXCHL=EXCML*RDZ
407 CXCHS=EXCMS*RDZ
408 !
409 BTGX=G/THLOW
410 ELFC=VKARMAN*BTGX
411 !
412 IF(PBLH>1000.)THEN
413 BTGH=BTGX*PBLH
414 ELSE
415 BTGH=BTGX*1000.
416 ENDIF
417 !
418 !----------------------------------------------------------------------
419 !
420 !*** SEA POINTS
421 !
422 !----------------------------------------------------------------------
423 !
424 IF(SEAMASK>0.5)THEN
425 !
426 !----------------------------------------------------------------------
427 DO ITR=1,ITRMX
428 !----------------------------------------------------------------------
429 Z0=MAX(USTFC*USTAR*USTAR,1.59E-5)
430 !
431 !*** VISCOUS SUBLAYER, JANJIC MWR 1994
432 !
433 !----------------------------------------------------------------------
434 IF(USTAR<USTC)THEN
435 !----------------------------------------------------------------------
436 !
437 IF(USTAR<USTR)THEN
438 !
439 #if ( NMM_CORE == 1 )
440 if(NTSD+1.eq.1) then
441 #else
442 IF(NTSD==1)THEN
443 #endif
444 !tgs IF(NTSD==1)THEN
445 AKMS=CXCHS
446 AKHS=CXCHS
447 QS=QLOW
448 ENDIF
449 !
450 ZU=FZU1*SQRT(SQRT(Z0*USTAR*RVISC))/USTAR
451 WGHT=AKMS*ZU*RVISC
452 RWGH=WGHT/(WGHT+1.)
453 UZ0=(ULOW*RWGH+UZ0)*0.5
454 VZ0=(VLOW*RWGH+VZ0)*0.5
455 !
456 ZT=FZT1*ZU
457 ZQ=FZQ1*ZT
458 WGHTT=AKHS*ZT*RTVISC
459 WGHTQ=AKHS*ZQ*RQVISC
460 !
461 #if ( NMM_CORE == 1 )
462 if(NTSD+1>1) then
463 #else
464 IF(NTSD>1)THEN
465 #endif
466 !tgs IF(NTSD>1)THEN
467 THZ0=((WGHTT*THLOW+THS)/(WGHTT+1.)+THZ0)*0.5
468 QZ0=((WGHTQ*QLOW+QS)/(WGHTQ+1.)+QZ0)*0.5
469 ELSE
470 THZ0=(WGHTT*THLOW+THS)/(WGHTT+1.)
471 QZ0=(WGHTQ*QLOW+QS)/(WGHTQ+1.)
472 ENDIF
473 !
474 ENDIF
475 !
476 IF(USTAR>=USTR.AND.USTAR<USTC)THEN
477 ZU=Z0
478 UZ0=0.
479 VZ0=0.
480 !
481 ZT=FZT2*SQRT(SQRT(Z0*USTAR*RVISC))/USTAR
482 ZQ=FZQ2*ZT
483 WGHTT=AKHS*ZT*RTVISC
484 WGHTQ=AKHS*ZQ*RQVISC
485 !
486 #if ( NMM_CORE == 1 )
487 if(NTSD+1>1) then
488 #else
489 IF(NTSD>1)THEN
490 #endif
491
492 !tgs IF(NTSD>1)THEN
493 THZ0=((WGHTT*THLOW+THS)/(WGHTT+1.)+THZ0)*0.5
494 QZ0=((WGHTQ*QLOW+QS)/(WGHTQ+1.)+QZ0)*0.5
495 ELSE
496 THZ0=(WGHTT*THLOW+THS)/(WGHTT+1.)
497 QZ0=(WGHTQ*QLOW+QS)/(WGHTQ+1.)
498 ENDIF
499 !
500 ENDIF
501 !----------------------------------------------------------------------
502 ELSE
503 !----------------------------------------------------------------------
504 ZU=Z0
505 UZ0=0.
506 VZ0=0.
507 !
508 ZT=Z0
509 THZ0=THS
510 !
511 ZQ=Z0
512 QZ0=QS
513 !----------------------------------------------------------------------
514 ENDIF
515 !----------------------------------------------------------------------
516 TEM=(TLOW+TZ0)*0.5
517 THM=(THELOW+THZ0)*0.5
518 !
519 A=THM*P608
520 B=(ELOCP/TEM-1.-P608)*THM
521 !
522 DTHV=((THELOW-THZ0)*((QLOW+QZ0+CWMLOW)*(0.5*P608)+1.) &
523 & +(QLOW-QZ0+CWMLOW)*A+CWMLOW*B)
524 !
525 DU2=MAX((ULOW-UZ0)**2+(VLOW-VZ0)**2,EPSU2)
526 RIB=BTGX*DTHV*ZSL/DU2
527 !----------------------------------------------------------------------
528 ! IF(RIB>=RIC)THEN
529 !----------------------------------------------------------------------
530 ! AKMS=MAX( VISC*RDZ,CXCHS)
531 ! AKHS=MAX(TVISC*RDZ,CXCHS)
532 !----------------------------------------------------------------------
533 ! ELSE ! turbulent branch
534 !----------------------------------------------------------------------
535 ZSLU=ZSL+ZU
536 ZSLT=ZSL+ZT
537 !
538 RZSU=ZSLU/ZU
539 RZST=ZSLT/ZT
540 !
541 RLOGU=LOG(RZSU)
542 RLOGT=LOG(RZST)
543 !
544 !----------------------------------------------------------------------
545 !*** 1./MONIN-OBUKHOV LENGTH
546 !----------------------------------------------------------------------
547 !
548 RLMO=ELFC*AKHS*DTHV/USTAR**3
549 !
550 ZETALU=ZSLU*RLMO
551 ZETALT=ZSLT*RLMO
552 ZETAU=ZU*RLMO
553 ZETAT=ZT*RLMO
554 !
555 ZETALU=MIN(MAX(ZETALU,ZTMIN1),ZTMAX1)
556 ZETALT=MIN(MAX(ZETALT,ZTMIN1),ZTMAX1)
557 ZETAU=MIN(MAX(ZETAU,ZTMIN1/RZSU),ZTMAX1/RZSU)
558 ZETAT=MIN(MAX(ZETAT,ZTMIN1/RZST),ZTMAX1/RZST)
559 !
560 !----------------------------------------------------------------------
561 !*** WATER FUNCTIONS
562 !----------------------------------------------------------------------
563 !
564 RZ=(ZETAU-ZTMIN1)/DZETA1
565 K=INT(RZ)
566 RDZT=RZ-REAL(K)
567 K=MIN(K,KZTM2)
568 K=MAX(K,0)
569 PSMZ=(PSIM1(K+2)-PSIM1(K+1))*RDZT+PSIM1(K+1)
570 !
571 RZ=(ZETALU-ZTMIN1)/DZETA1
572 K=INT(RZ)
573 RDZT=RZ-REAL(K)
574 K=MIN(K,KZTM2)
575 K=MAX(K,0)
576 PSMZL=(PSIM1(K+2)-PSIM1(K+1))*RDZT+PSIM1(K+1)
577 !
578 SIMM=PSMZL-PSMZ+RLOGU
579 !
580 RZ=(ZETAT-ZTMIN1)/DZETA1
581 K=INT(RZ)
582 RDZT=RZ-REAL(K)
583 K=MIN(K,KZTM2)
584 K=MAX(K,0)
585 PSHZ=(PSIH1(K+2)-PSIH1(K+1))*RDZT+PSIH1(K+1)
586 !
587 RZ=(ZETALT-ZTMIN1)/DZETA1
588 K=INT(RZ)
589 RDZT=RZ-REAL(K)
590 K=MIN(K,KZTM2)
591 K=MAX(K,0)
592 PSHZL=(PSIH1(K+2)-PSIH1(K+1))*RDZT+PSIH1(K+1)
593 !
594 SIMH=(PSHZL-PSHZ+RLOGT)*FH01
595 !----------------------------------------------------------------------
596 USTARK=USTAR*VKARMAN
597 AKMS=MAX(USTARK/SIMM,CXCHS)
598 AKHS=MAX(USTARK/SIMH,CXCHS)
599 !
600 !----------------------------------------------------------------------
601 !*** BELJAARS CORRECTION FOR USTAR
602 !----------------------------------------------------------------------
603 !
604 IF(DTHV<=0.)THEN !zj
605 WSTAR2=WWST2*ABS(BTGH*AKHS*DTHV)**(2./3.) !zj
606 ELSE !zj
607 WSTAR2=0. !zj
608 ENDIF !zj
609 USTAR=MAX(SQRT(AKMS*SQRT(DU2+WSTAR2)),EPSUST)
610 !
611 !----------------------------------------------------------------------
612 ! ENDIF ! End of turbulent branch
613 !----------------------------------------------------------------------
614 !
615 ENDDO ! End of the iteration loop over sea points
616 !
617 !----------------------------------------------------------------------
618 !
619 !*** LAND POINTS
620 !
621 !----------------------------------------------------------------------
622 !
623 ELSE
624 !
625 !----------------------------------------------------------------------
626 IF(NTSD==1)THEN
627 QS=QLOW
628 ENDIF
629 !
630 ZU=Z0
631 UZ0=0.
632 VZ0=0.
633 !
634 ZT=ZU*ZTFC
635 THZ0=THS
636 !
637 ZQ=ZT
638 QZ0=QS
639 !----------------------------------------------------------------------
640 TEM=(TLOW+TZ0)*0.5
641 THM=(THELOW+THZ0)*0.5
642 !
643 A=THM*P608
644 B=(ELOCP/TEM-1.-P608)*THM
645 !
646 DTHV=((THELOW-THZ0)*((QLOW+QZ0+CWMLOW)*(0.5*P608)+1.) &
647 & +(QLOW-QZ0+CWMLOW)*A+CWMLOW*B)
648 !
649 DU2=MAX((ULOW-UZ0)**2+(VLOW-VZ0)**2,EPSU2)
650 RIB=BTGX*DTHV*ZSL/DU2
651 !----------------------------------------------------------------------
652 ! IF(RIB>=RIC)THEN
653 ! AKMS=MAX( VISC*RDZ,CXCHL)
654 ! AKHS=MAX(TVISC*RDZ,CXCHL)
655 !----------------------------------------------------------------------
656 ! ELSE ! Turbulent branch
657 !----------------------------------------------------------------------
658 ZSLU=ZSL+ZU
659 !
660 RZSU=ZSLU/ZU
661 !
662 RLOGU=LOG(RZSU)
663
664 ZSLT=ZSL+ZU ! u,v and t are at the same level
665 !----------------------------------------------------------------------
666 !
667 !mp Topo modification of ZILFC term
668 !
669 TOPOTERM=TOPOFAC*ZSFC**2.
670 TOPOTERM=MAX(TOPOTERM,3.0)
671 !
672 IF(DTHV>0.)THEN
673 ZZIL=ZILFC*TOPOTERM
674 ELSE
675 ZZIL=ZILFC
676 ENDIF
677 !
678 !----------------------------------------------------------------------
679 !
680 land_point_iteration: DO ITR=1,ITRMX
681 !
682 !----------------------------------------------------------------------
683 !*** ZILITINKEVITCH FIX FOR ZT
684 !----------------------------------------------------------------------
685 !
686 ! oldform ZT=MAX(EXP(ZZIL*SQRT(USTAR*ZU))*ZU,EPSZT)
687 ZT=MAX(EXP(ZZIL*SQRT(USTAR*Z0BASE))*Z0BASE,EPSZT)
688 !
689 RZST=ZSLT/ZT
690 RLOGT=LOG(RZST)
691 !
692 !----------------------------------------------------------------------
693 !*** 1./MONIN-OBUKHOV LENGTH-SCALE
694 !----------------------------------------------------------------------
695 !
696 RLMO=ELFC*AKHS*DTHV/USTAR**3
697 ZETALU=ZSLU*RLMO
698 ZETALT=ZSLT*RLMO
699 ZETAU=ZU*RLMO
700 ZETAT=ZT*RLMO
701 !
702 ZETALU=MIN(MAX(ZETALU,ZTMIN2),ZTMAX2)
703 ZETALT=MIN(MAX(ZETALT,ZTMIN2),ZTMAX2)
704 ZETAU=MIN(MAX(ZETAU,ZTMIN2/RZSU),ZTMAX2/RZSU)
705 ZETAT=MIN(MAX(ZETAT,ZTMIN2/RZST),ZTMAX2/RZST)
706 !
707 !----------------------------------------------------------------------
708 !*** LAND FUNCTIONS
709 !----------------------------------------------------------------------
710 !
711 RZ=(ZETAU-ZTMIN2)/DZETA2
712 K=INT(RZ)
713 RDZT=RZ-REAL(K)
714 K=MIN(K,KZTM2)
715 K=MAX(K,0)
716 PSMZ=(PSIM2(K+2)-PSIM2(K+1))*RDZT+PSIM2(K+1)
717 !
718 RZ=(ZETALU-ZTMIN2)/DZETA2
719 K=INT(RZ)
720 RDZT=RZ-REAL(K)
721 K=MIN(K,KZTM2)
722 K=MAX(K,0)
723 PSMZL=(PSIM2(K+2)-PSIM2(K+1))*RDZT+PSIM2(K+1)
724 !
725 SIMM=PSMZL-PSMZ+RLOGU
726 !
727 RZ=(ZETAT-ZTMIN2)/DZETA2
728 K=INT(RZ)
729 RDZT=RZ-REAL(K)
730 K=MIN(K,KZTM2)
731 K=MAX(K,0)
732 PSHZ=(PSIH2(K+2)-PSIH2(K+1))*RDZT+PSIH2(K+1)
733 !
734 RZ=(ZETALT-ZTMIN2)/DZETA2
735 K=INT(RZ)
736 RDZT=RZ-REAL(K)
737 K=MIN(K,KZTM2)
738 K=MAX(K,0)
739 PSHZL=(PSIH2(K+2)-PSIH2(K+1))*RDZT+PSIH2(K+1)
740 !
741 SIMH=(PSHZL-PSHZ+RLOGT)*FH02
742 !----------------------------------------------------------------------
743 USTARK=USTAR*VKARMAN
744 AKMS=MAX(USTARK/SIMM,CXCHL)
745 AKHS=MAX(USTARK/SIMH,CXCHL)
746 !
747 !----------------------------------------------------------------------
748 !*** BELJAARS CORRECTION FOR USTAR
749 !----------------------------------------------------------------------
750 !
751 IF(DTHV<=0.)THEN !zj
752 WSTAR2=WWST2*ABS(BTGH*AKHS*DTHV)**(2./3.) !zj
753 ELSE !zj
754 WSTAR2=0. !zj
755 ENDIF !zj
756 !
757 USTAR=MAX(SQRT(AKMS*SQRT(DU2+WSTAR2)),EPSUST)
758 !
759 !----------------------------------------------------------------------
760 ENDDO land_point_iteration
761 !----------------------------------------------------------------------
762 !
763 ! ENDIF ! End of turbulant branch over land
764 !
765 !----------------------------------------------------------------------
766 !
767 ENDIF ! End of land/sea branch
768 !
769 !----------------------------------------------------------------------
770 !*** COUNTERGRADIENT FIX
771 !----------------------------------------------------------------------
772 !
773 ! HV=-AKHS*DTHV
774 ! IF(HV>0.)THEN
775 ! FCT=-10.*(BTGX)**(-1./3.)
776 ! CT=FCT*(HV/(PBLH*PBLH))**(2./3.)
777 ! ELSE
778 CT=0.
779 ! ENDIF
780 !
781 !----------------------------------------------------------------------
782 !----------------------------------------------------------------------
783 !*** THE FOLLOWING DIAGNOSTIC BLOCK PRODUCES 2-m and 10-m VALUES
784 !*** FOR TEMPERATURE, MOISTURE, AND WINDS. IT IS DONE HERE SINCE
785 !*** THE VARIOUS QUANTITIES NEEDED FOR THE COMPUTATION ARE LOST
786 !*** UPON EXIT FROM THE ROTUINE.
787 !----------------------------------------------------------------------
788 !----------------------------------------------------------------------
789 !
790 WSTAR=SQRT(WSTAR2)/WWST
791 !
792 UMFLX=AKMS*(ULOW -UZ0 )
793 VMFLX=AKMS*(VLOW -VZ0 )
794 HSFLX=AKHS*(THLOW-THZ0)
795 HLFLX=AKHS*(QLOW -QZ0 )
796 !----------------------------------------------------------------------
797 ! IF(RIB>=RIC)THEN
798 !----------------------------------------------------------------------
799 ! IF(SEAMASK>0.5)THEN
800 ! AKMS10=MAX( VISC/10.,CXCHS)
801 ! AKHS02=MAX(TVISC/02.,CXCHS)
802 ! AKHS10=MAX(TVISC/10.,CXCHS)
803 ! ELSE
804 ! AKMS10=MAX( VISC/10.,CXCHL)
805 ! AKHS02=MAX(TVISC/02.,CXCHL)
806 ! AKHS10=MAX(TVISC/10.,CXCHL)
807 ! ENDIF
808 !----------------------------------------------------------------------
809 ! ELSE
810 !----------------------------------------------------------------------
811 ZU10=ZU+10.
812 ZT02=ZT+02.
813 ZT10=ZT+10.
814 !
815 RLNU10=LOG(ZU10/ZU)
816 RLNT02=LOG(ZT02/ZT)
817 RLNT10=LOG(ZT10/ZT)
818 !
819 ZTAU10=ZU10*RLMO
820 ZTAT02=ZT02*RLMO
821 ZTAT10=ZT10*RLMO
822 !
823 !----------------------------------------------------------------------
824 !*** SEA
825 !----------------------------------------------------------------------
826 !
827 IF(SEAMASK>0.5)THEN
828 !
829 !----------------------------------------------------------------------
830 ZTAU10=MIN(MAX(ZTAU10,ZTMIN1),ZTMAX1)
831 ZTAT02=MIN(MAX(ZTAT02,ZTMIN1),ZTMAX1)
832 ZTAT10=MIN(MAX(ZTAT10,ZTMIN1),ZTMAX1)
833 !----------------------------------------------------------------------
834 RZ=(ZTAU10-ZTMIN1)/DZETA1
835 K=INT(RZ)
836 RDZT=RZ-REAL(K)
837 K=MIN(K,KZTM2)
838 K=MAX(K,0)
839 PSM10=(PSIM1(K+2)-PSIM1(K+1))*RDZT+PSIM1(K+1)
840 !
841 SIMM10=PSM10-PSMZ+RLNU10
842 !
843 RZ=(ZTAT02-ZTMIN1)/DZETA1
844 K=INT(RZ)
845 RDZT=RZ-REAL(K)
846 K=MIN(K,KZTM2)
847 K=MAX(K,0)
848 PSH02=(PSIH1(K+2)-PSIH1(K+1))*RDZT+PSIH1(K+1)
849 !
850 SIMH02=(PSH02-PSHZ+RLNT02)*FH01
851 !
852 RZ=(ZTAT10-ZTMIN1)/DZETA1
853 K=INT(RZ)
854 RDZT=RZ-REAL(K)
855 K=MIN(K,KZTM2)
856 K=MAX(K,0)
857 PSH10=(PSIH1(K+2)-PSIH1(K+1))*RDZT+PSIH1(K+1)
858 !
859 SIMH10=(PSH10-PSHZ+RLNT10)*FH01
860 !
861 AKMS10=MAX(USTARK/SIMM10,CXCHS)
862 AKHS02=MAX(USTARK/SIMH02,CXCHS)
863 AKHS10=MAX(USTARK/SIMH10,CXCHS)
864 !
865 !----------------------------------------------------------------------
866 !*** LAND
867 !----------------------------------------------------------------------
868 !
869 ELSE
870 !
871 !----------------------------------------------------------------------
872 ZTAU10=MIN(MAX(ZTAU10,ZTMIN2),ZTMAX2)
873 ZTAT02=MIN(MAX(ZTAT02,ZTMIN2),ZTMAX2)
874 ZTAT10=MIN(MAX(ZTAT10,ZTMIN2),ZTMAX2)
875 !----------------------------------------------------------------------
876 RZ=(ZTAU10-ZTMIN2)/DZETA2
877 K=INT(RZ)
878 RDZT=RZ-REAL(K)
879 K=MIN(K,KZTM2)
880 K=MAX(K,0)
881 PSM10=(PSIM2(K+2)-PSIM2(K+1))*RDZT+PSIM2(K+1)
882 !
883 SIMM10=PSM10-PSMZ+RLNU10
884 !
885 RZ=(ZTAT02-ZTMIN2)/DZETA2
886 K=INT(RZ)
887 RDZT=RZ-REAL(K)
888 K=MIN(K,KZTM2)
889 K=MAX(K,0)
890 PSH02=(PSIH2(K+2)-PSIH2(K+1))*RDZT+PSIH2(K+1)
891 !
892 SIMH02=(PSH02-PSHZ+RLNT02)*FH02
893 !
894 RZ=(ZTAT10-ZTMIN2)/DZETA2
895 K=INT(RZ)
896 RDZT=RZ-REAL(K)
897 K=MIN(K,KZTM2)
898 K=MAX(K,0)
899 PSH10=(PSIH2(K+2)-PSIH2(K+1))*RDZT+PSIH2(K+1)
900 !
901 SIMH10=(PSH10-PSHZ+RLNT10)*FH02
902 !
903 AKMS10=MAX(USTARK/SIMM10,CXCHL)
904 AKHS02=MAX(USTARK/SIMH02,CXCHL)
905 AKHS10=MAX(USTARK/SIMH10,CXCHL)
906 !----------------------------------------------------------------------
907 ENDIF
908 !----------------------------------------------------------------------
909 ! ENDIF
910 !----------------------------------------------------------------------
911 U10 =UMFLX/AKMS10+UZ0
912 V10 =VMFLX/AKMS10+VZ0
913 TH02=HSFLX/AKHS02+THZ0
914 !
915 !*** BE CERTAIN THAT THE 2-M THETA AND 10-M THETA ARE BRACKETED BY
916 !*** THE VALUES OF THZ0 AND THLOW.
917 !
918 IF(THLOW>THZ0.AND.(TH02<THZ0.OR.TH02>THLOW).OR. &
919 THLOW<THZ0.AND.(TH02>THZ0.OR.TH02<THLOW))THEN
920 TH02=THZ0+2.*RDZ*(THLOW-THZ0)
921 ENDIF
922 !
923 TH10=HSFLX/AKHS10+THZ0
924 !
925 IF(THLOW>THZ0.AND.(TH10<THZ0.OR.TH10>THLOW).OR. &
926 THLOW<THZ0.AND.(TH10>THZ0.OR.TH10<THLOW))THEN
927 TH10=THZ0+10.*RDZ*(THLOW-THZ0)
928 ENDIF
929 !
930 Q02 =HLFLX/AKHS02+QZ0
931 Q10 =HLFLX/AKHS10+QZ0
932 TERM1=-0.068283/TLOW
933 PSHLTR=PSFC*EXP(TERM1)
934 !
935 !----------------------------------------------------------------------
936 !*** COMPUTE "EQUIVALENT" Z0 TO APPROXIMATE LOCAL SHELTER READINGS.
937 !----------------------------------------------------------------------
938 !
939 U10E=U10
940 V10E=V10
941 !
942 IF(SEAMASK<0.5)THEN
943
944 !1st ZUUZ=MIN(0.5*ZU,0.1)
945 !1st ZU=MAX(0.1*ZU,ZUUZ)
946 !tst ZUUZ=amin1(ZU*0.50,0.3)
947 !tst ZU=amax1(ZU*0.3,ZUUZ)
948
949 ZUUZ=AMIN1(ZU*0.50,0.18)
950 ZU=AMAX1(ZU*0.35,ZUUZ)
951 !
952 ZU10=ZU+10.
953 RZSU=ZU10/ZU
954 RLNU10=LOG(RZSU)
955
956 ZETAU=ZU*RLMO
957 ZTAU10=ZU10*RLMO
958
959 ZTAU10=MIN(MAX(ZTAU10,ZTMIN2),ZTMAX2)
960 ZETAU=MIN(MAX(ZETAU,ZTMIN2/RZSU),ZTMAX2/RZSU)
961
962 RZ=(ZTAU10-ZTMIN2)/DZETA2
963 K=INT(RZ)
964 RDZT=RZ-REAL(K)
965 K=MIN(K,KZTM2)
966 K=MAX(K,0)
967 PSM10=(PSIM2(K+2)-PSIM2(K+1))*RDZT+PSIM2(K+1)
968 SIMM10=PSM10-PSMZ+RLNU10
969 EKMS10=MAX(USTARK/SIMM10,CXCHL)
970
971 U10E=UMFLX/EKMS10+UZ0
972 V10E=VMFLX/EKMS10+VZ0
973
974 ENDIF
975 !
976 U10=U10E
977 V10=V10E
978 !
979 !----------------------------------------------------------------------
980 !*** SET OTHER WRF DRIVER ARRAYS
981 !----------------------------------------------------------------------
982 !
983 RLOW=PLOW/(R_D*TLOW)
984 CHS=AKHS
985 CHS2=AKHS02
986 CQS2=AKHS02
987 HFX=-RLOW*CP*HSFLX
988 QFX=-RLOW*HLFLX*WETM
989 FLX_LH=XLV*QFX
990 FLHC=RLOW*CP*AKHS
991 FLQC=RLOW*AKHS*WETM
992 !!! QGH=PQ0/PSHLTR*EXP(A2S*(TSK-A3S)/(TSK-A4S))
993 QGH=((1.-SEAMASK)*PQ0+SEAMASK*PQ0SEA) &
994 & /PLOW*EXP(A2S*(TLOW-A3S)/(TLOW-A4S))
995 QGH=QGH/(1.-QGH) !Convert to mixing ratio
996 CPM=CP*(1.+0.8*QLOW)
997 !
998 !*** DO NOT COMPUTE QS OVER LAND POINTS HERE SINCE IT IS
999 !*** A PROGNOSTIC VARIABLE THERE. IT IS OKAY TO USE IT
1000 !*** AS A DIAGNOSTIC OVER WATER SINCE IT WILL CAUSE NO
1001 !*** INTERFERENCE BEFORE BEING RECOMPUTED IN MYJPBL.
1002 !
1003 IF(SEAMASK>0.5)THEN
1004 QS=QLOW+QFX/(RLOW*AKHS)
1005 QS=QS/(1.-QS)
1006 ENDIF
1007 !----------------------------------------------------------------------
1008 !
1009 END SUBROUTINE SFCDIF
1010 !
1011 !----------------------------------------------------------------------
1012 SUBROUTINE MYJSFCINIT(LOWLYR,USTAR,Z0 &
1013 & ,SEAMASK,XICE,IVGTYP,RESTART &
1014 & ,ALLOWED_TO_READ &
1015 & ,IDS,IDE,JDS,JDE,KDS,KDE &
1016 & ,IMS,IME,JMS,JME,KMS,KME &
1017 & ,ITS,ITE,JTS,JTE,KTS,KTE)
1018 !----------------------------------------------------------------------
1019 IMPLICIT NONE
1020 !----------------------------------------------------------------------
1021 LOGICAL,INTENT(IN) :: RESTART,ALLOWED_TO_READ
1022 !
1023 INTEGER,INTENT(IN) :: IDS,IDE,JDS,JDE,KDS,KDE &
1024 & ,IMS,IME,JMS,JME,KMS,KME &
1025 & ,ITS,ITE,JTS,JTE,KTS,KTE
1026 !
1027 INTEGER,DIMENSION(IMS:IME,JMS:JME),INTENT(IN) :: IVGTYP
1028 !
1029 INTEGER,DIMENSION(IMS:IME,JMS:JME),INTENT(INOUT) :: LOWLYR
1030 !
1031 REAL,DIMENSION(IMS:IME,JMS:JME),INTENT(IN) :: SEAMASK,XICE
1032 !
1033 REAL,DIMENSION(IMS:IME,JMS:JME),INTENT(INOUT) :: USTAR,Z0
1034 !
1035 REAL,DIMENSION(0:30) :: VZ0TBL
1036 REAL,DIMENSION(0:30) :: VZ0TBL_24
1037 !
1038 INTEGER :: I,IDUM,IRECV,J,JDUM,K,ITF,JTF,KTF,MAXGBL_IVGTYP &
1039 &, MAXLOC_IVGTYP,MPI_COMM_COMP
1040 !
1041 ! INTEGER :: MPI_INTEGER,MPI_MAX
1042 !
1043 REAL :: SM,X,ZETA1,ZETA2,ZRNG1,ZRNG2
1044 !
1045 REAL :: PIHF=3.1415926/2.,EPS=1.E-6
1046 !----------------------------------------------------------------------
1047 VZ0TBL= &
1048 & (/0., &
1049 & 2.653,0.826,0.563,1.089,0.854,0.856,0.035,0.238,0.065,0.076 &
1050 & ,0.011,0.035,0.011,0.000,0.000,0.000,0.000,0.000,0.000,0.000 &
1051 & ,0.000,0.000,0.000,0.000,0.000,0.000,0.000,0.000,0.000,0.000/)
1052
1053 VZ0TBL_24= (/0., &
1054 & 1.00, 0.07, 0.07, 0.07, 0.07, 0.15, &
1055 & 0.08, 0.03, 0.05, 0.86, 0.80, 0.85, &
1056 & 2.65, 1.09, 0.80, 0.001, 0.04, 0.05, &
1057 & 0.01, 0.04, 0.06, 0.05, 0.03, 0.001, &
1058 & 0.000, 0.000, 0.000, 0.000, 0.000, 0.000 /)
1059
1060 !----------------------------------------------------------------------
1061 !
1062 JTF=MIN0(JTE,JDE-1)
1063 KTF=MIN0(KTE,KDE-1)
1064 ITF=MIN0(ITE,IDE-1)
1065 !
1066 !
1067 !*** FOR NOW, ASSUME SIGMA MODE FOR LOWEST MODEL LAYER
1068 !
1069 DO J=JTS,JTF
1070 DO I=ITS,ITF
1071 LOWLYR(I,J)=1
1072 ! USTAR(I,J)=EPSUST
1073 ENDDO
1074 ENDDO
1075 !----------------------------------------------------------------------
1076 #if (NMM_CORE == 1)
1077 !
1078 IF(.NOT.RESTART)THEN
1079 ! CALL WRF_GET_DM_COMMUNICATOR(MPI_COMM_COMP)
1080 ! MAXLOC_IVGTYP=MAXVAL(IVGTYP)
1081 ! CALL MPI_ALLREDUCE(MAXLOC_IVGTYP,MAXGBL_IVGTYP,1,MPI_INTEGER &
1082 ! &, MPI_MAX,MPI_COMM_COMP,IRECV)
1083 MAXGBL_IVGTYP=MAXVAL(IVGTYP)
1084 CALL WRF_DM_MAXVAL(MAXGBL_IVGTYP,IDUM,JDUM)
1085 !
1086 IF (MAXGBL_IVGTYP<13) THEN
1087 DO J=JTS,JTE
1088 DO I=ITS,ITE
1089 SM=SEAMASK(I,J)-1.
1090 IF(SM+XICE(I,J)<0.5)THEN
1091 Z0(I,J)=VZ0TBL(IVGTYP(I,J))
1092 ENDIF
1093 ENDDO
1094 ENDDO
1095 !
1096 ELSE
1097 !
1098 DO J=JTS,JTE
1099 DO I=ITS,ITE
1100 SM=SEAMASK(I,J)-1.
1101 IF(SM+XICE(I,J)<0.5)THEN
1102 Z0(I,J)=VZ0TBL_24(IVGTYP(I,J))
1103 ENDIF
1104 ENDDO
1105 ENDDO
1106 !
1107 ENDIF ! Vegtype check
1108 !
1109 ENDIF ! Restart check
1110
1111 #endif
1112 !----------------------------------------------------------------------
1113 IF(.NOT.RESTART)THEN
1114 DO J=JTS,JTE
1115 DO I=ITS,ITF
1116 USTAR(I,J)=0.1
1117 ENDDO
1118 ENDDO
1119 ENDIF
1120 !----------------------------------------------------------------------
1121 !
1122 !*** COMPUTE SURFACE LAYER INTEGRAL FUNCTIONS
1123 !
1124 !----------------------------------------------------------------------
1125 FH01=1.
1126 FH02=1.
1127 !
1128 ! ZTMIN1=-10.0
1129 ! ZTMAX1=2.0
1130 ! ZTMIN2=-10.0
1131 ! ZTMAX2=2.0
1132 ZTMIN1=-5.0
1133 ZTMAX1=1.0
1134 ZTMIN2=-5.0
1135 ZTMAX2=1.0
1136 !
1137 ZRNG1=ZTMAX1-ZTMIN1
1138 ZRNG2=ZTMAX2-ZTMIN2
1139 !
1140 DZETA1=ZRNG1/(KZTM-1)
1141 DZETA2=ZRNG2/(KZTM-1)
1142 !
1143 !----------------------------------------------------------------------
1144 !*** FUNCTION DEFINITION LOOP
1145 !----------------------------------------------------------------------
1146 !
1147 ZETA1=ZTMIN1
1148 ZETA2=ZTMIN2
1149 !
1150 DO K=1,KZTM
1151 !
1152 !----------------------------------------------------------------------
1153 !*** UNSTABLE RANGE
1154 !----------------------------------------------------------------------
1155 !
1156 IF(ZETA1<0.)THEN
1157 !
1158 !----------------------------------------------------------------------
1159 !*** PAULSON 1970 FUNCTIONS
1160 !----------------------------------------------------------------------
1161 X=SQRT(SQRT(1.-16.*ZETA1))
1162 !
1163 PSIM1(K)=-2.*LOG((X+1.)/2.)-LOG((X*X+1.)/2.)+2.*ATAN(X)-PIHF
1164 PSIH1(K)=-2.*LOG((X*X+1.)/2.)
1165 !
1166 !----------------------------------------------------------------------
1167 !*** STABLE RANGE
1168 !----------------------------------------------------------------------
1169 !
1170 ELSE
1171 !
1172 !----------------------------------------------------------------------
1173 !*** PAULSON 1970 FUNCTIONS
1174 !----------------------------------------------------------------------
1175 !
1176 ! PSIM1(K)=5.*ZETA1
1177 ! PSIH1(K)=5.*ZETA1
1178 !----------------------------------------------------------------------
1179 !*** HOLTSLAG AND DE BRUIN 1988
1180 !----------------------------------------------------------------------
1181 !
1182 PSIM1(K)=0.7*ZETA1+0.75*ZETA1*(6.-0.35*ZETA1)*EXP(-0.35*ZETA1)
1183 PSIH1(K)=0.7*ZETA1+0.75*ZETA1*(6.-0.35*ZETA1)*EXP(-0.35*ZETA1)
1184 !----------------------------------------------------------------------
1185 !
1186 ENDIF
1187 !
1188 !----------------------------------------------------------------------
1189 !*** UNSTABLE RANGE
1190 !----------------------------------------------------------------------
1191 !
1192 IF(ZETA2<0.)THEN
1193 !
1194 !----------------------------------------------------------------------
1195 !*** PAULSON 1970 FUNCTIONS
1196 !----------------------------------------------------------------------
1197 !
1198 X=SQRT(SQRT(1.-16.*ZETA2))
1199 !
1200 PSIM2(K)=-2.*LOG((X+1.)/2.)-LOG((X*X+1.)/2.)+2.*ATAN(X)-PIHF
1201 PSIH2(K)=-2.*LOG((X*X+1.)/2.)
1202 !----------------------------------------------------------------------
1203 !*** STABLE RANGE
1204 !----------------------------------------------------------------------
1205 !
1206 ELSE
1207 !
1208 !----------------------------------------------------------------------
1209 !*** PAULSON 1970 FUNCTIONS
1210 !----------------------------------------------------------------------
1211 !
1212 ! PSIM2(K)=5.*ZETA2
1213 ! PSIH2(K)=5.*ZETA2
1214 !
1215 !----------------------------------------------------------------------
1216 !*** HOLTSLAG AND DE BRUIN 1988
1217 !----------------------------------------------------------------------
1218 !
1219 PSIM2(K)=0.7*ZETA2+0.75*ZETA2*(6.-0.35*ZETA2)*EXP(-0.35*ZETA2)
1220 PSIH2(K)=0.7*ZETA2+0.75*ZETA2*(6.-0.35*ZETA2)*EXP(-0.35*ZETA2)
1221 !----------------------------------------------------------------------
1222 !
1223 ENDIF
1224 !
1225 !----------------------------------------------------------------------
1226 IF(K==KZTM)THEN
1227 ZTMAX1=ZETA1
1228 ZTMAX2=ZETA2
1229 ENDIF
1230 !
1231 ZETA1=ZETA1+DZETA1
1232 ZETA2=ZETA2+DZETA2
1233 !----------------------------------------------------------------------
1234 ENDDO
1235 !----------------------------------------------------------------------
1236 ZTMAX1=ZTMAX1-EPS
1237 ZTMAX2=ZTMAX2-EPS
1238 !----------------------------------------------------------------------
1239 !
1240 END SUBROUTINE MYJSFCINIT
1241 !
1242 !----------------------------------------------------------------------
1243 !
1244 END MODULE MODULE_SF_MYJSFC
1245 !
1246 !----------------------------------------------------------------------