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