wrftladj/module_em_tl.F

   1 !WRF+/TL:MODEL_LAYER:DYNAMICS
   2 !
   3
   4 MODULE g_module_em
   5
   6
   7    USE module_model_constants
   8
   9    USE g_module_advect_em, only: g_advect_u, g_advect_v, g_advect_w, g_advect_scalar, g_advect_scalar_pd, g_advect_scalar_mono, &
  10         g_advect_weno_u, g_advect_weno_v, g_advect_weno_w, g_advect_scalar_weno,  g_advect_scalar_wenopd
  11
  12    USE g_module_big_step_utilities_em, only: grid_config_rec_type, g_calculate_full, g_couple_momentum, g_calc_mu_uv, g_calc_ww_cp, g_calc_cq, g_calc_alt, g_calc_php, g_set_tend, g_rhs_ph, &
  13         g_horizontal_pressure_gradient, g_pg_buoy_w, g_w_damp, g_perturbation_coriolis, g_coriolis, g_curvature, g_horizontal_diffusion, g_horizontal_diffusion_3dmp, g_vertical_diffusion_u, &
  14         g_vertical_diffusion_v, g_vertical_diffusion, g_vertical_diffusion_3dmp, g_sixth_order_diffusion, g_rk_rayleigh_damp, g_theta_relaxation, g_vertical_diffusion_mp, g_zero_tend, g_zero_tend2d
  15
  16    USE module_state_description, only: param_first_scalar, p_qr, p_qv, p_qc, p_qg, p_qi, p_qs, tiedtkescheme, ntiedtkescheme, heldsuarez, positivedef, &
  17        gdscheme, g3scheme, kfetascheme, mskfscheme, monotonic, wenopd_scalar, weno_scalar, weno_mom
  18
  19    !USE module_damping_em, only: held_suarez_damp
  20
  21 CONTAINS
  22
  23 !------------------------------------------------------------------------
  24
  25  SUBROUTINE g_rk_step_prep(config_flags,rk_step,u,g_u,v,g_v,w,g_w,t, &
  26  g_t,ph,g_ph,mu,g_mu,moist,g_moist,ru,g_ru,rv,g_rv,rw,g_rw,ww, &
  27  g_ww,php,g_php,alt,g_alt,muu,g_muu,muv,g_muv,mub,mut,g_mut,phb,pb, &
  28  p,g_p,al,g_al,alb,cqu,g_cqu,cqv,g_cqv,cqw,g_cqw,msfux,msfuy,msfvx, &
  29  msfvx_inv,msfvy,msftx,msfty,fnm,fnp,dnw,rdx,rdy,n_moist,ids,ide,jds,jde,kds,kde,ims, &
  30  ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
  31
  32  IMPLICIT NONE
  33
  34  REAL :: Tmpv1,g_Tmpv1
  35  TYPE(grid_config_rec_type) :: config_flags
  36  INTEGER :: ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte
  37  INTEGER :: n_moist,rk_step
  38  REAL :: rdx,rdy
  39  REAL,DIMENSION(ims:ime,kms:kme,jms:jme) :: u,g_u,v,g_v,w,g_w,t,g_t,ph, &
  40  g_ph,phb,pb,al,g_al,alb
  41  REAL,DIMENSION(ims:ime,kms:kme,jms:jme) :: ru,g_ru,rv,g_rv,rw,g_rw,ww, &
  42  g_ww,php,g_php,cqu,g_cqu,cqv,g_cqv,cqw,g_cqw,alt,g_alt
  43  REAL,DIMENSION(ims:ime,kms:kme,jms:jme) :: p,g_p
  44  REAL,DIMENSION(ims:ime,kms:kme,jms:jme,n_moist) :: moist,g_moist
  45  REAL,DIMENSION(ims:ime,jms:jme) :: msftx,msfty,msfux,msfuy,msfvx,msfvx_inv,msfvy,mu, &
  46  g_mu,mub
  47  REAL,DIMENSION(ims:ime,jms:jme) :: muu,g_muu,muv,g_muv,mut,g_mut
  48  REAL,DIMENSION(kms:kme) :: fnm,fnp,dnw
  49  integer :: k
  50
  51  CALL g_calculate_full(mut,g_mut,mub,mu,g_mu,ids,ide,jds,jde,1,2,ims,ime,jms, &
  52  jme,1,1,its,ite,jts,jte,1,1)
  53
  54  CALL g_calc_mu_uv(config_flags,mu,g_mu,mub,muu,g_muu,muv,g_muv,ids,ide, &
  55  jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
  56
  57  CALL g_couple_momentum(muu,g_muu,ru,g_ru,u,g_u,msfuy,muv,g_muv,rv, &
  58  g_rv,v,g_v,msfvx,msfvx_inv,mut,g_mut,rw,g_rw,w,g_w,msfty,ids,ide,jds, &
  59  jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
  60
  61  CALL g_calc_ww_cp(u,g_u,v,g_v,mu,g_mu,mub,ww,g_ww,rdx,rdy,msftx,msfty, &
  62  msfux,msfuy,msfvx,msfvx_inv,msfvy,dnw,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms, &
  63  kme,its,ite,jts,jte,kts,kte)
  64
  65  CALL g_calc_cq(moist,g_moist,cqu,g_cqu,cqv,g_cqv,cqw,g_cqw,n_moist, &
  66  ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
  67
  68  CALL g_calc_alt(alt,g_alt,al,g_al,alb,ids,ide,jds,jde,kds,kde,ims,ime,jms, &
  69  jme,kms,kme,its,ite,jts,jte,kts,kte)
  70
  71  CALL g_calc_php(php,g_php,ph,g_ph,phb,ids,ide,jds,jde,kds,kde,ims,ime,jms, &
  72  jme,kms,kme,its,ite,jts,jte,kts,kte)
  73
  74  END SUBROUTINE g_rk_step_prep
  75
  76 !-------------------------------------------------------------------------------
  77
  78  SUBROUTINE g_rk_tendency(config_flags,rk_step,ru_tend,g_ru_tend,rv_tend, &
  79  g_rv_tend,rw_tend,g_rw_tend,ph_tend,g_ph_tend,t_tend,g_t_tend,ru_tendf, &
  80  g_ru_tendf,rv_tendf,g_rv_tendf,rw_tendf,g_rw_tendf,ph_tendf,g_ph_tendf, &
  81  t_tendf,g_t_tendf,mu_tend,g_mu_tend,u_save,g_u_save,v_save,g_v_save, &
  82  w_save,g_w_save,ph_save,g_ph_save,t_save,g_t_save,mu_save,g_mu_save, &
  83  RTHFTEN,g_RTHFTEN,ru,g_ru,rv,g_rv,rw,g_rw,ww,g_ww,u,g_u,v,g_v,w, &
  84  g_w,t,g_t,ph,g_ph,u_old,g_u_old,v_old,g_v_old,w_old,g_w_old,t_old, &
  85 ! Revised by Ning Pan, 2010-07-29
  86 ! g_t_old,ph_old,g_ph_old,h_diabatic,g_h_diabatic,phb,t_init,g_t_init,mu, &
  87  g_t_old,ph_old,g_ph_old,h_diabatic,g_h_diabatic,phb,t_init,mu, &
  88  g_mu,mut,g_mut,muu,g_muu,muv,g_muv,mub,al,g_al,alt,g_alt,p,g_p, &
  89  pb,php,g_php,cqu,g_cqu,cqv,g_cqv,cqw,g_cqw,u_base,v_base,t_base,qv_base, &
  90 ! Revised by Ning Pan, 2010-07-29
  91 ! z_base,msfux,msfuy,msfvx,msfvx_inv,msfvy,msftx,msfty,xlat,g_xlat,f,e,sina,cosa, &
  92 ! fnm,fnp,rdn,rdnw,dt,rdx,rdy,khdif,kvdif,xkmhd,g_xkmhd,xkhh,g_xkhh,diff_6th_opt, &
  93 ! diff_6th_factor,g_diff_6th_factor,dampcoef,g_dampcoef,zdamp,g_zdamp, &
  94 ! damp_opt,cf1,cf2,cf3,cfn,cfn1,n_moist,non_hydrostatic,top_lid,u_frame,g_u_frame, &
  95 ! v_frame,g_v_frame,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte, &
  96 ! kts,kte,max_vert_cfl,g_max_vert_cfl,max_horiz_cfl,g_max_horiz_cfl)
  97  z_base,msfux,msfuy,msfvx,msfvx_inv,msfvy,msftx,msfty,xlat,f,e,sina,cosa, &
  98  fnm,fnp,rdn,rdnw,dt,rdx,rdy,khdif,kvdif,xkmhd,g_xkmhd,xkhh,g_xkhh,diff_6th_opt, &
  99  diff_6th_factor,adv_opt,dampcoef,zdamp, &
 100  damp_opt,rad_nudge,cf1,cf2,cf3,cfn,cfn1,n_moist,non_hydrostatic,top_lid,u_frame, &
 101  v_frame,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte, &
 102  kts,kte,max_vert_cfl,max_horiz_cfl)
 103
 104  IMPLICIT NONE
 105
 106  REAL :: Tmpv1,g_Tmpv1
 107  TYPE(grid_config_rec_type) :: config_flags
 108  INTEGER :: ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte
 109  LOGICAL :: non_hydrostatic,top_lid
 110  INTEGER :: n_moist,rk_step
 111  REAL,DIMENSION(ims:ime,kms:kme,jms:jme) :: ru,g_ru,rv,g_rv,rw,g_rw,ww, &
 112  g_ww,u,g_u,v,g_v,w,g_w,t,g_t,ph,g_ph,u_old,g_u_old,v_old, &
 113  g_v_old,w_old,g_w_old,t_old,g_t_old,ph_old,g_ph_old,phb,al,g_al,alt, &
 114 ! Revised by Ning Pan, 2010-07-29
 115 ! g_alt,p,g_p,pb,php,g_php,cqu,g_cqu,cqv,g_cqv,t_init,g_t_init,xkmhd, &
 116  g_alt,p,g_p,pb,php,g_php,cqu,g_cqu,cqv,g_cqv,t_init,xkmhd, &
 117  g_xkmhd,xkhh,g_xkhh,h_diabatic,g_h_diabatic
 118  REAL,DIMENSION(ims:ime,kms:kme,jms:jme) :: ru_tend,g_ru_tend,rv_tend,g_rv_tend, &
 119  rw_tend,g_rw_tend,t_tend,g_t_tend,ph_tend,g_ph_tend,RTHFTEN,g_RTHFTEN, &
 120  u_save,g_u_save,v_save,g_v_save,w_save,g_w_save,ph_save,g_ph_save,t_save, &
 121  g_t_save
 122  REAL,DIMENSION(ims:ime,kms:kme,jms:jme) :: ru_tendf,g_ru_tendf,rv_tendf, &
 123  g_rv_tendf,rw_tendf,g_rw_tendf,t_tendf,g_t_tendf,ph_tendf,g_ph_tendf,cqw,g_cqw
 124  REAL,DIMENSION(ims:ime,jms:jme) :: mu_tend,g_mu_tend,mu_save,g_mu_save
 125  REAL,DIMENSION(ims:ime,jms:jme) :: msfux,msfuy,msfvx,msfvx_inv,msfvy,msftx,msfty, &
 126 ! Revised by Ning Pan, 2010-07-30
 127 ! xlat,g_xlat,f,e,sina,cosa,mu,g_mu,mut,g_mut,mub,muu,g_muu,muv,g_muv
 128  xlat,f,e,sina,cosa,mu,g_mu,mut,g_mut,mub,muu,g_muu,muv,g_muv
 129  REAL,DIMENSION(kms:kme) :: fnm,fnp,rdn,rdnw,u_base,v_base,t_base,qv_base,z_base
 130 ! Revised by Ning Pan, 2010-07-29
 131 ! REAL :: rdx,rdy,dt,u_frame,g_u_frame,v_frame,g_v_frame,khdif,kvdif
 132  REAL :: rdx,rdy,dt,u_frame,v_frame,khdif,kvdif
 133  INTEGER :: diff_6th_opt
 134 ! Revised by Ning Pan, 2010-07-29
 135 ! REAL :: diff_6th_factor,g_diff_6th_factor
 136  REAL :: diff_6th_factor
 137  INTEGER :: adv_opt
 138  INTEGER :: damp_opt,rad_nudge
 139 ! Revised by Ning Pan, 2010-07-29
 140 ! REAL :: zdamp,g_zdamp,dampcoef,g_dampcoef
 141 ! REAL :: max_horiz_cfl,g_max_horiz_cfl
 142 ! REAL :: max_vert_cfl,g_max_vert_cfl
 143  REAL :: zdamp,dampcoef
 144  REAL :: max_horiz_cfl
 145  REAL :: max_vert_cfl
 146 ! Revised by Ning Pan, 2010-07-29
 147 ! REAL :: kdift,g_kdift,khdq,g_khdq,kvdq,g_kvdq,cfn,cfn1,cf1,cf2,cf3
 148  REAL :: kdift,khdq,kvdq,cfn,cfn1,cf1,cf2,cf3
 149  INTEGER :: i,j,k
 150  INTEGER :: time_step
 151
 152  CALL g_zero_tend(ru_tend,g_ru_tend,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms, &
 153  kme,its,ite,jts,jte,kts,kte)
 154
 155  CALL g_zero_tend(rv_tend,g_rv_tend,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms, &
 156  kme,its,ite,jts,jte,kts,kte)
 157
 158  CALL g_zero_tend(rw_tend,g_rw_tend,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms, &
 159  kme,its,ite,jts,jte,kts,kte)
 160
 161  CALL g_zero_tend(t_tend,g_t_tend,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms, &
 162  kme,its,ite,jts,jte,kts,kte)
 163
 164  CALL g_zero_tend(ph_tend,g_ph_tend,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms, &
 165  kme,its,ite,jts,jte,kts,kte)
 166
 167  CALL g_zero_tend(u_save,g_u_save,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms, &
 168  kme,its,ite,jts,jte,kts,kte)
 169
 170  CALL g_zero_tend(v_save,g_v_save,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms, &
 171  kme,its,ite,jts,jte,kts,kte)
 172
 173  CALL g_zero_tend(w_save,g_w_save,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms, &
 174  kme,its,ite,jts,jte,kts,kte)
 175
 176  CALL g_zero_tend(ph_save,g_ph_save,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms, &
 177  kme,its,ite,jts,jte,kts,kte)
 178
 179  CALL g_zero_tend(t_save,g_t_save,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms, &
 180  kme,its,ite,jts,jte,kts,kte)
 181
 182  CALL g_zero_tend2d(mu_tend,g_mu_tend,ids,ide,jds,jde,1,1,ims,ime,jms,jme,1,1,its, &
 183  ite,jts,jte,1,1)
 184
 185  CALL g_zero_tend2d(mu_save,g_mu_save,ids,ide,jds,jde,1,1,ims,ime,jms,jme,1,1,its, &
 186  ite,jts,jte,1,1)
 187
 188 !This line is fail to be recognized
 189     CALL nl_get_time_step ( 1, time_step )
 190
 191  IF( (rk_step == 3) .and. ( adv_opt == WENO_MOM ) ) THEN
 192
 193  CALL g_advect_weno_u ( u, g_u, u, g_u, ru_tend, &
 194                 g_ru_tend, ru, g_ru, rv, g_rv, ww, g_ww, &
 195                 mut, g_mut, time_step, config_flags, &
 196                 msfux, msfuy, msfvx, msfvy,   &
 197                 msftx, msfty,                 &
 198                 fnm, fnp, rdx, rdy, rdnw,     &
 199                 ids, ide, jds, jde, kds, kde, &
 200                 ims, ime, jms, jme, kms, kme, &
 201                 its, ite, jts, jte, kts, kte )
 202
 203  ELSE
 204
 205  CALL g_advect_u(u,g_u,u,g_u,ru_tend,g_ru_tend,ru,g_ru,rv,g_rv,ww, &
 206  g_ww,mut,g_mut,time_step,config_flags,msfux,msfuy,msfvx,msfvy,msftx,msfty,fnm, &
 207  fnp,rdx,rdy,rdnw,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
 208
 209  ENDIF
 210
 211  IF( (rk_step == 3) .and. ( adv_opt == WENO_MOM ) ) THEN
 212
 213  CALL g_advect_weno_v ( v, g_v, v, g_v, rv_tend, g_rv_tend, &
 214                 ru, g_ru, rv, g_rv, ww, g_ww, &
 215                 mut, g_mut, time_step, config_flags, &
 216                 msfux, msfuy, msfvx, msfvy,   &
 217                 msftx, msfty,                 &
 218                 fnm, fnp, rdx, rdy, rdnw,     &
 219                 ids, ide, jds, jde, kds, kde, &
 220                 ims, ime, jms, jme, kms, kme, &
 221                 its, ite, jts, jte, kts, kte )
 222
 223  ELSE
 224
 225  CALL g_advect_v(v,g_v,v,g_v,rv_tend,g_rv_tend,ru,g_ru,rv,g_rv,ww, &
 226  g_ww,mut,g_mut,time_step,config_flags,msfux,msfuy,msfvx,msfvy,msftx,msfty,fnm, &
 227  fnp,rdx,rdy,rdnw,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
 228
 229  ENDIF
 230
 231  IF(non_hydrostatic) THEN
 232  IF( (rk_step == 3) .and. ( adv_opt == WENO_MOM ) ) THEN
 233  CALL g_advect_weno_w ( w, g_w, w, g_w, rw_tend, g_rw_tend, &
 234                  ru, g_ru, rv, g_rv, ww, g_ww, &
 235                  mut, time_step, config_flags, &
 236                  msfux, msfuy, msfvx, msfvy,   &
 237                  msftx, msfty,                 &
 238                  fnm, fnp, rdx, rdy, rdn,      &
 239                  ids, ide, jds, jde, kds, kde, &
 240                  ims, ime, jms, jme, kms, kme, &
 241                  its, ite, jts, jte, kts, kte )
 242
 243  ELSE
 244
 245  CALL g_advect_w(w,g_w,w,g_w,rw_tend,g_rw_tend,ru, &
 246  g_ru,rv,g_rv,ww,g_ww,mut,time_step,config_flags,msfux,msfuy,msfvx, &
 247  msfvy,msftx,msfty,fnm,fnp,rdx,rdy,rdn,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms, &
 248  kme,its,ite,jts,jte,kts,kte)
 249
 250  ENDIF
 251  ENDIF
 252
 253 !  theta flux divergence
 254 !hcl 11/2016 ERM: Use WENO for theta flux on 3rd RK step if using WENO_SCALAR or WENOPD_SCALAR
 255 ! to be consistent with other scalar fluxes
 256  IF(  ( config_flags%scalar_adv_opt == WENO_SCALAR           &
 257          .or. config_flags%scalar_adv_opt == WENOPD_SCALAR   &
 258          .or. config_flags%moist_adv_opt == WENO_SCALAR      &
 259          .or. config_flags%moist_adv_opt == WENOPD_SCALAR    &
 260                        )  .and. (rk_step == 3) ) THEN
 261
 262  CALL g_advect_scalar_weno(t,g_t,t,g_t,t_tend,g_t_tend,ru,g_ru,rv,g_rv, &
 263  ww,g_ww,mut,time_step,config_flags,msfux,msfuy,msfvx,msfvy,msftx,msfty, &
 264  fnm,fnp,rdx,rdy,rdnw,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
 265
 266  ELSE
 267
 268  CALL g_advect_scalar(t,g_t,t,g_t,t_tend,g_t_tend,ru,g_ru,rv,g_rv, &
 269  ww,g_ww,mut,time_step,config_flags,msfux,msfuy,msfvx,msfvy,msftx,msfty, &
 270  fnm,fnp,rdx,rdy,rdnw,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
 271
 272  ENDIF
 273
 274  IF( config_flags%cu_physics == GDSCHEME  .OR.       &
 275      config_flags%cu_physics == G3SCHEME  .OR.       &
 276      config_flags%cu_physics == NTIEDTKESCHEME ) THEN
 277
 278  CALL g_set_tend(RTHFTEN,g_RTHFTEN,t_tend,g_t_tend,msfty,ids,ide,jds,jde,kds, &
 279  kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
 280  END IF
 281
 282  CALL g_rhs_ph(ph_tend,g_ph_tend,u,g_u,v,g_v,ww,g_ww,ph,g_ph,ph, &
 283  g_ph,phb,w,g_w,mut,g_mut,muu,g_muu,muv,g_muv,fnm,fnp,rdnw,cfn,cfn1, &
 284  rdx,rdy,msfux,msfuy,msfvx,msfvx_inv,msfvy,msftx,msfty,non_hydrostatic,config_flags, &
 285  ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
 286
 287  CALL g_horizontal_pressure_gradient(ru_tend,g_ru_tend,rv_tend,g_rv_tend,ph, &
 288  g_ph,alt,g_alt,p,g_p,pb,al,g_al,php,g_php,cqu,g_cqu,cqv,g_cqv, &
 289  muu,g_muu,muv,g_muv,mu,g_mu,fnm,fnp,rdnw,cf1,cf2,cf3,cfn,cfn1,rdx,rdy,msfux,msfuy, &
 290  msfvx,msfvy,msftx,msfty,config_flags,non_hydrostatic,top_lid,ids,ide,jds,jde,kds,kde, &
 291  ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
 292
 293  IF(non_hydrostatic) THEN
 294
 295  CALL g_pg_buoy_w(rw_tend,g_rw_tend,p,g_p,cqw,g_cqw,mu,g_mu,mub,rdnw, &
 296  rdn,g,msftx,msfty,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
 297  ENDIF
 298
 299 ! Revised by Ning Pan, 2010-07-30
 300 ! CALL g_w_damp(rw_tend,g_rw_tend,max_vert_cfl,g_max_vert_cfl,max_horiz_cfl, &
 301 ! g_max_horiz_cfl,u,g_u,v,g_v,ww,g_ww,w,g_w,mut,g_mut,rdnw,rdx,rdy, &
 302  CALL g_w_damp(rw_tend,g_rw_tend,max_vert_cfl,max_horiz_cfl, &
 303  u,g_u,v,g_v,ww,g_ww,w,g_w,mut,g_mut,rdnw,rdx,rdy, &
 304  msfux,msfuy,msfvx,msfvy,dt,config_flags,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms, &
 305  kme,its,ite,jts,jte,kts,kte)
 306
 307  IF(config_flags%pert_coriolis) THEN
 308
 309  CALL g_perturbation_coriolis(ru,g_ru,rv,g_rv,rw,g_rw,ru_tend, &
 310  g_ru_tend,rv_tend,g_rv_tend,rw_tend,g_rw_tend,config_flags,u_base,v_base, &
 311  z_base,muu,g_muu,muv,g_muv,phb,ph,g_ph,msftx,msfty,msfux,msfuy,msfvx,msfvy, &
 312  f,e,sina,cosa,fnm,fnp,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite,jts, &
 313  jte,kts,kte)
 314  ELSE
 315
 316  CALL g_coriolis(ru,g_ru,rv,g_rv,rw,g_rw,ru_tend,g_ru_tend,rv_tend, &
 317  g_rv_tend,rw_tend,g_rw_tend,config_flags,msftx,msfty,msfux,msfuy,msfvx,msfvy,f, &
 318  e,sina,cosa,fnm,fnp,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
 319  END IF
 320
 321  CALL g_curvature(ru,g_ru,rv,g_rv,rw,g_rw,u,g_u,v,g_v,w, &
 322  ru_tend,g_ru_tend,rv_tend,g_rv_tend,rw_tend,g_rw_tend,config_flags,msfux, &
 323 ! Revised by Ning Pan, 2010-07-30: xlat is a constant array
 324 ! msfuy,msfvx,msfvy,msftx,msfty,xlat,g_xlat,fnm,fnp,rdx,rdy,ids,ide,jds,jde,kds,kde, &
 325  msfuy,msfvx,msfvy,msftx,msfty,xlat,fnm,fnp,rdx,rdy,ids,ide,jds,jde,kds,kde, &
 326  ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
 327
 328  IF(config_flags%ra_lw_physics == HELDSUAREZ) THEN
 329 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
 330 !!! Reamarked by Ning Pan, 2010-07-30 : JUST FOR DEBUGGING DYNAMICS OF WRF+   !!!
 331 !!!              REMARK SHOULD BE REMOVED WHEN CONSTRUCTING PHYSICS OF WRF+   !!!
 332 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
 333 ! CALL g_held_suarez_damp(ru_tend,g_ru_tend,rv_tend,g_rv_tend,ru,g_ru,rv, &
 334 ! g_rv,p,g_p,pb,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
 335  END IF
 336
 337  IF( rk_step == 1 ) THEN
 338
 339  IF(config_flags%diff_opt .eq. 1) THEN
 340
 341  CALL g_horizontal_diffusion('u',u,g_u,ru_tendf,g_ru_tendf,mut,g_mut, &
 342  config_flags,msfux,msfuy,msfvx,msfvx_inv,msfvy,msftx,msfty,khdif,xkmhd,g_xkmhd, &
 343  rdx,rdy,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
 344
 345  CALL g_horizontal_diffusion('v',v,g_v,rv_tendf,g_rv_tendf,mut,g_mut, &
 346  config_flags,msfux,msfuy,msfvx,msfvx_inv,msfvy,msftx,msfty,khdif,xkmhd,g_xkmhd, &
 347  rdx,rdy,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
 348
 349  CALL g_horizontal_diffusion('w',w,g_w,rw_tendf,g_rw_tendf,mut,g_mut, &
 350  config_flags,msfux,msfuy,msfvx,msfvx_inv,msfvy,msftx,msfty,khdif,xkmhd,g_xkmhd, &
 351  rdx,rdy,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
 352
 353 ! g_khdq =0.0  ! Remarked by Ning Pan, 2010-07-30
 354  khdq =3. *khdif
 355
 356  CALL g_horizontal_diffusion_3dmp('m',t,g_t,t_tendf,g_t_tendf,mut,g_mut, &
 357 ! Revised by Ning Pan, 2010-07-30
 358 ! config_flags,t_init,g_t_init,msfux,msfuy,msfvx,msfvx_inv,msfvy,msftx,msfty,khdq, &
 359 ! g_khdq,xkhh,g_xkhh,rdx,rdy,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its, &
 360  config_flags,t_init,msfux,msfuy,msfvx,msfvx_inv,msfvy,msftx,msfty,khdq, &
 361  xkhh,g_xkhh,rdx,rdy,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its, &
 362  ite,jts,jte,kts,kte)
 363
 364  IF(config_flags%bl_pbl_physics .eq. 0) THEN
 365
 366  CALL g_vertical_diffusion_u(u,g_u,ru_tendf,g_ru_tendf,config_flags,u_base, &
 367  alt,g_alt,muu,g_muu,rdn,rdnw,kvdif,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms, &
 368  kme,its,ite,jts,jte,kts,kte)
 369
 370  CALL g_vertical_diffusion_v(v,g_v,rv_tendf,g_rv_tendf,config_flags,v_base, &
 371  alt,g_alt,muv,g_muv,rdn,rdnw,kvdif,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms, &
 372  kme,its,ite,jts,jte,kts,kte)
 373
 374  IF(non_hydrostatic) CALL g_vertical_diffusion('w',w,g_w,rw_tendf,g_rw_tendf, &
 375  config_flags,alt,g_alt,mut,g_mut,rdn,rdnw,kvdif,ids,ide,jds,jde,kds,kde,ims, &
 376  ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
 377
 378 ! g_kvdq =0.0  ! Remarked by Ning Pan, 2010-07-30
 379  kvdq =3. *kvdif
 380
 381  CALL g_vertical_diffusion_3dmp(t,g_t,t_tendf,g_t_tendf,config_flags,t_init, &
 382 ! Revised by Ning Pan, 2010-07-30
 383 ! g_t_init,alt,g_alt,mut,g_mut,rdn,rdnw,kvdq,g_kvdq,ids,ide,jds,jde,kds, &
 384  alt,g_alt,mut,g_mut,rdn,rdnw,kvdq,ids,ide,jds,jde,kds, &
 385  kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
 386  ENDIF
 387
 388  END IF
 389
 390  IF( diff_6th_opt .NE. 0 ) THEN
 391
 392  CALL g_sixth_order_diffusion('u',u,g_u,ru_tendf,g_ru_tendf,mut,g_mut,dt, &
 393 ! Revised by Ning Pan, 2010-07-30
 394 ! config_flags,diff_6th_opt,diff_6th_factor,g_diff_6th_factor,ids,ide,jds,jde,kds, &
 395  config_flags,diff_6th_opt,diff_6th_factor,ids,ide,jds,jde,kds, &
 396  kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
 397
 398  CALL g_sixth_order_diffusion('v',v,g_v,rv_tendf,g_rv_tendf,mut,g_mut,dt, &
 399 ! Revised by Ning Pan, 2010-07-30
 400 ! config_flags,diff_6th_opt,diff_6th_factor,g_diff_6th_factor,ids,ide,jds,jde,kds, &
 401  config_flags,diff_6th_opt,diff_6th_factor,ids,ide,jds,jde,kds, &
 402  kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
 403
 404  IF(non_hydrostatic) CALL g_sixth_order_diffusion('w',w,g_w,rw_tendf, &
 405  g_rw_tendf,mut,g_mut,dt,config_flags,diff_6th_opt,diff_6th_factor, &
 406 ! Revised by Ning Pan, 2010-07-30
 407 ! g_diff_6th_factor,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
 408  ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
 409
 410  CALL g_sixth_order_diffusion('m',t,g_t,t_tendf,g_t_tendf,mut,g_mut,dt, &
 411 ! Revised by Ning Pan, 2010-07-30
 412 ! config_flags,diff_6th_opt,diff_6th_factor,g_diff_6th_factor,ids,ide,jds,jde,kds, &
 413  config_flags,diff_6th_opt,diff_6th_factor,ids,ide,jds,jde,kds, &
 414  kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
 415  ENDIF
 416
 417  IF( damp_opt .eq. 2 ) CALL g_rk_rayleigh_damp(ru_tendf,g_ru_tendf,rv_tendf, &
 418  g_rv_tendf,rw_tendf,g_rw_tendf,t_tendf,g_t_tendf,u,g_u,v,g_v,w,g_w, &
 419 ! Revised by Ning Pan, 2010-07-23
 420 ! t,g_t,t_init,g_t_init,mut,g_mut,muu,g_muu,muv,g_muv,ph,g_ph,phb, &
 421 ! u_base,v_base,t_base,z_base,dampcoef,g_dampcoef,zdamp,g_zdamp,ids,ide,jds,jde, &
 422  t,g_t,t_init,mut,g_mut,muu,g_muu,muv,g_muv,ph,g_ph,phb, &
 423  u_base,v_base,t_base,z_base,dampcoef,zdamp,ids,ide,jds,jde, &
 424  kds,kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
 425
 426     IF( rad_nudge .eq. 1 )                                     &
 427        CALL g_theta_relaxation( t_tendf, g_t_tendf, t, g_t, t_init, &
 428                               mut, g_mut, ph, g_ph, phb,       &
 429                               t_base, z_base,                  &
 430                               ids, ide, jds, jde, kds, kde,    &
 431                               ims, ime, jms, jme, kms, kme,    &
 432                               its, ite, jts, jte, kts, kte   )
 433
 434   END IF
 435
 436  END SUBROUTINE g_rk_tendency
 437
 438 !-------------------------------------------------------------------------------
 439
 440 !        Generated by TAPENADE     (INRIA, Tropics team)
 441 !  Tapenade 3.5 (r3785) - 22 Mar 2011 18:35
 442 !
 443 !  Differentiation of rk_addtend_dry in forward (tangent) mode:
 444 !   variations   of useful results: ph_tendf rw_tendf mu_tend rv_tendf
 445 !                ru_tend rw_tend ru_tendf t_tend rv_tend t_tendf
 446 !                ph_tend
 447 !   with respect to varying inputs: ph_tendf rw_tendf u_save ph_save
 448 !                w_save mu_tend rv_tendf ru_tend rw_tend h_diabatic
 449 !                ru_tendf t_tend mu_tendf t_save v_save rv_tend
 450 !                t_tendf mut ph_tend
 451 !   RW status of diff variables: ph_tendf:in-out rw_tendf:in-out
 452 !                u_save:in ph_save:in w_save:in mu_tend:in-out
 453 !                rv_tendf:in-out ru_tend:in-out rw_tend:in-out
 454 !                h_diabatic:in ru_tendf:in-out t_tend:in-out mu_tendf:in
 455 !                t_save:in v_save:in rv_tend:in-out t_tendf:in-out
 456 !                mut:in ph_tend:in-out
 457 SUBROUTINE G_RK_ADDTEND_DRY(ru_tend, ru_tendd, rv_tend, rv_tendd, &
 458 &  rw_tend, rw_tendd, ph_tend, ph_tendd, t_tend, t_tendd, ru_tendf, &
 459 &  ru_tendfd, rv_tendf, rv_tendfd, rw_tendf, rw_tendfd, ph_tendf, &
 460 &  ph_tendfd, t_tendf, t_tendfd, u_save, u_saved, v_save, v_saved, w_save&
 461 &  , w_saved, ph_save, ph_saved, t_save, t_saved, mu_tend, mu_tendd, &
 462 &  mu_tendf, mu_tendfd, rk_step, h_diabatic, h_diabaticd, mut, mutd, &
 463 &  msftx, msfty, msfux, msfuy, msfvx, msfvx_inv, msfvy, ids, ide, jds, &
 464 &  jde, kds, kde, ims, ime, jms, jme, kms, kme, ips, ipe, jps, jpe, kps, &
 465 &  kpe, its, ite, jts, jte, kts, kte)
 466   IMPLICIT NONE
 467 !  Input data.
 468   INTEGER, INTENT(IN) :: ids, ide, jds, jde, kds, kde, ims, ime, jms, &
 469 &  jme, kms, kme, ips, ipe, jps, jpe, kps, kpe, its, ite, jts, jte, kts, &
 470 &  kte
 471   INTEGER, INTENT(IN) :: rk_step
 472   REAL, DIMENSION(ims:ime, kms:kme, jms:jme), INTENT(INOUT) :: ru_tend, &
 473 &  rv_tend, rw_tend, ph_tend, t_tend, ru_tendf, rv_tendf, rw_tendf, &
 474 &  ph_tendf, t_tendf
 475   REAL, DIMENSION(ims:ime, kms:kme, jms:jme), INTENT(INOUT) :: ru_tendd&
 476 &  , rv_tendd, rw_tendd, ph_tendd, t_tendd, ru_tendfd, rv_tendfd, &
 477 &  rw_tendfd, ph_tendfd, t_tendfd
 478   REAL, DIMENSION(ims:ime, jms:jme), INTENT(INOUT) :: mu_tend, mu_tendf
 479   REAL, DIMENSION(ims:ime, jms:jme), INTENT(INOUT) :: mu_tendd, &
 480 &  mu_tendfd
 481   REAL, DIMENSION(ims:ime, kms:kme, jms:jme), INTENT(IN) :: u_save, &
 482 &  v_save, w_save, ph_save, t_save, h_diabatic
 483   REAL, DIMENSION(ims:ime, kms:kme, jms:jme), INTENT(IN) :: u_saved, &
 484 &  v_saved, w_saved, ph_saved, t_saved, h_diabaticd
 485   REAL, DIMENSION(ims:ime, jms:jme), INTENT(IN) :: mut, msftx, msfty, &
 486 &  msfux, msfuy, msfvx, msfvx_inv, msfvy
 487   REAL, DIMENSION(ims:ime, jms:jme), INTENT(IN) :: mutd
 488 ! Local
 489   INTEGER :: i, j, k
 490   INTEGER :: min8
 491   INTEGER :: min7
 492   INTEGER :: min6
 493   INTEGER :: min5
 494   INTEGER :: min4
 495   INTEGER :: min3
 496   INTEGER :: min2
 497   INTEGER :: min1
 498   IF (jte .GT. jde - 1) THEN
 499     min1 = jde - 1
 500   ELSE
 501     min1 = jte
 502   END IF
 503 !<DESCRIPTION>
 504 !
 505 ! rk_addtend_dry constructs the full large-timestep tendency terms for
 506 ! momentum (u,v,w), theta and geopotential equations.   This is accomplished
 507 ! by combining the physics tendencies (in *tendf; these are computed
 508 ! the first RK substep, held fixed thereafter) with the RK tendencies
 509 ! (in *tend, these include advection, pressure gradient, etc;
 510 ! these change each rk substep).  Output is in *tend.
 511 !
 512 !</DESCRIPTION>
 513 !  Finally, add the forward-step tendency to the rk_tendency
 514 ! u/v/w/save contain bc tendency that needs to be multiplied by msf
 515 ! (u by msfuy, v by msfvx)
 516 !  before adding it to physics tendency (*tendf)
 517 ! For momentum we need the final tendency to include an inverse msf
 518 ! physics/bc tendency needs to be divided, advection tendency already has it
 519 ! For scalars we need the final tendency to include an inverse msf (msfty)
 520 ! advection tendency is OK, physics/bc tendency needs to be divided by msf
 521   DO j=jts,min1
 522     DO k=kts,kte-1
 523       DO i=its,ite
 524 ! multiply by my to uncouple u
 525         IF (rk_step .EQ. 1) THEN
 526           ru_tendfd(i, k, j) = ru_tendfd(i, k, j) + msfuy(i, j)*u_saved(&
 527 &            i, k, j)
 528           ru_tendf(i, k, j) = ru_tendf(i, k, j) + u_save(i, k, j)*msfuy(&
 529 &            i, j)
 530         END IF
 531 ! divide by my to couple u
 532         ru_tendd(i, k, j) = ru_tendd(i, k, j) + ru_tendfd(i, k, j)/msfuy&
 533 &          (i, j)
 534         ru_tend(i, k, j) = ru_tend(i, k, j) + ru_tendf(i, k, j)/msfuy(i&
 535 &          , j)
 536       END DO
 537     END DO
 538   END DO
 539   DO j=jts,jte
 540     DO k=kts,kte-1
 541       IF (ite .GT. ide - 1) THEN
 542         min2 = ide - 1
 543       ELSE
 544         min2 = ite
 545       END IF
 546       DO i=its,min2
 547 ! multiply by mx to uncouple v
 548         IF (rk_step .EQ. 1) THEN
 549           rv_tendfd(i, k, j) = rv_tendfd(i, k, j) + msfvx(i, j)*v_saved(&
 550 &            i, k, j)
 551           rv_tendf(i, k, j) = rv_tendf(i, k, j) + v_save(i, k, j)*msfvx(&
 552 &            i, j)
 553         END IF
 554 ! divide by mx to couple v
 555         rv_tendd(i, k, j) = rv_tendd(i, k, j) + msfvx_inv(i, j)*&
 556 &          rv_tendfd(i, k, j)
 557         rv_tend(i, k, j) = rv_tend(i, k, j) + rv_tendf(i, k, j)*&
 558 &          msfvx_inv(i, j)
 559       END DO
 560     END DO
 561   END DO
 562   IF (jte .GT. jde - 1) THEN
 563     min3 = jde - 1
 564   ELSE
 565     min3 = jte
 566   END IF
 567   DO j=jts,min3
 568     DO k=kts,kte
 569       IF (ite .GT. ide - 1) THEN
 570         min4 = ide - 1
 571       ELSE
 572         min4 = ite
 573       END IF
 574       DO i=its,min4
 575 ! multiply by my to uncouple w
 576         IF (rk_step .EQ. 1) THEN
 577           rw_tendfd(i, k, j) = rw_tendfd(i, k, j) + msfty(i, j)*w_saved(&
 578 &            i, k, j)
 579           rw_tendf(i, k, j) = rw_tendf(i, k, j) + w_save(i, k, j)*msfty(&
 580 &            i, j)
 581         END IF
 582 ! divide by my to couple w
 583         rw_tendd(i, k, j) = rw_tendd(i, k, j) + rw_tendfd(i, k, j)/msfty&
 584 &          (i, j)
 585         rw_tend(i, k, j) = rw_tend(i, k, j) + rw_tendf(i, k, j)/msfty(i&
 586 &          , j)
 587         IF (rk_step .EQ. 1) THEN
 588           ph_tendfd(i, k, j) = ph_tendfd(i, k, j) + ph_saved(i, k, j)
 589           ph_tendf(i, k, j) = ph_tendf(i, k, j) + ph_save(i, k, j)
 590         END IF
 591 ! divide by my to couple scalar
 592         ph_tendd(i, k, j) = ph_tendd(i, k, j) + ph_tendfd(i, k, j)/msfty&
 593 &          (i, j)
 594         ph_tend(i, k, j) = ph_tend(i, k, j) + ph_tendf(i, k, j)/msfty(i&
 595 &          , j)
 596       END DO
 597     END DO
 598   END DO
 599   IF (jte .GT. jde - 1) THEN
 600     min5 = jde - 1
 601   ELSE
 602     min5 = jte
 603   END IF
 604   DO j=jts,min5
 605     DO k=kts,kte-1
 606       IF (ite .GT. ide - 1) THEN
 607         min6 = ide - 1
 608       ELSE
 609         min6 = ite
 610       END IF
 611       DO i=its,min6
 612         IF (rk_step .EQ. 1) THEN
 613           t_tendfd(i, k, j) = t_tendfd(i, k, j) + t_saved(i, k, j)
 614           t_tendf(i, k, j) = t_tendf(i, k, j) + t_save(i, k, j)
 615         END IF
 616 ! divide by my to couple theta
 617         t_tendd(i, k, j) = t_tendd(i, k, j) + t_tendfd(i, k, j)/msfty(i&
 618 &          , j) + (mutd(i, j)*h_diabatic(i, k, j)+mut(i, j)*h_diabaticd(i&
 619 &          , k, j))/msfty(i, j)
 620         t_tend(i, k, j) = t_tend(i, k, j) + t_tendf(i, k, j)/msfty(i, j)&
 621 &          + mut(i, j)*h_diabatic(i, k, j)/msfty(i, j)
 622       END DO
 623     END DO
 624   END DO
 625   IF (jte .GT. jde - 1) THEN
 626     min7 = jde - 1
 627   ELSE
 628     min7 = jte
 629   END IF
 630 ! divide by my to couple heating
 631   DO j=jts,min7
 632     IF (ite .GT. ide - 1) THEN
 633       min8 = ide - 1
 634     ELSE
 635       min8 = ite
 636     END IF
 637     DO i=its,min8
 638 ! mu tendencies not coupled with 1/msf
 639       mu_tendd(i, j) = mu_tendd(i, j) + mu_tendfd(i, j)
 640       mu_tend(i, j) = mu_tend(i, j) + mu_tendf(i, j)
 641     END DO
 642   END DO
 643 END SUBROUTINE G_RK_ADDTEND_DRY
 644 !-------------------------------------------------------------------------------
 645
 646 ! Revised by Ning Pan, 2010-08-02
 647 ! SUBROUTINE g_rk_scalar_tend(scs,sce,config_flags,rk_step,dt,g_dt,ru,g_ru,rv, &
 648  SUBROUTINE g_rk_scalar_tend(scs,sce,config_flags,tenddec,rk_step,dt,ru,g_ru,rv, &
 649  g_rv,ww,g_ww,mut,g_mut,mub,mu_old,g_mu_old,alt,g_alt,scalar_old, &
 650  g_scalar_old,scalar,g_scalar,scalar_tends,g_scalar_tends,advect_tend, &
 651  g_advect_tend,h_tendency,g_h_tendency,z_tendency,g_z_tendency,RQVFTEN,g_RQVFTEN,base,moist_step,fnm,fnp,msfux,msfuy,msfvx, &
 652  msfvx_inv,msfvy,msftx,msfty,rdx,rdy,rdn,rdnw,khdif,kvdif,xkmhd,g_xkmhd, &
 653 ! Revised by Ning Pan, 2010-08-02
 654 ! diff_6th_opt,diff_6th_factor,g_diff_6th_factor,adv_opt,ids,ide,jds,jde,kds,kde, &
 655  diff_6th_opt,diff_6th_factor,adv_opt,ids,ide,jds,jde,kds,kde, &
 656  ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
 657
 658  IMPLICIT NONE
 659
 660  REAL :: Tmpv1,g_Tmpv1
 661  TYPE(grid_config_rec_type) :: config_flags
 662  LOGICAL :: tenddec
 663  INTEGER :: rk_step,scs,sce
 664  INTEGER :: ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte
 665  LOGICAL :: moist_step
 666  REAL,DIMENSION(ims:ime,kms:kme,jms:jme,scs:sce) :: scalar,g_scalar,scalar_old, &
 667  g_scalar_old
 668  REAL,DIMENSION(ims:ime,kms:kme,jms:jme,scs:sce) :: scalar_tends,g_scalar_tends
 669  REAL,DIMENSION(ims:ime,kms:kme,jms:jme) :: advect_tend,g_advect_tend
 670  REAL,DIMENSION(ims:ime,kms:kme,jms:jme) :: h_tendency, z_tendency
 671  REAL,DIMENSION(ims:ime,kms:kme,jms:jme) :: g_h_tendency, g_z_tendency
 672  REAL,DIMENSION(ims:ime,kms:kme,jms:jme) :: RQVFTEN,g_RQVFTEN
 673  REAL,DIMENSION(ims:ime,kms:kme,jms:jme) :: ru,g_ru,rv,g_rv,ww,g_ww,xkmhd, &
 674  g_xkmhd,alt,g_alt
 675  REAL,DIMENSION(kms:kme) :: fnm,fnp,rdn,rdnw,base
 676  REAL,DIMENSION(ims:ime,jms:jme) :: msfux,msfuy,msfvx,msfvx_inv,msfvy,msftx,msfty,mub, &
 677  mut,g_mut,mu_old,g_mu_old
 678  REAL :: rdx,rdy,khdif,kvdif
 679  INTEGER :: diff_6th_opt
 680 ! Revised by Ning Pan, 2010-08-02
 681 ! REAL :: diff_6th_factor,g_diff_6th_factor
 682 ! REAL :: dt,g_dt
 683  REAL :: diff_6th_factor
 684  REAL :: dt
 685  INTEGER :: adv_opt
 686
 687  INTEGER :: im,i,j,k
 688  INTEGER :: time_step
 689  REAL :: khdq,kvdq,tendency,g_tendency
 690
 691  khdq =khdif/prandtl
 692
 693  kvdq =kvdif/prandtl
 694
 695  DO im =scs,sce
 696
 697  CALL g_zero_tend(advect_tend(ims,kms,jms),g_advect_tend(ims,kms,jms) &
 698 ,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
 699  CALL g_zero_tend(h_tendency(ims,kms,jms),g_h_tendency(ims,kms,jms) &
 700 ,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
 701  CALL g_zero_tend(z_tendency(ims,kms,jms),g_z_tendency(ims,kms,jms) &
 702 ,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
 703 !This line is fail to be recognized
 704       CALL nl_get_time_step ( 1, time_step )
 705
 706  IF( (rk_step == 3) .and. (adv_opt == POSITIVEDEF) ) THEN
 707
 708  CALL g_advect_scalar_pd(scalar(ims,kms,jms,im),g_scalar(ims,kms,jms,im) &
 709  ,scalar_old(ims,kms,jms,im),g_scalar_old(ims,kms,jms,im),advect_tend(ims,kms,jms) &
 710  ,g_advect_tend(ims,kms,jms),h_tendency(ims,kms,jms),g_h_tendency(ims,kms,jms),z_tendency(ims,kms,jms),g_z_tendency(ims,kms,jms) &
 711  ,ru,g_ru,rv,g_rv,ww,g_ww,mut,g_mut,mub, &
 712  mu_old,g_mu_old,time_step,config_flags,tenddec,msfux,msfuy,msfvx,msfvy,msftx,msfty,fnm, &
 713 ! Revised by Ning Pan, 2010-08-02
 714 ! fnp,rdx,rdy,rdnw,dt,g_dt,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite, &
 715  fnp,rdx,rdy,rdnw,dt,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite, &
 716  jts,jte,kts,kte)
 717  ELSE IF( (rk_step == 3) .and. (adv_opt == MONOTONIC) ) THEN
 718
 719  CALL g_advect_scalar_mono(scalar(ims,kms,jms,im),g_scalar(ims,kms,jms,im) &
 720  ,scalar_old(ims,kms,jms,im),g_scalar_old(ims,kms,jms,im),advect_tend(ims,kms,jms) &
 721  ,g_advect_tend(ims,kms,jms),h_tendency(ims,kms,jms),g_h_tendency(ims,kms,jms),z_tendency(ims,kms,jms),g_z_tendency(ims,kms,jms) &
 722  ,ru,g_ru,rv,g_rv,ww,g_ww,mut,g_mut,mub, &
 723  mu_old,g_mu_old,config_flags,tenddec,msfux,msfuy,msfvx,msfvy,msftx,msfty,fnm,fnp,rdx,rdy, &
 724 ! Revised by Ning Pan, 2010-08-02
 725 ! rdnw,dt,g_dt,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
 726  rdnw,dt,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
 727
 728  ELSE IF( (rk_step == 3) .and. (adv_opt == WENO_SCALAR) ) THEN
 729
 730    CALL g_advect_scalar_weno ( scalar(ims,kms,jms,im),     &
 731                             g_scalar(ims,kms,jms,im),      &
 732                             scalar(ims,kms,jms,im),        &
 733                             g_scalar(ims,kms,jms,im),      &
 734                             advect_tend(ims,kms,jms),      &
 735                             g_advect_tend(ims,kms,jms),    &
 736                             ru, g_ru, rv, g_rv, ww, g_ww,  &
 737                             mut, time_step,    &
 738                             config_flags,                  &
 739                             msfux, msfuy, msfvx, msfvy,    &
 740                             msftx, msfty, fnm, fnp,        &
 741                             rdx, rdy, rdnw,                &
 742                             ids, ide, jds, jde, kds, kde,  &
 743                             ims, ime, jms, jme, kms, kme,  &
 744                             its, ite, jts, jte, kts, kte  )
 745
 746  ELSEIF( (rk_step == 3) .and. (adv_opt == WENOPD_SCALAR) ) THEN
 747
 748    CALL g_advect_scalar_wenopd   ( scalar(ims,kms,jms,im),             &
 749                                  g_scalar(ims,kms,jms,im),           &
 750                                  scalar_old(ims,kms,jms,im),         &
 751                                  g_scalar_old(ims,kms,jms,im),       &
 752                                  advect_tend(ims,kms,jms),           &
 753                                  g_advect_tend(ims,kms,jms),         &
 754                                  ru, g_ru, rv, g_rv, ww, g_ww,       &
 755                                  mut, g_mut, mub, mu_old, g_mu_old,  &
 756                                  time_step, config_flags,            &
 757                                  msfux, msfuy, msfvx, msfvy,         &
 758                                  msftx, msfty, fnm, fnp,             &
 759                                  rdx, rdy, rdnw,dt,                  &
 760                                  ids, ide, jds, jde, kds, kde,       &
 761                                  ims, ime, jms, jme, kms, kme,       &
 762                                  its, ite, jts, jte, kts, kte     )
 763
 764  ELSE
 765
 766  CALL g_advect_scalar(scalar(ims,kms,jms,im),g_scalar(ims,kms,jms,im) &
 767  ,scalar(ims,kms,jms,im),g_scalar(ims,kms,jms,im),advect_tend(ims,kms,jms) &
 768  ,g_advect_tend(ims,kms,jms),ru,g_ru,rv,g_rv,ww,g_ww,mut, &
 769  time_step,config_flags,msfux,msfuy,msfvx,msfvy,msftx,msfty,fnm,fnp,rdx,rdy,rdnw,ids, &
 770  ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
 771  END IF
 772
 773  IF((config_flags%cu_physics == GDSCHEME .OR. config_flags%cu_physics == G3SCHEME .OR. &
 774      config_flags%cu_physics == KFETASCHEME .OR. &      ! new trigger in KF
 775      config_flags%cu_physics == MSKFSCHEME  .OR. &
 776      config_flags%cu_physics == TIEDTKESCHEME .OR. &    ! Tiedtke
 777      config_flags%cu_physics == NTIEDTKESCHEME)    &    ! new Tiedtke
 778      .and. moist_step .and. ( im == P_QV) ) THEN
 779
 780  CALL g_set_tend(RQVFTEN,g_RQVFTEN,advect_tend,g_advect_tend,msfty,ids,ide, &
 781  jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
 782  ENDIF
 783
 784  IF( rk_step == 1 ) THEN
 785
 786  IF(config_flags%diff_opt .eq. 1) THEN
 787
 788  CALL g_horizontal_diffusion('m',scalar(ims,kms,jms,im),g_scalar(ims,kms,jms,im) &
 789  ,scalar_tends(ims,kms,jms,im),g_scalar_tends(ims,kms,jms,im),mut,g_mut, &
 790  config_flags,msfux,msfuy,msfvx,msfvx_inv,msfvy,msftx,msfty,khdq,xkmhd,g_xkmhd,rdx, &
 791  rdy,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
 792
 793  IF(config_flags%bl_pbl_physics .eq. 0) THEN
 794
 795  IF( (moist_step) .and. ( im == P_QV)) THEN
 796
 797  CALL g_vertical_diffusion_mp(scalar(ims,kms,jms,im),g_scalar(ims,kms,jms,im) &
 798  ,scalar_tends(ims,kms,jms,im),g_scalar_tends(ims,kms,jms,im),config_flags,base, &
 799  alt,g_alt,mut,g_mut,rdn,rdnw,kvdq,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms, &
 800  kme,its,ite,jts,jte,kts,kte)
 801  ELSE
 802
 803  CALL g_vertical_diffusion('m',scalar(ims,kms,jms,im),g_scalar(ims,kms,jms,im) &
 804  ,scalar_tends(ims,kms,jms,im),g_scalar_tends(ims,kms,jms,im),config_flags,alt, &
 805  g_alt,mut,g_mut,rdn,rdnw,kvdq,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme, &
 806  its,ite,jts,jte,kts,kte)
 807  END IF
 808  ENDIF
 809  ENDIF
 810
 811  IF( diff_6th_opt .NE. 0 ) CALL g_sixth_order_diffusion('m',scalar(ims,kms,jms,im) &
 812  ,g_scalar(ims,kms,jms,im),scalar_tends(ims,kms,jms,im),g_scalar_tends(ims,kms, &
 813 ! Revised by Ning Pan, 2010-08-02
 814 ! jms,im),mut,g_mut,dt,g_dt,config_flags,diff_6th_opt,diff_6th_factor, &
 815 ! g_diff_6th_factor,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
 816  jms,im),mut,g_mut,dt,config_flags,diff_6th_opt,diff_6th_factor, &
 817  ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
 818   ENDIF
 819  ENDDO
 820
 821  END SUBROUTINE g_rk_scalar_tend
 822
 823 !-------------------------------------------------------------------------------
 824 !        Generated by TAPENADE     (INRIA, Tropics team)
 825 !  Tapenade 3.10 (r5498) - 20 Jan 2015 09:48
 826 !
 827 !  Differentiation of q_diabatic_add in forward (tangent) mode:
 828 !   variations   of useful results: scalar_tends
 829 !   with respect to varying inputs: qc_diabatic qv_diabatic scalar_tends
 830 !                mu
 831 !   RW status of diff variables: qc_diabatic:in qv_diabatic:in
 832 !                scalar_tends:in-out mu:in
 833 SUBROUTINE g_Q_DIABATIC_ADD(scs, sce, dt, mu, mud, qv_diabatic, &
 834 & qv_diabaticd, qc_diabatic, qc_diabaticd, scalar_tends, scalar_tendsd, &
 835 & ids, ide, jds, jde, kds, kde, ims, ime, jms, jme, kms, kme, its, ite, &
 836 & jts, jte, kts, kte)
 837   IMPLICIT NONE
 838 !  Input data.
 839   INTEGER, INTENT(IN) :: scs, sce
 840   INTEGER, INTENT(IN) :: ids, ide, jds, jde, kds, kde, ims, ime, jms, &
 841 & jme, kms, kme, its, ite, jts, jte, kts, kte
 842   REAL, DIMENSION(ims:ime, jms:jme), INTENT(IN) :: mu
 843   REAL, DIMENSION(ims:ime, jms:jme), INTENT(IN) :: mud
 844   REAL, DIMENSION(ims:ime, kms:kme, jms:jme), INTENT(IN) :: qv_diabatic&
 845 & , qc_diabatic
 846   REAL, DIMENSION(ims:ime, kms:kme, jms:jme), INTENT(IN) :: qv_diabaticd&
 847 & , qc_diabaticd
 848   REAL, DIMENSION(ims:ime, kms:kme, jms:jme, scs:sce), INTENT(INOUT) :: &
 849 & scalar_tends
 850   REAL, DIMENSION(ims:ime, kms:kme, jms:jme, scs:sce), INTENT(INOUT) :: &
 851 & scalar_tendsd
 852   REAL, INTENT(IN) :: dt
 853 ! Local data
 854   INTEGER :: im, i, j, k
 855   INTEGER :: min4
 856   INTEGER :: min3
 857   INTEGER :: min2
 858   INTEGER :: min1
 859 scalar_loop:DO im=scs,sce
 860     IF (im .EQ. p_qv) THEN
 861       IF (jte .GT. jde - 1) THEN
 862         min1 = jde - 1
 863       ELSE
 864         min1 = jte
 865       END IF
 866       DO j=jts,min1
 867         DO k=kts,kte-1
 868           IF (ite .GT. ide - 1) THEN
 869             min2 = ide - 1
 870           ELSE
 871             min2 = ite
 872           END IF
 873           DO i=its,min2
 874             scalar_tendsd(i,k,j,im) = scalar_tendsd(i,k,j,im) +     &
 875                                       qv_diabaticd(i,k,j)*mu(i,j) + &
 876                                       qv_diabatic(i,k,j)*mud(i,j)
 877             scalar_tends(i,k,j,im) = scalar_tends(i,k,j,im) +  &
 878                                      qv_diabatic(i,k,j)*mu(i,j)
 879           END DO
 880         END DO
 881       END DO
 882     END IF
 883     IF (im .EQ. p_qc) THEN
 884       IF (jte .GT. jde - 1) THEN
 885         min3 = jde - 1
 886       ELSE
 887         min3 = jte
 888       END IF
 889       DO j=jts,min3
 890         DO k=kts,kte-1
 891           IF (ite .GT. ide - 1) THEN
 892             min4 = ide - 1
 893           ELSE
 894             min4 = ite
 895           END IF
 896           DO i=its,min4
 897             scalar_tendsd(i,k,j,im) = scalar_tendsd(i,k,j,im) +     &
 898                                       qc_diabaticd(i,k,j)*mu(i,j) + &
 899                                       qc_diabatic(i,k,j)*mud(i,j)
 900             scalar_tends(i,k,j,im) = scalar_tends(i,k,j,im) + &
 901                                      qc_diabatic(i,k,j)*mu(i,j)
 902           END DO
 903         END DO
 904       END DO
 905     END IF
 906   END DO scalar_loop
 907 END SUBROUTINE g_Q_DIABATIC_ADD
 908
 909 !-------------------------------------------------------------------------------
 910 !        Generated by TAPENADE     (INRIA, Tropics team)
 911 !  Tapenade 3.10 (r5498) - 20 Jan 2015 09:48
 912 !
 913 !  Differentiation of q_diabatic_subtr in forward (tangent) mode:
 914 !   variations   of useful results: scalar
 915 !   with respect to varying inputs: qc_diabatic qv_diabatic scalar
 916 !   RW status of diff variables: qc_diabatic:in qv_diabatic:in
 917 !                scalar:in-out
 918 SUBROUTINE g_Q_DIABATIC_SUBTR(scs, sce, dt, qv_diabatic, qv_diabaticd, &
 919 & qc_diabatic, qc_diabaticd, scalar, scalard, ids, ide, jds, jde, kds, &
 920 & kde, ims, ime, jms, jme, kms, kme, its, ite, jts, jte, kts, kte)
 921   IMPLICIT NONE
 922 !  Input data.
 923   INTEGER, INTENT(IN) :: scs, sce
 924   INTEGER, INTENT(IN) :: ids, ide, jds, jde, kds, kde, ims, ime, jms, &
 925 & jme, kms, kme, its, ite, jts, jte, kts, kte
 926   REAL, DIMENSION(ims:ime, kms:kme, jms:jme), INTENT(IN) :: qv_diabatic&
 927 & , qc_diabatic
 928   REAL, DIMENSION(ims:ime, kms:kme, jms:jme), INTENT(IN) :: qv_diabaticd&
 929 & , qc_diabaticd
 930   REAL, DIMENSION(ims:ime, kms:kme, jms:jme, scs:sce), INTENT(INOUT) :: &
 931 & scalar
 932   REAL, DIMENSION(ims:ime, kms:kme, jms:jme, scs:sce), INTENT(INOUT) :: &
 933 & scalard
 934   REAL, INTENT(IN) :: dt
 935 ! Local data
 936   INTEGER :: im, i, j, k
 937   INTEGER :: min4
 938   INTEGER :: min3
 939   INTEGER :: min2
 940   INTEGER :: min1
 941 scalar_loop:DO im=scs,sce
 942     IF (im .EQ. p_qv) THEN
 943       IF (jte .GT. jde - 1) THEN
 944         min1 = jde - 1
 945       ELSE
 946         min1 = jte
 947       END IF
 948       DO j=jts,min1
 949         DO k=kts,kte-1
 950           IF (ite .GT. ide - 1) THEN
 951             min2 = ide - 1
 952           ELSE
 953             min2 = ite
 954           END IF
 955           DO i=its,min2
 956             scalard(i,k,j,im) = scalard(i,k,j,im) - dt*qv_diabaticd(i,k,j)
 957             scalar(i,k,j,im) = scalar(i,k,j,im) - dt*qv_diabatic(i,k,j)
 958           END DO
 959         END DO
 960       END DO
 961     END IF
 962     IF (im .EQ. p_qc) THEN
 963       IF (jte .GT. jde - 1) THEN
 964         min3 = jde - 1
 965       ELSE
 966         min3 = jte
 967       END IF
 968       DO j=jts,min3
 969         DO k=kts,kte-1
 970           IF (ite .GT. ide - 1) THEN
 971             min4 = ide - 1
 972           ELSE
 973             min4 = ite
 974           END IF
 975           DO i=its,min4
 976             scalard(i,k,j,im) = scalard(i,k,j,im) - dt*qc_diabaticd(i,k,j)
 977             scalar(i,k,j,im) = scalar(i,k,j,im) - dt*qc_diabatic(i,k,j)
 978           END DO
 979         END DO
 980       END DO
 981     END IF
 982   END DO scalar_loop
 983 END SUBROUTINE g_Q_DIABATIC_SUBTR
 984
 985 !-------------------------------------------------------------------------------
 986
 987 SUBROUTINE g_rk_update_scalar ( scs, sce,                      &
 988                                 scalar_1, g_scalar_1, scalar_2, g_scalar_2, sc_tend, g_sc_tend,  &
 989                                 advh_t, g_advh_t, advz_t, g_advz_t,              &
 990                                 advect_tend, g_advect_tend,    &
 991                                 h_tendency, g_h_tendency,      &
 992                                 z_tendency, g_z_tendency,      &
 993                                 msftx, msfty,     &
 994                                 mu_old, g_mu_old, mu_new, g_mu_new, mu_base,  &
 995                                 rk_step, dt, spec_zone,        &
 996                                 config_flags,                  &
 997                                 tenddec,                       &
 998                                 ids, ide, jds, jde, kds, kde,  &
 999                                 ims, ime, jms, jme, kms, kme,  &
1000                                 its, ite, jts, jte, kts, kte  )
1001
1002    IMPLICIT NONE
1003
1004    !  Input data.
1005
1006    TYPE(grid_config_rec_type), INTENT(IN) :: config_flags
1007    LOGICAL :: tenddec
1008
1009    INTEGER, INTENT(IN) :: scs, sce, rk_step, spec_zone
1010    INTEGER, INTENT(IN) :: ids, ide, jds, jde, kds, kde, &
1011                           ims, ime, jms, jme, kms, kme, &
1012                           its, ite, jts, jte, kts, kte
1013
1014    REAL,    INTENT(IN) :: dt
1015
1016    REAL, DIMENSION(ims:ime, kms:kme, jms:jme , scs:sce),                &
1017          INTENT(INOUT)                                  :: g_scalar_1,  &
1018                                                            g_scalar_2
1019    REAL, DIMENSION(ims:ime, kms:kme, jms:jme , scs:sce),                &
1020          INTENT(INOUT)                                  :: scalar_1,    &
1021                                                            scalar_2
1022
1023    REAL, DIMENSION(ims:ime, kms:kme, jms:jme , scs:sce),                &
1024          INTENT(IN)                                     :: g_sc_tend
1025    REAL, DIMENSION(ims:ime, kms:kme, jms:jme , scs:sce),                &
1026          INTENT(IN)                                     :: sc_tend
1027
1028    REAL, DIMENSION(ims:ime, kms:kme, jms:jme ),                &
1029          INTENT(IN)                                  :: g_advect_tend
1030    REAL, DIMENSION(ims:ime, kms:kme, jms:jme ),                &
1031          INTENT(IN)                                  :: advect_tend
1032
1033    REAL, DIMENSION(ims:ime, kms:kme, jms:jme ), OPTIONAL :: advh_t,  advz_t ! accumulating for output
1034    REAL, DIMENSION(ims:ime, kms:kme, jms:jme ), OPTIONAL :: g_advh_t,  g_advz_t ! accumulating for output
1035    REAL, DIMENSION(ims:ime, kms:kme, jms:jme ) :: h_tendency, z_tendency ! from rk_scalar_tend
1036    REAL, DIMENSION(ims:ime, kms:kme, jms:jme ) :: g_h_tendency, g_z_tendency ! from rk_scalar_tend
1037
1038    REAL, DIMENSION(ims:ime, jms:jme  ), INTENT(IN   ) ::  g_mu_old,  &
1039                                                           g_mu_new
1040    REAL, DIMENSION(ims:ime, jms:jme  ), INTENT(IN   ) ::  mu_old,  &
1041                                                           mu_new,  &
1042                                                           mu_base, &
1043                                                           msftx,   &
1044                                                           msfty
1045
1046    INTEGER :: i,j,k,im
1047    REAL    :: sc_middle, msfsq
1048    REAL, DIMENSION(its:ite) :: g_muold, g_r_munew
1049    REAL, DIMENSION(its:ite) :: muold, r_munew
1050
1051    REAL, DIMENSION(its:ite, kts:kte, jts:jte) :: g_tendency
1052    REAL, DIMENSION(its:ite, kts:kte, jts:jte) :: tendency
1053
1054    INTEGER :: i_start,i_end,j_start,j_end,k_start,k_end
1055    INTEGER :: i_start_spc,i_end_spc,j_start_spc,j_end_spc,k_start_spc,k_end_spc
1056
1057 !<DESCRIPTION>
1058 !
1059 !  rk_scalar_update advances the scalar equation given the time t value
1060 !  of the scalar and the scalar tendency.
1061 !
1062 !</DESCRIPTION>
1063
1064
1065 !
1066 !  set loop limits.
1067
1068       i_start = its
1069       i_end   = min(ite,ide-1)
1070       j_start = jts
1071       j_end   = min(jte,jde-1)
1072       k_start = kts
1073       k_end   = kte-1
1074
1075       i_start_spc = i_start
1076       i_end_spc   = i_end
1077       j_start_spc = j_start
1078       j_end_spc   = j_end
1079       k_start_spc = k_start
1080       k_end_spc   = k_end
1081
1082     IF( config_flags%nested .or. config_flags%specified ) THEN
1083       IF( .NOT. config_flags%periodic_x)i_start = max( its,ids+spec_zone )
1084       IF( .NOT. config_flags%periodic_x)i_end   = min( ite,ide-spec_zone-1 )
1085       j_start = max( jts,jds+spec_zone )
1086       j_end   = min( jte,jde-spec_zone-1 )
1087       k_start = kts
1088       k_end   = min( kte, kde-1 )
1089     ENDIF
1090
1091     IF ( rk_step == 1 ) THEN
1092
1093       !  replace t-dt values (in scalar_1) with t values scalar_2,
1094       !  then compute new values by adding tendency to values at t
1095
1096       DO  im = scs,sce
1097
1098        DO  j = jts, min(jte,jde-1)
1099        DO  k = kts, min(kte,kde-1)
1100        DO  i = its, min(ite,ide-1)
1101            g_tendency(i,k,j) = 0.
1102            tendency(i,k,j) = 0.
1103        ENDDO
1104        ENDDO
1105        ENDDO
1106
1107        DO  j = j_start,j_end
1108        DO  k = k_start,k_end
1109        DO  i = i_start,i_end
1110           ! scalar was coupled with my
1111            g_tendency(i,k,j) = g_advect_tend(i,k,j) * msfty(i,j)
1112            tendency(i,k,j) = advect_tend(i,k,j) * msfty(i,j)
1113        ENDDO
1114        ENDDO
1115        ENDDO
1116
1117        DO  j = j_start_spc,j_end_spc
1118        DO  k = k_start_spc,k_end_spc
1119        DO  i = i_start_spc,i_end_spc
1120            g_tendency(i,k,j) = g_tendency(i,k,j) + g_sc_tend(i,k,j,im)
1121            tendency(i,k,j) = tendency(i,k,j) + sc_tend(i,k,j,im)
1122        ENDDO
1123        ENDDO
1124        ENDDO
1125
1126       DO  j = jts, min(jte,jde-1)
1127
1128       DO  i = its, min(ite,ide-1)
1129         g_muold(i) = g_mu_old(i,j)
1130         muold(i) = mu_old(i,j) + mu_base(i,j)
1131         g_r_munew(i) = -g_mu_new(i,j) / ((mu_new(i,j)+mu_base(i,j))*(mu_new(i,j)+mu_base(i,j)))
1132         r_munew(i) = 1./(mu_new(i,j) + mu_base(i,j))
1133       ENDDO
1134
1135       DO  k = kts, min(kte,kde-1)
1136       DO  i = its, min(ite,ide-1)
1137
1138         g_scalar_1(i,k,j,im) = g_scalar_2(i,k,j,im)
1139         scalar_1(i,k,j,im) = scalar_2(i,k,j,im)
1140         g_scalar_2(i,k,j,im) = (muold(i)*g_scalar_1(i,k,j,im)+g_muold(i)*scalar_1(i,k,j,im)   &
1141                              + dt*g_tendency(i,k,j))*r_munew(i) &
1142                              + (muold(i)*scalar_1(i,k,j,im)   &
1143                              + dt*tendency(i,k,j))*g_r_munew(i)
1144         scalar_2(i,k,j,im) = (muold(i)*scalar_1(i,k,j,im)   &
1145                              + dt*tendency(i,k,j))*r_munew(i)
1146
1147       ENDDO
1148       ENDDO
1149       ENDDO
1150
1151       ENDDO
1152
1153     ELSE
1154
1155       !  just compute new values, scalar_1 already at time t.
1156
1157       DO  im = scs, sce
1158
1159        DO  j = jts, min(jte,jde-1)
1160        DO  k = kts, min(kte,kde-1)
1161        DO  i = its, min(ite,ide-1)
1162            g_tendency(i,k,j) = 0.
1163            tendency(i,k,j) = 0.
1164        ENDDO
1165        ENDDO
1166        ENDDO
1167
1168        DO  j = j_start,j_end
1169        DO  k = k_start,k_end
1170        DO  i = i_start,i_end
1171            ! scalar was coupled with my
1172            g_tendency(i,k,j) = g_advect_tend(i,k,j) * msfty(i,j)
1173            tendency(i,k,j) = advect_tend(i,k,j) * msfty(i,j)
1174        ENDDO
1175        ENDDO
1176        ENDDO
1177
1178        DO  j = j_start_spc,j_end_spc
1179        DO  k = k_start_spc,k_end_spc
1180        DO  i = i_start_spc,i_end_spc
1181            g_tendency(i,k,j) = g_tendency(i,k,j) + g_sc_tend(i,k,j,im)
1182            tendency(i,k,j) = tendency(i,k,j) + sc_tend(i,k,j,im)
1183        ENDDO
1184        ENDDO
1185        ENDDO
1186
1187       DO  j = jts, min(jte,jde-1)
1188
1189       DO  i = its, min(ite,ide-1)
1190         g_muold(i) = g_mu_old(i,j)
1191         muold(i) = mu_old(i,j) + mu_base(i,j)
1192         g_r_munew(i) = -g_mu_new(i,j) / ((mu_new(i,j)+mu_base(i,j))*(mu_new(i,j)+mu_base(i,j)))
1193         r_munew(i) = 1./(mu_new(i,j) + mu_base(i,j))
1194       ENDDO
1195
1196       DO  k = kts, min(kte,kde-1)
1197       DO  i = its, min(ite,ide-1)
1198
1199         g_scalar_2(i,k,j,im) = (muold(i)*g_scalar_1(i,k,j,im)+g_muold(i)*scalar_1(i,k,j,im)   &
1200                              + dt*g_tendency(i,k,j))*r_munew(i)  &
1201                              + (muold(i)*scalar_1(i,k,j,im)      &
1202                              + dt*tendency(i,k,j))*g_r_munew(i)
1203         scalar_2(i,k,j,im) = (muold(i)*scalar_1(i,k,j,im)   &
1204                              + dt*tendency(i,k,j))*r_munew(i)
1205
1206       ENDDO
1207       ENDDO
1208
1209       ! This is separated from the k/i-loop above for better performance
1210       IF ( PRESENT(advh_t) .AND. PRESENT(advz_t) .AND. PRESENT(g_advh_t) .AND. PRESENT(g_advz_t) ) THEN
1211          IF (tenddec.and.rk_step.eq.config_flags%rk_ord) THEN
1212             DO k = kts, min(kte,kde-1)
1213             DO i = its, min(ite,ide-1)
1214
1215                g_advh_t(i,k,j) = g_advh_t(i,k,j) + (dt*g_h_tendency(i,k,j)* msfty(i,j))*r_munew(i) &
1216                                             + (dt*h_tendency(i,k,j)* msfty(i,j))*g_r_munew(i)
1217                advh_t(i,k,j) = advh_t(i,k,j) + (dt*h_tendency(i,k,j)* msfty(i,j))*r_munew(i)
1218                g_advz_t(i,k,j) = g_advz_t(i,k,j) + (dt*g_z_tendency(i,k,j)* msfty(i,j))*r_munew(i) &
1219                                             + (dt*z_tendency(i,k,j)* msfty(i,j))*g_r_munew(i)
1220                advz_t(i,k,j) = advz_t(i,k,j) + (dt*z_tendency(i,k,j)* msfty(i,j))*r_munew(i)
1221
1222             ENDDO
1223             ENDDO
1224          END IF
1225       END IF
1226
1227       ENDDO
1228
1229       ENDDO
1230
1231     END IF
1232
1233 END SUBROUTINE g_rk_update_scalar
1234
1235 !-------------------------------------------------------------------------------
1236
1237  SUBROUTINE g_rk_update_scalar_pd(scs,sce,scalar,g_scalar,sc_tend,g_sc_tend, &
1238 ! Revised by Ning Pan, 2010-08-03
1239 ! mu_old,g_mu_old,mu_new,g_mu_new,mu_base,g_mu_base,rk_step,dt,spec_zone, &
1240  mu_old,g_mu_old,mu_new,g_mu_new,mu_base,rk_step,dt,spec_zone, &
1241  config_flags,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
1242
1243  IMPLICIT NONE
1244
1245  REAL :: Tmpv1,g_Tmpv1,Tmpv2,g_Tmpv2,Tmpv3,g_Tmpv3
1246  TYPE(grid_config_rec_type) :: config_flags
1247  INTEGER :: scs,sce,rk_step,spec_zone
1248  INTEGER :: ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte
1249  REAL :: dt
1250  REAL,DIMENSION(ims:ime,kms:kme,jms:jme,scs:sce) :: scalar,g_scalar,sc_tend,g_sc_tend
1251 ! Revised by Ning Pan, 2010-08-03
1252 ! REAL,DIMENSION(ims:ime,jms:jme) :: mu_old,g_mu_old,mu_new,g_mu_new,mu_base,g_mu_base
1253  REAL,DIMENSION(ims:ime,jms:jme) :: mu_old,g_mu_old,mu_new,g_mu_new,mu_base
1254  INTEGER :: i,j,k,im
1255  REAL :: sc_middle,g_sc_middle,msfsq,g_msfsq
1256  REAL,DIMENSION(its:ite) :: muold,g_muold,r_munew,g_r_munew
1257  REAL,DIMENSION(its:ite,kts:kte,jts:jte) :: tendency,g_tendency
1258  INTEGER :: i_start,i_end,j_start,j_end,k_start,k_end
1259  INTEGER :: i_start_spc,i_end_spc,j_start_spc,j_end_spc,k_start_spc,k_end_spc
1260
1261  i_start =its
1262
1263  i_end =min(ite,ide-1)
1264
1265  j_start =jts
1266
1267  j_end =min(jte,jde-1)
1268
1269  k_start =kts
1270
1271  k_end =kte-1
1272
1273  i_start_spc =i_start
1274
1275  i_end_spc =i_end
1276
1277  j_start_spc =j_start
1278
1279  j_end_spc =j_end
1280
1281  k_start_spc =k_start
1282
1283  k_end_spc =k_end
1284
1285  IF( config_flags%nested .or. config_flags%specified ) THEN
1286
1287  IF( .NOT. config_flags%periodic_x) i_start =max(its,ids+spec_zone)
1288
1289  IF( .NOT. config_flags%periodic_x) i_end =min(ite,ide-spec_zone-1)
1290
1291  j_start =max(jts,jds+spec_zone)
1292
1293  j_end =min(jte,jde-spec_zone-1)
1294
1295  k_start =kts
1296
1297  k_end =min(kte,kde-1)
1298  ENDIF
1299
1300  DO im =scs,sce
1301  DO j =jts,min(jte,jde-1)
1302  DO k =kts,min(kte,kde-1)
1303  DO i =its,min(ite,ide-1)
1304
1305  g_tendency(i,k,j) =0.0
1306  tendency(i,k,j) =0.
1307
1308  ENDDO
1309  ENDDO
1310  ENDDO
1311
1312  DO j =j_start_spc,j_end_spc
1313  DO k =k_start_spc,k_end_spc
1314  DO i =i_start_spc,i_end_spc
1315
1316  g_tendency(i,k,j) =g_tendency(i,k,j) +g_sc_tend(i,k,j,im)
1317  tendency(i,k,j) =tendency(i,k,j) +sc_tend(i,k,j,im)
1318
1319  g_sc_tend(i,k,j,im) =0.0
1320  sc_tend(i,k,j,im) =0.
1321
1322  ENDDO
1323  ENDDO
1324  ENDDO
1325
1326  DO j =jts,min(jte,jde-1)
1327  DO i =its,min(ite,ide-1)
1328
1329 ! Revised by Ning Pan, 2010-08-03
1330 ! g_muold(i) =g_mu_old(i,j) +g_mu_base(i,j)
1331  g_muold(i) =g_mu_old(i,j)
1332  muold(i) =mu_old(i,j) +mu_base(i,j)
1333
1334 ! Revised by Ning Pan, 2010-08-03
1335 ! g_r_munew(i) =-1.*(g_mu_new(i,j) +g_mu_base(i,j))/((mu_new(i,j) &
1336 ! +mu_base(i,j))*(mu_new(i,j) +mu_base(i,j)))
1337  g_r_munew(i) =-1.*(g_mu_new(i,j))/((mu_new(i,j) &
1338  +mu_base(i,j))*(mu_new(i,j) +mu_base(i,j)))
1339  r_munew(i) =1./(mu_new(i,j) +mu_base(i,j))
1340
1341  ENDDO
1342
1343  DO k =kts,min(kte,kde-1)
1344  DO i =its,min(ite,ide-1)
1345
1346  g_Tmpv1 =muold(i)*g_scalar(i,k,j,im) +g_muold(i)*scalar(i,k,j,im)
1347  Tmpv1 =muold(i)*scalar(i,k,j,im)
1348
1349  g_Tmpv2 =(Tmpv1 +dt*tendency(i,k,j))*g_r_munew(i) +(g_Tmpv1 +dt* &
1350  g_tendency(i,k,j))*r_munew(i)
1351  Tmpv2 =(Tmpv1 +dt*tendency(i,k,j))*r_munew(i)
1352
1353  g_scalar(i,k,j,im) =g_Tmpv2
1354  scalar(i,k,j,im) =Tmpv2
1355
1356  ENDDO
1357  ENDDO
1358  ENDDO
1359  ENDDO
1360
1361  END SUBROUTINE g_rk_update_scalar_pd
1362
1363 !------------------------------------------------------------
1364
1365 SUBROUTINE g_init_zero_tendency(ru_tendf,  g_ru_tendf, &
1366                                 rv_tendf,  g_rv_tendf, &
1367                                 rw_tendf,  g_rw_tendf, &
1368                                 ph_tendf,  g_ph_tendf, &
1369                                 t_tendf,   g_t_tendf,  &
1370                                 tke_tendf, g_tke_tendf, &
1371                                 mu_tendf,  g_mu_tendf,  &
1372                                 moist_tendf, g_moist_tendf,  &
1373 !  NPan - 05/26/10 {
1374 !  Uncomment the corresponding args when chem is needed.
1375 !                                chem_tendf,  g_chem_tendf,   &
1376                                 scalar_tendf,g_scalar_tendf, &
1377                                  tracer_tendf,g_tracer_tendf, &
1378 !  NPan }
1379                                 n_tracer,                   &
1380                                 n_moist,n_chem,n_scalar,rk_step,         &
1381                                 ids, ide, jds, jde, kds, kde,            &
1382                                 ims, ime, jms, jme, kms, kme,            &
1383                                 its, ite, jts, jte, kts, kte             )
1384 !-----------------------------------------------------------------------
1385    IMPLICIT NONE
1386 !-----------------------------------------------------------------------
1387
1388    INTEGER ,       INTENT(IN   ) :: ids, ide, jds, jde, kds, kde, &
1389                                     ims, ime, jms, jme, kms, kme, &
1390                                     its, ite, jts, jte, kts, kte
1391
1392    INTEGER ,       INTENT(IN   ) :: n_moist,n_chem,n_scalar,n_tracer,rk_step
1393
1394    REAL , DIMENSION( ims:ime , kms:kme, jms:jme  ) , INTENT(INOUT) ::  &
1395                                                              g_ru_tendf, &
1396                                                              g_rv_tendf, &
1397                                                              g_rw_tendf, &
1398                                                              g_ph_tendf, &
1399                                                               g_t_tendf, &
1400                                                             g_tke_tendf
1401    REAL , DIMENSION( ims:ime , kms:kme, jms:jme  ) , INTENT(INOUT) ::  &
1402                                                              ru_tendf, &
1403                                                              rv_tendf, &
1404                                                              rw_tendf, &
1405                                                              ph_tendf, &
1406                                                               t_tendf, &
1407                                                             tke_tendf
1408
1409    REAL , DIMENSION( ims:ime , jms:jme  ) , INTENT(INOUT) ::  g_mu_tendf
1410    REAL , DIMENSION( ims:ime , jms:jme  ) , INTENT(INOUT) ::  mu_tendf
1411
1412    REAL , DIMENSION(ims:ime, kms:kme, jms:jme, n_moist),INTENT(INOUT)::&
1413                                                           g_moist_tendf
1414    REAL , DIMENSION(ims:ime, kms:kme, jms:jme, n_moist),INTENT(INOUT)::&
1415                                                           moist_tendf
1416
1417 !  NPan - 05/26/10 {
1418 !  Uncomment the corresponding definations when chem is needed.
1419 !   REAL , DIMENSION(ims:ime, kms:kme, jms:jme, n_chem ),INTENT(INOUT)::&
1420 !                                                          g_chem_tendf
1421 !   REAL , DIMENSION(ims:ime, kms:kme, jms:jme, n_chem ),INTENT(INOUT)::&
1422 !                                                          chem_tendf
1423     REAL , DIMENSION(ims:ime, kms:kme, jms:jme, n_tracer ),INTENT(INOUT)::&
1424                                                            g_tracer_tendf
1425     REAL , DIMENSION(ims:ime, kms:kme, jms:jme, n_tracer ),INTENT(INOUT)::&
1426                                                            tracer_tendf
1427
1428    REAL , DIMENSION(ims:ime, kms:kme, jms:jme, n_scalar ),INTENT(INOUT)::&
1429                                                           g_scalar_tendf
1430    REAL , DIMENSION(ims:ime, kms:kme, jms:jme, n_scalar ),INTENT(INOUT)::&
1431                                                           scalar_tendf
1432 !  NPan }
1433
1434 ! LOCAL VARS
1435
1436    INTEGER :: im, ic, is
1437
1438 !<DESCRIPTION>
1439 !
1440 ! init_zero_tendency
1441 ! sets tendency arrays to zero for all prognostic variables.
1442 !
1443 !</DESCRIPTION>
1444
1445
1446    CALL g_zero_tend ( ru_tendf, g_ru_tendf,            &
1447                       ids, ide, jds, jde, kds, kde,    &
1448                       ims, ime, jms, jme, kms, kme,    &
1449                       its, ite, jts, jte, kts, kte     )
1450
1451    CALL g_zero_tend ( rv_tendf, g_rv_tendf,            &
1452                       ids, ide, jds, jde, kds, kde,    &
1453                       ims, ime, jms, jme, kms, kme,    &
1454                       its, ite, jts, jte, kts, kte     )
1455
1456    CALL g_zero_tend ( rw_tendf, g_rw_tendf,            &
1457                       ids, ide, jds, jde, kds, kde,    &
1458                       ims, ime, jms, jme, kms, kme,    &
1459                       its, ite, jts, jte, kts, kte     )
1460
1461    CALL g_zero_tend ( ph_tendf, g_ph_tendf,            &
1462                       ids, ide, jds, jde, kds, kde,    &
1463                       ims, ime, jms, jme, kms, kme,    &
1464                       its, ite, jts, jte, kts, kte     )
1465
1466    CALL g_zero_tend ( t_tendf, g_t_tendf,            &
1467                       ids, ide, jds, jde, kds, kde,    &
1468                       ims, ime, jms, jme, kms, kme,    &
1469                       its, ite, jts, jte, kts, kte     )
1470
1471    CALL g_zero_tend ( tke_tendf, g_tke_tendf,            &
1472                       ids, ide, jds, jde, kds, kde,    &
1473                       ims, ime, jms, jme, kms, kme,    &
1474                       its, ite, jts, jte, kts, kte     )
1475
1476    CALL g_zero_tend2d ( mu_tendf, g_mu_tendf,            &
1477                       ids, ide, jds, jde, kds, kds,    &
1478                       ims, ime, jms, jme, kms, kms,    &
1479                       its, ite, jts, jte, kts, kts     )
1480
1481 !   DO im=PARAM_FIRST_SCALAR,n_moist
1482    DO im=1,n_moist                      ! make sure first one is zero too
1483       CALL g_zero_tend ( moist_tendf(ims,kms,jms,im), g_moist_tendf(ims,kms,jms,im), &
1484                          ids, ide, jds, jde, kds, kde, &
1485                          ims, ime, jms, jme, kms, kme, &
1486                          its, ite, jts, jte, kts, kte  )
1487    ENDDO
1488
1489 !  NPan - 05/26/10 {
1490 !  Uncomment the corresponding statements when chem is needed.
1491 !!   DO ic=PARAM_FIRST_SCALAR,n_chem
1492 !   DO ic=1,n_chem                       !! make sure first one is zero too
1493 !      CALL g_zero_tend ( chem_tendf(ims,kms,jms,ic), g_chem_tendf(ims,kms,jms,ic), &
1494 !                         ids, ide, jds, jde, kds, kde, &
1495 !                         ims, ime, jms, jme, kms, kme, &
1496 !                         its, ite, jts, jte, kts, kte  )
1497 !   ENDDO
1498
1499 !   DO ic=PARAM_FIRST_SCALAR,n_tracer
1500    DO ic=1,n_tracer                     !! make sure first one is zero too
1501       CALL g_zero_tend ( tracer_tendf(ims,kms,jms,ic), g_tracer_tendf(ims,kms,jms,ic), &
1502                          ids, ide, jds, jde, kds, kde, &
1503                          ims, ime, jms, jme, kms, kme, &
1504                          its, ite, jts, jte, kts, kte  )
1505    ENDDO
1506
1507 !   DO ic=PARAM_FIRST_SCALAR,n_scalar
1508    DO ic=1,n_scalar                       ! make sure first one is zero too
1509       CALL g_zero_tend ( scalar_tendf(ims,kms,jms,ic), g_scalar_tendf(ims,kms,jms,ic), &
1510                          ids, ide, jds, jde, kds, kde, &
1511                          ims, ime, jms, jme, kms, kme, &
1512                          its, ite, jts, jte, kts, kte  )
1513    ENDDO
1514 !  NPan }
1515
1516 END SUBROUTINE g_init_zero_tendency
1517
1518 !        Generated by TAPENADE     (INRIA, Tropics team)
1519 !  Tapenade 3.5 (r3785) - 22 Mar 2011 18:35
1520 !
1521 !  Differentiation of calculate_phy_tend in forward (tangent) mode:
1522 !   variations   of useful results: rthndgdten rublten rqvndgdten
1523 !                rthraten rqccuten rthcuten rqicuten rvndgdten
1524 !                rqscuten rqrshten rqvshten rucuten rvshten rqvblten
1525 !                rvblten rqcshten rthshten rqgshten rqishten rqcblten
1526 !                rthblten rqrcuten rqiblten rqsshten rqvcuten rvcuten
1527 !                rushten rundgdten
1528 !   with respect to varying inputs: rthndgdten rublten rqvndgdten
1529 !                rthraten rqccuten rthcuten rqicuten rvndgdten
1530 !                rqscuten rqrshten rqvshten rucuten rvshten rqvblten
1531 !                rvblten rqcshten rthshten rqgshten rqishten rqcblten
1532 !                rthblten rqrcuten rqiblten rqsshten rqvcuten rvcuten
1533 !                rushten muu muv rundgdten mu
1534 !   RW status of diff variables: rthndgdten:in-out rublten:in-out
1535 !                rqvndgdten:in-out rthraten:in-out rqccuten:in-out
1536 !                rthcuten:in-out rqicuten:in-out rvndgdten:in-out
1537 !                rqscuten:in-out rqrshten:in-out rqvshten:in-out
1538 !                rucuten:in-out rvshten:in-out rqvblten:in-out
1539 !                rvblten:in-out rqcshten:in-out rthshten:in-out
1540 !                rqgshten:in-out rqishten:in-out rqcblten:in-out
1541 !                rthblten:in-out rqrcuten:in-out rqiblten:in-out
1542 !                rqsshten:in-out rqvcuten:in-out rvcuten:in-out
1543 !                rushten:in-out muu:in muv:in rundgdten:in-out
1544 !                mu:in
1545 SUBROUTINE G_CALCULATE_PHY_TEND(config_flags, mu, mud, muu, muud, muv, &
1546 &  muvd, pi3d, rthraten, rthratend, rublten, rubltend, rvblten, rvbltend&
1547 &  , rthblten, rthbltend, rqvblten, rqvbltend, rqcblten, rqcbltend, &
1548 &  rqiblten, rqibltend, rucuten, rucutend, rvcuten, rvcutend, rthcuten, &
1549 &  rthcutend, rqvcuten, rqvcutend, rqccuten, rqccutend, rqrcuten, &
1550 &  rqrcutend, rqicuten, rqicutend, rqscuten, rqscutend, rushten, rushtend&
1551 &  , rvshten, rvshtend, rthshten, rthshtend, rqvshten, rqvshtend, &
1552 &  rqcshten, rqcshtend, rqrshten, rqrshtend, rqishten, rqishtend, &
1553 &  rqsshten, rqsshtend, rqgshten, rqgshtend, rundgdten, rundgdtend, &
1554 &  rvndgdten, rvndgdtend, rthndgdten, rthndgdtend, rqvndgdten, &
1555 &  rqvndgdtend, rmundgdten, ids, ide, jds, jde, kds, kde, ims, ime, jms, &
1556 &  jme, kms, kme, its, ite, jts, jte, kts, kte)
1557   IMPLICIT NONE
1558   TYPE(GRID_CONFIG_REC_TYPE), INTENT(IN) :: config_flags
1559   INTEGER, INTENT(IN) :: ids, ide, jds, jde, kds, kde, ims, ime, jms, &
1560 &  jme, kms, kme, its, ite, jts, jte, kts, kte
1561   REAL, DIMENSION(ims:ime, kms:kme, jms:jme), INTENT(IN) :: pi3d
1562   REAL, DIMENSION(ims:ime, jms:jme), INTENT(IN) :: mu, muu, muv
1563   REAL, DIMENSION(ims:ime, jms:jme), INTENT(IN) :: mud, muud, muvd
1564 ! radiation
1565   REAL, DIMENSION(ims:ime, kms:kme, jms:jme), INTENT(INOUT) :: rthraten
1566   REAL, DIMENSION(ims:ime, kms:kme, jms:jme), INTENT(INOUT) :: rthratend
1567 ! cumulus
1568   REAL, DIMENSION(ims:ime, kms:kme, jms:jme), INTENT(INOUT) :: rucuten, &
1569 &  rvcuten, rthcuten, rqvcuten, rqccuten, rqrcuten, rqicuten, rqscuten, &
1570 &  rushten, rvshten, rthshten, rqvshten, rqcshten, rqrshten, rqishten, &
1571 &  rqsshten, rqgshten
1572   REAL, DIMENSION(ims:ime, kms:kme, jms:jme), INTENT(INOUT) :: rucutend&
1573 &  , rvcutend, rthcutend, rqvcutend, rqccutend, rqrcutend, rqicutend, &
1574 &  rqscutend, rushtend, rvshtend, rthshtend, rqvshtend, rqcshtend, &
1575 &  rqrshtend, rqishtend, rqsshtend, rqgshtend
1576 ! pbl
1577   REAL, DIMENSION(ims:ime, kms:kme, jms:jme), INTENT(INOUT) :: rublten, &
1578 &  rvblten, rthblten, rqvblten, rqcblten, rqiblten
1579   REAL, DIMENSION(ims:ime, kms:kme, jms:jme), INTENT(INOUT) :: rubltend&
1580 &  , rvbltend, rthbltend, rqvbltend, rqcbltend, rqibltend
1581 ! fdda
1582   REAL, DIMENSION(ims:ime, kms:kme, jms:jme), INTENT(INOUT) :: rundgdten&
1583 &  , rvndgdten, rthndgdten, rqvndgdten
1584   REAL, DIMENSION(ims:ime, kms:kme, jms:jme), INTENT(INOUT) :: &
1585 &  rundgdtend, rvndgdtend, rthndgdtend, rqvndgdtend
1586   REAL, DIMENSION(ims:ime, jms:jme), INTENT(INOUT) :: rmundgdten
1587   INTEGER :: i, k, j
1588   INTEGER :: itf, ktf, jtf, itsu, jtsv
1589   IF (ite .GT. ide - 1) THEN
1590     itf = ide - 1
1591   ELSE
1592     itf = ite
1593   END IF
1594   IF (jte .GT. jde - 1) THEN
1595     jtf = jde - 1
1596   ELSE
1597     jtf = jte
1598   END IF
1599   IF (kte .GT. kde - 1) THEN
1600     ktf = kde - 1
1601   ELSE
1602     ktf = kte
1603   END IF
1604   IF (its .LT. ids + 1) THEN
1605     itsu = ids + 1
1606   ELSE
1607     itsu = its
1608   END IF
1609   IF (jts .LT. jds + 1) THEN
1610     jtsv = jds + 1
1611   ELSE
1612     jtsv = jts
1613   END IF
1614 ! radiation
1615   IF (config_flags%ra_lw_physics .GT. 0 .OR. config_flags%ra_sw_physics &
1616 &      .GT. 0) THEN
1617     DO j=jts,jtf
1618       DO k=kts,ktf
1619         DO i=its,itf
1620           rthratend(i, k, j) = mud(i, j)*rthraten(i, k, j) + mu(i, j)*&
1621 &            rthratend(i, k, j)
1622           rthraten(i, k, j) = mu(i, j)*rthraten(i, k, j)
1623         END DO
1624       END DO
1625     END DO
1626   END IF
1627 ! cumulus
1628   IF (config_flags%cu_physics .GT. 0) THEN
1629     DO j=jts,jtf
1630       DO i=its,itf
1631         DO k=kts,ktf
1632           rucutend(i, k, j) = mud(i, j)*rucuten(i, k, j) + mu(i, j)*&
1633 &            rucutend(i, k, j)
1634           rucuten(i, k, j) = mu(i, j)*rucuten(i, k, j)
1635           rvcutend(i, k, j) = mud(i, j)*rvcuten(i, k, j) + mu(i, j)*&
1636 &            rvcutend(i, k, j)
1637           rvcuten(i, k, j) = mu(i, j)*rvcuten(i, k, j)
1638           rthcutend(i, k, j) = mud(i, j)*rthcuten(i, k, j) + mu(i, j)*&
1639 &            rthcutend(i, k, j)
1640           rthcuten(i, k, j) = mu(i, j)*rthcuten(i, k, j)
1641           rqvcutend(i, k, j) = mud(i, j)*rqvcuten(i, k, j) + mu(i, j)*&
1642 &            rqvcutend(i, k, j)
1643           rqvcuten(i, k, j) = mu(i, j)*rqvcuten(i, k, j)
1644         END DO
1645       END DO
1646     END DO
1647     IF (p_qc .GE. param_first_scalar) THEN
1648       DO j=jts,jtf
1649         DO i=its,itf
1650           DO k=kts,ktf
1651             rqccutend(i, k, j) = mud(i, j)*rqccuten(i, k, j) + mu(i, j)*&
1652 &              rqccutend(i, k, j)
1653             rqccuten(i, k, j) = mu(i, j)*rqccuten(i, k, j)
1654           END DO
1655         END DO
1656       END DO
1657     END IF
1658     IF (p_qr .GE. param_first_scalar) THEN
1659       DO j=jts,jtf
1660         DO i=its,itf
1661           DO k=kts,ktf
1662             rqrcutend(i, k, j) = mud(i, j)*rqrcuten(i, k, j) + mu(i, j)*&
1663 &              rqrcutend(i, k, j)
1664             rqrcuten(i, k, j) = mu(i, j)*rqrcuten(i, k, j)
1665           END DO
1666         END DO
1667       END DO
1668     END IF
1669     IF (p_qi .GE. param_first_scalar) THEN
1670       DO j=jts,jtf
1671         DO i=its,itf
1672           DO k=kts,ktf
1673             rqicutend(i, k, j) = mud(i, j)*rqicuten(i, k, j) + mu(i, j)*&
1674 &              rqicutend(i, k, j)
1675             rqicuten(i, k, j) = mu(i, j)*rqicuten(i, k, j)
1676           END DO
1677         END DO
1678       END DO
1679     END IF
1680     IF (p_qs .GE. param_first_scalar) THEN
1681       DO j=jts,jtf
1682         DO i=its,itf
1683           DO k=kts,ktf
1684             rqscutend(i, k, j) = mud(i, j)*rqscuten(i, k, j) + mu(i, j)*&
1685 &              rqscutend(i, k, j)
1686             rqscuten(i, k, j) = mu(i, j)*rqscuten(i, k, j)
1687           END DO
1688         END DO
1689       END DO
1690     END IF
1691   END IF
1692 ! shallow cumulus
1693   IF (config_flags%shcu_physics .GT. 0) THEN
1694     DO j=jts,jtf
1695       DO i=its,itf
1696         DO k=kts,ktf
1697           rushtend(i, k, j) = mud(i, j)*rushten(i, k, j) + mu(i, j)*&
1698 &            rushtend(i, k, j)
1699           rushten(i, k, j) = mu(i, j)*rushten(i, k, j)
1700           rvshtend(i, k, j) = mud(i, j)*rvshten(i, k, j) + mu(i, j)*&
1701 &            rvshtend(i, k, j)
1702           rvshten(i, k, j) = mu(i, j)*rvshten(i, k, j)
1703           rthshtend(i, k, j) = mud(i, j)*rthshten(i, k, j) + mu(i, j)*&
1704 &            rthshtend(i, k, j)
1705           rthshten(i, k, j) = mu(i, j)*rthshten(i, k, j)
1706           rqvshtend(i, k, j) = mud(i, j)*rqvshten(i, k, j) + mu(i, j)*&
1707 &            rqvshtend(i, k, j)
1708           rqvshten(i, k, j) = mu(i, j)*rqvshten(i, k, j)
1709         END DO
1710       END DO
1711     END DO
1712     IF (p_qc .GE. param_first_scalar) THEN
1713       DO j=jts,jtf
1714         DO i=its,itf
1715           DO k=kts,ktf
1716             rqcshtend(i, k, j) = mud(i, j)*rqcshten(i, k, j) + mu(i, j)*&
1717 &              rqcshtend(i, k, j)
1718             rqcshten(i, k, j) = mu(i, j)*rqcshten(i, k, j)
1719           END DO
1720         END DO
1721       END DO
1722     END IF
1723     IF (p_qr .GE. param_first_scalar) THEN
1724       DO j=jts,jtf
1725         DO i=its,itf
1726           DO k=kts,ktf
1727             rqrshtend(i, k, j) = mud(i, j)*rqrshten(i, k, j) + mu(i, j)*&
1728 &              rqrshtend(i, k, j)
1729             rqrshten(i, k, j) = mu(i, j)*rqrshten(i, k, j)
1730           END DO
1731         END DO
1732       END DO
1733     END IF
1734     IF (p_qi .GE. param_first_scalar) THEN
1735       DO j=jts,jtf
1736         DO i=its,itf
1737           DO k=kts,ktf
1738             rqishtend(i, k, j) = mud(i, j)*rqishten(i, k, j) + mu(i, j)*&
1739 &              rqishtend(i, k, j)
1740             rqishten(i, k, j) = mu(i, j)*rqishten(i, k, j)
1741           END DO
1742         END DO
1743       END DO
1744     END IF
1745     IF (p_qs .GE. param_first_scalar) THEN
1746       DO j=jts,jtf
1747         DO i=its,itf
1748           DO k=kts,ktf
1749             rqsshtend(i, k, j) = mud(i, j)*rqsshten(i, k, j) + mu(i, j)*&
1750 &              rqsshtend(i, k, j)
1751             rqsshten(i, k, j) = mu(i, j)*rqsshten(i, k, j)
1752           END DO
1753         END DO
1754       END DO
1755     END IF
1756     IF (p_qg .GE. param_first_scalar) THEN
1757       DO j=jts,jtf
1758         DO i=its,itf
1759           DO k=kts,ktf
1760             rqgshtend(i, k, j) = mud(i, j)*rqgshten(i, k, j) + mu(i, j)*&
1761 &              rqgshtend(i, k, j)
1762             rqgshten(i, k, j) = mu(i, j)*rqgshten(i, k, j)
1763           END DO
1764         END DO
1765       END DO
1766     END IF
1767   END IF
1768 ! pbl
1769   IF (config_flags%bl_pbl_physics .GT. 0) THEN
1770     DO j=jts,jtf
1771       DO k=kts,ktf
1772         DO i=its,itf
1773           rubltend(i, k, j) = mud(i, j)*rublten(i, k, j) + mu(i, j)*&
1774 &            rubltend(i, k, j)
1775           rublten(i, k, j) = mu(i, j)*rublten(i, k, j)
1776           rvbltend(i, k, j) = mud(i, j)*rvblten(i, k, j) + mu(i, j)*&
1777 &            rvbltend(i, k, j)
1778           rvblten(i, k, j) = mu(i, j)*rvblten(i, k, j)
1779           rthbltend(i, k, j) = mud(i, j)*rthblten(i, k, j) + mu(i, j)*&
1780 &            rthbltend(i, k, j)
1781           rthblten(i, k, j) = mu(i, j)*rthblten(i, k, j)
1782         END DO
1783       END DO
1784     END DO
1785     IF (p_qv .GE. param_first_scalar) THEN
1786       DO j=jts,jtf
1787         DO k=kts,ktf
1788           DO i=its,itf
1789             rqvbltend(i, k, j) = mud(i, j)*rqvblten(i, k, j) + mu(i, j)*&
1790 &              rqvbltend(i, k, j)
1791             rqvblten(i, k, j) = mu(i, j)*rqvblten(i, k, j)
1792           END DO
1793         END DO
1794       END DO
1795     END IF
1796     IF (p_qc .GE. param_first_scalar) THEN
1797       DO j=jts,jtf
1798         DO k=kts,ktf
1799           DO i=its,itf
1800             rqcbltend(i, k, j) = mud(i, j)*rqcblten(i, k, j) + mu(i, j)*&
1801 &              rqcbltend(i, k, j)
1802             rqcblten(i, k, j) = mu(i, j)*rqcblten(i, k, j)
1803           END DO
1804         END DO
1805       END DO
1806     END IF
1807     IF (p_qi .GE. param_first_scalar) THEN
1808       DO j=jts,jtf
1809         DO k=kts,ktf
1810           DO i=its,itf
1811             rqibltend(i, k, j) = mud(i, j)*rqiblten(i, k, j) + mu(i, j)*&
1812 &              rqibltend(i, k, j)
1813             rqiblten(i, k, j) = mu(i, j)*rqiblten(i, k, j)
1814           END DO
1815         END DO
1816       END DO
1817     END IF
1818   END IF
1819 ! fdda
1820 ! note fdda u and v tendencies are staggered, also only interior points have muu/muv,
1821 !   so only couple those
1822   IF (config_flags%grid_fdda .GT. 0) THEN
1823     DO j=jts,jtf
1824       DO k=kts,ktf
1825         DO i=itsu,itf
1826 !     if( i == itf/2 .AND. j == jtf/2 .AND. k == ktf/2 ) &
1827 !     write(*,'(a,3i6,e15.5)') 'u_ten before=',i,k,j, RUNDGDTEN(i,k,j)
1828           rundgdtend(i, k, j) = muud(i, j)*rundgdten(i, k, j) + muu(i, j&
1829 &            )*rundgdtend(i, k, j)
1830           rundgdten(i, k, j) = muu(i, j)*rundgdten(i, k, j)
1831         END DO
1832       END DO
1833     END DO
1834 !        if( i == itf/2 .AND. j == jtf/2 .AND. k==ktf/2 ) &
1835 !          write(*,'(a,2f15.5)') 'mu, muu=',mu(i,j), muu(i,j)
1836 !     if( i == itf/2 .AND. j == jtf/2 .AND. k == ktf/2 ) &
1837 !     write(*,'(a,3i6,e15.5)') 'u_ten after=',i,k,j, RUNDGDTEN(i,k,j)
1838 !     if( RUNDGDTEN(i,k,j) > 30.0 ) write(*,*) 'IKJ=',i,k,j
1839 !     write(*,'(a,e15.5)') 'u_ten MAXIMUM after=', maxval(RUNDGDTEN)
1840     DO j=jtsv,jtf
1841       DO k=kts,ktf
1842         DO i=its,itf
1843           rvndgdtend(i, k, j) = muvd(i, j)*rvndgdten(i, k, j) + muv(i, j&
1844 &            )*rvndgdtend(i, k, j)
1845           rvndgdten(i, k, j) = muv(i, j)*rvndgdten(i, k, j)
1846         END DO
1847       END DO
1848     END DO
1849     DO j=jts,jtf
1850       DO k=kts,ktf
1851         DO i=its,itf
1852 !     if( i == itf/2 .AND. j == jtf/2 .AND. k == ktf/2 ) &
1853 !     write(*,'(a,3i6,e15.5)') 'th before=',i,k,j, RTHNDGDTEN(I,K,J)
1854           rthndgdtend(i, k, j) = mud(i, j)*rthndgdten(i, k, j) + mu(i, j&
1855 &            )*rthndgdtend(i, k, j)
1856           rthndgdten(i, k, j) = mu(i, j)*rthndgdten(i, k, j)
1857         END DO
1858       END DO
1859     END DO
1860 !        RMUNDGDTEN(I,J) - no coupling
1861 !     if( i == itf/2 .AND. j == jtf/2 .AND. k == ktf/2 ) &
1862 !     write(*,'(a,3i6,e15.5)') 'th after=',i,k,j, RTHNDGDTEN(I,K,J)
1863     IF (p_qv .GE. param_first_scalar) THEN
1864       DO j=jts,jtf
1865         DO k=kts,ktf
1866           DO i=its,itf
1867             rqvndgdtend(i, k, j) = mud(i, j)*rqvndgdten(i, k, j) + mu(i&
1868 &              , j)*rqvndgdtend(i, k, j)
1869             rqvndgdten(i, k, j) = mu(i, j)*rqvndgdten(i, k, j)
1870           END DO
1871         END DO
1872       END DO
1873     END IF
1874   END IF
1875 END SUBROUTINE G_CALCULATE_PHY_TEND
1876
1877 !-----------------------------------------------------------------------
1878
1879 ! Revised by Ning Pan, 2010-08-03
1880 ! SUBROUTINE g_bound_tke(tke,g_tke,tke_upper_bound,g_tke_upper_bound,ids,ide, &
1881  SUBROUTINE g_bound_tke(tke,g_tke,tke_upper_bound,ids,ide, &
1882  jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
1883
1884  IMPLICIT NONE
1885
1886  REAL :: Tmpv1,g_Tmpv1
1887  INTEGER :: ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte
1888  REAL,DIMENSION(ims:ime,kms:kme,jms:jme) :: tke,g_tke
1889 ! Revised by Ning Pan, 2010-08-03
1890 ! REAL :: tke_upper_bound,g_tke_upper_bound
1891  REAL :: tke_upper_bound
1892  INTEGER :: i,k,j
1893
1894  DO j =jts,min(jte,jde-1)
1895  DO k =kts,kte-1
1896  DO i =its,min(ite,ide-1)
1897
1898 ! g_tke(i,k,j) =(g_tke_upper_bound +(g_tke(i,k,j) +0.0 +(g_tke(i,k,j) -0.0) &
1899 ! *sign(1.0, tke(i,k,j) -(0.)))*0.5 -(g_tke_upper_bound -(g_tke(i,k,j) &
1900 ! +0.0 +(g_tke(i,k,j) -0.0)*sign(1.0, tke(i,k,j) -(0.)))*0.5)*sign(1.0, &
1901 ! tke_upper_bound -(max(tke(i,k,j),0.))))*0.5
1902  g_tke(i,k,j) =((g_tke(i,k,j) +0.0 +(g_tke(i,k,j) -0.0) &
1903  *sign(1.0, tke(i,k,j) -(0.)))*0.5 -(-(g_tke(i,k,j) &
1904  +0.0 +(g_tke(i,k,j) -0.0)*sign(1.0, tke(i,k,j) -(0.)))*0.5)*sign(1.0, &
1905  tke_upper_bound -(max(tke(i,k,j),0.))))*0.5
1906  tke(i,k,j) =min(tke_upper_bound,max(tke(i,k,j),0.))
1907
1908  ENDDO
1909  ENDDO
1910  ENDDO
1911
1912  END SUBROUTINE g_bound_tke
1913
1914
1915 END MODULE g_module_em