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