Context Navigation

← Previous Revision
Latest Revision
Next Revision →
Normal
Revision Log

module_bl_myjpbl.F @ 3026

Last change on this file since 3026 was 2759, checked in by aslmd, 3 years ago
adding unmodified code from WRFV3.0.1.1, expurged from useless data +1M size
Property svn:executable set to ``*
File size: 58.5 KB

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

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

Original Format