source: trunk/WRF.COMMON/WRFV2/dyn_em/solve_em.F @ 3026

Last change on this file since 3026 was 2211, checked in by aslmd, 5 years ago

MESOSCALE: adapted lmd_driver for changes in r2123 (LES was done in r2123, not needed)

File size: 170.0 KB
Line 
1!WRF:MEDIATION_LAYER:SOLVER
2
3SUBROUTINE solve_em ( grid , config_flags  &
4! Actual arguments generated from Registry
5#include "em_dummy_new_args.inc"
6!
7                    )
8
9! Driver layer modules
10   USE module_domain
11   USE module_configure
12   USE module_driver_constants
13   USE module_machine
14   USE module_tiles
15   USE module_dm
16! Mediation layer modules
17! Model layer modules
18   USE module_model_constants
19   USE module_small_step_em
20   USE module_em
21   USE module_big_step_utilities_em
22   USE module_bc
23   USE module_bc_em
24   USE module_solvedebug_em
25   USE module_physics_addtendc
26   USE module_diffusion_em
27! Registry generated module
28   USE module_state_description
29!!!******MARS MARS
30!!!******MARS MARS
31!!!******MARS MARS
32!   USE module_radiation_driver
33!   USE module_surface_driver
34!   USE module_cumulus_driver
35!   USE module_microphysics_driver
36!   USE module_microphysics_zero_out
37!   USE module_pbl_driver
38!   USE module_fddagd_driver
39!   USE module_fddaobs_driver
40!   USE module_diagnostics
41!****MARS
42! in the dyn_em folder physics driver are only used here
43   USE module_lmd_driver
44!****MARS
45#ifdef WRF_CHEM
46   USE module_input_chem_data
47   USE module_chem_utilities
48#endif
49
50   IMPLICIT NONE
51
52   !  Input data.
53
54   TYPE(domain) , TARGET          :: grid
55
56   !  Definitions of dummy arguments to this routine (generated from Registry).
57#include "em_dummy_new_decl.inc"
58
59   !  Structure that contains run-time configuration (namelist) data for domain
60   TYPE (grid_config_rec_type) , INTENT(IN)          :: config_flags
61
62   ! Local data
63
64   INTEGER                         :: k_start , k_end, its, ite, jts, jte
65   INTEGER                         :: ids , ide , jds , jde , kds , kde , &
66                                      ims , ime , jms , jme , kms , kme , &
67                                      ips , ipe , jps , jpe , kps , kpe
68   INTEGER                         :: ij , iteration
69   INTEGER                         :: im , num_3d_m , ic , num_3d_c , is , num_3d_s
70   INTEGER                         :: loop
71   INTEGER                         :: ijds, ijde
72   INTEGER                         :: sz
73
74   LOGICAL                         :: specified_bdy, channel_bdy
75
76
77! storage for tendencies and decoupled state (generated from Registry)
78
79#include <em_i1_decl.inc>
80! Previous time level of tracer arrays now defined as i1 variables;
81! the state 4d arrays now redefined as 1-time level arrays in Registry.
82! Benefit: save memory in nested runs, since only 1 domain is active at a
83! time.  Potential problem on stack-limited architectures: increases
84! amount of data on program stack by making these automatic arrays.
85
86   INTEGER :: rc
87   INTEGER :: number_of_small_timesteps, rk_step
88   INTEGER :: klevel,ijm,ijp,i,j,k,size1,size2    ! for prints/plots only
89   INTEGER :: idum1, idum2, dynamics_option
90
91   INTEGER :: rk_order, iwmax, jwmax, kwmax
92   REAL :: dt_rk, dts_rk, dtm, wmax
93   INTEGER :: l,kte,kk
94
95! urban related variables
96   INTEGER :: NUM_ROOF_LAYERS, NUM_WALL_LAYERS, NUM_ROAD_LAYERS   ! urban
97
98! Define benchmarking timers if -DBENCH is compiled
99#include <bench_solve_em_def.h>
100
101!----------------------
102! Executable statements
103!----------------------
104
105! Needed by some comm layers, grid%e.g. RSL. If needed, nmm_data_calls.inc is
106! generated from the registry.  The definition of REGISTER_I1 allows
107! I1 data to be communicated in this routine if necessary.
108#ifdef DM_PARALLEL
109#    define REGISTER_I1
110#      include "em_data_calls.inc"
111#endif
112
113!<DESCRIPTION>
114!<pre>
115! solve_em is the main driver for advancing a grid a single timestep.
116! It is a mediation-layer routine -> DM and SM calls are made where
117! needed for parallel processing. 
118!
119! solve_em can integrate the equations using 3 time-integration methods
120!     
121!    - 3rd order Runge-Kutta time integration (recommended)
122!     
123!    - 2nd order Runge-Kutta time integration
124!     
125! The main sections of solve_em are
126!     
127! (1) Runge-Kutta (RK) loop
128!     
129! (2) Non-timesplit physics (i.grid%e., tendencies computed for updating
130!     model state variables during the first RK sub-step (loop)
131!     
132! (3) Small (acoustic, sound) timestep loop - within the RK sub-steps
133!     
134! (4) scalar advance for moist and chem scalar variables (and TKE)
135!     within the RK sub-steps.
136!     
137! (5) time-split physics (after the RK step), currently this includes
138!     only microphyics
139!
140! A more detailed description of these sections follows.
141!</pre>
142!</DESCRIPTION>
143
144! Initialize timers if compiled with -DBENCH
145#include <bench_solve_em_init.h>
146
147!  set runge-kutta solver (2nd or 3rd order)
148
149   dynamics_option = config_flags%rk_ord
150
151!  Obtain dimension information stored in the grid data structure.
152  CALL get_ijk_from_grid (  grid ,                   &
153                            ids, ide, jds, jde, kds, kde,    &
154                            ims, ime, jms, jme, kms, kme,    &
155                            ips, ipe, jps, jpe, kps, kpe    )
156
157  k_start         = kps
158  k_end           = kpe
159
160  ijds = min(ids, jds)
161  ijde = max(ide, jde)
162
163  num_3d_m        = num_moist
164  num_3d_c        = num_chem
165  num_3d_s        = num_scalar
166
167!  Compute these starting and stopping locations for each tile and number of tiles.
168!  See: http://www.mmm.ucar.edu/wrf/WG2/topics/settiles
169  CALL set_tiles ( grid , ids , ide , jds , jde , ips , ipe , jps , jpe )
170!print *,grid%julian,grid%julday,' grid%julian,grid%julday in solve'
171
172  grid%itimestep = grid%itimestep + 1
173
174!**********************************************************************
175!
176!  LET US BEGIN.......
177!
178!<DESCRIPTION>
179!<pre>
180! (1) RK integration loop is named the "Runge_Kutta_loop:"
181!
182!   Predictor-corrector type time integration.
183!   Advection terms are evaluated at time t for the predictor step,
184!   and advection is re-evaluated with the latest predicted value for
185!   each succeeding time corrector step
186!
187!   2nd order Runge Kutta (rk_order = 2):
188!   Step 1 is taken to the midpoint predictor, step 2 is the full step.
189!
190!   3rd order Runge Kutta (rk_order = 3):
191!   Step 1 is taken to from t to dt/3, step 2 is from t to dt/2,
192!   and step 3 is from t to dt.
193!
194!   non-timesplit physics are evaluated during first RK step and
195!   these physics tendencies are stored for use in each RK pass.
196!</pre>
197!</DESCRIPTION>
198!**********************************************************************
199
200#ifdef WRF_CHEM
201!
202!    prepare chem aerosols for advection before communication
203!
204
205   kte=min(k_end,kde-1)
206# ifdef DM_PARALLEL
207   if ( num_chem >= PARAM_FIRST_SCALAR ) then
208!-----------------------------------------------------------------------
209! see above
210!--------------------------------------------------------------
211     CALL wrf_debug ( 200 , ' call HALO_RK_CHEM' )
212     IF      ( config_flags%h_mom_adv_order <= 4 ) THEN
213#      include "HALO_EM_CHEM_E_3.inc"
214     ELSE IF ( config_flags%h_mom_adv_order <= 6 ) THEN
215#      include "HALO_EM_CHEM_E_5.inc"
216     ELSE
217       WRITE(wrf_err_message,*)'solve_em: invalid h_mom_adv_order = ',config_flags%h_mom_adv_order
218       CALL wrf_error_fatal(TRIM(wrf_err_message))
219     ENDIF
220   endif
221# endif
222!--------------------------------------------------------------
223#endif
224
225 rk_order = config_flags%rk_ord
226 IF (grid%time_step_sound == 0) THEN
227! auto-set option
228! This function will give 4 for 6*dx and 6 for 10*dx and returns even numbers only
229   grid%time_step_sound = max ( 2 * ( INT (300.*grid%dt/grid%dx-0.01) + 1 ), 4 )
230   WRITE(wrf_err_message,*)'dx, dt, time_step_sound=',grid%dx,grid%dt,grid%time_step_sound
231   CALL wrf_debug ( 50 , wrf_err_message )
232 ENDIF
233
234 grid%dts = grid%dt/float(grid%time_step_sound)
235
236 Runge_Kutta_loop:  DO rk_step = 1, rk_order
237
238   !  Set the step size and number of small timesteps for
239   !  each part of the timestep
240
241   dtm = grid%dt
242   IF ( rk_order == 1 ) THEN   
243
244      write(wrf_err_message,*)' leapfrog removed, error exit for dynamics_option = ',dynamics_option
245      CALL wrf_error_fatal( wrf_err_message )
246
247   ELSE IF ( rk_order == 2 ) THEN   ! 2nd order Runge-Kutta timestep
248
249       IF ( rk_step == 1) THEN
250             dt_rk  = 0.5*grid%dt
251             dts_rk = grid%dts
252             number_of_small_timesteps = grid%time_step_sound/2
253       ELSE
254             dt_rk = grid%dt
255             dts_rk = grid%dts
256             number_of_small_timesteps = grid%time_step_sound
257       ENDIF
258
259   ELSE IF ( rk_order == 3 ) THEN ! third order Runge-Kutta
260
261       IF ( rk_step == 1) THEN
262            dt_rk = grid%dt/3.
263            dts_rk = dt_rk
264            number_of_small_timesteps = 1
265       ELSE IF (rk_step == 2) THEN
266            dt_rk  = 0.5*grid%dt
267            dts_rk = grid%dts
268            number_of_small_timesteps = grid%time_step_sound/2
269       ELSE
270            dt_rk = grid%dt
271            dts_rk = grid%dts
272            number_of_small_timesteps = grid%time_step_sound
273       ENDIF
274
275   ELSE
276
277      write(wrf_err_message,*)' unknown solver, error exit for dynamics_option = ',dynamics_option
278      CALL wrf_error_fatal( wrf_err_message )
279
280   END IF
281
282!
283!  Time level t is in the *_2 variable in the first part
284!  of the step, and in the *_1 variable after the predictor.
285!  the latest predicted values are stored in the *_2 variables.
286!
287   CALL wrf_debug ( 200 , ' call rk_step_prep ' )
288
289BENCH_START(step_prep_tim)
290   !$OMP PARALLEL DO   &
291   !$OMP PRIVATE ( ij )
292
293   DO ij = 1 , grid%num_tiles
294
295      CALL rk_step_prep  ( config_flags, rk_step,            &
296                           grid%em_u_2, grid%em_v_2, grid%em_w_2, grid%em_t_2, grid%em_ph_2, grid%em_mu_2,   &
297                           moist,                            &
298                           grid%em_ru, grid%em_rv, grid%em_rw, grid%em_ww, grid%em_php, grid%em_alt, grid%em_muu, grid%em_muv,   &
299                           grid%em_mub, grid%em_mut, grid%em_phb, grid%em_pb, grid%em_p, grid%em_al, grid%em_alb,    &
300                           cqu, cqv, cqw,                    &
301                           grid%msfu, grid%msfv, grid%msft,                 &
302                           grid%em_fnm, grid%em_fnp, grid%em_dnw, grid%rdx, grid%rdy,          &
303                           num_3d_m,                         &
304                           ids, ide, jds, jde, kds, kde,     &
305                           ims, ime, jms, jme, kms, kme,     &
306                           grid%i_start(ij), grid%i_end(ij), &
307                           grid%j_start(ij), grid%j_end(ij), &
308                           k_start, k_end                   )
309
310   END DO
311   !$OMP END PARALLEL DO
312BENCH_END(step_prep_tim)
313
314#ifdef DM_PARALLEL
315!-----------------------------------------------------------------------
316!  Stencils for patch communications  (WCS, 29 June 2001)
317!  Note:  the small size of this halo exchange reflects the
318!         fact that we are carrying the uncoupled variables
319!         as state variables in the mass coordinate model, as
320!         opposed to the coupled variables as in the height
321!         coordinate model.
322!
323!                           * * * * *
324!         *        * * *    * * * * *
325!       * + *      * + *    * * + * *
326!         *        * * *    * * * * *
327!                           * * * * *
328!
329!  3D variables - note staggering!  grid%em_ru(X), grid%em_rv(Y), grid%em_ww(grid%em_z), grid%em_php(grid%em_z)
330!
331!j grid%em_ru     x
332!j grid%em_rv     x
333!j grid%em_ww     x
334!j grid%em_php    x
335!j grid%em_alt    x
336!j grid%em_ph_2   x
337!j grid%em_phb    x
338!
339!  the following are 2D (xy) variables
340!
341!j grid%em_muu    x
342!j grid%em_muv    x
343!j grid%em_mut    x
344!--------------------------------------------------------------
345#    include "HALO_EM_A.inc"
346#endif
347
348! set boundary conditions on variables
349! from big_step_prep for use in big_step_proc
350
351#ifdef DM_PARALLEL
352#  include "PERIOD_BDY_EM_A.inc"
353#endif
354
355!   CALL set_tiles ( grid , ids , ide , jds , jde , ips-1 , ipe+1 , jps-1 , jpe+1 )
356
357BENCH_START(set_phys_bc_tim)
358   !$OMP PARALLEL DO   &
359   !$OMP PRIVATE ( ij )
360
361   DO ij = 1 , grid%num_tiles
362
363       CALL wrf_debug ( 200 , ' call rk_phys_bc_dry_1' )
364
365        CALL rk_phys_bc_dry_1( config_flags, grid%em_ru, grid%em_rv, grid%em_rw, grid%em_ww,      &
366                               grid%em_muu, grid%em_muv, grid%em_mut, grid%em_php, grid%em_alt, grid%em_p,        &
367                               ids, ide, jds, jde, kds, kde,      &
368                               ims, ime, jms, jme, kms, kme,      &
369                               ips, ipe, jps, jpe, kps, kpe,      &
370                               grid%i_start(ij), grid%i_end(ij),  &
371                               grid%j_start(ij), grid%j_end(ij),  &
372                               k_start, k_end                )
373       !TBH:  need this 2nd timestep and later
374       CALL set_physical_bc3d( grid%em_al, 'p', config_flags,            &
375                               ids, ide, jds, jde, kds, kde,     &
376                               ims, ime, jms, jme, kms, kme,     &
377                               ips, ipe, jps, jpe, kps, kpe,     &
378                               grid%i_start(ij), grid%i_end(ij), &
379                               grid%j_start(ij), grid%j_end(ij), &
380                               k_start    , k_end               )
381       CALL set_physical_bc3d( grid%em_ph_2, 'w', config_flags,            &
382                                 ids, ide, jds, jde, kds, kde, &
383                                 ims, ime, jms, jme, kms, kme, &
384                                 ips, ipe, jps, jpe, kps, kpe, &
385                               grid%i_start(ij), grid%i_end(ij),        &
386                               grid%j_start(ij), grid%j_end(ij),        &
387                               k_start, k_end                )
388
389   END DO
390   !$OMP END PARALLEL DO
391BENCH_END(set_phys_bc_tim)
392
393    rk_step_is_one : IF (rk_step == 1) THEN ! only need to initialize diffusion tendencies
394
395 ! initialize all tendencies to zero in order to update physics
396 ! tendencies first (separate from dry dynamics).
397 
398BENCH_START(init_zero_tend_tim)
399     !$OMP PARALLEL DO   &
400     !$OMP PRIVATE ( ij )
401
402     DO ij = 1 , grid%num_tiles
403
404        CALL wrf_debug ( 200 , ' call init_zero_tendency' )
405        CALL init_zero_tendency ( ru_tendf, rv_tendf, rw_tendf,     &
406                                  ph_tendf, t_tendf, tke_tend,      &
407                                  mu_tendf,                         &
408                                  moist_tend,chem_tend,scalar_tend, &
409                                  num_3d_m,num_3d_c,num_3d_s,       &
410                                  rk_step,                          &
411                                  ids, ide, jds, jde, kds, kde,     &
412                                  ims, ime, jms, jme, kms, kme,     &
413                                  grid%i_start(ij), grid%i_end(ij), &
414                                  grid%j_start(ij), grid%j_end(ij), &
415                                  k_start, k_end                   )
416
417     END DO
418   !$OMP END PARALLEL DO
419BENCH_END(init_zero_tend_tim)
420
421#ifdef DM_PARALLEL
422#     include "HALO_EM_PHYS_A.inc"
423#endif
424
425!<DESCRIPTION>
426!<pre>
427!(2) The non-timesplit physics begins with a call to "phy_prep"
428!    (which computes some diagnostic variables such as temperature,
429!    pressure, u and v at grid%em_p points, etc).  This is followed by
430!    calls to the physics drivers:
431!
432!              radiation,
433!              surface,
434!              pbl,
435!              cumulus,
436!              3D TKE and mixing.
437!<pre>
438!</DESCRIPTION>
439
440
441BENCH_START(phy_prep_tim)
442      !$OMP PARALLEL DO   &
443      !$OMP PRIVATE ( ij )
444      DO ij = 1 , grid%num_tiles
445
446         CALL wrf_debug ( 200 , ' call phy_prep' )
447         CALL phy_prep ( config_flags,                           &
448                         grid%em_mut, grid%em_muu, grid%em_muv, grid%em_u_2, &
449                         grid%em_v_2, grid%em_w_2, grid%em_p, grid%em_pb, grid%em_alt,              &
450                         grid%em_ph_2, grid%em_phb, grid%em_t_2, grid%tsk, moist, num_3d_m,   &
451                         mu_3d, rho,                             &
452                         th_phy, p_phy, pi_phy, u_phy, v_phy, w_phy,    &
453                         p8w, t_phy, t8w, grid%em_z, z_at_w,             &
454                         dz8w, grid%em_fnm, grid%em_fnp,                         &   
455                         grid%rthraten,                               &
456                         grid%rthblten, grid%rublten, grid%rvblten,             &
457                         grid%rqvblten, grid%rqcblten, grid%rqiblten,           &
458                         grid%rthcuten, grid%rqvcuten, grid%rqccuten,           &
459                         grid%rqrcuten, grid%rqicuten, grid%rqscuten,           &
460                         grid%rthften,  grid%rqvften,                      &
461                         grid%RUNDGDTEN, grid%RVNDGDTEN, grid%RTHNDGDTEN,       &
462                         grid%RQVNDGDTEN, grid%RMUNDGDTEN,                 &
463                         ids, ide, jds, jde, kds, kde,           &
464                         ims, ime, jms, jme, kms, kme,           &
465                         grid%i_start(ij), grid%i_end(ij),       &
466                         grid%j_start(ij), grid%j_end(ij),       &
467                         k_start, k_end                         )
468      ENDDO
469      !$OMP END PARALLEL DO
470
471
472BENCH_END(phy_prep_tim)
473
474!  physics to implement
475
476!      CALL set_tiles ( grid , ids , ide-1 , jds , jde-1 ips , ipe , jps , jpe )
477
478! Open MP loops are in physics drivers
479! radiation
480
481!-----------------------------------------------------------------
482! urban related variable are added to arguments of radiation_driver
483!-----------------------------------------------------------------
484
485!!!******MARS MARS
486!!!******MARS MARS
487!!!******MARS MARS
488
489
490!         CALL wrf_debug ( 200 , ' call radiation_driver' )
491!BENCH_START(rad_driver_tim)
492!
493!         CALL radiation_driver(                                           &
494!     &         ACFRCV=grid%acfrcv      ,ACFRST=grid%acfrst      ,ALBEDO=grid%albedo      &
495!     &        ,CFRACH=grid%cfrach      ,CFRACL=grid%cfracl      ,CFRACM=grid%cfracm      &
496!     &        ,CUPPT=grid%cuppt        ,CZMEAN=grid%czmean      ,DT=grid%dt              &
497!     &        ,DZ8W=dz8w          ,EMISS=grid%emiss        ,GLW=grid%glw            &
498!     &        ,GMT=grid%gmt            ,GSW=grid%gsw            ,HBOT=grid%hbot          &
499!     &        ,HTOP=grid%htop ,HBOTR=grid%hbotr, HTOPR=grid%htopr ,ICLOUD=config_flags%icloud &
500!     &        ,ITIMESTEP=grid%itimestep,JULDAY=grid%julday, JULIAN=grid%julian      &
501!     &        ,JULYR=grid%julyr        ,LW_PHYSICS=config_flags%ra_lw_physics  &
502!     &        ,NCFRCV=grid%ncfrcv      ,NCFRST=grid%ncfrst      ,NPHS=1             &
503!     &        ,P8W=p8w            ,P=p_phy            ,PI=pi_phy          &
504!     &        ,RADT=grid%radt     ,RA_CALL_OFFSET=grid%ra_call_offset     &
505!     &        ,RHO=rho            ,RLWTOA=grid%rlwtoa                          &
506!     &        ,RSWTOA=grid%rswtoa      ,RTHRATEN=grid%rthraten                      &
507!     &        ,RTHRATENLW=grid%rthratenlw                                      &
508!     &        ,RTHRATENSW=grid%rthratensw                  ,SNOW=grid%snow          &
509!     &        ,STEPRA=grid%stepra      ,SWDOWN=grid%swdown      ,SWDOWNC=grid%swdownc    &
510!     &        ,SW_PHYSICS=config_flags%ra_sw_physics  ,T8W=t8w            &
511!     &        ,T=t_phy            ,TAUCLDC=grid%taucldc    ,TAUCLDI=grid%taucldi    &
512!     &        ,TSK=grid%tsk            ,VEGFRA=grid%vegfra     ,WARM_RAIN=grid%warm_rain &
513!     &        ,XICE=grid%xice                                                  &
514!     &        ,XLAND=grid%xland        ,XLAT=grid%xlat          ,XLONG=grid%xlong        &
515!!Optional urban
516!     &        ,DECLIN_URB=grid%declin_urb        ,COSZ_URB2D=grid%cosz_urb2d        &
517!     &        ,OMG_URB2D=grid%omg_urb2d                                        &
518!!
519!     &        ,Z=grid%em_z                                                        &
520!     &        ,LEVSIZ=grid%levsiz, N_OZMIXM=num_ozmixm                    &
521!     &        ,N_AEROSOLC=num_aerosolc                                    &
522!     &        ,PAERLEV=grid%paerlev                                       &
523!     &        ,CAM_ABS_DIM1=grid%cam_abs_dim1, CAM_ABS_DIM2=grid%cam_abs_dim2 &
524!     &        ,CAM_ABS_FREQ_S=grid%cam_abs_freq_s                         &
525!     &        ,XTIME=grid%xtime                                                &
526!            ! indexes
527!     &        ,IDS=ids,IDE=ide, JDS=jds,JDE=jde, KDS=kds,KDE=kde          &
528!     &        ,IMS=ims,IME=ime, JMS=jms,JME=jme, KMS=kms,KME=kme          &
529!     &        ,i_start=grid%i_start,i_end=min(grid%i_end, ide-1)          &
530!     &        ,j_start=grid%j_start,j_end=min(grid%j_end, jde-1)          &
531!     &        ,kts=k_start, kte=min(k_end,kde-1)                          &
532!     &        ,num_tiles=grid%num_tiles                                   &
533!            ! Optional                         
534!     &        , CLDFRA=grid%cldfra                                             &
535!     &        , PB=grid%em_pb                                                     &
536!     &        , F_ICE_PHY=grid%f_ice_phy,F_RAIN_PHY=grid%f_rain_phy                 &
537!     &        , QV=moist(ims,kms,jms,P_QV), F_QV=F_QV                     &
538!     &        , QC=moist(ims,kms,jms,P_QC), F_QC=F_QC                                        &
539!     &        , QR=moist(ims,kms,jms,P_QR), F_QR=F_QR                     &
540!     &        , QI=moist(ims,kms,jms,P_QI), F_QI=F_QI                                        &
541!     &        , QS=moist(ims,kms,jms,P_QS), F_QS=F_QS                     &
542!     &        , QG=moist(ims,kms,jms,P_QG), F_QG=F_QG                     &
543!#ifdef ACFLUX
544!     &        ,ACSWUPT=acswupt    ,ACSWUPTC=acswuptc                      &
545!     &        ,ACSWDNT=acswdnt    ,ACSWDNTC=acswdntc                      &
546!     &        ,ACSWUPB=acswupb    ,ACSWUPBC=acswupbc                      &
547!     &        ,ACSWDNB=acswdnb    ,ACSWDNBC=acswdnbc                      &
548!     &        ,ACLWUPT=aclwupt    ,ACLWUPTC=aclwuptc                      &
549!     &        ,ACLWDNT=aclwdnt    ,ACLWDNTC=aclwdntc                      &
550!     &        ,ACLWUPB=aclwupb    ,ACLWUPBC=aclwupbc                      &
551!     &        ,ACLWDNB=aclwdnb    ,ACLWDNBC=aclwdnbc                      &
552!     &        ,SWUPT=swupt    ,SWUPTC=swuptc                              &
553!     &        ,SWDNT=swdnt    ,SWDNTC=swdntc                              &
554!     &        ,SWUPB=swupb    ,SWUPBC=swupbc                              &
555!     &        ,SWDNB=swdnb    ,SWDNBC=swdnbc                              &
556!     &        ,LWUPT=lwupt    ,LWUPTC=lwuptc                              &
557!     &        ,LWDNT=lwdnt    ,LWDNTC=lwdntc                              &
558!     &        ,LWUPB=lwupb    ,LWUPBC=lwupbc                              &
559!     &        ,LWDNB=lwdnb    ,LWDNBC=lwdnbc                              &
560!#endif
561!     &        ,LWCF=grid%lwcf                                                  &
562!     &        ,SWCF=grid%swcf                                                  &
563!     &        ,OLR=grid%olr                                                    &
564!     &        ,OZMIXM=grid%ozmixm, PIN=grid%pin                                     &
565!     &        ,M_PS_1=grid%m_ps_1, M_PS_2=grid%m_ps_2, AEROSOLC_1=grid%aerosolc_1        &
566!     &        ,AEROSOLC_2=grid%aerosolc_2, M_HYBI0=grid%m_hybi                      &
567!     &        ,ABSTOT=grid%abstot, ABSNXT=grid%absnxt, EMSTOT=grid%emstot                &
568!#ifdef WRF_CHEM
569!     &        ,QC_ADJUST=grid%GD_CLOUD , QI_ADJUST=grid%GD_CLOUD2         &
570!     &        ,PM2_5_DRY=grid%pm2_5_dry, PM2_5_WATER=grid%pm2_5_water               &
571!     &        ,PM2_5_DRY_EC=grid%pm2_5_dry_ec                                  &
572!     &        ,TAUAER300=grid%tauaer1, TAUAER400=grid%tauaer2 & ! jcb
573!     &        ,TAUAER600=grid%tauaer3, TAUAER999=grid%tauaer4 & ! jcb
574!     &        ,GAER300=grid%gaer1, GAER400=grid%gaer2, GAER600=grid%gaer3, GAER999=grid%gaer4 & ! jcb
575!     &        ,WAER300=grid%waer1, WAER400=grid%waer2, WAER600=grid%waer3, WAER999=grid%waer4 & ! jcb
576!#endif
577!     &                                                              )
578!
579!BENCH_END(rad_driver_tim)
580!
581!!********* Surface driver
582!! surface
583!
584!BENCH_START(surf_driver_tim)
585!
586!!-----------------------------------------------------------------
587!! urban related variable are added to arguments of surface_driver
588!!-----------------------------------------------------------------
589!      NUM_ROOF_LAYERS = grid%num_soil_layers !urban
590!      NUM_WALL_LAYERS = grid%num_soil_layers !urban
591!      NUM_ROAD_LAYERS = grid%num_soil_layers !urban
592!
593!      CALL wrf_debug ( 200 , ' call surface_driver' )
594!      CALL surface_driver(                                                &
595!     &         ACSNOM=grid%acsnom      ,ACSNOW=grid%acsnow      ,AKHS=grid%akhs          &
596!     &        ,AKMS=grid%akms          ,ALBBCK=grid%albbck      ,ALBEDO=grid%albedo      &
597!     &        ,BR=br              ,CANWAT=grid%canwat      ,CHKLOWQ=chklowq    &
598!     &        ,CT=grid%ct              ,DT=grid%dt         ,DX=grid%dx         &
599!     &        ,DZ8W=dz8w          ,DZS=grid%dzs            ,FLHC=grid%flhc          &
600!     &        ,FLQC=grid%flqc          ,GLW=grid%glw            ,GRDFLX=grid%grdflx      &
601!     &        ,GSW=grid%gsw    ,SWDOWN=grid%swdown        ,GZ1OZ0=gz1oz0      ,HFX=grid%hfx              &
602!     &        ,HT=grid%ht              ,IFSNOW=config_flags%ifsnow      ,ISFFLX=config_flags%isfflx      &
603!     &        ,ISLTYP=grid%isltyp      ,ITIMESTEP=grid%itimestep                    &
604!     &        ,IVGTYP=grid%ivgtyp      ,LH=grid%lh              ,LOWLYR=grid%lowlyr      &
605!     &        ,MAVAIL=grid%mavail      ,NUM_SOIL_LAYERS=config_flags%num_soil_layers        &
606!     &        ,P8W=p8w            ,PBLH=grid%pblh          ,PI_PHY=pi_phy      &
607!     &        ,PSFC=grid%psfc          ,PSHLTR=pshltr      ,PSIH=psih          &
608!     &        ,PSIM=psim          ,P_PHY=p_phy        ,Q10=q10            &
609!     &        ,Q2=grid%q2              ,QFX=grid%qfx            ,QSFC=grid%qsfc          &
610!     &        ,QSHLTR=qshltr      ,QZ0=grid%qz0            ,RAINCV=grid%raincv      &
611!     &        ,RA_LW_PHYSICS=config_flags%ra_lw_physics            ,RHO=rho            &
612!     &        ,RMOL=grid%rmol          ,SFCEVP=grid%sfcevp      ,SFCEXC=grid%sfcexc      &
613!     &        ,SFCRUNOFF=grid%sfcrunoff                                        &
614!     &        ,SF_SFCLAY_PHYSICS=config_flags%sf_sfclay_physics                        &
615!     &        ,SF_SURFACE_PHYSICS=config_flags%sf_surface_physics  ,SH2O=grid%sh2o          &
616!     &        ,SHDMAX=grid%shdmax      ,SHDMIN=grid%shdmin      ,SMOIS=grid%smois        &
617!     &        ,SMSTAV=grid%smstav      ,SMSTOT=grid%smstot      ,SNOALB=grid%snoalb      &
618!     &        ,SNOW=grid%snow          ,SNOWC=grid%snowc        ,SNOWH=grid%snowh        &
619!     &        ,SST=grid%sst            ,SST_UPDATE=grid%sst_update                  &
620!     &        ,STEPBL=grid%stepbl      ,TH10=th10          ,TH2=grid%th2            &
621!     &        ,THZ0=grid%thz0          ,TH_PHY=th_phy      ,TKE_MYJ=grid%tke_myj    &
622!     &        ,TMN=grid%tmn            ,TSHLTR=tshltr      ,TSK=grid%tsk            &
623!     &        ,TSLB=grid%tslb          ,T_PHY=t_phy        ,U10=grid%u10            &
624!     &        ,URATX=grid%uratx        ,VRATX=grid%vratx   ,TRATX=grid%tratx        &
625!     &        ,UDRUNOFF=grid%udrunoff  ,UST=grid%ust       ,UZ0=grid%uz0            &
626!     &        ,U_FRAME=grid%u_frame    ,U_PHY=u_phy        ,V10=grid%v10            &
627!     &        ,VEGFRA=grid%vegfra      ,VZ0=grid%vz0       ,V_FRAME=grid%v_frame    &
628!     &        ,V_PHY=v_phy             ,WARM_RAIN=grid%warm_rain                    &
629!     &        ,WSPD=wspd               ,XICE=grid%xice     ,XLAND=grid%xland        &
630!     &        ,Z0=grid%z0              ,Z=grid%em_z        ,ZNT=grid%znt            &
631!     &        ,ZS=grid%zs                                                           &
632!     &        ,DECLIN_URB=grid%declin_urb  ,COSZ_URB2D=grid%cosz_urb2d    & !I urban
633!     &        ,OMG_URB2D=grid%omg_urb2d    ,xlat_urb2d=grid%XLAT          & !I urban
634!     &        ,NUM_ROOF_LAYERS=num_roof_layers                            & !I urban
635!     &        ,NUM_WALL_LAYERS=num_wall_layers                            & !I urban
636!     &        ,NUM_ROAD_LAYERS=num_road_layers                            &
637!     &        ,DZR=grid%dzr ,DZB=grid%dzb ,DZG=grid%dzg                   & !I urban
638!     &        ,TR_URB2D=grid%tr_urb2d ,TB_URB2D=grid%tb_urb2d             &
639!     &        ,TG_URB2D=grid%tg_urb2d                                     & !H urban
640!     &        ,TC_URB2D=grid%tc_urb2d ,QC_URB2D=grid%qc_urb2d             & !H urban
641!     &        ,UC_URB2D=grid%uc_urb2d                                     & !H urban
642!     &        ,XXXR_URB2D=grid%xxxr_urb2d                                 &
643!     &        ,XXXB_URB2D=grid%xxxb_urb2d                                 & !H urban
644!     &        ,XXXG_URB2D=grid%xxxg_urb2d                                 &
645!     &        ,XXXC_URB2D=grid%xxxc_urb2d                                 & !H urban
646!     &        ,TRL_URB3D=grid%trl_urb3d   ,TBL_URB3D=grid%tbl_urb3d       & !H urban
647!     &        ,TGL_URB3D=grid%tgl_urb3d                                   & !H urban
648!     &        ,SH_URB2D=grid%sh_urb2d     ,LH_URB2D=grid%lh_urb2d         &
649!     &        ,G_URB2D=grid%g_urb2d                                       & !H urban
650!     &        ,RN_URB2D=grid%rn_urb2d     , TS_URB2D=grid%ts_urb2d        & !H urban
651!     &        ,FRC_URB2D=grid%frc_urb2d                                   & !H urban
652!     &        ,UTYPE_URB2D=grid%utype_urb2d                               & !H urban
653!     &        ,ucmcall=grid%ucmcall                                       & !H urban
654!           ! Indexes
655!     &        ,IDS=ids,IDE=ide, JDS=jds,JDE=jde, KDS=kds,KDE=kde          &
656!     &        ,IMS=ims,IME=ime, JMS=jms,JME=jme, KMS=kms,KME=kme          &
657!     &        , I_START=grid%i_start,I_END=min(grid%i_end, ide-1)         &
658!     &        , J_START=grid%j_start,J_END=min(grid%j_end, jde-1)         &
659!     &        , KTS=k_start, KTE=min(k_end,kde-1)                         &
660!     &        , NUM_TILES=grid%num_tiles                                  &
661!           ! Optional
662!     &        ,QV_CURR=moist(ims,kms,jms,P_QV), F_QV=F_QV                 &
663!     &        ,QC_CURR=moist(ims,kms,jms,P_QC), F_QC=F_QC                 &
664!     &        ,QR_CURR=moist(ims,kms,jms,P_QR), F_QR=F_QR                 &
665!     &        ,QI_CURR=moist(ims,kms,jms,P_QI), F_QI=F_QI                 &
666!     &        ,QS_CURR=moist(ims,kms,jms,P_QS), F_QS=F_QS                 &
667!     &        ,QG_CURR=moist(ims,kms,jms,P_QG), F_QG=F_QG                 &
668!     &        ,CAPG=grid%capg, EMISS=grid%emiss, HOL=hol,MOL=grid%mol                    &
669!     &        ,RAINBL=grid%rainbl,SR=grid%em_sr                                              &
670!     &        ,RAINNCV=grid%rainncv,REGIME=regime,T2=grid%t2,THC=grid%thc                &
671!     &        ,QSG=grid%qsg,QVG=grid%qvg,QCG=grid%qcg,SOILT1=grid%soilt1,TSNAV=grid%tsnav          & ! ruc lsm
672!     &        ,SMFR3D=grid%smfr3d,KEEPFR3DFLAG=grid%keepfr3dflag                    & ! ruc lsm
673!     &        ,POTEVP=grid%em_POTEVP, SNOPCX=grid%em_SNOPCX, SOILTB=grid%em_SOILTB                & ! ruc lsm
674!     &                                                              )
675!BENCH_END(surf_driver_tim)
676!
677!!*********
678!! pbl
679!
680!      CALL wrf_debug ( 200 , ' call pbl_driver' )
681!BENCH_START(pbl_driver_tim)
682!      CALL pbl_driver(                                                    &
683!     &         AKHS=grid%akhs          ,AKMS=grid%akms                              &
684!     &        ,BL_PBL_PHYSICS=config_flags%bl_pbl_physics                 &
685!     &        ,BR=br              ,CHKLOWQ=chklowq    ,CT=grid%ct              &
686!     &        ,DT=grid%dt              ,DX=grid%dx              ,DZ8W=dz8w          &
687!     &        ,EL_MYJ=grid%el_myj      ,EXCH_H=grid%exch_h      ,GRDFLX=grid%grdflx      &
688!     &        ,GZ1OZ0=gz1oz0      ,HFX=grid%hfx            ,HT=grid%ht              &
689!     &        ,ITIMESTEP=grid%itimestep                    ,KPBL=grid%kpbl          &
690!     &        ,LH=grid%lh              ,LOWLYR=grid%lowlyr      ,P8W=p8w            &
691!     &        ,PBLH=grid%pblh          ,PI_PHY=pi_phy      ,PSIH=psih          &
692!     &        ,PSIM=psim          ,P_PHY=p_phy        ,QFX=grid%qfx            &
693!     &        ,QSFC=grid%qsfc          ,QZ0=grid%qz0                                &
694!     &        ,RA_LW_PHYSICS=config_flags%ra_lw_physics                   &
695!     &        ,RHO=rho            ,RQCBLTEN=grid%rqcblten  ,RQIBLTEN=grid%rqiblten  &
696!     &        ,RQVBLTEN=grid%rqvblten  ,RTHBLTEN=grid%rthblten  ,RUBLTEN=grid%rublten    &
697!     &        ,RVBLTEN=grid%rvblten    ,SNOW=grid%snow          ,STEPBL=grid%stepbl      &
698!     &        ,THZ0=grid%thz0          ,TH_PHY=th_phy      ,TKE_MYJ=grid%tke_myj    &
699!     &        ,TSK=grid%tsk            ,T_PHY=t_phy        ,UST=grid%ust            &
700!     &        ,UZ0=grid%uz0            ,U_FRAME=grid%u_frame    ,U_PHY=u_phy        &
701!     &        ,VZ0=grid%vz0            ,V_FRAME=grid%v_frame    ,V_PHY=v_phy        &
702!     &        ,WARM_RAIN=grid%warm_rain                    ,WSPD=wspd          &
703!     &        ,XICE=grid%xice          ,XLAND=grid%xland        ,Z=grid%em_z                &
704!     &        ,ZNT=grid%znt                                                    &
705!     &        ,IDS=ids,IDE=ide, JDS=jds,JDE=jde, KDS=kds,KDE=kde          &
706!     &        ,IMS=ims,IME=ime, JMS=jms,JME=jme, KMS=kms,KME=kme          &
707!     &        ,I_START=grid%i_start,I_END=min(grid%i_end, ide-1)          &
708!     &        ,J_START=grid%j_start,J_END=min(grid%j_end, jde-1)          &
709!     &        ,KTS=k_start, KTE=min(k_end,kde-1)                          &
710!     &        ,NUM_TILES=grid%num_tiles                                   &
711!          ! optional
712!     &        ,QV_CURR=moist(ims,kms,jms,P_QV), F_QV=F_QV                 &
713!     &        ,QC_CURR=moist(ims,kms,jms,P_QC), F_QC=F_QC                 &
714!     &        ,QR_CURR=moist(ims,kms,jms,P_QR), F_QR=F_QR                 &
715!     &        ,QI_CURR=moist(ims,kms,jms,P_QI), F_QI=F_QI                 &
716!     &        ,QS_CURR=moist(ims,kms,jms,P_QS), F_QS=F_QS                 &
717!     &        ,QG_CURR=moist(ims,kms,jms,P_QG), F_QG=F_QG                 &
718!     &        ,HOL=HOL, MOL=grid%mol, REGIME=REGIME                            &
719!     &                                                          )
720!
721!BENCH_END(pbl_driver_tim)
722!
723!! cumulus para.
724!
725!          CALL wrf_debug ( 200 , ' call cumulus_driver' )
726!
727!BENCH_START(cu_driver_tim)
728!         CALL cumulus_driver(                                             &
729!                 ! Prognostic variables
730!     &              U=u_phy   ,V=v_phy   ,TH=th_phy  ,T=t_phy             &
731!     &             ,W=grid%em_w_2     ,P=p_phy   ,PI=pi_phy  ,RHO=rho             &
732!                 ! Other arguments
733!     &             ,ITIMESTEP=grid%itimestep ,DT=grid%dt      ,DX=grid%dx                &
734!     &             ,RAINC=grid%rainc   ,RAINCV=grid%raincv   ,NCA=grid%nca               &
735!     &             ,HTOP=grid%cutop    ,HBOT=grid%cubot      ,KPBL=grid%kpbl             &
736!     &             ,DZ8W=dz8w     ,P8W=p8w                                &
737!     &             ,W0AVG=grid%w0avg   ,STEPCU=grid%stepcu                          &
738!     &             ,CLDEFI=grid%cldefi ,LOWLYR=grid%lowlyr ,XLAND=grid%xland             &
739!     &             ,APR_GR=grid%apr_gr ,APR_W=grid%apr_w   ,APR_MC=grid%apr_mc           &
740!     &             ,APR_ST=grid%apr_st ,APR_AS=grid%apr_as ,APR_CAPMA=grid%apr_capma     &
741!     &             ,APR_CAPME=grid%apr_capme          ,APR_CAPMI=grid%apr_capmi     &
742!     &             ,MASS_FLUX=grid%mass_flux          ,XF_ENS=grid%xf_ens           &
743!     &             ,PR_ENS=grid%pr_ens ,HT=grid%ht                                  &
744!     &             ,ENSDIM=config_flags%ensdim ,MAXIENS=config_flags%maxiens ,MAXENS=config_flags%maxens         &
745!     &             ,MAXENS2=config_flags%maxens2                ,MAXENS3=config_flags%maxens3       &
746!     &             ,CU_ACT_FLAG=cu_act_flag   ,WARM_RAIN=grid%warm_rain        &
747!     &             ,GSW=grid%gsw                                               &
748!                 ! Selection flag
749!     &             ,CU_PHYSICS=config_flags%cu_physics                    &
750!                 ! Dimension arguments
751!     &             ,IDS=ids,IDE=ide, JDS=jds,JDE=jde, KDS=kds,KDE=kde     &
752!     &             ,IMS=ims,IME=ime, JMS=jms,JME=jme, KMS=kms,KME=kme     &
753!     &             ,I_START=grid%i_start,I_END=min(grid%i_end, ide-1)     &
754!     &             ,J_START=grid%j_start,J_END=min(grid%j_end, jde-1)     &
755!     &             ,KTS=k_start, KTE=min(k_end,kde-1)                     &
756!     &             ,NUM_TILES=grid%num_tiles                              &
757!                 ! Moisture tendency arguments
758!     &             ,RQVCUTEN=grid%rqvcuten , RQCCUTEN=grid%rqccuten                 &
759!     &             ,RQSCUTEN=grid%rqscuten , RQICUTEN=grid%rqicuten                 &
760!     &             ,RQRCUTEN=grid%rqrcuten , RQVBLTEN=grid%rqvblten                 &
761!     &             ,RQVFTEN=grid%rqvften                                       &
762!                 ! Other tendency arguments
763!     &             ,RTHRATEN=grid%rthraten , RTHBLTEN=grid%rthblten                 &
764!     &             ,RTHCUTEN=grid%rthcuten , RTHFTEN=grid%rthften                   &
765!                 ! Moisture tracer arguments
766!     &             ,QV_CURR=moist(ims,kms,jms,P_QV), F_QV=F_QV            &
767!     &             ,QC_CURR=moist(ims,kms,jms,P_QC), F_QC=F_QC            &
768!     &             ,QR_CURR=moist(ims,kms,jms,P_QR), F_QR=F_QR            &
769!     &             ,QI_CURR=moist(ims,kms,jms,P_QI), F_QI=F_QI            &
770!     &             ,QS_CURR=moist(ims,kms,jms,P_QS), F_QS=F_QS            &
771!     &             ,QG_CURR=moist(ims,kms,jms,P_QG), F_QG=F_QG            &
772!#ifdef WRF_CHEM
773!     &             ,GD_CLOUD=grid%GD_CLOUD,GD_CLOUD2=grid%GD_CLOUD2                          &
774!#endif
775!     &                                                          )
776!BENCH_END(cu_driver_tim)
777!
778!! fdda
779!
780!          CALL wrf_debug ( 200 , ' call fddagd_driver' )
781!
782!BENCH_START(fdda_driver_tim)
783!   CALL fddagd_driver(itimestep=grid%itimestep,dt=grid%dt,xtime=grid%XTIME,         &
784!                  id=grid%id,      &
785!                  RUNDGDTEN=grid%rundgdten,RVNDGDTEN=grid%rvndgdten,                &
786!                  RTHNDGDTEN=grid%rthndgdten,RQVNDGDTEN=grid%rqvndgdten,            &
787!                  RMUNDGDTEN=grid%rmundgdten,                                  &
788!                  u_ndg_old=fdda3d(ims,kms,jms,P_u_ndg_old),              &
789!                  v_ndg_old=fdda3d(ims,kms,jms,P_v_ndg_old),              &
790!                  t_ndg_old=fdda3d(ims,kms,jms,P_t_ndg_old),              &
791!                  q_ndg_old=fdda3d(ims,kms,jms,P_q_ndg_old),              &
792!                  mu_ndg_old=fdda2d(ims,1,jms,P_mu_ndg_old),              &
793!                  u_ndg_new=fdda3d(ims,kms,jms,P_u_ndg_new),              &
794!                  v_ndg_new=fdda3d(ims,kms,jms,P_v_ndg_new),              &
795!                  t_ndg_new=fdda3d(ims,kms,jms,P_t_ndg_new),              &
796!                  q_ndg_new=fdda3d(ims,kms,jms,P_q_ndg_new),              &
797!                  mu_ndg_new=fdda2d(ims,1,jms,P_mu_ndg_new),              &
798!                  u3d=grid%em_u_2,v3d=grid%em_v_2,th_phy=th_phy,rho=rho,moist=moist,      &
799!                  p_phy=p_phy,pi_phy=pi_phy,p8w=p8w,t_phy=t_phy,          &
800!                  dz8w=dz8w,z=grid%em_z,z_at_w=z_at_w,                            &
801!                  config_flags=config_flags,dx=grid%DX,n_moist=num_3d_m,  &
802!                  STEPFG=grid%STEPFG,                                          &
803!                  pblh=grid%pblh,ht=grid%ht,                                        &
804!                    IDS=ids,IDE=ide, JDS=jds,JDE=jde, KDS=kds,KDE=kde     &
805!                   ,IMS=ims,IME=ime, JMS=jms,JME=jme, KMS=kms,KME=kme     &
806!                   ,I_START=grid%i_start,I_END=min(grid%i_end, ide-1)     &
807!                   ,J_START=grid%j_start,J_END=min(grid%j_end, jde-1)     &
808!                   ,KTS=k_start, KTE=min(k_end,kde-1)                     &
809!                   , num_tiles=grid%num_tiles                             )
810!BENCH_END(fdda_driver_tim)
811
812
813!****MARS
814IF (config_flags%modif_wrf) THEN
815!!!!!!!!!!!!!!!!!!!!!!!
816! call to LMD physics !
817!!!!!!!!!!!!!!!!!!!!!!!
818         CALL wrf_debug ( 200 , ' call lmd_driver' )
819         CALL lmd_driver(                                                 &
820           ! structure
821     &         id=grid%id,max_dom=grid%max_dom      &
822     &        ,DT=grid%dt                   &
823     &        ,ITIMESTEP=grid%itimestep                                   &
824           ! position
825     &        ,XLAT=grid%xlat,XLONG=grid%xlong                            &
826     &        ,DX=grid%dx         ,DY=grid%dy                             &
827     &        ,MSFT=grid%msft,MSFU=grid%msfu,MSFV=grid%msfv               &
828           ! indexes
829     &        ,IDS=ids,IDE=ide, JDS=jds,JDE=jde, KDS=kds,KDE=kde          &
830     &        ,IMS=ims,IME=ime, JMS=jms,JME=jme, KMS=kms,KME=kme          &
831     &        ,i_start=grid%i_start,i_end=min(grid%i_end, ide-1)          &
832     &        ,j_start=grid%j_start,j_end=min(grid%j_end, jde-1)          &
833     &        ,kts=k_start, kte=min(k_end,kde-1)                          &
834     &        ,num_tiles=grid%num_tiles                                   &
835           ! time variables
836     &        ,GMT=grid%gmt       ,JULYR=grid%julyr   ,JULDAY=grid%julday &
837           ! output from phy_prep
838     &        ,P8W=p8w            ,DZ8W=dz8w          ,T8W=t8w            &
839           ! Prognostic variables at p points
840     &        ,Z=grid%em_z                                                &
841                        ! z is the geopotential height at p points
842                        ! ... (cf phy_prep in module_big_step_utilities_em)
843     &        ,HT=grid%ht                                                 &
844     &        ,U=u_phy            ,V=v_phy             ,W=w_phy           &
845     &        ,TH=th_phy          ,T=t_phy                                &
846     &        ,P=p_phy            ,EXNER=pi_phy               ,RHO=rho    &
847!!!!!ADDITION CAS IDEALISES
848     &        ,PTOP=grid%p_top                                            &
849           ! physics calls
850     &        ,RADT=grid%radt                          &
851     &        ,CUDT=grid%cudt                          &
852           ! surface temperature and surface pressure
853     &        ,TSK=grid%tsk,PSFC=grid%psfc                                &
854           ! tendencies for the dynamics
855     &        ,RTHBLTEN=grid%rthblten  ,RUBLTEN=grid%rublten,RVBLTEN=grid%rvblten    &
856           ! scalars
857     &        ,NUM_3D_S=num_3d_s,SCALAR=grid%scalar  &
858     &        ,NUM_3D_M=num_moist,MOIST=grid%moist   &
859           ! mode
860     &        ,MARS_MODE=config_flags%mars   &
861#ifdef NEWPHYS
862           ! planet
863     &        ,planet_type=config_flags%planet   &
864#endif
865           ! added variables
866     &        ,M_ALBEDO=grid%m_albedo,M_TI=grid%m_ti  &
867     &        ,M_CO2ICE=grid%m_co2ice,M_EMISS=grid%m_emiss   & 
868     &        ,M_H2OICE=grid%m_h2oice     &
869     &        ,M_TSOIL=grid%m_tsoil       &
870     &        ,M_Q2=grid%m_q2             &
871     &        ,M_TSURF=grid%m_tsurf       &
872#ifdef NEWPHYS
873     &        ,M_FLUXRAD=grid%m_fluxrad   &
874     &        ,M_WSTAR=grid%m_wstar       &
875     &        ,M_ISOIL=grid%m_isoil       &
876     &        ,M_DSOIL=grid%m_dsoil       &
877     &        ,M_Z0=grid%m_z0             &
878     &        ,CST_Z0=config_flags%init_Z0      &
879#endif
880     &        ,M_GW=grid%m_gw  & 
881     &        ,NUM_SOIL_LAYERS=config_flags%num_soil_layers    &
882           ! cst fields
883     &        ,CST_AL=config_flags%init_AL  &
884     &        ,CST_TI=config_flags%init_TI  &
885     &        ,ISFFLX=config_flags%isfflx   &
886     &        ,DIFF_OPT=config_flags%diff_opt     &
887     &        ,KM_OPT=config_flags%km_opt         &
888     &        ,HISTORY_INTERVAL=model_config_rec%history_interval(1)  &
889          !------------------!
890          ! OUTPUT VARIABLES !
891          !------------------!
892#ifdef NEWPHYS
893     &       ,HR_SW=grid%HR_SW   &
894     &       ,HR_LW=grid%HR_LW   &
895     &       ,HR_DYN=grid%HR_DYN   &
896     &       ,DDT=grid%DDT  &
897     &       ,DT_RAD=grid%DT_RAD  &
898     &       ,DT_VDF=grid%DT_VDF  &
899     &       ,DT_AJS=grid%DT_AJS  &
900     &       ,CLOUDFRAC=grid%CLOUDFRAC   &
901     &       ,TOTCLOUDFRAC=grid%TOTCLOUDFRAC   &
902     &       ,GRAIN=grid%GRAIN   &
903     &       ,GSNOW=grid%GSNOW   &
904     &       ,REEVAP=grid%REEVAP  &
905     &       ,SURFRAIN=grid%SURFRAIN  &
906     &       ,ALBEQ=grid%ALBEQ   &
907     &       ,FLUXTOP_DN=grid%FLUXTOP_DN   &
908     &       ,FLUXABS_SW=grid%FLUXABS_SW   &
909     &       ,FLUXTOP_LW=grid%FLUXTOP_LW   &
910     &       ,FLUXSURF_SW=grid%FLUXSURF_SW   &
911     &       ,FLUXSURF_LW=grid%FLUXSURF_LW   &
912     &       ,FLXGRD=grid%FLXGRD   &
913     &       ,LSCEZ=grid%LSCEZ   &
914     &       ,H2OICE_REFF=grid%H2OICE_REFF  &
915     &       ,LATENT_HF=grid%LATENT_HF  &
916     &       ,SWDOWNZ=grid%SWDOWNZ   &
917     &       ,TAU_DUST=grid%TAU_DUST   &
918     &       ,RDUST=grid%RDUST   &
919     &       ,QSURFDUST=grid%QSURFDUST   &
920     &       ,MTOT=grid%MTOT   &
921     &       ,ICETOT=grid%ICETOT   &
922     &       ,VMR_ICE=grid%VMR_ICE   &
923     &       ,TAU_ICE=grid%TAU_ICE   &
924     &       ,RICE=grid%RICE   &
925     &       ,ZMAX=grid%ZMAX   &
926     &       ,HFMAX=grid%HFMAX &
927     &       ,USTM=grid%USTM   &
928     &       ,HFX=grid%HFX &
929#else
930#include "module_lmd_driver_output4.inc"
931#endif
932     &        ,SLPX=grid%slpx,SLPY=grid%slpy,RESTART=config_flags%restart)
933ENDIF
934
935!!!!!!!!!!!!!!!!!!!!!!!
936! call to LMD physics !
937!!!!!!!!!!!!!!!!!!!!!!!
938!****MARS
939
940
941
942! calculate_phy_tend
943
944BENCH_START(cal_phy_tend)
945      !$OMP PARALLEL DO   &
946      !$OMP PRIVATE ( ij )
947
948      DO ij = 1 , grid%num_tiles
949
950          CALL wrf_debug ( 200 , ' call calculate_phy_tend' )
951          CALL calculate_phy_tend (config_flags,grid%em_mut,grid%em_muu,grid%em_muv,pi_phy,            &
952                     grid%rthraten,                                         &
953                     grid%rublten,grid%rvblten,grid%rthblten,                         &
954                     grid%rqvblten,grid%rqcblten,grid%rqiblten,                       &
955                     grid%rthcuten,grid%rqvcuten,grid%rqccuten,grid%rqrcuten,              &
956                     grid%rqicuten,grid%rqscuten,                                &
957                     grid%RUNDGDTEN,grid%RVNDGDTEN,grid%RTHNDGDTEN,grid%RQVNDGDTEN,        &
958                     grid%RMUNDGDTEN,                                       &
959                     ids,ide, jds,jde, kds,kde,                        &
960                     ims,ime, jms,jme, kms,kme,                        &
961                     grid%i_start(ij), min(grid%i_end(ij),ide-1),      &
962                     grid%j_start(ij), min(grid%j_end(ij),jde-1),      &
963                     k_start    , min(k_end,kde-1)                     )
964
965      ENDDO
966      !$OMP END PARALLEL DO
967BENCH_END(cal_phy_tend)
968
969! tke diffusion
970
971     IF(config_flags%diff_opt .eq. 2 .OR. config_flags%diff_opt .eq. 1) THEN
972
973BENCH_START(comp_diff_metrics_tim)
974       !$OMP PARALLEL DO   &
975       !$OMP PRIVATE ( ij )
976
977       DO ij = 1 , grid%num_tiles
978
979          CALL wrf_debug ( 200 , ' call compute_diff_metrics ' )
980          CALL compute_diff_metrics ( config_flags, grid%em_ph_2, grid%em_phb, grid%em_z, grid%em_rdz, grid%em_rdzw, &
981                                      grid%em_zx, grid%em_zy, grid%rdx, grid%rdy,                      &
982                                      ids, ide, jds, jde, kds, kde,          &
983                                      ims, ime, jms, jme, kms, kme,          &
984                                      grid%i_start(ij), grid%i_end(ij),      &
985                                      grid%j_start(ij), grid%j_end(ij),      &
986                                      k_start    , k_end                    )
987       ENDDO
988       !$OMP END PARALLEL DO
989BENCH_END(comp_diff_metrics_tim)
990
991#ifdef DM_PARALLEL
992#  include "PERIOD_BDY_EM_A1.inc"
993#endif
994
995BENCH_START(tke_diff_bc_tim)
996       DO ij = 1 , grid%num_tiles
997
998          CALL wrf_debug ( 200 , ' call bc for diffusion_metrics ' )
999          CALL set_physical_bc3d( grid%em_rdzw , 'w', config_flags,           &
1000                                  ids, ide, jds, jde, kds, kde,       &
1001                                  ims, ime, jms, jme, kms, kme,       &
1002                                  ips, ipe, jps, jpe, kps, kpe,       &
1003                                  grid%i_start(ij), grid%i_end(ij),   &
1004                                  grid%j_start(ij), grid%j_end(ij),   &
1005                                  k_start    , k_end                 )
1006          CALL set_physical_bc3d( grid%em_rdz , 'w', config_flags,            &
1007                                  ids, ide, jds, jde, kds, kde,       &
1008                                  ims, ime, jms, jme, kms, kme,       &
1009                                  ips, ipe, jps, jpe, kps, kpe,       &
1010                                  grid%i_start(ij), grid%i_end(ij),   &
1011                                  grid%j_start(ij), grid%j_end(ij),   &
1012                                  k_start    , k_end                 )
1013          CALL set_physical_bc3d( grid%em_z , 'w', config_flags,              &
1014                                  ids, ide, jds, jde, kds, kde,       &
1015                                  ims, ime, jms, jme, kms, kme,       &
1016                                  ips, ipe, jps, jpe, kps, kpe,       &
1017                                  grid%i_start(ij), grid%i_end(ij),   &
1018                                  grid%j_start(ij), grid%j_end(ij),   &
1019                                  k_start    , k_end                 )
1020          CALL set_physical_bc3d( grid%em_zx , 'w', config_flags,             &
1021                                  ids, ide, jds, jde, kds, kde,       &
1022                                  ims, ime, jms, jme, kms, kme,       &
1023                                  ips, ipe, jps, jpe, kps, kpe,       &
1024                                  grid%i_start(ij), grid%i_end(ij),   &
1025                                  grid%j_start(ij), grid%j_end(ij),   &
1026                                  k_start    , k_end                 )
1027          CALL set_physical_bc3d( grid%em_zy , 'w', config_flags,             &
1028                                  ids, ide, jds, jde, kds, kde,       &
1029                                  ims, ime, jms, jme, kms, kme,       &
1030                                  ips, ipe, jps, jpe, kps, kpe,       &
1031                                  grid%i_start(ij), grid%i_end(ij),   &
1032                                  grid%j_start(ij), grid%j_end(ij),   &
1033                                  k_start    , k_end                 )
1034
1035       ENDDO
1036BENCH_END(tke_diff_bc_tim)
1037
1038#ifdef DM_PARALLEL
1039#     include "HALO_EM_TKE_C.inc"
1040#endif
1041
1042BENCH_START(deform_div_tim)
1043
1044       !$OMP PARALLEL DO   &
1045       !$OMP PRIVATE ( ij )
1046
1047       DO ij = 1 , grid%num_tiles
1048
1049          CALL wrf_debug ( 200 , ' call cal_deform_and_div' )
1050          CALL cal_deform_and_div ( config_flags,grid%em_u_2,grid%em_v_2,grid%em_w_2,grid%div,        &
1051                                    grid%defor11,grid%defor22,grid%defor33,grid%defor12,     &
1052                                    grid%defor13,grid%defor23,                     &
1053                                    grid%u_base, grid%v_base,grid%msfu,grid%msfv,grid%msft,       &
1054                                    grid%rdx, grid%rdy, grid%em_dn, grid%em_dnw, grid%em_rdz, grid%em_rdzw,        &
1055                                    grid%em_fnm,grid%em_fnp,grid%cf1,grid%cf2,grid%cf3,grid%em_zx,grid%em_zy,           &
1056                                    ids, ide, jds, jde, kds, kde,        &
1057                                    ims, ime, jms, jme, kms, kme,        &
1058                                    grid%i_start(ij), grid%i_end(ij),    &
1059                                    grid%j_start(ij), grid%j_end(ij),    &
1060                                    k_start    , k_end                  )
1061       ENDDO
1062       !$OMP END PARALLEL DO
1063BENCH_END(deform_div_tim)
1064
1065
1066#ifdef DM_PARALLEL
1067#     include "HALO_EM_TKE_D.inc"
1068#endif
1069
1070
1071! calculate tke, kmh, and kmv
1072
1073BENCH_START(calc_tke_tim)
1074       !$OMP PARALLEL DO   &
1075       !$OMP PRIVATE ( ij )
1076
1077       DO ij = 1 , grid%num_tiles
1078
1079          CALL wrf_debug ( 200 , ' call calculate_km_kh' )
1080          CALL calculate_km_kh( config_flags,grid%dt,grid%dampcoef,grid%zdamp,config_flags%damp_opt,     &
1081                                grid%xkmh,grid%xkmhd,grid%xkmv,grid%xkhh,grid%xkhv,grid%bn2,               &
1082                                grid%khdif,grid%kvdif,grid%div,                             &
1083                                grid%defor11,grid%defor22,grid%defor33,grid%defor12,             &
1084                                grid%defor13,grid%defor23,                             &
1085                                grid%em_tke_2,p8w,t8w,th_phy,           &
1086                                t_phy,p_phy,moist,grid%em_dn,grid%em_dnw,                    &
1087                                grid%dx,grid%dy,grid%em_rdz,grid%em_rdzw,config_flags%mix_cr_len,num_3d_m,          &
1088                                grid%cf1, grid%cf2, grid%cf3, grid%warm_rain,                    &
1089                                grid%kh_tke_upper_bound, grid%kv_tke_upper_bound,      &
1090                                ids,ide, jds,jde, kds,kde,                   &
1091                                ims,ime, jms,jme, kms,kme,                   &
1092                                grid%i_start(ij), grid%i_end(ij),            &
1093                                grid%j_start(ij), grid%j_end(ij),            &
1094                                k_start    , k_end                          )
1095       ENDDO
1096       !$OMP END PARALLEL DO
1097BENCH_END(calc_tke_tim)
1098
1099#ifdef DM_PARALLEL
1100#     include "HALO_EM_TKE_E.inc"
1101#endif
1102
1103     ENDIF
1104
1105#ifdef DM_PARALLEL
1106#      include "PERIOD_BDY_EM_PHY_BC.inc"
1107      IF ( config_flags%grid_fdda .eq. 1) THEN
1108#      include "PERIOD_BDY_EM_FDDA_BC.inc"
1109      ENDIF
1110#      include "PERIOD_BDY_EM_CHEM.inc"
1111#endif
1112
1113BENCH_START(phy_bc_tim)
1114     !$OMP PARALLEL DO   &
1115     !$OMP PRIVATE ( ij )
1116
1117!
1118!****MARS: TODO: check that this is OK, because physical tendency is an input
1119!
1120     DO ij = 1 , grid%num_tiles
1121
1122       CALL wrf_debug ( 200 , ' call phy_bc' )
1123       CALL phy_bc (config_flags,grid%div,grid%defor11,grid%defor22,grid%defor33,            &
1124                            grid%defor12,grid%defor13,grid%defor23,                     &
1125                            grid%xkmh,grid%xkmhd,grid%xkmv,grid%xkhh,grid%xkhv,                   &
1126                            grid%em_tke_2,                          &
1127                            grid%rublten, grid%rvblten,                            &
1128                            ids, ide, jds, jde, kds, kde,                &
1129                            ims, ime, jms, jme, kms, kme,                &
1130                            ips, ipe, jps, jpe, kps, kpe,                &
1131                            grid%i_start(ij), grid%i_end(ij),            &
1132                            grid%j_start(ij), grid%j_end(ij),            &
1133                            k_start    , k_end                           )
1134     ENDDO
1135     !$OMP END PARALLEL DO
1136BENCH_END(phy_bc_tim)
1137
1138#ifdef DM_PARALLEL
1139!-----------------------------------------------------------------------
1140!
1141! MPP for some physics tendency, km, kh, deformation, and divergence
1142!
1143!               *                     *
1144!             * + *      * + *        +
1145!               *                     *
1146!
1147! (for PBL)
1148! grid%rublten                  x
1149! grid%rvblten                             x
1150!
1151! (for diff_opt >= 1)
1152! grid%defor11                  x
1153! grid%defor22                             x
1154! grid%defor12       x
1155! grid%defor13                  x
1156! grid%defor23                             x
1157! grid%div           x
1158! grid%xkmv          x
1159! grid%xkmh          x
1160! grid%xkmhd         x
1161! grid%xkhv          x
1162! grid%xkhh          x
1163! tke           x
1164!
1165!-----------------------------------------------------------------------
1166
1167!!****MARS: always include this HALO for Mars version ...
1168      IF ( ( config_flags%bl_pbl_physics .ge. 1 ) &
1169            .OR. ( config_flags%modif_wrf ) ) THEN
1170#      include "HALO_EM_PHYS_PBL.inc"
1171      ENDIF
1172
1173
1174      IF ( config_flags%grid_fdda .eq. 1) THEN
1175#      include "HALO_EM_FDDA.inc"
1176      ENDIF
1177      IF ( config_flags%diff_opt .ge. 1 ) THEN
1178#      include "HALO_EM_PHYS_DIFFUSION.inc"
1179      ENDIF
1180
1181      IF      ( config_flags%h_mom_adv_order <= 4 ) THEN
1182#       include "HALO_EM_TKE_3.inc"
1183      ELSE IF ( config_flags%h_mom_adv_order <= 6 ) THEN
1184#       include "HALO_EM_TKE_5.inc"
1185      ELSE
1186        WRITE(wrf_err_message,*)'solve_em: invalid h_mom_adv_order = ',config_flags%h_mom_adv_order
1187        CALL wrf_error_fatal(TRIM(wrf_err_message))
1188      ENDIF
1189#endif
1190
1191BENCH_START(update_phy_ten_tim)
1192      !$OMP PARALLEL DO   &
1193      !$OMP PRIVATE ( ij )
1194
1195      DO ij = 1 , grid%num_tiles
1196
1197          CALL wrf_debug ( 200 , ' call update_phy_ten' )
1198        CALL update_phy_ten(t_tendf, ru_tendf, rv_tendf,moist_tend,        &
1199                          mu_tendf,                                        &
1200                          grid%rthraten,grid%rthblten,grid%rthcuten,grid%rublten,grid%rvblten,      &
1201                          grid%rqvblten,grid%rqcblten,grid%rqiblten,                      &
1202                          grid%rqvcuten,grid%rqccuten,grid%rqrcuten,grid%rqicuten,grid%rqscuten,    &
1203                          grid%RUNDGDTEN,grid%RVNDGDTEN,grid%RTHNDGDTEN,grid%RQVNDGDTEN,       &
1204                          grid%RMUNDGDTEN,                                      &
1205                          num_3d_m,config_flags,rk_step,                   &
1206                          ids, ide, jds, jde, kds, kde,                    &
1207                          ims, ime, jms, jme, kms, kme,                    &
1208                          grid%i_start(ij), grid%i_end(ij),                &
1209                          grid%j_start(ij), grid%j_end(ij),                &
1210                          k_start, k_end                               )
1211
1212      END DO
1213      !$OMP END PARALLEL DO
1214BENCH_END(update_phy_ten_tim)
1215
1216     IF( config_flags%diff_opt .eq. 2 .and. config_flags%km_opt .eq. 2 ) THEN
1217
1218BENCH_START(tke_rhs_tim)
1219       !$OMP PARALLEL DO   &
1220       !$OMP PRIVATE ( ij )
1221
1222       DO ij = 1 , grid%num_tiles
1223
1224          CALL tke_rhs  ( tke_tend,grid%bn2,                               &
1225                          config_flags,grid%defor11,grid%defor22,grid%defor33,       &
1226                          grid%defor12,grid%defor13,grid%defor23,grid%em_u_2,grid%em_v_2,grid%em_w_2,grid%div,    &
1227                          grid%em_tke_2,grid%em_mut,                     &
1228                          th_phy,p_phy,p8w,t8w,grid%em_z,grid%em_fnm,grid%em_fnp,             &
1229                          grid%cf1,grid%cf2,grid%cf3,grid%msft,grid%xkmh,grid%xkmv,grid%xkhv,grid%rdx,grid%rdy,    &
1230                          grid%dx,grid%dy,grid%dt,grid%em_zx,grid%em_zy,grid%em_rdz,grid%em_rdzw,grid%em_dn,       &
1231                          grid%em_dnw,config_flags%mix_cr_len,  &
1232                          ids, ide, jds, jde, kds, kde,               &
1233                          ims, ime, jms, jme, kms, kme,               &
1234                          grid%i_start(ij), grid%i_end(ij),           &
1235                          grid%j_start(ij), grid%j_end(ij),           &
1236                          k_start    , k_end                         )
1237
1238       ENDDO
1239       !$OMP END PARALLEL DO
1240BENCH_END(tke_rhs_tim)
1241
1242     ENDIF
1243
1244! calculate vertical diffusion first and then horizontal
1245! (keep this order)
1246
1247     IF(config_flags%diff_opt .eq. 2) THEN
1248
1249!!!!****MARS: vertical diffusion is done by the physics
1250       IF ( (config_flags%bl_pbl_physics .eq. 0) &
1251             .AND. (.not. config_flags%modif_wrf ) ) THEN
1252
1253BENCH_START(vert_diff_tim)
1254         !$OMP PARALLEL DO   &
1255         !$OMP PRIVATE ( ij )
1256         DO ij = 1 , grid%num_tiles
1257
1258           CALL wrf_debug ( 200 , ' call vertical_diffusion_2 ' )
1259           CALL vertical_diffusion_2( ru_tendf, rv_tendf, rw_tendf,              &
1260                                      t_tendf, tke_tend,                         &
1261                                      moist_tend, num_3d_m,                      &
1262                                      chem_tend, num_3d_c,                       &
1263                                      scalar_tend, num_3d_s,                     &
1264                                      grid%em_u_2, grid%em_v_2,                                  &
1265                                      grid%em_t_2,grid%u_base,grid%v_base,grid%em_t_base,grid%qv_base,          &
1266                                      grid%em_mut,grid%em_tke_2,config_flags,                    &
1267                                      grid%defor13,grid%defor23,grid%defor33,                   &
1268                                      grid%div, moist, chem, scalar,                  &
1269                                      grid%xkmv, grid%xkhv, config_flags%km_opt,                        &
1270                                      grid%em_fnm, grid%em_fnp, grid%em_dn, grid%em_dnw, grid%em_rdz, grid%em_rdzw,              &
1271                                      ids, ide, jds, jde, kds, kde,              &
1272                                      ims, ime, jms, jme, kms, kme,              &
1273                                      grid%i_start(ij), grid%i_end(ij),          &
1274                                      grid%j_start(ij), grid%j_end(ij),          &
1275                                      k_start    , k_end                        )
1276
1277         ENDDO
1278         !$OMP END PARALLEL DO
1279BENCH_END(vert_diff_tim)
1280       ENDIF
1281
1282!
1283BENCH_START(hor_diff_tim)
1284       !$OMP PARALLEL DO   &
1285       !$OMP PRIVATE ( ij )
1286       DO ij = 1 , grid%num_tiles
1287
1288         CALL wrf_debug ( 200 , ' call horizontal_diffusion_2' )
1289         CALL horizontal_diffusion_2( t_tendf, ru_tendf, rv_tendf, rw_tendf, &
1290                                      tke_tend,                              &
1291                                      moist_tend, num_3d_m,                  &
1292                                      chem_tend, num_3d_c,                   &
1293                                      scalar_tend, num_3d_s,                 &
1294                                      grid%em_t_2, th_phy,                           &
1295                                      grid%em_mut, grid%em_tke_2, config_flags,              &
1296                                      grid%defor11, grid%defor22, grid%defor12,             &
1297                                      grid%defor13, grid%defor23, grid%div,                 &
1298                                      moist, chem, scalar,                   &
1299                                      grid%msfu, grid%msfv, grid%msft, grid%xkmhd, grid%xkhh, config_flags%km_opt, &
1300                                      grid%rdx, grid%rdy, grid%em_rdz, grid%em_rdzw,                   &
1301                                      grid%em_fnm, grid%em_fnp, grid%cf1, grid%cf2, grid%cf3,               &
1302                                      grid%em_zx, grid%em_zy, grid%em_dn, grid%em_dnw,                       &
1303                                      ids, ide, jds, jde, kds, kde,          &
1304                                      ims, ime, jms, jme, kms, kme,          &
1305                                      grid%i_start(ij), grid%i_end(ij),      &
1306                                      grid%j_start(ij), grid%j_end(ij),      &
1307                                      k_start    , k_end                    )
1308       ENDDO
1309       !$OMP END PARALLEL DO
1310BENCH_END(hor_diff_tim)
1311
1312     ENDIF
1313
1314!# ifdef DM_PARALLEL
1315!#     include "HALO_OBS_NUDGE.inc"
1316!#endif
1317!!***********************************************************************
1318!! This section for obs nudging
1319!      !$OMP PARALLEL DO   &
1320!      !$OMP PRIVATE ( ij )
1321!
1322!      DO ij = 1 , grid%num_tiles
1323!
1324!         CALL fddaobs_driver (grid%grid_id, model_config_rec%grid_id,  &
1325!                  model_config_rec%parent_id, config_flags%restart,    &
1326!                  grid%obs_nudge_opt,                                  &
1327!                  grid%obs_ipf_errob,                                  &
1328!                  grid%obs_ipf_nudob,                                  &
1329!                  grid%fdda_start,                                     &
1330!                  grid%fdda_end,                                       &
1331!                  grid%obs_nudge_wind,                                 &
1332!                  grid%obs_nudge_temp,                                 &
1333!                  grid%obs_nudge_mois,                                 &
1334!                  grid%obs_nudge_pstr,                                 &
1335!                  grid%obs_coef_wind,                                  &
1336!                  grid%obs_coef_temp,                                  &
1337!                  grid%obs_coef_mois,                                  &
1338!                  grid%obs_coef_pstr,                                  &             
1339!                  grid%obs_rinxy,                                      &
1340!                  grid%obs_rinsig,                                     &
1341!                  grid%obs_twindo,                                     &
1342!                  grid%obs_npfi,                                       &
1343!                  grid%obs_ionf,                                       &
1344!                  grid%obs_idynin,                                     &
1345!                  grid%obs_dtramp,                                     &
1346!                  model_config_rec%cen_lat(1),                         &
1347!                  model_config_rec%cen_lon(1),                         &
1348!                  config_flags%truelat1,                               &
1349!                  config_flags%truelat2,                               &
1350!                  config_flags%map_proj,                               &
1351!                  model_config_rec%i_parent_start,                     &
1352!                  model_config_rec%j_parent_start,                     &
1353!                  grid%parent_grid_ratio,                              &
1354!                  grid%max_dom, grid%itimestep,                        &
1355!                  grid%dt, grid%gmt, grid%julday, grid%fdob,           &
1356!                  grid%max_obs,                                        &
1357!                  model_config_rec%nobs_ndg_vars,                      &
1358!                  model_config_rec%nobs_err_flds,                      &
1359!                  grid%fdob%nstat, grid%fdob%varobs, grid%fdob%errf,   &
1360!                  grid%dx, grid%KPBL,grid%HT,                          &
1361!                  grid%em_mut, grid%em_muu, grid%em_muv,               &
1362!                  grid%msft, grid%msfu, grid%msfv,                     &
1363!                  p_phy, t_tendf, t0,                                  &
1364!                  grid%em_u_2, grid%em_v_2, grid%em_t_2,               &
1365!                  moist(:,:,:,P_QV),                                   &
1366!                  grid%em_pb, grid%p_top, grid%em_p,                   &
1367!                  grid%uratx, grid%vratx, grid%tratx,                  &
1368!                  ru_tendf, rv_tendf,                                  &
1369!                  moist_tend(:,:,:,P_QV), grid%em_obs_savwt,           &
1370!                  ids,ide, jds,jde, kds,kde,                           &
1371!                  ims,ime, jms,jme, kms,kme,                           &
1372!                  grid%i_start(ij), min(grid%i_end(ij),ide-1),         &
1373!                  grid%j_start(ij), min(grid%j_end(ij),jde-1),         &
1374!                  k_start    , min(k_end,kde-1)                     )
1375!
1376!      ENDDO
1377!
1378!     !$OMP END PARALLEL DO
1379!!
1380!!***********************************************************************
1381
1382     END IF rk_step_is_one
1383
1384BENCH_START(rk_tend_tim)
1385   !$OMP PARALLEL DO   &
1386   !$OMP PRIVATE ( ij )
1387   DO ij = 1 , grid%num_tiles
1388
1389      CALL wrf_debug ( 200 , ' call rk_tendency' )
1390      CALL rk_tendency ( config_flags, rk_step,                           &
1391                         grid%em_ru_tend, grid%em_rv_tend, rw_tend, ph_tend, t_tend,      &
1392                         ru_tendf, rv_tendf, rw_tendf, ph_tendf, t_tendf, &
1393                         mu_tend, grid%em_u_save, grid%em_v_save, w_save, ph_save,        &
1394                         grid%em_t_save, mu_save, grid%rthften,                        &
1395                         grid%em_ru, grid%em_rv, grid%em_rw, grid%em_ww,                                  &
1396                         grid%em_u_2, grid%em_v_2, grid%em_w_2, grid%em_t_2, grid%em_ph_2,                        &
1397                         grid%em_u_1, grid%em_v_1, grid%em_w_1, grid%em_t_1, grid%em_ph_1,                        &
1398                         grid%h_diabatic, grid%em_phb, grid%em_t_init,                         &
1399                         grid%em_mu_2, grid%em_mut, grid%em_muu, grid%em_muv, grid%em_mub,                        &
1400                         grid%em_al, grid%em_alt, grid%em_p, grid%em_pb, grid%em_php, cqu, cqv, cqw,              &
1401                         grid%u_base, grid%v_base, grid%em_t_base, grid%qv_base, grid%z_base,         &
1402                         grid%msfu, grid%msfv, grid%msft, grid%f, grid%e, grid%sina, grid%cosa,              &
1403                         grid%em_fnm, grid%em_fnp, grid%em_rdn, grid%em_rdnw,                             &
1404                         grid%dt, grid%rdx, grid%rdy, grid%khdif, grid%kvdif, grid%xkmhd,               &
1405                         grid%diff_6th_opt, grid%diff_6th_factor,           &
1406                         grid%dampcoef,grid%zdamp,config_flags%damp_opt,                         &
1407                         grid%cf1, grid%cf2, grid%cf3, grid%cfn, grid%cfn1, num_3d_m,              &
1408                         config_flags%non_hydrostatic,                    &
1409                         ids, ide, jds, jde, kds, kde,                    &
1410                         ims, ime, jms, jme, kms, kme,                    &
1411                         grid%i_start(ij), grid%i_end(ij),                &
1412                         grid%j_start(ij), grid%j_end(ij),                &
1413                         k_start, k_end                                  )
1414   END DO
1415   !$OMP END PARALLEL DO
1416BENCH_END(rk_tend_tim)
1417
1418BENCH_START(relax_bdy_dry_tim)
1419   !$OMP PARALLEL DO   &
1420   !$OMP PRIVATE ( ij )
1421   DO ij = 1 , grid%num_tiles
1422
1423     IF( (config_flags%specified .or. config_flags%nested) .and. rk_step == 1 ) THEN
1424
1425       CALL relax_bdy_dry ( config_flags,                                &
1426                            grid%em_u_save, grid%em_v_save, ph_save, grid%em_t_save,             &
1427                            w_save, mu_tend,                             &
1428                            grid%em_ru, grid%em_rv, grid%em_ph_2, grid%em_t_2,                           &
1429                            grid%em_w_2, grid%em_mu_2, grid%em_mut,                              &
1430                            grid%em_u_b, grid%em_v_b, grid%em_ph_b, grid%em_t_b, grid%em_w_b,                    &
1431                            grid%em_mu_b,                                        &
1432                            grid%em_u_bt, grid%em_v_bt, grid%em_ph_bt, grid%em_t_bt,                     &
1433                            grid%em_w_bt, grid%em_mu_bt,                                 &
1434                            config_flags%spec_bdy_width, grid%spec_zone, grid%relax_zone,       &
1435                            grid%dtbc, grid%fcx, grid%gcx,                              &
1436                            ijds, ijde,                                  &
1437                            ids,ide, jds,jde, kds,kde,                   &
1438                            ims,ime, jms,jme, kms,kme,                   &
1439                            ips,ipe, jps,jpe, kps,kpe,                   &
1440                            grid%i_start(ij), grid%i_end(ij),            &
1441                            grid%j_start(ij), grid%j_end(ij),            &
1442                            k_start, k_end                              )
1443
1444
1445     ENDIF
1446
1447     CALL rk_addtend_dry( grid%em_ru_tend,  grid%em_rv_tend,  rw_tend,  ph_tend,  t_tend,  &
1448                          ru_tendf, rv_tendf, rw_tendf, ph_tendf, t_tendf, &
1449                          grid%em_u_save, grid%em_v_save, w_save, ph_save, grid%em_t_save, &
1450                          mu_tend, mu_tendf, rk_step,                      &
1451                          grid%h_diabatic, grid%em_mut, grid%msft, grid%msfu, grid%msfv,               &
1452                          ids,ide, jds,jde, kds,kde,                       &
1453                          ims,ime, jms,jme, kms,kme,                       &
1454                          ips,ipe, jps,jpe, kps,kpe,                       &
1455                          grid%i_start(ij), grid%i_end(ij),                &
1456                          grid%j_start(ij), grid%j_end(ij),                &
1457                          k_start, k_end                                  )
1458
1459     IF( config_flags%specified .or. config_flags%nested ) THEN
1460       CALL spec_bdy_dry ( config_flags,                                    &
1461                           grid%em_ru_tend, grid%em_rv_tend, ph_tend, t_tend,               &
1462                           rw_tend, mu_tend,                                &
1463                           grid%em_u_b, grid%em_v_b, grid%em_ph_b, grid%em_t_b,                             &
1464                           grid%em_w_b, grid%em_mu_b,                                       &
1465                           grid%em_u_bt, grid%em_v_bt, grid%em_ph_bt, grid%em_t_bt,                         &
1466                           grid%em_w_bt, grid%em_mu_bt,                                     &
1467                           config_flags%spec_bdy_width, grid%spec_zone,                       &
1468                           ijds, ijde,                 & ! min/max(id,jd)
1469                           ids,ide, jds,jde, kds,kde,  & ! domain dims
1470                           ims,ime, jms,jme, kms,kme,  & ! memory dims
1471                           ips,ipe, jps,jpe, kps,kpe,  & ! patch  dims
1472                           grid%i_start(ij), grid%i_end(ij),                &
1473                           grid%j_start(ij), grid%j_end(ij),                &
1474                           k_start, k_end                                  )
1475     
1476     ENDIF
1477
1478   END DO
1479   !$OMP END PARALLEL DO
1480BENCH_END(relax_bdy_dry_tim)
1481
1482!<DESCRIPTION>
1483!<pre>
1484! (3) Small (acoustic,sound) steps.
1485!
1486!    Several acoustic steps are taken each RK pass.  A small step
1487!    sequence begins with calculating perturbation variables
1488!    and coupling them to the column dry-air-mass mu
1489!    (call to small_step_prep).  This is followed by computing
1490!    coefficients for the vertically implicit part of the
1491!    small timestep (call to calc_coef_w). 
1492!
1493!    The small steps are taken
1494!    in the named loop "small_steps:".  In the small_steps loop, first
1495!    the horizontal momentum (u and v) are advanced (call to advance_uv),
1496!    next mu and theta are advanced (call to advance_mu_t) followed by
1497!    advancing w and the geopotential (call to advance_w).  Diagnostic
1498!    values for pressure and inverse density are updated at the end of
1499!    each small_step.
1500!
1501!    The small-step section ends with the change of the perturbation variables
1502!    back to full variables (call to small_step_finish).
1503!</pre>
1504!</DESCRIPTION>
1505
1506BENCH_START(small_step_prep_tim)
1507   !$OMP PARALLEL DO   &
1508   !$OMP PRIVATE ( ij )
1509
1510   DO ij = 1 , grid%num_tiles
1511
1512    ! Calculate coefficients for the vertically implicit acoustic/gravity wave
1513    ! integration.  We only need calculate these for the first pass through -
1514    ! the predictor step.  They are reused as is for the corrector step.
1515    ! For third-order RK, we need to recompute these after the first
1516    ! predictor because we may have changed the small timestep -> grid%dts.
1517
1518      CALL wrf_debug ( 200 , ' call calc_coef_w' )
1519
1520      CALL small_step_prep( grid%em_u_1,grid%em_u_2,grid%em_v_1,grid%em_v_2,grid%em_w_1,grid%em_w_2,          &
1521                            grid%em_t_1,grid%em_t_2,grid%em_ph_1,grid%em_ph_2,                &
1522                            grid%em_mub, grid%em_mu_1, grid%em_mu_2,                  &
1523                            grid%em_muu, muus, grid%em_muv, muvs,             &
1524                            grid%em_mut, grid%em_muts, grid%em_mudf,                  &
1525                            grid%em_u_save, grid%em_v_save, w_save,           &
1526                            grid%em_t_save, ph_save, mu_save,         &
1527                            grid%em_ww, ww1,                          &
1528                            grid%em_dnw, c2a, grid%em_pb, grid%em_p, grid%em_alt,             &
1529                            grid%msfu, grid%msfv, grid%msft,                 &
1530                            rk_step,                          &
1531                            ids, ide, jds, jde, kds, kde,     &
1532                            ims, ime, jms, jme, kms, kme,     &
1533                            grid%i_start(ij), grid%i_end(ij), &
1534                            grid%j_start(ij), grid%j_end(ij), &
1535                            k_start    , k_end               )
1536      CALL calc_p_rho( grid%em_al, grid%em_p, grid%em_ph_2,                      &
1537                       grid%em_alt, grid%em_t_2, grid%em_t_save, c2a, pm1,       &
1538                       grid%em_mu_2, grid%em_muts, grid%em_znu, t0,              &
1539                       grid%em_rdnw, grid%em_dnw, grid%smdiv,                 &
1540                       config_flags%non_hydrostatic, 0,               &
1541                       ids, ide, jds, jde, kds, kde,     &
1542                       ims, ime, jms, jme, kms, kme,     &
1543                       grid%i_start(ij), grid%i_end(ij), &
1544                       grid%j_start(ij), grid%j_end(ij), &
1545                       k_start    , k_end               )
1546
1547      IF (config_flags%non_hydrostatic)                                &
1548      CALL calc_coef_w( a,alpha,gamma,                    &
1549                        grid%em_mut, cqw,                         &
1550                        grid%em_rdn, grid%em_rdnw, c2a,                   &
1551                        dts_rk, g, grid%epssm,                 &
1552                        ids, ide, jds, jde, kds, kde,     &
1553                        ims, ime, jms, jme, kms, kme,     &
1554                        grid%i_start(ij), grid%i_end(ij), &
1555                        grid%j_start(ij), grid%j_end(ij), &
1556                        k_start    , k_end               )
1557
1558
1559   ENDDO
1560   !$OMP END PARALLEL DO
1561BENCH_END(small_step_prep_tim)
1562
1563
1564#ifdef DM_PARALLEL
1565!-----------------------------------------------------------------------
1566!  Stencils for patch communications  (WCS, 29 June 2001)
1567!  Note:  the small size of this halo exchange reflects the
1568!         fact that we are carrying the uncoupled variables
1569!         as state variables in the mass coordinate model, as
1570!         opposed to the coupled variables as in the height
1571!         coordinate model.
1572!
1573!                              * * * * *
1574!            *        * * *    * * * * *
1575!          * + *      * + *    * * + * *
1576!            *        * * *    * * * * *
1577!                              * * * * *
1578!
1579!  3D variables - note staggering!  grid%em_ph_2(grid%em_z), grid%em_u_save(X), grid%em_v_save(Y)
1580!
1581!j grid%em_ph_2      x
1582!j grid%em_al        x
1583!j grid%em_p         x
1584!j grid%em_t_1       x
1585!j grid%em_t_save    x
1586!j grid%em_u_save    x
1587!j grid%em_v_save    x
1588!
1589!  the following are 2D (xy) variables
1590!
1591!j grid%em_mu_1      x
1592!j grid%em_mu_2      x
1593!j grid%em_mudf      x
1594!--------------------------------------------------------------
1595#      include "HALO_EM_B.inc"
1596#      include "PERIOD_BDY_EM_B.inc"
1597#endif
1598
1599BENCH_START(set_phys_bc2_tim)
1600   !$OMP PARALLEL DO   &
1601   !$OMP PRIVATE ( ij )
1602
1603   DO ij = 1 , grid%num_tiles
1604
1605         CALL set_physical_bc3d( grid%em_ru_tend, 'u', config_flags,          &
1606                                 ids, ide, jds, jde, kds, kde, &
1607                                 ims, ime, jms, jme, kms, kme, &
1608                                 ips, ipe, jps, jpe, kps, kpe, &
1609                           grid%i_start(ij), grid%i_end(ij),                 &
1610                           grid%j_start(ij), grid%j_end(ij),                 &
1611                           k_start    , k_end                     )
1612
1613         CALL set_physical_bc3d( grid%em_rv_tend, 'v', config_flags,            &
1614                                 ids, ide, jds, jde, kds, kde, &
1615                                 ims, ime, jms, jme, kms, kme, &
1616                                 ips, ipe, jps, jpe, kps, kpe, &
1617                           grid%i_start(ij), grid%i_end(ij),                 &
1618                           grid%j_start(ij), grid%j_end(ij),                 &
1619                           k_start    , k_end                     )
1620
1621         CALL set_physical_bc3d( grid%em_ph_2, 'w', config_flags,          &
1622                                 ids, ide, jds, jde, kds, kde, &
1623                                 ims, ime, jms, jme, kms, kme, &
1624                                 ips, ipe, jps, jpe, kps, kpe, &
1625                           grid%i_start(ij), grid%i_end(ij),                 &
1626                           grid%j_start(ij), grid%j_end(ij),                 &
1627                           k_start    , k_end                     )
1628
1629         CALL set_physical_bc3d( grid%em_al, 'p', config_flags,            &
1630                                 ids, ide, jds, jde, kds, kde, &
1631                                 ims, ime, jms, jme, kms, kme, &
1632                                 ips, ipe, jps, jpe, kps, kpe, &
1633                           grid%i_start(ij), grid%i_end(ij),                 &
1634                           grid%j_start(ij), grid%j_end(ij),                 &
1635                           k_start    , k_end                     )
1636
1637         CALL set_physical_bc3d( grid%em_p, 'p', config_flags,             &
1638                                 ids, ide, jds, jde, kds, kde, &
1639                                 ims, ime, jms, jme, kms, kme, &
1640                                 ips, ipe, jps, jpe, kps, kpe, &
1641                           grid%i_start(ij), grid%i_end(ij),                 &
1642                           grid%j_start(ij), grid%j_end(ij),                 &
1643                           k_start    , k_end                     )
1644
1645         CALL set_physical_bc3d( grid%em_t_1, 'p', config_flags,             &
1646                                 ids, ide, jds, jde, kds, kde, &
1647                                 ims, ime, jms, jme, kms, kme, &
1648                                 ips, ipe, jps, jpe, kps, kpe, &
1649                           grid%i_start(ij), grid%i_end(ij),                 &
1650                           grid%j_start(ij), grid%j_end(ij),                 &
1651                           k_start    , k_end                     )
1652
1653         CALL set_physical_bc3d( grid%em_t_save, 't', config_flags,             &
1654                                 ids, ide, jds, jde, kds, kde, &
1655                                 ims, ime, jms, jme, kms, kme, &
1656                                 ips, ipe, jps, jpe, kps, kpe, &
1657                           grid%i_start(ij), grid%i_end(ij),                 &
1658                           grid%j_start(ij), grid%j_end(ij),                 &
1659                           k_start    , k_end                     )
1660
1661         CALL set_physical_bc2d( grid%em_mu_1, 't', config_flags,          &
1662                                 ids, ide, jds, jde,               &
1663                                 ims, ime, jms, jme,               &
1664                                 ips, ipe, jps, jpe,               &
1665                                 grid%i_start(ij), grid%i_end(ij), &
1666                                 grid%j_start(ij), grid%j_end(ij) )
1667
1668         CALL set_physical_bc2d( grid%em_mu_2, 't', config_flags,          &
1669                                 ids, ide, jds, jde,               &
1670                                 ims, ime, jms, jme,               &
1671                                 ips, ipe, jps, jpe,               &
1672                                 grid%i_start(ij), grid%i_end(ij), &
1673                                 grid%j_start(ij), grid%j_end(ij) )
1674
1675         CALL set_physical_bc2d( grid%em_mudf, 't', config_flags,          &
1676                                 ids, ide, jds, jde,               &
1677                                 ims, ime, jms, jme,               &
1678                                 ips, ipe, jps, jpe,               &
1679                                 grid%i_start(ij), grid%i_end(ij), &
1680                                 grid%j_start(ij), grid%j_end(ij) )
1681
1682    END DO
1683    !$OMP END PARALLEL DO
1684BENCH_END(set_phys_bc2_tim)
1685
1686   small_steps : DO iteration = 1 , number_of_small_timesteps
1687
1688   ! Boundary condition time (or communication time). 
1689
1690#ifdef DM_PARALLEL
1691#      include "PERIOD_BDY_EM_B.inc"
1692#endif
1693
1694
1695      !$OMP PARALLEL DO   &
1696      !$OMP PRIVATE ( ij )
1697
1698      DO ij = 1 , grid%num_tiles
1699
1700BENCH_START(advance_uv_tim)
1701         CALL advance_uv ( grid%em_u_2, grid%em_ru_tend, grid%em_v_2, grid%em_rv_tend,       &
1702                           grid%em_p, grid%em_pb,                            &
1703                           grid%em_ph_2, grid%em_php, grid%em_alt, grid%em_al, grid%em_mu_2,         &
1704                           grid%em_muu, cqu, grid%em_muv, cqv, grid%em_mudf,         &
1705                           grid%rdx, grid%rdy, dts_rk,                 &
1706                           grid%cf1, grid%cf2, grid%cf3, grid%em_fnm, grid%em_fnp,          &
1707                           grid%emdiv,                            &
1708                           grid%em_rdnw, config_flags,grid%spec_zone,     &
1709                           config_flags%non_hydrostatic,                  &
1710                           ids, ide, jds, jde, kds, kde,     &
1711                           ims, ime, jms, jme, kms, kme,     &
1712                           grid%i_start(ij), grid%i_end(ij), &
1713                           grid%j_start(ij), grid%j_end(ij), &
1714                           k_start    , k_end               )
1715BENCH_END(advance_uv_tim)
1716
1717BENCH_START(spec_bdy_uv_tim)
1718         IF( config_flags%specified .or. config_flags%nested ) THEN
1719           CALL spec_bdyupdate(grid%em_u_2, grid%em_ru_tend, dts_rk,      &
1720                               'u'         , config_flags, &
1721                               grid%spec_zone,                  &
1722                               ids,ide, jds,jde, kds,kde,  & ! domain dims
1723                               ims,ime, jms,jme, kms,kme,  & ! memory dims
1724                               ips,ipe, jps,jpe, kps,kpe,  & ! patch  dims
1725                               grid%i_start(ij), grid%i_end(ij),         &
1726                               grid%j_start(ij), grid%j_end(ij),         &
1727                               k_start    , k_end             )
1728
1729           CALL spec_bdyupdate(grid%em_v_2, grid%em_rv_tend, dts_rk,      &
1730                               'v'         , config_flags, &
1731                               grid%spec_zone,                  &
1732                               ids,ide, jds,jde, kds,kde,  & ! domain dims
1733                               ims,ime, jms,jme, kms,kme,  & ! memory dims
1734                               ips,ipe, jps,jpe, kps,kpe,  & ! patch  dims
1735                               grid%i_start(ij), grid%i_end(ij),         &
1736                               grid%j_start(ij), grid%j_end(ij),         &
1737                               k_start    , k_end             )
1738
1739         ENDIF
1740BENCH_END(spec_bdy_uv_tim)
1741
1742      END DO
1743      !$OMP END PARALLEL DO
1744
1745#ifdef DM_PARALLEL
1746!
1747!  Stencils for patch communications  (WCS, 29 June 2001)
1748!
1749!         *                     *
1750!       * + *      * + *        +
1751!         *                     *
1752!
1753!  grid%em_u_2               x
1754!  grid%em_v_2                          x
1755!
1756#     include "HALO_EM_C.inc"
1757#endif
1758
1759      !$OMP PARALLEL DO   &
1760      !$OMP PRIVATE ( ij )
1761
1762      DO ij = 1 , grid%num_tiles
1763
1764        !  advance the mass in the column, theta, and calculate grid%em_ww
1765
1766BENCH_START(advance_mu_t_tim)
1767        CALL advance_mu_t( grid%em_ww, ww1, grid%em_u_2, grid%em_u_save, grid%em_v_2, grid%em_v_save, &
1768                           grid%em_mu_2, grid%em_mut, muave, grid%em_muts, grid%em_muu, grid%em_muv,  &
1769                           grid%em_mudf, grid%em_ru_m, grid%em_rv_m, grid%em_ww_m,                       &
1770                           grid%em_t_2, grid%em_t_save, t_2save, t_tend,                              &
1771                           mu_tend,                                                                   &
1772                           grid%rdx, grid%rdy, dts_rk, grid%epssm,                                    &
1773                           grid%em_dnw, grid%em_fnm, grid%em_fnp, grid%em_rdnw,                       &
1774                           grid%msfu, grid%msfv, grid%msft,                                           &
1775                           iteration, config_flags,                                                   &
1776                           ids, ide, jds, jde, kds, kde,      &
1777                           ims, ime, jms, jme, kms, kme,      &
1778                           grid%i_start(ij), grid%i_end(ij),  &
1779                           grid%j_start(ij), grid%j_end(ij),  &
1780                           k_start    , k_end                )
1781BENCH_END(advance_mu_t_tim)
1782
1783BENCH_START(spec_bdy_t_tim)
1784         IF( config_flags%specified .or. config_flags%nested ) THEN
1785
1786           CALL spec_bdyupdate(grid%em_t_2, t_tend, dts_rk,      &
1787                               't'         , config_flags, &
1788                               grid%spec_zone,                  &
1789                               ids,ide, jds,jde, kds,kde,  & ! domain dims
1790                               ims,ime, jms,jme, kms,kme,  & ! memory dims
1791                               ips,ipe, jps,jpe, kps,kpe,  & ! patch  dims
1792                               grid%i_start(ij), grid%i_end(ij),         &
1793                               grid%j_start(ij), grid%j_end(ij),         &
1794                               k_start    , k_end             )
1795
1796           CALL spec_bdyupdate(grid%em_mu_2, mu_tend, dts_rk,      &
1797                               'm'         , config_flags, &
1798                               grid%spec_zone,                  &
1799                               ids,ide, jds,jde, 1  ,1  ,  & ! domain dims
1800                               ims,ime, jms,jme, 1  ,1  ,  & ! memory dims
1801                               ips,ipe, jps,jpe, 1  ,1  ,  & ! patch  dims
1802                               grid%i_start(ij), grid%i_end(ij),         &
1803                               grid%j_start(ij), grid%j_end(ij),         &
1804                               1    , 1             )
1805
1806           CALL spec_bdyupdate(grid%em_muts, mu_tend, dts_rk,      &
1807                               'm'         , config_flags, &
1808                               grid%spec_zone,                  &
1809                               ids,ide, jds,jde, 1  ,1  ,  & ! domain dims
1810                               ims,ime, jms,jme, 1  ,1  ,  & ! memory dims
1811                               ips,ipe, jps,jpe, 1  ,1  ,  & ! patch  dims
1812                               grid%i_start(ij), grid%i_end(ij),         &
1813                               grid%j_start(ij), grid%j_end(ij),         &
1814                               1    , 1             )
1815         ENDIF
1816BENCH_END(spec_bdy_t_tim)
1817
1818         ! sumflux accumulates the time-averged mass flux
1819         ! (time averaged over the acoustic steps) for use
1820         ! in the scalar advection (flux divergence).  Using
1821         ! time averaged values gives us exact scalar conservation.
1822
1823BENCH_START(sumflux_tim)
1824         CALL sumflux ( grid%em_u_2, grid%em_v_2, grid%em_ww,                         &
1825                        grid%em_u_save, grid%em_v_save, ww1,                  &
1826                        grid%em_muu, grid%em_muv,                             &
1827                        grid%em_ru_m, grid%em_rv_m, grid%em_ww_m, grid%epssm,              &
1828                        grid%msfu, grid%msfv,                           &
1829                        iteration, number_of_small_timesteps, &
1830                        ids, ide, jds, jde, kds, kde,         &
1831                        ims, ime, jms, jme, kms, kme,         &
1832                        grid%i_start(ij), grid%i_end(ij),     &
1833                        grid%j_start(ij), grid%j_end(ij),     &
1834                        k_start    , k_end                   )
1835BENCH_END(sumflux_tim)
1836
1837         ! small (acoustic) step for the vertical momentum,
1838         ! density and coupled potential temperature.
1839
1840
1841BENCH_START(advance_w_tim)
1842        IF ( config_flags%non_hydrostatic ) THEN
1843          CALL advance_w( grid%em_w_2, rw_tend, grid%em_ww, grid%em_u_2, grid%em_v_2,       &
1844                          grid%em_mu_2, grid%em_mut, muave, grid%em_muts,           &
1845                          t_2save, grid%em_t_2, grid%em_t_save,             &
1846                          grid%em_ph_2, ph_save, grid%em_phb, ph_tend,      &
1847                          grid%ht, c2a, cqw, grid%em_alt, grid%em_alb,           &
1848                          a, alpha, gamma,                  &
1849                          grid%rdx, grid%rdy, dts_rk, t0, grid%epssm,      &
1850                          grid%em_dnw, grid%em_fnm, grid%em_fnp, grid%em_rdnw, grid%em_rdn,         &
1851                          grid%cf1, grid%cf2, grid%cf3, grid%msft,              &
1852                          config_flags,                     &
1853                          ids,ide, jds,jde, kds,kde,        & ! domain dims
1854                          ims,ime, jms,jme, kms,kme,        & ! memory dims
1855                          grid%i_start(ij), grid%i_end(ij), &
1856                          grid%j_start(ij), grid%j_end(ij), &
1857                          k_start    , k_end               )
1858        ENDIF
1859BENCH_END(advance_w_tim)
1860
1861        IF( config_flags%specified .or. config_flags%nested ) THEN
1862
1863BENCH_START(spec_bdynhyd_tim)
1864           IF (config_flags%non_hydrostatic)  THEN
1865             CALL spec_bdyupdate_ph( ph_save, grid%em_ph_2, ph_tend, mu_tend, grid%em_muts, dts_rk, &
1866                                     'h'         , config_flags, &
1867                                     grid%spec_zone,                  &
1868                                     ids,ide, jds,jde, kds,kde,  & ! domain dims
1869                                     ims,ime, jms,jme, kms,kme,  & ! memory dims
1870                                     ips,ipe, jps,jpe, kps,kpe,  & ! patch  dims
1871                                     grid%i_start(ij), grid%i_end(ij),         &
1872                                     grid%j_start(ij), grid%j_end(ij),         &
1873                                     k_start    , k_end             )
1874
1875             IF( config_flags%specified ) THEN
1876               CALL zero_grad_bdy ( grid%em_w_2,                        &
1877                                    'w'         , config_flags, &
1878                                    grid%spec_zone,                  &
1879                                    ids,ide, jds,jde, kds,kde,  & ! domain dims
1880                                    ims,ime, jms,jme, kms,kme,  & ! memory dims
1881                                    ips,ipe, jps,jpe, kps,kpe,  & ! patch  dims
1882                                    grid%i_start(ij), grid%i_end(ij),         &
1883                                    grid%j_start(ij), grid%j_end(ij),         &
1884                                    k_start    , k_end             )
1885             ELSE
1886               CALL spec_bdyupdate   ( grid%em_w_2, rw_tend, dts_rk,       &
1887                                       'h'         , config_flags, &
1888                                       grid%spec_zone,                  &
1889                                       ids,ide, jds,jde, kds,kde,  & ! domain dims
1890                                       ims,ime, jms,jme, kms,kme,  & ! memory dims
1891                                       ips,ipe, jps,jpe, kps,kpe,  & ! patch  dims
1892                                       grid%i_start(ij), grid%i_end(ij),         &
1893                                       grid%j_start(ij), grid%j_end(ij),         &
1894                                       k_start    , k_end             )
1895             ENDIF
1896          ENDIF
1897BENCH_END(spec_bdynhyd_tim)
1898        ENDIF
1899
1900BENCH_START(cald_p_rho_tim)
1901        CALL calc_p_rho( grid%em_al, grid%em_p, grid%em_ph_2,                      &
1902                         grid%em_alt, grid%em_t_2, grid%em_t_save, c2a, pm1,       &
1903                         grid%em_mu_2, grid%em_muts, grid%em_znu, t0,              &
1904                         grid%em_rdnw, grid%em_dnw, grid%smdiv,                 &
1905                         config_flags%non_hydrostatic, iteration,       &
1906                         ids, ide, jds, jde, kds, kde,     &
1907                         ims, ime, jms, jme, kms, kme,     &
1908                         grid%i_start(ij), grid%i_end(ij), &
1909                         grid%j_start(ij), grid%j_end(ij), &
1910                         k_start    , k_end               )
1911BENCH_END(cald_p_rho_tim)
1912
1913   ENDDO
1914   !$OMP END PARALLEL DO
1915
1916#ifdef DM_PARALLEL
1917!
1918!  Stencils for patch communications  (WCS, 29 June 2001)
1919!
1920!         *                     *
1921!       * + *      * + *        +
1922!         *                     *
1923!
1924!  grid%em_ph_2   x
1925!  grid%em_al     x
1926!  grid%em_p      x
1927!
1928!  2D variables (x,y)
1929!
1930!  grid%em_mu_2   x
1931!  grid%em_muts   x
1932!  grid%em_mudf   x
1933
1934#      include "HALO_EM_C2.inc"
1935#      include "PERIOD_BDY_EM_B3.inc"
1936#endif
1937
1938BENCH_START(phys_bc_tim)
1939      !$OMP PARALLEL DO   &
1940      !$OMP PRIVATE ( ij )
1941
1942      DO ij = 1 , grid%num_tiles
1943
1944        ! boundary condition set for next small timestep
1945
1946         CALL set_physical_bc3d( grid%em_ph_2, 'w', config_flags,          &
1947                                 ids, ide, jds, jde, kds, kde,     &
1948                                 ims, ime, jms, jme, kms, kme,     &
1949                                 ips, ipe, jps, jpe, kps, kpe,     &
1950                                 grid%i_start(ij), grid%i_end(ij), &
1951                                 grid%j_start(ij), grid%j_end(ij), &
1952                                 k_start    , k_end               )
1953
1954         CALL set_physical_bc3d( grid%em_al, 'p', config_flags,            &
1955                                 ids, ide, jds, jde, kds, kde,     &
1956                                 ims, ime, jms, jme, kms, kme,     &
1957                                 ips, ipe, jps, jpe, kps, kpe,     &
1958                                 grid%i_start(ij), grid%i_end(ij), &
1959                                 grid%j_start(ij), grid%j_end(ij), &
1960                                 k_start    , k_end               )
1961
1962         CALL set_physical_bc3d( grid%em_p, 'p', config_flags,             &
1963                                 ids, ide, jds, jde, kds, kde,     &
1964                                 ims, ime, jms, jme, kms, kme,     &
1965                                 ips, ipe, jps, jpe, kps, kpe,     &
1966                                 grid%i_start(ij), grid%i_end(ij), &
1967                                 grid%j_start(ij), grid%j_end(ij), &
1968                                 k_start    , k_end               )
1969
1970         CALL set_physical_bc2d( grid%em_muts, 't', config_flags,          &
1971                                 ids, ide, jds, jde,               &
1972                                 ims, ime, jms, jme,               &
1973                                 ips, ipe, jps, jpe,               &
1974                                 grid%i_start(ij), grid%i_end(ij), &
1975                                 grid%j_start(ij), grid%j_end(ij) )
1976
1977         CALL set_physical_bc2d( grid%em_mu_2, 't', config_flags,          &
1978                                 ids, ide, jds, jde,               &
1979                                 ims, ime, jms, jme,               &
1980                                 ips, ipe, jps, jpe,               &
1981                                 grid%i_start(ij), grid%i_end(ij), &
1982                                 grid%j_start(ij), grid%j_end(ij) )
1983
1984         CALL set_physical_bc2d( grid%em_mudf, 't', config_flags,          &
1985                                 ids, ide, jds, jde,               &
1986                                 ims, ime, jms, jme,               &
1987                                 ips, ipe, jps, jpe,               &
1988                                 grid%i_start(ij), grid%i_end(ij), &
1989                                 grid%j_start(ij), grid%j_end(ij) )
1990
1991      END DO
1992      !$OMP END PARALLEL DO
1993BENCH_END(phys_bc_tim)
1994
1995   END DO small_steps
1996
1997   !$OMP PARALLEL DO   &
1998   !$OMP PRIVATE ( ij )
1999
2000   DO ij = 1 , grid%num_tiles
2001
2002      CALL wrf_debug ( 200 , ' call rk_small_finish' )
2003
2004      ! change time-perturbation variables back to
2005      ! full perturbation variables.
2006      ! first get updated mu at u and v points
2007
2008BENCH_START(calc_mu_uv_tim)
2009      CALL calc_mu_uv_1 ( config_flags,                     &
2010                          grid%em_muts, muus, muvs,                 &
2011                          ids, ide, jds, jde, kds, kde,     &
2012                          ims, ime, jms, jme, kms, kme,     &
2013                          grid%i_start(ij), grid%i_end(ij), &
2014                          grid%j_start(ij), grid%j_end(ij), &
2015                          k_start    , k_end               )
2016BENCH_END(calc_mu_uv_tim)
2017BENCH_START(small_step_finish_tim)
2018      CALL small_step_finish( grid%em_u_2, grid%em_u_1, grid%em_v_2, grid%em_v_1, grid%em_w_2, grid%em_w_1,     &
2019                              grid%em_t_2, grid%em_t_1, grid%em_ph_2, grid%em_ph_1, grid%em_ww, ww1,    &
2020                              grid%em_mu_2, grid%em_mu_1,                       &
2021                              grid%em_mut, grid%em_muts, grid%em_muu, muus, grid%em_muv, muvs,  &
2022                              grid%em_u_save, grid%em_v_save, w_save,           &
2023                              grid%em_t_save, ph_save, mu_save,         &
2024                              grid%msfu, grid%msfv, grid%msft,                 &
2025                              grid%h_diabatic,                       &
2026                              number_of_small_timesteps,dts_rk, &
2027                              rk_step, rk_order,                &
2028                              ids, ide, jds, jde, kds, kde,     &
2029                              ims, ime, jms, jme, kms, kme,     &
2030                              grid%i_start(ij), grid%i_end(ij), &
2031                              grid%j_start(ij), grid%j_end(ij), &
2032                              k_start    , k_end               )
2033!  call  to set ru_m, rv_m and ww_m b.c's for PD advection
2034
2035         IF (rk_step == 3) THEN
2036
2037           CALL set_physical_bc3d( grid%em_ru_m, 'u', config_flags,   &
2038                                   ids, ide, jds, jde, kds, kde,      &
2039                                   ims, ime, jms, jme, kms, kme,      &
2040                                   ips, ipe, jps, jpe, kps, kpe,      &
2041                                   grid%i_start(ij), grid%i_end(ij),  &
2042                                   grid%j_start(ij), grid%j_end(ij),  &
2043                                   k_start    , k_end                )
2044
2045           CALL set_physical_bc3d( grid%em_rv_m, 'v', config_flags,   &
2046                                   ids, ide, jds, jde, kds, kde,      &
2047                                   ims, ime, jms, jme, kms, kme,      &
2048                                   ips, ipe, jps, jpe, kps, kpe,      &
2049                                   grid%i_start(ij), grid%i_end(ij),  &
2050                                   grid%j_start(ij), grid%j_end(ij),  &
2051                                   k_start    , k_end                )
2052
2053           CALL set_physical_bc3d( grid%em_ww_m, 'w', config_flags,   &
2054                                   ids, ide, jds, jde, kds, kde,      &
2055                                   ims, ime, jms, jme, kms, kme,      &
2056                                   ips, ipe, jps, jpe, kps, kpe,      &
2057                                   grid%i_start(ij), grid%i_end(ij),  &
2058                                   grid%j_start(ij), grid%j_end(ij),  &
2059                                   k_start    , k_end                )
2060
2061           CALL set_physical_bc2d( grid%em_mut, 't', config_flags,   &
2062                                   ids, ide, jds, jde,               &
2063                                   ims, ime, jms, jme,                &
2064                                   ips, ipe, jps, jpe,                &
2065                                   grid%i_start(ij), grid%i_end(ij),  &
2066                                   grid%j_start(ij), grid%j_end(ij) )
2067
2068          END IF
2069
2070BENCH_END(small_step_finish_tim)
2071
2072   END DO
2073   !$OMP END PARALLEL DO
2074
2075!-----------------------------------------------------------------------
2076!  add in physics tendency first if positive definite advection is used.
2077!  pd advection applies advective flux limiter on last runge-kutta step
2078!-----------------------------------------------------------------------
2079! first moisture
2080
2081  IF (config_flags%pd_moist .and. (rk_step == rk_order)) THEN
2082
2083   !$OMP PARALLEL DO   &
2084   !$OMP PRIVATE ( ij )
2085   DO ij = 1 , grid%num_tiles
2086       CALL wrf_debug ( 200 , ' call rk_update_scalar_pd' )
2087       do im = PARAM_FIRST_SCALAR, num_3d_m
2088       CALL rk_update_scalar_pd( im, im,                                   &
2089                                 moist_old(ims,kms,jms,im),                &
2090                                 moist_tend(ims,kms,jms,im),               &
2091                                 grid%em_mu_1, grid%em_mu_1, grid%em_mub,  &
2092                                 rk_step, dt_rk, grid%spec_zone,           &
2093                                 config_flags,                             &
2094                                 ids, ide, jds, jde, kds, kde,             &
2095                                 ims, ime, jms, jme, kms, kme,             &
2096                                 grid%i_start(ij), grid%i_end(ij),         &
2097                                 grid%j_start(ij), grid%j_end(ij),         &
2098                                 k_start    , k_end                       )
2099       ENDDO
2100   END DO
2101   !$OMP END PARALLEL DO
2102
2103!---------------------- positive definite bc call
2104
2105   !$OMP PARALLEL DO   &
2106   !$OMP PRIVATE ( ij )
2107   DO ij = 1 , grid%num_tiles
2108      IF (num_3d_m >= PARAM_FIRST_SCALAR) THEN
2109        DO im = PARAM_FIRST_SCALAR , num_3d_m
2110          CALL set_physical_bc3d( moist_old(ims,kms,jms,im), 'p', config_flags,   &
2111                                   ids, ide, jds, jde, kds, kde,                  &
2112                                   ims, ime, jms, jme, kms, kme,                  &
2113                                   ips, ipe, jps, jpe, kps, kpe,                  &
2114                                   grid%i_start(ij), grid%i_end(ij),              &
2115                                   grid%j_start(ij), grid%j_end(ij),              &
2116                                   k_start    , k_end                            )
2117         END DO
2118      ENDIF
2119   END DO
2120   !$OMP END PARALLEL DO
2121
2122#ifdef DM_PARALLEL
2123   if(config_flags%pd_moist) then
2124#ifndef RSL
2125     IF      ( config_flags%h_sca_adv_order <= 4 ) THEN
2126#     include "HALO_EM_MOIST_OLD_E_5.inc"
2127     ELSE IF ( config_flags%h_sca_adv_order <= 6 ) THEN
2128#     include "HALO_EM_MOIST_OLD_E_7.inc"
2129     ELSE
2130       WRITE(wrf_err_message,*)'solve_em: invalid h_sca_adv_order = ',config_flags%h_sca_adv_order
2131       CALL wrf_error_fatal(TRIM(wrf_err_message))
2132     ENDIF
2133#else
2134       WRITE(wrf_err_message,*)'cannot use pd scheme with RSL - use RSL-LITE'
2135       CALL wrf_error_fatal(TRIM(wrf_err_message))
2136#endif     
2137  endif
2138#endif
2139
2140   END IF  ! end if for pd_moist
2141
2142! scalars
2143
2144  IF (config_flags%pd_scalar .and. (rk_step == rk_order)) THEN
2145
2146   !$OMP PARALLEL DO   &
2147   !$OMP PRIVATE ( ij )
2148
2149   DO ij = 1 , grid%num_tiles
2150       CALL wrf_debug ( 200 , ' call rk_update_scalar_pd' )
2151       do im = PARAM_FIRST_SCALAR, num_3d_s
2152       CALL rk_update_scalar_pd( im, im,                                  &
2153                                 scalar_old(ims,kms,jms,im),              &
2154                                 scalar_tend(ims,kms,jms,im),             &
2155                                 grid%em_mu_1, grid%em_mu_1, grid%em_mub, &
2156                                 rk_step, dt_rk, grid%spec_zone,          &
2157                                 config_flags,                            &
2158                                 ids, ide, jds, jde, kds, kde,            &
2159                                 ims, ime, jms, jme, kms, kme,            &
2160                                 grid%i_start(ij), grid%i_end(ij),        &
2161                                 grid%j_start(ij), grid%j_end(ij),        &
2162                                 k_start    , k_end                      )
2163       ENDDO
2164   ENDDO
2165   !$OMP END PARALLEL DO
2166
2167!---------------------- positive definite bc call
2168
2169   !$OMP PARALLEL DO   &
2170   !$OMP PRIVATE ( ij )
2171
2172   DO ij = 1 , grid%num_tiles
2173      IF (num_3d_m >= PARAM_FIRST_SCALAR) THEN
2174        DO im = PARAM_FIRST_SCALAR , num_3d_s
2175            CALL set_physical_bc3d(  scalar_old(ims,kms,jms,im), 'p', config_flags, &
2176                                   ids, ide, jds, jde, kds, kde,                    &
2177                                   ims, ime, jms, jme, kms, kme,                    &
2178                                   ips, ipe, jps, jpe, kps, kpe,                    &
2179                                   grid%i_start(ij), grid%i_end(ij),                &
2180                                   grid%j_start(ij), grid%j_end(ij),                &
2181                                   k_start    , k_end                              )
2182         END DO
2183      ENDIF
2184   END DO
2185   !$OMP END PARALLEL DO
2186
2187#ifdef DM_PARALLEL
2188   if(config_flags%pd_scalar) then
2189#ifndef RSL
2190     IF      ( config_flags%h_sca_adv_order <= 4 ) THEN
2191#     include "HALO_EM_SCALAR_OLD_E_5.inc"
2192     ELSE IF ( config_flags%h_sca_adv_order <= 6 ) THEN
2193#     include "HALO_EM_SCALAR_OLD_E_7.inc"
2194     ELSE
2195       WRITE(wrf_err_message,*)'solve_em: invalid h_sca_adv_order = ',config_flags%h_sca_adv_order
2196       CALL wrf_error_fatal(TRIM(wrf_err_message))
2197     ENDIF
2198#else
2199       WRITE(wrf_err_message,*)'cannot use pd scheme with RSL - use RSL-LITE'
2200       CALL wrf_error_fatal(TRIM(wrf_err_message))
2201#endif     
2202  endif
2203#endif
2204
2205   END IF  ! end if for pd_scalar
2206
2207! chem
2208
2209  IF (config_flags%pd_chem .and. (rk_step == rk_order)) THEN
2210
2211   write(6,*) ' pd advection for chem '
2212
2213   !$OMP PARALLEL DO   &
2214   !$OMP PRIVATE ( ij )
2215   DO ij = 1 , grid%num_tiles
2216       CALL wrf_debug ( 200 , ' call rk_update_scalar_pd' )
2217       do im = PARAM_FIRST_SCALAR, num_3d_c
2218       CALL rk_update_scalar_pd( im, im,                                  &
2219                                 chem_old(ims,kms,jms,im),                &
2220                                 chem_tend(ims,kms,jms,im),               &
2221                                 grid%em_mu_1, grid%em_mu_1, grid%em_mub, &
2222                                 rk_step, dt_rk, grid%spec_zone,          &
2223                                 config_flags,                            &
2224                                 ids, ide, jds, jde, kds, kde,            &
2225                                 ims, ime, jms, jme, kms, kme,            &
2226                                 grid%i_start(ij), grid%i_end(ij),        &
2227                                 grid%j_start(ij), grid%j_end(ij),        &
2228                                 k_start    , k_end                      )
2229       ENDDO
2230   END DO
2231   !$OMP END PARALLEL DO
2232
2233!---------------------- positive definite bc call
2234
2235   !$OMP PARALLEL DO   &
2236   !$OMP PRIVATE ( ij )
2237   DO ij = 1 , grid%num_tiles
2238      IF (num_3d_m >= PARAM_FIRST_SCALAR) THEN
2239        DO im = PARAM_FIRST_SCALAR , num_3d_c
2240          CALL set_physical_bc3d(  chem_old(ims,kms,jms,im), 'p', config_flags,     &
2241                                   ids, ide, jds, jde, kds, kde,                    &
2242                                   ims, ime, jms, jme, kms, kme,                    &
2243                                   ips, ipe, jps, jpe, kps, kpe,                    &
2244                                   grid%i_start(ij), grid%i_end(ij),                &
2245                                   grid%j_start(ij), grid%j_end(ij),                &
2246                                   k_start    , k_end                              )
2247         END DO
2248      ENDIF
2249   END DO
2250   !$OMP END PARALLEL DO
2251
2252
2253#ifdef DM_PARALLEL
2254   if(config_flags%pd_chem) then
2255#ifndef RSL
2256     IF      ( config_flags%h_sca_adv_order <= 4 ) THEN
2257#     include "HALO_EM_CHEM_OLD_E_5.inc"
2258     ELSE IF ( config_flags%h_sca_adv_order <= 6 ) THEN
2259#     include "HALO_EM_CHEM_OLD_E_7.inc"
2260     ELSE
2261       WRITE(wrf_err_message,*)'solve_em: invalid h_sca_adv_order = ',config_flags%h_sca_adv_order
2262       CALL wrf_error_fatal(TRIM(wrf_err_message))
2263     ENDIF
2264#else
2265       WRITE(wrf_err_message,*)'cannot use pd scheme with RSL - use RSL-LITE'
2266       CALL wrf_error_fatal(TRIM(wrf_err_message))
2267#endif     
2268  endif
2269#endif
2270
2271  END IF  ! end if for pd_chem
2272
2273! tke
2274
2275  IF (config_flags%pd_tke .and. (rk_step == rk_order) &
2276      .and. (config_flags%km_opt .eq. 2)                ) THEN
2277
2278   write(6,*) ' pd advection for tke '
2279
2280   !$OMP PARALLEL DO   &
2281   !$OMP PRIVATE ( ij )
2282
2283   DO ij = 1 , grid%num_tiles
2284       CALL wrf_debug ( 200 , ' call rk_update_scalar_pd' )
2285       CALL rk_update_scalar_pd( 1, 1,                                    &
2286                                 grid%em_tke_1,                           &
2287                                 tke_tend(ims,kms,jms),                   &
2288                                 grid%em_mu_1, grid%em_mu_1, grid%em_mub, &
2289                                 rk_step, dt_rk, grid%spec_zone,          &
2290                                 config_flags,                            &
2291                                 ids, ide, jds, jde, kds, kde,            &
2292                                 ims, ime, jms, jme, kms, kme,            &
2293                                 grid%i_start(ij), grid%i_end(ij),        &
2294                                 grid%j_start(ij), grid%j_end(ij),        &
2295                                 k_start    , k_end                      )
2296   END DO
2297   !$OMP END PARALLEL DO
2298
2299!---------------------- positive definite bc call
2300
2301   !$OMP PARALLEL DO   &
2302   !$OMP PRIVATE ( ij )
2303   DO ij = 1 , grid%num_tiles
2304          CALL set_physical_bc3d(  grid%em_tke_1, 'p', config_flags,  &
2305                                   ids, ide, jds, jde, kds, kde,      &
2306                                   ims, ime, jms, jme, kms, kme,      &
2307                                   ips, ipe, jps, jpe, kps, kpe,      &
2308                                   grid%i_start(ij), grid%i_end(ij),  &
2309                                   grid%j_start(ij), grid%j_end(ij),  &
2310                                   k_start    , k_end                )
2311   END DO
2312   !$OMP END PARALLEL DO
2313
2314!---  end of positive definite physics tendency update
2315
2316#ifdef DM_PARALLEL
2317     IF      ( config_flags%h_sca_adv_order <= 4 ) THEN
2318#     include "HALO_EM_TKE_OLD_E_5.inc"
2319     ELSE IF ( config_flags%h_sca_adv_order <= 6 ) THEN
2320#     include "HALO_EM_TKE_OLD_E_7.inc"
2321     ELSE
2322       WRITE(wrf_err_message,*)'solve_em: invalid h_sca_adv_order = ',config_flags%h_sca_adv_order
2323       CALL wrf_error_fatal(TRIM(wrf_err_message))
2324     ENDIF
2325#endif
2326
2327   END IF  ! end if for pd_tke
2328
2329#ifdef DM_PARALLEL
2330!
2331!  Stencils for patch communications  (WCS, 29 June 2001)
2332!
2333!
2334! grid%em_ru_m      x
2335! grid%em_rv_m      x
2336! grid%em_ww_m      x
2337! grid%em_mut       x
2338!
2339!--------------------------------------------------------------
2340
2341#  include "HALO_EM_D.inc"
2342#endif
2343
2344!<DESCRIPTION>
2345!<pre>
2346! (4) Still within the RK loop, the scalar variables are advanced.
2347!
2348!    For the moist and chem variables, each one is advanced
2349!    individually, using named loops "moist_variable_loop:"
2350!    and "chem_variable_loop:".  Each RK substep begins by
2351!    calculating the advective tendency, and, for the first RK step,
2352!    3D mixing (calling rk_scalar_tend) followed by an update
2353!    of the scalar (calling rk_scalar_update).
2354!</pre>
2355!</DESCRIPTION>
2356
2357
2358  moist_scalar_advance: IF (num_3d_m >= PARAM_FIRST_SCALAR )  THEN
2359
2360   moist_variable_loop: do im = PARAM_FIRST_SCALAR, num_3d_m
2361
2362   if (grid%adv_moist_cond .or. im==p_qv ) then
2363
2364   !$OMP PARALLEL DO   &
2365   !$OMP PRIVATE ( ij )
2366
2367   moist_tile_loop_1: DO ij = 1 , grid%num_tiles
2368
2369       CALL wrf_debug ( 200 , ' call rk_scalar_tend' )
2370
2371BENCH_START(rk_scalar_tend_tim)
2372
2373       CALL rk_scalar_tend (  im, im, config_flags,                                      &
2374                              rk_step, dt_rk,                                            &
2375                              grid%em_ru_m, grid%em_rv_m, grid%em_ww_m,                  &
2376                              grid%em_mut, grid%em_mub, grid%em_mu_1,                    &
2377                              grid%em_alt,                                               &
2378                              moist_old(ims,kms,jms,im),                                 &
2379                              moist(ims,kms,jms,im),                                     &
2380                              moist_tend(ims,kms,jms,im),                                &
2381                              advect_tend,grid%rqvften,                                  &
2382                              grid%qv_base, .true., grid%em_fnm, grid%em_fnp,            &
2383                              grid%msfu, grid%msfv, grid%msft,                           &
2384                              grid%rdx, grid%rdy, grid%em_rdn, grid%em_rdnw, grid%khdif, &
2385                              grid%kvdif, grid%xkmhd,                                    &
2386                              grid%diff_6th_opt, grid%diff_6th_factor,                   &
2387                              config_flags%pd_moist,            &
2388                              ids, ide, jds, jde, kds, kde,     &
2389                              ims, ime, jms, jme, kms, kme,     &
2390                              grid%i_start(ij), grid%i_end(ij), &
2391                              grid%j_start(ij), grid%j_end(ij), &
2392                              k_start    , k_end               )
2393
2394BENCH_END(rk_scalar_tend_tim)
2395
2396BENCH_START(rlx_bdy_scalar_tim)
2397     IF( ( config_flags%specified .or. config_flags%nested ) .and. rk_step == 1 ) THEN
2398         IF ( im .EQ. P_QV .OR. config_flags%nested ) THEN
2399            CALL relax_bdy_scalar ( moist_tend(ims,kms,jms,im),            &
2400                                    moist(ims,kms,jms,im),  grid%em_mut,         &
2401                                    moist_b(1,1,1,1,im),                   &
2402                                    moist_bt(1,1,1,1,im),                  &
2403                                    config_flags%spec_bdy_width, grid%spec_zone, grid%relax_zone, &
2404                                    grid%dtbc, grid%fcx, grid%gcx,             &
2405                                    config_flags,               &
2406                                    ijds, ijde,                 & ! min/max(id,jd)
2407                                    ids,ide, jds,jde, kds,kde,  & ! domain dims
2408                                    ims,ime, jms,jme, kms,kme,  & ! memory dims
2409                                    ips,ipe, jps,jpe, kps,kpe,  & ! patch  dims
2410                                    grid%i_start(ij), grid%i_end(ij),      &
2411                                    grid%j_start(ij), grid%j_end(ij),      &
2412                                    k_start, k_end                        )
2413
2414            CALL spec_bdy_scalar  ( moist_tend(ims,kms,jms,im),                &
2415                                    moist_b(1,1,1,1,im),                     &
2416                                    moist_bt(1,1,1,1,im),                    &
2417                                    config_flags%spec_bdy_width, grid%spec_zone,                 &
2418                                    config_flags,               &
2419                                    ijds, ijde,                 & ! min/max(id,jd)
2420                                    ids,ide, jds,jde, kds,kde,  & ! domain dims
2421                                    ims,ime, jms,jme, kms,kme,  & ! memory dims
2422                                    ips,ipe, jps,jpe, kps,kpe,  & ! patch  dims
2423                                    grid%i_start(ij), grid%i_end(ij),          &
2424                                    grid%j_start(ij), grid%j_end(ij),          &
2425                                    k_start, k_end                               )
2426          ENDIF
2427     ENDIF
2428BENCH_END(rlx_bdy_scalar_tim)
2429
2430   ENDDO moist_tile_loop_1
2431   !$OMP END PARALLEL DO
2432
2433   !$OMP PARALLEL DO   &
2434   !$OMP PRIVATE ( ij )
2435
2436   moist_tile_loop_2: DO ij = 1 , grid%num_tiles
2437
2438       CALL wrf_debug ( 200 , ' call rk_update_scalar' )
2439
2440BENCH_START(update_scal_tim)
2441       CALL rk_update_scalar( im, im,                           &
2442                              moist_old(ims,kms,jms,im),        &
2443                              moist(ims,kms,jms,im),            &
2444                              moist_tend(ims,kms,jms,im),       &
2445                              advect_tend, grid%msft,                &
2446                              grid%em_mu_1, grid%em_mu_2, grid%em_mub,                  &
2447                              rk_step, dt_rk, grid%spec_zone,        &
2448                              config_flags,     &
2449                              ids, ide, jds, jde, kds, kde,     &
2450                              ims, ime, jms, jme, kms, kme,     &
2451                              grid%i_start(ij), grid%i_end(ij), &
2452                              grid%j_start(ij), grid%j_end(ij), &
2453                              k_start    , k_end               )
2454BENCH_END(update_scal_tim)
2455
2456BENCH_START(flow_depbdy_tim)
2457       IF( config_flags%specified ) THEN
2458         IF(im .ne. P_QV)THEN
2459           CALL flow_dep_bdy  (  moist(ims,kms,jms,im),                     &
2460                               grid%em_ru_m, grid%em_rv_m, config_flags, &
2461                               grid%spec_zone,                  &
2462                               ids,ide, jds,jde, kds,kde,  & ! domain dims
2463                               ims,ime, jms,jme, kms,kme,  & ! memory dims
2464                               ips,ipe, jps,jpe, kps,kpe,  & ! patch  dims
2465                               grid%i_start(ij), grid%i_end(ij),                      &
2466                               grid%j_start(ij), grid%j_end(ij),                      &
2467                               k_start, k_end                               )
2468         ENDIF
2469       ENDIF
2470BENCH_END(flow_depbdy_tim)
2471
2472   ENDDO moist_tile_loop_2
2473   !$OMP END PARALLEL DO
2474
2475     ENDIF  !-- if (grid%adv_moist_cond .or. im==p_qv ) then
2476
2477   ENDDO moist_variable_loop
2478
2479 ENDIF moist_scalar_advance
2480
2481BENCH_START(tke_adv_tim)
2482 TKE_advance: IF (config_flags%km_opt .eq. 2) then
2483
2484#ifdef DM_PARALLEL
2485      IF      ( config_flags%h_mom_adv_order <= 4 ) THEN
2486#       include "HALO_EM_TKE_ADVECT_3.inc"
2487      ELSE IF ( config_flags%h_mom_adv_order <= 6 ) THEN
2488#       include "HALO_EM_TKE_ADVECT_5.inc"
2489      ELSE
2490        WRITE(wrf_err_message,*)'solve_em: invalid h_mom_adv_order = ',config_flags%h_mom_adv_order
2491        CALL wrf_error_fatal(TRIM(wrf_err_message))
2492      ENDIF
2493#endif
2494
2495   !$OMP PARALLEL DO   &
2496   !$OMP PRIVATE ( ij )
2497
2498   tke_tile_loop_1: DO ij = 1 , grid%num_tiles
2499
2500     CALL wrf_debug ( 200 , ' call rk_scalar_tend for tke' )
2501     CALL rk_scalar_tend ( 1, 1, config_flags,                                        &
2502                           rk_step, dt_rk,                                            &
2503                           grid%em_ru_m, grid%em_rv_m, grid%em_ww_m,                  &
2504                           grid%em_mut, grid%em_mub, grid%em_mu_1,                    &
2505                           grid%em_alt,                                               &
2506                           grid%em_tke_1,                                             &
2507                           grid%em_tke_2,                                             &
2508                           tke_tend(ims,kms,jms),                                     &
2509                           advect_tend,grid%rqvften,                                  &
2510                           grid%qv_base, .false., grid%em_fnm, grid%em_fnp,           &
2511                           grid%msfu, grid%msfv, grid%msft,                           &
2512                           grid%rdx, grid%rdy, grid%em_rdn, grid%em_rdnw, grid%khdif, &
2513                           grid%kvdif, grid%xkmhd,                                    &
2514                           grid%diff_6th_opt, grid%diff_6th_factor,                   &
2515                           config_flags%pd_tke,              &
2516                           ids, ide, jds, jde, kds, kde,     &
2517                           ims, ime, jms, jme, kms, kme,     &
2518                           grid%i_start(ij), grid%i_end(ij), &
2519                           grid%j_start(ij), grid%j_end(ij), &
2520                           k_start    , k_end               )
2521
2522   ENDDO tke_tile_loop_1
2523   !$OMP END PARALLEL DO
2524
2525   !$OMP PARALLEL DO   &
2526   !$OMP PRIVATE ( ij )
2527
2528   tke_tile_loop_2: DO ij = 1 , grid%num_tiles
2529
2530     CALL wrf_debug ( 200 , ' call rk_update_scalar' )
2531     CALL rk_update_scalar( 1, 1,                             &
2532                            grid%em_tke_1,               &
2533                            grid%em_tke_2,               &
2534                            tke_tend(ims,kms,jms),            &
2535                            advect_tend,grid%msft,                 &
2536                            grid%em_mu_1, grid%em_mu_2, grid%em_mub,                  &
2537                            rk_step, dt_rk, grid%spec_zone,        &
2538                            config_flags,     &
2539                            ids, ide, jds, jde, kds, kde,     &
2540                            ims, ime, jms, jme, kms, kme,     &
2541                            grid%i_start(ij), grid%i_end(ij), &
2542                            grid%j_start(ij), grid%j_end(ij), &
2543                            k_start    , k_end               )
2544
2545! bound the tke (greater than 0, less than tke_upper_bound)
2546
2547     CALL bound_tke( grid%em_tke_2, grid%tke_upper_bound, &
2548                     ids, ide, jds, jde, kds, kde,        &
2549                     ims, ime, jms, jme, kms, kme,        &
2550                     grid%i_start(ij), grid%i_end(ij),    &
2551                     grid%j_start(ij), grid%j_end(ij),    &
2552                     k_start    , k_end                  )
2553
2554     IF( config_flags%specified .or. config_flags%nested ) THEN
2555         CALL flow_dep_bdy (  grid%em_tke_2,                     &
2556                              grid%em_ru_m, grid%em_rv_m, config_flags,               &
2557                              grid%spec_zone,                              &
2558                              ids,ide, jds,jde, kds,kde,  & ! domain dims
2559                              ims,ime, jms,jme, kms,kme,  & ! memory dims
2560                              ips,ipe, jps,jpe, kps,kpe,  & ! patch  dims
2561                              grid%i_start(ij), grid%i_end(ij),       &
2562                              grid%j_start(ij), grid%j_end(ij),       &
2563                              k_start, k_end                               )
2564     ENDIF
2565   ENDDO tke_tile_loop_2
2566   !$OMP END PARALLEL DO
2567
2568   END IF TKE_advance
2569BENCH_END(tke_adv_tim)
2570
2571#ifdef WRF_CHEM
2572!  next the chemical species
2573BENCH_START(chem_adv_tim)
2574  chem_scalar_advance: IF (num_3d_c >= PARAM_FIRST_SCALAR)  THEN
2575
2576   chem_variable_loop: do ic = PARAM_FIRST_SCALAR, num_3d_c
2577
2578   !$OMP PARALLEL DO   &
2579   !$OMP PRIVATE ( ij )
2580
2581   chem_tile_loop_1: DO ij = 1 , grid%num_tiles
2582
2583       CALL wrf_debug ( 200 , ' call rk_scalar_tend in chem_tile_loop_1' )
2584       CALL rk_scalar_tend ( ic, ic, config_flags,                            &
2585                             rk_step, dt_rk,                                  &
2586                             grid%em_ru_m, grid%em_rv_m, grid%em_ww_m,        &
2587                             grid%em_mut, grid%em_mub, grid%em_mu_1,          &
2588                             grid%em_alt,                                     &
2589                             chem_old(ims,kms,jms,ic),                        &
2590                             chem(ims,kms,jms,ic),                            &
2591                             chem_tend(ims,kms,jms,ic),                       &
2592                             advect_tend,grid%rqvften,                        &
2593                             grid%qv_base, .false., grid%em_fnm, grid%em_fnp, &
2594                             grid%msfu, grid%msfv, grid%msft,                 &
2595                             grid%rdx, grid%rdy, grid%em_rdn, grid%em_rdnw,   &
2596                             grid%khdif, grid%kvdif, grid%xkmhd,              &
2597                             grid%diff_6th_opt, grid%diff_6th_factor,         &
2598                             config_flags%pd_chem,             &
2599                             ids, ide, jds, jde, kds, kde,     &
2600                             ims, ime, jms, jme, kms, kme,     &
2601                             grid%i_start(ij), grid%i_end(ij), &
2602                             grid%j_start(ij), grid%j_end(ij), &
2603                             k_start    , k_end               )
2604!
2605! Currently, chemistry species with specified boundaries (i.grid%e. the mother
2606! domain)  are being over written by flow_dep_bdy_chem. So, relax_bdy and
2607! spec_bdy are only called for nests. For boundary conditions from global model or larger domain,
2608! chem is uncoupled, and only used for one row/column on inflow (if have_bcs_chem=.true.)
2609!
2610     IF( ( config_flags%nested ) .and. rk_step == 1 ) THEN
2611       if(ic.eq.1)CALL wrf_debug ( 10 , ' have_bcs_chem' )
2612
2613         CALL relax_bdy_scalar ( chem_tend(ims,kms,jms,ic),             &
2614                                 chem(ims,kms,jms,ic),  grid%em_mut,            &
2615                                 chem_b(1,1,1,1,ic),                    &
2616                                 chem_bt(1,1,1,1,ic),                   &
2617                                 config_flags%spec_bdy_width, grid%spec_zone, grid%relax_zone, &
2618                                 grid%dtbc, grid%fcx, grid%gcx,             &
2619                                 config_flags,               &
2620                                 ijds, ijde,                 & ! min/max(id,jd)
2621                                 ids,ide, jds,jde, kds,kde,  & ! domain dims
2622                                 ims,ime, jms,jme, kms,kme,  & ! memory dims
2623                                 ips,ipe, jps,jpe, kps,kpe,  & ! patch  dims
2624                                 grid%i_start(ij), grid%i_end(ij),      &
2625                                 grid%j_start(ij), grid%j_end(ij),      &
2626                                 k_start, k_end                         )
2627         CALL spec_bdy_scalar  ( chem_tend(ims,kms,jms,ic),                 &
2628                                 chem_b(1,1,1,1,ic),                        &
2629                                 chem_bt(1,1,1,1,ic),                       &
2630                                 config_flags%spec_bdy_width, grid%spec_zone,                 &
2631                                 config_flags,               &
2632                                 ijds, ijde,                 & ! min/max(id,jd)
2633                                 ids,ide, jds,jde, kds,kde,  & ! domain dims
2634                                 ims,ime, jms,jme, kms,kme,  & ! memory dims
2635                                 ips,ipe, jps,jpe, kps,kpe,  & ! patch  dims
2636                                 grid%i_start(ij), grid%i_end(ij),          &
2637                                 grid%j_start(ij), grid%j_end(ij),          &
2638                                 k_start, k_end                             )
2639     ENDIF
2640
2641   ENDDO chem_tile_loop_1
2642
2643   !$OMP END PARALLEL DO
2644
2645
2646   !$OMP PARALLEL DO   &
2647   !$OMP PRIVATE ( ij )
2648
2649   chem_tile_loop_2: DO ij = 1 , grid%num_tiles
2650
2651       CALL wrf_debug ( 200 , ' call rk_update_scalar' )
2652       CALL rk_update_scalar( ic, ic,                           &
2653                              chem_old(ims,kms,jms,ic),         &  ! was chem_1
2654                              chem(ims,kms,jms,ic),             &
2655                              chem_tend(ims,kms,jms,ic),        &
2656                              advect_tend, grid%msft,                &
2657                              grid%em_mu_1, grid%em_mu_2, grid%em_mub,                  &
2658                              rk_step, dt_rk, grid%spec_zone,        &
2659                              config_flags,     &
2660                              ids, ide, jds, jde, kds, kde,     &
2661                              ims, ime, jms, jme, kms, kme,     &
2662                              grid%i_start(ij), grid%i_end(ij), &
2663                              grid%j_start(ij), grid%j_end(ij), &
2664                              k_start    , k_end               )
2665
2666
2667       IF( config_flags%specified  ) THEN
2668! come back to this and figure out why two different routines are needed. JM 20041203
2669!#ifndef WRF_CHEM
2670!!$           CALL flow_dep_bdy  ( chem(ims,kms,jms,ic),       &
2671!!$                                grid%em_ru_m, grid%em_rv_m, config_flags,   &
2672!!$                                grid%spec_zone,             &
2673!!$                                ids,ide, jds,jde, kds,kde,  & ! domain dims
2674!!$                                ims,ime, jms,jme, kms,kme,  & ! memory dims
2675!!$                                ips,ipe, jps,jpe, kps,kpe,  & ! patch  dims
2676!!$                                grid%i_start(ij), grid%i_end(ij),  &
2677!!$                                grid%j_start(ij), grid%j_end(ij),  &
2678!!$                                k_start, k_end                    )
2679!#else
2680           CALL flow_dep_bdy_chem( chem(ims,kms,jms,ic),            &
2681                                chem_b(1,1,1,1,ic),                 &
2682                                chem_bt(1,1,1,1,ic), dt_rk,         &
2683                                config_flags%spec_bdy_width,grid%em_z,      &
2684                                ijds, ijde,grid%have_bcs_chem,      &
2685                                grid%em_ru_m, grid%em_rv_m, config_flags,grid%em_alt,       &
2686                                grid%em_t_1,grid%em_pb,grid%em_p,t0,p1000mb,rcp,grid%em_ph_2,grid%em_phb,g, &
2687                                grid%spec_zone,ic,                  &
2688                                ids,ide, jds,jde, kds,kde,  & ! domain dims
2689                                ims,ime, jms,jme, kms,kme,  & ! memory dims
2690                                ips,ipe, jps,jpe, kps,kpe,  & ! patch  dims
2691                                grid%i_start(ij), grid%i_end(ij),   &
2692                                grid%j_start(ij), grid%j_end(ij),   &
2693                                k_start, k_end                      )
2694!#endif
2695       ENDIF
2696
2697
2698   ENDDO chem_tile_loop_2
2699   !$OMP END PARALLEL DO
2700
2701   ENDDO chem_variable_loop
2702
2703 ENDIF chem_scalar_advance
2704BENCH_END(chem_adv_tim)
2705#endif
2706
2707!  next the other scalar species
2708  other_scalar_advance: IF (num_3d_s >= PARAM_FIRST_SCALAR)  THEN
2709
2710   scalar_variable_loop: do is = PARAM_FIRST_SCALAR, num_3d_s
2711   !$OMP PARALLEL DO   &
2712   !$OMP PRIVATE ( ij )
2713
2714   scalar_tile_loop_1: DO ij = 1 , grid%num_tiles
2715
2716       CALL wrf_debug ( 200 , ' call rk_scalar_tend' )
2717       CALL rk_scalar_tend ( is, is, config_flags,                            &
2718                             rk_step, dt_rk,                                  &
2719                             grid%em_ru_m, grid%em_rv_m, grid%em_ww_m,        &
2720                             grid%em_mut, grid%em_mub, grid%em_mu_1,          &
2721                             grid%em_alt,                                     &
2722                             scalar_old(ims,kms,jms,is),                      &
2723                             scalar(ims,kms,jms,is),                          &
2724                             scalar_tend(ims,kms,jms,is),                     &
2725                             advect_tend,grid%rqvften,                        &
2726                             grid%qv_base, .false., grid%em_fnm, grid%em_fnp, &
2727                             grid%msfu, grid%msfv, grid%msft,                 &
2728                             grid%rdx, grid%rdy, grid%em_rdn, grid%em_rdnw,   &
2729                             grid%khdif, grid%kvdif, grid%xkmhd,              &
2730                             grid%diff_6th_opt, grid%diff_6th_factor,         &
2731                             config_flags%pd_scalar,           &
2732                             ids, ide, jds, jde, kds, kde,     &
2733                             ims, ime, jms, jme, kms, kme,     &
2734                             grid%i_start(ij), grid%i_end(ij), &
2735                             grid%j_start(ij), grid%j_end(ij), &
2736                             k_start    , k_end               )
2737
2738!!!!****MARS
2739!     IF( config_flags%nested .and. (rk_step == 1) ) THEN
2740!       IF (is .eq. P_QNI) THEN
2741     IF ( ( config_flags%specified .or. config_flags%nested ) .and. rk_step == 1) THEN
2742
2743!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2744!!!!****MARS MARS MARS: what to do with tracers at the boundaries ??? !!!!!
2745!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2746!IF (  ( (config_flags%mars .eq.  1) .AND. (is .EQ. 2) ) &    !! 1   pass only water vapor
2747IF  ( (  config_flags%mars .eq.  1                    ) &    !! 1   pass water vapor AND ice
2748 .OR. (  config_flags%mars .eq.  2                    ) &    !! 2   pass the dust tracer
2749 .OR. (  config_flags%mars .eq.  3                    ) &    !! 3   pass dust Q and N
2750 .OR. (  config_flags%mars .eq.  10                   ) &    !! 10  pass the co2 tracer
2751! .OR. ( (config_flags%mars .eq. 11) .AND. (is .EQ. 2) ) &    !! 11  pass only water vapor
2752 .OR. (  config_flags%mars .eq. 11                    ) &    !! 11  pass EVERYTHING
2753 .OR. (  config_flags%mars .eq. 12                    ) &    !! 12  pass EVERYTHING
2754 .OR. (  config_flags%mars .eq. 34                    ) &
2755 .OR. config_flags%nested ) THEN                             !! *   pass all tracers if nested
2756         CALL relax_bdy_scalar ( scalar_tend(ims,kms,jms,is),            &
2757                                 scalar(ims,kms,jms,is),  grid%em_mut,         &
2758                                 scalar_b(1,1,1,1,is),                   &
2759                                 scalar_bt(1,1,1,1,is),                  &
2760                                 config_flags%spec_bdy_width, grid%spec_zone, grid%relax_zone, &
2761                                 grid%dtbc, grid%fcx, grid%gcx,             &
2762                                 config_flags,               &
2763                                 ijds, ijde,                 & ! min/max(id,jd)
2764                                 ids,ide, jds,jde, kds,kde,  & ! domain dims
2765                                 ims,ime, jms,jme, kms,kme,  & ! memory dims
2766                                 ips,ipe, jps,jpe, kps,kpe,  & ! patch  dims
2767                                 grid%i_start(ij), grid%i_end(ij),      &
2768                                 grid%j_start(ij), grid%j_end(ij),      &
2769                                 k_start, k_end                        )
2770
2771         CALL spec_bdy_scalar  ( scalar_tend(ims,kms,jms,is),                &
2772                                 scalar_b(1,1,1,1,is),                   &
2773                                 scalar_bt(1,1,1,1,is),                  &
2774                                 config_flags%spec_bdy_width, grid%spec_zone, &
2775                                 config_flags,               &
2776                                 ijds, ijde,                 & ! min/max(id,jd)
2777                                 ids,ide, jds,jde, kds,kde,  & ! domain dims
2778                                 ims,ime, jms,jme, kms,kme,  & ! memory dims
2779                                 ips,ipe, jps,jpe, kps,kpe,  & ! patch  dims
2780                                 grid%i_start(ij), grid%i_end(ij),          &
2781                                 grid%j_start(ij), grid%j_end(ij),          &
2782                                 k_start, k_end                               )
2783ENDIF
2784
2785      ENDIF
2786
2787   ENDDO scalar_tile_loop_1
2788   !$OMP END PARALLEL DO
2789
2790
2791   !$OMP PARALLEL DO   &
2792   !$OMP PRIVATE ( ij )
2793
2794   scalar_tile_loop_2: DO ij = 1 , grid%num_tiles
2795
2796       CALL wrf_debug ( 200 , ' call rk_update_scalar' )
2797       CALL rk_update_scalar( is, is,                           &
2798                              scalar_old(ims,kms,jms,is),       &  ! was scalar_1
2799                              scalar(ims,kms,jms,is),           &
2800                              scalar_tend(ims,kms,jms,is),      &
2801                              advect_tend, grid%msft,                &
2802                              grid%em_mu_1, grid%em_mu_2, grid%em_mub,                  &
2803                              rk_step, dt_rk, grid%spec_zone,        &
2804                              config_flags,     &
2805                              ids, ide, jds, jde, kds, kde,     &
2806                              ims, ime, jms, jme, kms, kme,     &
2807                              grid%i_start(ij), grid%i_end(ij), &
2808                              grid%j_start(ij), grid%j_end(ij), &
2809                              k_start    , k_end               )
2810
2811!!!!!****MARS: if parent domain and any other tracer than water vapor, no bdy and simple zero flux bdy condition
2812IF ( config_flags%specified ) THEN
2813!IF (      ((config_flags%mars .eq. 1 ) .and. (is .ne. 2))  &   
2814!     .OR. ((config_flags%mars .eq. 11) .and. (is .ne. 2))  ) THEN
2815IF ( config_flags%mars .gt. 50 ) THEN
2816
2817     !!! YOU HAVE TO ADD A CONDITION HERE IF YOU ADD A mars OPTION in REGISTRY  !!!
2818     !!! (this is OK if you are using LES with periodic boundary conditions)   
2819
2820        CALL flow_dep_bdy  ( scalar(ims,kms,jms,is),     &
2821                                grid%em_ru_m, grid%em_rv_m, config_flags,   &
2822                                grid%spec_zone,                  &
2823                                ids,ide, jds,jde, kds,kde,  & ! domain dims
2824                                ims,ime, jms,jme, kms,kme,  & ! memory dims
2825                                ips,ipe, jps,jpe, kps,kpe,  & ! patch  dims
2826                                grid%i_start(ij), grid%i_end(ij),  &
2827                                grid%j_start(ij), grid%j_end(ij),  &
2828                                k_start, k_end                    )
2829
2830ENDIF
2831ENDIF
2832
2833   ENDDO scalar_tile_loop_2
2834   !$OMP END PARALLEL DO
2835
2836   ENDDO scalar_variable_loop
2837
2838 ENDIF other_scalar_advance
2839
2840 !  update the pressure and density at the new time level
2841
2842   !$OMP PARALLEL DO   &
2843   !$OMP PRIVATE ( ij )
2844   DO ij = 1 , grid%num_tiles
2845
2846BENCH_START(calc_p_rho_tim)
2847
2848     CALL calc_p_rho_phi( moist, num_3d_m,                &
2849                          grid%em_al, grid%em_alb, grid%em_mu_2, grid%em_muts,              &
2850                          grid%em_ph_2, grid%em_p, grid%em_pb, grid%em_t_2,                 &
2851                          p0, t0, grid%em_znu, grid%em_dnw, grid%em_rdnw,           &
2852                          grid%em_rdn, config_flags%non_hydrostatic,             &
2853                          ids, ide, jds, jde, kds, kde,     &
2854                          ims, ime, jms, jme, kms, kme,     &
2855                          grid%i_start(ij), grid%i_end(ij), &
2856                          grid%j_start(ij), grid%j_end(ij), &
2857                          k_start    , k_end               )
2858
2859BENCH_END(calc_p_rho_tim)
2860
2861   ENDDO
2862   !$OMP END PARALLEL DO
2863
2864!  Reset the boundary conditions if there is another corrector step.
2865!  (rk_step < rk_order), else we'll handle it at the end of everything
2866!  (after the split physics, before exiting the timestep).
2867
2868   rk_step_1_check: IF ( rk_step < rk_order ) THEN
2869
2870!-----------------------------------------------------------
2871!  Stencils for patch communications  (WCS, 29 June 2001)
2872!
2873!  here's where we need a wide comm stencil - these are the
2874!  uncoupled variables so are used for high order calc in
2875!  advection and mixong routines.
2876!
2877!                              * * * * *
2878!            *        * * *    * * * * *
2879!          * + *      * + *    * * + * *
2880!            *        * * *    * * * * *
2881!                              * * * * *
2882!
2883!
2884! grid%em_u_2                              x
2885! grid%em_v_2                              x
2886! grid%em_w_2                              x
2887! grid%em_t_2                              x
2888! grid%em_ph_2                             x
2889! grid%em_al         x
2890!
2891!  2D variable
2892! grid%em_mu_2       x
2893!
2894!  4D variable
2895! moist               x
2896! chem                x
2897!scalar               x
2898
2899#ifdef DM_PARALLEL
2900   IF      ( config_flags%h_mom_adv_order <= 4 ) THEN
2901#    include "HALO_EM_D2_3.inc"
2902   ELSE IF ( config_flags%h_mom_adv_order <= 6 ) THEN
2903#    include "HALO_EM_D2_5.inc"
2904   ELSE
2905     WRITE(wrf_err_message,*)'solve_em: invalid h_mom_adv_order = ',config_flags%h_mom_adv_order
2906     CALL wrf_error_fatal(TRIM(wrf_err_message))
2907   ENDIF
2908#  include "PERIOD_BDY_EM_D.inc"
2909#  include "PERIOD_BDY_EM_MOIST2.inc"
2910#  include "PERIOD_BDY_EM_CHEM2.inc"
2911#  include "PERIOD_BDY_EM_SCALAR2.inc"
2912#endif
2913
2914BENCH_START(bc_end_tim)
2915   !$OMP PARALLEL DO   &
2916   !$OMP PRIVATE ( ij )
2917
2918    tile_bc_loop_1: DO ij = 1 , grid%num_tiles
2919
2920
2921      CALL wrf_debug ( 200 , ' call rk_phys_bc_dry_2' )
2922
2923      CALL rk_phys_bc_dry_2( config_flags,                         &
2924                             grid%em_u_2, grid%em_v_2, grid%em_w_2,                    &
2925                             grid%em_t_2, grid%em_ph_2, grid%em_mu_2,                  &
2926                             ids, ide, jds, jde, kds, kde,     &
2927                             ims, ime, jms, jme, kms, kme,     &
2928                             ips, ipe, jps, jpe, kps, kpe,     &
2929                             grid%i_start(ij), grid%i_end(ij), &
2930                             grid%j_start(ij), grid%j_end(ij), &
2931                             k_start    , k_end               )
2932
2933BENCH_START(diag_w_tim)
2934     IF (.not. config_flags%non_hydrostatic) THEN
2935     CALL diagnose_w( ph_tend, grid%em_ph_2, grid%em_ph_1, grid%em_w_2, grid%em_muts, dt_rk,  &
2936                      grid%em_u_2, grid%em_v_2, grid%ht,                           &
2937                      grid%cf1, grid%cf2, grid%cf3, grid%rdx, grid%rdy, grid%msft,          &
2938                      ids, ide, jds, jde, kds, kde,           &
2939                      ims, ime, jms, jme, kms, kme,           &
2940                      grid%i_start(ij), grid%i_end(ij),       &
2941                      grid%j_start(ij), grid%j_end(ij),       &
2942                      k_start    , k_end                     )
2943     ENDIF
2944BENCH_END(diag_w_tim)
2945
2946      IF (num_3d_m >= PARAM_FIRST_SCALAR) THEN
2947
2948        moisture_loop_bdy_1 : DO im = PARAM_FIRST_SCALAR , num_3d_m
2949 
2950          CALL set_physical_bc3d( moist(ims,kms,jms,im), 'p', config_flags,   &
2951                                   ids, ide, jds, jde, kds, kde,             &
2952                                   ims, ime, jms, jme, kms, kme,             &
2953                                   ips, ipe, jps, jpe, kps, kpe,             &
2954                                   grid%i_start(ij), grid%i_end(ij),                   &
2955                                   grid%j_start(ij), grid%j_end(ij),                   &
2956                                   k_start    , k_end                       )
2957         END DO moisture_loop_bdy_1
2958
2959      ENDIF
2960
2961      IF (num_3d_c >= PARAM_FIRST_SCALAR) THEN
2962
2963        chem_species_bdy_loop_1 : DO ic = PARAM_FIRST_SCALAR , num_3d_c
2964
2965          CALL set_physical_bc3d( chem(ims,kms,jms,ic), 'p', config_flags,   &
2966                                  ids, ide, jds, jde, kds, kde,            &
2967                                  ims, ime, jms, jme, kms, kme,            &
2968                                  ips, ipe, jps, jpe, kps, kpe,            &
2969                                  grid%i_start(ij), grid%i_end(ij),                  &
2970                                  grid%j_start(ij), grid%j_end(ij),                  &
2971                                  k_start    , k_end-1                    )
2972
2973        END DO chem_species_bdy_loop_1
2974
2975      END IF
2976
2977      IF (num_3d_s >= PARAM_FIRST_SCALAR) THEN
2978
2979        scalar_species_bdy_loop_1 : DO is = PARAM_FIRST_SCALAR , num_3d_s
2980
2981          CALL set_physical_bc3d( scalar(ims,kms,jms,is), 'p', config_flags,   &
2982                                  ids, ide, jds, jde, kds, kde,            &
2983                                  ims, ime, jms, jme, kms, kme,            &
2984                                  ips, ipe, jps, jpe, kps, kpe,            &
2985                                  grid%i_start(ij), grid%i_end(ij),                  &
2986                                  grid%j_start(ij), grid%j_end(ij),                  &
2987                                  k_start    , k_end-1                    )
2988
2989        END DO scalar_species_bdy_loop_1
2990
2991      END IF
2992
2993      IF (config_flags%km_opt .eq. 2) THEN
2994
2995        CALL set_physical_bc3d( grid%em_tke_2 , 'p', config_flags,  &
2996                                ids, ide, jds, jde, kds, kde,            &
2997                                ims, ime, jms, jme, kms, kme,            &
2998                                ips, ipe, jps, jpe, kps, kpe,            &
2999                                grid%i_start(ij), grid%i_end(ij),        &
3000                                grid%j_start(ij), grid%j_end(ij),        &
3001                                k_start    , k_end                      )
3002      END IF
3003
3004    END DO tile_bc_loop_1
3005   !$OMP END PARALLEL DO
3006BENCH_END(bc_end_tim)
3007
3008
3009#ifdef DM_PARALLEL
3010
3011!                           * * * * *
3012!         *        * * *    * * * * *
3013!       * + *      * + *    * * + * *
3014!         *        * * *    * * * * *
3015!                           * * * * *
3016
3017! moist, chem, scalar, tke      x
3018
3019
3020      IF      ( config_flags%h_mom_adv_order <= 4 ) THEN
3021        IF ( (config_flags%pd_tke) .and. (rk_step == rk_order-1) ) THEN
3022#         include "HALO_EM_TKE_5.inc"
3023        ELSE
3024#         include "HALO_EM_TKE_3.inc"
3025        ENDIF
3026      ELSE IF ( config_flags%h_mom_adv_order <= 6 ) THEN
3027        IF ( (config_flags%pd_tke) .and. (rk_step == rk_order-1) ) THEN
3028#         include "HALO_EM_TKE_7.inc"
3029        ELSE
3030#         include "HALO_EM_TKE_5.inc"
3031        ENDIF
3032      ELSE
3033        WRITE(wrf_err_message,*)'solve_em: invalid h_sca_adv_order = ',config_flags%h_sca_adv_order
3034        CALL wrf_error_fatal(TRIM(wrf_err_message))
3035      ENDIF
3036
3037#if 0
3038   IF (config_flags%km_opt .eq. 2) THEN
3039#      include  "HALO_EM_TKE_F.inc"
3040   ENDIF
3041#endif
3042
3043   if ( num_moist .ge. PARAM_FIRST_SCALAR ) then
3044     IF      ( config_flags%h_sca_adv_order <= 4 ) THEN
3045       IF ( (config_flags%pd_moist) .and. (rk_step == rk_order-1) ) THEN
3046#        include "HALO_EM_MOIST_E_5.inc"
3047       ELSE
3048#        include "HALO_EM_MOIST_E_3.inc"
3049       END IF
3050     ELSE IF ( config_flags%h_sca_adv_order <= 6 ) THEN
3051       IF ( (config_flags%pd_moist) .and. (rk_step == rk_order-1) ) THEN
3052#        include "HALO_EM_MOIST_E_7.inc"
3053       ELSE
3054#        include "HALO_EM_MOIST_E_5.inc"
3055       END IF
3056     ELSE
3057       WRITE(wrf_err_message,*)'solve_em: invalid h_sca_adv_order = ',config_flags%h_sca_adv_order
3058       CALL wrf_error_fatal(TRIM(wrf_err_message))
3059     ENDIF
3060   endif
3061   if ( num_chem >= PARAM_FIRST_SCALAR ) then
3062     IF      ( config_flags%h_sca_adv_order <= 4 ) THEN
3063       IF ( (config_flags%pd_chem) .and. (rk_step == rk_order-1) ) THEN
3064#        include "HALO_EM_CHEM_E_5.inc"
3065       ELSE
3066#        include "HALO_EM_CHEM_E_3.inc"
3067       END IF
3068     ELSE IF ( config_flags%h_sca_adv_order <= 6 ) THEN
3069       IF ( (config_flags%pd_chem) .and. (rk_step == rk_order-1) ) THEN
3070#        include "HALO_EM_CHEM_E_7.inc"
3071       ELSE
3072#        include "HALO_EM_CHEM_E_5.inc"
3073       END IF
3074     ELSE
3075       WRITE(wrf_err_message,*)'solve_em: invalid h_sca_adv_order = ',config_flags%h_sca_adv_order
3076       CALL wrf_error_fatal(TRIM(wrf_err_message))
3077     ENDIF
3078   endif
3079   if ( num_scalar >= PARAM_FIRST_SCALAR ) then
3080     IF      ( config_flags%h_sca_adv_order <= 4 ) THEN
3081       IF ( (config_flags%pd_scalar) .and. (rk_step == rk_order-1) ) THEN
3082#        include "HALO_EM_SCALAR_E_5.inc"
3083       ELSE
3084#        include "HALO_EM_SCALAR_E_3.inc"
3085       END IF
3086     ELSE IF ( config_flags%h_sca_adv_order <= 6 ) THEN
3087       IF ( (config_flags%pd_scalar) .and. (rk_step == rk_order-1) ) THEN
3088#        include "HALO_EM_SCALAR_E_7.inc"
3089       ELSE
3090#        include "HALO_EM_SCALAR_E_5.inc"
3091       END IF
3092     ELSE
3093       WRITE(wrf_err_message,*)'solve_em: invalid h_sca_adv_order = ',config_flags%h_sca_adv_order
3094       CALL wrf_error_fatal(TRIM(wrf_err_message))
3095     ENDIF
3096   endif
3097#endif
3098
3099   ENDIF rk_step_1_check
3100
3101
3102!**********************************************************
3103!
3104!  end of RK predictor-corrector loop
3105!
3106!**********************************************************
3107
3108 END DO Runge_Kutta_loop
3109
3110   !$OMP PARALLEL DO   &
3111   !$OMP PRIVATE ( ij )
3112
3113   DO ij = 1 , grid%num_tiles
3114
3115BENCH_START(advance_ppt_tim)
3116      CALL wrf_debug ( 200 , ' call advance_ppt' )
3117      CALL advance_ppt(grid%rthcuten,grid%rqvcuten,grid%rqccuten,grid%rqrcuten, &
3118                     grid%rqicuten,grid%rqscuten,grid%rainc,grid%raincv,grid%nca,    &
3119                     grid%htop,grid%hbot,grid%cutop,grid%cubot,                 &
3120                     grid%cuppt, config_flags,                   &
3121                     ids,ide, jds,jde, kds,kde,             &
3122                     ims,ime, jms,jme, kms,kme,             &
3123                     grid%i_start(ij), grid%i_end(ij),      &
3124                     grid%j_start(ij), grid%j_end(ij),      &
3125                     k_start    , k_end                    )
3126BENCH_END(advance_ppt_tim)
3127
3128   ENDDO
3129   !$OMP END PARALLEL DO
3130
3131!<DESCRIPTION>
3132!<pre>
3133! (5) time-split physics.
3134!
3135!     Microphysics are the only time  split physics in the WRF model
3136!     at this time.  Split-physics begins with the calculation of
3137!     needed diagnostic quantities (pressure, temperature, etc.)
3138!     followed by a call to the microphysics driver,
3139!     and finishes with a clean-up, storing off of a diabatic tendency
3140!     from the moist physics, and a re-calulation of the  diagnostic
3141!     quantities pressure and density.
3142!</pre>
3143!</DESCRIPTION>
3144
3145!!!******MARS MARS
3146!!!******MARS MARS
3147!!!******MARS MARS
3148
3149!  IF (config_flags%mp_physics /= 0)  then
3150!
3151!   IF( config_flags%specified .or. config_flags%nested ) THEN
3152!     sz = grid%spec_zone
3153!   ELSE
3154!     sz = 0
3155!   ENDIF
3156!
3157!   !$OMP PARALLEL DO   &
3158!   !$OMP PRIVATE ( ij, its, ite, jts, jte )
3159!
3160!   scalar_tile_loop_1a: DO ij = 1 , grid%num_tiles
3161!
3162!       IF ( config_flags%periodic_x ) THEN
3163!         its = max(grid%i_start(ij),ids)
3164!         ite = min(grid%i_end(ij),ide-1)
3165!       ELSE
3166!         its = max(grid%i_start(ij),ids+sz)
3167!         ite = min(grid%i_end(ij),ide-1-sz)
3168!       ENDIF
3169!       jts = max(grid%j_start(ij),jds+sz)
3170!       jte = min(grid%j_end(ij),jde-1-sz)
3171!
3172!       CALL wrf_debug ( 200 , ' call moist_physics_prep' )
3173!BENCH_START(moist_physics_prep_tim)
3174!       CALL moist_physics_prep_em( grid%em_t_2, grid%em_t_1, t0, rho,                &
3175!                                   grid%em_al, grid%em_alb, grid%em_p, p8w, p0, grid%em_pb,          &
3176!                                   grid%em_ph_2, grid%em_phb, th_phy, pi_phy, p_phy, &
3177!                                   grid%em_z, z_at_w, dz8w,                  &
3178!                                   dtm, grid%h_diabatic,                  &
3179!                                   config_flags,grid%em_fnm, grid%em_fnp,            &
3180!                                   ids, ide, jds, jde, kds, kde,     &
3181!                                   ims, ime, jms, jme, kms, kme,     &
3182!                                   its, ite, jts, jte,               &
3183!                                   k_start    , k_end               )
3184!BENCH_END(moist_physics_prep_tim)
3185!   END DO scalar_tile_loop_1a
3186!   !$OMP END PARALLEL DO
3187!
3188!       CALL wrf_debug ( 200 , ' call microphysics_driver' )
3189!
3190!       grid%em_sr = 0.
3191!       specified_bdy = config_flags%specified .OR. config_flags%nested
3192!       channel_bdy = config_flags%specified .AND. config_flags%periodic_x
3193!
3194!#if 0
3195!BENCH_START(microswap_1)
3196!! for load balancing; communication to redistribute the points
3197!   IF ( config_flags%mp_physics .EQ. ETAMPNEW ) THEN
3198!#include "SWAP_ETAMP_NEW.inc"
3199!   ELSE IF ( config_flags%mp_physics .EQ. WSM3SCHEME ) THEN
3200!#include "SWAP_WSM3.inc"
3201!   ENDIF
3202!BENCH_END(microswap_1)
3203!#endif
3204!
3205!BENCH_START(micro_driver_tim)
3206!
3207!       CALL microphysics_driver(                                          &
3208!     &         DT=dtm             ,DX=grid%dx              ,DY=grid%dy              &
3209!     &        ,DZ8W=dz8w          ,F_ICE_PHY=grid%f_ice_phy                    &
3210!     &        ,ITIMESTEP=grid%itimestep                    ,LOWLYR=grid%lowlyr      &
3211!     &        ,P8W=p8w            ,P=p_phy            ,PI_PHY=pi_phy      &
3212!     &        ,RHO=rho            ,SPEC_ZONE=grid%spec_zone                    &
3213!     &        ,SR=grid%em_sr              ,TH=th_phy                              &
3214!     &        ,WARM_RAIN=grid%warm_rain                    ,XLAND=grid%xland        &
3215!     &        ,SPECIFIED=specified_bdy, CHANNEL_SWITCH=channel_bdy        &
3216!     &        ,F_RAIN_PHY=grid%f_rain_phy                                      &
3217!     &        ,F_RIMEF_PHY=grid%f_rimef_phy                                    &
3218!     &        ,MP_PHYSICS=config_flags%mp_physics                         &
3219!     &        ,IDS=ids,IDE=ide, JDS=jds,JDE=jde, KDS=kds,KDE=kde          &
3220!     &        ,IMS=ims,IME=ime, JMS=jms,JME=jme, KMS=kms,KME=kme          &
3221!     &        ,I_START=grid%i_start,I_END=min(grid%i_end, ide-1)          &
3222!     &        ,J_START=grid%j_start,J_END=min(grid%j_end, jde-1)          &
3223!     &        ,KTS=k_start, KTE=min(k_end,kde-1)                          &
3224!     &        ,NUM_TILES=grid%num_tiles                                   &
3225!                 ! Optional
3226!     &        , RAINNC=grid%rainnc, RAINNCV=grid%rainncv                            &
3227!     &        , SNOWNC=grid%snownc, SNOWNCV=grid%snowncv                            &
3228!     &        , GRAUPELNC=grid%graupelnc, GRAUPELNCV=grid%graupelncv                &
3229!     &        , W=grid%em_w_2, Z=grid%em_z, HT=grid%ht                                         &
3230!     &        , MP_RESTART_STATE=grid%mp_restart_state                         &
3231!     &        , TBPVS_STATE=grid%tbpvs_state                                   & ! etampnew
3232!     &        , TBPVS0_STATE=grid%tbpvs0_state                                 & ! etampnew
3233!     &        , QV_CURR=moist(ims,kms,jms,P_QV), F_QV=F_QV              &
3234!     &        , QC_CURR=moist(ims,kms,jms,P_QC), F_QC=F_QC              &
3235!     &        , QR_CURR=moist(ims,kms,jms,P_QR), F_QR=F_QR              &
3236!     &        , QI_CURR=moist(ims,kms,jms,P_QI), F_QI=F_QI              &
3237!     &        , QS_CURR=moist(ims,kms,jms,P_QS), F_QS=F_QS              &
3238!     &        , QG_CURR=moist(ims,kms,jms,P_QG), F_QG=F_QG              &
3239!     &        , QNI_CURR=scalar(ims,kms,jms,P_QNI), F_QNI=F_QNI         &
3240!     &        , QT_CURR=scalar(ims,kms,jms,P_QT), F_QT=F_QT             &
3241!                                                                          )
3242!BENCH_END(micro_driver_tim)
3243!
3244!#if 0
3245!BENCH_START(microswap_2)
3246!! for load balancing; communication to redistribute the points
3247!   IF ( config_flags%mp_physics .EQ. ETAMPNEW ) THEN
3248!#include "SWAP_ETAMP_NEW.inc"
3249!   ELSE IF ( config_flags%mp_physics .EQ. WSM3SCHEME ) THEN
3250!#include "SWAP_WSM3.inc"
3251!   ENDIF
3252!BENCH_END(microswap_2)
3253!#endif
3254!
3255!       CALL wrf_debug ( 200 , ' call moist_physics_finish' )
3256!BENCH_START(moist_phys_end_tim)
3257!   !$OMP PARALLEL DO   &
3258!   !$OMP PRIVATE ( ij, its, ite, jts, jte )
3259!
3260!   scalar_tile_loop_1b: DO ij = 1 , grid%num_tiles
3261!
3262!       IF ( config_flags%periodic_x ) THEN
3263!         its = max(grid%i_start(ij),ids)
3264!         ite = min(grid%i_end(ij),ide-1)
3265!       ELSE
3266!         its = max(grid%i_start(ij),ids+sz)
3267!         ite = min(grid%i_end(ij),ide-1-sz)
3268!       ENDIF
3269!       jts = max(grid%j_start(ij),jds+sz)
3270!       jte = min(grid%j_end(ij),jde-1-sz)
3271!
3272!       CALL microphysics_zero_out (                                    &
3273!                     moist , num_moist , config_flags ,              &
3274!                     ids, ide, jds, jde, kds, kde,                     &
3275!                     ims, ime, jms, jme, kms, kme,                     &
3276!                     its, ite, jts, jte,                               &
3277!                     k_start    , k_end                                )
3278!
3279!       CALL moist_physics_finish_em( grid%em_t_2, grid%em_t_1, t0, grid%em_muts, th_phy,       &
3280!                                     grid%h_diabatic, dtm, config_flags,    &
3281!                                     ids, ide, jds, jde, kds, kde,     &
3282!                                     ims, ime, jms, jme, kms, kme,     &
3283!                                     its, ite, jts, jte,               &
3284!                                     k_start    , k_end               )
3285!
3286!       CALL calc_p_rho_phi( moist, num_3d_m,                &
3287!                            grid%em_al, grid%em_alb, grid%em_mu_2, grid%em_muts,              &
3288!                            grid%em_ph_2, grid%em_p, grid%em_pb, grid%em_t_2,                 &
3289!                            p0, t0, grid%em_znu, grid%em_dnw, grid%em_rdnw,           &
3290!                            grid%em_rdn, config_flags%non_hydrostatic,             &
3291!                            ids, ide, jds, jde, kds, kde,     &
3292!                            ims, ime, jms, jme, kms, kme,     &
3293!                            its, ite, jts, jte,               &
3294!                            k_start    , k_end               )
3295!
3296!   END DO scalar_tile_loop_1b
3297!   !$OMP END PARALLEL DO
3298!BENCH_END(moist_phys_end_tim)
3299!
3300!  ENDIF
3301
3302
3303  IF (.not. config_flags%non_hydrostatic) THEN
3304   !$OMP PARALLEL DO   &
3305   !$OMP PRIVATE ( ij )
3306   DO ij = 1 , grid%num_tiles
3307     CALL diagnose_w( ph_tend, grid%em_ph_2, grid%em_ph_1, grid%em_w_2, grid%em_muts, dt_rk,  &
3308                      grid%em_u_2, grid%em_v_2, grid%ht,                           &
3309                      grid%cf1, grid%cf2, grid%cf3, grid%rdx, grid%rdy, grid%msft,          &
3310                      ids, ide, jds, jde, kds, kde,           &
3311                      ims, ime, jms, jme, kms, kme,           &
3312                      grid%i_start(ij), grid%i_end(ij),       &
3313                      grid%j_start(ij), grid%j_end(ij),       &
3314                      k_start    , k_end                     )
3315
3316   END DO
3317   !$OMP END PARALLEL DO
3318  ENDIF
3319
3320   chem_tile_loop_3: DO ij = 1 , grid%num_tiles
3321
3322          CALL wrf_debug ( 200 , ' call scalar_tile_loop_2' )
3323
3324     IF ( num_3d_c >= PARAM_FIRST_SCALAR ) then
3325
3326!
3327!  tiled chemistry not here, it is called from solve_interface, and found in chem_driver
3328!
3329
3330     END IF
3331
3332   END DO chem_tile_loop_3
3333
3334
3335   !  We're finished except for boundary condition (and patch) update
3336
3337   ! Boundary condition time (or communication time).  At this time, we have
3338   ! implemented periodic and symmetric physical boundary conditions.
3339
3340   ! b.c. routine for data within patch.
3341
3342   ! we need to do both time levels of
3343   ! data because the time filter only works in the physical solution space.
3344
3345   ! First, do patch communications for boundary conditions (periodicity)
3346
3347!-----------------------------------------------------------
3348!  Stencils for patch communications  (WCS, 29 June 2001)
3349!
3350!  here's where we need a wide comm stencil - these are the
3351!  uncoupled variables so are used for high order calc in
3352!  advection and mixong routines.
3353!
3354!                              * * * * *
3355!            *        * * *    * * * * *
3356!          * + *      * + *    * * + * *
3357!            *        * * *    * * * * *
3358!                              * * * * *
3359!
3360!   grid%em_u_1                            x
3361!   grid%em_u_2                            x
3362!   grid%em_v_1                            x
3363!   grid%em_v_2                            x
3364!   grid%em_w_1                            x
3365!   grid%em_w_2                            x
3366!   grid%em_t_1                            x
3367!   grid%em_t_2                            x
3368!  grid%em_ph_1                            x
3369!  grid%em_ph_2                            x
3370!  grid%em_tke_1                           x
3371!  grid%em_tke_2                           x
3372!
3373!    2D variables
3374!  grid%em_mu_1     x
3375!  grid%em_mu_2     x
3376!
3377!    4D variables
3378!  moist                         x
3379!   chem                         x
3380! scalar                         x
3381!----------------------------------------------------------
3382
3383
3384#ifdef DM_PARALLEL
3385   IF      ( config_flags%h_mom_adv_order <= 4 ) THEN
3386#    include "HALO_EM_D3_3.inc"
3387   ELSE IF ( config_flags%h_mom_adv_order <= 6 ) THEN
3388#    include "HALO_EM_D3_5.inc"
3389   ELSE
3390     WRITE(wrf_err_message,*)'solve_em: invalid h_mom_adv_order = ',config_flags%h_mom_adv_order
3391     CALL wrf_error_fatal(TRIM(wrf_err_message))
3392   ENDIF
3393#  include "PERIOD_BDY_EM_D3.inc"
3394#  include "PERIOD_BDY_EM_MOIST.inc"
3395#  include "PERIOD_BDY_EM_CHEM.inc"
3396#  include "PERIOD_BDY_EM_SCALAR.inc"
3397#endif
3398
3399!  now set physical b.c on a patch
3400
3401BENCH_START(bc_2d_tim)
3402   !$OMP PARALLEL DO   &
3403   !$OMP PRIVATE ( ij )
3404
3405   tile_bc_loop_2: DO ij = 1 , grid%num_tiles
3406
3407
3408     CALL wrf_debug ( 200 , ' call set_phys_bc_dry_2' )
3409
3410     CALL set_phys_bc_dry_2( config_flags,                           &
3411                             grid%em_u_1, grid%em_u_2, grid%em_v_1, grid%em_v_2, grid%em_w_1, grid%em_w_2,           &
3412                             grid%em_t_1, grid%em_t_2, grid%em_ph_1, grid%em_ph_2, grid%em_mu_1, grid%em_mu_2,       &
3413                             ids, ide, jds, jde, kds, kde,           &
3414                             ims, ime, jms, jme, kms, kme,           &
3415                             ips, ipe, jps, jpe, kps, kpe,           &
3416                             grid%i_start(ij), grid%i_end(ij),       &
3417                             grid%j_start(ij), grid%j_end(ij),       &
3418                             k_start    , k_end                     )
3419
3420     CALL set_physical_bc3d( grid%em_tke_1, 'p', config_flags,   &
3421                             ids, ide, jds, jde, kds, kde,            &
3422                             ims, ime, jms, jme, kms, kme,            &
3423                             ips, ipe, jps, jpe, kps, kpe,            &
3424                             grid%i_start(ij), grid%i_end(ij),        &
3425                             grid%j_start(ij), grid%j_end(ij),        &
3426                             k_start    , k_end-1                    )
3427     CALL set_physical_bc3d( grid%em_tke_2 , 'p', config_flags,  &
3428                             ids, ide, jds, jde, kds, kde,            &
3429                             ims, ime, jms, jme, kms, kme,            &
3430                             ips, ipe, jps, jpe, kps, kpe,            &
3431                             grid%i_start(ij), grid%i_end(ij),        &
3432                             grid%j_start(ij), grid%j_end(ij),        &
3433                             k_start    , k_end                      )
3434
3435     moisture_loop_bdy_2 : DO im = PARAM_FIRST_SCALAR , num_3d_m
3436
3437       CALL set_physical_bc3d( moist(ims,kms,jms,im), 'p',           &
3438                               config_flags,                           &
3439                               ids, ide, jds, jde, kds, kde,           &
3440                               ims, ime, jms, jme, kms, kme,           &
3441                               ips, ipe, jps, jpe, kps, kpe,           &
3442                               grid%i_start(ij), grid%i_end(ij),       &
3443                               grid%j_start(ij), grid%j_end(ij),       &
3444                               k_start    , k_end                     )
3445
3446     END DO moisture_loop_bdy_2
3447
3448     chem_species_bdy_loop_2 : DO ic = PARAM_FIRST_SCALAR , num_3d_c
3449
3450       CALL set_physical_bc3d( chem(ims,kms,jms,ic) , 'p', config_flags,  &
3451                               ids, ide, jds, jde, kds, kde,            &
3452                               ims, ime, jms, jme, kms, kme,            &
3453                               ips, ipe, jps, jpe, kps, kpe,            &
3454                               grid%i_start(ij), grid%i_end(ij),                  &
3455                               grid%j_start(ij), grid%j_end(ij),                  &
3456                               k_start    , k_end                      )
3457
3458     END DO chem_species_bdy_loop_2
3459
3460     scalar_species_bdy_loop_2 : DO is = PARAM_FIRST_SCALAR , num_3d_s
3461
3462       CALL set_physical_bc3d( scalar(ims,kms,jms,is) , 'p', config_flags,  &
3463                               ids, ide, jds, jde, kds, kde,            &
3464                               ims, ime, jms, jme, kms, kme,            &
3465                               ips, ipe, jps, jpe, kps, kpe,            &
3466                               grid%i_start(ij), grid%i_end(ij),                  &
3467                               grid%j_start(ij), grid%j_end(ij),                  &
3468                               k_start    , k_end                      )
3469
3470     END DO scalar_species_bdy_loop_2
3471
3472   END DO tile_bc_loop_2
3473   !$OMP END PARALLEL DO
3474BENCH_END(bc_2d_tim)
3475
3476   IF( config_flags%specified .or. config_flags%nested ) THEN
3477     grid%dtbc = grid%dtbc + grid%dt
3478   ENDIF
3479
3480! calculate some model diagnostics.
3481
3482!!!******MARS MARS
3483!!!******MARS MARS
3484
3485!         CALL wrf_debug ( 200 , ' call diagnostic_driver' )
3486!   
3487!         CALL diagnostic_output_calc(                                     &
3488!     &              DPSDT=grid%dpsdt   ,DMUDT=grid%dmudt                  &
3489!     &             ,P_PHY=p_phy   ,PK1M=grid%pk1m                         &
3490!     &             ,MU_2=grid%em_mu_2  ,MU_2M=grid%mu_2m                  &
3491!     &             ,U=grid%em_u_2    ,V=grid%em_v_2                       &
3492!     &             ,RAINCV=grid%raincv    ,RAINNCV=grid%rainncv           &
3493!     &             ,RAINC=grid%rainc    ,RAINNC=grid%rainnc               &
3494!     &             ,HFX=grid%hfx   ,SFCEVP=grid%sfcevp    ,LH=grid%lh     &
3495!     &             ,DT=grid%dt      ,SBW=config_flags%spec_bdy_width      &
3496!     &             ,XTIME=grid%xtime                                      &
3497!                 ! Selection flag
3498!     &             ,DIAG_PRINT=config_flags%diag_print                    &
3499!                 ! Dimension arguments
3500!     &             ,IDS=ids,IDE=ide, JDS=jds,JDE=jde, KDS=kds,KDE=kde     &
3501!     &             ,IMS=ims,IME=ime, JMS=jms,JME=jme, KMS=kms,KME=kme     &
3502!     &             ,IPS=ips,IPE=ipe, JPS=jps,JPE=jpe, KPS=kps,KPE=kpe     &
3503!     &             ,I_START=grid%i_start,I_END=min(grid%i_end, ide-1)     &
3504!     &             ,J_START=grid%j_start,J_END=min(grid%j_end, jde-1)     &
3505!     &             ,KTS=k_start, KTE=min(k_end,kde-1)                     &
3506!     &             ,NUM_TILES=grid%num_tiles                              &
3507!     &                                                          )
3508
3509#ifdef DM_PARALLEL
3510!-----------------------------------------------------------------------
3511! see above
3512!--------------------------------------------------------------
3513   CALL wrf_debug ( 200 , ' call HALO_RK_E' )
3514   IF      ( config_flags%h_mom_adv_order <= 4 ) THEN
3515#    include "HALO_EM_E_3.inc"
3516   ELSE IF ( config_flags%h_mom_adv_order <= 6 ) THEN
3517#    include "HALO_EM_E_5.inc"
3518   ELSE
3519     WRITE(wrf_err_message,*)'solve_em: invalid h_mom_adv_order = ',config_flags%h_mom_adv_order
3520     CALL wrf_error_fatal(TRIM(wrf_err_message))
3521   ENDIF
3522#endif
3523
3524#ifdef DM_PARALLEL
3525   if ( num_moist >= PARAM_FIRST_SCALAR  ) then
3526!-----------------------------------------------------------------------
3527! see above
3528!--------------------------------------------------------------
3529     CALL wrf_debug ( 200 , ' call HALO_RK_MOIST' )
3530     IF      ( config_flags%h_mom_adv_order <= 4 ) THEN
3531#      include "HALO_EM_MOIST_E_3.inc"
3532     ELSE IF ( config_flags%h_mom_adv_order <= 6 ) THEN
3533#      include "HALO_EM_MOIST_E_5.inc"
3534     ELSE
3535       WRITE(wrf_err_message,*)'solve_em: invalid h_mom_adv_order = ',config_flags%h_mom_adv_order
3536       CALL wrf_error_fatal(TRIM(wrf_err_message))
3537     ENDIF
3538   endif
3539   if ( num_chem >= PARAM_FIRST_SCALAR ) then
3540!-----------------------------------------------------------------------
3541! see above
3542!--------------------------------------------------------------
3543     CALL wrf_debug ( 200 , ' call HALO_RK_CHEM' )
3544     IF      ( config_flags%h_mom_adv_order <= 4 ) THEN
3545#      include "HALO_EM_CHEM_E_3.inc"
3546     ELSE IF ( config_flags%h_mom_adv_order <= 6 ) THEN
3547#      include "HALO_EM_CHEM_E_5.inc"
3548     ELSE
3549       WRITE(wrf_err_message,*)'solve_em: invalid h_mom_adv_order = ',config_flags%h_mom_adv_order
3550       CALL wrf_error_fatal(TRIM(wrf_err_message))
3551     ENDIF
3552   endif
3553   if ( num_scalar >= PARAM_FIRST_SCALAR ) then
3554!-----------------------------------------------------------------------
3555! see above
3556!--------------------------------------------------------------
3557     CALL wrf_debug ( 200 , ' call HALO_RK_SCALAR' )
3558     IF      ( config_flags%h_mom_adv_order <= 4 ) THEN
3559#      include "HALO_EM_SCALAR_E_3.inc"
3560     ELSE IF ( config_flags%h_mom_adv_order <= 6 ) THEN
3561#      include "HALO_EM_SCALAR_E_5.inc"
3562     ELSE
3563       WRITE(wrf_err_message,*)'solve_em: invalid h_mom_adv_order = ',config_flags%h_mom_adv_order
3564       CALL wrf_error_fatal(TRIM(wrf_err_message))
3565     ENDIF
3566   endif
3567#endif
3568
3569
3570!!!!!!!!!!!!!MARS
3571!!!!!!!!!!!!!MARS
3572!!
3573!! to get lighter output files, output the sum of the constant value phb (pb) with ph (p)
3574!! --- these variables were added to the Registry
3575!!
3576!!pressure
3577     grid%em_ptot  = grid%em_p + grid%em_pb
3578!!geopotential: already in php
3579     grid%em_phtot = grid%em_php
3580!!mass
3581     !grid%em_mutot = grid%em_mub + grid%em_mu_2
3582!!
3583!!!!!!!!!!!!!MARS
3584!!!!!!!!!!!!!MARS
3585
3586
3587
3588
3589   CALL wrf_debug ( 200 , ' call end of solve_em' )
3590
3591! Finish timers if compiled with -DBENCH.
3592#include <bench_solve_em_end.h>
3593
3594   RETURN
3595
3596END SUBROUTINE solve_em
Note: See TracBrowser for help on using the repository browser.