SCALE-RM
Functions/Subroutines
scale_atmos_dyn_tstep_short_fvm_hevi Module Reference

module Atmosphere / Dynamics RK More...

Functions/Subroutines

subroutine, public atmos_dyn_tstep_short_fvm_hevi_regist (ATMOS_DYN_TYPE, VA_out, VAR_NAME, VAR_DESC, VAR_UNIT)
 Register. More...
 
subroutine, public atmos_dyn_tstep_short_fvm_hevi_setup
 Setup. More...
 
subroutine, public atmos_dyn_tstep_short_fvm_hevi (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)
 
subroutine solve_direct (C, F1, F2, F3)
 

Detailed Description

module Atmosphere / Dynamics RK

Description
HEVI FVM scheme for Atmospheric dynamical process
Author
Team SCALE
History
  • 2013-06-18 (S.Nishizawa) [new] newly impremented
  • 2014-04-04 (S.Nishizawa) [mod] support terrain-following coordinate
NAMELIST
  • No namelist group
History Output
namedescriptionunitvariable
DENS_t_advch tendency of density (horiz. advection) kg/m3/s advch_t
DENS_t_advcv tendency of density (vert. advection) kg/m3/s advcv_t
MOMX_t_advch tendency of momentum x (horiz. advection) kg/m2/s2 advch_t
MOMX_t_advcv tendency of momentum x (vert. advection) kg/m2/s2 advcv_t
MOMX_t_cf tendency of momentum x (coliolis force) kg/m2/s2 cf_t
MOMX_t_ddiv tendency of momentum x (divergence damping) kg/m2/s2 ddiv_t
MOMX_t_pg tendency of momentum x (pressure gradient) kg/m2/s2 pg_t
MOMY_t_advch tendency of momentum y (horiz. advection) kg/m2/s2 advch_t
MOMY_t_advcv tendency of momentum y (vert. advection) kg/m2/s2 advcv_t
MOMY_t_cf tendency of momentum y (coliolis force) kg/m2/s2 cf_t
MOMY_t_ddiv tendency of momentum y (divergence damping) kg/m2/s2 ddiv_t
MOMY_t_pg tendency of momentum y (pressure gradient) kg/m2/s2 pg_t
MOMZ_t_advch tendency of momentum z (horiz. advection) kg/m2/s2 advch_t
MOMZ_t_advcv tendency of momentum z (vert. advection) kg/m2/s2 advcv_t
MOMZ_t_ddiv tendency of momentum z (divergence damping) kg/m2/s2 ddiv_t
MOMZ_t_pg tendency of momentum z (pressure gradient) kg/m2/s2 pg_t
MOMZ_t_wdamp tendency of momentum z (Rayleigh damping) kg/m2/s2 wdmp_t
RHOT_t_advch tendency of rho*theta (horiz. advection) K kg/m3/s advch_t
RHOT_t_advcv tendency of rho*theta (vert. advection) K kg/m3/s advcv_t

Function/Subroutine Documentation

◆ atmos_dyn_tstep_short_fvm_hevi_regist()

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

Register.

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

Definition at line 90 of file scale_atmos_dyn_tstep_short_fvm_hevi.F90.

References scale_stdio::io_fid_log, scale_stdio::io_l, and scale_process::prc_mpistop().

Referenced by scale_atmos_dyn_tstep_short::atmos_dyn_tstep_short_regist().

90  use scale_process, only: &
91  prc_mpistop
92  implicit none
93 
94  character(len=*), intent(in) :: ATMOS_DYN_TYPE
95  integer, intent(out) :: VA_out
96  character(len=H_SHORT), intent(out) :: VAR_NAME(:)
97  character(len=H_MID), intent(out) :: VAR_DESC(:)
98  character(len=H_SHORT), intent(out) :: VAR_UNIT(:)
99  !---------------------------------------------------------------------------
100 
101  if ( atmos_dyn_type /= 'FVM-HEVI' .AND. atmos_dyn_type /= 'HEVI' ) then
102  write(*,*) 'xxx ATMOS_DYN_TYPE is not FVM-HEVI. Check!'
103  call prc_mpistop
104  endif
105 
106  va_out = va_fvm_hevi
107  var_name(:) = ""
108  var_desc(:) = ""
109  var_unit(:) = ""
110 
111  if( io_l ) write(io_fid_log,*)
112  if( io_l ) write(io_fid_log,*) '*** Register additional prognostic variables (HEVI)'
113  if ( va_out < 1 ) then
114  if( io_l ) write(io_fid_log,*) '*** => nothing.'
115  endif
116 
117  return
scale_process::prc_mpistop
subroutine, public prc_mpistop
Abort MPI.
Definition: scale_process.F90:267
scale_process
module PROCESS
Definition: scale_process.F90:14
Here is the call graph for this function:
Here is the caller graph for this function:

◆ atmos_dyn_tstep_short_fvm_hevi_setup()

subroutine, public scale_atmos_dyn_tstep_short_fvm_hevi::atmos_dyn_tstep_short_fvm_hevi_setup ( )

Setup.

Definition at line 123 of file scale_atmos_dyn_tstep_short_fvm_hevi.F90.

References scale_stdio::io_fid_log, and scale_stdio::io_l.

Referenced by scale_atmos_dyn_tstep_short::atmos_dyn_tstep_short_regist().

123  implicit none
124  !---------------------------------------------------------------------------
125 
126  if( io_l ) write(io_fid_log,*) '*** HEVI Setup'
127 
128  return
Here is the caller graph for this function:

◆ atmos_dyn_tstep_short_fvm_hevi()

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

Definition at line 150 of file scale_atmos_dyn_tstep_short_fvm_hevi.F90.

References scale_atmos_dyn_common::atmos_dyn_fct(), scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_flux_valuew_z, scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxj13_uyz, scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxj13_xvz, scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxj13_xyw, scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxj23_uyz, scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxj23_xvz, scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxj23_xyw, scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxx_uyz, scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxx_xvz, scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxx_xyw, scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxx_xyz, scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxy_uyz, scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxy_xvz, scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxy_xyw, scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxy_xyz, scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxz_uyz, scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxz_xvz, scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxz_xyw, scale_atmos_dyn_fvm_flux::atmos_dyn_fvm_fluxz_xyz, scale_debug::check(), scale_index::i_dens, scale_index::i_momx, scale_index::i_momy, scale_index::i_momz, scale_index::i_rhot, scale_gridtrans::i_uv, scale_gridtrans::i_uvz, scale_gridtrans::i_uy, scale_gridtrans::i_uyw, scale_gridtrans::i_uyz, scale_gridtrans::i_xv, scale_gridtrans::i_xvw, scale_gridtrans::i_xvz, scale_gridtrans::i_xy, scale_gridtrans::i_xyw, scale_gridtrans::i_xyz, scale_grid_index::iblock, scale_grid_index::ie, scale_grid_index::ieh, scale_grid_index::ihalo, scale_grid_index::is, scale_grid_index::jblock, scale_grid_index::je, scale_grid_index::jeh, scale_grid_index::jhalo, scale_grid_index::js, scale_grid_index::ka, scale_grid_index::ke, scale_grid_index::ks, solve_direct(), scale_grid_index::xdir, scale_grid_index::ydir, and scale_grid_index::zdir.

Referenced by scale_atmos_dyn_tstep_short::atmos_dyn_tstep_short_regist().

