Context Navigation

solve_em.F @ 114

Last change on this file since 114 was 6, checked in by lfita, 11 years ago

Fixing problem with the solar radiation:

1.- Correct values for the julian and the hour (0,1) to the LMDZ physiq

Adding NaN checking for 'PSFC' and 'T_2'

File size: 224.6 KB

Line
1	!WRF:MEDIATION_LAYER:SOLVER
2
3	SUBROUTINE solve_em ( grid , config_flags &
4	! Arguments generated from Registry
5	#include "dummy_new_args.inc"
6	!
7	)
8	! Driver layer modules
9	USE module_state_description
10	USE module_domain, ONLY : &
11	domain, get_ijk_from_grid, get_ijk_from_subgrid &
12	,domain_get_current_time, domain_get_start_time &
13	,domain_get_sim_start_time, domain_clock_get
14	USE module_domain_type, ONLY : history_alarm
15	USE module_configure, ONLY : grid_config_rec_type
16	USE module_driver_constants
17	USE module_machine
18	USE module_tiles, ONLY : set_tiles
19	#ifdef DM_PARALLEL
20	USE module_dm, ONLY : &
21	local_communicator, mytask, ntasks, ntasks_x, ntasks_y &
22	,local_communicator_periodic, wrf_dm_maxval
23	USE module_comm_dm, ONLY : &
24	halo_em_a_sub,halo_em_b_sub,halo_em_c2_sub,halo_em_chem_e_3_sub &
25	,halo_em_chem_e_5_sub,halo_em_chem_e_7_sub,halo_em_chem_old_e_5_sub &
26	,halo_em_chem_old_e_7_sub,halo_em_c_sub,halo_em_d2_3_sub &
27	,halo_em_d2_5_sub,halo_em_d3_3_sub,halo_em_d3_5_sub,halo_em_d_sub &
28	,halo_em_e_3_sub,halo_em_e_5_sub,halo_em_hydro_uv_sub &
29	,halo_em_moist_e_3_sub,halo_em_moist_e_5_sub,halo_em_moist_e_7_sub &
30	,halo_em_moist_old_e_5_sub,halo_em_moist_old_e_7_sub &
31	,halo_em_scalar_e_3_sub,halo_em_scalar_e_5_sub,halo_em_scalar_e_7_sub &
32	,halo_em_scalar_old_e_5_sub,halo_em_scalar_old_e_7_sub,halo_em_tke_3_sub &
33	,halo_em_tke_5_sub,halo_em_tke_7_sub,halo_em_tke_advect_3_sub &
34	,halo_em_tke_advect_5_sub,halo_em_tke_old_e_5_sub &
35	,halo_em_tke_old_e_7_sub,halo_em_tracer_e_3_sub,halo_em_tracer_e_5_sub &
36	,halo_em_tracer_e_7_sub,halo_em_tracer_old_e_5_sub &
37	,halo_em_tracer_old_e_7_sub,period_bdy_em_a_sub &
38	,period_bdy_em_b3_sub,period_bdy_em_b_sub,period_bdy_em_chem2_sub &
39	,period_bdy_em_chem_old_sub,period_bdy_em_chem_sub,period_bdy_em_d3_sub &
40	,period_bdy_em_d_sub,period_bdy_em_e_sub,period_bdy_em_moist2_sub &
41	,period_bdy_em_moist_old_sub,period_bdy_em_moist_sub &
42	,period_bdy_em_scalar2_sub,period_bdy_em_scalar_old_sub &
43	,period_bdy_em_scalar_sub,period_bdy_em_tke_old_sub &
44	,period_bdy_em_tracer2_sub,period_bdy_em_tracer_old_sub &
45	,period_bdy_em_tracer_sub,period_em_da_sub,period_em_hydro_uv_sub
46	#endif
47	USE module_utility
48	! Mediation layer modules
49	! Model layer modules
50	USE module_model_constants
51	USE module_small_step_em
52	USE module_em
53	USE module_big_step_utilities_em
54	USE module_bc
55	USE module_bc_em
56	USE module_solvedebug_em
57	USE module_physics_addtendc
58	USE module_diffusion_em
59	USE module_polarfft
60	USE module_microphysics_driver
61	USE module_microphysics_zero_out
62	USE module_fddaobs_driver
63	USE module_diagnostics
64	#ifdef WRF_CHEM
65	USE module_input_chem_data
66	USE module_input_tracer
67	USE module_chem_utilities
68	#endif
69	USE module_first_rk_step_part1
70	USE module_first_rk_step_part2
71	USE module_llxy, ONLY : proj_cassini
72	USE module_avgflx_em, ONLY : zero_avgflx, upd_avgflx
73	#ifdef LMDZ
74	USE module_lmdz_phys
75	USE module_domain_type, ONLY: RESTART_ALARM
76	USE module_domain, ONLY: domain_get_time_since_sim_start
77	USE module_streams
78
79	#endif
80
81	IMPLICIT NONE
82
83	! Input data.
84
85	TYPE(domain) , TARGET :: grid
86
87	! Definitions of dummy arguments to this routine (generated from Registry).
88	#include "dummy_new_decl.inc"
89
90	! Structure that contains run-time configuration (namelist) data for domain
91	TYPE (grid_config_rec_type) , INTENT(IN) :: config_flags
92
93	! Local data
94
95	INTEGER :: k_start , k_end, its, ite, jts, jte
96	INTEGER :: ids , ide , jds , jde , kds , kde , &
97	ims , ime , jms , jme , kms , kme , &
98	ips , ipe , jps , jpe , kps , kpe
99
100	INTEGER :: sids , side , sjds , sjde , skds , skde , &
101	sims , sime , sjms , sjme , skms , skme , &
102	sips , sipe , sjps , sjpe , skps , skpe
103
104
105	INTEGER :: imsx, imex, jmsx, jmex, kmsx, kmex, &
106	ipsx, ipex, jpsx, jpex, kpsx, kpex, &
107	imsy, imey, jmsy, jmey, kmsy, kmey, &
108	ipsy, ipey, jpsy, jpey, kpsy, kpey
109
110	INTEGER :: ij , iteration
111	INTEGER :: im , num_3d_m , ic , num_3d_c , is , num_3d_s
112	INTEGER :: loop
113	INTEGER :: sz
114	INTEGER :: iswater
115
116	LOGICAL :: specified_bdy, channel_bdy
117
118	REAL :: t_new
119
120	! Changes in tendency at this timestep
121	real ,DIMENSION(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33) :: h_tendency, &
122	z_tendency
123
124	! Whether advection should produce decoupled horizontal and vertical advective tendency outputs
125	LOGICAL :: tenddec
126
127	#ifdef WRF_CHEM
128	! Index cross-referencing array for tendency accumulation
129	INTEGER, DIMENSION( num_chem ) :: adv_ct_indices
130	#endif
131
132	! storage for tendencies and decoupled state (generated from Registry)
133
134	#include <i1_decl.inc>
135	! Previous time level of tracer arrays now defined as i1 variables;
136	! the state 4d arrays now redefined as 1-time level arrays in Registry.
137	! Benefit: save memory in nested runs, since only 1 domain is active at a
138	! time. Potential problem on stack-limited architectures: increases
139	! amount of data on program stack by making these automatic arrays.
140
141	INTEGER :: rc
142	INTEGER :: number_of_small_timesteps, rk_step
143	INTEGER :: klevel,ijm,ijp,i,j,k,size1,size2 ! for prints/plots only
144	INTEGER :: idum1, idum2, dynamics_option
145
146	INTEGER :: rk_order, iwmax, jwmax, kwmax
147	REAL :: dt_rk, dts_rk, dts, dtm, wmax
148	REAL , ALLOCATABLE , DIMENSION(:) :: max_vert_cfl_tmp, max_horiz_cfl_tmp
149	LOGICAL :: leapfrog
150	INTEGER :: l,kte,kk
151	LOGICAL :: f_flux ! flag for computing averaged fluxes in cu_gd
152	REAL :: curr_secs
153	INTEGER :: num_sound_steps
154	INTEGER :: idex, jdex
155	REAL :: max_msft
156	REAL :: spacing
157
158	INTEGER :: ii, jj !kk is above after l,kte
159	REAL :: dclat
160	INTEGER :: debug_level
161
162	! urban related variables
163	INTEGER :: NUM_ROOF_LAYERS, NUM_WALL_LAYERS, NUM_ROAD_LAYERS ! urban
164
165	TYPE(WRFU_TimeInterval) :: tmpTimeInterval
166	REAL :: real_time
167	LOGICAL :: adapt_step_flag
168	LOGICAL :: fill_w_flag
169
170	! variables for flux-averaging code 20091223
171	CHARACTER*256 :: message, message2
172	REAL :: old_dt
173	TYPE(WRFU_Time) :: temp_time, CurrTime
174	INTEGER, PARAMETER :: precision = 100
175	INTEGER :: num, den
176	TYPE(WRFU_TimeInterval) :: dtInterval
177
178	#ifdef LMDZ
179	INTEGER :: im2, jm2, km2
180	INTEGER :: ix,iy,iz
181	INTEGER :: days, seconds, Sn, Sd
182	TYPE(WRFU_TimeInterval) :: timeSinceSimStart
183	TYPE(WRFU_Time) :: initime, simtime
184	LOGICAL :: wrftestrst, wrftestin
185	REAL :: minSinceSimStart
186	CHARACTER(LEN=256) :: mminlu, mminsl
187	CHARACTER(LEN=50) :: errmsg
188	INTEGER :: hr, minute, sec, ms, julyr, &
189	julday
190	REAL :: gmt
191
192
193	errmsg = 'ERROR -- error -- ERROR -- error'
194
195	im2 = config_flags%i_check_point
196	jm2 = config_flags%j_check_point
197	km2 = config_flags%k_check_point
198	#endif
199
200	! Define benchmarking timers if -DBENCH is compiled
201	#include <bench_solve_em_def.h>
202
203	!----------------------
204	! Executable statements
205	!----------------------
206
207	!<DESCRIPTION>
208	!<pre>
209	! solve_em is the main driver for advancing a grid a single timestep.
210	! It is a mediation-layer routine -> DM and SM calls are made where
211	! needed for parallel processing.
212	!
213	! solve_em can integrate the equations using 3 time-integration methods
214	!
215	! - 3rd order Runge-Kutta time integration (recommended)
216	!
217	! - 2nd order Runge-Kutta time integration
218	!
219	! The main sections of solve_em are
220	!
221	! (1) Runge-Kutta (RK) loop
222	!
223	! (2) Non-timesplit physics (i.e., tendencies computed for updating
224	! model state variables during the first RK sub-step (loop)
225	!
226	! (3) Small (acoustic, sound) timestep loop - within the RK sub-steps
227	!
228	! (4) scalar advance for moist and chem scalar variables (and TKE)
229	! within the RK sub-steps.
230	!
231	! (5) time-split physics (after the RK step), currently this includes
232	! only microphyics
233	!
234	! A more detailed description of these sections follows.
235	!</pre>
236	!</DESCRIPTION>
237
238	! Initialize timers if compiled with -DBENCH
239	#include <bench_solve_em_init.h>
240
241	! set runge-kutta solver (2nd or 3rd order)
242
243	dynamics_option = config_flags%rk_ord
244
245	! Obtain dimension information stored in the grid data structure.
246
247	CALL get_ijk_from_grid ( grid , &
248	ids, ide, jds, jde, kds, kde, &
249	ims, ime, jms, jme, kms, kme, &
250	ips, ipe, jps, jpe, kps, kpe, &
251	imsx, imex, jmsx, jmex, kmsx, kmex, &
252	ipsx, ipex, jpsx, jpex, kpsx, kpex, &
253	imsy, imey, jmsy, jmey, kmsy, kmey, &
254	ipsy, ipey, jpsy, jpey, kpsy, kpey )
255
256	CALL get_ijk_from_subgrid ( grid , &
257	sids, side, sjds, sjde, skds, skde, &
258	sims, sime, sjms, sjme, skms, skme, &
259	sips, sipe, sjps, sjpe, skps, skpe )
260	k_start = kps
261	k_end = kpe
262
263	num_3d_m = num_moist
264	num_3d_c = num_chem
265	num_3d_s = num_scalar
266
267	f_flux = config_flags%do_avgflx_cugd .EQ. 1
268
269	! Compute these starting and stopping locations for each tile and number of tiles.
270	! See: http://www.mmm.ucar.edu/wrf/WG2/topics/settiles
271	CALL set_tiles ( grid , ids , ide , jds , jde , ips , ipe , jps , jpe )
272
273	! Max values of CFL for adaptive time step scheme
274
275	ALLOCATE (max_vert_cfl_tmp(grid%num_tiles))
276	ALLOCATE (max_horiz_cfl_tmp(grid%num_tiles))
277
278	!
279	! Calculate current time in seconds since beginning of model run.
280	! Unfortunately, ESMF does not seem to have a way to return
281	! floating point seconds based on a TimeInterval. So, we will
282	! calculate it here--but, this is not clean!!
283	!
284	tmpTimeInterval = domain_get_current_time ( grid ) - domain_get_sim_start_time ( grid )
285	curr_secs = real_time(tmpTimeInterval)
286
287	old_dt = grid%dt ! store old time step for flux averaging code at end of RK loop
288	!-----------------------------------------------------------------------------
289	! Adaptive time step: Added by T. Hutchinson, WSI 3/5/07
290	! In this call, we do the time-step adaptation and set time-dependent lateral
291	! boundary condition nudging weights.
292	!
293	IF ( (config_flags%use_adaptive_time_step) .and. &
294	( (.not. grid%nested) .or. &
295	( (grid%nested) .and. (abs(grid%dtbc) < 0.0001) ) ) )THEN
296	CALL adapt_timestep(grid, config_flags)
297	adapt_step_flag = .TRUE.
298	ELSE
299	adapt_step_flag = .FALSE.
300	ENDIF
301	! End of adaptive time step modifications
302	!-----------------------------------------------------------------------------
303
304	grid%itimestep = grid%itimestep + 1
305
306	IF (config_flags%polar) dclat = 90./REAL(jde-jds) !(0.5 * 180/ny)
307
308	#ifdef WRF_CHEM
309
310	kte=min(k_end,kde-1)
311	# ifdef DM_PARALLEL
312	if ( num_chem >= PARAM_FIRST_SCALAR ) then
313	!-----------------------------------------------------------------------
314	! see matching halo calls below for stencils
315	!--------------------------------------------------------------
316	CALL wrf_debug ( 200 , ' call HALO_RK_CHEM' )
317	IF ( config_flags%h_mom_adv_order <= 4 ) THEN
318	# include "HALO_EM_CHEM_E_3.inc"
319	IF( config_flags%progn > 0 ) THEN
320	# include "HALO_EM_SCALAR_E_3.inc"
321	ENDIF
322	ELSE IF ( config_flags%h_mom_adv_order <= 6 ) THEN
323	# include "HALO_EM_CHEM_E_5.inc"
324	IF( config_flags%progn > 0 ) THEN
325	# include "HALO_EM_SCALAR_E_5.inc"
326	ENDIF
327	ELSE
328	WRITE(wrf_err_message,*)'solve_em: invalid h_mom_adv_order = ',config_flags%h_mom_adv_order
329	CALL wrf_error_fatal(TRIM(wrf_err_message))
330	ENDIF
331	ENDIF
332	if ( num_tracer >= PARAM_FIRST_SCALAR ) then
333	!-----------------------------------------------------------------------
334	! see matching halo calls below for stencils
335	!--------------------------------------------------------------
336	CALL wrf_debug ( 200 , ' call HALO_RK_tracer' )
337	IF ( config_flags%h_mom_adv_order <= 4 ) THEN
338	# include "HALO_EM_TRACER_E_3.inc"
339	ELSE IF ( config_flags%h_mom_adv_order <= 6 ) THEN
340	# include "HALO_EM_TRACER_E_5.inc"
341	ELSE
342	WRITE(wrf_err_message,*)'solve_em: invalid h_mom_adv_order = ',config_flags%h_mom_adv_order
343	CALL wrf_error_fatal(TRIM(wrf_err_message))
344	ENDIF
345	ENDIF
346	# endif
347	!--------------------------------------------------------------
348	adv_ct_indices( : ) = 1
349	IF ( config_flags%chemdiag == USECHEMDIAG ) THEN
350	! modify tendency list here
351	! note that the referencing direction here is opposite of that in chem_driver
352	adv_ct_indices(p_co ) = p_advh_co
353	adv_ct_indices(p_o3 ) = p_advh_o3
354	adv_ct_indices(p_no ) = p_advh_no
355	adv_ct_indices(p_no2 ) = p_advh_no2
356	adv_ct_indices(p_hno3) = p_advh_hno3
357	adv_ct_indices(p_iso ) = p_advh_iso
358	adv_ct_indices(p_ho ) = p_advh_ho
359	adv_ct_indices(p_ho2 ) = p_advh_ho2
360	END IF
361	#endif
362
363	rk_order = config_flags%rk_ord
364
365	IF ( grid%time_step_sound == 0 ) THEN
366	! This function will give 4 for 6dx and 6 for 10dx and returns even numbers only
367	spacing = min(grid%dx, grid%dy)
368	IF ( ( config_flags%use_adaptive_time_step ) .AND. ( config_flags%map_proj == PROJ_CASSINI ) ) THEN
369	max_msft=MIN ( MAX(grid%max_msftx, grid%max_msfty) , &
370	1.0/COS(config_flags%fft_filter_lat*degrad) )
371	num_sound_steps = max ( 2 * ( INT (300. * grid%dt / (spacing / max_msft) - 0.01 ) + 1 ), 4 )
372	ELSE IF ( config_flags%use_adaptive_time_step ) THEN
373	max_msft= MAX(grid%max_msftx, grid%max_msfty)
374	num_sound_steps = max ( 2 * ( INT (300. * grid%dt / (spacing / max_msft) - 0.01 ) + 1 ), 4 )
375	ELSE
376	num_sound_steps = max ( 2 * ( INT (300. * grid%dt / spacing - 0.01 ) + 1 ), 4 )
377	END IF
378	WRITE(wrf_err_message,*)'grid spacing, dt, time_step_sound=',spacing,grid%dt,num_sound_steps
379	CALL wrf_debug ( 50 , wrf_err_message )
380	ELSE
381	num_sound_steps = grid%time_step_sound
382	ENDIF
383
384	dts = grid%dt/float(num_sound_steps)
385
386	IF (config_flags%use_adaptive_time_step) THEN
387
388	CALL get_wrf_debug_level( debug_level )
389	IF ((config_flags%time_step < 0) .AND. (debug_level.GE.50)) THEN
390	#ifdef DM_PARALLEL
391	CALL wrf_dm_maxval(grid%max_vert_cfl, idex, jdex)
392	#endif
393	WRITE(wrf_err_message,*)'variable dt, max horiz cfl, max vert cfl: ',&
394	grid%dt, grid%max_horiz_cfl, grid%max_vert_cfl
395	CALL wrf_debug ( 0 , wrf_err_message )
396	ENDIF
397
398	grid%max_cfl_val = 0
399	grid%max_horiz_cfl = 0
400	grid%max_vert_cfl = 0
401	ENDIF
402
403	! setting bdy tendencies to zero for DFI if constant_bc = true
404
405	!$OMP PARALLEL DO &
406	!$OMP PRIVATE ( ij )
407	DO ij = 1 , grid%num_tiles
408
409	! IF( config_flags%specified .AND. grid%dfi_opt .NE. DFI_NODFI &
410	! .AND. config_flags%constant_bc .AND. (grid%dfi_stage .EQ. DFI_BCK .OR. grid%dfi_stage .EQ. DFI_FWD) ) THEN
411	IF( config_flags%specified .AND. config_flags%constant_bc ) THEN
412
413	CALL zero_bdytend (grid%u_btxs,grid%u_btxe,grid%u_btys,grid%u_btye, &
414	grid%v_btxs,grid%v_btxe,grid%v_btys,grid%v_btye, &
415	grid%ph_btxs,grid%ph_btxe,grid%ph_btys,grid%ph_btye, &
416	grid%t_btxs,grid%t_btxe,grid%t_btys,grid%t_btye, &
417	grid%w_btxs,grid%w_btxe,grid%w_btys,grid%w_btye, &
418	grid%mu_btxs,grid%mu_btxe,grid%mu_btys,grid%mu_btye, &
419	moist_btxs,moist_btxe, &
420	moist_btys,moist_btye, &
421	grid%spec_bdy_width,num_3d_m, &
422	ids,ide, jds,jde, kds,kde, &
423	ims,ime, jms,jme, kms,kme, &
424	ips,ipe, jps,jpe, kps,kpe, &
425	grid%i_start(ij), grid%i_end(ij), &
426	grid%j_start(ij), grid%j_end(ij), &
427	k_start, k_end )
428
429	ENDIF
430	ENDDO
431	!$OMP END PARALLEL DO
432
433	!**********************************************************************
434	!
435	! LET US BEGIN.......
436	!
437	!<DESCRIPTION>
438	!<pre>
439	! (1) RK integration loop is named the "Runge_Kutta_loop:"
440	!
441	! Predictor-corrector type time integration.
442	! Advection terms are evaluated at time t for the predictor step,
443	! and advection is re-evaluated with the latest predicted value for
444	! each succeeding time corrector step
445	!
446	! 2nd order Runge Kutta (rk_order = 2):
447	! Step 1 is taken to the midpoint predictor, step 2 is the full step.
448	!
449	! 3rd order Runge Kutta (rk_order = 3):
450	! Step 1 is taken to from t to dt/3, step 2 is from t to dt/2,
451	! and step 3 is from t to dt.
452	!
453	! non-timesplit physics are evaluated during first RK step and
454	! these physics tendencies are stored for use in each RK pass.
455	!</pre>
456	!</DESCRIPTION>
457	!**********************************************************************
458
459	Runge_Kutta_loop: DO rk_step = 1, rk_order
460
461	! Set the step size and number of small timesteps for
462	! each part of the timestep
463
464	dtm = grid%dt
465	IF ( rk_order == 1 ) THEN
466
467	write(wrf_err_message,*)' leapfrog removed, error exit for dynamics_option = ',dynamics_option
468	CALL wrf_error_fatal( wrf_err_message )
469
470	ELSE IF ( rk_order == 2 ) THEN ! 2nd order Runge-Kutta timestep
471
472	IF ( rk_step == 1) THEN
473	dt_rk = 0.5*grid%dt
474	dts_rk = dts
475	number_of_small_timesteps = num_sound_steps/2
476	ELSE
477	dt_rk = grid%dt
478	dts_rk = dts
479	number_of_small_timesteps = num_sound_steps
480	ENDIF
481
482	ELSE IF ( rk_order == 3 ) THEN ! third order Runge-Kutta
483
484	IF ( rk_step == 1) THEN
485	dt_rk = grid%dt/3.
486	dts_rk = dt_rk
487	number_of_small_timesteps = 1
488	ELSE IF (rk_step == 2) THEN
489	dt_rk = 0.5*grid%dt
490	dts_rk = dts
491	number_of_small_timesteps = num_sound_steps/2
492	ELSE
493	dt_rk = grid%dt
494	dts_rk = dts
495	number_of_small_timesteps = num_sound_steps
496	ENDIF
497
498	ELSE
499
500	write(wrf_err_message,*)' unknown solver, error exit for dynamics_option = ',dynamics_option
501	CALL wrf_error_fatal( wrf_err_message )
502
503	END IF
504
505	! Ensure that polar meridional velocity is zero
506	IF (config_flags%polar) THEN
507	!$OMP PARALLEL DO &
508	!$OMP PRIVATE ( ij )
509	DO ij = 1 , grid%num_tiles
510	CALL zero_pole ( grid%v_1, &
511	ids, ide, jds, jde, kds, kde, &
512	ims, ime, jms, jme, kms, kme, &
513	grid%i_start(ij), grid%i_end(ij), &
514	grid%j_start(ij), grid%j_end(ij), &
515	k_start, k_end )
516	CALL zero_pole ( grid%v_2, &
517	ids, ide, jds, jde, kds, kde, &
518	ims, ime, jms, jme, kms, kme, &
519	grid%i_start(ij), grid%i_end(ij), &
520	grid%j_start(ij), grid%j_end(ij), &
521	k_start, k_end )
522	END DO
523	!$OMP END PARALLEL DO
524	END IF
525	!
526	! Time level t is in the *_2 variable in the first part
527	! of the step, and in the *_1 variable after the predictor.
528	! the latest predicted values are stored in the *_2 variables.
529	!
530	#ifdef LMDZ1
531	WRITE(message, *)' dyn_em: before rk_step_prep rk_step=', rk_step
532	CALL wrf_debug(200, message)
533	WRITE(message, *)' t_tend: ',t_tendf(im2,km2,jm2), &
534	' u_tend: ', ru_tendf(im2,1,jm2)
535	CALL wrf_debug(200, message)
536	WRITE(message,*)' psfc_tend: ',grid%dpsdt(im2,jm2), &
537	'p sfc: ',p8w(im2,kms,jm2)
538	CALL wrf_debug(200, message)
539	WRITE(message,*)' p 1: ',grid%p(im2,kms,jm2), ' ph 1: ',grid%ph_2(im2,kms,jm2)
540	CALL wrf_debug(200, message)
541	#endif
542
543	CALL wrf_debug ( 200 , ' call rk_step_prep ' )
544
545	BENCH_START(step_prep_tim)
546	!$OMP PARALLEL DO &
547	!$OMP PRIVATE ( ij )
548
549	DO ij = 1 , grid%num_tiles
550
551	CALL rk_step_prep ( config_flags, rk_step, &
552	grid%u_2, grid%v_2, grid%w_2, grid%t_2, grid%ph_2, grid%mu_2, &
553	moist, &
554	grid%ru, grid%rv, grid%rw, grid%ww, grid%php, grid%alt, grid%muu, grid%muv, &
555	grid%mub, grid%mut, grid%phb, grid%pb, grid%p, grid%al, grid%alb, &
556	cqu, cqv, cqw, &
557	grid%msfux, grid%msfuy, grid%msfvx, grid%msfvx_inv, &
558	grid%msfvy, grid%msftx, grid%msfty, &
559	grid%fnm, grid%fnp, grid%dnw, grid%rdx, grid%rdy, &
560	num_3d_m, &
561	ids, ide, jds, jde, kds, kde, &
562	ims, ime, jms, jme, kms, kme, &
563	grid%i_start(ij), grid%i_end(ij), &
564	grid%j_start(ij), grid%j_end(ij), &
565	k_start, k_end )
566
567	END DO
568	!$OMP END PARALLEL DO
569	BENCH_END(step_prep_tim)
570
571	#ifdef LMDZ1
572	WRITE(message, *)' dyn_em: after rk_step_prep'
573	CALL wrf_debug(200, message)
574	WRITE(message, *)' t_tend: ',t_tendf(im2,km2,jm2), &
575	' u_tend: ', ru_tendf(im2,1,jm2)
576	CALL wrf_debug(200, message)
577	WRITE(message,*)' psfc_tend: ',grid%dpsdt(im2,jm2), &
578	'p sfc: ',p8w(im2,kms,jm2)
579	CALL wrf_debug(200, message)
580	WRITE(message,*)' p 1: ',grid%p(im2,kms,jm2), ' ph 1: ',grid%ph_2(im2,kms,jm2)
581	CALL wrf_debug(200, message)
582	#endif
583
584	#ifdef DM_PARALLEL
585	!-----------------------------------------------------------------------
586	! Stencils for patch communications (WCS, 29 June 2001)
587	! Note: the small size of this halo exchange reflects the
588	! fact that we are carrying the uncoupled variables
589	! as state variables in the mass coordinate model, as
590	! opposed to the coupled variables as in the height
591	! coordinate model.
592	!
593	! * * * * *
594	! * * * * * * * * *
595	! * + * * + * * * + * *
596	! * * * * * * * * *
597	! * * * * *
598	!
599	! 3D variables - note staggering! ru(X), rv(Y), ww(Z), php(Z)
600	!
601	! ru x
602	! rv x
603	! ww x
604	! php x
605	! alt x
606	! ph_2 x
607	! phb x
608	!
609	! the following are 2D (xy) variables
610	!
611	! muu x
612	! muv x
613	! mut x
614	!--------------------------------------------------------------
615	# include "HALO_EM_A.inc"
616	#endif
617
618	! set boundary conditions on variables
619	! from big_step_prep for use in big_step_proc
620
621	#ifdef DM_PARALLEL
622	# include "PERIOD_BDY_EM_A.inc"
623	#endif
624
625	BENCH_START(set_phys_bc_tim)
626	!$OMP PARALLEL DO &
627	!$OMP PRIVATE ( ij, ii, jj, kk )
628
629	DO ij = 1 , grid%num_tiles
630
631	CALL wrf_debug ( 200 , ' call rk_phys_bc_dry_1' )
632
633	CALL rk_phys_bc_dry_1( config_flags, grid%ru, grid%rv, grid%rw, grid%ww, &
634	grid%muu, grid%muv, grid%mut, grid%php, grid%alt, grid%p, &
635	ids, ide, jds, jde, kds, kde, &
636	ims, ime, jms, jme, kms, kme, &
637	ips, ipe, jps, jpe, kps, kpe, &
638	grid%i_start(ij), grid%i_end(ij), &
639	grid%j_start(ij), grid%j_end(ij), &
640	k_start, k_end )
641	CALL set_physical_bc3d( grid%al, 'p', config_flags, &
642	ids, ide, jds, jde, kds, kde, &
643	ims, ime, jms, jme, kms, kme, &
644	ips, ipe, jps, jpe, kps, kpe, &
645	grid%i_start(ij), grid%i_end(ij), &
646	grid%j_start(ij), grid%j_end(ij), &
647	k_start , k_end )
648	CALL set_physical_bc3d( grid%ph_2, 'w', config_flags, &
649	ids, ide, jds, jde, kds, kde, &
650	ims, ime, jms, jme, kms, kme, &
651	ips, ipe, jps, jpe, kps, kpe, &
652	grid%i_start(ij), grid%i_end(ij), &
653	grid%j_start(ij), grid%j_end(ij), &
654	k_start, k_end )
655
656	#ifdef LMDZ1
657	WRITE(message, *)' dyn_em: before polar'
658	CALL wrf_debug(200, message)
659	WRITE(message, *)' t_tend: ',t_tendf(im2,km2,jm2), &
660	' u_tend: ', ru_tendf(im2,1,jm2)
661	CALL wrf_debug(200, message)
662	WRITE(message,*)' psfc_tend: ',grid%dpsdt(im2,jm2), &
663	'p sfc: ',p8w(im2,kms,jm2)
664	CALL wrf_debug(200, message)
665	WRITE(message,*)' p 1: ',grid%p(im2,kms,jm2), ' ph 1: ',grid%ph_2(im2,kms,jm2)
666	CALL wrf_debug(200, message)
667	#endif
668
669	IF (config_flags%polar) THEN
670
671	!-------------------------------------------------------
672	! lat-lon grid pole-point (v) specification (extrapolate v, rv to the pole)
673	!-------------------------------------------------------
674
675	CALL pole_point_bc ( grid%v_1, &
676	ids, ide, jds, jde, kds, kde, &
677	ims, ime, jms, jme, kms, kme, &
678	grid%i_start(ij), grid%i_end(ij), &
679	grid%j_start(ij), grid%j_end(ij), &
680	k_start, k_end )
681
682	CALL pole_point_bc ( grid%v_2, &
683	ids, ide, jds, jde, kds, kde, &
684	ims, ime, jms, jme, kms, kme, &
685	grid%i_start(ij), grid%i_end(ij), &
686	grid%j_start(ij), grid%j_end(ij), &
687	k_start, k_end )
688
689	!-------------------------------------------------------
690	! end lat-lon grid pole-point (v) specification
691	!-------------------------------------------------------
692
693	ENDIF
694	END DO
695	!$OMP END PARALLEL DO
696
697	BENCH_END(set_phys_bc_tim)
698
699	rk_step_is_one : IF (rk_step == 1) THEN ! only need to initialize diffusion tendencies
700
701	!<DESCRIPTION>
702	!<pre>
703	!(2) The non-timesplit physics begins with a call to "phy_prep"
704	! (which computes some diagnostic variables such as temperature,
705	! pressure, u and v at p points, etc). This is followed by
706	! calls to the physics drivers:
707	!
708	! radiation,
709	! surface,
710	! pbl,
711	! cumulus,
712	! fddagd,
713	! 3D TKE and mixing.
714	!<pre>
715	!</DESCRIPTION>
716
717	#ifdef LMDZ
718	WRITE(message, *)' dyn_em: pre step1'
719	CALL wrf_debug(200, message)
720	WRITE(message, *)' t_tend: ',t_tendf(im2,km2,jm2), &
721	' u_tend: ', ru_tendf(im2,1,jm2)
722	CALL wrf_debug(200, message)
723	IF (config_flags%lmdz_physics) THEN
724	IF (config_flags%mp_physics + config_flags%ra_lw_physics + &
725	config_flags%ra_sw_physics + config_flags%sf_sfclay_physics + &
726	config_flags%bl_pbl_physics + config_flags%cu_physics /= 0) THEN
727	PRINT *,TRIM(errmsg)
728	PRINT *,' LMDZ physics are selected from namelist. lmdz_physics= ', &
729	config_flags%lmdz_physics
730	PRINT *,' Which requires no WRF physics schemes [0 value] and tey are not:'
731	PRINT '(2(A12,1x,I2,8x))',' mp =', config_flags%mp_physics,'ra_lw =', &
732	config_flags%ra_lw_physics
733	PRINT '(2(A12,1x,I2,8x))',' ra_sw =', config_flags%ra_sw_physics, &
734	'sf_sfclay =',config_flags%sf_sfclay_physics
735	PRINT '(2(A12,1x,I2,8x))',' bl_pbl =', config_flags%bl_pbl_physics, &
736	'cu =', config_flags%cu_physics
737	message = 'WRONG namelist ste-up'
738	CALL wrf_error_fatal(TRIM(message))
739	END IF
740
741	END IF
742	WRITE(message, *)' dyn_em: before step1'
743	CALL wrf_debug(200, message)
744	WRITE(message, *)' t_tend: ',t_tendf(im2,km2,jm2), &
745	' u_tend: ', ru_tendf(im2,1,jm2)
746	CALL wrf_debug(200, message)
747	WRITE(message,*)' psfc_tend: ',grid%dpsdt(im2,jm2), &
748	'p sfc: ',p8w(im2,kms,jm2)
749	CALL wrf_debug(200, message)
750	WRITE(message,*)' p 1: ',grid%p(im2,kms,jm2), ' ph 1: ',grid%ph_2(im2,kms,jm2)
751	CALL wrf_debug(200, message)
752
753	#endif
754
755	CALL first_rk_step_part1 ( grid, config_flags &
756	, moist , moist_tend &
757	, chem , chem_tend &
758	, tracer, tracer_tend &
759	, scalar , scalar_tend &
760	, fdda3d, fdda2d &
761	, ru_tendf, rv_tendf &
762	, rw_tendf, t_tendf &
763	, ph_tendf, mu_tendf &
764	, tke_tend &
765	, adapt_step_flag , curr_secs &
766	, psim , psih , wspd , gz1oz0 &
767	, br , chklowq &
768	, cu_act_flag , hol , th_phy &
769	, pi_phy , p_phy , grid%t_phy &
770	, u_phy , v_phy &
771	, dz8w , p8w , t8w , rho_phy , rho &
772	, ids, ide, jds, jde, kds, kde &
773	, ims, ime, jms, jme, kms, kme &
774	, ips, ipe, jps, jpe, kps, kpe &
775	, imsx, imex, jmsx, jmex, kmsx, kmex &
776	, ipsx, ipex, jpsx, jpex, kpsx, kpex &
777	, imsy, imey, jmsy, jmey, kmsy, kmey &
778	, ipsy, ipey, jpsy, jpey, kpsy, kpey &
779	, k_start , k_end &
780	, f_flux &
781	)
782
783	#ifdef LMDZ1
784	WRITE(message, *)' dyn_em: post step1 pre lmdz'
785	CALL wrf_debug(200, message)
786	WRITE(message, *)' t_tend: ',t_tendf(im2,km2,jm2), &
787	' u_tend: ', ru_tendf(im2,1,jm2)
788	CALL wrf_debug(200, message)
789	WRITE(message,*)' psfc_tend: ',grid%dpsdt(im2,jm2), &
790	'p sfc: ',p8w(im2,kms,jm2)
791	CALL wrf_debug(200, message)
792	WRITE(message,*)' p 1: ',grid%p(im2,kms,jm2), ' ph 1: ',grid%ph_2(im2,kms,jm2)
793	CALL wrf_debug(200, message)
794	#endif
795
796	#ifdef DM_PARALLEL
797	IF ( config_flags%bl_pbl_physics == MYNNPBLSCHEME2 .OR. &
798	config_flags%bl_pbl_physics == MYNNPBLSCHEME3 ) THEN
799	# include "HALO_EM_SCALAR_E_5.inc"
800	ENDIF
801	#endif
802
803	#ifdef LMDZ
804	!!! Using grid%clock & domain_get_sim_start_time
805	!!!!
806	initime = domain_get_start_time(grid)
807	CALL domain_clock_get(grid, current_time = simtime)
808	!! timeSinceSimStart = domain_get_time_since_sim_start( grid )
809	timeSinceSimStart = simtime - initime
810	CALL WRFU_TimeIntervalGet( timeSinceSimStart, &
811	D=days, S=seconds, Sn=Sn, Sd=Sd, rc=rc )
812	IF ( rc /= WRFU_SUCCESS ) THEN
813	CALL wrf_error_fatal ( &
814	'domain_clock_get: WRFU_TimeIntervalGet() failed' )
815	ENDIF
816	! get rid of hard-coded constants
817	minSinceSimStart = ( REAL( days ) * 24. * 60. ) + &
818	( REAL( seconds ) / 60. )
819	IF ( Sd /= 0 ) THEN
820	minSinceSimStart = minSinceSimStart + &
821	( ( REAL( Sn ) / REAL( Sd ) ) / 60. )
822	ENDIF
823	IF (minSinceSimStart == 0. .AND. grid%id == 1) THEN
824	PRINT *,' WRF+LMDZ: simulation is starting or comes from a restart!'
825	PRINT *,' since ', minSinceSimStart,' minutes has passed since it started'
826	!! wrftestrst = .TRUE.
827	wrftestrst = .TRUE.
828	ELSE
829	wrftestrst = .FALSE.
830	END IF
831	PRINT *,' Lluis: minSinceSimStart: ', minSinceSimStart
832	! Checking for time for input from auxiliar input 4
833	IF( WRFU_AlarmIsRinging( grid%alarms( first_auxinput + 3 ), rc=rc ) ) THEN
834	PRINT *,' WRF lowbdy time-step!!!!!'
835	wrftestin = .TRUE.
836	CALL WRFU_AlarmRingerOff( grid%alarms( first_auxinput + 3 ), rc=rc )
837	ELSE
838	wrftestin = .FALSE.
839	END IF
840
841	PRINT *,' grid id: ',grid%id
842
843	CALL nl_get_mminlu ( grid%id, mminlu )
844	IF (config_flags%sf_surface_physics == RUCLSMSCHEME) THEN
845	mminsl = 'STAS-RUC'
846	ELSE
847	mminsl = 'STAS'
848	END IF
849
850	grid%qv_2 = moist(:,:,:,P_QV)
851
852	! L. Fita, LMD. July 2014. Getting hour of the day
853	CALL domain_clock_get( grid, current_time=CurrTime )
854	CALL WRFU_TimeGet( CurrTime, YY= julyr, dayOfYear=julday, H=hr, M=minute, S=sec, MS=ms, rc=rc)
855	! Julian day hour (0, 1) !!
856	gmt=(hr+real(minute)/60.+real(sec)/3600.+real(ms)/(1000*3600))/24.
857
858	! Checking for NaNs (should not be necessary but....)
859	im2 = ims + (ime - ims) / 2
860	jm2 = jms + (jme - jms) / 2
861
862	IF (grid%t_2(im2,iz,jm2) /= grid%t_2(im2,iz,jm2) .OR. ABS(grid%t_2(im2,iz,jm2)) > 10000. ) THEN
863	PRINT errmsg
864	WRITE(wrF_err_message,*)'solve_em: wrong T value=', &
865	grid%t_2(im2,iz,jm2),' at: ', im2,', ', iz,', ', jm2,' !!!'
866	#ifdef DM_PARALLEL
867	CALL wrf_error_fatal(TRIM(wrf_err_message))
868	#else
869	PRINT *,TRIM(wrf_err_message)
870	STOP
871	#endif
872	END IF
873	! Checking for NaNs (should not be necessary but....)
874	IF (grid%psfc(im2,jm2) /= grid%psfc(im2,jm2) .OR. ABS(grid%psfc(im2,jm2)) > 1000000. ) THEN
875	PRINT errmsg
876	WRITE(wrF_err_message,*)'solve_em: wrong PSFC value=', &
877	grid%psfc(im2,jm2),' at: ', im2 ,', ', jm2, ' !!!'
878	#ifdef DM_PARALLEL
879	CALL wrf_error_fatal(TRIM(wrf_err_message))
880	#else
881	PRINT *,TRIM(wrf_err_message)
882	STOP
883	#endif
884	END IF
885
886	IF (config_flags%lmdz_physics) THEN
887
888	CALL call_lmdz_phys( &
889	& WRF_GRID=grid, WRF_XTIME=grid%xtime, &
890	& WRF_RESTART_ALARM=WRFU_AlarmIsRinging( grid%alarms( RESTART_ALARM ), rc=rc)&
891	& ,WRF_LON = grid%xlong, WRF_LAT=grid%xlat, &
892	& WRF_T=grid%t_2, WRF_U=grid%u_2, WRF_V=grid%v_2, &
893	& WRF_PERP=grid%P, WRF_BASEP=grid%PB, &
894	& WBDYW=config_flags%spec_bdy_width, WRF_ISRESTART=wrftestrst, &
895	& WRF_ISLOWBDYIN=wrftestin &
896	! Dimension arguments
897	& ,WIDS=ids,WIDE=ide, WJDS=jds,WJDE=jde, WKDS=kds,WKDE=kde &
898	& ,WIMS=ims,WIME=ime, WJMS=jms,WJME=jme, WKMS=kms,WKME=kme &
899	& ,WIPS=ips,WIPE=ipe, WJPS=jps,WJPE=jpe, WKPS=kps,WKPE=kpe &
900	& ,WI_START=grid%i_start,WI_END=MIN(grid%i_end, ide-1) &
901	& ,WJ_START=grid%j_start,WJ_END=MIN(grid%j_end, jde-1) &
902	& ,WKTS=k_start, WKTE=MIN(k_end,kde-1) &
903	& ,WNUM_TILES=grid%num_tiles &
904	& ,WNUM3DM=num_3d_m, WPARFIRSTSCAL=PARAM_FIRST_SCALAR, &
905	& WNX=config_flags%e_we, WNY=config_flags%e_sn, &
906	& WNZ=config_flags%e_vert, WJULDAY=FLOAT(julday), WGMT=gmt, &
907	& WTIME_STEP=REAL(config_flags%time_step), &
908	& WRF_FULLETA=grid%znw, WRF_HALFETA=grid%znu, WRF_DFULLETA=grid%dnw, &
909	& WRF_FULLPRES=p8w, WRF_PERGEOPOT=grid%ph_2, &
910	& WRF_BASEGEOPOT=grid%phb, &
911	& WRF_MOIST=grid%moist, WRF_W=grid%w_2, &
912	& WRF_PTOP=config_flags%p_top_requested, &
913	& WRF_PERMASS=grid%mu_1, WRF_BASEMASS=grid%mub, &
914	& WRF_MUT=grid%mut, WRF_MUU=grid%muu, WRF_MUV=grid%muv, &
915	!! & WRF_UTEND=grid%ru_tend, WRF_VTEND=grid%rv_tend, &
916	!! & WRF_TTEND=t_tend, &
917	& WRF_UTEND=ru_tendf, WRF_VTEND=rv_tendf, &
918	& WRF_TTEND=t_tendf, &
919	& WRF_MOISTTEND=moist_tend, WRF_PSFCTEND=grid%dpsdt, &
920	& WRF_QVID=P_QV, WRF_QCID=P_QC, WRF_QRID=P_QR, &
921	& WRF_QSID=P_QS, WRF_QIID=P_QI, WRF_QHID=P_QH, WRF_QGID=P_QG, &
922	! L. Fita. July 2013. Now defined as local dummy variables
923	! & WRF_DUDYN=???????????????, WRF_PVTHETA=????????????????????????, &
924	! & WRF_CLESPHY=?????????????, WRF_PRESNIVS=???????????????????????, &
925	& WRF_MAPFT=grid%msft, WRF_MAPFU=grid%msfu, WRF_MAPFV=grid%msfv, &
926	& WRF_DX=grid%dx, WRF_DY=grid%dy, &
927	& WRF_DBG=model_config_rec%debug_level, LANDCAT=mminlu, &
928	& SOILCAT=mminsl, &
929	& WRF_L_PBL=config_flags%lmdz_iflag_pbl, &
930	& WRF_L_CON=config_flags%lmdz_iflag_con, &
931	& WRF_L_THERMALS=config_flags%lmdz_iflag_thermals, &
932	& WRF_L_WAKE=config_flags%lmdz_iflag_wake, &
933	& wrf_nsoillayers=config_flags%num_soil_layers, &
934	& ICHECK_P=config_flags%i_check_point, &
935	& JCHECK_P=config_flags%j_check_point, &
936	& KCHECK_P=config_flags%k_check_point &
937	& )
938	END IF
939	#endif
940
941	#ifdef LMDZ1
942	WRITE(message, *)' dyn_em: post lmdz pre step2 rk_step:', rk_step
943	CALL wrf_debug(200, message)
944	WRITE(message, *)' t_tend: ',t_tendf(im2,km2,jm2), &
945	' u_tend: ', ru_tendf(im2,1,jm2)
946	CALL wrf_debug(200, message)
947	WRITE(message,*)' psfc_tend: ',grid%dpsdt(im2,jm2), &
948	'p sfc: ',p8w(im2,kms,jm2)
949	CALL wrf_debug(200, message)
950	WRITE(message,*)' p 1: ',grid%p(im2,kms,jm2), ' ph 1: ',grid%ph_2(im2,kms,jm2)
951	CALL wrf_debug(200, message)
952	#endif
953
954	CALL first_rk_step_part2 ( grid, config_flags &
955	, moist , moist_tend &
956	, chem , chem_tend &
957	, tracer, tracer_tend &
958	, scalar , scalar_tend &
959	, fdda3d, fdda2d &
960	, ru_tendf, rv_tendf &
961	, rw_tendf, t_tendf &
962	, ph_tendf, mu_tendf &
963	, tke_tend &
964	, adapt_step_flag , curr_secs &
965	, psim , psih , wspd , gz1oz0 &
966	, br , chklowq &
967	, cu_act_flag , hol , th_phy &
968	, pi_phy , p_phy , grid%t_phy &
969	, u_phy , v_phy &
970	, dz8w , p8w , t8w , rho_phy , rho &
971	, nba_mij, num_nba_mij & !JDM
972	, nba_rij, num_nba_rij & !JDM
973	, ids, ide, jds, jde, kds, kde &
974	, ims, ime, jms, jme, kms, kme &
975	, ips, ipe, jps, jpe, kps, kpe &
976	, imsx, imex, jmsx, jmex, kmsx, kmex &
977	, ipsx, ipex, jpsx, jpex, kpsx, kpex &
978	, imsy, imey, jmsy, jmey, kmsy, kmey &
979	, ipsy, ipey, jpsy, jpey, kpsy, kpey &
980	, k_start , k_end &
981	)
982
983	END IF rk_step_is_one
984
985	#ifdef LMDZ1
986	WRITE(message, *)' dyn_em: post step2 pre rk_tendency rk_step: ',rk_step
987	CALL wrf_debug(200, message)
988	WRITE(message, *)' t_tend: ',t_tendf(im2,km2,jm2), &
989	' u_tend: ', ru_tendf(im2,1,jm2)
990	CALL wrf_debug(200, message)
991	WRITE(message,*)' psfc_tend: ',grid%dpsdt(im2,jm2), &
992	'p sfc: ',p8w(im2,kms,jm2)
993	CALL wrf_debug(200, message)
994	WRITE(message,*)' p 1: ',grid%p(im2,kms,jm2), ' ph 1: ',grid%ph_2(im2,kms,jm2)
995	CALL wrf_debug(200, message)
996	#endif
997
998	BENCH_START(rk_tend_tim)
999	!$OMP PARALLEL DO &
1000	!$OMP PRIVATE ( ij )
1001	DO ij = 1 , grid%num_tiles
1002
1003	CALL wrf_debug ( 200 , ' call rk_tendency' )
1004	CALL rk_tendency ( config_flags, rk_step &
1005	,grid%ru_tend, grid%rv_tend, rw_tend, ph_tend, t_tend &
1006	,ru_tendf, rv_tendf, rw_tendf, ph_tendf, t_tendf &
1007	,mu_tend, grid%u_save, grid%v_save, w_save, ph_save &
1008	,grid%t_save, mu_save, grid%rthften &
1009	,grid%ru, grid%rv, grid%rw, grid%ww &
1010	,grid%u_2, grid%v_2, grid%w_2, grid%t_2, grid%ph_2 &
1011	,grid%u_1, grid%v_1, grid%w_1, grid%t_1, grid%ph_1 &
1012	,grid%h_diabatic, grid%phb, grid%t_init &
1013	,grid%mu_2, grid%mut, grid%muu, grid%muv, grid%mub &
1014	,grid%al, grid%alt, grid%p, grid%pb, grid%php, cqu, cqv, cqw &
1015	,grid%u_base, grid%v_base, grid%t_base, grid%qv_base, grid%z_base &
1016	,grid%msfux,grid%msfuy, grid%msfvx, grid%msfvx_inv &
1017	,grid%msfvy, grid%msftx,grid%msfty, grid%clat, grid%f, grid%e, grid%sina, grid%cosa &
1018	,grid%fnm, grid%fnp, grid%rdn, grid%rdnw &
1019	,grid%dt, grid%rdx, grid%rdy, grid%khdif, grid%kvdif, grid%xkmh, grid%xkhh &
1020	,grid%diff_6th_opt, grid%diff_6th_factor &
1021	,grid%dampcoef,grid%zdamp,config_flags%damp_opt,config_flags%rad_nudge &
1022	,grid%cf1, grid%cf2, grid%cf3, grid%cfn, grid%cfn1, num_3d_m &
1023	,config_flags%non_hydrostatic, config_flags%top_lid &
1024	,grid%u_frame, grid%v_frame &
1025	,ids, ide, jds, jde, kds, kde &
1026	,ims, ime, jms, jme, kms, kme &
1027	,grid%i_start(ij), grid%i_end(ij) &
1028	,grid%j_start(ij), grid%j_end(ij) &
1029	,k_start, k_end &
1030	,max_vert_cfl_tmp(ij), max_horiz_cfl_tmp(ij) )
1031
1032	END DO
1033	#ifdef LMDZ1
1034	WRITE(message, *)' dyn_em: post rk_tendency'
1035	CALL wrf_debug(200, message)
1036	WRITE(message, *)' t_tend: ',t_tendf(im2,km2,jm2), &
1037	' u_tend: ', ru_tendf(im2,1,jm2)
1038	CALL wrf_debug(200, message)
1039	WRITE(message,*)' psfc_tend: ',grid%dpsdt(im2,jm2), &
1040	'p sfc: ',p8w(im2,kms,jm2)
1041	CALL wrf_debug(200, message)
1042	WRITE(message,*)' p 1: ',grid%p(im2,kms,jm2), ' ph 1: ',grid%ph_2(im2,kms,jm2)
1043	CALL wrf_debug(200, message)
1044	#endif
1045	!$OMP END PARALLEL DO
1046	BENCH_END(rk_tend_tim)
1047
1048	IF (config_flags%use_adaptive_time_step) THEN
1049	DO ij = 1 , grid%num_tiles
1050	IF (max_horiz_cfl_tmp(ij) .GT. grid%max_horiz_cfl) THEN
1051	grid%max_horiz_cfl = max_horiz_cfl_tmp(ij)
1052	ENDIF
1053	IF (max_vert_cfl_tmp(ij) .GT. grid%max_vert_cfl) THEN
1054	grid%max_vert_cfl = max_vert_cfl_tmp(ij)
1055	ENDIF
1056	END DO
1057
1058	IF (grid%max_horiz_cfl .GT. grid%max_cfl_val) THEN
1059	grid%max_cfl_val = grid%max_horiz_cfl
1060	ENDIF
1061	IF (grid%max_vert_cfl .GT. grid%max_cfl_val) THEN
1062	grid%max_cfl_val = grid%max_vert_cfl
1063	ENDIF
1064	ENDIF
1065
1066	BENCH_START(relax_bdy_dry_tim)
1067	!$OMP PARALLEL DO &
1068	!$OMP PRIVATE ( ij )
1069	DO ij = 1 , grid%num_tiles
1070
1071	IF( (config_flags%specified .or. config_flags%nested) .and. rk_step == 1 ) THEN
1072
1073	CALL relax_bdy_dry ( config_flags, &
1074	grid%u_save, grid%v_save, ph_save, grid%t_save, &
1075	w_save, mu_tend, &
1076	grid%ru, grid%rv, grid%ph_2, grid%t_2, &
1077	grid%w_2, grid%mu_2, grid%mut, &
1078	grid%u_bxs,grid%u_bxe,grid%u_bys,grid%u_bye, &
1079	grid%v_bxs,grid%v_bxe,grid%v_bys,grid%v_bye, &
1080	grid%ph_bxs,grid%ph_bxe,grid%ph_bys,grid%ph_bye, &
1081	grid%t_bxs,grid%t_bxe,grid%t_bys,grid%t_bye, &
1082	grid%w_bxs,grid%w_bxe,grid%w_bys,grid%w_bye, &
1083	grid%mu_bxs,grid%mu_bxe,grid%mu_bys,grid%mu_bye, &
1084	grid%u_btxs,grid%u_btxe,grid%u_btys,grid%u_btye, &
1085	grid%v_btxs,grid%v_btxe,grid%v_btys,grid%v_btye, &
1086	grid%ph_btxs,grid%ph_btxe,grid%ph_btys,grid%ph_btye, &
1087	grid%t_btxs,grid%t_btxe,grid%t_btys,grid%t_btye, &
1088	grid%w_btxs,grid%w_btxe,grid%w_btys,grid%w_btye, &
1089	grid%mu_btxs,grid%mu_btxe,grid%mu_btys,grid%mu_btye, &
1090	config_flags%spec_bdy_width, grid%spec_zone, grid%relax_zone, &
1091	grid%dtbc, grid%fcx, grid%gcx, &
1092	ids,ide, jds,jde, kds,kde, &
1093	ims,ime, jms,jme, kms,kme, &
1094	ips,ipe, jps,jpe, kps,kpe, &
1095	grid%i_start(ij), grid%i_end(ij), &
1096	grid%j_start(ij), grid%j_end(ij), &
1097	k_start, k_end )
1098
1099	ENDIF
1100
1101	CALL rk_addtend_dry( grid%ru_tend, grid%rv_tend, rw_tend, ph_tend, t_tend, &
1102	ru_tendf, rv_tendf, rw_tendf, ph_tendf, t_tendf, &
1103	grid%u_save, grid%v_save, w_save, ph_save, grid%t_save, &
1104	mu_tend, mu_tendf, rk_step, &
1105	grid%h_diabatic, grid%mut, grid%msftx, &
1106	grid%msfty, grid%msfux,grid%msfuy, &
1107	grid%msfvx, grid%msfvx_inv, grid%msfvy, &
1108	ids,ide, jds,jde, kds,kde, &
1109	ims,ime, jms,jme, kms,kme, &
1110	ips,ipe, jps,jpe, kps,kpe, &
1111	grid%i_start(ij), grid%i_end(ij), &
1112	grid%j_start(ij), grid%j_end(ij), &
1113	k_start, k_end )
1114
1115	IF( config_flags%specified .or. config_flags%nested ) THEN
1116	CALL spec_bdy_dry ( config_flags, &
1117	grid%ru_tend, grid%rv_tend, ph_tend, t_tend, &
1118	rw_tend, mu_tend, &
1119	grid%u_bxs,grid%u_bxe,grid%u_bys,grid%u_bye, &
1120	grid%v_bxs,grid%v_bxe,grid%v_bys,grid%v_bye, &
1121	grid%ph_bxs,grid%ph_bxe,grid%ph_bys,grid%ph_bye, &
1122	grid%t_bxs,grid%t_bxe,grid%t_bys,grid%t_bye, &
1123	grid%w_bxs,grid%w_bxe,grid%w_bys,grid%w_bye, &
1124	grid%mu_bxs,grid%mu_bxe,grid%mu_bys,grid%mu_bye, &
1125	grid%u_btxs,grid%u_btxe,grid%u_btys,grid%u_btye, &
1126	grid%v_btxs,grid%v_btxe,grid%v_btys,grid%v_btye, &
1127	grid%ph_btxs,grid%ph_btxe,grid%ph_btys,grid%ph_btye, &
1128	grid%t_btxs,grid%t_btxe,grid%t_btys,grid%t_btye, &
1129	grid%w_btxs,grid%w_btxe,grid%w_btys,grid%w_btye, &
1130	grid%mu_btxs,grid%mu_btxe,grid%mu_btys,grid%mu_btye, &
1131	config_flags%spec_bdy_width, grid%spec_zone, &
1132	ids,ide, jds,jde, kds,kde, & ! domain dims
1133	ims,ime, jms,jme, kms,kme, & ! memory dims
1134	ips,ipe, jps,jpe, kps,kpe, & ! patch dims
1135	grid%i_start(ij), grid%i_end(ij), &
1136	grid%j_start(ij), grid%j_end(ij), &
1137	k_start, k_end )
1138
1139	ENDIF
1140
1141	END DO
1142	!$OMP END PARALLEL DO
1143	BENCH_END(relax_bdy_dry_tim)
1144
1145	!<DESCRIPTION>
1146	!<pre>
1147	! (3) Small (acoustic,sound) steps.
1148	!
1149	! Several acoustic steps are taken each RK pass. A small step
1150	! sequence begins with calculating perturbation variables
1151	! and coupling them to the column dry-air-mass mu
1152	! (call to small_step_prep). This is followed by computing
1153	! coefficients for the vertically implicit part of the
1154	! small timestep (call to calc_coef_w).
1155	!
1156	! The small steps are taken
1157	! in the named loop "small_steps:". In the small_steps loop, first
1158	! the horizontal momentum (u and v) are advanced (call to advance_uv),
1159	! next mu and theta are advanced (call to advance_mu_t) followed by
1160	! advancing w and the geopotential (call to advance_w). Diagnostic
1161	! values for pressure and inverse density are updated at the end of
1162	! each small_step.
1163	!
1164	! The small-step section ends with the change of the perturbation variables
1165	! back to full variables (call to small_step_finish).
1166	!</pre>
1167	!</DESCRIPTION>
1168
1169	BENCH_START(small_step_prep_tim)
1170	!$OMP PARALLEL DO &
1171	!$OMP PRIVATE ( ij )
1172	DO ij = 1 , grid%num_tiles
1173
1174	! Calculate coefficients for the vertically implicit acoustic/gravity wave
1175	! integration. We only need calculate these for the first pass through -
1176	! the predictor step. They are reused as is for the corrector step.
1177	! For third-order RK, we need to recompute these after the first
1178	! predictor because we may have changed the small timestep -> grid%dts.
1179
1180	CALL wrf_debug ( 200 , ' call small_step_prep ' )
1181
1182	#ifdef LMDZ1
1183	WRITE(message, *)' dyn_em: before small_step_prep'
1184	CALL wrf_debug(200, message)
1185	WRITE(message, *)' t_tend: ',t_tendf(im2,km2,jm2), &
1186	' u_tend: ', ru_tendf(im2,1,jm2)
1187	CALL wrf_debug(200, message)
1188	WRITE(message,*)' psfc_tend: ',grid%dpsdt(im2,jm2), &
1189	'p sfc: ',p8w(im2,kms,jm2)
1190	CALL wrf_debug(200, message)
1191	WRITE(message,*)' p 1: ',grid%p(im2,kms,jm2), ' ph 1: ',grid%ph_2(im2,kms,jm2)
1192	CALL wrf_debug(200, message)
1193	#endif
1194	CALL small_step_prep( grid%u_1,grid%u_2,grid%v_1,grid%v_2,grid%w_1,grid%w_2, &
1195	grid%t_1,grid%t_2,grid%ph_1,grid%ph_2, &
1196	grid%mub, grid%mu_1, grid%mu_2, &
1197	grid%muu, muus, grid%muv, muvs, &
1198	grid%mut, grid%muts, grid%mudf, &
1199	grid%u_save, grid%v_save, w_save, &
1200	grid%t_save, ph_save, mu_save, &
1201	grid%ww, ww1, &
1202	grid%dnw, c2a, grid%pb, grid%p, grid%alt, &
1203	grid%msfux,grid%msfuy, grid%msfvx, grid%msfvx_inv, &
1204	grid%msfvy, grid%msftx,grid%msfty, &
1205	grid%rdx, grid%rdy, rk_step, &
1206	ids, ide, jds, jde, kds, kde, &
1207	ims, ime, jms, jme, kms, kme, &
1208	grid%i_start(ij), grid%i_end(ij), &
1209	grid%j_start(ij), grid%j_end(ij), &
1210	k_start , k_end )
1211
1212	#ifdef LMDZ1
1213	WRITE(message, *)' dyn_em: post small_step_prep'
1214	CALL wrf_debug(200, message)
1215	WRITE(message, *)' t_tend: ',t_tendf(im2,km2,jm2), &
1216	' u_tend: ', ru_tendf(im2,1,jm2)
1217	CALL wrf_debug(200, message)
1218	WRITE(message,*)' psfc_tend: ',grid%dpsdt(im2,jm2), &
1219	'p sfc: ',p8w(im2,kms,jm2)
1220	CALL wrf_debug(200, message)
1221	WRITE(message,*)' p 1: ',grid%p(im2,kms,jm2), ' ph 1: ',grid%ph_2(im2,kms,jm2)
1222	CALL wrf_debug(200, message)
1223	#endif
1224	CALL calc_p_rho( grid%al, grid%p, grid%ph_2, &
1225	grid%alt, grid%t_2, grid%t_save, c2a, pm1, &
1226	grid%mu_2, grid%muts, grid%znu, t0, &
1227	grid%rdnw, grid%dnw, grid%smdiv, &
1228	config_flags%non_hydrostatic, 0, &
1229	ids, ide, jds, jde, kds, kde, &
1230	ims, ime, jms, jme, kms, kme, &
1231	grid%i_start(ij), grid%i_end(ij), &
1232	grid%j_start(ij), grid%j_end(ij), &
1233	k_start , k_end )
1234	#ifdef LMDZ1
1235	WRITE(message, *)' dyn_em: post calc_p_rho'
1236	CALL wrf_debug(200, message)
1237	WRITE(message, *)' t_tend: ',t_tendf(im2,km2,jm2), &
1238	' u_tend: ', ru_tendf(im2,1,jm2)
1239	CALL wrf_debug(200, message)
1240	WRITE(message,*)' psfc_tend: ',grid%dpsdt(im2,jm2), &
1241	'p sfc: ',p8w(im2,kms,jm2)
1242	CALL wrf_debug(200, message)
1243	WRITE(message,*)' p 1: ',grid%p(im2,kms,jm2), ' ph 1: ',grid%ph_2(im2,kms,jm2)
1244	CALL wrf_debug(200, message)
1245	#endif
1246
1247	IF (config_flags%non_hydrostatic) THEN
1248	CALL calc_coef_w( a,alpha,gamma, &
1249	grid%mut, cqw, &
1250	grid%rdn, grid%rdnw, c2a, &
1251	dts_rk, g, grid%epssm, &
1252	config_flags%top_lid, &
1253	ids, ide, jds, jde, kds, kde, &
1254	ims, ime, jms, jme, kms, kme, &
1255	grid%i_start(ij), grid%i_end(ij), &
1256	grid%j_start(ij), grid%j_end(ij), &
1257	k_start , k_end )
1258	ENDIF
1259	#ifdef LMDZ1
1260	WRITE(message, *)' dyn_em: post calc_coef_w'
1261	CALL wrf_debug(200, message)
1262	WRITE(message, *)' t_tend: ',t_tendf(im2,km2,jm2), &
1263	' u_tend: ', ru_tendf(im2,1,jm2)
1264	CALL wrf_debug(200, message)
1265	WRITE(message,*)' psfc_tend: ',grid%dpsdt(im2,jm2), &
1266	'p sfc: ',p8w(im2,kms,jm2)
1267	CALL wrf_debug(200, message)
1268	WRITE(message,*)' p 1: ',grid%p(im2,kms,jm2), ' ph 1: ',grid%ph_2(im2,kms,jm2)
1269	CALL wrf_debug(200, message)
1270	WRITE(message,*)' al: ',grid%al(im2,km2,jm2), ' p: ', grid%p(im2,1,jm2), &
1271	' ph: ',grid%ph_2(im2,1,jm2)
1272	CALL wrf_debug(200, message)
1273	WRITE(message,*)' mu: ',grid%mu_2(im2,jm2), 'alt: ', grid%alt(im2,1,jm2), &
1274	' mu: ',grid%mu_2(im2,jm2), ' znu: ', grid%znu(1),' ph 1: ',grid%ph_2(im2,2,jm2)
1275	CALL wrf_debug(200, message)
1276	PRINT *,' c2a: ',c2a(im2,1,jm2), 't_2: ', grid%t_2(im2,1,jm2), &
1277	' t_save: ',grid%t_save(im2,1,jm2),' pm1: ',pm1(im2,1,jm2)
1278	CALL wrf_debug(200, message)
1279	#endif
1280	ENDDO
1281	!$OMP END PARALLEL DO
1282	BENCH_END(small_step_prep_tim)
1283
1284	#ifdef DM_PARALLEL
1285	!-----------------------------------------------------------------------
1286	! Stencils for patch communications (WCS, 29 June 2001)
1287	! Note: the small size of this halo exchange reflects the
1288	! fact that we are carrying the uncoupled variables
1289	! as state variables in the mass coordinate model, as
1290	! opposed to the coupled variables as in the height
1291	! coordinate model.
1292	!
1293	! * * * * *
1294	! * * * * * * * * *
1295	! * + * * + * * * + * *
1296	! * * * * * * * * *
1297	! * * * * *
1298	!
1299	! 3D variables - note staggering! ph_2(Z), u_save(X), v_save(Y)
1300	!
1301	! ph_2 x
1302	! al x
1303	! p x
1304	! t_1 x
1305	! t_save x
1306	! u_save x
1307	! v_save x
1308	!
1309	! the following are 2D (xy) variables
1310	!
1311	! mu_1 x
1312	! mu_2 x
1313	! mudf x
1314	! php x
1315	! alt x
1316	! pb x
1317	!--------------------------------------------------------------
1318	# include "HALO_EM_B.inc"
1319	# include "PERIOD_BDY_EM_B.inc"
1320	#endif
1321
1322	BENCH_START(set_phys_bc2_tim)
1323	!$OMP PARALLEL DO &
1324	!$OMP PRIVATE ( ij )
1325
1326	DO ij = 1 , grid%num_tiles
1327
1328	CALL set_physical_bc3d( grid%ru_tend, 'u', config_flags, &
1329	ids, ide, jds, jde, kds, kde, &
1330	ims, ime, jms, jme, kms, kme, &
1331	ips, ipe, jps, jpe, kps, kpe, &
1332	grid%i_start(ij), grid%i_end(ij), &
1333	grid%j_start(ij), grid%j_end(ij), &
1334	k_start , k_end )
1335
1336	CALL set_physical_bc3d( grid%rv_tend, 'v', config_flags, &
1337	ids, ide, jds, jde, kds, kde, &
1338	ims, ime, jms, jme, kms, kme, &
1339	ips, ipe, jps, jpe, kps, kpe, &
1340	grid%i_start(ij), grid%i_end(ij), &
1341	grid%j_start(ij), grid%j_end(ij), &
1342	k_start , k_end )
1343
1344	CALL set_physical_bc3d( grid%ph_2, 'w', config_flags, &
1345	ids, ide, jds, jde, kds, kde, &
1346	ims, ime, jms, jme, kms, kme, &
1347	ips, ipe, jps, jpe, kps, kpe, &
1348	grid%i_start(ij), grid%i_end(ij), &
1349	grid%j_start(ij), grid%j_end(ij), &
1350	k_start , k_end )
1351
1352	CALL set_physical_bc3d( grid%al, 'p', config_flags, &
1353	ids, ide, jds, jde, kds, kde, &
1354	ims, ime, jms, jme, kms, kme, &
1355	ips, ipe, jps, jpe, kps, kpe, &
1356	grid%i_start(ij), grid%i_end(ij), &
1357	grid%j_start(ij), grid%j_end(ij), &
1358	k_start , k_end )
1359
1360	CALL set_physical_bc3d( grid%p, 'p', config_flags, &
1361	ids, ide, jds, jde, kds, kde, &
1362	ims, ime, jms, jme, kms, kme, &
1363	ips, ipe, jps, jpe, kps, kpe, &
1364	grid%i_start(ij), grid%i_end(ij), &
1365	grid%j_start(ij), grid%j_end(ij), &
1366	k_start , k_end )
1367
1368	CALL set_physical_bc3d( grid%t_1, 'p', config_flags, &
1369	ids, ide, jds, jde, kds, kde, &
1370	ims, ime, jms, jme, kms, kme, &
1371	ips, ipe, jps, jpe, kps, kpe, &
1372	grid%i_start(ij), grid%i_end(ij), &
1373	grid%j_start(ij), grid%j_end(ij), &
1374	k_start , k_end )
1375
1376	CALL set_physical_bc3d( grid%t_save, 't', config_flags, &
1377	ids, ide, jds, jde, kds, kde, &
1378	ims, ime, jms, jme, kms, kme, &
1379	ips, ipe, jps, jpe, kps, kpe, &
1380	grid%i_start(ij), grid%i_end(ij), &
1381	grid%j_start(ij), grid%j_end(ij), &
1382	k_start , k_end )
1383
1384	CALL set_physical_bc2d( grid%mu_1, 't', config_flags, &
1385	ids, ide, jds, jde, &
1386	ims, ime, jms, jme, &
1387	ips, ipe, jps, jpe, &
1388	grid%i_start(ij), grid%i_end(ij), &
1389	grid%j_start(ij), grid%j_end(ij) )
1390
1391	CALL set_physical_bc2d( grid%mu_2, 't', config_flags, &
1392	ids, ide, jds, jde, &
1393	ims, ime, jms, jme, &
1394	ips, ipe, jps, jpe, &
1395	grid%i_start(ij), grid%i_end(ij), &
1396	grid%j_start(ij), grid%j_end(ij) )
1397
1398	CALL set_physical_bc2d( grid%mudf, 't', config_flags, &
1399	ids, ide, jds, jde, &
1400	ims, ime, jms, jme, &
1401	ips, ipe, jps, jpe, &
1402	grid%i_start(ij), grid%i_end(ij), &
1403	grid%j_start(ij), grid%j_end(ij) )
1404
1405	END DO
1406	!$OMP END PARALLEL DO
1407	BENCH_END(set_phys_bc2_tim)
1408	small_steps : DO iteration = 1 , number_of_small_timesteps
1409
1410	! Boundary condition time (or communication time).
1411	#ifdef DM_PARALLEL
1412	# include "PERIOD_BDY_EM_B.inc"
1413	#endif
1414
1415	!$OMP PARALLEL DO &
1416	!$OMP PRIVATE ( ij )
1417
1418	DO ij = 1 , grid%num_tiles
1419
1420	BENCH_START(advance_uv_tim)
1421	#ifdef LMDZ1
1422	WRITE(message, *)' dyn_em: before advance_uv'
1423	CALL wrf_debug(200, message)
1424	WRITE(message, *)' t_tend: ',t_tendf(im2,km2,jm2), &
1425	' u_tend: ', ru_tendf(im2,1,jm2)
1426	CALL wrf_debug(200, message)
1427	WRITE(message,*)' psfc_tend: ',grid%dpsdt(im2,jm2), &
1428	'p sfc: ',p8w(im2,kms,jm2)
1429	CALL wrf_debug(200, message)
1430	WRITE(message,*)' p 1: ',grid%p(im2,kms,jm2), ' ph 1: ',grid%ph_2(im2,kms,jm2)
1431	CALL wrf_debug(200, message)
1432	WRITE(message,*)' al: ',grid%al(im2,km2,jm2), ' p: ', grid%p(im2,1,jm2), &
1433	' ph: ',grid%ph_2(im2,1,jm2)
1434	CALL wrf_debug(200, message)
1435	WRITE(message,*)' mu: ',grid%mu_2(im2,jm2), 'alt: ', grid%alt(im2,1,jm2), &
1436	' mu: ',grid%mu_2(im2,jm2), ' znu: ', grid%znu(1),' ph 1: ',grid%ph_2(im2,2,jm2)
1437	CALL wrf_debug(200, message)
1438	PRINT *,' c2a: ',c2a(im2,1,jm2), 't_2: ', grid%t_2(im2,1,jm2), &
1439	' t_save: ',grid%t_save(im2,1,jm2),' pm1: ',pm1(im2,1,jm2)
1440	CALL wrf_debug(200, message)
1441	#endif
1442
1443	CALL advance_uv ( grid%u_2, grid%ru_tend, grid%v_2, grid%rv_tend, &
1444	grid%p, grid%pb, &
1445	grid%ph_2, grid%php, grid%alt, grid%al, &
1446	grid%mu_2, &
1447	grid%muu, cqu, grid%muv, cqv, grid%mudf, &
1448	grid%msfux, grid%msfuy, grid%msfvx, &
1449	grid%msfvx_inv, grid%msfvy, &
1450	grid%rdx, grid%rdy, dts_rk, &
1451	grid%cf1, grid%cf2, grid%cf3, grid%fnm, grid%fnp, &
1452	grid%emdiv, &
1453	grid%rdnw, config_flags,grid%spec_zone, &
1454	config_flags%non_hydrostatic, config_flags%top_lid, &
1455	ids, ide, jds, jde, kds, kde, &
1456	ims, ime, jms, jme, kms, kme, &
1457	grid%i_start(ij), grid%i_end(ij), &
1458	grid%j_start(ij), grid%j_end(ij), &
1459	k_start , k_end )
1460	BENCH_END(advance_uv_tim)
1461
1462	END DO
1463	!$OMP END PARALLEL DO
1464	#ifdef LMDZ1
1465	WRITE(message, *)' dyn_em: after advance_uv'
1466	CALL wrf_debug(200, message)
1467	WRITE(message, *)' t_tend: ',t_tendf(im2,km2,jm2), &
1468	' u_tend: ', ru_tendf(im2,1,jm2)
1469	CALL wrf_debug(200, message)
1470	WRITE(message,*)' psfc_tend: ',grid%dpsdt(im2,jm2), &
1471	'p sfc: ',p8w(im2,kms,jm2)
1472	CALL wrf_debug(200, message)
1473	WRITE(message,*)' p 1: ',grid%p(im2,kms,jm2), ' ph 1: ',grid%ph_2(im2,kms,jm2)
1474	CALL wrf_debug(200, message)
1475	WRITE(message,*)' al: ',grid%al(im2,km2,jm2), ' p: ', grid%p(im2,1,jm2), &
1476	' ph: ',grid%ph_2(im2,1,jm2)
1477	CALL wrf_debug(200, message)
1478	WRITE(message,*)' mu: ',grid%mu_2(im2,jm2), 'alt: ', grid%alt(im2,1,jm2), &
1479	' mu: ',grid%mu_2(im2,jm2), ' znu: ', grid%znu(1),' ph 1: ',grid%ph_2(im2,2,jm2)
1480	CALL wrf_debug(200, message)
1481	PRINT *,' c2a: ',c2a(im2,1,jm2), 't_2: ', grid%t_2(im2,1,jm2), &
1482	' t_save: ',grid%t_save(im2,1,jm2),' pm1: ',pm1(im2,1,jm2)
1483	CALL wrf_debug(200, message)
1484	#endif
1485
1486	!-----------------------------------------------------------
1487	! acoustic integration polar filter for smallstep u, v
1488	!-----------------------------------------------------------
1489
1490	IF (config_flags%polar) THEN
1491
1492	CALL pxft ( grid=grid &
1493	,lineno=__LINE__ &
1494	,flag_uv = 1 &
1495	,flag_rurv = 0 &
1496	,flag_wph = 0 &
1497	,flag_ww = 0 &
1498	,flag_t = 0 &
1499	,flag_mu = 0 &
1500	,flag_mut = 0 &
1501	,flag_moist = 0 &
1502	,flag_chem = 0 &
1503	,flag_tracer = 0 &
1504	,flag_scalar = 0 &
1505	,positive_definite = .FALSE. &
1506	,moist=moist,chem=chem,tracer=tracer,scalar=scalar &
1507	,fft_filter_lat = config_flags%fft_filter_lat &
1508	,dclat = dclat &
1509	,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
1510	,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
1511	,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe &
1512	,imsx=imsx,imex=imex,jmsx=jmsx,jmex=jmex,kmsx=kmsx,kmex=kmex &
1513	,ipsx=ipsx,ipex=ipex,jpsx=jmsx,jpex=jpex,kpsx=kpsx,kpex=kpex )
1514
1515	END IF
1516
1517	!-----------------------------------------------------------
1518	! end acoustic integration polar filter for smallstep u, v
1519	!-----------------------------------------------------------
1520
1521	!$OMP PARALLEL DO &
1522	!$OMP PRIVATE ( ij )
1523	DO ij = 1 , grid%num_tiles
1524
1525	BENCH_START(spec_bdy_uv_tim)
1526	IF( config_flags%specified .or. config_flags%nested ) THEN
1527	CALL spec_bdyupdate(grid%u_2, grid%ru_tend, dts_rk, &
1528	'u' , config_flags, &
1529	grid%spec_zone, &
1530	ids,ide, jds,jde, kds,kde, & ! domain dims
1531	ims,ime, jms,jme, kms,kme, & ! memory dims
1532	ips,ipe, jps,jpe, kps,kpe, & ! patch dims
1533	grid%i_start(ij), grid%i_end(ij), &
1534	grid%j_start(ij), grid%j_end(ij), &
1535	k_start , k_end )
1536
1537	CALL spec_bdyupdate(grid%v_2, grid%rv_tend, dts_rk, &
1538	'v' , config_flags, &
1539	grid%spec_zone, &
1540	ids,ide, jds,jde, kds,kde, & ! domain dims
1541	ims,ime, jms,jme, kms,kme, & ! memory dims
1542	ips,ipe, jps,jpe, kps,kpe, & ! patch dims
1543	grid%i_start(ij), grid%i_end(ij), &
1544	grid%j_start(ij), grid%j_end(ij), &
1545	k_start , k_end )
1546
1547	ENDIF
1548	BENCH_END(spec_bdy_uv_tim)
1549
1550	END DO
1551	!$OMP END PARALLEL DO
1552
1553	#ifdef DM_PARALLEL
1554	!
1555	! Stencils for patch communications (WCS, 29 June 2001)
1556	!
1557	! * *
1558	! * + * * + * +
1559	! * *
1560	!
1561	! u_2 x
1562	! v_2 x
1563	!
1564	# include "HALO_EM_C.inc"
1565	#endif
1566
1567	!$OMP PARALLEL DO &
1568	!$OMP PRIVATE ( ij )
1569	DO ij = 1 , grid%num_tiles
1570
1571	! advance the mass in the column, theta, and calculate ww
1572
1573	BENCH_START(advance_mu_t_tim)
1574	CALL advance_mu_t( grid%ww, ww1, grid%u_2, grid%u_save, grid%v_2, grid%v_save, &
1575	grid%mu_2, grid%mut, muave, grid%muts, grid%muu, grid%muv, &
1576	grid%mudf, grid%ru_m, grid%rv_m, grid%ww_m, &
1577	grid%t_2, grid%t_save, t_2save, t_tend, &
1578	mu_tend, &
1579	grid%rdx, grid%rdy, dts_rk, grid%epssm, &
1580	grid%dnw, grid%fnm, grid%fnp, grid%rdnw, &
1581	grid%msfux,grid%msfuy, grid%msfvx, grid%msfvx_inv, &
1582	grid%msfvy, grid%msftx,grid%msfty, &
1583	iteration, config_flags, &
1584	ids, ide, jds, jde, kds, kde, &
1585	ims, ime, jms, jme, kms, kme, &
1586	grid%i_start(ij), grid%i_end(ij), &
1587	grid%j_start(ij), grid%j_end(ij), &
1588	k_start , k_end )
1589	BENCH_END(advance_mu_t_tim)
1590	ENDDO
1591	!$OMP END PARALLEL DO
1592	#ifdef LMDZ1
1593	WRITE(message, *)' dyn_em: after advance_mut'
1594	CALL wrf_debug(200, message)
1595	WRITE(message, *)' t_tend: ',t_tendf(im2,km2,jm2), &
1596	' u_tend: ', ru_tendf(im2,1,jm2)
1597	CALL wrf_debug(200, message)
1598	WRITE(message,*)' psfc_tend: ',grid%dpsdt(im2,jm2), &
1599	'p sfc: ',p8w(im2,kms,jm2)
1600	CALL wrf_debug(200, message)
1601	WRITE(message,*)' p 1: ',grid%p(im2,kms,jm2), ' ph 1: ',grid%ph_2(im2,kms,jm2)
1602	CALL wrf_debug(200, message)
1603	WRITE(message,*)' al: ',grid%al(im2,km2,jm2), ' p: ', grid%p(im2,1,jm2), &
1604	' ph: ',grid%ph_2(im2,1,jm2)
1605	CALL wrf_debug(200, message)
1606	WRITE(message,*)' mu: ',grid%mu_2(im2,jm2), 'alt: ', grid%alt(im2,1,jm2), &
1607	' mu: ',grid%mu_2(im2,jm2), ' znu: ', grid%znu(1),' ph 1: ',grid%ph_2(im2,2,jm2)
1608	CALL wrf_debug(200, message)
1609	PRINT *,' c2a: ',c2a(im2,1,jm2), 't_2: ', grid%t_2(im2,1,jm2), &
1610	' t_save: ',grid%t_save(im2,1,jm2),' pm1: ',pm1(im2,1,jm2)
1611	CALL wrf_debug(200, message)
1612	#endif
1613
1614	!-----------------------------------------------------------
1615	! acoustic integration polar filter for smallstep mu, t
1616	!-----------------------------------------------------------
1617
1618	IF ( (config_flags%polar) ) THEN
1619
1620	CALL pxft ( grid=grid &
1621	,lineno=__LINE__ &
1622	,flag_uv = 0 &
1623	,flag_rurv = 0 &
1624	,flag_wph = 0 &
1625	,flag_ww = 0 &
1626	,flag_t = 1 &
1627	,flag_mu = 1 &
1628	,flag_mut = 0 &
1629	,flag_moist = 0 &
1630	,flag_chem = 0 &
1631	,flag_tracer = 0 &
1632	,flag_scalar = 0 &
1633	,positive_definite = .FALSE. &
1634	,moist=moist,chem=chem,tracer=tracer,scalar=scalar &
1635	,fft_filter_lat = config_flags%fft_filter_lat &
1636	,dclat = dclat &
1637	,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
1638	,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
1639	,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe &
1640	,imsx=imsx,imex=imex,jmsx=jmsx,jmex=jmex,kmsx=kmsx,kmex=kmex &
1641	,ipsx=ipsx,ipex=ipex,jpsx=jmsx,jpex=jpex,kpsx=kpsx,kpex=kpex )
1642
1643	grid%muts = grid%mut + grid%mu_2 ! reset muts using filtered mu_2
1644
1645	END IF
1646
1647	!-----------------------------------------------------------
1648	! end acoustic integration polar filter for smallstep mu, t
1649	!-----------------------------------------------------------
1650
1651	BENCH_START(spec_bdy_t_tim)
1652
1653	!$OMP PARALLEL DO &
1654	!$OMP PRIVATE ( ij )
1655	DO ij = 1 , grid%num_tiles
1656
1657	IF( config_flags%specified .or. config_flags%nested ) THEN
1658
1659	CALL spec_bdyupdate(grid%t_2, t_tend, dts_rk, &
1660	't' , config_flags, &
1661	grid%spec_zone, &
1662	ids,ide, jds,jde, kds,kde, &
1663	ims,ime, jms,jme, kms,kme, &
1664	ips,ipe, jps,jpe, kps,kpe, &
1665	grid%i_start(ij), grid%i_end(ij),&
1666	grid%j_start(ij), grid%j_end(ij),&
1667	k_start , k_end )
1668
1669	CALL spec_bdyupdate(grid%mu_2, mu_tend, dts_rk, &
1670	'm' , config_flags, &
1671	grid%spec_zone, &
1672	ids,ide, jds,jde, 1 ,1 , &
1673	ims,ime, jms,jme, 1 ,1 , &
1674	ips,ipe, jps,jpe, 1 ,1 , &
1675	grid%i_start(ij), grid%i_end(ij),&
1676	grid%j_start(ij), grid%j_end(ij),&
1677	1 , 1 )
1678
1679	CALL spec_bdyupdate(grid%muts, mu_tend, dts_rk, &
1680	'm' , config_flags, &
1681	grid%spec_zone, &
1682	ids,ide, jds,jde, 1 ,1 , & ! domain dims
1683	ims,ime, jms,jme, 1 ,1 , & ! memory dims
1684	ips,ipe, jps,jpe, 1 ,1 , & ! patch dims
1685	grid%i_start(ij), grid%i_end(ij), &
1686	grid%j_start(ij), grid%j_end(ij), &
1687	1 , 1 )
1688	ENDIF
1689	BENCH_END(spec_bdy_t_tim)
1690
1691	! small (acoustic) step for the vertical momentum,
1692	! density and coupled potential temperature.
1693
1694
1695	BENCH_START(advance_w_tim)
1696	IF ( config_flags%non_hydrostatic ) THEN
1697	CALL advance_w( grid%w_2, rw_tend, grid%ww, w_save, &
1698	grid%u_2, grid%v_2, &
1699	grid%mu_2, grid%mut, muave, grid%muts, &
1700	t_2save, grid%t_2, grid%t_save, &
1701	grid%ph_2, ph_save, grid%phb, ph_tend, &
1702	grid%ht, c2a, cqw, grid%alt, grid%alb, &
1703	a, alpha, gamma, &
1704	grid%rdx, grid%rdy, dts_rk, t0, grid%epssm, &
1705	grid%dnw, grid%fnm, grid%fnp, grid%rdnw, &
1706	grid%rdn, grid%cf1, grid%cf2, grid%cf3, &
1707	grid%msftx, grid%msfty, &
1708	config_flags, config_flags%top_lid, &
1709	ids,ide, jds,jde, kds,kde, &
1710	ims,ime, jms,jme, kms,kme, &
1711	grid%i_start(ij), grid%i_end(ij), &
1712	grid%j_start(ij), grid%j_end(ij), &
1713	k_start , k_end )
1714	ENDIF
1715	BENCH_END(advance_w_tim)
1716
1717	ENDDO
1718	!$OMP END PARALLEL DO
1719	#ifdef LMDZ1
1720	WRITE(message, *)' dyn_em: after advance_w'
1721	CALL wrf_debug(200, message)
1722	WRITE(message, *)' t_tend: ',t_tendf(im2,km2,jm2), &
1723	' u_tend: ', ru_tendf(im2,1,jm2)
1724	CALL wrf_debug(200, message)
1725	WRITE(message,*)' psfc_tend: ',grid%dpsdt(im2,jm2), &
1726	'p sfc: ',p8w(im2,kms,jm2)
1727	CALL wrf_debug(200, message)
1728	WRITE(message,*)' p 1: ',grid%p(im2,kms,jm2), ' ph 1: ',grid%ph_2(im2,kms,jm2)
1729	CALL wrf_debug(200, message)
1730	WRITE(message,*)' al: ',grid%al(im2,km2,jm2), ' p: ', grid%p(im2,1,jm2), &
1731	' ph: ',grid%ph_2(im2,1,jm2)
1732	CALL wrf_debug(200, message)
1733	WRITE(message,*)' mu: ',grid%mu_2(im2,jm2), 'alt: ', grid%alt(im2,1,jm2), &
1734	' mu: ',grid%mu_2(im2,jm2), ' znu: ', grid%znu(1),' ph 1: ',grid%ph_2(im2,2,jm2)
1735	CALL wrf_debug(200, message)
1736	PRINT *,' c2a: ',c2a(im2,1,jm2), 't_2: ', grid%t_2(im2,1,jm2), &
1737	' t_save: ',grid%t_save(im2,1,jm2),' pm1: ',pm1(im2,1,jm2)
1738	CALL wrf_debug(200, message)
1739	#endif
1740
1741	!-----------------------------------------------------------
1742	! acoustic integration polar filter for smallstep w, geopotential
1743	!-----------------------------------------------------------
1744
1745	IF ( (config_flags%polar) .AND. (config_flags%non_hydrostatic) ) THEN
1746
1747	CALL pxft ( grid=grid &
1748	,lineno=__LINE__ &
1749	,flag_uv = 0 &
1750	,flag_rurv = 0 &
1751	,flag_wph = 1 &
1752	,flag_ww = 0 &
1753	,flag_t = 0 &
1754	,flag_mu = 0 &
1755	,flag_mut = 0 &
1756	,flag_moist = 0 &
1757	,flag_chem = 0 &
1758	,flag_tracer = 0 &
1759	,flag_scalar = 0 &
1760	,positive_definite = .FALSE. &
1761	,moist=moist,chem=chem,tracer=tracer,scalar=scalar &
1762	,fft_filter_lat = config_flags%fft_filter_lat &
1763	,dclat = dclat &
1764	,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
1765	,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
1766	,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe &
1767	,imsx=imsx,imex=imex,jmsx=jmsx,jmex=jmex,kmsx=kmsx,kmex=kmex &
1768	,ipsx=ipsx,ipex=ipex,jpsx=jmsx,jpex=jpex,kpsx=kpsx,kpex=kpex )
1769
1770	END IF
1771
1772	!-----------------------------------------------------------
1773	! end acoustic integration polar filter for smallstep w, geopotential
1774	!-----------------------------------------------------------
1775
1776	!$OMP PARALLEL DO &
1777	!$OMP PRIVATE ( ij )
1778	DO ij = 1 , grid%num_tiles
1779
1780	BENCH_START(sumflux_tim)
1781	CALL sumflux ( grid%u_2, grid%v_2, grid%ww, &
1782	grid%u_save, grid%v_save, ww1, &
1783	grid%muu, grid%muv, &
1784	grid%ru_m, grid%rv_m, grid%ww_m, grid%epssm, &
1785	grid%msfux, grid% msfuy, grid%msfvx, &
1786	grid%msfvx_inv, grid%msfvy, &
1787	iteration, number_of_small_timesteps, &
1788	ids, ide, jds, jde, kds, kde, &
1789	ims, ime, jms, jme, kms, kme, &
1790	grid%i_start(ij), grid%i_end(ij), &
1791	grid%j_start(ij), grid%j_end(ij), &
1792	k_start , k_end )
1793	BENCH_END(sumflux_tim)
1794	#ifdef LMDZ1
1795	WRITE(message, *)' dyn_em: after sumflux'
1796	CALL wrf_debug(200, message)
1797	WRITE(message, *)' t_tend: ',t_tendf(im2,km2,jm2), &
1798	' u_tend: ', ru_tendf(im2,1,jm2)
1799	CALL wrf_debug(200, message)
1800	WRITE(message,*)' psfc_tend: ',grid%dpsdt(im2,jm2), &
1801	'p sfc: ',p8w(im2,kms,jm2)
1802	CALL wrf_debug(200, message)
1803	WRITE(message,*)' p 1: ',grid%p(im2,kms,jm2), ' ph 1: ',grid%ph_2(im2,kms,jm2)
1804	CALL wrf_debug(200, message)
1805	WRITE(message,*)' al: ',grid%al(im2,km2,jm2), ' p: ', grid%p(im2,1,jm2), &
1806	' ph: ',grid%ph_2(im2,1,jm2)
1807	CALL wrf_debug(200, message)
1808	WRITE(message,*)' mu: ',grid%mu_2(im2,jm2), 'alt: ', grid%alt(im2,1,jm2), &
1809	' mu: ',grid%mu_2(im2,jm2), ' znu: ', grid%znu(1),' ph 1: ',grid%ph_2(im2,2,jm2)
1810	CALL wrf_debug(200, message)
1811	PRINT *,' c2a: ',c2a(im2,1,jm2), 't_2: ', grid%t_2(im2,1,jm2), &
1812	' t_save: ',grid%t_save(im2,1,jm2),' pm1: ',pm1(im2,1,jm2)
1813	CALL wrf_debug(200, message)
1814	#endif
1815
1816	IF( config_flags%specified .or. config_flags%nested ) THEN
1817
1818	BENCH_START(spec_bdynhyd_tim)
1819	IF (config_flags%non_hydrostatic) THEN
1820	CALL spec_bdyupdate_ph( ph_save, grid%ph_2, ph_tend, &
1821	mu_tend, grid%muts, dts_rk, &
1822	'h' , config_flags, &
1823	grid%spec_zone, &
1824	ids,ide, jds,jde, kds,kde, &
1825	ims,ime, jms,jme, kms,kme, &
1826	ips,ipe, jps,jpe, kps,kpe, &
1827	grid%i_start(ij), grid%i_end(ij),&
1828	grid%j_start(ij), grid%j_end(ij),&
1829	k_start , k_end )
1830	#ifdef LMDZ1
1831	WRITE(message, *)' dyn_em: after spec_bdynhyd_ph'
1832	CALL wrf_debug(200, message)
1833	WRITE(message, *)' t_tend: ',t_tendf(im2,km2,jm2), &
1834	' u_tend: ', ru_tendf(im2,1,jm2)
1835	CALL wrf_debug(200, message)
1836	WRITE(message,*)' psfc_tend: ',grid%dpsdt(im2,jm2), &
1837	'p sfc: ',p8w(im2,kms,jm2)
1838	CALL wrf_debug(200, message)
1839	WRITE(message,*)' p 1: ',grid%p(im2,kms,jm2), ' ph 1: ',grid%ph_2(im2,kms,jm2)
1840	CALL wrf_debug(200, message)
1841	WRITE(message,*)' al: ',grid%al(im2,km2,jm2), ' p: ', grid%p(im2,1,jm2), &
1842	' ph: ',grid%ph_2(im2,1,jm2)
1843	CALL wrf_debug(200, message)
1844	WRITE(message,*)' mu: ',grid%mu_2(im2,jm2), 'alt: ', grid%alt(im2,1,jm2), &
1845	' mu: ',grid%mu_2(im2,jm2), ' znu: ', grid%znu(1),' ph 1: ',grid%ph_2(im2,2,jm2)
1846	CALL wrf_debug(200, message)
1847	PRINT *,' c2a: ',c2a(im2,1,jm2), 't_2: ', grid%t_2(im2,1,jm2), &
1848	' t_save: ',grid%t_save(im2,1,jm2),' pm1: ',pm1(im2,1,jm2)
1849	CALL wrf_debug(200, message)
1850	#endif
1851	IF( config_flags%specified ) THEN
1852	CALL zero_grad_bdy ( grid%w_2, &
1853	'w' , config_flags, &
1854	grid%spec_zone, &
1855	ids,ide, jds,jde, kds,kde, &
1856	ims,ime, jms,jme, kms,kme, &
1857	ips,ipe, jps,jpe, kps,kpe, &
1858	grid%i_start(ij), grid%i_end(ij), &
1859	grid%j_start(ij), grid%j_end(ij), &
1860	k_start , k_end )
1861	ELSE
1862	CALL spec_bdyupdate ( grid%w_2, rw_tend, dts_rk, &
1863	'h' , config_flags, &
1864	grid%spec_zone, &
1865	ids,ide, jds,jde, kds,kde, &
1866	ims,ime, jms,jme, kms,kme, &
1867	ips,ipe, jps,jpe, kps,kpe, &
1868	grid%i_start(ij), grid%i_end(ij),&
1869	grid%j_start(ij), grid%j_end(ij),&
1870	k_start , k_end )
1871	ENDIF
1872	ENDIF
1873	BENCH_END(spec_bdynhyd_tim)
1874	ENDIF
1875	#ifdef LMDZ1
1876	WRITE(message, *)' dyn_em: after spec_bdynhyd_tim'
1877	CALL wrf_debug(200, message)
1878	WRITE(message, *)' t_tend: ',t_tendf(im2,km2,jm2), &
1879	' u_tend: ', ru_tendf(im2,1,jm2)
1880	CALL wrf_debug(200, message)
1881	WRITE(message,*)' psfc_tend: ',grid%dpsdt(im2,jm2), &
1882	'p sfc: ',p8w(im2,kms,jm2)
1883	CALL wrf_debug(200, message)
1884	WRITE(message,*)' p 1: ',grid%p(im2,kms,jm2), ' ph 1: ',grid%ph_2(im2,kms,jm2)
1885	CALL wrf_debug(200, message)
1886	WRITE(message,*)' al: ',grid%al(im2,km2,jm2), ' p: ', grid%p(im2,1,jm2), &
1887	' ph: ',grid%ph_2(im2,1,jm2)
1888	CALL wrf_debug(200, message)
1889	WRITE(message,*)' mu: ',grid%mu_2(im2,jm2), 'alt: ', grid%alt(im2,1,jm2), &
1890	' mu: ',grid%mu_2(im2,jm2), ' znu: ', grid%znu(1),' ph 1: ',grid%ph_2(im2,2,jm2)
1891	CALL wrf_debug(200, message)
1892	PRINT *,' c2a: ',c2a(im2,1,jm2), 't_2: ', grid%t_2(im2,1,jm2), &
1893	' t_save: ',grid%t_save(im2,1,jm2),' pm1: ',pm1(im2,1,jm2)
1894	CALL wrf_debug(200, message)
1895	#endif
1896
1897	BENCH_START(cald_p_rho_tim)
1898	CALL calc_p_rho( grid%al, grid%p, grid%ph_2, &
1899	grid%alt, grid%t_2, grid%t_save, c2a, pm1, &
1900	grid%mu_2, grid%muts, grid%znu, t0, &
1901	grid%rdnw, grid%dnw, grid%smdiv, &
1902	config_flags%non_hydrostatic, iteration, &
1903	ids, ide, jds, jde, kds, kde, &
1904	ims, ime, jms, jme, kms, kme, &
1905	grid%i_start(ij), grid%i_end(ij), &
1906	grid%j_start(ij), grid%j_end(ij), &
1907	k_start , k_end )
1908	BENCH_END(cald_p_rho_tim)
1909	#ifdef LMDZ1
1910	WRITE(message, *)' dyn_em: after calc_p_rho'
1911	CALL wrf_debug(200, message)
1912	WRITE(message, *)' t_tend: ',t_tendf(im2,km2,jm2), &
1913	' u_tend: ', ru_tendf(im2,1,jm2)
1914	CALL wrf_debug(200, message)
1915	WRITE(message,*)' psfc_tend: ',grid%dpsdt(im2,jm2), &
1916	'p sfc: ',p8w(im2,kms,jm2)
1917	CALL wrf_debug(200, message)
1918	WRITE(message,*)' p 1: ',grid%p(im2,kms,jm2), ' ph 1: ',grid%ph_2(im2,kms,jm2)
1919	CALL wrf_debug(200, message)
1920	#endif
1921
1922	ENDDO
1923	!$OMP END PARALLEL DO
1924	#ifdef LMDZ1
1925	WRITE(message, *)' dyn_em: after geopotential'
1926	CALL wrf_debug(200, message)
1927	WRITE(message, *)' t_tend: ',t_tendf(im2,km2,jm2), &
1928	' u_tend: ', ru_tendf(im2,1,jm2)
1929	CALL wrf_debug(200, message)
1930	WRITE(message,*)' psfc_tend: ',grid%dpsdt(im2,jm2), &
1931	'p sfc: ',p8w(im2,kms,jm2)
1932	CALL wrf_debug(200, message)
1933	WRITE(message,*)' p 1: ',grid%p(im2,kms,jm2), ' ph 1: ',grid%ph_2(im2,kms,jm2)
1934	CALL wrf_debug(200, message)
1935	#endif
1936
1937	#ifdef DM_PARALLEL
1938	!
1939	! Stencils for patch communications (WCS, 29 June 2001)
1940	!
1941	! * *
1942	! * + * * + * +
1943	! * *
1944	!
1945	! ph_2 x
1946	! al x
1947	! p x
1948	!
1949	! 2D variables (x,y)
1950	!
1951	! mu_2 x
1952	! muts x
1953	! mudf x
1954
1955	# include "HALO_EM_C2.inc"
1956	# include "PERIOD_BDY_EM_B3.inc"
1957	#endif
1958
1959	BENCH_START(phys_bc_tim)
1960	!$OMP PARALLEL DO &
1961	!$OMP PRIVATE ( ij )
1962	DO ij = 1 , grid%num_tiles
1963
1964	! boundary condition set for next small timestep
1965
1966	CALL set_physical_bc3d( grid%ph_2, 'w', config_flags, &
1967	ids, ide, jds, jde, kds, kde, &
1968	ims, ime, jms, jme, kms, kme, &
1969	ips, ipe, jps, jpe, kps, kpe, &
1970	grid%i_start(ij), grid%i_end(ij), &
1971	grid%j_start(ij), grid%j_end(ij), &
1972	k_start , k_end )
1973
1974	CALL set_physical_bc3d( grid%al, 'p', config_flags, &
1975	ids, ide, jds, jde, kds, kde, &
1976	ims, ime, jms, jme, kms, kme, &
1977	ips, ipe, jps, jpe, kps, kpe, &
1978	grid%i_start(ij), grid%i_end(ij), &
1979	grid%j_start(ij), grid%j_end(ij), &
1980	k_start , k_end )
1981
1982	CALL set_physical_bc3d( grid%p, 'p', config_flags, &
1983	ids, ide, jds, jde, kds, kde, &
1984	ims, ime, jms, jme, kms, kme, &
1985	ips, ipe, jps, jpe, kps, kpe, &
1986	grid%i_start(ij), grid%i_end(ij), &
1987	grid%j_start(ij), grid%j_end(ij), &
1988	k_start , k_end )
1989
1990	CALL set_physical_bc2d( grid%muts, 't', config_flags, &
1991	ids, ide, jds, jde, &
1992	ims, ime, jms, jme, &
1993	ips, ipe, jps, jpe, &
1994	grid%i_start(ij), grid%i_end(ij), &
1995	grid%j_start(ij), grid%j_end(ij) )
1996
1997	CALL set_physical_bc2d( grid%mu_2, 't', config_flags, &
1998	ids, ide, jds, jde, &
1999	ims, ime, jms, jme, &
2000	ips, ipe, jps, jpe, &
2001	grid%i_start(ij), grid%i_end(ij), &
2002	grid%j_start(ij), grid%j_end(ij) )
2003
2004	CALL set_physical_bc2d( grid%mudf, 't', config_flags, &
2005	ids, ide, jds, jde, &
2006	ims, ime, jms, jme, &
2007	ips, ipe, jps, jpe, &
2008	grid%i_start(ij), grid%i_end(ij), &
2009	grid%j_start(ij), grid%j_end(ij) )
2010
2011	END DO
2012	!$OMP END PARALLEL DO
2013	BENCH_END(phys_bc_tim)
2014
2015	END DO small_steps
2016
2017	!$OMP PARALLEL DO &
2018	!$OMP PRIVATE ( ij )
2019	DO ij = 1 , grid%num_tiles
2020
2021	CALL wrf_debug ( 200 , ' call rk_small_finish' )
2022
2023	! change time-perturbation variables back to
2024	! full perturbation variables.
2025	! first get updated mu at u and v points
2026
2027	BENCH_START(calc_mu_uv_tim)
2028	CALL calc_mu_uv_1 ( config_flags, &
2029	grid%muts, muus, muvs, &
2030	ids, ide, jds, jde, kds, kde, &
2031	ims, ime, jms, jme, kms, kme, &
2032	grid%i_start(ij), grid%i_end(ij), &
2033	grid%j_start(ij), grid%j_end(ij), &
2034	k_start , k_end )
2035	BENCH_END(calc_mu_uv_tim)
2036	BENCH_START(small_step_finish_tim)
2037	CALL small_step_finish( grid%u_2, grid%u_1, grid%v_2, grid%v_1, grid%w_2, grid%w_1, &
2038	grid%t_2, grid%t_1, grid%ph_2, grid%ph_1, grid%ww, ww1, &
2039	grid%mu_2, grid%mu_1, &
2040	grid%mut, grid%muts, grid%muu, muus, grid%muv, muvs, &
2041	grid%u_save, grid%v_save, w_save, &
2042	grid%t_save, ph_save, mu_save, &
2043	grid%msfux,grid%msfuy, grid%msfvx,grid%msfvy, grid%msftx,grid%msfty, &
2044	grid%h_diabatic, &
2045	number_of_small_timesteps,dts_rk, &
2046	rk_step, rk_order, &
2047	ids, ide, jds, jde, kds, kde, &
2048	ims, ime, jms, jme, kms, kme, &
2049	grid%i_start(ij), grid%i_end(ij), &
2050	grid%j_start(ij), grid%j_end(ij), &
2051	k_start , k_end )
2052	! call to set ru_m, rv_m and ww_m b.c's for PD advection
2053
2054	IF (rk_step == rk_order) THEN
2055
2056	CALL set_physical_bc3d( grid%ru_m, 'u', config_flags, &
2057	ids, ide, jds, jde, kds, kde, &
2058	ims, ime, jms, jme, kms, kme, &
2059	ips, ipe, jps, jpe, kps, kpe, &
2060	grid%i_start(ij), grid%i_end(ij), &
2061	grid%j_start(ij), grid%j_end(ij), &
2062	k_start , k_end )
2063
2064	CALL set_physical_bc3d( grid%rv_m, 'v', config_flags, &
2065	ids, ide, jds, jde, kds, kde, &
2066	ims, ime, jms, jme, kms, kme, &
2067	ips, ipe, jps, jpe, kps, kpe, &
2068	grid%i_start(ij), grid%i_end(ij), &
2069	grid%j_start(ij), grid%j_end(ij), &
2070	k_start , k_end )
2071
2072	CALL set_physical_bc3d( grid%ww_m, 'w', config_flags, &
2073	ids, ide, jds, jde, kds, kde, &
2074	ims, ime, jms, jme, kms, kme, &
2075	ips, ipe, jps, jpe, kps, kpe, &
2076	grid%i_start(ij), grid%i_end(ij), &
2077	grid%j_start(ij), grid%j_end(ij), &
2078	k_start , k_end )
2079
2080	CALL set_physical_bc2d( grid%mut, 't', config_flags, &
2081	ids, ide, jds, jde, &
2082	ims, ime, jms, jme, &
2083	ips, ipe, jps, jpe, &
2084	grid%i_start(ij), grid%i_end(ij), &
2085	grid%j_start(ij), grid%j_end(ij) )
2086
2087	CALL set_physical_bc2d( grid%muts, 't', config_flags, &
2088	ids, ide, jds, jde, &
2089	ims, ime, jms, jme, &
2090	ips, ipe, jps, jpe, &
2091	grid%i_start(ij), grid%i_end(ij), &
2092	grid%j_start(ij), grid%j_end(ij) )
2093
2094	END IF
2095
2096	BENCH_END(small_step_finish_tim)
2097
2098	END DO
2099	!$OMP END PARALLEL DO
2100	#ifdef LMDZ1
2101	WRITE(message, *)' dyn_em: after rk_small_finish'
2102	CALL wrf_debug(200, message)
2103	WRITE(message, *)' t_tend: ',t_tendf(im2,km2,jm2), &
2104	' u_tend: ', ru_tendf(im2,1,jm2)
2105	CALL wrf_debug(200, message)
2106	WRITE(message,*)' psfc_tend: ',grid%dpsdt(im2,jm2), &
2107	'p sfc: ',p8w(im2,kms,jm2)
2108	CALL wrf_debug(200, message)
2109	WRITE(message,*)' p 1: ',grid%p(im2,kms,jm2), ' ph 1: ',grid%ph_2(im2,kms,jm2)
2110	CALL wrf_debug(200, message)
2111	#endif
2112
2113	!-----------------------------------------------------------
2114	! polar filter for full dynamics variables and time-averaged mass fluxes
2115	!-----------------------------------------------------------
2116
2117	IF (config_flags%polar) THEN
2118
2119	CALL pxft ( grid=grid &
2120	,lineno=__LINE__ &
2121	,flag_uv = 1 &
2122	,flag_rurv = 1 &
2123	,flag_wph = 1 &
2124	,flag_ww = 1 &
2125	,flag_t = 1 &
2126	,flag_mu = 1 &
2127	,flag_mut = 1 &
2128	,flag_moist = 0 &
2129	,flag_chem = 0 &
2130	,flag_tracer = 0 &
2131	,flag_scalar = 0 &
2132	,positive_definite = .FALSE. &
2133	,moist=moist,chem=chem,tracer=tracer,scalar=scalar &
2134	,fft_filter_lat = config_flags%fft_filter_lat &
2135	,dclat = dclat &
2136	,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
2137	,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
2138	,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe &
2139	,imsx=imsx,imex=imex,jmsx=jmsx,jmex=jmex,kmsx=kmsx,kmex=kmex &
2140	,ipsx=ipsx,ipex=ipex,jpsx=jmsx,jpex=jpex,kpsx=kpsx,kpex=kpex )
2141
2142	END IF
2143
2144	!-----------------------------------------------------------
2145	! end polar filter for full dynamics variables and time-averaged mass fluxes
2146	!-----------------------------------------------------------
2147
2148	!-----------------------------------------------------------------------
2149	! add in physics tendency first if positive definite advection is used.
2150	! pd advection applies advective flux limiter on last runge-kutta step
2151	!-----------------------------------------------------------------------
2152	! first moisture
2153
2154	IF ((config_flags%moist_adv_opt /= ORIGINAL) .and. (rk_step == rk_order)) THEN
2155
2156	!$OMP PARALLEL DO &
2157	!$OMP PRIVATE ( ij )
2158	DO ij = 1 , grid%num_tiles
2159	CALL wrf_debug ( 200 , ' call rk_update_scalar_pd' )
2160	DO im = PARAM_FIRST_SCALAR, num_3d_m
2161	CALL rk_update_scalar_pd( im, im, &
2162	moist_old(ims,kms,jms,im), &
2163	moist_tend(ims,kms,jms,im), &
2164	grid%mu_1, grid%mu_1, grid%mub, &
2165	rk_step, dt_rk, grid%spec_zone, &
2166	config_flags, &
2167	ids, ide, jds, jde, kds, kde, &
2168	ims, ime, jms, jme, kms, kme, &
2169	grid%i_start(ij), grid%i_end(ij), &
2170	grid%j_start(ij), grid%j_end(ij), &
2171	k_start , k_end )
2172	ENDDO
2173	END DO
2174	!$OMP END PARALLEL DO
2175
2176	!---------------------- positive definite bc call
2177	#ifdef DM_PARALLEL
2178	IF (config_flags%moist_adv_opt /= ORIGINAL) THEN
2179	IF ( config_flags%h_sca_adv_order <= 4 ) THEN
2180	# include "HALO_EM_MOIST_OLD_E_5.inc"
2181	ELSE IF ( config_flags%h_sca_adv_order <= 6 ) THEN
2182	# include "HALO_EM_MOIST_OLD_E_7.inc"
2183	ELSE
2184	WRITE(wrf_err_message,*)'solve_em: invalid h_sca_adv_order = ',config_flags%h_sca_adv_order
2185	CALL wrf_error_fatal(TRIM(wrf_err_message))
2186	ENDIF
2187	ENDIF
2188	#endif
2189
2190	#ifdef DM_PARALLEL
2191	# include "PERIOD_BDY_EM_MOIST_OLD.inc"
2192	#endif
2193
2194	!$OMP PARALLEL DO &
2195	!$OMP PRIVATE ( ij )
2196	DO ij = 1 , grid%num_tiles
2197	IF (num_3d_m >= PARAM_FIRST_SCALAR) THEN
2198	DO im = PARAM_FIRST_SCALAR , num_3d_m
2199	CALL set_physical_bc3d( moist_old(ims,kms,jms,im), 'p', config_flags, &
2200	ids, ide, jds, jde, kds, kde, &
2201	ims, ime, jms, jme, kms, kme, &
2202	ips, ipe, jps, jpe, kps, kpe, &
2203	grid%i_start(ij), grid%i_end(ij), &
2204	grid%j_start(ij), grid%j_end(ij), &
2205	k_start , k_end )
2206	END DO
2207	ENDIF
2208	END DO
2209	!$OMP END PARALLEL DO
2210
2211	END IF ! end if for moist_adv_opt
2212
2213	! scalars
2214
2215	IF ((config_flags%scalar_adv_opt /= ORIGINAL) .and. (rk_step == rk_order)) THEN
2216
2217	!$OMP PARALLEL DO &
2218	!$OMP PRIVATE ( ij )
2219	DO ij = 1 , grid%num_tiles
2220	CALL wrf_debug ( 200 , ' call rk_update_scalar_pd' )
2221	DO im = PARAM_FIRST_SCALAR, num_3d_s
2222	CALL rk_update_scalar_pd( im, im, &
2223	scalar_old(ims,kms,jms,im), &
2224	scalar_tend(ims,kms,jms,im), &
2225	grid%mu_1, grid%mu_1, grid%mub, &
2226	rk_step, dt_rk, grid%spec_zone, &
2227	config_flags, &
2228	ids, ide, jds, jde, kds, kde, &
2229	ims, ime, jms, jme, kms, kme, &
2230	grid%i_start(ij), grid%i_end(ij), &
2231	grid%j_start(ij), grid%j_end(ij), &
2232	k_start , k_end )
2233	ENDDO
2234	ENDDO
2235	!$OMP END PARALLEL DO
2236
2237	!---------------------- positive definite bc call
2238	#ifdef DM_PARALLEL
2239	IF (config_flags%scalar_adv_opt /= ORIGINAL) THEN
2240	#ifndef RSL
2241	IF ( config_flags%h_sca_adv_order <= 4 ) THEN
2242	# include "HALO_EM_SCALAR_OLD_E_5.inc"
2243	ELSE IF ( config_flags%h_sca_adv_order <= 6 ) THEN
2244	# include "HALO_EM_SCALAR_OLD_E_7.inc"
2245	ELSE
2246	WRITE(wrf_err_message,*)'solve_em: invalid h_sca_adv_order = ',config_flags%h_sca_adv_order
2247	CALL wrf_error_fatal(TRIM(wrf_err_message))
2248	ENDIF
2249	#else
2250	WRITE(wrf_err_message,*)'cannot use pd scheme with RSL - use RSL-LITE'
2251	CALL wrf_error_fatal(TRIM(wrf_err_message))
2252	#endif
2253	endif
2254	#endif
2255
2256	#ifdef DM_PARALLEL
2257	# include "PERIOD_BDY_EM_SCALAR_OLD.inc"
2258	#endif
2259	!$OMP PARALLEL DO &
2260	!$OMP PRIVATE ( ij )
2261
2262	DO ij = 1 , grid%num_tiles
2263	IF (num_3d_s >= PARAM_FIRST_SCALAR) THEN
2264	DO im = PARAM_FIRST_SCALAR , num_3d_s
2265	CALL set_physical_bc3d( scalar_old(ims,kms,jms,im), 'p', config_flags, &
2266	ids, ide, jds, jde, kds, kde, &
2267	ims, ime, jms, jme, kms, kme, &
2268	ips, ipe, jps, jpe, kps, kpe, &
2269	grid%i_start(ij), grid%i_end(ij), &
2270	grid%j_start(ij), grid%j_end(ij), &
2271	k_start , k_end )
2272	END DO
2273	ENDIF
2274	END DO
2275	!$OMP END PARALLEL DO
2276
2277	END IF ! end if for scalar_adv_opt
2278
2279	! chem
2280
2281	IF ((config_flags%chem_adv_opt /= ORIGINAL) .and. (rk_step == rk_order)) THEN
2282
2283	!$OMP PARALLEL DO &
2284	!$OMP PRIVATE ( ij )
2285	DO ij = 1 , grid%num_tiles
2286	CALL wrf_debug ( 200 , ' call rk_update_scalar_pd' )
2287	DO im = PARAM_FIRST_SCALAR, num_3d_c
2288	CALL rk_update_scalar_pd( im, im, &
2289	chem_old(ims,kms,jms,im), &
2290	chem_tend(ims,kms,jms,im), &
2291	grid%mu_1, grid%mu_1, grid%mub, &
2292	rk_step, dt_rk, grid%spec_zone, &
2293	config_flags, &
2294	ids, ide, jds, jde, kds, kde, &
2295	ims, ime, jms, jme, kms, kme, &
2296	grid%i_start(ij), grid%i_end(ij), &
2297	grid%j_start(ij), grid%j_end(ij), &
2298	k_start , k_end )
2299	ENDDO
2300	END DO
2301	!$OMP END PARALLEL DO
2302
2303	!---------------------- positive definite bc call
2304	#ifdef DM_PARALLEL
2305	IF (config_flags%chem_adv_opt /= ORIGINAL) THEN
2306	IF ( config_flags%h_sca_adv_order <= 4 ) THEN
2307	# include "HALO_EM_CHEM_OLD_E_5.inc"
2308	ELSE IF ( config_flags%h_sca_adv_order <= 6 ) THEN
2309	# include "HALO_EM_CHEM_OLD_E_7.inc"
2310	ELSE
2311	WRITE(wrf_err_message,*)'solve_em: invalid h_sca_adv_order = ',config_flags%h_sca_adv_order
2312	CALL wrf_error_fatal(TRIM(wrf_err_message))
2313	ENDIF
2314	ENDIF
2315	#endif
2316
2317	#ifdef DM_PARALLEL
2318	# include "PERIOD_BDY_EM_CHEM_OLD.inc"
2319	#endif
2320
2321	!$OMP PARALLEL DO &
2322	!$OMP PRIVATE ( ij )
2323	DO ij = 1 , grid%num_tiles
2324	IF (num_3d_c >= PARAM_FIRST_SCALAR) THEN
2325	DO im = PARAM_FIRST_SCALAR , num_3d_c
2326	CALL set_physical_bc3d( chem_old(ims,kms,jms,im), 'p', config_flags, &
2327	ids, ide, jds, jde, kds, kde, &
2328	ims, ime, jms, jme, kms, kme, &
2329	ips, ipe, jps, jpe, kps, kpe, &
2330	grid%i_start(ij), grid%i_end(ij), &
2331	grid%j_start(ij), grid%j_end(ij), &
2332	k_start , k_end )
2333	END DO
2334	ENDIF
2335	END DO
2336	!$OMP END PARALLEL DO
2337
2338	ENDIF ! end if for chem_adv_opt
2339
2340	! tracer
2341
2342	IF ((config_flags%tracer_adv_opt /= ORIGINAL) .and. (rk_step == rk_order)) THEN
2343
2344	!$OMP PARALLEL DO &
2345	!$OMP PRIVATE ( ij )
2346	DO ij = 1 , grid%num_tiles
2347	CALL wrf_debug ( 200 , ' call rk_update_scalar_pd' )
2348	DO im = PARAM_FIRST_SCALAR, num_tracer
2349	CALL rk_update_scalar_pd( im, im, &
2350	tracer_old(ims,kms,jms,im), &
2351	tracer_tend(ims,kms,jms,im), &
2352	grid%mu_1, grid%mu_1, grid%mub, &
2353	rk_step, dt_rk, grid%spec_zone, &
2354	config_flags, &
2355	ids, ide, jds, jde, kds, kde, &
2356	ims, ime, jms, jme, kms, kme, &
2357	grid%i_start(ij), grid%i_end(ij), &
2358	grid%j_start(ij), grid%j_end(ij), &
2359	k_start , k_end )
2360	ENDDO
2361	END DO
2362	!$OMP END PARALLEL DO
2363
2364	!---------------------- positive definite bc call
2365	#ifdef DM_PARALLEL
2366	IF (config_flags%tracer_adv_opt /= ORIGINAL) THEN
2367	IF ( config_flags%h_sca_adv_order <= 4 ) THEN
2368	# include "HALO_EM_TRACER_OLD_E_5.inc"
2369	ELSE IF ( config_flags%h_sca_adv_order <= 6 ) THEN
2370	# include "HALO_EM_TRACER_OLD_E_7.inc"
2371	ELSE
2372	WRITE(wrf_err_message,*)'solve_em: invalid h_sca_adv_order = ',config_flags%h_sca_adv_order
2373	CALL wrf_error_fatal(TRIM(wrf_err_message))
2374	ENDIF
2375	ENDIF
2376	#endif
2377
2378	#ifdef DM_PARALLEL
2379	# include "PERIOD_BDY_EM_TRACER_OLD.inc"
2380	#endif
2381
2382	!$OMP PARALLEL DO &
2383	!$OMP PRIVATE ( ij )
2384	DO ij = 1 , grid%num_tiles
2385	IF (num_tracer >= PARAM_FIRST_SCALAR) THEN
2386	DO im = PARAM_FIRST_SCALAR , num_tracer
2387	CALL set_physical_bc3d( tracer_old(ims,kms,jms,im), 'p', config_flags, &
2388	ids, ide, jds, jde, kds, kde, &
2389	ims, ime, jms, jme, kms, kme, &
2390	ips, ipe, jps, jpe, kps, kpe, &
2391	grid%i_start(ij), grid%i_end(ij), &
2392	grid%j_start(ij), grid%j_end(ij), &
2393	k_start , k_end )
2394	END DO
2395	ENDIF
2396	END DO
2397	!$OMP END PARALLEL DO
2398
2399	ENDIF ! end if for tracer_adv_opt
2400
2401	! tke
2402
2403	IF ((config_flags%tke_adv_opt /= ORIGINAL) .and. (rk_step == rk_order) &
2404	.and. (config_flags%km_opt .eq. 2) ) THEN
2405
2406	!$OMP PARALLEL DO &
2407	!$OMP PRIVATE ( ij )
2408	DO ij = 1 , grid%num_tiles
2409	CALL wrf_debug ( 200 , ' call rk_update_scalar_pd' )
2410	CALL rk_update_scalar_pd( 1, 1, &
2411	grid%tke_1, &
2412	tke_tend(ims,kms,jms), &
2413	grid%mu_1, grid%mu_1, grid%mub, &
2414	rk_step, dt_rk, grid%spec_zone, &
2415	config_flags, &
2416	ids, ide, jds, jde, kds, kde, &
2417	ims, ime, jms, jme, kms, kme, &
2418	grid%i_start(ij), grid%i_end(ij), &
2419	grid%j_start(ij), grid%j_end(ij), &
2420	k_start , k_end )
2421	ENDDO
2422	!$OMP END PARALLEL DO
2423
2424	!---------------------- positive definite bc call
2425	#ifdef DM_PARALLEL
2426	IF (config_flags%tke_adv_opt /= ORIGINAL) THEN
2427	IF ( config_flags%h_sca_adv_order <= 4 ) THEN
2428	# include "HALO_EM_TKE_OLD_E_5.inc"
2429	ELSE IF ( config_flags%h_sca_adv_order <= 6 ) THEN
2430	# include "HALO_EM_TKE_OLD_E_7.inc"
2431	ELSE
2432	WRITE(wrf_err_message,*)'solve_em: invalid h_sca_adv_order = ',config_flags%h_sca_adv_order
2433	CALL wrf_error_fatal(TRIM(wrf_err_message))
2434	ENDIF
2435	ENDIF
2436	#endif
2437
2438	#ifdef DM_PARALLEL
2439	# include "PERIOD_BDY_EM_TKE_OLD.inc"
2440	#endif
2441
2442	!$OMP PARALLEL DO &
2443	!$OMP PRIVATE ( ij )
2444	DO ij = 1 , grid%num_tiles
2445	CALL set_physical_bc3d( grid%tke_1, 'p', config_flags, &
2446	ids, ide, jds, jde, kds, kde, &
2447	ims, ime, jms, jme, kms, kme, &
2448	ips, ipe, jps, jpe, kps, kpe, &
2449	grid%i_start(ij), grid%i_end(ij), &
2450	grid%j_start(ij), grid%j_end(ij), &
2451	k_start , k_end )
2452	END DO
2453	!$OMP END PARALLEL DO
2454
2455	!--- end of positive definite physics tendency update
2456
2457	END IF ! end if for tke_adv_opt
2458
2459	#ifdef DM_PARALLEL
2460	!
2461	! Stencils for patch communications (WCS, 29 June 2001)
2462	!
2463	! * * * * *
2464	! * * * * *
2465	! * * + * *
2466	! * * * * *
2467	! * * * * *
2468	!
2469	! ru_m x
2470	! rv_m x
2471	! ww_m x
2472	! mut x
2473	!
2474	!--------------------------------------------------------------
2475
2476	# include "HALO_EM_D.inc"
2477	! WCS addition 11/19/08
2478	# include "PERIOD_EM_DA.inc"
2479	#endif
2480
2481	!<DESCRIPTION>
2482	!<pre>
2483	! (4) Still within the RK loop, the scalar variables are advanced.
2484	!
2485	! For the moist and chem variables, each one is advanced
2486	! individually, using named loops "moist_variable_loop:"
2487	! and "chem_variable_loop:". Each RK substep begins by
2488	! calculating the advective tendency, and, for the first RK step,
2489	! 3D mixing (calling rk_scalar_tend) followed by an update
2490	! of the scalar (calling rk_update_scalar).
2491	!</pre>
2492	!</DESCRIPTION>
2493
2494
2495	moist_scalar_advance: IF (num_3d_m >= PARAM_FIRST_SCALAR ) THEN
2496
2497	moist_variable_loop: DO im = PARAM_FIRST_SCALAR, num_3d_m
2498
2499	! adv_moist_cond is set in module_physics_init based on mp_physics choice
2500	! true except for Ferrier scheme
2501
2502	IF (grid%adv_moist_cond .or. im==p_qv ) THEN
2503
2504	!$OMP PARALLEL DO &
2505	!$OMP PRIVATE ( ij )
2506	moist_tile_loop_1: DO ij = 1 , grid%num_tiles
2507
2508	CALL wrf_debug ( 200 , ' call rk_scalar_tend' )
2509	tenddec = .false.
2510
2511	BENCH_START(rk_scalar_tend_tim)
2512	CALL rk_scalar_tend ( im, im, config_flags, tenddec, &
2513	rk_step, dt_rk, &
2514	grid%ru_m, grid%rv_m, grid%ww_m, &
2515	grid%muts, grid%mub, grid%mu_1, &
2516	grid%alt, &
2517	moist_old(ims,kms,jms,im), &
2518	moist(ims,kms,jms,im), &
2519	moist_tend(ims,kms,jms,im), &
2520	advect_tend,h_tendency,z_tendency,grid%rqvften, &
2521	grid%qv_base, .true., grid%fnm, grid%fnp, &
2522	grid%msfux,grid%msfuy, grid%msfvx, grid%msfvx_inv,&
2523	grid%msfvy, grid%msftx,grid%msfty, &
2524	grid%rdx, grid%rdy, grid%rdn, grid%rdnw, grid%khdif, &
2525	grid%kvdif, grid%xkhh, &
2526	grid%diff_6th_opt, grid%diff_6th_factor, &
2527	config_flags%moist_adv_opt, &
2528	ids, ide, jds, jde, kds, kde, &
2529	ims, ime, jms, jme, kms, kme, &
2530	grid%i_start(ij), grid%i_end(ij), &
2531	grid%j_start(ij), grid%j_end(ij), &
2532	k_start , k_end )
2533
2534	BENCH_END(rk_scalar_tend_tim)
2535
2536	BENCH_START(rlx_bdy_scalar_tim)
2537	IF( ( config_flags%specified .or. config_flags%nested ) .and. rk_step == 1 ) THEN
2538	IF ( im .EQ. P_QV .OR. config_flags%nested ) THEN
2539	CALL relax_bdy_scalar ( moist_tend(ims,kms,jms,im), &
2540	moist(ims,kms,jms,im), grid%mut, &
2541	moist_bxs(jms,kms,1,im),moist_bxe(jms,kms,1,im), &
2542	moist_bys(ims,kms,1,im),moist_bye(ims,kms,1,im), &
2543	moist_btxs(jms,kms,1,im),moist_btxe(jms,kms,1,im), &
2544	moist_btys(ims,kms,1,im),moist_btye(ims,kms,1,im), &
2545	config_flags%spec_bdy_width, grid%spec_zone, grid%relax_zone, &
2546	grid%dtbc, grid%fcx, grid%gcx, &
2547	config_flags, &
2548	ids,ide, jds,jde, kds,kde, & ! domain dims
2549	ims,ime, jms,jme, kms,kme, & ! memory dims
2550	ips,ipe, jps,jpe, kps,kpe, & ! patch dims
2551	grid%i_start(ij), grid%i_end(ij), &
2552	grid%j_start(ij), grid%j_end(ij), &
2553	k_start, k_end )
2554
2555	CALL spec_bdy_scalar ( moist_tend(ims,kms,jms,im), &
2556	moist_bxs(jms,kms,1,im),moist_bxe(jms,kms,1,im), &
2557	moist_bys(ims,kms,1,im),moist_bye(ims,kms,1,im), &
2558	moist_btxs(jms,kms,1,im),moist_btxe(jms,kms,1,im), &
2559	moist_btys(ims,kms,1,im),moist_btye(ims,kms,1,im), &
2560	config_flags%spec_bdy_width, grid%spec_zone, &
2561	config_flags, &
2562	ids,ide, jds,jde, kds,kde, & ! domain dims
2563	ims,ime, jms,jme, kms,kme, & ! memory dims
2564	ips,ipe, jps,jpe, kps,kpe, & ! patch dims
2565	grid%i_start(ij), grid%i_end(ij), &
2566	grid%j_start(ij), grid%j_end(ij), &
2567	k_start, k_end )
2568	ENDIF
2569	ENDIF
2570	BENCH_END(rlx_bdy_scalar_tim)
2571
2572	ENDDO moist_tile_loop_1
2573	!$OMP END PARALLEL DO
2574
2575	!$OMP PARALLEL DO &
2576	!$OMP PRIVATE ( ij )
2577	moist_tile_loop_2: DO ij = 1 , grid%num_tiles
2578
2579	CALL wrf_debug ( 200 , ' call rk_update_scalar' )
2580	tenddec = .false.
2581
2582	BENCH_START(update_scal_tim)
2583	CALL rk_update_scalar( scs=im, sce=im, &
2584	scalar_1=moist_old(ims,kms,jms,im), &
2585	scalar_2=moist(ims,kms,jms,im), &
2586	sc_tend=moist_tend(ims,kms,jms,im), &
2587	! advh_t=advh_t(ims,kms,jms,1), &
2588	! advz_t=advz_t(ims,kms,jms,1), &
2589	advect_tend=advect_tend, &
2590	h_tendency=h_tendency, z_tendency=z_tendency, &
2591	msftx=grid%msftx,msfty=grid%msfty, &
2592	mu_old=grid%mu_1, mu_new=grid%mu_2, mu_base=grid%mub, &
2593	rk_step=rk_step, dt=dt_rk, spec_zone=grid%spec_zone, &
2594	config_flags=config_flags, tenddec=tenddec, &
2595	ids=ids, ide=ide, jds=jds, jde=jde, kds=kds, kde=kde, &
2596	ims=ims, ime=ime, jms=jms, jme=jme, kms=kms, kme=kme, &
2597	its=grid%i_start(ij), ite=grid%i_end(ij), &
2598	jts=grid%j_start(ij), jte=grid%j_end(ij), &
2599	kts=k_start , kte=k_end )
2600	BENCH_END(update_scal_tim)
2601
2602	BENCH_START(flow_depbdy_tim)
2603	IF( config_flags%specified ) THEN
2604	IF(im .ne. P_QV)THEN
2605	CALL flow_dep_bdy ( moist(ims,kms,jms,im), &
2606	grid%ru_m, grid%rv_m, config_flags, &
2607	grid%spec_zone, &
2608	ids,ide, jds,jde, kds,kde, &
2609	ims,ime, jms,jme, kms,kme, &
2610	ips,ipe, jps,jpe, kps,kpe, &
2611	grid%i_start(ij), grid%i_end(ij), &
2612	grid%j_start(ij), grid%j_end(ij), &
2613	k_start, k_end )
2614	ENDIF
2615	ENDIF
2616	BENCH_END(flow_depbdy_tim)
2617
2618	ENDDO moist_tile_loop_2
2619	!$OMP END PARALLEL DO
2620
2621	ENDIF !-- if (grid%adv_moist_cond .or. im==p_qv ) then
2622
2623	ENDDO moist_variable_loop
2624
2625	ENDIF moist_scalar_advance
2626
2627	BENCH_START(tke_adv_tim)
2628	TKE_advance: IF (config_flags%km_opt .eq. 2) then
2629	#ifdef DM_PARALLEL
2630	IF ( config_flags%h_mom_adv_order <= 4 ) THEN
2631	# include "HALO_EM_TKE_ADVECT_3.inc"
2632	ELSE IF ( config_flags%h_mom_adv_order <= 6 ) THEN
2633	# include "HALO_EM_TKE_ADVECT_5.inc"
2634	ELSE
2635	WRITE(wrf_err_message,*)'solve_em: invalid h_mom_adv_order = ',config_flags%h_mom_adv_order
2636	CALL wrf_error_fatal(TRIM(wrf_err_message))
2637	ENDIF
2638	#endif
2639	!$OMP PARALLEL DO &
2640	!$OMP PRIVATE ( ij )
2641	tke_tile_loop_1: DO ij = 1 , grid%num_tiles
2642
2643	CALL wrf_debug ( 200 , ' call rk_scalar_tend for tke' )
2644	tenddec = .false.
2645	CALL rk_scalar_tend ( 1, 1, config_flags, tenddec, &
2646	rk_step, dt_rk, &
2647	grid%ru_m, grid%rv_m, grid%ww_m, &
2648	grid%muts, grid%mub, grid%mu_1, &
2649	grid%alt, &
2650	grid%tke_1, &
2651	grid%tke_2, &
2652	tke_tend(ims,kms,jms), &
2653	advect_tend,h_tendency,z_tendency,grid%rqvften, &
2654	grid%qv_base, .false., grid%fnm, grid%fnp, &
2655	grid%msfux,grid%msfuy, grid%msfvx, grid%msfvx_inv, &
2656	grid%msfvy, grid%msftx,grid%msfty, &
2657	grid%rdx, grid%rdy, grid%rdn, grid%rdnw, grid%khdif, &
2658	grid%kvdif, grid%xkhh, &
2659	grid%diff_6th_opt, grid%diff_6th_factor, &
2660	config_flags%tke_adv_opt, &
2661	ids, ide, jds, jde, kds, kde, &
2662	ims, ime, jms, jme, kms, kme, &
2663	grid%i_start(ij), grid%i_end(ij), &
2664	grid%j_start(ij), grid%j_end(ij), &
2665	k_start , k_end )
2666
2667	ENDDO tke_tile_loop_1
2668	!$OMP END PARALLEL DO
2669
2670	!$OMP PARALLEL DO &
2671	!$OMP PRIVATE ( ij )
2672	tke_tile_loop_2: DO ij = 1 , grid%num_tiles
2673
2674	CALL wrf_debug ( 200 , ' call rk_update_scalar' )
2675	tenddec = .false.
2676	CALL rk_update_scalar( scs=1, sce=1, &
2677	scalar_1=grid%tke_1, &
2678	scalar_2=grid%tke_2, &
2679	sc_tend=tke_tend(ims,kms,jms), &
2680	! advh_t=advh_t(ims,kms,jms,1), &
2681	! advz_t=advz_t(ims,kms,jms,1), &
2682	advect_tend=advect_tend, &
2683	h_tendency=h_tendency, z_tendency=z_tendency, &
2684	msftx=grid%msftx,msfty=grid%msfty, &
2685	mu_old=grid%mu_1, mu_new=grid%mu_2, mu_base=grid%mub, &
2686	rk_step=rk_step, dt=dt_rk, spec_zone=grid%spec_zone, &
2687	config_flags=config_flags, tenddec=tenddec, &
2688	ids=ids, ide=ide, jds=jds, jde=jde, kds=kds, kde=kde, &
2689	ims=ims, ime=ime, jms=jms, jme=jme, kms=kms, kme=kme, &
2690	its=grid%i_start(ij), ite=grid%i_end(ij), &
2691	jts=grid%j_start(ij), jte=grid%j_end(ij), &
2692	kts=k_start , kte=k_end )
2693
2694	! bound the tke (greater than 0, less than tke_upper_bound)
2695
2696	CALL bound_tke( grid%tke_2, grid%tke_upper_bound, &
2697	ids, ide, jds, jde, kds, kde, &
2698	ims, ime, jms, jme, kms, kme, &
2699	grid%i_start(ij), grid%i_end(ij), &
2700	grid%j_start(ij), grid%j_end(ij), &
2701	k_start , k_end )
2702
2703	IF( config_flags%specified .or. config_flags%nested ) THEN
2704	CALL flow_dep_bdy ( grid%tke_2, &
2705	grid%ru_m, grid%rv_m, config_flags, &
2706	grid%spec_zone, &
2707	ids,ide, jds,jde, kds,kde, & ! domain dims
2708	ims,ime, jms,jme, kms,kme, & ! memory dims
2709	ips,ipe, jps,jpe, kps,kpe, & ! patch dims
2710	grid%i_start(ij), grid%i_end(ij), &
2711	grid%j_start(ij), grid%j_end(ij), &
2712	k_start, k_end )
2713	ENDIF
2714	ENDDO tke_tile_loop_2
2715	!$OMP END PARALLEL DO
2716
2717	ENDIF TKE_advance
2718	BENCH_END(tke_adv_tim)
2719
2720	#ifdef WRF_CHEM
2721	! next the chemical species
2722	BENCH_START(chem_adv_tim)
2723	chem_scalar_advance: IF (num_3d_c >= PARAM_FIRST_SCALAR) THEN
2724
2725	chem_variable_loop: DO ic = PARAM_FIRST_SCALAR, num_3d_c
2726
2727	!$OMP PARALLEL DO &
2728	!$OMP PRIVATE ( ij )
2729	chem_tile_loop_1: DO ij = 1 , grid%num_tiles
2730
2731	CALL wrf_debug ( 200 , ' call rk_scalar_tend in chem_tile_loop_1' )
2732	tenddec = (( config_flags%chemdiag == USECHEMDIAG ) .and. &
2733	( adv_ct_indices(ic) >= PARAM_FIRST_SCALAR ))
2734	CALL rk_scalar_tend ( ic, ic, config_flags, tenddec, &
2735	rk_step, dt_rk, &
2736	grid%ru_m, grid%rv_m, grid%ww_m, &
2737	grid%muts, grid%mub, grid%mu_1, &
2738	grid%alt, &
2739	chem_old(ims,kms,jms,ic), &
2740	chem(ims,kms,jms,ic), &
2741	chem_tend(ims,kms,jms,ic), &
2742	advect_tend,h_tendency,z_tendency,grid%rqvften, &
2743	grid%qv_base, .false., grid%fnm, grid%fnp, &
2744	grid%msfux,grid%msfuy, grid%msfvx, grid%msfvx_inv, &
2745	grid%msfvy, grid%msftx,grid%msfty, &
2746	grid%rdx, grid%rdy, grid%rdn, grid%rdnw, &
2747	grid%khdif, grid%kvdif, grid%xkhh, &
2748	grid%diff_6th_opt, grid%diff_6th_factor, &
2749	config_flags%chem_adv_opt, &
2750	ids, ide, jds, jde, kds, kde, &
2751	ims, ime, jms, jme, kms, kme, &
2752	grid%i_start(ij), grid%i_end(ij), &
2753	grid%j_start(ij), grid%j_end(ij), &
2754	k_start , k_end )
2755	!
2756	! Currently, chemistry species with specified boundaries (i.e. the mother
2757	! domain) are being over written by flow_dep_bdy_chem. So, relax_bdy and
2758	! spec_bdy are only called for nests. For boundary conditions from global model or larger domain,
2759	! chem is uncoupled, and only used for one row/column on inflow (if have_bcs_chem=.true.)
2760	!
2761	IF( ( config_flags%nested ) .and. rk_step == 1 ) THEN
2762	IF(ic.eq.1)CALL wrf_debug ( 10 , ' have_bcs_chem' )
2763	CALL relax_bdy_scalar ( chem_tend(ims,kms,jms,ic), &
2764	chem(ims,kms,jms,ic), grid%mut, &
2765	chem_bxs(jms,kms,1,ic),chem_bxe(jms,kms,1,ic), &
2766	chem_bys(ims,kms,1,ic),chem_bye(ims,kms,1,ic), &
2767	chem_btxs(jms,kms,1,ic),chem_btxe(jms,kms,1,ic), &
2768	chem_btys(ims,kms,1,ic),chem_btye(ims,kms,1,ic), &
2769	config_flags%spec_bdy_width, grid%spec_zone, grid%relax_zone, &
2770	grid%dtbc, grid%fcx, grid%gcx, &
2771	config_flags, &
2772	ids,ide, jds,jde, kds,kde, &
2773	ims,ime, jms,jme, kms,kme, &
2774	ips,ipe, jps,jpe, kps,kpe, &
2775	grid%i_start(ij), grid%i_end(ij), &
2776	grid%j_start(ij), grid%j_end(ij), &
2777	k_start, k_end )
2778	CALL spec_bdy_scalar ( chem_tend(ims,kms,jms,ic), &
2779	chem_bxs(jms,kms,1,ic),chem_bxe(jms,kms,1,ic), &
2780	chem_bys(ims,kms,1,ic),chem_bye(ims,kms,1,ic), &
2781	chem_btxs(jms,kms,1,ic),chem_btxe(jms,kms,1,ic), &
2782	chem_btys(ims,kms,1,ic),chem_btye(ims,kms,1,ic), &
2783	config_flags%spec_bdy_width, grid%spec_zone, &
2784	config_flags, &
2785	ids,ide, jds,jde, kds,kde, &
2786	ims,ime, jms,jme, kms,kme, &
2787	ips,ipe, jps,jpe, kps,kpe, &
2788	grid%i_start(ij), grid%i_end(ij), &
2789	grid%j_start(ij), grid%j_end(ij), &
2790	k_start, k_end )
2791	ENDIF
2792
2793	ENDDO chem_tile_loop_1
2794	!$OMP END PARALLEL DO
2795
2796	!$OMP PARALLEL DO &
2797	!$OMP PRIVATE ( ij )
2798
2799	chem_tile_loop_2: DO ij = 1 , grid%num_tiles
2800
2801	CALL wrf_debug ( 200 , ' call rk_update_scalar' )
2802	tenddec = (( config_flags%chemdiag == USECHEMDIAG ) .and. &
2803	( adv_ct_indices(ic) >= PARAM_FIRST_SCALAR ))
2804	CALL rk_update_scalar( scs=ic, sce=ic, &
2805	scalar_1=chem_old(ims,kms,jms,ic), &
2806	scalar_2=chem(ims,kms,jms,ic), &
2807	sc_tend=chem_tend(ims,kms,jms,ic), &
2808	advh_t=advh_ct(ims,kms,jms,adv_ct_indices(ic)), &
2809	advz_t=advz_ct(ims,kms,jms,adv_ct_indices(ic)), &
2810	advect_tend=advect_tend, &
2811	h_tendency=h_tendency, z_tendency=z_tendency, &
2812	msftx=grid%msftx,msfty=grid%msfty, &
2813	mu_old=grid%mu_1, mu_new=grid%mu_2, mu_base=grid%mub, &
2814	rk_step=rk_step, dt=dt_rk, spec_zone=grid%spec_zone, &
2815	config_flags=config_flags, tenddec=tenddec, &
2816	ids=ids, ide=ide, jds=jds, jde=jde, kds=kds, kde=kde, &
2817	ims=ims, ime=ime, jms=jms, jme=jme, kms=kms, kme=kme, &
2818	its=grid%i_start(ij), ite=grid%i_end(ij), &
2819	jts=grid%j_start(ij), jte=grid%j_end(ij), &
2820	kts=k_start , kte=k_end )
2821
2822	IF( config_flags%specified ) THEN
2823	CALL flow_dep_bdy_chem( chem(ims,kms,jms,ic), &
2824	chem_bxs(jms,kms,1,ic), chem_btxs(jms,kms,1,ic), &
2825	chem_bxe(jms,kms,1,ic), chem_btxe(jms,kms,1,ic), &
2826	chem_bys(ims,kms,1,ic), chem_btys(ims,kms,1,ic), &
2827	chem_bye(ims,kms,1,ic), chem_btye(ims,kms,1,ic), &
2828	dt_rk+grid%dtbc, &
2829	config_flags%spec_bdy_width,grid%z, &
2830	grid%have_bcs_chem, &
2831	grid%ru_m, grid%rv_m, config_flags,grid%alt, &
2832	grid%t_1,grid%pb,grid%p,t0,p1000mb,rcp,grid%ph_2,grid%phb,g, &
2833	grid%spec_zone,ic, &
2834	ids,ide, jds,jde, kds,kde, & ! domain dims
2835	ims,ime, jms,jme, kms,kme, & ! memory dims
2836	ips,ipe, jps,jpe, kps,kpe, & ! patch dims
2837	grid%i_start(ij), grid%i_end(ij), &
2838	grid%j_start(ij), grid%j_end(ij), &
2839	k_start, k_end )
2840	ENDIF
2841	ENDDO chem_tile_loop_2
2842	!$OMP END PARALLEL DO
2843
2844	ENDDO chem_variable_loop
2845	ENDIF chem_scalar_advance
2846	BENCH_END(chem_adv_tim)
2847	#endif
2848	! next the chemical species
2849	BENCH_START(tracer_adv_tim)
2850	tracer_advance: IF (num_tracer >= PARAM_FIRST_SCALAR) THEN
2851
2852	tracer_variable_loop: DO ic = PARAM_FIRST_SCALAR, num_tracer
2853
2854	!$OMP PARALLEL DO &
2855	!$OMP PRIVATE ( ij )
2856	tracer_tile_loop_1: DO ij = 1 , grid%num_tiles
2857
2858	CALL wrf_debug ( 15 , ' call rk_scalar_tend in tracer_tile_loop_1' )
2859	tenddec = .false.
2860	CALL rk_scalar_tend ( ic, ic, config_flags, tenddec, &
2861	rk_step, dt_rk, &
2862	grid%ru_m, grid%rv_m, grid%ww_m, &
2863	grid%muts, grid%mub, grid%mu_1, &
2864	grid%alt, &
2865	tracer_old(ims,kms,jms,ic), &
2866	tracer(ims,kms,jms,ic), &
2867	tracer_tend(ims,kms,jms,ic), &
2868	advect_tend,h_tendency,z_tendency,grid%rqvften, &
2869	grid%qv_base, .false., grid%fnm, grid%fnp, &
2870	grid%msfux,grid%msfuy, grid%msfvx, grid%msfvx_inv, &
2871	grid%msfvy, grid%msftx,grid%msfty, &
2872	grid%rdx, grid%rdy, grid%rdn, grid%rdnw, &
2873	grid%khdif, grid%kvdif, grid%xkhh, &
2874	grid%diff_6th_opt, grid%diff_6th_factor, &
2875	config_flags%tracer_adv_opt, &
2876	ids, ide, jds, jde, kds, kde, &
2877	ims, ime, jms, jme, kms, kme, &
2878	grid%i_start(ij), grid%i_end(ij), &
2879	grid%j_start(ij), grid%j_end(ij), &
2880	k_start , k_end )
2881	!
2882	! Currently, chemistry species with specified boundaries (i.e. the mother
2883	! domain) are being over written by flow_dep_bdy_chem. So, relax_bdy and
2884	! spec_bdy are only called for nests. For boundary conditions from global model or larger domain,
2885	! chem is uncoupled, and only used for one row/column on inflow (if have_bcs_chem=.true.)
2886	!
2887	IF( ( config_flags%nested ) .and. rk_step == 1 ) THEN
2888	IF(ic.eq.1)CALL wrf_debug ( 10 , ' have_bcs_tracer' )
2889	CALL relax_bdy_scalar ( tracer_tend(ims,kms,jms,ic), &
2890	tracer(ims,kms,jms,ic), grid%mut, &
2891	tracer_bxs(jms,kms,1,ic),tracer_bxe(jms,kms,1,ic), &
2892	tracer_bys(ims,kms,1,ic),tracer_bye(ims,kms,1,ic), &
2893	tracer_btxs(jms,kms,1,ic),tracer_btxe(jms,kms,1,ic), &
2894	tracer_btys(ims,kms,1,ic),tracer_btye(ims,kms,1,ic), &
2895	config_flags%spec_bdy_width, grid%spec_zone, grid%relax_zone, &
2896	grid%dtbc, grid%fcx, grid%gcx, &
2897	config_flags, &
2898	ids,ide, jds,jde, kds,kde, &
2899	ims,ime, jms,jme, kms,kme, &
2900	ips,ipe, jps,jpe, kps,kpe, &
2901	grid%i_start(ij), grid%i_end(ij), &
2902	grid%j_start(ij), grid%j_end(ij), &
2903	k_start, k_end )
2904	CALL spec_bdy_scalar ( tracer_tend(ims,kms,jms,ic), &
2905	tracer_bxs(jms,kms,1,ic),tracer_bxe(jms,kms,1,ic), &
2906	tracer_bys(ims,kms,1,ic),tracer_bye(ims,kms,1,ic), &
2907	tracer_btxs(jms,kms,1,ic),tracer_btxe(jms,kms,1,ic), &
2908	tracer_btys(ims,kms,1,ic),tracer_btye(ims,kms,1,ic), &
2909	config_flags%spec_bdy_width, grid%spec_zone, &
2910	config_flags, &
2911	ids,ide, jds,jde, kds,kde, &
2912	ims,ime, jms,jme, kms,kme, &
2913	ips,ipe, jps,jpe, kps,kpe, &
2914	grid%i_start(ij), grid%i_end(ij), &
2915	grid%j_start(ij), grid%j_end(ij), &
2916	k_start, k_end )
2917	ENDIF
2918
2919	ENDDO tracer_tile_loop_1
2920	!$OMP END PARALLEL DO
2921
2922	!$OMP PARALLEL DO &
2923	!$OMP PRIVATE ( ij )
2924
2925	tracer_tile_loop_2: DO ij = 1 , grid%num_tiles
2926
2927	CALL wrf_debug ( 200 , ' call rk_update_scalar' )
2928	tenddec = .false.
2929	CALL rk_update_scalar( scs=ic, sce=ic, &
2930	scalar_1=tracer_old(ims,kms,jms,ic), &
2931	scalar_2=tracer(ims,kms,jms,ic), &
2932	sc_tend=tracer_tend(ims,kms,jms,ic), &
2933	! advh_t=advh_t(ims,kms,jms,1), &
2934	! advz_t=advz_t(ims,kms,jms,1), &
2935	advect_tend=advect_tend, &
2936	h_tendency=h_tendency, z_tendency=z_tendency, &
2937	msftx=grid%msftx,msfty=grid%msfty, &
2938	mu_old=grid%mu_1, mu_new=grid%mu_2, mu_base=grid%mub, &
2939	rk_step=rk_step, dt=dt_rk, spec_zone=grid%spec_zone, &
2940	config_flags=config_flags, tenddec=tenddec, &
2941	ids=ids, ide=ide, jds=jds, jde=jde, kds=kds, kde=kde, &
2942	ims=ims, ime=ime, jms=jms, jme=jme, kms=kms, kme=kme, &
2943	its=grid%i_start(ij), ite=grid%i_end(ij), &
2944	jts=grid%j_start(ij), jte=grid%j_end(ij), &
2945	kts=k_start , kte=k_end )
2946
2947	IF( config_flags%specified ) THEN
2948	#ifdef WRF_CHEM
2949	CALL flow_dep_bdy_tracer( tracer(ims,kms,jms,ic), &
2950	tracer_bxs(jms,kms,1,ic), tracer_btxs(jms,kms,1,ic), &
2951	tracer_bxe(jms,kms,1,ic), tracer_btxe(jms,kms,1,ic), &
2952	tracer_bys(ims,kms,1,ic), tracer_btys(ims,kms,1,ic), &
2953	tracer_bye(ims,kms,1,ic), tracer_btye(ims,kms,1,ic), &
2954	dt_rk+grid%dtbc, &
2955	config_flags%spec_bdy_width,grid%z, &
2956	grid%have_bcs_tracer, &
2957	grid%ru_m, grid%rv_m, config_flags%tracer_opt,grid%alt, &
2958	grid%t_1,grid%pb,grid%p,t0,p1000mb,rcp,grid%ph_2,grid%phb,g, &
2959	grid%spec_zone,ic, &
2960	ids,ide, jds,jde, kds,kde, & ! domain dims
2961	ims,ime, jms,jme, kms,kme, & ! memory dims
2962	ips,ipe, jps,jpe, kps,kpe, & ! patch dims
2963	grid%i_start(ij), grid%i_end(ij), &
2964	grid%j_start(ij), grid%j_end(ij), &
2965	k_start, k_end )
2966	#else
2967	CALL flow_dep_bdy ( tracer(ims,kms,jms,ic), &
2968	grid%ru_m, grid%rv_m, config_flags, &
2969	grid%spec_zone, &
2970	ids,ide, jds,jde, kds,kde, & ! domain dims
2971	ims,ime, jms,jme, kms,kme, & ! memory dims
2972	ips,ipe, jps,jpe, kps,kpe, & ! patch dims
2973	grid%i_start(ij), grid%i_end(ij), &
2974	grid%j_start(ij), grid%j_end(ij), &
2975	k_start, k_end )
2976	#endif
2977	ENDIF
2978	ENDDO tracer_tile_loop_2
2979	!$OMP END PARALLEL DO
2980
2981	ENDDO tracer_variable_loop
2982	ENDIF tracer_advance
2983	BENCH_END(tracer_adv_tim)
2984
2985	! next the other scalar species
2986	other_scalar_advance: IF (num_3d_s >= PARAM_FIRST_SCALAR) THEN
2987
2988	scalar_variable_loop: do is = PARAM_FIRST_SCALAR, num_3d_s
2989	!$OMP PARALLEL DO &
2990	!$OMP PRIVATE ( ij )
2991	scalar_tile_loop_1: DO ij = 1 , grid%num_tiles
2992
2993	CALL wrf_debug ( 200 , ' call rk_scalar_tend' )
2994	tenddec = .false.
2995	CALL rk_scalar_tend ( is, is, config_flags, tenddec, &
2996	rk_step, dt_rk, &
2997	grid%ru_m, grid%rv_m, grid%ww_m, &
2998	grid%muts, grid%mub, grid%mu_1, &
2999	grid%alt, &
3000	scalar_old(ims,kms,jms,is), &
3001	scalar(ims,kms,jms,is), &
3002	scalar_tend(ims,kms,jms,is), &
3003	advect_tend,h_tendency,z_tendency,grid%rqvften, &
3004	grid%qv_base, .false., grid%fnm, grid%fnp, &
3005	grid%msfux,grid%msfuy, grid%msfvx, grid%msfvx_inv, &
3006	grid%msfvy, grid%msftx,grid%msfty, &
3007	grid%rdx, grid%rdy, grid%rdn, grid%rdnw, &
3008	grid%khdif, grid%kvdif, grid%xkhh, &
3009	grid%diff_6th_opt, grid%diff_6th_factor, &
3010	config_flags%scalar_adv_opt, &
3011	ids, ide, jds, jde, kds, kde, &
3012	ims, ime, jms, jme, kms, kme, &
3013	grid%i_start(ij), grid%i_end(ij), &
3014	grid%j_start(ij), grid%j_end(ij), &
3015	k_start , k_end )
3016
3017	IF( config_flags%nested .and. (rk_step == 1) ) THEN
3018
3019	CALL relax_bdy_scalar ( scalar_tend(ims,kms,jms,is), &
3020	scalar(ims,kms,jms,is), grid%mut, &
3021	scalar_bxs(jms,kms,1,is),scalar_bxe(jms,kms,1,is), &
3022	scalar_bys(ims,kms,1,is),scalar_bye(ims,kms,1,is), &
3023	scalar_btxs(jms,kms,1,is),scalar_btxe(jms,kms,1,is), &
3024	scalar_btys(ims,kms,1,is),scalar_btye(ims,kms,1,is), &
3025	config_flags%spec_bdy_width, grid%spec_zone, grid%relax_zone, &
3026	grid%dtbc, grid%fcx, grid%gcx, &
3027	config_flags, &
3028	ids,ide, jds,jde, kds,kde, &
3029	ims,ime, jms,jme, kms,kme, &
3030	ips,ipe, jps,jpe, kps,kpe, &
3031	grid%i_start(ij), grid%i_end(ij), &
3032	grid%j_start(ij), grid%j_end(ij), &
3033	k_start, k_end )
3034
3035	CALL spec_bdy_scalar ( scalar_tend(ims,kms,jms,is), &
3036	scalar_bxs(jms,kms,1,is),scalar_bxe(jms,kms,1,is), &
3037	scalar_bys(ims,kms,1,is),scalar_bye(ims,kms,1,is), &
3038	scalar_btxs(jms,kms,1,is),scalar_btxe(jms,kms,1,is), &
3039	scalar_btys(ims,kms,1,is),scalar_btye(ims,kms,1,is), &
3040	config_flags%spec_bdy_width, grid%spec_zone, &
3041	config_flags, &
3042	ids,ide, jds,jde, kds,kde, &
3043	ims,ime, jms,jme, kms,kme, &
3044	ips,ipe, jps,jpe, kps,kpe, &
3045	grid%i_start(ij), grid%i_end(ij), &
3046	grid%j_start(ij), grid%j_end(ij), &
3047	k_start, k_end )
3048
3049	ENDIF ! b.c test for chem nested boundary condition
3050
3051	ENDDO scalar_tile_loop_1
3052	!$OMP END PARALLEL DO
3053
3054	!$OMP PARALLEL DO &
3055	!$OMP PRIVATE ( ij )
3056	scalar_tile_loop_2: DO ij = 1 , grid%num_tiles
3057
3058	CALL wrf_debug ( 200 , ' call rk_update_scalar' )
3059	tenddec = .false.
3060	CALL rk_update_scalar( scs=is, sce=is, &
3061	scalar_1=scalar_old(ims,kms,jms,is), &
3062	scalar_2=scalar(ims,kms,jms,is), &
3063	sc_tend=scalar_tend(ims,kms,jms,is), &
3064	! advh_t=advh_t(ims,kms,jms,1), &
3065	! advz_t=advz_t(ims,kms,jms,1), &
3066	advect_tend=advect_tend, &
3067	h_tendency=h_tendency, z_tendency=z_tendency, &
3068	msftx=grid%msftx,msfty=grid%msfty, &
3069	mu_old=grid%mu_1, mu_new=grid%mu_2, mu_base=grid%mub, &
3070	rk_step=rk_step, dt=dt_rk, spec_zone=grid%spec_zone, &
3071	config_flags=config_flags, tenddec=tenddec, &
3072	ids=ids, ide=ide, jds=jds, jde=jde, kds=kds, kde=kde, &
3073	ims=ims, ime=ime, jms=jms, jme=jme, kms=kms, kme=kme, &
3074	its=grid%i_start(ij), ite=grid%i_end(ij), &
3075	jts=grid%j_start(ij), jte=grid%j_end(ij), &
3076	kts=k_start , kte=k_end )
3077
3078	IF( config_flags%specified ) THEN
3079
3080	CALL flow_dep_bdy ( scalar(ims,kms,jms,is), &
3081	grid%ru_m, grid%rv_m, config_flags, &
3082	grid%spec_zone, &
3083	ids,ide, jds,jde, kds,kde, & ! domain dims
3084	ims,ime, jms,jme, kms,kme, & ! memory dims
3085	ips,ipe, jps,jpe, kps,kpe, & ! patch dims
3086	grid%i_start(ij), grid%i_end(ij), &
3087	grid%j_start(ij), grid%j_end(ij), &
3088	k_start, k_end )
3089	ENDIF
3090
3091	ENDDO scalar_tile_loop_2
3092	!$OMP END PARALLEL DO
3093
3094	ENDDO scalar_variable_loop
3095
3096	ENDIF other_scalar_advance
3097
3098	! update the pressure and density at the new time level
3099
3100	!$OMP PARALLEL DO &
3101	!$OMP PRIVATE ( ij )
3102	DO ij = 1 , grid%num_tiles
3103
3104	BENCH_START(calc_p_rho_tim)
3105
3106	CALL calc_p_rho_phi( moist, num_3d_m, &
3107	grid%al, grid%alb, grid%mu_2, grid%muts, &
3108	grid%ph_2, grid%p, grid%pb, grid%t_2, &
3109	p0, t0, grid%znu, grid%dnw, grid%rdnw, &
3110	grid%rdn, config_flags%non_hydrostatic, &
3111	ids, ide, jds, jde, kds, kde, &
3112	ims, ime, jms, jme, kms, kme, &
3113	grid%i_start(ij), grid%i_end(ij), &
3114	grid%j_start(ij), grid%j_end(ij), &
3115	k_start , k_end )
3116
3117	BENCH_END(calc_p_rho_tim)
3118
3119	ENDDO
3120	!$OMP END PARALLEL DO
3121
3122	! Reset the boundary conditions if there is another corrector step.
3123	! (rk_step < rk_order), else we'll handle it at the end of everything
3124	! (after the split physics, before exiting the timestep).
3125
3126	rk_step_1_check: IF ( rk_step < rk_order ) THEN
3127
3128	!-----------------------------------------------------------
3129	! rk3 substep polar filter for scalars (moist,chem,scalar)
3130	!-----------------------------------------------------------
3131
3132	IF (config_flags%polar) THEN
3133	IF ( num_3d_m >= PARAM_FIRST_SCALAR ) THEN
3134	CALL wrf_debug ( 200 , ' call filter moist ' )
3135	DO im = PARAM_FIRST_SCALAR, num_3d_m
3136	CALL couple_scalars_for_filter ( FIELD=moist(ims,kms,jms,im) &
3137	,MU=grid%mu_2 , MUB=grid%mub &
3138	,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
3139	,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
3140	,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe )
3141	CALL pxft ( grid=grid &
3142	,lineno=__LINE__ &
3143	,flag_uv = 0 &
3144	,flag_rurv = 0 &
3145	,flag_wph = 0 &
3146	,flag_ww = 0 &
3147	,flag_t = 0 &
3148	,flag_mu = 0 &
3149	,flag_mut = 0 &
3150	,flag_moist = im &
3151	,flag_chem = 0 &
3152	,flag_scalar = 0 &
3153	,flag_tracer = 0 &
3154	,positive_definite=.FALSE. &
3155	,moist=moist,chem=chem,tracer=tracer,scalar=scalar &
3156	,fft_filter_lat = config_flags%fft_filter_lat &
3157	,dclat = dclat &
3158	,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
3159	,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
3160	,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe &
3161	,imsx=imsx,imex=imex,jmsx=jmsx,jmex=jmex,kmsx=kmsx,kmex=kmex &
3162	,ipsx=ipsx,ipex=ipex,jpsx=jmsx,jpex=jpex,kpsx=kpsx,kpex=kpex )
3163	CALL uncouple_scalars_for_filter ( FIELD=moist(ims,kms,jms,im) &
3164	,MU=grid%mu_2 , MUB=grid%mub &
3165	,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
3166	,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
3167	,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe )
3168	END DO
3169	END IF
3170
3171	IF ( num_3d_c >= PARAM_FIRST_SCALAR ) THEN
3172	CALL wrf_debug ( 200 , ' call filter chem ' )
3173	DO im = PARAM_FIRST_SCALAR, num_3d_c
3174	CALL couple_scalars_for_filter ( FIELD=chem(ims,kms,jms,im) &
3175	,MU=grid%mu_2 , MUB=grid%mub &
3176	,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
3177	,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
3178	,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe )
3179	CALL pxft ( grid=grid &
3180	,lineno=__LINE__ &
3181	,flag_uv = 0 &
3182	,flag_rurv = 0 &
3183	,flag_wph = 0 &
3184	,flag_ww = 0 &
3185	,flag_t = 0 &
3186	,flag_mu = 0 &
3187	,flag_mut = 0 &
3188	,flag_moist = 0 &
3189	,flag_chem = im &
3190	,flag_tracer = 0 &
3191	,flag_scalar = 0 &
3192	,positive_definite=.FALSE. &
3193	,moist=moist,chem=chem,tracer=tracer,scalar=scalar &
3194	,fft_filter_lat = config_flags%fft_filter_lat &
3195	,dclat = dclat &
3196	,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
3197	,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
3198	,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe &
3199	,imsx=imsx,imex=imex,jmsx=jmsx,jmex=jmex,kmsx=kmsx,kmex=kmex &
3200	,ipsx=ipsx,ipex=ipex,jpsx=jmsx,jpex=jpex,kpsx=kpsx,kpex=kpex )
3201	CALL uncouple_scalars_for_filter ( FIELD=chem(ims,kms,jms,im) &
3202	,MU=grid%mu_2 , MUB=grid%mub &
3203	,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
3204	,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
3205	,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe )
3206	END DO
3207	END IF
3208	IF ( num_tracer >= PARAM_FIRST_SCALAR ) THEN
3209	CALL wrf_debug ( 200 , ' call filter tracer ' )
3210	DO im = PARAM_FIRST_SCALAR, num_tracer
3211	CALL couple_scalars_for_filter ( FIELD=tracer(ims,kms,jms,im) &
3212	,MU=grid%mu_2 , MUB=grid%mub &
3213	,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
3214	,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
3215	,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe )
3216	CALL pxft ( grid=grid &
3217	,lineno=__LINE__ &
3218	,flag_uv = 0 &
3219	,flag_rurv = 0 &
3220	,flag_wph = 0 &
3221	,flag_ww = 0 &
3222	,flag_t = 0 &
3223	,flag_mu = 0 &
3224	,flag_mut = 0 &
3225	,flag_moist = 0 &
3226	,flag_chem = 0 &
3227	,flag_tracer = im &
3228	,flag_scalar = 0 &
3229	,positive_definite=.FALSE. &
3230	,moist=moist,chem=chem,tracer=tracer,scalar=scalar &
3231	,fft_filter_lat = config_flags%fft_filter_lat &
3232	,dclat = dclat &
3233	,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
3234	,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
3235	,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe &
3236	,imsx=imsx,imex=imex,jmsx=jmsx,jmex=jmex,kmsx=kmsx,kmex=kmex &
3237	,ipsx=ipsx,ipex=ipex,jpsx=jmsx,jpex=jpex,kpsx=kpsx,kpex=kpex )
3238	CALL uncouple_scalars_for_filter ( FIELD=tracer(ims,kms,jms,im) &
3239	,MU=grid%mu_2 , MUB=grid%mub &
3240	,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
3241	,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
3242	,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe )
3243	END DO
3244	END IF
3245
3246	IF ( num_3d_s >= PARAM_FIRST_SCALAR ) THEN
3247	CALL wrf_debug ( 200 , ' call filter scalar ' )
3248	DO im = PARAM_FIRST_SCALAR, num_3d_s
3249	CALL couple_scalars_for_filter ( FIELD=scalar(ims,kms,jms,im) &
3250	,MU=grid%mu_2 , MUB=grid%mub &
3251	,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
3252	,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
3253	,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe )
3254	CALL pxft ( grid=grid &
3255	,lineno=__LINE__ &
3256	,flag_uv = 0 &
3257	,flag_rurv = 0 &
3258	,flag_wph = 0 &
3259	,flag_ww = 0 &
3260	,flag_t = 0 &
3261	,flag_mu = 0 &
3262	,flag_mut = 0 &
3263	,flag_moist = 0 &
3264	,flag_chem = 0 &
3265	,flag_tracer = 0 &
3266	,flag_scalar = im &
3267	,positive_definite=.FALSE. &
3268	,moist=moist,chem=chem,tracer=tracer,scalar=scalar &
3269	,fft_filter_lat = config_flags%fft_filter_lat &
3270	,dclat = dclat &
3271	,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
3272	,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
3273	,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe &
3274	,imsx=imsx,imex=imex,jmsx=jmsx,jmex=jmex,kmsx=kmsx,kmex=kmex &
3275	,ipsx=ipsx,ipex=ipex,jpsx=jmsx,jpex=jpex,kpsx=kpsx,kpex=kpex )
3276	CALL uncouple_scalars_for_filter ( FIELD=scalar(ims,kms,jms,im) &
3277	,MU=grid%mu_2 , MUB=grid%mub &
3278	,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
3279	,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
3280	,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe )
3281	END DO
3282	END IF
3283	END IF ! polar filter test
3284
3285	!-----------------------------------------------------------
3286	! END rk3 substep polar filter for scalars (moist,chem,scalar)
3287	!-----------------------------------------------------------
3288
3289	!-----------------------------------------------------------
3290	! Stencils for patch communications (WCS, 29 June 2001)
3291	!
3292	! here's where we need a wide comm stencil - these are the
3293	! uncoupled variables so are used for high order calc in
3294	! advection and mixong routines.
3295	!
3296	!
3297	! * * * * * * *
3298	! * * * * * * * * * * * *
3299	! * * * * * * * * * * * * *
3300	! * + * * * + * * * * * + * * *
3301	! * * * * * * * * * * * * *
3302	! * * * * * * * * * * * *
3303	! * * * * * * *
3304	!
3305	! al x
3306	!
3307	! 2D variable
3308	! mu_2 x
3309	!
3310	! (adv order <=4)
3311	! u_2 x
3312	! v_2 x
3313	! w_2 x
3314	! t_2 x
3315	! ph_2 x
3316	!
3317	! (adv order <=6)
3318	! u_2 x
3319	! v_2 x
3320	! w_2 x
3321	! t_2 x
3322	! ph_2 x
3323	!
3324	! 4D variable
3325	! moist x
3326	! chem x
3327	! scalar x
3328
3329	#ifdef DM_PARALLEL
3330	IF ( config_flags%h_mom_adv_order <= 4 ) THEN
3331	# include "HALO_EM_D2_3.inc"
3332	ELSE IF ( config_flags%h_mom_adv_order <= 6 ) THEN
3333	# include "HALO_EM_D2_5.inc"
3334	ELSE
3335	WRITE(wrf_err_message,*)'solve_em: invalid h_mom_adv_order = ',config_flags%h_mom_adv_order
3336	CALL wrf_error_fatal(TRIM(wrf_err_message))
3337	ENDIF
3338	# include "PERIOD_BDY_EM_D.inc"
3339	# include "PERIOD_BDY_EM_MOIST2.inc"
3340	# include "PERIOD_BDY_EM_CHEM2.inc"
3341	# include "PERIOD_BDY_EM_TRACER2.inc"
3342	# include "PERIOD_BDY_EM_SCALAR2.inc"
3343	#endif
3344
3345	BENCH_START(bc_end_tim)
3346	!$OMP PARALLEL DO &
3347	!$OMP PRIVATE ( ij )
3348	tile_bc_loop_1: DO ij = 1 , grid%num_tiles
3349	CALL wrf_debug ( 200 , ' call rk_phys_bc_dry_2' )
3350
3351	CALL rk_phys_bc_dry_2( config_flags, &
3352	grid%u_2, grid%v_2, grid%w_2, &
3353	grid%t_2, grid%ph_2, grid%mu_2, &
3354	ids, ide, jds, jde, kds, kde, &
3355	ims, ime, jms, jme, kms, kme, &
3356	ips, ipe, jps, jpe, kps, kpe, &
3357	grid%i_start(ij), grid%i_end(ij), &
3358	grid%j_start(ij), grid%j_end(ij), &
3359	k_start , k_end )
3360
3361	BENCH_START(diag_w_tim)
3362	IF (.not. config_flags%non_hydrostatic) THEN
3363	CALL diagnose_w( ph_tend, grid%ph_2, grid%ph_1, grid%w_2, grid%muts, dt_rk, &
3364	grid%u_2, grid%v_2, grid%ht, &
3365	grid%cf1, grid%cf2, grid%cf3, grid%rdx, grid%rdy, grid%msftx, grid%msfty, &
3366	ids, ide, jds, jde, kds, kde, &
3367	ims, ime, jms, jme, kms, kme, &
3368	grid%i_start(ij), grid%i_end(ij), &
3369	grid%j_start(ij), grid%j_end(ij), &
3370	k_start , k_end )
3371	ENDIF
3372	BENCH_END(diag_w_tim)
3373
3374	IF (num_3d_m >= PARAM_FIRST_SCALAR) THEN
3375
3376	moisture_loop_bdy_1 : DO im = PARAM_FIRST_SCALAR , num_3d_m
3377
3378	CALL set_physical_bc3d( moist(ims,kms,jms,im), 'p', config_flags, &
3379	ids, ide, jds, jde, kds, kde, &
3380	ims, ime, jms, jme, kms, kme, &
3381	ips, ipe, jps, jpe, kps, kpe, &
3382	grid%i_start(ij), grid%i_end(ij), &
3383	grid%j_start(ij), grid%j_end(ij), &
3384	k_start , k_end )
3385	END DO moisture_loop_bdy_1
3386
3387	ENDIF
3388
3389	IF (num_3d_c >= PARAM_FIRST_SCALAR) THEN
3390
3391	chem_species_bdy_loop_1 : DO ic = PARAM_FIRST_SCALAR , num_3d_c
3392
3393	CALL set_physical_bc3d( chem(ims,kms,jms,ic), 'p', config_flags, &
3394	ids, ide, jds, jde, kds, kde, &
3395	ims, ime, jms, jme, kms, kme, &
3396	ips, ipe, jps, jpe, kps, kpe, &
3397	grid%i_start(ij), grid%i_end(ij), &
3398	grid%j_start(ij), grid%j_end(ij), &
3399	k_start , k_end-1 )
3400
3401	END DO chem_species_bdy_loop_1
3402
3403	END IF
3404
3405	IF (num_tracer >= PARAM_FIRST_SCALAR) THEN
3406
3407	tracer_species_bdy_loop_1 : DO ic = PARAM_FIRST_SCALAR , num_tracer
3408
3409	CALL set_physical_bc3d( tracer(ims,kms,jms,ic), 'p', config_flags, &
3410	ids, ide, jds, jde, kds, kde, &
3411	ims, ime, jms, jme, kms, kme, &
3412	ips, ipe, jps, jpe, kps, kpe, &
3413	grid%i_start(ij), grid%i_end(ij), &
3414	grid%j_start(ij), grid%j_end(ij), &
3415	k_start , k_end-1 )
3416
3417	END DO tracer_species_bdy_loop_1
3418
3419	END IF
3420
3421	IF (num_3d_s >= PARAM_FIRST_SCALAR) THEN
3422
3423	scalar_species_bdy_loop_1 : DO is = PARAM_FIRST_SCALAR , num_3d_s
3424
3425	CALL set_physical_bc3d( scalar(ims,kms,jms,is), 'p', config_flags, &
3426	ids, ide, jds, jde, kds, kde, &
3427	ims, ime, jms, jme, kms, kme, &
3428	ips, ipe, jps, jpe, kps, kpe, &
3429	grid%i_start(ij), grid%i_end(ij), &
3430	grid%j_start(ij), grid%j_end(ij), &
3431	k_start , k_end-1 )
3432
3433	END DO scalar_species_bdy_loop_1
3434
3435	END IF
3436
3437	IF (config_flags%km_opt .eq. 2) THEN
3438
3439	CALL set_physical_bc3d( grid%tke_2 , 'p', config_flags, &
3440	ids, ide, jds, jde, kds, kde, &
3441	ims, ime, jms, jme, kms, kme, &
3442	ips, ipe, jps, jpe, kps, kpe, &
3443	grid%i_start(ij), grid%i_end(ij), &
3444	grid%j_start(ij), grid%j_end(ij), &
3445	k_start , k_end )
3446	END IF
3447
3448	END DO tile_bc_loop_1
3449	!$OMP END PARALLEL DO
3450	BENCH_END(bc_end_tim)
3451
3452
3453	#ifdef DM_PARALLEL
3454
3455	! * * * * *
3456	! * * * * * * * * *
3457	! * + * * + * * * + * *
3458	! * * * * * * * * *
3459	! * * * * *
3460
3461	! moist, chem, scalar, tke x
3462
3463
3464	IF ( config_flags%h_mom_adv_order <= 4 ) THEN
3465	IF ( (config_flags%tke_adv_opt /= ORIGINAL) .and. (rk_step == rk_order-1) ) THEN
3466	# include "HALO_EM_TKE_5.inc"
3467	ELSE
3468	# include "HALO_EM_TKE_3.inc"
3469	ENDIF
3470	ELSE IF ( config_flags%h_mom_adv_order <= 6 ) THEN
3471	IF ( (config_flags%tke_adv_opt /= ORIGINAL) .and. (rk_step == rk_order-1) ) THEN
3472	# include "HALO_EM_TKE_7.inc"
3473	ELSE
3474	# include "HALO_EM_TKE_5.inc"
3475	ENDIF
3476	ELSE
3477	WRITE(wrf_err_message,*)'solve_em: invalid h_sca_adv_order = ',config_flags%h_sca_adv_order
3478	CALL wrf_error_fatal(TRIM(wrf_err_message))
3479	ENDIF
3480
3481	IF ( num_moist .GE. PARAM_FIRST_SCALAR ) THEN
3482	IF ( config_flags%h_sca_adv_order <= 4 ) THEN
3483	IF ( (config_flags%moist_adv_opt /= ORIGINAL) .and. (rk_step == rk_order-1) ) THEN
3484	# include "HALO_EM_MOIST_E_5.inc"
3485	ELSE
3486	# include "HALO_EM_MOIST_E_3.inc"
3487	END IF
3488	ELSE IF ( config_flags%h_sca_adv_order <= 6 ) THEN
3489	IF ( (config_flags%moist_adv_opt /= ORIGINAL) .and. (rk_step == rk_order-1) ) THEN
3490	# include "HALO_EM_MOIST_E_7.inc"
3491	ELSE
3492	# include "HALO_EM_MOIST_E_5.inc"
3493	END IF
3494	ELSE
3495	WRITE(wrf_err_message,*)'solve_em: invalid h_sca_adv_order = ',config_flags%h_sca_adv_order
3496	CALL wrf_error_fatal(TRIM(wrf_err_message))
3497	ENDIF
3498	ENDIF
3499	IF ( num_chem >= PARAM_FIRST_SCALAR ) THEN
3500	IF ( config_flags%h_sca_adv_order <= 4 ) THEN
3501	IF ( (config_flags%chem_adv_opt /= ORIGINAL) .and. (rk_step == rk_order-1) ) THEN
3502	# include "HALO_EM_CHEM_E_5.inc"
3503	ELSE
3504	# include "HALO_EM_CHEM_E_3.inc"
3505	ENDIF
3506	ELSE IF ( config_flags%h_sca_adv_order <= 6 ) THEN
3507	IF ( (config_flags%chem_adv_opt /= ORIGINAL) .and. (rk_step == rk_order-1) ) THEN
3508	# include "HALO_EM_CHEM_E_7.inc"
3509	ELSE
3510	# include "HALO_EM_CHEM_E_5.inc"
3511	ENDIF
3512	ELSE
3513	WRITE(wrf_err_message,*)'solve_em: invalid h_sca_adv_order = ',config_flags%h_sca_adv_order
3514	CALL wrf_error_fatal(TRIM(wrf_err_message))
3515	ENDIF
3516	ENDIF
3517	IF ( num_tracer >= PARAM_FIRST_SCALAR ) THEN
3518	IF ( config_flags%h_sca_adv_order <= 4 ) THEN
3519	IF ( (config_flags%tracer_adv_opt /= ORIGINAL) .and. (rk_step == rk_order-1) ) THEN
3520	# include "HALO_EM_TRACER_E_5.inc"
3521	ELSE
3522	# include "HALO_EM_TRACER_E_3.inc"
3523	ENDIF
3524	ELSE IF ( config_flags%h_sca_adv_order <= 6 ) THEN
3525	IF ( (config_flags%tracer_adv_opt /= ORIGINAL) .and. (rk_step == rk_order-1) ) THEN
3526	# include "HALO_EM_TRACER_E_7.inc"
3527	ELSE
3528	# include "HALO_EM_TRACER_E_5.inc"
3529	ENDIF
3530	ELSE
3531	WRITE(wrf_err_message,*)'solve_em: invalid h_sca_adv_order = ',config_flags%h_sca_adv_order
3532	CALL wrf_error_fatal(TRIM(wrf_err_message))
3533	ENDIF
3534	ENDIF
3535	IF ( num_scalar >= PARAM_FIRST_SCALAR ) THEN
3536	IF ( config_flags%h_sca_adv_order <= 4 ) THEN
3537	IF ( (config_flags%scalar_adv_opt /= ORIGINAL) .and. (rk_step == rk_order-1) ) THEN
3538	# include "HALO_EM_SCALAR_E_5.inc"
3539	ELSE
3540	# include "HALO_EM_SCALAR_E_3.inc"
3541	ENDIF
3542	ELSE IF ( config_flags%h_sca_adv_order <= 6 ) THEN
3543	IF ( (config_flags%scalar_adv_opt /= ORIGINAL) .and. (rk_step == rk_order-1) ) THEN
3544	# include "HALO_EM_SCALAR_E_7.inc"
3545	ELSE
3546	# include "HALO_EM_SCALAR_E_5.inc"
3547	ENDIF
3548	ELSE
3549	WRITE(wrf_err_message,*)'solve_em: invalid h_sca_adv_order = ',config_flags%h_sca_adv_order
3550	CALL wrf_error_fatal(TRIM(wrf_err_message))
3551	ENDIF
3552	ENDIF
3553	#endif
3554
3555	ENDIF rk_step_1_check
3556
3557
3558	!**********************************************************
3559	!
3560	! end of RK predictor-corrector loop
3561	!
3562	!**********************************************************
3563
3564	END DO Runge_Kutta_loop
3565	#ifdef LMDZ1
3566	WRITE(message, *)' dyn_em: after runge_kutta_loop'
3567	CALL wrf_debug(200, message)
3568	WRITE(message, *)' t_tend: ',t_tendf(im2,km2,jm2), &
3569	' u_tend: ', ru_tendf(im2,1,jm2)
3570	CALL wrf_debug(200, message)
3571	WRITE(message,*)' psfc_tend: ',grid%dpsdt(im2,jm2), &
3572	'p sfc: ',p8w(im2,kms,jm2)
3573	CALL wrf_debug(200, message)
3574	WRITE(message,*)' p 1: ',grid%p(im2,kms,jm2), ' ph 1: ',grid%ph_2(im2,kms,jm2)
3575	CALL wrf_debug(200, message)
3576	#endif
3577
3578	IF (config_flags%do_avgflx_em .EQ. 1) THEN
3579	! Reinitialize time-averaged fluxes if history output was written after the previous time step:
3580	CALL WRFU_ALARMGET(grid%alarms( HISTORY_ALARM ),prevringtime=temp_time)
3581	CALL domain_clock_get ( grid, current_time=CurrTime, &
3582	current_timestr=message2 )
3583	! use overloaded -, .LT. operator to check whether to initialize avgflx:
3584	! reinitialize after each history output (detect this here by comparing current time
3585	! against last history time and time step - this code follows what's done in adapt_timestep_em):
3586	WRITE ( message , FMT = '("solve_em: old_dt =",g15.6,", dt=",g15.6," on domain ",I3)' ) &
3587	& old_dt,grid%dt,grid%id
3588	CALL wrf_debug(200,message)
3589	old_dt=min(old_dt,grid%dt)
3590	num = INT(old_dt * precision)
3591	den = precision
3592	CALL WRFU_TimeIntervalSet(dtInterval, Sn=num, Sd=den)
3593	IF (CurrTime .lt. temp_time + dtInterval) THEN
3594	WRITE ( message , FMT = '("solve_em: initializing avgflx at time ",A," on domain ",I3)' ) &
3595	& TRIM(message2), grid%id
3596	CALL wrf_message(trim(message))
3597	grid%avgflx_count = 0
3598	!tile-loop for zero_avgflx
3599	!$OMP PARALLEL DO &
3600	!$OMP PRIVATE ( ij )
3601
3602	DO ij = 1 , grid%num_tiles
3603	CALL wrf_debug(200,'In solve_em, before zero_avgflx call')
3604	CALL zero_avgflx(grid%avgflx_rum,grid%avgflx_rvm,grid%avgflx_wwm, &
3605	& ids, ide, jds, jde, kds, kde, &
3606	& ims, ime, jms, jme, kms, kme, &
3607	& grid%i_start(ij), grid%i_end(ij), grid%j_start(ij), grid%j_end(ij), &
3608	& k_start , k_end, f_flux, &
3609	& grid%avgflx_cfu1,grid%avgflx_cfd1,grid%avgflx_dfu1, &
3610	& grid%avgflx_efu1,grid%avgflx_dfd1,grid%avgflx_efd1 )
3611	CALL wrf_debug(200,'In solve_em, after zero_avgflx call')
3612	ENDDO
3613	ENDIF
3614
3615	! Update avgflx quantities
3616	!tile-loop for upd_avgflx
3617	!$OMP PARALLEL DO &
3618	!$OMP PRIVATE ( ij )
3619
3620	DO ij = 1 , grid%num_tiles
3621	CALL wrf_debug(200,'In solve_em, before upd_avgflx call')
3622	CALL upd_avgflx(grid%avgflx_count,grid%avgflx_rum,grid%avgflx_rvm,grid%avgflx_wwm, &
3623	& grid%ru_m, grid%rv_m, grid%ww_m, &
3624	& ids, ide, jds, jde, kds, kde, &
3625	& ims, ime, jms, jme, kms, kme, &
3626	& grid%i_start(ij), grid%i_end(ij), grid%j_start(ij), grid%j_end(ij), &
3627	& k_start , k_end, f_flux, &
3628	& grid%cfu1,grid%cfd1,grid%dfu1,grid%efu1,grid%dfd1,grid%efd1, &
3629	& grid%avgflx_cfu1,grid%avgflx_cfd1,grid%avgflx_dfu1, &
3630	& grid%avgflx_efu1,grid%avgflx_dfd1,grid%avgflx_efd1 )
3631	CALL wrf_debug(200,'In solve_em, after upd_avgflx call')
3632
3633	ENDDO
3634	grid%avgflx_count = grid%avgflx_count + 1
3635	ENDIF
3636	!
3637	!$OMP PARALLEL DO &
3638	!$OMP PRIVATE ( ij )
3639	DO ij = 1 , grid%num_tiles
3640
3641	BENCH_START(advance_ppt_tim)
3642	CALL wrf_debug ( 200 , ' call advance_ppt' )
3643	CALL advance_ppt(grid%rthcuten,grid%rqvcuten,grid%rqccuten,grid%rqrcuten, &
3644	grid%rqicuten,grid%rqscuten, &
3645	grid%rainc,grid%raincv,grid%rainsh,grid%pratec,grid%pratesh, &
3646	grid%nca,grid%htop,grid%hbot,grid%cutop,grid%cubot, &
3647	grid%cuppt, grid%dt, config_flags, &
3648	ids,ide, jds,jde, kds,kde, &
3649	ims,ime, jms,jme, kms,kme, &
3650	grid%i_start(ij), grid%i_end(ij), &
3651	grid%j_start(ij), grid%j_end(ij), &
3652	k_start , k_end )
3653	BENCH_END(advance_ppt_tim)
3654
3655	ENDDO
3656	!$OMP END PARALLEL DO
3657
3658	!<DESCRIPTION>
3659	!<pre>
3660	! (5) time-split physics.
3661	!
3662	! Microphysics are the only time split physics in the WRF model
3663	! at this time. Split-physics begins with the calculation of
3664	! needed diagnostic quantities (pressure, temperature, etc.)
3665	! followed by a call to the microphysics driver,
3666	! and finishes with a clean-up, storing off of a diabatic tendency
3667	! from the moist physics, and a re-calulation of the diagnostic
3668	! quantities pressure and density.
3669	!</pre>
3670	!</DESCRIPTION>
3671
3672	IF( config_flags%specified .or. config_flags%nested ) THEN
3673	sz = grid%spec_zone
3674	ELSE
3675	sz = 0
3676	ENDIF
3677
3678	IF (config_flags%mp_physics /= 0) then
3679
3680	!$OMP PARALLEL DO &
3681	!$OMP PRIVATE ( ij, its, ite, jts, jte )
3682
3683	scalar_tile_loop_1a: DO ij = 1 , grid%num_tiles
3684
3685	IF ( config_flags%periodic_x ) THEN
3686	its = max(grid%i_start(ij),ids)
3687	ite = min(grid%i_end(ij),ide-1)
3688	ELSE
3689	its = max(grid%i_start(ij),ids+sz)
3690	ite = min(grid%i_end(ij),ide-1-sz)
3691	ENDIF
3692	jts = max(grid%j_start(ij),jds+sz)
3693	jte = min(grid%j_end(ij),jde-1-sz)
3694
3695	CALL wrf_debug ( 200 , ' call moist_physics_prep' )
3696	BENCH_START(moist_physics_prep_tim)
3697	CALL moist_physics_prep_em( grid%t_2, grid%t_1, t0, rho, &
3698	grid%al, grid%alb, grid%p, p8w, p0, grid%pb, &
3699	grid%ph_2, grid%phb, th_phy, pi_phy, p_phy, &
3700	grid%z, grid%z_at_w, dz8w, &
3701	dtm, grid%h_diabatic, &
3702	config_flags,grid%fnm, grid%fnp, &
3703	ids, ide, jds, jde, kds, kde, &
3704	ims, ime, jms, jme, kms, kme, &
3705	its, ite, jts, jte, &
3706	k_start , k_end )
3707	BENCH_END(moist_physics_prep_tim)
3708	END DO scalar_tile_loop_1a
3709	!$OMP END PARALLEL DO
3710
3711	CALL wrf_debug ( 200 , ' call microphysics_driver' )
3712
3713	grid%sr = 0.
3714	specified_bdy = config_flags%specified .OR. config_flags%nested
3715	channel_bdy = config_flags%specified .AND. config_flags%periodic_x
3716
3717	BENCH_START(micro_driver_tim)
3718
3719	CALL microphysics_driver( &
3720	& DT=dtm ,DX=grid%dx ,DY=grid%dy &
3721	& ,DZ8W=dz8w ,F_ICE_PHY=grid%f_ice_phy &
3722	& ,ITIMESTEP=grid%itimestep ,LOWLYR=grid%lowlyr &
3723	& ,P8W=p8w ,P=p_phy ,PI_PHY=pi_phy &
3724	& ,RHO=rho ,SPEC_ZONE=grid%spec_zone &
3725	& ,SR=grid%sr ,TH=th_phy &
3726	& ,refl_10cm=grid%refl_10cm & ! hm, 9/22/09 for refl
3727	& ,WARM_RAIN=grid%warm_rain &
3728	& ,T8W=t8w &
3729	& ,CLDFRA=grid%cldfra, EXCH_H=grid%exch_h &
3730	& ,NSOURCE=grid%qndropsource &
3731	#ifdef WRF_CHEM
3732	& ,QLSINK=grid%qlsink,CLDFRA_OLD=grid%cldfra_old &
3733	& ,PRECR=grid%precr, PRECI=grid%preci, PRECS=grid%precs, PRECG=grid%precg &
3734	& ,CHEM_OPT=config_flags%chem_opt, PROGN=config_flags%progn &
3735	#endif
3736	& ,XLAND=grid%xland &
3737	& ,SPECIFIED=specified_bdy, CHANNEL_SWITCH=channel_bdy &
3738	& ,F_RAIN_PHY=grid%f_rain_phy &
3739	& ,F_RIMEF_PHY=grid%f_rimef_phy &
3740	& ,MP_PHYSICS=config_flags%mp_physics &
3741	& ,ID=grid%id &
3742	& ,IDS=ids,IDE=ide, JDS=jds,JDE=jde, KDS=kds,KDE=kde &
3743	& ,IMS=ims,IME=ime, JMS=jms,JME=jme, KMS=kms,KME=kme &
3744	#ifdef RUN_ON_GPU
3745	& ,IPS=ips,IPE=ipe, JPS=jps,JPE=jpe, KPS=kps,KPE=kpe &
3746	#endif
3747	& ,I_START=grid%i_start,I_END=min(grid%i_end, ide-1) &
3748	& ,J_START=grid%j_start,J_END=min(grid%j_end, jde-1) &
3749	& ,KTS=k_start, KTE=min(k_end,kde-1) &
3750	& ,NUM_TILES=grid%num_tiles &
3751	& ,NAER=grid%naer &
3752	! Optional
3753	& , RAINNC=grid%rainnc, RAINNCV=grid%rainncv &
3754	& , SNOWNC=grid%snownc, SNOWNCV=grid%snowncv &
3755	& , GRAUPELNC=grid%graupelnc, GRAUPELNCV=grid%graupelncv & ! for milbrandt2mom
3756	& , HAILNC=grid%hailnc, HAILNCV=grid%hailncv &
3757	& , W=grid%w_2, Z=grid%z, HT=grid%ht &
3758	& , MP_RESTART_STATE=grid%mp_restart_state &
3759	& , TBPVS_STATE=grid%tbpvs_state & ! etampnew
3760	& , TBPVS0_STATE=grid%tbpvs0_state & ! etampnew
3761	& , QV_CURR=moist(ims,kms,jms,P_QV), F_QV=F_QV &
3762	& , QC_CURR=moist(ims,kms,jms,P_QC), F_QC=F_QC &
3763	& , QR_CURR=moist(ims,kms,jms,P_QR), F_QR=F_QR &
3764	& , QI_CURR=moist(ims,kms,jms,P_QI), F_QI=F_QI &
3765	& , QS_CURR=moist(ims,kms,jms,P_QS), F_QS=F_QS &
3766	& , QG_CURR=moist(ims,kms,jms,P_QG), F_QG=F_QG &
3767	& , QH_CURR=moist(ims,kms,jms,P_QH), F_QH=F_QH & ! for milbrandt2mom
3768	& , QNDROP_CURR=scalar(ims,kms,jms,P_QNDROP), F_QNDROP=F_QNDROP &
3769	& , QT_CURR=scalar(ims,kms,jms,P_QT), F_QT=F_QT &
3770	& , QNN_CURR=scalar(ims,kms,jms,P_QNN), F_QNN=F_QNN &
3771	& , QNI_CURR=scalar(ims,kms,jms,P_QNI), F_QNI=F_QNI &
3772	& , QNC_CURR=scalar(ims,kms,jms,P_QNC), F_QNC=F_QNC &
3773	& , QNR_CURR=scalar(ims,kms,jms,P_QNR), F_QNR=F_QNR &
3774	& , QNS_CURR=scalar(ims,kms,jms,P_QNS), F_QNS=F_QNS &
3775	& , QNG_CURR=scalar(ims,kms,jms,P_QNG), F_QNG=F_QNG &
3776	& , QNH_CURR=scalar(ims,kms,jms,P_QNH), F_QNH=F_QNH & ! for milbrandt2mom
3777	! & , QZR_CURR=scalar(ims,kms,jms,P_QZR), F_QZR=F_QZR & ! for milbrandt3mom
3778	! & , QZI_CURR=scalar(ims,kms,jms,P_QZI), F_QZI=F_QZI & ! "
3779	! & , QZS_CURR=scalar(ims,kms,jms,P_QZS), F_QZS=F_QZS & ! "
3780	! & , QZG_CURR=scalar(ims,kms,jms,P_QZG), F_QZG=F_QZG & ! "
3781	! & , QZH_CURR=scalar(ims,kms,jms,P_QZH), F_QZH=F_QZH & ! "
3782	& , qrcuten=grid%rqrcuten, qscuten=grid%rqscuten &
3783	& , qicuten=grid%rqicuten,mu=grid%mut &
3784	& , HAIL=config_flags%gsfcgce_hail & ! for gsfcgce
3785	& , ICE2=config_flags%gsfcgce_2ice & ! for gsfcgce
3786	! & , ccntype=config_flags%milbrandt_ccntype & ! for milbrandt (2mom)
3787	! YLIN
3788	! RI_CURR INPUT
3789	& , RI_CURR=grid%rimi &
3790	)
3791	BENCH_END(micro_driver_tim)
3792
3793	#ifdef LMDZ
3794	grid%h_diabatic = 0.
3795	#endif
3796
3797	#if 0
3798	BENCH_START(microswap_2)
3799	! for load balancing; communication to redistribute the points
3800	IF ( config_flags%mp_physics .EQ. ETAMPNEW ) THEN
3801	#include "SWAP_ETAMP_NEW.inc"
3802	ELSE IF ( config_flags%mp_physics .EQ. WSM3SCHEME ) THEN
3803	#include "SWAP_WSM3.inc"
3804	ENDIF
3805	BENCH_END(microswap_2)
3806	#endif
3807
3808	CALL wrf_debug ( 200 , ' call moist_physics_finish' )
3809	BENCH_START(moist_phys_end_tim)
3810
3811	!$OMP PARALLEL DO &
3812	!$OMP PRIVATE ( ij, its, ite, jts, jte, im, ii, jj, kk )
3813
3814	DO ij = 1 , grid%num_tiles
3815
3816	its = max(grid%i_start(ij),ids)
3817	ite = min(grid%i_end(ij),ide-1)
3818	jts = max(grid%j_start(ij),jds)
3819	jte = min(grid%j_end(ij),jde-1)
3820
3821	CALL microphysics_zero_outb ( &
3822	moist , num_moist , config_flags , &
3823	ids, ide, jds, jde, kds, kde, &
3824	ims, ime, jms, jme, kms, kme, &
3825	its, ite, jts, jte, &
3826	k_start , k_end )
3827
3828	CALL microphysics_zero_outb ( &
3829	scalar , num_scalar , config_flags , &
3830	ids, ide, jds, jde, kds, kde, &
3831	ims, ime, jms, jme, kms, kme, &
3832	its, ite, jts, jte, &
3833	k_start , k_end )
3834
3835	CALL microphysics_zero_outb ( &
3836	chem , num_chem , config_flags , &
3837	ids, ide, jds, jde, kds, kde, &
3838	ims, ime, jms, jme, kms, kme, &
3839	its, ite, jts, jte, &
3840	k_start , k_end )
3841	CALL microphysics_zero_outb ( &
3842	tracer , num_tracer , config_flags , &
3843	ids, ide, jds, jde, kds, kde, &
3844	ims, ime, jms, jme, kms, kme, &
3845	its, ite, jts, jte, &
3846	k_start , k_end )
3847
3848	IF ( config_flags%periodic_x ) THEN
3849	its = max(grid%i_start(ij),ids)
3850	ite = min(grid%i_end(ij),ide-1)
3851	ELSE
3852	its = max(grid%i_start(ij),ids+sz)
3853	ite = min(grid%i_end(ij),ide-1-sz)
3854	ENDIF
3855	jts = max(grid%j_start(ij),jds+sz)
3856	jte = min(grid%j_end(ij),jde-1-sz)
3857
3858	CALL microphysics_zero_outa ( &
3859	moist , num_moist , config_flags , &
3860	ids, ide, jds, jde, kds, kde, &
3861	ims, ime, jms, jme, kms, kme, &
3862	its, ite, jts, jte, &
3863	k_start , k_end )
3864
3865	CALL microphysics_zero_outa ( &
3866	scalar , num_scalar , config_flags , &
3867	ids, ide, jds, jde, kds, kde, &
3868	ims, ime, jms, jme, kms, kme, &
3869	its, ite, jts, jte, &
3870	k_start , k_end )
3871
3872	CALL microphysics_zero_outa ( &
3873	chem , num_chem , config_flags , &
3874	ids, ide, jds, jde, kds, kde, &
3875	ims, ime, jms, jme, kms, kme, &
3876	its, ite, jts, jte, &
3877	k_start , k_end )
3878
3879	CALL microphysics_zero_outa ( &
3880	tracer , num_tracer , config_flags , &
3881	ids, ide, jds, jde, kds, kde, &
3882	ims, ime, jms, jme, kms, kme, &
3883	its, ite, jts, jte, &
3884	k_start , k_end )
3885
3886	CALL moist_physics_finish_em( grid%t_2, grid%t_1, t0, grid%muts, th_phy, &
3887	grid%h_diabatic, dtm, config_flags, &
3888	#if ( WRF_DFI_RADAR == 1 )
3889	grid%dfi_tten_rad,grid%dfi_stage, &
3890	#endif
3891	ids, ide, jds, jde, kds, kde, &
3892	ims, ime, jms, jme, kms, kme, &
3893	its, ite, jts, jte, &
3894	k_start , k_end )
3895
3896	END DO
3897	!$OMP END PARALLEL DO
3898
3899	ENDIF ! microphysics test
3900
3901	!-----------------------------------------------------------
3902	! filter for moist variables post-microphysics and end of timestep
3903	!-----------------------------------------------------------
3904
3905	IF (config_flags%polar) THEN
3906	IF ( num_3d_m >= PARAM_FIRST_SCALAR ) THEN
3907	CALL wrf_debug ( 200 , ' call filter moist' )
3908	DO im = PARAM_FIRST_SCALAR, num_3d_m
3909	DO jj = jps, MIN(jpe,jde-1)
3910	DO kk = kps, MIN(kpe,kde-1)
3911	DO ii = ips, MIN(ipe,ide-1)
3912	moist(ii,kk,jj,im)=moist(ii,kk,jj,im)*(grid%mu_2(ii,jj)+grid%mub(ii,jj))
3913	ENDDO
3914	ENDDO
3915	ENDDO
3916
3917	CALL pxft ( grid=grid &
3918	,lineno=__LINE__ &
3919	,flag_uv = 0 &
3920	,flag_rurv = 0 &
3921	,flag_wph = 0 &
3922	,flag_ww = 0 &
3923	,flag_t = 0 &
3924	,flag_mu = 0 &
3925	,flag_mut = 0 &
3926	,flag_moist = im &
3927	,flag_chem = 0 &
3928	,flag_tracer = 0 &
3929	,flag_scalar = 0 &
3930	,positive_definite=.FALSE. &
3931	,moist=moist,chem=chem,tracer=tracer,scalar=scalar &
3932	,fft_filter_lat = config_flags%fft_filter_lat &
3933	,dclat = dclat &
3934	,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
3935	,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
3936	,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe &
3937	,imsx=imsx,imex=imex,jmsx=jmsx,jmex=jmex,kmsx=kmsx,kmex=kmex &
3938	,ipsx=ipsx,ipex=ipex,jpsx=jmsx,jpex=jpex,kpsx=kpsx,kpex=kpex )
3939
3940	DO jj = jps, MIN(jpe,jde-1)
3941	DO kk = kps, MIN(kpe,kde-1)
3942	DO ii = ips, MIN(ipe,ide-1)
3943	moist(ii,kk,jj,im)=moist(ii,kk,jj,im)/(grid%mu_2(ii,jj)+grid%mub(ii,jj))
3944	ENDDO
3945	ENDDO
3946	ENDDO
3947	ENDDO
3948	ENDIF
3949	ENDIF
3950
3951	!-----------------------------------------------------------
3952	! end filter for moist variables post-microphysics and end of timestep
3953	!-----------------------------------------------------------
3954
3955	!$OMP PARALLEL DO &
3956	!$OMP PRIVATE ( ij, its, ite, jts, jte, im, ii, jj, kk )
3957	scalar_tile_loop_1ba: DO ij = 1 , grid%num_tiles
3958
3959	IF ( config_flags%periodic_x ) THEN
3960	its = max(grid%i_start(ij),ids)
3961	ite = min(grid%i_end(ij),ide-1)
3962	ELSE
3963	its = max(grid%i_start(ij),ids+sz)
3964	ite = min(grid%i_end(ij),ide-1-sz)
3965	ENDIF
3966	jts = max(grid%j_start(ij),jds+sz)
3967	jte = min(grid%j_end(ij),jde-1-sz)
3968
3969	CALL calc_p_rho_phi( moist, num_3d_m, &
3970	grid%al, grid%alb, grid%mu_2, grid%muts, &
3971	grid%ph_2, grid%p, grid%pb, grid%t_2, &
3972	p0, t0, grid%znu, grid%dnw, grid%rdnw, &
3973	grid%rdn, config_flags%non_hydrostatic, &
3974	ids, ide, jds, jde, kds, kde, &
3975	ims, ime, jms, jme, kms, kme, &
3976	its, ite, jts, jte, &
3977	k_start , k_end )
3978
3979	END DO scalar_tile_loop_1ba
3980	!$OMP END PARALLEL DO
3981	BENCH_END(moist_phys_end_tim)
3982
3983	IF (.not. config_flags%non_hydrostatic) THEN
3984	#ifdef DM_PARALLEL
3985	# include "HALO_EM_HYDRO_UV.inc"
3986	# include "PERIOD_EM_HYDRO_UV.inc"
3987	#endif
3988	!$OMP PARALLEL DO &
3989	!$OMP PRIVATE ( ij )
3990	DO ij = 1 , grid%num_tiles
3991	CALL diagnose_w( ph_tend, grid%ph_2, grid%ph_1, grid%w_2, grid%muts, dt_rk, &
3992	grid%u_2, grid%v_2, grid%ht, &
3993	grid%cf1, grid%cf2, grid%cf3, grid%rdx, grid%rdy, grid%msftx, grid%msfty, &
3994	ids, ide, jds, jde, kds, kde, &
3995	ims, ime, jms, jme, kms, kme, &
3996	grid%i_start(ij), grid%i_end(ij), &
3997	grid%j_start(ij), grid%j_end(ij), &
3998	k_start , k_end )
3999
4000	END DO
4001	!$OMP END PARALLEL DO
4002
4003	END IF
4004
4005	CALL wrf_debug ( 200 , ' call chem polar filter ' )
4006
4007	!-----------------------------------------------------------
4008	! filter for chem and scalar variables at end of timestep
4009	!-----------------------------------------------------------
4010
4011	IF (config_flags%polar) THEN
4012
4013	IF ( num_3d_c >= PARAM_FIRST_SCALAR ) then
4014	chem_filter_loop: DO im = PARAM_FIRST_SCALAR, num_3d_c
4015	DO jj = jps, MIN(jpe,jde-1)
4016	DO kk = kps, MIN(kpe,kde-1)
4017	DO ii = ips, MIN(ipe,ide-1)
4018	chem(ii,kk,jj,im)=chem(ii,kk,jj,im)*(grid%mu_2(ii,jj)+grid%mub(ii,jj))
4019	ENDDO
4020	ENDDO
4021	ENDDO
4022
4023	CALL pxft ( grid=grid &
4024	,lineno=__LINE__ &
4025	,flag_uv = 0 &
4026	,flag_rurv = 0 &
4027	,flag_wph = 0 &
4028	,flag_ww = 0 &
4029	,flag_t = 0 &
4030	,flag_mu = 0 &
4031	,flag_mut = 0 &
4032	,flag_moist = 0 &
4033	,flag_chem = im &
4034	,flag_tracer = 0 &
4035	,flag_scalar = 0 &
4036	,positive_definite=.FALSE. &
4037	,moist=moist,chem=chem,tracer=tracer,scalar=scalar &
4038	,fft_filter_lat = config_flags%fft_filter_lat &
4039	,dclat = dclat &
4040	,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
4041	,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
4042	,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe &
4043	,imsx=imsx,imex=imex,jmsx=jmsx,jmex=jmex,kmsx=kmsx,kmex=kmex &
4044	,ipsx=ipsx,ipex=ipex,jpsx=jmsx,jpex=jpex,kpsx=kpsx,kpex=kpex )
4045
4046	DO jj = jps, MIN(jpe,jde-1)
4047	DO kk = kps, MIN(kpe,kde-1)
4048	DO ii = ips, MIN(ipe,ide-1)
4049	chem(ii,kk,jj,im)=chem(ii,kk,jj,im)/(grid%mu_2(ii,jj)+grid%mub(ii,jj))
4050	ENDDO
4051	ENDDO
4052	ENDDO
4053	ENDDO chem_filter_loop
4054	ENDIF
4055	IF ( num_tracer >= PARAM_FIRST_SCALAR ) then
4056	tracer_filter_loop: DO im = PARAM_FIRST_SCALAR, num_tracer
4057	DO jj = jps, MIN(jpe,jde-1)
4058	DO kk = kps, MIN(kpe,kde-1)
4059	DO ii = ips, MIN(ipe,ide-1)
4060	tracer(ii,kk,jj,im)=tracer(ii,kk,jj,im)*(grid%mu_2(ii,jj)+grid%mub(ii,jj))
4061	ENDDO
4062	ENDDO
4063	ENDDO
4064
4065	CALL pxft ( grid=grid &
4066	,lineno=__LINE__ &
4067	,flag_uv = 0 &
4068	,flag_rurv = 0 &
4069	,flag_wph = 0 &
4070	,flag_ww = 0 &
4071	,flag_t = 0 &
4072	,flag_mu = 0 &
4073	,flag_mut = 0 &
4074	,flag_moist = 0 &
4075	,flag_chem = 0 &
4076	,flag_tracer = im &
4077	,flag_scalar = 0 &
4078	,positive_definite=.FALSE. &
4079	,moist=moist,chem=chem,tracer=tracer,scalar=scalar &
4080	,fft_filter_lat = config_flags%fft_filter_lat &
4081	,dclat = dclat &
4082	,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
4083	,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
4084	,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe &
4085	,imsx=imsx,imex=imex,jmsx=jmsx,jmex=jmex,kmsx=kmsx,kmex=kmex &
4086	,ipsx=ipsx,ipex=ipex,jpsx=jmsx,jpex=jpex,kpsx=kpsx,kpex=kpex )
4087
4088	DO jj = jps, MIN(jpe,jde-1)
4089	DO kk = kps, MIN(kpe,kde-1)
4090	DO ii = ips, MIN(ipe,ide-1)
4091	tracer(ii,kk,jj,im)=tracer(ii,kk,jj,im)/(grid%mu_2(ii,jj)+grid%mub(ii,jj))
4092	ENDDO
4093	ENDDO
4094	ENDDO
4095	ENDDO tracer_filter_loop
4096	ENDIF
4097
4098	IF ( num_3d_s >= PARAM_FIRST_SCALAR ) then
4099	scalar_filter_loop: DO im = PARAM_FIRST_SCALAR, num_3d_s
4100	DO jj = jps, MIN(jpe,jde-1)
4101	DO kk = kps, MIN(kpe,kde-1)
4102	DO ii = ips, MIN(ipe,ide-1)
4103	scalar(ii,kk,jj,im)=scalar(ii,kk,jj,im)*(grid%mu_2(ii,jj)+grid%mub(ii,jj))
4104	ENDDO
4105	ENDDO
4106	ENDDO
4107
4108	CALL pxft ( grid=grid &
4109	,lineno=__LINE__ &
4110	,flag_uv = 0 &
4111	,flag_rurv = 0 &
4112	,flag_wph = 0 &
4113	,flag_ww = 0 &
4114	,flag_t = 0 &
4115	,flag_mu = 0 &
4116	,flag_mut = 0 &
4117	,flag_moist = 0 &
4118	,flag_chem = 0 &
4119	,flag_tracer = 0 &
4120	,flag_scalar = im &
4121	,positive_definite=.FALSE. &
4122	,moist=moist,chem=chem,tracer=tracer,scalar=scalar &
4123	,fft_filter_lat = config_flags%fft_filter_lat &
4124	,dclat = dclat &
4125	,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
4126	,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
4127	,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe &
4128	,imsx=imsx,imex=imex,jmsx=jmsx,jmex=jmex,kmsx=kmsx,kmex=kmex &
4129	,ipsx=ipsx,ipex=ipex,jpsx=jmsx,jpex=jpex,kpsx=kpsx,kpex=kpex )
4130
4131	DO jj = jps, MIN(jpe,jde-1)
4132	DO kk = kps, MIN(kpe,kde-1)
4133	DO ii = ips, MIN(ipe,ide-1)
4134	scalar(ii,kk,jj,im)=scalar(ii,kk,jj,im)/(grid%mu_2(ii,jj)+grid%mub(ii,jj))
4135	ENDDO
4136	ENDDO
4137	ENDDO
4138	ENDDO scalar_filter_loop
4139	ENDIF
4140	ENDIF
4141
4142	!-----------------------------------------------------------
4143	! end filter for chem and scalar variables at end of timestep
4144	!-----------------------------------------------------------
4145
4146	! We're finished except for boundary condition (and patch) update
4147
4148	! Boundary condition time (or communication time). At this time, we have
4149	! implemented periodic and symmetric physical boundary conditions.
4150
4151	! b.c. routine for data within patch.
4152
4153	! we need to do both time levels of
4154	! data because the time filter only works in the physical solution space.
4155
4156	! First, do patch communications for boundary conditions (periodicity)
4157
4158	!-----------------------------------------------------------
4159	! Stencils for patch communications (WCS, 29 June 2001)
4160	!
4161	! here's where we need a wide comm stencil - these are the
4162	! uncoupled variables so are used for high order calc in
4163	! advection and mixong routines.
4164	!
4165	! * * * * *
4166	! * * * * * * * * *
4167	! * + * * + * * * + * *
4168	! * * * * * * * * *
4169	! * * * * *
4170	!
4171	! grid%u_1 x
4172	! grid%u_2 x
4173	! grid%v_1 x
4174	! grid%v_2 x
4175	! grid%w_1 x
4176	! grid%w_2 x
4177	! grid%t_1 x
4178	! grid%t_2 x
4179	! grid%ph_1 x
4180	! grid%ph_2 x
4181	! grid%tke_1 x
4182	! grid%tke_2 x
4183	!
4184	! 2D variables
4185	! grid%mu_1 x
4186	! grid%mu_2 x
4187	!
4188	! 4D variables
4189	! moist x
4190	! chem x
4191	! scalar x
4192	!----------------------------------------------------------
4193
4194
4195	#ifdef DM_PARALLEL
4196	IF ( config_flags%h_mom_adv_order <= 4 ) THEN
4197	# include "HALO_EM_D3_3.inc"
4198	ELSE IF ( config_flags%h_mom_adv_order <= 6 ) THEN
4199	# include "HALO_EM_D3_5.inc"
4200	ELSE
4201	WRITE(wrf_err_message,*)'solve_em: invalid h_mom_adv_order = ',config_flags%h_mom_adv_order
4202	CALL wrf_error_fatal(TRIM(wrf_err_message))
4203	ENDIF
4204	# include "PERIOD_BDY_EM_D3.inc"
4205	# include "PERIOD_BDY_EM_MOIST.inc"
4206	# include "PERIOD_BDY_EM_CHEM.inc"
4207	# include "PERIOD_BDY_EM_TRACER.inc"
4208	# include "PERIOD_BDY_EM_SCALAR.inc"
4209	#endif
4210
4211	! now set physical b.c on a patch
4212
4213	BENCH_START(bc_2d_tim)
4214	!$OMP PARALLEL DO &
4215	!$OMP PRIVATE ( ij )
4216	tile_bc_loop_2: DO ij = 1 , grid%num_tiles
4217
4218	CALL wrf_debug ( 200 , ' call set_phys_bc_dry_2' )
4219
4220	CALL set_phys_bc_dry_2( config_flags, &
4221	grid%u_1, grid%u_2, grid%v_1, grid%v_2, grid%w_1, grid%w_2, &
4222	grid%t_1, grid%t_2, grid%ph_1, grid%ph_2, grid%mu_1, grid%mu_2, &
4223	ids, ide, jds, jde, kds, kde, &
4224	ims, ime, jms, jme, kms, kme, &
4225	ips, ipe, jps, jpe, kps, kpe, &
4226	grid%i_start(ij), grid%i_end(ij), &
4227	grid%j_start(ij), grid%j_end(ij), &
4228	k_start , k_end )
4229
4230	CALL set_physical_bc3d( grid%tke_1, 'p', config_flags, &
4231	ids, ide, jds, jde, kds, kde, &
4232	ims, ime, jms, jme, kms, kme, &
4233	ips, ipe, jps, jpe, kps, kpe, &
4234	grid%i_start(ij), grid%i_end(ij), &
4235	grid%j_start(ij), grid%j_end(ij), &
4236	k_start , k_end-1 )
4237
4238	CALL set_physical_bc3d( grid%tke_2 , 'p', config_flags, &
4239	ids, ide, jds, jde, kds, kde, &
4240	ims, ime, jms, jme, kms, kme, &
4241	ips, ipe, jps, jpe, kps, kpe, &
4242	grid%i_start(ij), grid%i_end(ij), &
4243	grid%j_start(ij), grid%j_end(ij), &
4244	k_start , k_end )
4245
4246	moisture_loop_bdy_2 : DO im = PARAM_FIRST_SCALAR , num_3d_m
4247
4248	CALL set_physical_bc3d( moist(ims,kms,jms,im), 'p', &
4249	config_flags, &
4250	ids, ide, jds, jde, kds, kde, &
4251	ims, ime, jms, jme, kms, kme, &
4252	ips, ipe, jps, jpe, kps, kpe, &
4253	grid%i_start(ij), grid%i_end(ij), &
4254	grid%j_start(ij), grid%j_end(ij), &
4255	k_start , k_end )
4256
4257	END DO moisture_loop_bdy_2
4258
4259	chem_species_bdy_loop_2 : DO ic = PARAM_FIRST_SCALAR , num_3d_c
4260
4261	CALL set_physical_bc3d( chem(ims,kms,jms,ic) , 'p', config_flags, &
4262	ids, ide, jds, jde, kds, kde, &
4263	ims, ime, jms, jme, kms, kme, &
4264	ips, ipe, jps, jpe, kps, kpe, &
4265	grid%i_start(ij), grid%i_end(ij), &
4266	grid%j_start(ij), grid%j_end(ij), &
4267	k_start , k_end )
4268
4269	END DO chem_species_bdy_loop_2
4270
4271	tracer_species_bdy_loop_2 : DO ic = PARAM_FIRST_SCALAR , num_tracer
4272
4273	CALL set_physical_bc3d( tracer(ims,kms,jms,ic) , 'p', config_flags, &
4274	ids, ide, jds, jde, kds, kde, &
4275	ims, ime, jms, jme, kms, kme, &
4276	ips, ipe, jps, jpe, kps, kpe, &
4277	grid%i_start(ij), grid%i_end(ij), &
4278	grid%j_start(ij), grid%j_end(ij), &
4279	k_start , k_end )
4280
4281	END DO tracer_species_bdy_loop_2
4282
4283	scalar_species_bdy_loop_2 : DO is = PARAM_FIRST_SCALAR , num_3d_s
4284
4285	CALL set_physical_bc3d( scalar(ims,kms,jms,is) , 'p', config_flags, &
4286	ids, ide, jds, jde, kds, kde, &
4287	ims, ime, jms, jme, kms, kme, &
4288	ips, ipe, jps, jpe, kps, kpe, &
4289	grid%i_start(ij), grid%i_end(ij), &
4290	grid%j_start(ij), grid%j_end(ij), &
4291	k_start , k_end )
4292
4293	END DO scalar_species_bdy_loop_2
4294
4295	END DO tile_bc_loop_2
4296	!$OMP END PARALLEL DO
4297	BENCH_END(bc_2d_tim)
4298
4299	IF( config_flags%specified .or. config_flags%nested ) THEN
4300	grid%dtbc = grid%dtbc + grid%dt
4301	ENDIF
4302
4303	! reset surface w for consistency
4304
4305	#ifdef DM_PARALLEL
4306	# include "HALO_EM_C.inc"
4307	# include "PERIOD_BDY_EM_E.inc"
4308	#endif
4309
4310	CALL wrf_debug ( 10 , ' call set_w_surface' )
4311	fill_w_flag = .false.
4312
4313	!$OMP PARALLEL DO &
4314	!$OMP PRIVATE ( ij )
4315	DO ij = 1 , grid%num_tiles
4316	CALL set_w_surface( config_flags, grid%znw, fill_w_flag, &
4317	grid%w_2, grid%ht, grid%u_2, grid%v_2, &
4318	grid%cf1, grid%cf2, grid%cf3, grid%rdx, grid%rdy,&
4319	grid%msftx, grid%msfty, &
4320	ids, ide, jds, jde, kds, kde, &
4321	ims, ime, jms, jme, kms, kme, &
4322	grid%i_start(ij), grid%i_end(ij), &
4323	grid%j_start(ij), grid%j_end(ij), &
4324	k_start, k_end )
4325	! its, ite, jts, jte, k_start, min(k_end,kde-1), &
4326
4327	END DO
4328	!$OMP END PARALLEL DO
4329
4330	! calculate some model diagnostics.
4331
4332	#ifdef LMDZ1
4333	WRITE(message, *)' dyn_em: before diagnostics'
4334	CALL wrf_debug(200, message)
4335	WRITE(message, *)' t_tend: ',t_tendf(im2,km2,jm2), &
4336	' u_tend: ', ru_tendf(im2,1,jm2)
4337	CALL wrf_debug(200, message)
4338	WRITE(message,*)' psfc_tend: ',grid%dpsdt(im2,jm2), &
4339	'p sfc: ',p8w(im2,kms,jm2)
4340	CALL wrf_debug(200, message)
4341	WRITE(message,*)' p 1: ',grid%p(im2,kms,jm2), ' ph 1: ',grid%ph_2(im2,kms,jm2)
4342	CALL wrf_debug(200, message)
4343	#endif
4344	CALL wrf_debug ( 200 , ' call diagnostic_driver' )
4345
4346	CALL diagnostic_output_calc( &
4347	& DPSDT=grid%dpsdt ,DMUDT=grid%dmudt &
4348	& ,P8W=p8w ,PK1M=grid%pk1m &
4349	& ,MU_2=grid%mu_2 ,MU_2M=grid%mu_2m &
4350	& ,U=grid%u_2 ,V=grid%v_2 &
4351	& ,RAINCV=grid%raincv ,RAINNCV=grid%rainncv &
4352	& ,RAINC=grid%rainc ,RAINNC=grid%rainnc &
4353	& ,I_RAINC=grid%i_rainc ,I_RAINNC=grid%i_rainnc &
4354	& ,HFX=grid%hfx ,SFCEVP=grid%sfcevp ,LH=grid%lh &
4355	& ,DT=grid%dt ,SBW=config_flags%spec_bdy_width &
4356	& ,XTIME=grid%xtime ,T2=grid%t2 &
4357	& ,ACSWUPT=grid%acswupt ,ACSWUPTC=grid%acswuptc &
4358	& ,ACSWDNT=grid%acswdnt ,ACSWDNTC=grid%acswdntc &
4359	& ,ACSWUPB=grid%acswupb ,ACSWUPBC=grid%acswupbc &
4360	& ,ACSWDNB=grid%acswdnb ,ACSWDNBC=grid%acswdnbc &
4361	& ,ACLWUPT=grid%aclwupt ,ACLWUPTC=grid%aclwuptc &
4362	& ,ACLWDNT=grid%aclwdnt ,ACLWDNTC=grid%aclwdntc &
4363	& ,ACLWUPB=grid%aclwupb ,ACLWUPBC=grid%aclwupbc &
4364	& ,ACLWDNB=grid%aclwdnb ,ACLWDNBC=grid%aclwdnbc &
4365	& ,I_ACSWUPT=grid%i_acswupt ,I_ACSWUPTC=grid%i_acswuptc &
4366	& ,I_ACSWDNT=grid%i_acswdnt ,I_ACSWDNTC=grid%i_acswdntc &
4367	& ,I_ACSWUPB=grid%i_acswupb ,I_ACSWUPBC=grid%i_acswupbc &
4368	& ,I_ACSWDNB=grid%i_acswdnb ,I_ACSWDNBC=grid%i_acswdnbc &
4369	& ,I_ACLWUPT=grid%i_aclwupt ,I_ACLWUPTC=grid%i_aclwuptc &
4370	& ,I_ACLWDNT=grid%i_aclwdnt ,I_ACLWDNTC=grid%i_aclwdntc &
4371	& ,I_ACLWUPB=grid%i_aclwupb ,I_ACLWUPBC=grid%i_aclwupbc &
4372	& ,I_ACLWDNB=grid%i_aclwdnb ,I_ACLWDNBC=grid%i_aclwdnbc &
4373	! Selection flag
4374	& ,DIAG_PRINT=config_flags%diag_print &
4375	& ,BUCKET_MM=config_flags%bucket_mm &
4376	& ,BUCKET_J =config_flags%bucket_J &
4377	& ,SNOWNCV=grid%snowncv, SNOW_ACC_NC=grid%snow_acc_nc &
4378	& ,PREC_ACC_C=grid%prec_acc_c &
4379	& ,PREC_ACC_NC=grid%prec_acc_nc &
4380	& ,PREC_ACC_DT=config_flags%prec_acc_dt &
4381	& ,CURR_SECS=curr_secs &
4382	! Dimension arguments
4383	& ,IDS=ids,IDE=ide, JDS=jds,JDE=jde, KDS=kds,KDE=kde &
4384	& ,IMS=ims,IME=ime, JMS=jms,JME=jme, KMS=kms,KME=kme &
4385	& ,IPS=ips,IPE=ipe, JPS=jps,JPE=jpe, KPS=kps,KPE=kpe &
4386	& ,I_START=grid%i_start,I_END=min(grid%i_end, ide-1) &
4387	& ,J_START=grid%j_start,J_END=min(grid%j_end, jde-1) &
4388	& ,KTS=k_start, KTE=min(k_end,kde-1) &
4389	& ,NUM_TILES=grid%num_tiles &
4390	& )
4391
4392	#ifdef DM_PARALLEL
4393	!-----------------------------------------------------------------------
4394	! see above
4395	!--------------------------------------------------------------
4396	CALL wrf_debug ( 200 , ' call HALO_RK_E' )
4397	IF ( config_flags%h_mom_adv_order <= 4 ) THEN
4398	# include "HALO_EM_E_3.inc"
4399	ELSE IF ( config_flags%h_mom_adv_order <= 6 ) THEN
4400	# include "HALO_EM_E_5.inc"
4401	ELSE
4402	WRITE(wrf_err_message,*)'solve_em: invalid h_mom_adv_order = ',config_flags%h_mom_adv_order
4403	CALL wrf_error_fatal(TRIM(wrf_err_message))
4404	ENDIF
4405	#endif
4406
4407	#ifdef DM_PARALLEL
4408	IF ( num_moist >= PARAM_FIRST_SCALAR ) THEN
4409	!-----------------------------------------------------------------------
4410	! see above
4411	!--------------------------------------------------------------
4412	CALL wrf_debug ( 200 , ' call HALO_RK_MOIST' )
4413	IF ( config_flags%h_mom_adv_order <= 4 ) THEN
4414	# include "HALO_EM_MOIST_E_3.inc"
4415	ELSE IF ( config_flags%h_mom_adv_order <= 6 ) THEN
4416	# include "HALO_EM_MOIST_E_5.inc"
4417	ELSE
4418	WRITE(wrf_err_message,*)'solve_em: invalid h_mom_adv_order = ',config_flags%h_mom_adv_order
4419	CALL wrf_error_fatal(TRIM(wrf_err_message))
4420	ENDIF
4421	ENDIF
4422	IF ( num_chem >= PARAM_FIRST_SCALAR ) THEN
4423	!-----------------------------------------------------------------------
4424	! see above
4425	!--------------------------------------------------------------
4426	CALL wrf_debug ( 200 , ' call HALO_RK_CHEM' )
4427	IF ( config_flags%h_mom_adv_order <= 4 ) THEN
4428	# include "HALO_EM_CHEM_E_3.inc"
4429	ELSE IF ( config_flags%h_mom_adv_order <= 6 ) THEN
4430	# include "HALO_EM_CHEM_E_5.inc"
4431	ELSE
4432	WRITE(wrf_err_message,*)'solve_em: invalid h_mom_adv_order = ',config_flags%h_mom_adv_order
4433	CALL wrf_error_fatal(TRIM(wrf_err_message))
4434	ENDIF
4435	ENDIF
4436	IF ( num_tracer >= PARAM_FIRST_SCALAR ) THEN
4437	!-----------------------------------------------------------------------
4438	! see above
4439	!--------------------------------------------------------------
4440	CALL wrf_debug ( 200 , ' call HALO_RK_TRACER' )
4441	IF ( config_flags%h_mom_adv_order <= 4 ) THEN
4442	# include "HALO_EM_TRACER_E_3.inc"
4443	ELSE IF ( config_flags%h_mom_adv_order <= 6 ) THEN
4444	# include "HALO_EM_TRACER_E_5.inc"
4445	ELSE
4446	WRITE(wrf_err_message,*)'solve_em: invalid h_mom_adv_order = ',config_flags%h_mom_adv_order
4447	CALL wrf_error_fatal(TRIM(wrf_err_message))
4448	ENDIF
4449	ENDIF
4450	IF ( num_scalar >= PARAM_FIRST_SCALAR ) THEN
4451	!-----------------------------------------------------------------------
4452	! see above
4453	!--------------------------------------------------------------
4454	CALL wrf_debug ( 200 , ' call HALO_RK_SCALAR' )
4455	IF ( config_flags%h_mom_adv_order <= 4 ) THEN
4456	# include "HALO_EM_SCALAR_E_3.inc"
4457	ELSE IF ( config_flags%h_mom_adv_order <= 6 ) THEN
4458	# include "HALO_EM_SCALAR_E_5.inc"
4459	ELSE
4460	WRITE(wrf_err_message,*)'solve_em: invalid h_mom_adv_order = ',config_flags%h_mom_adv_order
4461	CALL wrf_error_fatal(TRIM(wrf_err_message))
4462	ENDIF
4463	ENDIF
4464	#endif
4465
4466	! Max values of CFL for adaptive time step scheme
4467
4468	DEALLOCATE(max_vert_cfl_tmp)
4469	DEALLOCATE(max_horiz_cfl_tmp)
4470
4471	CALL wrf_debug ( 200 , ' call end of solve_em' )
4472
4473	! Finish timers if compiled with -DBENCH.
4474	#include <bench_solve_em_end.h>
4475
4476	RETURN
4477
4478	END SUBROUTINE solve_em

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

Download in other formats:

Original Format