Context Navigation

module_ADVECTION.F

Last change on this file was 1, checked in by lfita, 10 years ago

-- --- Opening of the WRF+LMDZ coupling repository --- -- -

WRF: version v3.3
LMDZ: version v1818

More details in:

File size: 129.0 KB

Line
1	!----------------------------------------------------------------------
2	!#define BIT_FOR_BIT
3	!----------------------------------------------------------------------
4	#include "nmm_loop_basemacros.h"
5	#include "nmm_loop_macros.h"
6	!----------------------------------------------------------------------
7	!
8	!NCEP_MESO:MODEL_LAYER: HORIZONTAL AND VERTICAL ADVECTION
9	!
10	!----------------------------------------------------------------------
11	!
12	MODULE MODULE_ADVECTION
13	!
14	!----------------------------------------------------------------------
15	USE MODULE_MODEL_CONSTANTS
16	USE MODULE_EXT_INTERNAL
17	!----------------------------------------------------------------------
18	#if defined(DM_PARALLEL) && !defined(STUBMPI)
19	INCLUDE "mpif.h"
20	#endif
21	!----------------------------------------------------------------------
22	!
23	REAL,PARAMETER :: FF2=-0.64813,FF3=0.24520,FF4=-0.12189
24	REAL,PARAMETER :: FFC=1.533,FBC=1.-FFC
25	REAL :: CONSERVE_MIN=0.9,CONSERVE_MAX=1.1
26	!
27	!----------------------------------------------------------------------
28	!*** CRANK-NICHOLSON OFF-CENTER WEIGHTS FOR CURRENT AND FUTURE
29	!*** TIME LEVELS.
30	!-----------------------------------------------------------------------
31	!
32	REAL,PARAMETER :: WGT1=0.90
33	REAL,PARAMETER :: WGT2=2.-WGT1
34	!
35	!*** FOR CRANK_NICHOLSON CHECK ONLY.
36	!
37	INTEGER :: ITEST=47,JTEST=70
38	REAL :: ADTP,ADUP,ADVP,TTLO,TTUP,TULO,TUUP,TVLO,TVUP
39	!
40	!----------------------------------------------------------------------
41	CONTAINS
42	!
43	!***********************************************************************
44	SUBROUTINE ADVE(NTSD,DT,DETA1,DETA2,PDTOP &
45	& ,CURV,F,FAD,F4D,EM_LOC,EMT_LOC,EN,ENT,DX,DY &
46	& ,HBM2,VBM2 &
47	& ,T,U,V,PDSLO,TOLD,UOLD,VOLD &
48	& ,PETDT,UPSTRM &
49	& ,FEW,FNS,FNE,FSE &
50	& ,ADT,ADU,ADV &
51	& ,N_IUP_H,N_IUP_V &
52	& ,N_IUP_ADH,N_IUP_ADV &
53	& ,IUP_H,IUP_V,IUP_ADH,IUP_ADV &
54	& ,IHE,IHW,IVE,IVW &
55	& ,IDS,IDE,JDS,JDE,KDS,KDE &
56	& ,IMS,IME,JMS,JME,KMS,KME &
57	& ,ITS,ITE,JTS,JTE,KTS,KTE)
58	!***********************************************************************
59	!$$$ SUBPROGRAM DOCUMENTATION BLOCK
60	! . . .
61	! SUBPROGRAM: ADVE HORIZONTAL AND VERTICAL ADVECTION
62	! PRGRMMR: JANJIC ORG: W/NP22 DATE: 93-10-28
63	!
64	! ABSTRACT:
65	! ADVE CALCULATES THE CONTRIBUTION OF THE HORIZONTAL AND VERTICAL
66	! ADVECTION TO THE TENDENCIES OF TEMPERATURE AND WIND AND THEN
67	! UPDATES THOSE VARIABLES.
68	! THE JANJIC ADVECTION SCHEME FOR THE ARAKAWA E GRID IS USED
69	! FOR ALL VARIABLES INSIDE THE FIFTH ROW. AN UPSTREAM SCHEME
70	! IS USED ON ALL VARIABLES IN THE THIRD, FOURTH, AND FIFTH
71	! OUTERMOST ROWS. THE ADAMS-BASHFORTH TIME SCHEME IS USED.
72	!
73	! PROGRAM HISTORY LOG:
74	! 87-06-?? JANJIC - ORIGINATOR
75	! 95-03-25 BLACK - CONVERSION FROM 1-D TO 2-D IN HORIZONTAL
76	! 96-03-28 BLACK - ADDED EXTERNAL EDGE
77	! 98-10-30 BLACK - MODIFIED FOR DISTRIBUTED MEMORY
78	! 99-07- JANJIC - CONVERTED TO ADAMS-BASHFORTH SCHEME
79	! COMBINING HORIZONTAL AND VERTICAL ADVECTION
80	! 02-02-04 BLACK - ADDED VERTICAL CFL CHECK
81	! 02-02-05 BLACK - CONVERTED TO WRF FORMAT
82	! 02-08-29 MICHALAKES - CONDITIONAL COMPILATION OF MPI
83	! CONVERT TO GLOBAL INDEXING
84	! 02-09-06 WOLFE - MORE CONVERSION TO GLOBAL INDEXING
85	! 04-05-29 JANJIC,BLACK - CRANK-NICHOLSON VERTICAL ADVECTION
86	! 04-11-23 BLACK - THREADED
87	! 05-12-14 BLACK - CONVERTED FROM IKJ TO IJK
88	!
89	! USAGE: CALL ADVE FROM SUBROUTINE SOLVE_NMM
90	! INPUT ARGUMENT LIST:
91	!
92	! OUTPUT ARGUMENT LIST:
93	!
94	! OUTPUT FILES:
95	! NONE
96	!
97	! SUBPROGRAMS CALLED:
98	!
99	! UNIQUE: NONE
100	!
101	! LIBRARY: NONE
102	!
103	! ATTRIBUTES:
104	! LANGUAGE: FORTRAN 90
105	! MACHINE : IBM SP
106	!$$$
107	!***********************************************************************
108	!-----------------------------------------------------------------------
109	!
110	IMPLICIT NONE
111	!
112	!-----------------------------------------------------------------------
113	!
114	INTEGER,INTENT(IN) :: IDS,IDE,JDS,JDE,KDS,KDE &
115	& ,IMS,IME,JMS,JME,KMS,KME &
116	& ,ITS,ITE,JTS,JTE,KTS,KTE
117	!
118	INTEGER, DIMENSION(JMS:JME),INTENT(IN) :: IHE,IHW,IVE,IVW &
119	,N_IUP_H,N_IUP_V &
120	& ,N_IUP_ADH,N_IUP_ADV
121	!
122	INTEGER, DIMENSION(IMS:IME,JMS:JME),INTENT(IN) :: IUP_H,IUP_V &
123	& ,IUP_ADH,IUP_ADV
124	!
125	INTEGER,INTENT(IN) :: NTSD
126	!
127	REAL,INTENT(IN) :: DT,DY,EN,ENT,F4D,PDTOP
128	!
129	REAL,DIMENSION(NMM_MAX_DIM),INTENT(IN) :: EM_LOC,EMT_LOC
130	!
131	REAL,DIMENSION(KMS:KME),INTENT(IN) :: DETA1,DETA2
132	!
133	REAL,DIMENSION(IMS:IME,JMS:JME),INTENT(IN) :: CURV,DX,F,FAD,HBM2 &
134	& ,PDSLO,VBM2
135	!
136	REAL,DIMENSION(IMS:IME,JMS:JME),INTENT(OUT) :: ADT,ADU,ADV
137	!
138	REAL,DIMENSION(IMS:IME,JMS:JME,KMS:KME),INTENT(IN) :: PETDT
139	!
140	REAL,DIMENSION(IMS:IME,JMS:JME,KMS:KME),INTENT(INOUT) :: T,TOLD &
141	& ,U,UOLD &
142	& ,V,VOLD
143	!
144	REAL,DIMENSION(IMS:IME,JMS:JME,KMS:KME),INTENT(OUT) :: FEW,FNE &
145	& ,FNS,FSE
146	!
147	!-----------------------------------------------------------------------
148	!*** LOCAL VARIABLES
149	!-----------------------------------------------------------------------
150	!
151	LOGICAL :: UPSTRM
152	!
153	INTEGER :: I,IEND,IFP,IFQ,II,IPQ,ISP,ISQ,ISTART &
154	& ,IUP_ADH_J,IVH,IVL &
155	& ,J,J1,JA,JAK,JEND,JGLOBAL,JJ,JKNT,JP2,JSTART &
156	& ,K,KNTI_ADH,KSTART,KSTOP &
157	& ,N,N_IUPH_J,N_IUPADH_J,N_IUPADV_J
158	!
159	INTEGER :: MY_IS_GLB,MY_IE_GLB,MY_JS_GLB,MY_JE_GLB
160	!
161	INTEGER,DIMENSION(ITS-5:ITE+5,JTS-5:JTE+5) :: ISPA,ISQA
162	!
163	REAL :: ADPDX,ADPDY,ARRAY3_X,CFL,CFT,CFU,CFV,CMT,CMU,CMV &
164	& ,DTE,DTQ,F0,F1,F2,F3,FEWP,FNEP,FNSP,FPP,FSEP,HM &
165	& ,PDOP,PDOPU,PDOPV,PP &
166	& ,PVVLO,PVVLOU,PVVLOV,PVVUP,PVVUPU,PVVUPV &
167	& ,QP,RDP,RDPU,RDPV &
168	& ,TEMPA,TEMPB,TTA,TTB,UDY &
169	& ,VDX,VM,VVLO,VVLOU,VVLOV,VVUP,VVUPU,VVUPV
170	!
171	REAL,DIMENSION(ITS-5:ITE+5,JTS-5:JTE+5) :: ARRAY0,ARRAY1 &
172	& ,ARRAY2,ARRAY3 &
173	& ,DPDE,RDPD,RDPDX,RDPDY &
174	& ,TEW,TNE,TNS,TSE,TST &
175	& ,UNE,UNED,UEW,UNS,USE &
176	& ,USED,UST &
177	& ,VEW,VNE,VNS,VSE &
178	& ,VST
179	!
180	REAL,DIMENSION(ITS-5:ITE+5,JTS-5:JTE+5,KTS:KTE) :: VAD_TEND_T &
181	& ,VAD_TEND_U &
182	& ,VAD_TEND_V
183	!
184	REAL,DIMENSION(KTS:KTE) :: CRT,CRU,CRV,DETA1_PDTOP &
185	& ,RCMT,RCMU,RCMV,RSTT,RSTU,RSTV &
186	& ,T_K,TN,U_K,UN,V_K,VN
187	!
188	!-----------------------------------------------------------------------
189	!***********************************************************************
190	!
191	! DPDE ----- 3
192	! \| J Increasing
193	! \|
194	! \| ^
195	! FNS ----- 2 \|
196	! \| \|
197	! \| \|
198	! \| \|
199	! VNS ----- 1 \|
200	! \|
201	! \|
202	! \|
203	! ADV ----- 0 ------> Current J
204	! \|
205	! \|
206	! \|
207	! VNS ----- -1
208	! \|
209	! \|
210	! \|
211	! FNS ----- -2
212	! \|
213	! \|
214	! \|
215	! DPDE ----- -3
216	!
217	!***********************************************************************
218	!-----------------------------------------------------------------------
219	DO J=JTS-5,JTE+5
220	DO I=ITS-5,ITE+5
221	ARRAY0(I,J)=0.0
222	ARRAY1(I,J)=0.0
223	ARRAY2(I,J)=0.0
224	ARRAY3(I,J)=0.0
225	DPDE(I,J)=0.0
226	RDPD(I,J)=0.0
227	RDPDX(I,J)=0.0
228	RDPDY(I,J)=0.0
229	TEW(I,J)=0.0
230	TNE(I,J)=0.0
231	TNS(I,J)=0.0
232	TSE(I,J)=0.0
233	TST(I,J)=0.0
234	UNE(I,J)=0.0
235	UNED(I,J)=0.0
236	UEW(I,J)=0.0
237	UNS(I,J)=0.0
238	USE(I,J)=0.0
239	USED(I,J)=0.0
240	UST(I,J)=0.0
241	VEW(I,J)=0.0
242	VNE(I,J)=0.0
243	VNS(I,J)=0.0
244	VSE(I,J)=0.0
245	VST(I,J)=0.0
246	ENDDO
247	ENDDO
248	!-----------------------------------------------------------------------
249	!
250	DTQ=DT*0.25
251	DTE=DT(0.50.25)
252	!
253	!-----------------------------------------------------------------------
254	!***
255	!*** PRECOMPUTE DETA1 TIMES PDTOP.
256	!***
257	!-----------------------------------------------------------------------
258	!
259	DO K=KTS,KTE
260	DETA1_PDTOP(K)=DETA1(K)*PDTOP
261	ENDDO
262	!
263	!-----------------------------------------------------------------------
264	!***
265	!*** INITIALIZE SOME WORKING ARRAYS TO ZERO
266	!***
267	!
268	!-----------------------------------------------------------------------
269	!-----------------------------------------------------------------------
270	!
271	!*** COMPUTE VERTICAL ADVECTION TENDENCIES USING CRANK-NICHOLSON.
272	!
273	!-----------------------------------------------------------------------
274	!-----------------------------------------------------------------------
275	!
276	!-----------------------------------------------------------------------
277	!*** FIRST THE TEMPERATURE
278	!-----------------------------------------------------------------------
279	!$omp parallel do &
280	!$omp& private(cft,cfu,cfv,cmt,cmu,cmv,crt,cru,crv,i,k &
281	!$omp& ,pdop,pdopu,pdopv,pvvlo,pvvlou,pvvlov,pvvup,pvvupu,pvvupv &
282	!$omp& ,rcmt,rcmu,rcmv,rdp,rdpu,rdpv,rstt,rstu,rstv,t_k,tn &
283	!$omp& ,u_k,un,v_k,vn,vvlo,vvlou,vvlov,vvup,vvupu,vvupv)
284	!!$omp& private(adtp,adup,advp,ttlo,ttup,tulo,tuup,tvlo,tvup)
285	!-----------------------------------------------------------------------
286	!
287	main_vertical: DO J=MYJS2,MYJE2
288	!
289	!-----------------------------------------------------------------------
290	!
291	iloop_for_t: DO I=MYIS1,MYIE1
292	!
293	!-----------------------------------------------------------------------
294	!*** EXTRACT T FROM THE COLUMN
295	!-----------------------------------------------------------------------
296	!
297	DO K=KTS,KTE
298	T_K(K)=T(I,J,K)
299	ENDDO
300	!
301	!-----------------------------------------------------------------------
302	!
303	PDOP=PDSLO(I,J)
304	PVVLO=PETDT(I,J,KTE-1)*DTQ
305	VVLO=PVVLO/(DETA1_PDTOP(KTE)+DETA2(KTE)*PDOP)
306	CMT=-VVLO*WGT2+1.
307	RCMT(KTE)=1./CMT
308	CRT(KTE)=VVLO*WGT2
309	RSTT(KTE)=-VVLOWGT1(T_K(KTE-1)-T_K(KTE))+T_K(KTE)
310	!
311	!-----------------------------------------------------------------------
312	!
313	DO K=KTE-1,KTS+1,-1
314	RDP=1./(DETA1_PDTOP(K)+DETA2(K)*PDOP)
315	PVVUP=PVVLO
316	PVVLO=PETDT(I,J,K-1)*DTQ
317	VVUP=PVVUP*RDP
318	VVLO=PVVLO*RDP
319	CFT=-VVUPWGT2RCMT(K+1)
320	CMT=-CRT(K+1)CFT+((VVUP-VVLO)WGT2+1.)
321	RCMT(K)=1./CMT
322	CRT(K)=VVLO*WGT2
323	RSTT(K)=-RSTT(K+1)*CFT+T_K(K) &
324	& -(T_K(K)-T_K(K+1))VVUPWGT1 &
325	& -(T_K(K-1)-T_K(K))VVLOWGT1
326	ENDDO
327	!
328	!-----------------------------------------------------------------------
329	!
330	PVVUP=PVVLO
331	VVUP=PVVUP/(DETA1_PDTOP(KTS)+DETA2(KTS)*PDOP)
332	CFT=-VVUPWGT2RCMT(KTS+1)
333	CMT=-CRT(KTS+1)CFT+VVUPWGT2+1.
334	CRT(KTS)=0.
335	RSTT(KTS)=-(T_K(KTS)-T_K(KTS+1))VVUPWGT1 &
336	& -RSTT(KTS+1)*CFT+T_K(KTS)
337	TN(KTS)=RSTT(KTS)/CMT
338	VAD_TEND_T(I,J,KTS)=TN(KTS)-T_K(KTS)
339	!
340	DO K=KTS+1,KTE
341	TN(K)=(-CRT(K)TN(K-1)+RSTT(K))RCMT(K)
342	VAD_TEND_T(I,J,K)=TN(K)-T_K(K)
343	ENDDO
344	!
345	!-----------------------------------------------------------------------
346	!*** The following section is only for checking the implicit solution
347	!*** using back-substitution. Remove this section otherwise.
348	!-----------------------------------------------------------------------
349	! if(ntsd<=10.or.ntsd>=6000)then
350	! IF(I==ITEST.AND.J==JTEST)THEN
351	!!
352	! PVVLO=PETDT(I,J,KTE-1)DT0.25
353	! VVLO=PVVLO/(DETA1_PDTOP(KTE)+DETA2(KTE)*PDOP)
354	! TTLO=VVLO*(T(I,J,KTE-1)-T(I,J,KTE) &
355	! & +TN(KTE-1)-TN(KTE))
356	! ADTP=TTLO+TN(KTE)-T(I,J,KTE)
357	! WRITE(0,*)' NTSD=',NTSD,' I=',ITEST,' J=',JTEST,' K=',KTE &
358	! &, ' ADTP=',ADTP
359	! WRITE(0,*)' T=',T(I,J,KTE),' TN=',TN(KTE) &
360	! &, ' VAD_TEND_T=',VAD_TEND_T(I,J,KTE)
361	! WRITE(0,*)' '
362	!!
363	! DO K=KTE-1,KTS+1,-1
364	! RDP=1./(DETA1_PDTOP(K)+DETA2(K)*PDOP)
365	! PVVUP=PVVLO
366	! PVVLO=PETDT(I,J,K-1)DT0.25
367	! VVUP=PVVUP*RDP
368	! VVLO=PVVLO*RDP
369	! TTUP=VVUP*(T(I,J,K)-T(I,J,K+1)+TN(K)-TN(K+1))
370	! TTLO=VVLO*(T(I,J,K-1)-T(I,J,K)+TN(K-1)-TN(K))
371	! ADTP=TTLO+TTUP+TN(K)-T(I,J,K)
372	! WRITE(0,*)' NTSD=',NTSD,' I=',I,' J=',J,' K=',K &
373	! &, ' ADTP=',ADTP
374	! WRITE(0,*)' T=',T(I,J,K),' TN=',TN(K) &
375	! &, ' VAD_TEND_T=',VAD_TEND_T(I,J,K)
376	! WRITE(0,*)' '
377	! ENDDO
378	!!
379	! PVVUP=PVVLO
380	! VVUP=PVVUP/(DETA1_PDTOP(KTS)+DETA2(KTS)*PDOP)
381	! TTUP=VVUP*(T(I,J,KTS)-T(I,J,KTS+1)+TN(KTS)-TN(KTS+1))
382	! ADTP=TTUP+TN(KTS)-T(I,J,KTS)
383	! WRITE(0,*)' NTSD=',NTSD,' I=',I,' J=',J,' K=',KTS &
384	! &, ' ADTP=',ADTP
385	! WRITE(0,*)' T=',T(I,J,KTS),' TN=',TN(KTS) &
386	! &, ' VAD_TEND_T=',VAD_TEND_T(I,J,KTS)
387	! WRITE(0,*)' '
388	! ENDIF
389	! endif
390	!
391	!-----------------------------------------------------------------------
392	!*** End of check.
393	!-----------------------------------------------------------------------
394	!
395	ENDDO iloop_for_t
396	!
397	!-----------------------------------------------------------------------
398	!
399	!*** NOW VERTICAL ADVECTION OF WIND COMPONENTS
400	!
401	!-----------------------------------------------------------------------
402	!
403	iloop_for_uv: DO I=MYIS1,MYIE1
404	!
405	!-----------------------------------------------------------------------
406	!*** EXTRACT U AND V FROM THE COLUMN
407	!-----------------------------------------------------------------------
408	!
409	DO K=KTS,KTE
410	U_K(K)=U(I,J,K)
411	V_K(K)=V(I,J,K)
412	ENDDO
413	!
414	!-----------------------------------------------------------------------
415	!
416	PDOPU=(PDSLO(I+IVW(J),J)+PDSLO(I+IVE(J),J))*0.5
417	PDOPV=(PDSLO(I,J-1)+PDSLO(I,J+1))*0.5
418	PVVLOU=(PETDT(I+IVW(J),J,KTE-1)+PETDT(I+IVE(J),J,KTE-1))*DTE
419	PVVLOV=(PETDT(I,J-1,KTE-1)+PETDT(I,J+1,KTE-1))*DTE
420	VVLOU=PVVLOU/(DETA1_PDTOP(KTE)+DETA2(KTE)*PDOPU)
421	VVLOV=PVVLOV/(DETA1_PDTOP(KTE)+DETA2(KTE)*PDOPV)
422	CMU=-VVLOU*WGT2+1.
423	CMV=-VVLOV*WGT2+1.
424	RCMU(KTE)=1./CMU
425	RCMV(KTE)=1./CMV
426	CRU(KTE)=VVLOU*WGT2
427	CRV(KTE)=VVLOV*WGT2
428	RSTU(KTE)=-VVLOUWGT1(U_K(KTE-1)-U_K(KTE))+U_K(KTE)
429	RSTV(KTE)=-VVLOVWGT1(V_K(KTE-1)-V_K(KTE))+V_K(KTE)
430	!
431	!-----------------------------------------------------------------------
432	!
433	DO K=KTE-1,KTS+1,-1
434	RDPU=1./(DETA1_PDTOP(K)+DETA2(K)*PDOPU)
435	RDPV=1./(DETA1_PDTOP(K)+DETA2(K)*PDOPV)
436	PVVUPU=PVVLOU
437	PVVUPV=PVVLOV
438	PVVLOU=(PETDT(I+IVW(J),J,K-1)+PETDT(I+IVE(J),J,K-1))*DTE
439	PVVLOV=(PETDT(I,J-1,K-1)+PETDT(I,J+1,K-1))*DTE
440	VVUPU=PVVUPU*RDPU
441	VVUPV=PVVUPV*RDPV
442	VVLOU=PVVLOU*RDPU
443	VVLOV=PVVLOV*RDPV
444	CFU=-VVUPUWGT2RCMU(K+1)
445	CFV=-VVUPVWGT2RCMV(K+1)
446	CMU=-CRU(K+1)CFU+(VVUPU-VVLOU)WGT2+1.
447	CMV=-CRV(K+1)CFV+(VVUPV-VVLOV)WGT2+1.
448	RCMU(K)=1./CMU
449	RCMV(K)=1./CMV
450	CRU(K)=VVLOU*WGT2
451	CRV(K)=VVLOV*WGT2
452	RSTU(K)=-RSTU(K+1)*CFU+U_K(K) &
453	& -(U_K(K)-U_K(K+1))VVUPUWGT1 &
454	& -(U_K(K-1)-U_K(K))VVLOUWGT1
455	RSTV(K)=-RSTV(K+1)*CFV+V_K(K) &
456	& -(V_K(K)-V_K(K+1))VVUPVWGT1 &
457	& -(V_K(K-1)-V_K(K))VVLOVWGT1
458	ENDDO
459	!
460	!-----------------------------------------------------------------------
461	!
462	RDPU=1./(DETA1_PDTOP(KTS)+DETA2(KTS)*PDOPU)
463	RDPV=1./(DETA1_PDTOP(KTS)+DETA2(KTS)*PDOPV)
464	PVVUPU=PVVLOU
465	PVVUPV=PVVLOV
466	VVUPU=PVVUPU*RDPU
467	VVUPV=PVVUPV*RDPV
468	CFU=-VVUPUWGT2RCMU(KTS+1)
469	CFV=-VVUPVWGT2RCMV(KTS+1)
470	CMU=-CRU(KTS+1)CFU+VVUPUWGT2+1.
471	CMV=-CRV(KTS+1)CFV+VVUPVWGT2+1.
472	CRU(KTS)=0.
473	CRV(KTS)=0.
474	RSTU(KTS)=-(U_K(KTS)-U_K(KTS+1))VVUPUWGT1 &
475	& -RSTU(KTS+1)*CFU+U_K(KTS)
476	RSTV(KTS)=-(V_K(KTS)-V_K(KTS+1))VVUPVWGT1 &
477	& -RSTV(KTS+1)*CFV+V_K(KTS)
478	UN(KTS)=RSTU(KTS)/CMU
479	VN(KTS)=RSTV(KTS)/CMV
480	VAD_TEND_U(I,J,KTS)=UN(KTS)-U_K(KTS)
481	VAD_TEND_V(I,J,KTS)=VN(KTS)-V_K(KTS)
482	!
483	DO K=KTS+1,KTE
484	UN(K)=(-CRU(K)UN(K-1)+RSTU(K))RCMU(K)
485	VN(K)=(-CRV(K)VN(K-1)+RSTV(K))RCMV(K)
486	VAD_TEND_U(I,J,K)=UN(K)-U_K(K)
487	VAD_TEND_V(I,J,K)=VN(K)-V_K(K)
488	ENDDO
489	!
490	!-----------------------------------------------------------------------
491	!*** The following section is only for checking the implicit solution
492	!*** using back-substitution. Remove this section otherwise.
493	!-----------------------------------------------------------------------
494	!
495	! if(ntsd<=10.or.ntsd>=6000)then
496	! IF(I==ITEST.AND.J==JTEST)THEN
497	!!
498	! PDOPU=(PDSLO(I+IVW(J),J)+PDSLO(I+IVE(J),J))*0.5
499	! PDOPV=(PDSLO(I,J-1)+PDSLO(I,J+1))*0.5
500	! PVVLOU=(PETDT(I+IVW(J),J,KTE-1) &
501	! & +PETDT(I+IVE(J),J,KTE-1))*DTE
502	! PVVLOV=(PETDT(I,J-1,KTE-1) &
503	! & +PETDT(I,J+1,KTE-1))*DTE
504	! VVLOU=PVVLOU/(DETA1_PDTOP(KTE)+DETA2(KTE)*PDOPU)
505	! VVLOV=PVVLOV/(DETA1_PDTOP(KTE)+DETA2(KTE)*PDOPV)
506	! TULO=VVLOU*(U(I,J,KTE-1)-U(I,J,KTE)+UN(KTE-1)-UN(KTE))
507	! TVLO=VVLOV*(V(I,J,KTE-1)-V(I,J,KTE)+VN(KTE-1)-VN(KTE))
508	! ADUP=TULO+UN(KTE)-U(I,J,KTE)
509	! ADVP=TVLO+VN(KTE)-V(I,J,KTE)
510	! WRITE(0,*)' NTSD=',NTSD,' I=',I,' J=',J,' K=',KTE &
511	! &, ' ADUP=',ADUP,' ADVP=',ADVP
512	! WRITE(0,*)' U=',U(I,J,KTE),' UN=',UN(KTE) &
513	! &, ' VAD_TEND_U=',VAD_TEND_U(I,KTE) &
514	! &, ' V=',V(I,J,KTE),' VN=',VN(KTE) &
515	! &, ' VAD_TEND_V=',VAD_TEND_V(I,KTE)
516	! WRITE(0,*)' '
517	!!
518	! DO K=KTE-1,KTS+1,-1
519	! RDPU=1./(DETA1_PDTOP(K)+DETA2(K)*PDOPU)
520	! RDPV=1./(DETA1_PDTOP(K)+DETA2(K)*PDOPV)
521	! PVVUPU=PVVLOU
522	! PVVUPV=PVVLOV
523	! PVVLOU=(PETDT(I+IVW(J),J,K-1) &
524	! & +PETDT(I+IVE(J),J,K-1))*DTE
525	! PVVLOV=(PETDT(I,J-1,K-1)+PETDT(I,J+1,K-1))*DTE
526	! VVUPU=PVVUPU*RDPU
527	! VVUPV=PVVUPV*RDPV
528	! VVLOU=PVVLOU*RDPU
529	! VVLOV=PVVLOV*RDPV
530	! TUUP=VVUPU*(U(I,J,K)-U(I,J,K+1)+UN(K)-UN(K+1))
531	! TVUP=VVUPV*(V(I,J,K)-V(I,J,K+1)+VN(K)-VN(K+1))
532	! TULO=VVLOU*(U(I,J,K-1)-U(I,J,K)+UN(K-1)-UN(K))
533	! TVLO=VVLOV*(V(I,J,K-1)-V(I,J,K)+VN(K-1)-VN(K))
534	! ADUP=TUUP+TULO+UN(K)-U(I,J,K)
535	! ADVP=TVUP+TVLO+VN(K)-V(I,J,K)
536	! WRITE(0,*)' NTSD=',NTSD,' I=',ITEST,' J=',JTEST,' K=',K &
537	! &, ' ADUP=',ADUP,' ADVP=',ADVP
538	! WRITE(0,*)' U=',U(I,J,K),' UN=',UN(K) &
539	! &, ' VAD_TEND_U=',VAD_TEND_U(I,K) &
540	! &, ' V=',V(I,J,K),' VN=',VN(K) &
541	! &, ' VAD_TEND_V=',VAD_TEND_V(I,K)
542	! WRITE(0,*)' '
543	! ENDDO
544	!!
545	! PVVUPU=PVVLOU
546	! PVVUPV=PVVLOV
547	! VVUPU=PVVUPU/(DETA1_PDTOP(KTS)+DETA2(KTS)*PDOPU)
548	! VVUPV=PVVUPV/(DETA1_PDTOP(KTS)+DETA2(KTS)*PDOPV)
549	! TUUP=VVUPU*(U(I,J,KTS)-U(I,J,KTS+1)+UN(KTS)-UN(KTS+1))
550	! TVUP=VVUPV*(V(I,J,KTS)-V(I,J,KTS+1)+VN(KTS)-VN(KTS+1))
551	! ADUP=TUUP+UN(KTS)-U(I,J,KTS)
552	! ADVP=TVUP+VN(KTS)-V(I,J,KTS)
553	! WRITE(0,*)' NTSD=',NTSD,' I=',ITEST,' J=',JTEST,' K=',KTS &
554	! &, ' ADUP=',ADUP,' ADVP=',ADVP
555	! WRITE(0,*)' U=',U(I,J,KTS),' UN=',UN(KTS) &
556	! &, ' VAD_TEND_U=',VAD_TEND_U(I,KTS) &
557	! &, ' V=',V(I,J,KTS),' VN=',VN(KTS) &
558	! &, ' VAD_TEND_V=',VAD_TEND_V(I,KTS)
559	! WRITE(0,*)' '
560	! ENDIF
561	! endif
562	!
563	!-----------------------------------------------------------------------
564	!*** End of check.
565	!-----------------------------------------------------------------------
566	!
567	ENDDO iloop_for_uv
568	!
569	!-----------------------------------------------------------------------
570	!
571	ENDDO main_vertical
572	!
573	!-----------------------------------------------------------------------
574	!-----------------------------------------------------------------------
575	!
576	!*** COMPUTE HORIZONTAL ADVECTION TENDENCIES.
577	!
578	!-----------------------------------------------------------------------
579	!-----------------------------------------------------------------------
580	!$omp parallel do &
581	!$omp& private(adpdx,adpdy,adt,adu,adv,array0,array1,array2,array3 &
582	!$omp& ,array3_x,dpde,f0,f1,f2,f3,fewp,fnep,fnsp,fpp,fsep,hm &
583	!$omp& ,i,ifp,ifq,ii,ipq,isp,ispa,isq,isqa,iup_adh_j,j,k &
584	!$omp& ,knti_adh,n_iupadh_j,n_iupadv_j,n_iuph_j,pp,qp &
585	!$omp& ,rdpd,rdpdx,rdpdy,tew,tne,tns,tse,tst,tta,ttb &
586	!$omp& ,uew,udy,une,uned,uns,use,used,ust &
587	!$omp& ,vdx,vew,vm,vne,vns,vse,vst)
588	!-----------------------------------------------------------------------
589	!
590	main_horizontal: DO K=KTS,KTE
591	!
592	!-----------------------------------------------------------------------
593	!
594	DO J=MYJS_P4,MYJE_P4
595	DO I=MYIS_P4,MYIE_P4
596	DPDE(I,J)=DETA1_PDTOP(K)+DETA2(K)*PDSLO(I,J)
597	RDPD(I,J)=1./DPDE(I,J)
598	TST(I,J)=T(I,J,K)FFC+TOLD(I,J,K)FBC
599	UST(I,J)=U(I,J,K)FFC+UOLD(I,J,K)FBC
600	VST(I,J)=V(I,J,K)FFC+VOLD(I,J,K)FBC
601	ENDDO
602	ENDDO
603	!
604	!-----------------------------------------------------------------------
605	!*** MASS FLUXES AND MASS POINT ADVECTION COMPONENTS
606	!*** THE NS AND EW FLUXES IN THE FOLLOWING LOOP ARE ON V POINTS
607	!*** FOR T.
608	!-----------------------------------------------------------------------
609	!
610	DO J=MYJS1_P3,MYJE1_P3
611	DO I=MYIS_P3,MYIE_P3
612	!
613	ADPDX=DPDE(I+IVW(J),J)+DPDE(I+IVE(J),J)
614	ADPDY=DPDE(I,J-1)+DPDE(I,J+1)
615	RDPDX(I,J)=1./ADPDX
616	RDPDY(I,J)=1./ADPDY
617	!
618	UDY=U(I,J,K)*DY
619	VDX=V(I,J,K)*DX(I,J)
620	!
621	FEWP=UDY*ADPDX
622	FNSP=VDX*ADPDY
623	!
624	FEW(I,J,K)=FEWP
625	FNS(I,J,K)=FNSP
626	!
627	TEW(I,J)=FEWP*(TST(I+IVE(J),J)-TST(I+IVW(J),J))
628	TNS(I,J)=FNSP*(TST(I,J+1)-TST(I,J-1))
629	!
630	UNED(I,J)=UDY+VDX
631	USED(I,J)=UDY-VDX
632	!
633	ENDDO
634	ENDDO
635	!
636	!-----------------------------------------------------------------------
637	!*** DIAGONAL FLUXES AND DIAGONALLY AVERAGED WIND
638	!*** THE NE AND SE FLUXES ARE ASSOCIATED WITH H POINTS
639	!*** (ACTUALLY JUST TO THE NE AND SE OF EACH H POINT).
640	!-----------------------------------------------------------------------
641	!
642	DO J=MYJS1_P2,MYJE2_P2
643	DO I=MYIS_P2,MYIE_P2
644	FNEP=(UNED(I+IHE(J),J)+UNED(I ,J+1)) &
645	& *(DPDE(I ,J)+DPDE(I+IHE(J),J+1))
646	FNE(I,J,K)=FNEP
647	TNE(I,J)=FNEP*(TST(I+IHE(J),J+1)-TST(I,J))
648	ENDDO
649	ENDDO
650	!
651	DO J=MYJS2_P2,MYJE1_P2
652	DO I=MYIS_P2,MYIE_P2
653	FSEP=(USED(I+IHE(J),J)+USED(I ,J-1)) &
654	& *(DPDE(I ,J)+DPDE(I+IHE(J),J-1))
655	FSE(I,J,K)=FSEP
656	TSE(I,J)=FSEP*(TST(I+IHE(J),J-1)-TST(I,J))
657	!
658	ENDDO
659	ENDDO
660	!
661	!-----------------------------------------------------------------------
662	!*** HORIZONTAL T ADVECTION TENDENCY ADT IS ON H POINTS OF COURSE.
663	!-----------------------------------------------------------------------
664	!
665	DO J=MYJS5,MYJE5
666	DO I=MYIS2,MYIE2
667	ADT(I,J)=(TEW(I+IHW(J),J)+TEW(I+IHE(J),J) &
668	& +TNS(I,J-1)+TNS(I,J+1) &
669	& +TNE(I+IHW(J),J-1)+TNE(I,J) &
670	& +TSE(I,J)+TSE(I+IHW(J),J+1)) &
671	& RDPD(I,J)FAD(I,J)
672	ENDDO
673	ENDDO
674	!
675	!
676	!-----------------------------------------------------------------------
677	!*** CALCULATION OF MOMENTUM ADVECTION COMPONENTS.
678	!-----------------------------------------------------------------------
679	!
680	DO J=MYJS4_P1,MYJE4_P1
681	DO I=MYIS_P1,MYIE_P1
682	!
683	!-----------------------------------------------------------------------
684	!*** THE NS AND EW FLUXES ARE ON H POINTS FOR U AND V.
685	!-----------------------------------------------------------------------
686	!
687	UEW(I,J)=(FEW(I+IHW(J),J,K)+FEW(I+IHE(J),J,K)) &
688	& *(UST(I+IHE(J),J)-UST(I+IHW(J),J))
689	UNS(I,J)=(FNS(I+IHW(J),J,K)+FNS(I+IHE(J),J,K)) &
690	& *(UST(I,J+1)-UST(I,J-1))
691	VEW(I,J)=(FEW(I,J-1,K)+FEW(I,J+1,K)) &
692	& *(VST(I+IHE(J),J)-VST(I+IHW(J),J))
693	VNS(I,J)=(FNS(I,J-1,K)+FNS(I,J+1,K)) &
694	& *(VST(I,J+1)-VST(I,J-1))
695	!
696	!-----------------------------------------------------------------------
697	!*** THE FOLLOWING NE AND SE FLUXES ARE TIED TO V POINTS AND ARE
698	!*** LOCATED JUST TO THE NE AND SE OF THE GIVEN I,J.
699	!-----------------------------------------------------------------------
700	!
701	UNE(I,J)=(FNE(I+IVW(J),J,K)+FNE(I+IVE(J),J,K)) &
702	& *(UST(I+IVE(J),J+1)-UST(I,J))
703	USE(I,J)=(FSE(I+IVW(J),J,K)+FSE(I+IVE(J),J,K)) &
704	& *(UST(I+IVE(J),J-1)-UST(I,J))
705	VNE(I,J)=(FNE(I,J-1,K)+FNE(I,J+1,K)) &
706	& *(VST(I+IVE(J),J+1)-VST(I,J))
707	VSE(I,J)=(FSE(I,J-1,K)+FSE(I,J+1,K)) &
708	& *(VST(I+IVE(J),J-1)-VST(I,J))
709	!
710	!-----------------------------------------------------------------------
711	!
712	ENDDO
713	ENDDO
714	!
715	!-----------------------------------------------------------------------
716	!*** COMPUTE THE ADVECTION TENDENCIES FOR U AND V.
717	!*** THE AD ARRAYS ARE ON THE VELOCITY POINTS.
718	!-----------------------------------------------------------------------
719	!
720	DO J=MYJS5,MYJE5
721	DO I=MYIS2,MYIE2
722	ADU(I,J)=(UEW(I+IVW(J),J)+UEW(I+IVE(J),J) &
723	& +UNS(I,J-1)+UNS(I,J+1) &
724	& +UNE(I+IVW(J),J-1)+UNE(I,J) &
725	& +USE(I,J)+USE(I+IVW(J),J+1)) &
726	& RDPDX(I,J)FAD(I+IVW(J),J)
727	!
728	ADV(I,J)=(VEW(I+IVW(J),J)+VEW(I+IVE(J),J) &
729	& +VNS(I,J-1)+VNS(I,J+1) &
730	& +VNE(I+IVW(J),J-1)+VNE(I,J) &
731	& +VSE(I,J)+VSE(I+IVW(J),J+1)) &
732	& RDPDY(I,J)FAD(I+IVW(J),J)
733	ENDDO
734	ENDDO
735	!
736	!-----------------------------------------------------------------------
737	!
738	!*** END OF JANJIC HORIZONTAL ADVECTION
739	!
740	!-----------------------------------------------------------------------
741	!
742	!*** UPSTREAM ADVECTION OF T
743	!
744	!-----------------------------------------------------------------------
745	!
746	upstream: IF(UPSTRM)THEN
747	!
748	!-----------------------------------------------------------------------
749	!***
750	!*** COMPUTE UPSTREAM COMPUTATIONS ON THIS TASK'S ROWS.
751	!***
752	!-----------------------------------------------------------------------
753	!
754	jloop_upstream: DO J=MYJS2,MYJE2
755	!
756	N_IUPH_J=N_IUP_H(J) ! See explanation in START_DOMAIN_NMM
757	DO II=0,N_IUPH_J-1
758	!
759	I=IUP_H(IMS+II,J)
760	TTA=EMT_LOC(J)*(UST(I,J-1)+UST(I+IHW(J),J) &
761	& +UST(I+IHE(J),J)+UST(I,J+1))
762	TTB=ENT *(VST(I,J-1)+VST(I+IHW(J),J) &
763	& +VST(I+IHE(J),J)+VST(I,J+1))
764	PP=-TTA-TTB
765	QP= TTA-TTB
766	!
767	IF(PP<0.)THEN
768	ISPA(I,J)=-1
769	ELSE
770	ISPA(I,J)= 1
771	ENDIF
772	!
773	IF(QP<0.)THEN
774	ISQA(I,J)=-1
775	ELSE
776	ISQA(I,J)= 1
777	ENDIF
778	!
779	PP=ABS(PP)
780	QP=ABS(QP)
781	ARRAY3_X=PP*QP
782	ARRAY0(I,J)=ARRAY3_X-PP-QP
783	ARRAY1(I,J)=PP-ARRAY3_X
784	ARRAY2(I,J)=QP-ARRAY3_X
785	ARRAY3(I,J)=ARRAY3_X
786	ENDDO
787	!
788	!-----------------------------------------------------------------------
789	!
790	N_IUPADH_J=N_IUP_ADH(J)
791	KNTI_ADH=1
792	IUP_ADH_J=IUP_ADH(IMS,J)
793	!
794	iloop_T: DO II=0,N_IUPH_J-1
795	!
796	I=IUP_H(IMS+II,J)
797	!
798	ISP=ISPA(I,J)
799	ISQ=ISQA(I,J)
800	IFP=(ISP-1)/2
801	IFQ=(-ISQ-1)/2
802	IPQ=(ISP-ISQ)/2
803	!
804	!-----------------------------------------------------------------------
805	!
806	IF(I==IUP_ADH_J)THEN ! Upstream advection T tendencies
807	!
808	ISP=ISPA(I,J)
809	ISQ=ISQA(I,J)
810	IFP=(ISP-1)/2
811	IFQ=(-ISQ-1)/2
812	IPQ=(ISP-ISQ)/2
813	!
814	F0=ARRAY0(I,J)
815	F1=ARRAY1(I,J)
816	F2=ARRAY2(I,J)
817	F3=ARRAY3(I,J)
818	!
819	ADT(I,J)=F0*T(I,J,K) &
820	& +F1*T(I+IHE(J)+IFP,J+ISP,K) &
821	& +F2*T(I+IHE(J)+IFQ,J+ISQ,K) &
822	+F3*T(I+IPQ,J+ISP+ISQ,K)
823	!
824	!-----------------------------------------------------------------------
825	!
826	IF(KNTI_ADH<N_IUPADH_J)THEN
827	IUP_ADH_J=IUP_ADH(IMS+KNTI_ADH,J)
828	KNTI_ADH=KNTI_ADH+1
829	ENDIF
830	!
831	ENDIF ! End of upstream advection T tendency IF block
832	!
833	ENDDO iloop_T
834	!
835	!-----------------------------------------------------------------------
836	!
837	!*** UPSTREAM ADVECTION OF VELOCITY COMPONENTS
838	!
839	!-----------------------------------------------------------------------
840	!
841	N_IUPADV_J=N_IUP_ADV(J)
842	!
843	DO II=0,N_IUPADV_J-1
844	I=IUP_ADV(IMS+II,J)
845	!
846	TTA=EM_LOC(J)*UST(I,J)
847	TTB=EN *VST(I,J)
848	PP=-TTA-TTB
849	QP=TTA-TTB
850	!
851	IF(PP<0.)THEN
852	ISP=-1
853	ELSE
854	ISP= 1
855	ENDIF
856	!
857	IF(QP<0.)THEN
858	ISQ=-1
859	ELSE
860	ISQ= 1
861	ENDIF
862	!
863	IFP=(ISP-1)/2
864	IFQ=(-ISQ-1)/2
865	IPQ=(ISP-ISQ)/2
866	PP=ABS(PP)
867	QP=ABS(QP)
868	F3=PP*QP
869	F0=F3-PP-QP
870	F1=PP-F3
871	F2=QP-F3
872	!
873	ADU(I,J)=F0*U(I,J,K) &
874	& +F1*U(I+IVE(J)+IFP,J+ISP,K) &
875	& +F2*U(I+IVE(J)+IFQ,J+ISQ,K) &
876	& +F3*U(I+IPQ,J+ISP+ISQ,K)
877	!
878	ADV(I,J)=F0*V(I,J,K) &
879	& +F1*V(I+IVE(J)+IFP,J+ISP,K) &
880	& +F2*V(I+IVE(J)+IFQ,J+ISQ,K) &
881	& +F3*V(I+IPQ,J+ISP+ISQ,K)
882	!
883	ENDDO
884	!
885	ENDDO jloop_upstream
886	!
887	!-----------------------------------------------------------------------
888	!
889	ENDIF upstream
890	!
891	!-----------------------------------------------------------------------
892	!
893	!*** END OF HORIZONTAL ADVECTION
894	!
895	!-----------------------------------------------------------------------
896	!
897	!*** NOW SUM THE VERTICAL AND HORIZONTAL TENDENCIES,
898	!*** CURVATURE AND CORIOLIS TERMS.
899	!
900	!-----------------------------------------------------------------------
901	!
902	DO J=MYJS2,MYJE2
903	DO I=MYIS1,MYIE1
904	HM=HBM2(I,J)
905	VM=VBM2(I,J)
906	ADT(I,J)=(VAD_TEND_T(I,J,K)+2.ADT(I,J))HM
907	!
908	FPP=CURV(I,J)2.UST(I,J)+F(I,J)*2.
909	ADU(I,J)=(VAD_TEND_U(I,J,K)+2.ADU(I,J)+VST(I,J)FPP)*VM
910	ADV(I,J)=(VAD_TEND_V(I,J,K)+2.ADV(I,J)-UST(I,J)FPP)*VM
911	ENDDO
912	ENDDO
913	!
914	!-----------------------------------------------------------------------
915	!*** SAVE THE OLD VALUES FOR TIMESTEPPING
916	!-----------------------------------------------------------------------
917	!
918	DO J=MYJS_P4,MYJE_P4
919	DO I=MYIS_P4,MYIE_P4
920	TOLD(I,J,K)=T(I,J,K)
921	UOLD(I,J,K)=U(I,J,K)
922	VOLD(I,J,K)=V(I,J,K)
923	ENDDO
924	ENDDO
925	!
926	!-----------------------------------------------------------------------
927	!*** FINALLY UPDATE THE PROGNOSTIC VARIABLES
928	!-----------------------------------------------------------------------
929	!
930	DO J=MYJS2,MYJE2
931	DO I=MYIS1,MYIE1
932	T(I,J,K)=ADT(I,J)+T(I,J,K)
933	U(I,J,K)=ADU(I,J)+U(I,J,K)
934	V(I,J,K)=ADV(I,J)+V(I,J,K)
935	ENDDO
936	ENDDO
937	!
938	!-----------------------------------------------------------------------
939	!
940	ENDDO main_horizontal
941	!
942	!-----------------------------------------------------------------------
943	!
944	END SUBROUTINE ADVE
945	!
946	!-----------------------------------------------------------------------
947	!
948	!***********************************************************************
949	SUBROUTINE VAD2(NTSD,DT,IDTAD,DX,DY &
950	& ,AETA1,AETA2,DETA1,DETA2,PDSL,PDTOP,HBM2 &
951	& ,Q,Q2,CWM,PETDT &
952	& ,N_IUP_H,N_IUP_V &
953	& ,N_IUP_ADH,N_IUP_ADV &
954	& ,IUP_H,IUP_V,IUP_ADH,IUP_ADV &
955	& ,IHE,IHW,IVE,IVW &
956	& ,IDS,IDE,JDS,JDE,KDS,KDE &
957	& ,IMS,IME,JMS,JME,KMS,KME &
958	& ,ITS,ITE,JTS,JTE,KTS,KTE)
959	!***********************************************************************
960	!$$$ SUBPROGRAM DOCUMENTATION BLOCK
961	! . . .
962	! SUBPROGRAM: VAD2 VERTICAL ADVECTION OF H2O SUBSTANCE AND TKE
963	! PRGRMMR: JANJIC ORG: W/NP22 DATE: 96-07-19
964	!
965	! ABSTRACT:
966	! VAD2 CALCULATES THE CONTRIBUTION OF THE VERTICAL ADVECTION
967	! TO THE TENDENCIES OF WATER SUBSTANCE AND TKE AND THEN UPDATES
968	! THOSE VARIABLES. AN ANTI-FILTERING TECHNIQUE IS USED.
969	!
970	! PROGRAM HISTORY LOG:
971	! 96-07-19 JANJIC - ORIGINATOR
972	! 98-11-02 BLACK - MODIFIED FOR DISTRIBUTED MEMORY
973	! 99-03-17 TUCCILLO - INCORPORATED MPI_ALLREDUCE FOR GLOBAL SUM
974	! 02-02-06 BLACK - CONVERTED TO WRF FORMAT
975	! 02-09-06 WOLFE - MORE CONVERSION TO GLOBAL INDEXING
976	! 04-11-23 BLACK - THREADED
977	! 05-12-14 BLACK - CONVERTED FROM IKJ TO IJK
978	! 07-08-14 janjic - bc & no conservation in the advection step
979	!
980	! USAGE: CALL VAD2 FROM SUBROUTINE SOLVE_NMM
981	! INPUT ARGUMENT LIST:
982	!
983	! OUTPUT ARGUMENT LIST
984	!
985	! OUTPUT FILES:
986	! NONE
987	! SUBPROGRAMS CALLED:
988	!
989	! UNIQUE: NONE
990	!
991	! LIBRARY: NONE
992	!
993	! ATTRIBUTES:
994	! LANGUAGE: FORTRAN 90
995	! MACHINE : IBM SP
996	!$$$
997	!***********************************************************************
998	!----------------------------------------------------------------------
999	!
1000	IMPLICIT NONE
1001	!
1002	!----------------------------------------------------------------------
1003	!
1004	INTEGER,INTENT(IN) :: IDS,IDE,JDS,JDE,KDS,KDE &
1005	& ,IMS,IME,JMS,JME,KMS,KME &
1006	,ITS,ITE,JTS,JTE,KTS,KTE
1007	!
1008	INTEGER,DIMENSION(JMS:JME),INTENT(IN) :: IHE,IHW,IVE,IVW
1009	INTEGER,DIMENSION(JMS:JME),INTENT(IN) :: N_IUP_H,N_IUP_V &
1010	& ,N_IUP_ADH,N_IUP_ADV
1011	INTEGER,DIMENSION(IMS:IME,JMS:JME),INTENT(IN) :: IUP_H,IUP_V &
1012	& ,IUP_ADH,IUP_ADV
1013	!
1014	INTEGER,INTENT(IN) :: IDTAD,NTSD
1015	!
1016	REAL,INTENT(IN) :: DT,DY,PDTOP
1017	!
1018	REAL,DIMENSION(KMS:KME),INTENT(IN) :: AETA1,AETA2,DETA1,DETA2
1019	!
1020	REAL,DIMENSION(IMS:IME,JMS:JME),INTENT(IN) :: DX,HBM2,PDSL
1021	!
1022	REAL,DIMENSION(IMS:IME,JMS:JME,KMS:KME),INTENT(IN) :: PETDT
1023	!
1024	REAL,DIMENSION(IMS:IME,JMS:JME,KMS:KME),INTENT(INOUT) :: CWM,Q,Q2
1025	!
1026	!*** LOCAL VARIABLES
1027	!----------------------------------------------------------------------
1028	!
1029	REAL,PARAMETER :: FF1=0.500
1030	!
1031	LOGICAL,SAVE :: TRADITIONAL=.TRUE.
1032	!
1033	INTEGER :: I,IRECV,J,JFP,JFQ,K,LAP,LLAP
1034	!
1035	INTEGER,DIMENSION(KTS:KTE) :: LA
1036	!
1037	REAL*8 :: ADDT,AFRP,D2PQE,D2PQQ,D2PQW,DEP,DETAP,DQP &
1038	& ,DWP,E00,E4P,EP,EP0,HADDT,HBM2IJ &
1039	& ,Q00,Q4P,QP,QP0 &
1040	& ,rdpdn,rdpup,sfacek,sfacqk,sfacwk,RFC,RR &
1041	& ,SUMNE,SUMNQ,SUMNW,SUMPE,SUMPQ,SUMPW &
1042	& ,W00,W4P,WP,WP0
1043	!
1044	REAL,DIMENSION(KTS:KTE) :: AFR,DEL,DQL,DWL,E3,E4,PETDTK &
1045	& ,RFACE,RFACQ,RFACW,Q3,Q4,W3,W4
1046	!
1047	!-----------------------------------------------------------------------
1048	!***********************************************************************
1049	!-----------------------------------------------------------------------
1050	!
1051	ADDT=REAL(IDTAD)*DT
1052	!
1053	!-----------------------------------------------------------------------
1054	!$omp parallel do &
1055	!$omp& private(afr,afrp,bot,d2pqe,d2pqq,d2pqw,del,dep,detap,dpdn,dpup &
1056	!$omp& ,dql,dqp,dwl,dwp,e00,e3,e4,e4p,ep,ep0,haddt,i,j,k &
1057	!$omp& ,la,lap,llap,petdtk,q00,q3,q4,q4p,qp,qp0,rfacek,rfacqk &
1058	!$omp& ,rfacwk,rfc,rr,sumne,sumnq,sumnw,sumpe,sumpq,sumpw,top &
1059	!$omp& ,w00,w3,w4,w4p,wp,wp0)
1060	!-----------------------------------------------------------------------
1061	!
1062	main_integration: DO J=MYJS2,MYJE2
1063	!
1064	!-----------------------------------------------------------------------
1065	!
1066	main_iloop: DO I=MYIS1_P1,MYIE1_P1
1067	!
1068	!-----------------------------------------------------------------------
1069	!
1070	E3(KTE)=Q2(I,J,KTE)*0.5
1071	!
1072	DO K=KTE-1,KTS,-1
1073	E3(K)=MAX((Q2(I,J,K+1)+Q2(I,J,K))*0.5,EPSQ2)
1074	ENDDO
1075	!
1076	DO K=KTS,KTE
1077	Q3(K)=MAX(Q(I,J,K),EPSQ)
1078	W3(K)=MAX(CWM(I,J,K),CLIMIT)
1079	E4(K)=E3(K)
1080	Q4(K)=Q3(K)
1081	W4(K)=W3(K)
1082	ENDDO
1083	!
1084	IF(TRADITIONAL)THEN
1085	PETDTK(KTE)=PETDT(I,J,KTE-1)*0.5
1086	!
1087	DO K=KTE-1,KTS+1,-1
1088	PETDTK(K)=(PETDT(I,J,K)+PETDT(I,J,K-1))*0.5
1089	ENDDO
1090	!
1091	PETDTK(KTS)=PETDT(I,J,KTS)*0.5
1092	!
1093	ELSE
1094	!
1095	!-----------------------------------------------------------------------
1096	!*** PERFORM HORIZONTAL AVERAGING OF VERTICAL VELOCITY
1097	!-----------------------------------------------------------------------
1098	!
1099	PETDTK(KTE)=(PETDT(I+IHW(J-1),J-1,KTE-1) &
1100	& +PETDT(I+IHE(J-1),J-1,KTE-1) &
1101	& +PETDT(I+IHW(J+1),J+1,KTE-1) &
1102	& +PETDT(I+IHE(J+1),J+1,KTE-1) &
1103	& +PETDT(I,J,KTE-1)4. )0.0625
1104	!
1105	DO K=KTE-1,KTS+1,-1
1106	PETDTK(K)=(PETDT(I+IHW(J-1),J-1,K-1) &
1107	+PETDT(I+IHE(J-1),J-1,K-1) &
1108	& +PETDT(I+IHW(J+1),J+1,K-1) &
1109	& +PETDT(I+IHE(J+1),J+1,K-1) &
1110	& +PETDT(I+IHW(J-1),J-1,K ) &
1111	& +PETDT(I+IHE(J-1),J-1,K ) &
1112	& +PETDT(I+IHW(J+1),J+1,K ) &
1113	& +PETDT(I+IHE(J+1),J+1,K ) &
1114	& +(PETDT(I,J,K-1)+PETDT(I,J,K))*4. &
1115	& )*0.0625
1116	ENDDO
1117	!
1118	PETDTK(KTS)=(PETDT(I+IHW(J-1),J-1,KTS) &
1119	& +PETDT(I+IHE(J-1),J-1,KTS) &
1120	& +PETDT(I+IHW(J+1),J+1,KTS) &
1121	& +PETDT(I+IHE(J+1),J+1,KTS) &
1122	& +PETDT(I,J,KTS)4. )0.0625
1123
1124	ENDIF
1125	!
1126	!-----------------------------------------------------------------------
1127	!
1128	HADDT=-ADDT*HBM2(I,J)
1129	!
1130	DO K=KTE,KTS,-1
1131	RR=PETDTK(K)*HADDT
1132	!
1133	IF(RR<0.)THEN
1134	LAP=1
1135	ELSE
1136	LAP=-1
1137	ENDIF
1138	!
1139	LA(K)=LAP
1140	LLAP=K+LAP
1141	!
1142	if(llap.gt.kts-1.and.llap.lt.kte+1) then ! internal and outflow pts.
1143	rr=abs(rr &
1144	& /((aeta1(llap)-aeta1(k))*pdtop &
1145	& +(aeta2(llap)-aeta2(k))*pdsl(i,j)))
1146	if(rr.gt.0.999) rr=0.999
1147	!
1148	AFR(K)=(((FF4RR+FF3)RR+FF2)RR+FF1)RR
1149	dql(k)=(q3(llap)-q3(k))*rr
1150	dwl(k)=(w3(llap)-w3(k))*rr
1151	del(k)=(e3(llap)-e3(k))*rr
1152	elseif(llap.eq.kts-1) then
1153	!
1154	!chem rr=abs(rr &
1155	!chem /((1.-aeta2(kts))*pdsl(i,j)))
1156	!chem afr(kts)=0.
1157	!chem dql(kts)=(epsq -q3(kts))*rr
1158	!chem dwl(kts)=(climit-w3(kts))*rr
1159	!chem del(kts)=(epsq2 -e3(kts))*rr
1160	!
1161	rr=0.
1162	afr(kts)=0.
1163	dql(kts)=0.
1164	dwl(kts)=0.
1165	del(kts)=0.
1166	else
1167	rr=abs(rr &
1168	/(aeta1(kte)*pdtop))
1169	afr(kte)=0.
1170	dql(kte)=(epsq -q3(kte))*rr
1171	dwl(kte)=(climit-w3(kte))*rr
1172	del(kte)=(epsq2 -e3(kte))*rr
1173	endif
1174	ENDDO
1175	!
1176	!-----------------------------------------------------------------------
1177	!
1178	DO K=KTS,KTE
1179	Q4(K)=Q3(K)+DQL(K)
1180	W4(K)=W3(K)+DWL(K)
1181	E4(K)=E3(K)+DEL(K)
1182	ENDDO
1183	!
1184	!-----------------------------------------------------------------------
1185	!*** ANTI-FILTERING STEP
1186	!-----------------------------------------------------------------------
1187	!
1188	SUMPQ=0.
1189	SUMNQ=0.
1190	SUMPW=0.
1191	SUMNW=0.
1192	SUMPE=0.
1193	SUMNE=0.
1194	!
1195	!*** ANTI-FILTERING LIMITERS
1196	!
1197	antifilter: DO K=KTE-1,KTS+1,-1
1198	!
1199	DETAP=DETA1(K)PDTOP+DETA2(K)PDSL(I,J)
1200	!
1201	DQL(K)=0.
1202	DWL(K)=0.
1203	DEL(K)=0.
1204	!
1205	Q4P=Q4(K)
1206	W4P=W4(K)
1207	E4P=E4(K)
1208	!
1209	LAP=LA(K)
1210	!
1211	if(lap.ne.0)then
1212	rdpdn=1./((aeta1(k+lap)-aeta1(k))*pdtop &
1213	& +(aeta2(k+lap)-aeta2(k))*pdsl(i,j))
1214	rdpup=1./((aeta1(k)-aeta1(k-lap))*pdtop &
1215	& +(aeta2(k)-aeta2(k-lap))*pdsl(i,j))
1216	!
1217	afrp=afr(k)*detap
1218	!
1219	d2pqq=((q4(k+lap)-q4p)*rdpdn &
1220	& -(q4p-q4(k-lap))rdpup)afrp
1221	d2pqw=((w4(k+lap)-w4p)*rdpdn &
1222	& -(w4p-w4(k-lap))rdpup)afrp
1223	d2pqe=((e4(k+lap)-e4p)*rdpdn &
1224	& -(e4p-e4(k-lap))rdpup)afrp
1225	ELSE
1226	D2PQQ=0.
1227	D2PQW=0.
1228	D2PQE=0.
1229	ENDIF
1230	!
1231	QP=Q4P-D2PQQ
1232	WP=W4P-D2PQW
1233	EP=E4P-D2PQE
1234	!
1235	Q00=Q3(K)
1236	QP0=Q3(K+LAP)
1237	!
1238	W00=W3(K)
1239	WP0=W3(K+LAP)
1240	!
1241	E00=E3(K)
1242	EP0=E3(K+LAP)
1243	!
1244	IF(LAP/=0)THEN
1245	QP=MAX(QP,MIN(Q00,QP0))
1246	QP=MIN(QP,MAX(Q00,QP0))
1247	WP=MAX(WP,MIN(W00,WP0))
1248	WP=MIN(WP,MAX(W00,WP0))
1249	EP=MAX(EP,MIN(E00,EP0))
1250	EP=MIN(EP,MAX(E00,EP0))
1251	ENDIF
1252	!
1253	dqp=qp-q4p
1254	dwp=wp-w4p
1255	dep=ep-e4p
1256	!
1257	DQL(K)=DQP
1258	DWL(K)=DWP
1259	DEL(K)=DEP
1260	!
1261	DQP=DQP*DETAP
1262	DWP=DWP*DETAP
1263	DEP=DEP*DETAP
1264	!
1265	IF(DQP>0.)THEN
1266	SUMPQ=SUMPQ+DQP
1267	ELSE
1268	SUMNQ=SUMNQ+DQP
1269	ENDIF
1270	!
1271	IF(DWP>0.)THEN
1272	SUMPW=SUMPW+DWP
1273	ELSE
1274	SUMNW=SUMNW+DWP
1275	ENDIF
1276	!
1277	IF(DEP>0.)THEN
1278	SUMPE=SUMPE+DEP
1279	ELSE
1280	SUMNE=SUMNE+DEP
1281	ENDIF
1282	!
1283	ENDDO antifilter
1284	!
1285	!-----------------------------------------------------------------------
1286	!
1287	DQL(KTS)=0.
1288	DWL(KTS)=0.
1289	DEL(KTS)=0.
1290	!
1291	DQL(KTE)=0.
1292	DWL(KTE)=0.
1293	DEL(KTE)=0.
1294	!
1295	!-----------------------------------------------------------------------
1296	!*** FIRST MOMENT CONSERVING FACTOR
1297	!-----------------------------------------------------------------------
1298	!
1299	if(sumpq*(-sumnq).gt.1.e-9) then
1300	sfacqk=-sumnq/sumpq
1301	else
1302	sfacqk=0.
1303	endif
1304	!
1305	if(sumpw*(-sumnw).gt.1.e-9) then
1306	sfacwk=-sumnw/sumpw
1307	else
1308	sfacwk=0.
1309	endif
1310	!
1311	if(sumpe*(-sumne).gt.1.e-9) then
1312	sfacek=-sumne/sumpe
1313	else
1314	sfacek=0.
1315	endif
1316	!
1317	!-----------------------------------------------------------------------
1318	!*** IMPOSE CONSERVATION ON ANTI-FILTERING
1319	!-----------------------------------------------------------------------
1320	!
1321	DO K=KTE,KTS,-1
1322	!
1323	dqp=dql(k)
1324	if(sfacqk.gt.0.) then
1325	if(sfacqk.ge.1.) then
1326	if(dqp.lt.0.) dqp=dqp/sfacqk
1327	else
1328	if(dqp.gt.0.) dqp=dqp*sfacqk
1329	endif
1330	else
1331	dqp=0.
1332	endif
1333	q (i,j,k)=q4(k)+dqp
1334	!
1335	dwp=dwl(k)
1336	if(sfacwk.gt.0.) then
1337	if(sfacwk.ge.1.) then
1338	if(dwp.lt.0.) dwp=dwp/sfacwk
1339	else
1340	if(dwp.gt.0.) dwp=dwp*sfacwk
1341	endif
1342	else
1343	dwp=0.
1344	endif
1345	cwm(i,j,k)=w4(k)+dwp
1346	!
1347	dep=del(k)
1348	if(sfacek.gt.0.) then
1349	if(sfacek.ge.1.) then
1350	if(dep.lt.0.) dep=dep/sfacek
1351	else
1352	if(dep.gt.0.) dep=dep*sfacek
1353	endif
1354	else
1355	dep=0.
1356	endif
1357	e3 ( k)=e4(k)+dep
1358	!
1359	ENDDO
1360	!-----------------------------------------------------------------------
1361	HBM2IJ=HBM2(I,J)
1362	Q2(I,J,KTE)=MAX(E3(KTE)+E3(KTE)-EPSQ2,EPSQ2)*HBM2IJ &
1363	& +Q2(I,J,KTE)*(1.-HBM2IJ)
1364	DO K=KTE-1,KTS+1,-1
1365	Q2(I,J,K)=MAX(E3(K)+E3(K)-Q2(I,J,K+1),EPSQ2)*HBM2IJ &
1366	& +Q2(I,J,K)*(1.-HBM2IJ)
1367	ENDDO
1368	!-----------------------------------------------------------------------
1369	!
1370	ENDDO main_iloop
1371	!
1372	!-----------------------------------------------------------------------
1373	!
1374	ENDDO main_integration
1375	!
1376	!-----------------------------------------------------------------------
1377	!
1378	END SUBROUTINE VAD2
1379	!
1380
1381	!-----------------------------------------------------------------------
1382	!
1383	!-----------------------------------------------------------------------
1384	!***********************************************************************
1385	SUBROUTINE HAD2( &
1386	#if defined(DM_PARALLEL)
1387	& domdesc , &
1388	#endif
1389	& NTSD,DT,IDTAD,DX,DY &
1390	& ,AETA1,AETA2,DETA1,DETA2,PDSL,PDTOP &
1391	& ,HBM2,HBM3 &
1392	& ,Q,Q2,CWM,U,V,Z,HYDRO &
1393	& ,N_IUP_H,N_IUP_V &
1394	& ,N_IUP_ADH,N_IUP_ADV &
1395	& ,IUP_H,IUP_V,IUP_ADH,IUP_ADV &
1396	& ,IHE,IHW,IVE,IVW &
1397	& ,IDS,IDE,JDS,JDE,KDS,KDE &
1398	& ,IMS,IME,JMS,JME,KMS,KME &
1399	& ,ITS,ITE,JTS,JTE,KTS,KTE)
1400	!***********************************************************************
1401	!$$$ SUBPROGRAM DOCUMENTATION BLOCK
1402	! . . .
1403	! SUBPROGRAM: HAD2 HORIZONTAL ADVECTION OF H2O AND TKE
1404	! PRGRMMR: JANJIC ORG: W/NP22 DATE: 96-07-19
1405	!
1406	! ABSTRACT:
1407	! HAD2 CALCULATES THE CONTRIBUTION OF THE HORIZONTAL ADVECTION
1408	! TO THE TENDENCIES OF WATER SUBSTANCE AND TKE AND THEN
1409	! UPDATES THOSE VARIABLES. AN ANTI-FILTERING TECHNIQUE IS USED.
1410	!
1411	! PROGRAM HISTORY LOG:
1412	! 96-07-19 JANJIC - ORIGINATOR
1413	! 98-11-02 BLACK - MODIFIED FOR DISTRIBUTED MEMORY
1414	! 99-03-17 TUCCILLO - INCORPORATED MPI_ALLREDUCE FOR GLOBAL SUM
1415	! 02-02-06 BLACK - CONVERTED TO WRF FORMAT
1416	! 02-09-06 WOLFE - MORE CONVERSION TO GLOBAL INDEXING
1417	! 03-05-23 JANJIC - ADDED SLOPE FACTOR
1418	! 04-11-23 BLACK - THREADED
1419	! 05-12-14 BLACK - CONVERTED FROM IKJ TO IJK
1420	! 07-08-14 janjic - no conservation in advection step
1421	!
1422	! USAGE: CALL HAD2 FROM SUBROUTINE SOLVE_NMM
1423	! INPUT ARGUMENT LIST:
1424	!
1425	! OUTPUT ARGUMENT LIST
1426	!
1427	! OUTPUT FILES:
1428	! NONE
1429	! SUBPROGRAMS CALLED:
1430	!
1431	! UNIQUE: NONE
1432	!
1433	! LIBRARY: NONE
1434	!
1435	! ATTRIBUTES:
1436	! LANGUAGE: FORTRAN 90
1437	! MACHINE : IBM SP
1438	!$$$
1439	!***********************************************************************
1440	!-----------------------------------------------------------------------
1441	!
1442	IMPLICIT NONE
1443	!
1444	!-----------------------------------------------------------------------
1445	!
1446	INTEGER,INTENT(IN) :: IDS,IDE,JDS,JDE,KDS,KDE &
1447	& ,IMS,IME,JMS,JME,KMS,KME &
1448	& ,ITS,ITE,JTS,JTE,KTS,KTE
1449	!
1450	INTEGER,DIMENSION(JMS:JME),INTENT(IN) :: IHE,IHW,IVE,IVW
1451	INTEGER,DIMENSION(JMS:JME),INTENT(IN) :: N_IUP_H,N_IUP_V &
1452	& ,N_IUP_ADH,N_IUP_ADV
1453	INTEGER,DIMENSION(IMS:IME,JMS:JME),INTENT(IN) :: IUP_H,IUP_V &
1454	& ,IUP_ADH,IUP_ADV
1455	!
1456	!-----------------------------------------------------------------------
1457	!
1458	INTEGER,INTENT(IN) :: IDTAD,NTSD
1459	!
1460	REAL,INTENT(IN) :: DT,DY,PDTOP
1461	!
1462	REAL,DIMENSION(KMS:KME),INTENT(IN) :: AETA1,AETA2,DETA1,DETA2
1463	!
1464	REAL,DIMENSION(IMS:IME,JMS:JME),INTENT(IN) :: DX,HBM2,HBM3,PDSL
1465	!
1466	REAL,DIMENSION(IMS:IME,JMS:JME,KMS:KME),INTENT(IN) :: U,V,Z
1467	!
1468	REAL,DIMENSION(IMS:IME,JMS:JME,KMS:KME),INTENT(INOUT) :: CWM,Q,Q2
1469	!
1470	LOGICAL,INTENT(IN) :: HYDRO
1471	!
1472	!-----------------------------------------------------------------------
1473	!*** LOCAL VARIABLES
1474	!-----------------------------------------------------------------------
1475	!
1476	REAL,PARAMETER :: FF1=0.530
1477	!
1478	#ifdef DM_PARALLEL
1479	INTEGER :: DOMDESC
1480	#endif
1481	!
1482	#if defined(BIT_FOR_BIT) && defined(DM_PARALLEL) && !defined(STUBMPI)
1483	LOGICAL,EXTERNAL :: WRF_DM_ON_MONITOR
1484	REAL,DIMENSION(IMS:IME,JMS:JME,KMS:KME,6) :: XSUMS_L
1485	REAL,DIMENSION(IDS:IDE,JDS:JDE,KDS:KDE,6) :: XSUMS_G
1486	#endif
1487	!
1488	LOGICAL :: BOT,TOP
1489	!
1490	INTEGER :: I,IRECV,J,JFP,JFQ,K,LAP,LLAP,MPI_COMM_COMP
1491	INTEGER :: N
1492	!
1493	INTEGER,DIMENSION(ITS-5:ITE+5,JTS-5:JTE+5,KTS:KTE) :: IFPA,IFPF &
1494	& ,IFQA,IFQF &
1495	& ,JFPA,JFPF &
1496	& ,JFQA,JFQF
1497	!
1498	REAL :: ADDT,AFRP,CRIT,D2PQE,D2PQQ,D2PQW,DEP,DESTIJ,DQP,DQSTIJ &
1499	& ,DVOLP,DWP,DWSTIJ,DZA,DZB,E00,E0Q,E1X,E2IJ,E4P,ENH,EP,EP0 &
1500	& ,ESTIJ,FPQ,HAFP,HAFQ,HBM2IJ,HM,PP,PPQ00,Q00,Q0Q &
1501	& ,Q1IJ,Q4P,QP,QP0,QSTIJ,RDY,RFACE,RFACQ,RFACW,RFC &
1502	& ,RFEIJ,RFQIJ,RFWIJ,RR,SLOPAC,SPP,SQP,SSA,SSB,SUMNE,SUMNQ &
1503	& ,SUMNW,SUMPE,SUMPQ,SUMPW,TTA,TTB,W00,W0Q,W1IJ,W4P,WP,WP0 &
1504	& ,WSTIJ
1505	!
1506	DOUBLE PRECISION,DIMENSION(6,KTS:KTE) :: GSUMS,XSUMS
1507	!
1508	REAL,DIMENSION(KTS:KTE) :: AFR,DEL,DQL,DWL,E3,E4 &
1509	& ,Q3,Q4,W3,W4
1510	!
1511	REAL,DIMENSION(ITS-5:ITE+5,JTS-5:JTE+5) :: DARE,EMH
1512	!
1513	REAL,DIMENSION(ITS-5:ITE+5,JTS-5:JTE+5,KTS:KTE) :: AFP,AFQ,DEST &
1514	& ,DQST,DVOL,DWST &
1515	& ,E1,E2,Q1,W1
1516	!-----------------------------------------------------------------------
1517	integer :: nunit,ier
1518	save nunit
1519	!-----------------------------------------------------------------------
1520	!***********************************************************************
1521	!-----------------------------------------------------------------------
1522	!
1523	RDY=1./DY
1524	SLOPAC=SLOPHTSQRT(2.)0.5*50.
1525	CRIT=SLOPACREAL(IDTAD)DTRDY1000.
1526	!
1527	ADDT=REAL(IDTAD)*DT
1528	ENH=ADDT/(08.*DY)
1529	!
1530	!-----------------------------------------------------------------------
1531	!$omp parallel do &
1532	!$omp& private(i,j)
1533	DO J=MYJS_P3,MYJE_P3
1534	DO I=MYIS_P2,MYIE_P2
1535	EMH (I,J)=ADDT/(08.*DX(I,J))
1536	DARE(I,J)=HBM3(I,J)DX(I,J)DY
1537	E1(I,J,KTE)=MAX(Q2(I,J,KTE)*0.5,EPSQ2)
1538	E2(I,J,KTE)=E1(I,J,KTE)
1539	ENDDO
1540	ENDDO
1541	!-----------------------------------------------------------------------
1542	!$omp parallel do &
1543	!$omp& private(dza,dzb,e1x,fpq,hm,i,j,jfp,jfq,k,pp,qp,ssa,ssb,spp,sqp &
1544	!$omp& ,tta,ttb)
1545	!-----------------------------------------------------------------------
1546	!
1547	vertical_1: DO K=KTS,KTE
1548	!
1549	!-----------------------------------------------------------------------
1550	!
1551	DO J=MYJS_P3,MYJE_P3
1552	DO I=MYIS_P2,MYIE_P2
1553	DVOL(I,J,K)=DARE(I,J)(DETA1(K)PDTOP+DETA2(K)*PDSL(I,J))
1554	Q (I,J,K)=MAX(Q (I,J,K),EPSQ)
1555	CWM(I,J,K)=MAX(CWM(I,J,K),CLIMIT)
1556	Q1 (I,J,K)=Q (I,J,K)
1557	W1 (I,J,K)=CWM(I,J,K)
1558	ENDDO
1559	ENDDO
1560	!
1561	IF(K<KTE)THEN
1562	DO J=MYJS_P3,MYJE_P3
1563	DO I=MYIS_P2,MYIE_P2
1564	E1X=(Q2(I,J,K+1)+Q2(I,J,K))*0.5
1565	E1(I,J,K)=MAX(E1X,EPSQ2)
1566	E2(I,J,K)=E1(I,J,K)
1567	ENDDO
1568	ENDDO
1569	ENDIF
1570	!
1571	!-----------------------------------------------------------------------
1572	!
1573	DO J=MYJS2_P1,MYJE2_P1
1574	DO I=MYIS1_P1,MYIE1_P1
1575	!
1576	HM=HBM2(I,J)
1577	TTA=(U(I,J-1,K)+U(I+IHW(J),J,K)+U(I+IHE(J),J,K)+U(I,J+1,K)) &
1578	& EMH(I,J)HM
1579	TTB=(V(I,J-1,K)+V(I+IHW(J),J,K)+V(I+IHE(J),J,K)+V(I,J+1,K)) &
1580	& ENHHBM2(I,J)
1581	!
1582	SPP=-TTA-TTB
1583	SQP= TTA-TTB
1584	!
1585	IF(SPP<0.)THEN
1586	JFP=-1
1587	ELSE
1588	JFP=1
1589	ENDIF
1590	IF(SQP<0.)THEN
1591	JFQ=-1
1592	ELSE
1593	JFQ=1
1594	ENDIF
1595	!
1596	IFPA(I,J,K)=IHE(J)+I+( JFP-1)/2
1597	IFQA(I,J,K)=IHE(J)+I+(-JFQ-1)/2
1598	!
1599	JFPA(I,J,K)=J+JFP
1600	JFQA(I,J,K)=J+JFQ
1601	!
1602	IFPF(I,J,K)=IHE(J)+I+(-JFP-1)/2
1603	IFQF(I,J,K)=IHE(J)+I+( JFQ-1)/2
1604	!
1605	JFPF(I,J,K)=J-JFP
1606	JFQF(I,J,K)=J-JFQ
1607	! if(i==111.and.j==438.and.k==1)then
1608	! endif
1609	!
1610	!-----------------------------------------------------------------------
1611	IF(.NOT.HYDRO)THEN ! z currently not available for hydro=.true.
1612	DZA=(Z(IFPA(I,J,K),JFPA(I,J,K),K)-Z(I,J,K))*RDY
1613	DZB=(Z(IFQA(I,J,K),JFQA(I,J,K),K)-Z(I,J,K))*RDY
1614	!
1615	IF(ABS(DZA)>SLOPAC)THEN
1616	SSA=DZA*SPP
1617	IF(SSA>CRIT)THEN
1618	SPP=0. !spp*.1
1619	ENDIF
1620	ENDIF
1621	!
1622	IF(ABS(DZB)>SLOPAC)THEN
1623	SSB=DZB*SQP
1624	IF(SSB>CRIT)THEN
1625	SQP=0. !sqp*.1
1626	ENDIF
1627	ENDIF
1628	!
1629	ENDIF
1630	!
1631	!-----------------------------------------------------------------------
1632	!
1633	FPQ=SPPSQP0.25
1634	PP=ABS(SPP)
1635	QP=ABS(SQP)
1636	!
1637	AFP(I,J,K)=(((FF4PP+FF3)PP+FF2)PP+FF1)PP
1638	AFQ(I,J,K)=(((FF4QP+FF3)QP+FF2)QP+FF1)QP
1639	!
1640	Q1(I,J,K)=(Q (IFPA(I,J,K),JFPA(I,J,K),K)-Q (I,J,K))*PP &
1641	& +(Q (IFQA(I,J,K),JFQA(I,J,K),K)-Q (I,J,K))*QP &
1642	& +(Q (I,J-2,K)+Q (I,J+2,K) &
1643	& -Q (I-1,J,K)-Q (I+1,J,K))*FPQ &
1644	& +Q(I,J,K)
1645	!
1646	W1(I,J,K)=(CWM(IFPA(I,J,K),JFPA(I,J,K),K)-CWM(I,J,K))*PP &
1647	& +(CWM(IFQA(I,J,K),JFQA(I,J,K),K)-CWM(I,J,K))*QP &
1648	& +(CWM(I,J-2,K)+CWM(I,J+2,K) &
1649	& -CWM(I-1,J,K)-CWM(I+1,J,K))*FPQ &
1650	& +CWM(I,J,K)
1651	!
1652	E2(I,J,K)=(E1 (IFPA(I,J,K),JFPA(I,J,K),K)-E1 (I,J,K))*PP &
1653	& +(E1 (IFQA(I,J,K),JFQA(I,J,K),K)-E1 (I,J,K))*QP &
1654	& +(E1 (I,J-2,K)+E1 (I,J+2,K) &
1655	& -E1 (I-1,J,K)-E1 (I+1,J,K))*FPQ &
1656	& +E1(I,J,K)
1657	!
1658	ENDDO
1659	ENDDO
1660	!
1661	!-----------------------------------------------------------------------
1662	!
1663	ENDDO vertical_1
1664	!
1665	!-----------------------------------------------------------------------
1666	!*** ANTI-FILTERING STEP
1667	!-----------------------------------------------------------------------
1668	!
1669	DO K=KTS,KTE
1670	XSUMS(1,K)=0.
1671	XSUMS(2,K)=0.
1672	XSUMS(3,K)=0.
1673	XSUMS(4,K)=0.
1674	XSUMS(5,K)=0.
1675	XSUMS(6,K)=0.
1676	ENDDO
1677	!-----------------------------------------------------------------------
1678	!
1679	!*** ANTI-FILTERING LIMITERS
1680	!
1681	!-----------------------------------------------------------------------
1682	#if defined(BIT_FOR_BIT) && defined(DM_PARALLEL) && !defined(STUBMPI)
1683	DO N=1,6
1684	!
1685	!$omp parallel do &
1686	!$omp& private(i,j,k)
1687	DO K=KMS,KME
1688	DO J=JMS,JME
1689	DO I=IMS,IME
1690	XSUMS_L(I,J,K,N)=0.
1691	ENDDO
1692	ENDDO
1693	ENDDO
1694	!
1695	!$omp parallel do &
1696	!$omp& private(i,j,k)
1697	DO K=KDS,KDE
1698	DO J=JDS,JDE
1699	DO I=IDS,IDE
1700	XSUMS_G(I,J,K,N)=0.
1701	ENDDO
1702	ENDDO
1703	ENDDO
1704	!
1705	ENDDO
1706	!
1707	#endif
1708	!-----------------------------------------------------------------------
1709	!$omp parallel do &
1710	!$omp& private(d2pqe,d2pqq,d2pqw,destij,dqstij,dvolp,dwstij &
1711	!$omp& ,e00,e0q,e2ij,ep0,estij,hafp,hafq,i,j,k &
1712	!$omp& ,q00,q0q,q1ij,qp0,qstij,w00,w0q,w1ij,wp0,wstij)
1713	!-----------------------------------------------------------------------
1714	!
1715	vertical_2: DO K=KTS,KTE
1716	!
1717	!-----------------------------------------------------------------------
1718	!
1719	DO J=MYJS2,MYJE2
1720	DO I=MYIS1,MYIE1
1721	!
1722	DVOLP=DVOL(I,J,K)
1723	Q1IJ =Q1(I,J,K)
1724	W1IJ =W1(I,J,K)
1725	E2IJ =E2(I,J,K)
1726	!
1727	HAFP=AFP(I,J,K)
1728	HAFQ=AFQ(I,J,K)
1729	!
1730	D2PQQ=(Q1(IFPA(I,J,K),JFPA(I,J,K),K)-Q1IJ &
1731	& -Q1IJ+Q1(IFPF(I,J,K),JFPF(I,J,K),K)) &
1732	& *HAFP &
1733	& +(Q1(IFQA(I,J,K),JFQA(I,J,K),K)-Q1IJ &
1734	& -Q1IJ+Q1(IFQF(I,J,K),JFQF(I,J,K),K)) &
1735	& *HAFQ
1736	!
1737	D2PQW=(W1(IFPA(I,J,K),JFPA(I,J,K),K)-W1IJ &
1738	& -W1IJ+W1(IFPF(I,J,K),JFPF(I,J,K),K)) &
1739	& *HAFP &
1740	& +(W1(IFQA(I,J,K),JFQA(I,J,K),K)-W1IJ &
1741	& -W1IJ+W1(IFQF(I,J,K),JFQF(I,J,K),K)) &
1742	& *HAFQ
1743	!
1744	D2PQE=(E2(IFPA(I,J,K),JFPA(I,J,K),K)-E2IJ &
1745	& -E2IJ+E2(IFPF(I,J,K),JFPF(I,J,K),K)) &
1746	& *HAFP &
1747	& +(E2(IFQA(I,J,K),JFQA(I,J,K),K)-E2IJ &
1748	& -E2IJ+E2(IFQF(I,J,K),JFQF(I,J,K),K)) &
1749	& *HAFQ
1750	!
1751	QSTIJ=Q1IJ-D2PQQ
1752	WSTIJ=W1IJ-D2PQW
1753	ESTIJ=E2IJ-D2PQE
1754	!
1755	Q00=Q (I ,J ,K)
1756	QP0=Q (IFPA(I,J,K),JFPA(I,J,K),K)
1757	Q0Q=Q (IFQA(I,J,K),JFQA(I,J,K),K)
1758	!
1759	W00=CWM(I ,J ,K)
1760	WP0=CWM(IFPA(I,J,K),JFPA(I,J,K),K)
1761	W0Q=CWM(IFQA(I,J,K),JFQA(I,J,K),K)
1762	!
1763	E00=E1 (I ,J ,K)
1764	EP0=E1 (IFPA(I,J,K),JFPA(I,J,K),K)
1765	E0Q=E1 (IFQA(I,J,K),JFQA(I,J,K),K)
1766	!
1767	QSTIJ=MAX(QSTIJ,MIN(Q00,QP0,Q0Q))
1768	QSTIJ=MIN(QSTIJ,MAX(Q00,QP0,Q0Q))
1769	WSTIJ=MAX(WSTIJ,MIN(W00,WP0,W0Q))
1770	WSTIJ=MIN(WSTIJ,MAX(W00,WP0,W0Q))
1771	ESTIJ=MAX(ESTIJ,MIN(E00,EP0,E0Q))
1772	ESTIJ=MIN(ESTIJ,MAX(E00,EP0,E0Q))
1773	!
1774	! DQSTIJ=QSTIJ-Q(I,J,K)
1775	! DWSTIJ=WSTIJ-CWM(I,J,K)
1776	! DESTIJ=ESTIJ-E1(I,J,K)
1777	!
1778	dqstij=qstij-q1(i,j,k)
1779	dwstij=wstij-w1(i,j,k)
1780	destij=estij-e2(i,j,k)
1781	!
1782	DQST(I,J,K)=DQSTIJ
1783	DWST(I,J,K)=DWSTIJ
1784	DEST(I,J,K)=DESTIJ
1785	!
1786	DQSTIJ=DQSTIJ*DVOLP
1787	DWSTIJ=DWSTIJ*DVOLP
1788	DESTIJ=DESTIJ*DVOLP
1789	!
1790	!-----------------------------------------------------------------------
1791	#if defined(BIT_FOR_BIT) && defined(DM_PARALLEL) && !defined(STUBMPI)
1792	!-----------------------------------------------------------------------
1793	DO N=1,6
1794	XSUMS_L(I,J,K,N)=0.
1795	ENDDO
1796	!
1797	IF(DQSTIJ>0.)THEN
1798	XSUMS_L(I,J,K,1)=DQSTIJ
1799	ELSE
1800	XSUMS_L(I,J,K,2)=DQSTIJ
1801	ENDIF
1802	!
1803	IF(DWSTIJ>0.)THEN
1804	XSUMS_L(I,J,K,3)=DWSTIJ
1805	ELSE
1806	XSUMS_L(I,J,K,4)=DWSTIJ
1807	ENDIF
1808	!
1809	IF(DESTIJ>0.)THEN
1810	XSUMS_L(I,J,K,5)=DESTIJ
1811	ELSE
1812	XSUMS_L(I,J,K,6)=DESTIJ
1813	ENDIF
1814	!-----------------------------------------------------------------------
1815	#else
1816	!-----------------------------------------------------------------------
1817	IF(DQSTIJ>0.)THEN
1818	XSUMS(1,K)=XSUMS(1,K)+DQSTIJ
1819	ELSE
1820	XSUMS(2,K)=XSUMS(2,K)+DQSTIJ
1821	ENDIF
1822	!
1823	IF(DWSTIJ>0.)THEN
1824	XSUMS(3,K)=XSUMS(3,K)+DWSTIJ
1825	ELSE
1826	XSUMS(4,K)=XSUMS(4,K)+DWSTIJ
1827	ENDIF
1828	!
1829	IF(DESTIJ>0.)THEN
1830	XSUMS(5,K)=XSUMS(5,K)+DESTIJ
1831	ELSE
1832	XSUMS(6,K)=XSUMS(6,K)+DESTIJ
1833	ENDIF
1834	!-----------------------------------------------------------------------
1835	#endif
1836	!-----------------------------------------------------------------------
1837	!
1838	ENDDO
1839	ENDDO
1840	!
1841	!-----------------------------------------------------------------------
1842	!
1843	ENDDO vertical_2
1844	!
1845	!-----------------------------------------------------------------------
1846	#if defined(BIT_FOR_BIT) && defined(DM_PARALLEL) && !defined(STUBMPI)
1847	!-----------------------------------------------------------------------
1848	DO N=1,6
1849	CALL WRF_PATCH_TO_GLOBAL_REAL(XSUMS_L(IMS,JMS,KMS,N) &
1850	&, XSUMS_G(1,1,1,N),DOMDESC &
1851	&, 'xyz','xyz' &
1852	&, IDS,IDE,JDS,JDE,KDS,KDE &
1853	&, IMS,IME,JMS,JME,KMS,KME &
1854	&, ITS,ITE,JTS,JTE,KTS,KTE)
1855	ENDDO
1856	!
1857	DO K=KTS,KTE
1858	DO N=1,6
1859	GSUMS(N,K)=0.
1860	ENDDO
1861	ENDDO
1862	!
1863	IF(WRF_DM_ON_MONITOR())THEN
1864	DO N=1,6
1865	!$omp parallel do &
1866	!$omp& private(i,j,k)
1867	DO K=KTS,KTE
1868	DO J=JDS,JDE
1869	DO I=IDS,IDE
1870	GSUMS(N,K)=GSUMS(N,K)+XSUMS_G(I,J,K,N)
1871	ENDDO
1872	ENDDO
1873	ENDDO
1874	ENDDO
1875	ENDIF
1876
1877	CALL WRF_DM_BCAST_BYTES(GSUMS,2RWORDSIZE6*(KTE-KTS+1) )
1878
1879	!-----------------------------------------------------------------------
1880	#else
1881	!-----------------------------------------------------------------------
1882	!
1883	!-----------------------------------------------------------------------
1884	!*** GLOBAL REDUCTION
1885	!-----------------------------------------------------------------------
1886	!
1887	# if defined(DM_PARALLEL) && !defined(STUBMPI)
1888	CALL WRF_GET_DM_COMMUNICATOR(MPI_COMM_COMP)
1889	CALL MPI_ALLREDUCE(XSUMS,GSUMS,6*(KTE-KTS+1) &
1890	& ,MPI_DOUBLE_PRECISION,MPI_SUM &
1891	& ,MPI_COMM_COMP,IRECV)
1892	# else
1893	DO K=KTS,KTE
1894	DO N=1,6
1895	GSUMS(N,K)=XSUMS(N,K)
1896	ENDDO
1897	ENDDO
1898	# endif
1899	!
1900	!-----------------------------------------------------------------------
1901	#endif
1902	!-----------------------------------------------------------------------
1903	!
1904	!-----------------------------------------------------------------------
1905	!*** END OF GLOBAL REDUCTION
1906	!-----------------------------------------------------------------------
1907	!
1908	! if(mype==0)then
1909	!!! if(ntsd==0)then
1910	!!! call int_get_fresh_handle(nunit)
1911	!!! close(nunit)
1912	! nunit=56
1913	!!! open(unit=nunit,file='gsums',form='unformatted',iostat=ier)
1914	!!! endif
1915	! endif
1916	!-----------------------------------------------------------------------
1917	!$omp parallel do &
1918	!$omp& private(destij,dqstij,dwstij,i,j,k,rface,rfacq,rfacw &
1919	!$omp& ,rfeij,rfqij,rfwij,sumne,sumnq,sumnw,sumpe,sumpq,sumpw)
1920	!-----------------------------------------------------------------------
1921	!
1922	vertical_3: DO K=KTS,KTE
1923	!
1924	!-----------------------------------------------------------------------
1925	! if(mype==0)then
1926	! write(nunit)(gsums(i,k),i=1,6)
1927	! endif
1928	!!! read(nunit)(gsums(i,k),i=1,6)
1929	!-----------------------------------------------------------------------
1930	!
1931	SUMPQ=GSUMS(1,K)
1932	SUMNQ=GSUMS(2,K)
1933	SUMPW=GSUMS(3,K)
1934	SUMNW=GSUMS(4,K)
1935	SUMPE=GSUMS(5,K)
1936	SUMNE=GSUMS(6,K)
1937	!
1938	!-----------------------------------------------------------------------
1939	!*** FIRST MOMENT CONSERVING FACTOR
1940	!-----------------------------------------------------------------------
1941	!
1942	IF(SUMPQ>1.)THEN
1943	RFACQ=-SUMNQ/SUMPQ
1944	ELSE
1945	RFACQ=1.
1946	ENDIF
1947	!
1948	IF(SUMPW>1.)THEN
1949	RFACW=-SUMNW/SUMPW
1950	ELSE
1951	RFACW=1.
1952	ENDIF
1953	!
1954	IF(SUMPE>1.)THEN
1955	RFACE=-SUMNE/SUMPE
1956	ELSE
1957	RFACE=1.
1958	ENDIF
1959	!
1960	IF(RFACQ<CONSERVE_MIN.OR.RFACQ>CONSERVE_MAX)RFACQ=1.
1961	IF(RFACW<CONSERVE_MIN.OR.RFACW>CONSERVE_MAX)RFACW=1.
1962	IF(RFACE<CONSERVE_MIN.OR.RFACE>CONSERVE_MAX)RFACE=1.
1963	!
1964	!-----------------------------------------------------------------------
1965	! if(mype==0.and.ntsd==181)close(nunit)
1966	!-----------------------------------------------------------------------
1967	!
1968	!-----------------------------------------------------------------------
1969	!*** IMPOSE CONSERVATION ON ANTI-FILTERING
1970	!-----------------------------------------------------------------------
1971	!
1972	if(rfacq<1.)then
1973	do j=MYJS2,MYJE2
1974	DO I=MYIS1,MYIE1
1975	DQSTIJ=DQST(I,J,K)
1976	RFQIJ=HBM2(I,J)*(RFACQ-1.)+1.
1977	IF(DQSTIJ>=0.)DQSTIJ=DQSTIJ*RFQIJ
1978	q(i,j,k)=q1(i,j,k)+dqstij
1979	ENDDO
1980	enddo
1981	else
1982	do j=MYJS2,MYJE2
1983	DO I=MYIS1,MYIE1
1984	DQSTIJ=DQST(I,J,K)
1985	RFQIJ=HBM2(I,J)*(RFACQ-1.)+1.
1986	IF(DQSTIJ<0.)DQSTIJ=DQSTIJ/RFQIJ
1987	q(i,j,k)=q1(i,j,k)+dqstij
1988	ENDDO
1989	enddo
1990	endif
1991	!
1992	!-----------------------------------------------------------------------
1993	!
1994	if(rfacw<1.)then
1995	do j=MYJS2,MYJE2
1996	DO I=MYIS1,MYIE1
1997	DWSTIJ=DWST(I,J,K)
1998	RFWIJ=HBM2(I,J)*(RFACW-1.)+1.
1999	IF(DWSTIJ>=0.)DWSTIJ=DWSTIJ*RFWIJ
2000	cwm(i,j,k)=w1(i,j,k)+dwstij
2001	ENDDO
2002	enddo
2003	else
2004	do j=MYJS2,MYJE2
2005	DO I=MYIS1,MYIE1
2006	DWSTIJ=DWST(I,J,K)
2007	RFWIJ=HBM2(I,J)*(RFACW-1.)+1.
2008	IF(DWSTIJ<0.)DWSTIJ=DWSTIJ/RFWIJ
2009	cwm(i,j,k)=w1(i,j,k)+dwstij
2010	ENDDO
2011	enddo
2012	endif
2013
2014	!-----------------------------------------------------------------------
2015	!
2016	if(rface<1.)then
2017	do j=MYJS2,MYJE2
2018	DO I=MYIS1,MYIE1
2019	DESTIJ=DEST(I,J,K)
2020	RFEIJ=HBM2(I,J)*(RFACE-1.)+1.
2021	IF(DESTIJ>=0.)DESTIJ=DESTIJ*RFEIJ
2022	e1(i,j,k)=e2(i,j,k)+destij
2023	ENDDO
2024	enddo
2025	else
2026	do j=MYJS2,MYJE2
2027	DO I=MYIS1,MYIE1
2028	DESTIJ=DEST(I,J,K)
2029	RFEIJ=HBM2(I,J)*(RFACE-1.)+1.
2030	IF(DESTIJ<0.)DESTIJ=DESTIJ/RFEIJ
2031	e1(i,j,k)=e2(i,j,k)+destij
2032	ENDDO
2033	enddo
2034	endif
2035	!
2036	!-----------------------------------------------------------------------
2037	!
2038	DO J=MYJS,MYJE
2039	DO I=MYIS,MYIE
2040	Q (I,J,K)=MAX(Q (I,J,K),EPSQ)
2041	CWM(I,J,K)=MAX(CWM(I,J,K),CLIMIT)
2042	ENDDO
2043	ENDDO
2044	!
2045	!-----------------------------------------------------------------------
2046	!
2047	ENDDO vertical_3
2048	!
2049	!-----------------------------------------------------------------------
2050	!
2051	!$omp parallel do &
2052	!$omp& private(i,j)
2053	DO J=MYJS,MYJE
2054	DO I=MYIS,MYIE
2055	Q2(I,J,KTE)=MAX(E1(I,J,KTE)+E1(I,J,KTE)-EPSQ2,EPSQ2)
2056	ENDDO
2057	ENDDO
2058	!
2059	!-----------------------------------------------------------------------
2060	!
2061	DO K=KTE-1,KTS+1,-1
2062	!$omp parallel do &
2063	!$omp& private(i,j)
2064	DO J=MYJS,MYJE
2065	DO I=MYIS,MYIE
2066	IF(K>KTS)THEN
2067	Q2(I,J,K)=MAX(E1(I,J,K)+E1(I,J,K)-Q2(I,J,K+1),EPSQ2)
2068	ELSE
2069	Q2(I,J,K)=Q2(I,J,K+1)
2070	ENDIF
2071	ENDDO
2072	ENDDO
2073	ENDDO
2074	!-----------------------------------------------------------------------
2075	!
2076	END SUBROUTINE HAD2
2077	!
2078	!-----------------------------------------------------------------------
2079	!***********************************************************************
2080	!-----------------------------------------------------------------------
2081	! New routines added by Georg Grell to handle advection more like ARW
2082	! core. Instead of VAD2/HAD2 that advect TKE, specific humidity, and
2083	! condensed water species all in one routine, we call VAD2/HAD2_SCAL
2084	! with multidimensioned arrays to advect each variable. For purposes
2085	! here, solve_nmm.F calls this routine once for TKE, then again for
2086	! all the species held in the moist array (qv, qc, qi, qr, qs, qg),
2087	! then call again for number concentrations held in scalar array (qni).
2088	! The dummy argument lstart is the starting index of the multidimensioned
2089	! array for starting the advection since the 1st index of moist and
2090	! scalar are actually empty placeholders (and the 2nd element is vapor,
2091	! then qc, etc.) When calling with single 3D array (like TKE), just
2092	! set NUM_SCAL=1 and lstart=1. The variable to advect is called SCAL
2093	! herein.
2094	!***********************************************************************
2095	SUBROUTINE VAD2_SCAL(NTSD,DT,IDTAD,DX,DY &
2096	& ,AETA1,AETA2,DETA1,DETA2,PDSL,PDTOP &
2097	& ,HBM2 &
2098	& ,SCAL,PETDT &
2099	& ,N_IUP_H,N_IUP_V &
2100	& ,N_IUP_ADH,N_IUP_ADV &
2101	& ,IUP_H,IUP_V,IUP_ADH,IUP_ADV &
2102	& ,IHE,IHW,IVE,IVW &
2103	& ,NUM_SCAL,LSTART &
2104	& ,IDS,IDE,JDS,JDE,KDS,KDE &
2105	& ,IMS,IME,JMS,JME,KMS,KME &
2106	& ,ITS,ITE,JTS,JTE,KTS,KTE)
2107	!***********************************************************************
2108	!$$$ SUBPROGRAM DOCUMENTATION BLOCK
2109	! . . .
2110	! SUBPROGRAM: VAD2_SCAL VERTICAL ADVECTION OF SCALARS
2111	!
2112	! PRGRMMR: JANJIC ORG: W/NP22 DATE: 96-07-19
2113	! GRELL,PECKHAM ORG: NOAA/FSL DATE: 05-02-03
2114	!
2115	! ABSTRACT:
2116	! VAD2_SCAL CALCULATES THE CONTRIBUTION OF THE VERTICAL ADVECTION
2117	! TO THE TENDENCIES OF SCALAR SUBSTANCES AND THEN UPDATES
2118	! THOSE VARIABLES. AN ANTI-FILTERING TECHNIQUE IS USED.
2119	!
2120	! PROGRAM HISTORY LOG:
2121	! 96-07-19 JANJIC - ORIGINATOR
2122	! 05-02-03 GRELL,PECKHAM - MODIFIED FOR SCALARS
2123	!
2124	! USAGE: CALL VAD2_SCAL FROM SUBROUTINE SOLVE_NMM
2125	! INPUT ARGUMENT LIST:
2126	!
2127	! OUTPUT ARGUMENT LIST
2128	!
2129	! OUTPUT FILES:
2130	! NONE
2131	! SUBPROGRAMS CALLED:
2132	!
2133	! UNIQUE: NONE
2134	!
2135	! LIBRARY: NONE
2136	!
2137	! ATTRIBUTES:
2138	! LANGUAGE: FORTRAN 90
2139	! MACHINE : IBM
2140	!$$$
2141	!***********************************************************************
2142	!----------------------------------------------------------------------
2143	!
2144	IMPLICIT NONE
2145	!
2146	!----------------------------------------------------------------------
2147	!
2148	INTEGER,INTENT(IN) :: IDS,IDE,JDS,JDE,KDS,KDE &
2149	& ,IMS,IME,JMS,JME,KMS,KME &
2150	,ITS,ITE,JTS,JTE,KTS,KTE
2151	!
2152	INTEGER,INTENT(IN) :: LSTART,NUM_SCAL
2153	!
2154	INTEGER,DIMENSION(JMS:JME),INTENT(IN) :: IHE,IHW,IVE,IVW
2155	INTEGER,DIMENSION(JMS:JME),INTENT(IN) :: N_IUP_H,N_IUP_V &
2156	& ,N_IUP_ADH,N_IUP_ADV
2157	INTEGER,DIMENSION(IMS:IME,JMS:JME),INTENT(IN) :: IUP_H,IUP_V &
2158	& ,IUP_ADH,IUP_ADV
2159	!
2160	INTEGER,INTENT(IN) :: IDTAD,NTSD
2161	!
2162	REAL,INTENT(IN) :: DT,DY,PDTOP
2163	!
2164	REAL,DIMENSION(KMS:KME),INTENT(IN) :: AETA1,AETA2,DETA1,DETA2
2165	!
2166	REAL,DIMENSION(IMS:IME,JMS:JME),INTENT(IN) :: DX,HBM2,PDSL
2167	!
2168	REAL,DIMENSION(IMS:IME,JMS:JME,KMS:KME),INTENT(IN) :: PETDT
2169	!
2170	REAL,DIMENSION(IMS:IME,JMS:JME,KMS:KME,1:NUM_SCAL) &
2171	,INTENT(INOUT) :: SCAL
2172	!
2173	!----------------------------------------------------------------------
2174	!*** LOCAL VARIABLES
2175	!----------------------------------------------------------------------
2176	!
2177	REAL,PARAMETER :: FF1=0.500
2178	!
2179	LOGICAL,SAVE :: TRADITIONAL=.TRUE.
2180	!
2181	INTEGER :: I,IRECV,J,JFP,JFQ,K,L,LAP,LLAP
2182	!
2183	INTEGER,DIMENSION(KTS:KTE) :: LA
2184	!
2185	REAL*8 :: ADDT,AFRP,D2PQQ,DETAP,DPDN,DPUP,DQP &
2186	& ,HADDT,HBM2IJ &
2187	& ,Q00,Q4P,QP,QP0 &
2188	& ,RFACQK,RFC,RR &
2189	& ,SUMNQ,SUMPQ
2190	!
2191	REAL :: SFACQK
2192	!
2193	REAL,DIMENSION(KTS:KTE) :: AFR,DEL,DQL,DWL,E3,E4,PETDTK &
2194	& ,RFACE,RFACQ,RFACW,Q3,Q4,W3,W4
2195	!
2196	!-----------------------------------------------------------------------
2197	!***********************************************************************
2198	!-----------------------------------------------------------------------
2199	!
2200	ADDT=REAL(IDTAD)*DT
2201	!
2202	!-----------------------------------------------------------------------
2203	!$omp parallel do &
2204	!$omp& private(afr,afrp,d2pqq,detap,dpdn,dpup &
2205	!$omp& ,dql,dqp,haddt,i,j,k &
2206	!$omp& ,la,lap,llap,petdtk,q00,q3,q4,q4p,qp,qp0,rfacqk &
2207	!$omp& ,rfc,rr,sfacqk,sumnq,sumpq)
2208	!-----------------------------------------------------------------------
2209	!
2210	scalar_loop: DO L=LSTART,NUM_SCAL
2211	!
2212	main_integration: DO J=MYJS2,MYJE2
2213	!
2214	!-----------------------------------------------------------------------
2215	!
2216	main_iloop: DO I=MYIS1_P1,MYIE1_P1
2217	!
2218	!-----------------------------------------------------------------------
2219	!
2220	DO K=KTS,KTE
2221	Q3(K)=SCAL(I,J,K,L)
2222	Q4(K)=Q3(K)
2223	ENDDO
2224	!
2225	IF(TRADITIONAL)THEN
2226	PETDTK(KTE)=PETDT(I,J,KTE-1)*0.5
2227	!
2228	DO K=KTE-1,KTS+1,-1
2229	PETDTK(K)=(PETDT(I,J,K)+PETDT(I,J,K-1))*0.5
2230	ENDDO
2231	!
2232	PETDTK(KTS)=PETDT(I,J,KTS)*0.5
2233	!
2234	ELSE
2235	!
2236	!-----------------------------------------------------------------------
2237	!*** PERFORM HORIZONTAL AVERAGING OF VERTICAL VELOCITY
2238	!-----------------------------------------------------------------------
2239	!
2240	PETDTK(KTE)=(PETDT(I+IHW(J-1),J-1,KTE-1) &
2241	& +PETDT(I+IHE(J-1),J-1,KTE-1) &
2242	& +PETDT(I+IHW(J+1),J+1,KTE-1) &
2243	& +PETDT(I+IHE(J+1),J+1,KTE-1) &
2244	& +PETDT(I,J,KTE-1)4. )0.0625
2245	!
2246	DO K=KTE-1,KTS+1,-1
2247	PETDTK(K)=(PETDT(I+IHW(J-1),J-1,K-1) &
2248	+PETDT(I+IHE(J-1),J-1,K-1) &
2249	& +PETDT(I+IHW(J+1),J+1,K-1) &
2250	& +PETDT(I+IHE(J+1),J+1,K-1) &
2251	& +PETDT(I+IHW(J-1),J-1,K ) &
2252	& +PETDT(I+IHE(J-1),J-1,K ) &
2253	& +PETDT(I+IHW(J+1),J+1,K ) &
2254	& +PETDT(I+IHE(J+1),J+1,K ) &
2255	& +(PETDT(I,J,K-1)+PETDT(I,J,K))*4. &
2256	& )*0.0625
2257	ENDDO
2258	!
2259	PETDTK(KTS)=(PETDT(I+IHW(J-1),J-1,KTS) &
2260	& +PETDT(I+IHE(J-1),J-1,KTS) &
2261	& +PETDT(I+IHW(J+1),J+1,KTS) &
2262	& +PETDT(I+IHE(J+1),J+1,KTS) &
2263	& +PETDT(I,J,KTS)4. )0.0625
2264
2265	ENDIF
2266	!
2267	!-----------------------------------------------------------------------
2268	!
2269	HADDT=-ADDT*HBM2(I,J)
2270	!
2271	DO K=KTE,KTS,-1
2272	RR=PETDTK(K)*HADDT
2273	!
2274	IF(RR<0.)THEN
2275	LAP=1
2276	ELSE
2277	LAP=-1
2278	ENDIF
2279	!
2280	LA(K)=LAP
2281	LLAP=K+LAP
2282	!
2283	IF(LLAP>KTS-1.AND.LLAP<KTE+1)THEN
2284	RR=ABS(RR/((AETA1(LLAP)-AETA1(K))*PDTOP &
2285	& +(AETA2(LLAP)-AETA2(K))*PDSL(I,J)))
2286	IF(RR>0.9)RR=0.9
2287	!
2288	AFR(K)=(((FF4RR+FF3)RR+FF2)RR+FF1)RR
2289	DQP=(Q3(LLAP)-Q3(K))*RR
2290	DQL(K)=DQP
2291	ELSE
2292	RR=0.
2293	AFR(K)=0.
2294	DQL(K)=0.
2295	ENDIF
2296	ENDDO
2297	!
2298	!-----------------------------------------------------------------------
2299	!
2300	IF(LA(KTE-1)>0)THEN
2301	RFC=(DETA1(KTE-1)PDTOP+DETA2(KTE-1)PDSL(I,J)) &
2302	& /(DETA1(KTE )PDTOP+DETA2(KTE )PDSL(I,J))
2303	DQL(KTE)=-DQL(KTE-1)*RFC
2304	ENDIF
2305	!
2306	IF(LA(KTS+1)<0)THEN
2307	RFC=(DETA1(KTS+1)PDTOP+DETA2(KTS+1)PDSL(I,J)) &
2308	& /(DETA1(KTS )PDTOP+DETA2(KTS )PDSL(I,J))
2309	DQL(KTS)=-DQL(KTS+1)*RFC
2310	ENDIF
2311	!
2312	DO K=KTS,KTE
2313	Q4(K)=Q3(K)+DQL(K)
2314	ENDDO
2315	!
2316	!-----------------------------------------------------------------------
2317	!*** ANTI-FILTERING STEP
2318	!-----------------------------------------------------------------------
2319	!
2320	SUMPQ=0.
2321	SUMNQ=0.
2322	!
2323	!*** ANTI-FILTERING LIMITERS
2324	!
2325	antifilter: DO K=KTE-1,KTS+1,-1
2326	!
2327	DETAP=DETA1(K)PDTOP+DETA2(K)PDSL(I,J)
2328	!
2329	Q4P=Q4(K)
2330	!
2331	LAP=LA(K)
2332	!
2333	DPDN=(AETA1(K+LAP)-AETA1(K))*PDTOP &
2334	& +(AETA2(K+LAP)-AETA2(K))*PDSL(I,J)
2335	DPUP=(AETA1(K)-AETA1(K-LAP))*PDTOP &
2336	& +(AETA2(K)-AETA2(K-LAP))*PDSL(I,J)
2337	!
2338	AFRP=2.AFR(K)DPDN*DPDN/(DPDN+DPUP)
2339	D2PQQ=((Q4(K+LAP)-Q4P)/DPDN &
2340	& -(Q4P-Q4(K-LAP))/DPUP)*AFRP
2341	!
2342	QP=Q4P-D2PQQ
2343	!
2344	Q00=Q3(K)
2345	QP0=Q3(K+LAP)
2346	!
2347	QP=MAX(QP,MIN(Q00,QP0))
2348	QP=MIN(QP,MAX(Q00,QP0))
2349	!
2350	DQP=QP-Q00
2351	!
2352	DQL(K)=DQP
2353	!
2354	ENDDO antifilter
2355	!
2356	!-----------------------------------------------------------------------
2357	!
2358	IF(LA(KTE-1)>0)THEN
2359	RFC=(DETA1(KTE-1)PDTOP+DETA2(KTE-1)PDSL(I,J)) &
2360	& /(DETA1(KTE )PDTOP+DETA2(KTE )PDSL(I,J))
2361	DQL(KTE)=-DQL(KTE-1)*RFC+DQL(KTE)
2362	ENDIF
2363	!
2364	IF(LA(KTS+1)<0)THEN
2365	RFC=(DETA1(KTS+1)PDTOP+DETA2(KTS+1)PDSL(I,J)) &
2366	& /(DETA1(KTS )PDTOP+DETA2(KTS )PDSL(I,J))
2367	DQL(KTS)=-DQL(KTS+1)*RFC+DQL(KTS)
2368	ENDIF
2369	!
2370	DO K=KTS,KTE
2371	DETAP=DETA1(K)PDTOP+DETA2(K)PDSL(I,J)
2372	DQP=DQL(K)*DETAP
2373	!
2374	IF(DQP>0.)THEN
2375	SUMPQ=SUMPQ+DQP
2376	ELSE
2377	SUMNQ=SUMNQ+DQP
2378	ENDIF
2379	ENDDO
2380	!
2381	!-----------------------------------------------------------------------
2382	!*** FIRST MOMENT CONSERVING FACTOR
2383	!-----------------------------------------------------------------------
2384	!
2385	IF(SUMPQ>1.E-9)THEN
2386	SFACQK=-SUMNQ/SUMPQ
2387	ELSE
2388	SFACQK=1.
2389	ENDIF
2390	!
2391	IF(SFACQK<CONSERVE_MIN.OR.SFACQK>CONSERVE_MAX)SFACQK=1.
2392	!
2393	RFACQK=1./SFACQK
2394	!
2395	!-----------------------------------------------------------------------
2396	!*** IMPOSE CONSERVATION ON ANTI-FILTERING
2397	!-----------------------------------------------------------------------
2398	!
2399	DO K=KTE,KTS,-1
2400	DQP=DQL(K)
2401	IF(SFACQK>=1.)THEN
2402	IF(DQP<0.)DQP=DQP*RFACQK
2403	ELSE
2404	IF(DQP>0.)DQP=DQP*SFACQK
2405	ENDIF
2406	SCAL(I,J,K,L)=Q3(K)+DQP
2407	ENDDO
2408	!
2409	!-----------------------------------------------------------------------
2410	!
2411	ENDDO main_iloop
2412	!
2413	!-----------------------------------------------------------------------
2414	!
2415	ENDDO main_integration
2416	!
2417	!-----------------------------------------------------------------------
2418	!
2419	ENDDO scalar_loop
2420	!
2421	!-----------------------------------------------------------------------
2422	!
2423	END SUBROUTINE VAD2_SCAL
2424	!
2425	!-----------------------------------------------------------------------
2426	!
2427	!***********************************************************************
2428	SUBROUTINE HAD2_SCAL( &
2429	#if defined(DM_PARALLEL)
2430	& DOMDESC , &
2431	#endif
2432	& NTSD,DT,IDTAD,DX,DY &
2433	& ,AETA1,AETA2,DETA1,DETA2,PDSL,PDTOP &
2434	& ,HBM2,HBM3 &
2435	& ,SCAL,U,V,Z,HYDRO &
2436	& ,N_IUP_H,N_IUP_V &
2437	& ,N_IUP_ADH,N_IUP_ADV &
2438	& ,IUP_H,IUP_V,IUP_ADH,IUP_ADV &
2439	& ,IHE,IHW,IVE,IVW &
2440	& ,NUM_SCAL,LSTART &
2441	& ,IDS,IDE,JDS,JDE,KDS,KDE &
2442	& ,IMS,IME,JMS,JME,KMS,KME &
2443	& ,ITS,ITE,JTS,JTE,KTS,KTE)
2444	!***********************************************************************
2445	!$$$ SUBPROGRAM DOCUMENTATION BLOCK
2446	! . . .
2447	! SUBPROGRAM: HAD2_SCAL HORIZONTAL ADVECTION OF SCALAR
2448	! PRGRMMR: JANJIC ORG: W/NP22 DATE: 96-07-19
2449	! GRELL,PECKHAM ORG: NOAA/FSL DATE: 05-02-03
2450	!
2451	! ABSTRACT:
2452	! HAD2_SCAL CALCULATES THE CONTRIBUTION OF THE HORIZONTAL ADVECTION
2453	! TO THE TENDENCIES OF SCALAR SUBSTANCES AND THEN
2454	! UPDATES THOSE VARIABLES. AN ANTI-FILTERING TECHNIQUE IS USED.
2455	!
2456	! PROGRAM HISTORY LOG:
2457	! 96-07-19 JANJIC - ORIGINATOR
2458	! 05-01-03 GRELL,PECKHAM - MODIFIED FOR SCALAR
2459	!
2460	! USAGE: CALL HAD2_SCAL FROM SUBROUTINE SOLVE_NMM
2461	! INPUT ARGUMENT LIST:
2462	!
2463	! OUTPUT ARGUMENT LIST
2464	!
2465	! OUTPUT FILES:
2466	! NONE
2467	! SUBPROGRAMS CALLED:
2468	!
2469	! UNIQUE: NONE
2470	!
2471	! LIBRARY: NONE
2472	!
2473	! ATTRIBUTES:
2474	! LANGUAGE: FORTRAN 90
2475	! MACHINE : IBM
2476	!$$$
2477	!***********************************************************************
2478	!-----------------------------------------------------------------------
2479	!
2480	IMPLICIT NONE
2481	!
2482	!-----------------------------------------------------------------------
2483	!
2484	INTEGER,INTENT(IN) :: IDS,IDE,JDS,JDE,KDS,KDE &
2485	& ,IMS,IME,JMS,JME,KMS,KME &
2486	& ,ITS,ITE,JTS,JTE,KTS,KTE
2487	!
2488	INTEGER,DIMENSION(JMS:JME),INTENT(IN) :: IHE,IHW,IVE,IVW
2489	INTEGER,DIMENSION(JMS:JME),INTENT(IN) :: N_IUP_H,N_IUP_V &
2490	& ,N_IUP_ADH,N_IUP_ADV
2491	INTEGER,DIMENSION(IMS:IME,JMS:JME),INTENT(IN) :: IUP_H,IUP_V &
2492	& ,IUP_ADH,IUP_ADV
2493	!
2494	!-----------------------------------------------------------------------
2495	!
2496	INTEGER,INTENT(IN) :: IDTAD,LSTART,NTSD,NUM_SCAL
2497	!
2498	REAL,INTENT(IN) :: DT,DY,PDTOP
2499	!
2500	REAL,DIMENSION(KMS:KME),INTENT(IN) :: AETA1,AETA2,DETA1,DETA2
2501	!
2502	REAL,DIMENSION(IMS:IME,JMS:JME),INTENT(IN) :: DX,HBM2,HBM3,PDSL
2503	!
2504	REAL,DIMENSION(IMS:IME,JMS:JME,KMS:KME),INTENT(IN) :: U,V,Z
2505	!
2506	REAL,DIMENSION(IMS:IME,JMS:JME,KMS:KME,1:NUM_SCAL) &
2507	,INTENT(INOUT) :: SCAL
2508	!
2509	LOGICAL,INTENT(IN) :: HYDRO
2510	!
2511	!-----------------------------------------------------------------------
2512	!*** LOCAL VARIABLES
2513	!-----------------------------------------------------------------------
2514	!
2515	REAL,PARAMETER :: FF1=0.530
2516	!
2517	#ifdef DM_PARALLEL
2518	INTEGER :: DOMDESC
2519	#endif
2520	!
2521	#if defined(BIT_FOR_BIT) && defined(DM_PARALLEL) && !defined(STUBMPI)
2522	LOGICAL,EXTERNAL :: WRF_DM_ON_MONITOR
2523	REAL,DIMENSION(IMS:IME,JMS:JME,KMS:KME,2) :: XSUMS_L
2524	REAL,DIMENSION(IDS:IDE,JDS:JDE,KDS:KDE,2) :: XSUMS_G
2525	#endif
2526	!
2527	LOGICAL :: BOT,TOP
2528	!
2529	INTEGER :: I,IRECV,J,JFP,JFQ,K,L,LAP,LLAP,MPI_COMM_COMP
2530	INTEGER :: N
2531	!
2532	INTEGER,DIMENSION(ITS-5:ITE+5,JTS-5:JTE+5,KTS:KTE) :: IFPA,IFPF &
2533	& ,IFQA,IFQF &
2534	& ,JFPA,JFPF &
2535	& ,JFQA,JFQF
2536	!
2537	REAL :: ADDT,AFRP,CRIT,D2PQE,D2PQQ,D2PQW,DEP,DESTIJ,DQP,DQSTIJ &
2538	& ,DVOLP,DWP,DWSTIJ,DZA,DZB,E00,E0Q,E1X,E2IJ,E4P,ENH,EP,EP0 &
2539	& ,ESTIJ,FPQ,HAFP,HAFQ,HBM2IJ,HM,PP,PPQ00,Q00,Q0Q &
2540	& ,Q1IJ,Q4P,QP,QP0,QSTIJ,RDY,RFACE,RFACQ,RFACW,RFC &
2541	& ,RFEIJ,RFQIJ,RFWIJ,RR,SLOPAC,SPP,SQP,SSA,SSB,SUMNQ,SUMPQ &
2542	& ,TTA,TTB,W00,W0Q,W1IJ,W4P,WP,WP0,WSTIJ
2543	!
2544	DOUBLE PRECISION,DIMENSION(2,KTS:KTE) :: GSUMS,XSUMS
2545	!
2546	REAL,DIMENSION(KTS:KTE) :: AFR,DEL,DQL,DWL,Q3,Q4
2547	!
2548	REAL,DIMENSION(ITS-5:ITE+5,JTS-5:JTE+5) :: DARE,EMH
2549	!
2550	REAL,DIMENSION(ITS-5:ITE+5,JTS-5:JTE+5,KTS:KTE) :: AFP,AFQ,DEST &
2551	& ,DQST,DVOL,DWST &
2552	& ,Q1
2553	!
2554	REAL,DIMENSION(IMS:IME,JMS:JME,KMS:KME) :: Q
2555	!
2556	!-----------------------------------------------------------------------
2557	integer :: nunit,ier
2558	save nunit
2559	!-----------------------------------------------------------------------
2560	!***********************************************************************
2561	!-----------------------------------------------------------------------
2562	!
2563	RDY=1./DY
2564	SLOPAC=SLOPHTSQRT(2.)0.5*50.
2565	CRIT=SLOPACREAL(IDTAD)DTRDY1000.
2566	!
2567	ADDT=REAL(IDTAD)*DT
2568	ENH=ADDT/(08.*DY)
2569	!
2570	!-----------------------------------------------------------------------
2571	!$omp parallel do &
2572	!$omp& private(i,j)
2573	DO J=MYJS_P3,MYJE_P3
2574	DO I=MYIS_P2,MYIE_P2
2575	EMH (I,J)=ADDT/(08.*DX(I,J))
2576	DARE(I,J)=HBM3(I,J)DX(I,J)DY
2577	ENDDO
2578	ENDDO
2579	!-----------------------------------------------------------------------
2580	!
2581	scalar_loop: DO L=LSTART,NUM_SCAL
2582	!
2583	!-----------------------------------------------------------------------
2584	!$omp parallel do &
2585	!$omp& private(dza,dzb,e1x,fpq,hm,i,j,jfp,jfq,k,pp,qp,ssa,ssb,spp,sqp &
2586	!$omp& ,tta,ttb)
2587	!-----------------------------------------------------------------------
2588	!
2589	vertical_1: DO K=KTS,KTE
2590	!
2591	!-----------------------------------------------------------------------
2592	!
2593	DO J=MYJS_P3,MYJE_P3
2594	DO I=MYIS_P2,MYIE_P2
2595	DVOL(I,J,K)=DARE(I,J)(DETA1(K)PDTOP+DETA2(K)*PDSL(I,J))
2596	Q (I,J,K)=SCAL(I,J,K,L)
2597	Q1(I,J,K)=Q(I,J,K)
2598	ENDDO
2599	ENDDO
2600	!
2601	!-----------------------------------------------------------------------
2602	!
2603	DO J=MYJS2_P1,MYJE2_P1
2604	DO I=MYIS1_P1,MYIE1_P1
2605	!
2606	HM=HBM2(I,J)
2607	TTA=(U(I,J-1,K)+U(I+IHW(J),J,K)+U(I+IHE(J),J,K)+U(I,J+1,K)) &
2608	& EMH(I,J)HM
2609	TTB=(V(I,J-1,K)+V(I+IHW(J),J,K)+V(I+IHE(J),J,K)+V(I,J+1,K)) &
2610	& ENHHBM2(I,J)
2611	!
2612	SPP=-TTA-TTB
2613	SQP= TTA-TTB
2614	!
2615	IF(SPP<0.)THEN
2616	JFP=-1
2617	ELSE
2618	JFP=1
2619	ENDIF
2620	IF(SQP<0.)THEN
2621	JFQ=-1
2622	ELSE
2623	JFQ=1
2624	ENDIF
2625	!
2626	IFPA(I,J,K)=IHE(J)+I+( JFP-1)/2
2627	IFQA(I,J,K)=IHE(J)+I+(-JFQ-1)/2
2628	!
2629	JFPA(I,J,K)=J+JFP
2630	JFQA(I,J,K)=J+JFQ
2631	!
2632	IFPF(I,J,K)=IHE(J)+I+(-JFP-1)/2
2633	IFQF(I,J,K)=IHE(J)+I+( JFQ-1)/2
2634	!
2635	JFPF(I,J,K)=J-JFP
2636	JFQF(I,J,K)=J-JFQ
2637	!
2638	!-----------------------------------------------------------------------
2639	IF(.NOT.HYDRO)THEN ! z currently not available for hydro=.true.
2640	DZA=(Z(IFPA(I,J,K),JFPA(I,J,K),K)-Z(I,J,K))*RDY
2641	DZB=(Z(IFQA(I,J,K),JFQA(I,J,K),K)-Z(I,J,K))*RDY
2642	!
2643	IF(ABS(DZA)>SLOPAC)THEN
2644	SSA=DZA*SPP
2645	IF(SSA>CRIT)THEN
2646	SPP=0. !spp*.1
2647	ENDIF
2648	ENDIF
2649	!
2650	IF(ABS(DZB)>SLOPAC)THEN
2651	SSB=DZB*SQP
2652	IF(SSB>CRIT)THEN
2653	SQP=0. !sqp*.1
2654	ENDIF
2655	ENDIF
2656	!
2657	ENDIF
2658	!
2659	!-----------------------------------------------------------------------
2660	!
2661	FPQ=SPPSQP0.25
2662	PP=ABS(SPP)
2663	QP=ABS(SQP)
2664	!
2665	AFP(I,J,K)=(((FF4PP+FF3)PP+FF2)PP+FF1)PP
2666	AFQ(I,J,K)=(((FF4QP+FF3)QP+FF2)QP+FF1)QP
2667	!
2668	Q1(I,J,K)=(Q (IFPA(I,J,K),JFPA(I,J,K),K)-Q (I,J,K))*PP &
2669	& +(Q (IFQA(I,J,K),JFQA(I,J,K),K)-Q (I,J,K))*QP &
2670	& +(Q (I,J-2,K)+Q (I,J+2,K) &
2671	& -Q (I-1,J,K)-Q (I+1,J,K))*FPQ &
2672	& +Q(I,J,K)
2673	!
2674	ENDDO
2675	ENDDO
2676	!
2677	!-----------------------------------------------------------------------
2678	!
2679	ENDDO vertical_1
2680	!
2681	!-----------------------------------------------------------------------
2682	!*** ANTI-FILTERING STEP
2683	!-----------------------------------------------------------------------
2684	!
2685	DO K=KTS,KTE
2686	XSUMS(1,K)=0.
2687	XSUMS(2,K)=0.
2688	ENDDO
2689	!
2690	!-----------------------------------------------------------------------
2691	!
2692	!*** ANTI-FILTERING LIMITERS
2693	!
2694	!-----------------------------------------------------------------------
2695	#if defined(BIT_FOR_BIT) && defined(DM_PARALLEL) && !defined(STUBMPI)
2696	DO N=1,2
2697	!
2698	!$omp parallel do &
2699	!$omp& private(i,j,k)
2700	DO K=KMS,KME
2701	DO J=JMS,JME
2702	DO I=IMS,IME
2703	XSUMS_L(I,J,K,N)=0.
2704	ENDDO
2705	ENDDO
2706	ENDDO
2707	!
2708	!$omp parallel do &
2709	!$omp& private(i,j,k)
2710	DO K=KDS,KDE
2711	DO J=JDS,JDE
2712	DO I=IDS,IDE
2713	XSUMS_G(I,J,K,N)=0.
2714	ENDDO
2715	ENDDO
2716	ENDDO
2717	!
2718	ENDDO
2719	!
2720	#endif
2721	!-----------------------------------------------------------------------
2722	!$omp parallel do &
2723	!$omp& private(d2pqe,d2pqq,d2pqw,destij,dqstij,dvolp,dwstij &
2724	!$omp& ,e00,e0q,ep0,estij,hafp,hafq,i,j,k &
2725	!$omp& ,q00,q0q,q1ij,qp0,qstij,w00,w0q,wp0,wstij)
2726	!-----------------------------------------------------------------------
2727	!
2728	vertical_2: DO K=KTS,KTE
2729	!
2730	!-----------------------------------------------------------------------
2731	!
2732	DO J=MYJS2,MYJE2
2733	DO I=MYIS1,MYIE1
2734	!
2735	DVOLP=DVOL(I,J,K)
2736	Q1IJ =Q1(I,J,K)
2737	!
2738	HAFP=AFP(I,J,K)
2739	HAFQ=AFQ(I,J,K)
2740	!
2741	D2PQQ=(Q1(IFPA(I,J,K),JFPA(I,J,K),K)-Q1IJ &
2742	& -Q1IJ+Q1(IFPF(I,J,K),JFPF(I,J,K),K)) &
2743	& *HAFP &
2744	& +(Q1(IFQA(I,J,K),JFQA(I,J,K),K)-Q1IJ &
2745	& -Q1IJ+Q1(IFQF(I,J,K),JFQF(I,J,K),K)) &
2746	& *HAFQ
2747	!
2748	QSTIJ=Q1IJ-D2PQQ
2749	!
2750	Q00=Q (I ,J ,K)
2751	QP0=Q (IFPA(I,J,K),JFPA(I,J,K),K)
2752	Q0Q=Q (IFQA(I,J,K),JFQA(I,J,K),K)
2753	!
2754	QSTIJ=MAX(QSTIJ,MIN(Q00,QP0,Q0Q))
2755	QSTIJ=MIN(QSTIJ,MAX(Q00,QP0,Q0Q))
2756	!
2757	DQSTIJ=QSTIJ-Q(I,J,K)
2758	!
2759	DQST(I,J,K)=DQSTIJ
2760	!
2761	DQSTIJ=DQSTIJ*DVOLP
2762	!
2763	!-----------------------------------------------------------------------
2764	#if defined(BIT_FOR_BIT) && defined(DM_PARALLEL) && !defined(STUBMPI)
2765	!-----------------------------------------------------------------------
2766	DO N=1,2
2767	XSUMS_L(I,J,K,N)=0.
2768	ENDDO
2769	!
2770	IF(DQSTIJ>0.)THEN
2771	XSUMS_L(I,J,K,1)=DQSTIJ
2772	ELSE
2773	XSUMS_L(I,J,K,2)=DQSTIJ
2774	ENDIF
2775	!
2776	!-----------------------------------------------------------------------
2777	#else
2778	!-----------------------------------------------------------------------
2779	IF(DQSTIJ>0.)THEN
2780	XSUMS(1,K)=XSUMS(1,K)+DQSTIJ
2781	ELSE
2782	XSUMS(2,K)=XSUMS(2,K)+DQSTIJ
2783	ENDIF
2784	!
2785	!-----------------------------------------------------------------------
2786	#endif
2787	!-----------------------------------------------------------------------
2788	!
2789	ENDDO
2790	ENDDO
2791	!
2792	!-----------------------------------------------------------------------
2793	!
2794	ENDDO vertical_2
2795	!
2796	!-----------------------------------------------------------------------
2797	#if defined(BIT_FOR_BIT) && defined(DM_PARALLEL) && !defined(STUBMPI)
2798	!-----------------------------------------------------------------------
2799	DO N=1,2
2800	CALL WRF_PATCH_TO_GLOBAL_REAL(XSUMS_L(IMS,JMS,KMS,N) &
2801	&, XSUMS_G(1,1,1,N),DOMDESC &
2802	&, 'xyz','xzy' &
2803	&, IDS,IDE,KDS,KDE,JDS,JDE &
2804	&, IMS,IME,KMS,KME,JMS,JME &
2805	&, ITS,ITE,KTS,KTE,JTS,JTE)
2806	ENDDO
2807	!
2808	DO K=KTS,KTE
2809	DO N=1,2
2810	GSUMS(N,K)=0.
2811	ENDDO
2812	ENDDO
2813	!
2814	IF(WRF_DM_ON_MONITOR())THEN
2815	DO N=1,2
2816	!$omp parallel do &
2817	!$omp& private(i,j,k)
2818	DO K=KDS,KDE
2819	DO J=JDS,JDE
2820	DO I=IDS,IDE
2821	GSUMS(N,K)=GSUMS(N,K)+XSUMS_G(I,J,K,N)
2822	ENDDO
2823	ENDDO
2824	ENDDO
2825	ENDDO
2826	ENDIF
2827
2828	CALL WRF_DM_BCAST_BYTES(GSUMS,2RWORDSIZE2*(KDE-KDS+1) )
2829
2830	!-----------------------------------------------------------------------
2831	#else
2832	!-----------------------------------------------------------------------
2833	!
2834	!-----------------------------------------------------------------------
2835	!*** GLOBAL REDUCTION
2836	!-----------------------------------------------------------------------
2837	!
2838	# if defined(DM_PARALLEL) && !defined(STUBMPI)
2839	CALL WRF_GET_DM_COMMUNICATOR(MPI_COMM_COMP)
2840	CALL MPI_ALLREDUCE(XSUMS,GSUMS,2*(KTE-KTS+1) &
2841	& ,MPI_DOUBLE_PRECISION,MPI_SUM &
2842	& ,MPI_COMM_COMP,IRECV)
2843	# else
2844	DO K=KTS,KTE
2845	DO N=1,2
2846	GSUMS(N,K)=XSUMS(N,K)
2847	ENDDO
2848	ENDDO
2849	# endif
2850	!
2851	!-----------------------------------------------------------------------
2852	#endif
2853	!-----------------------------------------------------------------------
2854	!
2855	!-----------------------------------------------------------------------
2856	!*** END OF GLOBAL REDUCTION
2857	!-----------------------------------------------------------------------
2858	!
2859	! if(mype==0)then
2860	!!! if(ntsd==0)then
2861	!!! call int_get_fresh_handle(nunit)
2862	!!! close(nunit)
2863	! nunit=56
2864	!!! open(unit=nunit,file='gsums',form='unformatted',iostat=ier)
2865	!!! endif
2866	! endif
2867	!-----------------------------------------------------------------------
2868	!$omp parallel do &
2869	!$omp& private(destij,dqstij,dwstij,i,j,k,rface,rfacq,rfacw &
2870	!$omp& ,rfeij,rfqij,rfwij,sumnq,sumpq)
2871	!-----------------------------------------------------------------------
2872	!
2873	vertical_3: DO K=KTS,KTE
2874	!
2875	!-----------------------------------------------------------------------
2876	! if(mype==0)then
2877	! write(nunit)(gsums(i,k),i=1,6)
2878	! endif
2879	!!! read(nunit)(gsums(i,k),i=1,6)
2880	!-----------------------------------------------------------------------
2881	!
2882	SUMPQ=GSUMS(1,K)
2883	SUMNQ=GSUMS(2,K)
2884	!
2885	!-----------------------------------------------------------------------
2886	!*** FIRST MOMENT CONSERVING FACTOR
2887	!-----------------------------------------------------------------------
2888	!
2889	IF(SUMPQ>1.)THEN
2890	RFACQ=-SUMNQ/SUMPQ
2891	ELSE
2892	RFACQ=1.
2893	ENDIF
2894	!
2895	IF(RFACQ<CONSERVE_MIN.OR.RFACQ>CONSERVE_MAX)RFACQ=1.
2896	!
2897	!-----------------------------------------------------------------------
2898	! if(mype==0.and.ntsd==181)close(nunit)
2899	!-----------------------------------------------------------------------
2900	!
2901	!-----------------------------------------------------------------------
2902	!*** IMPOSE CONSERVATION ON ANTI-FILTERING
2903	!-----------------------------------------------------------------------
2904	!
2905	DO J=MYJS2,MYJE2
2906	IF(RFACQ<1.)THEN
2907	DO I=MYIS1,MYIE1
2908	DQSTIJ=DQST(I,J,K)
2909	RFQIJ=HBM2(I,J)*(RFACQ-1.)+1.
2910	IF(DQSTIJ>=0.)DQSTIJ=DQSTIJ*RFQIJ
2911	Q(I,J,K)=Q(I,J,K)+DQSTIJ
2912	ENDDO
2913	ELSE
2914	DO I=MYIS1,MYIE1
2915	DQSTIJ=DQST(I,J,K)
2916	RFQIJ=HBM2(I,J)*(RFACQ-1.)+1.
2917	IF(DQSTIJ<0.)DQSTIJ=DQSTIJ/RFQIJ
2918	Q(I,J,K)=Q(I,J,K)+DQSTIJ
2919	ENDDO
2920	ENDIF
2921	ENDDO
2922	!
2923	!-----------------------------------------------------------------------
2924	!
2925	DO J=MYJS,MYJE
2926	DO I=MYIS,MYIE
2927	SCAL(I,J,K,L)=Q(I,J,K)
2928	ENDDO
2929	ENDDO
2930	!
2931	!-----------------------------------------------------------------------
2932	!
2933	ENDDO vertical_3
2934	!
2935	!-----------------------------------------------------------------------
2936	!
2937	ENDDO scalar_loop
2938	!
2939	!-----------------------------------------------------------------------
2940	!
2941	END SUBROUTINE HAD2_SCAL
2942	!
2943	!-----------------------------------------------------------------------
2944	!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2945	subroutine adv2 &
2946	(UPSTRM &
2947	,mype,kss,kse &
2948	,ids,ide,jds,jde,kds,kde &
2949	,ims,ime,jms,jme,kms,kme &
2950	,its,ite,jts,jte,kts,kte &
2951	,N_IUP_H &
2952	,N_IUP_ADH &
2953	,IUP_H,IUP_ADH &
2954	,ENT &
2955	,idtad &
2956	,dt,pdtop &
2957	,ihe,ihw,ive,ivw &
2958	,deta1,deta2 &
2959	,EMT_LOC &
2960	,fad,hbm2,pdsl,pdslo &
2961	,petdt &
2962	,UOLD,VOLD &
2963	,s,sp &
2964	!---temporary arguments-------------------------------------------------
2965	,fne,fse,few,fns,s1,tcs)
2966	!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2967	implicit none
2968	!-----------------------------------------------------------------------
2969	real,parameter:: &
2970	cfc=1.533 & ! adams-bashforth positioning in time
2971	,bfc=1.-cfc & ! adams bashforth positioning in time
2972	,cflc=9.005 & !
2973	,epsq=1.e-20 & ! floor value for specific humidity
2974	,epsq2=0.2 & ! floor value for 2tke
2975	,epscm=2.e-6 & ! a floor value (not used)
2976	,w1=1.0 & ! crank-nicholson uncentering
2977	!,w1=-1.00 & ! crank-nicholson uncentering
2978	,w2=2.-w1 ! crank-nicholson uncentering
2979
2980	logical,intent(in):: &
2981	upstrm
2982
2983	integer,intent(in):: &
2984	idtad & ! time step multiplier
2985	,kse & ! terminal species index
2986	,kss & ! initial species index
2987	,mype & !
2988	,ids,ide,jds,jde,kds,kde &
2989	,ims,ime,jms,jme,kms,kme &
2990	,its,ite,jts,jte,kts,kte
2991
2992	real,intent(in):: &
2993	ent & !
2994	,dt & ! dynamics time step
2995	,pdtop !
2996
2997	real,dimension(kts:kte),intent(in):: &
2998	deta1 & ! delta sigmas
2999	,deta2 ! delta pressures
3000
3001	integer,dimension(jms:jme),intent(in):: &
3002	ihe,ihw,ive,ivw &
3003	,n_iup_adh,n_iup_h
3004
3005	integer,dimension(ims:ime,jms:jme),intent(in):: &
3006	iup_h,iup_adh
3007
3008	real,dimension(2600),intent(in):: & !!!zj see nmm_max_dim in adve !!!zj
3009	emt_loc
3010
3011	real,dimension(ims:ime,jms:jme),intent(in):: &
3012	fad & !
3013	,hbm2 & !
3014	,pdsl & ! sigma range pressure difference
3015	,pdslo ! sigma range pressure difference
3016
3017	real,dimension(ims:ime,jms:jme,kms:kme),intent(in):: &
3018	petdt & ! vertical mass flux
3019	,uold,vold
3020
3021	real,dimension(ims:ime,jms:jme,kms:kme,kss:kse),intent(inout):: &
3022	s ! tracers
3023
3024	real,dimension(ims:ime,jms:jme,kms:kme,kss:kse),intent(inout):: &
3025	sp ! s at previous time level
3026
3027	!---temporary arguments-------------------------------------------------
3028	real,dimension(ims:ime,jms:jme,kms:kme),intent(in):: &
3029	fne & ! mass flux, ne direction
3030	,fse & ! mass flux, se direction
3031	,few & ! mass flux, x direction
3032	,fns ! mass flux, y direction
3033
3034	real,dimension(ims:ime,jms:jme,kms:kme,kss:kse),intent(inout):: &
3035	s1 & ! intermediate value of s
3036	,tcs ! timechange of s
3037
3038	!--local variables------------------------------------------------------
3039	integer:: &
3040	i & !
3041	,j & !
3042	,k & !
3043	,ks !
3044
3045	INTEGER :: IEND,IFP,IFQ,II,IPQ,ISP,ISQ,ISTART &
3046	& ,IUP_ADH_J &
3047	& ,J1,JA,JAK,JEND,JGLOBAL,JJ,JKNT,JP2,JSTART &
3048	& ,KNTI_ADH,KSTART,KSTOP &
3049	& ,N,N_IUPH_J,N_IUPADH_J,N_IUPADV_J
3050
3051	INTEGER :: MY_IS_GLB,MY_IE_GLB,MY_JS_GLB,MY_JE_GLB
3052
3053	real:: &
3054	cf & ! temporary
3055	,cms & ! temporary
3056	,dtq & ! dt/4
3057	,fahp & ! temporary grid factor
3058	,sn & !
3059	,rdp & ! 1/deltap
3060	,vvlo & ! vertical velocity, lower interface
3061	,vvup & ! vertical velocity, upper interface
3062	,pvvup ! vertical mass flux, upper interface
3063
3064	REAL :: ARRAY3_X &
3065	& ,F0,F1,F2,F3 &
3066	& ,PP &
3067	& ,QP &
3068	& ,TEMPA,TEMPB,TTA,TTB
3069
3070	real,dimension(kts:kte):: &
3071	deta1_pdtop !
3072
3073	INTEGER,DIMENSION(ITS-5:ITE+5,JTS-5:JTE+5) :: ISPA,ISQA
3074
3075	real,dimension(its-5:ite+5,jts-5:jte+5):: &
3076	pdop & ! hydrostatic pressure difference at h points
3077	,pvvlo & ! vertical mass flux, lower interface
3078	,ss1 & ! extrapolated species between time levels
3079	,ssne & ! flux, ne direction
3080	,ssse & ! flux, nw direction
3081	,ssx & ! flux, x direction
3082	,ssy ! flux, y direction
3083
3084	REAL,DIMENSION(ITS-5:ITE+5,JTS-5:JTE+5) :: ARRAY0,ARRAY1 &
3085	& ,ARRAY2,ARRAY3
3086
3087	real,dimension(its-5:ite+5,jts-5:jte+5,kts:kte):: &
3088	crs & ! vertical advection temporary
3089	,rcms ! vertical advection temporary
3090
3091	real,dimension(its-5:ite+5,jts-5:jte+5,kts:kte,kss:kse):: &
3092	rsts ! vertical advection temporary
3093	!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
3094	DO J=JTS-5,JTE+5
3095	DO I=ITS-5,ITE+5
3096	pdop (i,j)=0.
3097	pvvlo(i,j)=0.
3098	ss1 (i,j)=0.
3099	ssne (i,j)=0.
3100	ssse (i,j)=0.
3101	enddo
3102	enddo
3103	!
3104	DO K=KTS,KTE
3105	DO J=JTS-5,JTE+5
3106	DO I=ITS-5,ITE+5
3107	crs (i,j,k)=0.
3108	rcms(i,j,k)=0.
3109	enddo
3110	enddo
3111	enddo
3112	!
3113	do ks=kss,kse
3114	DO K=KTS,KTE
3115	DO J=JTS-5,JTE+5
3116	DO I=ITS-5,ITE+5
3117	rsts(i,j,k,ks)=0.
3118	s1 (i,j,k,ks)=0.
3119	enddo
3120	enddo
3121	enddo
3122	enddo
3123	!
3124	do ks=kss,kse
3125	DO K=KMS,KME
3126	DO J=JMS,JME
3127	DO I=IMS,IME
3128	s1 (i,j,k,ks)=0.
3129	tcs(i,j,k,ks)=0.
3130	enddo
3131	enddo
3132	enddo
3133	enddo
3134	!-----------------------------------------------------------------------
3135	do k=kts,kte
3136	deta1_pdtop(k)=deta1(k)*pdtop
3137	enddo
3138	!-----------------------------------------------------------------------
3139	do ks=kss,kse ! loop by species
3140	!-----------------------------------------------------------------------
3141	DO K=KTS,KTE
3142	DO J=MYJS_P4,MYJE_P4
3143	DO I=MYIS_P4,MYIE_P4
3144	s1(i,j,k,ks)=sqrt(s(i,j,k,ks))
3145	enddo
3146	enddo
3147	enddo
3148	!-----------------------------------------------------------------------
3149	enddo ! end of the loop by species
3150	!-----------------------------------------------------------------------
3151	DO J=MYJS_P4,MYJE_P4
3152	DO I=MYIS_P4,MYIE_P4
3153	pdop(i,j)=(pdslo(i,j)+pdsl(i,j))*0.5
3154	enddo
3155	enddo
3156	!---crank-nicholson vertical advection----------------------------------
3157	dtq=dtidtad0.25
3158
3159	DO J=MYJS2,MYJE2
3160	DO I=MYIS1,MYIE1
3161	pvvlo(i,j)=petdt(i,j,kte-1)*dtq
3162	vvlo=pvvlo(i,j)/(deta2(kte)*pdop(i,j)+deta1_pdtop(kte))
3163	!
3164	cms=-vvlo*w2+1.
3165	rcms(i,j,kte)=1./cms
3166	crs(i,j,kte)=vvlo*w2
3167	!
3168	do ks=kss,kse
3169	rsts(i,j,kte,ks)=(-vvlo*w1) &
3170	*(s1(i,j,kte-1,ks)-s1(i,j,kte,ks)) &
3171	+s1(i,j,kte,ks)
3172	enddo
3173	enddo
3174	enddo
3175	DO K=KTE-1,KTS+1,-1
3176	DO J=MYJS2,MYJE2
3177	DO I=MYIS1,MYIE1
3178	rdp=1./(deta2(k)*pdop(i,j)+deta1_pdtop(k))
3179	pvvup=pvvlo(i,j)
3180	pvvlo(i,j)=petdt(i,j,k-1)*dtq
3181	!
3182	vvup=pvvup*rdp
3183	vvlo=pvvlo(i,j)*rdp
3184	!
3185	! if(abs(vvlo).gt.cflc) then
3186	! if(vvlo.lt.0.) then
3187	! vvlo=-cflc
3188	! else
3189	! vvlo= cflc
3190	! endif
3191	! endif
3192	!
3193	cf=-vvupw2rcms(i,j,k+1)
3194	cms=-crs(i,j,k+1)cf+((vvup-vvlo)w2+1.)
3195	rcms(i,j,k)=1./cms
3196	crs(i,j,k)=vvlo*w2
3197	!
3198	do ks=kss,kse
3199	rsts(i,j,k,ks)=-rsts(i,j,k+1,ks)*cf+s1(i,j,k,ks) &
3200	-(s1(i,j,k ,ks)-s1(i,j,k+1,ks))vvupw1 &
3201	-(s1(i,j,k-1,ks)-s1(i,j,k ,ks))vvlow1
3202	enddo
3203	enddo
3204	enddo
3205	enddo
3206	DO J=MYJS2,MYJE2
3207	DO I=MYIS1,MYIE1
3208	pvvup=pvvlo(i,j)
3209	vvup=pvvup/(deta2(kts)*pdop(i,j)+deta1_pdtop(kts))
3210	!
3211	cf=-vvupw2rcms(i,j,kts+1)
3212	cms=-crs(i,j,kts+1)cf+(vvupw2+1.)
3213	rcms(i,j,kts)=1./cms
3214	crs(i,j,kts)=0.
3215	!
3216	do ks=kss,kse
3217	rsts(i,j,kts,ks)=-rsts(i,j,kts+1,ks)*cf+s1(i,j,kts,ks) &
3218	-(s1(i,j,kts,ks)-s1(i,j,kts+1,ks))vvupw1
3219	!
3220	tcs(i,j,kts,ks)=rsts(i,j,kts,ks)*rcms(i,j,kts)-s1(i,j,kts,ks)
3221	enddo
3222	enddo
3223	enddo
3224	do ks=kss,kse
3225	DO K=KTS+1,KTE
3226	DO J=MYJS2,MYJE2
3227	DO I=MYIS1,MYIE1
3228	tcs(i,j,k,ks)=(-crs(i,j,k)*(s1(i,j,k-1,ks)+tcs(i,j,k-1,ks)) &
3229	+rsts(i,j,k,ks)) &
3230	*rcms(i,j,k)-s1(i,j,k,ks)
3231	enddo
3232	enddo
3233	enddo
3234	enddo
3235	!-----------------------------------------------------------------------
3236	do ks=kss,kse ! loop by species
3237	!-----------------------------------------------------------------------
3238	DO K=KTS,KTE
3239	DO J=MYJS_P5,MYJE_P5
3240	DO I=MYIS_P5,MYIE_P5
3241	ss1(i,j)=s1(i,j,k,ks)cfc+sp(i,j,k,ks)bfc
3242	sp(i,j,k,ks)=s1(i,j,k,ks)
3243	enddo
3244	enddo
3245	!---fluxes--------------------------------------------------------------
3246	DO J=MYJS1_P2,MYJE1_P2
3247	DO I=MYIS_P2,MYIE_P3
3248	ssx(i,j)=(ss1(i+ive(j),j )-ss1(i+ivw(j),j ))*few(i,j,k) &
3249	*hbm2(i,j)
3250	ssy(i,j)=(ss1(i ,j+1)-ss1(i ,j-1))*fns(i,j,k) &
3251	*hbm2(i,j)
3252	enddo
3253	enddo
3254	DO J=MYJS1_P2,MYJE2_P2
3255	DO I=MYIS_P2,MYIE_P2
3256	ssne(i,j)=(ss1(i+ihe(j),j+1)-ss1(i,j))fne(i,j,k)hbm2(i,j)
3257	enddo
3258	enddo
3259	DO J=MYJS2_P2,MYJE1_P2
3260	DO I=MYIS_P2,MYIE_P2
3261	ssse(i,j)=(ss1(i+ihe(j),j-1)-ss1(i,j))fse(i,j,k)hbm2(i,j)
3262	enddo
3263	enddo
3264	!---advection of species------------------------------------------------
3265	DO J=MYJS5,MYJE5
3266	DO I=MYIS2,MYIE2
3267	tcs(i,j,k,ks)=((ssx (i+ihw(j),j )+ssx (i+ihe(j),j ) &
3268	+ssy (i ,j-1)+ssy (i ,j+1) &
3269	+ssne(i+ihw(j),j-1)+ssne(i ,j ) &
3270	+ssse(i ,j )+ssse(i+ihw(j),j+1)) &
3271	fad(i,j)2.0*idtad & !! 2.0 compensates for fad
3272	/(deta2(k)*pdop(i,j)+deta1_pdtop(k)) &
3273	+tcs(i,j,k,ks))*hbm2(i,j)
3274	enddo
3275	enddo
3276	!-----------------------------------------------------------------------
3277	!
3278	!*** upstream advection
3279	!
3280	!-----------------------------------------------------------------------
3281	!
3282	upstream: IF(UPSTRM)THEN
3283	!
3284	!-----------------------------------------------------------------------
3285	!***
3286	!*** COMPUTE UPSTREAM COMPUTATIONS ON THIS TASK'S ROWS.
3287	!***
3288	!-----------------------------------------------------------------------
3289	!
3290	jloop_upstream: DO J=MYJS2,MYJE2
3291	!
3292	N_IUPH_J=N_IUP_H(J) ! See explanation in START_DOMAIN_NMM
3293	DO II=0,N_IUPH_J-1
3294	!
3295	I=IUP_H(IMS+II,J)
3296	tta=emt_loc(j) &
3297	*(uold(i ,j-1,k)+uold(i+ihw(j),j ,k) &
3298	+uold(i+ihe(j),j ,k)+uold(i ,j+1,k))
3299	ttb=ent &
3300	*(vold(i ,j-1,k)+vold(i+ihw(j),j ,k) &
3301	+vold(i+ihe(j),j ,k)+vold(i ,j+1,k))
3302	PP=-TTA-TTB
3303	QP= TTA-TTB
3304	!
3305	IF(PP<0.)THEN
3306	ISPA(I,J)=-1
3307	ELSE
3308	ISPA(I,J)= 1
3309	ENDIF
3310	!
3311	IF(QP<0.)THEN
3312	ISQA(I,J)=-1
3313	ELSE
3314	ISQA(I,J)= 1
3315	ENDIF
3316	!
3317	PP=ABS(PP)
3318	QP=ABS(QP)
3319	ARRAY3_X=PP*QP
3320	ARRAY0(I,J)=ARRAY3_X-PP-QP
3321	ARRAY1(I,J)=PP-ARRAY3_X
3322	ARRAY2(I,J)=QP-ARRAY3_X
3323	ARRAY3(I,J)=ARRAY3_X
3324	ENDDO
3325	!
3326	!-----------------------------------------------------------------------
3327	!
3328	N_IUPADH_J=N_IUP_ADH(J)
3329	KNTI_ADH=1
3330	IUP_ADH_J=IUP_ADH(IMS,J)
3331	!
3332	iloop_T: DO II=0,N_IUPH_J-1
3333	!
3334	I=IUP_H(IMS+II,J)
3335	!
3336	ISP=ISPA(I,J)
3337	ISQ=ISQA(I,J)
3338	IFP=(ISP-1)/2
3339	IFQ=(-ISQ-1)/2
3340	IPQ=(ISP-ISQ)/2
3341	!
3342	!-----------------------------------------------------------------------
3343	!
3344	IF(I==IUP_ADH_J)THEN ! Upstream advection T tendencies
3345	!
3346	ISP=ISPA(I,J)
3347	ISQ=ISQA(I,J)
3348	IFP=(ISP-1)/2
3349	IFQ=(-ISQ-1)/2
3350	IPQ=(ISP-ISQ)/2
3351	!
3352	F0=ARRAY0(I,J)
3353	F1=ARRAY1(I,J)
3354	F2=ARRAY2(I,J)
3355	F3=ARRAY3(I,J)
3356	!
3357	tcs(i,j,k,ks)=(f0*s1(i,j,k,ks) &
3358	& +f1*s1(i+ihe(j)+ifp,j+isp,k,ks) &
3359	& +f2*s1(i+ihe(j)+ifq,j+isq,k,ks) &
3360	& +f3s1(i+ipq,j+isp+isq,k,ks))2.0 &
3361	& *idtad &
3362	& +tcs(i,j,k,ks)*hbm2(i,j)
3363	!
3364	!-----------------------------------------------------------------------
3365	!
3366	IF(KNTI_ADH<N_IUPADH_J)THEN
3367	IUP_ADH_J=IUP_ADH(IMS+KNTI_ADH,J)
3368	KNTI_ADH=KNTI_ADH+1
3369	ENDIF
3370	!
3371	ENDIF ! End of upstream advection tendency IF block
3372	!
3373	ENDDO iloop_T
3374	!
3375	!-----------------------------------------------------------------------
3376	enddo jloop_upstream
3377	!-----------------------------------------------------------------------
3378	endif upstream
3379	!-----------------------------------------------------------------------
3380	enddo ! kts,kte
3381	!-----------------------------------------------------------------------
3382	enddo ! end of the loop by the species
3383	!-----------------------------------------------------------------------
3384	endsubroutine adv2
3385	!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
3386	subroutine mono &
3387	( &
3388	#if defined(DM_PARALLEL)
3389	domdesc, &
3390	#endif
3391	mype,ntsd,hours,kss,kse &
3392	,ids,ide,jds,jde,kds,kde &
3393	,ims,ime,jms,jme,kms,kme &
3394	,its,ite,jts,jte,kts,kte &
3395	,idtad &
3396	,dy,pdtop &
3397	,sumdrrw &
3398	,ihe,ihw &
3399	,deta1,deta2 &
3400	,dx,hbm2,pd &
3401	,s &
3402	!---temporary arguments-------------------------------------------------
3403	,s1,tcs)
3404	!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
3405	implicit none
3406	!-----------------------------------------------------------------------
3407	real,parameter:: &
3408	epsq=1.e-20 & ! floor value for specific humidity
3409	,epsq2=0.02 ! floor value for 2tke
3410
3411	#ifdef DM_PARALLEL
3412	INTEGER :: DOMDESC
3413	#endif
3414
3415	integer,intent(in):: &
3416	idtad & ! time step multiplier
3417	,kse & ! terminal species index
3418	,kss & ! initial species index
3419	,mype & !
3420	,ntsd & !
3421	,ids,ide,jds,jde,kds,kde &
3422	,ims,ime,jms,jme,kms,kme &
3423	,its,ite,jts,jte,kts,kte
3424
3425	real,intent(in):: &
3426	dy & !
3427	,pdtop & !
3428	,hours !
3429
3430	real,intent(inout):: &
3431	sumdrrw
3432
3433	integer,dimension(jms:jme),intent(in):: &
3434	ihe,ihw
3435
3436	real,dimension(kts:kte),intent(in):: &
3437	deta1 & ! delta sigmas
3438	,deta2 ! delta sigmas
3439
3440	real,dimension(ims:ime,jms:jme),intent(in):: &
3441	dx & !
3442	,hbm2 & !
3443	,pd ! sigma range pressure difference
3444
3445	real,dimension(ims:ime,jms:jme,kms:kme,kss:kse),intent(in):: &
3446	s ! s
3447	!---temporary arguments-------------------------------------------------
3448	real,dimension(ims:ime,jms:jme,kms:kme,kss:kse),intent(inout):: &
3449	s1 & ! intermediate value of s
3450	,tcs ! timechange of s
3451	!--local variables------------------------------------------------------
3452	integer:: &
3453	i & !
3454	,ierr & !
3455	,irecv & !
3456	,j & !
3457	,k & !
3458	,ks & !
3459	,lngth & !
3460	,mpi_comm_comp !
3461
3462	real:: &
3463	dsks & !
3464	,dvolp & !
3465	,rfacs & !
3466	,s1p & !
3467	,sfacs & !
3468	,smax & ! local maximum
3469	,smin & ! local minimum
3470	,smaxh & ! horizontal local maximum
3471	,sminh & ! horizontal local minimum
3472	,smaxv & ! vertical local maximum
3473	,sminv & ! vertical local minimum
3474	,sn & !
3475	,sumns & !
3476	,sumps & !
3477	,steep
3478
3479	double precision:: &
3480	xsmp,gsmp
3481
3482	double precision,dimension(2kss-1:2kse):: &
3483	vgsms
3484
3485	real,dimension(kts:kte):: &
3486	deta1_pdtop !
3487
3488	real,dimension(2kss-1:2kse,kts:kte):: &
3489	gsms_single !
3490
3491	double precision,dimension(2kss-1:2kse,kts:kte):: &
3492	gsms & ! sum of neg/pos changes all global fields
3493	,xsms ! sum of neg/pos changes all local fields
3494
3495	double precision,dimension(2kss-1:2kse):: &
3496	vgsums !
3497
3498	real,dimension(its-5:ite+5,jts-5:jte+5,kts:kte):: &
3499	dvol & ! grid box volume
3500	,rdvol ! 1./grid box volume
3501
3502	logical, save :: first=.true.
3503	real, save :: sumrrw_first, summass_first
3504	real :: summass
3505	double precision, save :: gsmp_first
3506	!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
3507
3508	! steep=1.-0.040*idtad
3509	steep=1.
3510
3511	DO K=KTS,KTE
3512	deta1_pdtop(k)=deta1(k)*pdtop
3513	enddo
3514	!
3515	DO K=KTS,KTE
3516	DO J=MYJS2,MYJE2
3517	DO I=MYIS1,MYIE1
3518	dvolp=(deta2(k)pd(i,j)+deta1_pdtop(k))dx(i,j)*dy
3519	rdvol(i,j,k)=hbm2(i,j)/dvolp
3520	dvol (i,j,k)=hbm2(i,j)*dvolp
3521	enddo
3522	enddo
3523	enddo
3524
3525	! go to 109
3526	!---monotonization------------------------------------------------------
3527	do ks=kss,kse ! loop by species
3528	!-----------------------------------------------------------------------
3529	DO K=KTS,KTE
3530	xsms(2*ks-1,k)=0.
3531	xsms(2*ks ,k)=0.
3532	gsms(2*ks-1,k)=0.
3533	gsms(2*ks ,k)=0.
3534	!
3535	DO J=MYJS2,MYJE2
3536	DO I=MYIS1,MYIE1
3537	!
3538	s1p=(s1(i,j,k,ks)+tcs(i,j,k,ks))**2
3539	tcs(i,j,k,ks)=s1p-s(i,j,k,ks)
3540	!
3541	sminh=min(s(i ,j-2,k,ks) &
3542	,s(i+ihw(j),j-1,k,ks) &
3543	,s(i+ihe(j),j-1,k,ks) &
3544	,s(i-1 ,j ,k,ks) &
3545	,s(i ,j ,k,ks) &
3546	,s(i+1 ,j ,k,ks) &
3547	,s(i+ihw(j),j+1,k,ks) &
3548	,s(i+ihe(j),j+1,k,ks) &
3549	,s(i ,j+2,k,ks))
3550	smaxh=max(s(i ,j-2,k,ks) &
3551	,s(i+ihw(j),j-1,k,ks) &
3552	,s(i+ihe(j),j-1,k,ks) &
3553	,s(i-1 ,j ,k,ks) &
3554	,s(i ,j ,k,ks) &
3555	,s(i+1 ,j ,k,ks) &
3556	,s(i+ihw(j),j+1,k,ks) &
3557	,s(i+ihe(j),j+1,k,ks) &
3558	,s(i ,j+2,k,ks))
3559	!
3560	if(k.gt.kts.and.k.lt.kte) then
3561	sminv=min(s(i,j,k-1,ks),s(i,j,k ,ks),s(i,j,k+1,ks))
3562	smaxv=max(s(i,j,k-1,ks),s(i,j,k ,ks),s(i,j,k+1,ks))
3563	elseif(k.eq.kts) then
3564	sminv=min(s(i,j,k ,ks),s(i,j,k+1,ks))
3565	smaxv=max(s(i,j,k ,ks),s(i,j,k+1,ks))
3566	elseif(k.eq.kte) then
3567	sminv=min(s(i,j,k-1,ks),s(i,j,k ,ks))
3568	smaxv=max(s(i,j,k-1,ks),s(i,j,k ,ks))
3569	endif
3570	!
3571	smin=min(sminh,sminv)
3572	smax=max(smaxh,smaxv)
3573	!
3574	sn=s1p
3575	if(sn.gt.steep*smax) sn=smax
3576	if(sn.lt. smin) sn=smin
3577	!
3578	dsks=(sn-s1p)*dvol(i,j,k)
3579	s1(i,j,k,ks)=dsks
3580	if(dsks.gt.0.) then
3581	xsms(2ks-1,k)=xsms(2ks-1,k)+dsks
3582	else
3583	xsms(2ks ,k)=xsms(2ks ,k)+dsks
3584	endif
3585	!
3586	enddo
3587	enddo
3588	enddo
3589	!-----------------------------------------------------------------------
3590	enddo ! end of the loop by species
3591	!-----------------------------------------------------------------------
3592	!-----------------------------------------------------------------------
3593	!*** GLOBAL REDUCTION
3594	!-----------------------------------------------------------------------
3595	!
3596	# if defined(DM_PARALLEL) && !defined(STUBMPI)
3597	CALL WRF_GET_DM_COMMUNICATOR(MPI_COMM_COMP)
3598	lngth=(2kse-2kss+2)*(KTE-KTS+1)
3599	CALL MPI_ALLREDUCE(XSMS,GSMS,lngth &
3600	& ,MPI_DOUBLE_PRECISION,MPI_SUM &
3601	& ,MPI_COMM_COMP,IRECV)
3602	# else
3603	DO K=KTS,KTE
3604	do ks=kss,kse
3605	gsms(2ks-1,k)=xsms(2ks-1,k)
3606	gsms(2ks ,k)=xsms(2ks ,k)
3607	enddo
3608	enddo
3609	# endif
3610	!
3611	DO K=KTS,KTE
3612	do ks=kss,kse
3613	gsms_single(2ks-1,k)=gsms(2ks-1,k)
3614	gsms_single(2ks ,k)=gsms(2ks ,k)
3615	enddo
3616	enddo
3617	!
3618	!-----------------------------------------------------------------------
3619	!*** END OF GLOBAL REDUCTION
3620	!-----------------------------------------------------------------------
3621	!
3622	!dusan!---forced conservation after monotonization----------------------------
3623	!dusan do ks=kss,kse
3624	!dusan DO K=KTS,KTE
3625	!dusan sumps=gsms_single(2*ks-1,k)
3626	!dusan sumns=gsms_single(2*ks ,k)
3627	!dusan!
3628	!dusan if(sumps*(-sumns).gt.1.) then
3629	!dusan sfacs=-sumns/sumps
3630	!dusan rfacs=1./sfacs
3631	!dusan else
3632	!dusan sfacs=0.
3633	!dusan rfacs=0.
3634	!dusan endif
3635	!dusan!
3636	!dusan DO J=MYJS2,MYJE2
3637	!dusan DO I=MYIS1,MYIE1
3638	!dusan dsks=s1(i,j,k,ks)*rdvol(i,j,k)
3639	!dusan if(sfacs.gt.1.) then
3640	!dusan if(dsks.lt.0.) dsks=dsks*rfacs
3641	!dusan!zjtest if(dsks.lt.0.) dsks=dsks
3642	!dusan else
3643	!dusan if(dsks.ge.0.) dsks=dsks*sfacs
3644	!dusan endif
3645	!dusan tcs(i,j,k,ks)=tcs(i,j,k,ks)+dsks
3646	!dusan enddo
3647	!dusan enddo
3648	!dusan enddo
3649	!dusan!-----------------------------------------------------------------------
3650	!dusan enddo ! end of the loop by species
3651
3652	!---forced conservation after monotonization----------------------------
3653	do ks=kss,kse
3654	vgsums(2*ks-1)=0.
3655	vgsums(2*ks )=0.
3656	DO K=KTS,KTE
3657	vgsums(2ks-1)=gsms(2ks-1,k)+vgsums(2*ks-1)
3658	vgsums(2ks )=gsms(2ks ,k)+vgsums(2*ks )
3659	enddo
3660	enddo
3661
3662	do ks=kss,kse
3663	sumps=vgsums(2*ks-1)
3664	sumns=vgsums(2*ks )
3665	!
3666	if(sumps*(-sumns).gt.1.) then
3667	sfacs=-sumns/sumps
3668	rfacs=1./sfacs
3669	else
3670	sfacs=0.
3671	rfacs=0.
3672	endif
3673	!
3674	DO K=KTS,KTE
3675	DO J=MYJS2,MYJE2
3676	DO I=MYIS1,MYIE1
3677	dsks=s1(i,j,k,ks)*rdvol(i,j,k)
3678	if(sfacs.lt.1.) then
3679	if(dsks.gt.0.) dsks=dsks*sfacs
3680	endif
3681	tcs(i,j,k,ks)=tcs(i,j,k,ks)+dsks
3682	enddo
3683	enddo
3684	enddo
3685	!-----------------------------------------------------------------------
3686	enddo ! end of the loop by species
3687	109 continue
3688	!-----------------------------------------------------------------------
3689	!zjwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwwww
3690	!-----------------------------------------------------------------------
3691	do ks=kss,kse ! loop by species
3692	DO K=KTS,KTE
3693	xsms(2*ks-1,k)=0.
3694	xsms(2*ks ,k)=0.
3695	gsms(2*ks-1,k)=0.
3696	gsms(2*ks ,k)=0.
3697	!
3698	DO J=MYJS2,MYJE2
3699	DO I=MYIS1,MYIE1
3700	xsms(2ks-1,k)=xsms(2ks-1,k)+s (i,j,k,ks)*dvol(i,j,k)
3701	xsms(2ks ,k)=xsms(2ks ,k)+tcs(i,j,k,ks)*dvol(i,j,k)
3702	enddo
3703	enddo
3704	enddo
3705	enddo
3706	!
3707	xsmp=0.
3708	DO J=MYJS2,MYJE2
3709	DO I=MYIS1,MYIE1
3710	xsmp=pd(i,j)dx(i,j)dy*hbm2(i,j)+xsmp
3711	enddo
3712	enddo
3713	!
3714	!-----------------------------------------------------------------------
3715	!*** GLOBAL REDUCTION
3716	!-----------------------------------------------------------------------
3717	!
3718	# if defined(DM_PARALLEL) && !defined(STUBMPI)
3719	CALL WRF_GET_DM_COMMUNICATOR(MPI_COMM_COMP)
3720	lngth=1
3721	CALL MPI_ALLREDUCE(xsmp,gsmp,lngth &
3722	& ,MPI_DOUBLE_PRECISION,MPI_SUM &
3723	& ,MPI_COMM_COMP,IRECV)
3724	# else
3725	gsmp=xsmp
3726	# endif
3727	!
3728	!-----------------------------------------------------------------------
3729	!*** END OF GLOBAL REDUCTION
3730	!-----------------------------------------------------------------------
3731	!
3732	!
3733	!-----------------------------------------------------------------------
3734	!*** GLOBAL REDUCTION
3735	!-----------------------------------------------------------------------
3736	!
3737	# if defined(DM_PARALLEL) && !defined(STUBMPI)
3738	CALL WRF_GET_DM_COMMUNICATOR(MPI_COMM_COMP)
3739	lngth=(2kse-2kss+2)*(KTE-KTS+1)
3740	CALL MPI_ALLREDUCE(XSMS,GSMS,lngth &
3741	& ,MPI_DOUBLE_PRECISION,MPI_SUM &
3742	& ,MPI_COMM_COMP,IRECV)
3743	# else
3744	DO K=KTS,KTE
3745	do ks=kss,kse
3746	gsms(2ks-1,k)=xsms(2ks-1,k)
3747	gsms(2ks ,k)=xsms(2ks ,k)
3748	enddo
3749	enddo
3750	# endif
3751	!
3752	DO K=KTS,KTE
3753	do ks=kss,kse
3754	gsms_single(2ks-1,k)=gsms(2ks-1,k)
3755	gsms_single(2ks ,k)=gsms(2ks ,k)
3756	enddo
3757	enddo
3758	!
3759	!-----------------------------------------------------------------------
3760	!*** END OF GLOBAL REDUCTION
3761	!-----------------------------------------------------------------------
3762	!
3763
3764	do ks=kss,kse
3765	vgsms(2*ks-1)=0.
3766	vgsms(2*ks )=0.
3767	enddo
3768	!
3769	do ks=kss,kse
3770	DO K=KTS,KTE
3771	vgsms(2ks-1)=gsms(2ks-1,k)+vgsms(2*ks-1)
3772	vgsms(2ks )=gsms(2ks ,k)+vgsms(2*ks )
3773	enddo
3774	enddo
3775	!
3776	sumdrrw=vgsms(6)+sumdrrw
3777	!
3778	summass=0.0
3779	DO K=KTS,KTE
3780	DO J=MYJS2,MYJE2
3781	DO I=MYIS1,MYIE1
3782	summass=summass+dvol(i,j,k)
3783	ENDDO
3784	ENDDO
3785	ENDDO
3786	!
3787	if (first) then
3788	sumrrw_first = vgsms(5)
3789	gsmp_first = gsmp
3790	summass_first = summass
3791	first=.false.
3792	end if
3793	!
3794	! write(0,1000) ntsd,hours, & ! 4,5
3795	! (vgsms(ks),ks=2kss-1,2kse), & ! 6-13
3796	! gsmp,gsmp_first,gsmp/gsmp_first, & ! 14-16
3797	! sumdrrw, & ! 17
3798	! vgsms(5)/sumrrw_first,sumrrw_first, & ! 18-19
3799	! summass,summass_first,summass/summass_first ! 20-22
3800	1000 format('global vol sums ',i6,f8.3,30d13.5)
3801	!zjmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmm
3802	endsubroutine mono
3803	!-----------------------------------------------------------------------
3804	!-----------------------------------------------------------------------
3805	!
3806	END MODULE MODULE_ADVECTION
3807	!
3808	!-----------------------------------------------------------------------

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

Download in other formats:

Original Format