1 | ! |
---|
2 | MODULE MODULE_BL_MYJURB |
---|
3 | ! |
---|
4 | !----------------------------------------------------------------------- |
---|
5 | ! |
---|
6 | USE MODULE_MODEL_CONSTANTS |
---|
7 | ! |
---|
8 | !----------------------------------------------------------------------- |
---|
9 | ! |
---|
10 | ! REFERENCES: Janjic (2002), NCEP Office Note 437 |
---|
11 | ! Mellor and Yamada (1982), Rev. Geophys. Space Phys. |
---|
12 | ! |
---|
13 | ! ABSTRACT: |
---|
14 | ! MYJ UPDATES THE TURBULENT KINETIC ENERGY WITH THE PRODUCTION/ |
---|
15 | ! DISSIPATION TERM AND THE VERTICAL DIFFUSION TERM |
---|
16 | ! (USING AN IMPLICIT FORMULATION) FROM MELLOR-YAMADA |
---|
17 | ! LEVEL 2.5 AS EXTENDED BY JANJIC. EXCHANGE COEFFICIENTS FOR |
---|
18 | ! THE SURFACE AND FOR ALL LAYER INTERFACES ARE COMPUTED FROM |
---|
19 | ! MONIN-OBUKHOV THEORY. |
---|
20 | ! THE TURBULENT VERTICAL EXCHANGE IS THEN EXECUTED. |
---|
21 | ! |
---|
22 | !----------------------------------------------------------------------- |
---|
23 | ! |
---|
24 | INTEGER :: ITRMX=5 ! Iteration count for mixing length computation |
---|
25 | ! |
---|
26 | ! REAL,PARAMETER :: G=9.81,PI=3.1415926,R_D=287.04,R_V=461.6 & |
---|
27 | ! & ,VKARMAN=0.4 |
---|
28 | REAL,PARAMETER :: PI=3.1415926,VKARMAN=0.4 |
---|
29 | ! REAL,PARAMETER :: CP=7.*R_D/2. |
---|
30 | REAL,PARAMETER :: CAPA=R_D/CP |
---|
31 | REAL,PARAMETER :: RLIVWV=XLS/XLV,ELOCP=2.72E6/CP |
---|
32 | REAL,PARAMETER :: EPS1=1.E-12,EPS2=0. |
---|
33 | REAL,PARAMETER :: EPSL=0.32,EPSRU=1.E-7,EPSRS=1.E-7 & |
---|
34 | & ,EPSTRB=1.E-24 |
---|
35 | REAL,PARAMETER :: EPSA=1.E-8,EPSIT=1.E-4,EPSU2=1.E-4,EPSUST=0.07 & |
---|
36 | & ,FH=1.01 |
---|
37 | REAL,PARAMETER :: ALPH=0.30,BETA=1./273.,EL0MAX=1000.,EL0MIN=1. & |
---|
38 | & ,ELFC=0.23*0.5,GAM1=0.2222222222222222222 & |
---|
39 | & ,PRT=1. |
---|
40 | REAL,PARAMETER :: A1=0.659888514560862645 & |
---|
41 | & ,A2x=0.6574209922667784586 & |
---|
42 | & ,B1=11.87799326209552761 & |
---|
43 | & ,B2=7.226971804046074028 & |
---|
44 | & ,C1=0.000830955950095854396 |
---|
45 | REAL,PARAMETER :: A2S=17.2693882,A3S=273.16,A4S=35.86 |
---|
46 | REAL,PARAMETER :: ELZ0=0.,ESQ=5.0,EXCM=0.001 & |
---|
47 | & ,FHNEU=0.8,GLKBR=10.,GLKBS=30. & |
---|
48 | & ,QVISC=2.1E-5,RFC=0.191,RIC=0.505,SMALL=0.35 & |
---|
49 | & ,SQPR=0.84,SQSC=0.84,SQVISC=258.2,TVISC=2.1E-5 & |
---|
50 | & ,USTC=0.7,USTR=0.225,VISC=1.5E-5 & |
---|
51 | & ,WOLD=0.15,WWST=1.2,ZTMAX=1.,ZTFC=1.,ZTMIN=-5. |
---|
52 | ! |
---|
53 | REAL,PARAMETER :: SEAFC=0.98,PQ0SEA=PQ0*SEAFC |
---|
54 | ! |
---|
55 | REAL,PARAMETER :: BTG=BETA*G,CZIV=SMALL*GLKBS & |
---|
56 | ! & ,EP_1=R_V/R_D-1.,ESQHF=0.5*5.0,GRRS=GLKBR/GLKBS & |
---|
57 | & ,ESQHF=0.5*5.0,GRRS=GLKBR/GLKBS & |
---|
58 | & ,RB1=1./B1,RTVISC=1./TVISC,RVISC=1./VISC & |
---|
59 | & ,ZQRZT=SQSC/SQPR |
---|
60 | ! |
---|
61 | REAL,PARAMETER :: ADNH= 9.*A1*A2x*A2x*(12.*A1+3.*B2)*BTG*BTG & |
---|
62 | & ,ADNM=18.*A1*A1*A2x*(B2-3.*A2x)*BTG & |
---|
63 | & ,ANMH=-9.*A1*A2x*A2x*BTG*BTG & |
---|
64 | & ,ANMM=-3.*A1*A2x*(3.*A2x+3.*B2*C1+18.*A1*C1-B2) & |
---|
65 | & *BTG & |
---|
66 | & ,BDNH= 3.*A2x*(7.*A1+B2)*BTG & |
---|
67 | & ,BDNM= 6.*A1*A1 & |
---|
68 | & ,BEQH= A2x*B1*BTG+3.*A2x*(7.*A1+B2)*BTG & |
---|
69 | & ,BEQM=-A1*B1*(1.-3.*C1)+6.*A1*A1 & |
---|
70 | & ,BNMH=-A2x*BTG & |
---|
71 | & ,BNMM=A1*(1.-3.*C1) & |
---|
72 | & ,BSHH=9.*A1*A2x*A2x*BTG & |
---|
73 | & ,BSHM=18.*A1*A1*A2x*C1 & |
---|
74 | & ,BSMH=-3.*A1*A2x*(3.*A2x+3.*B2*C1+12.*A1*C1-B2) & |
---|
75 | & *BTG & |
---|
76 | & ,CESH=A2x & |
---|
77 | & ,CESM=A1*(1.-3.*C1) & |
---|
78 | & ,CNV=EP_1*G/BTG & |
---|
79 | & ,ELFCS=VKARMAN*BTG & |
---|
80 | & ,FZQ1=RTVISC*QVISC*ZQRZT & |
---|
81 | & ,FZQ2=RTVISC*QVISC*ZQRZT & |
---|
82 | & ,FZT1=RVISC *TVISC*SQPR & |
---|
83 | & ,FZT2=CZIV*GRRS*TVISC*SQPR & |
---|
84 | & ,FZU1=CZIV*VISC & |
---|
85 | & ,PIHF=0.5*PI & |
---|
86 | & ,RFAC=RIC/(FHNEU*RFC*RFC) & |
---|
87 | & ,RQVISC=1./QVISC & |
---|
88 | & ,RRIC=1./RIC & |
---|
89 | & ,USTFC=0.018/G & |
---|
90 | & ,WNEW=1.-WOLD & |
---|
91 | & ,WWST2=WWST*WWST |
---|
92 | ! |
---|
93 | !----------------------------------------------------------------------- |
---|
94 | !*** FREE TERM IN THE EQUILIBRIUM EQUATION FOR (L/Q)**2 |
---|
95 | !----------------------------------------------------------------------- |
---|
96 | ! |
---|
97 | REAL,PARAMETER :: AEQH=9.*A1*A2x*A2x*B1*BTG*BTG & |
---|
98 | & +9.*A1*A2x*A2x*(12.*A1+3.*B2)*BTG*BTG & |
---|
99 | & ,AEQM=3.*A1*A2x*B1*(3.*A2x+3.*B2*C1+18.*A1*C1-B2)& |
---|
100 | & *BTG+18.*A1*A1*A2x*(B2-3.*A2x)*BTG |
---|
101 | ! |
---|
102 | !----------------------------------------------------------------------- |
---|
103 | !*** FORBIDDEN TURBULENCE AREA |
---|
104 | !----------------------------------------------------------------------- |
---|
105 | ! |
---|
106 | REAL,PARAMETER :: REQU=-AEQH/AEQM & |
---|
107 | & ,EPSGH=1.E-9,EPSGM=REQU*EPSGH |
---|
108 | ! |
---|
109 | !----------------------------------------------------------------------- |
---|
110 | !*** NEAR ISOTROPY FOR SHEAR TURBULENCE, WW/Q2 LOWER LIMIT |
---|
111 | !----------------------------------------------------------------------- |
---|
112 | ! |
---|
113 | REAL,PARAMETER :: UBRYL=(18.*REQU*A1*A1*A2x*B2*C1*BTG & |
---|
114 | & +9.*A1*A2x*A2x*B2*BTG*BTG) & |
---|
115 | & /(REQU*ADNM+ADNH) & |
---|
116 | & ,UBRY=(1.+EPSRS)*UBRYL,UBRY3=3.*UBRY |
---|
117 | ! |
---|
118 | REAL,PARAMETER :: AUBH=27.*A1*A2x*A2x*B2*BTG*BTG-ADNH*UBRY3 & |
---|
119 | & ,AUBM=54.*A1*A1*A2x*B2*C1*BTG -ADNM*UBRY3 & |
---|
120 | & ,BUBH=(9.*A1*A2x+3.*A2x*B2)*BTG-BDNH*UBRY3 & |
---|
121 | & ,BUBM=18.*A1*A1*C1 -BDNM*UBRY3 & |
---|
122 | & ,CUBR=1. - UBRY3 & |
---|
123 | & ,RCUBR=1./CUBR |
---|
124 | ! |
---|
125 | !----------------------------------------------------------------------- |
---|
126 | ! |
---|
127 | CONTAINS |
---|
128 | ! |
---|
129 | !---------------------------------------------------------------------- |
---|
130 | SUBROUTINE MYJURB(IDIFF,FLAG_BEP,DT,STEPBL,HT,DZ & |
---|
131 | & ,PMID,PINT,TH,T,EXNER,QV,CWM,U,V,RHO & |
---|
132 | & ,TSK,QSFC,CHKLOWQ,THZ0,QZ0,UZ0,VZ0 & |
---|
133 | & ,LOWLYR,XLAND,SICE,SNOW & |
---|
134 | & ,TKE_MYJ,EXCH_H,EXCH_M,USTAR,ZNT,EL_MYJ,PBLH,KPBL,CT & |
---|
135 | & ,AKHS,AKMS,ELFLX & |
---|
136 | & ,RUBLTEN,RVBLTEN,RTHBLTEN,RQVBLTEN,RQCBLTEN & |
---|
137 | ! Variables added for multi-layer UCM |
---|
138 | & ,FRC_URB2D & |
---|
139 | & ,A_U_BEP,A_V_BEP,A_T_BEP,A_Q_BEP & |
---|
140 | & ,A_E_BEP,B_U_BEP,B_V_BEP & |
---|
141 | & ,B_T_BEP,B_Q_BEP,B_E_BEP,DLG_BEP & |
---|
142 | & ,DL_U_BEP,SF_BEP,VL_BEP & |
---|
143 | ! Multi-layer UCM end |
---|
144 | & ,IDS,IDE,JDS,JDE,KDS,KDE & |
---|
145 | & ,IMS,IME,JMS,JME,KMS,KME & |
---|
146 | & ,ITS,ITE,JTS,JTE,KTS,KTE) |
---|
147 | !---------------------------------------------------------------------- |
---|
148 | ! |
---|
149 | IMPLICIT NONE |
---|
150 | ! |
---|
151 | !---------------------------------------------------------------------- |
---|
152 | INTEGER,INTENT(IN) :: IDIFF |
---|
153 | INTEGER,INTENT(IN) :: IDS,IDE,JDS,JDE,KDS,KDE & |
---|
154 | & ,IMS,IME,JMS,JME,KMS,KME & |
---|
155 | & ,ITS,ITE,JTS,JTE,KTS,KTE |
---|
156 | ! |
---|
157 | INTEGER,INTENT(IN) :: STEPBL |
---|
158 | |
---|
159 | INTEGER,DIMENSION(IMS:IME,JMS:JME),INTENT(IN) :: LOWLYR |
---|
160 | ! |
---|
161 | INTEGER,DIMENSION(IMS:IME,JMS:JME),INTENT(OUT) :: KPBL |
---|
162 | ! |
---|
163 | REAL,INTENT(IN) :: DT |
---|
164 | ! |
---|
165 | REAL,DIMENSION(IMS:IME,JMS:JME),INTENT(IN) :: HT,SICE,SNOW & |
---|
166 | & ,TSK,XLAND |
---|
167 | ! |
---|
168 | REAL,DIMENSION(IMS:IME,KMS:KME,JMS:JME),INTENT(IN) :: CWM,DZ & |
---|
169 | & ,EXNER & |
---|
170 | & ,PMID,PINT & |
---|
171 | & ,QV,RHO & |
---|
172 | & ,T,TH,U,V |
---|
173 | ! |
---|
174 | REAL,DIMENSION(IMS:IME,JMS:JME),INTENT(OUT) :: PBLH |
---|
175 | ! |
---|
176 | REAL,DIMENSION(IMS:IME,JMS:JME),INTENT(INOUT) :: AKHS,AKMS |
---|
177 | ! |
---|
178 | REAL,DIMENSION(IMS:IME,KMS:KME,JMS:JME) & |
---|
179 | & ,INTENT(OUT) :: EL_MYJ & |
---|
180 | & ,RQCBLTEN,RQVBLTEN & |
---|
181 | & ,RTHBLTEN & |
---|
182 | & ,RUBLTEN,RVBLTEN |
---|
183 | ! Variables added for multi-layer UCM |
---|
184 | REAL,DIMENSION(IMS:IME,KMS:KME,JMS:JME),INTENT(IN) :: A_U_BEP & |
---|
185 | & ,A_V_BEP,A_T_BEP & |
---|
186 | & ,A_E_BEP,B_U_BEP & |
---|
187 | & ,A_Q_BEP,B_Q_BEP & |
---|
188 | & ,B_V_BEP,B_T_BEP & |
---|
189 | & ,B_E_BEP,DLG_BEP & |
---|
190 | & ,DL_U_BEP & |
---|
191 | & ,VL_BEP,SF_BEP |
---|
192 | REAL,DIMENSION(IMS:IME,JMS:JME),INTENT(IN) :: FRC_URB2D |
---|
193 | REAL,DIMENSION(IMS:IME,KMS:KME,JMS:JME) & |
---|
194 | & ,INTENT(INOUT) :: EXCH_H,EXCH_M,TKE_MYJ |
---|
195 | ! UCM end |
---|
196 | |
---|
197 | REAL,DIMENSION(IMS:IME,JMS:JME),INTENT(INOUT) :: CT,QSFC,QZ0 & |
---|
198 | & ,THZ0,USTAR & |
---|
199 | & ,UZ0,VZ0,ZNT |
---|
200 | ! |
---|
201 | REAL,DIMENSION(IMS:IME,JMS:JME),INTENT(IN) :: CHKLOWQ,ELFLX |
---|
202 | LOGICAL, INTENT(IN ) :: FLAG_BEP |
---|
203 | ! |
---|
204 | !---------------------------------------------------------------------- |
---|
205 | !*** |
---|
206 | !*** LOCAL VARIABLES |
---|
207 | !*** |
---|
208 | INTEGER :: I,J,K,KFLIP,LLOW,LMH,LMXL |
---|
209 | ! Variables added for multi-layer UCM |
---|
210 | REAL,DIMENSION(KTS:KTE) :: A_U,A_V,A_T,A_E,B_U,B_V,B_T,B_E !Vertically flipped fields |
---|
211 | REAL,DIMENSION(KTS:KTE) :: A_Q,B_Q,VLK,DL_U,B_E_FACE,DZ_KFLIP !Vertically flipped fields |
---|
212 | REAL,DIMENSION(KTS:KTE) :: SFK !Vertically flipped SF,VL |
---|
213 | INTEGER :: IURB |
---|
214 | REAL :: SLOPE,ICEPT |
---|
215 | ! UCM end |
---|
216 | |
---|
217 | ! |
---|
218 | INTEGER,DIMENSION(ITS:ITE,JTS:JTE) :: LPBL |
---|
219 | ! |
---|
220 | REAL :: AKHS_DENS,AKMS_DENS,APEX,DCDT,DELTAZ,DQDT,DTDIF,DTDT & |
---|
221 | & ,DTTURBL,DUDT,DVDT,EXNSFC,PSFC,PTOP,QFC1,QLOW,QOLD & |
---|
222 | & ,RATIOMX,RDTTURBL,RG,RWMSK,SEAMASK,THNEW,THOLD,TX & |
---|
223 | & ,ULOW,VLOW,WMSK |
---|
224 | ! |
---|
225 | REAL,DIMENSION(KTS:KTE) :: CWMK,PK,Q2K,QK,THEK,TK,UK,VK |
---|
226 | ! |
---|
227 | REAL,DIMENSION(KTS:KTE-1) :: AKHK,AKMK,EL,GH,GM |
---|
228 | ! |
---|
229 | REAL,DIMENSION(KTS:KTE+1) :: ZHK |
---|
230 | ! |
---|
231 | REAL,DIMENSION(ITS:ITE,JTS:JTE) :: THSK |
---|
232 | ! |
---|
233 | REAL,DIMENSION(KTS:KTE,ITS:ITE) :: RHOK |
---|
234 | ! |
---|
235 | REAL,DIMENSION(ITS:ITE,KTS:KTE,JTS:JTE) :: APE,THE |
---|
236 | ! |
---|
237 | REAL,DIMENSION(ITS:ITE,KTS:KTE-1,JTS:JTE) :: AKH,AKM |
---|
238 | ! |
---|
239 | REAL,DIMENSION(ITS:ITE,KTS:KTE+1,JTS:JTE) :: ZINT |
---|
240 | ! |
---|
241 | !*** Begin debugging |
---|
242 | REAL :: ZSL_DIAG |
---|
243 | INTEGER :: IMD,JMD,PRINT_DIAG |
---|
244 | !*** End debugging |
---|
245 | ! |
---|
246 | !---------------------------------------------------------------------- |
---|
247 | !********************************************************************** |
---|
248 | !---------------------------------------------------------------------- |
---|
249 | ! |
---|
250 | !*** Begin debugging |
---|
251 | IMD=(IMS+IME)/2 |
---|
252 | JMD=(JMS+JME)/2 |
---|
253 | !*** End debugging |
---|
254 | ! |
---|
255 | !*** MAKE PREPARATIONS |
---|
256 | ! |
---|
257 | |
---|
258 | !---------------------------------------------------------------------- |
---|
259 | DTTURBL=DT*STEPBL |
---|
260 | RDTTURBL=1./DTTURBL |
---|
261 | DTDIF=DTTURBL |
---|
262 | RG=1./G |
---|
263 | ! |
---|
264 | DO J=JTS,JTE |
---|
265 | DO K=KTS,KTE-1 |
---|
266 | DO I=ITS,ITE |
---|
267 | AKM(I,K,J)=0. |
---|
268 | ENDDO |
---|
269 | ENDDO |
---|
270 | ENDDO |
---|
271 | ! |
---|
272 | DO J=JTS,JTE |
---|
273 | DO K=KTS,KTE+1 |
---|
274 | DO I=ITS,ITE |
---|
275 | ZINT(I,K,J)=0. |
---|
276 | ENDDO |
---|
277 | ENDDO |
---|
278 | ENDDO |
---|
279 | ! |
---|
280 | DO J=JTS,JTE |
---|
281 | DO I=ITS,ITE |
---|
282 | ZINT(I,KTE+1,J)=HT(I,J) ! Z at bottom of lowest sigma layer |
---|
283 | ! |
---|
284 | !!!!!!!!! |
---|
285 | !!!!!! UNCOMMENT THESE LINES IF USING ETA COORDINATES |
---|
286 | !!!!!!!!! |
---|
287 | !!!!!! ZINT(I,KTE+1,J)=1.E-4 ! Z of bottom of lowest eta layer |
---|
288 | !!!!!! ZHK(KTE+1)=1.E-4 ! Z of bottom of lowest eta layer |
---|
289 | ! |
---|
290 | ENDDO |
---|
291 | ENDDO |
---|
292 | ! |
---|
293 | DO J=JTS,JTE |
---|
294 | DO K=KTE,KTS,-1 |
---|
295 | KFLIP=KTE+1-K |
---|
296 | DO I=ITS,ITE |
---|
297 | ZINT(I,K,J)=ZINT(I,K+1,J)+DZ(I,KFLIP,J) |
---|
298 | APEX=1./EXNER(I,K,J) |
---|
299 | APE(I,K,J)=APEX |
---|
300 | TX=T(I,K,J) |
---|
301 | THE(I,K,J)=(CWM(I,K,J)*(-ELOCP/TX)+1.)*TH(I,K,J) |
---|
302 | ENDDO |
---|
303 | ENDDO |
---|
304 | ENDDO |
---|
305 | EL_MYJ(its:ite,:,jts:jte) = 0. |
---|
306 | ! |
---|
307 | !---------------------------------------------------------------------- |
---|
308 | setup_integration: DO J=JTS,JTE |
---|
309 | !---------------------------------------------------------------------- |
---|
310 | ! |
---|
311 | DO I=ITS,ITE |
---|
312 | ! |
---|
313 | !*** LOWEST LAYER ABOVE GROUND MUST BE FLIPPED |
---|
314 | ! |
---|
315 | LMH=KTE-LOWLYR(I,J)+1 |
---|
316 | ! |
---|
317 | PTOP=PINT(I,KTE+1,J) ! KTE+1=KME |
---|
318 | PSFC=PINT(I,LOWLYR(I,J),J) |
---|
319 | ! |
---|
320 | !*** CONVERT LAND MASK (1 FOR SEA; 0 FOR LAND) |
---|
321 | ! |
---|
322 | SEAMASK=XLAND(I,J)-1. |
---|
323 | ! |
---|
324 | !*** FILL 1-D VERTICAL ARRAYS |
---|
325 | !*** AND FLIP DIRECTION SINCE MYJ SCHEME |
---|
326 | !*** COUNTS DOWNWARD FROM THE DOMAIN'S TOP |
---|
327 | ! |
---|
328 | DO K=KTE,KTS,-1 |
---|
329 | KFLIP=KTE+1-K |
---|
330 | TK(K)=T(I,KFLIP,J) |
---|
331 | THEK(K)=THE(I,KFLIP,J) |
---|
332 | RATIOMX=QV(I,KFLIP,J) |
---|
333 | QK(K)=RATIOMX/(1.+RATIOMX) |
---|
334 | CWMK(K)=CWM(I,KFLIP,J) |
---|
335 | PK(K)=PMID(I,KFLIP,J) |
---|
336 | UK(K)=U(I,KFLIP,J) |
---|
337 | VK(K)=V(I,KFLIP,J) |
---|
338 | |
---|
339 | DZ_KFLIP(K)=DZ(I,KFLIP,J) |
---|
340 | |
---|
341 | |
---|
342 | |
---|
343 | ! |
---|
344 | !*** TKE=0.5*(q**2) ==> q**2=2.*TKE |
---|
345 | ! |
---|
346 | Q2K(K)=2.*TKE_MYJ(I,KFLIP,J) |
---|
347 | ! |
---|
348 | !*** COMPUTE THE HEIGHTS OF THE LAYER INTERFACES |
---|
349 | ! |
---|
350 | ZHK(K)=ZINT(I,K,J) |
---|
351 | ENDDO |
---|
352 | |
---|
353 | ZHK(KTE+1)=HT(I,J) ! Z at bottom of lowest sigma layer |
---|
354 | ! Multi-layer UCM |
---|
355 | if(flag_bep)then |
---|
356 | DO K=KTE,KTS,-1 |
---|
357 | KFLIP=KTE+1-K |
---|
358 | B_E(K)=B_E_BEP(I,KFLIP,J) |
---|
359 | A_E(K)=A_E_BEP(I,KFLIP,J) |
---|
360 | SFK(K)=SF_BEP(I,KFLIP,J) |
---|
361 | VLK(K)=VL_BEP(I,KFLIP,J) |
---|
362 | B_E_FACE(K)=0. |
---|
363 | DL_U(K)=DL_U_BEP(I,KFLIP,J) |
---|
364 | enddo |
---|
365 | else |
---|
366 | DO K=KTE,KTS,-1 |
---|
367 | B_E(K)=0. |
---|
368 | A_E(K)=0. |
---|
369 | SFK(K)=1. |
---|
370 | VLK(K)=1. |
---|
371 | B_E_FACE(K)=0. |
---|
372 | DL_U(K)=0. |
---|
373 | enddo |
---|
374 | endif |
---|
375 | !UCM end |
---|
376 | ! Multi-layer UCM - caculate B_E at grid faces |
---|
377 | ! SFK(1)=1 |
---|
378 | IF (FLAG_BEP) THEN |
---|
379 | DO K=KTE,KTS,-1 |
---|
380 | IF (K==KTE) THEN |
---|
381 | B_E_FACE(K)=0.0 |
---|
382 | ELSEIF ( K == KTS ) THEN |
---|
383 | B_E_FACE(K) = 0.0 |
---|
384 | ELSE |
---|
385 | SLOPE=(B_E(K-1)-B_E(K))/0.5/(DZ_KFLIP(K-1)+DZ_KFLIP(K)) |
---|
386 | ICEPT=B_E(K)-SLOPE*(ZHK(K)-0.5*DZ_KFLIP(K)) |
---|
387 | B_E_FACE(K)=2.*(SLOPE*ZHK(K)+ICEPT) !*2 accounts for Q2K=2.*TKE_MYJ |
---|
388 | ENDIF |
---|
389 | ENDDO |
---|
390 | ENDIF |
---|
391 | ! |
---|
392 | !*** Begin debugging |
---|
393 | ! IF(I==IMD.AND.J==JMD)THEN |
---|
394 | ! PRINT_DIAG=1 |
---|
395 | ! ELSE |
---|
396 | ! PRINT_DIAG=0 |
---|
397 | ! ENDIF |
---|
398 | ! IF(I==227.AND.J==363)PRINT_DIAG=2 |
---|
399 | !*** End debugging |
---|
400 | ! |
---|
401 | !---------------------------------------------------------------------- |
---|
402 | !*** |
---|
403 | !*** FIND THE MIXING LENGTH |
---|
404 | !*** |
---|
405 | CALL MIXLEN(LMH,UK,VK,TK,THEK,QK,CWMK & |
---|
406 | & ,Q2K,ZHK,GM,GH,EL & |
---|
407 | & ,PBLH(I,J),LPBL(I,J),LMXL,CT(I,J) & |
---|
408 | ! Multi-layer UCM |
---|
409 | ,DL_U & |
---|
410 | ! UCM end |
---|
411 | & ,IDS,IDE,JDS,JDE,KDS,KDE & |
---|
412 | & ,IMS,IME,JMS,JME,KMS,KME & |
---|
413 | & ,ITS,ITE,JTS,JTE,KTS,KTE) |
---|
414 | ! |
---|
415 | !---------------------------------------------------------------------- |
---|
416 | !*** |
---|
417 | !*** SOLVE FOR THE PRODUCTION/DISSIPATION OF |
---|
418 | !*** THE TURBULENT KINETIC ENERGY |
---|
419 | !*** |
---|
420 | ! |
---|
421 | CALL PRODQ2(LMH,DTTURBL,USTAR(I,J),GM,GH,EL,Q2K & |
---|
422 | & ,IDS,IDE,JDS,JDE,KDS,KDE & |
---|
423 | & ,IMS,IME,JMS,JME,KMS,KME & |
---|
424 | & ,ITS,ITE,JTS,JTE,KTS,KTE) |
---|
425 | ! |
---|
426 | ! Multi-layer UCM |
---|
427 | IF (FLAG_BEP) THEN |
---|
428 | Q2K(LMH)=Q2K(LMH)*(1-FRC_URB2D(I,J))+(B_E_FACE(LMH)* & |
---|
429 | & 0.5*EL(LMH-1)/11.788/2.)**(2./3.) |
---|
430 | ENDIF |
---|
431 | !UCM end |
---|
432 | ! |
---|
433 | !---------------------------------------------------------------------- |
---|
434 | !*** THE MODEL LAYER (COUNTING UPWARD) CONTAINING THE TOP OF THE PBL |
---|
435 | !---------------------------------------------------------------------- |
---|
436 | ! |
---|
437 | KPBL(I,J)=KTE-LPBL(I,J)+1 |
---|
438 | ! |
---|
439 | !---------------------------------------------------------------------- |
---|
440 | !*** |
---|
441 | !*** FIND THE EXCHANGE COEFFICIENTS IN THE FREE ATMOSPHERE |
---|
442 | !*** |
---|
443 | CALL DIFCOF(LMH,LMXL,GM,GH,EL,TK,Q2K,ZHK,AKMK,AKHK & |
---|
444 | & ,IDS,IDE,JDS,JDE,KDS,KDE & |
---|
445 | & ,IMS,IME,JMS,JME,KMS,KME & |
---|
446 | & ,ITS,ITE,JTS,JTE,KTS,KTE,PRINT_DIAG) ! debug |
---|
447 | ! |
---|
448 | !*** COUNTING DOWNWARD FROM THE TOP, THE EXCHANGE COEFFICIENTS AKH |
---|
449 | !*** ARE DEFINED ON THE BOTTOMS OF THE LAYERS KTS TO KTE-1. COUNTING |
---|
450 | !*** COUNTING UPWARD FROM THE BOTTOM, THOSE SAME COEFFICIENTS EXCH_H |
---|
451 | !*** ARE DEFINED ON THE TOPS OF THE LAYERS KTS TO KTE-1. |
---|
452 | ! |
---|
453 | DO K=KTS,KTE-1 |
---|
454 | KFLIP=KTE-K |
---|
455 | AKH(I,K,J)=AKHK(K) |
---|
456 | AKM(I,K,J)=AKMK(K) |
---|
457 | DELTAZ=0.5*(ZHK(KFLIP)-ZHK(KFLIP+2)) |
---|
458 | ! Multi-layer UCM |
---|
459 | EXCH_H(I,K+1,J)=AKHK(KFLIP)*DELTAZ |
---|
460 | EXCH_M(I,K+1,J)=AKMK(KFLIP)*DELTAZ |
---|
461 | ! UCM end |
---|
462 | ENDDO |
---|
463 | !---------------------------------------------------------------------- |
---|
464 | !*** |
---|
465 | !*** CARRY OUT THE VERTICAL DIFFUSION OF |
---|
466 | !*** TURBULENT KINETIC ENERGY |
---|
467 | !*** |
---|
468 | ! |
---|
469 | |
---|
470 | CALL VDIFQ(LMH,DTDIF,Q2K,EL,ZHK & |
---|
471 | ! Multi-layer UCM |
---|
472 | & ,A_E,B_E_FACE,SFK,VLK & |
---|
473 | ! UCM end |
---|
474 | & ,IDS,IDE,JDS,JDE,KDS,KDE & |
---|
475 | & ,IMS,IME,JMS,JME,KMS,KME & |
---|
476 | & ,ITS,ITE,JTS,JTE,KTS,KTE) |
---|
477 | ! |
---|
478 | !*** SAVE THE NEW TKE AND MIXING LENGTH. |
---|
479 | ! |
---|
480 | DO K=KTS,KTE |
---|
481 | KFLIP=KTE+1-K |
---|
482 | Q2K(KFLIP)=AMAX1(Q2K(KFLIP),EPSQ2) |
---|
483 | TKE_MYJ(I,K,J)=0.5*Q2K(KFLIP) |
---|
484 | IF(KFLIP/=KTE)EL_MYJ(I,K,J)=EL(KFLIP) ! EL IS NOT DEFINED AT KTE (ground surface) |
---|
485 | ENDDO |
---|
486 | ! |
---|
487 | ENDDO |
---|
488 | ! |
---|
489 | !---------------------------------------------------------------------- |
---|
490 | ENDDO setup_integration |
---|
491 | !---------------------------------------------------------------------- |
---|
492 | ! Multi-layer UCM - if idiff=1 then integration in PBL driver |
---|
493 | IF(IDIFF.NE.1)then |
---|
494 | ! UCM end |
---|
495 | ! |
---|
496 | !*** CONVERT SURFACE SENSIBLE TEMPERATURE TO POTENTIAL TEMPERATURE. |
---|
497 | ! |
---|
498 | DO J=JTS,JTE |
---|
499 | DO I=ITS,ITE |
---|
500 | PSFC=PINT(I,LOWLYR(I,J),J) |
---|
501 | THSK(I,J)=TSK(I,J)*(1.E5/PSFC)**CAPA |
---|
502 | ENDDO |
---|
503 | ENDDO |
---|
504 | ! |
---|
505 | !---------------------------------------------------------------------- |
---|
506 | ! |
---|
507 | !---------------------------------------------------------------------- |
---|
508 | main_integration: DO J=JTS,JTE |
---|
509 | !---------------------------------------------------------------------- |
---|
510 | ! |
---|
511 | DO I=ITS,ITE |
---|
512 | ! |
---|
513 | !*** FILL 1-D VERTICAL ARRAYS |
---|
514 | !*** AND FLIP DIRECTION SINCE MYJ SCHEME |
---|
515 | !*** COUNTS DOWNWARD FROM THE DOMAIN'S TOP |
---|
516 | ! |
---|
517 | DO K=KTE,KTS,-1 |
---|
518 | KFLIP=KTE+1-K |
---|
519 | THEK(K)=THE(I,KFLIP,J) |
---|
520 | RATIOMX=QV(I,KFLIP,J) !mixing ratio |
---|
521 | QK(K)=RATIOMX/(1.+RATIOMX) !specific humidity |
---|
522 | CWMK(K)=CWM(I,KFLIP,J) |
---|
523 | ZHK(K)=ZINT(I,K,J) |
---|
524 | RHOK(K,I)=PMID(I,KFLIP,J)/(R_D*T(I,KFLIP,J)* & |
---|
525 | & (1.+P608*QK(K)-CWMK(K))) |
---|
526 | |
---|
527 | ENDDO |
---|
528 | ! SFK(1)=1. |
---|
529 | if(flag_bep)then |
---|
530 | ! Multi-layer UCM |
---|
531 | |
---|
532 | DO K=KTE,KTS,-1 |
---|
533 | KFLIP=KTE+1-K |
---|
534 | A_Q(K)=A_Q_BEP(I,KFLIP,J) |
---|
535 | B_Q(K)=B_Q_BEP(I,KFLIP,J) |
---|
536 | A_T(K)=A_T_BEP(I,KFLIP,J) |
---|
537 | B_T(K)=B_T_BEP(I,KFLIP,J) |
---|
538 | SFK(K)=SF_BEP(I,KFLIP,J) |
---|
539 | VLK(K)=VL_BEP(I,KFLIP,J) |
---|
540 | ENDDO |
---|
541 | ! UCM end |
---|
542 | else |
---|
543 | DO K=KTE,KTS,-1 |
---|
544 | A_Q(K)=0. |
---|
545 | B_Q(K)=0. |
---|
546 | A_T(K)=0. |
---|
547 | B_T(K)=0. |
---|
548 | SFK(K)=1. |
---|
549 | VLK(K)=1. |
---|
550 | ENDDO |
---|
551 | endif |
---|
552 | ! |
---|
553 | !*** COUNTING DOWNWARD FROM THE TOP, THE EXCHANGE COEFFICIENTS AKH |
---|
554 | !*** ARE DEFINED ON THE BOTTOMS OF THE LAYERS KTS TO KTE-1. THESE COEFFICIENTS |
---|
555 | !*** ARE ALSO MULTIPLIED BY THE DENSITY AT THE BOTTOM INTERFACE LEVEL. |
---|
556 | ! |
---|
557 | DO K=KTS,KTE-1 |
---|
558 | AKHK(K)=AKH(I,K,J)*0.5*(RHOK(K,I)+RHOK(K+1,I)) |
---|
559 | ENDDO |
---|
560 | ! |
---|
561 | ZHK(KTE+1)=ZINT(I,KTE+1,J) |
---|
562 | ! |
---|
563 | SEAMASK=XLAND(I,J)-1. |
---|
564 | THZ0(I,J)=(1.-SEAMASK)*THSK(I,J)+SEAMASK*THZ0(I,J) |
---|
565 | ! |
---|
566 | LLOW=LOWLYR(I,J) |
---|
567 | AKHS_DENS=AKHS(I,J)*RHOK(KTE+1-LLOW,I) |
---|
568 | ! |
---|
569 | IF(SEAMASK<0.5)THEN |
---|
570 | QFC1=XLV*CHKLOWQ(I,J)*AKHS_DENS |
---|
571 | ! |
---|
572 | IF(SNOW(I,J)>0..OR.SICE(I,J)>0.5)THEN |
---|
573 | QFC1=QFC1*RLIVWV |
---|
574 | ENDIF |
---|
575 | ! |
---|
576 | IF(QFC1>0.)THEN |
---|
577 | QLOW=QK(KTE+1-LLOW) |
---|
578 | QSFC(I,J)=QLOW+ELFLX(I,J)/QFC1 |
---|
579 | ENDIF |
---|
580 | ! |
---|
581 | ELSE |
---|
582 | PSFC=PINT(I,LOWLYR(I,J),J) |
---|
583 | EXNSFC=(1.E5/PSFC)**CAPA |
---|
584 | QSFC(I,J)=PQ0SEA/PSFC & |
---|
585 | & *EXP(A2*(THSK(I,J)-A3*EXNSFC)/(THSK(I,J)-A4*EXNSFC)) |
---|
586 | ENDIF |
---|
587 | ! |
---|
588 | QZ0 (I,J)=(1.-SEAMASK)*QSFC(I,J)+SEAMASK*QZ0 (I,J) |
---|
589 | ! |
---|
590 | !*** LOWEST LAYER ABOVE GROUND MUST BE FLIPPED |
---|
591 | ! |
---|
592 | LMH=KTE-LOWLYR(I,J)+1 |
---|
593 | ! |
---|
594 | !---------------------------------------------------------------------- |
---|
595 | !*** CARRY OUT THE VERTICAL DIFFUSION OF |
---|
596 | !*** TEMPERATURE AND WATER VAPOR |
---|
597 | !---------------------------------------------------------------------- |
---|
598 | ! |
---|
599 | |
---|
600 | CALL VDIFH(DTDIF,LMH,THZ0(I,J),QZ0(I,J) & |
---|
601 | & ,AKHS_DENS,CHKLOWQ(I,J),CT(I,J) & |
---|
602 | & ,THEK,QK,CWMK,AKHK,ZHK,RHOK(KTS,I) & |
---|
603 | ! Multi-layer UCM |
---|
604 | & ,A_T,B_T,A_Q,B_Q,SFK,VLK & |
---|
605 | ! UCM end |
---|
606 | & ,IDS,IDE,JDS,JDE,KDS,KDE & |
---|
607 | & ,IMS,IME,JMS,JME,KMS,KME & |
---|
608 | & ,ITS,ITE,JTS,JTE,KTS,KTE,I,J) |
---|
609 | !---------------------------------------------------------------------- |
---|
610 | !*** |
---|
611 | !*** COMPUTE PRIMARY VARIABLE TENDENCIES |
---|
612 | !*** |
---|
613 | DO K=KTS,KTE |
---|
614 | KFLIP=KTE+1-K |
---|
615 | THOLD=TH(I,K,J) |
---|
616 | THNEW=THEK(KFLIP)+CWMK(KFLIP)*ELOCP*APE(I,K,J) |
---|
617 | DTDT=(THNEW-THOLD)*RDTTURBL |
---|
618 | QOLD=QV(I,K,J)/(1.+QV(I,K,J)) |
---|
619 | DQDT=(QK(KFLIP)-QOLD)*RDTTURBL |
---|
620 | DCDT=(CWMK(KFLIP)-CWM(I,K,J))*RDTTURBL |
---|
621 | ! |
---|
622 | RTHBLTEN(I,K,J)=DTDT |
---|
623 | RQVBLTEN(I,K,J)=DQDT/(1.-QK(KFLIP))**2 |
---|
624 | RQCBLTEN(I,K,J)=DCDT |
---|
625 | ENDDO |
---|
626 | ! |
---|
627 | |
---|
628 | !*** Begin debugging |
---|
629 | ! IF(I==IMD.AND.J==JMD)THEN |
---|
630 | ! PRINT_DIAG=0 |
---|
631 | ! ELSE |
---|
632 | ! PRINT_DIAG=0 |
---|
633 | ! ENDIF |
---|
634 | ! IF(I==227.AND.J==363)PRINT_DIAG=0 |
---|
635 | !*** End debugging |
---|
636 | ! |
---|
637 | PSFC=.01*PINT(I,LOWLYR(I,J),J) |
---|
638 | ZSL_DIAG=0.5*DZ(I,1,J) |
---|
639 | ! |
---|
640 | !*** Begin debugging |
---|
641 | ! IF(PRINT_DIAG==1)THEN |
---|
642 | ! |
---|
643 | ! write(6,"(a, 2i5, 2i3, 2f8.2, f6.2, 2f8.2)") & |
---|
644 | ! '{turb4 i,j, Kpbl, Kmxl, Psfc, Zsfc, Zsl, Zpbl, Zmxl = ' & |
---|
645 | ! , i, j, KPBL(i,j), KTE-LMXL+1, PSFC, ZHK(LMH+1), ZSL_diag & |
---|
646 | ! , PBLH(i,j), ZHK(LMXL)-ZHK(LMH+1) |
---|
647 | ! write(6,"(a, 2f7.2, f7.3, 3e11.4)") & |
---|
648 | ! '{turb4 tsk, thsk, qz0, q**2_0, akhs, exch_0 = ' & |
---|
649 | ! , tsk(i,j)-273.15, thsk(i,j), 1000.*qz0(i,j) & |
---|
650 | ! , 2.*tke_myj(i,1,j), akhs(i,j), akhs(i,j)*ZSL_diag |
---|
651 | ! write(6,"(a)") & |
---|
652 | ! '{turb5 k, Pmid, Pint_1, Tc, TH, DTH, GH, GM, EL, Q**2, Akh, EXCH_h, Dz, Dp' |
---|
653 | ! do k=kts,kte/2 |
---|
654 | ! KFLIP=KTE-K !-- Includes the KFLIP-1 in earlier versions |
---|
655 | ! write(6,"(a,i3, 2f8.2, 2f8.3, 3e12.4, 4e11.4, f7.2, f6.2)") & |
---|
656 | ! '{turb5 ', k, .01*pmid(i,k,j),.01*pint(i,k,j), T(i,k,j)-273.15 & |
---|
657 | ! , th(i,k,j), DTTURBL*rthblten(i,k,j), GH(KFLIP), GM(KFLIP) & |
---|
658 | ! , el_myj(i,KFLIP,j), 2.*tke_myj(i,k+1,j), akh(i,KFLIP,j) & |
---|
659 | ! , exch_h(i,k,j), dz(i,k,j), .01*(pint(i,k,j)-pint(i,k+1,j)) |
---|
660 | ! enddo |
---|
661 | ! |
---|
662 | ! ELSEIF(PRINT_DIAG==2)THEN |
---|
663 | ! |
---|
664 | ! write(6,"(a, 2i5, 2i3, 2f8.2, f6.2, 2f8.2)") & |
---|
665 | ! '}turb4 i,j, Kpbl, Kmxl, Psfc, Zsfc, Zsl, Zpbl, Zmxl = ' & |
---|
666 | ! , i, j, KPBL(i,j), KTE-LMXL+1, PSFC, ZHK(LMH+1), ZSL_diag & |
---|
667 | ! , PBLH(i,j), ZHK(LMXL)-ZHK(LMH+1) |
---|
668 | ! write(6,"(a, 2f7.2, f7.3, 3e11.4)") & |
---|
669 | ! '}turb4 tsk, thsk, qz0, q**2_0, akhs, exch_0 = ' & |
---|
670 | ! , tsk(i,j)-273.15, thsk(i,j), 1000.*qz0(i,j) & |
---|
671 | ! , 2.*tke_myj(i,1,j), akhs(i,j), akhs(i,j)*ZSL_diag |
---|
672 | ! write(6,"(a)") & |
---|
673 | ! '}turb5 k, Pmid, Pint_1, Tc, TH, DTH, GH, GM, EL, Q**2, Akh, EXCH_h, Dz, Dp' |
---|
674 | ! do k=kts,kte/2 |
---|
675 | ! KFLIP=KTE-K !-- Includes the KFLIP-1 in earlier versions |
---|
676 | ! write(6,"(a,i3, 2f8.2, 2f8.3, 3e12.4, 4e11.4, f7.2, f6.2)") & |
---|
677 | ! '}turb5 ', k, .01*pmid(i,k,j),.01*pint(i,k,j), T(i,k,j)-273.15 & |
---|
678 | ! , th(i,k,j), DTTURBL*rthblten(i,k,j), GH(KFLIP), GM(KFLIP) & |
---|
679 | ! , el_myj(i,KFLIP,j), 2.*tke_myj(i,k+1,j), akh(i,KFLIP,j) & |
---|
680 | ! , exch_h(i,k,j), dz(i,k,j), .01*(pint(i,k,j)-pint(i,k+1,j)) |
---|
681 | ! enddo |
---|
682 | ! ENDIF |
---|
683 | !*** End debugging |
---|
684 | ! |
---|
685 | !---------------------------------------------------------------------- |
---|
686 | ENDDO |
---|
687 | !---------------------------------------------------------------------- |
---|
688 | DO I=ITS,ITE |
---|
689 | ! |
---|
690 | !*** FILL 1-D VERTICAL ARRAYS |
---|
691 | !*** AND FLIP DIRECTION SINCE MYJ SCHEME |
---|
692 | !*** COUNTS DOWNWARD FROM THE DOMAIN'S TOP |
---|
693 | ! |
---|
694 | |
---|
695 | DO K=KTS,KTE-1 |
---|
696 | AKMK(K)=AKM(I,K,J) |
---|
697 | AKMK(K)=AKMK(K)*(RHOK(K,I)+RHOK(K+1,I))*0.5 |
---|
698 | ENDDO |
---|
699 | ! |
---|
700 | LLOW=LOWLYR(I,J) |
---|
701 | AKMS_DENS=AKMS(I,J)*RHOK(KTE+1-LLOW,I) |
---|
702 | ! |
---|
703 | |
---|
704 | DO K=KTE,KTS,-1 |
---|
705 | KFLIP=KTE+1-K |
---|
706 | UK(K)=U(I,KFLIP,J) |
---|
707 | VK(K)=V(I,KFLIP,J) |
---|
708 | ZHK(K)=ZINT(I,K,J) |
---|
709 | |
---|
710 | ENDDO |
---|
711 | ! SFK(1)=1. |
---|
712 | if(flag_bep)then |
---|
713 | ! Multi-layer UCM |
---|
714 | DO K=KTE,KTS,-1 |
---|
715 | KFLIP=KTE+1-K |
---|
716 | A_U(K)=A_U_BEP(I,KFLIP,J) |
---|
717 | A_V(K)=A_V_BEP(I,KFLIP,J) |
---|
718 | B_U(K)=B_U_BEP(I,KFLIP,J) |
---|
719 | B_V(K)=B_V_BEP(I,KFLIP,J) |
---|
720 | SFK(K)=SF_BEP(I,KFLIP,J) |
---|
721 | VLK(K)=VL_BEP(I,KFLIP,J) |
---|
722 | ENDDO |
---|
723 | ! UCM end |
---|
724 | else |
---|
725 | DO K=KTE,KTS,-1 |
---|
726 | A_U(K)=0. |
---|
727 | A_V(K)=0. |
---|
728 | B_U(K)=0. |
---|
729 | B_V(K)=0. |
---|
730 | SFK(K)=1. |
---|
731 | VLK(K)=1. |
---|
732 | ENDDO |
---|
733 | endif |
---|
734 | |
---|
735 | ZHK(KTE+1)=ZINT(I,KTE+1,J) |
---|
736 | ! |
---|
737 | !---------------------------------------------------------------------- |
---|
738 | !*** CARRY OUT THE VERTICAL DIFFUSION OF |
---|
739 | !*** VELOCITY COMPONENTS |
---|
740 | !---------------------------------------------------------------------- |
---|
741 | ! |
---|
742 | CALL VDIFV(LMH,DTDIF,UZ0(I,J),VZ0(I,J) & |
---|
743 | & ,AKMS_DENS,UK,VK,AKMK,ZHK,RHOK(KTS,I) & |
---|
744 | ! Multi-layer UCM |
---|
745 | & ,A_U,A_V,B_U,B_V,SFK,VLK & |
---|
746 | ! UCM end |
---|
747 | & ,IDS,IDE,JDS,JDE,KDS,KDE & |
---|
748 | & ,IMS,IME,JMS,JME,KMS,KME & |
---|
749 | & ,ITS,ITE,JTS,JTE,KTS,KTE,I,J) |
---|
750 | ! |
---|
751 | !---------------------------------------------------------------------- |
---|
752 | !*** |
---|
753 | !*** COMPUTE PRIMARY VARIABLE TENDENCIES |
---|
754 | !*** |
---|
755 | DO K=KTS,KTE |
---|
756 | KFLIP=KTE+1-K |
---|
757 | DUDT=(UK(KFLIP)-U(I,K,J))*RDTTURBL |
---|
758 | DVDT=(VK(KFLIP)-V(I,K,J))*RDTTURBL |
---|
759 | RUBLTEN(I,K,J)=DUDT |
---|
760 | RVBLTEN(I,K,J)=DVDT |
---|
761 | ENDDO |
---|
762 | ! |
---|
763 | ENDDO |
---|
764 | !---------------------------------------------------------------------- |
---|
765 | ! |
---|
766 | ENDDO main_integration |
---|
767 | ! |
---|
768 | ENDIF |
---|
769 | !---------------------------------------------------------------------- |
---|
770 | ! |
---|
771 | END SUBROUTINE MYJURB |
---|
772 | ! |
---|
773 | !---------------------------------------------------------------------- |
---|
774 | !XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX |
---|
775 | !---------------------------------------------------------------------- |
---|
776 | SUBROUTINE MIXLEN & |
---|
777 | !---------------------------------------------------------------------- |
---|
778 | ! ****************************************************************** |
---|
779 | ! * * |
---|
780 | ! * LEVEL 2.5 MIXING LENGTH * |
---|
781 | ! * * |
---|
782 | ! ****************************************************************** |
---|
783 | ! |
---|
784 | &(LMH,U,V,T,THE,Q,CWM,Q2,Z,GM,GH,EL,PBLH,LPBL,LMXL,CT & |
---|
785 | ! Multi-layer UCM |
---|
786 | &,DL_U & |
---|
787 | ! UCM end |
---|
788 | &,IDS,IDE,JDS,JDE,KDS,KDE & |
---|
789 | &,IMS,IME,JMS,JME,KMS,KME & |
---|
790 | &,ITS,ITE,JTS,JTE,KTS,KTE) |
---|
791 | !---------------------------------------------------------------------- |
---|
792 | ! |
---|
793 | IMPLICIT NONE |
---|
794 | ! |
---|
795 | !---------------------------------------------------------------------- |
---|
796 | INTEGER,INTENT(IN) :: IDS,IDE,JDS,JDE,KDS,KDE & |
---|
797 | & ,IMS,IME,JMS,JME,KMS,KME & |
---|
798 | & ,ITS,ITE,JTS,JTE,KTS,KTE |
---|
799 | ! |
---|
800 | INTEGER,INTENT(IN) :: LMH |
---|
801 | ! |
---|
802 | INTEGER,INTENT(OUT) :: LMXL,LPBL |
---|
803 | ! |
---|
804 | REAL,DIMENSION(KTS:KTE),INTENT(IN) :: CWM,Q,Q2,T,THE,U,V |
---|
805 | ! Multi-layer UCM |
---|
806 | REAL,DIMENSION(KTS:KTE),INTENT(IN) :: DL_U |
---|
807 | ! UCM end |
---|
808 | ! |
---|
809 | REAL,DIMENSION(KTS:KTE+1),INTENT(IN) :: Z |
---|
810 | ! |
---|
811 | REAL,INTENT(OUT) :: PBLH |
---|
812 | ! |
---|
813 | REAL,DIMENSION(KTS:KTE-1),INTENT(OUT) :: EL,GH,GM |
---|
814 | ! |
---|
815 | REAL,INTENT(INOUT) :: CT |
---|
816 | !---------------------------------------------------------------------- |
---|
817 | !*** |
---|
818 | !*** LOCAL VARIABLES |
---|
819 | !*** |
---|
820 | INTEGER :: K,LPBLM |
---|
821 | ! |
---|
822 | REAL :: A,ADEN,B,BDEN,AUBR,BUBR,BLMX,EL0,ELOQ2X,GHL,GML & |
---|
823 | & ,QOL2ST,QOL2UN,QDZL,RDZ,SQ,SREL,SZQ,TEM,THM,VKRMZ |
---|
824 | ! |
---|
825 | REAL,DIMENSION(KTS:KTE) :: Q1 |
---|
826 | ! |
---|
827 | REAL,DIMENSION(KTS:KTE-1) :: DTH,ELM,REL |
---|
828 | ! |
---|
829 | !---------------------------------------------------------------------- |
---|
830 | !********************************************************************** |
---|
831 | !--------------FIND THE HEIGHT OF THE PBL------------------------------- |
---|
832 | LPBL=LMH |
---|
833 | ! |
---|
834 | DO K=LMH-1,1,-1 |
---|
835 | IF(Q2(K)<=EPSQ2*FH)THEN |
---|
836 | LPBL=K |
---|
837 | GO TO 110 |
---|
838 | ENDIF |
---|
839 | ENDDO |
---|
840 | ! |
---|
841 | LPBL=1 |
---|
842 | ! |
---|
843 | !--------------THE HEIGHT OF THE PBL------------------------------------ |
---|
844 | ! |
---|
845 | 110 PBLH=Z(LPBL+1)-Z(LMH+1) |
---|
846 | ! |
---|
847 | !----------------------------------------------------------------------- |
---|
848 | DO K=KTS,LMH |
---|
849 | Q1(K)=0. |
---|
850 | ENDDO |
---|
851 | ! |
---|
852 | DO K=1,LMH-1 |
---|
853 | DTH(K)=THE(K)-THE(K+1) |
---|
854 | ENDDO |
---|
855 | ! |
---|
856 | DO K=LMH-2,1,-1 |
---|
857 | IF(DTH(K)>0..AND.DTH(K+1)<=0.)THEN |
---|
858 | DTH(K)=DTH(K)+CT |
---|
859 | EXIT |
---|
860 | ENDIF |
---|
861 | ENDDO |
---|
862 | ! |
---|
863 | CT=0. |
---|
864 | !---------------------------------------------------------------------- |
---|
865 | DO K=KTS,LMH-1 |
---|
866 | RDZ=2./(Z(K)-Z(K+2)) |
---|
867 | GML=((U(K)-U(K+1))**2+(V(K)-V(K+1))**2)*RDZ*RDZ |
---|
868 | GM(K)=MAX(GML,EPSGM) |
---|
869 | ! |
---|
870 | TEM=(T(K)+T(K+1))*0.5 |
---|
871 | THM=(THE(K)+THE(K+1))*0.5 |
---|
872 | ! |
---|
873 | A=THM*P608 |
---|
874 | B=(ELOCP/TEM-1.-P608)*THM |
---|
875 | ! |
---|
876 | GHL=(DTH(K)*((Q(K)+Q(K+1)+CWM(K)+CWM(K+1))*(0.5*P608)+1.) & |
---|
877 | & +(Q(K)-Q(K+1)+CWM(K)-CWM(K+1))*A & |
---|
878 | & +(CWM(K)-CWM(K+1))*B)*RDZ |
---|
879 | ! |
---|
880 | IF(ABS(GHL)<=EPSGH)GHL=EPSGH |
---|
881 | GH(K)=GHL |
---|
882 | ENDDO |
---|
883 | ! |
---|
884 | !---------------------------------------------------------------------- |
---|
885 | !*** FIND MAXIMUM MIXING LENGTHS AND THE LEVEL OF THE PBL TOP |
---|
886 | !---------------------------------------------------------------------- |
---|
887 | ! |
---|
888 | LMXL=LMH |
---|
889 | ! |
---|
890 | DO K=KTS,LMH-1 |
---|
891 | GML=GM(K) |
---|
892 | GHL=GH(K) |
---|
893 | ! |
---|
894 | IF(GHL>=EPSGH)THEN |
---|
895 | IF(GML/GHL<=REQU)THEN |
---|
896 | ELM(K)=EPSL |
---|
897 | LMXL=K |
---|
898 | ELSE |
---|
899 | AUBR=(AUBM*GML+AUBH*GHL)*GHL |
---|
900 | BUBR= BUBM*GML+BUBH*GHL |
---|
901 | QOL2ST=(-0.5*BUBR+SQRT(BUBR*BUBR*0.25-AUBR*CUBR))*RCUBR |
---|
902 | ELOQ2X=1./QOL2ST |
---|
903 | ELM(K)=MAX(SQRT(ELOQ2X*Q2(K)),EPSL) |
---|
904 | ENDIF |
---|
905 | ELSE |
---|
906 | ADEN=(ADNM*GML+ADNH*GHL)*GHL |
---|
907 | BDEN= BDNM*GML+BDNH*GHL |
---|
908 | QOL2UN=-0.5*BDEN+SQRT(BDEN*BDEN*0.25-ADEN) |
---|
909 | ELOQ2X=1./(QOL2UN+EPSRU) ! repsr1/qol2un |
---|
910 | ELM(K)=MAX(SQRT(ELOQ2X*Q2(K)),EPSL) |
---|
911 | ENDIF |
---|
912 | ENDDO |
---|
913 | ! |
---|
914 | IF(ELM(LMH-1)==EPSL)LMXL=LMH |
---|
915 | ! |
---|
916 | !---------------------------------------------------------------------- |
---|
917 | !*** THE HEIGHT OF THE MIXED LAYER |
---|
918 | !---------------------------------------------------------------------- |
---|
919 | ! |
---|
920 | BLMX=Z(LMXL)-Z(LMH+1) |
---|
921 | ! |
---|
922 | !---------------------------------------------------------------------- |
---|
923 | DO K=LPBL,LMH |
---|
924 | Q1(K)=SQRT(Q2(K)) |
---|
925 | ENDDO |
---|
926 | !---------------------------------------------------------------------- |
---|
927 | SZQ=0. |
---|
928 | SQ =0. |
---|
929 | ! |
---|
930 | DO K=KTS,LMH-1 |
---|
931 | QDZL=(Q1(K)+Q1(K+1))*(Z(K+1)-Z(K+2)) |
---|
932 | SZQ=(Z(K+1)+Z(K+2)-Z(LMH+1)-Z(LMH+1))*QDZL+SZQ |
---|
933 | SQ=QDZL+SQ |
---|
934 | ENDDO |
---|
935 | ! |
---|
936 | !---------------------------------------------------------------------- |
---|
937 | !*** COMPUTATION OF ASYMPTOTIC L IN BLACKADAR FORMULA |
---|
938 | !---------------------------------------------------------------------- |
---|
939 | ! |
---|
940 | EL0=MIN(ALPH*SZQ*0.5/SQ,EL0MAX) |
---|
941 | EL0=MAX(EL0 ,EL0MIN) |
---|
942 | ! |
---|
943 | !---------------------------------------------------------------------- |
---|
944 | !*** ABOVE THE PBL TOP |
---|
945 | !---------------------------------------------------------------------- |
---|
946 | ! |
---|
947 | LPBLM=MAX(LPBL-1,1) |
---|
948 | ! |
---|
949 | DO K=KTS,LPBLM |
---|
950 | EL(K)=MIN((Z(K)-Z(K+2))*ELFC,ELM(K)) |
---|
951 | REL(K)=EL(K)/ELM(K) |
---|
952 | ENDDO |
---|
953 | ! |
---|
954 | !---------------------------------------------------------------------- |
---|
955 | !*** INSIDE THE PBL |
---|
956 | !---------------------------------------------------------------------- |
---|
957 | ! |
---|
958 | IF(LPBL<LMH)THEN |
---|
959 | DO K=LPBL,LMH-1 |
---|
960 | VKRMZ=(Z(K+1)-Z(LMH+1))*VKARMAN |
---|
961 | ! Multi-layer UCM |
---|
962 | IF (DL_U(K).GT.0.0) THEN |
---|
963 | EL(K)=MIN(VKRMZ/(VKRMZ/EL0+1.),ELM(K)) |
---|
964 | EL(K)=1./(1./EL(K)+1./DL_U(K)) |
---|
965 | ELSE |
---|
966 | EL(K)=MIN(VKRMZ/(VKRMZ/EL0+1.),ELM(K)) |
---|
967 | ENDIF |
---|
968 | ! UCM end |
---|
969 | REL(K)=EL(K)/ELM(K) |
---|
970 | ENDDO |
---|
971 | ENDIF |
---|
972 | ! |
---|
973 | DO K=LPBL+1,LMH-2 |
---|
974 | SREL=MIN(((REL(K-1)+REL(K+1))*0.5+REL(K))*0.5,REL(K)) |
---|
975 | EL(K)=MAX(SREL*ELM(K),EPSL) |
---|
976 | ENDDO |
---|
977 | ! |
---|
978 | !---------------------------------------------------------------------- |
---|
979 | END SUBROUTINE MIXLEN |
---|
980 | !---------------------------------------------------------------------- |
---|
981 | !XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX |
---|
982 | !---------------------------------------------------------------------- |
---|
983 | SUBROUTINE PRODQ2 & |
---|
984 | !---------------------------------------------------------------------- |
---|
985 | ! ****************************************************************** |
---|
986 | ! * * |
---|
987 | ! * LEVEL 2.5 Q2 PRODUCTION/DISSIPATION * |
---|
988 | ! * * |
---|
989 | ! ****************************************************************** |
---|
990 | ! |
---|
991 | &(LMH,DTTURBL,USTAR,GM,GH,EL,Q2 & |
---|
992 | &,IDS,IDE,JDS,JDE,KDS,KDE & |
---|
993 | &,IMS,IME,JMS,JME,KMS,KME & |
---|
994 | &,ITS,ITE,JTS,JTE,KTS,KTE) |
---|
995 | !---------------------------------------------------------------------- |
---|
996 | ! |
---|
997 | IMPLICIT NONE |
---|
998 | ! |
---|
999 | !---------------------------------------------------------------------- |
---|
1000 | INTEGER,INTENT(IN) :: IDS,IDE,JDS,JDE,KDS,KDE & |
---|
1001 | & ,IMS,IME,JMS,JME,KMS,KME & |
---|
1002 | & ,ITS,ITE,JTS,JTE,KTS,KTE |
---|
1003 | ! |
---|
1004 | INTEGER,INTENT(IN) :: LMH |
---|
1005 | ! |
---|
1006 | REAL,INTENT(IN) :: DTTURBL,USTAR |
---|
1007 | ! |
---|
1008 | REAL,DIMENSION(KTS:KTE-1),INTENT(IN) :: GH,GM |
---|
1009 | REAL,DIMENSION(KTS:KTE-1),INTENT(INOUT) :: EL |
---|
1010 | ! |
---|
1011 | REAL,DIMENSION(KTS:KTE),INTENT(INOUT) :: Q2 |
---|
1012 | !---------------------------------------------------------------------- |
---|
1013 | !*** |
---|
1014 | !*** LOCAL VARIABLES |
---|
1015 | !*** |
---|
1016 | INTEGER :: K |
---|
1017 | ! |
---|
1018 | REAL :: ADEN,AEQU,ANUM,ARHS,BDEN,BEQU,BNUM,BRHS,CDEN,CRHS & |
---|
1019 | & ,DLOQ1,ELOQ11,ELOQ12,ELOQ13,ELOQ21,ELOQ22,ELOQ31,ELOQ32 & |
---|
1020 | & ,ELOQ41,ELOQ42,ELOQ51,ELOQ52,ELOQN,EQOL2,GHL,GML & |
---|
1021 | & ,RDEN1,RDEN2,RHS2,RHSP1,RHSP2,RHST2 |
---|
1022 | ! |
---|
1023 | !---------------------------------------------------------------------- |
---|
1024 | !********************************************************************** |
---|
1025 | !---------------------------------------------------------------------- |
---|
1026 | ! |
---|
1027 | main_integration: DO K=1,LMH-1 |
---|
1028 | GML=GM(K) |
---|
1029 | GHL=GH(K) |
---|
1030 | ! |
---|
1031 | !---------------------------------------------------------------------- |
---|
1032 | !*** COEFFICIENTS OF THE EQUILIBRIUM EQUATION |
---|
1033 | !---------------------------------------------------------------------- |
---|
1034 | ! |
---|
1035 | AEQU=(AEQM*GML+AEQH*GHL)*GHL |
---|
1036 | BEQU= BEQM*GML+BEQH*GHL |
---|
1037 | ! |
---|
1038 | !---------------------------------------------------------------------- |
---|
1039 | !*** EQUILIBRIUM SOLUTION FOR L/Q |
---|
1040 | !---------------------------------------------------------------------- |
---|
1041 | ! |
---|
1042 | EQOL2=-0.5*BEQU+SQRT(BEQU*BEQU*0.25-AEQU) |
---|
1043 | ! |
---|
1044 | !---------------------------------------------------------------------- |
---|
1045 | !*** IS THERE PRODUCTION/DISSIPATION ? |
---|
1046 | !---------------------------------------------------------------------- |
---|
1047 | ! |
---|
1048 | IF((GML+GHL*GHL<=EPSTRB) & |
---|
1049 | & .OR.(GHL>=EPSGH.AND.GML/GHL<=REQU) & |
---|
1050 | & .OR.(EQOL2<=EPS2))THEN |
---|
1051 | ! |
---|
1052 | !---------------------------------------------------------------------- |
---|
1053 | !*** NO TURBULENCE |
---|
1054 | !---------------------------------------------------------------------- |
---|
1055 | ! |
---|
1056 | Q2(K)=EPSQ2 |
---|
1057 | EL(K)=EPSL |
---|
1058 | !---------------------------------------------------------------------- |
---|
1059 | ! |
---|
1060 | ELSE |
---|
1061 | ! |
---|
1062 | !---------------------------------------------------------------------- |
---|
1063 | !*** TURBULENCE |
---|
1064 | !---------------------------------------------------------------------- |
---|
1065 | !---------------------------------------------------------------------- |
---|
1066 | !*** COEFFICIENTS OF THE TERMS IN THE NUMERATOR |
---|
1067 | !---------------------------------------------------------------------- |
---|
1068 | ! |
---|
1069 | ANUM=(ANMM*GML+ANMH*GHL)*GHL |
---|
1070 | BNUM= BNMM*GML+BNMH*GHL |
---|
1071 | ! |
---|
1072 | !---------------------------------------------------------------------- |
---|
1073 | !*** COEFFICIENTS OF THE TERMS IN THE DENOMINATOR |
---|
1074 | !---------------------------------------------------------------------- |
---|
1075 | ! |
---|
1076 | ADEN=(ADNM*GML+ADNH*GHL)*GHL |
---|
1077 | BDEN= BDNM*GML+BDNH*GHL |
---|
1078 | CDEN= 1. |
---|
1079 | ! |
---|
1080 | !---------------------------------------------------------------------- |
---|
1081 | !*** COEFFICIENTS OF THE NUMERATOR OF THE LINEARIZED EQ. |
---|
1082 | !---------------------------------------------------------------------- |
---|
1083 | ! |
---|
1084 | ARHS=-(ANUM*BDEN-BNUM*ADEN)*2. |
---|
1085 | BRHS=- ANUM*4. |
---|
1086 | CRHS=- BNUM*2. |
---|
1087 | ! |
---|
1088 | !---------------------------------------------------------------------- |
---|
1089 | !*** INITIAL VALUE OF L/Q |
---|
1090 | !---------------------------------------------------------------------- |
---|
1091 | ! |
---|
1092 | DLOQ1=EL(K)/SQRT(Q2(K)) |
---|
1093 | ! |
---|
1094 | !---------------------------------------------------------------------- |
---|
1095 | !*** FIRST ITERATION FOR L/Q, RHS=0 |
---|
1096 | !---------------------------------------------------------------------- |
---|
1097 | ! |
---|
1098 | ELOQ21=1./EQOL2 |
---|
1099 | ELOQ11=SQRT(ELOQ21) |
---|
1100 | ELOQ31=ELOQ21*ELOQ11 |
---|
1101 | ELOQ41=ELOQ21*ELOQ21 |
---|
1102 | ELOQ51=ELOQ21*ELOQ31 |
---|
1103 | ! |
---|
1104 | !---------------------------------------------------------------------- |
---|
1105 | !*** 1./DENOMINATOR |
---|
1106 | !---------------------------------------------------------------------- |
---|
1107 | ! |
---|
1108 | RDEN1=1./(ADEN*ELOQ41+BDEN*ELOQ21+CDEN) |
---|
1109 | ! |
---|
1110 | !---------------------------------------------------------------------- |
---|
1111 | !*** D(RHS)/D(L/Q) |
---|
1112 | !---------------------------------------------------------------------- |
---|
1113 | ! |
---|
1114 | RHSP1=(ARHS*ELOQ51+BRHS*ELOQ31+CRHS*ELOQ11)*RDEN1*RDEN1 |
---|
1115 | ! |
---|
1116 | !---------------------------------------------------------------------- |
---|
1117 | !*** FIRST-GUESS SOLUTION |
---|
1118 | !---------------------------------------------------------------------- |
---|
1119 | ! |
---|
1120 | ELOQ12=ELOQ11+(DLOQ1-ELOQ11)*EXP(RHSP1*DTTURBL) |
---|
1121 | ELOQ12=MAX(ELOQ12,EPS1) |
---|
1122 | ! |
---|
1123 | !---------------------------------------------------------------------- |
---|
1124 | !*** SECOND ITERATION FOR L/Q |
---|
1125 | !---------------------------------------------------------------------- |
---|
1126 | ! |
---|
1127 | ELOQ22=ELOQ12*ELOQ12 |
---|
1128 | ELOQ32=ELOQ22*ELOQ12 |
---|
1129 | ELOQ42=ELOQ22*ELOQ22 |
---|
1130 | ELOQ52=ELOQ22*ELOQ32 |
---|
1131 | ! |
---|
1132 | !---------------------------------------------------------------------- |
---|
1133 | !*** 1./DENOMINATOR |
---|
1134 | !---------------------------------------------------------------------- |
---|
1135 | ! |
---|
1136 | RDEN2=1./(ADEN*ELOQ42+BDEN*ELOQ22+CDEN) |
---|
1137 | RHS2 =-(ANUM*ELOQ42+BNUM*ELOQ22)*RDEN2+RB1 |
---|
1138 | RHSP2= (ARHS*ELOQ52+BRHS*ELOQ32+CRHS*ELOQ12)*RDEN2*RDEN2 |
---|
1139 | RHST2=RHS2/RHSP2 |
---|
1140 | ! |
---|
1141 | !---------------------------------------------------------------------- |
---|
1142 | !*** CORRECTED SOLUTION |
---|
1143 | !---------------------------------------------------------------------- |
---|
1144 | ! |
---|
1145 | ELOQ13=ELOQ12-RHST2+(RHST2+DLOQ1-ELOQ12)*EXP(RHSP2*DTTURBL) |
---|
1146 | ELOQ13=AMAX1(ELOQ13,EPS1) |
---|
1147 | ! |
---|
1148 | !---------------------------------------------------------------------- |
---|
1149 | !*** TWO ITERATIONS IS ENOUGH IN MOST CASES ... |
---|
1150 | !---------------------------------------------------------------------- |
---|
1151 | ! |
---|
1152 | ELOQN=ELOQ13 |
---|
1153 | ! |
---|
1154 | IF(ELOQN>EPS1)THEN |
---|
1155 | Q2(K)=EL(K)*EL(K)/(ELOQN*ELOQN) |
---|
1156 | Q2(K)=AMAX1(Q2(K),EPSQ2) |
---|
1157 | IF(Q2(K)==EPSQ2)THEN |
---|
1158 | EL(K)=EPSL |
---|
1159 | ENDIF |
---|
1160 | ELSE |
---|
1161 | Q2(K)=EPSQ2 |
---|
1162 | EL(K)=EPSL |
---|
1163 | ENDIF |
---|
1164 | ! |
---|
1165 | !---------------------------------------------------------------------- |
---|
1166 | !*** END OF TURBULENT BRANCH |
---|
1167 | !---------------------------------------------------------------------- |
---|
1168 | ! |
---|
1169 | ENDIF |
---|
1170 | !---------------------------------------------------------------------- |
---|
1171 | !*** END OF PRODUCTION/DISSIPATION LOOP |
---|
1172 | !---------------------------------------------------------------------- |
---|
1173 | ! |
---|
1174 | ENDDO main_integration |
---|
1175 | ! |
---|
1176 | !---------------------------------------------------------------------- |
---|
1177 | !*** LOWER BOUNDARY CONDITION FOR Q2 |
---|
1178 | !---------------------------------------------------------------------- |
---|
1179 | ! |
---|
1180 | Q2(LMH)=AMAX1(B1**(2./3.)*USTAR*USTAR,EPSQ2) |
---|
1181 | !---------------------------------------------------------------------- |
---|
1182 | ! |
---|
1183 | END SUBROUTINE PRODQ2 |
---|
1184 | ! |
---|
1185 | !---------------------------------------------------------------------- |
---|
1186 | !XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX |
---|
1187 | !---------------------------------------------------------------------- |
---|
1188 | SUBROUTINE DIFCOF & |
---|
1189 | ! ****************************************************************** |
---|
1190 | ! * * |
---|
1191 | ! * LEVEL 2.5 DIFFUSION COEFFICIENTS * |
---|
1192 | ! * * |
---|
1193 | ! ****************************************************************** |
---|
1194 | &(LMH,LMXL,GM,GH,EL,T,Q2,Z,AKM,AKH & |
---|
1195 | &,IDS,IDE,JDS,JDE,KDS,KDE & |
---|
1196 | &,IMS,IME,JMS,JME,KMS,KME & |
---|
1197 | &,ITS,ITE,JTS,JTE,KTS,KTE,PRINT_DIAG) ! debug |
---|
1198 | !---------------------------------------------------------------------- |
---|
1199 | ! |
---|
1200 | IMPLICIT NONE |
---|
1201 | ! |
---|
1202 | !---------------------------------------------------------------------- |
---|
1203 | INTEGER,INTENT(IN) :: IDS,IDE,JDS,JDE,KDS,KDE & |
---|
1204 | & ,IMS,IME,JMS,JME,KMS,KME & |
---|
1205 | & ,ITS,ITE,JTS,JTE,KTS,KTE |
---|
1206 | ! |
---|
1207 | INTEGER,INTENT(IN) :: LMH,LMXL |
---|
1208 | ! |
---|
1209 | REAL,DIMENSION(KTS:KTE),INTENT(IN) :: Q2,T |
---|
1210 | REAL,DIMENSION(KTS:KTE-1),INTENT(IN) :: EL,GH,GM |
---|
1211 | REAL,DIMENSION(KTS:KTE+1),INTENT(IN) :: Z |
---|
1212 | ! |
---|
1213 | REAL,DIMENSION(KTS:KTE-1),INTENT(OUT) :: AKH,AKM |
---|
1214 | !---------------------------------------------------------------------- |
---|
1215 | !*** |
---|
1216 | !*** LOCAL VARIABLES |
---|
1217 | !*** |
---|
1218 | INTEGER :: K,KINV |
---|
1219 | ! |
---|
1220 | REAL :: ADEN,AKMIN,BDEN,BESH,BESM,CDEN,D2T,ELL,ELOQ2,ELOQ4,ELQDZ & |
---|
1221 | & ,ESH,ESM,GHL,GML,Q1L,RDEN,RDZ |
---|
1222 | ! |
---|
1223 | !*** Begin debugging |
---|
1224 | INTEGER,INTENT(IN) :: PRINT_DIAG |
---|
1225 | ! REAL :: D2Tmin |
---|
1226 | !*** End debugging |
---|
1227 | ! |
---|
1228 | !---------------------------------------------------------------------- |
---|
1229 | !********************************************************************** |
---|
1230 | !---------------------------------------------------------------------- |
---|
1231 | ! |
---|
1232 | DO K=1,LMH-1 |
---|
1233 | ELL=EL(K) |
---|
1234 | ! |
---|
1235 | ELOQ2=ELL*ELL/Q2(K) |
---|
1236 | ELOQ4=ELOQ2*ELOQ2 |
---|
1237 | ! |
---|
1238 | GML=GM(K) |
---|
1239 | GHL=GH(K) |
---|
1240 | ! |
---|
1241 | !---------------------------------------------------------------------- |
---|
1242 | !*** COEFFICIENTS OF THE TERMS IN THE DENOMINATOR |
---|
1243 | !---------------------------------------------------------------------- |
---|
1244 | ! |
---|
1245 | ADEN=(ADNM*GML+ADNH*GHL)*GHL |
---|
1246 | BDEN= BDNM*GML+BDNH*GHL |
---|
1247 | CDEN= 1. |
---|
1248 | ! |
---|
1249 | !---------------------------------------------------------------------- |
---|
1250 | !*** COEFFICIENTS FOR THE SM DETERMINANT |
---|
1251 | !---------------------------------------------------------------------- |
---|
1252 | ! |
---|
1253 | BESM=BSMH*GHL |
---|
1254 | ! |
---|
1255 | !---------------------------------------------------------------------- |
---|
1256 | !*** COEFFICIENTS FOR THE SH DETERMINANT |
---|
1257 | !---------------------------------------------------------------------- |
---|
1258 | ! |
---|
1259 | BESH=BSHM*GML+BSHH*GHL |
---|
1260 | ! |
---|
1261 | !---------------------------------------------------------------------- |
---|
1262 | !*** 1./DENOMINATOR |
---|
1263 | !---------------------------------------------------------------------- |
---|
1264 | ! |
---|
1265 | RDEN=1./(ADEN*ELOQ4+BDEN*ELOQ2+CDEN) |
---|
1266 | ! |
---|
1267 | !---------------------------------------------------------------------- |
---|
1268 | !*** SM AND SH |
---|
1269 | !---------------------------------------------------------------------- |
---|
1270 | ! |
---|
1271 | ESM=(BESM*ELOQ2+CESM)*RDEN |
---|
1272 | ESH=(BESH*ELOQ2+CESH)*RDEN |
---|
1273 | ! |
---|
1274 | !---------------------------------------------------------------------- |
---|
1275 | !*** DIFFUSION COEFFICIENTS |
---|
1276 | !---------------------------------------------------------------------- |
---|
1277 | ! |
---|
1278 | RDZ=2./(Z(K)-Z(K+2)) |
---|
1279 | Q1L=SQRT(Q2(K)) |
---|
1280 | ELQDZ=ELL*Q1L*RDZ |
---|
1281 | AKM(K)=ELQDZ*ESM |
---|
1282 | AKH(K)=ELQDZ*ESH |
---|
1283 | !---------------------------------------------------------------------- |
---|
1284 | ENDDO |
---|
1285 | !---------------------------------------------------------------------- |
---|
1286 | ! |
---|
1287 | !---------------------------------------------------------------------- |
---|
1288 | !*** INVERSIONS |
---|
1289 | !---------------------------------------------------------------------- |
---|
1290 | ! |
---|
1291 | ! IF(LMXL==LMH)THEN |
---|
1292 | ! KINV=LMH |
---|
1293 | ! D2Tmin=0. |
---|
1294 | ! |
---|
1295 | ! DO K=LMH/2,LMH-1 |
---|
1296 | ! D2T=T(K-1)-2.*T(K)+T(K+1) |
---|
1297 | ! IF(D2T<D2Tmin)THEN |
---|
1298 | ! D2Tmin=D2T |
---|
1299 | ! IF(D2T<0)KINV=K |
---|
1300 | ! ENDIF |
---|
1301 | ! ENDDO |
---|
1302 | ! |
---|
1303 | ! IF(KINV<LMH)THEN |
---|
1304 | ! DO K=KINV-1,LMH-1 |
---|
1305 | ! RDZ=2./(Z(K)-Z(K+2)) |
---|
1306 | ! AKMIN=0.5*RDZ |
---|
1307 | ! AKM(K)=MAX(AKM(K),AKMIN) |
---|
1308 | ! AKH(K)=MAX(AKH(K),AKMIN) |
---|
1309 | ! ENDDO |
---|
1310 | ! |
---|
1311 | !*** Begin debugging |
---|
1312 | ! IF(PRINT_DIAG>0)THEN |
---|
1313 | ! write(6,"(a,3i3)") '{turb1 lmxl,lmh,kinv=',lmxl,lmh,kinv |
---|
1314 | ! write(6,"(a,3i3)") '}turb1 lmxl,lmh,kinv=',lmxl,lmh,kinv |
---|
1315 | ! IF(PRINT_DIAG==1)THEN |
---|
1316 | ! write(6,"(a)") & |
---|
1317 | ! '{turb3 k, t, d2t, rdz, z(k), z(k+2), akmin, akh ' |
---|
1318 | ! ELSE |
---|
1319 | ! write(6,"(a)") & |
---|
1320 | ! '}turb3 k, t, d2t, rdz, z(k), z(k+2), akmin, akh ' |
---|
1321 | ! ENDIF |
---|
1322 | ! DO K=LMH-1,KINV-1,-1 |
---|
1323 | ! D2T=T(K-1)-2.*T(K)+T(K+1) |
---|
1324 | ! RDZ=2./(Z(K)-Z(K+2)) |
---|
1325 | ! AKMIN=0.5*RDZ |
---|
1326 | ! IF(PRINT_DIAG==1)THEN |
---|
1327 | ! write(6,"(a,i3,f8.3,2e12.5,2f9.2,2e12.5)") '{turb3 ' & |
---|
1328 | ! ,k,t(k)-273.15,d2t,rdz,z(k),z(k+2),akmin,akh(k) |
---|
1329 | ! ELSE |
---|
1330 | ! write(6,"(a,i3,f8.3,2e12.5,2f9.2,2e12.5)") '}turb3 ' & |
---|
1331 | ! ,k,t(k)-273.15,d2t,rdz,z(k),z(k+2),akmin,akh(k) |
---|
1332 | ! ENDIF |
---|
1333 | ! ENDDO |
---|
1334 | ! ENDIF !- IF (print_diag > 0) THEN |
---|
1335 | ! ENDIF !- IF(KINV<LMH)THEN |
---|
1336 | !*** End debugging |
---|
1337 | ! |
---|
1338 | ! ENDIF |
---|
1339 | !---------------------------------------------------------------------- |
---|
1340 | ! |
---|
1341 | END SUBROUTINE DIFCOF |
---|
1342 | ! |
---|
1343 | !---------------------------------------------------------------------- |
---|
1344 | !XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX |
---|
1345 | !---------------------------------------------------------------------- |
---|
1346 | SUBROUTINE VDIFQ & |
---|
1347 | ! ****************************************************************** |
---|
1348 | ! * * |
---|
1349 | ! * VERTICAL DIFFUSION OF Q2 (TKE) * |
---|
1350 | ! * * |
---|
1351 | ! ****************************************************************** |
---|
1352 | &(LMH,DTDIF,Q2,EL,Z & |
---|
1353 | ! Multi-layer UCM |
---|
1354 | &,A_E,B_E,SFK,VLK & |
---|
1355 | !UCM end |
---|
1356 | &,IDS,IDE,JDS,JDE,KDS,KDE & |
---|
1357 | &,IMS,IME,JMS,JME,KMS,KME & |
---|
1358 | &,ITS,ITE,JTS,JTE,KTS,KTE) |
---|
1359 | !---------------------------------------------------------------------- |
---|
1360 | ! |
---|
1361 | IMPLICIT NONE |
---|
1362 | ! |
---|
1363 | !---------------------------------------------------------------------- |
---|
1364 | INTEGER,INTENT(IN) :: IDS,IDE,JDS,JDE,KDS,KDE & |
---|
1365 | & ,IMS,IME,JMS,JME,KMS,KME & |
---|
1366 | & ,ITS,ITE,JTS,JTE,KTS,KTE |
---|
1367 | ! |
---|
1368 | INTEGER,INTENT(IN) :: LMH |
---|
1369 | ! |
---|
1370 | REAL,INTENT(IN) :: DTDIF |
---|
1371 | ! |
---|
1372 | REAL,DIMENSION(KTS:KTE-1),INTENT(IN) :: EL |
---|
1373 | REAL,DIMENSION(KTS:KTE+1),INTENT(IN) :: Z |
---|
1374 | ! |
---|
1375 | REAL,DIMENSION(KTS:KTE),INTENT(INOUT) :: Q2 |
---|
1376 | ! Multi-layer UCM |
---|
1377 | REAL,DIMENSION(KTS:KTE),INTENT(IN) :: A_E,B_E,VLK |
---|
1378 | REAL,DIMENSION(KTS:KTE),INTENT(IN) :: SFK |
---|
1379 | ! UCM end |
---|
1380 | !---------------------------------------------------------------------- |
---|
1381 | !*** |
---|
1382 | !*** LOCAL VARIABLES |
---|
1383 | !*** |
---|
1384 | INTEGER :: K |
---|
1385 | ! |
---|
1386 | REAL :: ADEN,AKQS,BDEN,BESH,BESM,CDEN,CF,DTOZS,ELL,ELOQ2,ELOQ4 & |
---|
1387 | & ,ELQDZ,ESH,ESM,ESQHF,GHL,GML,Q1L,RDEN,RDZ |
---|
1388 | ! |
---|
1389 | REAL,DIMENSION(KTS:KTE-2) :: AKQ,CM,CR,DTOZ,RSQ2 |
---|
1390 | !---------------------------------------------------------------------- |
---|
1391 | !********************************************************************** |
---|
1392 | !---------------------------------------------------------------------- |
---|
1393 | !*** |
---|
1394 | !*** VERTICAL TURBULENT DIFFUSION |
---|
1395 | !*** |
---|
1396 | !---------------------------------------------------------------------- |
---|
1397 | ESQHF=0.5*ESQ |
---|
1398 | ! |
---|
1399 | DO K=KTS,LMH-2 |
---|
1400 | DTOZ(K)=(DTDIF+DTDIF)/(Z(K)-Z(K+2))/VLK(K) !Bep |
---|
1401 | AKQ(K)=SFK(K)*SQRT((Q2(K)+Q2(K+1))*0.5)*(EL(K)+EL(K+1))*ESQHF & |
---|
1402 | & /(Z(K+1)-Z(K+2)) !Bep |
---|
1403 | CR(K)=-DTOZ(K)*AKQ(K) |
---|
1404 | ENDDO |
---|
1405 | ! |
---|
1406 | CM(1)=DTOZ(1)*AKQ(1)+1. |
---|
1407 | RSQ2(1)=Q2(1)+DTDIF*B_E(1) !Bep |
---|
1408 | ! |
---|
1409 | DO K=KTS+1,LMH-2 |
---|
1410 | CF=-DTOZ(K)*AKQ(K-1)/CM(K-1) |
---|
1411 | CM(K)=-CR(K-1)*CF+(AKQ(K-1)+AKQ(K))*DTOZ(K)+1. |
---|
1412 | CM(K)=-CR(K-1)*CF+(AKQ(K-1)+AKQ(K))*DTOZ(K)+1. |
---|
1413 | RSQ2(K)=-RSQ2(K-1)*CF+Q2(K)+DTDIF*B_E(K) !Bep |
---|
1414 | ENDDO |
---|
1415 | ! |
---|
1416 | DTOZS=(DTDIF+DTDIF)/(Z(LMH-1)-Z(LMH+1))/VLK(K) |
---|
1417 | AKQS=SFK(K)*SQRT((Q2(LMH-1)+Q2(LMH))*0.5)*(EL(LMH-1)+ELZ0)*ESQHF& |
---|
1418 | & /(Z(LMH)-Z(LMH+1)) |
---|
1419 | ! |
---|
1420 | CF=-DTOZS*AKQ(LMH-2)/CM(LMH-2) |
---|
1421 | ! |
---|
1422 | Q2(LMH-1)=(DTOZS*AKQS*Q2(LMH)-RSQ2(LMH-2)*CF+Q2(LMH-1)& |
---|
1423 | & +DTDIF*B_E(LMH-1)) & !Bep |
---|
1424 | & /((AKQ(LMH-2)+AKQS)*DTOZS-CR(LMH-2)*CF+1.) |
---|
1425 | ! |
---|
1426 | DO K=LMH-2,KTS,-1 |
---|
1427 | Q2(K)=(-CR(K)*Q2(K+1)+RSQ2(K))/CM(K) |
---|
1428 | ENDDO |
---|
1429 | !---------------------------------------------------------------------- |
---|
1430 | ! |
---|
1431 | END SUBROUTINE VDIFQ |
---|
1432 | ! |
---|
1433 | !---------------------------------------------------------------------- |
---|
1434 | !XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX |
---|
1435 | !--------------------------------------------------------------------- |
---|
1436 | SUBROUTINE VDIFH(DTDIF,LMH,THZ0,QZ0,RKHS,CHKLOWQ,CT & |
---|
1437 | & ,THE,Q,CWM,RKH_IN,Z,RHO & |
---|
1438 | ! Multi-layer UCM |
---|
1439 | & ,A_T,B_T,A_Q,B_Q,SFK,VLK & |
---|
1440 | ! UCM end |
---|
1441 | & ,IDS,IDE,JDS,JDE,KDS,KDE & |
---|
1442 | & ,IMS,IME,JMS,JME,KMS,KME & |
---|
1443 | & ,ITS,ITE,JTS,JTE,KTS,KTE,I,J) |
---|
1444 | ! *************************************************************** |
---|
1445 | ! * * |
---|
1446 | ! * VERTICAL DIFFUSION OF MASS VARIABLES * |
---|
1447 | ! * * |
---|
1448 | ! *************************************************************** |
---|
1449 | !--------------------------------------------------------------------- |
---|
1450 | ! |
---|
1451 | IMPLICIT NONE |
---|
1452 | ! |
---|
1453 | !--------------------------------------------------------------------- |
---|
1454 | INTEGER,INTENT(IN) :: IDS,IDE,JDS,JDE,KDS,KDE & |
---|
1455 | & ,IMS,IME,JMS,JME,KMS,KME & |
---|
1456 | & ,ITS,ITE,JTS,JTE,KTS,KTE,I,J |
---|
1457 | ! |
---|
1458 | INTEGER,INTENT(IN) :: LMH |
---|
1459 | ! |
---|
1460 | REAL,INTENT(IN) :: CHKLOWQ,CT,DTDIF,QZ0,RKHS,THZ0 |
---|
1461 | ! |
---|
1462 | REAL,DIMENSION(KTS:KTE-1),INTENT(IN) :: RKH_IN |
---|
1463 | REAL,DIMENSION(KTS:KTE),INTENT(IN) :: RHO |
---|
1464 | REAL,DIMENSION(KTS:KTE+1),INTENT(IN) :: Z |
---|
1465 | REAL,DIMENSION(KTS:KTE),INTENT(INOUT) :: CWM,Q,THE |
---|
1466 | ! Variables added for Bep |
---|
1467 | REAL,DIMENSION(KTS:KTE),INTENT(IN) :: A_T,B_T,A_Q,B_Q,VLK |
---|
1468 | REAL,DIMENSION(KTS:KTE),INTENT(IN) :: SFK |
---|
1469 | ! Bep end |
---|
1470 | ! |
---|
1471 | !---------------------------------------------------------------------- |
---|
1472 | !*** |
---|
1473 | !*** LOCAL VARIABLES |
---|
1474 | !*** |
---|
1475 | INTEGER :: K |
---|
1476 | ! |
---|
1477 | REAL :: CF,CMB,CMCB,CMQB,CMTB,CTHF,DTOZL,DTOZS & |
---|
1478 | & ,RCML,RKHH,RKQS,RSCB,RSQB,RSTB |
---|
1479 | ! Variables added for Bep |
---|
1480 | REAL :: CF_T,CF_Q,RCML_T,RCML_Q,CMB_T,CMB_Q |
---|
1481 | ! Bep end |
---|
1482 | ! |
---|
1483 | REAL,DIMENSION(KTS:KTE-1) :: CM,CR,DTOZ,RKCT,RSC,RSQ,RST |
---|
1484 | ! Variables added for Bep |
---|
1485 | REAL,DIMENSION(KTS:KTE-1) :: CM_T,CM_Q |
---|
1486 | REAL,DIMENSION(KTS:KTE-1) :: RKH |
---|
1487 | ! Bep end |
---|
1488 | ! |
---|
1489 | !---------------------------------------------------------------------- |
---|
1490 | !********************************************************************** |
---|
1491 | !---------------------------------------------------------------------- |
---|
1492 | CTHF=0.5*CT |
---|
1493 | ! |
---|
1494 | DO K=KTS,LMH-1 |
---|
1495 | RKH(K)=RKH_IN(K)*SFK(K) !Bep |
---|
1496 | DTOZ(K)=DTDIF/(Z(K)-Z(K+1))/VLK(K) !Bep |
---|
1497 | CR(K)=-DTOZ(K)*RKH(K) |
---|
1498 | RKCT(K)=RKH(K)*(Z(K)-Z(K+2))*CTHF |
---|
1499 | ENDDO |
---|
1500 | ! |
---|
1501 | CM(KTS)=DTOZ(KTS)*RKH(KTS)+RHO(KTS) |
---|
1502 | CM_T(KTS)=DTOZ(KTS)*RKH(KTS)+RHO(KTS)*(1.-DTDIF*A_T(KTS)) !Bep |
---|
1503 | CM_Q(KTS)=DTOZ(KTS)*RKH(KTS)+RHO(KTS)*(1.-DTDIF*A_Q(KTS)) !Bep |
---|
1504 | !---------------------------------------------------------------------- |
---|
1505 | RST(KTS) = -RKCT(KTS)*DTOZ(KTS)+& |
---|
1506 | &RHO(KTS)*(THE(KTS)+DTDIF*B_T(KTS)) !Bep |
---|
1507 | RSQ(KTS) = RHO(KTS)*( Q(KTS)+DTDIF*B_Q(KTS)) !Bep |
---|
1508 | RSC(KTS) = RHO(KTS)* CWM(KTS) |
---|
1509 | !---------------------------------------------------------------------- |
---|
1510 | DO K=KTS+1,LMH-1 |
---|
1511 | DTOZL=DTOZ(K) |
---|
1512 | CF = -DTOZL*RKH(K-1)/CM(K-1) |
---|
1513 | CF_T = -DTOZL*RKH(K-1)/CM_T(K-1) !Bep |
---|
1514 | CF_Q = -DTOZL*RKH(K-1)/CM_Q(K-1) !Bep |
---|
1515 | CM(K) = -CR(K-1)*CF+ (RKH(K-1)+RKH(K))*DTOZL+RHO(K) |
---|
1516 | CM_T(K)= -CR(K-1)*CF_T+(RKH(K-1)+RKH(K))*DTOZL+RHO(K)*& |
---|
1517 | &(1.-DTDIF*A_T(K)) |
---|
1518 | CM_Q(K)= -CR(K-1)*CF_Q+(RKH(K-1)+RKH(K))*DTOZL+RHO(K)*& |
---|
1519 | &(1.-DTDIF*A_Q(K)) !Bep |
---|
1520 | RST(K)=-RST(K-1)*CF_T+(RKCT(K-1)-RKCT(K))*DTOZL+RHO(K)*& |
---|
1521 | &(THE(K)+DTDIF*B_T(K)) !Bep |
---|
1522 | RSQ(K)=-RSQ(K-1)*CF_Q+RHO(K)*( Q(K)+DTDIF*B_Q(K)) !Bep |
---|
1523 | RSC(K)=-RSC(K-1)*CF+CWM(K)*RHO(K) |
---|
1524 | ENDDO |
---|
1525 | ! |
---|
1526 | DTOZS=DTDIF/(Z(LMH)-Z(LMH+1))/VLK(LMH) !Bep |
---|
1527 | RKHH=RKH(LMH-1) |
---|
1528 | ! |
---|
1529 | CF=-DTOZS*RKHH/CM(LMH-1) |
---|
1530 | CF_T=-DTOZS*RKHH/CM_T(LMH-1) !Bep |
---|
1531 | CF_Q=-DTOZS*RKHH/CM_Q(LMH-1) !Bep |
---|
1532 | RKQS=RKHS*CHKLOWQ |
---|
1533 | ! |
---|
1534 | CMB=CR(LMH-1)*CF |
---|
1535 | CMB_T=CR(LMH-1)*CF_T |
---|
1536 | CMB_Q=CR(LMH-1)*CF_Q |
---|
1537 | CMTB=-CMB_T+RKHH*DTOZS+RHO(LMH)*(1.-DTDIF*A_T(LMH)) !Bep |
---|
1538 | CMQB=-CMB_Q+RKHH*DTOZS+RHO(LMH)*(1.-DTDIF*A_Q(LMH)) !Bep |
---|
1539 | CMCB=-CMB+(RKHH )*DTOZS+RHO(LMH) |
---|
1540 | ! |
---|
1541 | RSTB=-RST(LMH-1)*CF_T+RHO(LMH)*(THE(LMH)+DTDIF*B_T(LMH))+& |
---|
1542 | &RKCT(LMH-1)*DTOZS !Bep |
---|
1543 | RSQB=-RSQ(LMH-1)*CF_Q+RHO(LMH)*( Q(LMH)+DTDIF*B_Q(LMH)) !Bep |
---|
1544 | RSCB=-RSC(LMH-1)*CF+CWM(LMH)*RHO(LMH) |
---|
1545 | !---------------------------------------------------------------------- |
---|
1546 | THE(LMH)= RSTB/CMTB |
---|
1547 | Q(LMH) = RSQB/CMQB |
---|
1548 | CWM(LMH)= RSCB/CMCB |
---|
1549 | !---------------------------------------------------------------------- |
---|
1550 | DO K=LMH-1,KTS,-1 |
---|
1551 | RCML=1./CM(K) |
---|
1552 | RCML_T=1./CM_T(K) !Bep |
---|
1553 | RCML_Q=1./CM_Q(K) !Bep |
---|
1554 | THE(K)=(-CR(K)*THE(K+1)+RST(K))*RCML_T |
---|
1555 | Q(K) =(-CR(K)* Q(K+1)+RSQ(K))*RCML_Q |
---|
1556 | CWM(K)=(-CR(K)*CWM(K+1)+RSC(K))*RCML |
---|
1557 | ENDDO |
---|
1558 | !---------------------------------------------------------------------- |
---|
1559 | ! |
---|
1560 | |
---|
1561 | |
---|
1562 | END SUBROUTINE VDIFH |
---|
1563 | ! |
---|
1564 | !---------------------------------------------------------------------- |
---|
1565 | !XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX |
---|
1566 | !--------------------------------------------------------------------- |
---|
1567 | SUBROUTINE VDIFX(DTDIF,LMH,RKHS,CT & |
---|
1568 | ,DUST1,DUST2,RKH,Z,RHO & |
---|
1569 | ,IDS,IDE,JDS,JDE,KDS,KDE & |
---|
1570 | ,IMS,IME,JMS,JME,KMS,KME & |
---|
1571 | ,ITS,ITE,JTS,JTE,KTS,KTE) |
---|
1572 | ! *************************************************************** |
---|
1573 | ! * * |
---|
1574 | ! * VERTICAL DIFFUSION OF MASS VARIABLES * |
---|
1575 | ! * * |
---|
1576 | ! *************************************************************** |
---|
1577 | !--------------------------------------------------------------------- |
---|
1578 | ! |
---|
1579 | IMPLICIT NONE |
---|
1580 | ! |
---|
1581 | !--------------------------------------------------------------------- |
---|
1582 | INTEGER,INTENT(IN) :: IDS,IDE,JDS,JDE,KDS,KDE & |
---|
1583 | ,IMS,IME,JMS,JME,KMS,KME & |
---|
1584 | ,ITS,ITE,JTS,JTE,KTS,KTE |
---|
1585 | ! |
---|
1586 | INTEGER,INTENT(IN) :: LMH |
---|
1587 | ! |
---|
1588 | REAL,INTENT(IN) :: CT,DTDIF,RKHS |
---|
1589 | ! |
---|
1590 | REAL,DIMENSION(KTS:KTE-1),INTENT(IN) :: RKH |
---|
1591 | REAL,DIMENSION(KTS:KTE),INTENT(IN) :: RHO |
---|
1592 | REAL,DIMENSION(KTS:KTE+1),INTENT(IN) :: Z |
---|
1593 | REAL,DIMENSION(KTS:KTE),INTENT(INOUT) :: DUST1,DUST2 |
---|
1594 | ! |
---|
1595 | !---------------------------------------------------------------------- |
---|
1596 | !*** |
---|
1597 | !*** LOCAL VARIABLES |
---|
1598 | !*** |
---|
1599 | INTEGER :: K |
---|
1600 | ! |
---|
1601 | REAL :: CF,CMB,CMDB,CTHF,DTOZL,DTOZS & |
---|
1602 | ,RCML,RKHH,RSD1B,RSD2B |
---|
1603 | ! |
---|
1604 | REAL,DIMENSION(KTS:KTE-1) :: CM,CR,DTOZ,RKCT,RSD1,RSD2 |
---|
1605 | ! |
---|
1606 | !---------------------------------------------------------------------- |
---|
1607 | !********************************************************************** |
---|
1608 | !---------------------------------------------------------------------- |
---|
1609 | CTHF=0.5*CT |
---|
1610 | ! |
---|
1611 | DO K=KTS,LMH-1 |
---|
1612 | DTOZ(K)=DTDIF/(Z(K)-Z(K+1)) |
---|
1613 | CR(K)=-DTOZ(K)*RKH(K) |
---|
1614 | RKCT(K)=RKH(K)*(Z(K)-Z(K+2))*CTHF |
---|
1615 | ENDDO |
---|
1616 | ! |
---|
1617 | CM(KTS)=DTOZ(KTS)*RKH(KTS)+RHO(KTS) |
---|
1618 | !---------------------------------------------------------------------- |
---|
1619 | RSD1(KTS)=DUST1(KTS)*RHO(KTS) |
---|
1620 | RSD2(KTS)=DUST2(KTS)*RHO(KTS) |
---|
1621 | !---------------------------------------------------------------------- |
---|
1622 | DO K=KTS+1,LMH-1 |
---|
1623 | DTOZL=DTOZ(K) |
---|
1624 | CF=-DTOZL*RKH(K-1)/CM(K-1) |
---|
1625 | CM(K)=-CR(K-1)*CF+(RKH(K-1)+RKH(K))*DTOZL+RHO(K) |
---|
1626 | RSD1(K)=-RSD1(K-1)*CF+DUST1(K)*RHO(K) |
---|
1627 | RSD2(K)=-RSD2(K-1)*CF+DUST2(K)*RHO(K) |
---|
1628 | ENDDO |
---|
1629 | ! |
---|
1630 | DTOZS=DTDIF/(Z(LMH)-Z(LMH+1)) |
---|
1631 | RKHH=RKH(LMH-1) |
---|
1632 | ! |
---|
1633 | CF=-DTOZS*RKHH/CM(LMH-1) |
---|
1634 | ! |
---|
1635 | CMB=CR(LMH-1)*CF |
---|
1636 | CMDB=-CMB+RKHH*DTOZS+RHO(LMH) |
---|
1637 | ! |
---|
1638 | RSD1B=-RSD1(LMH-1)*CF+DUST1(LMH)*RHO(LMH) |
---|
1639 | RSD2B=-RSD2(LMH-1)*CF+DUST2(LMH)*RHO(LMH) |
---|
1640 | !---------------------------------------------------------------------- |
---|
1641 | DUST1(LMH)=RSD1B/CMDB |
---|
1642 | DUST2(LMH)=RSD2B/CMDB |
---|
1643 | !---------------------------------------------------------------------- |
---|
1644 | DO K=LMH-1,KTS,-1 |
---|
1645 | RCML=1./CM(K) |
---|
1646 | DUST1(K)=(-CR(K)*DUST1(K+1)+RSD1(K))*RCML |
---|
1647 | DUST2(K)=(-CR(K)*DUST2(K+1)+RSD2(K))*RCML |
---|
1648 | ENDDO |
---|
1649 | !---------------------------------------------------------------------- |
---|
1650 | ! |
---|
1651 | END SUBROUTINE VDIFX |
---|
1652 | |
---|
1653 | !--------------------------------------------------------------------- |
---|
1654 | !XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX |
---|
1655 | !--------------------------------------------------------------------- |
---|
1656 | SUBROUTINE VDIFV(LMH,DTDIF,UZ0,VZ0,RKMS,U,V,RKM_IN,Z,RHO & |
---|
1657 | ! Variables added for Bep |
---|
1658 | & ,A_UK,A_VK,B_UK,B_VK,SFK,VLK & |
---|
1659 | ! Bep end |
---|
1660 | & ,IDS,IDE,JDS,JDE,KDS,KDE & |
---|
1661 | & ,IMS,IME,JMS,JME,KMS,KME & |
---|
1662 | ,ITS,ITE,JTS,JTE,KTS,KTE,I,J) |
---|
1663 | ! *************************************************************** |
---|
1664 | ! * * |
---|
1665 | ! * VERTICAL DIFFUSION OF VELOCITY COMPONENTS * |
---|
1666 | ! * * |
---|
1667 | ! *************************************************************** |
---|
1668 | !--------------------------------------------------------------------- |
---|
1669 | ! |
---|
1670 | IMPLICIT NONE |
---|
1671 | ! |
---|
1672 | !--------------------------------------------------------------------- |
---|
1673 | INTEGER,INTENT(IN) :: IDS,IDE,JDS,JDE,KDS,KDE & |
---|
1674 | & ,IMS,IME,JMS,JME,KMS,KME & |
---|
1675 | & ,ITS,ITE,JTS,JTE,KTS,KTE,I,J |
---|
1676 | ! |
---|
1677 | INTEGER,INTENT(IN) :: LMH |
---|
1678 | ! |
---|
1679 | REAL,INTENT(IN) :: RKMS,DTDIF,UZ0,VZ0 |
---|
1680 | ! |
---|
1681 | REAL,DIMENSION(KTS:KTE-1),INTENT(IN) :: RKM_IN !Bep |
---|
1682 | REAL,DIMENSION(KTS:KTE),INTENT(IN) :: RHO |
---|
1683 | REAL,DIMENSION(KTS:KTE+1),INTENT(IN) :: Z |
---|
1684 | ! |
---|
1685 | REAL,DIMENSION(KTS:KTE),INTENT(INOUT) :: U,V |
---|
1686 | REAL,DIMENSION(KTS:KTE),INTENT(IN) :: A_UK,A_VK,B_UK,B_VK,VLK !Bep |
---|
1687 | REAL,DIMENSION(KTS:KTE),INTENT(IN) :: SFK !Bep |
---|
1688 | !---------------------------------------------------------------------- |
---|
1689 | !*** |
---|
1690 | !*** LOCAL VARIABLES |
---|
1691 | !*** |
---|
1692 | INTEGER :: K |
---|
1693 | ! |
---|
1694 | REAL :: CFU,RCMLU,RCMVBU |
---|
1695 | REAL :: CFV,RCMLV,RCMVBV |
---|
1696 | REAL :: DTOZAK,DTOZL,DTOZS,RHOK,RKMH |
---|
1697 | ! |
---|
1698 | REAL,DIMENSION(KTS:KTE-1) :: CMU,CMV,CR,DTOZ,RSU,RSV |
---|
1699 | REAL,DIMENSION(KTS:KTE-1) :: RKM !Bep |
---|
1700 | !---------------------------------------------------------------------- |
---|
1701 | !********************************************************************** |
---|
1702 | !---------------------------------------------------------------------- |
---|
1703 | DO K=1,LMH-1 |
---|
1704 | RKM(K)=RKM_IN(K)*SFK(K) |
---|
1705 | DTOZ(K)=DTDIF/(Z(K)-Z(K+1))/VLK(K) |
---|
1706 | CR(K)=-DTOZ(K)*RKM(K) !alpha, gamma |
---|
1707 | ENDDO |
---|
1708 | ! |
---|
1709 | RHOK=RHO(1) |
---|
1710 | CMU(1)=DTOZ(1)*RKM(1)+RHOK*(1.-DTDIF*A_UK(1)) !Bep change beta |
---|
1711 | CMV(1)=DTOZ(1)*RKM(1)+RHOK*(1.-DTDIF*A_VK(1)) !Bep change beta |
---|
1712 | RSU(1)=RHOK*(U(1)+DTDIF*B_UK(1)) !Bep change W |
---|
1713 | RSV(1)=RHOK*(V(1)+DTDIF*B_VK(1)) !Bep change W |
---|
1714 | !---------------------------------------------------------------------- |
---|
1715 | DO K=2,LMH-1 |
---|
1716 | DTOZL=DTOZ(K) |
---|
1717 | CFU=-DTOZL*RKM(K-1)/CMU(K-1) |
---|
1718 | CFV=-DTOZL*RKM(K-1)/CMV(K-1) |
---|
1719 | RHOK=RHO(K) |
---|
1720 | CMU(K)=-CR(K-1)*CFU+(RKM(K-1)+RKM(K))& |
---|
1721 | &*DTOZL+RHOK*(1.-DTDIF*A_UK(K)) !Bep |
---|
1722 | CMV(K)=-CR(K-1)*CFV+(RKM(K-1)+RKM(K))& |
---|
1723 | &*DTOZL+RHOK*(1.-DTDIF*A_VK(K)) !Bep |
---|
1724 | RSU(K)=-RSU(K-1)*CFU+RHOK*(U(K)+DTDIF*B_UK(K)) !Bep |
---|
1725 | RSV(K)=-RSV(K-1)*CFV+RHOK*(V(K)+DTDIF*B_VK(K)) !Bep |
---|
1726 | ENDDO |
---|
1727 | !---------------------------------------------------------------------- |
---|
1728 | ! Bep: change lower BC from value UZ0 at LMH+1 to flux derived BC: Z(LMH+1)=0. |
---|
1729 | DTOZS=DTDIF/(Z(LMH)-Z(LMH+1))/VLK(LMH) !Bep |
---|
1730 | RKMH=RKM(LMH-1) !Bep |
---|
1731 | ! |
---|
1732 | CFU=-DTOZS*RKMH/CMU(LMH-1) |
---|
1733 | CFV=-DTOZS*RKMH/CMV(LMH-1) |
---|
1734 | RHOK=RHO(LMH) |
---|
1735 | RCMVBU=1./(RKMH*DTOZS-CR(LMH-1)*CFU+RHOK*(1.-DTDIF*A_UK(LMH))) !Bep |
---|
1736 | RCMVBV=1./(RKMH*DTOZS-CR(LMH-1)*CFV+RHOK*(1.-DTDIF*A_VK(LMH))) !Bep |
---|
1737 | |
---|
1738 | !---------------------------------------------------------------------- |
---|
1739 | U(LMH)=(-RSU(LMH-1)*CFU+RHOK*(U(LMH)+DTDIF*B_UK(LMH)))*RCMVBU !Bep |
---|
1740 | V(LMH)=(-RSV(LMH-1)*CFV+RHOK*(V(LMH)+DTDIF*B_VK(LMH)))*RCMVBV !Bep |
---|
1741 | !---------------------------------------------------------------------- |
---|
1742 | DO K=LMH-1,1,-1 |
---|
1743 | RCMLU=1./CMU(K) !Bep |
---|
1744 | RCMLV=1./CMV(K) !Bep |
---|
1745 | U(K)=(-CR(K)*U(K+1)+RSU(K))*RCMLU !Bep |
---|
1746 | V(K)=(-CR(K)*V(K+1)+RSV(K))*RCMLV !Bep |
---|
1747 | ENDDO |
---|
1748 | !---------------------------------------------------------------------- |
---|
1749 | ! |
---|
1750 | |
---|
1751 | END SUBROUTINE VDIFV |
---|
1752 | ! |
---|
1753 | !----------------------------------------------------------------------- |
---|
1754 | !XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX |
---|
1755 | !======================================================================= |
---|
1756 | SUBROUTINE MYJURBINIT(RUBLTEN,RVBLTEN,RTHBLTEN,RQVBLTEN, & |
---|
1757 | & TKE_MYJ,EXCH_H,RESTART,ALLOWED_TO_READ, & |
---|
1758 | & IDS,IDE,JDS,JDE,KDS,KDE, & |
---|
1759 | & IMS,IME,JMS,JME,KMS,KME, & |
---|
1760 | & ITS,ITE,JTS,JTE,KTS,KTE ) |
---|
1761 | !----------------------------------------------------------------------- |
---|
1762 | IMPLICIT NONE |
---|
1763 | !----------------------------------------------------------------------- |
---|
1764 | LOGICAL,INTENT(IN) :: ALLOWED_TO_READ,RESTART |
---|
1765 | INTEGER,INTENT(IN) :: IDS,IDE,JDS,JDE,KDS,KDE, & |
---|
1766 | & IMS,IME,JMS,JME,KMS,KME, & |
---|
1767 | & ITS,ITE,JTS,JTE,KTS,KTE |
---|
1768 | |
---|
1769 | REAL,DIMENSION(IMS:IME,KMS:KME,JMS:JME),INTENT(OUT) :: EXCH_H, & |
---|
1770 | & RUBLTEN, & |
---|
1771 | & RVBLTEN, & |
---|
1772 | & RTHBLTEN, & |
---|
1773 | & RQVBLTEN, & |
---|
1774 | & TKE_MYJ |
---|
1775 | INTEGER :: I,J,K,ITF,JTF,KTF |
---|
1776 | !----------------------------------------------------------------------- |
---|
1777 | !----------------------------------------------------------------------- |
---|
1778 | |
---|
1779 | JTF=MIN0(JTE,JDE-1) |
---|
1780 | KTF=MIN0(KTE,KDE-1) |
---|
1781 | ITF=MIN0(ITE,IDE-1) |
---|
1782 | |
---|
1783 | IF(.NOT.RESTART)THEN |
---|
1784 | DO J=JTS,JTF |
---|
1785 | DO K=KTS,KTF |
---|
1786 | DO I=ITS,ITF |
---|
1787 | TKE_MYJ(I,K,J)=EPSQ2 |
---|
1788 | RUBLTEN(I,K,J)=0. |
---|
1789 | RVBLTEN(I,K,J)=0. |
---|
1790 | RTHBLTEN(I,K,J)=0. |
---|
1791 | RQVBLTEN(I,K,J)=0. |
---|
1792 | EXCH_H(I,K,J)=0. |
---|
1793 | ENDDO |
---|
1794 | ENDDO |
---|
1795 | ENDDO |
---|
1796 | ENDIF |
---|
1797 | |
---|
1798 | END SUBROUTINE MYJURBINIT |
---|
1799 | !----------------------------------------------------------------------- |
---|
1800 | ! |
---|
1801 | END MODULE MODULE_BL_MYJURB |
---|
1802 | ! |
---|
1803 | !----------------------------------------------------------------------- |
---|