SCALE-RM
scale_atmos_dyn_tstep_short_fvm_hivi.F90
Go to the documentation of this file.
1 !-------------------------------------------------------------------------------
13 !-------------------------------------------------------------------------------
14 #define HIVI_BICGSTAB 1
15 #include "scalelib.h"
16 module scale_atmos_dyn_tstep_short_fvm_hivi
17  !-----------------------------------------------------------------------------
18  !
19  !++ used modules
20  !
21  use mpi
22  use scale_precision
23  use scale_io
24  use scale_prof
25  use scale_atmos_grid_cartesc_index
26  use scale_index
27  use scale_tracer
28 #if defined DEBUG || defined QUICKDEBUG
29  use scale_debug, only: &
30  check
31  use scale_const, only: &
32  undef => const_undef, &
33  iundef => const_undef2
34 #endif
35  !-----------------------------------------------------------------------------
36  implicit none
37  private
38  !-----------------------------------------------------------------------------
39  !
40  !++ Public procedure
41  !
42  public :: atmos_dyn_tstep_short_fvm_hivi_regist
43  public :: atmos_dyn_tstep_short_fvm_hivi_setup
44  public :: atmos_dyn_tstep_short_fvm_hivi
45 
46  !-----------------------------------------------------------------------------
47  !
48  !++ Public parameters & variables
49  !
50  !-----------------------------------------------------------------------------
51  !
52  !++ Private procedure
53  !
54  !-----------------------------------------------------------------------------
55  !
56  !++ Private parameters & variables
57  !
58  integer, private, parameter :: va_fvm_hivi = 0
59 
60  integer, private :: itmax
61  real(RP), private :: epsilon
62 
63  integer, private :: mtype ! MPI DATATYPE
64 
65  ! tentative
66  real(RP), private, parameter :: fact_n = 7.0_rp / 12.0_rp
67  real(RP), private, parameter :: fact_f = -1.0_rp / 12.0_rp
68 
69  !-----------------------------------------------------------------------------
70 contains
71  !-----------------------------------------------------------------------------
73  subroutine atmos_dyn_tstep_short_fvm_hivi_regist( &
74  ATMOS_DYN_TYPE, &
75  VA_out, &
76  VAR_NAME, &
77  VAR_DESC, &
78  VAR_UNIT )
79  use scale_prc, only: &
80  prc_abort
81  implicit none
82 
83  character(len=*), intent(in) :: ATMOS_DYN_TYPE
84  integer, intent(out) :: VA_out
85  character(len=H_SHORT), intent(out) :: VAR_NAME(:)
86  character(len=H_MID), intent(out) :: VAR_DESC(:)
87  character(len=H_SHORT), intent(out) :: VAR_UNIT(:)
88  !---------------------------------------------------------------------------
89 
90  if ( atmos_dyn_type /= 'FVM-HIVI' .AND. atmos_dyn_type /= 'HIVI' ) then
91  log_error("ATMOS_DYN_Tstep_short_fvm_hivi_regist",*) 'ATMOS_DYN_TYPE is not FVM-HIVI. Check!'
92  call prc_abort
93  endif
94 
95  va_out = va_fvm_hivi
96  var_name(:) = ""
97  var_desc(:) = ""
98  var_unit(:) = ""
99 
100  log_newline
101  log_info("ATMOS_DYN_Tstep_short_fvm_hivi_regist",*) 'Register additional prognostic variables (HIVI)'
102  if ( va_out < 1 ) then
103  log_info_cont(*) '=> nothing.'
104  endif
105 
106  return
107  end subroutine atmos_dyn_tstep_short_fvm_hivi_regist
108 
109  !-----------------------------------------------------------------------------
111  subroutine atmos_dyn_tstep_short_fvm_hivi_setup
112  use scale_prc, only: &
113  prc_abort
114  implicit none
115 
116  namelist / param_atmos_dyn_tstep_fvm_hivi / &
117  itmax, &
118  epsilon
119 
120  integer :: ierr
121  !---------------------------------------------------------------------------
122 
123  log_info("ATMOS_DYN_Tstep_short_fvm_hivi_setup",*) 'HIVI Setup'
124 #ifdef HIVI_BICGSTAB
125  log_info("ATMOS_DYN_Tstep_short_fvm_hivi_setup",*) 'USING Bi-CGSTAB'
126 #else
127  log_info("ATMOS_DYN_Tstep_short_fvm_hivi_setup",*) 'USING Multi-Grid'
128  log_error("ATMOS_DYN_Tstep_short_fvm_hivi_setup",*) 'Not Implemented yet'
129  call prc_abort
130 #endif
131 
132  ! currently, vertical difference scheme for potential temperature is the CD4
133  itmax = 100
134  epsilon = 0.1_rp ** (rp*2)
135 
136  !--- read namelist
137  rewind(io_fid_conf)
138  read(io_fid_conf,nml=param_atmos_dyn_tstep_fvm_hivi,iostat=ierr)
139 
140  if( ierr < 0 ) then !--- missing
141  log_info("ATMOS_DYN_Tstep_short_fvm_hivi_setup",*) 'Not found namelist. Default used.'
142  elseif( ierr > 0 ) then !--- fatal error
143  log_error("ATMOS_DYN_Tstep_short_fvm_hivi_setup",*) 'Not appropriate names in namelist PARAM_ATMOS_DYN_TSTEP_FVM_HIVI. Check!'
144  call prc_abort
145  endif
146  log_nml(param_atmos_dyn_tstep_fvm_hivi)
147 
148  if ( rp == dp ) then
149  mtype = mpi_double_precision
150  elseif( rp == sp ) then
151  mtype = mpi_real
152  else
153  log_error("ATMOS_DYN_Tstep_short_fvm_hivi_setup",*) 'Unsupported precision'
154  call prc_abort
155  endif
156 
157  return
158  end subroutine atmos_dyn_tstep_short_fvm_hivi_setup
159 
160  !-----------------------------------------------------------------------------
161  subroutine atmos_dyn_tstep_short_fvm_hivi( &
162  DENS_RK, MOMZ_RK, MOMX_RK, MOMY_RK, RHOT_RK, &
163  PROG_RK, &
164  mflx_hi, tflx_hi, &
165  DENS0, MOMZ0, MOMX0, MOMY0, RHOT0, &
166  DENS, MOMZ, MOMX, MOMY, RHOT, &
167  DENS_t, MOMZ_t, MOMX_t, MOMY_t, RHOT_t, &
168  PROG0, PROG, &
169  DPRES0, RT2P, CORIOLI, &
170  num_diff, wdamp_coef, divdmp_coef, DDIV, &
171  FLAG_FCT_MOMENTUM, FLAG_FCT_T, &
172  FLAG_FCT_ALONG_STREAM, &
173  CDZ, FDZ, FDX, FDY, &
174  RCDZ, RCDX, RCDY, RFDZ, RFDX, RFDY, &
175  PHI, GSQRT, J13G, J23G, J33G, MAPF, &
176  REF_dens, REF_rhot, &
177  BND_W, BND_E, BND_S, BND_N, &
178  dtrk, last )
179  use scale_atmos_grid_cartesc_index
180  use scale_const, only: &
181  grav => const_grav, &
182  p00 => const_pre00
183  use scale_comm_cartesc, only: &
184  comm_vars8, &
185  comm_wait
186  use scale_atmos_dyn_common, only: &
187  atmos_dyn_fct
188  use scale_atmos_dyn_fvm_flux, only: &
189  atmos_dyn_fvm_flux_valuew_z, &
190  atmos_dyn_fvm_fluxz_xyz, &
191  atmos_dyn_fvm_fluxx_xyz, &
192  atmos_dyn_fvm_fluxy_xyz, &
193  atmos_dyn_fvm_fluxz_xyw, &
194  atmos_dyn_fvm_fluxj13_xyw, &
195  atmos_dyn_fvm_fluxj23_xyw, &
196  atmos_dyn_fvm_fluxx_xyw, &
197  atmos_dyn_fvm_fluxy_xyw, &
198  atmos_dyn_fvm_fluxz_uyz, &
199  atmos_dyn_fvm_fluxj13_uyz, &
200  atmos_dyn_fvm_fluxj23_uyz, &
201  atmos_dyn_fvm_fluxx_uyz, &
202  atmos_dyn_fvm_fluxy_uyz, &
203  atmos_dyn_fvm_fluxz_xvz, &
204  atmos_dyn_fvm_fluxj13_xvz, &
205  atmos_dyn_fvm_fluxj23_xvz, &
206  atmos_dyn_fvm_fluxx_xvz, &
207  atmos_dyn_fvm_fluxy_xvz
208  implicit none
209 
210  real(RP), intent(out) :: DENS_RK(ka,ia,ja) ! prognostic variables
211  real(RP), intent(out) :: MOMZ_RK(ka,ia,ja) !
212  real(RP), intent(out) :: MOMX_RK(ka,ia,ja) !
213  real(RP), intent(out) :: MOMY_RK(ka,ia,ja) !
214  real(RP), intent(out) :: RHOT_RK(ka,ia,ja) !
215 
216  real(RP), intent(out) :: PROG_RK(ka,ia,ja,va) !
217 
218  real(RP), intent(inout) :: mflx_hi(ka,ia,ja,3) ! rho * vel(x,y,z)
219  real(RP), intent(out) :: tflx_hi(ka,ia,ja,3) ! rho * theta * vel(x,y,z)
220 
221  real(RP), intent(in),target :: DENS0(ka,ia,ja) ! prognostic variables
222  real(RP), intent(in),target :: MOMZ0(ka,ia,ja) ! at previous dynamical time step
223  real(RP), intent(in),target :: MOMX0(ka,ia,ja) !
224  real(RP), intent(in),target :: MOMY0(ka,ia,ja) !
225  real(RP), intent(in),target :: RHOT0(ka,ia,ja) !
226 
227  real(RP), intent(in) :: DENS(ka,ia,ja) ! prognostic variables
228  real(RP), intent(in) :: MOMZ(ka,ia,ja) ! at previous RK step
229  real(RP), intent(in) :: MOMX(ka,ia,ja) !
230  real(RP), intent(in) :: MOMY(ka,ia,ja) !
231  real(RP), intent(in) :: RHOT(ka,ia,ja) !
232 
233  real(RP), intent(in) :: DENS_t(ka,ia,ja)
234  real(RP), intent(in) :: MOMZ_t(ka,ia,ja)
235  real(RP), intent(in) :: MOMX_t(ka,ia,ja)
236  real(RP), intent(in) :: MOMY_t(ka,ia,ja)
237  real(RP), intent(in) :: RHOT_t(ka,ia,ja)
238 
239  real(RP), intent(in) :: PROG0(ka,ia,ja,va)
240  real(RP), intent(in) :: PROG (ka,ia,ja,va)
241 
242  real(RP), intent(in) :: DPRES0(ka,ia,ja)
243  real(RP), intent(in) :: RT2P(ka,ia,ja)
244  real(RP), intent(in) :: CORIOLI(1,ia,ja)
245  real(RP), intent(in) :: num_diff(ka,ia,ja,5,3)
246  real(RP), intent(in) :: wdamp_coef(ka)
247  real(RP), intent(in) :: divdmp_coef
248  real(RP), intent(in) :: DDIV(ka,ia,ja)
249 
250  logical, intent(in) :: FLAG_FCT_MOMENTUM
251  logical, intent(in) :: FLAG_FCT_T
252  logical, intent(in) :: FLAG_FCT_ALONG_STREAM
253 
254  real(RP), intent(in) :: CDZ(ka)
255  real(RP), intent(in) :: FDZ(ka-1)
256  real(RP), intent(in) :: FDX(ia-1)
257  real(RP), intent(in) :: FDY(ja-1)
258  real(RP), intent(in) :: RCDZ(ka)
259  real(RP), intent(in) :: RCDX(ia)
260  real(RP), intent(in) :: RCDY(ja)
261  real(RP), intent(in) :: RFDZ(ka-1)
262  real(RP), intent(in) :: RFDX(ia-1)
263  real(RP), intent(in) :: RFDY(ja-1)
264 
265  real(RP), intent(in) :: PHI (ka,ia,ja)
266  real(RP), intent(in) :: GSQRT (ka,ia,ja,7)
267  real(RP), intent(in) :: J13G (ka,ia,ja,7)
268  real(RP), intent(in) :: J23G (ka,ia,ja,7)
269  real(RP), intent(in) :: J33G
270  real(RP), intent(in) :: MAPF (ia,ja,2,4)
271  real(RP), intent(in) :: REF_dens(ka,ia,ja)
272  real(RP), intent(in) :: REF_rhot(ka,ia,ja)
273 
274  logical, intent(in) :: BND_W
275  logical, intent(in) :: BND_E
276  logical, intent(in) :: BND_S
277  logical, intent(in) :: BND_N
278 
279  real(RP), intent(in) :: dtrk
280  logical, intent(in) :: last
281 
282 
283  ! diagnostic variables (work space)
284  real(RP) :: POTT(ka,ia,ja) ! potential temperature [K]
285  real(RP) :: DDENS(ka,ia,ja) ! density deviation from reference
286  real(RP) :: DPRES(ka,ia,ja) ! pressure deviation from reference
287  real(RP) :: DPRES_N(ka,ia,ja) ! pressure deviation at t=n+1 [Pa]
288 
289  real(RP) :: Sr(ka,ia,ja)
290  real(RP) :: Sw(ka,ia,ja)
291  real(RP) :: Su(ka,ia,ja)
292  real(RP) :: Sv(ka,ia,ja)
293  real(RP) :: St(ka,ia,ja)
294  real(RP) :: RCs2T(ka,ia,ja)
295  real(RP) :: B(ka,ia,ja)
296 
297  real(RP) :: duvw
298 
299  real(RP) :: qflx_hi (ka,ia,ja,3)
300  real(RP) :: qflx_J13(ka,ia,ja)
301  real(RP) :: qflx_J23(ka,ia,ja)
302  real(RP) :: mflx_hi2(ka,ia,ja,3)
303 
304  ! for implicit solver
305  real(RP) :: M(7,ka,ia,ja)
306  real(RP) :: r(ka,ia,ja)
307  real(RP) :: p(ka,ia,ja)
308 
309  real(RP) :: zero(ka,ia,ja)
310 
311  real(RP) :: rdt
312 
313  integer :: IIS, IIE, JJS, JJE
314  integer :: k, i, j
315 
316  zero = 0.0_rp
317 
318 #ifdef DEBUG
319  dpres(:,:,:) = undef
320  pott(:,:,:) = undef
321 
322  mflx_hi(:,:,:,:) = undef
323  mflx_hi(ks-1,:,:,zdir) = 0.0_rp
324 
325  qflx_hi(:,:,:,:) = undef
326  tflx_hi(:,:,:,:) = undef
327  qflx_j13(:,:,:) = undef
328  qflx_j23(:,:,:) = undef
329  mflx_hi2(:,:,:,:) = undef
330 
331  sr(:,:,:) = undef
332  sw(:,:,:) = undef
333  su(:,:,:) = undef
334  sv(:,:,:) = undef
335  st(:,:,:) = undef
336  rcs2t(:,:,:) = undef
337 
338  b(:,:,:) = undef
339  r(:,:,:) = undef
340  p(:,:,:) = undef
341 #endif
342 
343 #if defined DEBUG || defined QUICKDEBUG
344  dens_rk( 1:ks-1,:,:) = undef
345  dens_rk(ke+1:ka ,:,:) = undef
346  momz_rk( 1:ks-1,:,:) = undef
347  momz_rk(ke+1:ka ,:,:) = undef
348  momx_rk( 1:ks-1,:,:) = undef
349  momx_rk(ke+1:ka ,:,:) = undef
350  momy_rk( 1:ks-1,:,:) = undef
351  momy_rk(ke+1:ka ,:,:) = undef
352  rhot_rk( 1:ks-1,:,:) = undef
353  rhot_rk(ke+1:ka ,:,:) = undef
354  prog_rk( 1:ks-1,:,:,:) = undef
355  prog_rk(ke+1:ka ,:,:,:) = undef
356 #endif
357 
358  rdt = 1.0_rp / dtrk
359 
360  do jjs = js, je, jblock
361  jje = jjs+jblock-1
362  do iis = is, ie, iblock
363  iie = iis+iblock-1
364 
365  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
366  do j = jjs-1, jje+1
367  do i = iis-1, iie+1
368  do k = ks, ke
369 #ifdef DEBUG
370  call check( __line__, dpres0(k,i,j) )
371  call check( __line__, rt2p(k,i,j) )
372  call check( __line__, rhot(k,i,j) )
373  call check( __line__, ref_rhot(k,i,j) )
374 #endif
375  dpres(k,i,j) = dpres0(k,i,j) + rt2p(k,i,j) * ( rhot(k,i,j) - ref_rhot(k,i,j) )
376  enddo
377  dpres(ks-1,i,j) = dpres0(ks-1,i,j) - dens(ks,i,j) * ( phi(ks-1,i,j) - phi(ks+1,i,j) )
378  dpres(ke+1,i,j) = dpres0(ke+1,i,j) - dens(ke,i,j) * ( phi(ke+1,i,j) - phi(ke-1,i,j) )
379  enddo
380  enddo
381 
382  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
383  do j = jjs-1, jje+1
384  do i = iis-1, iie+1
385  do k = ks, ke
386 #ifdef DEBUG
387  call check( __line__, dens(k,i,j) )
388  call check( __line__, ref_dens(k,i,j) )
389 #endif
390  ddens(k,i,j) = dens(k,i,j) - ref_dens(k,i,j)
391  enddo
392  enddo
393  enddo
394 
395  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
396  do j = jjs-2, jje+2
397  do i = iis-2, iie+2
398  do k = ks, ke
399 #ifdef DEBUG
400  call check( __line__, rhot(k,i,j) )
401  call check( __line__, dens(k,i,j) )
402 #endif
403  pott(k,i,j) = rhot(k,i,j) / dens(k,i,j)
404  enddo
405  enddo
406  enddo
407 #ifdef DEBUG
408  k = iundef; i = iundef; j = iundef
409 #endif
410 
411  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
412  do j = jjs, jje
413  do i = iis, iie
414  do k = ks, ke
415 #ifdef DEBUG
416  call check( __line__, rt2p(k,i,j) )
417 #endif
418  rcs2t(k,i,j) = 1.0_rp / rt2p(k,i,j)
419  enddo
420  enddo
421  enddo
422 #ifdef DEBUG
423  k = iundef; i = iundef; j = iundef
424 #endif
425 
426  !##### continuity equation #####
427 
428  ! at (x, y, w)
429  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
430  do j = jjs, jje
431  do i = iis-1, iie
432  do k = ks, ke-1
433 #ifdef DEBUG
434  call check( __line__, momx(k+1,i ,j) )
435  call check( __line__, momx(k+1,i-1,j) )
436  call check( __line__, momx(k ,i ,j) )
437  call check( __line__, momx(k ,i+1,j) )
438  call check( __line__, momy(k+1,i,j) )
439  call check( __line__, momy(k+1,i,j-1) )
440  call check( __line__, momy(k ,i,j) )
441  call check( __line__, momy(k ,i,j-1) )
442  call check( __line__, num_diff(k,i,j,i_dens,zdir) )
443  call check( __line__, gsqrt(k,i,j,i_xyw) )
444 #endif
445  mflx_hi2(k,i,j,zdir) = j13g(k,i,j,i_xyw) * 0.25_rp * ( momx(k+1,i,j)+momx(k+1,i-1,j) &
446  + momx(k ,i,j)+momx(k ,i-1,j) ) &
447  + j23g(k,i,j,i_xyw) * 0.25_rp * ( momy(k+1,i,j)+momy(k+1,i,j-1) &
448  + momy(k ,i,j)+momy(k ,i,j-1) ) &
449  + gsqrt(k,i,j,i_xyw) * num_diff(k,i,j,i_dens,zdir)
450  enddo
451  enddo
452  enddo
453 #ifdef DEBUG
454  k = iundef; i = iundef; j = iundef
455 #endif
456 
457  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
458  do j = jjs, jje
459  do i = iis, iie
460  mflx_hi2(ks-1,i,j,zdir) = 0.0_rp
461  mflx_hi2(ke ,i,j,zdir) = 0.0_rp
462  enddo
463  enddo
464 
465  ! at (u, y, z)
466  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
467  do j = jjs , jje
468  do i = iis-1, iie
469  do k = ks, ke
470 #ifdef DEBUG
471  call check( __line__, gsqrt(k,i,j,i_uyz) )
472  call check( __line__, num_diff(k,i,j,i_dens,xdir) )
473 #endif
474  mflx_hi2(k,i,j,xdir) = gsqrt(k,i,j,i_uyz) * num_diff(k,i,j,i_dens,xdir)
475  enddo
476  enddo
477  enddo
478 #ifdef DEBUG
479  k = iundef; i = iundef; j = iundef
480 #endif
481  ! at (x, v, z)
482  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
483  do j = jjs-1, jje
484  do i = iis , iie
485  do k = ks, ke
486 #ifdef DEBUG
487  call check( __line__, gsqrt(k,i,j,i_xvz) )
488  call check( __line__, num_diff(k,i,j,i_dens,ydir) )
489 #endif
490  mflx_hi2(k,i,j,ydir) = gsqrt(k,i,j,i_xvz) * num_diff(k,i,j,i_dens,ydir)
491  enddo
492  enddo
493  enddo
494 #ifdef DEBUG
495  k = iundef; i = iundef; j = iundef
496 #endif
497 
498  !##### momentum equation (z) #####
499 
500  ! at (x, y, z)
501  ! note than z-index is added by -1
502  call atmos_dyn_fvm_fluxz_xyw( qflx_hi(:,:,:,zdir), & ! (out)
503  momz, momz, dens, & ! (in)
504  gsqrt(:,:,:,i_xyz), j33g, & ! (in)
505  num_diff(:,:,:,i_momz,zdir), & ! (in)
506  cdz, fdz, dtrk, &
507  iis, iie, jjs, jje ) ! (in)
508  call atmos_dyn_fvm_fluxj13_xyw( qflx_j13, & ! (out)
509  momx, momz, dens, & ! (in)
510  gsqrt(:,:,:,i_xyz), j13g(:,:,:,i_xyz), mapf(:,:,:,i_xy), & ! (in)
511  cdz, &
512  iis, iie, jjs, jje ) ! (in)
513  call atmos_dyn_fvm_fluxj23_xyw( qflx_j23, & ! (out)
514  momy, momz, dens, & ! (in)
515  gsqrt(:,:,:,i_xyz), j23g(:,:,:,i_xyz), mapf(:,:,:,i_xy), & ! (in)
516  cdz, &
517  iis, iie, jjs, jje ) ! (in)
518 
519  call atmos_dyn_fvm_fluxx_xyw( qflx_hi(:,:,:,xdir), & ! (out)
520  momx, momz, dens, & ! (in)
521  gsqrt(:,:,:,i_uyw), mapf(:,:,:,i_uy), & ! (in)
522  num_diff(:,:,:,i_momz,xdir), & ! (in)
523  cdz, & ! (in)
524  iis, iie, jjs, jje ) ! (in)
525 
526  ! at (x, v, w)
527  call atmos_dyn_fvm_fluxy_xyw( qflx_hi(:,:,:,ydir), & ! (out)
528  momy, momz, dens, & ! (in)
529  gsqrt(:,:,:,i_xvw), mapf(:,:,:,i_xv), & ! (in)
530  num_diff(:,:,:,i_momz,ydir), & ! (in)
531  cdz, & ! (in)
532  iis, iie, jjs, jje ) ! (in)
533 
534  !--- update momentum(z)
535  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
536  do j = jjs, jje
537  do i = iis, iie
538  do k = ks, ke-1
539 #ifdef DEBUG
540  call check( __line__, qflx_hi(k ,i ,j ,zdir) )
541  call check( __line__, qflx_hi(k-1,i ,j ,zdir) )
542  call check( __line__, qflx_j13(k ,i ,j ) )
543  call check( __line__, qflx_j13(k-1,i ,j ) )
544  call check( __line__, qflx_j23(k ,i ,j ) )
545  call check( __line__, qflx_j23(k-1,i ,j ) )
546  call check( __line__, qflx_hi(k ,i ,j ,xdir) )
547  call check( __line__, qflx_hi(k ,i-1,j ,xdir) )
548  call check( __line__, qflx_hi(k ,i ,j ,ydir) )
549  call check( __line__, qflx_hi(k ,i ,j-1,ydir) )
550  call check( __line__, ddiv(k ,i,j) )
551  call check( __line__, ddiv(k+1,i,j) )
552  call check( __line__, momz0(k,i,j) )
553  call check( __line__, momz_t(k,i,j) )
554 #endif
555  sw(k,i,j) = &
556  - ( ( qflx_hi(k,i,j,zdir) - qflx_hi(k-1,i ,j ,zdir) &
557  + qflx_j13(k,i,j) - qflx_j13(k-1,i ,j ) &
558  + qflx_j23(k,i,j) - qflx_j23(k-1,i ,j ) ) * rfdz(k) &
559  + ( qflx_hi(k,i,j,xdir) - qflx_hi(k ,i-1,j ,xdir) ) * rcdx(i) &
560  + ( qflx_hi(k,i,j,ydir) - qflx_hi(k ,i ,j-1,ydir) ) * rcdy(j) &
561  ) / gsqrt(k,i,j,i_xyw) &
562  - wdamp_coef(k) * momz0(k,i,j) & ! Rayleigh damping
563  + divdmp_coef * rdt * ( ddiv(k+1,i,j)-ddiv(k,i,j) ) * fdz(k) & ! divergence damping
564  + momz_t(k,i,j)
565  enddo
566  enddo
567  enddo
568 #ifdef DEBUG
569  k = iundef; i = iundef; j = iundef
570 #endif
571  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
572  do j = jjs, jje
573  do i = iis, iie
574  sw(ks-1,i,j) = 0.0_rp
575  sw(ke ,i,j) = 0.0_rp
576  enddo
577  enddo
578 #ifdef DEBUG
579  k = iundef; i = iundef; j = iundef
580 #endif
581 
582 
583  !##### momentum equation (x) #####
584 
585  ! at (u, y, w)
586  call atmos_dyn_fvm_fluxz_uyz( qflx_hi(:,:,:,zdir), & ! (out)
587  momz, momx, dens, & ! (in)
588  gsqrt(:,:,:,i_uyw), j33g, & ! (in)
589  num_diff(:,:,:,i_momx,zdir), & ! (in)
590  cdz, & ! (in)
591  iis, iie, jjs, jje ) ! (in)
592  call atmos_dyn_fvm_fluxj13_uyz( qflx_j13, & ! (out)
593  momx, momx, dens, & ! (in)
594  gsqrt(:,:,:,i_uyz), j13g(:,:,:,i_uyw), mapf(:,:,:,i_uy), & ! (in)
595  cdz, & ! (in)
596  iis, iie, jjs, jje ) ! (in)
597  call atmos_dyn_fvm_fluxj23_uyz( qflx_j23, & ! (out)
598  momy, momx, dens, & ! (in)
599  gsqrt(:,:,:,i_uyz), j23g(:,:,:,i_uyw), mapf(:,:,:,i_uy), & ! (in)
600  cdz, & ! (in)
601  iis, iie, jjs, jje ) ! (in)
602 
603  ! at (x, y, z)
604  ! note that x-index is added by -1
605  call atmos_dyn_fvm_fluxx_uyz( qflx_hi(:,:,:,xdir), & ! (out)
606  momx, momx, dens, & ! (in)
607  gsqrt(:,:,:,i_xyz), mapf(:,:,:,i_xy), & ! (in)
608  num_diff(:,:,:,i_momx,xdir), & ! (in)
609  cdz, & ! (in)
610  iis, iie, jjs, jje ) ! (in)
611 
612  ! at (u, v, z)
613  call atmos_dyn_fvm_fluxy_uyz( qflx_hi(:,:,:,ydir), & ! (out)
614  momy, momx, dens, & ! (in)
615  gsqrt(:,:,:,i_uvz), mapf(:,:,:,i_uv), & ! (in)
616  num_diff(:,:,:,i_momx,ydir), & ! (in)
617  cdz, & ! (in)
618  iis, iie, jjs, jje ) ! (in)
619 
620  !--- update momentum(x)
621  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
622  do j = jjs, jje
623  do i = iis, iie
624  do k = ks, ke
625 #ifdef DEBUG
626  call check( __line__, qflx_hi(k ,i ,j ,zdir) )
627  call check( __line__, qflx_hi(k-1,i ,j ,zdir) )
628  call check( __line__, qflx_hi(k ,i ,j ,xdir) )
629  call check( __line__, qflx_hi(k ,i-1,j ,xdir) )
630  call check( __line__, qflx_hi(k ,i ,j ,ydir) )
631  call check( __line__, qflx_hi(k ,i ,j-1,ydir) )
632  call check( __line__, dpres(k+1,i+1,j) )
633  call check( __line__, dpres(k+1,i ,j) )
634  call check( __line__, dpres(k-1,i+1,j) )
635  call check( __line__, dpres(k-1,i ,j) )
636  call check( __line__, corioli(1,i ,j) )
637  call check( __line__, corioli(1,i+1,j) )
638  call check( __line__, momy(k,i ,j ) )
639  call check( __line__, momy(k,i+1,j ) )
640  call check( __line__, momy(k,i ,j-1) )
641  call check( __line__, momy(k,i+1,j-1) )
642  call check( __line__, ddiv(k,i+1,j) )
643  call check( __line__, ddiv(k,i ,j) )
644  call check( __line__, momx0(k,i,j) )
645 #endif
646  su(k,i,j) = ( &
647  - ( ( qflx_hi(k,i,j,zdir) - qflx_hi(k-1,i ,j ,zdir) &
648  + qflx_j13(k,i,j) - qflx_j13(k-1,i ,j) &
649  + qflx_j23(k,i,j) - qflx_j23(k-1,i ,j) ) * rcdz(k) &
650  + ( qflx_hi(k,i,j,xdir) - qflx_hi(k ,i-1,j ,xdir) ) * rfdx(i) &
651  + ( qflx_hi(k,i,j,ydir) - qflx_hi(k ,i ,j-1,ydir) ) * rcdy(j) ) &
652  - ( j13g(k+1,i,j,i_uyz) * ( dpres(k+1,i+1,j)+dpres(k+1,i,j) ) &
653  - j13g(k-1,i,j,i_uyz) * ( dpres(k-1,i+1,j)+dpres(k-1,i,j) ) ) &
654  * 0.5_rp / ( fdz(k+1)+fdz(k) ) &
655  ) / gsqrt(k,i,j,i_uyz) &
656  + 0.125_rp * ( corioli(1,i,j)+corioli(1,i+1,j) ) &
657  * ( momy(k,i,j)+momy(k,i+1,j)+momy(k,i,j-1)+momy(k,i+1,j-1) ) & ! coriolis force
658  + divdmp_coef * rdt * ( ddiv(k,i+1,j)-ddiv(k,i,j) ) * fdx(i) & ! divergence damping
659  + momx_t(k,i,j)
660  enddo
661  enddo
662  enddo
663 #ifdef DEBUG
664  k = iundef; i = iundef; j = iundef
665 #endif
666 
667  !##### momentum equation (y) #####
668 
669  ! at (x, v, w)
670  call atmos_dyn_fvm_fluxz_xvz( qflx_hi(:,:,:,zdir), & ! (out)
671  momz, momy, dens, & ! (in)
672  gsqrt(:,:,:,i_xvw), j33g, & ! (in)
673  num_diff(:,:,:,i_momy,zdir), & ! (in)
674  cdz, & ! (in)
675  iis, iie, jjs, jje ) ! (in)
676  call atmos_dyn_fvm_fluxj13_xvz( qflx_j13, & ! (out)
677  momx, momy, dens, & ! (in)
678  gsqrt(:,:,:,i_xvz), j13g(:,:,:,i_xvw), mapf(:,:,:,i_xv), & ! (in)
679  cdz, & ! (in)
680  iis, iie, jjs, jje ) ! (in)
681  call atmos_dyn_fvm_fluxj23_xvz( qflx_j23, & ! (out)
682  momy, momy, dens, & ! (in)
683  gsqrt(:,:,:,i_xvz), j23g(:,:,:,i_xvw), mapf(:,:,:,i_xv), & ! (in)
684  cdz, & ! (in)
685  iis, iie, jjs, jje ) ! (in)
686 
687  ! at (u, v, z)
688  call atmos_dyn_fvm_fluxx_xvz( qflx_hi(:,:,:,xdir), & ! (out)
689  momx, momy, dens, & ! (in)
690  gsqrt(:,:,:,i_uvz), mapf(:,:,:,i_uv), & ! (in)
691  num_diff(:,:,:,i_momy,xdir), & ! (in)
692  cdz, & ! (in)
693  iis, iie, jjs, jje ) ! (in)
694 
695  ! at (x, y, z)
696  ! note that y-index is added by -1
697  call atmos_dyn_fvm_fluxy_xvz( qflx_hi(:,:,:,ydir), & ! (out)
698  momy, momy, dens, & ! (in)
699  gsqrt(:,:,:,i_xyz), mapf(:,:,:,i_xy), & ! (in)
700  num_diff(:,:,:,i_momy,ydir), & ! (in
701  cdz, & ! (in)
702  iis, iie, jjs, jje ) ! (in)
703 
704  !--- update momentum(y)
705  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
706  do j = jjs, jje
707  do i = iis, iie
708  do k = ks, ke
709 #ifdef DEBUG
710  call check( __line__, qflx_hi(k ,i ,j ,zdir) )
711  call check( __line__, qflx_hi(k-1,i ,j ,zdir) )
712  call check( __line__, qflx_hi(k ,i ,j ,xdir) )
713  call check( __line__, qflx_hi(k ,i-1,j ,xdir) )
714  call check( __line__, qflx_hi(k ,i ,j ,ydir) )
715  call check( __line__, qflx_hi(k ,i ,j-1,ydir) )
716  call check( __line__, dpres(k+1,i,j ) )
717  call check( __line__, dpres(k+1,i,j+1) )
718  call check( __line__, dpres(k-1,i,j ) )
719  call check( __line__, dpres(k-1,i,j+1) )
720  call check( __line__, corioli(1,i,j ) )
721  call check( __line__, corioli(1,i,j+1) )
722  call check( __line__, momx(k,i ,j ) )
723  call check( __line__, momx(k,i ,j+1) )
724  call check( __line__, momx(k,i-1,j ) )
725  call check( __line__, momx(k,i-1,j+1) )
726  call check( __line__, ddiv(k,i,j+1) )
727  call check( __line__, ddiv(k,i,j ) )
728  call check( __line__, momy_t(k,i,j) )
729 #endif
730  sv(k,i,j) = &
731  ( - ( ( qflx_hi(k,i,j,zdir) - qflx_hi(k-1,i ,j ,zdir) &
732  + qflx_j13(k,i,j) - qflx_j13(k-1,i ,j ) &
733  + qflx_j23(k,i,j) - qflx_j23(k-1,i ,j ) ) * rcdz(k) &
734  + ( qflx_hi(k,i,j,xdir) - qflx_hi(k ,i-1,j ,xdir) ) * rcdx(i) &
735  + ( qflx_hi(k,i,j,ydir) - qflx_hi(k ,i ,j-1,ydir) ) * rfdy(j) ) &
736  - ( j23g(k+1,i,j,i_xvz) * ( dpres(k+1,i,j+1)+dpres(k+1,i,j) ) &
737  - j23g(k-1,i,j,i_xvz) * ( dpres(k-1,i,j+1)+dpres(k-1,i,j) ) ) &
738  * 0.5_rp / ( fdz(k+1)+fdz(k) ) & ! pressure gradient force
739  ) / gsqrt(k,i,j,i_xvz) &
740  - 0.125_rp * ( corioli(1,i ,j+1)+corioli(1,i ,j) ) &
741  * ( momx(k,i,j+1)+momx(k,i,j)+momx(k,i-1,j+1)+momx(k,i-1,j) ) & ! coriolis force
742  + divdmp_coef * rdt * ( ddiv(k,i,j+1)-ddiv(k,i,j) ) * fdy(j) & ! divergence damping
743  + momy_t(k,i,j)
744  enddo
745  enddo
746  enddo
747 #ifdef DEBUG
748  k = iundef; i = iundef; j = iundef
749 #endif
750 
751  !##### Thermodynamic Equation #####
752 
753 
754  ! at (x, y, w)
755  call atmos_dyn_fvm_fluxz_xyz( tflx_hi(:,:,:,zdir), & ! (out)
756  mflx_hi2(:,:,:,zdir), pott, gsqrt(:,:,:,i_xyw), & ! (in)
757  num_diff(:,:,:,i_rhot,zdir), & ! (in)
758  cdz, & ! (in)
759  iis, iie, jjs, jje ) ! (in)
760 
761  ! at (u, y, z)
762  call atmos_dyn_fvm_fluxx_xyz( tflx_hi(:,:,:,xdir), & ! (out)
763  mflx_hi2(:,:,:,xdir), pott, gsqrt(:,:,:,i_uyz), & ! (in)
764  num_diff(:,:,:,i_rhot,xdir), & ! (in)
765  cdz, & ! (in)
766  iis, iie, jjs, jje ) ! (in)
767 
768  ! at (x, v, z)
769  call atmos_dyn_fvm_fluxy_xyz( tflx_hi(:,:,:,ydir), & ! (out)
770  mflx_hi2(:,:,:,ydir), pott, gsqrt(:,:,:,i_xvz), & ! (in)
771  num_diff(:,:,:,i_rhot,ydir), & ! (in)
772  cdz, & ! (in)
773  iis, iie, jjs, jje ) ! (in)
774 
775  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
776  do j = jjs, jje
777  do i = iis, iie
778  do k = ks, ke
779 #ifdef DEBUG
780  call check( __line__, tflx_hi(k ,i ,j ,zdir) )
781  call check( __line__, tflx_hi(k-1,i ,j ,zdir) )
782  call check( __line__, tflx_hi(k ,i ,j ,xdir) )
783  call check( __line__, tflx_hi(k ,i-1,j ,xdir) )
784  call check( __line__, tflx_hi(k ,i ,j ,ydir) )
785  call check( __line__, tflx_hi(k ,i ,j-1,ydir) )
786  call check( __line__, rhot_t(k,i,j) )
787 #endif
788  st(k,i,j) = &
789  - ( ( tflx_hi(k,i,j,zdir) - tflx_hi(k-1,i ,j ,zdir) ) * rcdz(k) &
790  + ( tflx_hi(k,i,j,xdir) - tflx_hi(k ,i-1,j ,xdir) ) * rcdx(i) &
791  + ( tflx_hi(k,i,j,ydir) - tflx_hi(k ,i ,j-1,ydir) ) * rcdy(j) &
792  ) / gsqrt(k,i,j,i_xyz) &
793  + rhot_t(k,i,j)
794  enddo
795  enddo
796  enddo
797 #ifdef DEBUG
798  k = iundef; i = iundef; j = iundef
799 #endif
800 
801  enddo
802  enddo
803 
804  ! implicit solver
805 #ifdef HIVI_BICGSTAB
806 
807 #ifdef DEBUG
808  m(:,:,:,:) = undef
809 #endif
810 
811  call comm_vars8( su, 1 )
812  call comm_vars8( sv, 2 )
813  call comm_wait ( su, 1 )
814  call comm_wait ( sv, 2 )
815 
816  do jjs = js, je, jblock
817  jje = jjs+jblock-1
818  do iis = is, ie, iblock
819  iie = iis+iblock-1
820  ! r = b - M x0
821  do j = jjs, jje
822  do i = iis, iie
823  do k = ks+2, ke-2
824 #ifdef DEBUG
825  call check( __line__, momz(k-1,i,j) )
826  call check( __line__, momz(k ,i,j) )
827  call check( __line__, momx(k,i-1,j) )
828  call check( __line__, momx(k,i ,j) )
829  call check( __line__, momy(k,i,j-1) )
830  call check( __line__, momy(k,i,j ) )
831  call check( __line__, pott(k-2,i,j) )
832  call check( __line__, pott(k-1,i,j) )
833  call check( __line__, pott(k ,i,j) )
834  call check( __line__, pott(k+1,i,j) )
835  call check( __line__, pott(k+2,i,j) )
836  call check( __line__, pott(k,i-2,j) )
837  call check( __line__, pott(k,i-1,j) )
838  call check( __line__, pott(k,i ,j) )
839  call check( __line__, pott(k,i+1,j) )
840  call check( __line__, pott(k,i+2,j) )
841  call check( __line__, pott(k,i,j-2) )
842  call check( __line__, pott(k,i,j-1) )
843  call check( __line__, pott(k,i,j ) )
844  call check( __line__, pott(k,i,j+1) )
845  call check( __line__, pott(k,i,j+2) )
846  call check( __line__, sw(k-1,i,j) )
847  call check( __line__, sw(k ,i,j) )
848  call check( __line__, su(k,i-1,j) )
849  call check( __line__, su(k,i ,j) )
850  call check( __line__, sv(k,i,j-1) )
851  call check( __line__, sv(k,i,j ) )
852  call check( __line__, gsqrt(k,i ,j,i_uyz) )
853  call check( __line__, gsqrt(k,i-1,j,i_uyz) )
854  call check( __line__, gsqrt(k,i,j-1,i_xvz) )
855  call check( __line__, gsqrt(k,i,j ,i_xvz) )
856  call check( __line__, gsqrt(k,i,j,i_xyz) )
857  call check( __line__, st(k,i,j) )
858  call check( __line__, dpres(k-1,i,j) )
859  call check( __line__, dpres(k ,i,j) )
860  call check( __line__, dpres(k+1,i,j) )
861  call check( __line__, rt2p(k-1,i,j) )
862  call check( __line__, rt2p(k ,i,j) )
863  call check( __line__, rt2p(k+1,i,j) )
864  call check( __line__, rhot(k-1,i,j) )
865  call check( __line__, rhot(k+1,i,j) )
866  call check( __line__, ddens(k-1,i,j) )
867  call check( __line__, ddens(k+1,i,j) )
868 #endif
869  b(k,i,j) = ( &
870  ( j33g * ( momz(k ,i,j) + dtrk*sw(k ,i,j) ) &
871  * ( fact_n*(pott(k+1,i,j)+pott(k ,i,j)) &
872  + fact_f*(pott(k+2,i,j)+pott(k-1,i,j)) ) &
873  - j33g * ( momz(k-1,i,j) + dtrk*sw(k-1,i,j) ) &
874  * ( fact_n*(pott(k ,i,j)+pott(k-1,i,j)) &
875  + fact_f*(pott(k+1,i,j)+pott(k-2,i,j)) ) ) * rcdz(k) &
876  + ( gsqrt(k,i ,j,i_uyz) * ( momx(k,i ,j) + dtrk*su(k,i ,j) ) &
877  * ( fact_n*(pott(k,i+1,j)+pott(k,i ,j)) &
878  + fact_f*(pott(k,i+2,j)+pott(k,i-1,j)) ) &
879  - gsqrt(k,i-1,j,i_uyz) * ( momx(k,i-1,j) + dtrk*su(k,i-1,j) ) &
880  * ( fact_n*(pott(k,i ,j)+pott(k,i-1,j)) &
881  + fact_f*(pott(k,i+1,j)+pott(k,i-2,j)) ) ) * rcdx(i) &
882  + ( gsqrt(k,i,j ,i_xvz) * ( momy(k,i,j ) + dtrk*sv(k,i,j ) ) &
883  * ( fact_n*(pott(k,i,j+1)+pott(k,i,j )) &
884  + fact_f*(pott(k,i,j+2)+pott(k,i,j-1)) ) &
885  - gsqrt(k,i,j-1,i_xvz) * ( momy(k,i,j-1) + dtrk*sv(k,i,j-1) ) &
886  * ( fact_n*(pott(k,i,j )+pott(k,i,j-1)) &
887  + fact_f*(pott(k,i,j+1)+pott(k,i,j-2)) ) ) * rcdy(j) &
888  + gsqrt(k,i,j,i_xyz) * ( st(k,i,j) - dpres(k,i,j) * rcs2t(k,i,j) * rdt ) &
889  ) * rdt &
890  + grav * j33g * ( ( dpres(k+1,i,j)*rcs2t(k+1,i,j) &
891  - dpres(k-1,i,j)*rcs2t(k-1,i,j) ) &
892  - ( rhot(k+1,i,j)*ddens(k+1,i,j)/dens(k+1,i,j) &
893  - rhot(k-1,i,j)*ddens(k-1,i,j)/dens(k-1,i,j) ) &
894  ) / ( fdz(k) + fdz(k-1) )
895  enddo
896  enddo
897  enddo
898 #ifdef DEBUG
899  k = iundef; i = iundef; j = iundef
900 #endif
901  do j = jjs, jje
902  do i = iis, iie
903 #ifdef DEBUG
904  call check( __line__, momz(ks,i,j) )
905  call check( __line__, momx(ks,i,j) )
906  call check( __line__, momy(ks,i,j) )
907  call check( __line__, sw(ks,i,j) )
908  call check( __line__, su(ks,i-1,j) )
909  call check( __line__, su(ks,i ,j) )
910  call check( __line__, sv(ks,i,j-1) )
911  call check( __line__, sv(ks,i,j ) )
912  call check( __line__, st(ks,i,j) )
913  call check( __line__, pott(ks ,i,j) )
914  call check( __line__, pott(ks+1,i,j) )
915  call check( __line__, pott(ks,i-2,j) )
916  call check( __line__, pott(ks,i-1,j) )
917  call check( __line__, pott(ks,i ,j) )
918  call check( __line__, pott(ks,i+1,j) )
919  call check( __line__, pott(ks,i+2,j) )
920  call check( __line__, pott(ks,i,j-2) )
921  call check( __line__, pott(ks,i,j-1) )
922  call check( __line__, pott(ks,i,j ) )
923  call check( __line__, pott(ks,i,j+1) )
924  call check( __line__, pott(ks,i,j+2) )
925  call check( __line__, gsqrt(ks,i,j,i_uyz) )
926  call check( __line__, gsqrt(ks,i,j,i_xvz) )
927  call check( __line__, gsqrt(ks,i,j,i_xyz) )
928  call check( __line__, dpres(ks ,i,j) )
929  call check( __line__, dpres(ks+1,i,j) )
930  call check( __line__, rt2p(ks ,i,j) )
931  call check( __line__, rt2p(ks+1,i,j) )
932  call check( __line__, rhot(ks+1,i,j) )
933  call check( __line__, ddens(ks+1,i,j) )
934 #endif
935  b(ks,i,j) = ( &
936  ( j33g * ( momz(ks ,i,j) + dtrk*sw(ks ,i,j) ) &
937  * 0.5_rp*(pott(ks+1,i,j)+pott(ks ,i,j)) ) * rcdz(ks) &
938  + ( gsqrt(ks,i ,j,i_uyz) * ( momx(ks,i ,j) + dtrk*su(ks,i ,j) ) &
939  * ( fact_n*(pott(ks,i+1,j)+pott(ks,i ,j)) &
940  + fact_f*(pott(ks,i+2,j)+pott(ks,i-1,j)) ) &
941  - gsqrt(ks,i-1,j,i_uyz) * ( momx(ks,i-1,j) + dtrk*su(ks,i-1,j) ) &
942  * ( fact_n*(pott(ks,i ,j)+pott(ks,i-1,j)) &
943  + fact_f*(pott(ks,i+1,j)+pott(ks,i-2,j)) ) ) * rcdx(i) &
944  + ( gsqrt(ks,i,j ,i_xvz) * ( momy(ks,i,j ) + dtrk*sv(ks,i,j ) ) &
945  * ( fact_n*(pott(ks,i,j+1)+pott(ks,i,j )) &
946  + fact_f*(pott(ks,i,j+2)+pott(ks,i,j-1)) ) &
947  - gsqrt(ks,i,j-1,i_xvz) * ( momy(ks,i,j-1) + dtrk*sv(ks,i,j-1) ) &
948  * ( fact_n*(pott(ks,i,j )+pott(ks,i,j-1)) &
949  + fact_f*(pott(ks,i,j+1)+pott(ks,i,j-2)) ) ) * rcdy(j) &
950  + gsqrt(ks,i,j,i_xyz) * ( st(ks,i,j) - dpres(ks,i,j) * rcs2t(ks,i,j) * rdt ) &
951  ) * rdt &
952  + grav * j33g * 0.5_rp * ( ( dpres(ks,i,j)+dpres(ks+1,i,j) ) * rcs2t(ks,i,j) &
953  - ( ddens(ks,i,j)+ddens(ks+1,i,j) ) * pott(ks,i,j) ) * rcdz(ks)
954 #ifdef DEBUG
955  call check( __line__, momz(ks ,i,j) )
956  call check( __line__, momz(ks+1 ,i,j) )
957  call check( __line__, momx(ks+1,i-1,j) )
958  call check( __line__, momx(ks+1,i ,j) )
959  call check( __line__, momy(ks+1,i,j-1) )
960  call check( __line__, momy(ks+1,i,j ) )
961  call check( __line__, pott(ks ,i,j) )
962  call check( __line__, pott(ks+1 ,i,j) )
963  call check( __line__, pott(ks+1+1,i,j) )
964  call check( __line__, pott(ks+1+2,i,j) )
965  call check( __line__, pott(ks+1,i-2,j) )
966  call check( __line__, pott(ks+1,i-1,j) )
967  call check( __line__, pott(ks+1,i ,j) )
968  call check( __line__, pott(ks+1,i+1,j) )
969  call check( __line__, pott(ks+1,i+2,j) )
970  call check( __line__, pott(ks+1,i,j-2) )
971  call check( __line__, pott(ks+1,i,j-1) )
972  call check( __line__, pott(ks+1,i,j ) )
973  call check( __line__, pott(ks+1,i,j+1) )
974  call check( __line__, pott(ks+1,i,j+2) )
975  call check( __line__, sw(ks ,i,j) )
976  call check( __line__, sw(ks+1,i,j) )
977  call check( __line__, su(ks+1,i-1,j) )
978  call check( __line__, su(ks+1,i ,j) )
979  call check( __line__, sv(ks+1,i,j-1) )
980  call check( __line__, sv(ks+1,i,j ) )
981  call check( __line__, gsqrt(ks+1,i ,j,i_uyz) )
982  call check( __line__, gsqrt(ks+1,i-1,j,i_uyz) )
983  call check( __line__, gsqrt(ks+1,i,j-1,i_xvz) )
984  call check( __line__, gsqrt(ks+1,i,j ,i_xvz) )
985  call check( __line__, gsqrt(ks+1,i,j,i_xyz) )
986  call check( __line__, st(ks+1,i,j) )
987  call check( __line__, dpres(ks ,i,j) )
988  call check( __line__, dpres(ks+1,i,j) )
989  call check( __line__, dpres(ks+2,i,j) )
990  call check( __line__, rt2p(ks ,i,j) )
991  call check( __line__, rt2p(ks+1,i,j) )
992  call check( __line__, rt2p(ks+2,i,j) )
993  call check( __line__, rhot(ks ,i,j) )
994  call check( __line__, rhot(ks+2,i,j) )
995  call check( __line__, ddens(ks ,i,j) )
996  call check( __line__, ddens(ks+2,i,j) )
997 #endif
998  b(ks+1,i,j) = ( &
999  ( j33g * ( momz(ks+1,i,j) + dtrk*sw(ks+1,i,j) ) &
1000  * ( fact_n*(pott(ks+2,i,j)+pott(ks+1,i,j)) &
1001  + fact_f*(pott(ks+3,i,j)+pott(ks ,i,j)) ) &
1002  - j33g * ( momz(ks+1-1,i,j) + dtrk*sw(ks+1-1,i,j) ) &
1003  * ( 0.5_rp*(pott(ks+1,i,j)+pott(ks,i,j)) ) ) * rcdz(ks+1) &
1004  + ( gsqrt(ks+1,i ,j,i_uyz) * ( momx(ks+1,i ,j) + dtrk*su(ks+1,i ,j) ) &
1005  * ( fact_n*(pott(ks+1,i+1,j)+pott(ks+1,i ,j)) &
1006  + fact_f*(pott(ks+1,i+2,j)+pott(ks+1,i-1,j)) ) &
1007  - gsqrt(ks+1,i-1,j,i_uyz) * ( momx(ks+1,i-1,j) + dtrk*su(ks+1,i-1,j) ) &
1008  * ( fact_n*(pott(ks+1,i ,j)+pott(ks+1,i-1,j)) &
1009  + fact_f*(pott(ks+1,i+1,j)+pott(ks+1,i-2,j)) ) ) * rcdx(i) &
1010  + ( gsqrt(ks+1,i,j ,i_xvz) * ( momy(ks+1,i,j ) + dtrk*sv(ks+1,i,j ) ) &
1011  * ( fact_n*(pott(ks+1,i,j+1)+pott(ks+1,i,j )) &
1012  + fact_f*(pott(ks+1,i,j+2)+pott(ks+1,i,j-1)) ) &
1013  - gsqrt(ks+1,i,j-1,i_xvz) * ( momy(ks+1,i,j-1) + dtrk*sv(ks+1,i,j-1) ) &
1014  * ( fact_n*(pott(ks+1,i,j )+pott(ks+1,i,j-1)) &
1015  + fact_f*(pott(ks+1,i,j+1)+pott(ks+1,i,j-2)) ) ) * rcdy(j) &
1016  + gsqrt(ks+1,i,j,i_xyz) * ( st(ks+1,i,j) - dpres(ks+1,i,j) * rcs2t(ks+1,i,j) * rdt ) &
1017  ) * rdt &
1018  + grav * j33g * ( ( dpres(ks+2,i,j)*rcs2t(ks+2,i,j) &
1019  - dpres(ks ,i,j)*rcs2t(ks ,i,j) ) &
1020  - ( rhot(ks+2,i,j)*ddens(ks+2,i,j)/dens(ks+2,i,j) &
1021  - rhot(ks ,i,j)*ddens(ks ,i,j)/dens(ks ,i,j) ) &
1022  ) / ( fdz(ks+1) + fdz(ks) )
1023 #ifdef DEBUG
1024  call check( __line__, momz(ke-2,i,j) )
1025  call check( __line__, momz(ke-1,i,j) )
1026  call check( __line__, momx(ke-1,i-1,j) )
1027  call check( __line__, momx(ke-1,i ,j) )
1028  call check( __line__, momy(ke-1,i,j-1) )
1029  call check( __line__, momy(ke-1,i,j ) )
1030  call check( __line__, pott(ke-3,i,j) )
1031  call check( __line__, pott(ke-2,i,j) )
1032  call check( __line__, pott(ke-1,i,j) )
1033  call check( __line__, pott(ke ,i,j) )
1034  call check( __line__, pott(ke-1,i-2,j) )
1035  call check( __line__, pott(ke-1,i-1,j) )
1036  call check( __line__, pott(ke-1,i ,j) )
1037  call check( __line__, pott(ke-1,i+1,j) )
1038  call check( __line__, pott(ke-1,i+2,j) )
1039  call check( __line__, pott(ke-1,i,j-2) )
1040  call check( __line__, pott(ke-1,i,j-1) )
1041  call check( __line__, pott(ke-1,i,j ) )
1042  call check( __line__, pott(ke-1,i,j+1) )
1043  call check( __line__, pott(ke-1,i,j+2) )
1044  call check( __line__, sw(ke-2,i,j) )
1045  call check( __line__, sw(ke-1,i,j) )
1046  call check( __line__, su(ke-1,i-1,j) )
1047  call check( __line__, su(ke-1,i ,j) )
1048  call check( __line__, sv(ke-1,i,j-1) )
1049  call check( __line__, sv(ke-1,i,j ) )
1050  call check( __line__, gsqrt(ke-1,i ,j,i_uyz) )
1051  call check( __line__, gsqrt(ke-1,i-1,j,i_uyz) )
1052  call check( __line__, gsqrt(ke-1,i,j-1,i_xvz) )
1053  call check( __line__, gsqrt(ke-1,i,j ,i_xvz) )
1054  call check( __line__, gsqrt(ke-1,i,j,i_xyz) )
1055  call check( __line__, st(ke-1,i,j) )
1056  call check( __line__, dpres(ke-2,i,j) )
1057  call check( __line__, dpres(ke-1,i,j) )
1058  call check( __line__, dpres(ke ,i,j) )
1059  call check( __line__, rt2p(ke-2,i,j) )
1060  call check( __line__, rt2p(ke-1,i,j) )
1061  call check( __line__, rt2p(ke ,i,j) )
1062  call check( __line__, rhot(ke-2,i,j) )
1063  call check( __line__, rhot(ke,i,j) )
1064  call check( __line__, ddens(ke-2,i,j) )
1065  call check( __line__, ddens(ke,i,j) )
1066 #endif
1067  b(ke-1,i,j) = ( &
1068  ( j33g * ( momz(ke-1,i,j) + dtrk*sw(ke-1,i,j) ) &
1069  * ( 0.5_rp*(pott(ke ,i,j)+pott(ke-1,i,j)) ) &
1070  - j33g * ( momz(ke-2,i,j) + dtrk*sw(ke-2,i,j) ) &
1071  * ( fact_n*(pott(ke-1,i,j)+pott(ke-2,i,j)) &
1072  + fact_f*(pott(ke ,i,j)+pott(ke-3,i,j)) ) ) * rcdz(ke-1) &
1073  + ( gsqrt(ke-1,i ,j,i_uyz) * ( momx(ke-1,i ,j) + dtrk*su(ke-1,i ,j) ) &
1074  * ( fact_n*(pott(ke-1,i+1,j)+pott(ke-1,i ,j)) &
1075  + fact_f*(pott(ke-1,i+2,j)+pott(ke-1,i-1,j)) ) &
1076  - gsqrt(ke-1,i-1,j,i_uyz) * ( momx(ke-1,i-1,j) + dtrk*su(ke-1,i-1,j) ) &
1077  * ( fact_n*(pott(ke-1,i ,j)+pott(ke-1,i-1,j)) &
1078  + fact_f*(pott(ke-1,i+1,j)+pott(ke-1,i-2,j)) ) ) * rcdx(i) &
1079  + ( gsqrt(ke-1,i,j ,i_xvz) * ( momy(ke-1,i,j ) + dtrk*sv(ke-1,i,j ) ) &
1080  * ( fact_n*(pott(ke-1,i,j+1)+pott(ke-1,i,j )) &
1081  + fact_f*(pott(ke-1,i,j+2)+pott(ke-1,i,j-1)) ) &
1082  - gsqrt(ke-1,i,j-1,i_xvz) * ( momy(ke-1,i,j-1) + dtrk*sv(ke-1,i,j-1) ) &
1083  * ( fact_n*(pott(ke-1,i,j )+pott(ke-1,i,j-1)) &
1084  + fact_f*(pott(ke-1,i,j+1)+pott(ke-1,i,j-2)) ) ) * rcdy(j) &
1085  + gsqrt(ke-1,i,j,i_xyz) * ( st(ke-1,i,j) - dpres(ke-1,i,j) * rcs2t(ke-1,i,j) * rdt ) &
1086  ) * rdt &
1087  + grav * j33g * ( ( dpres(ke ,i,j)*rcs2t(ke ,i,j) &
1088  - dpres(ke-2,i,j)*rcs2t(ke-2,i,j) ) &
1089  - ( rhot(ke ,i,j)*ddens(ke ,i,j)/dens(ke ,i,j) &
1090  - rhot(ke-2,i,j)*ddens(ke-2,i,j)/dens(ke-2,i,j) )&
1091  ) / ( fdz(ke-1) + fdz(ke-1-1) )
1092 #ifdef DEBUG
1093  call check( __line__, momz(ke-1,i,j) )
1094  call check( __line__, momx(ke,i-1,j) )
1095  call check( __line__, momx(ke,i ,j) )
1096  call check( __line__, momy(ke,i,j-1) )
1097  call check( __line__, momy(ke,i,j ) )
1098  call check( __line__, sw(ke-1,i,j) )
1099  call check( __line__, su(ke,i-1,j) )
1100  call check( __line__, su(ke,i ,j) )
1101  call check( __line__, sv(ke,i,j-1) )
1102  call check( __line__, sv(ke,i,j ) )
1103  call check( __line__, gsqrt(ke,i-1,j,i_uyz) )
1104  call check( __line__, pott(ke-1,i,j) )
1105  call check( __line__, pott(ke ,i,j) )
1106  call check( __line__, pott(ke,i-2,j) )
1107  call check( __line__, pott(ke,i-1,j) )
1108  call check( __line__, pott(ke,i ,j) )
1109  call check( __line__, pott(ke,i+1,j) )
1110  call check( __line__, pott(ke,i+2,j) )
1111  call check( __line__, pott(ke,i,j-2) )
1112  call check( __line__, pott(ke,i,j-1) )
1113  call check( __line__, pott(ke,i,j ) )
1114  call check( __line__, pott(ke,i,j+1) )
1115  call check( __line__, pott(ke,i,j+2) )
1116  call check( __line__, gsqrt(ke,i-1,j,i_uyz) )
1117  call check( __line__, gsqrt(ke,i ,j,i_uyz) )
1118  call check( __line__, gsqrt(ke,i,j-1,i_xvz) )
1119  call check( __line__, gsqrt(ke,i,j ,i_xvz) )
1120  call check( __line__, gsqrt(ke,i,j,i_xyz) )
1121  call check( __line__, st(ke,i,j) )
1122  call check( __line__, dpres(ke-1,i,j) )
1123  call check( __line__, dpres(ke ,i,j) )
1124  call check( __line__, rt2p(ke-1,i,j) )
1125  call check( __line__, rt2p(ke ,i,j) )
1126  call check( __line__, rhot(ke-1,i,j) )
1127  call check( __line__, ddens(ke-1,i,j) )
1128 #endif
1129  b(ke,i,j) = ( &
1130  ( &
1131  - j33g * ( momz(ke-1,i,j) + dtrk*sw(ke-1,i,j) ) &
1132  * 0.5_rp*(pott(ke ,i,j)+pott(ke-1,i,j)) ) * rcdz(ke) &
1133  + ( gsqrt(ke,i ,j,i_uyz) * ( momx(ke,i ,j) + dtrk*su(ke,i ,j) ) &
1134  * ( fact_n*(pott(ke,i+1,j)+pott(ke,i ,j)) &
1135  + fact_f*(pott(ke,i+2,j)+pott(ke,i-1,j)) ) &
1136  - gsqrt(ke,i-1,j,i_uyz) * ( momx(ke,i-1,j) + dtrk*su(ke,i-1,j) ) &
1137  * ( fact_n*(pott(ke,i ,j)+pott(ke,i-1,j)) &
1138  + fact_f*(pott(ke,i+1,j)+pott(ke,i-2,j)) ) ) * rcdx(i) &
1139  + ( gsqrt(ke,i,j ,i_xvz) * ( momy(ke,i,j ) + dtrk*sv(ke,i,j ) ) &
1140  * ( fact_n*(pott(ke,i,j+1)+pott(ke,i,j )) &
1141  + fact_f*(pott(ke,i,j+2)+pott(ke,i,j-1)) ) &
1142  - gsqrt(ke,i,j-1,i_xvz) * ( momy(ke,i,j-1) + dtrk*sv(ke,i,j-1) ) &
1143  * ( fact_n*(pott(ke,i,j )+pott(ke,i,j-1)) &
1144  + fact_f*(pott(ke,i,j+1)+pott(ke,i,j-2)) ) ) * rcdy(j) &
1145  + gsqrt(ke,i,j,i_xyz) * ( st(ke,i,j) - dpres(ke,i,j) * rcs2t(ke,i,j) * rdt ) &
1146  ) * rdt &
1147  + grav * j33g * 0.5_rp * ( - ( dpres(ke,i,j)+dpres(ke-1,i,j) ) * rcs2t(ke,i,j) &
1148  + ( ddens(ke,i,j)+ddens(ke-1,i,j) ) * pott(ke,i,j) ) * rcdz(ke)
1149  enddo
1150  enddo
1151 #ifdef DEBUG
1152  k = iundef; i = iundef; j = iundef
1153 #endif
1154 
1155  call make_matrix( &
1156  m, & ! (inout)
1157  pott, rcs2t, grav, &
1158  gsqrt, j33g, &
1159  rcdz, rfdz, rcdx, rfdx, rcdy,rfdy, fdz, &
1160  rdt, &
1161  fact_n, fact_f, &
1162  i_xyz, i_xyw, &
1163  iis, iie, jjs, jje )
1164 
1165  enddo
1166  enddo
1167 
1168  call solve_bicgstab( &
1169  dpres_n, & ! (out)
1170  dpres, & ! (in)
1171  m, b ) ! (in)
1172 
1173  call comm_vars8( dpres_n, 1 )
1174  call comm_wait ( dpres_n, 1 )
1175 
1176 #ifdef DEBUG
1177  call check_solver( dpres_n, m, b )
1178 #endif
1179 
1180 #else
1181 
1182  call solve_multigrid
1183 
1184 #endif
1185 
1186  do jjs = js, je, jblock
1187  jje = jjs+jblock-1
1188  do iis = is, ie, iblock
1189  iie = iis+iblock-1
1190 
1191  !##### momentum equation (z) #####
1192 
1193  !--- update momentum(z)
1194  !$omp parallel do private(i,j,k,duvw) OMP_SCHEDULE_ collapse(2)
1195  do j = jjs, jje
1196  do i = iis, iie
1197  do k = ks, ke-1
1198 #ifdef DEBUG
1199  call check( __line__, dpres_n(k+1,i,j) )
1200  call check( __line__, dpres_n(k ,i,j) )
1201  call check( __line__, dpres(k+1,i,j) )
1202  call check( __line__, dpres(k ,i,j) )
1203  call check( __line__, ddens(k+1,i,j) )
1204  call check( __line__, ddens(k ,i,j) )
1205  call check( __line__, ref_dens(k+1,i,j) )
1206  call check( __line__, ref_dens(k,i,j) )
1207  call check( __line__, momz0(k,i,j) )
1208 #endif
1209  duvw = dtrk * ( ( &
1210  - j33g * ( dpres_n(k+1,i,j) - dpres_n(k,i,j) ) * rfdz(k) &
1211  ) / gsqrt(k,i,j,i_uyz) &
1212  - 0.5_rp * grav &
1213  * ( ddens(k+1,i,j) &
1214  + ( dpres_n(k+1,i,j) - dpres(k+1,i,j) ) &
1215  / ( pott(k+1,i,j) * rt2p(k+1,i,j) ) &
1216  + ddens(k ,i,j) &
1217  + ( dpres_n(k ,i,j) - dpres(k ,i,j) ) &
1218  / ( pott(k ,i,j) * rt2p(k ,i,j) ) ) &
1219  + sw(k,i,j) )
1220  momz_rk(k,i,j) = momz0(k,i,j) + duvw
1221  mflx_hi(k,i,j,zdir) = j33g * ( momz(k,i,j) + duvw )
1222  enddo
1223  enddo
1224  enddo
1225 #ifdef DEBUG
1226  k = iundef; i = iundef; j = iundef
1227 #endif
1228  do j = jjs, jje
1229  do i = iis, iie
1230  mflx_hi(ks-1,i,j,zdir) = 0.0_rp
1231  mflx_hi(ke ,i,j,zdir) = 0.0_rp
1232  enddo
1233  enddo
1234 #ifdef DEBUG
1235  k = iundef; i = iundef; j = iundef
1236 #endif
1237 
1238 
1239  !##### momentum equation (x) #####
1240 
1241  !--- update momentum(x)
1242  !$omp parallel do private(i,j,k,duvw) OMP_SCHEDULE_ collapse(2)
1243  do j = jjs , jje
1244  do i = iis-1, iie
1245  do k = ks, ke
1246 #ifdef DEBUG
1247  call check( __line__, dpres_n(k,i+1,j) )
1248  call check( __line__, dpres_n(k,i ,j) )
1249  call check( __line__, gsqrt(k,i+1,j,i_xyz) )
1250  call check( __line__, gsqrt(k,i ,j,i_xyz) )
1251  call check( __line__, gsqrt(k,i,j,i_uyz) )
1252  call check( __line__, su(k,i,j) )
1253  call check( __line__, momx0(k,i,j) )
1254 #endif
1255  duvw = dtrk * ( ( &
1256  - ( gsqrt(k,i+1,j,i_xyz) * dpres_n(k,i+1,j) &
1257  - gsqrt(k,i ,j,i_xyz) * dpres_n(k,i ,j) ) * rfdx(i) & ! pressure gradient force
1258  ) / gsqrt(k,i,j,i_uyz) &
1259  + su(k,i,j) )
1260 
1261  momx_rk(k,i,j) = momx0(k,i,j) + duvw
1262  mflx_hi(k,i,j,xdir) = gsqrt(k,i,j,i_uyz) * ( momx(k,i,j) + duvw )
1263  enddo
1264  enddo
1265  enddo
1266 #ifdef DEBUG
1267  k = iundef; i = iundef; j = iundef
1268 #endif
1269 
1270  !##### momentum equation (y) #####
1271 
1272  !--- update momentum(y)
1273  !$omp parallel do private(i,j,k,duvw) OMP_SCHEDULE_ collapse(2)
1274  do j = jjs-1, jje
1275  do i = iis, iie
1276  do k = ks, ke
1277 #ifdef DEBUG
1278  call check( __line__, dpres_n(k,i,j ) )
1279  call check( __line__, dpres_n(k,i,j+1) )
1280  call check( __line__, gsqrt(k,i,j ,i_xyz) )
1281  call check( __line__, gsqrt(k,i,j+1,i_xyz) )
1282  call check( __line__, gsqrt(k,i,j,i_xvz) )
1283  call check( __line__, sv(k,i,j) )
1284  call check( __line__, momy0(k,i,j) )
1285 #endif
1286  duvw = dtrk * ( ( &
1287  - ( gsqrt(k,i,j+1,i_xyz) * dpres_n(k,i,j+1) &
1288  - gsqrt(k,i,j ,i_xyz) * dpres_n(k,i,j ) ) * rfdy(j) &
1289  ) / gsqrt(k,i,j,i_xvz) &
1290  + sv(k,i,j) )
1291  momy_rk(k,i,j) = momy0(k,i,j) + duvw
1292  mflx_hi(k,i,j,ydir) = gsqrt(k,i,j,i_xvz) * ( momy(k,i,j) + duvw )
1293  enddo
1294  enddo
1295  enddo
1296 #ifdef DEBUG
1297  k = iundef; i = iundef; j = iundef
1298 #endif
1299 
1300 
1301  !##### Thermodynamic Equation #####
1302 
1303  ! at (x, y, w)
1304  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
1305  do j = jjs, jje
1306  do i = iis, iie
1307  do k = ks+1, ke-2
1308 #ifdef DEBUG
1309  call check( __line__, mflx_hi(k,i,j,zdir) )
1310  call check( __line__, pott(k,i-1,j) )
1311  call check( __line__, pott(k,i ,j) )
1312  call check( __line__, pott(k,i+1,j) )
1313  call check( __line__, pott(k,i+1,j) )
1314  call check( __line__, num_diff(k,i,j,i_rhot,xdir) )
1315 #endif
1316  tflx_hi(k,i,j,zdir) = mflx_hi(k,i,j,zdir) &
1317  * ( fact_n * ( pott(k+1,i,j) + pott(k ,i,j) ) &
1318  + fact_f * ( pott(k+2,i,j) + pott(k-1,i,j) ) )
1319  enddo
1320  enddo
1321  enddo
1322 #ifdef DEBUG
1323  k = iundef; i = iundef; j = iundef
1324 #endif
1325  !$omp parallel do private(i,j) OMP_SCHEDULE_ collapse(2)
1326  do j = jjs, jje
1327  do i = iis, iie
1328  tflx_hi(ks-1,i,j,zdir) = 0.0_rp
1329  tflx_hi(ks ,i,j,zdir) = mflx_hi(ks ,i,j,zdir) * 0.5_rp * ( pott(ks+1,i,j) + pott(ks ,i,j) )
1330  tflx_hi(ke-1,i,j,zdir) = mflx_hi(ke-1,i,j,zdir) * 0.5_rp * ( pott(ke ,i,j) + pott(ke-1,i,j) )
1331  tflx_hi(ke ,i,j,zdir) = 0.0_rp
1332  enddo
1333  enddo
1334 #ifdef DEBUG
1335  k = iundef; i = iundef; j = iundef
1336 #endif
1337  ! at (u, y, z)
1338  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
1339  do j = jjs, jje
1340  do i = iis-1, iie
1341  do k = ks, ke
1342 #ifdef DEBUG
1343  call check( __line__, mflx_hi(k,i,j,xdir) )
1344  call check( __line__, pott(k,i-1,j) )
1345  call check( __line__, pott(k,i ,j) )
1346  call check( __line__, pott(k,i+1,j) )
1347  call check( __line__, pott(k,i+1,j) )
1348  call check( __line__, num_diff(k,i,j,i_rhot,xdir) )
1349 #endif
1350  tflx_hi(k,i,j,xdir) = mflx_hi(k,i,j,xdir) &
1351  * ( fact_n * ( pott(k,i+1,j)+pott(k,i ,j) ) &
1352  + fact_f * ( pott(k,i+2,j)+pott(k,i-1,j) ) )
1353  enddo
1354  enddo
1355  enddo
1356 #ifdef DEBUG
1357  k = iundef; i = iundef; j = iundef
1358 #endif
1359  ! at (x, v, z)
1360  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
1361  do j = jjs-1, jje
1362  do i = iis, iie
1363  do k = ks, ke
1364 #ifdef DEBUG
1365  call check( __line__, mflx_hi(k,i,j,ydir) )
1366  call check( __line__, pott(k,i,j-1) )
1367  call check( __line__, pott(k,i,j ) )
1368  call check( __line__, pott(k,i,j+1) )
1369  call check( __line__, pott(k,i,j+2) )
1370  call check( __line__, num_diff(k,i,j,i_rhot,ydir) )
1371 #endif
1372  tflx_hi(k,i,j,ydir) = mflx_hi(k,i,j,ydir) &
1373  * ( fact_n * ( pott(k,i,j+1)+pott(k,i,j ) ) &
1374  + fact_f * ( pott(k,i,j+2)+pott(k,i,j-1) ) )
1375  enddo
1376  enddo
1377  enddo
1378 #ifdef DEBUG
1379  k = iundef; i = iundef; j = iundef
1380 #endif
1381 
1382  ! rho*theta
1383  do j = jjs, jje
1384  do i = iis, iie
1385  do k = ks, ke
1386  rhot_rk(k,i,j) = rhot0(k,i,j) &
1387  + dtrk * ( - ( ( tflx_hi(k,i,j,zdir) - tflx_hi(k-1,i ,j ,zdir) ) * rcdz(k) &
1388  + ( tflx_hi(k,i,j,xdir) - tflx_hi(k ,i-1,j ,xdir) ) * rcdx(i) &
1389  + ( tflx_hi(k,i,j,ydir) - tflx_hi(k ,i ,j-1,ydir) ) * rcdy(j) ) &
1390  / gsqrt(k,i,j,i_xyz) &
1391  + st(k,i,j) )
1392  enddo
1393  enddo
1394  enddo
1395 
1396  !##### continuous equation #####
1397 
1398  ! total momentum flux
1399  do j = jjs, jje
1400  do i = iis, iie
1401  do k = ks, ke-1
1402  mflx_hi(k,i,j,zdir) = mflx_hi(k,i,j,zdir) + mflx_hi2(k,i,j,zdir)
1403  enddo
1404  enddo
1405  enddo
1406  do j = jjs , jje
1407  do i = iis-1, iie
1408  do k = ks, ke
1409  mflx_hi(k,i,j,xdir) = mflx_hi(k,i,j,xdir) + mflx_hi2(k,i,j,xdir)
1410  enddo
1411  enddo
1412  enddo
1413  do j = jjs-1, jje
1414  do i = iis , iie
1415  do k = ks, ke
1416  mflx_hi(k,i,j,ydir) = mflx_hi(k,i,j,ydir) + mflx_hi2(k,i,j,ydir)
1417  enddo
1418  enddo
1419  enddo
1420 
1421  !--- update density
1422  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
1423  do j = jjs, jje
1424  do i = iis, iie
1425  do k = ks, ke
1426 #ifdef DEBUG
1427  call check( __line__, dens0(k,i,j) )
1428  call check( __line__, mflx_hi(k ,i ,j ,xdir) )
1429  call check( __line__, mflx_hi(k ,i-1,j ,xdir) )
1430  call check( __line__, mflx_hi(k ,i ,j ,ydir) )
1431  call check( __line__, mflx_hi(k ,i ,j-1,ydir) )
1432  call check( __line__, dens_t(k,i,j) )
1433 #endif
1434  dens_rk(k,i,j) = dens0(k,i,j) &
1435  + dtrk * ( - ( ( mflx_hi(k,i,j,zdir)-mflx_hi(k-1,i ,j, zdir) ) * rcdz(k) &
1436  + ( mflx_hi(k,i,j,xdir)-mflx_hi(k ,i-1,j, xdir) ) * rcdx(i) &
1437  + ( mflx_hi(k,i,j,ydir)-mflx_hi(k ,i, j-1,ydir) ) * rcdy(j) ) &
1438  / gsqrt(k,i,j,i_xyz) & ! divergence
1439  + dens_t(k,i,j) )
1440  enddo
1441  enddo
1442  enddo
1443 #ifdef DEBUG
1444  k = iundef; i = iundef; j = iundef
1445 #endif
1446 
1447 
1448 #ifdef DEBUG
1449  call check_pres( &
1450  dpres_n, dpres, rhot_rk, rhot, dens_rk, dens, b,&
1451  rt2p )
1452 #endif
1453 
1454  enddo
1455  enddo
1456 
1457  return
1458  end subroutine atmos_dyn_tstep_short_fvm_hivi
1459 
1460 #ifdef HIVI_BICGSTAB
1461  subroutine solve_bicgstab( &
1462  DPRES_N, &
1463  DPRES, &
1464  M, B )
1465  use scale_prc, only: &
1466  prc_abort
1467  use scale_comm_cartesc, only: &
1468  comm_world, &
1469  comm_vars8, &
1470  comm_wait
1471  implicit none
1472  real(RP), intent(out) :: DPRES_N(ka,ia,ja)
1473  real(RP), intent(in) :: DPRES(ka,ia,ja)
1474  real(RP), intent(in) :: M(7,ka,ia,ja)
1475  real(RP), intent(in) :: B(ka,ia,ja)
1476 
1477  real(RP) :: r0(ka,ia,ja)
1478 
1479  real(RP) :: p(ka,ia,ja)
1480  real(RP) :: Mp(ka,ia,ja)
1481  real(RP) :: s(ka,ia,ja)
1482  real(RP) :: Ms(ka,ia,ja)
1483  real(RP) :: al, be, w
1484 
1485  real(RP), pointer :: r(:,:,:)
1486  real(RP), pointer :: rn(:,:,:)
1487  real(RP), pointer :: swap(:,:,:)
1488  real(RP), target :: v0(ka,ia,ja)
1489  real(RP), target :: v1(ka,ia,ja)
1490  real(RP) :: r0r
1491  real(RP) :: norm, error, error2
1492 
1493  real(RP) :: iprod(2)
1494  real(RP) :: buf(2)
1495 
1496  integer :: k, i, j
1497  integer :: iis, iie, jjs, jje
1498  integer :: iter
1499  integer :: ierror
1500 
1501 #ifdef DEBUG
1502  r0(:,:,:) = undef
1503  p(:,:,:) = undef
1504  mp(:,:,:) = undef
1505  s(:,:,:) = undef
1506  ms(:,:,:) = undef
1507 
1508  v0(:,:,:) = undef
1509  v1(:,:,:) = undef
1510 #endif
1511 
1512  r => v0
1513  rn => v1
1514 
1515  call mul_matrix( v1, m, dpres ) ! v1 = M x0
1516 
1517  norm = 0.0_rp
1518  r0r = 0.0_rp
1519 
1520  do jjs = js, je, jblock
1521  jje = jjs+jblock-1
1522  do iis = is, ie, iblock
1523  iie = iis+iblock-1
1524 
1525  do j = jjs, jje
1526  do i = iis, iie
1527  do k = ks, ke
1528  norm = norm + b(k,i,j)**2
1529  enddo
1530  enddo
1531  enddo
1532 #ifdef DEBUG
1533  k = iundef; i = iundef; j = iundef
1534 #endif
1535 
1536  ! r = b - M x0
1537  do j = jjs, jje
1538  do i = iis, iie
1539  do k = ks, ke
1540 #ifdef DEBUG
1541  call check( __line__, b(k,i,j) )
1542  call check( __line__, v1(k,i,j) )
1543 #endif
1544  r(k,i,j) = b(k,i,j) - v1(k,i,j)
1545  enddo
1546  enddo
1547  enddo
1548 #ifdef DEBUG
1549  k = iundef; i = iundef; j = iundef
1550 #endif
1551 
1552  do j = jjs, jje
1553  do i = iis, iie
1554  do k = ks, ke
1555  r0(k,i,j) = r(k,i,j)
1556  p(k,i,j) = r(k,i,j)
1557  enddo
1558  enddo
1559  enddo
1560 #ifdef DEBUG
1561  k = iundef; i = iundef; j = iundef
1562 #endif
1563 
1564  do j = jjs, jje
1565  do i = iis, iie
1566  do k = ks, ke
1567 #ifdef DEBUG
1568  call check( __line__, r(k,i,j) )
1569  call check( __line__, r0(k,i,j) )
1570 #endif
1571  r0r = r0r + r0(k,i,j) * r(k,i,j)
1572  enddo
1573  enddo
1574  enddo
1575 #ifdef DEBUG
1576  k = iundef; i = iundef; j = iundef
1577 #endif
1578 
1579  enddo
1580  enddo
1581 
1582  do j = js-1, je+1
1583  do i = is-1, ie+1
1584  do k = ks, ke
1585  dpres_n(k,i,j) = dpres(k,i,j)
1586  end do
1587  end do
1588  end do
1589 
1590 
1591  iprod(1) = r0r
1592  iprod(2) = norm
1593  call mpi_allreduce(iprod, buf, 2, mtype, mpi_sum, comm_world, ierror)
1594  r0r = buf(1)
1595  norm = buf(2)
1596 
1597  error2 = norm
1598 
1599  do iter = 1, itmax
1600 
1601  error = 0.0_rp
1602  do j = js, je
1603  do i = is, ie
1604  do k = ks, ke
1605 #ifdef DEBUG
1606  call check( __line__, r(k,i,j) )
1607 #endif
1608  error = error + r(k,i,j)**2
1609  enddo
1610  enddo
1611  enddo
1612 #ifdef DEBUG
1613  k = iundef; i = iundef; j = iundef
1614 #endif
1615  call mpi_allreduce(error, buf, 1, mtype, mpi_sum, comm_world, ierror)
1616  error = buf(1)
1617 
1618 #ifdef DEBUG
1619  log_info("solve_bicgstab",*) iter, error/norm
1620 #endif
1621  if ( sqrt(error/norm) < epsilon .OR. error > error2 ) then
1622 #ifdef DEBUG
1623  log_info("solve_bicgstab",*) "Bi-CGSTAB converged:", iter
1624 #endif
1625  exit
1626  endif
1627  error2 = error
1628 
1629  call comm_vars8( p, 1 )
1630  call comm_wait ( p, 1 )
1631  call mul_matrix( mp, m, p )
1632 
1633  iprod(1) = 0.0_rp
1634  do j = js, je
1635  do i = is, ie
1636  do k = ks, ke
1637 #ifdef DEBUG
1638  call check( __line__, r0(k,i,j) )
1639  call check( __line__, mp(k,i,j) )
1640 #endif
1641  iprod(1) = iprod(1) + r0(k,i,j) * mp(k,i,j)
1642  enddo
1643  enddo
1644  enddo
1645 #ifdef DEBUG
1646  k = iundef; i = iundef; j = iundef
1647 #endif
1648  call mpi_allreduce(iprod, buf, 1, mtype, mpi_sum, comm_world, ierror)
1649  al = r0r / buf(1) ! (r0,r) / (r0,Mp)
1650 
1651  do j = js, je
1652  do i = is, ie
1653  do k = ks, ke
1654 #ifdef DEBUG
1655  call check( __line__, r(k,i,j) )
1656  call check( __line__, mp(k,i,j) )
1657 #endif
1658  s(k,i,j) = r(k,i,j) - al*mp(k,i,j)
1659  enddo
1660  enddo
1661  enddo
1662 #ifdef DEBUG
1663  k = iundef; i = iundef; j = iundef
1664 #endif
1665 
1666  call comm_vars8( s, 1 )
1667  call comm_wait ( s, 1 )
1668  call mul_matrix( ms, m, s )
1669  iprod(1) = 0.0_rp
1670  iprod(2) = 0.0_rp
1671  do j = js, je
1672  do i = is, ie
1673  do k = ks, ke
1674 #ifdef DEBUG
1675  call check( __line__, ms(k,i,j) )
1676  call check( __line__, s(k,i,j) )
1677 #endif
1678  iprod(1) = iprod(1) + ms(k,i,j) * s(k,i,j)
1679  iprod(2) = iprod(2) + ms(k,i,j) * ms(k,i,j)
1680  enddo
1681  enddo
1682  enddo
1683 #ifdef DEBUG
1684  k = iundef; i = iundef; j = iundef
1685 #endif
1686  call mpi_allreduce(iprod, buf, 2, mtype, mpi_sum, comm_world, ierror)
1687  w = buf(1) / buf(2) ! (Ms,s) / (Ms,Ms)
1688 
1689  iprod(1) = 0.0_rp
1690 
1691  do jjs = js, je, jblock
1692  jje = jjs+jblock-1
1693  do iis = is, ie, iblock
1694  iie = iis+iblock-1
1695 
1696  do j = jjs, jje
1697  do i = iis, iie
1698  do k = ks, ke
1699 #ifdef DEBUG
1700  call check( __line__, dpres_n(k,i,j) )
1701  call check( __line__, p(k,i,j) )
1702  call check( __line__, s(k,i,j) )
1703 #endif
1704  dpres_n(k,i,j) = dpres_n(k,i,j) + al*p(k,i,j) + w*s(k,i,j)
1705  enddo
1706  enddo
1707  enddo
1708 #ifdef DEBUG
1709  k = iundef; i = iundef; j = iundef
1710 #endif
1711 
1712  do j = jjs, jje
1713  do i = iis, iie
1714  do k = ks, ke
1715 #ifdef DEBUG
1716  call check( __line__, s(k,i,j) )
1717  call check( __line__, ms(k,i,j) )
1718 #endif
1719  rn(k,i,j) = s(k,i,j) - w*ms(k,i,j)
1720  enddo
1721  enddo
1722  enddo
1723 #ifdef DEBUG
1724  k = iundef; i = iundef; j = iundef
1725 #endif
1726 
1727  do j = jjs, jje
1728  do i = iis, iie
1729  do k = ks, ke
1730 #ifdef DEBUG
1731  call check( __line__, r0(k,i,j) )
1732  call check( __line__, rn(k,i,j) )
1733 #endif
1734  iprod(1) = iprod(1) + r0(k,i,j) * rn(k,i,j)
1735  enddo
1736  enddo
1737  enddo
1738 #ifdef DEBUG
1739  k = iundef; i = iundef; j = iundef
1740 #endif
1741 
1742  enddo
1743  enddo
1744 
1745  be = al/w / r0r
1746 
1747  call mpi_allreduce(iprod, r0r, 1, mtype, mpi_sum, comm_world, ierror)
1748 
1749  be = be * r0r ! al/w * (r0,rn)/(r0,r)
1750  do j = js, je
1751  do i = is, ie
1752  do k = ks, ke
1753 #ifdef DEBUG
1754  call check( __line__, rn(k,i,j) )
1755  call check( __line__, p(k,i,j) )
1756  call check( __line__, mp(k,i,j) )
1757 #endif
1758  p(k,i,j) = rn(k,i,j) + be * ( p(k,i,j) - w*mp(k,i,j) )
1759  enddo
1760  enddo
1761  enddo
1762 #ifdef DEBUG
1763  k = iundef; i = iundef; j = iundef
1764 #endif
1765 
1766  swap => rn
1767  rn => r
1768  r => swap
1769 #ifdef DEBUG
1770  rn(:,:,:) = undef
1771 #endif
1772  enddo
1773 
1774  if ( iter >= itmax ) then
1775  log_error("solve_bicgstab",*) 'not converged', error, norm
1776  call prc_abort
1777  endif
1778 
1779  return
1780  end subroutine solve_bicgstab
1781 
1782  subroutine make_matrix(&
1783  M, &
1784  POTT, RCs2T, GRAV, &
1785  G, J33G, &
1786  RCDZ, RFDZ, RCDX, RFDX, RCDY, RFDY, FDZ, &
1787  rdt, &
1788  FACT_N, FACT_F, &
1789  I_XYZ, I_XYW, &
1790  IIS, IIE, JJS, JJE )
1791  implicit none
1792  real(RP), intent(inout) :: M(7,ka,ia,ja)
1793  real(RP), intent(in) :: POTT(ka,ia,ja)
1794  real(RP), intent(in) :: RCs2T(ka,ia,ja)
1795  real(RP), intent(in) :: GRAV
1796  real(RP), intent(in) :: G(ka,ia,ja,7)
1797  real(RP), intent(in) :: J33G
1798  real(RP), intent(in) :: RCDZ(ka)
1799  real(RP), intent(in) :: RFDZ(ka-1)
1800  real(RP), intent(in) :: RFDX(ia-1)
1801  real(RP), intent(in) :: RCDX(ia)
1802  real(RP), intent(in) :: RCDY(ja)
1803  real(RP), intent(in) :: RFDY(ja-1)
1804  real(RP), intent(in) :: FDZ(ka-1)
1805  real(RP), intent(in) :: rdt
1806  real(RP), intent(in) :: FACT_N
1807  real(RP), intent(in) :: FACT_F
1808  integer, intent(in) :: I_XYZ
1809  integer, intent(in) :: I_XYW
1810  integer, intent(in) :: IIS
1811  integer, intent(in) :: IIE
1812  integer, intent(in) :: JJS
1813  integer, intent(in) :: JJE
1814 
1815  integer :: k, i, j
1816 
1817  do j = jjs, jje
1818  do i = iis, iie
1819  do k = ks+2, ke-2
1820 #ifdef DEBUG
1821  call check( __line__, pott(k-2,i,j) )
1822  call check( __line__, pott(k-1,i,j) )
1823  call check( __line__, pott(k ,i,j) )
1824  call check( __line__, pott(k+1,i,j) )
1825  call check( __line__, pott(k+2,i,j) )
1826  call check( __line__, pott(k,i-2,j) )
1827  call check( __line__, pott(k,i-1,j) )
1828  call check( __line__, pott(k,i ,j) )
1829  call check( __line__, pott(k,i+1,j) )
1830  call check( __line__, pott(k,i+2,j) )
1831  call check( __line__, pott(k,i,j-2) )
1832  call check( __line__, pott(k,i,j-1) )
1833  call check( __line__, pott(k,i,j ) )
1834  call check( __line__, pott(k,i,j+1) )
1835  call check( __line__, pott(k,i,j+2) )
1836  call check( __line__, g(k-1,i,j,i_xyw) )
1837  call check( __line__, g(k+1,i,j,i_xyw) )
1838  call check( __line__, g(k,i,j,i_xyz) )
1839  call check( __line__, rcs2t(k-1,i,j) )
1840  call check( __line__, rcs2t(k ,i,j) )
1841  call check( __line__, rcs2t(k+1,i,j) )
1842 #endif
1843  ! k,i,j
1844  m(1,k,i,j) = &
1845  - ( ( fact_n * (pott(k+1,i,j)+pott(k ,i,j)) &
1846  + fact_f * (pott(k+2,i,j)+pott(k-1,i,j)) ) * rfdz(k ) / g(k ,i,j,i_xyw) &
1847  + ( fact_n * (pott(k ,i,j)+pott(k-1,i,j)) &
1848  + fact_f * (pott(k+1,i,j)+pott(k-2,i,j)) ) * rfdz(k-1) / g(k-1,i,j,i_xyw) &
1849  ) * j33g * j33g * rfdz(k) &
1850  - ( ( fact_n * (pott(k,i+1,j)+pott(k,i ,j)) &
1851  + fact_f * (pott(k,i+2,j)+pott(k,i-1,j)) ) * rfdx(i ) &
1852  + ( fact_n * (pott(k,i ,j)+pott(k,i-1,j)) &
1853  + fact_f * (pott(k,i+1,j)+pott(k,i-2,j)) ) * rfdx(i-1) &
1854  ) * g(k,i,j,i_xyz) * rfdx(i) &
1855  - ( ( fact_n * (pott(k,i,j+1)+pott(k,i,j )) &
1856  + fact_f * (pott(k,i,j+2)+pott(k,i,j-1)) ) * rfdy(j ) &
1857  + ( fact_n * (pott(k,i,j )+pott(k,i,j-1)) &
1858  + fact_f * (pott(k,i,j-1)+pott(k,i,j-2)) ) * rfdy(j-1) &
1859  ) * g(k,i,j,i_xyz) * rfdy(j) &
1860  - g(k,i,j,i_xyz) * rcs2t(k,i,j) * rdt * rdt
1861  ! k-1
1862  m(2,k,i,j) = j33g * j33g / g(k-1,i,j,i_xyw) &
1863  * ( fact_n * (pott(k ,i,j)+pott(k-1,i,j)) &
1864  + fact_f * (pott(k+1,i,j)+pott(k-2,i,j)) ) &
1865  * rfdz(k-1) * rcdz(k) &
1866  - grav * j33g * rcs2t(k-1,i,j) / ( fdz(k)+fdz(k-1) )
1867  ! k+1
1868  m(3,k,i,j) = j33g * j33g / g(k+1,i,j,i_xyw) &
1869  * ( fact_n * (pott(k+1,i,j)+pott(k ,i,j)) &
1870  + fact_f * (pott(k+2,i,j)+pott(k-1,i,j)) ) &
1871  * rfdz(k ) * rcdz(k) &
1872  + grav * j33g * rcs2t(k+1,i,j) / ( fdz(k)+fdz(k-1) )
1873  enddo
1874  enddo
1875  enddo
1876 #ifdef DEBUG
1877  k = iundef; i = iundef; j = iundef
1878 #endif
1879  do j = jjs, jje
1880  do i = iis, iie
1881 #ifdef DEBUG
1882  call check( __line__, pott(ks ,i,j) )
1883  call check( __line__, pott(ks+1,i,j) )
1884  call check( __line__, pott(ks,i-2,j) )
1885  call check( __line__, pott(ks,i-1,j) )
1886  call check( __line__, pott(ks,i ,j) )
1887  call check( __line__, pott(ks,i+1,j) )
1888  call check( __line__, pott(ks,i+2,j) )
1889  call check( __line__, pott(ks,i,j-2) )
1890  call check( __line__, pott(ks,i,j-1) )
1891  call check( __line__, pott(ks,i,j ) )
1892  call check( __line__, pott(ks,i,j+1) )
1893  call check( __line__, pott(ks,i,j+2) )
1894  call check( __line__, g(ks,i,j,i_xyz) )
1895  call check( __line__, rcs2t(ks,i,j) )
1896  call check( __line__, pott(ks ,i,j) )
1897  call check( __line__, pott(ks+1,i,j) )
1898  call check( __line__, g(ks+1,i,j,i_xyw) )
1899  call check( __line__, rcs2t(ks+1,i,j) )
1900 #endif
1901  ! k,i,j
1902  m(1,ks,i,j) = &
1903  - ( 0.5_rp * (pott(ks+1,i,j)+pott(ks ,i,j)) * rfdz(ks ) &
1904  ) * j33g * j33g / g(ks ,i,j,i_xyw) * rfdz(ks) &
1905  - ( ( fact_n * (pott(ks,i+1,j)+pott(ks,i ,j)) &
1906  + fact_f * (pott(ks,i+2,j)+pott(ks,i-1,j)) ) * rfdx(i ) &
1907  + ( fact_n * (pott(ks,i ,j)+pott(ks,i-1,j)) &
1908  + fact_f * (pott(ks,i+1,j)+pott(ks,i-2,j)) ) * rfdx(i-1) &
1909  ) * g(ks,i,j,i_xyz) * rfdx(i) &
1910  - ( ( fact_n * (pott(ks,i,j+1)+pott(ks,i,j )) &
1911  + fact_f * (pott(ks,i,j+2)+pott(ks,i,j-1)) ) * rfdy(j ) &
1912  + ( fact_n * (pott(ks,i,j )+pott(ks,i,j-1)) &
1913  + fact_f * (pott(ks,i,j-1)+pott(ks,i,j-2)) ) * rfdy(j-1) &
1914  ) * g(ks,i,j,i_xyz) * rfdy(j) &
1915  - g(ks,i,j,i_xyz) * rcs2t(ks,i,j) * rdt * rdt &
1916  + grav * j33g * 0.5_rp * rcs2t(ks,i,j) * rcdz(ks)
1917  ! k+1
1918  m(3,ks,i,j) = j33g * j33g / g(ks+1,i,j,i_xyw) &
1919  * 0.5_rp * (pott(ks+1,i,j)+pott(ks ,i,j)) &
1920  * rfdz(ks ) * rcdz(ks) &
1921  + grav * j33g * 0.5_rp * rcs2t(ks+1,i,j) * rcdz(ks)
1922 #ifdef DEBUG
1923  call check( __line__, pott(ks ,i,j) )
1924  call check( __line__, pott(ks+1,i,j) )
1925  call check( __line__, pott(ks+2,i,j) )
1926  call check( __line__, pott(ks+3,i,j) )
1927  call check( __line__, pott(ks+1,i-2,j) )
1928  call check( __line__, pott(ks+1,i-1,j) )
1929  call check( __line__, pott(ks+1,i ,j) )
1930  call check( __line__, pott(ks+1,i+1,j) )
1931  call check( __line__, pott(ks+1,i+2,j) )
1932  call check( __line__, pott(ks+1,i,j-2) )
1933  call check( __line__, pott(ks+1,i,j-1) )
1934  call check( __line__, pott(ks+1,i,j ) )
1935  call check( __line__, pott(ks+1,i,j+1) )
1936  call check( __line__, pott(ks+1,i,j+2) )
1937  call check( __line__, g(ks ,i,j,i_xyw) )
1938  call check( __line__, g(ks+2,i,j,i_xyw) )
1939  call check( __line__, g(ks+1,i,j,i_xyz) )
1940  call check( __line__, rcs2t(ks ,i,j) )
1941  call check( __line__, rcs2t(ks+1 ,i,j) )
1942  call check( __line__, rcs2t(ks+2,i,j) )
1943 #endif
1944  ! k,i,j
1945  m(1,ks+1,i,j) = &
1946  - ( ( fact_n * (pott(ks+2,i,j)+pott(ks+1,i,j)) &
1947  + fact_f * (pott(ks+3,i,j)+pott(ks ,i,j)) ) * rfdz(ks+1) / g(ks+1,i,j,i_xyw) &
1948  + 0.5_rp * (pott(ks+1,i,j)+pott(ks ,i,j)) * rfdz(ks ) / g(ks ,i,j,i_xyw) &
1949  ) * j33g * j33g * rfdz(ks+1) &
1950  - ( ( fact_n * (pott(ks+1,i+1,j)+pott(ks+1,i ,j)) &
1951  + fact_f * (pott(ks+1,i+2,j)+pott(ks+1,i-1,j)) ) * rfdx(i ) &
1952  + ( fact_n * (pott(ks+1,i ,j)+pott(ks+1,i-1,j)) &
1953  + fact_f * (pott(ks+1,i+1,j)+pott(ks+1,i-2,j)) ) * rfdx(i-1) &
1954  ) * g(ks+1,i,j,i_xyz) * rfdx(i) &
1955  - ( ( fact_n * (pott(ks+1,i,j+1)+pott(ks+1,i,j )) &
1956  + fact_f * (pott(ks+1,i,j+2)+pott(ks+1,i,j-1)) ) * rfdy(j ) &
1957  + ( fact_n * (pott(ks+1,i,j )+pott(ks+1,i,j-1)) &
1958  + fact_f * (pott(ks+1,i,j-1)+pott(ks+1,i,j-2)) ) * rfdy(j-1) &
1959  ) * g(ks+1,i,j,i_xyz) * rfdy(j) &
1960  - g(ks+1,i,j,i_xyz) * rcs2t(ks+1,i,j) * rdt * rdt
1961  ! k-1
1962  m(2,ks+1,i,j) = j33g * j33g / g(ks,i,j,i_xyw) &
1963  * 0.5_rp * (pott(ks+1,i,j)+pott(ks ,i,j)) &
1964  * rfdz(ks ) * rcdz(ks+1) &
1965  - grav * j33g * rcs2t(ks ,i,j) / ( fdz(ks+1)+fdz(ks) )
1966  ! k+1
1967  m(3,ks+1,i,j) = j33g * j33g / g(ks+2,i,j,i_xyw) &
1968  * ( fact_n * (pott(ks+2,i,j)+pott(ks+1,i,j)) &
1969  + fact_f * (pott(ks+3,i,j)+pott(ks ,i,j)) ) &
1970  * rfdz(ks+1) * rcdz(ks+1) &
1971  + grav * j33g * rcs2t(ks+2,i,j) / ( fdz(ks+1)+fdz(ks) )
1972 #ifdef DEBUG
1973  call check( __line__, pott(ke-3,i,j) )
1974  call check( __line__, pott(ke-2,i,j) )
1975  call check( __line__, pott(ke-1,i,j) )
1976  call check( __line__, pott(ke ,i,j) )
1977  call check( __line__, pott(ke-1,i-2,j) )
1978  call check( __line__, pott(ke-1,i-1,j) )
1979  call check( __line__, pott(ke-1,i ,j) )
1980  call check( __line__, pott(ke-1,i+1,j) )
1981  call check( __line__, pott(ke-1,i+2,j) )
1982  call check( __line__, pott(ke-1,i,j-2) )
1983  call check( __line__, pott(ke-1,i,j-1) )
1984  call check( __line__, pott(ke-1,i,j ) )
1985  call check( __line__, pott(ke-1,i,j+1) )
1986  call check( __line__, pott(ke-1,i,j+2) )
1987  call check( __line__, g(ke-2,i,j,i_xyw) )
1988  call check( __line__, g(ke ,i,j,i_xyw) )
1989  call check( __line__, g(ke-1,i,j,i_xyz) )
1990  call check( __line__, rcs2t(ke-2,i,j) )
1991  call check( __line__, rcs2t(ke-1,i,j) )
1992  call check( __line__, rcs2t(ke ,i,j) )
1993 #endif
1994  ! k,i,j
1995  m(1,ke-1,i,j) = &
1996  - ( 0.5_rp * (pott(ke ,i,j)+pott(ke-1,i,j)) * rfdz(ke-1) / g(ke-1,i,j,i_xyw) &
1997  + ( fact_n * (pott(ke-1,i,j)+pott(ke-2,i,j)) &
1998  + fact_f * (pott(ke ,i,j)+pott(ke-3,i,j)) ) * rfdz(ke-2) / g(ke-2,i,j,i_xyw) &
1999  ) * j33g * j33g * rfdz(ke-1) &
2000  - ( ( fact_n * (pott(ke-1,i+1,j)+pott(ke-1,i ,j)) &
2001  + fact_f * (pott(ke-1,i+2,j)+pott(ke-1,i-1,j)) ) * rfdx(i ) &
2002  + ( fact_n * (pott(ke-1,i ,j)+pott(ke-1,i-1,j)) &
2003  + fact_f * (pott(ke-1,i+1,j)+pott(ke-1,i-2,j)) ) * rfdx(i-1) &
2004  ) * g(ke-1,i,j,i_xyz) * rfdx(i) &
2005  - ( ( fact_n * (pott(ke-1,i,j+1)+pott(ke-1,i,j )) &
2006  + fact_f * (pott(ke-1,i,j+2)+pott(ke-1,i,j-1)) ) * rfdy(j ) &
2007  + ( fact_n * (pott(ke-1,i,j )+pott(ke-1,i,j-1)) &
2008  + fact_f * (pott(ke-1,i,j-1)+pott(ke-1,i,j-2)) ) * rfdy(j-1) &
2009  ) * g(ke-1,i,j,i_xyz) * rfdy(j) &
2010  - g(ke-1,i,j,i_xyz) * rcs2t(ke-1,i,j) * rdt * rdt
2011  ! k-1
2012  m(2,ke-1,i,j) = j33g * j33g / g(ke-2,i,j,i_xyw) &
2013  * ( fact_n * (pott(ke-1,i,j)+pott(ke-2,i,j)) &
2014  + fact_f * (pott(ke ,i,j)+pott(ke-3,i,j)) ) &
2015  * rfdz(ke-2) * rcdz(ke-1) &
2016  - grav * j33g * rcs2t(ke-2,i,j) / ( fdz(ke-1)+fdz(ke-2) )
2017  ! k+1
2018  m(3,ke-1,i,j) = j33g * j33g / g(ke ,i,j,i_xyw) &
2019  * 0.5_rp * (pott(ke ,i,j)+pott(ke-1,i,j)) &
2020  * rfdz(ke-1) * rcdz(ke-1) &
2021  + grav * j33g * rcs2t(ke ,i,j) / ( fdz(ke-1)+fdz(ke-2) )
2022 #ifdef DEBUG
2023  call check( __line__, pott(ke-1,i,j) )
2024  call check( __line__, pott(ke ,i,j) )
2025  call check( __line__, pott(ke,i-2,j) )
2026  call check( __line__, pott(ke,i-1,j) )
2027  call check( __line__, pott(ke,i ,j) )
2028  call check( __line__, pott(ke,i+1,j) )
2029  call check( __line__, pott(ke,i+2,j) )
2030  call check( __line__, pott(ke,i,j-2) )
2031  call check( __line__, pott(ke,i,j-1) )
2032  call check( __line__, pott(ke,i,j ) )
2033  call check( __line__, pott(ke,i,j+1) )
2034  call check( __line__, pott(ke,i,j+2) )
2035  call check( __line__, g(ke-1,i,j,i_xyw) )
2036  call check( __line__, g(ke,i,j,i_xyz) )
2037  call check( __line__, rcs2t(ke-1,i,j) )
2038  call check( __line__, rcs2t(ke,i,j) )
2039 #endif
2040  ! k,i,j
2041  m(1,ke,i,j) = &
2042  - ( &
2043  + 0.5_rp * (pott(ke ,i,j)+pott(ke-1,i,j)) * rfdz(ke-1) / g(ke-1,i,j,i_xyw) &
2044  ) * j33g * j33g * rfdz(ke) &
2045  - ( ( fact_n * (pott(ke,i+1,j)+pott(ke,i ,j)) &
2046  + fact_f * (pott(ke,i+2,j)+pott(ke,i-1,j)) ) * rfdx(i ) &
2047  + ( fact_n * (pott(ke,i ,j)+pott(ke,i-1,j)) &
2048  + fact_f * (pott(ke,i+1,j)+pott(ke,i-2,j)) ) * rfdx(i-1) &
2049  ) * g(ke,i,j,i_xyz) * rfdx(i) &
2050  - ( ( fact_n * (pott(ke,i,j+1)+pott(ke,i,j )) &
2051  + fact_f * (pott(ke,i,j+2)+pott(ke,i,j-1)) ) * rfdy(j ) &
2052  + ( fact_n * (pott(ke,i,j )+pott(ke,i,j-1)) &
2053  + fact_f * (pott(ke,i,j-1)+pott(ke,i,j-2)) ) * rfdy(j-1) &
2054  ) * g(ke,i,j,i_xyz) * rfdy(j) &
2055  - g(ke,i,j,i_xyz) * rcs2t(ke,i,j) * rdt * rdt &
2056  - grav * j33g * 0.5_rp * rcs2t(ke,i,j) * rcdz(ke)
2057  ! k-1
2058  m(2,ke,i,j) = j33g * j33g / g(ke-1,i,j,i_xyw) &
2059  * 0.5_rp * (pott(ke ,i,j)+pott(ke-1,i,j)) &
2060  * rfdz(ke-1) * rcdz(ke) &
2061  - grav * j33g * 0.5_rp * rcs2t(ke,i,j) * rcdz(ke)
2062  enddo
2063  enddo
2064 #ifdef DEBUG
2065  k = iundef; i = iundef; j = iundef
2066 #endif
2067 
2068  do j = jjs, jje
2069  do i = iis, iie
2070  do k = ks, ke
2071 #ifdef DEBUG
2072  call check( __line__, g(k,i-1,j,i_xyz) )
2073  call check( __line__, g(k,i+1,j,i_xyz) )
2074  call check( __line__, g(k,i,j-1,i_xyz) )
2075  call check( __line__, g(k,i,j+1,i_xyz) )
2076  call check( __line__, pott(k,i-2,j ) )
2077  call check( __line__, pott(k,i-1,j ) )
2078  call check( __line__, pott(k,i ,j ) )
2079  call check( __line__, pott(k,i+1,j ) )
2080  call check( __line__, pott(k,i+2,j ) )
2081  call check( __line__, pott(k,i ,j-2) )
2082  call check( __line__, pott(k,i ,j-1) )
2083  call check( __line__, pott(k,i ,j ) )
2084  call check( __line__, pott(k,i ,j+1) )
2085  call check( __line__, pott(k,i ,j+2) )
2086 #endif
2087  ! i-1
2088  m(4,k,i,j) = g(k,i-1,j,i_xyz) &
2089  * ( fact_n * (pott(k,i ,j)+pott(k,i-1,j)) &
2090  + fact_f * (pott(k,i+1,j)+pott(k,i-2,j)) ) &
2091  * rfdx(i-1) * rcdx(i)
2092  ! i+1
2093  m(5,k,i,j) = g(k,i+1,j,i_xyz) &
2094  * ( fact_n * (pott(k,i+1,j)+pott(k,i ,j)) &
2095  + fact_f * (pott(k,i+2,j)+pott(k,i-1,j)) ) &
2096  * rfdx(i ) * rcdx(i)
2097  ! j-1
2098  m(6,k,i,j) = g(k,i,j-1,i_xyz) &
2099  * ( fact_n * (pott(k,i,j )+pott(k,i,j-1)) &
2100  + fact_f * (pott(k,i,j+1)+pott(k,i,j-2)) ) &
2101  * rfdy(j-1) * rcdy(j)
2102  ! j+1
2103  m(7,k,i,j) = g(k,i,j+1,i_xyz) &
2104  * ( fact_n * (pott(k,i,j+1)+pott(k,i,j )) &
2105  + fact_f * (pott(k,i,j+2)+pott(k,i,j-1)) ) &
2106  * rfdy(j ) * rcdy(j)
2107  enddo
2108  enddo
2109  enddo
2110 #ifdef DEBUG
2111  k = iundef; i = iundef; j = iundef
2112 #endif
2113 
2114  return
2115  end subroutine make_matrix
2116 
2117  subroutine mul_matrix(V, M, C)
2118  implicit none
2119  real(RP), intent(out) :: V(ka,ia,ja)
2120  real(RP), intent(in) :: M(7,ka,ia,ja)
2121  real(RP), intent(in) :: C(ka,ia,ja)
2122 
2123  integer :: k, i, j
2124 
2125 #ifdef DEBUG
2126  v(:,:,:) = undef
2127 #endif
2128 
2129  do j = js, je
2130  do i = is, ie
2131  do k = ks+1, ke-1
2132 #ifdef DEBUG
2133  call check( __line__, m(1,k,i,j) )
2134  call check( __line__, m(2,k,i,j) )
2135  call check( __line__, m(3,k,i,j) )
2136  call check( __line__, m(4,k,i,j) )
2137  call check( __line__, m(5,k,i,j) )
2138  call check( __line__, m(6,k,i,j) )
2139  call check( __line__, m(7,k,i,j) )
2140  call check( __line__, c(k ,i ,j ) )
2141  call check( __line__, c(k-1,i ,j ) )
2142  call check( __line__, c(k+1,i ,j ) )
2143  call check( __line__, c(k ,i-1,j ) )
2144  call check( __line__, c(k ,i+1,j ) )
2145  call check( __line__, c(k ,i ,j-1) )
2146  call check( __line__, c(k ,i ,j+1) )
2147 #endif
2148  v(k,i,j) = m(1,k,i,j) * c(k ,i ,j ) &
2149  + m(2,k,i,j) * c(k-1,i ,j ) &
2150  + m(3,k,i,j) * c(k+1,i ,j ) &
2151  + m(4,k,i,j) * c(k ,i-1,j ) &
2152  + m(5,k,i,j) * c(k ,i+1,j ) &
2153  + m(6,k,i,j) * c(k ,i ,j-1) &
2154  + m(7,k,i,j) * c(k ,i ,j+1)
2155  enddo
2156  enddo
2157  enddo
2158 #ifdef DEBUG
2159  k = iundef; i = iundef; j = iundef
2160 #endif
2161  do j = js, je
2162  do i = is, ie
2163 #ifdef DEBUG
2164  call check( __line__, m(1,ks,i,j) )
2165  call check( __line__, m(3,ks,i,j) )
2166  call check( __line__, m(4,ks,i,j) )
2167  call check( __line__, m(5,ks,i,j) )
2168  call check( __line__, m(6,ks,i,j) )
2169  call check( __line__, m(7,ks,i,j) )
2170  call check( __line__, c(ks ,i ,j ) )
2171  call check( __line__, c(ks+1,i ,j ) )
2172  call check( __line__, c(ks ,i-1,j ) )
2173  call check( __line__, c(ks ,i+1,j ) )
2174  call check( __line__, c(ks ,i ,j-1) )
2175  call check( __line__, c(ks ,i ,j+1) )
2176  call check( __line__, m(1,ke,i,j) )
2177  call check( __line__, m(2,ke,i,j) )
2178  call check( __line__, m(4,ke,i,j) )
2179  call check( __line__, m(5,ke,i,j) )
2180  call check( __line__, m(6,ke,i,j) )
2181  call check( __line__, m(7,ke,i,j) )
2182  call check( __line__, c(ke ,i ,j ) )
2183  call check( __line__, c(ke-1,i ,j ) )
2184  call check( __line__, c(ke ,i-1,j ) )
2185  call check( __line__, c(ke ,i+1,j ) )
2186  call check( __line__, c(ke ,i ,j-1) )
2187  call check( __line__, c(ke ,i ,j+1) )
2188 #endif
2189  v(ks,i,j) = m(1,ks,i,j) * c(ks ,i ,j ) &
2190  + m(3,ks,i,j) * c(ks+1,i ,j ) &
2191  + m(4,ks,i,j) * c(ks ,i-1,j ) &
2192  + m(5,ks,i,j) * c(ks ,i+1,j ) &
2193  + m(6,ks,i,j) * c(ks ,i ,j-1) &
2194  + m(7,ks,i,j) * c(ks ,i ,j+1)
2195  v(ke,i,j) = m(1,ke,i,j) * c(ke ,i ,j ) &
2196  + m(2,ke,i,j) * c(ke-1,i ,j ) &
2197  + m(4,ke,i,j) * c(ke ,i-1,j ) &
2198  + m(5,ke,i,j) * c(ke ,i+1,j ) &
2199  + m(6,ke,i,j) * c(ke ,i ,j-1) &
2200  + m(7,ke,i,j) * c(ke ,i ,j+1)
2201  enddo
2202  enddo
2203 #ifdef DEBUG
2204  k = iundef; i = iundef; j = iundef
2205 #endif
2206 
2207  return
2208  end subroutine mul_matrix
2209 
2210 #else
2211  subroutine solve_multigrid
2212  end subroutine solve_multigrid
2213 #endif
2214 
2215 #ifdef DEBUG
2216  subroutine check_solver( &
2217  DPRES, M, B )
2218  use scale_prc, only: &
2219  prc_abort
2220  real(RP), intent(in) :: DPRES(ka,ia,ja)
2221  real(RP), intent(in) :: M(7,ka,ia,ja)
2222  real(RP), intent(in) :: B(ka,ia,ja)
2223 
2224  real(RP) :: B2(ka,ia,ja)
2225  integer :: k, i, j
2226  real(RP) :: err
2227 
2228  call mul_matrix(b2, m, dpres)
2229 
2230  do j = js, je
2231  do i = is, ie
2232  do k = ks, ke
2233  err = abs( b2(k,i,j) - b(k,i,j) )
2234  if ( err > 1.e-5_rp .and. abs( err / b(k,i,j) ) > 1.e-5_rp ) then
2235  log_error("check_solver",*) "solver error is too large: ", k,i,j, b(k,i,j), b2(k,i,j)
2236  call prc_abort
2237  endif
2238  enddo
2239  enddo
2240  enddo
2241 
2242  end subroutine check_solver
2243 
2244  subroutine check_pres( &
2245  DPRES_N, DPRES, &
2246  RHOT_RK, RHOT, &
2247  DENS_RK, DENS, &
2248  B, &
2249  RT2P )
2250  use scale_prc, only: &
2251  prc_abort
2252  real(RP), intent(in) :: DPRES_N(ka,ia,ja)
2253  real(RP), intent(in) :: DPRES(ka,ia,ja)
2254  real(RP), intent(in) :: RHOT_RK(ka,ia,ja)
2255  real(RP), intent(in) :: RHOT(ka,ia,ja)
2256  real(RP), intent(in) :: DENS_RK(ka,ia,ja)
2257  real(RP), intent(in) :: DENS(ka,ia,ja)
2258  real(RP), intent(in) :: B(ka,ia,ja)
2259  real(RP), intent(in) :: RT2P(ka,ia,ja)
2260 
2261  real(RP) :: lhs, rhs
2262  integer :: k,i,j
2263 
2264  do j = js, je
2265  do i = is, ie
2266  do k = ks, ke
2267  lhs = dpres_n(k,i,j) - dpres(k,i,j)
2268  rhs = rt2p(k,i,j) * ( rhot_rk(k,i,j) - rhot(k,i,j) )
2269  if ( abs( (lhs - rhs) / lhs ) > 1e-15 ) then
2270  log_error("check_pres",*) "error is too large: ", k,i,j, lhs, rhs, &
2271  dpres_n(k,i,j),dpres(k,i,j),rhot_rk(k,i,j),rhot(k,i,j),dens_rk(k,i,j),dens(k,i,j),b(k,i,j)
2272  call prc_abort
2273  endif
2274  enddo
2275  enddo
2276  enddo
2277  end subroutine check_pres
2278 #endif
2279 
2280 end module scale_atmos_dyn_tstep_short_fvm_hivi
scale_index::i_rhot
integer, parameter, public i_rhot
Definition: scale_index.F90:32
scale_debug
module DEBUG
Definition: scale_debug.F90:11
scale_atmos_grid_cartesc_index::i_xyz
integer, public i_xyz
Definition: scale_atmos_grid_cartesC_index.F90:90
scale_comm_cartesc::comm_world
integer, public comm_world
communication world ID
Definition: scale_comm_cartesC.F90:94
scale_index::va
integer, public va
Definition: scale_index.F90:35
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:236
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:225
scale_atmos_grid_cartesc_index::ia
integer, public ia
of whole cells: x, local, with HALO
Definition: scale_atmos_grid_cartesC_index.F90:48
scale_index::i_momx
integer, parameter, public i_momx
Definition: scale_index.F90:30
scale_index::i_momz
integer, parameter, public i_momz
Definition: scale_index.F90:29
scale_atmos_grid_cartesc_index::i_xy
integer, public i_xy
Definition: scale_atmos_grid_cartesC_index.F90:98
scale_atmos_grid_cartesc_index::i_xvw
integer, public i_xvw
Definition: scale_atmos_grid_cartesC_index.F90:93
scale_atmos_grid_cartesc_index::i_xvz
integer, public i_xvz
Definition: scale_atmos_grid_cartesC_index.F90:95
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: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:165
scale_atmos_dyn_tstep_short_fvm_hivi
module Atmosphere / Dynamics RK
Definition: scale_atmos_dyn_tstep_short_fvm_hivi.F90:16
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_fluxj13_xyw
procedure(flux_j), pointer, public atmos_dyn_fvm_fluxj13_xyw
Definition: scale_atmos_dyn_fvm_flux.F90:193
scale_atmos_grid_cartesc_index::ja
integer, public ja
of whole cells: y, local, with HALO
Definition: scale_atmos_grid_cartesC_index.F90:49
scale_io::io_fid_conf
integer, public io_fid_conf
Config file ID.
Definition: scale_io.F90:55
scale_atmos_dyn_tstep_short_fvm_hivi::atmos_dyn_tstep_short_fvm_hivi_regist
subroutine, public atmos_dyn_tstep_short_fvm_hivi_regist(ATMOS_DYN_TYPE, VA_out, VAR_NAME, VAR_DESC, VAR_UNIT)
Register.
Definition: scale_atmos_dyn_tstep_short_fvm_hivi.F90:79
scale_debug::check
subroutine, public check(current_line, v)
Undefined value checker.
Definition: scale_debug.F90:56
scale_index::i_dens
integer, parameter, public i_dens
Definition: scale_index.F90:28
scale_atmos_dyn_tstep_short_fvm_hivi::atmos_dyn_tstep_short_fvm_hivi
subroutine, public atmos_dyn_tstep_short_fvm_hivi(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, dtrk, last)
Definition: scale_atmos_dyn_tstep_short_fvm_hivi.F90:179
scale_index::i_momy
integer, parameter, public i_momy
Definition: scale_index.F90:31
scale_const::const_undef
real(rp), public const_undef
Definition: scale_const.F90:41
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:182
scale_comm_cartesc
module COMMUNICATION
Definition: scale_comm_cartesC.F90:11
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_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:187
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_tracer
module TRACER
Definition: scale_tracer.F90:16
scale_index
module Index
Definition: scale_index.F90:11
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:231
scale_atmos_grid_cartesc_index::ydir
integer, parameter, public ydir
Definition: scale_atmos_grid_cartesC_index.F90:34
scale_atmos_grid_cartesc_index::i_uvz
integer, public i_uvz
Definition: scale_atmos_grid_cartesC_index.F90:96
scale_atmos_grid_cartesc_index
module atmosphere / grid / cartesC index
Definition: scale_atmos_grid_cartesC_index.F90:12
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_tstep_short_fvm_hivi::solve_bicgstab
subroutine solve_bicgstab(DPRES_N, DPRES, M, B)
Definition: scale_atmos_dyn_tstep_short_fvm_hivi.F90:1465
scale_const::const_pre00
real(rp), public const_pre00
pressure reference [Pa]
Definition: scale_const.F90:88
scale_prc
module PROCESS
Definition: scale_prc.F90:11
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_tstep_short_fvm_hivi::atmos_dyn_tstep_short_fvm_hivi_setup
subroutine, public atmos_dyn_tstep_short_fvm_hivi_setup
Setup.
Definition: scale_atmos_dyn_tstep_short_fvm_hivi.F90:112
scale_const::const_grav
real(rp), public const_grav
standard acceleration of gravity [m/s2]
Definition: scale_const.F90:46
scale_const::const_undef2
integer, parameter, public const_undef2
undefined value (INT2)
Definition: scale_const.F90:38
scale_atmos_dyn_common
module Atmosphere / Dynamics common
Definition: scale_atmos_dyn_common.F90:12
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_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_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:338
scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_flux_valuew_z
procedure(valuew), pointer, public atmos_dyn_fvm_flux_valuew_z
Definition: scale_atmos_dyn_fvm_flux.F90:156
scale_const
module CONSTANT
Definition: scale_const.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::xdir
integer, parameter, public xdir
Definition: scale_atmos_grid_cartesC_index.F90:33
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:247
scale_atmos_grid_cartesc_index::i_uyw
integer, public i_uyw
Definition: scale_atmos_grid_cartesC_index.F90:92
scale_atmos_grid_cartesc_index::i_uy
integer, public i_uy
Definition: scale_atmos_grid_cartesC_index.F90:99
scale_atmos_grid_cartesc_index::i_uyz
integer, public i_uyz
Definition: scale_atmos_grid_cartesC_index.F90:94
scale_prof
module profiler
Definition: scale_prof.F90:11
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:204
scale_atmos_dyn_fvm_flux
module scale_atmos_dyn_fvm_flux
Definition: scale_atmos_dyn_fvm_flux.F90:13
scale_precision
module PRECISION
Definition: scale_precision.F90:14
scale_atmos_grid_cartesc_index::ka
integer, public ka
of whole cells: z, local, with HALO
Definition: scale_atmos_grid_cartesC_index.F90:47
scale_atmos_dyn_common::atmos_dyn_fct
subroutine, public atmos_dyn_fct(qflx_anti, phi_in, DENS0, DENS, qflx_hi, qflx_lo, mflx_hi, rdz, rdx, rdy, GSQRT, MAPF, dt, flag_vect)
Flux Correction Transport Limiter.
Definition: scale_atmos_dyn_common.F90:1930
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:220
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:252
scale_atmos_grid_cartesc_index::i_uv
integer, public i_uv
Definition: scale_atmos_grid_cartesC_index.F90:101
scale_precision::sp
integer, parameter, public sp
Definition: scale_precision.F90:31
scale_io
module STDIO
Definition: scale_io.F90:10
scale_atmos_dyn_tstep_short_fvm_hivi::make_matrix
subroutine make_matrix(M, POTT, RCs2T, GRAV, G, J33G, RCDZ, RFDZ, RCDX, RFDX, RCDY, RFDY, FDZ, rdt, FACT_N, FACT_F, I_XYZ, I_XYW, IIS, IIE, JJS, JJE)
Definition: scale_atmos_dyn_tstep_short_fvm_hivi.F90:1791
scale_precision::rp
integer, parameter, public rp
Definition: scale_precision.F90:41
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:160
scale_atmos_grid_cartesc_index::zdir
integer, parameter, public zdir
Definition: scale_atmos_grid_cartesC_index.F90:32
scale_precision::dp
integer, parameter, public dp
Definition: scale_precision.F90:32
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:209
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:170
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:214
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:177
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:241
scale_atmos_grid_cartesc_index::i_xyw
integer, public i_xyw
Definition: scale_atmos_grid_cartesC_index.F90:91
Generated by   doxygen 1.8.13