SCALE-RM
Functions/Subroutines
scale_atmos_dyn_tstep_short_fvm_heve Module Reference

module Atmosphere / Dynamics RK More...

Functions/Subroutines

subroutine, public atmos_dyn_tstep_short_fvm_heve_regist (ATMOS_DYN_TYPE, VA_out, VAR_NAME, VAR_DESC, VAR_UNIT)
 Register. More...
 
subroutine, public atmos_dyn_tstep_short_fvm_heve_setup
 Setup. More...
 
subroutine, public atmos_dyn_tstep_short_fvm_heve (DENS_RK, MOMZ_RK, MOMX_RK, MOMY_RK, RHOT_RK, PROG_RK, mflx_hi, tflx_hi, DENS0, MOMZ0, MOMX0, MOMY0, RHOT0, DENS, MOMZ, MOMX, MOMY, RHOT, DENS_t, MOMZ_t, MOMX_t, MOMY_t, RHOT_t, PROG0, PROG, DPRES0, RT2P, CORIOLI, num_diff, wdamp_coef, divdmp_coef, DDIV, FLAG_FCT_MOMENTUM, FLAG_FCT_T, FLAG_FCT_ALONG_STREAM, CDZ, FDZ, FDX, FDY, RCDZ, RCDX, RCDY, RFDZ, RFDX, RFDY, PHI, GSQRT, J13G, J23G, J33G, MAPF, REF_dens, REF_rhot, BND_W, BND_E, BND_S, BND_N, TwoD, dtrk, last)
 

Detailed Description

module Atmosphere / Dynamics RK

Description
HEVE FVM scheme for Atmospheric dynamical process
Author
Team SCALE
NAMELIST
  • No namelist group
History Output
namedescriptionunitvariable
DENS_t_advch tendency of density (horiz. advection) (w/ HIST_TEND) kg/m3/s advch_t
DENS_t_advcv tendency of density (vert. advection) (w/ HIST_TEND) kg/m3/s advcv_t
MOMX_t_advch tendency of momentum x (horiz. advection) (w/ HIST_TEND) kg/m2/s2 advch_t
MOMX_t_advcv tendency of momentum x (vert. advection) (w/ HIST_TEND) kg/m2/s2 advcv_t
MOMX_t_cf tendency of momentum x (coliolis force) (w/ HIST_TEND) kg/m2/s2 cf_t
MOMX_t_ddiv tendency of momentum x (divergence damping) (w/ HIST_TEND) kg/m2/s2 ddiv_t
MOMX_t_pg tendency of momentum x (pressure gradient) (w/ HIST_TEND) kg/m2/s2 pg_t
MOMY_t_advch tendency of momentum y (horiz. advection) (w/ HIST_TEND) kg/m2/s2 advch_t
MOMY_t_advcv tendency of momentum y (vert. advection) (w/ HIST_TEND) kg/m2/s2 advcv_t
MOMY_t_cf tendency of momentum y (coliolis force) (w/ HIST_TEND) kg/m2/s2 cf_t
MOMY_t_ddiv tendency of momentum y (divergence damping) (w/ HIST_TEND) kg/m2/s2 ddiv_t
MOMY_t_pg tendency of momentum y (pressure gradient) (w/ HIST_TEND) kg/m2/s2 pg_t
MOMZ_t_advch tendency of momentum z (horiz. advection) (w/ HIST_TEND) kg/m2/s2 advch_t
MOMZ_t_advcv tendency of momentum z (vert. advection) (w/ HIST_TEND) kg/m2/s2 advcv_t
MOMZ_t_ddiv tendency of momentum z (divergence damping) (w/ HIST_TEND) kg/m2/s2 ddiv_t
MOMZ_t_pg tendency of momentum z (pressure gradient) (w/ HIST_TEND) kg/m2/s2 pg_t
MOMZ_t_wdamp tendency of momentum z (Rayleigh damping) (w/ HIST_TEND) kg/m2/s2 wdmp_t
RHOT_t_advch tendency of rho*theta (horiz. advection) (w/ HIST_TEND) K kg/m3/s advch_t
RHOT_t_advcv tendency of rho*theta (vert. advection) (w/ HIST_TEND) K kg/m3/s advcv_t

Function/Subroutine Documentation

◆ atmos_dyn_tstep_short_fvm_heve_regist()

subroutine, public scale_atmos_dyn_tstep_short_fvm_heve::atmos_dyn_tstep_short_fvm_heve_regist ( character(len=*), intent(in)  ATMOS_DYN_TYPE,
integer, intent(out)  VA_out,
character(len=h_short), dimension(:), intent(out)  VAR_NAME,
character(len=h_mid), dimension(:), intent(out)  VAR_DESC,
character(len=h_short), dimension(:), intent(out)  VAR_UNIT 
)

Register.

Parameters
[out]va_outnumber of prognostic variables
[out]var_namename of the variables
[out]var_descdesc. of the variables
[out]var_unitunit of the variables

Definition at line 72 of file scale_atmos_dyn_tstep_short_fvm_heve.F90.

72  use scale_prc, only: &
73  prc_abort
74  implicit none
75 
76  character(len=*), intent(in) :: ATMOS_DYN_TYPE
77  integer, intent(out) :: VA_out
78  character(len=H_SHORT), intent(out) :: VAR_NAME(:)
79  character(len=H_MID), intent(out) :: VAR_DESC(:)
80  character(len=H_SHORT), intent(out) :: VAR_UNIT(:)
81  !---------------------------------------------------------------------------
82 
83  if ( atmos_dyn_type /= 'FVM-HEVE' .AND. atmos_dyn_type /= 'HEVE' ) then
84  log_error("ATMOS_DYN_Tstep_short_fvm_heve_regist",*) 'ATMOS_DYN_TYPE is not FVM-HEVE. Check!'
85  call prc_abort
86  endif
87 
88  va_out = va_fvm_heve
89  var_name(:) = ""
90  var_desc(:) = ""
91  var_unit(:) = ""
92 
93  log_newline
94  log_info("ATMOS_DYN_Tstep_short_fvm_heve_regist",*) 'Register additional prognostic variables (HEVE)'
95  if ( va_out < 1 ) then
96  log_info_cont(*) '=> nothing.'
97  endif
98 
99  return

References scale_prc::prc_abort().

Referenced by scale_atmos_dyn_tstep_short::atmos_dyn_tstep_short_regist().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ atmos_dyn_tstep_short_fvm_heve_setup()

subroutine, public scale_atmos_dyn_tstep_short_fvm_heve::atmos_dyn_tstep_short_fvm_heve_setup

Setup.

Definition at line 105 of file scale_atmos_dyn_tstep_short_fvm_heve.F90.

105  return

Referenced by scale_atmos_dyn_tstep_short::atmos_dyn_tstep_short_regist().

Here is the caller graph for this function:

◆ atmos_dyn_tstep_short_fvm_heve()

subroutine, public scale_atmos_dyn_tstep_short_fvm_heve::atmos_dyn_tstep_short_fvm_heve ( real(rp), dimension (ka,ia,ja), intent(out)  DENS_RK,
real(rp), dimension (ka,ia,ja), intent(out)  MOMZ_RK,
real(rp), dimension (ka,ia,ja), intent(out)  MOMX_RK,
real(rp), dimension (ka,ia,ja), intent(out)  MOMY_RK,
real(rp), dimension (ka,ia,ja), intent(out)  RHOT_RK,
real(rp), dimension (ka,ia,ja,va), intent(out)  PROG_RK,
real(rp), dimension (ka,ia,ja,3), intent(inout)  mflx_hi,
real(rp), dimension (ka,ia,ja,3), intent(out)  tflx_hi,
real(rp), dimension (ka,ia,ja), intent(in), target  DENS0,
real(rp), dimension (ka,ia,ja), intent(in), target  MOMZ0,
real(rp), dimension (ka,ia,ja), intent(in), target  MOMX0,
real(rp), dimension (ka,ia,ja), intent(in), target  MOMY0,
real(rp), dimension (ka,ia,ja), intent(in), target  RHOT0,
real(rp), dimension (ka,ia,ja), intent(in)  DENS,
real(rp), dimension (ka,ia,ja), intent(in)  MOMZ,
real(rp), dimension (ka,ia,ja), intent(in)  MOMX,
real(rp), dimension (ka,ia,ja), intent(in)  MOMY,
real(rp), dimension (ka,ia,ja), intent(in)  RHOT,
real(rp), dimension (ka,ia,ja), intent(in)  DENS_t,
real(rp), dimension (ka,ia,ja), intent(in)  MOMZ_t,
real(rp), dimension (ka,ia,ja), intent(in)  MOMX_t,
real(rp), dimension (ka,ia,ja), intent(in)  MOMY_t,
real(rp), dimension (ka,ia,ja), intent(in)  RHOT_t,
real(rp), dimension (ka,ia,ja,va), intent(in)  PROG0,
real(rp), dimension (ka,ia,ja,va), intent(in)  PROG,
real(rp), dimension (ka,ia,ja), intent(in)  DPRES0,
real(rp), dimension (ka,ia,ja), intent(in)  RT2P,
real(rp), dimension ( ia,ja), intent(in)  CORIOLI,
real(rp), dimension(ka,ia,ja,5,3), intent(in)  num_diff,
real(rp), dimension(ka), intent(in)  wdamp_coef,
real(rp), intent(in)  divdmp_coef,
real(rp), dimension (ka,ia,ja), intent(in)  DDIV,
logical, intent(in)  FLAG_FCT_MOMENTUM,
logical, intent(in)  FLAG_FCT_T,
logical, intent(in)  FLAG_FCT_ALONG_STREAM,
real(rp), dimension (ka), intent(in)  CDZ,
real(rp), dimension (ka-1), intent(in)  FDZ,
real(rp), dimension (ia-1), intent(in)  FDX,
real(rp), dimension (ja-1), intent(in)  FDY,
real(rp), dimension(ka), intent(in)  RCDZ,
real(rp), dimension(ia), intent(in)  RCDX,
real(rp), dimension(ja), intent(in)  RCDY,
real(rp), dimension(ka-1), intent(in)  RFDZ,
real(rp), dimension(ia-1), intent(in)  RFDX,
real(rp), dimension(ja-1), intent(in)  RFDY,
real(rp), dimension (ka,ia,ja), intent(in)  PHI,
real(rp), dimension (ka,ia,ja,7), intent(in)  GSQRT,
real(rp), dimension (ka,ia,ja,7), intent(in)  J13G,
real(rp), dimension (ka,ia,ja,7), intent(in)  J23G,
real(rp), intent(in)  J33G,
real(rp), dimension (ia,ja,2,4), intent(in)  MAPF,
real(rp), dimension(ka,ia,ja), intent(in)  REF_dens,
real(rp), dimension(ka,ia,ja), intent(in)  REF_rhot,
logical, intent(in)  BND_W,
logical, intent(in)  BND_E,
logical, intent(in)  BND_S,
logical, intent(in)  BND_N,
logical, intent(in)  TwoD,
real(rp), intent(in)  dtrk,
logical, intent(in)  last 
)
Parameters
[in]phigeopotential
[in]gsqrtvertical metrics {G}^1/2
[in]j13g(1,3) element of Jacobian matrix
[in]j23g(2,3) element of Jacobian matrix
[in]j33g(3,3) element of Jacobian matrix
[in]mapfmap factor
[in]ref_densreference density

Definition at line 127 of file scale_atmos_dyn_tstep_short_fvm_heve.F90.

