Context Navigation

solve_em.F @ 1

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

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

WRF: version v3.3
LMDZ: version v1818

More details in:

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

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

Download in other formats:

Original Format