150  use scale_grid_index
151  use scale_const, only: &
152 #ifdef dry
153  rdry => const_rdry, &
154  cvdry => const_cvdry, &
155  cpdry => const_cpdry, &
156 #endif
157  eps => const_eps, &
158  grav => const_grav, &
159  p00 => const_pre00
160  use scale_atmos_dyn_common, only: &
161  atmos_dyn_fct
162  use scale_atmos_dyn_fvm_flux, only: &
163  atmos_dyn_fvm_flux_valuew_z, &
164  atmos_dyn_fvm_fluxz_xyz, &
165  atmos_dyn_fvm_fluxx_xyz, &
166  atmos_dyn_fvm_fluxy_xyz, &
167  atmos_dyn_fvm_fluxz_xyw, &
168  atmos_dyn_fvm_fluxj13_xyw, &
169  atmos_dyn_fvm_fluxj23_xyw, &
170  atmos_dyn_fvm_fluxx_xyw, &
171  atmos_dyn_fvm_fluxy_xyw, &
172  atmos_dyn_fvm_fluxz_uyz, &
173  atmos_dyn_fvm_fluxj13_uyz, &
174  atmos_dyn_fvm_fluxj23_uyz, &
175  atmos_dyn_fvm_fluxx_uyz, &
176  atmos_dyn_fvm_fluxy_uyz, &
177  atmos_dyn_fvm_fluxz_xvz, &
178  atmos_dyn_fvm_fluxj13_xvz, &
179  atmos_dyn_fvm_fluxj23_xvz, &
180  atmos_dyn_fvm_fluxx_xvz, &
181  atmos_dyn_fvm_fluxy_xvz
182  use scale_gridtrans, only: &
183  i_xyz, &
184  i_xyw, &
185  i_uyw, &
186  i_xvw, &
187  i_uyz, &
188  i_xvz, &
189  i_uvz, &
190  i_xy , &
191  i_uy , &
192  i_xv , &
193  i_uv
194 #ifdef HIST_TEND
195  use scale_history, only: &
196  hist_in
197 #endif
198  implicit none
199 
200  real(RP), intent(out) :: DENS_RK(KA,IA,JA) ! prognostic variables
201  real(RP), intent(out) :: MOMZ_RK(KA,IA,JA) !
202  real(RP), intent(out) :: MOMX_RK(KA,IA,JA) !
203  real(RP), intent(out) :: MOMY_RK(KA,IA,JA) !
204  real(RP), intent(out) :: RHOT_RK(KA,IA,JA) !
205 
206  real(RP), intent(out) :: PROG_RK(KA,IA,JA,VA) !
207 
208  real(RP), intent(inout) :: mflx_hi(KA,IA,JA,3) ! rho * vel(x,y,z)
209  real(RP), intent(out) :: tflx_hi(KA,IA,JA,3) ! rho * theta * vel(x,y,z)
210 
211  real(RP), intent(in),target :: DENS0(KA,IA,JA) ! prognostic variables
212  real(RP), intent(in),target :: MOMZ0(KA,IA,JA) ! at previous dynamical time step
213  real(RP), intent(in),target :: MOMX0(KA,IA,JA) !
214  real(RP), intent(in),target :: MOMY0(KA,IA,JA) !
215  real(RP), intent(in),target :: RHOT0(KA,IA,JA) !
216 
217  real(RP), intent(in) :: DENS(KA,IA,JA) ! prognostic variables
218  real(RP), intent(in) :: MOMZ(KA,IA,JA) ! at previous RK step
219  real(RP), intent(in) :: MOMX(KA,IA,JA) !
220  real(RP), intent(in) :: MOMY(KA,IA,JA) !
221  real(RP), intent(in) :: RHOT(KA,IA,JA) !
222 
223  real(RP), intent(in) :: DENS_t(KA,IA,JA)
224  real(RP), intent(in) :: MOMZ_t(KA,IA,JA)
225  real(RP), intent(in) :: MOMX_t(KA,IA,JA)
226  real(RP), intent(in) :: MOMY_t(KA,IA,JA)
227  real(RP), intent(in) :: RHOT_t(KA,IA,JA)
228 
229  real(RP), intent(in) :: PROG0(KA,IA,JA,VA)
230  real(RP), intent(in) :: PROG (KA,IA,JA,VA)
231 
232  real(RP), intent(in) :: DPRES0(KA,IA,JA)
233  real(RP), intent(in) :: RT2P(KA,IA,JA)
234  real(RP), intent(in) :: CORIOLI(1,IA,JA)
235  real(RP), intent(in) :: num_diff(KA,IA,JA,5,3)
236  real(RP), intent(in) :: wdamp_coef(KA)
237  real(RP), intent(in) :: divdmp_coef
238  real(RP), intent(in) :: DDIV(KA,IA,JA)
239 
240  logical, intent(in) :: FLAG_FCT_MOMENTUM
241  logical, intent(in) :: FLAG_FCT_T
242  logical, intent(in) :: FLAG_FCT_ALONG_STREAM
243 
244  real(RP), intent(in) :: CDZ(KA)
245  real(RP), intent(in) :: FDZ(KA-1)
246  real(RP), intent(in) :: FDX(IA-1)
247  real(RP), intent(in) :: FDY(JA-1)
248  real(RP), intent(in) :: RCDZ(KA)
249  real(RP), intent(in) :: RCDX(IA)
250  real(RP), intent(in) :: RCDY(JA)
251  real(RP), intent(in) :: RFDZ(KA-1)
252  real(RP), intent(in) :: RFDX(IA-1)
253  real(RP), intent(in) :: RFDY(JA-1)
254 
255  real(RP), intent(in) :: PHI (KA,IA,JA)
256  real(RP), intent(in) :: GSQRT (KA,IA,JA,7)
257  real(RP), intent(in) :: J13G (KA,IA,JA,7)
258  real(RP), intent(in) :: J23G (KA,IA,JA,7)
259  real(RP), intent(in) :: J33G
260  real(RP), intent(in) :: MAPF (IA,JA,2,4)
261  real(RP), intent(in) :: REF_dens(KA,IA,JA)
262  real(RP), intent(in) :: REF_rhot(KA,IA,JA)
263 
264  logical, intent(in) :: BND_W
265  logical, intent(in) :: BND_E
266  logical, intent(in) :: BND_S
267  logical, intent(in) :: BND_N
268 
269  real(RP), intent(in) :: dtrk
270  logical, intent(in) :: last
271 
272 
273  ! diagnostic variables (work space)
274  real(RP) :: POTT(KA,IA,JA) ! potential temperature [K]
275  real(RP) :: DPRES(KA,IA,JA) ! pressure deviation from reference pressure
276 
277  real(RP) :: qflx_hi (KA,IA,JA,3)
278  real(RP) :: qflx_J13(KA,IA,JA)
279  real(RP) :: qflx_J23(KA,IA,JA)
280 
281  real(RP) :: advch ! horizontal advection
282  real(RP) :: advcv ! vertical advection
283  real(RP) :: wdmp ! Raileight damping
284  real(RP) :: div ! divergence damping
285  real(RP) :: pg ! pressure gradient force
286  real(RP) :: cf ! colioris force
287 #ifdef HIST_TEND
288  real(RP) :: advch_t(KA,IA,JA,5)
289  real(RP) :: advcv_t(KA,IA,JA,5)
290  real(RP) :: wdmp_t(KA,IA,JA)
291  real(RP) :: ddiv_t(KA,IA,JA,3)
292  real(RP) :: pg_t(KA,IA,JA,3)
293  real(RP) :: cf_t(KA,IA,JA,2)
294  logical :: lhist
295 #endif
296 
297  ! for implicit solver
298  real(RP) :: A(KA,IA,JA)
299  real(RP) :: B
300  real(RP) :: Sr(KA,IA,JA)
301  real(RP) :: Sw(KA,IA,JA)
302  real(RP) :: St(KA,IA,JA)
303  real(RP) :: PT(KA,IA,JA)
304  real(RP) :: C(KMAX-1,IA,JA)
305 
306  real(RP) :: F1(KA,IA,JA)
307  real(RP) :: F2(KA,IA,JA)
308  real(RP) :: F3(KA,IA,JA)
309 
310  integer :: IIS, IIE, JJS, JJE
311  integer :: k, i, j
312  integer :: iss, iee
313 
314 #ifdef DEBUG
315  pott(:,:,:) = undef
316  dpres(:,:,:) = undef
317 
318  pt(:,:,:) = undef
319 
320  qflx_hi(:,:,:,:) = undef
321  qflx_j13(:,:,:) = undef
322  qflx_j23(:,:,:) = undef
323 #endif
324 
325 #if defined DEBUG || defined QUICKDEBUG
326  dens_rk( 1:ks-1,:,:) = undef
327  dens_rk(ke+1:ka ,:,:) = undef
328  momz_rk( 1:ks-1,:,:) = undef
329  momz_rk(ke+1:ka ,:,:) = undef
330  momx_rk( 1:ks-1,:,:) = undef
331  momx_rk(ke+1:ka ,:,:) = undef
332  momy_rk( 1:ks-1,:,:) = undef
333  momy_rk(ke+1:ka ,:,:) = undef
334  rhot_rk( 1:ks-1,:,:) = undef
335  rhot_rk(ke+1:ka ,:,:) = undef
336  prog_rk( 1:ks-1,:,:,:) = undef
337  prog_rk(ke+1:ka ,:,:,:) = undef
338 #endif
339 
340 #ifdef HIST_TEND
341  lhist = last
342 #endif
343 
344 #ifdef PROFILE_FIPP
345  call fipp_start()
346 #endif
347 
348  ifs_off = 1
349  jfs_off = 1
350  if ( bnd_w ) ifs_off = 0
351  if ( bnd_s ) jfs_off = 0
352 
353 
354  do jjs = js, je, jblock
355  jje = jjs+jblock-1
356  do iis = is, ie, iblock
357  iie = iis+iblock-1
358 
359  profile_start("hevi_pres")
360  !$omp parallel do default(none) private(i,j,k) OMP_SCHEDULE_ collapse(2) &
361  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,DPRES0,RT2P,RHOT,REF_rhot,DPRES,DENS,PHI)
362  do j = jjs, jje+1
363  do i = iis, iie+1
364  do k = ks, ke
365 #ifdef DEBUG
366  call check( __line__, dpres0(k,i,j) )
367  call check( __line__, rt2p(k,i,j) )
368  call check( __line__, rhot(k,i,j) )
369  call check( __line__, ref_rhot(k,i,j) )
370 #endif
371  dpres(k,i,j) = dpres0(k,i,j) + rt2p(k,i,j) * ( rhot(k,i,j) - ref_rhot(k,i,j) )
372  enddo
373  dpres(ks-1,i,j) = dpres0(ks-1,i,j) - dens(ks,i,j) * ( phi(ks-1,i,j) - phi(ks+1,i,j) )
374  dpres(ke+1,i,j) = dpres0(ke+1,i,j) - dens(ke,i,j) * ( phi(ke+1,i,j) - phi(ke-1,i,j) )
375  enddo
376  enddo
377  profile_stop("hevi_pres")
378 
379  !##### continuity equation #####
380 
381  profile_start("hevi_mflx_z")
382  ! at (x, y, w)
383  !$omp parallel do default(none) private(i,j,k) OMP_SCHEDULE_ collapse(2) &
384  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,GSQRT,I_XYW,MOMX,MOMY,J13G,J23G,mflx_hi,MAPF,I_XY,num_diff)
385  do j = jjs, jje
386  do i = iis-1, iie
387  mflx_hi(ks-1,i,j,zdir) = 0.0_rp
388  do k = ks, ke-1
389 #ifdef DEBUG
390  call check( __line__, momx(k+1,i ,j) )
391  call check( __line__, momx(k+1,i-1,j) )
392  call check( __line__, momx(k ,i ,j) )
393  call check( __line__, momx(k ,i+1,j) )
394  call check( __line__, momy(k+1,i,j) )
395  call check( __line__, momy(k+1,i,j-1) )
396  call check( __line__, momy(k ,i,j) )
397  call check( __line__, momy(k ,i,j-1) )
398 #endif
399  mflx_hi(k,i,j,zdir) = j13g(k,i,j,i_xyw) * 0.25_rp / mapf(i,j,2,i_xy) &
400  * ( momx(k+1,i,j)+momx(k+1,i-1,j) &
401  + momx(k ,i,j)+momx(k ,i-1,j) ) &
402  + j23g(k,i,j,i_xyw) * 0.25_rp / mapf(i,j,1,i_xy) &
403  * ( momy(k+1,i,j)+momy(k+1,i,j-1) &
404  + momy(k ,i,j)+momy(k ,i,j-1) ) &
405  + gsqrt(k,i,j,i_xyw) / ( mapf(i,j,1,i_xy)*mapf(i,j,2,i_xy) ) * num_diff(k,i,j,i_dens,zdir)
406  enddo
407  mflx_hi(ke ,i,j,zdir) = 0.0_rp
408  enddo
409  enddo
410 #ifdef DEBUG
411  k = iundef; i = iundef; j = iundef
412 #endif
413  profile_stop("hevi_mflx_z")
414 
415  profile_start("hevi_mflx_x")
416  iss = max(iis-1,is-ifs_off)
417  iee = min(iie,ieh)
418  ! at (u, y, z)
419  !$omp parallel do default(none) private(i,j,k) OMP_SCHEDULE_ collapse(2) &
420  !$omp shared(JJS,JJE,iss,iee,KS,KE,GSQRT,I_UYZ,MOMX,num_diff,mflx_hi,MAPF,I_UY)
421  do j = jjs, jje
422  do i = iss, iee
423  do k = ks, ke
424 #ifdef DEBUG
425  call check( __line__, gsqrt(k,i,j,i_uyz) )
426  call check( __line__, momx(k,i,j) )
427  call check( __line__, num_diff(k,i,j,i_dens,xdir) )
428 #endif
429  mflx_hi(k,i,j,xdir) = gsqrt(k,i,j,i_uyz) / mapf(i,j,2,i_uy) &
430  * ( momx(k,i,j) + num_diff(k,i,j,i_dens,xdir) )
431  enddo
432  enddo
433  enddo
434 #ifdef DEBUG
435  k = iundef; i = iundef; j = iundef
436 #endif
437  profile_stop("hevi_mflx_x")
438 
439  ! at (x, v, z)
440  !$omp parallel do default(none) private(i,j,k) OMP_SCHEDULE_ collapse(2) &
441  !$omp shared(JJS,JS,JFS_OFF,JJE,JEH,IIS,IIE,KS,KE,GSQRT,I_XVZ,MOMY,num_diff,mflx_hi,MAPF,I_XV)
442  do j = max(jjs-1,js-jfs_off), min(jje,jeh)
443  do i = iis, iie
444  do k = ks, ke
445 #ifdef DEBUG
446  call check( __line__, gsqrt(k,i,j,i_xvz) )
447  call check( __line__, momy(k,i,j) )
448  call check( __line__, num_diff(k,i,j,i_dens,ydir) )
449 #endif
450  mflx_hi(k,i,j,ydir) = gsqrt(k,i,j,i_xvz) / mapf(i,j,1,i_xv) &
451  * ( momy(k,i,j) + num_diff(k,i,j,i_dens,ydir) )
452  enddo
453  enddo
454  enddo
455 #ifdef DEBUG
456  k = iundef; i = iundef; j = iundef
457 #endif
458 
459  !--- update density
460  profile_start("hevi_sr")
461  !$omp parallel do default(none) OMP_SCHEDULE_ collapse(2) &
462  !$omp private(i,j,k,advcv,advch) &
463 #ifdef HIST_TEND
464  !$omp shared(advcv_t,advch_t,lhist) &
465 #endif
466  !$omp shared(JJS,JJE,IIS,IIE,KS,KE) &
467  !$omp shared(DENS0,Sr,mflx_hi,DENS_t) &
468  !$omp shared(RCDZ,RCDX,RCDY,MAPF,GSQRT,I_XY,I_XYZ)
469  do j = jjs, jje
470  do i = iis, iie
471  do k = ks, ke
472 #ifdef DEBUG
473  call check( __line__, dens0(k,i,j) )
474  call check( __line__, mflx_hi(k ,i ,j ,xdir) )
475  call check( __line__, mflx_hi(k ,i-1,j ,xdir) )
476  call check( __line__, mflx_hi(k ,i ,j ,ydir) )
477  call check( __line__, mflx_hi(k ,i ,j-1,ydir) )
478  call check( __line__, dens_t(k,i,j) )
479 #endif
480  advcv = - ( mflx_hi(k,i,j,zdir)-mflx_hi(k-1,i ,j, zdir) ) * rcdz(k)
481  advch = - ( ( mflx_hi(k,i,j,xdir)-mflx_hi(k ,i-1,j, xdir) ) * rcdx(i) &
482  + ( mflx_hi(k,i,j,ydir)-mflx_hi(k ,i, j-1,ydir) ) * rcdy(j) ) &
483  * mapf(i,j,1,i_xy) * mapf(i,j,2,i_xy)
484  sr(k,i,j) = ( advcv + advch ) / gsqrt(k,i,j,i_xyz) + dens_t(k,i,j)
485 #ifdef HIST_TEND
486  if ( lhist ) then
487  advcv_t(k,i,j,i_dens) = advcv / gsqrt(k,i,j,i_xyz)
488  advch_t(k,i,j,i_dens) = advch / gsqrt(k,i,j,i_xyz)
489  end if
490 #endif
491  enddo
492  enddo
493  enddo
494 #ifdef DEBUG
495  k = iundef; i = iundef; j = iundef
496 #endif
497  profile_stop("hevi_sr")
498 
499  !##### momentum equation (z) #####
500 
501  ! at (x, y, z)
502  ! not that z-index is added by -1
503  profile_start("hevi_momz_qflxhi_z")
504  call atmos_dyn_fvm_fluxz_xyw( qflx_hi(:,:,:,zdir), & ! (out)
505  momz, momz, dens, & ! (in)
506  gsqrt(:,:,:,i_xyz), j33g, & ! (in)
507  num_diff(:,:,:,i_momz,zdir), & ! (in)
508  cdz, fdz, dtrk, &
509  iis, iie, jjs, jje ) ! (in)
510  profile_stop("hevi_momz_qflxhi_z")
511 
512  profile_start("hevi_momz_qflxj")
513  call atmos_dyn_fvm_fluxj13_xyw( qflx_j13, & ! (out)
514  momx, momz, dens, & ! (in)
515  gsqrt(:,:,:,i_xyz), j13g(:,:,:,i_xyz), mapf(:,:,:,i_xy), & ! (in)
516  cdz, &
517  iis, iie, jjs, jje ) ! (in)
518  call atmos_dyn_fvm_fluxj23_xyw( qflx_j23, & ! (out)
519  momy, momz, dens, & ! (in)
520  gsqrt(:,:,:,i_xyz), j23g(:,:,:,i_xyz), mapf(:,:,:,i_xy), & ! (in)
521  cdz, &
522  iis, iie, jjs, jje ) ! (in)
523  profile_stop("hevi_momz_qflxj")
524 
525  ! at (u, y, w)
526  profile_start("hevi_momz_qflxhi_x")
527  call atmos_dyn_fvm_fluxx_xyw( qflx_hi(:,:,:,xdir), & ! (out)
528  momx, momz, dens, & ! (in)
529  gsqrt(:,:,:,i_uyw), mapf(:,:,:,i_uy), & ! (in)
530  num_diff(:,:,:,i_momz,xdir), & ! (in)
531  cdz, & ! (in)
532  iis, iie, jjs, jje ) ! (in)
533  profile_stop("hevi_momz_qflxhi_x")
534 
535  ! at (x, v, w)
536  profile_start("hevi_momz_qflxhi_y")
537  call atmos_dyn_fvm_fluxy_xyw( qflx_hi(:,:,:,ydir), & ! (out)
538  momy, momz, dens, & ! (in)
539  gsqrt(:,:,:,i_xvw), mapf(:,:,:,i_xv), & ! (in)
540  num_diff(:,:,:,i_momz,ydir), & ! (in)
541  cdz, & ! (in)
542  iis, iie, jjs, jje ) ! (in)
543  profile_stop("hevi_momz_qflxhi_y")
544 
545  !--- update momentum(z)
546  profile_start("hevi_sw")
547  !$omp parallel do default(none) OMP_SCHEDULE_ collapse(2) &
548  !$omp private(i,j,k,advcv,advch,cf,wdmp,div) &
549 #ifdef HIST_TEND
550  !$omp shared(lhist,advcv_t,advch_t,wdmp_t,ddiv_t) &
551 #endif
552  !$omp shared(JJS,JJE,IIS,IIE,KS,KE) &
553  !$omp shared(qflx_hi,qflx_J13,qflx_J23,DDIV,MOMZ0,MOMZ_t,Sw) &
554  !$omp shared(RFDZ,RCDX,RCDY,FDZ,dtrk,wdamp_coef,divdmp_coef) &
555  !$omp shared(MAPF,GSQRT,I_XY,I_XYW)
556  do j = jjs, jje
557  do i = iis, iie
558  do k = ks, ke-1
559 #ifdef DEBUG
560  call check( __line__, qflx_hi(k ,i ,j ,zdir) )
561  call check( __line__, qflx_hi(k-1,i ,j ,zdir) )
562  call check( __line__, qflx_j13(k ,i ,j) )
563  call check( __line__, qflx_j13(k-1,i ,j) )
564  call check( __line__, qflx_j23(k ,i ,j) )
565  call check( __line__, qflx_j23(k-1,i ,j) )
566  call check( __line__, qflx_hi(k ,i ,j ,xdir) )
567  call check( __line__, qflx_hi(k ,i-1,j ,xdir) )
568  call check( __line__, qflx_hi(k ,i ,j ,ydir) )
569  call check( __line__, qflx_hi(k ,i ,j-1,ydir) )
570  call check( __line__, ddiv(k ,i,j) )
571  call check( __line__, ddiv(k+1,i,j) )
572  call check( __line__, momz0(k,i,j) )
573  call check( __line__, momz_t(k,i,j) )
574 #endif
575  advcv = - ( qflx_hi(k,i,j,zdir) - qflx_hi(k-1,i ,j ,zdir) &
576  + qflx_j13(k,i,j) - qflx_j13(k-1,i ,j ) &
577  + qflx_j23(k,i,j) - qflx_j23(k-1,i ,j ) ) * rfdz(k)
578  advch = - ( ( qflx_hi(k,i,j,xdir) - qflx_hi(k,i-1,j ,xdir) ) * rcdx(i) &
579  + ( qflx_hi(k,i,j,ydir) - qflx_hi(k,i ,j-1,ydir) ) * rcdy(j) ) &
580  * mapf(i,j,1,i_xy) * mapf(i,j,2,i_xy)
581  wdmp = - wdamp_coef(k) * momz0(k,i,j)
582  div = divdmp_coef / dtrk * ( ddiv(k+1,i,j)-ddiv(k,i,j) ) * fdz(k) ! divergence damping
583  sw(k,i,j) = ( advcv + advch ) / gsqrt(k,i,j,i_xyw) &
584  + wdmp + div + momz_t(k,i,j)
585 #ifdef HIST_TEND
586  if ( lhist ) then
587  advcv_t(k,i,j,i_momz) = advcv / gsqrt(k,i,j,i_xyw)
588  advch_t(k,i,j,i_momz) = advch / gsqrt(k,i,j,i_xyw)
589  wdmp_t(k,i,j) = wdmp
590  ddiv_t(k,i,j,1) = div
591  endif
592 #endif
593  enddo
594  enddo
595  enddo
596  profile_stop("hevi_sw")
597 #ifdef DEBUG
598  k = iundef; i = iundef; j = iundef
599 #endif
600 
601 
602  !##### Thermodynamic Equation #####
603 
604  !$omp parallel do default(none) private(i,j,k) OMP_SCHEDULE_ collapse(2) &
605  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,JHALO,IHALO,RHOT,DENS,POTT)
606  do j = jjs-jhalo, jje+jhalo
607  do i = iis-ihalo, iie+ihalo
608  do k = ks, ke
609 #ifdef DEBUG
610  call check( __line__, rhot(k,i,j) )
611  call check( __line__, dens(k,i,j) )
612 #endif
613  pott(k,i,j) = rhot(k,i,j) / dens(k,i,j)
614  enddo
615  enddo
616  enddo
617 #ifdef DEBUG
618  k = iundef; i = iundef; j = iundef
619 #endif
620 
621  ! at (x, y, w)
622  call atmos_dyn_fvm_fluxz_xyz( tflx_hi(:,:,:,zdir), & ! (out)
623  mflx_hi(:,:,:,zdir), pott, gsqrt(:,:,:,i_xyw), & ! (in)
624  num_diff(:,:,:,i_rhot,zdir), & ! (in)
625  cdz, & ! (in)
626  iis, iie, jjs, jje ) ! (in)
627 
628  ! at (u, y, z)
629  call atmos_dyn_fvm_fluxx_xyz( tflx_hi(:,:,:,xdir), & ! (out)
630  mflx_hi(:,:,:,xdir), pott, gsqrt(:,:,:,i_uyz), & ! (in)
631  num_diff(:,:,:,i_rhot,xdir), & ! (in)
632  cdz, & ! (in)
633  iis, iie, jjs, jje ) ! (in)
634 
635  ! at (x, v, z)
636  call atmos_dyn_fvm_fluxy_xyz( tflx_hi(:,:,:,ydir), & ! (out)
637  mflx_hi(:,:,:,ydir), pott, gsqrt(:,:,:,i_xvz), & ! (in)
638  num_diff(:,:,:,i_rhot,ydir), & ! (in)
639  cdz, & ! (in)
640  iis, iie, jjs, jje ) ! (in)
641 
642 
643  profile_start("hevi_st")
644  !$omp parallel do default(none) OMP_SCHEDULE_ collapse(2) &
645  !$omp private(i,j,k,advcv,advch) &
646 #ifdef HIST_TEND
647  !$omp shared(lhist,advcv_t,advch_t) &
648 #endif
649  !$omp shared(JJS,JJE,IIS,IIE,KS,KE) &
650  !$omp shared(tflx_hi,RHOT_t,St,RCDZ,RCDX,RCDY) &
651  !$omp shared(MAPF,GSQRT,I_XY,I_XYZ)
652  do j = jjs, jje
653  do i = iis, iie
654  do k = ks, ke
655 #ifdef DEBUG
656  call check( __line__, tflx_hi(k ,i ,j ,zdir) )
657  call check( __line__, tflx_hi(k-1,i ,j ,zdir) )
658  call check( __line__, tflx_hi(k ,i ,j ,xdir) )
659  call check( __line__, tflx_hi(k ,i-1,j ,xdir) )
660  call check( __line__, tflx_hi(k ,i ,j ,ydir) )
661  call check( __line__, tflx_hi(k ,i ,j-1,ydir) )
662  call check( __line__, rhot_t(k,i,j) )
663 #endif
664  advcv = - ( tflx_hi(k,i,j,zdir) - tflx_hi(k-1,i ,j ,zdir) ) * rcdz(k)
665  advch = - ( ( tflx_hi(k,i,j,xdir) - tflx_hi(k ,i-1,j ,xdir) ) * rcdx(i) &
666  + ( tflx_hi(k,i,j,ydir) - tflx_hi(k ,i ,j-1,ydir) ) * rcdy(j) ) &
667  * mapf(i,j,1,i_xy) * mapf(i,j,2,i_xy)
668  st(k,i,j) = ( advcv + advch ) / gsqrt(k,i,j,i_xyz) + rhot_t(k,i,j)
669 #ifdef HIST_TEND
670  if ( lhist ) then
671  advcv_t(k,i,j,i_rhot) = advcv / gsqrt(k,i,j,i_xyz)
672  advch_t(k,i,j,i_rhot) = advch / gsqrt(k,i,j,i_xyz)
673  endif
674 #endif
675  enddo
676  enddo
677  enddo
678 #ifdef DEBUG
679  k = iundef; i = iundef; j = iundef
680 #endif
681  profile_stop("hevi_st")
682 
683  ! implicit solver
684 
685  profile_start("hevi_solver")
686 
687 !OCL INDEPENDENT
688 !OCL PREFETCH_SEQUENTIAL(SOFT)
689 #ifndef __GFORTRAN__
690  !$omp parallel do default(none) OMP_SCHEDULE_ collapse(2) &
691  !$omp private(i,j,B,pg,advcv) &
692 #ifdef HIST_TEND
693  !$omp shared(lhist,pg_t,advcv_t) &
694 #endif
695 #ifdef DEBUG
696  !$omp shared(RHOT) &
697 #endif
698  !$omp shared(JJS,JJE,IIS,IIE,KS,KE) &
699  !$omp shared(mflx_hi,MOMZ_RK,MOMZ0) &
700  !$omp shared(DENS_RK,RHOT_RK,DENS0,RHOT0,DENS,MOMZ,PT,POTT,DPRES) &
701  !$omp shared(GRAV,dtrk,A,C,F1,F2,F3,REF_dens,Sr,Sw,St,RT2P) &
702  !$omp shared(ATMOS_DYN_FVM_flux_valueW_Z) &
703  !$omp shared(MAPF,GSQRT,J33G,I_XY,I_XYZ,I_XYW,CDZ,RCDZ,RFDZ)
704 #else
705  !$omp parallel do default(shared) private(i,j,k) OMP_SCHEDULE_ collapse(2) &
706  !$omp private(B,pg,advcv)
707 #endif
708  do j = jjs, jje
709  do i = iis, iie
710 
711  call atmos_dyn_fvm_flux_valuew_z( pt(:,i,j), & ! (out)
712  momz(:,i,j), pott(:,i,j), gsqrt(:,i,j,i_xyz), & ! (in)
713  cdz )
714 
715  do k = ks, ke
716  a(k,i,j) = dtrk**2 * j33g * rcdz(k) * rt2p(k,i,j) * j33g / gsqrt(k,i,j,i_xyz)
717  enddo
718  b = grav * dtrk**2 * j33g / ( cdz(ks+1) + cdz(ks) )
719  f1(ks,i,j) = - ( pt(ks+1,i,j) * rfdz(ks) * a(ks+1,i,j) + b ) / gsqrt(ks,i,j,i_xyw)
720  f2(ks,i,j) = 1.0_rp + ( pt(ks ,i,j) * rfdz(ks) * ( a(ks+1,i,j)+a(ks,i,j) ) ) / gsqrt(ks,i,j,i_xyw)
721  do k = ks+1, ke-2
722  b = grav * dtrk**2 * j33g / ( cdz(k+1) + cdz(k) )
723  f1(k,i,j) = - ( pt(k+1,i,j) * rfdz(k) * a(k+1,i,j) + b ) / gsqrt(k,i,j,i_xyw)
724  f2(k,i,j) = 1.0_rp + ( pt(k ,i,j) * rfdz(k) * ( a(k+1,i,j)+a(k,i,j) ) ) / gsqrt(k,i,j,i_xyw)
725  f3(k,i,j) = - ( pt(k-1,i,j) * rfdz(k) * a(k,i,j) - b ) / gsqrt(k,i,j,i_xyw)
726  enddo
727  b = grav * dtrk**2 * j33g / ( cdz(ke) + cdz(ke-1) )
728  f2(ke-1,i,j) = 1.0_rp + ( pt(ke-1,i,j) * rfdz(ke-1) * ( a(ke,i,j)+a(ke-1,i,j) ) ) / gsqrt(ke-1,i,j,i_xyw)
729  f3(ke-1,i,j) = - ( pt(ke-2,i,j) * rfdz(ke-1) * a(ke-1,i,j) - b ) / gsqrt(ke-1,i,j,i_xyw)
730  do k = ks, ke-1
731  pg = - ( dpres(k+1,i,j) + rt2p(k+1,i,j)*dtrk*st(k+1,i,j) &
732  - dpres(k ,i,j) - rt2p(k ,i,j)*dtrk*st(k ,i,j) ) &
733  * rfdz(k) * j33g / gsqrt(k,i,j,i_xyw) &
734  - grav &
735  * ( f2h(k,1,i_xyz) * ( dens(k+1,i,j) - ref_dens(k+1,i,j) + sr(k+1,i,j) * dtrk ) &
736  + f2h(k,2,i_xyz) * ( dens(k ,i,j) - ref_dens(k ,i,j) + sr(k ,i,j) * dtrk ) )
737  c(k-ks+1,i,j) = momz(k,i,j) + dtrk * ( pg + sw(k,i,j) )
738 #ifdef HIST_TEND
739  if ( lhist ) pg_t(k,i,j,1) = pg
740 #endif
741  enddo
742 
743  call solve_direct( &
744  c(:,i,j), & ! (inout)
745  f1(:,i,j), f2(:,i,j), f3(:,i,j) ) ! (in)
746 
747  do k = ks, ke-1
748 #ifdef DEBUG_HEVI2HEVE
749  ! for debug (change to explicit integration)
750  c(k-ks+1,i,j) = momz(k,i,j)
751  mflx_hi(k,i,j,zdir) = mflx_hi(k,i,j,zdir) &
752  + j33g * momz(k,i,j) / ( mapf(i,j,1,i_xy) * mapf(i,j,2,i_xy) )
753  momz_rk(k,i,j) = momz0(k,i,j) &
754  + dtrk*( &
755  - j33g * ( dpres(k+1,i,j)-dpres(k,i,j) ) * rfdz(k) / gsqrt(k,i,j,i_xyw) &
756  - grav * ( f2h(k,2,i_xyz)*(dens(k,i,j)-ref_dens(k,i,j))+f2h(k,1,i_xyz)*(dens(k+1,i,j)-ref_dens(k+1,i,j)) ) &
757  + sw(k,i,j) )
758 #else
759  ! z-momentum flux
760  mflx_hi(k,i,j,zdir) = mflx_hi(k,i,j,zdir) &
761  + j33g * c(k-ks+1,i,j) / ( mapf(i,j,1,i_xy) * mapf(i,j,2,i_xy) )
762  ! z-momentum
763  momz_rk(k,i,j) = momz0(k,i,j) &
764  + ( c(k-ks+1,i,j) - momz(k,i,j) )
765 #endif
766  enddo
767  momz_rk(ks-1,i,j) = 0.0_rp
768  momz_rk(ke ,i,j) = 0.0_rp
769 
770  ! density and rho*theta
771  advcv = - c(1,i,j) * j33g * rcdz(ks) / gsqrt(ks,i,j,i_xyz) ! C(0) = 0
772  dens_rk(ks,i,j) = dens0(ks,i,j) + dtrk * ( advcv + sr(ks,i,j) )
773 #ifdef HIST_TEND
774  if ( lhist ) advcv_t(ks,i,j,i_dens) = advcv
775 #endif
776  advcv = - c(1,i,j)*pt(ks,i,j) * j33g * rcdz(ks) / gsqrt(ks,i,j,i_xyz) ! C(0) = 0
777  rhot_rk(ks,i,j) = rhot0(ks,i,j) + dtrk * ( advcv + st(ks,i,j) )
778 #ifdef HIST_TEND
779  if ( lhist ) advcv_t(ks,i,j,i_rhot) = advcv
780 #endif
781  do k = ks+1, ke-1
782  advcv = - ( c(k-ks+1,i,j) - c(k-ks,i,j) ) &
783  * j33g * rcdz(k) / gsqrt(k,i,j,i_xyz)
784  dens_rk(k,i,j) = dens0(k,i,j) + dtrk * ( advcv + sr(k,i,j) )
785 #ifdef HIST_TEND
786  if ( lhist ) advcv_t(k,i,j,i_dens) = advcv
787 #endif
788  advcv = - ( c(k-ks+1,i,j)*pt(k,i,j) - c(k-ks,i,j)*pt(k-1,i,j) ) &
789  * j33g * rcdz(k) / gsqrt(k,i,j,i_xyz)
790  rhot_rk(k,i,j) = rhot0(k,i,j) + dtrk * ( advcv + st(k,i,j) )
791 #ifdef HIST_TEND
792  if ( lhist ) advcv_t(k,i,j,i_rhot) = advcv
793 #endif
794  enddo
795  advcv = c(ke-ks,i,j) * j33g * rcdz(ke) / gsqrt(ke,i,j,i_xyz) ! C(KE-KS+1) = 0
796  dens_rk(ke,i,j) = dens0(ke,i,j) + dtrk * ( advcv + sr(ke,i,j) )
797 #ifdef HIST_TEND
798  if ( lhist ) advcv_t(ke,i,j,i_dens) = advcv
799 #endif
800  advcv = c(ke-ks,i,j) * pt(ke-1,i,j) * j33g * rcdz(ke) / gsqrt(ke,i,j,i_xyz) ! C(KE-KS+1) = 0
801  rhot_rk(ke,i,j) = rhot0(ke,i,j) + dtrk * ( advcv + st(ke,i,j) )
802 #ifdef HIST_TEND
803  if ( lhist ) advcv_t(ke,i,j,i_rhot) = advcv
804 #endif
805 
806 #ifdef DEBUG
807  call check_equation( &
808  c(:,i,j), &
809  dens(:,i,j), momz(:,i,j), rhot(:,i,j), dpres(:,i,j), &
810  ref_dens(:,i,j), &
811  sr(:,i,j), sw(:,i,j), st(:,i,j), &
812  j33g, gsqrt(:,i,j,:), &
813  rt2p(:,i,j), &
814  dtrk, i, j )
815 #endif
816 
817  enddo
818  enddo
819 #ifdef DEBUG
820  k = iundef; i = iundef; j = iundef
821 #endif
822 
823  profile_stop("hevi_solver")
824 
825 
826  !##### momentum equation (x) #####
827 
828  profile_start("hevi_momx")
829  ! at (u, y, w)
830  call atmos_dyn_fvm_fluxz_uyz( qflx_hi(:,:,:,zdir), & ! (out)
831  momz, momx, dens, & ! (in)
832  gsqrt(:,:,:,i_uyw), j33g, & ! (in)
833  num_diff(:,:,:,i_momx,zdir), & ! (in)
834  cdz, & ! (in)
835  iis, iie, jjs, jje ) ! (in)
836  call atmos_dyn_fvm_fluxj13_uyz( qflx_j13, & ! (out)
837  momx, momx, dens, & ! (in)
838  gsqrt(:,:,:,i_uyz), j13g(:,:,:,i_uyw), mapf(:,:,:,i_uy), & ! (in)
839  cdz, & ! (in)
840  iis, iie, jjs, jje ) ! (in)
841  call atmos_dyn_fvm_fluxj23_uyz( qflx_j23, & ! (out)
842  momy, momx, dens, & ! (in)
843  gsqrt(:,:,:,i_uyz), j23g(:,:,:,i_uyw), mapf(:,:,:,i_uy), & ! (in)
844  cdz, & ! (in)
845  iis, iie, jjs, jje ) ! (in)
846 
847  ! at (x, y, z)
848  ! note that x-index is added by -1
849  call atmos_dyn_fvm_fluxx_uyz( qflx_hi(:,:,:,xdir), & ! (out)
850  momx, momx, dens, & ! (in)
851  gsqrt(:,:,:,i_xyz), mapf(:,:,:,i_xy), & ! (in)
852  num_diff(:,:,:,i_momx,xdir), & ! (in)
853  cdz, & ! (in)
854  iis, iie, jjs, jje ) ! (in)
855 
856  ! at (u, v, z)
857  call atmos_dyn_fvm_fluxy_uyz( qflx_hi(:,:,:,ydir), & ! (out)
858  momy, momx, dens, & ! (in)
859  gsqrt(:,:,:,i_uvz), mapf(:,:,1,i_uv), & ! (in)
860  num_diff(:,:,:,i_momx,ydir), & ! (in)
861  cdz, & ! (in)
862  iis, iie, jjs, jje ) ! (in)
863 
864  !--- update momentum(x)
865  iee = min(iie,ieh)
866  !$omp parallel do default(none) OMP_SCHEDULE_ collapse(2) &
867  !$omp private(i,j,k,advch,advcv,pg,cf,div) &
868 #ifdef HIST_TEND
869  !$omp shared(lhist,advch_t,advcv_t,pg_t,cf_t,ddiv_t) &
870 #endif
871  !$omp shared(JJS,JJE,IIS,iee,KS,KE) &
872  !$omp shared(MOMX_RK,DPRES,DENS,MOMX,MOMY,DDIV,MOMX0,MOMX_t) &
873  !$omp shared(qflx_hi,qflx_J13,qflx_J23) &
874  !$omp shared(RCDZ,RCDY,RFDX,CDZ,FDX) &
875  !$omp shared(MAPF,GSQRT,J13G,I_XY,I_UY,I_UV,I_XYZ,I_UYW,I_UYZ) &
876  !$omp shared(dtrk,CORIOLI,divdmp_coef)
877  do j = jjs, jje
878  do i = iis, iee
879  do k = ks, ke
880 #ifdef DEBUG
881  call check( __line__, qflx_hi(k ,i ,j ,zdir) )
882  call check( __line__, qflx_hi(k-1,i ,j ,zdir) )
883  call check( __line__, qflx_hi(k ,i ,j ,xdir) )
884  call check( __line__, qflx_hi(k ,i-1,j ,xdir) )
885  call check( __line__, qflx_hi(k ,i ,j ,ydir) )
886  call check( __line__, qflx_hi(k ,i ,j-1,ydir) )
887  call check( __line__, dpres(k,i+1,j) )
888  call check( __line__, dpres(k,i ,j) )
889  call check( __line__, corioli(1,i ,j) )
890  call check( __line__, corioli(1,i+1,j) )
891  call check( __line__, momy(k,i ,j ) )
892  call check( __line__, momy(k,i+1,j ) )
893  call check( __line__, momy(k,i ,j-1) )
894  call check( __line__, momy(k,i+1,j-1) )
895  call check( __line__, ddiv(k,i+1,j) )
896  call check( __line__, ddiv(k,i ,j) )
897  call check( __line__, momx0(k,i,j) )
898 #endif
899  advcv = - ( qflx_hi(k,i,j,zdir) - qflx_hi(k-1,i ,j ,zdir) &
900  + qflx_j13(k,i,j) - qflx_j13(k-1,i ,j ) &
901  + qflx_j23(k,i,j) - qflx_j23(k-1,i ,j ) ) * rcdz(k)
902  advch = - ( ( qflx_hi(k,i,j,xdir) - qflx_hi(k ,i-1,j ,xdir) ) * rfdx(i) &
903  + ( qflx_hi(k,i,j,ydir) - qflx_hi(k ,i ,j-1,ydir) ) * rcdy(j) ) &
904  * mapf(i,j,1,i_uy) * mapf(i,j,2,i_uy)
905  pg = ( ( gsqrt(k,i+1,j,i_xyz) * dpres(k,i+1,j) & ! [x,y,z]
906  - gsqrt(k,i ,j,i_xyz) * dpres(k,i ,j) & ! [x,y,z]
907  ) * rfdx(i) &
908  + ( j13g(k ,i,j,i_uyw) &
909  * 0.5_rp * ( f2h(k,1,i_uyz) * ( dpres(k+1,i+1,j)+dpres(k+1,i,j) ) &
910  + f2h(k,2,i_uyz) * ( dpres(k ,i+1,j)+dpres(k ,i,j) ) ) & ! [x,y,z->u,y,w]
911  - j13g(k-1,i,j,i_uyw) &
912  * 0.5_rp * ( f2h(k,1,i_uyz) * ( dpres(k ,i+1,j)+dpres(k ,i,j) ) &
913  + f2h(k,2,i_uyz) * ( dpres(k-1,i+1,j)+dpres(k-1,i,j) ) ) & ! [x,y,z->u,y,w]
914  ) * rcdz(k) ) &
915  * mapf(i,j,1,i_uy)
916  cf = 0.125_rp * ( corioli(1,i+1,j )+corioli(1,i,j ) ) & ! [x,y,z->u,y,z]
917  * ( momy(k,i+1,j )+momy(k,i,j ) &
918  + momy(k,i+1,j-1)+momy(k,i,j-1) ) & ! [x,v,z->u,y,z]
919  + 0.25_rp * mapf(i,j,1,i_uy) * mapf(i,j,2,i_uy) &
920  * ( momy(k,i,j) + momy(k,i,j-1) + momy(k,i+1,j) + momy(k,i+1,j-1) ) &
921  * ( ( momy(k,i,j) + momy(k,i,j-1) + momy(k,i+1,j) + momy(k,i+1,j-1) ) * 0.25_rp &
922  * ( 1.0_rp/mapf(i+1,j,2,i_xy) - 1.0_rp/mapf(i,j,2,i_xy) ) * rfdx(i) &
923  - momx(k,i,j) &
924  * ( 1.0_rp/mapf(i,j,1,i_uv) - 1.0_rp/mapf(i,j-1,1,i_uv) ) * rcdy(j) ) &
925  * 2.0_rp / ( dens(k,i+1,j) + dens(k,i,j) ) ! metric term
926  div = divdmp_coef / dtrk * ( ddiv(k,i+1,j)/mapf(i+1,j,2,i_xy) - ddiv(k,i,j)/mapf(i,j,1,i_xy) ) &
927  * mapf(i,j,1,i_uy) * mapf(i,j,2,i_uy) * fdx(i) ! divergence damping
928  momx_rk(k,i,j) = momx0(k,i,j) &
929  + dtrk * ( ( advcv + advch - pg ) / gsqrt(k,i,j,i_uyz) + cf + div + momx_t(k,i,j) )
930 #ifdef HIST_TEND
931  if ( lhist ) then
932  advcv_t(k,i,j,i_momx) = advcv / gsqrt(k,i,j,i_uyz)
933  advch_t(k,i,j,i_momx) = advch / gsqrt(k,i,j,i_uyz)
934  pg_t(k,i,j,2) = - pg / gsqrt(k,i,j,i_uyz)
935  cf_t(k,i,j,1) = cf
936  ddiv_t(k,i,j,2) = div
937  endif
938 #endif
939  enddo
940  enddo
941  enddo
942  profile_stop("hevi_momx")
943 #ifdef DEBUG
944  k = iundef; i = iundef; j = iundef
945 #endif
946 
947  !##### momentum equation (y) #####
948  profile_start("hevi_momy")
949  ! at (x, v, w)
950  call atmos_dyn_fvm_fluxz_xvz( qflx_hi(:,:,:,zdir), & ! (out)
951  momz, momy, dens, & ! (in)
952  gsqrt(:,:,:,i_xvw), j33g, & ! (in)
953  num_diff(:,:,:,i_momy,zdir), & ! (in)
954  cdz, & ! (in)
955  iis, iie, jjs, jje ) ! (in)
956  call atmos_dyn_fvm_fluxj13_xvz( qflx_j13, & ! (out)
957  momx, momy, dens, & ! (in)
958  gsqrt(:,:,:,i_xvz), j13g(:,:,:,i_xvw), mapf(:,:,:,i_xv), & ! (in)
959  cdz, & ! (in)
960  iis, iie, jjs, jje ) ! (in)
961  call atmos_dyn_fvm_fluxj23_xvz( qflx_j23, & ! (out)
962  momy, momy, dens, & ! (in)
963  gsqrt(:,:,:,i_xvz), j23g(:,:,:,i_xvw), mapf(:,:,:,i_xv), & ! (in)
964  cdz, & ! (in)
965  iis, iie, jjs, jje ) ! (in)
966 
967  ! at (u, v, z)
968  call atmos_dyn_fvm_fluxx_xvz( qflx_hi(:,:,:,xdir), & ! (out)
969  momx, momy, dens, & ! (in)
970  gsqrt(:,:,:,i_uvz), mapf(:,:,:,i_uv), & ! (in)
971  num_diff(:,:,:,i_momy,xdir), & ! (in)
972  cdz, & ! (in)
973  iis, iie, jjs, jje ) ! (in)
974 
975  ! at (x, y, z)
976  ! note that y-index is added by -1
977  call atmos_dyn_fvm_fluxy_xvz( qflx_hi(:,:,:,ydir), & ! (out)
978  momy, momy, dens, & ! (in)
979  gsqrt(:,:,:,i_xyz), mapf(:,:,:,i_xy), & ! (in)
980  num_diff(:,:,:,i_momy,ydir), & ! (in
981  cdz, & ! (in)
982  iis, iie, jjs, jje ) ! (in)
983 
984  !--- update momentum(y)
985  !$omp parallel do default(none) OMP_SCHEDULE_ collapse(2) &
986  !$omp private(i,j,k,advch,advcv,pg,cf,div) &
987 #ifdef HIST_TEND
988  !$omp shared(lhist,advch_t,advcv_t,pg_t,cf_t,ddiv_t) &
989 #endif
990  !$omp shared(JJS,JJE,JEH,IIS,IIE,KS,KE) &
991  !$omp shared(MOMY_RK,DPRES,DENS,MOMX,MOMY,DDIV,MOMY0,MOMY_t) &
992  !$omp shared(qflx_hi,qflx_J13,qflx_J23) &
993  !$omp shared(RCDZ,RCDX,RFDY,CDZ,FDY) &
994  !$omp shared(MAPF,GSQRT,J23G,I_XY,I_XV,I_UV,I_XYZ,I_XVW,I_XVZ) &
995  !$omp shared(dtrk,CORIOLI,divdmp_coef)
996  do j = jjs, min(jje,jeh)
997  do i = iis, iie
998  do k = ks, ke
999 #ifdef DEBUG
1000  call check( __line__, qflx_hi(k ,i ,j ,zdir) )
1001  call check( __line__, qflx_hi(k-1,i ,j ,zdir) )
1002  call check( __line__, qflx_hi(k ,i ,j ,xdir) )
1003  call check( __line__, qflx_hi(k ,i-1,j ,xdir) )
1004  call check( __line__, qflx_hi(k ,i ,j ,ydir) )
1005  call check( __line__, qflx_hi(k ,i ,j-1,ydir) )
1006  call check( __line__, dpres(k,i,j ) )
1007  call check( __line__, dpres(k,i,j+1) )
1008  call check( __line__, corioli(1,i,j ) )
1009  call check( __line__, corioli(1,i,j+1) )
1010  call check( __line__, momx(k,i ,j ) )
1011  call check( __line__, momx(k,i ,j+1) )
1012  call check( __line__, momx(k,i-1,j ) )
1013  call check( __line__, momx(k,i-1,j+1) )
1014  call check( __line__, ddiv(k,i,j+1) )
1015  call check( __line__, ddiv(k,i,j ) )
1016  call check( __line__, momy_t(k,i,j) )
1017  call check( __line__, momy0(k,i,j) )
1018 #endif
1019  advcv = - ( qflx_hi(k,i,j,zdir) - qflx_hi(k-1,i ,j ,zdir) &
1020  + qflx_j13(k,i,j) - qflx_j13(k-1,i ,j ) &
1021  + qflx_j23(k,i,j) - qflx_j23(k-1,i ,j ) ) * rcdz(k)
1022  advch = - ( ( qflx_hi(k,i,j,xdir) - qflx_hi(k ,i-1,j ,xdir) ) * rcdx(i) &
1023  + ( qflx_hi(k,i,j,ydir) - qflx_hi(k ,i ,j-1,ydir) ) * rfdy(j) ) &
1024  * mapf(i,j,1,i_xv) * mapf(i,j,2,i_xv)
1025  pg = ( ( gsqrt(k,i,j+1,i_xyz) * dpres(k,i,j+1) & ! [x,y,z]
1026  - gsqrt(k,i,j ,i_xyz) * dpres(k,i,j ) & ! [x,y,z]
1027  ) * rfdy(j) &
1028  + ( j23g(k ,i,j,i_xvw) &
1029  * 0.5_rp * ( f2h(k ,1,i_xvz) * ( dpres(k+1,i,j+1)+dpres(k+1,i,j) ) &
1030  + f2h(k ,2,i_xvz) * ( dpres(k ,i,j+1)+dpres(k ,i,j) ) ) & ! [x,y,z->x,v,w]
1031  - j23g(k-1,i,j,i_xvw) &
1032  * 0.5_rp * ( f2h(k-1,1,i_xvz) * ( dpres(k ,i,j+1)+dpres(k ,i,j) ) &
1033  + f2h(k-1,2,i_xvz) * ( dpres(k-1,i,j+1)+dpres(k-1,i,j) ) ) & ! [x,y,z->x,v,w]
1034  ) * rcdz(k) ) &
1035  * mapf(i,j,2,i_xv)
1036  cf = - 0.125_rp * ( corioli(1,i ,j+1)+corioli(1,i ,j) ) & ! [x,y,z->x,v,z]
1037  * ( momx(k,i ,j+1)+momx(k,i ,j) &
1038  + momx(k,i-1,j+1)+momx(k,i-1,j) ) & ! [u,y,z->x,v,z]
1039  - 0.25_rp * mapf(i,j,1,i_xv) * mapf(i,j,2,i_xv) &
1040  * ( momx(k,i,j) + momx(k,i-1,j) + momx(k,i,j+1) + momx(k,i-1,j+1) ) &
1041  * ( momy(k,i,j) &
1042  * ( 1.0_rp/mapf(i,j,2,i_uv) - 1.0_rp/mapf(i-1,j,2,i_uv) ) * rcdx(i) &
1043  - 0.25_rp * ( momx(k,i,j)+momx(k,i-1,j)+momx(k,i,j+1)+momx(k,i-1,j+1) ) &
1044  * ( 1.0_rp/mapf(i,j+1,1,i_xy) - 1.0_rp/mapf(i,j,1,i_xy) ) * rfdy(j) ) &
1045  * 2.0_rp / ( dens(k,i,j+1) + dens(k,i,j) ) ! metoric term
1046  div = divdmp_coef / dtrk * ( ddiv(k,i,j+1)/mapf(i,j+1,1,i_xy) - ddiv(k,i,j)/mapf(i,j,1,i_xy) ) &
1047  * mapf(i,j,1,i_xv) * mapf(i,j,2,i_xv) * fdy(j) ! divergence damping
1048  momy_rk(k,i,j) = momy0(k,i,j) &
1049  + dtrk * ( ( advcv + advch - pg ) / gsqrt(k,i,j,i_xvz) + cf + div + momy_t(k,i,j) )
1050 #ifdef HIST_TEND
1051  if ( lhist ) then
1052  advcv_t(k,i,j,i_momy) = advcv / gsqrt(k,i,j,i_xvz)
1053  advch_t(k,i,j,i_momy) = advch / gsqrt(k,i,j,i_xvz)
1054  pg_t(k,i,j,3) = - pg / gsqrt(k,i,j,i_xvz)
1055  cf_t(k,i,j,2) = cf
1056  ddiv_t(k,i,j,3) = div
1057  endif
1058 #endif
1059  enddo
1060  enddo
1061  enddo
1062  profile_stop("hevi_momy")
1063 #ifdef DEBUG
1064  k = iundef; i = iundef; j = iundef
1065 #endif
1066 
1067  enddo
1068  enddo
1069 
1070 #ifdef PROFILE_FIPP
1071  call fipp_stop()
1072 #endif
1073 
1074 #ifdef HIST_TEND
1075  if ( lhist ) then
1076  call hist_in(advcv_t(:,:,:,i_dens), 'DENS_t_advcv', 'tendency of density (vert. advection)', 'kg/m3/s' )
1077  call hist_in(advcv_t(:,:,:,i_momz), 'MOMZ_t_advcv', 'tendency of momentum z (vert. advection)', 'kg/m2/s2', zdim='half')
1078  call hist_in(advcv_t(:,:,:,i_momx), 'MOMX_t_advcv', 'tendency of momentum x (vert. advection)', 'kg/m2/s2', xdim='half')
1079  call hist_in(advcv_t(:,:,:,i_momy), 'MOMY_t_advcv', 'tendency of momentum y (vert. advection)', 'kg/m2/s2', ydim='half')
1080  call hist_in(advcv_t(:,:,:,i_rhot), 'RHOT_t_advcv', 'tendency of rho*theta (vert. advection)', 'K kg/m3/s' )
1081 
1082  call hist_in(advch_t(:,:,:,i_dens), 'DENS_t_advch', 'tendency of density (horiz. advection)', 'kg/m3/s' )
1083  call hist_in(advch_t(:,:,:,i_momz), 'MOMZ_t_advch', 'tendency of momentum z (horiz. advection)', 'kg/m2/s2', zdim='half')
1084  call hist_in(advch_t(:,:,:,i_momx), 'MOMX_t_advch', 'tendency of momentum x (horiz. advection)', 'kg/m2/s2', xdim='half')
1085  call hist_in(advch_t(:,:,:,i_momy), 'MOMY_t_advch', 'tendency of momentum y (horiz. advection)', 'kg/m2/s2', ydim='half')
1086  call hist_in(advch_t(:,:,:,i_rhot), 'RHOT_t_advch', 'tendency of rho*theta (horiz. advection)', 'K kg/m3/s' )
1087 
1088  call hist_in(pg_t(:,:,:,1), 'MOMZ_t_pg', 'tendency of momentum z (pressure gradient)', 'kg/m2/s2', zdim='half')
1089  call hist_in(pg_t(:,:,:,2), 'MOMX_t_pg', 'tendency of momentum x (pressure gradient)', 'kg/m2/s2', xdim='half')
1090  call hist_in(pg_t(:,:,:,3), 'MOMY_t_pg', 'tendency of momentum y (pressure gradient)', 'kg/m2/s2', ydim='half')
1091 
1092  call hist_in(wdmp_t(:,:,:), 'MOMZ_t_wdamp', 'tendency of momentum z (Rayleigh damping)', 'kg/m2/s2', zdim='half')
1093 
1094  call hist_in(ddiv_t(:,:,:,1), 'MOMZ_t_ddiv', 'tendency of momentum z (divergence damping)', 'kg/m2/s2', zdim='half')
1095  call hist_in(ddiv_t(:,:,:,2), 'MOMX_t_ddiv', 'tendency of momentum x (divergence damping)', 'kg/m2/s2', xdim='half')
1096  call hist_in(ddiv_t(:,:,:,3), 'MOMY_t_ddiv', 'tendency of momentum y (divergence damping)', 'kg/m2/s2', ydim='half')
1097 
1098  call hist_in(cf_t(:,:,:,1), 'MOMX_t_cf', 'tendency of momentum x (coliolis force)', 'kg/m2/s2', xdim='half')
1099  call hist_in(cf_t(:,:,:,2), 'MOMY_t_cf', 'tendency of momentum y (coliolis force)', 'kg/m2/s2', ydim='half')
1100  endif
1101 #endif
1102 
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_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_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_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_grid_index::ydir
integer, parameter, public ydir
Definition: scale_grid_index.F90:37
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_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_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_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_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_grid_index::jhalo
integer, public jhalo
of halo cells: y
Definition: scale_grid_index.F90:52
scale_grid_index::js
integer, public js
start point of inner domain: y, local
Definition: scale_grid_index.F90:62
scale_atmos_dyn_common
module Atmosphere / Dynamics common
Definition: scale_atmos_dyn_common.F90:15
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_grid_index::jeh
integer, public jeh
end point of inner domain: y, local (half level)
Definition: scale_grid_index.F90:76
scale_grid_index::ieh
integer, public ieh
end point of inner domain: x, local (half level)
Definition: scale_grid_index.F90:75
scale_gridtrans::i_xv
integer, public i_xv
Definition: scale_gridtrans.F90:59
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_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_history
module HISTORY
Definition: scale_history.F90:16
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_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_gridtrans::i_uv
integer, public i_uv
Definition: scale_gridtrans.F90:60
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_grid_index::ihalo
integer, public ihalo
of halo cells: x
Definition: scale_grid_index.F90:51
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
Here is the call graph for this function:
Here is the caller graph for this function:

◆ solve_direct()

subroutine scale_atmos_dyn_tstep_short_fvm_hevi::solve_direct ( real(rp), dimension(kmax-1), intent(inout)  C,
real(rp), dimension(ka), intent(in)  F1,
real(rp), dimension(ka), intent(in)  F2,
real(rp), dimension(ka), intent(in)  F3 
)

Definition at line 1109 of file scale_atmos_dyn_tstep_short_fvm_hevi.F90.

References scale_const::const_eps, scale_const::const_grav, scale_grid::grid_rcdz, scale_grid::grid_rfdz, scale_gridtrans::i_xyw, scale_gridtrans::i_xyz, scale_grid_index::ke, scale_grid_index::ks, and scale_process::prc_mpistop().

Referenced by atmos_dyn_tstep_short_fvm_hevi().

1109 
1110  use scale_process, only: &
1111  prc_mpistop
1112  implicit none
1113  real(RP), intent(inout) :: C(KMAX-1)
1114  real(RP), intent(in) :: F1(KA)
1115  real(RP), intent(in) :: F2(KA)
1116  real(RP), intent(in) :: F3(KA)
1117 
1118  real(RP) :: e(KMAX-2)
1119  real(RP) :: f(KMAX-2)
1120 
1121  real(RP) :: rdenom
1122 
1123  integer :: k
1124 
1125  rdenom = 1.0_rp / f2(ks)
1126  e(1) = - f1(ks) * rdenom
1127  f(1) = c(1) * rdenom
1128  do k = 2, kmax-2
1129  rdenom = 1.0_rp / ( f2(k+ks-1) + f3(k+ks-1) * e(k-1) )
1130  e(k) = - f1(k+ks-1) * rdenom
1131  f(k) = ( c(k) - f3(k+ks-1) * f(k-1) ) * rdenom
1132  enddo
1133 
1134  ! C = \rho w
1135  c(kmax-1) = ( c(kmax-1) - f3(ke-1) * f(kmax-2) ) &
1136  / ( f2(ke-1) + f3(ke-1) * e(kmax-2) ) ! C(KMAX-1) = f(KMAX-1)
1137  do k = kmax-2, 1, -1
1138  c(k) = e(k) * c(k+1) + f(k)
1139  enddo
1140 
1141  return
scale_process::prc_mpistop
subroutine, public prc_mpistop
Abort MPI.
Definition: scale_process.F90:267
scale_grid_index::ke
integer, public ke
end point of inner domain: z, local
Definition: scale_grid_index.F90:59
scale_grid_index::kmax
integer, public kmax
of computational cells: z, local
Definition: scale_grid_index.F90:43
scale_process
module PROCESS
Definition: scale_process.F90:14
scale_grid_index::ks
integer, public ks
start point of inner domain: z, local
Definition: scale_grid_index.F90:58
Here is the call graph for this function:
Here is the caller graph for this function:
Generated by   doxygen 1.8.14