SCALE-RM
scale_atmos_dyn_tstep_short_fvm_heve.F90
Go to the documentation of this file.
1 !-------------------------------------------------------------------------------
10 !-------------------------------------------------------------------------------
11 #include "scalelib.h"
12 module scale_atmos_dyn_tstep_short_fvm_heve
13  !-----------------------------------------------------------------------------
14  !
15  !++ used modules
16  !
17  use scale_precision
18  use scale_io
19  use scale_prof
20  use scale_atmos_grid_cartesc_index
21  use scale_index
22  use scale_tracer
23  use scale_prc
24 #if defined DEBUG || defined QUICKDEBUG
25  use scale_debug, only: &
26  check
27  use scale_const, only: &
28  undef => const_undef, &
29  iundef => const_undef2
30 #endif
31  !-----------------------------------------------------------------------------
32  implicit none
33  private
34  !-----------------------------------------------------------------------------
35  !
36  !++ Public procedure
37  !
38  public :: atmos_dyn_tstep_short_fvm_heve_regist
39  public :: atmos_dyn_tstep_short_fvm_heve_setup
40  public :: atmos_dyn_tstep_short_fvm_heve
41 
42  !-----------------------------------------------------------------------------
43  !
44  !++ Public parameters & variables
45  !
46  !-----------------------------------------------------------------------------
47  !
48  !++ Private procedure
49  !
50 #if 1
51 #define F2H(k,p,idx) (CDZ(k+p-1)*GSQRT(k+p-1,i,j,idx)/(CDZ(k)*GSQRT(k,i,j,idx)+CDZ(k+1)*GSQRT(k+1,i,j,idx)))
52 #else
53 #define F2H(k,p,idx) 0.5_RP
54 #endif
55 
56  !-----------------------------------------------------------------------------
57  !
58  !++ Private parameters & variables
59  !
60  integer, parameter :: VA_FVM_HEVE = 0
61 
62  !-----------------------------------------------------------------------------
63 contains
64  !-----------------------------------------------------------------------------
66  subroutine atmos_dyn_tstep_short_fvm_heve_regist( &
67  ATMOS_DYN_TYPE, &
68  VA_out, &
69  VAR_NAME, &
70  VAR_DESC, &
71  VAR_UNIT )
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
100  end subroutine atmos_dyn_tstep_short_fvm_heve_regist
101 
102  !-----------------------------------------------------------------------------
104  subroutine atmos_dyn_tstep_short_fvm_heve_setup
105  return
106  end subroutine atmos_dyn_tstep_short_fvm_heve_setup
107 
108  !-----------------------------------------------------------------------------
109  subroutine atmos_dyn_tstep_short_fvm_heve( &
110  DENS_RK, MOMZ_RK, MOMX_RK, MOMY_RK, RHOT_RK, &
111  PROG_RK, &
112  mflx_hi, tflx_hi, &
113  DENS0, MOMZ0, MOMX0, MOMY0, RHOT0, &
114  DENS, MOMZ, MOMX, MOMY, RHOT, &
115  DENS_t, MOMZ_t, MOMX_t, MOMY_t, RHOT_t, &
116  PROG0, PROG, &
117  DPRES0, RT2P, CORIOLI, &
118  num_diff, wdamp_coef, divdmp_coef, DDIV, &
119  FLAG_FCT_MOMENTUM, FLAG_FCT_T, &
120  FLAG_FCT_ALONG_STREAM, &
121  CDZ, FDZ, FDX, FDY, &
122  RCDZ, RCDX, RCDY, RFDZ, RFDX, RFDY, &
123  PHI, GSQRT, J13G, J23G, J33G, MAPF, &
124  REF_dens, REF_rhot, &
125  BND_W, BND_E, BND_S, BND_N, TwoD, &
126  dtrk, last )
127  use scale_atmos_grid_cartesc_index
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: &
136  atmos_dyn_fvm_fct
137  use scale_atmos_dyn_fvm_flux_ud1, only: &
138  atmos_dyn_fvm_fluxz_xyw_ud1, &
139  atmos_dyn_fvm_fluxx_xyw_ud1, &
140  atmos_dyn_fvm_fluxy_xyw_ud1, &
141  atmos_dyn_fvm_fluxz_uyz_ud1, &
142  atmos_dyn_fvm_fluxx_uyz_ud1, &
143  atmos_dyn_fvm_fluxy_uyz_ud1, &
144  atmos_dyn_fvm_fluxz_xvz_ud1, &
145  atmos_dyn_fvm_fluxx_xvz_ud1, &
146  atmos_dyn_fvm_fluxy_xvz_ud1, &
147  atmos_dyn_fvm_fluxz_xyz_ud1, &
148  atmos_dyn_fvm_fluxx_xyz_ud1, &
149  atmos_dyn_fvm_fluxy_xyz_ud1
150  use scale_atmos_dyn_fvm_flux, only: &
151  atmos_dyn_fvm_fluxz_xyz, &
152  atmos_dyn_fvm_fluxx_xyz, &
153  atmos_dyn_fvm_fluxy_xyz, &
154  atmos_dyn_fvm_fluxz_xyw, &
155  atmos_dyn_fvm_fluxj13_xyw, &
156  atmos_dyn_fvm_fluxj23_xyw, &
157  atmos_dyn_fvm_fluxx_xyw, &
158  atmos_dyn_fvm_fluxy_xyw, &
159  atmos_dyn_fvm_fluxz_uyz, &
160  atmos_dyn_fvm_fluxj13_uyz, &
161  atmos_dyn_fvm_fluxj23_uyz, &
162  atmos_dyn_fvm_fluxx_uyz, &
163  atmos_dyn_fvm_fluxy_uyz, &
164  atmos_dyn_fvm_fluxz_xvz, &
165  atmos_dyn_fvm_fluxj13_xvz, &
166  atmos_dyn_fvm_fluxj23_xvz, &
167  atmos_dyn_fvm_fluxx_xvz, &
168  atmos_dyn_fvm_fluxy_xvz
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
1901  end subroutine atmos_dyn_tstep_short_fvm_heve
1902 
1903 end module scale_atmos_dyn_tstep_short_fvm_heve
scale_const::const_grav
real(rp), public const_grav
standard acceleration of gravity [m/s2]
Definition: scale_const.F90:46
scale_atmos_dyn_tstep_short_fvm_heve::atmos_dyn_tstep_short_fvm_heve_setup
subroutine, public atmos_dyn_tstep_short_fvm_heve_setup
Setup.
Definition: scale_atmos_dyn_tstep_short_fvm_heve.F90:105
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:1150
scale_atmos_grid_cartesc_index::i_uy
integer, public i_uy
Definition: scale_atmos_grid_cartesC_index.F90:99
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:541
scale_atmos_grid_cartesc_index::i_xv
integer, public i_xv
Definition: scale_atmos_grid_cartesC_index.F90:100
scale_prc::prc_abort
subroutine, public prc_abort
Abort Process.
Definition: scale_prc.F90:342
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:212
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:173
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:201
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_index
module Index
Definition: scale_index.F90:11
scale_atmos_grid_cartesc_index::ihalo
integer, public ihalo
Definition: scale_atmos_grid_cartesC_index.F90:44
scale_const::const_undef2
integer, parameter, public const_undef2
undefined value (INT2)
Definition: scale_const.F90:38
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:228
scale_precision
module PRECISION
Definition: scale_precision.F90:14
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:101
scale_const::const_eps
real(rp), public const_eps
small number
Definition: scale_const.F90:33
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:196
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:168
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:239
scale_atmos_dyn_tstep_short_fvm_heve::atmos_dyn_tstep_short_fvm_heve_regist
subroutine, public atmos_dyn_tstep_short_fvm_heve_regist(ATMOS_DYN_TYPE, VA_out, VAR_NAME, VAR_DESC, VAR_UNIT)
Register.
Definition: scale_atmos_dyn_tstep_short_fvm_heve.F90:72
scale_atmos_dyn_tstep_short_fvm_heve::atmos_dyn_tstep_short_fvm_heve
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)
Definition: scale_atmos_dyn_tstep_short_fvm_heve.F90:127
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:68
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:212
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:90
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:190
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:260
scale_index::va
integer, public va
Definition: scale_index.F90:35
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:185
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:1000
scale_atmos_grid_cartesc_index::i_xy
integer, public i_xy
Definition: scale_atmos_grid_cartesC_index.F90:98
scale_precision::rp
integer, parameter, public rp
Definition: scale_precision.F90:41
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:94
scale_io
module STDIO
Definition: scale_io.F90:10
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:142
scale_atmos_grid_cartesc_index::i_xvw
integer, public i_xvw
Definition: scale_atmos_grid_cartesC_index.F90:93
scale_tracer::k
real(rp), public k
Definition: scale_tracer.F90:44
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_debug::check
subroutine, public check(current_line, v)
Undefined value checker.
Definition: scale_debug.F90:56
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:95
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:1431
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:694
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:67
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:234
scale_prof
module profiler
Definition: scale_prof.F90:11
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:92
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_grid_cartesc_index::ja
integer, public ja
Definition: scale_atmos_grid_cartesC_index.F90:49
scale_tracer
module TRACER
Definition: scale_tracer.F90:12
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:311
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:207
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_debug
module DEBUG
Definition: scale_debug.F90:11
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:255
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:250
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:163
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:223
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:617
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:1373
scale_atmos_grid_cartesc_index::i_uvz
integer, public i_uvz
Definition: scale_atmos_grid_cartesC_index.F90:96
scale_const::const_undef
real(rp), public const_undef
Definition: scale_const.F90:41
scale_atmos_dyn_fvm_fct
module Atmosphere / Dynamics common
Definition: scale_atmos_dyn_fvm_fct.F90:12
scale_atmos_dyn_tstep_short_fvm_heve
module Atmosphere / Dynamics RK
Definition: scale_atmos_dyn_tstep_short_fvm_heve.F90:12
scale_const::const_pre00
real(rp), public const_pre00
pressure reference [Pa]
Definition: scale_const.F90:88
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:1060
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:180
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:244
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:217
scale_atmos_grid_cartesc_index::i_xyw
integer, public i_xyw
Definition: scale_atmos_grid_cartesC_index.F90:91
Generated by   doxygen 1.8.17