128  use scale_const, only: &
129  eps => const_eps, &
130  grav => const_grav, &
131  p00 => const_pre00
132  use scale_comm_cartesc, only: &
133  comm_vars8, &
134  comm_wait
135  use scale_atmos_dyn_fvm_fct, only: &
137  use scale_atmos_dyn_fvm_flux_ud1, only: &
150  use scale_atmos_dyn_fvm_flux, only: &
169 #ifdef HIST_TEND
170  use scale_file_history, only: &
171  file_history_in
172 #endif
173  implicit none
174 
175  real(RP), intent(out) :: DENS_RK (KA,IA,JA) ! prognostic variables
176  real(RP), intent(out) :: MOMZ_RK (KA,IA,JA) !
177  real(RP), intent(out) :: MOMX_RK (KA,IA,JA) !
178  real(RP), intent(out) :: MOMY_RK (KA,IA,JA) !
179  real(RP), intent(out) :: RHOT_RK (KA,IA,JA) !
180  real(RP), intent(out) :: PROG_RK (KA,IA,JA,VA) !
181 
182  real(RP), intent(inout) :: mflx_hi (KA,IA,JA,3) ! mass flux
183  real(RP), intent(out) :: tflx_hi (KA,IA,JA,3) ! internal energy flux
184 
185  real(RP), intent(in), target :: DENS0 (KA,IA,JA) ! prognostic variables at previous dynamical time step
186  real(RP), intent(in), target :: MOMZ0 (KA,IA,JA) !
187  real(RP), intent(in), target :: MOMX0 (KA,IA,JA) !
188  real(RP), intent(in), target :: MOMY0 (KA,IA,JA) !
189  real(RP), intent(in), target :: RHOT0 (KA,IA,JA) !
190  real(RP), intent(in) :: PROG0 (KA,IA,JA,VA)
191 
192  real(RP), intent(in) :: DENS (KA,IA,JA) ! prognostic variables at previous RK step
193  real(RP), intent(in) :: MOMZ (KA,IA,JA) !
194  real(RP), intent(in) :: MOMX (KA,IA,JA) !
195  real(RP), intent(in) :: MOMY (KA,IA,JA) !
196  real(RP), intent(in) :: RHOT (KA,IA,JA) !
197  real(RP), intent(in) :: PROG (KA,IA,JA,VA)
198 
199  real(RP), intent(in) :: DENS_t (KA,IA,JA) ! tendency
200  real(RP), intent(in) :: MOMZ_t (KA,IA,JA) !
201  real(RP), intent(in) :: MOMX_t (KA,IA,JA) !
202  real(RP), intent(in) :: MOMY_t (KA,IA,JA) !
203  real(RP), intent(in) :: RHOT_t (KA,IA,JA) !
204 
205  real(RP), intent(in) :: DPRES0 (KA,IA,JA)
206  real(RP), intent(in) :: RT2P (KA,IA,JA)
207  real(RP), intent(in) :: CORIOLI ( IA,JA)
208  real(RP), intent(in) :: num_diff(KA,IA,JA,5,3)
209  real(RP), intent(in) :: wdamp_coef(KA)
210  real(RP), intent(in) :: divdmp_coef
211  real(RP), intent(in) :: DDIV (KA,IA,JA)
212 
213  logical, intent(in) :: FLAG_FCT_MOMENTUM
214  logical, intent(in) :: FLAG_FCT_T
215  logical, intent(in) :: FLAG_FCT_ALONG_STREAM
216 
217  real(RP), intent(in) :: CDZ (KA)
218  real(RP), intent(in) :: FDZ (KA-1)
219  real(RP), intent(in) :: FDX (IA-1)
220  real(RP), intent(in) :: FDY (JA-1)
221  real(RP), intent(in) :: RCDZ(KA)
222  real(RP), intent(in) :: RCDX(IA)
223  real(RP), intent(in) :: RCDY(JA)
224  real(RP), intent(in) :: RFDZ(KA-1)
225  real(RP), intent(in) :: RFDX(IA-1)
226  real(RP), intent(in) :: RFDY(JA-1)
227 
228  real(RP), intent(in) :: PHI (KA,IA,JA)
229  real(RP), intent(in) :: GSQRT (KA,IA,JA,7)
230  real(RP), intent(in) :: J13G (KA,IA,JA,7)
231  real(RP), intent(in) :: J23G (KA,IA,JA,7)
232  real(RP), intent(in) :: J33G
233  real(RP), intent(in) :: MAPF (IA,JA,2,4)
234  real(RP), intent(in) :: REF_dens(KA,IA,JA)
235  real(RP), intent(in) :: REF_rhot(KA,IA,JA)
236 
237  logical, intent(in) :: BND_W
238  logical, intent(in) :: BND_E
239  logical, intent(in) :: BND_S
240  logical, intent(in) :: BND_N
241  logical, intent(in) :: TwoD
242 
243  real(RP), intent(in) :: dtrk
244  logical, intent(in) :: last
245 
246  ! diagnostic variables
247  real(RP) :: VELZ (KA,IA,JA) ! velocity w [m/s]
248  real(RP) :: VELX (KA,IA,JA) ! velocity u [m/s]
249  real(RP) :: VELY (KA,IA,JA) ! velocity v [m/s]
250  real(RP) :: POTT (KA,IA,JA) ! potential temperature [K]
251  real(RP) :: DPRES(KA,IA,JA) ! pressure - reference pressure
252 
253  real(RP) :: qflx_J13(KA,IA,JA)
254  real(RP) :: qflx_J23(KA,IA,JA)
255  real(RP) :: pgf (KA,IA,JA) ! pressure gradient force
256  real(RP) :: buoy (KA,IA,JA) ! buoyancy force
257  real(RP) :: cor (KA,IA,JA) ! Coriolis force
258 
259  ! flux
260  real(RP) :: qflx_hi (KA,IA,JA,3)
261 #ifndef NO_FCT_DYN
262  real(RP) :: qflx_lo (KA,IA,JA,3)
263  real(RP) :: qflx_anti(KA,IA,JA,3)
264  real(RP) :: tflx_lo (KA,IA,JA,3)
265  real(RP) :: tflx_anti(KA,IA,JA,3)
266  real(RP) :: DENS0_uvw(KA,IA,JA)
267  real(RP) :: DENS_uvw (KA,IA,JA)
268 #endif
269  real(RP) :: advch ! horizontal advection
270  real(RP) :: advcv ! vertical advection
271  real(RP) :: wdamp ! rayleight damping for W
272  real(RP) :: div ! divergence damping
273 #ifdef HIST_TEND
274  real(RP) :: advch_t(KA,IA,JA,5)
275  real(RP) :: advcv_t(KA,IA,JA,5)
276  real(RP) :: wdmp_t(KA,IA,JA)
277  real(RP) :: ddiv_t(KA,IA,JA,3)
278  real(RP) :: pg_t(KA,IA,JA,3)
279  real(RP) :: cf_t(KA,IA,JA,2)
280  logical :: lhist
281 #endif
282 
283  integer :: IIS, IIE
284  integer :: JJS, JJE
285  integer :: IFS_OFF, JFS_OFF
286  integer :: k, i, j
287  !---------------------------------------------------------------------------
288 
289 #ifdef DEBUG
290  velz(:,:,:) = undef
291  velx(:,:,:) = undef
292  vely(:,:,:) = undef
293  pott(:,:,:) = undef
294 
295  dpres(:,:,:) = undef
296 
297 #ifndef NO_FCT_DYN
298  qflx_lo(:,:,:,:) = undef
299  qflx_anti(:,:,:,:) = undef
300  tflx_lo(:,:,:,:) = undef
301  tflx_anti(:,:,:,:) = undef
302 #endif
303 #endif
304 
305 #if defined DEBUG || defined QUICKDEBUG
306  dens_rk( 1:ks-1,:,:) = undef
307  dens_rk(ke+1:ka ,:,:) = undef
308  momz_rk( 1:ks-1,:,:) = undef
309  momz_rk(ke+1:ka ,:,:) = undef
310  momx_rk( 1:ks-1,:,:) = undef
311  momx_rk(ke+1:ka ,:,:) = undef
312  momy_rk( 1:ks-1,:,:) = undef
313  momy_rk(ke+1:ka ,:,:) = undef
314  rhot_rk( 1:ks-1,:,:) = undef
315  rhot_rk(ke+1:ka ,:,:) = undef
316  prog_rk( 1:ks-1,:,:,:) = undef
317  prog_rk(ke+1:ka ,:,:,:) = undef
318 #endif
319 
320 #ifdef HIST_TEND
321  advch_t = 0.0_rp
322  advcv_t = 0.0_rp
323  wdmp_t = 0.0_rp
324  ddiv_t = 0.0_rp
325  pg_t = 0.0_rp
326  cf_t = 0.0_rp
327 
328  lhist = last
329 #endif
330 
331  ifs_off = 1
332  jfs_off = 1
333  if ( bnd_w ) ifs_off = 0
334  if ( bnd_s ) jfs_off = 0
335 
336 
337  do jjs = js, je, jblock
338  jje = jjs+jblock-1
339  do iis = is, ie, iblock
340  iie = iis+iblock-1
341 
342  ! pressure, pott. temp.
343 
344  !$omp parallel private(i,j,k)
345 
346  !$omp do OMP_SCHEDULE_ collapse(2)
347  do j = jjs-1, jje+1
348  do i = max(iis-1,1), min(iie+1,ia)
349  do k = ks, ke
350 #ifdef DEBUG
351  call check( __line__, dpres0(k,i,j) )
352  call check( __line__, rt2p(k,i,j) )
353  call check( __line__, rhot(k,i,j) )
354  call check( __line__, ref_rhot(k,i,j) )
355 #endif
356  dpres(k,i,j) = dpres0(k,i,j) + rt2p(k,i,j) * ( rhot(k,i,j) - ref_rhot(k,i,j) )
357  enddo
358  dpres(ks-1,i,j) = dpres0(ks-1,i,j) - dens(ks,i,j) * ( phi(ks-1,i,j) - phi(ks+1,i,j) )
359  dpres(ke+1,i,j) = dpres0(ke+1,i,j) - dens(ke,i,j) * ( phi(ke+1,i,j) - phi(ke-1,i,j) )
360  enddo
361  enddo
362  !$omp end do nowait
363 #ifdef DEBUG
364  k = iundef; i = iundef; j = iundef
365 #endif
366  !$omp do OMP_SCHEDULE_ collapse(2)
367  do j = jjs-jhalo, jje+jhalo
368  do i = iis-ihalo, iie+ihalo
369  do k = ks, ke
370 #ifdef DEBUG
371  call check( __line__, rhot(k,i,j) )
372  call check( __line__, dens(k,i,j) )
373 #endif
374  pott(k,i,j) = rhot(k,i,j) / dens(k,i,j)
375  enddo
376  enddo
377  enddo
378  !$omp end do nowait
379 
380  !$omp end parallel
381 #ifdef DEBUG
382  k = iundef; i = iundef; j = iundef
383 #endif
384 
385 #ifndef NO_FCT_DYN
386  enddo
387  enddo
388 
389  if ( flag_fct_momentum ) then
390  ! momentum -> velocity
391  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
392  do j = js-1, je+2
393  do i = max(is-1,1), min(ie+2,ia)
394  do k = ks, ke-1
395 #ifdef DEBUG
396  call check( __line__, momz0(k,i,j) )
397  call check( __line__, dens0(k ,i,j) )
398  call check( __line__, dens0(k+1,i,j) )
399 #endif
400  velz(k,i,j) = 2.0_rp * momz0(k,i,j) / ( dens0(k+1,i,j)+dens0(k,i,j) )
401  enddo
402  enddo
403  enddo
404 #ifdef DEBUG
405  k = iundef; i = iundef; j = iundef
406 #endif
407  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
408  do j = js-1, je+2
409  do i = max(is-1,1), min(ie+2,ia)
410  velz(ke,i,j) = 0.0_rp
411  enddo
412  enddo
413 #ifdef DEBUG
414  k = iundef; i = iundef; j = iundef
415 #endif
416 
417  if ( twod ) then
418  !$omp parallel do private(j,k) OMP_SCHEDULE_
419  do j = js-1, je+2
420  do k = ks, ke
421 #ifdef DEBUG
422  call check( __line__, momx0(k,is,j) )
423  call check( __line__, dens0(k,is,j) )
424 #endif
425  velx(k,is,j) = momx0(k,is,j) / dens0(k,is,j)
426  enddo
427  enddo
428  else
429  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
430  do j = js-1, je+2
431  do i = is-2, ie+1
432  do k = ks, ke
433 #ifdef DEBUG
434  call check( __line__, momx0(k,i,j) )
435  call check( __line__, dens0(k,i ,j) )
436  call check( __line__, dens0(k,i+1,j) )
437 #endif
438  velx(k,i,j) = 2.0_rp * momx0(k,i,j) / ( dens0(k,i+1,j)+dens0(k,i,j) )
439  enddo
440  enddo
441  enddo
442  end if
443 #ifdef DEBUG
444  k = iundef; i = iundef; j = iundef
445 #endif
446 
447  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
448  do j = js-2, je+1
449  do i = max(is-1,1), min(ie+2,ia)
450  do k = ks, ke
451 #ifdef DEBUG
452  call check( __line__, momy0(k,i,j) )
453  call check( __line__, dens0(k,i,j ) )
454  call check( __line__, dens0(k,i,j+1) )
455 #endif
456  vely(k,i,j) = 2.0_rp * momy0(k,i,j) / ( dens0(k,i,j+1)+dens0(k,i,j) )
457  enddo
458  enddo
459  enddo
460 #ifdef DEBUG
461  k = iundef; i = iundef; j = iundef
462 #endif
463  call comm_vars8( velz(:,:,:), 4 )
464  call comm_vars8( velx(:,:,:), 5 )
465  call comm_vars8( vely(:,:,:), 6 )
466  endif
467 
468  do jjs = js, je, jblock
469  jje = jjs+jblock-1
470  do iis = is, ie, iblock
471  iie = iis+iblock-1
472 #endif
473 
474  !########################################################################
475  ! continuity equation (total rho)
476  !########################################################################
477 
478  !$omp parallel private(i,j,k,advcv,advch)
479 
480  !-----< high order flux >-----
481 
482  ! at (x, y, w)
483 
484  if ( twod ) then
485  !$omp do OMP_SCHEDULE_ collapse(2)
486  do j = jjs, jje
487  do i = iis, iie
488  do k = ks, ke-1
489 #ifdef DEBUG
490  call check( __line__, momz(k+1,i,j) )
491  call check( __line__, momz(k ,i,j) )
492  call check( __line__, momz(k-1,i,j) )
493  call check( __line__, num_diff(k,i,j,i_dens,zdir) )
494 #endif
495  mflx_hi(k,i,j,zdir) = j33g * momz(k,i,j) / ( mapf(i,j,1,i_xy)*mapf(i,j,2,i_xy) ) &
496  + j23g(k,i,j,i_xyw) * 0.25_rp * ( momy(k+1,i,j)+momy(k+1,i,j-1) &
497  + momy(k ,i,j)+momy(k ,i,j-1) ) &
498  / mapf(i,j,1,i_xy) & ! [{x,v,z->x,y,w}]
499  + gsqrt(k,i,j,i_xyw) * num_diff(k,i,j,i_dens,zdir) / ( mapf(i,j,1,i_xy)*mapf(i,j,2,i_xy) )
500  enddo
501  enddo
502  enddo
503  !$omp end do nowait
504  else
505  !$omp do OMP_SCHEDULE_ collapse(2)
506  do j = jjs, jje
507  do i = iis, iie
508  do k = ks, ke-1
509 #ifdef DEBUG
510  call check( __line__, momz(k+1,i,j) )
511  call check( __line__, momz(k ,i,j) )
512  call check( __line__, momz(k-1,i,j) )
513  call check( __line__, num_diff(k,i,j,i_dens,zdir) )
514 #endif
515  mflx_hi(k,i,j,zdir) = j33g * momz(k,i,j) / ( mapf(i,j,1,i_xy)*mapf(i,j,2,i_xy) ) &
516  + j13g(k,i,j,i_xyw) * 0.25_rp * ( momx(k+1,i,j)+momx(k+1,i-1,j) &
517  + momx(k ,i,j)+momx(k ,i-1,j) ) &
518  / mapf(i,j,2,i_xy) & ! [{u,y,z->x,y,w}]
519  + j23g(k,i,j,i_xyw) * 0.25_rp * ( momy(k+1,i,j)+momy(k+1,i,j-1) &
520  + momy(k ,i,j)+momy(k ,i,j-1) ) &
521  / mapf(i,j,1,i_xy) & ! [{x,v,z->x,y,w}]
522  + gsqrt(k,i,j,i_xyw) * num_diff(k,i,j,i_dens,zdir) / ( mapf(i,j,1,i_xy)*mapf(i,j,2,i_xy) )
523  enddo
524  enddo
525  enddo
526  !$omp end do nowait
527  end if
528 #ifdef DEBUG
529  k = iundef; i = iundef; j = iundef
530 #endif
531 
532  !$omp do OMP_SCHEDULE_ collapse(2)
533  do j = jjs, jje
534  do i = iis, iie
535  mflx_hi(ks-1,i,j,zdir) = 0.0_rp
536  mflx_hi(ke ,i,j,zdir) = 0.0_rp
537  enddo
538  enddo
539  !$omp end do nowait
540 #ifdef DEBUG
541  k = iundef; i = iundef; j = iundef
542 #endif
543 
544  ! at (u, y, z)
545 
546  if ( .not. twod ) then
547  !$omp do OMP_SCHEDULE_ collapse(2)
548  do j = jjs, jje
549  do i = iis-ifs_off, min(iie,ieh)
550  do k = ks, ke
551 #ifdef DEBUG
552  call check( __line__, momx(k,i+1,j) )
553  call check( __line__, momx(k,i ,j) )
554  call check( __line__, momx(k,i-1,j) )
555  call check( __line__, num_diff(k,i,j,i_dens,xdir) )
556 #endif
557  mflx_hi(k,i,j,xdir) = gsqrt(k,i,j,i_uyz) / mapf(i,j,2,i_uy) &
558  * ( momx(k,i,j) + num_diff(k,i,j,i_dens,xdir) )
559  enddo
560  enddo
561  enddo
562  !$omp end do nowait
563  end if
564 #ifdef DEBUG
565  k = iundef; i = iundef; j = iundef
566 #endif
567 
568  ! at (x, v, z)
569 
570  !$omp do OMP_SCHEDULE_ collapse(2)
571  do j = jjs-jfs_off, min(jje,jeh)
572  do i = iis , iie
573  do k = ks, ke
574 #ifdef DEBUG
575  call check( __line__, momy(k,i,j+1) )
576  call check( __line__, momy(k,i,j ) )
577  call check( __line__, momy(k,i,j-1) )
578  call check( __line__, num_diff(k,i,j,i_dens,ydir) )
579 #endif
580  mflx_hi(k,i,j,ydir) = gsqrt(k,i,j,i_xvz) / mapf(i,j,1,i_xv) &
581  * ( momy(k,i,j) + num_diff(k,i,j,i_dens,ydir) )
582  enddo
583  enddo
584  enddo
585  !$omp end do
586 
587 #ifdef DEBUG
588  k = iundef; i = iundef; j = iundef
589 #endif
590 
591  !-----< update density >-----
592 
593  if ( twod ) then
594  !$omp do OMP_SCHEDULE_
595  do j = jjs, jje
596  do k = ks, ke
597 #ifdef DEBUG
598  call check( __line__, dens0(k,is,j) )
599  call check( __line__, mflx_hi(k ,is,j ,zdir) )
600  call check( __line__, mflx_hi(k-1,is,j ,zdir) )
601  call check( __line__, mflx_hi(k ,is,j ,ydir) )
602  call check( __line__, mflx_hi(k ,is,j-1,ydir) )
603  call check( __line__, dens_t(k,is,j) )
604 #endif
605  advcv = - ( mflx_hi(k,is,j,zdir)-mflx_hi(k-1,is,j, zdir) ) * rcdz(k)
606  advch = - ( mflx_hi(k,is,j,ydir)-mflx_hi(k ,is,j-1,ydir) ) * rcdy(j)
607  dens_rk(k,is,j) = dens0(k,is,j) &
608  + dtrk * ( ( advcv + advch ) * mapf(is,j,2,i_xy) / gsqrt(k,is,j,i_xyz) &
609  + dens_t(k,is,j) )
610 #ifdef HIST_TEND
611  if ( lhist ) then
612  advcv_t(k,is,j,i_dens) = advcv * mapf(is,j,2,i_xy) / gsqrt(k,is,j,i_xyz)
613  advch_t(k,is,j,i_dens) = advch * mapf(is,j,2,i_xy) / gsqrt(k,is,j,i_xyz)
614  endif
615 #endif
616  enddo
617  enddo
618  !$omp end do nowait
619  else
620  !$omp do OMP_SCHEDULE_ collapse(2)
621  do j = jjs, jje
622  do i = iis, iie
623  do k = ks, ke
624 #ifdef DEBUG
625  call check( __line__, dens0(k,i,j) )
626  call check( __line__, mflx_hi(k ,i ,j ,zdir) )
627  call check( __line__, mflx_hi(k-1,i ,j ,zdir) )
628  call check( __line__, mflx_hi(k ,i ,j ,xdir) )
629  call check( __line__, mflx_hi(k ,i-1,j ,xdir) )
630  call check( __line__, mflx_hi(k ,i ,j ,ydir) )
631  call check( __line__, mflx_hi(k ,i ,j-1,ydir) )
632  call check( __line__, dens_t(k,i,j) )
633 #endif
634  advcv = - ( mflx_hi(k,i,j,zdir)-mflx_hi(k-1,i, j, zdir) ) * rcdz(k)
635  advch = - ( mflx_hi(k,i,j,xdir)-mflx_hi(k ,i-1,j, xdir) ) * rcdx(i) &
636  - ( mflx_hi(k,i,j,ydir)-mflx_hi(k ,i, j-1,ydir) ) * rcdy(j)
637  dens_rk(k,i,j) = dens0(k,i,j) &
638  + dtrk * ( ( advcv + advch ) * mapf(i,j,1,i_xy) * mapf(i,j,2,i_xy) / gsqrt(k,i,j,i_xyz) &
639  + dens_t(k,i,j) )
640 #ifdef HIST_TEND
641  if ( lhist ) then
642  advcv_t(k,i,j,i_dens) = advcv * mapf(i,j,1,i_xy) * mapf(i,j,2,i_xy) / gsqrt(k,i,j,i_xyz)
643  advch_t(k,i,j,i_dens) = advch * mapf(i,j,1,i_xy) * mapf(i,j,2,i_xy) / gsqrt(k,i,j,i_xyz)
644  endif
645 #endif
646  enddo
647  enddo
648  enddo
649  !$omp end do nowait
650  end if
651 
652  !$omp end parallel
653 
654 #ifdef DEBUG
655  k = iundef; i = iundef; j = iundef
656 #endif
657 
658  !########################################################################
659  ! momentum equation (z)
660  !########################################################################
661 
662  !-----< high order flux >-----
663 
664  ! at (x, y, z)
665  ! note than z-index is added by -1
666  call atmos_dyn_fvm_fluxz_xyw( qflx_hi(:,:,:,zdir), & ! (out)
667  momz, momz, dens, & ! (in)
668  gsqrt(:,:,:,i_xyz), j33g, & ! (in)
669  num_diff(:,:,:,i_momz,zdir), & ! (in)
670  cdz, fdz, dtrk, &
671  iis, iie, jjs, jje ) ! (in)
672  if ( .not. twod ) &
673  call atmos_dyn_fvm_fluxj13_xyw( qflx_j13, & ! (out)
674  momx, momz, dens, & ! (in)
675  gsqrt(:,:,:,i_xyz), j13g(:,:,:,i_xyz), mapf(:,:,:,i_xy), & ! (in)
676  cdz, twod, &
677  iis, iie, jjs, jje ) ! (in)
678  call atmos_dyn_fvm_fluxj23_xyw( qflx_j23, & ! (out)
679  momy, momz, dens, & ! (in)
680  gsqrt(:,:,:,i_xyz), j23g(:,:,:,i_xyz), mapf(:,:,:,i_xy), & ! (in)
681  cdz, twod, &
682  iis, iie, jjs, jje ) ! (in)
683 
684  ! at (u, y, w)
685  if ( .not. twod ) &
686  call atmos_dyn_fvm_fluxx_xyw( qflx_hi(:,:,:,xdir), & ! (out)
687  momx, momz, dens, & ! (in)
688  gsqrt(:,:,:,i_uyw), mapf(:,:,:,i_uy), & ! (in)
689  num_diff(:,:,:,i_momz,xdir), & ! (in)
690  cdz, twod, & ! (in)
691  iis, iie, jjs, jje ) ! (in)
692 
693  ! at (x, v, w)
694  call atmos_dyn_fvm_fluxy_xyw( qflx_hi(:,:,:,ydir), & ! (out)
695  momy, momz, dens, & ! (in)
696  gsqrt(:,:,:,i_xvw), mapf(:,:,:,i_xv), & ! (in)
697  num_diff(:,:,:,i_momz,ydir), & ! (in)
698  cdz, twod, & ! (in)
699  iis, iie, jjs, jje ) ! (in)
700 
701 
702  !$omp parallel private(i,j,k,advcv,advch,wdamp,div)
703 
704  ! pressure gradient force at (x, y, w)
705 
706  !$omp do OMP_SCHEDULE_ collapse(2)
707  do j = jjs, jje
708  do i = iis, iie
709  do k = ks, ke-1
710  pgf(k,i,j) = j33g * ( dpres(k+1,i,j)-dpres(k,i,j) ) * rfdz(k) ! [x,y,z]
711  enddo
712  enddo
713  enddo
714  !$omp end do nowait
715 
716  ! buoyancy force at (x, y, w)
717 
718  !$omp do OMP_SCHEDULE_ collapse(2)
719  do j = jjs, jje
720  do i = iis, iie
721  do k = ks, ke-1
722  ! use not density at the half level but mean density between CZ(k) and C(k+1)
723  buoy(k,i,j) = grav * gsqrt(k,i,j,i_xyw) &
724  * 0.5_rp * ( ( dens(k+1,i,j)-ref_dens(k+1,i,j) ) &
725  + ( dens(k ,i,j)-ref_dens(k ,i,j) ) )
726  enddo
727  enddo
728  enddo
729  !$omp end do nowait
730 
731  !-----< update momentum (z) -----
732 
733  if ( twod ) then
734  !$omp do OMP_SCHEDULE_
735  do j = jjs, jje
736  do k = ks, ke-1
737 #ifdef DEBUG
738  call check( __line__, qflx_hi(k ,is,j ,zdir) )
739  call check( __line__, qflx_hi(k-1,is,j ,zdir) )
740  call check( __line__, qflx_hi(k ,is,j ,ydir) )
741  call check( __line__, qflx_hi(k ,is,j-1,ydir) )
742  call check( __line__, ddiv(k ,is,j) )
743  call check( __line__, ddiv(k+1,is,j) )
744  call check( __line__, momz0(k,is,j) )
745  call check( __line__, momz_t(k,is,j) )
746 #endif
747  advcv = - ( qflx_hi(k,is,j,zdir) - qflx_hi(k-1,is,j ,zdir) &
748  + qflx_j23(k,is,j) - qflx_j23(k-1,is,j) ) * rfdz(k)
749  advch = - ( qflx_hi(k,is,j,ydir) - qflx_hi(k,is,j-1,ydir) ) * rcdy(j) &
750  * mapf(is,j,2,i_xy)
751  wdamp = - wdamp_coef(k) * momz0(k,is,j)
752  div = divdmp_coef / dtrk * fdz(k) * ( ddiv(k+1,is,j)-ddiv(k,is,j) ) ! divergence damping
753  momz_rk(k,is,j) = momz0(k,is,j) &
754  + dtrk * ( ( advcv + advch &
755  - pgf(k,is,j) & ! pressure gradient force
756  - buoy(k,is,j) & ! buoyancy force
757  ) / gsqrt(k,is,j,i_xyw) &
758  + wdamp & ! Rayleigh damping
759  + div &
760  + momz_t(k,is,j) ) ! physics tendency
761 #ifdef HIST_TEND
762  if ( lhist ) then
763  advcv_t(k,is,j,i_momz) = advcv / gsqrt(k,is,j,i_xyw)
764  advch_t(k,is,j,i_momz) = advch / gsqrt(k,is,j,i_xyw)
765  pg_t(k,is,j,1) = ( - pgf(k,is,j) - buoy(k,is,j) ) / gsqrt(k,is,j,i_xyw)
766  wdmp_t(k,is,j) = wdamp
767  ddiv_t(k,is,j,1) = div
768  endif
769 #endif
770  enddo
771  enddo
772  !$omp end do nowait
773  else
774  !$omp do OMP_SCHEDULE_ collapse(2)
775  do j = jjs, jje
776  do i = iis, iie
777  do k = ks, ke-1
778 #ifdef DEBUG
779  call check( __line__, qflx_hi(k ,i ,j ,zdir) )
780  call check( __line__, qflx_hi(k-1,i ,j ,zdir) )
781  call check( __line__, qflx_hi(k ,i ,j ,xdir) )
782  call check( __line__, qflx_hi(k ,i-1,j ,xdir) )
783  call check( __line__, qflx_hi(k ,i ,j ,ydir) )
784  call check( __line__, qflx_hi(k ,i ,j-1,ydir) )
785  call check( __line__, ddiv(k ,i,j) )
786  call check( __line__, ddiv(k+1,i,j) )
787  call check( __line__, momz0(k,i,j) )
788  call check( __line__, momz_t(k,i,j) )
789 #endif
790  advcv = - ( qflx_hi(k,i,j,zdir) - qflx_hi(k-1,i ,j ,zdir) &
791  + qflx_j13(k,i,j) - qflx_j13(k-1,i,j) &
792  + qflx_j23(k,i,j) - qflx_j23(k-1,i,j) ) * rfdz(k)
793  advch = - ( ( qflx_hi(k,i,j,xdir) - qflx_hi(k,i-1,j,xdir) ) * rcdx(i) &
794  + ( qflx_hi(k,i,j,ydir) - qflx_hi(k,i,j-1,ydir) ) * rcdy(j) ) &
795  * mapf(i,j,1,i_xy) * mapf(i,j,2,i_xy)
796  wdamp = - wdamp_coef(k) * momz0(k,i,j)
797  div = divdmp_coef / dtrk * fdz(k) * ( ddiv(k+1,i,j)-ddiv(k,i,j) ) ! divergence damping
798  momz_rk(k,i,j) = momz0(k,i,j) &
799  + dtrk * ( ( advcv + advch &
800  - pgf(k,i,j) & ! pressure gradient force
801  - buoy(k,i,j) & ! buoyancy force
802  ) / gsqrt(k,i,j,i_xyw) &
803  + wdamp & ! Rayleigh damping
804  + div &
805  + momz_t(k,i,j) ) ! physics tendency
806 #ifdef HIST_TEND
807  if ( lhist ) then
808  advcv_t(k,i,j,i_momz) = advcv / gsqrt(k,i,j,i_xyw)
809  advch_t(k,i,j,i_momz) = advch / gsqrt(k,i,j,i_xyw)
810  pg_t(k,i,j,1) = ( - pgf(k,i,j) - buoy(k,i,j) ) / gsqrt(k,i,j,i_xyw)
811  wdmp_t(k,i,j) = wdamp
812  ddiv_t(k,i,j,1) = div
813  endif
814 #endif
815  enddo
816  enddo
817  enddo
818  !$omp end do nowait
819  end if
820 #ifdef DEBUG
821  k = iundef; i = iundef; j = iundef
822 #endif
823 
824  !$omp do OMP_SCHEDULE_ collapse(2)
825  do j = jjs, jje
826  do i = iis, iie
827  momz_rk(ks-1,i,j) = 0.0_rp
828  momz_rk(ke ,i,j) = 0.0_rp
829 #ifdef HIST_TEND
830  if ( lhist ) then
831  advcv_t(ke,i,j,i_momz) = 0.0_rp
832  advch_t(ke,i,j,i_momz) = 0.0_rp
833  pg_t(ke,i,j,1) = 0.0_rp
834  wdmp_t(ke,i,j) = 0.0_rp
835  ddiv_t(ke,i,j,1) = 0.0_rp
836  endif
837 #endif
838  enddo
839  enddo
840  !$omp end do nowait
841 
842  !$omp end parallel
843 #ifdef DEBUG
844  k = iundef; i = iundef; j = iundef
845 #endif
846 
847 #ifndef NO_FCT_DYN
848  if ( flag_fct_momentum ) then
849 
850  ! monotonic flux
851  ! at (x, y, layer)
852  ! note than z-index is added by -1
853  call atmos_dyn_fvm_fluxz_xyw_ud1( qflx_lo(:,:,:,zdir), & ! (out)
854  momz, momz0, dens, & ! (in)
855  gsqrt(:,:,:,i_xyz), j33g, & ! (in)
856  num_diff(:,:,:,i_momz,zdir), & ! (in)
857  cdz, fdz, dtrk, & ! (in)
858  iis-1, iie+1, jjs-1, jje+1 ) ! (in)
859 
860  if ( .not. twod ) &
861  call atmos_dyn_fvm_fluxx_xyw_ud1( qflx_lo(:,:,:,xdir), & ! (out)
862  momx, momz0, dens, & ! (in)
863  gsqrt(:,:,:,i_uyz), mapf(:,:,:,i_uy), & ! (in)
864  num_diff(:,:,:,i_momz,xdir), & ! (in)
865  cdz, twod, & ! (in)
866  iis-1, iie+1, jjs-1, jje+1 ) ! (in)
867 
868  call atmos_dyn_fvm_fluxy_xyw_ud1( qflx_lo(:,:,:,ydir), & ! (out)
869  momy, momz0, dens, & ! (in)
870  gsqrt(:,:,:,i_xvz), mapf(:,:,:,i_xv), & ! (in)
871  num_diff(:,:,:,i_momz,ydir), & ! (in)
872  cdz, twod, & ! (in)
873  iis-1, iie+1, jjs-1, jje+1 ) ! (in)
874 
875  endif
876  enddo
877  enddo
878 
879  if ( flag_fct_momentum ) then
880 
881  call comm_vars8( dens_rk, 1 )
882  call comm_wait ( dens_rk, 1, .false. )
883 
884  do j = js, je
885  do i = is, ie
886  do k = ks, ke
887  qflx_hi(k,i,j,zdir) = qflx_hi(k,i,j,zdir) / ( mapf(i,j,1,i_xy) * mapf(i,j,2,i_xy) ) &
888  + qflx_j13(k,i,j) + qflx_j23(k,i,j)
889  enddo
890  enddo
891  enddo
892 
893  do j = js-1, je+1
894  do i = is-1, ie+1
895  do k = ks, ke-1
896  dens0_uvw(k,i,j) = 0.5_rp * ( dens0(k,i,j) + dens0(k+1,i,j) )
897  dens_uvw(k,i,j) = 0.5_rp * ( dens_rk(k,i,j) + dens_rk(k+1,i,j) )
898  enddo
899  enddo
900  enddo
901 
902  do j = js-1, je+1
903  do i = is-1, ie+1
904  dens_uvw(ke,i,j) = dens_uvw(ke-1,i,j)
905  dens0_uvw(ke,i,j) = dens0_uvw(ke-1,i,j)
906  enddo
907  enddo
908 
909  call comm_wait ( velz(:,:,:), 4 )
910 
911  call atmos_dyn_fvm_fct( qflx_anti, & ! (out)
912  velz, dens0_uvw, dens_uvw, & ! (in)
913  qflx_hi, qflx_lo, & ! (in)
914  mflx_hi, & ! (in)
915  rfdz, rcdx, rcdy, & ! (in)
916  gsqrt(:,:,:,i_xyw), & ! (in)
917  mapf(:,:,:,i_xy), & ! (in)
918  twod, dtrk, & ! (in)
919  flag_fct_along_stream ) ! (in)
920 
921  do jjs = js, je, jblock
922  jje = jjs+jblock-1
923  do iis = is, ie, iblock
924  iie = iis+iblock-1
925 
926  !--- update momentum(z)
927  if ( twod ) then
928  !$omp parallel do private(i,j,k) OMP_SCHEDULE_
929  do j = jjs, jje
930  do k = ks, ke-1
931  momz_rk(k,is,j) = momz_rk(k,is,j) &
932  + dtrk * ( ( qflx_anti(k,is,j,zdir) - qflx_anti(k-1,is,j ,zdir) ) * rfdz(k) &
933  + ( qflx_anti(k,is,j,ydir) - qflx_anti(k ,is,j-1,ydir) ) * rcdy(j) &
934  * mapf(is,j,2,i_xy) ) &
935  / gsqrt(k,is,j,i_xyw)
936  enddo
937  enddo
938  else
939  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
940  do j = jjs, jje
941  do i = iis, iie
942  do k = ks, ke-1
943  momz_rk(k,i,j) = momz_rk(k,i,j) &
944  + dtrk * ( ( qflx_anti(k,i,j,zdir) - qflx_anti(k-1,i ,j ,zdir) ) * rfdz(k) &
945  + ( ( qflx_anti(k,i,j,xdir) - qflx_anti(k ,i-1,j ,xdir) ) * rcdx(i) &
946  + ( qflx_anti(k,i,j,ydir) - qflx_anti(k ,i ,j-1,ydir) ) * rcdy(j) ) &
947  * mapf(i,j,1,i_xy) * mapf(i,j,2,i_xy) ) &
948  / gsqrt(k,i,j,i_xyw)
949  enddo
950  enddo
951  enddo
952  end if
953 
954  enddo
955  enddo
956 
957  endif ! FLAG_FCT_MOMENTUM
958 
959 #ifdef DEBUG
960  qflx_hi(:,:,:,:) = undef
961 #endif
962 
963  do jjs = js, je, jblock
964  jje = jjs+jblock-1
965  do iis = is, ie, iblock
966  iie = iis+iblock-1
967 #endif
968 
969  !########################################################################
970  ! momentum equation (x)
971  !########################################################################
972 
973  !-----< high order flux >-----
974 
975  ! at (u, y, w)
976  call atmos_dyn_fvm_fluxz_uyz( qflx_hi(:,:,:,zdir), & ! (out)
977  momz, momx, dens, & ! (in)
978  gsqrt(:,:,:,i_uyw), j33g, & ! (in)
979  num_diff(:,:,:,i_momx,zdir), & ! (in)
980  cdz, twod, & ! (in)
981  iis, iie, jjs, jje ) ! (in)
982  if ( .not. twod ) &
983  call atmos_dyn_fvm_fluxj13_uyz( qflx_j13, & ! (out)
984  momx, momx, dens, & ! (in)
985  gsqrt(:,:,:,i_uyz), j13g(:,:,:,i_uyw), mapf(:,:,:,i_uy), & ! (in)
986  cdz, twod, & ! (in)
987  iis, iie, jjs, jje ) ! (in)
988  call atmos_dyn_fvm_fluxj23_uyz( qflx_j23, & ! (out)
989  momy, momx, dens, & ! (in)
990  gsqrt(:,:,:,i_uyz), j23g(:,:,:,i_uyw), mapf(:,:,:,i_uy), & ! (in)
991  cdz, twod, & ! (in)
992  iis, iie, jjs, jje ) ! (in)
993 
994  ! at (x, y, z)
995  ! note that x-index is added by -1
996  if ( .not. twod ) &
997  call atmos_dyn_fvm_fluxx_uyz( qflx_hi(:,:,:,xdir), & ! (out)
998  momx, momx, dens, & ! (in)
999  gsqrt(:,:,:,i_xyz), mapf(:,:,:,i_xy), & ! (in)
1000  num_diff(:,:,:,i_momx,xdir), & ! (in)
1001  cdz, twod, & ! (in)
1002  iis, iie, jjs, jje ) ! (in)
1003 
1004  ! at (u, v, z)
1005  call atmos_dyn_fvm_fluxy_uyz( qflx_hi(:,:,:,ydir), & ! (out)
1006  momy, momx, dens, & ! (in)
1007  gsqrt(:,:,:,i_uvz), mapf(:,:,1,i_uv), & ! (in)
1008  num_diff(:,:,:,i_momx,ydir), & ! (in)
1009  cdz, twod, & ! (in)
1010  iis, iie, jjs, jje ) ! (in)
1011 
1012 
1013  !$omp parallel private(i,j,k,advcv,advch,div)
1014 
1015  ! pressure gradient force at (u, y, z)
1016 
1017  if ( .not. twod ) then
1018  !$omp do OMP_SCHEDULE_ collapse(2)
1019  do j = jjs, jje
1020  do i = iis, iie
1021  do k = ks, ke
1022  pgf(k,i,j) = ( ( gsqrt(k,i+1,j,i_xyz) * dpres(k,i+1,j) & ! [x,y,z]
1023  - gsqrt(k,i ,j,i_xyz) * dpres(k,i ,j) & ! [x,y,z]
1024  ) * rfdx(i) &
1025  + ( j13g(k ,i,j,i_uyw) &
1026  * 0.5_rp * ( f2h(k,1,i_uyz) * ( dpres(k+1,i+1,j)+dpres(k+1,i,j) ) &
1027  + f2h(k,2,i_uyz) * ( dpres(k ,i+1,j)+dpres(k ,i,j) ) ) & ! [x,y,z->u,y,w]
1028  - j13g(k-1,i,j,i_uyw) &
1029  * 0.5_rp * ( f2h(k,1,i_uyz) * ( dpres(k ,i+1,j)+dpres(k ,i,j) ) &
1030  + f2h(k,2,i_uyz) * ( dpres(k-1,i+1,j)+dpres(k-1,i,j) ) ) & ! [x,y,z->u,y,w]
1031  ) * rcdz(k) ) &
1032  * mapf(i,j,1,i_uy)
1033  enddo
1034  enddo
1035  enddo
1036  !$omp end do nowait
1037  end if
1038 
1039  ! coriolis force at (u, y, z)
1040 
1041  if ( twod ) then
1042  !$omp do OMP_SCHEDULE_
1043  do j = jjs, jje
1044  do k = ks, ke
1045 #ifdef DEBUG
1046  call check( __line__, momy(k,is,j ) )
1047  call check( __line__, momy(k,is,j-1) )
1048 #endif
1049  cor(k,is,j) = 0.5_rp * corioli(is,j) * ( momy(k,is,j) + momy(k,is,j-1) )
1050  enddo
1051  enddo
1052  !$omp end do
1053  else
1054  !$omp do OMP_SCHEDULE_ collapse(2)
1055  do j = jjs, jje
1056  do i = iis, iie
1057  do k = ks, ke
1058 #ifdef DEBUG
1059  call check( __line__, momy(k,i ,j ) )
1060  call check( __line__, momy(k,i+1,j ) )
1061  call check( __line__, momy(k,i ,j-1) )
1062  call check( __line__, momy(k,i+1,j-1) )
1063 #endif
1064  cor(k,i,j) = 0.125_rp * ( corioli( i+1,j )+corioli( i,j ) ) & ! [x,y,z->u,y,z]
1065  * ( momy(k,i+1,j )+momy(k,i,j ) &
1066  + momy(k,i+1,j-1)+momy(k,i,j-1) ) & ! [x,v,z->u,y,z]
1067  + 0.25_rp * mapf(i,j,1,i_uy) * mapf(i,j,2,i_uy) &
1068  * ( momy(k,i,j) + momy(k,i,j-1) + momy(k,i+1,j) + momy(k,i+1,j-1) ) &
1069  * ( ( momy(k,i,j) + momy(k,i,j-1) + momy(k,i+1,j) + momy(k,i+1,j-1) ) * 0.25_rp &
1070  * ( 1.0_rp/mapf(i+1,j,2,i_xy) - 1.0_rp/mapf(i,j,2,i_xy) ) * rfdx(i) &
1071  - momx(k,i,j) &
1072  * ( 1.0_rp/mapf(i,j,1,i_uv) - 1.0_rp/mapf(i,j-1,1,i_uv) ) * rcdy(j) ) &
1073  * 2.0_rp / ( dens(k,i+1,j) + dens(k,i,j) ) ! metric term
1074  enddo
1075  enddo
1076  enddo
1077  !$omp end do
1078  end if
1079 
1080  !-----< update momentum (x) >-----
1081 
1082  if ( twod ) then
1083  !$omp do OMP_SCHEDULE_
1084  do j = jjs, jje
1085  do k = ks, ke
1086 #ifdef DEBUG
1087  call check( __line__, qflx_hi(k ,is,j ,zdir) )
1088  call check( __line__, qflx_hi(k-1,is,j ,zdir) )
1089  call check( __line__, qflx_hi(k ,is,j ,ydir) )
1090  call check( __line__, qflx_hi(k ,is,j-1,ydir) )
1091  call check( __line__, momx0(k,is,j) )
1092 #endif
1093  ! advection
1094  advcv = - ( qflx_hi(k,is,j,zdir) - qflx_hi(k-1,is,j ,zdir) &
1095  + qflx_j23(k,is,j) - qflx_j23(k-1,is,j) ) * rcdz(k)
1096  advch = - ( qflx_hi(k,is,j,ydir) - qflx_hi(k ,is,j-1,ydir) ) * rcdy(j) &
1097  * mapf(is,j,2,i_uy)
1098  momx_rk(k,is,j) = momx0(k,is,j) &
1099  + dtrk * ( ( advcv + advch & ! advection
1100  ) / gsqrt(k,is,j,i_uyz) &
1101  + cor(k,is,j) & ! coriolis force
1102  + div & ! divergence damping
1103  + momx_t(k,is,j) ) ! physics tendency
1104 #ifdef HIST_TEND
1105  if ( lhist ) then
1106  advcv_t(k,is,j,i_momx) = advcv / gsqrt(k,is,j,i_uyz)
1107  advch_t(k,is,j,i_momx) = advch / gsqrt(k,is,j,i_uyz)
1108  pg_t(k,is,j,2) = 0.0_rp
1109  cf_t(k,is,j,1) = cor(k,i,j)
1110  ddiv_t(k,is,j,2) = 0.0_rp
1111  endif
1112 #endif
1113  enddo
1114  enddo
1115  !$omp end do nowait
1116  else
1117  !$omp do OMP_SCHEDULE_ collapse(2)
1118  do j = jjs, jje
1119  do i = iis, min(iie, ieh)
1120  do k = ks, ke
1121 #ifdef DEBUG
1122  call check( __line__, qflx_hi(k ,i ,j ,zdir) )
1123  call check( __line__, qflx_hi(k-1,i ,j ,zdir) )
1124  call check( __line__, qflx_hi(k ,i ,j ,xdir) )
1125  call check( __line__, qflx_hi(k ,i-1,j ,xdir) )
1126  call check( __line__, qflx_hi(k ,i ,j ,ydir) )
1127  call check( __line__, qflx_hi(k ,i ,j-1,ydir) )
1128  call check( __line__, ddiv(k,i+1,j) )
1129  call check( __line__, ddiv(k,i ,j) )
1130  call check( __line__, momx0(k,i,j) )
1131 #endif
1132  ! advection
1133  advcv = - ( qflx_hi(k,i,j,zdir) - qflx_hi(k-1,i ,j ,zdir) &
1134  + qflx_j13(k,i,j) - qflx_j13(k-1,i,j) &
1135  + qflx_j23(k,i,j) - qflx_j23(k-1,i,j) ) * rcdz(k)
1136  advch = - ( ( qflx_hi(k,i,j,xdir) - qflx_hi(k ,i-1,j ,xdir) ) * rfdx(i) &
1137  + ( qflx_hi(k,i,j,ydir) - qflx_hi(k ,i ,j-1,ydir) ) * rcdy(j) ) &
1138  * mapf(i,j,1,i_uy) * mapf(i,j,2,i_uy)
1139  div = divdmp_coef / dtrk * fdx(i) * ( ddiv(k,i+1,j)-ddiv(k,i,j) ) ! divergence damping
1140  momx_rk(k,i,j) = momx0(k,i,j) &
1141  + dtrk * ( ( advcv + advch & ! advection
1142  - pgf(k,i,j) & ! pressure gradient force
1143  ) / gsqrt(k,i,j,i_uyz) &
1144  + cor(k,i,j) & ! coriolis force
1145  + div & ! divergence damping
1146  + momx_t(k,i,j) ) ! physics tendency
1147 #ifdef HIST_TEND
1148  if ( lhist ) then
1149  advcv_t(k,i,j,i_momx) = advcv / gsqrt(k,i,j,i_uyz)
1150  advch_t(k,i,j,i_momx) = advch / gsqrt(k,i,j,i_uyz)
1151  pg_t(k,i,j,2) = - pgf(k,i,j) / gsqrt(k,i,j,i_uyz)
1152  cf_t(k,i,j,1) = cor(k,i,j)
1153  ddiv_t(k,i,j,2) = div
1154  endif
1155 #endif
1156  enddo
1157  enddo
1158  enddo
1159  !$omp end do nowait
1160  end if
1161 
1162  !$omp end parallel
1163 #ifdef DEBUG
1164  k = iundef; i = iundef; j = iundef
1165 #endif
1166 
1167 #ifndef NO_FCT_DYN
1168  if ( flag_fct_momentum ) then
1169 
1170  call atmos_dyn_fvm_fluxz_uyz_ud1( qflx_lo(:,:,:,zdir), & ! (out)
1171  momz, momx0, dens, & ! (in)
1172  gsqrt(:,:,:,i_uyw), j33g, & ! (in)
1173  num_diff(:,:,:,i_momx,zdir), & ! (in)
1174  cdz, twod, & ! (in)
1175  iis-1, iie+1, jjs-1, jje+1 ) ! (in)
1176 
1177  ! note that x-index is added by -1
1178  if ( .not. twod ) &
1179  call atmos_dyn_fvm_fluxx_uyz_ud1( qflx_lo(:,:,:,xdir), & ! (out)
1180  momx, momx0, dens, & ! (in)
1181  gsqrt(:,:,:,i_xyz), mapf(:,:,:,i_uy), & ! (in)
1182  num_diff(:,:,:,i_momx,xdir), & ! (in)
1183  cdz, twod, & ! (in)
1184  iis-1, iie+1, jjs-1, jje+1 ) ! (in)
1185 
1186  call atmos_dyn_fvm_fluxy_uyz_ud1( qflx_lo(:,:,:,ydir), & ! (out)
1187  momy, momx0, dens, & ! (in)
1188  gsqrt(:,:,:,i_uvz), mapf(:,:,:,i_xv), & ! (in)
1189  num_diff(:,:,:,i_momx,ydir), & ! (in)
1190  cdz, twod, & ! (in)
1191  iis-1, iie+1, jjs-1, jje+1 ) ! (in)
1192 
1193  endif
1194  enddo
1195  enddo
1196 
1197  if ( flag_fct_momentum ) then
1198 
1199  !$omp parallel do
1200  do j = js, je
1201  do i = is, ie
1202  do k = ks, ke
1203  qflx_hi(k,i,j,zdir) = qflx_hi(k,i,j,zdir) / ( mapf(i,j,1,i_uy) * mapf(i,j,2,i_uy) )&
1204  + qflx_j13(k,i,j) + qflx_j23(k,i,j)
1205  enddo
1206  enddo
1207  enddo
1208 
1209  if ( twod ) then
1210  !$omp parallel do
1211  do j = js-1, je+1
1212  do k = ks, ke
1213  dens0_uvw(k,is,j) = dens0(k,is,j)
1214  dens_uvw(k,is,j) = dens_rk(k,is,j)
1215  enddo
1216  enddo
1217  else
1218  !$omp parallel do
1219  do j = js-1, je+1
1220  do i = is-1, ie+1
1221  do k = ks, ke
1222  dens0_uvw(k,i,j) = 0.5_rp * ( dens0(k,i,j) + dens0(k,i+1,j) )
1223  dens_uvw(k,i,j) = 0.5_rp * ( dens_rk(k,i,j) + dens_rk(k,i+1,j) )
1224  enddo
1225  enddo
1226  enddo
1227  end if
1228 
1229  call comm_wait ( velx(:,:,:), 5 )
1230 
1231  call atmos_dyn_fvm_fct( qflx_anti, & ! (out)
1232  velx, dens0_uvw, dens_uvw, & ! (in)
1233  qflx_hi, qflx_lo, & ! (in)
1234  mflx_hi, & ! (in)
1235  rcdz, rfdx, rcdy, & ! (in)
1236  gsqrt(:,:,:,i_uyz), & ! (in)
1237  mapf(:,:,:,i_uy), & ! (in)
1238  twod, dtrk, & ! (in)
1239  flag_fct_along_stream ) ! (in)
1240 
1241  do jjs = js, je, jblock
1242  jje = jjs+jblock-1
1243  do iis = is, ie, iblock
1244  iie = iis+iblock-1
1245 
1246  !--- update momentum(x)
1247  if ( twod ) then
1248  !$omp parallel do private(j,k) OMP_SCHEDULE_
1249  do j = jjs, jje
1250  do k = ks, ke
1251 #ifdef DEBUG
1252  call check( __line__, momx_rk(k,is,j) )
1253  call check( __line__, qflx_anti(k ,is,j ,zdir) )
1254  call check( __line__, qflx_anti(k-1,is,j ,zdir) )
1255  call check( __line__, qflx_anti(k ,is,j ,ydir) )
1256  call check( __line__, qflx_anti(k ,is,j-1,ydir) )
1257 #endif
1258  momx_rk(k,is,j) = momx_rk(k,is,j) &
1259  + dtrk * ( ( ( qflx_anti(k,is,j,zdir) - qflx_anti(k-1,is,j ,zdir) ) * rcdz(k) &
1260  + ( qflx_anti(k,is,j,ydir) - qflx_anti(k ,is,j-1,ydir) ) * rcdy(j) ) ) &
1261  * mapf(is,j,2,i_uy) &
1262  / gsqrt(k,is,j,i_uyz)
1263  enddo
1264  enddo
1265  else
1266  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
1267  do j = jjs, jje
1268  do i = iis, min(iie,ieh)
1269  do k = ks, ke
1270 #ifdef DEBUG
1271  call check( __line__, momx_rk(k,i,j) )
1272  call check( __line__, qflx_anti(k ,i ,j ,zdir) )
1273  call check( __line__, qflx_anti(k-1,i ,j ,zdir) )
1274  call check( __line__, qflx_anti(k ,i ,j ,xdir) )
1275  call check( __line__, qflx_anti(k ,i-1,j ,xdir) )
1276  call check( __line__, qflx_anti(k ,i ,j ,ydir) )
1277  call check( __line__, qflx_anti(k ,i ,j-1,ydir) )
1278 #endif
1279  momx_rk(k,i,j) = momx_rk(k,i,j) &
1280  + dtrk * ( ( ( qflx_anti(k,i,j,zdir) - qflx_anti(k-1,i ,j ,zdir) ) * rcdz(k) &
1281  + ( qflx_anti(k,i,j,xdir) - qflx_anti(k ,i-1,j ,xdir) ) * rfdx(i) &
1282  + ( qflx_anti(k,i,j,ydir) - qflx_anti(k ,i ,j-1,ydir) ) * rcdy(j) ) ) &
1283  * mapf(i,j,1,i_uy) * mapf(i,j,2,i_uy) &
1284  / gsqrt(k,i,j,i_uyz)
1285  enddo
1286  enddo
1287  enddo
1288  end if
1289 #ifdef DEBUG
1290  k = iundef; i = iundef; j = iundef
1291 #endif
1292 
1293  enddo
1294  enddo
1295 
1296 #ifdef DEBUG
1297  qflx_lo(:,:,:,:) = undef
1298  qflx_anti(:,:,:,:) = undef
1299 #endif
1300 
1301  endif ! FLAG_FCT_MOMENTUM
1302 
1303 #ifdef DEBUG
1304  qflx_hi(:,:,:,:) = undef
1305 #endif
1306 
1307  do jjs = js, je, jblock
1308  jje = jjs+jblock-1
1309  do iis = is, ie, iblock
1310  iie = iis+iblock-1
1311 #endif
1312 
1313  !########################################################################
1314  ! momentum equation (y)
1315  !########################################################################
1316 
1317  !-----< high order flux >-----
1318 
1319  ! at (x, v, w)
1320  call atmos_dyn_fvm_fluxz_xvz( qflx_hi(:,:,:,zdir), & ! (out)
1321  momz, momy, dens, & ! (in)
1322  gsqrt(:,:,:,i_xvw), j33g, & ! (in)
1323  num_diff(:,:,:,i_momy,zdir), & ! (in)
1324  cdz, twod, & ! (in)
1325  iis, iie, jjs, jje ) ! (in)
1326  if ( .not. twod ) &
1327  call atmos_dyn_fvm_fluxj13_xvz( qflx_j13, & ! (out)
1328  momx, momy, dens, & ! (in)
1329  gsqrt(:,:,:,i_xvz), j13g(:,:,:,i_xvw), mapf(:,:,:,i_xv), & ! (in)
1330  cdz, twod, & ! (in)
1331  iis, iie, jjs, jje ) ! (in)
1332  call atmos_dyn_fvm_fluxj23_xvz( qflx_j23, & ! (out)
1333  momy, momy, dens, & ! (in)
1334  gsqrt(:,:,:,i_xvz), j23g(:,:,:,i_xvw), mapf(:,:,:,i_xv), & ! (in)
1335  cdz, twod, & ! (in)
1336  iis, iie, jjs, jje ) ! (in)
1337 
1338  ! at (u, v, z)
1339  if ( .not. twod ) &
1340  call atmos_dyn_fvm_fluxx_xvz( qflx_hi(:,:,:,xdir), & ! (out)
1341  momx, momy, dens, & ! (in)
1342  gsqrt(:,:,:,i_uvz), mapf(:,:,:,i_uv), & ! (in)
1343  num_diff(:,:,:,i_momy,xdir), & ! (in)
1344  cdz, twod, & ! (in)
1345  iis, iie, jjs, jje ) ! (in)
1346 
1347  ! at (x, y, z)
1348  ! note that y-index is added by -1
1349  call atmos_dyn_fvm_fluxy_xvz( qflx_hi(:,:,:,ydir), & ! (out)
1350  momy, momy, dens, & ! (in)
1351  gsqrt(:,:,:,i_xyz), mapf(:,:,:,i_xy), & ! (in)
1352  num_diff(:,:,:,i_momy,ydir), & ! (in
1353  cdz, twod, & ! (in)
1354  iis, iie, jjs, jje ) ! (in)
1355 
1356 
1357  !$omp parallel private(i,j,k,advcv,advch,div)
1358 
1359  ! pressure gradient force at (x, v, z)
1360 
1361  !$omp do OMP_SCHEDULE_ collapse(2)
1362  do j = jjs, jje
1363  do i = iis, iie
1364  do k = ks, ke
1365  pgf(k,i,j) = ( ( gsqrt(k,i,j+1,i_xyz) * dpres(k,i,j+1) & ! [x,y,z]
1366  - gsqrt(k,i,j ,i_xyz) * dpres(k,i,j ) & ! [x,y,z]
1367  ) * rfdy(j) &
1368  + ( j23g(k ,i,j,i_xvw) &
1369  * 0.5_rp * ( f2h(k ,1,i_xvz) * ( dpres(k+1,i,j+1)+dpres(k+1,i,j) ) &
1370  + f2h(k ,2,i_xvz) * ( dpres(k ,i,j+1)+dpres(k ,i,j) ) ) & ! [x,y,z->x,v,w]
1371  - j23g(k-1,i,j,i_xvw) &
1372  * 0.5_rp * ( f2h(k-1,1,i_xvz) * ( dpres(k ,i,j+1)+dpres(k ,i,j) ) &
1373  + f2h(k-1,2,i_xvz) * ( dpres(k-1,i,j+1)+dpres(k-1,i,j) ) ) & ! [x,y,z->x,v,w]
1374  ) * rcdz(k) ) &
1375  * mapf(i,j,2,i_xv)
1376  enddo
1377  enddo
1378  enddo
1379  !$omp end do nowait
1380 
1381  ! coriolis force at (x, v, z)
1382 
1383  if ( twod ) then
1384  !$omp do OMP_SCHEDULE_
1385  do j = jjs, jje
1386  do k = ks, ke
1387 #ifdef DEBUG
1388  call check( __line__, momx(k,is,j ) )
1389  call check( __line__, momx(k,is,j+1) )
1390 #endif
1391  cor(k,is,j) = - 0.25_rp * ( corioli( is,j+1)+corioli( is,j) ) & ! [x,y,z->x,v,z]
1392  * ( momx(k,is,j+1)+momx(k,is,j) )
1393  enddo
1394  enddo
1395  !$omp end do
1396  else
1397  !$omp do OMP_SCHEDULE_ collapse(2)
1398  do j = jjs, jje
1399  do i = iis, iie
1400  do k = ks, ke
1401 #ifdef DEBUG
1402  call check( __line__, momx(k,i ,j ) )
1403  call check( __line__, momx(k,i ,j+1) )
1404  call check( __line__, momx(k,i-1,j ) )
1405  call check( __line__, momx(k,i-1,j+1) )
1406 #endif
1407  cor(k,i,j) = - 0.125_rp * ( corioli( i ,j+1)+corioli( i ,j) ) & ! [x,y,z->x,v,z]
1408  * ( momx(k,i ,j+1)+momx(k,i ,j) &
1409  + momx(k,i-1,j+1)+momx(k,i-1,j) ) & ! [u,y,z->x,v,z]
1410  - 0.25_rp * mapf(i,j,1,i_xv) * mapf(i,j,2,i_xv) &
1411  * ( momx(k,i,j) + momx(k,i-1,j) + momx(k,i,j+1) + momx(k,i-1,j+1) )&
1412  * ( momy(k,i,j) &
1413  * ( 1.0_rp/mapf(i,j,2,i_uv) - 1.0_rp/mapf(i-1,j,2,i_uv) ) * rcdx(i) &
1414  - 0.25_rp * ( momx(k,i,j)+momx(k,i-1,j)+momx(k,i,j+1)+momx(k,i-1,j+1) ) &
1415  * ( 1.0_rp/mapf(i,j+1,1,i_xy) - 1.0_rp/mapf(i,j,1,i_xy) ) * rfdy(j) ) &
1416  * 2.0_rp / ( dens(k,i,j) + dens(k,i,j+1) ) ! metoric term
1417  enddo
1418  enddo
1419  enddo
1420  !$omp end do
1421  end if
1422 
1423  !-----< update momentum (y) >-----
1424 
1425  if ( twod ) then
1426  !$omp do OMP_SCHEDULE_
1427  do j = jjs, min(jje, jeh)
1428  do k = ks, ke
1429 #ifdef DEBUG
1430  call check( __line__, qflx_hi(k ,is,j ,zdir) )
1431  call check( __line__, qflx_hi(k-1,is,j ,zdir) )
1432  call check( __line__, qflx_hi(k ,is,j ,ydir) )
1433  call check( __line__, qflx_hi(k ,is,j-1,ydir) )
1434  call check( __line__, ddiv(k,is,j+1) )
1435  call check( __line__, ddiv(k,is,j ) )
1436  call check( __line__, momy_t(k,is,j) )
1437  call check( __line__, momy0(k,is,j) )
1438 #endif
1439  advcv = - ( qflx_hi(k,is,j,zdir) - qflx_hi(k-1,is,j ,zdir) &
1440  + qflx_j23(k,is,j) - qflx_j23(k-1,is,j) ) * rcdz(k)
1441  advch = - ( qflx_hi(k,is,j,ydir) - qflx_hi(k ,is,j-1,ydir) ) * rfdy(j) &
1442  * mapf(is,j,2,i_xv)
1443  div = divdmp_coef / dtrk * fdy(j) * ( ddiv(k,is,j+1)-ddiv(k,is,j) )
1444  momy_rk(k,is,j) = momy0(k,is,j) &
1445  + dtrk * ( ( advcv + advch & ! advection
1446  - pgf(k,is,j) & ! pressure gradient force
1447  ) / gsqrt(k,is,j,i_xvz) &
1448  + cor(k,is,j) & ! coriolis force
1449  + div & ! divergence damping
1450  + momy_t(k,is,j) ) ! physics tendency
1451 #ifdef HIST_TEND
1452  if ( lhist ) then
1453  advcv_t(k,is,j,i_momy) = advcv / gsqrt(k,is,j,i_uyz)
1454  advch_t(k,is,j,i_momy) = advch / gsqrt(k,is,j,i_uyz)
1455  pg_t(k,is,j,3) = - pgf(k,i,j) / gsqrt(k,is,j,i_uyz)
1456  cf_t(k,is,j,2) = cor(k,is,j)
1457  ddiv_t(k,is,j,3) = div
1458  endif
1459 #endif
1460  enddo
1461  enddo
1462  !$omp end do nowait
1463  else
1464  !$omp do OMP_SCHEDULE_ collapse(2)
1465  do j = jjs, min(jje, jeh)
1466  do i = iis, iie
1467  do k = ks, ke
1468 #ifdef DEBUG
1469  call check( __line__, qflx_hi(k ,i ,j ,zdir) )
1470  call check( __line__, qflx_hi(k-1,i ,j ,zdir) )
1471  call check( __line__, qflx_hi(k ,i ,j ,xdir) )
1472  call check( __line__, qflx_hi(k ,i-1,j ,xdir) )
1473  call check( __line__, qflx_hi(k ,i ,j ,ydir) )
1474  call check( __line__, qflx_hi(k ,i ,j-1,ydir) )
1475  call check( __line__, ddiv(k,i,j+1) )
1476  call check( __line__, ddiv(k,i,j ) )
1477  call check( __line__, momy_t(k,i,j) )
1478  call check( __line__, momy0(k,i,j) )
1479 #endif
1480  advcv = - ( qflx_hi(k,i,j,zdir) - qflx_hi(k-1,i ,j ,zdir) &
1481  + qflx_j13(k,i,j) - qflx_j13(k-1,i,j) &
1482  + qflx_j23(k,i,j) - qflx_j23(k-1,i,j) ) * rcdz(k)
1483  advch = - ( ( qflx_hi(k,i,j,xdir) - qflx_hi(k ,i-1,j ,xdir) ) * rcdx(i) &
1484  + ( qflx_hi(k,i,j,ydir) - qflx_hi(k ,i ,j-1,ydir) ) * rfdy(j) ) &
1485  * mapf(i,j,1,i_xv) * mapf(i,j,2,i_xv)
1486  div = divdmp_coef / dtrk * fdy(j) * ( ddiv(k,i,j+1)-ddiv(k,i,j) )
1487  momy_rk(k,i,j) = momy0(k,i,j) &
1488  + dtrk * ( ( advcv + advch & ! advection
1489  - pgf(k,i,j) & ! pressure gradient force
1490  ) / gsqrt(k,i,j,i_xvz) &
1491  + cor(k,i,j) & ! coriolis force
1492  + div & ! divergence damping
1493  + momy_t(k,i,j) ) ! physics tendency
1494 #ifdef HIST_TEND
1495  if ( lhist ) then
1496  advcv_t(k,i,j,i_momy) = advcv / gsqrt(k,i,j,i_uyz)
1497  advch_t(k,i,j,i_momy) = advch / gsqrt(k,i,j,i_uyz)
1498  pg_t(k,i,j,3) = - pgf(k,i,j) / gsqrt(k,i,j,i_uyz)
1499  cf_t(k,i,j,2) = cor(k,i,j)
1500  ddiv_t(k,i,j,3) = div
1501  endif
1502 #endif
1503  enddo
1504  enddo
1505  enddo
1506  !$omp end do nowait
1507  end if
1508 
1509  !$omp end parallel
1510 #ifdef DEBUG
1511  k = iundef; i = iundef; j = iundef
1512 #endif
1513 
1514 #ifndef NO_FCT_DYN
1515  if ( flag_fct_momentum ) then
1516 
1517  ! monotonic flux
1518  ! at (x, v, interface)
1519  call atmos_dyn_fvm_fluxz_xvz_ud1( qflx_lo(:,:,:,zdir), & ! (out)
1520  momz, momy0, dens, & ! (in)
1521  gsqrt(:,:,:,i_xvz), j33g, & ! (in)
1522  num_diff(:,:,:,i_momy,zdir), & ! (in)
1523  cdz, twod, & ! (in)
1524  iis-1, iie+1, jjs-1, jje+1 ) ! (in)
1525 
1526  ! at (u, v, layer)
1527  if ( .not. twod ) &
1528  call atmos_dyn_fvm_fluxx_xvz_ud1( qflx_lo(:,:,:,xdir), & ! (out)
1529  momx, momy0, dens, & ! (in)
1530  gsqrt(:,:,:,i_uvz), mapf(:,:,:,i_xy), & ! (in)
1531  num_diff(:,:,:,i_momy,xdir), & ! (in)
1532  cdz, twod, & ! (in)
1533  iis-1, iie+1, jjs-1, jje+1 ) ! (in)
1534 
1535  ! at (x, y, layer)
1536  ! note that y-index is added by -1
1537  call atmos_dyn_fvm_fluxy_xvz_ud1( qflx_lo(:,:,:,ydir), & ! (out)
1538  momy, momy0, dens, & ! (in)
1539  gsqrt(:,:,:,i_xyz), mapf(:,:,:,i_xy), & ! (in)
1540  num_diff(:,:,:,i_momy,ydir), & ! (in)
1541  cdz, twod, & ! (in)
1542  iis-1, iie+1, jjs-1, jje+1 ) ! (in)
1543 
1544  endif
1545  enddo
1546  enddo
1547 
1548  if ( flag_fct_momentum ) then
1549 
1550  do j = js, je
1551  do i = is, ie
1552  do k = ks, ke
1553  qflx_hi(k,i,j,zdir) = qflx_hi(k,i,j,zdir) / ( mapf(i,j,1,i_xv) * mapf(i,j,2,i_xv) ) &
1554  + qflx_j13(k,i,j) + qflx_j23(k,i,j)
1555  enddo
1556  enddo
1557  enddo
1558 
1559  do j = js-1, je+1
1560  do i = is-1, ie+1
1561  do k = ks, ke
1562  dens0_uvw(k,i,j) = 0.5_rp * ( dens0(k,i,j) + dens0(k,i,j+1) )
1563  dens_uvw(k,i,j) = 0.5_rp * ( dens_rk(k,i,j) + dens_rk(k,i,j+1) )
1564  enddo
1565  enddo
1566  enddo
1567 
1568  call comm_wait ( vely(:,:,:), 6 )
1569 
1570  call atmos_dyn_fvm_fct( qflx_anti, & ! (out)
1571  vely, dens0_uvw, dens_uvw, & ! (in)
1572  qflx_hi, qflx_lo, & ! (in)
1573  mflx_hi, & ! (in)
1574  rcdz, rcdx, rfdy, & ! (in)
1575  gsqrt(:,:,:,i_xvz), & ! (in)
1576  mapf(:,:,:,i_xv), & ! (in)
1577  twod, dtrk, & ! (in)
1578  flag_fct_along_stream ) ! (in)
1579 
1580  do jjs = js, je, jblock
1581  jje = jjs+jblock-1
1582  do iis = is, ie, iblock
1583  iie = iis+iblock-1
1584 
1585  if ( twod ) then
1586  !--- update momentum(y)
1587  !$omp parallel do private(i,j,k) OMP_SCHEDULE_
1588  do j = jjs, min(jje, jeh)
1589  do k = ks, ke
1590 #ifdef DEBUG
1591  call check( __line__, momy_rk(k,is,j) )
1592  call check( __line__, qflx_anti(k ,is,j ,zdir) )
1593  call check( __line__, qflx_anti(k-1,is,j ,zdir) )
1594  call check( __line__, qflx_anti(k ,is,j ,ydir) )
1595  call check( __line__, qflx_anti(k ,is,j-1,ydir) )
1596 #endif
1597  momy_rk(k,is,j) = momy_rk(k,is,j) &
1598  + dtrk * ( ( ( qflx_anti(k,is,j,zdir) - qflx_anti(k-1,is,j ,zdir) ) * rcdz(k) &
1599  + ( qflx_anti(k,is,j,ydir) - qflx_anti(k ,is,j-1,ydir) ) * rfdy(j) ) ) &
1600  * mapf(is,j,2,i_xv) &
1601  / gsqrt(k,is,j,i_xvz)
1602  enddo
1603  enddo
1604  else
1605  !--- update momentum(y)
1606  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
1607  do j = jjs, min(jje, jeh)
1608  do i = iis, iie
1609  do k = ks, ke
1610 #ifdef DEBUG
1611  call check( __line__, momy_rk(k,i,j) )
1612  call check( __line__, qflx_anti(k ,i ,j ,zdir) )
1613  call check( __line__, qflx_anti(k-1,i ,j ,zdir) )
1614  call check( __line__, qflx_anti(k ,i ,j ,xdir) )
1615  call check( __line__, qflx_anti(k ,i-1,j ,xdir) )
1616  call check( __line__, qflx_anti(k ,i ,j ,ydir) )
1617  call check( __line__, qflx_anti(k ,i ,j-1,ydir) )
1618 #endif
1619  momy_rk(k,i,j) = momy_rk(k,i,j) &
1620  + dtrk * ( ( ( qflx_anti(k,i,j,zdir) - qflx_anti(k-1,i ,j ,zdir) ) * rcdz(k) &
1621  + ( qflx_anti(k,i,j,xdir) - qflx_anti(k ,i-1,j ,xdir) ) * rcdx(i) &
1622  + ( qflx_anti(k,i,j,ydir) - qflx_anti(k ,i ,j-1,ydir) ) * rfdy(j) ) ) &
1623  * mapf(i,j,1,i_xv) * mapf(i,j,2,i_xv) &
1624  / gsqrt(k,i,j,i_xvz)
1625  enddo
1626  enddo
1627  enddo
1628  end if
1629 #ifdef DEBUG
1630  k = iundef; i = iundef; j = iundef
1631 #endif
1632 
1633  enddo
1634  enddo
1635 
1636 #ifdef DEBUG
1637  qflx_lo(:,:,:,:) = undef
1638  qflx_anti(:,:,:,:) = undef
1639 #endif
1640  endif ! FLAG_FCT_MOMENTUM
1641 
1642 #ifdef DEBUG
1643  qflx_hi(ks:,:,:,:) = undef
1644 #endif
1645 
1646  do jjs = js, je, jblock
1647  jje = jjs+jblock-1
1648  do iis = is, ie, iblock
1649  iie = iis+iblock-1
1650 #endif
1651 
1652  !########################################################################
1653  ! Thermodynamic equation
1654  !########################################################################
1655 
1656  !-----< high order flux >-----
1657 
1658  ! at (x, y, w)
1659  call atmos_dyn_fvm_fluxz_xyz( tflx_hi(:,:,:,zdir), & ! (out)
1660  mflx_hi(:,:,:,zdir), pott, gsqrt(:,:,:,i_xyw), & ! (in)
1661  num_diff(:,:,:,i_rhot,zdir), & ! (in)
1662  cdz, & ! (in)
1663  iis, iie, jjs, jje ) ! (in)
1664 
1665  ! at (u, y, z)
1666  if ( .not. twod ) &
1667  call atmos_dyn_fvm_fluxx_xyz( tflx_hi(:,:,:,xdir), & ! (out)
1668  mflx_hi(:,:,:,xdir), pott, gsqrt(:,:,:,i_uyz), & ! (in)
1669  num_diff(:,:,:,i_rhot,xdir), & ! (in)
1670  cdz, & ! (in)
1671  iis, iie, jjs, jje ) ! (in)
1672 
1673  ! at (x, v, z)
1674  call atmos_dyn_fvm_fluxy_xyz( tflx_hi(:,:,:,ydir), & ! (out)
1675  mflx_hi(:,:,:,ydir), pott, gsqrt(:,:,:,i_xvz), & ! (in)
1676  num_diff(:,:,:,i_rhot,ydir), & ! (in)
1677  cdz, & ! (in)
1678  iis, iie, jjs, jje ) ! (in)
1679 
1680  !-----< update rho*theta >-----
1681 
1682  if ( twod ) then
1683  !$omp parallel do private(i,j,k,advcv,advch) OMP_SCHEDULE_
1684  do j = jjs, jje
1685  do k = ks, ke
1686 #ifdef DEBUG
1687  call check( __line__, tflx_hi(k ,is,j ,zdir) )
1688  call check( __line__, tflx_hi(k-1,is,j ,zdir) )
1689  call check( __line__, tflx_hi(k ,is,j ,ydir) )
1690  call check( __line__, tflx_hi(k ,is,j-1,ydir) )
1691  call check( __line__, rhot_t(k,is,j) )
1692  call check( __line__, rhot0(k,is,j) )
1693 #endif
1694  advcv = - ( tflx_hi(k,is,j,zdir) - tflx_hi(k-1,is,j ,zdir) ) * rcdz(k)
1695  advch = - ( tflx_hi(k,is,j,ydir) - tflx_hi(k ,is,j-1,ydir) ) * rcdy(j)
1696  rhot_rk(k,is,j) = rhot0(k,is,j) &
1697  + dtrk * ( ( advcv + advch ) * mapf(is,j,2,i_xy) / gsqrt(k,is,j,i_xyz) &
1698  + rhot_t(k,is,j) )
1699 #ifdef HIST_TEND
1700  if ( lhist ) then
1701  advcv_t(k,is,j,i_rhot) = advcv * mapf(is,j,2,i_xy)/ gsqrt(k,is,j,i_xyz)
1702  advch_t(k,is,j,i_rhot) = advch * mapf(is,j,2,i_xy)/ gsqrt(k,is,j,i_xyz)
1703  endif
1704 #endif
1705  enddo
1706  enddo
1707  else
1708  !$omp parallel do private(i,j,k,advcv,advch) OMP_SCHEDULE_ collapse(2)
1709  do j = jjs, jje
1710  do i = iis, iie
1711  do k = ks, ke
1712 #ifdef DEBUG
1713  call check( __line__, tflx_hi(k ,i ,j ,zdir) )
1714  call check( __line__, tflx_hi(k-1,i ,j ,zdir) )
1715  call check( __line__, tflx_hi(k ,i ,j ,xdir) )
1716  call check( __line__, tflx_hi(k ,i-1,j ,xdir) )
1717  call check( __line__, tflx_hi(k ,i ,j ,ydir) )
1718  call check( __line__, tflx_hi(k ,i ,j-1,ydir) )
1719  call check( __line__, rhot_t(k,i,j) )
1720  call check( __line__, rhot0(k,i,j) )
1721 #endif
1722  advcv = - ( tflx_hi(k,i,j,zdir) - tflx_hi(k-1,i ,j ,zdir) ) * rcdz(k)
1723  advch = - ( tflx_hi(k,i,j,xdir) - tflx_hi(k ,i-1,j ,xdir) ) * rcdx(i) &
1724  - ( tflx_hi(k,i,j,ydir) - tflx_hi(k ,i ,j-1,ydir) ) * rcdy(j)
1725  rhot_rk(k,i,j) = rhot0(k,i,j) &
1726  + dtrk * ( ( advcv + advch ) * mapf(i,j,1,i_xy) * mapf(i,j,2,i_xy) / gsqrt(k,i,j,i_xyz) &
1727  + rhot_t(k,i,j) )
1728 #ifdef HIST_TEND
1729  if ( lhist ) then
1730  advcv_t(k,i,j,i_rhot) = advcv * mapf(i,j,1,i_xy) * mapf(i,j,2,i_xy)/ gsqrt(k,i,j,i_xyz)
1731  advch_t(k,i,j,i_rhot) = advch * mapf(i,j,1,i_xy) * mapf(i,j,2,i_xy)/ gsqrt(k,i,j,i_xyz)
1732  endif
1733 #endif
1734  enddo
1735  enddo
1736  enddo
1737  end if
1738 #ifdef DEBUG
1739  k = iundef; i = iundef; j = iundef
1740 #endif
1741 
1742  enddo
1743  enddo
1744 
1745 #ifndef NO_FCT_DYN
1746 
1747  if ( flag_fct_t ) then
1748 
1749  call comm_vars8( mflx_hi(:,:,:,zdir), 1 )
1750  call comm_vars8( mflx_hi(:,:,:,xdir), 2 )
1751  call comm_vars8( mflx_hi(:,:,:,ydir), 3 )
1752  call comm_wait ( mflx_hi(:,:,:,zdir), 1, .false. )
1753  call comm_wait ( mflx_hi(:,:,:,xdir), 2, .false. )
1754  call comm_wait ( mflx_hi(:,:,:,ydir), 3, .false. )
1755 
1756  if ( .NOT. flag_fct_momentum ) then
1757  call comm_vars8( dens_rk, 1 )
1758  call comm_wait ( dens_rk, 1, .false. )
1759  endif
1760 
1761  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
1762  do j = js-1, je+2
1763  do i = is-1, ie+2
1764  do k = ks, ke
1765  pott(k,i,j) = rhot0(k,i,j) / dens0(k,i,j)
1766  enddo
1767  enddo
1768  enddo
1769 
1770  do jjs = js, je, jblock
1771  jje = jjs+jblock-1
1772  do iis = is, ie, iblock
1773  iie = iis+iblock-1
1774 
1775  ! monotonic flux
1776  ! at (x, y, interface)
1777  call atmos_dyn_fvm_fluxz_xyz_ud1( tflx_lo(:,:,:,zdir), & ! (out)
1778  mflx_hi(:,:,:,zdir), pott, gsqrt(:,:,:,i_xyz), & ! (in)
1779  num_diff(:,:,:,i_rhot,zdir), & ! (in)
1780  cdz, & ! (in)
1781  iis-1, iie+1, jjs-1, jje+1 ) ! (in)
1782 
1783  ! at (u, y, layer)
1784  if ( .not. twod ) &
1785  call atmos_dyn_fvm_fluxx_xyz_ud1( tflx_lo(:,:,:,xdir), & ! (out)
1786  mflx_hi(:,:,:,xdir), pott, gsqrt(:,:,:,i_uyz), & ! (in)
1787  num_diff(:,:,:,i_rhot,xdir), & ! (in)
1788  cdz, & ! (in)
1789  iis-1, iie+1, jjs-1, jje+1 ) ! (in)
1790 
1791  ! at (x, v, layer)
1792  call atmos_dyn_fvm_fluxy_xyz_ud1( tflx_lo(:,:,:,ydir), & ! (out)
1793  mflx_hi(:,:,:,ydir), pott, gsqrt(:,:,:,i_xvz), & ! (in)
1794  num_diff(:,:,:,i_rhot,ydir), & ! (in)
1795  cdz, & ! (in)
1796  iis-1, iie+1, jjs-1, jje+1 ) ! (in)
1797 
1798  enddo
1799  enddo
1800 
1801  call atmos_dyn_fvm_fct( tflx_anti, & ! (out)
1802  pott, dens0, dens_rk, & ! (out)
1803  tflx_hi, tflx_lo, & ! (in)
1804  mflx_hi, & ! (in)
1805  rcdz, rcdx, rcdy, & ! (in)
1806  gsqrt(:,:,:,i_xyz), & ! (in)
1807  mapf(:,:,:,i_xy), & ! (in)
1808  twod, dtrk, & ! (in)
1809  flag_fct_along_stream ) ! (in)
1810 
1811  do jjs = js, je, jblock
1812  jje = jjs+jblock-1
1813  do iis = is, ie, iblock
1814  iie = iis+iblock-1
1815 
1816  !--- update rho*theta
1817  if ( twod ) then
1818  !$omp parallel do private(i,j,k) OMP_SCHEDULE_
1819  do j = jjs, jje
1820  do i = iis, iie
1821  do k = ks, ke
1822 #ifdef DEBUG
1823  call check( __line__, rhot_rk(k,is,j) )
1824  call check( __line__, tflx_anti(k ,is,j ,zdir) )
1825  call check( __line__, tflx_anti(k-1,is,j ,zdir) )
1826  call check( __line__, tflx_anti(k ,is,j ,ydir) )
1827  call check( __line__, tflx_anti(k ,is,j-1,ydir) )
1828 #endif
1829  rhot_rk(k,is,j) = rhot_rk(k,is,j) &
1830  + dtrk * ( ( tflx_anti(k,is,j,zdir) - tflx_anti(k-1,is,j ,zdir) ) * rcdz(k) &
1831  + ( tflx_anti(k,is,j,ydir) - tflx_anti(k ,is,j-1,ydir) ) * rcdy(j) ) &
1832  * mapf(is,j,2,i_xy) &
1833  / gsqrt(k,is,j,i_xyz)
1834  enddo
1835  enddo
1836  enddo
1837  else
1838  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
1839  do j = jjs, jje
1840  do i = iis, iie
1841  do k = ks, ke
1842 #ifdef DEBUG
1843  call check( __line__, rhot_rk(k,i,j) )
1844  call check( __line__, tflx_anti(k ,i ,j ,zdir) )
1845  call check( __line__, tflx_anti(k-1,i ,j ,zdir) )
1846  call check( __line__, tflx_anti(k ,i ,j ,xdir) )
1847  call check( __line__, tflx_anti(k ,i-1,j ,xdir) )
1848  call check( __line__, tflx_anti(k ,i ,j ,ydir) )
1849  call check( __line__, tflx_anti(k ,i ,j-1,ydir) )
1850 #endif
1851  rhot_rk(k,i,j) = rhot_rk(k,i,j) &
1852  + dtrk * ( ( tflx_anti(k,i,j,zdir) - tflx_anti(k-1,i ,j ,zdir) ) * rcdz(k) &
1853  + ( tflx_anti(k,i,j,xdir) - tflx_anti(k ,i-1,j ,xdir) ) * rcdx(i) &
1854  + ( tflx_anti(k,i,j,ydir) - tflx_anti(k ,i ,j-1,ydir) ) * rcdy(j) ) &
1855  * mapf(i,j,1,i_xy) * mapf(i,j,2,i_xy) &
1856  / gsqrt(k,i,j,i_xyz)
1857  enddo
1858  enddo
1859  enddo
1860  end if
1861 #ifdef DEBUG
1862  k = iundef; i = iundef; j = iundef
1863 #endif
1864  enddo
1865  enddo
1866  endif
1867 #endif
1868 
1869 #ifdef HIST_TEND
1870  if ( lhist ) then
1871  call file_history_in(advcv_t(:,:,:,i_dens), 'DENS_t_advcv', 'tendency of density (vert. advection) (w/ HIST_TEND)', 'kg/m3/s' )
1872  call file_history_in(advcv_t(:,:,:,i_momz), 'MOMZ_t_advcv', 'tendency of momentum z (vert. advection) (w/ HIST_TEND)', 'kg/m2/s2', dim_type='ZHXY')
1873  call file_history_in(advcv_t(:,:,:,i_momx), 'MOMX_t_advcv', 'tendency of momentum x (vert. advection) (w/ HIST_TEND)', 'kg/m2/s2', dim_type='ZXHY')
1874  call file_history_in(advcv_t(:,:,:,i_momy), 'MOMY_t_advcv', 'tendency of momentum y (vert. advection) (w/ HIST_TEND)', 'kg/m2/s2', dim_type='ZXYH')
1875  call file_history_in(advcv_t(:,:,:,i_rhot), 'RHOT_t_advcv', 'tendency of rho*theta (vert. advection) (w/ HIST_TEND)', 'K kg/m3/s' )
1876 
1877  call file_history_in(advch_t(:,:,:,i_dens), 'DENS_t_advch', 'tendency of density (horiz. advection) (w/ HIST_TEND)', 'kg/m3/s' )
1878  call file_history_in(advch_t(:,:,:,i_momz), 'MOMZ_t_advch', 'tendency of momentum z (horiz. advection) (w/ HIST_TEND)', 'kg/m2/s2', dim_type='ZHXY')
1879  call file_history_in(advch_t(:,:,:,i_momx), 'MOMX_t_advch', 'tendency of momentum x (horiz. advection) (w/ HIST_TEND)', 'kg/m2/s2', dim_type='ZXHY')
1880  call file_history_in(advch_t(:,:,:,i_momy), 'MOMY_t_advch', 'tendency of momentum y (horiz. advection) (w/ HIST_TEND)', 'kg/m2/s2', dim_type='ZXYH')
1881  call file_history_in(advch_t(:,:,:,i_rhot), 'RHOT_t_advch', 'tendency of rho*theta (horiz. advection) (w/ HIST_TEND)', 'K kg/m3/s' )
1882 
1883  call file_history_in(pg_t(:,:,:,1), 'MOMZ_t_pg', 'tendency of momentum z (pressure gradient) (w/ HIST_TEND)', 'kg/m2/s2', dim_type='ZHXY')
1884  if ( .not. twod ) &
1885  call file_history_in(pg_t(:,:,:,2), 'MOMX_t_pg', 'tendency of momentum x (pressure gradient) (w/ HIST_TEND)', 'kg/m2/s2', dim_type='ZXHY')
1886  call file_history_in(pg_t(:,:,:,3), 'MOMY_t_pg', 'tendency of momentum y (pressure gradient) (w/ HIST_TEND)', 'kg/m2/s2', dim_type='ZXYH')
1887 
1888  call file_history_in(wdmp_t(:,:,:), 'MOMZ_t_wdamp', 'tendency of momentum z (Rayleigh damping) (w/ HIST_TEND)', 'kg/m2/s2', dim_type='ZHXY')
1889 
1890  call file_history_in(ddiv_t(:,:,:,1), 'MOMZ_t_ddiv', 'tendency of momentum z (divergence damping) (w/ HIST_TEND)', 'kg/m2/s2', dim_type='ZHXY')
1891  if ( .not. twod ) &
1892  call file_history_in(ddiv_t(:,:,:,2), 'MOMX_t_ddiv', 'tendency of momentum x (divergence damping) (w/ HIST_TEND)', 'kg/m2/s2', dim_type='ZXHY')
1893  call file_history_in(ddiv_t(:,:,:,3), 'MOMY_t_ddiv', 'tendency of momentum y (divergence damping) (w/ HIST_TEND)', 'kg/m2/s2', dim_type='ZXYH')
1894 
1895  call file_history_in(cf_t(:,:,:,1), 'MOMX_t_cf', 'tendency of momentum x (coliolis force) (w/ HIST_TEND)', 'kg/m2/s2', dim_type='ZXHY')
1896  call file_history_in(cf_t(:,:,:,2), 'MOMY_t_cf', 'tendency of momentum y (coliolis force) (w/ HIST_TEND)', 'kg/m2/s2', dim_type='ZXYH')
1897  endif
1898 #endif
1899 
1900  return

References scale_atmos_dyn_fvm_fct::atmos_dyn_fvm_fct(), scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxj13_uyz, scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxj13_xvz, scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxj13_xyw, scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxj23_uyz, scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxj23_xvz, scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxj23_xyw, scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxx_uyz, scale_atmos_dyn_fvm_flux_ud1::atmos_dyn_fvm_fluxx_uyz_ud1(), scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxx_xvz, scale_atmos_dyn_fvm_flux_ud1::atmos_dyn_fvm_fluxx_xvz_ud1(), scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxx_xyw, scale_atmos_dyn_fvm_flux_ud1::atmos_dyn_fvm_fluxx_xyw_ud1(), scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxx_xyz, scale_atmos_dyn_fvm_flux_ud1::atmos_dyn_fvm_fluxx_xyz_ud1(), scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxy_uyz, scale_atmos_dyn_fvm_flux_ud1::atmos_dyn_fvm_fluxy_uyz_ud1(), scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxy_xvz, scale_atmos_dyn_fvm_flux_ud1::atmos_dyn_fvm_fluxy_xvz_ud1(), scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxy_xyw, scale_atmos_dyn_fvm_flux_ud1::atmos_dyn_fvm_fluxy_xyw_ud1(), scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxy_xyz, scale_atmos_dyn_fvm_flux_ud1::atmos_dyn_fvm_fluxy_xyz_ud1(), scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxz_uyz, scale_atmos_dyn_fvm_flux_ud1::atmos_dyn_fvm_fluxz_uyz_ud1(), scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxz_xvz, scale_atmos_dyn_fvm_flux_ud1::atmos_dyn_fvm_fluxz_xvz_ud1(), scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxz_xyw, scale_atmos_dyn_fvm_flux_ud1::atmos_dyn_fvm_fluxz_xyw_ud1(), scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxz_xyz, scale_atmos_dyn_fvm_flux_ud1::atmos_dyn_fvm_fluxz_xyz_ud1(), scale_debug::check(), scale_const::const_eps, scale_const::const_grav, scale_const::const_pre00, scale_index::i_dens, scale_index::i_momx, scale_index::i_momy, scale_index::i_momz, scale_index::i_rhot, scale_atmos_grid_cartesc_index::i_uv, scale_atmos_grid_cartesc_index::i_uvz, scale_atmos_grid_cartesc_index::i_uy, scale_atmos_grid_cartesc_index::i_uyw, scale_atmos_grid_cartesc_index::i_uyz, scale_atmos_grid_cartesc_index::i_xv, scale_atmos_grid_cartesc_index::i_xvw, scale_atmos_grid_cartesc_index::i_xvz, scale_atmos_grid_cartesc_index::i_xy, scale_atmos_grid_cartesc_index::i_xyw, scale_atmos_grid_cartesc_index::i_xyz, scale_atmos_grid_cartesc_index::ia, scale_atmos_grid_cartesc_index::iblock, scale_atmos_grid_cartesc_index::ie, scale_atmos_grid_cartesc_index::ieh, scale_atmos_grid_cartesc_index::ihalo, scale_atmos_grid_cartesc_index::is, scale_atmos_grid_cartesc_index::jblock, scale_atmos_grid_cartesc_index::je, scale_atmos_grid_cartesc_index::jeh, scale_atmos_grid_cartesc_index::jhalo, scale_atmos_grid_cartesc_index::js, scale_tracer::k, scale_atmos_grid_cartesc_index::ka, scale_atmos_grid_cartesc_index::ke, scale_atmos_grid_cartesc_index::ks, scale_atmos_grid_cartesc_index::xdir, scale_atmos_grid_cartesc_index::ydir, and scale_atmos_grid_cartesc_index::zdir.

Referenced by scale_atmos_dyn_tstep_short::atmos_dyn_tstep_short_regist().

Here is the call graph for this function:
Here is the caller graph for this function:
scale_const::const_grav
real(rp), public const_grav
standard acceleration of gravity [m/s2]
Definition: scale_const.F90:49
scale_atmos_dyn_fvm_flux_ud1::atmos_dyn_fvm_fluxz_xvz_ud1
subroutine, public atmos_dyn_fvm_fluxz_xvz_ud1(flux, mom, val, DENS, GSQRT, J33G, num_diff, CDZ, TwoD, IIS, IIE, JJS, JJE)
calculation z-flux at XV
Definition: scale_atmos_dyn_fvm_flux_ud1.F90:1241
scale_atmos_grid_cartesc_index::i_uy
integer, public i_uy
Definition: scale_atmos_grid_cartesC_index.F90:100
scale_atmos_grid_cartesc_index::ke
integer, public ke
end point of inner domain: z, local
Definition: scale_atmos_grid_cartesC_index.F90:52
scale_atmos_dyn_fvm_flux_ud1::atmos_dyn_fvm_fluxx_xyw_ud1
subroutine, public atmos_dyn_fvm_fluxx_xyw_ud1(flux, mom, val, DENS, GSQRT, MAPF, num_diff, CDZ, TwoD, IIS, IIE, JJS, JJE)
calculation X-flux at XYW
Definition: scale_atmos_dyn_fvm_flux_ud1.F90:578
scale_atmos_grid_cartesc_index::i_xv
integer, public i_xv
Definition: scale_atmos_grid_cartesC_index.F90:101
scale_prc::prc_abort
subroutine, public prc_abort
Abort Process.
Definition: scale_prc.F90:350
scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxx_uyz
procedure(flux_mom), pointer, public atmos_dyn_fvm_fluxx_uyz
Definition: scale_atmos_dyn_fvm_flux.F90:214
scale_index::i_momz
integer, parameter, public i_momz
Definition: scale_index.F90:29
scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxy_xyz
procedure(flux_phi), pointer, public atmos_dyn_fvm_fluxy_xyz
Definition: scale_atmos_dyn_fvm_flux.F90:175
scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxj23_xyw
procedure(flux_j), pointer, public atmos_dyn_fvm_fluxj23_xyw
Definition: scale_atmos_dyn_fvm_flux.F90:203
scale_atmos_dyn_fvm_fct::atmos_dyn_fvm_fct
subroutine, public atmos_dyn_fvm_fct(qflx_anti, phi_in, DENS0, DENS, qflx_hi, qflx_lo, mflx_hi, rdz, rdx, rdy, GSQRT, MAPF, TwoD, dt, flag_vect)
Setup.
Definition: scale_atmos_dyn_fvm_fct.F90:72
scale_index::i_momx
integer, parameter, public i_momx
Definition: scale_index.F90:30
scale_atmos_grid_cartesc_index::ihalo
integer, public ihalo
Definition: scale_atmos_grid_cartesC_index.F90:44
scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxj23_uyz
procedure(flux_j), pointer, public atmos_dyn_fvm_fluxj23_uyz
Definition: scale_atmos_dyn_fvm_flux.F90:230
scale_atmos_grid_cartesc_index::ka
integer, public ka
Definition: scale_atmos_grid_cartesC_index.F90:47
scale_atmos_dyn_fvm_flux
module scale_atmos_dyn_fvm_flux
Definition: scale_atmos_dyn_fvm_flux.F90:13
scale_atmos_grid_cartesc_index::jblock
integer, public jblock
block size for cache blocking: y
Definition: scale_atmos_grid_cartesC_index.F90:41
scale_atmos_grid_cartesc_index::i_uv
integer, public i_uv
Definition: scale_atmos_grid_cartesC_index.F90:102
scale_const::const_eps
real(rp), public const_eps
small number
Definition: scale_const.F90:35
scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxj13_xyw
procedure(flux_j), pointer, public atmos_dyn_fvm_fluxj13_xyw
Definition: scale_atmos_dyn_fvm_flux.F90:198
scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxx_xyz
procedure(flux_phi), pointer, public atmos_dyn_fvm_fluxx_xyz
Definition: scale_atmos_dyn_fvm_flux.F90:170
scale_atmos_grid_cartesc_index::iblock
integer, public iblock
block size for cache blocking: x
Definition: scale_atmos_grid_cartesC_index.F90:40
scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxx_xvz
procedure(flux_mom), pointer, public atmos_dyn_fvm_fluxx_xvz
Definition: scale_atmos_dyn_fvm_flux.F90:241
scale_atmos_grid_cartesc_index::jeh
integer, public jeh
end point of inner domain: y, local (half level)
Definition: scale_atmos_grid_cartesC_index.F90:69
scale_index::i_dens
integer, parameter, public i_dens
Definition: scale_index.F90:28
scale_atmos_dyn_fvm_flux_ud1::atmos_dyn_fvm_fluxx_xyz_ud1
subroutine, public atmos_dyn_fvm_fluxx_xyz_ud1(flux, mflx, val, GSQRT, num_diff, CDZ, IIS, IIE, JJS, JJE)
calculation X-flux at XYZ
Definition: scale_atmos_dyn_fvm_flux_ud1.F90:221
scale_file_history
module file_history
Definition: scale_file_history.F90:15
scale_index::i_momy
integer, parameter, public i_momy
Definition: scale_index.F90:31
scale_atmos_grid_cartesc_index::i_xyz
integer, public i_xyz
Definition: scale_atmos_grid_cartesC_index.F90:91
scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxy_xyw
procedure(flux_mom), pointer, public atmos_dyn_fvm_fluxy_xyw
Definition: scale_atmos_dyn_fvm_flux.F90:192
scale_atmos_dyn_fvm_flux_ud1::atmos_dyn_fvm_fluxy_xyz_ud1
subroutine, public atmos_dyn_fvm_fluxy_xyz_ud1(flux, mflx, val, GSQRT, num_diff, CDZ, IIS, IIE, JJS, JJE)
calculation Y-flux at XYZ
Definition: scale_atmos_dyn_fvm_flux_ud1.F90:271
scale_prc
module PROCESS
Definition: scale_prc.F90:11
scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxx_xyw
procedure(flux_mom), pointer, public atmos_dyn_fvm_fluxx_xyw
Definition: scale_atmos_dyn_fvm_flux.F90:187
scale_atmos_dyn_fvm_flux_ud1::atmos_dyn_fvm_fluxx_uyz_ud1
subroutine, public atmos_dyn_fvm_fluxx_uyz_ud1(flux, mom, val, DENS, GSQRT, MAPF, num_diff, CDZ, TwoD, IIS, IIE, JJS, JJE)
calculation X-flux at UY
Definition: scale_atmos_dyn_fvm_flux_ud1.F90:1085
scale_atmos_grid_cartesc_index::i_xy
integer, public i_xy
Definition: scale_atmos_grid_cartesC_index.F90:99
scale_atmos_grid_cartesc_index::ie
integer, public ie
end point of inner domain: x, local
Definition: scale_atmos_grid_cartesC_index.F90:54
scale_atmos_grid_cartesc_index::i_uyz
integer, public i_uyz
Definition: scale_atmos_grid_cartesC_index.F90:95
scale_atmos_dyn_fvm_flux_ud1::atmos_dyn_fvm_fluxz_xyz_ud1
subroutine, public atmos_dyn_fvm_fluxz_xyz_ud1(flux, mflx, val, GSQRT, num_diff, CDZ, IIS, IIE, JJS, JJE)
calculation z-flux at XYZ
Definition: scale_atmos_dyn_fvm_flux_ud1.F90:143
scale_atmos_grid_cartesc_index::i_xvw
integer, public i_xvw
Definition: scale_atmos_grid_cartesC_index.F90:94
scale_atmos_grid_cartesc_index
module atmosphere / grid / cartesC index
Definition: scale_atmos_grid_cartesC_index.F90:12
scale_const
module CONSTANT
Definition: scale_const.F90:11
scale_atmos_grid_cartesc_index::ia
integer, public ia
Definition: scale_atmos_grid_cartesC_index.F90:48
scale_atmos_dyn_fvm_flux_ud1
module scale_atmos_dyn_fvm_flux_ud1
Definition: scale_atmos_dyn_fvm_flux_ud1.F90:16
scale_atmos_grid_cartesc_index::zdir
integer, parameter, public zdir
Definition: scale_atmos_grid_cartesC_index.F90:32
scale_atmos_grid_cartesc_index::i_xvz
integer, public i_xvz
Definition: scale_atmos_grid_cartesC_index.F90:96
scale_atmos_dyn_fvm_flux_ud1::atmos_dyn_fvm_fluxy_xvz_ud1
subroutine, public atmos_dyn_fvm_fluxy_xvz_ud1(flux, mom, val, DENS, GSQRT, MAPF, num_diff, CDZ, TwoD, IIS, IIE, JJS, JJE)
calculation Y-flux at XV
Definition: scale_atmos_dyn_fvm_flux_ud1.F90:1548
scale_atmos_dyn_fvm_flux_ud1::atmos_dyn_fvm_fluxz_uyz_ud1
subroutine, public atmos_dyn_fvm_fluxz_uyz_ud1(flux, mom, val, DENS, GSQRT, J33G, num_diff, CDZ, TwoD, IIS, IIE, JJS, JJE)
calculation z-flux at UY
Definition: scale_atmos_dyn_fvm_flux_ud1.F90:747
scale_atmos_grid_cartesc_index::ieh
integer, public ieh
end point of inner domain: x, local (half level)
Definition: scale_atmos_grid_cartesC_index.F90:68
scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxz_xvz
procedure(flux_z), pointer, public atmos_dyn_fvm_fluxz_xvz
Definition: scale_atmos_dyn_fvm_flux.F90:236
scale_index::i_rhot
integer, parameter, public i_rhot
Definition: scale_index.F90:32
scale_atmos_grid_cartesc_index::i_uyw
integer, public i_uyw
Definition: scale_atmos_grid_cartesC_index.F90:93
scale_atmos_grid_cartesc_index::is
integer, public is
start point of inner domain: x, local
Definition: scale_atmos_grid_cartesC_index.F90:53
scale_atmos_grid_cartesc_index::jhalo
integer, public jhalo
Definition: scale_atmos_grid_cartesC_index.F90:45
scale_atmos_dyn_fvm_flux_ud1::atmos_dyn_fvm_fluxz_xyw_ud1
subroutine, public atmos_dyn_fvm_fluxz_xyw_ud1(flux, mom, val, DENS, GSQRT, J33G, num_diff, CDZ, FDZ, dtrk, IIS, IIE, JJS, JJE)
calculation z-flux at XYW
Definition: scale_atmos_dyn_fvm_flux_ud1.F90:324
scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxz_uyz
procedure(flux_z), pointer, public atmos_dyn_fvm_fluxz_uyz
Definition: scale_atmos_dyn_fvm_flux.F90:209
scale_atmos_grid_cartesc_index::xdir
integer, parameter, public xdir
Definition: scale_atmos_grid_cartesC_index.F90:33
scale_atmos_grid_cartesc_index::ks
integer, public ks
start point of inner domain: z, local
Definition: scale_atmos_grid_cartesC_index.F90:51
scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxj23_xvz
procedure(flux_j), pointer, public atmos_dyn_fvm_fluxj23_xvz
Definition: scale_atmos_dyn_fvm_flux.F90:257
scale_comm_cartesc
module COMMUNICATION
Definition: scale_comm_cartesC.F90:11
scale_atmos_grid_cartesc_index::js
integer, public js
start point of inner domain: y, local
Definition: scale_atmos_grid_cartesC_index.F90:55
scale_atmos_grid_cartesc_index::ydir
integer, parameter, public ydir
Definition: scale_atmos_grid_cartesC_index.F90:34
scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxj13_xvz
procedure(flux_j), pointer, public atmos_dyn_fvm_fluxj13_xvz
Definition: scale_atmos_dyn_fvm_flux.F90:252
scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxz_xyz
procedure(flux_phi), pointer, public atmos_dyn_fvm_fluxz_xyz
Definition: scale_atmos_dyn_fvm_flux.F90:165
scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxj13_uyz
procedure(flux_j), pointer, public atmos_dyn_fvm_fluxj13_uyz
Definition: scale_atmos_dyn_fvm_flux.F90:225
scale_atmos_dyn_fvm_flux_ud1::atmos_dyn_fvm_fluxy_xyw_ud1
subroutine, public atmos_dyn_fvm_fluxy_xyw_ud1(flux, mom, val, DENS, GSQRT, MAPF, num_diff, CDZ, TwoD, IIS, IIE, JJS, JJE)
calculation Y-flux at XYW
Definition: scale_atmos_dyn_fvm_flux_ud1.F90:662
scale_atmos_dyn_fvm_flux_ud1::atmos_dyn_fvm_fluxx_xvz_ud1
subroutine, public atmos_dyn_fvm_fluxx_xvz_ud1(flux, mom, val, DENS, GSQRT, MAPF, num_diff, CDZ, TwoD, IIS, IIE, JJS, JJE)
calculation X-flux at XV
Definition: scale_atmos_dyn_fvm_flux_ud1.F90:1488
scale_atmos_grid_cartesc_index::i_uvz
integer, public i_uvz
Definition: scale_atmos_grid_cartesC_index.F90:97
scale_atmos_dyn_fvm_fct
module Atmosphere / Dynamics common
Definition: scale_atmos_dyn_fvm_fct.F90:12
scale_const::const_pre00
real(rp), public const_pre00
pressure reference [Pa]
Definition: scale_const.F90:97
scale_atmos_grid_cartesc_index::je
integer, public je
end point of inner domain: y, local
Definition: scale_atmos_grid_cartesC_index.F90:56
scale_atmos_dyn_fvm_flux_ud1::atmos_dyn_fvm_fluxy_uyz_ud1
subroutine, public atmos_dyn_fvm_fluxy_uyz_ud1(flux, mom, val, DENS, GSQRT, MAPF, num_diff, CDZ, TwoD, IIS, IIE, JJS, JJE)
calculation Y-flux at UY
Definition: scale_atmos_dyn_fvm_flux_ud1.F90:1147
scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxz_xyw
procedure(flux_wz), pointer, public atmos_dyn_fvm_fluxz_xyw
Definition: scale_atmos_dyn_fvm_flux.F90:182
scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxy_xvz
procedure(flux_mom), pointer, public atmos_dyn_fvm_fluxy_xvz
Definition: scale_atmos_dyn_fvm_flux.F90:246
scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxy_uyz
procedure(flux_mom), pointer, public atmos_dyn_fvm_fluxy_uyz
Definition: scale_atmos_dyn_fvm_flux.F90:219
scale_atmos_grid_cartesc_index::i_xyw
integer, public i_xyw
Definition: scale_atmos_grid_cartesC_index.F90:92