1 !WRF:MEDIATION_LAYER:SOLVER
3 SUBROUTINE solve_em ( grid , config_flags &
4 ! Arguments generated from Registry
5 #include "dummy_new_args.inc"
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
18 USE module_tiles, ONLY : set_tiles
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
48 ! Mediation layer modules
50 USE module_model_constants
51 USE module_small_step_em
53 USE module_big_step_utilities_em
56 USE module_solvedebug_em
57 USE module_physics_addtendc
58 USE module_diffusion_em
60 USE module_microphysics_driver
61 USE module_microphysics_zero_out
62 USE module_fddaobs_driver
63 USE module_diagnostics
65 USE module_input_chem_data
66 USE module_input_tracer
67 USE module_chem_utilities
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
78 TYPE(domain) , TARGET :: grid
80 ! Definitions of dummy arguments to this routine (generated from Registry).
81 #include "dummy_new_decl.inc"
83 ! Structure that contains run-time configuration (namelist) data for domain
84 TYPE (grid_config_rec_type) , INTENT(IN) :: config_flags
88 INTEGER :: k_start , k_end, its, ite, jts, jte
89 INTEGER :: ids , ide , jds , jde , kds , kde , &
90 ims , ime , jms , jme , kms , kme , &
91 ips , ipe , jps , jpe , kps , kpe
93 INTEGER :: sids , side , sjds , sjde , skds , skde , &
94 sims , sime , sjms , sjme , skms , skme , &
95 sips , sipe , sjps , sjpe , skps , skpe
98 INTEGER :: imsx, imex, jmsx, jmex, kmsx, kmex, &
99 ipsx, ipex, jpsx, jpex, kpsx, kpex, &
100 imsy, imey, jmsy, jmey, kmsy, kmey, &
101 ipsy, ipey, jpsy, jpey, kpsy, kpey
103 INTEGER :: ij , iteration
104 INTEGER :: im , num_3d_m , ic , num_3d_c , is , num_3d_s
109 LOGICAL :: specified_bdy, channel_bdy
113 ! Changes in tendency at this timestep
114 real ,DIMENSION(grid%sm31:grid%em31,grid%sm32:grid%em32,grid%sm33:grid%em33) :: h_tendency, &
117 ! Whether advection should produce decoupled horizontal and vertical advective tendency outputs
121 ! Index cross-referencing array for tendency accumulation
122 INTEGER, DIMENSION( num_chem ) :: adv_ct_indices
125 ! storage for tendencies and decoupled state (generated from Registry)
127 #include <i1_decl.inc>
128 ! Previous time level of tracer arrays now defined as i1 variables;
129 ! the state 4d arrays now redefined as 1-time level arrays in Registry.
130 ! Benefit: save memory in nested runs, since only 1 domain is active at a
131 ! time. Potential problem on stack-limited architectures: increases
132 ! amount of data on program stack by making these automatic arrays.
135 INTEGER :: number_of_small_timesteps, rk_step
136 INTEGER :: klevel,ijm,ijp,i,j,k,size1,size2 ! for prints/plots only
137 INTEGER :: idum1, idum2, dynamics_option
139 INTEGER :: rk_order, iwmax, jwmax, kwmax
140 REAL :: dt_rk, dts_rk, dts, dtm, wmax
141 REAL , ALLOCATABLE , DIMENSION(:) :: max_vert_cfl_tmp, max_horiz_cfl_tmp
144 LOGICAL :: f_flux ! flag for computing averaged fluxes in cu_gd
146 INTEGER :: num_sound_steps
147 INTEGER :: idex, jdex
151 INTEGER :: ii, jj !kk is above after l,kte
153 INTEGER :: debug_level
155 ! urban related variables
156 INTEGER :: NUM_ROOF_LAYERS, NUM_WALL_LAYERS, NUM_ROAD_LAYERS ! urban
158 TYPE(WRFU_TimeInterval) :: tmpTimeInterval
160 LOGICAL :: adapt_step_flag
161 LOGICAL :: fill_w_flag
163 ! variables for flux-averaging code 20091223
164 CHARACTER*256 :: message, message2
166 TYPE(WRFU_Time) :: temp_time, CurrTime
167 INTEGER, PARAMETER :: precision = 100
169 TYPE(WRFU_TimeInterval) :: dtInterval
171 ! Define benchmarking timers if -DBENCH is compiled
172 #include <bench_solve_em_def.h>
174 !----------------------
175 ! Executable statements
176 !----------------------
180 ! solve_em is the main driver for advancing a grid a single timestep.
181 ! It is a mediation-layer routine -> DM and SM calls are made where
182 ! needed for parallel processing.
184 ! solve_em can integrate the equations using 3 time-integration methods
186 ! - 3rd order Runge-Kutta time integration (recommended)
188 ! - 2nd order Runge-Kutta time integration
190 ! The main sections of solve_em are
192 ! (1) Runge-Kutta (RK) loop
194 ! (2) Non-timesplit physics (i.e., tendencies computed for updating
195 ! model state variables during the first RK sub-step (loop)
197 ! (3) Small (acoustic, sound) timestep loop - within the RK sub-steps
199 ! (4) scalar advance for moist and chem scalar variables (and TKE)
200 ! within the RK sub-steps.
202 ! (5) time-split physics (after the RK step), currently this includes
205 ! A more detailed description of these sections follows.
209 ! Initialize timers if compiled with -DBENCH
210 #include <bench_solve_em_init.h>
212 ! set runge-kutta solver (2nd or 3rd order)
214 dynamics_option = config_flags%rk_ord
216 ! Obtain dimension information stored in the grid data structure.
218 CALL get_ijk_from_grid ( grid , &
219 ids, ide, jds, jde, kds, kde, &
220 ims, ime, jms, jme, kms, kme, &
221 ips, ipe, jps, jpe, kps, kpe, &
222 imsx, imex, jmsx, jmex, kmsx, kmex, &
223 ipsx, ipex, jpsx, jpex, kpsx, kpex, &
224 imsy, imey, jmsy, jmey, kmsy, kmey, &
225 ipsy, ipey, jpsy, jpey, kpsy, kpey )
227 CALL get_ijk_from_subgrid ( grid , &
228 sids, side, sjds, sjde, skds, skde, &
229 sims, sime, sjms, sjme, skms, skme, &
230 sips, sipe, sjps, sjpe, skps, skpe )
236 num_3d_s = num_scalar
238 f_flux = config_flags%do_avgflx_cugd .EQ. 1
240 ! Compute these starting and stopping locations for each tile and number of tiles.
241 ! See: http://www.mmm.ucar.edu/wrf/WG2/topics/settiles
242 CALL set_tiles ( grid , ids , ide , jds , jde , ips , ipe , jps , jpe )
244 ! Max values of CFL for adaptive time step scheme
246 ALLOCATE (max_vert_cfl_tmp(grid%num_tiles))
247 ALLOCATE (max_horiz_cfl_tmp(grid%num_tiles))
250 ! Calculate current time in seconds since beginning of model run.
251 ! Unfortunately, ESMF does not seem to have a way to return
252 ! floating point seconds based on a TimeInterval. So, we will
253 ! calculate it here--but, this is not clean!!
255 tmpTimeInterval = domain_get_current_time ( grid ) - domain_get_sim_start_time ( grid )
256 curr_secs = real_time(tmpTimeInterval)
258 old_dt = grid%dt ! store old time step for flux averaging code at end of RK loop
259 !-----------------------------------------------------------------------------
260 ! Adaptive time step: Added by T. Hutchinson, WSI 3/5/07
261 ! In this call, we do the time-step adaptation and set time-dependent lateral
262 ! boundary condition nudging weights.
264 IF ( (config_flags%use_adaptive_time_step) .and. &
265 ( (.not. grid%nested) .or. &
266 ( (grid%nested) .and. (abs(grid%dtbc) < 0.0001) ) ) )THEN
267 CALL adapt_timestep(grid, config_flags)
268 adapt_step_flag = .TRUE.
270 adapt_step_flag = .FALSE.
272 ! End of adaptive time step modifications
273 !-----------------------------------------------------------------------------
275 grid%itimestep = grid%itimestep + 1
277 IF (config_flags%polar) dclat = 90./REAL(jde-jds) !(0.5 * 180/ny)
283 if ( num_chem >= PARAM_FIRST_SCALAR ) then
284 !-----------------------------------------------------------------------
285 ! see matching halo calls below for stencils
286 !--------------------------------------------------------------
287 CALL wrf_debug ( 200 , ' call HALO_RK_CHEM' )
288 IF ( config_flags%h_mom_adv_order <= 4 ) THEN
289 # include "HALO_EM_CHEM_E_3.inc"
290 IF( config_flags%progn > 0 ) THEN
291 # include "HALO_EM_SCALAR_E_3.inc"
293 ELSE IF ( config_flags%h_mom_adv_order <= 6 ) THEN
294 # include "HALO_EM_CHEM_E_5.inc"
295 IF( config_flags%progn > 0 ) THEN
296 # include "HALO_EM_SCALAR_E_5.inc"
299 WRITE(wrf_err_message,*)'solve_em: invalid h_mom_adv_order = ',config_flags%h_mom_adv_order
300 CALL wrf_error_fatal(TRIM(wrf_err_message))
303 if ( num_tracer >= PARAM_FIRST_SCALAR ) then
304 !-----------------------------------------------------------------------
305 ! see matching halo calls below for stencils
306 !--------------------------------------------------------------
307 CALL wrf_debug ( 200 , ' call HALO_RK_tracer' )
308 IF ( config_flags%h_mom_adv_order <= 4 ) THEN
309 # include "HALO_EM_TRACER_E_3.inc"
310 ELSE IF ( config_flags%h_mom_adv_order <= 6 ) THEN
311 # include "HALO_EM_TRACER_E_5.inc"
313 WRITE(wrf_err_message,*)'solve_em: invalid h_mom_adv_order = ',config_flags%h_mom_adv_order
314 CALL wrf_error_fatal(TRIM(wrf_err_message))
318 !--------------------------------------------------------------
319 adv_ct_indices( : ) = 1
320 IF ( config_flags%chemdiag == USECHEMDIAG ) THEN
321 ! modify tendency list here
322 ! note that the referencing direction here is opposite of that in chem_driver
323 adv_ct_indices(p_co ) = p_advh_co
324 adv_ct_indices(p_o3 ) = p_advh_o3
325 adv_ct_indices(p_no ) = p_advh_no
326 adv_ct_indices(p_no2 ) = p_advh_no2
327 adv_ct_indices(p_hno3) = p_advh_hno3
328 adv_ct_indices(p_iso ) = p_advh_iso
329 adv_ct_indices(p_ho ) = p_advh_ho
330 adv_ct_indices(p_ho2 ) = p_advh_ho2
334 rk_order = config_flags%rk_ord
336 IF ( grid%time_step_sound == 0 ) THEN
337 ! This function will give 4 for 6*dx and 6 for 10*dx and returns even numbers only
338 spacing = min(grid%dx, grid%dy)
339 IF ( ( config_flags%use_adaptive_time_step ) .AND. ( config_flags%map_proj == PROJ_CASSINI ) ) THEN
340 max_msft=MIN ( MAX(grid%max_msftx, grid%max_msfty) , &
341 1.0/COS(config_flags%fft_filter_lat*degrad) )
342 num_sound_steps = max ( 2 * ( INT (300. * grid%dt / (spacing / max_msft) - 0.01 ) + 1 ), 4 )
343 ELSE IF ( config_flags%use_adaptive_time_step ) THEN
344 max_msft= MAX(grid%max_msftx, grid%max_msfty)
345 num_sound_steps = max ( 2 * ( INT (300. * grid%dt / (spacing / max_msft) - 0.01 ) + 1 ), 4 )
347 num_sound_steps = max ( 2 * ( INT (300. * grid%dt / spacing - 0.01 ) + 1 ), 4 )
349 WRITE(wrf_err_message,*)'grid spacing, dt, time_step_sound=',spacing,grid%dt,num_sound_steps
350 CALL wrf_debug ( 50 , wrf_err_message )
352 num_sound_steps = grid%time_step_sound
355 dts = grid%dt/float(num_sound_steps)
357 IF (config_flags%use_adaptive_time_step) THEN
359 CALL get_wrf_debug_level( debug_level )
360 IF ((config_flags%time_step < 0) .AND. (debug_level.GE.50)) THEN
362 CALL wrf_dm_maxval(grid%max_vert_cfl, idex, jdex)
364 WRITE(wrf_err_message,*)'variable dt, max horiz cfl, max vert cfl: ',&
365 grid%dt, grid%max_horiz_cfl, grid%max_vert_cfl
366 CALL wrf_debug ( 0 , wrf_err_message )
370 grid%max_horiz_cfl = 0
371 grid%max_vert_cfl = 0
374 ! setting bdy tendencies to zero for DFI if constant_bc = true
378 DO ij = 1 , grid%num_tiles
380 ! IF( config_flags%specified .AND. grid%dfi_opt .NE. DFI_NODFI &
381 ! .AND. config_flags%constant_bc .AND. (grid%dfi_stage .EQ. DFI_BCK .OR. grid%dfi_stage .EQ. DFI_FWD) ) THEN
382 IF( config_flags%specified .AND. config_flags%constant_bc ) THEN
384 CALL zero_bdytend (grid%u_btxs,grid%u_btxe,grid%u_btys,grid%u_btye, &
385 grid%v_btxs,grid%v_btxe,grid%v_btys,grid%v_btye, &
386 grid%ph_btxs,grid%ph_btxe,grid%ph_btys,grid%ph_btye, &
387 grid%t_btxs,grid%t_btxe,grid%t_btys,grid%t_btye, &
388 grid%w_btxs,grid%w_btxe,grid%w_btys,grid%w_btye, &
389 grid%mu_btxs,grid%mu_btxe,grid%mu_btys,grid%mu_btye, &
390 moist_btxs,moist_btxe, &
391 moist_btys,moist_btye, &
392 grid%spec_bdy_width,num_3d_m, &
393 ids,ide, jds,jde, kds,kde, &
394 ims,ime, jms,jme, kms,kme, &
395 ips,ipe, jps,jpe, kps,kpe, &
396 grid%i_start(ij), grid%i_end(ij), &
397 grid%j_start(ij), grid%j_end(ij), &
402 !$OMP END PARALLEL DO
404 !**********************************************************************
406 ! LET US BEGIN.......
410 ! (1) RK integration loop is named the "Runge_Kutta_loop:"
412 ! Predictor-corrector type time integration.
413 ! Advection terms are evaluated at time t for the predictor step,
414 ! and advection is re-evaluated with the latest predicted value for
415 ! each succeeding time corrector step
417 ! 2nd order Runge Kutta (rk_order = 2):
418 ! Step 1 is taken to the midpoint predictor, step 2 is the full step.
420 ! 3rd order Runge Kutta (rk_order = 3):
421 ! Step 1 is taken to from t to dt/3, step 2 is from t to dt/2,
422 ! and step 3 is from t to dt.
424 ! non-timesplit physics are evaluated during first RK step and
425 ! these physics tendencies are stored for use in each RK pass.
428 !**********************************************************************
430 Runge_Kutta_loop: DO rk_step = 1, rk_order
432 ! Set the step size and number of small timesteps for
433 ! each part of the timestep
436 IF ( rk_order == 1 ) THEN
438 write(wrf_err_message,*)' leapfrog removed, error exit for dynamics_option = ',dynamics_option
439 CALL wrf_error_fatal( wrf_err_message )
441 ELSE IF ( rk_order == 2 ) THEN ! 2nd order Runge-Kutta timestep
443 IF ( rk_step == 1) THEN
446 number_of_small_timesteps = num_sound_steps/2
450 number_of_small_timesteps = num_sound_steps
453 ELSE IF ( rk_order == 3 ) THEN ! third order Runge-Kutta
455 IF ( rk_step == 1) THEN
458 number_of_small_timesteps = 1
459 ELSE IF (rk_step == 2) THEN
462 number_of_small_timesteps = num_sound_steps/2
466 number_of_small_timesteps = num_sound_steps
471 write(wrf_err_message,*)' unknown solver, error exit for dynamics_option = ',dynamics_option
472 CALL wrf_error_fatal( wrf_err_message )
476 ! Ensure that polar meridional velocity is zero
477 IF (config_flags%polar) THEN
480 DO ij = 1 , grid%num_tiles
481 CALL zero_pole ( grid%v_1, &
482 ids, ide, jds, jde, kds, kde, &
483 ims, ime, jms, jme, kms, kme, &
484 grid%i_start(ij), grid%i_end(ij), &
485 grid%j_start(ij), grid%j_end(ij), &
487 CALL zero_pole ( grid%v_2, &
488 ids, ide, jds, jde, kds, kde, &
489 ims, ime, jms, jme, kms, kme, &
490 grid%i_start(ij), grid%i_end(ij), &
491 grid%j_start(ij), grid%j_end(ij), &
494 !$OMP END PARALLEL DO
497 ! Time level t is in the *_2 variable in the first part
498 ! of the step, and in the *_1 variable after the predictor.
499 ! the latest predicted values are stored in the *_2 variables.
501 CALL wrf_debug ( 200 , ' call rk_step_prep ' )
503 BENCH_START(step_prep_tim)
507 DO ij = 1 , grid%num_tiles
509 CALL rk_step_prep ( config_flags, rk_step, &
510 grid%u_2, grid%v_2, grid%w_2, grid%t_2, grid%ph_2, grid%mu_2, &
512 grid%ru, grid%rv, grid%rw, grid%ww, grid%php, grid%alt, grid%muu, grid%muv, &
513 grid%mub, grid%mut, grid%phb, grid%pb, grid%p, grid%al, grid%alb, &
515 grid%msfux, grid%msfuy, grid%msfvx, grid%msfvx_inv, &
516 grid%msfvy, grid%msftx, grid%msfty, &
517 grid%fnm, grid%fnp, grid%dnw, grid%rdx, grid%rdy, &
519 ids, ide, jds, jde, kds, kde, &
520 ims, ime, jms, jme, kms, kme, &
521 grid%i_start(ij), grid%i_end(ij), &
522 grid%j_start(ij), grid%j_end(ij), &
526 !$OMP END PARALLEL DO
527 BENCH_END(step_prep_tim)
530 !-----------------------------------------------------------------------
531 ! Stencils for patch communications (WCS, 29 June 2001)
532 ! Note: the small size of this halo exchange reflects the
533 ! fact that we are carrying the uncoupled variables
534 ! as state variables in the mass coordinate model, as
535 ! opposed to the coupled variables as in the height
540 ! * + * * + * * * + * *
544 ! 3D variables - note staggering! ru(X), rv(Y), ww(Z), php(Z)
554 ! the following are 2D (xy) variables
559 !--------------------------------------------------------------
560 # include "HALO_EM_A.inc"
563 ! set boundary conditions on variables
564 ! from big_step_prep for use in big_step_proc
567 # include "PERIOD_BDY_EM_A.inc"
570 BENCH_START(set_phys_bc_tim)
572 !$OMP PRIVATE ( ij, ii, jj, kk )
574 DO ij = 1 , grid%num_tiles
576 CALL wrf_debug ( 200 , ' call rk_phys_bc_dry_1' )
578 CALL rk_phys_bc_dry_1( config_flags, grid%ru, grid%rv, grid%rw, grid%ww, &
579 grid%muu, grid%muv, grid%mut, grid%php, grid%alt, grid%p, &
580 ids, ide, jds, jde, kds, kde, &
581 ims, ime, jms, jme, kms, kme, &
582 ips, ipe, jps, jpe, kps, kpe, &
583 grid%i_start(ij), grid%i_end(ij), &
584 grid%j_start(ij), grid%j_end(ij), &
586 CALL set_physical_bc3d( grid%al, 'p', config_flags, &
587 ids, ide, jds, jde, kds, kde, &
588 ims, ime, jms, jme, kms, kme, &
589 ips, ipe, jps, jpe, kps, kpe, &
590 grid%i_start(ij), grid%i_end(ij), &
591 grid%j_start(ij), grid%j_end(ij), &
593 CALL set_physical_bc3d( grid%ph_2, 'w', config_flags, &
594 ids, ide, jds, jde, kds, kde, &
595 ims, ime, jms, jme, kms, kme, &
596 ips, ipe, jps, jpe, kps, kpe, &
597 grid%i_start(ij), grid%i_end(ij), &
598 grid%j_start(ij), grid%j_end(ij), &
601 IF (config_flags%polar) THEN
603 !-------------------------------------------------------
604 ! lat-lon grid pole-point (v) specification (extrapolate v, rv to the pole)
605 !-------------------------------------------------------
607 CALL pole_point_bc ( grid%v_1, &
608 ids, ide, jds, jde, kds, kde, &
609 ims, ime, jms, jme, kms, kme, &
610 grid%i_start(ij), grid%i_end(ij), &
611 grid%j_start(ij), grid%j_end(ij), &
614 CALL pole_point_bc ( grid%v_2, &
615 ids, ide, jds, jde, kds, kde, &
616 ims, ime, jms, jme, kms, kme, &
617 grid%i_start(ij), grid%i_end(ij), &
618 grid%j_start(ij), grid%j_end(ij), &
621 !-------------------------------------------------------
622 ! end lat-lon grid pole-point (v) specification
623 !-------------------------------------------------------
627 !$OMP END PARALLEL DO
628 BENCH_END(set_phys_bc_tim)
630 rk_step_is_one : IF (rk_step == 1) THEN ! only need to initialize diffusion tendencies
634 !(2) The non-timesplit physics begins with a call to "phy_prep"
635 ! (which computes some diagnostic variables such as temperature,
636 ! pressure, u and v at p points, etc). This is followed by
637 ! calls to the physics drivers:
648 CALL first_rk_step_part1 ( grid, config_flags &
649 , moist , moist_tend &
651 , tracer, tracer_tend &
652 , scalar , scalar_tend &
654 , ru_tendf, rv_tendf &
655 , rw_tendf, t_tendf &
656 , ph_tendf, mu_tendf &
658 , config_flags%use_adaptive_time_step &
660 , psim , psih , wspd , gz1oz0 &
662 , cu_act_flag , hol , th_phy &
663 , pi_phy , p_phy , grid%t_phy &
665 , dz8w , p8w , t8w , rho_phy , rho &
666 , ids, ide, jds, jde, kds, kde &
667 , ims, ime, jms, jme, kms, kme &
668 , ips, ipe, jps, jpe, kps, kpe &
669 , imsx, imex, jmsx, jmex, kmsx, kmex &
670 , ipsx, ipex, jpsx, jpex, kpsx, kpex &
671 , imsy, imey, jmsy, jmey, kmsy, kmey &
672 , ipsy, ipey, jpsy, jpey, kpsy, kpey &
678 IF ( config_flags%bl_pbl_physics == MYNNPBLSCHEME2 .OR. &
679 config_flags%bl_pbl_physics == MYNNPBLSCHEME3 ) THEN
680 # include "HALO_EM_SCALAR_E_5.inc"
684 CALL first_rk_step_part2 ( grid, config_flags &
685 , moist , moist_tend &
687 , tracer, tracer_tend &
688 , scalar , scalar_tend &
690 , ru_tendf, rv_tendf &
691 , rw_tendf, t_tendf &
692 , ph_tendf, mu_tendf &
694 , adapt_step_flag , curr_secs &
695 , psim , psih , wspd , gz1oz0 &
697 , cu_act_flag , hol , th_phy &
698 , pi_phy , p_phy , grid%t_phy &
700 , dz8w , p8w , t8w , rho_phy , rho &
701 , nba_mij, num_nba_mij & !JDM
702 , nba_rij, num_nba_rij & !JDM
703 , ids, ide, jds, jde, kds, kde &
704 , ims, ime, jms, jme, kms, kme &
705 , ips, ipe, jps, jpe, kps, kpe &
706 , imsx, imex, jmsx, jmex, kmsx, kmex &
707 , ipsx, ipex, jpsx, jpex, kpsx, kpex &
708 , imsy, imey, jmsy, jmey, kmsy, kmey &
709 , ipsy, ipey, jpsy, jpey, kpsy, kpey &
713 END IF rk_step_is_one
715 BENCH_START(rk_tend_tim)
718 DO ij = 1 , grid%num_tiles
720 CALL wrf_debug ( 200 , ' call rk_tendency' )
721 CALL rk_tendency ( config_flags, rk_step &
722 ,grid%ru_tend, grid%rv_tend, rw_tend, ph_tend, t_tend &
723 ,ru_tendf, rv_tendf, rw_tendf, ph_tendf, t_tendf &
724 ,mu_tend, grid%u_save, grid%v_save, w_save, ph_save &
725 ,grid%t_save, mu_save, grid%rthften &
726 ,grid%ru, grid%rv, grid%rw, grid%ww &
727 ,grid%u_2, grid%v_2, grid%w_2, grid%t_2, grid%ph_2 &
728 ,grid%u_1, grid%v_1, grid%w_1, grid%t_1, grid%ph_1 &
729 ,grid%h_diabatic, grid%phb, grid%t_init &
730 ,grid%mu_2, grid%mut, grid%muu, grid%muv, grid%mub &
731 ,grid%al, grid%alt, grid%p, grid%pb, grid%php, cqu, cqv, cqw &
732 ,grid%u_base, grid%v_base, grid%t_base, grid%qv_base, grid%z_base &
733 ,grid%msfux,grid%msfuy, grid%msfvx, grid%msfvx_inv &
734 ,grid%msfvy, grid%msftx,grid%msfty, grid%clat, grid%f, grid%e, grid%sina, grid%cosa &
735 ,grid%fnm, grid%fnp, grid%rdn, grid%rdnw &
736 ,grid%dt, grid%rdx, grid%rdy, grid%khdif, grid%kvdif, grid%xkmh, grid%xkhh &
737 ,grid%diff_6th_opt, grid%diff_6th_factor &
738 ,grid%dampcoef,grid%zdamp,config_flags%damp_opt,config_flags%rad_nudge &
739 ,grid%cf1, grid%cf2, grid%cf3, grid%cfn, grid%cfn1, num_3d_m &
740 ,config_flags%non_hydrostatic, config_flags%top_lid &
741 ,grid%u_frame, grid%v_frame &
742 ,ids, ide, jds, jde, kds, kde &
743 ,ims, ime, jms, jme, kms, kme &
744 ,grid%i_start(ij), grid%i_end(ij) &
745 ,grid%j_start(ij), grid%j_end(ij) &
747 ,max_vert_cfl_tmp(ij), max_horiz_cfl_tmp(ij) )
749 !$OMP END PARALLEL DO
750 BENCH_END(rk_tend_tim)
752 IF (config_flags%use_adaptive_time_step) THEN
753 DO ij = 1 , grid%num_tiles
754 IF (max_horiz_cfl_tmp(ij) .GT. grid%max_horiz_cfl) THEN
755 grid%max_horiz_cfl = max_horiz_cfl_tmp(ij)
757 IF (max_vert_cfl_tmp(ij) .GT. grid%max_vert_cfl) THEN
758 grid%max_vert_cfl = max_vert_cfl_tmp(ij)
762 IF (grid%max_horiz_cfl .GT. grid%max_cfl_val) THEN
763 grid%max_cfl_val = grid%max_horiz_cfl
765 IF (grid%max_vert_cfl .GT. grid%max_cfl_val) THEN
766 grid%max_cfl_val = grid%max_vert_cfl
770 BENCH_START(relax_bdy_dry_tim)
773 DO ij = 1 , grid%num_tiles
775 IF( (config_flags%specified .or. config_flags%nested) .and. rk_step == 1 ) THEN
777 CALL relax_bdy_dry ( config_flags, &
778 grid%u_save, grid%v_save, ph_save, grid%t_save, &
780 grid%ru, grid%rv, grid%ph_2, grid%t_2, &
781 grid%w_2, grid%mu_2, grid%mut, &
782 grid%u_bxs,grid%u_bxe,grid%u_bys,grid%u_bye, &
783 grid%v_bxs,grid%v_bxe,grid%v_bys,grid%v_bye, &
784 grid%ph_bxs,grid%ph_bxe,grid%ph_bys,grid%ph_bye, &
785 grid%t_bxs,grid%t_bxe,grid%t_bys,grid%t_bye, &
786 grid%w_bxs,grid%w_bxe,grid%w_bys,grid%w_bye, &
787 grid%mu_bxs,grid%mu_bxe,grid%mu_bys,grid%mu_bye, &
788 grid%u_btxs,grid%u_btxe,grid%u_btys,grid%u_btye, &
789 grid%v_btxs,grid%v_btxe,grid%v_btys,grid%v_btye, &
790 grid%ph_btxs,grid%ph_btxe,grid%ph_btys,grid%ph_btye, &
791 grid%t_btxs,grid%t_btxe,grid%t_btys,grid%t_btye, &
792 grid%w_btxs,grid%w_btxe,grid%w_btys,grid%w_btye, &
793 grid%mu_btxs,grid%mu_btxe,grid%mu_btys,grid%mu_btye, &
794 config_flags%spec_bdy_width, grid%spec_zone, grid%relax_zone, &
795 grid%dtbc, grid%fcx, grid%gcx, &
796 ids,ide, jds,jde, kds,kde, &
797 ims,ime, jms,jme, kms,kme, &
798 ips,ipe, jps,jpe, kps,kpe, &
799 grid%i_start(ij), grid%i_end(ij), &
800 grid%j_start(ij), grid%j_end(ij), &
805 CALL rk_addtend_dry( grid%ru_tend, grid%rv_tend, rw_tend, ph_tend, t_tend, &
806 ru_tendf, rv_tendf, rw_tendf, ph_tendf, t_tendf, &
807 grid%u_save, grid%v_save, w_save, ph_save, grid%t_save, &
808 mu_tend, mu_tendf, rk_step, &
809 grid%h_diabatic, grid%mut, grid%msftx, &
810 grid%msfty, grid%msfux,grid%msfuy, &
811 grid%msfvx, grid%msfvx_inv, grid%msfvy, &
812 ids,ide, jds,jde, kds,kde, &
813 ims,ime, jms,jme, kms,kme, &
814 ips,ipe, jps,jpe, kps,kpe, &
815 grid%i_start(ij), grid%i_end(ij), &
816 grid%j_start(ij), grid%j_end(ij), &
819 IF( config_flags%specified .or. config_flags%nested ) THEN
820 CALL spec_bdy_dry ( config_flags, &
821 grid%ru_tend, grid%rv_tend, ph_tend, t_tend, &
823 grid%u_bxs,grid%u_bxe,grid%u_bys,grid%u_bye, &
824 grid%v_bxs,grid%v_bxe,grid%v_bys,grid%v_bye, &
825 grid%ph_bxs,grid%ph_bxe,grid%ph_bys,grid%ph_bye, &
826 grid%t_bxs,grid%t_bxe,grid%t_bys,grid%t_bye, &
827 grid%w_bxs,grid%w_bxe,grid%w_bys,grid%w_bye, &
828 grid%mu_bxs,grid%mu_bxe,grid%mu_bys,grid%mu_bye, &
829 grid%u_btxs,grid%u_btxe,grid%u_btys,grid%u_btye, &
830 grid%v_btxs,grid%v_btxe,grid%v_btys,grid%v_btye, &
831 grid%ph_btxs,grid%ph_btxe,grid%ph_btys,grid%ph_btye, &
832 grid%t_btxs,grid%t_btxe,grid%t_btys,grid%t_btye, &
833 grid%w_btxs,grid%w_btxe,grid%w_btys,grid%w_btye, &
834 grid%mu_btxs,grid%mu_btxe,grid%mu_btys,grid%mu_btye, &
835 config_flags%spec_bdy_width, grid%spec_zone, &
836 ids,ide, jds,jde, kds,kde, & ! domain dims
837 ims,ime, jms,jme, kms,kme, & ! memory dims
838 ips,ipe, jps,jpe, kps,kpe, & ! patch dims
839 grid%i_start(ij), grid%i_end(ij), &
840 grid%j_start(ij), grid%j_end(ij), &
846 !$OMP END PARALLEL DO
847 BENCH_END(relax_bdy_dry_tim)
851 ! (3) Small (acoustic,sound) steps.
853 ! Several acoustic steps are taken each RK pass. A small step
854 ! sequence begins with calculating perturbation variables
855 ! and coupling them to the column dry-air-mass mu
856 ! (call to small_step_prep). This is followed by computing
857 ! coefficients for the vertically implicit part of the
858 ! small timestep (call to calc_coef_w).
860 ! The small steps are taken
861 ! in the named loop "small_steps:". In the small_steps loop, first
862 ! the horizontal momentum (u and v) are advanced (call to advance_uv),
863 ! next mu and theta are advanced (call to advance_mu_t) followed by
864 ! advancing w and the geopotential (call to advance_w). Diagnostic
865 ! values for pressure and inverse density are updated at the end of
868 ! The small-step section ends with the change of the perturbation variables
869 ! back to full variables (call to small_step_finish).
873 BENCH_START(small_step_prep_tim)
876 DO ij = 1 , grid%num_tiles
878 ! Calculate coefficients for the vertically implicit acoustic/gravity wave
879 ! integration. We only need calculate these for the first pass through -
880 ! the predictor step. They are reused as is for the corrector step.
881 ! For third-order RK, we need to recompute these after the first
882 ! predictor because we may have changed the small timestep -> grid%dts.
884 CALL wrf_debug ( 200 , ' call small_step_prep ' )
886 CALL small_step_prep( grid%u_1,grid%u_2,grid%v_1,grid%v_2,grid%w_1,grid%w_2, &
887 grid%t_1,grid%t_2,grid%ph_1,grid%ph_2, &
888 grid%mub, grid%mu_1, grid%mu_2, &
889 grid%muu, muus, grid%muv, muvs, &
890 grid%mut, grid%muts, grid%mudf, &
891 grid%u_save, grid%v_save, w_save, &
892 grid%t_save, ph_save, mu_save, &
894 grid%dnw, c2a, grid%pb, grid%p, grid%alt, &
895 grid%msfux,grid%msfuy, grid%msfvx, grid%msfvx_inv, &
896 grid%msfvy, grid%msftx,grid%msfty, &
897 grid%rdx, grid%rdy, rk_step, &
898 ids, ide, jds, jde, kds, kde, &
899 ims, ime, jms, jme, kms, kme, &
900 grid%i_start(ij), grid%i_end(ij), &
901 grid%j_start(ij), grid%j_end(ij), &
904 CALL calc_p_rho( grid%al, grid%p, grid%ph_2, &
905 grid%alt, grid%t_2, grid%t_save, c2a, pm1, &
906 grid%mu_2, grid%muts, grid%znu, t0, &
907 grid%rdnw, grid%dnw, grid%smdiv, &
908 config_flags%non_hydrostatic, 0, &
909 ids, ide, jds, jde, kds, kde, &
910 ims, ime, jms, jme, kms, kme, &
911 grid%i_start(ij), grid%i_end(ij), &
912 grid%j_start(ij), grid%j_end(ij), &
915 IF (config_flags%non_hydrostatic) THEN
916 CALL calc_coef_w( a,alpha,gamma, &
918 grid%rdn, grid%rdnw, c2a, &
919 dts_rk, g, grid%epssm, &
920 config_flags%top_lid, &
921 ids, ide, jds, jde, kds, kde, &
922 ims, ime, jms, jme, kms, kme, &
923 grid%i_start(ij), grid%i_end(ij), &
924 grid%j_start(ij), grid%j_end(ij), &
929 !$OMP END PARALLEL DO
930 BENCH_END(small_step_prep_tim)
933 !-----------------------------------------------------------------------
934 ! Stencils for patch communications (WCS, 29 June 2001)
935 ! Note: the small size of this halo exchange reflects the
936 ! fact that we are carrying the uncoupled variables
937 ! as state variables in the mass coordinate model, as
938 ! opposed to the coupled variables as in the height
943 ! * + * * + * * * + * *
947 ! 3D variables - note staggering! ph_2(Z), u_save(X), v_save(Y)
957 ! the following are 2D (xy) variables
965 !--------------------------------------------------------------
966 # include "HALO_EM_B.inc"
967 # include "PERIOD_BDY_EM_B.inc"
970 BENCH_START(set_phys_bc2_tim)
974 DO ij = 1 , grid%num_tiles
976 CALL set_physical_bc3d( grid%ru_tend, 'u', config_flags, &
977 ids, ide, jds, jde, kds, kde, &
978 ims, ime, jms, jme, kms, kme, &
979 ips, ipe, jps, jpe, kps, kpe, &
980 grid%i_start(ij), grid%i_end(ij), &
981 grid%j_start(ij), grid%j_end(ij), &
984 CALL set_physical_bc3d( grid%rv_tend, 'v', config_flags, &
985 ids, ide, jds, jde, kds, kde, &
986 ims, ime, jms, jme, kms, kme, &
987 ips, ipe, jps, jpe, kps, kpe, &
988 grid%i_start(ij), grid%i_end(ij), &
989 grid%j_start(ij), grid%j_end(ij), &
992 CALL set_physical_bc3d( grid%ph_2, 'w', config_flags, &
993 ids, ide, jds, jde, kds, kde, &
994 ims, ime, jms, jme, kms, kme, &
995 ips, ipe, jps, jpe, kps, kpe, &
996 grid%i_start(ij), grid%i_end(ij), &
997 grid%j_start(ij), grid%j_end(ij), &
1000 CALL set_physical_bc3d( grid%al, 'p', config_flags, &
1001 ids, ide, jds, jde, kds, kde, &
1002 ims, ime, jms, jme, kms, kme, &
1003 ips, ipe, jps, jpe, kps, kpe, &
1004 grid%i_start(ij), grid%i_end(ij), &
1005 grid%j_start(ij), grid%j_end(ij), &
1008 CALL set_physical_bc3d( grid%p, 'p', config_flags, &
1009 ids, ide, jds, jde, kds, kde, &
1010 ims, ime, jms, jme, kms, kme, &
1011 ips, ipe, jps, jpe, kps, kpe, &
1012 grid%i_start(ij), grid%i_end(ij), &
1013 grid%j_start(ij), grid%j_end(ij), &
1016 CALL set_physical_bc3d( grid%t_1, 'p', config_flags, &
1017 ids, ide, jds, jde, kds, kde, &
1018 ims, ime, jms, jme, kms, kme, &
1019 ips, ipe, jps, jpe, kps, kpe, &
1020 grid%i_start(ij), grid%i_end(ij), &
1021 grid%j_start(ij), grid%j_end(ij), &
1024 CALL set_physical_bc3d( grid%t_save, 't', config_flags, &
1025 ids, ide, jds, jde, kds, kde, &
1026 ims, ime, jms, jme, kms, kme, &
1027 ips, ipe, jps, jpe, kps, kpe, &
1028 grid%i_start(ij), grid%i_end(ij), &
1029 grid%j_start(ij), grid%j_end(ij), &
1032 CALL set_physical_bc2d( grid%mu_1, 't', config_flags, &
1033 ids, ide, jds, jde, &
1034 ims, ime, jms, jme, &
1035 ips, ipe, jps, jpe, &
1036 grid%i_start(ij), grid%i_end(ij), &
1037 grid%j_start(ij), grid%j_end(ij) )
1039 CALL set_physical_bc2d( grid%mu_2, 't', config_flags, &
1040 ids, ide, jds, jde, &
1041 ims, ime, jms, jme, &
1042 ips, ipe, jps, jpe, &
1043 grid%i_start(ij), grid%i_end(ij), &
1044 grid%j_start(ij), grid%j_end(ij) )
1046 CALL set_physical_bc2d( grid%mudf, 't', config_flags, &
1047 ids, ide, jds, jde, &
1048 ims, ime, jms, jme, &
1049 ips, ipe, jps, jpe, &
1050 grid%i_start(ij), grid%i_end(ij), &
1051 grid%j_start(ij), grid%j_end(ij) )
1054 !$OMP END PARALLEL DO
1055 BENCH_END(set_phys_bc2_tim)
1056 small_steps : DO iteration = 1 , number_of_small_timesteps
1058 ! Boundary condition time (or communication time).
1060 # include "PERIOD_BDY_EM_B.inc"
1064 !$OMP PRIVATE ( ij )
1066 DO ij = 1 , grid%num_tiles
1068 BENCH_START(advance_uv_tim)
1069 CALL advance_uv ( grid%u_2, grid%ru_tend, grid%v_2, grid%rv_tend, &
1071 grid%ph_2, grid%php, grid%alt, grid%al, &
1073 grid%muu, cqu, grid%muv, cqv, grid%mudf, &
1074 grid%msfux, grid%msfuy, grid%msfvx, &
1075 grid%msfvx_inv, grid%msfvy, &
1076 grid%rdx, grid%rdy, dts_rk, &
1077 grid%cf1, grid%cf2, grid%cf3, grid%fnm, grid%fnp, &
1079 grid%rdnw, config_flags,grid%spec_zone, &
1080 config_flags%non_hydrostatic, config_flags%top_lid, &
1081 ids, ide, jds, jde, kds, kde, &
1082 ims, ime, jms, jme, kms, kme, &
1083 grid%i_start(ij), grid%i_end(ij), &
1084 grid%j_start(ij), grid%j_end(ij), &
1086 BENCH_END(advance_uv_tim)
1089 !$OMP END PARALLEL DO
1091 !-----------------------------------------------------------
1092 ! acoustic integration polar filter for smallstep u, v
1093 !-----------------------------------------------------------
1095 IF (config_flags%polar) THEN
1097 CALL pxft ( grid=grid &
1110 ,positive_definite = .FALSE. &
1111 ,moist=moist,chem=chem,tracer=tracer,scalar=scalar &
1112 ,fft_filter_lat = config_flags%fft_filter_lat &
1114 ,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
1115 ,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
1116 ,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe &
1117 ,imsx=imsx,imex=imex,jmsx=jmsx,jmex=jmex,kmsx=kmsx,kmex=kmex &
1118 ,ipsx=ipsx,ipex=ipex,jpsx=jmsx,jpex=jpex,kpsx=kpsx,kpex=kpex )
1122 !-----------------------------------------------------------
1123 ! end acoustic integration polar filter for smallstep u, v
1124 !-----------------------------------------------------------
1127 !$OMP PRIVATE ( ij )
1128 DO ij = 1 , grid%num_tiles
1130 BENCH_START(spec_bdy_uv_tim)
1131 IF( config_flags%specified .or. config_flags%nested ) THEN
1132 CALL spec_bdyupdate(grid%u_2, grid%ru_tend, dts_rk, &
1133 'u' , config_flags, &
1135 ids,ide, jds,jde, kds,kde, & ! domain dims
1136 ims,ime, jms,jme, kms,kme, & ! memory dims
1137 ips,ipe, jps,jpe, kps,kpe, & ! patch dims
1138 grid%i_start(ij), grid%i_end(ij), &
1139 grid%j_start(ij), grid%j_end(ij), &
1142 CALL spec_bdyupdate(grid%v_2, grid%rv_tend, dts_rk, &
1143 'v' , config_flags, &
1145 ids,ide, jds,jde, kds,kde, & ! domain dims
1146 ims,ime, jms,jme, kms,kme, & ! memory dims
1147 ips,ipe, jps,jpe, kps,kpe, & ! patch dims
1148 grid%i_start(ij), grid%i_end(ij), &
1149 grid%j_start(ij), grid%j_end(ij), &
1153 BENCH_END(spec_bdy_uv_tim)
1156 !$OMP END PARALLEL DO
1160 ! Stencils for patch communications (WCS, 29 June 2001)
1169 # include "HALO_EM_C.inc"
1173 !$OMP PRIVATE ( ij )
1174 DO ij = 1 , grid%num_tiles
1176 ! advance the mass in the column, theta, and calculate ww
1178 BENCH_START(advance_mu_t_tim)
1179 CALL advance_mu_t( grid%ww, ww1, grid%u_2, grid%u_save, grid%v_2, grid%v_save, &
1180 grid%mu_2, grid%mut, muave, grid%muts, grid%muu, grid%muv, &
1181 grid%mudf, grid%ru_m, grid%rv_m, grid%ww_m, &
1182 grid%t_2, grid%t_save, t_2save, t_tend, &
1184 grid%rdx, grid%rdy, dts_rk, grid%epssm, &
1185 grid%dnw, grid%fnm, grid%fnp, grid%rdnw, &
1186 grid%msfux,grid%msfuy, grid%msfvx, grid%msfvx_inv, &
1187 grid%msfvy, grid%msftx,grid%msfty, &
1188 iteration, config_flags, &
1189 ids, ide, jds, jde, kds, kde, &
1190 ims, ime, jms, jme, kms, kme, &
1191 grid%i_start(ij), grid%i_end(ij), &
1192 grid%j_start(ij), grid%j_end(ij), &
1194 BENCH_END(advance_mu_t_tim)
1196 !$OMP END PARALLEL DO
1198 !-----------------------------------------------------------
1199 ! acoustic integration polar filter for smallstep mu, t
1200 !-----------------------------------------------------------
1202 IF ( (config_flags%polar) ) THEN
1204 CALL pxft ( grid=grid &
1217 ,positive_definite = .FALSE. &
1218 ,moist=moist,chem=chem,tracer=tracer,scalar=scalar &
1219 ,fft_filter_lat = config_flags%fft_filter_lat &
1221 ,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
1222 ,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
1223 ,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe &
1224 ,imsx=imsx,imex=imex,jmsx=jmsx,jmex=jmex,kmsx=kmsx,kmex=kmex &
1225 ,ipsx=ipsx,ipex=ipex,jpsx=jmsx,jpex=jpex,kpsx=kpsx,kpex=kpex )
1227 grid%muts = grid%mut + grid%mu_2 ! reset muts using filtered mu_2
1231 !-----------------------------------------------------------
1232 ! end acoustic integration polar filter for smallstep mu, t
1233 !-----------------------------------------------------------
1235 BENCH_START(spec_bdy_t_tim)
1238 !$OMP PRIVATE ( ij )
1239 DO ij = 1 , grid%num_tiles
1241 IF( config_flags%specified .or. config_flags%nested ) THEN
1243 CALL spec_bdyupdate(grid%t_2, t_tend, dts_rk, &
1244 't' , config_flags, &
1246 ids,ide, jds,jde, kds,kde, &
1247 ims,ime, jms,jme, kms,kme, &
1248 ips,ipe, jps,jpe, kps,kpe, &
1249 grid%i_start(ij), grid%i_end(ij),&
1250 grid%j_start(ij), grid%j_end(ij),&
1253 CALL spec_bdyupdate(grid%mu_2, mu_tend, dts_rk, &
1254 'm' , config_flags, &
1256 ids,ide, jds,jde, 1 ,1 , &
1257 ims,ime, jms,jme, 1 ,1 , &
1258 ips,ipe, jps,jpe, 1 ,1 , &
1259 grid%i_start(ij), grid%i_end(ij),&
1260 grid%j_start(ij), grid%j_end(ij),&
1263 CALL spec_bdyupdate(grid%muts, mu_tend, dts_rk, &
1264 'm' , config_flags, &
1266 ids,ide, jds,jde, 1 ,1 , & ! domain dims
1267 ims,ime, jms,jme, 1 ,1 , & ! memory dims
1268 ips,ipe, jps,jpe, 1 ,1 , & ! patch dims
1269 grid%i_start(ij), grid%i_end(ij), &
1270 grid%j_start(ij), grid%j_end(ij), &
1273 BENCH_END(spec_bdy_t_tim)
1275 ! small (acoustic) step for the vertical momentum,
1276 ! density and coupled potential temperature.
1279 BENCH_START(advance_w_tim)
1280 IF ( config_flags%non_hydrostatic ) THEN
1281 CALL advance_w( grid%w_2, rw_tend, grid%ww, w_save, &
1282 grid%u_2, grid%v_2, &
1283 grid%mu_2, grid%mut, muave, grid%muts, &
1284 t_2save, grid%t_2, grid%t_save, &
1285 grid%ph_2, ph_save, grid%phb, ph_tend, &
1286 grid%ht, c2a, cqw, grid%alt, grid%alb, &
1288 grid%rdx, grid%rdy, dts_rk, t0, grid%epssm, &
1289 grid%dnw, grid%fnm, grid%fnp, grid%rdnw, &
1290 grid%rdn, grid%cf1, grid%cf2, grid%cf3, &
1291 grid%msftx, grid%msfty, &
1292 config_flags, config_flags%top_lid, &
1293 ids,ide, jds,jde, kds,kde, &
1294 ims,ime, jms,jme, kms,kme, &
1295 grid%i_start(ij), grid%i_end(ij), &
1296 grid%j_start(ij), grid%j_end(ij), &
1299 BENCH_END(advance_w_tim)
1302 !$OMP END PARALLEL DO
1304 !-----------------------------------------------------------
1305 ! acoustic integration polar filter for smallstep w, geopotential
1306 !-----------------------------------------------------------
1308 IF ( (config_flags%polar) .AND. (config_flags%non_hydrostatic) ) THEN
1310 CALL pxft ( grid=grid &
1323 ,positive_definite = .FALSE. &
1324 ,moist=moist,chem=chem,tracer=tracer,scalar=scalar &
1325 ,fft_filter_lat = config_flags%fft_filter_lat &
1327 ,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
1328 ,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
1329 ,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe &
1330 ,imsx=imsx,imex=imex,jmsx=jmsx,jmex=jmex,kmsx=kmsx,kmex=kmex &
1331 ,ipsx=ipsx,ipex=ipex,jpsx=jmsx,jpex=jpex,kpsx=kpsx,kpex=kpex )
1335 !-----------------------------------------------------------
1336 ! end acoustic integration polar filter for smallstep w, geopotential
1337 !-----------------------------------------------------------
1340 !$OMP PRIVATE ( ij )
1341 DO ij = 1 , grid%num_tiles
1343 BENCH_START(sumflux_tim)
1344 CALL sumflux ( grid%u_2, grid%v_2, grid%ww, &
1345 grid%u_save, grid%v_save, ww1, &
1346 grid%muu, grid%muv, &
1347 grid%ru_m, grid%rv_m, grid%ww_m, grid%epssm, &
1348 grid%msfux, grid% msfuy, grid%msfvx, &
1349 grid%msfvx_inv, grid%msfvy, &
1350 iteration, number_of_small_timesteps, &
1351 ids, ide, jds, jde, kds, kde, &
1352 ims, ime, jms, jme, kms, kme, &
1353 grid%i_start(ij), grid%i_end(ij), &
1354 grid%j_start(ij), grid%j_end(ij), &
1356 BENCH_END(sumflux_tim)
1358 IF( config_flags%specified .or. config_flags%nested ) THEN
1360 BENCH_START(spec_bdynhyd_tim)
1361 IF (config_flags%non_hydrostatic) THEN
1362 CALL spec_bdyupdate_ph( ph_save, grid%ph_2, ph_tend, &
1363 mu_tend, grid%muts, dts_rk, &
1364 'h' , config_flags, &
1366 ids,ide, jds,jde, kds,kde, &
1367 ims,ime, jms,jme, kms,kme, &
1368 ips,ipe, jps,jpe, kps,kpe, &
1369 grid%i_start(ij), grid%i_end(ij),&
1370 grid%j_start(ij), grid%j_end(ij),&
1372 IF( config_flags%specified ) THEN
1373 CALL zero_grad_bdy ( grid%w_2, &
1374 'w' , config_flags, &
1376 ids,ide, jds,jde, kds,kde, &
1377 ims,ime, jms,jme, kms,kme, &
1378 ips,ipe, jps,jpe, kps,kpe, &
1379 grid%i_start(ij), grid%i_end(ij), &
1380 grid%j_start(ij), grid%j_end(ij), &
1383 CALL spec_bdyupdate ( grid%w_2, rw_tend, dts_rk, &
1384 'h' , config_flags, &
1386 ids,ide, jds,jde, kds,kde, &
1387 ims,ime, jms,jme, kms,kme, &
1388 ips,ipe, jps,jpe, kps,kpe, &
1389 grid%i_start(ij), grid%i_end(ij),&
1390 grid%j_start(ij), grid%j_end(ij),&
1394 BENCH_END(spec_bdynhyd_tim)
1397 BENCH_START(cald_p_rho_tim)
1398 CALL calc_p_rho( grid%al, grid%p, grid%ph_2, &
1399 grid%alt, grid%t_2, grid%t_save, c2a, pm1, &
1400 grid%mu_2, grid%muts, grid%znu, t0, &
1401 grid%rdnw, grid%dnw, grid%smdiv, &
1402 config_flags%non_hydrostatic, iteration, &
1403 ids, ide, jds, jde, kds, kde, &
1404 ims, ime, jms, jme, kms, kme, &
1405 grid%i_start(ij), grid%i_end(ij), &
1406 grid%j_start(ij), grid%j_end(ij), &
1408 BENCH_END(cald_p_rho_tim)
1411 !$OMP END PARALLEL DO
1415 ! Stencils for patch communications (WCS, 29 June 2001)
1425 ! 2D variables (x,y)
1431 # include "HALO_EM_C2.inc"
1432 # include "PERIOD_BDY_EM_B3.inc"
1435 BENCH_START(phys_bc_tim)
1437 !$OMP PRIVATE ( ij )
1438 DO ij = 1 , grid%num_tiles
1440 ! boundary condition set for next small timestep
1442 CALL set_physical_bc3d( grid%ph_2, 'w', config_flags, &
1443 ids, ide, jds, jde, kds, kde, &
1444 ims, ime, jms, jme, kms, kme, &
1445 ips, ipe, jps, jpe, kps, kpe, &
1446 grid%i_start(ij), grid%i_end(ij), &
1447 grid%j_start(ij), grid%j_end(ij), &
1450 CALL set_physical_bc3d( grid%al, 'p', config_flags, &
1451 ids, ide, jds, jde, kds, kde, &
1452 ims, ime, jms, jme, kms, kme, &
1453 ips, ipe, jps, jpe, kps, kpe, &
1454 grid%i_start(ij), grid%i_end(ij), &
1455 grid%j_start(ij), grid%j_end(ij), &
1458 CALL set_physical_bc3d( grid%p, 'p', config_flags, &
1459 ids, ide, jds, jde, kds, kde, &
1460 ims, ime, jms, jme, kms, kme, &
1461 ips, ipe, jps, jpe, kps, kpe, &
1462 grid%i_start(ij), grid%i_end(ij), &
1463 grid%j_start(ij), grid%j_end(ij), &
1466 CALL set_physical_bc2d( grid%muts, 't', config_flags, &
1467 ids, ide, jds, jde, &
1468 ims, ime, jms, jme, &
1469 ips, ipe, jps, jpe, &
1470 grid%i_start(ij), grid%i_end(ij), &
1471 grid%j_start(ij), grid%j_end(ij) )
1473 CALL set_physical_bc2d( grid%mu_2, 't', config_flags, &
1474 ids, ide, jds, jde, &
1475 ims, ime, jms, jme, &
1476 ips, ipe, jps, jpe, &
1477 grid%i_start(ij), grid%i_end(ij), &
1478 grid%j_start(ij), grid%j_end(ij) )
1480 CALL set_physical_bc2d( grid%mudf, 't', config_flags, &
1481 ids, ide, jds, jde, &
1482 ims, ime, jms, jme, &
1483 ips, ipe, jps, jpe, &
1484 grid%i_start(ij), grid%i_end(ij), &
1485 grid%j_start(ij), grid%j_end(ij) )
1488 !$OMP END PARALLEL DO
1489 BENCH_END(phys_bc_tim)
1494 !$OMP PRIVATE ( ij )
1495 DO ij = 1 , grid%num_tiles
1497 CALL wrf_debug ( 200 , ' call rk_small_finish' )
1499 ! change time-perturbation variables back to
1500 ! full perturbation variables.
1501 ! first get updated mu at u and v points
1503 BENCH_START(calc_mu_uv_tim)
1504 CALL calc_mu_uv_1 ( config_flags, &
1505 grid%muts, muus, muvs, &
1506 ids, ide, jds, jde, kds, kde, &
1507 ims, ime, jms, jme, kms, kme, &
1508 grid%i_start(ij), grid%i_end(ij), &
1509 grid%j_start(ij), grid%j_end(ij), &
1511 BENCH_END(calc_mu_uv_tim)
1512 BENCH_START(small_step_finish_tim)
1513 CALL small_step_finish( grid%u_2, grid%u_1, grid%v_2, grid%v_1, grid%w_2, grid%w_1, &
1514 grid%t_2, grid%t_1, grid%ph_2, grid%ph_1, grid%ww, ww1, &
1515 grid%mu_2, grid%mu_1, &
1516 grid%mut, grid%muts, grid%muu, muus, grid%muv, muvs, &
1517 grid%u_save, grid%v_save, w_save, &
1518 grid%t_save, ph_save, mu_save, &
1519 grid%msfux,grid%msfuy, grid%msfvx,grid%msfvy, grid%msftx,grid%msfty, &
1521 number_of_small_timesteps,dts_rk, &
1522 rk_step, rk_order, &
1523 ids, ide, jds, jde, kds, kde, &
1524 ims, ime, jms, jme, kms, kme, &
1525 grid%i_start(ij), grid%i_end(ij), &
1526 grid%j_start(ij), grid%j_end(ij), &
1528 ! call to set ru_m, rv_m and ww_m b.c's for PD advection
1530 IF (rk_step == rk_order) THEN
1532 CALL set_physical_bc3d( grid%ru_m, 'u', config_flags, &
1533 ids, ide, jds, jde, kds, kde, &
1534 ims, ime, jms, jme, kms, kme, &
1535 ips, ipe, jps, jpe, kps, kpe, &
1536 grid%i_start(ij), grid%i_end(ij), &
1537 grid%j_start(ij), grid%j_end(ij), &
1540 CALL set_physical_bc3d( grid%rv_m, 'v', config_flags, &
1541 ids, ide, jds, jde, kds, kde, &
1542 ims, ime, jms, jme, kms, kme, &
1543 ips, ipe, jps, jpe, kps, kpe, &
1544 grid%i_start(ij), grid%i_end(ij), &
1545 grid%j_start(ij), grid%j_end(ij), &
1548 CALL set_physical_bc3d( grid%ww_m, 'w', config_flags, &
1549 ids, ide, jds, jde, kds, kde, &
1550 ims, ime, jms, jme, kms, kme, &
1551 ips, ipe, jps, jpe, kps, kpe, &
1552 grid%i_start(ij), grid%i_end(ij), &
1553 grid%j_start(ij), grid%j_end(ij), &
1556 CALL set_physical_bc2d( grid%mut, 't', config_flags, &
1557 ids, ide, jds, jde, &
1558 ims, ime, jms, jme, &
1559 ips, ipe, jps, jpe, &
1560 grid%i_start(ij), grid%i_end(ij), &
1561 grid%j_start(ij), grid%j_end(ij) )
1563 CALL set_physical_bc2d( grid%muts, 't', config_flags, &
1564 ids, ide, jds, jde, &
1565 ims, ime, jms, jme, &
1566 ips, ipe, jps, jpe, &
1567 grid%i_start(ij), grid%i_end(ij), &
1568 grid%j_start(ij), grid%j_end(ij) )
1572 BENCH_END(small_step_finish_tim)
1575 !$OMP END PARALLEL DO
1577 !-----------------------------------------------------------
1578 ! polar filter for full dynamics variables and time-averaged mass fluxes
1579 !-----------------------------------------------------------
1581 IF (config_flags%polar) THEN
1583 CALL pxft ( grid=grid &
1596 ,positive_definite = .FALSE. &
1597 ,moist=moist,chem=chem,tracer=tracer,scalar=scalar &
1598 ,fft_filter_lat = config_flags%fft_filter_lat &
1600 ,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
1601 ,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
1602 ,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe &
1603 ,imsx=imsx,imex=imex,jmsx=jmsx,jmex=jmex,kmsx=kmsx,kmex=kmex &
1604 ,ipsx=ipsx,ipex=ipex,jpsx=jmsx,jpex=jpex,kpsx=kpsx,kpex=kpex )
1608 !-----------------------------------------------------------
1609 ! end polar filter for full dynamics variables and time-averaged mass fluxes
1610 !-----------------------------------------------------------
1612 !-----------------------------------------------------------------------
1613 ! add in physics tendency first if positive definite advection is used.
1614 ! pd advection applies advective flux limiter on last runge-kutta step
1615 !-----------------------------------------------------------------------
1618 IF ((config_flags%moist_adv_opt /= ORIGINAL) .and. (rk_step == rk_order)) THEN
1621 !$OMP PRIVATE ( ij )
1622 DO ij = 1 , grid%num_tiles
1623 CALL wrf_debug ( 200 , ' call rk_update_scalar_pd' )
1624 DO im = PARAM_FIRST_SCALAR, num_3d_m
1625 CALL rk_update_scalar_pd( im, im, &
1626 moist_old(ims,kms,jms,im), &
1627 moist_tend(ims,kms,jms,im), &
1628 grid%mu_1, grid%mu_1, grid%mub, &
1629 rk_step, dt_rk, grid%spec_zone, &
1631 ids, ide, jds, jde, kds, kde, &
1632 ims, ime, jms, jme, kms, kme, &
1633 grid%i_start(ij), grid%i_end(ij), &
1634 grid%j_start(ij), grid%j_end(ij), &
1638 !$OMP END PARALLEL DO
1640 !---------------------- positive definite bc call
1642 IF (config_flags%moist_adv_opt /= ORIGINAL) THEN
1643 IF ( config_flags%h_sca_adv_order <= 4 ) THEN
1644 # include "HALO_EM_MOIST_OLD_E_5.inc"
1645 ELSE IF ( config_flags%h_sca_adv_order <= 6 ) THEN
1646 # include "HALO_EM_MOIST_OLD_E_7.inc"
1648 WRITE(wrf_err_message,*)'solve_em: invalid h_sca_adv_order = ',config_flags%h_sca_adv_order
1649 CALL wrf_error_fatal(TRIM(wrf_err_message))
1655 # include "PERIOD_BDY_EM_MOIST_OLD.inc"
1659 !$OMP PRIVATE ( ij )
1660 DO ij = 1 , grid%num_tiles
1661 IF (num_3d_m >= PARAM_FIRST_SCALAR) THEN
1662 DO im = PARAM_FIRST_SCALAR , num_3d_m
1663 CALL set_physical_bc3d( moist_old(ims,kms,jms,im), 'p', config_flags, &
1664 ids, ide, jds, jde, kds, kde, &
1665 ims, ime, jms, jme, kms, kme, &
1666 ips, ipe, jps, jpe, kps, kpe, &
1667 grid%i_start(ij), grid%i_end(ij), &
1668 grid%j_start(ij), grid%j_end(ij), &
1673 !$OMP END PARALLEL DO
1675 END IF ! end if for moist_adv_opt
1679 IF ((config_flags%scalar_adv_opt /= ORIGINAL) .and. (rk_step == rk_order)) THEN
1682 !$OMP PRIVATE ( ij )
1683 DO ij = 1 , grid%num_tiles
1684 CALL wrf_debug ( 200 , ' call rk_update_scalar_pd' )
1685 DO im = PARAM_FIRST_SCALAR, num_3d_s
1686 CALL rk_update_scalar_pd( im, im, &
1687 scalar_old(ims,kms,jms,im), &
1688 scalar_tend(ims,kms,jms,im), &
1689 grid%mu_1, grid%mu_1, grid%mub, &
1690 rk_step, dt_rk, grid%spec_zone, &
1692 ids, ide, jds, jde, kds, kde, &
1693 ims, ime, jms, jme, kms, kme, &
1694 grid%i_start(ij), grid%i_end(ij), &
1695 grid%j_start(ij), grid%j_end(ij), &
1699 !$OMP END PARALLEL DO
1701 !---------------------- positive definite bc call
1703 IF (config_flags%scalar_adv_opt /= ORIGINAL) THEN
1705 IF ( config_flags%h_sca_adv_order <= 4 ) THEN
1706 # include "HALO_EM_SCALAR_OLD_E_5.inc"
1707 ELSE IF ( config_flags%h_sca_adv_order <= 6 ) THEN
1708 # include "HALO_EM_SCALAR_OLD_E_7.inc"
1710 WRITE(wrf_err_message,*)'solve_em: invalid h_sca_adv_order = ',config_flags%h_sca_adv_order
1711 CALL wrf_error_fatal(TRIM(wrf_err_message))
1714 WRITE(wrf_err_message,*)'cannot use pd scheme with RSL - use RSL-LITE'
1715 CALL wrf_error_fatal(TRIM(wrf_err_message))
1721 # include "PERIOD_BDY_EM_SCALAR_OLD.inc"
1724 !$OMP PRIVATE ( ij )
1726 DO ij = 1 , grid%num_tiles
1727 IF (num_3d_s >= PARAM_FIRST_SCALAR) THEN
1728 DO im = PARAM_FIRST_SCALAR , num_3d_s
1729 CALL set_physical_bc3d( scalar_old(ims,kms,jms,im), 'p', config_flags, &
1730 ids, ide, jds, jde, kds, kde, &
1731 ims, ime, jms, jme, kms, kme, &
1732 ips, ipe, jps, jpe, kps, kpe, &
1733 grid%i_start(ij), grid%i_end(ij), &
1734 grid%j_start(ij), grid%j_end(ij), &
1739 !$OMP END PARALLEL DO
1741 END IF ! end if for scalar_adv_opt
1745 IF ((config_flags%chem_adv_opt /= ORIGINAL) .and. (rk_step == rk_order)) THEN
1748 !$OMP PRIVATE ( ij )
1749 DO ij = 1 , grid%num_tiles
1750 CALL wrf_debug ( 200 , ' call rk_update_scalar_pd' )
1751 DO im = PARAM_FIRST_SCALAR, num_3d_c
1752 CALL rk_update_scalar_pd( im, im, &
1753 chem_old(ims,kms,jms,im), &
1754 chem_tend(ims,kms,jms,im), &
1755 grid%mu_1, grid%mu_1, grid%mub, &
1756 rk_step, dt_rk, grid%spec_zone, &
1758 ids, ide, jds, jde, kds, kde, &
1759 ims, ime, jms, jme, kms, kme, &
1760 grid%i_start(ij), grid%i_end(ij), &
1761 grid%j_start(ij), grid%j_end(ij), &
1765 !$OMP END PARALLEL DO
1767 !---------------------- positive definite bc call
1769 IF (config_flags%chem_adv_opt /= ORIGINAL) THEN
1770 IF ( config_flags%h_sca_adv_order <= 4 ) THEN
1771 # include "HALO_EM_CHEM_OLD_E_5.inc"
1772 ELSE IF ( config_flags%h_sca_adv_order <= 6 ) THEN
1773 # include "HALO_EM_CHEM_OLD_E_7.inc"
1775 WRITE(wrf_err_message,*)'solve_em: invalid h_sca_adv_order = ',config_flags%h_sca_adv_order
1776 CALL wrf_error_fatal(TRIM(wrf_err_message))
1782 # include "PERIOD_BDY_EM_CHEM_OLD.inc"
1786 !$OMP PRIVATE ( ij )
1787 DO ij = 1 , grid%num_tiles
1788 IF (num_3d_c >= PARAM_FIRST_SCALAR) THEN
1789 DO im = PARAM_FIRST_SCALAR , num_3d_c
1790 CALL set_physical_bc3d( chem_old(ims,kms,jms,im), 'p', config_flags, &
1791 ids, ide, jds, jde, kds, kde, &
1792 ims, ime, jms, jme, kms, kme, &
1793 ips, ipe, jps, jpe, kps, kpe, &
1794 grid%i_start(ij), grid%i_end(ij), &
1795 grid%j_start(ij), grid%j_end(ij), &
1800 !$OMP END PARALLEL DO
1802 ENDIF ! end if for chem_adv_opt
1806 IF ((config_flags%tracer_adv_opt /= ORIGINAL) .and. (rk_step == rk_order)) THEN
1809 !$OMP PRIVATE ( ij )
1810 DO ij = 1 , grid%num_tiles
1811 CALL wrf_debug ( 200 , ' call rk_update_scalar_pd' )
1812 DO im = PARAM_FIRST_SCALAR, num_tracer
1813 CALL rk_update_scalar_pd( im, im, &
1814 tracer_old(ims,kms,jms,im), &
1815 tracer_tend(ims,kms,jms,im), &
1816 grid%mu_1, grid%mu_1, grid%mub, &
1817 rk_step, dt_rk, grid%spec_zone, &
1819 ids, ide, jds, jde, kds, kde, &
1820 ims, ime, jms, jme, kms, kme, &
1821 grid%i_start(ij), grid%i_end(ij), &
1822 grid%j_start(ij), grid%j_end(ij), &
1826 !$OMP END PARALLEL DO
1828 !---------------------- positive definite bc call
1830 IF (config_flags%tracer_adv_opt /= ORIGINAL) THEN
1831 IF ( config_flags%h_sca_adv_order <= 4 ) THEN
1832 # include "HALO_EM_TRACER_OLD_E_5.inc"
1833 ELSE IF ( config_flags%h_sca_adv_order <= 6 ) THEN
1834 # include "HALO_EM_TRACER_OLD_E_7.inc"
1836 WRITE(wrf_err_message,*)'solve_em: invalid h_sca_adv_order = ',config_flags%h_sca_adv_order
1837 CALL wrf_error_fatal(TRIM(wrf_err_message))
1843 # include "PERIOD_BDY_EM_TRACER_OLD.inc"
1847 !$OMP PRIVATE ( ij )
1848 DO ij = 1 , grid%num_tiles
1849 IF (num_tracer >= PARAM_FIRST_SCALAR) THEN
1850 DO im = PARAM_FIRST_SCALAR , num_tracer
1851 CALL set_physical_bc3d( tracer_old(ims,kms,jms,im), 'p', config_flags, &
1852 ids, ide, jds, jde, kds, kde, &
1853 ims, ime, jms, jme, kms, kme, &
1854 ips, ipe, jps, jpe, kps, kpe, &
1855 grid%i_start(ij), grid%i_end(ij), &
1856 grid%j_start(ij), grid%j_end(ij), &
1861 !$OMP END PARALLEL DO
1863 ENDIF ! end if for tracer_adv_opt
1867 IF ((config_flags%tke_adv_opt /= ORIGINAL) .and. (rk_step == rk_order) &
1868 .and. (config_flags%km_opt .eq. 2) ) THEN
1871 !$OMP PRIVATE ( ij )
1872 DO ij = 1 , grid%num_tiles
1873 CALL wrf_debug ( 200 , ' call rk_update_scalar_pd' )
1874 CALL rk_update_scalar_pd( 1, 1, &
1876 tke_tend(ims,kms,jms), &
1877 grid%mu_1, grid%mu_1, grid%mub, &
1878 rk_step, dt_rk, grid%spec_zone, &
1880 ids, ide, jds, jde, kds, kde, &
1881 ims, ime, jms, jme, kms, kme, &
1882 grid%i_start(ij), grid%i_end(ij), &
1883 grid%j_start(ij), grid%j_end(ij), &
1886 !$OMP END PARALLEL DO
1888 !---------------------- positive definite bc call
1890 IF (config_flags%tke_adv_opt /= ORIGINAL) THEN
1891 IF ( config_flags%h_sca_adv_order <= 4 ) THEN
1892 # include "HALO_EM_TKE_OLD_E_5.inc"
1893 ELSE IF ( config_flags%h_sca_adv_order <= 6 ) THEN
1894 # include "HALO_EM_TKE_OLD_E_7.inc"
1896 WRITE(wrf_err_message,*)'solve_em: invalid h_sca_adv_order = ',config_flags%h_sca_adv_order
1897 CALL wrf_error_fatal(TRIM(wrf_err_message))
1903 # include "PERIOD_BDY_EM_TKE_OLD.inc"
1907 !$OMP PRIVATE ( ij )
1908 DO ij = 1 , grid%num_tiles
1909 CALL set_physical_bc3d( grid%tke_1, 'p', config_flags, &
1910 ids, ide, jds, jde, kds, kde, &
1911 ims, ime, jms, jme, kms, kme, &
1912 ips, ipe, jps, jpe, kps, kpe, &
1913 grid%i_start(ij), grid%i_end(ij), &
1914 grid%j_start(ij), grid%j_end(ij), &
1917 !$OMP END PARALLEL DO
1919 !--- end of positive definite physics tendency update
1921 END IF ! end if for tke_adv_opt
1925 ! Stencils for patch communications (WCS, 29 June 2001)
1938 !--------------------------------------------------------------
1940 # include "HALO_EM_D.inc"
1941 ! WCS addition 11/19/08
1942 # include "PERIOD_EM_DA.inc"
1947 ! (4) Still within the RK loop, the scalar variables are advanced.
1949 ! For the moist and chem variables, each one is advanced
1950 ! individually, using named loops "moist_variable_loop:"
1951 ! and "chem_variable_loop:". Each RK substep begins by
1952 ! calculating the advective tendency, and, for the first RK step,
1953 ! 3D mixing (calling rk_scalar_tend) followed by an update
1954 ! of the scalar (calling rk_update_scalar).
1959 moist_scalar_advance: IF (num_3d_m >= PARAM_FIRST_SCALAR ) THEN
1961 moist_variable_loop: DO im = PARAM_FIRST_SCALAR, num_3d_m
1963 ! adv_moist_cond is set in module_physics_init based on mp_physics choice
1964 ! true except for Ferrier scheme
1966 IF (grid%adv_moist_cond .or. im==p_qv ) THEN
1969 !$OMP PRIVATE ( ij )
1970 moist_tile_loop_1: DO ij = 1 , grid%num_tiles
1972 CALL wrf_debug ( 200 , ' call rk_scalar_tend' )
1975 BENCH_START(rk_scalar_tend_tim)
1976 CALL rk_scalar_tend ( im, im, config_flags, tenddec, &
1978 grid%ru_m, grid%rv_m, grid%ww_m, &
1979 grid%muts, grid%mub, grid%mu_1, &
1981 moist_old(ims,kms,jms,im), &
1982 moist(ims,kms,jms,im), &
1983 moist_tend(ims,kms,jms,im), &
1984 advect_tend,h_tendency,z_tendency,grid%rqvften, &
1985 grid%qv_base, .true., grid%fnm, grid%fnp, &
1986 grid%msfux,grid%msfuy, grid%msfvx, grid%msfvx_inv,&
1987 grid%msfvy, grid%msftx,grid%msfty, &
1988 grid%rdx, grid%rdy, grid%rdn, grid%rdnw, grid%khdif, &
1989 grid%kvdif, grid%xkhh, &
1990 grid%diff_6th_opt, grid%diff_6th_factor, &
1991 config_flags%moist_adv_opt, &
1992 ids, ide, jds, jde, kds, kde, &
1993 ims, ime, jms, jme, kms, kme, &
1994 grid%i_start(ij), grid%i_end(ij), &
1995 grid%j_start(ij), grid%j_end(ij), &
1998 BENCH_END(rk_scalar_tend_tim)
2000 BENCH_START(rlx_bdy_scalar_tim)
2001 IF( ( config_flags%specified .or. config_flags%nested ) .and. rk_step == 1 ) THEN
2002 IF ( im .EQ. P_QV .OR. config_flags%nested ) THEN
2003 CALL relax_bdy_scalar ( moist_tend(ims,kms,jms,im), &
2004 moist(ims,kms,jms,im), grid%mut, &
2005 moist_bxs(jms,kms,1,im),moist_bxe(jms,kms,1,im), &
2006 moist_bys(ims,kms,1,im),moist_bye(ims,kms,1,im), &
2007 moist_btxs(jms,kms,1,im),moist_btxe(jms,kms,1,im), &
2008 moist_btys(ims,kms,1,im),moist_btye(ims,kms,1,im), &
2009 config_flags%spec_bdy_width, grid%spec_zone, grid%relax_zone, &
2010 grid%dtbc, grid%fcx, grid%gcx, &
2012 ids,ide, jds,jde, kds,kde, & ! domain dims
2013 ims,ime, jms,jme, kms,kme, & ! memory dims
2014 ips,ipe, jps,jpe, kps,kpe, & ! patch dims
2015 grid%i_start(ij), grid%i_end(ij), &
2016 grid%j_start(ij), grid%j_end(ij), &
2019 CALL spec_bdy_scalar ( moist_tend(ims,kms,jms,im), &
2020 moist_bxs(jms,kms,1,im),moist_bxe(jms,kms,1,im), &
2021 moist_bys(ims,kms,1,im),moist_bye(ims,kms,1,im), &
2022 moist_btxs(jms,kms,1,im),moist_btxe(jms,kms,1,im), &
2023 moist_btys(ims,kms,1,im),moist_btye(ims,kms,1,im), &
2024 config_flags%spec_bdy_width, grid%spec_zone, &
2026 ids,ide, jds,jde, kds,kde, & ! domain dims
2027 ims,ime, jms,jme, kms,kme, & ! memory dims
2028 ips,ipe, jps,jpe, kps,kpe, & ! patch dims
2029 grid%i_start(ij), grid%i_end(ij), &
2030 grid%j_start(ij), grid%j_end(ij), &
2034 BENCH_END(rlx_bdy_scalar_tim)
2036 ENDDO moist_tile_loop_1
2037 !$OMP END PARALLEL DO
2040 !$OMP PRIVATE ( ij )
2041 moist_tile_loop_2: DO ij = 1 , grid%num_tiles
2043 CALL wrf_debug ( 200 , ' call rk_update_scalar' )
2046 BENCH_START(update_scal_tim)
2047 CALL rk_update_scalar( scs=im, sce=im, &
2048 scalar_1=moist_old(ims,kms,jms,im), &
2049 scalar_2=moist(ims,kms,jms,im), &
2050 sc_tend=moist_tend(ims,kms,jms,im), &
2051 ! advh_t=advh_t(ims,kms,jms,1), &
2052 ! advz_t=advz_t(ims,kms,jms,1), &
2053 advect_tend=advect_tend, &
2054 h_tendency=h_tendency, z_tendency=z_tendency, &
2055 msftx=grid%msftx,msfty=grid%msfty, &
2056 mu_old=grid%mu_1, mu_new=grid%mu_2, mu_base=grid%mub, &
2057 rk_step=rk_step, dt=dt_rk, spec_zone=grid%spec_zone, &
2058 config_flags=config_flags, tenddec=tenddec, &
2059 ids=ids, ide=ide, jds=jds, jde=jde, kds=kds, kde=kde, &
2060 ims=ims, ime=ime, jms=jms, jme=jme, kms=kms, kme=kme, &
2061 its=grid%i_start(ij), ite=grid%i_end(ij), &
2062 jts=grid%j_start(ij), jte=grid%j_end(ij), &
2063 kts=k_start , kte=k_end )
2064 BENCH_END(update_scal_tim)
2066 BENCH_START(flow_depbdy_tim)
2067 IF( config_flags%specified ) THEN
2068 IF(im .ne. P_QV)THEN
2069 CALL flow_dep_bdy ( moist(ims,kms,jms,im), &
2070 grid%ru_m, grid%rv_m, config_flags, &
2072 ids,ide, jds,jde, kds,kde, &
2073 ims,ime, jms,jme, kms,kme, &
2074 ips,ipe, jps,jpe, kps,kpe, &
2075 grid%i_start(ij), grid%i_end(ij), &
2076 grid%j_start(ij), grid%j_end(ij), &
2080 BENCH_END(flow_depbdy_tim)
2082 ENDDO moist_tile_loop_2
2083 !$OMP END PARALLEL DO
2085 ENDIF !-- if (grid%adv_moist_cond .or. im==p_qv ) then
2087 ENDDO moist_variable_loop
2089 ENDIF moist_scalar_advance
2091 BENCH_START(tke_adv_tim)
2092 TKE_advance: IF (config_flags%km_opt .eq. 2) then
2094 IF ( config_flags%h_mom_adv_order <= 4 ) THEN
2095 # include "HALO_EM_TKE_ADVECT_3.inc"
2096 ELSE IF ( config_flags%h_mom_adv_order <= 6 ) THEN
2097 # include "HALO_EM_TKE_ADVECT_5.inc"
2099 WRITE(wrf_err_message,*)'solve_em: invalid h_mom_adv_order = ',config_flags%h_mom_adv_order
2100 CALL wrf_error_fatal(TRIM(wrf_err_message))
2104 !$OMP PRIVATE ( ij )
2105 tke_tile_loop_1: DO ij = 1 , grid%num_tiles
2107 CALL wrf_debug ( 200 , ' call rk_scalar_tend for tke' )
2109 CALL rk_scalar_tend ( 1, 1, config_flags, tenddec, &
2111 grid%ru_m, grid%rv_m, grid%ww_m, &
2112 grid%muts, grid%mub, grid%mu_1, &
2116 tke_tend(ims,kms,jms), &
2117 advect_tend,h_tendency,z_tendency,grid%rqvften, &
2118 grid%qv_base, .false., grid%fnm, grid%fnp, &
2119 grid%msfux,grid%msfuy, grid%msfvx, grid%msfvx_inv, &
2120 grid%msfvy, grid%msftx,grid%msfty, &
2121 grid%rdx, grid%rdy, grid%rdn, grid%rdnw, grid%khdif, &
2122 grid%kvdif, grid%xkhh, &
2123 grid%diff_6th_opt, grid%diff_6th_factor, &
2124 config_flags%tke_adv_opt, &
2125 ids, ide, jds, jde, kds, kde, &
2126 ims, ime, jms, jme, kms, kme, &
2127 grid%i_start(ij), grid%i_end(ij), &
2128 grid%j_start(ij), grid%j_end(ij), &
2131 ENDDO tke_tile_loop_1
2132 !$OMP END PARALLEL DO
2135 !$OMP PRIVATE ( ij )
2136 tke_tile_loop_2: DO ij = 1 , grid%num_tiles
2138 CALL wrf_debug ( 200 , ' call rk_update_scalar' )
2140 CALL rk_update_scalar( scs=1, sce=1, &
2141 scalar_1=grid%tke_1, &
2142 scalar_2=grid%tke_2, &
2143 sc_tend=tke_tend(ims,kms,jms), &
2144 ! advh_t=advh_t(ims,kms,jms,1), &
2145 ! advz_t=advz_t(ims,kms,jms,1), &
2146 advect_tend=advect_tend, &
2147 h_tendency=h_tendency, z_tendency=z_tendency, &
2148 msftx=grid%msftx,msfty=grid%msfty, &
2149 mu_old=grid%mu_1, mu_new=grid%mu_2, mu_base=grid%mub, &
2150 rk_step=rk_step, dt=dt_rk, spec_zone=grid%spec_zone, &
2151 config_flags=config_flags, tenddec=tenddec, &
2152 ids=ids, ide=ide, jds=jds, jde=jde, kds=kds, kde=kde, &
2153 ims=ims, ime=ime, jms=jms, jme=jme, kms=kms, kme=kme, &
2154 its=grid%i_start(ij), ite=grid%i_end(ij), &
2155 jts=grid%j_start(ij), jte=grid%j_end(ij), &
2156 kts=k_start , kte=k_end )
2158 ! bound the tke (greater than 0, less than tke_upper_bound)
2160 CALL bound_tke( grid%tke_2, grid%tke_upper_bound, &
2161 ids, ide, jds, jde, kds, kde, &
2162 ims, ime, jms, jme, kms, kme, &
2163 grid%i_start(ij), grid%i_end(ij), &
2164 grid%j_start(ij), grid%j_end(ij), &
2167 IF( config_flags%specified .or. config_flags%nested ) THEN
2168 CALL flow_dep_bdy ( grid%tke_2, &
2169 grid%ru_m, grid%rv_m, config_flags, &
2171 ids,ide, jds,jde, kds,kde, & ! domain dims
2172 ims,ime, jms,jme, kms,kme, & ! memory dims
2173 ips,ipe, jps,jpe, kps,kpe, & ! patch dims
2174 grid%i_start(ij), grid%i_end(ij), &
2175 grid%j_start(ij), grid%j_end(ij), &
2178 ENDDO tke_tile_loop_2
2179 !$OMP END PARALLEL DO
2182 BENCH_END(tke_adv_tim)
2185 ! next the chemical species
2186 BENCH_START(chem_adv_tim)
2187 chem_scalar_advance: IF (num_3d_c >= PARAM_FIRST_SCALAR) THEN
2189 chem_variable_loop: DO ic = PARAM_FIRST_SCALAR, num_3d_c
2192 !$OMP PRIVATE ( ij )
2193 chem_tile_loop_1: DO ij = 1 , grid%num_tiles
2195 CALL wrf_debug ( 200 , ' call rk_scalar_tend in chem_tile_loop_1' )
2196 tenddec = (( config_flags%chemdiag == USECHEMDIAG ) .and. &
2197 ( adv_ct_indices(ic) >= PARAM_FIRST_SCALAR ))
2198 CALL rk_scalar_tend ( ic, ic, config_flags, tenddec, &
2200 grid%ru_m, grid%rv_m, grid%ww_m, &
2201 grid%muts, grid%mub, grid%mu_1, &
2203 chem_old(ims,kms,jms,ic), &
2204 chem(ims,kms,jms,ic), &
2205 chem_tend(ims,kms,jms,ic), &
2206 advect_tend,h_tendency,z_tendency,grid%rqvften, &
2207 grid%qv_base, .false., grid%fnm, grid%fnp, &
2208 grid%msfux,grid%msfuy, grid%msfvx, grid%msfvx_inv, &
2209 grid%msfvy, grid%msftx,grid%msfty, &
2210 grid%rdx, grid%rdy, grid%rdn, grid%rdnw, &
2211 grid%khdif, grid%kvdif, grid%xkhh, &
2212 grid%diff_6th_opt, grid%diff_6th_factor, &
2213 config_flags%chem_adv_opt, &
2214 ids, ide, jds, jde, kds, kde, &
2215 ims, ime, jms, jme, kms, kme, &
2216 grid%i_start(ij), grid%i_end(ij), &
2217 grid%j_start(ij), grid%j_end(ij), &
2220 ! Currently, chemistry species with specified boundaries (i.e. the mother
2221 ! domain) are being over written by flow_dep_bdy_chem. So, relax_bdy and
2222 ! spec_bdy are only called for nests. For boundary conditions from global model or larger domain,
2223 ! chem is uncoupled, and only used for one row/column on inflow (if have_bcs_chem=.true.)
2225 IF( ( config_flags%nested ) .and. rk_step == 1 ) THEN
2226 IF(ic.eq.1)CALL wrf_debug ( 10 , ' have_bcs_chem' )
2227 CALL relax_bdy_scalar ( chem_tend(ims,kms,jms,ic), &
2228 chem(ims,kms,jms,ic), grid%mut, &
2229 chem_bxs(jms,kms,1,ic),chem_bxe(jms,kms,1,ic), &
2230 chem_bys(ims,kms,1,ic),chem_bye(ims,kms,1,ic), &
2231 chem_btxs(jms,kms,1,ic),chem_btxe(jms,kms,1,ic), &
2232 chem_btys(ims,kms,1,ic),chem_btye(ims,kms,1,ic), &
2233 config_flags%spec_bdy_width, grid%spec_zone, grid%relax_zone, &
2234 grid%dtbc, grid%fcx, grid%gcx, &
2236 ids,ide, jds,jde, kds,kde, &
2237 ims,ime, jms,jme, kms,kme, &
2238 ips,ipe, jps,jpe, kps,kpe, &
2239 grid%i_start(ij), grid%i_end(ij), &
2240 grid%j_start(ij), grid%j_end(ij), &
2242 CALL spec_bdy_scalar ( chem_tend(ims,kms,jms,ic), &
2243 chem_bxs(jms,kms,1,ic),chem_bxe(jms,kms,1,ic), &
2244 chem_bys(ims,kms,1,ic),chem_bye(ims,kms,1,ic), &
2245 chem_btxs(jms,kms,1,ic),chem_btxe(jms,kms,1,ic), &
2246 chem_btys(ims,kms,1,ic),chem_btye(ims,kms,1,ic), &
2247 config_flags%spec_bdy_width, grid%spec_zone, &
2249 ids,ide, jds,jde, kds,kde, &
2250 ims,ime, jms,jme, kms,kme, &
2251 ips,ipe, jps,jpe, kps,kpe, &
2252 grid%i_start(ij), grid%i_end(ij), &
2253 grid%j_start(ij), grid%j_end(ij), &
2257 ENDDO chem_tile_loop_1
2258 !$OMP END PARALLEL DO
2261 !$OMP PRIVATE ( ij )
2263 chem_tile_loop_2: DO ij = 1 , grid%num_tiles
2265 CALL wrf_debug ( 200 , ' call rk_update_scalar' )
2266 tenddec = (( config_flags%chemdiag == USECHEMDIAG ) .and. &
2267 ( adv_ct_indices(ic) >= PARAM_FIRST_SCALAR ))
2268 CALL rk_update_scalar( scs=ic, sce=ic, &
2269 scalar_1=chem_old(ims,kms,jms,ic), &
2270 scalar_2=chem(ims,kms,jms,ic), &
2271 sc_tend=chem_tend(ims,kms,jms,ic), &
2272 advh_t=advh_ct(ims,kms,jms,adv_ct_indices(ic)), &
2273 advz_t=advz_ct(ims,kms,jms,adv_ct_indices(ic)), &
2274 advect_tend=advect_tend, &
2275 h_tendency=h_tendency, z_tendency=z_tendency, &
2276 msftx=grid%msftx,msfty=grid%msfty, &
2277 mu_old=grid%mu_1, mu_new=grid%mu_2, mu_base=grid%mub, &
2278 rk_step=rk_step, dt=dt_rk, spec_zone=grid%spec_zone, &
2279 config_flags=config_flags, tenddec=tenddec, &
2280 ids=ids, ide=ide, jds=jds, jde=jde, kds=kds, kde=kde, &
2281 ims=ims, ime=ime, jms=jms, jme=jme, kms=kms, kme=kme, &
2282 its=grid%i_start(ij), ite=grid%i_end(ij), &
2283 jts=grid%j_start(ij), jte=grid%j_end(ij), &
2284 kts=k_start , kte=k_end )
2286 IF( config_flags%specified ) THEN
2287 CALL flow_dep_bdy_chem( chem(ims,kms,jms,ic), &
2288 chem_bxs(jms,kms,1,ic), chem_btxs(jms,kms,1,ic), &
2289 chem_bxe(jms,kms,1,ic), chem_btxe(jms,kms,1,ic), &
2290 chem_bys(ims,kms,1,ic), chem_btys(ims,kms,1,ic), &
2291 chem_bye(ims,kms,1,ic), chem_btye(ims,kms,1,ic), &
2293 config_flags%spec_bdy_width,grid%z, &
2294 grid%have_bcs_chem, &
2295 grid%ru_m, grid%rv_m, config_flags,grid%alt, &
2296 grid%t_1,grid%pb,grid%p,t0,p1000mb,rcp,grid%ph_2,grid%phb,g, &
2297 grid%spec_zone,ic, &
2298 ids,ide, jds,jde, kds,kde, & ! domain dims
2299 ims,ime, jms,jme, kms,kme, & ! memory dims
2300 ips,ipe, jps,jpe, kps,kpe, & ! patch dims
2301 grid%i_start(ij), grid%i_end(ij), &
2302 grid%j_start(ij), grid%j_end(ij), &
2305 ENDDO chem_tile_loop_2
2306 !$OMP END PARALLEL DO
2308 ENDDO chem_variable_loop
2309 ENDIF chem_scalar_advance
2310 BENCH_END(chem_adv_tim)
2312 ! next the chemical species
2313 BENCH_START(tracer_adv_tim)
2314 tracer_advance: IF (num_tracer >= PARAM_FIRST_SCALAR) THEN
2316 tracer_variable_loop: DO ic = PARAM_FIRST_SCALAR, num_tracer
2319 !$OMP PRIVATE ( ij )
2320 tracer_tile_loop_1: DO ij = 1 , grid%num_tiles
2322 CALL wrf_debug ( 15 , ' call rk_scalar_tend in tracer_tile_loop_1' )
2324 CALL rk_scalar_tend ( ic, ic, config_flags, tenddec, &
2326 grid%ru_m, grid%rv_m, grid%ww_m, &
2327 grid%muts, grid%mub, grid%mu_1, &
2329 tracer_old(ims,kms,jms,ic), &
2330 tracer(ims,kms,jms,ic), &
2331 tracer_tend(ims,kms,jms,ic), &
2332 advect_tend,h_tendency,z_tendency,grid%rqvften, &
2333 grid%qv_base, .false., grid%fnm, grid%fnp, &
2334 grid%msfux,grid%msfuy, grid%msfvx, grid%msfvx_inv, &
2335 grid%msfvy, grid%msftx,grid%msfty, &
2336 grid%rdx, grid%rdy, grid%rdn, grid%rdnw, &
2337 grid%khdif, grid%kvdif, grid%xkhh, &
2338 grid%diff_6th_opt, grid%diff_6th_factor, &
2339 config_flags%tracer_adv_opt, &
2340 ids, ide, jds, jde, kds, kde, &
2341 ims, ime, jms, jme, kms, kme, &
2342 grid%i_start(ij), grid%i_end(ij), &
2343 grid%j_start(ij), grid%j_end(ij), &
2346 ! Currently, chemistry species with specified boundaries (i.e. the mother
2347 ! domain) are being over written by flow_dep_bdy_chem. So, relax_bdy and
2348 ! spec_bdy are only called for nests. For boundary conditions from global model or larger domain,
2349 ! chem is uncoupled, and only used for one row/column on inflow (if have_bcs_chem=.true.)
2351 IF( ( config_flags%nested ) .and. rk_step == 1 ) THEN
2352 IF(ic.eq.1)CALL wrf_debug ( 10 , ' have_bcs_tracer' )
2353 CALL relax_bdy_scalar ( tracer_tend(ims,kms,jms,ic), &
2354 tracer(ims,kms,jms,ic), grid%mut, &
2355 tracer_bxs(jms,kms,1,ic),tracer_bxe(jms,kms,1,ic), &
2356 tracer_bys(ims,kms,1,ic),tracer_bye(ims,kms,1,ic), &
2357 tracer_btxs(jms,kms,1,ic),tracer_btxe(jms,kms,1,ic), &
2358 tracer_btys(ims,kms,1,ic),tracer_btye(ims,kms,1,ic), &
2359 config_flags%spec_bdy_width, grid%spec_zone, grid%relax_zone, &
2360 grid%dtbc, grid%fcx, grid%gcx, &
2362 ids,ide, jds,jde, kds,kde, &
2363 ims,ime, jms,jme, kms,kme, &
2364 ips,ipe, jps,jpe, kps,kpe, &
2365 grid%i_start(ij), grid%i_end(ij), &
2366 grid%j_start(ij), grid%j_end(ij), &
2368 CALL spec_bdy_scalar ( tracer_tend(ims,kms,jms,ic), &
2369 tracer_bxs(jms,kms,1,ic),tracer_bxe(jms,kms,1,ic), &
2370 tracer_bys(ims,kms,1,ic),tracer_bye(ims,kms,1,ic), &
2371 tracer_btxs(jms,kms,1,ic),tracer_btxe(jms,kms,1,ic), &
2372 tracer_btys(ims,kms,1,ic),tracer_btye(ims,kms,1,ic), &
2373 config_flags%spec_bdy_width, grid%spec_zone, &
2375 ids,ide, jds,jde, kds,kde, &
2376 ims,ime, jms,jme, kms,kme, &
2377 ips,ipe, jps,jpe, kps,kpe, &
2378 grid%i_start(ij), grid%i_end(ij), &
2379 grid%j_start(ij), grid%j_end(ij), &
2383 ENDDO tracer_tile_loop_1
2384 !$OMP END PARALLEL DO
2387 !$OMP PRIVATE ( ij )
2389 tracer_tile_loop_2: DO ij = 1 , grid%num_tiles
2391 CALL wrf_debug ( 200 , ' call rk_update_scalar' )
2393 CALL rk_update_scalar( scs=ic, sce=ic, &
2394 scalar_1=tracer_old(ims,kms,jms,ic), &
2395 scalar_2=tracer(ims,kms,jms,ic), &
2396 sc_tend=tracer_tend(ims,kms,jms,ic), &
2397 ! advh_t=advh_t(ims,kms,jms,1), &
2398 ! advz_t=advz_t(ims,kms,jms,1), &
2399 advect_tend=advect_tend, &
2400 h_tendency=h_tendency, z_tendency=z_tendency, &
2401 msftx=grid%msftx,msfty=grid%msfty, &
2402 mu_old=grid%mu_1, mu_new=grid%mu_2, mu_base=grid%mub, &
2403 rk_step=rk_step, dt=dt_rk, spec_zone=grid%spec_zone, &
2404 config_flags=config_flags, tenddec=tenddec, &
2405 ids=ids, ide=ide, jds=jds, jde=jde, kds=kds, kde=kde, &
2406 ims=ims, ime=ime, jms=jms, jme=jme, kms=kms, kme=kme, &
2407 its=grid%i_start(ij), ite=grid%i_end(ij), &
2408 jts=grid%j_start(ij), jte=grid%j_end(ij), &
2409 kts=k_start , kte=k_end )
2411 IF( config_flags%specified ) THEN
2413 CALL flow_dep_bdy_tracer( tracer(ims,kms,jms,ic), &
2414 tracer_bxs(jms,kms,1,ic), tracer_btxs(jms,kms,1,ic), &
2415 tracer_bxe(jms,kms,1,ic), tracer_btxe(jms,kms,1,ic), &
2416 tracer_bys(ims,kms,1,ic), tracer_btys(ims,kms,1,ic), &
2417 tracer_bye(ims,kms,1,ic), tracer_btye(ims,kms,1,ic), &
2419 config_flags%spec_bdy_width,grid%z, &
2420 grid%have_bcs_tracer, &
2421 grid%ru_m, grid%rv_m, config_flags%tracer_opt,grid%alt, &
2422 grid%t_1,grid%pb,grid%p,t0,p1000mb,rcp,grid%ph_2,grid%phb,g, &
2423 grid%spec_zone,ic, &
2424 ids,ide, jds,jde, kds,kde, & ! domain dims
2425 ims,ime, jms,jme, kms,kme, & ! memory dims
2426 ips,ipe, jps,jpe, kps,kpe, & ! patch dims
2427 grid%i_start(ij), grid%i_end(ij), &
2428 grid%j_start(ij), grid%j_end(ij), &
2431 CALL flow_dep_bdy ( tracer(ims,kms,jms,ic), &
2432 grid%ru_m, grid%rv_m, config_flags, &
2434 ids,ide, jds,jde, kds,kde, & ! domain dims
2435 ims,ime, jms,jme, kms,kme, & ! memory dims
2436 ips,ipe, jps,jpe, kps,kpe, & ! patch dims
2437 grid%i_start(ij), grid%i_end(ij), &
2438 grid%j_start(ij), grid%j_end(ij), &
2442 ENDDO tracer_tile_loop_2
2443 !$OMP END PARALLEL DO
2445 ENDDO tracer_variable_loop
2446 ENDIF tracer_advance
2447 BENCH_END(tracer_adv_tim)
2449 ! next the other scalar species
2450 other_scalar_advance: IF (num_3d_s >= PARAM_FIRST_SCALAR) THEN
2452 scalar_variable_loop: do is = PARAM_FIRST_SCALAR, num_3d_s
2454 !$OMP PRIVATE ( ij )
2455 scalar_tile_loop_1: DO ij = 1 , grid%num_tiles
2457 CALL wrf_debug ( 200 , ' call rk_scalar_tend' )
2459 CALL rk_scalar_tend ( is, is, config_flags, tenddec, &
2461 grid%ru_m, grid%rv_m, grid%ww_m, &
2462 grid%muts, grid%mub, grid%mu_1, &
2464 scalar_old(ims,kms,jms,is), &
2465 scalar(ims,kms,jms,is), &
2466 scalar_tend(ims,kms,jms,is), &
2467 advect_tend,h_tendency,z_tendency,grid%rqvften, &
2468 grid%qv_base, .false., grid%fnm, grid%fnp, &
2469 grid%msfux,grid%msfuy, grid%msfvx, grid%msfvx_inv, &
2470 grid%msfvy, grid%msftx,grid%msfty, &
2471 grid%rdx, grid%rdy, grid%rdn, grid%rdnw, &
2472 grid%khdif, grid%kvdif, grid%xkhh, &
2473 grid%diff_6th_opt, grid%diff_6th_factor, &
2474 config_flags%scalar_adv_opt, &
2475 ids, ide, jds, jde, kds, kde, &
2476 ims, ime, jms, jme, kms, kme, &
2477 grid%i_start(ij), grid%i_end(ij), &
2478 grid%j_start(ij), grid%j_end(ij), &
2481 IF( config_flags%nested .and. (rk_step == 1) ) THEN
2483 CALL relax_bdy_scalar ( scalar_tend(ims,kms,jms,is), &
2484 scalar(ims,kms,jms,is), grid%mut, &
2485 scalar_bxs(jms,kms,1,is),scalar_bxe(jms,kms,1,is), &
2486 scalar_bys(ims,kms,1,is),scalar_bye(ims,kms,1,is), &
2487 scalar_btxs(jms,kms,1,is),scalar_btxe(jms,kms,1,is), &
2488 scalar_btys(ims,kms,1,is),scalar_btye(ims,kms,1,is), &
2489 config_flags%spec_bdy_width, grid%spec_zone, grid%relax_zone, &
2490 grid%dtbc, grid%fcx, grid%gcx, &
2492 ids,ide, jds,jde, kds,kde, &
2493 ims,ime, jms,jme, kms,kme, &
2494 ips,ipe, jps,jpe, kps,kpe, &
2495 grid%i_start(ij), grid%i_end(ij), &
2496 grid%j_start(ij), grid%j_end(ij), &
2499 CALL spec_bdy_scalar ( scalar_tend(ims,kms,jms,is), &
2500 scalar_bxs(jms,kms,1,is),scalar_bxe(jms,kms,1,is), &
2501 scalar_bys(ims,kms,1,is),scalar_bye(ims,kms,1,is), &
2502 scalar_btxs(jms,kms,1,is),scalar_btxe(jms,kms,1,is), &
2503 scalar_btys(ims,kms,1,is),scalar_btye(ims,kms,1,is), &
2504 config_flags%spec_bdy_width, grid%spec_zone, &
2506 ids,ide, jds,jde, kds,kde, &
2507 ims,ime, jms,jme, kms,kme, &
2508 ips,ipe, jps,jpe, kps,kpe, &
2509 grid%i_start(ij), grid%i_end(ij), &
2510 grid%j_start(ij), grid%j_end(ij), &
2513 ENDIF ! b.c test for chem nested boundary condition
2515 ENDDO scalar_tile_loop_1
2516 !$OMP END PARALLEL DO
2519 !$OMP PRIVATE ( ij )
2520 scalar_tile_loop_2: DO ij = 1 , grid%num_tiles
2522 CALL wrf_debug ( 200 , ' call rk_update_scalar' )
2524 CALL rk_update_scalar( scs=is, sce=is, &
2525 scalar_1=scalar_old(ims,kms,jms,is), &
2526 scalar_2=scalar(ims,kms,jms,is), &
2527 sc_tend=scalar_tend(ims,kms,jms,is), &
2528 ! advh_t=advh_t(ims,kms,jms,1), &
2529 ! advz_t=advz_t(ims,kms,jms,1), &
2530 advect_tend=advect_tend, &
2531 h_tendency=h_tendency, z_tendency=z_tendency, &
2532 msftx=grid%msftx,msfty=grid%msfty, &
2533 mu_old=grid%mu_1, mu_new=grid%mu_2, mu_base=grid%mub, &
2534 rk_step=rk_step, dt=dt_rk, spec_zone=grid%spec_zone, &
2535 config_flags=config_flags, tenddec=tenddec, &
2536 ids=ids, ide=ide, jds=jds, jde=jde, kds=kds, kde=kde, &
2537 ims=ims, ime=ime, jms=jms, jme=jme, kms=kms, kme=kme, &
2538 its=grid%i_start(ij), ite=grid%i_end(ij), &
2539 jts=grid%j_start(ij), jte=grid%j_end(ij), &
2540 kts=k_start , kte=k_end )
2542 IF( config_flags%specified ) THEN
2544 CALL flow_dep_bdy ( scalar(ims,kms,jms,is), &
2545 grid%ru_m, grid%rv_m, config_flags, &
2547 ids,ide, jds,jde, kds,kde, & ! domain dims
2548 ims,ime, jms,jme, kms,kme, & ! memory dims
2549 ips,ipe, jps,jpe, kps,kpe, & ! patch dims
2550 grid%i_start(ij), grid%i_end(ij), &
2551 grid%j_start(ij), grid%j_end(ij), &
2555 ENDDO scalar_tile_loop_2
2556 !$OMP END PARALLEL DO
2558 ENDDO scalar_variable_loop
2560 ENDIF other_scalar_advance
2562 ! update the pressure and density at the new time level
2565 !$OMP PRIVATE ( ij )
2566 DO ij = 1 , grid%num_tiles
2568 BENCH_START(calc_p_rho_tim)
2570 CALL calc_p_rho_phi( moist, num_3d_m, config_flags%hypsometric_opt, &
2571 grid%al, grid%alb, grid%mu_2, grid%muts, &
2572 grid%ph_2, grid%phb, grid%p, grid%pb, grid%t_2, &
2573 p0, t0, grid%p_top, grid%znu, grid%znw, grid%dnw, grid%rdnw, &
2574 grid%rdn, config_flags%non_hydrostatic, &
2575 ids, ide, jds, jde, kds, kde, &
2576 ims, ime, jms, jme, kms, kme, &
2577 grid%i_start(ij), grid%i_end(ij), &
2578 grid%j_start(ij), grid%j_end(ij), &
2581 BENCH_END(calc_p_rho_tim)
2584 !$OMP END PARALLEL DO
2586 ! Reset the boundary conditions if there is another corrector step.
2587 ! (rk_step < rk_order), else we'll handle it at the end of everything
2588 ! (after the split physics, before exiting the timestep).
2590 rk_step_1_check: IF ( rk_step < rk_order ) THEN
2592 !-----------------------------------------------------------
2593 ! rk3 substep polar filter for scalars (moist,chem,scalar)
2594 !-----------------------------------------------------------
2596 IF (config_flags%polar) THEN
2597 IF ( num_3d_m >= PARAM_FIRST_SCALAR ) THEN
2598 CALL wrf_debug ( 200 , ' call filter moist ' )
2599 DO im = PARAM_FIRST_SCALAR, num_3d_m
2600 CALL couple_scalars_for_filter ( FIELD=moist(ims,kms,jms,im) &
2601 ,MU=grid%mu_2 , MUB=grid%mub &
2602 ,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
2603 ,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
2604 ,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe )
2605 CALL pxft ( grid=grid &
2618 ,positive_definite=.FALSE. &
2619 ,moist=moist,chem=chem,tracer=tracer,scalar=scalar &
2620 ,fft_filter_lat = config_flags%fft_filter_lat &
2622 ,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
2623 ,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
2624 ,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe &
2625 ,imsx=imsx,imex=imex,jmsx=jmsx,jmex=jmex,kmsx=kmsx,kmex=kmex &
2626 ,ipsx=ipsx,ipex=ipex,jpsx=jmsx,jpex=jpex,kpsx=kpsx,kpex=kpex )
2627 CALL uncouple_scalars_for_filter ( FIELD=moist(ims,kms,jms,im) &
2628 ,MU=grid%mu_2 , MUB=grid%mub &
2629 ,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
2630 ,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
2631 ,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe )
2635 IF ( num_3d_c >= PARAM_FIRST_SCALAR ) THEN
2636 CALL wrf_debug ( 200 , ' call filter chem ' )
2637 DO im = PARAM_FIRST_SCALAR, num_3d_c
2638 CALL couple_scalars_for_filter ( FIELD=chem(ims,kms,jms,im) &
2639 ,MU=grid%mu_2 , MUB=grid%mub &
2640 ,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
2641 ,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
2642 ,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe )
2643 CALL pxft ( grid=grid &
2656 ,positive_definite=.FALSE. &
2657 ,moist=moist,chem=chem,tracer=tracer,scalar=scalar &
2658 ,fft_filter_lat = config_flags%fft_filter_lat &
2660 ,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
2661 ,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
2662 ,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe &
2663 ,imsx=imsx,imex=imex,jmsx=jmsx,jmex=jmex,kmsx=kmsx,kmex=kmex &
2664 ,ipsx=ipsx,ipex=ipex,jpsx=jmsx,jpex=jpex,kpsx=kpsx,kpex=kpex )
2665 CALL uncouple_scalars_for_filter ( FIELD=chem(ims,kms,jms,im) &
2666 ,MU=grid%mu_2 , MUB=grid%mub &
2667 ,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
2668 ,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
2669 ,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe )
2672 IF ( num_tracer >= PARAM_FIRST_SCALAR ) THEN
2673 CALL wrf_debug ( 200 , ' call filter tracer ' )
2674 DO im = PARAM_FIRST_SCALAR, num_tracer
2675 CALL couple_scalars_for_filter ( FIELD=tracer(ims,kms,jms,im) &
2676 ,MU=grid%mu_2 , MUB=grid%mub &
2677 ,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
2678 ,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
2679 ,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe )
2680 CALL pxft ( grid=grid &
2693 ,positive_definite=.FALSE. &
2694 ,moist=moist,chem=chem,tracer=tracer,scalar=scalar &
2695 ,fft_filter_lat = config_flags%fft_filter_lat &
2697 ,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
2698 ,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
2699 ,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe &
2700 ,imsx=imsx,imex=imex,jmsx=jmsx,jmex=jmex,kmsx=kmsx,kmex=kmex &
2701 ,ipsx=ipsx,ipex=ipex,jpsx=jmsx,jpex=jpex,kpsx=kpsx,kpex=kpex )
2702 CALL uncouple_scalars_for_filter ( FIELD=tracer(ims,kms,jms,im) &
2703 ,MU=grid%mu_2 , MUB=grid%mub &
2704 ,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
2705 ,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
2706 ,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe )
2710 IF ( num_3d_s >= PARAM_FIRST_SCALAR ) THEN
2711 CALL wrf_debug ( 200 , ' call filter scalar ' )
2712 DO im = PARAM_FIRST_SCALAR, num_3d_s
2713 CALL couple_scalars_for_filter ( FIELD=scalar(ims,kms,jms,im) &
2714 ,MU=grid%mu_2 , MUB=grid%mub &
2715 ,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
2716 ,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
2717 ,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe )
2718 CALL pxft ( grid=grid &
2731 ,positive_definite=.FALSE. &
2732 ,moist=moist,chem=chem,tracer=tracer,scalar=scalar &
2733 ,fft_filter_lat = config_flags%fft_filter_lat &
2735 ,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
2736 ,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
2737 ,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe &
2738 ,imsx=imsx,imex=imex,jmsx=jmsx,jmex=jmex,kmsx=kmsx,kmex=kmex &
2739 ,ipsx=ipsx,ipex=ipex,jpsx=jmsx,jpex=jpex,kpsx=kpsx,kpex=kpex )
2740 CALL uncouple_scalars_for_filter ( FIELD=scalar(ims,kms,jms,im) &
2741 ,MU=grid%mu_2 , MUB=grid%mub &
2742 ,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
2743 ,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
2744 ,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe )
2747 END IF ! polar filter test
2749 !-----------------------------------------------------------
2750 ! END rk3 substep polar filter for scalars (moist,chem,scalar)
2751 !-----------------------------------------------------------
2753 !-----------------------------------------------------------
2754 ! Stencils for patch communications (WCS, 29 June 2001)
2756 ! here's where we need a wide comm stencil - these are the
2757 ! uncoupled variables so are used for high order calc in
2758 ! advection and mixong routines.
2762 ! * * * * * * * * * * * *
2763 ! * * * * * * * * * * * * *
2764 ! * + * * * + * * * * * + * * *
2765 ! * * * * * * * * * * * * *
2766 ! * * * * * * * * * * * *
2794 IF ( config_flags%h_mom_adv_order <= 4 ) THEN
2795 # include "HALO_EM_D2_3.inc"
2796 ELSE IF ( config_flags%h_mom_adv_order <= 6 ) THEN
2797 # include "HALO_EM_D2_5.inc"
2799 WRITE(wrf_err_message,*)'solve_em: invalid h_mom_adv_order = ',config_flags%h_mom_adv_order
2800 CALL wrf_error_fatal(TRIM(wrf_err_message))
2802 # include "PERIOD_BDY_EM_D.inc"
2803 # include "PERIOD_BDY_EM_MOIST2.inc"
2804 # include "PERIOD_BDY_EM_CHEM2.inc"
2805 # include "PERIOD_BDY_EM_TRACER2.inc"
2806 # include "PERIOD_BDY_EM_SCALAR2.inc"
2809 BENCH_START(bc_end_tim)
2811 !$OMP PRIVATE ( ij )
2812 tile_bc_loop_1: DO ij = 1 , grid%num_tiles
2813 CALL wrf_debug ( 200 , ' call rk_phys_bc_dry_2' )
2815 CALL rk_phys_bc_dry_2( config_flags, &
2816 grid%u_2, grid%v_2, grid%w_2, &
2817 grid%t_2, grid%ph_2, grid%mu_2, &
2818 ids, ide, jds, jde, kds, kde, &
2819 ims, ime, jms, jme, kms, kme, &
2820 ips, ipe, jps, jpe, kps, kpe, &
2821 grid%i_start(ij), grid%i_end(ij), &
2822 grid%j_start(ij), grid%j_end(ij), &
2825 BENCH_START(diag_w_tim)
2826 IF (.not. config_flags%non_hydrostatic) THEN
2827 CALL diagnose_w( ph_tend, grid%ph_2, grid%ph_1, grid%w_2, grid%muts, dt_rk, &
2828 grid%u_2, grid%v_2, grid%ht, &
2829 grid%cf1, grid%cf2, grid%cf3, grid%rdx, grid%rdy, grid%msftx, grid%msfty, &
2830 ids, ide, jds, jde, kds, kde, &
2831 ims, ime, jms, jme, kms, kme, &
2832 grid%i_start(ij), grid%i_end(ij), &
2833 grid%j_start(ij), grid%j_end(ij), &
2836 BENCH_END(diag_w_tim)
2838 IF (num_3d_m >= PARAM_FIRST_SCALAR) THEN
2840 moisture_loop_bdy_1 : DO im = PARAM_FIRST_SCALAR , num_3d_m
2842 CALL set_physical_bc3d( moist(ims,kms,jms,im), 'p', config_flags, &
2843 ids, ide, jds, jde, kds, kde, &
2844 ims, ime, jms, jme, kms, kme, &
2845 ips, ipe, jps, jpe, kps, kpe, &
2846 grid%i_start(ij), grid%i_end(ij), &
2847 grid%j_start(ij), grid%j_end(ij), &
2849 END DO moisture_loop_bdy_1
2853 IF (num_3d_c >= PARAM_FIRST_SCALAR) THEN
2855 chem_species_bdy_loop_1 : DO ic = PARAM_FIRST_SCALAR , num_3d_c
2857 CALL set_physical_bc3d( chem(ims,kms,jms,ic), 'p', config_flags, &
2858 ids, ide, jds, jde, kds, kde, &
2859 ims, ime, jms, jme, kms, kme, &
2860 ips, ipe, jps, jpe, kps, kpe, &
2861 grid%i_start(ij), grid%i_end(ij), &
2862 grid%j_start(ij), grid%j_end(ij), &
2865 END DO chem_species_bdy_loop_1
2869 IF (num_tracer >= PARAM_FIRST_SCALAR) THEN
2871 tracer_species_bdy_loop_1 : DO ic = PARAM_FIRST_SCALAR , num_tracer
2873 CALL set_physical_bc3d( tracer(ims,kms,jms,ic), 'p', config_flags, &
2874 ids, ide, jds, jde, kds, kde, &
2875 ims, ime, jms, jme, kms, kme, &
2876 ips, ipe, jps, jpe, kps, kpe, &
2877 grid%i_start(ij), grid%i_end(ij), &
2878 grid%j_start(ij), grid%j_end(ij), &
2881 END DO tracer_species_bdy_loop_1
2885 IF (num_3d_s >= PARAM_FIRST_SCALAR) THEN
2887 scalar_species_bdy_loop_1 : DO is = PARAM_FIRST_SCALAR , num_3d_s
2889 CALL set_physical_bc3d( scalar(ims,kms,jms,is), 'p', config_flags, &
2890 ids, ide, jds, jde, kds, kde, &
2891 ims, ime, jms, jme, kms, kme, &
2892 ips, ipe, jps, jpe, kps, kpe, &
2893 grid%i_start(ij), grid%i_end(ij), &
2894 grid%j_start(ij), grid%j_end(ij), &
2897 END DO scalar_species_bdy_loop_1
2901 IF (config_flags%km_opt .eq. 2) THEN
2903 CALL set_physical_bc3d( grid%tke_2 , 'p', config_flags, &
2904 ids, ide, jds, jde, kds, kde, &
2905 ims, ime, jms, jme, kms, kme, &
2906 ips, ipe, jps, jpe, kps, kpe, &
2907 grid%i_start(ij), grid%i_end(ij), &
2908 grid%j_start(ij), grid%j_end(ij), &
2912 END DO tile_bc_loop_1
2913 !$OMP END PARALLEL DO
2914 BENCH_END(bc_end_tim)
2921 ! * + * * + * * * + * *
2925 ! moist, chem, scalar, tke x
2928 IF ( config_flags%h_mom_adv_order <= 4 ) THEN
2929 IF ( (config_flags%tke_adv_opt /= ORIGINAL) .and. (rk_step == rk_order-1) ) THEN
2930 # include "HALO_EM_TKE_5.inc"
2932 # include "HALO_EM_TKE_3.inc"
2934 ELSE IF ( config_flags%h_mom_adv_order <= 6 ) THEN
2935 IF ( (config_flags%tke_adv_opt /= ORIGINAL) .and. (rk_step == rk_order-1) ) THEN
2936 # include "HALO_EM_TKE_7.inc"
2938 # include "HALO_EM_TKE_5.inc"
2941 WRITE(wrf_err_message,*)'solve_em: invalid h_sca_adv_order = ',config_flags%h_sca_adv_order
2942 CALL wrf_error_fatal(TRIM(wrf_err_message))
2945 IF ( num_moist .GE. PARAM_FIRST_SCALAR ) THEN
2946 IF ( config_flags%h_sca_adv_order <= 4 ) THEN
2947 IF ( (config_flags%moist_adv_opt /= ORIGINAL) .and. (rk_step == rk_order-1) ) THEN
2948 # include "HALO_EM_MOIST_E_5.inc"
2950 # include "HALO_EM_MOIST_E_3.inc"
2952 ELSE IF ( config_flags%h_sca_adv_order <= 6 ) THEN
2953 IF ( (config_flags%moist_adv_opt /= ORIGINAL) .and. (rk_step == rk_order-1) ) THEN
2954 # include "HALO_EM_MOIST_E_7.inc"
2956 # include "HALO_EM_MOIST_E_5.inc"
2959 WRITE(wrf_err_message,*)'solve_em: invalid h_sca_adv_order = ',config_flags%h_sca_adv_order
2960 CALL wrf_error_fatal(TRIM(wrf_err_message))
2963 IF ( num_chem >= PARAM_FIRST_SCALAR ) THEN
2964 IF ( config_flags%h_sca_adv_order <= 4 ) THEN
2965 IF ( (config_flags%chem_adv_opt /= ORIGINAL) .and. (rk_step == rk_order-1) ) THEN
2966 # include "HALO_EM_CHEM_E_5.inc"
2968 # include "HALO_EM_CHEM_E_3.inc"
2970 ELSE IF ( config_flags%h_sca_adv_order <= 6 ) THEN
2971 IF ( (config_flags%chem_adv_opt /= ORIGINAL) .and. (rk_step == rk_order-1) ) THEN
2972 # include "HALO_EM_CHEM_E_7.inc"
2974 # include "HALO_EM_CHEM_E_5.inc"
2977 WRITE(wrf_err_message,*)'solve_em: invalid h_sca_adv_order = ',config_flags%h_sca_adv_order
2978 CALL wrf_error_fatal(TRIM(wrf_err_message))
2981 IF ( num_tracer >= PARAM_FIRST_SCALAR ) THEN
2982 IF ( config_flags%h_sca_adv_order <= 4 ) THEN
2983 IF ( (config_flags%tracer_adv_opt /= ORIGINAL) .and. (rk_step == rk_order-1) ) THEN
2984 # include "HALO_EM_TRACER_E_5.inc"
2986 # include "HALO_EM_TRACER_E_3.inc"
2988 ELSE IF ( config_flags%h_sca_adv_order <= 6 ) THEN
2989 IF ( (config_flags%tracer_adv_opt /= ORIGINAL) .and. (rk_step == rk_order-1) ) THEN
2990 # include "HALO_EM_TRACER_E_7.inc"
2992 # include "HALO_EM_TRACER_E_5.inc"
2995 WRITE(wrf_err_message,*)'solve_em: invalid h_sca_adv_order = ',config_flags%h_sca_adv_order
2996 CALL wrf_error_fatal(TRIM(wrf_err_message))
2999 IF ( num_scalar >= PARAM_FIRST_SCALAR ) THEN
3000 IF ( config_flags%h_sca_adv_order <= 4 ) THEN
3001 IF ( (config_flags%scalar_adv_opt /= ORIGINAL) .and. (rk_step == rk_order-1) ) THEN
3002 # include "HALO_EM_SCALAR_E_5.inc"
3004 # include "HALO_EM_SCALAR_E_3.inc"
3006 ELSE IF ( config_flags%h_sca_adv_order <= 6 ) THEN
3007 IF ( (config_flags%scalar_adv_opt /= ORIGINAL) .and. (rk_step == rk_order-1) ) THEN
3008 # include "HALO_EM_SCALAR_E_7.inc"
3010 # include "HALO_EM_SCALAR_E_5.inc"
3013 WRITE(wrf_err_message,*)'solve_em: invalid h_sca_adv_order = ',config_flags%h_sca_adv_order
3014 CALL wrf_error_fatal(TRIM(wrf_err_message))
3019 ENDIF rk_step_1_check
3022 !**********************************************************
3024 ! end of RK predictor-corrector loop
3026 !**********************************************************
3028 END DO Runge_Kutta_loop
3030 IF (config_flags%do_avgflx_em .EQ. 1) THEN
3031 ! Reinitialize time-averaged fluxes if history output was written after the previous time step:
3032 CALL WRFU_ALARMGET(grid%alarms( HISTORY_ALARM ),prevringtime=temp_time)
3033 CALL domain_clock_get ( grid, current_time=CurrTime, &
3034 current_timestr=message2 )
3035 ! use overloaded -, .LT. operator to check whether to initialize avgflx:
3036 ! reinitialize after each history output (detect this here by comparing current time
3037 ! against last history time and time step - this code follows what's done in adapt_timestep_em):
3038 WRITE ( message , FMT = '("solve_em: old_dt =",g15.6,", dt=",g15.6," on domain ",I3)' ) &
3039 & old_dt,grid%dt,grid%id
3040 CALL wrf_debug(200,message)
3041 old_dt=min(old_dt,grid%dt)
3042 num = INT(old_dt * precision)
3044 CALL WRFU_TimeIntervalSet(dtInterval, Sn=num, Sd=den)
3045 IF (CurrTime .lt. temp_time + dtInterval) THEN
3046 WRITE ( message , FMT = '("solve_em: initializing avgflx at time ",A," on domain ",I3)' ) &
3047 & TRIM(message2), grid%id
3048 CALL wrf_message(trim(message))
3049 grid%avgflx_count = 0
3050 !tile-loop for zero_avgflx
3052 !$OMP PRIVATE ( ij )
3054 DO ij = 1 , grid%num_tiles
3055 CALL wrf_debug(200,'In solve_em, before zero_avgflx call')
3056 CALL zero_avgflx(grid%avgflx_rum,grid%avgflx_rvm,grid%avgflx_wwm, &
3057 & ids, ide, jds, jde, kds, kde, &
3058 & ims, ime, jms, jme, kms, kme, &
3059 & grid%i_start(ij), grid%i_end(ij), grid%j_start(ij), grid%j_end(ij), &
3060 & k_start , k_end, f_flux, &
3061 & grid%avgflx_cfu1,grid%avgflx_cfd1,grid%avgflx_dfu1, &
3062 & grid%avgflx_efu1,grid%avgflx_dfd1,grid%avgflx_efd1 )
3063 CALL wrf_debug(200,'In solve_em, after zero_avgflx call')
3067 ! Update avgflx quantities
3068 !tile-loop for upd_avgflx
3070 !$OMP PRIVATE ( ij )
3072 DO ij = 1 , grid%num_tiles
3073 CALL wrf_debug(200,'In solve_em, before upd_avgflx call')
3074 CALL upd_avgflx(grid%avgflx_count,grid%avgflx_rum,grid%avgflx_rvm,grid%avgflx_wwm, &
3075 & grid%ru_m, grid%rv_m, grid%ww_m, &
3076 & ids, ide, jds, jde, kds, kde, &
3077 & ims, ime, jms, jme, kms, kme, &
3078 & grid%i_start(ij), grid%i_end(ij), grid%j_start(ij), grid%j_end(ij), &
3079 & k_start , k_end, f_flux, &
3080 & grid%cfu1,grid%cfd1,grid%dfu1,grid%efu1,grid%dfd1,grid%efd1, &
3081 & grid%avgflx_cfu1,grid%avgflx_cfd1,grid%avgflx_dfu1, &
3082 & grid%avgflx_efu1,grid%avgflx_dfd1,grid%avgflx_efd1 )
3083 CALL wrf_debug(200,'In solve_em, after upd_avgflx call')
3086 grid%avgflx_count = grid%avgflx_count + 1
3090 !$OMP PRIVATE ( ij )
3091 DO ij = 1 , grid%num_tiles
3093 BENCH_START(advance_ppt_tim)
3094 CALL wrf_debug ( 200 , ' call advance_ppt' )
3095 CALL advance_ppt(grid%rthcuten,grid%rqvcuten,grid%rqccuten,grid%rqrcuten, &
3096 grid%rqicuten,grid%rqscuten, &
3097 grid%rainc,grid%raincv,grid%rainsh,grid%pratec,grid%pratesh, &
3098 grid%nca,grid%htop,grid%hbot,grid%cutop,grid%cubot, &
3099 grid%cuppt, grid%dt, config_flags, &
3100 ids,ide, jds,jde, kds,kde, &
3101 ims,ime, jms,jme, kms,kme, &
3102 grid%i_start(ij), grid%i_end(ij), &
3103 grid%j_start(ij), grid%j_end(ij), &
3105 BENCH_END(advance_ppt_tim)
3108 !$OMP END PARALLEL DO
3112 ! (5) time-split physics.
3114 ! Microphysics are the only time split physics in the WRF model
3115 ! at this time. Split-physics begins with the calculation of
3116 ! needed diagnostic quantities (pressure, temperature, etc.)
3117 ! followed by a call to the microphysics driver,
3118 ! and finishes with a clean-up, storing off of a diabatic tendency
3119 ! from the moist physics, and a re-calulation of the diagnostic
3120 ! quantities pressure and density.
3124 IF( config_flags%specified .or. config_flags%nested ) THEN
3130 IF (config_flags%mp_physics /= 0) then
3133 !$OMP PRIVATE ( ij, its, ite, jts, jte )
3135 scalar_tile_loop_1a: DO ij = 1 , grid%num_tiles
3137 IF ( config_flags%periodic_x ) THEN
3138 its = max(grid%i_start(ij),ids)
3139 ite = min(grid%i_end(ij),ide-1)
3141 its = max(grid%i_start(ij),ids+sz)
3142 ite = min(grid%i_end(ij),ide-1-sz)
3144 jts = max(grid%j_start(ij),jds+sz)
3145 jte = min(grid%j_end(ij),jde-1-sz)
3147 CALL wrf_debug ( 200 , ' call moist_physics_prep' )
3148 BENCH_START(moist_physics_prep_tim)
3149 CALL moist_physics_prep_em( grid%t_2, grid%t_1, t0, rho, &
3150 grid%al, grid%alb, grid%p, p8w, p0, grid%pb, &
3151 grid%ph_2, grid%phb, th_phy, pi_phy, p_phy, &
3152 grid%z, grid%z_at_w, dz8w, &
3153 dtm, grid%h_diabatic, &
3154 config_flags,grid%fnm, grid%fnp, &
3155 ids, ide, jds, jde, kds, kde, &
3156 ims, ime, jms, jme, kms, kme, &
3157 its, ite, jts, jte, &
3159 BENCH_END(moist_physics_prep_tim)
3160 END DO scalar_tile_loop_1a
3161 !$OMP END PARALLEL DO
3163 CALL wrf_debug ( 200 , ' call microphysics_driver' )
3166 specified_bdy = config_flags%specified .OR. config_flags%nested
3167 channel_bdy = config_flags%specified .AND. config_flags%periodic_x
3169 BENCH_START(micro_driver_tim)
3171 CALL microphysics_driver( &
3172 & DT=dtm ,DX=grid%dx ,DY=grid%dy &
3173 & ,DZ8W=dz8w ,F_ICE_PHY=grid%f_ice_phy &
3174 & ,ITIMESTEP=grid%itimestep ,LOWLYR=grid%lowlyr &
3175 & ,P8W=p8w ,P=p_phy ,PI_PHY=pi_phy &
3176 & ,RHO=rho ,SPEC_ZONE=grid%spec_zone &
3177 & ,SR=grid%sr ,TH=th_phy &
3178 & ,refl_10cm=grid%refl_10cm & ! hm, 9/22/09 for refl
3179 & ,WARM_RAIN=grid%warm_rain &
3181 & ,CLDFRA=grid%cldfra, EXCH_H=grid%exch_h &
3182 & ,NSOURCE=grid%qndropsource &
3184 & ,QLSINK=grid%qlsink,CLDFRA_OLD=grid%cldfra_old &
3185 & ,PRECR=grid%precr, PRECI=grid%preci, PRECS=grid%precs, PRECG=grid%precg &
3186 & ,CHEM_OPT=config_flags%chem_opt, PROGN=config_flags%progn &
3188 & ,XLAND=grid%xland &
3189 & ,SPECIFIED=specified_bdy, CHANNEL_SWITCH=channel_bdy &
3190 & ,F_RAIN_PHY=grid%f_rain_phy &
3191 & ,F_RIMEF_PHY=grid%f_rimef_phy &
3192 & ,MP_PHYSICS=config_flags%mp_physics &
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 &
3197 & ,IPS=ips,IPE=ipe, JPS=jps,JPE=jpe, KPS=kps,KPE=kpe &
3199 & ,I_START=grid%i_start,I_END=min(grid%i_end, ide-1) &
3200 & ,J_START=grid%j_start,J_END=min(grid%j_end, jde-1) &
3201 & ,KTS=k_start, KTE=min(k_end,kde-1) &
3202 & ,NUM_TILES=grid%num_tiles &
3205 & , RAINNC=grid%rainnc, RAINNCV=grid%rainncv &
3206 & , SNOWNC=grid%snownc, SNOWNCV=grid%snowncv &
3207 & , GRAUPELNC=grid%graupelnc, GRAUPELNCV=grid%graupelncv & ! for milbrandt2mom
3208 & , HAILNC=grid%hailnc, HAILNCV=grid%hailncv &
3209 & , W=grid%w_2, Z=grid%z, HT=grid%ht &
3210 & , MP_RESTART_STATE=grid%mp_restart_state &
3211 & , TBPVS_STATE=grid%tbpvs_state & ! etampnew
3212 & , TBPVS0_STATE=grid%tbpvs0_state & ! etampnew
3213 & , QV_CURR=moist(ims,kms,jms,P_QV), F_QV=F_QV &
3214 & , QC_CURR=moist(ims,kms,jms,P_QC), F_QC=F_QC &
3215 & , QR_CURR=moist(ims,kms,jms,P_QR), F_QR=F_QR &
3216 & , QI_CURR=moist(ims,kms,jms,P_QI), F_QI=F_QI &
3217 & , QS_CURR=moist(ims,kms,jms,P_QS), F_QS=F_QS &
3218 & , QG_CURR=moist(ims,kms,jms,P_QG), F_QG=F_QG &
3219 & , QH_CURR=moist(ims,kms,jms,P_QH), F_QH=F_QH & ! for milbrandt2mom
3220 & , QNDROP_CURR=scalar(ims,kms,jms,P_QNDROP), F_QNDROP=F_QNDROP &
3221 & , QT_CURR=scalar(ims,kms,jms,P_QT), F_QT=F_QT &
3222 & , QNN_CURR=scalar(ims,kms,jms,P_QNN), F_QNN=F_QNN &
3223 & , QNI_CURR=scalar(ims,kms,jms,P_QNI), F_QNI=F_QNI &
3224 & , QNC_CURR=scalar(ims,kms,jms,P_QNC), F_QNC=F_QNC &
3225 & , QNR_CURR=scalar(ims,kms,jms,P_QNR), F_QNR=F_QNR &
3226 & , QNS_CURR=scalar(ims,kms,jms,P_QNS), F_QNS=F_QNS &
3227 & , QNG_CURR=scalar(ims,kms,jms,P_QNG), F_QNG=F_QNG &
3228 & , QNH_CURR=scalar(ims,kms,jms,P_QNH), F_QNH=F_QNH & ! for milbrandt2mom
3229 ! & , QZR_CURR=scalar(ims,kms,jms,P_QZR), F_QZR=F_QZR & ! for milbrandt3mom
3230 ! & , QZI_CURR=scalar(ims,kms,jms,P_QZI), F_QZI=F_QZI & ! "
3231 ! & , QZS_CURR=scalar(ims,kms,jms,P_QZS), F_QZS=F_QZS & ! "
3232 ! & , QZG_CURR=scalar(ims,kms,jms,P_QZG), F_QZG=F_QZG & ! "
3233 ! & , QZH_CURR=scalar(ims,kms,jms,P_QZH), F_QZH=F_QZH & ! "
3234 & , qrcuten=grid%rqrcuten, qscuten=grid%rqscuten &
3235 & , qicuten=grid%rqicuten,mu=grid%mut &
3236 & , HAIL=config_flags%gsfcgce_hail & ! for gsfcgce
3237 & , ICE2=config_flags%gsfcgce_2ice & ! for gsfcgce
3238 ! & , ccntype=config_flags%milbrandt_ccntype & ! for milbrandt (2mom)
3241 & , RI_CURR=grid%rimi &
3243 BENCH_END(micro_driver_tim)
3246 BENCH_START(microswap_2)
3247 ! for load balancing; communication to redistribute the points
3248 IF ( config_flags%mp_physics .EQ. ETAMPNEW ) THEN
3249 #include "SWAP_ETAMP_NEW.inc"
3250 ELSE IF ( config_flags%mp_physics .EQ. WSM3SCHEME ) THEN
3251 #include "SWAP_WSM3.inc"
3253 BENCH_END(microswap_2)
3256 CALL wrf_debug ( 200 , ' call moist_physics_finish' )
3257 BENCH_START(moist_phys_end_tim)
3260 !$OMP PRIVATE ( ij, its, ite, jts, jte, im, ii, jj, kk )
3262 DO ij = 1 , grid%num_tiles
3264 its = max(grid%i_start(ij),ids)
3265 ite = min(grid%i_end(ij),ide-1)
3266 jts = max(grid%j_start(ij),jds)
3267 jte = min(grid%j_end(ij),jde-1)
3269 CALL microphysics_zero_outb ( &
3270 moist , num_moist , config_flags , &
3271 ids, ide, jds, jde, kds, kde, &
3272 ims, ime, jms, jme, kms, kme, &
3273 its, ite, jts, jte, &
3276 CALL microphysics_zero_outb ( &
3277 scalar , num_scalar , config_flags , &
3278 ids, ide, jds, jde, kds, kde, &
3279 ims, ime, jms, jme, kms, kme, &
3280 its, ite, jts, jte, &
3283 CALL microphysics_zero_outb ( &
3284 chem , num_chem , config_flags , &
3285 ids, ide, jds, jde, kds, kde, &
3286 ims, ime, jms, jme, kms, kme, &
3287 its, ite, jts, jte, &
3289 CALL microphysics_zero_outb ( &
3290 tracer , num_tracer , config_flags , &
3291 ids, ide, jds, jde, kds, kde, &
3292 ims, ime, jms, jme, kms, kme, &
3293 its, ite, jts, jte, &
3296 IF ( config_flags%periodic_x ) THEN
3297 its = max(grid%i_start(ij),ids)
3298 ite = min(grid%i_end(ij),ide-1)
3300 its = max(grid%i_start(ij),ids+sz)
3301 ite = min(grid%i_end(ij),ide-1-sz)
3303 jts = max(grid%j_start(ij),jds+sz)
3304 jte = min(grid%j_end(ij),jde-1-sz)
3306 CALL microphysics_zero_outa ( &
3307 moist , num_moist , config_flags , &
3308 ids, ide, jds, jde, kds, kde, &
3309 ims, ime, jms, jme, kms, kme, &
3310 its, ite, jts, jte, &
3313 CALL microphysics_zero_outa ( &
3314 scalar , num_scalar , config_flags , &
3315 ids, ide, jds, jde, kds, kde, &
3316 ims, ime, jms, jme, kms, kme, &
3317 its, ite, jts, jte, &
3320 CALL microphysics_zero_outa ( &
3321 chem , num_chem , config_flags , &
3322 ids, ide, jds, jde, kds, kde, &
3323 ims, ime, jms, jme, kms, kme, &
3324 its, ite, jts, jte, &
3327 CALL microphysics_zero_outa ( &
3328 tracer , num_tracer , config_flags , &
3329 ids, ide, jds, jde, kds, kde, &
3330 ims, ime, jms, jme, kms, kme, &
3331 its, ite, jts, jte, &
3334 CALL moist_physics_finish_em( grid%t_2, grid%t_1, t0, grid%muts, th_phy, &
3335 grid%h_diabatic, dtm, config_flags, &
3336 #if ( WRF_DFI_RADAR == 1 )
3337 grid%dfi_tten_rad,grid%dfi_stage, &
3339 ids, ide, jds, jde, kds, kde, &
3340 ims, ime, jms, jme, kms, kme, &
3341 its, ite, jts, jte, &
3345 !$OMP END PARALLEL DO
3347 ENDIF ! microphysics test
3349 !-----------------------------------------------------------
3350 ! filter for moist variables post-microphysics and end of timestep
3351 !-----------------------------------------------------------
3353 IF (config_flags%polar) THEN
3354 IF ( num_3d_m >= PARAM_FIRST_SCALAR ) THEN
3355 CALL wrf_debug ( 200 , ' call filter moist' )
3356 DO im = PARAM_FIRST_SCALAR, num_3d_m
3357 DO jj = jps, MIN(jpe,jde-1)
3358 DO kk = kps, MIN(kpe,kde-1)
3359 DO ii = ips, MIN(ipe,ide-1)
3360 moist(ii,kk,jj,im)=moist(ii,kk,jj,im)*(grid%mu_2(ii,jj)+grid%mub(ii,jj))
3365 CALL pxft ( grid=grid &
3378 ,positive_definite=.FALSE. &
3379 ,moist=moist,chem=chem,tracer=tracer,scalar=scalar &
3380 ,fft_filter_lat = config_flags%fft_filter_lat &
3382 ,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
3383 ,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
3384 ,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe &
3385 ,imsx=imsx,imex=imex,jmsx=jmsx,jmex=jmex,kmsx=kmsx,kmex=kmex &
3386 ,ipsx=ipsx,ipex=ipex,jpsx=jmsx,jpex=jpex,kpsx=kpsx,kpex=kpex )
3388 DO jj = jps, MIN(jpe,jde-1)
3389 DO kk = kps, MIN(kpe,kde-1)
3390 DO ii = ips, MIN(ipe,ide-1)
3391 moist(ii,kk,jj,im)=moist(ii,kk,jj,im)/(grid%mu_2(ii,jj)+grid%mub(ii,jj))
3399 !-----------------------------------------------------------
3400 ! end filter for moist variables post-microphysics and end of timestep
3401 !-----------------------------------------------------------
3404 !$OMP PRIVATE ( ij, its, ite, jts, jte, im, ii, jj, kk )
3405 scalar_tile_loop_1ba: DO ij = 1 , grid%num_tiles
3407 IF ( config_flags%periodic_x ) THEN
3408 its = max(grid%i_start(ij),ids)
3409 ite = min(grid%i_end(ij),ide-1)
3411 its = max(grid%i_start(ij),ids+sz)
3412 ite = min(grid%i_end(ij),ide-1-sz)
3414 jts = max(grid%j_start(ij),jds+sz)
3415 jte = min(grid%j_end(ij),jde-1-sz)
3417 CALL calc_p_rho_phi( moist, num_3d_m, config_flags%hypsometric_opt, &
3418 grid%al, grid%alb, grid%mu_2, grid%muts, &
3419 grid%ph_2, grid%phb, grid%p, grid%pb, grid%t_2, &
3420 p0, t0, grid%p_top, grid%znu, grid%znw, grid%dnw, grid%rdnw, &
3421 grid%rdn, config_flags%non_hydrostatic, &
3422 ids, ide, jds, jde, kds, kde, &
3423 ims, ime, jms, jme, kms, kme, &
3424 its, ite, jts, jte, &
3427 END DO scalar_tile_loop_1ba
3428 !$OMP END PARALLEL DO
3429 BENCH_END(moist_phys_end_tim)
3431 IF (.not. config_flags%non_hydrostatic) THEN
3433 # include "HALO_EM_HYDRO_UV.inc"
3434 # include "PERIOD_EM_HYDRO_UV.inc"
3437 !$OMP PRIVATE ( ij )
3438 DO ij = 1 , grid%num_tiles
3439 CALL diagnose_w( ph_tend, grid%ph_2, grid%ph_1, grid%w_2, grid%muts, dt_rk, &
3440 grid%u_2, grid%v_2, grid%ht, &
3441 grid%cf1, grid%cf2, grid%cf3, grid%rdx, grid%rdy, grid%msftx, grid%msfty, &
3442 ids, ide, jds, jde, kds, kde, &
3443 ims, ime, jms, jme, kms, kme, &
3444 grid%i_start(ij), grid%i_end(ij), &
3445 grid%j_start(ij), grid%j_end(ij), &
3449 !$OMP END PARALLEL DO
3453 CALL wrf_debug ( 200 , ' call chem polar filter ' )
3455 !-----------------------------------------------------------
3456 ! filter for chem and scalar variables at end of timestep
3457 !-----------------------------------------------------------
3459 IF (config_flags%polar) THEN
3461 IF ( num_3d_c >= PARAM_FIRST_SCALAR ) then
3462 chem_filter_loop: DO im = PARAM_FIRST_SCALAR, num_3d_c
3463 DO jj = jps, MIN(jpe,jde-1)
3464 DO kk = kps, MIN(kpe,kde-1)
3465 DO ii = ips, MIN(ipe,ide-1)
3466 chem(ii,kk,jj,im)=chem(ii,kk,jj,im)*(grid%mu_2(ii,jj)+grid%mub(ii,jj))
3471 CALL pxft ( grid=grid &
3484 ,positive_definite=.FALSE. &
3485 ,moist=moist,chem=chem,tracer=tracer,scalar=scalar &
3486 ,fft_filter_lat = config_flags%fft_filter_lat &
3488 ,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
3489 ,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
3490 ,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe &
3491 ,imsx=imsx,imex=imex,jmsx=jmsx,jmex=jmex,kmsx=kmsx,kmex=kmex &
3492 ,ipsx=ipsx,ipex=ipex,jpsx=jmsx,jpex=jpex,kpsx=kpsx,kpex=kpex )
3494 DO jj = jps, MIN(jpe,jde-1)
3495 DO kk = kps, MIN(kpe,kde-1)
3496 DO ii = ips, MIN(ipe,ide-1)
3497 chem(ii,kk,jj,im)=chem(ii,kk,jj,im)/(grid%mu_2(ii,jj)+grid%mub(ii,jj))
3501 ENDDO chem_filter_loop
3503 IF ( num_tracer >= PARAM_FIRST_SCALAR ) then
3504 tracer_filter_loop: DO im = PARAM_FIRST_SCALAR, num_tracer
3505 DO jj = jps, MIN(jpe,jde-1)
3506 DO kk = kps, MIN(kpe,kde-1)
3507 DO ii = ips, MIN(ipe,ide-1)
3508 tracer(ii,kk,jj,im)=tracer(ii,kk,jj,im)*(grid%mu_2(ii,jj)+grid%mub(ii,jj))
3513 CALL pxft ( grid=grid &
3526 ,positive_definite=.FALSE. &
3527 ,moist=moist,chem=chem,tracer=tracer,scalar=scalar &
3528 ,fft_filter_lat = config_flags%fft_filter_lat &
3530 ,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
3531 ,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
3532 ,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe &
3533 ,imsx=imsx,imex=imex,jmsx=jmsx,jmex=jmex,kmsx=kmsx,kmex=kmex &
3534 ,ipsx=ipsx,ipex=ipex,jpsx=jmsx,jpex=jpex,kpsx=kpsx,kpex=kpex )
3536 DO jj = jps, MIN(jpe,jde-1)
3537 DO kk = kps, MIN(kpe,kde-1)
3538 DO ii = ips, MIN(ipe,ide-1)
3539 tracer(ii,kk,jj,im)=tracer(ii,kk,jj,im)/(grid%mu_2(ii,jj)+grid%mub(ii,jj))
3543 ENDDO tracer_filter_loop
3546 IF ( num_3d_s >= PARAM_FIRST_SCALAR ) then
3547 scalar_filter_loop: DO im = PARAM_FIRST_SCALAR, num_3d_s
3548 DO jj = jps, MIN(jpe,jde-1)
3549 DO kk = kps, MIN(kpe,kde-1)
3550 DO ii = ips, MIN(ipe,ide-1)
3551 scalar(ii,kk,jj,im)=scalar(ii,kk,jj,im)*(grid%mu_2(ii,jj)+grid%mub(ii,jj))
3556 CALL pxft ( grid=grid &
3569 ,positive_definite=.FALSE. &
3570 ,moist=moist,chem=chem,tracer=tracer,scalar=scalar &
3571 ,fft_filter_lat = config_flags%fft_filter_lat &
3573 ,ids=ids,ide=ide,jds=jds,jde=jde,kds=kds,kde=kde &
3574 ,ims=ims,ime=ime,jms=jms,jme=jme,kms=kms,kme=kme &
3575 ,ips=ips,ipe=ipe,jps=jps,jpe=jpe,kps=kps,kpe=kpe &
3576 ,imsx=imsx,imex=imex,jmsx=jmsx,jmex=jmex,kmsx=kmsx,kmex=kmex &
3577 ,ipsx=ipsx,ipex=ipex,jpsx=jmsx,jpex=jpex,kpsx=kpsx,kpex=kpex )
3579 DO jj = jps, MIN(jpe,jde-1)
3580 DO kk = kps, MIN(kpe,kde-1)
3581 DO ii = ips, MIN(ipe,ide-1)
3582 scalar(ii,kk,jj,im)=scalar(ii,kk,jj,im)/(grid%mu_2(ii,jj)+grid%mub(ii,jj))
3586 ENDDO scalar_filter_loop
3590 !-----------------------------------------------------------
3591 ! end filter for chem and scalar variables at end of timestep
3592 !-----------------------------------------------------------
3594 ! We're finished except for boundary condition (and patch) update
3596 ! Boundary condition time (or communication time). At this time, we have
3597 ! implemented periodic and symmetric physical boundary conditions.
3599 ! b.c. routine for data within patch.
3601 ! we need to do both time levels of
3602 ! data because the time filter only works in the physical solution space.
3604 ! First, do patch communications for boundary conditions (periodicity)
3606 !-----------------------------------------------------------
3607 ! Stencils for patch communications (WCS, 29 June 2001)
3609 ! here's where we need a wide comm stencil - these are the
3610 ! uncoupled variables so are used for high order calc in
3611 ! advection and mixong routines.
3615 ! * + * * + * * * + * *
3640 !----------------------------------------------------------
3644 IF ( config_flags%h_mom_adv_order <= 4 ) THEN
3645 # include "HALO_EM_D3_3.inc"
3646 ELSE IF ( config_flags%h_mom_adv_order <= 6 ) THEN
3647 # include "HALO_EM_D3_5.inc"
3649 WRITE(wrf_err_message,*)'solve_em: invalid h_mom_adv_order = ',config_flags%h_mom_adv_order
3650 CALL wrf_error_fatal(TRIM(wrf_err_message))
3652 # include "PERIOD_BDY_EM_D3.inc"
3653 # include "PERIOD_BDY_EM_MOIST.inc"
3654 # include "PERIOD_BDY_EM_CHEM.inc"
3655 # include "PERIOD_BDY_EM_TRACER.inc"
3656 # include "PERIOD_BDY_EM_SCALAR.inc"
3659 ! now set physical b.c on a patch
3661 BENCH_START(bc_2d_tim)
3663 !$OMP PRIVATE ( ij )
3664 tile_bc_loop_2: DO ij = 1 , grid%num_tiles
3666 CALL wrf_debug ( 200 , ' call set_phys_bc_dry_2' )
3668 CALL set_phys_bc_dry_2( config_flags, &
3669 grid%u_1, grid%u_2, grid%v_1, grid%v_2, grid%w_1, grid%w_2, &
3670 grid%t_1, grid%t_2, grid%ph_1, grid%ph_2, grid%mu_1, grid%mu_2, &
3671 ids, ide, jds, jde, kds, kde, &
3672 ims, ime, jms, jme, kms, kme, &
3673 ips, ipe, jps, jpe, kps, kpe, &
3674 grid%i_start(ij), grid%i_end(ij), &
3675 grid%j_start(ij), grid%j_end(ij), &
3678 CALL set_physical_bc3d( grid%tke_1, 'p', config_flags, &
3679 ids, ide, jds, jde, kds, kde, &
3680 ims, ime, jms, jme, kms, kme, &
3681 ips, ipe, jps, jpe, kps, kpe, &
3682 grid%i_start(ij), grid%i_end(ij), &
3683 grid%j_start(ij), grid%j_end(ij), &
3686 CALL set_physical_bc3d( grid%tke_2 , 'p', config_flags, &
3687 ids, ide, jds, jde, kds, kde, &
3688 ims, ime, jms, jme, kms, kme, &
3689 ips, ipe, jps, jpe, kps, kpe, &
3690 grid%i_start(ij), grid%i_end(ij), &
3691 grid%j_start(ij), grid%j_end(ij), &
3694 moisture_loop_bdy_2 : DO im = PARAM_FIRST_SCALAR , num_3d_m
3696 CALL set_physical_bc3d( moist(ims,kms,jms,im), 'p', &
3698 ids, ide, jds, jde, kds, kde, &
3699 ims, ime, jms, jme, kms, kme, &
3700 ips, ipe, jps, jpe, kps, kpe, &
3701 grid%i_start(ij), grid%i_end(ij), &
3702 grid%j_start(ij), grid%j_end(ij), &
3705 END DO moisture_loop_bdy_2
3707 chem_species_bdy_loop_2 : DO ic = PARAM_FIRST_SCALAR , num_3d_c
3709 CALL set_physical_bc3d( chem(ims,kms,jms,ic) , 'p', config_flags, &
3710 ids, ide, jds, jde, kds, kde, &
3711 ims, ime, jms, jme, kms, kme, &
3712 ips, ipe, jps, jpe, kps, kpe, &
3713 grid%i_start(ij), grid%i_end(ij), &
3714 grid%j_start(ij), grid%j_end(ij), &
3717 END DO chem_species_bdy_loop_2
3719 tracer_species_bdy_loop_2 : DO ic = PARAM_FIRST_SCALAR , num_tracer
3721 CALL set_physical_bc3d( tracer(ims,kms,jms,ic) , 'p', config_flags, &
3722 ids, ide, jds, jde, kds, kde, &
3723 ims, ime, jms, jme, kms, kme, &
3724 ips, ipe, jps, jpe, kps, kpe, &
3725 grid%i_start(ij), grid%i_end(ij), &
3726 grid%j_start(ij), grid%j_end(ij), &
3729 END DO tracer_species_bdy_loop_2
3731 scalar_species_bdy_loop_2 : DO is = PARAM_FIRST_SCALAR , num_3d_s
3733 CALL set_physical_bc3d( scalar(ims,kms,jms,is) , 'p', config_flags, &
3734 ids, ide, jds, jde, kds, kde, &
3735 ims, ime, jms, jme, kms, kme, &
3736 ips, ipe, jps, jpe, kps, kpe, &
3737 grid%i_start(ij), grid%i_end(ij), &
3738 grid%j_start(ij), grid%j_end(ij), &
3741 END DO scalar_species_bdy_loop_2
3743 END DO tile_bc_loop_2
3744 !$OMP END PARALLEL DO
3745 BENCH_END(bc_2d_tim)
3747 IF( config_flags%specified .or. config_flags%nested ) THEN
3748 grid%dtbc = grid%dtbc + grid%dt
3751 ! reset surface w for consistency
3754 # include "HALO_EM_C.inc"
3755 # include "PERIOD_BDY_EM_E.inc"
3758 CALL wrf_debug ( 10 , ' call set_w_surface' )
3759 fill_w_flag = .false.
3762 !$OMP PRIVATE ( ij )
3763 DO ij = 1 , grid%num_tiles
3764 CALL set_w_surface( config_flags, grid%znw, fill_w_flag, &
3765 grid%w_2, grid%ht, grid%u_2, grid%v_2, &
3766 grid%cf1, grid%cf2, grid%cf3, grid%rdx, grid%rdy,&
3767 grid%msftx, grid%msfty, &
3768 ids, ide, jds, jde, kds, kde, &
3769 ims, ime, jms, jme, kms, kme, &
3770 grid%i_start(ij), grid%i_end(ij), &
3771 grid%j_start(ij), grid%j_end(ij), &
3773 ! its, ite, jts, jte, k_start, min(k_end,kde-1), &
3776 !$OMP END PARALLEL DO
3778 ! calculate some model diagnostics.
3780 CALL wrf_debug ( 200 , ' call diagnostic_driver' )
3782 CALL diagnostic_output_calc( &
3783 & DPSDT=grid%dpsdt ,DMUDT=grid%dmudt &
3784 & ,P8W=p8w ,PK1M=grid%pk1m &
3785 & ,MU_2=grid%mu_2 ,MU_2M=grid%mu_2m &
3786 & ,U=grid%u_2 ,V=grid%v_2 &
3787 & ,RAINCV=grid%raincv ,RAINNCV=grid%rainncv &
3788 & ,RAINC=grid%rainc ,RAINNC=grid%rainnc &
3789 & ,I_RAINC=grid%i_rainc ,I_RAINNC=grid%i_rainnc &
3790 & ,HFX=grid%hfx ,SFCEVP=grid%sfcevp ,LH=grid%lh &
3791 & ,DT=grid%dt ,SBW=config_flags%spec_bdy_width &
3792 & ,XTIME=grid%xtime ,T2=grid%t2 &
3793 & ,ACSWUPT=grid%acswupt ,ACSWUPTC=grid%acswuptc &
3794 & ,ACSWDNT=grid%acswdnt ,ACSWDNTC=grid%acswdntc &
3795 & ,ACSWUPB=grid%acswupb ,ACSWUPBC=grid%acswupbc &
3796 & ,ACSWDNB=grid%acswdnb ,ACSWDNBC=grid%acswdnbc &
3797 & ,ACLWUPT=grid%aclwupt ,ACLWUPTC=grid%aclwuptc &
3798 & ,ACLWDNT=grid%aclwdnt ,ACLWDNTC=grid%aclwdntc &
3799 & ,ACLWUPB=grid%aclwupb ,ACLWUPBC=grid%aclwupbc &
3800 & ,ACLWDNB=grid%aclwdnb ,ACLWDNBC=grid%aclwdnbc &
3801 & ,I_ACSWUPT=grid%i_acswupt ,I_ACSWUPTC=grid%i_acswuptc &
3802 & ,I_ACSWDNT=grid%i_acswdnt ,I_ACSWDNTC=grid%i_acswdntc &
3803 & ,I_ACSWUPB=grid%i_acswupb ,I_ACSWUPBC=grid%i_acswupbc &
3804 & ,I_ACSWDNB=grid%i_acswdnb ,I_ACSWDNBC=grid%i_acswdnbc &
3805 & ,I_ACLWUPT=grid%i_aclwupt ,I_ACLWUPTC=grid%i_aclwuptc &
3806 & ,I_ACLWDNT=grid%i_aclwdnt ,I_ACLWDNTC=grid%i_aclwdntc &
3807 & ,I_ACLWUPB=grid%i_aclwupb ,I_ACLWUPBC=grid%i_aclwupbc &
3808 & ,I_ACLWDNB=grid%i_aclwdnb ,I_ACLWDNBC=grid%i_aclwdnbc &
3810 & ,DIAG_PRINT=config_flags%diag_print &
3811 & ,BUCKET_MM=config_flags%bucket_mm &
3812 & ,BUCKET_J =config_flags%bucket_J &
3813 & ,SNOWNCV=grid%snowncv, SNOW_ACC_NC=grid%snow_acc_nc &
3814 & ,PREC_ACC_C=grid%prec_acc_c &
3815 & ,PREC_ACC_NC=grid%prec_acc_nc &
3816 & ,PREC_ACC_DT=config_flags%prec_acc_dt &
3817 & ,CURR_SECS=curr_secs &
3818 ! Dimension arguments
3819 & ,IDS=ids,IDE=ide, JDS=jds,JDE=jde, KDS=kds,KDE=kde &
3820 & ,IMS=ims,IME=ime, JMS=jms,JME=jme, KMS=kms,KME=kme &
3821 & ,IPS=ips,IPE=ipe, JPS=jps,JPE=jpe, KPS=kps,KPE=kpe &
3822 & ,I_START=grid%i_start,I_END=min(grid%i_end, ide-1) &
3823 & ,J_START=grid%j_start,J_END=min(grid%j_end, jde-1) &
3824 & ,KTS=k_start, KTE=min(k_end,kde-1) &
3825 & ,NUM_TILES=grid%num_tiles &
3828 IF (config_flags%output_diagnostics == 1) THEN
3829 IF ((config_flags%auxhist3_interval == 0 ) ) THEN
3830 WRITE (wrf_err_message , * )"CLWRF: ERROR -- error -- ERROR -- error : NO 'auxhist3_interval' has been defined in 'namelist.input'"
3831 CALL wrf_error_fatal ( TRIM(wrf_err_message) )
3833 CALL wrf_debug ( 200 , ' CLWRF: call diagnostic_calc' )
3834 CALL clwrf_output_calc( &
3835 & DPSDT=grid%dpsdt ,DMUDT=grid%dmudt &
3836 & ,P8W=p8w ,PK1M=grid%pk1m &
3837 & ,MU_2=grid%mu_2 ,MU_2M=grid%mu_2m &
3838 & ,U=grid%u_2 ,V=grid%v_2 &
3839 & ,is_restart=config_flags%restart &
3840 & ,clwrfH=config_flags%auxhist3_interval &
3841 & ,T2=grid%t2, Q2=grid%q2, U10=grid%u10, V10=grid%v10 &
3842 & ,SKINTEMP=grid%tsk &
3843 & ,T2CLMIN=grid%t2min, T2CLMAX=grid%t2max &
3844 & ,TT2CLMIN=grid%tt2min, TT2CLMAX=grid%tt2max &
3845 & ,T2CLMEAN=grid%t2mean, T2CLSTD=grid%t2std &
3846 & ,Q2CLMIN=grid%q2min, Q2CLMAX=grid%q2max &
3847 & ,TQ2CLMIN=grid%tq2min, TQ2CLMAX=grid%tq2max &
3848 & ,Q2CLMEAN=grid%q2mean, Q2CLSTD=grid%q2std &
3849 & ,U10CLMAX=grid%u10max, V10CLMAX=grid%v10max &
3850 & ,SPDUV10CLMAX=grid%spduv10max &
3851 & ,TSPDUV10CLMAX=grid%tspduv10max &
3852 & ,U10CLMEAN=grid%u10mean, V10CLMEAN=grid%v10mean &
3853 & ,SPDUV10CLMEAN=grid%spduv10mean &
3854 & ,U10CLSTD=grid%u10std, V10CLSTD=grid%v10std &
3855 & ,SPDUV10CLSTD=grid%spduv10std &
3856 & ,RAINCCLMAX=grid%raincvmax &
3857 & ,RAINNCCLMAX=grid%rainncvmax &
3858 & ,TRAINCCLMAX=grid%traincvmax &
3859 & ,TRAINNCCLMAX=grid%trainncvmax &
3860 & ,RAINCCLMEAN=grid%raincvmean &
3861 & ,RAINNCCLMEAN=grid%rainncvmean &
3862 & ,RAINCCLSTD=grid%raincvstd &
3863 & ,RAINNCCLSTD=grid%rainncvstd &
3864 & ,SKINTEMPCLMIN=grid%skintempmin &
3865 & ,SKINTEMPCLMAX=grid%skintempmax &
3866 & ,TSKINTEMPCLMIN=grid%tskintempmin &
3867 & ,TSKINTEMPCLMAX=grid%tskintempmax &
3868 & ,SKINTEMPCLMEAN=grid%skintempmean &
3869 & ,SKINTEMPCLSTD=grid%skintempstd &
3870 & ,RAINCV=grid%raincv ,RAINNCV=grid%rainncv &
3871 & ,RAINC=grid%rainc ,RAINNC=grid%rainnc &
3872 & ,I_RAINC=grid%i_rainc ,I_RAINNC=grid%i_rainnc &
3873 & ,HFX=grid%hfx ,SFCEVP=grid%sfcevp ,LH=grid%lh &
3874 & ,DT=grid%dt ,SBW=config_flags%spec_bdy_width &
3875 & ,XTIME=grid%xtime &
3877 & ,DIAG_PRINT=config_flags%diag_print &
3878 & ,BUCKET_MM=config_flags%bucket_mm &
3879 & ,BUCKET_J =config_flags%bucket_J &
3880 ! Dimension arguments
3881 & ,IDS=ids,IDE=ide, JDS=jds,JDE=jde, KDS=kds,KDE=kde &
3882 & ,IMS=ims,IME=ime, JMS=jms,JME=jme, KMS=kms,KME=kme &
3883 & ,IPS=ips,IPE=ipe, JPS=jps,JPE=jpe, KPS=kps,KPE=kpe &
3884 & ,I_START=grid%i_start,I_END=min(grid%i_end, ide-1) &
3885 & ,J_START=grid%j_start,J_END=min(grid%j_end, jde-1) &
3886 & ,KTS=k_start, KTE=min(k_end,kde-1) &
3887 & ,NUM_TILES=grid%num_tiles &
3892 !-----------------------------------------------------------------------
3894 !--------------------------------------------------------------
3895 CALL wrf_debug ( 200 , ' call HALO_RK_E' )
3896 IF ( config_flags%h_mom_adv_order <= 4 ) THEN
3897 # include "HALO_EM_E_3.inc"
3898 ELSE IF ( config_flags%h_mom_adv_order <= 6 ) THEN
3899 # include "HALO_EM_E_5.inc"
3901 WRITE(wrf_err_message,*)'solve_em: invalid h_mom_adv_order = ',config_flags%h_mom_adv_order
3902 CALL wrf_error_fatal(TRIM(wrf_err_message))
3907 IF ( num_moist >= PARAM_FIRST_SCALAR ) THEN
3908 !-----------------------------------------------------------------------
3910 !--------------------------------------------------------------
3911 CALL wrf_debug ( 200 , ' call HALO_RK_MOIST' )
3912 IF ( config_flags%h_mom_adv_order <= 4 ) THEN
3913 # include "HALO_EM_MOIST_E_3.inc"
3914 ELSE IF ( config_flags%h_mom_adv_order <= 6 ) THEN
3915 # include "HALO_EM_MOIST_E_5.inc"
3917 WRITE(wrf_err_message,*)'solve_em: invalid h_mom_adv_order = ',config_flags%h_mom_adv_order
3918 CALL wrf_error_fatal(TRIM(wrf_err_message))
3921 IF ( num_chem >= PARAM_FIRST_SCALAR ) THEN
3922 !-----------------------------------------------------------------------
3924 !--------------------------------------------------------------
3925 CALL wrf_debug ( 200 , ' call HALO_RK_CHEM' )
3926 IF ( config_flags%h_mom_adv_order <= 4 ) THEN
3927 # include "HALO_EM_CHEM_E_3.inc"
3928 ELSE IF ( config_flags%h_mom_adv_order <= 6 ) THEN
3929 # include "HALO_EM_CHEM_E_5.inc"
3931 WRITE(wrf_err_message,*)'solve_em: invalid h_mom_adv_order = ',config_flags%h_mom_adv_order
3932 CALL wrf_error_fatal(TRIM(wrf_err_message))
3935 IF ( num_tracer >= PARAM_FIRST_SCALAR ) THEN
3936 !-----------------------------------------------------------------------
3938 !--------------------------------------------------------------
3939 CALL wrf_debug ( 200 , ' call HALO_RK_TRACER' )
3940 IF ( config_flags%h_mom_adv_order <= 4 ) THEN
3941 # include "HALO_EM_TRACER_E_3.inc"
3942 ELSE IF ( config_flags%h_mom_adv_order <= 6 ) THEN
3943 # include "HALO_EM_TRACER_E_5.inc"
3945 WRITE(wrf_err_message,*)'solve_em: invalid h_mom_adv_order = ',config_flags%h_mom_adv_order
3946 CALL wrf_error_fatal(TRIM(wrf_err_message))
3949 IF ( num_scalar >= PARAM_FIRST_SCALAR ) THEN
3950 !-----------------------------------------------------------------------
3952 !--------------------------------------------------------------
3953 CALL wrf_debug ( 200 , ' call HALO_RK_SCALAR' )
3954 IF ( config_flags%h_mom_adv_order <= 4 ) THEN
3955 # include "HALO_EM_SCALAR_E_3.inc"
3956 ELSE IF ( config_flags%h_mom_adv_order <= 6 ) THEN
3957 # include "HALO_EM_SCALAR_E_5.inc"
3959 WRITE(wrf_err_message,*)'solve_em: invalid h_mom_adv_order = ',config_flags%h_mom_adv_order
3960 CALL wrf_error_fatal(TRIM(wrf_err_message))
3965 ! Max values of CFL for adaptive time step scheme
3967 DEALLOCATE(max_vert_cfl_tmp)
3968 DEALLOCATE(max_horiz_cfl_tmp)
3970 CALL wrf_debug ( 200 , ' call end of solve_em' )
3972 ! Finish timers if compiled with -DBENCH.
3973 #include <bench_solve_em_end.h>
3977 END SUBROUTINE solve_em