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